Myers, R.H., Montgomery, D.C. & Anderson-Cook, C.M. (2009)”Response Surface Methodology,” Third Edition.
Chapter 7. Experimantal Designs for Fitting Response Surfaces - I.


たとえば「おいしいグレープフルーツジュースをつくるためには、酸味と甘味をどうすればよいか？」というような問題があるとき、それを調べるアプローチには、大きく分けてふたつあります。ひとつは、いろんな
グレープフルーツジュースの酸味と甘味とおいしさを、手当たり次第に調べる方法です。もうひとつは、酸味と甘味のレベルを変えながら、いくつかのグレープフルーツジュースのサンプル品をつくり、それらのおいしさを調べる方法です。
後者のアプローチ、すなわち実験による製品最適化のための強力な武器となるのが、応答曲面法 (response

surface methodology) です。応答曲面法では、まず応答変数 (グレープフルーツジュースのおいしさ) を規定
する少数の独立変数 (酸味と甘味) を突き止めます。そのうえで、それらのレベルを変えたいくつかの刺激 (サンプル品) を用意して応答変数の値を調べ、応答変数を最大化 (ないし最小化) する組み合わせを統計的に推測します。なお、応答変数は 1 個、独立変数は 2〜3 個、いずれも連続的な変数であることがふつうです。
応答曲面法を理解するためには、回帰分析・分散分析の基礎知識と、実験計画法の基礎知識の両方が必要です。というわけで、定評ある教科書である Myers-Montgomery 本の最新版を読む勉強会を行ってまいりました。最終回の今回は、7 章の内容をご紹介します。この章では、「応答変数を最適化する組み合わせをうまく推測するためには、どんな刺激を用意すればよいか」という問題について、その基本的な考え方が紹介されています。
では、お楽しみください . . .



2013/04/09 RSM 読書会 (誤字訂正版). 報告: 小野 (株式会社インサイト・ファクトリー)

1

1 応答曲面計画が持つことが望ましい性質
応答曲面計画が持つことが望まれる、いくつかの重要な性質を以下に挙げる:

1. 結果としてモデルがデータによく当てはまること
2. 当てはまりの悪さの検定を可能にするだけの十分な情報を与えること
3. より高次のモデルを順に作ることができること
4.「純粋な」実験上のエラーの推定値を提供すること
5. データのなかに存在する外れ値に対して鈍感 (ロバスト) であること
6. 計画水準の制御におけるエラーに対してロバストであること
7. コストの面で効率的であること
8. ブロック化した実験が可能であること *1
9. 分散等質性という仮定が正しいかどうかチェックできること
10. V ar[ŷ(x)]/σ 2 の分布がよいこと *2
それぞれの性質の重要さは場合によって異なる。またこれらの性質は、往々にしてトレードオフの関係にある。

*1 本章 4.10 節を参照。
*2 本章 3.4 節を参照。

2

2 操作可能領域、関心領域、モデル不適切性6 章で示したように、応答変数 y と計画変数 x1 , x2 , . . . , xk との関係を示す真の関数は未知である。そこで私たち
は、E(y) = f (x, θ) である関数 f が、領域 R(x) において、一次なり二次なりの多項式で近似できると仮定する。
この R を関心領域 (region of interest) という。R は実験ごとに異なる。
いっぽう、装置や薬品といったものが働き、理論的には実験が可能であるような範囲 O(x) を考えることもできる。これを操作可能領域 (region of operability) という。
通常、R と O がどこかという知識を私たちは持っていない。Figure 7.1 を参照。ここで R は最適点の場所についての current best guess であり、手元の ŷ(x) が f (x, θ) のよい予測となっていると私たちが感じている領域である。

2.1 モデル不適切性とモデルバイアス
不適切な応答曲面モデルを指定したら何が起きるか。
私たちがE(y) = X1 β1と仮定し、しかしほんとうは

E(y) = X1 β1 + X2 β2と仮定した方がよかった、という場合について考えよう。このとき、β1 の推定量 b1 は次のように歪む:E(b1 ) = β1 + Aβ2
ただし A = (X1′ X1 )−1 X1′ X2 。この A を別名行列 (alias matrix) という。
モデルを当てはめた結果は次のように歪む。

E(ŷ) − E(y) = Rβ23ただし R = X1 A − X2 。

■バイアスと適合度検定との関係

