2. Estimation of small area means

Isabel Molina, J.N.K. Rao and Gauri Sankar Datta

Previous | Next

Consider a population partitioned into $m$  areas and let ${\theta }_i$  be the mean of the variable of interest for area $i,i=1,\dots ,m.$  We assume that a sample is drawn independently from each area. Let $y_i$  be a design-unbiased direct estimator of ${\theta }_i$  obtained using survey data from the sampled area $i.$  Direct estimators are very inefficient for areas with small sample sizes. We study small area estimation under an area level model, in which the values of area level covariates are available for all areas. The basic model of this type is the Fay-Herriot model, introduced by Fay and Herriot (1979), to estimate per capita income for small places in the United States. This model consists of two parts. The first part assumes that direct estimators, $y_i,$  of small area means, ${\theta }_i,$  are design unbiased, satisfying

$$y_i={\theta }_i+e_i\mathrm{,}{e}_i\stackrel{\text{ind}}{\sim }N\left(\mathrm{0,}{D}_i\right)\mathrm{,}i=1,\dots \mathrm{,}m\mathrm{.}(2.1)$$

Here, the sampling variance $D_i\text{=Var}\left(y_i|{\theta }_i\right)$  is assumed to be known for all areas $i=1,\dots \mathrm{,}m.$  In practice, the $D_i’\text{s}$  are ascertained from external sources or by smoothing the estimated sampling variances using a generalized variance function method (Fay and Herriot 1979).

In the second part, the Fay-Herriot model treats ${\theta }_i$  as random and assumes that a $p-$vector of area level covariates, $x_i,$  linearly related to ${\theta }_i,$  is available for each area $i,$  i.e.,

$${\theta }_i={x^{\prime }}_i\beta +v_i\mathrm{,}{v}_i\stackrel{\text{iid}}{\sim }N\left(\mathrm{0,}A\right)\mathrm{,}i=1,\dots \mathrm{,}m\mathrm{,}(2.2)$$

where $v_i$  is the random effect of area $i,$  assumed to be independent of $e_i$  and $A\ge 0$  is the variance of the random effects. Observe that marginally,

$$y_i\stackrel{\text{ind}}{\sim }N\left({x^{\prime }}_i\beta ,D_i+A\right)\mathrm{,}i=1,\dots \mathrm{,}m.(2.3)$$

Letting $y={\left(y_1\mathrm{,}\dots \mathrm{,}{y}_m\right)}^{\prime },X={\left(x_1\mathrm{,}\dots \mathrm{,}{x}_m\right)}^{\prime }$  and $D=\text{diag}\left(D_1\mathrm{,}\dots \mathrm{,}{D}_m\right),$  model (2.3) may be expressed in matrix notation as $y\sim N\left\{X\beta \mathrm{,}\Sigma \left(A\right)\right\}$  with $\Sigma \left(A\right)=D+A{I}_m,$  where $I_m$  denotes the $m\times m$  identity matrix. If $A$  is known, the componentwise best linear unbiased predictor (BLUP) of $\theta ={\left({\theta }_1\mathrm{,}\dots \mathrm{,}{\theta }_m\right)}^{\prime }$  is given by

$$\tilde{\theta }\left(A\right)={\left({\tilde{\theta }}_1\left(A\right)\mathrm{,}\dots \mathrm{,}{\tilde{\theta }}_m\left(A\right)\right)}^{\prime }=X\tilde{\beta }\left(A\right)+A{\Sigma }^{-1}\left(A\right)\left\{y-X\tilde{\beta }\left(A\right)\right\}\mathrm{,}(2.4)$$

where

$$\begin{array}{lll}\tilde{\beta }\left(A\right)\hfill & =\hfill & {\left\{X^{\prime }{\Sigma }^{-1}\left(A\right)X\right\}}^{-1}{X}^{\prime }{\Sigma }^{-1}\left(A\right)y\hfill \\ \hfill & =\hfill & {\left\{{\displaystyle \sum _{i=1}^m}{\left(A+D_i\right)}^{-1}{x}_i{x^{\prime }}_i\right\}}^{-1}{\displaystyle \sum _{i=1}^m}{\left(A+D_i\right)}^{-1}{x}_i{y}_i(2.5)\hfill \end{array}$$

