I have two random variables: x1 and x2. Each random variable has the density function: 2*Sqrt[1-x^2]/Pi. Both random variables are semi-circles on the -1 <= x <= 1 interval. The “area” under each random variable equals unity.
When I try to find the density function of Y = Abs[x1 - x2], I get stuck. Does anyone know of any tricks to find the density function for Y? My code is below.
I suppose this is the correct syntax:
dist = ProbabilityDistribution[2*Sqrt[1 - x^2]/Pi, {x, -1, 1}];
trasf = TransformedDistribution[
Abs[a - b], {a \[Distributed] dist, b \[Distributed] dist}]
NProbability[x > 1, x \[Distributed] trasf]
Plot[NProbability[x > y, x \[Distributed] trasf], {y, 0, 2}]
Unfortunately, it seems that the calculation of the PDF is very slow.
A numerical method to calculate the PDF:
pdf[z_] :=
If[0 < z < 2,
8/π^2 NIntegrate[
Sqrt[(1 - x^2) (1 - (x - z)^2)], {x, Max[z - 1, -1], Min[z + 1, 1]}],
0]
Maybe there’s some kind of trig substitution that would make
solvable?
The documentation page for PDF used to point to the probability distribution function. With 14.3 it points to Import and Export format.
Using David Trimas’ integral but with slightly simpler integration limits Mathematica 14.3 (Windows 11) results in the following:
pdf2 = Integrate[8 Sqrt[(1 - x^2) (1 - (d + x)^2)]/\[Pi]^2, {x, -1, 1 - d},
Assumptions -> 0 < d < 2]
(* (4 z ((4 + d^2) EllipticE[1 - 4/d^2] - 8 EllipticK[1 - 4/d^2]))/(3 \[Pi]^2) *)
Note that I used d rather than z as in David Trimas’ example. I only get the input repeated back using z and even using y. I have quit and restarted the kernel, closed and reopened Mathematica and still the same issue.
Here is a “plot” of the two implementations:
pdf[z_] := If[0 < z < 2, 8/\[Pi]^2 NIntegrate[Sqrt[(1 - x^2) (1 - (x - z)^2)], {x, Max[z - 1, -1], Min[z + 1, 1]}], 0]
Plot[{pdf2, pdf[d]}, {d, 0, 2}, PlotStyle -> {Thickness[0.02], {Thickness[0.005], Yellow}}]
