How do you find shape parameters of a Johnson distribution from skew and kurtosis? I'm trying to generate random numbers from a distribution with given moments, in my field the Johnson family of distributions is typically used but I can't find anywhere how to set it's shape parameters to give specific population moments?
If anyone knows about any code that does it, that would be ideal but I would be happy with a book recommendation.
Thanks
 A: There is this paper by Shenton and Bowman (1975) titled Johnson's $S_U$ and the Skewness and Kurtosis Statistics. In it, the authors provide the third and fourth moments of Johnson's $S_U$ distribution.
Following the second section therein... 
Define the density as follows $$ p(y) = \frac{\delta}{\sqrt{2\pi}} \frac{\exp(-\frac{1}{2} z^2)}{\sqrt{1+y^2}} \, , $$ with $ z = \gamma + \delta \ln \left(y + \sqrt{1+y^2}\right) $.
Further, define $\ln \omega = 1 /\delta^2$ and $\Omega = \gamma / \delta$, then the third and fourth moments are 


*

*Skewness. $$ \mu_3 = -(\omega -1)^2 \sqrt{\omega}\left[\omega(\omega + 2)\sinh(3\Omega)\, + \,3\sinh(\Omega) \right] /4 \,.  $$

*Kurtosis. $$\mu_4 = (\omega -1)^2 \left[d_4 \cosh(4\Omega)\, + \, d_2 \cosh(2\Omega) + d_0   \right], $$ with $d_4 = \omega^2(\omega^4  +  2\omega^3 + 3\omega^2 - 3)/8  ,\, d_2 = \frac{1}{2}\omega^2(\omega + 2)  ,\, d_0 = \frac{3}{8}(2\omega + 1)$.


With regard to generating a random sample, these R packages may be useful:


*

*https://www.rdocumentation.org/packages/SuppDists/versions/1.1-9.4/topics/Johnson

*https://www.rdocumentation.org/packages/ExtDist/versions/0.6-3/topics/JohnsonSU
