R で二元配置分散分析を実行する方法

二元配置 ANOVA (「分散分析」) は、2 つの因子に分割された 3 つ以上の独立したグループの平均間に統計的に有意な差があるかどうかを判断するために使用されます。

このチュートリアルでは、R で二元配置分散分析を実行する方法について説明します。

例: R での二元配置分散分析

運動強度と性別が減量に影響を与えるかどうかを判断したいとします。この場合、注目している 2 つの要因は運動量と性別であり、応答変数はポンド単位で測定される体重減少です。

二元配置分散分析を実行して、運動と性別が体重減少に影響を与えるかどうか、また運動と性別の間に体重減少に関する相互作用があるかどうかを判断できます。

実験に参加する男性 30 人、女性 30 人を募集します。この実験では、各 10 人が無作為に割り当てられ、運動なし、軽い運動、または激しい運動プログラムを 1 か月間実行します。

次のコードは、作業するデータ フレームを作成します。

 #make this example reproducible
set.seed(10)

#create data frame
data <- data.frame(gender = rep(c("Male", "Female"), each = 30),
                   exercise = rep(c("None", "Light", "Intense"), each = 10, times = 2),
                   weight_loss = c(runif(10, -3, 3), runif(10, 0, 5), runif(10, 5, 9),
                                   runif(10, -4, 2), runif(10, 0, 3), runif(10, 3, 8)))

#view first six rows of data frame
head(data)

# gender exercise weight_loss
#1 Male None 0.04486922
#2 Male None -1.15938896
#3 Male None -0.43855400
#4 Male None 1.15861249
#5 Male None -2.48918419
#6 Male None -1.64738030

#see how many participants are in each group
table(data$gender, data$exercise)

# Intense Light None
# Female 10 10 10
# Male 10 10 10

データを探索する

二元配置 ANOVA モデルをフィッティングする前に、 dplyrパッケージを使用して 6 つの治療グループのそれぞれの体重減少の平均と標準偏差を見つけることで、データをよりよく理解できます。

 #load dplyr package
library(dplyr)

#find mean and standard deviation of weight loss for each treatment group
data %>%
  group_by (gender, exercise) %>%
  summarize (mean = mean(weight_loss),
            sd = sd(weight_loss))

# A tibble: 6 x 4
# Groups: gender [2]
# gender exercise means sd
#          
#1 Female Intense 5.31 1.02 
#2 Female Light 0.920 0.835
#3 Female None -0.501 1.77 
#4 Male Intense 7.37 0.928
#5 Male Light 2.13 1.22 
#6 Male None -0.698 1.12 

6 つの治療グループごとに箱ひげ図を作成して、各グループの体重減少の分布を視覚化することもできます。

 #set margins so that axis labels on boxplot don't get cut off
by(mar=c(8, 4.1, 4.1, 2.1))

#create boxplots
boxplot(weight_loss ~ gender:exercise,
data = data,
main = "Weight Loss Distribution by Group",
xlab = "Group",
ylab = "Weight Loss",
col = "steelblue",
border = "black", 
las = 2 #make x-axis labels perpendicular
)

激しい運動に参加した 2 つのグループの方が、体重減少値が高いことがすぐにわかります。また、激しい運動グループと軽い運動グループの両方で、男性の方が女性よりも体重減少値が高い傾向があることもわかります。

次に、二元配置分散分析モデルをデータに当てはめて、これらの視覚的な違いが実際に統計的に有意であるかどうかを確認します。

二元配置分散分析モデルのフィッティング

R で二元配置 ANOVA モデルを近似するための一般的な構文は次のとおりです。

aov(応答変数 ~predictor_variable1 *predictor_variable2, data = データセット)

2 つの予測子変数間の*は、2 つの予測子変数間の交互作用効果もテストすることを示していることに注意してください。

この例では、次のコードを使用して、 weight_loss を応答変数として、性別と運動量を 2 つの予測変数として使用して、二元配置分散分析モデルを近似することができます。

次に、 summary()関数を使用して、モデルの結果を表示します。

 #fit the two-way ANOVA model
model <- aov(weight_loss ~ gender * exercise, data = data)

#view the model output
summary(model)

# Df Sum Sq Mean Sq F value Pr(>F)    
#gender 1 15.8 15.80 11.197 0.0015 ** 
#exercise 2 505.6 252.78 179.087 <2e-16 ***
#gender:exercise 2 13.0 6.51 4.615 0.0141 *  
#Residuals 54 76.2 1.41                   
#---
#Significant. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

モデルの結果から、性別、運動、および 2 つの変数間の相互作用はすべて、有意水準 0.05 で統計的に有意であることがわかります。

モデルの仮定を確認する

次に進む前に、モデルの結果が信頼できるように、モデルの仮定が満たされていることを確認する必要があります。特に、二元配置分散分析では次のことを前提としています。

1. 独立性– 各グループの観察は互いに独立している必要があります。ランダム化された設計を使用したため、この仮定は満たされるはずなので、あまり心配する必要はありません。

