1

我有一个数据集,我正在使用 lme4 拟合混合模型回归。

dat <- structure(list(dv280 = c(41L, 68L, 0L, 6L, 20L, 30L, 8L, 1L, 
15L, NA, 59L, 5L, 21L, 41L, 11L, 14L, -2L, 20L, 25L, 33L, 32L, 
30L, 68L, 16L, 11L, -1L, 8L, 0L), v0 = c(55L, 90L, 30L, 23L, 
74L, 48L, 25L, 25L, 46L, NA, 60L, 69L, 55L, 41L, 34L, 41L, 53L, 
76L, 72L, 64L, 34L, 37L, 75L, 21L, 26L, 14L, 24L, 19L), treatment = structure(c(2L, 
1L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 2L, 2L, 1L, 1L, 
2L, 1L, 2L, 2L, 1L, 2L, 1L, 1L, 2L, 2L, 1L), .Label = c("hc", 
"nhc"), class = "factor"), cse = structure(c(2, 2, 6, 6, -4, 
-4, 5, 5, NA, NA, -4, -4, -3, -3, -2, -2, 3, 3, 2, 2, -4, -4, 
-7, -7, 4, 4, 2, 2), .Dim = c(28L, 1L)), pp = structure(c(1L, 
1L, 2L, 2L, 3L, 3L, 4L, 4L, 5L, 5L, 6L, 6L, 7L, 7L, 8L, 8L, 9L, 
9L, 10L, 10L, 11L, 11L, 12L, 12L, 13L, 13L, 14L, 14L), .Label = c("j.1", 
"j.3", "j.6", "j.11", "j.13", "j.16", "j.17", "j.18", "j.19", 
"j.22", "j.24", "j.30", "j.32", "j.36"), class = "factor")), .Names = c("dv280", 
"v0", "treatment", "cse", "pp"), row.names = c(NA, 28L), class = "data.frame")
head(dat)
require(lme4)
m <- lmer(dv280 ~ 1 + v0:treatment + cse + (0 + v0 | pp), data=dat, REML=TRUE)
summary(m)

接下来我想绘制结果拟合,因此使用来自http://glmm.wikidot.com/faq的示例代码,我为原始值生成了预测。

newdat <- data.frame(
  v0=dat$v0,
  treatment=dat$treatment, 
  cse=dat$cse,
  dv280=0)
newdat <- newdat[-c(9,10),]
mm <- model.matrix(terms(m), newdat)
newdat$dv280 <- mm %*% fixef(m)
pvar1 <- diag(mm %*% tcrossprod(vcov(m), mm))
tvar1 <- pvar1 + VarCorr(m)$pp[1]
newdat <- data.frame(newdat, plo=newdat$dv280 - 2 * sqrt(pvar1), 
                     phi=newdat$dv280 + 2 * sqrt(pvar1), 
                     tlo=newdat$dv280 - 2 * sqrt(tvar1), 
                     thi=newdat$dv280 + 2 * sqrt(tvar1)) 

这很容易用 ggplot 绘制:

p <- ggplot(data=newdat, mapping=aes(x=v0, y=dv280, colour=treatment)) +
  geom_point() +
  geom_smooth(method='lm', se=TRUE) +
  scale_colour_discrete(guide=guide_legend(title.position='left', title.hjust=1))
p + .mytheme + coord_cartesian(xlim=c(-20,100)) +     
  geom_hline(yintercept=0, colour='gray35', linetype='dashed') +
  geom_vline(xintercept=0, colour='gray35', linetype='dashed') 

但如图所示,ggplot 绘制的回归线的截距与估计的回归截距相比是关闭的;至少 Nhc 线的情况是这样,它显然是负截距,而常见的估计截距是 0.54。

阴谋

我的错误很容易(我认为):使用 geom_smooth 导致的拟合与回归所说的不同,因为 geom_smooth 每次处理都单独获取数据,并且在没有拟合模型的上下文的情况下以面值计算。但我不知道如何正确绘制线条。

4

1 回答 1

1

您需要考虑cse. 如果我在 中设置cse=0newdat则这些行将通过截距:

p <- ggplot(data=newdat, mapping=aes(x=v0, y=dv280, colour=treatment)) +
  geom_point() +
  geom_smooth(method='lm', se=TRUE, fullrange=TRUE) +
  scale_colour_discrete(guide=guide_legend(title.position='left', title.hjust=1))
p + theme_bw() + coord_cartesian(xlim=c(-1,1), ylim=c(0.4,0.65)) +     
  geom_hline(yintercept=0, colour='gray35', linetype='dashed') +
  geom_vline(xintercept=0, colour='gray35', linetype='dashed') 

在此处输入图像描述

请注意,通常不建议在模型中包含没有相应主效应的交互作用。

编辑:

根据您的评论,您可能想要这样的东西(使用您的原件newdat):

p <- ggplot(data=newdat, mapping=aes(x=v0, y=dv280, ymin=tlo, ymax=thi, colour=treatment, fill=treatment)) +
  geom_ribbon(alpha=0.2, aes(colour=NULL)) +
  geom_point() +
  geom_line() +
  scale_colour_discrete(guide=guide_legend(title.position='left', title.hjust=1))
p + theme_bw() +   
  geom_hline(yintercept=0, colour='gray35', linetype='dashed') +
  geom_vline(xintercept=0, colour='gray35', linetype='dashed') 

在此处输入图像描述

然而,这个图似乎不是很有用(尤其是点之间的线),因为它只描绘了 dv280-v0-plane 上的投影,您甚至无法以这种方式看到线性关系。

于 2013-12-13T11:05:27.960 回答