【統計】二元配置分散分析(Two-way ANOVA)

本記事では二元配置分散分析(Two-way ANOVA)について記載します。

【目次】

 

前回の 一元配置分散分析 は1要因に対する分析手法でした。
2要因の場合は、二元配置分散分析を用います。
※本記事の二元配置分散分析は、対応なしの被験者間実験計画(ABSタイプ)です。
 

帰無仮説、対立仮設


 要因Aの主効果

  • 帰無仮説  H_{0}: 要因Aの水準間の主効果が同じ(平均が同じ)
  • 対立仮設  H_{1}: 要因Aの水準間のうち、少なくとも一つの水準で主効果に差がある(平均に差がある)

要因Bの主効果

  • 帰無仮説  H_{0}: 要因Bの水準間の主効果が同じ(平均が同じ)
  • 対立仮設  H_{1}: 要因Bの水準間のうち、少なくとも一つの水準で主効果に差がある(平均に差がある)

要因Aと要因Bの交互作用

  • 帰無仮説  H_{0}: 要因Aと要因Bに交互作用がある
  • 対立仮設  H_{1}: 要因Aと要因Bに交互作用がない


二元配置分散分析では、要因Aの水準間の平均の違い、要因Bの水準間の平均の違い、要因AとBの交互作用の有無を評価できます。

また、要因AとBに交互作用が無いと分かっている場合は交互作用を含めないモデルで分析可能です。交互作用が無いとは、例えば「要因A1の時の要因B1とB2の差」と「要因A2の時の要因B1とB2の差」に違いがないと言うことです。

  

計算式等


※本計算は各水準のサンプルサイズが等しいことを前提とします。また、対応なしの被験者間実験計画(ABSタイプ)です。

交互作用項を含めた分析
■ 平方和
\quad S_{A} = \sum _{k=1} ^{a} n_{k} \cdot \left( \overline{y_{k}} - \bar{\bar{y}} \right)^{2}
\quad S_{B} = \sum _{m=1} ^{b} n_{m} \cdot \left( \overline{y_{m}} - \bar{\bar{y}} \right)^{2}
\quad S_{AB} = \sum _{k=1} ^{a} \sum _{m=1} ^{b} n_{km} \cdot \left( \overline{y_{km}} - \overline{y_{k}} - \overline{y_{m}} + \bar{\bar{y}} \right)^{2}
\quad S_{e} = \sum _{k=1} ^{a} \sum _{m=1} ^{b} \sum _{i=1} ^{n_{km}} \left( y_{km_{i}} - \overline{y_{km}} \right)^{2}
\quad S_{T} = \sum _{i=1} ^{n} \left( y_{i} - \bar{\bar{y}} \right)^{2}

■ 分散分析表

変動因 平方和 自由度 平均平方 F値(F比)
水準間A S_{A} f_{A}=a-1 V_{A}=\dfrac{S_{A}}{f_{A}} F_{A}=\dfrac{V_{A}}{V_{e}}
水準間B S_{B} f_{B}=b-1 V_{B}=\dfrac{S_{B}}{f_{B}} F_{B}=\dfrac{V_{B}}{V_{e}}
交互作用AB S_{AB} f_{AB}=a\cdot b-a-b+1 V_{AB}=\dfrac{S_{AB}}{f_{AB}} F_{AB}=\dfrac{V_{AB}}{V_{e}}
残差e S_{e} f_{e}=n-a \cdot b V_{e}=\dfrac{S_{e}}{f_{e}}  
全体T S_{T} f_{T}=n-1    

※ k=要因Aの特定水準、m=要因Bの特定水準、i=特定の観測値、\overline{y_{k}}=要因Aの特定水準の平均値、\overline{y_{m}}=要因Bの特定水準の平均値、\bar{\bar{y}}=全体の平均値(総平均)、a=要因Aの水準数、b=要因Bの水準数

ここで平方和と自由度には次の関係があります。
\quad S_{T} = S_{A} + S_{B} + S_{AB} + S_{e}
\quad f_{T} = f_{A} + f_{B} + f_{AB} + f_{e}

 