この式は重要である。まず 6 章で述べた適合度検定 (lack-of-fit test) との関係2を示そう。残差の二乗平均 s についてE(s2 ) = σ 2 +β2′ R′ Rβ2p1

つまり、s2 のバイアスは、当てはめた値におけるバイアス Rβ2 の平方和に比例する。適合度検定とは、非ゼロのβ2 が存在するかどうかを検出する検定である。

■バイアスと実験計画の選定との関係

さきほどの式は実験計画の選択とも関わっている。それぞれの値 yj について、私たちは損失 E[ŷj − f (xj , β)] を最小化したいと考えているわけである。ここで2E[ŷj − f (xj , β)]2 = V ar(ŷj ) − [Bias(ŷj )]2である。つまり、損失は分散誤差 V ar(ŷj ) とバイアス誤差 [Bias(ŷj )]2 に分かれる。さきほどの式とは∑β2′ R′ Rβ2=[Bias(ŷj )]2
p1
j=1という関係がある。

■バイアスを小さくすることも大事だ

統計学の本は実験計画の選択にあたってたいてい V ar(ŷj ) のことばかり述べている。本書の以下の説明もそうである。これは指定したモデルが正しい場合について考えているからだ。でもほんとは [Bias(ŷj )]2 を小さくすることも大事である。

4

3 一次モデルのための実験計画
計画変数 x1 , x2 , . . . , xk 、応答 y 、ランの数を N とする。各変数の水準は −1 から +1 までの範囲内だとする。以下のモデルを当てはめる:ŷi = b0 + b1 xi1 + b2 xi2 + · · · + bk xik
推定される係数の分散は小さいほうが望ましい。この観点に基づく実験計画の選択基準として一番わかりやすいのは V ar(bi ) であろう。これを最小化する計画を選択するために重要になるのが、直交性という概念である。

3.1 一次直交計画
X ′ X が対角行列になる計画のことを一次直交計画という。つまり、X の列がお互いに直交しているような計画である *3 。次の性質がある:

一次モデルと固定された標本 N があるとき、もし j = 1, 2, . . . , k に関して Xj ∈ [−1, +1] ならば、j =1, 2, . . . , k において V ar(bi )/σ 2 が最小化されるのは、計画が直交しており、かつ計画に含まれるすべての xiの水準が i = 1, 2, . . . , k において ±1 のときである。



例. 計画変数が 3 つ、ランが 8 つのとき、23−1 一部実施要因計画 (分解能 III) *4 が用いられる。この計画の X
は



1
 1

 1

 1
X=
 1

 1

 1
1

1
1
−1
−1
−1
−1
1
1


−1 −1
−1 −1 1 −1 
1 −1 

−1
1 

−1
1 

1
1 
1
1

列は直交しており、すべての水準が ±1 である。従ってこの計画は、モデル

ŷ = b0 + b1 x1 + b2 x2 + b3 x3
に対して分散最適である。つまり、観察事例数が同じなら、係数の分散は 23 完全実施要因計画と等価であり、これより小さくなる計画はない。

3.2 交互作用を含むモデルのための直交計画
一次モデルの場合と同じく、直交計画が分散最適である。
*3*4

X の列のうちどの異なる 2 列を取りだしても、積の合計が 0 になる、ということ。
実験計画の分野での独特の書き方で、2 水準の変数 3 つのすべての組み合わせ (23 個) について調べる実験計画を「23 計画」、そのうち半分だけをうまいこと選んで調べる実験計画を「23−1 計画」と書きます。その選び方が持っている性質を分解能という考え方で表現します。ここではあまり気にしなくていいと思います。

5

例. 計画変数が 5 つあり、モデルが下式であるとしよう *5 :yi = β0 +5∑

βi xi +

i=1

5∑∑

βij xi xj + ϵ

i<j=2

25 完全実施要因計画、ないし、25−1 部分実施要因計画 (分解能 V) がここでの直交計画になる。後者の X を示す:

■飽和計画

上の例では、モデルの項は 16 個、計画点も 16 個である。すなわち、残差の自由度がゼロになっている。こういう計画のことを飽和している (saturated) という。
このとき、適合度検定はできない。つまり、モデルが不適切かどうかを検証するための情報がない。データに対する R2 は 100% になる *6 。

