División de los Datos

Como tenemos los datos de cada estado así como el nacional, vamos a separar los datos para obtener solo los datos a nivel Nacional

dataSetPecuario2020 <- subset(Pecuario_2020_Clean,  Pecuario_2020_Clean$Entidad %in% "Nacional")

dataSetPecuario2021 <- subset(Pecuario_2021_Clean,  Pecuario_2021_Clean$Entidad %in% "Nacional")

dataSetLluvia2020 <- subset(Lluvia_2020_Clean,  Lluvia_2020_Clean$Entidad %in% "Nacional")
dataSetLluvia2021 <- subset(Lluvia_2021_Clean, Lluvia_2021_Clean$Entidad %in% "Nacional")

newDataPecuario = rbind(dataSetPecuario2020,dataSetPecuario2021)

newDataLluvia = rbind(dataSetLluvia2020,dataSetLluvia2021)

Matriz de Diagramas de dispersión

Con le fin de conprender mejor el tema, una matriz de diagramas de dispersión, es basicamente una cuadricula la cual contiene varios gráficos, los cuales nos permiten conocer las relaciones existenes entre las variables. En dicha matríz podemos observar gráficos de dispersion de hasta cinco variables de procedencia numérica.

Una características de las matrices de diagramas de dispersión es que utilizan simbolos únicos para identificarse en caso de que no se utilice la variable “Color por”, por otro lado, si se llega asignar la variable “Color por”, se requiere por medio de valores únicos.

DatosNacional <- subset(Datos, Datos$Entidad %in% "Nacional")
pairs(DatosNacional)

Como podemos ver, los datos no parece haber una relación entre el Índice de Volumen Físico y los milimetros de lluvia, pero tal parece que puede existir una relación lineal entre la Fecha y el Índice de Volumen Físico.

Recta de los mínimos cuadrados

La recta de los mínimos cuadrados se puede definir como un método de regresión, el cual es una estrategia que se utiliza para modelar el comportamiento de todos los datos a travéz de una recta que cruza todos ellos. Esta recta se utiliza basicamente para observar y estudiar como se relacionan dos variables naturalmente.

Para encontrar la recta de los mínimos cuadrados es necesario seguir una serie de pasos:

1. Se necesita calcular la media de los valores de las variables X y Y

2. Sumar los cuadrados de la variable x

3. Se requiere sumar cada valor de x, este multiplicado por el valor y que le corresponde

4. Calcular la pendiente de la recta (Debemos utilizar la siguiente fórmula);

,donde n es el número total de puntos de los datos

5. Debemos de calcular la intercepción de la recta (Utilizar la siguiente fórmula):

,donde y y x corresponde a las medias de as coordenadas de los puntos de datos correspondientes.

6. Se utiliza la pendiente e interceptación para en y para fórmular la ecuación de la recta.

Ahora vamos a ver la correlación entre nuestra variable de Indice de Volumen Físico en función de la cantidad de lluvia

regresion <- lm(Datos$IVF ~ Datos$MM, data = Datos)
summary(regresion)

## 
## Call:
## lm(formula = Datos$IVF ~ Datos$MM, data = Datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -41.213 -13.489  -2.639   8.902 103.713 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 1.157e+02  9.776e-01 118.313   <2e-16 ***
## Datos$MM    7.046e-03  7.520e-03   0.937    0.349    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 22.19 on 790 degrees of freedom
## Multiple R-squared:  0.00111,    Adjusted R-squared:  -0.0001542 
## F-statistic: 0.8781 on 1 and 790 DF,  p-value: 0.349

Con el modelo utilizando las tres variables obtenemos una R^2 ajustada de -0.0001542, lo que este modelo no es capaz de explicar en función de la lluvia, la producción pecuaria, por lo que tendremos que utilizar otro método que nos permita estudiar mejor la relación entre ambos datos.

Buscando datos Faltantes

anyNA(newDataLluvia)

## [1] FALSE

Esto signigica que no tenemos ningún valor faltante.

Descomposición de la serie de tiempo

