100本ノック / 確率統計 / 確率・統計理論100本ノック

故障・品質・保全時間をどう見積もるか――製造業で使う連続分布10選

2020-04-13

故障・品質・保全時間をどう見積もるか――製造業で使う連続分布10選

概要

製造現場で観測する寸法、待ち時間、修理時間、寿命は、同じ「連続量」でも分布の形が異なります。本記事では、架空の精密部品工場を題材に、一様分布、正規分布、指数分布、ガンマ分布、ベータ分布、カイ二乗分布、t分布、F分布、対数正規分布、ワイブル分布をPythonで確認します。

目的は分布名を暗記することではありません。規格内率、予防保全周期、修理要員、少数データの不確実性、設備間のばらつき差といった意思決定に、どの分布をどの前提で使うかを整理することです。

[!NOTE] 本資料は、数理工房 (もしくは代表である和山個人) が過去に企業研修において使用した notebook を企業様の許可を得て再構成・編集のうえ公開しています。
掲載データはすべて架空のものであり、実在する企業・工場・数値とは一切関係ありません。

はじめに：この記事で扱う製造業の実務課題

今回の舞台は、精密シャフトを加工する架空工場です。品質保証部門は外径20.00 mm、規格19.95〜20.05 mmの適合率を管理し、保全部門は突発故障と修理時間から保全計画を立てています。管理者が答えたい問いは次の通りです。

工程能力から、規格外がどの程度発生しそうか
「前回の故障から長く動いた設備」は、直ちに壊れやすいのか
修理完了までの時間と、部品寿命の長い裾をどう見積もるか
小標本しかないとき、平均や分散をどこまで信頼できるか
設備Aと設備Bのばらつきに、実務上意味のある差があるか

連続分布は、これらの問いを「平均だけ」ではなく、確率、分位点、区間、リスクとして表現する共通言語です。

現場でよくある状況

月次会議では、平均外径、平均修理時間、平均故障間隔が一つの表に並びます。しかし、平均が同じでも意思決定は同じになりません。

外径は目標値を中心にほぼ左右対称だが、修理時間は右に長い裾を持つ
故障間隔は偶発故障期と摩耗故障期で、経過時間に対する故障率の動きが違う
新設備は標本数が少なく、推定値そのものの不確実性が大きい
全数検査できないため、検査位置や検査時刻の無作為化が必要になる

同じ平均・標準偏差を報告するだけでは、欠品、停止、品質流出の確率を比較できません。

なぜこの問題は判断が難しいのか

第一に、分布の選択には発生メカニズムの仮定が含まれます。多数の微小要因が加算される寸法は正規分布が候補になりますが、正の値だけを取り、複数工程の所要時間を合計する完了時間にはガンマ分布が自然です。

第二に、理論分布への適合と意思決定への有用性は同じではありません。ヒストグラムが似ていても、裾の確率を過小評価すれば保全部品や要員を不足させます。設備、製品、故障モードを混ぜると、単一分布の前提自体が崩れます。

第三に、推定値には標本誤差があります。平均の不確実性にはt分布、分散の不確実性にはカイ二乗分布、2群の分散比にはF分布が現れます。点推定だけで設備差を断定せず、区間と前提を合わせて判断する必要があります。

今回扱うノックの全体像

No.	連続分布	製造業での主な問い
051	一様分布	検査位置を偏りなく選べているか
052	正規分布	寸法の規格内率をどう見積もるか
053	指数分布	偶発故障の待ち時間をどう扱うか
054	ガンマ分布	複数作業の合計時間をどう見積もるか
055	ベータ分布	不良率の不確実性をどう表すか
056	カイ二乗分布	工程分散の推定誤差をどう評価するか
057	t分布	少数測定で母平均をどう推定するか
058	F分布	2設備の分散をどう比較するか
059	対数正規分布	右に裾を引く修理時間をどう計画へ入れるか
060	ワイブル分布	摩耗を含む寿命と保全周期をどう設計するか

各分布について、確率密度関数の形だけでなく、前提、Pythonによる確認、出力を意思決定へ変換する際の注意を説明します。

Python 環境の準備

NumPyで再現可能な架空データを生成し、pandasで集計、SciPyで確率分布と区間を計算し、matplotlibで可視化します。seabornや外部データは使用しません。乱数seedは固定します。

import sys

import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import scipy
from IPython.display import display
from scipy import stats

import japanize_matplotlib  # matplotlibの日本語表示を有効化

SEED = 20260711
rng = np.random.default_rng(SEED)

pd.set_option("display.max_columns", 20)
pd.set_option("display.float_format", lambda x: f"{x:,.4f}")

print(f"Python     : {sys.version.split()[0]}")
print(f"NumPy      : {np.__version__}")
print(f"pandas     : {pd.__version__}")
print(f"SciPy      : {scipy.__version__}")
print(f"matplotlib : {matplotlib.__version__}")
print(f"random seed: {SEED}")

Python     : 3.13.1
NumPy      : 2.5.1
pandas     : 3.0.3
SciPy      : 1.18.0
matplotlib : 3.11.0
random seed: 20260711

架空データの作成

分析用に、次の4種類のデータを作ります。

品質測定：設備A・Bで加工した精密シャフト各300本の外径
故障間隔：偶発故障を想定した300件の故障間隔
復旧作業：診断・部品交換・調整・確認という4工程の合計時間250件
寿命・修理：摩耗を含む部品寿命400件と、右裾の長い修理時間250件

