  Survey data are often collected under multistage sampling designs where units
are binned to clusters that are sampled in a first stage. The unit-indexed
population variables of interest are typically dependent within cluster. We
propose a Fully Bayesian method that constructs an exact likelihood for the
observed sample to incorporate unit-level marginal sampling weights for
performing unbiased inference for population parameters while simultaneously
accounting for the dependence induced by sampling clusters of units to produce
correct uncertainty quantification. Our approach parameterizes cluster-indexed
random effects in both a marginal model for the response and a conditional
model for published, unit-level sampling weights. We compare our method to
plug-in Bayesian and frequentist alternatives in a simulation study and
demonstrate that our method most closely achieves correct uncertainty
quantification for model parameters, including the generating variances for
cluster-indexed random effects. We demonstrate our method in two applications
with NHANES data.