La descomposición de series de tiempo es el proceso que identifica y calcula las diversas componentes existentes en una serie temporal, así como la forma en que estas se relacionan entre sí. Por ende, se procederá hacer una descomposición de serie de tiempo de las precipitaciones mensuales en MM de México.

datosLluvia<-ts(newDataLluvia$MM , frequency = 2 ,start = c(2020, 1),end = c(2021,12))
datos_dcLluvia<- datosLluvia %>%
  decompose(type="multiplicative") %>%
  autoplot()
datos_dcLluvia

1. Usando un modelo de suavizado exponencial

• Validación cruzada

La validación cruzada es una técnica para evaluar modelos de ML mediante el entrenamiento de varios modelos de ML en subconjuntos de los datos de entrada disponibles y evaluarlos con el subconjunto complementario de los datos. Por lo tanto, los datos de entrenamiento y test se procederá a dividirlos en un 20 y 80 porciento respectivamente.

test_datosLluvia <- tail(datosLluvia, 3) 
train_datosLluvia <- head(datosLluvia, length(datosLluvia)-3) 

• Holt-Winters

El método de Holt Winters es utilizado para realizar pronósticos del comportamiento de una serie temporal a partir de los datos obtenidos anteriormente. El método se basa en un algoritmo iterativo que a cada tiempo realiza un pronóstico sobre el comportamiento de la serie en base a promedios debidamente ponderados de los datos anteriores.

datos_sesLluvia <- HoltWinters(train_datosLluvia, seasonal = "multiplicative", )
datos_sesLluvia

## Holt-Winters exponential smoothing with trend and multiplicative seasonal component.
## 
## Call:
## HoltWinters(x = train_datosLluvia, seasonal = "multiplicative")
## 
## Smoothing parameters:
##  alpha: 1
##  beta : 0
##  gamma: 1
## 
## Coefficients:
##         [,1]
## a  35.990891
## b  -2.325000
## s1  0.983076
## s2  1.016924

• Forecasting

Forecasting es el proceso de hacer predicciones sobre el futuro mediante el análisis estadístico de tendencias observadas en datos históricos.

datos_forecastLluvia <- forecast(object = datos_sesLluvia, h = 3)

• Visualizando resultados

Por ende una vez aplicado el proceso de hacer predicciones con Forecasting, se procederá hacer la visualización del resultado en forma gráfica, viendo los datos de entrenamiento vs datos de prueba vs datos reales.

datosLluvia %>%
  autoplot(series = "Datos Reales") +
  autolayer(datos_forecastLluvia$fitted, series = "Datos Entrenamiento") +
  autolayer(datos_forecastLluvia$mean , series = "Datos de Prueba") +
  theme_minimal()

## Warning: Removed 2 row(s) containing missing values (geom_path).

Evaluar la precisión del modelo

Una vez realizada la prediccion se procedera a evluar la precision del modelo.

eval_sesLluvias <- accuracy(datos_forecastLluvia, test_datosLluvia)
eval_sesLluvias

##                      ME     RMSE      MAE         MPE      MAPE      MASE
## Training set   3.696165 33.56470 23.69206   -6.077623  39.43366 0.5564418
## Test set     -14.964087 15.82444 14.96409 -130.078319 130.07832 0.3514530
##                    ACF1 Theil's U
## Training set  0.3785737        NA
## Test set     -0.5806394  1.283706

Deacuerdo con la prueba el modelo tiene un error medio del 130% Por lo que nuestro modelo no es capaz de predecir el comportamiento de las precipitaciones de las lluvias en Mexico a nivel NAcional.

Utilizando el modelo SARIMA o ARIMA

Los modelos SARIMA captan el comportamiento puramente estacional de una serie, en forma similar, como hemos visto, se realiza para la componente regular o no estacional. Una serie con influencia solamente por la componente estacional puede ser descrito por un modelo SARIMA (P,D,Q)

Empezaremos haciendo una prueba al conjunto de datos.

adf.test(datosLluvia)

## 
##  Augmented Dickey-Fuller Test
## 
## data:  datosLluvia
## Dickey-Fuller = -1.5945, Lag order = 2, p-value = 0.7269
## alternative hypothesis: stationary

