【R講座】分散分析(ANOVA)

R講座

研究室に配属されたばかりの新入生や、これからRで統計分析を始めたいと思っている方へ向けて、【R講座】では、RとRStudioの基本的な使い方から統計手法の選び方、基本的なデータ分析方法を解説しています。特にRが初めての方でも安心して学べるように、RStudioのクリック操作も紹介していきます。実際のコード例を交えながら進めるので、これからの研究やデータ分析に、役立てていただけたら嬉しいです。

みなさん、こんにちは!

ここでは、この回で紹介した分散分析(ANOVA)について解説しています。

今回の内容
  • 分散分析の概要
  • 分散分析の方法
  • 結果の見方

分散分析(ANOVA)とは

分散分析では、データのばらつき(分散)を「要因による変動」と「誤差(個体差)による変動」に分解し、要因が影響を与えているかを検定します。

例えば、植物の成長 (growth) に対して、光 (light) の影響を調べる場合、成長のばらつきは次の2つの要因に分けられます。

  1. 要因による変動(処理変動): 光の違いによる成長量の変化
  2. 誤差変動(個体差や測定誤差): 同じ光条件でも成長量が異なる部分

この2つの変動を比べ、要因による変動が誤差変動よりも大きいかどうかをF検定で評価します。

分散分析の種類

分散分析には、いくつかの種類があります。ここでは代表的な次の分散分析を紹介します。

  • 一元配置分散分析:1つの独立変数(因子)を持つ分散分析。
  • 二元配置分散分析:2つの独立変数(因子)を持つ分散分析。
  • 対応のある分散分析:同じ被験者を複数回測定したデータに対する分散分析。

これらの他にも、要因の水準数や組み合わせによってさまざまな分散分析の方法があります。分散分析の種類は、研究の目的に応じて適切なものを選択します。

前提条件

分散分析を行う前に、以下の条件を確認しましょう。

  • 対応の有無
    • 対応のあるデータ:同じ対象の前後比較(例:同じ人のダイエット前後の体重)
    • 対応のないデータ:異なる対象の比較(例:異なる人々のグループ間比較)
  • 正規性:データが正規分布に従っているか。
  • 等分散性:各群の分散が等しいか。
  • 群数:3群以上。
  • データ尺度:データが間隔尺度または比例尺度であること。

仮説の設定

分散分析を行う前に、以下の仮説を設定します。

