Power analysis of HLM LMM
1.0.0
심리학의 재현성이 점점 더 중요 해짐에 따라 통계 분석에 대한 전력 계산에 대한 수요도 증가하고 있습니다. 대부분의 기본 통계 방법 ( t test , ANOVA 등)의 전력 계산은 일부 소프트웨어 (예 : jamovi )에서 완료 될 수 있습니다. 그러나 고급 통계 방법으로서, 선형 혼합 모델 (LMM)의 전력 계산은 시각적 인터페이스를 갖춘 jamovi 와 같은 소프트웨어에서 수행 할 수 없으며 R 의 관련 데이터 패킷 및 함수는 여전히 필요합니다. 이 기사는 주로 데이터 소개 및 기능 소개, 적용 가능 및 제한된 시나리오 및 특정 예제 데모를 포함하여 R 에서 LMM의 전력 계산을 수행하는 방법을 소개합니다.
simr 데이터 패킷은 R 에서 LMM 전력 분석에 더 편리한 패키지입니다. 분석에서 모델 데이터는 Monte Carlo 시뮬레이션에 의해 시뮬레이션되며 상당한 비율은 설정된 시뮬레이션 수에서 계산됩니다. 주요 함수는 powerSim() 입니다.
사후 효과의 힘, 즉 기존 데이터를 기반으로 구축 된 모델의 효과를 기반으로 한 전력을 계산합니다.
선험적 효과의 힘, 즉 기존 연구에서 효과의 크기에 따라 힘을 계산합니다.
전력 규모는 피험자 또는 항목의 수에 의해 영향을 받기 때문에 피험자/항목의 수를 추정하십시오. 특정 모델은 숫자가 많을수록 전력이 클수록, 실험 전의 효과에 기초하여 공식적인 실험에 필요한 피험자/항목의 수를 추정하는 기능을 제공한다는 것입니다.
여기서 우리는 심리학 실험 (데이터 다운로드 주소)의 데이터를 사용하여 입증합니다. 이 실험은인지 과제 ( DV )를 수행 할 때 10 명의 피험자 ( subj )의 반응을 수집했습니다. 실험에는 두 가지 조건 ( CondA 및 CondB )이 있었는데, 이는 2*2 피험자와 프로젝트 내에서 설계되었습니다 ( item ). 즉, 두 가지 임의의 요인이 동시에 고려됩니다.
우리는 여기서 고정 효과 및 비 랜덤 효과에 중점을두기 때문에 모델링시 무작위 효과 부분에서 랜덤 절편 만 고려하고 고정 효과 부분은 두 가지 주요 효과와 그 상호 작용을 검사합니다.
모델링은 다음과 같습니다.
> Model = lmer(data = DemoData, DV~CondA*CondB+(1|subj)+(1|item))
모델의 고정 효과 부분은 다음과 같습니다.
> summary(Model)$coef
# Estimate Std. Error df t value Pr(>|t|)
# (Intercept) 279.43090 23.37537 11.71977 11.9540740 6.422210e-08
# CondAA2 24.29565 13.49694 498.43142 1.8000866 7.245150e-02
# CondBB2 12.18484 13.36445 509.31328 0.9117351 3.623395e-01
# CondAA2:CondBB2 -32.82881 19.29629 509.76684 -1.7013020 8.949610e-02
여기서 우리는 예를 들어 두 사람 ( CondAA2:CondBB2 ) 사이의 상호 작용을 검사하는 고정 효과를 취합니다.
> PowerAB_ttest = simr::powerSim(fit = Model, # 要考察的模型
test = fixed('CondAA2:CondBB2', # 要考察的固定效应的名称
method = 't'), # 选取检验方法，因为固定效应为t检验，因此method设置为t
nsim=50) # 设置模拟次数，建议设置为500 (此时可以获取到较稳定的power)
> PowerAB_ttest
# Power for predictor 'CondAA2:CondBB2', (95% confidence interval):
# 46.00% (31.81, 60.68)
#
# Test: t-test with Satterthwaite degrees of freedom (package lmerTest)
# Effect size for CondAA2:CondBB2 is -33.
#
# Based on 50 simulations, (2 warnings, 0 errors)
# alpha = 0.05, nrow = 553
#
# Time elapsed: 0 h 0 m 8 s
#
# nb: result might be an observed power calculation
상호 작용의 고정 효과의 힘은 0.46 이었다 (신뢰 구간 0.318 – 0.607 ).
마찬가지로 CondA 의 고정 효과를 계산하는 전력은 다음 코드를 통해 구현할 수 있습니다.
> PowerA_ttest = simr::powerSim(fit = Model, test = fixed('CondAA2', method = 't'), nsim=50)
모델의 주요 효과/상호 작용 부분은 다음과 같습니다.
> anova(Model)
# Type III Analysis of Variance Table with Satterthwaite's method
# Sum Sq Mean Sq NumDF DenDF F value Pr(>F)
# CondA 8546 8546 1 489.18 0.6694 0.4137
# CondB 2453 2453 1 509.78 0.1921 0.6613
# CondA:CondB 36951 36951 1 509.77 2.8944 0.0895 .
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
여기서 우리는 CondA 의 주요 효과를 예를 들어 검토합니다.
> PowerA_Ftest = simr::powerSim(fit = Model, # 要考察的模型
test = fixed('CondA', # 要考察的主效应/交互作用的名称
method = 'f'), # 选取检验方法，因为主效应为f检验，因此method设置为f
nsim=50) # 设置模拟次数，建议设置为500 (此时可以获取到较稳定的power)
> PowerA_Ftest
# Power for predictor 'CondA', (95% confidence interval):
# 12.00% ( 4.53, 24.31)
#
# Test: Type-II F-test (package car)
#
# Based on 50 simulations, (50 warnings, 0 errors)
# alpha = 0.05, nrow = 553
#
# Time elapsed: 0 h 0 m 58 s
#
# nb: result might be an observed power calculation
CondA 의 주요 효과의 힘은 0.12 입니다 (신뢰 구간은 0.045 – 0.243 ).
예를 들어, 이전 연구에서 비교적 안정적인 효과 크기가 얻어 졌다면, CondA 의 고정 효과는 일반적으로 약 50 인 것으로 밝혀졌으며, 이전 효과에 따라 전력을 계산할 수 있습니다. 단계는 다음과 같습니다.
> fixef(Model) # 查看模型固定效应
# (Intercept) CondAA2 CondBB2 CondAA2:CondBB2
# 279.43090 24.29565 12.18484 -32.82881
> fixef(Model)[2] = 50 # 修改CondA的效应
> fixef(Model) # 查看修改后的固定效应
# (Intercept) CondAA2 CondBB2 CondAA2:CondBB2
# 279.43090 50.00000 12.18484 -32.82881
> PowerA_ttest = simr::powerSim(fit = Model, test = fixed('CondAA2', method = 't'), nsim=50)
> PowerA_ttest
# Power for predictor 'CondAA2', (95% confidence interval):
# 88.00% (75.69, 95.47)
#
# Test: t-test with Satterthwaite degrees of freedom (package lmerTest)
# Effect size for CondAA2 is 50.
#
# Based on 50 simulations, (3 warnings, 0 errors)
# alpha = 0.05, nrow = 553
#
# Time elapsed: 0 h 0 m 9 s
현재 CondA 의 고정 효과 전력은 0.88 입니다.
> PowerA_Ftest = simr::powerSim(fit = Model, test = fixed('CondA', method = 'f'), nsim=50)
> PowerA_Ftest
# Power for predictor 'CondA', (95% confidence interval):
# 94.00% (83.45, 98.75)
#
# Test: Type-II F-test (package car)
#
# Based on 50 simulations, (50 warnings, 0 errors)
# alpha = 0.05, nrow = 553
#
# Time elapsed: 0 h 1 m 4 s
CondA 의 주요 효과는 0.94 입니다.
> Model = lmer(data = DemoData, DV~CondA*CondB+(1|subj)+(1|item)) # 这里仍以后验power为例
상호 작용의 고정 효과의 힘은 0.46 이며, 우리는 피험자의 수를 테스트하려고합니다. 전력은 0.8 도달 할 수 있습니다. 단계는 다음과 같습니다.
> Model2 = extend(object = Model, along = 'subj', n = 30)
> PowerA_ttest2 = powerSim(Model2, fixed('CondAA2', method = 't'), nsim=50)
> PowerA_ttest2
# Power for predictor 'CondAA2', (95% confidence interval):
# 84.00% (70.89, 92.83)
#
# Test: t-test with Satterthwaite degrees of freedom (package lmerTest)
# Effect size for CondAA2 is 24.
#
# Based on 50 simulations, (6 warnings, 0 errors)
# alpha = 0.05, nrow = 1659
#
# Time elapsed: 0 h 0 m 11 s
#
# nb: result might be an observed power calculation
현재 전력은 0.84 로 상승했습니다.
powerCurve() 함수를 사용할 수 있습니다. > Pcurve = powerCurve(fit = Model2, test = fixed('CondAA2', method = 't'), along = 'subj', nsim=50)
> Pcurve
# Power for predictor 'CondAA2', (95% confidence interval), by largest value of subj:
# 11: 22.00% (11.53, 35.96) - 180 rows
# 14: 36.00% (22.92, 50.81) - 325 rows
# 17: 50.00% (35.53, 64.47) - 507 rows
# 2: 58.00% (43.21, 71.81) - 667 rows
# 22: 64.00% (49.19, 77.08) - 837 rows
# 25: 58.00% (43.21, 71.81) - 989 rows
# 28: 74.00% (59.66, 85.37) - 1166 rows
# 30: 82.00% (68.56, 91.42) - 1318 rows
# 6: 86.00% (73.26, 94.18) - 1489 rows
# 9: 90.00% (78.19, 96.67) - 1659 rows
#
# Time elapsed: 0 h 1 m 39 s
현재 각 전력 앞의 값은 현재 부정확 하며 실제 값은 다음과 같습니다.
> Pcurve$nlevels
# [1] 3 6 9 12 15 18 21 24 27 30
약 24 명의 피험자가 권력에 도달하면 전력이 0.8 이상에 도달 할 수 있음을 알 수 있습니다.
샘플 크기 범위는 수동으로 설정할 수 있습니다.
> Pcurve = powerCurve(fit = Model2, test = fixed('CondAA2', method = 't'), along = 'subj', nsim=50,
breaks=c(21, 22, 23, 24, 25, 30)) # 考察样本量为21、22、23、24、25、30时的power
> Pcurve
# Power for predictor 'CondAA2', (95% confidence interval), by largest value of subj:
# 28: 72.00% (57.51, 83.77) - 1166 rows
# 29: 74.00% (59.66, 85.37) - 1220 rows
# 3: 70.00% (55.39, 82.14) - 1262 rows
# 30: 74.00% (59.66, 85.37) - 1318 rows
# 4: 76.00% (61.83, 86.94) - 1369 rows
# 9: 80.00% (66.28, 89.97) - 1659 rows
#
# Time elapsed: 0 h 1 m 12 s
> Pcurve$nlevels
# [1] 21 22 23 24 25 30
일반화 된 모델에서 효과의 전력 계산은 위와 유사하며 주요 차이점은 다음과 같습니다.
fixed() 의 method 매개 변수는 'z' 로 설정해야합니다.fixed() 의 method 매개 변수는 'chisq' 로 설정해야합니다. 일반화 된 모델의 주요 효과는 카이-제곱 테스트이기 때문입니다.
