Bayesian spatial models for estimating means of sampled and non-sampled small areas
Section 6. Conclusions

• Previous

In this paper, we followed a Bayesian approach to investigate four spatial random-effects models as alternatives to the independent Fay-Herriot model to estimate small area means. In particular, we considered four spatial models with different autocorrelation structures. We further extended the spatial models to allow multiple small areas without any direct estimates in predicting small area means for all the areas. For a class of noninformative priors, we established the propriety of posterior densities of the proposed models for both setups.

A simulation study in Section 4 illustrates that prediction accuracy can be greatly improved by considering spatial models when effective covariates are unavailable. Datta, Hall and Mandal (2011) noted that the prediction accuracy of small area estimation models largely depends on the availability of good covariates. In other words, when suitable covariates are unavailable, the independent Fay-Herriot model may not provide a significant advantage over direct estimates. The simulation results indicated that, in such cases, the spatial models considerably increase the prediction accuracy by exploiting information from adjacent areas.

We applied the proposed spatial random-effects models to estimate four-person family median incomes. Even when a good covariate exists, the spatial models exhibited noticeable improvements in terms of mean squared deviation and average posterior standard deviations. When a good covariate is unavailable, spatial models provided significantly more accurate median income predictions with much smaller variability, which agrees with the simulation results. Furthermore, the SAR and LCAR models provide more precise small area estimates when direct estimates of some states are excluded in model fitting.

In summary, the spatial models considered in this paper outperform the independent Fay-Herriot model. A significant improvement can be expected when effective covariates are unavailable. Since useful covariates are not always available, the utility of the proposed models in small area estimation can be substantial. Our simulation study and real data analysis demonstrate no clear winner among the proposed models. Nonetheless, the SAR and LCAR models show better performance compared with other spatial models. Also, the LCAR model performs robustly well with simulated data from the SAR model and real data with unknown spatial dependence. Thus, in the context of real applications where true dependency is unknown, we recommend the LCAR model.

 This work assumes that all areas have at least one neighborhood. In real applications, however, there are many situations that data contain small areas with no neighborhood (stand-alone areas). Although the proposed models can accommodate stand-alone areas by adjusting the diagonal entries of the precision matrices as in Brown, Datta and Lazar (2017), we find that this approach results in a counterintuitive prior, where stand-alone areas have smaller prior random effect variances than areas with neighborhoods. Also, we find that this prior can considerably deteriorate predictions of stand-alone areas. This is a practically important problem as many countries have islands, and this will possibly be our future research to pursue.

Disclaimer and acknowledgements

This report is released to inform interested parties of ongoing research and to encourage discussion. The views expressed on statistical, methodological, technical, or operational issues are those of the authors and not those of the U.S. Census Bureau, or the University of Georgia. The authors are grateful to Dr. William R. Bell for his insightful comments on an earlier version of this work that led to an improved manuscript.

Appendix

A.   Proof of the propriety of the posterior pdf

Proof of Theorem 1. For convenience of notation, we denote ${\Omega }_k\left(\rho \right)$ by ${\Omega }_k,$ and, for a given square matrix $A,$ the determinant of $A$ is denoted by $\left|A\right|.$ We use $K$ to denote a generic positive constant, not depending on the variables we are integrating out.

Let $m_1\ge 0$ be the number of small areas with no direct estimates and let $m_2=m-m_1.$ Also, let $Y_{\left(2\right)}$ be the $m_2\times 1$ vector with direct estimates of the sampled small areas. Without loss of generality, we assume that ${\theta }_1\mathrm{,}\dots \mathrm{,}{\theta }_m$ are arranged so that $\theta ={\left({\theta }_{\left(1\right)}^T\mathrm{,}{\theta }_{\left(2\right)}^T\right)}^T.$ Let $D_{\left(2\right)}={\left\{D_i\right\}}_{i\mathrm{<mo>=</mo>}{m}_1+1}^m$ be the diagonal matrix with sampling variances corresponding to the components of $Y_{\left(2\right)}$ and $\delta ={{\displaystyle \max }}_{m_1<i\le m}{D}_i<\infty .$

The joint pdf of $Y_{\left(2\right)}\mathrm{,}\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2$ and $\rho$ is given by

$$f\left(y_{\left(2\right)}\mathrm{,}\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2\mathrm{,}\rho \right)=N_{m_2}\left(y_{\left(2\right)}|{\theta }_{\left(2\right)}\mathrm{,}{D}_{\left(2\right)}\right)N_m\left(\theta |X\beta \mathrm{,}{\sigma }_v^2{\Omega }^{-1}\right)g\left({\sigma }_v^2\right)h\left(\rho \right)\mathrm{,}(\text{A}.1)$$

