Series de Tiempo

Divorcios en la CDMX

Meza Pale Pablo

1. Breve referencia de la serie de tiempo elegida, así como frecuencia, unidad de medida, fuente, etc.

Nuestra serie es sobre el número de divorcios mensuales, que ocurrieron en la Ciudad de México en el periodo que abarca del año 2000 al año 2016. La unidad de medida es por divorcio ocurrido. FUENTE: INEGI. Estadísticas de Registros administrativos; Estadísticas de nupcialidad: Divorcios.

2. Graficar la serie de tiempo, analizar patrones y observaciones atípicas. Apoye su análisis con una descomposición clásica.

base1<-read.csv("C:/Users/user/Documents/economia/Sexto semestre/Series de tiempo/Seriedivorcios.csv")
divorcios<-ts(base1$divorcios,start = c(2000,1), frequency = 12)
plot(divorcios, xlab="Periodo", main= "Divorcios")

summary(divorcios)
##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
##    1.00   49.75   92.50   90.03  130.00  177.00

Analizando ésta gráfica, observamos que nuestra serie no presenta tantas complicaciones sin embargo, si presenta observaciones atípicas; si vemos el gráfico, hubo un gran incremento en el número de divorcios a partir de 2008, esto debido a una reforma en el código civil de la CDMX que permitía decretar el divorcio de manera expedita, lo que facilitó este trámite en gran medida. Podemos intuir que nuestra serie contiene patrones de “estacionalidad” y ligera “tendencia”.

library(ggplot2)
library(ggfortify)
library(forecast)
fit<-decompose(divorcios, type = "multiplicative")
autoplot(fit, ts.colour = "red", ts.linetype = "dashed")

Utilizando la descomposición clásica podemos ver que existe la estacionalidad y una tendencia, en particular en los últimos años de la serie, asimismo notamos que no hay errores tan graves en nuestra serie.

3. Si es necesario, utilizar transformación Box-Cox / logaritmos para estabilizar la varianza.

Presentar gráfica.

library(TSA)
## Loading required package: leaps
## Loading required package: locfit
## locfit 1.5-9.1    2013-03-22
## Loading required package: mgcv
## Loading required package: nlme
## 
## Attaching package: 'nlme'
## The following object is masked from 'package:forecast':
## 
##     getResponse
## This is mgcv 1.8-23. For overview type 'help("mgcv-package")'.
## Loading required package: tseries
## 
## Attaching package: 'TSA'
## The following objects are masked from 'package:stats':
## 
##     acf, arima
## The following object is masked from 'package:utils':
## 
##     tar
BoxCox.ar(divorcios)

Observamos que el valor de ?? se acerca mucho a 1, lo que sugiere que no requiere una transformación logarítmica.

4. Modelar la serie de tiempo elegida como función de una tendencia, tendencua cuadrática y/o medias estacionales. Presente el diagnóstico del modelo elegido.

modelo1<-lm(divorcios~time(divorcios))
summary(modelo1)
## 
## Call:
## lm(formula = divorcios ~ time(divorcios))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -84.379 -44.588  -1.667  37.109  90.279 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)
## (Intercept)     2309.5866  1420.8038   1.626    0.106
## time(divorcios)   -1.1051     0.7074  -1.562    0.120
## 
## Residual standard error: 49.58 on 202 degrees of freedom
## Multiple R-squared:  0.01194,    Adjusted R-squared:  0.007046 
## F-statistic:  2.44 on 1 and 202 DF,  p-value: 0.1198

En esta primer regresión podemos ver que nuestro intercepto es de 2309.58 y nuestra variable independiente es de -1.10 lo que significa que cada año que pasa disminuyen en 1.10 los divorcios en la CDMX. De igual manera podemos observar que el tiempo no es significativo. Vemos que nuestro error estándar es bastante grande y nuestra R cuadrada muy pequeña por lo que nuestra estimación no es eficiente.

A continuación se hizo un modelo con un ajuste cuadrático para mejorar el ajuste.

modelo2<-lm(divorcios~time(divorcios)+ I(time(divorcios)^2))
summary(modelo2)
## 
## Call:
## lm(formula = divorcios ~ time(divorcios) + I(time(divorcios)^2))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -90.831 -41.748  -3.356  39.118  93.238 
## 
## Coefficients:
##                        Estimate Std. Error t value Pr(>|t|)  
## (Intercept)          -1.299e+06  6.453e+05  -2.013   0.0455 *
## time(divorcios)       1.294e+03  6.426e+02   2.015   0.0453 *
## I(time(divorcios)^2) -3.225e-01  1.600e-01  -2.016   0.0451 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 49.21 on 201 degrees of freedom
## Multiple R-squared:  0.03153,    Adjusted R-squared:  0.02189 
## F-statistic: 3.272 on 2 and 201 DF,  p-value: 0.03998

