[Insight-developers] faster PCA for the shape levelset

[Carl-Fredrik Westin]

This is a comment to the attendees of yesterday's Pittsburgh meeting

As I mentioned yesterday there is a faster way to calculate the
eigen modes (if I understood correctly the way it is done now).
Here are the details:

Doing PCA of the sum of the outer products of the distance transform
example images is extremely expensive and is not possible for larger
images, better to do eigen decomposition on the matrix based on the
inner products of the example images, and then transform these vectors
to the image domain.

Example:
   3 images of size NxN

   Stack images as column vectors in a matrix A = (im1 im2 im3), giving
   an (NxN) by 3 matrix.

   - Doing PCA on AA^T, an (NxN) by (NxN)  matrix, gives:  U*S*V = AA^T,
   where the three first column vectors in U are the eigen modes. This is
   only doable if N is reasonably small.

   - Better to do eigenvector decomposition (or PCA) of A^TA, in this
   example a 3x3 matrix, E Lam E = A^TA, where E is a 3x3 matrix with
   the stacked eigenvectors, and then transform these vectors using A:
   U2 = A*E. This should give the same vectors as in U (need to play
   with the sign, and normalize the vectors to compare them).

I don't know how Leventon did in his implementation, but I would assume
he used this trick.

I played with coding this example on the flight home yesterday, and it
seems to work as expected. Let me know if you want the example code (not
itk though).

C-F

-- 
Dr. Carl-Fredrik Westin
Director, Laboratory of Mathematics in Imaging (LMI)
Assistant Professor of Radiology, Harvard Medical School
Brigham and Women's Hospital, 75 Francis St., Boston, MA 02115

Phone: (+1) 617-278-0639, Fax: (+1) 617-264-6887
Email:  westin at bwh . harvard . edu, Office: Thorn Building 323
http://www . spl . harvard . edu/~westin