[R] truncated normal
Duncan Murdoch
murdoch at stats.uwo.ca
Wed Jul 23 22:34:06 CEST 2008
On 7/23/2008 4:22 PM, Duncan Murdoch wrote:
> On 7/23/2008 3:41 PM, cindy Guo wrote:
>> Hi, I want to generate random samples from truncated normal say
>> Normal(0,1)Indicator((0,1),(2,4)). It has more than one intervals. In the
>> library msm, it seems to me that the 'lower' and 'upper' arguments can only
>> be a number. I tried rtnorm(1,mean=0,sd=1, lower=c(0,2),upper=c(1,4)) and it
>> didn't work. Can you tell me how I can do truncated normal at more than one
>> intervals?
>
> The inverse CDF method will work. For example, this untested code:
>
> rtruncnorm <- function(n, mean=0, sd=1, lower=-Inf, upper=Inf) {
> plower <- pnorm(lower, mean, sd) # get the probability values for the
> # lower limits
> pupper <- pnorm(upper, mean, sd) # ditto for the upper limits
> p <- plower + (pupper - plower)*runif(n) # get random values between
> # those
> qnorm(p, mean, sd) # invert the CDFs
> }
>
> As I said, this is untested, but it should work if all of mean, sd,
> lower, and upper are the same length, and in some cases where they aren't.
>
> Duncan Murdoch
One case where the code above will *not* work is if your truncation is
too far out in the tails, e.g. a standard normal, truncated to be in the
interval [100, 101]. It is possible to do those in R, but it's tricky
to avoid rounding error: in the example above, both plower and pupper
will be exactly 1 in this case. You need to do calculations on a log
scale, and work with upper tails instead of lower ones. It all gets
kind of messy.
Duncan Murdoch
More information about the R-help
mailing list