source("https://riseki.cloudfree.jp/?plugin=attach&refer=ANOVA%E5%90%9B&openfile=anovakun_489.txt")9 多群の平均値差の検定
心理学実験においては，古典的に分散分析モデルが多用されてきた。平均値差を見ることで実験の効果，因果関係を明らかにできるように，巧妙に実験デザインが組み立てられる。その精緻さは理論的一貫性という意味である種の美しさを持ち，多くの研究者が魅了されてきた。
分散分析にのせることを目的にした実験計画であり，実験デザインの不自由さ(とにかく分散分析をしなければならない!)が批判的に論じられることもあるが，心理学が測定している対象が平均値差以上の精度で議論できる性質でないという反論もあるだろう。
今や分散分析を超えたより高度な統計モデルがあり，現在の研究においては分散分析はもはや過去のものにすぎないかもしれないが，以後のモデルも分散分析を基本としたその発展系であるので，改めて基本を押さえておくことも重要である。
9.1 分散分析の基礎
分散分析は「分散」の分析であるかのような名称であるが，平均値差を検定するためのものである。なぜ「分散」を冠するかといえば，効果量のところで見たように，平均値差の判断には群内分散の情報が必要だからである。
多群の平均値の差，その散らばりを群間分散といい，群に含まれるデータの散らばりを群内分散とよぶ。分散分析は群内分散に対する群間分散の比が十分に大きいと考えられる場合，群間に統計的な有意差があると判断する。分散の比を表す確率分布はF分布と呼ばれる。F分布は群間・群内それぞれの自由度を母数にもつ。
また，実験計画はBetweenデザインとWithinデザインに区分される。t検定でみたような，対応のない独立した群を対象にしたデザインがBetween，群間に相関が想定される対応のある群を対象にしたデザインがWithinである。Withinデザインは同じ個体から複数回の反応を得る(ex. period 1-2-3…)ため，反復測定デザインRepeated measured designともよばれることがある。この場合，群内分散から個人内の分散すなわち個人差を取り出すことができるため，これが分離できないBetweenデザインよりも基本的にWithinデザインのほうが目的となる変動を捉えやすい。ただし，反復測定による個体への負担を考えると，毎回Withinデザインでいいというわけにもいかないところが難点である。
\[Between Design: \text{全変動} = \text{群間変動}+ \text{群内変動(誤差)}\] \[Within Design: \text{全変動} = \text{群間変動}+ \text{個人差変動} + \text{誤差}\]
分散分析は要因が複数ある場合も考えられるから，要因AがBetween，要因BがWithinといった場合は混合計画と呼ばれることがある。慣例的に，要因Factorとその要因に含まれる水準Levelを同時に表現し，\(\text{間}2 \times \text{間}3\)の分散分析(二要因の分散分析で，いずれもBetweenデザインであり，水準数がそれぞれ2と3)，といった言い方をすることがある。
9.2 分散分析のステップ
t検定において等分散性の仮定が成立するかどうかが事前に問題になったように，Betweenデザインにおいても分散の等質性は仮定されており，Leveneの検定などで事前に検証しておくべきである。またWithinデザインにおいては，データの組成に関わる分散共分散行列の非対角要素が全て等しいことが望ましいが，実践的にはそこまでの仮定が成立しているとは考えにくい。ただし分散分析としては，等分散性の仮定よりも，より緩やかな球面性の仮定が成立していればよいとされており，これを事前に検定することが一般的である。Welchの補正のように，球面性の仮定が成立していない場合は，自由度を補正することで検定の精度が維持される。
分散分析は多要因・多水準の平均値差の検定である。各水準ごとにt検定を繰り返せば良いのではないか，というアイデアは誰しも思いつくことであろうが，この方法は検定の目的である\(\alpha\)水準の制御ができなくなるという問題を含む。そこで多水準の場合は分散分析を行うことで，すべての要因・水準の母平均が同じであるという帰無仮説を検定し，効果の有無をまず明確にする。この帰無仮説が棄却されたらどこかに差があるわけだから，以後は慎重に\(\alpha\)水準を制御しつつ事後的な検定にすすむ。
水準間の差をみるための事後的な検定は，下位検定とも呼ばれる。その方法は多岐に渡り，ゴールドスタンダードは存在せず，往々にして分析者が利用しているソフトウェアが対応する手法が選択される。要因・水準が多くなると検証すべき組み合わせも多くなり，下位検定の手続きも非常に煩雑になる。統計ソフトウェアはそれこそ機械的に，幾重にも細かく分散分析表を分解して下位検定をつづけていってくれるが，いくら制御されているとはいえ検定を繰り返していることに変わりはないし，各下位検定の結果を一貫した総合的解釈をするのは困難である。実験計画はシンプルであるほうが望ましいし，複雑なモデルになるようであれば分散分析を超えた，階層線形モデルやベイジアンアプローチなどを取る方が良いだろう。
9.3 ANOVA君を使う
分散分析をRで実行するには，基本関数であるaovやcarパッケージなどを用いることができる。 もっとも，その出力は必ずしも親切ではないし，下位検定や効果量の算出などは別のパッケージ，別の関数を用いる必要がある。
筆者がお勧めするのは，大正大学の井関龍太が開発したanovakunである。パッケージ化されていないので，リンク先からソースコードを読み込んでanovakun関数を実行する必要があるが，さまざまな実験デザインに対応し，また下位検定や効果量，球面性の補正などおよそ分散分析で必要な手法は網羅されている。以下ではこれを用いた実践を行う。
anovakunの読み込みは，ソースコードをプロジェクトフォルダにダウンロードしてsource関数で読み込むか，インターネットに繋がっている状態でリンク先から直接ソースファイル(anovakun_489.txt)1をsource関数で読み込むといいだろう。
読み込みが終わるとEnvironタブにanovakun関数が含まれていることを確認しよう。
9.3.1 ANOVA君の入力とデータ
ANOVA君は伝統的にワイド型データから読み込むようになっている。すなわち，一行に1オブザベーション入っている形式である。Between計画の場合は，データの前に水準数を表すインデックスと最終的な従属変数の形に整形したデータが必要である。Within計画の場合は1行に1Obs.なのだから，反復した水準の数だけ右にデータを入れていく形に整形する。
しかしChapter 3.7 で述べたように，昨今は計算機にとって優しい型，ロング型での入力もおおく，ANOVA君もversion 4.4.0からロング型での入力も許すようになった。その場合はオプションlong=TRUE とロング型であることを明記する必要がある。
ANOVA君を使う時は，関数anovakunに，データ，要因計画の型，各要因の水準の順で入力する。ここで要因計画の型とは，文字列でBetween/Withinの違いを明示することになる。被験者のラベルを表す小文字のsを挟んで，左側に間(Between)要因，右側に内(Within)要因を入れる。例えば一要因Between計画の場合は"As"，二要因Within計画の場合は"sAB"，間1内2の混合計画であれば"AsBC"のようにする。
続いて入力する水準数は，要因の数だけ必要である。ただし，ロング型で入力した場合は自動的に水準数が計算されるので入力の必要がない。
このテキストでは，データの持ち替えについてすでに触れているので，色々扱いやすいロング型に整形して利用していくものとする。
9.4 Betweenデザイン
9.4.1 1way-ANOVA
もっとも単純な一要因3水準，Between計画の例から始めよう。仮想データの生成を行うことで，分散分析のメカニズムと共に見ていくことにする。
set.seed(123)
# 各群のサンプルサイズ
n1 <- 5
n2 <- 4
n3 <- 6
# 母平均，効果量，母SD
mu <- 10
delta <- 1
sigma <- 3
# 群平均
mu1 <- mu - (delta * sigma)
mu2 <- mu
mu3 <- mu + (delta * sigma)
# データセット
X1 <- rnorm(n1, mu1, sigma)
X2 <- rnorm(n2, mu2, sigma)
X3 <- rnorm(n3, mu3, sigma)
## 組み上げる
dat <- data.frame(
ID = 1:(n1 + n2 + n3),
group = as.factor(rep(LETTERS[1:3], c(n1, n2, n3))),
value = c(X1, X2, X3)
)
## データの確認
dat ID group value
1 1 A 5.318573
2 2 A 6.309468
3 3 A 11.676125
4 4 A 7.211525
5 5 A 7.387863
6 6 B 15.145195
7 7 B 11.382749
8 8 B 6.204816
9 9 B 7.939441
10 10 C 11.663014
11 11 C 16.672245
12 12 C 14.079441
13 13 C 14.202314
14 14 C 13.332048
15 15 C 11.332477### 実行
anovakun(dat, "As", long = TRUE, peta = TRUE)
[ As-Type Design ]
This output was generated by anovakun 4.8.9 under R version 4.4.1.
It was executed on Mon Aug 5 15:41:21 2024.
<< DESCRIPTIVE STATISTICS >>
------------------------------
group n Mean S.D.
------------------------------
A 5 7.5807 2.4331
B 4 10.1681 3.9548
C 6 13.5469 1.9483
------------------------------
<< ANOVA TABLE >>
== This data is UNBALANCED!! ==
== Type III SS is applied. ==
--------------------------------------------------------------
Source SS df MS F-ratio p-value p.eta^2
--------------------------------------------------------------
group 98.3840 2 49.1920 6.5897 0.0117 * 0.5234
Error 89.5804 12 7.4650
--------------------------------------------------------------
Total 187.9644 14 13.4260
+p < .10, *p < .05, **p < .01, ***p < .001
<< POST ANALYSES >>
< MULTIPLE COMPARISON for "group" >
== Shaffer's Modified Sequentially Rejective Bonferroni Procedure ==
== The factor < group > is analysed as independent means. ==
== Alpha level is 0.05. ==
------------------------------
group n Mean S.D.
------------------------------
A 5 7.5807 2.4331
B 4 10.1681 3.9548
C 6 13.5469 1.9483
------------------------------
-------------------------------------------------------
Pair Diff t-value df p adj.p
-------------------------------------------------------
A-C -5.9662 3.6062 12 0.0036 0.0108 A < C *
B-C -3.3789 1.9159 12 0.0795 0.0795 B = C
A-B -2.5873 1.4117 12 0.1834 0.1834 A = B
-------------------------------------------------------
output is over --------------------///出力結果は大きく分けて記述統計<< DESCRIPTIVE STATISTICS >>と，分散分析表<< ANOVA TABLE >>，下位検定<< POST ANALYSES >>に分けられる。記述統計はデータが正しく読み込めているかどうかのチェックに使おう。
一番のメインは分散分析表であり，平方和sum of squaresを自由度dfで割った，1自由度あたりのデータの散らばりを，群間と群内(誤差)との比で検証しているのが見て取れる。群間平方和が\(98.38\)，群内平方和が\(89.58\)であり，それぞれ自由度\(2\)(\(3\)水準\(-1\))と\(12\)(\(\sum_{j=1}^3 n_j-1\))から生じているので，平均平方Mean Squaresがそれぞれ\(49.19\)と\(7.47\)である。この比が\(6.5897\)で，自由度\(F(2,12)\)のF分布においてこの値以上の極端な数字が出る確率が5%を下回っている(実に\(p=0.0117\)である)ため，統計的に有意であると判断できる。 分散分析表のTotalのところで，全体のSSが群間SS+群内SSに一致していること，自由度も全体df=群間df+群内dfになっていることを確認しておこう。
また，anovakun関数の引数としてpeta = TRUEを指定したが，これは偏\(\eta^2\)(partial eta)と呼ばれる効果量を出力するためのオプションである。
今回は分散分析の時点で統計的な有意差が認められたため(\(F(2,12)=6.59, p < 0.05, \eta^2=0.52\))，続いて下位検定が表示されている。ANOVA君は下位検定についても複数のオプションを持っているが，デフォルトではShafferの修正Bonferroni検定が行われる。詳しくは専門書(永田 and 吉田 1997)を参照してほしいが，概略を説明すると，検証すべき仮説の数で有意水準を分割するというBonferroniの方法を，競合する仮説の数も考慮して分母を調整するというものである。
この計算の結果，A群とC群の間にのみ統計的な有意差が確認された(\(t(12)=3.61,p<0.05\))と言える。
9.4.2 2way-ANOVA
二要因の場合も見ておこう。ANOVA君の表記方法は要因計画の型が変わるだけで大きな変更はないが，交互作用interactionを考える必要があるところがポイントである。これも仮想データの組成を見ることでその意義がわかりやすくなるだろう。間2\(\times\)間2の実験デザインを例に，まずは各水準の理論的平均値がどのようにつくられるかをみておこう。
set.seed(123)
# 各群のサンプルサイズ
n <- 10
# 全体平均，効果量，母SD
mu <- 10
delta1 <- 1
delta2 <- 0 # ここではあえて要因Bの効果を0にしている
delta3 <- 2
sigma <- 3
# 効果の計算
effectA <- delta1 * sigma # Factor A
effectB <- delta2 * sigma # Factor B
effectAB <- delta3 * sigma # interaction
# 各群の平均
mu11 <- mu + effectA + effectB + effectAB
mu12 <- mu + effectA - effectB - effectAB
mu21 <- mu - effectA + effectB - effectAB
mu22 <- mu - effectA - effectB + effectAB効果の現れ方は相対的だから，要因Aが第一水準に+effectAの形で現れたら，第二水準には-effectA の形で現れる。要因Bについても同様である。交互作用については組み合わせにおいて生じるから，要因Aの第一水準と要因Bの第一水準の組み合わせのところに+effectABを充てる。ここでも効果は相対的に現れるという条件を守るために，要因Aの第一水準の中で+effectABの効果を相殺するために，要因Aの第一水準と要因Bの第二水準の組み合わせの符号が反転する。同様に，要因Bの第一水準の中で相殺するために要因Aの第二水準と要因Bの第一水準には-effectABが加わる。
このようにして考えられる理論的平均値に対して，外乱要因である誤差が生じて実現値が得られる。 組み上げて得られたデータを確認しておこう。
X11 <- rnorm(n, mean = mu11, sd = sigma)
X12 <- rnorm(n, mean = mu12, sd = sigma)
X21 <- rnorm(n, mean = mu21, sd = sigma)
X22 <- rnorm(n, mean = mu22, sd = sigma)
dat <- data.frame(
ID = 1:(n * 4),
FactorA = rep(1:2, each = n * 2),
FactorB = rep(rep(1:2, each = n), 2),
value = c(X11, X12, X21, X22)
)
dat ID FactorA FactorB value
1 1 1 1 17.3185731
2 2 1 1 18.3094675
3 3 1 1 23.6761249
4 4 1 1 19.2115252
5 5 1 1 19.3878632
6 6 1 1 24.1451950
7 7 1 1 20.3827486
8 8 1 1 15.2048163
9 9 1 1 16.9394414
10 10 1 1 17.6630141
11 11 1 2 10.6722454
12 12 1 2 8.0794415
13 13 1 2 8.2023144
14 14 1 2 7.3320481
15 15 1 2 5.3324766
16 16 1 2 12.3607394
17 17 1 2 8.4935514
18 18 1 2 1.1001485
19 19 1 2 9.1040677
20 20 1 2 5.5816258
21 21 2 1 -2.2034711
22 22 2 1 0.3460753
23 23 2 1 -2.0780133
24 24 2 1 -1.1866737
25 25 2 1 -0.8751178
26 26 2 1 -4.0600799
27 27 2 1 3.5133611
28 28 2 1 1.4601194
29 29 2 1 -2.4144108
30 30 2 1 4.7614448
31 31 2 2 14.2793927
32 32 2 2 12.1147856
33 33 2 2 15.6853770
34 34 2 2 15.6344005
35 35 2 2 15.4647432
36 36 2 2 15.0659208
37 37 2 2 14.6617530
38 38 2 2 12.8142649
39 39 2 2 12.0821120
40 40 2 2 11.8585870もちろん実際には，計画に応じたデータセットが得られているはずであり，各群のサンプルサイズが異なるなどの事情もあるだろう。しかしこうして，理論的にデータの組成を見ておくことで，サンプルサイズを変えたり効果量を変えたりしながら，どのように結果が変わってくるかを確認しながら進めることができる2。
それではこの仮想データを分析してみよう。
anovakun(dat, "ABs", long = TRUE, peta = TRUE)
[ ABs-Type Design ]
This output was generated by anovakun 4.8.9 under R version 4.4.1.
It was executed on Mon Aug 5 15:41:21 2024.
<< DESCRIPTIVE STATISTICS >>
-----------------------------------------
FactorA FactorB n Mean S.D.
-----------------------------------------
1 1 10 19.2239 2.8614
1 2 10 7.6259 3.1142
2 1 10 -0.2737 2.7924
2 2 10 13.9661 1.5819
-----------------------------------------
<< ANOVA TABLE >>
---------------------------------------------------------------------------
Source SS df MS F-ratio p-value p.eta^2
---------------------------------------------------------------------------
FactorA 432.7854 1 432.7854 61.4190 0.0000 *** 0.6305
FactorB 17.4478 1 17.4478 2.4761 0.1243 ns 0.0644
FactorA x FactorB 1668.9825 1 1668.9825 236.8545 0.0000 *** 0.8681
Error 253.6721 36 7.0464
---------------------------------------------------------------------------
Total 2372.8878 39 60.8433
+p < .10, *p < .05, **p < .01, ***p < .001
<< POST ANALYSES >>
< SIMPLE EFFECTS for "FactorA x FactorB" INTERACTION >
----------------------------------------------------------------------
Source SS df MS F-ratio p-value p.eta^2
----------------------------------------------------------------------
FactorA at 1 1900.7730 1 1900.7730 269.7492 0.0000 *** 0.8823
FactorA at 2 200.9950 1 200.9950 28.5243 0.0000 *** 0.4421
FactorB at 1 672.5693 1 672.5693 95.4480 0.0000 *** 0.7261
FactorB at 2 1013.8610 1 1013.8610 143.8826 0.0000 *** 0.7999
Error 253.6721 36 7.0464
----------------------------------------------------------------------
+p < .10, *p < .05, **p < .01, ***p < .001
output is over --------------------///基本的な結果の見方については，一要因のときと同じである。今回は要因Aと交互作用の効果を作り，正しく検出されている。下位検定については，要因Aが2水準であったためこちらの主効果の検証は必要なく(記述統計を見て群平均比較をすればよい)，交互作用についての単純効果の検証が行われている。
9.5 Withinデザイン
Withinデザインは対応のあるt検定の時と同じように，多次元正規分布からの生成として考えよう。すなわち各個体から得られるデータが相関しているという仮定をおくのである。以下のサンプルコードを読んで，データ生成過程を確認しよう。なお共分散は\(\rho_{xy}=\frac{s_{xy}}{s_xs_y}\)より\(s_{xy}=\rho_{xy}s_xs_y\)として整形している。
library(tidyverse)── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
✔ dplyr 1.1.4 ✔ readr 2.1.5
✔ forcats 1.0.0 ✔ stringr 1.5.1
✔ ggplot2 3.5.1 ✔ tibble 3.2.1
✔ lubridate 1.9.3 ✔ tidyr 1.3.1
✔ purrr 1.0.2
── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
✖ dplyr::filter() masks stats::filter()
✖ dplyr::lag() masks stats::lag()
ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errorslibrary(MASS)
次のパッケージを付け加えます: 'MASS'
以下のオブジェクトは 'package:dplyr' からマスクされています:
selectset.seed(42)
# 各群のサンプルサイズ
n <- 10
# 全体平均，効果量，母SD
mu <- 10
delta <- 1
s1 <- s2 <- s3 <- 1
rho12 <- 0.1
rho13 <- 0.3
rho23 <- 0.8
mus <- c(mu, mu + s1 * delta, mu - s1 * delta)
# 共分散行列の生成
Sigma <- matrix(NA, ncol = 3, nrow = 3)
Sigma[1, 1] <- s1^2
Sigma[2, 2] <- s2^2
Sigma[3, 3] <- s3^2
Sigma[1, 2] <- Sigma[2, 1] <- rho12 * s1 * s2
Sigma[1, 3] <- Sigma[3, 1] <- rho13 * s1 * s3
Sigma[2, 3] <- Sigma[3, 2] <- rho23 * s2 * s3
# データの生成
X <- mvrnorm(n, mus, Sigma) %>% as.data.frame()
# データの確認
X V1 V2 V3
1 10.625304 9.418518 7.493325
2 12.437964 11.249993 8.806719
3 8.604481 11.182418 8.722798
4 9.390742 10.181310 8.786312
5 9.567609 10.147592 9.194809
6 10.651739 11.005419 8.917299
7 9.125913 9.805634 7.511082
8 7.770294 12.462671 8.790231
9 6.909722 9.863405 7.429485
10 11.267590 10.798088 8.754522# Long型に整形
X <- X %>%
rowid_to_column("ID") %>%
pivot_longer(-ID) %>%
print()# A tibble: 30 × 3
ID name value
<int> <chr> <dbl>
1 1 V1 10.6
2 1 V2 9.42
3 1 V3 7.49
4 2 V1 12.4
5 2 V2 11.2
6 2 V3 8.81
7 3 V1 8.60
8 3 V2 11.2
9 3 V3 8.72
10 4 V1 9.39
# ℹ 20 more rows# 分析の実行
anovakun(X, "sA", long = TRUE, peta = TRUE, GG = TRUE)
[ sA-Type Design ]
This output was generated by anovakun 4.8.9 under R version 4.4.1.
It was executed on Mon Aug 5 15:41:21 2024.
<< DESCRIPTIVE STATISTICS >>
-----------------------------
name n Mean S.D.
-----------------------------
V1 10 9.6351 1.6609
V2 10 10.6115 0.9057
V3 10 8.4407 0.6777
-----------------------------
<< SPHERICITY INDICES >>
== Mendoza's Multisample Sphericity Test and Epsilons ==
-------------------------------------------------------------------------
Effect Lambda approx.Chi df p LB GG HF CM
-------------------------------------------------------------------------
name 0.0068 8.8720 2 0.0118 * 0.5000 0.5988 0.6392 0.5547
-------------------------------------------------------------------------
LB = lower.bound, GG = Greenhouse-Geisser
HF = Huynh-Feldt-Lecoutre, CM = Chi-Muller
<< ANOVA TABLE >>
--------------------------------------------------------------
Source SS df MS F-ratio p-value p.eta^2
--------------------------------------------------------------
s 16.4609 9 1.8290
--------------------------------------------------------------
name 23.6422 2 11.8211 10.7022 0.0009 *** 0.5432
s x name 19.8819 18 1.1045
--------------------------------------------------------------
Total 59.9849 29 2.0684
+p < .10, *p < .05, **p < .01, ***p < .001
<< POST ANALYSES >>
< MULTIPLE COMPARISON for "name" >
== Shaffer's Modified Sequentially Rejective Bonferroni Procedure ==
== The factor < name > is analysed as dependent means. ==
== Alpha level is 0.05. ==
-----------------------------
name n Mean S.D.
-----------------------------
V1 10 9.6351 1.6609
V2 10 10.6115 0.9057
V3 10 8.4407 0.6777
-----------------------------
----------------------------------------------------------
Pair Diff t-value df p adj.p
----------------------------------------------------------
V2-V3 2.1708 9.5342 9 0.0000 0.0000 V2 > V3 *
V1-V3 1.1945 2.3896 9 0.0406 0.0406 V1 > V3 *
V1-V2 -0.9764 1.6250 9 0.1386 0.1386 V1 = V2
----------------------------------------------------------
output is over --------------------///上記コードについていくつか解説をしておこう。 今回，各群の分散は同じにしつつ，変数間相関に大きな違いを持たせた。あえて球面性の仮定が成立しないような例を得たかったからで，出力の<< SPHERICITY INDICES >>をみると統計量\(\lambda\)のあとの\(p\)値が5%を下回っており，「球面性が成立している」という帰無仮説が棄却されていることがわかる。この時いくつかの補正法があるが，今回はGreenhouse-Geisserの補正を当てることにしている。それがanovakun関数のなかのGG=TRUEの箇所である。
これを踏まえて分散分析表が示されている。ここでも全体平方和SS，自由度dfが各行の要素の和になっていることが確認できるが，その因子名のところにsが含まれていることが確認できる。これが個体ごとの変動を表しており，誤差から個人差を取り除いて効果の検証ができていることがわかる。
分散分析は加法的，線形的な分解であるからわかりやすく，要因が複雑に組み合わさることがあっても基本的に今回のパーツを組み上げることで理解できる。言い換えると，データが先にある実践的な場合には平方和をひとつひとつ丁寧に紐解いていくことで理解できる。実にanovakunの前進であるanova4では4要因，anovakunでは26要因までのデザインを分析することが可能である。もっとも4要因計画にもなると3次の交互作用まで考えられ，主効果と合わせてこれらの交互作用効果を解釈するのは困難である。anoakunは2次以上の交互作用が見られた場合，自動的に下位検定を行ってくれないので，要因の水準ごとにデータを分割して，分散分析表を解体しつつ分析する必要がある。3
しかしもちろん，これには検定の多重性の問題が関わってくるから，あまり推奨される手法ではない。ごく少ない要因で，主効果の有無を検証することを主眼においた丁寧な実験デザインを組み立てることを試みるべきである。
またここでは，仮想データを作ることで，得られたデータの背後にある生成メカニズムに注目した。「与えられたデータを分解する」のが分散分析であるのに対し，リバースエンジニアリングからアプローチしたのである。こうすることで，分散分析の見えない仮定に注意が向くことを期待している。簡便のために，いくつかのパラメータを均質化するなどしたが，実践的には群ごとのサンプルサイズが異なることも少なくないだろうし，群間の分散や共分散が均質であることを前提とするのは難しいだろう。これを考慮した細かい作り込みも，リバースエンジニアリングによって生成メカニズムがわかっている場合には応用が可能である。さらに，どこの水準間にどのような効果があると仮定されるか，といった精緻な仮説があるのなら，そこだけをターゲットにした分析を行うことも可能である。
分散分析は，あくまでも大雑把な全体的傾向を見るためのものであることに留意しよう。心理学のデータがより精緻な仮定に耐えうる精度を持つものになれば，分散分析は過去の遺物となるかもしれない。
9.6 課題
- 以下のデータセットは一要因4水準Between計画で得られたものです。分散分析を行って，要因の効果があるかどうか，水準間に差があるとすればどこに見られるかを報告してください。なおこのデータセットはこちらex_anova1.csvからダウンロード可能です。
ID group value
1 1 A 14.37
2 2 A 15.11
3 3 A 16.11
4 4 A 11.17
5 5 A 14.51
6 6 A 7.85
7 7 A 10.65
8 8 B 16.45
9 9 B 11.76
10 10 B 19.11
11 11 B 19.62
12 12 C 2.92
13 13 C 6.27
14 14 C 1.82
15 15 C -0.10
16 16 C 5.30
17 17 C 1.57
18 18 D 8.33
19 19 D 2.71
20 20 D 5.97
21 21 D 4.97
22 22 D 1.65
23 23 D 8.73
24 24 D 5.93
25 25 D 4.27- 以下のデータセットは一要因4水準Within計画で得られたものです。分散分析を行って，要因の効果があるかどうか，水準間に差があるとすればどこに見られるかを報告してください。なおこのデータセットはこちらex_anova2.csvからダウンロード可能です。
V1 V2 V3 V4
1 11.32 12.99 9.34 -0.14
2 10.77 13.84 14.74 3.52
3 9.86 12.26 12.56 2.60
4 8.74 11.59 14.27 0.68
5 11.12 12.93 12.92 1.13
6 9.65 16.55 12.60 2.32
7 9.72 14.64 9.69 -1.34
8 12.02 11.18 14.43 2.64
9 10.00 10.79 9.19 -1.09
10 10.04 15.53 13.38 1.82
11 10.20 11.56 11.02 -0.05
12 7.81 9.29 12.20 -3.25- 間(3)\(\times\)間(3)の分散分析モデルの仮想データセットを作りましょう。そのデータに分散分析を適用し，仮定した要因の効果がみられるか(あるいは効果がないと仮定した場合に正しく検出されないか)を確認しましょう。
- 【発展課題】二要因混合計画分散分析(間\(\times\)内)の仮想データセットを作りましょう。そのデータに分散分析を適用し，仮定した要因の効果がみられるか(あるいは効果がないと仮定した場合に正しく検出されないか)を確認しましょう。
2024.06.12時点での最新バージョンが4.8.9である。リンク先URLは，公式サイトからソースファイルのリンクをコピーして貼り付けると良い。↩︎
かつては分散分析は手計算でできる分析モデルであり，得られたデータを平方和に分解していくプロセスをたどりながら分散分析のメカニズムが体得されるという教育が多く見られた。ただしその方法は計算に時間がかかること，ミスが混在しやすいことに加え，手元のデータが唯一無二のものであるという印象を強くすることが懸念される。推測統計学においては，手元のデータはあくまでも実現値に過ぎないと考えるのであり，乱数を生成して幾つでも自在に作り出せる経験を得た方が教育効果として良いのではないか，と筆者は考えている。↩︎
anovakunの補助関数anovatanを用いることで，注目したい要因ごとにデータを分割してくれる。詳しくは公式サイトのマニュアルを参照してほしい。↩︎