2021/09/02 (更新日: 2023/12/08)

実験計画法で欠測値を推定する方法がわかる【重要】

実験計画法

こういう疑問に答えます。

本記事のテーマ

欠測値を推定する方法がわかる【重要】

• ①欠測値の処理方法がわかる
• ②二元配置実験(繰返し無し)の場合の欠測値の推定方法がわかる
• ③二元配置実験(繰返し有り)の場合の欠測値の推定方法がわかる

記事の信頼性

記事を書いている私は、QC検定®１級合格した後、さらに実験計画法に磨きをかけています。とはいえ、QC検定®1級合格前の1.5年前までは、実験計画法すら知りませんでした。実験計画法を初めて勉強して3ヶ月後にQC検定®2級を合格しました。実験計画法はまったく理解できていませんでしたが、計算方法だけ暗記して点数を稼ぐレベルでした。

本記事は、実験計画法を学び始めるときに、なぜ？と不思議に思う内容をわかりやすく解説します。すぐ読めます！

①欠測値の処理方法がわかる

データを捨てるか推定して入れるか

データを用意します。因子A(3水準)と因子B(4水準)のデータがあるが、AB22が欠測していまい、11個のデータしかないとします。このままでは分散分析ができません。さて、どうしますか？

	A1	A2	A3
B1	30	40	50
B2	40	×	60
B3	50	70	80
B4	100	120	130

欠測値があると、繰返し数が異なります。一元配置実験ならば、繰返し数の異なる場合として分散分析できますが、二元以上の配置実験になると分散分析ができません。なぜなら平方和が分解できないからです。

繰返し数が異なる二元配置実験は分散分析ができない理由については、関連記事繰返し数が異なる場合は一元配置実験だけである理由がわかるに詳細な計算過程を解説しています。

欠測値の処理方法は、次の2つです。

1. データを捨てる
2. 欠測値を推定して追加する

データを捨てる(推奨しない)

下表のように、データを削って分散分析する方法です。

	A1	A3
B1	30	50
B3	50	80
B4	100	130

でも、推奨しません。。
理由は
・データを削るのはもったいない
・データの特性が変化する
からです。データ分析結果が変わってしまうため、推奨しません。

欠測値を推定して追加する(推奨)

欠測といっても、1個か2個程度なので、下表のように、推定値を作って分散分析します。
これなら、測定した全データを活かせるので、分析結果の精度を低下することはありません。

	A1	A2	A3
B1	30	40	50
B2	40	推定値	60
B3	50	70	80
B4	100	120	130

でも、ここで１つ疑問が湧きます。

データの推定方法

データの推定方法

唯一無二な方法ではなく、都度推定値を考えて算出すればよいですが、
ある変数が最小値になる条件(微分=0)を使います。
誤差が大きくならないように気をつけて推定したいことから、残差平方和Seが最小となる推定値を決めます。

以下、二元配置実験での欠測値の推定方法を解説します。一元配置実験は、関連記事一元配置実験の平方和の分解ができる【初心者必見】に詳細な計算過程を解説しています。

②二元配置実験(繰返し無し)の場合の欠測値の推定方法がわかる

欠測値の数によらず、解析方法は同じですが、途中経過の式は欠測値の数によって異なるので、個別に解説します。

㉑欠測値1つの場合

以下のような2因子実験を考えます。AB11だけ欠測だったとします。

	A1	A2	・・・	Aa	計
B1	–	x21	・・・	xa1	\(T’_{・1}\)
B2	x12	x22	・・・	xa2	\(T_{・2}\)
・・・	・・・	・・・	・・・	・・・	・・・
Bb	x1b	x2b	・・・	xab	\(T_{・b}\)
計	\(T’_{1・}\)	\(T_{2・}\)	・・・	\(T_{a・}\)	T’

AB11に推定値を入れます。未知数なのでxとして代入します。

	A1	A2	・・・	Aa	計
B1	x	x21	・・・	xa1	x+\(T’_{・1}\)
B2	x12	x22	・・・	xa2	\(T_{・2}\)
・・・	・・・	・・・	・・・	・・・	・・・
Bb	x1b	x2b	・・・	xab	\(T_{・b}\)
計	x+\(T’_{1・}\)	\(T_{2・}\)	・・・	\(T_{a・}\)	x+T’

