R で一元配置分散分析を実行する方法
一元配置分散分析は、 3 つ以上の独立したグループの平均間に統計的に有意な差があるかどうかを判断するために使用されます。
このタイプの検定は、応答変数に対する予測変数の影響を分析するため、一元配置分散分析と呼ばれます。
注: 応答変数に対する 2 つの予測変数の影響に興味がある場合は、二元配置分散分析を実行できます。
R で一元配置分散分析を実行する方法
次の例は、R で一元配置分散分析を実行する方法を示しています。
背景
3 つの異なる運動プログラムが減量に異なる影響を与えるかどうかを判断したいとします。私たちが研究する予測変数は運動プログラムであり、応答変数はポンド単位で測定される体重減少です。
一元配置分散分析を実行して、3 つのプログラムによる体重減少の間に統計的に有意な差があるかどうかを判断できます。
実験に参加する人を 90 人募集し、30 人をランダムに割り当てて、プログラム A、プログラム B、またはプログラム C のいずれかを 1 か月間実行します。
次のコードは、作業するデータ フレームを作成します。
#make this example reproducible
set.seed(0)
#create data frame
data <- data.frame(program = rep(c("A", "B", "C"), each = 30),
weight_loss = c(runif(30, 0, 3),
runif(30, 0, 5),
runif(30, 1, 7)))
#view first six rows of data frame
head(data)
# program weight_loss
#1 A 2.6900916
#2 A 0.7965260
#3 A 1.1163717
#4 A 1.7185601
#5 A 2.7246234
#6 A 0.6050458
データ フレームの最初の列は、その人が 1 か月間参加したプログラムを示し、2 番目の列は、その人がプログラム終了時に経験した総体重減少をポンド単位で示します。
データを探索する
一元配置 ANOVA モデルを当てはめる前に、 dplyrパッケージを使用して 3 つのプログラムそれぞれの体重減少の平均と標準偏差を見つけることで、データをよりよく理解できます。
#load dplyr package library (dplyr) #find mean and standard deviation of weight loss for each treatment group data %>% group_by (program) %>% summarize (mean = mean(weight_loss), sd = sd(weight_loss)) # A tibble: 3 x 3 # program mean sd # #1 A 1.58 0.905 #2 B 2.56 1.24 #3 C 4.13 1.57
3 つのプログラムごとに箱ひげ図を作成して、各プログラムの体重減少の分布を視覚化することもできます。
#create boxplots
boxplot(weight_loss ~ program,
data = data,
main = "Weight Loss Distribution by Program",
xlab = "Program",
ylab = "Weight Loss",
col = "steelblue",
border = "black")
これらの箱ひげ図から、プログラム C の参加者の平均体重減少が最も高く、プログラム A の参加者の平均体重減少が最も低いことがわかります。
また、体重減少の標準偏差 (箱ひげ図の「長さ」) が、プログラム C では他の 2 つのプログラムに比べて少し高いこともわかります。
次に、一元配置分散分析モデルをデータに当てはめて、これらの視覚的な違いが実際に統計的に有意であるかどうかを確認します。
一元配置分散分析モデル フィッティング
R で一元配置分散分析モデルを近似するための一般的な構文は次のとおりです。
aov(応答変数 ~ 予測変数、データ = データセット)
この例では、次のコードを使用して、 weight_loss を応答変数として使用し、プログラムを予測変数として使用して、一元配置 ANOVA モデルを近似することができます。次に、 summary()関数を使用して、モデルの結果を表示します。
#fit the one-way ANOVA model model <- aov(weight_loss ~ program, data = data) #view the model output summary(model) # Df Sum Sq Mean Sq F value Pr(>F) #program 2 98.93 49.46 30.83 7.55e-11 *** #Residuals 87 139.57 1.60 #--- #Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
モデルの結果から、予測変数のプログラムは0.05 有意水準で統計的に有意であることがわかります。
言い換えれば、3 つのプログラムによる平均体重減少の間には統計的に有意な差があるということです。
モデルの仮定を確認する
次に進む前に、モデルの結果が信頼できるように、モデルの仮定が満たされていることを確認する必要があります。特に、一元配置分散分析では次のことを前提としています。
1. 独立性– 各グループの観察は互いに独立している必要があります。ランダム化された設計を使用したため(つまり、参加者を運動プログラムにランダムに割り当てました)、この仮定は満たされるはずなので、あまり心配する必要はありません。
2. 正規性– 従属変数は、予測変数の各レベルに対してほぼ正規分布を持つ必要があります。
3. 等分散– 各グループの分散は等しいかほぼ等しい。
正規性と等分散の仮定をチェックする 1 つの方法は、 plot()関数を使用することです。この関数は、4 つのモデル検査プロットを生成します。特に、次の 2 つのプロットに特に関心があります。
- 残差 vs.近似– このグラフは、残差と近似値の関係を示します。このグラフを使用して、グループ間の分散がほぼ等しいかどうかを大まかに評価できます。
- QQ プロット– このプロットには、理論的な分位数に対する標準化された残差が表示されます。このグラフを使用して、正規性の仮定が満たされているかどうかを大まかに評価できます。
次のコードを使用して、これらのモデル検査プロットを作成できます。
plot(model)
上のQQ グラフを使用すると、正規性の仮定を検証できます。理想的には、標準化された残差はプロットの直線の対角線に沿って配置されます。ただし、上のグラフでは、最初と最後に向かって残差が線から少しずれていることがわかります。これは、正規性の仮定が違反される可能性があることを示しています。
ザ・レジデュアルズ vs.上の調整されたグラフにより、分散が等しいという仮定を検証できます。理想的には、近似値の各レベルで残差が均等に分散されることが望まれます。
近似値が高くなるほど残差がさらに広がっていることがわかり、分散が等しいという仮定に違反している可能性があることがわかります。
等分散を正式にテストするには、 carパッケージを使用して Levene テストを実行できます。
#load car package library (car) #conduct Levene's Test for equality of variances leveneTest(weight_loss ~ program, data = data) #Levene's Test for Homogeneity of Variance (center = median) # Df F value Pr(>F) #group 2 4.1716 0.01862 * #87 #--- #Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
検定の p 値は0.01862です。有意水準 0.05 を使用すると、分散が 3 つのプログラム間で等しいという帰無仮説が棄却されます。ただし、有意水準 0.01 を使用すると、帰無仮説は棄却されません。
正規性と分散の等価性の仮定が確実に満たされるようにデータの変換を試みることもできますが、今のところはあまり心配しません。
治療法の違いを分析する
モデルの仮定が満たされている (または合理的に満たされている) ことを確認したら、事後テストを実行して、どの治療グループが互いに異なるかを正確に判断できます。
事後テストでは、 TukeyHSD()関数を使用して、多重比較の Tukey テストを実行します。
#perform Tukey's Test for multiple comparisons
TukeyHSD(model, conf.level=.95)
#Tukey multiple comparisons of means
# 95% family-wise confidence level
#
#Fit: aov(formula = weight_loss ~ program, data = data)
#
#$program
# diff lwr upr p adj
#BA 0.9777414 0.1979466 1.757536 0.0100545
#CA 2.5454024 1.7656076 3.325197 0.0000000
#CB 1.5676610 0.7878662 2.347456 0.0000199
p 値は、各プログラム間に統計的に有意な差があるかどうかを示します。結果は、各プログラムの平均体重減少の間に 0.05 有意水準で統計的に有意な差があることを示しています。
R のLot(TukeyHSD())関数を使用して、Tukey 検定の結果として得られる 95% 信頼区間を視覚化することもできます。
#create confidence interval for each comparison
plot(TukeyHSD(model, conf.level=.95), las = 2)
信頼区間の結果は仮説検定の結果と一致します。
特に、プログラム間の平均体重減少の信頼区間には値0が含まれていないことがわかります。これは、3 つのプログラム間の平均体重減少に統計的に有意な差があることを示しています。
これは、仮説検定のすべてのp 値が0.05 未満であることと一致しています。
一元配置分散分析結果のレポート
最後に、結果を要約した方法で一元配置分散分析の結果をレポートできます。
運動プログラムの効果を調べるために一元配置分散分析が実行されました。 体重減少(ポンド単位で測定)。体重減少に対する 3 つのプログラムの効果の間には、統計的に有意な差がありました (F(2, 87) = 30.83、p = 7.55e-11)。事後 Tukey の HSD テストを実行しました。
プログラム C の参加者の平均体重減少は、プログラム B の参加者の平均体重減少よりも大幅に大きくなっています (p < 0.0001)。
プログラム C の参加者の平均体重減少は、プログラム A の参加者の平均体重減少よりも大幅に大きくなっています (p < 0.0001)。
さらに、プログラム B の参加者の平均体重減少は、プログラム A の参加者の平均体重減少よりも有意に大きかった (p = 0.01)。
追加リソース
次のチュートリアルでは、一元配置分散分析に関する追加情報を提供します。