R/RcppExports.R

In shabbychef/fromo: Fast Robust Moments

# Generated by using Rcpp::compileAttributes() -> do not edit by hand
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393

#' @title
#' Compute first K moments
#' @description
#' Compute the (standardized) 2nd through kth moments, the mean, and the number of elements.
#' 
#' @param v a vector
#' @param na_rm whether to remove NA, false by default.
#' @param max_order the maximum order of the centered moment to be computed.
#' @param used_df the number of degrees of freedom consumed, used in the denominator
#' of the centered moments computation. These are subtracted from the number of
#' observations. 
#' @param sg_df the number of degrees of freedom consumed in the computation of
#' the variance or standard deviation. This defaults to 1 to match the 
#' \sQuote{Bessel correction}.
#'
#' @details
#'
#' Computes the number of elements, the mean, and the 2nd through kth
#' centered standardized moment, for \eqn{k=2,3,4}{k=2,3,4}. These
#' are computed via the numerically robust one-pass method of Bennett \emph{et. al.}
#' In general they will \emph{not} match exactly with the 'standard'
#' implementations, due to differences in roundoff.
#'
#' These methods are reasonably fast, on par with the 'standard' implementations.
#' However, they will usually be faster than calling the various standard implementations
#' more than once.
#'
#' Moments are computed as follows, given some values \eqn{x_i} and optional weights \eqn{w_i},
#' defaulting to 1, the weighted mean is computed as
#' \deqn{\mu = \frac{\sum_i x_i w_i}{\sum w_i}.}
#' The weighted kth central sum is computed as
#' \deqn{\mu = \sum_i \left(x_i - \mu\right)^k w_i.}
#' Let \eqn{n = \sum_i w_i} be the sum of weights (or number of observations in the unweighted case).
#' Then the weighted kth central moment is computed as that weighted sum divided by the
#' adjusted sum weights:
#' \deqn{\mu_k = \frac{\sum_i \left(x_i - \mu\right)^k w_i}{n - \nu},}
#' where \eqn{\nu} is the \sQuote{used df}, provided by the user to adjust the denominator.
#' (Typical values are 0 or 1.)
#' The weighted kth standardized moment is the central moment divided by the second central moment
#' to the \eqn{k/2} power:
#' \deqn{\tilde{\mu}_k = \frac{\mu_k}{\mu_2^{k/2}}.}
#' The (centered) rth cumulant, for \eqn{r \ge 2} is then computed using the formula of Willink, namely
#' \deqn{\kappa_r = \mu_r - \sum_{j=0}^{r - 2} {r - 1 \choose j} \mu_j \kappa {r-j}.}
#' The standardized rth cumulant is the rth centered cumulant divided by \eqn{\mu_2^{r/2}}.
#'
#' @return a vector, filled out as follows:
#' \describe{
#' \item{sd3}{A vector of the (sample) standard devation, mean, and number of elements (or the total weight when \code{wts}
#' are given).}
#' \item{skew4}{A vector of the (sample) skewness, standard devation, mean, and number of elements (or the total weight when 
#' \code{wts} are given).}
#' \item{kurt5}{A vector of the (sample) excess kurtosis, skewness, standard devation, mean, and number of elements (or the
#' total weight when \code{wts} are given).}
#' \item{cent_moments}{A vector of the (sample) \eqn{k}th centered moment, then \eqn{k-1}th centered moment, ..., 
#'  then the \emph{variance}, the mean, and number of elements (total weight when \code{wts} are given).}
#' \item{std_moments}{A vector of the (sample) \eqn{k}th standardized (and centered) moment, then 
#'  \eqn{k-1}th, ..., then standard devation, mean, and number of elements (total weight).}
#' \item{cent_cumulants}{A vector of the (sample) \eqn{k}th (centered, but this is redundant) cumulant, then the \eqn{k-1}th, ...,
#'  then the \emph{variance} (which is the second cumulant), then \emph{the mean}, then the number of elements (total weight).}
#' \item{std_cumulants}{A vector of the (sample) \eqn{k}th standardized (and centered, but this is redundant) cumulant, then the \eqn{k-1}th, ...,
#'  down to the third, then \emph{the variance}, \emph{the mean}, then the number of elements (total weight).}
#' }
#'
#' @note
#' The first centered (and standardized) moment is often defined to be identically 0. Instead \code{cent_moments}
#' and \code{std_moments} returns the mean. 
#' Similarly, the second standardized moments defined to be identically 1; \code{std_moments} instead returns the standard
#' deviation. The reason is that a user can always decide to ignore the results and fill in a 0 or 1 as they need, but 
#' could not efficiently compute the mean and standard deviation from scratch if we discard it.
#' The antepenultimate element of the output of \code{std_cumulants} is not a one, even though that \sQuote{should} be
#' the standardized second cumulant.
#' 
#' @note
#' The antepenultimate element of the output of \code{cent_moments}, \code{cent_cumulants} and \code{std_cumulants} is the \emph{variance},
#' not the standard deviation. All other code return the standard deviation in that place.
#'
#' @note
#' The kurtosis is \emph{excess kurtosis}, with a 3 subtracted, and should be nearly zero
#' for Gaussian input.
#'
#' @note
#' The term 'centered cumulants' is redundant. The intent was to avoid possible collision with existing code named 'cumulants'.
#'
#' @examples
#' x <- rnorm(1e5)
#' sd3(x)[1] - sd(x)
#' skew4(x)[4] - length(x)
#' skew4(x)[3] - mean(x)
#' skew4(x)[2] - sd(x)
#' if (require(moments)) {
#'   skew4(x)[1] - skewness(x)
#' }
#'
#'
#' # check 'robustness'; only the mean should change:
#' kurt5(x + 1e12) - kurt5(x)
#' # check speed
#' if (require(microbenchmark) && require(moments)) {
#'   dumbk <- function(x) { c(kurtosis(x) - 3.0,skewness(x),sd(x),mean(x),length(x)) }
#'   set.seed(1234)
#'   x <- rnorm(1e6)
#'   print(kurt5(x) - dumbk(x))
#'   microbenchmark(dumbk(x),kurt5(x),times=10L)
#' }
#' y <- std_moments(x,6)
#' cml <- cent_cumulants(x,6)
#' std <- std_cumulants(x,6)
#' 
#' # check that skew matches moments::skewness
#' if (require(moments)) {
#'     set.seed(1234)
#'     x <- rnorm(1000)
#'     resu <- fromo::skew4(x)
#' 
#'     msku <- moments::skewness(x)
#'     stopifnot(abs(msku - resu[1]) < 1e-14)
#' }
#' 
#' # check skew vs e1071 skewness, which has a different denominator
#' if (require(e1071)) {
#'     set.seed(1234)
#'     x <- rnorm(1000)
#'     resu <- fromo::skew4(x)
#' 
#'     esku <- e1071::skewness(x,type=3)
#'     nobs <- resu[4]
#'     stopifnot(abs(esku - resu[1] * ((nobs-1)/nobs)^(3/2)) < 1e-14)
#'
#'     # similarly:
#'     resu <- fromo::std_moments(x,max_order=3,used_df=0)
#'     stopifnot(abs(esku - resu[1] * ((nobs-1)/nobs)^(3/2)) < 1e-14)
#' }
#'
#' @template etc
#' @template ref-romo
#' @template ref-cumulants
#' @template param-wts
#' @template note-wts
#' @rdname firstmoments
#' @export
sd3 <- function(v, na_rm = FALSE, wts = NULL, sg_df = 1.0, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_sd3', PACKAGE = 'fromo', v, na_rm, wts, sg_df, check_wts, normalize_wts)
}