xだけ増加したことがわかります。

次に各平方和を求めます。
CT=\(\frac{1}{ab}(x+T’)^2\)
S_T=\(x^2+x_{21}^2+…+x_{ab}^2\)
S_A=\(\frac{1}{a}((x+T’_{1･})^2+T_{2･}^2+…+T_{a･}^2)\)-CT
S_B=\(\frac{1}{b}((x+T’_{･1})^2+T_{･2}^2+…+T_{･b}^2)\)-CT
S_e= S_T– S_A– S_B
=\(x^2+x_{21}^2+…+x_{ab}^2\)
-\(\frac{1}{a}((x+T’_{1･})^2+T_{2･}^2+…+T_{a･}^2)\)
-\(\frac{1}{b}((x+T’_{･1})^2+T_{･2}^2+…+T_{･b}^2)\)
+\(\frac{1}{ab}(x+T’)^2\)
と式で書けます。

\(\frac{dSe}{dx}\)
=2x-\(\frac{2}{a}( x+T’_{1･})\)-\(\frac{2}{b}( x+T’_{･1})\)+\(\frac{2}{ab}(x+T’)\)
=0
整理すると、
x=\(\frac{b T’_{1･}+a b T’_{･1}-T’}{(a-1)(b-1)}\)
となります。1をi,jに変えて一般化します。
x=\(\frac{b T’_{i･}+a T’_{･j}-T’}{(a-1)(b-1)}\)

下表の場合、AB11を推定しましょう。

	A1	A2	A3	A4	計
B1	×	20	50	60	130
B2	50	70	90	110	320
B3	70	90	110	130	400
B4	150	170	180	200	700
B5	210	220	240	260	930
計	480	570	670	760	2480

x=\(\frac{b T’_{i･}+a b T’_{･j}-T’}{(a-1)(b-1)}\)
に代入します。上表の黄色枠から
a=4,b=5, \(T’_{i･}\)=130,\( T’_{･j}\)=480,T=2480を代入します。
x=(5×130+4×480-2480)/((4-1)(5-1))=7.5

グラフで、推定値xと残差平方和Seの関係をグラフにすると、
確かにx=7.5付近で残差平方和Seが最小になっているのがわかります。

㉒欠測値2つの場合

実は、欠測値1つの場合が理解できれば、簡単です。
変数を2つ(x,y)にしたらよいからです。

欠測値が1つ場合の残差平方和Seは、
S_e= S_T– S_A– S_B
=\(x^2+x_{21}^2+…+x_{ab}^2\)
-\(\frac{1}{a}((x+T’_{1･})^2+T_{2･}^2+…+T_{a･}^2)\)
-\(\frac{1}{b}((x+T’_{･1})^2+T_{･2}^2+…+T_{･b}^2)\)
+\(\frac{1}{ab}(x+T’)^2\)
ですね。

欠測値が2つ場合の残差平方和Seは、
S_e= S_T– S_A– S_B
=\(x^2+ y^2+x_{22}^2+…+x_{ab}^2\)
-\(\frac{1}{a}((x+T’_{i･})^2+(y+T’_{i’･})^2+T_{2･}^2+…+T_{a･}^2)\)
-\(\frac{1}{b}((x+T’_{･j})^2+(y+T’_{･j’})^2+T_{･2}^2+…+T_{･b}^2)\)
+\(\frac{1}{ab}(x+y+T’)^2\)
と、yを増やせば式が作れます。

\(\frac{\partial Se}{\partial x}\)
=2x-\(\frac{2}{a}( x+T’_{i･})\)-\(\frac{2}{b}( x+T’_{･j})\)+\(\frac{2}{ab}(x+y+T’)\)
=0
\(\frac{\partial Se}{\partial y}\)
=2y-\(\frac{2}{a}( y+T’_{i’･})\)-\(\frac{2}{b}( y+T’_{･j’})\)+\(\frac{2}{ab}(x+y+T’)\)
=0

整理すると,
(a-1)(b-1)x+y=b\(T’_{i･}\)+a\(T’_{･j}\)-T’
(a-1)(b-1)y+x=b\(T’_{i’･}\)+a\(T’_{･j’}\)-T’
を満たす、x,yとすればよいです。

