平方和と自由度の計算・分散分析表の作成

実験計画法　第2章 §2.2　★★★ $S_T=S_A+S_e$ の分解・自由度・期待平均平方・$F$ 検定。1級で毎回問われる中心の手続き

平方和の分解

データの総平均からのずれを、水準平均を経由して2つに分ける。

$$x_{ij}-\bar{\bar x}=(\bar x_{i\cdot}-\bar{\bar x})+(x_{ij}-\bar x_{i\cdot}).\tag{2.5}$$

右辺の第1項は「どの水準でとったか」によるずれ（主効果の推定値 $\hat a_i$）、第2項は「同じ水準の中でのばらつき」（残差 $e_{ij}$）である。両辺を2乗して $i,j$ について和をとる。

証明

(2.5) を2乗して和をとると

$$S_T=S_A+S_e+2\sum_{i=1}^a\sum_{j=1}^r(\bar x_{i\cdot}-\bar{\bar x})(x_{ij}-\bar x_{i\cdot}).$$

交差項では、因子 $\bar x_{i\cdot}-\bar{\bar x}$ が $j$ によらないので、内側の和の外に出せる。

$$2\sum_{i=1}^a(\bar x_{i\cdot}-\bar{\bar x})\sum_{j=1}^r(x_{ij}-\bar x_{i\cdot})$$

内側の和は残差の水準内の和で、§2.1 で見たとおり0である。よって交差項は消える。

分解の意味は幾何で見ると分かりやすい。データ $N$ 個を並べたベクトルを $\bm{x}$ とする。

この分解は、総平均を引いた $\bm{x}-\bar{\bar x}\bm{1}$ を、「水準ごとに一定の値をとるベクトル」の成分と「水準ごとに和が0のベクトル」の成分に分けている。前者と後者は内積が0（直交）なので、Pythagoras の定理によって長さの2乗が足し算になる。$S_T=S_A+S_e$ はこの直交分解の式である。

計算式と修正項

定義式 (2.6) のまま計算すると、平均が割り切れないときに端数の2乗を何度も扱うことになる。手計算では次の形を使う。

証明

$S_T=\sum\sum(x_{ij}^2-2x_{ij}\bar{\bar x}+\bar{\bar x}^2)=\sum\sum x_{ij}^2-2\bar{\bar x}T+N\bar{\bar x}^2$。$\bar{\bar x}=T/N$ を代入すると $-2T^2/N+T^2/N=-T^2/N=-CT$ となる。$S_A$ も同様に展開する。

$$\begin{aligned}S_A&=r\sum_i(\bar x_{i\cdot}^2-2\bar x_{i\cdot}\bar{\bar x}+\bar{\bar x}^2) =r\sum_i\bar x_{i\cdot}^2-2\bar{\bar x}\cdot r\sum_i\bar x_{i\cdot}+ra\bar{\bar x}^2\\ &=\sum_i\frac{T_{i\cdot}^2}{r}-CT.\end{aligned}$$

最後の行では $r\bar x_{i\cdot}^2=T_{i\cdot}^2/r$ と $r\sum_i\bar x_{i\cdot}=T$ を使った。総数は $ra=N$ である。$S_e$ は定理 2.4 から従う。

「修正項」の名前は、原点のまわりの2乗和 $\sum\sum x_{ij}^2$ を、総平均のまわりの2乗和へ直すために引く量であることに由来する。覚え方の要点は「合計を2乗して、その合計に含まれるデータの個数で割る」ことで、$CT=T^2/N$ も $T_{i\cdot}^2/r$ もこの形をしている。後の章の2元配置や直交配列表でも、平方和はすべて「（合計）$^2$／（個数）の和 $-CT$」の形で書ける。

自由度

平方和は足し合わせる項の数が多いほど大きくなるので、平方和の大小をそのまま比べても意味が無い。平方和を「独立に動ける成分の数」で割って比べる。この数が自由度である。

各値の根拠を数える。