■中心ランの効果 直交計画に中心ラン (全計画変数の値が 0 のラン) を追加しても、直交性は失われないが、分散
最適とはいえなくなる。
主効果や交互作用が線形だと考えている限り、中心ランを加えても得られる情報は何もない。ただし、二次効果の検出はとても大事だし、中心ランの追加はそのための役に立つ。
*5

記号がわかりにくいけど、

yi = β0 + β1 x1 + (3 個省略) + β5 x5 + β12 x1 x2 + (8 個省略) + β45 x4 + x5 + ϵ*6という意味です。下の行列 X の列見出しを参照。
この本を通じて、計画点あたり 1 個しかデータを収集しない状況が想定されていることに注意してください。

6

■どんな計画を用いるべきか?

分散最適な計画が常に最良とは限らない。計画変数が 3 つあって、次の一次モデル

をあてはめる場合について考えよう:

yi = β0 + β1 xi1 + β2 xi2 + β3 xi3 + ϵi(i = 1, 2, . . . , 8)

Figure 7.2 の 3 つの計画を比較してみよう。計画 1 と計画 3 は分散最適である。回帰係数の分散は σ 2 /8 となる。
いっぽう計画 2 では、切片項の分散は σ 2 /8, ほかの係数の分散は σ 2 /4 である。

• 計画 1 では適合度検定ができない。この計画は、モデルが適切であるという強い先行情報があるときに適切である。

• 計画 2 では一次効果の検定ならびに二次効果の検定ができる。いっぽう、一次項の推定は計画 1・計画 3 と比べて不正確である。この計画は、二次項があるかもしれないが交互作用は絶対にないとわかっているときに適切である。

• 計画 3 ではモデルの適切性はわからない。一次効果の検定はできるが交互作用効果や二次効果の検定はできない。二次効果や交互作用効果が絶対にないとわかっているときに適切。

3.3 そのほかの一次直交計画 - シンプレクス計画
飽和している一次直交計画のなかに、シンプレクス計画という特殊なクラスがある。
シンプレクス計画では、k 個の変数を持つ一次モデルをあてはめるために N = k + 1 個の観察事例を必要とする。
計画点は、k 次元空間で原点からみて 2 つの計画点がなす角度 θ がcosθ = −1kとなるように置く。k = 2, N = 3 なら Figure 7.3 のように正三角形になる。k = 3, N = 4 なら Figure 7.4 になる。これは 23−1 一部実施実験計画 (分解能 III) でもある。

■シンプレクスの回転 シンプレクス計画を原点を中心にぐるぐる回転させても一次直交であることはかわらない。
つまり、所与の k についてシンプレクス計画は無数にある。

7

3.4 もうひとつの分散特性 - 予測分散
ここまでは係数の分散に注目してきた。実験計画についてかんがえる際のもうひとつの重要な基準は、予測分散である。

■一次モデルの予測分散

予測値の分散は下式であらわされる:′P V (x) = V ar[ŷ(x)] = σ 2 x(m) (X ′ X)−1 x(m)ただし、x(m) はモデルによって決まる値である。一次モデルなら′

x(1) = [1, x1 , x2 , . . . , xk ]交互作用を含む k = 2 のモデルなら′

x(1) = [1, x1 , x2 , x1 x2 ]である。P V (x) が計画空間の場所によって変わる点に注意。

8

■尺度化された予測分散計画を比較する際には、予測分散を次のように尺度化しておくと都合がよい。

SP V (x) =′N V ar[ŷ(x)]= N x(m) (X ′ X)−1 x(m)σ2

分散最適な一次モデルでは、x と原点との距離を ρx として、

SP V (x) = 1 + ρ2xとなる。このように、原点から離れるにつれて SPV は大きくなる。分散最適でない計画ならばなおさらそうである。

9

4 二次モデルをあてはめるための計画
応答曲面法にとって欠かせないフェイズのひとつが変数スクリーニングである (1 章で詳述した)。そこでは 2 水準の要因計画と部分実施 (fractions) が主要な役割を演じる。その次の領域探索のフェイズ *7 も一次の計画でおこなわれる。
どこかの時点で、リサーチャーは計画変数に二次の応答曲面をあてはめることに関心を持つだろう。ここで必要なのは、ユーザが下式の二次モデル *8 をあてはめることを可能にしてくれるような実験計画である。

y = β0 +k∑i=1βi xi +k∑

βii x2i +

i=1

