运行 将 R 中的方差分析与两个因子内和一个因子间混合时出错
Error when running mixed anova in R with two within factor and one between factor
我想要一些帮助来使用 R 来分析一个实验，在这个实验中，三组参与者分别被展示了两种类型的刺激物三次。因变量是连续测量。这是数据集的外观示例。
SubjectID Group Trial StimType Measure
1 1 group1 trial3_stimA A 0.55908866
2 2 group1 trial3_stimA A 0.98884446
3 3 group2 trial3_stimA A 0.00000000
4 4 group2 trial3_stimA A 0.27067991
5 5 group3 trial3_stimA A 0.37169285
6 6 group3 trial3_stimA A 0.42113984
7 1 group1 trial3_stimB B 0.00000000
8 2 group1 trial3_stimB B 0.49892807
9 3 group2 trial3_stimB B 0.14602589
10 4 group2 trial3_stimB B 0.50946555
11 5 group3 trial3_stimB B 0.25572820
12 6 group3 trial3_stimB B 0.22932966
13 1 group1 trial1_stimA A 0.42207604
14 2 group1 trial1_stimA A 0.85599588
15 3 group2 trial1_stimA A 0.36428381
16 4 group2 trial1_stimA A 0.46679336
17 5 group3 trial1_stimA A 0.69379734
18 6 group3 trial1_stimA A 0.55607716
19 1 group1 trial1_stimB B 0.24261465
20 2 group1 trial1_stimB B 0.35176384
21 3 group2 trial1_stimB B 0.21116215
22 4 group2 trial1_stimB B 0.33112544
23 5 group3 trial1_stimB B 0.00000000
24 6 group3 trial1_stimB B 0.00000000
25 1 group1 trial2_stimA A 0.05506943
26 2 group1 trial2_stimA A 0.22537470
27 3 group2 trial2_stimA A 0.00000000
28 4 group2 trial2_stimA A 0.18511144
29 5 group3 trial2_stimA A 0.15586156
30 6 group3 trial2_stimA A 0.04467100
31 1 group1 trial2_stimB B 0.03890585
32 2 group1 trial2_stimB B 0.29787709
33 3 group2 trial2_stimB B 0.00000000
34 4 group2 trial2_stimB B 0.28971992
35 5 group3 trial2_stimB B 0.12993238
36 6 group3 trial2_stimB B 0.05066011
这是我的数据结构
'data.frame': 36 obs. of 5 variables:
$ SubjectID: Factor w/ 6 levels "1","2","3","4",..: 1 2 3 4 5 6 1 2 3 4 ...
$ Group : Factor w/ 3 levels "group1","group2",..: 1 1 2 2 3 3 1 1 2 2 ...
$ Trial : Factor w/ 6 levels "trial1_stimA",..: 5 5 5 5 5 5 6 6 6 6 ...
$ StimType : Factor w/ 2 levels "A","B": 1 1 1 1 1 1 2 2 2 2 ...
$ Measure : num 0.559 0.989 0 0.271 0.372 ...
我需要 运行 混合方差分析，其中组作为受试者因素和试验之间的因素，刺激类型作为受试者因素内的因素。我已经尝试使用三个不同的 R 包和不同的语法，但 R 要么 returns 一条错误消息，要么输出缺少交互组 x 试验 x stim 类型。
例如，当我使用 rstatix 包中的 anova_test() 时
#Trying mixed ANOVA with rstatix
>
> mixed.anova <- anova_test(
+ data = prepared_data, dv = Measure, wid = SubjectID,
+ between = Group, within = c(Trial,StimType)
+ )
Error in check.imatrix(X.design) :
Terms in the intra-subject model matrix are not orthogonal.
> get_anova_table(all_subjects)
Error in is.data.frame(x) : object 'all_subjects' not found
当我使用 afex 包中的 aov 时
> #using afex
>
>
> mixed.anova2 <- aov_car(Measure ~ Group*Trial*StimType + Error(1|SubjectID/(Trial*StimType)), prepared_data)
Error: Empty cells in within-subjects design (i.e., bad data structure).
table(data[c("Trial", "StimType")])
# StimType
# Trial A B
# trial1_stimA 6 0
# trial1_stimB 0 6
# trial2_stimA 6 0
# trial2_stimB 0 6
# trial3_stimA 6 0
# trial3_stimB 0 6
>
> aov.bww
Call:
aov(formula = SCR ~ Group * Trial * CSType + Error(SubjectID) +
Group, data = sixPhasesAbs2)
Grand Mean: 397.1325
Stratum 1: SubjectID
Terms:
Group Residuals
Sum of Squares 187283464 6399838881
Deg. of Freedom 2 81
Residual standard error: 8888.777
18 out of 20 effects not estimable
Estimated effects may be unbalanced
Stratum 2: Within
Terms:
Trial Group:Trial Residuals
Sum of Squares 158766098 374583941 12799639740
Deg. of Freedom 5 10 405
Residual standard error: 5621.748
12 out of 27 effects not estimable
Estimated effects may be unbalanced
最后，我尝试使用 lmer
> #using lmer
> model = lmer(Measure ~Group*Trial*StimType +(1|SubjectID), data = prepared_data)
fixed-effect model matrix is rank deficient so dropping 18 columns / coefficients
>
> summary(model)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Measure ~ Group * Trial * StimType + (1 | SubjectID)
Data: prepared_data
REML criterion at convergence: -16.8
Scaled residuals:
Min 1Q Median 3Q Max
-1.1854 -0.5424 0.0000 0.5424 1.1854
Random effects:
Groups Name Variance Std.Dev.
SubjectID (Intercept) 0.024226 0.15565
Residual 0.006861 0.08283
Number of obs: 36, groups: SubjectID, 6
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.639036 0.124672 4.459241 5.126 0.005095 **
Groupgroup2 -0.223497 0.176313 4.459241 -1.268 0.267088
Groupgroup3 -0.014099 0.176313 4.459241 -0.080 0.939728
Trialtrial1_stimB -0.341847 0.082829 15.000000 -4.127 0.000896 ***
Trialtrial2_stimA -0.498814 0.082829 15.000000 -6.022 2.34e-05 ***
Trialtrial2_stimB -0.470644 0.082829 15.000000 -5.682 4.35e-05 ***
Trialtrial3_stimA 0.134931 0.082829 15.000000 1.629 0.124124
Trialtrial3_stimB -0.389572 0.082829 15.000000 -4.703 0.000283 ***
Groupgroup2:Trialtrial1_stimB 0.197452 0.117138 15.000000 1.686 0.112550
Groupgroup3:Trialtrial1_stimB -0.283091 0.117138 15.000000 -2.417 0.028865 *
Groupgroup2:Trialtrial2_stimA 0.175831 0.117138 15.000000 1.501 0.154095
Groupgroup3:Trialtrial2_stimA -0.025857 0.117138 15.000000 -0.221 0.828271
Groupgroup2:Trialtrial2_stimB 0.199966 0.117138 15.000000 1.707 0.108413
Groupgroup3:Trialtrial2_stimB -0.063997 0.117138 15.000000 -0.546 0.592870
Groupgroup2:Trialtrial3_stimA -0.415129 0.117138 15.000000 -3.544 0.002946 **
Groupgroup3:Trialtrial3_stimA -0.363452 0.117138 15.000000 -3.103 0.007276 **
Groupgroup2:Trialtrial3_stimB 0.301779 0.117138 15.000000 2.576 0.021070 *
Groupgroup3:Trialtrial3_stimB 0.007164 0.117138 15.000000 0.061 0.952043
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation matrix not shown by default, as p = 18 > 12.
Use print(x, correlation=TRUE) or
vcov(x) if you need it
fit warnings:
fixed-effect model matrix is rank deficient so dropping 18 columns / coefficients
>
> anova(model)
Missing cells for: Trialtrial1_stimB:StimTypeA, Trialtrial2_stimB:StimTypeA, Trialtrial3_stimB:StimTypeA, Trialtrial1_stimA:StimTypeB, Trialtrial2_stimA:StimTypeB, Trialtrial3_stimA:StimTypeB, Groupgroup1:Trialtrial1_stimB:StimTypeA, Groupgroup2:Trialtrial1_stimB:StimTypeA, Groupgroup3:Trialtrial1_stimB:StimTypeA, Groupgroup1:Trialtrial2_stimB:StimTypeA, Groupgroup2:Trialtrial2_stimB:StimTypeA, Groupgroup3:Trialtrial2_stimB:StimTypeA, Groupgroup1:Trialtrial3_stimB:StimTypeA, Groupgroup2:Trialtrial3_stimB:StimTypeA, Groupgroup3:Trialtrial3_stimB:StimTypeA, Groupgroup1:Trialtrial1_stimA:StimTypeB, Groupgroup2:Trialtrial1_stimA:StimTypeB, Groupgroup3:Trialtrial1_stimA:StimTypeB, Groupgroup1:Trialtrial2_stimA:StimTypeB, Groupgroup2:Trialtrial2_stimA:StimTypeB, Groupgroup3:Trialtrial2_stimA:StimTypeB, Groupgroup1:Trialtrial3_stimA:StimTypeB, Groupgroup2:Trialtrial3_stimA:StimTypeB, Groupgroup3:Trialtrial3_stimA:StimTypeB.
Interpret type III hypotheses with care.
Type III Analysis of Variance Table with Satterthwaite's method
Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
Group 0.00723 0.003614 2 3 0.5267 0.6367151
Trial 0.96171 0.192343 5 15 28.0356 4.172e-07 ***
Group:Trial 0.44102 0.044102 10 15 6.4283 0.0007411 ***
StimType
Group:StimType
Trial:StimType
Group:Trial:StimType
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
我搜索过类似的问题，但没有找到任何可以帮助我解决这个问题的答案。我怀疑（鉴于错误消息）我可能必须更改数据集的结构，但我不知道该怎么做，因为我是 R 的初学者。我如何 运行 混合方差分析数据集？
由于 Trial 包含有关 Trial 和 StimType 的信息，因此固定效应模型矩阵秩亏。从原始 post 中的 afex::aov_car() 输出可以看出这一点，它说明了空单元格是 Trial 的效果，其中包含来自两个不同变量的信息。
> mixed.anova2 <- aov_car(Measure ~ Group*Trial*StimType + Error(1|SubjectID/(Trial*StimType)), prepared_data)
Error: Empty cells in within-subjects design (i.e., bad data structure).
table(data[c("Trial", "StimType")])
# StimType
# Trial A B
# trial1_stimA 6 0
# trial1_stimB 0 6
# trial2_stimA 6 0
# trial2_stimB 0 6
# trial3_stimA 6 0
# trial3_stimB 0 6
>
> aov.bww
我们可以通过将 Trial 变量拆分为两个变量来纠正排名不足。
使用 lme4 包和 tidyr::separate() 将试验值与刺激值分开，分析如下所示：
textFile <- "rowId SubjectID Group Trial StimType Measure
1 1 group1 trial3_stimA A 0.55908866
2 2 group1 trial3_stimA A 0.98884446
3 3 group2 trial3_stimA A 0.00000000
4 4 group2 trial3_stimA A 0.27067991
5 5 group3 trial3_stimA A 0.37169285
6 6 group3 trial3_stimA A 0.42113984
7 1 group1 trial3_stimB B 0.00000000
8 2 group1 trial3_stimB B 0.49892807
9 3 group2 trial3_stimB B 0.14602589
10 4 group2 trial3_stimB B 0.50946555
11 5 group3 trial3_stimB B 0.25572820
12 6 group3 trial3_stimB B 0.22932966
13 1 group1 trial1_stimA A 0.42207604
14 2 group1 trial1_stimA A 0.85599588
15 3 group2 trial1_stimA A 0.36428381
16 4 group2 trial1_stimA A 0.46679336
17 5 group3 trial1_stimA A 0.69379734
18 6 group3 trial1_stimA A 0.55607716
19 1 group1 trial1_stimB B 0.24261465
20 2 group1 trial1_stimB B 0.35176384
21 3 group2 trial1_stimB B 0.21116215
22 4 group2 trial1_stimB B 0.33112544
23 5 group3 trial1_stimB B 0.00000000
24 6 group3 trial1_stimB B 0.00000000
25 1 group1 trial2_stimA A 0.05506943
26 2 group1 trial2_stimA A 0.22537470
27 3 group2 trial2_stimA A 0.00000000
28 4 group2 trial2_stimA A 0.18511144
29 5 group3 trial2_stimA A 0.15586156
30 6 group3 trial2_stimA A 0.04467100
31 1 group1 trial2_stimB B 0.03890585
32 2 group1 trial2_stimB B 0.29787709
33 3 group2 trial2_stimB B 0.00000000
34 4 group2 trial2_stimB B 0.28971992
35 5 group3 trial2_stimB B 0.12993238
36 6 group3 trial2_stimB B 0.05066011
"
data <- read.table(text = textFile,header = TRUE)
library(dplyr)
library(tidyr)
# split trial from Stim
data %>% separate(Trial,into = c("trialName","stimName")) -> data
library(lme4)
model = lmer(Measure ~Group*trialName*StimType +(1|SubjectID), data = data)
summary(model)
...输出：
> summary(model)
Linear mixed model fit by REML ['lmerMod']
Formula: Measure ~ Group * trialName * StimType + (1 | SubjectID)
Data: data
REML criterion at convergence: -16.8
Scaled residuals:
Min 1Q Median 3Q Max
-1.1854 -0.5424 0.0000 0.5424 1.1854
Random effects:
Groups Name Variance Std.Dev.
SubjectID (Intercept) 0.024226 0.15565
Residual 0.006861 0.08283
Number of obs: 36, groups: SubjectID, 6
Fixed effects:
Estimate Std. Error t value
(Intercept) 0.63904 0.12467 5.126
Groupgroup2 -0.22350 0.17631 -1.268
Groupgroup3 -0.01410 0.17631 -0.080
trialNametrial2 -0.49881 0.08283 -6.022
trialNametrial3 0.13493 0.08283 1.629
StimTypeB -0.34185 0.08283 -4.127
Groupgroup2:trialNametrial2 0.17583 0.11714 1.501
Groupgroup3:trialNametrial2 -0.02586 0.11714 -0.221
Groupgroup2:trialNametrial3 -0.41513 0.11714 -3.544
Groupgroup3:trialNametrial3 -0.36345 0.11714 -3.103
Groupgroup2:StimTypeB 0.19745 0.11714 1.686
Groupgroup3:StimTypeB -0.28309 0.11714 -2.417
trialNametrial2:StimTypeB 0.37002 0.11714 3.159
trialNametrial3:StimTypeB -0.18266 0.11714 -1.559
Groupgroup2:trialNametrial2:StimTypeB -0.17332 0.16566 -1.046
Groupgroup3:trialNametrial2:StimTypeB 0.24495 0.16566 1.479
Groupgroup2:trialNametrial3:StimTypeB 0.51946 0.16566 3.136
Groupgroup3:trialNametrial3:StimTypeB 0.65371 0.16566 3.946
Correlation matrix not shown by default, as p = 18 > 12.
Use print(x, correlation=TRUE) or
vcov(x) if you need it
>
更新：解析效果名称
在对我的回答的评论中，发问者解释说，用于描述试验和刺激类型的数据与 post 随问题一起提供的数据中表示的数据不同。鉴于评论的内容，这里有一个 dplyr 的方法来解析试验条件和刺激类型。
# solving the data cleaning problem in comments
df <- data.frame(treatmentName = c("MeanCondPlusLate",
"MeanCondMinusLate",
"MeanExtPlusLate",
"MeanExtMinusLate",
"TestPlus1",
"TestMinus1"))
library(dplyr)
df %>% mutate(trial = case_when(grepl("Cond",treatmentName) == TRUE ~ 1,
grepl("Ext",treatmentName) == TRUE ~ 2,
grepl("Test",treatmentName) == TRUE ~ 3),
stimLevel = case_when(grepl("Plus",treatmentName) == TRUE ~ "A",
grepl("Minus",treatmentName) == TRUE ~ "B"))
...输出：
treatmentName trial stimLevel
1 MeanCondPlusLate 1 A
2 MeanCondMinusLate 1 B
3 MeanExtPlusLate 2 A
4 MeanExtMinusLate 2 B
5 TestPlus1 3 A
6 TestMinus1 3 B
>
我想要一些帮助来使用 R 来分析一个实验，在这个实验中，三组参与者分别被展示了两种类型的刺激物三次。因变量是连续测量。这是数据集的外观示例。
SubjectID Group Trial StimType Measure
1 1 group1 trial3_stimA A 0.55908866
2 2 group1 trial3_stimA A 0.98884446
3 3 group2 trial3_stimA A 0.00000000
4 4 group2 trial3_stimA A 0.27067991
5 5 group3 trial3_stimA A 0.37169285
6 6 group3 trial3_stimA A 0.42113984
7 1 group1 trial3_stimB B 0.00000000
8 2 group1 trial3_stimB B 0.49892807
9 3 group2 trial3_stimB B 0.14602589
10 4 group2 trial3_stimB B 0.50946555
11 5 group3 trial3_stimB B 0.25572820
12 6 group3 trial3_stimB B 0.22932966
13 1 group1 trial1_stimA A 0.42207604
14 2 group1 trial1_stimA A 0.85599588
15 3 group2 trial1_stimA A 0.36428381
16 4 group2 trial1_stimA A 0.46679336
17 5 group3 trial1_stimA A 0.69379734
18 6 group3 trial1_stimA A 0.55607716
19 1 group1 trial1_stimB B 0.24261465
20 2 group1 trial1_stimB B 0.35176384
21 3 group2 trial1_stimB B 0.21116215
22 4 group2 trial1_stimB B 0.33112544
23 5 group3 trial1_stimB B 0.00000000
24 6 group3 trial1_stimB B 0.00000000
25 1 group1 trial2_stimA A 0.05506943
26 2 group1 trial2_stimA A 0.22537470
27 3 group2 trial2_stimA A 0.00000000
28 4 group2 trial2_stimA A 0.18511144
29 5 group3 trial2_stimA A 0.15586156
30 6 group3 trial2_stimA A 0.04467100
31 1 group1 trial2_stimB B 0.03890585
32 2 group1 trial2_stimB B 0.29787709
33 3 group2 trial2_stimB B 0.00000000
34 4 group2 trial2_stimB B 0.28971992
35 5 group3 trial2_stimB B 0.12993238
36 6 group3 trial2_stimB B 0.05066011
这是我的数据结构
'data.frame': 36 obs. of 5 variables:
$ SubjectID: Factor w/ 6 levels "1","2","3","4",..: 1 2 3 4 5 6 1 2 3 4 ...
$ Group : Factor w/ 3 levels "group1","group2",..: 1 1 2 2 3 3 1 1 2 2 ...
$ Trial : Factor w/ 6 levels "trial1_stimA",..: 5 5 5 5 5 5 6 6 6 6 ...
$ StimType : Factor w/ 2 levels "A","B": 1 1 1 1 1 1 2 2 2 2 ...
$ Measure : num 0.559 0.989 0 0.271 0.372 ...
我需要 运行 混合方差分析，其中组作为受试者因素和试验之间的因素，刺激类型作为受试者因素内的因素。我已经尝试使用三个不同的 R 包和不同的语法，但 R 要么 returns 一条错误消息，要么输出缺少交互组 x 试验 x stim 类型。
例如，当我使用 rstatix 包中的 anova_test() 时
#Trying mixed ANOVA with rstatix
>
> mixed.anova <- anova_test(
+ data = prepared_data, dv = Measure, wid = SubjectID,
+ between = Group, within = c(Trial,StimType)
+ )
Error in check.imatrix(X.design) :
Terms in the intra-subject model matrix are not orthogonal.
> get_anova_table(all_subjects)
Error in is.data.frame(x) : object 'all_subjects' not found
当我使用 afex 包中的 aov 时
> #using afex
>
>
> mixed.anova2 <- aov_car(Measure ~ Group*Trial*StimType + Error(1|SubjectID/(Trial*StimType)), prepared_data)
Error: Empty cells in within-subjects design (i.e., bad data structure).
table(data[c("Trial", "StimType")])
# StimType
# Trial A B
# trial1_stimA 6 0
# trial1_stimB 0 6
# trial2_stimA 6 0
# trial2_stimB 0 6
# trial3_stimA 6 0
# trial3_stimB 0 6
>
> aov.bww
Call:
aov(formula = SCR ~ Group * Trial * CSType + Error(SubjectID) +
Group, data = sixPhasesAbs2)
Grand Mean: 397.1325
Stratum 1: SubjectID
Terms:
Group Residuals
Sum of Squares 187283464 6399838881
Deg. of Freedom 2 81
Residual standard error: 8888.777
18 out of 20 effects not estimable
Estimated effects may be unbalanced
Stratum 2: Within
Terms:
Trial Group:Trial Residuals
Sum of Squares 158766098 374583941 12799639740
Deg. of Freedom 5 10 405
Residual standard error: 5621.748
12 out of 27 effects not estimable
Estimated effects may be unbalanced
最后，我尝试使用 lmer
> #using lmer
> model = lmer(Measure ~Group*Trial*StimType +(1|SubjectID), data = prepared_data)
fixed-effect model matrix is rank deficient so dropping 18 columns / coefficients
>
> summary(model)
Linear mixed model fit by REML. t-tests use Satterthwaite's method ['lmerModLmerTest']
Formula: Measure ~ Group * Trial * StimType + (1 | SubjectID)
Data: prepared_data
REML criterion at convergence: -16.8
Scaled residuals:
Min 1Q Median 3Q Max
-1.1854 -0.5424 0.0000 0.5424 1.1854
Random effects:
Groups Name Variance Std.Dev.
SubjectID (Intercept) 0.024226 0.15565
Residual 0.006861 0.08283
Number of obs: 36, groups: SubjectID, 6
Fixed effects:
Estimate Std. Error df t value Pr(>|t|)
(Intercept) 0.639036 0.124672 4.459241 5.126 0.005095 **
Groupgroup2 -0.223497 0.176313 4.459241 -1.268 0.267088
Groupgroup3 -0.014099 0.176313 4.459241 -0.080 0.939728
Trialtrial1_stimB -0.341847 0.082829 15.000000 -4.127 0.000896 ***
Trialtrial2_stimA -0.498814 0.082829 15.000000 -6.022 2.34e-05 ***
Trialtrial2_stimB -0.470644 0.082829 15.000000 -5.682 4.35e-05 ***
Trialtrial3_stimA 0.134931 0.082829 15.000000 1.629 0.124124
Trialtrial3_stimB -0.389572 0.082829 15.000000 -4.703 0.000283 ***
Groupgroup2:Trialtrial1_stimB 0.197452 0.117138 15.000000 1.686 0.112550
Groupgroup3:Trialtrial1_stimB -0.283091 0.117138 15.000000 -2.417 0.028865 *
Groupgroup2:Trialtrial2_stimA 0.175831 0.117138 15.000000 1.501 0.154095
Groupgroup3:Trialtrial2_stimA -0.025857 0.117138 15.000000 -0.221 0.828271
Groupgroup2:Trialtrial2_stimB 0.199966 0.117138 15.000000 1.707 0.108413
Groupgroup3:Trialtrial2_stimB -0.063997 0.117138 15.000000 -0.546 0.592870
Groupgroup2:Trialtrial3_stimA -0.415129 0.117138 15.000000 -3.544 0.002946 **
Groupgroup3:Trialtrial3_stimA -0.363452 0.117138 15.000000 -3.103 0.007276 **
Groupgroup2:Trialtrial3_stimB 0.301779 0.117138 15.000000 2.576 0.021070 *
Groupgroup3:Trialtrial3_stimB 0.007164 0.117138 15.000000 0.061 0.952043
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
Correlation matrix not shown by default, as p = 18 > 12.
Use print(x, correlation=TRUE) or
vcov(x) if you need it
fit warnings:
fixed-effect model matrix is rank deficient so dropping 18 columns / coefficients
>
> anova(model)
Missing cells for: Trialtrial1_stimB:StimTypeA, Trialtrial2_stimB:StimTypeA, Trialtrial3_stimB:StimTypeA, Trialtrial1_stimA:StimTypeB, Trialtrial2_stimA:StimTypeB, Trialtrial3_stimA:StimTypeB, Groupgroup1:Trialtrial1_stimB:StimTypeA, Groupgroup2:Trialtrial1_stimB:StimTypeA, Groupgroup3:Trialtrial1_stimB:StimTypeA, Groupgroup1:Trialtrial2_stimB:StimTypeA, Groupgroup2:Trialtrial2_stimB:StimTypeA, Groupgroup3:Trialtrial2_stimB:StimTypeA, Groupgroup1:Trialtrial3_stimB:StimTypeA, Groupgroup2:Trialtrial3_stimB:StimTypeA, Groupgroup3:Trialtrial3_stimB:StimTypeA, Groupgroup1:Trialtrial1_stimA:StimTypeB, Groupgroup2:Trialtrial1_stimA:StimTypeB, Groupgroup3:Trialtrial1_stimA:StimTypeB, Groupgroup1:Trialtrial2_stimA:StimTypeB, Groupgroup2:Trialtrial2_stimA:StimTypeB, Groupgroup3:Trialtrial2_stimA:StimTypeB, Groupgroup1:Trialtrial3_stimA:StimTypeB, Groupgroup2:Trialtrial3_stimA:StimTypeB, Groupgroup3:Trialtrial3_stimA:StimTypeB.
Interpret type III hypotheses with care.
Type III Analysis of Variance Table with Satterthwaite's method
Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
Group 0.00723 0.003614 2 3 0.5267 0.6367151
Trial 0.96171 0.192343 5 15 28.0356 4.172e-07 ***
Group:Trial 0.44102 0.044102 10 15 6.4283 0.0007411 ***
StimType
Group:StimType
Trial:StimType
Group:Trial:StimType
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
我搜索过类似的问题，但没有找到任何可以帮助我解决这个问题的答案。我怀疑（鉴于错误消息）我可能必须更改数据集的结构，但我不知道该怎么做，因为我是 R 的初学者。我如何 运行 混合方差分析数据集？
由于 Trial 包含有关 Trial 和 StimType 的信息，因此固定效应模型矩阵秩亏。从原始 post 中的 afex::aov_car() 输出可以看出这一点，它说明了空单元格是 Trial 的效果，其中包含来自两个不同变量的信息。
> mixed.anova2 <- aov_car(Measure ~ Group*Trial*StimType + Error(1|SubjectID/(Trial*StimType)), prepared_data)
Error: Empty cells in within-subjects design (i.e., bad data structure).
table(data[c("Trial", "StimType")])
# StimType
# Trial A B
# trial1_stimA 6 0
# trial1_stimB 0 6
# trial2_stimA 6 0
# trial2_stimB 0 6
# trial3_stimA 6 0
# trial3_stimB 0 6
>
> aov.bww
我们可以通过将 Trial 变量拆分为两个变量来纠正排名不足。
使用 lme4 包和 tidyr::separate() 将试验值与刺激值分开，分析如下所示：
textFile <- "rowId SubjectID Group Trial StimType Measure
1 1 group1 trial3_stimA A 0.55908866
2 2 group1 trial3_stimA A 0.98884446
3 3 group2 trial3_stimA A 0.00000000
4 4 group2 trial3_stimA A 0.27067991
5 5 group3 trial3_stimA A 0.37169285
6 6 group3 trial3_stimA A 0.42113984
7 1 group1 trial3_stimB B 0.00000000
8 2 group1 trial3_stimB B 0.49892807
9 3 group2 trial3_stimB B 0.14602589
10 4 group2 trial3_stimB B 0.50946555
11 5 group3 trial3_stimB B 0.25572820
12 6 group3 trial3_stimB B 0.22932966
13 1 group1 trial1_stimA A 0.42207604
14 2 group1 trial1_stimA A 0.85599588
15 3 group2 trial1_stimA A 0.36428381
16 4 group2 trial1_stimA A 0.46679336
17 5 group3 trial1_stimA A 0.69379734
18 6 group3 trial1_stimA A 0.55607716
19 1 group1 trial1_stimB B 0.24261465
20 2 group1 trial1_stimB B 0.35176384
21 3 group2 trial1_stimB B 0.21116215
22 4 group2 trial1_stimB B 0.33112544
23 5 group3 trial1_stimB B 0.00000000
24 6 group3 trial1_stimB B 0.00000000
25 1 group1 trial2_stimA A 0.05506943
26 2 group1 trial2_stimA A 0.22537470
27 3 group2 trial2_stimA A 0.00000000
28 4 group2 trial2_stimA A 0.18511144
29 5 group3 trial2_stimA A 0.15586156
30 6 group3 trial2_stimA A 0.04467100
31 1 group1 trial2_stimB B 0.03890585
32 2 group1 trial2_stimB B 0.29787709
33 3 group2 trial2_stimB B 0.00000000
34 4 group2 trial2_stimB B 0.28971992
35 5 group3 trial2_stimB B 0.12993238
36 6 group3 trial2_stimB B 0.05066011
"
data <- read.table(text = textFile,header = TRUE)
library(dplyr)
library(tidyr)
# split trial from Stim
data %>% separate(Trial,into = c("trialName","stimName")) -> data
library(lme4)
model = lmer(Measure ~Group*trialName*StimType +(1|SubjectID), data = data)
summary(model)
...输出：
> summary(model)
Linear mixed model fit by REML ['lmerMod']
Formula: Measure ~ Group * trialName * StimType + (1 | SubjectID)
Data: data
REML criterion at convergence: -16.8
Scaled residuals:
Min 1Q Median 3Q Max
-1.1854 -0.5424 0.0000 0.5424 1.1854
Random effects:
Groups Name Variance Std.Dev.
SubjectID (Intercept) 0.024226 0.15565
Residual 0.006861 0.08283
Number of obs: 36, groups: SubjectID, 6
Fixed effects:
Estimate Std. Error t value
(Intercept) 0.63904 0.12467 5.126
Groupgroup2 -0.22350 0.17631 -1.268
Groupgroup3 -0.01410 0.17631 -0.080
trialNametrial2 -0.49881 0.08283 -6.022
trialNametrial3 0.13493 0.08283 1.629
StimTypeB -0.34185 0.08283 -4.127
Groupgroup2:trialNametrial2 0.17583 0.11714 1.501
Groupgroup3:trialNametrial2 -0.02586 0.11714 -0.221
Groupgroup2:trialNametrial3 -0.41513 0.11714 -3.544
Groupgroup3:trialNametrial3 -0.36345 0.11714 -3.103
Groupgroup2:StimTypeB 0.19745 0.11714 1.686
Groupgroup3:StimTypeB -0.28309 0.11714 -2.417
trialNametrial2:StimTypeB 0.37002 0.11714 3.159
trialNametrial3:StimTypeB -0.18266 0.11714 -1.559
Groupgroup2:trialNametrial2:StimTypeB -0.17332 0.16566 -1.046
Groupgroup3:trialNametrial2:StimTypeB 0.24495 0.16566 1.479
Groupgroup2:trialNametrial3:StimTypeB 0.51946 0.16566 3.136
Groupgroup3:trialNametrial3:StimTypeB 0.65371 0.16566 3.946
Correlation matrix not shown by default, as p = 18 > 12.
Use print(x, correlation=TRUE) or
vcov(x) if you need it
>
更新：解析效果名称
在对我的回答的评论中，发问者解释说，用于描述试验和刺激类型的数据与 post 随问题一起提供的数据中表示的数据不同。鉴于评论的内容，这里有一个 dplyr 的方法来解析试验条件和刺激类型。
# solving the data cleaning problem in comments
df <- data.frame(treatmentName = c("MeanCondPlusLate",
"MeanCondMinusLate",
"MeanExtPlusLate",
"MeanExtMinusLate",
"TestPlus1",
"TestMinus1"))
library(dplyr)
df %>% mutate(trial = case_when(grepl("Cond",treatmentName) == TRUE ~ 1,
grepl("Ext",treatmentName) == TRUE ~ 2,
grepl("Test",treatmentName) == TRUE ~ 3),
stimLevel = case_when(grepl("Plus",treatmentName) == TRUE ~ "A",
grepl("Minus",treatmentName) == TRUE ~ "B"))
...输出：
treatmentName trial stimLevel
1 MeanCondPlusLate 1 A
2 MeanCondMinusLate 1 B
3 MeanExtPlusLate 2 A
4 MeanExtMinusLate 2 B
5 TestPlus1 3 A
6 TestMinus1 3 B
>