perf: BLAS projection, per-width caching, chunked fits, and accuracy fixes for ExtremeLearningMachine

[CyrilJl]

Summary

Rewrites ExtremeLearningMachine for speed and memory, fixes a correctness bug that corrupted its output when used as the imputer regressor, and fixes a second, separate accuracy issue that made it systematically underperform (and sometimes actively hurt) compared to FastRidge in typical imputation workloads:

• BLAS projection: the naive numba triple-loop matmul kernel is replaced with np.matmul + in-place bias/ReLU on float32 (numba dependency removed from elm.py).
• Per-width projection cache (bug fix): fit previously initialized the projection matrix once and reused it forever. Inside MultivariateImputer, the same instance is refit on column subsets of varying widths (one per missingness pattern); refits with more input features read past the projection matrix (numba, no bounds check) and produced non-finite imputations, while narrower refits silently used a truncated projection. Projections are now sampled lazily per input width and cached, seeded identically — reproducible across instances, and identical values to before for the first width. set_params invalidates the cache when n_features/random_state change.
• Chunked large fits/predicts: above 65,536 rows, fit accumulates the augmented Gram matrix [H, y, 1] chunk by chunk and solves via the existing fit_ridge_from_gram; predict projects into a reusable chunk buffer. The full (n_samples, n_features) hidden matrix is never materialized.
• Fan-in scaling + sample-aware hidden width (accuracy fix): MultivariateImputer refits a fresh model per missingness pattern on a small, pattern-specific row subset (often tens of rows). Two compounding issues made ExtremeLearningMachine a poor fit for that regime: (1) projection weights/bias were drawn ~N(0,1) with no fan-in scaling, so pre-activation variance — and therefore the effective strength of alpha — grew with the number of predictor columns per pattern; (2) the hidden width (n_features, default 100) was fixed regardless of how many samples a given pattern's fit had, turning the internal ridge fit into a severely underdetermined (parameters > samples) problem. Weights/bias are now scaled by 1/sqrt(n_input_features), and a new min_samples_per_feature parameter (default 5) caps the hidden width to n_samples // min_samples_per_feature at fit time (set to 0 to opt out and keep the old fixed-width behavior).

Measurements

Same machine, warm runs (local perf/elm_benchmark.py):

Scenario	Before	After
Standalone fit, 20k rows × 35 → 100 components	8.7 ms	3.3 ms
Standalone predict, 5k rows	1.9 ms	0.7 ms
MultivariateImputer(regressor=ELM), 10k × 60 panel, 3 targets	0.67 s, non-finite output	0.46 s, finite output
Fit + predict, 500k rows × 35	0.32 s / 204 MB peak	0.18 s / 29 MB peak

Accuracy: ExtremeLearningMachine vs FastRidge as the imputer regressor

Benchmarked via MultivariateImputer(regressor=...) across 10 seeds, 10% MAR missingness (mean r² on numeric targets):

Dataset	FastRidge	ELM before accuracy fix	ELM after accuracy fix
Diabetes	0.541	0.295	0.490
Titanic (mixed types)	0.447	-0.092	0.532
SyntheticMixed	0.748	0.725	0.737
Nonlinear-relationship synthetic*	0.618	0.918	0.918

*constructed so every column is a nonlinear transform of shared latent factors — the one case ELM should structurally beat a linear model, kept to confirm the fix doesn't just regress ELM toward FastRidge everywhere.

Root cause on Titanic: per-pattern fits use a median of 38 rows against a 100-unit hidden layer (parameters > samples by 2.6x), which is what produced the negative r² before the fix. On typical tabular/time-series data FastRidge remains the stronger default; ExtremeLearningMachine is best reserved for imputation contexts with known nonlinear structure and enough samples per pattern.

Tests

• Two regression tests from the original perf work: varying-width refits stay finite and match a fresh instance; chunked fit matches single-shot at float32 rounding level.
• Five new tests for the accuracy fix: hidden-width capping by sample count, disabling the cap via min_samples_per_feature=0, predict/fit width consistency after capping, and fan-in variance stability across input widths.
• Full suite passes (46 tests).

Docs: docs/api.rst autodocs ExtremeLearningMachine via :members:, so the updated docstring (documenting fan-in scaling and min_samples_per_feature) is already reflected on the API page — verified with a clean local sphinx-build (no warnings).

Note: imputation results with an ELM regressor will differ from previous releases whenever pattern widths varied — the old outputs were garbage (out-of-bounds reads) or badly overfit (fixed hidden width vs. pattern sample size), so no parity is intended.

🤖 Generated with Claude Code