I am going through the theory behind mixing distributions, and I came across a portion in an essay where it talks about mixing two gamma distributions, it does not show how to get the resulting distribution. But I am interested if you could somehow.
So the question goes like this:
We have $S|\Lambda \sim \Gamma(\alpha,\Lambda) $ and $\Lambda\sim\Gamma(\beta, \mu) $, we want to find $S$ such that the pdf of S is: $$ f_S(s) = \int_0^\infty f_{S|\Lambda=\lambda}(s) f_\Lambda(\lambda)d\lambda $$ So far I have that: $$ f_S(s) = \int_0^\infty ( \frac {s^{\alpha-1}e^{-s/\lambda}}{\lambda^{\alpha}\Gamma(\alpha)} ) ( \frac {\lambda^{\beta-1}e^{-\lambda/\mu}}{\mu^{\beta}\Gamma(\beta)} )d\lambda = \frac{s^{\alpha-1}}{\mu^\beta\Gamma(\alpha)\Gamma(\beta)} \int_0^\infty \lambda^{\beta-\alpha-1}e^{-(s/\lambda+\lambda/\mu)}d\lambda $$ I tried various substitutions and manipulations etc, can't figure how to evaluate this integral.
If you don't mind special functions: $$ f_S(s)=\frac{2} {\Gamma (\alpha ) \Gamma (\beta )s}\left(\frac{s}{\mu}\right)^{\frac{\alpha+\beta}{2}} K_{\alpha -\beta }\left(2 \sqrt{\frac{s}{\mu }}\right), $$ where $K$ is a modified Bessel function of the second kind.