$\newcommand{\Reals}{\mathbf{R}}$The "missing link" in my earlier sketch was the principle that if $R$ is a smooth function with $R(x) = O(|x|^{N})$ near $0$, and if we write $R(x) = x^{N} r(x)$, then $r$ itself is smooth. The simplest argument I could construct was more delicate than I'd originally expected, so I've included it in some detail. It's entirely possible there is a more elegant argument.
Lemma: If $R$ is a smooth, real-valued function on $\Reals$ such that $R(x) = O(|x|)$ near $0$, and if $R(x) = xr(x)$ for all real $x$, then $r$ is smooth.
Proof: Because $R$ is smooth, $r$ is smooth except possibly at $0$. Let $D$ denote the derivative operator and write, e.g., $D^{k}(x^{n})$ to denote the $k$ th derivative of the $n$th power function. If $n \geq 0$ is an integer, a straightforward induction (compare the identity $R^{(n)}(x) = xr^{(n)}(x) + nr^{(n-1)}(x)$, or $r^{(n)}(x) = \frac{1}{x}\bigl[R^{(n)}(x) - nr^{(n-1)}(x)\bigr]$, which holds for $x \neq 0$), gives the formula
$$
r^{(n)}(x) = \frac{1}{x^{n+1}} \sum_{k=0}^{n} (-1)^{k} D^{k}(x^{n}) D^{n-k}R(x).
\tag{1}
$$
For example,
\begin{align*}
r'(x) &= \frac{xR'(x) - R(x)}{x^{2}}, \\
r''(x) &= \frac{x^{2} R''(x) - 2x R'(x) + 2R(x)}{x^{3}}, \\
r'''(x) &= \frac{x^{3} R'''(x) - 3x^{2} R''(x) + 6x R'(x) - 6R(x)}{x^{4}}.
\end{align*}
The numerator (i.e., the sum) in (1) is a smooth function, and clearly vanishes at $0$. (The terms containing a positive power of $x$ vanish because $R$ is smooth, while $R(0) = 0$ by hypothesis.)
Upon differentiating this numerator and shifting the index of summation, everything cancels except one term:
\begin{align*}
D\sum_{k=0}^{n} (-1)^{k} D^{k}(x^{n}) D^{n-k}R(x)
&= \sum_{k=0}^{n} (-1)^{k} \bigl[D^{k+1}(x^{n}) D^{n-k}R(x) + D^{k}(x^{n}) D^{n-k+1}R(x)\bigr] \\
&= -\sum_{k=1}^{n+1} (-1)^{k} D^{k}(x^{n}) D^{n-k+1}R(x)
+ \sum_{k=0}^{n} (-1)^{k} D^{k}(x^{n}) D^{n-k+1}R(x) \\
&= x^{n} R^{(n+1)}(0).
\end{align*}
Consequently, applying l'Hôpital's rule to (1) gives
$$
\lim_{x \to 0} r^{(n)}(x)
= \lim_{x \to 0} \frac{1}{x^{n+1}} \sum_{k=0}^{n} (-1)^{k} D^{k}(x^{n}) D^{n-k}R(x)
= \lim_{x \to 0} \frac{x^{n} R^{(n+1)}(0)}{(n+1)x^{n}}
= \frac{R^{(n+1)}(0)}{n+1},\quad n\geq 0.
$$
To show $r$ is of class $C^{1}$, expand $R$ to third order at $0$:
$$
R(x) = R'(0)x + \frac{R''(0)}{2!}x^{2} + O(|x|^{3}),\qquad
r(x) = R'(0) + \frac{R''(0)}{2!}x + O(|x|^{2}).
$$
It follows at once that $r'(0)$ exists and is equal to $\frac{1}{2}R''(0) = \lim(r', 0)$.
Using this fact as a base case, assume inductively that $r$ is of class $C^{n}$ for some $n \geq 0$. The function $r^{(n)}$ is continuous on $\Reals$, differentiable except possibly at $0$, and $\lim(r^{(n+1)}, 0)$ exists. It follows that $r^{(n)}$ is differentiable at $0$, and that
$$
r^{(n+1)}(0) = \lim_{x \to 0} r^{(n+1)}(x).
$$
This establishes that $r$ is of class $C^{n+1}$ at $0$. By induction on $n$, $r$ is smooth.
Theorem: If $R$ is a smooth function such that $R(x) = O(|x|^{N})$ near $0$ for some integer $N \geq 1$, and if $R(x) = x^{N}r(x)$, then $r$ is smooth.
Proof: Immediate from the lemma. (We have $R(x) = O(|x|)$, so $x^{N-1}r(x) = O(|x|^{N-1})$ is smooth, etc.)
Armed with the theorem (which I'd blithely been assuming), we turn to the question at hand.
Suppose $f$ is smooth, and that $f^{(n)}(0) = 0$ when $n$ is not a multiple of $3$. For each positive integer $N$, $f$ has an order-$3N$ Taylor expansion
$$
f(x) = \sum_{k=0}^{N} \frac{f^{(3k)}(0) x^{3k}}{(3k)!} + R_{3N}(x)
\tag{2}
$$
with $R_{3N}$ smooth, and $R_{3N}(x) = O(|x|^{3N + 3})$ at $0$. Write $R_{3N}(x) = x^{3N + 3} r(x)$. By the theorem, $r$ is smooth.
Define $g$ by $g(x) = f(x^{1/3})$ for all real $x$. Fix an integer $N \geq 1$, and use (2) to write
$$
g(x) = \sum_{k=0}^{N} \frac{f^{(3k)}(0) x^{k}}{(3k)!} + R_{3N}(x^{1/3})
= \sum_{k=0}^{N} \frac{f^{(3k)}(0) x^{k}}{(3k)!} + x^{N+1} r(x^{1/3}).
$$
To show $g$ is smooth, it suffices to show the remainder $R(x) := R_{3N}(x^{1/3}) = x^{N+1} r(x^{1/3})$ is of class $C^{N}$ (since $N \geq 1$ is arbitrary).
Since $r$ is smooth, $R$ is smooth except possibly at $0$. For $x \neq 0$,
\begin{align*}
R'(x)
&= (N + 1)x^{N} r(x^{1/3}) + x^{N+1} r'(x^{1/3}) \tfrac{1}{3} x^{-2/3} \\
&= x^{N}\bigl[(N + 1) r(x^{1/3}) + \tfrac{1}{3} x^{1/3} r'(x^{1/3})\bigr] \\
&= x^{N} r_{1}(x^{1/3}),
\end{align*}
with the final equality serving to define $r_{1}$. Clearly, $r_{1}$ is smooth, which (i) proves $R$ is of class $C^{1}$, and (ii) expresses $R'(x) = x^{N} r_{1}(x^{1/3})$, allowing us to continue differentiating inductively, with the result that $R^{(N)}(x) = x r_{N}(x^{1/3})$ for some smooth function $r_{N}$. Particularly, $R$ is of class $C^{N}$.
