Generalized bootstrapping as an alternative to classical weighting approaches

Generalized bootstrapping as an alternative to classical weighting approaches Felix Bittmann Leibniz Institute for Educational Trajectories Bamberg, Germany [email protected] 0000-0003-0802-5854 November 2025 (v1.0) Abstract. Some Stata commands do not accept (sampling) weights, which limits their applicability, particularly when working with survey data. As an alternative that can be used with any command, we propose generalized bootstrapping. Unlike the conventional bootstrap, in which each observation has an equal probability of being selected into a resample, generalized bootstrapping incorporates sampling weights, resulting in unequal probability sampling (UPS). We introduce the new command gbs, which implements generalized bootstrapping in Stata. Several examples illustrate its use, and a simulation study evaluates how well the command performs in comparison to standard sampling approaches. Keywords: bootstrapping, weighting, resampling, survey data, simulation, Stata 1 Introduction Sampling weights are of central importance in many statistical analyses. When working with survey data, sampling weights are often provided to obtain estimates that are representative of a target population. In addition, synthetic sampling weights can be used to address causal questions and create balanced samples (e.g., through entropy balancing). Statistics of interest are then computed by specifying the desired command with the [pweight=varname] option. However, some Stata commands do not permit the use of weights, such as spearman,ranksum, or alpha, which can limit their applicability. One solution to this problem is generalized bootstrapping. In contrast to regular bootstrapping (Efron and Tibshirani 1986; Bittmann 2021), where each observation has an equal probability of being included in a bootstrap resample, generalized bootstrapping incorporates sampling weights (Rao and Wu 1988). This results in unequal probability sampling (UPS), meaning that some observations have a higher probability of being selected. Instead of applying the sampling weight directly within the command of interest, one can generate a synthetic bootstrap resample and then apply the command to this newly generated sample. By repeating this process many times, a bootstrap distribution is created. Unlike the regular bootstrap, where the point estimate is derived from the original sample and resampling is performed solely for inference, in the generalized bootstrap © 2025 F. Bittmann DOI: 10.5281/zenodo.17581256 2Generalized bootstrap the point estimate is typically taken as either the mean or median of the bootstrap distribution. Inference can then be conducted by constructing confidence intervals from this distribution. Although the generalized bootstrap appears conceptually similar to the regular bootstrap, it presents a key limitation regarding inference. It has been proven that, in the regular bootstrap, both normal and percentile confidence intervals provide valid inference; however, this does not hold for the generalized bootstrap. With unequal probability sampling, some observations are more likely to be selected into a resample. While this can be advantageous—since the inclusion probability of extreme outliers may be reduced, depending on how the weights are constructed—it also means that each resample exhibits less variability than the original sample. Consequently, confidence intervals derived from a generalized bootstrap distribution are usually too narrow, resulting in overly liberal intervals. In practice, this makes it “too easy” to reject the null hypothesis, potentially leading to incorrect conclusions in favor of the alternative. The extent of this bias depends on the sampling weights used. One possible remedy is to forgo inference entirely and use the generalized bootstrap only to obtain point estimates, which remain valid. However, if inference is desired, several considerations apply. Because the resulting intervals tend to be too narrow, a non-significant effect can be interpreted with confidence, as even this liberal method fails to detect a significant result.1Conversely, a significant finding should be interpreted with caution. As an additional adjustment, one might consider applying a correction factor to the confidence intervals to account for the effect of the sampling weights. Multiplying the interval limits by a factor greater than one widens the intervals while leaving point estimates unchanged. However, no general rule exists for determining this factor, as the extent of interval shrinkage depends on the specific data and variables. A simple heuristic is to base the factor on the coefficient of variation (CV) of the sampling weights—the ratio of the standard deviation to the mean. A weight variable assigning similar values to all observations yields a smaller CV than one with widely varying weights. Yet, even this approach provides only limited guidance on how strongly to adjust the intervals. Users should therefore remain cautious and consider this limitation when interpreting results. 2 The gbs command 2.1 Installation and dependencies To install gbs, type: net install gbs , replace from(https://raw.githubusercontent.com/fbittmann/gbs/main/). Make sure that the dependencies, gsample2,moremata3, and parallel4(Vega Yon and Quistorff 2019), are also installed. 1. Of course, researchers should always think about the statistical power of their analysis as well. 2. https://github.com/benjann/gsample 3. https://github.com/benjann/moremata 4. https://github.com/gvegayon/parallel F. Bittmann 3 2.2 Syntax gbs if  in , expression(string) weight(varname)reps(#) seed(#) level(#) graph format(%fmt) parallel(#) cluster(varname) idcluster(varname) strata(varname) saving(filename,...) rround noadjust noisily : command 2.3 Options •expression(string)specifies the list of coefficients or statistics to bootstrap. To determine how other Stata commands return these values, type return list or ereturn list after executing the command of interest. Estimation results are often stored in a matrix called r(table). Specifying expression(string) is mandatory. •weight(varname)specifies the variable containing the sampling weights for each observation. Specifying weight(varname) is mandatory. •reps(#)specifies the number of bootstrap resamples to be drawn. The higher this value, the more precise the results will be. However, if a command requires substantial computation time, setting a very high number may be prohibitive. Using multiple threads (see parallel below) can reduce computation time. For testing purposes, as few as 100 replications may suffice, while typically at least 500 or 1,000 resamples are recommended. Inspecting the distribution of the bootstrap resamples (see the graph option below) can help assess whether a sufficient number of samples has been generated. For further guidance, see the simulation study below. •seed(#)specifies the random-number seed. Because bootstrapping depends on random resampling, estimation results will vary even with the same data. To obtain reproducible results, set a seed. When using more than one thread (see parallel below), the number of seeds specified must match the number of threads. •level(#)specifies the confidence level for the intervals. The default is level(95). •graph displays, for each statistic of interest, a histogram of the final bootstrap distribution. Ideally, these distributions should be smooth and approximately normal. If this is not the case, consider increasing the number of resamples. If the problem persists, bootstrapping may not be an appropriate method for the data at hand. •parallel(#)specifies the number of threads to be used for computation. Increasing this number can speed up processing, but it should not exceed the number of available system threads. Setting a value higher than your system supports may cause instability or crashes. For further details, refer to parallel. 4Generalized bootstrap •cluster(varname)specifies the variables identifying sampling clusters (i.e., primary sampling units). If cluster() is specified, the sample drawn consists of clusters. The cluster variables may be numeric or string. The weights must be constant within each cluster. •idcluster(varname)creates a new variable containing a unique identifier for each resampled cluster. This option requires that cluster() also be specified. •strata(varname)specifies the variables identifying strata. If strata() is specified, samples are selected within each stratum. The strata variables may be numeric or string. •saving(filename)creates a Stata data file (.dta) containing, for each statistic in expression(), a variable with the bootstrap replicates. •rround applies random rounding to non-integer sample sizes across strata. For details, see the documentation of gsample. •noisily displays additional output during execution. •noadjust specifies that the reported confidence intervals should not be adjusted, even though they are likely to be too narrow. In generalized bootstrapping with unequal sampling probabilities, observations with larger weights have a higher probability of being included in a resample. Consequently, the resamples tend to exhibit less variation than the original sample. While this does not affect point estimates substantially, it underestimates variances and produces overly narrow (liberal) confidence intervals, making it “easy” to reject the null hypothesis. Because this effect can be severe, depending on the weight variable, a crude adjustment step is implemented by default. This adjustment widens the confidence intervals based on the coefficient of variation of the weight variable. The noadjust option suppresses this correction. 3 Examples In our first example, we use gbs with spearman, a command for computing Spearman’s rank correlation coefficient ρ, which does not support sampling weights. We aim to examine whether the number of individuals living in a household is associated with self-reported health. First, we open the dataset and remove missing values from one variable. Next, we reverse this variable so that higher numerical values represent better health. Because spearman does not set the estimation sample, we manually remove missing values from all variables included in the analysis. Finally, we call gbs with 500 resamples. Setting the seed ensures reproducible results. . webuse nhanes2, clear . replace hlthstat = .a if hlthstat == 8 (14 real changes made, 14 to missing) . generate goodhealth = 6 - hlthstat F. Bittmann 5 (16 missing values generated) . drop if missing(goodhealth, finalwgt, houssiz) (16 observations deleted) . gbs, expression(rho=r(rho)) weight(finalwgt) reps(500) seed(123): /// > spearman goodhealth houssiz warning: the specified command does not set e(sample), so no observations will be excluded from the resampling because of missing values or other reasons. To exclude observations, press Break, save the data, drop any observations that are to be excluded, and rerun gbs. stat mean median lower upper type rho .0564884 .0563265 .0341004 .0825615 perc. rho . . .0342025 .081346 norm. Note: Normal CIs are based on the mean. CI level: 95% based on 500 resamples At first, gbs issues a warning message because spearman does not set e(sample). Since we have manually removed missing values from all relevant variables, this message can safely be ignored. Estimation commands such as regress automatically define the estimation sample, so less manual preparation is required in those cases. gbs then produces a concise output table. In this example, only a single statistic is computed. The table displays both the mean and the median across all 500 bootstrap resamples. As these values are typically very similar, the choice of which to report is usually inconsequential. However, if they differ substantially, this may indicate unreliable bootstrapping results. In such cases, consider increasing the number of bootstrap replications and using the graph option to inspect the distribution. The estimated coefficients are positive, indicating that self-reported health increases with the number of individuals living in a household. Alongside the point estimates, gbs also produces confidence intervals—both normal-based and percentile-based—which are generally quite similar. In our second example, we compute the reliability measure Cronbach’s alpha using sampling weights. For demonstration purposes, we open the auto dataset and generate an artificial sampling weight. We then call gbs, this time specifying the graph option to visualize the bootstrap distribution (see Figure 1). . sysuse auto, clear (1978 automobile data) . set seed 123 . generate w = runiform() . gbs, expression(alpha=r(alpha)) weight(w) reps(1500) seed(123) graph: /// > alpha headroom length turn trunk warning: the specified command does not set e(sample), so no observations will be excluded from the resampling because of missing values or other reasons. To exclude observations, press Break, save the data, drop any observations that are to be excluded, and rerun gbs. 6Generalized bootstrap stat mean median lower upper type alpha .549632 .5501717 .4780048 .6232775 perc. alpha . . .4795421 .6244531 norm. Note: Normal CIs are based on the mean. CI level: 95% based on 1500 resamples 0 5 10 15 20 .45 .5 .55 .6 alpha Figure 1: Bootstrap distribution As the distribution appears relatively smooth and approximately normal—even though the dataset contains only 74 observations—we conclude that there are no major issues with bootstrapping Cronbach’s alpha in this case. In the third and final example, we demonstrate the use of gbs with clustered data. We employ panel data and estimate a random-effects regression model. We begin by opening the dataset, removing a few survey waves to speed up computation, and generating a sampling weight. When working with clustered data, it is required that weights remain constant within each cluster (in this case, within each individual). We then create a new ID variable identical to the original one. Finally, we call gbs. To improve efficiency, we use two threads (which also requires specifying two seeds). The statistic of interest is a regression coefficient, extracted from r(table). Additionally, we retrieve the coefficient for the within-R2, stored in e(r2_w). . webuse nlswork, clear (National Longitudinal Survey of Young Women, 14-24 years old in 1968) . keep if year < 77 (16,575 observations deleted) . set seed 123 . generate temp = runiform() . bysort idcode: egen w = max(temp) . generate newid = idcode . xtset newid year F. Bittmann 7 Panel variable: newid (unbalanced) Time variable: year, 68 to 75, but with gaps Delta: 1 unit . gbs, expression(union=r(table)[1,1] r2_w=e(r2_w)) weight(w) reps(500) seed(15 55) /// > parallel(2) cluster(idcode) idcluster(newid): /// > xtreg ln_wage union south hours tenure, re stat mean median lower upper type union .0887347 .0870876 .059476 .1244409 perc. union . . .0562352 .1225394 norm. r2_w .1072725 .1057711 .070034 .1549264 perc. r2_w . . .065465 .1507651 norm. Note: Normal CIs are based on the mean. CI level: 95% based on 500 resamples 4 Simulation study To evaluate the statistical validity of generalized bootstrapping, we conducted a simulation study. Because gbs can also be used with commands that support sampling weights, the results of both methods can be directly compared. Treating the standard approach as the gold standard, these results serve as a benchmark for assessing the performance of gbs. For this purpose, a series of simulations were carried out using linear (OLS) and binary logistic regression models. The NHANES2 dataset was employed, and various aspects were systematically varied, including sample size, dependent and independent variables, and synthetically created sampling weights to ensure a high degree of randomness and variation. For each sample and setting, both the standard approach (with robust standard errors) and the generalized bootstrap approach using gbs were applied. Point estimates and confidence intervals were stored and compared. A total of 2,500 simulations per model (regress/logit) were conducted.5For gbs, 750 resamples were specified to achieve high precision, while all other parameters were kept at their default values. After conducting the simulations, we computed standardized differences as follows:     coefstandard −coefgbs coefstandard     . Multiplying this statistic by 100 yields the deviation between methods in percent. This transformation is necessary because the variables are measured on different scales, and a common metric facilitates comparison. For the confidence intervals, we first computed their widths for each method and then formed the ratio widthstandard widthgbs . Values greater than 1 indicate that the standard intervals are wider; values below 1 indicate that the bootstrap intervals are wider. To evaluate the coefficients, median 5. Final NLogit = 2,415, as a few simulations did not converge. 8Generalized bootstrap (least absolute deviations) regressions were estimated, and the widths of the confidence bands were visualized using violin plots. Results for the coefficients are shown in Figure 2, created with coefplot (Jann 2014). Figure 2: Differences between methods’ coefficients by model. 95% confidence intervals included. As the results show, the differences between the approaches are quite small. Although statistically significant at the 5% level, the maximum difference is only about 1.9%. The deviations are slightly smaller for the mean of the bootstrap distributions and somewhat larger for logistic than for OLS regression coefficients. Because errors are likely to decrease further when using more than 750 resamples (a feasible choice for applied work but often impractical in large-scale simulations), the results appear satisfactory. Figure 3 shows the ratio of confidence interval widths. As normal-based and percentile-based intervals are virtually identical, only the normal-based ones are displayed. For both models, the median ratio exceeds 1, indicating that the gbs confidence intervals remain too narrow even after adjustment. The interquartile range spans from 0.81 to 1.61 for regress and from 1.08 to 1.81 for logit. Because the intervals are still liberal on average, they should be interpreted with caution, as the null hypothesis may be rejected too frequently. Next, we examine how the quality of the results depends on sample size. To this end, we divided the sample sizes into quintiles (the median sample size is 2,992, ranging from 432 to 5,605; by design, the sample sizes are approximately uniformly distributed across this range). We then estimated median regression models for both the difference from the standard approach and the ratio of confidence interval widths. Results are shown in Figure 4. The differences between the standard and generalized bootstrap approaches diminish as sample size increases. However, since the largest quintile contains on average more than 5,000 observations, the differences—particularly in the CI width ratios— remain noticeable. F. Bittmann 9 Figure 3: Confidence interval width ratio by model. Values greater than 1 indicate that regular CIs are wider. The median ratios are 1.20 and 1.38. Finally, we assess how many bootstrap resamples are required to obtain stable and precise point estimates. Especially when the command of interest is computationally intensive, a very large number of resamples may be infeasible. While more resamples generally improve precision, computation time can become a limiting factor. As an empirical test, we repeated the linear regression estimations 300 times. For each run, 15,000 resamples were generated, and the estimation results were stored (Hesterberg 2015). The regression coefficient obtained with the standard sampling approach served as the benchmark. We then examined how quickly the bootstrapping results converged toward this benchmark. As in the previous simulations, we used a standardized difference metric expressed in percent. Factors such as sample size and the dependent and independent variables were randomly varied. Unlike in Figure 2, here we report arithmetic means with 95% confidence intervals rather than medians.6Results are presented in Figure 5. When the number of resamples is very low (50), the average deviation from the benchmark is approximately 10%. Once the number of resamples increases to around 200, the deviation drops to about 5%. Overall, the results clearly show that more resamples reduce error, though the rate of improvement slows (both axes in Figure 5 are logarithmic). Beyond roughly 5,000 resamples, further gains in precision become minimal. 5 Discussion Generalized bootstrapping is a flexible approach that is particularly useful when a Stata command does not accept (sampling) weights. Instead of applying weights directly 6. When medians are used, the results are very similar to those shown in Figure 2 (at 750 resamples).