How to solve ergodic density/distribution using R

May I have a question on how to solve the following problem by R code?

Mainly we want to solve the equation show in the attached image.  The equation is a continuous version of Markov process.  

In the equation, we have been able to achieve two things using R code:
[1]  From year-2009 sample data, we can estimate the marginal density ?f(x ; 2009)? by using R function ?density()?
[2]  From appropriately grouped all sample data, we can estimate the conditional density ?g(z|x ; 1)? by using R function ?cde()? from hdrcde package

Now we are searching and reading reference docs on how to solve for the long-run (ergodic) distribution ?f(x)?.  Much appreciated if any suggestions on solution steps using R; and we will be happy to provide additional references.  Thanks a lot!

Regards,

Hong Yu