Política monetaria en Costa Rica: efectos y transmisión en el periodo 2011-2022¶

Preliminar¶

Importación de paquetes¶

In [ ]:
# Paquete bccr, del profesor Randall Romero (https://randall-romero.github.io/bccr/)
try:
    import bccr
    from bccr import SW
except:
    print('Paquete bccr no instalado, se procede con la instalación')
    !pip install bccr
    from bccr import SW


# Paquete statsmodels, para modelos econométricos y pruebas (https://www.statsmodels.org/stable/index.html)
## General
try:
    import statsmodels.api as sm
except:
    print('Paquete statsmodels no instalado, se procede con la instalación')
    !pip install statsmodels
    import statsmodels.api as sm
## Desestacionalización
import statsmodels.tsa.seasonal
## Pruebas Econométricas
from statsmodels.tsa.stattools import adfuller, kpss
from statsmodels.stats.stattools import durbin_watson
from statsmodels.tsa.stattools import zivot_andrews
from scipy.stats import normaltest, jarque_bera
## VAR
from statsmodels.tsa.vector_ar.vecm import *



# Paquetes estándar
## Análisis de Datos
import pandas as pd
from functools import reduce
import numpy as np
import datetime
## Graficación
import matplotlib.pyplot as plt
import matplotlib.pylab as pylab
import seaborn as sns

### Valores generales para graficación
plt.rcParams["figure.dpi"] = 300

# Diagnóstico (no se usan todos, REVISAR...)
from statsmodels.stats.diagnostic import acorr_breusch_godfrey, acorr_ljungbox, acorr_lm, het_breuschpagan

#Shortcuts
g=globals()

Entorno de carpetas y graficación¶

In [ ]:
!mkdir Datos
!mkdir PDF
!mkdir PNG

mkdir: cannot create directory ‘Datos’: File exists
mkdir: cannot create directory ‘PDF’: File exists
mkdir: cannot create directory ‘PNG’: File exists
In [ ]:
plt.style.use('seaborn-v0_8')
params = {'figure.figsize':(5,3),
        'legend.frameon':True,
        'axes.labelsize': '11',
        'legend.fontsize':'9',
        'xtick.labelsize':'11',
        'ytick.labelsize':'11',
        'savefig.pad_inches':'0.1',
        'savefig.bbox': 'tight',
        'axes.edgecolor': 'w',
        'legend.fancybox':True,
        'legend.frameon':True}

pylab.rcParams.update(params)

Conformación de la base de datos¶

NOTA: EL IMAE COMO ES DESESTACIONALIZADO, CUANDO SE AGREGAN OBSERVACIONES NUEVAS CAMBIAN LAS VIEJAS. Por lo tanto, si se vuelve a importar las series cambian ligeramente los resultados (t-estadísticos y p-values en las pruebas del IMAE y del Modelo Conjunto). Para esto al final incluimos la base con la que corrimos todo y la usamos para correr el modelo. Aún así, incluímos el código que empleamos para formar la base originalmente.

Esta base es ligeramente distinta a la que usamos para la Entrega 2, pues la guardamos más tarde. Esto significa que los p-values y t-estadísticos de la entrega final *son distintos a los de la Entrega 2. Sin embargo, en ningún caso afecta la conclusión (significancia o no significancia) pues el cambio es casi despreciable.

Hicimos todos los resultados del documento final consistentes con esta base y con la versión actual de este código, para que sean verificables y replicables.

Datos del Banco Central de Costa Rica¶

Se empiezan importando aquellas series del BCCR que no requieren ningún formato adicional. Estas son todas las series que el BCCR reporta de manera mensual.

In [ ]:
#Mensuales
base_bccr_temp = SW(imae=87853, ipc=89635, cred=1415, epi = 84990,
  FechaInicio='2011-07-01', FechaFinal='2022-12-31')

# Se especifica el índice como un formato de fecha
base_bccr_temp.index = pd.PeriodIndex(base_bccr_temp.index)

Ahora se empiezan a importar aquellas series que sí requieren de un formato adicional. Estas son: la Tasa de Política Monetaria (TPM) que se reporta en forma diaria, y el tipo de cambio nominal que se reporta también de forma diaria. En el caso de la TPM, lo que nos interesa es ver los cambios mes a mes, por ende tomamos la observación de final de mes.

In [ ]:
# Se importa la base en frecuencia mensual, especificando que se tome la última observación de cada mes con func = 'last'
base_bccr_temp['R'] = SW(3541,
  FechaInicio='2011-07-01', FechaFinal='2022-12-31',
  func='last', freq='M')

En el caso del tipo de cambio nominal, tomamos el promedio del tipo de cambio transado en MONEX en ese mes; el promedio se calcula eliminando aquellas observaciones en el que el tipo de cambio transado en MONEX es 0, pues eso no representa un valor del tipo de cambio sino que un día en el que no se transó en el mercado.

In [ ]:
# Se importa la base en frecuencia diaria
tipo_de_cambio = SW(E=3323, FechaInicio='2011-07-01', FechaFinal='2022-12-31')

# Se eliminan las observaciones de días en los que no se transó en MONEX
tipo_de_cambio = tipo_de_cambio.loc[tipo_de_cambio['E'] != 0.00]

# Se calcula el promedio mensual
tipo_de_cambio = tipo_de_cambio.groupby(pd.Grouper(freq = 'M')).mean()

# Se une la base con base_bccr_temp
base_bccr_temp = base_bccr_temp.merge(tipo_de_cambio, how = 'inner', left_index = True, right_index = True)

