cumSeg (version 1.3)

jumpoints: Change point estimation in piecewise constant models

Description

Estimation of change points and model selection via generalized BIC and other criteria

Usage

jumpoints(y, x, k = min(30, round(length(y)/10)), output = "2", 
    psi = NULL, round = TRUE, control = fit.control(), 
    selection = sel.control(), ...)

Arguments

y

the observed (genomic) sequence supposed to have a piecewise constant mean function.

x

the `segmented' variable, e.g. the genomic location. If missing simple indices 1,2,... n (length of y) are assumed.

k

the starting number of changepoints. It should be quite larger than the supposed number of (true) changepoints. This argument is ignored if starting values of the changepoints are specified via psi.

output

which output should be produced? Possible values are "1", "2", or "3"; see Details

psi

numeric vector to indicate the starting values for the changepoints. When psi=NULL (default), k quantiles are assumed.

round

logical; should the values of the changepoints be rounded?

control

a list returned by fit.control.

selection

a list returned by sel.control.

…

additional arguments.

Value

A list including several components depending on the value of output

If output="1" the most relevant components are

fitted.values

the fitted values

n.psi

the estimated number of changepoints

est.means

the estimated means

psi

the estimated changepoints

If output="2" the most relevant components are

fitted.values

the fitted values

n.psi

the estimated number of changepoints

criterion

the values of the selection criterion

psi

the estimated changepoints

est.means

the estimated means

psi0

the estimated changepoints at ouput 1 (before applying the selection criterion)

est.means0

the estimated means at ouput 1 (before applying the selection criterion)

If output="3" the most relevant components are those of output 2 but

psi0

the estimated changepoints at ouput 1

psi1

the estimated changepoints at ouput 2

psi

the estimated changepoints at ouput 3 (after applying again the segmented algorithm).

Details

The algorithm works by suitably transforming the observed responses to fit a continuous piecewise linear model. Starting from k changepoints, a large number of changepoints is first estimated. This number will be (typically slightly) lower than k since some changepoints will be discarded during the iterative steps when taking non admissible values. If output="1", jumpoints returns them which typically will be more than the actual ones. If output="2" the appropriate number of changepoints is selected via the criterion specified in argument selection via sel.control (e.g. BIC, MDL, ..). Finally if output="3", the segmented algorithm is run again to try to improve the changepoint estimates returned by the previous step.

References

Muggeo, V.M.R., Adelfio, G., Efficient change point detection for genomic sequences of continuous measurements, Bioinformatics 27, 161-166.

See Also

lars, sel.control, fit.control.

Examples

# NOT RUN {
# }
# NOT RUN {
n<-100
x<-1:n/n
 
lp<-I(x>.1) -1*I(x>.15)+.585*I(x>.45)-.585*I(x>.6) -I(x>.9)
e<-rnorm(n,0,.154)
y<-lp+e #data

#fit the model without selecting the changepoints
o1<-jumpoints(y,output="1")
plot(o1, typeL="l")
lines(lp, col=2) #true regression function
legend("topright", c("true","fit with output=1"),bty="n", col=c(2,1), lty=1)

#fit model and select the changepoints
o2<-jumpoints(y,output="2")
par(mfrow=c(1,2))
plot(o2, what="c")
plot(o2, typeL="s")
lines(lp, col=3) #true regression function
legend("topright", c("true","fit with output=2"),bty="n", col=c(3,1), lty=1)
# }
# NOT RUN {
    
# }