乱塊法とは

実験計画法　第11章 §11.1　★★★ ブロック因子を変量として構造式に入れ、期待平均平方・検定・推定・分散成分の推定を導く

実験を何日かに分けて行うと、日によって室温や原料の状態が変わり、同じ条件でも測定値の水準がずれる。このずれを放っておくと、すべて誤差に混ざって誤差分散が大きくなり、本当は差のある要因が「有意でない」と判定されることになる。乱塊法（randomized block design）は、このずれを生む単位（日、原料のロット、装置、圃場の区画など）をブロックとしてまとめ、各ブロックの中で要因のすべての水準を1回ずつ無作為な順に実施する計画である。

解析では、ブロックを1つの因子として分散分析表に1行立て、ブロック間の変動を誤差から分けて取り出す。

この章で作る道具は4つある。

1. ブロックを変量因子とみなしたデータの構造式と、その下での期待平均平方。
2. 要因の検定が、繰返しのない2元配置（第6章）と同じ計算で済む理由。
3. 分散分析後の推定の非対称性。水準の母平均そのものの信頼区間にはブロックの分散成分が入るのに、水準間の差の信頼区間には入らない。
4. ブロック因子の分散成分 $\sigma_R^2$ の点推定と区間推定。Satterthwaite の等価自由度による近似区間と、分散比 $\sigma_R^2/\sigma^2$ の $F$ 分布による正確な区間の2通りを作る。

統計応用（理工学）の実験計画法の大問では、「データから分散分析表を作る → 期待平均平方を書く → 検定する →最適水準の母平均と差を区間推定する → ブロックの分散成分を推定する」という連鎖で、この章の道具がそのまま使われる。

局所管理という考え方

Fisher は実験の計画に3つの原則を置いた。

1. 反復：同じ条件を複数回行い、誤差の大きさを推定できるようにする。
2. 無作為化：実験の順序や割り付けを乱数で決め、系統的なずれが特定の水準に偏ってかからないようにする。
3. 局所管理：実験の場を均一な小さいまとまりに分け、まとまりの中で比較する。

乱塊法は局所管理をもっとも素直に実現した計画で、「乱塊」の「塊」がブロックを、「乱」がブロックの中での無作為化を指す。

例として、触媒 $A$ の4種類 $A_1,\dots,A_4$ を比べたいが、1日に4回しか反応を行えないとする。5日かけて合計20回の実験をする。

20回の順番を完全に無作為に決める（完全無作為化法の1元配置）と、日と水準の組合せが偏る。たとえばある日に $A_1$ が2回、別の日に $A_4$ が0回、ということが起きる。日による水準のずれが大きいと、その偏りがそのまま水準間の差に混ざる。

乱塊法では、各日（ブロック）に4種類を1回ずつ割り当て、日の中での順番だけを無作為に決める（図 11.1）。すると、どの水準も各日に1回ずつ現れるので、日のずれはすべての水準に同じだけ加わり、水準どうしを比べる段階で打ち消し合う。

ブロックの中で順番を無作為にする理由は、ブロックの中にも時間とともに変わる要因（装置の温まり方、作業者の疲れ）があるからである。毎日 $A_1\to A_2\to A_3\to A_4$ の順に実施すると、「朝に行った」という条件が $A_1$ にだけ重なる。すると朝と夕方の差が、$A_1$ と $A_4$ の差に混ざる。ブロック内の無作為化は、ブロック化で取り除けなかった残りのずれを、どの水準にも偏らない誤差として扱えるようにする手続きである。

データの構造式

因子 $A$ の水準数を $a$ とおく。ブロックの数は $r$ とする。ブロック $R_j$ の中で水準 $A_i$ を実施した測定値を、次のように書く。

$$x_{ij}\qquad(i=1,\dots,a,\quad j=1,\dots,r).$$

総データ数は $N=ar$ である。

