Difference of Gamma random variables Given two independent random variables $X\sim \mathrm{Gamma}(\alpha_X,\beta_X)$ and $Y\sim \mathrm{Gamma}(\alpha_Y,\beta_Y)$, what is the distribution of the difference, i.e. $D=X-Y$?
If the result is not well-known, how would I go about deriving the result?      
 A: I will outline how the problem can be approached and state
what I think the end result will be for the special case
when the shape parameters are integers, but not fill in the
details.


*

*First, note that $X-Y$ takes on values in $(-\infty,\infty)$
and so $f_{X-Y}(z)$ has support $(-\infty,\infty)$.

*Second, from the standard results that the 
density of the sum of two independent continuous random variables is the
convolution of their densities, that is,
$$f_{X+Y}(z) = \int_{-\infty}^\infty f_X(x)f_Y(z-x)\,\mathrm dx$$
and that the density of the random variable $-Y$ is
$f_{-Y}(\alpha) = f_Y(-\alpha)$, deduce that
$$f_{X-Y}(z) = f_{X+(-Y)}(z) = \int_{-\infty}^\infty f_X(x)f_{-Y}(z-x)\,\mathrm dx
= \int_{-\infty}^\infty f_X(x)f_Y(x-z)\,\mathrm dx.$$

*Third, for non-negative random variables $X$ and $Y$, note that the
above expression simplifies to
$$f_{X-Y}(z) = \begin{cases}
\int_0^\infty f_X(x)f_Y(x-z)\,\mathrm dx, & z < 0,\\
\int_{0}^\infty f_X(y+z)f_Y(y)\,\mathrm dy, & z > 0.
\end{cases}$$

*Finally, using parametrization $\Gamma(s,\lambda)$ to mean a
random variable with density 
$\lambda\frac{(\lambda x)^{s-1}}{\Gamma(s)}\exp(-\lambda x)\mathbf 1_{x>0}(x)$,
and with
$X \sim \Gamma(s,\lambda)$ and $Y \sim \Gamma(t,\mu)$  random variables, 
we have for $z > 0$ that
$$\begin{align*}f_{X-Y}(z) &=  \int_{0}^\infty 
\lambda\frac{(\lambda (y+z))^{s-1}}{\Gamma(s)}\exp(-\lambda (y+z))
\mu\frac{(\mu y)^{t-1}}{\Gamma(t)}\exp(-\mu y)\,\mathrm dy\\
&= \exp(-\lambda z) \int_0^\infty p(y,z)\exp(-(\lambda+\mu)y)\,\mathrm dy.\tag{1}
\end{align*}$$
Similarly, for $z < 0$,
$$\begin{align*}f_{X-Y}(z) &=  \int_{0}^\infty 
\lambda\frac{(\lambda x)^{s-1}}{\Gamma(s)}\exp(-\lambda x)
\mu\frac{(\mu (x-z))^{t-1}}{\Gamma(t)}\exp(-\mu (x-z))\,\mathrm dx\\
&= \exp(\mu z) \int_0^\infty q(x,z)\exp(-(\lambda+\mu)x)\,\mathrm dx.\tag{2}
\end{align*}$$

These integrals are not easy to evaluate but for the special case
$s = t$, Gradshteyn and Ryzhik, Tables of Integrals, Series, and Products,
Section 3.383, lists the value of
$$\int_0^\infty x^{s-1}(x+\beta)^{s-1}\exp(-\nu x)\,\mathrm dx$$
in terms of polynomial, exponential and Bessel functions of $\beta$
and this can be used to write down explicit expressions for $f_{X-Y}(z)$.

From here on, we assume that $s$ and $t$ are integers so
that $p(y,z)$ is a polynomial in $y$ and $z$ of degree $(s+t-2, s-1)$
and $q(x,z)$ is a polynomial in $x$ and $z$ of degree $(s+t-2,t-1)$.


*

*For $z > 0$, the integral $(1)$ 
is the sum of $s$ Gamma integrals with respect to $y$ with coefficients
$1, z, z^2, \ldots z^{s-1}$. It follows that the density of
$X-Y$ is proportional to a mixture density of 
$\Gamma(1,\lambda), \Gamma(2,\lambda), \cdots, \Gamma(s,\lambda)$
random variables for $z > 0$. Note that this result
will hold even if $t$ is not an integer.

*Similarly, for $z < 0$,
the density of
$X-Y$ is proportional to a mixture density of 
$\Gamma(1,\mu), \Gamma(2,\mu), \cdots, \Gamma(t,\mu)$
random variables flipped over, that is,
it will have terms such as $(\mu|z|)^{k-1}\exp(\mu z)$
instead of the usual $(\mu z)^{k-1}\exp(-\mu z)$.
Also, this result will hold even if $s$ is not an integer.
A: To my knowledge the distribution of the difference of two independent gamma r.v.’s was first studied by Mathai in 1993. He derived a closed form solution. I will not reproduce his work here. Instead I will point you to the original source. The closed form solution can be found on page 241 as theorem 2.1 in his paper On non-central generalized Laplacianness of quadratic forms in normal variables.
A: the difference of two independent or correlated Gamma random variables are special cases of McKay distribution.
The exact and complete answer can be find in:
Sum and difference of two squared correlated Nakagami variates in connection with the McKay distribution
Holm, H., Alouini, M.-S.
IEEE Transactions on Communications, 2004 Vol. 52; Iss. 8
