From b8e57d89ebb5c94a066c52f5a54cc858b864cf84 Mon Sep 17 00:00:00 2001 From: Chris Nicholas <4948774+cnicholas@users.noreply.github.com> Date: Mon, 12 Oct 2026 19:38:07 -0400 Subject: [PATCH 1/4] perf(stratified): linear in the number of strata, not strata x rows Stratified X/mR looked each stratum up with a scan of the whole table, so with every row its own condition an X/mR pair took 48 s at 30,165 rows and an estimated 80 minutes at 307,683. It now takes 0.2 s and 1.9 s. - One groupby gives every stratum's row positions (_stratum_row_positions / _stratum_rows); the X/mR, Xbar, S and lane-boundary loops use it instead of df[df[col] == stratum]. - Per-stratum and per-phase X/mR limits use calculate_limits_vectorized (_mr_family_limits) instead of one calculate_limits call per row through DataFrame.apply. - The per-stratum statistics dict is built in one pass (_mr_stratum_statistics), keeping the numpy float64 values the row lookup returned. - Strata lists are searched as sets; a one-row stratum skips the lane-boundary call; focus_on reads result.strata once. Results are unchanged. 66 stratified scenarios (seven T100 responses; X/mR by one and both factors, phased, R2, mR alone; Xbar/S by time and n_mode='average'; a 30,165-row gene file) give identical data, statistics, metadata, strata and focus() output before and after, including value types. Full suite 2,511 passed; e2e gate 1511 pass, 0 fail, 43 pending, docs/reference/validation.md unchanged. scripts/benchmark.py adds a one-row-strata scenario (SDS-2) and tests/test_performance.py checks it scales linearly (per-row ratio 0.94 from 10K to 100K rows). test_stratified_insufficiency.py checks the one-row-strata result: nothing published, every stratum reported. Stratified Xbar/S is still about 2.5 ms per stratum (linear, not yet vectorised). --- CHANGELOG.md | 11 ++ processbehavior/analysis.py | 211 +++++++++++++------------ processbehavior/focused_result.py | 8 +- processbehavior/spc_constants.py | 11 +- scripts/benchmark.py | 20 ++- tests/test_performance.py | 8 + tests/test_stratified_insufficiency.py | 25 +++ 7 files changed, 178 insertions(+), 116 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 58ce42f..86f4b00 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -101,6 +101,17 @@ Equation numbers below are that manual's. PPU 5.74; ȳ alone gives 5.73). It equals ȳ in ADS 1. ### Fixed +- **Stratified charts take time in proportion to the data, however many strata there are.** + The stratified X/mR path looked each stratum up with a scan of the whole table, so its cost + grew with strata × rows. With every row its own condition, an X/mR pair took 48 s at 30,165 + rows, and an estimated 80 minutes at 307,683 (a gene-expression study, 6,033 genes × 51 + subjects). It now takes 0.2 s and 1.9 s. The same scan sat in the stratified Xbar/S loops and + the lane-boundary loop, the per-stratum limits were computed one `pd.Series` at a time, and + two lists of strata were searched as lists. Results are unchanged: every chart, statistic, + lane boundary and stratum list was compared before and after on the reference data and a + 30,165-row gene-effects file. `scripts/benchmark.py` adds a one-row-strata scenario, and + `tests/test_performance.py` checks that it scales linearly. Stratified Xbar/S still costs + about 2.5 ms per stratum (15 s at 6,033): linear, but not yet vectorised. - **Xbar and S charts with subgroups of unequal size follow the manual (Eqs 11-5, 11-16, 11-17).** sigma-hat is the average of each subgroup's S_r / c4(N_r) (was S-bar / c4(N_r)), and each subgroup's limits use its own N_r, so they step with subgroup size, as VAS draws them. N_r now diff --git a/processbehavior/analysis.py b/processbehavior/analysis.py index 050d687..8db5514 100644 --- a/processbehavior/analysis.py +++ b/processbehavior/analysis.py @@ -197,6 +197,65 @@ def _split_strata_by_sufficiency( return insufficient, published +def _stratum_row_positions(df: pd.DataFrame, stratify_col: str) -> dict: + """Row positions of every stratum in `df`, from one grouping pass. + + The per-stratum loops used to select with ``df[df[stratify_col] == stratum]``, which + compares the whole column for every stratum: O(strata x rows). With one row per + stratum that is quadratic in the data: 48s for an X/mR pair at 30,165 rows, an + estimated 80 minutes at 307,683. Look rows up here instead with :func:`_stratum_rows`. + + NaN strata are absent (groupby drops them), so they select no rows, exactly as + ``col == NaN`` did. + """ + return df.groupby(stratify_col, sort=False, observed=True).indices + + +def _stratum_rows(df: pd.DataFrame, positions: dict, stratum) -> pd.DataFrame: + """The rows of `df` in `stratum`, in `df`'s order (a view; copy before writing).""" + pos = positions.get(stratum) + return df.iloc[pos] if pos is not None else df.iloc[0:0] + + +def _mr_family_limits(mr_spec: _MRChartSpec, *, mean, mR) -> pd.DataFrame: + """Limits for one MR-family chart per row of a per-stratum or per-phase table. + + The array form of the scalar ``calculate_limits`` these sites used to call through + ``DataFrame.apply(axis=1)``: the same formulae, 10-12 manual Eqs 12-10/12-11 (X) and + 12-5 (mR). The scalar call built a ``pd.Series`` per stratum, about 0.1 ms each, so + it cost 30s per chart at 300,000 strata. ``mean`` is ignored for the mR chart. + """ + if mr_spec.limits_type == 'XmR': + return calculate_limits_vectorized('XmR', mean=mean, mR=mR) + return calculate_limits_vectorized(mr_spec.limits_type, mR=mR) + + +def _mr_stratum_statistics(mr_spec: _MRChartSpec, grouped: pd.DataFrame, stratify_col: str, + lpl_col: str, upl_col: str, round_to: int) -> dict: + """The per-stratum statistics dict of a stratified MR-family chart, in one pass. + + `grouped` has one row per stratum. Replaces a loop that found each stratum's row + with ``grouped[grouped[stratify_col] == stratum]``: O(strata^2). The values stay + numpy floats (``to_numpy``, not ``tolist``), as that loop's row lookup returned them. + """ + n = _mr_family_subgroup_n(mr_spec.chart_type) + return { + stratum: { + 'N': n, + 'center': round(center, round_to), + 'lpl': round(lpl, round_to), + 'upl': round(upl, round_to), + } + for stratum, center, lpl, upl in zip( + grouped[stratify_col].tolist(), + grouped['center'].to_numpy(), + grouped[lpl_col].to_numpy(), + grouped[upl_col].to_numpy(), + strict=True, + ) + } + + def _no_replication_alternative(chart: str) -> str: """What to reach for instead, when a subgroup chart has nothing to subgroup.""" return ( @@ -278,7 +337,9 @@ def _build_stratified_mr_payload( published_strata = list(strata) else: metadata['insufficient_strata'] = insufficient_strata or None - published_strata = [s for s in strata if s not in (insufficient_strata or [])] + # A set: against the list this was O(strata^2) when most strata are insufficient. + insufficient_set = set(insufficient_strata or ()) + published_strata = [s for s in strata if s not in insufficient_set] return { mr_spec.chart_type: { 'data': chart_out, @@ -951,6 +1012,28 @@ def _resolve_by_grouping(self, value_col: str) -> tuple[list[str], str | None, l groupby_cols = list(by) return groupby_cols, None, [] + def _stratum_lane_boundaries(self, out: pd.DataFrame, stratify_col: str, strata: list, + collapsed_factors: list[str]) -> dict: + """Lane boundaries for each stratum of a stratified chart, keyed by stratum. + + Strata with no boundary are left out. A stratum of fewer than two rows can't have + one, so it is skipped without the per-call cost, which matters when nearly every + row is its own stratum. + """ + lane_boundaries: dict = {} + if not collapsed_factors: + return lane_boundaries + positions = _stratum_row_positions(out, stratify_col) + for stratum in strata: + pos = positions.get(stratum) + if pos is None or len(pos) < 2: + continue + stratum_data = out.iloc[pos].reset_index(drop=True) + boundaries = self._calculate_lane_boundaries(stratum_data, collapsed_factors) + if boundaries: + lane_boundaries[stratum] = boundaries + return lane_boundaries + def _calculate_lane_boundaries(self, df: pd.DataFrame, collapsed_vars: list[str]) -> list[dict]: """ Calculate lane boundaries where collapsed factors change. @@ -1368,6 +1451,7 @@ def _calculate_xbar_stratified( _limits_col = self._resolve_limits_column(value_col, df) strata = df[stratify_col].unique().tolist() + positions = _stratum_row_positions(df, stratify_col) all_xbar_frames = [] chart_statistics = {} intermediates_per_stratum = {} @@ -1378,7 +1462,7 @@ def _calculate_xbar_stratified( insufficient_strata: list = [] for stratum in strata: - sdf = df[df[stratify_col] == stratum].copy() + sdf = _stratum_rows(df, positions, stratum).copy() # Group by time within this stratum # Note: we always compute mean directly from filtered data, NOT using @@ -1457,7 +1541,8 @@ def _calculate_xbar_stratified( chart_out = pd.concat(all_xbar_frames, ignore_index=True) # Publish only strata a caller can actually focus() into. - published_strata = [s for s in strata if s not in set(insufficient_strata)] + insufficient_set = set(insufficient_strata) + published_strata = [s for s in strata if s not in insufficient_set] # Select output columns group_col = groupby_cols[0] if groupby_cols else None @@ -1541,10 +1626,12 @@ def _calculate_s_stratified( df['_stratify_key'] = encode_rsg_series(df, stratify_by) stratify_col = '_stratify_key' strata = df[stratify_col].unique().tolist() + positions = _stratum_row_positions(df, stratify_col) per_stratum = None insufficient_strata = [] - all_strata = list(strata) + [s for s in insufficient_strata if s not in strata] + strata_set = set(strata) + all_strata = list(strata) + [s for s in insufficient_strata if s not in strata_set] all_s_frames = [] chart_statistics = {} @@ -1562,7 +1649,7 @@ def _calculate_s_stratified( n_max = inter['n_max'] n_avg = inter['n_avg'] else: - sdf = df[df[stratify_col] == stratum].copy() + sdf = _stratum_rows(df, positions, stratum).copy() _limits_col = self._resolve_limits_column(value_col, sdf) agg_dict = { 's': pd.NamedAgg(column=_limits_col, aggfunc='std'), @@ -1619,7 +1706,8 @@ def _calculate_s_stratified( _raise_no_replication('S', sds, insufficient_strata, len(all_strata)) chart_out = pd.concat(all_s_frames, ignore_index=True) - published_strata = [s for s in all_strata if s not in set(insufficient_strata)] + insufficient_set = set(insufficient_strata) + published_strata = [s for s in all_strata if s not in insufficient_set] group_col = groupby_cols[0] if groupby_cols else None cols_to_keep = [stratify_col, group_col, 's', 'center', 'lpl', 'upl', 'beyond_limits'] @@ -1983,17 +2071,7 @@ def _calculate_mr_chart_from_precomputed( # Recompute R-specific limits per phase within each stratum r_phase_stats = phase_stats[[stratify_col, '_phase_id', 'mR']].copy() - r_lims = r_phase_stats.apply( - lambda row: calculate_limits( - mean=0, - sd=0, - N=0, - mR=row['mR'], - limits_type=mr_spec.limits_type, - round_to=spec.round_to, - ), - axis=1, - ) + r_lims = _mr_family_limits(mr_spec, mean=0, mR=r_phase_stats['mR']) r_phase_stats = _join_limits(r_phase_stats, r_lims) r_phase_stats['center'] = r_phase_stats['mR'] @@ -2048,17 +2126,7 @@ def _calculate_mr_chart_from_precomputed( # grouped has columns: stratify_col, center (=mean), mR, lpl, upl (X limits) r_grouped = grouped.rename(columns={'center': 'xmr_center', 'mR': 'center'}) - r_lims = r_grouped.apply( - lambda row: calculate_limits( - mean=0, - sd=0, - N=0, - mR=row['center'], - limits_type=mr_spec.limits_type, - round_to=spec.round_to, - ), - axis=1, - ) + r_lims = _mr_family_limits(mr_spec, mean=0, mR=r_grouped['center']) r_grouped = _join_limits(r_grouped, r_lims, rsuffix='_r') @@ -2084,15 +2152,9 @@ def _calculate_mr_chart_from_precomputed( r_out = self._add_beyond_limits_flag(r_out, value_col=plot_col) # Build statistics - chart_statistics = {} - for stratum in strata: - row = r_grouped[r_grouped[stratify_col] == stratum].iloc[0] - chart_statistics[stratum] = { - 'N': _mr_family_subgroup_n(mr_spec.chart_type), - 'center': round(row['center'], spec.round_to), - 'lpl': round(row[r_lpl_col], spec.round_to), - 'upl': round(row[r_upl_col], spec.round_to), - } + chart_statistics = _mr_stratum_statistics( + mr_spec, r_grouped, stratify_col, r_lpl_col, r_upl_col, spec.round_to, + ) # Format output chart_out = self._build_output_columns( @@ -2186,30 +2248,10 @@ def _calculate_mr_chart_stratified( # noqa: C901 # 5. Compute per-phase limits if mr_spec.center_source == 'mean': phase_stats['center'] = phase_stats['_mean'] - phase_lims = phase_stats.apply( - lambda row: calculate_limits( - mean=row['_mean'], - sd=0, - N=0, - mR=row['mR'], - limits_type=mr_spec.limits_type, - round_to=spec.round_to, - ), - axis=1, - ) + phase_lims = _mr_family_limits(mr_spec, mean=phase_stats['_mean'], mR=phase_stats['mR']) else: # center_source == 'mR' (R chart) phase_stats['center'] = phase_stats['mR'] - phase_lims = phase_stats.apply( - lambda row: calculate_limits( - mean=0, - sd=0, - N=0, - mR=row['mR'], - limits_type=mr_spec.limits_type, - round_to=spec.round_to, - ), - axis=1, - ) + phase_lims = _mr_family_limits(mr_spec, mean=0, mR=phase_stats['mR']) phase_stats = _join_limits(phase_stats, phase_lims) # 6. Merge per-phase limits back to row-level data @@ -2228,13 +2270,7 @@ def _calculate_mr_chart_stratified( # noqa: C901 out = self._add_beyond_limits_flag(out, value_col=plot_col) # 9. Lane boundaries per stratum - lane_boundaries = {} - if collapsed_factors: - for stratum in strata: - stratum_data = out[out[stratify_col] == stratum].reset_index(drop=True) - boundaries = self._calculate_lane_boundaries(stratum_data, collapsed_factors) - if boundaries: - lane_boundaries[stratum] = boundaries + lane_boundaries = self._stratum_lane_boundaries(out, stratify_col, strata, collapsed_factors) # 10. Statistics per stratum: limits vary per phase chart_statistics = { @@ -2311,29 +2347,10 @@ def _calculate_mr_chart_stratified( # noqa: C901 # For X: center=mean, limits based on mean ± E2*mR # For mR: center=mR, limits based on 0 to D4*mR if mr_spec.center_source == 'mean': - lims = grouped.apply( - lambda row: calculate_limits( - mean=row['mean'], - sd=0, - N=0, - mR=(row['mR'] if pd.notna(row['mR']) else 0.0), - limits_type=mr_spec.limits_type, - round_to=spec.round_to, - ), - axis=1, - ) + # A one-row stratum has no moving range: its X limits collapse to the mean. + lims = _mr_family_limits(mr_spec, mean=grouped['mean'], mR=grouped['mR'].fillna(0.0)) else: # center_source == 'mR' - lims = grouped.apply( - lambda row: calculate_limits( - mean=0, - sd=0, - N=0, - mR=row['mR'], - limits_type=mr_spec.limits_type, - round_to=spec.round_to, - ), - axis=1, - ) + lims = _mr_family_limits(mr_spec, mean=0, mR=grouped['mR']) # Normalize/join lpl/upl grouped = _join_limits(grouped, lims) @@ -2363,24 +2380,10 @@ def _calculate_mr_chart_stratified( # noqa: C901 strata = grouped[stratify_col].tolist() # Lane boundaries per stratum - lane_boundaries = {} - if collapsed_factors: - for stratum in strata: - stratum_data = out[out[stratify_col] == stratum].reset_index(drop=True) - boundaries = self._calculate_lane_boundaries(stratum_data, collapsed_factors) - if boundaries: - lane_boundaries[stratum] = boundaries + lane_boundaries = self._stratum_lane_boundaries(out, stratify_col, strata, collapsed_factors) # Build statistics nested by stratum - chart_statistics = {} - for stratum in strata: - row = grouped[grouped[stratify_col] == stratum].iloc[0] - chart_statistics[stratum] = { - 'N': _mr_family_subgroup_n(mr_spec.chart_type), - 'center': round(row['center'], spec.round_to), - 'lpl': round(row['lpl'], spec.round_to), - 'upl': round(row['upl'], spec.round_to), - } + chart_statistics = _mr_stratum_statistics(mr_spec, grouped, stratify_col, 'lpl', 'upl', spec.round_to) # Format output extra_cols = [c for c in stratify_by if c != spec.rsg_var_name and c != spec.time_var] diff --git a/processbehavior/focused_result.py b/processbehavior/focused_result.py index 84085d2..1bc4982 100644 --- a/processbehavior/focused_result.py +++ b/processbehavior/focused_result.py @@ -40,14 +40,16 @@ def focus_on(result: AnalysisResult, stratum: str) -> AnalysisResult: ProcessBehaviorError If the chart lacks both stratify_col metadata and an rsg column. """ - if not result.strata: + # `strata` is a property that rebuilds the list on every read: read it once. + strata = result.strata + if not strata: raise ValidationError( "Cannot focus: this result is not stratified. " "Use result.strata to check available subgroups." ) - if stratum not in result.strata: + if stratum not in strata: raise ValidationError( - f"Stratum '{stratum}' not found. Available strata: {result.strata}" + f"Stratum '{stratum}' not found. Available strata: {strata}" ) focused_charts: dict[str, dict[str, Any]] = {} diff --git a/processbehavior/spc_constants.py b/processbehavior/spc_constants.py index 459e03f..1b4323a 100644 --- a/processbehavior/spc_constants.py +++ b/processbehavior/spc_constants.py @@ -389,11 +389,12 @@ def calculate_limits_vectorized( Notes ----- **Conversion is deliberately partial.** The four primary Xbar/S limit computations in - ``analysis.py`` use this; the phase-segment and per-stratum sites still call the scalar - form. Those iterate over phases or strata (a handful of rows), not subgroups, so they - were not hot when measured — converting them would carry risk for no gain. If a future - profile shows otherwise, they are the next candidates; do not assume the absence of a - call here means a site was checked and rejected on correctness grounds. + ``analysis.py`` use this, and so do the stratified X/mR paths (through + ``_mr_family_limits``): those once iterated over a handful of strata, until a study + with one stratum per row (300,000 of them) made the scalar calls cost 30s per chart. + The ungrouped phase-segment sites still call the scalar form; they iterate over phases, + not subgroups or strata. Do not assume the absence of a call here means a site was + checked and rejected on correctness grounds. """ index = next((arg.index for arg in (mean, sd, N, mR) if isinstance(arg, pd.Series)), None) diff --git a/scripts/benchmark.py b/scripts/benchmark.py index a0b0e50..7c138f8 100644 --- a/scripts/benchmark.py +++ b/scripts/benchmark.py @@ -2,8 +2,9 @@ """Performance benchmark + drift check for processbehavior. Runs a fixed set of scenarios — design {SDS-1 exact-R2, SDS-2 MA2 worst-case} × size -{10K, 100K, 1M} × op {init, formulate, execute, execute-stratified} — prints a report table, writes -``benchmarks/last_run.json``, and compares against the committed ``benchmarks/baseline.json``. +{10K, 100K, 1M} × op {init, formulate, execute, execute-stratified, and for SDS-2 one-row +strata} — prints a report table, writes ``benchmarks/last_run.json``, and compares against the +committed ``benchmarks/baseline.json``. The baseline is *same-machine* (absolute times are hardware-specific), so a regression flag is meaningful when you re-run on the same box. Cross-machine, rely on the shape checks in @@ -109,13 +110,24 @@ def run_scenario(sds: int, target_rows: int) -> list[dict]: peak = _peak_mb(lambda: ProcessBehavior(df).formulate(**spec)) common = {"sds": sds, "size": target_rows, "rows": rows} - return [ + results = [ {**common, "op": "init", "seconds": t_init}, {**common, "op": "formulate", "seconds": t_form}, {**common, "op": f"execute_{exec_chart.lower()}", "seconds": t_exec}, {**common, "op": "execute_x_stratified", "seconds": t_strat}, - {**common, "op": "formulate_peak", "seconds": None, "peak_mb": peak, "input_mb": input_mb}, ] + if sds == 2: + # One row per stratum: the time column taken as a factor makes every row its own + # condition. The scenario above has at most a few thousand strata; this one has as + # many as rows. The stratified X/mR path was quadratic in strata until 2026-10-12: + # about 80 minutes per chart at 307,683 one-row strata, unnoticed because nothing + # measured this shape. + one_row = [*DESIGNS[sds]["factors"], "time"] + one_row_study = pdf.formulate(response="y", factors=one_row) + t_one_row, _ = _time(lambda: one_row_study.execute(chart="X", by=one_row, companion=True), repeat) + results.append({**common, "op": "execute_x_one_row_strata", "seconds": t_one_row}) + results.append({**common, "op": "formulate_peak", "seconds": None, "peak_mb": peak, "input_mb": input_mb}) + return results def run_benchmarks(quick: bool = False) -> list[dict]: diff --git a/tests/test_performance.py b/tests/test_performance.py index c02e87f..ad113f0 100644 --- a/tests/test_performance.py +++ b/tests/test_performance.py @@ -592,6 +592,14 @@ def test_ma2_scales_near_linearly(self, benchmark_results): print(f'\nSDS-2 MA2 per-row: {small:.2f}us (10K) -> {large:.2f}us (100K), ratio {large / small:.2f}') assert large <= small * self.MAX_RATIO, f'MA2 per-row grew {large / small:.1f}x (>{self.MAX_RATIO}x)' + def test_one_row_strata_scale_linearly(self, benchmark_results): + """Stratified X/mR with every row its own stratum. A per-stratum table scan makes + this quadratic: 10x the rows took ~100x the time before 2026-10-12.""" + small = _per_row_us(benchmark_results, 2, 'execute_x_one_row_strata', 10_000) + large = _per_row_us(benchmark_results, 2, 'execute_x_one_row_strata', 100_000) + print(f'\nOne-row strata per-row: {small:.2f}us (10K) -> {large:.2f}us (100K), ratio {large / small:.2f}') + assert large <= small * self.MAX_RATIO, f'one-row strata per-row grew {large / small:.1f}x (>{self.MAX_RATIO}x)' + @pytest.mark.slow @pytest.mark.benchmark diff --git a/tests/test_stratified_insufficiency.py b/tests/test_stratified_insufficiency.py index d9691b2..4e9c33e 100644 --- a/tests/test_stratified_insufficiency.py +++ b/tests/test_stratified_insufficiency.py @@ -184,3 +184,28 @@ def test_sufficient_data_publishes_every_stratum(df, response, chart): assert _meta(result, chart)['insufficient_strata'] is None for stratum in result.strata: assert len(result.focus(stratum).get_chart(chart)) > 0 + + +# --------------------------------------------------------------------------- +# Every row its own stratum (X/mR) +# --------------------------------------------------------------------------- + + +def test_one_row_strata_publish_nothing_and_report_every_stratum(df): + """Every row its own condition, as in a 307,683-row study that took 80 minutes per chart. + + Each stratum holds one reading, so none has a moving range and none is published, but + every one is reported. The time column taken as a factor does it on the reference data. + The cost of this shape is guarded in tests/test_performance.py (one-row strata). + """ + by = [*FACTORS, TIME] + study = pb.formulate(df, response='PM SDS 2', factors=by) + result = study.execute(chart='X', by=by, companion=True) + + rows = int(df['PM SDS 2'].notna().sum()) + assert result.strata == [] + for chart in ('X', 'mR'): + meta = _meta(result, chart) + assert len(meta['insufficient_strata']) == rows + assert len(set(meta['insufficient_strata'])) == rows + assert meta['lane_boundaries'] is None From c672584831059f11ba1d45923ef076f4e7500bc8 Mon Sep 17 00:00:00 2001 From: Chris Nicholas <4948774+cnicholas@users.noreply.github.com> Date: Mon, 12 Oct 2026 20:10:04 -0400 Subject: [PATCH 2/4] perf(plotting): add lane dividers in one assignment; none above 1,000 Each lane divider was added with its own fig.add_shape call, and plotly re-validates every shape already on the figure, so drawing them cost O(lanes^2): 200 dividers 1.3 s, 400 5.7 s, 800 27 s. An X chart of a 6,033-gene study (by=[], one lane per gene) would have taken about 25 minutes; the app's first Primary Analyses fold draws exactly that chart. - The first shape or annotation goes through add_shape / add_annotation so plotly resolves the subplot's axis references for row/col; the rest copy them and are appended in one assignment. The figure is identical to the call-by-call one (to_dict equal for single and faceted charts, with and without labels, two themes); 800 dividers now take 0.03 s. - Above MAX_LANE_DIVIDERS (1,000) none are drawn: they would merge into a solid band at any chart width, and the browser still draws every one. Every Primary Analyses fold on a 30,165-row, 6,033-gene file now renders (figure, table, signals) in at most 2.4 s, 8 s for all of them; at 307,683 rows the slowest fold is about 20 s. Suite 2,517 passed; e2e 1511 pass, 0 fail, 43 pending, validation.md unchanged. --- CHANGELOG.md | 7 +++ processbehavior/plotting/lane_boundaries.py | 56 +++++++++++++-------- tests/test_plotting_decorations.py | 43 +++++++++++++++- 3 files changed, 85 insertions(+), 21 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 86f4b00..7a39165 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -112,6 +112,13 @@ Equation numbers below are that manual's. 30,165-row gene-effects file. `scripts/benchmark.py` adds a one-row-strata scenario, and `tests/test_performance.py` checks that it scales linearly. Stratified Xbar/S still costs about 2.5 ms per stratum (15 s at 6,033): linear, but not yet vectorised. +- **Charts with many lanes draw in a fraction of a second.** Each lane divider was added to the + figure with its own `add_shape` call, and plotly re-validates every shape already there, so the + cost grew with the square of the number of lanes: 800 dividers took 27 s, and an X chart of a + 6,033-gene study (one lane per gene) would have taken about 25 minutes. They are now added in + one assignment, which builds the same figure (0.03 s for 800). Above 1,000 + (`plotting.lane_boundaries.MAX_LANE_DIVIDERS`) none are drawn: the dividers would merge into a + solid band at any chart width. - **Xbar and S charts with subgroups of unequal size follow the manual (Eqs 11-5, 11-16, 11-17).** sigma-hat is the average of each subgroup's S_r / c4(N_r) (was S-bar / c4(N_r)), and each subgroup's limits use its own N_r, so they step with subgroup size, as VAS draws them. N_r now diff --git a/processbehavior/plotting/lane_boundaries.py b/processbehavior/plotting/lane_boundaries.py index 0d0d5d6..9edef09 100644 --- a/processbehavior/plotting/lane_boundaries.py +++ b/processbehavior/plotting/lane_boundaries.py @@ -10,6 +10,12 @@ from .themes import ChartTheme +# More dividers than this can't be told apart at any chart width: they merge into a solid band, +# and the browser still has to draw every one. Above it none are drawn. A chart with that many +# lanes (one per condition: 6,032 for a 6,033-gene study) is read from its points and limits. +MAX_LANE_DIVIDERS = 1_000 + + def add_lane_boundaries( fig: go.Figure, lane_boundaries: list[dict] | None, @@ -25,6 +31,7 @@ def add_lane_boundaries( Unified function for both single charts and faceted subplots. Lane boundaries show where collapsed factors change within the chart, helping distinguish groups of observations from different factor levels. + Nothing is drawn above ``MAX_LANE_DIVIDERS`` boundaries. Parameters ---------- @@ -43,7 +50,7 @@ def add_lane_boundaries( show_labels : bool, default True Whether to show factor labels at boundary positions """ - if not lane_boundaries: + if not lane_boundaries or len(lane_boundaries) > MAX_LANE_DIVIDERS: return y_min, y_max = y_range @@ -52,27 +59,36 @@ def add_lane_boundaries( if row is not None and col is not None: subplot_kwargs = {'row': row, 'col': col} + line = dict(color=theme.lane_boundary_color, dash=theme.lane_boundary_dash, width=theme.lane_boundary_width) + font = dict(size=theme.lane_boundary_annotation_size, color=theme.lane_boundary_color) + shapes, annotations = [], [] for boundary in lane_boundaries: x_pos = boundary['position'] label = boundary.get('label', '') - - fig.add_shape( - type='line', - x0=x_pos, - x1=x_pos, - y0=y_min, - y1=y_max, - line=dict(color=theme.lane_boundary_color, dash=theme.lane_boundary_dash, width=theme.lane_boundary_width), - **subplot_kwargs, - ) - + shapes.append(dict(type='line', x0=x_pos, x1=x_pos, y0=y_min, y1=y_max, line=line)) if show_labels and label: - fig.add_annotation( - x=x_pos, - y=y_max, - text=label, - showarrow=False, - yanchor='bottom', - font=dict(size=theme.lane_boundary_annotation_size, color=theme.lane_boundary_color), - **subplot_kwargs, + annotations.append( + dict(x=x_pos, y=y_max, text=label, showarrow=False, yanchor='bottom', font=font) ) + + _add_all(fig, 'shapes', fig.add_shape, shapes, subplot_kwargs) + _add_all(fig, 'annotations', fig.add_annotation, annotations, subplot_kwargs) + + +def _add_all(fig: go.Figure, prop: str, add_one, items: list[dict], subplot_kwargs: dict) -> None: + """Add layout shapes or annotations in one assignment instead of one call each. + + ``fig.add_shape`` re-validates every shape already on the figure, so adding N one at a + time costs O(N^2): 800 dividers took 27s, and 6,032 would take about 25 minutes. The first + item goes through ``add_one`` so plotly resolves the subplot's axis references for + ``row``/``col``. The rest copy those references and are appended in one assignment. The + figure is the same as one built call by call. + """ + if not items: + return + add_one(**items[0], **subplot_kwargs) + if len(items) == 1: + return + first = fig.layout[prop][-1] + refs = {'xref': first.xref, 'yref': first.yref} + fig.layout[prop] = fig.layout[prop] + tuple({**item, **refs} for item in items[1:]) diff --git a/tests/test_plotting_decorations.py b/tests/test_plotting_decorations.py index 68364b6..0583b5b 100644 --- a/tests/test_plotting_decorations.py +++ b/tests/test_plotting_decorations.py @@ -12,7 +12,7 @@ from plotly.subplots import make_subplots from processbehavior.plotting import get_theme -from processbehavior.plotting.lane_boundaries import add_lane_boundaries +from processbehavior.plotting.lane_boundaries import MAX_LANE_DIVIDERS, add_lane_boundaries from processbehavior.plotting.limits import ( add_stepped_limit_line, build_stepped_coordinates, @@ -224,6 +224,47 @@ def test_boundary_without_label(self, theme): assert len(fig.layout.shapes) == 1 assert len(fig.layout.annotations) == 0 + @staticmethod + def _one_call_each(fig, boundaries, y_range, theme, row=None, col=None, show_labels=True): + """The figure add_lane_boundaries built before batching: one add_shape per divider.""" + y_min, y_max = y_range + kw = {'row': row, 'col': col} if row is not None and col is not None else {} + for b in boundaries: + fig.add_shape(type='line', x0=b['position'], x1=b['position'], y0=y_min, y1=y_max, + line=dict(color=theme.lane_boundary_color, dash=theme.lane_boundary_dash, + width=theme.lane_boundary_width), **kw) + if show_labels and b.get('label', ''): + fig.add_annotation(x=b['position'], y=y_max, text=b['label'], showarrow=False, yanchor='bottom', + font=dict(size=theme.lane_boundary_annotation_size, + color=theme.lane_boundary_color), **kw) + + @pytest.mark.parametrize(('make', 'subplot'), [ + (go.Figure, {}), + (lambda: make_subplots(rows=2, cols=2), {'row': 2, 'col': 1}), + ]) + @pytest.mark.parametrize('show_labels', [True, False]) + def test_batched_figure_matches_one_call_each(self, theme, make, subplot, show_labels): + """Adding the dividers in one assignment builds the same figure, references included.""" + boundaries = [{'position': i * 3, 'label': f'L{i}' if i % 4 else ''} for i in range(1, 40)] + expected, actual = make(), make() + for fig in (expected, actual): + fig.add_shape(type='rect', x0=0, x1=1, y0=0, y1=1) # a shape already on the figure + self._one_call_each(expected, boundaries, (-1, 5), theme, show_labels=show_labels, **subplot) + add_lane_boundaries(actual, boundaries, (-1, 5), theme, show_labels=show_labels, **subplot) + assert actual.to_dict() == expected.to_dict() + + def test_draws_up_to_the_divider_limit(self, theme): + fig = go.Figure() + add_lane_boundaries(fig, [{'position': i} for i in range(MAX_LANE_DIVIDERS)], (0, 1), theme) + assert len(fig.layout.shapes) == MAX_LANE_DIVIDERS + + def test_draws_none_above_the_divider_limit(self, theme): + """Past the limit dividers merge into a solid band: one lane per gene of a 6,033-gene study.""" + fig = go.Figure() + add_lane_boundaries(fig, [{'position': i, 'label': 'g'} for i in range(6_032)], (0, 1), theme) + assert len(fig.layout.shapes) == 0 + assert len(fig.layout.annotations) == 0 + # ========================================================================= # Stats Box From 3ac08da3f33d23ff2f268386e965e5528de49619 Mon Sep 17 00:00:00 2001 From: Chris Nicholas <4948774+cnicholas@users.noreply.github.com> Date: Tue, 13 Oct 2026 23:44:54 -0400 Subject: [PATCH 3/4] fix(signals): detect_signals() works on phased charts (#141) detect_signals() raised TypeError on every phased=True chart. Each phase has its own centre, so the statistics carry center=None, and the per-row zone calculation subtracted that scalar. The chart marks run rules not applicable (metadata run_rules_applicable=False, 70503be), but only the renderer read the flag. - result_signals passes run_rules_applicable to SignalDetector.detect. When it is False and the rules are the chart-type default, the default is Rule 1 only: rules 2-8 are not applicable rather than skipped, so the result is not marked partial on their account. - An explicit rule list runs as given (SignalConfig's "explicit rules, no filtering" contract). The per-row zones, the run rule (rule 4) and the violation records read each row's own centre when the scalar centre is None. The limits and the points the chart flags were always correct; Rule 1 now reports the same points (flagged_observations equals the rows with beyond_limits != 0). Unphased charts are unchanged: X still evaluates rules 1-8, Xbar/S Rule 1. tests/test_phased_limits.py TestPhasedSignals: combined, companion and stratified (via focus) phased charts; the violation records carry the phase's centre and limits; explicit rules 2 and 4 judged against each phase's centre; unphased X keeps every rule. --- CHANGELOG.md | 6 +++ processbehavior/analysis_result.py | 4 ++ processbehavior/result_signals.py | 2 + processbehavior/signals/detector.py | 33 ++++++++++++- processbehavior/signals/detectors.py | 3 ++ tests/test_phased_limits.py | 74 ++++++++++++++++++++++++++++ 6 files changed, 120 insertions(+), 2 deletions(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 7a39165..a95766f 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -119,6 +119,12 @@ Equation numbers below are that manual's. one assignment, which builds the same figure (0.03 s for 800). Above 1,000 (`plotting.lane_boundaries.MAX_LANE_DIVIDERS`) none are drawn: the dividers would merge into a solid band at any chart width. +- **`detect_signals()` works on phased charts.** It raised `TypeError` on every `phased=True` + chart: each phase has its own centre, so the statistics carry `center=None`, and the zone + calculation read that scalar. The chart already marks run rules not applicable, but only the + plot read the flag. Phased charts now get Rule 1 by default, judged against each phase's + limits, with the same points the chart flags; an explicit rule list runs as given, against each + row's own centre. The limits and the flagged points on the chart were always correct. - **Xbar and S charts with subgroups of unequal size follow the manual (Eqs 11-5, 11-16, 11-17).** sigma-hat is the average of each subgroup's S_r / c4(N_r) (was S-bar / c4(N_r)), and each subgroup's limits use its own N_r, so they step with subgroup size, as VAS draws them. N_r now diff --git a/processbehavior/analysis_result.py b/processbehavior/analysis_result.py index 9bf8619..4a86b39 100644 --- a/processbehavior/analysis_result.py +++ b/processbehavior/analysis_result.py @@ -1039,6 +1039,10 @@ def detect_signals( ``signals.is_partial`` is True and ``signals.summary`` says which rules were evaluated. + On a phased chart (``phased=True``) each phase has its own centre line, so + the default rules are Rule 1 only, judged against each phase's limits. An + explicit rule list still runs as given, against each row's own centre. + Access violations: >>> signals.by_rule['rule_2'] # Rule 2 violations diff --git a/processbehavior/result_signals.py b/processbehavior/result_signals.py index b561ab1..37f247d 100644 --- a/processbehavior/result_signals.py +++ b/processbehavior/result_signals.py @@ -146,4 +146,6 @@ def _detect_for_chart( value_col=value_col, chart_name=chart_name, chart_type=chart_type, + # Phased charts mark run rules not applicable; until now only the renderer read it. + run_rules_applicable=metadata.get("run_rules_applicable", True), ) diff --git a/processbehavior/signals/detector.py b/processbehavior/signals/detector.py index 414f056..8bb2753 100644 --- a/processbehavior/signals/detector.py +++ b/processbehavior/signals/detector.py @@ -91,6 +91,7 @@ def detect( value_col: str = 'mean', chart_name: str = 'Chart', chart_type: str = 'Xbar', + run_rules_applicable: bool = True, ) -> SignalResult: """ Detect signals in control chart data. @@ -109,6 +110,11 @@ def detect( Name of chart for reporting chart_type : str, default 'Xbar' Chart type ('Xbar', 'S', 'X', 'mR') - determines applicable rules + run_rules_applicable : bool, default True + False for charts whose centre moves along the series (phased limits, + ``metadata['run_rules_applicable']``). The chart-type default rules + narrow to Rule 1; an explicit rule list runs as given, against each + row's own centre. Returns ------- @@ -126,6 +132,13 @@ def detect( # Get applicable rules for this chart type applicable_rules = config.get_rules_for_chart(chart_type) + # Phased limits give each phase its own centre, so the chart marks run rules + # 2-8 not applicable: they are not part of the default set, rather than + # skipped (a skip would mark the evaluation partial). An explicit rule list is + # the analyst's choice and runs as given. + if not run_rules_applicable and config.enabled_rules == 'default': + applicable_rules = [r for r in applicable_rules if r == 'rule_1'] + # Log rule selection if logger.isEnabledFor(logging.INFO): all_rules = [f'rule_{i}' for i in range(1, 9)] @@ -291,7 +304,7 @@ def _build_result( 'rule_number': int(rule_name.split('_')[1]), 'description': self.RULE_DESCRIPTIONS[rule_name], 'value': data.loc[idx, value_col], - 'center': stats['center'], + 'center': self._row_center(data, stats, idx), 'upl': ucl, 'lpl': lcl, } @@ -309,6 +322,22 @@ def _limits_vary(self, stats: dict) -> bool: """Check if control limits vary (per-row limits).""" return bool(stats.get('limits_vary')) + @staticmethod + def _center(data: pd.DataFrame, stats: dict): + """The centre line: the scalar, or the per-row column when there is none. + + Phased charts give each phase its own centre, so their statistics carry + ``center=None`` and the per-phase value lives on each row. + """ + if stats.get('center') is None and 'center' in data.columns: + return data['center'] + return stats['center'] + + def _row_center(self, data: pd.DataFrame, stats: dict, idx): + """The centre line at one row.""" + center = self._center(data, stats) + return center.loc[idx] if isinstance(center, pd.Series) else center + def _calculate_per_row_zones(self, data: pd.DataFrame, stats: dict, config: SignalConfig) -> pd.DataFrame: """ Calculate per-row zone boundaries for varying process limits. @@ -334,7 +363,7 @@ def _calculate_per_row_zones(self, data: pd.DataFrame, stats: dict, config: Sign if 'upl' not in data.columns or 'lpl' not in data.columns: raise ValueError("Cannot calculate per-row zones: data missing 'upl' or 'lpl' columns") - center = stats['center'] + center = self._center(data, stats) zone_def = config.zone_definition # Calculate sigma per row diff --git a/processbehavior/signals/detectors.py b/processbehavior/signals/detectors.py index 38e18e0..2c2431b 100644 --- a/processbehavior/signals/detectors.py +++ b/processbehavior/signals/detectors.py @@ -173,6 +173,9 @@ def detect_run(data: pd.DataFrame, stats: dict, value_col: str, length: int = 8) """ values = data[value_col] center = stats.get('center', stats.get('Mean')) # Try 'center' first, fallback to 'Mean' + if center is None and 'center' in data.columns: + # Phased charts: each phase has its own centre, carried on the rows. + center = data['center'] # Above or below center above = (values > center).astype(int) diff --git a/tests/test_phased_limits.py b/tests/test_phased_limits.py index ff451cf..d4c5cf6 100644 --- a/tests/test_phased_limits.py +++ b/tests/test_phased_limits.py @@ -412,3 +412,77 @@ def test_phased_residual_limits_vary(self, sds2_study): result = sds2_study.execute(chart='X', by=[], phased=True, value='R2') data = result.get_chart('X') assert data['center'].nunique() > 1 or data['upl'].nunique() > 1 + + +# ============================================================================= +# Phased + signal detection +# ============================================================================= + + +class TestPhasedSignals: + """detect_signals() on phased charts: Rule 1 against each phase's own limits. + + Phased statistics carry ``center=None`` (each phase has its own centre) and the + chart marks run rules not applicable. detect_signals() used to read the scalar + centre and raise TypeError on every phased chart. + """ + + def test_phased_x_rule_1_against_phase_limits(self, sds1_study): + result = sds1_study.execute(chart='X', by=[], phased=True) + signals = result.detect_signals(chart='X') + assert signals.rules_evaluated == ['rule_1'] + assert signals.rules_skipped == {} + assert not signals.is_partial + data = result.get_chart('X') + beyond = set(data.index[data['beyond_limits'] != 0]) + assert beyond, 'fixture should have a point beyond its phase limits' + assert signals.flagged_observations == beyond + + def test_violation_records_the_phase_centre(self, sds1_study): + result = sds1_study.execute(chart='X', by=[], phased=True) + data = result.get_chart('X') + violations = result.detect_signals(chart='X').to_dataframe() + assert violations['center'].notna().all() + for _, row in violations.iterrows(): + assert row['center'] == data.loc[row['obs_id'], 'center'] + assert row['upl'] == data.loc[row['obs_id'], 'upl'] + + def test_phased_companion_both_charts(self, sds1_study): + result = sds1_study.execute(chart='X', by=[], phased=True, companion=True) + signals = result.detect_signals() + assert set(signals) == {'X', 'mR'} + for name, sig in signals.items(): + data = result.get_chart(name) + assert sig.rules_evaluated == ['rule_1'] + assert sig.flagged_observations == set(data.index[data['beyond_limits'] != 0]) + + def test_phased_stratified_focus(self, sds1_study): + result = sds1_study.execute(chart='X', by=['factor 1'], phased=True, companion=True) + for stratum in result.strata: + signals = result.focus(stratum).detect_signals() + assert {sig.rules_evaluated[0] for sig in signals.values()} == {'rule_1'} + + def test_explicit_rules_run_against_each_phase_centre(self): + """An explicit rule list runs as given; run and zone rules read the row's centre.""" + from processbehavior.signals import SignalConfig, SignalDetector + + # Phase 1 sits above its centre 0, phase 2 below its centre 10: one run each. + # Phase 2's last two points fall in its lower zone A (7 to 8, sigma 1). + data = pd.DataFrame({ + 'x': [1.0] * 8 + [9.0] * 6 + [7.5, 7.5], + 'center': [0.0] * 8 + [10.0] * 8, + 'upl': [3.0] * 8 + [13.0] * 8, + 'lpl': [-3.0] * 8 + [7.0] * 8, + }) + stats = {'N': None, 'center': None, 'lpl': None, 'upl': None, 'limits_vary': True} + signals = SignalDetector().detect( + data, stats, SignalConfig(enabled_rules=['rule_2', 'rule_4'], min_observations=1), + value_col='x', chart_type='X', run_rules_applicable=False, + ) + assert signals.rules_evaluated == ['rule_2', 'rule_4'] + assert set(signals.by_rule['rule_4']['obs_id']) == {7, 15} + assert set(signals.by_rule['rule_2']['obs_id']) == {15} + + def test_unphased_x_keeps_every_rule(self, sds1_study): + signals = sds1_study.execute(chart='X', by=[]).detect_signals(chart='X') + assert signals.rules_evaluated == [f'rule_{i}' for i in range(1, 9)] From 5947ec5681a7aa9b9d76a879e31f6c20ee12a789 Mon Sep 17 00:00:00 2001 From: Chris Nicholas <4948774+cnicholas@users.noreply.github.com> Date: Tue, 13 Oct 2026 23:46:08 -0400 Subject: [PATCH 4/4] fix(study): value='response' charts the response (#142) execute(value='response') is in the signature (Literal['response', 'R1', ...]) and the docstring, but every chart refused it: "Valid charts for RESPONSE: None". _validate_execute_request tells residuals apart by value.upper().startswith('R'), which 'RESPONSE' also satisfies, so the value went down the residual path and failed the pairing check. why_not() and supports_calibration() refused it the same way. _normalise_value turns 'response' (any case) into None, the response's own spelling, at the start of execute(), supports_calibration() and why_not(); _resolve_value_column already treats None as the response. recentered=True with the response is still refused, and the message now says why: the response is already on the measurement scale. tests/test_execute_validation.py TestValueResponse: 'response', 'Response' and 'RESPONSE' give the same charts, statistics and metadata as value=None on Xbar, S, Xbar+S, X, mR and Histogram; why_not and supports_calibration agree; recentered stays refused. --- CHANGELOG.md | 4 +++ processbehavior/study.py | 21 ++++++++++++- tests/test_execute_validation.py | 51 ++++++++++++++++++++++++++++++++ 3 files changed, 75 insertions(+), 1 deletion(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index a95766f..71cf7ca 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -125,6 +125,10 @@ Equation numbers below are that manual's. plot read the flag. Phased charts now get Rule 1 by default, judged against each phase's limits, with the same points the chart flags; an explicit rule list runs as given, against each row's own centre. The limits and the flagged points on the chart were always correct. +- **`value='response'` charts the response**, as the `execute()` signature and docstring say. + Every chart refused it ("Valid charts for RESPONSE: None"): residuals are told apart by a leading + `R`, which `RESPONSE` also has. `execute()`, `supports_calibration()` and `why_not()` now treat + it, in any case, exactly as `value=None`. - **Xbar and S charts with subgroups of unequal size follow the manual (Eqs 11-5, 11-16, 11-17).** sigma-hat is the average of each subgroup's S_r / c4(N_r) (was S-bar / c4(N_r)), and each subgroup's limits use its own N_r, so they step with subgroup size, as VAS draws them. N_r now diff --git a/processbehavior/study.py b/processbehavior/study.py index fda6421..2f4a935 100644 --- a/processbehavior/study.py +++ b/processbehavior/study.py @@ -75,6 +75,17 @@ def _base_residual_code(value: str) -> str: return upper[2:] if upper.startswith('RCR') else upper +def _normalise_value(value: str | None) -> str | None: + """'response' (any case) -> None, the response's own spelling; residuals unchanged. + + Residuals are told apart by their leading 'R', which 'RESPONSE' also has, so the + documented ``value='response'`` must become None before anything tests the prefix. + """ + if isinstance(value, str) and value.lower() == 'response': + return None + return value + + @dataclass class DesignReport: """ @@ -1521,6 +1532,8 @@ def why_not(self, chart: str, value: str | None = None) -> str: >>> study.why_not('Xbar', value='R2') "'Xbar' with value='R2' is not valid for ADS 1.\\nValid charts for R2: S, X" """ + value = _normalise_value(value) + # Handle old fused syntax: why_not('R5_Xbar') -> why_not('Xbar', value='R5') if '_' in chart and value is None: parts = chart.split('_', 1) @@ -2018,6 +2031,8 @@ def execute( ) by = stratify_by + value = _normalise_value(value) + # 0. Resolve calibration (label string -> attached Calibration object). calibration_obj = self._resolve_calibration(calibration) @@ -2114,6 +2129,7 @@ def supports_calibration( 'parameter.' ) by = stratify_by + value = _normalise_value(value) from .analysis import Analysis from .formulation_spec import ChartRequest @@ -2396,7 +2412,10 @@ def _validate_recentered(self, value: str | None) -> None: ) recenterable = set(RESIDUAL_CODES) if value is None: - raise ValidationError('recentered=True requires a residual value (R1-R6), got value=None') + raise ValidationError( + 'recentered=True requires a residual value (R1-R6); the response ' + '(value=None or "response") is already on the measurement scale.' + ) if value.upper() not in recenterable: raise ValidationError(f"recentered=True requires a residual value (R1-R6), got '{value}'") diff --git a/tests/test_execute_validation.py b/tests/test_execute_validation.py index 054d113..a54591c 100644 --- a/tests/test_execute_validation.py +++ b/tests/test_execute_validation.py @@ -218,3 +218,54 @@ def test_execute_with_accessor(self, sds1_study): """End-to-end: study.execute(chart=..., value=study.residuals.R5).""" result = sds1_study.execute(chart=sds1_study.charts.Xbar, value=sds1_study.residuals.R5) assert result is not None + + +# ============================================================================= +# value='response' is the response (#142) +# ============================================================================= + + +class TestValueResponse: + """The documented value='response' charts the response, like value=None. + + Residuals are told apart by a leading 'R', which 'RESPONSE' also has, so the + documented spelling was routed down the residual path and refused. + """ + + @pytest.mark.parametrize( + 'kwargs', + [ + {'chart': 'Xbar'}, + {'chart': 'S'}, + {'chart': 'Xbar', 'companion': True}, + {'chart': 'X', 'by': []}, + {'chart': 'mR', 'by': []}, + {'chart': 'Histogram'}, + ], + ids=['Xbar', 'S', 'Xbar+S', 'X', 'mR', 'Histogram'], + ) + @pytest.mark.parametrize('spelling', ['response', 'Response', 'RESPONSE']) + def test_same_result_as_none(self, sds1_study, kwargs, spelling): + from pandas.testing import assert_frame_equal + + expected = sds1_study.execute(**kwargs) + actual = sds1_study.execute(value=spelling, **kwargs) + assert set(actual.charts) == set(expected.charts) + for name in expected.charts: + assert_frame_equal(actual.get_chart(name), expected.get_chart(name)) + assert actual.get_statistics(name) == expected.get_statistics(name) + assert actual.charts[name]['metadata'] == expected.charts[name]['metadata'] + + def test_why_not_agrees(self, sds1_study): + for chart in ('Xbar', 'S', 'X', 'mR'): + assert sds1_study.why_not(chart, value='response') == sds1_study.why_not(chart) + + def test_supports_calibration_agrees(self, sds1_study): + for kwargs in ({'chart': 'X', 'by': []}, {'chart': 'Xbar'}, {'chart': 'X', 'by': ['factor 1']}): + assert sds1_study.supports_calibration(value='response', **kwargs) == ( + sds1_study.supports_calibration(**kwargs) + ) + + def test_recentered_response_still_refused(self, sds1_study): + with pytest.raises(ValidationError, match='recentered.*requires.*residual'): + sds1_study.execute(chart='Xbar', value='response', recentered=True)