Simulation of multivariate normal rectangle probabilities and their derivatives - Theoretical and computational results

An extensive literature in econometrics and in numerical analysis has considered the problem of evaluating the multiple integral P(B; mu, Omega) = integral(a)(b)n(upsilon - mu, Omega)d upsilon = E(V)1(V is an element of B). where V is a m-dimensional normal vector with mean mu, covariance matrix Omega, and density n(upsilon - mu, Omega), and 1(V is an element of B) is an indicator for the event B = {V\a < V < b). A leading case of such an integral is the negative orthant probability, where B = {V\V < 0). The problem is computationally difficult except in very special cases. The multinomial probit (MNP) model used in econometrics and biometrics has cell probabilities that are negative orthant probabilities, with mu and Omega depending on unknown parameters (and, in general, on covariates). Estimation of this model requires. for each trial parameter vector and each observation in a sample, evaluation of P(B; mu, Omega) and of its derivatives with respect to mu and Omega. This paper surveys Monte Carlo techniques that have been developed for approximations of P(B; mu, Omega) and its linear and logarithmic derivatives, that limit computation while possessing properties that facilitate their use in iterative calculations for statistical inference: the Crude Frequency Simulator (CFS), Normal Importance Sampling (NIS), a Kernel-Smoothed Frequency Simulator (KFS), Stern's Decomposition Simulator (SDS), the Geweke-Hajivassiliou-Keane Simulator (GHK), a Parabolic Cylinder Function Simulator (PCF), Deak's Chi-squared Simulator (DCS), an Acceptance/Rejection Simulator JARS), the Gibbs Sampler Simulator (GSS), a Sequentially Unbiased Simulator (SUS), and an Approximately Unbiased Simulator (AUS). We also discuss Gauss and FORTRAN implementations of these algorithms and present our computational experience with them. We iind that GHK is overall the most reliable method.