Adicionalmente, se crea un índice de expectativas de precios a 12 meses. Esto se hace considerando el valor del índice de precios al consumidor en $t$ y aplicándole la inflación esperada a 12 meses. Por ejemplo, mediante el IPC de Enero 2014 y la Inflación Esperada a 12 meses de Enero 2014 se genera un índice de precios que los consumidores esperan observar doce meses después de la observación. Para esto se parte de la fórmula de inflación esperada:

$$ E(\pi_{t+12})_t = \frac{E(IPC_{t+12})_t - IPC_{t}}{IPC_{t}}$$

Y se despeja el valor de $E(IPC_{t+12})$ tal que:

$$E(IPC_{t+12})_t = IPC_{t} (1 + E(\pi_{t+12})_t)$$

In [ ]:
# Se importan datos para crear un índice de Expectativas de Inflación
expectativas_inflacion = SW(epi = 84990, ipc = 89635, FechaInicio = '2011-07-01', FechaFinal = '2022-12-31')
expectativas_inflacion.index = pd.PeriodIndex(expectativas_inflacion.index)

# Los datos inicialmente se reportan en %, lo dividmos entre 100 para poder hacer el cálculo
expectativas_inflacion['epi'] = expectativas_inflacion['epi']/100

# Definimos el ipe12m de acuerdo a la ecuación descrita anteriormente
expectativas_inflacion['ipe12m'] = expectativas_inflacion['ipc']* (1 + expectativas_inflacion['epi'])

# Hacemos el merge con la base para ya obtener la base final del bccr, el merge se hace por el índice (la fecha)
ipe12m = expectativas_inflacion['ipe12m']
base_bccr = base_bccr_temp.merge(ipe12m, how = 'inner', left_index = True, right_index = True)

Datos de SECMCA¶

La base de datos de SECMCA está en un drive accesible a cualquier persona con el link. Bajamos los datos de esta base.

In [ ]:
sheet_id = '152_Ehaf--o-na2Trodj0-Vr4yWIOhCtx'

i = pd.read_excel(f'https://docs.google.com/spreadsheets/d/{sheet_id}/export?format=xlsx')
i = i['i']
i.index = pd.period_range('2011-07', '2022-12', freq='M')
i.index.names = ['fecha']

Índice de Herfindahl-Hirschman¶

In [ ]:
hhi = pd.read_excel('https://drive.google.com/uc?id={}'.format('1H5mOmoZh5yR9ZFvmY-NmzhHrNOHwA9ji'),
  usecols=['HHI_CC','Mes'],
  index_col='Mes')
hhi.index = pd.PeriodIndex(hhi.index, freq='M')
hhi.index.names = ['fecha']

Base Final¶

In [ ]:
# Se consolidan los datos en una sola base
base = reduce(lambda left,right:pd.merge(left,right,on='fecha'), [base_bccr, hhi, i])
base.index = base.index.to_timestamp()

Añadimos a la base el número efectivo de competidores, definido como el inverso del índice de Herfindahl-Hirschman, así como transformaciones logarítmicas de ipc, imae, cred, E.

In [ ]:
base['lnipc'] = np.log(base['ipc'])
base['lnimae'] = np.log(base['imae'])
base['lncred'] = np.log(base['cred'])
base['lnE'] = np.log(base['E'])
base['ncomp'] = 1/base['HHI_CC']

Por último, realizamos un reordenamiento de las columnas y descargamos la base como Excel.

In [ ]:
# Reordenamos las variables de la base
#base = base[['imae','lnimae','ipc','lnipc','epi','ipe12m','i','cred','lncred','R','E','lnE','ncomp']]

# Visualizamos algunas observaciones de nuestra base consolidada
base.tail(5)
Out[ ]:
imae ipc cred epi R E ipe12m HHI_CC i lnipc lnimae lncred lnE ncomp
fecha
2022-08-01 111.995071 113.062770 1.768658e+07 2.80968 7.5 659.366667 116.239472 0.117422 11.68 4.727943 4.718455 16.688317 6.491280 8.516307
2022-09-01 111.697783 111.984726 1.751371e+07 3.91994 8.5 642.351429 116.374460 0.117594 11.42 4.718362 4.715797 16.678495 6.465136 8.503849
2022-10-01 113.270731 111.131513 1.754836e+07 2.36522 9.0 626.962381 113.760017 0.117780 12.64 4.710714 4.729781 16.680471 6.440887 8.490421
2022-11-01 113.901807 111.288755 1.748661e+07 2.79630 9.0 613.374091 114.400722 0.118082 13.30 4.712128 4.735337 16.676946 6.418975 8.468681
2022-12-01 115.072753 111.435922 1.756736e+07 3.44901 9.0 593.758000 115.279358 0.118639 12.43 4.713450 4.745565 16.681553 6.386472 8.428942
In [ ]:
# Exportamos los datos a Excel para poder desestacionalizar variables
base.to_excel('./Datos/Base.xlsx')

LO ANTERIOR FUE EL PROCEDIMIENTO QUE HICIMOS PARA CONFORMAR LA BASE ORIGINALMENTE, AHORA IMPORTAMOS LA BASE QUE USAMOS. SE PUEDE VERIFICAR QUE TODO ES IDENTICO SALVO EL IMAE y LNIMAE, PUES CAMBIARON POR LA NUEVA DESESTACIONALIZACION. LOS DATOS ESTÁN EN LA MISMA CARPETA QUE EL COLAB.

In [ ]:
sheet_id = '1aawDz8YugG2qOqONt8M5h8oluBU_uzLs'