k
∑∑

βij xi xj + ϵi<j=2

この式は 1 + 2k + k(k − 1)/2 個のパラメータを持っている。

6 章で述べたように、二次の応答曲面のあてはめを行うためには、実験計画は次の特性を持っていなければならない:

1. それぞれの計画変数について、少なくとも 3 水準ある。
2. k 個の計画変数に対して、少なくとも 1 + 2k + k(k − 1)/2 個の異なる計画点がある *9 。
ただしこれは最低限の条件である。本章の最初に挙げた 10 個の性質を忘れないこと。
本章の以下の部分と 8 章では、二次モデルの重要な性質と、二次計画のいくつかのクラスについて述べる。二次計画の場合、計画の直交性はあまり大事な話ではなくなり、個々の係数の推定よりも、尺度化された予測分散のほうがもっと大事になる。これは、どの変数をモデルに含めるかということよりも、y の予測 ŷ(x) の良さに関心が向けられるからである。

4.1 中心複合計画
二次の計画のなかでいちばん人気があるのは中心複合計画 (central composite designs; CCDs) である。この計画は、次の 3 種類の計画点からなる。計画変数の数を k 個として、• 2 水準要因計画から得られる計画点。完全実施計画、ないし分解能 V の部分実施計画を用いる *10 。一次の項と 2 要因の交互作用効果の推定に寄与する。中心複合計画は、これらの項に関しては分散最適となる。

• 軸上の 2k 個の計画点。ある点は、ある計画変数について α ないし −α、他の計画変数について 0 とする。二*7 5 章で述べられています。
*8

記号がまたわかりにくいけど、計画変数が 2 つだとして、

yi = β0 + β1 x1 + β2 x2 + β11 x21 + +β22 x22 + β12 x1 x2 + ϵという意味です。
*9 計画変数が 2 つだったら 6 個, 3 つだったら 10 個ってことですね。
*10

ここで「分解能 V の部分実施計画」といっているのは「各要因の主効果および 2 要因交互作用がすべて推定可能なら部分実施計画でもいいよ」という意味です。k = 2 ないし k = 3 の場合には完全実施計画でなければならないです。

10

次の項の推定に寄与する。

• nc 個の中心ラン。pure error の推定 *11 と二次の項の推定に寄与する。
Figure 7.5 に k = 2, α =√√2 の例を、Figure 7.6 に k = 3, α = 3 の例を示す。中心複合計画では α と nc の選択が重要である。

4.2 計画積率と回転可能性■計画積率

実験計画の性質の多くは、実験領域における計画点の散らばりかたによって決まる。計画点の分布を表すために、計画積率 (design moments) という考え方を用いる。



1
 1

X= .
 ..

x11
x12
..
.

1 x1N
*11


x21 . . . xk1
x22 . . . xk1 .. . .
.. .
.
. 
x2N . . . xkN

2 章に説明がありましたね。測定誤差など、モデルと無関係な誤差のこと。

11

について、以下のように定義する *12 :• 一次積率 (first moments): [i] = N1∑N

u=1 xiu∑N• 二次純積率 (second pure moments): [ii] = N1 u=1 x2iu
∑N
• 二次混合積率 (second mixed moments): [ij] = N1 u=1 xiu xju ただし i ̸= j

上から順に、標本平均、標本分散、標本共分散に相当する。
少なくともひとつの変数が偶数次である積率を奇積率 (odd moments)、そうでない積率を偶積率 (even moments)と呼ぶ。一次積率と二次混合積率は奇積率、二次純積率は偶積率である。

■積率行列 下式であらわされる行列を積率行列 (moment matrix) と呼ぶ。

M=

X ′X
N

一次モデルのための 2k 完全実施要因計画ないし部分実施要因計画は、直交であるから、積率行列はもちろんこう
なる:



1 0

′
XX
 0 1
M=
= .
N
 ..
0

... 0
..

.




 = Ik

1

2k 完全実施要因計画ないし部分実施要因計画では、M は二次積率までを含んでいる。中心複合計画の場合、Mは四次積率までを含む (下記)。

*12

要するに、ある計画変数についてすべてのランを見通したとき、値の平均を一次積率、値の二乗の平均を二次純積率、ほかの計画変数の値との積の平均を二次混合積率と呼んでいます。

12