交互作用項を含めない分析
■ 平方和 ※ S_{e} 以外同じ。交互作用項を残差項にプールしています。
\quad S_{A} = \sum _{k=1} ^{a} n_{k} \cdot \left( \overline{y_{k}} - \bar{\bar{y}} \right)^{2}
\quad S_{B} = \sum _{m=1} ^{b} n_{m} \cdot \left( \overline{y_{m}} - \bar{\bar{y}} \right)^{2}
\quad S_{e} = \sum _{k=1} ^{a} \sum _{m=1} ^{b} \sum _{i=1} ^{n_{km}} \left( y_{km_{i}} - \overline{y_{k}} - \overline{y_{m}} + \bar{\bar{y}} \right)^{2}
\quad S_{T} = \sum _{i=1} ^{n} \left( y_{i} - \bar{\bar{y}} \right)^{2}

■ 分散分析表 ※ S_{e} 以外同じ。

変動因 平方和 自由度 平均平方 F値(F比)
水準間A S_{A} f_{A}=a-1 V_{A}=\dfrac{S_{A}}{f_{A}} F_{A}=\dfrac{V_{A}}{V_{e}}
水準間B S_{B} f_{B}=b-1 V_{B}=\dfrac{S_{B}}{f_{B}} F_{B}=\dfrac{V_{B}}{V_{e}}
残差e S_{e} f_{e}=n-a-b+1 V_{e}=\dfrac{S_{e}}{f_{e}}  
全体T S_{T} f_{T}=n-1    

※ k=要因Aの特定水準、m=要因Bの特定水準、i=特定の観測値、\overline{y_{k}}=要因Aの特定水準の平均値、\overline{y_{m}}=要因Bの特定水準の平均値、\bar{\bar{y}}=全体の平均値(総平均)、a=要因Aの水準数、b=要因Bの水準数

ここで平方和と自由度には次の関係があります。
\quad S_{T} = S_{A} + S_{B} + S_{e}
\quad f_{T} = f_{A} + f_{B} + f_{e}

 

また、質的変数をダミー変数(0と1の数量変数)にして重回帰分析を行う(正規方程式を解く)ことで、各水準の条件付き期待値を算出できます(数量化理論1類)。
こちらは 一元配置分散分析 の記事を参照ください。

 

ノート


簡単なデータで計算方法を確認します。
今回は薬剤(A薬、B薬、C薬の3水準)と前処置(有、無の2水準)の2要因において、評価指標Xの結果に差があるかを見ていきます(架空データです)。交互作用項を含めた分析を行います。なお、評価指標Xは特定のものではありませんが、低いほど良いとしておきます。

 A薬+前処置無 {20,21,22,20,21}
 A薬+前処置有 {18,16,17,16,17}
 B薬+前処置無 {19,18,17,19,19}
 B薬+前処置有 {19,18,18,17,17}
 C薬+前処置無 {20,18,16,19,18}
 C薬+前処置有 {15,14,17,16,18}

 

まず、全体・薬剤別・前処置別で平均を求め、総平均からのズレ(偏差平方和)を求めます。

薬剤 前処置 N 平均 要因:薬剤によるズレ(偏差平方和)と
要因:前処置によるズレ(偏差平方和)
全体 全体 30 18.0 ---
A薬 全体 10 18.8  \left(18.8-18.0\right)^{2}\cdot 10=6.4
B薬 全体 10 18.1  \left(18.1-18.0\right)^{2}\cdot 10=0.1
C薬 全体 10 17.1  \left(17.1-18.0\right)^{2}\cdot 10=8.1
全体 15 19.1  \left(19.1-18.0\right)^{2}\cdot 10≒19.267
全体 15 16.9  \left(16.9-18.0\right)^{2}\cdot 10≒19.267

同じく、薬剤×前処置別で平均を求め、総平均・各要因の平均からのズレ(交互作用に関する変動)を求めます。

