6. Discussion

Benmei Liu, Partha Lahiri and Graham Kalton

Previous

In this paper, we report the results of a simulation study from a real finite population to evaluate the credible intervals obtained from four different hierarchical models in terms of their interval lengths and their design-based coverage properties. To the best of our knowledge, such a design-based evaluation of small area credible (or confidence) intervals has not previously been performed in the evaluation of small area estimates. 

In the simulation study, we have compared the design-based coverage properties of credible intervals resulting from different hierarchical Bayes models for estimating small area proportions from a stratified simple random sample design. Overall, none of the models emerges as a clear winner and so we are not in a position to recommend any of the models studied. 

The hierarchical Bayes version of the well-known Fay-Herriot model appears to produce overly conservative credible intervals. The non-normality of both the sampling and the linking models is a possible source of this problem. The credible intervals for the beta-logistic hierarchical model achieve almost the nominal coverage for the finite population proportions and the bias property for this model is the best among the four models being compared when the sample sizes are small. However, since one of the full conditionals for the beta-logistic model involves the survey-weighted proportions, there is a problem with the MCMC whenever the survey-weighted proportion is zero. The credible intervals for this model are also wider than those for the other two models with a logistic linking model. It may be possible to reduce the width of the credible interval for the beta-logistic model by modifying the model in some way, such as by employing a suitable two-part mixture random effect model that will avoid the problem of survey-weighted proportions of zero. Further investigation is needed. Also consideration could usefully be given to other possible models, possibly a discrete probability model for Level 1, to improve on interval estimation of small proportions for small areas.

The simulation study found that the coverage of the Bayesian credible intervals for the finite population proportions was far from the nominal 95 percent level for all four models, and a similar finding was also obtained for the design-based coverage of the widely-used Fay-Herriot model. In the light of these findings we carried out a number of further analyses in a search for an explanation. These analyses included: adding predictor variables to the models; using a uniform prior distribution for ${\sigma }_{\nu }^2$ (based on arguments made by Gelman 2006); the use of empirical best prediction approach for the M1 model; increasing the sample size in states with few births to a minimum of 50; and applying the methods to estimate the proportion of births in each state below the national median birthweight. Although there are some differences in the coverage properties for the state finite population proportions, none of these analyses produced coverage rates close to the nominal rates.  The only case where the nominal rates coincided with the actual coverage rates was for a simulated dataset constructed under model M1 for the state proportions below the national median birthweight; the average coverage rates were 5.1 and 5.2 percent for the EBP and HB approaches, respectively. 

This simulation study was restricted to a single stage sample design. In addition, for simplicity no auxiliary variables were included in the linking models in our main analyses, whereas in practice the inclusion of such variables is routine and almost essential. Further simulation studies are needed to cover different sample designs, different sample sizes, and to incorporate some auxiliary variables in the linking models. We hope that our study will encourage others to conduct similar design-based simulations to evaluate small area estimation methods. Based on our limited results, users of small area estimates need to be cautioned about the interpretation of the credible intervals associated with the estimates.

Acknowledgements

The authors wish to thank the associate editor and two referees for a number of constructive suggestions which led to substantial improvement of the original manuscript. The second author's research was supported by the National Science Foundation SES-085100.

Appendix

Appendix A

A1. Full conditional distributions for the parameters of each model

Let $\overrightarrow{p}={(p_{1w},\mathrm{...},p_{mw})}^t$ and $r_i=\frac{{\psi }_i}{{\psi }_i+{\sigma }_v^2}$.

The full conditional distributions for the Fay-Herriot model (M1) are given as follows:

i) ${\theta }_i|\mu ,{\sigma }_v^2,\overrightarrow{p}~N((1-r_i)p_{iw}+r_i\mu ,{\psi }_i(1-r_i))$;

ii) $\mu |{\theta }_i,{\sigma }_v^2,p~N\left(\frac{1}{m}{\displaystyle \sum _{i=1}^m{\theta }_i},\frac{{\sigma }_v^2}{m}\right)$;

iii) ${\sigma }_v^2|\mu ,{\theta }_i,\overrightarrow{p}~ING\left(a+\frac{1}{2}m,b+\frac{1}{2}{\displaystyle \sum _{i=1}^m({\theta }_i}-\mu {)}^2\right)$.

