I'm struggling to find the density function of $Y=X^3$ where $X \sim \mathcal{N}(0,1)$ with density function $\phi(\cdot)$.
Since $g(X) = X^3$ is monotonic I suppose that the "change of variables" formula $$f_Y(y) = f_X(g^{-1}(y)) \cdot \biggr\lvert\frac{d}{dy}g^{-1}(y) \biggr\rvert$$ is applicable here? With $g^{-1}(Y)=Y^{1/3}$ and $f_X(x) = \phi(x)$ this yields $$f_Y(y) = \phi(y^{1/3}) \cdot \biggr\lvert\frac{y^{-2/3}}{3}\biggr\rvert = \frac{\phi(y^{1/3})}{3y^{2/3}}$$ which is only defined for $y>0$ and does not integrate to unity so it can't be correct.
What am I doing wrong?
One can directly plug in the formula with the absolute value of Jacobian $|J|$ only when $J$ is finite. At $y = 0$ the Jacobian blows up, corresponding to the horizontal slope of $Y = X^3$, which is infinite slope in inverse transform $X = g^{-1}(Y)$.
In this kind of situation where you start out knowing $g(X)$ being monotonic (thus invertible), you should just split the domain (support).
Formally, re-define the (forward) transformation such that \begin{align} Y &\equiv \begin{cases} \hphantom{8pt}X^3 & \text{for}~~ X > 0 \\ -(-X)^3 & \text{for}~~ X < 0 \end{cases} &&\implies &X &\equiv \begin{cases} \hphantom{-}Y^{\frac13} & \text{for}~~ Y > 0 \\ -(-Y)^{\frac13} & \text{for}~~ Y < 0 \end{cases} \end{align} For the missing point at zero, note that at the natural limit $X \to 0$ (either from the left $0^{-}$ or from the right $0^{+}$), or equivalently at $Y \to 0$, the transformation is the identity function $Y = X$. Hence one will define $Y = g(X) \equiv X$ at zero as such, and have the Jacobian be unity by construct. The blowing-up issue is explicitly resolved.
Now you can plug in the formula (as correctly stated in the question post) for each of the two pieces respectively and (necessarily) conclude that the density $f_Y$ is symmetric with respect to zero.