Process Capability Analysis (Cp, Cpk, Pp, Ppk): Automated Capability Reporting in Python

Process capability analysis quantifies the alignment between a manufacturing process's inherent variation and its engineering specification limits. Where a control chart tells you whether a process is stable, capability indices tell you whether a stable process is good enough — whether its natural spread fits inside the tolerance band the customer specified. Within the broader SPC Fundamentals & Control Chart Taxonomy, capability metrics are the bridge between real-time monitoring and the quality numbers that end up on a PPAP submission, a supplier scorecard, or an audit report. Automating them correctly means far more than plugging four numbers into a formula: it means estimating two different sigmas without confusing them, gating the whole calculation on process stability and normality, and producing figures that reconcile to the legacy quality system to the third decimal place.

What Breaks in Production Without a Disciplined Capability Engine

The single most common capability defect in automated pipelines is reporting Cpk on an unstable process. Cp and Cpk assume the process is in statistical control — that the only variation present is common cause. If the X-Bar R chart is signalling, the within-subgroup sigma no longer represents the true process, and the capability number is a fiction that will not predict field defect rates. A pipeline that computes Cpk on live data without first confirming control ships confident-looking numbers that are quietly meaningless.

The second failure mode is confusing the two sigma families. Cp/Cpk use a within-subgroup sigma estimated from control-chart statistics; Pp/Ppk use an overall sample sigma across every reading. Mixing the estimators — for example computing "Cpk" from the total sample standard deviation — inflates or deflates the index and breaks parity with AIAG-conformant software. The gap between the two families is not noise to be smoothed away; it is a diagnostic signal, and a pipeline that collapses them throws away the most valuable information capability analysis produces.

The third trap is silent non-normality. Every textbook capability index assumes a normal distribution. Real manufacturing data — flatness, runout, concentricity, pull-force, one-sided geometric characteristics — is frequently skewed or bounded at zero. Feeding skewed data into the normal-theory formulas produces a Ppk that badly misestimates the true parts-per-million defect rate, often by an order of magnitude. A production engine must test normality, warn loudly, and route non-normal characteristics to a transformation or percentile method rather than reporting a wrong number without comment. These three failures share one outcome: capability reports that pass internally but do not survive a customer audit or predict a warranty return.

Within-subgroup versus overall capability against the tolerance band Two normal curves sit inside a tolerance band bounded by a lower and upper specification limit. A tall narrow curve, offset right of the tolerance midpoint, is the within-subgroup distribution that drives Cp and Cpk; its full width spans six within-subgroup sigma. A shorter, wider curve behind it is the overall distribution that drives Pp and Ppk, absorbing between-subgroup drift. Cpk is set by the shorter distance from the process mean to the nearer upper limit. The small area of the overall curve past the upper specification limit is shaded to mark the defect parts-per-million. LSL USL tolerance band (USL − LSL) μ 6 σ̂ within → Cp = (USL − LSL) / 6 σ̂within USL − μ (3 σ̂ ⇒ Cpk) within-subgroup short-term σ̂ → Cp / Cpk overall long-term σ̂ → Pp / Ppk defect ppm

Statistical Specification

Capability indices come in two families that differ only in which sigma goes into the denominator. The potential indices (Cp, Pp) measure spread against the tolerance width and ignore centering; the performance indices (Cpk, Ppk) penalise any offset of the process mean from the middle of the tolerance.

Within-subgroup family (Cp, Cpk) — uses $\hat{\sigma}_{within}$, the short-term common-cause sigma:

$$C_p = \frac{USL - LSL}{6\,\hat{\sigma}_{within}}$$

$$C_{pk} = \min\!\left(\frac{USL - \mu}{3\,\hat{\sigma}_{within}},\; \frac{\mu - LSL}{3\,\hat{\sigma}_{within}}\right)$$

Overall family (Pp, Ppk) — uses $\hat{\sigma}_{overall}$, the long-term sample sigma:

$$P_p = \frac{USL - LSL}{6\,\hat{\sigma}_{overall}}$$

