Comparison of unit level and area level small area estimators 3. Area level modelComparison of unit level and area level small area estimators 3. Area level model

• Previous
• Next

The Fay-Herriot model (Fay and Herriot 1979) is a basic area level model widely used in small area estimation to improve the direct survey estimates. The Fay-Herriot model has two components, namely, a sampling model for the direct survey estimates and a linking model for the small area parameters of interest. The sampling model assumes that given the area-specific sample size $n_i>1,$ there exists a direct survey estimator ${\widehat{\theta }}_i^{\text{DIR}}.$ The direct survey estimator is design unbiased for the small area parameter ${\theta }_i.$ The sampling model is given by

$${\widehat{\theta }}_i^{\text{DIR}}={\theta }_i+e_i,i=1,\dots ,m,(3.1)$$

where the $e_i$ is the sampling error associated with the direct estimator ${\widehat{\theta }}_i^{\text{DIR}}$ and $m$ is the number of small areas. It is customary in practice to assume that the $e_i’\text{s}$ are independently normal random variables with mean $E\left(e_i\right)=0$ and sampling variance $\mathrm{var}\left(e_i\right)={\sigma }_i^2.$ The linking model is obtained by assuming that the small area parameter of interest ${\theta }_i$ is related to area level auxiliary variables $z_i=\left(z_{i1},\dots ,z_{ip}\right)\prime$ through the following linear regression model

$${\theta }_i=z_i\beta +v_i,i=1,\dots ,m,(3.2)$$

where $\beta =\left({\beta }_1,\dots ,{\beta }_p\right)\prime$ is a $p\times 1$ vector of regression coefficients, and the $v_i’\text{s}$ are area-specific random effects assumed to be $i.i.d.$ with $E\left(v_i\right)=0$ and $\mathrm{var}\left(v_i\right)={\sigma }_v^2.$ The assumption of normality is generally also made, even though it is more difficult to justify the assumption. This assumption is needed to obtain the MSE estimation. The model variance ${\sigma }_v^2$ is unknown and needs to be estimated from the data. The area level random effect $v_i$ capture the unstructured heterogeneity among areas that is not explained by the sampling variances. Combining models (3.1) and (3.2) leads to a linear mixed area level model given by

$${\widehat{\theta }}_i^{\text{DIR}}=z_i\beta +v_i+e_i.(3.3)$$

Model (3.3) involves both design-based random errors $e_i$ and model-based random effects $v_i.$ For the Fay-Herriot model, the sampling variance ${\sigma }_i^2$ is assumed to be known in model (3.3). This is a very strong assumption. Generally smoothed estimators of the sampling variances are used in the Fay-Herriot model and then ${\sigma }_i^2’\text{s}$ are treated as known. However, if direct estimators of sampling variances are used in the Fay-Herriot model, an extra term needs to be added to the MSE estimator to account for the extra variation (Wang and Fuller 2003).

Assuming that the model variance ${\sigma }_v^2$ is known, the best linear unbiased predictor (BLUP) of the small area parameter ${\theta }_i$ can be obtained as

$${\tilde{\theta }}_i={\gamma }_i{\widehat{\theta }}_i^{\text{DIR}}+\left(1-{\gamma }_i\right)z_i{\tilde{\beta }}_{\text{WLS}},(3.4)$$

where ${\gamma }_i={\sigma }_v^2/\left({\sigma }_v^2+{\sigma }_i^2\right),$ and ${\tilde{\beta }}_{\text{WLS}}$ is the weighted least squared (WLS) estimator of $\beta$ given by

