Suppose that $X$ is gamma distributed with scale $s$ and shape $\alpha$, and that $P$ is a polynomial of degree $n$ : $$P(x)=\sum_{k=0}^n p_k x^k.$$
When it exists (e.g., when the coefficients $p_k$ are all negative), how could I approximate the expectation $$\mathbb E(e^{P(X)}) ?$$
Context: I am trying to catch the non-assymptotic distribution of the sample moments of a gamma random variable. Denoting $Y = (1,X,X^2,...X^n)$, this expectation is the moment generating function of $Y$, and i really want the distribution of $1/N \sum_{j=1}^N Y_j$ for i.i.d such $Y_j$'s, which led me to this question.
Trials: I tried to use a Taylor series expansion for the exponential, and then the bare moments of the gamma distributions, but this is not really effective as it is very combinatorial.
Hopes: the multivariate distribution of $Y$ is singular but simply defined, thus I hope we could get an efficient expansion for its mgf… but maybe I'm wrong :)