解説/Review

特集：品質・信頼性・安全性確保に関わる横幹的体系

田口玄一の統計学に対する貢献
椿広計∗1

Gen-ichi Taguchi’s Contribution to Statistics
Hiroe TSUBAKI∗1
Abstract– The author introduces and appreciates original statistical methods proposed by Gen-ichi
Taguchi from 1957 to 1977 such as his shrunken estimators, nonparametric inferences to clarify
significant differences of probability distributions and unique design of randomized experiments.
The author also reviews Dr. Taguchi’s classification of input variables to clarify their roles in thede-signed experiments.
Keywords– Classification of Factors, Robust Parameter Design, Ridge Regression, Semiparametric
Poisson Modelling, Randomized Design

1. はじめに：工学を横断する実験統計学

一方，田口は，ロバスト・パラメータ設計を確率モデル
が前提となる推測統計学的方法の範疇には入らないと考筆者は，1975 年東京大学教養学部で田口玄一の「統えた．そのため，ロバスト・パラメータ設計を Fisher[1]計学」受講した．筆者の講義ノートによれば，田口は，の統計的実験計画法の優れた拡張として解釈しようとす1975 年 4 月 22 日「工学を横断する工学がある」と講義した．情報収集効率を最大化する Fisher の実験統計学，
情報伝達効率を最大化する Shannon の通信理論，その複合技術としての Wiener のサイバネティックスをその例とした．その上で，１年間実験統計学の考え方を紹介するとした．
田口は，Fisher が創生した「実験計画法」[1] を新製
品・新技術開発の効率を大きく向上させる工業実験計画法に進化させた．そもそも産業界における実験計画法応用の目的は製品の特性をデータに基づいて効率的に最適化することである．田口は，1950 年代から，直交表や線点図を用いた多因子部分実施実験を日本の産業界に投入し，新製品開発実験の生産性を欧米で常態化していた一因子実験の組み合わせの数十倍にした．1980 年代以降，
海外でタグチメソッド，あるいはロバスト・パラメータ設計と呼ばれる，独創的設計科学体系を考案し，国際的に実践した．これを通じて，欧米で最もよく知られた日本の工学者となり，本田宗一郎，豊田英二についで，日本人で 3 人目のアメリカ自動車殿堂入りを 1997 年に果たした．

る国際的統計コミュニティと，田口同様それを統計的方∗1 情報・システム研究機構

法の枠組みを超えたものとする国際的工学コミュニティの間でも大きな論争を巻き起こした [2]．現在，ロバス
ト・パラメータ設計には，ISO 規格も存在する [3]．
統計学コミュニティは，田口の統計的方法をより洗練したものとしようと試みたのに対して，機械工学などのコミュニティは，統計学が，田口が示した工学原理を十分理解していないことを批判した．1992 年，筆者はミズーリ大学設計生産性センターを訪問した. センター長であり米国機械工学会最適化工学分野の重鎮だった

Kenneth Ragsdell が，「田口のエンジニアリング哲学は，
自分の生涯で最も偉大なもの」と述べた [4]．
現在，日本の学術コミュニティでは，ロバスト・パラ
メータ設計を創始し，推測統計学から決別する前に，田口がどのような研究構想を統計学に対して提起したかは，殆ど知られていない．その理由は，田口が国際的に著名になり，ロバスト・パラメータ設計に転向した 1980
年代以降，自身が提案したユニークな統計的方法の進化に対して関心を失ったこと，田口の 1960 年代を中心とした統計的方法の提案が，当時の技術者が計算可能にするために二次形式の分解を中心とした近似的方法だった

統計数理研究所

東京都立川市緑町 10-3
Institute of Statistical Mathematics, 10-3 Midori-cho,Tachikawa-shi, Tokyo∗1 Theこと，このため，田口の著書が 1980 年代に翻訳され海外に紹介されたときに，海外の統計学者が，田口の統計的方法を誤りとし，統計的方法の細部に対する貢献を否

Received: 21 July 2019, Accepted: 16 August 2019.

110

横幹 第 13 巻 第 2 号

Gen-ichi Taguchi’s Contribution to Statistics

定したことが挙げられる．

示因子，ブロック因子，補助因子（共変量）を，古典的しかし，今日統計科学や統計的機械学習の中で用い実験計画法に対する重要な貢献としている．その後，田られる基本原理の幾つかは，田口が 1950 年代から 1970

口は，信号因子，調整因子，座標因子（ω （オメガ）因年代末まで，出版・改訂を続け，英訳もされた「実験計子）といった概念も導入した．

