tests/test_clubSandwich_interoperability.R

In plm: Linear Models for Panel Data

## Test of interoperability of vcov generated by package 'clubSandwich' with pwaldtest
## test of detection of those vcovs and test of translation function trans_clubSandwich_vcov
# library(plm)
# library(lmtest)
# 
# if (!requireNamespace("clubSandwich", quietly = TRUE)) {
# 	print("package 'clubSandwich not available locally, skipping associated test")
# } else {
# 	library(clubSandwich)
# 	data("MortalityRates", package = "clubSandwich")
# 	
# 	# partly from clubSandwich's vignette
# 	# subset for deaths in motor vehicle accidents, 1970-1983
# 	MV_deaths <- subset(MortalityRates, cause=="Motor Vehicle" & 
# 												year <= 1983 & !is.na(beertaxa))
# 	
# 	plm_unweighted <- plm(mrate ~ legal + beertaxa, data = MV_deaths, 
# 												effect = "individual", index = c("state","year"))
# 	
# 	plm_unweighted_re <- plm(mrate ~ legal + beertaxa, data = MV_deaths, model = "random",
# 													 effect = "twoways", index = c("state","year"), random.method = "walhus")
# 	
# 	coef_test(plm_unweighted, vcov = "CR1", cluster = "individual", test = "naive-t")
# 	coef_test(plm_unweighted, vcov = "CR1", cluster = "individual", test = "Satterthwaite")
# 	
# 	
# 	vcov_club_id   <- vcovCR(plm_unweighted, type = "CR0", cluster = "individual")
# 	vcov_club_time <- vcovCR(plm_unweighted, type = "CR0", cluster = "time")
# 	vcov_plm_id    <- vcovHC(plm_unweighted, method="arellano", type="HC0", cluster = "group")
# 	vcov_plm_time  <- vcovHC(plm_unweighted, method="arellano", type="HC0", cluster = "time")
# 	
# 	if (!isTRUE(all.equal(as.matrix(vcov_club_id),   vcov_plm_id,   check.attributes = FALSE))) stop("vcov_club_id and vcov_plm_id differ")
# 	if (!isTRUE(all.equal(as.matrix(vcov_club_time), vcov_plm_time, check.attributes = FALSE))) stop("vcov_club_time and vcov_plm_time differ")
# 	
# 	
# 	# some features:
# 	is.matrix(vcov_club_id) # TRUE
# 	class(vcov_club_id) # c("vcovCR", "clubSandwich")
# 	str(vcov_club_id)
# 	
# 	## clubSandwich's vcov holds the actual variable in attribute "cluster", 
# 	## not a character string like plm's vcov
# 	
# 	# Translation of attribute "cluster"
# 	attr(plm:::trans_clubSandwich_vcov(vcov_club_id,   attr(model.frame(plm_unweighted), "index")), "cluster")
# 	attr(plm:::trans_clubSandwich_vcov(vcov_club_time, attr(model.frame(plm_unweighted), "index")), "cluster")
# 	
# 	# Test if df adjustment in pwaldtest works (needs plm ver >= 1.5-34)
# 	df2_id <- pwaldtest(plm_unweighted, vcov = vcov_plm_id,  test = "F")[["parameter"]][2]
# 	if (!pwaldtest(plm_unweighted, vcov = vcov_club_id, test = "F")[["parameter"]][2] == df2_id) stop("wrong df2")
# 	
# 	df2_time <- pwaldtest(plm_unweighted, vcov = vcov_plm_time,  test = "F")[["parameter"]][2]
# 	if (!pwaldtest(plm_unweighted, vcov = vcov_club_time, test = "F")[["parameter"]][2] == df2_time) stop("wrong df2")
# 	
# 	attr(vcov_club_id, "cluster") <- "asdf"
# 	if (!is.null(attr(plm:::trans_clubSandwich_vcov(vcov_club_id, attr(model.frame(plm_unweighted), "index")), "cluster")))
# 		stop("attr 'cluster' other than known values found but not set to NULL")
# 	
# 	# should give warning
# 	attr(vcov_club_id, "cluster") <- NULL
# 	attr(plm:::trans_clubSandwich_vcov(vcov_club_id, attr(model.frame(plm_unweighted), "index")), "cluster")
# 	
# 	vcov_club_id   <- vcovCR(plm_unweighted, type = "CR0", cluster = "individual") # get "fresh" vcov
# 	
# 	# clubSandwich's Wald tests
# 	Wald_test(plm_unweighted, constraints = c("legal","beertaxa"), vcov = "CR1", test = "chi-sq")
# 	Wald_test(plm_unweighted, constraints = c("legal","beertaxa"), vcov = "CR0", test = "Naive-F")
# 	Wald_test(plm_unweighted, constraints = c("legal","beertaxa"), vcov = "CR0", test = "All")
# 	
# 	
# 	summary(plm_unweighted, vcov = vcov_plm_id)
# 	summary(plm_unweighted, vcov = vcov_club_id)  # just a run test
# 	summary(plm_unweighted, vcov = vcov_plm_time)
# 	summary(plm_unweighted, vcov = vcov_club_time) # just a run test
# 	
# 	# run test, gave error pre rev. 637 (in case an intercept is in the model):
# 	pwaldtest(plm_unweighted_re, vcov = vcovCR(plm_unweighted, type = "CR2"), test = "F")
# 	
# 	clubSandwich::coef_test(plm_unweighted, vcov = "CR0")
# 	clubSandwich::coef_test(plm_unweighted, vcov = "CR0", test = "naive-t")
# 	clubSandwich::coef_test(plm_unweighted, vcov = "CR0", test = "saddlepoint")
# 	
# 	clubSandwich::coef_test(plm_unweighted, vcov = vcov_club_id, test = "naive-t")
# 	clubSandwich::coef_test(plm_unweighted, vcov = "CR0", test = "All")
# 	lmtest::coeftest(plm_unweighted, vcov. = vcov_plm_id)
# 	lmtest::coeftest(plm_unweighted, vcov. = vcov_club_id)
# 	lmtest::coeftest(plm_unweighted, vcov. = function(x) vcovCR(x, type = "CR0", cluster = "individual"))
# 	# lmtest::coeftest(plm_unweighted, vcov. = vcovCR)# does not work due to missing default values
# }

### multiwayvcov
# library(multiwayvcov)
# lm_mod <- lm(mrate ~ legal + beertaxa + factor(state), data = MV_deaths)
# summary(lm_mod)
# vcov_multiwayvcov <- cluster.vcov(lm_mod, ~ year)
# lmtest::coeftest(lm_mod, vcov_multiwayvcov)
# str(vcov_multiwayvcov) # no additional information to infer vcov comes from multiwayvcov or about cluster variable
# class(vcov_multiwayvcov)
# on plm model
# Error in UseMethod("estfun") : 
#   no applicable method for 'estfun' applied to an object of class "c('plm', 'panelmodel')"

##### lfe
## robust
# library(lfe)
# data("Grunfeld", package = "plm")
# mod_felm_firm_robust <- felm(inv ~ value + capital | firm | 0 | firm, data = Grunfeld)
# str(mod_felm_firm_robust)
# 
# mod_felm_firm_robust$clustervar # in est. object: clustervar holds the actual variable, not just a character
# summary(mod_felm_firm_robust) # gives df adjusted F test for projected model
# lfe:::vcov.felm(mod_felm_firm_robust) 
# attributes(lfe:::vcov.felm(mod_felm_firm_robust)) # no special attributes -> no information about clustering in vcov object
# str(lfe:::vcov.felm(mod_felm_firm_robust))