10
$\begingroup$

I understand the relation between Principal Component Analysis and Singular Value Decomposition at an algebraic/exact level. My question is about the scikit-learn implementation.

The documentation says: "[TruncatedSVD] is very similar to PCA, but operates on sample vectors directly, instead of on a covariance matrix.", which would reflect the algebraic difference between both approaches. However, it later says: "This estimator [TruncatedSVD] supports two algorithm: a fast randomized SVD solver, and a “naive” algorithm that uses ARPACK as an eigensolver on (X * X.T) or (X.T * X), whichever is more efficient.". Regarding PCA, it says: "Linear dimensionality reduction using Singular Value Decomposition of the data to project it ...". And PCA implementation supports the same two algorithms (randomized and ARPACK) solvers plus another one, LAPACK. Looking into the code I can see that both ARPACK and LAPACK in both PCA and TruncatedSVD do svd on sample data X, ARPACK being able to deal with sparse matrices (using svds).

So, aside from different attributes and methods and that PCA can additionally do exact full singular value decomposition using LAPACK, PCA and TruncatedSVD scikit-learn implementations seem to be exactly the same algorithm. First question: Is this correct?

Second question: even though LAPACK and ARPACK use scipy.linalg.svd(X) and scipy.linalg.svds(X), being X the sample matrix, they compute the singular value decomposition or eigen-decomposition of $X^T*X$ or $X*X^T$ internally. While the "randomized" solver doesn't need to compute the product. (This is relevant in connection with numerical stability, see Why PCA of data by means of SVD of the data?). Is this correct?

Relevant code: PCA line 415. TruncatedSVD line 137.

$\endgroup$
  • 1
    $\begingroup$ could you add a link to the code $\endgroup$ – seanv507 Oct 10 '16 at 21:49
  • 1
    $\begingroup$ drake - I think I agree with you on first Q. don't understand the second. what do you mean 'they compute the singular value decomposition or eigen-decomposition of XT∗XXT∗X or X∗XTX∗XT internally' - you have just shown the code where it is all done using SVD on X? - numerical issues refer to first computing covariance matrix (call it C) then finding eigenvectors of C $\endgroup$ – seanv507 Oct 11 '16 at 0:24
  • $\begingroup$ @seanv507 Regarding 2nd question -- I guess that scipy.linalg.svd(X) computes the svd by doing the eigen-decomposition of $X^T * X$ or/and $X * X^T$. Same thing for linalg.svds(X). Quoting: "a fast randomized SVD solver, and a “naive” algorithm that uses ARPACK as an eigensolver on (X * X.T) or (X.T * X)". See also last line in docs.scipy.org/doc/scipy/reference/generated/… . The only way I can understand the first quote is that the randomized algorithm is the only one that doesn't compute the covariance/gram matrix $\endgroup$ – drake Oct 11 '16 at 0:49
  • 1
    $\begingroup$ I would guess the ARPACK approach has to do with something like Arnoldi iteration, so it only has to do matrix-vector products. (In principle these sort of iterative methods do not even an explicit $X$, just a pair of routines Xtimes() and Xt_times(). This is common for large sparse matrices in PDE solvers, for example.) $\endgroup$ – GeoMatt22 Oct 11 '16 at 1:01
  • $\begingroup$ @GeoMatt22 Can you elaborate on your comment? Do you mean that ARPACK or LAPACK approaches don't suffer from numerical instabilities because they don't need to compute the covariance matrix? $\endgroup$ – drake Oct 11 '16 at 1:10
9
$\begingroup$

PCA and TruncatedSVD scikit-learn implementations seem to be exactly the same algorithm.

No: PCA is (truncated) SVD on centered data (by per-feature mean substraction). If the data is already centered, those two classes will do the same.

In practice TruncatedSVD is useful on large sparse datasets which cannot be centered without making the memory usage explode.

  • numpy.linalg.svd and scipy.linalg.svd both rely on LAPACK _GESDD described here: http://www.netlib.org/lapack/lug/node32.html (divide and conquer driver)

  • scipy.sparse.linalg.svds relies on ARPACK to do a eigen value decomposition of XT . X or X . X.T (depending on the shape of the data) via the Arnoldi iteration method. The HTML user guide of ARPACK has a broken formatting which hides the computational details but the Arnoldi iteration is well described on wikipedia: https://en.wikipedia.org/wiki/Arnoldi_iteration

Here is the code for the ARPACK-based SVD in scipy:

https://github.com/scipy/scipy/blob/master/scipy/sparse/linalg/eigen/arpack/arpack.py#L1642 (search for the string for "def svds" in case of line change in the source code).

$\endgroup$
  • 2
    $\begingroup$ One can support sparse data efficiently (TruncatedSVD), the other cannot (PCA). This is why we have 2 classes. $\endgroup$ – ogrisel Oct 12 '16 at 14:43
  • 1
    $\begingroup$ If that's the reason, then I'd call them SVD and SparseSVD (or similar) to avoid confusion. $\endgroup$ – drake Oct 12 '16 at 22:37
  • 2
    $\begingroup$ But people want PCA and they might not know that PCA is just SVD on centered data. $\endgroup$ – ogrisel Oct 13 '16 at 7:51
  • 3
    $\begingroup$ @drake I disagree that "the procedures are different (PCA uses covariance matrix and SVD uses data matrix)". PCA is a name for the type of analysis. One can use different algorithms and implementations to carry it out. EIG of cov matrix is one method, SVD of centered data matrix is another method, and then EIG and SVD can be performed by various methods as well. Does not matter - all of that is PCA. $\endgroup$ – amoeba Oct 18 '16 at 12:10
  • 1
    $\begingroup$ @amoeba Thanks for the clarification/correction on the terminology. What you say makes more sense to me, given that SVD and EIG are algebraic theorems/methods with a wider scope than PCA $\endgroup$ – drake Oct 18 '16 at 15:56

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.