資料來源

讀進資料

dta <- read.csv("rt_dummy_data.csv")
head(dta,20)
##    PID   RT   modality   stim
## 1  301 1024 Audio-only   gown
## 2  301  838 Audio-only  might
## 3  301 1060 Audio-only   fern
## 4  301  882 Audio-only   vane
## 5  301  971 Audio-only    pup
## 6  301 1064 Audio-only   rise
## 7  301 1283 Audio-only   hill
## 8  301  982 Audio-only caught
## 9  301 1144 Audio-only    leg
## 10 301  743 Audio-only   hurt
## 11 301  639 Audio-only  shack
## 12 301  853 Audio-only   tell
## 13 301  649 Audio-only  worth
## 14 301  757 Audio-only  shook
## 15 301  958 Audio-only  douse
## 16 301  942 Audio-only   jail
## 17 301  748 Audio-only    bog
## 18 301  994 Audio-only thatch
## 19 301  810 Audio-only    bit
## 20 301 1241 Audio-only    pop

把需要套件載入

library(pacman)
## Warning: 套件 'pacman' 是用 R 版本 4.2.3 來建造的
p_load(lme4,tidyverse,afex)

一般來說,會做 ANOVA

dta_anova <- aggregate(dta$RT,list(dta$PID,dta$modality),mean)
names(dta_anova) <- c('PID','modality','RTmean')
dta_anova$PID <- as.factor(dta_anova$PID)
head(dta_anova,10)
##    PID   modality    RTmean
## 1  301 Audio-only 1027.4106
## 2  302 Audio-only 1047.3455
## 3  303 Audio-only  882.6735
## 4  304 Audio-only 1238.2754
## 5  306 Audio-only 1043.3231
## 6  307 Audio-only 1115.6829
## 7  308 Audio-only 1252.6053
## 8  309 Audio-only  795.7773
## 9  310 Audio-only 1175.5823
## 10 311 Audio-only 1015.0566
rst_anova <- aov(formula = RTmean ~ modality + Error(PID/modality),data=dta_anova)
summary(rst_anova)
## 
## Error: PID
##           Df  Sum Sq Mean Sq F value Pr(>F)
## Residuals 52 2928806   56323               
## 
## Error: PID:modality
##           Df Sum Sq Mean Sq F value   Pr(>F)    
## modality   1 183264  183264   43.72 2.05e-08 ***
## Residuals 52 217966    4192                     
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

看一下兩個情境差異

MNS <- tapply(dta_anova$RTmean,dta_anova$modality,mean)
MNS
##  Audio-only Audiovisual 
##    1044.073    1127.233

-兩組差 83 ms

看一下效果量

p_load(effectsize)
eta_squared(rst_anova, partial = FALSE)
## # Effect Size for ANOVA (Type I)
## 
## Group        | Parameter | Eta2 |       95% CI
## ----------------------------------------------
## PID:modality |  modality | 0.06 | [0.00, 1.00]
## 
## - One-sided CIs: upper bound fixed at [1.00].

以 multilevel 方式分析

完整模型

dta$modality <- ifelse(dta$modality == "Audio-only", 0, 1)
rst <- lmer(RT ~ 1 + modality + (1 + modality|PID) + (1 + modality|stim), 
    data = dta,
    control = lmerControl(optimizer = "bobyqa"))
summary(rst)
## Linear mixed model fit by REML. t-tests use Satterthwaite's method [
## lmerModLmerTest]
## Formula: RT ~ 1 + modality + (1 + modality | PID) + (1 + modality | stim)
##    Data: dta
## Control: lmerControl(optimizer = "bobyqa")
## 
## REML criterion at convergence: 302385.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.3646 -0.6964 -0.0140  0.5886  5.0003 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev. Corr 
##  stim     (Intercept)   303.9   17.43        
##           modality      216.7   14.72   0.16 
##  PID      (Intercept) 28552.7  168.98        
##           modality     7709.8   87.81   -0.17
##  Residual             65258.8  255.46        
## Number of obs: 21679, groups:  stim, 543; PID, 53
## 
## Fixed effects:
##             Estimate Std. Error      df t value Pr(>|t|)    
## (Intercept)  1044.14      23.36   52.14  44.704  < 2e-16 ***
## modality       83.18      12.58   52.09   6.615 2.02e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Correlation of Fixed Effects:
##          (Intr)
## modality -0.178
  • 這邊也是看到 83,顯著