# Analytically integrating bivariate Gaussian PDF

I have this somewhat complicated expression in terms of the following 8 variables (note: actually, $$X$$ and $$Y$$ are variables, and $$a,b,c,d,e,f$$ are longer expressions as a funtion two other variables, say $$Z$$ and $$T$$, which I need to keep intact) of the form:

$$F(X,Y,a,b,c,d,e,f)=aXYe^{bX+cY+dXY+eX^2+fY^2}$$

And I need to perform the double integral over the two dimensions X and Y:

$$G(a,b,c,d,e,f)=\int\limits^1_0\int\limits^\infty_0 F(X,Y,a,b,c,d,e,f) dYdX$$

where the parameters a, b, c, d, e, and f all contain other variables which I need to keep intact. Because of this, Numerical Integration won't quite work. Ultimately I need to find the analytical form of $$G(a,b,c,d,e,f)$$ in Mathematica, even though it'll be quite complicated.

I've tried simply using the Integrate function:

$Assumptions = {{a, b, c, d, e, f, X, Y} \[Element] Reals} F = a X Y Exp[b X + c Y + d X Y + e X^2 + f Y^2] G=Integrate[F,{Y,0,1},{X,0,\[Infinity]}]  but no luck unfortunately, Mathematica works for quite awhile and then finally gives up, returning the input. I can indeed perform one of the integrals like so: $Assumptions = {{a, b, c, d, e, f, X, Y} \[Element] Reals,f<0}
F = a X Y Exp[b X + c Y + d X Y + e X^2 + f Y^2]
G1=Integrate[F,{Y,0,\[Infinity]}]//Simplify


but the result is complex enough that Mathematica still can't solve the remaining integral separately. Simplifying the expression gives me three terms that could be integrated over X separately:

int1 = (I a)/(4 f^(3/2)) X (Exp[G1[[3, 2]] + G1[[6, 1, 2, 2]]] // Simplify) (2 Sqrt[-f]) // Simplify
int2 = (I a)/(4 f^(3/2)) X G1[[3]] (Sqrt[\[Pi]] (c + d X)) //Simplify
int3 = (I a)/(4 f^(3/2))X G1[[3]] (Sqrt[\[Pi]] (c + d X) G1[[6, 2, 3, 2]]) // Simplify


Luckily, the first two integrals can be done fairly easily:

result1=Integrate[int1, {X, 0, 1}]
result2=Integrate[int2, {X, 0, 1}]


Which Mathematica can indeed perform. However, the third portion, int3, seemingly can't be solved as it has the following form:

$$\int\limits^1_0 j X(1+kX) \text{Erf} \lbrack m+nX \rbrack e^{oX+pX^2} dX$$

Where Erf[] is the error function, and $$j,k,m,n,o,p$$ are new parameters in terms of $$a,b,c,d,e,f$$. So ultimately, it boils down to that integral, which I have no idea how to solve! Does anyone have any ides which could help?

• I think that that does not have finite closed-form expression in terms of very large class of special functions.Try NIntegrate. – Mariusz Iwaniuk Jul 28 at 15:21