一元配置分散分析
  • 帰無仮説(H0:全ての群の平均値は等しい。
  • 対立仮説(H1:少なくとも1つの群の平均値が他と異なる。
二元配置分散分析

要因1と要因2にそれぞれ次の仮説を設定します。

  • 帰無仮説(H0:全ての群の平均値は等しい。
  • 対立仮説(H1:少なくとも1つの群の平均値が他と異なる。

また、要因1と要因2の交互作用について次の仮説を設置します。

  • 帰無仮説(H0:要因1と要因2は互いに影響しあわない。
  • 対立仮説(H1:要因1と要因2は互いに影響する。

関数の構造と引数オプション

Rで分散分析を行う方法はいくつかありますが、今回は、分散分析にaov()関数を使用する方法を紹介します。

関数の構造

aov(formula, data = NULL, projections = FALSE, qr = TRUE, contrasts = NULL, ...)

引数オプション
  • formula: 分析モデルの式(例:目的変数 ~ 説明変数
  • data: データフレーム
  • projections: 投影行列を含むかどうか
  • qr: QR分解を行うかどうか
  • contrasts: 対比の設定

まとめると、aov()関数は次のように使います。

aov(モデル式)

この分散分析モデルをsummary()関数を組み合わせることで分散分析を実行できます。

分析の実践

次のステップでデータを分析していきます。

Rで分散分析の手順
  • STEP 1
    データの読み込み

    csvファイルからデータを読み込みます。

  • STEP 2
    データの型変換(キャスト)

    読み込んだデータの要因データを、Factor型に変更します。

  • STEP 3
    分散分析の実行

    summary()関数とaov()関数で分散分析を実行します。

  • STEP 4
    結果の出力

    計算された結果がコンソールペインに出力されます。

使用するデータ

この講座では、説明のために同じ CSV データを使い回しています。
実際には、データの性質(分布・尺度・サンプル数など)に合わせて、適切な統計検定を選びましょう。

一元配置分散分析

一元配置分散分析は、次のcsvファイルを使用して説明します。

csvファイル

このcsvファイルには次のデータが含まれています。

Temperature (°C)Height (cm)Weight (g)
2511.40856335683434.032819857028
259.4582396378988954.9096737259706
2510.278664724106454.2160336538475
259.806027255329768.739038985953
2511.576158181218760.3451432394433
258.5244523647394850.8181031035401
259.8553917926894649.1747623798284
258.9249898091818356.060734308621
2510.406542731944941.125798546829
2512.229262201640951.0542139019377
3013.4855029917497103.528744733185
3014.938292577988105.503933584551
3014.85272920996188.6566903148316
3016.541593068827114.623515387464
3014.0181443311961107.021167106676
3015.4965781726617125.071111484834
3016.696947880723181.0997285637598
3014.739263691431994.1018720980883
3014.294071414332382.8549770315418
3014.838821493827795.790021021833
3515.501321827723778.1014137650469
3513.986460329505292.025705860145
3516.614752235468170.5661519563607
3515.005641984852563.0140291668191
3512.095100939654471.926190857745
3513.892835181031381.2105420273818
3516.547566932618376.8190219140427
3514.023169649653388.1840093119247
3514.898496552368372.0109068656934
3515.042650249796758.5177825784117
対応のある一元配置分散分析

対応のある一元配置分散分析は、次のcsvファイルを使用して説明します。

csvファイル

このcsvファイルには次のデータが含まれています。

IDTemperature (°C)Height (cm)Weight (g)
12511.4134.03
2259.4654.91
32510.2854.22
4259.8168.74
52511.5860.35
6258.5250.82
7259.8649.17
8258.9256.06
92510.4141.13
102512.2351.05
13013.49103.53
23014.94105.50
33014.8588.66
43016.54114.62
53014.02107.02
63015.50125.07
73016.7081.10
83014.7494.10
93014.2982.86
103014.8495.79
13515.5078.10
23513.9992.03
33516.6170.57
43515.0163.01
53512.1071.93
63513.8981.21
73516.5576.82
83514.0288.18
93514.9072.01
103515.0458.52
二元配置分散分析

二元配置分散分析は、次のcsvファイルを使用して説明します。

csvファイル

このcsvファイルには次のデータが含まれています。

Period (day)Temperature (°C)Height (cm)Weight (g)
302510.9548.76
30258.8964.67
302510.6256.74
302510.5169.56
302510.3747.31
602531.7287.55
602529.7996.04
602528.69100.97
602530.0697.62
602529.7795.88
303020.6459.23
303021.6367.03
303018.1964.04
303019.7978.08
303020.0778.45
603040.55165.40
603039.30146.71
603040.39159.48
603040.38145.21
603039.99134.85
対応のある二元配置分散分析

対応のある二元配置分散分析は、次のcsvファイルを使用して説明します。

csvファイル

このcsvファイルには次のデータが含まれています。

IDPeriod (day)Temperature (°C)Height (cm)Weight (g)
1302510.9548.76
230258.8964.67
3302510.6256.74
4302510.5169.56
5302510.3747.31
6602531.7287.55
7602529.7996.04
8602528.69100.97
9602530.0697.62
10602529.7795.88
1303020.6459.23
2303021.6367.03
3303018.1964.04
4303019.7978.08
5303020.0778.45
6603040.55165.40
7603039.30146.71
8603040.39159.48
9603040.38145.21
10603039.99134.85

データの読み込み

まずは、次のコードを使って、オブジェクト「data」にread.csv()関数でcsvファイルのデータを代入します。

# データの読み込み
data <- read.csv(file.choose(),
                 check.names = F)

データの前処理

次に、データの因子変換をします。

Rでは説明変数をas.factor()関数でFactor型のデータとします。

Factor型へ変換

# 説明変数をデータ型をFactor型に変更
data$列の名前 <- as.factor(data$列の名前)

例:一元配置分散分析の場合
# データ型の変換
data$`Temperature (°C)` <- as.factor(data$`Temperature (°C)`)

二元配置分散分析はもう一つの説明変数を同じようにFactor型へ変換します。また、対応のある分散分析はIDの部分をFactor型に変換します。

実行結果は、Environmentタブで確認できます。Temperature (°C)のデータ型がinteger型(左図)からFactor型(右図)に変わります。

バッククォート「`」について、少しだけ補足します。

R言語においてバッククォート「 ` 」で文字列を挟むと、通常は予約語や空白、特殊文字を含む名前を使うことができます。これにより、通常は無効な識別子として扱われる文字列を有効なオブジェクト名として扱うことができます。

ヘッダー名であるHeight (cm)は、スペース「 」や括弧「()」が含まれているため、そのままオブジェクト名として使用できませんが、バッククォートを使うことでオブジェクト名として使用できます。

分散分析と結果の見方

有意水準・p値・信頼区間
項目説明
有意水準 (α)帰無仮説を棄却する基準。通常 0.05(5%) や 0.01(1%) を使用する。
例: α = 0.05なら、5%未満の確率で偶然起こる差を「有意」と判断する。
p値検定統計量が観測された値以上になる確率。
p値 < 有意水準 (α) なら、統計的に有意と判断し、帰無仮説を棄却する。
統計的に有意とは?「偶然の変動では説明できない差がある」と判断すること。ただし「実験的に重要」や「因果関係がある」とは限らない。
信頼区間 (Confidence Interval, CI)母集団の真の値(2つのグループの平均値の差)が含まれる範囲を示す。例えば 95%信頼区間 は、繰り返し実験したときに95%の確率で真の値(平均値の差)を含む。
信頼区間と有意性の関係もし信頼区間が ゼロ(または比較対象の値)を含まなければ、統計的に有意と判断できる。
例:平均差の95%信頼区間が (0.5, 2.3) なら、有意水準5%で有意。
注意点統計的有意でも「効果の大きさ(実用的な意味)」とは異なる。
p値が大きくても「差がない」とは言えない(サンプル数が少ない可能性)。
一元配置分散分析

summary()aov()を組み込むことで結果が出力されます。

summary(aov(目的変数 ~ 説明変数))

# 一元配置分散分析
# データの読み込み
data <- read.csv(file.choose(), check.names = F)

# データ型の変換
data$`Temperature (°C)` <- as.factor(data$`Temperature (°C)`)

# aov関数を使用して分散分析を実行
summary(aov(data$`Height (cm)` ~ data$`Temperature (°C)`))
                        Df Sum Sq Mean Sq F value   Pr(>F)    
data$`Temperature (°C)`  2 143.09   71.55   50.21 7.99e-10 ***
Residuals               27  38.47    1.42                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
結果の見方

この結果では、p値( Pr(>F) )が0.05より小さいので、有意差があると判断します。

oneway.test()関数を使うことで等分散を仮定できないときの一元配置分散分析を行えます。

summary()aov()を組み合わせた計算結果と同じになります。

# 一元配置分散分析(等分散を仮定)
oneway.test(data$`Height (cm)` ~ data$`Temperature (°C)`, var.equal = TRUE)
    One-way analysis of means

data:  data$`Height (cm)` and data$`Temperature (°C)`
F = 50.211, num df = 2, denom df = 27, p-value = 7.992e-10
# 一元配置分散分析(等分散を仮定しない)
oneway.test(data$`Height (cm)` ~ data$`Temperature (°C)`, var.equal = FALSE)
    One-way analysis of means (not assuming equal variances)

data:  data$`Height (cm)` and data$`Temperature (°C)`
F = 50.411, num df = 2.000, denom df = 17.769, p-value = 4.741e-08
二元配置分散分析

一元配置分散分析モデルに二つ目の説明変数を組み合わせます。

summary(aov(目的変数 ~ 説明変数1 * 説明変数2))

# データ読み込み
data <- read.csv(file.choose(), check.names = FALSE)

# データの整理
data$`Temperature (°C)` <- as.factor(data$`Temperature (°C)`)
data$`Period (day)` <- as.factor(data$`Period (day)`)

# aov関数を使用して分散分析を実行
summary(aov(data$`Height (cm)` ~ data$`Temperature (°C)` * data$`Period (day)`))
                                            Df Sum Sq Mean Sq  F value  Pr(>F)    
data$`Temperature (°C)`                      1  495.6   495.6  537.476 9.7e-14 ***
data$`Period (day)`                          1 1979.9  1979.9 2147.256 < 2e-16 ***
data$`Temperature (°C)`:data$`Period (day)`  1    0.1     0.1    0.138   0.715    
Residuals                                   16   14.8     0.9                     
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
結果の見方

二元配置分散分析では、2の要因に対するp値と要因の交互作用のp値が出力されます。

有意水準を5%とする場合、p値が0.05未満であれば、対立仮説を採択し、少なくとも1つの群の平均値が他と異なると判断します。

この結果では、温度(data$`Temperature (°C)`)と期間(data$`Period (day)`)のp値( Pr(>F) )が0.05より小さいので、「有意な差がある」と判断します。また、交互作用(data$`Temperature (°C)`:data$`Period (day)` )は、p値( Pr(>F) )が0.05より大きいので、有意差はないと判断します。

対応のある一元配置分散分析

一元配置分散分析モデルにIDの説明変数を組み合わせます。

Error(ID / 説明変数)を指定すると、同じ個体の中での変動を誤差として扱いながら、処理効果を評価することになります。

summary(aov(目的変数 ~ 説明変数 + Error(ID / 説明変数)))

# データ読み込み
data <- read.csv(file.choose(), check.names = FALSE)

# データの整理
data$`Temperature (°C)` <- as.factor(data$`Temperature (°C)`)
data$ID <- as.factor(data$ID)

# aov関数を使用して分散分析を実行
summary(aov(data$`Height (cm)` ~ data$`Temperature (°C)` + Error(data$ID/data$`Temperature (°C)`)))
Error: data$ID
          Df Sum Sq Mean Sq F value Pr(>F)
Residuals  9   12.3   1.367                

Error: data$ID:data$`Temperature (°C)`
                        Df Sum Sq Mean Sq F value   Pr(>F)    
data$`Temperature (°C)`  2 143.09   71.55   49.21 5.05e-08 ***
Residuals               18  26.17    1.45                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
結果の見方

この結果では、p値( Pr(>F) )が0.05より小さいので、有意差があると判断します。

対応のある二元配置分散分析

二元配置分散分析モデルにIDの説明変数を組み合わせます。

Error(ID / 説明変数)を指定すると、同じ個体の中での変動を誤差として扱いながら、処理効果を評価することになります。

summary(aov(目的変数 ~ 説明変数1 * 説明変数2 + Error(ID / 説明変数1 * 説明変数2)))

# データ読み込み
data <- read.csv(file.choose(), check.names = FALSE)

# データの整理
data$`Temperature (°C)` <- as.factor(data$`Temperature (°C)`)
data$`Period (day)` <- as.factor(data$`Period (day)`)
data$ID <- as.factor(data$ID)

# aov関数を使用して分散分析を実行
summary(aov(data$`Height (cm)` ~ data$`Temperature (°C)` * data$`Period (day)` + Error(data$ID/(data$`Temperature (°C)` * data$`Period (day)`))))
Error: data$ID
          Df Sum Sq Mean Sq F value Pr(>F)
Residuals  4  4.763   1.191                

Error: data$ID:data$`Temperature (°C)`
                        Df Sum Sq Mean Sq F value   Pr(>F)    
data$`Temperature (°C)`  1  495.6   495.6   994.9 6.02e-06 ***
Residuals                4    2.0     0.5                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Error: data$ID:data$`Period (day)`
                    Df Sum Sq Mean Sq F value   Pr(>F)    
data$`Period (day)`  1 1979.9  1979.9   11082 4.88e-08 ***
Residuals            4    0.7     0.2                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Error: data$ID:data$`Temperature (°C)`:data$`Period (day)`
                                            Df Sum Sq Mean Sq F value Pr(>F)
data$`Temperature (°C)`:data$`Period (day)`  1  0.128  0.1277    0.07  0.804
Residuals                                    4  7.283  1.8207           
結果の見方

二元配置分散分析では、2の要因に対するp値と要因の交互作用のp値が出力されます。

有意水準を5%とする場合、p値が0.05未満であれば、対立仮説を採択し、少なくとも1つの群の平均値が他と異なると判断します。

この結果では、温度(data$`Temperature (°C)`)と期間(data$`Period (day)`)のp値( Pr(>F) )が0.05より小さいので、「有意な差がある」と判断します。また、交互作用(data$`Temperature (°C)`:data$`Period (day)` )は、p値( Pr(>F) )が0.05より大きいので、有意差はないと判断します。

Error()関数内の説明変数1 * 説明変数2を括弧で囲み忘れに注意しましょう。

Error(ID/(説明変数1 * 説明変数2))

Error(ID/説明変数1 * 説明変数2)

まとめ

まとめ
  • 分散分析の概要
    • 正規分布・3群以上・間隔尺度または比例尺度
  • 分散分析の方法
    • summary(aov())が基本の形
  • 結果の見方
    • 一元配置分散分析は要因に対するp値
    • 二元配置分散分析は要因に対するp値と交互作用に対するp値

ここでは、Rを使用して分散分析を行う方法について解説しました。

配布したcsvデータにはWeight (g)のデータもあるので、コードを書き換えて練習してみてください。

Rでは、コードの文法が正しければどんな統計手法でも計算が実行されてしまうため、結果が出たからといって安心するのは危険です。統計解析では、データの種類や前提条件に適した手法を選ぶことが重要であり、誤った検定を適用すると意味のない結果を得てしまう可能性があります。

そのため、単にコードを動かすだけでなく、

  • データの分布(正規分布かどうか)
  • 変数の尺度(カテゴリーデータか、連続データか)
  • サンプルサイズ(母数が少ないと適さない検定もある)

などを考慮し、統計の知識をもとに適切な解析方法を選択する必要があります。

タイトルとURLをコピーしました