mcse.multi {mcmcse} | R Documentation |
Multivariate Monte Carlo standard errors for expectations.
Description
Function returns the estimate of the covariance matrix in the Markov Chain CLT using batch means or spectral variance methods (with different lag windows). The function also returns the Monte Carlo estimate.
Usage
mcse.multi(x, method = "bm", r = 3, size = NULL, g = NULL,
adjust = TRUE, blather = FALSE)
Arguments
x |
A matrix or data frame of Markov chain output. Number of rows is the Monte Carlo sample size. |
method |
Any of “ |
r |
The lugsail parameters ( |
size |
Represents the batch size in “ |
g |
A function that represents features of interest. |
adjust |
Defaults to |
blather |
If |
Value
A list is returned with the following components,
cov |
a covariance matrix estimate. |
est |
estimate of g(x). |
nsim |
number of rows of the input x. |
eigen_values |
eigen values of the estimate cov. |
method |
method used to calculate matrix cov. |
size |
value of size used to calculate cov. |
Adjustment_Used |
whether an adjustment was used to calculate cov. |
References
Vats, D., Flegal, J. M., and, Jones, G. L Multivariate output analysis for Markov chain Monte Carlo, Biometrika, 106, 321–-337.
Vats, D., Flegal, J. M., and, Jones, G. L. (2018) Strong Consistency of multivariate spectral variance estimators for Markov chain Monte Carlo, Bernoulli, 24, 1860–-1909.
See Also
batchSize
, which computes an optimal batch size.
mcse.initseq
, which computes an initial sequence estimator.
mcse
, which acts on a vector.
mcse.mat
, which applies mcse to each column of a matrix or data frame.
mcse.q
and mcse.q.mat
, which compute standard
errors for quantiles.
Examples
## Bivariate Normal with mean (mu1, mu2) and covariance sigma
n <- 1e3
mu <- c(2, 50)
sigma <- matrix(c(1, 0.5, 0.5, 1), nrow = 2)
out <- BVN_Gibbs(n, mu, sigma)
mcse.bm <- mcse.multi(x = out)
mcse.tuk <- mcse.multi(x = out, method = "tukey")
# If we are only estimating the mean of the first component,
# and the second moment of the second component
g <- function(x) return(c(x[1], x[2]^2))
mcse <- mcse.multi(x = out, g = g)