4
$\begingroup$

I'm playing around with the standard deviation of the samples taken from the normal distribution. When calculating the integral

$$\int _{-\infty }^{\infty }\int _{-\infty }^{\infty }f(x) f(y) \sqrt{x^2-\frac{1}{2} (x+y)^2+y^2}dydx$$

Where $$f(x)=\text{PDF}[\text{NormalDistribution}[],x]=\frac{e^{-\frac{x^2}{2}}}{\sqrt{2 \pi }}$$

Mathematica returns 0:

f = PDF[NormalDistribution[]];
Integrate[Sqrt[(x^2 + y^2 - (x + y)^2/2)] f[x] f[y], {x, -Infinity, Infinity}, {y, -Infinity, Infinity}]
(*Outputs 0*)

I know that it can't be 0 because for all $x$, $f(x)>0$, and for all $(x,y)\in\mathbb{R}^2$, $\sqrt{x^2-\frac{1}{2} (x+y)^2+y^2} \geq 0$.

If we replace Integrate with NIntegrate, the output is 0.797884, which seems reasonable, but why does Integrate break?

Improve this question
$\endgroup$
1
3
$\begingroup$

Unfortunately, the problematic part of your integrand is hidden too well to make Mathematica see it. The sqrt-term that resembles a distance from the diagonal is the place where it goes wrong. If you look at it closely and assuming that x and y are real-valued, we get the following simplification

FullSimplify[Sqrt[x^2 + y^2 - 1/2 (x + y)^2], Element[{x, y}, Reals]]
(* Abs[x - y]/Sqrt[2] *)

This term needs to be treated with care. One way, as pointed out by Jim is to use the PrincipalValue option. However, just knowing the fact where the problem lies let's us solve this on our own.

Firstly, we can simplify your integrand and throw it into Mathematica's face that it needs to be careful (my version is 11.2 btw):

f = PDF[NormalDistribution[]];
g = Sqrt[(x^2 + y^2 - (x + y)^2/2)] f[x] f[y];

Integrate[FullSimplify[g, Element[{x, y}, Reals]],
  {x, -Infinity, Infinity}, {y, -Infinity, Infinity}]
(* Sqrt[2/π] *)

Secondly and highly underappreciated is to see the symmetry and integrate only in the valid half of the plane. This result needs than to be multiplied by two:

2*Integrate[g, {x, -Infinity, Infinity}, {y, x, Infinity}] // Simplify
(* Sqrt[2/π] *)

As Geoge pointed out, there is a way of doing the above but calculating both sides:

Integrate[g, {x, -Infinity, Infinity}, {y, -Infinity, x, Infinity}]
Improve this answer
$\endgroup$
4
  • $\begingroup$ MMA Version 8.0 gives Sqrt[2/π] to the original integral directly and without problem. So this seems to be a setback in higher versions. $\endgroup$ – Akku14 Feb 27 '18 at 7:40
  • $\begingroup$ this works also: Integrate[Sqrt[(x^2 + y^2 - (x + y)^2/2)] f[x] f[y], {x, -Infinity, Infinity}, {y, -Infinity, x, Infinity}] $\endgroup$ – george2079 Feb 27 '18 at 17:20
  • $\begingroup$ @george2079 That is nice. I did not know that you can give something like {y, -Infinity, x, Infinity} as integration bounds. Is this mentioned anywhere in the documentation? $\endgroup$ – halirutan♦ Feb 27 '18 at 17:32
  • $\begingroup$ that syntax is documented for NIntegrate ( "test for singularities at the intermediate points" ). I don't see it in the Integrate docs. $\endgroup$ – george2079 Feb 27 '18 at 18:11

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.