jmv::ANOVA & summary(lm) output

iainjgallagher · Post by **iainjgallagher** » Tue Feb 02, 2021 2:34 pm

Hi

I'm writing a tutorial for some students to illustrate the ANOVA can be conceptualised as a linear model. Using jmv::ANOVA seems to mask the ability of the lm() function to properly label factors in the summary(lm) output if you try to run the lm() after jmv::ANOVA.

Code: Select all

library(jmv)

# make fake data
bia <- rnorm(10, 17.4, 6.6)
DW <- rnorm(10, 23, 5.9)
Gth <- rnorm(10, 21, 6.4)
Jck <- rnorm(10, 14.9, 7.1)

data_in <- data.frame(method = c( rep('bia', 10), rep('DW', 10), rep('Gth', 10), rep('Jck', 10)), value = c(bia, DW, Gth, Jck))
mod1 <- lm(value ~ method, data = data_in)
summary(mod1)

The output is:

Code: Select all

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   18.903      2.169   8.715 2.14e-10 ***
methodDW       1.672      3.068   0.545   0.5890    
methodGth      5.741      3.068   1.872   0.0694 .  
methodJck     -4.543      3.068  -1.481   0.1473

Easy to see which category each coefficient refers to.

Then:

Code: Select all

ANOVA(formula = value ~ method, data = data_in)

 ANOVA                                                                         
 ───────────────────────────────────────────────────────────────────────────── 
                Sum of Squares    df    Mean Square    F           p           
 ───────────────────────────────────────────────────────────────────────────── 
   method             512.7074     3      170.90245    4.896498    0.0058978   
   Residuals         1256.5079    36       34.90300                            
 ─────────────────────────────────────────────────────────────────────────────

Great!

Then:

Code: Select all

mod2 <- lm(value ~ method, data = data_in)
summary(mod2)

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)  18.8055     0.9341  20.132  < 2e-16 ***
method1      -2.9132     1.6179  -1.801  0.08016 .  
method2       5.4571     1.6179   3.373  0.00179 ** 
method3       0.9381     1.6179   0.580  0.56565

The 'method' labels (DW, Gth, Jck) vanish from the summary(lm) output when you run a new/subsequent lm(). I've tried this in RStudio and the console with the same result.

Code: Select all

sessionInfo()

R version 3.5.1 (2018-07-02)
Platform: x86_64-apple-darwin15.6.0 (64-bit)
Running under: macOS  10.14.6

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_GB.UTF-8/en_GB.UTF-8/en_GB.UTF-8/C/en_GB.UTF-8/en_GB.UTF-8

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

other attached packages:
[1] jmv_1.2.5

loaded via a namespace (and not attached):
 [1] zip_2.0.4         Rcpp_1.0.4.6      pillar_1.4.4      compiler_3.5.1    cellranger_1.1.0  forcats_0.5.0     base64enc_0.1-3  
 [8] tools_3.5.1       jsonlite_1.6.1    lifecycle_0.2.0   tibble_3.0.1      gtable_0.3.0      jmvcore_1.2.19    pkgconfig_2.0.3  
[15] rlang_0.4.8       openxlsx_4.1.5    rstudioapi_0.11   curl_4.3          haven_2.2.0       rio_0.5.16        dplyr_1.0.2      
[22] generics_0.0.2    vctrs_0.3.6       hms_0.5.3         grid_3.5.1        tidyselect_1.1.0  glue_1.4.1        data.table_1.12.8
[29] R6_2.4.1          readxl_1.3.1      foreign_0.8-76    carData_3.0-3     ggplot2_3.3.0     purrr_0.3.4       car_3.0-7        
[36] magrittr_1.5      scales_1.1.1      ellipsis_0.3.1    abind_1.4-5       colorspace_1.4-1  stringi_1.4.6     munsell_0.5.0    
[43] crayon_1.3.4

Can you help?

Best,

Iain

Post by **jonathon** » Wed Feb 03, 2021 2:22 am

the only thing i can think of that might be playing a role is that ANOVA() changes the contrasts:

https://github.com/jamovi/jmv/blob/master/R/ancova.b.R#L79

jonathon

iainjgallagher · Post by **iainjgallagher** » Wed Feb 03, 2021 9:55 am

Yes, that's it. Resetting the default contrasts before the second lm() fixes things.

Code: Select all

base::options(contrasts = c("contr.treatment","contr.poly"))

mod2 <- lm(value ~ method, data = data_in)
summary(mod2)

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept)   15.501      1.962   7.902 2.24e-09 ***
methodDW       5.872      2.774   2.116   0.0413 *  
methodGth      3.627      2.774   1.307   0.1994    
methodJck     -1.285      2.774  -0.463   0.6462

Thanks... I can work with that.

best,

i

jamovi