■積率行列はなぜ大事か 実験者はどこに最適点があるのか知らない。だから、将来の応答の値についての予測は、関心領域のどこでもだいたい同じ質を持っていないと困る。つまり、尺度化された予測分散が安定していないと困る。従って、二次の計画においては、尺度化された予測分散が実験計画領域を通じて安定した分布を持っていることが大事である。

3.4 節で、尺度化された予測分散を下式のように定義した:′

SP V (x) = N x(m) (X ′ X)−1 x(m)
ここから、
(m)′

SP V (x) = x

(

X ′X
N

)−1

′

x(m) = x(m) M −1 x(m)である。このように、積率行列は尺度化された予測分散に大きな影響を与える。

■回転可能性 次の性質のことを、実験の回転可能性 (design rotatability) という。





計画の中心から等距離な任意の 2 点において N V ar[ŷ(x)]/σ 2 が等しい計画を、回転可能な計画と呼ぶ。回転可能な計画では、N V ar[ŷ(x)]/σ 2 が球の上で一定である。



この性質それ自体は予測分散が安定していることを意味しないけれども、この性質は計画のパラメータ (たとえば、
中心複合計画の α と nc ) を選ぶ際のガイドラインを提供してくれる。それに、回転可能性 (ないし、それに近づくこと) はふつう、他の重要な性質を犠牲にすることなく、簡単に達成できる。
というわけで、以下では計画の回転可能性の必要十分条件を決定するための基礎について述べる。詳細はAppendix 1. で述べることにし、ここでは結論だけ紹介する。

■回転可能な一次計画の積率行列

まず一次計画の場合について考えよう。計画が回転可能になる必要十分条件は、「二次までの奇積率が 0 で、二次純積率がすべて等しい」である。言い換えると、• [i] = 0(i = 1, 2, . . . , k)

• [ij] = 0(i ̸= j, i, j = 1, 2, . . . , k)

• [ii] = λ2(i = 1, 2, . . . , k)

λ2 は計画の尺度化で決まる値で、分解能 III 以上の 2 水準計画でいえば、すべての水準が ±1 であれば λ2 = 1 である。
一次の計画の場合、計画が回転可能であるための条件は、計画が分散最適であるための条件と等しい。あまり驚いてはいけない。分散最適な一次モデルでは x と原点との距離を ρx として

SP V (x) = 1 + ρ2xとなることを思い出してほしい。この式は計画が回転可能であることを示している。

13

■回転可能な二次計画の積率行列 二次計画の場合、回転可能性の必要十分条件は、「四次までのすべての奇積率が0 で、[iiii]/[iijj] (i ̸= j) が 3」である。ややこしい条件だが、達成するのは比較的に容易である。

4.3 回転可能性と中心複合計画
中心複合計画の場合、四次までのすべての奇積率は 0 なので、問題は [iiii]/[iijj] だけである。要因計画から来た計画点の数を F とすると (完全実施要因計画なら F = 2k )、F + 2α4[iiii]=[iijj]
Fとなるので、これを 3 にするためには α =k = 2 のとき α = 22/4 =√4F とすればよい。中心ランの数とは無関係である点に注意。

√
2 = 1.414, k = 3 のとき α = 23/4 = 1.682 である *13 。k = 4 までであれば、要するに、(中心ランを別にして) だいたい球形に計画点をおけばよいことになる。

■回転可能な中心複合計画における中心ランの数回転可能な中心複合計画を用いるとき、中心ランの数を少なくするのはあまりよくない。Figure 7.8, Figure 7.9 を参照。k = 2 の中心複合計画で、前者は nc = 1 のとき、後者
は nc = 5 のときの、尺度化された予測分散を示している。N は後者のほうが大きいのに、N V ar[ŷ(x)]/σ 2 は概して後者のほうが小さい。中心のあたりでは約 3.5 倍も違う。
この性質は、k の値が大きくても変わらない。詳細は 8 章で述べるが、(中心ラン以外の) 計画点がだいたい球形の上に載っている計画では、中心ランは 3 個から 5 個程度必要である。

■回転可能性は重要か

二次計画において、計画が厳密に回転可能であることは必須ではない。実のところ、次の性質が成り立つ。



計画領域として球状の領域が望まれているとき、分散の観点から見てもっとも効率的な中心複合計画は、

α=√k, 中心ランが 3〜5 個の計画である。