画法」[5, 6, 7] に暗示されている．1950 年代の初版 [5]因子の分類を通じて，製品・技術開発で，どのようには，比較的 Fisher の実験計画法に則しつつ，随所に工学

統計的モデルを構築すべきかの指針が与えられた．竹内的センスの投入が開始されたもの，1960 年代の第 2 版[10] は，田口が「天上の実験計画法を地上に引きずりおろした」という逆説的評価も行ったが，1990 年代米国
で，田口が Fisher の後継者とも評価されたのは，種々批判を受けた実験計画法の分析技術的側面ではなく，戦略的側面に対する独創的貢献からと考える．
以下，簡単に田口の様々な因子を再定義する．

（新版）[6] は，田口のオリジナルな統計的方法にあふれた内容で比較的数理的内容も解説されているもの，そして，最後の版となり，英訳もされた 1970 年代の第 3 版[7] は，統計理論から離陸し，ロバスト・パラメータ設計に向かうアイデアが示されたものとなっている．
本解説は，筆者が重要と考える田口の統計学に関する提案の一部を紹介したい．なお，品質工学の基幹となる
ロバスト・パラメータ設計の基本原理に関する筆者の見解は，椿，河村 [8] に詳述したが，本特集でも，立林和夫氏が紹介しているので，参考にされたい．

b) 制御因子 (Control Factor)：意図的にその水準（変数値）を制御可能で，その効果を推定すべき因子．
b-1) 信号因子 (Signal Factor)：その水準の変化が，比例性，線形性など出力変数と特定の強い関係性にあることが既知の制御因子．

2. 出力最適化戦略を導いた因子の分類
2.1 古典的実験計画法における因子の分類「因子 (Factor)」とは，「出力変数（被予測変数あるいは応答変数）と呼ばれる制御対象となる変数の変動に，系統的影響を与える可能性があると考えられるために観測する入力変数」の実験計画法分野での呼称である．

Fisher[1] の実験計画法は，要因効果を評価したい制御因
子と，要因効果には関心がないが，出力変数の変動に影響を与えるブロック因子を区別した．例えば，農業実験における品種は制御因子，どの畑で栽培したかはブロック因子であり，第 i 品種 (i = 1, . . . , a) を第 j 畑 ( j = 1, . . . , b)
で栽培した時の収量を Yi j とすれば，「一元配置乱塊法
実験（One way Randomized Block Design）」の分散分析
(ANOVA, analysis of variance) モデルと呼ばれる加法モ
デルは，(2.1) となる．
Yi j = µ + αi + β j + εi ja) 因子：その水準（変数値）が出力変数の期待値やバラつきに対して影響を与える可能性のある入力変
数．

(2.1)

ここで，µ は総平均（general mean，通常，一般平均と呼称）を示す未知母数，αi ，β j は，それぞれ制御因子，
ブロック因子の「主効果 (main effect)」を示す未知母数
で，識別可能とする線形制約条件が導入される．また，

εi j は誤差変動であり，実験順序を無作為化することで，確率変数と見なす．単純な実験繰り返しを排除し，系統的変動要因であるブロック因子を通じて実験を行うのが，Fisher の実験計画法の基本原理である「局所管理の原則」である．

b-2) 調整因子 (Adjustment Factor)：出力変数のバラツ
キには影響を与えないが期待値には影響を与える制
御因子．
c) 標示因子 (Descriptive Factor)：観測可能で，意図的に実験で取り上げられるが，その水準を制御することは不可能な因子．
c-1) 座標因子（ω 因子）：出力変数ないしはその変数変換の観測区間を表現するだけの形式的標示因子．
d) 誤差因子（あるいはノイズ因子, Noise Factor）：実験
室内ではその水準を制御可能だが，実利用の現場では，観測も制御も不可能な因子と考えられ，その要因効果を可能な限り最小化することが必要な因子．
標示因子を誤差因子とする実験を計画する場合もある．
d-1) ブロック因子 (Block Factor)：制御因子との交互作用を持たない，誤差因子．
e) 補助因子（あるいは共変量, Covariate）：意図的に実験で取り上げられていないが，観測は可能な因子．
制御因子，標示因子，誤差因子の影響を受けて変動する場合もあれば，標示因子や実利用における誤差
因子の変動に影響を与える場合もある．

2.3 因子の分類が導く出力変数最適化戦略
因子の分類が重要なのは，d-1) の定義に関わるブロッ

2.2 田口の実験計画法における因子の分類

ク因子と制御因子との「交互作用 (Interaction)」, すなわ

