Skip to content

MBPCA / MBPLS hardcode default_rng(0) in fit() and expose no random_state #503

Description

@kgdunn

Summary

MBPCA.fit() and MBPLS.fit() seed their NIPALS super-score initialisation from a hard-coded np.random.default_rng(0), and neither estimator accepts a random_state. This violates the package's own written contract in docs/development/reproducibility.rst, which forbids literal seeds in production code and requires every public function touching an RNG to accept random_state.

Where

Line numbers as of ac52aab (branch claude/repo-error-audit-si62bq, PR #500).

  • src/process_improve/multivariate/_mbpca.py:311 - rng = np.random.default_rng(0)
  • src/process_improve/multivariate/_mbpls.py:380 - rng = np.random.default_rng(0)
  • Contract: docs/development/reproducibility.rst
  • Helper to use: src/process_improve/_random.py (check_random_state), already used correctly by _resampling.py

Evidence

# src/process_improve/multivariate/_mbpca.py:311
rng = np.random.default_rng(0)
# src/process_improve/multivariate/_mbpls.py:380
rng = np.random.default_rng(0)

Neither MBPCA.__init__ nor MBPLS.__init__ takes random_state (grep random_state src/process_improve/multivariate/_mbpca.py returns nothing).

Why it is wrong

  1. It breaks the documented contract. docs/development/reproducibility.rst states that hard-coded literal seeds in production code are forbidden, and carves out an exception only for cosmetic things like plot jitter. These are fit paths, not jitter.
  2. It cannot be audited. Because the seed is fixed and unreachable, a user cannot re-run a fit from a different start to check whether NIPALS landed in the same fixed point. On near-degenerate blocks a different start genuinely can converge elsewhere - that is exactly the check the contract exists to make possible.
  3. A random start is a poor initialiser anyway. The single-block PLS path deliberately seeds from the highest-variance Y column instead (start_col = int(np.argmax(start_SSY_col)) in _pls.py), because it is closer to the leading component, needs fewer iterations, and is deterministic without a seed at all. See Multivariate bug fixes and numerical refinements #195 for that change.

Suggested fix

Preferred: remove the RNG entirely and seed the super-score deterministically the way _pls.py does - from the highest sum-of-squares column/block. That satisfies the contract by having no RNG at all, converges faster, and needs no new constructor parameter.

If a random start is genuinely wanted, add random_state: int | np.random.Generator | None = None to both constructors, resolve it once at the top of fit() via check_random_state from process_improve._random, and document it per the reproducibility contract.

Either way, the sign convention already applied after convergence makes the fitted signs independent of the start.

Acceptance criteria

  • No literal seed remains in _mbpca.py / _mbpls.py fit paths.
  • Either no RNG is used, or random_state is a constructor parameter resolved through process_improve._random.check_random_state.
  • Test: two fits on the same data give identical loadings/scores (determinism), and - if random_state was added - two different seeds give the same converged component up to sign on well-conditioned data.
  • docs/development/reproducibility.rst needs no exception added for these estimators.

Scope / non-goals

_mbpls.py:1050 (rng = np.random.default_rng(seed) in the randomization test) already takes a seed from the caller - leave it alone.

References

Found during the repo-wide audit in PR #500, deferred there because it is an API addition rather than a bug fix.

Metadata

Metadata

Assignees

No one assigned

    Labels

    No labels
    No labels

    Projects

    No projects

    Milestone

    No milestone

    Relationships

    None yet

    Development

    No branches or pull requests

    Issue actions