O(N) Hierarchical algorithm for computing the expectations of truncated multi-variate normal distributions in N dimensions

In this paper, we study the $N$-dimensional integral $\phi(a,b; A) = \int_{a}^{b} H(x) f(x | A) \text{d} x$ representing the expectation of a function $H(X)$ where $f(x | A)$ is the truncated multi-variate normal (TMVN) distribution with zero mean, $x$ is the vector of integration variables for the $N$-dimensional random vector $X$, $A$ is the inverse of the covariance matrix $\Sigma$, and $a$ and $b$ are constant vectors. We present a new hierarchical algorithm which can evaluate $\phi(a,b; A)$ using asymptotically optimal $O(N)$ operations when $A$ has "low-rank" blocks with "low-dimensional" features and $H(x)$ is "low-rank". We demonstrate the divide-and-conquer idea when $A$ is a symmetric positive definite tridiagonal matrix, and present the necessary building blocks and rigorous potential theory based algorithm analysis when $A$ is given by the exponential covariance model. Numerical results are presented to demonstrate the algorithm accuracy and efficiency for these two cases. We also briefly discuss how the algorithm can be generalized to a wider class of covariance models and its limitations.