where $N_{m_2}\left(y_{\left(2\right)}|{\theta }_{\left(2\right)}\mathrm{,}{D}_{\left(2\right)}\right)$ is the normal pdf with mean ${\theta }_{\left(2\right)}$ and covariance matrix $D_{\left(2\right)}.$ The posterior pdf $\pi \left(\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2\mathrm{,}\rho |y_{\left(2\right)}\right)$ will be proper if and only if the function $f\left(y_{\left(2\right)}\mathrm{,}\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2\mathrm{,}\rho \right)$ is integrable with respect to $\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2$ and $\rho .$ Since

$$N_{m_2}\left(y_{\left(2\right)}|{\theta }_{\left(2\right)}\mathrm{,}{D}_{\left(2\right)}\right)\le K\exp \left\{-\frac{1}{2\delta }{\left(y_{\left(2\right)}-{\theta }_{\left(2\right)}\right)}^T\left(y_{\left(2\right)}-{\theta }_{\left(2\right)}\right)\right\}\mathrm{,}$$

we have from (A.1)

$$\begin{array}{ll}f\left(y_{\left(2\right)}\mathrm{,}\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2\mathrm{,}\rho \right)\hfill & \le K\exp \left\{-\frac{1}{2\delta }{\left(y_{\left(2\right)}-{\theta }_{\left(2\right)}\right)}^T\left(y_{\left(2\right)}-{\theta }_{\left(2\right)}\right)\right\}{N}_m\left(\theta |X\beta \mathrm{,}{\sigma }_v^2{\Omega }^{-1}\right)g\left({\sigma }_v^2\right)h\left(\rho \right)\hfill \\ \hfill & =K{\displaystyle \int }\exp \left\{-\frac{1}{2\delta }{\left(y-\theta \right)}^T\left(y-\theta \right)\right\}d{y}_{\left(1\right)}{N}_m\left(\theta |X\beta \mathrm{,}{\sigma }_v^2{\Omega }^{-1}\right)g\left({\sigma }_v^2\right)h\left(\rho \right)\mathrm{,}(\text{A}.2)\hfill \end{array}$$

where $y={\left(y_{\left(1\right)}^T\mathrm{,}{y}_{\left(2\right)}^T\right)}^T.$ By integrating both sides of (A.2) with respect to $\theta ,$ we get

$${\displaystyle \int }f\left(y_{\left(2\right)}\mathrm{,}\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2\mathrm{,}\rho \right)d\theta \le Kg\left({\sigma }_v^2\right)h\left(\rho \right){\displaystyle \int }{N}_m\left(y|X\beta \mathrm{,}\delta I_m+{\sigma }_v^2{\Omega }^{-1}\right)d{y}_{\left(1\right)}\mathrm{.}(\text{A}.3)$$

Partition $X$ as $X={\left[X_1\mathrm{,}{X}_2\right]}^T,$ where $X_1^T$ is $m_1\times p$ and $X_2^T$ is $m_2\times p.$ We assume that rank $\left(X_2\right)=p.$ Let $d={\left(0_{m_1}^T\mathrm{,}{y}_{\left(2\right)}^T\right)}^T,$ $\varphi ={\left(y_{\left(1\right)}^T\mathrm{,}{\beta }^T\right)}^T$ and

$$G=\left[\begin{array}{ll}-I_{m_1}\hfill & X_1^T\hfill \\ 0_{m_2\mathrm{,}{m}_1}\hfill & X_2^T\hfill \end{array}\right]\mathrm{.}$$

Then, we can write

$$y-X\beta =d-G\varphi \mathrm{,}$$

where $G$ is $m\times \left(m_1+p\right),$ $\varphi$ is $\left(m_1+p\right)\times 1.$ Hence, (A.3) can be written as

$${\displaystyle \int }f\left(y_{\left(2\right)}\mathrm{,}\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2\mathrm{,}\rho \right)d\theta \le Kg\left({\sigma }_v^2\right)h\left(\rho \right){\displaystyle \int }{N}_m\left(d|G\varphi \mathrm{,}\delta I_m+{\sigma }_v^2{\Omega }^{-1}\right)d{y}_{\left(1\right)}\mathrm{.}(\text{A}.4)$$

By integrating both sides of (A.4) with respect to $\beta ,$ we get