base = pd.read_excel(f'https://docs.google.com/spreadsheets/d/{sheet_id}/export?format=xlsx')


base.index = base['fecha']

base = base[['imae','ipc', 'cred', 'epi','ipe12m','R','E','ipe12m', 'i','lnipc','lnimae','lncred','lnE','ncomp']]

base.tail(5)
Out[ ]:
imae ipc cred epi ipe12m R E ipe12m i lnipc lnimae lncred lnE ncomp
fecha
2022-08-01 112.694931 113.062770 1.768658e+07 2.80968 116.239472 7.5 659.366667 116.239472 11.68 4.727943 4.724684 16.688317 6.491280 8.516307
2022-09-01 112.144385 111.984726 1.751371e+07 3.91994 116.374460 8.5 642.351429 116.374460 11.42 4.718362 4.719787 16.678495 6.465136 8.503849
2022-10-01 113.663680 111.131513 1.754836e+07 2.36522 113.760017 9.0 626.962381 113.760017 12.64 4.710714 4.733244 16.680471 6.440887 8.490421
2022-11-01 114.291444 111.288755 1.748661e+07 2.79630 114.400722 9.0 613.374091 114.400722 13.30 4.712128 4.738752 16.676946 6.418975 8.468681
2022-12-01 115.471478 111.435922 1.756736e+07 3.44901 115.279358 9.0 593.758000 115.279358 12.43 4.713450 4.749024 16.681553 6.386472 8.428942

Análisis descriptivo¶

IMAE¶

Se importa la serie no desestacionalizada para comparar.

In [ ]:
imae_no_des = SW(imae_nd = 87703, FechaInicio = '2011-07-01', FechaFinal = '2022-12-31')
imae_no_des.index = pd.PeriodIndex(imae_no_des.index)
In [ ]:
fig, ax = plt.subplots()

imae_no_des['imae_nd'].plot(ax = ax, label = 'IMAE')
base['imae'].plot(ax = ax, label = 'IMAE Desestacionalizado')

plt.legend()

ax.set_title("Índice Mensual de Actividad Económica")

ax.set_ylabel('Índice')
ax.set_xlabel('Fecha')

plt.show()
fig.savefig('./PNG/imae.png')
fig.savefig('./PDF/imae.pdf')
No description has been provided for this image

IPC¶

In [ ]:
fig, ax = plt.subplots()
base['ipc'].plot(ax = ax, label = 'IPC')

ax.set_title("Índice de Precios al Consumidor")
ax.set_ylabel('Índice')
ax.set_xlabel('Fecha')

plt.show()
fig.savefig('./PNG/ipc.png')
fig.savefig('./PDF/ipc.pdf')
No description has been provided for this image

Expectativas de inflación¶

In [ ]:
fig, ax = plt.subplots()
base['epi'].plot(ax = ax, label = 'Expectativas de Inflación')

ax.set_title("Expectativas de Inflación a 12 meses ")
ax.set_ylabel('Porcentaje (%)')
ax.set_xlabel('Fecha')

plt.show()
fig.savefig('./PNG/epi.png')
fig.savefig('./PDF/epi.pdf')
No description has been provided for this image

Tasa de interés¶

In [ ]:
fig, ax = plt.subplots()
(base['i']).plot(ax = ax, label='')

ax.set_title("Tasa de interés activa promedio")
ax.set_ylabel('Porcentaje (%)')
ax.set_xlabel('Fecha')

plt.show()
fig.savefig('./PNG/i.png')
fig.savefig('./PDF/i.pdf')
No description has been provided for this image

Crédito¶

In [ ]:
fig, ax = plt.subplots()
(base['cred']/1000000).plot(ax = ax, label = 'Crédito')

ax.set_title("Crédito al Sector Privado")
ax.set_ylabel('Billones de Colones')
ax.set_xlabel('Fecha')

plt.show()
fig.savefig('./PNG/cred.png')
fig.savefig('./PDF/cred.pdf')
No description has been provided for this image

Tipo de Cambio¶

In [ ]:
fig, ax = plt.subplots()
base['E'].plot(ax = ax, label = 'Tipo de Cambio')

ax.set_title("Tipo de Cambio Promedio en Monex")
ax.set_ylabel('Colones')
ax.set_xlabel('Fecha')

plt.show()
fig.savefig('./PNG/E.png')
fig.savefig('./PDF/E.pdf')
No description has been provided for this image

TPM¶

In [ ]:
fig, ax = plt.subplots()
base['R'].plot(ax = ax, label='')

ax.set_title("Tasa de Política Monetaria")
ax.set_ylabel('Porcentaje (%)')
ax.set_xlabel('Fecha')

plt.show()
fig.savefig('./PNG/tpm.png')
fig.savefig('./PDF/tpm.pdf')
No description has been provided for this image

Competidores efectivos¶

In [ ]:
fig, ax = plt.subplots()
base['ncomp'].plot(ax = ax, label = 'Número de Competidores Efectivos')

ax.set_title("Competidores efectivos en el mercado de crédito")
ax.set_ylabel('Cantidad')
ax.set_xlabel('Fecha')

plt.show()
fig.savefig('./PNG/ncomp.png')
fig.savefig('./PDF/ncomp.pdf')
No description has been provided for this image

Análisis estadístico¶

Prueba ADF¶

In [ ]:
def ADFTable(variables, lag_criteria):
    prueba_dickey_fuller = pd.DataFrame.from_dict({'Serie': [], 'Sin constante sin tendencia': [],
                                                   'Con constante sin tendencia': [], 'Con constante con tendencia': [],
                                                   'Primera Diferencia':[]})
    for variable in variables:
        row_data = []
        for test_type in ['n', 'c', 'ct']:
            results = adfuller(base[variable], autolag = lag_criteria, regression = test_type)
            test_stat = results[0]
            crit_1 = list(results[4].values())[0]
            crit_5 = list(results[4].values())[1]
            crit_10 = list(results[4].values())[2]
            if test_stat <= crit_1:
                test_result = str(round(test_stat,2)) + '***'
            elif test_stat <= crit_5:
                test_result = str(round(test_stat,2)) + '**'
            elif test_stat <= crit_10:
                test_result = str(round(test_stat,2)) + '*'
            else:
                test_result = str(round(test_stat,2))
            row_data.append(test_result)
        first_difference = adfuller(base[variable].diff(1).dropna(), autolag = lag_criteria, regression = 'n')
        test_stat_dif = first_difference[0]
        crit_1 = list(first_difference[4].values())[0]
        crit_5 = list(first_difference[4].values())[1]
        crit_10 = list(first_difference[4].values())[2]
        if test_stat_dif <= crit_1:
            test_result_dif = str(round(test_stat_dif,2)) + '***'
        elif test_stat_dif <= crit_5:
            test_result_dif = str(round(test_stat_dif,2)) + '**'
        elif test_stat_dif <= crit_10:
            test_result_dif = str(round(test_stat_dif,2)) + '*'
        else:
            test_result_dif = str(round(test_stat_dif,2))
        row_data.append(test_result_dif)
        row = pd.DataFrame.from_dict({'Serie': [variable], 'Sin constante sin tendencia': [row_data[0]],
                                      'Con constante sin tendencia': [row_data[1]], 'Con constante con tendencia': [row_data[2]],
                                      'Primera Diferencia':[row_data[3]]})
        prueba_dickey_fuller = pd.concat([prueba_dickey_fuller, row], ignore_index = True)

    return prueba_dickey_fuller
In [ ]:
ADFTable(['lnimae', 'lnipc', 'epi', 'i', 'lncred','lnE', 'R', 'ncomp'], 'AIC')
Out[ ]:
Serie Sin constante sin tendencia Con constante sin tendencia Con constante con tendencia Primera Diferencia
0 lnimae 2.55 -0.89 -2.3 -13.61***
1 lnipc 3.32 -0.44 -2.19 -3.04***
2 epi -0.69 -2.98** -3.26* -5.68***
3 i -1.13 -1.12 -1.91 -10.08***
4 lncred 3.0 -2.73* -0.88 -3.09***
5 lnE 1.0 -1.72 -4.49*** -4.4***
6 R -0.79 -3.15** -3.05 -3.38***
7 ncomp 1.58 -1.63 0.06 -2.25**

Prueba KPSS¶

Seguimos el procedimiento de Castrillo y Rodríguez (2009) por lo que no empleamos la prueba KPSS.

In [ ]:
prueba_kpss = pd.DataFrame.from_dict({'Serie':[], 'Niveles':[], 'Tendencia':[]})

for varname in ['lnimae','lnipc', 'epi', 'i','lncred', 'R' ,'lnE', 'ncomp']:
  prueba_c = kpss(base[varname], regression = 'c', nlags = 'auto')

  if prueba_c[0] >= list(prueba_c[3].values())[-1]:
    significancia_c = '***'
  elif prueba_c[0] >= list(prueba_c[3].values())[-3]:
    significancia_c = '**'
  elif prueba_c[0] >= list(prueba_c[3].values())[-4]:
    significancia_c = '*'
  else:
    significancia_c = ''

  estadistico_c = str(np.round(prueba_c[0], 2)) + significancia_c

  prueba_ct = kpss(base[varname], regression = 'ct', nlags = 'auto')

  if prueba_ct[0] >= list(prueba_ct[3].values())[-1]:
    significancia_ct = '***'
  elif prueba_ct[0] >= list(prueba_ct[3].values())[-3]:
    significancia_ct = '**'
  elif prueba_ct[0] >= list(prueba_ct[3].values())[-4]:
    significancia_ct = '*'
  else:
    significancia_ct = ''

  estadistico_ct = str(np.round(prueba_ct[0], 3)) + significancia_ct

  row = pd.DataFrame.from_dict({'Serie':[varname], 'Niveles': [estadistico_c], 'Tendencia': [estadistico_ct]})

  prueba_kpss = pd.concat([prueba_kpss, row], ignore_index = True)

prueba_kpss

<ipython-input-265-d07243ddd272>:4: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is smaller than the p-value returned.


<ipython-input-265-d07243ddd272>:17: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is smaller than the p-value returned.


<ipython-input-265-d07243ddd272>:4: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is smaller than the p-value returned.


<ipython-input-265-d07243ddd272>:17: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is greater than the p-value returned.


<ipython-input-265-d07243ddd272>:4: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is smaller than the p-value returned.


<ipython-input-265-d07243ddd272>:17: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is greater than the p-value returned.


<ipython-input-265-d07243ddd272>:4: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is smaller than the p-value returned.


<ipython-input-265-d07243ddd272>:17: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is smaller than the p-value returned.


<ipython-input-265-d07243ddd272>:4: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is greater than the p-value returned.


<ipython-input-265-d07243ddd272>:17: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is greater than the p-value returned.


<ipython-input-265-d07243ddd272>:4: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is smaller than the p-value returned.


<ipython-input-265-d07243ddd272>:17: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is greater than the p-value returned.


