Solve \begin {align*} y^{\prime } & =4\left ( yx\right ) ^{\frac {1}{3}}\\ y^{\prime } & =\omega \left ( x,y\right ) \end {align*}
The symmetry condition results in the pde\begin {equation} \eta _{x}+\omega \left ( \eta _{y}-\xi _{x}\right ) -\omega ^{2}\xi _{y}-\omega _{x}\xi -\omega _{y}\eta =0 \tag {1} \end {equation} Trying polynomial anstaz \begin {align*} \xi & =a_{0}+a_{1}x\\ \eta & =b_{0}+b_{1}y \end {align*}
And substituting these into (1) and simplifying gives\[ \left ( -16a_{1}+8b_{1}\right ) yx-4xb_{0}-4ya_{0}=0 \] Setting all coefficients to zero gives\begin {align*} -16a_{1}+8b_{1} & =0\\ b_{0} & =0\\ a_{0} & =0 \end {align*}
Hence \(a_{1}=\frac {1}{2}b_{1}\). Letting \(b_{1}=1\) then \(a_{1}=\frac {1}{2}\) and the infinitesimals are \begin {align*} \xi & =\frac {1}{2}x\\ \eta & =y \end {align*}
The integrating factor is therefore \begin {align*} \mu \left ( x,y\right ) & =\frac {1}{\eta -\xi \omega }\\ & =\frac {1}{y-\frac {1}{2}x\left ( 4\left ( yx\right ) ^{\frac {1}{3}}\right ) }\\ & =\frac {1}{y-2x\left ( xy\right ) ^{\frac {1}{3}}} \end {align*}
The next step is to determine the canonical coordinates \(R,S\). This is done by using the standard characteristic equation by writing\begin {equation} \frac {dx}{\xi }=\frac {dy}{\eta }=dS\nonumber \end {equation} The first pair of equations gives\[ \frac {dy}{dx}=\frac {\eta }{\xi }=\frac {2y}{x}\] Solving gives\[ y=c_{1}x^{2}\] Hence \begin {equation} R=c_{1}=\frac {y}{x^{2}} \tag {2} \end {equation} And \(S\) is found from \[ dS=\frac {dx}{\xi }=2\frac {dx}{x}\] Integrating gives\begin {align*} S & =2\ln x+c_{1}\\ & =2\ln x \end {align*}
By choosing \(c_{1}=0\). Now the ODE \(\frac {dS}{dR}=F\left ( R\right ) \) is found from\begin {align*} \frac {dS}{dR} & =\frac {\frac {dS}{dx}+\omega \left ( x,y\right ) \frac {dS}{dy}}{\frac {dR}{dx}+\omega \left ( x,y\right ) \frac {dR}{dy}}\\ & =\frac {S_{x}+\omega \left ( x,y\right ) S_{y}}{R_{x}+\omega \left ( x,y\right ) R_{y}} \end {align*}
But \(S_{x}=\frac {2}{x},R_{x}=-2\frac {y}{x^{3}},S_{y}=0,R_{y}=\frac {2}{x^{2}}\). Substituting these into the above and simplifying gives\begin {align*} \frac {dS}{dR} & =\frac {x^{2}}{2x\left ( yx\right ) ^{\frac {1}{3}}-y}\\ & =\frac {1}{2\frac {1}{x}\left ( yx\right ) ^{\frac {1}{3}}-\frac {y}{x^{2}}}\\ & =\frac {1}{2y^{\frac {1}{3}}x^{-\frac {2}{3}}-\frac {y}{x^{2}}}\\ & =\frac {1}{2\left ( \frac {y}{x^{2}}\right ) ^{\frac {1}{3}}-\frac {y}{x^{2}}}\\ & =\frac {1}{2\left ( R\right ) ^{\frac {1}{3}}-R} \end {align*}
Hence\[ \frac {dS}{dR}=\frac {1}{2R^{\frac {1}{3}}-R}\] Which is a quadrature. Solving gives\begin {align*} \int dS & =\int \frac {1}{2R^{\frac {1}{3}}-R}dR\\ S & =-\frac {3}{2}\ln \left ( -2+R^{\frac {2}{3}}\right ) +c_{1} \end {align*}
Converting back to \(x,y\) gives\[ 2\ln x=-\frac {3}{2}\ln \left ( -2+\left ( \frac {y}{x^{2}}\right ) ^{\frac {2}{3}}\right ) +c_{1}\] The above can be simplified more if needed to solve for \(y\left ( x\right ) \) explicitly.