# modwtvar

Multiscale variance of maximal overlap discrete wavelet transform

## Description

example

wvar = modwtvar(w) returns unbiased estimates of the wavelet variance by scale for the maximal overlap discrete wavelet transform (MODWT). The default wavelet type is sym4.

example

wvar = modwtvar(w,wname) uses the wavelet wname to determine the number of boundary coefficients by level for unbiased estimates.

example

[wvar,wvarci] = modwtvar(___) returns the 95% confidence intervals for the variance estimates by scale.

example

[___] = modwtvar(w,wname,___,conflevel) uses conflevel for the coverage probability of the confidence interval.

example

[___] = modwtvar(w,wname,___,Name,Value,) returns wavelet variance with additional options specified by one or more Name,Value pair arguments.

example

[wvar,wvarci,nj] = modwtvar(w,wname,___) returns the number of coefficients used to form the variance and confidence intervals by level.

example

wvartable = modwtvar(w,wname,'table'), where 'table' returns a MATLAB® table, wvartable, containing the number of MODWT coefficients by level, the confidence boundaries, and the variance estimates. You can place 'table' anywhere after input w, except between the name and value of another Name,Value pair.

modwtvar(___) with no output arguments plots the wavelet variances by scale with lower and upper confidence bounds. The scaling variance is not included in the plot because the scaling variance can be much larger than the wavelet variances.

## Examples

collapse all

Obtain the MODWT of the Southern Oscillation Index data using the default symlets wavelet with 4 vanishing moments. Compute the unbiased estimates of the wavelet variance by scale.

wsoi = modwt(soi);
wvar = modwtvar(wsoi)
wvar = 10×1

0.3568
0.9026
1.1576
1.0952
0.9678
0.5478
0.6353
1.9570
0.8398
0.8247

Obtain the MODWT of the Southern Oscillation Index data using the Daubechies wavelet with 2 vanishing moments ('db2'). Compute the unbiased estimates of the wavelet variance by scale.

wsoi = modwt(soi,'db2');
wvar = modwtvar(wsoi,'db2')
wvar = 12×1

0.4296
0.9204
1.1370
1.0847
0.9255
0.5932
0.7630
1.6672
0.8048
0.7555
⋮

Obtain the MODWT of the Nile River minimum level data using the Fejer- Korovkin wavelet with eight coefficients down to level five. Use modwtvar to obtain and plot the variance estimates and 95% confidence intervals.

wtnile = modwt(nileriverminima,'fk8',5);
[wnilevar,wvarci] = modwtvar(wtnile,'fk8');

errlower = (wnilevar-wvarci(:,1));
errupper = (wvarci(:,2)-wnilevar);
errorbar(1:5,wnilevar(1:5),errlower(1:5),...
errupper(1:5),'ko','markerfacecolor','k')
hold on
title('Wavelet Variance by Scale of Nile River Levels','fontsize',14);
ylabel('Variance');
xlabel('Time (in years)');
ax = gca;
ax.XTick = [1:5];
ax.XTickLabel = {'2','4','8','16','32'};
hold off

Show how different confidence level values affect the width of the confidence intervals. An increased confidence level value increases the confidence interval width.

Obtain the MODWT of the Southern Oscillation Index data using the Fejer-Korovkin wavelet with eight coefficients.

wsoi = modwt(soi,'fk8');

Obtain the width of the.90, .95, and .99 confidence intervals for each level.

[~,wvarci90] = modwtvar(wsoi,'fk8',0.90);
w90 = wvarci90(:,2)-wvarci90(:,1);
[~,wvarci95] = modwtvar(wsoi,'fk8',0.95);
w95 = wvarci95(:,2)-wvarci95(:,1);
[~,wvarci99] = modwtvar(wsoi,'fk8',0.99);
w99 = wvarci99(:,2)-wvarci99(:,1);

Compare the three columns. The first column shows the .90 confidence level values, the second the .95 values, and the third the .99 values. Each row is the width of the interval at each wavelet scale. You can see that the width of the confidence interval increases with larger confidence level values.

