人工知能を作りたい!

統計的仮説検定

くり返し数が等しい二元配置分散分析(対応なし)




はじめに

分散分析において、母平均に差をもたらす原因を因子といいます。また、因子の中に含まれる個々の条件を水準といいます。そして、因子を2つ配置した場合の分散分析を二元配置分散分析といいます。

因子Aの水準と因子Bの水準のどの組み合わせでも測定回数が等しい場合、くり返し数が等しいといいます。因子Aの水準と因子Bの水準のどの組み合わせでも測定回数が1回のみの場合、くり返しがないといいます。

くり返しのない二元配置分酸分析の場合、数式の構造上、交互作用効果を計算することはできません。ですので、もし、2つの因子に交互作用がありそうなら、くり返し数が等しいデータを得た後、くり返し数が等しい二元配置分散分析を行った方が良いです。
全て異なる被験者から得たデータ(被験者間のデータ)は対応がないデータといいます。同じ被験者から異なる水準で得たデータ(被験者内のデータ)は対応があるデータといいます。




例題

健康な人の血糖値は、空腹時は約70 mg/dL、食後2時間後は約130 mg/dLです。食後のインスリン分泌に遅れのある2型糖尿病患者の血糖値は、空腹時は約70 mg/dLで正常値ですが、食後2時間後は約200 mg/dL以上になります。

新しく開発したインスリン分泌促進薬の効果を確かめるため、インスリン分泌に遅れのある2型糖尿病患者群から18人の被験者を無作為に抽出しました。

それぞれの被験者に、新薬を食前10分前服用してもらい、食後2時間後に血糖値を測定しました。そして、「新薬を服用しない時の食後血糖値」と「新薬を服用した時の食後血糖値」の差を血糖降下度 (mg/dL) としてデータ化しました。

因子とその水準を以下の通り設定しました。

因子A:種類の違い(水準:新薬A、新薬B)
因子B:服用量の違い(水準:100μg、500μg、1000μg)

6通りある因子A水準と因子B水準の組み合わせそれぞれに、18人の被験者から3人ずつを割り当てました。つまり、以下のデータは、繰り返し数3回で対応のないデータです。

   100μg  500μg  1000μg
 薬A  30
 25
 28
 49
 59
 68
 80
 76
 82
 薬B  80
 65
 59
 79
 89
 78
 70
 86
 65

種類と服用量の違いによって、血糖降下度の母平均は異なるといえるでしょうか。




主効果と交互作用

くり返し数が等しい二元配置分酸分析では、二つの因子それぞれの主効果と、因子同士がお互いに影響を与え合って生じる交互作用を考える必要があります。

因子Aと因子Bに交互作用がない場合のデータを以下に示します。


   水準B1  水準B2  水準B3
 水準A1  3
 4
 5
 7
 8
 9
 4
 5
 6
 水準A2  5
 6
 7
 9
 10
 11
 6
 7
 8

「因子Aの水準と因子Bの水準の組み合わせ」ごとの平均値を求めると以下のようになります。


   水準B1  水準B2  水準B3
 水準A1  4  8  5
 水準A2  6  10  7

関数 interaction.plot() を使って、平均値のプロットを描いてみましょう。

関数 interaction.plot の書式は
interaction.plot(横軸に置く一つ目の因子, 軸外に置く2つ目の因子, 平均値を求めるデータ)
です。

ます、関数 interaction.plot() に入力する変数を作成します。


> データ <- c(3, 4, 5, 5, 6, 7, 7, 8, 9, 9, 10, 11, 4, 5,6, 6, 7, 8)


平均していない数値を、変数「データ」に格納します。


> 文字型の因子A <- rep(c(rep("水準A1", 3), rep("水準A2", 3)), 3) > 因子A <- factor(文字型の因子A)
> 因子A  [1] 水準A1 水準A1 水準A1 水準A2 水準A2 水準A2 水準A1 水準A1 水準A1 水準A2 [11] 水準A2 水準A2 水準A1 水準A1 水準A1 水準A2 水準A2 水準A2 Levels: 水準A1 水準A2


関数 rep(値, 回数) は、指定した値を指定した回数分、繰り返したデータを作成します。関数 rep() を利用して、数値の順番に対応するように水準A1と水準A2を並べ、変数「文字型の因子A」に格納します。次に、変数「文字型の因子A」を文字型データから因子型データに変換してから、変数「因子A」に格納します。最後に確認のため出力しました。


