Adaptive Rejection Sampling Algorithm

Description

rARS generates a sequence of random numbers using the adaptive rejection sampling algorithm.

Usage

rARS(n, formula, min = -Inf, max = Inf, sp)

Arguments

Author(s)

Dong Zhang <dzhang0716@126.com>

Examples

# Example 1: Standard normal distribution
x1 <- rARS(100,"exp(-x^2/2)",-Inf,Inf,c(-2,2))

# Example 2: Truncated normal distribution
x2 <- rARS(100,"exp(-x^2/2)",-2.1,2.1,c(-2,2))

# Example 3: Normal distribution with mean=2 and sd=2
x3 <- rARS(100,"exp(-(x-2)^2/(2*4))",-Inf,Inf,c(-3,3))

# Example 4: Exponential distribution with rate=3
x4 <- rARS(100,"exp(-3*x)",0,Inf,c(2,3,100))

# Example 5: Beta distribution with alpha=3 and beta=4
x5 <- rARS(100,"x^2*(1-x)^3",0,1,c(0.4,0.6))

# Example 6: Gamma distribution with alpha=5 and lambda=2
x6 <- rARS(100,"x^(5-1)*exp(-2*x)",0,Inf,c(1,10))

# Example 7: Student distribution with df=10
x7 <- rARS(100,"(1+x^2/10)^(-(10+1)/2)",-Inf,Inf,c(-10,2))

# Example 8: F distribution with m=10 and n=5
x8 <- rARS(100,"x^(10/2-1)/(1+10/5*x)^(15/2)",0,Inf,c(3,10))

# Example 9:Cauchy distribution
x9 <- rARS(100,"1/(1+(x-1)^2)",-Inf,Inf,c(-2,2,10))

# Example 10:Rayleigh distribution with lambda=1
x10 <- rARS(100,"2*x*exp(-x^2)",0,Inf,c(0.01,10))

Adaptive Slice Sampling Algorithm With Stepping-Out Procedures

Description

rASS generates a sequence of random numbers by the adaptive slice sampling algorithm with stepping-out procedures.

Usage

rASS(n, x0 = 0, formula, w = 3)

Arguments

Author(s)

Dong Zhang <dzhang0716@126.com>

References

Neal R M. Slice sampling - Rejoinder[J]. Annals of Statistics, 2003, 31(3):758-767.

Examples

# Example 1: Sampling from exponential distribution with bounded domain
x<-rASS(100,-1,"1.114283*exp(-(4-x^2)^2)",3)
plot(density(x))

Concave-Convex Adaptive Rejection Sampling Algorithm

Description

rCCARS generates a sequence of random numbers by the concave-convex adaptive rejection sampling algorithm from target distributions with bounded domain.

Usage

rCCARS(n, cvformula, ccformula, min, max, sp)

Arguments

Details

Strictly speaking, the concave-convex adaptive rejection sampling algorithm can generate samples from target distributions who have bounded domains. For distributions with unbounded domain, rCCARS can also be used for sampling approximately. For example, if we want draw a sequence from N(0,1) by the concave-convex adaptive rejection sampling algorithm. We know that X~N(0,1) has a so small probability in two tails that we can ingore the parts at both ends. Pr(X>20)=P(X<-20)=2.753624e-89, therefore we can get random numbers approximately from N(0,1) with the bound [-20,20]. Also, you can make this bound large enough to reduce sampling error.

Author(s)

Dong Zhang <dzhang0716@126.com>

References

Teh Y W. Concave-Convex Adaptive Rejection Sampling[J]. Journal of Computational & Graphical Statistics, 2011, 20(3):670-691.

Examples

# Example 1: Generalized inverse bounded gaussian distribution with lambda=-1 and a=b=2
x<-rCCARS(100,"x+x^-1","2*log(x)",0.001,100,1)
hist(x,breaks=20,probability =TRUE);lines(density(x,bw=0.1),col="red",lwd=2,lty=2)
f <- function(x) {x^(-2)*exp(-x-x^(-1))/0.2797318}
lines(seq(0,5,0.01),f(seq(0,5,0.01)),lwd=2,lty=3,col="blue")

