Don't know if this will be of help, but if you change variables to polar coordinates, we get an ODE DSolve can handle:
subs = Simplify@NestList[ (* change of variables *) D[First@#, x] -> Dt[Last@#]/Dt[r[t] Cos[t]] &, y[x] -> r[t] Sin[t], 2]~Join~ {x -> r[t] Cos[t]}; (* new ode *) ode = -x*Derivative[1][y][x] - x*Derivative[1][y][x]^3 + y[x]*(1 + Derivative[1][y][x]^2) + y[x]^2*y''[x] + (x^2 - R*Sqrt[x^2 + y[x]^2])*y''[x] /. subs // Together // Numerator // Simplify[#, r[t] > 0] &; odes = FactorList[ode][[2 ;;, 1]] == 0 // Thread
(* the process yielded a spurious factor r[t] {r[t] == 0, r[t] (r[t] - R) r''[t] + (2 R - r[t]) r'[t]^2 + R r[t]^2 == 0} *)
dsol = DSolve[Last@odes, r, t]
(* {{r -> Function[{t}, InverseFunction[-I (ArcTanh[(R - #1)/Sqrt[ R^2 - 2 R #1 - C[1] #1^2]] - ArcTan[(-R - C[1] #1)/( Sqrt[C[1]] Sqrt[R^2 - 2 R #1 - C[1] #1^2])]/Sqrt[C[1]]) & ][t + C[2]]]}, {r -> Function[{t}, InverseFunction[ I (ArcTanh[(R - #1)/Sqrt[R^2 - 2 R #1 - C[1] #1^2]] - ArcTan[(-R - C[1] #1)/( Sqrt[C[1]] Sqrt[R^2 - 2 R #1 - C[1] #1^2])]/Sqrt[C[1]]) & ][t + C[2]]]}} *)
(To connect with my symmetry comments, the change of variables shows that the system is invariant under t -> t + t0, which is a rotation in polar coordinates.)
DSolvean ODE that is too hard to solve. (That's what the result indicates.) $\endgroup$Rterm by eitherxory[x], then the ODE is homogenenous, andDSolvecan find a first integral and return an implicit solution. I'm just pointing out it's close to an ODE that is not completely impossible, in case you made a mistake in entering the equation. $\endgroup$ParametricNDSolve[]). You can also useAsymptoticDSolveValue[], to get a series expansion about a point of interest, which may or may not be useful to you. $\endgroup$