runehaubo / lmerTestR

Repository for the R-package lmerTest
48 stars 9 forks source link

update function reverting to lme4 #21

Closed ahmjebb1 closed 5 years ago

ahmjebb1 commented 5 years ago

Hi

I'm having problems when I try update on an lmer model to remove a fixed effect term. Somehow I am changing the summary output I obtain for the second model (after the update).

Can you help me? :) I included an example of code and output below!

Alexandra`


>m.juv.mass.speed <- lmerTest::lmer(speed ~ poly(mass,2)*sex + incline  + substrate + valley + age.days + lhf + (1|yrborn) + (1|colony) + (1|uid),
                                   data = juv.speed, REML = FALSE)

>summary(m.juv.mass.speed)

Linear mixed model fit by maximum likelihood t-tests use Satterthwaite
  approximations to degrees of freedom [lmerMod]
Formula: speed ~ poly(mass, 2) * sex + incline + substrate + valley +  
    age.days + lhf + (1 | yrborn) + (1 | colony) + (1 | uid)
   Data: juv.speed

     AIC      BIC   logLik deviance df.resid 
   819.7    882.5   -392.8    785.7      281 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-4.1871 -0.5894 -0.0609  0.5093  4.4341 

Random effects:
 Groups   Name        Variance Std.Dev.
 uid      (Intercept) 0.00000  0.0000  
 colony   (Intercept) 0.04342  0.2084  
 yrborn   (Intercept) 0.17366  0.4167  
 Residual             0.74590  0.8637  
Number of obs: 298, groups:  uid, 211; colony, 11; yrborn, 11

Fixed effects:
                      Estimate Std. Error         df t value Pr(>|t|)    
(Intercept)           3.227890   1.197598 291.661962   2.695 0.007441 ** 
poly(mass, 2)1        5.405079   2.215408 295.941403   2.440 0.015284 *  
poly(mass, 2)2       -2.985483   1.399798 290.811526  -2.133 0.033779 *  
sexM                 -0.070374   0.112435 291.234124  -0.626 0.531867    
incline              -0.006782   0.003945 280.954871  -1.719 0.086709 .  
substrateHV          -0.710809   0.190512 297.130001  -3.731 0.000228 ***
substrateLV          -0.131812   0.163409 293.408435  -0.807 0.420527    
substrateS           -0.739791   0.229598 287.276441  -3.222 0.001419 ** 
valleyuv             -0.188394   0.198440   8.166936  -0.949 0.369674    
age.days             -0.008113   0.008212 288.568398  -0.988 0.324002    
lhf                   0.020053   0.021262 297.675935   0.943 0.346375    
poly(mass, 2)1:sexM  -2.919722   2.069648 288.136440  -1.411 0.159402    
poly(mass, 2)2:sexM  -0.549925   1.973113 289.022641  -0.279 0.780667    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

>m.juv.mass.speed1 <- update(m.juv.mass.speed, ~. - incline)

>summary(m.juv.mass.speed1)

Linear mixed model fit by maximum likelihood  ['lmerMod']
Formula: speed ~ poly(mass, 2) + sex + substrate + valley + age.days +  
    lhf + (1 | yrborn) + (1 | colony) + (1 | uid) + poly(mass,      2):sex
   Data: juv.speed

     AIC      BIC   logLik deviance df.resid 
   834.9    894.4   -401.4    802.9      289 

Scaled residuals: 
    Min      1Q  Median      3Q     Max 
-4.1281 -0.6326 -0.0810  0.5138  4.4430 

Random effects:
 Groups   Name        Variance Std.Dev.
 uid      (Intercept) 0.00000  0.0000  
 colony   (Intercept) 0.03411  0.1847  
 yrborn   (Intercept) 0.15650  0.3956  
 Residual             0.74859  0.8652  
Number of obs: 305, groups:  uid, 215; colony, 11; yrborn, 11

Fixed effects:
                     Estimate Std. Error t value
(Intercept)          3.528505   1.167180   3.023
poly(mass, 2)1       5.986848   2.174571   2.753
poly(mass, 2)2      -3.004301   1.393610  -2.156
sexM                -0.078828   0.111552  -0.707
substrateHV         -0.724048   0.189908  -3.813
substrateLV         -0.142861   0.162130  -0.881
substrateS          -0.732200   0.227649  -3.216
valleyuv            -0.135783   0.185037  -0.734
age.days            -0.011076   0.008007  -1.383
lhf                  0.016731   0.020766   0.806
poly(mass, 2)1:sexM -2.701502   2.048046  -1.319
poly(mass, 2)2:sexM -0.590334   1.953080  -0.302

Correlation of Fixed Effects:
            (Intr) pl(,2)1 pl(,2)2 sexM   sbstHV sbstLV sbstrS vallyv ag.dys
ply(mss,2)1  0.576                                                          
ply(mss,2)2 -0.237 -0.215                                                   
sexM         0.161  0.029   0.006                                           
substrateHV -0.137 -0.074   0.111  -0.043                                   
substrateLV -0.090  0.018   0.101   0.045  0.593                            
substrateS  -0.148 -0.050   0.119   0.044  0.461  0.483                     
valleyuv    -0.035  0.104  -0.045  -0.060 -0.200 -0.187 -0.233              
age.days    -0.569 -0.613   0.217   0.106  0.098  0.000  0.110 -0.032       
lhf         -0.936 -0.464   0.184  -0.288  0.049  0.004  0.069 -0.004  0.300
ply(m,2)1:M  0.046 -0.488   0.027  -0.078  0.003 -0.001  0.016 -0.029 -0.020
ply(m,2)2:M  0.047  0.056  -0.658  -0.012 -0.075 -0.119 -0.087  0.041 -0.067
            lhf    p(,2)1:
ply(mss,2)1               
ply(mss,2)2               
sexM                      
substrateHV               
substrateLV               
substrateS                
valleyuv                  
age.days                  
lhf                       
ply(m,2)1:M -0.042        
ply(m,2)2:M -0.017 -0.079 
runehaubo commented 5 years ago

Hi,

I cannot reproduce with current version of lmerTest and since you appear to be using an older version, please update, try again and report back here using an lmerTest-dataset if the problem persist. Also please include the results of sessionInfo(). Here is what I get:

> data("TVbo", package="lmerTest")
> fm <- lmerTest::lmer(Coloursaturation ~ TVset + Picture + (1|Assessor:TVset) + (1|Assessor),
+            data=TVbo)
> # summary(fm)
> class(fm)
[1] "lmerModLmerTest"
attr(,"package")
[1] "lmerTest"
> fm2 <- update(fm, ~.-TVset)
> class(fm2)
[1] "lmerModLmerTest"
attr(,"package")
[1] "lmerTest"
> # summary(fm2)
> class(as(fm2, "lmerMod"))
[1] "lmerMod"
attr(,"package")
[1] "lme4"
> sessionInfo()
R version 3.5.1 (2018-07-02)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS  10.14

Matrix products: default
BLAS: /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib
LAPACK: /Library/Frameworks/R.framework/Versions/3.5/Resources/lib/libRlapack.dylib

locale:
[1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
[1] fortunes_1.5-4

loaded via a namespace (and not attached):
 [1] zip_1.0.0           Rcpp_0.12.18        cellranger_1.1.0    pillar_1.3.0       
 [5] compiler_3.5.1      nloptr_1.0.4        plyr_1.8.4          bindr_0.1.1        
 [9] forcats_0.3.0       tools_3.5.1         lme4_1.1-18-1       tibble_1.4.2       
[13] gtable_0.2.0        nlme_3.1-137        lattice_0.20-35     pkgconfig_2.0.2    
[17] rlang_0.2.2         openxlsx_4.1.0      Matrix_1.2-14       rstudioapi_0.7     
[21] curl_3.2            yaml_2.2.0          haven_1.1.2         rio_0.5.10         
[25] bindrcpp_0.2.2      dplyr_0.7.6         hms_0.4.2           lmerTest_3.0-1.9003
[29] grid_3.5.1          tidyselect_0.2.4    data.table_1.11.4   glue_1.3.0         
[33] R6_2.2.2            readxl_1.1.0        foreign_0.8-71      carData_3.0-1      
[37] minqa_1.2.4         ggplot2_3.0.0       purrr_0.2.5         car_3.0-2          
[41] magrittr_1.5        scales_1.0.0        MASS_7.3-50         splines_3.5.1      
[45] abind_1.4-5         assertthat_0.2.0    colorspace_1.3-2    numDeriv_2016.8-1  
[49] lazyeval_0.2.1      munsell_0.5.0       crayon_1.3.4       

You use lmerTest::lmer which suggests that you may not have called library(lmerTest) first but that does not make any difference in my example above.

I hope this helps, Rune

PS: using REML=FALSE and at the same time asking for a Satterthwaite's approximation makes no sense.