[Math] Integrating a special skew normal — the CDF of a convolution of a normal with a truncated normal

definite integralsimproper-integralsintegrationprobabilityprobability distributions

I am having a little trouble trying to compute an integral. In short, I wish to solve the following:
$$F(x) = \int_{-\infty}^x \phi(au-b)\,\Phi(au+b)\,du $$

My intuition is that this might be impossible. So, failing that, I am trying to solve a special case:
$$ F(0) = \int_{-\infty}^0 \phi(au-b)\,\Phi(au+b)\,du $$

Edit: Where $$\phi(x) = \frac{1}{\sqrt{2\pi}}e^{-\frac{x^2}{2}} $$
And: $$ \Phi(x) = \int_{-\infty}^x \phi(t) \, dt $$

My solution method (so far) has been trying to use integration by parts (several times) but this is resulting in a rather large mess. I was curious if anyone had a (perhaps elegant) way of solving this problem?

I have spent considerable time using the following resources:

  • Wikipedia's list of integrals of Gaussian functions. These are very helpful.

  • This math overflow answer solves for the integral over the entire domain. I am looking for either the indefinite integral or the solution given the limits $(-\infty,0)$.

  • Here is a paper that includes many identities that I have found helpful.

  • I think this is a special case of (some sort of) Skew Normal distribution: another way of framing the question is to look for the CDF of the convolution of a normal pdf with a truncated normal pdf. This paper by Azzalini might be of passing interest.

Addition: I have noticed the following identity from Owen:
$$\int_{-\infty}^{0}\phi(ax)\Phi(bx)dx = \frac{1}{2\pi a} \arctan\left( \frac{a}{b} \right) $$

Best Answer

For every positive $a$, the change of variable $v=au$ yields $F(0)=G(b)/a$ with

$$ G(b)=\int_{-\infty}^0\varphi(v-b)\Phi(v+b)\mathrm dv. $$

First, some easy remarks: $G(0)=\frac12\left.\Phi(v)^2\right|_{-\infty}^0=\frac18$, and, since $\Phi\leqslant1$, $G(b)\leqslant\Phi(-b)$, in particular $G(+\infty)=0$.

Now, the computation of $G$. Note that $\Phi'=\varphi$ and $\varphi'(u)=-u\varphi(u)$, hence $G'(b)=H(b)+K(b)$ with $$ H(b)=\int_{-\infty}^0(v-b)\varphi(v-b)\Phi(v+b)\mathrm dv, $$ and $$ K(b)=\int_{-\infty}^0\varphi(v-b)\varphi(v+b)\mathrm dv. $$ Integrating by parts $H(b)$ yields $$ H(b)=\left.-\varphi(v-b)\Phi(v+b)\right|_{-\infty}^0+\int_{-\infty}^0\varphi(v-b)\varphi(v+b)\mathrm dv, $$ that is, $$ H(b)=-\varphi(b)\Phi(b)+K(b). $$ Now, $\varphi(v-b)\varphi(v+b)=\varphi(b\sqrt2)\varphi(v\sqrt2)$, hence $$ \sqrt2K(b)=\varphi(b\sqrt2)\left.\Phi(v\sqrt2)\right|_{-\infty}^0=\frac1{2}\varphi(b\sqrt2). $$ Thus, $$ G'(b)=-\varphi(b)\Phi(b)+2K(b)=-\varphi(b)\Phi(b)+\frac1{\sqrt2}\varphi(b\sqrt2). $$ Together with the limit values computed above, this shows that, for every positive $a$ and for every $b$,

$$ \int_{-\infty}^0\varphi(au-b)\Phi(au+b)\mathrm du=\frac{\Phi(b\sqrt2)-\Phi(b)^2}{2a}. $$

The same technique seems to yield, not explicit formulas for $G_x(b)$, but some funny-looking duality relations between $G_x(b)$ and $G_b(x)$, where $$ G_x(b)=\int_{-\infty}^x\varphi(v-b)\Phi(v+b)\mathrm dv, $$ namely, one might have (but this should be checked)

$$ G_x(b)+G_b(x)=\Phi(b\sqrt2)\Phi(x\sqrt2). $$

Related Question