$$P_{pk} = \min\!\left(\frac{USL - \mu}{3\,\hat{\sigma}_{overall}},\; \frac{\mu - LSL}{3\,\hat{\sigma}_{overall}}\right)$$

The two sigma estimates are calculated by fundamentally different routes:

$$\hat{\sigma}_{within} = \frac{\overline{R}}{d_2} \quad (n \le 8) \qquad\text{or}\qquad \hat{\sigma}_{within} = \frac{\overline{S}}{c_4} \quad (n \ge 9)$$

$$\hat{\sigma}_{overall} = \sqrt{\frac{1}{N-1}\sum_{i=1}^{N}(x_i - \bar{x})^2}$$

The within-subgroup sigma pools only the variation inside rational subgroups, so it captures common cause alone. The overall sigma is the ordinary sample standard deviation of every individual reading (with ddof=1), so it absorbs common cause and any drift, tool wear, or lot shift that occurred across the study. That is why $\hat{\sigma}_{overall} \ge \hat{\sigma}_{within}$ whenever between-subgroup variation exists, and therefore why Ppk ≤ Cpk on a drifting process.

The bias-correction constants must come from the same table your X-Bar R chart and X-Bar S chart use — carry them to four decimals, because rounding $d_2$ or $c_4$ propagates directly into every index.

Subgroup size $n$ $d_2$ (for $\overline{R}/d_2$) $c_4$ (for $\overline{S}/c_4$) Preferred estimator
2 1.128 0.7979 $\overline{R}/d_2$
3 1.693 0.8862 $\overline{R}/d_2$
4 2.059 0.9213 $\overline{R}/d_2$
5 2.326 0.9400 $\overline{R}/d_2$
6 2.534 0.9515 $\overline{R}/d_2$
7 2.704 0.9594 $\overline{R}/d_2$
8 2.847 0.9650 $\overline{R}/d_2$
9 2.970 0.9693 $\overline{S}/c_4$
10 3.078 0.9727 $\overline{S}/c_4$

When to Use Cpk vs Ppk vs the Chart Alternatives

Report both Cpk and Ppk whenever you can form rational subgroups, and read the gap between them as a diagnosis. When Cpk and Ppk are close, the process is stable and the within-subgroup sigma is a fair picture of long-term behaviour. When Cpk materially exceeds Ppk (a ratio above roughly 1.3), the process carries between-subgroup instability — the very drift the X-Bar R chart is designed to catch — and that instability must be resolved before the capability number means anything to a customer.

Choose the sigma estimator to match how the data was collected. For consistently sized rational subgroups of two to nine, the $\overline{R}/d_2$ route feeds Cpk directly from your X-Bar R baseline. Once subgroups regularly exceed nine the range loses efficiency, so switch to the standard-deviation route of the X-Bar S chart for large subgroups and estimate $\hat{\sigma}_{within}$ as $\overline{S}/c_4$. When you can only measure one unit at a time and no rational subgroup exists — high-mix low-volume machining, destructive tests, continuous chemical streams — there is no within-subgroup sigma to speak of; use the moving-range estimate from the Individual Moving Range (I-MR) chart, where $\hat{\sigma} = \overline{MR}/1.128$, and in that case Cpk and Ppk collapse toward the same value.

Capability indices apply to variable measurements against two-sided or one-sided tolerances. For pass/fail, defect-count, or classification data there are no Cp/Cpk indices at all — monitor those characteristics with attribute control charts (p, np, c, u) and report defect rates instead. And before any of this, the incoming series must be clean: missing values must be resolved and genuine outliers separated from real shifts, because a single dropped-sensor spike will corrupt the overall sigma and drag Ppk down for the whole run.

Production-Ready Python Implementation

The engine below computes all four indices from a subgrouped DataFrame, estimates the two sigmas by the correct route for the observed subgroup size, and refuses to return numbers when its preconditions are violated. It validates input structure, checks normality, guards against zero variation, and returns a structured dictionary ready for a dashboard or a compliance export.