Deacuerdo con el valor de P, los datos no son estacionarios

Análisis de ajuste

En esta seccion se representan los datos de una manera grafica con el apoyo de autoplot().

datos_diffLluvia <- datosLluvia

datos_diffLluvia %>%
  autoplot()

Analizando la serie de tiempo por medio de un modelo ARIMA

datos_autoLluvia <- auto.arima(y = train_datosLluvia, seasonal = T)
summary(datos_autoLluvia)

## Series: train_datosLluvia 
## ARIMA(0,1,0) 
## 
## sigma^2 = 1048:  log likelihood = -48.96
## AIC=99.92   AICc=100.42   BIC=100.22
## 
## Training set error measures:
##                  ME     RMSE      MAE       MPE     MAPE      MASE      ACF1
## Training set 1.3656 30.86077 19.60196 -8.681991 35.08629 0.4603801 0.3355203

Deacuerdo con la prueba del modeloARima tiene un error medio del 35% Por lo que nuestro modelo es capaz de predecir el comportamiento de las precipitaciones de las lluvias en Mexico a nivel Nacional, sin embargo, no sera tan exacta la precision ya que un 35.08% es un gran numero.

Visualizacion grafica de los datos reales versus el modelo SARIMA

Para observar mejor los resultados de la evaluacion anterior se procedera a visualizarlos de una forma grafica.

train_datosLluvia %>%
  autoplot(series= "Datos Reales") +
  autolayer (datos_autoLluvia$fitted, series = "SARIMA auto") +
  theme_minimal()

Analisis de los cambios diferenciales en la serie de tiempo

En esta seccion se muestra los cambios diferenciales en la serie de tiempo de tiempo evaluada.

tsdisplay(datos_diffLluvia)

Buscando datos Faltantes

anyNA(newDataPecuario)

## [1] FALSE

Esto signigica que no tenemos ningún valor faltante.

Descomposición de la serie de tiempo

La descomposición de series de tiempo es el proceso que identifica y calcula las diversas componentes existentes en una serie temporal, así como la forma en que estas se relacionan entre sí. Por ende, se procederá hacer una descomposición de serie de tiempo de las precipitaciones mensuales en MM de México.

datosPecuario<-ts(newDataPecuario$IVF , frequency = 2 ,start = c(2020, 1),end = c(2021,12))
datos_dc <- datosLluvia %>%
  decompose(type="multiplicative") %>%
  autoplot()
datos_dc

1. Usando un modelo de suavizado exponencial

• Validación cruzada

La validación cruzada es una técnica para evaluar modelos de ML mediante el entrenamiento de varios modelos de ML en subconjuntos de los datos de entrada disponibles y evaluarlos con el subconjunto complementario de los datos. Por lo tanto, los datos de entrenamiento y test se procederá a dividirlos en un 20 y 80 porciento respectivamente.

test_datosPecuario <- tail(datosPecuario, 3) # 20 porciento de los datos para pruebas
train_datosPecuario <- head(datosPecuario, length(datosPecuario)-3) # 80 porciento de los datos para entrenamiento

• Holt-Winters

El método de Holt Winters es utilizado para realizar pronósticos del comportamiento de una serie temporal a partir de los datos obtenidos anteriormente. El método se basa en un algoritmo iterativo que a cada tiempo realiza un pronóstico sobre el comportamiento de la serie en base a promedios debidamente ponderados de los datos anteriores.

datos_sesPecuario <- HoltWinters(train_datosPecuario, seasonal = "multiplicative", )
datos_sesPecuario

## Holt-Winters exponential smoothing with trend and multiplicative seasonal component.
## 
## Call:
## HoltWinters(x = train_datosPecuario, seasonal = "multiplicative")
## 
## Smoothing parameters:
##  alpha: 1
##  beta : 0.02141579
##  gamma: 1
## 
## Coefficients:
##           [,1]
## a  125.8119237
## b    1.2242424
## s1   0.9977103
## s2   1.0022897

• Forecasting

