| Title: | Multiple Testing of Local Extrema for Detection of Change Points |
|---|---|
| Description: | A new approach to detect change points based on smoothing and multiple testing, which is for long data sequence modeled as piecewise constant functions plus stationary Gaussian noise, see Dan Cheng and Armin Schwartzman (2015) <arXiv:1504.06384>. |
| Authors: | Zhibing He and Dan Cheng |
| Maintainer: | Zhibing He <[email protected]> |
| License: | GPL-3 |
| Version: | 1.0-1 |
| Built: | 2024-12-13 06:30:57 UTC |
| Source: | CRAN |
Estimate s2,lambda2,lambda4,Delta
ch.est(nu, gamma, size, B = 100)ch.est(nu, gamma, size, B = 100)
nu |
bandwidth of Gaussian kernal applied to White-noise, Whitenoise error if |
gamma |
bandwidth of nonparameter smoothing |
size |
sample size |
B |
Montelarlo iteration times |
a list of s2,lambda2,lambda4,Delta
Multiple Testing of Local Extrema for Detection of Change Points https://arxiv.org/abs/1504.06384
ch.est(nu=2,gamma=4,size=1000,B=100)ch.est(nu=2,gamma=4,size=1000,B=100)
Compute convolution function using FFT, similar to the function 'conv' in matlab
conv(u, v, shape = c("same", "full"))conv(u, v, shape = c("same", "full"))
u |
vector |
v |
vector |
shape |
if 'same', return central part of the convolution,the same size as u; ortherwise return the whole sequence with size lenth(u)+length(v)-1 |
a vector of convolution, as specified by shape.
Matlab document on 'conv' https://www.mathworks.com/help/matlab/ref/conv.html
u = c(-1,2,3,-2,0,1,2) v = c(2,4,-1,1) w = conv(u,v,'same')u = c(-1,2,3,-2,0,1,2) v = c(2,4,-1,1) w = conv(u,v,'same')
Evaluate performance of estimated change points
Fdr(uh, b, th)Fdr(uh, b, th)
uh |
a vector of estimated change points locations |
b |
a scalar of location tolerance, specified by user |
th |
a vector of true change points locations |
a list of vector of FDR and Power
FDR |
a scalar of fdr (false discovery rate) |
Power |
a scalar of power (true positive rate) |
Fdr(uh=c(7,15,32,47),b=4,th=c(10,20,30,40,50))Fdr(uh=c(7,15,32,47),b=4,th=c(10,20,30,40,50))
gamma and nu
Parallel computing fdr and power of change points estimation for different gamma and nu
fdr.gam(c, mu, Gamma, Nu, b, th, B = 100, level = 0.1, iter = 100)fdr.gam(c, mu, Gamma, Nu, b, th, B = 100, level = 0.1, iter = 100)
c |
number of cpu cores used for parallel computing |
mu |
a vector of piecewise constant |
Gamma |
a vector of different |
Nu |
a vector of different |
b |
a scalar of location tolerance, specified by user |
th |
a vector of true change points locations |
B |
Montelarlo iteration times |
level |
FDR control level |
iter |
iteration times for each combination of |
a list of matrix with the same length as Nu, FDR and Power for different Gamma are displayed within each matrix
size=12000 a = 1 A = a*(1:119) H = seq(100,11900,100) mu = GenMu(A,H,size=size) z = GenZ(nu=2,size=size) Gamma = seq(1,5,1) Nu = seq(0,2,0.5) model = fdr.gam(2,mu,Gamma,Nu,8,H,iter=100)size=12000 a = 1 A = a*(1:119) H = seq(100,11900,100) mu = GenMu(A,H,size=size) z = GenZ(nu=2,size=size) Gamma = seq(1,5,1) Nu = seq(0,2,0.5) model = fdr.gam(2,mu,Gamma,Nu,8,H,iter=100)
FDR threshold based on the Benjamini-Hochberg algorithm
fdrBH(p, q)fdrBH(p, q)
p |
a vector of p-values |
q |
False Discovery Rate level |
p-value threshold based on independence or positive dependence
fdrBH(seq(0.01,0.1,0.01),q=0.1)fdrBH(seq(0.01,0.1,0.01),q=0.1)
Generate first-order differential of a smoothed sequence Y
GenDY(mu, z, gamma)GenDY(mu, z, gamma)
mu |
a vector of piecewise constant |
z |
a vector of stationary Gaussian random error |
gamma |
bandwidth of nonparameter smoothing |
a vector of the differential of Y
mu = GenMu(x=1:10,pos=seq(10,100,10),size=150) z = GenZ(nu=2,size=150) GenDY(mu=mu,z=z,gamma=4)mu = GenMu(x=1:10,pos=seq(10,100,10),size=150) z = GenZ(nu=2,size=150) GenDY(mu=mu,z=z,gamma=4)
Generate a piecewise constant sequence starting from 0
GenMu(x, pos, size)GenMu(x, pos, size)
x |
a vector containing all values of change points |
pos |
positions of change points, corresponding to |
size |
sample size |
a piecewise constant sequence
GenMu(x=1:10,pos=seq(10,100,10),size=150)GenMu(x=1:10,pos=seq(10,100,10),size=150)
Generate Gaussian autocorrelated random error sequence based on White-noise and Gaussian kernal
GenZ(nu, size)GenZ(nu, size)
nu |
bandwidth of Gaussian kernal applied to White-noise, Whitenoise error if |
size |
sample size |
a vector of random error
GenZ(nu=2,size=1000)GenZ(nu=2,size=1000)
Illustration plot of the procedure t0 detect change points
illu.plot(mu, z, gamma, whichcp, b, Tmax, Tmin)illu.plot(mu, z, gamma, whichcp, b, Tmax, Tmin)
mu |
a vector of piecewise constant |
z |
a vector of stationary Gaussian random error |
gamma |
bandwidth of nonparameter smoothing |
whichcp |
output of the function |
b |
a scalar of location tolerance, specified by user |
Tmax |
a vector of true peak locations |
Tmin |
a vector true valley locations |
a figure plot showing detection of change points
set.seed(2019) L = 1200 A = c(2.8,0,-2.4,0,-3,0.5,3,5,2,0)/1.5 Tmax = c(150,410,680,770,980) Tmin = c(250,320,550,1000,1100) H = c(150,250,320,410,550,680,770,980,1000,1100) mu = GenMu(A,H,L); z = GenZ(nu=2,L) y1 = GenDY(mu=mu,z=z,gamma=6) chest = ch.est(nu=2,gamma=6,size=L,B=100) chp= which.cp(y1,chest,level=0.1) illu.plot(mu,z,gamma=6,chp,b=5,Tmax,Tmin)set.seed(2019) L = 1200 A = c(2.8,0,-2.4,0,-3,0.5,3,5,2,0)/1.5 Tmax = c(150,410,680,770,980) Tmin = c(250,320,550,1000,1100) H = c(150,250,320,410,550,680,770,980,1000,1100) mu = GenMu(A,H,L); z = GenZ(nu=2,L) y1 = GenDY(mu=mu,z=z,gamma=6) chest = ch.est(nu=2,gamma=6,size=L,B=100) chp= which.cp(y1,chest,level=0.1) illu.plot(mu,z,gamma=6,chp,b=5,Tmax,Tmin)
Find locations of change points
which.cp(y1, chest, level = 0.1)which.cp(y1, chest, level = 0.1)
y1 |
a vector of the differential of sequence Y |
chest |
output of function |
level |
FDR control level |
a list of components
peak |
a vector of peaks location |
vall |
a vector of valleys location |
pval |
a scalar of adjusted p-value based on FDR control |
thresh |
a scalar of threshold for |
mu = GenMu(x=1:10,pos=seq(10,100,10),size=150) z = GenZ(nu=2,size=150) y1 = GenDY(mu,z,gamma=4) chest = ch.est(nu=2,gamma=8,size=150,B=100) which.cp(y1,chest,level=0.1)mu = GenMu(x=1:10,pos=seq(10,100,10),size=150) z = GenZ(nu=2,size=150) y1 = GenDY(mu,z,gamma=4) chest = ch.est(nu=2,gamma=8,size=150,B=100) which.cp(y1,chest,level=0.1)
Find local maxima and minima in a sequence
which.peaks(x, partial = FALSE, decreasing = FALSE)which.peaks(x, partial = FALSE, decreasing = FALSE)
x |
a vector with maxima or minima |
partial |
endpoints will be considered if 'true' |
decreasing |
find local minima if 'true', ortherwise local maxima |
a vector of positions of local maxima or minima
a = 100:1 which.peaks(a*sin(a/3))a = 100:1 which.peaks(a*sin(a/3))