竹内 [9] は，田口が「実験計画法（上）」2.1.2「因子のち因子の組み合わせ効果をどう扱うべきかに関わるから種類」で導入した因子の分類，すなわち，制御因子，標である．仮に，ブロック因子を想定したにも関わらず，

Oukan Vol.13, No.2

111

Tsubaki, H.

制御因子とブロック因子とに交互作用が存在するなら
ば，モデル (2.1) に αi と β j の交互作用効果 (αβi j ) を加えた (2.2) を用いなければならない．

φ (x)V (µ (xc )) で，Y の制御因子を所与とした際の条件付分散を表現したこととなる．φ (x) が，田口の望目特性の SN 比を一般化した概念となる．ここで，Y が望目特性とは，目標値に近いほど望ましいと考えられることを

Yi j = µ + αi + β j + αβi j + εi j

(2.2)

田口は，制御因子と標示因子及び誤差因子との交互作
用が全て推定可能な直積実験を主として用いてきた．制
御因子と標示因子との直積実験とは，制御因子の水準組
み合わせに関わる大きさ N の実験計画と標示因子の水
準組み合わせに関する大きさ M の実験計画を独立に行い，両実験計画の直積集合（大きさ MN ）の実験を行うことである．直積実験を行うことで，標示因子の水準毎に制御因子を最適化することが可能となる．
田口が標示因子を可能な限り誤差因子として扱うこと意味する．田口は，望目特性に対する 2 段階設計を提唱した．第一段階は目標値に依存せず，バラつきの尺度である SN 比を最適化し，第 2 段階では SN 比に影響を与えない調整因子を用いて，目標値に Y を合わせ込むのである．
推測統計学では Y の変動は，特定の確率分布を想定
し，期待値 µ と分散との関数関係が定まる．例えば，二項分布では，V (µ ) = µ (1 − µ ) である．

Box による，田口のロバスト・パラメータ設計の解釈
は，次の 2 段階設計である．

や制御因子間の交互作用の無視を推奨したロバスト・パ

STEP 1: xopt = arg minx φ (x)ラメータ設計に転換した 1980 年代以降，誤差因子の水⊤⊤
STEP 2: xadjust,opt = {xadjust | µ ((x⊤opt , xadjust ) ) = µtargrt }準によらず出力変数のバラつきを最小化することが目的に代わった．制御因子の最適化によるバラつき最小化後に，調整因子を用いて，期待値を目標値に合わせ込む 2
段階設計を行うのが基本戦略となった．

3.2 確率モデルとの訣別：誤差因子とその調合
田口は，1980 年代以降，誤差に対する確率モデルを想定する既存の統計科学と訣別し，ロバスト・パラメータ設計に移行した．ロバスト・パラメータ設計では，制御因子が所与の状況で，Y の系統変動を引き起こすのは，誤3. 敵対的学習の先駆：SN 比最適化

差因子 n の変動である．従って，出力変数 Y の期待値・3.1 位置・尺度同時要因解析と２段階最適化

分散は，n を定められた閉領域内で実験的に n1 , . . . , nm ,2 章でも触れたように田口が，統計モデル論に与えた
影響は，出力変数のバラツキの尺度としての SN 比の導入と，その最適化である．これは，1980 年代以降，米国統計学の中心にあった G.E.P. Box の影響下，ベル研究所に集っていた Jeff Wu, Vijar Nair ら少壮統計学者による主要研究対象となった．
Fisher の実験計画法が期待値に関する要因分析を行うのに対して，彼らはバラつきに対する要因を同時に分析することに関心を寄せた．このような位置（Location）と尺度（Scale）との同時要因解析の提案がなされたが，特に Box[11] は，田口の SN 比を Nelder ら [12] の一般化
線形モデルの「散らばり母数 (Dispersion Parameter)」に一般化し，散らばり母数こそ最小化すべき出力変数の変動に関する母数とした．Box の 2 段階設計は，次のようなことである．出力変数 Y ，入力制御因子ベクトルを調⊤
整因子でない制御因子と調整因子とに x⊤ = (x⊤c , xadjust )と分解し，次のようなモデルを考える．

と変動させた際の，算術平均ないしは標本分散に過ぎない．しかし，出力変数が二項変数に類似した 0-1 変数な
らば，平均が 0 あるいは 1 のとき，必ず標本分散は 0 となり，分散関数が想定可能というだけである．
統計学者には評判が悪いのだが，田口による誤差因子の調合（水準の選択）は興味深い．調合とは，出力変動を最大化する，すなわち SN 比を最小化する最も不利な誤差因子の水準設定での実験実施である．これは，一般化 SN 比を誤差因子の関数として φ (x, n) と表現した上で，STEP 1 を，次のように行えと提唱したことと解釈できる．