• 総平方和：偏差 $x_{ij}-\bar{\bar x}$ は $N$ 個ある。偏差の総和は0なので、$N-1$ 個を決めれば残り1個は決まる。よって $\phi_T=N-1$。
• 級間平方和：偏差 $\bar x_{i\cdot}-\bar{\bar x}$ は $a$ 個ある。その和は $\sum_i\bar x_{i\cdot}-a\bar{\bar x}=0$ なので、$\phi_A=a-1$。
• 誤差平方和：残差 $e_{ij}$ は $N$ 個ある。制約は水準ごとに1本ずつ（$\sum_j e_{ij}=0$）あるので、$\phi_e=N-a$。

自由度は幾何の言葉では部分空間の次元である。平方和の分解のところで述べた「水準ごとに一定のベクトルで、成分の総和が0のもの」の全体は$a-1$ 次元、「水準ごとに和が0のベクトル」の全体は $a(r-1)$ 次元である。$S_A$ と $S_e$ は、データをそれぞれの部分空間へ直交射影したベクトルの長さの2乗であり、自由度はその部分空間の次元に等しい。この見方は第4章で検定の理論を述べるとき、また数理統計の線形モデルの章とつなげるときに使う。

なぜ $n-1$ で割るのか

自由度で割る操作は、不偏分散

$$S^2=\frac{1}{n-1}\sum_{i=1}^n(X_i-\bar X)^2$$

が $n$ でなく $n-1$ で割る理由と同じである。

証明

$\sum_k(X_k-\bar X)^2=\sum_k X_k^2-n\bar X^2$ と書ける（命題 2.5 の $S_T$ と同じ展開）。$E[X_k^2]=V[X_k]+(E[X_k])^2=\sigma^2+\mu^2$、$E[\bar X^2]=V[\bar X]+\mu^2=\sigma^2/n+\mu^2$ を代入すると

$$E\Bigl[\sum_k(X_k-\bar X)^2\Bigr]=n(\sigma^2+\mu^2)-n\Bigl(\frac{\sigma^2}{n}+\mu^2\Bigr)=(n-1)\sigma^2. \quad\square$$

偏差 $X_k-\bar X$ の基準は、真の平均でなく、データから計算した $\bar X$ である。2乗和 $\sum(X_k-c)^2$ を最小にする定数は $c=\bar X$ である。そのため次の不等式が成り立ち、左辺の期待値は右辺の期待値 $n\sigma^2$ より小さくなる。

$$\sum_k(X_k-\bar X)^2\le\sum_k(X_k-\mu)^2$$

減る量の期待値はちょうど $\sigma^2$ の1個分（$nV[\bar X]$）である。平均を1つ推定するごとに自由度が1つ減るという言い方は、このことを指す。2乗和を自由度で割れば、平均の推定で失われた分を補正した $\sigma^2$ の不偏推定量になる。

期待平均平方

平方和を自由度で割ったものを平均平方と呼ぶ。

$$V_A=\frac{S_A}{\phi_A},\qquad V_e=\frac{S_e}{\phi_e}.\tag{2.9}$$

平均平方の期待値（期待平均平方 $E(V)$）を求めると、$F$ 検定が何をしているかが分かる。

証明

$\mu+a_i$ は引き算で消えるので、残差は誤差だけで書ける。

$$x_{ij}-\bar x_{i\cdot}=\varepsilon_{ij}-\bar\varepsilon_{i\cdot}$$

よって水準ごとの部分は、独立な誤差 $r$ 個の偏差平方和である。各誤差は期待値0・分散 $\sigma^2$ をもつ。命題 2.7 でその期待値は $(r-1)\sigma^2$。$a$ 水準分を足して $E(S_e)=a(r-1)\sigma^2$。

$S_A$ については、(2.3) と $\sum a_i=0$ から次のように書ける。

$$\bar x_{i\cdot}-\bar{\bar x}=a_i+(\bar\varepsilon_{i\cdot}-\bar{\bar\varepsilon})$$

2乗して和をとると

$$S_A=r\sum_i a_i^2+2r\sum_i a_i(\bar\varepsilon_{i\cdot}-\bar{\bar\varepsilon})+r\sum_i(\bar\varepsilon_{i\cdot}-\bar{\bar\varepsilon})^2.$$

