$\newcommand\R{\mathbb R}\newcommand\C{\mathbb C}$You were almost there.
Indeed, suppose that $$f(x)=\sum_{j=0}^\infty a_j e^{-jx}\in\R \quad \forall x\in(0,\infty). \tag{1}\label{1}$$ Then $$g_f(s):=f(-\ln s)=\sum_{j=0}^\infty a_j s^j\in\R \quad \forall s\in(0,1).$$ So, the radius of convergence of the power series $\sum_{j=0}^\infty a_j z^j$ is $1$. So, $g_f$ can be extended to the analytic function $G_f$ on the open unit disc $D$ in $\C$ given by the formula $G_f(z)=\sum_{j=0}^\infty a_j z^j$ for $z\in D$.
Vice versa, if the function $g_f$ can be extended to an analytic function $G_f$ on $D$, then $G_f(z)=\sum_{j=0}^\infty a_j z^j$ with $a_j=G_f^{(j)}/j!$ and hence \eqref{1} holds.
Thus, \eqref{1} holds for some $a_j$'s if and only if the function $g_f$ can be extended to an analytic function on $D$.
Similarly, \eqref{1} holds for some $a_j$'s if and only if the function $g_f$ can be extended to a real-analytic function on the interval $(-1,1)$.