xopt = arg minx arg maxn φ (x, n)
人工知能分野での強化学習は，統計学的には Wald の逐次推論 [13] を実験計画法に適用した逐次実験計画法と密接な関係がある．しかし，田口は制御因子に敵対する誤差因子を提唱し，その最適調合を提案した [14]．この田口の構想は，近年の敵対的学習 [15] の先駆的事例

E[Y ] = µ (xc )と考えることができる．

Var[Y ] = φ (x)V (µ (xc ))ただし，µ (xc ) は，Y の期待値に関する構造を示し，古4. リッジ回帰の先駆：割引係数法

典的には回帰関数と呼ばれる．既知関数 V は，一般化線形モデルでは，出力変数の期待値と分散との関数関係

田口 [6] の 37.1「割引係数」では，ベイズモデルによを表現し，「分散関数 (Variance Function)」と呼ばれる．

る予測と同一の結果を導く「割引係数法」が，提案され112

横幹 第 13 巻 第 2 号

Gen-ichi Taguchi’s Contribution to Statisticsている．

し，指標関数 I( · ) は引数である論理式が真のときは 1,

田口の直交表実験を通常の線形モデルで表現すると(4.1) のようになる．

偽の時は 0 をとる関数である．

Y = X β + ε , E[ε ] = 0, Cov[ε ] = σ 2 I(4.1)

出力変数ベクトルおよび誤差ベクトル Y，ε は，それぞれ n 次元ベクトル，未知母数ベクトル β は，p 次元ベ
クトル，従って計画行列（Design Matrix） X は，n × p
行列，情報行列 X ⊤ X は，実験計画の直交性より，(4.2)
を満たす．

X ⊤ X = diag(δ j )(4.2)

ここで，diag( · ) は，引数を対角要素に持つ対角行列を示す．このとき，β = (β j ) の線形最良不偏推定量（BLUE,

Best Linear Unbiased Estimator） b = (b j ) は，⊤

標関数 I(x ≤ t) を 2 値データととらえ，時点 t の適当な範囲に m 点 t1 ,t2 , . . . ,tm を設定し，確率変数 x を m 変量

0, 1 確率変数
Ci = I(x ≤ t), i = 1, . . . , mに変換することである．Ci を観測変量 x の代わりに推論する方が，情報量が大きくなると考えたのである．田口は，Mi = 1 −Ci を x に対する「精密累積データ」と呼び，出力変数とした．さらに，分布形 Fx (t) を推論の対象とするために，ti を座標因子 (ω 因子) と呼び，要因分
析の対象とした．こうして，分布形 Fx (t) への制御因子
や標示因子等の影響を座標因子との交互作用を検討する
ことで可能とする「精密累積法」が提唱され，マグネッb = diag(1/δ j )X Y, Cov[b] = σ diag(1/δ j )2となる．偏りはあるが，平均二乗誤差が最も小さくなる
線形最良不変推定量は，

diag(∆ j /(1 + ∆ j ))b

田口が意図した分布関数に基づく一般的方法とは，指トワイヤーの寿命試験に適用された．

5.2 精密累積法の統計モデルによる解釈(4.3)

椿は，田口の追悼論文 [19] の中で，精密累積法のアとなる．ただし，∆ j = δ j /σ 2 であり，非心度と呼ばれイデアが， Lindsay[20] の強度乗法モデルに継承されてる．しかし，非心度 ∆ j は未知なので，b j と σ 2 の不偏推定量 s2 = Y⊤ (I − X(X ⊤ X)−1 X ⊤ )Y/(n − p) を用いることが考えられる．田口は，その推定量を max(Fj − 1, 0),

Fj = δ j b2j /s2 とし，この推定方式を割引係数法と呼んだ．
(4.3) は，Ridge Regression の名で提唱された推定量 [16]であり，現在機械学習分野でも用いられている．この推定量は，経験ベイズ推定量 [17] と見なすこともできる．
日本品質管理学会では 1980 年代に縮小推定量と呼ばれたこの種の方法の応用可能性を集中的に検討した [18]．
一方，1960 年代に遡ると田口 [6] は，最尤法と条件付分布に基づく，より数理的割引係数を導入している．そこでは，(4.1) が「不偏推定であることは，データが出る
前にのみ正しい」，しかしひとたびデータを見た後，つまり非心度が観測されたという条件の下では，それを条件付けた β の最尤推定量を求めるのが良いというアイデアを基に，小標本での割引係数を求めている．

