mean.tnorm<-function(mu,sd,lower,upper){
##return the expectation of a truncated normal distribution
lower.std=(lower-mu)/sd
upper.std=(upper-mu)/sd
mean=mu+sd*(dnorm(lower.std)-dnorm(upper.std))/
(pnorm(upper.std)-pnorm(lower.std))
return(mean)
}
var.tnorm<-function(mu,sd,lower,upper){
##return the variance of a truncated normal distribution
lower.std=(lower-mu)/sd
upper.std=(upper-mu)/sd
variance=sd^2*(1+(lower.std*dnorm(lower.std)-upper.std*dnorm(upper.std))/
(pnorm(upper.std)-pnorm(lower.std))-((dnorm(lower.std)-dnorm(upper.std))/
(pnorm(upper.std)-pnorm(lower.std)))^2)
return(variance)
}
###Testing
> library(msm)
> a=rtnorm(1000000,-5,2,1,3)
> paste(mean(a),var(a))
[1] "1.52135857341077 0.197281057170982"
> paste(mean.tnorm(-5,2,1,3),var.tnorm(-5,2,1,3))
[1] "1.52090857118 0.197111175109889"
Showing posts with label random numbers. Show all posts
Showing posts with label random numbers. Show all posts
Thursday, September 3, 2009
Truncated Normal Distribution
Many distributions may be used to describe patterns that are non-negative; however, there are not as many choices when an upper bound is also needed (although the beta distribution is very flexible). For various reasons, truncated distributions are sometimes preferred, and the truncated normal is particularly popular. While R has a package that includes the standard functions for this distribution (see rtnorm, dtnorm, etc. in the msm pacakge), the true expectation and variance of the distribution may be of interest. It turns out that the first two moments of the truncated normal are not too hard to calculate (but worth writing functions for):
Thursday, May 17, 2007
Random number generation
Luc Devroye from McGill University's School of Computer Science has provided pdf copies of his book, "Non-Uniform Random Variate Generation," on-line.
http://cg.scs.carleton.ca/~luc/rnbookindex.html
Most common distributions are already built into R, but sometimes you need to build your own -- Luc provides the algorithms to do it.
http://cg.scs.carleton.ca/~luc/rnbookindex.html
Most common distributions are already built into R, but sometimes you need to build your own -- Luc provides the algorithms to do it.
Subscribe to:
Posts (Atom)