将鼠标(Multiple Equations by Chained Equations)软件包更新到版本> 3后,我的代码停止工作。我希望从多重插补数据集上的线性回归中检索估计的方差-协方差矩阵。在2.46.0版中,使用pool函数可以轻松访问此数量(鼠标称为t)。在版本大于3.0的小鼠中,池函数不再返回完整的方差-协方差矩阵,它仅返回方差-协方差矩阵的对角元素。
这是一个工作示例:
首先创建一些缺少值的数据集:
set.seed(243)
iris$Sepal.Length[sample(length(iris$Sepal.Length), size = 5)] <- NA
iris$Sepal.Width[sample(length(iris$Sepal.Width), size = 5)] <- NA
iris$Petal.Length[sample(length(iris$Petal.Length), size = 5)] <- NA
iris$Species[sample(length(iris$Species), size = 5)] <- NA
第二次乘法估算丢失的数据
iris.mi <- mice(iris, 5)
第三,对每个乘插补数据集执行线性回归,将结果存储在mira对象中
mira.out <- with(iris.mi, lm(Sepal.Width ~ Sepal.Length + Petal.Length + Petal.Width + Species))
第四,使用鲁宾规则汇总这些分析的结果。这是通过小鼠的库功能实现的。
pool.out <- pool(analyses)
在mouse pack的2.46.0版中,可以通过输入来检索完整方差协方差矩阵t
pool.out$t
在mouses包的较新版本(> 3.0)中,不存在pool.out $ t对象。所有可以做的就是通过键入来检索差异
pool.out$pooled
并选择标记为t的列。似乎没有办法访问完整的方差-协方差矩阵。所有人都可以访问矩阵的对角元素,这些元素存储在pool.out $ pooled data.frame的t列中。
我想访问完整的方差协方差矩阵,因为我需要在具有多重估算数据的线性回归中计算交互项的边际效应和置信区间。这些置信区间可以通过仅使用方差-协方差矩阵的对角元素来近似,但是使用完整的方差-协方差矩阵会更有意义。
我想知道为什么在mouses包中实现了此更改,以及我如何能够在较新版本中访问方差-协方差矩阵。
感谢您的帮助。
最佳答案
至于为什么不再可用,我认为它与mice > 3.0
有关,具体取决于broom
功能来收集模型输出。它的小插图(https://cran.r-project.org/web/packages/broom/vignettes/broom.html)指出:“扫帚包获取R中内置函数(例如lm,nls或t.test)的混乱输出,并将它们转换为整洁的数据帧。”尽管确实尝试过,但这并不能轻易地容纳导出完整矩阵vcov
而不是当前显示的对角t
的对角线所必需的各个模型的t
矩阵。
作为一种变通方法,您可以使用以下代码(包括带有少量重命名的示例代码)自己派生它(例如,参见Dong和Peng 2013 http://springerplus.springeropen.com/articles/10.1186/2193-1801-2-222):
set.seed(243)
iris$Sepal.Length[sample(length(iris$Sepal.Length), size = 5)] <- NA
iris$Sepal.Width[sample(length(iris$Sepal.Width), size = 5)] <- NA
iris$Petal.Length[sample(length(iris$Petal.Length), size = 5)] <- NA
iris$Species[sample(length(iris$Species), size = 5)] <- NA
iris.mi <- mice(iris, 5)
fit.mi <- with(iris.mi, lm(Sepal.Width ~ Sepal.Length + Petal.Length + Petal.Width + Species))
fil.pooled <- pool(fit.mi)
# get the full matrix ubar (instead of only the diagonal)
m <- fil.pooled$m
ubar <- Reduce("+", lapply(fit.mi$analyses, vcov)) / (m)
b <- fil.pooled$pooled$b # this one is still provided by mice
# # or by hand as well
# qbar <- getqbar(fil.pooled) # pooled estimates
# b <- 1 / (m-1) * rowSums((sapply(fit.mi$analyses, coef) - qbar)^2)
t <- ubar + (1 + 1 / (m)) * b # this is t as it used to be
# check versus the diagonal of t that is still provided
all.equal(as.numeric(diag(t)), fil.pooled$pooled$t) # check
关于r - R中的鼠标包,更新到鼠标3.0后,mipo对象不再返回方差协方差矩阵,我们在Stack Overflow上找到一个类似的问题:https://stackoverflow.com/questions/51224252/