Fix riemann ASR eigpair mismatch in nonlinear_eigenspace

Author: nbara

meegkit/utils/covariances.py

@@ -503,5 +503,10 @@ def ehess(X, U):
     problem = Problem(manifold=manifold, cost=cost, euclidean_gradient=egrad,
                       euclidean_hessian=ehess)
     Xsol = solver.run(problem, initial_point=U0)
+    X = Xsol.point
 
-    return S0, Xsol.point
+    # Keep eigpairs consistent: compute each optimized vector's Rayleigh
+    # quotient against L instead of returning the random-init eigenvalues S0.
+    S = np.real(np.diag(X.T @ L @ X))
+
+    return S, X

tests/test_cov.py

@@ -1,4 +1,5 @@
 import numpy as np
+import pytest
 from numpy.testing import assert_almost_equal
 
 from meegkit.utils import convmtx, tscov, tsxcov
@@ -88,6 +89,25 @@ def test_convmtx():
                   ])
     )
 
+
+def test_nonlinear_eigenspace_consistent_eigpairs():
+    """Returned eigenvalues must match returned eigenvectors."""
+    pytest.importorskip("pymanopt")
+    from meegkit.utils.covariances import nonlinear_eigenspace
+
+    A = rng.standard_normal((6, 6))
+    L = A.T @ A + np.eye(6) * 1e-6
+
+    S, X = nonlinear_eigenspace(L, 6)
+    rayleigh = np.real(np.diag(X.T @ L @ X))
+
+    np.testing.assert_allclose(
+        np.sort(np.real(S)),
+        np.sort(rayleigh),
+        rtol=1e-3,
+        atol=1e-6,
+    )
+
 if __name__ == "__main__":
     import pytest
     pytest.main([__file__])