Add Kernel PCA decomposition - #337
Merged
Merged
Conversation
Implements Scholar.Decomposition.KernelPCA, the kernel extension of PCA that runs the decomposition in a reproducing kernel Hilbert space, capturing non-linear structure that ordinary PCA cannot. Supports the linear, poly, rbf, sigmoid and cosine kernels, with fit/2, transform/2 and fit_transform/2 mirroring the existing PCA API. The kernel matrix is double-centered, decomposed with Nx.LinAlg.eigh, and the eigenvectors are sign-flipped for deterministic output. Results match scikit-learn's KernelPCA on all five kernels (covered by the tests). Addresses the Kernel PCA item of elixir-nx#246.
Three issues found while reviewing the initial implementation:
* transform/2 now validates the input has the same number of features used
to fit the model, matching PCA's own validation, instead of raising an
Nx-internal shape error.
* The cosine kernel no longer divides by zero for all-zero rows; it treats
them as zero similarity, matching scikit-learn's normalize().
* Centering the kernel matrix can leave it asymmetric by a few ULPs, which
Nx.LinAlg.eigh rejects outright since it requires exact symmetry (this
could crash fit/2 depending on kernel/gamma/dataset). The centered
kernel is now explicitly symmetrized before the decomposition.
Also expands test coverage: transform/2 and eigenvalues for poly/sigmoid/
cosine (previously only checked via fit_transform/2), custom gamma and
degree, the num_components == num_samples boundary, default kernel/gamma,
eigenvalue ordering, and regression tests for the three fixes above.
Verifying every output against scikit-learn surfaced two more issues:
* With num_components close to num_samples, the smallest eigenvalue of the
centered kernel is zero up to floating-point noise (the centering forces
it), and could come out slightly negative, crashing Nx.sqrt in
fit_transform/2. Indefinite kernels such as sigmoid can also produce
genuinely negative eigenvalues. They are now clipped to zero and their
components project to zero, as in scikit-learn.
* The eps previously passed to Nx.LinAlg.eigh (1.0e-8) stopped the
iteration too early for f64 inputs, leaving eigenpair residuals around
1.0e-4 where LAPACK reaches 1.0e-16, while tighter tolerances made f32
accumulate rounding noise past convergence. The tolerance now matches
the input precision (1.0e-8 for f32, 1.0e-11 for f64), which brings the
f64 projections within 1.0e-5 of the SciPy reference in the
worst-conditioned case measured (500x better than before).
Also validates gamma as a positive number and accepts integers for gamma
and coef0, following the option style used elsewhere in Scholar.
Verifying against scikit-learn on a larger dataset (25 samples) exposed that Nx.LinAlg.eigh can stop before the eigenvalues converge: the leading eigenvalue was exact, but the following ones came out up to a few percent too small (their eigenpair residuals were around 1.0e-2 where the leading one was at 1.0e-6), throwing fit_transform/2 off by the same margin. The 6-sample data used by the other tests was too small to show this. The eigenvectors themselves were accurate, so the eigenvalues are now recomputed as the Rayleigh quotient of the renormalized eigenvectors, whose error is quadratic in the eigenvector error. On the 25-sample data this brings the eigenvalues from an absolute error of 2.7 (poly kernel) to 1.0e-10 against scikit-learn, and fit_transform from 0.3 to 1.0e-5. Adds a regression test with the 25-sample dataset and refreshes the doctest values, which moved by a few ULPs.
The component order comes from sorting the raw eigh eigenvalues, but the final values are the refined ones, so two eigenvalues that eigh resolved within its error could end up out of decreasing order, breaking the documented invariant. They are sorted again after the refinement. Also documents the kernel option fields stored in the struct.
josevalim
approved these changes
Jul 12, 2026
Contributor
|
💚 💙 💜 💛 ❤️ |
2 tasks
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.
Adds
Scholar.Decomposition.KernelPCA, the kernel extension of PCA that runsthe decomposition in a reproducing kernel Hilbert space, capturing non-linear
structure that ordinary PCA cannot. Addresses the Kernel PCA item of #246.
Supports the
:linear,:poly,:rbf,:sigmoidand:cosinekernels,with
fit/2,transform/2andfit_transform/2mirroring the existing PCAAPI. The kernel matrix is double-centered, decomposed with
Nx.LinAlg.eigh,and the eigenvectors are sign-flipped for deterministic output.
A few numerical decisions worth calling out:
since centering can leave asymmetries of a few ULPs and
eighrequiresexact symmetry.
eighcan stop before the eigenvalues fully converge while theeigenvectors are already accurate, so the eigenvalues are recomputed as
the Rayleigh quotient of the renormalized eigenvectors. On a 25-sample
dataset this brings them from a few percent off to within 1.0e-10 of
scikit-learn (covered by a regression test).
eightolerance follows the input precision (f32/f64): a tolerancebelow the floating-point resolution makes the iteration accumulate
rounding noise past convergence.
kernels such as sigmoid) are clipped and their components project to
zero, as in scikit-learn.
Verification
Results match
sklearn.decomposition.KernelPCAon all five kernels foreigenvalues,fit_transformandtransform, including customgamma/degree, zero-norm rows with the cosine kernel, and
num_components == num_samples. f32 agrees to roughly 1.0e-4 and f64 toroughly 1.0e-5 or better. The reference values in the tests were generated
with scikit-learn 1.6.1.