いることを指摘した．以下，その要旨である．
正値確率変数 x に対応する精密累積データ Ck が，第k 単位時刻 (k = 1, . . . , M) に観測されたとし，それらが平均 λk のポアソン分布に独立に従うと仮定する．ここで，座標因子で指定された区間幅は十分小さく，単位時
刻に 2 回以上の事象が観測される確率が無視可能なほど小さいと仮定する．
このとき，第 m 単位時刻 (m ≤ M) に事象が生起し，第 m 単位時刻以降，精密累積データが 1 となったとする．すなわち，C1 = . . . = Cm−1 = 0, Cm = . . . = CM = 1の場合を考える．この仮定の下で，精密累積データのポアソン尤度関数は (5.1) となる．
exp(−

∑

λk ) (∏λk ) exp(−

k=m,...,M

k=1,...,m−1∑λk ) (5.1)

k=m,...,M

一方，原変数 x を生成する正値分布のハザード関数（瞬間故障率関数）を λ (t) とすれば，第 m 単位時刻に対応する時刻 tm に事象が発生した場合の尤度関数への寄

5. 確率分布モデルに対する推論

与は，(5.2) となる．

5.1 精密累積データの提案

λ (tm ) exp(−Λ(tm ))田口 [6] 32.2 では，「どんな特性値もすべて “ある，なし” という 2 組の分類によって表現される」とした．更

ここで，Λ(t) =に，任意の確率変数 x は，それが従う分布関数 Fx (t) の推(5.2) から，定値として用いられることを指摘した．田口 [6] には，

∫ t
0

(5.2)

λ (s)ds は，累積ハザード関数である．

exp(−

∑

λi ) λm(5.3)

数理的記述はないが，指標関数 I(x ≤ t) が Fx (t) の不偏

k=1,...,m

推定量であることを利用すれば，どのような確率変数もが，背後の正値確率変数 x の分布の厳密な尤度関数を離散統一的に扱えることに注目したものと考えられる．ただ近似したものとなることがわかる．これが Lindsay[20]

Oukan Vol.13, No.2

113

Tsubaki, H.

Table 1: Modified Accurate Minute Data C in case of X = 5

は事象の生起を促進する要因を抽出しようとするのに対し，逆ハザード関数のモデリングは事象の生起を抑制す
る要因を抽出するものと考えられる [22]．この種のモデリング戦略も十分な議論が必要である．

6. 因果推論の萌芽としての確率対応法のアイデアである．以上から，田口の精密累積データ
は，Table １のように修正しなければならないことがわ
かる．k = 1, . . . , m までの修正精密累積データが観測されたものとすれば，分布形に依存せず近似尤度関数に基づく漸近有効な推論がポアソン回帰モデルを通じて可能となる．
しかし，Lindsay は，田口が考えた座標（時点）因子と制御因子などとの交互作用を検討可能なモデリングを提案していない．また，ポアソン変量の座標因子水準毎の強度関数λをノンパラメトリックに分析することも扱っていない．椿 [19] が示したように一般化加法モデル (Generalized Additive Model)[21] を用いることで，田口の意図に不完全ながらも近づくことは可能である．田口は，制御因子と座標因子との交互作用の検討を分割法（Split Plot Design）を用いて行っている．これは，確率的誤差を制御因子の水準設定に付随して起きる環境誤差（一次誤差）と，座標因子の水準において生じる偶然変動（二次誤差）との２つを考えるものである．今後，時間依存性共変量の導入時に本質的役割を果たすアイデア

6.1 無作為化計画としての確率対応法
Propensity Score[23] に基づく統計的因果推論は，統計的方法における無作為割付や無作為抽出の意味を明確にした．既に田口 [5] は，1950 年代に Fisher の無作為割付が要因効果の不偏推定を構成する操作であることを明確に指摘している．これを積極的に利用した技術が「確率対応法」であり，統計モデル論と Neyman[24] の無作為抽出を結合した方法論である．この方法を，田口は 1955
年に電気通信研究所で行われた発泡 PVC の実験で用いた．田口の確率対応法は，通常の実験回数では実現できない多因子多水準実験に対して，因子を 2 分し,第一群の因子に対する実験は行うという条件の下で，第一群の因子に直交する第 2 群の因子に関する複数の実験計画
を確率的（無作為）に選択するという方法となっている．
今日に至るまで色あせていない独創的アイデアである．
椿 [25] を要約し，田口のアイデアを単純化しよう．先ず線形モデルから発生する実験データに対して大きさ Nの有限母集団 (Y1 , x1 ), . . . ., (YN , xN ) を考える．ただし，だが，階層的誤差モデルは Lindsay[20] も椿 [19] も想定

