Let $\Delta$ denote the forward difference operator on a sequence: $$\Delta s_n = s_{n+1} - s_n,$$ and $\Delta^m$ denote the forward difference of the order $m$: $$\Delta^0 s_n = s_n, \quad \Delta^{m+1} s_n = \Delta\left(\Delta^m s_n\right).$$ We say that a sequence $s_n$ is completely monotone iff $$(-1)^m \Delta^{m} s_n > 0,$$ i.e. the sequence $s_n$ itself is positive and decreasing, its first differences $\Delta s_n$ are negative and increasing (that is, decreasing in absolute values), its second differences $\Delta^2 s_n$ are positive and decreasing (like the sequence $s_n$ itself), and so on. So, the sequence of second differences of a completely monotone sequence is also completely monotone (and so is the sequence of its differences of any even order).
Let $a$ be a real number in the interval $1<a<e^{1/e}$. Let ${^n a}$ denote tetration: $${^0 a}=1, \quad {^{(n+1)} a} = a^{\left({^n a}\right)}.$$ The sequence ${^n a}$ is strictly increasing and converges to the limit $${^\infty a} = \lim_{n\to\infty}\left({^n a}\right) = \frac{W(-\ln a)}{-\ln a},$$ where $W(x)$ is the Lambert $W$-function. Let $s_n = {^\infty a} - {^n a}$. It appears that the sequence $s_n$ is completely monotone. How can we prove it?
So from my comments above, $\lambda$ is the multiplier at the fixed point. $a^{L+x}=L+\lambda x + O x^2\;\;\ln_a(L+x) = L +\frac{x}{\lambda} + O x^2 $ $$L= \frac{W(-\ln a)}{-\ln a}\;\;\;\lambda=\ln(L)$$
For the bases in question $1<L<e,\;\;0<\lambda<1$. Then an approximate solution for the superfunction for large values of z is $$S_1(z)=L-\lambda^z=L-e^{z \ln(\ln(L))}$$
$S_1(x)$ is a completely monotonic analytic function at the real axis, where the odd derivatives are all >0, and the even derivatives are all <0, which we will use as the definition of a completely monotonic analytic function. This is because $\ln(\ln(L))<0$ for all bases in question; and we substitute that into the Taylor series for exp(z). $\ln_a(x)$ is also completely monotonic from the Taylor series of $\ln(x)$, for the Op's bases, in the domain $0<x<L$, which is the domain we will use in the iterations below. From that initial $S_1(z)$ estimate, we generate a sequence of better and better estimates, which converges to the Schröder function solution by lemma1 which is described below.
$$S_n(z)=\log_a(S_{n-1}(z+1))$$ We can show that $S_n(x)$ is completely monotonic if $Re(x+1)>S^{-1}_{n-1}(0)>$; by lemma2, which is that the composition of two completely monotonic analytic functions is completely monotonic in the appropriate range. In the limit, $S(x)$ is completely monotonic at the real axis if $Re(x+1)>S^{-1}(0)$
$$S(z) = \lim_{n\to\infty} S_n(z)\;\;\;S(z+1)=a^{S(z)} \;\;\; \lim_{x\to\infty}S(x)=L$$
lemma 1, the sequence of $S_n(z)$ functions; $S(z)$ converges to the Superfunction generated from the Schröder equation. I would prove this by generating error estimates from the Superfunction generated from the formal inverse Schröder equation which is periodic; $p=\frac{2\pi i}{\ln{\lambda}}$ so $S(z)$ can be expressed as an analytic Fourier series $S(z)=L+\sum a_n(-1)^n\lambda^{n z}$ where $a_1=1$ and $a_n$ are the formal inverse Schröder function Taylor series coefficients and $S(z)$ converges if $\Re(z+1)>S^{-1}(0)\;\;$ end lemma 1
Then we shift by a constant to get the desired sexp function. $$\text{sexp}_a(z)=S(z+k)\;\;\;k=S^{-1}(0)+1\;\;\;S(k)=1$$
Finally, $\text{sexp}_a(z)$ is completely monotonic at the real axis if z>-2. Then the sequence $^n a=\text{sexp}_a(n)$ is also a completely monotonic sequence by lemma3 that if an analytic function is completely monotonic over a given range, then a sequence of equally spaced samples of that function in that range is also completely monotonic.
So, if we prove lemma(1,2,3), then we have the Ops result. Lemma1 would come from a good book on complex dynamics, Lennart Carlson's book, would be the place to start, but I have seen it mentioned before. The other two Lemmas should be well known in the theory of monotonic functions. Real valued Schröder superfunction solutions will often be fully monotonic; so presumably this is known.
There is a possible alternative approach to lemma1, that $S(x)$ is completely monotonic. As noted above, the formal superfunction solution can be expressed as a Fourier series: $S(z)=L+\sum a_n(-1)^n\lambda^{n z}$ Then $a_1=1;$ and by observation $a_n(-1)^n<0 \;\forall n$. This sequence leads immediately to a completely monotonic solution for $S(z)$ at the real axis although I haven't tried proving the observation.