Bayesian predictive inference of small area proportions under selection bias
Section 4. Numerical studies

• Previous
• Next

In Section 4.1, we describe an example on severe activity limitation. We present the results of the IS, HoS and HeS models for the comparison. In Section 4.2, we describe a simulation study to assess the performance of the three models under two kinds of assumptions of the distribution of the sample selection probabilities. That is, data are generated from either the homogeneous nonignorable or the heterogeneous nonignorable selection model and all three models are fit.

4.1   Illustrative example

In our application, we use data from the 1995 National Health Interview Survey (NHIS95). These data were first used by Nandram, Bai and Choi (2011) to estimate change point in activity limitation. NBBS constructed synthetic data arising from a segment of NHIS95 for this study. For adults, 30-80 years old, NBBS analyzed data on severe activity limitation (SAL) for a single area (no pooling), where $y\mathrm{<mo>=</mo>}1$ if an adult has SAL and $y\mathrm{<mo>=</mo>}0$ otherwise. CNK used the data from NBBS to perform small area estimation on SAL. We use the data in a manner similar to NBBS but we fit a more general model, the key contribution of this paper.

NBBS formed twelve domains (small areas) by crossing education, sex and race. They have categorized education into three levels (pre-college: L, college: M and post-college: H). Sex (male: M and female: F) and Race (white: W and nonwhite: B) have two levels each. We will continue to call these domains LMW, LMB, LFW, LFB, MMW, MMB, MFW, MFB, HMW, HMB, HFW, HFB (e.g., LMW: white males with pre-college education, LMB: black males with pre-college education, etc.).

NHIS95 used a multistage sample design to draw samples from the population of the United States. Therefore, it is necessary to use an adult’s survey weight for accurate analysis. See NBBS and CNK for a discussion of the survey weights. Like NBBS and CNK, we considered the reciprocal of a survey weight as the “selection” probability of each adult. Selection probabilities are a major part of the survey weights. It would have been better if we had the selection probabilities. But this is an approximation we make in our data analysis. For convenience, we have presented the data again in Table 4.1.

We reduce the sample size from the original data set to increase the effect of small area model. For the HeS model, we order the selection probabilities from smallest to largest ${\pi }_{\left(\text{}1\text{}\right)}\mathrm{,}\dots ,{\pi }_{\left(n_i\right)}$ within each small area. Let the quantiles be $t_{i1}\mathrm{<mo>=</mo>}{\pi }_{\left(0.20n_i\right)}\mathrm{,}{t}_{i2}\mathrm{<mo>=</mo>}{\pi }_{\left(0.40n_i\right)}\mathrm{,}{t}_{i3}\mathrm{<mo>=</mo>}{\pi }_{\left(0.60n_i\right)}\mathrm{,}{t}_{i4}\mathrm{<mo>=</mo>}{\pi }_{\left(0.80n_i\right)}$ and let $t_{i0}\mathrm{<mo>=</mo>}{\pi }_{\left(\text{}1\text{}\right)}$ and $t_{i5}\mathrm{<mo>=</mo>}{\pi }_{\left(n_i\right)}$ for $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,l.$ We define ${\pi }_{iu}^{*}\mathrm{<mo>=</mo>}(t_{i(u-1)}+t_{iu})/2,$ $u\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,U_i,$ $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,l$ (i.e., the mid point of each quantile within each small area). Note that ${\theta }_{iuy}$ is the proportion of sampled units in the $u^{\text{th}}$ quantile conditional on $y\mathrm{<mo>=</mo>}\mathrm{0,}1$ for $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,l.$ If the ${\theta }_{iu0}$ are considerably different from ${\theta }_{iu1}$ in some areas, there is strong evidence that the sampled values are biased. To specify ${\theta }_{i0}^{\left(\text{}0\text{}\right)}$ and ${\theta }_{i1}^{\left(\text{}0\text{}\right)}$ in the prior distributions, we take ${\theta }_{iy}^{\left(\text{}0\text{}\right)}\mathrm{<mo>=</mo>}{\widehat{\theta }}_{iy},$ the maximum likelihood estimator of ${\theta }_{iy},$ $y\mathrm{<mo>=</mo>}\mathrm{0,}1,$ $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,l,$ which we call a “MLE” prior (see NBBS).

