computeStat: Computation of the multiscale statistic

Description

Computes the multiscale vector of penalised statistics, (3.7) in the vignette, or the penalised multiscale statistic, (3.6) in the vignette, for given signal.

Usage

computeStat(y, signal = 0, family = NULL, intervalSystem = NULL, lengths = NULL,
            penalty = NULL, nq = length(y),
            output = c("list", "vector", "maximum"), ...)

Value

If output == list a list containing in maximum the penalised multiscale statistic, i.e. the maximum over all test statistics, in stat the multiscale vector of penalised statistics, i.e. a vector of length lengths giving the maximum over all tests of that length, and in lengths the vector of lengths. If output == vector a numeric vector giving the multiscale vector of penalised statistics. If output == maximum a single numeric giving the penalised multiscale statistic. -Inf is returned for lengths for which on all intervals of that length contained in the set of intervals the signal is not constant and, hence, no test statistic can be computed. This behaves similar to max(numeric(0)).

Arguments

y: a numeric vector containing the observations
signal: the given signal, either a single numeric for a constant function equal to the given value or an object of class stepfit. More precisely, a list containing an integer vector leftIndex, an integer vector rightIndex and a numeric vector value, all of the same length, e.g. a data.frame, specifying a step function is enough
family: a string specifying the assumed parametric family, for more details see parametricFamily, currently "gauss", "hsmuce" and "mDependentPS" are supported. By default (NULL) "gauss" is assumed
intervalSystem: a string giving the used interval system, either "all" for all intervals, "dyaLen" for all intervals of dyadic length or "dyaPar" for the dyadic partition, for more details see intervalSystem. By default (NULL) the default interval system of the specified parametric family will be used, which one this will be is described in parametricFamily
lengths: an integer vector giving the set of lengths, i.e. only intervals of these lengths will be considered. Note that not all lengths are possible for all interval systems and for all parametric families, see intervalSystem and parametricFamily, respectively, to see which ones are allowed. By default (NULL) all lengths that are possible for the specified intervalSystem and for the specified parametric family will be used
penalty: a string specifying how the statistics will be penalised, either "sqrt", "log" or "none", see penalty and section 3.2 in the vignette for more details. By default (NULL) the default penalty of the specified parametric family will be used, which one this will be is described in parametricFamily
nq: a single integer larger than or equal to length(y) giving the number of observations used in the penalty term, see penalty for more details. The possibility to use a number larger than length(y) is given for comparisons, since a (slightly) larger number can be chosen in critVal and monteCarloSimulation to avoid extensive recomputations for (slightly) varying number of observations. For more details see also the Section Simulating, saving and loading of Monte-Carlo simulations in critVal
output: a string specifying the output, see Value
...: further parameters of the parametric family. Depending on argument family some might be required, but others might be optional, please see parametricFamily for more details

References

Frick, K., Munk, A., Sieling, H. (2014) Multiscale change-point inference. With discussion and rejoinder by the authors. Journal of the Royal Statistical Society, Series B 76(3), 495--580.

Pein, F., Sieling, H., Munk, A. (2017) Heterogeneous change point inference. Journal of the Royal Statistical Society, Series B, 79(4), 1207--1227.

Examples

Run this code

# \dontshow{
savePathRcache <- R.cache::getCacheRootPath()

R.cache::setCacheRootPath(path = file.path(R.cache::getCacheRootPath(), "test"))
# }

y <- rnorm(100)
# for the default signal = 0 a signal constant 0 is assumed
identical(computeStat(y), computeStat(y,
            signal = list(leftIndex = 1L, rightIndex = 100L, value = 0)))

# different constant value
ret <- computeStat(y, signal = 1)
# penalised multiscale statistic
identical(ret$maximum, computeStat(y, signal = 1, output = "maximum"))
# multiscale vector of penalised statistics
identical(ret$stat, computeStat(y, signal = 1, output = "vector"))

y <- c(rnorm(50), rnorm(50, 1))
# true signal
computeStat(y, signal = list(leftIndex = c(1L, 51L), rightIndex = c(50L, 100L),
                             value = c(0, 1)))

# fit satisfies the multiscale contraint, i.e.
# the penalised multiscale statistic is not larger than the used global quantile 1
computeStat(y, signal = stepFit(y, q = 1), output = "maximum") <= 1

# different interval system, lengths, penalty, given parameter sd
# and computed for an increased number of observations nq
computeStat(y, signal = list(leftIndex = c(1L, 51L), rightIndex = c(50L, 100L),
                             value = c(0, 1)), nq = 128, sd = 0.5,
            intervalSystem = "dyaLen", lengths = c(1L, 2L, 4L, 8L), penalty = "none")

# family "hsmuce"
computeStat(y, signal = mean(y), family = "hsmuce")

# family "mDependentPS"
signal <- list(leftIndex = c(1L, 13L), rightIndex = c(12L, 17L), value = c(0, -1))
y <- c(rep(0, 13), rep(-1, 4)) + 
  as.numeric(arima.sim(n = 17, list(ar = c(), ma = c(0.8, 0.5, 0.3)), sd = 1))
covariances <- as.numeric(ARMAacf(ar = c(), ma = c(0.8, 0.5, 0.3), lag.max = 3))
computeStat(y, signal = signal, family = "mDependentPS", covariances = covariances)

# \dontshow{
unlink(R.cache::getCacheRootPath(), force = TRUE, recursive = TRUE)

R.cache::setCacheRootPath(savePathRcache)
# }