[R] Pseudo-random numbers between two numbers
Berwin A Turlach
berwin at maths.uwa.edu.au
Thu Mar 12 14:52:01 CET 2009
On Wed, 11 Mar 2009 06:57:02 -0400
Duncan Murdoch <murdoch at stats.uwo.ca> wrote:
> guox at ucalgary.ca wrote:
> > Please forget the last email I sent with the same subject.
> > =================
> > I would like to generate pseudo-random numbers between two numbers
> > using R, up to a given distribution,
> > for instance, norm.
> > That is something like
> > rnorm(HowMany,Min,Max,mean,sd) over rnorm(HowMany,mean,sd).
> > I am wondering if
> >
> > qnorm(runif(HowMany, pnorm(Min,mean,sd), pnorm(Max,mean, sd)),
> > mean, sd)
> >
> > is good.
This is the approach taken in
http://www.jstatsoft.org/v16/c02
http://www.jstatsoft.org/v16/c02/paper
which you may want to consult. They give general routines to truncate
arbitrary distributions.
> It depends on what Min and Max are. If you get far out in the tails,
> rounding error will kill you. For example, pnorm(x) is exactly 1 for
> x bigger than 10 or so, so this approach would fail if Min and Max
> were both bigger than 10.
I agree with this assessment and guess that similar problem will arise
with other truncated distributions created with the code from the paper
above.
And, note, that if sd is small, you can easily be in the situation
where you evaluate pnorm() at extreme values.
> The solution is to switch to lower=FALSE in the upper tail, and
> possibly switch to a log scale if you want to be really extreme.
Or use package truncnorm. Though their c code seem to sample with
rejection from a normal proposal which would be quite inefficient for
"extreme" truncations, e.g., standard normal truncated to [8,10].
Or use package msm, whose rtnorm() implements a more efficient
algorithm for simulating from a truncated normal distribution.
Cheers,
Berwin
More information about the R-help
mailing list