[r-cran-coda] 07/60: Imported Upstream version 0.10-2
Andreas Tille
tille at debian.org
Fri Dec 16 12:11:22 UTC 2016
This is an automated email from the git hooks/post-receive script.
tille pushed a commit to branch master
in repository r-cran-coda.
commit 8bba452b893a3c4367db7e59beab216e77be3996
Author: Andreas Tille <tille at debian.org>
Date: Fri Dec 16 12:09:54 2016 +0100
Imported Upstream version 0.10-2
---
CHANGELOG | 14 ++
DESCRIPTION | 12 +-
INDEX | 13 +-
NAMESPACE | 24 +++
R.Rcheck/00check.log | 3 -
R/HPDinterval.R | 20 ++
R/cumuplot.R | 9 +-
R/gelman.R | 2 +-
R/geweke.R | 2 +-
R/output.R | 20 +-
R/trellisplots.R | 546 +++++++++++++++++++++++++++++++++++++++++++++++++++
man/HPDinterval.Rd | 44 +++++
man/autocorr.diag.Rd | 2 +-
man/autocorr.plot.Rd | 3 +-
man/cumuplot.Rd | 4 +-
man/gelman.plot.Rd | 4 +-
man/geweke.plot.Rd | 2 +-
man/plot.mcmc.Rd | 2 +-
man/trellisplots.Rd | 195 ++++++++++++++++++
19 files changed, 881 insertions(+), 40 deletions(-)
diff --git a/CHANGELOG b/CHANGELOG
index 22775f3..ac6facb 100644
--- a/CHANGELOG
+++ b/CHANGELOG
@@ -1,3 +1,17 @@
+0.10-2
+- Import generics from lattice
+- In plotting functions, "ask" now defaults to dev.interactive(),
+ and not the default value in par(), as introduced in 0.9-2
+
+0.10-1
+
+- Added "mcmc" and "mcmc.list" methods for several lattice functions
+ (xyplot, qqmath, densityplot and acfplot (with generic here as
+ well)). There's a levelplot method too, which is currently
+ experimental. These are Trellis analogs of existing coda functions,
+ and may in future replace them.
+
+
0.9-5
- Fixed bug in summary.mcmc (safespec0 misspelled) which affected
diff --git a/DESCRIPTION b/DESCRIPTION
index f53fe1e..5569570 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,11 +1,11 @@
Package: coda
-Version: 0.9-5
-Date: 2005/5/4
+Version: 0.10-2
+Date: 2005-11-10
Title: Output analysis and diagnostics for MCMC
Author: Martyn Plummer, Nicky Best, Kate Cowles, Karen Vines
Maintainer: Martyn Plummer <plummer at iarc.fr>
-Depends: R (>= 1.9)
-Description: Output analysis and diagnostics for Markov Chain Monte Carlo simulations.
+Depends: R (>= 2.2.0), lattice
+Description: Output analysis and diagnostics for Markov Chain Monte
+ Carlo simulations.
License: GPL Version 2 or later.
-URL: http://www-fis.iarc.fr/coda/
-Packaged: Wed May 4 14:51:50 2005; martyn
+Packaged: Thu Nov 10 18:19:02 2005; martyn
diff --git a/INDEX b/INDEX
index 377df45..0c61bd9 100644
--- a/INDEX
+++ b/INDEX
@@ -10,15 +10,16 @@ coda.options Options settings for the codamenu driver
crosscorr Cross correlations for MCMC output
crosscorr.plot Plot image of correlation matrix
cumuplot Cumulative quantile plot
-densplot Probability density function estimate from
- MCMC output
+densityplot.mcmc Trellis plots for mcmc objects
+densplot Probability density function estimate from MCMC
+ output
effectiveSize Effective sample size for estimating the mean
gelman.diag Gelman and Rubin's convergence diagnostic
gelman.plot Gelman-Rubin-Brooks plot
geweke.diag Geweke's convergence diagnostic
geweke.plot Geweke-Brooks plot
-heidel.diag Heidelberger and Welch's convergence
- diagnostic
+heidel.diag Heidelberger and Welch's convergence diagnostic
+HPDinterval Highest Posterior Density intervals
line Simple linear regression example
mcmc Markov Chain Monte Carlo Objects
[.mcmc Extract or replace parts of MCMC objects
@@ -38,8 +39,8 @@ read.openbugs Read CODA output files produced by OpenBUGS
rejectionRate Rejection Rate for Metropolis-Hastings chains
spectrum0 Estimate spectral density at zero
spectrum0.ar Estimate spectral density at zero
-summary.mcmc Summary statistics for Markov Chain Monte
- Carlo chains
+summary.mcmc Summary statistics for Markov Chain Monte Carlo
+ chains
thin Thinning interval
time.mcmc Time attributes for mcmc objects
traceplot Trace plot of MCMC output
diff --git a/NAMESPACE b/NAMESPACE
index 7dc871d..f7d687b 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -1,5 +1,7 @@
## export functions
export(
+ acfplot,
+ HPDinterval,
as.array.mcmc.list,
as.mcmc,
as.ts.mcmc,
@@ -82,6 +84,8 @@ S3method(autocorr.diag, mcmc)
S3method(rejectionRate, mcmc)
S3method(batchSE,mcmc)
S3method(as.data.frame,mcmc)
+##DMB
+S3method(HPDinterval, mcmc)
## mcmc.list
S3method("[", mcmc.list)
@@ -98,6 +102,8 @@ S3method(window, mcmc.list)
S3method(autocorr.diag, mcmc.list)
S3method(rejectionRate, mcmc.list)
S3method(batchSE,mcmc.list)
+##DMB
+S3method(HPDinterval, mcmc.list)
## misc
S3method(print, gelman.diag)
@@ -108,3 +114,21 @@ S3method(print, summary.mcmc)
## Imports
importFrom(stats, start, end, frequency, time, window)
+importFrom(lattice, xyplot, qqmath, densityplot, levelplot)
+
+
+## lattice methods
+
+S3method(xyplot, mcmc)
+S3method(qqmath, mcmc)
+S3method(densityplot, mcmc)
+S3method(acfplot, mcmc)
+S3method(levelplot, mcmc)
+
+S3method(xyplot, mcmc.list)
+S3method(qqmath, mcmc.list)
+S3method(densityplot, mcmc.list)
+S3method(acfplot, mcmc.list)
+
+
+
diff --git a/R.Rcheck/00check.log b/R.Rcheck/00check.log
deleted file mode 100644
index 4bc706b..0000000
--- a/R.Rcheck/00check.log
+++ /dev/null
@@ -1,3 +0,0 @@
-* using log directory '/home/martyn/bayes/diagnostics/coda/0.9/coda/R.Rcheck'
-* using R version 2.1.0, 2005-04-18
-* checking for file 'R/DESCRIPTION' ... NO
diff --git a/R/HPDinterval.R b/R/HPDinterval.R
new file mode 100644
index 0000000..eaab4a7
--- /dev/null
+++ b/R/HPDinterval.R
@@ -0,0 +1,20 @@
+HPDinterval <- function(obj, prob = 0.95, ...) UseMethod("HPDinterval")
+
+HPDinterval.mcmc <- function(obj, prob = 0.95, ...)
+{
+ vals <- apply(obj, 2, sort)
+ if (!is.matrix(vals)) stop("obj must have nsamp > 1")
+ nsamp <- nrow(vals)
+ npar <- ncol(vals)
+ gap <- max(1, min(nsamp - 1, round(nsamp * prob)))
+ init <- 1:(nsamp - gap)
+ inds <- apply(vals[init + gap, ] - vals[init, ], 2, which.min)
+ ans <- cbind(vals[cbind(inds, 1:npar)],
+ vals[cbind(inds + gap, 1:npar)])
+ dimnames(ans) <- list(colnames(obj), c("lower", "upper"))
+ attr(ans, "Probability") <- gap/nsamp
+ ans
+}
+
+HPDinterval.mcmc.list <- function(obj, prob = 0.95, ...)
+ lapply(obj, HPDinterval, prob)
diff --git a/R/cumuplot.R b/R/cumuplot.R
index 4b16ff3..2e55ef1 100644
--- a/R/cumuplot.R
+++ b/R/cumuplot.R
@@ -1,6 +1,6 @@
cumuplot <- function(x, probs=c(0.025,0.5,0.975), ylab="", lty=c(2,1),
- lwd=c(1,2), type="l", ask=TRUE, auto.layout=TRUE,
- col=1, ...)
+ lwd=c(1,2), type="l", ask=dev.interactive(),
+ auto.layout=TRUE, col=1, ...)
{
cquantile <- function(z, probs)
{
@@ -14,12 +14,11 @@ cumuplot <- function(x, probs=c(0.025,0.5,0.975), ylab="", lty=c(2,1),
return(cquant)
}
- oldpar <- par(ask = ask)
+ oldpar <- NULL
on.exit(par(oldpar))
if (auto.layout) {
oldpar <- par(mfrow = set.mfrow(Nchains = nchain(x),
Nparms = nvar(x)))
- oldpar <- c(oldpar, par(ask = ask))
}
if (!is.mcmc.list(x))
@@ -33,6 +32,8 @@ cumuplot <- function(x, probs=c(0.025,0.5,0.975), ylab="", lty=c(2,1),
col=col, ...)
title(paste(varnames(x)[j], ifelse(is.null(chanames(x)),
"", ":"), chanames(x)[i], sep = ""))
+ if (i == 1 & j == 1)
+ oldpar <- c(oldpar, par(ask=ask))
}
}
}
diff --git a/R/gelman.R b/R/gelman.R
index d3d91d0..f78c0fd 100644
--- a/R/gelman.R
+++ b/R/gelman.R
@@ -144,7 +144,7 @@
"gelman.plot" <-
function (x, bin.width = 10, max.bins = 50, confidence = 0.95,
- transform = FALSE, auto.layout = TRUE, ask = TRUE,
+ transform = FALSE, auto.layout = TRUE, ask = dev.interactive(),
col = 1:2, lty = 1:2, xlab = "last iteration in chain",
ylab = "shrink factor", type = "l", ...)
{
diff --git a/R/geweke.R b/R/geweke.R
index d5f52fb..2ef4680 100644
--- a/R/geweke.R
+++ b/R/geweke.R
@@ -21,7 +21,7 @@
"geweke.plot" <-
function (x, frac1 = 0.1, frac2 = 0.5, nbins = 20,
- pvalue = 0.05, auto.layout = TRUE, ask = TRUE, ...)
+ pvalue = 0.05, auto.layout = TRUE, ask = dev.interactive(), ...)
{
x <- as.mcmc.list(x)
oldpar <- NULL
diff --git a/R/output.R b/R/output.R
index ba11f2b..e3c6e85 100644
--- a/R/output.R
+++ b/R/output.R
@@ -21,14 +21,13 @@ function (x, lags = c(0, 1, 5, 10, 50), relative = TRUE)
}
"autocorr.plot" <-
-function (x, lag.max, auto.layout = TRUE, ask = par("ask"), ...)
+function (x, lag.max, auto.layout = TRUE, ask = dev.interactive(), ...)
{
-## RGA fixed to use default ask value.
- oldpar <- NULL
+ oldpar <- NULL
on.exit(par(oldpar))
if (auto.layout)
oldpar <- par(mfrow = set.mfrow(Nchains = nchain(x),
- Nparms = nvar(x)))
+ Nparms = nvar(x)))
if (!is.mcmc.list(x))
x <- mcmc.list(as.mcmc(x))
for (i in 1:nchain(x)) {
@@ -37,12 +36,12 @@ function (x, lag.max, auto.layout = TRUE, ask = par("ask"), ...)
else acf(as.ts.mcmc(x[[i]]), lag.max = lag.max, plot = FALSE)
for (j in 1:nvar(x)) {
plot(xacf$lag[, j, j], xacf$acf[, j, j], type = "h",
- ylab = "Autocorrelation", xlab = "Lag", ylim = c(-1,
- 1), ...)
- title(paste(varnames(x)[j], ifelse(is.null(chanames(x)),
- "", ":"), chanames(x)[i], sep = ""))
+ ylab = "Autocorrelation", xlab = "Lag", ylim = c(-1, 1), ...)
+ title(paste(varnames(x)[j],
+ ifelse(is.null(chanames(x)), "", ":"),
+ chanames(x)[i], sep = ""))
if (i==1 && j==1)
- oldpar <- c(oldpar, par(ask = ask))
+ oldpar <- c(oldpar, par(ask = ask))
}
}
invisible(x)
@@ -320,9 +319,8 @@ function (x, smooth = TRUE, col = 1:6, type = "l", ylab = "",
}
"plot.mcmc" <- function (x, trace = TRUE, density = TRUE, smooth = TRUE, bwf,
- auto.layout = TRUE, ask = par("ask"), ...)
+ auto.layout = TRUE, ask = dev.interactive(), ...)
{
-## RGA fixed to use default ask value.
oldpar <- NULL
on.exit(par(oldpar))
if (auto.layout) {
diff --git a/R/trellisplots.R b/R/trellisplots.R
new file mode 100644
index 0000000..c17a822
--- /dev/null
+++ b/R/trellisplots.R
@@ -0,0 +1,546 @@
+
+
+### Copyright (C) 2005 Deepayan Sarkar
+### <Deepayan.Sarkar at R-project.org>, Douglas Bates
+### <Douglas.Bates at R-project.org>. See file COPYING for license
+### terms.
+
+
+
+
+
+### unexported helper function to obtain a valid subset argument.
+### Mostly to ensure that we don't up with a huge vector of integer
+### indices in the trivial cases.
+
+
+
+thinned.indices <- function(object, n = nrow(object), start = 1, thin = 1)
+{
+ if (is.mcmc(object) &&
+ (start * thin != 1) &&
+ !all(mcpar(object)[-2] == 1))
+ warning("mcmc object is already thinned")
+ if (start < 1) stop("Invalid start")
+ else if (start == 1)
+ {
+ if (thin < 1) stop("Invalid thin")
+ else if (thin == 1) TRUE
+ else rep(c(TRUE, FALSE), c(1, thin-1))
+ }
+ else
+ {
+ if (thin < 1) stop("Invalid thin")
+ else if (thin == 1) -seq(length = start-1)
+ else start + thin * (0:(floor(n - start) / thin))
+ }
+}
+
+
+
+## most functions will have methods for mcmc and mcmclist objects. In
+## the second case, another grouping variable is added. By default,
+## this will be used for ``grouped displays'' (when possible), but
+## could also be used for conditioning by setting groups=FALSE
+
+
+
+
+
+## levelplot, analog of plot.crosscorr. Defaults changed to match
+## those of plot.crosscorr as much as possible.
+
+
+levelplot.mcmc <-
+ function(x, main = attr(x, "title"),
+ start = 1, thin = 1,
+ ...,
+ xlab = "", ylab = "",
+ cuts = 10,
+ at = do.breaks(c(-1.001, 1.001), cuts),
+ col.regions = topo.colors(100),
+ subset = thinned.indices(x, start = start, thin = thin))
+{
+ cormat <- cor(x[subset, ])
+ cormat <- cormat[, rev(seq(length = ncol(cormat)))]
+ levelplot(cormat,
+ main = main, ...,
+ cuts = cuts, at = at,
+ col.regions = col.regions,
+ xlab = xlab, ylab = ylab)
+}
+
+
+## the mcmc.list method wouldn't do any grouping (so should have
+## outer=TRUE by default). It hasn't been written yet because
+
+
+
+
+## The splom (FIXME: not yet written) method may be more useful
+
+## in progress, unexported. Planning to make it be like levelplot on
+## the lower diagonal (maybe ellipses instead of plain boxes) and
+## normal splom on the upper diagonal. Not much point in having tick
+## marks. Names in the middle save space, unlike in levelplot which
+## stupidly shows correlation=1 on the diagonal.
+
+
+splom.mcmc <-
+ function(x, main = attr(x, "title"),
+ start = 1, thin = 1,
+ as.matrix = TRUE,
+ xlab = "", ylab = "",
+ cuts = 10,
+ at = do.breaks(c(-1.001, 1.001), cuts),
+ col.regions = topo.colors(100),
+ ...,
+ pscales = 0,
+ subset = thinned.indices(x, start = start, thin = thin))
+{
+## cormat <- cor(x[subset, ])
+## cormat <- cormat[, rev(seq(length = ncol(cormat)))]
+ splom(as.data.frame(x[subset, ]),
+ as.matrix = as.matrix,
+ main = main, ...,
+ pscales = pscales,
+ cuts = cuts, at = at,
+ lower.panel = function(x, y, ...) {
+ corval <- cor(x, y)
+ grid::grid.text(lab = round(corval, 2))
+ },
+ col.regions = col.regions,
+ xlab = xlab, ylab = ylab)
+}
+
+
+
+
+
+
+
+
+### methods for densityplot (mcmc and mcmc.list)
+
+
+
+densityplot.mcmc <-
+ function(x, outer, aspect = "xy",
+ default.scales = list(relation = "free"),
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ xlab = "",
+ plot.points = "rug",
+ ...,
+ subset = thinned.indices(x, start = start, thin = thin))
+{
+ if (!missing(outer)) warning("specification of outer ignored")
+ data <- as.data.frame(x)
+ form <-
+ as.formula(paste("~",
+ paste(lapply(names(data), as.name),
+ collapse = "+")))
+
+
+### The following is one possible approach, but it does not generalize
+### to mcmclist objects.
+
+
+## densityplot(form, data = data,
+## outer = outer,
+## aspect = aspect,
+## default.scales = default.scales,
+## main = main,
+## xlab = xlab,
+## plot.points = plot.points,
+## subset = eval(subset),
+## ...)
+
+
+### This one does, with the only downside I can think of being that
+### subscripts, if used, will give indices in subsetted data, not
+### original. But that's true even if the original mcmc object was
+### itself already thinned.
+
+ densityplot(form, data = data[subset, ],
+ outer = TRUE,
+ aspect = aspect,
+ default.scales = default.scales,
+ main = main,
+ xlab = xlab,
+ plot.points = plot.points,
+ ...)
+}
+
+
+densityplot.mcmc.list <-
+ function(x, outer = FALSE, groups = !outer,
+ aspect = "xy",
+ default.scales = list(relation = "free"),
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ xlab = "",
+ plot.points = "rug",
+ ...,
+ subset = thinned.indices(x[[1]], start = start, thin = thin))
+{
+ if (groups && outer) warning("'groups=TRUE' ignored when 'outer=TRUE'")
+ datalist <- lapply(x, function(x) as.data.frame(x)[subset, ])
+ data <- do.call("rbind", datalist)
+ form <-
+ if (outer)
+ as.formula(paste("~",
+ paste(lapply(names(data), as.name),
+ collapse = "+"),
+ "| .run"))
+## as.formula(paste("~",
+## paste(names(data),
+## collapse = "+"),
+## "| .run"))
+ else
+ as.formula(paste("~",
+ paste(lapply(names(data), as.name),
+ collapse = "+")))
+## as.formula(paste("~",
+## paste(names(data),
+## collapse = "+")))
+ ##data[["index"]] <- seq(length = nrow(x[[1]]))[subset]
+ data[[".run"]] <- gl(length(datalist), nrow(datalist[[1]]))
+ if (groups && !outer)
+ densityplot(form, data = data,
+ outer = TRUE,
+ groups = .run,
+ aspect = aspect,
+ default.scales = default.scales,
+ main = main,
+ xlab = xlab,
+ plot.points = plot.points,
+ ...)
+ else
+ densityplot(form, data = data,
+ outer = TRUE,
+ aspect = aspect,
+ default.scales = default.scales,
+ main = main,
+ xlab = xlab,
+ plot.points = plot.points,
+ ...)
+}
+
+
+### methods for qqmath (mcmc and mcmc.list)
+
+
+
+
+qqmath.mcmc <-
+ function(x, outer, aspect = "xy",
+ default.scales = list(y = list(relation = "free")),
+ prepanel = prepanel.qqmathline,
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ ylab = "",
+ ...,
+ subset = thinned.indices(x, start = start, thin = thin))
+{
+ if (!missing(outer)) warning("specification of outer ignored")
+ data <- as.data.frame(x)
+ form <-
+ as.formula(paste("~",
+ paste(lapply(names(data), as.name),
+ collapse = "+")))
+ qqmath(form, data = data[subset, ],
+ outer = TRUE,
+ aspect = aspect,
+ prepanel = prepanel,
+ default.scales = default.scales,
+ main = main,
+ ylab = ylab,
+ ...)
+}
+
+
+qqmath.mcmc.list <-
+ function(x, outer = FALSE, groups = !outer,
+ aspect = "xy",
+ default.scales = list(y = list(relation = "free")),
+ prepanel = prepanel.qqmathline,
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ ylab = "",
+ ...,
+ subset = thinned.indices(x[[1]], start = start, thin = thin))
+{
+ if (groups && outer) warning("'groups=TRUE' ignored when 'outer=TRUE'")
+ datalist <- lapply(x, function(x) as.data.frame(x)[subset, ])
+ data <- do.call("rbind", datalist)
+ form <-
+ if (outer)
+ as.formula(paste("~",
+ paste(lapply(names(data), as.name),
+ collapse = "+"),
+ "| .run"))
+ else
+ as.formula(paste("~",
+ paste(lapply(names(data), as.name),
+ collapse = "+")))
+ ##data[["index"]] <- seq(length = nrow(x[[1]]))[subset]
+ data[[".run"]] <- gl(length(datalist), nrow(datalist[[1]]))
+ if (groups && !outer)
+ qqmath(form, data = data,
+ outer = TRUE,
+ groups = .run,
+ aspect = aspect,
+ prepanel = prepanel,
+ default.scales = default.scales,
+ main = main,
+ ylab = ylab,
+ ...)
+ else
+ qqmath(form, data = data,
+ outer = TRUE,
+ aspect = aspect,
+ default.scales = default.scales,
+ main = main,
+ ylab = ylab,
+ ...)
+}
+
+
+
+### methods for xyplot (mcmc and mcmc.list)
+
+
+
+xyplot.mcmc <-
+ function(x, outer, layout = c(1, ncol(x)),
+ default.scales = list(y = list(relation = "free")),
+ type = 'l',
+ start = 1, thin = 1,
+ ylab = "",
+ xlab = "Iteration number",
+ main = attr(x, "title"),
+ ...,
+ subset = thinned.indices(x, start = start, thin = thin))
+{
+ if (!missing(outer)) warning("specification of outer ignored")
+ data <- as.data.frame(x)
+ form <- eval(parse(text = paste(paste(lapply(names(data), as.name),
+ collapse = "+"), "~.index")))
+ data[[".index"]] <- seq(length = nrow(data))
+ xyplot(form, data = data[subset, ],
+ outer = TRUE,
+ layout = layout,
+ default.scales = default.scales,
+ type = type,
+ xlab = xlab,
+ ylab = ylab,
+ main = main,
+ ...)
+}
+
+
+
+xyplot.mcmc.list <-
+ function(x, outer = FALSE, groups = !outer,
+ aspect = "xy", layout = c(1, ncol(x[[1]])),
+ default.scales = list(y = list(relation = "free")),
+ type = 'l',
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ ylab = "",
+ ...,
+ subset = thinned.indices(x[[1]], start = start, thin = thin))
+{
+ if (groups && outer) warning("'groups=TRUE' ignored when 'outer=TRUE'")
+ datalist <- lapply(x, function(x) as.data.frame(x)[subset, ])
+ data <- do.call("rbind", datalist)
+ form <-
+ if (outer)
+ eval(parse(text = paste(paste(lapply(names(data), as.name),
+ collapse = "+"), "~.index | .run")))
+ else
+ eval(parse(text = paste(paste(lapply(names(data), as.name),
+ collapse = "+"), "~.index")))
+## form <-
+## if (outer)
+## as.formula(paste(paste(names(data),
+## collapse = "+"),
+## "~ index | .run"))
+## else
+## as.formula(paste(paste(names(data),
+## collapse = "+"),
+## "~ index"))
+ data[[".index"]] <- seq(length = nrow(datalist[[1]])) ## repeated
+ data[[".run"]] <- gl(length(datalist), nrow(datalist[[1]]))
+ if (groups && !outer)
+ xyplot(form, data = data,
+ outer = TRUE,
+ layout = layout,
+ groups = .run,
+ default.scales = default.scales,
+ type = type,
+ main = main,
+ ylab = ylab,
+ ...)
+ else
+ xyplot(form, data = data,
+ outer = TRUE,
+ layout = layout,
+ default.scales = default.scales,
+ type = type,
+ main = main,
+ ylab = ylab,
+ ...)
+}
+
+
+
+
+
+
+
+
+### methods for acfplot (mcmc and mcmc.list)
+
+
+
+panel.acfplot <-
+ function(..., groups = NULL)
+{
+ reference.line <- trellis.par.get("reference.line")
+ panel.abline(h = 0,
+ col = reference.line$col,
+ lty = reference.line$lty,
+ lwd = reference.line$lwd,
+ alpha = reference.line$alpha)
+ if (is.null(groups)) panel.xyplot(...)
+ else panel.superpose(..., groups = groups)
+}
+
+
+
+
+acfplot <- function(x, ...)
+ UseMethod("acfplot")
+
+
+
+acfplot.mcmc <-
+ function(x, outer,
+ prepanel = function(x, y, ...) list(ylim= c(-1, 1) * max(abs(y[-1]))),
+ panel = panel.acfplot,
+ type = "h",
+ aspect = "xy",
+ start = 1, thin = 1,
+ lag.max = NULL,
+ ylab = "Autocorrelation",
+ xlab = "Lag",
+ main = attr(x, "title"),
+ ...,
+ subset = thinned.indices(x, start = start, thin = thin))
+{
+ if (!missing(outer)) warning("specification of outer ignored")
+ getAcf <- function(x, lag.max)
+ {
+ as.vector(acf(x, lag.max = lag.max, plot = FALSE)$acf)
+ }
+ data <- as.data.frame(apply(x[subset, ], 2, getAcf, lag.max = lag.max))
+ form <-
+ eval(parse(text = paste(paste(lapply(names(data), as.name),
+ collapse = "+"), "~.lag")))
+ data[[".lag"]] <- seq(length = nrow(data))
+ xyplot(form, data = data,
+ outer = TRUE,
+ prepanel = prepanel,
+ panel = panel,
+ type = type,
+ aspect = aspect,
+ xlab = xlab,
+ ylab = ylab,
+ main = main,
+ ...)
+}
+
+
+
+
+acfplot.mcmc.list <-
+ function(x, outer = FALSE, groups = !outer,
+ prepanel = function(x, y, ..., groups = NULL, subscripts) {
+ if (is.null(groups)) list(ylim= c(-1, 1) * max(abs(y[-1])))
+ else list(ylim = c(-1, 1) * max(sapply(split(y, groups[subscripts]),
+ function(x)
+ max(abs(x[-1]), na.rm = TRUE ))))
+ },
+ panel = panel.acfplot,
+ type = if (groups) 'b' else 'h',
+ aspect = "xy",
+ start = 1, thin = 1,
+ lag.max = NULL,
+ ylab = "Autocorrelation",
+ xlab = "Lag",
+ main = attr(x, "title"),
+ ...,
+ subset = thinned.indices(x[[1]], start = start, thin = thin))
+{
+ if (groups && outer) warning("'groups=TRUE' ignored when 'outer=TRUE'")
+ getAcf <- function(x, lag.max)
+ {
+ as.vector(acf(x, lag.max = lag.max, plot = FALSE)$acf)
+ }
+ if (groups || outer)
+ {
+ datalist <-
+ lapply(x, function(x)
+ as.data.frame(apply(x[subset, ], 2,
+ getAcf, lag.max = lag.max)))
+ data <- do.call("rbind", datalist)
+ }
+ else
+ {
+ ## this is not quite valid, as we are combining multiple
+ ## series, but shouldn't be too bad (FIXME: should we warn?)
+
+ datalist <- lapply(x, function(x) x[subset, ])
+ data <-
+ as.data.frame(apply(do.call("rbind", datalist),
+ 2, getAcf, lag.max = lag.max))
+ }
+ form <-
+ if (outer)
+ as.formula(paste(paste(lapply(names(data), as.name),
+ collapse = "+"),
+ "~ .lag | .run"))
+ else
+ as.formula(paste(paste(lapply(names(data), as.name),
+ collapse = "+"),
+ "~ .lag"))
+ data[[".lag"]] <- seq(length = nrow(datalist[[1]])) ## repeated
+ data[[".run"]] <- gl(length(datalist), nrow(datalist[[1]]))
+ if (groups && !outer)
+ xyplot(form, data = data,
+ outer = TRUE,
+ groups = .run,
+ prepanel = prepanel,
+ panel = panel,
+ type = type,
+ aspect = aspect,
+ xlab = xlab,
+ ylab = ylab,
+ main = main,
+ ...)
+ else
+ xyplot(form, data = data,
+ outer = TRUE,
+ prepanel = prepanel,
+ panel = panel,
+ type = type,
+ aspect = aspect,
+ xlab = xlab,
+ ylab = ylab,
+ main = main,
+ ...)
+}
+
+
diff --git a/man/HPDinterval.Rd b/man/HPDinterval.Rd
new file mode 100644
index 0000000..34d3bcc
--- /dev/null
+++ b/man/HPDinterval.Rd
@@ -0,0 +1,44 @@
+\name{HPDinterval}
+\alias{HPDinterval}
+\alias{HPDinterval.mcmc}
+\alias{HPDinterval.mcmc.list}
+\title{Highest Posterior Density intervals}
+\description{
+ Create Highest Posterior Density (HPD) intervals for the parameters in
+ an MCMC sample.
+}
+\usage{
+HPDinterval(obj, prob = 0.95, ...)
+\method{HPDinterval}{mcmc}(obj, prob = 0.95, ...)
+\method{HPDinterval}{mcmc.list}(obj, prob = 0.95, ...)
+}
+\arguments{
+ \item{obj}{The object containing the MCMC sample - usually of class
+ \code{"mcmc"} or \code{"mcmc.list"}}.
+ \item{prob}{A numeric scalar in the interval (0,1) giving the target
+ probability content of the intervals. The nominal probability
+ content of the intervals is the multiple of \code{1/nrow(obj)}
+ nearest to \code{prob}.}
+ \item{\dots}{Optional additional arguments for methods. None are used
+ at present.}
+}
+\details{
+ For each parameter the interval is constructed from the empirical cdf
+ of the sample as the shortest interval for which the difference in
+ the ecdf values of the endpoints is the nominal probability. Assuming
+ that the distribution is not severely multimodal, this is the HPD interval.
+}
+\value{
+ For an \code{"mcmc"} object, a matrix with columns \code{"lower"} and
+ \code{"upper"} and rows corresponding to the parameters. The
+ attribute \code{"Probability"} is the nominal probability content of
+ the intervals. A list of such matrices is returned for an
+ \code{"mcmc.list"} object.
+}
+\author{Douglas Bates}
+\examples{
+data(line)
+HPDinterval(line)
+}
+\keyword{univar}
+\keyword{htest}
diff --git a/man/autocorr.diag.Rd b/man/autocorr.diag.Rd
index 3aec818..0441cff 100644
--- a/man/autocorr.diag.Rd
+++ b/man/autocorr.diag.Rd
@@ -31,6 +31,6 @@ A vector containing the autocorrelations.
\author{Russell Almond}
\seealso{
-\code{\link{autocorr}, \link{acf}}, \code{\link{autocorr.plot}}.}
+\code{\link{autocorr}}, \code{\link{acf}}, \code{\link{autocorr.plot}}.}
}
\keyword{ts}
diff --git a/man/autocorr.plot.Rd b/man/autocorr.plot.Rd
index 0912f79..61124a2 100644
--- a/man/autocorr.plot.Rd
+++ b/man/autocorr.plot.Rd
@@ -2,7 +2,8 @@
\alias{autocorr.plot}
\title{Plot autocorrelations for Markov Chains}
-\usage{autocorr.plot(x, lag.max, auto.layout = TRUE, ask=par("ask"), ...)}
+\usage{autocorr.plot(x, lag.max, auto.layout = TRUE,
+ask=dev.interactive(), ...)}
\arguments{
\item{x}{A Markov Chain}
diff --git a/man/cumuplot.Rd b/man/cumuplot.Rd
index 0e60a4a..248eed5 100644
--- a/man/cumuplot.Rd
+++ b/man/cumuplot.Rd
@@ -4,8 +4,8 @@
\usage{
cumuplot(x, probs=c(0.025,0.5,0.975), ylab="",
- lty=c(2,1), lwd=c(1,2), type="l", ask=TRUE, auto.layout=TRUE,
- col=1, ...)
+ lty=c(2,1), lwd=c(1,2), type="l", ask=dev.interactive(),
+ auto.layout=TRUE, col=1, ...)
}
\arguments{
diff --git a/man/gelman.plot.Rd b/man/gelman.plot.Rd
index 2c5d40d..aa1c257 100644
--- a/man/gelman.plot.Rd
+++ b/man/gelman.plot.Rd
@@ -5,8 +5,8 @@
\usage{
gelman.plot(x, bin.width = 10, max.bins = 50,
-confidence = 0.95, transform = FALSE, auto.layout = TRUE, ask = TRUE,
-col, lty, xlab, ylab, type, ...)
+confidence = 0.95, transform = FALSE, auto.layout = TRUE,
+ask = dev.interactive(), col, lty, xlab, ylab, type, ...)
}
\arguments{
diff --git a/man/geweke.plot.Rd b/man/geweke.plot.Rd
index bd54ee2..5c85ffd 100644
--- a/man/geweke.plot.Rd
+++ b/man/geweke.plot.Rd
@@ -4,7 +4,7 @@
\usage{
geweke.plot(x, frac1 = 0.1, frac2 = 0.5, nbins = 20,
- pvalue = 0.05, auto.layout = TRUE, ask = TRUE, ...)
+ pvalue = 0.05, auto.layout = TRUE, ask = dev.interactive(), ...)
}
\arguments{
diff --git a/man/plot.mcmc.Rd b/man/plot.mcmc.Rd
index a835b60..6bdd65b 100644
--- a/man/plot.mcmc.Rd
+++ b/man/plot.mcmc.Rd
@@ -4,7 +4,7 @@
\usage{
\method{plot}{mcmc}(x, trace = TRUE, density = TRUE, smooth = TRUE, bwf,
- auto.layout = TRUE, ask = par("ask"), ...)
+ auto.layout = TRUE, ask = dev.interactive(), ...)
}
\arguments{
diff --git a/man/trellisplots.Rd b/man/trellisplots.Rd
new file mode 100644
index 0000000..9b86395
--- /dev/null
+++ b/man/trellisplots.Rd
@@ -0,0 +1,195 @@
+\name{trellisplots}
+\title{Trellis plots for mcmc objects}
+\alias{densityplot.mcmc}
+\alias{levelplot.mcmc}
+\alias{qqmath.mcmc}
+\alias{xyplot.mcmc}
+\alias{densityplot.mcmc.list}
+\alias{levelplot.mcmc.list}
+\alias{qqmath.mcmc.list}
+\alias{xyplot.mcmc.list}
+\alias{acfplot}
+\alias{acfplot.mcmc}
+\alias{acfplot.mcmc.list}
+\usage{
+\method{densityplot}{mcmc}(x,
+ outer, aspect = "xy",
+ default.scales = list(relation = "free"),
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ xlab = "",
+ plot.points = "rug",
+ \dots,
+ subset)
+\method{densityplot}{mcmc.list}(x,
+ outer = FALSE, groups = !outer,
+ aspect = "xy",
+ default.scales = list(relation = "free"),
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ xlab = "",
+ plot.points = "rug",
+ \dots,
+ subset)
+\method{levelplot}{mcmc}(x, main = attr(x, "title"),
+ start = 1, thin = 1,
+ \dots,
+ xlab = "", ylab = "",
+ cuts = 10, at,
+ col.regions = topo.colors(100),
+ subset)
+\method{qqmath}{mcmc}(x,
+ outer, aspect = "xy",
+ default.scales = list(y = list(relation = "free")),
+ prepanel = prepanel.qqmathline,
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ ylab = "",
+ \dots,
+ subset)
+\method{qqmath}{mcmc.list}(x,
+ outer = FALSE, groups = !outer,
+ aspect = "xy",
+ default.scales = list(y = list(relation = "free")),
+ prepanel = prepanel.qqmathline,
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ ylab = "",
+ \dots,
+ subset)
+\method{xyplot}{mcmc}(x,
+ outer, layout = c(1, ncol(x)),
+ default.scales = list(y = list(relation = "free")),
+ type = 'l',
+ start = 1, thin = 1,
+ ylab = "",
+ xlab = "Iteration number",
+ main = attr(x, "title"),
+ \dots,
+ subset)
+\method{xyplot}{mcmc.list}(x, outer = FALSE, groups = !outer,
+ aspect = "xy", layout = c(1, ncol(x[[1]])),
+ default.scales = list(y = list(relation = "free")),
+ type = 'l',
+ start = 1, thin = 1,
+ main = attr(x, "title"),
+ ylab = "",
+ \dots,
+ subset)
+acfplot(x, \dots)
+\method{acfplot}{mcmc}(x, outer,
+ prepanel, panel,
+ type = 'h',
+ aspect = "xy",
+ start = 1, thin = 1,
+ lag.max = NULL,
+ ylab = "Autocorrelation",
+ xlab = "Lag",
+ main = attr(x, "title"),
+ ...,
+ subset)
+\method{acfplot}{mcmc.list}(x, outer = FALSE, groups = !outer,
+ prepanel, panel,
+ type = if (groups) 'b' else 'h',
+ aspect = "xy",
+ start = 1, thin = 1,
+ lag.max = NULL,
+ ylab = "Autocorrelation",
+ xlab = "Lag",
+ main = attr(x, "title"),
+ ...,
+ subset)
+}
+\description{
+ These methods use the Trellis framework as implemented in the
+ \code{lattice} package to produce space-conserving diagnostic plots
+ from \code{"mcmc"} and \code{"mcmc.list"} objects. The \code{xyplot}
+ methods produce trace plots. The \code{densityplot} methods and
+ \code{qqmath} methods produce empirical density and probability
+ plots. The \code{levelplot} method depicts the correlation of the
+ series. The \code{acfplot} methods plot the auto-correlation in the
+ series.
+}
+\arguments{
+ \item{x}{ an \code{"mcmc"} or \code{"mcmc.list"} object. }
+ \item{outer}{ for the \code{"mcmc.list"} methods, a logical flag to
+ control whether multiple runs of a series are displayed in the same
+ panel (they are if \code{FALSE}, not if \code{TRUE}). If specified
+ in the \code{"mcmc"} methods, this argument is ignored with a
+ warning. }
+ \item{groups}{ for the \code{"mcmc.list"} methods, a logical flag to
+ control whether the underlying \code{lattice} call will be supplied
+ a \code{groups} arguments indicating which run a data point
+ originated from. The panel function is responsible for handling
+ such an argument, and will usually differentiate runs within a panel
+ by using different graphical parameters. When \code{outer=FALSE},
+ the default of \code{groups} is \code{TRUE} if the corresponding
+ default panel function is able to make use of such information.
+ When \code{outer=FALSE}, \code{groups=TRUE} will be ignored with a
+ warning. }
+ \item{aspect}{ controls the physical aspect ratio of the panel. See
+ \code{\link[lattice:xyplot]{xyplot}} for details. The default for
+ these methods is chosen carefully - check what the default plot
+ looks like before changing this parameter.}
+ \item{default.scales}{ this parameter provides a reasonable default
+ value of the \code{scales} parameter for the method. It is unlikely
+ that a user will wish to change this parameter. Pass a value for
+ \code{scales} (see \code{\link[lattice:xyplot]{xyplot}}) instead,
+ which will override values specified here. }
+ \item{type}{ a character vector that determines if lines, points,
+ etc. are drawn on the panel. The default values for the methods are
+ carefully chosen. See
+ \code{\link[lattice:panel.xyplot]{panel.xyplot}} for possible
+ values. }
+ \item{thin}{ an optional thinning interval that is applied before the
+ plot is drawn.}
+ \item{start}{ an optional value for the starting point within the
+ series. Values before the starting point are considered part of the
+ "burn-in" of the series and dropped.}
+ \item{plot.points}{ character argument giving the style in which
+ points are added to the plot. See
+ \code{\link[lattice:panel.densityplot]{panel.densityplot}} for
+ details. }
+ \item{layout}{a method-specific default for the \code{layout} argument
+ to the lattice functions.}
+ \item{xlab,ylab,main}{Used to provide default axis annotations and
+ plot labels.}
+ \item{cuts, at}{ defines number and location of values where colors
+ change }
+ \item{col.regions}{ color palette used }
+ \item{lag.max}{ maximum lag for which autocorrelation is computed. By
+ default, the value chosen by \code{\link{acf}} is used }
+ \item{prepanel,panel}{ suitable prepanel and panel functions for
+ \code{acfplot}. The prepanel function omits the lag-0
+ auto-correlation (which is always 1) from the range calculations. }
+ \item{\dots}{other arguments, passed to the lattice function.
+ Documentation of the corresponding generics in the \code{lattice}
+ package should be consulted for possible arguments. }
+ \item{subset}{indices of the subset of the series to plot. The
+ default is constructed from the \code{start} and \code{thin}
+ arguments.}
+}
+\value{
+
+ An object of class \code{"trellis"}. The relevant
+ \code{\link[lattice:update.trellis]{update}} method can be used to
+ update components of the object and the
+ \code{\link[lattice:print.trellis]{print}} method (usually called by
+ default) will plot it on an appropriate plotting device.
+
+}
+\seealso{
+ \code{\link[lattice:Lattice]{Lattice}} for a brief introduction to
+ lattice displays and links to further documentation.
+}
+\author{ Deepayan Sarkar \email{Deepayan.Sarkar at R-project.org}}
+\examples{
+data(line)
+xyplot(line)
+xyplot(line[[1]], start = 10)
+densityplot(line, start = 10)
+qqmath(line, start = 10)
+levelplot(line[[2]])
+acfplot(line, outer = TRUE)
+}
+\keyword{hplot}
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-science/packages/r-cran-coda.git
More information about the debian-science-commits
mailing list