In models of single field inflation the Friedmann equation reads
$$ H(t)^2 = \frac{8\pi G}{3} \big(V(\varphi) + \frac{\dot{\varphi}^2}{2}\big) \tag{1}\label{eq1}. $$
In the slow-roll approximation one neglects the kinetic term, $\dot{\varphi}\ll V(\varphi)$, and then slow-roll parameters are usually defined in terms of derivatives of the potential, i.e., $\epsilon_{V} \propto (V'(\varphi)/V(\varphi))$ and so on. Then, the number of e-folds from a pivot time $t_\ast$ until the end of inflation reads
$$N_\ast =\int_{t_\ast}^{t_{\textrm{end}}} H(t)\, \textrm{d} t \approx 8\pi G \int_{\varphi_\ast}^{\varphi_{\textrm{end}}} \frac{V(\varphi)}{V'(\varphi)} \, \mathrm{d}\varphi \tag{2}\label{eq2}.$$
However, I'm interested in finding an analytical self-consistent expansion order-by-order of $N_\ast$, in terms of the Hubble-flow parameters
$$ \epsilon_{1}(t) \equiv - \frac{\dot{H}(t)}{H(t)^2}, \quad \epsilon_{n}(t) \equiv - \frac{\dot{\epsilon}_{n-1}(t)}{H(t) \epsilon_{n-1}(t)}\tag{3}\label{eq3}.$$
In other words, I'm looking for something like $N_\ast^{-1} = (\ldots)\, \epsilon_{1\ast} +(\ldots)\, \epsilon_{1\ast}^2+ \mathcal{O}(\epsilon_{1}^3)$. Why? Because I'd like to invert the expansion and find something like
$$\epsilon_{1\ast} = (\ldots) \, N^{-1}_\ast +(\ldots) \,N^{-2}_\ast + \mathcal{O}(N_{\ast}^{-3}) \tag{4}\label{eq4} ,$$
where $\epsilon_{1\ast}=\epsilon_{1}(t_\ast)$. By finding the coefficients of Eq. \eqref{eq4} one might be able to plug these results in other expressions of interest, e.g., higher corrections to the primordial power spectrum, that explicitly depend on $\epsilon_{1\ast}, \epsilon_{2\ast}$, etc (but not on $\epsilon_{V})$.
Is it possible -in general- to find an expansion of the form Eq. \eqref{eq4} for a given potential $V(\varphi)$? And if yes, does it depend on solving the Klein-Gordon equation $\ddot{\varphi}+3H(t)\dot{\varphi} + V'(\varphi)=0$, or is there any shortcut?