$${\displaystyle \int }f\left(y_{\left(2\right)}\mathrm{,}\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2\mathrm{,}\rho \right)d\theta d\beta \le Kg\left({\sigma }_v^2\right)h\left(\rho \right){\displaystyle \int }{N}_m\left(d|G\varphi \mathrm{,}\delta I_m+{\sigma }_v^2{\Omega }^{-1}\right)d\varphi \mathrm{.}(\text{A}.5)$$

Since rank $\left(X_2\right)=p,$ we immediately get that rank $\left(G\right)=m_1+p.$ Thus $G$ has full column rank. We denote $m_1+p$ by $q.$ For $k=\mathrm{2,}\dots \mathrm{,}5,$ we now derive upper bounds for

$${\left|\delta I_m+{\sigma }_v^2{\Omega }_k^{-1}\right|}^{-1/2}{\displaystyle \int }\exp \left\{-\frac{1}{2}{\left(d-G\varphi \right)}^T{\left(\delta I_m+{\sigma }_v^2{\Omega }_k^{-1}\right)}^{-1}\left(d-G\varphi \right)\right\}d\varphi$$

which will be integrable with respect to ${\sigma }_v^2$ and $\rho$ as follows.

A.1   Details for the SCAR Model

We first consider the CAR model where $k=2.$ Let $P_W$ be an orthogonal matrix such that $P_W^{T}W{P}_W=\text{diag}{\left\{{\lambda }_i\right\}}_{i\mathrm{<mo>=</mo>}1}^m=\Lambda .$ Then ${\Omega }_2{\left(\rho \right)}^{-1}=P_W{\left\{I-\rho \Lambda \right\}}^{-1}{P}_W^T,$ and hence,

$$\begin{array}{ll}{\left(d-G\varphi \right)}^T{\left(\delta I_m+{\sigma }_v^2{\Omega }_2^{-1}\right)}^{-1}\left(d-G\varphi \right)\hfill & \mathrm{<mo>=</mo>}{\left(P_W^{T}d-P_W^{T}G\varphi \right)}^T{\left(\delta I_m+{\sigma }_v^2{\mathrm{<mo>\{</mo>}I-\rho \Lambda \mathrm{<mo>\}</mo>}}^{-1}\right)}^{-1}\left(P_W^{T}d-P_W^{T}G\varphi \right)\hfill \\ \hfill & \mathrm{<mo>=</mo>}{\left(d_{*}-G_{*}\varphi \right)}^T{\left(\delta I_m+{\sigma }_v^2{\left\{I-\rho \Lambda \right\}}^{-1}\right)}^{-1}\left(d_{*}-G_{*}\varphi \right)\mathrm{,}\hfill \end{array}$$

where $d_{*}=P_W^{T}d$ and $G_{*}=P_W^{T}G.$ Suppose the rows of $G_{*}$ corresponding to distinct $q$ indices $\left\{i_1\mathrm{,}\dots \mathrm{,}{i}_q\right\}\subseteq \left\{\mathrm{1,}\dots \mathrm{,}m\right\}$ are linearly independent. We denote these rows by $g_{i_k^{*}}^T,$ $k=\mathrm{1,}\dots \mathrm{,}q.$ Let $A$ be the $q\times q$ non-singular matrix ${\left[g_{i_1^{*}}\mathrm{,}\dots \mathrm{,}{g}_{i_q^{*}}\right]}^T$ and $\eta ={\left({\eta }_1\mathrm{,}\dots \mathrm{,}{\eta }_q\right)}^T=A\varphi .$ Note that

$${\left(d-G\varphi \right)}^T{\left(\delta I_m+{\sigma }_v^2{\Omega }_2^{-1}\right)}^{-1}\left(d-G\varphi \right)\ge {\displaystyle \sum _{k\mathrm{<mo>=</mo>}1}^q}\frac{{\left(d_{i_k^{*}}-{\eta }_{i_k}\right)}^2}{\delta +{\sigma }_v^2{\left(1-\rho {\lambda }_{i_k}\right)}^{-1}}\mathrm{<mi>\; </mi>.}$$

From this, we get that

