生物医学統計を学ぶ！: 分散分析（analysis of variance, ANOVA）

2群の差の検定と異なり、3群以上を対象に差の検定を行います。

一般に、3群以上の群は多群とも呼べます。
分散分析は、データのばらつきを因子によるものと誤差によるものに分解し、その大きさを比較することで群間に差があるかどうかを分析する手法です。
また、分散分析からさらに詳細な分析を行うための多重比較法に移行するパイプラインの役割を果たします。

多群の検定では、「要因」を軸にして分析できるようになります。
順に説明していきます。

一元配置分散分析

検定の目的

1元配置分散分析（One-factor ANOVA, One-way ANOVA）とは、3群以上のデータ（変数は群につき1つ）に対して行う差の検定です。
この方法はすべての群の母平均（各群の平均）が等しいという帰無仮説を検定します。
一元配置分散分析の結果が有意であったときの結論は、「平均値の異なる群が少なくとも1つある」ということになります。
一元配置分散分析はパラメトリックな検定に属します
これとは別に、一元配置分散分析でいうところの代表値を中央値の差を用いて検定するノンパラメトリックな検定にはクラスカル・ワリス検定などがあります。

データの構造

データの要因を把握する必要があります。

一元配置分散分析に利用するデータは、1要因・多群のデータです。

例えば、健常群、疾患A群、疾患B群の3群があり、それぞれ空腹時血糖値のデータがあるとします。

この例では、要因は空腹時血糖値（変数）です。

フローチャート

フローチャートは次に示す通りです。

フローチャート中のウェルチの補正は、一元配置分散分析にウェルチの補正を加えた手法です。

フローチャート中、多重比較法という見慣れない手法もありますが、別途ご紹介致します。

実践

1因子（1変数）について多群間での違いを検討するために、一元配置分散分析を適用します。

ある臨床研究で、A薬、B薬、C薬、 D薬について、対応のないデータそれぞれで10回分の血糖値検査結果を得ました。

帰無仮説は「これらの薬剤間で、薬剤が血糖値に及ぼす影響に差はない」です。

```
  # pingouinパッケージを利用します
  !pip install --upgrade pingouin
  
  # データをダウンロードします
  # https://drive.google.com/file/d/19Q9lujUGMskIQ10febSIFWj04nQ5F17B/view?usp=sharing
  !curl -L -o "One-Way_ANOVA.csv" "https://drive.google.com/uc?export=download&id=19Q9lujUGMskIQ10febSIFWj04nQ5F17B"
  
  # データを読み込みます
  import pandas as pd
  data = pd.read_csv(filepath_or_buffer="One-Way_ANOVA.csv", names=["group","values"], encoding="utf-8")
  print(data.head())
  data.describe()
  
  # 標準化
  bs = (data["values"] - data["values"].values.mean())/np.std(data["values"]) # blood sugar
  data = pd.DataFrame(np.stack([bs,data["group"].values],axis=1),columns=["values","group"])
  data["values"] = data["values"].astype(float)
  
  # 正規性を確認します。
  from scipy import stats
  # データを薬剤ごとに分けて、検定します。
  for name in ["A","B","C","D"]:
    one_group = data[data["group"] == name]
    stat,p = stats.shapiro(one_group["values"])
    print("norm dist test",name,stat,p) # 正規性が仮定できることが分かります
  
  # 次に等分散性を調べます。
  stat_l,p_l = stats.levene(data[data["group"] == "A"]["values"],data[data["group"] == "B"]["values"],data[data["group"] == "C"]["values"],data[data["group"] == "D"]["values"], center="mean")
  print("levene test",stat_l,p_l) # 等分散が仮定できることが分かります
  
  # 正規性+、等分散性+なので、一元配置分散分析を行います
  import pingouin as pg
  aov = pg.anova(dv='values', between='group', data=data, detailed=True)
  print(aov)
  # 0行目は一元配置モデルから算出した値
  # 1行目は誤差です。
  # 'Source' : 因子
  # 'SS' : Sums of squares：平方和
  # 'DF' : Degrees of freedom：自由度(n-1)
  # 'MS' : Mean squares：平均平方（平方和を自由度で除したもの）
  # 'F' : F-values：F統計量
  # 'p-unc' : uncorrected p-values：非補正p値
  # 'np2' : Partial eta-square effect sizes

# p-uncは非常に小さい値(0.0001未満)になりましたので、有意水準5％で有意差が認められました。
  # 何かの群が同等でないことがわかったということです。
  #ここからさらにどの群がどの程度違うのかを調べていくには、多重比較法であるテュ―キーの方法などを用います。
  
  # 補足１
  # 先の例では、正規性が認められましたが、認められなかったパターンも手順のみ確認しておきます。
  # （このデータでは）実際には行いませんが、同じデータで試してみましょう。
  # フローチャートから、正規性が仮定できない場合は、「クラスカル・ウォリス検定」を使います。
  from scipy import stats
  stat,p = stats.kruskal(data[data["group"] == "A"]["values"],data[data["group"] == "B"]["values"],data[data["group"] == "C"]["values"],data[data["group"] == "D"]["values"])
  print("kruskal test",stat,p)
  
  # 補足２
  # このデータは正規性（+）、等分散（+）ですが、
  # 等分散が仮定できなかったパターンを手順の確認の意味で実行します。
  # 等分散が仮定できない場合は、ウェルチ補正ありの一元配置分散分析を用います。
  from pingouin import welch_anova
  waov = welch_anova(dv='values', between='group', data=data)
  print(waov)
  # 'Source' : 因子
  # 'ddof1' : 自由度1
  # 'ddof2' : 自由度2
  # 'F' : F-values：F統計量
  # 'p-unc' : uncorrected p-values：非補正p値
```

