2 EBLUPs and WFQ estimators

Yong You, J.N.K. Rao and Mike Hidiroglou

Previous | Next

Suppose that we have $m$ small areas with design-unbiased direct estimators, $y_i,$ of the area means ${\theta }_i,i=1,\dots ,m.$ The FH model refers to $(y_i,{\theta }_i)$ and associated area level auxiliary variables $x_i=(x_{i1},\dots ,x_{ip}{)}^{\prime }$ with $x_{i1}=1.$ Assuming independent sampling across areas, the FH model may be written as a linear mixed model given by

$y_i={\theta }_i+e_i={x^{\prime }}_i\beta +v_i+e_i,i=1,\dots ,m,$(2.1)

where ${\theta }_i={x^{\prime }}_i\beta +v_i$ is the linking model, and $e_i$ is the sampling error with mean 0 and known variance ${\sigma }_i^2,$ which is independent of the area specific random effect $v_i.$ Sampling is independent across areas and the $v_i’\text{s}$ are assumed to be independent and identically distributed with mean 0 and variance ${\sigma }_v^2.$

The best linear unbiased predictor (BLUP) of ${\theta }_i$ under the "true� model (2.1) is given by

${\tilde{\theta }}_i={\gamma }_i{y}_i+(1-{\gamma }_i){x^{\prime }}_i\tilde{\beta },$(2.2)

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

$\tilde{\beta }={\left({\displaystyle \sum _{i=1}^m{\gamma }_i}{x}_i{x^{\prime }}_i\right)}^{-1}\left({\displaystyle \sum _{i=1}^m{\gamma }_i{x}_i{y}_i}\right),$(2.3)

see Rao (2003, page 116). The estimator ${\tilde{\theta }}_i$ depends on the unknown model variance ${\sigma }_v^2,$ and replacing ${\sigma }_v^2$ in (2.2) by a suitable estimator ${\widehat{\sigma }}_v^2,$ we get the EBLUP:

${\widehat{\theta }}_i^{\text{EBLUP}}={\widehat{\gamma }}_i{y}_i+(1-{\widehat{\gamma }}_i){x^{\prime }}_i\widehat{\beta },$(2.4)

where ${\widehat{\gamma }}_i={\widehat{\sigma }}_v^2/({\widehat{\sigma }}_v^2+{\sigma }_i^2)$ and $\widehat{\beta }$ is obtained from (2.3) by replacing ${\gamma }_i$ by ${\widehat{\gamma }}_i.$ In this paper, we use the restricted maximum likelihood (REML) estimator of ${\sigma }_v^2,$ assuming normality of $v_i$ and $e_i.$ The weighted sum ${\displaystyle {\sum }_{i=1}^m{w}_i{\widehat{\theta }}_i^{\text{EBLUP}}}$ of the EBLUPs (2.4) does not necessarily agree with the corresponding weighted direct estimator ${\displaystyle {\sum }_{i=1}^m{w}_i{y}_i}$ of the aggregate, where the $w_i’\text{s}$ are pre-specified weights such that ${\displaystyle {\sum }_{i=1}^m{w}_i{y}_i}$ is a design-consistent estimator of the aggregate (total or mean). If the gap between ${\displaystyle {\sum }_{i=1}^m{w}_i}{\widehat{\theta }}_i^{\text{EBLUP}}$ and ${\displaystyle {\sum }_{i=1}^m{w}_i{y}_i}$ is large, it may indicate some model failure that should be taken care of before proceeding to benchmarking, as noted by the Associate Editor.

An estimator of the mean squared prediction error $\text{MSPE}({\widehat{\theta }}_i^{\text{EBLUP}})=\text{E}{({\widehat{\theta }}_i^{\text{EBLUP}}-{\theta }_i)}^2$ correct to second-order terms, under REML estimation, is given by

$\text{mspe}({\widehat{\theta }}_i^{\text{EBLUP}})=g_{1i}+g_{2i}+2g_{3i},$(2.5)