import numpy as np
import pandas as pd
from scipy import stats
from typing import Dict
import warnings


class CapabilityCalculator:
    """
    Modular process-capability engine with explicit error handling
    and factory-grade data validation.

    Computes Cp, Cpk (within-subgroup sigma) and Pp, Ppk (overall sigma)
    from a subgrouped DataFrame, selecting the R-bar/d2 or S-bar/c4
    estimator automatically from the observed subgroup size.
    """

    # Bias-correction constants for unbiased within-subgroup sigma.
    # Carry four decimals — rounding propagates into every index.
    D2_LOOKUP = {
        2: 1.128, 3: 1.693, 4: 2.059, 5: 2.326,
        6: 2.534, 7: 2.704, 8: 2.847, 9: 2.970,
    }
    C4_LOOKUP = {
        2: 0.7979, 3: 0.8862, 4: 0.9213, 5: 0.9400,
        6: 0.9515, 7: 0.9594, 8: 0.9650, 9: 0.9693,
        10: 0.9727,
    }

    def __init__(
        self,
        df: pd.DataFrame,
        measurement_col: str,
        subgroup_col: str,
        usl: float,
        lsl: float,
        normality_alpha: float = 0.05,
    ) -> None:
        self.df = df.dropna(subset=[measurement_col, subgroup_col]).copy()
        self.measurement_col = measurement_col
        self.subgroup_col = subgroup_col
        self.usl = usl
        self.lsl = lsl
        self.alpha = normality_alpha

    def _validate_data(self) -> None:
        """Enforce factory-floor preconditions before any index is computed."""
        if self.df.empty:
            raise ValueError("Empty dataset after NaN removal.")
        if self.usl <= self.lsl:
            raise ValueError("USL must be strictly greater than LSL.")
        if self.df[self.subgroup_col].nunique() < 2:
            raise ValueError(
                "At least two subgroups are required for within-sigma estimation."
            )

    def _estimate_sigma_within(self) -> float:
        """R-bar/d2 for n <= 8, S-bar/c4 for n >= 9 (short-term sigma)."""
        grouped = self.df.groupby(self.subgroup_col)[self.measurement_col]
        # Use the modal subgroup size; inconsistent sizes are flagged below.
        n = int(grouped.count().mode().iloc[0])

        if n <= 8:
            r_bar = (grouped.max() - grouped.min()).mean()
            d2 = self.D2_LOOKUP.get(n)
            if d2 is None:
                raise ValueError(f"Unsupported subgroup size n={n} for R-bar/d2.")
            return r_bar / d2

        s_bar = grouped.std(ddof=1).mean()
        c4 = self.C4_LOOKUP.get(n)
        if c4 is None:
            raise ValueError(
                f"Unsupported subgroup size n={n} for S-bar/c4. "
                "Compute c4 dynamically via the gamma function for non-tabulated n."
            )
        return s_bar / c4

    def _estimate_sigma_overall(self) -> float:
        """Ordinary sample standard deviation of every reading (long-term sigma)."""
        return self.df[self.measurement_col].std(ddof=1)

    def _check_normality(self) -> bool:
        """Shapiro-Wilk normality test. Warns but does not block non-normal data."""
        series = self.df[self.measurement_col]
        n = len(series)
        if n > 5000:
            warnings.warn(
                "Shapiro-Wilk is unreliable for N > 5000; prefer Anderson-Darling."
            )
        _, p_val = stats.shapiro(series.sample(min(n, 5000), random_state=0))
        if p_val < self.alpha:
            warnings.warn(
                f"Data fails the normality test (p={p_val:.4f}). "
                "Apply a Box-Cox transform or a percentile-based capability method "
                "before trusting Ppk as a ppm estimate."
            )
            return False
        return True

    def compute(self) -> Dict[str, object]:
        self._validate_data()
        is_normal = self._check_normality()

        sigma_within = self._estimate_sigma_within()
        sigma_overall = self._estimate_sigma_overall()

        if sigma_within == 0 or sigma_overall == 0:
            raise ValueError(
                "Zero variation detected. Verify sensor resolution and data quality."
            )

        mu = self.df[self.measurement_col].mean()

        cp = (self.usl - self.lsl) / (6 * sigma_within)
        cpk = min(
            (self.usl - mu) / (3 * sigma_within),
            (mu - self.lsl) / (3 * sigma_within),
        )
        pp = (self.usl - self.lsl) / (6 * sigma_overall)
        ppk = min(
            (self.usl - mu) / (3 * sigma_overall),
            (mu - self.lsl) / (3 * sigma_overall),
        )

        return {
            "mu": round(mu, 4),
            "sigma_within": round(sigma_within, 4),
            "sigma_overall": round(sigma_overall, 4),
            "Cp": round(cp, 3),
            "Cpk": round(cpk, 3),
            "Pp": round(pp, 3),
            "Ppk": round(ppk, 3),
            "cpk_ppk_ratio": round(cpk / ppk, 3) if ppk else None,
            "is_normal": is_normal,
        }

