User manual

Why Multitaper?

The multitaper method is a widely used set of methods for nonparametric estimation of frequency domain quantities relevant to time series such as the power spectrum, coherence, and transfer function. For an accessible introduction to multitaper methods, we recommend (Park, Lindberg, and Vernon, 1987).

Structs

Multitaper.jl outputs spectra, coherences, transfer functions, etc. in the form of Julia structs. The main reason for this is in the ease of Plotting.

The basic frequency domain structs are

where MT signifies that these are multitaper estimates. Under the hood, one will find the above structs may contain the following

Finally, one may be interested in computing the following time domain quantities

Below we expand on these basic structs.

MTSpectrum

The MTSpectrum struct contains a multitaper spectrum estimate which consists of the following fields:

  • f is the range of frequencies
  • S contains the spectrum estimate
  • phase contains the vector of phases, if this is a cross-spectrum, and otherwise nothing.
  • params is a MTParameters struct which stores the choice of bandwidth and other quantities for the multitaper.
  • coef is a EigenCoefficient struct which stores the intermediate values of the calculation, if desired. These may be useful for subsequent calculations.
  • Fpval is the vector of F-test p-values at every frequency, see

(Thomson, 1982).

  • jkvar is the vector of jackknifed variance estimates at every frequency.
  • Tsq_pval contains the result of a $T^2$ test for multiple lines, if requested. This field otherwise contains nothing.

The MTSpectrum struct is generated by the function multispec when called on a time series

julia> using Multitaper, Plots

julia> x = rand(100); # A time series consisting of 100 random points

julia> S = multispec(x); # Compute the multitaper spectrum of the data using NW = 4

julia> typeof(S)
MTSpectrum{Nothing,Nothing,Nothing}

julia> plot(S)
Plot{Plots.GRBackend() n=3}

There is a plots recipe, see Plotting, for plotting structs of this type.

There is an IJulia.jl notebook in Multitaper.jl/Examples/01_basic_multitaper.ipynb showing an example taken from the excellent text (Percival and Walden, 1993). To run this notebook on a cloned copy of the repo, follow the instructions in Examples.

The multispec function, when run on a single time series, additionally allows one to compute the following quantities:

  • F-test p-value, by setting the kwarg Fpval = true. This tests for a spectral line at every frequency.
  • T-squared test p-value (Thomson, 2011)

An example of the former are available in the Multitaper.jl/Examples/01_basic_multitaper.ipynb notebook.

The main advantage of the multitaper spectrum is that one can simultaneously control the bias and variance of the estimator by varying the time-bandwidth product, $NW$. This parameter is used to compute the multiple orthogonal tapers, see dpss_tapers, on which the method depends. When one computes the multitaper spectrum, this parameter, along with a number of other quantities is given in the output in the form of a MTParameters struct.

MTParameters

This struct holds the following quantities

comments on optimal choice of $NW$.

  • K which is the number of tapers. One typically chooses $K <= 2NW$ as the eigenvalue problem that gives the tapers has only $2NW$ large eigenvalues. For details on discrete prolate spheroidal sequences, which are used as data tapers, consult dpss_tapers and (Slepian, 1978).
  • N which is the number of data points in the time series.
  • dt is the sampling interval, typically in seconds.
  • M is the length of the output spectrum, usually different from $N$ because of zero-padding.
  • nsegments is the number of data blocks that have been made from the time series. When one calls multispec, nsegments is one, but if one has called welch, one can control the number of segments using the keyword argument of the same name.
  • overlap is the proportion of data by which the segments overlap, which is set to nothing if there is only one segment, but by default is 0.5 if one uses the welch method.

These quantities are important when plotting, as they give indicators of the resolution of the estimate, as indicated by a small cross-shape in the plotted spectrum. See notebooks in Examples and Plotting for details.

Finally, if one selects the value true for the keyword argument guts, then the intermediate values of the multitaper calculation are returned as well. These intermediate values are eigencoefficients.

EigenCoefficient

The struct containing eigencoefficients has two fields

  • coef the complex-valued tapered Fourier transformed data, and
  • wts which are the result of an adaptive weighting calculation, which are only computed if a_weight is set to true in a call to multispec.

(See (Thomson, 1982) for the optional adaptive weighting procedure).

For details of the additional keyword arguments relevant to multispec, see multispec entry in index.

MTCoherence

The MTCoherence struct contains a multitaper spectrum estimate which consists of the following fields:

  • f is the range of frequencies
  • coh contains the coherence estimate, in transformed units, i.e. transformed by the function atanhtransf. One can convert back to original squared coherence using the tanhtrans function.
  • phase contains the vector of phases. Phase is by default unwrapped, which allows values outside of [-180,180) degrees.
  • params is a MTParameters struct which stores the choice of bandwidth and other quantities for the multitaper.
  • coef is a EigenCoefficient struct which stores the intermediate values of the calculation, if desired. These may be useful for subsequent calculations.
  • jkvar is the vector of jackknifed variance estimates at every frequency.
  • Tsq_pval contains the result of a $T^2$ test for multiple lines, if requested. This field otherwise contains nothing.

