# constructors of the waveletVar class -----------------------------------------------
# Build a waveletVar object from a numeric vector (private function)
#
# @param vpr Numeric vector containing the wavelet coefficients' variances on
# each wavelet resolution level (resolution levels are assumed to be arranged
# in ascending order).
# @param family Specifies the family of wavelets used in the computation of the
# wavelet transform. See \code{\link[wavethresh]{wd}}.
# @param filter_number An integer identifying the concrete filter used within the
# \code{family} of wavelets selected. \code{filter_number} is related with the
# smoothness of the wavelet. See \code{\link[wavethresh]{wd}}.
buildWaveletVar = function(vpr, family, filter_number) {
names(vpr) = 0:(length(vpr) - 1)
class(vpr) = c("waveletVar", class(vpr))
attr(vpr, "wt_info") = list(family = family,
filter_number = filter_number)
vpr
}
#' Wavelet coefficients' variances.
#'
#' \code{waveletVar} computes the wavelet coefficients' variances of an object
#' representing a wavelet transform in each of its resolution levels.
#'
#' In addition to the specific methods for the \code{waveletVar} class, all the
#' methods for numeric vectors can be used with a \code{waveletVar} object.
#'
#' @param x A \code{wd} object (see \code{\link[wavethresh]{wd}}). The wavelet
#' coefficients' variances are computed in each of the resolution levels
#' of \code{wd}.
#' @return A S3 \code{waveletVar} object that stores the
#' Wavelet coefficients' variances depending on the resolution level.
#' @export
#' @examples
#' use_H = 0.3
#' fbm = fbmSim(n = 2 ^ 10, H = use_H)
#' w_fbm = wd(fbm, bc = "symmetric")
#' vpr = waveletVar(w_fbm)
#' plot(vpr)
#' # the theoreticalWaveletVar also returns a 'waveletVar' object:
#' theo_vpr = theoreticalWaveletVar(H = use_H, sigma2 = 1,
#' family = wtInfo(vpr)[["family"]],
#' filter_number = wtInfo(vpr)[["filter_number"]],
#' nlevels = length(vpr))
#' points(theo_vpr, col = "red")
#' @seealso \code{\link{wtInfo}}, \code{\link{resolutionLevels}},
#' \code{\link{theoreticalWaveletVar}}, \code{\link{fracdet}},
#' \code{\link{estimateFbmPars}}
#' @import wavethresh
#' @export
waveletVar = function(x){
kMIN_POINTS = 5
if (!inherits(x, "wd")) {
stop("x is not an \"wd\" object")
}
x_nlevels = wavethresh::nlevelsWT(x)
vpr = rep(0.0, x_nlevels)
filter_len = length(
wavethresh::filter.select(
filter.number = x$filter$filter.number,
family = x$filter$family
)$H
)
# first index unaffected by the circularity assumption of the wavelet transform
unaffected_index = ceiling( (filter_len - 2) * (1 - 1 / 2 ^ (x_nlevels:1)) )
# estimate the variance of the resolution level 0 ...
vpr[[1]] = wavethresh::accessD(x, 0) ^ 2
# ... and loop for the remainder levels
resolution_levels = 0:(x_nlevels - 1)
for (i in 2:x_nlevels) {
available_points = (2 ^ resolution_levels[[i]] - 2 * unaffected_index[[i]])
if (available_points > kMIN_POINTS) {
segment_index = unaffected_index[[i]]:(2 ^ resolution_levels[[i]] -
unaffected_index[[i]])
vpr[[i]] = var(wavethresh::accessD(x,
resolution_levels[[i]])[segment_index])
}else {
# not enough points for a variance-estimation free from the circularity
# assumption: we use all points of the level
vpr[[i]] = var(wavethresh::accessD(x, resolution_levels[[i]]))
}
}
buildWaveletVar(vpr,
family = x$filter$family,
filter_number = x$filter$filter.number)
}
# Accesing attributes of the waveletVar class ------------------------------------
#' Get the resolution levels at which the wavelet coefficients' variances was
#' computed.
#' @param x A \code{waveletVar} object.
#' @return A numerical vector of the same length as \code{waveletVar}.
#' @export
resolutionLevels = function(x){
UseMethod("resolutionLevels")
}
#' @export
resolutionLevels.waveletVar = function(x){
0:(length(x) - 1)
}
#' Get the family and the filter-number used in the wavelet transform asociated
#' with the \code{waveletVar} object.
#' @inheritParams resolutionLevels
#' @return An R list containing the family (\code{family} field) and the
#' filter-number (\code{filter_number}).
#' @export
wtInfo = function(x){
UseMethod("wtInfo")
}
#' @export
wtInfo.waveletVar = function(x){
attr(x, "wt_info")
}
# Plotting functions for waveletVar class ----------------------------------------
#' @export
plot.waveletVar = function(x, main = "Variance per Resolution Level",
xlab = "Wavelet Resolution Level",
ylab = "Variance", log = "y", ...){
plot(resolutionLevels(x), as.numeric(x), main = main, xlab = xlab, ylab = ylab,
log = log, ...)
}
#' @export
points.waveletVar = function(x, ...) {
points(resolutionLevels(x), as.numeric(x), ...)
}
#' @export
lines.waveletVar = function(x, ...){
lines(resolutionLevels(x), as.numeric(x), ...)
}
# Theoretical expression of waveletVar -----------------------------------------------
# private function
# get wavelet filters
getWaveletFilters = function(family, filter_number){
h = wavethresh::filter.select(filter.number = filter_number,
family = family)$H
len = length(h)
# Definition of g[n]:
# g[n] = (-1) ^ (1 - n) * h[1 - n]
index = 1 - (len - 1):0
g = rev(h) * (-1) ^ (1 - index)
list(h = h, g = g)
}
#' Wavelet coefficients' variances of a fractional Brownian motion
#'
#' \code{theoreticalWaveletVar} calculates the theoretical wavelet coefficients'
#' variance of a fractional Brownian motion (fBm) with known parameters.
#'
#' @param H The Hurst exponent of the fBm.
#' @param sigma2 Variance of the fractional Gaussian noise that results after
#' differentiating the fBm.
#' @param family Specifies the family of wavelets used in the computation of the
#' wavelet transform. See \code{\link[wavethresh]{wd}}.
#' @param filter_number An integer identifying the concrete filter used within the
#' \code{family} of wavelets selected. \code{filter_number} is related with the
#' smoothness of the wavelet. See \code{\link[wavethresh]{wd}}.
#' @param nlevels Number of resolution levels to compute. Thus, the resulting
#' wavelet coefficients' variances are the expected ones from a fBm of length
#' \eqn{2 ^ n}.
#' @return A \code{waveletVar} object with the theoretical wavelet variance.
#' @examples
#' use_H = 0.3
#' fbm = fbmSim(n = 2 ^ 10, H = use_H)
#' w_fbm = wd(fbm, bc = "symmetric")
#' vpr = waveletVar(w_fbm)
#' plot(vpr)
#' # the theoreticalWaveletVar also returns a 'waveletVar' object:
#' theo_vpr = theoreticalWaveletVar(H = use_H, sigma2 = 1,
#' family = wtInfo(vpr)[["family"]],
#' filter_number = wtInfo(vpr)[["filter_number"]],
#' nlevels = length(vpr))
#' points(theo_vpr, col = "red")
#' @seealso \code{\link{waveletVar}}, \code{\link{estimateFbmPars}}
#' @export
#' @import Rcpp
#' @useDynLib fracdet
theoreticalWaveletVar = function(H, sigma2, family, filter_number, nlevels) {
filters = getWaveletFilters(family = family,
filter_number = filter_number)
buildWaveletVar(.Call('fracdet_theoreticalWaveletVarCpp', PACKAGE = 'fracdet',
filters$g, filters$h, H, sigma2, nlevels),
family = family,
filter_number = filter_number)
}
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.