Punto 1

Descripcion de la serie elegida, asi como su frecuencia, unidad de medida,fuente, etc.

Caracteristicas

Ruta temática: Energía > Producción de petroquímicos por entidad federativa y principales productos > Total

Periocidad mensual
Unidad de medida toneladadas
Fuente PEMEX. Dirección Corporativa de Planeación Estratégica.
Fecha inicial 2010/01
Fecha final 2018/01
Última actualización 2018/02/28

Fuente de información: Enlace a Banco de Información Financiera

Punto 2

Graficar los datos, analizar patrones y observaciones atipicas.

Grafica Original

library(urca)
library(forecast)
library(TSA)
## 
## Attaching package: 'TSA'
## The following objects are masked from 'package:stats':
## 
##     acf, arima
## The following object is masked from 'package:utils':
## 
##     tar
library(ggplot2)


prodp <- read.csv('E:/basest.csv')
prodp<-ts(prodp$Total, frequency = 12, start= c(2010,1))
autoplot(prodp, colors="darkblue", main= "Producción de petroquímicos por entidad federativa y principales productos. ", ylab="toneladas", xlab= "años")

En esta grafica se muestran los datos originales, comprobando la estacionalidad de serie de tiempo y los ciclos que muestra.

Descomposicion

library(forecast)
Fit<-decompose(prodp)
autoplot(Fit, col= "darkblue")

-Componentes de la serie de tiempo

Tendencia - Ciclo

El comportamineto de la producción petroquimica en el pais a través del tiempo, presenta una tendencia NEGATIVA.
Actualmente se tienen plantas fuera de operación por falta de competitividad, de mercado o falta de materia prima, las cuales siguen formando parte de los activos de Pemex Petroquímica.
El cierre obedece a que su operación no genera ingresos suficientes para cubrir ni siquiera los costos variables, lo que repercute en una menor produccion petroquimica
El volumen nacional de producción de petroquímicos cayó 10.7% en un año, afectado por la caída en el suministro de materia prima de Petróleos Mexicanos (Pemex), la falta de estos insumos precursores ha provocado que las fabricas privadas de petroquímicos en el país operen a una tercera parte de su capacidad en los últimos años, según diversos industriales. 
plot(Fit$trend, col= "turquoise2")

Estacionalidad:

Como la estacionalidad siempre es fija y conocida, muestra como a través de los años disminuye la producción de dicho bien.

plot(Fit$seasonal, col= "darkmagenta")

Observaciones atipicas

ggseasonplot(prodp, main= ("Producción de petroquímicos por entidad federativa y principales productos."), polar="TRUE")

Esta gráfica se puede observar mejor la estacionalidad, mostrando que en se tiene una menor producción petroquimica en en el país como van pasando los años, y el principal factor fue que la producción de petroquímicos cayó 10.7% en un año, afectado por la caída en el suministro de materia prima de Pemex, que descendió 11.6%, según el último anuario de la Asociación Nacional de la Industria Química (ANIQ).

La produccion petroquimica hilo en el 2015 dos años a la baja.

Enlace a “Mexico: Produccion petroqimica” 2016

Punto 3

Si es necesario, utilizar transformacion Box-Cox / logaritmos para estabilizar la varianza.

Tranformacion Box-Cox

BoxCox.ar(prodp)
## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

## Warning in arima0(x, order = c(i, 0L, 0L), include.mean = demean): possible
## convergence problem: optim gave code = 1

El intervalo a 95% para \(\?? \??\), contiene un valor de \(\??=0\??=0\) muy cerca de su centro y sugiere una transformación logarítmica.

Tranformacion logaritmos

ts.prodp <- (log(prodp))
plot(log(ts.prodp), main="Producción de petroquímicos por entidad federativa y principales productos.", col="chartreuse2", ylab="toneladas", xlab="años", lwd=2, type="o", pch=1)

Punto 4

En caso de estacionalidad, aplicar diferencias estacionales

Primera Diferencia

difpp <- diff(ts.prodp)
plot(diff(log(prodp)), col="deeppink2")

Segunda Diferencia

dif2pp <- diff(prodp,difference=2)
plot(diff(prodp,difference=2),col="cyan3")

Media Estable

library(gridExtra)
library(forecast)

logpp<- autoplot (log(prodp)) 
pd<- autoplot (diff(log(prodp)))
sd<- autoplot (diff(prodp, difference=2))

grid.arrange(logpp, pd, sd)

Punto 5

Usar prueba Dickey-Fuller para evaluar el orden de integracion de la serie.

library(urca)

