CODEWITHSUNDEEP
rstatisticsanova
user9935785
user9935785

Reputation: 13

Repeated measures ANOVA and link to mixed-effect models in R

I have a problem when performing a two-way rm ANOVA in R on the following data (link : https://drive.google.com/open?id=1nIlFfijUm4Ib6TJoHUUNeEJnZnnNzO29):

subjectnbr is the id of the subject and blockType and linesTTL are the independent variables. RT2 is the dependent variable

I first performed the rm ANOVA through using ezANOVA with the following code: 

ANOVA_RTS <- ezANOVA(
    data=castRTs
    , dv=RT2
    , wid=subjectnbr
    , within = .(blockType,linesTTL)
    , type = 2
    , detailed = TRUE
    , return_aov = FALSE
)
ANOVA_RTS

The result is correct (I double-checked using statistica).

However, when I perform the rm ANOVA using the lme function, I do not get the same answer and I have no clue why.

There is my code: 

lmeRTs <- lme(
      RT2 ~ blockType*linesTTL, 
      random = ~1|subjectnbr/blockType/linesTTL, 
      data=castRTs)

anova(lmeRTs)

Here are the outputs of both ezANOVA and lme.

I hope I have been clear enough and have given you all the information needed.

I'm looking forward for your help as I am trying to figure it out for at least 4 hours!

Thanks in advance.

Upvotes: 1

Views: 610

Answers (1)

Maurits Evers
Maurits Evers

Reputation: 50728

Here is a step-by-step example on how to reproduce ezANOVA results with nlme::lme.

The data

We read in the data and ensure that all categorical variables are factors.

# Read in data
library(tidyverse);
df <- read.csv("castRTs.csv");
df <- df %>%
    mutate(
        blockType = factor(blockType),
        linesTTL = factor(linesTTL));

Results from ezANOVA

As a check, we reproduce the ez::ezANOVA results.

## ANOVA using ez::ezANOVA
library(ez);
model1 <- ezANOVA(
    data = df,
    dv = RT2,
    wid = subjectnbr,
    within = .(blockType, linesTTL),
    type = 2,
    detailed = TRUE,
    return_aov = FALSE);
model1;
#        $ANOVA
#              Effect DFn DFd          SSn       SSd           F            p
#1        (Intercept)   1  13 2047405.6654 34886.767 762.9332235 6.260010e-13
#2          blockType   1  13     236.5412  5011.442   0.6136028 4.474711e-01
#3           linesTTL   1  13    6584.7222  7294.620  11.7348665 4.514589e-03
#4 blockType:linesTTL   1  13    1019.1854  2521.860   5.2538251 3.922784e-02
#  p<.05         ges
#1     * 0.976293831
#2       0.004735442
#3     * 0.116958989
#4     * 0.020088855

Results from nlme::lme

We now run nlme::lme

## ANOVA using nlme::lme
library(nlme);
model2 <- anova(lme(
    RT2 ~ blockType * linesTTL,
    random = list(subjectnbr = pdBlocked(list(~1, pdIdent(~blockType - 1), pdIdent(~linesTTL - 1)))),
    data = df))
model2;
#                   numDF denDF  F-value p-value
#(Intercept)            1    39 762.9332  <.0001
#blockType              1    39   0.6136  0.4382
#linesTTL               1    39  11.7349  0.0015
#blockType:linesTTL     1    39   5.2538  0.0274

Results/conclusion

We can see that the F test results from both methods are identical. The somewhat complicated structure of the random effect definition in lme arises from the fact that you have two crossed random effects. Here "crossed" means that for every combination of blockType and linesTTL there exists an observation for every subjectnbr.

Some additional (optional) details

To understand the role of pdBlocked and pdIdent we need to take a look at the corresponding two-level mixed effect model

enter image description here

The predictor variables enter image description here are your categorical variables blockType and linesTTL, which are generally encoded using dummy variables.

The variance-covariance matrix for the random effects enter image description here can take different forms, depending on the underlying correlation structure of your random effect coefficients. To be consistent with the assumptions of a two-level repeated measure ANOVA, we must specify a block-diagonal variance-covariance matrix pdBlocked, where we create diagonal blocks for the offset ~1, and for the categorical predictor variables blockType pdIdent(~blockType - 1) and linesTTL pdIdent(~linesTTL - 1), respectively. Note that we need to subtract the offset from the last two blocks (since we've already accounted for the offset).

Some relevant/interesting resources

Pinheiro and Bates, Mixed-Effects Models in S and S-PLUS, Springer (2000)

Potvin and Schutz, Statistical power for the two-factor repeated measures ANOVA, Behavior Research Methods, Instruments & Computers, 32, 347-356 (2000)

Deming Mi, How to understand and apply mixed-effect models, Department of Biostatistics, Vanderbilt university

Upvotes: 0

Related Questions