#' @rdname firstmoments
#' @export
skew4 <- function(v, na_rm = FALSE, wts = NULL, sg_df = 1.0, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_skew4', PACKAGE = 'fromo', v, na_rm, wts, sg_df, check_wts, normalize_wts)
}

#' @rdname firstmoments
#' @export
kurt5 <- function(v, na_rm = FALSE, wts = NULL, sg_df = 1.0, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_kurt5', PACKAGE = 'fromo', v, na_rm, wts, sg_df, check_wts, normalize_wts)
}

#' @rdname firstmoments
#' @export
cent_moments <- function(v, max_order = 5L, used_df = 0L, na_rm = FALSE, wts = NULL, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_cent_moments', PACKAGE = 'fromo', v, max_order, used_df, na_rm, wts, check_wts, normalize_wts)
}

#' @rdname firstmoments
#' @export
std_moments <- function(v, max_order = 5L, used_df = 0L, na_rm = FALSE, wts = NULL, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_std_moments', PACKAGE = 'fromo', v, max_order, used_df, na_rm, wts, check_wts, normalize_wts)
}

#' @rdname firstmoments
#' @export
cent_cumulants <- function(v, max_order = 5L, used_df = 0L, na_rm = FALSE, wts = NULL, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_cent_cumulants', PACKAGE = 'fromo', v, max_order, used_df, na_rm, wts, check_wts, normalize_wts)
}

#' @rdname firstmoments
#' @export
std_cumulants <- function(v, max_order = 5L, used_df = 0L, na_rm = FALSE, wts = NULL, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_std_cumulants', PACKAGE = 'fromo', v, max_order, used_df, na_rm, wts, check_wts, normalize_wts)
}

#' @title
#' Centered sums; join and unjoined.
#'
#' @description
#'
#' Compute, join, or unjoin centered sums.
#'
#' @param ret1 an \eqn{ord+1} vector as output by \code{\link{cent_sums}} consisting of
#' the count, the mean, then the k through ordth centered sum of some observations.
#' @param ret2 an \eqn{ord+1} vector as output by \code{\link{cent_sums}} consisting of
#' the count, the mean, then the k through ordth centered sum of some observations.
#' @param ret3 an \eqn{ord+1} vector as output by \code{\link{cent_sums}} consisting of
#' the count, the mean, then the k through ordth centered sum of some observations.
#' @inheritParams cent_moments
#'
#' @return a vector the same size as the input consisting of the adjusted version of the input.
#' When there are not sufficient (non-nan) elements for the computation, \code{NaN} are returned.
#'
#' @examples
#'
#'  set.seed(1234)
#'  x1 <- rnorm(1e3,mean=1)
#'  x2 <- rnorm(1e3,mean=1)
#'  max_ord <- 6L
#'  rs1 <- cent_sums(x1,max_ord)
#'  rs2 <- cent_sums(x2,max_ord)
#'  rs3 <- cent_sums(c(x1,x2),max_ord)
#'  rs3alt <- join_cent_sums(rs1,rs2)
#'  stopifnot(max(abs(rs3 - rs3alt)) < 1e-7)
#'  rs1alt <- unjoin_cent_sums(rs3,rs2)
#'  rs2alt <- unjoin_cent_sums(rs3,rs1)
#'  stopifnot(max(abs(rs1 - rs1alt)) < 1e-7)
#'  stopifnot(max(abs(rs2 - rs2alt)) < 1e-7)
#'
#' @template etc
#' @template ref-romo
#' @template param-wts
#' @rdname centsums 
#' @export
cent_sums <- function(v, max_order = 5L, na_rm = FALSE, wts = NULL, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_cent_sums', PACKAGE = 'fromo', v, max_order, na_rm, wts, check_wts, normalize_wts)
}

#' @rdname centsums 
#' @export
join_cent_sums <- function(ret1, ret2) {
    .Call('_fromo_join_cent_sums', PACKAGE = 'fromo', ret1, ret2)
}

