Introduction to Small Area Estimation

Rajesh Majumder

07-12-2025

Definition

Small Area Estimation (SAE) refers to statistical methods used to produce reliable estimates for small geographic areas (e.g., districts, blocks, villages) or small population subgroups (e.g., by age, sex, occupation) where sample sizes from surveys are too small to yield direct reliable estimates.

       A small area does not always mean geographically small; it can also mean small subpopulations.

       The main purpose is to borrow strength from related areas or auxiliary information (like census data, administrative records, or satellite data) to improve estimation.

Why SAE is Needed?

       National surveys are typically designed to give accurate estimates at national or state/provincial level.

       For local planning and policy decisions, estimates at the district, sub-district, or demographic subgroup level are required.

       Direct estimates (based only on local sample data) are unreliable because of high variance (small n).

       SAE provides more precise estimates by combining survey data with auxiliary sources and statistical models.

Types of Small Area Estimation Techniques

SAE techniques are generally divided into two categories:

1.    Direct Estimation (Design-based)

       Based only on sample data from the target area.

       Example: Using weighted means or proportions for the district sample.

       Limitation: If sample size is small, estimates are unstable (large standard error).

       Direct methods are rarely useful for very small domains.

2.    Indirect (Model-based) Estimation

These methods borrow strength across areas by linking small areas using models and auxiliary data. There is various type of modelling strategies:

A.   Synthetic Estimation

Idea: Apply relationship observed at a larger area (e.g., state level) to smaller areas.

Example: If state-level survey shows 20% prevalence of diabetes, assume each district has the same prevalence.

Limitation: Ignores local variations → may be biased if small areas differ substantially.

B.   Regression Synthetic Estimation

Uses regression models linking survey outcome to auxiliary variables (from census/administrative data).

Predicted values for small areas are taken as estimates.

Limitation: Doesn t account for small area random effects, so still may be biased.

C.   Composite Estimation

Combines direct survey estimate and synthetic estimate with weights.

Example:   

where weight depends on sample size/variance.

Advantage: Balances between unbiased direct estimate and stable synthetic estimate.

D.  Model-Based SAE (Hierarchical / Mixed Models)

This is the most widely used category.

a.     Area-level Models (Fay-Herriot Model)

       Developed by Fay & Herriot (1979) for estimating per-capita income in small places in the U.S.

       It is an area-level model, meaning it works with direct survey estimates at the area level (not unit-level data).

       Widely used in official statistics because many national surveys release only area-level direct estimates and their variances.

Suppose we want to estimate a parameter (e.g., mean income, disease prevalence) for small areas

From a sample survey, we have:

       Direct estimate for area

      Its known sampling variance:

We also have auxiliary information (from census, admin records, etc.):

       Covariates: (vector of predictors for area )

Fay Herriot Model is a two-level model (hierarchical / mixed model):

Level 1: Sampling model (measurement error model)

       : direct survey estimate for area

       : true small area mean (unknown)

       : sampling error, variance (assumed known from survey design)

Level 2: Linking model (area-level regression)

      : regression component (fixed effect of covariates)

       : area-specific random effect, accounts for unexplained area heterogeneity

       : variance of random effects (unknown, to be estimated)

Combined Model

This is a linear mixed model with two sources of error:

       Sampling error ()

       Area effect ()

Estimation

                                  I.         Parameters ()

       Estimated by Maximum Likelihood (ML) or Restricted Maximum Likelihood (REML).

       Some applications use Method of Moments for .

                          II.        Small Area Predictions

The Empirical Best Linear Unbiased Predictor (EBLUP) of is:

 

Where

       is a shrinkage factor:

       If area sample size is large → small → → estimate relies mostly on direct survey data.

       If area sample size is small → large → → estimate relies mostly on model predictions.

Thus, EBLUP is a weighted compromise between direct and synthetic estimates.

Mean Squared Error (MSE)

       MSE of is estimated using analytical formulas or bootstrap.

       Confidence intervals can be constructed accordingly.

 

b.   Unit-level Models (Battese-Harter-Fuller Model)

Why Unit-Level Models?

Unit-level models are used when:

       We have micro-data (individual or household records).

       We want small area means/totals, but sample sizes in each area are small.

       Auxiliary variables are available for every unit in the population, or at least their area-level means are known.

Unlike Fay Herriot (which uses area-level aggregates), BHF directly models the unit-level responses.

The Model Structure

For each sampled unit j in area i, the model is:

Where:

