#Bibliotecas
#-----------------------
library(vars)
## Warning: package 'vars' was built under R version 4.3.3
## Loading required package: MASS
## Loading required package: strucchange
## Warning: package 'strucchange' was built under R version 4.3.3
## 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: urca
## Loading required package: lmtest
library(tseries)
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(forecast)
library(urca)
library(highcharter)
## Warning: package 'highcharter' was built under R version 4.3.3
##
## Attaching package: 'highcharter'
## The following object is masked from 'package:lmtest':
##
## unemployment
# Datos
#------------------------------------
data(Canada)
Canada = as.data.frame(Canada)
Canada
## e prod rw U
## 1 929.6105 405.3665 386.1361 7.53
## 2 929.8040 404.6398 388.1358 7.70
## 3 930.3184 403.8149 390.5401 7.47
## 4 931.4277 404.2158 393.9638 7.27
## 5 932.6620 405.0467 396.7647 7.37
## 6 933.5509 404.4167 400.0217 7.13
## 7 933.5315 402.8191 400.7515 7.40
## 8 933.0769 401.9773 405.7335 8.33
## 9 932.1238 402.0897 409.0504 8.83
## 10 930.6359 401.3067 411.3984 10.43
## 11 929.0971 401.6302 413.0194 12.20
## 12 928.5633 401.5638 415.1670 12.77
## 13 929.0694 402.8157 414.6621 12.43
## 14 930.2655 403.1421 415.7319 12.23
## 15 931.6770 403.0786 416.2315 11.70
## 16 932.1390 403.7188 418.1439 11.20
## 17 932.2767 404.8668 419.7352 11.27
## 18 932.8328 405.6362 420.4842 11.47
## 19 933.7334 405.1363 420.9309 11.30
## 20 934.1772 406.0246 422.1124 11.17
## 21 934.5928 406.4123 423.6278 11.00
## 22 935.6067 406.3009 423.9887 10.63
## 23 936.5111 406.3354 424.1902 10.27
## 24 937.4201 406.7737 426.1270 10.20
## 25 938.4159 405.1525 426.8578 9.67
## 26 938.9992 404.9298 426.7457 9.60
## 27 939.2354 404.5765 426.8858 9.60
## 28 939.6795 404.1995 428.8403 9.50
## 29 940.2497 405.9499 430.1223 9.50
## 30 941.4358 405.8221 430.2307 9.03
## 31 942.2981 406.4463 430.3930 8.70
## 32 943.5322 407.0512 432.0284 8.13
## 33 944.3490 407.9460 433.3886 7.87
## 34 944.8215 408.1796 433.9641 7.67
## 35 945.0671 408.5998 434.4844 7.80
## 36 945.8067 409.0906 436.1569 7.73
## 37 946.8697 408.7042 438.2651 7.57
## 38 946.8766 408.9803 438.7636 7.57
## 39 947.2497 408.3287 439.9498 7.33
## 40 947.6513 407.8857 441.8359 7.57
## 41 948.1840 407.2605 443.1769 7.63
## 42 948.3492 406.7752 444.3592 7.60
## 43 948.0322 406.1794 444.5236 8.17
## 44 947.1065 405.4398 446.9694 9.20
## 45 946.0796 403.2800 450.1586 10.17
## 46 946.1838 403.3649 451.5464 10.33
## 47 946.2258 403.3807 452.2984 10.40
## 48 945.9978 404.0032 453.1201 10.37
## 49 945.5183 404.4774 453.9991 10.60
## 50 945.3514 404.7868 454.9552 11.00
## 51 945.2918 405.2710 455.4824 11.40
## 52 945.4008 405.3830 456.1009 11.73
## 53 945.9058 405.1564 457.2027 11.07
## 54 945.9035 406.4700 457.3886 11.67
## 55 946.3190 406.2293 457.7799 11.47
## 56 946.5796 406.7265 457.5535 11.30
## 57 946.7800 408.5785 458.8024 10.97
## 58 947.6283 409.6767 459.0564 10.63
## 59 948.6221 410.3858 459.1578 10.10
## 60 949.3992 410.5395 459.7037 9.67
## 61 949.9481 410.4453 459.7037 9.53
## 62 949.7945 410.6256 460.0258 9.47
## 63 949.9534 410.8672 461.0257 9.50
## 64 950.2502 411.2359 461.3039 9.27
## 65 950.5380 410.6637 461.4031 9.50
## 66 950.7871 410.8085 462.9277 9.43
## 67 950.8695 412.1160 464.6888 9.70
## 68 950.9281 412.9994 465.0717 9.90
## 69 951.8457 412.9551 464.2851 9.43
## 70 952.6005 412.8241 464.0344 9.30
## 71 953.5976 413.0489 463.4535 8.87
## 72 954.1434 413.6110 465.0717 8.77
## 73 954.5426 413.6048 466.0889 8.60
## 74 955.2631 412.9684 466.6171 8.33
## 75 956.0561 412.2659 465.7478 8.17
## 76 956.7966 412.9106 465.8995 8.03
## 77 957.3865 413.8294 466.4099 7.90
## 78 958.0634 414.2242 466.9552 7.87
## 79 958.7166 415.1678 467.6281 7.53
## 80 959.4881 415.7016 467.7026 6.93
## 81 960.3625 416.8674 469.1348 6.80
## 82 960.7834 417.6104 469.3364 6.70
## 83 961.0290 418.0030 470.0117 6.93
## 84 961.7657 417.2667 469.6472 6.87
# Plot de las series temporales
#--------------------------------
layout(matrix(1:4, nrow = 2, ncol = 2))
plot.ts(Canada$e, main = "Employment", ylab ="", xlab="")
plot.ts(Canada$prod, main = "Productivity", ylab ="", xlab="")
plot.ts(Canada$rw, main = "Real Wage", ylab ="", xlab="")
plot.ts(Canada$U, main = "Unemployment Rate", ylab ="", xlab="")
var_canada <- ts(Canada[,c(1,2,3,4)],frequency = 12)
plot.ts(var_canada)
Estacionariedad
#Empleo
#----------------------------
pp.test(var_canada[,1])
##
## Phillips-Perron Unit Root Test
##
## data: var_canada[, 1]
## Dickey-Fuller Z(alpha) = -5.8701, Truncation lag parameter = 3, p-value
## = 0.7735
## alternative hypothesis: stationary
P valor = 0.77 > 0.05 –> Acepto H0 –> Raíz Unitaria (No estacionaria)
#Productividad
#-----------------------------
pp.test(var_canada[,2])
##
## Phillips-Perron Unit Root Test
##
## data: var_canada[, 2]
## Dickey-Fuller Z(alpha) = -7.4247, Truncation lag parameter = 3, p-value
## = 0.6816
## alternative hypothesis: stationary
P valor = 0.68 > 0.05 –> Acepto H0 –> Raíz Unitaria (No estacionaria)
# Salario real
#---------------------
pp.test(var_canada[,3])
##
## Phillips-Perron Unit Root Test
##
## data: var_canada[, 3]
## Dickey-Fuller Z(alpha) = -4.6296, Truncation lag parameter = 3, p-value
## = 0.8468
## alternative hypothesis: stationary
P valor = 0.84 > 0.05 –> Acepto H0 –> Raíz Unitaria (No estacionaria)
# Tasa de desempleo
#---------------------
pp.test(var_canada[,4])
##
## Phillips-Perron Unit Root Test
##
## data: var_canada[, 4]
## Dickey-Fuller Z(alpha) = -7.1861, Truncation lag parameter = 3, p-value
## = 0.6957
## alternative hypothesis: stationary
P valor = 0.69 > 0.05 –> Acepto H0 –> Raíz Unitaria (No estacionaria)
# ADF empleo
#-------------------------------
adf.test(var_canada[,1])
##
## Augmented Dickey-Fuller Test
##
## data: var_canada[, 1]
## Dickey-Fuller = -2.148, Lag order = 4, p-value = 0.5152
## alternative hypothesis: stationary
P valor = 0.51 > 0.05 –> Acepto H0 –> Raíz Unitaria (No estacionaria)
# ADF productividad
#---------------------------------------
adf.test(var_canada[,2])
##
## Augmented Dickey-Fuller Test
##
## data: var_canada[, 2]
## Dickey-Fuller = -2.8725, Lag order = 4, p-value = 0.218
## alternative hypothesis: stationary
P valor = 0.218 > 0.05 –> Acepto H0 –> Raíz Unitaria (No estacionaria)
# ADF salario real
#---------------------------------------
adf.test(var_canada[,3])
##
## Augmented Dickey-Fuller Test
##
## data: var_canada[, 3]
## Dickey-Fuller = -2.0558, Lag order = 4, p-value = 0.553
## alternative hypothesis: stationary
P valor = 0.553 > 0.05 –> Acepto H0 –> Raíz Unitaria (No estacionaria)
# ADF tasa de desempleo
#---------------------------------------
adf.test(var_canada[,4])
##
## Augmented Dickey-Fuller Test
##
## data: var_canada[, 4]
## Dickey-Fuller = -2.5988, Lag order = 4, p-value = 0.3303
## alternative hypothesis: stationary
P valor = 0.3303 > 0.05 –> Acepto H0 –> Raíz Unitaria (No estacionaria)
Las series temporales en niveles no son estaiconarias y debemos estudiar su primera diferencia.
# PP primera diferencia del empleo
#---------------------------------------
pp.test(diff(var_canada[,1]))
##
## Phillips-Perron Unit Root Test
##
## data: diff(var_canada[, 1])
## Dickey-Fuller Z(alpha) = -26.227, Truncation lag parameter = 3, p-value
## = 0.01233
## alternative hypothesis: stationary
P valor = 0.01233 > 0.05 –> Rechazo H0 –> NO Raíz Unitaria (Sí estacionaria)
# PP primera diferencia de la productividad
#---------------------------------------
pp.test(diff(var_canada[,2]))
## Warning in pp.test(diff(var_canada[, 2])): p-value smaller than printed p-value
##
## Phillips-Perron Unit Root Test
##
## data: diff(var_canada[, 2])
## Dickey-Fuller Z(alpha) = -60.867, Truncation lag parameter = 3, p-value
## = 0.01
## alternative hypothesis: stationary
P valor = 0.01 > 0.05 –> Rechazo H0 –> NO Raíz Unitaria (Sí estacionaria)
# PP primera diferencia del salario real
#---------------------------------------
pp.test(diff(var_canada[,3]))
## Warning in pp.test(diff(var_canada[, 3])): p-value smaller than printed p-value
##
## Phillips-Perron Unit Root Test
##
## data: diff(var_canada[, 3])
## Dickey-Fuller Z(alpha) = -61.576, Truncation lag parameter = 3, p-value
## = 0.01
## alternative hypothesis: stationary
P valor = 0.01 > 0.05 –> Rechazo H0 –> NO Raíz Unitaria (Sí estacionaria)
# PP primera diferencia de la tasa de desempleo
#---------------------------------------
pp.test(diff(var_canada[,4]))
## Warning in pp.test(diff(var_canada[, 4])): p-value smaller than printed p-value
##
## Phillips-Perron Unit Root Test
##
## data: diff(var_canada[, 4])
## Dickey-Fuller Z(alpha) = -38.946, Truncation lag parameter = 3, p-value
## = 0.01
## alternative hypothesis: stationary
P valor = 0.01 > 0.05 –> Rechazo H0 –> NO Raíz Unitaria (Sí estacionaria)
# ADF de la primera diferencia del empleo
#--------------------------------------------
adf.test(diff(var_canada[,1]))
##
## Augmented Dickey-Fuller Test
##
## data: diff(var_canada[, 1])
## Dickey-Fuller = -3.2668, Lag order = 4, p-value = 0.08274
## alternative hypothesis: stationary
Al 10% sí es estacionaria, pero no al 5%.
# ADF de la primera diferencia de la productividad
#--------------------------------------------
adf.test(diff(var_canada[,2]))
##
## Augmented Dickey-Fuller Test
##
## data: diff(var_canada[, 2])
## Dickey-Fuller = -3.2361, Lag order = 4, p-value = 0.08775
## alternative hypothesis: stationary
Al 10% sí es estacionaria, pero no al 5%.
# ADF de la primera diferencia del salario real
#--------------------------------------------
adf.test(diff(var_canada[,3]))
##
## Augmented Dickey-Fuller Test
##
## data: diff(var_canada[, 3])
## Dickey-Fuller = -3.2352, Lag order = 4, p-value = 0.08789
## alternative hypothesis: stationary
Al 10% sí es estacionaria, pero no al 5%.
# ADF de la primera diferencia de la tasa de desempleo
#--------------------------------------------
adf.test(diff(var_canada[,4]))
##
## Augmented Dickey-Fuller Test
##
## data: diff(var_canada[, 4])
## Dickey-Fuller = -3.7325, Lag order = 4, p-value = 0.02698
## alternative hypothesis: stationary
Al 10% y al 5% sí es estacionaria.
Con PP: son estacionarias todas las series al 5% Con ADF: Todas las series son estacionarias al 10%.
Diferencia de las variables
# Primera diferencia de todas las variables
#-------------------------------------
var_canada_dif <- diff(var_canada[,c(1,2,3,4)])
head(var_canada_dif)
## e prod rw U
## Feb 1 0.19347066 -0.7266325 1.9996500 0.17
## Mar 1 0.51440302 -0.8249509 2.4043538 -0.23
## Apr 1 1.10929985 0.4008901 3.4237043 -0.20
## May 1 1.23431817 0.8309404 2.8008737 0.10
## Jun 1 0.88893413 -0.6299749 3.2570107 -0.24
## Jul 1 -0.01941353 -1.5976113 0.7297971 0.27
#Gráfico de las primeras diferencias
#-----------------------------------------
plot.ts(var_canada_dif)
hchart(var_canada_dif)
## Warning: Deprecated function. Use the `create_axis` function.
Una vez las series son estacionarias, se procede a la estimación
# Selección del orden del VAR
#------------------------------------
VARselect(var_canada_dif, lag.max=8, type = "const")
## $selection
## AIC(n) HQ(n) SC(n) FPE(n)
## 2 1 1 2
##
## $criteria
## 1 2 3 4 5
## AIC(n) -6.270032208 -6.285458074 -6.03818873 -6.035665365 -5.808725091
## HQ(n) -6.023272860 -5.841291247 -5.39661443 -5.196683581 -4.772335828
## SC(n) -5.652035378 -5.173063780 -4.43139698 -3.934476142 -3.213138403
## FPE(n) 0.001893667 0.001871884 0.00241986 0.002469804 0.003191486
## 6 7 8
## AIC(n) -5.558322503 -5.370524857 -5.33131938
## HQ(n) -4.324525761 -3.939320637 -3.70270768
## SC(n) -2.468338351 -1.786143241 -1.25254030
## FPE(n) 0.004286004 0.005511821 0.00626064
$criteria 2
AIC(n) -6.285458074 HQ(n) -5.841291247 SC(n) -5.173063780 FPE(n)
0.001871884
Usaremos dos retardos (que son los que minimizan el AIC)
#Modelo VAR p = 2
#----------------------------
modelo_var <- VAR(Canada, p=2)
summary(modelo_var)
##
## VAR Estimation Results:
## =========================
## Endogenous variables: e, prod, rw, U
## Deterministic variables: const
## Sample size: 82
## Log Likelihood: -175.819
## Roots of the characteristic polynomial:
## 0.995 0.9081 0.9081 0.7381 0.7381 0.1856 0.1429 0.1429
## Call:
## VAR(y = Canada, p = 2)
##
##
## Estimation results for equation e:
## ==================================
## e = e.l1 + prod.l1 + rw.l1 + U.l1 + e.l2 + prod.l2 + rw.l2 + U.l2 + const
##
## Estimate Std. Error t value Pr(>|t|)
## e.l1 1.638e+00 1.500e-01 10.918 < 2e-16 ***
## prod.l1 1.673e-01 6.114e-02 2.736 0.00780 **
## rw.l1 -6.312e-02 5.524e-02 -1.143 0.25692
## U.l1 2.656e-01 2.028e-01 1.310 0.19444
## e.l2 -4.971e-01 1.595e-01 -3.116 0.00262 **
## prod.l2 -1.017e-01 6.607e-02 -1.539 0.12824
## rw.l2 3.844e-03 5.552e-02 0.069 0.94499
## U.l2 1.327e-01 2.073e-01 0.640 0.52418
## const -1.370e+02 5.585e+01 -2.453 0.01655 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.3628 on 73 degrees of freedom
## Multiple R-Squared: 0.9985, Adjusted R-squared: 0.9984
## F-statistic: 6189 on 8 and 73 DF, p-value: < 2.2e-16
##
##
## Estimation results for equation prod:
## =====================================
## prod = e.l1 + prod.l1 + rw.l1 + U.l1 + e.l2 + prod.l2 + rw.l2 + U.l2 + const
##
## Estimate Std. Error t value Pr(>|t|)
## e.l1 -0.17277 0.26977 -0.640 0.52390
## prod.l1 1.15043 0.10995 10.464 3.57e-16 ***
## rw.l1 0.05130 0.09934 0.516 0.60710
## U.l1 -0.47850 0.36470 -1.312 0.19362
## e.l2 0.38526 0.28688 1.343 0.18346
## prod.l2 -0.17241 0.11881 -1.451 0.15104
## rw.l2 -0.11885 0.09985 -1.190 0.23778
## U.l2 1.01592 0.37285 2.725 0.00805 **
## const -166.77552 100.43388 -1.661 0.10109
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.6525 on 73 degrees of freedom
## Multiple R-Squared: 0.9787, Adjusted R-squared: 0.9764
## F-statistic: 419.3 on 8 and 73 DF, p-value: < 2.2e-16
##
##
## Estimation results for equation rw:
## ===================================
## rw = e.l1 + prod.l1 + rw.l1 + U.l1 + e.l2 + prod.l2 + rw.l2 + U.l2 + const
##
## Estimate Std. Error t value Pr(>|t|)
## e.l1 -0.268833 0.322619 -0.833 0.407
## prod.l1 -0.081065 0.131487 -0.617 0.539
## rw.l1 0.895478 0.118800 7.538 1.04e-10 ***
## U.l1 0.012130 0.436149 0.028 0.978
## e.l2 0.367849 0.343087 1.072 0.287
## prod.l2 -0.005181 0.142093 -0.036 0.971
## rw.l2 0.052677 0.119410 0.441 0.660
## U.l2 -0.127708 0.445892 -0.286 0.775
## const -33.188339 120.110525 -0.276 0.783
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.7803 on 73 degrees of freedom
## Multiple R-Squared: 0.9989, Adjusted R-squared: 0.9987
## F-statistic: 8009 on 8 and 73 DF, p-value: < 2.2e-16
##
##
## Estimation results for equation U:
## ==================================
## U = e.l1 + prod.l1 + rw.l1 + U.l1 + e.l2 + prod.l2 + rw.l2 + U.l2 + const
##
## Estimate Std. Error t value Pr(>|t|)
## e.l1 -0.58076 0.11563 -5.023 3.49e-06 ***
## prod.l1 -0.07812 0.04713 -1.658 0.101682
## rw.l1 0.01866 0.04258 0.438 0.662463
## U.l1 0.61893 0.15632 3.959 0.000173 ***
## e.l2 0.40982 0.12296 3.333 0.001352 **
## prod.l2 0.05212 0.05093 1.023 0.309513
## rw.l2 0.04180 0.04280 0.977 0.331928
## U.l2 -0.07117 0.15981 -0.445 0.657395
## const 149.78056 43.04810 3.479 0.000851 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
##
## Residual standard error: 0.2797 on 73 degrees of freedom
## Multiple R-Squared: 0.9726, Adjusted R-squared: 0.9696
## F-statistic: 324 on 8 and 73 DF, p-value: < 2.2e-16
##
##
##
## Covariance matrix of residuals:
## e prod rw U
## e 0.131635 -0.007469 -0.04210 -0.06909
## prod -0.007469 0.425711 0.06461 0.01392
## rw -0.042099 0.064613 0.60886 0.03422
## U -0.069087 0.013923 0.03422 0.07821
##
## Correlation matrix of residuals:
## e prod rw U
## e 1.00000 -0.03155 -0.1487 -0.6809
## prod -0.03155 1.00000 0.1269 0.0763
## rw -0.14870 0.12691 1.0000 0.1568
## U -0.68090 0.07630 0.1568 1.0000
Diagnosis
roots(modelo_var)
## [1] 0.9950338 0.9081062 0.9081062 0.7380565 0.7380565 0.1856381 0.1428889
## [8] 0.1428889
Al ser inferiores a uno todas las raíces del polinomio, el VAR es estacionario.
Contraste de autocorrelación de los errores
serial.test(modelo_var,lags.pt=16, type ="PT.adjusted")
##
## Portmanteau Test (adjusted)
##
## data: Residuals of VAR object modelo_var
## Chi-squared = 231.59, df = 224, p-value = 0.3497
El p valor es 0.3497 > 0.05. Aceptamos H0, es decir, que los errores no están correlacionados.
Prueba de Normaliad
normality.test(modelo_var, multivariate.only = FALSE)
## $e
##
## JB-Test (univariate)
##
## data: Residual of e equation
## Chi-squared = 0.15347, df = 2, p-value = 0.9261
##
##
## $prod
##
## JB-Test (univariate)
##
## data: Residual of prod equation
## Chi-squared = 4.2651, df = 2, p-value = 0.1185
##
##
## $rw
##
## JB-Test (univariate)
##
## data: Residual of rw equation
## Chi-squared = 0.33482, df = 2, p-value = 0.8459
##
##
## $U
##
## JB-Test (univariate)
##
## data: Residual of U equation
## Chi-squared = 0.56643, df = 2, p-value = 0.7534
##
##
## $JB
##
## JB-Test (multivariate)
##
## data: Residuals of VAR object modelo_var
## Chi-squared = 5.094, df = 8, p-value = 0.7475
##
##
## $Skewness
##
## Skewness only (multivariate)
##
## data: Residuals of VAR object modelo_var
## Chi-squared = 1.7761, df = 4, p-value = 0.7769
##
##
## $Kurtosis
##
## Kurtosis only (multivariate)
##
## data: Residuals of VAR object modelo_var
## Chi-squared = 3.3179, df = 4, p-value = 0.5061
Todos los errores son normales, al tener un p valor mayor a un 5% en todos los casos.