$$\begin{array}{ll}{\displaystyle \int }\exp \left\{-\frac{1}{2}{\left(d-G\varphi \right)}^T{\left(\delta I_m+{\sigma }_v^2{\Omega }_2^{-1}\right)}^{-1}\left(d-G\varphi \right)\right\}d\varphi \hfill & \le {\displaystyle \int }\exp \left\{-\frac{1}{2}{\displaystyle \sum _{k\mathrm{<mo>=</mo>}1}^q}\frac{{\left(d_{i_k^{*}}-{\eta }_{i_k}\right)}^2}{\delta +{\sigma }_v^2{\left(1-\rho {\lambda }_{i_k}\right)}^{-1}}\right\}d\varphi \hfill \\ \hfill & \mathrm{<mo>=</mo>}{\displaystyle \int }\exp \left\{-\frac{1}{2}{\displaystyle \sum _{k\mathrm{<mo>=</mo>1}}^q}\frac{{\left(d_{i_k^{*}}-{\eta }_{i_k}\right)}^2}{\delta +{\sigma }_v^2{\left(1-\rho {\lambda }_{i_k}\right)}^{-1}}\right\}d\eta {\left|A^T\text{}A\right|}^{-1/2}\hfill \\ \hfill & \mathrm{<mo>=</mo>}K{\displaystyle \prod _{k\mathrm{<mo>=</mo>1}}^q}{\left\{\delta +{\sigma }_v^2{\left(1-\rho {\lambda }_{i_k}\right)}^{-1}\right\}}^{1/2}\mathrm{,}(\text{A}.6)\hfill \end{array}$$

where $K>0$ is a finite generic constant. Also, we know that

$${\left|\delta I_m+{\sigma }_v^2{\Omega }_2^{-1}\right|}^{-1/2}={\displaystyle \prod _{i\mathrm{<mo>=</mo>}1}^m}{\left\{\delta +{\sigma }_v^2{\left(1-\rho {\lambda }_i\right)}^{-1}\right\}}^{-1/2}\mathrm{.}(\text{A}.7)$$

By (A.6) and (A.7), we get

$$\begin{array}{l}{\left|\delta I_m+{\sigma }_v^2{\Omega }_2^{-1}\right|}^{-1/2}{\displaystyle \int }\exp \left\{-\frac{1}{2}{\left(d-G\varphi \right)}^T{\left(\delta I_m+{\sigma }_v^2{\Omega }_2^{-1}\right)}^{-1}\left(d-G\varphi \right)\right\}d\varphi \hfill \\ \le K{\displaystyle \prod _{i\notin \mathrm{<mo>\{</mo>}{i}_1\mathrm{,}\dots \mathrm{,}{i}_q\mathrm{<mo>\}</mo>}}}{\left\{\delta +{\sigma }_v^2{\left(1-\rho {\lambda }_i\right)}^{-1}\right\}}^{-1/2}\hfill \\ \le K\left\{I\left({\sigma }_v^2<N\right)+{\left({\sigma }_v^2\right)}^{-\left(m-q\right)/2}{\displaystyle \prod _{i\notin \left\{i_1\mathrm{,}\dots \mathrm{,}{i}_q\right\}}}{\left(1-\rho {\lambda }_i\right)}^{1/2}I\left({\sigma }_v^2>N\right)\right\}(\text{A}.8)\hfill \end{array}$$

for any positive number $N.$ Recall that ${\lambda }_m^{-1}<\rho \mathrm{<mo><</mo>}{\lambda }_1^{-1}.$ We know $1-\rho {\lambda }_i$ is an eigenvalue of ${\Omega }_2.$ Thus, for ${\lambda }_m^{-1}<\rho <{\lambda }_1^{-1},$ for $i=\mathrm{1,}\dots ,m,$ $1-\rho {\lambda }_i>0.$ Also, ${\displaystyle {\sum }_{i\mathrm{<mo>=</mo>}1}^m}\left(1-\rho {\lambda }_i\right)=m.$ These imply that $0<1-\rho {\lambda }_i<m.$ Then from (A.8), we get

$$\begin{array}{l}{\left|\delta I_m+{\sigma }_v^2{\Omega }_2^{-1}\right|}^{-1/2}{\displaystyle \int }\exp \left\{-\frac{1}{2}{\left(d-G\varphi \right)}^T{\left(\delta I_m+{\sigma }_v^2{\Omega }_2^{-1}\right)}^{-1}\left(d-G\varphi \right)\right\}d\varphi \hfill \\ \le K\left\{I\left({\sigma }_v^2<N\right)+{\left({\sigma }_v^2\right)}^{-\left(m-q\right)/2}I\left({\sigma }_v^2>N\right)\right\}\mathrm{.}(\text{A}.9)\hfill \end{array}$$

From (A.5) and (A.9), it follows under the conditions of the theorem that the desired integral ${\displaystyle \int }f\left(y_{\left(2\right)}\mathrm{,}\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2\mathrm{,}\rho \right)d\theta d\beta d{\sigma }_v^{2}d\rho$ is finite.

A.2   Details for the SAR Model

We now consider $k=3$ for the SAR model. With $W_{*}=L^{-1/2}W{L}^{-1/2},$ we have