$${\tilde{\beta }}_{\text{WLS}}={\left[{\displaystyle \sum _{i=1}^m{\left({\sigma }_i^2+{\sigma }_v^2\right)}^{-1}}{z}_i{z}_i{}^{\prime }\right]}^{-1}\left[{\displaystyle \sum _{i=1}^m{\left({\sigma }_i^2+{\sigma }_v^2\right)}^{-1}{z}_i{y}_i}\right]={\left[{\displaystyle \sum _{i=1}^m{\gamma }_i}{z}_i{z}_i{}^{\prime }\right]}^{-1}\left[{\displaystyle \sum _{i=1}^m{\gamma }_i{z}_i{y}_i}\right].$$

There are several methods available to estimate the unknown model variance ${\sigma }_v^2;$ You (2010) provides a review of these methods. We chose the restricted maximum likelihood (REML) obtained by Cressie (1992) to estimate the model variance under the Fay-Herriot model. Using the scoring algorithm, the REML estimator ${\widehat{\sigma }}_v^2$ is obtained as

$${\sigma }_v^{2\left(k+1\right)}={\sigma }_v^{2\left(k\right)}+{\left[I_R\left({\sigma }_v^{2\left(k\right)}\right)\right]}^{-1}{S}_R\left({\sigma }_v^{2\left(k\right)}\right),\text{ for }k=1,2,\dots ,$$

where $I_R\left({\sigma }_v^2\right)=1/2\text{tr}\left[PP\right],$ and $S_R\left({\sigma }_v^2\right)=1/2y^{\prime }PPy-1/2\text{tr}\left[P\right],$ and $P=V^{-1}-V^{-1}Z{\left(Z^{\prime }{V}^{-1}Z\right)}^{-1}{Z}^{\prime }{V}^{-1}.$ Using a guessing value for ${\sigma }_v^{2\left(1\right)}$ as the starting value, the algorithm converges very fast.

Replacing ${\sigma }_v^2$ in equation (3.4) by the REML estimator ${\widehat{\sigma }}_v^2,$ we obtain the EBLUP of the small area parameter ${\theta }_i$ based on the Fay-Herriot model as

$${\widehat{\theta }}_i^{\text{FH}}={\widehat{\gamma }}_i{\widehat{\theta }}_i^{\text{DIR}}+\left(1-{\widehat{\gamma }}_i\right)z_i{\widehat{\beta }}_{\text{WLS}},(3.5)$$

where ${\widehat{\gamma }}_i={\widehat{\sigma }}_v^2/\left({\widehat{\sigma }}_v^2+{\sigma }_i^2\right).$ The MSE estimator of ${\widehat{\theta }}_i^{\text{FH}}$ is given by (see Rao 2003)

$$\text{mse}\left({\widehat{\theta }}_i^{\text{FH}}\right)=g_{1i}+g_{2i}+2g_{3i},(3.6)$$

where $g_{1i}$ is the leading term, $g_{2i}$ accounts for the variability due to estimation of the regression parameter $\beta ,$ and $g_{3i}$ is due to the estimation of the model variance. These $g-$ terms are defined as follow:

$$g_{1i}={\widehat{\gamma }}_i{\sigma }_i^2,g_{2i}={\left(1-{\widehat{\gamma }}_i\right)}^2{z}_i\text{var}\left({\widehat{\beta }}_{\text{WLS}}\right)z_i={\widehat{\sigma }}_v^2{\left(1-{\widehat{\gamma }}_i\right)}^2{z}_i{\left({\displaystyle \sum _{i=1}^m{\widehat{\gamma }}_i{z}_i{z}_i}\right)}^{-1}{z}_i$$

and $g_{3i}={\left({\sigma }_i^2\right)}^2{\left({\widehat{\sigma }}_v^2+{\sigma }_i^2\right)}^{-3}var\left({\widehat{\sigma }}_v^2\right).$

The estimated variance of ${\widehat{\sigma }}_v^2$ is given by $var\left({\widehat{\sigma }}_v^2\right)=2{\left({\displaystyle {\sum }_{i=1}^m{\left({\widehat{\sigma }}_v^2+{\sigma }_i^2\right)}^{-2}}\right)}^{-1};$ see Datta and Lahiri (2000).