第2項の期待値は0である。第3項について考える。$\bar\varepsilon_{1\cdot},\dots,\bar\varepsilon_{a\cdot}$ は独立で、期待値0・分散 $\sigma^2/r$ をもつ。その平均が $\bar{\bar\varepsilon}$ なので、命題 2.7 を $n=a$ として使えば次を得る。

$$E\Bigl[\sum_i(\bar\varepsilon_{i\cdot}-\bar{\bar\varepsilon})^2\Bigr]=\frac{(a-1)\sigma^2}{r}$$

$r$ 倍して $(a-1)\sigma^2$ を得る。自由度で割れば $E(V_e)$、$E(V_A)$ の式になる。

$\sigma_A^2$ は主効果の大きさを分散の形で表した量で、母数模型（水準を固定して選んだ模型）では確率変数でなく定数である。

定理 2.8 から、2つの平均平方は次の量の不偏推定量である。

• $V_e$ は誤差分散 $\sigma^2$ を推定する。
• $V_A$ は、誤差の分に主効果の分を足した $\sigma^2+r\sigma_A^2$ を推定する。

帰無仮説は次のとおりである。

$$H_0:\ a_1=\dots=a_a=0\quad(\text{すなわち}\ \sigma_A^2=0)$$

$H_0$ のもとでは2つの平均平方がどちらも誤差分散を推定するので、比 $V_A/V_e$ は1の近くに来る。主効果があれば $V_A$ の期待値は $r\sigma_A^2$ だけ大きくなり、比は1より大きくなる。これが分散比で検定する理由である。また、この式から $\sigma_A^2$ の不偏推定量 $\hat\sigma_A^2=(V_A-V_e)/r$ も得られる。

分散比の分布と $F$ 検定

検定に使うには、比の分布が要る。正規母集団の標本分布の結果（数理統計の第5章 §5.2）として、$X_1,\dots,X_n$ が独立に $\mathrm{N}(\mu,\sigma^2)$ に従うとき

$$\frac{1}{\sigma^2}\sum_{k=1}^n(X_k-\bar X)^2\sim\chi^2_{n-1},\qquad \bar X\ \text{と}\ \sum_{k=1}^n(X_k-\bar X)^2\ \text{は独立}\tag{2.12}$$

が成り立つ。これを借りると次が言える。

証明

水準 $i$ のデータは、独立に $\mathrm{N}(\mu+a_i,\sigma^2)$ に従う。(2.12) より、次の量は $\chi^2_{r-1}$ に従い、$\bar x_{i\cdot}$ と独立である。

$$\frac{1}{\sigma^2}\sum_{j=1}^r(x_{ij}-\bar x_{i\cdot})^2$$

異なる水準のデータは独立である。そのため $S_e/\sigma^2$ は、独立な $\chi^2_{r-1}$ を水準の数だけ足したものである。カイ二乗分布の再生性により、これは $\chi^2_{a(r-1)}$ に従う。また $S_e$ は、水準平均の組全体と独立である。

帰無仮説のもとでは、水準平均 $\bar x_{1\cdot},\dots,\bar x_{a\cdot}$ は独立に同じ分布 $\mathrm{N}(\mu,\sigma^2/r)$ に従う。再び (2.12) を $n=a$ として使うと次を得る。

$$\frac{\sum_i(\bar x_{i\cdot}-\bar{\bar x})^2}{\sigma^2/r}=\frac{S_A}{\sigma^2}\sim\chi^2_{a-1}$$

$S_A$ は水準平均だけの関数なので $S_e$ と独立である。独立なカイ二乗変数をそれぞれの自由度で割った比は $F$ 分布に従うので、$F_0=\dfrac{(S_A/\sigma^2)/\phi_A}{(S_e/\sigma^2)/\phi_e}\sim F_{\phi_A,\phi_e}$。$\sigma^2$ は約分で消える。

$F_0$ が大きいほど主効果がある証拠になるので、棄却域は上側にとる。有意水準 $\alpha$ の検定は

$$F_0\ge F_{\phi_A,\phi_e}(\alpha)\ \text{ならば}\ H_0\ \text{を棄却する}\tag{2.14}$$

