RustyML - v0.13.0

RustyML v0.13.0 is a large performance-and-correctness release. The matrix-product backend moves to the pure-Rust gemm crate, model persistence switches to a compact binary format, several estimators gain new capabilities, and a broad numerics pass aligns the crate with scikit-learn / NumPy / PyTorch conventions.

Heads-up: this release contains breaking API and on-disk-format changes.
Read Breaking Changes before upgrading.

Highlights

• New matmul backend — all GEMM/GEMV now runs on the gemm crate's runtime-dispatched SIMD kernels, with shape-aware parallelism.
• Compact binary model persistence — models and weights serialize via postcard instead of JSON (about 5x smaller on a fitted KMeans).
• Richer estimators — data-dependent kernel gamma, LDA probabilities, a closed-form linear-regression solver, Isolation Forest labels, and a squared-hinge LinearSVC.
• Runtime-tunable parallelism — the new public tuning module overrides the serial/parallel thresholds per machine, with no recompile.
• scikit-learn-aligned numerics — convergence criteria, regularization scaling, and stability constants now match the reference implementations.

Breaking Changes

Model persistence: JSON to binary (postcard)

• save_to_path / load_from_path on every classical-ML model, PCA / KernelPCA, and Sequential now write and read postcard binary instead of serde_json.
• Old .json model files can no longer be loaded — re-save any persisted models. The method signatures are unchanged; only the byte format differs (a .bin extension is now suggested).
• IoError::Json(serde_json::Error) is renamed to IoError::Serialization(postcard::Error); code matching on IoError::Json must switch to IoError::Serialization.
• The serialized format also changed for LinearRegression (new solver field) and LinearSVC (new loss and learning_rate_decay fields), so previously saved instances of those two models must also be re-saved.

Matrix-product backend: matrixmultiply to gemm

• The crate-internal GEMM/GEMV path now calls the gemm crate; matrixmultiply is dropped.
• The public matrix-product API in math::matmul — gemm, gemv, gemm_par, and gemv_par — is removed, with no public replacement. Out-of-crate callers should use ndarray's dot or the gemm crate directly.
• Public model signatures (Dense, PCA, KMeans, SVM, LDA, kernels) are unchanged; only the underlying kernel and its numerics differ (see Behavior Changes).

Kernel gamma is now a Gamma enum

• The gamma field of KernelType::Poly / RBF / Sigmoid changes from f64 to the new Gamma enum.
• Update every construction site, e.g. KernelType::RBF { gamma: 0.5 } becomes KernelType::RBF { gamma: Gamma::Value(0.5) }.

standardize drops its epsilon parameter

• The signature changes from standardize(data, axis, epsilon) to standardize(&data, axis).
• Drop the third argument, e.g. standardize(&data, StandardizationAxis::Column, 1e-8) becomes standardize(&data, StandardizationAxis::Column).

LDA::get_n_components returns Option<usize>

• The default n_components is now None (auto-resolved at fit time to min(n_classes - 1, n_features)) instead of 2.
• Callers reading get_n_components() must handle Option<usize>.

New Features

Estimators

• Data-dependent kernel gamma — the new Gamma enum adds Gamma::Scale (1 / (n_features * X.var())) and Gamma::Auto (1 / n_features) alongside Gamma::Value(f64); SVC and KernelPCA resolve Scale / Auto against the training data at fit.
• LDA::decision_function and LDA::predict_proba — per-class discriminant scores and a numerically-stable softmax posterior, plus automatic n_components selection.
• LinearRegression closed-form solver and score — the new Solver::Normal (selectable via with_solver) solves the ridge least-squares system by SVD, and score(x, y) returns the coefficient of determination (R²).
• IsolationForest::predict_labels(x, contamination) — classifies samples as inlier (+1) or outlier (-1) by flagging the top ceil(contamination * n_samples) scores.
• LinearSVC squared-hinge loss and learning-rate decay — the new Loss enum adds Loss::SquaredHinge, and with_learning_rate_decay applies inverse-scaling decay per epoch.
• t-SNE min_grad_norm early stopping — TSNE::with_min_grad_norm stops optimization once the largest gradient entry drops below the threshold (default 1e-7); pass 0.0 to keep the previous full-n_iter behavior.

tuning module

• The new public tuning module exposes a flat set_* / get_* facade (grouped into matmul, elementwise, reduction, tree, conv, pool, norm, and metrics) for overriding the serial-vs-parallel and GEMM-strategy thresholds at runtime, per machine, without recompiling.
• Defaults are unchanged, so results stay reproducible across runs on the same machine.

Performance

• gemm-crate GEMM/GEMV backend with shape-aware parallelism: serial below a per-dtype FLOP gate, column-parallel for wide outputs, and a rayon row-split for thin / column-starved outputs and matvecs.
• Bounded nested GEMM parallelism — three hotspots that re-forked rayon inside an already-parallel loop now serialize the inner GEMM once the outer axis fills the pool: conv backward ~13% faster (batch=32), and many-class LDA::fit and high-k KMeans fit ~6% faster each.
• Matvecs routed through the gemm backend — svc_predict -5%, knn_predict -22%, mean_shift_fit -47%.
• silhouette_score symmetric fill — each unordered pair is evaluated once instead of twice: ~33% faster for Euclidean / Manhattan, ~45% for Minkowski(3).
• generate_polynomial_features gating — small-input expansion -94.5%, while the large first-order path stays parallel.
• Power-iteration eigensolver — one matvec per step (was two) and an in-place rank-1 Hotelling deflation, feeding the PCA / KernelPCA PowerIteration solvers.
• Every performance change above is numerically unchanged and stays reproducible across runs on the same machine.

Behavior Changes

• Matrix products stay reproducible across runs on the same machine, but are no longer bit-for-bit identical to the old matrixmultiply backend (or identical across different machines); tests or pipelines asserting bit-exact matmul output must switch to an approximate comparison.
• KMeans declares convergence on centroid shift (<= tol * mean(var(X))) rather than inertia change, so the iteration count and final centroids can change for a given tol.
• LogisticRegression no longer divides the L1/L2 penalty by the sample count (matching scikit-learn's SGD convention), so fitted coefficients change for any regularized model at a fixed alpha.
• DecisionTree scales min_impurity_decrease by the node's sample fraction and enforces min_samples_leaf during the split search, which can change the resulting tree.
• PCA flips principal-axis signs deterministically, so all SVD solvers and repeated runs agree on orientation; reconstructions are unaffected, but component signs may flip relative to the previous release.
• NaN and stability alignment — math::sigmoid drops its ±500 input clamp, the softmax layer drops its 1e-8 denominator floor (NaN now propagates), utils::normalize leaves near-zero lanes unchanged, and log_loss, mean_absolute_percentage_error, and t-SNE align their stability epsilons to the reference values. Finite, well-formed inputs are numerically unchanged.

Bug Fixes

• LinearSVC::fit now rejects labels other than 0.0 / 1.0 instead of silently mishandling multi-class input.
• MeanShift keeps the current center on a zero-weight window instead of teleporting it to the origin (which injected a spurious cluster).
• DecisionTree no longer discards a whole categorical split because of a single under-sized rare branch.
• IsolationForest::predict handles non-contiguous (sliced or transposed) input rows instead of panicking.
• The lightweight metrics-only build compiles again.

Dependencies

• Added gemm and postcard; dropped matrixmultiply and serde_json.
• All other dependencies are at their latest stable versions.

Full Changelog

See CHANGELOG.md for the complete, commit-level list of changes.