生成時の分布を既知とするのは、理論値と標本から得た値を照合するためです。実務では、故障モードの層別、分布適合の診断、打切りデータの扱いを別途行います。

# 品質：設備ごとに平均とばらつきが異なる外径データ
n_per_machine = 300
quality_df = pd.DataFrame({
    "設備": np.repeat(["設備A", "設備B"], n_per_machine),
    "外径_mm": np.concatenate([
        rng.normal(20.000, 0.018, n_per_machine),
        rng.normal(20.006, 0.027, n_per_machine),
    ]),
})
quality_df["規格内"] = quality_df["外径_mm"].between(19.95, 20.05)

# 検査・保全：位置、故障間隔、復旧4工程、修理時間、摩耗寿命
inspection_df = pd.DataFrame({"コイル上の検査位置_m": rng.uniform(0, 100, 500)})
failure_df = pd.DataFrame({"偶発故障間隔_h": rng.exponential(96, 300)})
task_times = rng.exponential(0.8, size=(250, 4))
recovery_df = pd.DataFrame(task_times, columns=["診断_h", "交換_h", "調整_h", "確認_h"])
recovery_df["復旧合計_h"] = recovery_df.sum(axis=1)
repair_df = pd.DataFrame({"修理時間_h": rng.lognormal(np.log(2.5), 0.65, 250)})
life_df = pd.DataFrame({"部品寿命_h": 1200 * rng.weibull(1.8, 400)})

quality_summary = quality_df.groupby("設備").agg(
    測定数=("外径_mm", "size"),
    平均外径_mm=("外径_mm", "mean"),
    標準偏差_mm=("外径_mm", "std"),
    規格内率=("規格内", "mean"),
)
display(quality_summary)

data_summary = pd.DataFrame({
    "データ": ["検査位置", "偶発故障間隔", "復旧合計時間", "修理時間", "部品寿命"],
    "件数": [len(inspection_df), len(failure_df), len(recovery_df), len(repair_df), len(life_df)],
    "平均": [inspection_df.iloc[:, 0].mean(), failure_df.iloc[:, 0].mean(),
           recovery_df["復旧合計_h"].mean(), repair_df.iloc[:, 0].mean(), life_df.iloc[:, 0].mean()],
    "中央値": [inspection_df.iloc[:, 0].median(), failure_df.iloc[:, 0].median(),
            recovery_df["復旧合計_h"].median(), repair_df.iloc[:, 0].median(), life_df.iloc[:, 0].median()],
    "単位": ["m", "h", "h", "h", "h"],
})
display(data_summary.round(3))

	測定数	平均外径_mm	標準偏差_mm	規格内率
設備
設備A	300	19.9999	0.0189	0.9933
設備B	300	20.0069	0.0260	0.9600

	データ	件数	平均	中央値	単位
0	検査位置	500	47.4910	46.2350	m
1	偶発故障間隔	300	101.9570	69.6000	h
2	復旧合計時間	250	3.1670	2.8150	h
3	修理時間	250	3.2540	2.4820	h
4	部品寿命	400	1,030.8990	941.3880	h

No.051：一様分布

実務での意味

長さ100 mのコイルから検査位置を無作為に1点選ぶとき、すべての位置が同じ選択機会を持つ設計が一様分布です。先頭付近だけを測る運用では、巻き終わり側の異常を見落とす可能性があります。一様乱数は「品質値が一様」という仮定ではなく、検査機会を空間や時間へ公平に配る仕組みとして使います。

分析・モデル化の考え方

区間 $a\le x\le b$ の連続一様分布の確率密度は

f(x)=\frac{1}{b-a}, \qquad E[X]=\frac{a+b}{2}, \qquad \mathrm{Var}(X)=\frac{(b-a)^2}{12}

です。20〜40 mを選ぶ確率は区間長の比 $(40-20)/(100-0)=0.2$ です。端点ちょうどの確率は0ですが、区間には正の確率があります。

Pythonで確認する

positions = inspection_df["コイル上の検査位置_m"]
segment_rate = positions.between(20, 40).mean()
uniform_result = pd.DataFrame({
    "指標": ["標本平均(m)", "理論平均(m)", "20〜40mの標本比率", "20〜40mの理論確率"],
    "値": [positions.mean(), 50.0, segment_rate, 0.2],
})
display(uniform_result.round(4))

fig, ax = plt.subplots(figsize=(8, 4))
ax.hist(positions, bins=10, range=(0, 100), density=True,
        color="steelblue", edgecolor="white", alpha=0.8, label="検査位置")