薬剤 前処置 N 平均 総平均・各要因の平均からのズレ
(交互作用に関する変動)
A薬 5 20.8  \left(20.8-18.8-19.1+18.0\right)^{2}\cdot 5≒3.756
A薬 5 16.8  \left(16.8-18.8-16.9+18.0\right)^{2}\cdot 5≒3.756
B薬 5 18.4  \left(18.4-18.1-19.1+18.0\right)^{2}\cdot 5≒3.472
B薬 5 17.8  \left(17.8-18.1-16.9+18.0\right)^{2}\cdot 5≒3.472
C薬 5 18.2  \left(18.2-17.1-19.1+18.0\right)^{2}\cdot 5≒0.006
C薬 5 16.0  \left(16.0-17.1-16.9+18.0\right)^{2}\cdot 5≒0.006

各表の最右列のズレを足し合わせて、各要因の水準間平方和と交互作用の平方和を計算します。
\quad S_{A}=6.4+0.1+8.1=14.6
\quad S_{B}=19.267+19.267≒38.533
\quad S_{AB}=3.756+3.756+3.472+3.472+0.006+0.006=14.467

また、今回は薬剤が3水準(A薬、B薬、C薬)、前処置が2水準(有、無)のため自由度は次のとおりです。
\quad f_{A}=3-1=2
\quad f_{B}=2-1=1
\quad f_{AB}=3\times2-3-2+1=2

 

次に、個々の観測値について、各水準の平均からのズレ(群内の偏差平方)と総平均からのズレ(全体の偏差平方)を求めます。
※  S_{T} = S_{A} + S_{B} + S_{AB} + S_{e} なので  S_{T} S_{e} の1つ分かれば十分ですが、計算の理解のため両方求めます。

薬剤 前処置 評価指標 各水準の平均からのズレ
(群内の偏差平方)
総平均からのズレ
全体の偏差平方
A薬 20 (20 - 20.8)^{2}=0.64 (20-18)^{2}=4
A薬 21 (21 - 20.8)^{2}=0.04 (21-18)^{2}=9
A薬 22 (22 - 20.8)^{2}=1.44 (22-18)^{2}=16
A薬 20 (20 - 20.8)^{2}=0.64 (20-18)^{2}=4
A薬 21 (21 - 20.8)^{2}=0.04 (21-18)^{2}=9
A薬 18 (18 - 16.8)^{2}=1.44 (18-18)^{2}=0
A薬 16 (16 - 16.8)^{2}=0.64 (16-18)^{2}=4
A薬 17 (17 - 16.8)^{2}=0.04 (17-18)^{2}=1
A薬 16 (16 - 16.8)^{2}=0.64 (16-18)^{2}=4
A薬 17 (17 - 16.8)^{2}=0.04 (17-18)^{2}=1
B薬 19 (19 - 18.4)^{2}=0.36 (19-18)^{2}=1
B薬 18 (18 - 18.4)^{2}=0.16 (18-18)^{2}=0
B薬 17 (17 - 18.4)^{2}=1.96 (17-18)^{2}=1
B薬 19 (19 - 18.4)^{2}=0.36 (19-18)^{2}=1
B薬 19 (19 - 18.4)^{2}=0.36 (19-18)^{2}=1
B薬 19 (19 - 17.8)^{2}=1.44 (19-18)^{2}=1
B薬 18 (18 - 17.8)^{2}=0.04 (18-18)^{2}=0
B薬 18 (18 - 17.8)^{2}=0.04 (18-18)^{2}=0
B薬 17 (17 - 17.8)^{2}=0.64 (17-18)^{2}=1
B薬 17 (17 - 17.8)^{2}=0.64 (17-18)^{2}=1
C薬 20 (20 - 18.2)^{2}=3.24 (20-18)^{2}=4
C薬 18 (18 - 18.2)^{2}=0.04 (18-18)^{2}=0
C薬 16 (16 - 18.2)^{2}=4.84 (16-18)^{2}=4
C薬 19 (19 - 18.2)^{2}=0.64 (19-18)^{2}=1
C薬 18 (18 - 18.2)^{2}=0.04 (18-18)^{2}=0
C薬 15 (15 - 16)^{2}=1 (15-18)^{2}=9
C薬 14 (14 - 16)^{2}=4 (14-18)^{2}=16
C薬 17 (17 - 16)^{2}=1 (17-18)^{2}=1
C薬 16 (16 - 16)^{2}=0 (16-18)^{2}=4
C薬 18 (18 - 16)^{2}=4 (18-18)^{2}=0