Component	Meaning
	Observed unit-level response (income, yield, etc.)
	Unit-level auxiliary covariates
	Fixed-effect regression parameters
	Area-level random effect
	Unit-level residual error

 

Distributional Assumptions

       Area effects (between-area variation):

       Unit errors (within-area variation):

Both independent of each other.

What Does the Model Capture?

✔ Between-area variation

captures differences among areas (e.g., districts, villages).

✔ Within-area individual variation

captures variation among individuals inside a given area.

✔ Uses complete micro-data

This allows more precision than area-level models.

The Target: Predict the Area Mean

We want to estimate the small area mean:

But we only observe part of the area (small samples).

The EBLUP Under the BHF Model

The Best Linear Unbiased Predictor (BLUP) for the area mean is:

Where:

       is the population mean of auxiliary variables (known or estimated).

      is the fixed-effect estimate (from lmer).

       is the predicted random effect:

Here:

       = sample mean in area

       = sample mean of covariates in area

       = unit-level shrinkage weight

Shrinkage Weight in BHF Model

Where:

       = sample size in area

Interpretation:

Case	Meaning
Large	Shrinkage low → direct estimator dominates
Small	Shrinkage high → model-based prediction dominates
Large	High between-area variation → random effect more important
Large	High within-area noise → shrink more toward synthetic

Connection to Composite Estimator

The BHF predictor is equivalent to:

So:

✔ It is a composite of direct and model-based predictions

✔ Automatically weights them based on model variances

✔ No need to manually compute composite weights (unlike FH)

Advantages of Unit-Level (BHF) Model

       Uses full micro-data, not just area means

Much more efficient than FH.

       Handles individual covariates

Allows rich modeling (age, education, land size, etc.).

       Natural mixed-model structure

Estimation via REML is straightforward.

       Provides EBLUP automatically

No manual composite estimator required.

What You Need to Apply BHF

You need:

       Unit-level survey data (sample)

       Population totals or means of auxiliary variables

       Unit- and area-level covariates

       Mixed model software (like lmer())

 

c.     Hierarchical Bayesian (HB) Models

       Extends Fay-Herriot and unit-level models using Bayesian framework.

       Prior distributions are assumed for parameters, and posterior distributions are used for inference.

       Advantages:

o   Produces full probability distributions (credible intervals).

o   Flexible to include spatial correlations, non-linear effects.

E.   Empirical Bayes Methods

       A frequentist version of Bayesian approach.

       Area estimates are shrunk towards regression predictions depending on reliability of survey data.

       Example: If a district has large sample size, direct estimate is given more weight; if small sample size, synthetic estimate dominates.

 

F.    Spatial Small Area Estimation

       Incorporates spatial correlation among neighboring areas.

       Example: Disease prevalence in one district may resemble neighboring districts.

       Spatial Fay-Herriot models and Conditional Autoregressive (CAR) priors are often used.

 

G.  Advanced Methods (AI/ML)

       M-quantile Models: Robust alternative to mixed models, less sensitive to outliers.

       Machine Learning SAE: Random Forest, Gradient Boosting applied with auxiliary data.

       Model-assisted approaches: Combine survey design weights with regression models.

Applications of SAE

• Public Health: Estimating district-level prevalence of diseases, malnutrition, risk factors (WHO, CDC often use SAE).
• Poverty Mapping: Estimating poverty at district/block level using survey + census data (World Bank).
• Agriculture: Crop yield estimation at small administrative levels.
• Education: Literacy or school attendance at local levels.
• Epidemiology: Rare disease prevalence estimation at local level.

References for Further Reading

1. Rao, J.N.K., & Molina, I. (2015). Small Area Estimation (2nd ed.). Wiley.
2. Fay, R.E., & Herriot, R.A. (1979). Estimates of income for small places: An application of James-Stein procedures to census data. Journal of the American Statistical Association.
3. Pfeffermann, D. (2013). New important developments in small area estimation. Statistical Science, 28(1), 40 68.

Packages for R and Python