[w90,w95,w99]
ans = 10×3

0.0195    0.0233    0.0306
0.0739    0.0880    0.1158
0.1347    0.1606    0.2113
0.1798    0.2145    0.2826
0.2304    0.2751    0.3634
0.1825    0.2184    0.2900
0.2858    0.3435    0.4613
1.5445    1.8757    2.5837
1.0625    1.3262    1.9551
2.8460    3.9883    7.8724

Specify non-default confidence methods using name-value pairs to compare the width of their confidence levels. Note that for Gaussian confidence level intervals, it is possible to obtain negative lower confidence bounds.

Obtain the MODWT of the Southern Oscillation Index data using the Fejer-Korovkin wavelet with eight coefficients.

wsoi = modwt(soi,'fk8');

Use the Chi2Eta and Gaussian confidence methods to obtain the variances and confidence interval bounds for each method.

[wvar_c,wvarci_c] = modwtvar(wsoi,'fk8',[],'ConfidenceMethod','chi2eta1');
[wvar_g,wvarci_g] = modwtvar(wsoi,'fk8',[],'ConfidenceMethod','gaussian');

Compute the upper and lower errors for each confidence interval and plot the results. Note that the Gaussian intervals are slightly shifted to enable better visualization.

errlower_c = wvar_c-wvarci_c(:,1);
errupper_c = wvarci_c(:,2)-wvar_c;

errlower_g = wvar_g(:,1)-wvarci_g(:,1);
errupper_g = wvarci_g(:,2)-wvar_g;

errorbar(1:10,wvar_c(1:10),errlower_c(1:10),...
errupper_c(1:10),'ko','markerfacecolor','b')
hold on;
xoffset = (1.3:10.3);
errorbar(xoffset,wvar_g(1:10),errlower_g(1:10),...
errupper_g(1:10),'ro','markerfacecolor','r')

title('Wavelet Chi2Eta2 vs. Gaussian Confidence Intervals','fontsize',14);
ylabel('Variance');
xlabel('Level')
ax = gca;
ax.XTick = [1:10];
legend('Chi2Eta','Gaussian','Location','northwest');
hold off

Compare the number of coefficients for unbiased and biased wavelet variance estimates. For the unbiased (default) estimates, the number of nonboundary coefficients decreases by scale. For biased estimates, the number of coefficients matches the number of input rows and is constant for every scale.

Obtain the MODWT of the Southern Oscillation Index data using the Fejer-Korovkin wavelet with eight coefficients. Compute the unbiased and biased estimates of the wavelet variance down to level ten. The number of coefficients used in the unbiased estimates decrease by scale.

wsoi = modwt(soi,'fk8');
[wvar_unb,wvarci_unb,nj_unb] = modwtvar(wsoi,'fk8');
[wvar_b,wvarci_b,nj_b] = modwtvar(wsoi,'fk8',[],'EstimatorType','biased');
[nj_unb(1:10),nj_b(1:10)]
ans = 10×2

12991       12998
12977       12998
12949       12998
12893       12998
12781       12998
12557       12998
12109       12998
11213       12998
9421       12998
5837       12998

Compute the MODWT of the Southern Oscillation Index data using the Fejer- Korovkin wavelet with eight coefficients. Compute a variance table for the data.

wsoi = modwt(soi,'fk8');
[wvartable] = modwtvar(wsoi,'fk8',0.90,'ConfidenceMethod','Gaussian',...
'table')
wvartable=10×4 table
NJ       Lower     Variance     Upper
_____    _______    ________    _______

D1     12991     0.3291    0.33848     0.34786
D2     12977    0.87172     0.9034     0.93508
D3     12949     1.1041     1.1628      1.2216
D4     12893     1.0204     1.0933      1.1662
D5     12781     0.8833    0.98255      1.0818
D6     12557    0.47178    0.54152     0.61125
D7     12109    0.41916    0.57934     0.73951
D8     11213    0.33639      2.055      3.7736
D9      9421     0.4752    0.83369      1.1922
D10     5837    0.37485    0.84386      1.3129

