R Programming/Random Number Generation
From Wikibooks, the open-content textbooks collection
Contents |
[edit] Intro
Computer are deterministic and therefore are not able to generate purely random numbers. We generally use pseudo random number generator. The default algorithm in R is Mersenne-Twister but a long list of methods is available. Look at the help of RNGkind() to learn about random number generator.
?RNGkind
[edit] Seed
A pseudo random number generator is an algorithm based on a starting point called "seed". If you want to perform an exact replication of your program, you have to specify the seed using the function set.seed(). The argument of set.seed has to be an integer.
> set.seed(1) > runif(1) [1] 0.2655087 > set.seed(1) > runif(1) [1] 0.2655087
[edit] Sampling in a vector
Toss 10 coins
> sample(0:1,10,replace=T) [1] 1 0 0 0 1 0 0 1 1 1
Roll 10 dice
> sample(1:6,10,replace=T) [1] 4 1 5 3 2 5 5 6 3 2
play lottery (6 random numbers out of 49 without replacement)
> sample(1:49,6,replace=F) [1] 18 35 29 1 33 11
[edit] Sampling in a standard univariate distribution
You can use rnorm, rt, etc.
[edit] Misspecified argument
Note that if you put as argument of rnorm a vector instead of a number, R takes by default the length of the vector instead of returning an error. Here is an example :
x <- rnorm(10) # Sample a normal random vector set.seed(1) # use the seed z <- rnorm(x) # put a vector instead of a number as an argument of rnorm set.seed(1) # initalize the seed again z2 <- rnorm(length(x)) # sample in a vector with the same length as x plot(z2,z) # check that z and z2 are the same
[edit] Inverse CDF method
- If you know the inverse CDF (quantile function), you can generate the random variable by sampling in the standard uniform distribution and transforming using the CDF.
[edit] Importance sampling
- See Jeff Gill's routine
[edit] Metropolis algorithm
[edit] Gibbs algorithm
[edit] Quasi random numbers
Sometimes you need to generate quasi random sequences. The randtoolbox library provides several quasi random number generators.
> library(randtoolbox)
> halton(10, dim = 2, init = TRUE, normal = FALSE, usetime = FALSE)
[,1] [,2]
[1,] 0.5000 0.33333333
[2,] 0.2500 0.66666667
[3,] 0.7500 0.11111111
[4,] 0.1250 0.44444444
[5,] 0.6250 0.77777778
[6,] 0.3750 0.22222222
[7,] 0.8750 0.55555556
[8,] 0.0625 0.88888889
[9,] 0.5625 0.03703704
[10,] 0.3125 0.37037037
You can compare Halton draws with the standard R (pseudo) random number generator. Halton draws are much more systematic
>random <- cbind(runif(1000),runif(1000)) >halton <- halton(1000, dim = 2, init = TRUE, normal = FALSE, usetime = FALSE) >par(mfrow=c(2,2)) >plot(halton[,1],halton[,2]) >plot(random[,1],random[,2])
[edit] Resources
- Revolution Computing entry on pseudo random numbers
- Christian Robert's thought on random number generators
- Statistical Distributions Module Wessa.net it is online application which generates random numbers using R. You can have access to the R code and use it in your own programs.
