Thursday, March 09, 2006

Re: st: Producing random numbers from a kernel density

Another popular method is by acceptance-rejection algorithm: if your density is h(x) dominated by a functions like Kf(x) where K>1 is a relatively well known constant, and f(x) is an exactly known distribution from which it is easy to generate, then A-R algorithm proceeds as follows: (1) generate x from f(x) (2) generate an independent U~U[0,1] (3) check if K*U*f(x) is smaller than g(x) (note that you are using your density of interest here). If it is, then keep it, and this will be a new number drawn from your distribution; if not, discard it and go for a new pair of numbers in (1) and (2).

With this method, you may have to re-estimate your kernel density for each x. That's going to be a mess.

You might want to check Monahan, Numerical Methods of Statistics, and/or Liu, Monte Carlo Strategies in Scientific Computing.

I think yet another option you might have is to approximate your density by a mixture of known densities, and then generate easily from that mixture. I wrote -denormix- package for normal mixtures in the ice age of Stata 6, if not Stata 5, but it is only intended for estimation, not for Monte Carlo generation.

On 3/9/06, Marcello Pagano <pagano@hsph.harvard.edu> wrote: > Depends how many numbers you want to generate, how much accuracy you > want etc., but a quick and (not so) dirty method is to (1) make sure you > have a decent density estimate (such as eqprhistogram ;-) ) that is > always non-negative and integrates to one, then (2) calculate the > distribution function of this density at as fine a grid as you want, > then (3) generate a uniform random number, and (4) use the probability > integral transform (which says that your distribution function is > uniformly distributed) to get your variate Voila. > > m.p. > > > Peter Willemé wrote: > > > Dear all, > > I am trying to draw random numbers from an arbitrary distribution. One > > idea that occurred to me was to use a kernel density estimator > > (-kdensity-), but I have no idea how to produce random numbers from > > the result. All suggestions are welcome. > > Best, > > Peter > > > > > * > * For searches and help try: > * http://www.stata.com/support/faqs/res/findit.html > * http://www.stata.com/support/statalist/faq > * http://www.ats.ucla.edu/stat/stata/ >

-- Stas Kolenikov http://stas.kolenikov.name

* * For searches and help try: * http://www.stata.com/support/faqs/res/findit.html * http://www.stata.com/support/statalist/faq * http://www.ats.ucla.edu/stat/stata/


Tag:


Links to this post:

Create a Link



<< Home

This page is powered by Blogger. Isn't yours?