context("Check that results agree with sandwich package")
# generate data
n = 300
x1 = runif(n)
x2 = rnorm(n)
z = sample(1:6, n, replace = T)
y = 2 + runif(1, -1, 1) * x1 + runif(1, -1, 1) * x2 + rnorm(n)
# fit model
m = lm(y ~ x1 + x2)
for (a in paste0("HC", 0:4)) {
test_that(paste0("Het: estimators agree on ", a), {
expect_equal(
jars::.robustSE(model.matrix(m), resid(m), type = a),
sandwich::vcovHC(m, type = a),
check.attributes = F
)
expect_equal(
vcov(summary(m, robust = TRUE, type = a)),
sandwich::vcovHC(m, type = a),
check.attributes = F
)
}
)
}
for (a in paste0("HC", 0:1)) {
test_that(paste0("Clust: estimators agree on ", a), {
expect_equal(
jars::.clusterSE(model.matrix(m), resid(m), clust = z, type = a),
sandwich::vcovCL(m, cluster = z, type = a),
check.attributes = F
)
expect_equal(
vcov(summary(m, cluster = z, type = a)),
sandwich::vcovCL(m, cluster = z, type = a),
check.attributes = F
)
}
)
}
context("Behavior When Missing Values are Present")
# generate some missing values
dat = data.frame(y, x1, x2)
n_miss = 20L
which_var = sample(1:3, n_miss, replace = T)
which_row = sample.int(n, n_miss, replace = F)
for (i in 1:20) dat[which_row[i], which_var[i]] = NA
# refit model
m = lm(y ~ x1 + x2, data = dat)
for (a in paste0("HC", 0:4)) {
test_that(paste0("Het: estimators agree on ",
a, " in the presence of NA"), {
expect_equal(
vcov(summary(m, robust = TRUE, type = a)),
sandwich::vcovHC(m, type = a),
check.attributes = F
)
}
)
}
for (a in paste0("HC", 0:1)) {
test_that(paste0("Clust: estimators agree on ",
a, " in the presence of NA (in vars)"), {
expect_equal(
vcov(summary(m, cluster = z, type = a)),
sandwich::vcovCL(m, cluster = z, type = a),
check.attributes = F
)
}
)
}
# add missing values to cluster var
z[sample.int(n, 5, replace = F)] = NA
for (a in paste0("HC", 0:1)) {
test_that(paste0("Clust: estimators agree on ",
a, " in the presence of NA (in vars & clust)"), {
expect_false(
isTRUE(all.equal(
vcov(summary(m, cluster = z, type = a)),
sandwich::vcovCL(m, cluster = z, type = a),
check.attributes = F
))
)
expect_warning(
vcov(summary(m, cluster = z, type = a))
)
}
)
}
# recode missing values in clusters as own cluster
z[is.na(z)] = 99
for (a in paste0("HC", 0:1)) {
test_that(
paste0("Clust: estimators agree on ",
a, " in the presence of NA (in vars & clust)"), {
expect_equal(
vcov(summary(m, cluster = z, type = a)),
sandwich::vcovCL(m, cluster = z, type = a),
check.attributes = F
)
}
)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.