PCA? about nd4js HOT 3 CLOSED

Just for heads-up: I've compared my nd4js implementation with the dedicated ml-pca library, which also uses SVD under the hood (as it should). There were a few things I had not properly done in my notebook (for example, the "variance explained" by each PC should have been the normalized squared eigenvalues), but now I have it fully aligned. I'll do a proper write-up … later, because… time.

from nd4js.

Hi @Fil

Short Answer

In case of the SVD, both U and V are orthonormal matrices, which means the inverse of U is simply its transpose U.T. My knowledge on the PCA is a little blurry, but here's a little example (which is hopefully correct). The data is from Josh Starmer's YouTube PCA Tutorial.

Edit: I've made a mistake in example: You have to project by V.T instead of U.T in the second-to-last line.

const X_raw = nd.array(
  //Gene1, Gene2, Gene3
  [[   10,   6  ,  12   ], // Mouse 1
   [   11,   4  ,   9   ], // Mouse 2
   [    8,   5  ,  10   ], // Mouse 3
   [    3,   3  ,   2.5 ], // Mouse 4
   [    2,   2.8,   1.3 ], // Mouse 5
   [    1,   1  ,   2   ]] // Mouse 6
);

const [M,N] = X_raw.shape;

const X_mean = X_raw.reduceElems(/*along_axis=*/0, (x,y) => x+y).mapElems(sum => sum/M);

const X_centered = nd.zip_elems([X_raw, X_mean], (raw,mean) => raw-mean );

const [U,sv,V] = nd.la.svd_decomp(X_centered);

const X_projected = nd.la.matmul(X_centered, V.T);
console.log('X_projected:\n' + X_projected);

Longer Answer

There is no matrix_invert in ND4JS and there probably never will be. Computing the inverse of a matrix directly is almost always a bad idea, as explained in this blog post. Instead of inv(A) @ y, You can use nd.la.solve(A,y) in ND4JS. Alternatively, You can use one of the many matrix factorization directly. To use the QR decomposition for example, You could use: nd.la.qr_solve(nd.la.qr_decomp(A),y).

Hope this helps
Dirk

from nd4js.

Ah yes, transposition! I knew there was a trick—my maths are quite a bit rusty… thanks a lot for the detailed explanation too. (I had noticed that V was the unit matrix, not U.)

It works quite well: https://observablehq.com/d/0dd301be2ad20e63 (not published yet since it's in dear need of explanations… but the 3D chart below (made with the Wine UCI dataset) appears to be correct.

Interestingly, the Penguins dataset doesn't "work" with 2 PC, you really need 3 to cleanly separate the species.

PS: I still like the "explorability" of matrix.inverse, but now I understand why it's preferable to avoid it in an applied library.

from nd4js.