反復測定による分散分析

検定の目的

反復測定一元配置分散分析（One factor repeated measures ANOVA, One way Repeated Measure ANOVA, RM-ANOVA）は対応のある3群以上の群間の差を検定する手法です。

繰り返し測定する項目が要因（変数）となります。

この方法はすべての群の母平均（群の平均）が等しいという帰無仮説を検定します。

検定の結果が有意であったときの結論は、「平均値の異なる群が少なくとも1つある」ということになります。

データの構造

例えば、Aさん、Bさん、Cさん3人の対象者に対して、一週間おきに問診を3回繰り返して取得したデータなどを利用します。

フローチャート

フローチャートは次のようになります。

実践

あるプラセボ対照試験において、被験者ごとに0週、1週、3週、6週の検査値（1変数）が得られている。

実薬群において、時点間の差を考慮したい。

要因は経過時間ごとの検査値です。

帰無仮説は「時点間の差はない（時点間によって差があるものはない）」です。

```
  # データ読み込み
  # https://drive.google.com/file/d/1BoiFdcooyDnr5a4BdATUGSMHrsoHQtD6/view?usp=sharing
  !curl -L -o "RM_ANOVA.csv" "https://drive.google.com/uc?export=download&id=1BoiFdcooyDnr5a4BdATUGSMHrsoHQtD6"
  import pandas as pd
  data = pd.read_csv(filepath_or_buffer="RM_ANOVA.csv", header=0,encoding="utf-8")[["ID","Group","W0","W1","W3","W6"]]
  print(data.head())
  
  # フローチャートに従って進めます。
  # まずは、正規性の検定です
  # このデータのIDとGroupはカテゴリ変数です。（解析時にはデータから外します）
  # Groupの1は実薬群、0はプラセボです。
  # Repeated measure ANOVAは、対応のあるデータを用いるので、今回はGroup1のみを利用します。
  # データをGroup1のみにします。
  control = data[data["Group"] == 1]
  control = control.reset_index().drop(["index"],axis=1)
  print("データを実薬群にしました")
  print(control.head())

# 各変数の正規性を確認します
  from scipy import stats
  for s in ["W0","W1","W3","W6"]:
    X = control[s]
    stat,p = stats.shapiro(X)
    print(s,stat,p)
  
  # このデータでは正規性が仮定できませんでした。
  # よって、フリードマンの検定を行います。
  stat, p = stats.friedmanchisquare(control["W0"],control["W1"],control["W3"],control["W6"])
  print("フリードマンの検定",stat,p) # 515.0118881118887 2.662434646149274e-111
  # 今回はp=0.05未満で、有意差（有意水準5%）が認められました。
  # どこかの時点で時点全体との差が大きいようです。

# ここからさらに詳しく調べるには、多重比較法（スティール・デュアスの方法など）を行います。
  # スティール・デュアスの方法は、多重比較法で取り扱います。
    
  # このデータでは正規性が仮定できませんでしたが、手順の確認の意味で、
  # 同じデータを使って正規性が仮定できた場合の操作を確認しておきます。

# 正規性（+）、モークリーの球面性検定
  import pingouin as pg
  s,w,chi,ddof,p = pg.sphericity(control[["W0","W1","W3","W6"]])
  print(s,w,chi,ddof,p)
  # spher       boolean True if data have the sphericity property.
  # Wfloat      Test statistic.
  # chi_sqfloat Chi-square statistic.
  # ddofint     Degrees of freedom.
  # pfloat      p-value  
  
  # 球面性が仮定できなかったため、Repeated Measure ANOVAをグリーンハウスカイザーのε修正を付けて実行します。
  import pingouin as pg
  # オプションにcorrection=Booleanで指定すれば、グリーンハウスカイザーのε修正ありの計算になります
  pg.rm_anova(data=control[["W0","W1","W3","W6"]],correction=True)

# 結果、p値は2つあります。
  # p-uncはRM-ANOVAのp値です。球面性が仮定できない場合でもRM-ANOVAのp値（p-unc）で補正のないp値を確認できます。
  # p-GG-corrはRM-ANOVAをグリーンハウスカイザーのε修正ありで補正したp値です。
  # 今回は後者のp値ですから、p=0.05未満で「差がある」という結果になりました。
  ```

