u_analytics/capability/
indices.rs

1//! Process capability indices (Cp, Cpk, Pp, Ppk, Cpm).
2//!
3//! Process capability indices quantify how well a process output fits within
4//! specification limits. Short-term indices (Cp, Cpk) use within-group
5//! variation, while long-term indices (Pp, Ppk) use overall variation.
6//!
7//! # References
8//!
9//! - Montgomery (2019), *Introduction to Statistical Quality Control*, 8th ed.,
10//!   Chapter 8.
11//! - Kane (1986), "Process Capability Indices", *Journal of Quality Technology*
12//!   18(1), pp. 41--52.
13//! - Chan, Cheng & Spiring (1988), "A New Measure of Process Capability: Cpm",
14//!   *Journal of Quality Technology* 20(3), pp. 162--175.
15
16use u_numflow::stats;
17
18/// Input specification for process capability analysis.
19///
20/// Defines the upper and/or lower specification limits and an optional target
21/// value. At least one specification limit must be provided.
22///
23/// # Examples
24///
25/// ```
26/// use u_analytics::capability::ProcessCapability;
27///
28/// // Two-sided specification: LSL = 9.0, USL = 11.0
29/// let spec = ProcessCapability::new(Some(11.0), Some(9.0)).unwrap();
30///
31/// let data = [9.5, 10.0, 10.2, 9.8, 10.1, 10.3, 9.9, 10.0];
32/// let indices = spec.compute(&data, 0.15).unwrap();
33/// assert!(indices.cp.is_some());
34/// assert!(indices.cpk.is_some());
35/// ```
36pub struct ProcessCapability {
37    usl: Option<f64>,
38    lsl: Option<f64>,
39    target: Option<f64>,
40}
41
42/// Computed capability indices.
43///
44/// Fields are `Option<f64>` because not all indices can be computed for
45/// one-sided specifications. For example, Cp requires both USL and LSL.
46///
47/// # Index interpretation
48///
49/// | Index | Value | Interpretation |
50/// |-------|-------|----------------|
51/// | Cp/Pp | >= 1.33 | Process is capable |
52/// | Cpk/Ppk | >= 1.33 | Process is capable and centered |
53/// | Cpm | >= 1.33 | Process meets Taguchi loss criterion |
54///
55/// Reference: Montgomery (2019), Chapter 8, Table 8.5.
56#[derive(Debug, Clone)]
57pub struct CapabilityIndices {
58    /// Cp = (USL - LSL) / (6 * sigma_within). Requires both limits.
59    pub cp: Option<f64>,
60    /// Cpk = min(Cpu, Cpl). Requires at least one limit.
61    pub cpk: Option<f64>,
62    /// Cpu = (USL - mean) / (3 * sigma_within). Requires USL.
63    pub cpu: Option<f64>,
64    /// Cpl = (mean - LSL) / (3 * sigma_within). Requires LSL.
65    pub cpl: Option<f64>,
66    /// Pp = (USL - LSL) / (6 * sigma_overall). Requires both limits.
67    pub pp: Option<f64>,
68    /// Ppk = min(Ppu, Ppl). Requires at least one limit.
69    pub ppk: Option<f64>,
70    /// Ppu = (USL - mean) / (3 * sigma_overall). Requires USL.
71    pub ppu: Option<f64>,
72    /// Ppl = (mean - LSL) / (3 * sigma_overall). Requires LSL.
73    pub ppl: Option<f64>,
74    /// Cpm = Cp / sqrt(1 + ((mean - target) / sigma_within)^2).
75    /// Requires both limits and a target.
76    ///
77    /// Reference: Chan, Cheng & Spiring (1988).
78    pub cpm: Option<f64>,
79    /// Sample mean of the data.
80    pub mean: f64,
81    /// Short-term (within-group) standard deviation.
82    pub std_dev_within: f64,
83    /// Long-term (overall) standard deviation.
84    pub std_dev_overall: f64,
85}
86
87impl ProcessCapability {
88    /// Creates a new process capability specification.
89    ///
90    /// At least one of `usl` or `lsl` must be `Some`. If both are provided,
91    /// `usl` must be greater than `lsl`.
92    ///
93    /// # Errors
94    ///
95    /// Returns an error string if:
96    /// - Both `usl` and `lsl` are `None`
97    /// - `usl <= lsl` when both are provided
98    /// - Either limit is non-finite (NaN or infinity)
99    ///
100    /// # Examples
101    ///
102    /// ```
103    /// use u_analytics::capability::ProcessCapability;
104    ///
105    /// // Two-sided
106    /// let spec = ProcessCapability::new(Some(10.0), Some(5.0)).unwrap();
107    ///
108    /// // Upper limit only
109    /// let spec = ProcessCapability::new(Some(10.0), None).unwrap();
110    ///
111    /// // Lower limit only
112    /// let spec = ProcessCapability::new(None, Some(5.0)).unwrap();
113    ///
114    /// // Error: no limits
115    /// assert!(ProcessCapability::new(None, None).is_err());
116    ///
117    /// // Error: USL <= LSL
118    /// assert!(ProcessCapability::new(Some(5.0), Some(10.0)).is_err());
119    /// ```
120    pub fn new(usl: Option<f64>, lsl: Option<f64>) -> Result<Self, &'static str> {
121        if usl.is_none() && lsl.is_none() {
122            return Err("at least one specification limit (USL or LSL) is required");
123        }
124        if let Some(u) = usl {
125            if !u.is_finite() {
126                return Err("USL must be finite");
127            }
128        }
129        if let Some(l) = lsl {
130            if !l.is_finite() {
131                return Err("LSL must be finite");
132            }
133        }
134        if let (Some(u), Some(l)) = (usl, lsl) {
135            if u <= l {
136                return Err("USL must be greater than LSL");
137            }
138        }
139        Ok(Self {
140            usl,
141            lsl,
142            target: None,
143        })
144    }
145
146    /// Sets the target value for Cpm calculation.
147    ///
148    /// If not set, the target defaults to the midpoint `(USL + LSL) / 2`
149    /// when both limits are available. Cpm is not computed for one-sided
150    /// specifications without an explicit target.
151    ///
152    /// # Examples
153    ///
154    /// ```
155    /// use u_analytics::capability::ProcessCapability;
156    ///
157    /// let spec = ProcessCapability::new(Some(11.0), Some(9.0))
158    ///     .unwrap()
159    ///     .with_target(10.0);
160    /// ```
161    pub fn with_target(mut self, target: f64) -> Self {
162        self.target = Some(target);
163        self
164    }
165
166    /// Computes all capability indices using the provided within-group sigma.
167    ///
168    /// The `sigma_within` parameter represents the short-term (within-group)
169    /// standard deviation, typically estimated from a control chart as R-bar/d2
170    /// or S-bar/c4.
171    ///
172    /// The overall (long-term) standard deviation is computed from the data
173    /// using the sample standard deviation.
174    ///
175    /// # Returns
176    ///
177    /// `None` if:
178    /// - `data` has fewer than 2 elements
179    /// - `sigma_within` is not positive or not finite
180    /// - `data` contains NaN or infinity values
181    ///
182    /// # Examples
183    ///
184    /// ```
185    /// use u_analytics::capability::ProcessCapability;
186    ///
187    /// let spec = ProcessCapability::new(Some(11.0), Some(9.0)).unwrap();
188    /// let data = [9.5, 10.0, 10.2, 9.8, 10.1, 10.3, 9.9, 10.0];
189    /// let sigma_within = 0.15; // from control chart
190    ///
191    /// let indices = spec.compute(&data, sigma_within).unwrap();
192    /// assert!(indices.cp.unwrap() > 0.0);
193    /// assert!(indices.cpk.unwrap() > 0.0);
194    /// assert!(indices.pp.unwrap() > 0.0);
195    /// assert!(indices.ppk.unwrap() > 0.0);
196    /// ```
197    pub fn compute(&self, data: &[f64], sigma_within: f64) -> Option<CapabilityIndices> {
198        if !sigma_within.is_finite() || sigma_within <= 0.0 {
199            return None;
200        }
201        let x_bar = stats::mean(data)?;
202        let sigma_overall = stats::std_dev(data)?;
203
204        Some(self.compute_indices(x_bar, sigma_within, sigma_overall))
205    }
206
207    /// Computes capability indices using overall sigma for both short-term
208    /// and long-term estimates.
209    ///
210    /// Use this when no within-group sigma estimate is available (e.g., no
211    /// rational subgrouping). Both Cp/Cpk and Pp/Ppk will use the same
212    /// sigma, so Cp == Pp and Cpk == Ppk.
213    ///
214    /// # Returns
215    ///
216    /// `None` if:
217    /// - `data` has fewer than 2 elements
218    /// - `data` contains NaN or infinity values
219    ///
220    /// # Examples
221    ///
222    /// ```
223    /// use u_analytics::capability::ProcessCapability;
224    ///
225    /// let spec = ProcessCapability::new(Some(11.0), Some(9.0)).unwrap();
226    /// let data = [9.5, 10.0, 10.2, 9.8, 10.1, 10.3, 9.9, 10.0];
227    ///
228    /// let indices = spec.compute_overall(&data).unwrap();
229    /// // When using overall sigma for both, Cp == Pp
230    /// assert!((indices.cp.unwrap() - indices.pp.unwrap()).abs() < 1e-15);
231    /// ```
232    pub fn compute_overall(&self, data: &[f64]) -> Option<CapabilityIndices> {
233        let x_bar = stats::mean(data)?;
234        let sigma_overall = stats::std_dev(data)?;
235
236        Some(self.compute_indices(x_bar, sigma_overall, sigma_overall))
237    }
238
239    /// Internal computation of all indices given mean and sigma values.
240    fn compute_indices(
241        &self,
242        x_bar: f64,
243        sigma_within: f64,
244        sigma_overall: f64,
245    ) -> CapabilityIndices {
246        // Short-term indices (within-group sigma)
247        let cpu = self.usl.map(|u| (u - x_bar) / (3.0 * sigma_within));
248        let cpl = self.lsl.map(|l| (x_bar - l) / (3.0 * sigma_within));
249        let cp = match (self.usl, self.lsl) {
250            (Some(u), Some(l)) => Some((u - l) / (6.0 * sigma_within)),
251            _ => None,
252        };
253        let cpk = match (cpu, cpl) {
254            (Some(u), Some(l)) => Some(u.min(l)),
255            (Some(u), None) => Some(u),
256            (None, Some(l)) => Some(l),
257            (None, None) => None,
258        };
259
260        // Long-term indices (overall sigma)
261        let ppu = self.usl.map(|u| (u - x_bar) / (3.0 * sigma_overall));
262        let ppl = self.lsl.map(|l| (x_bar - l) / (3.0 * sigma_overall));
263        let pp = match (self.usl, self.lsl) {
264            (Some(u), Some(l)) => Some((u - l) / (6.0 * sigma_overall)),
265            _ => None,
266        };
267        let ppk = match (ppu, ppl) {
268            (Some(u), Some(l)) => Some(u.min(l)),
269            (Some(u), None) => Some(u),
270            (None, Some(l)) => Some(l),
271            (None, None) => None,
272        };
273
274        // Taguchi Cpm index
275        let cpm = cp.and_then(|cp_val| {
276            let target = self.target.or_else(|| match (self.usl, self.lsl) {
277                (Some(u), Some(l)) => Some((u + l) / 2.0),
278                _ => None,
279            })?;
280            let deviation_ratio = (x_bar - target) / sigma_within;
281            Some(cp_val / (1.0 + deviation_ratio * deviation_ratio).sqrt())
282        });
283
284        CapabilityIndices {
285            cp,
286            cpk,
287            cpu,
288            cpl,
289            pp,
290            ppk,
291            ppu,
292            ppl,
293            cpm,
294            mean: x_bar,
295            std_dev_within: sigma_within,
296            std_dev_overall: sigma_overall,
297        }
298    }
299}
300
301#[cfg(test)]
302mod tests {
303    use super::*;
304
305    // -----------------------------------------------------------------------
306    // Construction tests
307    // -----------------------------------------------------------------------
308
309    #[test]
310    fn new_requires_at_least_one_limit() {
311        assert!(ProcessCapability::new(None, None).is_err());
312    }
313
314    #[test]
315    fn new_rejects_usl_leq_lsl() {
316        assert!(ProcessCapability::new(Some(5.0), Some(10.0)).is_err());
317        assert!(ProcessCapability::new(Some(5.0), Some(5.0)).is_err());
318    }
319
320    #[test]
321    fn new_rejects_non_finite() {
322        assert!(ProcessCapability::new(Some(f64::NAN), Some(1.0)).is_err());
323        assert!(ProcessCapability::new(Some(10.0), Some(f64::INFINITY)).is_err());
324    }
325
326    #[test]
327    fn new_accepts_valid_two_sided() {
328        assert!(ProcessCapability::new(Some(10.0), Some(5.0)).is_ok());
329    }
330
331    #[test]
332    fn new_accepts_usl_only() {
333        assert!(ProcessCapability::new(Some(10.0), None).is_ok());
334    }
335
336    #[test]
337    fn new_accepts_lsl_only() {
338        assert!(ProcessCapability::new(None, Some(5.0)).is_ok());
339    }
340
341    // -----------------------------------------------------------------------
342    // Computation -- textbook example
343    // -----------------------------------------------------------------------
344
345    /// Textbook example: Montgomery (2019), Example 8.1
346    ///
347    /// LSL = 200, USL = 220, target = 210 (midpoint)
348    /// Process mean ~ 210, sigma_within = 2.0
349    ///
350    /// Cp = (220 - 200) / (6 * 2) = 20/12 = 1.6667
351    #[test]
352    fn textbook_centered_process() {
353        let spec = ProcessCapability::new(Some(220.0), Some(200.0)).unwrap();
354
355        let data = [
356            208.0, 209.0, 210.0, 211.0, 212.0, 208.5, 209.5, 210.5, 211.5, 210.0, 209.0, 211.0,
357            210.0, 209.5, 210.5, 210.0, 210.0, 210.0, 209.0, 211.0,
358        ];
359
360        let sigma_within = 2.0;
361        let indices = spec.compute(&data, sigma_within).unwrap();
362
363        // Cp = (220 - 200) / (6 * 2) = 1.6667
364        let cp = indices.cp.unwrap();
365        assert!(
366            (cp - 1.6667).abs() < 0.001,
367            "expected Cp ~ 1.6667, got {cp}"
368        );
369
370        let cpk = indices.cpk.unwrap();
371        assert!(cpk > 0.0, "Cpk should be positive");
372
373        let cpm = indices.cpm.unwrap();
374        assert!(cpm > 0.0, "Cpm should be positive for centered process");
375    }
376
377    /// Off-center process: mean shifted toward USL.
378    ///
379    /// LSL = 200, USL = 220, mean ~ 215, sigma_within = 2.0
380    /// Cp = (220 - 200) / (6 * 2) = 1.6667
381    /// Cpu = (220 - 215) / (3 * 2) = 0.8333
382    /// Cpl = (215 - 200) / (3 * 2) = 2.5
383    /// Cpk = min(0.8333, 2.5) = 0.8333
384    #[test]
385    fn off_center_process() {
386        let spec = ProcessCapability::new(Some(220.0), Some(200.0)).unwrap();
387
388        let data = [
389            213.0, 214.0, 215.0, 216.0, 217.0, 213.5, 214.5, 215.5, 216.5, 215.0, 214.0, 216.0,
390            215.0, 214.5, 215.5, 215.0, 215.0, 215.0, 214.0, 216.0,
391        ];
392
393        let sigma_within = 2.0;
394        let indices = spec.compute(&data, sigma_within).unwrap();
395
396        let cp = indices.cp.unwrap();
397        assert!(
398            (cp - 1.6667).abs() < 0.001,
399            "expected Cp ~ 1.6667, got {cp}"
400        );
401
402        let cpu = indices.cpu.unwrap();
403        assert!(
404            (cpu - 0.8333).abs() < 0.05,
405            "expected Cpu ~ 0.8333, got {cpu}"
406        );
407
408        let cpl = indices.cpl.unwrap();
409        assert!((cpl - 2.5).abs() < 0.05, "expected Cpl ~ 2.5, got {cpl}");
410
411        let cpk = indices.cpk.unwrap();
412        assert!((cpk - cpu).abs() < 1e-15, "Cpk should equal min(Cpu, Cpl)");
413    }
414
415    // -----------------------------------------------------------------------
416    // One-sided specifications
417    // -----------------------------------------------------------------------
418
419    #[test]
420    fn usl_only_computes_cpu_not_cpl() {
421        let spec = ProcessCapability::new(Some(10.0), None).unwrap();
422        let data = [7.0, 8.0, 9.0, 7.5, 8.5, 8.0, 7.0, 9.0, 8.0, 8.5];
423        let indices = spec.compute(&data, 0.5).unwrap();
424
425        assert!(indices.cpu.is_some());
426        assert!(indices.cpl.is_none());
427        assert!(indices.cp.is_none(), "Cp requires both limits");
428        assert!(indices.cpk.is_some(), "Cpk should be Cpu for USL-only");
429        assert!(
430            (indices.cpk.unwrap() - indices.cpu.unwrap()).abs() < 1e-15,
431            "Cpk == Cpu when only USL is set"
432        );
433    }
434
435    #[test]
436    fn lsl_only_computes_cpl_not_cpu() {
437        let spec = ProcessCapability::new(None, Some(5.0)).unwrap();
438        let data = [7.0, 8.0, 9.0, 7.5, 8.5, 8.0, 7.0, 9.0, 8.0, 8.5];
439        let indices = spec.compute(&data, 0.5).unwrap();
440
441        assert!(indices.cpl.is_some());
442        assert!(indices.cpu.is_none());
443        assert!(indices.cp.is_none());
444        assert!(indices.cpk.is_some(), "Cpk should be Cpl for LSL-only");
445        assert!(
446            (indices.cpk.unwrap() - indices.cpl.unwrap()).abs() < 1e-15,
447            "Cpk == Cpl when only LSL is set"
448        );
449    }
450
451    // -----------------------------------------------------------------------
452    // Overall-only computation
453    // -----------------------------------------------------------------------
454
455    #[test]
456    fn compute_overall_matches_pp_equals_cp() {
457        let spec = ProcessCapability::new(Some(220.0), Some(200.0)).unwrap();
458        let data = [
459            208.0, 209.0, 210.0, 211.0, 212.0, 208.5, 209.5, 210.5, 211.5, 210.0,
460        ];
461        let indices = spec.compute_overall(&data).unwrap();
462
463        let cp = indices.cp.unwrap();
464        let pp = indices.pp.unwrap();
465        assert!(
466            (cp - pp).abs() < 1e-15,
467            "Cp should equal Pp in compute_overall"
468        );
469
470        let cpk = indices.cpk.unwrap();
471        let ppk = indices.ppk.unwrap();
472        assert!(
473            (cpk - ppk).abs() < 1e-15,
474            "Cpk should equal Ppk in compute_overall"
475        );
476    }
477
478    // -----------------------------------------------------------------------
479    // Cpm with explicit target
480    // -----------------------------------------------------------------------
481
482    #[test]
483    fn cpm_with_explicit_target() {
484        let spec = ProcessCapability::new(Some(220.0), Some(200.0))
485            .unwrap()
486            .with_target(212.0);
487
488        let data = [
489            208.0, 209.0, 210.0, 211.0, 212.0, 208.5, 209.5, 210.5, 211.5, 210.0,
490        ];
491
492        let sigma_within = 2.0;
493        let indices = spec.compute(&data, sigma_within).unwrap();
494
495        let cpm = indices.cpm.unwrap();
496        let cp = indices.cp.unwrap();
497
498        assert!(
499            cpm < cp,
500            "Cpm ({cpm}) should be less than Cp ({cp}) when mean != target"
501        );
502    }
503
504    #[test]
505    fn cpm_equals_cp_when_on_target() {
506        let spec = ProcessCapability::new(Some(220.0), Some(200.0)).unwrap();
507
508        // Data perfectly at midpoint (target = 210)
509        let data = [210.0; 20];
510
511        let sigma_within = 2.0;
512        let indices = spec.compute(&data, sigma_within).unwrap();
513
514        let cpm = indices.cpm.unwrap();
515        let cp = indices.cp.unwrap();
516
517        assert!(
518            (cpm - cp).abs() < 1e-10,
519            "Cpm ({cpm}) should equal Cp ({cp}) when mean == target"
520        );
521    }
522
523    // -----------------------------------------------------------------------
524    // Edge cases
525    // -----------------------------------------------------------------------
526
527    #[test]
528    fn compute_returns_none_for_insufficient_data() {
529        let spec = ProcessCapability::new(Some(10.0), Some(0.0)).unwrap();
530        assert!(spec.compute(&[5.0], 1.0).is_none());
531        assert!(spec.compute(&[], 1.0).is_none());
532    }
533
534    #[test]
535    fn compute_returns_none_for_invalid_sigma() {
536        let spec = ProcessCapability::new(Some(10.0), Some(0.0)).unwrap();
537        let data = [4.0, 5.0, 6.0, 5.0, 5.5];
538        assert!(spec.compute(&data, 0.0).is_none());
539        assert!(spec.compute(&data, -1.0).is_none());
540        assert!(spec.compute(&data, f64::NAN).is_none());
541        assert!(spec.compute(&data, f64::INFINITY).is_none());
542    }
543
544    #[test]
545    fn compute_returns_none_for_nan_data() {
546        let spec = ProcessCapability::new(Some(10.0), Some(0.0)).unwrap();
547        let data = [4.0, f64::NAN, 6.0];
548        assert!(spec.compute(&data, 1.0).is_none());
549    }
550
551    #[test]
552    fn compute_overall_returns_none_for_insufficient_data() {
553        let spec = ProcessCapability::new(Some(10.0), Some(0.0)).unwrap();
554        assert!(spec.compute_overall(&[5.0]).is_none());
555        assert!(spec.compute_overall(&[]).is_none());
556    }
557
558    // -----------------------------------------------------------------------
559    // Pp/Ppk differ from Cp/Cpk when sigmas differ
560    // -----------------------------------------------------------------------
561
562    #[test]
563    fn pp_differs_from_cp_when_sigmas_differ() {
564        let spec = ProcessCapability::new(Some(220.0), Some(200.0)).unwrap();
565
566        let data = [
567            205.0, 207.0, 210.0, 213.0, 215.0, 206.0, 208.0, 212.0, 214.0, 210.0, 204.0, 216.0,
568            209.0, 211.0, 210.0,
569        ];
570
571        let sigma_within = 1.5;
572        let indices = spec.compute(&data, sigma_within).unwrap();
573
574        let cp = indices.cp.unwrap();
575        let pp = indices.pp.unwrap();
576
577        assert!(
578            pp < cp,
579            "Pp ({pp}) should be < Cp ({cp}) when sigma_overall > sigma_within"
580        );
581    }
582
583    // -----------------------------------------------------------------------
584    // Known numerical example
585    // -----------------------------------------------------------------------
586
587    /// Exact numerical verification.
588    ///
589    /// USL = 10, LSL = 0, sigma_within = 1.0
590    /// Cp = (10 - 0) / (6 * 1) = 1.6667
591    #[test]
592    fn exact_numerical_verification() {
593        let spec = ProcessCapability::new(Some(10.0), Some(0.0)).unwrap();
594
595        let data = [4.0, 4.5, 5.0, 5.5, 6.0, 4.0, 5.0, 6.0, 5.0, 5.0];
596        let x_bar = stats::mean(&data).unwrap();
597
598        let sigma_within = 1.0;
599        let indices = spec.compute(&data, sigma_within).unwrap();
600
601        let expected_cp = 10.0 / 6.0;
602        let expected_cpu = (10.0 - x_bar) / 3.0;
603        let expected_cpl = (x_bar - 0.0) / 3.0;
604
605        assert!(
606            (indices.cp.unwrap() - expected_cp).abs() < 1e-10,
607            "Cp: expected {expected_cp}, got {}",
608            indices.cp.unwrap()
609        );
610        assert!(
611            (indices.cpu.unwrap() - expected_cpu).abs() < 1e-10,
612            "Cpu: expected {expected_cpu}, got {}",
613            indices.cpu.unwrap()
614        );
615        assert!(
616            (indices.cpl.unwrap() - expected_cpl).abs() < 1e-10,
617            "Cpl: expected {expected_cpl}, got {}",
618            indices.cpl.unwrap()
619        );
620        assert!(
621            (indices.cpk.unwrap() - expected_cpu.min(expected_cpl)).abs() < 1e-10,
622            "Cpk: expected {}, got {}",
623            expected_cpu.min(expected_cpl),
624            indices.cpk.unwrap()
625        );
626    }
627}
u_analytics/capability/indices.rs

u_analytics/capability/
indices.rs
