8
$\begingroup$

I'm using R to provide bootstrap (percentile and t methods) of estimated population totals, using data from a complex survey. It is a stratified survey of tourists expenditure that is weighted to population (the known number of tourists from customs information). I want a confidence interval for total expenditure of various combinations of tourist types (eg "Australian business travellers"). Those tourist characteristics are part of the poststratification weighting scheme, but not the original stratification. Stratification is by departing flight where the interview took place. The weights are very complex and definitely not independent of the variable of interest.

My procedure is to produce a stratified resample; take the cases in this sample's weights from the original sample and scale them up/down so they add again to the correct population totals; and calculate my statistic. Repeat r times. If I don't do the reweighting procedure in the middle my bootstrapped estimates are on average significantly higher than the estimates from the original sample and hence useless. My reweighting procedure does not exactly duplicate the original poststratification weighting (it is a simpler version) but is callibrated to produce the same marginal population totals for combinations of tourists' country of residence, purpose of visit, and age.

I have written code in R to do this but am interested in if there is an existing function I can compare my results with. I've looked at both the boot and survey packages and can't find anything. While boot allows bootstrapping from a stratified sample, I can't see a way to perform the reweighting to population of each replicate resample. My first question is - any pointers to a prebuilt function in R that does this reweighting to population as part of its bootstrap?

I also haven't found anything on exactly this problem in my cursory examination of the literature. However, it must be a common challenge for surveys that have been weighted to populations. My second question is - any pointers to discussion in the literature (not just on R) about the merits or otherwise of this reweighting procedure in the middle of a bootstrap?

Improve this question
$\endgroup$
4
  • $\begingroup$ Don't know about bootstrap, but the Australian Buruea of Statistics uses a jacknife method for variance estimation that sounds very similar to what your bootstrap is doing. The algorithm is called GREGWT (generalised regression weighting), and the jacknife replicates are also "benchmarked" same as the original weights. This may be useful to compare. Implementation is in SAS though... $\endgroup$
    – probabilityislogic
    Mar 14, 2012 at 20:20
  • 3
    $\begingroup$ Thanks. In the end Lumley's survey package did what I wanted once I got rid of the idea of doing a bootstrap in a loop myself and instead generated replicate weights and raked them to the original marginal total weights. It also does a variety of jacknife replicate weights. I also found the EVER library immensely useful - it has an implementation of Kott's delete-a-group jacknife which is a method trialled by Statistics New Zealand (in SAS) and in fact used in our own surveys (again in SAS). I might look into ABS' GREGWT. $\endgroup$
    – Peter Ellis
    Mar 14, 2012 at 22:00
  • $\begingroup$ Canty Davison (1999) have a good paper on this topic jstor.org/stable/2681000 $\endgroup$
    – Estimate the estimators
    Jan 24 at 7:16
  • $\begingroup$ The 'svrep' package implements a few different survey bootstrap methods, and includes this vignette to help choose which is most appropriate: cran.r-project.org/web/packages/svrep/vignettes/…. $\endgroup$
    – bschneidr
    Feb 15 at 0:31

2 Answers 2

Reset to default
7
$\begingroup$

In R i would tell you to see if the functions related with "bootweights" in the survey package suit you in any way. But since you have already gone over that package I don't think you will find many alternatives ... I also looked for a similiar thing a couple of weeks ago and ended up implementing my own code.

For the discussion of bootstrapping and survey weights in general you can find some references in this presentation which also contains references to an implementation of a bootstrapping procedure for complex survey designs in STATA.

Improve this answer
$\endgroup$
1
1
$\begingroup$

References on the reweighting method

When you reweight your data to match known population totals (using raking, post-stratification, or some other form of calibration), it has long been common practice to repeat the reweighting procedure for each replicate sample. This practice and its justification are described clearly in the following classical references on variance estimation for surveys:

Packages in R

The 'survey' and 'svrep' packages both provide a few different methods for bootstrapping with survey data. This vignette from the 'svrep' package provides guidance on how to choose an appropriate bootstrap method and number of bootstrap replicates:

https://cran.r-project.org/web/packages/svrep/vignettes/bootstrap-replicates.html.

When you use the 'survey' package's functions (rake(), postStratify(), or calibrate()) to reweight data with bootstrap replicate weights, the package will automatically repeat the reweighting procedure separately for each bootstrap replicate.

Below is example R code for how to implement this:

library(survey)
library(svrep)

# Load example survey data ----
data('lou_vax_survey', package = 'svrep')

# Create bootstrap weights ----
  vax_survey_design <- svydesign(data = lou_vax_survey,
                                 ids = ~ 1,
                                 prob = ~ SAMPLING_WEIGHT)
  
  boot_design <- as_bootstrap_design(
    design = vax_survey_design,
    replicates = 500
  )

# Define control totals (i.e. known population values) ----
  
  control_totals <- list(
    'RACE_ETHNICITY' = data.frame(
      'RACE_ETHNICITY' = c(
        "Black or African American alone, not Hispanic or Latino", 
        "Hispanic or Latino", "Other Race, not Hispanic or Latino",
        "White alone, not Hispanic or Latino"),
      'TOTAL' = c(119041, 27001, 27633, 423027)
    ),
    'SEX' = data.frame(
      'SEX' = c("Male", "Female"),
      'TOTAL' = c(283688, 313014)
    )
  )
  
# Reweight the data to match control totals, using raking ----
  raked_boot_design <- rake(
    design = boot_design,
    sample.margins = list(~ RACE_ETHNICITY, ~ SEX),
    population.margins = control_totals,
    control = list(maxit = 20, epsilon = 0.0001)
  )
  
# Check the resulting estimates ----
  
  raked_boot_design |> svytable(formula = ~ RACE_ETHNICITY)
#> RACE_ETHNICITY
#> Black or African American alone, not Hispanic or Latino 
#>                                                  119041 
#>                                      Hispanic or Latino 
#>                                                   27001 
#>                      Other Race, not Hispanic or Latino 
#>                                                   27633 
#>                     White alone, not Hispanic or Latino 
#>                                                  423027
  raked_boot_design |> svytable(formula = ~ SEX)
#> SEX
#> Female   Male 
#> 313014 283688
Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.