In [1]:
# rm(list=ls()) 
options(OutDec = ",")
#==============================================================================
# Descricao
#
# Exemplo da normal bivariada. Analise da estrutura de dependencia.
#
#==============================================================================
silence <- suppressPackageStartupMessages # Para omitir mensagens de alertas
silence(library(mvtnorm))                 # uma biblioteca da normal multivariada
In [2]:
#==============================================================================
# Criando estrutura de grade bidimensional de pontos para fazer graficos
#==============================================================================
nx  <- 100                              # numero de pontos no eixo x
ny  <- 100                              # numero de pontos no eixo y
x   <- seq(-4,4,length=nx)              # sequencia de ponto em x
y   <- seq(-4,4,length=ny)              # sequencia de ponto em x
z   <- matrix(NA,nx,ny)                 # z = f(x,y)
mu  <- c(0,0)                           # vetor de medias
S   <- matrix(c(1, .7, .7, 1),2,2)      # matriz de covariancias 
In [3]:
#==============================================================================
# Avaliando f(x,y) para cada x e cada y na grade de valores
#==============================================================================
for(i in 1:nx)
    for(j in 1:ny)
        z[i,j] <- dmvnorm(c(x[i],y[j]), mean=mu,sigma=S)
In [4]:
#==============================================================================
# Grafico da funcao de densidade de probabilidade (fdp) da normal bivariada
#==============================================================================
persp(x,y,z, phi = 45, theta = 30,col="darkorange")
No description has been provided for this image
In [5]:
#==============================================================================
# Grafico de contorno da fdp da normal bivariada
#==============================================================================
contour(x,y,z,col="red",lwd=2,bty="n")
No description has been provided for this image
In [6]:
#==============================================================================
# Gerando dados da normal bivariada
#==============================================================================
set.seed(5485)                    # fixando a semente pseudo-aleatorio

w <- rmvnorm(100000,mean=mu,sigma=S)
In [7]:
#==============================================================================
# Grafico de dispersao
#==============================================================================
plot(w,pch=".",main="",xlab="x",ylab="y",bty="n")
contour(x,y,z,col="red",lwd=3,add=T,drawlabels=F)
No description has been provided for this image
In [8]:
#==============================================================================
# Histograma das marginais
#==============================================================================
hist(w[,1],prob=T,nclass=50,main="",xlab="x",ylab="Densidade")
lines(x,dnorm(x),lwd=2,col="red")
hist(w[,2],prob=T,nclass=50,main="",xlab="y",ylab="Densidade")
lines(y,dnorm(y),lwd=2,col="red")
No description has been provided for this image
No description has been provided for this image
In [9]:
print("==== Médias ====")
print(apply(w,2,"mean"))  # vetor (linhas) de medias amostral

print("==== Covariâncias ====")
print(var(w))             # matriz de covariancias amostral

print("==== Correlações ====")
print(cor(w))             # matriz de correlacoes amostral
[1] "==== Médias ===="
[1] -0,002417111  0,002322957
[1] "==== Covariâncias ===="
          [,1]      [,2]
[1,] 0,9960662 0,6987780
[2,] 0,6987780 0,9988886
[1] "==== Correlações ===="
          [,1]      [,2]
[1,] 1,0000000 0,7005459
[2,] 0,7005459 1,0000000