2. 正規性– 従属変数は、2 つの因子のグループの組み合わせごとにほぼ正規分布を持つ必要があります。

この仮定をテストする 1 つの方法は、モデル残差のヒストグラムを作成することです。残差がほぼ正規分布している場合、この仮定は満たされるはずです。

 #define model residuals
reside <- model$residuals

#create histogram of residuals
hist(resid, main = "Histogram of Residuals", xlab = "Residuals", col = "steelblue")

残差はほぼ正規分布しているため、正規性の仮定が満たされていると想定できます。

3. 等分散– 各グループの分散は等しいかほぼ等しい。

この仮定を確認する 1 つの方法は、 carパッケージを使用して分散の等価性について Levene 検定を実行することです。

 #load car package
library(car)

#conduct Levene's Test for equality of variances
leveneTest(weight_loss ~ gender * exercise, data = data)

#Levene's Test for Homogeneity of Variance (center = median)
# Df F value Pr(>F)
#group 5 1.8547 0.1177
#54  

検定の p 値は有意水準 0.05 より大きいため、グループ間の分散が等しいという仮定が満たされていると考えることができます。

治療法の違いを分析する

モデルの仮定が満たされていることを確認したら、事後テストを実行して、どの治療グループが互いに異なるかを正確に判断できます。

事後テストでは、 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 ~ gender * exercise, data = data)
#
#$gender
# diff lwr upr p adj
#Male-Female 1.026456 0.4114451 1.641467 0.0014967
#
#$exercise
# diff lwr upr p adj
#Light-Intense -4.813064 -5.718493 -3.907635 0.0e+00
#None-Intense -6.938966 -7.844395 -6.033537 0.0e+00
#None-Light -2.125902 -3.031331 -1.220473 1.8e-06
#
#$`gender:exercise`
# diff lwr upr p adj
#Male:Intense-Female:Intense 2.0628297 0.4930588 3.63260067 0.0036746
#Female:Light-Female:Intense -4.3883563 -5.9581272 -2.81858535 0.0000000
#Male:Light-Female:Intense -3.1749419 -4.7447128 -1.60517092 0.0000027
#Female:None-Female:Intense -5.8091131 -7.3788841 -4.23934219 0.0000000
#Male:None-Female:Intense -6.0059891 -7.5757600 -4.43621813 0.0000000
#Female:Light-Male:Intense -6.4511860 -8.0209570 -4.88141508 0.0000000
#Male:Light-Male:Intense -5.2377716 -6.8075425 -3.66800066 0.0000000
#Female:None-Male:Intense -7.8719429 -9.4417138 -6.30217192 0.0000000
#Male:None-Male:Intense -8.0688188 -9.6385897 -6.49904786 0.0000000
#Male:Light-Female:Light 1.2134144 -0.3563565 2.78318536 0.2185439
#Female:None-Female:Light -1.4207568 -2.9905278 0.14901410 0.0974193
#Male:None-Female:Light -1.6176328 -3.1874037 -0.04786184 0.0398106
#Female:None-Male:Light -2.6341713 -4.2039422 -1.06440032 0.0001050
#Male:None-Male:Light -2.8310472 -4.4008181 -1.26127627 0.0000284
#Male:None-Female:None -0.1968759 -1.7666469 1.37289500 0.9990364

p 値は、各グループ間に統計的に有意な差があるかどうかを示します。

たとえば、上の最後の行では、男性の運動なしグループは、女性の運動なしグループと比較して、体重減少において統計的に有意な差が見られなかったことがわかります (p 値: 0.990364)。

R のLot()関数を使用して、Tukey 検定の結果として得られる 95% 信頼区間を視覚化することもできます。

 #set axis margins so labels don't get cut off
by(mar=c(4.1, 13, 4.1, 2.1))

#create confidence interval for each comparison
plot(TukeyHSD(model, conf.level=.95), las = 2)

二元配置分散分析結果のレポート

最後に、結果を要約した方法で二元配置分散分析の結果をレポートできます。

二元配置分散分析を実行して、性別 (男性、女性)と運動プログラム(なし、軽い、激しい)が体重減少(ポンドで測定) に及ぼす影響を調べました。体重減少に対する性別と運動の効果の間には、統計的に有意な相互作用がありました (F(2, 54) = 4.615、p = 0.0141)。事後 Tukey の HSD テストを実行しました。

男性の場合、激しい運動プログラムは、軽いプログラム (p < 0.0001) または運動プログラムなし(p < 0.0001) よりも大幅に大きな体重減少をもたらしました。さらに、男性では、軽い食事療法は、運動療法を行わない場合よりも有意に大きな体重減少をもたらしました（p < 0.0001）。

女性の場合、激しい運動プログラムは、軽いプログラム (p < 0.0001) または運動プログラムなし(p < 0.0001) よりも大幅に大きな体重減少をもたらしました。

正規性チェックと Levene 検定を実行して、ANOVA の仮定が満たされていることを確認しました。