Learn R Programming

trunmnt (version 1.0.0)

mtrunmnt: Creating an S3 object for computing the product moment

Description

mtrunmnt creates an S3 object designed to compute the product moment for a truncated multivariate normal distribution, utilizing the algorithm by Lee (2020). Key attributes of this object are the nodes and weights of the multivariate Gaussian quadrature and the probability of the truncation interval

Usage

mtrunmnt(
  mu,
  lower = -Inf,
  upper = Inf,
  Sigma = 1,
  Sigmae = 1,
  Z = matrix(rep(1, length(mu)), ncol = 1),
  D = matrix(1, ncol = 1, nrow = 1),
  nGH = 35
)

Value

A mtrunmnt object

Arguments

mu

**(Required)** Mean vector of the parent multivariate normal distribution.

lower

vector of lower limits. If the lower limits are the same, a scalar value can be given. Defaults to -Inf.

upper

Vector of upper limits. If the upper limits are the same, a scalar value can be given. Defaults to Inf.

Sigma

The variance-covariance matrix of the parent multivariate normal distribution. It must be given a symmetric positive definite matrix, if Sigmae, D and Z are not specified.

Sigmae

Vector of variances of error terms. If the variances are the same, a scalar value can be given. Defaults to 1.

Z

Design matrix for the random components. Defaults to \(n \times 1\) matrix of 1's where \(n\) is the dimension of mu. It must be specified carefully with the argument D. ncol(Z) determines the dimension of D.

D

Variance-covariance matrix of \(u\). See Details. If the random components are independent, you can specify either a vector of variances or a scalar value. A scalar value means that the random components have the same variance. Defaults to 1.

nGH

Number of quadrature points. Defaults to 35.

Details

Assume the parent multivariate normal distribution comes from a mixed-effects linear model: $$\mathbf{Y} = \mathbf{X}\boldsymbol{\beta} + \mathbf{Zu} + \boldsymbol{\epsilon}$$ where \(\mathbf{X}\) and \(\mathbf{Z}\) are design matrices corresponding to \(\boldsymbol{\beta}\) and \(\mathbf{u}\) representing fixed and random effects, respectively, and \(\boldsymbol{\epsilon}\) is the vector of errors. It is assumed that the random effects \(\mathbf{u}\) follows a multivariate normal distribution with mean \(\mathbf{0}\), and symmetric positive definite variance-covariance matrix \(\mathbf{D}\). As usual, the distribution of \(\boldsymbol{\epsilon}\) is assumed to be a multivariate normal with mean \(\mathbf{0}\) and variance-covariance matrix \(\sigma^2_{\epsilon} \mathbf{I}\), but for more flexibility, it can be assumed that the error terms are independent, but do not have equal variance. That is, **we assume** \(\boldsymbol{\epsilon} \sim N(\mathbf{0}, \mathbf{E})\) where \(\mathbf{E}\) is a diagonal matrix. Then, $$\mathbf{Y} \sim N(\boldsymbol{\mu}, \mathbf{\Sigma})$$ where \(\mathbf{\Sigma} = \mathbf{Z}\mathbf{D}\mathbf{Z}' + \mathbf{E}\). The variance-covariance structure in mtrunmnt can thus be specified either by providing the individual components \(\mathbf{D}, \mathbf{Z}\), and \(\mathbf{E}\), or by directly supplying the resulting overall variance-covariance matrix \(boldsymbol{\Sigma}\).

References

Lee, S.-C. (2020). Moments calculation for truncated multivariate normal in nonlinear generalized mixed models. Communications for Statistical Applications and Methods, Vol. 27, No. 3, 377–383.

Examples

Run this code
### Create a mtrunmnt objective ###

set.seed(123)
sigma2e <- 1
sigma2a <- 2
n <- 5
mu <- seq(-1,1, length.out = n)
y <- mu + rnorm(1, sd = sqrt(sigma2a)) + rnorm(n, sd = sqrt(sigma2e))
S <- matrix(sigma2a, ncol = n, nrow = n) + diag(sigma2e, n)
a  <- rep(-Inf, n)
b  <- rep(Inf, n)
a[y >= 0] <- 0
b[y <  0] <- 0
obj1 <- mtrunmnt(mu, lower = a, upper = b, Sigmae = sigma2e, D = sigma2a) 
obj2 <- mtrunmnt(mu, lower = a, upper = b, Sigma = S) 
probntrun(obj1)
probntrun(obj2)
prodmnt(obj1, c(2,2,0,0,0))
prodmnt(obj2, c(2,2,0,0,0))
meanvar(obj1)
meanvar(obj2)

Run the code above in your browser using DataLab