I would like to reproduce the standard output of a repeated measures anova in SPSS with R. For a balanced design (i.e. equal group sizes) I can exactly reproduce the SPSS output of the within-subjects contrasts. But when using an unbalanced (i.e. different group sizes) the results of the within-subjects contrasts are not exactly the same, especially for the main effect of the within-subject factor.
Create a test dataset: (un)comment the line cond to create (un)balanced dataset
set.seed(1234) t <- 4 idmax <- 40 myData <- data.frame(PID = rep(seq(from = 1, to = idmax, by = 1), t), #cond = rep(c(rep(c("A"), 20), rep(c("B"), 20)), t), #BALANCED cond = rep(c(rep(c("A"), 22), rep(c("B"), 18)), t), #UNBALANCED time = rep(x = 1:t, each = idmax), acc = sample(x = 1:100, size = idmax*t, replace = TRUE) ) myData <- within(myData, { PID <- factor(PID) cond <- factor(cond) time <- factor(time) acc <- acc/100 }) Using the ezAnova function of the ez package I can exactly reproduce the standard SPSS output of the Between-Subject Effects and the Within-Subject Effects for both balanced and unbalanced designs.
library(ez) ezANOVA( data = myData, dv = acc, wid = PID, within = c("time"), between = cond, type = 3) But I cannot reproduce the Within-subjects Contrasts table of the SPSS output for unbalanced designs.
For a balanced design I can reproduce the results using the aov function.
time.L <- C(myData$time, poly, 1) time.Q <- C(myData$time, poly, 2) time.C <- C(myData$time, poly, 3) mdl_balanced <- aov(acc ~ cond*time + Error(factor(PID)/(time.L + time.Q + time.C)), data = myData) summary(mdl_balanced) R Output:
Error: factor(PID) Df Sum Sq Mean Sq F value Pr(>F) cond 1 0.128 0.12826 1.524 0.225 Residuals 38 3.198 0.08416 Error: factor(PID):time.L Df Sum Sq Mean Sq F value Pr(>F) time 1 0.2999 0.29993 4.213 0.0471 * <--1 cond:time 1 0.0000 0.00003 0.000 0.9842 Residuals 38 2.7055 0.07120 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Error: factor(PID):time.Q Df Sum Sq Mean Sq F value Pr(>F) time 1 0.553 0.5534 6.210 0.0172 * <--2 cond:time 1 0.014 0.0135 0.152 0.6992 Residuals 38 3.386 0.0891 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Error: factor(PID):time.C Df Sum Sq Mean Sq F value Pr(>F) time 1 0.1208 0.12079 1.507 0.227 <--3 cond:time 1 0.0259 0.02588 0.323 0.573 Residuals 38 3.0454 0.08014 SPSS output
Tests of Within-Subjects Contrasts Measure: MEASURE_1 Source time Type III SS df Mean Square F Sig. time Linear ,300 1 ,300 4,213 ,047 <--1 Quadratic ,553 1 ,553 6,210 ,017 <--2 Cubic ,121 1 ,121 1,507 ,227 <--3 time * cond Linear 2,812E-005 1 2,812E-005 ,000 ,984 Quadratic ,014 1 ,014 ,152 ,699 Cubic ,026 1 ,026 ,323 ,573 Error(time) Linear 2,706 38 ,071 Quadratic 3,386 38 ,089 Cubic 3,045 38 ,080 But if I run the same aov code for the unbalanced design. The results in the R output for the within-subject factor time clearly differ from the SPSS output (indicated with <--). The R output for cond:time match the SPSS output.
Error: factor(PID) Df Sum Sq Mean Sq F value Pr(>F) cond 1 0.114 0.11398 1.348 0.253 Residuals 38 3.213 0.08454 Error: factor(PID):time.L Df Sum Sq Mean Sq F value Pr(>F) time 1 0.2999 0.29993 4.229 0.0466 * <--1 cond:time 1 0.0106 0.01063 0.150 0.7009 Residuals 38 2.6949 0.07092 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Error: factor(PID):time.Q Df Sum Sq Mean Sq F value Pr(>F) time 1 0.553 0.5534 6.202 0.0172 * <--2 cond:time 1 0.009 0.0093 0.104 0.7491 Residuals 38 3.391 0.0892 --- Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 Error: factor(PID):time.C Df Sum Sq Mean Sq F value Pr(>F) time 1 0.1208 0.12079 1.527 0.224 <--3 cond:time 1 0.0649 0.06486 0.820 0.371 Residuals 38 3.0064 0.07912 SPSS output
Tests of Within-Subjects Contrasts Measure: MEASURE_1 Source time Type III SS df Mean Square F Sig. time Linear ,286 1 ,286 4,030 ,052 <--1 Quadratic ,562 1 ,562 6,301 ,016 <--2 Cubic ,103 1 ,103 1,297 ,262 <--3 time * cond Linear ,011 1 ,011 ,150 ,701 Quadratic ,009 1 ,009 ,104 ,749 Cubic ,065 1 ,065 ,820 ,371 Error(time) Linear 2,695 38 ,071 Quadratic 3,391 38 ,089 Cubic 3,006 38 ,079 So, my question is: How can I compute the correct F-values for polynomial contrasts for the factor time in an unbalanced repeated measures anova design using R?
PS. I can use a linear-mixed effects model for the unbalanced design in combination with fit.contrast, but this produces t-values for the contrasts instead of F-values.
options(contrasts = c("contr.sum", "contr.poly")) library(nlme) model.lme <- lme(acc ~ cond*time, random = ~ 1 | PID, data = myData) anova(model.lme, type="marginal") library(gmodels) fit.contrast(model.lme, "time", t(contr.poly(4, c(1, 2, 3, 4)))) Estimate Std. Error t-value Pr(>|t|) time.L 0.08495364 0.0448776 1.893008 0.060892471 time.Q 0.11915404 0.0448776 2.655089 0.009063148 time.C 0.05090443 0.0448776 1.134295 0.259050482 Best Answer
I'm not familiar with the R package ez, but I can offer an alternative approach with the R package afex:
require(afex) fm <- aov_car(acc ~ cond+Error(PID/time), data = myData) Then you have the results that match up your SPSS output:
require(phia) testInteractions(fm$lm, custom=list(time=contr.poly(4, c(1, 2, 3, 4))[,1]), idata = fm$data[["idata"]]) Value Df test stat approx F num Df den Df Pr(>F) time1 0.084954 1 0.095882 4.0299 1 38 0.05185 . testInteractions(fm$lm, custom=list(time=contr.poly(4, c(1, 2, 3, 4))[,2]), idata = fm$data[["idata"]]) Value Df test stat approx F num Df den Df Pr(>F) time1 0.11915 1 0.14223 6.3011 1 38 0.01644 * testInteractions(fm$lm, custom=list(time=contr.poly(4, c(1, 2, 3, 4))[,3]), idata = fm$data[["idata"]]) Value Df test stat approx F num Df den Df Pr(>F) time1 0.050904 1 0.033005 1.297 1 38 0.2619 Similar Posts:
- Solved – Testing anova hypothesis with contrasts in R and SPSS
- Solved – Analyse unbalanced repeated measures 2x2x2x2 type II anova interactions
- Solved – Justification for default contr.poly() polynomial contrasts in R
- Solved – aov() versus ezANOVA F-values
- Solved – Why are results different when using aov_ez{afex} and Anova{car}, Type III SS in R