In [1]:
# rm(list=ls())
options(OutDec = ",") 
#==============================================================================
# Exemplo para gerar de uma normal truncada entre 0 e infinito pelo metodo da
# rejeicao utilizando proposta Gama(1,1)=Exp(1)
#==============================================================================
set.seed(54345) # Semente
n  <- 100000    # Tamanho da amostra
In [2]:
#==============================================================================
# Funcao auxiliar de padronizacao: TCL sobre U(0,1)
#==============================================================================
z.pad <- function(u){
	return(sqrt(12*length(u))*(mean(u)-1/2))
}
In [3]:
#==============================================================================
# O TCL aqui se aplica sobre o tamanho da amostra m
#==============================================================================
In [4]:
#==============================================================================
# Para m = 2
#==============================================================================
m    <- 2  # Tamanho da amostra das U(0,1) 
U    <- matrix(runif(n*m),n,m)
x1   <- apply(U,1,"z.pad")
In [5]:
#==============================================================================
# Histograma
#==============================================================================
par(mfrow=c(1,1),lwd=2,cex.lab=1.5,cex.axis=1.5,lab=c(10,5,0),
    mar=c(4.5,5,2,1),bty="n")
xseq <- seq(-5,5,length=1000)
yseq <- dnorm(xseq)
hist(x1,nclass=30,prob=T,main="",ylim=c(0,0.40),xlim=c(-5,5),col="darkgreen",
     ylab=expression(f(x[1])),xlab=expression(x[1]))
lines(xseq,yseq,col="red",lwd=3)
No description has been provided for this image
In [6]:
#==============================================================================
# Para m = 5
#==============================================================================
m    <- 5  # Tamanho da amostra das U(0,1) 
U    <- matrix(runif(n*m),n,m)
x2   <- apply(U,1,"z.pad")
In [7]:
#==============================================================================
# Histograma
#==============================================================================
par(mfrow=c(1,1),lwd=2,cex.lab=1.5,cex.axis=1.5,lab=c(10,5,0),
    mar=c(4.5,5,2,1),bty="n")
hist(x2,nclass=30,prob=T,main="",ylim=c(0,0.40),xlim=c(-5,5),col="darkgreen",
     ylab=expression(f(x[2])),xlab=expression(x[2]))
lines(xseq,yseq,col="red",lwd=3)
No description has been provided for this image
In [8]:
#==============================================================================
# Para m = 12
#==============================================================================
m    <- 12  # Tamanho da amostra das U(0,1) 
U    <- matrix(runif(n*m),n,m)
x3   <- apply(U,1,"z.pad")
In [9]:
#==============================================================================
# Histograma
#==============================================================================
par(mfrow=c(1,1),lwd=2,cex.lab=1.5,cex.axis=1.5,lab=c(10,5,0),
    mar=c(4.5,5,2,1),bty="n")
hist(x3,nclass=30,prob=T,main="",ylim=c(0,0.40),xlim=c(-5,5),col="darkgreen",
     ylab=expression(f(x[3])),xlab=expression(x[3]))
lines(xseq,yseq,col="red",lwd=3)
No description has been provided for this image
In [10]:
#==============================================================================
# Para m = 50
#==============================================================================
m    <- 50  # Tamanho da amostra das U(0,1) 
U    <- matrix(runif(n*m),n,m)
x4   <- apply(U,1,"z.pad")
In [11]:
#==============================================================================
# Histograma
#==============================================================================
par(mfrow=c(1,1),lwd=2,cex.lab=1.5,cex.axis=1.5,lab=c(10,5,0),
    mar=c(4.5,5,2,1),bty="n")
hist(x4,nclass=50,prob=T,main="",ylim=c(0,0.40),xlim=c(-5,5),col="darkgreen",
     ylab=expression(f(x[4])),xlab=expression(x[4]))
lines(xseq,yseq,col="red",lwd=3)
No description has been provided for this image
In [12]:
# rm(list=ls()) 
#==============================================================================
# graphics.off()
#==============================================================================
# Fim
#==============================================================================