For our data, mixing is slow, so we draw 120,000 samples and burn in 60,000. Then we take every ${30}^{\text{th}}$ iterate to obtain a sample of 2,000 iterates for inference. This burn-in period is sufficiently long to get random samples, which is based on the trace plots and Geweke test. We have enough samples since the effective sample size (ESS) of parameters sampled from an MCMC should be less than or equal to 2,000. The correlation is nonsignificant for all parameters. Also, stationarity of our sampler is demonstrated by Geweke test. The results of convergence diagnostic for hyper-parameters are shown in Table 4.8 and Figure 4.2. By taking ${\kappa }_0\mathrm{<mo>=</mo>}\left(\mathrm{7.3,}\mathrm{3,}\mathrm{11,}\mathrm{4.8,}\mathrm{8.5,}\mathrm{9,}\mathrm{11.5,}\mathrm{10,}\mathrm{10,}\mathrm{9.7,}\mathrm{10.2,}11\right)$ and ${\kappa }_1\mathrm{<mo>=</mo>}\left(\mathrm{3.2,}\mathrm{1.7,}\mathrm{2.2,}\mathrm{1.5,}\mathrm{5,}\mathrm{2.2,}\mathrm{3.3,}\mathrm{2.2,}\mathrm{2.5,}\mathrm{1.5,}\mathrm{1.5,}1.3\right),$ we obtained acceptance rates between 25% and 45% of Metropolis-Hastings samplers.

The summaries of data are shown in Table 4.1. The averages of the selection probabilities in each area are mostly similar for $y\mathrm{<mo>=</mo>}0$ and $y\mathrm{<mo>=</mo>}1.$ Since the sample sizes within the domains are not quite large, the sampling fractions are very small. Thus, the simultaneous analysis using a small area model is appropriate. In column 4 of Table 4.1 we have also presented the proportion of individuals with SAL, and we can see that this value is relatively large for low level education. We compared the counts in the two sets of bins from the histograms of the selection probabilities for $y\mathrm{<mo>=</mo>}\mathrm{0,}1$ in each domain. In fact, this is a test of independence in a $2\times U_i$ categorical table, and we use a chi-squared test and a Fisher’s exact test for equality of ${\theta }_{iu0}\mathrm{<mo>=</mo>}{\theta }_{iu1},$ $u\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,U_i$ $(U_i\mathrm{<mo>=</mo>}5$ cells) in each domain. As is evident from the $p$ -values which are presented in the last column of Table 4.1, selection bias should matter mostly in Domain MMW and perhaps in Domains LMW, MMB, HMW and HFW. As a measure of selection bias, we also calculated the biserial correlation between the $y_{ij}$ and ${\pi }_{ij}$ for all areas combined and got a value of 0.033 with a $p$ -value of 0.03.

In Table 4.2, we provide the results of the unpooled (individual analysis) of the finite population proportions under the ignorable selection model and the nonignorable selection model for each domain separately (the nonignorable selection model is the NBBS model). We compare them using the posterior means (PM), the posterior standard deviations (PSD) of PMs, the coefficient of variation (CV) and 95% highest posterior density (HPD) intervals. This is a repetition of the analysis under the NBBS model.

In Table 4.3, we compare summaries of the pooled estimators of the finite population proportions under the IS model, the HoS model and the HeS model for the 12 domains. The effect of the selection bias is seen because the PMs under the IS model and the other models are different. The estimators under the HoS model and the HeS model are similar in some domains because the selection bias is more severe for some domains. The PSDs under the HoS model and the HeS model are bigger than under the IS model in most domains, but the 95% HPD intervals overlap. The effect of the simultaneous analysis is also seen because the PMs for pooled estimators are smoothed relative to PMs for separately analyzed finite populations.

