2

在使用 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
4

1 回答 1

2

该模型说响应变量受两种随机变化的影响:由于块的变化,以及由于品种的变化。每个品种的手段都包括这些变异来源;但是这些方法的比较排除了块变化,因为品种是在同一个块上比较的。

您只有 6 个块,因此估计块的变化有 5 个自由度,这就解释了多样性均值的自由度。比较有更多的自由度,因为您不必考虑块变化。

这里要考虑的另一件事是对nlme包的支持使用包含方法来获得自由度。这基本上涉及查看每种效应的自由度的最坏情况。如果您改为使用lme4包和lmer函数来拟合模型,lsmeans将使用 Satterthwaite 或 Kendall-Roger 方法来获得自由度,这些结果可能会更大一些。但是,均值的自由度仍将大大低于比较的自由度。

附录:SAS 结果

这是一些具有相同数据和模型的 SAS 代码:

proc mixed data = Oats;
  class Variety Block;
  model yield = Variety / ddfm = satterth;
  random Block;
  lsmeans Variety / tdiff;

...和 ​​lsmeans 结果:

                           Least Squares Means

                                   Standard
Effect     Variety     Estimate       Error      DF    t Value    Pr > |t|
Variety    Golden_R      104.50      7.6809    8.87      13.61      <.0001
Variety    Marvello      109.79      7.6809    8.87      14.29      <.0001
Variety    Victory      97.6250      7.6809    8.87      12.71      <.0001

                      Differences of Least Squares Means
                                           Standard
Effect    Variety    _Variety   Estimate      Error     DF   t Value   Pr > |t|

Variety   Golden_R   Marvello    -5.2917     6.6853     64     -0.79     0.4316
Variety   Golden_R   Victory      6.8750     6.6853     64      1.03     0.3076
Variety   Marvello   Victory     12.1667     6.6853     64      1.82     0.0734

请注意,当 Satterthwaite 方法用于自由度时,SAS 显示比较的 df 为 64,但均值本身仅显示 8.87 df。

如果ddfm在语句中省略了该选项model,则默认为 df 的包含方法,并且在两个表中都列出了 64 df。但是,我认为 SAS 在实施遏制方面是不正确的;请参阅我之前在 CrossValidated 中关于此主题的帖子:https ://stats.stackexchange.com/questions/140156/degrees-of-freedom-using-containment-method

于 2017-03-16T00:50:00.797 回答