Con este ajuste vemos que está mejor nuestro modelo, podemos ver que el tiempo, y el tiempo al cuadrado son significativos al 1% y aumentó nuestra R cuadrada, de igual forma nuestro error estándar disminuyó considerablemente.

ggplot(divorcios, aes(y=divorcios, x=time(divorcios)))+geom_point()+geom_smooth(se=FALSE)
## Don't know how to automatically pick scale for object of type ts. Defaulting to continuous.
## Don't know how to automatically pick scale for object of type ts. Defaulting to continuous.
## `geom_smooth()` using method = 'loess'

autoplot(modelo1, which = 1:6, ncol = 2, label.size = 3)

En este gráfico podemos ver los residuales de nuestra serie, los cuales son muy diferentes con respecto a los ajustados, vemos que valores como el 143 afectan nuestra serie, vemos que nuestra gráfica quantil-quantil no es buena ya que hay algunos valores que afectan, analizando la distancia de Cook podemos observar que hay momentos de gran importancia casi al final de la serie.

Prueba Shapiro-Wilk

shapiro.test(rstudent(modelo1))
## 
##  Shapiro-Wilk normality test
## 
## data:  rstudent(modelo1)
## W = 0.96201, p-value = 2.731e-05

En esta primera prueba vemos que el componente estocástico no está normalmente distribuido por lo que se rechaza la hipótesis nula ya que el valor P es de: 2.731e-05

shapiro.test(rstudent(modelo2))
## 
##  Shapiro-Wilk normality test
## 
## data:  rstudent(modelo2)
## W = 0.96803, p-value = 0.0001373

En esta prueba también se rechaza la hipótesis nula, de igual manera, el componente estocástico no está normamlmente distribuido. Aunque el valor p, es ligeramente mayor que el anterior.

Correlograma (función de autocorrelación)

ggAcf(rstudent(modelo1))

Aquí vemos los errores estándar de nuestra serie, vemos que no son tantos y que la mayoría de los componentes se comportan de manera correcta.

ggAcf(rstudent(modelo2))

Este es el Correlograma de la segunda regresión y vemos que tiene mejores resultados con respecto a la anterior.

Usando la serie elegida como variable dependiente y dummys estacionales como independientes.

mes<-season(divorcios)
modelo3<-lm(divorcios~mes)
summary(modelo3)
## 
## Call:
## lm(formula = divorcios ~ mes)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -101.176  -31.500   -1.647   40.221   92.059 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)    92.235     11.922   7.736 5.65e-13 ***
## mesFebruary    11.941     16.861   0.708   0.4797    
## mesMarch       13.647     16.861   0.809   0.4193    
## mesApril        4.000     16.861   0.237   0.8127    
## mesMay        -10.824     16.861  -0.642   0.5217    
## mesJune         2.647     16.861   0.157   0.8754    
## mesJuly       -17.235     16.861  -1.022   0.3080    
## mesAugust      -3.588     16.861  -0.213   0.8317    
## mesSeptember    8.824     16.861   0.523   0.6014    
## mesOctober    -15.294     16.861  -0.907   0.3655    
## mesNovember    12.412     16.861   0.736   0.4625    
## mesDecember   -33.000     16.861  -1.957   0.0518 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 49.16 on 192 degrees of freedom
## Multiple R-squared:  0.07696,    Adjusted R-squared:  0.02408 
## F-statistic: 1.455 on 11 and 192 DF,  p-value: 0.1514

En este modelo observamos que ninguna de las variables es significativa más que el intercepto. Vemos que hay un error estándar grande y una R cuadrada muy pequeña.

autoplot(modelo3, which = 1:6, ncol = 2, label.size = 3)

En estos gráficos vemos que nuestros residuos no se comportan como un componente estocástico. La distancia de Cook muestra que hay momentos importantes en la seri(153, 182)

shapiro.test(rstudent(modelo3))
## 
##  Shapiro-Wilk normality test
## 
## data:  rstudent(modelo3)
## W = 0.96993, p-value = 0.0002356

Con esta prueba vemos que se rechaza la hipótesis nula de normalidad ya que el valor p es de: .0002356

ggAcf(rstudent(modelo3))

En este Correlograma podemos analizar los errores estándar de nuestra serie y vemos que aquí hay más errores que con los otros procedimientos.

