Singular Value Decomposition in 3 Dimensions [message #62289] |
Tue, 02 September 2008 09:33  |
tomandwilltamu08
Messages: 4 Registered: September 2008
|
Junior Member |
|
|
I am wondering how to do Singular Value Decomposition in 3 Dimensions
in IDL. All of the canned routines seem to work only on 2D arrays.
Specifically, I am trying to preform Principle Component Analysis on
stacks of 2D images.
For example, how can one preform an SVD on a 2048x2048xn array to get
2048x2048 principle components?
Thanks much,
-Will
|
|
|
Re: Singular Value Decomposition in 3 Dimensions [message #62367 is a reply to message #62289] |
Wed, 03 September 2008 04:56  |
Juggernaut
Messages: 83 Registered: June 2008
|
Member |
|
|
On Sep 3, 7:52 am, Bennett <juggernau...@gmail.com> wrote:
> On Sep 2, 12:33 pm, tomandwilltam...@gmail.com wrote:
>
>> I am wondering how to do Singular Value Decomposition in 3 Dimensions
>> in IDL. All of the canned routines seem to work only on 2D arrays.
>
>> Specifically, I am trying to preform Principle Component Analysis on
>> stacks of 2D images.
>
>> For example, how can one preform an SVD on a 2048x2048xn array to get
>> 2048x2048 principle components?
>
>> Thanks much,
>> -Will
>
> If you want the principal components for the
> 3D array you can do something like this
> sz = size(array, /dimensions)
> newArray = fltarr(sz[2], sz[1]*sz[0])
> FOR i=0, sz[2]-1 DO BEGIN
> newArray[i,*] = transpose(reform(array[*,*,i], sz[0]*sz[1]))
> ENDFOR
> result = pcomp(newArray, eigenvalues=evals, /standardize)
>
> pcomp() is IDLs built in for doing PCA and result will be
> an array of I believe the same dimensions of newArray which to
> get back into viewing form you could just reform it back like
>
> tv, reform(result[0,*],sz[0],sz[1])
>
> There may be better ways of doing it but I may as well give
> you a point to jump off of
By the way the for loop can be eliminated by just putting
transpose(reform(array, sz[0]*sz[1], sz[2])) into pcomp
> FOR i=0, sz[2]-1 DO BEGIN
> newArray[i,*] = transpose(reform(array[*,*,i], sz[0]*sz[1]))
> ENDFOR
> result = pcomp(newArray, eigenvalues=evals, /standardize)
becomes
result = pcomp(transpose(reform(array, sz[0]*sz[1], sz[2])), ...)
|
|
|
Re: Singular Value Decomposition in 3 Dimensions [message #62368 is a reply to message #62289] |
Wed, 03 September 2008 04:46  |
Mort Canty
Messages: 134 Registered: March 2003
|
Senior Member |
|
|
tomandwilltamu08@gmail.com schrieb:
> I am wondering how to do Singular Value Decomposition in 3 Dimensions
> in IDL. All of the canned routines seem to work only on 2D arrays.
>
> Specifically, I am trying to preform Principle Component Analysis on
> stacks of 2D images.
>
> For example, how can one preform an SVD on a 2048x2048xn array to get
> 2048x2048 principle components?
>
> Thanks much,
> -Will
SVD is a decomposition theorem for matrices (2D arrays). I think you may
have an incorrect understanding of what principal (not principle)
component analysis means. If you apply PCA to a stack of n 2048x2048
images, you will get n 2048x2048 principal component images. The n
images will be uncorrelated and stacked in the order of decreasing
variance.
Mort
|
|
|
Re: Singular Value Decomposition in 3 Dimensions [message #62369 is a reply to message #62289] |
Wed, 03 September 2008 04:52  |
Juggernaut
Messages: 83 Registered: June 2008
|
Member |
|
|
On Sep 2, 12:33 pm, tomandwilltam...@gmail.com wrote:
> I am wondering how to do Singular Value Decomposition in 3 Dimensions
> in IDL. All of the canned routines seem to work only on 2D arrays.
>
> Specifically, I am trying to preform Principle Component Analysis on
> stacks of 2D images.
>
> For example, how can one preform an SVD on a 2048x2048xn array to get
> 2048x2048 principle components?
>
> Thanks much,
> -Will
If you want the principal components for the
3D array you can do something like this
sz = size(array, /dimensions)
newArray = fltarr(sz[2], sz[1]*sz[0])
FOR i=0, sz[2]-1 DO BEGIN
newArray[i,*] = transpose(reform(array[*,*,i], sz[0]*sz[1]))
ENDFOR
result = pcomp(newArray, eigenvalues=evals, /standardize)
pcomp() is IDLs built in for doing PCA and result will be
an array of I believe the same dimensions of newArray which to
get back into viewing form you could just reform it back like
tv, reform(result[0,*],sz[0],sz[1])
There may be better ways of doing it but I may as well give
you a point to jump off of
|
|
|