mnorm: An R Package for Calculation and Differentiation of Conditional Multivariate Normal Densities and Probabilities

Technical Summary of "mnorm: An R Package for Calculation and Differentiation of Conditional Multivariate Normal Densities and Probabilities"

Problem Statement
The multivariate normal distribution (MVN) is fundamental to econometric modeling, particularly for maximum-likelihood estimation (MLE) and Bayesian inference. While existing R packages such as mvtnorm, pbivnorm, and mvnfast provide robust tools for calculating MVN probabilities and densities, they often lack specific functionalities required for advanced optimization tasks. Specifically, many practical applications require:

1. Analytical Gradients: The ability to differentiate MVN densities and probabilities with respect to various parameters (covariances, integration limits, and conditioned components) to accelerate numerical optimization routines.
2. Conditional Calculations: Efficient computation of conditional and truncated MVN probabilities.
3. Non-Normal Marginals: Handling cases where marginal distributions are non-normal via the Gaussian copula.
4. Performance: High-speed computation, particularly for bivariate and trivariate cases, and support for parallel processing.

Existing solutions often rely on numerical differentiation, which is computationally expensive, or lack the specific gradient capabilities needed for complex econometric models like sample selection and endogenous switching models.

Methodology
The mnorm package addresses these limitations through a combination of theoretical derivations and high-performance computing implementations:

• Analytical Derivations: The core contribution is the derivation of closed-form analytical formulas for the gradients of MVN densities and probabilities. These formulas cover:
  • Derivatives with respect to the argument vector ($x$).
  • Derivatives with respect to the covariance matrix ($\Sigma$), including handling symmetry constraints.
  • Derivatives with respect to integration limits (lower and upper bounds).
  • Derivatives with respect to conditioned (given) components in conditional distributions.
  • Derivatives with respect to parameters of non-normal marginal distributions (e.g., Student's $t$, logistic, and PGN distributions) when used within a Gaussian copula framework.
• Algorithmic Implementation:
  • Low Dimensions ($m \leq 3$): The package utilizes Gauss-Legendre quadrature with 30 sample points for bivariate and trivariate probabilities, implementing methods by Drezner and Wesolowsky (bivariate) and Drezner (trivariate).
  • High Dimensions ($m > 3$): The package employs the Geweke-Hajivassiliou-Keane (GHK) algorithm with a deterministic Richtmyer sequence and variable reordering.
  • Conditional Probabilities: The package calculates conditional densities and probabilities by deriving the conditional mean and covariance matrices and applying the standard density/probability functions to the reduced dimensionality.
• Computational Architecture:
  • The package is written in Rcpp and utilizes RcppArmadillo for linear algebra operations to ensure speed.
  • It supports Open Multi-Processing (OMP) for parallel (multi-core) computing.
  • It includes specific bounds for log-probabilities and their derivatives to assist in proving consistency and asymptotic normality for estimators.

Key Contributions

1. Differentiation Capabilities: The primary contribution is the provision of functions (dmnorm and pmnorm) that return not only densities and probabilities but also their analytical gradients with respect to $x$, $\Sigma$, integration limits, and conditioned values.
2. Conditional and Truncated Support: The package offers efficient routines for conditional MVN densities and probabilities, allowing users to specify indices of conditioned components.
3. Gaussian Copula Integration: It supports non-normal marginals by allowing users to specify distributions (PGN, Student's $t$, logistic, normal) and differentiating probabilities with respect to their parameters.
4. Performance Optimization: By implementing analytical gradients, the package avoids the computational overhead of numerical differentiation. It also optimizes speed for low-dimensional cases using specialized quadrature methods and for high-dimensional cases via the GHK algorithm.
5. Theoretical Bounds: The authors derive useful upper bounds on the derivatives of log-probabilities and lower bounds on log-probabilities, which are relevant for theoretical econometric analysis.

Results
The paper presents extensive efficiency analyses comparing mnorm against popular alternatives (mvtnorm, mvnfast, pbivnorm, TruncatedNormal, tlrmvnmvt):

• Gradient Speed: Analytical gradients provided by mnorm significantly outperform numerical gradients. In 5-dimensional cases, calculating analytical gradients for densities is comparable in time to calculating the density itself, whereas numerical gradients require multiple function calls ($m+1$ or more), making them 5–10 times slower. For probabilities, the speed advantage of analytical gradients over numerical ones is even more pronounced, especially when the number of unconditioned dimensions is low.
• Density Calculation: The dmnorm function is faster than dmvn (from mvnfast) and dmvnorm (from mvtnorm) in both single-core and 8-core parallel configurations for 5-dimensional and 100-dimensional cases.
• Probability Calculation:
  • Bivariate/Trivariate: pmnorm is notably faster than pbivnorm and the TVPACK algorithm from mvtnorm (which requires loops for vectorization), while maintaining comparable accuracy (maximum absolute differences $\approx 10^{-16}$).
  • High-Dimensional: In dimensions 16, 64, 128, and 512, pmnorm demonstrates speed and accuracy comparable to the tlrmvnmvt package and TruncatedNormal, often outperforming mvtnorm in speed.
• Conditional Efficiency: The calculation of conditional probabilities is highly efficient because the derivatives with respect to conditioned values are functions of derivatives with respect to other parameters, requiring minimal additional computation time.

Significance and Claims
The authors claim that the mnorm package is of "special interest for the implementation of the maximum-likelihood estimators in econometric models." Specifically, the ability to calculate analytical gradients substantially speeds up numeric optimization routines, which is critical for complex models such as:

• Multivariate and multinomial sample selection models.
• Endogenous switching econometric models.

The paper asserts that by providing fast, accurate, and differentiable routines for conditional MVN densities and probabilities, mnorm eases the implementation of these estimators. The authors remain modest regarding the package's novelty, positioning it as an extension of existing tools that fills specific gaps in differentiation and conditional handling, rather than replacing the fundamental algorithms of MVN calculation found in other packages. The provided bounds on log-probability derivatives are presented as a supplementary contribution to the literature on inverse Mills ratios and estimator consistency.