$$\begin{array}{ll}{\Omega }_3\hfill & ={\left(I_m-\rho \tilde{W}\right)}^T\left(I_m-\rho \tilde{W}\right)\hfill \\ \hfill & ={\left(L-\rho W\right)}^T{L}^{-2}\left(L-\rho W\right)\hfill \\ \hfill & =L^{1/2}\left(I_m-\rho W_{*}\right)L^{-1}\left(I_m-\rho W_{*}\right)L^{1/2}\mathrm{.}\hfill \end{array}$$

First, $\text{tr}{\Omega }_3=m+{\rho }^2{\displaystyle {\sum }_i}{\displaystyle {\sum }_j}{\tilde{w}}_{ij}^2\le m+{\rho }^2{\displaystyle {\sum }_i}{\displaystyle {\sum }_j}{\tilde{w}}_{ij}=m+{\rho }^{2}m<2m$ since $0\le {\tilde{w}}_{ij}\le 1,$ ${\displaystyle {\sum }_j}{\tilde{w}}_{ij}=1,$ and $-1<\rho <1.$

Note that the eigenvalues ${\nu }_1\mathrm{,}\dots \mathrm{,}{\nu }_m$ of $W_{*}$ are all real (since $W_{*}$ is symmetric). Also, $W_{*}$ and $\tilde{W}$ have identical eigenvalues. Being a stochastic matrix, $\tilde{W}$ has at least one eigenvalue which is one and the remaining eigenvalues are bounded above by 1, that is $\left|{\nu }_i\right|\le 1$ and ${{\displaystyle \max }}_i{\nu }_i=1.$ As $-1<\rho <1$ and $1-\rho {\nu }_i>0,$ $\left|{\Omega }_3\right|\mathrm{<mo>=</mo>}{\displaystyle {\prod }_{i\mathrm{<mo>=</mo>}1}^m}{\left(1-\rho {\nu }_i\right)}^2>0.$ Thus, the eigenvalues of ${\Omega }_3$ are positive and bounded above by $2m.$ Let $l_{\left(1\right)}=\min l_i$ and $l_{\left(m\right)}=\max l_i,$ where $L=\text{diag}{\left\{l_i\right\}}_{i\mathrm{<mo>=</mo>}1}^m.$ Then $l_{\left(1\right)}>0$ and $l_{\left(m\right)}$ is bounded above. By writing

$${\Sigma }_3\mathrm{<mo>=</mo>}\delta I_m+{\sigma }_v^2{\Omega }_3^{-1}=L^{-1/2}\left\{\delta L+{\sigma }_v^2{\left(I_m-\rho W_{*}\right)}^{-1}L{\left(I-\rho W_{*}\right)}^{-1}\right\}{L}^{-1/2}\mathrm{,}$$

we have

$$\begin{array}{ll}\left|{\Sigma }_3\right|={\left|L\right|}^{-1}\left|\delta L+{\sigma }_v^2{\left(I_m-\rho W_{*}\right)}^{-1}L{\left(I_m-\rho W_{*}\right)}^{-1}\right|\hfill & \ge {\left|L\right|}^{-1}{l}_{\left(1\right)}^m\left|\delta I_m+{\sigma }_v^2{\left(I_m-\rho W_{*}\right)}^{-2}\right|\mathrm{<mi>\; </mi>}\hfill \\ \hfill & ={\left|L\right|}^{-1}{l}_{\left(1\right)}^m\mathrm{<mi>\; </mi>}{\displaystyle \prod _{i\mathrm{<mo>=</mo>}1}^m}\left\{\delta +{\sigma }_v^2{\left(1-\rho {\nu }_i\right)}^{-2}\right\}\mathrm{,}(\text{A}.10)\hfill \end{array}$$

Letting $P_{W_{*}}$ be the matrix of eigenvectors of $W_{*}$ such that $P_{W_{*}}^T{W}_{*}{P}_{W_{*}}=\text{diag}{\left\{{\nu }_i\right\}}_{i\mathrm{<mo>=</mo>}1}^m=N_{*},$ we also have

