plot.bcp {bcp} | R Documentation |
`plot.bcp' produces summary plots of the results of the `bcp' function.
plot.bcp(x, ...)
x |
the result of a call to `bcp'. |
... |
additional arguments. |
`plot.bcp' produces the following plots:
Posterior Means: location in the sequence versus the posterior mean over the iterations.
Posterior Probability of a Change: location in the sequence versus the relative frequency of iterations which resulted in a change point.
Chandra Erdman
Bai J., Perron P. (2003), Computation and Analysis of Multiple Structural Change Models, Journal of Applied Econometrics, 18, 1-22. url: http: //qed.econ.queensu.ca/jae/2003-v18.1/bai-perron/.
Daniel Barry and J. A. Hartigan (1993), A Bayesian Analysis for Change Point Problems, Journal of The American Statistical Association, 88, 309-19.
Olshen, A. B., Venkatraman, E. S., Lucito, R., Wigler, M. (2004), Circular binary segmentation for the analysis of array-based DNA copy number data, Biostatistics, 5, 557-572. url: http://www.bioconductor.org/repository/release1.5/package/html/DNAcopy.html.
Martyn Plummer, Nicky Best, Kate Cowles, and Karen Vines (2006), The coda Package, version 0.10-7, CRAN: The Comprehensive R Network.
Snijders et al. (2001), Assembly of microarrays for genome-wide measurement of DNA copy number, Nature Genetics, 29, 263-264.
Achim Zeileis, Friedrich Leisch, Bruce Hansen, Kurt Hornik, Christian Kleiber (2006), The strucchange Package, version 1.3-1, CRAN: The Comprehensive R Network.
`bcp', `summary.bcp', and `print.bcp' for complete results and summary statistics.
##### A random sample from a few normal distributions ##### testdata <- c(rnorm(20), rnorm(10, 5, 1), rnorm(20)) bcp.0 <- bcp(testdata) summary.bcp(bcp.0) plot.bcp(bcp.0) ##### Coriell chromosome 11 ##### data(coriell) chrom11 <- na.omit(coriell$Coriell.05296[coriell$Chromosome==11]) n <- length(chrom11) bcp.11 <- bcp(chrom11[1:n]) summary.bcp(bcp.11) plot.bcp(bcp.11) # to see bcp and Circular Binary Segmentation results run: if(require("DNAcopy")) { bcp.11$rhos[,ncol(bcp.11$rhos)] <- 0 changepoint.freq <- apply(bcp.11$rhos[bcp.11$burnin:(bcp.11$burnin+bcp.11$burnin),1:dim(bcp.11$results)[2]],2,mean) cbs <- segment(CNA(chrom11[1:n], rep(1, n), 1:n), verbose = 0) cbs.ests <- rep(unlist(cbs$output[6]), unlist(cbs$output[5])) op <- par(mfrow=c(2,1),col.lab="black",col.main="black") plot(1:n, bcp.11$posterior.mean, type="l", xlab="Location", ylab="Posterior Mean", main="Posterior Means") lines(cbs.ests, col="red") points(chrom11) plot(1:n, changepoint.freq, type="l", xlab="Location", ylab="Posterior Probability of a Change", main="Change Point Locations") for(i in 1:(dim(cbs$output)[1]-1)) abline(v=cbs$output$loc.end[i], col="red") par(op) } else { cat("DNAcopy is not loaded") } ##### RealInt ##### data("RealInt") bcp.ri <- bcp(as.vector(RealInt)) summary.bcp(bcp.ri) plot.bcp(bcp.ri) # to see bcp and Bai and Perron results run: if(require("strucchange")) { bcp.ri$rhos[,ncol(bcp.ri$rhos)] <- 0 changepoint.freq <- apply(bcp.ri$rhos[bcp.ri$burnin:(bcp.ri$burnin+bcp.ri$burnin),1:dim(bcp.ri$results)[2]],2,mean) bp <- breakpoints(RealInt ~ 1, h = 2)$breakpoints rho <- rep(0, length(RealInt)) rho[bp] <- 1 b.num<-1 + c(0,cumsum(rho[1:(length(rho)-1)])) bp.mean <- unlist(lapply(split(RealInt,b.num),mean)) bp.ri <- rep(0,length(RealInt)) for(i in 1:length(bp.ri)) bp.ri[i] <- bp.mean[b.num[i]] op <- par(mfrow=c(2,1),col.lab="black",col.main="black") xax <- seq(1961, 1987, length=103) plot(xax, bcp.ri$posterior.mean, type="l", xlab="Time", ylab="Mean", main="Posterior Means") lines(xax, bp.ri, col="blue") points(RealInt) plot(xax, changepoint.freq, type="l", xlab="Time", ylab="Posterior Probability", main="Posterior Probability of a Change") for(i in 1:length(bp.ri)) abline(v=xax[bp[i]], col="blue") par(op) } else { cat("strucchange is not loaded") }