The posterior summaries of ${\theta }_{i0}$ and ${\theta }_{i1}$ for $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,l$ are shown in Tables 4.4. The ${\theta }_{iu0}$ are different from the ${\theta }_{iu1}$ in some areas. It is strong evidence to indicate the presence of selection bias. We also present the posterior summaries of ${\theta }_0$ and ${\theta }_1$ under the HoS model in Table 4.5. Note that in Tables 4.4-4.5, these ${\theta }_{i0}$ and ${\theta }_{i1}$ are not selection probabilities. The $\theta$ under the HoS model appears different from those under the HeS model. We present the posterior summaries of $a_{i0}$ and $a_{i1}$ for $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,l$ in Table 4.7, and $a_0$ and $a_1$ of the HoS model in Table 4.6. These values are very small, but its ratio is large in each area. It is also strong evidence of selection bias, albeit these are small sample sizes within areas.

We present the posterior densities of the finite population proportions in Figure 4.1. The solid line represents the density of $P$ under the IS model, the dotted line represents that under the HoS model, and the dashed line represents that under the HeS model for each domain. We can see that the plots of IS model are shifted due to the effect of the selection bias.

$$ Our example shows that it is important to include a component for the selection mechanism into a model when a biased sample is available, otherwise there is likely to be misleading estimates of the small area proportions. Because the distribution of the sample selection probabilities is a little different across domains, the heterogeneity assumption is more reasonable for our numerical example.

4.2   Simulation study

We perform a simulation study to assess the accuracy of our model. Specifically, we consider two situations, which are homogeneous or heterogeneous distributions of the sample selection probabilities to show how different the IS, HoS and HeS models can be. We also show that when the correlation between the binary responses and the selection probabilities is strong, the IS model can perform badly.

To perform a simulation study, we generate $l\mathrm{<mo>=</mo>}12$ finite populations with $i^{\text{th}}$ finite population having $N_i\mathrm{<mo>=</mo>}$ 1,000 units and selection probability ${\pi }_{ij}\mathrm{,}$ $j\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,$ 1,000, $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,12.$ Then samples of size $n_i\mathrm{<mo>=</mo>}$ 50, $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,12$ are taken from each area. We have generated 100 data sets. The outline of the data collection is as follows.

We control the biserial correlation between the binary responses and the selection probabilities by changing the $a$ values in Step 2. For example, if we set $a\mathrm{<mo>=</mo>}$ 0.975, 0.95, 0.9, 0.8, 0.7 in the simulation scheme, then the biserial correlations are respectably $\rho \mathrm{<mo>=</mo>}$ 0.05, 0.1, 0.2, 0.4, 0.7. We generate data having several biserial correlations, and categorized them as three levels (Low: $\rho \mathrm{<mo><</mo>}$ 0.3, Medium: 0.3  $\le \rho \mathrm{<mo><</mo>}$ 0.6, High: $\rho \ge$ 0.6). These correspond to weak, medium and strong selection bias.

We can also control the heterogeneity or homogeneity assumptions by changing a support of the sample selection probabilities. For homogeneity assumption, we generate ${\pi }_{ij}$ of Step 4 from similar supports for all areas. On the other hand, we set up a support of ${\pi }_{ij}$ from a different interval by area to make heterogeneous distributions of the sample selection probabilities.

