100本ノック / 確率統計 / 確率・統計マーケティング応用100本ノック

製造業の品質管理をPythonで実践｜不良率・Cp/Cpk・管理図10本ノック

2020-05-21

品質データを「判断」に変える：製造業の確率・統計10本ノック

本記事は、製造業の品質管理を題材に、不良率の把握から工程能力、管理図、検定、信頼性、将来リスクのシミュレーションまでを一つの架空データでつなぐ実践 notebook です。数式を計算すること自体ではなく、「いつ工程を止めるか」「どこへ改善資源を配分するか」「顧客リスクをどう説明するか」という意思決定に統計を結び付けます。

この100本ノックでは、品質管理、設備保全、需要予測、在庫最適化、生産計画、マーケティングサイエンス、機械学習、意思決定科学、シミュレーション、数理工房DI実践を順に扱います。最初の10本は、その土台となる品質管理です。Pythonで再現可能な分析を作り、現場・品質保証・経営が同じ数字を見て議論できる状態を目指します。

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

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

架空の精密部品工場では、3ラインで同一部品を生産しています。重要寸法の規格は 49.5〜50.5 mm、各ロットから50個を測定し、外観検査では各ロット200個を検査します。最近、顧客から「不良率は安定しているのか」「ライン間に差があるのか」「来月の損失はいくらになり得るのか」と説明を求められました。

ここで必要なのは単一のKPIではありません。平均不良率が低くても突発異常を見逃すことがあり、Cpが高くても中心が規格中央からずれていれば不良を生みます。本記事では、水準・ばらつき・時間変化・推定誤差・経済影響を分けて観察し、最後に統合します。

2. 現場でよくある状況

月次不良率だけが報告され、いつ・どのラインで悪化したか分からない
CpとCpkが混同され、「設備能力はあるが芯ずれしている」状態を見落とす
管理限界と規格限界を同じものとして扱う
少数サンプルの差を、そのまま工程差だと断定する
異常検知の閾値や抜取検査の条件が、顧客リスクと結び付いていない
点推定だけで計画し、推定誤差や翌月の変動幅を予算に反映していない

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

品質指標には異なる役割があります。規格限界は顧客・設計が求める許容範囲、管理限界は安定工程で通常起こる変動範囲です。また、不良率は有限標本からの推定値であり、真の不良率そのものではありません。したがって「基準を超えた／超えていない」の二択だけでなく、データ生成過程、標本サイズ、時間順序、誤判定コストを明示する必要があります。

4. 今回扱うノックの全体像

No.	テーマ	主な意思決定
001	不良率推定	改善対象の優先順位
002	工程能力指数 Cp	ばらつき削減の必要性
003	工程能力指数 Cpk	芯出しと設備改善の切り分け
004	管理図	通常変動と特殊原因の識別
005	異常検知	調査対象ロットの絞り込み
006	抜取検査	検査コストと流出リスクの設計
007	信頼区間による品質評価	推定の不確実性を含む説明
008	工程の比較検定	ライン差の統計的確認
009	故障率推定	保全・予備品計画
010	品質シミュレーション	翌月損失の確率的予算化

5. Python 環境の準備

numpy と pandas でデータを扱い、統計量は scipy、可視化は matplotlib を使います。乱数シードを固定するため、実行するたびに同じ結果を再現できます。グラフ内の表記は実行環境のフォント差を避けるため英語に統一します。

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

SEED = 20260711
rng = np.random.default_rng(SEED)
pd.set_option("display.precision", 4)

print(f"Python     : {sys.version.split()[0]}")
print(f"numpy      : {np.__version__}")
print(f"pandas     : {pd.__version__}")
print(f"matplotlib : {matplotlib.__version__}")

Python     : 3.13.1
numpy      : 2.5.1
pandas     : 3.0.3
matplotlib : 3.11.0

6. 架空データの作成

120ロット（3ライン×40ロット）について、寸法測定、外観検査、設備稼働時間、故障件数を生成します。ラインごとに平均・標準偏差・基礎不良確率・故障率をわずかに変え、さらに一部ロットへ特殊原因による寸法シフトを加えます。これは分析手法の挙動を確認するための架空設定です。

変数	意味
`dimension_mm`	重要寸法の実測値（各ロット50個）
`inspected` / `defects`	外観検査数 / 不良数
`exposure_hours` / `failures`	設備稼働時間 / 故障件数
`LSL`, `USL`	寸法規格下限・上限

