在使用 nlme 包构建的线性混合模型的情况下,我有一个关于 lsmeans 包使用的自由度的问题。
这是一个基于 Oats 数据集来说明我的问题的示例。我不是想讨论这个模型在给定数据集的情况下是否相关,我只是想重现我在另一个数据集上遇到的问题;-)。
Oats.lme <- lme(yield ~ Variety, random = ~1 | Block, data = Oats)
anova(Oats.lme)
使用方差分析,我获得了预期的 64 个自由度。
numDF denDF F-value p-value
(Intercept) 1 64 245.1409 <.0001
Variety 2 64 1.6654 0.1972
然后我使用 lsmeans 函数:
lsmeans(Oats.lme, list(poly ~ Variety))
我得到
$`lsmeans of Variety`
Variety lsmean SE df lower.CL upper.CL
Golden Rain 104.5000 7.680866 5 84.75571 124.2443
Marvellous 109.7917 7.680866 5 90.04737 129.5360
Victory 97.6250 7.680866 5 77.88071 117.3693
Confidence level used: 0.95
$`polynomial contrasts of contrast`
contrast estimate SE df t.ratio p.value
linear -6.87500 6.68529 64 -1.028 0.3076
quadratic -17.45833 11.57926 64 -1.508 0.1365
对于对比,我获得了相同的 64 df,但对于 lsmeans 本身,我只有 5 df。我也使用 SAS,对于相同类型的模型,lsmeans 和 contrasts 的 df 数量相同(当前示例为 64)。
我已经看到使用 lme4 包时可能会改变自由度,但是我的代码嵌入在基于 nlme 的内部开发工具中,所以我基本上坚持使用 nlme。
现在有人会为什么会发生这种情况以及是否可以改变它?还是我错过了什么?
更新 - 初始错误消息
我最初注意到在一个特定情况下 lsmeans 的自由度降低了,我的随机运行效果只有 2 个级别,并且当我对 Dunnett 的调整感兴趣时。由于我对对比比对 lsmeans 更感兴趣,现在我了解了它的来源,我仍然可以使用它,但我把它放在那里以防万一有人遇到同样的错误并想知道原因。
我用 Oats 数据示例在下面复制了它。我获得的错误发生在 lsmeans:::.qdunnx 函数中,是由于 lsmeans 的 df 为 1。
Oats.lme <- lme(yield ~ Variety, random = ~1 | Block, data = subset(Oats,Block %in% c("I","II")))
lsm <- lsmeans(Oats.lme, trt.vs.ctrl ~ Variety)
summary(lsm,adjust = "dunnettx", infer = c(T, T), level = 0.95)
这是结果
$lsmeans
Variety lsmean SE df lower.CL upper.CL
Golden Rain 123.250 15.88642 1 -78.60608 325.1061
Marvellous 125.500 15.88642 1 -76.35608 327.3561
Victory 115.125 15.88642 1 -86.73108 316.9811
Confidence level used: 0.95
$contrasts
contrast estimate SE df t.ratio p.value
Marvellous - Golden Rain 2.250 12.8697 20 0.175 0.9695
Victory - Golden Rain -8.125 12.8697 20 -0.631 0.7482
P value adjustment: dunnettx method for 2 tests
Error in if (abs(diff(r[1:2])) < 5e-04) return(r[1]) :
missing value where TRUE/FALSE needed
In addition: Warning message:
In qtukey(p, (1 + sqrt(1 + 8 * k))/2, df) : production de NaN