我是新手 R 用户。我安装了 Zelig 4.1-3 版和 Amelia II 1.7 版。我对如何使用 R 包和函数获得组合的多重插补数据的自由度、t 统计量和 f 值感到困惑。
首先,我加载了 Amelia 和 Zelig:
require(Amelia)
require(Zelig)
然后,我加载了 Amelia 附带的示例数据:
data(freetrade)
我使用 amelia 函数为此数据集创建了 5 个插补。
a.out <- amelia(freetrade, m = 5, ts = "year", cs = "country")
然后,为了结合插补,我使用了 zelig 函数:
z.out.imp <- zelig(tariff ~ polity + pop + gdp.pc + year + country,
data = a.out$imputations, model = "ls" )
但是,当我使用此代码时,我得到的系数似乎是单个插补的系数,而不是组合集的系数:
summary(z.out.imp)
它们如下:
Coefficients:
Value Std. Error t-stat p-value
(Intercept) 2.766176e+03 6.670110e+02 4.1471215 0.0003572868
polity 1.645011e-01 3.078134e-01 0.5344183 0.5938286336
pop -6.079963e-08 6.518429e-08 -0.9327345 0.3774275934
gdp.pc -4.246794e-04 1.945866e-03 -0.2182470 0.8319093062
year -1.335563e+00 3.519513e-01 -3.7947390 0.0009787456
countryIndonesia -7.000319e+01 4.646330e+01 -1.5066343 0.1700377061
countryKorea -8.643855e+01 4.671629e+01 -1.8502870 0.0926657863
countryMalaysia -8.815182e+01 5.389486e+01 -1.6356256 0.1393312364
countryNepal -8.215250e+01 5.475828e+01 -1.5002753 0.1702129176
countryPakistan -4.349869e+01 5.149729e+01 -0.8446791 0.4238033944
countryPhilippines -8.088975e+01 5.320694e+01 -1.5202857 0.1673234716
countrySriLanka -7.668840e+01 5.695485e+01 -1.3464771 0.2161986616
countryThailand -7.400481e+01 5.186395e+01 -1.4269026 0.1903428838
Amelia 手册显示了组合多重插补数据集的一些系数应该是什么,尽管没有解释如何使用 R 获得所有这些系数(参见http://cran.r-project.org/web/的第 46 页包/Amelia/vignettes/amelia.pdf)
Complete DF = 167
DF: min = 10.36
avg = 18.81
max = 37.62
F( 2, 10.4) = 15.50
Prob > F = 0.0008
Value Std. Error t-stat p-value
polity -0.206 0.39 -0.53 0.61
pop -3.21 e-08 8.72e-09 3.68 0.004
gdp.pc -0.0027 0.000644 -4.28 0.000
Intercept 32.7 2.66 12.29 0.000
因为 amelia 函数使用蒙特卡罗模拟,所以我们可以预期运行之间的微小差异。然而,截距的巨大差异表明 zelig 函数返回的回归统计数据不是组合集。
Amelia 手册提供了以下代码:
> b.out <-NULL
> se.out <-NULL
> for(i in 1:a.out$m){
+ ols.out <- lm(tariff ~ polity + pop + gdp.pc, data = a.out$imputations[[i]])
+ b.out <- rbind(b.out, ols.out$coef)
+ se.out <-rbind(se.out, coef(summary(ols.out))[,2])
+ }
> combined.results<-mi.meld(q=b.out, se = se.out)
> combined.results
我尝试使用它。返回的结果非常接近第 46 页上显示的值和标准错误:
$q.mi
(Intercept) polity pop gdp.pc
[1,] 33.17325 -0.1499587 2.967196e-08 -0.002724229
$se.mi
(Intercept) polity pop gdp.pc
[1,] 2.116721 0.276968 6.061993e-09 0.0006596203
但是,它们不包括 t 统计量、自由度或 f 值。
R 中是否有可用的开源包或函数,以便我无需手动计算即可获得自由度、t 统计量和 f 值?
谢谢。