#install.packages(“yuima”) #install.packages(“pdfetch”) library(pdfetch) library(tidyverse) library(yuima) #obtenemos los datos desde el sitio de yahoo NVCRdata <- pdfetch_YAHOO(“NVCR”,from = c(“2019-01-01”),to = c(“2020-01-01”), interval = ‘1d’) #Obtenemos la columna que es de nuestro interés Novocure <- NVCRdata[,4] length(Novocure) #Convertimos en una serie temporal tsNovoCure <- ts(Novocure, start = c(2019,1),frequency=365) #Graficamos plot(tsNovoCure) #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=“mux", diffusion="sigmax”, xinit=x0) 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)
#comparación coef(fit) sigma mu
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(tnsim), 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 valores_simulados <- simulate(gBm, true.parameter = list(mu=mu, sigma=sigma)) plot(valores_simulados)
#PROBAR CON FORMA 1 Y FORMA 2 nsim <- 1000 t <- 366 mu <- 1.517177 sigma <- 0.531907 S0 <- 32.72 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[253, ]) %>% ggplot(aes(x = price)) + geom_histogram(aes(y = ..density..), binwidth = 0.1) + geom_density() + ggtitle(‘terminal price distribution’)
D <- gbm[253, ] %>% density()
D\(x[which.max(D\)y)] #resultado del modelo
NVCRdata2 <- pdfetch_YAHOO(“NVCR”,from = c(“2019-01-01”),to = c(“2020-01-10”), interval = ‘1d’) Novocure2 <- NVCRdata2[,4] View(Novocure2)