$\DeclareMathOperator{\erf}{erf}$
I want to find the cumulative distribution function (CDF) of the random variable $Z = |X| + |Y|$ where $X$ and $Y$ are two independent random variables which are normally distributed with mean 0 but different variance $\sigma_X^2$ and $\sigma^2_{Y}$.
$|X|$ and $|Y|$ follow the so-called Half-Normal distribution, which is a special case of the Folded-Normal distribution, when $\mu = 0$. The probability density functions (PDF) of $|X|$ and $|Y|$ for $x > 0$ are defined as follows (wikipedia):
\begin{align}
\label{pdf}
f_{|X|}(x, \sigma_{X}) = \frac{\sqrt{2}}{\sigma_{X} \sqrt{\pi}} \exp \left(- \frac{x^2}{2\sigma_{X}^2}\right) && f_{|Y|}(y, \sigma_{Y}) = \frac{\sqrt{2}}{\sigma_{Y} \sqrt{\pi}} \exp\left(- \frac{y^2}{2\sigma_{Y}^2}\right)
\end{align}
And their CDF:
\begin{align}
\label{cdf}
F_{|X|}(x, \sigma_{X}) = \erf\left(\frac{x}{\sigma_{X} \sqrt{2}}\right) && F_{|Y|}(y, \sigma_{Y}) = \erf\left(\frac{y}{\sigma_{Y} \sqrt{2}}\right)
\end{align}
@DilipSarwate gave a very elegant solution for the case $\sigma_{Y} = \sigma_{X}$ here, but I am interested in the general case when this is not true, i.e. one of the variable has a higher weight in the final value (bigger scale).
So far, I have determined the PDF of $Z$ by doing the convolution of $f_{|X|}$ and $f_{|Y|}$. Since $f_{|X|}$ and $f_{|Y|}$ are defined only for $x > 0$, we can truncate the convolution to the interval $[0,z]$:
\begin{align}
f_Z(z) &= \int_0^z f_{|Y|}(z-x) f_{|X|}(x) \, dx \\
f_Z(z) &= \frac{2}{\sigma_X \sigma_Y \pi} \int_0^z \exp\left(\frac{-(z-x)^2}{2\sigma_Y^2}\right) \exp\left(\frac{-x^2}{2\sigma_X^2}\right) \, dx \\
f_Z(z) &= \frac{2}{\sigma_X \sigma_Y \pi} \int_0^z
\exp\left(\frac{-(x^2 – 2xz + z^2)}{2\sigma_Y^2}\right) \exp\left(\frac{-x^2}{2\sigma_X^2}\right) \, dx \\
f_Z(z) &= \frac{2}{\sigma_X \sigma_Y \pi} \exp\left(\frac{-z^2}{2\sigma_Y^2}\right) \int_0^z \exp\left(\frac{-x^2+2xz}{2\sigma_Y^2}\right) \exp\left(\frac{-x^2}{2\sigma_X^2}\right) \, dx \\
f_Z(z) &= \frac{2}{\sigma_X \sigma_Y \pi} \exp\left(\frac{-z^2}{2\sigma_Y^2}\right) \int_0^z \exp\left(\frac{-(\sigma_X^2 + \sigma_Y^2)x^2 + 2z\sigma_X^2 x}{2 \sigma_X^2 \sigma_Y^2}\right) \, dx
\end{align}
We define $\sigma_Z = \sqrt{\left(\sigma_X^2 + \sigma_Y^2\right)}$:
\begin{align}
f_{Z}(z) &= \frac{2}{\sigma_X \sigma_Y \pi} \exp\left(\frac{-z^2}{2\sigma_Y^2}\right) \int_{0}^{z} \exp\left(\frac{-\sigma_Z^2 x^2 + 2z\sigma_X^2 x}{2 \sigma_X^2 \sigma_Y^2}\right) dx
\end{align}
We complete the square $ax^2 + bx$ to the form $a(x-h)^2 + k$ with $h = -\frac{b}{2a}$ and $k = – \frac{b^2}{4a}$, i.e., $h = \frac{z\sigma_X^2}{\sigma_Z^2}$ and $k = \frac{(z \sigma_X^2)^2}{\sigma_Z^2}$. This results in:
\begin{align}
f_{Z}(z) &= \frac{2}{\sigma_X \sigma_Y \pi}
\exp\left(\frac{-z^2}{2\sigma_Y^2}\right)
\int_{0}^{z}
\exp\left(\frac{-\sigma_Z^2\left(x – \frac{z \sigma_X^2}{\sigma_Z^2}\right)^2
+ \frac{(z \sigma_X^2)^2}{\sigma_Z^2}}{2\sigma_X^2 \sigma_Y^2}\right) dx \\
f_{Z}(z) &= \frac{2}{\sigma_X \sigma_Y \pi}
\exp\left(\frac{-z^2}{2\sigma_Y^2}\right)
\exp\left(\frac{(z\sigma_X^2)^2}{2\sigma_X^2 \sigma_Y^2 \sigma_Z^2}\right) \int_{0}^{z}
\exp\left(\frac{-\sigma_Z^2 \left(x – \frac{z \sigma_X^2}{\sigma_Z^2}\right)^2}{2\sigma_X^2 \sigma_Y^2}\right) dx \\
f_{Z}(z) &= \frac{2}{\sigma_X \sigma_Y \pi} \exp\left(\frac{-z^2(\sigma_Z^2 – \sigma_X^2)}{2\sigma_Y^2 \sigma_Z^2}\right) \int_{0}^{z} \exp\left(\frac{-\sigma_Z^2 \left(x – \frac{z \sigma_X^2}{\sigma_Z^2}\right)^2}{2\sigma_X^2 \sigma_Y^2}\right) dx \\
f_{Z}(z) &= \frac{2}{\sigma_X \sigma_Y \pi} \exp\left(\frac{-z^2}{2 \sigma_Z^2}\right) \int_{0}^{z} \exp\left(\frac{-\sigma_Z^2 \left(x – \frac{z \sigma_X^2}{\sigma_Z^2} \right)^2}{2\sigma_X^2 \sigma_Y^2}\right) dx
\end{align}
\begin{align}
f_{Z}(z) &= \frac{2}{\sigma_X \sigma_Y \pi} \exp\left(\frac{-z^2}{2 \sigma_Z^2}\right) \int_{0}^{z} \exp\left(\frac{-\left(x – \frac{z \sigma_X^2}{\sigma_Z^2}\right)^2}{2\frac{\sigma_X^2 \sigma_Y^2}{\sigma_Z^2}}\right) dx \\
\end{align}
We add the term $\frac{\sqrt{2 \pi (\frac{\sigma_X \sigma_Y}{\sigma_Z})^2}}{\sqrt{2 \pi (\frac{\sigma_X \sigma_Y}{\sigma_Z})^2}}$:
\begin{align}
f_{Z}(z) &= \frac{2}{\sigma_X \sigma_Y \pi}
\exp\left(\frac{-z^2}{2 \sigma_Z^2}\right)
\sqrt{2 \pi \left(\frac{\sigma_X \sigma_Y}{\sigma_Z}\right)^2}
\int_{0}^{z}
\frac{1}{\sqrt{2 \pi \left(\frac{\sigma_X \sigma_Y}{\sigma_Z}\right)^2}}
\exp\left(\frac{-\left(x – \frac{z \sigma_X^2}{\sigma_Z^2}\right)^2}{2\left(\frac{\sigma_X \sigma_Y}{\sigma_Z}\right)^2}\right) dx \\
f_{Z}(z) &= \frac{2 \sqrt{2}}{\sigma_Z \sqrt{\pi}}
\exp\left(\frac{-z^2}{2 \sigma_Z^2}\right)
\int_{0}^{z}
\frac{1}{\sqrt{2 \pi \left(\frac{\sigma_X \sigma_Y}{\sigma_Z}\right)^2}}
\exp\left(\frac{-\left(x – \frac{z \sigma_X^2}{\sigma_Z^2}\right)^2}{2\left(\frac{\sigma_X \sigma_Y}{\sigma_Z}\right)^2}\right) dx
\end{align}
As we can see, the expression in the integral represents the density of a Gaussian distribution with mean $\frac{z \sigma_X^2}{\sigma_Z^2}$ and standard deviation $\frac{\sigma_X \sigma_Y}{\sigma_Z}$. As a results we can evaluate it as follows:
\begin{align}
f_{Z}(z) &= \frac{2 \sqrt{2}}{\sigma_Z \sqrt{\pi}} \exp\left(\frac{-z^2}{2 \sigma_Z^2}\right) \left[\Phi\left(z, \frac{z \sigma_X^2}{\sigma_Z^2}, \frac{\sigma_X \sigma_Y}{\sigma_Z}\right) – \Phi\left(0, \frac{z \sigma_X^2}{\sigma_Z^2}, \frac{\sigma_X \sigma_Y}{\sigma_Z}\right)\right]
\end{align}
where $\Phi(x,m,s)$ is the value of the cumulative distribution function at $x$ of a Gaussian distribution with mean $m$ and standard deviation $s$.
Now we can compute the cumulative distribution function $F_Z$:
\begin{align}
F_{Z}(z) &= \int_{0}^{z} f_{Z}(x) dx \\
F_{Z}(z) &= \frac{2 \sqrt{2}}{\sigma_Z \sqrt{\pi}} \int_{0}^{z} \exp\left(\frac{-x^2}{2 \sigma_Z^2}\right) \left[\Phi\left(x, \frac{x \sigma_X^2}{\sigma_Z^2}, \frac{\sigma_X \sigma_Y}{\sigma_Z}\right) – \Phi\left(0, \frac{x \sigma_X^2}{\sigma_Z^2}, \frac{\sigma_X \sigma_Y}{\sigma_Z}\right)\right] dx \\
F_{Z}(z) &= \frac{2 \sqrt{2}}{\sigma_Z \sqrt{\pi}} \int_{0}^{z} \exp\left(\frac{-x^2}{2 \sigma_Z^2}\right) \left[
\frac{1}{2} \left( 1 + \erf \left(\frac{x-\frac{x \sigma_X^2}{\sigma_Z^2}}{\frac{\sigma_X \sigma_Y}{\sigma_Z} \sqrt{2}} \right) \right) –
\frac{1}{2} \left( 1 + \erf \left( \frac{-\frac{x \sigma_X^2}{\sigma_Z^2}}{\frac{\sigma_X \sigma_Y}{\sigma_Z} \sqrt{2}} \right) \right)
\right] dx
\\
F_{Z}(z) &= \frac{\sqrt{2}}{\sigma_Z \sqrt{\pi}} \int_{0}^{z} \exp\left(\frac{-x^2}{2 \sigma_Z^2}\right) \left[
\erf \left( \frac{x-\frac{x \sigma_X^2}{\sigma_Z^2}}{\frac{\sigma_X \sigma_Y}{\sigma_Z} \sqrt{2}} \right) –
\erf \left( \frac{-\frac{x \sigma_X^2}{\sigma_Z^2}}{\frac{\sigma_X \sigma_Y}{\sigma_Z} \sqrt{2}} \right)
\right] dx
\end{align}
Since the integral of the sum is equal to the sum of the integrals:
\begin{align}
F_{Z}(z) &= \frac{\sqrt{2}}{\sigma_Z \sqrt{\pi}} \left[
\int_{0}^{z}
\exp\left(\frac{-x^2}{2 \sigma_Z^2}\right)
\erf \left( \frac{x-\frac{x \sigma_X^2}{\sigma_Z^2}}{\frac{\sigma_X \sigma_Y}{\sigma_Z} \sqrt{2}} \right) dx –
\int_{0}^{z}
\exp\left(\frac{-x^2}{2 \sigma_Z^2}\right)
\erf \left( \frac{-\frac{x \sigma_X^2}{\sigma_Z^2}}{\frac{\sigma_X \sigma_Y}{\sigma_Z} \sqrt{2}} \right)
dx \right] \\
F_{Z}(z) &= \frac{\sqrt{2}}{\sigma_Z \sqrt{\pi}} \left[
\int_{0}^{z}
\exp\left(\frac{-x^2}{2 \sigma_Z^2}\right)
\erf \left( \frac{x \left(1-\frac{\sigma_X^2}{\sigma_Z^2} \right)}{\frac{\sigma_X \sigma_Y}{\sigma_Z} \sqrt{2}} \right) dx –
\int_{0}^{z}
\exp\left(\frac{-x^2}{2 \sigma_Z^2}\right)
\erf \left( \frac{-\frac{x \sigma_X^2}{\sigma_Z^2}}{\frac{\sigma_X \sigma_Y}{\sigma_Z} \sqrt{2}} \right)
dx \right] \\
F_{Z}(z) &= \frac{\sqrt{2}}{\sigma_Z \sqrt{\pi}} \left[
\int_{0}^{z}
\exp\left(\frac{-x^2}{2 \sigma_Z^2}\right)
\erf \left( x \frac{\sigma_Y}{\sqrt{2} \sigma_Z \sigma_X } \right) dx –
\int_{0}^{z}
\exp\left(\frac{-x^2}{2 \sigma_Z^2}\right)
\erf \left(- x \frac{\sigma_X}{\sqrt{2} \sigma_Z \sigma_Y } \right)
dx \right]
\end{align}
We know that $\erf(x)$ is an odd function and that $\exp(x)$ is an even function, as a result, their product is odd, and thus:
\begin{align}
F_{Z}(z) &= \frac{\sqrt{2}}{\sigma_Z \sqrt{\pi}} \left[
\int_{0}^{z}
\exp\left(\frac{-x^2}{2 \sigma_Z^2}\right)
\erf \left( x \frac{\sigma_Y}{\sqrt{2} \sigma_Z \sigma_X } \right) dx +
\int_{0}^{z}
\exp\left(\frac{-x^2}{2 \sigma_Z^2}\right)
\erf \left(x \frac{\sigma_X}{\sqrt{2} \sigma_Z \sigma_Y } \right)
dx \right]
\end{align}
This is where I am stuck. It seems that the integrals are not elementary. What can I do?
Just to verify my results in the case of equal variance, I am going to assume that $\sigma_X = \sigma_Y$. From 3.3 (41) in this paper, We know that:
$$
\int \exp \left( -a^2 x^2 \right) \left[ \erf(ax)\right]^n \, dx = \frac{\sqrt{\pi}}{2a(n+1)} \left[\erf(ax)\right]^{n+1}
$$
Thus:
\begin{align}
F^{\sigma_X = \sigma_Y}_{Z}(z) &= \frac{2 \sqrt{2}}{\sigma_Z \sqrt{\pi}}
\int_{0}^{z}
\exp\left(\frac{-x^2}{2 \sigma_Z^2}\right)
\erf \left( \frac{x}{\sqrt{2} \sigma_Z } \right) dx \\
F^{\sigma_X = \sigma_Y}_{Z}(z) &= \frac{2 \sqrt{2}}{\sigma_Z \sqrt{\pi}}
\frac{\sqrt{\pi}}{4 \frac{1}{\sqrt{2} \sigma_Z}} \erf(\frac{z}{\sqrt{2} \sigma_Z})^2 \\
F^{\sigma_X = \sigma_Y}_{Z}(z) &= \erf(\frac{z}{\sqrt{2} \sigma_Z})^2
\end{align}
Which is consistent with the solution from here.
However, as I said, I don't know how to solve the problem when $\sigma_X \neq \sigma_Y$. Any suggestions? I think I have no choice but to use a numerical approximation method, such as Simpson's rule.
Best Answer
$\DeclareMathOperator{\erf}{erf}$In the last step (so the one where you got stuck) you found two integrals of the form $$\int_{0}^{z} e^{-a^{2} x^{2}} {\rm erf}(bx)\, dx \quad.$$ I asked about this integral on Mathoverflow via this question. It turns out it can be written by means of a special function called Owen's T: $$\int_{0}^{z} e^{-a^{2} x^{2}} {\rm erf}(bx)\, dx = \frac{\tan^{-1}(b/a)}{a \sqrt{\pi}} - 2 \frac{\sqrt{\pi}}{a} T(\sqrt{2}az,b/a) \quad.$$ When you plug in the values $a = \sqrt{\frac{1}{2\sigma_{Z}^{2}}} = \frac{1}{\sqrt{2}\sigma_{Z}}, b_{1} = \frac{\sigma_{Y}}{\sqrt{2}\sigma_{Z}\sigma_{X}}$ for the first integral, and $b_{2} = \frac{\sigma_{X}}{\sqrt{2}\sigma_{Z}\sigma_{Y}} $ for the second integral (with the same value for $a$ as in the first one), you will obtain the cumulative distribution function of the convolution of two independent half-normal distributions.
When $\sigma_{X} = \sigma_{Y} = \sigma$, you can use the property $T(h,1) = \frac{1}{2} \Phi(h) (1- \Phi(h)) $ to find $F_{Z} (z) = {\rm erf} \big{(} \frac{z}{\sqrt{2}\sigma} \big{)}^{2}$, which I think is similar but not quite the same as the solution you found already, because you have $\sigma_{Z}$ in stead of $\sigma$, and $\sigma_{Z} = \sqrt{2} \sigma$ when $\sigma_{X} = \sigma_{Y} = \sigma$. So I'm not sure what's going on over here. Perhaps you know?
Edit: I made a small mistake. The cumulative distribution function and the error function are related as follows: $\Phi(x) = 1/2 + 1/2 \erf(x/\sqrt{2})$. I did not take into account the $\frac{1}{\sqrt{2}}$ factor. So the expression you found in case $\sigma_{1} = \sigma_{2} = \sigma$ is correct.