I'm stuck with this convolution integral ($z \geq 0$)... \begin{equation} f_{Z}(z)=\int^{\infty}_{-\infty}f_{1}(x)f_{2}(z-x)dx = \mbox{ } ??? \end{equation} which represents the pdf of the sum $Z = X_1 + X_2$ of two random variables $X_1$ and $X_2$. Each $X_1$ and $X_2$ belongs to the same family (each of them is proportional to a non-central $\chi^2$ variable): \begin{equation} \begin{aligned} f_{1}(x) &= \alpha_1 e^{-\beta_1 x} (\gamma_1 x)^{\frac{1}{2} \nu} I_{\nu}(\delta_1 \sqrt{x}) \times \mathcal{I}_{(0,\infty)}(x)\\ f_{2}(x) &= \alpha_2 e^{-\beta_2 x} (\gamma_2 x)^{\frac{1}{2} \mu} I_{\mu}(\delta_2 \sqrt{x}) \times \mathcal{I}_{(0,\infty)}(x) \end{aligned} \end{equation} where $\mathcal{I}_{(0,\infty)}(\cdot)$ is the indicator function of the positive real axis and $I_{\nu}(\cdot)$ and $I_{\mu}(\cdot)$ are modified Bessel functions of the first kind of orders $\nu \geq 0$ and $\mu \geq 0$, respectively. Moreover, coefficients are all nonnegative: $\alpha_i \geq 0, \beta_i \geq 0, \gamma_i \geq 0, \delta_i \geq 0$.
My questions:
1) can anybody provide a closed-form expression for the integral? Off course I know I can easily use Fourier, but I'm really wondering whether such pdf can be explicitly write down.
2) Moreover, is it possible to write $f_{Z}(z)$ as $f_1$ and $f_2$, i.e. do (?) $\left\{\alpha_z, \beta_z, \gamma_z, \delta_z, q \right\}$ exist such that \begin{equation} f_{Z}(z) = \alpha_z e^{-\beta_z z} (\gamma_z z)^{\frac{1}{2} q} I_{q}(\delta_z \sqrt{z}) \end{equation}
3) I'm separately interested in the case: $\nu = \mu = 0$.
P.S. This question arises in this way: in a famous paper by Cox-Ingersoll-Ross (1985), there is a positive stochastic process $R(t)$ satisfying a square root SDE $$ dR = k(\theta - R)dt + \sigma \sqrt{R} dW $$ where $\alpha \geq 0$, $\theta \geq 0$, $\sigma \geq 0$ and $W_t$ is a standard Wiener process. Well, it can be shown that the conditional transition density (for time $s \geq t$) takes the same form above $$ f(R(s),s|R(t),t) = ce^{-u-v}\left(\frac{v}{u}\right)^{q/2}I_{q}(2\sqrt{uv}) $$ where \begin{equation} \begin{aligned} c &=\frac{2k}{\sigma^{2}\left(1-e^{-k(s-t)}\right)}\\ u &=cR(t)e^{-k(s-t)}\\ v &=cR(s)\\ q &=\frac{2k\theta}{\sigma^2}-1 \end{aligned} \end{equation} and $I_{q}(\cdot)$ is a modified Bessel function of the first kind of order $q$. Moreover it can be shown that: $2 c R(s) \sim \chi^2(2q + 2, 2u)$ i.e. the rescaled variable $2cR(s)$ distributes as a non central $\chi^2$ distribution with $2q+2$ degrees of freedom and parameter of non centrality $2u$. Note that, since the non central $\chi^2$ distribution is not closed w.r.t. the rescaling of the variable, $R(s)$ itself does not distribute according to a non central $\chi^2$. Therefore, my question is about the transitional distribution of the sum two independent processes $R_1$ and $R_2$ following \begin{equation} \begin{aligned} R &= R_1 + R_2 \\ dR_1 &= k_1(\theta_2 - R_1)dt + \sigma_1 \sqrt{R} dW_1 \\ dR_2 &= k_2(\theta_2 - R_2)dt + \sigma_2 \sqrt{R} dW_1 \\ \end{aligned} \end{equation} with $Corr[W_1, W_2] = 0$ and each having transitional density $f_{R_i}$. and therefore \begin{equation} f_{R}(r)=\int^{\infty}_{0}f_{R_1}(x)f_{R_2}(r-x)dx \end{equation} Off course $f_{R}(r)$ can be easily computed by numerical Fourier inversion, but as I said before, I would be courious about a closed form expression for the pdf.
Thanks for your attention.