I would like to generate pseudo-random numbers between two numbers using R, up to a given distribution, for instance, rnorm. That is something like rnorm(HowMany,Min,Max,mean,sd) over rnorm(HowMany,mean,sd). I am wondering if dnorm(runif(HowMany, Min, Max), mean, sd) is good. Any idea? Thanks. -james
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. Any idea? Thanks. -james
On 10-Mar-09 23:01:45, 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. Any idea? Thanks. > -jamesIt would certainly work as intended! For instance: truncnorm<-function(HowMany,Min,Max,mean=0,sd=1){ qnorm(runif(HowMany, pnorm(Min,mean,sd), pnorm(Max,mean, sd)), mean, sd)} Sample <- truncnorm(1000,-1,2.5) hist(Sample,breaks=100) Hoping this helps, Ted. -------------------------------------------------------------------- E-Mail: (Ted Harding) <Ted.Harding at manchester.ac.uk> Fax-to-email: +44 (0)870 094 0861 Date: 10-Mar-09 Time: 23:39:13 ------------------------------ XFMail ------------------------------
I have modified my example to make it more convincing! See at end. On 10-Mar-09 23:39:17, Ted Harding wrote:> On 10-Mar-09 23:01:45, 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. Any idea? Thanks. >> -james > > It would certainly work as intended! For instance: > > truncnorm<-function(HowMany,Min,Max,mean=0,sd=1){ > qnorm(runif(HowMany, pnorm(Min,mean,sd), pnorm(Max,mean, sd)), > mean, sd)} > Sample <- truncnorm(1000,-1,2.5) > hist(Sample,breaks=100) > > Hoping this helps, > Ted.truncnorm<-function(HowMany,Min,Max,mean=0,sd=1){ qnorm(runif(HowMany, pnorm(Min,mean,sd), pnorm(Max,mean, sd)), mean, sd)} Sample <- truncnorm(100000,-1,2.5) hist(Sample,breaks=100) u0<-(-0.975)+0.05*(0:69) yy<-dnorm(u0) du<-min(diff(H$breaks)) Y <- 100000*yy*du/(pnorm(2.5)-pnorm(-1.0)) lines(u0,Y) Ted. -------------------------------------------------------------------- E-Mail: (Ted Harding) <Ted.Harding at manchester.ac.uk> Fax-to-email: +44 (0)870 094 0861 Date: 11-Mar-09 Time: 00:05:23 ------------------------------ XFMail ------------------------------
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.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. 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. Duncan Murdoch
Because of a mistake I made in copying code into email, there has been confusion (expressed in some private emails). Hence the corrected version is appended at the end. On 11-Mar-09 00:05:26, Ted Harding wrote:> I have modified my example to make it more convincing! See at end. > > On 10-Mar-09 23:39:17, Ted Harding wrote: >> On 10-Mar-09 23:01:45, 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. Any idea? Thanks. >>> -james >> >> It would certainly work as intended! For instance: >> >> truncnorm<-function(HowMany,Min,Max,mean=0,sd=1){ >> qnorm(runif(HowMany, pnorm(Min,mean,sd), pnorm(Max,mean, sd)), >> mean, sd)} >> Sample <- truncnorm(1000,-1,2.5) >> hist(Sample,breaks=100) >> >> Hoping this helps, >> Ted. > > truncnorm<-function(HowMany,Min,Max,mean=0,sd=1){ > qnorm(runif(HowMany, pnorm(Min,mean,sd), pnorm(Max,mean, sd)), > mean, sd)} > Sample <- truncnorm(100000,-1,2.5) > hist(Sample,breaks=100) > u0<-(-0.975)+0.05*(0:69) > yy<-dnorm(u0) > du<-min(diff(H$breaks)) > Y <- 100000*yy*du/(pnorm(2.5)-pnorm(-1.0)) > lines(u0,Y) > > Ted.Corrected version: ------------------ truncnorm<-function(HowMany,Min,Max,mean=0,sd=1){ qnorm(runif(HowMany, pnorm(Min,mean,sd), pnorm(Max,mean, sd)), mean, sd)} Sample <- truncnorm(100000,-1,2.5) H <- hist(Sample,breaks=100) u0<-(-0.975)+0.05*(0:69) yy<-dnorm(u0) du<-min(diff(H$breaks)) Y <- 100000*yy*du/(pnorm(2.5)-pnorm(-1.0)) lines(u0,Y) Apologies for the error. Ted. -------------------------------------------------------------------- E-Mail: (Ted Harding) <Ted.Harding at manchester.ac.uk> Fax-to-email: +44 (0)870 094 0861 Date: 11-Mar-09 Time: 19:48:12 ------------------------------ XFMail ------------------------------
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