De acuerdo a todos los modelos que hicimos, las pruebas que se les aplicaron, el análisis de los residuales, los diagnósticos de los modelos, creemos que el mejor modelo para nuestra serie, sería el modelo de tendencia y estacionalidad, ya que es el que presenta un mejor ajuste a nuestra serie.

5. En caso de estacionalidad de la serie de tiempo, aplicar diferencias estacionales. Presentar gráfica.

STDIFF<-diff(divorcios, 12)
autoplot(STDIFF)

Se le aplicó a la serie una primera diferencia estacional.

6. Usar prueba Dickey-Fuller para evaluar el orden de integración de la serie. Diferenciar hasta que la serie sea estacionaria.

library(urca)
summary(ur.df(STDIFF, type = "none",selectlags = "AIC"))
## 
## ############################################### 
## # 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 
## -160.585  -20.766    2.263   21.804  192.752 
## 
## Coefficients:
##            Estimate Std. Error t value Pr(>|t|)    
## z.lag.1    -0.59797    0.08681  -6.888 8.24e-11 ***
## z.diff.lag -0.18231    0.07175  -2.541   0.0119 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 53.1 on 188 degrees of freedom
## Multiple R-squared:  0.3865, Adjusted R-squared:   0.38 
## F-statistic: 59.21 on 2 and 188 DF,  p-value: < 2.2e-16
## 
## 
## Value of test-statistic is: -6.8883 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.58 -1.95 -1.62

Con la prueba Dickey-Fuller comparamos el valor estadístico, el cual es de -6.8883 con los valores críticos, por lo tanto se rechaza la hipótesis nula de que no es estacionaria, por lo que ya podemos empezar a crear modelos.

7. Usar herramientas ACF, PACF, EACF y/o criterios de Akaike/ Bayes para construir propuestas de modelos.

ggtsdisplay(diff(divorcios,12), main='Dif estacional del no. de divorcios')

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

La Función de Autocorrelación sugiere un MA (3) con un MA (1) estacional. La Función de Autocorrelación Parcial sugiere un AR (1) con un AR (1) estacional. La Función de autocorrelación extendida sugiere un AR (1) MA (1)

Propuestas

fit1<-Arima(STDIFF, order = c(0,0,3), seasonal = c(0,0,1))
fit1
## Series: STDIFF 
## ARIMA(0,0,3)(0,0,1)[12] with non-zero mean 
## 
## Coefficients:
##          ma1     ma2     ma3     sma1     mean
##       0.1379  0.2183  0.1938  -0.5757  -2.2111
## s.e.  0.0722  0.0702  0.0769   0.1011   2.4779
## 
## sigma^2 estimated as 2285:  log likelihood=-1014.9
## AIC=2041.81   AICc=2042.26   BIC=2061.35

Criterio de Akaike: 2041.81

fit2<-Arima(STDIFF, order= c(1,0,0), seasonal = c(1,0,0))
fit2
## Series: STDIFF 
## ARIMA(1,0,0)(1,0,0)[12] with non-zero mean 
## 
## Coefficients:
##          ar1     sar1     mean
##       0.2591  -0.3081  -3.4668
## s.e.  0.0696   0.0683   3.8211
## 
## sigma^2 estimated as 2606:  log likelihood=-1026.64
## AIC=2061.27   AICc=2061.49   BIC=2074.3

Criterio de Akaike: 2061.27, es menos efectivo que la especificación anterior.

fit3<-Arima(STDIFF, order = c(1,0,1), seasonal = c(1,0,1))
fit3
## Series: STDIFF 
## ARIMA(1,0,1)(1,0,1)[12] with non-zero mean 
## 
## Coefficients:
##          ar1      ma1    sar1     sma1     mean
##       0.8486  -0.6584  0.3448  -0.8494  -1.5466
## s.e.  0.0755   0.1021  0.1132   0.0929   2.5025
## 
## sigma^2 estimated as 2139:  log likelihood=-1010.31
## AIC=2032.62   AICc=2033.07   BIC=2052.16

Criterio de Akaike: 2032.62, Ha sido el mejor de los tres modelos

8. Como parte del diagnóstico del modelo, analizar los residuos y aplicar la prueba Ljung-Box.

checkresiduals(fit1)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(0,0,3)(0,0,1)[12] with non-zero mean
## Q* = 36.151, df = 19, p-value = 0.01011
## 
## Model df: 5.   Total lags used: 24

La anterior especificación no pasa los tests básicos por lo que nos vamos a la opción dos.

checkresiduals(fit2)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,0,0)(1,0,0)[12] with non-zero mean
## Q* = 53.761, df = 21, p-value = 0.0001069
## 
## Model df: 3.   Total lags used: 24