The MTCoherence struct is generated by the function multispec when called on two time series with the default output option.

There is a plots recipe, see Plotting, for plotting structs of this type.

There is an IJulia.jl notebook in Multitaper.jl/Examples/02_multivariate.ipynb showing an example taken from the excellent (and freely available) text (Shumway and Stoffer). To run this notebook on a cloned copy of the repo, follow the instructions in Examples.

MTTransf

The MTTransf struct contains a multitaper transfer function estimate (Thomson and Chave, 1991) which consists of the following fields:

  • f is the range of frequencies
  • transf contains the transfer function estimate, in transformed units, i.e. transformed by the function atanhtransf. One can convert back to original squared coherence using the tanhtrans function.
  • phase contains the vector of phases. Phase is by default unwrapped, which allows values outside of [-180,180) degrees.
  • params is a MTParameters struct which stores the choice of bandwidth and other quantities for the multitaper.
  • coef is a EigenCoefficient struct which stores the intermediate values of the calculation, if desired. These may be useful for subsequent calculations.
  • jkvar is the vector of jackknifed variance estimates at every frequency.

The MTTransf struct is generated by the function multispec when called on two time series with the output option set to :transf.

Time Domain Structs

Finally, we come to the time domain structs, which are indexed by lag.

Demodulate

The Demodulate struct contains a multitaper complex demodulate estimate which consists of the following fields:

  • time is the range of times over which the time series is collected.
  • mag contains the vector of magnitudes of the complex demodulate.
  • phase contains the vector of phases. Phase is by default unwrapped, which allows values outside of [-180,180) degrees.

The Demodulate struct is generated by the function demodulate when called on the time series.

There is a plots recipe, see Plotting, for plotting structs of this type.

There is an IJulia.jl notebook in Multitaper.jl/Examples/04_Demodulation.ipynb showing an example taken from the paper (Thomson, 1995) and R multitaper. To run this notebook on a cloned copy of the repo, follow the instructions in Examples. The authors acknowledge the work of K. Rahim and W. Burr, authors of the aforementioned R package, for the function template and example in docs.

MTAutocorrelationFunction

The MTAutocorrelationFunction struct contains a multitaper estimate of autocorrelation which consists of the following fields:

  • lags is the range of times over which the time series is collected.
  • acf contains the vector of autocorrelations computed on the time series.
  • params contains the MTParameters struct containing the multitaper bandwidth and other parameters

The MTAutocorrelationFunction struct is generated by the function mt_acf when called on the time series.

There is a plots recipe, see Plotting, for plotting structs of this type.

MTAutocovarianceFunction

The MTAutocovarianceFunction struct contains a multitaper estimate of autocovariance which consists of the following fields:

  • lags is the range of times over which the time series is collected.
  • acvf contains the vector of autocovariances computed on the time series.
  • params contains the MTParameters struct containing the multitaper bandwidth and other parameters

The MTAutocorrelationFunction struct is generated by the function mt_acvf when called on the time series.

There is a plots recipe, see Plotting, for plotting structs of this type.

MTCepstrum

The MTCepstrum struct contains a multitaper cepstrum estimate which consists of the following fields:

  • time is the range of times over which the time series is collected.
  • ceps contains the vector of cepstrum coefficients.
  • phase contains the vector of phases. Phase is by default unwrapped, which allows values outside of [-180,180) degrees.

The MTCepstrum struct is generated by the function mt_cepstrum when called on the time series.

MtCrossCovarianceFunction

The MTCrossCovarianceFunction struct contains a multitaper estimate of cross covariance which consists of the following fields:

  • lags is the range of times over which the time series is collected.
  • ccvf contains the vector of cross covariances computed on the two time series.
  • params contains the MTParameters struct containing the multitaper bandwidth and other parameters

The MTCrossCovarianceFunction struct is generated by the function mt_ccvf when called on the two time series.

There is a plots recipe, see Plotting, for plotting structs of this type.

MTCrossCorrelationFunction

The MTCrossCorrelationFunction struct contains a multitaper estimate of cross correlation which consists of the following fields:

  • lags is the range of times over which the time series is collected.
  • ccf contains the vector of cross correlations computed on the two time series.
  • params contains the MTParameters struct containing the multitaper bandwidth and other parameters

The MTCrossCorrelationFunction struct is generated by the function mt_ccvf when called on the two time series.

There is a plots recipe, see Plotting, for plotting structs of this type.

Functions creating the structs

Now that the basic structs have been defined, we can begin to describe the functions that create them. See Function index for all of these.