Sensory-to-chemistry enhancements: centring/scaling traps, panel diagnostics, chemistry preprocessing, permutation nulls - #532
Merged
Conversation
Placeholder commit so the pull request exists before the work lands. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01HvMyPfqTG7bLKTm5zr2vPF
Codecov Report❌ Patch coverage is
📢 Thoughts on this report? Let us know! |
Two conditions in the PLS API produce a plausible wrong answer rather than an error, and both are now announced. A1: PLS(scale=False) centres nothing and fits no intercept, so a response left on its natural scale displaces every prediction by roughly the response mean. Adding a constant to the response of an otherwise perfect fit moves in-sample R2 from +0.99 to -4.12, which a caller reasonably reads as "no relationship in this data". fit() now raises a SpecificationWarning naming the offending columns, the block, and the symptom. It does not centre for you: scale=False means "touch nothing", and centring silently would move the numbers for every caller who already centres correctly (a breaking change, not a quiet fix). A3: select_n_components(scale_inside_folds=True) re-fits MCUVScaler inside every training fold, which overwrites a scaling the caller applied on purpose. Two deliberately different strategies then report RMSECV identical to six decimal places, so a comparison between them shows no difference for reasons unrelated to the data. The existing warning fires on scale_inside_folds=False, which is the wrong way round for pre-scaled input; the True branch now warns when X arrives already centred and unit-variance scaled. Docstrings, no behaviour change: - A2: center() returns the value it subtracted, scale() returns the multiplier it applied (the reciprocal), and scale() defaults to ddof=0 where MCUVScaler uses ddof=1. Both docstrings now state the asymmetry and cross-reference. - A4: spe_ and hotellings_t2_ are (n_samples, n_components), one column per component. np.asarray(m.spe_).ravel() is right at one component and silently yields n * A values above it; use .iloc[:, -1]. - A6: sum(VIP ** 2) is exactly the number of X variables, so a VIP-exceedance count describes the shape of the VIP distribution rather than the presence of a relationship, and has almost no power as a null statistic. - A7 does not reproduce (fits at, and beyond, as many components as there are variables all succeed), so nothing functional changes and no component cap is added. safe_inverse's message now lists the likely causes. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01HvMyPfqTG7bLKTm5zr2vPF
mixed_assessor_model on a panel with no rows returned frames with no columns at all, not merely no rows, so the caller met the problem as a KeyError on ftests["f_product_mam"] rather than as a diagnosable error at the source. An empty panel is reachable whenever an upstream filter removes every attribute, so it now raises a ValueError naming that condition and pointing at the filter. Missing schema columns are checked in the same place. New module process_improve.sensory.diagnostics, three preconditions worth testing before any attribute is modelled as an intensity: boundary_occupancy reports floor, ceiling and exact-zero occupancy per attribute, separately. An attribute pinned against a scale bound violates the Mixed Assessor Model's premise that assessors compress or expand a linear scale, because no scaling difference is expressible where everyone records the same value. The three are separate columns because they are separate questions: a panel whose convention is to record "not perceived" as a small positive number looks floor-pinned but has no exact zeros at all. detection_rate gives the response that suits such an attribute, a product-by- attribute probability of detection rather than an intensity. The docstring says plainly that it is not comparable with intensity scores. A product-attribute pair nobody assessed comes back NaN, not 0: never detected and never asked are different answers. assessor_variance_equality tests the precondition under which Grossmann et al. (2023, Food Quality and Preference, doi:10.1016/j.foodqual.2022.104792) show the Mixed Assessor Model misreads unequal assessor variance as a scaling effect, shifting its F-test so disagreement is understated. Residuals are taken against the product mean first, so genuine product effects are removed before Levene's test (median-centred) compares assessors. A spread ratio accompanies the p-value as an effect size. Tests: the A5 empty-panel KeyError, and constructed data with one deliberately noisy assessor where the check fires for the unequal attribute and stays quiet for a matched equal one. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01HvMyPfqTG7bLKTm5zr2vPF
New subpackage process_improve.chemistry, for getting a block of concentrations or peak areas ready to relate to a block of sensory attributes. The order is fixed throughout: trim, then transform, then centre, then scale. classify_zero_states records what a zero in each column is taken to mean, and defaults to "unknown". Calling a zero left-censored asserts a latent value below a detection limit; calling it structurally absent asserts the compound is not there. Those need opposite handling and the distinction cannot be recovered from an exported table, so the caller declares it. Passing a detection limit is that declaration, and classifies the compound as rounded; nothing else does. trim_by_prevalence splits on the non-zero count and returns a presence layer covering every compound, kept and dropped alike. A compound seen in three of forty products has no concentration worth modelling but its binary fingerprint may carry more than a column full of zeros ever could. A missing measurement stays NaN in that layer: it is not an absence. normalisation_check reports the row totals and the ones outside a fold-band around the median, which is how an unrecorded dilution or a dropped peak shows up before it reaches the transform. choose_transform reads the range ratio of the detected values, returning "log", "linear", or an honest "ambiguous" that apply_transform resolves with a caller-chosen default. apply_transform substitutes non-detects (half the declared limit, or half the smallest value seen) before taking logs, leaving detected values exact; a linear compound passes through untouched, since a zero on a linear scale is a usable number. center_and_scale offers autoscale and Pareto, and takes the presence layer as a required argument so detected_only cannot be switched on without the mask to hand. It defaults to False: on zeros that were never imputed, excluding them puts every non-detect many standard deviations below a centre estimated from a handful of values, the column becomes a large-magnitude binary, and PLS then tracks how sparse a variable is rather than how it relates to the response, so every attribute returns the same rare compounds at the top of its list. apply_fitted_transform and apply_fitted_center_scale replay tables computed elsewhere. They are what make honest nested cross-validation possible: without them the transform offsets and the scaling constants would both have seen the row being predicted. The constants table names its column "divisor" and is divided by, so the multiplier-versus-divisor confusion has nowhere to hide. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01HvMyPfqTG7bLKTm5zr2vPF
New module process_improve.multivariate._null, re-exported alongside vip, center and scale. permutation_q2 is the replacement for a VIP-exceedance null (A6). It permutes the response across products, refits through a caller-supplied fit_predict, and compares observed out-of-sample performance with what a random reassignment reaches. Out-of-sample is what makes it bite: the same test on in-sample R2 would mostly measure model capacity. Whole rows are permuted, so each product keeps its full response vector and the correlation structure among attributes survives; permuting column-wise would inflate the null. The p-value uses the (1 + count) / (n_perm + 1) form so it can never be exactly zero, and the docstring states the consequence, that the smallest attainable p-value is 1 / (n_perm + 1), where a caller choosing n_perm will read it. Two things the docstring makes the caller's own: the cross-validation scheme (leave-one-out is right when every product is precious and wasteful at hundreds of rows, and a permutation null multiplies that cost by n_perm), and re-deriving the response centring and scaling inside each fold, which is the half that is easy to forget. pipeline_null generalises it to a whole selection procedure: pass a callable that runs filtering, transformation, scaling and selection end to end, and it counts discoveries under a permuted response. The response-independent steps are deliberately not hoisted out of the loop, since hoisting them would be an assumption about the caller's code. What is asserted instead is determinism: a selector that disagrees with itself on two identical calls makes the empirical FDR meaningless, so that draws a SpecificationWarning rather than a silently contaminated number. class_enrichment tests whether a named class of compounds is over-represented at the top of a ranking, hypergeometrically. At small sample sizes this is frequently stronger evidence than R2 or Q2: recovering the chemically expected class for an attribute is structure noise does not produce, whereas a high R2 on few products with several components nearly is. Its docstring carries the one-component caveat, that a single-component model's coefficient matrix is outer(x_weights, y_loadings) and therefore orders identically for every attribute, so an enrichment that looks attribute-specific may be one ranking reported many times. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01HvMyPfqTG7bLKTm5zr2vPF
…ty selection New module process_improve.interactions. Its docstring says plainly, in a warning block, that all three functions are unvalidated on real data: the sample-size conditions they exist to test are rarely met, so this code may go unreached for a long time after shipping and its coverage is unit tests only. Shipping them as public API commits to maintaining code whose behaviour on real inputs is unverified, and a reader deserves to know that before quoting a result from it. pair_coverage asks whether the observations populate all four corners of a pair's plane, splitting each variable at its own median so only the joint distribution can fail. Two variables that co-vary occupy the agreeing corners only and fail the check; the returned detail carries their correlation, because that is nearly always the explanation. A failure here is the correct answer rather than a defect to work around: the interaction is not identifiable from those observations and no regularisation makes it so. interaction_terms builds the products in the one order that is defensible: transform, centre and scale, multiply, then re-centre and re-scale the products. A product of two standardised columns is not itself centred; for approximately bivariate normal parents its mean is the correlation r and its variance is 1 + r squared, so skipping the second pass leaks correlation into the intercept and inflates exactly the pairs whose interactions deserve least trust. Parents that do not look standardised draw a SpecificationWarning, since both that reasoning and the reported parent_correlation depend on it. stability_selection reports selection frequency across complementary half-samples: each split is used in both directions, so every product appears in exactly half the subsamples and the two runs of a split share no rows. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01HvMyPfqTG7bLKTm5zr2vPF
One page per new module, in the user-guide toctree: - chemistry.rst: the fixed order (trim, transform, centre, scale), why a zero is not self-describing, why a trimmed compound comes back as a presence layer, why detected_only is off by default, and how the apply_fitted_* pair make a nested fold honest. - sensory_diagnostics.rst: whether an attribute can be modelled as an intensity at all, the floor-pinned-but-not-zero case that only the exact_zero column separates, and the Grossmann et al. (2023) result that assessor_variance_ equality tests. - permutation_nulls.rst: why a high in-sample R2 and a VIP-exceedance count are both unusable as evidence, what permutation_q2 asks instead, the p-value floor of 1 / (n_perm + 1), and the one-component caveat on any per-attribute ranking. - interactions.rst: leads with the provisional warning, then pair coverage, the non-negotiable re-centring of product terms, and stability selection. Autodoc for the new modules lives on these pages rather than under docs/api, so each object is described exactly once and the strict build stays clean. Version bumped 1.75.2 -> 1.76.0 (MINOR: new modules and API additions), with CITATION.cff kept in step and the CHANGELOG section and link footer updated. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01HvMyPfqTG7bLKTm5zr2vPF
Codecov flagged 36 lines on the patch. Most were guards this PR documents in a Raises section but never exercised; two were real behaviour. The behaviour: permutation_q2 on a constant response column, and on predictions that are entirely missing, now has its NaN path pinned rather than inferred. _looks_prescaled, which the A3 warning is built on, gets direct tests for the four cases that matter: a standardised block, a block scaled with ddof=0 (off by sqrt(n / (n - 1)), inside the tolerance on purpose), a raw block, and an all-constant block, where no column carries evidence either way and the answer is False rather than a vacuous True. The guards: non-DataFrame input, a block with no compounds, a block of unreadable dtypes (distinguished from an all-NaN block, which is legitimate), a non-DataFrame detected mask, constants and applied tables missing a required column, mismatched row counts, and n_perm below 1. Every new module is now fully covered except one line: _log_offset's "nothing was ever detected" branch, which cannot be reached through apply_transform because the log rule requires two detected values to form its range ratio. Rather than contort a test to reach it, the comment there now says so, and says why the branch is kept. 2967 passed, coverage 94.45%. Co-Authored-By: Claude Opus 5 <noreply@anthropic.com> Claude-Session: https://claude.ai/code/session_01HvMyPfqTG7bLKTm5zr2vPF
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
2 participants
Add this suggestion to a batch that can be applied as a single commit.This suggestion is invalid because no changes were made to the code.Suggestions cannot be applied while the pull request is closed.Suggestions cannot be applied while viewing a subset of changes.Only one suggestion per line can be applied in a batch.Add this suggestion to a batch that can be applied as a single commit.Applying suggestions on deleted lines is not supported.You must change the existing code in this line in order to create a valid suggestion.Outdated suggestions cannot be applied.This suggestion has been applied or marked resolved.Suggestions cannot be applied from pending reviews.Suggestions cannot be applied on multi-line comments.Suggestions cannot be applied while the pull request is queued to merge.Suggestion cannot be applied right now. Please check back later.
Summary
fix(multivariate)- the two centring/scaling conditions that produce a plausible wrong answer rather than an error are now announced.PLS(scale=False)fits no intercept, so an un-centred response displaces every prediction by roughly the response mean (in-sample R2 goes from +0.99 to -4.12 on an otherwise perfect fit);select_n_components(scale_inside_folds=True)re-standardises inside every fold and erases a scaling the caller chose, so two deliberately different strategies report RMSECV identical to six decimal places. Neither is fixed silently:scale=Falsemeans "touch nothing", and auto-centring would move the numbers for every caller who already centres correctly. Docstrings also corrected for thecenter/scaleasymmetry (subtrahend vs multiplier,ddof=0vsddof=1), the per-component shape ofspe_/hotellings_t2_, and thesum(VIP**2) = Kidentity that makes a VIP-exceedance count useless as a null statistic. A7 does not reproduce, so nothing functional changed there and no component cap was added.feat(sensory)-mixed_assessor_modelon an empty panel raisedKeyError: 'f_product_mam'somewhere downstream; it now raises aValueErrornaming the condition. Newsensory.diagnostics:boundary_occupancy(floor, ceiling and exact-zero counted separately, so a panel recording "not perceived" as a small positive number is not mistaken for floor-pinned),detection_rate, andassessor_variance_equality(the Grossmann et al. 2023 precondition, Levene on product-mean residuals).feat(chemistry)- new subpackage: trim, transform, centre, scale, in that order, withapply_fitted_*partners so a held-out fold is preprocessed with training constants alone. Zero states default tounknown, never censored. Trimmed compounds return as a presence layer covering every compound.detected_onlydefaults toFalse.feat(multivariate)-permutation_q2,pipeline_null,class_enrichment: nulls that respond to signal. Whole response rows are permuted; the p-value floor of1 / (n_perm + 1)is stated where a caller choosesn_perm.feat(interactions)-pair_coverage,interaction_terms,stability_selection, landing last and marked provisional in the module docstring: unvalidated on real data, unit-test coverage only.Version bumped 1.75.2 -> 1.76.0 (MINOR),
CITATION.cffkept in step,CHANGELOG.mdupdated, and a user-guide page added per new module.Test plan
pytestpasses with the project's default options: 2967 passed, 6 skipped, coverage 94.45% against the 92% gate. The 6 skips are pre-existing (network-dependent openmv.net fetches and two known API gaps).python -O, thetest-under-dash-OCI gate.ruff check .andruff format --check .clean;mypy src/process_improveclean.sum(VIP**2) = 5.0000and the empty-panelKeyError). A7 did not reproduce: fits at 3, 4 and 5 components on a 3-variable block all succeeded, as didvip,cross_validateandselect_n_componentson the result. No component cap was added.KeyErrorin place of the expectedValueError.test_multivariate_centring_traps.py, 26 intest_sensory_diagnostics.py, 64 intest_chemistry.py, 37 intest_multivariate_null.py, 31 intest_interactions.py. Tests over 2 s carry theslowmarker. Every new module is fully covered bar one line:_log_offset's "nothing was ever detected" branch, unreachable throughapply_transformbecause the log rule requires two detected values to form its range ratio; the comment there says so, rather than a test contorting to reach it.-W, notebooks executed) is green with the four new pages. Locally I could only build with the notebook directories excluded (nopandocin the sandbox, and the intersphinx inventories are proxy-blocked); that partial build showed no warnings attributable to the new pages, and all four render their autodoc content (11, 3, 3 and 3 documented objects). CI confirms the full build.Known-good downstream regression (acceptance criterion 4)
Two analyses whose answers are established independently of this package were run at
origin/mainand at this branch, through a script printing 10 significant figures:datasets/multivariate/LDPE/). Scores, loadings, weights, T2 at A=6, back-transformed predictions, cumulative R2, RMSE, VIP, andselect_n_components.The output is byte-identical before and after: nothing moved. That is the expected result, since the change set adds warnings, docstrings and new modules and touches no fitting path. The reproduced agreement figures were, on both sides: published Wold loadings
P1 = [0.5410, 0.3493, 0.5410, 0.5410]and R2 per component[0.8307, 0.1693]; LDPE versus SIMCA, max absolute deviation 9.1e-05 on scores, 2.8e-05 on loadings, 2.9e-05 on weights, 3.1e-04 on T2, and 3.7e-05 relative on predictions.sum(VIP**2)is 14.0000000000 for K=14, pinning the A6 identity on real data as well as synthetic.Checklist
pyproject.toml(MINOR: new modules and API additions), withCITATION.cffset to the same version and date in the same commitruff check .passesCHANGELOG.mdupdated