Ctrl+K

Extras

Extras#

R coding conventions for interactions#

In fact we can write this code more compactly, as R will automatically include the main effects for the two variables as well as the interaction, if we use the * to denote which variables we want to model an interaction for. For example, we obtain the same output with the more compact coding here:

%%R
model.int<-lm(CognitionB ~ VisitNum*Sex, dat = cogDat)
summary(model.int)

Call:
lm(formula = CognitionB ~ VisitNum * Sex, data = cogDat)

Residuals:
    Min      1Q  Median      3Q     Max 
-12.937  -2.827   0.218   2.982  12.064 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)     7.7202     0.4353  17.737  < 2e-16 ***
VisitNum        0.3977     0.1196   3.325 0.000931 ***
SexM            1.2563     0.6497   1.934 0.053556 .  
VisitNum:SexM  -0.4590     0.1858  -2.471 0.013729 *  
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 4.527 on 693 degrees of freedom
Multiple R-squared:  0.01617,	Adjusted R-squared:  0.01192 
F-statistic: 3.798 on 3 and 693 DF,  p-value: 0.01017

If in fact, we want to include just the interaction without the main effect terms, we can use “:” instead.

For example:

%%R
model.int<-lm(CognitionB ~ VisitNum:Sex, dat = cogDat)
summary(model.int)

Call:
lm(formula = CognitionB ~ VisitNum:Sex, data = cogDat)

Residuals:
     Min       1Q   Median       3Q      Max 
-12.7225  -2.9448   0.2117   3.0424  11.6298 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)    8.28413    0.32378  25.585  < 2e-16 ***
VisitNum:SexF  0.26797    0.09922   2.701  0.00709 ** 
VisitNum:SexM  0.11395    0.10974   1.038  0.29946    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 4.536 on 694 degrees of freedom
Multiple R-squared:  0.01087,	Adjusted R-squared:  0.008016 
F-statistic: 3.812 on 2 and 694 DF,  p-value: 0.02257

Because we have omitted the main effect terms, we need two interaction terms to capture the sex specific effects (i.e. we need two regression coefficients to enable us to estimate a female-specific slope and a male-specific slope). When we have age as a main effect, the regression coefficient is equivalent to the “VisitNum:SexF” variable. As shown here

%%R
model.int<-lm(CognitionB ~ Sex + VisitNum:Sex, dat = cogDat)
summary(model.int)

Call:
lm(formula = CognitionB ~ Sex + VisitNum:Sex, data = cogDat)

Residuals:
    Min      1Q  Median      3Q     Max 
-12.937  -2.827   0.218   2.982  12.064 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept)    7.72023    0.43527  17.737  < 2e-16 ***
SexM           1.25633    0.64970   1.934 0.053556 .  
SexF:VisitNum  0.39774    0.11962   3.325 0.000931 ***
SexM:VisitNum -0.06122    0.14213  -0.431 0.666781    
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 4.527 on 693 degrees of freedom
Multiple R-squared:  0.01617,	Adjusted R-squared:  0.01192 
F-statistic: 3.798 on 3 and 693 DF,  p-value: 0.01017

In general though, it is advisable to have the main effects for each predictor variable as well as the interaction, to ensure that effects are correctly attributed to the right source.

More complex mixed effects models#

Take a look at the lme4 vignette for more details on how to specify more complex mixed effect models with this package.

Also this post: https://rstudio-pubs-static.s3.amazonaws.com/63556_e35cc7e2dfb54a5bb551f3fa4b3ec4ae.html

Notes on REML here: http://users.stat.umn.edu/~gary/classes/5303/handouts/REML.pdf

A common error message when using lmer() is Error in KhatriRao(sm, t(mm)) : (p \<- ncol(X)) == ncol(Y) is not TRUE

If you get this error, try removing observations with missing data. While lm() and glm() were good at automatically handling the presence of these lmer throws an arguably confusing error.

Contents