为什么 arm::sim 与 merTools::predictInterval 预测的置信区间不同?
Why are the coinfidence intervals predicted by arm::sim vs merTools::predictInterval different?
我正在比较 arm 的 sim()
函数生成的 confidence-interval (CI
) 和 merTools
生成的 predictInterval()
。
我使用来自 lme4
的 sleepstudy
数据集作为示例。
我期望这两种方法产生相同的结果,但事实并非如此。我所缺少的两种方法之间的根本区别是什么?
代码如下:
正在导入测试数据
sleepstudy <- as_tibble(sleepstudy) %>%
mutate(id = rep(1:18, each = 10)) %>%
dplyr::select(id, Days, Reaction) %>%
filter(id <= 16)
lme4 的多级模型
lmerfit <- lmer(Reaction ~ Days + (Days | id), data = sleepstudy)
生成预测
这是为了比较后面sim和preditInterval.
生成的中值
sleepstudy$predicted <- predict(lmerfit, newdata=sleepstudy, allow.new.levels=T)
CIs使用手臂:个人水平
sims <- sim(lmerfit, n.sims = 1000)
yhat <- fitted(sims, lmerfit)
sleepstudy$lower <- apply(yhat, 1, quantile, prob=0.025)
sleepstudy$median <- apply(yhat, 1, quantile, prob=0.5)
sleepstudy$upper <- apply(yhat, 1, quantile, prob=0.975)
CIs 使用 merTols
preds <- predictInterval(lmerfit,
newdata = sleepstudy,
n.sims = 1000,
include.resid.var=FALSE,
level=0.95,
stat="median")
sleepstudy <- cbind(sleepstudy, preds)
例如,我将第一个数据与两个不同的 CI 预测一起绘制。黑点是数据。红点是 lmerfit
的预测值。
黑线和黑色虚线分别是 arm::sim
的中位数和 95% CIs。
红线和虚线分别是 merTools::predictInterval
的中位数和 95% CIs。
预测值和模拟中值相同,但 CI 有很大差异。可能是什么原因?哪个准确?
ggplot(data = filter(sleepstudy, id == 1), aes(x=Days, y=Reaction)) +
geom_point() +
geom_point(aes(y=predicted), col = "red") +
geom_line(aes(y=median), col ="black" ) +
geom_line(aes(y=lower), col ="black", lty = 2) +
geom_line(aes(y=upper), col ="black", lty = 2) +
geom_line(aes(y=fit), col = "red") +
geom_line(aes(y=lwr), col = "red", lty = 2) +
geom_line(aes(y=upr), col = "red", lty = 2)
merTools CRAN 页面对此进行了介绍 (https://cran.r-project.org/web/packages/merTools/vignettes/Using_predictInterval.html),对 sim 和 predictInterval 进行了直接比较。基本上,我的理解是 sim 忽略了随机截距的不确定性,使用众数作为点估计。 predictInterval 的间隔更宽,因为它们考虑了这种额外的不确定性,因此可能更现实。
我正在比较 arm 的 sim()
函数生成的 confidence-interval (CI
) 和 merTools
生成的 predictInterval()
。
我使用来自 lme4
的 sleepstudy
数据集作为示例。
我期望这两种方法产生相同的结果,但事实并非如此。我所缺少的两种方法之间的根本区别是什么?
代码如下:
正在导入测试数据
sleepstudy <- as_tibble(sleepstudy) %>%
mutate(id = rep(1:18, each = 10)) %>%
dplyr::select(id, Days, Reaction) %>%
filter(id <= 16)
lme4 的多级模型
lmerfit <- lmer(Reaction ~ Days + (Days | id), data = sleepstudy)
生成预测
这是为了比较后面sim和preditInterval.
生成的中值sleepstudy$predicted <- predict(lmerfit, newdata=sleepstudy, allow.new.levels=T)
CIs使用手臂:个人水平
sims <- sim(lmerfit, n.sims = 1000)
yhat <- fitted(sims, lmerfit)
sleepstudy$lower <- apply(yhat, 1, quantile, prob=0.025)
sleepstudy$median <- apply(yhat, 1, quantile, prob=0.5)
sleepstudy$upper <- apply(yhat, 1, quantile, prob=0.975)
CIs 使用 merTols
preds <- predictInterval(lmerfit,
newdata = sleepstudy,
n.sims = 1000,
include.resid.var=FALSE,
level=0.95,
stat="median")
sleepstudy <- cbind(sleepstudy, preds)
例如,我将第一个数据与两个不同的 CI 预测一起绘制。黑点是数据。红点是 lmerfit
的预测值。
黑线和黑色虚线分别是 arm::sim
的中位数和 95% CIs。
红线和虚线分别是 merTools::predictInterval
的中位数和 95% CIs。
预测值和模拟中值相同,但 CI 有很大差异。可能是什么原因?哪个准确?
ggplot(data = filter(sleepstudy, id == 1), aes(x=Days, y=Reaction)) +
geom_point() +
geom_point(aes(y=predicted), col = "red") +
geom_line(aes(y=median), col ="black" ) +
geom_line(aes(y=lower), col ="black", lty = 2) +
geom_line(aes(y=upper), col ="black", lty = 2) +
geom_line(aes(y=fit), col = "red") +
geom_line(aes(y=lwr), col = "red", lty = 2) +
geom_line(aes(y=upr), col = "red", lty = 2)
merTools CRAN 页面对此进行了介绍 (https://cran.r-project.org/web/packages/merTools/vignettes/Using_predictInterval.html),对 sim 和 predictInterval 进行了直接比较。基本上,我的理解是 sim 忽略了随机截距的不确定性,使用众数作为点估计。 predictInterval 的间隔更宽,因为它们考虑了这种额外的不确定性,因此可能更现实。