㉓欠測値3つ以上の場合

2つの場合と同様にx,y,z,…と未知数を用意して残差平方和が最小となるように式を立てます。

S_e= S_T– S_A– S_B
=\(x^2+ y^2+z^2+x_{21}^2+…+x_{ab}^2\)
-\(\frac{1}{a}((x+T’_{i･})^2+(y+T’_{i’･})^2+(z+T’_{i”･})^2+T_{2･}^2+…+T_{a･}^2)\)
-\(\frac{1}{b}((x+T’_{･j})^2+(y+T’_{･j’})^2+(z+T’_{･j”})^2+T_{･2}^2+…+T_{･b}^2)\)
+\(\frac{1}{ab}(x+y+z+T’)^2\)
と、zを増やせば式が作れます。

\(\frac{\partial Se}{\partial x}\)
=2x-\(\frac{2}{a}( x+T’_{i･})\)-\(\frac{2}{b}( x+T’_{･j})\)+\(\frac{2}{ab}(x+y+z+T’)\)
=0
\(\frac{\partial Se}{\partial y}\)
=2y-\(\frac{2}{a}( y+T’_{i’･})\)-\(\frac{2}{b}( y+T’_{･j’})\)+\(\frac{2}{ab}( x+y+z+T’)\)
=0
\(\frac{\partial Se}{\partial z}\)
=2z-\(\frac{2}{a}( z+T’_{i”･})\)-\(\frac{2}{b}( z+T’_{･j”})\)+\(\frac{2}{ab}( x+y+z+T’)\)
=0

整理すると,
(a-1)(b-1)x+y+z=b\(T’_{i･}\)+a\(T’_{･j}\)-T’
(a-1)(b-1)y+x+z=b\(T’_{i’･}\)+a\(T’_{･j’}\)-T’
(a-1)(b-1)z+x+y=b\(T’_{i”･}\)+a\(T’_{･j”}\)-T’
を満たす、x,y,zとすればよいです。

欠測値の数が4以上となれば、同様の式を立てて連立方程式を作って解けばOKです。欠測値の求め方を理解するだけでOKです。

③二元配置実験(繰返し有り)の場合の欠測値の推定方法がわかる

欠測値の数によらず、解析方法は同じですが、途中経過の式は欠測値の数によって異なるので、個別に解説します。

㉛欠測値1つの場合

以下のような2因子実験を考えます。AB111だけ欠測だったとします。

	A1	A2	…	Aa	計
B1	–	x211	…	xa11	\(T’_{・1・}\)
	x112	…	…	…
	x11c	x21c	…	xa1c
B2	x121	x221	…	xa21	\(T’_{・2・}\)
	…	…	…	…
	x12c	x22c	…	xa2c
…	…	…	…	…	…
Bb	x1b1	x2b1	…	xab1	\(T’_{・b・}\)
	…	…	…	…
	x1bc	x2bc	…	xabc
計	\(T’_{1・・}\)	\(T’_{2・・}\)	…	\(T’_{a・・}\)	T’

AB111に推定値を入れます。未知数なのでxとして代入します。

	A1	A2	…	Aa	計
B1	x	x211	…	xa11	x+\(T’_{・1・}\)
	x112	…	…	…
	x11c	x21c	…	xa1c
B2	x121	x221	…	xa21	\(T’_{・2・}\)
	…	…	…	…
	x12c	x22c	…	xa2c
…	…	…	…	…	…
Bb	x1b1	x2b1	…	xab1	\(T’_{・b・}\)
	…	…	…	…
	x1bc	x2bc	…	xabc
計	x+\(T’_{1・・}\)	\(T’_{2・・}\)	…	\(T’_{a・・}\)	x+T’

xだけ増加したことがわかります。

次に各平方和を求めます。

S_A×B= S_AB– S_A– S_B
は,
S_AB= S_A+S_B+ S_A×B
より、
S_e= S_T– S_A– S_B– S_A×B
= S_T-( S_A+S_B+ S_A×B)
= S_T– S_AB
とシンプルな式になります。

