U2A7

setwd("~/EA1011")
library(markovchain)

## Package:  markovchain
## Version:  0.8.5-2
## Date:     2020-09-07
## BugReport: https://github.com/spedygiorgio/markovchain/issues

Para mayor información consultar el siguiente enlace: https://cran.r-project.org/web/packages/markovchain/vignettes/an_introduction_to_markovchain_package.pdf

-> Asignación U2A7, tomar una serie de datos y aplicar la cadena de markov con un planteamiento similar al de este ejercicio para hacer predicciones

Se tratara de determinar la probabilidad de que el club de beisbol de la Liga Mexicana del Pacifico, Yaquis de Cd. Obregon gane sus encuentros en base a los juegos ganados por serie.

library(readr)
gxs <- read_csv("gxserie_yaquis.csv")

## 
## -- Column specification --------------------------------------------------------
## cols(
##   Serie = col_double(),
##   Ganados = col_double()
## )

data(gxs)

## Warning in data(gxs): data set 'gxs' not found

str(gxs)

## tibble [10 x 2] (S3: spec_tbl_df/tbl_df/tbl/data.frame)
##  $ Serie  : num [1:10] 1 2 3 4 5 6 7 8 9 10
##  $ Ganados: num [1:10] 4 4 1 3 3 2 0 2 2 0
##  - attr(*, "spec")=
##   .. cols(
##   ..   Serie = col_double(),
##   ..   Ganados = col_double()
##   .. )

La variable gxs (ganados por serie) de los datos guarda el registro del número de veces que se gana de la forma, la tabla de frecuencias para gxs es:

table(gxs$Ganados)

## 
## 0 1 2 3 4 
## 2 1 3 2 2

Podemos observar los primeros registros de los datos de la siguiente forma:

head(gxs)

## # A tibble: 6 x 2
##   Serie Ganados
##   <dbl>   <dbl>
## 1     1       4
## 2     2       4
## 3     3       1
## 4     4       3
## 5     5       3
## 6     6       2

PG <- gxs$Ganados
head(PG)

## [1] 4 4 1 3 3 2

Como modelo se propone una cadena de markov finita de 3 estados para modelar los registros semanales del numero de veces que gana almacenados en gxs

P1 <- createSequenceMatrix(PG)
P1

##   0 1 2 3 4
## 0 0 0 1 0 0
## 1 0 0 0 1 0
## 2 2 0 1 0 0
## 3 0 0 1 1 0
## 4 0 1 0 0 1

Dicha función devuelve una lista con todos los resultados de la estimaciones, incluyendo un objeto markovchain que posee la matriz de transición. Para el ejemplo anterior se tienen los siguientes resultados.

Fit <-markovchainFit(data = PG, confidencelevel = 0.95)
Fit

## $estimate
## MLE Fit 
##  A  5 - dimensional discrete Markov Chain defined by the following states: 
##  0, 1, 2, 3, 4 
##  The transition matrix  (by rows)  is defined as follows: 
##           0   1         2   3   4
## 0 0.0000000 0.0 1.0000000 0.0 0.0
## 1 0.0000000 0.0 0.0000000 1.0 0.0
## 2 0.6666667 0.0 0.3333333 0.0 0.0
## 3 0.0000000 0.0 0.5000000 0.5 0.0
## 4 0.0000000 0.5 0.0000000 0.0 0.5
## 
## 
## $standardError
##           0   1         2   3   4
## 0 0.0000000 0.0 1.0000000 0.0 0.0
## 1 0.0000000 0.0 0.0000000 1.0 0.0
## 2 0.4714045 0.0 0.3333333 0.0 0.0
## 3 0.0000000 0.0 0.5000000 0.5 0.0
## 4 0.0000000 0.5 0.0000000 0.0 0.5
## 
## $confidenceLevel
## [1] 0.95
## 
## $lowerEndpointMatrix
##   0 1 2 3 4
## 0 0 0 0 0 0
## 1 0 0 0 0 0
## 2 0 0 0 0 0
## 3 0 0 0 0 0
## 4 0 0 0 0 0
## 
## $upperEndpointMatrix
##   0 1         2 3 4
## 0 0 0 1.0000000 0 0
## 1 0 0 0.0000000 1 0
## 2 1 0 0.9866548 0 0
## 3 0 0 1.0000000 1 0
## 4 0 1 0.0000000 0 1
## 
## $logLikelihood
## [1] -4.682131

mc <- Fit$estimate
summary(mc)

## MLE Fit  Markov chain that is composed by: 
## Closed classes: 
## 0 2 
## Recurrent classes: 
## {0,2}
## Transient classes: 
## {1},{3},{4}
## The Markov chain is not irreducible 
## The absorbing states are: NONE

El diagrama de la cadena markov

plot(mc)

# Y la distribución estacionaria de la matriz es:

steadyStates(mc)

##        0 1   2 3 4
## [1,] 0.4 0 0.6 0 0

Predicciones Una de las utilidades mas grandes de los Procesos Estocásticos, es su capacidad de realizar predicciones, esto se puede realizar mediante la función “predict()”, para eso observemos los ultimos registros del proceso

tail(PG)

## [1] 3 2 0 2 2 0

Se observa que en esta ultima serie se ganaron 0 juegos, como el modelo es una cadena de markov, para realizar una predicción de la siguiente serie, necesita conocer la serie actual, las prediciones para las siguientes* n = 3* series dado que en la utlima vez se registraron 0 ganados , son:

predict(mc, newdata="0",n.ahead=3)

## [1] "2" "0" "2"