Re: About Rotated EOF patterns of NCEP reanalysis 500mb monthly data

From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>
Date: Tue Dec 04 2012 - 21:26:45 MST

FYI: To test for orthogonality
      (1) Space: test via dot product
      (2) Time: test vis correlation
; =================================================================
; Example for say neof=3
;==================================================================
; correlation
;==================================================================
   r01 = escorc(eof_ts(0,:), eof_ts(1,:))
   r12 = escorc(eof_ts(1,:), eof_ts(2,:))
   r02 = escorc(eof_ts(0,:), eof_ts(2,:))
   print("r01="+r01+" r12="+r12+" r02="+r02)
;==================================================================
; dot product
;==================================================================
   d01 = sum(eof(0,:,:)*eof(1,:,:))
   d12 = sum(eof(1,:,:)*eof(2,:,:))
   d02 = sum(eof(0,:,:)*eof(2,:,:))
   print("d01="+d01+" d12="+d12+" d02="+d02)

Note: If the original data are type float, it is possible (likely)
       the resulting dot products and correlations will not be
       exactly zero. The reason is that NCL promotes the floats
       to double; does the computations in double precision;
       converts double to float for the return. Roundoff occurs.

You may see numbers like:

(0) r01= 7.46909e-10 r12= 1.70884e-08 r02=6.63595e-09
(0) d01=-1.79949e-07 d12=-2.32831e-08 d02=3.76895e-08

which are 'numerically 0.0'

