#source("c:/Work07/ShortCourse07/Programs/BetaAR-2.txt",print.eval=TRUE)#
#This gives an accept-reject algorithm for a beta rv#
a<-2.7; b<-6.3; c<-1.67;nsim<-25000;
#Generate u and v#
u<-runif(nsim);
v<-rbeta(nsim,2,6)				#beta candidate
#---------Generate Y, the beta random variable--------------#
test<-dbeta(v, a, b)/(c*dbeta(v, 2, 6));				#density ratio
Y<-v*(u<test)					#accepted values
Y<-Y[Y!=0]					#eliminate zeros
length(Y)/nsim					#percent accepted
#----------Plot---------------------------------------------#
par(mfrow=c(1,2))
hist(v,col="grey",freq=F,xlim=c(0,1),ylim=c(0,3))
par(new=T)
plot(function(x)(dbeta(x, 2, 6)),xlim=c(0,1),ylim=c(0,3),xlab="",ylab="")
hist(Y,col="grey",freq=F,xlim=c(0,1),ylim=c(0,3))
par(new=T)
plot(function(x)(dbeta(x, a, b)),xlim=c(0,1),ylim=c(0,3),xlab="",ylab="");
#------------------------------------------------------------