Disclamer:我不确定在哪里发布这个问题:CV 或 SO,但最终决定先在这里尝试
一位审稿人要求我为我论文的辅助分析中报告和,这是我所在领域的标准) 。
具体来说,我提到了调用的固定A:B:C:response
效果之间的重要 4 向交互(anova(lmer.model)
emmeans
pairs()
我知道至少存在关于计算混合效果模型的效果大小是否有意义的争论,但我想在没有任何进一步讨论的情况下让审稿人满意。
我遇到了一些在 LMEM 的上下文中提到效果大小的消息来源,但我觉得数学不够强大,无法理解它。
我的问题是双重的:
如果有任何可引用的方式来为 anovas测试、或 我该怎么做?一个包/函数/脚本会很有帮助。最终提示从对象中检索关键值以手动获取 eff.sizes 也会有所帮助
R
lmer()
例如。我知道并 给出但没有我不知道如何获得它
anova(lmer_obj)
R
- 如果没有办法计算它 - 在回答审稿人为什么我坚持跳过效果大小时,我可以引用什么论文(不是博客文章)。
请注意,我对随机效应本身不感兴趣——我定义了一个由随机效应的设计结构证明的最大值,只是为了控制更多的误差方差,我唯一感兴趣的是固定效应——即来自 anova 表的 Fs 和相应的边际均值。
一些示例性结果(对应于具有随机值的真实数据结构)如下:
A:B:C:response
交互作用在 F(1; 12082,1)=4,60 时显着 )
> anova(m0)
Type III Analysis of Variance Table with Satterthwaite's method
Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
A 1.6717 1.6717 1 134.7 0.4210 0.51755
B 0.1375 0.1375 1 11860.0 0.0346 0.85238
C 7.1708 7.1708 1 133.5 1.8058 0.18129
response 3.7775 3.7775 1 12070.3 0.9513 0.32940
A:B 1.1291 1.1291 1 11872.8 0.2844 0.59387
A:C 3.2427 3.2427 1 121.1 0.8166 0.36796
B:C 4.1048 4.1048 1 121.8 1.0337 0.31130
A:response 0.0000 0.0000 1 12080.8 0.0000 0.99828
B:response 4.2350 4.2350 1 12078.9 1.0665 0.30176
C:response 1.6567 1.6567 1 12082.3 0.4172 0.51834
A:B:C 8.5249 8.5249 1 131.4 2.1469 0.14525
A:B:response 0.8765 0.8765 1 132.0 0.2207 0.63926
A:C:response 1.5150 1.5150 1 119.5 0.3815 0.53797
B:C:response 0.5921 0.5921 1 122.5 0.1491 0.70005
A:B:C:response 18.2803 18.2803 1 12082.1 4.6036 0.03192 *
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
后续比较表明 的交互作用A:B:C
仅对response=0
: )
> joint_tests(m0, by="response")
Note: D.f. calculations have been disabled because the number of observations exceeds 3000.
To enable adjustments, set emm_options(pbkrtest.limit = 12292) or larger,
but be warned that this may result in large computation time and memory use.
response = no:
model term df1 df2 F.ratio p.value
A 1 858.01 0.161 0.6884
B 1 12177.54 0.530 0.4664
C 1 837.27 0.212 0.6453
A:B 1 785.54 0.365 0.5457
A:C 1 507.23 0.872 0.3509
B:C 1 476.41 0.145 0.7035
A:B:C 1 809.32 4.680 0.0308
response = yes:
model term df1 df2 F.ratio p.value
A 1 184.95 0.361 0.5489
B 1 11853.63 0.601 0.4381
C 1 182.47 3.237 0.0736
A:B 1 177.43 0.000 0.9874
A:C 1 130.20 0.072 0.7882
B:C 1 123.71 1.465 0.2285
A:B:C 1 178.96 0.235 0.6281
在我的真实数据集中,有更多的预测效果,但想法仍然相同。我在APA 6 ed s 或 s 。规则,我被要求为它们添加某种效果大小。
任何有关此事的帮助将不胜感激。
完全可重现的例子
###############################################################
#Simulate data replicating the structure of a real life example
###############################################################
library(AlgDesign) #for generating all levels a factorial design)
set.seed(2)
df <-gen.factorial(c(16,2,2,2,100), factors = "all",
varNames = c("rep", "A", "B", "C", "Subject"))
df$rep <- as.numeric(df$rep)
df$Subject <- as.numeric(df$Subject)
logRT <- rnorm(n=12800, m=7, sd=2)
trialno<- rep(1:128, times = 100)
response <- factor(sample(0:1, 12800, prob = c(0.3, 0.7), replace = T),
labels= c("no", "yes"))
#Simulate some values to drop (set as missings) due to extremly low latency
missingRTs<- as.logical(sample(0:1, 12800, prob = c(0.96, 0.04), replace = T))
logRT[missingRTs==1] <- NA
df <- cbind(df, logRT, trialno, response)
df <- df[complete.cases(df),]
##########################
#Setup model with afex
##########################
library(afex)
m0 <- mixed(logRT ~ A*B*C*response + (A*B*C*response||Subject),
data = df, return = 'merMod', method = "S", expand_re = TRUE)
##########################
#Get results for paper
##########################
anova(m0)
emm_options(lmerTest.limit = 12292)
em0<-emmeans(m0, ~A*B*C*response)
joint_tests(m0, by="response")
joint_tests(m0, by=c("response", "B"))