On 12/4/12 5:39 PM, Adam Phillips wrote:
> Hi Waqar,
> I do not have any experience with calculating rotated EOFs, but I know
> that the reason why one would regress anomalies onto the standardized
> PC-timeseries is so that you get a resulting pattern in terms of "hPa
> (ex. if your EOF pattern is based on SLP) per standard deviation change
> in the timeseries". It does not have anything to do with preserving
> orthogonality as far as I know. If you were to alternatively plot the
> EOF pattern then the units would simply be "hPa (again, for example)",
> but the values would quite possibly be much smaller than typical SLP
> values you are used to dealing with.
>
> Hopefully that answers your question. If not, please let ncl-talk know.
> Best regards,
> Adam
>
> On 12/4/12 4:16 PM, Waqar Younas wrote:
>> Hi i need to post my question again as in previous email, there font
>> size is too small to see. Sorry for that.
>> My question is related to orthogonality of the rotated EOFs
>> In NCL code eof_4.ncl, the following is to calculate normal EOFs
>> eof = eofunc_Wrap(wx, neof, optEOF)
>> eof = -1*eof ; *special* match
>> sign of CPC
>>
>> eof_ts = eofunc_ts_Wrap (wx, eof, optETS)
>>
>> printVarSummary( eof ) ; examine EOF variables
>> printVarSummary( eof_ts )
>>
>> eof_ts = dim_standardize_n( eof_ts, 0, 1) ; normalize
>>
>> ; =================================================================
>> ; Regress
>> ; =================================================================
>>
>> eof_regres = eof ; create an array w
>> meta data
>> do ne=0,neof-1
>> eof_regres(ne,:,:) = (/ regCoef(eof_ts(ne,:),
>> xAnom(lat|:,lon|:,time|:)) /)
>> end do
>>
>> Here, First EOFs were calculated and the the time series from these
>> normal EOFs is regressed on anomalies to get "new" EOFs. I am
>> wondering why this regression step is necessary? Is it to ensure
>>
>> to ensure the orthogonality? For the rotated EOFs, i just rotate the
>> normal EOFs and using these rotated EOFs, i got the rotated time
>> series as following.
>>
>> ev_rot = eofunc_varimax_Wrap(eof,1)
>>
>> eofunc_varimax_reorder(ev_rot)
>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>> ;;;;;;;;;;;;;;;;;;;;
>>
>> ev_rot1 = eofunc_ts_Wrap (wx, ev_rot, optETS)
>> ev_rot1 = dim_standardize_n( ev_rot1, 0, 1)
>> ; =================================================================
>> ; Regress
>> ; =================================================================
>>
>> ev_regres = ev_rot ; create an array w
>> meta data
>> do ne=0,neof-1
>> ev_regres(ne,:,:) = (/ regCoef(ev_rot1(ne,:),
>> xAnom(lat|:,lon|:,time|:)) /)
>> end do
>>
>> Now this ev_regres is my regressed rotated EOFs by regressing the
>> rotated time series on anomalies. Again is this step is to ensure the
>> orthogonality? I am little confused about this step of regressing the
>> time series on anomaly data in normal EOF and rotated EOFs. Normally
>> we project the EOF pattern on anomalies data, instead of projecting
>> time series on anomalies.
>> I would be very thankful for any explanation because, i used this
>> regressed pattern to get my time series again.
>>
>>
>>
>>
>>
>>
>>
>>
>> On Tue, Dec 4, 2012 at 3:11 PM, Waqar Younas <vickyqau@gmail.com
>> <mailto:vickyqau@gmail.com>> wrote:
>>
>> Hi
>> I have one more question regarding orthogonality of rotated EOFs.
>> In the ncl code eof_4.ncl,
>>
>> eof = eofunc_Wrap(wx, neof, optEOF)
>> eof = -1*eof ; *special* match sign of CPC
>>
>> eof_ts = eofunc_ts_Wrap (wx, eof, optETS)
>>
>> printVarSummary( eof ) ; examine EOF variables
>> printVarSummary( eof_ts )
>>
>> eof_ts = dim_standardize_n( eof_ts, 0, 1) ; normalize
>>
>> ; =================================================================
>> ; Regress
>> ; =================================================================
>>
>> eof_regres = eof ; create an array w meta data
>> do ne=0,neof-1
>> eof_regres(ne,:,:) = (/ regCoef(eof_ts(ne,:), xAnom(lat|:,lon|:,time|:)) /)
>> end do
>>
>> First EOFs were calculated and then the time series from these EOFs is regressed on anomalies to get "new" EOFs. I am wondering why this step is necessary in simple EOFs. Is it necessary to
>>
>>
>> to ensure the orthogonality? For the rotated EOFs, i just rotate
>> the normal EOFs and using these rotated EOFs, i got the rotated
>> time series as following.
>>
>> ev_rot = eofunc_varimax_Wrap(eof,1)
>>
>> eofunc_varimax_reorder(ev_rot)
>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
>>
>> ev_rot1 = eofunc_ts_Wrap (wx, ev_rot, optETS)
>> ev_rot1 = dim_standardize_n( ev_rot1, 0, 1)
>> ; =================================================================
>> ; Regress
>> ; =================================================================
>>
>> ev_regres = ev_rot ; create an
>> array w meta data
>> do ne=0,neof-1
>> ev_regres(ne,:,:) = (/ regCoef(ev_rot1(ne,:),
>> xAnom(lat|:,lon|:,time|:)) /)
>> end do
>>
>> Now this ev_regres is my regressed rotated EOFs by regressing the
>> rotated time series on anomalies. Again is this step is to ensure
>> the orthogonality? I am little confused about this step of
>> regressing the time series on anomaly data in normal EOF and
>> rotated EOFs. Normally we project the EOF pattern on anomalies
>> data, instead of projecting time series on anomalies.
>> I would be very thankful for any explanation because, i used this
>> regressed pattern to get my time series again.
>>
>>
>> On Fri, Nov 30, 2012 at 12:39 PM, Cara-Lyn Lappen
>> <lappen7@gmail.com <mailto:lappen7@gmail.com>> wrote:
>>
>> Yes, this is a correct assessment of rotated vs regular EOFs
>>
>>
>>
>>
>>> Hi
>>> Recently i came to know that the procedure to get the time
>>> series of rotated EOF is not simple as getting time series of
>>> simple EOFs. The reason is that rotated EOFs are not
>>> orthogonal to each other so to get its time series we need
>>> regress the rotated EOF patterns on anomaly data to get the
>>> time series. In simple EOFs, we can simple project (multiply)
>>> the EOF on anomaly data to get time series.
>>> I just want to confirm that my understanding is right?
>>> Thanks for your help in advance.
>>>
>>>
>>>
>>> On Thu, Nov 8, 2012 at 9:09 PM, Dennis Shea <shea@ucar.edu
>>> <mailto:shea@ucar.edu>> wrote:
>>>
>>> Did you read the documentation for
>>>
>>> http://www.ncl.ucar.edu/Document/Functions/Built-in/eofunc_varimax.shtml
>>> and
>>> http://www.ncl.ucar.edu/Document/Functions/Contributed/eofunc_varimax_reorder.shtml
>>>
>>> [snip]
>>> eof = eofunc_Wrap(x, neof, optEOF)
>>> eof_ts = eofunc_ts_Wrap (x, eof, optETS)
>>>
>>> printVarSummary( eof ) ; examine EOF
>>> variables
>>> printVarSummary( eof_ts )
>>> print("eof_ts: min="+min(eof_ts)+" max="+max(eof_ts) )
>>>
>>> ;
>>> =================================================================
>>> ; Perform varimax rotation
>>> ;
>>> =================================================================
>>> eof_rot = eofunc_varimax_Wrap( eof, 1 )
>>> printVarSummary( eof_rot )
>>> print("eof_rot: min="+min(eof_rot)+" max="+max(eof_rot) )
>>>
>>> ;
>>> =================================================================
>>> ; put rotated EOFs into descending order (% variance
>>> explained)
>>> ;
>>> =================================================================
>>> eofunc_varimax_reorder( eof_rot )
>>> printVarSummary( eof_rot )
>>> [snip]
>>>
>>> The percent variance explained after rotation is returned as
>>> the attribute 'pcvar_varimax'
>>>
>>>
>>>
>>>
>>> On 11/8/12 8:10 PM, Waqar Younas wrote:
>>>
>>> Hi
>>>
>>> I have one more question about rotated EOFs.
>>> When i use the following loop to plot my rotated EOF
>>> with variance
>>> explained, it is not working for variance explained
>>> of rotated EOF. For
>>> example,
>>> do n=0,2-1
>>> res@gsnLeftString = "EOF "+(n+1)
>>> *res@gsnRightString = sprintf("%5.1f", eof@pcvar(n))
>>> +"%"*
>>>
>>> plot(n)=gsn_csm_contour_map_polar(wks,ev_regres(n,:,:),res)
>>> end do
>>> gsn_panel(wks,plot,(/2,1/),resP) ; now draw as
>>> one plot
>>>
>>> This (res@gsnRightString = sprintf("%5.1f",
>>> eof@pcvar(n)) +"%") is
>>> variance explained by EOF instead of rotated EOF.
>>> When i modifiy this
>>> line by
>>> res@gsnRightString = sprintf("%5.1f",
>>> ev_rot@pcvar(n)) +"%", it gives me
>>> segmentation fault that pcvar is undefined variable
>>> or function. Here
>>> ev_rot are my rotated EOFs.
>>> Also, is there any difference between variance
>>> explained in EOF and
>>> rotated EOFs?
>>> Thanks in advance for help.
>>>
>>>
>>>
>>>
>>> On Tue, Nov 6, 2012 at 4:49 PM, Dennis Shea
>>> <shea@ucar.edu <mailto:shea@ucar.edu>
>>> <mailto:shea@ucar.edu <mailto:shea@ucar.edu>>> wrote:
>>>
>>> As noted in the eofunc_varimax documentation:
>>>
>>> "The results may be very dependent upon the user
>>> specified number of
>>> modes used in the rotation. The "best" number of
>>> modes to use may
>>> have to be determined by experiment."
>>>
>>> --
>>> Let's say the maximum possible number of
>>> eigenvalues/vectors is 1000.
>>> This means the 2nd EOF must be orthogonal to the
>>> 1st; the 3rd must
>>> be orthogonal to 1 and 2, etc. The orthogonality
>>> property impose
>>> constraints of the patterns.
>>>
>>> Even if you requested a subset of the eigen
>>> information, the LAPACK
>>> computes all the eigenvalues and returns just the
>>> subset but the
>>> full orthogonality is still in effect!
>>>
>>> If you use
>>> eof = eofunc(x,10, False)
>>> eofunc_varimax(eof)
>>>
>>> then the varimax imposes orthogonality on these
>>> 10 modes only.
>>>
>>> If you used only 2 modes .... then I'm not sure
>>> what you get
>>> but 2 modes in not 'right'
>>> eof = eofunc(x, 2, False)
>>> eofunc_varimax(eof)
>>>
>>>
>>>
>>>
>>>
>>>
>>> On 11/06/2012 03:00 PM, Waqar Younas wrote:
>>>
>>> Hi Phillips
>>> The data i used is NCEP/NCAP reanalysis of
>>> 500mb monthly
>>> geopotential
>>> height from 1950-2000.
>>> I got my results correct. The only confusion
>>> remaining is when i
>>> calculate only 2 rotational components the
>>> results are not
>>> exact, but
>>> when i calculate the 10 rotational components
>>> and plot only 2, the
>>> results are more like NAO (first component)
>>> and PNA (second
>>> component).
>>> Please find attached figure when i calculate
>>> 10 components and plot
>>> first 2 components. The other figure (when i
>>> calculate only 2
>>> components
>>> and plot them) is attached in the first
>>> message above. I do not
>>> know why
>>> there is such difference?
>>>
>>> Thanks.
>>>
>>>
>>>
>>> On Tue, Nov 6, 2012 at 1:29 PM, Adam Phillips
>>> <asphilli@ucar.edu <mailto:asphilli@ucar.edu>
>>> <mailto:asphilli@ucar.edu
>>> <mailto:asphilli@ucar.edu>>
>>> <mailto:asphilli@ucar.edu
>>> <mailto:asphilli@ucar.edu> <mailto:asphilli@ucar.edu
>>> <mailto:asphilli@ucar.edu>>>> wrote:
>>>
>>> Hi Waqar,
>>> Looking at your script, I don't see
>>> anything wrong. That
>>> leaves me
>>> trying to answer why the patterns are
>>> not what you expected
>>> them to
>>> be, and of course I cannot answer that.
>>> The patterns depend
>>> on the
>>> data you use. You can try running your
>>> script on other
>>> data, say,
>>> the NCEP/NCAR Reanalysis, to see if the
>>> results are what
>>> you expect.
>>> Good luck,
>>> Adam
>>>
>>>
>>> On 11/6/12 12:31 PM, Waqar Younas wrote:
>>>
>>> Hi
>>> It does not make any difference
>>> because even without
>>> reordering, i
>>> was getting the EOF patterns in
>>> order as first mode
>>> explained
>>> variance was 13.3% and second one is
>>> 10% and so on. The
>>> patterns
>>> did not change.
>>> I do not know where i am making
>>> mistake because
>>> especially for
>>> PNA, the most variance should be in
>>> Western side but my PNA
>>> pattern is capturing variability on
>>> Eastern side too.
>>> Now i am
>>> attaching the whole code which i
>>> modified using eof_4.ncl.
>>> Thanks.
>>>
>>>
>>>
>>>
>>>
>>> On Tue, Nov 6, 2012 at 9:05 AM,
>>> Dennis Shea
>>> <shea@ucar.edu <mailto:shea@ucar.edu>
>>> <mailto:shea@ucar.edu <mailto:shea@ucar.edu>>
>>> <mailto:shea@ucar.edu
>>> <mailto:shea@ucar.edu> <mailto:shea@ucar.edu
>>> <mailto:shea@ucar.edu>>>> wrote:
>>>
>>> The documentation for
>>> eofunc_varimax includes the
>>> following
>>>
>>> " The order of the returned
>>> varimax rotated EOFs is not
>>> necessarily in strict descending order as returned
>>> by eofunc
>>> or eofunc_varimax_Wrap. The user
>>> should use the
>>> *eofunc_varimax_reorder* procedure to place the
>>> returned
>>> varimax rotated EOFs into strict
>>> descending order."
>>>
>>> http://www.ncl.ucar.edu/__Document/Functions/__Contributed/eofunc_varimax___reorder.shtml
>>>
>>>
>>>
>>> <http://www.ncl.ucar.edu/Document/Functions/Contributed/eofunc_varimax_reorder.shtml>
>>>
>>> It seems to me you should have"
>>>
>>> ev_rot = eofunc_varimax_Wrap(eof,1)
>>> printVarSummary(ev_rot)
>>>
>>> eofunc_varimax_reorder(ev_rot)
>>> printVarSummary(ev_rot)
>>>
>>> ; ev_rot will be orthogonal.
>>> However, 'ev_rot1'
>>> ; will not tield principal
>>> components that are
>>> independent
>>>
>>>
>>>
>>> ev_rot1 = eofunc_ts_Wrap (wx, ev_rot, optETS)
>>>
>>>
>>> On 11/05/2012 06:23 PM, Waqar
>>> Younas wrote:
>>>
>>> Hi all,
>>> I need help in performing
>>> the rotated EOF of
>>> NCEP/NCAR
>>> reanalysis data
>>> using NCL. I used the code eof_4.ncl present on ncl
>>> functions page. This
>>> code use the same data i.e. 500mb monthly data
>>> set from
>>> NCEP/NCAR from
>>> 20S to 90S but this code only calculates the
>>> EOF. This code is
>>> originally written for calculating AAO
>>> oscillations.
>>> For my case, i change the region from 20-90 N
>>> and for
>>> rotated EOF, i
>>> added these lines in the code
>>>
>>> ;;;;;;; Rotation of EOF
>>>
>>> ev_rot = eofunc_varimax_Wrap(eof,1)
>>> printVarSummary(ev_rot)
>>>
>>> ;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;__;;;;;;;;;;;;;;;;;;;;
>>>
>>> ev_rot1 = eofunc_ts_Wrap (wx, ev_rot, optETS)
>>> ev_rot1 = dim_standardize_n( ev_rot1, 0, 1)
>>> ;
>>>
>>> ==============================__==============================__=====
>>> ; Regress
>>> ;
>>>
>>> ==============================__==============================__=====
>>>
>>>
>>> ev_regres = ev_rot
>>> ;
>>> create an array w
>>> meta data
>>> do ne=0,neof-1
>>> ev_regres(ne,:,:) = (/ regCoef(ev_rot1(ne,:),
>>> xAnom(lat|:,lon|:,time|:)) /)
>>> end do
>>> At the end, i plot this ev_regres which
>>> contains REOF1 and
>>> REOF2.
>>> According to many papers and definitions of Climate
>>> prediction centre,
>>> the first REOF should be NAO and second should
>>> be PNA. In
>>> my case, the
>>> first EOF is NAO but it is slightly different.
>>> The second
>>> EOF is
>>> capturing PNA but it is not exact. I am
>>> attaching the figure.
>>> For my research i need PNA pattern which is
>>> second REOF.
>>> It is capturing
>>> PNA but also some other variabilties too.
>>> Thanks for the help in advance.
>>>
>>>
>>>
>>>
>>>
>>> --
>>> Best Regards
>>>
>>> Waqar Younas
>>> PhD Candidate
>>> Natural Resources and Environmental Studies (NRES)
>>> University of Northern British Columbia (UNBC)
>>> Canada
>>>
>>>
>>>
>>> _________________________________________________
>>>
>>> ncl-talk mailing list
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/__mailman/listinfo/ncl-talk
>>>
>>>
>>> <http://mailman.ucar.edu/mailman/listinfo/ncl-talk>
>>>
>>>
>>>
>>>
>>>
>>> --
>>> Best Regards
>>>
>>> Waqar Younas
>>> PhD Candidate
>>> Natural Resources and Environmental
>>> Studies (NRES)
>>> University of Northern British
>>> Columbia (UNBC)
>>> Canada
>>>
>>>
>>>
>>> _________________________________________________
>>>
>>> ncl-talk mailing list
>>> List instructions, subscriber
>>> options, unsubscribe:
>>> http://mailman.ucar.edu/__mailman/listinfo/ncl-talk
>>>
>>> <http://mailman.ucar.edu/mailman/listinfo/ncl-talk>
>>>
>>>
>>>
>>> _________________________________________________
>>>
>>> ncl-talk mailing list
>>> List instructions, subscriber options,
>>> unsubscribe:
>>> http://mailman.ucar.edu/__mailman/listinfo/ncl-talk
>>>
>>>
>>> <http://mailman.ucar.edu/mailman/listinfo/ncl-talk>
>>>
>>>
>>>
>>>
>>> --
>>> Best Regards
>>>
>>> Waqar Younas
>>> PhD Candidate
>>> Natural Resources and Environmental Studies
>>> (NRES)
>>> University of Northern British Columbia (UNBC)
>>> Canada
>>>
>>>
>>>
>>> _________________________________________________
>>>
>>> ncl-talk mailing list
>>> List instructions, subscriber options,
>>> unsubscribe:
>>> http://mailman.ucar.edu/__mailman/listinfo/ncl-talk
>>>
>>>
>>> <http://mailman.ucar.edu/mailman/listinfo/ncl-talk>
>>>
>>>
>>>
>>>
>>>
>>> --
>>> Best Regards
>>>
>>> Waqar Younas
>>> PhD Candidate
>>> Natural Resources and Environmental Studies (NRES)
>>> University of Northern British Columbia (UNBC)
>>> Canada
>>>
>>>
>>>
>>>
>>> --
>>> Best Regards
>>>
>>> Waqar Younas
>>> PhD Candidate
>>> Natural Resources and Environmental Studies (NRES)
>>> University of Northern British Columbia (UNBC)
>>> Canada
>>>
>>> _______________________________________________
>>> ncl-talk mailing list
>>> List instructions, subscriber options, unsubscribe:
>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>>
>>
>>
>> --
>> Best Regards
>>
>> Waqar Younas
>> PhD Candidate
>> Natural Resources and Environmental Studies (NRES)
>> University of Northern British Columbia (UNBC)
>> Canada
>>
>>
>>
>>
>> --
>> Best Regards
>>
>> Waqar Younas
>> PhD Candidate
>> Natural Resources and Environmental Studies (NRES)
>> University of Northern British Columbia (UNBC)
>> Canada
>>
>>
>>
>> _______________________________________________
>> ncl-talk mailing list
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
>
>
> _______________________________________________
> ncl-talk mailing list
> List instructions, subscriber options, unsubscribe:
> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Tue Dec 4 21:27:03 2012

This archive was generated by hypermail 2.1.8 : Fri Dec 07 2012 - 13:30:06 MST