$$\begin{array}{ll}{\left(d-G\varphi \right)}^T{\Sigma }_3^{-1}\left(d-G\varphi \right)\hfill & \mathrm{<mo>=</mo>}{\left(L^{1/2}d-L^{1/2}G\varphi \right)}^T{\left\{\delta L+{\sigma }_v^2{\left(I_m-\rho W_{*}\right)}^{-1}L{\left(I_m-\rho W_{*}\right)}^{-1}\right\}}^{-1}\left(L^{1/2}d-L^{1/2}G\varphi \right)\hfill \\ \hfill & \mathrm{<mo>=</mo>}{\left(r-S\varphi \right)}^T{\left\{\delta L+{\sigma }_v^2{\left(I_m-\rho W_{*}\right)}^{-1}L{\left(I_m-\rho W_{*}\right)}^{-1}\right\}}^{-1}\left(r-S\varphi \right)\hfill \\ \hfill & \ge {\left(l_{\left(m\right)}^{-1/2}r-l_{\left(m\right)}^{-1/2}S\varphi \right)}^T{\left\{\delta I_m+{\sigma }_v^2{\left(I_m-\rho W_{*}\right)}^{-2}\right\}}^{-1}\left(l_{\left(m\right)}^{-1/2}r-l_{\left(m\right)}^{-1/2}S\varphi \right)\hfill \\ \hfill & \mathrm{<mo>=</mo>}{\left(\tilde{r}-\tilde{S}\varphi \right)}^T{\left\{\delta I_m+{\sigma }_v^2{\left(I_m-\rho N_{*}\right)}^{-2}\right\}}^{-1}\left(\tilde{r}-\tilde{S}\varphi \right)\hfill \\ \hfill & \ge {\displaystyle \sum _{k\mathrm{<mo>=</mo>1}}^q}\frac{{\left({\tilde{r}}_{i_k}-{\tilde{s}}_{i_k}^T\mathbf{\varphi }\right)}^2}{\delta +{\sigma }_v^2{\left(1-\rho {\nu }_{i_k}\right)}^{-2}}\mathrm{,}(\text{A}.11)\hfill \end{array}$$

where $r=L^{1/2}d,$ $S=L^{1/2}G,$ $\tilde{r}=l_{\left(m\right)}^{-1/2}{P}_{W_{*}}r,$ $\tilde{S}=l_{\left(m\right)}^{-1/2}{P}_{W_{*}}S,$ and $\left\{i_1\mathrm{,}\dots \mathrm{,}{i}_q\right\}$ is a subset of $\left\{\mathrm{1,}\dots \mathrm{,}m\right\}$ so that the $q\times q$ matrix ${\left[{\tilde{s}}_{i_1}\mathrm{,}\dots \mathrm{,}{\tilde{s}}_{i_q}\right]}^T={\tilde{S}}_1,$ a submatrix of $\tilde{S},$ is non-singular. Note that ${\tilde{S}}_1$ is determined by $W.$ Using (A.11) we get

$${\displaystyle \int }\exp \left\{-\frac{1}{2}{\left(d-G\varphi \right)}^T{\Sigma }_3^{-1}\left(d-G\varphi \right)\right\}d\varphi \le {\left(2\pi \right)}^{q/2}{\left|{\tilde{S}}_1^T{\tilde{S}}_1\right|}^{-1/2}{\displaystyle \prod _{k\mathrm{<mo>=</mo>}1}^q}{\left\{\delta +{\sigma }_v^2{\left(1-\rho {\nu }_{i_k}\right)}^{-2}\right\}}^{1/2}\mathrm{.}(\text{A}.12)$$

Based on (A.10) and (A.12), we get that

$$\begin{array}{l}{\displaystyle \int }{\left|{\Sigma }_3\right|}^{-1/2}\exp \left\{-\frac{1}{2}{\left(d-G\varphi \right)}^T{\Sigma }_3^{-1}\left(d-G\varphi \right)\right\}d\varphi \le K{\displaystyle \prod _{i\notin \left\{i_1\mathrm{,}\dots \mathrm{,}{i}_q\right\}}}{\left\{\delta +{\sigma }_v^2{\left(1-\rho {\nu }_i\right)}^{-2}\right\}}^{-1/2}\\ \le K\left\{I\left({\sigma }_v^2<N\right)+{\left({\sigma }_v^2\right)}^{-\left(m-q\right)/2}I\left({\sigma }_v^2>N\right){\displaystyle \prod _{i\notin \left\{i_1\mathrm{,}\dots \mathrm{,}{i}_q\right\}}}\left(1-\rho {\nu }_i\right)\right\}\\ \le K\left\{I\left({\sigma }_v^2<N\right)+{\left({\sigma }_v^2\right)}^{-\left(m-q\right)/2}I\left({\sigma }_v^2>N\right)\right\}\mathrm{,}(\text{A}.13)\end{array}$$

where we use the fact that $-1<\rho <1$ and $-1\le {\nu }_i\le 1$ to claim $0<1-\rho {\nu }_i<2.$ From (A.5) and (A.13), it follows by proceeding along the lines we did for the CAR model that the desired integral ${\displaystyle \int }f\left(y_{\left(2\right)}\mathrm{,}\theta \mathrm{,}\beta \mathrm{,}{\sigma }_v^2\mathrm{,}\rho \right)d\theta d\beta d{\sigma }_v^{2}d\rho$ is finite under the conditions of the theorem.

