I want to propose a (recursive) method to obtain 
$F_k(x)$ by exploiting a formula well known from
the theory of Borel re-summation. In the present
context it looks like this:

\begin{equation}
\begin{array}{ll}
\displaystyle \sum_{n \geq 0} \, \langle E_k \rangle_n \, z^n 
&\displaystyle = \
\int_0^\infty \, dt \,  e^{-t} \sum_{n \geq 0} \, 
{\langle E_k \rangle_n \over {n!}} \, (zt)^n  \\
&\displaystyle = \
\int_0^\infty \, dt \, e^{-t} \, F_k(zt) \\
&\displaystyle = \ \Big( {1 \over z} \Big) \cdot \,\text{$\frak{L}$} \big\{ F_k \big\} \Big({1 \over z} \Big) \quad \text{Laplace transform}
\end{array}
\end{equation}

I'll illustrate the approach for $k=1$. Let $E_1(x) := \sum_{n \geq 0} \, \langle E_1 \rangle_n \, x^n$. In this case the recurrence relation, as indicated the initial post above, is

\begin{equation}
\underbrace{n \langle E_1 \rangle_n}_{xE_1'(x) - x} 
\ = \  
\underbrace{\langle E_1 \rangle_{n-1}}_{xE_1(x)}
\ + \ 
\underbrace{(n-1) \langle E_1 \rangle_{n-2}}_{x^3E_1'(x) \ + \ x^2E_1(x)}
\ +  
\underbrace{\ n^2 \ }_{\sum_{n \geq 2} n^2 x^n}
\end{equation}

for $n \geq 2$ with the initial values $\langle E_1\rangle_0 = 0$ and $\langle E_1 \rangle_1 = 1$.
I've "underbraced" the terms of the recurrence by
their respective contributions to the $1$-st order
non-homogeneous ODE for the generating function $E_1(x)$
which is, upon simplification,


\begin{equation}
E´_1(x) \ + \ {1 \over {x-1}} \, E_1(x)
\ = \ {1 \over {(x-1)^4}}
\end{equation}

If I've done my calculus without any glitches
the general solution is of the form

\begin{equation}
\ {\text{const} \over {x-1}} 
\ - \ {1 \over 2} \, {1 \over {(x-1)^3}}
\end{equation}

whose constant term must be ${1 \over 2}$ since $\langle E_1\rangle_0 = 0$. Therefore

\begin{equation}
E_1(x) \ = \ {1 \over 2} \, {x(x-2) \over {(x-1)^3}}
\end{equation}

and after performing the requisite inverse transforms
we should (excuse me for mimicking Wikipedia) we get

\begin{equation}
\begin{array}{ll}
\displaystyle F_k(zt) 
&\displaystyle = \ \sum_{n \geq 0} \, {\langle E_1 \rangle_n \over {n!}} \, (zt)^n \\
&\displaystyle = \ {1 \over {2\pi}} \, \int_{-\pi}^{\pi} \,
d\theta \, E_1 \big(zt e^{-i\theta} \big) \exp \big( e^{i\theta} \big) \\
&\displaystyle \ \ \ \ \ \text{(Not sure what this is yet.)}
\end{array}
\end{equation}

ines.