Modeling Gaussian Distribution for 4 Variables

Question

I have 4 random variables A,B,C,D, that all have the same variance $\sigma^2$, but they are independent of each other. What I am trying to find, or to at least approximate, is the standard deviation of $\Delta$R, where $\Delta$R = $(A^2+B^2+C^2+D^2)^{1/2}$.

A,B,C,D have means $\mu_a$,$\mu_b$,$\mu_c$,$\mu_d$ respectively. The probability distribution function is:

E^(1/2 (-((x - Subscript[μ, a])^2/σ^2) - (y - 
 Subscript[μ, b])^2/σ^2 - (z - Subscript[μ, 
 c])^2/σ^2 - (u - Subscript[μ, 
 d])^2/σ^2))/(4 π^2 Sqrt[σ^8])

My first attempt to solve this was to convert to hyper-spherical coordinates, where

• a = $\rho$Cos($\psi$)
• b = $\rho$Sin($\psi$)Cos($\theta$)
• c = $\rho$Sin($\psi$)Sin($\theta$)Cos($\phi$)
• d = $\rho$Sin($\psi$)Sin($\theta$)Sin($\phi$)
• dV = $\rho^3$$Sin^2$($\psi$)Sin($\theta$)d$\psi$d$\theta$d$\phi$d$\rho$

and $\phi$ runs from 0 to $\pi$, $\theta$ runs from 0 to 2$\pi$, $\psi$ runs from 0 to 2$\pi$, and $\rho$ runs from 0 to r. Then $$p_r(r) = d/dr\int_0^r \int_0^{\pi} \int_0^{2\pi}\ \int_0^\pi\ PDF \rho^3Sin^2(\psi)Sin(\theta)\,d\psi\,d\theta\,d\phi\,d\rho$$

Mathematica can do the first integral and most of the second integral, but then it gets stuck.The result of the first integral is

1/3 Sin[θ] (2 B + 3 π BesselI[0, B] + 
(3 π BesselI[1, A])/A -(3 π BesselI[1, B])/B + 
(3 π StruveL[1, B])/B + 3 π StruveL[2, B])

Where A and B are constants, different from the ones earlier. (Sorry.) But the constants aren't super important.

A = 2$\mu_a$, and

B = 2$\mu_b$Cos($\theta$)+{2$\mu_c$$\rho$Cos($\phi$)+2$\mu_d$$\rho$Sin($\phi$)}Sin($\theta$)

The code for the first integral is below

temp1 = Sin[θ]*(Exp[A*Cos[ψ]] + 
Exp[B*Sin[ψ]])*(1/(2 I)*(Exp[I*ψ] - Exp[-I*ψ]))^2
Integrate[temp1, {ψ, 0, Pi}]

Now Mathematica can do the first 4 parts of the second integral, but it cannot find the integrals of the Struve function. When I simplify the fifth part of the second integral as much as i can to

\!\(
\*SubsuperscriptBox[\(∫\), \(0\), \(2\ π\)]\(
\*FractionBox[\(π\ StruveL[1, 
     m\ Cos[θ] + n\ Sin[θ]]\), \(n + 
    m\ Cot[θ]\)] ⅆθ\)\)

Mathematica cannot compute it. Here m = 2$\mu_b$$\rho$ and n = 2$\mu_c$$\rho$Cos($\phi$)+2$\mu_d$$\rho$Sin($\phi$) . If i was able to do the second and third integral, I could then find the first and second moment to calculate the standard deviation of $\Delta$R.

So my question is, is there a way to do this quadruple integral and find $p_r(r)$? If not, is there a way to approximate or numerically calculate the value I want? Or is there a simpler/better way to find what I want?

Thanks, Kyle.

Answer 1

Following @kglr the sum of squares $A^2+B^2+C^2+D^2$ divided by $\sigma^2$ has a noncentral chisquare distribution with 4 degrees of freedom and noncentrality parameter $\lambda=(\mu_A^2+\mu_B^2+\mu_C^2+\mu_D^2)/\sigma^2$. So the variance of the square root of the sum of squares can be obtained by finding the first and second raw moments:

m1 = Expectation[σ x^(1/2), x \[Distributed] NoncentralChiSquareDistribution[4, λ]];
m2 = Expectation[σ^2 x, x \[Distributed] NoncentralChiSquareDistribution[4, λ]];
variance = FullSimplify[m2 - m1^2]

(* (4 + λ) σ^2 - 1/8 E^(-λ/2) π ((3 + λ) σ BesselI[0, λ/4] + 
   (1 + λ) σ BesselI[1, λ/4])^2 *)

Answer 2

Update: Standard deviation of $(a^2 + b^2 + c^2 + d^2)^{1/2}$

With $ \lambda =( \mu _1^2+\mu _2^2+\mu _3^2+\mu _4^2 ) / \sigma^2$,

td = TransformedDistribution[Sqrt@t, Distributed[t, NoncentralChiSquareDistribution[4,λ]]];

FullSimplify[StandardDeviation[td]]

FullSimplify[Variance[td]]

Original post: Standard deviation of $a^2 + b^2 + c^2 + d^2$

StandardDeviation[NoncentralChiSquareDistribution[4, 
        Total[Table[Subscript[μ, i]^2, {i, {a, b, c, d}}]]/σ^2]]

$$\sqrt{\frac{4 \left(\mu _a^2+\mu _b^2+\mu _c^2+\mu _d^2\right)}{\sigma ^2}+8}$$

Why?

Let

aa = a / σ; bb = b / σ; cc = c  / σ; dd = d / σ;

Thus, aa is $N(\frac{\mu_a}{\sigma}, 1)$. Similarly, bb, cc and dd ar normal with variance 1.

Sum of squares of four Normal random variables each with variance 1 has a NoncentralChiSquareDistribution[ν, λ] with parameters $\nu = 4$ and $ \lambda =( \mu _1^2+\mu _2^2+\mu _3^2+\mu _4^2 ) / \sigma^2$:

tdist = TransformedDistribution[aa^2 + bb^2 + cc^2 + dd^2, 
          {Distributed[aa, NormalDistribution[Subscript[μ, a]/σ, 1], 
           Distributed[bb, NormalDistribution[Subscript[μ, b]/σ, 1],
           Distributed[cc, NormalDistribution[Subscript[μ, c]/σ, 1],
           Distributed[dd, NormalDistribution[Subscript[μ, d]/σ, 1]}]

and

StandardDeviation[tdist]