N元配置分散分析(N-way ANOVA)

検定の目的

多元配置分散分析（N-way ANOVA）とは、一元配置分散分析（1要因）をN要因に拡張した分散分析です。

一元配置分散分析との大きな違いは、一元配置分散分析では要因は1つでしたので、要因間の作用は考慮できませんでしたが、多元配置の場合は複数の要因を指定するため、要因の差を検定するだけでなく、要因間の交互作用も解析できます。

ただし、交互作用の分析は相関分析とは異なります。

データの構造

例えば、手術後1週で退院、2週で退院、3週で退院の手術後から退院までの期間のカテゴリカルな要因と、手術A、手術Bという手術の種類のカテゴリカルな要因の2要因（独立変数）、および、被験者ごとの値（例えば、炎症の定量値として扱えるCRPなど。この例における従属変数）が存在するようなデータを例に示します。

ID	手術種類	退院カテゴリ	術後3日CRP
1	A	早期退院	2.5
2	A	早期退院	0.5
3	A	満期退院	3.4
4	A	満期退院	0.6
5	B	満期退院	6.7
6	B	退院遅延	10.8

この場合、手術の種類と退院のカテゴリとで要因は2つとなり、2元配置分散分析が適用できます。

データが3要因になるときは3元配置分散分析、4要因の時は4元配置分散分析と、要因の数に沿って N元配置分散分析として扱われます。

N元配置分散分析では、要因を組み合わせて仮説を設定し、あるひとつの変数（従属変数）への影響（これを「効果」という）について検討します。

一般的に、N元配置分散分析では要因の有意差の検定（少なくとも1つの要因に有意差がある（主効果））と、交互作用（interaction）の検定が同時に計算されます。

主効果(main effect)とは、それぞれの要因（N元配置分散分析では変数を指定します）がそれぞれ「独自」に各変数へ与える単純効果のことです。

先の例で説明すると、N元配置分散分析で有意となった場合に、主効果は、手術の種類と退院のカテゴリがそれぞれCRPに与える影響を調べるということになります。

交互作用(interaction)とは、要因（変数）を組み合わせた場合の複合効果のことです。

ある要因Aの主効果とある要因Bの主効果だけでは説明できない、組み合わせにより生じる効果を意味します。

N要因の分散分析では、「主効果があるか」よりも先に、要因の交互作用を検証します。

要因間の交互作用が認められとき、主効果の検定結果を考察しに行くという手順になります。

例えば、要因Aと要因Bの交互作用が有意である時、要因Bのある水準での要因Aの主効果、要因Aのある水準での要因Bの主効果について確認します。

