I have to calculate expectations of the form $\int_0^\infty f(x)\frac{x^{k-1}e^{-\frac{x}{\theta}}}{\theta^k\Gamma(k)}dx$, basically $\mathbb{E}[f(x)]$ where $x\sim\text{Gamma}(k,\theta)$. Unfortunately many of these integrals don't have a closed form solution, so I reside to using a power series expansion for $f(x)=\sum_{n=0}^{\infty}a_nx^n$ and solve this integrals through moments of gamma distribution. However I want to write this results in terms of well known generalized functions like generalized hypergeometric function or Meijer-G function.
If we know the fact that for any $\alpha$, we have $\mathbb{E}[x^\alpha] = \frac{\theta^\alpha\Gamma(k+\alpha)}{\Gamma(k)}$ and knowing $k,\alpha\ge 0$, we can write
$$\mathbb{E}[f(x)] = \sum_{n=0}^{\infty} a_n\frac{\theta^n\Gamma(k+n)}{\Gamma(k)}$$
How can I write this result in terms of generalized hypergeometric function or Meijer-G function?
Furthermore maybe in general this is not possible, when it is possible to do this, for example if $f(x)=e^x$, $f(x)=\ln x$ or $\frac{ax+b}{cx+d}$ and similar functions.
Thank you!