Normal variance mixtures: Distribution, density and parameter estimation

TL;DR

This paper introduces efficient algorithms for computing the distribution, density, and parameter estimation of normal variance mixtures, including high-dimensional cases, using advanced Monte Carlo methods, with implementations in R.

Contribution

It proposes novel algorithms for joint distribution, density evaluation, and parameter estimation of normal variance mixtures, extending existing methods to high dimensions and general cases.

Findings

Algorithms are fast and accurate even in high dimensions (~1000)

Monte Carlo methods improve computation over existing techniques

Implementation available in R package nvmix

Abstract

Normal variance mixtures are a class of multivariate distributions that generalize the multivariate normal by randomizing (or mixing) the covariance matrix via multiplication by a non-negative random variable W. The multivariate t distribution is an example of such mixture, where W has an inverse-gamma distribution. Algorithms to compute the joint distribution function and perform parameter estimation for the multivariate normal and t (with integer degrees of freedom) can be found in the literature and are implemented in, e.g., the R package mvtnorm. In this paper, efficient algorithms to perform these tasks in the general case of a normal variance mixture are proposed. In addition to the above two tasks, the evaluation of the joint (logarithmic) density function of a general normal variance mixture is tackled as well, as it is needed for parameter estimation and does not always exist in closed form in this more general setup. For the evaluation of the joint distribution function, the proposed algorithms apply randomized quasi-Monte Carlo (RQMC) methods in a way that improves upon existing methods proposed for the multivariate normal and t distributions. An adaptive RQMC algorithm that similarly exploits the superior convergence properties of RQMC methods is presented for the task of evaluating the joint log-density function. This allows the parameter estimation task to be accomplished via an EM-like algorithm where all weights and log-densities are numerically estimated. It is demonstrated through numerical examples that the suggested algorithms are quite fast; even for high dimensions around 1000 the distribution function can be estimated with moderate accuracy using only a few seconds of run time. Even log-densities around -100 can be estimated accurately and quickly. An implementation of all algorithms presented in this work is available in the R package nvmix (version >= 0.0.4).