A.3   Details for the CAR Model

We now consider $k=4$ for the IAR model where

$${\Omega }_4=L-\rho W=L^{1/2}\left(I_m-\rho W_{*}\right)L^{1/2}\mathrm{.}$$

Let ${\Sigma }_4=\delta I_m+{\sigma }_v^2{\Omega }_4^{-1}=L^{-1/2}\left\{\delta L+{\sigma }_v^2{\left(I_m-\rho W_{*}\right)}^{-1}\right\}{L}^{-1/2}.$ Then

$$\left|{\Sigma }_4\right|\ge {\left|L\right|}^{-1}{k}_{*}^m{\displaystyle \prod _{i=1}^m}\left\{\delta +{\sigma }_v^2{\mathrm{<mo\; stretchy="false">(</mo>1}-\rho {\nu }_i\mathrm{<mo\; stretchy="false">)</mo>}}^{-1}\right\}\mathrm{,}(\text{A}.14)$$

where $k_{*}^m=\min \left\{l_{\left(1\right)}\mathrm{,}1\right\}.$ Proceeding along the same line as in (A.11), we get that

$${\left(d-G\varphi \right)}^T{\Sigma }_4^{-1}\left(d-G\varphi \right)\ge {\displaystyle \sum _{k=1}^q}\frac{{\left({\tilde{r}}_{i_k}-{\tilde{s}}_{i_k}^T\varphi \right)}^2}{\delta +{\sigma }_v^2{\left(1-\rho {\nu }_{i_k}\right)}^{-1}}\mathrm{.}(\text{A}.15)$$

Again, as we had for the two previous cases, we can use (A.14) and (A.15) to establish that the desired integral is finite under the conditions stated in the theorem.

A.4   Details for the LCAR Model

Finally, we consider $k=5,$ where for the LCAR case we have

$${\Omega }_5=\rho R+\left(1-\rho \right)I_m\mathrm{.}$$

Suppose $r_1\mathrm{,}\dots \mathrm{,}{r}_m$ are the eigenvalues of $R$ and $P_R$ is an orthogonal matrix such that $P_R^{T}R{P}_R=\text{diag}{\left\{r_i\right\}}_{i=1}^m.$ Since $R$ is a non-negative definite matrix, $r_i\ge 0,$ $i=\mathrm{1,}\dots \mathrm{,}m,$ and ${\displaystyle {\sum }_{i\mathrm{<mo>=</mo>1}}^m}{r}_i=\text{tr}R={\displaystyle {\sum }_{i\mathrm{<mo>=</mo>}1}^m}{l}_i,$ implying that $r_1\mathrm{,}\dots \mathrm{,}{r}_m$ are all bounded between 0 and $l={\displaystyle {\sum }_{i\mathrm{<mo>=</mo>}1}^m}{l}_i.$ Then we can write

$${\Omega }_5=P_R\left\{\text{diag}{\left\{\rho r_i+1-\rho \right\}}_{i\mathrm{<mo>=</mo>}1}^m\right\}{P}_R^T\mathrm{,}$$

and claim that for $0<\rho <1,$ the eigenvalues of ${\Omega }_5$ are all positive and bounded above by ${\displaystyle {\sum }_{i=1}^m}{r}_i+1=l+1.$ Then, with $\tilde{r}=P_R^{T}d,$ and $\tilde{S}=P_R^{T}G,$ we can establish an inequality similar to (A.11). Note that the nonsingular matrix ${\tilde{S}}_1$ is a submatrix of $\tilde{S}$ and is free from $\rho .$ Boundedness of the eigenvalues of ${\Omega }_5$ will lead to an inequality similar to (A.14). Finally, we get the desired integral is finite under the conditions of the theorem.

References

Banerjee, S., Carlin, B.P., and Gelfand, A.E. (2003). Hierarchical Modeling and Analysis for Spatial Data, Chapman and Hall/CRC.

Berger, J.O. (1985). Statistical Decision Theory and Bayesian Analysis, Springer Science & Business Media.

Besag, J., and Kooperberg, C. (1995). On conditional and intrinsic autoregressions. Biometrika, 82, 733-746.

Brown, D.A., Datta, G.S. and Lazar, N.A. (2017). A Bayesian generalized CAR model for correlated signal detection. Statistica Sinica, 27, 1125-1153.

Datta, G.S., and Lahiri, P. (2000). A unified measure of uncertainty of estimated best linear unbiased predictors in small area estimation problems. Statistica Sinica, 10, 613-627.