ブロックの効果を $r_j$ と書く本も多い。ここではブロック数 $r$ との混同を避けて、ギリシャ文字のガンマで書く。

構造式には $A\times R$ の交互作用の項が無い。これは「ブロックによって水準間の差が変わらない」という仮定である。

繰返しのない2元配置（第6章）と同じく、各セルのデータが1個しか無いので、交互作用を入れると誤差と区別できなくなる。そのため仮定として置くしかない。実際に交互作用があれば、それは誤差に含まれて推定される。

構造式 (11.1) からデータの共分散の構造が決まる。同じブロックの異なる水準のデータ $x_{ij}$ と $x_{i'j}$ を比べる（ただし水準は異なる）。2つは共通のブロック効果を含む。よって

$$V[x_{ij}]=\sigma_R^2+\sigma^2,\qquad \mathrm{Cov}(x_{ij},x_{i'j})=V[\gamma_j]=\sigma_R^2,\qquad \mathrm{Cov}(x_{ij},x_{i'j'})=0\quad(j\neq j')\tag{11.2}$$

となる。同じブロックの中のデータは正の相関を持つ。この相関

$$\frac{\sigma_R^2}{\sigma_R^2+\sigma^2}$$

をブロック内相関（級内相関）と呼ぶ。乱塊法の解析では、この正の相関を、同じブロックの中で差を取ることによって打ち消す。

平方和の分解と分散分析表

記号を次のように置く。

$$\begin{aligned}&\text{水準 }{A_i}\text{ の合計と平均}: && T_{i\cdot}=\sum_j x_{ij},\quad \bar x_{i\cdot}=T_{i\cdot}/r,\\ &\text{ブロック }{R_j}\text{ の合計と平均}: && T_{\cdot j}=\sum_i x_{ij},\quad \bar x_{\cdot j}=T_{\cdot j}/a,\\ &\text{総合計・総平均・修正項}: && T=\sum_{i,j}x_{ij},\quad \bar x=T/N,\quad CT=T^2/N .\end{aligned}$$

各データの総平均からの偏差を

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

と3つに分ける。右辺の3つの項は、$(i,j)$ について和を取ると互いの積の和がすべて 0 になる。確かめ方は次のとおりである。

• 第1項と第3項：積の和は次の形になる。

$$\sum_i(\bar x_{i\cdot}-\bar x)\sum_j(x_{ij}-\bar x_{i\cdot}-\bar x_{\cdot j}+\bar x).$$

内側の和は、$\sum_j\bar x_{\cdot j}=r\bar x$ を使うと

$$T_{i\cdot}-r\bar x_{i\cdot}-\sum_j\bar x_{\cdot j}+r\bar x=0-r\bar x+r\bar x=0 .$$

• 第2項と第3項：添字 $i$ と $j$ の役割を入れ替えた同じ計算で 0 になる。
• 第1項と第2項：$\sum_i(\bar x_{i\cdot}-\bar x)=0$ から 0 になる。

したがって両辺を2乗して和を取ると

$$S_T=S_A+S_R+S_e,\tag{11.4}$$

$$\begin{align*}S_T&=\sum_{i,j}(x_{ij}-\bar x)^2=\sum_{i,j}x_{ij}^2-CT, & S_A&=r\sum_i(\bar x_{i\cdot}-\bar x)^2=\sum_i\frac{T_{i\cdot}^2}{r}-CT,\\S_R&=a\sum_j(\bar x_{\cdot j}-\bar x)^2=\sum_j\frac{T_{\cdot j}^2}{a}-CT, & S_e&=\sum_{i,j}(x_{ij}-\bar x_{i\cdot}-\bar x_{\cdot j}+\bar x)^2=S_T-S_A-S_R\tag{11.5}\end{align*}$$

を得る。自由度は次のとおりで、$\phi_T=\phi_A+\phi_R+\phi_e$ が成り立つ。

$$\phi_T=N-1,\qquad \phi_A=a-1,\qquad \phi_R=r-1,\qquad \phi_e=(a-1)(r-1).$$

平方和と自由度の計算は、繰返しのない2元配置で $B$ をブロック因子 $R$ に読み替えたものとまったく同じである。乱塊法に固有の内容は、次に導く期待平均平方と分散分析後の推定である。

期待平均平方と分散分析表

証明

構造式を各平均に代入する。誤差とブロック効果の平均を次のように書く。

$$\bar\gamma=\frac1r\sum_j\gamma_j,\quad \bar\varepsilon_{i\cdot}=\frac1r\sum_j\varepsilon_{ij},\quad \bar\varepsilon_{\cdot j}=\frac1a\sum_i\varepsilon_{ij},\quad \bar\varepsilon=\frac1N\sum_{i,j}\varepsilon_{ij}.$$

$\sum_i a_i=0$ より

$$\bar x_{i\cdot}=\mu+a_i+\bar\gamma+\bar\varepsilon_{i\cdot},\qquad \bar x_{\cdot j}=\mu+\gamma_j+\bar\varepsilon_{\cdot j},\qquad \bar x=\mu+\bar\gamma+\bar\varepsilon .$$

したがって

$$\begin{aligned}\bar x_{i\cdot}-\bar x&=a_i+(\bar\varepsilon_{i\cdot}-\bar\varepsilon),\\ \bar x_{\cdot j}-\bar x&=(\gamma_j-\bar\gamma)+(\bar\varepsilon_{\cdot j}-\bar\varepsilon),\\ x_{ij}-\bar x_{i\cdot}-\bar x_{\cdot j}+\bar x&=\varepsilon_{ij}-\bar\varepsilon_{i\cdot}-\bar\varepsilon_{\cdot j}+\bar\varepsilon .\end{aligned}$$

1行目にブロックの効果が残っていないことに注意する。水準の平均 $\bar x_{i\cdot}$ はどれも全ブロックを1回ずつ通っているので、$\bar\gamma$ を同じだけ含み、総平均を引くと消える。3行目でも、ブロック効果は $\bar x_{\cdot j}$ に含まれる同じ項と打ち消し合う。

$S_A$ の期待値。誤差の水準平均 $\bar\varepsilon_{1\cdot},\dots,\bar\varepsilon_{a\cdot}$ は独立に $\mathrm{N}(0,\sigma^2/r)$ に従う。それらの平均が $\bar\varepsilon$ である。標本分散の性質から

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

交差項は期待値 0 なので

$$E(S_A)=r\sum_i a_i^2+r\cdot\frac{(a-1)\sigma^2}{r}=(a-1)(\sigma^2+r\sigma_A^2).$$

$S_R$ の期待値。$u_j=\gamma_j+\bar\varepsilon_{\cdot j}$ とおく。

これらは独立に $\mathrm{N}(0,\sigma_R^2+\sigma^2/a)$ に従う。2つの項は独立で、異なるブロックの間に共通の項が無いからである。また $\bar x_{\cdot j}-\bar x=u_j-\bar u$ である。よって

$$S_R=a\sum_j(u_j-\bar u)^2,\qquad E(S_R)=a(r-1)\Bigl(\sigma_R^2+\frac{\sigma^2}{a}\Bigr)=(r-1)(\sigma^2+a\sigma_R^2).\tag{11.7}$$

$S_e$ の期待値。残差 $e_{ij}=\varepsilon_{ij}-\bar\varepsilon_{i\cdot}-\bar\varepsilon_{\cdot j}+\bar\varepsilon$ の分散を直接計算する。これを誤差の1次結合 $\sum_{k,l}c_{kl}\varepsilon_{kl}$ と書くと、係数は次のように因数分解できる。

$$c_{kl}=\delta_{ik}\delta_{jl}-\frac{\delta_{ik}}{r}-\frac{\delta_{jl}}{a}+\frac1N=\Bigl(\delta_{ik}-\frac1a\Bigr)\Bigl(\delta_{jl}-\frac1r\Bigr).$$

$\sum_k(\delta_{ik}-1/a)^2=1-1/a$ を使うと

$$\sum_{k,l}c_{kl}^2=\Bigl(1-\frac1r\Bigr)\Bigl(1-\frac1a\Bigr)$$

となる。よって $E(e_{ij}^2)=\sigma^2(1-1/a)(1-1/r)$ で、$N$ 個を足して

$$E(S_e)=\frac{ar\sigma^2(a-1)(r-1)}{ar}=(a-1)(r-1)\sigma^2 .$$

各平方和を自由度で割れば (11.6) を得る。

分布と独立性。誤差ベクトル $\bm{\varepsilon}=(\varepsilon_{ij})\in\mathbb{R}^N$ の空間を、互いに直交する4つの部分空間に分ける。

1. 総平均の方向（1次元）。
2. 「$i$ だけによる成分の総平均からの偏差」の空間（$a-1$ 次元）。
3. 「$j$ だけによる成分の偏差」の空間（$r-1$ 次元）。
4. 残り（$(a-1)(r-1)$ 次元）。

誤差ベクトルを2番目・3番目・4番目の部分空間へ直交射影すると、その成分はそれぞれ次のようになる。

$$\bar\varepsilon_{i\cdot}-\bar\varepsilon,\qquad \bar\varepsilon_{\cdot j}-\bar\varepsilon,\qquad e_{ij}.$$

Cochran の定理によれば、$\bm{\varepsilon}\sim\mathrm{N}(\bm 0,\sigma^2 I_N)$ の互いに直交する部分空間への射影は独立である。しかも射影の長さの2乗を $\sigma^2$ で割ると、次元を自由度とするカイ二乗分布に従う。ここから次が出る。

• $S_e$ は4番目の射影だけから作られているので、$S_e/\sigma^2\sim\chi^2_{(a-1)(r-1)}$。
• 帰無仮説のもとで $S_A$ は2番目の射影だけから作られるので、$S_A/\sigma^2\sim\chi^2_{a-1}$。
• $S_R$ は3番目の射影と $\bm\gamma$ から作られ、どちらも2番目・4番目の射影と独立である。よって3つの平方和は独立である。

最後に、(11.7) で $u_j$ は独立同分布の正規変数なので

$$\frac{\sum_j(u_j-\bar u)^2}{\sigma_R^2+\sigma^2/a}\sim\chi^2_{r-1},\qquad\text{すなわち}\qquad \frac{S_R}{\sigma^2+a\sigma_R^2}\sim\chi^2_{r-1}$$

である。

定理 11.2 から、分散分析表は表 11.1 の形になる。

因子の検定 $H_0\colon\sigma_A^2=0$ には $F_0=V_A/V_e$ を使う。$F_0\ge F_{\phi_A,\phi_e}(\alpha)$ なら有意とする。この比が正確に $F_{\phi_A,\phi_e}$ 分布に従う理由は2つある。

• 期待値の差 $E(V_A)-E(V_e)$ はちょうど $r\sigma_A^2$ である。よって帰無仮説のもとで、分子と分母は同じ誤差分散を推定する。
• 分子と分母は独立なカイ二乗変数をそれぞれの自由度で割ったものである。

ブロックの検定 $H_0\colon\sigma_R^2=0$ にも $F_0=V_R/V_e$ を使う。ただしブロックは比べたい要因ではないので、この検定の目的は「ブロックを作ったことに意味があったか」の確認である。有意でなくても、計画の段階でブロックを作った以上、ブロックを誤差へプールせずに残すのが基本の扱いである。

分散分析後の推定

因子 $A$ が有意になったら、水準 $A_i$ の母平均と、2水準の母平均の差を推定する。ここで乱塊法に固有の注意がある。変量のブロックを考えるとき、水準 $A_i$ の母平均とは「ブロックの母集団全体で平均した値」$\mu(A_i)=\mu+a_i$ である。その点推定量は

$$\hat\mu(A_i)=\bar x_{i\cdot}=\mu+a_i+\bar\gamma+\bar\varepsilon_{i\cdot}$$

で、実験に使った $r$ 個のブロックの平均 $\bar\gamma$ を含む。よって

$$V[\bar x_{i\cdot}]=\frac{\sigma_R^2}{r}+\frac{\sigma^2}{r}=\frac{\sigma_R^2+\sigma^2}{r}.\tag{11.8}$$

一方、2水準の差ではブロック効果の平均が消える。

$$\bar x_{i\cdot}-\bar x_{i'\cdot}=a_i-a_{i'}+\bar\varepsilon_{i\cdot}-\bar\varepsilon_{i'\cdot}$$

なので

$$V[\bar x_{i\cdot}-\bar x_{i'\cdot}]=\frac{2\sigma^2}{r}.\tag{11.9}$$

Satterthwaite の等価自由度。

独立な平均平方 $V_1,\dots,V_k$ が、それぞれ $\phi_kV_k/E(V_k)\sim\chi^2_{\phi_k}$ を満たすとする。正の係数の1次結合 $\hat\theta=\sum_k c_kV_k$ の分布を、次のように近似する。

$$\frac{\phi^*\hat\theta}{E(\hat\theta)}\ \approx\ \chi^2_{\phi^*}.$$

自由度 $\phi^*$ は分散が一致するように決める。カイ二乗分布の分散は自由度の2倍なので、$V[V_k]=2E(V_k)^2/\phi_k$ であり、

$$V[\hat\theta]=\sum_k c_k^2\frac{2E(V_k)^2}{\phi_k} \quad\text{と}\quad V\Bigl[\frac{E(\hat\theta)}{\phi^*}\chi^2_{\phi^*}\Bigr]=\frac{2E(\hat\theta)^2}{\phi^*}$$

を等しいとおき、期待値を推定値で置き換えると

$$\phi^*=\frac{\bigl(\sum_k c_kV_k\bigr)^2}{\sum_k (c_kV_k)^2/\phi_k}\tag{11.12}$$

を得る。これを Satterthwaite の等価自由度と呼ぶ。数表を引くときは最も近い整数に丸める。

(11.11) では $c_1V_1=V_R/(ar)$、$c_2V_2=(a-1)V_e/(ar)$ である。共通の因子 $1/(ar)$ は分子と分母で消えるので

$$\phi^*=\frac{\{V_R+(a-1)V_e\}^2}{V_R^2/\phi_R+(a-1)^2V_e^2/\phi_e}$$

となる。分割法（第12章）でも同じ近似を使う。

ブロック因子の分散成分の推定

ブロック間のばらつき $\sigma_R^2$ は、それ自体が工程の情報になる（日による変動が大きいなら、日々の管理を見直す必要がある）。$E(V_R)-E(V_e)=a\sigma_R^2$ なので、不偏な点推定量は

$$\hat\sigma_R^2=\frac{V_R-V_e}{a}\tag{11.13}$$

である。$V_R<V_e$ のときこの値は負になる。分散は負にならないので、その場合は 0 と報告する（こうすると不偏性は失われる）。

近似区間（Satterthwaite）。

点推定量は2つの平均平方の差なので、係数に負のもの（$-V_e/a$）を含む。それでも近似の考え方は同じで、$\phi^*\hat\sigma_R^2/\sigma_R^2$ を近似的に $\chi^2_{\phi^*}$ とみなし、

$$\phi^*=\frac{(V_R-V_e)^2}{V_R^2/\phi_R+V_e^2/\phi_e},\qquad \frac{\phi^*\hat\sigma_R^2}{\chi^2_{\phi^*}(\alpha/2)}\le\sigma_R^2\le\frac{\phi^*\hat\sigma_R^2}{\chi^2_{\phi^*}(1-\alpha/2)}\tag{11.14}$$

を信頼率 $1-\alpha$ の近似信頼区間とする。等価自由度の式では分子に係数 $1/a^2$ が出るが、分母にも同じ係数が出るので消える。

負の係数を含む1次結合では近似の精度が落ち、$V_R$ が $V_e$ に近いほど悪くなる。逆にブロックの平均平方が誤差分散よりはるかに大きいときは、$\phi^*\approx\phi_R$ になる。このとき区間は、$\hat\sigma_R^2\approx V_R/a$ のカイ二乗区間に近づく。

分散比の正確な区間。

比 $\theta=\sigma_R^2/\sigma^2$ なら正確な区間が作れる。定理 11.2 から、次の2つは独立なカイ二乗変数である。

$$\frac{S_R}{\sigma^2+a\sigma_R^2},\qquad \frac{S_e}{\sigma^2}.$$

よって

$$\frac{V_R/(\sigma^2+a\sigma_R^2)}{V_e/\sigma^2}=\frac{V_R/V_e}{1+a\theta}\sim F_{\phi_R,\phi_e}.\tag{11.15}$$

$F_0=V_R/V_e$ と書き、次の不等式を $\theta$ について解く。

$$F_{\phi_R,\phi_e}(1-\alpha/2)\le \frac{F_0}{1+a\theta}\le F_{\phi_R,\phi_e}(\alpha/2).$$

すると

$$\frac{1}{a}\Bigl(\frac{F_0}{F_{\phi_R,\phi_e}(\alpha/2)}-1\Bigr)\le\theta\le\frac{1}{a}\Bigl(\frac{F_0}{F_{\phi_R,\phi_e}(1-\alpha/2)}-1\Bigr).\tag{11.16}$$

下限が負になったら 0 に置き換える。下側の点は $F_{m,n}(1-\alpha/2)=1/F_{n,m}(\alpha/2)$ で、上側の表から求める。

ブロック内相関は $\sigma_R^2/(\sigma_R^2+\sigma^2)=\theta/(1+\theta)$ と書ける。これは $\theta$ の単調増加関数なので、(11.16) の両端を変換すれば正確な区間になる。

確認テスト

オリジナル

1. 11.1-1 Fisher の3原則を挙げ、乱塊法でブロックを作ることとブロックの中で順番を乱数で決めることが、それぞれどの原則に当たるかを答えよ。
2. 11.1-2 構造式 (11.1) のもとで、$E(S_R)=(r-1)(\sigma^2+a\sigma_R^2)$ を導け。ブロック平均と総平均の差を構造式で表すところから始めよ。
3. 11.1-3 構造式 (11.1) のもとで、$\bar x_{1\cdot}$ と $\bar x_{2\cdot}$ の共分散を求めよ。それを使って次を示せ。

$$V[\bar x_{1\cdot}-\bar x_{2\cdot}]=2\sigma^2/r.$$ 4. 11.1-4 3水準・6ブロックの乱塊法で $V_R=20$、$V_e=2$ を得た。$\theta=\sigma_R^2/\sigma^2$ の点推定値と、信頼率 95% の正確な信頼区間を求めよ。数表の値 $F_{5,10}(0.025)=4.236$、$F_{10,5}(0.025)=6.619$ を用いてよい。5. 11.1-5 ブロックを無視して1元配置で解析したときの誤差分散 $V_e'=(S_R+S_e)/\{a(r-1)\}$ について、$E(V_e')=\sigma^2+\sigma_R^2$ を示せ。また、このとき $F_0'=V_A/V_e'$ を $F_{a-1,a(r-1)}$ 分布で判定する検定が正確でない理由を、データの共分散の構造から述べよ。