8.6 Random Number from Gamma Distribution

Rgamma
## function (n, alph, bet) 
## {
##     Res = vector(length = n)
##     for (i in 1:n) {
##         u = runif(2)
##         v = -log(u)
##         if (v[2] > (alph - 1) * (v[1] - log(v[1]) - 1)) 
##             Res[i] = bet * v[1]
##     }
##     return(Res)
## }
## <bytecode: 0x000001e2acc30308>
## <environment: namespace:math>

Test script

Rgamma(1, 1, 2)
## [1] 0.6656
Rgamma(2, 1, 2)
## [1] 5.471 0.683
x = Rgamma(10000, 1, 2)
summary(x)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##   0.001   0.583   1.393   2.008   2.775  19.604
sd(x)
## [1] 2.014
var(x)
## [1] 4.055