Yi = x⊤i β + εi , i = 1, . . . , N

していない．

(6.1)

田口は，産業界で自らの統計的方法を利活用した品の 成 立 を 仮 定 す る ．こ こ で ，xi は ，計 画 行 列 X =質工学研究グループに対して「精密累積法を考案した(x1 , x2 , . . . , xN ) の列ベクトルとし，母集団全てで既知，
β は未知パラメータベクトルで，いずれも p 次元ベクトルとする．また εi は期待値 0, 分散 σ 2 の正規確率変数とする．
更に，出力変数 Yi は，計画母集団から確率 pi で抽出され，実験されることで観測されるとする．また，確率的に抽出されない場合は，実験データが欠測したものと考える．Fisher の実験計画法の 3 原則の一つである無作為化の原則は，抽出された実実験に対して系統誤差を偶然誤差化するための無作為化である．しかし，ここでの無作為化は極めて大きな実験計画の空間からの部分実施
実験を無作為に選択することを意味する．
以上の準備の下で確率対応法を線形モデルの推論一般に拡張し，今日的なアルゴリズム風に記載すれば，以下のようになる．

ことで，統計学に対してやるべきことは全て行い，興味を失った」と語ったとされている．実際，精密累積法は生存時間分析に限らず，正値確率変数全般に適用可能である．また，一般の確率変数を正値確率変数に変換する
ことは容易なので，田口が指摘しているようにあらゆるデータの解析に利用可能である．
しかし，精密累積データを 0,1 データの解析法と考え
ると，田口の精密累積法は 0 と 1 との置換に対しては不変な推論となっていた．一方，ポアッソンモデル当ては
めは，事象生起を 1 としなければならないので，0 と 1との置換は不変ではない．田口が SN 比解析の望大特性解析，すなわち，特性が大きくなればなるほど良いと考えられる場合に用いた，データの逆数変換を正値確率変数に適用すれば，無限大が 0, 0 が無限大に変換される．
通常の時間軸座標で，事象が生起し寿命が観測された以降の事象を 0，寿命が観測された時点を 1 として，精密累積データを構成し，Lindsay の近似尤度を構成することも可能である．しかし，それは逆ハザードモデルと呼ばれる異なった結果を導く．ハザード関数のモデリング

114

STEP 0: (未知パラメータ推定値の初期化)
k = 0 とし, β の初期推定量 bk を適当に与える．
STEP 1: 予測値 Qki の算出
Yi の予測値を Qki = Xi⊤ bk とする．

横幹 第 13 巻 第 2 号

Gen-ichi Taguchi’s Contribution to Statistics

STEP 2: Yi の不偏推定量 yki の算出
Yi が抽出されない場合：yki = Qki
Yi が抽出された場合：yki = Qki + (Yi − Qki )/piy⊤k = (yk1 , . . . , ykN )るまでの分散の増大である．確率対応法の魅力と問題点
は，Propensity Score 法など因果推論にも全て引き継がれている．このため，統計的方法の中でも最も注目されているベイズ統計と因果推論とに齟齬が生じている．筆者
は，ベイズ統計学のように条件付き原理を認めれば，因STEP 3: β の推定量の OLS での更新
bk+1 = (X ⊤ X)−1 X ⊤ yk

果推論は本来許容されないと考えている [26]．しかし，多くの統計モデルの前提条件がデータの不偏性である以STEP 4: 収束判定
∥bk+1 − bk ∥ < δ → bk+1 を β の推定量とする
Otherwise → Step 1 に戻る
STEP 5: σ 2 の推定量 s2 の算出
e = yk − Xbk+1
s2 = e⊤ e/(N − p), Cov[bk+1 ] の推定量 = s2 (X ⊤ X)−1

上，データを編成して，ユーザーに提供する立場の統計家には，不偏性も必要な制約となる場合も多い．田口が1950 年代に今日的統計的方法を発想し，産業界で実践した先駆性は評価されるべきである．

7. おわりになお，田口の確率対応法はより巧妙で，計画空間から

本稿で概説した統計学の原理に関わる方法以外にも，の単純確率抽出ではない．推論効率の高い直交計画が一