The full conditional distributions for the Normal-Logistic model (M2) are given as follows:

i) ${\theta }_i|\mu ,{\sigma }_v^2,\overrightarrow{p}\propto \frac{1}{{\theta }_i(1-{\theta }_i){\sigma }_v\sqrt{{\psi }_i}}\exp \left(-\frac{{(p_{iw}-{\theta }_i)}^2}{2{\psi }_i}-\frac{{(\text{logit}({\theta }_i)-\mu )}^2}{2{\sigma }_v^2}\right)$;

ii) $\mu |{\theta }_i,{\sigma }_v^2,p~N\left(\frac{1}{m}{\displaystyle \sum _{i=1}^m\text{logit}({\theta }_i}),\frac{{\sigma }_v^2}{m}\right)$;

iii) ${\sigma }_v^2|\mu ,{\theta }_i,\overrightarrow{p}~ING\left(a+\frac{1}{2}m,b+\frac{1}{2}{\displaystyle \sum _{i=1}^m(\text{logit}({\theta }_i})-\mu {)}^2\right)$.

The full conditional distributions for the Normal-Logistic model with unknown variance (M3) are the same as those of M2 except that replacing ${\psi }_i$ by ${\theta }_i(1-{\theta }_i)def{f}_{iw}/n_i$ for the distribution of ${\theta }_i$ given other parameters.

Let ${\delta }_{iw}=\frac{n_i}{def{f}_{iw}}-1.$ The full conditional distributions for the Beta-Logistic model (M4) are given as follows:

i) ${\theta }_i|\mu ,{\sigma }_v^2,\overrightarrow{p}\propto \frac{1}{{\theta }_i(1-{\theta }_i){\sigma }_v}\frac{p_{iw}^{{\theta }_i{\delta }_{iw}-1}{(1-p_{iw})}^{(1-{\theta }_i){\delta }_{iw}-1}}{\Gamma ({\theta }_i{\delta }_{iw})\Gamma ((1-{\theta }_i){\delta }_{iw})}\exp \left(-\frac{{(\text{logit}({\theta }_i)-\mu )}^2}{2{\sigma }_v^2}\right);$

ii) $\mu |{\theta }_i,{\sigma }_v^2,p~N\left(\frac{1}{m}{\displaystyle \sum _{i=1}^m\text{logit}({\theta }_i}),\frac{{\sigma }_v^2}{m}\right)$;

iii) ${\sigma }_v^2|\mu ,{\theta }_i,\overrightarrow{p}~ING\left(a+\frac{1}{2}m,b+\frac{1}{2}{\displaystyle \sum _{i=1}^m(\text{logit}({\theta }_i})-\mu {)}^2\right)$.

Appendix B

WinBUGS code for Model 1:

model {
      for ( i in 1:N)  {
         pobs[i] ~ dnorm(theta[i], D[i])
         D[i] <- 1/varhat[i]
         theta[i]<-u+v[i]
         v[i]~dnorm(0, tau)
                          }
       u~dflat()
       tau~dgamma(0.001, 0.001)
       sigma_v2<-1/tau
          }

WinBUGS code for Model 2:

model {
      for ( i in 1:N)  {
         pobs[i] ~ dnorm(theta[i], D[i])
         D[i] <- 1/varhat[i]
         logit(theta[i])<-u+v[i]
         v[i]~dnorm(0, tau)
                           }
       u~dflat()
       tau~dgamma(0.001, 0.001)
       sigma_v2<-1/tau
          }

WinBUGS code for Model 3:

model {
      for ( i in 1:N)  {
         pobs[i] ~ dnorm(theta[i], E[i])
         E[i] <- SAMPn[i]/(theta[i]*(1-theta[i])*DEFF_kish[i])
         logit(theta[i])<-u+v[i]
         v[i]~dnorm(0, tau)
         D[i]<-1/E[i]
                          }
       u~dflat()
       tau~dgamma(0.001, 0.001)
       sigma_v2<-1/tau
          }

WinBUGS code for Model 4:

model {
      for ( i in 1:N)  {
         pobs[i] ~ dbeta(a[i], b[i])
         a[i] <- theta[i]*(theta[i]*(1-theta[i])/D[i]-1)
         b[i] <- (1-theta[i])*(theta[i]*(1-theta[i])/D[i]-1)
         logit(theta[i])<-u+v[i]
         v[i]~dnorm(0, tau)
         D[i]<-theta[i]*(1-theta[i])*DEFF_kish[i]/SAMPn[i]
                          }
       u~dflat()
       tau~dgamma(0.001, 0.001)
       sigma_v2<-1/tau
          }

References

Arora, V. and Lahiri, P. (1997). On the superiority of the Bayesian method over the BLUP in small area estimation problems. Statistica Sinica, 7, 1053-1063.

Chattopadhyay, M., Lahiri, P., Larsen, M., and Reimnitz, J. (1999). Composite estimation of drug prevalence for sub-state areas. Survey Methodology, 25, 81-86.

Chen, M., Shao, Q., and Ibraham, J.G. (2000). Monte Carlo Methods in Bayesian Computation. New York: Springer-Verlag.

Citro, C., and Kalton, G. (Eds.) (2000). Small-Area Income and Poverty Estimates: Priorities for 2000 and Beyond. Washington, DC: National Academy Press.

Farrell, P.J., MacGibbon, B., and Tomberlin, T.J. (1997). Empirical Bayes estimators of small area proportions in multistage designs. Statistical Sinica, 7, 1065-1083.

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.  (2006). Prior distributions for variance parameters in hierarchical models. Bayesian Analysis, 1, 515-533.

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

Jiang, J., and Lahiri, P. (2006a). Mixed model prediction and small area estimation. Test, 15, 111-999.

Jiang, J., and Lahiri, P. (2006b). Estimation of finite population domain means: A model-assisted empirical best prediction approach. Journal of the American Statistical Association, 101, 301-311.

Kish, L. (1965). Survey sampling. New York: John Wiley.

Lunn, D.J., Thomas, A., Best, N., and Spiegelhalter, D. (2000). WinBUGS - A Bayesian modeling framework: Concepts, structure, and extensibility. Statistics and Computing, 10, 325-337.

Malec, D., Davis, W., and Cao, X. (1999). Small area estimates of overweight prevalence using sample selection adjustment. Statistics in Medicine, 18, 3189-3200.

Malec, D., Sedransk, J., Moriarity, C.L., and LeClere, F.B. (1997). Small area inference for binary variables in the National Health Interview Survey. Journal of the American Statistical Association, 92, 815-826.

Maples, J., and Bell, W.R. (2005). Evaluation of school district poverty estimates: Predictive models using IRS income tax data. Proceedings of the Survey Research Methods Section, American Statistical Association, 1322-1329.

Mohadjer, L., Rao, J.N.K., Liu, B., Krenzke, T., and Van De Kerckhove, W. (2012). Hierarchical Bayes small area estimates of adult literacy using unmatched sampling and linking models. Journal of the Indian Society of Agricultural Statistics, 66 (1), 55-63.

Rao, J.N.K. (2003). Small area estimation. New York: John Wiley and Sons.

Robert, C.P., and Casella, G. (1999). Monte Carlo Statistical Methods. New York: Springer-Verlag.

Singh, A.C., Folsom, R.E. JR., and Vaish, A.K. (2005). Small area modeling for survey data with smoothed error covariance structure via generalized design effects. Federal Committee on Statistical Methods Statistical Working Paper, No. 39 (http://www.fcsm.gov/05papers/Singh_etal_IIIC.pdf).

U.S. National Center for Health Statistics (2009). National Vital Statistics System. Birth Data. (http://www.cdc.gov/nchs/births.htm)

Wright, D., Sathe, N., and Spagnola, K. (2007). State Estimates of Substance Use from the 2004-2005 National Surveys on Drug Use and Health. (DHHS Publication No. SMA 07-4235, NSDUH Series H-31). Rockville, MD: Substance Abuse and Mental Health Services Administration, Office of Applied Studies.

You, Y. (2008). An integrated modeling approach to unemployment rate estimation for subprovincial areas of Canada. Survey Methodology, 34, 19-27.

You, Y., and Rao, J.N.K. (2002). Small area estimation using unmatched sampling and linking models. Canadian Journal of Statistics, 30, 3-15.

Previous