Executive Summary
Whether Shaft Diameter mm is capable of meeting its specification
The short answer
No—the process is marginally capable at best. Cpk is 1.24, below the 1.33 standard for acceptability, and the upper specification limit is the binding constraint at only 3.72 within-process sigma away.
The detail
Cpk = 1.241 (95% CI: 1.09 to 1.39); Ppk = 1.146. The process mean is 10.04 mm with within-process sigma 0.043 mm and overall sigma 0.0465 mm. The upper limit binds: the mean sits 0.16 mm below 10.2 mm. Zero units fell outside the specification in 150 measurements (0 parts per million observed); the normal model predicts 291.9 parts per million long-term. Normality was not rejected, validating the indices. Cp is 1.551 and Pp is 1.433; the drop reflects between-subgroup drift over the study period.
What this can't tell you
The observed defect rate of zero cannot resolve capability below roughly 6,667 parts per million at this sample size. The modelled defect estimate relies on the normality assumption; although the test did not reject it, the prediction is only as good as that assumption holds in the tails.
Analysis Overview
Capability of Shaft Diameter mm against the specification 9.8 to 10.2, over 150 units.
The short answer
Capability indices measure how much room the process has inside the specification window. Cpk of 1.24 means the shaft diameter process has about 1.24 times the spread needed to fit safely inside the 9.8–10.2 mm limits, assuming the process stays centred. The process mean sits at 10.04 mm, closer to the upper limit (10.2) than the lower one (9.8), which is why the upper side binds and reduces the index from Cp 1.55 to Cpk 1.24.
The detail
The specification window is 0.4 mm wide. The process uses 0.2579 mm of it (six within-process sigma of 0.043). Cpk is 1.241 with a 95 percent confidence interval of 1.090 to 1.391. Cp is 1.551 and Pp is 1.433; the gap between them reflects drift in the long-term record. The centring index k is 0.2, meaning the mean sits 20 percent of the way from mid-spec toward the nearer limit. The upper limit binds: CPU is 1.2408 and CPL is 1.8612, so the upper side is the constraint.
What this can't tell you
The 95 percent confidence interval on Cpk (1.09 to 1.39) is wide enough that quoting the index to two decimals overstates precision at 150 units. The process has not yet produced a defect in this sample, so the observed rate of 0 parts per million cannot distinguish the process from anything better than roughly 6,667 ppm.
Data Quality
How the rows, the limits and the subgroups were read.
The short answer
All 150 measurements loaded without gaps and were usable. The specification limits came from the data columns as supplied: lower 9.8 mm and upper 10.2 mm. The gauge produced 150 distinct values, indicating no rounding or quantization artifacts that would inflate capability artificially. No subgroup structure was available, so capability indices use the moving range method.
The detail
150 rows loaded; 150 usable units. Zero measurements were missing. The lower limit 9.8 and upper limit 10.2 were read from the 'Lower Spec' and 'Upper Spec' columns. No target value was supplied, so Cpm is not reported. No usable subgroup column was mapped, so short-term sigma comes from the moving range between consecutive measurements. 150 distinct measured values appear in the data, confirming the gauge resolution is adequate for capability analysis.
What this can't tell you
Without subgroups, we cannot separate special causes that occur within a production run from the natural variation of the process. The moving range method is valid but assumes consecutive measurements are independent; if the process exhibits autocorrelation, short-term sigma may be understated.
Distribution Against the Specification
Where Shaft Diameter mm actually falls relative to the limits it is judged against.
The short answer
The shaft diameter distribution sits inside the specification window with no observed defects, but the mean is offset toward the upper limit. The process occupies 0.2579 of the 0.4 mm specification width, and the shape matches the normal curve assumption the indices rest on.
The detail
The specification window is 0.4 mm wide (9.8 to 10.2). The process occupies 0.2579 mm (six within-process sigma of 0.043). The mean at 10.04 mm sits 20 percent of the way from mid-spec to the nearer (upper) limit. No measurement fell outside the specification in this sample of 150 units. The distribution shape is consistent with normality (skewness -0.085, excess kurtosis -0.084), supporting the indices as reported.
What this can't tell you
At 150 units, even a perfect process cannot resolve a defect rate below roughly 6,667 ppm. The absence of defects in the sample does not confirm the process is defect-free, only that any defect rate is small enough not to appear in this sample size.
Capability Indices
Every index, what it is measuring, and which one binds.
| Index | Value | Interpretation |
|---|---|---|
| Cp | 1.551 | Spec width divided by six within-process sigma — the best Shaft Diameter mm could do if it were perfectly centred |
| CPU | 1.241 | Distance from the mean up to the upper limit, in three-sigma units |
| CPL | 1.861 | Distance from the mean down to the lower limit, in three-sigma units |
| Cpk | 1.241 | The smaller of the two sides — short-term capability as actually centred (the upper side binds) |
| Pp | 1.433 | Cp recomputed on total long-term variation instead of within-subgroup variation |
| Ppk | 1.146 | Cpk recomputed on total long-term variation — what the customer actually receives over time |
| k | 0.2 | Centring index: 0 means the mean sits exactly mid-spec, 1 means it sits on a limit |
The short answer
Cpk is 1.241, below the 1.33 pass mark. The upper specification limit binds: CPU is 1.2408 while CPL is 1.8612, so the process has less headroom on the high side. The gap from Cp 1.551 to Cpk 1.241 is partly centring (the mean is offset) and partly spread.
The detail
Cpk is 1.241 with a 95 percent confidence interval of 1.090 to 1.391. Cp is 1.551 and Pp is 1.433. Cpk is 80 percent of Cp, meaning 20 percent of the shortfall against a perfectly centred process comes from centring rather than spread. The centring index k is 0.2. The upper side binds: CPU 1.2408 versus CPL 1.8612. One measurement falls outside three within-process sigma of the centre. The 1.33 pass mark is a widely used convention corresponding to the nearer limit sitting four sigma from the mean.
What this can't tell you
The confidence interval on Cpk (1.09 to 1.39) is wide at 150 units, so the two-decimal precision in the quoted index overstates certainty. Whether 1.24 is acceptable depends on your tolerance; the index describes the process as it ran, not how it will perform after changes.
Defect Rate and Sigma Level
What fell outside the limits, and what the normal model predicts.
| Statistic | Value | Basis |
|---|---|---|
| Observed out of spec (PPM) | 0 | 0 of 150 units measured outside the specification |
| Observed below lower limit (PPM) | 0 | 0 units below 9.8 |
| Observed above upper limit (PPM) | 0 | 0 units above 10.2 |
| Expected out of spec, long-term normal model (PPM) | 291.9 | normal curve at mean 10.04 with the overall standard deviation 0.0465 |
| Expected out of spec, short-term normal model (PPM) | 98.65 | normal curve at mean 10.04 with the within-process sigma 0.043 |
| Observed yield (%) | 100 | share of the 150 measurements inside the specification |
| Process sigma, short-term (Z bench) | 3.722 | the normal deviate matching the short-term out-of-spec probability |
| Sigma level with the 1.5 shift convention | 5.222 | Z bench plus 1.5, the long-term drift allowance used in Six Sigma reporting |
The short answer
The process produced zero defects in 150 units measured, but this sample size is too small to detect defect rates below roughly 6,667 parts per million. The normal model predicts 291.9 parts per million long-term, a rate that depends on the normality assumption, which was not rejected.
The detail
Observed defect rate: 0 out of 150 units, or 0 parts per million. Expected defect rate (long-term normal model): 291.9 parts per million, based on the overall standard deviation of 0.0465 mm at mean 10.04 mm. Short-term model: 98.6513 parts per million. The short-term process sigma is 3.7224 (Z bench); with the 1.5 shift convention used in Six Sigma reporting, this becomes 5.2224. The observed yield is 100%. Normality was not rejected, so the modelled rate is a reasonable extrapolation into the tails the sample never reached.
What this can't tell you
Zero observed defects cannot confirm the process is truly defect-free. The sample is too small to resolve rates below roughly 6,667 parts per million; the modelled rate is the only meaningful estimate at this size, and it is only as good as the normality assumption.
Normality Check
Whether the assumption the indices rest on actually holds.
The short answer
The data follows a normal distribution closely enough that the capability indices can be trusted. Shapiro-Wilk and Anderson-Darling tests both fail to reject normality, and the Q-Q plot shows points tight to the diagonal with no departures at either tail.
The detail
Shapiro-Wilk p = 0.8588 and Anderson-Darling p = 0.6032; neither rejects normality at the 5 percent level. Skewness is -0.085 and excess kurtosis is -0.084, both negligible. The Q-Q plot shows a tight upward pattern along the diagonal with no outliers detaching from the bulk. The normal assumption behind Cp, Cpk, Pp, Ppk, and the modelled defect rate survives here.
What this can't tell you
Normality holds in the centre and observed range, but the indices depend entirely on the tails. The test cannot rule out heavy tails beyond the sample extremes; it only confirms the sample itself is consistent with a normal distribution. A larger sample or external process knowledge would be needed to confirm tail behaviour.
Where the Variation Comes From
Short-term spread versus drift — the difference between Cp and Pp.
The short answer
The process shows two sources of variation: within-subgroup spread (0.043 mm) and between-subgroup drift (0.0178 mm), which together make up the overall standard deviation of 0.0465 mm. The drift accounts for part of the gap between short-term capability Cp 1.551 and long-term capability Pp 1.433.
The detail
Within-subgroup (short-term) standard deviation is 0.043 mm, which drives Cp and Cpk. Between-subgroup (drift) standard deviation is 0.0178 mm, the gap between Cp and Pp. Overall (long-term) standard deviation is 0.0465 mm, which drives Pp and Ppk. The drift component is 14.6 percent of the total variance. No subgroup column was mapped, so the 0.0178 mm is inferred from the residual rather than directly attributed. Eliminating the drift entirely would move Pp from 1.433 toward Cp 1.551.
What this can't tell you
Without a subgroup structure (shift, batch, machine), the drift cannot be assigned to a specific source. Mapping such a column would allow the between-subgroup variance to be directly attributed and would clarify which operational factor drives the 0.0178 mm of variation.
Methodology
Statistical methodology and diagnostics for Process Capability — Cp/Cpk
Statistical Method
Standard-library analysis: a full Cp / Cpk capability study of one measured characteristic against its specification limits. Computes the short-term indices (Cp, Cpk, CPU, CPL) from within-subgroup variation and the long-term indices (Pp, Ppk, Cpm) from total variation, estimates the defect rate in parts per million both as observed and under the normal model, reports the process sigma level, draws the distribution against the spec limits, and tests normality — because the indices assume a normal distribution and mislead badly when it fails. Works on any dataset: map the measured column and supply the limits.
- The measurements come from a stable process — capability describes a process that is already in statistical control
- The measured characteristic is approximately normally distributed (tested explicitly, and reported when it fails)
- The specification limits supplied are the real limits the output is judged against
- Subgroups, when mapped, group measurements taken close together in time or condition
- Cp, Cpk, Pp and Ppk all assume normality; on skewed or heavy-tailed data the indices and the modelled defect rate are wrong, and the report says so rather than quietly reporting them
- Capability is only meaningful for a stable process — an out-of-control process has no single capability to measure
- The modelled defect rate extrapolates far into tails that the sample never observed, so it is a model result and not a measurement
- The 1.33 pass mark is an industry convention, not a statistical threshold
Analysis Code
Complete R source code for this analysis
Process Capability — Cp / Cpk
Is the process capable of meeting its specification limits? Given a measured characteristic and the spec limits it is judged against, the analysis computes the short-term indices (Cp, Cpk) from within-subgroup variation, the long-term indices (Pp, Ppk) from total variation, the estimated defect rate in parts per million, and the process sigma level.
Why This Method?
Control charts answer "is the process stable?" — the voice of the process. Capability answers the different question the customer actually asks: "does the output fit inside the window it has to fit in?" — the voice of the customer. Cp is the ratio of the spec width to the natural process width; Cpk penalises a process that is off-centre; Pp and Ppk repeat both using long-term variation, so the gap between them measures drift between subgroups.
What This Analysis Covers
- The distribution of the characteristic drawn against its spec limits
- Cp, Cpk, CPU, CPL, Pp, Ppk, Cpm and the centring index k
- Observed and normal-model defect rates (PPM) and the process sigma level
- A normality check (Shapiro-Wilk plus a hand-implemented Anderson-Darling),
because the indices assume normality and mislead badly when it fails
- Within- versus between-subgroup variation, which is exactly the Cp/Pp gap
Standard Library
Platform standard-library module (LAT-1441): runs on ANY dataset via the semantic mapping {measurement, lower_spec, upper_spec, subgroup, target}. All narrative is derived from the user's own column names and computed values.
suppressPackageStartupMessages(library(DT))
suppressPackageStartupMessages(library(htmlwidgets))
suppressPackageStartupMessages(library(arrow))
suppressPackageStartupMessages(library(knitr))
suppressPackageStartupMessages(library(rmarkdown))
suppressPackageStartupMessages(library(dplyr))
suppressPackageStartupMessages(library(tidyr))
suppressPackageStartupMessages(library(ggplot2))
suppressPackageStartupMessages(library(stringr))
suppressPackageStartupMessages(library(lubridate))
suppressPackageStartupMessages(library(broom))
suppressPackageStartupMessages(library(Matrix))
suppressPackageStartupMessages(library(cluster))
suppressPackageStartupMessages(library(data.table))Helpers
Core Analysis Pipeline
compute_shared <- function(df, params, col_map = list()) {
# === SHARED EXPORTS ===
# initial_rows/final_rows/rows_removed $ row accounting
# measurement_name / subgroup_name $ humanized user names
# lsl / usl / target / spec_sides $ resolved specification
# lsl_source / usl_source $ where each limit came from
# mu / sd_overall / sigma_within $ location and the two spreads
# sigma_source $ how within-sigma was estimated
# cp/cpu/cpl/cpk/pp/ppu/ppl/ppk/cpm/k $ capability indices (NA when undefined)
# cpk_ci_low / cpk_ci_high $ 95% interval on Cpk
# ppm_obs_* / ppm_exp_st / ppm_exp_lt $ defect rates
# z_bench_st / z_bench_lt / sigma_level $ sigma level
# shapiro_p / ad_p / skewness / kurtosis / normal_ok / normality_text
# n_groups / sigma_between / pct_between / anova_p $ variance sources
# n_unstable / stability_text
# dist_df / qq_df / indices_df / defects_df / variation_df
# metrics / json_output
# === /SHARED EXPORTS ===
initial_rows <- nrow(df)
measurement_name <- humanize_semantic("measurement", col_map)[1]
if (!("measurement" %in% names(df))) {
stop(sprintf("A numeric measurement column('%s') must be mapped — it is the characteristic whose capability is being judged.",
measurement_name))
}Step 1: Read the measurement (95% numeric rule) and drop missing rows
mv <- coerce_numeric_95(df$measurement)
if (is.null(mv)) {
stop(sprintf("The column '%s' is not numeric — process capability needs a measured numeric characteristic.",
measurement_name))
}
keep <- is.finite(mv)
n_dropped_na <- sum(!keep)
x <- mv[keep]
n <- length(x)
if (n < 20) {
stop(sprintf("Only %d usable %s of '%s' — a capability study needs at least 20 measurements to estimate the spread at all (30 or more is the usual practical floor).",
n, units_word(n), measurement_name))
}
final_rows <- n
rows_removed <- initial_rows - final_rows
mu <- mean(x)
sd_overall <- stats::sd(x)
if (!is.finite(sd_overall) || sd_overall <= 0) {
stop(sprintf("Every value of '%s' is identical, so the process has no measurable spread and no capability index can be computed.",
measurement_name))
}Step 2: Resolve the specification limits
Three routes, in priority order: module parameters, then a mapped spec-limit column, then nothing (in which case the analysis refuses).
lsl_h <- if ("lower_spec" %in% names(df)) humanize_semantic("lower_spec", col_map)[1] else "lower spec"
usl_h <- if ("upper_spec" %in% names(df)) humanize_semantic("upper_spec", col_map)[1] else "upper spec"
tgt_h <- if ("target" %in% names(df)) humanize_semantic("target", col_map)[1] else "target"
lsl_r <- resolve_spec_limit(params$lsl %||% params$LSL %||% params$lower_spec_limit %||%
params$lower_spec %||% params$lower_specification_limit,
if ("lower_spec" %in% names(df)) df$lower_spec[keep] else NULL,
"lsl", lsl_h)
usl_r <- resolve_spec_limit(params$usl %||% params$USL %||% params$upper_spec_limit %||%
params$upper_spec %||% params$upper_specification_limit,
if ("upper_spec" %in% names(df)) df$upper_spec[keep] else NULL,
"usl", usl_h)
tgt_r <- resolve_spec_limit(params$target %||% params$nominal %||% params$target_value,
if ("target" %in% names(df)) df$target[keep] else NULL,
"target", tgt_h)
lsl <- lsl_r$value; usl <- usl_r$value; target <- tgt_r$value
has_lsl <- is.finite(lsl); has_usl <- is.finite(usl)
if (!has_lsl && !has_usl) {
stop(sprintf("No specification limit was supplied for '%s', so Cp and Cpk cannot be computed — capability is a comparison against a spec, and inventing one would invent the answer. Supply a limit in any of three ways: map a column holding the limit to lower_spec or upper_spec, pass lsl or usl as a module parameter, or add a constant limit column to the file.",
measurement_name))
}
if (has_lsl && has_usl && !(usl > lsl)) {
stop(sprintf("The upper specification limit(%s) must be greater than the lower specification limit(%s).",
fmt_num(usl, 4), fmt_num(lsl, 4)))
}
spec_sides <- if (has_lsl && has_usl) "two-sided" else if (has_usl) "upper only" else "lower only"
if (is.finite(target) && has_lsl && has_usl && (target < lsl || target > usl)) {
target <- NA_real_ # a target outside the spec window is not usable for Cpm
}Step 3: Subgroups (optional) and the two sigma estimates
Within-subgroup (short-term) sigma drives Cp/Cpk; total (long-term) sigma drives Pp/Ppk. The difference between them IS the between-subgroup drift.
has_sub <- "subgroup" %in% names(df)
subgroup_name <- if (has_sub) humanize_semantic("subgroup", col_map)[1] else NA_character_
g <- if (has_sub) as.character(df$subgroup)[keep] else NULL
if (has_sub) {
g[is.na(g) | !nzchar(trimws(g))] <- "Missing"
}
n_groups <- 0L; sigma_between <- NA_real_; pct_between <- NA_real_
anova_f <- NA_real_; anova_p <- NA_real_
group_means <- NULL; group_sizes <- NULL
sigma_within <- NA_real_; sigma_source <- ""
if (has_sub) {
sizes <- table(g)
usable <- names(sizes)[sizes >= 2]
if (length(usable) >= 2) {
idx <- g %in% usable
xg <- x[idx]; gg <- g[idx]
lev <- unique(gg)
n_groups <- length(lev)
group_sizes <- as.integer(table(gg)[lev])
group_means <- vapply(lev, function(l) mean(xg[gg == l]), numeric(1))
ss_w <- sum(vapply(lev, function(l) {
v <- xg[gg == l]; sum((v - mean(v))^2)
}, numeric(1)))
df_w <- length(xg) - n_groups
sigma_within <- sqrt(ss_w / df_w)
sigma_source <- sprintf("pooled within-subgroup standard deviation across %d subgroups of '%s'",
n_groups, subgroup_name)One-way analysis of variance, computed directly from sums of squares.
grand <- mean(xg)
ss_b <- sum(group_sizes * (group_means - grand)^2)
df_b <- n_groups - 1
ms_b <- ss_b / df_b; ms_w <- ss_w / df_w
if (is.finite(ms_w) && ms_w > 0) {
anova_f <- ms_b / ms_w
anova_p <- stats::pf(anova_f, df_b, df_w, lower.tail = FALSE)
}
n0 <- (length(xg) - sum(group_sizes^2) / length(xg)) / df_b
if (is.finite(n0) && n0 > 0) {
var_b <- max(0, (ms_b - ms_w) / n0)
sigma_between <- sqrt(var_b)
denom <- var_b + sigma_within^2
if (is.finite(denom) && denom > 0) pct_between <- 100 * var_b / denom
}
}
}
if (!is.finite(sigma_within) || sigma_within <= 0) {No usable subgroups: fall back to the Individuals estimator, the average moving range between consecutive measurements divided by d2 = 1.128.
mr <- abs(diff(x))
mrbar <- mean(mr)
sigma_within <- mrbar / 1.128
sigma_source <- "average moving range between consecutive measurements divided by 1.128 (the Individuals estimator, used because no usable subgroups were mapped)"
n_groups <- 0L
if (!is.finite(sigma_within) || sigma_within <= 0) {
sigma_within <- sd_overall
sigma_source <- "the overall standard deviation(consecutive measurements repeat exactly, so no short-term estimate was possible)"
}
}
if (!is.finite(sigma_between) || is.na(sigma_between)) {
sigma_between <- sqrt(max(0, sd_overall^2 - sigma_within^2))
if (!is.finite(pct_between)) {
denom <- sigma_between^2 + sigma_within^2
pct_between <- if (denom > 0) 100 * sigma_between^2 / denom else 0
}
}Step 4: Capability indices
cp <- if (has_lsl && has_usl) (usl - lsl) / (6 * sigma_within) else NA_real_
cpu <- if (has_usl) (usl - mu) / (3 * sigma_within) else NA_real_
cpl <- if (has_lsl) (mu - lsl) / (3 * sigma_within) else NA_real_
cpk <- min(c(cpu, cpl), na.rm = TRUE)
pp <- if (has_lsl && has_usl) (usl - lsl) / (6 * sd_overall) else NA_real_
ppu <- if (has_usl) (usl - mu) / (3 * sd_overall) else NA_real_
ppl <- if (has_lsl) (mu - lsl) / (3 * sd_overall) else NA_real_
ppk <- min(c(ppu, ppl), na.rm = TRUE)
cpm <- if (has_lsl && has_usl && is.finite(target)) {
(usl - lsl) / (6 * sqrt(sd_overall^2 + (mu - target)^2))
} else NA_real_
k <- if (has_lsl && has_usl) abs(mu - (usl + lsl) / 2) / ((usl - lsl) / 2) else NA_real_
binding_side <- if (has_lsl && has_usl) {
if (is.finite(cpu) && is.finite(cpl) && cpu <= cpl) "upper" else "lower"
} else if (has_usl) "upper" else "lower"95% interval on Cpk (Bissell's normal approximation). It widens fast as n falls, which is the honest counterweight to quoting Cpk to two decimals.
cpk_se <- sqrt(1 / (9 * n) + cpk^2 / (2 * (n - 1)))
cpk_ci_low <- cpk - stats::qnorm(0.975) * cpk_se
cpk_ci_high <- cpk + stats::qnorm(0.975) * cpk_seStep 5: Defect rates — what was observed, and what the normal model says
n_below <- if (has_lsl) sum(x < lsl) else 0L
n_above <- if (has_usl) sum(x > usl) else 0L
n_out <- n_below + n_above
ppm_obs_below <- 1e6 * n_below / n
ppm_obs_above <- 1e6 * n_above / n
ppm_obs <- 1e6 * n_out / n
tail_p <- function(sig) {
lo <- if (has_lsl) stats::pnorm((lsl - mu) / sig) else 0
hi <- if (has_usl) stats::pnorm((usl - mu) / sig, lower.tail = FALSE) else 0
c(lo = lo, hi = hi, tot = lo + hi)
}
p_st <- tail_p(sigma_within)
p_lt <- tail_p(sd_overall)
ppm_exp_st <- 1e6 * p_st[["tot"]]
ppm_exp_lt <- 1e6 * p_lt[["tot"]]
ppm_exp_lt_below <- 1e6 * p_lt[["lo"]]
ppm_exp_lt_above <- 1e6 * p_lt[["hi"]]
z_of <- function(p_tot) {
p_tot <- min(max(p_tot, 1e-16), 1 - 1e-16)
stats::qnorm(p_tot, lower.tail = FALSE)
}
z_bench_st <- z_of(p_st[["tot"]])
z_bench_lt <- z_of(p_lt[["tot"]])
sigma_level <- z_bench_st + 1.5 # the industry "six sigma" convention
yield_obs <- 100 * (n - n_out) / nStep 6: Normality — the assumption the indices stand on
set.seed(42)
x_test <- if (n > 5000) sample(x, 5000) else x
shapiro_n <- length(x_test)
shapiro_p <- NA_real_
if (shapiro_n >= 3 && length(unique(x_test)) > 1) {
sw <- tryCatch(stats::shapiro.test(x_test), error = function(e) NULL)
if (!is.null(sw)) shapiro_p <- sw$p.value
}
ad <- ad_test_normal(x)
ad_p <- ad$p
m2 <- mean((x - mu)^2); m3 <- mean((x - mu)^3); m4 <- mean((x - mu)^4)
skewness <- if (m2 > 0) m3 / m2^1.5 else NA_real_
kurtosis <- if (m2 > 0) m4 / m2^2 - 3 else NA_real_
rejects <- c(if (is.finite(shapiro_p) && shapiro_p < 0.05) "Shapiro-Wilk",
if (is.finite(ad_p) && ad_p < 0.05) "Anderson-Darling")
normal_ok <- length(rejects) == 0 && (is.finite(shapiro_p) || is.finite(ad_p))Gauge discrimination: the smallest step the instrument actually recorded. When that step covers a meaningful fraction of the process spread, every normality test is partly reading rounding rather than process shape.
n_distinct <- length(unique(x))
ux <- sort(unique(x))
gauge_res <- if (length(ux) >= 2) min(diff(ux)) else NA_real_
gauge_coarse <- (is.finite(gauge_res) && gauge_res >= sd_overall / 4) || n_distinct < 10
shapiro_txt <- if (is.finite(shapiro_p)) {
sprintf("Shapiro-Wilk p = %s%s", fmt_num(shapiro_p, 4),
if (n > 5000) sprintf(" (computed on a random sample of %s of the %s measurements)",
format(shapiro_n, big.mark = ","), format(n, big.mark = ",")) else "")
} else "Shapiro-Wilk could not be computed"
ad_txt <- if (is.finite(ad_p)) sprintf("Anderson-Darling p = %s", fmt_num(ad_p, 4))
else "Anderson-Darling could not be computed"
shape_txt <- sprintf("skewness %s, excess kurtosis %s", fmt_num(skewness, 3), fmt_num(kurtosis, 3))
normality_text <- if (normal_ok) {
paste0(shapiro_txt, " and ", ad_txt,
" — neither test rejects normality at the 5 percent level(", shape_txt,
"). The capability indices rest on a normal distribution, and that assumption survives here.")
} else if (length(rejects) > 0) {
paste0(shapiro_txt, " and ", ad_txt, " — ",
paste(rejects, collapse = " and "),
if (length(rejects) == 1) " rejects" else " reject",
" normality at the 5 percent level(", shape_txt,
"). Cp, Cpk, Pp, Ppk and the normal-model defect rate all assume a normal distribution, so on this data they should not be trusted as stated",
if (is.finite(skewness) && abs(skewness) > 0.5) {
sprintf("; the distribution is %s-skewed, which means the %s tail is longer than a normal curve allows and the defect rate on that side is understated",
if (skewness > 0) "right" else "left",
if (skewness > 0) "upper" else "lower")
} else "", ".")
} else {
paste0("Normality could not be tested on this data(", shape_txt,
"), so the capability indices should be read as descriptive only.")
}
if (gauge_coarse) {
normality_text <- paste0(normality_text, " Only ", n_distinct,
" distinct values appear across ", format(n, big.mark = ","),
" measurements, in steps of ", fmt_num(gauge_res, 6),
" against an overall standard deviation of ", fmt_num(sd_overall, 4),
", so the gauge resolution is coarse relative to the process spread and every normality test is reading rounding as much as shape.")
}Step 7: Stability — capability only means something for a stable process
if (n_groups >= 2 && !is.null(group_means)) {
grand <- mean(x)
se_g <- sigma_within / sqrt(group_sizes)
n_unstable <- sum(group_means > grand + 3 * se_g | group_means < grand - 3 * se_g)
stability_unit <- "subgroup"
} else {
n_unstable <- sum(x > mu + 3 * sigma_within | x < mu - 3 * sigma_within)
stability_unit <- "measurement"
}
stability_text <- if (n_unstable == 0) {
sprintf("No %s falls outside three within-process sigma of the centre, so nothing in this data contradicts the stability that a capability study assumes.",
stability_unit)
} else {
sprintf("%d %s%s fall outside three within-process sigma of the centre. Capability describes a stable process; where special-cause variation is present the indices describe a moving target, and a control chart should settle the process first.",
n_unstable, stability_unit, if (n_unstable == 1) "" else "s")
}Step 9: Verdict wording
capability_word <- if (cpk >= 1.33) "capable" else if (cpk >= 1.0) "marginal" else "not capable"
verdict <- sprintf("%s(Cpk %s)", tools::toTitleCase(capability_word), fmt_idx(cpk, 2))
metrics <- list()
metrics[["Measurements"]] <- n
metrics[["Process Mean"]] <- round(mu, 4)
metrics[["Sigma(Within)"]] <- round(sigma_within, 4)
metrics[["Sigma(Overall)"]] <- round(sd_overall, 4)
if (is.finite(cp)) metrics[["Cp"]] <- round(cp, 3)
metrics[["Cpk"]] <- round(cpk, 3)
if (is.finite(pp)) metrics[["Pp"]] <- round(pp, 3)
metrics[["Ppk"]] <- round(ppk, 3)
metrics[["Expected PPM"]] <- round(ppm_exp_lt, 2)
metrics[["Observed PPM"]] <- round(ppm_obs, 1)
metrics[["Process Sigma Level"]] <- round(sigma_level, 2)
metrics[["Normality"]] <- if (normal_ok) "not rejected" else "rejected"
metrics[["Verdict"]] <- verdict
spec_phrase <- if (spec_sides == "two-sided") {
sprintf("specification %s to %s", fmt_num(lsl, 4), fmt_num(usl, 4))
} else if (spec_sides == "upper only") {
sprintf("upper specification limit %s", fmt_num(usl, 4))
} else {
sprintf("lower specification limit %s", fmt_num(lsl, 4))
}
json_output <- list(
answer = paste0(
"Capability of ", measurement_name, " against the ", spec_phrase,
" across ", format(n, big.mark = ","), " ", units_word(n), ": Cpk = ",
fmt_idx(cpk, 2), " (95 percent interval ", fmt_idx(cpk_ci_low, 2), " to ",
fmt_idx(cpk_ci_high, 2), "), Ppk = ", fmt_idx(ppk, 2),
if (is.finite(cp)) paste0(", Cp = ", fmt_idx(cp, 2)) else "",
" — the process is ", capability_word,
" against the usual 1.33 threshold. The normal model puts the long-term defect rate at ",
fmt_ppm(ppm_exp_lt), " parts per million(process sigma ",
fmt_idx(sigma_level, 2), " with the 1.5 shift convention), against ",
fmt_ppm(ppm_obs), " parts per million actually observed. ",
if (normal_ok) {
"Normality is not rejected, so those model-based numbers stand."
} else {
"Normality is rejected on this data, so the indices and the model-based defect rate should not be trusted as stated."
}
),
cards = lapply(
c("tldr", "overview", "preprocessing", "capability_histogram",
"capability_indices", "defect_estimate", "normality_check",
"variation_sources"),
function(cid) list(id = cid, metrics = metrics)
)
)
list(
initial_rows = initial_rows, final_rows = final_rows, rows_removed = rows_removed,
n = n, n_dropped_na = n_dropped_na, n_distinct = n_distinct,
measurement_name = measurement_name, subgroup_name = subgroup_name,
lsl = lsl, usl = usl, target = target, has_lsl = has_lsl, has_usl = has_usl,
spec_sides = spec_sides, spec_phrase = spec_phrase,
lsl_source = lsl_r$source, usl_source = usl_r$source, target_source = tgt_r$source,
mu = mu, sd_overall = sd_overall, sigma_within = sigma_within,
sigma_source = sigma_source,
cp = cp, cpu = cpu, cpl = cpl, cpk = cpk, pp = pp, ppu = ppu, ppl = ppl,
ppk = ppk, cpm = cpm, k = k, binding_side = binding_side,
cpk_ci_low = cpk_ci_low, cpk_ci_high = cpk_ci_high,
n_below = n_below, n_above = n_above, n_out = n_out,
ppm_obs = ppm_obs, ppm_obs_below = ppm_obs_below, ppm_obs_above = ppm_obs_above,
ppm_exp_st = ppm_exp_st, ppm_exp_lt = ppm_exp_lt,
ppm_exp_lt_below = ppm_exp_lt_below, ppm_exp_lt_above = ppm_exp_lt_above,
z_bench_st = z_bench_st, z_bench_lt = z_bench_lt, sigma_level = sigma_level,
yield_obs = yield_obs,
shapiro_p = shapiro_p, ad_p = ad_p, ad_stat = ad$A2,
skewness = skewness, kurtosis = kurtosis, normal_ok = normal_ok,
normality_text = normality_text,
n_groups = n_groups, sigma_between = sigma_between, pct_between = pct_between,
anova_f = anova_f, anova_p = anova_p,
n_unstable = n_unstable, stability_text = stability_text,
gauge_res = gauge_res, gauge_coarse = gauge_coarse,
capability_word = capability_word, verdict = verdict,
dist_df = dist_df, qq_df = qq_df, indices_df = indices_df,
defects_df = defects_df, variation_df = variation_df,
metrics = metrics, json_output = json_output
)
}