Multiplying two Gamma distributions over the same variable I am looking at a software library where there is a function that multiplies two Gamma distributions defined over the same random variable. So, it is basically multiplying two Gamma densities with shape and rate parameters.
So, the Gamma distribution parameterised by the shape and rate parameters is given as:
$$
D(x; \alpha, \beta) = \frac{\beta^{\alpha} x^{\alpha-1} e^{-\beta x}}{\Gamma(\alpha)}
$$ 
where $\alpha$ is the shape parameter and $\beta$ is the rate parameter and $\Gamma$ is the Gamma function defined as $\Gamma(\alpha) = \int_{0}^{\infty} x^{\alpha-1} e^{-x} dx$. 
Do, I want to get an expression for the product of $D(x; \alpha_0, \beta_0)$ and $D(x; \alpha_1, \beta_1)$. Apparently, this should be another Gamma distribution and according to this software implementation, the shape and rate parameters for this new Gamma distribution should be: $\alpha_{new} = \alpha_0 + \alpha_1 - 1$ and $\beta_{new} = \beta_0 + \beta_1$. 
Now, apparently this should be simple algebra but I have been unsuccessfully been trying to formulate this product in a way that I can equate the terms with the expression for the Gamma distribution that would arrive at this result.
 A: The product of two gamma densities on the same variable is proportional to another gamma density ... as long as the sum of the two shape parameters in the product is > 1. When that's not true, it's not the case that the product is proportional to a density.
consider: 
$$g(x) = f(x;\alpha_1,\beta_1)\cdot f(x;\alpha_2,\beta_2)$$
where $f(x;\alpha,\beta) = \frac{\beta^\alpha}{\Gamma(\alpha)} x^{\alpha-1} e^{- \beta x }\,,\quad x>0,\,\alpha,\beta>0$
$$g(x) = \frac{\beta_1^{\alpha_1}}{\Gamma(\alpha_1)} x^{\alpha_1-1} e^{- \beta_1 x }\frac{\beta_2^{\alpha_2}}{\Gamma(\alpha_2)} x^{\alpha_2-1} e^{- \beta_2 x }$$
$$\quad\quad = \frac{\beta_1^{\alpha_1}}{\Gamma(\alpha_1)}\frac{\beta_2^{\alpha_2}}{\Gamma(\alpha_2)}  x^{\alpha_1-1}x^{\alpha_2-1}  e^{- \beta_1 x }e^{- \beta_2 x }$$
$$\quad\quad = \frac{\beta_1^{\alpha_1}}{\Gamma(\alpha_1)}\frac{\beta_2^{\alpha_2}}{\Gamma(\alpha_2)}  x^{(\alpha_1+\alpha_2-1)-1}  e^{- (\beta_1+ \beta_2) x }$$
Clearly, for that to be proportional to a $\text{gamma}(\alpha_1+\alpha_2-1,\beta_1+ \beta_2)$ density, its parameters must obey the constraints on the parameters (i.e. $\alpha_1+\alpha_2-1>0$, or $\alpha_1+\alpha_2>1$).
What happens when it's not the case that $\alpha_1+\alpha_2-1>0$? 
To simplify things consider $\beta_1+ \beta_2=1$, and let $\delta=\alpha_1+\alpha_2-1$. The gamma pdf is:
$$\frac{1}{\Gamma(\delta)} x^{\delta-1} e^{- x }\,,\quad x>0\,.$$
Now in the case where $\delta\leq 0$, the integral that defines $\Gamma(\delta)$ doesn't converge (though the function is generally extended to negative values via analytic continuation, that's no longer the value of the integral).
Which is to say, we no longer have something proportional to a valid pdf.
When $\alpha_1+\alpha_2>1$, the multiplicative constant out the front of the gamma density is
$$\frac{\beta_1^{\alpha_1}\beta_2^{\alpha_2}}{(\beta_1+\beta_2)^{\alpha_1+\alpha_2-1}}\frac{\Gamma(\alpha_1+\alpha_2-1)}{\Gamma(\alpha_1)\Gamma(\alpha_2)}\,.$$
which can be simplified a little, but which I won't pursue further here as I don't think that's what you're after.
