In [1]:
options(OutDec = ",")
#==============================================================================
# Descricao
#==============================================================================
# Os dados a seguir correspondem aos indices diarios da BOVESPA de 3 de
# janeiro de 2011 a 10 de setembro de 2021.
#
# Os dados foram obtidos do sitio http://br.financas.yahoo.com/
#
# Concentraremos nos indices (precos) de fechamento ("Close")
#
#==============================================================================
# Carregando pacotes
#==============================================================================
silence <- suppressPackageStartupMessages # Para omitir mensagens de alertas
silence(library(tseries)) # Teste de Jarque-Bera: 'jarque.bera.test'
silence(library(psych)) # Estatistica descritiva: 'describe'
silence(library(latticeExtra)) # Para construir dois graficos em um so'
#==============================================================================
# Lendo os dados do arquivo BVSP.csv
#==============================================================================
BVSP <- read.csv("../dados/BVSP.csv",header=TRUE,sep=",")
print(head(BVSP))
print(tail(BVSP))
Date Open High Low Close Adj.Close Volume
1 2011-01-03 69310 70471 69305 69962 69962 1862400
2 2011-01-04 69962 70318 69560 70318 70318 2427200
3 2011-01-05 70311 71173 69802 71091 71091 2309200
4 2011-01-06 71093 71167 70469 70579 70579 2546000
5 2011-01-07 70580 70783 69718 70057 70057 1761000
6 2011-01-10 70056 70133 69666 70127 70127 1610800
Date Open High Low Close Adj.Close Volume
2609 2021-09-02 119394 119397 116534 116677 116677 9862800
2610 2021-09-03 116679 117396 115583 116933 116933 12517500
2611 2021-09-06 116926 117981 116156 117869 117869 5111300
2612 2021-09-08 117866 117866 113172 113413 113413 12254400
2613 2021-09-09 113413 116354 112435 115361 115361 13890600
2614 2021-09-10 115370 116896 114286 114286 114286 11155700
In [2]:
#==============================================================================
# Copiando os indices de fechamento como a serie de precos
#==============================================================================
preco <- BVSP$Close
dias <- as.Date(BVSP$Date) # Dias
n <- length(preco) # Numero de observacoes
In [3]:
#==============================================================================
# Grafico da serie de precos
#==============================================================================
par(mfrow=c(1,1),lwd=2.0,cex.lab=1.5,cex.axis=1.5,lab=c(15,5,5),
mar=c(5,5,2,2.5),xpd=T,cex.main=2.0,bty="n")
plot(dias,preco/1000,lwd=2,type="l",col="black",main="",xlab="Dias",
ylab="Preço (milhares)",ylim=c(20,140),
xlim=c(as.Date("2011-01-01"),as.Date("2021-12-31")))
In [4]:
#==============================================================================
# Calculando o retorno de capitalizacao continua, ou simplesmente log-retorno
#==============================================================================
rt <- diff(log(preco)) # Calcula a primeira diferenca
# Grafico da serie de retornos brutos: r_t = log(1+R_t) = log(P_t)-log(P_{t-1})
par(mfrow=c(1,1),lwd=2.0,cex.lab=1.5,cex.axis=1.4,lab=c(15,5,5),
mar=c(5,5,2,2.5),xpd=T,cex.main=2.0,bty="n")
plot(dias[2:n],rt,lwd=2,col="black",main="",xlab="Dias",ylab="Log Retorno",
type="l",xlim=c(as.Date("2011-01-01"),as.Date("2021-12-31")),
ylim=c(-0.16,0.16))
In [5]:
#==============================================================================
# A serie de precos e log retornos em um unico grafico
#==============================================================================
# Primeira opcao: depende de 'library(latticeExtra)'
#==============================================================================
obj1 <- xyplot(rt ~ dias,type="l",lwd=2,col="black",xlab="Dias",
ylab="Log Retorno")
obj2 <- xyplot(preco ~ dias,type="l",lwd=2,col="red",ylab="Preço")
doubleYScale(obj1,obj2,add.ylab2=TRUE,use.style=FALSE)
In [6]:
#==============================================================================
# Segunda opcao de graficos para analise
#==============================================================================
par(mfrow=c(2,1),lwd=2.0,cex.lab=1.5,cex.axis=1.5,lab=c(10,5,5),
mar=c(5,5,2,2.5),xpd=T,cex.main=2.0,bty="n")
plot(dias,preco/1000,lwd=2,col="black",main="",xlab="Dias",ylab="Preço",
type="l",xlim=c(as.Date("2011-01-01"),as.Date("2021-12-31")),
ylim=c(20,140))
plot(dias[2:n],rt,lwd=2,col="black",main="",xlab="Dias",ylab="Log Retorno",
type="l",xlim=c(as.Date("2011-01-01"),as.Date("2021-12-31")),
ylim=c(-0.16,0.16))
In [7]:
#==============================================================================
# Calculando varios momentos amostrais
#==============================================================================
options(digits=4)
r.describe <- describe(rt)
r.media <- r.describe[,3]
r.desvpad <- r.describe[,4]
r.mediana <- r.describe[,5]
r.assimetria <- r.describe[,11]
r.exc.curtose <- r.describe[,12]
STATS <- cbind(r.media,r.desvpad,r.mediana,r.assimetria,
r.exc.curtose)
print(STATS)
r.media r.desvpad r.mediana r.assimetria r.exc.curtose [1,] 0,0001878 0,01616 0,0003192 -0,8946 12,23
In [8]:
#==============================================================================
# Histograma dos log retornos
#==============================================================================
par(mfrow=c(1,1),lwd=2.0,cex.lab=1.5,cex.axis=1.5,lab=c(12,5,5),
mar=c(5,5,2,2.5),xpd=T,cex.main=2.0,bty="n")
hist(rt,lwd=2,col=0,main="",xlab="Log Retornos",ylab="Densidade",prob=T,
nclass=50,xlim=c(-0.17,0.17))
xseq <- seq(r.media-4*r.desvpad,r.media+4*r.desvpad,length=1000)
yseq <- dnorm(xseq,r.media,r.desvpad)
lines(xseq,yseq,lwd=3,col="red")
text(-0.09,25,"Certamente não normal",cex=1.2)
In [9]:
#==============================================================================
# Boxplot dos log retornos
#==============================================================================
par(mfrow=c(1,1),lwd=2.0,cex.lab=1.5,cex.axis=1.5,lab=c(10,5,5),
mar=c(5,5,2,2.5),xpd=T,cex.main=2.0,bty="n")
boxplot(rt,lwd=2,col="red",main="",ylab="Log Retornos",pch=16,
ylim=c(-0.2,0.2))
In [10]:
#==============================================================================
# Normal QQ-plot dos log retornos
#==============================================================================
par(mfrow=c(1,1),lwd=2.0,cex.lab=1.5,cex.axis=1.5,lab=c(10,5,5),
mar=c(5,5,2,2.5),xpd=T,cex.main=2.0,bty="n")
qqnorm(rt,lwd=2,col="red",main="",xlab="Quantis teóricos",
ylab="Quantis amostrais",pch=16,xlim=c(-4,4),ylim=c(-0.2,0.2))
qqline(rt,lwd=2,col="black",xlim=c(-4,4),xpd=F)
In [11]:
#==============================================================================
# Testando a hipotese de distribuicao normal para os log retornos
#==============================================================================
# Lembre-se de verificar de novo a assimetria e a curtose
print(STATS[,4:5])
print("====================")
# Nivel de significancia dos testes
alpha <- 0.05
m <- length(rt) # Tamanho da serie de log retornos
# Testando a assimetria zero
tSr <- STATS[1,4]/sqrt(6/m) # Razao t-Student
test1 <- abs(tSr) > qnorm(1-alpha/2) # (TRUE=Rejeitar, FALSE=Nao rejeitar)
valorpSr <- 2*(1-pnorm(abs(tSr))) # Valor p
print(test1)
print(paste("P-valor = ",valorpSr,sep=""))
print("====================")
# Testando a excesso de curtose zero
tKr <- STATS[1,5]/sqrt(24/m) # Razao t-Student
test2 <- abs(tKr) > qnorm(1-alpha/2) # (TRUE=Rejeitar, FALSE=Nao rejeitar)
valorpKr <- 2*(1-pnorm(abs(tKr))) # Valor p
print(test2)
print(paste("P-valor = ",valorpKr,sep=""))
print("====================")
# Testando a normalidade - teste de Jarque Bera
JB <- tSr^2 + tKr^2
names(JB) <- "Jarque-Bera"
test3 <- JB > qchisq(1-alpha,2) # (TRUE=Rejeitar, FALSE=Nao rejeitar)
print(JB)
print(test3)
print("====================")
print(jarque.bera.test(rt)) # Verificar o valor p (mesma conclusao)
print("====================")
r.assimetria r.exc.curtose
-0,8946 12,2347
[1] "===================="
r.assimetria
TRUE
[1] "P-valor = 0"
[1] "===================="
r.exc.curtose
TRUE
[1] "P-valor = 0"
[1] "===================="
Jarque-Bera
16646
Jarque-Bera
TRUE
[1] "===================="
Jarque Bera Test
data: rt
X-squared = 16677, df = 2, p-value <2e-16
[1] "===================="