Semana 1

Juan Zapata - USMP

2020-05-14

Install.packages

pkges <- c("WDI","tidyverse","tseries","forecast")
#install.packages("pkges")
lapply(pkges,library,character.only=T)
#> -- Attaching packages -------------- tidyverse 1.3.0 --
#> v ggplot2 3.3.0     v purrr   0.3.3
#> v tibble  2.1.3     v dplyr   0.8.5
#> v tidyr   1.0.2     v stringr 1.4.0
#> v readr   1.3.1     v forcats 0.5.0
#> -- Conflicts ----------------- tidyverse_conflicts() --
#> x dplyr::filter() masks stats::filter()
#> x dplyr::lag()    masks stats::lag()
#> Registered S3 method overwritten by 'quantmod':
#>   method            from
#>   as.zoo.data.frame zoo
#> [[1]]
#> [1] "WDI"       "stats"     "graphics"  "grDevices" "utils"     "datasets" 
#> [7] "methods"   "base"     
#> 
#> [[2]]
#>  [1] "forcats"   "stringr"   "dplyr"     "purrr"     "readr"     "tidyr"    
#>  [7] "tibble"    "ggplot2"   "tidyverse" "WDI"       "stats"     "graphics" 
#> [13] "grDevices" "utils"     "datasets"  "methods"   "base"     
#> 
#> [[3]]
#>  [1] "tseries"   "forcats"   "stringr"   "dplyr"     "purrr"     "readr"    
#>  [7] "tidyr"     "tibble"    "ggplot2"   "tidyverse" "WDI"       "stats"    
#> [13] "graphics"  "grDevices" "utils"     "datasets"  "methods"   "base"     
#> 
#> [[4]]
#>  [1] "forecast"  "tseries"   "forcats"   "stringr"   "dplyr"     "purrr"    
#>  [7] "readr"     "tidyr"     "tibble"    "ggplot2"   "tidyverse" "WDI"      
#> [13] "stats"     "graphics"  "grDevices" "utils"     "datasets"  "methods"  
#> [19] "base"
## -- Attaching packages ---------------- tidyverse 1.3.0 --
## v ggplot2 3.3.0     v purrr   0.3.4
## v tibble  3.0.0     v dplyr   0.8.5
## v tidyr   1.0.2     v stringr 1.4.0
## v readr   1.3.1     v forcats 0.5.0
## -- Conflicts ------------- tidyverse_conflicts() --
## x dplyr::filter() masks stats::filter()
## x dplyr::lag()    masks stats::lag()
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
## [[1]]
## [1] "WDI"       "stats"     "graphics"  "grDevices" "utils"     "datasets" 
## [7] "methods"   "base"     
## 
## [[2]]
##  [1] "forcats"   "stringr"   "dplyr"     "purrr"     "readr"     "tidyr"    
##  [7] "tibble"    "ggplot2"   "tidyverse" "WDI"       "stats"     "graphics" 
## [13] "grDevices" "utils"     "datasets"  "methods"   "base"     
## 
## [[3]]
##  [1] "tseries"   "forcats"   "stringr"   "dplyr"     "purrr"     "readr"    
##  [7] "tidyr"     "tibble"    "ggplot2"   "tidyverse" "WDI"       "stats"    
## [13] "graphics"  "grDevices" "utils"     "datasets"  "methods"   "base"     
## 
## [[4]]
##  [1] "forecast"  "tseries"   "forcats"   "stringr"   "dplyr"     "purrr"    
##  [7] "readr"     "tidyr"     "tibble"    "ggplot2"   "tidyverse" "WDI"      
## [13] "stats"     "graphics"  "grDevices" "utils"     "datasets"  "methods"  
## [19] "base"

Importación de datos

options(scipen=999)
datos <- WDI(indicator='NY.GDP.MKTP.KN', country=c('CL','PE','BR',"AR"), start=1960, end=2018)
summary(datos)
#>     iso2c             country          NY.GDP.MKTP.KN                 year     
#>  Length:236         Length:236         Min.   :    69946000000   Min.   :1960  
#>  Class :character   Class :character   1st Qu.:   292830863724   1st Qu.:1974  
#>  Mode  :character   Mode  :character   Median :   664997245255   Median :1989  
#>                                        Mean   : 15935318938400   Mean   :1989  
#>                                        3rd Qu.:  7182050744270   3rd Qu.:2004  
#>                                        Max.   :153758254309000   Max.   :2018

