Dismiss Notice
Join Physics Forums Today!
The friendliest, high quality science and math community on the planet! Everyone who loves science is here!

Evaluating an integral in Mathematica, analytical form

  1. Sep 28, 2012 #1

    fluidistic

    User Avatar
    Gold Member

    So I've several huge integrals to solve by hands and I would like to check my result with Mathematica. However I don't really know how to do so.
    Here is one of the several integrals I must solve:
    [tex]\int \frac{1}{2\pi \sigma _x \sigma _y \sqrt{1-\rho ^2 }} \exp \{ \frac{-1}{2(1-\rho ^2)} \left ( \left ( \frac{x-\mu _x}{\sigma _x} \right ) ^2 +\left ( \frac{y-\mu _y}{\sigma _y} \right ) ^2 -\frac{2\rho (x-\mu _x)(y-\mu _y)}{\sigma _x \sigma _y} \right ] \} dx[/tex].
    Solving by hands I reached [tex]\frac{1}{\sqrt{2\pi}\sigma _y} \exp \{ \left [ -\frac{1}{2 (1-\rho ^2)} \right ] \left [ \left ( \frac{y-\mu _y }{\sigma _y} \right ) ^2 +\left ( \frac{\mu _x }{\sigma _x} \right ) ^2 + \frac{2 \rho \mu _x y}{\sigma _x \sigma _y} +\frac{\mu _x}{2} +\frac{\rho \sigma _x (y-\mu _y)}{\sigma _y} \right ] \}[/tex].
    I want to confirm this result via mathematica, I've tried the following command:
    Code (Text):
    f[x_,y_]:=1/(2*\[Pi]*sigma_X*sigma_Y*Sqrt[1-rho^2])*E^(-1/(2(1-rho^2))*(((x-mu_X)/sigma_X)^2+((y-mu_Y)/sigma_Y)^2-(2*rho*(x-mu_X)*(y-mu_Y))/(sigma_X*sigma_Y)))
    which produced no error.
    Then
    Code (Text):
    \[Integral]f[x,y]\[DifferentialD]x
    returned [tex]\int \frac{e^{-\frac{\frac{(x-\text{mu$\_$X}){}^2}{\text{sigma$\_$X}{}^2}+\frac{(y-\text{mu$\_$Y}){}^2}{\text{sigma$\_$Y}{}^2}-\frac{2 \text{rho} (x-\text{mu$\_$X}) (y-\text{mu$\_$Y})}{\text{sigma$\_$X} \text{sigma$\_$Y}}}{2 \left(1-\text{rho}^2\right)}}}{2 \pi \sqrt{1-\text{rho}^2} \text{sigma$\_$X} \text{sigma$\_$Y}} \, dx[/tex]
    I've tried "Simplify", "FullSimplify" but the integral is never calculated.
    Do you have any idea how to do so? Thanks!
     
  2. jcsd
  3. Sep 28, 2012 #2
    Does this verify your manual calculation?

    In[1]:= Integrate[1/(2*π*sigmaX*sigmaY*Sqrt[1-rho^2])*E^(-1/(2(1-rho^2))*(((x-muX)/sigmaX)^2+((y-muY)/sigmaY)^2-(2*rho*(x-muX)*(y-muY))/(sigmaX*sigmaY))),x]

    Out[1]= (Sqrt[-1 + rho^2]*Erfi[(muY*rho*sigmaX - muX*sigmaY + sigmaY*x - rho*sigmaX*y)/(Sqrt[2]*Sqrt[-1 + rho^2]*sigmaX*sigmaY)])/ (2*E^((muY - y)^2/(2*sigmaY^2))*Sqrt[2*Pi]*Sqrt[1 - rho^2]*sigmaY)
     
  4. Sep 29, 2012 #3

    fluidistic

    User Avatar
    Gold Member

    Hmm nice thank you, that's a way to integrate that works. I forgot to mention I did the integral from [itex]-\infty[/itex] to [itex]\infty[/itex].
    So I tried
    Code (Text):
    Integrate[
     1/(2*\[Pi]*sigmaX*sigmaY*Sqrt[1 - rho^2])*
      E^(-1/(2 (1 - rho^2))*(((x - muX)/sigmaX)^2 + ((y - muY)/
              sigmaY)^2 - (2*rho*(x - muX)*(y - muY))/(sigmaX*
              sigmaY))), {x, -\[Infinity], \[Infinity]}]
    which, after a long time, returned
    Code (Text):
    If[Re[1/((-1 + rho^2) sigmaX^2)] < 0 &&
       Re[(-1 + rho^2) sigmaX^2] <= 0,
      E^(-((muY - y)^2/(2 sigmaY^2))) Sqrt[2 \[Pi]]
        Sqrt[-(-1 + rho^2) sigmaX^2],
      Integrate[
       E^(-(((muX - x)^2/sigmaX^2 - (2 rho (muX - x) (muY - y))/(
         sigmaX sigmaY) + (muY - y)^2/sigmaY^2)/(
        2 (1 - rho^2)))), {x, -\[Infinity], \[Infinity]},
       Assumptions ->
        Re[1/((-1 + rho^2) sigmaX^2)] >= 0 ||
         Re[sigmaX^2 - rho^2 sigmaX^2] < 0]]/(2 \[Pi] Sqrt[1 - rho^2]
        sigmaX sigmaY)
    . In all cases I don't think it confirms my result by hands so I guess I have made some errors.
     
  5. Sep 29, 2012 #4
    If I make some hopefully reasonable assumptions, that you should verify, I get

    In[1]:= Assuming[Element[rho,Reals]&& Element[sigmaX,Reals]&& -1<rho<1&&sigmaX≠0,
    Integrate[1/(2*Pi*sigmaX*sigmaY*Sqrt[1-rho^2])*E^(-1/(2 (1-rho^2))*(((x-muX)/sigmaX)^2+((y-muY)/sigmaY)^2-(2*rho*(x-muX)*(y-muY))/(sigmaX*sigmaY))),{x,-Infinity,Infinity}]]

    Out[1]= Abs[sigmaX]/(E^((muY - y)^2/(2*sigmaY^2))*Sqrt[2*Pi]*sigmaX*sigmaY)

    Warning: If you look carefully at my previous post you will notice there is "Integrate[1/(2*π" and you might wonder where in the world I came up with "n"?!?! That is a scrape-n-paste error where the desktop published Greek Pi got mangled into an n in the process of going from the forum screen to the clipboard to Mathematica to the clipboard to the forum screen. This leads me to my next sentence.

    One bit of general advice for using Mathematica, the more of the "I absolutely must desktop publish my every keystroke" Koolaid(tm) you have drunk, just like befor your first post, the more difficulties you will tend to have with Mathematica. In the beginning Mathematica was a calculation engine. Only much later did they try to glue on a veneer of desktop publishing and that almost never works. Again and again here we see posts "Why doesn't this work?" and the answer is often that they are confusing making the expression look like a textbook with making the expression something they can get a correct answer from. There are long lists of things that look like a textbook but don't work to get answers.

    My apologies for not catching that. Hopefully there are no more errors in this. Damn. There was another of exactly the same error in this post which I have now manually fixed.
     
    Last edited: Sep 29, 2012
Know someone interested in this topic? Share this thread via Reddit, Google+, Twitter, or Facebook




Similar Discussions: Evaluating an integral in Mathematica, analytical form
Loading...