This is a generic function, with specific methods defined for lm, plm, glm, gls, lme, robu, rma.uni, and rma.mv objects.

vcovCR returns a sandwich estimate of the variance-covariance matrix of a set of regression coefficient estimates.

## Usage

vcovCR(obj, cluster, type, target, inverse_var, form, ...)

# S3 method for default
vcovCR(
obj,
cluster,
type,
target = NULL,
inverse_var = FALSE,
form = "sandwich",
...
)

## Arguments

obj

Fitted model for which to calculate the variance-covariance matrix

cluster

Expression or vector indicating which observations belong to the same cluster. For some classes, the cluster will be detected automatically if not specified.

type

Character string specifying which small-sample adjustment should be used, with available options "CR0", "CR1", "CR1p", "CR1S", "CR2", or "CR3". See "Details" section of vcovCR for further information.

target

Optional matrix or vector describing the working variance-covariance model used to calculate the CR2 and CR4 adjustment matrices. If a vector, the target matrix is assumed to be diagonal. If not specified, vcovCR will attempt to infer a value.

inverse_var

Optional logical indicating whether the weights used in fitting the model are inverse-variance. If not specified, vcovCR will attempt to infer a value.

form

Controls the form of the returned matrix. The default "sandwich" will return the sandwich variance-covariance matrix. Alternately, setting form = "meat" will return only the meat of the sandwich and setting form = B, where B is a matrix of appropriate dimension, will return the sandwich variance-covariance matrix calculated using B as the bread. form = "estfun" will return the (appropriately scaled) estimating function, the transposed crossproduct of which is equal to the sandwich variance-covariance matrix.

...

Additional arguments available for some classes of objects.

## Value

An object of class c("vcovCR","clubSandwich"), which consists of a matrix of the estimated variance of and covariances between the regression coefficient estimates. The matrix has several attributes:

type

indicates which small-sample adjustment was used

cluster

contains the factor vector that defines independent clusters

v_scale

constant used in scaling the sandwich estimator

est_mats

contains a list of estimating matrices used to calculate the sandwich estimator

contains a list of adjustment matrices used to calculate the sandwich estimator

target

contains the working variance-covariance model used to calculate the adjustment matrices. This is needed for calculating small-sample corrections for Wald tests.

## Details

vcovCR returns a sandwich estimate of the variance-covariance matrix of a set of regression coefficient estimates.

Several different small sample corrections are available, which run parallel with the "HC" corrections for heteroskedasticity-consistent variance estimators, as implemented in vcovHC. The "CR2" adjustment is recommended (Pustejovsky & Tipton, 2017; Imbens & Kolesar, 2016). See Pustejovsky and Tipton (2017) and Cameron and Miller (2015) for further technical details. Available options include:

"CR0"

is the original form of the sandwich estimator (Liang & Zeger, 1986), which does not make any small-sample correction.

"CR1"

multiplies CR0 by m / (m - 1), where m is the number of clusters.

"CR1p"

multiplies CR0 by m / (m - p), where m is the number of clusters and p is the number of covariates.

"CR1S"

multiplies CR0 by (m (N-1)) / [(m - 1)(N - p)], where m is the number of clusters, N is the total number of observations, and p is the number of covariates. Some Stata commands use this correction by default.

"CR2"

is the "bias-reduced linearization" adjustment proposed by Bell and McCaffrey (2002) and further developed in Pustejovsky and Tipton (2017). The adjustment is chosen so that the variance-covariance estimator is exactly unbiased under a user-specified working model.

"CR3"

approximates the leave-one-cluster-out jackknife variance estimator (Bell & McCaffrey, 2002).

## References

Bell, R. M., & McCaffrey, D. F. (2002). Bias reduction in standard errors for linear regression with multi-stage samples. Survey Methodology, 28(2), 169-181.

Cameron, A. C., & Miller, D. L. (2015). A Practitioner's Guide to Cluster-Robust Inference. Journal of Human Resources, 50(2), 317-372. doi:10.3368/jhr.50.2.317

Imbens, G. W., & Kolesar, M. (2016). Robust standard errors in small samples: Some practical advice. Review of Economics and Statistics, 98(4), 701-712. doi:10.1162/rest_a_00552

Liang, K.-Y., & Zeger, S. L. (1986). Longitudinal data analysis using generalized linear models. Biometrika, 73(1), 13-22. doi:10.1093/biomet/73.1.13

Pustejovsky, J. E. & Tipton, E. (2018). Small sample methods for cluster-robust variance estimation and hypothesis testing in fixed effects models. Journal of Business and Economic Statistics, 36(4), 672-683. doi:10.1080/07350015.2016.1247004

vcovCR.lm, vcovCR.plm, vcovCR.glm, vcovCR.gls, vcovCR.lme, vcovCR.lmerMod, vcovCR.robu, vcovCR.rma.uni, vcovCR.rma.mv

## Examples


# simulate design with cluster-dependence
m <- 8
cluster <- factor(rep(LETTERS[1:m], 3 + rpois(m, 5)))
n <- length(cluster)
X <- matrix(rnorm(3 * n), n, 3)
nu <- rnorm(m)[cluster]
e <- rnorm(n)
y <- X %*% c(.4, .3, -.3) + nu + e
dat <- data.frame(y, X, cluster, row = 1:n)

# fit linear model
lm_fit <- lm(y ~ X1 + X2 + X3, data = dat)
vcov(lm_fit)
#>               (Intercept)            X1            X2           X3
#> (Intercept)  3.943130e-02  5.888022e-06  0.0093848910 -0.004472691
#> X1           5.888022e-06  3.684308e-02 -0.0009074324 -0.005102985
#> X2           9.384891e-03 -9.074324e-04  0.0450887035 -0.002456576
#> X3          -4.472691e-03 -5.102985e-03 -0.0024565764  0.041755459

# cluster-robust variance estimator with CR2 small-sample correction
vcovCR(lm_fit, cluster = dat$cluster, type = "CR2") #> (Intercept) X1 X2 X3 #> (Intercept) 0.209628936 0.059185662 -0.01929938 -0.007673005 #> X1 0.059185662 0.033754003 -0.01057876 -0.003128643 #> X2 -0.019299384 -0.010578760 0.03281238 -0.020750773 #> X3 -0.007673005 -0.003128643 -0.02075077 0.019237942 # compare small-sample adjustments CR_types <- paste0("CR",c("0","1","1S","2","3")) sapply(CR_types, function(type) sqrt(diag(vcovCR(lm_fit, cluster = dat$cluster, type = type))))
#>                   CR0       CR1      CR1S       CR2       CR3
#> (Intercept) 0.4162838 0.4450261 0.4567909 0.4578525 0.5048910
#> X1          0.1605749 0.1716618 0.1761999 0.1837226 0.2107155
#> X2          0.1571100 0.1679577 0.1723978 0.1811419 0.2103259
#> X3          0.1189471 0.1271598 0.1305215 0.1387009 0.1630821