R Packages

       Sae

       sae2

       hbsae

       saeRobust

       rstanarm / brms / INLA

       emdi (excellent for poverty SAE)

       survey (for direct estimators)

       mboost / ranger (ML-assisted SAE)

       Josae (Johannes Breidenbach). Unit and area-level EBLUP estimators including also heteroscedasticity.

       hbsae (Harm Jan Boonstra): Hierarchical Bayes methods for basic area-level and unit-level models (also REML fitting).

       BayesSAE (Chengchun Shi) Bayesian methods for a variety of models, including unmatched models and spatial models.

       saeeb (Rizki Ananda Fauziah, Ika Yuni Wulansari): EB estimators under small area models for counts.

       rsae (Tobias Schoch). Robust methods for area and unit-level models.

       saeRobust (Sebastian Warnhol): Robust EBLUP under area-level models, including models with spatial and temporal correlation.

       saeSim (Sebastian Warnholz, Timo Schmid): Simulation and model fitting in SAE.

       saery (M.D. Esteban, D. Morales, A. Perez): EBLUP for temporal area-level Rao-Yu model.

       mme (E. Lopez-Vizcaino, M.J. Lombardia and D. Morales). Multinomial area-level models for small area estimation of proportions, including models with temporal correlation.

       saeME (Muhammad Rifqi Mubarak, Azka Ubaidillah): SAE under measurement error of covariates.

       emdi (Sylvia Harmening, Ann-Kristin Kreutzmann, Soeren Pannier, Natalia Rojas-Perilla, Nicola Salvati, Timo Schmid, Matthias Templ, Nikos Tzavidis, Nora Wurz): Functions that support estimating, assessing and mapping regional disaggregated indicators.

Python

       PyMC (Bayesian SAE)

       Statsmodels (limited but usable)

Simulation example in R

rm(list=ls())
library(dplyr)

library(survey)

library(sae)

library(mgcv)

library(randomForest)

###############################################
# COMPLETE WORKING SMALL AREA ESTIMATION DEMO
###############################################
set.seed(123)
#
##############################
### 1. SIMULATE POPULATION
##############################
m= 20 # Number of small areas
Ni= sample(80:150, m, replace = TRUE) # Population size in each area

# Area-level covariate
x_area= rnorm(m, 50, 10)

# True area means (unknown in practice)
beta0= 5
beta1= 0.4
sigma_u= 4

u_area= rnorm(m, 0, sigma_u)
theta_true= beta0 + beta1 * x_area + u_area # True but unknown


##################################
### 2. GENERATE POPULATION DATA
##################################

population= data.frame(
area = rep(1:m, times = Ni),
x_area = rep(x_area, times = Ni)
)

# Auxiliary unit-level covariate
population$x_unit= rnorm(sum(Ni), 0, 1)

# Unit-level model: y = θ_i + β * x_unit + error
beta_unit1= 2
sigma_e= 5

population$y= theta_true[population$area] +
beta_unit1 * population$x_unit +
rnorm(sum(Ni), 0, sigma_e)

##########################################
### 3. DRAW SAMPLE (MIMICS A REAL SURVEY)
##########################################

sample_indices= unlist(
lapply(1:m, function(a) {
which(population$area == a) |>
sample(size = sample(8:15, 1))
})
)

sample_data= population[sample_indices, ]

###############################################
### 4. COMPUTE DIRECT ESTIMATES (θ_hat_D)
###############################################

library(dplyr)

direct_df= sample_data %>%
group_by(area) %>%
summarise(
theta_direct = mean(y),
D = var(y) / n(), # sampling variance of direct estimate
x_area = mean(x_area)
) %>% data.frame()

##################################################
### 5. FIT FAY HERRIOT MODEL (AREA LEVEL)
##################################################

library(sae)

fh_model= mseFH(theta_direct ~ x_area,
data = direct_df,
vardir = D, method="ML")

FH_EB= fh_model$est$eblup # Model-based EB estimator
###############################################
### 6. UNIT-LEVEL BHF MODEL (Battese-Harter)
###############################################

# population-level mean of x_unit per area
meanxpop= population %>%
group_by(area) %>%
summarise(mean_xunit_pop = mean(x_unit))

# population-size
popnsize= population %>%
group_by(area) %>%
summarise(Ni = n())

# Fit BHF model
bhf_model= pbmseBHF(
formula = y ~ x_unit,
dom = area,
data = sample_data,
meanxpop = meanxpop,
popnsize = popnsize,
method = "REML"
)

BHF_pred= bhf_model$est$eblup

###################################################################################
# 7. FIT FAY HERRIOT MODEL (AREA LEVEL) with simple Weighted Least Square Technique
###################################################################################

# Weighted regression = Fay-Herriot fixed effects
fh_lm= lm(theta_direct ~ x_area, data = direct_df, weights = 1/direct_df$D)
#summary(fh_lm)