である。実験計画法の慣習では、有意水準5%で有意なら $F_0$ の右に「$*$」を付ける。1%で有意なら「$**$」を付ける。検出力や非心 $F$ 分布を使った検定の意味づけは第4章 §4.1 で扱う。

分散分析表

ここまでの量を1つの表にまとめたものが分散分析表である（表 2.2）。

表の作り方の手順をまとめておく。

• (i) 水準合計 $T_{i\cdot}$、総合計 $T$、データの2乗和を計算する。
• (ii) $CT=T^2/N$ を求める。
• (iii) (2.7) で $S_T$ と $S_A$ を求める。残りを $S_e=S_T-S_A$ とする。
• (iv) 自由度を (2.8) で数え、$\phi_T=\phi_A+\phi_e$ を確かめる。
• (v) 平均平方と $F_0$ を計算し、$F$ 表の値と比べて判定する。

図 2.2 は、通し例で使う $F_{3,12}$ 分布の密度と棄却域である。

数値変換

測定値の桁が大きいときや小数点以下が長いときは、平方和の計算の前にデータを一次変換して扱いやすい数にする。定数 $x_0$ と $h>0$ を選び、次のようにおく。

$$X_{ij}=(x_{ij}-x_0)\,h$$

証明

平方和はすべて「データから何らかの平均を引いた偏差の2乗和」である。$X_{ij}=(x_{ij}-x_0)h$ のとき、どの平均も同じ変換を受ける（たとえば $\bar X_{i\cdot}=(\bar x_{i\cdot}-x_0)h$）。

そのため偏差では $x_0$ が引き算で消え、偏差は $h$ 倍になる。2乗和は $h^2$ 倍になる。$F_0$ は2つの平均平方の比なので、$h^2$ が約分で消える。

平行移動 $x_0$ は平方和に影響しないので、総平均に近い値を選んで数を小さくする。倍率 $h$ は小数点を消すために10の累乗を選ぶことが多い。平方和を元の単位に戻すときは $h$ でなく $h^2$ で割る。一方、信頼区間の幅は $\sqrt{V_e}$ に比例するので、$h$ で割って戻す。

寄与率

$F$ 検定は「効果があるか無いか」を判定するが、効果が全体のばらつきのどれだけを占めるかは表さない。それを表すのが寄与率である。

単純に $S_A/S_T$ をとると、効果が無くても0にならない。定理 2.8 のとおり、$S_A$ の期待値には誤差の分 $\phi_A\sigma^2$ が含まれるためである。そこで誤差の分を差し引く。

$E(S_A-\phi_A V_e)=\phi_A\sigma^2+r\sum a_i^2-\phi_A\sigma^2=r\sum a_i^2$ なので、純変動は主効果だけによる変動の不偏推定量である。差し引いた $\phi_A V_e$ は誤差の側へ戻すので、$S_e'$ は誤差平方和に同じ量を足したものになる。

確認テスト

オリジナル

1. 2.2-1 1元配置で誤差平方和の自由度が $a(r-1)$ になる理由を、残差が満たす制約を使って1〜2文で説明せよ。
2. 2.2-2 確認テスト 2.1-2 のデータについて、$CT$ と平方和 $S_T,S_A,S_e$ を計算式で求め、分散分析表を作って有意水準5%と1%で判定せよ。$F_{2,6}(0.05)=5.14$、$F_{2,6}(0.01)=10.92$ を用いてよい。
3. 2.2-3 確認テスト 2.1-2 のデータが、元の測定値 $x$ を $X=(x-2.30)\times100$ と変換したものであったとする。元の単位での平方和 $S_A,S_e$ と分散比 $F_0$ を求めよ。
4. 2.2-4 定義式 $S_A=r\sum_i(\bar x_{i\cdot}-\bar{\bar x})^2$ から、計算式 $S_A=\sum_i T_{i\cdot}^2/r-T^2/N$ を導け。
5. 2.2-5 確認テスト 2.2-2 の結果から、純変動 $S_A'$ と寄与率 $\rho_A,\rho_e$ を求めよ。また、$S_A/S_T$ を寄与率として使ってはいけない理由を期待値を使って述べよ。