Up to now we have assumed that the sampling variance ${\sigma }_i^2$ is assumed known in the Fay-Herriot model (3.3). This is a very strong assumption. Usually a direct survey estimator, say $s_i^2,$ of the sampling variance ${\sigma }_i^2$ is available. As these estimated variances can be quite variable, they are smoothed using external models and generalized variance functions: these smoothed variances are denoted as ${\tilde{s}}_i^2.$ The smoothed sampling variance estimates ${\tilde{s}}_i^2$ are used in the Fay-Herriot model and treated as known. The associated $\text{mse}({\widehat{\theta }}_i^{\text{FH}})$ is obtained by replacing ${\sigma }_i^2$ by ${\tilde{s}}_i^2$ in equation (3.6). Rivest and Vandal (2003) and Wang and Fuller (2003) considered the small area estimation using the Fay-Herriot model with the direct sampling variance estimates $s_i^2$ under the assumption that the estimators $s_i^2$ are independent of the direct survey estimators $y_i$ and $d_i{s}_i^2\sim {\sigma }_i^2{\chi }_{d_i}^2,$ where $d_i=n_i-1$ and $n_i$ is the sample size for the $i^{\text{th}}$ area. When the direct sampling variance estimate $s_i^2$ is used in the place of the true sampling variance ${\sigma }_i^2,$ an extra term accounts for the uncertainty of using $s_i^2$ is needed in the MSE estimator (3.6), and this term, denoted as $g_{4i},$ is given by

$$g_{4i}=\frac{4}{n_i-1}\frac{{\widehat{\sigma }}_v^4{s}_i^4}{{\left({\widehat{\sigma }}_v^2+s_i^2\right)}^3};$$

see Rivest and Vandal (2003) and Wang and Fuller (2003) for details.

To apply the Fay-Herriot model, we need to obtain area level direct estimates and the corresponding sampling variance estimates as input values for the Fay-Herriot model. We consider three area level direct estimators; namely, the direct sample mean estimator assuming simple random sampling (SRS), the Horvitz-Thompson estimator (HT), and the weighted Hájek estimator (HA). The weighted Hájek estimator is also used in the pseudo-EBLUP estimator for the unit level model denoted as ${\overline{y}}_{iw}$ in equation (2.7). Table 3.1 presents these three area level direct estimators and the corresponding sampling variance estimators.

These area level estimators are used as input values into the Fay-Herriot model. Correspondingly, the three area level model-based estimators are denoted as: FH-SRS, FH-HT, and FH-HA. That is, we replace ${\widehat{\theta }}_i^{\text{DIR}}$ by ${\widehat{\theta }}_i^{\text{SRS}},$ ${\widehat{\theta }}_i^{\text{HT}}$ or ${\widehat{\theta }}_i^{\text{HA}}$ in (3.5) and obtain the corresponding model-based estimator ${\widehat{\theta }}_i^{\text{FH-SRS}},$ ${\widehat{\theta }}_i^{\text{FH-HT}}$ and ${\widehat{\theta }}_i^{\text{FH-HA}}.$ The SRS direct estimator ${\widehat{\theta }}_i^{\mathrm{SRS}}$ ignores the sample design and is not design consistent, unless the sample design is based on simple random sampling. Note that ${\widehat{\theta }}_i^{\text{HT}}$ and ${\widehat{\theta }}_i^{\text{HA}}$ are design consistent estimators. It follows that the corresponding model-based estimators ${\widehat{\theta }}_i^{\text{FH-HT}}$ and ${\widehat{\theta }}_i^{\text{FH-HA}}$ are design consistent as the sample size increases. Furthermore, this means that these estimators are robust to model misspecification.

In the next section, we compare the unit level model with the Fay-Herriot model through a simulation study. The statistics used for these comparisons are bias, relative root MSE and confidence intervals of the model-based estimators.

• Previous
• Next