Esta especificación tampoco cumple con los requisitos necesarios. Por lo que vamos con la siguiente opción.

checkresiduals(fit3)

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

Esta última propuesta cumple con todos los requisitos, el valor p es de: 0.076 por lo que se acepta la Hipótesis nula= No correlación serial, además podemos notar que los rezagos se comportan de manera estable, solo se sale un rezago, de la zona de aceptación, por muy poco.

9. Usar la función auto.arima() y compare sus resultados con el inciso anterior.

Propuestaauto<-auto.arima(STDIFF)
Propuestaauto
## Series: STDIFF 
## ARIMA(1,0,1)(0,0,2)[12] with zero mean 
## 
## Coefficients:
##          ar1      ma1     sma1     sma2
##       0.8520  -0.6614  -0.5175  -0.1918
## s.e.  0.0755   0.1029   0.0738   0.0757
## 
## sigma^2 estimated as 2161:  log likelihood=-1011.18
## AIC=2032.36   AICc=2032.68   BIC=2048.65

Utilizamos dos funciones de auto.arima. La primera nos presentó el siguiente modelo. Nos sugiere un modelo ARIMA (1, 0, 1) (0, 0, 2) Podemos observar que el criterio de Akaike es de 2032.36

checkresiduals(Propuestaauto)

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,0,1)(0,0,2)[12] with zero mean
## Q* = 29.815, df = 20, p-value = 0.07291
## 
## Model df: 4.   Total lags used: 24

Analizando, notamos que solo se sale un rezago y que no hay correlación serial, ya que el valor “p” es de: 0.072, es ligeramente peor que el modelo que propusimos. Por este motivo usamos otro auto.arima pero con más especificaciones:

propuesta2<-auto.arima(STDIFF, stepwise = FALSE, approximation = FALSE)
propuesta2
## Series: STDIFF 
## ARIMA(1,0,1)(1,0,1)[12] with zero mean 
## 
## Coefficients:
##          ar1      ma1    sar1     sma1
##       0.8541  -0.6613  0.3433  -0.8468
## s.e.  0.0731   0.1000  0.1127   0.0920
## 
## sigma^2 estimated as 2133:  log likelihood=-1010.5
## AIC=2030.99   AICc=2031.31   BIC=2047.28

El modelo que propone es el mismo que nuestra tercera opción: un modelo ARIMA (1, 0, 1) (1, 0, 1) Criterio de Akaike: 2030.99

checkresiduals(propuesta2) 

## 
##  Ljung-Box test
## 
## data:  Residuals from ARIMA(1,0,1)(1,0,1)[12] with zero mean
## Q* = 28.24, df = 20, p-value = 0.1038
## 
## Model df: 4.   Total lags used: 24

Este modelo cumple con todos los requisitos, solo que el valor “p” cambia respecto a nuestro modelo, éste es de 0.1038.

10. Presente la ecuación final y describa.

ARIMA (1, 0, 1) (1, 0, 1) [12] Nuestro modelo presenta un pico en el primer rezago tanto en la ACF como en la PACF, asimismo muestra un pico significativo en el rezago 12 en ambas funciones de autocorrelación.

11. Realizar el pronóstico ARIMA para dos años.

fit5<-Arima(divorcios,order = c(1,1,0), seasonal = c(1,0,0))
fit5
## Series: divorcios 
## ARIMA(1,1,0)(1,0,0)[12] 
## 
## Coefficients:
##           ar1    sar1
##       -0.5341  0.3835
## s.e.   0.0602  0.0656
## 
## sigma^2 estimated as 2403:  log likelihood=-1078.27
## AIC=2162.54   AICc=2162.66   BIC=2172.48
autoplot(forecast(fit5))

12. Utilizando metodos de pronósticos simples y/o suavizamiento exponencial generar un pronóstico para dos años. Elegir el método que presente la raíz del error medio al cuadrado más pequeño.

Para poder realizar el pronostico de nuestra serie, se realizo un recorte de los ultimos dos años para pronosticar nuestra serie mediante los metodos mean(media), naive y naive estacional.

DIVPRO<-window(divorcios,start=2000,end=c(2014,12))

autoplot(DIVPRO) 

library(forecast)
library(ggplot2)
DIVPRO<-window(divorcios,start=2000,end=c(2014,12))
View(DIVPRO)
divorcio1 <- meanf(DIVPRO,h=24)
divorcio2 <- rwf(DIVPRO,h=24)
divorcio3 <- snaive(DIVPRO,h=24)
divorcio4 <- naive(DIVPRO, h=24)
plot(divorcio1)