## The classical Fay Herriot EB shrinkage factor is: B_i= A/(A+D_i)
## where A = variance of area random effects.

## A simple estimator (Prais Winsten type):

# estimate FH random-effect variance A
A_hat= max(0, var(residuals(fh_lm)) - mean(direct_df$D))
#A_hat

## Compute EBLUP (Manually)

# model predictions
theta_model= predict(fh_lm)

# shrinkage factor
B_i= A_hat / (A_hat + direct_df$D)

# EB estimators (same as Composite Estimates )
FH_EB_lm= B_i * direct_df$theta_direct + (1 - B_i) * theta_model

#############################################################################
### 8. UNIT-LEVEL BHF MODEL (Battese-Harter) with simple Mixed effect Model
#############################################################################

bhf_lmer= lmer(y ~ x_unit + (1 | area), data = sample_data)
# Predict summary(bhf_lmer)

## Predict Small-Area Means (Manual BHF EBLUP)
# population mean of x_unit in each area
meanx_pop= population %>% group_by(area) %>%
summarise(mx = mean(x_unit))

# fixed effects
beta_hat= fixef(bhf_lmer)

# random effects
u_hat= ranef(bhf_lmer)$area[,1]

# Extract variance components
sigma_u2 = as.numeric(VarCorr(bhf_lmer)$area[1])
sigma_e2 = attr(VarCorr(bhf_lmer), "sc")^2

# Sample size per area
n_i= table(sample_data$area)

## Synthetic estimator = Model-based prediction ignoring area sample mean
synthetic_unit= beta_hat[1] + beta_hat[2] * meanx_pop$mx

# BHF predictor
BHF_pred_lme4= synthetic_unit + u_hat

############################################################
### 9. Function for calculate composite etimats' weights

############################################################

get_composite_weights= function(
method = c("FH_area", "FH_unit"),
D = NULL,
A = NULL,
sigma_u2 = NULL,
sigma_e2 = NULL,
n_i = NULL
) {
method= match.arg(method)

if (method == "FH_area") {
if (is.null(D) | is.null(A))
stop("FH_area requires D (direct variances) and A (FH random-effect variance).")

# Fay Herriot shrinkage weights
w= A / (A + D)
return(w)
}

if (method == "FH_unit") {
if (is.null(sigma_u2) | is.null(sigma_e2) | is.null(n_i))
stop("FH_unit requires sigma_u2, sigma_e2, and n_i (sample sizes per area).")

# Unit-level shrinkage weights (BHF-type)
w= sigma_u2 / (sigma_u2 + sigma_e2 / n_i)
return(w)
}
}

## Composite etimats' weights for Area-Level Fay Herriot model
w_FH_area= get_composite_weights(method = "FH_area",D = direct_df$D,A = A_hat)

## Composite etimats' weights for Unit-Level BHF model
w_FH_unit= get_composite_weights(method = "FH_unit",sigma_u2 = sigma_u2,sigma_e2 = sigma_e2,n_i = as.numeric(n_i))

 

 

##########################################
### 10. Calculating Composite Estimates
##########################################

## Composite etimates for Area-Level Fay Herriot model
theta_composite= w_FH_area * direct_df$theta_direct +(1 - w_FH_area) * predict(fh_lm)

theta_composite_unit= w_FH_unit * direct_df$theta_direct + (1 - w_FH_unit) * synthetic_unit

###########################################
## 11. GAM (Generalized Additive Model)
##########################################
# Fit GAM on sampled units
# using smooth term for x_unit and linear for area-level covariate x_area
gam_fit= gam(y ~ s(x_unit) + x_area, data = sample_data, method = "REML")

# Predict at population unit level and then compute area synthetic means
population$pred_gam= predict(gam_fit, newdata = population)

synthetic_gam_area= population %>%
group_by(area) %>%
summarise(syn_gam = mean(pred_gam)) %>%
arrange(area)

## Composite Estimates for GAM Synthetic ##
direct_and_syn_gam= direct_df %>%
arrange(area) %>%
left_join(synthetic_gam_area, by = "area")

# Composite GAM estimate
composite_gam_unit= w_FH_unit * direct_and_syn_gam$theta_direct +
(1 - w_FH_unit) * direct_and_syn_gam$syn_gam

############################################################
## 12 RANDOM FOREST SYNTHETIC + COMPOSITE ESTIMATES (RF)
############################################################

library(randomForest)
set.seed(2026)