Datta, G.S., and Smith, D.D. (2003). On propriety of posterior distributions of variance components in small area estimation. Journal of Statistical Planning and Inference, 112, 175-183.

Datta, G.S., Hall, P. and Mandal, A. (2011). Model selection by testing for the presence of small-area effects, and application to area-level data. Journal of the American Statistical Association, 106, 362-374.

Datta, G.S., Rao, J.N.K. and Smith, D.D. (2005). On measuring the variability of small area estimators under a basic area level model. Biometrika, 92, 183-196.

Datta, G.S., Lahiri, P., Maiti, T. and Lu, K.L. (1999). Hierarchical Bayes estimation of unemployment rates for the states of the U.S. Journal of the American Statistical Association, 94, 1074-1082.

Fay, R.E., and 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, 74, 269-277.

Gelman, A., and Rubin, D.B. (1992). Inference from iterative simulation using multiple sequences. Statistical Science, 7, 457-472.

Gelman, A., Carlin, J.B., Stern, H.S., Dunson, D.B., Vehtari, A. and Rubin, D.B. (2013). Bayesian Data Analysis, CRC Press, 3^rd ed.

Ghosh, M. (1992). Hierarchical and empirical Bayes multivariate estimation. In Current Issues in Statistical Inference: Essays in Honor of D. Basu, (Eds., M. Ghosh and P.K. Pathak), Institute of Mathematical Statistics, 151-177.

Hodges, J.S. (2019). Richly Parameterized Linear Models: Additive, Time Series, and Spatial Models Using Random Effects, Chapman and Hall/CRC.

Leroux, B.G., Lei, X. and Breslow, N. (2000). Estimation of disease rates in small areas: A new mixed model for spatial dependence. In Statistical Models in Epidemiology, the Environment, and Clinical Trials, Springer, 179-191.

MacNab, Y.C. (2003). Hierarchical Bayesian spatial modelling of small-area rates of non-rare disease. Statistics in Medicine, 22, 1761-1773.

Opsomer, J.D., Claeskens, G., Ranalli, M.G., Kauermann, G. and Breidt, F. (2008). Non-parametric small area estimation using penalized spline regression. Journal of the Royal Statistical Society: Series B (Statistical Methodology), 70, 265-286.

Porter, A.T., Wikle, C.K. and Holan, S.H. (2015). Small area estimation via multivariate Fay-Herriot models with latent spatial dependence. Australian & New Zealand Journal of Statistics, 57, 15-29.

Porter, A.T., Holan, S.H., Wikle, C.K. and Cressie, N. (2014). Spatial Fay-Herriot models for small area estimation with functional covariates. Spatial Statistics, 10, 27-42.

Prasad, N.G.N., and Rao, J.N.K. (1990). The estimation of the mean squared error of small-area estimators. Journal of the American Statistical Association, 85, 163-171.

Rao, J.N.K., and Molina, I. (2015). Small Area Estimation, New York: John Wiley & Sons, Inc.

Rao, J.N.K., Sinha, S.K. and Dumitrescu, L. (2014). Robust small area estimation under semi-parametric mixed models. Canadian Journal of Statistics, 42, 126-141.

Speckman, P.L., and Sun, D. (2003). Fully Bayesian spline smoothing and intrinsic autoregressive priors. Biometrika, 90, 289-302.

Stan Development Team (2018). RStan: The R interface to Stan. R package version 2.17.3.

Sun, D., Tsutakawa, R.K. and Speckman, P.L. (1999). Posterior distribution of hierarchical models using CAR (1) distributions. Biometrika, 86, 341-350.

Torabi, M. (2012). Hierarchical Bayes estimation of spatial statistics for rates. Journal of Statistical Planning and Inference, 142, 358-365.

Trevisani, M., and Gelfand, A. (2013). Spatial misalignment models for small area estimation: A simulation study. In Advances in Theoretical and Applied Statistics, Springer, 269-279.

Watanabe, S., and Opper, M. (2010). Asymptotic equivalence of Bayes cross validation and widely applicable information criterion in singular learning theory. Journal of Machine Learning Research, 11.

Whittle, P. (1954). On stationary processes in the plane. Biometrika, 41, 434-449.

You, Y., and Zhou, Q.M. (2011). Hierarchical Bayes small area estimation under a spatial model with application to health survey data. Survey Methodology, 37, 1, 25-37. Paper available at https://www150.statcan.gc.ca/n1/en/pub/12-001-x/2011001/article/11445-eng.pdf.

• Previous
• Article main page