ax.axhline(1 / 100, color="darkred", linestyle="--", label="理論密度 1/100")
ax.set_title("一様乱数で選んだコイル検査位置")
ax.set_xlabel("コイル上の位置 (m)")
ax.set_ylabel("確率密度")
ax.grid(alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()

	指標	値
0	標本平均(m)	47.4910
1	理論平均(m)	50.0000
2	20〜40mの標本比率	0.1940
3	20〜40mの理論確率	0.2000

png

結果の読み取り

500点の検査位置は区間全体へおおむね均等に分布し、20〜40 mの標本比率も理論値0.2に近づきます。ただし有限標本なので各ビンの高さは完全には揃いません。実務では乱数を作るだけでなく、停止中・段取り中など選べない時間がないか、選ばれた位置が現場都合で差し替えられていないかを監査します。周期的な異常が疑われる場合は、単純無作為抽出だけでなく層別抽出も検討します。

No.052：正規分布

実務での意味

加工寸法は、温度、工具状態、材料、測定誤差など多数の小さな要因が加算的に作用すると、中心付近が厚く左右対称な正規分布で近似できることがあります。平均と標準偏差から規格外確率を計算できるため、工程能力や検査負荷の見積もりに有用です。

分析・モデル化の考え方

平均 $\mu$ 、標準偏差 $\sigma$ の正規分布の密度は

f(x)=\frac{1}{\sqrt{2\pi}\sigma} \exp\left\{-\frac{(x-\mu)^2}{2\sigma^2}\right\}.

規格下限LSL、上限USLに入る確率は、標準正規分布の累積分布関数 $\Phi$ を用いて

P(\mathrm{LSL}\le X\le\mathrm{USL}) =\Phi\!\left(\frac{\mathrm{USL}-\mu}{\sigma}\right) -\Phi\!\left(\frac{\mathrm{LSL}-\mu}{\sigma}\right)

と表せます。ここでは設備Bを対象に、標本平均・標本標準偏差を代入します。

Pythonで確認する

lsl, usl = 19.95, 20.05
x_b = quality_df.loc[quality_df["設備"] == "設備B", "外径_mm"]
mu_hat, sigma_hat = x_b.mean(), x_b.std(ddof=1)
predicted_yield = stats.norm.cdf(usl, mu_hat, sigma_hat) - stats.norm.cdf(lsl, mu_hat, sigma_hat)
observed_yield = x_b.between(lsl, usl).mean()

display(pd.DataFrame({
    "平均外径_mm": [mu_hat], "標準偏差_mm": [sigma_hat],
    "正規分布による規格内率": [predicted_yield], "観測規格内率": [observed_yield],
}).round(4))

grid = np.linspace(19.90, 20.10, 500)
fig, ax = plt.subplots(figsize=(8, 4))
ax.hist(x_b, bins=24, density=True, alpha=0.65, color="slateblue",
        edgecolor="white", label="設備Bの測定値")
ax.plot(grid, stats.norm.pdf(grid, mu_hat, sigma_hat), color="darkred", label="推定正規密度")
ax.axvline(lsl, color="black", linestyle="--", label="規格限界")
ax.axvline(usl, color="black", linestyle="--")
ax.set_title("設備Bの外径分布と規格限界")
ax.set_xlabel("外径 (mm)")
ax.set_ylabel("確率密度")
ax.grid(alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()

	平均外径_mm	標準偏差_mm	正規分布による規格内率	観測規格内率
0	20.0069	0.0260	0.9365	0.9600

png

結果の読み取り

設備Bは平均が規格中心より上側にあり、ばらつきもあるため、上限側の規格外に注意が必要です。推定正規分布による規格内率と観測率は近いものの、正規性は自動的に保証されません。量産データでは時系列のドリフト、ロット差、混合分布、外れ値を確認し、安定状態を確認してから将来の規格内率へ外挿します。平均を中心へ戻す施策と、標準偏差を小さくする施策は分けて管理します。

No.053：指数分布

実務での意味

偶発故障が一定の発生率で起こる期間では、次の故障までの待ち時間を指数分布で近似できます。保全要員の呼出頻度や、一定期間を無故障で運転できる確率の初期見積もりに使えます。

分析・モデル化の考え方

故障率 $\lambda>0$ の指数分布は

f(t)=\lambda e^{-\lambda t},\quad S(t)=P(T>t)=e^{-\lambda t},\quad E[T]=\frac{1}{\lambda}

です。重要な性質は無記憶性

P(T>s+t\mid T>s)=P(T>t)

です。これは「長く動いた設備ほど壊れやすい」という摩耗を表しません。摩耗がある場合はNo.060のワイブル分布などを検討します。

Pythonで確認する

intervals = failure_df["偶発故障間隔_h"]
mean_hat = intervals.mean()
lambda_hat = 1 / mean_hat
s, t = 72, 48
conditional_emp = (intervals > s + t).sum() / (intervals > s).sum()

display(pd.DataFrame({
    "平均故障間隔_h": [mean_hat],
    "推定故障率_1/h": [lambda_hat],
    "P(T>48)理論": [np.exp(-lambda_hat * t)],
    "P(T>120|T>72)標本": [conditional_emp],
}).round(4))

time_grid = np.linspace(0, 400, 500)
fig, ax = plt.subplots(figsize=(8, 4))
ax.hist(intervals, bins=25, density=True, alpha=0.7, color="teal",
        edgecolor="white", label="故障間隔")
ax.plot(time_grid, stats.expon.pdf(time_grid, scale=mean_hat), color="darkred", label="推定指数密度")
ax.set_title("偶発故障を想定した故障間隔")
ax.set_xlabel("故障間隔 (h)")
ax.set_ylabel("確率密度")
ax.grid(alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()

	平均故障間隔_h	推定故障率_1/h	P(T>48)理論	P(T>120\|T>72)標本
0	101.9568	0.0098	0.6245	0.6351

png

結果の読み取り

故障間隔は右に裾を持ち、標本平均の逆数から故障率を推定できます。無記憶性の標本確認は有限件数のため理論値と完全には一致しません。実務上の核心は、指数分布を採用すると経過時間にかかわらず瞬間故障率を一定と仮定する点です。故障モードを混ぜず、摩耗・初期故障・保全条件の変化がある場合は区間やモードを分けます。

No.054：ガンマ分布

実務での意味

復旧作業が診断、交換、調整、確認のような複数段階からなり、各段階の所要時間が正の値で変動するとき、合計時間はガンマ分布で表せる場合があります。平均だけでなく「4時間以内に復旧できる確率」や要員拘束時間の上側分位点を計画できます。

分析・モデル化の考え方

形状母数 $k$ 、尺度母数 $\theta$ のガンマ分布の密度は

f(x)=\frac{x^{k-1}e^{-x/\theta}}{\Gamma(k)\theta^k},\quad x>0,

E[X]=k\theta,\qquad \mathrm{Var}(X)=k\theta^2.

同じ尺度 $\theta$ の独立な指数分布を $k$ 個足すと、形状 $k$ のガンマ分布になります。ここでは4作業、各平均0.8時間なので、理論平均は3.2時間です。

Pythonで確認する

recovery = recovery_df["復旧合計_h"]
k, theta = 4, 0.8
gamma_result = pd.DataFrame({
    "指標": ["標本平均(h)", "理論平均(h)", "標本90%分位(h)", "理論90%分位(h)", "4時間以内の理論確率"],
    "値": [recovery.mean(), k * theta, recovery.quantile(0.9),
          stats.gamma.ppf(0.9, a=k, scale=theta), stats.gamma.cdf(4, a=k, scale=theta)],
})
display(gamma_result.round(4))

grid = np.linspace(0, recovery.max() * 1.05, 500)
fig, ax = plt.subplots(figsize=(8, 4))
ax.hist(recovery, bins=24, density=True, alpha=0.7, color="darkorange",
        edgecolor="white", label="復旧合計時間")
ax.plot(grid, stats.gamma.pdf(grid, a=k, scale=theta), color="navy", label="Gamma(4, 0.8)")
ax.axvline(stats.gamma.ppf(0.9, a=k, scale=theta), color="darkred", linestyle="--", label="理論90%分位")
ax.set_title("4作業の復旧合計時間とガンマ分布")
ax.set_xlabel("復旧合計時間 (h)")
ax.set_ylabel("確率密度")
ax.grid(alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()

	指標	値
0	標本平均(h)	3.1670
1	理論平均(h)	3.2000
2	標本90%分位(h)	5.3064
3	理論90%分位(h)	5.3446
4	4時間以内の理論確率	0.7350

png

結果の読み取り

復旧合計時間は0未満を取らず、右に裾を持ちます。標本平均と90%分位は理論値に近く、要員計画では平均3.2時間だけでなく、約9割が収まる時間も使えます。ただし、各作業が独立で同じ尺度を持つという設定は簡略化です。部品待ちや承認待ちが混ざる場合は、待ち時間を別KPIに分けるか、混合分布・シミュレーションで扱います。

No.055：ベータ分布

実務での意味

不良率や成功率は0から1の範囲にあります。新製品の立上げ直後のようにデータが少ない場合、「不良率は2%」という一点だけでなく、不良率そのものの不確実性をベータ分布で表すと、追加検査や出荷判断を確率で議論できます。

分析・モデル化の考え方

形状母数 $\alpha,\beta>0$ のベータ分布は

f(p)=\frac{p^{\alpha-1}(1-p)^{\beta-1}}{B(\alpha,\beta)},\quad 0<p<1,

E[P]=\frac{\alpha}{\alpha+\beta},\qquad \mathrm{Var}(P)=\frac{\alpha\beta}{(\alpha+\beta)^2(\alpha+\beta+1)}

です。二項データに対する共役事前分布でもあり、事前分布 $\mathrm{Beta}(\alpha_0,\beta_0)$ に、不良 $d$ 個・良品 $n-d$ 個を観測すると、事後分布は $\mathrm{Beta}(\alpha_0+d,\beta_0+n-d)$ になります。ここでは過去知見を $\mathrm{Beta}(2,98)$ 、新規検査を200個中4個不良とします。

Pythonで確認する

alpha0, beta0 = 2, 98
n, defects = 200, 4
alpha1, beta1 = alpha0 + defects, beta0 + n - defects
ci_low, ci_high = stats.beta.ppf([0.025, 0.975], alpha1, beta1)

beta_result = pd.DataFrame({
    "分布": ["事前", "事後"],
    "alpha": [alpha0, alpha1], "beta": [beta0, beta1],
    "平均不良率": [alpha0 / (alpha0 + beta0), alpha1 / (alpha1 + beta1)],
})
display(beta_result.round(4))
print(f"事後95%区間: {ci_low:.4%} 〜 {ci_high:.4%}")
print(f"不良率が3%を超える事後確率: {stats.beta.sf(0.03, alpha1, beta1):.2%}")

p_grid = np.linspace(0, 0.08, 500)
fig, ax = plt.subplots(figsize=(8, 4))
ax.plot(p_grid, stats.beta.pdf(p_grid, alpha0, beta0), label="事前 Beta(2, 98)")
ax.plot(p_grid, stats.beta.pdf(p_grid, alpha1, beta1), label="事後 Beta(6, 294)")
ax.axvspan(ci_low, ci_high, color="orange", alpha=0.2, label="事後95%区間")
ax.set_title("検査データ更新前後の不良率分布")
ax.set_xlabel("不良率")
ax.set_ylabel("確率密度")
ax.grid(alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()

	分布	alpha	beta	平均不良率
0	事前	2	98	0.0200
1	事後	6	294	0.0200

事後95%区間: 0.7399% 〜 3.8591%
不良率が3%を超える事後確率: 11.38%

png

結果の読み取り

観測不良率は4/200=2%ですが、真の不良率を一点に固定せず、95%区間と「3%を超える確率」を得られます。出荷判定へ使う場合は、3%という管理基準、誤判定の損失、事前分布の根拠を関係者と合意します。事前分布を変えた感度分析も必要です。ベータ分布は個々の製品品質値ではなく、0〜1に制約された未知の率を表している点に注意します。

No.056：カイ二乗分布

実務での意味

標準偏差は工程のばらつきを表しますが、少数標本から計算した標本分散そのものも変動します。測定15本だけで「標準偏差は0.02 mm」と断定せず、母分散が取り得る範囲を示すためにカイ二乗分布を使います。

分析・モデル化の考え方

母集団が正規分布で、標本サイズを $n$ 、不偏標本分散を $S^2$ 、母分散を $\sigma^2$ とすると

\frac{(n-1)S^2}{\sigma^2}\sim\chi^2_{n-1}.

したがって、信頼係数 $1-\alpha$ の母分散の信頼区間は

\left[ \frac{(n-1)S^2}{\chi^2_{1-\alpha/2,n-1}}, \frac{(n-1)S^2}{\chi^2_{\alpha/2,n-1}} \right]

です。区間は左右対称ではなく、標本が少ないほど広くなります。

Pythonで確認する

sigma_true = 0.020
n = 15
simulations = 5000
samples = rng.normal(20.0, sigma_true, size=(simulations, n))
sample_vars = samples.var(axis=1, ddof=1)
chi_stats = (n - 1) * sample_vars / sigma_true**2

one_sample = samples[0]
s2 = one_sample.var(ddof=1)
df_chi = n - 1
var_low = df_chi * s2 / stats.chi2.ppf(0.975, df_chi)
var_high = df_chi * s2 / stats.chi2.ppf(0.025, df_chi)

display(pd.DataFrame({
    "標本標準偏差_mm": [np.sqrt(s2)],
    "母標準偏差95%下限_mm": [np.sqrt(var_low)],
    "母標準偏差95%上限_mm": [np.sqrt(var_high)],
    "真の標準偏差_mm": [sigma_true],
}).round(5))

grid = np.linspace(0, stats.chi2.ppf(0.995, df_chi), 500)
fig, ax = plt.subplots(figsize=(8, 4))
ax.hist(chi_stats, bins=40, density=True, alpha=0.7, color="mediumpurple",
        edgecolor="white", label="シミュレーション統計量")
ax.plot(grid, stats.chi2.pdf(grid, df_chi), color="darkred", label=f"カイ二乗分布 (自由度{df_chi})")
ax.set_title("標本分散から作る統計量の分布")
ax.set_xlabel(r"$(n-1)S^2/\sigma^2$")
ax.set_ylabel("確率密度")
ax.grid(alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()

	標本標準偏差_mm	母標準偏差95%下限_mm	母標準偏差95%上限_mm	真の標準偏差_mm
0	0.0177	0.0130	0.0280	0.0200

png

結果の読み取り

15本の標本から得た標準偏差の95%区間は、点推定よりかなり広くなります。シミュレーションした統計量は自由度14のカイ二乗分布に沿い、理論式の意味を確認できます。この区間は母集団の正規性と独立性に依存します。工程のドリフト、自己相関、外れ値がある場合は、単純な区間を使う前に工程状態を見直します。工程能力指数の信頼性も標本数に左右されるため、点推定だけで合否を決めません。

No.057：t分布

実務での意味

新設備の立上げで測定が12本しかない場合、母標準偏差は未知です。標本平均を標本標準偏差で標準化すると、その不確実性は標準正規分布より裾の厚いt分布で表されます。少数データの平均値を過度に確実視しないための分布です。

分析・モデル化の考え方

正規母集団からの独立標本 $X_1,\ldots,X_n$ に対し

T=\frac{\bar X-\mu}{S/\sqrt{n}}\sim t_{n-1}.

母平均の $100(1-\alpha)\%$ 信頼区間は

\bar X\pm t_{1-\alpha/2,n-1}\frac{S}{\sqrt n}

です。自由度が増えるとt分布は標準正規分布へ近づきます。

Pythonで確認する

startup_sample = quality_df.loc[quality_df["設備"] == "設備B", "外径_mm"].iloc[:12].to_numpy()
n_t = len(startup_sample)
mean_t = startup_sample.mean()
se_t = startup_sample.std(ddof=1) / np.sqrt(n_t)
t_critical = stats.t.ppf(0.975, df=n_t - 1)
z_critical = stats.norm.ppf(0.975)
t_ci = (mean_t - t_critical * se_t, mean_t + t_critical * se_t)
z_ci = (mean_t - z_critical * se_t, mean_t + z_critical * se_t)

display(pd.DataFrame({
    "方法": ["t分布", "標準正規分布（比較用）"],
    "95%下限_mm": [t_ci[0], z_ci[0]],
    "95%上限_mm": [t_ci[1], z_ci[1]],
    "区間幅_mm": [t_ci[1] - t_ci[0], z_ci[1] - z_ci[0]],
}).round(5))

grid = np.linspace(-4.5, 4.5, 500)
fig, ax = plt.subplots(figsize=(8, 4))
ax.plot(grid, stats.t.pdf(grid, df=n_t - 1), label=f"t分布 (自由度{n_t - 1})")
ax.plot(grid, stats.norm.pdf(grid), linestyle="--", label="標準正規分布")
ax.set_title("少数標本で用いるt分布の厚い裾")
ax.set_xlabel("標準化した値")
ax.set_ylabel("確率密度")
ax.grid(alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()

	方法	95%下限_mm	95%上限_mm	区間幅_mm
0	t分布	19.9970	20.0302	0.0332
1	標準正規分布（比較用）	19.9988	20.0283	0.0295

png

結果の読み取り

t分布による区間は、同じ標準誤差へ標準正規分布を機械的に当てた区間より広くなります。これは少数標本で母標準偏差も推定している不確実性を反映します。区間が規格中心や管理上の許容差をまたぐなら、「差がない」と結論するのではなく、現時点では判断精度が足りないと読み、追加測定の価値を検討します。非正規性が強い少数標本ではt法も安定しないため、測定過程と分布形状の確認が先です。

No.058：F分布

実務での意味

2台の設備で平均寸法が同程度でも、ばらつきが違えば規格外リスクは変わります。正規母集団を仮定した2群の分散比はF分布に従い、設備間のばらつき差を区間と検定で評価できます。

分析・モデル化の考え方

独立な2標本の不偏分散を $S_1^2,S_2^2$ 、母分散を $\sigma_1^2,\sigma_2^2$ とすると

\frac{S_1^2/\sigma_1^2}{S_2^2/\sigma_2^2}\sim F_{n_1-1,n_2-1}.

母分散が等しいという帰無仮説では $F=S_1^2/S_2^2$ をF分布と比較します。ただしF検定は非正規性や外れ値に敏感です。実務では箱ひげ図、時系列、Levene検定などの頑健な方法も併用します。

Pythonで確認する

a = quality_df.loc[quality_df["設備"] == "設備A", "外径_mm"].iloc[:60].to_numpy()
b = quality_df.loc[quality_df["設備"] == "設備B", "外径_mm"].iloc[:60].to_numpy()
var_a, var_b = a.var(ddof=1), b.var(ddof=1)
f_stat = var_b / var_a
df1, df2 = len(b) - 1, len(a) - 1
p_value = 2 * min(stats.f.cdf(f_stat, df1, df2), stats.f.sf(f_stat, df1, df2))
ratio_low = f_stat / stats.f.ppf(0.975, df1, df2)
ratio_high = f_stat / stats.f.ppf(0.025, df1, df2)

display(pd.DataFrame({
    "設備B/Aの標本分散比": [f_stat],
    "母分散比95%下限": [ratio_low],
    "母分散比95%上限": [ratio_high],
    "両側p値": [min(p_value, 1.0)],
}).round(4))

grid = np.linspace(0, stats.f.ppf(0.995, df1, df2), 500)
fig, ax = plt.subplots(figsize=(8, 4))
ax.plot(grid, stats.f.pdf(grid, df1, df2), color="navy", label=f"F({df1}, {df2})")
ax.axvline(f_stat, color="darkred", linestyle="--", label=f"観測分散比 {f_stat:.2f}")
ax.set_title("等分散仮説のF分布と観測分散比")
ax.set_xlabel("設備B / 設備A の分散比")
ax.set_ylabel("確率密度")
ax.grid(alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()

	設備B/Aの標本分散比	母分散比95%下限	母分散比95%上限	両側p値
0	1.9165	1.1448	3.2085	0.0136

png

結果の読み取り

設備B/Aの分散比は1より大きく、95%区間とp値から、今回の設定では設備Bのばらつきが大きいことを示す材料が得られます。しかし統計的な差だけでは改善投資の優先順位は決まりません。分散差が規格外個数、手直し時間、顧客影響へ換算するとどれだけの損失になるかを評価します。また、設備別の製品構成や測定器が同じかを確認し、外れ値一件で結論が変わらないかも点検します。

No.059：対数正規分布

実務での意味

修理時間は0未満にならず、多くは短時間で終わる一方、診断難航や部品手配で非常に長くなることがあります。複数の倍率要因が積み重なる量は、対数を取ると正規分布に近づき、元の尺度では対数正規分布になる場合があります。

分析・モデル化の考え方

$Y=\log X$ が $N(\mu,\sigma^2)$ に従うとき、 $X$ は対数正規分布に従い、

f(x)=\frac{1}{x\sigma\sqrt{2\pi}} \exp\left\{-\frac{(\log x-\mu)^2}{2\sigma^2}\right\},\quad x>0,

\mathrm{median}(X)=e^\mu,\qquad E[X]=e^{\mu+\sigma^2/2}

です。右裾が長いため平均は中央値より大きくなります。典型時間には中央値、人員・停止損失には平均や上側分位点を使うなど、目的でKPIを分けます。

Pythonで確認する

repair = repair_df["修理時間_h"]
log_repair = np.log(repair)
mu_log, sigma_log = log_repair.mean(), log_repair.std(ddof=1)
mean_model = np.exp(mu_log + sigma_log**2 / 2)
median_model = np.exp(mu_log)
p95_model = stats.lognorm.ppf(0.95, s=sigma_log, scale=np.exp(mu_log))

display(pd.DataFrame({
    "指標": ["標本平均", "標本中央値", "モデル平均", "モデル中央値", "モデル95%分位"],
    "修理時間_h": [repair.mean(), repair.median(), mean_model, median_model, p95_model],
}).round(3))

grid = np.linspace(0.01, repair.quantile(0.995) * 1.2, 500)
fig, axes = plt.subplots(1, 2, figsize=(11, 4))
axes[0].hist(repair, bins=25, density=True, alpha=0.7, color="seagreen", edgecolor="white")
axes[0].plot(grid, stats.lognorm.pdf(grid, s=sigma_log, scale=np.exp(mu_log)), color="darkred")
axes[0].set_title("元尺度：右裾の長い修理時間")
axes[0].set_xlabel("修理時間 (h)")
axes[0].set_ylabel("確率密度")
axes[0].grid(alpha=0.3)
axes[1].hist(log_repair, bins=25, density=True, alpha=0.7, color="slateblue", edgecolor="white")
log_grid = np.linspace(log_repair.min(), log_repair.max(), 500)
axes[1].plot(log_grid, stats.norm.pdf(log_grid, mu_log, sigma_log), color="darkred")
axes[1].set_title("対数尺度：ほぼ左右対称")
axes[1].set_xlabel("log(修理時間)")
axes[1].set_ylabel("確率密度")
axes[1].grid(alpha=0.3)
plt.tight_layout()
plt.show()

	指標	修理時間_h
0	標本平均	3.2540
1	標本中央値	2.4820
2	モデル平均	3.2810
3	モデル中央値	2.6130
4	モデル95%分位	7.9280

png

結果の読み取り

元尺度では平均が中央値を上回り、少数の長時間修理が平均停止時間を押し上げます。中央値だけで要員や損失を見積もると長時間案件を過小評価し、最大値だけでは偶然に振り回されます。平均、中央値、90%・95%分位を役割別に報告するのが有効です。部品待ちで別の山が生じるなら、単一の対数正規分布に押し込まず、「実作業時間」と「待ち時間」を分けます。

No.060：ワイブル分布

実務での意味

軸受、シール、工具などの寿命では、経過時間とともに故障しやすさが変わります。ワイブル分布は形状母数により、初期故障、偶発故障、摩耗故障を表現でき、予防交換周期や保証期間の検討に広く使われます。

分析・モデル化の考え方

形状母数 $k$ 、尺度母数 $\lambda$ のワイブル分布は

f(t)=\frac{k}{\lambda}\left(\frac{t}{\lambda}\right)^{k-1} \exp\left\{-\left(\frac{t}{\lambda}\right)^k\right\},

S(t)=\exp\left\{-\left(\frac{t}{\lambda}\right)^k\right\},\qquad h(t)=\frac{k}{\lambda}\left(\frac{t}{\lambda}\right)^{k-1}.

$k<1$ は故障率低下、 $k=1$ は一定、 $k>1$ は故障率上昇を表します。尺度 $\lambda$ は生存率が $e^{-1}\approx36.8\%$ になる時点であり、平均寿命そのものではありません。

Pythonで確認する

life = life_df["部品寿命_h"]
# 打切りなしの架空データに対して、位置0固定で母数を最尤推定
shape_hat, loc_hat, scale_hat = stats.weibull_min.fit(life, floc=0)
replace_at_10pct_failure = scale_hat * (-np.log(0.90)) ** (1 / shape_hat)
median_life = stats.weibull_min.ppf(0.5, shape_hat, scale=scale_hat)

display(pd.DataFrame({
    "推定形状k": [shape_hat], "推定尺度lambda_h": [scale_hat],
    "推定中央値_h": [median_life],
    "累積故障10%となる時間_h": [replace_at_10pct_failure],
}).round(2))

time_grid = np.linspace(1, 2200, 500)
survival = stats.weibull_min.sf(time_grid, shape_hat, scale=scale_hat)
hazard = (shape_hat / scale_hat) * (time_grid / scale_hat) ** (shape_hat - 1)

fig, axes = plt.subplots(1, 2, figsize=(11, 4))
axes[0].plot(time_grid, survival, color="navy")
axes[0].axvline(replace_at_10pct_failure, color="darkred", linestyle="--", label="累積故障10%")
axes[0].set_title("推定ワイブル生存関数")
axes[0].set_xlabel("使用時間 (h)")
axes[0].set_ylabel("生存確率")
axes[0].grid(alpha=0.3)
axes[0].legend()
axes[1].plot(time_grid, hazard, color="darkorange")
axes[1].set_title("推定ハザード関数")
axes[1].set_xlabel("使用時間 (h)")
axes[1].set_ylabel("瞬間故障率 (1/h)")
axes[1].grid(alpha=0.3)
plt.tight_layout()
plt.show()

	推定形状k	推定尺度lambda_h	推定中央値_h	累積故障10%となる時間_h
0	1.7600	1,156.4200	939.4800	322.9300

png

結果の読み取り

推定形状母数は1より大きく、使用時間とともに瞬間故障率が上昇する摩耗型を示します。「累積故障を10%以内に抑える」というサービス水準なら、その条件に対応する時間が交換周期の候補になります。ただし、ここで計算した時点をそのまま採用するのではなく、計画停止費、突発停止費、交換後の初期故障、在庫、交換作業能力を含む総費用で最適化します。実務の寿命データには観測期間内に壊れなかった右打切りが多いため、単純に故障品だけへ分布を当てると寿命を短く見積もります。

対象ノックを通して見える実務上の示唆

No.051〜No.060を通して、連続分布は次の4つの役割に整理できます。

観測を設計する：一様分布で検査位置・時刻の選択機会を公平にする
現象を表現する：正規、指数、ガンマ、対数正規、ワイブル分布を発生メカニズムに合わせる
未知量の不確実性を表す：ベータ分布で不良率、カイ二乗分布で分散、t分布で少数標本の平均を扱う
工程を比較する：F分布で分散比を評価し、差を規格外・損失へ換算する

分布選択の順番は「ヒストグラムに似た曲線を探す」だけでは不十分です。対象量の範囲、生成過程、独立性、時間変化、打切り、意思決定で重視する裾を先に確認します。

判断対象	中心だけでは不足する理由	併記したい指標
寸法品質	平均が中心でもばらつきで規格外が出る	標準偏差、規格内率、時系列
復旧計画	右裾により平均超過が頻発する	中央値、90%・95%分位
新工程評価	少数標本では推定値が大きく動く	信頼区間、追加測定数
予防保全	故障率が使用時間で変化する	生存確率、ハザード、総費用

実務導入する場合に必要なこと

1. 分析単位とデータ定義を固定する

設備ID、製品、材料ロット、測定器、故障モード、停止開始・復旧時刻、計画停止か突発停止かを定義します。修理時間に部品待ちを含めるかなど、KPIの開始・終了条件を統一します。

2. 層別と時系列を先に確認する

異なる設備・製品・故障モードを混ぜると、単一分布に見えなくなります。管理図や時系列で工程の安定性を確認し、層別後の分布を評価します。相関した連続測定を独立標本として扱わないことも重要です。

3. 適合度だけでなく、意思決定する裾を検証する

Q-Qプロット、経験累積分布、残差、AICなどを使いながら、規格限界や95%分位付近の誤差を確認します。中心部への適合が良くても、停止損失を決める右裾が合わなければ目的に適しません。

4. 打切り・欠測・選択バイアスを管理する

寿命分析では未故障品を右打切りとして保持します。重大故障だけが記録される、短時間停止が欠ける、検査位置が現場判断で変わるといった記録過程もデータ化します。

5. 確率を費用と運用ルールへ接続する

規格外率、停止確率、交換周期を、廃棄費、納期遅延、予備品、保全要員、顧客影響へ変換します。モデルを定期更新し、予測確率と実績のずれを監視する責任者・頻度・判断基準を決めます。

まとめ

一様分布は検査位置などの無作為化、正規分布は安定した寸法ばらつきの近似に使える
指数分布は一定故障率、ガンマ分布は複数の正の待ち時間の合計を表現する
ベータ分布は0〜1の率の不確実性を表し、検査結果で更新できる
カイ二乗分布、t分布、F分布は、分散・少数標本の平均・分散比の推論を支える
対数正規分布は右裾の長い修理時間、ワイブル分布は時間で変化する故障率に向く
分布は発生メカニズムと意思決定目的から選び、層別、時系列、打切り、裾の適合を確認する

連続分布を実務で使う価値は、現場のばらつきを一つの平均へ潰さず、「どの事象が、いつまでに、どの程度の確率で起きるか」へ変換できることにあります。

法人向けのご相談

数理工房では、製造業の品質・生産・保全データを対象に、以下のような支援を行っています。

工程能力、規格外率、検査設計の可視化
故障・修理・寿命データの分布分析と予防保全設計
統計モデル、シミュレーション、KPI設計の伴走支援
製造業の実データに合わせたPython・統計研修
分析結果を現場運用や管理会議へ定着させる仕組みづくり

📩 お問い合わせ: surikobo.co.jp/contact
まずはお気軽にご相談ください。

故障・品質・保全時間をどう見積もるか――製造業で使う連続分布10選

故障・品質・保全時間をどう見積もるか――製造業で使う連続分布10選

概要

はじめに：この記事で扱う製造業の実務課題

現場でよくある状況

なぜこの問題は判断が難しいのか

今回扱うノックの全体像

Python 環境の準備

架空データの作成

No.051：一様分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.052：正規分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.053：指数分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.054：ガンマ分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.055：ベータ分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.056：カイ二乗分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.057：t分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.058：F分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.059：対数正規分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

No.060：ワイブル分布

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

対象ノックを通して見える実務上の示唆

実務導入する場合に必要なこと

1. 分析単位とデータ定義を固定する

2. 層別と時系列を先に確認する

3. 適合度だけでなく、意思決定する裾を検証する

4. 打切り・欠測・選択バイアスを管理する

5. 確率を費用と運用ルールへ接続する

まとめ

法人向けのご相談

Related Posts

工程変更を「効いた気がする」で終わらせない――統計的推論で品質投資を判断する

製造業の統計推定をPythonで学ぶ｜抜取検査からEMアルゴリズムまで

標本数を、品質判断の根拠に変える――製造KPIの極限定理10本ノック

複数KPIを同時に読む――製造品質・設備負荷・納期リスクの多変量確率

製造業の不良・設備停止を離散分布で予測する｜Python実践10本ノック