LSL, TARGET, USL = 49.5, 50.0, 50.5
n_lots, sample_per_lot, inspected_per_lot = 120, 50, 200
lines = np.tile(["A", "B", "C"], n_lots // 3)
dates = pd.date_range("2026-01-05", periods=n_lots, freq="D")

line_mean = {"A": 50.00, "B": 50.12, "C": 49.94}
line_sigma = {"A": 0.12, "B": 0.16, "C": 0.10}
base_defect_p = {"A": 0.007, "B": 0.015, "C": 0.009}
failure_rate = {"A": 1 / 1400, "B": 1 / 650, "C": 1 / 2200}
special_shift = {41: 0.32, 76: -0.28, 102: 0.35}  # 0始まりのロット位置

lot_rows, measurement_rows = [], []
for i, (date, line) in enumerate(zip(dates, lines)):
    mu = line_mean[line] + special_shift.get(i, 0.0)
    values = rng.normal(mu, line_sigma[line], sample_per_lot)
    p = min(base_defect_p[line] + 0.035 * abs(mu - TARGET) + (0.03 if i in special_shift else 0), 0.20)
    defects = rng.binomial(inspected_per_lot, p)
    exposure = int(rng.integers(650, 951))
    failures = rng.poisson(exposure * failure_rate[line])
    lot_id = f"L{i + 1:03d}"
    lot_rows.append((lot_id, date, line, inspected_per_lot, defects, exposure, failures,
                     values.mean(), values.std(ddof=1)))
    measurement_rows.extend((lot_id, date, line, x) for x in values)

lot_df = pd.DataFrame(lot_rows, columns=[
    "lot_id", "date", "line", "inspected", "defects", "exposure_hours",
    "failures", "dimension_mean", "dimension_std"
])
measurements = pd.DataFrame(measurement_rows, columns=["lot_id", "date", "line", "dimension_mm"])
lot_df["defect_rate"] = lot_df["defects"] / lot_df["inspected"]

print(f"ロット数: {len(lot_df):,} / 寸法測定数: {len(measurements):,}")
display(lot_df.head())

ロット数: 120 / 寸法測定数: 6,000

	lot_id	date	line	inspected	defects	exposure_hours	failures	dimension_mean	dimension_std	defect_rate
0	L001	2026-01-05	A	200	1	887	1	49.9851	0.1170	0.005
1	L002	2026-01-06	B	200	6	874	1	50.1248	0.1722	0.030
2	L003	2026-01-07	C	200	2	798	1	49.9519	0.1064	0.010
3	L004	2026-01-08	A	200	4	779	1	49.9892	0.1213	0.020
4	L005	2026-01-09	B	200	1	720	2	50.1023	0.1498	0.005

7. No.001：不良率推定

実務での意味

不良率は品質コスト、出荷判定、改善テーマ選定の入口です。ただし、ライン別不良率を単純平均すると、検査数が異なる場合に誤った値になります。基本は「不良総数 ÷ 検査総数」で集計します。

分析・モデル化の考え方

各製品が不良となる確率を $p$ 、検査数を $n$ 、不良数を $X$ とすると、独立かつ同一条件という近似の下で $X\sim\mathrm{Binomial}(n,p)$ とみなせます。最尤推定量は

$\hat{p}=\frac{X}{n}$

です。まずライン別の点推定と寄与不良数を見て、全体KPIの悪化要因を特定します。

Pythonで確認する

defect_summary = lot_df.groupby("line").agg(
    lots=("lot_id", "count"), inspected=("inspected", "sum"), defects=("defects", "sum")
)
defect_summary["defect_rate_pct"] = 100 * defect_summary["defects"] / defect_summary["inspected"]
defect_summary["defect_share_pct"] = 100 * defect_summary["defects"] / defect_summary["defects"].sum()
display(defect_summary)

ax = defect_summary["defect_rate_pct"].plot(kind="bar", color=["#4C78A8", "#F58518", "#54A24B"], legend=False)
ax.set_title("Estimated Defect Rate by Line")
ax.set_xlabel("Line")
ax.set_ylabel("Defect rate (%)")
ax.grid(axis="y", alpha=0.3)
plt.tight_layout()
plt.show()

	lots	inspected	defects	defect_rate_pct	defect_share_pct
line
A	40	8000	61	0.7625	18.8272
B	40	8000	171	2.1375	52.7778
C	40	8000	92	1.1500	28.3951

png

結果の読み取り

ラインBは不良率と不良数への寄与がともに大きく、改善の第一候補です。ただし、点推定だけでは偶然変動との区別はできません。No.007の信頼区間、No.008の検定と組み合わせて、差の確からしさを確認します。

8. No.002：工程能力指数 Cp

実務での意味

Cpは、工程が規格幅に対して十分小さいばらつきで加工できるかを示します。設備更新、治具改善、条件標準化など「ばらつきを縮める投資」の判断材料になります。

分析・モデル化の考え方

工程が安定し、寸法がおおむね正規分布に従うとき、

$C_p=\frac{USL-LSL}{6s}$

で評価します。 $6s$ は平均の周囲に広がる工程幅です。Cpは平均位置を考慮しないため、潜在能力の指標と解釈します。ここでは特殊原因3ロットを除いた基準データでライン別に計算します。

Pythonで確認する

special_lots = {"L042", "L077", "L103"}
baseline = measurements.loc[~measurements["lot_id"].isin(special_lots)].copy()
capability = baseline.groupby("line")["dimension_mm"].agg(mean="mean", std="std")
capability["Cp"] = (USL - LSL) / (6 * capability["std"])
display(capability)

ax = capability["Cp"].plot(kind="bar", color="#4C78A8", legend=False)
ax.axhline(1.33, color="#E45756", linestyle="--", label="Reference 1.33")
ax.set_title("Potential Process Capability (Cp)")
ax.set_xlabel("Line")
ax.set_ylabel("Cp")
ax.grid(axis="y", alpha=0.3)
ax.legend()
plt.tight_layout()
plt.show()

	mean	std	Cp
line
A	50.0002	0.1182	1.4100
B	50.1127	0.1612	1.0341
C	49.9389	0.0986	1.6895

png

結果の読み取り

ラインBは他ラインより標準偏差が大きく、潜在的な工程能力が低いことが分かります。一般に使われる1.33はあくまで参考線であり、採用基準は製品リスクや顧客要求に合わせて決めます。Cpだけでは中心ずれを評価できない点が次の課題です。

9. No.003：工程能力指数 Cpk

実務での意味

Cpkは、ばらつきに加えて工程平均が近い側の規格限界までどれだけ余裕を持つかを表します。Cpは良好なのにCpkが低いなら、設備能力そのものより、設定値・補正値・工具摩耗による芯ずれを優先して疑えます。

分析・モデル化の考え方

$C_{pk}=\min\left(\frac{USL-\bar{x}}{3s},\frac{\bar{x}-LSL}{3s}\right)$

CpとCpkの差は中心ずれの情報を持ちます。ただし両指標とも工程の安定性を前提とするため、管理図で特殊原因を除外してから定常期間を評価します。

Pythonで確認する

capability["Cpu"] = (USL - capability["mean"]) / (3 * capability["std"])
capability["Cpl"] = (capability["mean"] - LSL) / (3 * capability["std"])
capability["Cpk"] = capability[["Cpu", "Cpl"]].min(axis=1)
capability["Cp_minus_Cpk"] = capability["Cp"] - capability["Cpk"]
display(capability[["mean", "std", "Cp", "Cpk", "Cp_minus_Cpk"]])

ax = capability[["Cp", "Cpk"]].plot(kind="bar", color=["#4C78A8", "#F58518"])
ax.set_title("Cp and Cpk by Line")
ax.set_xlabel("Line")
ax.set_ylabel("Capability index")
ax.grid(axis="y", alpha=0.3)
plt.tight_layout()
plt.show()

	mean	std	Cp	Cpk	Cp_minus_Cpk
line
A	50.0002	0.1182	1.4100	1.4094	0.0006
B	50.1127	0.1612	1.0341	0.8010	0.2331
C	49.9389	0.0986	1.6895	1.4830	0.2065

png

結果の読み取り

ラインBではCpよりCpkがさらに低く、ばらつきと上側への中心ずれが同時に存在します。ラインCも中心が下側へ寄っているため、Cpだけを報告するとリスクを過小評価します。短期的には芯出し、恒久対策としてラインBのばらつき要因分析、という施策分解が可能です。

10. No.004：管理図

実務での意味

管理図は、不良率の時系列変化から「工程に内在する通常変動」と「調査すべき特殊原因」を分ける道具です。管理限界を超えたロットを起点に、材料、設備、作業、測定、環境の変更履歴を調べます。

分析・モデル化の考え方

検査数が一定のp管理図では、全期間のプール不良率を $\bar{p}$ として

製造業の品質管理をPythonで実践｜不良率・Cp/Cpk・管理図10本ノック

品質データを「判断」に変える：製造業の確率・統計10本ノック

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

2. 現場でよくある状況

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

4. 今回扱うノックの全体像

5. Python 環境の準備

6. 架空データの作成

7. No.001：不良率推定

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

8. No.002：工程能力指数 Cp

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

9. No.003：工程能力指数 Cpk

実務での意味

分析・モデル化の考え方

Pythonで確認する

結果の読み取り

10. No.004：管理図

実務での意味

分析・モデル化の考え方

Related Posts

製造業DIをPythonで実践｜KPI・知識グラフ・最適化・AIエージェント10本ノック

製造業シミュレーションをPythonで学ぶ｜モンテカルロからデジタルツインまで10本ノック

製造業の意思決定科学をPythonで学ぶ｜期待値・ベイズ・ロバスト最適化10本ノック

製造業の機械学習をPythonで実践｜予測・因果推論・最適化10本ノック

製造業の顧客分析をPythonで実践｜RFM・離反・CLV・Uplift 10本ノック