From f7c5b35e8f651b2496338b160a2ea0925ff44c01 Mon Sep 17 00:00:00 2001 From: Chris Nicholas <4948774+cnicholas@users.noreply.github.com> Date: Tue, 29 Sep 2026 21:52:38 -0400 Subject: [PATCH] fix(spc-constants): individuals and moving-range limits use the manual's constants unrounded ## Summary - **What:** X-chart multiplier 2.66 -> 3/1.128 (2.659574) and mR upper-limit multiplier 3.268 -> 1 + 3(0.8525)/1.128 (3.267287), built from d2 = 1.128 and d3 = 0.8525. - **Why:** Bishop's VAS manual defines the limits as 3*mR/1.128 (Eq 12.10-12.11) and (1 + 3*d3/d2)*mR (Eq 12.5); 2.66 and 3.268 are its printed roundings. His VAS software computes them unrounded, so PB's X/mR limits now match VAS digit for digit (Medicare first chart: +/-3838.32 and 4715.38) and validation needs no tolerance for this. - **Scope:** spc_constants (one definition site), maximum_information (drops a duplicate 1.128), MR_SIGMA_INTERVAL_80 regenerated, tests, golden masters, docs. ## Contract / Invariants (must remain true) - [x] **Default behavior unchanged:** except X/mR limits, which move by 0.016% / 0.022% - [x] **Residuals unaffected:** no residual/SDS recomputation; charting-constant change only - [x] **Row/index alignment preserved:** golden-master rows, order and signals unchanged - [x] **Lane boundaries preserved:** unchanged - [x] **Output schema compatible:** stats dict {N, center, lpl, upl} unchanged ## Behavior Changes (explicit) - **New API:** `processbehavior.spc_constants.D3_N2 = 0.8525`; `D2_N2` is now defined directly as 1.128 (it was derived as 3/2.66 = 1.1278). - **New semantics:** XMR_LIMIT_MULTIPLIER = 3/D2_N2; R_UPPER_LIMIT_MULTIPLIER = 1 + 3*D3_N2/D2_N2. Calibrated X/mR limits and the maximum-information chart use d2 = 1.128 exactly. ## Key Design Decisions - **Constants:** follow the manual's building blocks (d2 = 1.128, d3 = 0.8525) rather than the true closed forms (2/sqrt(pi)), because the manual and VAS use 1.128. - **Single definition site:** E2 and D4 derive from D2_N2/D3_N2; no literal copies left in processbehavior/ (maximum_information's _D2 removed). - **Series-length table:** regenerated with validation/short_series_bands.py; five p10/p90 entries move by 0.001. Design-report wording untouched. ## Risks & Mitigations - **P0 risks:** golden-master regeneration hiding an unintended change. - **Mitigations:** diffed old vs new snapshots: Xbar/S values identical (files restored), no centre or signal changes, every stored X/mR limit equals its own stored mR-bar times the new constants within 3-dp rounding (and matched the old constants before). ## Tests ### Added / Updated - [x] tests/test_vas_chart_constants.py - pins Bishop's VAS 9/29 chart values (Medicare first and per-organisation charts, PM SDS 2 first chart); 7 of 8 fail on the old constants - [x] tests/test_spc_constants.py::TestMovingRangeConstants - constants equal the manual's formulas and the VAS values; XmR/R expected limits use the exact expressions - [x] tests/fixtures/golden_masters (5 X/mR scenarios) - regenerated ### Regression Focus (must fail if broken) - [x] X half-width and mR upper limit against VAS (Medicare 3838.32 / 4715.38) - [x] Xbar/S golden masters unchanged ## Manual Verification - [x] `pytest tests/` (full suite green) - [x] `pytest tests/test_vas_chart_constants.py -v` - [x] `python validation/e2e_bishop_report.py` -> 280 passed, 0 failed - [x] `ruff check .`; mypy unchanged from main (17 pre-existing advisory errors) ## Notes - **Docs:** chart-types, wheeler-terminology, series-length table, two examples' printed outputs; CHANGELOG [Unreleased]. Stored notebook outputs refresh on the next docs rebuild. - **Follow-ups (not in this PR):** Release v0.3.3 commit; bump processbehavior-app pin. --- CHANGELOG.md | 12 +++ docs/appendix/wheeler-terminology.md | 4 +- docs/examples/process-behavior-chart.md | 2 +- docs/examples/xmr-chart-in-python.md | 4 +- docs/user-guide/chart-types.md | 4 +- docs/user-guide/series-length.md | 22 ++--- processbehavior/maximum_information.py | 11 +-- processbehavior/spc_constants.py | 41 +++++---- .../missing_values/X_data.parquet | Bin 6209 -> 6209 bytes .../missing_values/mR_data.parquet | Bin 6894 -> 6894 bytes .../missing_values/mR_statistics.json | 2 +- .../paired_xmr_r/X_data.parquet | Bin 7057 -> 7057 bytes .../paired_xmr_r/X_statistics.json | 2 +- .../paired_xmr_r/mR_data.parquet | Bin 7941 -> 7941 bytes .../pathological_ordering/X_data.parquet | Bin 6030 -> 6030 bytes .../pathological_ordering/X_statistics.json | 2 +- .../pathological_ordering/mR_data.parquet | Bin 6646 -> 6646 bytes .../pathological_ordering/mR_statistics.json | 2 +- .../single_obs_strata/X_data.parquet | Bin 6117 -> 6117 bytes .../single_obs_strata/X_statistics.json | 4 +- .../single_obs_strata/mR_data.parquet | Bin 6757 -> 6757 bytes .../single_obs_strata/mR_statistics.json | 2 +- .../unstratified_small/X_data.parquet | Bin 5774 -> 5774 bytes .../unstratified_small/X_statistics.json | 2 +- .../unstratified_small/mR_data.parquet | Bin 6334 -> 6334 bytes .../unstratified_small/mR_statistics.json | 2 +- tests/test_spc_constants.py | 40 +++++++- tests/test_vas_chart_constants.py | 86 ++++++++++++++++++ 28 files changed, 186 insertions(+), 58 deletions(-) create mode 100644 tests/test_vas_chart_constants.py diff --git a/CHANGELOG.md b/CHANGELOG.md index 7f95e29..266f386 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -7,6 +7,18 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0 ## [Unreleased] +### Changed +- **Individuals and moving-range limits now use the constants exactly as Bishop's VAS manual + defines them.** The individuals-chart multiplier is 3/1.128 = 2.659574 (previously the + rounded 2.66) and the moving-range upper limit is 1 + 3(0.8525)/1.128 = 3.267287 (previously + 3.268), built from the manual's d₂ = 1.128 and d₃ = 0.8525 (Eq 12.4, 12.5, 12.10, 12.11). The + manual prints the rounded values; Bishop's VAS software computes them without rounding, and PB's + X and mR limits now match it digit for digit (Medicare first chart: ±3838.32 and 4715.38). + Limits move by 0.016% and 0.02%, so a signal changes only for a point within that distance of a + limit. The d₂ used on the calibration path and by the maximum-information chart is now exactly + 1.128 (it was 3/2.66 = 1.1278), and `MR_SIGMA_INTERVAL_80` was regenerated: five entries move + by 0.001. New constant `D3_N2` in `processbehavior.spc_constants`. + ### Added - Docs: the chart-types guide now explains the three views of a study with factors and time (the combined ``by=[]`` chart, the ``phased=True`` view, and full stratification), diff --git a/docs/appendix/wheeler-terminology.md b/docs/appendix/wheeler-terminology.md index 87848ee..a932260 100644 --- a/docs/appendix/wheeler-terminology.md +++ b/docs/appendix/wheeler-terminology.md @@ -233,8 +233,8 @@ Wheeler uses the standard SPC constants from Shewhart's work: | c₄ | Unbiasing s | Related to gamma function | | A₃ | Xbar limits from s | 3 / (c₄√n) | | B₃, B₄ | S chart limits | Functions of c₄ and n | -| d₂ | Unbiasing range | Tabulated | -| D₃, D₄ | mR chart limits | Functions of d₂ | +| d₂ | Unbiasing range | Tabulated; 1.128 for moving ranges of two | +| D₃, D₄ | mR chart limits | Functions of d₂ and d₃; D₄ = 1 + 3(0.8525)/1.128 for moving ranges of two | --- diff --git a/docs/examples/process-behavior-chart.md b/docs/examples/process-behavior-chart.md index 513ceea..0470e5a 100644 --- a/docs/examples/process-behavior-chart.md +++ b/docs/examples/process-behavior-chart.md @@ -36,7 +36,7 @@ Output: ``` Center line: 147.38 -Natural process limits: (139.07, 155.69) +Natural process limits: (139.071, 155.689) Signals: 26 ``` diff --git a/docs/examples/xmr-chart-in-python.md b/docs/examples/xmr-chart-in-python.md index 142f80d..8b99c02 100644 --- a/docs/examples/xmr-chart-in-python.md +++ b/docs/examples/xmr-chart-in-python.md @@ -29,8 +29,8 @@ result.plot() # interactive plotly figure: X on top, mR below Output: ``` -{'N': 1, 'center': 91.883, 'lpl': 83.209, 'upl': 100.557} -{'N': 2, 'center': 3.261, 'lpl': 0.0, 'upl': 10.657} +{'N': 1, 'center': 91.883, 'lpl': 83.211, 'upl': 100.556} +{'N': 2, 'center': 3.261, 'lpl': 0.0, 'upl': 10.654} ``` ## Reading it diff --git a/docs/user-guide/chart-types.md b/docs/user-guide/chart-types.md index 2850a3f..3c350e4 100644 --- a/docs/user-guide/chart-types.md +++ b/docs/user-guide/chart-types.md @@ -242,7 +242,7 @@ Plots the **standard deviation** of each subgroup. Plots each individual observation. - **Centerline**: Average of all observations (X̅) -- **Control Limits**: X̅ ± 2.66 × R̅ (average moving range) +- **Control Limits**: X̅ ± (3/1.128) × R̅, about X̅ ± 2.66 × R̅ (R̅ is the average moving range) - **Interpretation**: Points beyond limits indicate special causes ### The mR (Moving Range) Chart @@ -250,7 +250,7 @@ Plots each individual observation. Plots the absolute difference between consecutive observations. - **Centerline**: Average moving range (R̅) -- **UCL**: 3.27 × R̅ +- **UCL**: (1 + 3 × 0.8525/1.128) × R̅, about 3.267 × R̅ - **LCL**: 0 (range cannot be negative) - **Interpretation**: Large ranges indicate sudden changes diff --git a/docs/user-guide/series-length.md b/docs/user-guide/series-length.md index 20e9db0..e5fe207 100644 --- a/docs/user-guide/series-length.md +++ b/docs/user-guide/series-length.md @@ -72,30 +72,30 @@ anyone who wants them in code; the design report prints only the sentence. | 4 | 0.430 | 0.928 | 1.676 | 3.90 | | 5 | 0.488 | 0.948 | 1.593 | 3.27 | | 6 | 0.528 | 0.951 | 1.526 | 2.89 | -| 7 | 0.565 | 0.962 | 1.484 | 2.63 | -| 8 | 0.595 | 0.971 | 1.453 | 2.44 | +| 7 | 0.565 | 0.961 | 1.484 | 2.63 | +| 8 | 0.595 | 0.970 | 1.453 | 2.44 | | 9 | 0.618 | 0.974 | 1.418 | 2.30 | | 10 | 0.637 | 0.974 | 1.396 | 2.19 | | 11 | 0.656 | 0.977 | 1.375 | 2.10 | -| 12 | 0.670 | 0.981 | 1.360 | 2.03 | +| 12 | 0.670 | 0.981 | 1.359 | 2.03 | | 13 | 0.680 | 0.983 | 1.344 | 1.98 | -| 14 | 0.694 | 0.984 | 1.330 | 1.92 | +| 14 | 0.694 | 0.983 | 1.329 | 1.92 | | 15 | 0.703 | 0.983 | 1.315 | 1.87 | | 16 | 0.714 | 0.987 | 1.310 | 1.83 | -| 17 | 0.722 | 0.986 | 1.300 | 1.80 | +| 17 | 0.721 | 0.986 | 1.300 | 1.80 | | 18 | 0.730 | 0.986 | 1.288 | 1.77 | | 19 | 0.739 | 0.987 | 1.280 | 1.73 | | 20 | 0.745 | 0.988 | 1.274 | 1.71 | | 21 | 0.749 | 0.989 | 1.266 | 1.69 | | 22 | 0.756 | 0.989 | 1.258 | 1.67 | -| 23 | 0.761 | 0.989 | 1.253 | 1.65 | +| 23 | 0.761 | 0.988 | 1.253 | 1.65 | | 24 | 0.766 | 0.990 | 1.249 | 1.63 | | 25 | 0.771 | 0.991 | 1.243 | 1.61 | -| 26 | 0.775 | 0.991 | 1.238 | 1.60 | -| 27 | 0.780 | 0.992 | 1.234 | 1.58 | -| 28 | 0.784 | 0.993 | 1.229 | 1.57 | -| 29 | 0.787 | 0.992 | 1.224 | 1.55 | -| 30 | 0.790 | 0.993 | 1.220 | 1.54 | +| 26 | 0.774 | 0.991 | 1.238 | 1.60 | +| 27 | 0.780 | 0.991 | 1.234 | 1.58 | +| 28 | 0.784 | 0.992 | 1.229 | 1.57 | +| 29 | 0.787 | 0.991 | 1.223 | 1.55 | +| 30 | 0.790 | 0.992 | 1.220 | 1.54 | Read a row as: at four points, the middle 80% of moving-range sigma estimates fall between 0.43 and 1.68 times the true sigma. Most of the improvement is spent by eight to ten diff --git a/processbehavior/maximum_information.py b/processbehavior/maximum_information.py index 3ff25b5..b817a68 100644 --- a/processbehavior/maximum_information.py +++ b/processbehavior/maximum_information.py @@ -21,7 +21,7 @@ import pandas as pd from .exceptions import ValidationError -from .spc_constants import R_UPPER_LIMIT_MULTIPLIER, XMR_LIMIT_MULTIPLIER +from .spc_constants import D2_N2, R_UPPER_LIMIT_MULTIPLIER, XMR_LIMIT_MULTIPLIER if TYPE_CHECKING: from .analysis_dataset import AnalysisDataSet @@ -48,9 +48,9 @@ class MaximumInformationResult: sigma_hat : float mR / d2 — noise floor sigma estimate. upl : float - Upper natural process limit (R2 mean + 2.66 * mR). + Upper natural process limit (R2 mean + E2 * mR, E2 = 3/1.128 ≈ 2.66). lpl : float - Lower natural process limit (R2 mean - 2.66 * mR). + Lower natural process limit (R2 mean - E2 * mR). n_signals : int Points beyond limits on XmR. round_to : int @@ -159,9 +159,6 @@ def __repr__(self) -> str: # Pure Function # ============================================================================ -# d2 for n=2 (moving range of consecutive pairs) -_D2 = 1.128 - def assess_maximum_information( ads: AnalysisDataSet, @@ -204,7 +201,7 @@ def assess_maximum_information( r2_mean = float(np.mean(r2_values)) mr_values = np.abs(np.diff(r2_values)) r2_mR = float(np.mean(mr_values)) - sigma_hat = r2_mR / _D2 + sigma_hat = r2_mR / D2_N2 # Limits: mean ± E2 * mR upl = r2_mean + XMR_LIMIT_MULTIPLIER * r2_mR diff --git a/processbehavior/spc_constants.py b/processbehavior/spc_constants.py index 9e9d38f..e96a4ab 100644 --- a/processbehavior/spc_constants.py +++ b/processbehavior/spc_constants.py @@ -29,15 +29,19 @@ # Control limit multiplier (3-sigma limits are standard in SPC) SIGMA_MULTIPLIER = 3 -# E2 constant for XmR charts (n=2, moving range of 2 consecutive observations) -# Used for calculating control limits on individual values -# E2 = d2 / d3 for n=2, where d2 = 1.128 and d3 = 0.8525 -XMR_LIMIT_MULTIPLIER = 2.66 +# Moving-range constants for n = 2 (consecutive pairs), as in Bishop's VAS manual: +# sigma = mR / d2 (Eq 12.4), X limits = X ± 3·mR/d2 (Eq 12.10-12.11), and the +# moving-range upper limit = (1 + 3·d3/d2)·mR (Eq 12.5). The manual prints the +# rounded values 2.66 and 3.268; Bishop's VAS software computes them from d2 and +# d3 without rounding, and so does this library. +D2_N2 = 1.128 +D3_N2 = 0.8525 -# D4 constant for R charts (n=2, range of 2 consecutive observations) -# Used for upper control limit on moving range -# D4 = 1 + 3(d3/d2) for n=2 -R_UPPER_LIMIT_MULTIPLIER = 3.268 +# E2 for XmR charts: X ± E2·mR, E2 = 3/d2 (≈ 2.66) +XMR_LIMIT_MULTIPLIER = SIGMA_MULTIPLIER / D2_N2 + +# D4 for the moving-range chart: upper limit = D4·mR, D4 = 1 + 3·d3/d2 (≈ 3.267) +R_UPPER_LIMIT_MULTIPLIER = 1 + SIGMA_MULTIPLIER * D3_N2 / D2_N2 # ============================================================================ @@ -260,8 +264,8 @@ def calculate_limits( ----- Xbar limits: X̄ ± (3 * Wd) / sqrt(n), where Wd = S / c4(n) S limits: S * b3(n) to S * b4(n) - XmR limits: X̄ ± (E2 * mR), where E2 = 2.66 - R limits: 0 to mR * D4, where D4 = 3.268 + XmR limits: X̄ ± (E2 * mR), where E2 = 3/d2 = 3/1.128 (≈ 2.66) + R limits: 0 to mR * D4, where D4 = 1 + 3·d3/d2 = 1 + 3(0.8525)/1.128 (≈ 3.267) References ---------- @@ -414,11 +418,8 @@ def _per_n(func, sizes, *args): return pd.DataFrame({'lpl': lpl, 'upl': upl}, index=index) -# d2 bias constant for the n=2 moving range, derived from the library's own E2 -# (E2 = sigma_multiplier / d2 at n=2, with the default 3-sigma multiplier). Used -# only on the calibration path so calibrated X/mR limits stay internally -# consistent with the data-path XmR/R constants. -D2_N2 = 3.0 / XMR_LIMIT_MULTIPLIER +# calibrated_limits uses D2_N2 (defined with the XmR constants above), so +# calibrated X/mR limits stay consistent with the data-path XmR/R constants. def calibrated_limits( @@ -714,12 +715,12 @@ def suggest_chart_name(name: str) -> str: 9: (0.618, 1.418), 10: (0.637, 1.396), 11: (0.656, 1.375), - 12: (0.670, 1.360), + 12: (0.670, 1.359), 13: (0.680, 1.344), - 14: (0.694, 1.330), + 14: (0.694, 1.329), 15: (0.703, 1.315), 16: (0.714, 1.310), - 17: (0.722, 1.300), + 17: (0.721, 1.300), 18: (0.730, 1.288), 19: (0.739, 1.280), 20: (0.745, 1.274), @@ -728,9 +729,9 @@ def suggest_chart_name(name: str) -> str: 23: (0.761, 1.253), 24: (0.766, 1.249), 25: (0.771, 1.243), - 26: (0.775, 1.238), + 26: (0.774, 1.238), 27: (0.780, 1.234), 28: (0.784, 1.229), - 29: (0.787, 1.224), + 29: (0.787, 1.223), 30: (0.790, 1.220), } diff --git a/tests/fixtures/golden_masters/missing_values/X_data.parquet b/tests/fixtures/golden_masters/missing_values/X_data.parquet index a1a9ee9a99d5563883a80cfa460d7b54ba1cacf7..d69d769ef7035b731a8629248020b9df124c7b6d 100644 GIT binary patch delta 40 rcmX?TaL{0b6(6gqo`Ig>VtYPzcE4bga-S^M&2FNXnIOWG8{IAAgemKcKuhUn&l%xhRR zh*QhtF!=+!$mDhGyqmXi=CQGw>KW)6Zq^itXJGXMYp delta 171 zcmaE7`p$GiBeTQP3re?t)-gGVGKdO^is^_(hzf|ZiL%K^aKK~`EHMTR4AIR8nb)vt z5T};OVe$udk;&`Wc{gw4%wuCU(lgLA*sLiK&(2=yljT_I>%6&N+=Q7OEG{``@&zf; I$@iuB0aQysFaQ7m diff --git a/tests/fixtures/golden_masters/missing_values/mR_statistics.json b/tests/fixtures/golden_masters/missing_values/mR_statistics.json index 443a3f7..83c4b0c 100644 --- a/tests/fixtures/golden_masters/missing_values/mR_statistics.json +++ b/tests/fixtures/golden_masters/missing_values/mR_statistics.json @@ -2,5 +2,5 @@ "N": 2, "center": 0.707, "lpl": 0.0, - "upl": 2.312 + "upl": 2.311 } \ No newline at end of file diff --git a/tests/fixtures/golden_masters/paired_xmr_r/X_data.parquet b/tests/fixtures/golden_masters/paired_xmr_r/X_data.parquet index 98a70ca72f24a76a6058e1c66283e7e776272560..49699c884d6f89a09d028eb9579c94e325173cff 100644 GIT binary patch delta 165 zcmbPeKGA%G5v!B+egEF+aqbSH45DJ9QaYj+7)1p{*+ki7BsgGF2$mRw2Bzp_epcnp x(^+d+)rr+`flFs|8n*-+tErxWp5f+ZK_+&#&{AJ#lZwsV68D%O!jto)l>uk7GP?i( delta 165 zcmbPeKGA%G5vx;#<(3m=aqbSH45DJ9QaYj+7)1p{*+ki7BsgGF2$mRw2Bzp_epcnp x(^+d+)rr+`flFs|8n*-+tC608p26m3K_+%K&r)A!lk&~o68D%O!jto)l>tE*GE4vf diff --git a/tests/fixtures/golden_masters/paired_xmr_r/X_statistics.json b/tests/fixtures/golden_masters/paired_xmr_r/X_statistics.json index f4e047a..54a1cd4 100644 --- a/tests/fixtures/golden_masters/paired_xmr_r/X_statistics.json +++ b/tests/fixtures/golden_masters/paired_xmr_r/X_statistics.json @@ -1,6 +1,6 @@ { "N": 1, "center": 49.085, - "lpl": 46.738, + "lpl": 46.739, "upl": 51.432 } \ No newline at end of file diff --git a/tests/fixtures/golden_masters/paired_xmr_r/mR_data.parquet b/tests/fixtures/golden_masters/paired_xmr_r/mR_data.parquet index 843fe255d6cdf4c7bba5ac1603240bbeba4f001c..35f445acb930ff480f79be0faec696af19b100e8 100644 GIT binary patch delta 41 scmZp*Yqi^8Ey8N5XP{@e*+C?po!!tU%dylqc=Ix86K06`8{IAAgemKcKuy69$pW-%5G zqLczHnJmaAGWk3!@8&#qNj6qfJp(<%#XY?2?0&%}l>S|nW9|;345C7!VmhJ@q5`69qHHn}955LKON>DSU34=)vlxp8 zQA&Z9OcrDlnS7p=cXJ-QBpa)do`Igh;vQahcHc~sa-S^6&67ngGqHn(C3_|ph>P;b QFfcGo6k=cqa11g80Lz~-RR910 diff --git a/tests/fixtures/golden_masters/pathological_ordering/X_statistics.json b/tests/fixtures/golden_masters/pathological_ordering/X_statistics.json index 870816f..ecfabb8 100644 --- a/tests/fixtures/golden_masters/pathological_ordering/X_statistics.json +++ b/tests/fixtures/golden_masters/pathological_ordering/X_statistics.json @@ -1,6 +1,6 @@ { "N": 1, "center": 49.774, - "lpl": 47.548, + "lpl": 47.549, "upl": 52.0 } \ No newline at end of file diff --git a/tests/fixtures/golden_masters/pathological_ordering/mR_data.parquet b/tests/fixtures/golden_masters/pathological_ordering/mR_data.parquet index 961fec46315e2d248d60ecdc9061033bf0e514dd..bd245ca42493df2d6e3a8c275b24818198dda630 100644 GIT binary patch delta 165 zcmexn{LOfS7_-9~nWz-&`>YP445C7!VmhLBq5`69qHHn}955LKON>DSLv(W(vplN? zacWr|CV$|Rovg#jv-ucD8XK#ro`Ig>W<&mXc6LLbEXPvc;LX#;Oqe0!lP^fB0RY=M BE;9fC delta 165 zcmexn{LOfS7_-B}Dv@>W4_O^V8AOFd#dJjNL8{IAAgemKcKuhUn%nW_eZ( z;?%M_O#Z+rJ6VU5XY(8{IAAgemKcKuy69v-W@YR3@VOAy1x`ph{>S)@57Gp6q!E5pbHl@k)IXE_( la;#)!HPtiFGhFuLfTBQI0 delta 272 zcmaE=|5SeiKeK&QiuHZ}xo!@k45C7!VmhK8q5`69qHHn}955LKON>DSU39V^v$B2f z^gN;F1TQRRz+@yiU@VXsFg8?l^9JTK%&Np$w~*ag9ZfsXVl0Lxcul^*rZjmz2ghbp lj+LyeMtTN%28$o_va|bUnw0xwIc|O>a+wJtJo%8gG5{iVQMUj9 diff --git a/tests/fixtures/golden_masters/single_obs_strata/X_statistics.json b/tests/fixtures/golden_masters/single_obs_strata/X_statistics.json index 83381a4..fea9e49 100644 --- a/tests/fixtures/golden_masters/single_obs_strata/X_statistics.json +++ b/tests/fixtures/golden_masters/single_obs_strata/X_statistics.json @@ -1,6 +1,6 @@ { "N": 1, "center": 48.991, - "lpl": 45.229, - "upl": 52.754 + "lpl": 45.23, + "upl": 52.753 } \ No newline at end of file diff --git a/tests/fixtures/golden_masters/single_obs_strata/mR_data.parquet b/tests/fixtures/golden_masters/single_obs_strata/mR_data.parquet index 5667b7867b8b95bd01052d0086372d5779fc3f39..216f9e478f61891c6235813ebe06a9547b46907d 100644 GIT binary patch delta 165 zcmaEA^3-I5FSElm=kA>6)j|%U45C7!VmhL3q5`69qHHn}955LKON>DSLv(X5b2zI8 zacYGeCJVAFPM*)nvDuU}jg8e*&p^*`^8)^Oc6LLbEXPvc;LU8}Cd?4=$@Wre00=`e AoB#j- delta 165 zcmaEA^3-I5FSEmjqKt`eYlIv`8AOFd#dJj7L8{IAAgemKcKuhUn&A=5SUG z;?xQ`OcrEUoIIbCW3wq|8XK#To`Igh<^}xm?Ch03S&pT?&YRi9O_(9#lkKI{04ZNJ At^fc4 diff --git a/tests/fixtures/golden_masters/single_obs_strata/mR_statistics.json b/tests/fixtures/golden_masters/single_obs_strata/mR_statistics.json index c253ce4..47c6333 100644 --- a/tests/fixtures/golden_masters/single_obs_strata/mR_statistics.json +++ b/tests/fixtures/golden_masters/single_obs_strata/mR_statistics.json @@ -2,5 +2,5 @@ "N": 2, "center": 1.414, "lpl": 0.0, - "upl": 4.622 + "upl": 4.621 } \ No newline at end of file diff --git a/tests/fixtures/golden_masters/unstratified_small/X_data.parquet b/tests/fixtures/golden_masters/unstratified_small/X_data.parquet index 45acbde80bc525c5bec02c9458ac1c4a75141725..ceed21874c1b11dd61df6fb703ad0f7e075dc4d8 100644 GIT binary patch delta 167 zcmeCv?bF@hz+^Aq`Z{ojp1Xr6gQ$?Gn2xBNsDLP&D4UD~2TTUR5@XOn7u{UXw2WDm zIHg`}wrs4XdIoxii+gz3+5Lh|%6+n2H%}J6%)|~Bmh72aASTKq!@$7MF2ukP;22~G E03tFk`Tzg` delta 167 zcmeCv?bF@hz+`WI-@kXdp1Xr6gQ$?Gn2xBNsDLP&D4UD~2TTUR5@XOn7u{UXw2WDm zIHg`}wrs3MdIow1i+gz3*?luj%6+mNH%}J6%)|~Bmh72aASTKq!@$7MF2ukP;22~G E00Igxpa1{> diff --git a/tests/fixtures/golden_masters/unstratified_small/X_statistics.json b/tests/fixtures/golden_masters/unstratified_small/X_statistics.json index c4af14e..db6493d 100644 --- a/tests/fixtures/golden_masters/unstratified_small/X_statistics.json +++ b/tests/fixtures/golden_masters/unstratified_small/X_statistics.json @@ -1,6 +1,6 @@ { "N": 1, "center": 49.213, - "lpl": 46.364, + "lpl": 46.365, "upl": 52.061 } \ No newline at end of file diff --git a/tests/fixtures/golden_masters/unstratified_small/mR_data.parquet b/tests/fixtures/golden_masters/unstratified_small/mR_data.parquet index b34898bb51f44e86ffba42258f2edd2822b801fd..d8303b3a8819cc9b15b388ef4e077a8e78c43ca6 100644 GIT binary patch delta 172 zcmdmIxX*CIbSB5s3thK;{>$wk${;EvDyAbUBPt-uCdwuw!2y#&u*4WNFhw`tU`l1t zAX2fzdB9&Tim~RBrIA3C delta 153 zcmdmIxX*CIbS4f4FyL{R%*ZUE$m1XZ7Z796fV1JEn{P0svS`5cz=a*KDt4HBkX2-| y9UJfFHugL=RwF$FJ%i0ZdE?pHD}AyYOMRU;$BUXUvxCJY=S*H8Av$@z1U~>1M;Q75 diff --git a/tests/fixtures/golden_masters/unstratified_small/mR_statistics.json b/tests/fixtures/golden_masters/unstratified_small/mR_statistics.json index 361d11e..58db956 100644 --- a/tests/fixtures/golden_masters/unstratified_small/mR_statistics.json +++ b/tests/fixtures/golden_masters/unstratified_small/mR_statistics.json @@ -2,5 +2,5 @@ "N": 2, "center": 1.071, "lpl": 0.0, - "upl": 3.5 + "upl": 3.499 } \ No newline at end of file diff --git a/tests/test_spc_constants.py b/tests/test_spc_constants.py index fdd3a25..79e888c 100644 --- a/tests/test_spc_constants.py +++ b/tests/test_spc_constants.py @@ -12,6 +12,10 @@ import pytest from processbehavior.spc_constants import ( + D2_N2, + D3_N2, + R_UPPER_LIMIT_MULTIPLIER, + XMR_LIMIT_MULTIPLIER, b3, b4, c4, @@ -121,6 +125,34 @@ def test_b3_b4_raises_on_invalid_n(func): func(1) +# ============================================================================ +# Test: moving-range constants (n = 2) +# ============================================================================ + + +class TestMovingRangeConstants: + """The XmR and mR constants follow Bishop's VAS manual (Eq 12.4, 12.5, 12.10, 12.11) + and are computed from d2 = 1.128 and d3 = 0.8525 without rounding, as VAS does.""" + + def test_building_blocks_are_the_manuals(self): + assert D2_N2 == 1.128 + assert D3_N2 == 0.8525 + + def test_multipliers_are_the_manuals_formulas(self): + assert XMR_LIMIT_MULTIPLIER == 3 / D2_N2 + assert R_UPPER_LIMIT_MULTIPLIER == 1 + 3 * D3_N2 / D2_N2 + + def test_multipliers_match_the_vas_software(self): + # Measured from Bishop's VAS charts: Medicare 3838.32 / 1443.21 and 4715.38 / 1443.21. + assert abs(XMR_LIMIT_MULTIPLIER - 2.659574) < 1e-6 + assert abs(R_UPPER_LIMIT_MULTIPLIER - 3.267287) < 1e-6 + + def test_multipliers_round_to_the_manuals_printed_values(self): + # The manual prints 2.66 and 3.268 (3.686/1.128, with 3.686 itself rounded). + assert round(XMR_LIMIT_MULTIPLIER, 2) == 2.66 + assert abs(R_UPPER_LIMIT_MULTIPLIER - 3.268) < 0.001 + + # ============================================================================ # Test: calculate_limits # ============================================================================ @@ -133,10 +165,10 @@ def test_b3_b4_raises_on_invalid_n(func): ('Xbar', dict(mean=10.0, sd=0.5, N=5), 9.286, 10.714, 0.01, 0.01), # S: sd=0.5, N=5 → LPL=0.5*b3(5)=0, UPL=0.5*b4(5)≈1.044 ('S', dict(sd=0.5, N=5), 0.0, 1.044, 0.001, 0.01), - # XmR: mean=10, mR=0.3 → 10±2.66*0.3=10±0.798 - ('XmR', dict(mean=10.0, mR=0.3), 10.0 - 2.66 * 0.3, 10.0 + 2.66 * 0.3, 0.001, 0.001), - # R: mR=0.3 → LPL=0, UPL=0.3*3.268 - ('R', dict(mR=0.3), 0.0, 0.3 * 3.268, 0.001, 0.001), + # XmR: mean=10, mR=0.3 → 10 ± (3/1.128)·0.3 = 10 ± 0.7979 + ('XmR', dict(mean=10.0, mR=0.3), 10.0 - 3 / 1.128 * 0.3, 10.0 + 3 / 1.128 * 0.3, 1e-12, 1e-12), + # R: mR=0.3 → LPL=0, UPL = (1 + 3·0.8525/1.128)·0.3 + ('R', dict(mR=0.3), 0.0, (1 + 3 * 0.8525 / 1.128) * 0.3, 1e-12, 1e-12), ], ids=['Xbar', 'S', 'XmR', 'R'], ) diff --git a/tests/test_vas_chart_constants.py b/tests/test_vas_chart_constants.py new file mode 100644 index 0000000..087bf9b --- /dev/null +++ b/tests/test_vas_chart_constants.py @@ -0,0 +1,86 @@ +"""Individuals and moving-range chart limits match Bishop's VAS software. + +The expected values are read from Bishop's VAS output decks of 29 September 2026 (Medicare and +PM SDS 2), at the precision the charts print. They pin the XmR constants E2 = 3/1.128 and +D4 = 1 + 3(0.8525)/1.128 from the manual (Eq 12.4, 12.5, 12.10, 12.11): with the rounded 2.66 and +3.268 the Medicare limits miss by 0.6 and 1.0, and the PM SDS 2 moving-range limit prints 2.84. +""" + +from pathlib import Path + +import pandas as pd +import pytest + +import processbehavior as pb + +VALIDATION = Path(__file__).resolve().parent.parent / 'validation' +MEDICARE_CSV = VALIDATION / 'aco_per_capita_expenditure.csv' +T100_CSV = VALIDATION / 'PBTESTDATABASE_T100.csv' + + +# precision=6 so the comparison rounds once, at the precision VAS prints. +@pytest.fixture(scope='module') +def medicare(): + df = pd.read_csv(MEDICARE_CSV) + return pb.formulate(df, response='PER CAPITA EXPENDITURE', factors=['ACO'], time='YEAR', precision=6) + + +@pytest.fixture(scope='module') +def pm_sds_2(): + df = pd.read_csv(T100_CSV, na_values=['*']) + return pb.formulate(df, response='PM SDS 2', factors=['FACTOR 1', 'FACTOR 2'], time='PRODUCTION TIME', precision=6) + + +class TestMedicareFirstChart: + """VAS MEDICARE ANALYSIS, slides 1-2: individuals chart ±3838.32, moving range 1443.21 / 4715.38.""" + + def test_individuals_half_width(self, medicare): + s = medicare.execute(chart='X', by=[]).get_statistics('X') + assert s['upl'] - s['center'] == pytest.approx(3838.32, abs=0.005) + assert s['center'] - s['lpl'] == pytest.approx(3838.32, abs=0.005) + + def test_moving_range_limits(self, medicare): + s = medicare.execute(chart='mR', by=[]).get_statistics('mR') + assert s['center'] == pytest.approx(1443.21, abs=0.005) + assert s['upl'] == pytest.approx(4715.38, abs=0.005) + + +class TestMedicarePerOrganisationCharts: + """VAS MEDICARE ANALYSIS, slides 3-6: the per-organisation individuals and moving-range charts. + + VAS prints six significant figures and appears to place the printed limits about a centre line + already rounded to that precision (ACO-002: 12637.4 + 2169.97 = 14807.37, printed 14807.4; the exact + limit is 14807.346), so one-decimal values allow 0.1. + """ + + @pytest.mark.parametrize( + 'aco, lpl, upl, lpl_tol, upl_tol', + [('ACO-001', 7503.91, 11991.3, 0.005, 0.05), ('ACO-002', 10467.4, 14807.4, 0.1, 0.1)], + ) + def test_individuals_limits(self, medicare, aco, lpl, upl, lpl_tol, upl_tol): + t = medicare.execute(chart='X', by=['ACO']).chart_table() + row = t[t['subgroup'] == aco].iloc[0] + assert row['lpl'] == pytest.approx(lpl, abs=lpl_tol) + assert row['upl'] == pytest.approx(upl, abs=upl_tol) + + @pytest.mark.parametrize('aco, center, upl', [('ACO-001', 843.62, 2756.36), ('ACO-002', 815.91, 2665.81)]) + def test_moving_range_limits(self, medicare, aco, center, upl): + t = medicare.execute(chart='mR', by=['ACO']).chart_table() + row = t[t['subgroup'] == aco].iloc[0] + assert row['center'] == pytest.approx(center, abs=0.005) + assert row['upl'] == pytest.approx(upl, abs=0.005) + + +class TestPmSds2FirstChart: + """VAS PM SDS 2 ANALYSIS, slides 1-2: 237.78 (235.47 / 240.09), moving range 0.87 / 2.83.""" + + def test_individuals_limits(self, pm_sds_2): + s = pm_sds_2.execute(chart='X', by=[]).get_statistics('X') + assert round(s['center'], 2) == 237.78 + assert round(s['lpl'], 2) == 235.47 + assert round(s['upl'], 2) == 240.09 + + def test_moving_range_limits(self, pm_sds_2): + s = pm_sds_2.execute(chart='mR', by=[]).get_statistics('mR') + assert round(s['center'], 2) == 0.87 + assert round(s['upl'], 2) == 2.83