knitr::opts_chunk$set(collapse = TRUE, comment = "#>", dev = "png", fig.width = 7, fig.height = 5, message = FALSE, warning = FALSE) if (!requireNamespace("sandwich", quietly = TRUE) || !requireNamespace("clubSandwich", quietly = TRUE)) { knitr::opts_chunk$set(eval = FALSE) }
This vignette demonstrate how to compute confidence intervals based on (cluster) robust variance-covariance matrices for standard errors.
First, we load the required packages and create a sample data set with a binomial and continuous variable as predictor as well as a group factor.
library(ggeffects) set.seed(123) # example taken from "?clubSandwich::vcovCR" 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(0.4, 0.3, -0.3) + nu + e dat <- data.frame(y, X, cluster, row = 1:n) # fit linear model model <- lm(y ~ X1 + X2 + X3, data = dat)
In this example, we use the normal standard errors, as returned by predict()
, to compute confidence intervals.
predict_response(model, "X1")
me <- predict_response(model, "X1") plot(me)
Now, we use sandwich::vcovHC()
to estimate heteroskedasticity-consistent standard errors. To do so, first the name of a related function must be supplied or the type of the HC-estimation as string.
E.g., to use the default sandwich::vcovHC()
function, set vcov = "HC"
, in which case the default type in sandwich::vcovHC()
is called. Setting vcov = "HC1"
is a convenient shortcut for vcov = "HC", vcov_args = list(type = "HC1")
, which would call sandwich::vcovHC(type = "HC1")
.
# short: predict_response(model, "X1", vcov = "HC0") # This is equivalent to the following: predict_response(model, "X1", vcov = "HC", vcov_args = list(type = "HC0"))
me <- predict_response(model, "X1", vcov = "HC", vcov_args = list(type = "HC0")) plot(me)
vcov
Instead of character strings, the vcov
argument also accepts a function that returns a variance-covariance matrix. Further arguments that need to be passed to that functions should be provided as list to the vcov_args
argument. Thus, we can rewrite the above code-chunk in the following way:
predict_response( model, "X1", vcov = sandwich::vcovHC, vcov_args = list(type = "HC0") )
The last example shows how to define cluster-robust standard errors. These are based on clubSandwich::vcovCR()
. Thus, vcov = "CR"
is always required when estimating cluster robust standard errors. clubSandwich::vcovCR()
has also different estimation types, which must be specified in vcov_args
. Furthermore, clubSandwich::vcovCR()
requires the cluster
-argument, which must also be specified in vcov_args
:
# short: # predict_response(model, "X1", vcov = "CR0", # vcov_args = list(cluster = dat$cluster)). # This is equivalent to the following: predict_response( model, "X1", vcov = "CR", vcov_args = list(type = "CR0", cluster = dat$cluster) )
me <- predict_response( model, "X1", vcov = "CR", vcov_args = list(type = "CR0", cluster = dat$cluster) ) plot(me)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.