## 1. FIT RANDOM FOREST ON SAMPLE DATA
rf_fit= randomForest(y ~ x_unit + x_area,
data = sample_data,
ntree = 500,
mtry = 2,
importance = TRUE)
#
## 2. SYNTHETIC PREDICTION FOR ENTIRE POPULATION
population$rf_pred= predict(rf_fit, newdata = population)
#
## 3. SYNTHETIC AREA-LEVEL MEANS (RF)
rf_synthetic_area= population %>%
group_by(area) %>%
summarise(RF_synthetic = mean(rf_pred)) %>%
arrange(area)
#
## 4. MERGE WITH DIRECT ESTIMATES
direct_rf= direct_df %>%
arrange(area) %>%
left_join(rf_synthetic_area, by = "area")
#
## 6. COMPOSITE ESTIMATES USING RF SYNTHETIC MEANS
Composite_RF= w_FH_unit * direct_rf$theta_direct +
(1 - w_FH_unit) * direct_rf$RF_synthetic

###############################################
### 13. COLLECT ALL RESULTS TOGETHER
###############################################

final_results= data.frame(area=1:m,
true_theta= theta_true,
direct= direct_df$theta_direct,
Fay_Herriot= FH_EB,
Fay_Herriot_area_by_WLS= FH_EB_lm,
BHF_unit_level= BHF_pred$eblup,
BHF_unit_level_by_GLMM= BHF_pred_lme4,
Composite_Estimate_Area= theta_composite,
Composite_Estimate_Unit= theta_composite_unit,
GAM_synthetic_area= direct_and_syn_gam$syn_gam,
Composite_GAM_Unit= composite_gam_unit,
RF_synthetic_area= direct_rf$RF_synthetic,
Composite_RF_Unit= Composite_RF)

####################
### 14. Plotting
####################
par(mfrow=c(2,2))
plot(final_results$area, final_results$true_theta, type = "b", pch = 1,cex=2,
ylim = range(final_results[,c("true_theta","direct","Fay_Herriot","Fay_Herriot_area_by_WLS","Composite_Estimate_Area")], na.rm=TRUE),
xlab = "Area", ylab = "Estimate", main = "Area level Fay Herriot")
lines(final_results$area, final_results$direct, type = "b", col = "blue", pch = 17)
lines(final_results$area, final_results$Fay_Herriot, type = "b", col = "purple", pch = 15)
lines(final_results$area, final_results$Fay_Herriot_area_by_WLS, type = "b", col = "green", pch = 13)
lines(final_results$area, final_results$Composite_Estimate_Area, type = "b", col = "red", pch = 11)
legend("bottomright",
legend = c("True","Direct","Fay Herriot(by sae package)","Fay Herriot(by WLS)","Area level Composite estimate"),
col = c("black","blue","purple","green","red"),
pch = c(1,17,15,13,11))

plot(final_results$area, final_results$true_theta, type = "b", pch = 1,cex=2,
ylim = range(final_results[,c("true_theta","direct","BHF_unit_level","BHF_unit_level_by_GLMM","Composite_Estimate_Unit")], na.rm=TRUE),
xlab = "Area", ylab = "Estimate", main = "Unit level BHF")
lines(final_results$area, final_results$direct, type = "b", col = "blue", pch = 17)
lines(final_results$area, final_results$BHF_unit_level, type = "b", col = "purple", pch = 15)
lines(final_results$area, final_results$BHF_unit_level_by_GLMM, type = "b", col = "green", pch = 13)
lines(final_results$area, final_results$Composite_Estimate_Unit, type = "b", col = "red", pch = 11)
legend("bottomright",
legend = c("True","Direct","BHF(by sae package)","BHF(by GLMM)","Unit level Composite estimate"),
col = c("black","blue","purple","green","red"),
pch = c(1,17,15,13,11))

plot(final_results$area, final_results$true_theta, type = "b", pch = 1,cex=2,
ylim = range(final_results[,c("true_theta","direct","GAM_synthetic_area","Composite_GAM_Unit","BHF_unit_level","Composite_Estimate_Unit")], na.rm=TRUE),
xlab = "Area", ylab = "Estimate", main = "Unit level GAM")
lines(final_results$area, final_results$direct, type = "b", col = "blue", pch = 17)
lines(final_results$area, final_results$GAM_synthetic_area, type = "b", col = "purple", pch = 15)
lines(final_results$area, final_results$Composite_GAM_Unit, type = "b", col = "green", pch = 13)
lines(final_results$area, final_results$BHF_unit_level, type = "b", col = "skyblue", pch = 11)
lines(final_results$area, final_results$Composite_Estimate_Unit, type = "b", col = "red", pch = 10)
legend("bottomright",
legend = c("True","Direct","GAM synthetic","Unit level GAM Composite Estimate","BHF (by sae package)","Unit level Composite estimate"),
col = c("black","blue","purple","green","skyblue","red"),
pch = c(1,17,15,13,11,10))