where $g_{1i}={\widehat{\gamma }}_i{\sigma }_i^2$ is the leading term of order $O(1),$ and $g_{2i}$ and $g_{3i}$ are lower order terms of order $O(m^{-1})$ accounting for the variability of $\widehat{\beta }$ and ${\widehat{\sigma }}_v^2$ respectively (Rao 2003, page 128). We have

$g_{2i}={(1-{\widehat{\gamma }}_i)}^2{x^{\prime }}_i\widehat{V}(\tilde{\beta })x_i={\widehat{\sigma }}_v^2{(1-{\widehat{\gamma }}_i)}^2{x^{\prime }}_i{\left({\displaystyle \sum _{i=1}^m{\widehat{\gamma }}_i{x}_i{x^{\prime }}_i}\right)}^{-1}{x}_i$

and

$g_{3i}=2{({\sigma }_i^2)}^2{({\widehat{\sigma }}_v^2+{\sigma }_i^2)}^{-3}{\left\{{\displaystyle \sum _{i=1}^m{({\widehat{\sigma }}_v^2+{\sigma }_i^2)}^{-2}}\right\}}^{-1},$

where $\widehat{V}(\tilde{\beta })$ is the estimator of $\text{V}(\tilde{\beta })={\sigma }_v^2{\left({\displaystyle {\sum }_{i=1}^m{\gamma }_i}{x}_i{x^{\prime }}_i\right)}^{-1}.$ The MSPE estimator (2.5) is nearly unbiased in the sense that

$\text{E{mspe}({\widehat{\theta }}_i^{\text{EBLUP}})\text{}}=\text{MSPE}({\widehat{\theta }}_i^{\text{EBLUP}})+O(m^{-2}).$

WFQ obtained an EBLUP estimator, ${\widehat{\eta }}_i^{\text{EBLUP}}\text{},$ under the following augmented FH model:

$y_i={x^{\prime }}_i\delta +w_{ei}\lambda +u_i+e_i\equiv {\eta }_i+e_i\text{},$(2.6)

where the random effects $u_i$ are independent with $\text{E}(u_i)=$ 0 and $\text{Var}(u_i)={\sigma }_u^2,$ and is independent of $e_i.$ The augmenting auxiliary variable is taken as $w_{ei}=w_i{\sigma }_i^2.$ WFQ showed that the EBLUP estimator of ${\eta }_i$ under the augmented model (2.6), ${\widehat{\eta }}_i^{\text{EBLUP}}\equiv {\widehat{\theta }}_i^{\text{WFQ}},$ is self-benchmarking in the sense of satisfying

${\displaystyle \sum _{i=1}^m{w}_i{\widehat{\theta }}_i^{\text{WFQ}}}={\displaystyle \sum _{i=1}^m{w}_i{y}_i}.$

The EBLUP ${\widehat{\theta }}_i^{\text{WFQ}}$ under the augmented model (2.6) is obtained from (2.4) by changing ${x^{\prime }}_i$ to $({x^{\prime }}_i,w_{ei}),{\widehat{\sigma }}_v^2$ to ${\widehat{\sigma }}_u^2$ and $\widehat{{\beta }^{\prime }}$ to $(\widehat{{\delta }^{\prime }},\widehat{\lambda }).$ The estimator of $\text{MSPE}({\widehat{\theta }}_i^{\text{WFQ}})$ is similarly obtained from (2.5). Under the true model (2.1), the MSPE of ${\widehat{\theta }}_i^{\text{WFQ}}$ will be larger than the MSPE of ${\widehat{\theta }}_i^{\text{EBLUP}},$ but its estimator, $\text{mspe}({\widehat{\theta }}_i^{\text{WFQ}}),$ will remain nearly unbiased under the true model, as noted by the Associate Editor and a referee, because the true model is a special case of the augmented model with $\lambda =0.$

Previous | Next