summary(ur.df(dif2pp))
## 
## ############################################### 
## # Augmented Dickey-Fuller Test Unit Root Test # 
## ############################################### 
## 
## Test regression none 
## 
## 
## Call:
## lm(formula = z.diff ~ z.lag.1 - 1 + z.diff.lag)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -177143  -59009    7507   44679  386752 
## 
## Coefficients:
##            Estimate Std. Error t value Pr(>|t|)    
## z.lag.1    -2.22668    0.17770 -12.531  < 2e-16 ***
## z.diff.lag  0.34183    0.09606   3.558 0.000595 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 89300 on 91 degrees of freedom
## Multiple R-squared:  0.854,  Adjusted R-squared:  0.8508 
## F-statistic: 266.2 on 2 and 91 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic is: -12.5306 
## 
## Critical values for test statistics: 
##      1pct  5pct 10pct
## tau1 -2.6 -1.95 -1.61

Para que la serie sea estacionaria: phi tiene que ser menor a 1.

\(H_0=\) No es estacionaria. \(H_1=\) Ya es estacionaria.

Al tomar la segunada diferencia, se muestra que la serie ya es estacinaria, ya que el valor estadistico es de -12.5306, es mayor a los valores críticos de Tau.

Punto 6

Usar herramientas ACF, PACF, EACF y criterios de Akaike/ Bayes para construirpropuestas de modelos

ACF

ggAcf(dif2pp)

PACF

ggPacf(dif2pp)

EACF

eacf(dif2pp)
## AR/MA
##   0 1 2 3 4 5 6 7 8 9 10 11 12 13
## 0 x x o o o o o o o o x  x  x  o 
## 1 x o o x o o o o o o o  x  x  o 
## 2 x o x o o o o o o o o  x  o  o 
## 3 x x x x o o o o o o o  x  o  o 
## 4 x o o x x o o o o o o  o  o  o 
## 5 x o x o x o o o o o o  o  o  o 
## 6 o o x o x o o o o o o  o  o  o 
## 7 o o x o x o x o o o o  x  o  o
ggtsdisplay(dif2pp)

Akaike/ Bayes

res <- armasubsets(y=dif2pp, nar=12, nma=12, y.name='test', ar.method='ols')
plot(res)

Punto 7

Como parte del diagnostico del modelo, analizar los residuos y aplicar la prueba Ljung-Box

Analisis de los residuos

propuesta1 <- arima(dif2pp, order=c(2,0,0),method='ML') # AR(1)
propuesta1
## 
## Call:
## arima(x = dif2pp, order = c(2, 0, 0), method = "ML")
## 
## Coefficients:
##           ar1      ar2  intercept
##       -0.9126  -0.3515  -184.6596
## s.e.   0.0961   0.0974  4089.6492
## 
## sigma^2 estimated as 8.02e+09:  log likelihood = -1218.48,  aic = 2442.97
propuesta2 <- arima(dif2pp, order=c(0,0,1),method='ML') # MA(1)
propuesta2
## 
## Call:
## arima(x = dif2pp, order = c(0, 0, 1), method = "ML")
## 
## Coefficients:
##           ma1  intercept
##       -1.0000  -100.8005
## s.e.   0.0266   289.6591
## 
## sigma^2 estimated as 6.171e+09:  log likelihood = -1207.88,  aic = 2419.77
propuesta3 <- arima(dif2pp, order=c(4,0,1),method='ML') # ARMA(4,1)
propuesta3
## 
## Call:
## arima(x = dif2pp, order = c(4, 0, 1), method = "ML")
## 
## Coefficients:
##           ar1      ar2      ar3      ar4      ma1  intercept
##       -0.4113  -0.0841  -0.2131  -0.2862  -1.0000  -118.0359
## s.e.   0.0979   0.1059   0.1047   0.0986   0.0304   133.3103
## 
## sigma^2 estimated as 4.781e+09:  log likelihood = -1196.69,  aic = 2405.39

Residuos estandarizados

library(forecast)
library(TSA)
library(ggplot2)



plot(rstandard(propuesta1),ylab ='Residuos estandarizados', type='o'); abline(h=0)

plot(rstandard(propuesta2),ylab ='Residuos estandarizados', type='o'); abline(h=0)

plot(rstandard(propuesta3),ylab ='Residuos estandarizados', type='o'); abline(h=0)

Prueba Ljung-Box

tsdiag(propuesta1,gof=8,omit.initial=F)

tsdiag(propuesta2,gof=8,omit.initial=F)

tsdiag(propuesta3,gof=8,omit.initial=F)

Revision de los residuales

checkresiduals(propuesta1)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(2,0,0) with non-zero mean
## Q* = 41.606, df = 16, p-value = 0.0004517
## 
## Model df: 3.   Total lags used: 19
checkresiduals(propuesta2)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(0,0,1) with non-zero mean
## Q* = 57.342, df = 17, p-value = 2.866e-06
## 
## Model df: 2.   Total lags used: 19
checkresiduals(propuesta3)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(4,0,1) with non-zero mean
## Q* = 28.351, df = 13, p-value = 0.008082
## 
## Model df: 6.   Total lags used: 19

