Power analysis of HLM LMMダウンロード - Power analysis of HLM LMMダウンロード

線形混合モデルの電力分析

心理学の再現性がますます評価されるにつれて、統計分析のための電力計算の需要も増加しています。ほとんどの基本的な統計的方法（ t test 、 ANOVAなど）の電力計算は、一部のソフトウェア（ jamoviなど）で完了できます。ただし、高度な統計的方法として、視覚インターフェイスを備えたjamoviなどのソフトウェアで線形混合モデル（LMM）の電力計算は実行できず、 Rの関連データパケットと関数がまだ必要です。この記事では、主にRでLMMの電力計算を実行する方法を紹介します。これには、データとその機能の紹介、適用可能なシナリオと制限されたシナリオ、具体的な例デモンストレーションが含まれます。

データパケットとその機能

simrデータパケットは、 RのLMM電力分析により便利なパッケージです。分析では、モデルデータはモンテカルロシミュレーションによってシミュレートされ、セット数のシミュレーションでかなりの割合が計算されます。主な機能はpowerSim()です。

適用可能なシナリオ

事後効果のパワーを計算します。つまり、既存のデータに基づいて構築されたモデルの効果に基づいて電力を計算します。
先験的効果の力、つまり既存の研究の効果のサイズに基づいて電力を計算します。
電力のサイズは被験者またはアイテムの数の影響を受けるため、被験者/アイテムの数を推定します。特定のモデルは、より多くの数が多いほど電力が大きいことであり、これにより、実験の影響に基づいた正式な実験に必要な主題/項目の数を推定する機能が提供されます。

制限されたシナリオ

パワー分析全体に長い時間がかかります（このデータパケットの著者もこの問題を解決しようとしており、ユーザーはマルチスレッドコンピューティングを使用してそれを軽減することもできます）
アプリオリ効果のパワーが必要であっても、実験前のデータが必要です。

例

後方効果のパワーを計算します

ここでは、心理学実験（データダウンロードアドレス）のデータを使用してデモを行います。この実験では、認知タスク（ DV ）を実行するときに、10人の被験者（ subj ）の反応を収集しました。実験には2つの条件（ CondAとCondB ）があり、 2*2被験者とプロジェクト内（ item ）のために設計されていました。つまり、2つのランダムな要因が同時に考慮されます。

ここでは、固定効果と非ランダム効果に焦点を当てているため、モデリング時にランダム効果部分のランダムインターセプトのみを検討し、固定効果部分は2つの主要な効果とそれらの相互作用を調べます。

モデリングは次のとおりです。

 > 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

ここでは、例として次の例として、2つ（ 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達する可能性があります。手順は次のとおりです。

たとえば、被験者を拡大した後にモデルを確立するために、30人の被験者に拡大するために

 > 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

一般化された線形混合モデルの効果の力を計算する

一般化されたモデルの効果の電力計算は、上記に似ています。主な違いは次のとおりです。

固定効果を調べる場合、一般化モデルの固定効果はz検定であるためfixed()のmethodパラメーターを'z'に設定する必要があります。
主効果を調べるとき、一般化モデルの主効果はカイ二乗検定であるため、 fixed()のmethodパラメーターを'chisq'に設定する必要があります。

拡大する