<ipython-input-265-d07243ddd272>:4: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is smaller than the p-value returned.


<ipython-input-265-d07243ddd272>:17: InterpolationWarning:

The test statistic is outside of the range of p-values available in the
look-up table. The actual p-value is smaller than the p-value returned.
Out[ ]:
Serie Niveles Tendencia
0 lnimae 1.9*** 0.234***
1 lnipc 1.87*** 0.181**
2 epi 0.59** 0.099
3 i 1.73*** 0.093
4 lncred 1.89*** 0.509***
5 R 0.29 0.092
6 lnE 1.86*** 0.065
7 ncomp 1.94*** 0.319***

Prueba Zivot y Andrews¶

In [ ]:
zivot_andrews(base['lnimae'], autolag = 'AIC', regression = 'c')
Out[ ]:
(-3.6520362953164525,
 0.5800620009610958,
 {'1%': -5.27644, '5%': -4.81067, '10%': -4.56618},
 1,
 100)
In [ ]:
zivot_andrews(base['lnimae'], autolag = 'AIC', regression = 't')
Out[ ]:
(-2.4461340182846674,
 0.9151003087126337,
 {'1%': -5.03421, '5%': -4.4058, '10%': -4.13678},
 1,
 49)
In [ ]:
zivot_andrews(base['lnimae'], autolag = 'AIC', regression = 'ct')
Out[ ]:
(-6.679868534254452,
 0.0009922414372378083,
 {'1%': -5.57556, '5%': -5.07332, '10%': -4.82668},
 1,
 103)

En este caso se rechaza la hipótesis a cualquier nivel de significancia

In [ ]:
base.iloc[103] # Para ver la fecha del quiebre
Out[ ]:
imae      1.058404e+02
ipc       9.950899e+01
cred      1.576620e+07
epi       1.978560e+00
ipe12m    1.014778e+02
R         2.250000e+00
E         5.722565e+02
ipe12m    1.014778e+02
i         1.367000e+01
lnipc     4.600248e+00
lnimae    8.684901e-02
lncred    1.657338e+01
lnE       6.349587e+00
ncomp     8.671826e+00
Name: 2020-02-01 00:00:00, dtype: float64
In [ ]:
zivot_andrews(base['lnipc'], autolag = 'AIC', regression = 'c')
Out[ ]:
(-2.791006675222394,
 0.956074370684779,
 {'1%': -5.27644, '5%': -4.81067, '10%': -4.56618},
 5,
 54)
In [ ]:
zivot_andrews(base['lnipc'], autolag = 'AIC', regression = 't')
Out[ ]:
(-3.1149968725443267,
 0.5856770198750145,
 {'1%': -5.03421, '5%': -4.4058, '10%': -4.13678},
 5,
 117)
In [ ]:
zivot_andrews(base['lnipc'], autolag = 'AIC', regression = 'ct')
Out[ ]:
(-2.9618618283897655,
 0.9698639662155544,
 {'1%': -5.57556, '5%': -5.07332, '10%': -4.82668},
 5,
 113)
In [ ]:
zivot_andrews(base['i'], autolag = 'AIC', regression = 'c')
Out[ ]:
(-2.3432065420038217,
 0.9910077890162107,
 {'1%': -5.27644, '5%': -4.81067, '10%': -4.56618},
 2,
 89)
In [ ]:
zivot_andrews(base['i'], autolag = 'AIC', regression = 't')
Out[ ]:
(-2.0005051691002245,
 0.9892509301346064,
 {'1%': -5.03421, '5%': -4.4058, '10%': -4.13678},
 2,
 117)
In [ ]:
zivot_andrews(base['i'], autolag = 'AIC', regression = 'ct')
Out[ ]:
(-2.916431875978606,
 0.9742319813531298,
 {'1%': -5.57556, '5%': -5.07332, '10%': -4.82668},
 2,
 103)
In [ ]:
zivot_andrews(base['lncred'], autolag = 'AIC', regression = 'c')
Out[ ]:
(-2.996509303946015,
 0.9125126666077129,
 {'1%': -5.27644, '5%': -4.81067, '10%': -4.56618},
 2,
 87)
In [ ]:
zivot_andrews(base['lncred'], autolag = 'AIC', regression = 't')
Out[ ]:
(-3.9090006717836747,
 0.16479603874873106,
 {'1%': -5.03421, '5%': -4.4058, '10%': -4.13678},
 2,
 65)
In [ ]:
zivot_andrews(base['lncred'], autolag = 'AIC', regression = 'ct')
Out[ ]:
(-3.819290394439536,
 0.6234384512724646,
 {'1%': -5.57556, '5%': -5.07332, '10%': -4.82668},
 2,
 57)
In [ ]:
zivot_andrews(base['ncomp'], autolag = 'AIC', regression = 'c')
Out[ ]:
(-1.699988096195281,
 0.9990000000000001,
 {'1%': -5.27644, '5%': -4.81067, '10%': -4.56618},
 3,
 117)
In [ ]:
zivot_andrews(base['ncomp'], autolag = 'AIC', regression = 't')
Out[ ]:
(-3.015432173426013,
 0.6474704589007232,
 {'1%': -5.03421, '5%': -4.4058, '10%': -4.13678},
 3,
 104)
In [ ]:
zivot_andrews(base['ncomp'], autolag = 'AIC', regression = 'ct')
Out[ ]:
(-2.9505582156994214,
 0.9709849461751264,
 {'1%': -5.57556, '5%': -5.07332, '10%': -4.82668},
 3,
 94)