この計画は必ずしも回転可能ではないが、ほぼ回転可能である。また、厳密に回転可能な計画よりも少し性質が良い。その理由については 8 章で述べる。

*13√

読書会でご指摘いただいて気がついたのですが、k = 3, α = 3 の中心複合計画が厳密には回転可能でないという点、不思議ですね。
中心ランを除くすべての計画点が球の表面に載っているのに...。
√よくわからないのですが、k = 3, α = 3 の場合、12 個の計画点は、正六面体の 6 枚の面がそれぞれピラミッド型に外側に張り出したような形をつくっていて、24 面体だけど「正 24 面体」ではない。中心から見ると、12 個の点は中心から等距離だけれど、球の上にまんべんなく広がっているのではなく、ちょっと濃淡がある。だから回転可能ではない。ということなのかなあ、と思うのですが...

14

4.4 予測分散についてのさらなる話題
3.4 節で述べたように、予測分散は下式で定義される:′

P V (x) = V ar[ŷ(x)] = σ 2 x(m) (X ′ X)−1 x(m)■計画選択の段階

計画を選択する段階では σ 2 についてはわからないし、わかったところで計画の良し悪しを考える際には役に立たない。むしろ、′

SP V (x) = N x(m) (X ′ X)−1 x(m)15を使う方がよい。4.2 節で述べたように、SPV は積率行列

M=

X ′X
N

で求められる (積率行列のことを情報行列ともいう)。ここで、N はコストに相当し、積率行列はコストあたりの結果の品質に相当する。
しかし、実験計画のコストは N に比例するとは限らない。たとえば、実験をやること自体が大変で、ランを増やすのはそんなに大変でない、という場合もある。こういう場合は、できるだけ大きな実験をやったほうが安上がりだということになる。こういう状況では、計画を比較する際に、SPV ではなく、尺度化されていない予測分散′

U P V (x) = V ar[ŷ(x)]/σ 2 = x(m) (X ′ X)−1 x(m)を使ったほうがよいかもしれない。どちらを使うべきかは議論が分かれるところである。

■実験データ解析の段階 実験が済んだあとでは、「この計画から得られたモデルはどのくらい良いモデルか」に焦点が当てられる。SPV や UPV にはもう意味がない。σ 2 の推定値 (MSE) がもうわかっているのだから、UPV に

MSE を掛けるべきである。これを推定された予測分散という。
′

EP V = M SE x(m) (X ′ X)−1 x(m)この平方根が、ご存じ「予測の標準誤差」である。

4.5 立方体状領域と面心立方
関心領域と操作可能領域が同一で、計画領域がぴったり立方体の形をとる場合がある。これを立方体状計画領域(cuboidal design region) という。
この場合には、回転可能性は無視して、計画点をぎりぎり外側にするのが効率的である。二次計画の場合は、

α = 1 の中心複合計画を用いるのが良い。これを面心立方 (face-centered cube, FCD) と呼ぶ。Figure 7.11 を参照。

16

この場合、予測分散の観点からは、中心ランは 1 個か 2 個あれば十分である。Figure 7.12 は 0 個の場合、Figure

7.13 は 1 個の場合、Figure 7.14 は 3 個の場合の SPV を示している。たいして変わらないことがわかる。ただし、pure error の推定という観点からは、やはり中心ランは複数個あったほうがよい。

4.6 計画領域が球状なのはどういうときか操作可能領域と関心領域は往々にして異なる。いま、立方体領域に関心があり、操作可能領域がその少し外側に出ているとしよう。このときは、中心複合計画の α を調整して、「実用的な α 値」にするという方法がある。その値は関心領域の外側に落ちる場合もある。
関心領域・操作可能領域が明確でない場合も多い。しかし、計画のタイプを混同してしまうことがあるからといって、実験計画を用いない言い訳ににはならない。立方体領域を使うのが自然なときに球状領域を使ってしまっ
たからといって、依然として重要な情報が手に入るし、今後の実験でもっと良い領域を選ぶ助けになる。

4.7 中心複合計画についてのまとめ
中心複合計画は効率のよい計画であり、逐次的な実験に最適である。典型例として次の 2 つがある。

• 球状計画 (α =√k)。各要因について 5 水準を使う。(ほぼ) 回転可能である。中心ランは 3 個から 5 個必要である。