The resulting table contains the number of nonboundary coefficients, the lower and upper confidence level bounds, and the variance estimate for each level.

## Input Arguments

collapse all

MODWT transform, specified as a matrix. w is the output of modwt.

Data Types: double

Wavelet, specified as a character vector or string scalar corresponding to a valid wavelet, or as a positive even scalar indicating the length of the wavelet and scaling filters. The wavelet filter length must match the length used in the MODWT of the input.

If you use Name,Value pairs arguments or the 'table' syntax and you did not specify a wname , you must use [] as the second argument.

Confidence level, specified as a real scalar value greater than 0 and less than 1. The confidence level determines the coverage probability of the confidence intervals. The confidence levels are also shown in wvartable, if you specify 'table' as an input.

### Name-Value Pair Arguments

Specify optional comma-separated pairs of Name,Value arguments. Name is the argument name and Value is the corresponding value. Name must appear inside quotes. You can specify several name and value pair arguments in any order as Name1,Value1,...,NameN,ValueN.

Example: 'EstimatorType','biased' specifies a biased estimator.

Type of estimator used for variance estimates and confidence bounds, specified as the comma-separated pair consisting of 'EstimatorType' and one of these values.

• 'unbiased' — Unbiased estimator, which identifies and removes boundary coefficients prior to computing the variance estimates and confidence bounds. Unbiased estimates are used more frequently for wavelet variance computations.

• 'biased' — Biased estimator, which uses all coefficients to compute the variance estimates and confidence bounds.

Confidence method used to compute the confidence intervals, specified as the comma-separated pair consisting of 'ConfidenceMethod' and one of these values:

 'chi2eta3' Chi-square probability density method three, which determines the degrees of freedom.[1]. 'chi2eta1' Chi-square probability density method one, which determines the degrees of freedom [1]. 'gaussian' Gaussian method [1] . This method can result in negative lower bounds.

See Algorithm for information on each of these confidence methods.

Boundary condition used to compute the variance estimates and confidence bounds, specified as the comma-separated pair consisting of 'Boundary' and one of these values:

 'periodic' Periodic boundary handling, which does not change the original signal before computing the MODWT. If modwt uses periodic boundary handling, you must specify 'Boundary','periodic' for modwtvar to obtain a correct estimate. 'reflection' Reflection boundary handling. If the MODWT uses reflection boundary handling, you must also specify 'Boundary','reflection' for modwtvar to obtain a correct unbiased estimate. The MODWT, with reflection boundary handling, extends the original signal symmetrically at the right boundary to twice the signal length. The MODWTVAR algorithm has to know about this extended signal to calculate the correct unbiased estimate. For biased estimators, all the coefficients are used to form the variance estimates and confidence intervals regardless of the boundary handling.

## Output Arguments

collapse all

Wavelet variance estimates, returned as vector. The number of elements in wvar depends on the number of scales in the input matrix and, for unbiased estimates, on the wavelet length. For the unbiased case, modwtvar returns estimates only where nonboundary coefficients exist. This condition is satisfied when the transform level is not greater than floor(log2(N/(L-1)+1)), where N is the input signal length and L is the length of the wavelet filter. The number of biased estimates equals the input signal length. If the final level has sufficient nonboundary coefficients, modwtvar returns the scaling variance in the final element of wvar.

Confidence bounds, expressed as upper and lower confidence bounds, for the variance estimates in wvar, returned as a matrix. The default is 95% confidence bounds, but you can use a different value using the conflevel input argument. The confidence bounds matrix is M-by-2, where M is the number of levels. For unbiased estimates, the number of levels is limited by the number of nonboundary coefficients. For biased estimates, all levels are used. The first column of the confidence interval matrix contains the lower confidence bound and the second column contains the upper confidence bound. By default, modwtvar calculates the confidence intervals using the chi-square probability density, with the equivalent degrees of freedom estimated using the 'Chi2Eta3' confidence method.

