(HCI研究) 2要因のデータに行う分散分析の正規性の仮定について

:::note alert
記事の内容に不適切な箇所がありました。訂正箇所とコメント欄をご参照ください。
:::

# この記事について
この記事はHCI分野の研究室に所属している私が学部の研究のときに疑問に思ったことを調べ、メモとして残した記事です。
決まり文句ですが、私は統計学の専門家ではありません。記事の内容の使用についてはご自身の責任でお願いします。また、記事について何か不適切な点がありましたらご指摘頂けると幸いです。

# HCI研究における統計手法
HCI(ヒューマンコンピュータインタラクション)の分野では行った実験で得られたデータに対し、実験計画とデータの性質に合わせた統計手法を正しく適用して分析する必要があります。

ネット上においては以下のサイトがよく統計手法をまとめてくださっており、私も参考にさせて頂いています、ありがとうございます。

https://zenn.dev/tmizuho/books/3d511e017bfd23/viewer/810182

また、Jacob O. Wobbrock 先生のPractical Statistics for HCIの資料も参考になるかと思います。

https://depts.washington.edu/acelab/proj/ps4hci/index.html

こういった資料に助けていただき、学部の研究では統計手法の選別を楽に行うことが出来ました。

しかし、それでも手法の選別に迷う場面がありました。この記事ではその場面の説明、そしてこの疑問点に対して私が調べて得た知見を記載しています。

# 2要因のデータの分析手順
2要因のデータが得られた場合、~~まずそのデータに対する正規性の検定が必要になります~~

:::note alert
検定手法の選択をシャピロ・ウィルク検定、コルモゴロフ・スミルノフ検定のような正規性の検定"のみ"で行うことは不適切です。
本記事はこの部分について正しく認識せず、「正規性検定を行なってからパラメトリック検定・ノンパラメトリック検定を選択する」という手順を"安易"に踏んでいますが、正規性検定の結果の信ぴょう性については様々な議論があることを認識していたくようにお願いします。
https://best-biostatistics.com/summary/is-normality-testing-necessary.html
https://biolab.sakura.ne.jp/normality-test-nonparametric.html

また、正規性の検定の結果の疑わしさとは別に、そもそも正規性検定の結果で帰無仮説を棄却できない場合「正規分布に従っているかどうか」は何とも言えないことにも注意してください。
https://qiita.com/goodclues/items/0c91f6ac6df6a081fe86

コメント欄もご参照ください。
:::

<details><summary>記事の続きを読む場合は以上のことを理解した上でお願いします。</summary>
そして得られたデータに正規性および等分散性が仮定できるならANOVAを、正規性が仮定できない場合は整列ランク変換(Aligned rank transform, ART)を施す必要があります。