Estimación VAR¶

In [ ]:
endogenas = ['lnimae', 'lnipc', 'epi', 'i', 'lncred', 'lnE', 'R']

model = VAR(endog = base[endogenas], exog = base[['ncomp']])

max_lags = 10

model.select_order(max_lags).summary()

/usr/local/lib/python3.10/dist-packages/statsmodels/tsa/base/tsa_model.py:473: ValueWarning:

No frequency information was provided, so inferred frequency MS will be used.
Out[ ]:
VAR Order Selection (* highlights the minimums)
AIC BIC FPE HQIC
0 -65.14 -64.82 5.151e-29 -65.01
1 -81.18 -79.78* 5.530e-36 -80.61
2 -81.71 -79.21 3.303e-36 -80.69*
3 -81.83 -78.24 2.970e-36 -80.37
4 -81.89* -77.21 2.915e-36* -79.98
5 -81.69 -75.92 3.741e-36 -79.35
6 -81.39 -74.53 5.504e-36 -78.60
7 -81.28 -73.33 6.930e-36 -78.05
8 -81.19 -72.14 9.028e-36 -77.51
9 -81.46 -71.32 8.703e-36 -77.34
10 -81.72 -70.49 9.142e-36 -77.15

Se selecciona el criterio AIC y FPE de 4 rezagos. Inicialmente, se estimó con el criterio BIC de 1 rezago pero los residuos estaban autocorrelacionados.

In [ ]:
# Para el test el número rezagos en especificación VECM es -1 a los rezagos del VAR
# Se incluyen únicamente las variables I(1)
# 10% significancia
rank_test = select_coint_rank(base[['lnipc', 'lncred', 'i' ,'ncomp']], det_order = 0, k_ar_diff = model.select_order().aic -1, method = 'trace', signif = 0.1)

rank_test.summary()
Out[ ]:
Johansen cointegration test using trace test statistic with 10% significance level
r_0 r_1 test statistic critical value
0 4 60.77 44.49
1 4 21.70 27.07
In [ ]:
# 5% significancia
rank_test = select_coint_rank(base[['lnipc', 'lncred', 'i' ,'ncomp']], det_order = 0, k_ar_diff = model.select_order().aic -1, method = 'trace', signif = 0.05)

rank_test.summary()
Out[ ]:
Johansen cointegration test using trace test statistic with 5% significance level
r_0 r_1 test statistic critical value
0 4 60.77 47.85
1 4 21.70 29.80
In [ ]:
# 1% significancia
rank_test = select_coint_rank(base[['lnipc', 'lncred', 'i' ,'ncomp']], det_order = 0, k_ar_diff = model.select_order().aic -1, method = 'trace', signif = 0.01)

rank_test.summary()
Out[ ]:
Johansen cointegration test using trace test statistic with 1% significance level
r_0 r_1 test statistic critical value
0 4 60.77 54.68
1 4 21.70 35.46
In [ ]:
rank_test = select_coint_rank(base[['lnipc', 'lncred', 'i', 'ncomp']], det_order = 0, k_ar_diff = model.select_order().aic -1, method = 'maxeig', signif = 0.1)

rank_test.summary()
Out[ ]:
Johansen cointegration test using maximum eigenvalue test statistic with 10% significance level
r_0 r_1 test statistic critical value
0 1 39.07 25.12
1 2 13.28 18.89
In [ ]:
rank_test = select_coint_rank(base[['lnipc', 'lncred', 'i', 'ncomp']], det_order = 0, k_ar_diff = model.select_order().aic -1, method = 'maxeig', signif = 0.05)

rank_test.summary()
Out[ ]:
Johansen cointegration test using maximum eigenvalue test statistic with 5% significance level
r_0 r_1 test statistic critical value
0 1 39.07 27.59
1 2 13.28 21.13
In [ ]:
rank_test = select_coint_rank(base[['lnipc', 'lncred', 'i', 'ncomp']], det_order = 0, k_ar_diff = model.select_order().aic -1, method = 'maxeig', signif = 0.01)

rank_test.summary()
Out[ ]:
Johansen cointegration test using maximum eigenvalue test statistic with 1% significance level
r_0 r_1 test statistic critical value
0 1 39.07 32.72
1 2 13.28 25.86

Las variables sí estan cointegradas, por ende se procede con la estimación en log-niveles siguiendo a Can et al (2020) y Castrillo et al (2008)

In [ ]:
results = model.fit(4)

En la siguiente sección obtendremos funciones impulso-respuesta a partir de impulsos de una desviación estándar, por lo que es necesario obtener los valores correspondientes.

In [ ]:
print('stdev(R)=', np.sqrt(results.sigma_u['R'].loc['R']) )
print('stdev(lnE)=', np.sqrt(results.sigma_u['lnE'].loc['lnE']) )
print('stdev(epi)=', np.sqrt(results.sigma_u['epi'].loc['epi']) )

stdev(R)= 0.2904250491122893
stdev(lnE)= 0.008612031258951842
stdev(epi)= 0.2572984603990874

IRF¶

Recuperamos las funciones impulso-respuesta para un plazo de 24 meses posteriores al choque inicial. En primer lugar, un impulso de una desviación estándar en la Tasa de Política Monetaria ('R'), con intervalos de confianza al 90% y 95%.

In [ ]:
irf = results.irf(24)
In [ ]:
endogenas_dict = {'lnimae':'lnimae',
                  'lnipc':'lnipc',
                  'epi':'epi (Porcentaje, %)',
                  'i':'i (Porcentaje, %)',
                  'lncred':'lncred',
                  'lnE':'lnE',
                  'R':'R (Porcentaje, %)'}