• 立方体状計画 (α = 1)。各要因について 3 水準を使う。中心ランは 1 個か 2 個でよい。
そのほかにもさまざまなバリエーションがある。
中心複合計画の競争相手は 3k 要因計画である。k = 2 の場合は面心立方であり、すなわち中心複合計画である。

k = 3 で立方体状領域の場合は、もし実験者が計画点を 27 個を作れるならば有用かもしれない。k > 3 の場合は現実的でない。

4.8 ボックス・ベンケン計画
Box & Behnken (1960) は、二次の応答局面を当てはめるための 3 水準計画をつくるための、きわめて創造的な方法を提案した。
この考え方は釣り合い型不完備ブロック計画 (bakanced incomplete block design, BIBD) に基づいている。たとえば、処理数 3, ブロック数 3 の BIBD は下表である。処理 1 と 2 がペアになっているブロック (ブロック 1) では、処理 3 はない。

17

18

これと同じように、計画変数 x1 と x2 をペアにして 22 要因計画をつくり (水準は ±1 に尺度化する)、そこでは計画変数 x3 は中心 (0) とする。これを計画変数 x1 と x2 のペア、x2 と x3 のペアについて繰り返す。こうして、

k = 3 のボックス・ベンケン計画 (Box-Behnken design) ができあがる。最後の行は中心ランのベクトルを表す (つまり、中心ランは複数個あってよい)。

同様のやりかたで、k = 4 のボックス・ベンケン計画をつくることができる。k = 2 の場合はつくれない。

k = 3, k = 4 におけるボックス・ベンケン計画の計画点の数は 12 + nc , 24 + nc である。中心複合計画では
14 + nc , 24 + nc だから、あまり変わらない。

■ボックス・ベンケン計画の特徴

ボックス・ベンケン計画は中心複合計画にならんで効率的な計画である。以下の特徴がある。

• 水準間の間隔が等しい。
• ほぼ回転可能性がある (k=4,7 の場合は厳密に回転可能である)。
• 実は球状計画である。Figure 7.17 を参照。
• 中心ランとして 3 個から 5 個が必要である。

19

4.9 そのほかの球状計画：等放射状計画
2 要因実験の特別なクラスとして、等放射状計画 (equiradial designs) がある。同心円上に等間隔に実験点を置
き、さらに中心ランを置く計画である。五角形計画、六角形計画、. . . がある。すべて回転可能である。

Figure 7.20 に六角形計画の例を示す。なお、八角形計画は中心複合計画である。

■等放射状計画の特徴

計画変数の数が 2 の場合でも、本来は中心複合計画を使ったほうが良いが、コスト上の観点から難しいときには、等放射状計画が候補になる。その場合は六角形計画が望ましい。七角形計画が有用な場合もある。なお、五角形計画は飽和デザインなので、よほどのことがない限り用いてはならない。

20

中心ランは 2 個から 4 個必要である。

4.10 二次計画における直交ブロック化■直交ブロック化とは応答曲面法の実用場面では、実験が大きくなりすぎてしまうことが多い。そんなときは、実験をブロック化できると都合がよい。つまり、係数への影響を最小にするかたちで、ランをブロックに分離できると都合が良い。そのためには、ブロックによる効果が係数と直交であるような形でブロック化する必要がある。
これを直交ブロック化という。

■直交ブロック化の例

例として、22 要因計画をつかって一次モデルを当てはめる場合について考える。I = AB

を定義関係として 2 つの 1/2 の部分をつくるとしよう *14 。モデルは

yi = β0 + β1 xi1 + β2 xi2 + δ1 zi1 + δ2 zi2 + ϵi
ここで δ1 , δ2 はブロック効果の係数である。zi1 , zi2 はダミー変数で、yi がブロック 1 にあるとき zi1 = 1, そうでないとき zi1 = 0 となる。X 行列は

これではうまく推定できないので、z2 を除去し、z1 を中心化する。モデルは

yi = β0 + β1 xi1 + β2 xi2 + δ1 (zi1 − z̄1 ) + ϵi
X 行列は

*14

これも実験計画の分野に特有の言い回しだと思います。4 つのランを、2 つの要因を掛けた値が +1 になるか −1 になるかで 2 つにわける、という意味です。

21この例で、ブロック効果が回帰係数と直交しているということは
4∑(zu1 − z̄1 )xu1 = 0

u=1
4
∑

(zu1 − z̄1 )xu2 = 0

