如果不包括所有固定效应水平 (lme4),则排除组?
Exclude group if not including all fixed effect levels (lme4)?
我不确定在以下情况下该怎么做,其中某些级别的固定效果缺失(在随机效果内)- 它们是不平衡的。
想象一个有 5,000 条鱼的水族馆。它们是 100 种不同物种的一部分。我想测试他们的体重(连续)与他们是否由 Alan 或 Susie 喂食(只有两名员工喂鱼)之间是否存在关系。物种是随机效应。
我的模型如下所示:weight ~ employee + (1 + employee | species)
:具有随机截距和斜率的混合模型 (lmer
)。
但对于某些物种,所有鱼都由同一名员工(Alan 或 Susie)喂食。我应该将这些观察结果保留在模型中,还是应该排除它们?有这方面的文献吗?
这样应该没问题。混合模型非常适合这种缺失,除非它真的很极端(例如,有 没有 个物种,或者很少,这两个员工都测量过)。下面是一个虚构的小例子。
员工 1 的测量值缺失的情况的置信区间略宽;员工 2 的测量值缺失的情况对员工 2 的影响有相当广泛的置信区间(不确定为什么这些 不完全 为零,但我的猜测是它与模拟的特定随机效应值 - 即随机效应总体上具有零均值,因此这些可能略微 >0 以使总体估计平衡......?)
n_emp <- 2
n_spp <- 10
n_rep <- 20
dd <- expand.grid(emp = factor(seq(n_emp)),
spp = factor(seq(n_spp)),
rep = seq(n_rep))
dd2 <- subset(dd,
!((emp=="1" & (spp %in% 1:2)) |
(emp=="2" & (spp %in% 3:4))))
library(lme4)
form <- weight ~ emp + (1 + emp | spp)
## BUG/edge case: form[-2] breaks?
dd2$weight <- simulate(~ emp + (1 + emp | spp),
seed = 101,
newdata = dd2,
newparams = list(beta = c(1,1),
theta = c(1,1,1),
sigma = 1),
family = gaussian)[[1]]
m <- lmer(form, data = dd2)
rr <- as.data.frame(ranef(m))
rr$miss <- with(rr,
ifelse(grp %in% 1:2, "miss1",
ifelse(grp %in% 3:4, "miss2",
"nomiss")))
library(ggplot2)
ggplot(rr, aes(y = grp, x = condval, xmin = condval-2*condsd,
xmax = condval + 2*condsd, colour = miss)) +
geom_pointrange() +
facet_wrap(~term)
我不确定在以下情况下该怎么做,其中某些级别的固定效果缺失(在随机效果内)- 它们是不平衡的。
想象一个有 5,000 条鱼的水族馆。它们是 100 种不同物种的一部分。我想测试他们的体重(连续)与他们是否由 Alan 或 Susie 喂食(只有两名员工喂鱼)之间是否存在关系。物种是随机效应。
我的模型如下所示:weight ~ employee + (1 + employee | species)
:具有随机截距和斜率的混合模型 (lmer
)。
但对于某些物种,所有鱼都由同一名员工(Alan 或 Susie)喂食。我应该将这些观察结果保留在模型中,还是应该排除它们?有这方面的文献吗?
这样应该没问题。混合模型非常适合这种缺失,除非它真的很极端(例如,有 没有 个物种,或者很少,这两个员工都测量过)。下面是一个虚构的小例子。
员工 1 的测量值缺失的情况的置信区间略宽;员工 2 的测量值缺失的情况对员工 2 的影响有相当广泛的置信区间(不确定为什么这些 不完全 为零,但我的猜测是它与模拟的特定随机效应值 - 即随机效应总体上具有零均值,因此这些可能略微 >0 以使总体估计平衡......?)
n_emp <- 2
n_spp <- 10
n_rep <- 20
dd <- expand.grid(emp = factor(seq(n_emp)),
spp = factor(seq(n_spp)),
rep = seq(n_rep))
dd2 <- subset(dd,
!((emp=="1" & (spp %in% 1:2)) |
(emp=="2" & (spp %in% 3:4))))
library(lme4)
form <- weight ~ emp + (1 + emp | spp)
## BUG/edge case: form[-2] breaks?
dd2$weight <- simulate(~ emp + (1 + emp | spp),
seed = 101,
newdata = dd2,
newparams = list(beta = c(1,1),
theta = c(1,1,1),
sigma = 1),
family = gaussian)[[1]]
m <- lmer(form, data = dd2)
rr <- as.data.frame(ranef(m))
rr$miss <- with(rr,
ifelse(grp %in% 1:2, "miss1",
ifelse(grp %in% 3:4, "miss2",
"nomiss")))
library(ggplot2)
ggplot(rr, aes(y = grp, x = condval, xmin = condval-2*condsd,
xmax = condval + 2*condsd, colour = miss)) +
geom_pointrange() +
facet_wrap(~term)