# Estimating the cumulative probability of a bivariate normal distribution

I have a quick question regarding working out the probabilities of a bivariate normal distribution. To my knowledge, there is no nice closed-form for a cumulative distribution function for the bivariate normal distribution (Botev, 2016) so instead, we must numerically integrate through the bivariate normal distribution's probability density function (I think?).

I am referring to this thing in particular: https://mathworld.wolfram.com/BivariateNormalDistribution.html

So how is this done? By using the trapezoid rule or Simpson's rule? Or would that be an incorrect approach and a more accurate one is called for?

EDIT: a user pointed out a paper by Drzner & Welowosky (2010) which discusses numerical methods to integrate over a univariate distribution. This user also pointed out to me that the c.d.f. for this kind of distribution is

$$\mathbb P(X

• What do you mean by "not yet discovered"? The multivariate cdf is well-defined as$$\mathbb P(X<x,Y<y)=\int_{-\infty}^x\mathbb P(Y<y|X=x)\varphi(x;\mu_X,\sigma_X)\text dx$$and as a univariate integral can be accurately evaluated by numerical methods. See, eg, Drzner&Welowosky (2010) Jul 24, 2022 at 8:19
• Oh thank you @Xi'an --- I said "to my knowledge," but I'm not familiar with that paper :P Jul 25, 2022 at 3:18

The following is a screen-shot from Mark Schervish's (1984) paper on the approximation of the multivariate Normal probability function, based on (recursive) quadrature formulae, where $$G_L(x_N,\ldots,X_{N-L+1})=\int_{A(L)}^{B(L)}\cdots\int_{A(1)}^{B(1)} f(x_1,\ldots,x_N)\,\text dx_1\cdots\text dx_N$$ is the $$L$$-th inner integral. Schervish (1984) also provides a formula for the choice of the errors $$\delta_L$$ $$(L=1,\ldots, N-2)$$ towards the overall numerical error being at most $$\epsilon$$.