田口は，第 2 世代の人工知能の CART(Classification andRegression Tree)[27] に相当するアイデアを「実験計画法（上）」9 章で，「分類によるデータの解析法（分類法，カード解析法）」，「新版実験計画法（下）」24 章「逐次分類法」として構想している．そこでの記述はこの分野の嚆矢とされている AID(Automatic Interaction Detectoe)[28]のカード分類の自動化と類似している．
現代統計学の標準的技法となった欠測値処理に関する
EM アルゴリズム [29] については「実験計画法（下）」22章，「新版実験計画法（下）」の「欠測値の処理」として提示している．
増山 [30] は，「田口の『実験計画法』は，知られていな
い様々なアイデアが含まれており」とし，若手研究者に積極的な研究を勧めた．増山は 1989 年，筆者に「田口が文献に書いているのは，氷山が水面に浮かんでいる部分のようなもの」と述べた．実際，本稿で紹介した田口の統計的方法は，学術論文として，数理的根拠を明示して提案されたものではない．田口が本来は実現したかった
ことは，「新版実験計画法」序文に記載された通り，尤度（エントロピー）に基づく統計学的厳密解だった．これは完成せず，公表されなかった．田口は，単純な分散分
析技術をミニマックス解に近い方法として実用に供し，膨大な事例を叙事詩の如く創生した．このため，内外の統計学者から様々な批判があった．しかし，田口の「実験計画法」は．新版以降，下巻を中心にほぼ全ての章に何らかの独創性が示されている．また，上巻に工学全般に関わる考え方が示され，1980 年代以降，ロバスト・パ
ラメータ設計を導いたアイデアも示されている．統計科
学に関心を失った後も，本特集で大久保豪人氏が示した
マハラノビス・タグチシステムなどを提唱している．
筆者は，田口の統計学に関わるアイデアは，厳密に考えれば，統計学的にも正当化できる，あるいは，すでに正当化された深い意味を持つと考える．

部の重要な因子では必ず実現する前提の上で，他の因子に関する直交計画をクラスター抽出となっている．一次元データ Yi が独立に確率抽出されるのではなく，計画で規定される一群の観測値ベクトル Y = (Yi ) が，抽出されるのである．また，初期推定量の作り方についても観測されたデータに基づき一部の因子を用いた最小二乗推定量を用いるなどの工夫がなされている．確率対応法が提示されてから 60 年以上たった今日でも，実験計画の無作為化によって，観測データだけでは推定不可能なパラメータが推定可能となることは，一般の工学者からすれば，Neyman の標本調査論 [24]，今日隆盛を極める欠測値理論 [23] が導く「手品」以外のなにものでもなかろう．

6.2 不偏データで得ること失うこと無から有が作れる方法が可能となるのは，無作為化に基づく Yi の不偏推定量 yi を観測データの代わりに擬似的に使うからである．一方，この擬似データ yi には，(6.1) 式で示した実験計画法の統計モデルの誤差変動 σ 2に加えて，標本抽出に伴う変動を考慮することが必要である．誤差の自由度を大幅に増やし，事実上大標本理論が適用できることの代替に，誤差分散が大きくなる．
本来推定可能ではない要因効果が推定できるようになったのは，無作為化だけが役割を果たした訳ではない．ベイズ統計の事前分布想定に相当するのが確率対応
法では，STEP 0 のパラメータ初期値の推定である．しかし，通常のベイズ推定を，頻度論的パラメータ推定と見なすと不偏推定は導かず，割引係数法のような偏りのある縮小推定量を導く．不偏性は推論に対する制約条件
であり，無作為化に基づき全ての未知パラメータを不偏推定可能とする制約に対しては何らかのペナルティが存在する．それが，通常の推定論では非許容的とみなされ

Oukan Vol.13, No.2

115

Tsubaki, H.