#' @rdname centsums 
#' @export
unjoin_cent_sums <- function(ret3, ret2) {
    .Call('_fromo_unjoin_cent_sums', PACKAGE = 'fromo', ret3, ret2)
}

#' @title
#' Multivariate centered sums; join and unjoined.
#'
#' @description
#'
#' Compute, join, or unjoin multivariate centered (co-) sums.
#'
#' @param v an \eqn{m} by \eqn{n} matrix, each row an independent observation of some
#' \eqn{n} variate variable.
#' @param max_order the maximum order of cosum to compute. For now this can only be
#' 2; in the future higher order cosums should be possible.
#' @param na_omit a boolean; if \code{TRUE}, then only rows of \code{v} with complete
#' observations will be used.
#' @param ret1 a multdimensional array as output by \code{\link{cent_cosums}}.
#' @param ret2 a multdimensional array as output by \code{\link{cent_cosums}}.
#' @param ret3 a multdimensional array as output by \code{\link{cent_cosums}}.
#' @param used_df the number of degrees of freedom consumed, used in the denominator
#' of the centered moments computation. These are subtracted from the number of
#' observations.
#'
#' @return a multidimensional arry of dimension \code{max_order}, each side of length
#' \eqn{1+n}. For the case currently implemented where \code{max_order} must be 2, the
#' output is a symmetric matrix, where the element in the \code{1,1} position is the count of 
#' complete) rows of \code{v}, the \code{2:(n+1),1} column is the mean, and the
#' \code{2:(n+1),2:(n+1)} is the co \emph{sums} matrix, which is the covariance up to scaling
#' by the count. \code{cent_comoments} performs this normalization for you.
#'
#' @seealso cent_sums
#'
#' @examples
#'
#'  set.seed(1234)
#'  x1 <- matrix(rnorm(1e3*5,mean=1),ncol=5)
#'  x2 <- matrix(rnorm(1e3*5,mean=1),ncol=5)
#'  max_ord <- 2L
#'  rs1 <- cent_cosums(x1,max_ord)
#'  rs2 <- cent_cosums(x2,max_ord)
#'  rs3 <- cent_cosums(rbind(x1,x2),max_ord)
#'  rs3alt <- join_cent_cosums(rs1,rs2)
#'  stopifnot(max(abs(rs3 - rs3alt)) < 1e-7)
#'  rs1alt <- unjoin_cent_cosums(rs3,rs2)
#'  rs2alt <- unjoin_cent_cosums(rs3,rs1)
#'  stopifnot(max(abs(rs1 - rs1alt)) < 1e-7)
#'  stopifnot(max(abs(rs2 - rs2alt)) < 1e-7)
#'
#' @template etc
#' @template ref-romo
#' @rdname centcosums 
#' @export
cent_cosums <- function(v, max_order = 2L, na_omit = FALSE) {
    .Call('_fromo_cent_cosums', PACKAGE = 'fromo', v, max_order, na_omit)
}

#' @rdname centcosums 
#' @export
cent_comoments <- function(v, max_order = 2L, used_df = 0L, na_omit = FALSE) {
    .Call('_fromo_cent_comoments', PACKAGE = 'fromo', v, max_order, used_df, na_omit)
}

#' @rdname centcosums 
#' @export
join_cent_cosums <- function(ret1, ret2) {
    .Call('_fromo_join_cent_cosums', PACKAGE = 'fromo', ret1, ret2)
}

#' @rdname centcosums 
#' @export
unjoin_cent_cosums <- function(ret3, ret2) {
    .Call('_fromo_unjoin_cent_cosums', PACKAGE = 'fromo', ret3, ret2)
}

#' @title
#' Convert between different types of moments, raw, central, standardized.
#' @description
#' Given raw or central or standardized moments, convert to another type.
#' 
#' @param input a vector of the count, then the mean, then the \code{2} through \code{k}
#' raw or central moments.
#'
#' @template etc
#' @rdname moment_conversions
#' @export
cent2raw <- function(input) {
    .Call('_fromo_cent2raw', PACKAGE = 'fromo', input)
}