plot(divorcio2)

plot(divorcio3)

plot(divorcio4)

Una vez comparado con los dos años siguientes de la serie, concluimos que el metodo de naive estacional presenta un comportamiento mas logico y razonable para los siguientes años. Cabe mencionar que se tuvo que hacer de esta manera porque hubo problemas con el codigo previamente empleado: autoplot(DIVPRO) + forecast::autolayer(meanf(DIVPRO, h=24), PI=FALSE, series=“Mean”) + forecast::autolayer(naive(DIVPRO, h=24), PI=FALSE, series=“Naïve”) + forecast::autolayer(snaive(DIVPRO, h=24), PI=FALSE, series=“Seasonal naïve”) + ggtitle(“Pronostico del N[umero de Embarazos en la CDMX”) + xlab(“Año”) + ylab(“Embarazos”) + guides(colour=guide_legend(title=“Pronostico”))… El error, fue el siguiente. Error: Invalid input: date_trans works with objects of class Date only

14Analizar un VAR, necesita incluir los puntos siguientes:

  1. Buscar otra variable (incluso más variables) y redactar la teoría que se encuentra detrás de la relaci ´on entre variables que propone.
  2. Evaluar la estacionalidad de la nueva variable.
  3. Proponga el orden de rezagos óptimo para el VAR y presentar la estimación.
  4. Probar la correlación serial y normalidad en los residuos.
  5. Resultados de la prueba de causalidad de Granger, y f ) Crear un pron´ostico para dos años usando el VAR especificado.
  6. Interpretar las funciones impulso-respuesta.

a) La variable que vamos a utilizar es la de “Desempleos” ya que consideramos que una de las causas por las que hay divorcios, es que un miembro de la pareja está desempleado(a) y esto provoca una inestabilidad financiera que conlleva al divorcio. Para realizar este ejercicio, recortaremos nuestra serie para que sea mas eficiente la estimación.

b)

library(fpp2)
## Loading required package: fma
## Loading required package: expsmooth
library(vars)
## Loading required package: MASS
## 
## Attaching package: 'MASS'
## The following objects are masked from 'package:fma':
## 
##     cement, housing, petrol
## Loading required package: strucchange
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
## Loading required package: sandwich
## Loading required package: lmtest
library(foreign)
library(zoo)
desempleo<- read.csv("C:/Users/user/Documents/economia/Sexto semestre/Series de tiempo/desempleo.csv")
plot.ts(desempleo)

ts.des<-ts(desempleo$desempleos, start = c(2005,1), frequency = 12)
plot(ts.des)

dts.des<-diff(diff(ts.des, 12))
autoplot(dts.des)

library(urca)
summary(ur.df(dts.des, type = "none",selectlags = "AIC"))
## 
## ############################################### 
## # 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 
## -395471 -141717  -45519  122387  405908 
## 
## Coefficients:
##            Estimate Std. Error t value Pr(>|t|)    
## z.lag.1    -0.98960    0.21600  -4.581 3.94e-05 ***
## z.diff.lag  0.01565    0.16067   0.097    0.923    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 198800 on 43 degrees of freedom
## Multiple R-squared:  0.4729, Adjusted R-squared:  0.4484 
## F-statistic: 19.29 on 2 and 43 DF,  p-value: 1.049e-06
## 
## 
## Value of test-statistic is: -4.5815 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.62 -1.95 -1.61
Divorciovar<-window(divorcios,start=2005,end=c(2009,12))
autoplot(Divorciovar)

diffdivorcio<- diff(Divorciovar, 12)
summary(ur.df(diffdivorcio, type = "none",selectlags = "AIC"))
## 
## ############################################### 
## # 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 
## -135.364   -9.815   11.295   31.358   55.849 
## 
## Coefficients:
##            Estimate Std. Error t value Pr(>|t|)    
## z.lag.1     -0.8876     0.2136  -4.156 0.000147 ***
## z.diff.lag  -0.1371     0.1490  -0.920 0.362733    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 35.13 on 44 degrees of freedom
## Multiple R-squared:  0.5263, Adjusted R-squared:  0.5047 
## F-statistic: 24.44 on 2 and 44 DF,  p-value: 7.272e-08
## 
## 
## Value of test-statistic is: -4.1565 
## 
## Critical values for test statistics: 
##       1pct  5pct 10pct
## tau1 -2.62 -1.95 -1.61
basecombinada<- cbind.zoo(des= dts.des, div=diffdivorcio)

15. Evaluar la potencial cointegraci´on usando la metodolog´ia de Engle-Granger y Johansen.