与lmerTest一起使用时anova()不显示p值 [英] anova() does not display p-value when used with lmerTest
问题描述
我正在尝试使用lmerTest为固定效果设置p值. 我有4种不同的随机截距,3种交叉截距和1种嵌套:
I'm trying to use lmerTest to have p-values for my fixed effects. I have 4 different random intercepts, 3 crossed and one nested :
test.reml <- lmerTest::lmer(y ~ s1 + min + cot + min:cot + ge
+ vis + dur + mo + nps + dist + st1 + st2 + di1 + s1:cot
+ s1:min + s1:cot:min + s1:ge + s1:vis + s1:dur + s1:mo
+ s1:nps + s1:dist + s1:st1 + s1:st2 + s1:di1 + (1|Unique_key)
+ (s1-1|object) + (ns1-1|object)
+ (1|region), bdr, REML=1)
对对象进行两次观察,并通过对Unique_key(区域j中对象i的唯一标识符)的随机效应来引入两个度量之间的相关性.每个物体都可以在任何区域被观察到. S1是一个二进制变量,如果在第一个时间段观察到观测值,则取值为1,否则取值为0.每个对象的第一个周期有一个随机截距,第二个周期有一个随机截距.实际上,ns1是一个二进制变量,它是每个观察值的s1和s1 + ns1 = 1的补码.
The objects are observed two times and the correlation between the two measures is introduced by the random effect on Unique_key, an unique identifier of the object i in a region j. Each object could be observed in any region. S1 is a binary variable that takes the value 1 if the observation is observed the first time period and 0 either. There is one random intercept for the first periode and one random intercept for the second period for each object. ns1 is in fact a binary variable that is the complement of s1 and s1 + ns1 = 1 for each observation.
我可以拟合模型并使用summary()获得估计值和p值:
I can fit the model and get the estimates and p-values with summary() :
summary( test.reml)
Linear mixed model fit by REML ['merModLmerTest']
Formula: y ~ s1 + min + cot + min:cot + ge
+ vis + dur + mo + nps + dist + st1 + st2 + di1 + s1:cot
+ s1:min + s1:cot:min + s1:ge + s1:vis + s1:dur + s1:mo
+ s1:nps + s1:dist + s1:st1 + s1:st2 + s1:di1 + (1|Unique_key)
+ (s1-1|object) + (ns1-1|object)
+ (1|region), bdr, REML=1)
Data: bdr
REML criterion at convergence: 204569.1
Random effects:
Groups Name Variance Std.Dev.
Unique_key (Intercept) 0.2023 0.4497
object s1 0.3528 0.5940
object.1 ns1 0.5954 0.7716
Region (Intercept) 0.7563 0.8697
Residual 0.1795 0.4237
Number of obs: 113396, groups: Unique_key , 58541; object, 1065; Region, 87
Fixed effects:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 6.7341569 0.2382673 28.263 < 2e-16 ***
s1 0.7391924 0.2004413 3.688 0.000233 ***
min -0.0067606 0.0171385 -0.394 0.694205
cot 0.1235093 0.0353693 3.492 0.000499 ***
ge2 -0.1535452 0.0800998 -1.917 0.055525 .
ge3 -0.2131246 0.0986559 -2.160 0.030982 *
ge4 -0.1032694 0.1115603 -0.926 0.354830
ge5 -0.1769347 0.1296558 -1.365 0.172663
ge6 0.0117401 0.1115897 0.105 0.916231
ge7 -0.2692483 0.1022565 -2.633 0.008589 **
vis2 -0.0928661 0.0607950 -1.528 0.126938
vis3 -0.3026112 0.1246595 -2.428 0.015375 *
dur2 0.1479195 0.0786369 1.881 0.060249 .
dur3 0.1406340 0.0809379 1.738 0.082590 .
dur4 0.2742243 0.0884301 3.101 0.001981 **
dur5 0.1946761 0.1065815 1.827 0.068059 .
mo2 -0.1168591 0.1256017 -0.930 0.352386
mo3 -0.0611162 0.1267657 -0.482 0.629824
mo4 -0.2725720 0.1263740 -2.157 0.031248 *
mo5 -0.6107000 0.1379264 -4.428 1.05e-05 ***
mo6 -0.3635142 0.1299799 -2.797 0.005260 **
mo7 -0.0899233 0.1275164 -0.705 0.480846
mo8 -0.2349548 0.1253422 -1.875 0.061140 .
mo9 -0.2624888 0.1263051 -2.078 0.037934 *
mo10 -0.2882749 0.1244404 -2.317 0.020724 *
mo11 -0.1702823 0.1356031 -1.256 0.209497
mo12 0.1989155 0.1322339 1.504 0.132819
nps 0.0278418 0.0010393 26.790 < 2e-16 ***
dist2 0.4065093 0.1118916 3.633 0.000294 ***
dist3 0.0155691 0.0906664 0.172 0.863693
dist4 -0.2910960 0.1595805 -1.824 0.068424 .
dist5 -0.1316553 0.0913394 -1.441 0.149782
dist6 0.0477956 0.0995679 0.480 0.631308
dist7 0.1383000 0.0981247 1.409 0.159011
dist8 -0.3985620 0.0886316 -4.497 7.69e-06 ***
dist9 -0.2036683 0.0799584 -2.547 0.011005 *
st11 -0.0258775 0.0591631 -0.437 0.661919
st21 0.0089230 0.0573352 0.156 0.876356
di11 -0.0910207 0.0838321 -1.086 0.277846
min:cot 0.0066210 0.0006195 10.688 < 2e-16 ***
s1:cot -0.1505670 0.0443186 -3.397 0.000694 ***
s1:min 0.0079478 0.0015051 5.280 1.29e-07 ***
s1:ge2 0.0329272 0.1007943 0.327 0.743948
s1:ge3 0.2150927 0.1241590 1.732 0.083367 .
s1:ge4 0.1786057 0.1404119 1.272 0.203526
s1:ge5 -0.0422380 0.1631757 -0.259 0.795780
s1:ge6 0.1372051 0.1404415 0.977 0.328717
s1:ge7 0.1343314 0.1287059 1.044 0.296755
s1:vis2 0.1354091 0.0765084 1.770 0.076913 .
s1:vis3 0.2449180 0.1568745 1.561 0.118637
s1:dur2 -0.0888179 0.0989573 -0.898 0.369547
s1:dur3 -0.0532473 0.1018481 -0.523 0.601167
s1:dur4 -0.1239068 0.1112907 -1.113 0.265696
s1:dur5 -0.1191069 0.1341435 -0.888 0.374705
s1:mo2 -0.1357615 0.1574365 -0.862 0.388618
s1:mo3 0.0130976 0.1588743 0.082 0.934306
s1:mo4 0.0343900 0.1579532 0.218 0.827669
s1:mo5 0.2257241 0.1732449 1.303 0.192761
s1:mo6 0.0500347 0.1628755 0.307 0.758728
s1:mo7 -0.0451271 0.1596277 -0.283 0.777435
s1:mo8 -0.0200467 0.1572383 -0.127 0.898564
s1:mo9 0.0394005 0.1584268 0.249 0.803620
s1:mo10 0.0641038 0.1562518 0.410 0.681662
s1:mo11 -0.3136235 0.1703456 -1.841 0.065764 .
s1:mo12 -0.7003775 0.1660455 -4.218 2.58e-05 ***
s1:nps -0.0095428 0.0013077 -7.297 4.31e-13 ***
s1:dist2 -0.3867962 0.1407463 -2.748 0.006050 **
s1:dist3 -0.0516400 0.1140519 -0.453 0.650762
s1:dist4 -0.0567491 0.2008542 -0.283 0.777562
s1:dist5 0.0025780 0.1147143 0.022 0.982073
s1:dist6 -0.1456445 0.1252219 -1.163 0.244940
s1:dist7 -0.0452712 0.1234110 -0.367 0.713785
s1:dist8 0.0546400 0.1114865 0.490 0.624117
s1:dist9 0.0540697 0.1000415 0.540 0.588934
s1:st11 0.0784027 0.0744677 1.053 0.292549
s1:st21 -0.0394419 0.0721720 -0.546 0.584788
s1:di11 0.0463040 0.1055326 0.439 0.660882
s1:min:cot -0.0012850 0.0006004 -2.140 0.032344 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
但是使用anova()我得到了:
but with anova() I get :
type3.bonmodele <- lmerTest::anova(test.reml, ddf="Satterthwaite")
Analysis of Variance Table
Df Sum Sq Mean Sq F value
s1 1 7.385 7.385 41.1448
min 1 0.081 0.081 0.4536
cot 1 29.384 29.384 163.7026
ge 6 25.198 4.200 23.3968
vis 2 0.464 0.232 1.2929
dur 4 22.763 5.691 31.7042
mo 11 15.581 1.416 7.8914
nps 1 234.535 234.535 1306.6487
dist 8 18.547 2.318 12.9162
st1 1 0.034 0.034 0.1879
st2 1 0.058 0.058 0.3220
di1 1 0.261 0.261 1.4549
min:cot 1 22.537 22.537 125.5611
s1:cot 1 9.146 9.146 50.9555
s1:min 1 18.383 18.383 102.4171
s1:ge 6 5.152 0.859 4.7843
s1:vis 2 1.698 0.849 4.7311
s1:dur 4 2.829 0.707 3.9404
s1:mo 11 8.157 0.742 4.1312
s1:nps 1 10.102 10.102 56.2803
s1:dist 8 2.233 0.279 1.5550
s1:st1 1 0.188 0.188 1.0481
s1:st2 1 0.046 0.046 0.2560
s1:di1 1 0.035 0.035 0.1927
s1:min:cot 1 0.822 0.822 4.5804
当我尝试删除三元交互时,anova()函数返回p值...我也尝试过拆分数据框并使模型适合一半的数据,而anova()可以很好地工作.
When I try to remove the triple interaction the anova() function returns the p-values...I have also tried to split my data frame and to fit the model on half the data and anova() works well to.
使用这些功能时没有警告,并且我也尝试更改ddf选项和方法,但似乎无济于事.
There is no warning when I use the functions and I have also tried to change the ddf option and the method but nothing seems to work.
这是我的会话信息:
R version 3.0.0 (2013-04-03)
Platform: x86_64-w64-mingw32/x64 (64-bit)
locale:
[1] LC_COLLATE=French_Canada.1252 LC_CTYPE=French_Canada.1252 LC_MONETARY=French_Canada.1252 LC_NUMERIC=C
[5] LC_TIME=French_Canada.1252
attached base packages:
[1] parallel splines stats graphics grDevices utils datasets methods base
other attached packages:
[1] ggplot2_0.9.3.1 snow_0.3-13 Snowball_0.0-10 xtable_1.7-1 lmerTest_2.0-0 pbkrtest_0.3-7 MASS_7.3-29
[8] papeR_0.3 gmodels_2.15.4 survival_2.37-4 nlme_3.1-111 car_2.0-19 lme4_1.1-1 Matrix_1.1-0
[15] lattice_0.20-15
loaded via a namespace (and not attached):
[1] bitops_1.0-6 caTools_1.16 cluster_1.14.4 colorspace_1.2-4 dichromat_2.0-0 digest_0.6.3
[7] gdata_2.13.2 gplots_2.12.1 grid_3.0.0 gtable_0.1.2 gtools_3.0.0 Hmisc_3.12-2
[13] KernSmooth_2.23-10 labeling_0.2 minqa_1.2.1 munsell_0.4.2 nnet_7.3-7 numDeriv_2012.9-1
[19] plyr_1.8 proto_0.3-10 RColorBrewer_1.0-5 RCurl_1.95-4.1 reshape2_1.2.2 rJava_0.9-4
[25] ROAuth_0.9.3 rpart_4.1-3 scales_0.2.3 stringr_0.6.2 tools_3.0.0 twitteR_1.1.7
我无法共享数据,但是如果需要,我可以添加更多信息! 我想对自由度使用Satterthwaite近似值,但是如果您有任何其他获取p值的建议,请分享! 非常感谢你!
I can't share the data but I can add more infos if needed! I wanted to use the Satterthwaite approximation for the degrees of freedom but if you have any other suggestions to get p-values please share! Thank you very much!
推荐答案
如果lmerTest中发生错误,则默认情况下将给出lme4中的方差分析.因此,在您的情况下,发生了一些错误,但是如果不对数据进行测试,很难说出什么.可能是由于grad函数的简单方法所致,这是默认设置.您可以尝试:anova(test.reml,method.grad ="Richardson").否则就像我不看例子一样很难说...
If some error occurs in lmerTest then by default the anova from lme4 will be given. So in your case some error occurred, but it is quite hard to say what without testing on the data. Probably it is due to the simple method for the grad function, which is the default. You may try: anova(test.reml, method.grad="Richardson"). Otherwise again like I said quite hard to say without looking at the example...
亚历山大·库兹涅佐娃(Alexandra Kuznetsova)
Alexandra Kuznetsova
这篇关于与lmerTest一起使用时anova()不显示p值的文章就介绍到这了,希望我们推荐的答案对大家有所帮助,也希望大家多多支持IT屋!