Forecasting es el proceso de hacer predicciones sobre el futuro mediante el análisis estadístico de tendencias observadas en datos históricos.

datos_forecastPecuario <- forecast(object = datos_sesPecuario, h = 3)

• Visualizando resultados

Por ende una vez aplicado el proceso de hacer predicciones con Forecasting, se procederá hacer la visualización del resultado en forma gráfica, viendo los datos de entrenamiento vs datos de prueba vs datos reales.

datosPecuario %>%
  autoplot(series = "Datos Reales") +
  autolayer(datos_forecastPecuario$fitted, series = "Datos Entrenamiento") +
  autolayer(datos_forecastPecuario$mean , series = "Datos de Prueba") +
  theme_minimal()

## Warning: Removed 2 row(s) containing missing values (geom_path).

Evaluar la precisión del modelo

Una vez realizada la prediccion se procedera a evluar la precision del modelo.

eval_sesPecuario <- accuracy(datos_forecastPecuario, test_datosPecuario)
eval_sesPecuario

##                      ME       RMSE        MAE        MPE      MAPE      MASE
## Training set  0.2559506  0.7664347  0.6791311  0.2177197 0.5659706 0.2589907
## Test set     -8.2958476 11.3969436 10.1323232 -7.3001261 8.7182539 3.8640216
##                    ACF1 Theil's U
## Training set  0.1483145        NA
## Test set     -0.1747509  1.394874

Deacuerdo con la prueba de precisión, el modelo tiene un error medio del 8.17% Por lo que nuestro modelo es capaz de predecir el comportamiento del indice del volumen fisico en Mexico a nivel Nacional.

Utilizando el modelo SARIMA o ARIMA

Los modelos SARIMA captan el comportamiento puramente estacional de una serie, en forma similar, como hemos visto, se realiza para la componente regular o no estacional. Una serie con influencia solamente por la componente estacional puede ser descrito por un modelo SARIMA (P,D,Q)

Empezaremos haciendo una prueba al conjunto de datos

adf.test(datosPecuario)

## Warning in adf.test(datosPecuario): p-value greater than printed p-value

## 
##  Augmented Dickey-Fuller Test
## 
## data:  datosPecuario
## Dickey-Fuller = 0.20028, Lag order = 2, p-value = 0.99
## alternative hypothesis: stationary

Deacuerdo con el valor de P, los datos no son estacionarios

Análisis de ajuste

En esta seccion se representan los datos de una manera grafica con el apoyo de autoplot().

datos_diffPecuario <- datosPecuario

datos_diffPecuario %>%
  autoplot()

Analizando la serie de tiempo por medio de un modelo ARIMA

datos_autoPecuario <- auto.arima(y = train_datosPecuario, seasonal = T)
summary(datos_autoPecuario)

## Series: train_datosPecuario 
## ARIMA(0,0,0)(0,1,0)[2] with drift 
## 
## Coefficients:
##        drift
##       1.3111
## s.e.  0.1687
## 
## sigma^2 = 1.155:  log likelihood = -12.88
## AIC=29.75   AICc=31.75   BIC=30.15
## 
## Training set error measures:
##                      ME      RMSE       MAE        MPE      MAPE      MASE
## Training set 0.02020605 0.9165279 0.6545495 0.01936489 0.5396498 0.2496163
##                   ACF1
## Training set 0.4263745

Deacuerdo con la prueba del modeloARima tiene un error medio del 0.5396% por lo que nuestro modelo es capaz de predecir el comportamiento del indice del volumen fisico en Mexico a nivel Nacional.

Visualizacion grafica de los datos reales versus el modelo SARIMA

Para observar mejor los resultados de la evaluacion anterior se procedera a visualizarlos de una forma grafica.

train_datosPecuario %>%
  autoplot(series= "Datos Reales") +
  autolayer (datos_autoPecuario$fitted, series = "SARIMA auto") +
  theme_minimal()

Analisis de los cambios diferenciales en la serie de tiempo

En esta seccion se muestra los cambios diferenciales en la serie de tiempo de tiempo evaluada.

tsdisplay(datos_diffPecuario)