#' @title
#' Compute first K moments over a sliding window
#' @description
#' Compute the (standardized) 2nd through kth moments, the mean, and the number of elements over
#' an infinite or finite sliding window, returning a matrix.
#' 
#' @param v a vector
#' @param window the window size. if given as finite integer or double, passed through.
#' If \code{NULL}, \code{NA_integer_}, \code{NA_real_} or \code{Inf} are given, equivalent
#' to an infinite window size. If negative, an error will be thrown.
#' @param restart_period the recompute period. because subtraction of elements can cause
#' loss of precision, the computation of moments is restarted periodically based on 
#' this parameter. Larger values mean fewer restarts and faster, though less accurate
#' results. 
#' @param na_rm whether to remove NA, false by default.
#' @param max_order the maximum order of the centered moment to be computed.
#' @param min_df the minimum df to return a value, otherwise \code{NaN} is returned.
#' This can be used to prevent moments from being computed on too few observations.
#' Defaults to zero, meaning no restriction.
#' @param used_df the number of degrees of freedom consumed, used in the denominator
#' of the centered moments computation. These are subtracted from the number of
#' observations.
#'
#' @details
#'
#' Computes the number of elements, the mean, and the 2nd through kth
#' centered (and typically standardized) moments, for \eqn{k=2,3,4}{k=2,3,4}. These
#' are computed via the numerically robust one-pass method of Bennett \emph{et. al.}
#'
#' Given the length \eqn{n} vector \eqn{x}, we output matrix \eqn{M} where
#' \eqn{M_{i,j}}{M_i,j} is the \eqn{order - j + 1} moment (\emph{i.e.}
#' excess kurtosis, skewness, standard deviation, mean or number of elements)
#' of \eqn{x_{i-window+1},x_{i-window+2},...,x_{i}}{x_(i-window+1),x_(i-window+2),...,x_i}.
#' Barring \code{NA} or \code{NaN}, this is over a window of size \code{window}.
#' During the 'burn-in' phase, we take fewer elements.
#'
#' @return Typically a matrix, where the first columns are the kth, k-1th through 2nd standardized, 
#' centered moments, then a column of the mean, then a column of the number of (non-nan) elements in the input,
#' with the following exceptions:
#' \describe{
#' \item{running_cent_moments}{Computes arbitrary order centered moments. When \code{max_order_only} is set,
#' only a column of the maximum order centered moment is returned.}
#' \item{running_std_moments}{Computes arbitrary order standardized moments, then the standard deviation, the mean,
#' and the count. There is not yet an option for \code{max_order_only}, but probably should be.}
#' \item{running_cumulants}{Computes arbitrary order cumulants, and returns the kth, k-1th, through the second 
#' (which is the variance) cumulant, then the mean, and the count.}
#' }
#'
#' @note
#' the kurtosis is \emph{excess kurtosis}, with a 3 subtracted, and should be nearly zero
#' for Gaussian input.
#'
#' @examples
#' x <- rnorm(1e5)
#' xs3 <- running_sd3(x,10)
#' xs4 <- running_skew4(x,10)
#'
#' if (require(moments)) {
#'     set.seed(123)
#'     x <- rnorm(5e1)
#'     window <- 10L
#'     kt5 <- running_kurt5(x,window=window)
#'     rm1 <- t(sapply(seq_len(length(x)),function(iii) { 
#'                 xrang <- x[max(1,iii-window+1):iii]
#'                 c(moments::kurtosis(xrang)-3.0,moments::skewness(xrang),
#'                 sd(xrang),mean(xrang),length(xrang)) },
#'              simplify=TRUE))
#'     stopifnot(max(abs(kt5 - rm1),na.rm=TRUE) < 1e-12)
#' }
#'
#' xc6 <- running_cent_moments(x,window=100L,max_order=6L)
#'
#' @template etc
#' @template ref-romo
#' @template param-wts
#' @template note-wts
#' @template note-heywood
#' @rdname runningmoments
#' @export
running_sd3 <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_sd3', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @rdname runningmoments
#' @export
running_skew4 <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_skew4', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @rdname runningmoments
#' @export
running_kurt5 <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_kurt5', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @rdname runningmoments
#' @export
running_sd <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_sd', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @rdname runningmoments
#' @export
running_skew <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_skew', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @rdname runningmoments
#' @export
running_kurt <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_kurt', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @param max_order_only for \code{running_cent_moments}, if this flag is set, only compute
#' the maximum order centered moment, and return in a vector.
#' @rdname runningmoments
#' @export
running_cent_moments <- function(v, window = NULL, wts = NULL, max_order = 5L, na_rm = FALSE, max_order_only = FALSE, min_df = 0L, used_df = 0.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_cent_moments', PACKAGE = 'fromo', v, window, wts, max_order, na_rm, max_order_only, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @rdname runningmoments
#' @export
running_std_moments <- function(v, window = NULL, wts = NULL, max_order = 5L, na_rm = FALSE, min_df = 0L, used_df = 0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_std_moments', PACKAGE = 'fromo', v, window, wts, max_order, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @rdname runningmoments
#' @export
running_cumulants <- function(v, window = NULL, wts = NULL, max_order = 5L, na_rm = FALSE, min_df = 0L, used_df = 0.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_cumulants', PACKAGE = 'fromo', v, window, wts, max_order, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @title
#' Compute approximate quantiles over a sliding window
#' @description
#' Computes cumulants up to some given order, then employs the Cornish-Fisher approximation
#' to compute approximate quantiles using a Gaussian basis.
#' 
#' @param p the probability points at which to compute the quantiles. Should be in the range (0,1).
#' @inheritParams running_cumulants
#'
#' @details
#'
#' Computes the cumulants, then approximates quantiles using AS269 of Lee & Lin.
#'
#' @return A matrix, with one row for each element of \code{x}, and one column for each element of \code{q}.
#'
#' @note
#' The current implementation is not as space-efficient as it could be, as it first computes
#' the cumulants for each row, then performs the Cornish-Fisher approximation on a row-by-row
#' basis. In the future, this computation may be moved earlier into the pipeline to be more
#' space efficient. File an issue if the memory footprint is an issue for you.
#'
#' @examples
#' x <- rnorm(1e5)
#' xq <- running_apx_quantiles(x,c(0.1,0.25,0.5,0.75,0.9))
#' xm <- running_apx_median(x)
#'
#' @seealso \code{\link{t_running_apx_quantiles}}, \code{\link{running_cumulants}}, \code{PDQutils::qapx_cf}, \code{PDQutils::AS269}.
#' @template etc
#' @template ref-cf
#' @template ref-romo
#' @template param-wts
#' @template note-wts
#' @rdname runningquantiles
#' @export
running_apx_quantiles <- function(v, p, window = NULL, wts = NULL, max_order = 5L, na_rm = FALSE, min_df = 0L, used_df = 0.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_apx_quantiles', PACKAGE = 'fromo', v, p, window, wts, max_order, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @rdname runningquantiles
#' @export
running_apx_median <- function(v, window = NULL, wts = NULL, max_order = 5L, na_rm = FALSE, min_df = 0L, used_df = 0.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_apx_median', PACKAGE = 'fromo', v, window, wts, max_order, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @title
#' Compare data to moments computed over a sliding window.
#' @description
#' Computes moments over a sliding window, then adjusts the data accordingly, centering, or scaling,
#' or z-scoring, and so on.
#' 
#' @inheritParams running_cent_moments
#' @param min_df the minimum df to return a value, otherwise \code{NaN} is returned.
#' This can be used to prevent \emph{e.g.} Z-scores from being computed on only 3
#' observations. Defaults to zero, meaning no restriction, which can result in 
#' infinite Z-scores during the burn-in period.
#' @param lookahead for some of the operations, the value is compared to 
#' mean and standard deviation possibly using 'future' or 'past' information
#' by means of a non-zero lookahead. Positive values mean data are taken from
#' the future.
#' @param compute_se for \code{running_sharpe}, return an extra column of the
#' standard error, as computed by Mertens' correction.
#'
#' @details
#'
#' Given the length \eqn{n} vector \eqn{x}, for
#' a given index \eqn{i}, define \eqn{x^{(i)}}{x^(i)}
#' as the vector of 
#' \eqn{x_{i-window+1},x_{i-window+2},...,x_{i}}{x_(i-window+1),x_(i-window+2),...,x_i},
#' where we do not run over the 'edge' of the vector. In code, this is essentially
#' \code{x[(max(1,i-window+1)):i]}. Then define \eqn{\mu_i}{mu_i}, \eqn{\sigma_i}{sigma_i}
#' and \eqn{n_i}{n_i} as, respectively, the sample mean, standard deviation and number of
#' non-NA elements in \eqn{x^{(i)}}{x^(i)}. 
#'
#' We compute output vector \eqn{m} the same size as \eqn{x}. 
#' For the 'centered' version of \eqn{x}, we have \eqn{m_i = x_i - \mu_i}{m_i = x_i - mu_i}.
#' For the 'scaled' version of \eqn{x}, we have \eqn{m_i = x_i / \sigma_i}{m_i = x_i / sigma_i}.
#' For the 'z-scored' version of \eqn{x}, we have \eqn{m_i = (x_i - \mu_i) / \sigma_i}{m_i = (x_i - mu_i) / sigma_i}.
#' For the 't-scored' version of \eqn{x}, we have \eqn{m_i = \sqrt{n_i} \mu_i / \sigma_i}{m_i = sqrt(n_i) mu_i / sigma_i}.
#'
#' We also allow a 'lookahead' for some of these operations.
#' If positive, the moments are computed using data from larger indices;
#' if negative, from smaller indices. Letting \eqn{j = i + lookahead}{j = i + lookahead}:
#' For the 'centered' version of \eqn{x}, we have \eqn{m_i = x_i - \mu_j}{m_i = x_i - mu_j}.
#' For the 'scaled' version of \eqn{x}, we have \eqn{m_i = x_i / \sigma_j}{m_i = x_i / sigma_j}.
#' For the 'z-scored' version of \eqn{x}, we have \eqn{m_i = (x_i - \mu_j) / \sigma_j}{m_i = (x_i - mu_j) / sigma_j}.
#'
#' @return a vector the same size as the input consisting of the adjusted version of the input.
#' When there are not sufficient (non-nan) elements for the computation, \code{NaN} are returned.
#'
#' @examples
#'
#' if (require(moments)) {
#'     set.seed(123)
#'     x <- rnorm(5e1)
#'     window <- 10L
#'     rm1 <- t(sapply(seq_len(length(x)),function(iii) { 
#'                   xrang <- x[max(1,iii-window+1):iii]
#'                   c(sd(xrang),mean(xrang),length(xrang)) },
#'                   simplify=TRUE))
#'     rcent <- running_centered(x,window=window)
#'     rscal <- running_scaled(x,window=window)
#'     rzsco <- running_zscored(x,window=window)
#'     rshrp <- running_sharpe(x,window=window)
#'     rtsco <- running_tstat(x,window=window)
#'     rsrse <- running_sharpe(x,window=window,compute_se=TRUE)
#'     stopifnot(max(abs(rcent - (x - rm1[,2])),na.rm=TRUE) < 1e-12)
#'     stopifnot(max(abs(rscal - (x / rm1[,1])),na.rm=TRUE) < 1e-12)
#'     stopifnot(max(abs(rzsco - ((x - rm1[,2]) / rm1[,1])),na.rm=TRUE) < 1e-12)
#'     stopifnot(max(abs(rshrp - (rm1[,2] / rm1[,1])),na.rm=TRUE) < 1e-12)
#'     stopifnot(max(abs(rtsco - ((sqrt(rm1[,3]) * rm1[,2]) / rm1[,1])),na.rm=TRUE) < 1e-12)
#'     stopifnot(max(abs(rsrse[,1] - rshrp),na.rm=TRUE) < 1e-12)
#'
#'     rm2 <- t(sapply(seq_len(length(x)),function(iii) { 
#'                   xrang <- x[max(1,iii-window+1):iii]
#'                   c(kurtosis(xrang)-3.0,skewness(xrang)) },
#'                   simplify=TRUE))
#'     mertens_se <- sqrt((1 + ((2 + rm2[,1])/4) * rshrp^2 - rm2[,2]*rshrp) / rm1[,3])
#'     stopifnot(max(abs(rsrse[,2] - mertens_se),na.rm=TRUE) < 1e-12)
#' }
#'
#' @seealso \code{\link{t_running_centered}}, \code{\link{scale}}
#' @template etc
#' @template ref-romo
#' @template note-wts
#' @rdname runningadjustments
#' @export
running_centered <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, lookahead = 0L, restart_period = 100L, check_wts = FALSE, normalize_wts = FALSE) {
    .Call('_fromo_running_centered', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, lookahead, restart_period, check_wts, normalize_wts)
}

#' @rdname runningadjustments
#' @export
running_scaled <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, lookahead = 0L, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_scaled', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, lookahead, restart_period, check_wts, normalize_wts)
}

#' @rdname runningadjustments
#' @export
running_zscored <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, lookahead = 0L, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_zscored', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, lookahead, restart_period, check_wts, normalize_wts)
}

#' @rdname runningadjustments
#' @export
running_sharpe <- function(v, window = NULL, wts = NULL, na_rm = FALSE, compute_se = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_sharpe', PACKAGE = 'fromo', v, window, wts, na_rm, compute_se, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @rdname runningadjustments
#' @export
running_tstat <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_running_tstat', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, used_df, restart_period, check_wts, normalize_wts)
}

#' @title
#' Compute sums or means over a sliding window.
#'
#' @description
#' Compute the mean or sum over 
#' an infinite or finite sliding window, returning a vector the same size as the input.
#' 
#' @param v a vector.
#' @param window the window size. if given as finite integer or double, passed through.
#' If \code{NULL}, \code{NA_integer_}, \code{NA_real_} or \code{Inf} are given, equivalent
#' to an infinite window size. If negative, an error will be thrown.
#' @param restart_period the recompute period. because subtraction of elements can cause
#' loss of precision, the computation of moments is restarted periodically based on 
#' this parameter. Larger values mean fewer restarts and faster, though potentially less 
#' accurate results. Unlike in the computation of even order moments, loss of precision
#' is unlikely to be disastrous, so the default value is rather large.
#' @param na_rm whether to remove NA, false by default.
#' @param min_df the minimum df to return a value, otherwise \code{NaN} is returned,
#' only for the means computation.
#' This can be used to prevent moments from being computed on too few observations.
#' Defaults to zero, meaning no restriction.
#'
#' @details
#'
#' Computes the mean or sum of the elements, using a Kahan's Compensated Summation Algorithm,
#' a numerically robust one-pass method.
#'
#' Given the length \eqn{n} vector \eqn{x}, we output matrix \eqn{M} where
#' \eqn{M_{i,1}}{M_i,1} is the sum or mean 
#' of \eqn{x_{i-window+1},x_{i-window+2},...,x_{i}}{x_(i-window+1),x_(i-window+2),...,x_i}.
#' Barring \code{NA} or \code{NaN}, this is over a window of size \code{window}.
#' During the 'burn-in' phase, we take fewer elements. If fewer than \code{min_df} for
#' \code{running_mean}, returns \code{NA}.
#'
#' @return A vector the same size as the input.
#' @examples
#' x <- rnorm(1e5)
#' xs <- running_sum(x,10)
#' xm <- running_mean(x,100)
#'
#' @template etc
#' @template ref-romo
#' @template ref-kahan
#' @inheritParams sd3
#' @rdname runningmean 
#' @export
running_sum <- function(v, window = NULL, wts = NULL, na_rm = FALSE, restart_period = 10000L, check_wts = FALSE) {
    .Call('_fromo_running_sum', PACKAGE = 'fromo', v, window, wts, na_rm, restart_period, check_wts)
}

#' @rdname runningmean
#' @export
running_mean <- function(v, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, restart_period = 10000L, check_wts = FALSE) {
    .Call('_fromo_running_mean', PACKAGE = 'fromo', v, window, wts, na_rm, min_df, restart_period, check_wts)
}

#' @title
#' Compute first K moments over a sliding time-based window
#' @description
#' Compute the (standardized) 2nd through kth moments, the mean, and the number of elements over
#' an infinite or finite sliding time based window, returning a matrix. 
#' 
#' @param v a vector of data.
#' @param time  an optional vector of the timestamps of \code{v}. If given, must be
#'  the same length as \code{v}. If not given, we try to infer it by summing the
#'  \code{time_deltas}.
#' @param time_deltas  an optional vector of the deltas of timestamps. If given, must be
#'  the same length as \code{v}. If not given, and \code{wts} are given and \code{wts_as_delta} is true,
#'  we take the \code{wts} as the time deltas.  The deltas must be positive. We sum them to arrive
#'  at the times.
#' @param window the window size, in time units. if given as finite integer or double, passed through.
#' If \code{NULL}, \code{NA_integer_}, \code{NA_real_} or \code{Inf} are given, 
#'  and \code{variable_win} is true, then we infer the window from the lookback times: the
#'  first window is infinite, but the remaining is the deltas between lookback times.
#'  If \code{variable_win} is false, then these undefined values are equivalent to an
#'  infinite window.
#'  If negative, an error will be thrown.
#' @param lb_time  a vector of the times from which lookback will be performed. The output should
#'  be the same size as this vector. If not given, defaults to \code{time}.
#' @param na_rm whether to remove NA, false by default.
#' @param max_order the maximum order of the centered moment to be computed.
#' @param min_df the minimum df to return a value, otherwise \code{NaN} is returned.
#' This can be used to prevent moments from being computed on too few observations.
#' Defaults to zero, meaning no restriction.
#' @param used_df the number of degrees of freedom consumed, used in the denominator
#' of the centered moments computation. These are subtracted from the number of
#' observations.
#' @param variable_win  if true, and the \code{window} is not a concrete number,
#'  the computation window becomes the time between lookback times.
#' @param wts_as_delta  if true and the \code{time} and \code{time_deltas} are not
#' given, but \code{wts} are given, we take \code{wts} as the \code{time_deltas}.
#'
#' @details
#'
#' Computes the number of elements, the mean, and the 2nd through kth
#' centered (and typically standardized) moments, for \eqn{k=2,3,4}{k=2,3,4}. These
#' are computed via the numerically robust one-pass method of Bennett \emph{et. al.}
#'
#' Given the length \eqn{n} vector \eqn{x}, we output matrix \eqn{M} where
#' \eqn{M_{i,j}}{M_i,j} is the \eqn{order - j + 1} moment (\emph{i.e.}
#' excess kurtosis, skewness, standard deviation, mean or number of elements)
#' of some elements \eqn{x_i} defined by the sliding time window.
#' Barring \code{NA} or \code{NaN}, this is over a window of time width \code{window}.
#'
#' @template sec-t-win
#'
#' @return Typically a matrix, where the first columns are the kth, k-1th through 2nd standardized, 
#' centered moments, then a column of the mean, then a column of the number of (non-nan) elements in the input,
#' with the following exceptions:
#' \describe{
#' \item{t_running_cent_moments}{Computes arbitrary order centered moments. When \code{max_order_only} is set,
#' only a column of the maximum order centered moment is returned.}
#' \item{t_running_std_moments}{Computes arbitrary order standardized moments, then the standard deviation, the mean,
#' and the count. There is not yet an option for \code{max_order_only}, but probably should be.}
#' \item{t_running_cumulants}{Computes arbitrary order cumulants, and returns the kth, k-1th, through the second 
#' (which is the variance) cumulant, then the mean, and the count.}
#' }
#'
#' @note
#' the kurtosis is \emph{excess kurtosis}, with a 3 subtracted, and should be nearly zero
#' for Gaussian input.
#'
#' @examples
#' x <- rnorm(1e5)
#' xs3 <- t_running_sd3(x,time=seq_along(x),window=10)
#' xs4 <- t_running_skew4(x,time=seq_along(x),window=10)
#' # but what if you only cared about some middle values?
#' xs4 <- t_running_skew4(x,time=seq_along(x),lb_time=(length(x) / 2) + 0:10,window=20)
#'
#' @inheritParams running_sd3
#' @template etc
#' @template ref-romo
#' @template param-wts
#' @template note-wts
#' @template note-heywood
#' @seealso \code{\link{running_sd3}}.
#' @rdname t_runningmoments
#' @export
t_running_sd3 <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_sd3', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningmoments
#' @export
t_running_skew4 <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_skew4', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningmoments
#' @export
t_running_kurt5 <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_kurt5', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningmoments
#' @export
t_running_sd <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_sd', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningmoments
#' @export
t_running_skew <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_skew', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningmoments
#' @export
t_running_kurt <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_kurt', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @param max_order_only for \code{running_cent_moments}, if this flag is set, only compute
#' the maximum order centered moment, and return in a vector.
#' @rdname t_runningmoments
#' @export
t_running_cent_moments <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, max_order = 5L, na_rm = FALSE, max_order_only = FALSE, min_df = 0L, used_df = 0.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_cent_moments', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, max_order, na_rm, max_order_only, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningmoments
#' @export
t_running_std_moments <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, max_order = 5L, na_rm = FALSE, min_df = 0L, used_df = 0.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_std_moments', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, max_order, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningmoments
#' @export
t_running_cumulants <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, max_order = 5L, na_rm = FALSE, min_df = 0L, used_df = 0.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_cumulants', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, max_order, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @title
#' Compute approximate quantiles over a sliding time window
#' @description
#' Computes cumulants up to some given order, then employs the Cornish-Fisher approximation
#' to compute approximate quantiles using a Gaussian basis.
#' 
#' @param p the probability points at which to compute the quantiles. Should be in the range (0,1).
#' @inheritParams t_running_cumulants
#'
#' @details
#'
#' Computes the cumulants, then approximates quantiles using AS269 of Lee & Lin.
#' @template sec-t-win
#'
#' @return A matrix, with one row for each element of \code{x}, and one column for each element of \code{q}.
#'
#' @note
#' The current implementation is not as space-efficient as it could be, as it first computes
#' the cumulants for each row, then performs the Cornish-Fisher approximation on a row-by-row
#' basis. In the future, this computation may be moved earlier into the pipeline to be more
#' space efficient. File an issue if the memory footprint is an issue for you.
#'
#' @examples
#' x <- rnorm(1e5)
#' xq <- t_running_apx_quantiles(x,c(0.1,0.25,0.5,0.75,0.9),
#'        time=seq_along(x),window=200,lb_time=c(100,200,400))
#'
#' xq <- t_running_apx_median(x,time=seq_along(x),window=200,lb_time=c(100,200,400))
#' xq <- t_running_apx_median(x,time=cumsum(runif(length(x),min=0.5,max=1.5)),
#'       window=200,lb_time=c(100,200,400))
#'
#' # weighted median?
#' wts <- runif(length(x),min=1,max=5)
#' xq <- t_running_apx_median(x,wts=wts,wts_as_delta=TRUE,window=1000,lb_time=seq(1000,10000,by=1000))
#'
#' # these should give the same answer:
#' xr <- running_apx_median(x,window=200);
#' xt <- t_running_apx_median(x,time=seq_along(x),window=199.99)
#'
#' @seealso \code{\link{running_apx_quantiles}}, \code{\link{t_running_cumulants}}, \code{PDQutils::qapx_cf}, \code{PDQutils::AS269}.
#' @template etc
#' @template ref-cf
#' @template ref-romo
#' @template param-wts
#' @template note-wts
#' @rdname t_runningquantiles
#' @export
t_running_apx_quantiles <- function(v, p, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, max_order = 5L, na_rm = FALSE, min_df = 0L, used_df = 0.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_apx_quantiles', PACKAGE = 'fromo', v, p, time, time_deltas, window, wts, lb_time, max_order, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningquantiles
#' @export
t_running_apx_median <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, max_order = 5L, na_rm = FALSE, min_df = 0L, used_df = 0.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_apx_median', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, max_order, na_rm, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @title
#' Compare data to moments computed over a time sliding window.
#' @description
#' Computes moments over a sliding window, then adjusts the data accordingly, centering, or scaling,
#' or z-scoring, and so on.
#' 
#' @inheritParams t_running_cent_moments
#' @param min_df the minimum df to return a value, otherwise \code{NaN} is returned.
#' This can be used to prevent \emph{e.g.} Z-scores from being computed on only 3
#' observations. Defaults to zero, meaning no restriction, which can result in 
#' infinite Z-scores during the burn-in period.
#' @param lookahead for some of the operations, the value is compared to 
#' mean and standard deviation possibly using 'future' or 'past' information
#' by means of a non-zero lookahead. Positive values mean data are taken from
#' the future. This is in time units, and so should be a real.
#' @param compute_se for \code{running_sharpe}, return an extra column of the
#' standard error, as computed by Mertens' correction.
#'
#' @details
#'
#' Given the length \eqn{n} vector \eqn{x}, for
#' a given index \eqn{i}, define \eqn{x^{(i)}}{x^(i)}
#' as the elements of \eqn{x} defined by the sliding time window (see the section
#' on time windowing).
#' Then define \eqn{\mu_i}{mu_i}, \eqn{\sigma_i}{sigma_i}
#' and \eqn{n_i}{n_i} as, respectively, the sample mean, standard deviation and number of
#' non-NA elements in \eqn{x^{(i)}}{x^(i)}. 
#'
#' We compute output vector \eqn{m} the same size as \eqn{x}. 
#' For the 'centered' version of \eqn{x}, we have \eqn{m_i = x_i - \mu_i}{m_i = x_i - mu_i}.
#' For the 'scaled' version of \eqn{x}, we have \eqn{m_i = x_i / \sigma_i}{m_i = x_i / sigma_i}.
#' For the 'z-scored' version of \eqn{x}, we have \eqn{m_i = (x_i - \mu_i) / \sigma_i}{m_i = (x_i - mu_i) / sigma_i}.
#' For the 't-scored' version of \eqn{x}, we have \eqn{m_i = \sqrt{n_i} \mu_i / \sigma_i}{m_i = sqrt(n_i) mu_i / sigma_i}.
#'
#' We also allow a 'lookahead' for some of these operations.
#' If positive, the moments are computed using data from larger indices;
#' if negative, from smaller indices. 
#'
#' @template sec-t-win
#'
#' @return a vector the same size as the input consisting of the adjusted version of the input.
#' When there are not sufficient (non-nan) elements for the computation, \code{NaN} are returned.
#'
#' @seealso \code{\link{running_centered}}, \code{\link{scale}}
#' @template etc
#' @template ref-romo
#' @template note-wts
#' @rdname t_runningadjustments
#' @export
t_running_centered <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, lookahead = 0.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_centered', PACKAGE = 'fromo', v, time, time_deltas, window, wts, na_rm, min_df, used_df, lookahead, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningadjustments
#' @export
t_running_scaled <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, lookahead = 0.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_scaled', PACKAGE = 'fromo', v, time, time_deltas, window, wts, na_rm, min_df, used_df, lookahead, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningadjustments
#' @export
t_running_zscored <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, na_rm = FALSE, min_df = 0L, used_df = 1.0, lookahead = 0.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_zscored', PACKAGE = 'fromo', v, time, time_deltas, window, wts, na_rm, min_df, used_df, lookahead, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningadjustments
#' @export
t_running_sharpe <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, compute_se = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_sharpe', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, compute_se, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @rdname t_runningadjustments
#' @export
t_running_tstat <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, compute_se = FALSE, min_df = 0L, used_df = 1.0, restart_period = 100L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE, normalize_wts = TRUE) {
    .Call('_fromo_t_running_tstat', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, compute_se, min_df, used_df, restart_period, variable_win, wts_as_delta, check_wts, normalize_wts)
}

#' @title
#' Compute sums or means over a sliding time window.
#'
#' @description
#' Compute the mean or sum over 
#' an infinite or finite sliding time window, returning a vector the same size as the lookback
#' times.
#' 
#' @param v a vector.
#' @inheritParams t_running_sd
#'
#' @param restart_period the recompute period. because subtraction of elements can cause
#' loss of precision, the computation of moments is restarted periodically based on 
#' this parameter. Larger values mean fewer restarts and faster, though potentially less 
#' accurate results. Unlike in the computation of even order moments, loss of precision
#' is unlikely to be disastrous, so the default value is rather large.
#' @param min_df the minimum df to return a value, otherwise \code{NaN} is returned,
#' only for the means computation.
#' This can be used to prevent moments from being computed on too few observations.
#' Defaults to zero, meaning no restriction.
#'
#' @details
#'
#' Computes the mean or sum of the elements, using a Kahan's Compensated Summation Algorithm,
#' a numerically robust one-pass method.
#'
#' Given the length \eqn{n} vector \eqn{x}, we output matrix \eqn{M} where
#' \eqn{M_{i,1}}{M_i,1} is the sum or mean 
#' of some elements \eqn{x_i} defined by the sliding time window.
#' Barring \code{NA} or \code{NaN}, this is over a window of time width \code{window}.
#'
#' @template sec-t-win
#'
#' @return A vector the same size as the lookback times.
#' @examples
#' x <- rnorm(1e5)
#' xs <- t_running_sum(x,time=seq_along(x),window=10)
#' xm <- t_running_mean(x,time=cumsum(runif(length(x))),window=7.3)
#'
#' @template etc
#' @template ref-romo
#' @template ref-kahan
#' @inheritParams sd3
#' @template note-wts
#' @rdname t_runningmean 
#' @export
t_running_sum <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, min_df = 0L, restart_period = 10000L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE) {
    .Call('_fromo_t_running_sum', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, min_df, restart_period, variable_win, wts_as_delta, check_wts)
}

#' @rdname t_runningmean
#' @export
t_running_mean <- function(v, time = NULL, time_deltas = NULL, window = NULL, wts = NULL, lb_time = NULL, na_rm = FALSE, min_df = 0L, restart_period = 10000L, variable_win = FALSE, wts_as_delta = TRUE, check_wts = FALSE) {
    .Call('_fromo_t_running_mean', PACKAGE = 'fromo', v, time, time_deltas, window, wts, lb_time, na_rm, min_df, restart_period, variable_win, wts_as_delta, check_wts)
}