To compare the performance of three models, we compute several frequentist measures. First, we calculate the finite population proportion $P_i^{(h)},$ the posterior mean ${\text{PM}}_i^{(h)}$ and the posterior standard deviation ${\text{PSD}}_i^{(h)},$ $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,l,$ $h\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,100.$ Then we compute the absolute bias ${\text{AB}}_i^{(h)}\mathrm{<mo>=</mo>}\left|{\text{PM}}_i^{(h)}-P_i^{(h)}\right|$ and the root mean squared error ${\text{RMSE}}_i^{(h)}\mathrm{<mo>=</mo>}\sqrt{{\text{PSD}}_i^{(h)2}+{\text{AB}}_i^{(h)2}},$ $i\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,l,$ $h\mathrm{<mo>=</mo>}\mathrm{1,}\dots ,100.$ Using these frequentist quantities we obtain ${\text{AB}}_i\mathrm{<mo>=</mo>}\frac{1}{100}{\displaystyle {\sum }_{h\mathrm{<mo>=</mo>1}}^{100}{\text{AB}}_i^{(h)}}$ and ${\text{RMSE}}_i\mathrm{<mo>=</mo>}\frac{1}{100}{\displaystyle {\sum }_{h\mathrm{<mo>=</mo>1}}^{100}{\text{RMSE}}_i^{(h)}}.$ Also we compute the 95% highest posterior density (HPD) interval for each of the 100 simulated runs in each area. Then we look at the width $(W_i^{(h)})$ and the HPD incidence $(I_i^{(h)}).$ Let $I_i^{(h)}\mathrm{<mo>=</mo>}1$ if the 95% HPD interval contains the true value $P_i$ and let $I_i^{(h)}\mathrm{<mo>=</mo>}0$ otherwise. Then we calculate the coverage $C_i\mathrm{<mo>=</mo>}{\displaystyle {\sum }_{h\mathrm{<mo>=</mo>1}}^{100}{I}_i^{(h)}/100}$ and $W_i\mathrm{<mo>=</mo>}{\displaystyle {\sum }_{h\mathrm{<mo>=</mo>1}}^{100}{W}_i^{(h)}/100}.$

The biserial correlation (i.e., Pearson correlation coefficient) $\rho$ between $y$ and $\pi$ is calculated since the variable $y$ is binary (e.g., see Cox, 1974). The summaries based on AB, RMSE, coverage and width of the 95% HPD intervals for three correlation cases under the homogeneity assumption are shown in Table 4.9. Under this assumption, the performances of the HoS model is better than the HeS model in some domains, but the difference is very small. Both of them are better than the IS model in the sense of the closeness to the true $P_i.$ In particular, we can see that as the correlation between $y$ and $\pi$ increases, there is greater discrepancy between the IS model and the others.

Table 4.10 shows the summaries for three correlation cases under the heterogeneity assumption. From this table, we find that the HeS model is well behaved in the sense that it is closer to true $P_i$ than the other models when the effect of the selection bias is moderate to strong. It has smaller bias, smaller mean square error and better coverage.

As the biserial correlation increases there are increased disparities between the IS model and the HeS model. We present some plots of the posterior distributions of the finite population proportions for low, medium, high correlation cases under the homogeneity assumption (Figures 4.3-4.5). It is worth noting that there are large differences of the distributions of the finite population means as $\rho$ increases. There are no differences in the posterior densities of the HoS model and the HeS model under the homogeneity assumption. The posterior distributions under the HeS model departs from the right of the distribution under the IS model, with little changes in their spreads for all domains.

Similarly, we present some plots of the posterior distributions of the finite population proportions for low, medium, high correlation cases under the heterogeneity assumption (Figures 4.6-4.8). The posterior densities under the HeS model are different from the others for all domains.

Also, we have selected two values of $l\mathrm{<mo>=</mo>}\mathrm{12,}24,$ the number of areas, to see changes by increasing $l.$ We compute AB, RMSE, coverage and width of 95% HPD intervals for each domain, then we average these values over all domains. We present the results under the homogeneity and the heterogeneity assumptions in respectively Tables 4.11 and 4.12. The results are consistent as $l$  increases.

The simulation study shows that as the biserial correlation between the binary responses and the selection probabilities increases, there is greater discrepancy among the IS model, the HoS moel and the HeS model. That is, as the correlation is stronger, the effect of the selection bias becomes larger. The HeS model is better than the other models when the each area has the different distribution for the sample selection probabilities under moderate to strong selection bias.

• Previous
• Article main page
• Next