Punto 8

Usar la funcion auto.arima()y compare sus resultados con el inciso anterior

propuesta.auto <- auto.arima(dif2pp)
propuesta.auto
## Series: dif2pp 
## ARIMA(5,0,0)(1,0,0)[12] with zero mean 
## 
## Coefficients:
##           ar1      ar2      ar3      ar4      ar5   sar1
##       -1.1011  -0.8857  -0.7526  -0.6011  -0.2857  0.389
## s.e.   0.0984   0.1394   0.1473   0.1348   0.0972  0.103
## 
## sigma^2 estimated as 5.864e+09:  log likelihood=-1201.95
## AIC=2417.89   AICc=2419.18   BIC=2435.77
plot(rstandard(propuesta.auto),ylab ='Residuos estandarizados', type='o'); abline(h=0)

qqnorm(residuals(propuesta.auto)); qqline(residuals(propuesta.auto))

ggAcf(residuals(propuesta.auto))

tsdiag(propuesta.auto, gof=12, omit.initial=F)

checkresiduals(propuesta.auto)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(5,0,0)(1,0,0)[12] with zero mean
## Q* = 21.25, df = 13, p-value = 0.06814
## 
## Model df: 6.   Total lags used: 19

Punto 9

Ecuacion final

\[ Yt= -1.1011e_1-0.8857e_{t-1}-0.7526e_{t-2}-0.6011e_{t-3}-0.2857e_{t-4}+0.389+et \]

Punto 10

Crear un pronostico a dos años.

pronostico <- plot(forecast(propuesta.auto,h=24))

pronostico
## $mean
##              Jan         Feb         Mar         Apr         May
## 2018             -105693.695  106909.577  -82951.732   49781.701
## 2019   24214.907  -40395.127   42831.298  -33152.830   19855.574
## 2020    9531.214                                                
##              Jun         Jul         Aug         Sep         Oct
## 2018   -5348.718  -21226.080  -13582.947  -17736.387   47328.410
## 2019   -2653.692   -8347.996   -4869.311   -6887.387   18303.547
## 2020                                                            
##              Nov         Dec
## 2018  -15550.820   -3450.926
## 2019   -6035.498   -1494.994
## 2020                        
## 
## $lower
##                80%       95%
## Feb 2018 -203827.4 -255776.2
## Mar 2018  -39055.8 -116325.2
## Apr 2018 -232397.4 -311509.1
## May 2018 -100269.0 -179701.0
## Jun 2018 -155654.9 -235222.1
## Jul 2018 -172284.9 -252250.6
## Aug 2018 -164664.4 -244642.0
## Sep 2018 -169900.7 -250451.7
## Oct 2018 -105133.1 -185841.3
## Nov 2018 -168065.2 -248801.4
## Dec 2018 -156121.3 -236940.2
## Jan 2019 -128456.1 -209275.3
## Feb 2019 -196820.3 -279626.7
## Mar 2019 -119500.9 -205434.4
## Apr 2019 -196112.0 -282377.3
## May 2019 -143185.3 -229494.0
## Jun 2019 -165774.2 -252125.0
## Jul 2019 -171532.4 -257917.0
## Aug 2019 -168064.1 -254454.2
## Sep 2019 -170224.9 -256690.6
## Oct 2019 -145086.0 -231579.3
## Nov 2019 -169430.3 -255926.3
## Dec 2019 -164912.6 -251420.7
## Jan 2020 -153886.7 -240394.9
## 
## $upper
##                 80%       95%
## Feb 2018  -7560.036  44388.78
## Mar 2018 252874.952 330144.35
## Apr 2018  66493.894 145605.63
## May 2018 199832.398 279264.44
## Jun 2018 144957.432 224524.70
## Jul 2018 129832.732 209798.44
## Aug 2018 137498.477 217476.15
## Sep 2018 134427.964 214978.90
## Oct 2018 199789.921 280498.17
## Nov 2018 136963.540 217699.76
## Dec 2018 149219.475 230038.30
## Jan 2019 176885.930 257705.09
## Feb 2019 116030.006 198836.47
## Mar 2019 205163.503 291096.98
## Apr 2019 129806.299 216071.65
## May 2019 182896.483 269205.13
## Jun 2019 160466.801 246817.57
## Jul 2019 154836.409 241221.01
## Aug 2019 158325.489 244715.60
## Sep 2019 156450.144 242915.81
## Oct 2019 181693.135 268186.36
## Nov 2019 157359.333 243855.33
## Dec 2019 161922.628 248430.69
## Jan 2020 172949.149 259457.38