u=1ということである *15 。この計画では実際にそうなっている。つまり、この例は直交ブロック化である。

■二次計画での直交ブロック化の条件

計画変数 k 個の二次モデルに b 個のブロック変数を組み込んでみよう。モデルはyu = β0 +k∑

βi xui +i=1k∑βii x2ui +

i=1

k
∑∑

b
∑

βij xui xuj +δm (zum − z̄m ), u = 1, 2, . . . , N

m=1

i<j=2

zmu は、u 番目の観察事例が m 番目のブロックにはいっているときに 1 となるダミー変数である。このままでは推定できないが、その問題は脇においておく。
ブロック効果が回帰係数と直交しているということはN∑

xui (zum − z̄m ) = 0, i = 1, 2, . . . , k, m = 1, 2, . . . , b

u=1
N
∑

xui xuj (zum − z̄m ) = 0, i ̸= j, m = 1, 2, . . . , b

u=1
N
∑

x2ui (zum − z̄m ) = 0, i = 1, 2, . . . , k, m = 1, 2, . . . , bu=1ということである。ここで、積率 [i] = 0, [ij] = 0 だとしよう (中心複合計画、ボックス・ベンケン計画、2 水準要因
計画のすべてでこれが成り立つ)。このとき、ひとつめとふたつめの式は、それぞれ「各ブロック内で「各ブロック内で∑∑xui = 0」

xui xuj = 0(i ̸= j)」ということを意味する。あわせると、「各ブロックが一次直交計画でなければならない」ということを意味する。

3 番目の式はどうか。この式は、「ブロック m において

∑N2
u=1 xui zum = ẑm∑N2
u=1 xui 」と書き換えられる。つまり、「それぞれの計画変数において、その平方和に対する各ブロックの寄与は、ブロックのサイズに比例していなければならない」ということを意味する。
まとめると、



二次モデルの直交ブロック化計画は、以下の条件を必要とする。

• それぞれのブロックが一次直交計画を構成していること• それぞれの計画変数において、その平方和に対する各ブロックの寄与が、ブロックのサイズに比例して
いること





*15 原文では x の添字が x1u , x2u となっていますが、誤記だと思います。

22

■中心複合計画における直交ブロック化 k = 2 の回転可能な中心複合計画、2 ブロック、中心ランが各ブロックにつき 2 個、という場合について考えてみよう。下の計画 (α =k = 3 で 3 ブロックなら、下の計画で α =√2) が直交ブロック化になっている。

√
3 とすれば直交ブロック化になる。

この 2 つの例からわかるように、次の原則が成り立つ。2 ブロックの場合、いっぽうを要因のブロック、もういっぽうを軸のブロックにする。3 ブロックの場合は、要因の部分を 2 ブロックにわけ (それぞれのブロックが直交計画になるようにする)、軸の部分を 1 ブロックにする。
中心複合計画で可能な直交ブロックの数は、2, 3, 5, 9, 17, . . . である。k = 2 のときは 2 ブロック, k = 3 のとき

3 ブロックだけが可能である。軸の部分はブロックに分けてはいけない。

■立方体状の中心複合計画における直交ブロック化

立方体状の中心複合計画 (α = 1) の場合は、かなりサイズのちがうブロックをつくらないと直交しない。k = 2 の場合は23

k = 3 の場合はもっと大変である。したがって、ブロック化する場合は球状計画を使ったほうが良い。

■ボックスベンケン計画における直交ブロック化

k = 4, k = 5 のボックス・ベンケン計画は直交ブロック化できる。k = 4 の場合を示す。

■等放射状計画における直交ブロック化

等放射状計画は直交ブロック化できる。六角形計画の場合、点 1 と 3 と5, 点 2 と 4 と 6 をそれぞれブロックにし、各ブロックに同じ数の中心ランを加えればよい。

■ブロック化実験の分析形式

さきに述べたモデルは、ブロックが加算的な効果しか持たないことを前提としている。応答曲面法におけるブロック化とは、ブロックが応答の切片にしか効果を持たない、という考え方の上に成り

24

立っている。もしブロックが他の項と交互作用するならば *16 、ブロックごとに別々に応答局面を当てはめる必要がある。解釈は難しくなる。

おしまい

*16

つまり、いずれかの回帰係数が、ブロックのあいだで異なる、と考えられる場合。

25