Gary wrote:
Stephen Walton wrote:
An application question for a change: I need to produce pseudorandom numbers drawn from a lognormal distribution (see http://www.itl.nist.gov/div898/handbook/eda/section3/eda3669.htm). Does anyone have existing code for this in Numeric or numarray?
It is interesting to me to see this question.
If I understand you correctly, I can say that Speakeasy has what you want. Actually, Speakeasy can produce random numbers from any user-defined distribution, even one that is not analytical. A colleague of mine once (in the early '90s) used that capability in modeling nuclear scintillation detectors. He took a *measured* pulse height distribution from a PMT and used it to generate random numbers that he input into a model of a scintillation crystal. His "Monte Carlo" simulation was 13 lines long. I was intrigued by the usefulness of this feature, so I always have my eye open to see if any other software package has that capability. I haven't searched exhaustively, but I haven't seen it in Matlab (or octave or scilab) or Mathematica or Mathcad. Or any of the usual python packages.
Well, it depends on how tricky the function is. If it's univariate, and you can write out the pdf or cdf as a function, then I believe you can subclass scipy.stats.rv_continuous, and it's rvs() method will numerically invert the cdf to generate it's random numbers. If the function is highly multivariate, you might need to do Markov-Chain Monte Carlo which is implemented by PyMC. If you have a bunch of data points from a continuous distribution, but no functional description, then you can make a kernel density estimate and draw random numbers from that. As of last night, that functionality is in scipy.stats.gaussian_kde. Resampling from discrete data is pretty trivial to handwrite. Is there anything else you need? -- Robert Kern rkern@ucsd.edu "In the fields of hell where the grass grows high Are the graves of dreams allowed to die." -- Richard Harter