Pilliat_test_calibrate {HDCD} | R Documentation |
Generates detection thresholds for the Pilliat algorithm for testing for a single change-point using Monte Carlo simulation
Description
R wrapper for function choosing detection thresholds for the Dense, Partial sum and Berk-Jones statistics in the multiple change-point detection algorithm of Pilliat et al. (2023) for single change-point testing using Monte Carlo simulation. When Bonferroni==TRUE
, the detection thresholds are chosen by simulating the leading constant in the theoretical detection thresholds given in Pilliat et al. (2023), similarly as described in Appendix B in Moen et al. (2023) for ESAC. When Bonferroni==TRUE
, the thresholds for the Berk-Jones statistic are theoretical and not chosen by Monte Carlo simulation.
Usage
Pilliat_test_calibrate(
n,
p,
N = 100,
tol = 1/100,
threshold_bj_const = 6,
bonferroni = TRUE,
rescale_variance = TRUE,
fast = FALSE,
debug = FALSE
)
Arguments
n |
Number of observations |
p |
Number time series |
N |
Number of Monte Carlo samples used |
tol |
False error probability tolerance |
threshold_bj_const |
Leading constant for |
bonferroni |
If |
rescale_variance |
If |
fast |
If |
debug |
If |
Value
A list containing
thresholds_partial |
vector of thresholds for the Partial Sum statistic (respectively for |
threshold_dense |
threshold for the dense statistic |
thresholds_bj |
vector of thresholds for the Berk-Jones static (respectively for |
Examples
library(HDCD)
n = 50
p = 50
set.seed(100)
thresholds_test_emp = Pilliat_test_calibrate(n,p, bonferroni=TRUE,N=100, tol=1/100)
set.seed(100)
thresholds_test_emp_without_bonferroni = Pilliat_test_calibrate(n,p,
bonferroni=FALSE,N=100, tol=1/100)
thresholds_test_emp # thresholds with bonferroni correction
thresholds_test_emp_without_bonferroni # thresholds without bonferroni correction
# Generating data
X = matrix(rnorm(n*p), ncol = n, nrow=p)
Y = matrix(rnorm(n*p), ncol = n, nrow=p)
# Adding a single sparse change-point to X (and not Y):
X[1:5, 25:50] = X[1:5, 25:50] +2
resX = Pilliat_test(X,
threshold_dense=thresholds_test_emp$threshold_dense,
thresholds_bj = thresholds_test_emp$thresholds_bj,
thresholds_partial = thresholds_test_emp$thresholds_partial
)
resX
resY = Pilliat_test(Y,
threshold_dense=thresholds_test_emp$threshold_dense,
thresholds_bj = thresholds_test_emp$thresholds_bj,
thresholds_partial = thresholds_test_emp$thresholds_partial
)
resY