pycwt icon indicating copy to clipboard operation
pycwt copied to clipboard

the dof for the global wavelet spectrum will change if scales are in days instead of years?

Open JeancarloFU opened this issue 3 years ago • 2 comments

According to the tutorial sample.py, for the global wavelet spectrum the input parameter "dof" for "wavelet.significance(...,dof=dof)" (lines 91-92) is calculated using dof = N - scales, where N is the size of the input data and the parameter scales has units of years, then this "dof" value is used as "na" in eq (23) of Torrence and Compo (1998). However, if I want to use scales in days for the same problem the dof will change and will be negative when scales > N (N=504 for Nino3 example), which is not consistent.

According to Torrence and Compo (1998) (page 73), they mention that for the global wavelet if the points going into the average of eq (22) are within the cone of influence, then the number of effective points "na" of eq (23) is reduced by approximately one-half of the number within the COI to reflect the decreased amplitude (and information) within that region. So it means that the effective number of points "na" of eq (23) should be: N - the half of the points inside the COI at each scale? or it might be better to use N - all the points inside the COI at each scale? If so, when we are at the limit of the larger scales, in which all the points are almost inside the COI, na=N-N/2=N/2 for the former case and na=N-N=0 and for the latter one.

In conclusion, what approximation should we use to get the most realistic "na" for eq (23) that does not depend on the units of the parameter scales?

Thanks. Jeancarlo Fajardo U.

JeancarloFU avatar Feb 27 '21 01:02 JeancarloFU

Hi Jeancarlo,

Similar as with FFT, you need to keep in mind that N and the scales do not have units. If you want to calculate the equivalent oscillation frequency (or period), you need to use the sampling frequency dt. The whole spectral analysis is scaled by the sampling frequency.

Let's suppose you have daily samples of a given quantity for let's say 10-years, your timeseries will have about 3650 data points (N). With more data points you can actually increase the spectral resolution.

But, if you just want to change the unit of your sampling frequency, for example from months to days, and keep the number of data points N, you will need to change dt: from dt=1 (months) to dt=30 (days).

regeirk avatar Feb 28 '21 23:02 regeirk

Dear Sebastian,

You are right, the wavelet analysis doesn't depend on the units of dt. In line (83) inside the function cwt of the file wavelet.py, scales is defined by sj = s0 * 2 ** (np.arange(0, J + 1) * dj), where s0=2*dt (so it seems that sj has units of dt?), then in all the wavelet analysis sj/dt is always used and as you mentioned it doesn't have units. Therefore, in line (91) of the file sample.py, where dof = N - scales, scales should be divided by dt because it depends on s0 and hence on dt? Then this dof, that is actually the parameter "number of points average = na" of equation (23) of Torrence and Compo (1998), is used in the line (278) inside the function significance of the file wavelet.py, which also has dependence on the dimensionless quantity sj/dt.

Luckily, the above issue only affects the return parameter signif of the function significance (in wavelet.py) for the case of the global spectrum.

Finally, I think that dof = na = N - scales might not be appropriated because at the limit of the larger scales, in which all the points average for the global wavelet are mostly inside the COI, na is still in the order of N.

Thanks. Jeancarlo

JeancarloFU avatar Mar 01 '21 09:03 JeancarloFU