# BCPA/PlotBCPA

(Difference between revisions)
Jump to: navigation, search
Eli (Talk | contribs)
(Created page with "<pre> PlotBCPA <- function (t, x, ws, pp, threshold=10) { plot(t, x, type = "n", main = "", ylim = c(min(x), max(x) * 1.3)) x.breaks <- ws\$Break x.model <- ...")
Newer edit →

## Revision as of 09:57, 19 July 2011

```PlotBCPA <- function (t, x, ws, pp, threshold=10)
{
plot(t, x, type = "n", main = "", ylim = c(min(x), max(x) *
1.3))
x.breaks <- ws\$Break
x.model <- ws\$Model
x.breaks <- x.breaks[x.model > 0]
x.model <- x.model[x.model > 0]
mids <- hist(x.breaks, breaks = t, plot = F)\$mid
freq <- hist(x.breaks, breaks = t, plot = F)\$count
goodbreaks <- mids[freq > threshold]
freq <- freq[freq > threshold]
break.cols <- heat.colors(max(freq))
break.cols <- break.cols[length(break.cols):1]
abline(v = goodbreaks, lwd = freq, col = "orange")
lines(t, x, col = "darkgrey")
lines(t, pp\$mu.hat, lwd = 2)
lines(t, pp\$mu.hat + pp\$s.hat, col = 2, lwd = 1.5)
lines(t, pp\$mu.hat - pp\$s.hat, col = 2, lwd = 1.5)
rho.hat <- pp\$rho.hat
rho.int <- round((rho.hat - min(rho.hat, na.rm = 1)) * 1000) +
1
rho.col <- topo.colors(max(rho.int, na.rm=TRUE))[rho.int]
rho.cex <- (rho.hat - min(rho.hat, na.rm=TRUE))/(diff(range(rho.hat))) *
1.5 + 0.5
points(t, x, col = rho.col, pch = 19)
rho.index <- quantile(0:max(rho.int))
rho.index[1] <- 1
legend.cols <- topo.colors(max(rho.int, na.rm=TRUE))[round(rho.index)]
legend.rhos <- round(quantile((min(rho.hat) * 1000):(max(rho.hat) *
1000))/1000, 2)
legend("bottomright", bg = "white", legend = c(expression(hat(rho)),
legend.rhos), pch = c(0, rep(19, 5)), ncol = 2, col = c(0,
legend.cols), xjust = 0.5, yjust = 0.5, cex = 1.2)
legend("topright", bg = "white", legend = c(expression(hat(mu)),
expression(hat(mu) %+-% hat(sigma))), lty = 1, lwd = 2:1,
col = 1:2, xjust = 0.5, yjust = 0.5, cex = 1.2)
}
```