Joint distribution of the sum and product of two i.i.d. centered normal random variables

If $X$ and $Y$ are independent normally distributed random variables

$$X,Y\sim\mathcal{N}(0,\sigma^2)$$

How are the sum and product, $X+Y$ and $XY$, co-distributed? You can write the moment generating function

$$V=X+Y\qquad Z=XY$$ \begin{align} M_{V,Z}(s,t)&=\int _{-\infty }^{\infty }\int _{-\infty }^{\infty }\exp\left[s (v+z)+t (vz)\right]\psi (v,0,\sigma^2) \psi (z,0,\sigma^2)\mathrm{d}v\mathrm{d}z\\ &=\frac{\exp\left[\dfrac{s^2 \sigma ^2}{1-\sigma ^2 t}\right]}{\sqrt{1-\sigma ^4 t^2}} \end{align}

But I do not know what probability density function this belongs to, besides that it looks like the MGF for a modified Bessel function of the second kind and normal distribution. Attached is a histogram generated from $10^7$ such pairs of random numbers.

EDIT: Turns out Mathematica can handle this thing, and spits this monster out $$f(v,z;\sigma^2)=\frac{\sqrt[4]{-1} \left[\sqrt{-i \left(v^2-4 z\right)}-i \sqrt{i \left(v^2-4 z\right)}\right] }{\left| v^2-4 z\right| }\frac{\exp\left(-\dfrac{v^2-2 z}{2 \sigma ^2}\right)}{2 \pi \sigma ^2}$$ Does anyone know if this related to some well-known result I could read around?

• Ummm ... Your question seeks the joint pdf of $(V, Z)$ where $V=X+Y$ and $Z=X*Y$, where presumably you intend $X$ and $Y$ are assumed independent Normals. You then post some formulae that you claim Mathematica spat out as the solution to your problem, namely $f(x,y) = blah$, which (a) is self-evidently not posed in the form of your question i.e. is not the joint pdf of $(V,Z)$, and (b) in any event, we already know the joint pdf of $(X,Y)$ which is bivariate Normal with zero correlation, because that is your starting assumption ... which is not the equation you post. – wolfies Apr 8 '14 at 17:18
• Many apologies, the additional equation introduces confusion by using x and y as the independent variables. It gives the joint PDF of (V,Z). – user5751 Apr 8 '14 at 17:53
• Although it just a symbol but my professor always says to me to not ever use symbol $Z$ when dealing with normal distribution but not standard normal distribution. It just a random comment. – Tunk-Fey Apr 8 '14 at 18:39
• Thanks for your help Wolfies. I've included a plot of that function. It appears to match the numerical results, although the asymptote at $v^2-4z=0$ is very poorly resolved. The function Mathematica spat at me appears to be the solution, I've simply never encountered this flavour of distribution before. – user5751 Apr 8 '14 at 18:42
• @user5751 Perhaps you have entered the Mathematica output incorrectly here? Also, the Mathematica output you have posted appears unnecessarily messy - and methinks is incomplete without giving some thought to defining the domain of support. – wolfies Apr 9 '14 at 5:07

The usual method works. Assume without loss of generality that $\sigma^2=1$, then, for every measurable bounded function $u$, $$E(u(V,Z))=E(u(X+Y,XY)),$$ hence $$E(u(V,Z))=\iint u(x+y,xy)\,\mathrm e^{-(x^2+y^2)/2}\,\frac{\mathrm dx\mathrm dy}{2\pi},$$ and the task is to transform the RHS into $$\iint u(v,z)\,g(v,z)\,\mathrm dv\mathrm dz,$$ since then the function $g$ shall be the density of $(V,Z)$.
The change of variable $(x,y)\to(v,z)=(x+y,xy)$ is two-to-one, it has Jacobian $\mathrm dv\mathrm dz=|x-y|\,\mathrm dx\mathrm dy$, and it is inverted through the formulas $x,y=\frac12(v\pm\sqrt{v^2-4z})$ for every $(v,z)$ such that $4z\leqslant v^2$. Finally, $|x-y|=\sqrt{v^2-4z}$ and $x^2+y^2=v^2-2z$ hence $$E(u(V,Z))=\iint u(v,z)\,2\,\mathrm e^{-(v^2-2z)/2}\,\mathbf 1_{v^2\gt4z}\,\frac{\mathrm dv\mathrm dz}{2\pi\sqrt{v^2-4z}},$$ and the density of $(V,Z)$ is $$g(v,z)=\frac1\pi\,\mathrm e^{-v^2/2}\,\mathrm e^{z}\,\frac{\mathbf 1_{v^2\gt4z}}{\sqrt{v^2-4z}}.$$ If $\sigma^2\ne1$, $(V,Z)$ is transformed into $(\sigma V,\sigma^2Z)$ hence the density becomes $$g(v,z)=\frac1\pi\,\frac1{\sigma^2}\,\mathrm e^{-v^2/(2\sigma^2)}\,\mathrm e^{z/\sigma^2}\,\frac{\mathbf 1_{v^2\gt4z}}{\sqrt{v^2-4z}}.$$ In other words, in the formula you suggest at the end of your post, the funny looking factor $$\frac{\sqrt[4]{-1} \left[\sqrt{-i \left(v^2-4 z\right)}-i \sqrt{i \left(v^2-4 z\right)}\right] }{\left| v^2-4 z\right| }$$ should read $$\frac{2\,\mathbf 1_{v^2\gt4z}}{\sqrt{v^2-4z}}.$$
• Did wrote: $$g(v,z)=\frac1\pi\,\frac1{\sigma^2}\,\mathrm e^{-v^2/(2\sigma^2)}\,\mathrm e^{2z/\sigma^2}\,\frac{\mathbf 1_{v^2\gt4z}}{\sqrt{v^2-4z}}.$$ Do you perhaps mean: $$g(v,z)=\frac1\pi\,\frac1{\sigma^2}\,\mathrm e^{-(v^2 - 2 z)/(2\sigma^2)}\,\frac{\mathbf 1_{v^2\gt4z}}{\sqrt{v^2-4z}}.$$ – wolfies Apr 8 '14 at 20:05