library(pdfetch)
## Warning: package 'pdfetch' was built under R version 4.2.2
library(tidyverse)
## Warning: package 'tidyverse' was built under R version 4.2.2
## ── Attaching packages ─────────────────────────────────────── tidyverse 1.3.2 ──
## ✔ ggplot2 3.3.6      ✔ purrr   1.0.1 
## ✔ tibble  3.1.8      ✔ dplyr   1.0.10
## ✔ tidyr   1.2.1      ✔ stringr 1.4.1 
## ✔ readr   2.1.3      ✔ forcats 1.0.0
## Warning: package 'tidyr' was built under R version 4.2.2
## Warning: package 'readr' was built under R version 4.2.2
## Warning: package 'purrr' was built under R version 4.2.2
## Warning: package 'dplyr' was built under R version 4.2.2
## Warning: package 'stringr' was built under R version 4.2.2
## Warning: package 'forcats' was built under R version 4.2.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag()    masks stats::lag()
library(yuima)
## Warning: package 'yuima' was built under R version 4.2.2
## Loading required package: zoo
## Warning: package 'zoo' was built under R version 4.2.2
## 
## Attaching package: 'zoo'
## 
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## 
## Loading required package: stats4
## Loading required package: expm
## Warning: package 'expm' was built under R version 4.2.2
## Loading required package: Matrix
## Warning: package 'Matrix' was built under R version 4.2.2
## 
## 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
## Warning: package 'cubature' was built under R version 4.2.2
## 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':
## 
##     simulate
#obtenemos los datos desde el sitio de yahoo
NVCRdata <- pdfetch_YAHOO("AMZN",from = c("2018-05-03"),to = c("2023-05-03"), interval = '1d')
#Obtenemos la columna que es de nuestro interés
Novocure <- NVCRdata[,4]
length(Novocure)
## [1] 1218
#Convertimos en una serie temporal
tsNovoCure <- ts(Novocure, start = c(2018,1),frequency=365)
#Graficamos
plot(tsNovoCure)

#En la gráfica podemos observar que del año 2019 al año 2020 tuvo un comportamiento a la alza, pero empezo a disminuir e ir a la baja durante 2021. 

#Ver si estacionaria
#.-....
#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.436174 0.008840632
## mu    0.147963 0.238869619
## 
## -2 log L: 5857.482
#comparación
coef(fit)
##    sigma       mu 
## 0.436174 0.147963
sigma
##            AMZN.close
## AMZN.close  0.4364539
mu
##            AMZN.close
## AMZN.close  0.1480943
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 parameters
valores_simulados <- simulate(gBm, true.parameter = list(mu=mu, sigma=sigma))
## Warning in yuima.warn("'delta' (re)defined."): 
## YUIMA: 'delta' (re)defined.
plot(valores_simulados)

#Amazon a mostrado varaiaciones durante el tiempo, estas variaciones van hacia la alza. 

#PROBAR CON FORMA 1 Y FORMA 2
nsim <- 1219
t <- 1218
mu <- 0.147963  
sigma <- 0.436174
S0 <- x0 
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[1218, ]) %>%
  ggplot(aes(x = price)) +
  geom_histogram(aes(y = ..density..), binwidth = 0.1) +
  geom_density() + 
  ggtitle('terminal price distribution')

D <- gbm[1218, ] %>%
  density()

D$x[which.max(D$y)]  #resultado del modelo
## [1] 54.54574
summary(cars)
##      speed           dist       
##  Min.   : 4.0   Min.   :  2.00  
##  1st Qu.:12.0   1st Qu.: 26.00  
##  Median :15.0   Median : 36.00  
##  Mean   :15.4   Mean   : 42.98  
##  3rd Qu.:19.0   3rd Qu.: 56.00  
##  Max.   :25.0   Max.   :120.00

Note that the echo = FALSE parameter was added to the code chunk to prevent printing of the R code that generated the plot.