I am trying to compute the probability density function of the forward rate implied by the SABR formula approximation in order to see how the density implied by the approximation has negative probabilities or explosive behaviour depending on the input parameters (especially for low strikes).
All I have is the stochastic SABR process and the closed form solutions with the first few terms of the asymptotic expansion of the term lognormal volatility.
How would I go about computing the PDF for a range of points? What I want to achieve is something like this: