Validation

Two suites assert the published claims of the source papers. That is what separates software that runs from software you can cite.

just repro     # or: pytest tests/test_reproduce_carrera2004.py tests/test_reproduce_tubau2014.py

Carrera et al. (2004)

Their Application 2 — three end-members, five species, noisy end-member analyses and precise mixtures — is regenerated from the published tables and swept over sample count and noise level. The assertions:

  • least-squares correlations fall in the published 0.97–0.98 band at low end-member noise and 0.93–0.94 at high noise;

  • least squares does not improve with sample count, because it estimates each sample independently — the paper’s central diagnostic claim;

  • maximum likelihood beats it at every sample count and noise level tested;

  • maximum likelihood is insensitive to end-member noise where least squares is not;

  • the improvement index lands in the published 3–8 range.

The implementation is also cross-checked against the paper’s own formulation: profile_objective() computes the eliminated objective of their equation 21, and the suite asserts it agrees with what the solver minimises to eight decimal places.

Tubau et al. (2014)

The published Besòs end-member table and standard-deviation scheme are bundled and checked, including that a species carried at 1000× σ genuinely does not move the ratios.

The suite also reproduces a qualitative claim as a number. The paper states that its two dry-period end-members are “distinguished mainly by high ammonium and low calcium and magnesium”. Run identifiability_report() on the four conservative species alone and those end-members show condition number 159 and relative separation 0.105 — a failure by the library’s own thresholds. Adding NH₄, Ca and Mg drops the condition number to 12. The diagnostics rediscover the paper’s finding without being told.

Note

One inconsistency in the source is recorded rather than smoothed over: the paper’s prose says the dry end-members dominate the wet one “4:1”, but its own percentages (74% against 26%) give 2.8:1. The library follows the percentages.

The solvers themselves

The two suites above check that the library reproduces published answers. They cannot check that the optimisers find their own optimum, because a solver that stops early agrees with itself. tests/test_solver_optimality.py checks that separately, against references that share no code with the implementation.

Constrained least squares is checked against an exhaustive enumeration of the faces of the simplex — wasteful, but exact for the single-digit ne this field uses — over random problems at 2 to 6 end-members, with the samples deliberately placed outside the mixing hull so the non-negativity constraints activate. Every returned solution is also required to satisfy the KKT conditions directly.

This found a real bug. The previous active set zeroed the most negative ratio and re-solved, never reconsidering. That is a greedy heuristic, not an active-set method: for a sample outside the hull it can pin the wrong end-member out and stop on the wrong face of the simplex. It affected roughly 0.3% of such samples at four end-members and 1.2% at six, and it was invisible to every one of the 298 pre-existing tests. The worst case found reported a sample as pure end-member 4 where the optimum was a three-way mixture fitting 37 times better. The solver now releases constraints as well as adding them, which for a convex objective makes the returned point the global optimum; the failing case is kept as an explicit regression test.

Maximum likelihood is checked two ways: SciPy’s SLSQP is run on the paper’s eliminated objective (equation 21) from 25 starts on small problems and must not find a better optimum than the alternation, and a finite-difference KKT certificate confirms the converged fixed point is genuinely stationary. Both pass — the block coordinate descent does reach the optimum. Note that the default max_iter=2000 stops short of the default ratio_tol on Carrera’s own Application 2 (relative KKT residual 6e-4, falling to 4e-9 by 20 000 iterations). The ratios are settled far below analytical precision long before that, and the run reports converged=False with final_ratio_change rather than claiming success — but read those, as the docstring says, and raise max_iter if you want the optimum rather than a good answer.

Note

Mixing ratios are dimensionless, so changing concentration units may not move them. The suite asserts that across twenty orders of magnitude. This caught an absolute tolerance in the non-detect handling — a comparison between two concentrations — which has been made relative. The ratio tolerances elsewhere are dimensionless and stay absolute.

Everything else

319 tests in total: exactness on noise-free data, the simplex constraints, invariance to species order and to unit scaling, ragged and censored input, the isotope nonlinearity, and the figures. just check runs lint and the whole suite; just pyversions runs it on 3.11, 3.12 and 3.14.