Nothing
R/RcppExports.R
In gmwm: Generalized Method of Wavelet Moments
# This file was generated by Rcpp::compileAttributes
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393
#' @title Expected value DR
#' @description This function computes the expected value of a drift process.
#' @param omega A \code{double} corresponding to variance of drift.
#' @param n_ts An \code{int} indicating the length of the time series.
#' @return A \code{vec} containing the expected value of the drift.
#' @keywords internal
#' @examples
#' e_drift(1,200)
e_drift <- function(omega, n_ts) {
.Call('gmwm_e_drift', PACKAGE = 'gmwm', omega, n_ts)
}
#' @title Second moment DR
#' @description This function computes the second moment of a drift process.
#' @param omega A \code{double} corresponding to variance of drift.
#' @param n_ts An \code{int} indicating the length of the time series.
#' @return A \code{vec} containing the second moment of the drift.
#' @keywords internal
#' @examples
#' m2_drift(1, 200)
m2_drift <- function(omega, n_ts) {
.Call('gmwm_m2_drift', PACKAGE = 'gmwm', omega, n_ts)
}
#' @title Variance DR
#' @description This function computes the variance of a drift process.
#' @param omega A \code{double} corresponding to variance of drift.
#' @param n_ts An \code{int} indicating the length of the time series.
#' @return A \code{vec} containing the variance of the drift.
#' @keywords internal
#' @examples
#' var_drift(1, 200)
var_drift <- function(omega, n_ts) {
.Call('gmwm_var_drift', PACKAGE = 'gmwm', omega, n_ts)
}
#' @title Compute Tau-Overlap Allan Variance
#' @description Computation of Tau-Overlap Allan Variance
#' @usage avar_to_cpp(x)
#' @param x A \code{vector} with dimensions N x 1.
#' @return av A \code{matrix} that contains:
#' \itemize{
#' \item{Col 1}{The size of the cluster}
#' \item{Col 2}{The Allan variance}
#' \item{Col 3}{The error associated with the variance estimation.}
#' }
#' @details
#' Given \eqn{N} equally spaced samples with averaging time \eqn{\tau = n\tau _0}{tau = n*tau_0},
#' where \eqn{n} is an integer such that \eqn{ 1 \le n \le \frac{N}{2}}{1<= n <= N/2}.
#' Therefore, \eqn{n} is able to be selected from \eqn{\left\{ {n|n < \left\lfloor {{{\log }_2}\left( N \right)} \right\rfloor } \right\}}{{n|n< floor(log2(N))}}
#' Then, a sampling of \eqn{m = \left\lfloor {\frac{{N - 1}}{n}} \right\rfloor - 1} samples exist.
#' The tau-overlap estimator is given by:
#'
#' where \eqn{ {{\bar y}_t}\left( \tau \right) = \frac{1}{\tau }\sum\limits_{i = 0}^{\tau - 1} {{{\bar y}_{t - i}}} }.
#'
#' @author JJB
#' @references Long-Memory Processes, the Allan Variance and Wavelets, D. B. Percival and P. Guttorp
#' @examples
#' set.seed(999)
#' # Simulate white noise (P 1) with sigma^2 = 4
#' N = 100000
#' white.noise = rnorm(N, 0, 2)
#' #plot(white.noise,ylab="Simulated white noise process",xlab="Time",type="o")
#' #Simulate random walk (P 4)
#' random.walk = cumsum(0.1*rnorm(N, 0, 2))
#' combined.ts = white.noise+random.walk
#' av_mat = avar_to_cpp(combined.ts)
#' @keywords internal
avar_to_cpp <- function(x) {
.Call('gmwm_avar_to_cpp', PACKAGE = 'gmwm', x)
}
#' @title Compute Maximal-Overlap Allan Variance using Means
#' @description Computation of Maximal-Overlap Allan Variance
#' @usage avar_mo_cpp(x)
#' @param x A \code{vector} with dimensions N x 1.
#' @return av A \code{list} that contains:
#' \itemize{
#' \item{"clusters"}{The size of the cluster}
#' \item{"allan"}{The Allan variance}
#' \item{"errors"}{The error associated with the variance estimation.}
#' }
#' @details
#' Given \eqn{N} equally spaced samples with averaging time \eqn{\tau = n\tau _0}{tau = n*tau_0},
#' where \eqn{n} is an integer such that \eqn{ 1 \le n \le \frac{N}{2}}{1<= n <= N/2}.
#' Therefore, \eqn{n} is able to be selected from \eqn{\left\{ {n|n < \left\lfloor {{{\log }_2}\left( N \right)} \right\rfloor } \right\}}{{n|n< floor(log2(N))}}
#' Then, \eqn{M = N - 2n} samples exist.
#' The Maximal-overlap estimator is given by:
#' \eqn{\frac{1}{{2\left( {N - 2k + 1} \right)}}\sum\limits_{t = 2k}^N {{{\left[ {{{\bar Y}_t}\left( k \right) - {{\bar Y}_{t - k}}\left( k \right)} \right]}^2}} }
#'
#' where \eqn{ {{\bar y}_t}\left( \tau \right) = \frac{1}{\tau }\sum\limits_{i = 0}^{\tau - 1} {{{\bar y}_{t - i}}} }.
#' @author JJB
#' @references Long-Memory Processes, the Allan Variance and Wavelets, D. B. Percival and P. Guttorp
#' @examples
#' set.seed(999)
#' # Simulate white noise (P 1) with sigma^2 = 4
#' N = 100000
#' white.noise = rnorm(N, 0, 2)
#' #plot(white.noise,ylab="Simulated white noise process",xlab="Time",type="o")
#' #Simulate random walk (P 4)
#' random.walk = cumsum(0.1*rnorm(N, 0, 2))
#' combined.ts = white.noise+random.walk
#' av_mat = avar_mo_cpp(combined.ts)
#' @keywords internal
avar_mo_cpp <- function(x) {
.Call('gmwm_avar_mo_cpp', PACKAGE = 'gmwm', x)
}
#' @title ARMA Adapter to ARMA to WV Process function
#' @description Molds the data so that it works with the arma_to_wv function.
#' @param theta A \code{vec} that contains all the parameter estimates.
#' @param p A \code{int} that indicates the number of AR coefficients
#' @param q A \code{int} that indicates the number of MA coefficients.
#' @param tau A \code{vec} that lists the scales of the process e.g. 2^(1:J)
#' @return A \code{vec} containing the ARMA to WV results
#' @keywords internal
#' @backref src/analytical_matrix_derivatives.cpp
#' @backref src/analytical_matrix_derivatives.h
arma_adapter <- function(theta, p, q, tau) {
.Call('gmwm_arma_adapter', PACKAGE = 'gmwm', theta, p, q, tau)
}
#' @title Calculates the Jacobian for the ARMA process
#' @description Figure out the Jacobian for an ARMA process.
#' @param theta A \code{vec} that contains all the parameter estimates.
#' @param p A \code{int} that indicates the number of AR coefficients
#' @param q A \code{int} that indicates the number of MA coefficients.
#' @param tau A \code{vec} that lists the scales of the process e.g. 2^(1:J)
#' @return A \code{mat} that returns the numerical jacobian of the ARMA process.
#' @keywords internal
#' @backref src/analytical_matrix_derivatives.cpp
#' @backref src/analytical_matrix_derivatives.h
jacobian_arma <- function(theta, p, q, tau) {
.Call('gmwm_jacobian_arma', PACKAGE = 'gmwm', theta, p, q, tau)
}
#' Analytic D matrix for AR(1) process
#' @param phi A \code{double} corresponding to the phi coefficient of an AR(1) process.
#' @param sig2 A \code{double} corresponding to the error term of an AR(1) process.
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{\phi ^2}{sigma^2} and the second column contains the partial derivative with respect to \eqn{\sigma ^2}{sigma^2}
#' @details
#' The haar wavelet variance is given as \eqn{\frac{{\left( {\frac{\tau }{2} - 3{\rho _0} - \frac{{\tau \rho _0^2}}{2} + 4\rho _0^{\frac{\tau }{2} + 1} - \rho _0^{\tau + 1}} \right)\nu _0^2}}{{\frac{{{\tau ^2}}}{8}{{\left( {1 - {\rho _0}} \right)}^2}\left( {1 - \rho _0^2} \right)}}}{See PDF Manual for equation}
#' Note: \eqn{\phi = \rho}{phi = rho} and \eqn{V _0^2 = \sigma _0^2}{V[0]^2 = sigma[0]^2}.
#' Due to length, the analytical derivations of the AR(1) haar wavelet variance are given in a supplied file within vignette.
#' @author JJB
#' @examples
#' deriv_ar1(.3, 1, 2^(1:5))
deriv_ar1 <- function(phi, sig2, tau) {
.Call('gmwm_deriv_ar1', PACKAGE = 'gmwm', phi, sig2, tau)
}
#' Analytic second derivative matrix for AR(1) process
#' @param phi A \code{double} corresponding to the phi coefficient of an AR(1) process.
#' @param sig2 A \code{double} corresponding to the error term of an AR(1) process.
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the second partial derivative with respect to \eqn{\phi ^2}{sigma^2} and the second column contains the second partial derivative with respect to \eqn{\sigma ^2}{sigma^2}
#' @details
#' The haar wavelet variance is given as \eqn{\frac{{\left( {\frac{\tau }{2} - 3{\rho _0} - \frac{{\tau \rho _0^2}}{2} + 4\rho _0^{\frac{\tau }{2} + 1} - \rho _0^{\tau + 1}} \right)\nu _0^2}}{{\frac{{{\tau ^2}}}{8}{{\left( {1 - {\rho _0}} \right)}^2}\left( {1 - \rho _0^2} \right)}}}{See PDF Manual for equation}
#' Note: \eqn{\phi = \rho}{phi = rho} and \eqn{V _0^2 = \sigma _0^2}{V[0]^2 = sigma[0]^2}.
#' Due to length, the analytical derivations of the AR(1) haar wavelet variance are given in a supplied file within vignette.
#' @author JJB
#' @examples
#' deriv_2nd_ar1(.3, 1, 2^(1:5))
deriv_2nd_ar1 <- function(phi, sig2, tau) {
.Call('gmwm_deriv_2nd_ar1', PACKAGE = 'gmwm', phi, sig2, tau)
}
#' Analytic D matrix for drift process
#' @param omega A \code{double} that is the slope of the drift.
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{\omega _0}{omega[0]}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{{\tau ^2}\omega _0^2}}{2}}{nu^2(tau) = tau^2 omega_0^2 / 2}.
#' Taking the derivative with respect to \eqn{\omega _0^2}{omega_0^2} yields: \eqn{\frac{\partial }{{\partial {\omega _0}}}{\nu ^2}\left( \tau \right) = {\tau ^2}{\omega _0}}{tau^2 * omega_0}
#' @author JJB
#' @examples
#' deriv_dr(5.3, 2^(1:5))
deriv_dr <- function(omega, tau) {
.Call('gmwm_deriv_dr', PACKAGE = 'gmwm', omega, tau)
}
#' Analytic second derivative matrix for drift process
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the second partial derivative with respect to \eqn{\omega _0}{omega[0]}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{{\tau ^2}\omega _0^2}}{2}}{nu^2(tau) = tau^2 omega_0^2 / 2}.
#' Taking the derivative with respect to \eqn{\omega _0^2}{omega_0^2} yields: \eqn{\frac{\partial }{{\partial {\omega _0}}}{\nu ^2}\left( \tau \right) = {\tau ^2}{\omega _0}}{tau^2 * omega_0}
#' Taking second derivative with respect to \eqn{\omega _0^2}{omega_0^2} yields: \eqn{\frac{{{\partial ^2}}}{{\partial \omega _0^2}}{\nu ^2}\left( \tau \right) = {\tau ^2}}{tau^2}
#' @author JJB
#' @examples
#' deriv_2nd_dr(2^(1:5))
deriv_2nd_dr <- function(tau) {
.Call('gmwm_deriv_2nd_dr', PACKAGE = 'gmwm', tau)
}
#' Analytic D matrix quantisation noise process
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{Q _0^2}{Q[0]^2}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{3Q_0^2}}{{2{\tau ^2}}}}{nu^2(tau) = 3*Q[0]^2 / 2*tau^2}.
#' Taking the derivative with respect to \eqn{Q _0^2}{Q[0]^2} yields: \deqn{\frac{\partial }{{\partial Q_0^2}}{\nu ^2}\left( \tau \right) = \frac{3}{{2{\tau ^2}}}}{3/(2*tau^2)}.
#' The second derivative derivative with respect to \eqn{Q _0^2}{Q[0]^2} is then: \deqn{\frac{{{\partial ^2}}}{{\partial Q_0^4}}{\nu ^2}\left( \tau \right) = 0}{0}.
#' @author JJB
#' @examples
#' deriv_qn(2^(1:5))
deriv_qn <- function(tau) {
.Call('gmwm_deriv_qn', PACKAGE = 'gmwm', tau)
}
#' Analytic D matrix random walk process
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{\gamma _0^2}{gamma[0]^2}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{\left( {2{\tau ^2} + 1} \right)\gamma _0^2}}{{24\tau }}}{nu^2(tau) = (2*tau^2+1)*gamma^2 / (24*tau)}.
#' Taking the first derivative with respect to \eqn{\gamma _0^2}{gamma_0^2} yields: \deqn{\frac{{{\partial ^2}}}{{\partial \gamma _0^4}}{\nu ^2}\left( \tau \right) = 0}{(2*tau^2+1) / (24*tau)}
#' The second derivative derivative with respect to \eqn{\gamma _0^2}{gamma[0]^2} is then: \deqn{\frac{{{\partial ^2}}}{{\partial \sigma_0^4}}{\nu ^2}\left( \tau \right) = 0}{0}.
#' @author JJB
#' @examples
#' deriv_rw(2^(1:5))
deriv_rw <- function(tau) {
.Call('gmwm_deriv_rw', PACKAGE = 'gmwm', tau)
}
#' Analytic D matrix white noise process
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{\sigma _0^2}{sigma[0]^2}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{\sigma _0^2}}{\tau }}{nu^2(tau) = sigma_0^2 / tau}.
#' Taking the derivative with respect to \eqn{\sigma _0^2}{sigma_0^2} yields: \eqn{\frac{\partial }{{\partial \sigma _0^2}}{\nu ^2}\left( \tau \right) = \frac{1}{\tau }}{1/tau}
#' @author JJB
#' @examples
#' deriv_wn(2^(1:5))
deriv_wn <- function(tau) {
.Call('gmwm_deriv_wn', PACKAGE = 'gmwm', tau)
}
#' Analytic D matrix of Processes
#' @description This function computes each process to WV (haar) in a given model.
#' @param theta A \code{vec} containing the list of estimated parameters.
#' @param desc A \code{vector<string>} containing a list of descriptors.
#' @param objdesc A \code{field<vec>} containing a list of object descriptors.
#' @param tau A \code{vec} containing the scales e.g. 2^(1:J)
#' @return A \code{matrix} with the process derivatives going down the column
#' @details
#' Function returns the matrix effectively known as "D"
#' @author JJB
#' @examples
#' mod = AR1(.4,1) + WN(.2) + DR(.005)
#' derivative_first_matrix(mod$theta, mod$desc, mod$obj.desc, 2^(1:9))
derivative_first_matrix <- function(theta, desc, objdesc, tau) {
.Call('gmwm_derivative_first_matrix', PACKAGE = 'gmwm', theta, desc, objdesc, tau)
}
#' Analytic D matrix of Processes
#' @description This function computes each process to WV (haar) in a given model.
#' @param theta A \code{vec} containing the list of estimated parameters.
#' @param desc A \code{vector<string>} containing a list of descriptors.
#' @param objdesc A \code{field<vec>} containing a list of object descriptors.
#' @param tau A \code{vec} containing the scales e.g. 2^(1:J)
#' @param omegadiff A \code{vec} that contains the result of Omega * (wv_empir - wv_theo)
#' @return A \code{matrix} with the process derivatives going down the column
#' @details
#' Function returns the matrix effectively known as "D"
#' @author JJB
#' @examples
#' #TBA
#' @keywords internal
D_matrix <- function(theta, desc, objdesc, tau, omegadiff) {
.Call('gmwm_D_matrix', PACKAGE = 'gmwm', theta, desc, objdesc, tau, omegadiff)
}
#' @title Hook into R's ARIMA function
#' @description Uses R's ARIMA function to obtain CSS values for starting condition
#' @param data A \code{vec} of data.
#' @param params A \code{vec} of the ARMA parameters
#' @return A \code{vec} containing the CSS of the ARMA parameters.
#' @keywords internal
Rcpp_ARIMA <- function(data, params) {
.Call('gmwm_Rcpp_ARIMA', PACKAGE = 'gmwm', data, params)
}
#' Sort Matrix by Column
#'
#' Sorts a given matrix by a specific column while retain the elements in each row.
#'
#' @param x A \code{matrix} to sort
#' @param col A \code{int} that indicates the column the matrix should sort by.
#' @details
#' The functional difference between armadillo's sort() and sort_mat() is straight forward.
#' sort() will sort each column without respect to the rows.
#' Using sort_matrix will sort only 1 column and retain the other elements to be in the same row.
#' @return The matrix sorted by values in the specified column.
#' @keywords internal
sort_mat <- function(x, col) {
.Call('gmwm_sort_mat', PACKAGE = 'gmwm', x, col)
}
#' @title Reverse Subset Column
#' @description
#' Subsets the column by going from high indices to low (the reverse of the supported practice)
#' @param x A \code{matrix} of dimensions M x N
#' @param start A \code{unsigned int} that indicates the starting column.
#' @param end A \code{unsigned int} that indicates the ending column.
#' @return x A \code{matrix} with matrix rows displayed in reverse order
#' @details Consider a vector x=[[1,2],[3,4]].
#' By setting \code{start=1} and \code{end=0}, the function would output x=[[2,1],[4,1]].
#' Start and end must be valid C++ matrix locations. (e.g. matrix cols start at 0 and not 1)
#' @author JJB
#' @examples
#' x = matrix(c(1,2,3,4), nrow = 2,byrow = TRUE)
#' rev_col_subset(x, 1, 0)
#' @keywords internal
rev_col_subset <- function(x, start, end) {
.Call('gmwm_rev_col_subset', PACKAGE = 'gmwm', x, start, end)
}
#' @title Reverse Subset Row
#' @description Subsets the row by going from high indices to low (the reverse of the supported practice)
#' @param x A \code{matrix} of dimensions M x N
#' @param start A \code{unsigned int} that indicates the starting row.
#' @param end A \code{unsigned int} that indicates the ending row.
#' @return x A \code{matrix} with matrix rows displayed in reversed order
#' @details Consider a vector x=[[1,2],[3,4]], the function would output x=[[3,4],[1,2]].
#' Start and end must be valid C++ matrix locations. (e.g. matrix rows start at 0 and not 1)
#' @author JJB
#' @examples
#' x = matrix(c(1,2,3,4), nrow=2,byrow=TRUE)
#' rev_row_subset(x, 1, 0)
#' @keywords internal
rev_row_subset <- function(x, start, end) {
.Call('gmwm_rev_row_subset', PACKAGE = 'gmwm', x, start, end)
}
#' @title Reverse Armadillo Vector
#' @description Reverses the order of an Armadillo Vector
#' @usage reverse_vec(x)
#' @param x A \code{column vector} of length N
#' @return x A \code{column vector} with its contents reversed.
#' @details Consider a vector x=[1,2,3,4,5], the function would output x=[5,4,3,2,1].
#' @author JJB
#' @examples
#' x = 1:5
#' reverse_vec(x)
#' @keywords internal
reverse_vec <- function(x) {
.Call('gmwm_reverse_vec', PACKAGE = 'gmwm', x)
}
#' @title Transform an Armadillo field<vec> to a matrix
#' @description Unlists vectors in a field and places them into a matrix
#' @param x A \code{field<vec>}.
#' @return A \code{mat} containing the field elements within a column.
#' @author JJB
#' @examples
#' x=rnorm(100)
#' @keywords internal
field_to_matrix <- function(x) {
.Call('gmwm_field_to_matrix', PACKAGE = 'gmwm', x)
}
#' @title Accumulation of Armadillo field<vec>
#' @description Sums vectors in a field into a single variable.
#' @param x A \code{field<vec>}.
#' @return An \code{mat} containing the field elements within a column.
#' @author JJB
#' @examples
#' x=rnorm(100)
#' @keywords internal
sum_field_vec <- function(x) {
.Call('gmwm_sum_field_vec', PACKAGE = 'gmwm', x)
}
#' @title Build List of Unique Models
#' @description Creates a set containing unique strings.
#' @param combs A \code{mat} that is a binary matrix (0,1) containing the combinations of different variables.
#' @param x A \code{vec<string>} that contains a list of model descriptors.
#' @return A \code{set<string>} that contains the list of unique models.
#' @keywords internal
build_model_set <- function(combs, x) {
.Call('gmwm_build_model_set', PACKAGE = 'gmwm', combs, x)
}
#' @title Conversion function of Vector to Set
#' @description Converts a vector into a set
#' @param x A \code{vec<vec<string>>} that contains a list of model descriptors.
#' @return A \code{set<vector<string>>} that contains the list of unique models.
#' @keywords internal
vector_to_set <- function(model_str) {
.Call('gmwm_vector_to_set', PACKAGE = 'gmwm', model_str)
}
#' @title Find the Common Denominator of the Models
#' @description Determines the common denominator among models
#' @param x A \code{vector< vector<string> >} that contains all possible models under consideration
#' @return A \code{vector<string>} that contains the terms of the common denominator of all models
#' @keywords internal
find_full_model <- function(x) {
.Call('gmwm_find_full_model', PACKAGE = 'gmwm', x)
}
#' @title Find the Rank Models result
#' @description Provides the core material to create an S3 object for rank.models
#' @param data A \code{vec} of data.
#' @param model_str A \code{vector<vector<string>>} that gives a list of models to test.
#' @param full_model A \code{vector<string>} that contains the largest / full model.
#' @param alpha A \code{double} that indicates the alpha level for CIs.
#' @param compute_v A \code{string} indicating the type of V matrix to generate
#' @param model_type A \code{string} that describes the model generation / transformation: 'ssm' or 'imu'
#' @param K A \code{int} that controls how many times the GMWM is run.
#' @param H A \code{int} that controls how many bootstraps occur.
#' @param G A \code{int} that controls how many guesses occur.
#' @param robust A \code{bool} that indicates whether to use classical or robust wavelet variance.
#' @param eff A \code{double} that indicates the efficiency to use.
#' @param bs_optimism A \code{bool} that indicates whether the model selection score should be calculated with bootstrap or asymptotics.
#' @return A \code{field<field<field<mat>>>} that contains the model score matrix and the best GMWM model object.
#' @keywords internal
rank_models <- function(data, model_str, full_model, alpha, compute_v, model_type, K, H, G, robust, eff, bs_optimism) {
.Call('gmwm_rank_models', PACKAGE = 'gmwm', data, model_str, full_model, alpha, compute_v, model_type, K, H, G, robust, eff, bs_optimism)
}
#' @title Find the auto imu result
#' @description Provides the core material to create an S3 object for auto.imu
#' @param data A \code{mat} containing multiple columns of independent data with the same number of observations.
#' @param model_str A \code{vector<vector<string>>} that gives a list of models to test.
#' @param full_model A \code{vector<string>} that contains the largest / full model.
#' @param alpha A \code{double} that indicates the alpha level for CIs.
#' @param compute_v A \code{string} indicating the type of V matrix to generate
#' @param model_type A \code{string} that describes the model generation / transformation: 'ssm' or 'imu'
#' @param K A \code{int} that controls how many times the GMWM is run.
#' @param H A \code{int} that controls how many bootstraps occur.
#' @param G A \code{int} that controls how many guesses occur.
#' @param robust A \code{bool} that indicates whether to use classical or robust wavelet variance.
#' @param eff A \code{double} that indicates the efficiency to use.
#' @param bs_optimism A \code{bool} that indicates whether the model selection score should be calculated with bootstrap or asymptotics.
#' @return A \code{field<field<field<mat>>>} that contains the model score matrix and the best GMWM model object.
#' @keywords internal
auto_imu <- function(data, combs, full_model, alpha, compute_v, model_type, K, H, G, robust, eff, bs_optimism) {
.Call('gmwm_auto_imu', PACKAGE = 'gmwm', data, combs, full_model, alpha, compute_v, model_type, K, H, G, robust, eff, bs_optimism)
}
#' @title Bootstrap for Matrix V
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
cov_bootstrapper <- function(theta, desc, objdesc, N, robust, eff, H, diagonal_matrix) {
.Call('gmwm_cov_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, N, robust, eff, H, diagonal_matrix)
}
#' @title Bootstrap for Optimism
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
optimism_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_optimism_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Bootstrap for Optimism and GoF
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @param scales A \code{vec} containing the scales of the process.
#' @param model_type A \code{string} containing the model type either: SSM or IMU
#' @param N A \code{int} indicating how long the integer is.
#' @param robust A \code{bool} indicating robust (T) or classical (F).
#' @param eff A \code{double} that handles efficiency.
#' @param H A \code{int} that indicates how many bootstraps should be obtained.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
opt_n_gof_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_opt_n_gof_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Bootstrap for Standard Deviations of Theta Estimates
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
gmwm_sd_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_gmwm_sd_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Generate the Confidence Interval for GOF Bootstrapped
#' @description yaya
#' @param obj A \code{double} containing the objective value of the solution.
#' @param obj_boot A \code{vec} containing the objective values obtained while bootstrapping under Theta_hat.
#' @param B An \code{int} indicating how many times the bootstrapper should be run.
#' @param alpha A \code{double} indicating the amount of confidence for CI
#' @return A \code{vec} that has the alpha/2.0 quantile and then the 1-alpha/2.0 quantile.
#' @keywords internal
boot_pval_gof <- function(obj, obj_boot, B = 1000L, alpha = 0.05) {
.Call('gmwm_boot_pval_gof', PACKAGE = 'gmwm', obj, obj_boot, B, alpha)
}
#' @title Bootstrap for Estimating Both Theta and Theta SD
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
gmwm_param_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_gmwm_param_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Bootstrap for Everything!
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
all_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_all_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Absolute Value or Modulus of a Complex Number Squared.
#' @description Computes the squared value of the Modulus.
#' @param x A \code{cx_vec}.
#' @return A \code{vec} containing the modulus squared for each element.
#' @details Consider a complex number defined as: \eqn{z = x + i y} with real \eqn{x} and \eqn{y},
#' The modulus is defined as: \eqn{r = Mod\left(z\right) = \sqrt{\left(x^2 + y^2\right)}}{r = Mod(z) = sqrt(x^2 + y^2)}
#' This function will return: \eqn{r^2 = Mod\left(z\right)^2 = x^2 + y^2}
#' @examples
#' Mod_squared_cpp(c(1+.5i, 2+1i, 5+9i))
#' @keywords internal
Mod_squared_cpp <- function(x) {
.Call('gmwm_Mod_squared_cpp', PACKAGE = 'gmwm', x)
}
#' @title Absolute Value or Modulus of a Complex Number.
#' @description Computes the value of the Modulus.
#' @param x A \code{cx_vec}.
#' @return A \code{vec} containing the modulus for each element.
#' @details Consider a complex number defined as: \eqn{z = x + i y} with real \eqn{x} and \eqn{y},
#' The modulus is defined as: \eqn{r = Mod(z) = \sqrt{(x^2 + y^2)}}
#' @examples
#' Mod_cpp(c(1+.5i, 2+1i, 5+9i))
#' @keywords internal
Mod_cpp <- function(x) {
.Call('gmwm_Mod_cpp', PACKAGE = 'gmwm', x)
}
#' @title Computes the (MODWT) wavelet covariance matrix
#' @description Calculates the (MODWT) wavelet covariance matrix
#' @param signal_modwt A \code{field<vec>} that contains the modwt decomposition.
#' @param nb_level A \code{integer} that contains the level of decomposition J.
#' @param compute_v A \code{string} that indicates what kind of matrix should be created. Possible options: "diag" or "none"
#' @param robust A \code{boolean} that triggers the use of the robust estimate.
#' @param eff A \code{double} that indicates the efficiency as it relates to an MLE.
#' @return A \code{field<mat>} containing the covariance matrix.
#' @keywords internal
#' @examples
#' \dontrun{
#' x=rnorm(100)
#' decomp = modwt(x)
#' V = compute_cov_cpp(decomp$data, decomp$nlevels, compute_v="diag", robust = TRUE, eff=0.6)
#' }
compute_cov_cpp <- function(signal_modwt, nb_level, compute_v = "diag", robust = TRUE, eff = 0.6) {
.Call('gmwm_compute_cov_cpp', PACKAGE = 'gmwm', signal_modwt, nb_level, compute_v, robust, eff)
}
#' @title Computes the (MODWT) wavelet covariance matrix using Chi-square confidence interval bounds
#' @description Calculates the (MODWT) wavelet covariance matrix using Chi-square confidence interval bounds
#' @param ci_hi A \code{vec} that contains the upper confidence interval points.
#' @param ci_lo A \code{vec} that contains the lower confidence interval points.
#' @return A diagonal matrix.
#' @keywords internal
#' @examples
#' \dontrun{
#' x=runif(100)
#' y=x+3
#' fast_cov_cpp(y,x)
#' }
fast_cov_cpp <- function(ci_hi, ci_lo) {
.Call('gmwm_fast_cov_cpp', PACKAGE = 'gmwm', ci_hi, ci_lo)
}
#' @title Discrete Wavelet Transform
#' @description Calculation of the coefficients for the discrete wavelet transformation.
#' @param x A \code{vector} with dimensions \eqn{N\times 1}{N x 1}.
#' @param filter_name A \code{string} indicating the filter.
#' @param nlevels An \code{integer}, \eqn{J}, indicating the level of the decomposition.
#' @param boundary A \code{string} indicating the type of boundary method to use. Either \code{boundary="periodic"} or \code{"reflection"}.
#' @param brickwall A \code{bool} indicating whether the a brick wall procedure should be applied to the coefficients.
#' @return y A \code{field<vec>} that contains the wavelet coefficients for each decomposition level
#' @details
#' Performs a level J decomposition of the time series using the pyramid algorithm
#' @author JJB
#' @keywords internal
#' @examples
#' set.seed(999)
#' x = rnorm(2^8)
#' dwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
dwt_cpp <- function(x, filter_name, nlevels, boundary, brickwall) {
.Call('gmwm_dwt_cpp', PACKAGE = 'gmwm', x, filter_name, nlevels, boundary, brickwall)
}
#' @title Maximum Overlap Discrete Wavelet Transform
#' @description
#' Calculation of the coefficients for the discrete wavelet transformation
#' @inheritParams dwt_cpp
#' @return y A \code{field<vec>} that contains the wavelet coefficients for each decomposition level
#' @keywords internal
#' @details
#' Performs a level J decomposition of the time series using the pyramid algorithm.
#' Use this implementation to supply custom parameters instead of modwt(x),
#' which serves as a wrapper function.
#' @author JJB
#' @keywords internal
#' @examples
#' set.seed(999)
#' x = rnorm(100)
#' modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
modwt_cpp <- function(x, filter_name, nlevels, boundary, brickwall) {
.Call('gmwm_modwt_cpp', PACKAGE = 'gmwm', x, filter_name, nlevels, boundary, brickwall)
}
#' @title Removal of Boundary Wavelet Coefficients
#' @description Removes the first n wavelet coefficients.
#' @param x A \code{field<vec>} that contains the nlevel decomposition using either modwt or dwt.
#' @param wave_filter A \code{field<vec>} containing filter information. Only "haar" is implemented.
#' @param method A \code{string} to describe the mode. Choose between "modwt" and "dwt"
#' @return A \code{field<vec>} with boundary modwt or dwt taken care of.
#' @keywords internal
#' @details
#' The vectors are truncated by removing the first n wavelet coefficients.
#' These vectors are then stored into the field that is returned.
#' Note: As a result, there are no NA's introduced and hence the na.omit is not needed.
#' @examples
#' x = rnorm(100)
#' me = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = FALSE)
#' brick_wall(me, select_filter("haar"), "modwt")
brick_wall <- function(x, wave_filter, method) {
.Call('gmwm_brick_wall', PACKAGE = 'gmwm', x, wave_filter, method)
}
#' @title Generate a white noise process
#' @description Generates a white noise process with variance parameter sigma.
#' @param N An \code{integer} for signal length.
#' @param sigma2 A \code{double} that contains process variance.
#' @return wn A \code{vec} containing the white noise.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_wn(10, 1.5)
gen_wn <- function(N, sigma2 = 1) {
.Call('gmwm_gen_wn', PACKAGE = 'gmwm', N, sigma2)
}
#' @title Generate a drift
#' @description Generates a drift sequence with a given slope.
#' @param N An \code{integer} for signal length.
#' @param slope A \code{double} that contains drift slope
#' @return gd A \code{vec} containing the drift.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_dr(10, 8.2)
gen_dr <- function(N, slope = 5) {
.Call('gmwm_gen_dr', PACKAGE = 'gmwm', N, slope)
}
#' @title Generate a Quantisation Noise (QN) sequence
#' @description Generate an QN sequence given q2
#' @param N An \code{integer} for signal length.
#' @param q2 A \code{double} that contains autocorrection.
#' @return A \code{vec} containing the QN process.
#' @keywords internal
#' @details
#' To generate the quantisation noise, we follow this recipe:
#' First, we generate using a random uniform distribution:
#' \deqn{U_k^*\sim U\left[ {0,1} \right]}{U_k^*~U[0,1]}
#'
#' Then, we multiple the sequence by \eqn{\sqrt{12}}{sqrt(12)} so:
#' \deqn{{U_k} = \sqrt{12} U_k^*}{U_k = sqrt(12)*U_k^*}
#'
#' Next, we find the derivative of \eqn{{U_k}}{U_k}
#' \deqn{{{\dot U}_k} = \frac{{{U_{k + \Delta t}} - {U_k}}}{{\Delta t}}}{U_k^. = (U_(k + (delta)t) - U_k)}
#'
#' In this case, we modify the derivative such that:
#' \eqn{{{\dot U}_k}\Delta t = {U_{k + \Delta t}} - {U_k}}{U_k^. * (delta)t = U_{k + (delta)*t} - U_k}
#'
#' Thus, we end up with:
#' \deqn{{x_k} = \sqrt Q {{\dot U}_k}\Delta t}{x_k = sqrt(Q)*U_k^.*(delta)t}
#' \deqn{{x_k} = \sqrt Q \left( {{U_{k + 1}} - {U_k}} \right)}{x_k = sqrt(Q)* (U_(k+1) - U_(k))}
#'
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @examples
#' gen_qn(10, 5)
gen_qn <- function(N, q2 = .1) {
.Call('gmwm_gen_qn', PACKAGE = 'gmwm', N, q2)
}
#' @title Generate an AR(1) sequence
#' @description Generate an AR sequence given phi and sig2.
#' @details This needs to be extended to AR(p) see \code{arima.sim} and \code{filter}.
#' @param N An \code{integer} for signal length.
#' @param phi A \code{double} that contains autocorrection.
#' @param sigma2 A \code{double} that contains process variance.
#' @return gm A \code{vec} containing the AR(1) process.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_ar1(10, 5, 1.2)
gen_ar1 <- function(N, phi = .3, sigma2 = 1) {
.Call('gmwm_gen_ar1', PACKAGE = 'gmwm', N, phi, sigma2)
}
#' @title Generate a random walk without drift
#' @description Generates a random walk without drift.
#' @param N An \code{integer} for signal length.
#' @param sigma2 A \code{double} that contains process variance.
#' @return grw A \code{vec} containing the random walk without drift.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_rw(10, 8.2)
gen_rw <- function(N, sigma2 = 1) {
.Call('gmwm_gen_rw', PACKAGE = 'gmwm', N, sigma2)
}
#' @title Generate ARMA
#' @description Generate observations for a supplied ARMA model.
#' @param N An \code{integer} for signal length.
#' @param ar A \code{vec} that contains the AR coefficients.
#' @param ma A \code{vec} that contains the MA coefficients.
#' @param sigma2 A \code{double} that contains process variance.
#' @param n_start An \code{unsigned int} that indicates the amount of observations to be used for the burn in period.
#' @details The innovations are generated from a normal distribution.
#' @return A \code{vec} that contains the generated observations.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_arma(100, c(.3,.5), c(.1), 1, 0)
gen_arma <- function(N, ar, ma, sigma2 = 1.5, n_start = 0L) {
.Call('gmwm_gen_arma', PACKAGE = 'gmwm', N, ar, ma, sigma2, n_start)
}
#' @title Generate Time Series based on Model (Internal)
#' @description Create a time series based on a supplied time series model.
#' @param N An \code{interger} containing the amount of observations for the time series.
#' @param theta A \code{vec} containing the parameters to use to generate the model
#' @param desc A \code{vector<string>} containing the different model types (AR1, WN, etc..)
#' @param objdesc A \code{field<vec>} contains the different model objects e.g. AR1 = c(1,1)
#' @return A \code{vec} that contains combined time series.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' # AR
#' set.seed(1336)
#' gen_model(1000, c(.9,1), "AR1", list(c(1,1)))
gen_model <- function(N, theta, desc, objdesc) {
.Call('gmwm_gen_model', PACKAGE = 'gmwm', N, theta, desc, objdesc)
}
#' @title Generate Latent Time Series based on Model (Internal)
#' @description Create a latent time series based on a supplied time series model.
#' @param N An \code{interger} containing the amount of observations for the time series.
#' @param theta A \code{vec} containing the parameters to use to generate the model.
#' @param desc A \code{vector<string>} containing the different model types (AR1, WN, etc..).
#' @param objdesc A \code{field<vec>} containing the different model objects e.g. AR1 = c(1,1)
#' @return A \code{mat} containing data for each decomposed and combined time series.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' # AR
#' set.seed(1336)
#' gen_lts(1000, c(.9,1), "AR1", list(c(1,1)))
gen_lts <- function(N, theta, desc, objdesc) {
.Call('gmwm_gen_lts', PACKAGE = 'gmwm', N, theta, desc, objdesc)
}
#' @title Optim loses NaN
#' @description This function takes numbers that are very small and sets them to the minimal tolerance for C++.
#' Doing this prevents NaN from entering the optim routine.
#' @param theta A \code{vec} that contains estimated GMWM theta values (untransformed).
#' @return A \code{vec} that contains safe theta values.
#' @keywords internal
code_zero <- function(theta) {
.Call('gmwm_code_zero', PACKAGE = 'gmwm', theta)
}
#' @title Engine for obtaining the GMWM Estimator
#' @description This function uses the Generalized Method of Wavelet Moments (GMWM) to estimate the parameters of a time series model.
#' @param theta A \code{vec} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} containing a list of parameters (e.g. AR(1) = c(1,1), ARMA(p,q) = c(p,q,1))
#' @param model_type A \code{string} that represents the model transformation
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance
#' @param omega A \code{mat} that represents the covariance matrix.
#' @param scales A \code{vec} that contains the scales or taus (2^(1:J))
#' @param starting A \code{bool} that indicates whether we guessed starting (T) or the user supplied estimates (F).
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' If type = "imu" or "ssm", then parameter vector should indicate the characters of the models that compose the latent or state-space model.
#' The model options are:
#' \itemize{
#' \item{"AR1"}{a first order autoregressive process with parameters \eqn{(\phi,\sigma^2)}{phi, sigma^2}}
#' \item{"ARMA"}{an autoregressive moving average process with parameters \eqn{(\phi _p, \theta _q, \sigma^2)}{phi[p], theta[q], sigma^2}}
#' \item{"DR"}{a drift with parameter \eqn{\omega}{omega}}
#' \item{"QN"}{a quantization noise process with parameter \eqn{Q}}
#' \item{"RW"}{a random walk process with parameter \eqn{\sigma^2}{sigma^2}}
#' \item{"WN"}{a white noise process with parameter \eqn{\sigma^2}{sigma^2}}
#' }
#' If model_type = "imu" or type = "ssm" then
#' starting values pass through an initial bootstrap and pseudo-optimization before being passed to the GMWM optimization.
#' If robust = TRUE the function takes the robust estimate of the wavelet variance to be used in the GMWM estimation procedure.
#'
#' @author JJB
#' @references Wavelet variance based estimation for composite stochastic processes, S. Guerrier and Robust Inference for Time Series Models: a Wavelet-Based Framework, S. Guerrier
#' @keywords internal
#' @backref src/gmwm_logic.cpp
#' @backref src/gmwm_logic.h
gmwm_engine <- function(theta, desc, objdesc, model_type, wv_empir, omega, scales, starting) {
.Call('gmwm_gmwm_engine', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, wv_empir, omega, scales, starting)
}
#' @title Update Wrapper for the GMWM Estimator
#' @description This function uses information obtained previously (e.g. WV covariance matrix) to re-estimate a different model parameterization
#' @param theta A \code{vec} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} containing a list of parameters (e.g. AR(1) = c(1,1), ARMA(p,q) = c(p,q,1))
#' @param model_type A \code{string} that represents the model transformation
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance
#' @param omega A \code{mat} that represents the covariance matrix.
#' @param scales A \code{vec} that contains the scales or taus (2^(1:J))
#' @param starting A \code{bool} that indicates whether we guessed starting (T) or the user supplied estimates (F).
#' @return A \code{field<mat>} that contains the parameter estimates from GMWM estimator.
#' @author JJB
#' @references Wavelet variance based estimation for composite stochastic processes, S. Guerrier and Robust Inference for Time Series Models: a Wavelet-Based Framework, S. Guerrier
#' @keywords internal
#' @backref src/gmwm_logic.cpp
#' @backref src/gmwm_logic.h
gmwm_update_cpp <- function(theta, desc, objdesc, model_type, N, expect_diff, ranged, orgV, scales, wv, starting, compute_v, K, H, G, robust, eff) {
.Call('gmwm_gmwm_update_cpp', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, N, expect_diff, ranged, orgV, scales, wv, starting, compute_v, K, H, G, robust, eff)
}
#' @title Master Wrapper for the GMWM Estimator
#' @description This function generates WV, GMWM Estimator, and an initial test estimate.
#' @param data A \code{vec} containing the data.
#' @param theta A \code{vec} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} containing a list of parameters (e.g. AR(1) = c(1,1), ARMA(p,q) = c(p,q,1))
#' @param model_type A \code{string} that represents the model transformation
#' @param starting A \code{bool} that indicates whether the supplied values are guessed (T) or are user-based (F).
#' @param alpha A \code{double} that handles the alpha level of the confidence interval (1-alpha)*100
#' @param compute_v A \code{string} that describes what kind of covariance matrix should be computed.
#' @param K An \code{int} that controls how many times theta is updated.
#' @param H An \code{int} that controls how many bootstrap replications are done.
#' @param G An \code{int} that controls how many guesses at different parameters are made.
#' @param robust A \code{bool} that indicates whether the estimation should be robust or not.
#' @param eff A \code{double} that specifies the amount of efficiency required by the robust estimator.
#' @return A \code{field<mat>} that contains a list of ever-changing estimates...
#' @author JJB
#' @references Wavelet variance based estimation for composite stochastic processes, S. Guerrier and Robust Inference for Time Series Models: a Wavelet-Based Framework, S. Guerrier
#' @keywords internal
#' @backref src/gmwm_logic.cpp
#' @backref src/gmwm_logic.h
gmwm_master_cpp <- function(data, theta, desc, objdesc, model_type, starting, alpha, compute_v, K, H, G, robust, eff) {
.Call('gmwm_gmwm_master_cpp', PACKAGE = 'gmwm', data, theta, desc, objdesc, model_type, starting, alpha, compute_v, K, H, G, robust, eff)
}
#' @title Randomly guess a starting parameter
#' @description Sets starting parameters for each of the given parameters.
#' @param desc A \code{vector<string>} that contains the model's components.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @param model_type A \code{string} that indicates whether it is an SSM or sensor.
#' @param num_param An \code{unsigned int} number of parameters in the model (e.g. # of thetas).
#' @param expect_diff A \code{double} that contains the mean of the first difference of the data
#' @param N A \code{integer} that contains the number of observations in the data.
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance.
#' @param tau A \code{vec} that contains the scales. (e.g. 2^(1:J))
#' @param double A \code{double} that contains the drift slope given by \eqn{\frac{max-min}{N}}{(Max-Min)/N}
#' @param G A \code{integer} that indicates how many random draws that should be performed.
#' @return A \code{vec} containing smart parameter starting guesses to be iterated over.
#' @keywords internal
#' @examples
#' #TBA
guess_initial <- function(desc, objdesc, model_type, num_param, expect_diff, N, wv, tau, ranged, G) {
.Call('gmwm_guess_initial', PACKAGE = 'gmwm', desc, objdesc, model_type, num_param, expect_diff, N, wv, tau, ranged, G)
}
#' @title Randomly guess starting parameters for AR1
#' @description Sets starting parameters for each of the given parameters.
#' @param draw_id An \code{unsigned int} that contains the draw principles.
#' @param last_phi A \code{double} containing the last guessed phi value.
#' @param sigma2_total A \code{double} that contains the sum of all WVs.
#' @param model_type A \code{string} that describes the model transformation.
#' @return A \code{vec} containing smart parameter starting guesses to be iterated over.
#' @keywords internal
#' @examples
#' #TBA
ar1_draw <- function(draw_id, last_phi, sigma2_total, model_type) {
.Call('gmwm_ar1_draw', PACKAGE = 'gmwm', draw_id, last_phi, sigma2_total, model_type)
}
#' @title Randomly guess starting parameters for ARMA
#' @description Sets starting parameters for each of the given parameters.
#' @param p An \code{unsigned int} that contains the amount of AR parameters to generate.
#' @param q An \code{unsigned int} that contains the amount of MA parameters to generate.
#' @param sigma2_total A \code{double} that contains the sum of all WVs.
#' @return A \code{vec} containing smart parameter starting guesses to be iterated over.
#' @keywords internal
#' @examples
#' #TBA
arma_draws <- function(p, q, sigma2_total) {
.Call('gmwm_arma_draws', PACKAGE = 'gmwm', p, q, sigma2_total)
}
#' @title Randomly guess a starting parameter
#' @description Sets starting parameters for each of the given parameters.
#' @param desc A \code{vector<string>} that contains the model's components.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @param model_type A \code{string} that indicates whether it is an SSM or sensor.
#' @param num_param An \code{unsigned int} number of parameters in the model (e.g. # of thetas).
#' @param expect_diff A \code{double} that contains the mean of the first difference of the data
#' @param N A \code{integer} that contains the number of observations in the data.
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance.
#' @param tau A \code{vec} that contains the scales. (e.g. 2^(1:J))
#' @param B A \code{integer} that indicates how many random draws that should be performed.
#' @return A \code{vec} containing smart parameter starting guesses to be iterated over.
#' @keywords internal
#' @examples
#' #TBA
guess_initial_old <- function(desc, objdesc, model_type, num_param, expect_diff, N, wv_empir, tau, B) {
.Call('gmwm_guess_initial_old', PACKAGE = 'gmwm', desc, objdesc, model_type, num_param, expect_diff, N, wv_empir, tau, B)
}
#' @title Indirect Inference for ARMA
#' @description Option for indirect inference
#' @param ar A \code{vec} that contains the coefficients of the AR process.
#' @param ma A \code{vec} that contains the coefficients of the MA process.
#' @param sigma2 A \code{double} that indicates the sigma2 parameter of the ARMA process.
#' @param N A \code{int} that indicates how long the time series is.
#' @param robust A \code{bool} that indicates whether the estimation should be robust or not.
#' @param eff A \code{double} that specifies the amount of efficiency required by the robust estimator.
#' @param H A \code{int} that indicates how many iterations should take place.
#' @return A \code{vec} with the indirect inference results.
#' @keywords internal
idf_arma <- function(ar, ma, sigma2, N, robust, eff, H) {
.Call('gmwm_idf_arma', PACKAGE = 'gmwm', ar, ma, sigma2, N, robust, eff, H)
}
#' @title Indirect Inference for ARMA
#' @description Option for indirect inference
#' @param ar A \code{vec} that contains the coefficients of the AR process.
#' @param ma A \code{vec} that contains the coefficients of the MA process.
#' @param sigma2 A \code{double} that indicates the sigma2 parameter of the ARMA process.
#' @param N A \code{int} that indicates how long the time series is.
#' @param robust A \code{bool} that indicates whether the estimation should be robust or not.
#' @param eff A \code{double} that specifies the amount of efficiency required by the robust estimator.
#' @param H A \code{int} that indicates how many iterations should take place.
#' @return A \code{vec} with the indirect inference results.
#' @keywords internal
idf_arma_total <- function(ar, ma, sigma2, N, robust, eff, H) {
.Call('gmwm_idf_arma_total', PACKAGE = 'gmwm', ar, ma, sigma2, N, robust, eff, H)
}
#' @title Calculate the Psi matrix
#' @description Computes the Psi matrix using supplied parameters
#' @param A first derivative matrix
#' @param v_hat bootstrapped V
#' @param omega original omega matrix
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @return A \code{mat} that has the first column
#' @keywords internal
calculate_psi_matrix <- function(A, v_hat, omega) {
.Call('gmwm_calculate_psi_matrix', PACKAGE = 'gmwm', A, v_hat, omega)
}
#' @title Format the Confidence Interval for Estimates
#' @description Creates hi and lo confidence based on SE and alpha.
#' @param theta A \code{vec} containing the estimates
#' @param se A \code{vec} containing the standard error
#' @param alpha A \code{double} that contains the confidence level.
#' @return A \code{mat} that has:
#' \itemize{
#' \item Column 1: Lo CI
#' \item Column 2: Hi CI
#' \item Column 3: SE
#' }
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @keywords internal
format_ci <- function(theta, se, alpha) {
.Call('gmwm_format_ci', PACKAGE = 'gmwm', theta, se, alpha)
}
#' @title Generate the Confidence Interval for Theta Estimates
#' @description Create an Asymptotic CI for the Theta Estimates.
#' @param theta A \code{vec} containing the estimates
#' @param A A \code{mat} that is the first derivative matrix.
#' @param v_hat A \code{mat} that is the bootstrapped V matrix
#' @param omega A \code{mat} that is the inverse of the diagonal V matrix.
#' @param alpha A \code{double} that contains the confidence level.
#' @return A \code{mat} that has the first column
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @keywords internal
theta_ci <- function(theta, A, v_hat, omega, alpha) {
.Call('gmwm_theta_ci', PACKAGE = 'gmwm', theta, A, v_hat, omega, alpha)
}
#' @title Compute the GOF Test
#' @description yaya
#' @template tsobj_cpp
#' @param model_type A \code{string} that contains the model type: \code{"imu"} or \code{"ssm"}
#' @param tau A \code{vec} containing the scales of a proccess.
#' @param v_hat A \code{mat} that contains the bootstrapped matrix.
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance.
#' @return A \code{vec} that has
#' \itemize{
#' \item Test Statistic
#' \item P-Value
#' \item DF
#' }
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @keywords internal
gof_test <- function(theta, desc, objdesc, model_type, tau, v_hat, wv_empir) {
.Call('gmwm_gof_test', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, tau, v_hat, wv_empir)
}
#' @title Compute the Bootstrapped GoF Test
#' @description Handles the bootstrap computation and the bootstrapped p-value.
#' @param obj_value A \code{double} that contains the optimized objective function value.
#' @param bs_obj_values A \code{vec} that contains the objective function values under bootstrap.
#' @param alpha A \code{double} that indicates the confidence.
#' @param bs_gof_p_ci A \code{bool} that indicates whether CIs should be included or not.
#' @return A \code{vec} that has
#' \itemize{
#' \item Test Statistic
#' \item Low CI
#' \item Upper CI - BS
#' }
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @keywords internal
bootstrap_gof_test <- function(obj_value, bs_obj_values, alpha, bs_gof_p_ci) {
.Call('gmwm_bootstrap_gof_test', PACKAGE = 'gmwm', obj_value, bs_obj_values, alpha, bs_gof_p_ci)
}
#' @title B Matrix
#' @description B Matrix
#' @param A A \code{mat} containing the first derivatives of the process.
#' @param at_omega A \code{mat} containing A^T * Omega
#' @return A \code{mat}
#' @backref src/model_selection.cpp
#' @backref src/model_selection.h
#' @keywords internal
B_matrix <- function(A, at_omega) {
.Call('gmwm_B_matrix', PACKAGE = 'gmwm', A, at_omega)
}
#' @title Model Score
#' @description Calculates the modeling score of a GMWM
#' @param A A \code{mat} that contains the first derivatives of the processes
#' @param At_j A \code{mat} that contains the second derivative of each process
#' @param omega A \code{mat} that contains the omega used when calculating the GMWM
#' @param v_hat A \code{mat} that contains the covariance matrix
#' @param diff A \code{vec} that is the difference of the WV empirical and WV theoretical
#' @return A \code{vec}
#' @keywords internal
#' @backref src/model_selection.cpp
#' @backref src/model_selection.h
#' @details
#' The equation is slightly different than that stated in the paper due to the bootstrap already incorporating in
#' N.
model_score <- function(A, D, omega, v_hat, obj_value) {
.Call('gmwm_model_score', PACKAGE = 'gmwm', A, D, omega, v_hat, obj_value)
}
#' @title Extract Object
#' @description Extracts the object information and returns it.
#' @param theta A \code{vec} containing the theta values.
#' @param objdesc A \code{vec} at the desc point.
#' @param cur_position An \code{integer} at the current position.
#' @return A \code{field<vec>} containing the breakdown of the object.
#' @keywords internal
obj_extract <- function(theta, objdesc, cur_position) {
.Call('gmwm_obj_extract', PACKAGE = 'gmwm', theta, objdesc, cur_position)
}
#' @title Retrieve GMWM starting value from Yannick's objective function
#' @description Obtains the GMWM starting value given by Yannick's objective function optimization
#' @template tsobj_cpp
#' @param model_type A \code{string} containing the model type. Either 'imu' or 'ssm'
#' @param wv_empir A \code{vec} containing the empirical wavelet variance.
#' @param tau A \code{vec} that contains the scales of 2^(1:J), where J is the number of scales created by the decomposition.
#' @return A \code{double} that is the value of the Objective function under Yannick's starting algorithm
#' @keywords internal
getObjFunStarting <- function(theta, desc, objdesc, model_type, wv_empir, tau) {
.Call('gmwm_getObjFunStarting', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, wv_empir, tau)
}
#' @title Retrieve GMWM starting value from Yannick's objective function
#' @description Obtains the GMWM starting value given by Yannick's objective function optimization
#' @template tsobj_cpp
#' @param model_type A \code{string} containing the model type. Either 'imu' or 'ssm'
#' @param omega A \code{mat} that is the inverse of the diagonal of the V matrix.
#' @param wv_empir A \code{vec} containing the empirical wavelet variance.
#' @param tau A \code{vec} that contains the scales of 2^(1:J), where J is the number of scales created by the decomposition.
#' @return A \code{double} that is the value of the Objective function under Yannick's starting algorithm
#' @keywords internal
getObjFun <- function(theta, desc, objdesc, model_type, omega, wv_empir, tau) {
.Call('gmwm_getObjFun', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, omega, wv_empir, tau)
}
#' @title Root Finding C++
#' @description Used to interface with Armadillo
#' @param z A \code{cx_vec} (complex vector) that has 1 in the beginning (e.g. c(1,3i,-3i))
#' @keywords internal
do_polyroot_arma <- function(z) {
.Call('gmwm_do_polyroot_arma', PACKAGE = 'gmwm', z)
}
#' @title Root Finding C++
#' @description Vroom Vroom
#' @param z A \code{vec<complex<double>} (complex vector) that has 1 in the beginning (e.g. c(1,3i,-3i))
#' @keywords internal
do_polyroot_cpp <- function(z) {
.Call('gmwm_do_polyroot_cpp', PACKAGE = 'gmwm', z)
}
#' @title ARMA process to WV
#' @description This function computes the (haar) WV of an ARMA process
#' @param ar A \code{vec} containing the coefficients of the AR process
#' @param ma A \code{vec} containing the coefficients of the MA process
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @param sigma A \code{double} containing the residual variance
#' @return A \code{vec} containing the wavelet variance of the ARMA process.
#' @examples
#' arma_to_wv(c(.23,.43), c(.34,.41,.59), 2^(1:9), 3)
#' @seealso \code{\link{ARMAtoMA_cpp}},\code{\link{ARMAacf_cpp}}
arma_to_wv <- function(ar, ma, tau, sigma) {
.Call('gmwm_arma_to_wv', PACKAGE = 'gmwm', ar, ma, tau, sigma)
}
#' @title Helper Function for ARMA to WV Approximation
#' @description Indicates where the minimum ARMAacf value is and returns that as an index.
#' @param ar A \code{vec} containing the coefficients of the AR process
#' @param ma A \code{vec} containing the coefficients of the MA process
#' @param last_tau An \code{int} the Jth scale of 2^(1:J)
#' @param alpha A \code{double} indicating the cutoff.
#' @return A \code{vec} containing the wavelet variance of the ARMA process.
#' @keywords internal
#' @seealso \code{\link{arma_to_wv_app}}
acf_sum <- function(ar, ma, last_tau, alpha = 0.99) {
.Call('gmwm_acf_sum', PACKAGE = 'gmwm', ar, ma, last_tau, alpha)
}
#' @title ARMA process to WV approximation
#' @description This function computes the (haar) WV of an ARMA process
#' @param ar A \code{vec} containing the coefficients of the AR process
#' @param ma A \code{vec} containing the coefficients of the MA process
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @param sigma A \code{double} containing the residual variance
#' @param alpha A \code{double} indicating the cutoff.
#' @return A \code{vec} containing the wavelet variance of the ARMA process.
#' @keywords internal
#' @examples
#' arma_to_wv_app(c(.23,.43), c(.34,.41,.59), 2^(1:9), 3, .9)
#' @seealso \code{\link{ARMAtoMA_cpp}},\code{\link{ARMAacf_cpp}}
arma_to_wv_app <- function(ar, ma, tau, sigma, alpha = 0.9999) {
.Call('gmwm_arma_to_wv_app', PACKAGE = 'gmwm', ar, ma, tau, sigma, alpha)
}
#' @title Quantisation Noise to WV
#' @description This function compute the WV (haar) of a Quantisation Noise (QN) process
#' @param q2 A \code{double} corresponding to variance of drift
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the QN.
#' @examples
#' x.sim = 1:1000
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = qn_to_wv(1, tau)
#' plot(tau, wv.theo, col = "red")
qn_to_wv <- function(q2, tau) {
.Call('gmwm_qn_to_wv', PACKAGE = 'gmwm', q2, tau)
}
#' @title White Noise to WV
#' @description This function compute the WV (haar) of a White Noise process
#' @param sig2 A \code{double} corresponding to variance of WN
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the white noise.
#' @examples
#' x.sim = cumsum(rnorm(100000))
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = wn_to_wv(1, tau)
#' plot(tau, wv.theo, col = "red")
wn_to_wv <- function(sig2, tau) {
.Call('gmwm_wn_to_wv', PACKAGE = 'gmwm', sig2, tau)
}
#' @title Random Walk to WV
#' @description This function compute the WV (haar) of a Random Walk process
#' @param sig2 A \code{double} corresponding to variance of RW
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the random walk.
#' @examples
#' x.sim = cumsum(rnorm(100000))
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = rw_to_wv(1,tau)
#' plot(tau, wv.theo, col = "red")
rw_to_wv <- function(sig2, tau) {
.Call('gmwm_rw_to_wv', PACKAGE = 'gmwm', sig2, tau)
}
#' @title Drift to WV
#' @description This function compute the WV (haar) of a Drift process
#' @param omega A \code{double} corresponding to variance of drift
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the drift.
#' @examples
#' x.sim = 1:1000
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = dr_to_wv(1, tau)
#' plot(tau, wv.theo, col = "red")
dr_to_wv <- function(omega, tau) {
.Call('gmwm_dr_to_wv', PACKAGE = 'gmwm', omega, tau)
}
#' @title AR1 process to WV
#' @description This function compute the WV (haar) of an AR(1) process
#' @param phi A \code{double} that is the phi term of the AR(1) process
#' @param sig2 A \code{double} corresponding to variance of AR(1) process
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the AR(1) process.
#' @examples
#' x.sim = gen_ar1( N = 10000, phi = 0.9, sigma2 = 4 )
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = ar1_to_wv(phi = 0.9, sig2 = 16, tau)
#' plot(tau, wv.theo, col = "red")
ar1_to_wv <- function(phi, sig2, tau) {
.Call('gmwm_ar1_to_wv', PACKAGE = 'gmwm', phi, sig2, tau)
}
#' @title Model Process to WV
#' @description This function computes the summation of all Processes to WV (haar) in a given model
#' @param theta A \code{vec} containing the list of estimated parameters.
#' @param desc A \code{vector<string>} containing a list of descriptors.
#' @param objdesc A \code{field<vec>} containing a list of object descriptors.
#' @param tau A \code{vec} containing the scales e.g. 2^(1:J)
#' @return A \code{vec} containing the wavelet variance of the model.
#' @examples
#' x.sim = gen_ar1( N = 10000, phi = 0.9, sigma2 = 4 )
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = ar1_to_wv(phi = 0.9, sig2 = 16, tau)
#' plot(tau, wv.theo, col = "red")
#' @keywords internal
theoretical_wv <- function(theta, desc, objdesc, tau) {
.Call('gmwm_theoretical_wv', PACKAGE = 'gmwm', theta, desc, objdesc, tau)
}
#' @title Each Models Process Decomposed to WV
#' @description This function computes each process to WV (haar) in a given model.
#' @param theta A \code{vec} containing the list of estimated parameters.
#' @param desc A \code{vector<string>} containing a list of descriptors.
#' @param objdesc A \code{field<vec>} containing a list of object descriptors.
#' @param tau A \code{vec} containing the scales e.g. 2^(1:J)
#' @return A \code{mat} containing the wavelet variance of each process in the model
#' @examples
#' x.sim = gen_ar1( N = 10000, phi = 0.9, sigma2 = 4 )
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = ar1_to_wv(phi = 0.9, sig2 = 16, tau)
#' plot(tau, wv.theo, col = "red")
#' @keywords internal
decomp_theoretical_wv <- function(theta, desc, objdesc, tau) {
.Call('gmwm_decomp_theoretical_wv', PACKAGE = 'gmwm', theta, desc, objdesc, tau)
}
#' @title Decomposed WV to Single WV
#' @description This function computes the combined processes to WV (haar) in a given model.
#' @param decomp A \code{mat} with scales as rows and processes as columns
#' @return A \code{vec} containing the wavelet variance of the process for the overall model
#' @examples
#' x.sim = gen_ar1( N = 10000, phi = 0.9, sigma2 = 4 )
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = ar1_to_wv(phi = 0.9, sig2 = 16, tau)
#' plot(tau, wv.theo, col = "red")
#' @keywords internal
decomp_to_theo_wv <- function(decomp) {
.Call('gmwm_decomp_to_theo_wv', PACKAGE = 'gmwm', decomp)
}
#' @title Read an IMU Binary File into R
#'
#' @description
#' The function will take a file location in addition to the type of sensor it
#' came from and read the data into R.
#'
#' @param file_path A \code{string} that contains the full file path.
#' @param imu_type A \code{string} that contains a supported IMU type given below.
#' @details
#' Currently supports the following IMUs:
#' \itemize{
#' \item IMAR
#' \item LN200
#' \item LN200IG
#' \item IXSEA
#' \item NAVCHIP_INT
#' \item NAVCHIP_FLT
#' }
#'
#' We hope to soon be able to support delimited files.
#' @return A matrix with dimensions N x 7, where the columns represent:
#' \describe{
#' \item{Col 0}{Time}
#' \item{Col 1}{Gyro 1}
#' \item{Col 2}{Gyro 2}
#' \item{Col 3}{Gyro 3}
#' \item{Col 4}{Accel 1}
#' \item{Col 5}{Accel 2}
#' \item{Col 6}{Accel 3}
#' }
#' @references
#' Thanks goes to Philipp Clausen of Labo TOPO, EPFL, Switzerland, topo.epfl.ch, Tel:+41(0)21 693 27 55
#' for providing a matlab function that reads in IMUs.
#' The function below is a heavily modified port of MATLAB code into Armadillo/C++.
#'
#' @examples
#' \dontrun{
#' read_imu(file_path = "F:/Desktop/short_test_data.imu", imu_type = "IXSEA")
#' }
#' @keywords internal
read_imu <- function(file_path, imu_type) {
.Call('gmwm_read_imu', PACKAGE = 'gmwm', file_path, imu_type)
}
#' @title Generate a sequence of values
#' @description Creates a vector containing a sequence of values starting at the initial point and going to the terminal point.
#' @param a An \code{int}, that denotes the starting point.
#' @param b An \code{int}, that denotes the ending point.
#' @return A \code{vector} containing values moving from a to b. There are no restrictions on A's range.
#' @author James J Balamuta
#' @keywords internal
#' @examples
#' #Call with the following data:
#' seq_cpp(3, 5)
#' seq_cpp(5, 3)
seq_cpp <- function(a, b) {
.Call('gmwm_seq_cpp', PACKAGE = 'gmwm', a, b)
}
#' @title Generate a sequence of values based on supplied number
#' @description Creates a vector containing a sequence of values starting at 1 and going to the terminal point.
#' @param n An \code{int} that denotes the length of the vector.
#' @return A \code{vector} containing values moving from 1 to n.
#' @author James J Balamuta
#' @keywords internal
#' @examples
#' #Call with the following data:
#' seq_len_cpp(5)
seq_len_cpp <- function(n) {
.Call('gmwm_seq_len_cpp', PACKAGE = 'gmwm', n)
}
#' @title Find Quantiles
#' @description Attempts to find quantiles
#' @param x A \code{vec} of data
#' @param probs A \code{vec} of the quantiles to find.
#' @return A \code{vector} containing the quantiles
#' @author James J Balamuta
#' @keywords internal
#' @examples
#' #Call with the following data:
#' quantile_cpp(c(1,2,3,4,5,6,7), c(.25,.5,.75))
#' quantile(c(1,2,3,4,5,6,7), c(.25,.5,.75))
quantile_cpp <- function(x, probs) {
.Call('gmwm_quantile_cpp', PACKAGE = 'gmwm', x, probs)
}
#' @title Lagged Differences in Armadillo
#' @description Returns the ith difference of a time series of rth lag.
#' @param x A \code{vec} that is the time series
#' @param lag A \code{unsigned int} that indicates the lag
#' @param differences A \code{dif} that indicates how many differences should be taken
#' @return A \code{vector} containing the differenced time series.
#' @author JJB
#' @keywords internal
#' @examples
#' x = rnorm(10000, 0, 1)
#' diff_cpp(x,1,1)
diff_cpp <- function(x, lag, differences) {
.Call('gmwm_diff_cpp', PACKAGE = 'gmwm', x, lag, differences)
}
#' @title Converting an ARMA Process to an Infinite MA Process
#' @description Takes an ARMA function and converts it to an infinite MA process.
#' @param ar A \code{column vector} of length p
#' @param ma A \code{column vector} of length q
#' @param lag_max A \code{int} of the largest MA(Inf) coefficient required.
#' @return A \code{column vector} containing coefficients
#' @details This function is a port of the base stats package's ARMAtoMA. There is no significant speed difference between the two.
#' @author R Core Team and JJB
#' @keywords internal
#' @examples
#' # ARMA(2,1)
#' ARMAtoMA_cpp(c(1.0, -0.25), 1.0, 10)
#' # ARMA(0,1)
#' ARMAtoMA_cpp(numeric(0), 1.0, 10)
ARMAtoMA_cpp <- function(ar, ma, lag_max) {
.Call('gmwm_ARMAtoMA_cpp', PACKAGE = 'gmwm', ar, ma, lag_max)
}
#' @title Time Series Convolution Filters
#' @description Applies a convolution filter to a univariate time series.
#' @param x A \code{column vector} of length T
#' @param filter A \code{column vector} of length f
#' @param sides An \code{int} that takes either 1:for using past values only or 2: filter coefficients are centered around lag 0.
#' @param circular A \code{bool} that indicates if the filter should be wrapped around the ends of the time series.
#' @return A \code{column vec} that contains the results of the filtering process.
#' @details This is a port of the cfilter function harnessed by the filter function in stats.
#' It is about 5-7 times faster than R's base function. The benchmark was done on iMac Late 2013 using vecLib as the BLAS.
#' @author R Core Team and JJB
#' @keywords internal
#' @examples
#' x = 1:100
#' #
#' cfilter(x, rep(1, 3), sides = 2, circular = FALSE)
#' # Using R's function
#' filter(x, rep(1, 3))
#' #
#' cfilter(x, rep(1, 3), sides = 1, circular = FALSE)
#' # Using R's function
#' filter(x, rep(1, 3), sides = 1)
#' #
#' cfilter(x, rep(1, 3), sides = 1, circular = TRUE)
#' # Using R's function
#' filter(x, rep(1, 3), sides = 1, circular = TRUE)
cfilter <- function(x, filter, sides, circular) {
.Call('gmwm_cfilter', PACKAGE = 'gmwm', x, filter, sides, circular)
}
#' @title Time Series Recursive Filters
#' @description Applies a recursive filter to a univariate time series.
#' @usage rfilter(x, filter, init)
#' @param x A \code{column vector} of length T
#' @param filter A \code{column vector} of length f
#' @param init A \code{column vector} of length f that contains the initial values of the time series in reverse.
#' @return x A \code{column vector} with its contents reversed.
#' @details Note: The length of 'init' must be equal to the length of 'filter'.
#' This is a port of the rfilter function harnessed by the filter function in stats.
#' It is about 6-7 times faster than R's base function. The benchmark was done on iMac Late 2013 using vecLib as the BLAS.
#' @author R Core Team and JJB
#' @keywords internal
#' @examples
#' x = 1:100
#' #
#' rfilter(x, rep(1, 3), rep(1, 3))
#' # Using R's function
#' filter(x, rep(1, 3), method="recursive", init=rep(1, 3))
rfilter <- function(x, filter, init) {
.Call('gmwm_rfilter', PACKAGE = 'gmwm', x, filter, init)
}
#' @title Compute Theoretical ACF for an ARMA Process
#' @description Compute the theoretical autocorrelation function for an ARMA process.
#' @usage ARMAacf_cpp(ar,ma,lag_max)
#' @param ar A \code{vector} of length p containing AR coefficients
#' @param ma A \code{vector} of length q containing MA coefficients
#' @param lag_max A \code{unsigned integer} indicating the maximum lag necessary
#' @return x A \code{matrix} listing values from 1...nx in one column and 1...1, 2...2,....,n...n, in the other
#' @details This is an implementaiton of the ARMAacf function in R. It is approximately 40x times faster. The benchmark was done on iMac Late 2013 using vecLib as the BLAS.
#' @author R Core Team and JJB
#' @keywords internal
#' @examples
#' # ARMA(2,1)
#' ARMAacf_cpp(c(1.0, -0.25), 1.0, lag_max = 10)
#' # ARMA(0,1)
#' ARMAacf_cpp(numeric(0), .35, lag_max = 10)
ARMAacf_cpp <- function(ar, ma, lag_max) {
.Call('gmwm_ARMAacf_cpp', PACKAGE = 'gmwm', ar, ma, lag_max)
}
#' @title Discrete Fourier Transformation for Autocovariance Function
#' @description Calculates the autovariance function (ACF) using Discrete Fourier Transformation.
#' @param x A \code{cx_vec}.
#' @return A \code{vec} containing the ACF.
#' @details
#' This implementation is 2x as slow as Rs.
#' Two issues: 1. memory resize and 2. unoptimized fft algorithm in arma.
#' Consider piping back into R and rewrapping the object. (Decrease of about 10 microseconds.)
#' @keywords internal
#' @examples
#' x=rnorm(100)
#' dft_acf(x)
dft_acf <- function(x) {
.Call('gmwm_dft_acf', PACKAGE = 'gmwm', x)
}
#' @title Mean of the First Difference of the Data
#' @description The mean of the first difference of the data
#' @param x A \code{vec} containing the data
#' @return A \code{double} that contains the mean of the first difference of the data.
#' @keywords internal
#' @examples
#' x=rnorm(100)
#' mean_diff(x)
mean_diff <- function(x) {
.Call('gmwm_mean_diff', PACKAGE = 'gmwm', x)
}
#' @title Routing function for summary info
#' @description Gets all the data for the summary.gmwm function.
#' @param theta A \code{vec} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} containing a list of parameters (e.g. AR(1) = c(1,1), ARMA(p,q) = c(p,q,1))
#' @param model_type A \code{string} that represents the model transformation
#' @param wv_empir A \code{vec} that
#' @param theo A \code{vec} that
#' @param scales A \code{vec} that
#' @param V A \code{mat} that contains the V matrix used to obtain the GMWM.
#' @param omega A \code{mat} that
#' @param obj_value A \code{double} that contains the objective function value at the optimized solution.
#' @param N A \code{int} that indicates how long the time series is.
#' @param alpha A \code{double} that handles the alpha level of the confidence interval (1-alpha)*100
#' @param robust A \code{bool} that indicates whether the estimation should be robust or not.
#' @param eff A \code{double} that specifies the amount of efficiency required by the robust estimator.
#' @param inference A \code{bool} that indicates whether inference (e.g. GoF) should be run.
#' @param fullV A \code{bool} that indicates whether the matrix has been fully bootstrapped.
#' @param bs_gof A \code{bool} indicating whether the GoF should be bootstrapped or done asymptotically.
#' @param bs_gof_p_ci A \code{bool} indicating whether a bootstrapped p-value should be generated during the bootstrapped GoF
#' @param bs_ci A \code{bool} that indicates whether a bootstrapped CI should be obtained or to use analytical derivatives.
#' @param B A \code{int} that indicates how many iterations should take place.
#' @return A \code{field<mat>} that contains bootstrapped / asymptotic GoF results as well as CIs.
#' @keywords internal
get_summary <- function(theta, desc, objdesc, model_type, wv_empir, theo, scales, V, omega, obj_value, N, alpha, robust, eff, inference, fullV, bs_gof, bs_gof_p_ci, bs_theta_est, bs_ci, B) {
.Call('gmwm_get_summary', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, wv_empir, theo, scales, V, omega, obj_value, N, alpha, robust, eff, inference, fullV, bs_gof, bs_gof_p_ci, bs_theta_est, bs_ci, B)
}
#' @title Pseudo Logit Inverse Function
#' @description This function computes the pseudo inverse of a logit transformation of the parameters in order to constrain them to a positive domain
#' @param x A \code{vec} containing real numbers.
#' @return A \code{vec} containing logit probabilities.
#' @keywords internal
#' @examples
#' x.sim = rnorm(100)
#' pseudo_logit_inv(x.sim)
pseudo_logit_inv <- function(x) {
.Call('gmwm_pseudo_logit_inv', PACKAGE = 'gmwm', x)
}
#' @title Logit Inverse Function
#' @description This function computes the inverse of a logit transformation of the parameters.
#' @param x A \code{vec} containing real numbers.
#' @return A \code{vec} containing logit probabilities.
#' @keywords internal
#' @examples
#' x.sim = rnorm(100)
#' logit_inv(x.sim)
logit_inv <- function(x) {
.Call('gmwm_logit_inv', PACKAGE = 'gmwm', x)
}
#' @title Pseudo Logit Function
#' @description This function compute the link function to constrain parameters to a positive domain.
#' @param x A \code{vec} containing probabilities (e.g. 0 <= x <= 1)
#' @return A \code{vec} containing logit terms.
#' @keywords internal
#' @examples
#' x.sim = runif(100)
#' pseudo_logit(x.sim)
pseudo_logit <- function(x) {
.Call('gmwm_pseudo_logit', PACKAGE = 'gmwm', x)
}
#' @title Logit Function
#' @description This function computes the logit link function.
#' @param x A \code{vec} containing probabilities (e.g. -1 <= x <= 1)
#' @return A \code{vec} containing logit terms.
#' @keywords internal
#' @examples
#' x.sim = runif(100)
#' logit(x.sim)
logit <- function(x) {
.Call('gmwm_logit', PACKAGE = 'gmwm', x)
}
#' @title Logit Function
#' @description This function computes the logit link function.
#' @param x A \code{vec} containing probabilities (e.g. -1 <= x <= 1)
#' @return A \code{vec} containing logit terms.
#' @keywords internal
#' @examples
#' x.sim = runif(100)
#' logit(x.sim)
logit2 <- function(x) {
.Call('gmwm_logit2', PACKAGE = 'gmwm', x)
}
#' @title Logit2 Inverse Function
#' @description This function computes the inverse of a logit transformation of the parameters.
#' @param x A \code{vec} containing real numbers.
#' @return A \code{vec} containing logit probabilities.
#' @keywords internal
#' @examples
#' x.sim = rnorm(100)
#' logit_inv(x.sim)
logit2_inv <- function(x) {
.Call('gmwm_logit2_inv', PACKAGE = 'gmwm', x)
}
#' @title Transform Values for Optimization
#' @description Transform parameter guesses prior to estimating with GMWM
#' @template tsobj_cpp
#' @param model_type A \code{string} that contains the model type: \code{"imu"} or \code{"ssm"}
#' @return A \code{vec} containing the transformed guesses.
#' @keywords internal
transform_values <- function(theta, desc, objdesc, model_type) {
.Call('gmwm_transform_values', PACKAGE = 'gmwm', theta, desc, objdesc, model_type)
}
#' @title Revert Transform Values for Display
#' @description Undo the previous transform of parameter guesses to obtain the GMWM estimates.
#' @template tsobj_cpp
#' @param model_type A \code{string} that contains the model type: \code{"imu"} or \code{"ssm"}
#' @return A \code{vec} containing the undone transformation of parameters.
#' @keywords internal
untransform_values <- function(theta, desc, objdesc, model_type) {
.Call('gmwm_untransform_values', PACKAGE = 'gmwm', theta, desc, objdesc, model_type)
}
#' @title Obtain the smallest polynomial root
#' @description Calculates all the roots of a polynomial and returns the root that is the smallest.
#' @param x A \code{cx_vec} that has a 1 appended before the coefficents. (e.g. c(1, x))
#' @return A \code{double} with the minimum root value.
#' @keywords internal
minroot <- function(x) {
.Call('gmwm_minroot', PACKAGE = 'gmwm', x)
}
#' @title Check Invertibility Conditions
#' @description Checks the invertiveness of series of coefficients.
#' @param x A \code{cx_vec} that has a 1 appended before the coefficents. (e.g. c(1, x))
#' @return True (if outside unit circle) || False (if inside unit circle)
#' @keywords internal
invert_check <- function(x) {
.Call('gmwm_invert_check', PACKAGE = 'gmwm', x)
}
#' @title Count Models
#' @description Count the amount of models that exist.
#' @param desc A \code{vector<string>} that contains the model's components.
#' @return A \code{map<string, int>} containing how frequent the model component appears.
#' @keywords internal
count_models <- function(desc) {
.Call('gmwm_count_models', PACKAGE = 'gmwm', desc)
}
#' @title Order AR1s by size of phi.
#' @description Changes the order of AR1s in a string by size.
#' @template tsobj_cpp
#' @return A \code{vec} that has AR1s shown in descending parameter value.
#' @keywords internal
order_AR1s <- function(theta, desc, objdesc) {
.Call('gmwm_order_AR1s', PACKAGE = 'gmwm', theta, desc, objdesc)
}
#' @title Transform AR1 to GM
#' @description
#' Takes AR1 values and transforms them to GM
#' @param theta A \code{vec} that contains AR1 values.
#' @param freq A \code{double} indicating the frequency of the data.
#' @return A \code{vec} containing GM values.
#' @details
#' The function takes a vector of AR1 values \eqn{\phi}{phi} and \eqn{\sigma ^2}{sigma ^2}
#' and transforms them to GM values \eqn{\beta}{beta} and \eqn{\sigma ^2_{gm}}{sigma ^2[gm]}
#' using the formulas:
#' \eqn{\beta = - \frac{{\ln \left( \phi \right)}}{{\Delta t}}}{beta = -ln(phi)/delta_t}
#' \eqn{\sigma _{gm}^2 = \frac{{{\sigma ^2}}}{{1 - {\phi ^2}}} }{sigma^2[gm] = sigma^2/(1-phi^2)}
#' @keywords internal
#' @author JJB
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
#' @examples
#' ar1_to_gm(c(0.3,1,0.6,.3), 2)
ar1_to_gm <- function(theta, freq) {
.Call('gmwm_ar1_to_gm', PACKAGE = 'gmwm', theta, freq)
}
#' @title Transform GM to AR1
#' @description Takes GM values and transforms them to AR1
#' @param theta A \code{vec} that contains AR1 values.
#' @param freq A \code{double} indicating the frequency of the data.
#' @return A \code{vec} containing GM values.
#' @keywords internal
#' @author JJB
#' The function takes a vector of GM values \eqn{\beta}{beta} and \eqn{\sigma ^2_{gm}}{sigma ^2[gm]}
#' and transforms them to AR1 values \eqn{\phi}{phi} and \eqn{\sigma ^2}{sigma ^2}
#' using the formulas:
#' \eqn{\phi = \exp \left( { - \beta \Delta t} \right)}{phi = exp(-beta * delta[t])}
#' \eqn{{\sigma ^2} = \sigma _{gm}^2\left( {1 - \exp \left( { - 2\beta \Delta t} \right)} \right)}{sigma^2 = sigma^2[gm]*(1-exp(-2*beta*delta[t]))}
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
#' @examples
#' gm_to_ar1(c(0.3,1,0.6,.3), 2)
gm_to_ar1 <- function(theta, freq) {
.Call('gmwm_gm_to_ar1', PACKAGE = 'gmwm', theta, freq)
}
#' @title Generate the ts model object description
#' @description Creates the ts.model's obj.desc value
#' @param desc A \code{vector<string>} that contains a list of the strings of each process.
#' @return A \code{field<vec>} that contains the object description of each process.
#' @details
#' This function currently does NOT support ARMA models.
#' That is, there is no support for ARMA, AR, or MA.
#' There is support for AR1, GM, WN, DR, QN, and RW.
#' @keywords internal
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
model_objdesc <- function(desc) {
.Call('gmwm_model_objdesc', PACKAGE = 'gmwm', desc)
}
#' @title Generate the ts model object's theta vector
#' @description Creates the ts.model's theta vector
#' @param desc A \code{vector<string>} that contains a list of the strings of each process.
#' @return A \code{vec} with values initialized at 0 that span the space of parameters to be estimated.
#' @details
#' This function currently does NOT support ARMA models.
#' That is, there is no support for ARMA, AR, or MA.
#' There is support for AR1, GM, WN, DR, QN, and RW.
#' @keywords internal
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
model_theta <- function(desc) {
.Call('gmwm_model_theta', PACKAGE = 'gmwm', desc)
}
#' @title Generate the ts model object's process desc
#' @description Creates the ts.model's process desc
#' @param desc A \code{vector<string>} that contains a list of the strings of each process.
#' @return A \code{vector<string>} with a list of descriptive values to label the estimate matrix with
#' @details
#' This function currently does NOT support ARMA models.
#' That is, there is no support for ARMA, AR, or MA.
#' There is support for AR1, GM, WN, DR, QN, and RW.
#' @keywords internal
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
model_process_desc <- function(desc) {
.Call('gmwm_model_process_desc', PACKAGE = 'gmwm', desc)
}
#' @title Generate eta3 confidence interval
#' @description Computes the eta3 CI
#' @param y A \code{vec} that computes the brickwalled modwt dot product of each wavelet coefficient divided by their length.
#' @param dims A \code{String} indicating the confidence interval being calculated.
#' @param alpha_ov_2 A \code{double} that indicates the \eqn{\left(1-p\right)*\alpha}{(1-p)*alpha} confidence level
#' @return A \code{matrix} with the structure:
#' \itemize{
#' \item{Column 1}{Wavelet Variance}
#' \item{Column 2}{Chi-squared Lower Bounds}
#' \item{Column 3}{Chi-squared Upper Bounds}
#' }
#' @keywords internal
#' @examples
#' x = rnorm(100)
#' # Uses the internal MODWT function not associated with an S3 class.
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' y = wave_variance(decomp)
#' ci_wave_variance(decomp, y, type = "eta3", alpha_ov_2 = 0.025)
ci_eta3 <- function(y, dims, alpha_ov_2) {
.Call('gmwm_ci_eta3', PACKAGE = 'gmwm', y, dims, alpha_ov_2)
}
#' @title Generate eta3 robust confidence interval
#' @description Computes the eta3 robust CI
#' @param wv_robust A \code{vec} that computes the brickwalled modwt dot product of each wavelet coefficient divided by their length.
#' @param wv_ci_class A \code{mat} that contains the CI mean, CI Lower, and CI Upper
#' @param alpha_ov_2 A \code{double} that indicates the \eqn{\left(1-p\right)*\alpha}{(1-p)*alpha} confidence level
#' @param eff A \code{double} that indicates the efficiency.
#' @return A \code{matrix} with the structure:
#' \itemize{
#' \item{Column 1}{Robust Wavelet Variance}
#' \item{Column 2}{Chi-squared Lower Bounds}
#' \item{Column 3}{Chi-squared Upper Bounds}
#' }
#' @details
#' Within this function we are scaling the classical
#' @keywords internal
#' @examples
#' x = rnorm(100)
#' # Uses the internal MODWT function not associated with an S3 class.
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' y = wave_variance(decomp, robust = TRUE, eff = 0.6)
#' ci_wave_variance(decomp, y, type = "eta3", alpha_ov_2 = 0.025, robust = TRUE, eff = 0.6)
ci_eta3_robust <- function(wv_robust, wv_ci_class, alpha_ov_2, eff) {
.Call('gmwm_ci_eta3_robust', PACKAGE = 'gmwm', wv_robust, wv_ci_class, alpha_ov_2, eff)
}
#' @title Generate a Confidence intervval for a Univariate Time Series
#' @description Computes an estimate of the multiscale variance and a chi-squared confidence interval
#' @param signal_modwt_bw A \code{field<vec>} that contains the brick walled modwt or dwt decomposition
#' @param wv A \code{vec} that contains the wave variance.
#' @param type A \code{String} indicating the confidence interval being calculated.
#' @param alpha_ov_2 A \code{double} that indicates the \eqn{\left(1-p\right)*\alpha}{(1-p)*alpha} confidence level.
#' @param robust A \code{boolean} to determine the type of wave estimation.
#' @param eff A \code{double} that indicates the efficiency.
#' @return A \code{matrix} with the structure:
#' \itemize{
#' \item{Column 1}{Wavelet Variance}
#' \item{Column 2}{Chi-squared Lower Bounds}
#' \item{Column 3}{Chi-squared Upper Bounds}
#' }
#' @keywords internal
#' @details
#' This function can be expanded to allow for other confidence interval calculations.
#' @examples
#' set.seed(1337)
#' x = rnorm(100)
#' # Uses the internal MODWT function not associated with an S3 class.
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' y = wave_variance(decomp)
#' ci_wave_variance(decomp, y, type = "eta3", alpha_ov_2 = 0.025)
ci_wave_variance <- function(signal_modwt_bw, wv, type = "eta3", alpha_ov_2 = 0.025, robust = FALSE, eff = 0.6) {
.Call('gmwm_ci_wave_variance', PACKAGE = 'gmwm', signal_modwt_bw, wv, type, alpha_ov_2, robust, eff)
}
#' @title Generate a Wave Variance for a Univariate Time Series
#' @description Computes an estimate of the wave variance
#' @param signal_modwt_bw A \code{field<vec>} that contains the brick walled modwt or dwt decomposition
#' @param robust A \code{boolean} to determine the type of wave estimation.
#' @param eff A \code{double} that indicates the efficiency.
#' @return A \code{vec} that contains the wave variance.
#' @keywords internal
#' @examples
#' set.seed(1337)
#' x = rnorm(100)
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' wave_variance(decomp)
#'
#' wave_variance(decomp, robust = TRUE, eff = 0.6)
wave_variance <- function(signal_modwt_bw, robust = FALSE, eff = 0.6) {
.Call('gmwm_wave_variance', PACKAGE = 'gmwm', signal_modwt_bw, robust, eff)
}
#' @title Computes the (MODWT) wavelet variance
#' @description Calculates the (MODWT) wavelet variance
#' @param signal_modwt_bw A \code{field<vec>} that contains the modwt decomposition after it has been brick walled.
#' @param robust A \code{boolean} that triggers the use of the robust estimate.
#' @param eff A \code{double} that indicates the efficiency as it relates to an MLE.
#' @param alpha A \code{double} that indicates the \eqn{\left(1-p\right)*\alpha}{(1-p)*alpha} confidence level
#' @param ci_type A \code{String} indicating the confidence interval being calculated. Valid value: "eta3"
#' @return A \code{mat} with the structure:
#' \itemize{
#' \item{"variance"}{Wavelet Variance}
#' \item{"low"}{Lower CI}
#' \item{"high"}{Upper CI}
#' }
#' @keywords internal
#' @details
#' This function does the heavy lifting with the signal_modwt_bw
#' @examples
#' x = rnorm(100)
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' wvar_cpp(decomp, robust=FALSE, eff=0.6, alpha = 0.05, ci_type="eta3")
wvar_cpp <- function(signal_modwt_bw, robust, eff, alpha, ci_type) {
.Call('gmwm_wvar_cpp', PACKAGE = 'gmwm', signal_modwt_bw, robust, eff, alpha, ci_type)
}
#' @title Computes the (MODWT) wavelet variance
#' @description Calculates the (MODWT) wavelet variance
#' @param signal A \code{vec} that contains the data.
#' @param robust A \code{boolean} that triggers the use of the robust estimate.
#' @param eff A \code{double} that indicates the efficiency as it relates to an MLE.
#' @param alpha A \code{double} that indicates the \eqn{\left(1-p\right)\times \alpha}{(1-p)*alpha} confidence level
#' @param ci_type A \code{string} indicating the confidence interval being calculated. Valid value: "eta3"
#' @param strWavelet A \code{string} indicating the type of wave filter to be applied. Must be "haar"
#' @param decomp A \code{string} indicating whether to use "modwt" or "dwt" decomp
#' @return A \code{mat} with the structure:
#' \itemize{
#' \item{"variance"}{Wavelet Variance}
#' \item{"low"}{Lower CI}
#' \item{"high"}{Upper CI}
#' }
#' @keywords internal
#' @details
#' This function powers the wvar object. It is also extendable...
#' @examples
#' x=rnorm(100)
#' modwt_wvar_cpp(x, nlevels=4, robust=FALSE, eff=0.6, alpha = 0.05,
#' ci_type="eta3", strWavelet="haar", decomp="modwt")
modwt_wvar_cpp <- function(signal, nlevels, robust, eff, alpha, ci_type, strWavelet, decomp) {
.Call('gmwm_modwt_wvar_cpp', PACKAGE = 'gmwm', signal, nlevels, robust, eff, alpha, ci_type, strWavelet, decomp)
}
#' @title Computes the MO/DWT wavelet variance for multiple processes
#' @description Calculates the MO/DWT wavelet variance
#' @param signal A \code{matrix} that contains the same number of observations per dataset
#' @param robust A \code{boolean} that triggers the use of the robust estimate.
#' @param eff A \code{double} that indicates the efficiency as it relates to an MLE.
#' @param alpha A \code{double} that indicates the \eqn{\left(1-p\right)\times \alpha}{(1-p)*alpha} confidence level
#' @param ci_type A \code{string} indicating the confidence interval being calculated. Valid value: "eta3"
#' @param strWavelet A \code{string} indicating the type of wave filter to be applied. Must be "haar"
#' @param decomp A \code{string} indicating whether to use "modwt" or "dwt" decomp
#' @return A \code{field<mat>} with the structure:
#' \itemize{
#' \item{"variance"}{Wavelet Variance}
#' \item{"low"}{Lower CI}
#' \item{"high"}{Upper CI}
#' }
#' @keywords internal
#' @details
#' This function processes the decomposition of multiple signals quickly
#' @examples
#' x = cbind(rnorm(100),rnorm(100))
#' batch_modwt_wvar_cpp(x, nlevels=4, robust=FALSE, eff=0.6,
#' alpha = 0.05, ci_type="eta3", strWavelet="haar",
#' decomp="modwt")
batch_modwt_wvar_cpp <- function(signal, nlevels, robust, eff, alpha, ci_type, strWavelet, decomp) {
.Call('gmwm_batch_modwt_wvar_cpp', PACKAGE = 'gmwm', signal, nlevels, robust, eff, alpha, ci_type, strWavelet, decomp)
}
#' @title Computes the MODWT scales
#' @description Calculates the MODWT scales
#' @param nb_level A \code{integer} that contains the level of decomposition J.
#' @return A \code{vec} that contains 2^1, ... , 2^J
#' @keywords internal
#' @details
#' Used in wvar object.
#' @examples
#' scales_cpp(5)
scales_cpp <- function(nb_level) {
.Call('gmwm_scales_cpp', PACKAGE = 'gmwm', nb_level)
}
#' @title Quadrature Mirror Filter
#' @description Calculate the series quadrature mirror filter (QMF). Requires a series of an even length.
#' @usage qmf(g, inverse)
#' @param g A \code{vector} that contains the filter constants.
#' @param inverse A \code{bool} that indicates whether the inverse quadrature mirror filter is computed.
#' By default, the inverse quadrature mirror is computed.
#' @return A \code{vector} that contains either the forward QMF (evalute in order) or the inverse QMF (reverse order).
#' @author JJB
#' @keywords internal
#' @examples
#' # Haar values
#' g = rep(1/sqrt(2),2)
#' qmf(g)
qmf <- function(g, inverse = TRUE) {
.Call('gmwm_qmf', PACKAGE = 'gmwm', g, inverse)
}
#' @title Haar filter construction
#' @description Creates the haar filter
#' @usage haar_filter()
#' @return A \code{field<vec>} that contains:
#' \itemize{
#' \item{"L"}{A \code{integer} specifying the length of the filter}
#' \item{"h"}{A \code{vector} containing the coefficients for the wavelet filter}
#' \item{"g"}{A \code{vector} containing the coefficients for the scaling filter}
#' }
#' @details
#' This template can be used to increase the amount of filters available for selection.
#' @author JJB
#' @keywords internal
#' @examples
#' haar_filter()
haar_filter <- function() {
.Call('gmwm_haar_filter', PACKAGE = 'gmwm')
}
#' @title Select the Wavelet Filter
#' @description Constructs the wavelet filter to be used.
#' @usage select_filter(filter_name)
#' @param filter_name A \code{String} that must receive: \code{"haar"}.
#' @return info A \code{field<vec>} that contains:
#' \itemize{
#' \item{"L"}{A \code{integer} specifying the length of the filter}
#' \item{"h"}{A \code{vector} containing the coefficients for the wavelet filter}
#' \item{"g"}{A \code{vector} containing the coefficients for the scaling filter}
#' }
#' @details
#' The package is oriented toward using only the haar filter. If the package extends at a later time, then the supporting infrastructure is there.
#' @author JJB
#' @keywords internal
#' @examples
#' select_filter("haar")
select_filter <- function(filter_name = "haar") {
.Call('gmwm_select_filter', PACKAGE = 'gmwm', filter_name)
}
Try the gmwm package in your browser
Any scripts or data that you put into this service are public.
gmwm documentation built on April 14, 2017, 4:38 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
# This file was generated by Rcpp::compileAttributes
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393
#' @title Expected value DR
#' @description This function computes the expected value of a drift process.
#' @param omega A \code{double} corresponding to variance of drift.
#' @param n_ts An \code{int} indicating the length of the time series.
#' @return A \code{vec} containing the expected value of the drift.
#' @keywords internal
#' @examples
#' e_drift(1,200)
e_drift <- function(omega, n_ts) {
.Call('gmwm_e_drift', PACKAGE = 'gmwm', omega, n_ts)
}
#' @title Second moment DR
#' @description This function computes the second moment of a drift process.
#' @param omega A \code{double} corresponding to variance of drift.
#' @param n_ts An \code{int} indicating the length of the time series.
#' @return A \code{vec} containing the second moment of the drift.
#' @keywords internal
#' @examples
#' m2_drift(1, 200)
m2_drift <- function(omega, n_ts) {
.Call('gmwm_m2_drift', PACKAGE = 'gmwm', omega, n_ts)
}
#' @title Variance DR
#' @description This function computes the variance of a drift process.
#' @param omega A \code{double} corresponding to variance of drift.
#' @param n_ts An \code{int} indicating the length of the time series.
#' @return A \code{vec} containing the variance of the drift.
#' @keywords internal
#' @examples
#' var_drift(1, 200)
var_drift <- function(omega, n_ts) {
.Call('gmwm_var_drift', PACKAGE = 'gmwm', omega, n_ts)
}
#' @title Compute Tau-Overlap Allan Variance
#' @description Computation of Tau-Overlap Allan Variance
#' @usage avar_to_cpp(x)
#' @param x A \code{vector} with dimensions N x 1.
#' @return av A \code{matrix} that contains:
#' \itemize{
#' \item{Col 1}{The size of the cluster}
#' \item{Col 2}{The Allan variance}
#' \item{Col 3}{The error associated with the variance estimation.}
#' }
#' @details
#' Given \eqn{N} equally spaced samples with averaging time \eqn{\tau = n\tau _0}{tau = n*tau_0},
#' where \eqn{n} is an integer such that \eqn{ 1 \le n \le \frac{N}{2}}{1<= n <= N/2}.
#' Therefore, \eqn{n} is able to be selected from \eqn{\left\{ {n|n < \left\lfloor {{{\log }_2}\left( N \right)} \right\rfloor } \right\}}{{n|n< floor(log2(N))}}
#' Then, a sampling of \eqn{m = \left\lfloor {\frac{{N - 1}}{n}} \right\rfloor - 1} samples exist.
#' The tau-overlap estimator is given by:
#'
#' where \eqn{ {{\bar y}_t}\left( \tau \right) = \frac{1}{\tau }\sum\limits_{i = 0}^{\tau - 1} {{{\bar y}_{t - i}}} }.
#'
#' @author JJB
#' @references Long-Memory Processes, the Allan Variance and Wavelets, D. B. Percival and P. Guttorp
#' @examples
#' set.seed(999)
#' # Simulate white noise (P 1) with sigma^2 = 4
#' N = 100000
#' white.noise = rnorm(N, 0, 2)
#' #plot(white.noise,ylab="Simulated white noise process",xlab="Time",type="o")
#' #Simulate random walk (P 4)
#' random.walk = cumsum(0.1*rnorm(N, 0, 2))
#' combined.ts = white.noise+random.walk
#' av_mat = avar_to_cpp(combined.ts)
#' @keywords internal
avar_to_cpp <- function(x) {
.Call('gmwm_avar_to_cpp', PACKAGE = 'gmwm', x)
}
#' @title Compute Maximal-Overlap Allan Variance using Means
#' @description Computation of Maximal-Overlap Allan Variance
#' @usage avar_mo_cpp(x)
#' @param x A \code{vector} with dimensions N x 1.
#' @return av A \code{list} that contains:
#' \itemize{
#' \item{"clusters"}{The size of the cluster}
#' \item{"allan"}{The Allan variance}
#' \item{"errors"}{The error associated with the variance estimation.}
#' }
#' @details
#' Given \eqn{N} equally spaced samples with averaging time \eqn{\tau = n\tau _0}{tau = n*tau_0},
#' where \eqn{n} is an integer such that \eqn{ 1 \le n \le \frac{N}{2}}{1<= n <= N/2}.
#' Therefore, \eqn{n} is able to be selected from \eqn{\left\{ {n|n < \left\lfloor {{{\log }_2}\left( N \right)} \right\rfloor } \right\}}{{n|n< floor(log2(N))}}
#' Then, \eqn{M = N - 2n} samples exist.
#' The Maximal-overlap estimator is given by:
#' \eqn{\frac{1}{{2\left( {N - 2k + 1} \right)}}\sum\limits_{t = 2k}^N {{{\left[ {{{\bar Y}_t}\left( k \right) - {{\bar Y}_{t - k}}\left( k \right)} \right]}^2}} }
#'
#' where \eqn{ {{\bar y}_t}\left( \tau \right) = \frac{1}{\tau }\sum\limits_{i = 0}^{\tau - 1} {{{\bar y}_{t - i}}} }.
#' @author JJB
#' @references Long-Memory Processes, the Allan Variance and Wavelets, D. B. Percival and P. Guttorp
#' @examples
#' set.seed(999)
#' # Simulate white noise (P 1) with sigma^2 = 4
#' N = 100000
#' white.noise = rnorm(N, 0, 2)
#' #plot(white.noise,ylab="Simulated white noise process",xlab="Time",type="o")
#' #Simulate random walk (P 4)
#' random.walk = cumsum(0.1*rnorm(N, 0, 2))
#' combined.ts = white.noise+random.walk
#' av_mat = avar_mo_cpp(combined.ts)
#' @keywords internal
avar_mo_cpp <- function(x) {
.Call('gmwm_avar_mo_cpp', PACKAGE = 'gmwm', x)
}
#' @title ARMA Adapter to ARMA to WV Process function
#' @description Molds the data so that it works with the arma_to_wv function.
#' @param theta A \code{vec} that contains all the parameter estimates.
#' @param p A \code{int} that indicates the number of AR coefficients
#' @param q A \code{int} that indicates the number of MA coefficients.
#' @param tau A \code{vec} that lists the scales of the process e.g. 2^(1:J)
#' @return A \code{vec} containing the ARMA to WV results
#' @keywords internal
#' @backref src/analytical_matrix_derivatives.cpp
#' @backref src/analytical_matrix_derivatives.h
arma_adapter <- function(theta, p, q, tau) {
.Call('gmwm_arma_adapter', PACKAGE = 'gmwm', theta, p, q, tau)
}
#' @title Calculates the Jacobian for the ARMA process
#' @description Figure out the Jacobian for an ARMA process.
#' @param theta A \code{vec} that contains all the parameter estimates.
#' @param p A \code{int} that indicates the number of AR coefficients
#' @param q A \code{int} that indicates the number of MA coefficients.
#' @param tau A \code{vec} that lists the scales of the process e.g. 2^(1:J)
#' @return A \code{mat} that returns the numerical jacobian of the ARMA process.
#' @keywords internal
#' @backref src/analytical_matrix_derivatives.cpp
#' @backref src/analytical_matrix_derivatives.h
jacobian_arma <- function(theta, p, q, tau) {
.Call('gmwm_jacobian_arma', PACKAGE = 'gmwm', theta, p, q, tau)
}
#' Analytic D matrix for AR(1) process
#' @param phi A \code{double} corresponding to the phi coefficient of an AR(1) process.
#' @param sig2 A \code{double} corresponding to the error term of an AR(1) process.
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{\phi ^2}{sigma^2} and the second column contains the partial derivative with respect to \eqn{\sigma ^2}{sigma^2}
#' @details
#' The haar wavelet variance is given as \eqn{\frac{{\left( {\frac{\tau }{2} - 3{\rho _0} - \frac{{\tau \rho _0^2}}{2} + 4\rho _0^{\frac{\tau }{2} + 1} - \rho _0^{\tau + 1}} \right)\nu _0^2}}{{\frac{{{\tau ^2}}}{8}{{\left( {1 - {\rho _0}} \right)}^2}\left( {1 - \rho _0^2} \right)}}}{See PDF Manual for equation}
#' Note: \eqn{\phi = \rho}{phi = rho} and \eqn{V _0^2 = \sigma _0^2}{V[0]^2 = sigma[0]^2}.
#' Due to length, the analytical derivations of the AR(1) haar wavelet variance are given in a supplied file within vignette.
#' @author JJB
#' @examples
#' deriv_ar1(.3, 1, 2^(1:5))
deriv_ar1 <- function(phi, sig2, tau) {
.Call('gmwm_deriv_ar1', PACKAGE = 'gmwm', phi, sig2, tau)
}
#' Analytic second derivative matrix for AR(1) process
#' @param phi A \code{double} corresponding to the phi coefficient of an AR(1) process.
#' @param sig2 A \code{double} corresponding to the error term of an AR(1) process.
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the second partial derivative with respect to \eqn{\phi ^2}{sigma^2} and the second column contains the second partial derivative with respect to \eqn{\sigma ^2}{sigma^2}
#' @details
#' The haar wavelet variance is given as \eqn{\frac{{\left( {\frac{\tau }{2} - 3{\rho _0} - \frac{{\tau \rho _0^2}}{2} + 4\rho _0^{\frac{\tau }{2} + 1} - \rho _0^{\tau + 1}} \right)\nu _0^2}}{{\frac{{{\tau ^2}}}{8}{{\left( {1 - {\rho _0}} \right)}^2}\left( {1 - \rho _0^2} \right)}}}{See PDF Manual for equation}
#' Note: \eqn{\phi = \rho}{phi = rho} and \eqn{V _0^2 = \sigma _0^2}{V[0]^2 = sigma[0]^2}.
#' Due to length, the analytical derivations of the AR(1) haar wavelet variance are given in a supplied file within vignette.
#' @author JJB
#' @examples
#' deriv_2nd_ar1(.3, 1, 2^(1:5))
deriv_2nd_ar1 <- function(phi, sig2, tau) {
.Call('gmwm_deriv_2nd_ar1', PACKAGE = 'gmwm', phi, sig2, tau)
}
#' Analytic D matrix for drift process
#' @param omega A \code{double} that is the slope of the drift.
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{\omega _0}{omega[0]}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{{\tau ^2}\omega _0^2}}{2}}{nu^2(tau) = tau^2 omega_0^2 / 2}.
#' Taking the derivative with respect to \eqn{\omega _0^2}{omega_0^2} yields: \eqn{\frac{\partial }{{\partial {\omega _0}}}{\nu ^2}\left( \tau \right) = {\tau ^2}{\omega _0}}{tau^2 * omega_0}
#' @author JJB
#' @examples
#' deriv_dr(5.3, 2^(1:5))
deriv_dr <- function(omega, tau) {
.Call('gmwm_deriv_dr', PACKAGE = 'gmwm', omega, tau)
}
#' Analytic second derivative matrix for drift process
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the second partial derivative with respect to \eqn{\omega _0}{omega[0]}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{{\tau ^2}\omega _0^2}}{2}}{nu^2(tau) = tau^2 omega_0^2 / 2}.
#' Taking the derivative with respect to \eqn{\omega _0^2}{omega_0^2} yields: \eqn{\frac{\partial }{{\partial {\omega _0}}}{\nu ^2}\left( \tau \right) = {\tau ^2}{\omega _0}}{tau^2 * omega_0}
#' Taking second derivative with respect to \eqn{\omega _0^2}{omega_0^2} yields: \eqn{\frac{{{\partial ^2}}}{{\partial \omega _0^2}}{\nu ^2}\left( \tau \right) = {\tau ^2}}{tau^2}
#' @author JJB
#' @examples
#' deriv_2nd_dr(2^(1:5))
deriv_2nd_dr <- function(tau) {
.Call('gmwm_deriv_2nd_dr', PACKAGE = 'gmwm', tau)
}
#' Analytic D matrix quantisation noise process
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{Q _0^2}{Q[0]^2}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{3Q_0^2}}{{2{\tau ^2}}}}{nu^2(tau) = 3*Q[0]^2 / 2*tau^2}.
#' Taking the derivative with respect to \eqn{Q _0^2}{Q[0]^2} yields: \deqn{\frac{\partial }{{\partial Q_0^2}}{\nu ^2}\left( \tau \right) = \frac{3}{{2{\tau ^2}}}}{3/(2*tau^2)}.
#' The second derivative derivative with respect to \eqn{Q _0^2}{Q[0]^2} is then: \deqn{\frac{{{\partial ^2}}}{{\partial Q_0^4}}{\nu ^2}\left( \tau \right) = 0}{0}.
#' @author JJB
#' @examples
#' deriv_qn(2^(1:5))
deriv_qn <- function(tau) {
.Call('gmwm_deriv_qn', PACKAGE = 'gmwm', tau)
}
#' Analytic D matrix random walk process
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{\gamma _0^2}{gamma[0]^2}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{\left( {2{\tau ^2} + 1} \right)\gamma _0^2}}{{24\tau }}}{nu^2(tau) = (2*tau^2+1)*gamma^2 / (24*tau)}.
#' Taking the first derivative with respect to \eqn{\gamma _0^2}{gamma_0^2} yields: \deqn{\frac{{{\partial ^2}}}{{\partial \gamma _0^4}}{\nu ^2}\left( \tau \right) = 0}{(2*tau^2+1) / (24*tau)}
#' The second derivative derivative with respect to \eqn{\gamma _0^2}{gamma[0]^2} is then: \deqn{\frac{{{\partial ^2}}}{{\partial \sigma_0^4}}{\nu ^2}\left( \tau \right) = 0}{0}.
#' @author JJB
#' @examples
#' deriv_rw(2^(1:5))
deriv_rw <- function(tau) {
.Call('gmwm_deriv_rw', PACKAGE = 'gmwm', tau)
}
#' Analytic D matrix white noise process
#' @param tau A \code{vec} that contains the scales to be processed (e.g. 2^(1:J))
#' @return A \code{matrix} with the first column containing the partial derivative with respect to \eqn{\sigma _0^2}{sigma[0]^2}.
#' @details
#' The haar wavelet variance is given as \eqn{{\nu ^2}\left( \tau \right) = \frac{{\sigma _0^2}}{\tau }}{nu^2(tau) = sigma_0^2 / tau}.
#' Taking the derivative with respect to \eqn{\sigma _0^2}{sigma_0^2} yields: \eqn{\frac{\partial }{{\partial \sigma _0^2}}{\nu ^2}\left( \tau \right) = \frac{1}{\tau }}{1/tau}
#' @author JJB
#' @examples
#' deriv_wn(2^(1:5))
deriv_wn <- function(tau) {
.Call('gmwm_deriv_wn', PACKAGE = 'gmwm', tau)
}
#' Analytic D matrix of Processes
#' @description This function computes each process to WV (haar) in a given model.
#' @param theta A \code{vec} containing the list of estimated parameters.
#' @param desc A \code{vector<string>} containing a list of descriptors.
#' @param objdesc A \code{field<vec>} containing a list of object descriptors.
#' @param tau A \code{vec} containing the scales e.g. 2^(1:J)
#' @return A \code{matrix} with the process derivatives going down the column
#' @details
#' Function returns the matrix effectively known as "D"
#' @author JJB
#' @examples
#' mod = AR1(.4,1) + WN(.2) + DR(.005)
#' derivative_first_matrix(mod$theta, mod$desc, mod$obj.desc, 2^(1:9))
derivative_first_matrix <- function(theta, desc, objdesc, tau) {
.Call('gmwm_derivative_first_matrix', PACKAGE = 'gmwm', theta, desc, objdesc, tau)
}
#' Analytic D matrix of Processes
#' @description This function computes each process to WV (haar) in a given model.
#' @param theta A \code{vec} containing the list of estimated parameters.
#' @param desc A \code{vector<string>} containing a list of descriptors.
#' @param objdesc A \code{field<vec>} containing a list of object descriptors.
#' @param tau A \code{vec} containing the scales e.g. 2^(1:J)
#' @param omegadiff A \code{vec} that contains the result of Omega * (wv_empir - wv_theo)
#' @return A \code{matrix} with the process derivatives going down the column
#' @details
#' Function returns the matrix effectively known as "D"
#' @author JJB
#' @examples
#' #TBA
#' @keywords internal
D_matrix <- function(theta, desc, objdesc, tau, omegadiff) {
.Call('gmwm_D_matrix', PACKAGE = 'gmwm', theta, desc, objdesc, tau, omegadiff)
}
#' @title Hook into R's ARIMA function
#' @description Uses R's ARIMA function to obtain CSS values for starting condition
#' @param data A \code{vec} of data.
#' @param params A \code{vec} of the ARMA parameters
#' @return A \code{vec} containing the CSS of the ARMA parameters.
#' @keywords internal
Rcpp_ARIMA <- function(data, params) {
.Call('gmwm_Rcpp_ARIMA', PACKAGE = 'gmwm', data, params)
}
#' Sort Matrix by Column
#'
#' Sorts a given matrix by a specific column while retain the elements in each row.
#'
#' @param x A \code{matrix} to sort
#' @param col A \code{int} that indicates the column the matrix should sort by.
#' @details
#' The functional difference between armadillo's sort() and sort_mat() is straight forward.
#' sort() will sort each column without respect to the rows.
#' Using sort_matrix will sort only 1 column and retain the other elements to be in the same row.
#' @return The matrix sorted by values in the specified column.
#' @keywords internal
sort_mat <- function(x, col) {
.Call('gmwm_sort_mat', PACKAGE = 'gmwm', x, col)
}
#' @title Reverse Subset Column
#' @description
#' Subsets the column by going from high indices to low (the reverse of the supported practice)
#' @param x A \code{matrix} of dimensions M x N
#' @param start A \code{unsigned int} that indicates the starting column.
#' @param end A \code{unsigned int} that indicates the ending column.
#' @return x A \code{matrix} with matrix rows displayed in reverse order
#' @details Consider a vector x=[[1,2],[3,4]].
#' By setting \code{start=1} and \code{end=0}, the function would output x=[[2,1],[4,1]].
#' Start and end must be valid C++ matrix locations. (e.g. matrix cols start at 0 and not 1)
#' @author JJB
#' @examples
#' x = matrix(c(1,2,3,4), nrow = 2,byrow = TRUE)
#' rev_col_subset(x, 1, 0)
#' @keywords internal
rev_col_subset <- function(x, start, end) {
.Call('gmwm_rev_col_subset', PACKAGE = 'gmwm', x, start, end)
}
#' @title Reverse Subset Row
#' @description Subsets the row by going from high indices to low (the reverse of the supported practice)
#' @param x A \code{matrix} of dimensions M x N
#' @param start A \code{unsigned int} that indicates the starting row.
#' @param end A \code{unsigned int} that indicates the ending row.
#' @return x A \code{matrix} with matrix rows displayed in reversed order
#' @details Consider a vector x=[[1,2],[3,4]], the function would output x=[[3,4],[1,2]].
#' Start and end must be valid C++ matrix locations. (e.g. matrix rows start at 0 and not 1)
#' @author JJB
#' @examples
#' x = matrix(c(1,2,3,4), nrow=2,byrow=TRUE)
#' rev_row_subset(x, 1, 0)
#' @keywords internal
rev_row_subset <- function(x, start, end) {
.Call('gmwm_rev_row_subset', PACKAGE = 'gmwm', x, start, end)
}
#' @title Reverse Armadillo Vector
#' @description Reverses the order of an Armadillo Vector
#' @usage reverse_vec(x)
#' @param x A \code{column vector} of length N
#' @return x A \code{column vector} with its contents reversed.
#' @details Consider a vector x=[1,2,3,4,5], the function would output x=[5,4,3,2,1].
#' @author JJB
#' @examples
#' x = 1:5
#' reverse_vec(x)
#' @keywords internal
reverse_vec <- function(x) {
.Call('gmwm_reverse_vec', PACKAGE = 'gmwm', x)
}
#' @title Transform an Armadillo field<vec> to a matrix
#' @description Unlists vectors in a field and places them into a matrix
#' @param x A \code{field<vec>}.
#' @return A \code{mat} containing the field elements within a column.
#' @author JJB
#' @examples
#' x=rnorm(100)
#' @keywords internal
field_to_matrix <- function(x) {
.Call('gmwm_field_to_matrix', PACKAGE = 'gmwm', x)
}
#' @title Accumulation of Armadillo field<vec>
#' @description Sums vectors in a field into a single variable.
#' @param x A \code{field<vec>}.
#' @return An \code{mat} containing the field elements within a column.
#' @author JJB
#' @examples
#' x=rnorm(100)
#' @keywords internal
sum_field_vec <- function(x) {
.Call('gmwm_sum_field_vec', PACKAGE = 'gmwm', x)
}
#' @title Build List of Unique Models
#' @description Creates a set containing unique strings.
#' @param combs A \code{mat} that is a binary matrix (0,1) containing the combinations of different variables.
#' @param x A \code{vec<string>} that contains a list of model descriptors.
#' @return A \code{set<string>} that contains the list of unique models.
#' @keywords internal
build_model_set <- function(combs, x) {
.Call('gmwm_build_model_set', PACKAGE = 'gmwm', combs, x)
}
#' @title Conversion function of Vector to Set
#' @description Converts a vector into a set
#' @param x A \code{vec<vec<string>>} that contains a list of model descriptors.
#' @return A \code{set<vector<string>>} that contains the list of unique models.
#' @keywords internal
vector_to_set <- function(model_str) {
.Call('gmwm_vector_to_set', PACKAGE = 'gmwm', model_str)
}
#' @title Find the Common Denominator of the Models
#' @description Determines the common denominator among models
#' @param x A \code{vector< vector<string> >} that contains all possible models under consideration
#' @return A \code{vector<string>} that contains the terms of the common denominator of all models
#' @keywords internal
find_full_model <- function(x) {
.Call('gmwm_find_full_model', PACKAGE = 'gmwm', x)
}
#' @title Find the Rank Models result
#' @description Provides the core material to create an S3 object for rank.models
#' @param data A \code{vec} of data.
#' @param model_str A \code{vector<vector<string>>} that gives a list of models to test.
#' @param full_model A \code{vector<string>} that contains the largest / full model.
#' @param alpha A \code{double} that indicates the alpha level for CIs.
#' @param compute_v A \code{string} indicating the type of V matrix to generate
#' @param model_type A \code{string} that describes the model generation / transformation: 'ssm' or 'imu'
#' @param K A \code{int} that controls how many times the GMWM is run.
#' @param H A \code{int} that controls how many bootstraps occur.
#' @param G A \code{int} that controls how many guesses occur.
#' @param robust A \code{bool} that indicates whether to use classical or robust wavelet variance.
#' @param eff A \code{double} that indicates the efficiency to use.
#' @param bs_optimism A \code{bool} that indicates whether the model selection score should be calculated with bootstrap or asymptotics.
#' @return A \code{field<field<field<mat>>>} that contains the model score matrix and the best GMWM model object.
#' @keywords internal
rank_models <- function(data, model_str, full_model, alpha, compute_v, model_type, K, H, G, robust, eff, bs_optimism) {
.Call('gmwm_rank_models', PACKAGE = 'gmwm', data, model_str, full_model, alpha, compute_v, model_type, K, H, G, robust, eff, bs_optimism)
}
#' @title Find the auto imu result
#' @description Provides the core material to create an S3 object for auto.imu
#' @param data A \code{mat} containing multiple columns of independent data with the same number of observations.
#' @param model_str A \code{vector<vector<string>>} that gives a list of models to test.
#' @param full_model A \code{vector<string>} that contains the largest / full model.
#' @param alpha A \code{double} that indicates the alpha level for CIs.
#' @param compute_v A \code{string} indicating the type of V matrix to generate
#' @param model_type A \code{string} that describes the model generation / transformation: 'ssm' or 'imu'
#' @param K A \code{int} that controls how many times the GMWM is run.
#' @param H A \code{int} that controls how many bootstraps occur.
#' @param G A \code{int} that controls how many guesses occur.
#' @param robust A \code{bool} that indicates whether to use classical or robust wavelet variance.
#' @param eff A \code{double} that indicates the efficiency to use.
#' @param bs_optimism A \code{bool} that indicates whether the model selection score should be calculated with bootstrap or asymptotics.
#' @return A \code{field<field<field<mat>>>} that contains the model score matrix and the best GMWM model object.
#' @keywords internal
auto_imu <- function(data, combs, full_model, alpha, compute_v, model_type, K, H, G, robust, eff, bs_optimism) {
.Call('gmwm_auto_imu', PACKAGE = 'gmwm', data, combs, full_model, alpha, compute_v, model_type, K, H, G, robust, eff, bs_optimism)
}
#' @title Bootstrap for Matrix V
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
cov_bootstrapper <- function(theta, desc, objdesc, N, robust, eff, H, diagonal_matrix) {
.Call('gmwm_cov_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, N, robust, eff, H, diagonal_matrix)
}
#' @title Bootstrap for Optimism
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
optimism_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_optimism_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Bootstrap for Optimism and GoF
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @param scales A \code{vec} containing the scales of the process.
#' @param model_type A \code{string} containing the model type either: SSM or IMU
#' @param N A \code{int} indicating how long the integer is.
#' @param robust A \code{bool} indicating robust (T) or classical (F).
#' @param eff A \code{double} that handles efficiency.
#' @param H A \code{int} that indicates how many bootstraps should be obtained.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
opt_n_gof_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_opt_n_gof_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Bootstrap for Standard Deviations of Theta Estimates
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
gmwm_sd_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_gmwm_sd_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Generate the Confidence Interval for GOF Bootstrapped
#' @description yaya
#' @param obj A \code{double} containing the objective value of the solution.
#' @param obj_boot A \code{vec} containing the objective values obtained while bootstrapping under Theta_hat.
#' @param B An \code{int} indicating how many times the bootstrapper should be run.
#' @param alpha A \code{double} indicating the amount of confidence for CI
#' @return A \code{vec} that has the alpha/2.0 quantile and then the 1-alpha/2.0 quantile.
#' @keywords internal
boot_pval_gof <- function(obj, obj_boot, B = 1000L, alpha = 0.05) {
.Call('gmwm_boot_pval_gof', PACKAGE = 'gmwm', obj, obj_boot, B, alpha)
}
#' @title Bootstrap for Estimating Both Theta and Theta SD
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
gmwm_param_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_gmwm_param_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Bootstrap for Everything!
#' @description Using the bootstrap approach, we simulate a model based on user supplied parameters, obtain the wavelet variance, and then V.
#' @param theta A \code{vector} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' Expand in detail...
#' @author JJB
#' @keywords internal
#' @examples
#' # Coming soon
all_bootstrapper <- function(theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H) {
.Call('gmwm_all_bootstrapper', PACKAGE = 'gmwm', theta, desc, objdesc, scales, model_type, N, robust, eff, alpha, H)
}
#' @title Absolute Value or Modulus of a Complex Number Squared.
#' @description Computes the squared value of the Modulus.
#' @param x A \code{cx_vec}.
#' @return A \code{vec} containing the modulus squared for each element.
#' @details Consider a complex number defined as: \eqn{z = x + i y} with real \eqn{x} and \eqn{y},
#' The modulus is defined as: \eqn{r = Mod\left(z\right) = \sqrt{\left(x^2 + y^2\right)}}{r = Mod(z) = sqrt(x^2 + y^2)}
#' This function will return: \eqn{r^2 = Mod\left(z\right)^2 = x^2 + y^2}
#' @examples
#' Mod_squared_cpp(c(1+.5i, 2+1i, 5+9i))
#' @keywords internal
Mod_squared_cpp <- function(x) {
.Call('gmwm_Mod_squared_cpp', PACKAGE = 'gmwm', x)
}
#' @title Absolute Value or Modulus of a Complex Number.
#' @description Computes the value of the Modulus.
#' @param x A \code{cx_vec}.
#' @return A \code{vec} containing the modulus for each element.
#' @details Consider a complex number defined as: \eqn{z = x + i y} with real \eqn{x} and \eqn{y},
#' The modulus is defined as: \eqn{r = Mod(z) = \sqrt{(x^2 + y^2)}}
#' @examples
#' Mod_cpp(c(1+.5i, 2+1i, 5+9i))
#' @keywords internal
Mod_cpp <- function(x) {
.Call('gmwm_Mod_cpp', PACKAGE = 'gmwm', x)
}
#' @title Computes the (MODWT) wavelet covariance matrix
#' @description Calculates the (MODWT) wavelet covariance matrix
#' @param signal_modwt A \code{field<vec>} that contains the modwt decomposition.
#' @param nb_level A \code{integer} that contains the level of decomposition J.
#' @param compute_v A \code{string} that indicates what kind of matrix should be created. Possible options: "diag" or "none"
#' @param robust A \code{boolean} that triggers the use of the robust estimate.
#' @param eff A \code{double} that indicates the efficiency as it relates to an MLE.
#' @return A \code{field<mat>} containing the covariance matrix.
#' @keywords internal
#' @examples
#' \dontrun{
#' x=rnorm(100)
#' decomp = modwt(x)
#' V = compute_cov_cpp(decomp$data, decomp$nlevels, compute_v="diag", robust = TRUE, eff=0.6)
#' }
compute_cov_cpp <- function(signal_modwt, nb_level, compute_v = "diag", robust = TRUE, eff = 0.6) {
.Call('gmwm_compute_cov_cpp', PACKAGE = 'gmwm', signal_modwt, nb_level, compute_v, robust, eff)
}
#' @title Computes the (MODWT) wavelet covariance matrix using Chi-square confidence interval bounds
#' @description Calculates the (MODWT) wavelet covariance matrix using Chi-square confidence interval bounds
#' @param ci_hi A \code{vec} that contains the upper confidence interval points.
#' @param ci_lo A \code{vec} that contains the lower confidence interval points.
#' @return A diagonal matrix.
#' @keywords internal
#' @examples
#' \dontrun{
#' x=runif(100)
#' y=x+3
#' fast_cov_cpp(y,x)
#' }
fast_cov_cpp <- function(ci_hi, ci_lo) {
.Call('gmwm_fast_cov_cpp', PACKAGE = 'gmwm', ci_hi, ci_lo)
}
#' @title Discrete Wavelet Transform
#' @description Calculation of the coefficients for the discrete wavelet transformation.
#' @param x A \code{vector} with dimensions \eqn{N\times 1}{N x 1}.
#' @param filter_name A \code{string} indicating the filter.
#' @param nlevels An \code{integer}, \eqn{J}, indicating the level of the decomposition.
#' @param boundary A \code{string} indicating the type of boundary method to use. Either \code{boundary="periodic"} or \code{"reflection"}.
#' @param brickwall A \code{bool} indicating whether the a brick wall procedure should be applied to the coefficients.
#' @return y A \code{field<vec>} that contains the wavelet coefficients for each decomposition level
#' @details
#' Performs a level J decomposition of the time series using the pyramid algorithm
#' @author JJB
#' @keywords internal
#' @examples
#' set.seed(999)
#' x = rnorm(2^8)
#' dwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
dwt_cpp <- function(x, filter_name, nlevels, boundary, brickwall) {
.Call('gmwm_dwt_cpp', PACKAGE = 'gmwm', x, filter_name, nlevels, boundary, brickwall)
}
#' @title Maximum Overlap Discrete Wavelet Transform
#' @description
#' Calculation of the coefficients for the discrete wavelet transformation
#' @inheritParams dwt_cpp
#' @return y A \code{field<vec>} that contains the wavelet coefficients for each decomposition level
#' @keywords internal
#' @details
#' Performs a level J decomposition of the time series using the pyramid algorithm.
#' Use this implementation to supply custom parameters instead of modwt(x),
#' which serves as a wrapper function.
#' @author JJB
#' @keywords internal
#' @examples
#' set.seed(999)
#' x = rnorm(100)
#' modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
modwt_cpp <- function(x, filter_name, nlevels, boundary, brickwall) {
.Call('gmwm_modwt_cpp', PACKAGE = 'gmwm', x, filter_name, nlevels, boundary, brickwall)
}
#' @title Removal of Boundary Wavelet Coefficients
#' @description Removes the first n wavelet coefficients.
#' @param x A \code{field<vec>} that contains the nlevel decomposition using either modwt or dwt.
#' @param wave_filter A \code{field<vec>} containing filter information. Only "haar" is implemented.
#' @param method A \code{string} to describe the mode. Choose between "modwt" and "dwt"
#' @return A \code{field<vec>} with boundary modwt or dwt taken care of.
#' @keywords internal
#' @details
#' The vectors are truncated by removing the first n wavelet coefficients.
#' These vectors are then stored into the field that is returned.
#' Note: As a result, there are no NA's introduced and hence the na.omit is not needed.
#' @examples
#' x = rnorm(100)
#' me = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = FALSE)
#' brick_wall(me, select_filter("haar"), "modwt")
brick_wall <- function(x, wave_filter, method) {
.Call('gmwm_brick_wall', PACKAGE = 'gmwm', x, wave_filter, method)
}
#' @title Generate a white noise process
#' @description Generates a white noise process with variance parameter sigma.
#' @param N An \code{integer} for signal length.
#' @param sigma2 A \code{double} that contains process variance.
#' @return wn A \code{vec} containing the white noise.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_wn(10, 1.5)
gen_wn <- function(N, sigma2 = 1) {
.Call('gmwm_gen_wn', PACKAGE = 'gmwm', N, sigma2)
}
#' @title Generate a drift
#' @description Generates a drift sequence with a given slope.
#' @param N An \code{integer} for signal length.
#' @param slope A \code{double} that contains drift slope
#' @return gd A \code{vec} containing the drift.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_dr(10, 8.2)
gen_dr <- function(N, slope = 5) {
.Call('gmwm_gen_dr', PACKAGE = 'gmwm', N, slope)
}
#' @title Generate a Quantisation Noise (QN) sequence
#' @description Generate an QN sequence given q2
#' @param N An \code{integer} for signal length.
#' @param q2 A \code{double} that contains autocorrection.
#' @return A \code{vec} containing the QN process.
#' @keywords internal
#' @details
#' To generate the quantisation noise, we follow this recipe:
#' First, we generate using a random uniform distribution:
#' \deqn{U_k^*\sim U\left[ {0,1} \right]}{U_k^*~U[0,1]}
#'
#' Then, we multiple the sequence by \eqn{\sqrt{12}}{sqrt(12)} so:
#' \deqn{{U_k} = \sqrt{12} U_k^*}{U_k = sqrt(12)*U_k^*}
#'
#' Next, we find the derivative of \eqn{{U_k}}{U_k}
#' \deqn{{{\dot U}_k} = \frac{{{U_{k + \Delta t}} - {U_k}}}{{\Delta t}}}{U_k^. = (U_(k + (delta)t) - U_k)}
#'
#' In this case, we modify the derivative such that:
#' \eqn{{{\dot U}_k}\Delta t = {U_{k + \Delta t}} - {U_k}}{U_k^. * (delta)t = U_{k + (delta)*t} - U_k}
#'
#' Thus, we end up with:
#' \deqn{{x_k} = \sqrt Q {{\dot U}_k}\Delta t}{x_k = sqrt(Q)*U_k^.*(delta)t}
#' \deqn{{x_k} = \sqrt Q \left( {{U_{k + 1}} - {U_k}} \right)}{x_k = sqrt(Q)* (U_(k+1) - U_(k))}
#'
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @examples
#' gen_qn(10, 5)
gen_qn <- function(N, q2 = .1) {
.Call('gmwm_gen_qn', PACKAGE = 'gmwm', N, q2)
}
#' @title Generate an AR(1) sequence
#' @description Generate an AR sequence given phi and sig2.
#' @details This needs to be extended to AR(p) see \code{arima.sim} and \code{filter}.
#' @param N An \code{integer} for signal length.
#' @param phi A \code{double} that contains autocorrection.
#' @param sigma2 A \code{double} that contains process variance.
#' @return gm A \code{vec} containing the AR(1) process.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_ar1(10, 5, 1.2)
gen_ar1 <- function(N, phi = .3, sigma2 = 1) {
.Call('gmwm_gen_ar1', PACKAGE = 'gmwm', N, phi, sigma2)
}
#' @title Generate a random walk without drift
#' @description Generates a random walk without drift.
#' @param N An \code{integer} for signal length.
#' @param sigma2 A \code{double} that contains process variance.
#' @return grw A \code{vec} containing the random walk without drift.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_rw(10, 8.2)
gen_rw <- function(N, sigma2 = 1) {
.Call('gmwm_gen_rw', PACKAGE = 'gmwm', N, sigma2)
}
#' @title Generate ARMA
#' @description Generate observations for a supplied ARMA model.
#' @param N An \code{integer} for signal length.
#' @param ar A \code{vec} that contains the AR coefficients.
#' @param ma A \code{vec} that contains the MA coefficients.
#' @param sigma2 A \code{double} that contains process variance.
#' @param n_start An \code{unsigned int} that indicates the amount of observations to be used for the burn in period.
#' @details The innovations are generated from a normal distribution.
#' @return A \code{vec} that contains the generated observations.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' gen_arma(100, c(.3,.5), c(.1), 1, 0)
gen_arma <- function(N, ar, ma, sigma2 = 1.5, n_start = 0L) {
.Call('gmwm_gen_arma', PACKAGE = 'gmwm', N, ar, ma, sigma2, n_start)
}
#' @title Generate Time Series based on Model (Internal)
#' @description Create a time series based on a supplied time series model.
#' @param N An \code{interger} containing the amount of observations for the time series.
#' @param theta A \code{vec} containing the parameters to use to generate the model
#' @param desc A \code{vector<string>} containing the different model types (AR1, WN, etc..)
#' @param objdesc A \code{field<vec>} contains the different model objects e.g. AR1 = c(1,1)
#' @return A \code{vec} that contains combined time series.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' # AR
#' set.seed(1336)
#' gen_model(1000, c(.9,1), "AR1", list(c(1,1)))
gen_model <- function(N, theta, desc, objdesc) {
.Call('gmwm_gen_model', PACKAGE = 'gmwm', N, theta, desc, objdesc)
}
#' @title Generate Latent Time Series based on Model (Internal)
#' @description Create a latent time series based on a supplied time series model.
#' @param N An \code{interger} containing the amount of observations for the time series.
#' @param theta A \code{vec} containing the parameters to use to generate the model.
#' @param desc A \code{vector<string>} containing the different model types (AR1, WN, etc..).
#' @param objdesc A \code{field<vec>} containing the different model objects e.g. AR1 = c(1,1)
#' @return A \code{mat} containing data for each decomposed and combined time series.
#' @backref src/gen_process.cpp
#' @backref src/gen_process.h
#' @keywords internal
#' @examples
#' # AR
#' set.seed(1336)
#' gen_lts(1000, c(.9,1), "AR1", list(c(1,1)))
gen_lts <- function(N, theta, desc, objdesc) {
.Call('gmwm_gen_lts', PACKAGE = 'gmwm', N, theta, desc, objdesc)
}
#' @title Optim loses NaN
#' @description This function takes numbers that are very small and sets them to the minimal tolerance for C++.
#' Doing this prevents NaN from entering the optim routine.
#' @param theta A \code{vec} that contains estimated GMWM theta values (untransformed).
#' @return A \code{vec} that contains safe theta values.
#' @keywords internal
code_zero <- function(theta) {
.Call('gmwm_code_zero', PACKAGE = 'gmwm', theta)
}
#' @title Engine for obtaining the GMWM Estimator
#' @description This function uses the Generalized Method of Wavelet Moments (GMWM) to estimate the parameters of a time series model.
#' @param theta A \code{vec} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} containing a list of parameters (e.g. AR(1) = c(1,1), ARMA(p,q) = c(p,q,1))
#' @param model_type A \code{string} that represents the model transformation
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance
#' @param omega A \code{mat} that represents the covariance matrix.
#' @param scales A \code{vec} that contains the scales or taus (2^(1:J))
#' @param starting A \code{bool} that indicates whether we guessed starting (T) or the user supplied estimates (F).
#' @return A \code{vec} that contains the parameter estimates from GMWM estimator.
#' @details
#' If type = "imu" or "ssm", then parameter vector should indicate the characters of the models that compose the latent or state-space model.
#' The model options are:
#' \itemize{
#' \item{"AR1"}{a first order autoregressive process with parameters \eqn{(\phi,\sigma^2)}{phi, sigma^2}}
#' \item{"ARMA"}{an autoregressive moving average process with parameters \eqn{(\phi _p, \theta _q, \sigma^2)}{phi[p], theta[q], sigma^2}}
#' \item{"DR"}{a drift with parameter \eqn{\omega}{omega}}
#' \item{"QN"}{a quantization noise process with parameter \eqn{Q}}
#' \item{"RW"}{a random walk process with parameter \eqn{\sigma^2}{sigma^2}}
#' \item{"WN"}{a white noise process with parameter \eqn{\sigma^2}{sigma^2}}
#' }
#' If model_type = "imu" or type = "ssm" then
#' starting values pass through an initial bootstrap and pseudo-optimization before being passed to the GMWM optimization.
#' If robust = TRUE the function takes the robust estimate of the wavelet variance to be used in the GMWM estimation procedure.
#'
#' @author JJB
#' @references Wavelet variance based estimation for composite stochastic processes, S. Guerrier and Robust Inference for Time Series Models: a Wavelet-Based Framework, S. Guerrier
#' @keywords internal
#' @backref src/gmwm_logic.cpp
#' @backref src/gmwm_logic.h
gmwm_engine <- function(theta, desc, objdesc, model_type, wv_empir, omega, scales, starting) {
.Call('gmwm_gmwm_engine', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, wv_empir, omega, scales, starting)
}
#' @title Update Wrapper for the GMWM Estimator
#' @description This function uses information obtained previously (e.g. WV covariance matrix) to re-estimate a different model parameterization
#' @param theta A \code{vec} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} containing a list of parameters (e.g. AR(1) = c(1,1), ARMA(p,q) = c(p,q,1))
#' @param model_type A \code{string} that represents the model transformation
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance
#' @param omega A \code{mat} that represents the covariance matrix.
#' @param scales A \code{vec} that contains the scales or taus (2^(1:J))
#' @param starting A \code{bool} that indicates whether we guessed starting (T) or the user supplied estimates (F).
#' @return A \code{field<mat>} that contains the parameter estimates from GMWM estimator.
#' @author JJB
#' @references Wavelet variance based estimation for composite stochastic processes, S. Guerrier and Robust Inference for Time Series Models: a Wavelet-Based Framework, S. Guerrier
#' @keywords internal
#' @backref src/gmwm_logic.cpp
#' @backref src/gmwm_logic.h
gmwm_update_cpp <- function(theta, desc, objdesc, model_type, N, expect_diff, ranged, orgV, scales, wv, starting, compute_v, K, H, G, robust, eff) {
.Call('gmwm_gmwm_update_cpp', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, N, expect_diff, ranged, orgV, scales, wv, starting, compute_v, K, H, G, robust, eff)
}
#' @title Master Wrapper for the GMWM Estimator
#' @description This function generates WV, GMWM Estimator, and an initial test estimate.
#' @param data A \code{vec} containing the data.
#' @param theta A \code{vec} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} containing a list of parameters (e.g. AR(1) = c(1,1), ARMA(p,q) = c(p,q,1))
#' @param model_type A \code{string} that represents the model transformation
#' @param starting A \code{bool} that indicates whether the supplied values are guessed (T) or are user-based (F).
#' @param alpha A \code{double} that handles the alpha level of the confidence interval (1-alpha)*100
#' @param compute_v A \code{string} that describes what kind of covariance matrix should be computed.
#' @param K An \code{int} that controls how many times theta is updated.
#' @param H An \code{int} that controls how many bootstrap replications are done.
#' @param G An \code{int} that controls how many guesses at different parameters are made.
#' @param robust A \code{bool} that indicates whether the estimation should be robust or not.
#' @param eff A \code{double} that specifies the amount of efficiency required by the robust estimator.
#' @return A \code{field<mat>} that contains a list of ever-changing estimates...
#' @author JJB
#' @references Wavelet variance based estimation for composite stochastic processes, S. Guerrier and Robust Inference for Time Series Models: a Wavelet-Based Framework, S. Guerrier
#' @keywords internal
#' @backref src/gmwm_logic.cpp
#' @backref src/gmwm_logic.h
gmwm_master_cpp <- function(data, theta, desc, objdesc, model_type, starting, alpha, compute_v, K, H, G, robust, eff) {
.Call('gmwm_gmwm_master_cpp', PACKAGE = 'gmwm', data, theta, desc, objdesc, model_type, starting, alpha, compute_v, K, H, G, robust, eff)
}
#' @title Randomly guess a starting parameter
#' @description Sets starting parameters for each of the given parameters.
#' @param desc A \code{vector<string>} that contains the model's components.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @param model_type A \code{string} that indicates whether it is an SSM or sensor.
#' @param num_param An \code{unsigned int} number of parameters in the model (e.g. # of thetas).
#' @param expect_diff A \code{double} that contains the mean of the first difference of the data
#' @param N A \code{integer} that contains the number of observations in the data.
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance.
#' @param tau A \code{vec} that contains the scales. (e.g. 2^(1:J))
#' @param double A \code{double} that contains the drift slope given by \eqn{\frac{max-min}{N}}{(Max-Min)/N}
#' @param G A \code{integer} that indicates how many random draws that should be performed.
#' @return A \code{vec} containing smart parameter starting guesses to be iterated over.
#' @keywords internal
#' @examples
#' #TBA
guess_initial <- function(desc, objdesc, model_type, num_param, expect_diff, N, wv, tau, ranged, G) {
.Call('gmwm_guess_initial', PACKAGE = 'gmwm', desc, objdesc, model_type, num_param, expect_diff, N, wv, tau, ranged, G)
}
#' @title Randomly guess starting parameters for AR1
#' @description Sets starting parameters for each of the given parameters.
#' @param draw_id An \code{unsigned int} that contains the draw principles.
#' @param last_phi A \code{double} containing the last guessed phi value.
#' @param sigma2_total A \code{double} that contains the sum of all WVs.
#' @param model_type A \code{string} that describes the model transformation.
#' @return A \code{vec} containing smart parameter starting guesses to be iterated over.
#' @keywords internal
#' @examples
#' #TBA
ar1_draw <- function(draw_id, last_phi, sigma2_total, model_type) {
.Call('gmwm_ar1_draw', PACKAGE = 'gmwm', draw_id, last_phi, sigma2_total, model_type)
}
#' @title Randomly guess starting parameters for ARMA
#' @description Sets starting parameters for each of the given parameters.
#' @param p An \code{unsigned int} that contains the amount of AR parameters to generate.
#' @param q An \code{unsigned int} that contains the amount of MA parameters to generate.
#' @param sigma2_total A \code{double} that contains the sum of all WVs.
#' @return A \code{vec} containing smart parameter starting guesses to be iterated over.
#' @keywords internal
#' @examples
#' #TBA
arma_draws <- function(p, q, sigma2_total) {
.Call('gmwm_arma_draws', PACKAGE = 'gmwm', p, q, sigma2_total)
}
#' @title Randomly guess a starting parameter
#' @description Sets starting parameters for each of the given parameters.
#' @param desc A \code{vector<string>} that contains the model's components.
#' @param objdesc A \code{field<vec>} that contains an object description (e.g. values) of the model.
#' @param model_type A \code{string} that indicates whether it is an SSM or sensor.
#' @param num_param An \code{unsigned int} number of parameters in the model (e.g. # of thetas).
#' @param expect_diff A \code{double} that contains the mean of the first difference of the data
#' @param N A \code{integer} that contains the number of observations in the data.
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance.
#' @param tau A \code{vec} that contains the scales. (e.g. 2^(1:J))
#' @param B A \code{integer} that indicates how many random draws that should be performed.
#' @return A \code{vec} containing smart parameter starting guesses to be iterated over.
#' @keywords internal
#' @examples
#' #TBA
guess_initial_old <- function(desc, objdesc, model_type, num_param, expect_diff, N, wv_empir, tau, B) {
.Call('gmwm_guess_initial_old', PACKAGE = 'gmwm', desc, objdesc, model_type, num_param, expect_diff, N, wv_empir, tau, B)
}
#' @title Indirect Inference for ARMA
#' @description Option for indirect inference
#' @param ar A \code{vec} that contains the coefficients of the AR process.
#' @param ma A \code{vec} that contains the coefficients of the MA process.
#' @param sigma2 A \code{double} that indicates the sigma2 parameter of the ARMA process.
#' @param N A \code{int} that indicates how long the time series is.
#' @param robust A \code{bool} that indicates whether the estimation should be robust or not.
#' @param eff A \code{double} that specifies the amount of efficiency required by the robust estimator.
#' @param H A \code{int} that indicates how many iterations should take place.
#' @return A \code{vec} with the indirect inference results.
#' @keywords internal
idf_arma <- function(ar, ma, sigma2, N, robust, eff, H) {
.Call('gmwm_idf_arma', PACKAGE = 'gmwm', ar, ma, sigma2, N, robust, eff, H)
}
#' @title Indirect Inference for ARMA
#' @description Option for indirect inference
#' @param ar A \code{vec} that contains the coefficients of the AR process.
#' @param ma A \code{vec} that contains the coefficients of the MA process.
#' @param sigma2 A \code{double} that indicates the sigma2 parameter of the ARMA process.
#' @param N A \code{int} that indicates how long the time series is.
#' @param robust A \code{bool} that indicates whether the estimation should be robust or not.
#' @param eff A \code{double} that specifies the amount of efficiency required by the robust estimator.
#' @param H A \code{int} that indicates how many iterations should take place.
#' @return A \code{vec} with the indirect inference results.
#' @keywords internal
idf_arma_total <- function(ar, ma, sigma2, N, robust, eff, H) {
.Call('gmwm_idf_arma_total', PACKAGE = 'gmwm', ar, ma, sigma2, N, robust, eff, H)
}
#' @title Calculate the Psi matrix
#' @description Computes the Psi matrix using supplied parameters
#' @param A first derivative matrix
#' @param v_hat bootstrapped V
#' @param omega original omega matrix
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @return A \code{mat} that has the first column
#' @keywords internal
calculate_psi_matrix <- function(A, v_hat, omega) {
.Call('gmwm_calculate_psi_matrix', PACKAGE = 'gmwm', A, v_hat, omega)
}
#' @title Format the Confidence Interval for Estimates
#' @description Creates hi and lo confidence based on SE and alpha.
#' @param theta A \code{vec} containing the estimates
#' @param se A \code{vec} containing the standard error
#' @param alpha A \code{double} that contains the confidence level.
#' @return A \code{mat} that has:
#' \itemize{
#' \item Column 1: Lo CI
#' \item Column 2: Hi CI
#' \item Column 3: SE
#' }
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @keywords internal
format_ci <- function(theta, se, alpha) {
.Call('gmwm_format_ci', PACKAGE = 'gmwm', theta, se, alpha)
}
#' @title Generate the Confidence Interval for Theta Estimates
#' @description Create an Asymptotic CI for the Theta Estimates.
#' @param theta A \code{vec} containing the estimates
#' @param A A \code{mat} that is the first derivative matrix.
#' @param v_hat A \code{mat} that is the bootstrapped V matrix
#' @param omega A \code{mat} that is the inverse of the diagonal V matrix.
#' @param alpha A \code{double} that contains the confidence level.
#' @return A \code{mat} that has the first column
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @keywords internal
theta_ci <- function(theta, A, v_hat, omega, alpha) {
.Call('gmwm_theta_ci', PACKAGE = 'gmwm', theta, A, v_hat, omega, alpha)
}
#' @title Compute the GOF Test
#' @description yaya
#' @template tsobj_cpp
#' @param model_type A \code{string} that contains the model type: \code{"imu"} or \code{"ssm"}
#' @param tau A \code{vec} containing the scales of a proccess.
#' @param v_hat A \code{mat} that contains the bootstrapped matrix.
#' @param wv_empir A \code{vec} that contains the empirical wavelet variance.
#' @return A \code{vec} that has
#' \itemize{
#' \item Test Statistic
#' \item P-Value
#' \item DF
#' }
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @keywords internal
gof_test <- function(theta, desc, objdesc, model_type, tau, v_hat, wv_empir) {
.Call('gmwm_gof_test', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, tau, v_hat, wv_empir)
}
#' @title Compute the Bootstrapped GoF Test
#' @description Handles the bootstrap computation and the bootstrapped p-value.
#' @param obj_value A \code{double} that contains the optimized objective function value.
#' @param bs_obj_values A \code{vec} that contains the objective function values under bootstrap.
#' @param alpha A \code{double} that indicates the confidence.
#' @param bs_gof_p_ci A \code{bool} that indicates whether CIs should be included or not.
#' @return A \code{vec} that has
#' \itemize{
#' \item Test Statistic
#' \item Low CI
#' \item Upper CI - BS
#' }
#' @backref src/inference.cpp
#' @backref src/inference.h
#' @keywords internal
bootstrap_gof_test <- function(obj_value, bs_obj_values, alpha, bs_gof_p_ci) {
.Call('gmwm_bootstrap_gof_test', PACKAGE = 'gmwm', obj_value, bs_obj_values, alpha, bs_gof_p_ci)
}
#' @title B Matrix
#' @description B Matrix
#' @param A A \code{mat} containing the first derivatives of the process.
#' @param at_omega A \code{mat} containing A^T * Omega
#' @return A \code{mat}
#' @backref src/model_selection.cpp
#' @backref src/model_selection.h
#' @keywords internal
B_matrix <- function(A, at_omega) {
.Call('gmwm_B_matrix', PACKAGE = 'gmwm', A, at_omega)
}
#' @title Model Score
#' @description Calculates the modeling score of a GMWM
#' @param A A \code{mat} that contains the first derivatives of the processes
#' @param At_j A \code{mat} that contains the second derivative of each process
#' @param omega A \code{mat} that contains the omega used when calculating the GMWM
#' @param v_hat A \code{mat} that contains the covariance matrix
#' @param diff A \code{vec} that is the difference of the WV empirical and WV theoretical
#' @return A \code{vec}
#' @keywords internal
#' @backref src/model_selection.cpp
#' @backref src/model_selection.h
#' @details
#' The equation is slightly different than that stated in the paper due to the bootstrap already incorporating in
#' N.
model_score <- function(A, D, omega, v_hat, obj_value) {
.Call('gmwm_model_score', PACKAGE = 'gmwm', A, D, omega, v_hat, obj_value)
}
#' @title Extract Object
#' @description Extracts the object information and returns it.
#' @param theta A \code{vec} containing the theta values.
#' @param objdesc A \code{vec} at the desc point.
#' @param cur_position An \code{integer} at the current position.
#' @return A \code{field<vec>} containing the breakdown of the object.
#' @keywords internal
obj_extract <- function(theta, objdesc, cur_position) {
.Call('gmwm_obj_extract', PACKAGE = 'gmwm', theta, objdesc, cur_position)
}
#' @title Retrieve GMWM starting value from Yannick's objective function
#' @description Obtains the GMWM starting value given by Yannick's objective function optimization
#' @template tsobj_cpp
#' @param model_type A \code{string} containing the model type. Either 'imu' or 'ssm'
#' @param wv_empir A \code{vec} containing the empirical wavelet variance.
#' @param tau A \code{vec} that contains the scales of 2^(1:J), where J is the number of scales created by the decomposition.
#' @return A \code{double} that is the value of the Objective function under Yannick's starting algorithm
#' @keywords internal
getObjFunStarting <- function(theta, desc, objdesc, model_type, wv_empir, tau) {
.Call('gmwm_getObjFunStarting', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, wv_empir, tau)
}
#' @title Retrieve GMWM starting value from Yannick's objective function
#' @description Obtains the GMWM starting value given by Yannick's objective function optimization
#' @template tsobj_cpp
#' @param model_type A \code{string} containing the model type. Either 'imu' or 'ssm'
#' @param omega A \code{mat} that is the inverse of the diagonal of the V matrix.
#' @param wv_empir A \code{vec} containing the empirical wavelet variance.
#' @param tau A \code{vec} that contains the scales of 2^(1:J), where J is the number of scales created by the decomposition.
#' @return A \code{double} that is the value of the Objective function under Yannick's starting algorithm
#' @keywords internal
getObjFun <- function(theta, desc, objdesc, model_type, omega, wv_empir, tau) {
.Call('gmwm_getObjFun', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, omega, wv_empir, tau)
}
#' @title Root Finding C++
#' @description Used to interface with Armadillo
#' @param z A \code{cx_vec} (complex vector) that has 1 in the beginning (e.g. c(1,3i,-3i))
#' @keywords internal
do_polyroot_arma <- function(z) {
.Call('gmwm_do_polyroot_arma', PACKAGE = 'gmwm', z)
}
#' @title Root Finding C++
#' @description Vroom Vroom
#' @param z A \code{vec<complex<double>} (complex vector) that has 1 in the beginning (e.g. c(1,3i,-3i))
#' @keywords internal
do_polyroot_cpp <- function(z) {
.Call('gmwm_do_polyroot_cpp', PACKAGE = 'gmwm', z)
}
#' @title ARMA process to WV
#' @description This function computes the (haar) WV of an ARMA process
#' @param ar A \code{vec} containing the coefficients of the AR process
#' @param ma A \code{vec} containing the coefficients of the MA process
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @param sigma A \code{double} containing the residual variance
#' @return A \code{vec} containing the wavelet variance of the ARMA process.
#' @examples
#' arma_to_wv(c(.23,.43), c(.34,.41,.59), 2^(1:9), 3)
#' @seealso \code{\link{ARMAtoMA_cpp}},\code{\link{ARMAacf_cpp}}
arma_to_wv <- function(ar, ma, tau, sigma) {
.Call('gmwm_arma_to_wv', PACKAGE = 'gmwm', ar, ma, tau, sigma)
}
#' @title Helper Function for ARMA to WV Approximation
#' @description Indicates where the minimum ARMAacf value is and returns that as an index.
#' @param ar A \code{vec} containing the coefficients of the AR process
#' @param ma A \code{vec} containing the coefficients of the MA process
#' @param last_tau An \code{int} the Jth scale of 2^(1:J)
#' @param alpha A \code{double} indicating the cutoff.
#' @return A \code{vec} containing the wavelet variance of the ARMA process.
#' @keywords internal
#' @seealso \code{\link{arma_to_wv_app}}
acf_sum <- function(ar, ma, last_tau, alpha = 0.99) {
.Call('gmwm_acf_sum', PACKAGE = 'gmwm', ar, ma, last_tau, alpha)
}
#' @title ARMA process to WV approximation
#' @description This function computes the (haar) WV of an ARMA process
#' @param ar A \code{vec} containing the coefficients of the AR process
#' @param ma A \code{vec} containing the coefficients of the MA process
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @param sigma A \code{double} containing the residual variance
#' @param alpha A \code{double} indicating the cutoff.
#' @return A \code{vec} containing the wavelet variance of the ARMA process.
#' @keywords internal
#' @examples
#' arma_to_wv_app(c(.23,.43), c(.34,.41,.59), 2^(1:9), 3, .9)
#' @seealso \code{\link{ARMAtoMA_cpp}},\code{\link{ARMAacf_cpp}}
arma_to_wv_app <- function(ar, ma, tau, sigma, alpha = 0.9999) {
.Call('gmwm_arma_to_wv_app', PACKAGE = 'gmwm', ar, ma, tau, sigma, alpha)
}
#' @title Quantisation Noise to WV
#' @description This function compute the WV (haar) of a Quantisation Noise (QN) process
#' @param q2 A \code{double} corresponding to variance of drift
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the QN.
#' @examples
#' x.sim = 1:1000
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = qn_to_wv(1, tau)
#' plot(tau, wv.theo, col = "red")
qn_to_wv <- function(q2, tau) {
.Call('gmwm_qn_to_wv', PACKAGE = 'gmwm', q2, tau)
}
#' @title White Noise to WV
#' @description This function compute the WV (haar) of a White Noise process
#' @param sig2 A \code{double} corresponding to variance of WN
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the white noise.
#' @examples
#' x.sim = cumsum(rnorm(100000))
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = wn_to_wv(1, tau)
#' plot(tau, wv.theo, col = "red")
wn_to_wv <- function(sig2, tau) {
.Call('gmwm_wn_to_wv', PACKAGE = 'gmwm', sig2, tau)
}
#' @title Random Walk to WV
#' @description This function compute the WV (haar) of a Random Walk process
#' @param sig2 A \code{double} corresponding to variance of RW
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the random walk.
#' @examples
#' x.sim = cumsum(rnorm(100000))
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = rw_to_wv(1,tau)
#' plot(tau, wv.theo, col = "red")
rw_to_wv <- function(sig2, tau) {
.Call('gmwm_rw_to_wv', PACKAGE = 'gmwm', sig2, tau)
}
#' @title Drift to WV
#' @description This function compute the WV (haar) of a Drift process
#' @param omega A \code{double} corresponding to variance of drift
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the drift.
#' @examples
#' x.sim = 1:1000
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = dr_to_wv(1, tau)
#' plot(tau, wv.theo, col = "red")
dr_to_wv <- function(omega, tau) {
.Call('gmwm_dr_to_wv', PACKAGE = 'gmwm', omega, tau)
}
#' @title AR1 process to WV
#' @description This function compute the WV (haar) of an AR(1) process
#' @param phi A \code{double} that is the phi term of the AR(1) process
#' @param sig2 A \code{double} corresponding to variance of AR(1) process
#' @param tau A \code{vec} containing the scales e.g. 2^tau
#' @return A \code{vec} containing the wavelet variance of the AR(1) process.
#' @examples
#' x.sim = gen_ar1( N = 10000, phi = 0.9, sigma2 = 4 )
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = ar1_to_wv(phi = 0.9, sig2 = 16, tau)
#' plot(tau, wv.theo, col = "red")
ar1_to_wv <- function(phi, sig2, tau) {
.Call('gmwm_ar1_to_wv', PACKAGE = 'gmwm', phi, sig2, tau)
}
#' @title Model Process to WV
#' @description This function computes the summation of all Processes to WV (haar) in a given model
#' @param theta A \code{vec} containing the list of estimated parameters.
#' @param desc A \code{vector<string>} containing a list of descriptors.
#' @param objdesc A \code{field<vec>} containing a list of object descriptors.
#' @param tau A \code{vec} containing the scales e.g. 2^(1:J)
#' @return A \code{vec} containing the wavelet variance of the model.
#' @examples
#' x.sim = gen_ar1( N = 10000, phi = 0.9, sigma2 = 4 )
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = ar1_to_wv(phi = 0.9, sig2 = 16, tau)
#' plot(tau, wv.theo, col = "red")
#' @keywords internal
theoretical_wv <- function(theta, desc, objdesc, tau) {
.Call('gmwm_theoretical_wv', PACKAGE = 'gmwm', theta, desc, objdesc, tau)
}
#' @title Each Models Process Decomposed to WV
#' @description This function computes each process to WV (haar) in a given model.
#' @param theta A \code{vec} containing the list of estimated parameters.
#' @param desc A \code{vector<string>} containing a list of descriptors.
#' @param objdesc A \code{field<vec>} containing a list of object descriptors.
#' @param tau A \code{vec} containing the scales e.g. 2^(1:J)
#' @return A \code{mat} containing the wavelet variance of each process in the model
#' @examples
#' x.sim = gen_ar1( N = 10000, phi = 0.9, sigma2 = 4 )
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = ar1_to_wv(phi = 0.9, sig2 = 16, tau)
#' plot(tau, wv.theo, col = "red")
#' @keywords internal
decomp_theoretical_wv <- function(theta, desc, objdesc, tau) {
.Call('gmwm_decomp_theoretical_wv', PACKAGE = 'gmwm', theta, desc, objdesc, tau)
}
#' @title Decomposed WV to Single WV
#' @description This function computes the combined processes to WV (haar) in a given model.
#' @param decomp A \code{mat} with scales as rows and processes as columns
#' @return A \code{vec} containing the wavelet variance of the process for the overall model
#' @examples
#' x.sim = gen_ar1( N = 10000, phi = 0.9, sigma2 = 4 )
#' ntau = floor(log(length(x.sim),2))
#' tau = 2^(1:ntau)
#' wv.theo = ar1_to_wv(phi = 0.9, sig2 = 16, tau)
#' plot(tau, wv.theo, col = "red")
#' @keywords internal
decomp_to_theo_wv <- function(decomp) {
.Call('gmwm_decomp_to_theo_wv', PACKAGE = 'gmwm', decomp)
}
#' @title Read an IMU Binary File into R
#'
#' @description
#' The function will take a file location in addition to the type of sensor it
#' came from and read the data into R.
#'
#' @param file_path A \code{string} that contains the full file path.
#' @param imu_type A \code{string} that contains a supported IMU type given below.
#' @details
#' Currently supports the following IMUs:
#' \itemize{
#' \item IMAR
#' \item LN200
#' \item LN200IG
#' \item IXSEA
#' \item NAVCHIP_INT
#' \item NAVCHIP_FLT
#' }
#'
#' We hope to soon be able to support delimited files.
#' @return A matrix with dimensions N x 7, where the columns represent:
#' \describe{
#' \item{Col 0}{Time}
#' \item{Col 1}{Gyro 1}
#' \item{Col 2}{Gyro 2}
#' \item{Col 3}{Gyro 3}
#' \item{Col 4}{Accel 1}
#' \item{Col 5}{Accel 2}
#' \item{Col 6}{Accel 3}
#' }
#' @references
#' Thanks goes to Philipp Clausen of Labo TOPO, EPFL, Switzerland, topo.epfl.ch, Tel:+41(0)21 693 27 55
#' for providing a matlab function that reads in IMUs.
#' The function below is a heavily modified port of MATLAB code into Armadillo/C++.
#'
#' @examples
#' \dontrun{
#' read_imu(file_path = "F:/Desktop/short_test_data.imu", imu_type = "IXSEA")
#' }
#' @keywords internal
read_imu <- function(file_path, imu_type) {
.Call('gmwm_read_imu', PACKAGE = 'gmwm', file_path, imu_type)
}
#' @title Generate a sequence of values
#' @description Creates a vector containing a sequence of values starting at the initial point and going to the terminal point.
#' @param a An \code{int}, that denotes the starting point.
#' @param b An \code{int}, that denotes the ending point.
#' @return A \code{vector} containing values moving from a to b. There are no restrictions on A's range.
#' @author James J Balamuta
#' @keywords internal
#' @examples
#' #Call with the following data:
#' seq_cpp(3, 5)
#' seq_cpp(5, 3)
seq_cpp <- function(a, b) {
.Call('gmwm_seq_cpp', PACKAGE = 'gmwm', a, b)
}
#' @title Generate a sequence of values based on supplied number
#' @description Creates a vector containing a sequence of values starting at 1 and going to the terminal point.
#' @param n An \code{int} that denotes the length of the vector.
#' @return A \code{vector} containing values moving from 1 to n.
#' @author James J Balamuta
#' @keywords internal
#' @examples
#' #Call with the following data:
#' seq_len_cpp(5)
seq_len_cpp <- function(n) {
.Call('gmwm_seq_len_cpp', PACKAGE = 'gmwm', n)
}
#' @title Find Quantiles
#' @description Attempts to find quantiles
#' @param x A \code{vec} of data
#' @param probs A \code{vec} of the quantiles to find.
#' @return A \code{vector} containing the quantiles
#' @author James J Balamuta
#' @keywords internal
#' @examples
#' #Call with the following data:
#' quantile_cpp(c(1,2,3,4,5,6,7), c(.25,.5,.75))
#' quantile(c(1,2,3,4,5,6,7), c(.25,.5,.75))
quantile_cpp <- function(x, probs) {
.Call('gmwm_quantile_cpp', PACKAGE = 'gmwm', x, probs)
}
#' @title Lagged Differences in Armadillo
#' @description Returns the ith difference of a time series of rth lag.
#' @param x A \code{vec} that is the time series
#' @param lag A \code{unsigned int} that indicates the lag
#' @param differences A \code{dif} that indicates how many differences should be taken
#' @return A \code{vector} containing the differenced time series.
#' @author JJB
#' @keywords internal
#' @examples
#' x = rnorm(10000, 0, 1)
#' diff_cpp(x,1,1)
diff_cpp <- function(x, lag, differences) {
.Call('gmwm_diff_cpp', PACKAGE = 'gmwm', x, lag, differences)
}
#' @title Converting an ARMA Process to an Infinite MA Process
#' @description Takes an ARMA function and converts it to an infinite MA process.
#' @param ar A \code{column vector} of length p
#' @param ma A \code{column vector} of length q
#' @param lag_max A \code{int} of the largest MA(Inf) coefficient required.
#' @return A \code{column vector} containing coefficients
#' @details This function is a port of the base stats package's ARMAtoMA. There is no significant speed difference between the two.
#' @author R Core Team and JJB
#' @keywords internal
#' @examples
#' # ARMA(2,1)
#' ARMAtoMA_cpp(c(1.0, -0.25), 1.0, 10)
#' # ARMA(0,1)
#' ARMAtoMA_cpp(numeric(0), 1.0, 10)
ARMAtoMA_cpp <- function(ar, ma, lag_max) {
.Call('gmwm_ARMAtoMA_cpp', PACKAGE = 'gmwm', ar, ma, lag_max)
}
#' @title Time Series Convolution Filters
#' @description Applies a convolution filter to a univariate time series.
#' @param x A \code{column vector} of length T
#' @param filter A \code{column vector} of length f
#' @param sides An \code{int} that takes either 1:for using past values only or 2: filter coefficients are centered around lag 0.
#' @param circular A \code{bool} that indicates if the filter should be wrapped around the ends of the time series.
#' @return A \code{column vec} that contains the results of the filtering process.
#' @details This is a port of the cfilter function harnessed by the filter function in stats.
#' It is about 5-7 times faster than R's base function. The benchmark was done on iMac Late 2013 using vecLib as the BLAS.
#' @author R Core Team and JJB
#' @keywords internal
#' @examples
#' x = 1:100
#' #
#' cfilter(x, rep(1, 3), sides = 2, circular = FALSE)
#' # Using R's function
#' filter(x, rep(1, 3))
#' #
#' cfilter(x, rep(1, 3), sides = 1, circular = FALSE)
#' # Using R's function
#' filter(x, rep(1, 3), sides = 1)
#' #
#' cfilter(x, rep(1, 3), sides = 1, circular = TRUE)
#' # Using R's function
#' filter(x, rep(1, 3), sides = 1, circular = TRUE)
cfilter <- function(x, filter, sides, circular) {
.Call('gmwm_cfilter', PACKAGE = 'gmwm', x, filter, sides, circular)
}
#' @title Time Series Recursive Filters
#' @description Applies a recursive filter to a univariate time series.
#' @usage rfilter(x, filter, init)
#' @param x A \code{column vector} of length T
#' @param filter A \code{column vector} of length f
#' @param init A \code{column vector} of length f that contains the initial values of the time series in reverse.
#' @return x A \code{column vector} with its contents reversed.
#' @details Note: The length of 'init' must be equal to the length of 'filter'.
#' This is a port of the rfilter function harnessed by the filter function in stats.
#' It is about 6-7 times faster than R's base function. The benchmark was done on iMac Late 2013 using vecLib as the BLAS.
#' @author R Core Team and JJB
#' @keywords internal
#' @examples
#' x = 1:100
#' #
#' rfilter(x, rep(1, 3), rep(1, 3))
#' # Using R's function
#' filter(x, rep(1, 3), method="recursive", init=rep(1, 3))
rfilter <- function(x, filter, init) {
.Call('gmwm_rfilter', PACKAGE = 'gmwm', x, filter, init)
}
#' @title Compute Theoretical ACF for an ARMA Process
#' @description Compute the theoretical autocorrelation function for an ARMA process.
#' @usage ARMAacf_cpp(ar,ma,lag_max)
#' @param ar A \code{vector} of length p containing AR coefficients
#' @param ma A \code{vector} of length q containing MA coefficients
#' @param lag_max A \code{unsigned integer} indicating the maximum lag necessary
#' @return x A \code{matrix} listing values from 1...nx in one column and 1...1, 2...2,....,n...n, in the other
#' @details This is an implementaiton of the ARMAacf function in R. It is approximately 40x times faster. The benchmark was done on iMac Late 2013 using vecLib as the BLAS.
#' @author R Core Team and JJB
#' @keywords internal
#' @examples
#' # ARMA(2,1)
#' ARMAacf_cpp(c(1.0, -0.25), 1.0, lag_max = 10)
#' # ARMA(0,1)
#' ARMAacf_cpp(numeric(0), .35, lag_max = 10)
ARMAacf_cpp <- function(ar, ma, lag_max) {
.Call('gmwm_ARMAacf_cpp', PACKAGE = 'gmwm', ar, ma, lag_max)
}
#' @title Discrete Fourier Transformation for Autocovariance Function
#' @description Calculates the autovariance function (ACF) using Discrete Fourier Transformation.
#' @param x A \code{cx_vec}.
#' @return A \code{vec} containing the ACF.
#' @details
#' This implementation is 2x as slow as Rs.
#' Two issues: 1. memory resize and 2. unoptimized fft algorithm in arma.
#' Consider piping back into R and rewrapping the object. (Decrease of about 10 microseconds.)
#' @keywords internal
#' @examples
#' x=rnorm(100)
#' dft_acf(x)
dft_acf <- function(x) {
.Call('gmwm_dft_acf', PACKAGE = 'gmwm', x)
}
#' @title Mean of the First Difference of the Data
#' @description The mean of the first difference of the data
#' @param x A \code{vec} containing the data
#' @return A \code{double} that contains the mean of the first difference of the data.
#' @keywords internal
#' @examples
#' x=rnorm(100)
#' mean_diff(x)
mean_diff <- function(x) {
.Call('gmwm_mean_diff', PACKAGE = 'gmwm', x)
}
#' @title Routing function for summary info
#' @description Gets all the data for the summary.gmwm function.
#' @param theta A \code{vec} with dimensions N x 1 that contains user-supplied initial values for parameters
#' @param desc A \code{vector<string>} indicating the models that should be considered.
#' @param objdesc A \code{field<vec>} containing a list of parameters (e.g. AR(1) = c(1,1), ARMA(p,q) = c(p,q,1))
#' @param model_type A \code{string} that represents the model transformation
#' @param wv_empir A \code{vec} that
#' @param theo A \code{vec} that
#' @param scales A \code{vec} that
#' @param V A \code{mat} that contains the V matrix used to obtain the GMWM.
#' @param omega A \code{mat} that
#' @param obj_value A \code{double} that contains the objective function value at the optimized solution.
#' @param N A \code{int} that indicates how long the time series is.
#' @param alpha A \code{double} that handles the alpha level of the confidence interval (1-alpha)*100
#' @param robust A \code{bool} that indicates whether the estimation should be robust or not.
#' @param eff A \code{double} that specifies the amount of efficiency required by the robust estimator.
#' @param inference A \code{bool} that indicates whether inference (e.g. GoF) should be run.
#' @param fullV A \code{bool} that indicates whether the matrix has been fully bootstrapped.
#' @param bs_gof A \code{bool} indicating whether the GoF should be bootstrapped or done asymptotically.
#' @param bs_gof_p_ci A \code{bool} indicating whether a bootstrapped p-value should be generated during the bootstrapped GoF
#' @param bs_ci A \code{bool} that indicates whether a bootstrapped CI should be obtained or to use analytical derivatives.
#' @param B A \code{int} that indicates how many iterations should take place.
#' @return A \code{field<mat>} that contains bootstrapped / asymptotic GoF results as well as CIs.
#' @keywords internal
get_summary <- function(theta, desc, objdesc, model_type, wv_empir, theo, scales, V, omega, obj_value, N, alpha, robust, eff, inference, fullV, bs_gof, bs_gof_p_ci, bs_theta_est, bs_ci, B) {
.Call('gmwm_get_summary', PACKAGE = 'gmwm', theta, desc, objdesc, model_type, wv_empir, theo, scales, V, omega, obj_value, N, alpha, robust, eff, inference, fullV, bs_gof, bs_gof_p_ci, bs_theta_est, bs_ci, B)
}
#' @title Pseudo Logit Inverse Function
#' @description This function computes the pseudo inverse of a logit transformation of the parameters in order to constrain them to a positive domain
#' @param x A \code{vec} containing real numbers.
#' @return A \code{vec} containing logit probabilities.
#' @keywords internal
#' @examples
#' x.sim = rnorm(100)
#' pseudo_logit_inv(x.sim)
pseudo_logit_inv <- function(x) {
.Call('gmwm_pseudo_logit_inv', PACKAGE = 'gmwm', x)
}
#' @title Logit Inverse Function
#' @description This function computes the inverse of a logit transformation of the parameters.
#' @param x A \code{vec} containing real numbers.
#' @return A \code{vec} containing logit probabilities.
#' @keywords internal
#' @examples
#' x.sim = rnorm(100)
#' logit_inv(x.sim)
logit_inv <- function(x) {
.Call('gmwm_logit_inv', PACKAGE = 'gmwm', x)
}
#' @title Pseudo Logit Function
#' @description This function compute the link function to constrain parameters to a positive domain.
#' @param x A \code{vec} containing probabilities (e.g. 0 <= x <= 1)
#' @return A \code{vec} containing logit terms.
#' @keywords internal
#' @examples
#' x.sim = runif(100)
#' pseudo_logit(x.sim)
pseudo_logit <- function(x) {
.Call('gmwm_pseudo_logit', PACKAGE = 'gmwm', x)
}
#' @title Logit Function
#' @description This function computes the logit link function.
#' @param x A \code{vec} containing probabilities (e.g. -1 <= x <= 1)
#' @return A \code{vec} containing logit terms.
#' @keywords internal
#' @examples
#' x.sim = runif(100)
#' logit(x.sim)
logit <- function(x) {
.Call('gmwm_logit', PACKAGE = 'gmwm', x)
}
#' @title Logit Function
#' @description This function computes the logit link function.
#' @param x A \code{vec} containing probabilities (e.g. -1 <= x <= 1)
#' @return A \code{vec} containing logit terms.
#' @keywords internal
#' @examples
#' x.sim = runif(100)
#' logit(x.sim)
logit2 <- function(x) {
.Call('gmwm_logit2', PACKAGE = 'gmwm', x)
}
#' @title Logit2 Inverse Function
#' @description This function computes the inverse of a logit transformation of the parameters.
#' @param x A \code{vec} containing real numbers.
#' @return A \code{vec} containing logit probabilities.
#' @keywords internal
#' @examples
#' x.sim = rnorm(100)
#' logit_inv(x.sim)
logit2_inv <- function(x) {
.Call('gmwm_logit2_inv', PACKAGE = 'gmwm', x)
}
#' @title Transform Values for Optimization
#' @description Transform parameter guesses prior to estimating with GMWM
#' @template tsobj_cpp
#' @param model_type A \code{string} that contains the model type: \code{"imu"} or \code{"ssm"}
#' @return A \code{vec} containing the transformed guesses.
#' @keywords internal
transform_values <- function(theta, desc, objdesc, model_type) {
.Call('gmwm_transform_values', PACKAGE = 'gmwm', theta, desc, objdesc, model_type)
}
#' @title Revert Transform Values for Display
#' @description Undo the previous transform of parameter guesses to obtain the GMWM estimates.
#' @template tsobj_cpp
#' @param model_type A \code{string} that contains the model type: \code{"imu"} or \code{"ssm"}
#' @return A \code{vec} containing the undone transformation of parameters.
#' @keywords internal
untransform_values <- function(theta, desc, objdesc, model_type) {
.Call('gmwm_untransform_values', PACKAGE = 'gmwm', theta, desc, objdesc, model_type)
}
#' @title Obtain the smallest polynomial root
#' @description Calculates all the roots of a polynomial and returns the root that is the smallest.
#' @param x A \code{cx_vec} that has a 1 appended before the coefficents. (e.g. c(1, x))
#' @return A \code{double} with the minimum root value.
#' @keywords internal
minroot <- function(x) {
.Call('gmwm_minroot', PACKAGE = 'gmwm', x)
}
#' @title Check Invertibility Conditions
#' @description Checks the invertiveness of series of coefficients.
#' @param x A \code{cx_vec} that has a 1 appended before the coefficents. (e.g. c(1, x))
#' @return True (if outside unit circle) || False (if inside unit circle)
#' @keywords internal
invert_check <- function(x) {
.Call('gmwm_invert_check', PACKAGE = 'gmwm', x)
}
#' @title Count Models
#' @description Count the amount of models that exist.
#' @param desc A \code{vector<string>} that contains the model's components.
#' @return A \code{map<string, int>} containing how frequent the model component appears.
#' @keywords internal
count_models <- function(desc) {
.Call('gmwm_count_models', PACKAGE = 'gmwm', desc)
}
#' @title Order AR1s by size of phi.
#' @description Changes the order of AR1s in a string by size.
#' @template tsobj_cpp
#' @return A \code{vec} that has AR1s shown in descending parameter value.
#' @keywords internal
order_AR1s <- function(theta, desc, objdesc) {
.Call('gmwm_order_AR1s', PACKAGE = 'gmwm', theta, desc, objdesc)
}
#' @title Transform AR1 to GM
#' @description
#' Takes AR1 values and transforms them to GM
#' @param theta A \code{vec} that contains AR1 values.
#' @param freq A \code{double} indicating the frequency of the data.
#' @return A \code{vec} containing GM values.
#' @details
#' The function takes a vector of AR1 values \eqn{\phi}{phi} and \eqn{\sigma ^2}{sigma ^2}
#' and transforms them to GM values \eqn{\beta}{beta} and \eqn{\sigma ^2_{gm}}{sigma ^2[gm]}
#' using the formulas:
#' \eqn{\beta = - \frac{{\ln \left( \phi \right)}}{{\Delta t}}}{beta = -ln(phi)/delta_t}
#' \eqn{\sigma _{gm}^2 = \frac{{{\sigma ^2}}}{{1 - {\phi ^2}}} }{sigma^2[gm] = sigma^2/(1-phi^2)}
#' @keywords internal
#' @author JJB
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
#' @examples
#' ar1_to_gm(c(0.3,1,0.6,.3), 2)
ar1_to_gm <- function(theta, freq) {
.Call('gmwm_ar1_to_gm', PACKAGE = 'gmwm', theta, freq)
}
#' @title Transform GM to AR1
#' @description Takes GM values and transforms them to AR1
#' @param theta A \code{vec} that contains AR1 values.
#' @param freq A \code{double} indicating the frequency of the data.
#' @return A \code{vec} containing GM values.
#' @keywords internal
#' @author JJB
#' The function takes a vector of GM values \eqn{\beta}{beta} and \eqn{\sigma ^2_{gm}}{sigma ^2[gm]}
#' and transforms them to AR1 values \eqn{\phi}{phi} and \eqn{\sigma ^2}{sigma ^2}
#' using the formulas:
#' \eqn{\phi = \exp \left( { - \beta \Delta t} \right)}{phi = exp(-beta * delta[t])}
#' \eqn{{\sigma ^2} = \sigma _{gm}^2\left( {1 - \exp \left( { - 2\beta \Delta t} \right)} \right)}{sigma^2 = sigma^2[gm]*(1-exp(-2*beta*delta[t]))}
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
#' @examples
#' gm_to_ar1(c(0.3,1,0.6,.3), 2)
gm_to_ar1 <- function(theta, freq) {
.Call('gmwm_gm_to_ar1', PACKAGE = 'gmwm', theta, freq)
}
#' @title Generate the ts model object description
#' @description Creates the ts.model's obj.desc value
#' @param desc A \code{vector<string>} that contains a list of the strings of each process.
#' @return A \code{field<vec>} that contains the object description of each process.
#' @details
#' This function currently does NOT support ARMA models.
#' That is, there is no support for ARMA, AR, or MA.
#' There is support for AR1, GM, WN, DR, QN, and RW.
#' @keywords internal
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
model_objdesc <- function(desc) {
.Call('gmwm_model_objdesc', PACKAGE = 'gmwm', desc)
}
#' @title Generate the ts model object's theta vector
#' @description Creates the ts.model's theta vector
#' @param desc A \code{vector<string>} that contains a list of the strings of each process.
#' @return A \code{vec} with values initialized at 0 that span the space of parameters to be estimated.
#' @details
#' This function currently does NOT support ARMA models.
#' That is, there is no support for ARMA, AR, or MA.
#' There is support for AR1, GM, WN, DR, QN, and RW.
#' @keywords internal
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
model_theta <- function(desc) {
.Call('gmwm_model_theta', PACKAGE = 'gmwm', desc)
}
#' @title Generate the ts model object's process desc
#' @description Creates the ts.model's process desc
#' @param desc A \code{vector<string>} that contains a list of the strings of each process.
#' @return A \code{vector<string>} with a list of descriptive values to label the estimate matrix with
#' @details
#' This function currently does NOT support ARMA models.
#' That is, there is no support for ARMA, AR, or MA.
#' There is support for AR1, GM, WN, DR, QN, and RW.
#' @keywords internal
#' @backref src/ts_model_cpp.cpp
#' @backref src/ts_model_cpp.h
model_process_desc <- function(desc) {
.Call('gmwm_model_process_desc', PACKAGE = 'gmwm', desc)
}
#' @title Generate eta3 confidence interval
#' @description Computes the eta3 CI
#' @param y A \code{vec} that computes the brickwalled modwt dot product of each wavelet coefficient divided by their length.
#' @param dims A \code{String} indicating the confidence interval being calculated.
#' @param alpha_ov_2 A \code{double} that indicates the \eqn{\left(1-p\right)*\alpha}{(1-p)*alpha} confidence level
#' @return A \code{matrix} with the structure:
#' \itemize{
#' \item{Column 1}{Wavelet Variance}
#' \item{Column 2}{Chi-squared Lower Bounds}
#' \item{Column 3}{Chi-squared Upper Bounds}
#' }
#' @keywords internal
#' @examples
#' x = rnorm(100)
#' # Uses the internal MODWT function not associated with an S3 class.
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' y = wave_variance(decomp)
#' ci_wave_variance(decomp, y, type = "eta3", alpha_ov_2 = 0.025)
ci_eta3 <- function(y, dims, alpha_ov_2) {
.Call('gmwm_ci_eta3', PACKAGE = 'gmwm', y, dims, alpha_ov_2)
}
#' @title Generate eta3 robust confidence interval
#' @description Computes the eta3 robust CI
#' @param wv_robust A \code{vec} that computes the brickwalled modwt dot product of each wavelet coefficient divided by their length.
#' @param wv_ci_class A \code{mat} that contains the CI mean, CI Lower, and CI Upper
#' @param alpha_ov_2 A \code{double} that indicates the \eqn{\left(1-p\right)*\alpha}{(1-p)*alpha} confidence level
#' @param eff A \code{double} that indicates the efficiency.
#' @return A \code{matrix} with the structure:
#' \itemize{
#' \item{Column 1}{Robust Wavelet Variance}
#' \item{Column 2}{Chi-squared Lower Bounds}
#' \item{Column 3}{Chi-squared Upper Bounds}
#' }
#' @details
#' Within this function we are scaling the classical
#' @keywords internal
#' @examples
#' x = rnorm(100)
#' # Uses the internal MODWT function not associated with an S3 class.
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' y = wave_variance(decomp, robust = TRUE, eff = 0.6)
#' ci_wave_variance(decomp, y, type = "eta3", alpha_ov_2 = 0.025, robust = TRUE, eff = 0.6)
ci_eta3_robust <- function(wv_robust, wv_ci_class, alpha_ov_2, eff) {
.Call('gmwm_ci_eta3_robust', PACKAGE = 'gmwm', wv_robust, wv_ci_class, alpha_ov_2, eff)
}
#' @title Generate a Confidence intervval for a Univariate Time Series
#' @description Computes an estimate of the multiscale variance and a chi-squared confidence interval
#' @param signal_modwt_bw A \code{field<vec>} that contains the brick walled modwt or dwt decomposition
#' @param wv A \code{vec} that contains the wave variance.
#' @param type A \code{String} indicating the confidence interval being calculated.
#' @param alpha_ov_2 A \code{double} that indicates the \eqn{\left(1-p\right)*\alpha}{(1-p)*alpha} confidence level.
#' @param robust A \code{boolean} to determine the type of wave estimation.
#' @param eff A \code{double} that indicates the efficiency.
#' @return A \code{matrix} with the structure:
#' \itemize{
#' \item{Column 1}{Wavelet Variance}
#' \item{Column 2}{Chi-squared Lower Bounds}
#' \item{Column 3}{Chi-squared Upper Bounds}
#' }
#' @keywords internal
#' @details
#' This function can be expanded to allow for other confidence interval calculations.
#' @examples
#' set.seed(1337)
#' x = rnorm(100)
#' # Uses the internal MODWT function not associated with an S3 class.
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' y = wave_variance(decomp)
#' ci_wave_variance(decomp, y, type = "eta3", alpha_ov_2 = 0.025)
ci_wave_variance <- function(signal_modwt_bw, wv, type = "eta3", alpha_ov_2 = 0.025, robust = FALSE, eff = 0.6) {
.Call('gmwm_ci_wave_variance', PACKAGE = 'gmwm', signal_modwt_bw, wv, type, alpha_ov_2, robust, eff)
}
#' @title Generate a Wave Variance for a Univariate Time Series
#' @description Computes an estimate of the wave variance
#' @param signal_modwt_bw A \code{field<vec>} that contains the brick walled modwt or dwt decomposition
#' @param robust A \code{boolean} to determine the type of wave estimation.
#' @param eff A \code{double} that indicates the efficiency.
#' @return A \code{vec} that contains the wave variance.
#' @keywords internal
#' @examples
#' set.seed(1337)
#' x = rnorm(100)
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' wave_variance(decomp)
#'
#' wave_variance(decomp, robust = TRUE, eff = 0.6)
wave_variance <- function(signal_modwt_bw, robust = FALSE, eff = 0.6) {
.Call('gmwm_wave_variance', PACKAGE = 'gmwm', signal_modwt_bw, robust, eff)
}
#' @title Computes the (MODWT) wavelet variance
#' @description Calculates the (MODWT) wavelet variance
#' @param signal_modwt_bw A \code{field<vec>} that contains the modwt decomposition after it has been brick walled.
#' @param robust A \code{boolean} that triggers the use of the robust estimate.
#' @param eff A \code{double} that indicates the efficiency as it relates to an MLE.
#' @param alpha A \code{double} that indicates the \eqn{\left(1-p\right)*\alpha}{(1-p)*alpha} confidence level
#' @param ci_type A \code{String} indicating the confidence interval being calculated. Valid value: "eta3"
#' @return A \code{mat} with the structure:
#' \itemize{
#' \item{"variance"}{Wavelet Variance}
#' \item{"low"}{Lower CI}
#' \item{"high"}{Upper CI}
#' }
#' @keywords internal
#' @details
#' This function does the heavy lifting with the signal_modwt_bw
#' @examples
#' x = rnorm(100)
#' decomp = modwt_cpp(x, filter_name = "haar", nlevels = 4, boundary = "periodic", brickwall = TRUE)
#' wvar_cpp(decomp, robust=FALSE, eff=0.6, alpha = 0.05, ci_type="eta3")
wvar_cpp <- function(signal_modwt_bw, robust, eff, alpha, ci_type) {
.Call('gmwm_wvar_cpp', PACKAGE = 'gmwm', signal_modwt_bw, robust, eff, alpha, ci_type)
}
#' @title Computes the (MODWT) wavelet variance
#' @description Calculates the (MODWT) wavelet variance
#' @param signal A \code{vec} that contains the data.
#' @param robust A \code{boolean} that triggers the use of the robust estimate.
#' @param eff A \code{double} that indicates the efficiency as it relates to an MLE.
#' @param alpha A \code{double} that indicates the \eqn{\left(1-p\right)\times \alpha}{(1-p)*alpha} confidence level
#' @param ci_type A \code{string} indicating the confidence interval being calculated. Valid value: "eta3"
#' @param strWavelet A \code{string} indicating the type of wave filter to be applied. Must be "haar"
#' @param decomp A \code{string} indicating whether to use "modwt" or "dwt" decomp
#' @return A \code{mat} with the structure:
#' \itemize{
#' \item{"variance"}{Wavelet Variance}
#' \item{"low"}{Lower CI}
#' \item{"high"}{Upper CI}
#' }
#' @keywords internal
#' @details
#' This function powers the wvar object. It is also extendable...
#' @examples
#' x=rnorm(100)
#' modwt_wvar_cpp(x, nlevels=4, robust=FALSE, eff=0.6, alpha = 0.05,
#' ci_type="eta3", strWavelet="haar", decomp="modwt")
modwt_wvar_cpp <- function(signal, nlevels, robust, eff, alpha, ci_type, strWavelet, decomp) {
.Call('gmwm_modwt_wvar_cpp', PACKAGE = 'gmwm', signal, nlevels, robust, eff, alpha, ci_type, strWavelet, decomp)
}
#' @title Computes the MO/DWT wavelet variance for multiple processes
#' @description Calculates the MO/DWT wavelet variance
#' @param signal A \code{matrix} that contains the same number of observations per dataset
#' @param robust A \code{boolean} that triggers the use of the robust estimate.
#' @param eff A \code{double} that indicates the efficiency as it relates to an MLE.
#' @param alpha A \code{double} that indicates the \eqn{\left(1-p\right)\times \alpha}{(1-p)*alpha} confidence level
#' @param ci_type A \code{string} indicating the confidence interval being calculated. Valid value: "eta3"
#' @param strWavelet A \code{string} indicating the type of wave filter to be applied. Must be "haar"
#' @param decomp A \code{string} indicating whether to use "modwt" or "dwt" decomp
#' @return A \code{field<mat>} with the structure:
#' \itemize{
#' \item{"variance"}{Wavelet Variance}
#' \item{"low"}{Lower CI}
#' \item{"high"}{Upper CI}
#' }
#' @keywords internal
#' @details
#' This function processes the decomposition of multiple signals quickly
#' @examples
#' x = cbind(rnorm(100),rnorm(100))
#' batch_modwt_wvar_cpp(x, nlevels=4, robust=FALSE, eff=0.6,
#' alpha = 0.05, ci_type="eta3", strWavelet="haar",
#' decomp="modwt")
batch_modwt_wvar_cpp <- function(signal, nlevels, robust, eff, alpha, ci_type, strWavelet, decomp) {
.Call('gmwm_batch_modwt_wvar_cpp', PACKAGE = 'gmwm', signal, nlevels, robust, eff, alpha, ci_type, strWavelet, decomp)
}
#' @title Computes the MODWT scales
#' @description Calculates the MODWT scales
#' @param nb_level A \code{integer} that contains the level of decomposition J.
#' @return A \code{vec} that contains 2^1, ... , 2^J
#' @keywords internal
#' @details
#' Used in wvar object.
#' @examples
#' scales_cpp(5)
scales_cpp <- function(nb_level) {
.Call('gmwm_scales_cpp', PACKAGE = 'gmwm', nb_level)
}
#' @title Quadrature Mirror Filter
#' @description Calculate the series quadrature mirror filter (QMF). Requires a series of an even length.
#' @usage qmf(g, inverse)
#' @param g A \code{vector} that contains the filter constants.
#' @param inverse A \code{bool} that indicates whether the inverse quadrature mirror filter is computed.
#' By default, the inverse quadrature mirror is computed.
#' @return A \code{vector} that contains either the forward QMF (evalute in order) or the inverse QMF (reverse order).
#' @author JJB
#' @keywords internal
#' @examples
#' # Haar values
#' g = rep(1/sqrt(2),2)
#' qmf(g)
qmf <- function(g, inverse = TRUE) {
.Call('gmwm_qmf', PACKAGE = 'gmwm', g, inverse)
}
#' @title Haar filter construction
#' @description Creates the haar filter
#' @usage haar_filter()
#' @return A \code{field<vec>} that contains:
#' \itemize{
#' \item{"L"}{A \code{integer} specifying the length of the filter}
#' \item{"h"}{A \code{vector} containing the coefficients for the wavelet filter}
#' \item{"g"}{A \code{vector} containing the coefficients for the scaling filter}
#' }
#' @details
#' This template can be used to increase the amount of filters available for selection.
#' @author JJB
#' @keywords internal
#' @examples
#' haar_filter()
haar_filter <- function() {
.Call('gmwm_haar_filter', PACKAGE = 'gmwm')
}
#' @title Select the Wavelet Filter
#' @description Constructs the wavelet filter to be used.
#' @usage select_filter(filter_name)
#' @param filter_name A \code{String} that must receive: \code{"haar"}.
#' @return info A \code{field<vec>} that contains:
#' \itemize{
#' \item{"L"}{A \code{integer} specifying the length of the filter}
#' \item{"h"}{A \code{vector} containing the coefficients for the wavelet filter}
#' \item{"g"}{A \code{vector} containing the coefficients for the scaling filter}
#' }
#' @details
#' The package is oriented toward using only the haar filter. If the package extends at a later time, then the supporting infrastructure is there.
#' @author JJB
#' @keywords internal
#' @examples
#' select_filter("haar")
select_filter <- function(filter_name = "haar") {
.Call('gmwm_select_filter', PACKAGE = 'gmwm', filter_name)
}
Try the gmwm package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.