diff --git a/CHANGELOG.md b/CHANGELOG.md index a4a9826..647db9a 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -11,6 +11,23 @@ those changes. ## [Unreleased] +## [1.59.0] - 2026-07-23 + +### Added + +- `process_improve.sensory.panel_consistency`: a per-panelist reliability test that + compares each panelist's product signal against a reshuffle of their own scores. For + every panelist the product labels are permuted within each attribute (keeping the + panelist's own values, only reassigning which product got which) and the mean + eta-squared product effect and the agreement with the leave-one-out panel consensus are + recomputed; the null is built entirely from the panelist's own data, no simulation. A + panelist who carries a real, reproducible signal beats almost every reshuffle (small + p-value), while a noise rater is indistinguishable from its own permutations. Returns a + `PanelConsistency` with per-panelist `discrimination` / `p_discrimination`, `agreement` + / `p_agreement`, and a `consistent` flag (gated on the discrimination p-value; agreement + is reported but not gated, since a reproducible rater may still rank the products + differently from the panel). + ## [1.58.0] - 2026-07-23 ### Added @@ -2599,7 +2616,8 @@ this entry records them together. - Reworked the README with a sharper value proposition and a "Why not scikit-learn?" comparison table. -[Unreleased]: https://github.com/kgdunn/process-improve/compare/v1.58.0...HEAD +[Unreleased]: https://github.com/kgdunn/process-improve/compare/v1.59.0...HEAD +[1.59.0]: https://github.com/kgdunn/process-improve/compare/v1.58.0...v1.59.0 [1.58.0]: https://github.com/kgdunn/process-improve/compare/v1.57.0...v1.58.0 [1.57.0]: https://github.com/kgdunn/process-improve/compare/v1.56.0...v1.57.0 [1.56.0]: https://github.com/kgdunn/process-improve/compare/v1.55.1...v1.56.0 diff --git a/CITATION.cff b/CITATION.cff index a544588..fca9fa5 100644 --- a/CITATION.cff +++ b/CITATION.cff @@ -12,7 +12,7 @@ authors: repository-code: "https://github.com/kgdunn/process-improve" url: "https://kgdunn.github.io/process-improve/" license: MIT -version: 1.58.0 +version: 1.59.0 date-released: "2026-07-23" keywords: - chemometrics diff --git a/pyproject.toml b/pyproject.toml index 3996999..7fc5fc3 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -1,6 +1,6 @@ [project] name = "process-improve" -version = "1.58.0" +version = "1.59.0" description = 'Designed Experiments; Latent Variables (PCA, PLS, multivariate methods with missing data); Process Monitoring; Batch data analysis.' readme = "README.md" license = "MIT" diff --git a/src/process_improve/sensory/__init__.py b/src/process_improve/sensory/__init__.py index 63443bb..4ed7272 100644 --- a/src/process_improve/sensory/__init__.py +++ b/src/process_improve/sensory/__init__.py @@ -27,7 +27,13 @@ ) from process_improve.sensory.ingest import reshape_to_long from process_improve.sensory.mam import MAMResult, align_scores, mixed_assessor_model -from process_improve.sensory.panel import PanelScorecard, apply_correction, panel_scorecard +from process_improve.sensory.panel import ( + PanelConsistency, + PanelScorecard, + apply_correction, + panel_consistency, + panel_scorecard, +) from process_improve.sensory.recipes import SENSORY_RECIPES from process_improve.sensory.validation import ( DESCRIPTIVE_LONG_COLUMNS, @@ -40,6 +46,7 @@ "SENSORY_RECIPES", "AnalysisResult", "MAMResult", + "PanelConsistency", "PanelScorecard", "ValidationResult", "aggregate_to_product", @@ -48,6 +55,7 @@ "apply_correction", "discriminate_observational", "mixed_assessor_model", + "panel_consistency", "panel_scorecard", "permutation_column_null", "product_means", diff --git a/src/process_improve/sensory/panel.py b/src/process_improve/sensory/panel.py index bbc1b78..d243833 100644 --- a/src/process_improve/sensory/panel.py +++ b/src/process_improve/sensory/panel.py @@ -197,6 +197,213 @@ def panel_scorecard(panel: pd.DataFrame) -> PanelScorecard: return PanelScorecard(table=table, flagged=sorted(reasons), reasons=reasons) +#: Minimum products before a panelist's own-permutation null is meaningful. +_MIN_PRODUCTS_FOR_PERMUTATION = 3 + + +@dataclass +class PanelConsistency: + """Outcome of :func:`panel_consistency`. + + Attributes + ---------- + table : pandas.DataFrame + One row per panelist, indexed by ``panelist_id``, with the columns + ``discrimination``, ``p_discrimination``, ``agreement``, ``p_agreement``, + and ``consistent``. + inconsistent : list of str + Panelist ids whose product signal is not distinguishable from a reshuffle of + their own scores (``consistent`` is ``False``). + n_permutations : int + Number of within-panelist permutations used for each null. + """ + + table: pd.DataFrame + inconsistent: list[str] + n_permutations: int + + +def _eta_squared_codes(scores: np.ndarray, codes: np.ndarray, n_groups: int) -> float: + """Return the eta-squared of an integer-``codes`` grouping of ``scores``.""" + grand = float(scores.mean()) + ss_total = float(((scores - grand) ** 2).sum()) + if ss_total <= 0.0: + return float("nan") + sums = np.bincount(codes, weights=scores, minlength=n_groups) + counts = np.bincount(codes, minlength=n_groups).astype(float) + means = np.divide(sums, counts, out=np.full(n_groups, grand), where=counts > 0) + ss_between = float((counts * (means - grand) ** 2).sum()) + return ss_between / ss_total + + +def _group_means(scores: np.ndarray, codes: np.ndarray, n_groups: int) -> np.ndarray: + """Return the per-group mean of ``scores`` (NaN for an empty group).""" + sums = np.bincount(codes, weights=scores, minlength=n_groups) + counts = np.bincount(codes, minlength=n_groups).astype(float) + return np.divide(sums, counts, out=np.full(n_groups, np.nan), where=counts > 0) + + +def _correlation(a: np.ndarray, b: np.ndarray) -> float: + """Return the Pearson correlation of two aligned vectors, or NaN if undefined.""" + if a.size < 2 or a.std() <= 0.0 or b.std() <= 0.0: + return float("nan") + return float(np.corrcoef(a, b)[0, 1]) + + +def _consistency_pvalues( # noqa: PLR0913 + blocks: list[tuple[str, np.ndarray, np.ndarray, int]], + locators: list[tuple[int, int]], + consensus: np.ndarray, + real: tuple[float, float], + rng: np.random.Generator, + n_permutations: int, +) -> tuple[float, float]: + """Return one-sided permutation p-values for discrimination and agreement. + + Each permutation reshuffles the product labels within every attribute (breaking the + panelist's product-to-score link while keeping their own score multiset), then + recomputes the mean eta-squared and the agreement with the fixed leave-one-out + consensus. The p-value is the add-one fraction of permutations reaching at least the + real statistic (higher is better for both). + """ + disc_real, agree_real = real + ge_disc = 0 + ge_agree = 0 + for _ in range(n_permutations): + # Permute the product labels once per attribute block and reuse that single + # reshuffle for both the eta-squared and the agreement statistic. + etas = [] + means_by_block = [] + for _attr, scores, codes, ng in blocks: + pcodes = rng.permutation(codes) + etas.append(_eta_squared_codes(scores, pcodes, ng)) + means_by_block.append(_group_means(scores, pcodes, ng)) + null_disc = float(np.nanmean(etas)) if etas else float("nan") + if not np.isnan(null_disc) and null_disc >= disc_real: + ge_disc += 1 + if not np.isnan(agree_real): + null_own = np.array([means_by_block[bi][gi] for bi, gi in locators]) + null_agree = _correlation(null_own, consensus) + if not np.isnan(null_agree) and null_agree >= agree_real: + ge_agree += 1 + p_disc = (ge_disc + 1.0) / (n_permutations + 1.0) + p_agree = (ge_agree + 1.0) / (n_permutations + 1.0) if not np.isnan(agree_real) else float("nan") + return p_disc, p_agree + + +def panel_consistency( + panel: pd.DataFrame, + *, + n_permutations: int = 299, + random_state: int = 0, + alpha: float = 0.05, +) -> PanelConsistency: + """Test each panelist's product signal against a reshuffle of their own scores. + + For every panelist the product labels are permuted within each attribute (keeping the + panelist's own score values, only reassigning which product got which), and the mean + eta-squared product effect (:func:`_eta_squared_codes`) and the agreement with the + leave-one-out panel consensus are recomputed. A panelist who carries a real, + reproducible product signal beats almost every reshuffle of their own scores (small + p-value); one whose scores are noise is indistinguishable from their own permutations + (p near 1). The null is built entirely from the panelist's own data - no simulated or + parametric reference. + + ``consistent`` gates on the discrimination p-value (does the panelist separate the + products beyond their own noise). The agreement p-value is reported alongside but not + gated: a panelist can discriminate reproducibly yet rank the products differently from + the panel, which is disagreement rather than inconsistency. + + Parameters + ---------- + panel : pandas.DataFrame + Validated ``descriptive_long`` panel data. + n_permutations : int + Within-panelist permutations per null. + random_state : int + Seed for the permutations. + alpha : float + Significance level for the ``consistent`` flag. + + Returns + ------- + PanelConsistency + See the class docstring. + + Examples + -------- + >>> result = panel_consistency(validated.normalized_df, n_permutations=299) + >>> result.inconsistent + ['P8'] + """ + if n_permutations < 1: + raise ValueError(f"n_permutations must be at least 1, got {n_permutations}.") + + cell = panel.groupby(["panelist_id", "product", "attribute"], observed=True)["score"].mean().reset_index() + rng = np.random.default_rng(random_state) + records: dict[str, dict[str, float | bool]] = {} + for pid in panel["panelist_id"].unique(): + pan = panel[panel["panelist_id"] == pid] + blocks: list[tuple[str, np.ndarray, np.ndarray, int]] = [] + block_labels: list[np.ndarray] = [] + etas_real: list[float] = [] + for attr, grp in pan.groupby("attribute", observed=True): + scores = grp["score"].to_numpy(dtype=float) + valid = ~np.isnan(scores) + labels, codes = np.unique(grp["product"].to_numpy()[valid], return_inverse=True) + if valid.sum() < 2 or labels.size < 2: + continue + blocks.append((str(attr), scores[valid], codes, labels.size)) + block_labels.append(labels) + etas_real.append(_eta_squared_codes(scores[valid], codes, labels.size)) + + # Agreement against the leave-one-out consensus, over the (product, attribute) + # cells this panelist and the rest of the panel share. + loo = cell[cell["panelist_id"] != pid].groupby(["product", "attribute"], observed=True)["score"].mean() + loo_map = {(str(p), str(a)): float(v) for (p, a), v in loo.items()} + locators: list[tuple[int, int]] = [] + consensus: list[float] = [] + real_own: list[float] = [] + for bi, (attr, scores, codes, ng) in enumerate(blocks): + means = _group_means(scores, codes, ng) + labels = block_labels[bi] + for gi in range(ng): + key = (str(labels[gi]), attr) + if key in loo_map and not np.isnan(means[gi]): + locators.append((bi, gi)) + consensus.append(loo_map[key]) + real_own.append(float(means[gi])) + finite_etas = [e for e in etas_real if not np.isnan(e)] + discrimination = float(np.mean(finite_etas)) if finite_etas else float("nan") + agreement = _correlation(np.array(real_own), np.array(consensus)) if len(real_own) >= 2 else float("nan") + + if not blocks or pan["product"].nunique() < _MIN_PRODUCTS_FOR_PERMUTATION or np.isnan(discrimination): + records[str(pid)] = { + "discrimination": discrimination, + "p_discrimination": float("nan"), + "agreement": agreement, + "p_agreement": float("nan"), + "consistent": False, + } + continue + + p_disc, p_agree = _consistency_pvalues( + blocks, locators, np.array(consensus), (discrimination, agreement), rng, n_permutations + ) + records[str(pid)] = { + "discrimination": discrimination, + "p_discrimination": p_disc, + "agreement": agreement, + "p_agreement": p_agree, + "consistent": bool(p_disc <= alpha), + } + + table = pd.DataFrame.from_dict(records, orient="index") + table.index.name = "panelist_id" + inconsistent = sorted(str(pid) for pid in table.index if not bool(table.loc[pid, "consistent"])) + return PanelConsistency(table=table, inconsistent=inconsistent, n_permutations=n_permutations) + + def apply_correction(panel: pd.DataFrame, drop: list[str]) -> pd.DataFrame: """Return ``panel`` with the listed panelists removed. diff --git a/tests/test_sensory.py b/tests/test_sensory.py index 35af93e..56a2776 100644 --- a/tests/test_sensory.py +++ b/tests/test_sensory.py @@ -15,6 +15,7 @@ align_scores, analyze_descriptive, mixed_assessor_model, + panel_consistency, panel_scorecard, validate_descriptive, ) @@ -271,6 +272,77 @@ def test_scorecard_clean_panel_has_no_flags(): assert card.flagged == [] +def test_consistency_flags_the_noise_rater(): + """A panelist scoring at random is indistinguishable from a reshuffle of its own scores.""" + result = panel_consistency(_panel(anomalous="P8"), n_permutations=199, random_state=0) + table = result.table + assert {"discrimination", "p_discrimination", "agreement", "p_agreement", "consistent"} == set(table.columns) + + # The random rater cannot beat its own permutations: high discrimination p-value, + # not consistent, and the least consistent of the whole panel. + assert "P8" in result.inconsistent + assert table.loc["P8", "p_discrimination"] > 0.5 + assert table["p_discrimination"].idxmax() == "P8" + assert table.loc["P8", "agreement"] < 0.5 + + # The strongest genuine raters clear their own-permutation null. + assert table.loc["P5", "consistent"] + assert table.loc["P4", "consistent"] + + +def _row(panelist: str, product: str, attribute: str, replicate: int, score: float) -> dict: + """Build a single descriptive_long row.""" + return { + "panelist_id": panelist, + "session": 1, + "product": product, + "attribute": attribute, + "replicate": replicate, + "score": score, + } + + +def test_consistency_handles_degenerate_panelists(): + """A flat rater and a single-product rater cannot be assessed and are not consistent.""" + products = list("WXYZ") + truth = {"W": -1, "X": 0, "Y": 1, "Z": 2} + rows = [_row("good", prod, "A", rep, 5.0 + 2.0 * truth[prod]) for prod in products for rep in (1, 2)] + rows += [_row("flat", prod, "A", rep, 5.0) for prod in products for rep in (1, 2)] # no variance + rows += [_row("one", "W", "A", rep, 6.0) for rep in (1, 2)] # only one product + result = panel_consistency(pd.DataFrame(rows), n_permutations=50, random_state=0) + table = result.table + + assert np.isnan(table.loc["flat", "p_discrimination"]) + assert np.isnan(table.loc["flat", "agreement"]) + assert not bool(table.loc["flat", "consistent"]) + assert not bool(table.loc["one", "consistent"]) + assert {"flat", "one"}.issubset(result.inconsistent) + + +def test_consistency_rater_without_panel_overlap_has_no_agreement(): + """A rater discriminating over products the rest of the panel never rated gets no agreement.""" + truth = {"W": -1, "X": 0, "Y": 1, "Z": 2} + rows = [_row("solo", p, "A", rep, 5.0 + 2.0 * truth[p]) for p in "WXYZ" for rep in (1, 2)] + # Two other raters cover only product W, so ``solo`` shares a single cell with the panel. + rows += [_row("b1", "W", "A", rep, 5.0) for rep in (1, 2)] + rows += [_row("b2", "W", "A", rep, 5.2) for rep in (1, 2)] + row = panel_consistency(pd.DataFrame(rows), n_permutations=30, random_state=0).table.loc["solo"] + assert not np.isnan(row["p_discrimination"]) # solo still gets a discrimination p-value + assert np.isnan(row["agreement"]) # but no agreement can be formed against the panel + assert np.isnan(row["p_agreement"]) + + +def test_consistency_is_deterministic_and_validates(): + """Same seed gives identical output; a non-positive permutation count is rejected.""" + panel = _panel(anomalous="P8") + first = panel_consistency(panel, n_permutations=99, random_state=1) + second = panel_consistency(panel, n_permutations=99, random_state=1) + assert first.table.equals(second.table) + assert first.inconsistent == second.inconsistent + with pytest.raises(ValueError, match="n_permutations must be at least 1"): + panel_consistency(panel, n_permutations=0) + + def test_dropping_panelist_changes_means(): validated = validate_descriptive(_panel(), _obs(), mode="observational") kept = analyze_descriptive(validated, drop_panelists=None, discriminator=False)