`vcovBS.Rd`

Object-oriented estimation of basic bootstrap covariances, using
simple (clustered) case-based resampling, plus more refined methods
for `lm`

and `glm`

models.

vcovBS(x, ...) # S3 method for default vcovBS(x, cluster = NULL, R = 250, start = FALSE, ..., fix = FALSE, use = "pairwise.complete.obs", applyfun = NULL, cores = NULL) # S3 method for lm vcovBS(x, cluster = NULL, R = 250, type = "xy", ..., fix = FALSE, use = "pairwise.complete.obs", applyfun = NULL, cores = NULL) # S3 method for glm vcovBS(x, cluster = NULL, R = 250, start = FALSE, ..., fix = FALSE, use = "pairwise.complete.obs", applyfun = NULL, cores = NULL)

x | a fitted model object. |
---|---|

cluster | a variable indicating the clustering of observations,
a |

R | integer. Number of bootstrap replications. |

start | logical. Should |

type | character or function. The character string specifies the type of
bootstrap to use: One of |

... | arguments passed to methods. For the default method, this is
passed to |

fix | logical. Should the covariance matrix be fixed to be positive semi-definite in case it is not? |

use | character. Specification passed to |

applyfun | an optional |

cores | numeric. If set to an integer the |

Clustered sandwich estimators are used to adjust inference when errors
are correlated within (but not between) clusters. See the documentation for `vcovCL`

for specifics about covariance clustering. This function allows
for clustering in arbitrary many cluster dimensions (e.g., firm, time, industry), given all
dimensions have enough clusters (for more details, see Cameron et al. 2011).
Unlike `vcovCL`

, `vcovBS`

uses a bootstrap rather than an asymptotic solution.

Basic (clustered) bootstrap covariance matrix estimation is provided by
the default `vcovBS`

method. It samples clusters (where each observation
is its own cluster by default), i.e., using case-based resampling. For obtaining
a covariance matrix estimate it is assumed that an `update`

of the model with the resampled `subset`

can be obtained, the `coef`

extracted, and finally the covariance computed with `cov`

.

The `update`

model is evaluated in the `environment(terms(x))`

(if available).
To speed up computations two further arguments can be leveraged.

Instead of

`lapply`

a parallelized function such as`parLapply`

or`mclapply`

can be specified to iterate over the bootstrap replications. For the latter, specifying`cores = ...`

is a convenience shortcut.When specifying

`start = TRUE`

, the`coef(x)`

are passed to`update`

as`start = coef(x)`

. This may not be supported by all model fitting functions and is hence not turned on by default.

The ``xy'' or ``pairs'' bootstrap is consistent for heteroscedasticity and clustered errors,
and converges to the asymptotic solution used in `vcovCL`

as `R`

, \(n\), and \(g\) become large (\(n\) and \(g\) are the number of
observations and the number of clusters, respectively; see Efron 1979, or Mammen 1992, for a
discussion of bootstrap asymptotics). For small \(g\)--particularly under 30 groups--the
bootstrap will converge to a slightly different value than the asymptotic method, due to
the limited number of distinct bootstrap replications possible (see Webb 2014 for a discussion
of this phenomonon). The bootstrap will not necessarily converge to an asymptotic estimate
that has been corrected for small samples.

The xy approach to bootstrapping is generally only of interest to the
practitioner when the asymptotic solution is unavailable (this can happen when using
estimators that have no `estfun`

function, for example). The residual bootstrap,
by contrast, is rarely of practical interest, because while it provides consistent
inference for for clustered standard errors, it is not robust to heteroscedasticity.
More generally, bootstrapping is useful when the bootstrap makes different assumptions than the asymptotic
estimator, in particular when the number of clusters is small and large \(n\) or
\(g\) assumptions are unreasonable. Bootstrapping is also often effective for nonlinear models,
particularly in smaller samples, where asymptotic approaches often perform relatively poorly.
See Cameron and Miller (2015) for further discussion of bootstrap techniques in practical applications,
and Zeileis et al. (2020) show simulations comparing `vcovBS`

to `vcovCL`

in several
settings.

The `glm`

method works essentially like the default method but call
`glm.fit`

instead of codeupdate.

The `lm`

method provides additional bootstrapping `type`

s
and computes the bootstrapped coefficient estimates somewhat more efficiently using
`lm.fit`

(for case-based resampling) or `qr.coef`

rather than `update`

. The default `type`

is case-based resampling
(`type = "xy"`

) as in the default method. Alternative `type`

specifications are:

`"residual"`