このうち、各水準の平均からのズレ(群内の偏差平方)を足し合わせたものが残差平方和で、 総平均からのズレ(全体の偏差平方)を足し合わせたものが全体平方和です。
 \quad S_{e}≒0.64 + 0.04 + ... + 0 + 4≒30.4
 \quad S_{T}≒4+9+...+4+0≒98

また、全体のN数が30、薬剤が3水準(A薬、B薬、C薬)、前処置が2水準(有、無)のため自由度は次のとおりです。
\quad f_{e}=n-a\cdot b=30-3\cdot2=24
\quad f_{T}=n-1=30-1=29

 

最後に、これら平方和・自由度を分散分析表にまとめて F値 と p値 を算出します。

変動因 平方和 自由度 平均平方 F値(F比) p値
薬剤 S_{A}=14.6 f_{A}=2 V_{A}=7.3 F_{A}≒5.7632 0.0090
前処置 S_{B}≒38.533 f_{B}=1 V_{B}≒38.533 F_{B}≒30.4211 0.0000
交互作用 S_{AB}≒14.467 f_{AB}=2 V_{AB}≒7.233 F_{AB}≒5.7105 0.0094
残差e S_{e}=30.4 f_{e}=24 V_{e}≒1.267    
全体T S_{T}=98 f_{T}=29      

分散分析表のとおり、今回の架空データでは薬剤と前処置の交互作用が p<0.01 で有意であり、これら2要因の間には交互作用があるようです(何れかの組み合わせで相乗効果か相殺効果がある)。

なお、薬剤の主効果は p<0.01、前処置の主効果は p<0.001 でともに統計学的に有意でしたが、交互作用があるため参考程度にとどめておきます(独立していないため、薬剤と前処置を切り離して考えない方が良い)。

交互作用プロットを見ると、薬剤毎の前処置有無の平均の差が一定ではないため、交互作用が有意になったと考えられます。特にB薬は前処置に関する効果が薄いように見えます。 

f:id:cochineal19:20210220141315p:plain

※変数名:TRT01AN=薬剤、水準:1=A薬、2=B薬、3=C薬
 変数名:PRCFL=前処置有無、水準:0=無、1=有
 変数名:AVAL=評価指標X

 

二元配置分散分析の結果の解釈と進め方

今回の架空データでは交互作用ありでしたが、二元配置分散分析の結果によって次のアクションに違いがあるため、備忘録がてらメモします。
※多重比較を事後検定と認識するのはどうだ、交互作用があったら単純主効果検定など色々な見解を見かけましたが、こちらのサイト様 が分かりやすかったので、参考にさせていただきました。

1.交互作用が認められる
 ⇒ 薬剤と前処置が互いに影響を与えていると結論。
   必要であれば薬剤×前処置の6つの組み合わせで多重比較。
   各要因を切り離して考えない。

2.交互作用が認められない&薬剤の主効果が有意
 ⇒ 薬剤に有意な変動があったと結論。
   薬剤が3水準なので、必要であれば多重比較。

3. 交互作用が認められない&前処置の主効果が有意

 ⇒ 前処置に有意な変動があったと結論。
   2水準なので多重比較する必要なし。

4.交互作用が認められない&薬剤・前処置ともに主効果が有意

 ⇒ 薬剤・前処置にともに有意な変動があったと結論。
   加えて、薬剤×前処置に交互作用がなく互いに影響を与えていないと判断。
   必要であれば、薬剤が3水準なので、薬剤について多重比較。
   前処置は2水準なので完結。

 

Rでの実行:

