The main theme here is the attempt to generalize the Isserlis-Wick Theorem for Gaussian integrals to the case of exponentials of polynomials of higher degree.
Let
$$
P(x)=a_n x^n+\cdots+a_1 x+a_0
$$
be a polynomial of even degree $n$ with $\Re\ a_n>0$. Then the integral
$$
J=\int_{-\infty}^{\infty}e^{-P(x)}\ dx
$$
converges.
One can expand $\exp(-\sum_{i=1}^{n-1}a_i x^i)$ inside the integral and factor out the constant $e^{-a_0}$
$$
J=e^{-a_0}\int_{-\infty}^{\infty}e^{-a_n x^n}\sum_{r_1,\ldots,r_{n-1}=0}^{\infty}
\prod_{i=1}^{n-1}\frac{(-a_ix^i)^{r_i}}{r_i!}
\ \ dx
$$
and also take the sum out
$$
J=e^{-a_0}\sum_{r_1,\ldots,r_{n-1}=0}^{\infty}
\prod_{i=1}^{n-1}\frac{(-a_i)^{r_i}}{r_i!}
\int_{-\infty}^{\infty}e^{-a_n x^n} x^{r_1+\cdots+r_{n-1}}\ dx
$$
using absolute convergence.
Finally one is reduced to integrals of the form
$$
\int_{-\infty}^{\infty}e^{-a_nx^n} x^r\ dx
$$
which vanish if $r$ is odd and otherwise are equal to
$$
2\int_{0}^{\infty} e^{-a_nx^n}x^{r+1}\ \frac{dx}{x}=\frac{2}{n}a_n^{-\frac{r+1}{n}}\int_{0}^{\infty}
e^{-t}t^{\frac{r+1}{n}}\frac{dt}{t}=\frac{2}{n}a_n^{-\frac{r+1}{n}}
\Gamma(\frac{r+1}{n})
$$
after the change of variables $t=a_nx^n$.
For a no-go result about higher degree Isserlis-Wick theorems see the article "When is the Fourier transform of an elementary function elementary?" by Etingof, Kazhdan, and Polishchuk in Selecta Math 2002.
Also, things work well if one expands the interaction (in QFT lingo) lower degree terms using the leading monomial as the free reference measure. If one does it the other way as suggested in Igor's comment then exchanging summation and integration is not allowed. Nevertheless, one should still be able to recover the integral from the series using Borel-Leroy summation, see the article "Loop vertex expansion for $\Phi^{2k}$ theory in zero dimension" by Rivasseau and Wang and the follow up "Note on the intermediate field representation of $\phi^{2k}$ theory in zero dimension" by Lionni and Rivasseau.