> 文字型の因子B <- c(rep("水準B1", 6), rep("水準B2", 6), rep("水準B3", 6)) > 因子B <- factor(文字型の因子B)
> 因子B  [1] 水準B1 水準B1 水準B1 水準B1 水準B1 水準B1 水準B2 水準B2 水準B2 水準B2 [11] 水準B2 水準B2 水準B3 水準B3 水準B3 水準B3 水準B3 水準B3 Levels: 水準B1 水準B2 水準B3


数値の順番に対応するように水準B1と水準B2と水準B3を並べてから、変数「文字型の因子B」に格納します。次に変数「文字型の因子B」を文字型データから因子型データに変換してから、変数「因子B」に格納します。最後に確認のため出力しました。


変数の下準備ができたので、関数 interaction.plot() を実行します。


> interaction.plot(因子A, 因子B, データ)


因子Aを横軸に取ります。平均値を結んだ3つの直線はすべて平行になっています。


> interaction.plot(因子B, 因子A, データ)


因子Bを横軸に取ります。平均値を結んだ2つの折れ線はすべて平行になっています。



因子Aと因子Bに交互作用がある場合のデータを以下に示します。


   水準B1  水準B2  水準B3
 水準A1  6
 7
 8
 7
 8
 9
 8
 9
 10
 水準A2  5
 6
 7
 10
 11
 12
 13
 14
 15


「因子Aの水準と因子Bの水準の組み合わせ」ごとの平均値を求めると以下のようになります。


   水準B1  水準B2  水準B3
 水準A1  7  8  9
 水準A2  6  11  14

関数 interaction.plot() を使って、平均値のプロットを描いてみましょう。まず、変数の下準備をします。


> データ <- c(6, 7, 8, 5, 6, 7, 7, 8, 9, 10, 11, 12, 8, 9,10, 13, 14, 15)


> 因子A <- factor(rep(c(rep("水準A1", 3), rep("水準A2", 3)), 3))
> 因子A
 [1] 水準A1 水準A1 水準A1 水準A2 水準A2 水準A2 水準A1 水準A1 水準A1 水準A2
[11] 水準A2 水準A2 水準A1 水準A1 水準A1 水準A2 水準A2 水準A2
Levels: 水準A1 水準A2


> 因子B <- factor(c(rep("水準B1", 6), rep("水準B2", 6), rep("水準B3", 6)))
> 因子B
 [1] 水準B1 水準B1 水準B1 水準B1 水準B1 水準B1 水準B2 水準B2 水準B2 水準B2
[11] 水準B2 水準B2 水準B3 水準B3 水準B3 水準B3 水準B3 水準B3
Levels: 水準B1 水準B2 水準B3


変数の下準備ができたので、関数 interaction.plot() を実行します。


> interaction.plot(因子A, 因子B, データ)


因子Aを横軸に取ります。平均値を結んだ3つの直線は平行になっていません。



> interaction.plot(因子B, 因子A, データ)


因子Bを横軸に取ります。平均値を結んだ2つの折れ線はすべて平行になっていません。


平均値を結んだ直線や折れ線が平行である場合、グラフで見た通りに交互作用効果はないと言えます。

しかし、平行でない場合、交互作用があるように見えますが、分散分析表を出力してみると、交互作用の効果はそれほどではなく、「有意ではない」と判定されることもあります。

ですので、交互作用効果があるかどうかは、きちんと分散分析表を出力して、最終確認する必要があります。




帰無仮説と対立仮説の設定

因子Aの主効果
帰無仮説 H0:因子Aの水準が違っても母平均は等しい(因子Aに主効果はない)
対立仮説 H1:因子Aの水準の違いにより母平均に差がある(因子Aには主効果がある)

因子Bの主効果
帰無仮説 H0:因子Bの水準が違っても母平均は等しい(因子Bに主効果はない)
対立仮説 H1:因子Bの水準の違いにより母平均に差がある(因子Bには主効果がある)

因子Aと因子Bの交互作用効果
帰無仮説 H0:因子Aの水準と因子Bの水準は互いに影響を及ぼしあっていない(因子Aと因子Bに交互作用はない)
対立仮説 H1:因子Aの水準と因子Bの水準は互いに影響を及ぼしあっている(因子Aと因子Bには交互作用がある)




Rによる統計解析

関数 aov() を使って、くり返し数の等しい二元配置分散分析を実行してみましょう。

関数 aov() の書式は
aov(データ~一つ目の因子*二つ目の因子)
です。「一つ目の因子*二つ目の因子」は、「一つ目の因子の主効果」と「二つ目の因子の主効果」と「一つ目の因子と二つ目の因子の交互作用」の3つすべてを含めるという意味です。

