公開日2021年10月16日 最終更新日 2022年12月1日
みなさんこんにちは、michiです。
前回の記事「【実践 5】OC曲線の考え方・作り方」に引き続き、今回はQC検定実践編として、Excelでワイブル分布の作り方について学んでいきます。
ワイブル分布の詳細は、記事「ワイブル分布 とは」をご参照ください。
この記事を読めば、QC検定1級の論述対策になる!…カモ
キーワード:「ワイブル分布確率紙」「検査の意味」
目次
①データ表の作成
ワイブル分布の解析において、まずは取得したデータを順番に並べます。
例えば、 n=15 のサンプルで寿命試験をしたとします。(図1)
サンプルはすべて同じ時間に壊れることはなく、壊れる時間に差があります。
寿命のばらつきです。
計算をするため、壊れた順番に並び変えます。(図2)
\[\]
②累積故障頻度\(F(t(i))\) の計算
寿命試験のデータを整理できたら、累積故障頻度\(F(t(i))\)を計算します。
計算方法には、平均ランク法とメジアンランク法があります。
- 平均ランク法:\(F(t(i))=\frac{i}{n+1}×100(\%)\)
- メジアンランク法:\(F(t(i))=\frac{i-0.3}{n+0.4}×100(\%)\)
※\(i\) は「順番」になります。
データ数が30個以下の場合はメジアンランク法を、30個よりも多い場合は平均ランク法を使うと良いでしょう。
\[\]
今回はデータ数が30個以下なので、メジアンランク法で計算します。(図3)
計算式は、「=(“順番”-0.3)/(COUNT(“順番の範囲”)+0.4)」です。
計算式では、n数(サンプル数)を「COUNT関数」を使って求めています。
サンプル数によって数字が変わります。
この計算を累積故障頻度\(F(t(i))\) すべてに対して実行します。(図4)
図4を見ると、15番目(最後)の寿命試験が終わっても累積故障頻度\(F(t(i))\)は100%にはなっていません!
これは、ワイブル分布の計算では「試験データよりも寿命が長い製品もあるだろう」と仮定しているためです。
\[\]
③ワイブル分布縦軸 \(Y\) の計算
累積故障頻度\(F(t(i))\) を計算できたので、ワイブル分布確率紙の縦軸 \(Y\) を計算します。
\[Y=lnln\frac{1}{1-F(t(i))}\]
を計算します。(図5)
計算式は「=LN(LN(1/(1-“対象のF(t(i))”)))」です。
\[\]
ほかのセルにも同様に計算します。(図6)
累積故障頻度\(F(t(i))\) をパーセント表示に変更しておきました。
G列に次の計算用にセルを用意します。
\[\]
④散布図を描く
「③ワイブル分布縦軸 \(Y\) の計算」ができたら、散布図を描きます。
散布図は横軸を「時間」に、縦軸を「F(t)」で描きます。
\[\]
プロットした点が直線上にあれば大丈夫です。
「直線が途中で折れる」場合、二つ(以上)の故障モードが混在している可能性があります。
二つ(以上)の故障モードに層別して、それぞれの故障モード毎にパラメータを決めていきましょう。
\[\]
「直線ではなく、曲線が引ける」場合、適当な位置パラメータ\(\gamma\) を設定し、プロットした点が直線上に並ぶように調整します。
\[\]
具体的には「調整後時間」の列を追加し、任意の位置パラメータ\(\gamma\) の値だけ「時間」から引きます。
\[↓↓↓\]
\[\]
⑤y切片(絶対値)を計算する
次にy切片(絶対値)を計算します。
※以下、位置パラメータ\(\gamma\) の調整がないものとして進めます。
\[\]
y切片(絶対値)の計算が必要だなんて、記事「ワイブル分布 とは」では説明していません。
(# ゚Д゚)
しかし、エクセルではワイブル分布確率紙の計算をするときに便利なので、計算します。(図7)
計算式は「=ABS(INTERCEPT(“Yの範囲”,LN(“故障時間の範囲”)))」 となります。
\[\]
⑥形状パラメータ\(m\) の計算
次に形状パラメータ\(m\) を計算します。(図8)
計算式は、「=SLOPE(“Yの範囲”,LN(“時間の範囲”))」です。
\[\]
SLOPE関数で二つのパラメータの傾きを計算します。
縦軸は”Yの範囲”で、横軸は”時間の範囲”⇒”\(X=lnt\)” となります。
\[\]
⑦尺度パラメータ\(\eta\) の計算
次に尺度パラメータ\(\eta\) を計算します。(図9)
計算式は、「=EXP(“Y切片(絶対値)”/”形状パラメータ \(m\)”)」です。
\[\]
この計算は少しわかりづらいので、解説します。
そもそも尺度パラメータ\(\eta\) とは何かというと、「\(Y=0\) の時の\(t\) の値」です。
\(Y=0\) の時の直線との交点を\(X\) とします。
\(X=lnt\) で表せるので、\(exp(X)=t\) となります。
\[\]
このとき、「\(X=\)(“Y切片(絶対値)”/”形状パラメータ \(m\)”)」で表せます。
形状パラメータ\(m\) は、直線の傾きを表します。
言い換えると、\(X\) が「1」変化した時の高さが「\(m\)」です。
では、\(X\) した時の高さはどうなるのでしょうか?
(。´・ω・)?
答えは、「Y切片(絶対値)」だけ変化します。
この時のイメージ図は下のようになります。
比で表すと、「\(1:m=X:Y\)」
つまり、「\(X=\frac{Y}{m}\)」 となるわけです。
\[\]
⑧予測値を計算する
形状パラメータ\(m\) と尺度パラメータ\(\eta\) を求められました。
後は予測値を計算すれば目標達成です。
\[\]
最初に累積故障確率を求めます。(図11)
計算式は、「=WEIBULL(“時間(予測値)”,”形状パラメータ \(m\)”,”尺度パラメータ\(\eta\)”,”TRUE”)」です。
\[\]
実際に予測時間を「780」として計算してみると、51.5% となります。(図12)
このように、WEIBULL関数を使うことで、任意の時間経過後の全体の故障率を見積もることができます。
わかりやすくパーセント表示にしてみます。(図13)
\[\]
・・・ (。´・ω・)?アレ?
左の表(実測データ)と比較すると、780時間の累積故障確率が違います!
この違いは、「実験データ+1を全体として考えた累積故障確率」と、「ワイブル分布から全体の分布を想定して見積もった累積故障確率」の違いになります。
\[\]
つぎに、累積故障確率から故障時間を見積もってみます。(図14)
計算式は、「=EXP((LN(LN(1/(1-“累積故障確率(予測値)”)))+”y切片(絶対値)”)/”形状パラメータ \(m\)”)」です。
\[\]
この計算をすることで、「全体の〇%が故障するまでの時間」を見積もることができるようになります。
実際に累積故障確率50%、すなわち半分が故障する時間を計算すると、765.51時間となります。(図15)
\[\]
(;´・ω・)
先ほどの「時間ー累積故障確率」の関係と同様に、実測データと予測値は異なることに気をつけてください。
\[\]
まとめ
①データ表を作成する
②累積故障頻度\(F(t(i))\) を計算する
③ワイブル分布縦軸\(Y\) を計算する
④散布図を描いて直線上になるか確認する
⑤y切片(絶対値)を計算する
⑥形状パラメータ\(m\) を計算する
⑦尺度パラメータ\(\eta\) を計算する
⑧時間⇒累積故障確率、累積故障確率⇒時間 の予測値を計算する
\[\]
今回は長めの記事となってしましました。
ワイブル分布は製品の寿命予測によく使われますが、日常生活でも活かせる場面があるかもしれません。
例えば、お店の繁盛具合なんかも分析できるかもしれませんね。
是非いろいろ試してみてください!
とても分かりやすい解説だと思います。さっそく使わせて(?)いただきました。ありがとうございます。
解説文の「計算式は、=EXP(LN...」の左括弧が一つ足りないので、エラーを解消するために式をよーく見ることができて、理解が深まりました。
ご指摘ありがとうございます。
計算式の修正をしました。
引き続き当ブログをご愛顧賜りますようよろしくお願いいたします。