Analyzing the decay rate of Taylor series coefficients when high-order derivatives are intractable
Complex analysis can help. The rate of Taylor coefficients is determined by:
a) the radius of convergence, which is equal to the radius of the largest disk $|z|<r$ where your function is analytic. This radius is responsible for the exponential asymptotics, and
b) the nature of singularities on the circle $|z|=r$.
Your function $f$ can have only square root singularities, since $g$ has only square root singularities. Since the singularities of $g$ are $\pm1$, to determine the radius of convergence, you have to show that equation $g(z)=\pm1$ has no solutions in $|z|<1$. This will justify your asymptotics.
I have not done the calculation, perhaps you can do it yourself.
Ref: P. Flajolet and R. Sedgewick, Analytic combinatorics, Chap. VI.
Edit: Conrad gave a simple argument in the comment below below which shows that $f$ has no other singularities in the closed unit disk, except at $z=\pm1$, so your conjecture on the asymptotics is correct.
Note that $g(1)=g'(1)=1$ and for real $x\in(-1,1)$ \begin{equation*} g''(x)=\frac1{\pi\sqrt{1-x^2}}. \end{equation*} The map $z\mapsto1-z^2$ maps the set \begin{equation*} R:=\mathbb C\setminus(-\infty,-1]\setminus[1,\infty) \end{equation*} onto $\mathbb C\setminus(-\infty,0]$. So, the map \begin{equation*} R\ni z\mapsto h(z):=\frac1{\pi\sqrt{1-z^2}} \end{equation*} is analytic, and hence $g$ can be continued analytically to $R$ by the Taylor formula \begin{equation*} R\ni z\mapsto g(z):=g(1)+g'(1)(z-1)+\int_1^z(z-u)h(u)\,du \\ =z+\int_1^z(z-u)h(u)\,du. \tag{0} \end{equation*}
Take a real $r>1$ and let \begin{equation*} D_r:=\{z\in R\colon|z|<r,|\arg(z-1)|>\pi/4\}. \end{equation*} The main difficulty is to show that $g(D_r)\subseteq R$ for some $r>1$.
First here, take any real $t>0$. Then there is some real $u_t>0$ such that for all complex $z$ with $|z|\le1$ and $|\arg z|\ge t$ we have $|\frac1\pi+\frac z2|\le\frac1\pi+\frac12-u_t$, whence \begin{equation*} |g(z)|\le\Big|\frac1\pi+\frac z2\Big|+ \sum_{n=1}^\infty \frac{(2n-3)!!}{(2n-1)n!2^n \pi}|z|^{2n} \le\frac1\pi+\frac12-u_t+ \sum_{n=1}^\infty \frac{(2n-3)!!}{(2n-1)n!2^n \pi} =g(1)-u_t=1-u_t. \end{equation*} Since $g$ is analytic on $R$, it is uniformly continuous on any compact subset of $R$. So, there is some real $r_t>1$ such that $|g(z)|<1$ for all complex $z\ne1$ with $|z|\le r_t$ and $|\arg z|\ge t$. Also, it follows that $|g(z)|<1$ for all complex $z$ with $|z|<1$.
Thus, to prove that $g(D_r)\subseteq R$ for some $r>1$, it suffices to show that $\Im g(z)\ne0$ for all $z\ne1$ with $|z|\ge1$ and $|\arg(z-1)|>\pi/4$ that are close enough to $1$.
To see this, note that $h(u)\sim\frac1{\pi\sqrt2\,\sqrt{1-u}}$ as $u\to1$. Then (0) yields \begin{equation*} g(z)=z+c_1\cdot(z-1)^{3/2}, \tag{1} \end{equation*} where $c_1=c_1(z)$ converges to a nonzero complex number as $z\to1$. So, the conclusion that $\Im g(z)\ne0$ for all $z\ne1$ with $|z|\ge1$ and $|\arg(z-1)|>\pi/4$ that are close enough to $1$ follows, which does show that $g(D_r)\subseteq R$ for some $r>1$.
So, $f=g\circ g$ is analytic on $D_r$ for such an $r$.
Moreover, (1) implies
\begin{equation*}
f(z)=g(g(z))=g(z+c_1\cdot(z-1)^{3/2})\\
=z+c_1\cdot(z-1)^{3/2}+\hat c_1\cdot(z-1+c_1\cdot(z-1)^{3/2})^{3/2} \\
=z+\tilde c_1\cdot(z-1)^{3/2}, \tag{2}
\end{equation*}
where $\hat c_1=\hat c_1(z)\sim c_1(z)$ and $\tilde c_1=\tilde c_1(z)\sim2c_1(z)$ as $z\to1$.
Similarly, since $h(u)\sim\frac1{\pi\sqrt2\,\sqrt{1+u}}$ as $u\to-1$, the Taylor formula for $g(z)$ at $z=-1$ (together with the observation $g(-1)=g'(-1)=0$) yields $g(z)=c_{-1}\cdot(z+1)^{3/2}$, where $c_{-1}=c_{-1}(z)$ converges to a complex number as $z\to-1$. Therefore and because $g$ is analytic at $0$, we have
\begin{equation*}
f(z)=g(g(z))=c_0+\hat c_0\cdot g(z)
=c_0+\tilde c_0\cdot(z+1)^{3/2}, \tag{3}
\end{equation*}
where $c_0:=g(0)$, $\hat c_0=\hat c_0(z)=O(1)$, and $\tilde c_0=\tilde c_0(z)=O(1)$ as $z\to-1$.
Now we are finally in a position to use (with thanks to Alexandre Eremenko) Theorem VI.5 (with $\alpha=-3/2$, $\beta=0$, $\rho=1$, $r=2$, $\zeta_1=1$, $\zeta_2=-1$, $\mathbf D=D_r$, $\sigma_1(z)=z$, $\sigma_2(z)=c_0$), which yields the $n$th coefficient for $f$: \begin{equation*} [z^n]f(z)=O(n^{\alpha-1})=O(n^{-5/2}), \end{equation*} as you conjectured.
For an illustration, here is the set $\{g(z)\colon z\in R,|z|<2\}$:
