<s>[INST] Density of sum of truncated normal and normal distribution Suppose that $\varepsilon\sim N(0, \sigma_\varepsilon)$ and $\delta\sim N^+(0, \sigma_\delta)$. What is the density function for $X = \varepsilon - \delta$?
This proof apparently appeared in a Query by M.A. Weinstein in Technometrics 6 in 1964, which stated that the density of $X$ is given by
$$f_X(x) = \frac{2}{\sigma} \phi\left(\frac{x}{\sigma}\right) \left(1 - \Phi\left(\frac{x\lambda}{\sigma}\right)\right),$$
where $\sigma^2 = \sigma_\varepsilon^2 + \sigma_\delta^2$ and $\lambda = \sigma_\delta / \sigma_\varepsilon$ and $\phi$ and $\Phi$ are the standard normal density and distribution functions, respectively. However, that paper is very difficult to find online. What is the proof that the density of $X$ takes the above form?
 [/INST] Ultimately I needed to work through the algebra a bit more to arrive at the specified form. For posterity, the full proof is given below.
Proof
First consider the distribution function of $X$, which is given by
    $$F(x) = \Pr(X \leq x) = \Pr(\varepsilon - \delta \leq x)$$
$$= \int_{\varepsilon - \delta \leq x} f_\varepsilon(\varepsilon) f_\delta(\delta) d\delta d\varepsilon$$
$$= \int_{\delta\in\mathbb{R}^+} f_\delta(\delta) \int_{\varepsilon\in (-\infty, x + \delta]} f_\varepsilon(\varepsilon) d\varepsilon d\delta.$$
    Substituting in known density functions yields
    $$\int_0^\infty 2\phi(\delta | 0, \sigma_\delta) \int_{-\infty}^{x + \delta} \phi(\varepsilon | 0, \sigma_\varepsilon) d\varepsilon d\delta$$
$$= 2\int_0^\infty \phi(\delta | 0, \sigma_\delta) \Phi(x + \delta | 0, \sigma_\varepsilon) d\delta.$$
    The density of $X$ is then given by
    $$f(x) = \frac{dF}{dx} = 2\int_0^\infty \phi(\delta | 0, \sigma_\delta) \phi(x + \delta | 0, \sigma_\varepsilon) d\delta.$$
    Using Sage to perform this integration, the result is given by
    $$f(x) = -\frac{{\left(\operatorname{erf}\left(\frac{\sigma_{\delta} x}{2 \, \sqrt{\frac{1}{2} \, \sigma_{\delta}^{2} + \frac{1}{2} \, \sigma_{\varepsilon}^{2}} \sigma_{\varepsilon}}\right) e^{\left(\frac{\sigma_{\delta}^{2} x^{2}}{2 \, {\left(\sigma_{\delta}^{2} \sigma_{\varepsilon}^{2} + \sigma_{\varepsilon}^{4}\right)}}\right)} - e^{\left(\frac{\sigma_{\delta}^{2} x^{2}}{2 \, {\left(\sigma_{\delta}^{2} \sigma_{\varepsilon}^{2} + \sigma_{\varepsilon}^{4}\right)}}\right)}\right)} e^{\left(-\frac{x^{2}}{2 \, \sigma_{\varepsilon}^{2}}\right)}}{2 \, \sqrt{\pi} \sqrt{\frac{1}{2} \, \sigma_{\delta}^{2} + \frac{1}{2} \, \sigma_{\varepsilon}^{2}}}.$$
    Defining $\lambda = \sigma_\delta / \sigma_\varepsilon$ and $\sigma^2 = \sigma_\varepsilon^2 + \sigma_\delta^2$, the following can be simplified:
    $$\frac{\sigma_{\delta} x}{2 \, \sqrt{\frac{1}{2} \, \sigma_{\delta}^{2} + \frac{1}{2} \, \sigma_{\varepsilon}^{2}} \sigma_{\varepsilon}} = \frac{\lambda x}{\sigma\sqrt{2}} = \frac{x}{(\sigma / \lambda) \sqrt{2}},$$
$$\frac{\sigma_{\delta}^{2} x^{2}}{2 \, {\left(\sigma_{\delta}^{2} \sigma_{\varepsilon}^{2} + \sigma_{\varepsilon}^{4}\right)}} = \frac{\lambda^2 x^2}{2\sigma^2} = \frac{x^2}{2(\sigma / \lambda)^2}.$$
    Thus,
    $$f(x) = -\frac{\exp\left(\frac{x^2}{2(\sigma / \lambda)^2}\right)\left(\operatorname{erf}\left(\frac{x}{(\sigma / \lambda) \sqrt{2}}\right) - 1\right)\exp\left(-\frac{x^2}{2\sigma_\varepsilon^2}\right)}{\sigma\sqrt{2\pi}}$$
$$= -\frac{\left(\operatorname{erf}\left(\frac{x}{(\sigma / \lambda) \sqrt{2}}\right) - 1\right) \exp\left(-x^2\left(\frac{1}{2\sigma_\varepsilon^2} - \frac{1}{2(\sigma / \lambda)^2}\right)\right)}{\sigma\sqrt{2\pi}}.$$
    Now,
    $$\operatorname{erf}\left(\frac{x}{(\sigma / \lambda) \sqrt{2}}\right) - 1 = \left(1 + \operatorname{erf}\left(\frac{x}{(\sigma / \lambda) \sqrt{2}}\right)\right) - 2 = 2\left(\frac{1}{2}\left(1 + \operatorname{erf}\left(\frac{x}{(\sigma / \lambda) \sqrt{2}}\right)\right) - 1\right)$$
$$= 2\left(\Phi\left(\frac{x\lambda}{\sigma}\right) - 1\right) = -2\left(1 - \Phi\left(\frac{x\lambda}{\sigma}\right)\right).$$
    Also,
    $$\frac{1}{2\sigma_\varepsilon^2} - \frac{1}{2(\sigma / \lambda)^2} = \frac{1}{2\sigma_\varepsilon^2} - \frac{\sigma_\delta^2}{2\sigma_\varepsilon^2(\sigma_\delta^2 + \sigma_\varepsilon^2)} = \frac{\sigma_\delta^2 + \sigma_\varepsilon^2 - \sigma_\delta^2}{2\sigma_\varepsilon^2(\sigma_\delta^2 + \sigma_\varepsilon^2)} = \frac{\sigma_\varepsilon^2}{2\sigma_\varepsilon^2(\sigma_\delta^2 + \sigma_\varepsilon^2)} = \frac{1}{2\sigma^2}.$$
    So,
    $$f(x) = 2\left(1 - \Phi\left(\frac{x\lambda}{\sigma}\right)\right)\frac{\exp\left(-\frac{x^2}{2\sigma^2}\right)}{\sigma\sqrt{2\pi}} = 2\left(1 - \Phi\left(\frac{x\lambda}{\sigma}\right)\right) \phi(x | 0, \sigma) = \frac{2}{\sigma}\phi\left(\frac{x}{\sigma}\right) \left(1 - \Phi\left(\frac{x\lambda}{\sigma}\right)\right).$$
</s>