The cpk_ppk_ratio in the returned dictionary is the diagnostic to alert on: a ratio above 1.3 flags between-subgroup instability that must be cleared before the capability figure is reported.

Validation and Testing

Capability numbers are only trustworthy when three preconditions hold, and each should be an explicit gate rather than an assumption:

  • Measurement-system capability first. Run a Gage R&R and confirm the measurement system consumes an acceptable share of the tolerance (AIAG MSA: %GRR under 10% ideal, under 30% conditional) before any capability study. A gauge that eats a third of the tolerance will masquerade as process variation and depress Cpk artificially.
  • Statistical control before capability. Confirm the process is stable on its control chart — no out-of-control signals, and limits that are not silently drifting — before computing Cp/Cpk. If baseline limits move as data arrives, lock them through a disciplined rolling-window recalibration step rather than on autopilot. An unstable process has no single "process sigma" to be capable of.
  • Normality (or an explicit alternative). The engine's Shapiro-Wilk gate warns on non-normal data; for skewed or bounded characteristics, apply a Box-Cox / Yeo-Johnson transform or a percentile-based (ISO 22514) method and document the choice.

For a numerical smoke test, generate data with a known sigma and offset and assert the recovered index. A centred process with tolerance width equal to eight sigma should return Cp ≈ 1.333, and Cpk should equal Cp when the mean sits exactly at nominal:

def test_centered_process_cp_equals_cpk():
    rng = np.random.default_rng(0)
    n_sub, size, sigma = 40, 5, 1.0
    rows = []
    for g in range(n_sub):
        for x in rng.normal(loc=10.0, scale=sigma, size=size):
            rows.append({"sg": g, "meas": x})
    df = pd.DataFrame(rows)
    # Tolerance = +/- 4 sigma around nominal 10 -> Cp target 1.333
    res = CapabilityCalculator(df, "meas", "sg", usl=14.0, lsl=6.0).compute()
    assert abs(res["Cp"] - 1.333) < 0.1
    assert abs(res["Cp"] - res["Cpk"]) < 0.1   # centered -> Cpk == Cp

Require at least 20–25 subgroups (≥ 100 individual readings is the customary PPAP floor) before freezing a capability figure; below that the sigma estimate is too noisy for the third decimal to be meaningful.

Failure Modes and Edge Cases

Symptom Root cause Fix
Cpk far exceeds Ppk Between-subgroup drift, tool wear, or lot shift inflating overall sigma Investigate the process; do not report until the control chart is stable
"Cpk" disagrees with the legacy quality system Total sample sigma used where within-subgroup sigma was required Estimate $\hat{\sigma}_{within}$ via $\overline{R}/d_2$ or $\overline{S}/c_4$, never the pooled sample std
Ppk implies far more defects than are actually seen Non-normal (skewed/bounded) characteristic run through normal-theory formula Box-Cox transform or percentile (ISO 22514) method; never report the raw index
Index jumps between software packages ddof mismatch, outlier trimming, or rounded constants Standardise on ddof=1, no silent trimming, four-decimal $d_2$/$c_4$
ValueError: Zero variation detected Sensor resolution coarser than process spread (quantised readings) Verify gauge resolution; the measurement system cannot resolve the variation
Wildly optimistic Cpk after a data merge Two process streams pooled into one column, halving apparent sigma Split streams; a subgroup must contain one homogeneous condition

