For a group of random variables $X_{i:\,0\leq i < k}$, I am interested in finding the liklihood $f(x\in X_i)=\text{P}(x|\text{model}$) given that $X_i \sim \mathcal{N}(\mu_i+\Theta_i, e^{\beta_0}(\mu_i+\Theta_i)^{\beta_1})$, $\Theta \sim \text{Exp}(\phi)$. Or, in other words, an exponentially modified gaussian with a mean-variance relationship modeled with two parameters $\beta_0,\beta_1$ s.t. $\log\text{Var}(X) = \beta_0+\beta_1\log\mathbb{E}(X)$. This gives the marginal liklihood:$$f(x) = \frac{1}{\phi\sqrt{2\pi}}\int_0^\infty e^{-\frac{1}{2}\beta_0}(\mu_i+\Theta_i)^{-\frac{1}{2}\beta_1}\exp\left(-\frac{(x-(\mu_i+\Theta_i))^2}{e^{\beta_0}(\mu_i+\Theta_i)^{\beta_1}}-\frac{1}{\phi}\Theta_i \right)d\Theta_i $$ Currently I use numerical integration, but I am wondering if there is a canonical form that could speed up the integration (such as in terms of hypergeometric functions)? Any suggestions, including suggestions on algorithms that avoid integration all together, are welcome.
EDIT: While writing this, I found this post about the integral of a polynomial exponential, which I adapted to my particular rational function to yield: $$\frac{1}{\phi\sqrt{2\pi}}e^{-\frac{\mu_i}{\phi}}\sum_k^\infty \left[\sum_{\substack{p \\2n+m=k}}^\infty\frac{\phi^N\Gamma(N)}{(-2)^{n+p}e^{(n+m+p)\beta_0}n!m!p!}\right]x^k$$ where $N=m+2p+1-(n+m+p+\frac{1}{2})\beta_1$. Not sure if this actually converges, or if it is any better than numerical integration, though if the series converges for at least $|x| < 1,\phi < 1,$ and is rapid enough, it might be. It would also be nice to have a closed form for the coefficients. They seem to converge only when $\beta_1 > 1$, even though the integral should converge for all values of $\beta_1$. Using $2-\beta_1$ as the highest power when evaluating the integral might fix that.
