I'm trying to perform a Wald test on an interaction term in a survey-weighted model with imputed data - see below for a toy reprex.
When I run the last three lines of code with regTermTest
using syntax modeled directly off the function documentation I get the following error: Error in terms.default(model) : no terms component nor attribute
.
A quick Google search seems to suggest this error means I'm passing an unsupported object type to the function; I read this as perhaps regTermTest
does not support passing an MIResult
to model
. However, it seems that as of version 3.6 of the survey
package, regTermTest
supports MIResult
models? This page also seems to suggest that as well.
Appreciate any guidance on what I'm doing wrong. Alternatively I could be happy if someone knows how to get p-values on individual model terms from an MIResult
object (e.g. like is shown in this post for a regular model
object).
# load packages
library(tidyverse)
library(survey)
library(mi)
library(mitools)
# load data on school performance included in survey package
# documentation available here: https://r-survey.r-forge.r-project.org/survey/html/api.html
data(api)
# remove problematic variables that are unnecessary for this example
apisub <- apiclus1 %>% select(-c("name", "sname", "dname", "cname", "flag",
"acs.46", "acs.core"))
# create and update variable types in missing_data.frame
mdf <- missing_data.frame(apisub)
mdf <- change(mdf, "cds", what = "type", to = "irrelevant")
mdf <- change(mdf, "stype", what = "type", to = "irrelevant")
mdf <- change(mdf, "snum", what = "type", to = "irrelevant")
mdf <- change(mdf, "dnum", what = "type", to = "irrelevant")
mdf <- change(mdf, "cnum", what = "type", to = "irrelevant")
mdf <- change(mdf, "fpc", what = "type", to = "irrelevant")
mdf <- change(mdf, "pw", what = "type", to = "irrelevant")
# summarize the missing_data.frame
show(mdf)
# impute missing data
imputations <- mi(mdf)
# create imputation list to pass to svydesign
imp_list <- complete(imputations, m = 5)
# create complex survey design using imputed data
dsn <- svydesign(id = ~dnum,
weights = ~pw,
data = imputationList(imp_list),
fpc = ~fpc)
# subset the survey design to remove schools that did not meet both targets
# just as an example of subsetting
dsn_sub <- subset(dsn, both == "No")
# specify analytic model
anl <- with(dsn_sub,
svyglm(api99 ~ enroll + meals + avg.ed*ell,
family = gaussian(),
design = dsn
)
)
# combine results into a single output
res <- MIcombine(anl)
# perform wald test for main and ixn terms
regTermTest(res, ~meals)
regTermTest(res, ~avg.ed:ell)
regTermTest(res, ~avg.ed*ell)
Didn't figure out how to do this with regTermTest
but I learned that you can run the following to run a global test for interaction instead:
library(aod)
> aod::wald.test(Sigma = vcov(res),
+ b = coef(res),
+ Terms = 6)
Wald test:
----------
Chi-squared test:
X2 = 0.031, df = 1, P(> X2) = 0.86