plot(final_results$area, final_results$true_theta, type = "b", pch = 1,cex=2,
ylim = range(final_results[,c("true_theta","GAM_synthetic_area","direct","RF_synthetic_area","Composite_RF_Unit","BHF_unit_level","Composite_Estimate_Unit")], na.rm=TRUE),
xlab = "Area", ylab = "Estimate", main = "Unit level Random Forest")
lines(final_results$area, final_results$direct, type = "b", col = "blue", pch = 17)
lines(final_results$area, final_results$RF_synthetic_area, type = "b", col = "purple", pch = 15)
lines(final_results$area, final_results$Composite_RF_Unit, type = "b", col = "green", pch = 13)
lines(final_results$area, final_results$BHF_unit_level, type = "b", col = "skyblue", pch = 11)
lines(final_results$area, final_results$Composite_Estimate_Unit, type = "b", col = "red", pch = 10)
legend("bottomright",
legend = c("True","Direct","Random Forest synthetic","Unit level Random Forest Composite Estimate","BHF (by sae package)","Unit level Composite estimate"),
col = c("black","blue","purple","green","skyblue","red"),
pch = c(1,17,15,13,11,10))

par(mfrow=c(1,1))

Output:

Area	True theta	Direct	Fay Herriot	Fay Herriot area by WLS	BHF unit level	BHF unit level by GLMM	Composite Estimate Area	Composite Estimate Unit	GAM Synthetic Area	Composite GAM Unit	RF Synthetic Area	Composite RF Unit
1	23.30	24.96	24.95	24.94	25.02	25.03	24.94	25.02	24.71	24.94	25.42	25.01
2	15.83	17.24	18.07	18.13	17.25	17.33	18.13	17.85	21.57	17.55	17.60	17.26
3	30.76	29.42	28.41	28.47	29.37	29.32	28.47	29.11	22.63	28.87	28.66	29.36
4	27.33	26.87	26.55	26.58	27.74	27.72	26.58	26.75	23.22	26.45	26.05	26.78
5	13.76	14.34	15.02	15.11	15.10	15.19	15.11	15.19	18.54	14.66	14.90	14.39
6	26.74	27.91	28.06	27.97	27.98	27.96	27.97	27.75	29.16	28.00	28.79	27.97
7	23.75	23.44	24.07	24.01	24.20	24.21	24.01	23.75	26.54	23.83	25.65	23.72
8	23.57	23.51	23.26	23.32	23.81	23.83	23.32	23.68	20.78	23.27	22.33	23.41
9	29.68	28.58	28.77	28.69	28.60	28.57	28.69	28.35	31.10	28.80	29.80	28.69
10	27.72	27.95	27.78	27.67	26.95	26.94	27.67	27.70	27.40	27.89	25.28	27.67
11	23.71	26.29	26.02	26.03	25.47	25.47	26.03	26.22	24.41	26.11	24.98	26.16
12	28.41	29.81	29.69	29.57	29.79	29.74	29.57	29.50	29.45	29.78	30.17	29.83
13	33.99	34.56	33.02	32.87	33.26	33.18	32.87	33.66	29.54	34.03	32.09	34.30
14	27.38	26.22	26.90	26.74	26.83	26.82	26.74	26.18	29.24	26.60	28.49	26.51
15	33.82	32.83	31.57	31.45	32.69	32.63	31.45	32.29	28.58	32.50	31.42	32.72
16	21.02	20.98	21.84	21.76	21.34	21.40	21.76	21.32	28.06	21.49	22.73	21.10
17	27.09	27.81	27.59	27.58	27.75	27.73	27.58	27.62	25.55	27.60	26.43	27.68
18	24.27	24.89	24.77	24.78	24.90	24.90	24.78	24.94	24.24	24.83	24.16	24.83
19	24.34	22.19	22.45	22.45	21.90	21.94	22.45	22.53	23.93	22.37	23.86	22.37
20	23.74	23.93	23.73	23.79	23.77	23.79	23.79	24.07	22.59	23.80	24.79	24.01