Re: svd_lapack vs eofunc/eofunc_ts

From: Gil Compo <Gilbert.P.Compo_at_nyahnyahspammersnyahnyah>
Date: Tue Nov 29 2011 - 14:30:07 MST

Kerrie,

They should be different by a factor of 1/(n-1), since eofunc is
computing the covariance matrix itself. In your case, that factor is
exactly 1/3.

best wishes,

gil

On 11/29/11 2:14 PM, Kerrie Geil wrote:
> Hi All,
>
> Can anyone tell me why the singular value decomposition (svd_lapack)
> below yields different eigenvalues than the eigenanalysis (eofunc,
> eofunc_ts), while both output the exact same EOFs and PCs? Can't see
> what I'm doing wrong. Code and output are below.
>
> Thanks!
> Kerrie
>
>
> ;######################################################################################
> ;### data array A=(/ntime,nspace/) ###
> A=(/ (/-2.,1./),(/7.,-8./),(/-3.,0./),(/5.,-9./) /)
> A!0="time"
> A!1="space"
> ntime=4 ;### number of rows, time dimension 0 ###
> nspace=2 ;### number of columns, space dimension 1 ###
>
> Aprime = dim_rmvmean_n(A,0)
> ;### singular value decomposition of Aprime ###
> ;### Aprime=(u)(s)(vtranspose)-----> (4x2)=(4x2)(2x2)(2x2) ###
> u = new ( (/ntime,nspace/) , float ) ;###
> columns of u are PC eigenvectors ###
> vtranspose = new ( (/nspace,nspace/) , float )
> ;### rows of vtranspose are EOF eigenvectors ###
> s = svd_lapack (Aprime, "S" , "S", 0, vtranspose, u)
> ;### s^2=eigenvalues
> ;### s is returned with only the non-zero (diagonal) values
>
> print("for SVD")
> print("U:")
> write_matrix (u, "2f10.4", False)
> print("V:")
> write_matrix (vtranspose, "2f10.4", False)
> print("lambda:")
> print(s^2)
>
>
> At=transpose(A) ;### make time the rightmost dimension ###
> eof=eofunc(At,2,False) ;### rows of eof are EOFs ###
> eoftranspose=transpose(eof)
> pc=eofunc_ts_Wrap(At,eof,False) ;### rows of pc are PCs ###
> pctranspose=transpose(pc) ;### for similarity to SVD
> matrices, now columns of pc are PCs ###
> pc_colsum=dim_sum_n(pctranspose^2,0) ;### square each
> component and sum the columns ###
> pctranspose(:,0)=pctranspose(:,0)/sqrt(pc_colsum(0)) ;###
> normalize PC1
> pctranspose(:,1)=pctranspose(:,1)/sqrt(pc_colsum(1)) ;###
> normalize PC2
>
> print("for eofunc")
> print("U:")
> write_matrix (pctranspose, "2f10.4", False)
> print("V:")
> write_matrix (eof, "2f10.4", False)
> print("lambda:")
> print(eof@eval)
> ;######################################################################################
>
>
> OUPUT IS:
>
> (0) for SVD
> (0) U:
>
> 0.5010 0.4050
> -0.5261 0.5745
> 0.4983 -0.3748
> -0.4732 -0.6046
>
> (0) V:
>
> -0.6898 0.7240
> 0.7240 0.6898
>
> (0) lambda:
> (0) 153.4626
> (1) 3.287447
> (0)
> (0)
> (0) for eofunc
> (0) U:
>
> 0.5010 0.4050
> -0.5261 0.5745
> 0.4983 -0.3748
> -0.4732 -0.6046
>
> (0) V:
>
> -0.6898 0.7240
> 0.7240 0.6898
>
> (0) lambda:
> (0) 51.15419
> (1) 1.095816
>
>
>
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk

-- 
+++++++++++++++++++++++++++++++++++++++++++++++++++
Gil Compo, Research Scientist, CIRES
University of Colorado
Mail : CIRES Climate Diagnostics Center
NOAA Physical Sciences Division
Earth System Research Laboratory
325 Broadway R/PSD1, Boulder, CO 80305-3328
Email: compo@colorado.edu
Phone: (303) 497-6115 Fax: (303) 497-6449
http://www.esrl.noaa.gov/psd/people/gilbert.p.compo
+++++++++++++++++++++++++++++++++++++++++++++++++++
"Stop and consider the wondrous works of God."
  Job 37:34
---------------------------------------------------
The contents of this message are mine personally
and do not necessarily reflect any position of
NOAA or the University of Colorado.

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Tue Nov 29 14:30:20 2011

This archive was generated by hypermail 2.1.8 : Wed Nov 30 2011 - 19:52:47 MST