先の例で言い換えると、退院のそれぞれのカテゴリと「手術の種類」の主効果、手術のそれぞれの種類と「退院のカテゴリ」の主効果を確認するという流れです。

交互作用がすべての要因間の組み合わせで認められない場合も、主効果の解析結果を参照します。

主効果が有意である場合には、さらに詳細を調べるために必要に応じて多重比較を行います。

まとめると解析のパターンは次のようになります。

交互作用の有意差あり

主効果の有意差に関わらず全組み合わせの多重比較検定のフローを検討する

交互作用の有意差なし

主効果の有意差あり

要因ごとに多重比較検定

主効果の有意差なし

解析を終了する（あるいは、念のため要因ごとに多重比較検定を行う）

あくまで、一元配置分散分析が基本となっていますから、主効果の有意差は「要因間の差がない」を帰無仮説とします。

主効果の対立仮説は少なくとも1つの要因に有意差があるです。

交互作用はある要因とある要因と相性のような意味もあります。

例えば、要因Aと要因Bがそれぞれ全く同じようなデータだった場合、あるいは、要因Aと要因Bの間のなんの法則性もない場合は交互作用はないと考えられます

一方で、要因Aが要因Bの相乗効果を生じさせる、その反対に、相殺効果を生じさせるようなケースでは交互作用があると言えます

多重比較法以外で交互作用のさらなる考察のためには、データの要約による背景因子の把握、相関分析、交絡を考慮した共分散分析などが手段として挙げられます。

実践

pingoinパッケージのサンプルデータ「anova3.csv」を使って、利用方法を確認します。

anova3は、sex、risk、drug、cholesterolのデータセットです。

sex、risk、drugはカテゴリカル変数、cholesterolは実数の変数です。

上から順に主効果、交互効果となっています。

さきに交互効果のp値を確認します。

交互効果はすべてp≧0.05で、有意差は認められません。

次に、主効果のp値を確認します。

主効果のうち、「Risk」のみ有意差（p<0.05）が認められます。

この結果から、より詳細な検定は「Risk」についてまとめたデータについて行えるということがわかります。

反復測定による二元配置分散分析（2要因）

検定の目的

反復測定2元配置分散分析（Two- factor repeated measure ANOVA, Two- way repeated measure ANOVA）は繰り返し取得された対応のあるデータを対象に、2要因（カテゴリカルな独立変数）で多群を検定する手法です。

解析の手順は前述した反復測定による1要因の分散分析と同様です。

ノンパラメトリック検定は存在しないのでそのまま実行することができます。

1要因、2要因に関わらず、反復測定の分散分析では対象者（あるいは何かのサンプル）が繰り返しデータを取られます。

データの構造

反復測定による2要因の分散分析は次のような対応のあるデータです。

ID	手術形式	性別	術後W0	術後W1	術後W3...
100	A	M	10	9	8
200	A	F	15	14	13
100	B	F	30	29	28
200	B	M	35	34	33

対応のあるデータですので、同一個体に対して、繰り返し取得しているデータです。

この場合、要因は手術形式と性別です。

N元配置分散分析と同様に要因間の差を検定します

帰無仮説は「要因間に差がない」です。

実践

pingouinパッケージのサンプルデータ「rm_anova2」を利用します。

例として、要因は「性別」と「恐ろしさ」、従属変数は犯罪の起こしやすさの定量値を指定します。

分割プロットデザインによる分散分析

検定の目的

分割プロットデザインによる分散分析（ Split Plot Design ANOVA）は、対応のない要因と、反復測定要因（対応のある要因）の2要因以上を組み合わせ、同時に差の検定を行いたいときに適用する分散分析です。

帰無仮説は「要因間に差がない」です。

データの構造

例えば、対応のないデータ（要因）として手術形式（独立変数）、対応のある要因として反復的に術後に調べた何かの検査値（従属変数）がある場合、このようなデータになります。

ID	手術形式	時期	検査値
1	A	W0	9
1	A	W1	8
2	B	W0	29
2	B	W1	34

実践

あるプラセボ対照試験において、被験者ごとに0週、1週、3週、6週の検査値（1変数）が得られている。

各群において、群間差と時点間の差を考慮したいとします。

Visionary Imaging Services, Inc.