SMAC-Group / wv

:alarm_clock: This R package provides the tools to perform standard and robust wavelet variance analysis for time series (signal processing). Among others, aside from computing the wavelet variance and cross-covariance (classic and robust), the package provides inference tools (e.g. confidence intervals) and plotting tools allowing to perform some visual analysis and assess the properties of the underlying time series.
https://smac-group.github.io/wv/
15 stars 10 forks source link

bug theoretical wavelet variance AR1() #11

Open lionelvoirol opened 4 years ago

lionelvoirol commented 4 years ago

Issue 1)⚠️ When computing the theoretical wavelet variance of an AR1() process, one can see that if the phi parameter is smaller than 0 (but still contained in the bounds [-1,1]), the wavelet variance is undefined (NaN) for the odd scales.

Reproduce this problem with the following code: wv::ar1_to_wv(phi = -.5, sigma2 = 1, tau = 1:20)

Issue 2) ⚠️The theoretical wavelet variance of the AR1 process is not calculated correctly. One can indeed observe that the theoretical and empirical wavelet variance do not match when the correct parameters are given as input to wv::ar1_to_wv.

Reproduce this problem with the following code: y = gen_gts(model = AR1(phi = .7, sigma2 = 2), n = 100000)

mywv = wv::wvar(y)

plot(x = mywv$scales, mywv$variance, log = "xy", type ="b", col = "blue4")

wv_theo = wv::ar1_to_wv(phi = .7, sigma2 = 2, tau = 1:15)

lines(y= wv_theo, x = mywv$scales, col= "darkorange", type ="b")

Compare this to the following correct function where the wavelet is computing according to the SMAC Group website

perso_ar1_to_wv = function(phi, sigma2, tau){
  all_tau = 2^tau
  (sigma2*((phi^2-1)*all_tau + 2*phi*(phi^all_tau-4*phi^(all_tau/2)+3)) ) / ((phi-1)^3*(phi+1)*all_tau^2)

}