#The following examples are also available;
#But it may take a few minutes to run them.

# Example 2: Expontional bounded distribution
# x<-rCCARS(1000,"x^4","-8*x^2+16",-3,4,c(-2,1))
# hist(x,breaks=30,probability=TRUE);lines(density(x,bw=0.05),col="blue",lwd=2,lty=2)
# f <- function(x) exp(-(x^2-4)^2)/ 0.8974381
# lines(seq(-3,4,0.01),f(seq(-3,4,0.01)),col="red",lty=3,lwd=2)

# Example 3: Makeham bounded distribution
# x<-rCCARS(1000,"x+1/log(2)*(2^x-1)","-log(1+2^x)",0,5,c(1,2,3))
# hist(x,breaks=30,probability=TRUE);lines(density(x,bw=0.05),col="blue",lwd=2,lty=2)
# f <- function(x){(1+2^x)*exp(-x-1/log(2)*(2^x-1))}
# lines(seq(0,5,0.01),f(seq(0,5,0.01)),col="red",lty=3,lwd=2,type="l")

Modified Adaptive Rejection Sampling Algorithm

Description

rMARS generates a sequence of random numbers using the modified adaptive rejection sampling algorithm.

Usage

rMARS(n, formula, min = -Inf, max = Inf, sp, infp, m = 10^(-4))

Arguments

Author(s)

Dong Zhang <dzhang0716@126.com>

References

Martino L, Miguez J. A generalization of the adaptive rejection sampling algorithm[J]. Statistics & Computing, 2011, 21(4):633-647.

Examples

# Example 1: Exponential distribution
x <- rMARS(100,"exp(-(4-x^2)^2)",-Inf,Inf, c(-2.5,0,2.5),c(-2/sqrt(3),2/sqrt(3)))
hist(x,probability=TRUE,xlim=c(-3,3),ylim=c(0,1.2),breaks=20)
lines(density(x,bw=0.05),col="blue")
f <- function(x)(exp(-(4-x^2)^2))
lines(seq(-3,3,0.01),f(seq(-3,3,0.01))/integrate(f,-3,3)[[1]],lwd=2,lty=2,col="red")

#The following examples are also available;
#But it may take a few minutes to run them.

# Example 2: Distribution with bounded domain
# x <- rMARS(1000,"exp(-(x^2-x^3))",-3,2,c(-1,1),1/3)
# hist(x,probability=TRUE,xlim=c(-3,2.5),ylim=c(0,1.2),breaks=20)
# lines(density(x,bw=0.2),col="blue")
# f <- function(x) exp(-(x^2-x^3))
# lines(seq(-3,2,0.01),f(seq(-3,2,0.01))/integrate(f,-3,2)[[1]],lwd=2,lty=2,col="red",type="l")


# Example 3: Weibull distribution with k=3 and lambda=1
# x <- rMARS(100,"3*x^2*exp(-x^3)",10^-15,Inf,c(0.01,1),(1/3)^(1/3),m=10^-4)
# hist(x,probability=TRUE,breaks=20,xlim=c(0,2))
# lines(density(x,bw=0.15),col="blue")
# f <- function(x) 3*x^2*exp(-x^3)
# lines(seq(0,2,0.01),f(seq(0,2,0.01)),lwd=2,lty=2,col="red",type="l")


# Example 4: Mixed normal distribution with p=0.3,m1=2,m2=8,sigma1=1,sigma2=2
# x <- rMARS(100,"0.3/sqrt(2*pi)*exp(-(x-2)^2/2)+(1-0.3)/sqrt(2*pi)/2*exp(-(x-8)^2/8)",-Inf,Inf,
# c(-6,-4,0,3,6,15),c(-5.120801,-3.357761,3.357761,5.120801),m=10^-8)
# hist(x,breaks=20,probability=TRUE);lines(density(x,bw=0.45),col="blue",lwd=2)
# f <- function(x)0.3/sqrt(2*pi)*exp(-(x-2)^2/2)+(1-0.3)/sqrt(2*pi)/2*exp(-(x-8)^2/8)
# lines(seq(0,14,0.01),f(seq(0,14,0.01)),lty=3,col="red",lwd=2 )