sowas tutorial - Matlab version


sowas tutorial - Matlab version Next: , Previous: (dir), Up: (dir)

This is a short tutorial for the sowas package.


Next: , Previous: Top, Up: Top

Introduction

This tutorial should give you some easy and illustrative examples to get used to the sowas package.

The sowas package provides functions for the estimation (i.e. analysis) of wavelet spectra, cross spectra and coherence. Furthermore, it provides an algorithm to simulate surrogate data (i.e. synthesis) for every given wavelet spectrum.

The user should initially confine himself to using the functions wsp() for estimating the wavelet spectrum, wco() for estimating the coherence, plotwt() for plotting the results and surrwave() for creating surrogate data. All the other functions are called by the former functions and should be used with great care.


Next: , Previous: Introduction, Up: Top

References

When using this software, please cite the following two papers, the first for the areawise test or the surrogate data, the second for the pitfalls in general:

D. Maraun, J. Kurths and M. Holschneider: Nonstationary Gaussian Processes in Wavelet Domain: Synthesis, Estimation and Significance Testing. Phys. Rev. E 75, 016707, 2007,pdf, 738 kB

D. Maraun and J. Kurths: Cross Wavelet Analysis. Significance Testing and Pitfalls, Nonlin. Proc. Geoph. 11(4), 505-514, 2004, pdf, 282 kB, clarifications


Next: , Previous: References, Up: Top

Data Format

sowas in Matlab uses univariate time series as input, i.e. a two column matrix with times in the first column and values in the second column. If only one column is provided, this will automatically be interpreted as a vector of values with sampling time 1 and starting time 1.

Matlab provides the function load()to read a two-column ascii-file into such a matrix. Three data sets are provided with the package.


Next: , Previous: Data Format, Up: Top

Wavelet Sample Spectrum

To get a feeling, how sowas works, load the NINO3 time series as an example:

load('data/nino3.dat')

Then type

n3wsp=wsp(nino3,0.5,-1,-1,-1,-1,-1,-1,-1,-1,0)

The values set to -1 are set to the default value. This call calculates the wavelet periodogram (i.e. the non-averaged wavelet sample spectrum) of the NINO3-index. By setting nreal=0, no significance test is performed. Try now

n3wsp=wsp(nino3,0.5,-1,-1,-1,-1,-1,-1,-1,-1,100,0)

Now the pointwise significance test is performed additionally. You can see a lot of small patches which are likely to be false positive artefacts. Apply the areawise test to identify them:

n3wsp=wsp(nino3,0.5,-1,-1,-1,-1,-1,-1,-1,-1,100)

Higlight the El Nino frequency band by plotting the result with additional lines:

plotwt(n3wsp,-1,[2,7])

To understand, why the pointwise significance test is producing false positive results, create a white noise time series

wgn=randn(1920,2)
wgn(1:1920,1)=1850:1/12:(2010-1/12)

The time series has no structure:

plot(wgn(1:length(wgn),1),wgn(1:length(wgn),2))

Now apply both the pointwise test and the areawise test to the wavelet sample spectrum

wgnwsp=wsp(wgn,0.5,-1,-1,-1,-1,-1,-1,-1,-1,100)

Even though the data have no preferred frequencies, i.e. no structure, false positive patches occur, that are mostly identified by the areawise test.


Next: , Previous: Wavelet Sample Spectrum, Up: Top

Wavelet Sample Cross Spectrum

Now, as an example for a cross analysis, additionally load the AIR index:

load('data/air.dat')

Estimate the wavelet cross spectrum:

n3airwcs=wcs(nino3,air,0.5)

Now, as a comparison try the sample cross spectrum between NINO3 and Gaussian white noise. Type

n3wgnwcs=wcs(nino3,wgn,0.5)

Also here, many peaks related to El Nino events should occur (in case not, the white noise realization might be quite untypical), even though both processes are completely independent. Thus, the wavelet sample cross spectrum provides no reliable measure for a relationship between two processes. The reason for this is that the cross spectrum (like the covariance) is not normalized. The estimation might exhibit large values simply because of strong power in one of the (single) wavelet sample spectra.


Next: , Previous: Wavelet Sample Cross Spectrum, Up: Top

Wavelet Sample Coherence

To avoid the misleading results of the wavelet sample cross spectrum, use the wavelet sample coherency:

n3airwco=wco(nino3,air,0.5,-1,-1,-1,0,0)

This plot should be completely red, this result is for sure wrong. The reason is that coherency needs to be averaged over time and scale, otherwise you would try to estimate similarities of periodicities from single points, which is impossible. Thus, try

nino3airwco=wco(nino3,air,0.5)

This estimates the coherence by smoothing one octave (sw=0.5 in each direction) in scale direction and 3 periods in time direction on every scale (tw=1.5 periods in every direction). The result includes both significance tests. You still see relations between both processes around 1900-1920 and 1960-1980. For a comparison, try now the sample coherence between NINO3 and white noise:

nino3wgnwco=wco(nino3,wgn,0.5)

This result should show no relations between both processes (unless the realization of the white noise is very weird). Thus, wavelet coherence proves to be a much more reliable measure for the inference of interrelations.


Next: , Previous: Wavelet Sample Coherence, Up: Top

Functions

see the online documentation of sowas