. The residual cluster bootstrap resamples the residuals (as above, by cluster) which are subsequently added to the fitted values to obtain the bootstrapped response variable: \(y^{*} = \hat{y} + e^{*}\). Coefficients can then be estimated using`qr.coef()`

, reusing the QR decomposition from the original fit. As Cameron et al. (2008) point out, the residual cluster bootstrap is not well-defined when the clusters are unbalanced as residuals from one cluster cannot be easily assigned to another cluster with different size. Hence a warning is issued in that case.`"wild"`

(or equivalently`"wild-rademacher"`

or`"rademacher"`

). The wild cluster bootstrap does not actually resample the residuals but instead reforms the dependent variable by multiplying the residual by a randomly drawn value and adding the result to the fitted value: \(y^{*} = \hat{y} + e \cdot w\) (see Cameron et al. 2008). By default, the factors are drawn from the Rademacher distribution:`function(n) sample(c(-1, 1), n, replace = TRUE)`

.`"mammen"`

(or`"wild-mammen"`

). This draws the wild bootstrap factors as suggested by Mammen (1993):`sample(c(-1, 1) * (sqrt(5) + c(-1, 1))/2, n, replace = TRUE, prob = (sqrt(5) + c(1, -1))/(2 * sqrt(5)))`

.`"webb"`

(or`"wild-webb"`

). This implements the six-point distribution suggested by Webb (2014), which may improve inference when the number of clusters is small:`sample(c(-sqrt((3:1)/2), sqrt((1:3)/2)), n, replace = TRUE)`

.`"norm"`

(or`"wild-norm"`

). The standard normal/Gaussian distribution is used for drawing the wild bootstrap factors:`function(n) rnorm(n)`

.User-defined function. This needs of the form as above, i.e., a

`function(n)`

returning a vector of random wild bootstrap factors of corresponding length.

A matrix containing the covariance matrix estimate.

Cameron AC, Gelbach JB, Miller DL (2008).
“Bootstrap-Based Improvements for Inference with Clustered Errors”,
*The Review of Economics and Statistics*, **90**(3), 414--427.
doi: 10.3386/t0344

Cameron AC, Gelbach JB, Miller DL (2011).
“Robust Inference with Multiway Clustering”,
*Journal of Business & Economic Statistics*, **29**(2), 238--249.
doi: 10.1198/jbes.2010.07136

Cameron AC, Miller DL (2015).
“A Practitioner's Guide to Cluster-Robust Inference”,
*Journal of Human Resources*, **50**(2), 317--372.
doi: 10.3368/jhr.50.2.317

Efron B (1979).
“Bootstrap Methods: Another Look at the Jackknife”,
*The Annals of Statistics*, **7**(1), 1--26.
doi: 10.1214/aos/1176344552

Mammen, E (1992).
“When Does Bootstrap Work?: Asymptotic Results and Simulations”,
*Lecture Notes in Statistics*, **77**.
Springer Science & Business Media.

Mammen, E (1993).
“Bootstrap and Wild Bootstrap for High Dimensional Linear Models”,
*The Annals of Statistics*, **21**(1), 255--285.
doi: 10.1214/aos/1176349025

Webb, MD (2014).
“Reworking Wild Bootstrap Based Inference for Clustered Errors”,
Working Paper 1315, *Queen's Economics Department.*
https://www.econ.queensu.ca/sites/econ.queensu.ca/files/qed_wp_1315.pdf.

Zeileis A, Köll S, Graham N (2020).
“Various Versatile Variances: An Object-Oriented Implementation of Clustered Covariances in R.”
*Journal of Statistical Software*, **95**(1), 1--36.
doi: 10.18637/jss.v095.i01

## Petersen's data data("PetersenCL", package = "sandwich") m <- lm(y ~ x, data = PetersenCL) ## comparison of different standard errors suppressWarnings(RNGversion("3.5.0")) set.seed(1) cbind( "classical" = sqrt(diag(vcov(m))), "HC-cluster" = sqrt(diag(vcovCL(m, cluster = ~ firm))), "BS-cluster" = sqrt(diag(vcovBS(m, cluster = ~ firm))) )#> classical HC-cluster BS-cluster #> (Intercept) 0.02835932 0.06701270 0.07067533 #> x 0.02858329 0.05059573 0.04878784## two-way wild cluster bootstrap with Mammen distribution vcovBS(m, cluster = ~ firm + year, type = "wild-mammen")#> (Intercept) x #> (Intercept) 0.003535619 0.000488708 #> x 0.000488708 0.002589785