lpmvnorm  R Documentation 
Computes the loglikelihood (contributions) of multiple exact or intervalcensored observations (or a mix thereof) from multivariate normal distributions and evaluates corresponding score functions.
lpmvnorm(lower, upper, mean = 0, center = NULL, chol, invchol, logLik = TRUE,
M = NULL, w = NULL, seed = NULL, tol = .Machine$double.eps, fast = FALSE)
slpmvnorm(lower, upper, mean = 0, center = NULL, chol, invchol, logLik = TRUE,
M = NULL, w = NULL, seed = NULL, tol = .Machine$double.eps, fast = FALSE)
ldmvnorm(obs, mean = 0, chol, invchol, logLik = TRUE)
sldmvnorm(obs, mean = 0, chol, invchol, logLik = TRUE)
ldpmvnorm(obs, lower, upper, mean = 0, chol, invchol, logLik = TRUE, ...)
sldpmvnorm(obs, lower, upper, mean = 0, chol, invchol, logLik = TRUE, ...)
lower 
matrix of lower limits (one column for each observation, 
upper 
matrix of upper limits (one column for each observation, 
obs 
matrix of exact observations (one column for each observation, 
mean 
matrix of means (one column for each observation, length is
recycled to length of 
center 
matrix of negative rescaled means (one column for each observation, length is
recycled to length of 
chol 
Cholesky factors of covariance matrices as

invchol 
Cholesky factors of precision matrices as

logLik 
logical, if 
M 
number of iterations, early stopping based on estimated errors is NOT implemented. 
w 
an optional matrix of weights with 
seed 
an object specifying if and how the random number generator
should be initialized, see 
tol 
tolerance limit, values smaller than 
fast 
logical, if 
... 
additional arguments to 
Evaluates the multivariate normal loglikelihood defined by means
and
chol
over boxes defined by lower
and upper
or for
exact observations obs
.
MonteCarlo (Genz, 1992, the default) and quasiMonteCarlo (Genz & Bretz, 2002)
integration is implemented, the latter with weights obtained, for example,
from packages qrng or randtoolbox. It is the responsibility of
the user to ensure a meaningful lattice is used. In case of doubt, use
plain MonteCarlo (w = NULL
) or pmvnorm
.
slpmvnorm
computes both the individual loglikelihood contributions
and the corresponding score matrix (of dimension J \times (J + 1) / 2 \times N
) if
chol
contains diagonal elements. Otherwise, the dimension is J
\times (J  1) / 2 \times N
. The scores for exact or mixed exactinterval
observations are computed by sldmvnorm
and sldpmvnorm
,
respectively.
More details can be found in the lmvnorm_src
package vignette.
The loglikelihood (logLik = TRUE
) or the individual contributions to the loglikelihood.
slpmvnorm
, sldmvnorm
, and sldpmvnorm
return the score
matrices and, optionally (logLik = TRUE
), the individual loglikelihood contributions
as well as scores for obs
, lower
, upper
, and
mean
.
Genz, A. (1992). Numerical computation of multivariate normal probabilities. Journal of Computational and Graphical Statistics, 1, 141–150.
Genz, A. and Bretz, F. (2002), Methods for the computation of multivariate tprobabilities. Journal of Computational and Graphical Statistics, 11, 950–971.
dmvnorm
, vignette("lmvnorm_src", package = "mvtnorm")
### five observations
N < 5L
### dimension
J < 4L
### lower and upper bounds, ie intervalcensoring
lwr < matrix(runif(N * J), nrow = J)
upr < matrix(runif(N * J), nrow = J)
### Cholesky factor
(C < ltMatrices(runif(J * (J + 1) / 2), diag = TRUE))
### corresponding covariance matrix
(S < as.array(Tcrossprod(C))[,,1])
### plain MonteCarlo (Genz, 1992)
w < NULL
M < 25000
### quasiMonteCarlo (Genz & Bretz, 2002, but with different weights)
if (require("qrng")) w < t(ghalton(M * N, J  1))
### loglikelihood
lpmvnorm(lower = lwr, upper = upr, chol = C, w = w, M = M)
### compare with pmvnorm
exp(lpmvnorm(lower = lwr, upper = upr, chol = C, logLik = FALSE, w = w, M = M))
sapply(1:N, function(i) pmvnorm(lower = lwr[,i], upper = upr[,i], sigma = S))
### loglik contributions and score matrix
slpmvnorm(lower = lwr, upper = upr, chol = C, w = w, M = M, logLik = TRUE)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.