Laboratorio No. 5
Ricardo Cruz 1048220
2023-03-18
library(urca)
library(vars)## Loading required package: MASS## Loading required package: strucchange## Loading required package: zoo##
## Attaching package: 'zoo'## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric## Loading required package: sandwich## Loading required package: lmtestlibrary(mFilter)## Warning: package 'mFilter' was built under R version 4.2.3library(tseries)## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoolibrary(TSstudio)## Warning: package 'TSstudio' was built under R version 4.2.3library(forecast)
library(fpp2)## ── Attaching packages ────────────────────────────────────────────── fpp2 2.5 ──## ✔ ggplot2 3.4.0 ✔ expsmooth 2.3
## ✔ fma 2.5## library(TSA)## Registered S3 methods overwritten by 'TSA':
## method from
## fitted.Arima forecast
## plot.Arima forecast##
## Attaching package: 'TSA'## The following objects are masked from 'package:stats':
##
## acf, arima## The following object is masked from 'package:utils':
##
## tarlibrary(tidyverse)## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.0.10 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ lubridate 1.9.2 ✔ tibble 3.1.8
## ✔ purrr 1.0.1 ✔ tidyr 1.3.0## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ stringr::boundary() masks strucchange::boundary()
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ✖ dplyr::select() masks MASS::select()
## ✖ readr::spec() masks TSA::spec()
## ℹ Use the ]8;;http://conflicted.r-lib.org/conflicted package]8;; to force all conflicts to become errorslibrary(pdfetch)
library(tidyverse)
library(yuima)## Loading required package: stats4
## Loading required package: expm
## Loading required package: Matrix
##
## Attaching package: 'Matrix'
##
## The following objects are masked from 'package:tidyr':
##
## expand, pack, unpack
##
##
## Attaching package: 'expm'
##
## The following object is masked from 'package:Matrix':
##
## expm
##
## Loading required package: cubature
## Loading required package: mvtnorm
## ########################################
## This is YUIMA Project package v.1.15.22
## Why don't you try yuimaGUI package?
## Visit: http://www.yuima-project.com
## ########################################
##
## Attaching package: 'yuima'
##
## The following object is masked from 'package:tibble':
##
## char
##
## The following object is masked from 'package:stats':
##
## simulateSegundo ejercicio Yahoo Finance
Obtenemos los datos desde el sitio de yahoo, buscamos el acronimo o ticker de una empresa y lo colocamos, en este caso escogí Walmart que tiene el ticker WMT con fecha desde el 01/01/2019 hasta el 1/03/2023
NVCRdata <- pdfetch_YAHOO("WMT",from = c("2019-01-01"),to = c("2023-03-1"), interval = '1d') Obtenemos la columna que es de nuestro interés que en este caso es la cuarta
Novocure <- NVCRdata[,4]
length(Novocure)## [1] 1047Convertimos en los datos en una serie de tiempo
tsNovoCure <- ts(Novocure, start = c(2019,1),frequency=365)Gráficamos la serie de tiempo
plot(tsNovoCure)
Se calculan las diferencias de la serie de datos con logaritmo
l_NovoCure<-diff(log(tsNovoCure))
plot(l_NovoCure)
Calculo parámetros iniciales manera 1
Delta <- 1/365
alpha <- mean(l_NovoCure)/Delta
sigma <- sqrt(var(l_NovoCure)/Delta)
mu <- alpha +0.5*sigma^2
x0<-tsNovoCure[1]Calculo parámetros iniciales manera 2
x <- tsNovoCure
gBm <- setModel(drift="mu*x", diffusion="sigma*x", xinit=x0)## Warning in yuima.warn("Solution variable (lhs) not specified. Trying to use state variables."):
## YUIMA: Solution variable (lhs) not specified. Trying to use state variables.mod <- setYuima(model=gBm, data=setData(tsNovoCure, delta=Delta))
set.seed(123)
fit <- qmle(mod, start=list(mu=0, sigma=1),
lower=list(mu=0.1, sigma=0.1),
upper=list(mu=100, sigma=10))
summary(fit)## Quasi-Maximum likelihood estimation
##
## Call:
## qmle(yuima = mod, start = list(mu = 0, sigma = 1), lower = list(mu = 0.1,
## sigma = 0.1), upper = list(mu = 100, sigma = 10))
##
## Coefficients:
## Estimate Std. Error
## sigma 0.2797374 0.006117664
## mu 0.1858239 0.165246062
##
## -2 log L: 4296.521#comparación
coef(fit)## sigma mu
## 0.2797374 0.1858239sigma## WMT.close
## WMT.close 0.2796809mu## WMT.close
## WMT.close 0.1858412# nsim
# t
gbm_vec <- function(nsim = 10000, t = 25, mu = 0, sigma = 0.1, S0 = 100, dt = 1./365) {
# matrix of random draws - one for each day for each simulation
epsilon <- matrix(rnorm(t*nsim), ncol = nsim, nrow = t)
# get GBM and convert to price paths
gbm <- exp((mu - sigma * sigma / 2) * dt + sigma * epsilon * sqrt(dt))
gbm <- apply(rbind(rep(S0, nsim), gbm), 2, cumprod)
return(gbm)
}
gBm##
## Diffusion process
## Number of equations: 1
## Number of Wiener noises: 1
## Parametric model with 2 parametersvalores_simulados <- simulate(gBm, true.parameter = list(mu=mu, sigma=sigma))## Warning in yuima.warn("'delta' (re)defined."):
## YUIMA: 'delta' (re)defined.plot(valores_simulados)
PROBAR CON FORMA 1 Y FORMA 2. Además fijamos los datos como
cantidad de simulaciones, Días a predecir, mu, sigma, xo y dt. Por último mandamos la información a las simulaciones
nsim <- 1000
t <- 1074
mu <- 0.1858412
sigma <- 0.2796809
S0 <- 93.339996
dt = 1/365
gbm <- gbm_vec(nsim, t, mu, sigma, S0, dt)
gbm_df <- as.data.frame(gbm) %>%
mutate(ix = 1:nrow(gbm)) %>%
pivot_longer(-ix, names_to = 'sim', values_to = 'price')
gbm_df %>%
ggplot(aes(x=ix, y=price, color=sim)) +
geom_line() +
theme(legend.position = 'none')
data.frame(price = gbm[1075,]) %>%
ggplot(aes(x = price)) +
geom_histogram(aes(y = ..density..), binwidth = 0.1) +
geom_density() +
ggtitle('terminal price distribution')## Warning: The dot-dot notation (`..density..`) was deprecated in ggplot2 3.4.0.
## ℹ Please use `after_stat(density)` instead.
Predecimos el precio en el que estará la acción en los siguientes 27 días
D <- gbm[1048, ] %>%
density()
D$x[which.max(D$y)]## [1] 112.8646