私の研究においてはこの手順を実施することになりました。
(実は私の場合3要因なのですがややこしいので無視します。ちなみに、今考えると無闇に3要因以上の実験計画を立てたのはあまりよろしくなかったと思います[^6]。）

# 2要因のデータに対して正規性の検定をどう行えばいいのか
しかし、いざやってみると手が止まってしまいました。というのも、「正規性の検定を行う」という部分について詳細な手順の情報がなかったのです。例えば全体の分布の正規性だけを見ればいいのか、各要因ごとの分布それぞれについて確かめるのか、各要因の水準の組み合わせごとの分布について確かめればいいのか…という詳細について確たる情報がありませんでした。

一応、「各群(各要因の水準の組み合わせ)ごとに正規性の検定を行う」という手順が正しいのではという意見が見られましたが、これも明確な根拠があるわけではありませんでした。

というわけで、記事の最初に述べさせていただいたように、この点について自分で調べてました。

この点について、少ない情報の中でも以下のような指摘が見られます。

> 要するに，多要因，多群の分散分析の正規性は，その残差の正規性を調べれば良いのだが，意外と知らない人が多いのである。
(中略)
では，群ごとや，要因の水準ごとに正規性を検定すると，何が問題なのか？ごく簡単なことで，検定の多重性が問題になるのである。

https://biolab.sakura.ne.jp/normality-test-for-anova.html

こちらの記事では群ごとの正規性の検定は多重性の問題から適切ではないこと、その代わりに残差に対して正規性を検定すべきだと述べられています。

実際に、ANOVAの正規性の仮定は残差の仮定だという文言は確認できます[^1][^2][^3][^4][^7]。

ただ、当時はHCIの分野でこの残差の正規性の検定についての言及が見られず、確信を得られませんでした。
そこで「各群ごとに正規性の検定を行う」と「残差について正規性の検定を行う」という2つの手順の妥当性について、自分でシミュレーションして確かめてみました。上記の記事においてもシミュレーションがなされていますが、少しデータ数を増やし、できるだけHCIの実際の実験データに寄せてシミュレーションを行いました。その結果を下記に記載します。

# 2要因参加者間デザインでのシミュレーション
記事を参考に2要因参加者間デザイン(ABsデザイン)、各要因の水準をA1B1, A1B2, A2B1, A2B2としてシミュレーションを行います。ただし、それぞれの水準の大きさは100とし、全体でn=400としました。これはHCIの実験では1つの独立変数の組み合わせ当り100データぐらいはとってくるだろうという私の感覚による調整です。

以下のRのコードによってシミュレーションを行いまいた。ほぼ参考にさせていただいた記事そのままです。
標準正規母集団からデータを取り出し、各群ごとにShapiro-Wilk 正規性検定、残差についてShapiro-Wilk 正規性検定をおこなうという操作を100000回繰り返し、どのぐらいの割合で正規分布でないとの判定をしてしまうかを出力しています。

```R
np <- npres<- 0
k <- 100000
for (i in 1:k) {  

  participant <- factor(rep(1:10, each=40))
  A <- factor(rep(1:2, each=20, times=10))
  B <- factor(rep(1:2, each=10, times=20))
  time <- rnorm(400)

  data <- data.frame(participant, A, B, time)


 # 分散分析
  mod<- aov(time ~ A * B)

 # 残差
  resid<- mod$residuals

  # 各要因の水準組み合わせの正規性検定
  p.A1B1<- shapiro.test(subset(data, A==1 & B==1)$time)$p.value
  p.A1B2<- shapiro.test(subset(data, A==1 & B==2)$time)$p.value
  p.A2B1<- shapiro.test(subset(data, A==2 & B==1)$time)$p.value
  p.A2B2<- shapiro.test(subset(data, A==1 & B==2)$time)$p.value

  if (p.A1B1 < 0.05 || p.A1B2 < 0.05 || p.A2B1 < 0.05 || p.A2B2 < 0.05) {
   np <- np + 1
  }

  # 残差の正規性検定
  p.res<- shapiro.test(resid)$p.value

  if (p.res < 0.05) {
   npres <- npres + 1
  }

}

# 5% 水準での棄却率
# 要因別の正規性検定
npr <- np/k

# 分散分析の残差検定
npresr <- npres/k


sprintf("npr= %f", npr)
sprintf("npresr = %f", npresr)
```
出力結果
```
[1] "npr= 0.140980"
[1] "npresr = 0.050650"
```
結果より、100データを集めたとしても2つの正規性の検定の精度に差があり、残差の正規性の検定の精度の方が高いことが確認できました。

# 2要因参加者内デザインでのシミュレーション
同様のシミュレーションを2要因参加者内デザイン(sABデザイン)でも行いました。

```R
np <- npres<- 0
k <- 100000
for (i in 1:k) {  

  y<- rnorm(40)

  participant <- factor(rep(1:10, each=40))
  A <- factor(rep(1:2, each=20, times=10))
  B <- factor(rep(1:2, each=10, times=20))
  time <- rnorm(400)

  data <- data.frame(participant, A, B, time)


 # 分散分析
  mod <- aov(time~A*B +Error(participant+participant:A+participant:B+participant:A:B), data)

 # 残差
  resid<- mod$Within$residuals

  # 各要因の水準組み合わせの正規性検定
  p.A1B1<- shapiro.test(subset(data, A==1 & B==1)$time)$p.value
  p.A1B2<- shapiro.test(subset(data, A==1 & B==2)$time)$p.value
  p.A2B1<- shapiro.test(subset(data, A==2 & B==1)$time)$p.value
  p.A2B2<- shapiro.test(subset(data, A==1 & B==2)$time)$p.value

  if (p.A1B1 < 0.05 || p.A1B2 < 0.05 || p.A2B1 < 0.05 || p.A2B2 < 0.05) {
   np <- np + 1
  }

  # 残差の正規性検定
  p.res<- shapiro.test(resid)$p.value

  if (p.res < 0.05) {
   npres <- npres + 1
  }

}

# 5% 水準での棄却率
# 要因別の正規性検定
npr <- np/k

# 分散分析の残差検定
npresr <- npres/k

sprintf("npr= %f", npr)
sprintf("npresr = %f", npresr)
```
---------------------------------------
コードの下記部分について、
```R
 # 残差
  resid<- mod$Within$residuals
```
対応のある二元配置分散分析だとaovの結果は因子ごと、および各因子の組合せごと
今回の場合、
* Stratum 1: participant
* Stratum 2: participant:A
* Stratum 3: participant:B
* Stratum 4: participant:A:B
* Stratum 5: Within

に分けて出力されます。このうちどの残差が正規性を検定すべき残差なのか、確かな情報を得ることはできませんでした。

ある記事がone-way repeated-measure ANOVAの出力結果におけるWithinの項目を errors (Residuals)と述べていたこと[^5]、また下記のpythonのコードのResidualの値に一致するのが「Stratum 5: Within」の値であったことから、今回は「Stratum 5: Within」のResidualに対し正規性検定を行っています。

<details><summary>pythonコード</summary>

```python
# ANOVAモデルを設定
model = ols("""time ~ C(A) + C(B) +
               C(A):C(B) + C(participant):C(A):C(B)""", data=data).fit()

# ANOVAを実行
anova_table = sm.stats.anova_lm(model, typ=1)
print("ANOVA Table:")
print(anova_table)
```
</details>


また下記記事を参考にすると、aovの出力結果は順に
* 全体のズレ
* 要因1によるズレ
* 要因2によるズレ
* 交互作用によるズレ
* 残りのズレ（残差）

と出力しているような気がします。

https://kogolab.chillout.jp/elearn/hamburger/chap7/sec1.html



根拠が曖昧で申し訳ありませんが、この点について間違いがあればご指摘して頂けると幸いです。

---------------------------------------
Rの出力結果
```
[1] "npr= 0.142270"
[1] "npresr = 0.050990"
```
結果は参加者内計画においてもやはり2つの正規性の検定の精度に差があり、残差の正規性の検定の精度の方が高いことが確認できます。

~~以上の結果から、HCI研究における2要因のデータの正規性の検定を行う際には、残差についての正規性検定を行うことが有用である可能性が高いと分かりました。~~

# まとめ
2要因のデータの正規性を検定について、各群ごとに正規性の検定を行うことは不適切であり、残差に対して正規性の検定を行うべきであるという指摘に対し、それをHCI分野の実験に近いデータでそのことを検証しました。

その結果、やはり残差に対する正規性の検定の方が精度の高い結果を示し、HCI研究においてもこの手順を用いることの適切さが示唆されました。この結果は2要因のデータに限らず、多群・多要因のデータに関わることだと考えられます[^8]。

しかし、だからといって今後の研究で残差の正規性だけを確認しようとするのは、~~恐らく正しくはあるのですが~~、良い考え方ではないのかなとも思います。
最初に紹介させて頂いたJacob O. Wobbrock 先生のPractical Statistics for HCIでは以下のように述べられていました。
> Many researchers use two (or more) tools when analyzing data. Doing so helps to compare results to ensure correctness.

つまり、一流の研究者でも色々な統計手法を試してみて分析結果の妥当性を確認しているのだと思います。

このような考え方も大事にして、~~残差に対する正規性の検定をあくまで確からしい統計手順の1つだと認識し、様々な確からしい統計を行いながら~~、今後の研究の分析結果の妥当性を高めいこうと今は考えています。


</details>

# 備考
以下の記事・資料もHCIの統計手法の理解に役立つと思います。

https://yatani.jp/teaching/doku.php?id=hcistats:start_jp

https://www.stat.cmu.edu/~hseltman/309/Book/

また私が見逃しているだけで、今回の多群・多要因のデータの~~正規性検定~~正規性の確認の手順についてHCIの分野で何か言及があった場合はご教示いただけると幸いです。



[^1]:[Wikipedia - Analysis of variance, Assumptions](https://en.wikipedia.org/wiki/Analysis_of_variance)
[^2]:[Wikipedia - Two-way analysis of variance, Assumptions](https://en.wikipedia.org/wiki/Two-way_analysis_of_variance)
[^3]:[Practical Statistics for HCI, ps4hci.key, p81](https://depts.washington.edu/acelab/proj/ps4hci/index.html)
[^4]:[Two-Way ANOVA - Full Course(Youtube)](https://youtu.be/7OY0DdFyJBg?si=zcMOj233-3ncaDov)
[^5]:[HCI研究者向け統計手法 - Koji Yatani / 矢谷浩司, Analysis of Variance (ANOVA) for comparing the means](https://yatani.jp/teaching/doku.php?id=hcistats:anova)
[^6]:[*So, I recommend you to design your experiment so that you don't have to use a more complicated test than two-way ANOVA.* (HCI研究者向け統計手法 - Koji Yatani / 矢谷浩司)](https://yatani.jp/teaching/doku.php?id=hcistats:anova)
[^7]:[一元配置分散分析（ANOVA）の残差プロット - Support - Minitab](https://support.minitab.com/ja-jp/minitab/help-and-how-to/statistical-modeling/anova/how-to/one-way-anova/before-you-start/data-considerations/)
[^8]:[t 検定の正規性も残差を調べる：検定の多重性問題](https://biolab.sakura.ne.jp/normality-for-t-test.html)