1.Promedio de vector de datos (variable)

df <- datos %>%  rename(GDP=NY.GDP.MKTP.KN) %>% select(-country) %>% mutate(GDP=GDP/(10^9))
datos.rshp<-reshape(data = df,timevar = "iso2c",idvar = "year",v.names = "GDP",direction = "wide") %>%
  rename(AR=GDP.AR,
         BR=GDP.BR,
         CL=GDP.CL,
         PE=GDP.PE)

datos.rshp<- datos.rshp %>% arrange(year) %>% select(-year)
ts <- ts(datos.rshp$AR,start = c(1960,1),end = c(2018,1))
plot(ts,main="PBI real MN (1960-2018)",xlab="",ylab="Miles de millones")

l_ts<-diff(log(ts))
plot(l_ts,main="Variación % PBI",col=2,xlab="",ylab="%")

Estimación de Parámetros

mean(l_ts)
#> [1] 0.02331694
var(l_ts)
#> [1] 0.002741129
sd(l_ts)
#> [1] 0.05235579
adf.test(l_ts)
#> 
#>  Augmented Dickey-Fuller Test
#> 
#> data:  l_ts
#> Dickey-Fuller = -3.7345, Lag order = 3, p-value = 0.02969
#> alternative hypothesis: stationary
par(mfrow=c(3,1))
acf(ts)

acf(l_ts)

pacf(l_ts)

fit <- arima(x = log(ts), 
              c(1, 3, 0));fit
#> 
#> Call:
#> arima(x = log(ts), order = c(1, 3, 0))
#> 
#> Coefficients:
#>           ar1
#>       -0.5448
#> s.e.   0.1106
#> 
#> sigma^2 estimated as 0.00925:  log likelihood = 51.49,  aic = -98.98

pred <- predict(fit, n.ahead = 5)
par(mfrow=c(1,1))
ts.plot(ts,exp(pred$pred))

#simulacion de montecarlo

days <- 200

changes <- rnorm(200,mean=1.001,sd=0.005)

plot(cumprod(c(20,changes)),type='l',
     ylab="Price",
     xlab="day",
     main="BAYZ closing price (sample path)")

runs <- 100000

#simulates future movements and returns the #closing price on day 200

generate.path <- function(){
  days <- 200
  changes <- rnorm(200,mean=1.001,sd=0.008)
  sample.path <- cumprod(c(20,changes))
  closing.price <- sample.path[days+1] #+1 because we add the opening price
  return(closing.price)
}

old.closing <- cumprod(c(20,changes))
mc.closing <- replicate(runs,generate.path())
mean(old.closing)
#> [1] 22.70613
mean(mc.closing)
#> [1] 24.43217
quantile(mc.closing,0.95)
#>      95% 
#> 29.25689
quantile(mc.closing,0.05)
#>       5% 
#> 20.14765

#REMUESTREO

x<-c(58,67,74,74,80,89,95,97,98,107)
est<-median(x)
B<-100
n<-length(x)
res<-matrix(0,B,n)
for(h in 1:B){
  res[h,]<-sample(x,n,replace=T)
}
mediana.B<-apply(res,1,median)
s.B<-sd(mediana.B)
s.B #Error Estandar bootstrap
#> [1] 7.928868

Metodo percentiles

quantile(mediana.B,probs=c(0.025,0.975))
#>  2.5% 97.5% 
#>  70.5  97.0

#2.5% 97.5% 74 97

#JACKNIFE

x<-c(58,67,74,74,80,89,95,97,98,107)
est<-median(x)
n<-length(x)
median.j<-numeric()
for(h in 1:n){
  median.j[h]<-median(x[-h])
}
est.j<-mean(median.j)
ts.plot(est.j,exp(pred$pred))

sesgo.j<-(n-1)*(est.j-est)
ts.plot(sesgo.j,exp(pred$pred))

se.j<- sqrt((n-1)*sum((median.j-est.j)^2)/n)
ts.plot(se.j,exp(pred$pred))