sig = [0.05, 0.10]

for s in sig:
  for e in endogenas_dict:
    irf.plot(orth=True, impulse='R', response=e, signif=s, figsize=(4.5,2.5), subplot_params={'ytick.major.formatter': '{:.2f}'})
    plt.xlim(0,24)
    plt.suptitle('')
    plt.ylabel(endogenas_dict[e])
    plt.savefig(f'./PDF/IRF_R_{e}_{s}.pdf')
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image

Con el objetivo de analizar el impacto de la flexibilización cambiaria, recuperamos las funciones impulso-respuesta para un impulso de una desviación estándar en el logaritmo del tipo de cambio ('lnE'), sobre los logaritmos del nivel de precios ('lnipc') y actividad económica ('lnimae'). Nivel de confianza = 95%.

In [ ]:
sig = 0.05

for e in ['lnimae','lnipc', 'lnE']:
  irf.plot(orth=True, impulse='lnE', response=e, signif=sig, figsize=(4.5,2.5), subplot_params={'ytick.major.formatter': '{:.2f}'})
  plt.xlim(0,24)
  plt.suptitle('')
  plt.ylabel(e)
  plt.savefig(f'./PDF/IRF_lnE_{e}_{sig}.pdf')
No description has been provided for this image
No description has been provided for this image
No description has been provided for this image

De igual manera, un impacto de las expectativas de inflación ('epi') sobre el logaritmo del índice de precios al consumidor ('lnipc').

In [ ]:
irf.plot(orth=True, impulse='epi', response='lnipc', signif=0.05, figsize=(4.5,2.5), subplot_params={'ytick.major.formatter': '{:.2f}'})
plt.xlim(0,24)
plt.suptitle('')
plt.ylabel('lnipc')
plt.savefig(f'./PDF/IRF_epi_lnipc_0.05.pdf')
No description has been provided for this image

Residuos¶

Modelo¶

In [ ]:
# Modelo VAR(4)
results.test_whiteness(nlags = 12).summary()
Out[ ]:
Portmanteau-test for residual autocorrelation. H_0: residual autocorrelation up to lag 12 is zero. Conclusion: fail to reject H_0 at 5% significance level.
Test statistic Critical value p-value df
435.3 439.2 0.065 392
In [ ]:
# Distintas Especificaciones del VAR
results_alternative = model.fit(10)

results_alternative.test_whiteness(nlags = 12).summary()
Out[ ]:
Portmanteau-test for residual autocorrelation. H_0: residual autocorrelation up to lag 12 is zero. Conclusion: reject H_0 at 5% significance level.
Test statistic Critical value p-value df
433.9 122.1 0.000 98

No se rechaza que los residuos sean ruido blanco a un 5%. Sí se rechaza a un 10%.

In [ ]:
results.test_normality().summary()
Out[ ]:
normality (skew and kurtosis) test. H_0: data generated by normally-distributed process. Conclusion: reject H_0 at 5% significance level.
Test statistic Critical value p-value df
280.4 23.68 0.000 14

Se rechaza a cualquier nivel de significancia que los residuos tengan una distribución normal conjunta.

Residuos Ecuación IMAE¶

In [ ]:
residuos_imae = results.resid['lnimae']
In [ ]:
residuos_imae.plot()
Out[ ]:
<Axes: xlabel='fecha'>
No description has been provided for this image
In [ ]:
fig, ax = plt.subplots()

residuos_imae.hist(bins = 20, color = 'tab:blue')

ax.set_title('Residuos de Ecuación del IMAE')
Out[ ]:
Text(0.5, 1.0, 'Residuos de Ecuación del IMAE')
No description has been provided for this image
In [ ]:
# Autocorrelación de los residuos
acorr_ljungbox(residuos_imae, lags = 12,return_df = True)
Out[ ]:
lb_stat lb_pvalue
1 0.614104 0.433247
2 0.637965 0.726888
3 1.023335 0.795606
4 1.391953 0.845592
5 3.337680 0.648078
6 3.708232 0.716090
7 3.708728 0.812648
8 3.900062 0.866026
9 3.961380 0.913936
10 4.498294 0.922082
11 4.500379 0.952935
12 7.975988 0.787004

En ningún caso se rechaza la hipótesis nula: los residuos no están autocorrelacionados.

In [ ]:
# P-value del test de normalidad
jarque_bera(residuos_imae)
Out[ ]:
SignificanceResult(statistic=27.41527003307288, pvalue=1.11390910122266e-06)

Se rechaza la hipótesis nula para cualquier nivel de significancia: los residuos no siguen una distribución normal

Residuos Ecuación IPC¶

In [ ]:
residuos_ipc = results.resid['lnipc']
In [ ]:
residuos_ipc.plot()
Out[ ]:
<Axes: xlabel='fecha'>
No description has been provided for this image
In [ ]:
fig, ax = plt.subplots()

residuos_ipc.hist(bins = 20)

ax.set_title('Residuos de Ecuación del IPC')
Out[ ]:
Text(0.5, 1.0, 'Residuos de Ecuación del IPC')
No description has been provided for this image
In [ ]:
# Autocorrelación de los residuos
acorr_ljungbox(residuos_ipc, lags = 12,return_df = True)
Out[ ]:
lb_stat lb_pvalue
1 0.273345 0.601097
2 0.275764 0.871201
3 0.289733 0.961948
4 5.431166 0.245847
5 6.009259 0.305319
6 6.154021 0.406160
7 7.819684 0.348762
8 8.276892 0.406901
9 9.157329 0.422880
10 10.109897 0.430905
11 10.415316 0.493461
12 10.698448 0.554923
In [ ]:
jarque_bera(residuos_ipc)
Out[ ]:
SignificanceResult(statistic=0.3585091459844308, pvalue=0.8358930765070464)