is the weighted least squares (WLS) estimator of $\beta .$  In practice, however, $A$  is not known. Substituting a consistent estimator $\widehat{A}$  for $A$  in the BLUP (2.4), we get the EBLUP given by

$$\widehat{\theta }={\left({\widehat{\theta }}_1\mathrm{,}\dots \mathrm{,}{\widehat{\theta }}_m\right)}^{\prime }=X\widehat{\beta }+\widehat{A}{\widehat{\Sigma }}^{-1}\left(y-X\widehat{\beta }\right)\mathrm{,}(2.6)$$

where $\widehat{\beta }=\tilde{\beta }\left(\widehat{A}\right)$  and $\widehat{\Sigma }=D+\widehat{A}{I}_m.$  For the $i^{\text{th}}$  area, the EBLUP of ${\theta }_i$  can be expressed as a convex linear combination of the regression-synthetic estimator ${x^{\prime }}_i\widehat{\beta }$  and the direct estimator $y_i,$  as

$${\widehat{\theta }}_i=B_i\left(\widehat{A}\right){x^{\prime }}_i\widehat{\beta }+\left\{1-B_i\left(\widehat{A}\right)\right\}{y}_i\mathrm{,}(2.7)$$

where the weight attached to the regression-synthetic estimator ${x^{\prime }}_i\widehat{\beta }$  is given by $B_i\left(\widehat{A}\right),$  where $B_i\left(A\right)=D_i/\left(A+D_i\right).$  Observe that the weight increases with the sampling variance $D_i.$  Thus, when the direct estimator is not reliable, i.e., $D_i$  is large as compared with the total variance $\widehat{A}+D_i,$  more weight is attached to the regression-synthetic estimator ${x^{\prime }}_i\widehat{\beta }.$  On the other hand, when the direct estimator is efficient, $D_i$  is small relative to $\widehat{A}+D_i,$  and then more weight is given to the direct estimator $y_i.$

Several estimators of $A$  have been proposed in the literature including moment estimators without normality assumption, ML estimator and restricted (or residual) ML estimator (REML) estimator. The ML estimator of $A$  is ${\widehat{A}}_{\text{ML}}=\text{m}\text{a}\text{x}\left(\mathrm{0,}{\widehat{A}}_{\text{ML}}^{*}\right),$  where ${\widehat{A}}_{\text{ML}}^{*}$  can be obtained by maximizing the profile likelihood function given by

$$L_P\left(A\right)=c{\left|\Sigma \left(A\right)\right|}^{-1/2}\text{e}\text{x}\text{p}\left\{-\frac{1}{2}{y}^{\prime }P\left(A\right)y\right\}\mathrm{,}$$

where $c$  denotes a generic constant and

$$P\left(A\right)={\Sigma }^{-1}\left(A\right)-{\Sigma }^{-1}\left(A\right)X{\left\{X^{\prime }{\Sigma }^{-1}\left(A\right)X\right\}}^{-1}{X}^{\prime }{\Sigma }^{-1}\left(A\right)\mathrm{.}$$

The REML estimator of $A$  is ${\widehat{A}}_{\text{RE}}=\text{m}\text{a}\text{x}\left(\mathrm{0,}{\widehat{A}}_{\text{RE}}^{*}\right),$  where ${\widehat{A}}_{\text{RE}}^{*}$  is obtained by maximizing the restricted/residual likelihood, given by

$$L_{\text{RE}}\left(A\right)=c{\left|X^{\prime }{\Sigma }^{-1}\left(A\right)X\right|}^{-1/2}{\left|\Sigma \left(A\right)\right|}^{-1/2}\text{e}\text{x}\text{p}\left\{-\frac{1}{2}{y}^{\prime }P\left(A\right)y\right\}\mathrm{.}$$

In this paper, we focus on the REML estimator ${\widehat{A}}_{\text{RE}}$  which is frequently used in practice, and we denote by ${\widehat{\theta }}_{\text{RE}}={\left({\widehat{\theta }}_{\text{RE}\mathrm{,1}}\mathrm{,}\dots \mathrm{,}{\widehat{\theta }}_{\text{RE}\mathrm{,}m}\right)}^{\prime }$  the EBLUP given in (2.6) obtained with $\widehat{A}={\widehat{A}}_{\text{RE}}.$

Previous | Next