Number of nonboundary coefficients by scale, returned as a vector. For unbiased estimates, nj is the number of nonboundary coefficients and decreases by level. For biased estimates, nj is a vector of constants equal to the number of columns in the input matrix.

Variance table, returned as a MATLAB table. The four variables in the table are:

• NJ — Number of MODWT coefficients by level. For biased estimates, NJ is the number of coefficients in the MODWT. For unbiased estimates, NJ is the number of nonboundary coefficients.

• Lower — Lower confidence bound for the variance estimate.

• Variance — Variance estimate by level.

• Upper — Upper confidence bound for the variance estimate.

The row names of wvartable indicate the type and level of each estimate. For example, D1 indicates that the row corresponds to a wavelet or detail estimate at level 1. S6 indicates that the row corresponds to the scaling estimate at level 6. The scaling variance is computed for the final level of the MODWT. For unbiased estimates, modwtvar computes the scaling variance only when nonboundary scaling coefficients exist.

## Algorithms

The following expressions define the variance and confidence methods used in the MODWTVAR. The variables are

• Nj — Number of coefficients at level j

• v2 — Variance

• j — Level

• Wj,t — Wavelet coefficients

The variance estimate is

${\stackrel{^}{v}}_{j}^{2}=\frac{1}{{N}_{j}}\sum _{t=0}^{{N}_{j}-1}{W}_{j,t}^{2}$

The degrees of freedom for the Chi2Eta1 (chi2eta1) method are defined as

${\eta }_{1}=\frac{{N}_{j}{\stackrel{^}{v}}_{j}^{4}}{{\stackrel{^}{A}}_{j}}$

where

${\stackrel{^}{A}}_{j}=\frac{1}{2}{\underset{-1/2}{\overset{1/2}{\int }}\left[{\stackrel{^}{S}}_{j}^{\left(p\right)}\left(f\right)\right]}^{2}df$

In this equation, ${\stackrel{^}{S}}_{j}^{\left(p\right)}$ is the spectral density function estimate of the wavelet coefficients at level j.

The chi-square statistic is

$\frac{{\eta }_{1}{N}_{j}{\stackrel{^}{v}}_{j}^{2}}{{v}_{j}^{2}}~{Χ}_{{\eta }_{1}}^{2}$

The degrees of freedom for the Chi2Eta3 (chi2eta3) method are defined as

${\eta }_{3}=\text{max}\left(\frac{{N}_{j}}{{2}^{j}},1\right)$

The chi-square statistic is

$\frac{{\eta }_{3}{N}_{j}{\stackrel{^}{v}}_{j}^{2}}{{v}_{j}^{2}}~{Χ}_{{\eta }_{3}}^{2}$

For the Gaussian method, the statistic

$\frac{{N}_{j}^{1/2}\left(\left({\stackrel{^}{v}}_{j}^{2}-{v}_{j}^{2}\right)\right)}{{\left(2{\stackrel{^}{A}}_{j}\right)}^{1/2}}$

is distributed as N(0,1). The variable ${\stackrel{^}{A}}_{j}$ is as described for chi2eta1.

## References

[1] Percival, D. B., and A. T. Walden. Wavelet Methods for Time Series Analysis. Cambridge, UK: Cambridge University Press, 2000.

[2] Percival, D. B., D. Mondal, "A Wavelet Variance Primer." Handbook of Statistics, Volume. 300, Time Series Analysis: Methods and Applications, (T. S. Rao, S. S. Rao, and C. R. Rao, eds.). Oxford, UK: Elsevier, 2012, pp. 623–658.

[3] Cornish, C. R., C. S. Bretherton, and D. B. Percival. "Maximal Overlap Wavelet Statistical Analysis with Application to Atmospheric Turbulence." Boundary-Layer Meteorology. Vol. 119, Number 2, 2005, pp. 339–374.