Let us consider the following second-order stationary Fokker-Planck type PDE on $\mathbb R^2$ $$ \sigma\,\partial^2_xR(x,z)-\alpha \,x^{2q-1}\partial_xR(x,z)+\partial_z\left((z-x)\,R(x,z)\right)=0, $$ where $\sigma, \alpha > 0$ and $q=1,2,3,...$. I am looking for an analytical solution satisfying $$ \int_{\mathbb R}R(x,z) \, dz = 1, \quad R(x,z)>0 ,\qquad \qquad{(1)} $$ i.e., $R(x,z)$ is a pdf on $\mathbb R$ with respect to $z$ for all $x \in \mathbb R$. In case $q=1$, it is possible to verify that the solution is given by $$ R(x,z)=C\exp\left(-\frac{1}{2\sigma}\left(x-(1+\alpha )z\right)^2\right), $$ where $C$ is such that $R(x,z)$ satisfies (1). In order to find this, one can use the change of variables $(x,z)\to(y,z)$ where $y=x -(1+\alpha)z$, and then try for separable solutions, which works as the solution is only a function of $y$.
Is there a way to solve the equation for a general $q \geq 1$? I first tried tried to find a change of variables $(x,z)\to(x,y)$ such that the solution only depend on $y$, or is separable wrt $x$ and $y$, but failed. Then I tried to employ the Fourier transform, but failed nonetheless, as the term $x^{2q-1}$ introduces higher order derivatives wrt the Fourier variable.
Is there any approach which could be helpful?
EDIT: On this question a similar problem is treated, but the nonlinear coefficient seems to break down the method.