CT=\(\frac{1}{abc}(x+x_{112}+…+x_{abc})^2\)
S_T=(\(x^2+x_{112}^2+…+x_{abc}^2\))-CT
S_AB=\(\frac{1}{c}((x+…+x_{11c})^2\)+\((x_{211}+…+x_{21c})^2\)+…+\((x_{ab1}+…+x_{abc})^2)\)-CT
より、
S_e
=(\(x^2+x_{112}^2+…+x_{abc}^2\))-\(\frac{1}{c}((x+…+x_{11c})^2\)+\((x_{211}+…+x_{21c})^2\)+…+\((x_{ab1}+…+x_{abc})^2)\)

\(\frac{dSe}{dx}\)
=2x-\(\frac{2}{c}(x+…+x_{11c})\)=0
x=\(x_{112}+…+x_{11c}\)
x=\(\frac{1}{c-1}(x_{112}+…+x_{11c})\)

x=x111でしたが、これを一般化すると、
x=\(\frac{1}{c-1}(x_{111}+…+x_{11c})\)
(ただし、右辺は欠測値x_abcを除く)
となります。

よく見ると、欠測値以外の合計値をその個数で割った平均値になっていますね。

㉜欠測値2つの場合

実は、欠測値1つの場合が理解できれば、簡単です。
変数を2つ(x,y)にしたらよいからです。

欠測値が2つ場合の残差平方和Seは、
S_e
=(\(x^2+y^2+…+x_{abc}^2\))-\(\frac{1}{c}((x+y+…+x_{11c})^2\)+\((x_{211}+…+x_{21c})^2\)+…+\((x_{ab1}+…+x_{abc})^2)\)
と、yを増やせば式が作れます。

\(\frac{\partial Se}{\partial x}\)
=2x-\(\frac{2}{c}(x+y+…+x_{11c})\)=0
\(\frac{\partial Se}{\partial y}\)
=2y-\(\frac{2}{c}(x+y+…+x_{11c})\)=0
を整理すると、
(c-1)x=y+…+\(x_{11c}\)
(c-1)y=x+…+\(x_{11c}\)
一般化すると、
(c-1)x=y+…+\(x_{ijc}\)
(c-1)y=x+…+\(x_{ijc}\)
の連立方程式からx,yを求めることができます。

㉝欠測値3つ以上の場合

2つの場合と同様にx,y,z,…と未知数を用意して残差平方和が最小となるように式を立てます。

S_e
=\((x^2+y^2+z^2+…+x_{abc}^2)\)-\(\frac{1}{c}(x+y+z+…+x_{11c}^2)\)+\((x_{211}+…+x_{21c})^2\)+…+\((x_{ab1}+…+x_{abc}^2)\)
と、zを増やせば式が作れます。

\(\frac{\partial Se}{\partial x}\)
=2x-\(\frac{2}{c} (x+y+z+…+x_{11c})\)=0
\(\frac{\partial Se}{\partial y}\)
=2y-\(\frac{2}{c} (x+y+z+…+x_{11c})\)=0
\(\frac{\partial Se}{\partial z}\)
=2y-\(\frac{2}{c} (x+y+z+…+x_{11c})\)=0

を整理すると、
(c-1)x=y+z+…+\(x_{11c}\)
(c-1)y=x+z+…+\(x_{11c}\)
(c-1)z=y+z+…+\(x_{11c}\)
一般化すると、
(c-1)x=y+z+…+\(x_{ijc}\)
(c-1)y=x+z+…+\(x_{ijc}\)
(c-1)z=y+z+…+\(x_{ijc}\)
の連立方程式からx,y,zを求めることができます。

欠測値の数が4以上となれば、同様の式を立てて連立方程式を作って解けばOKです。欠測値の求め方を理解するだけでOKです。

まとめ

欠測値がある場合の推定方法を解説しました。

• ①欠測値の処理方法がわかる
• ②二元配置実験(繰返し無し)の場合の欠測値の推定方法がわかる
• ③二元配置実験(繰返し有り)の場合の欠測値の推定方法がわかる

実験計画法

Warning: count(): Parameter must be an array or an object that implements Countable in /home/qcplanets/qcplanets.com/public_html/wp-content/themes/m_theme/sns.php on line 119

Warning: Invalid argument supplied for foreach() in /home/qcplanets/qcplanets.com/public_html/wp-content/themes/m_theme/sns.php on line 122