まず、関数 aov() に入力する3つの変数を作成しましょう。


> 血糖降下度 <- c(30, 25, 28,49, 59, 68, 80, 76, 82, 80, 65, 59, 79, 89, 78, 70, 86, 65)
> 血糖降下度
  [1] 30 25 28 49 59 68 80 76 82 80 65 59 79 89 78 70 86 65
> 種類 <- factor(c(rep("A", 9), rep("B", 9)))
> 種類
 [1] A A A A A A A A A B B B B B B B B B
Levels: A B
> 服用量 <- factor(rep(c(rep("100", 3), rep("500", 3), rep("1000", 3)), 2))
> 服用量
 [1] 100 100 100 500 500 500 1000 1000 1000 100 100 100 500 500
[15] 500 1000 1000 1000
Levels: 100 1000 500


関数 aov() を実行する準備が整ったので、くり返し数の等しい二元分散分析を実行してみましょう。


> summary(aov(血糖降下度~種類*服用量))

            Df Sum Sq Mean Sq F value Pr(>F)
種類 1 1682.0 1682.0 26.54 0.000241 ***
服用量 2 2732.1 1366.1 21.55 0.000107 ***
種類:服用量 2 1623.0 811.5 12.80 0.001056 **
Residuals 12 760.7 63.4
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


分散分析表のp値を見ると以下の結果が導かれます。

●種類:5%水準で有意な効果がある (p = 0.000241)
つまり、「対立仮説 H1:因子Aの水準の違いにより母平均に差がある(因子Aには主効果がある)」が採択されました。

●服用量:5%水準で有意な効果がある (p = 0.000107)
つまり、「対立仮説 H1:因子Bの水準の違いにより母平均に差がある(因子Bには主効果がある)」が採択されました。

●種類と服用量の交互作用:5%水準で有意な効果がある (p = 0.001056)
つまり、「対立仮説 H1:因子Aの水準と因子Bの水準は互いに影響を及ぼしあっている(因子Aと因子Bには交互作用がある)」が採択されました。


関数 interaction.plot() を使って、平均値のプロットを描いてみましょう。


> interaction.plot(種類, 服用量, 血糖降下度)



> interaction.plot(服用量, 種類, 血糖降下度)




上記の「主効果と交互作用」の章で載せた、交互作用がないデータと交互作用があるデータは、「因子Aの水準と因子Bの水準の組み合わせごとの平均値」を折れ線グラフにしただけで、本当のところ交互作用があるのかないのかはわからないままです。

それぞれのデータの分散分析表を出力してみましょう。



交互作用がないデータの分散分析表

> データ <- c(3, 4, 5, 5, 6, 7, 7, 8, 9, 9, 10, 11, 4, 5,6, 6, 7, 8)
> 因子A <- factor(rep(c(rep("水準A1", 3), rep("水準A2", 3)), 3))
> 因子B <- factor(c(rep("水準B1", 6), rep("水準B2", 6), rep("水準B3", 6)))
> summary(aov(データ~因子A*因子B))

            Df Sum Sq Mean Sq F value Pr(>F)
因子A 1 18 18 18 0.00114 **
因子B 2 52 26 26 4.35e-05 ***
因子A:因子B 2 0 0 0 1.00000
Residuals 12 12 1
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1


因子Aと因子Bの交互作用を表す「因子A:因子B 」の行のp値は1ですので、5%水準で有意な効果はありません。つまり、「帰無仮説 H0:因子Aと因子Bに交互作用はない」が採択されました。



交互作用があるデータの分散分析表

> データ <- c(6, 7, 8, 5, 6, 7, 7, 8, 9, 10, 11, 12, 8, 9,10, 13, 14, 15)
> 因子A <- factor(rep(c(rep("水準A1", 3), rep("水準A2", 3)), 3))
> 因子B <- factor(c(rep("水準B1", 6), rep("水準B2", 6), rep("水準B3", 6)))
> summary(aov(データ~因子A*因子B))

            Df Sum Sq Mean Sq F value Pr(>F)
因子A 1 24.5 24.5 24.5 0.000337 ***
因子B 2 76.0 38.0 38.0 6.43e-06 ***
因子A:因子B 2 28.0 14.0 14.0 0.000729 ***
Residuals 12 12.0 1.0
---
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’


因子Aと因子Bの交互作用を表す「因子A:因子B 」の行のp値は0.000729ですので、5%水準で有意な効果があります。つまり、「帰無仮説 H0:因子Aと因子Bに交互作用はない」は棄却され、「対立仮説 H1:因子Aと因子Bに交互作用がある」が採択されました。