参考文献[1] Fisher, R. A., The Design of Experiments, Macmillan(1935)[2] Nair, V. N., Abraham, B., MacKay, J., Box, G., Kacker, R.
N., Lorenzen, T. J., Lucas, J. M., Myers, R. H., Vining, G.
G., Nelder, J. A., Phadke, M.
S., Sacks, J., Welch, W. J.,
Shoemaker, A. C., Tsui, K. L. , Taguchi, S. and Wu, C. F. J.,
Taguchi’s Parameter Design: A Panel Discussion, Technometrics,Vol. 34(2), pp. 127-161 (1992)[3] ISO, Applications of statistical and related methods to newtechnology and product development process – Robust parameter design (RPD)，ISO 16336:2014, International Organization for Standardization (2014)[4] 椿 広計，アメリカのタグチメソッド関連シンポジウム，特集「新製品開発に関するパラメータ設計」，品質管理，
Vol.43(12), pp.2063-2068 (1992)[5] 田口玄一，実験計画法 (上)，(下)，丸善 (1957, 1958)[6] 田口玄一，新版実験計画法 (上)，(下)，丸善 (1962)[7] 田口玄一，第 3 版実験計画法 (上)，(下)，丸善 (1976, 1977)[8] 椿 広計，河村敏彦，設計科学におけるタグチメソッド‐
パラメータ設計の体系化と新たな SN 比解析，日科技連
出版 (2008)[9] 竹内 啓, 要因実験における因子の種類と推測の方法, 京都
大学数理解析研究所講究録 526「実験データ解析の理論的背景」, 1/12 (1984)[10] 竹内 啓, 田口の方法をめぐって，品質, 16-2, 54/67 (1986)[11] Box, G. E. P.，Signal-to-Noise Rations, Performance Criteria, and Transformations, Technometrics, Vol. 30(1), pp.
1-17 (1988)[12] Nelder, J. A. and Wedderburn, R.W.M., Generalized Linear
Models, Journal of the Royal Statistical Society, Series A.,
Vol. 135(3), pp. 370-384 (1972)[13] Wald, A., Sequential Analysis, John Wiley (1947)[14] 田口玄一，三次元機構設計の誤差解析，精密工学会誌，Vol.52(8), pp. 1285-1288 (1986)[15] Gooffwllow, I. J., Pouget-Abadie, J., Mirza, M., Xu, B.,
Warde-Farley, D., Ozair, S., Courvillem A, Bengio, Y., Generative Adversarial Network, arXiv: 1406, 2661 (2014)[16] Hoerl, A. E., Kennard, R. W., Ridge regression: Biased estimation for nonorthogonal problems, Technometrics, Vol,12(1), pp. 55-67 (1970)[17] Efron, B. and Morris, C., Stein’s Estimation Rule and , Its
Competitors–An Empirical Bayes Approach，Journal of theAmerican Statistical Association，68, 117/130 (1973)[18] 鷲尾泰俊，篠崎信雄．椿 広計，張 元宗，縮小推定量とその応用‐統計手法研究会中間報告，品質，Vol.12(2), pp.
127-139 (1982)

116[19] 椿 広計, 田口の精密累積法のセミパラメトリックポアソン
モデルによる再定式化，応用統計学，Vol.42(3), pp.145-159(2013)[20] Lindsay, J. K., Parametric multiplicative intensities models
fitted to bus motor failure data, Applied Statistics, Vol.46(2),pp. 245-252 (1997)[21] Hastie, T. and Tibshirani, R., Generalized Additive Models;Some Application, Journal of the American Statistical Association, Vol.82, pp. 371-386 (1987)[22] 椿 広計，大野忠士，定量的リスク評価と定性的リスク評価との架橋：定量的リスク評価モデル当てはめにおける質的選択モデルの役割，計量生物学，Vol.29，pp. s133-s141(2008)[23] Rosenbaum, P. R. and Rubin, D. B., The central role of thepropensity score in observational studies for causal effects,
Biometrika, 70-1, 41/55 (1983)[24] Neyman, J., On the two different aspects of the representative method: the method of stratified sampling and themethod of purposive selection, Journal of the Royal Statistical Society, 97-4, 558/625 with discussion (1934)[25] 椿 広計，田口玄一の実験計画法―統計学に対する貢献―，第 9 回横幹コンファレンス予稿集，https://doi.
org/10.11487/oukan.2018.0_B-1-4[26] 椿 広計，CTI 作成における偏り補正の数理，特集 消費者動向指数 (CTI) の開発と応用，ESTRELA, No.292, pp.11-18(2018)[27] Breiman, L., Classification and Regression Trees, Taylor &Francis (1984)[28] Morgan, J. S. and Sonquist, J.A., Problems on the Analysis
of Survey Data, and a Proposal, Journal of the American
Statistical Association, Vol. 58, pp. 415-434 (1963)[29] Dempster, A. P., Laird, N.M. and Rubin, D. B., Maximum
Likelihood from Incomplete Data Via the EM Algorithm,
Jouranal of the Royal Statistical Society, Ser. B., Vol.39(1),pp. 1-38 (1977)[30] 増山元三郎，実験計画法第 2 版，岩波書店 (1972)

椿 広計

横幹 第 13 巻 第 2 号

1956 年 10 月 20 日生．82 年東京大学大学院工学研
究科修士課程計数工学専攻修了．東京大学助手，慶應
義塾大学講師，筑波大学助教授・教授，統計数理研究
所教授，（独）統計センター理事長を経て，2019 年情報・システム研究機構理事・統計数理研究所長，現在に至る．統計家として実務・研究に従事．工学博士．
日本品質管理学会，応用統計学会などの会員