Inconsistent subgroup sizes deserve special care: the engine takes the modal size for its constant lookup, so a study with mostly-five but occasional-four subgroups silently uses the $n=5$ constant. Enforce a single subgroup size in the ingestion layer, or compute a size-weighted pooled sigma, rather than letting the mode paper over the inconsistency.

Compliance Notes

Cite the four indices and their sigma estimators to the AIAG SPC Reference Manual (2nd ed.), which defines Cp/Cpk from within-subgroup variation and Pp/Ppk from overall variation and tabulates $d_2$ and $c_4$. The Gage R&R prerequisite follows AIAG MSA (4th ed.); demonstrate an acceptable measurement system before any capability study. For non-normal and non-standard distributions, follow ISO 22514-2 (process-performance and capability statistics for time-dependent process models) rather than forcing the normal-theory formulas. Under ISO 9001:2015 clause 9.1.1, retain the study dataset, the frozen control limits, the normality decision, and the engine version as documented evidence — serialize the full result dictionary alongside the raw data. The NIST/SEMATECH e-Handbook of Statistical Methods, section 6.1.6 provides the derivation of the capability indices and the ppm-to-index correspondence used to sanity-check reported defect rates.

Frequently Asked Questions

What is the practical difference between Cpk and Ppk?

Both penalise off-centre processes; they differ only in the sigma. Cpk uses the within-subgroup (short-term) sigma from control-chart statistics, so it describes what the process is capable of when only common-cause variation is present. Ppk uses the overall sample sigma across every reading, so it also absorbs drift, tool wear, and lot shifts that happened during the study. On a perfectly stable process the two are nearly equal; the gap between them measures instability.

Why is my Cpk much higher than my Ppk?

Because the process is drifting between subgroups. The within-subgroup sigma behind Cpk sees only the tight short-term spread, while the overall sigma behind Ppk also sees the between-subgroup movement, making it larger and Ppk smaller. A Cpk/Ppk ratio above roughly 1.3 means real instability — the process is not in statistical control — and you must resolve that on the control chart before the capability figure is meaningful.

Can I compute Cpk on data that fails a normality test?

Not with the standard formula. The normal-theory indices tie a Cpk value to a specific ppm defect rate, and that mapping is wrong for skewed or bounded data — often by an order of magnitude. Apply a Box-Cox or Yeo-Johnson transform and compute the index on the transformed scale, or use a percentile-based (ISO 22514) method that estimates the 0.135th and 99.865th percentiles directly. The engine warns rather than blocks so you make that decision explicitly.

Which within-subgroup sigma estimator should I use, R-bar/d2 or S-bar/c4?

Match it to the subgroup size, exactly as you would when choosing a chart. For rational subgroups of two to eight use $\overline{R}/d_2$; once subgroups reach nine or more the range loses efficiency, so switch to $\overline{S}/c_4$ from the X-Bar S chart. Whichever you pick, take $d_2$ and $c_4$ to four decimals — rounding them shifts the index enough to flip a borderline PPAP result.

What Cpk do I actually need?

The common automotive floor is 1.33 for ongoing production and 1.67 for new-process qualification, but the binding number is whatever the customer's control plan specifies — always defer to it. A Cpk of 1.33 corresponds to roughly 63 ppm defective on a centred normal process; 1.67 corresponds to roughly 0.6 ppm. If both Cpk and Ppk sit below target, the process needs either centring (reduce the mean offset) or spread reduction (reduce $\hat{\sigma}_{within}$), not simply tighter monitoring.

For chart selection criteria across every variable and attribute chart, see SPC Fundamentals & Control Chart Taxonomy.