Residuos Expectativas de Inflación¶

In [ ]:
residuos_epi = results.resid['epi']
In [ ]:
residuos_epi.plot()
Out[ ]:
<Axes: xlabel='fecha'>
No description has been provided for this image
In [ ]:
residuos_epi.hist(bins = 20)
Out[ ]:
<Axes: >
No description has been provided for this image
In [ ]:
# Autocorrelación de los residuos
acorr_ljungbox(residuos_epi, lags = 12,return_df = True)
Out[ ]:
lb_stat lb_pvalue
1 0.372371 0.541715
2 0.872509 0.646453
3 0.873032 0.831931
4 1.341394 0.854314
5 1.348119 0.929903
6 2.622619 0.854500
7 2.860814 0.897573
8 2.921438 0.939190
9 5.944589 0.745451
10 7.459204 0.681498
11 9.658783 0.561316
12 9.658897 0.645857
In [ ]:
jarque_bera(residuos_epi)
Out[ ]:
SignificanceResult(statistic=31.107470305078976, pvalue=1.7583229727799838e-07)

Residuos Tasas de Interés¶

In [ ]:
residuos_i = results.resid['i']
In [ ]:
residuos_i.plot()
Out[ ]:
<Axes: xlabel='fecha'>
No description has been provided for this image
In [ ]:
residuos_i.hist(bins = 20)
Out[ ]:
<Axes: >
No description has been provided for this image
In [ ]:
# Autocorrelación de los residuos
acorr_ljungbox(residuos_i, lags = 12,return_df = True)
Out[ ]:
lb_stat lb_pvalue
1 0.023244 0.878823
2 0.235342 0.888988
3 0.240647 0.970775
4 0.593482 0.963785
5 1.155629 0.949063
6 3.632939 0.726207
7 3.769260 0.805938
8 4.352657 0.823987
9 5.300366 0.807378
10 8.010494 0.627812
11 8.127114 0.701874
12 8.580429 0.738288
In [ ]:
jarque_bera(residuos_i)
Out[ ]:
SignificanceResult(statistic=70.78310277318377, pvalue=4.2623050321815174e-16)

Residuos Crédito¶

In [ ]:
residuos_cred = results.resid['lncred']
In [ ]:
residuos_cred.plot()
Out[ ]:
<Axes: xlabel='fecha'>
No description has been provided for this image
In [ ]:
residuos_cred.hist()
Out[ ]:
<Axes: >
No description has been provided for this image
In [ ]:
acorr_ljungbox(residuos_cred, lags = 12, return_df = True)
Out[ ]:
lb_stat lb_pvalue
1 0.025737 0.872545
2 0.172583 0.917327
3 1.435457 0.697245
4 2.818872 0.588580
5 2.999291 0.700095
6 3.959112 0.682210
7 5.410305 0.610024
8 5.518463 0.700994
9 6.362904 0.703119
10 6.369278 0.783343
11 6.415700 0.844239
12 11.190948 0.512629
In [ ]:
jarque_bera(residuos_cred)
Out[ ]:
SignificanceResult(statistic=35.33712390487484, pvalue=2.1214902535431453e-08)

Residuos Tipo de Cambio¶

In [ ]:
residuos_E= results.resid['lnE']
In [ ]:
residuos_E.plot()
Out[ ]:
<Axes: xlabel='fecha'>
No description has been provided for this image
In [ ]:
residuos_E.hist()
Out[ ]:
<Axes: >
No description has been provided for this image
In [ ]:
acorr_ljungbox(residuos_E, lags = 12, return_df = True)
Out[ ]:
lb_stat lb_pvalue
1 0.354251 0.551716
2 0.441675 0.801847
3 0.958277 0.811346
4 2.925881 0.570305
5 3.313030 0.651847
6 3.313985 0.768531
7 3.314135 0.854502
8 4.036864 0.853782
9 4.036997 0.908961
10 5.619999 0.846115
11 7.139016 0.787717
12 7.213967 0.843156
In [ ]:
jarque_bera(residuos_E)
Out[ ]:
SignificanceResult(statistic=47.3027201021324, pvalue=5.349889849772032e-11)

Residuos R¶

In [ ]:
residuos_R = results.resid['R']
In [ ]:
residuos_R.plot()
Out[ ]:
<Axes: xlabel='fecha'>
No description has been provided for this image
In [ ]:
residuos_R.hist()
Out[ ]:
<Axes: >
No description has been provided for this image
In [ ]:
acorr_ljungbox(residuos_R, lags = 12, return_df = True)
Out[ ]:
lb_stat lb_pvalue
1 0.137963 0.710315
2 0.237510 0.888025
3 1.153760 0.764115
4 2.663287 0.615654
5 2.857224 0.721986
6 4.174342 0.653096
7 5.000287 0.659928
8 5.194178 0.736635
9 7.903221 0.543928
10 8.034184 0.625498
11 9.576863 0.568801
12 10.106641 0.606606
In [ ]:
jarque_bera(residuos_R)
Out[ ]:
SignificanceResult(statistic=99.15365794117055, pvalue=2.944806396704625e-22)

En ningún caso se encuentra autocorrelación de los residuos.