> library(car) 
> ADS <- data.frame(
+   TRT01AN=c(0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1,2,2,2,2,2,2,2,2,2,2)
+   ,PRCFL=c(0,0,0,0,0,1,1,1,1,1,0,0,0,0,0,1,1,1,1,1,0,0,0,0,0,1,1,1,1,1)
+   ,AVAL=c(20,21,22,20,21,18,16,17,16,17,19,18,17,19,19,19,18,18,17,17,20,18,16,19,18,15,14,17,16,18))
> ADS$前処置 <- factor(ifelse(ADS$PRCFL==0,"無","有"),levels=c("無","有"))
> ADS$投与群 <- factor(ifelse(ADS$TRT01AN==0,"A薬",ifelse(ADS$TRT01AN==1,"B薬","C薬")),levels=c("A薬","B薬","C薬"))
> fit1 <- lm(AVAL ~ 投与群 + 前処置 + 投与群:前処置, data=ADS)
> car::Anova(fit1, Type="II")
Anova Table (Type II tests)

Response: AVAL
              Sum Sq Df F value    Pr(>F)    
投与群        14.600  2  5.7632  0.009035 ** 
前処置        38.533  1 30.4211 1.134e-05 ***
投与群:前処置 14.467  2  5.7105  0.009363 ** 
Residuals     30.400 24                      
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

SASでの実行:

data ADS;
do TRT01AN = 1; do PRCFL = 0; do AVAL = 20,21,22,20,21; output; end; end; end;
do TRT01AN = 1; do PRCFL = 1; do AVAL = 18,16,17,16,17; output; end; end; end;
do TRT01AN = 2; do PRCFL = 0; do AVAL = 19,18,17,19,19; output; end; end; end;
do TRT01AN = 2; do PRCFL = 1; do AVAL = 19,18,18,17,17; output; end; end; end;
do TRT01AN = 3; do PRCFL = 0; do AVAL = 20,18,16,19,18; output; end; end; end;
do TRT01AN = 3; do PRCFL = 1; do AVAL = 15,14,17,16,18; output; end; end; end;
run;
proc glm data=ADS;
class TRT01AN PRCFL;
model AVAL = TRT01AN PRCFL TRT01AN*PRCFL / ss1 ss2 ss3 e1 e2 e3 solution;
run;

f:id:cochineal19:20210220140739p:plain

f:id:cochineal19:20210220140807p:plain

プログラムコード 


Rのコード

library(car)
fit1 <- lm(AVAL ~ FACTOR1 + FACTOR2 + FACTOR1:FACTOR2, data=ADS) #-- FACTOR1:FACTOR2は交互作用項
car::Anova(fit1, Type="II") #-- Type2平方和

SASのコード

proc glm data=ADS1 outstat=OUTDS;
  class FACTOR1; /*class指定必須*/
  model AVAL = FACTOR1 FACTOR2 FACTOR1*FACTOR2 / ss1 ss2 ss3 e1 e2 e3 solution;
run;
/* option: ss1=Type1平方和, ss2=Type2平方和, ss3=Type3平方和 */
/* option: e#= Type#平方和のEstimable Functions(推定可能関数)を図示してくれる */
/* option: solution=回帰モデルの推定値を出力してくれる */

Pythonのコード

整備中

 

雑談


交互作用有無のイメージです。

f:id:cochineal19:20210220145954p:plain

 ⇒ 前処置の効果がA~C薬どこでも一定であれば交互作用なし。

 

f:id:cochineal19:20210220152536p:plain

 ⇒ 前処置の効果が一定でなければ交互作用あり。
   薬剤と前処置の組み合わせでより効果が出たり(相乗効果)、
   効果が薄まったり(相殺効果)する。

 

参考


統計学基礎(出版:東京図書), 日本統計学会編

統計学(出版:東京図書), 日本統計学会編

改訂増補版:統計検定を理解せずに使っている人のためにIII

  


サイトマップ

cochineal19.hatenablog.com

本ブログは個人メモです。 本ブログの内容によって生じた損害等の一切の責任を負いかねますのでご了承ください。