Skip to content

Fix pca_numpy conditioning: SVD the centered data matrix, not the covariance matrix#1524

Open
mengxihex wants to merge 2 commits into
compas-dev:mainfrom
mengxihex:pca-numpy-conditioning
Open

Fix pca_numpy conditioning: SVD the centered data matrix, not the covariance matrix#1524
mengxihex wants to merge 2 commits into
compas-dev:mainfrom
mengxihex:pca-numpy-conditioning

Conversation

@mengxihex

Copy link
Copy Markdown

Fixes #1522.

What

pca_numpy formed the covariance matrix C = Y.T @ Y / (n - 1) and SVD'd it. Forming C squares the condition number, so for near-degenerate inputs (an almost collinear point cloud whose smallest extent is ~1e-8 of its largest — still an exactly planar, perfectly posed fit) the smallest principal direction drowns at machine epsilon and comes back wrong: bestfit_plane_numpy normals were off by 11.6° (1e-8 aspect) to 37.3° (1e-9 aspect) on the reproduction in #1522.

This PR runs the SVD on the centered data matrix Y directly. The right-singular vectors of Y are exactly the eigenvectors of C, and the singular values are rescaled (s**2 / (n - 1)) so the returned eigenvalues keep their meaning (variance along each principal direction). Well-conditioned inputs return the same results as before; only the ill-conditioned ones improve. Cost is equivalent (both routes are O(n·dim²)).

Checks

  • New tests/compas/geometry/test_pca_numpy.py: (1) well-conditioned parity with the covariance route (eigenvalues match eigvalsh(C), eigenvectors diagonalize C), (2) the near-collinear cloud recovers its smallest direction to <1e-4°, (3) bestfit_plane_numpy regression on the same cloud (fails on the previous implementation with ~11.6° error).
  • Full tests/compas/geometry suite: 519 passed (includes the pca_numpy consumers bestfit_numpy, bbox_numpy, icp_numpy).
  • CHANGELOG entry added under Unreleased → Changed.

One behavioural note for review: SVD sign indeterminacy means individual principal directions may flip sign relative to the old implementation on some inputs (the sign was already unspecified in both routes).

mengxihex added 2 commits July 6, 2026 19:21
… matrix

Forming C = Y.T @ Y / (n - 1) squares the condition number, so near-
degenerate inputs (almost collinear point clouds) lose their smallest
principal direction to floating-point rounding - bestfit_plane_numpy
normals came back 11-37 degrees off on exactly planar sliver clouds.
The right-singular vectors of Y are the same eigenvectors; eigenvalues
are rescaled (s**2 / (n - 1)) to keep their variance meaning, so well-
conditioned results are unchanged.

Fixes compas-dev#1522
@mengxihex

Copy link
Copy Markdown
Author

Fixed the windows-ironpython failure: the new test file used the @ matmul operator, which is a parse-time SyntaxError under IronPython 2.7. It now follows the test_bbox.py idiom — py2-parseable, if compas.IPY: return guards, numpy imported inside the guarded test bodies. Verified locally: 3/3 new tests + the full tests/compas/geometry suite pass on CPython. Workflows on the new commit are awaiting approval.

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment

Labels

None yet

Projects

None yet

Development

Successfully merging this pull request may close these issues.

bestfit_plane_numpy returns wrong normals on near-collinear point clouds (covariance step squares the conditioning)

1 participant