Test Data for the Cumulative Bivariate Normal Distribution

Using Java I implemented the double precision algorithm to compute the cumulative bivariate normal distribution found in A.Genz, “Numerical computation of rectangular bivariate and trivariate normal and t probabilities”, Statistics and Computing, 14, (3), 2004.

$$M(a,b, \rho)=\frac{1}{2\pi\sqrt{1-\rho^2}}\int_{-\infty}^a\int_{-\infty}^b \exp\left(-\frac{x^2-2\rho xy +y^2}{2(1-\rho^2)}\right)dxdy$$

To validate the implementation I had hoped to find an easily accessible table of values to check against, alas even my trusted copy of Abramowitz and Stegun: Handbook of Mathematical Functions did not have tablulated values, so I considered alternatives. Firstly, I considered the analytic result
$$M(0,0,\rho)=\frac{1}{4}+\frac{\arcsin(\rho)}{2\pi}$$

which leads quite naturally to the simple test values

Elegant as they are, these special cases only cover only a small part of the algorithm, however randomly picking \(\rho\in(-1,1)\) and comparing to a computed value of \(1/4+\arcsin(\rho)/2\pi\), covers much of the code and is quite accessible since a double precision implementation of \(\arcsin\) is available in most languages.

Finally, I generated values from the statistical package R (library mvtnorm) at specific points which ensure 100% coverage of the algorithm.

There is a table of computed values in Haug’s book, “The complete guide to option pricing formulas“, Second Edition, however it has only one value in the interesting area \(|\rho| \in(0.925,1)\), and that test value is an extreme.