Exercicio de machine learning para prever termos yield curve brasileira e americada

Question

Em geral, os modelos de Diebold & Li (2006) ou alguma variacao de Nelson-Siegel (1987) sao usados para prever os termos da curva de juros.

Nesse exercicio, vamos abordar uma outra perspectivo dado que iremos usar tecnicas de forecasting de machine learning.

A base de dados da curva de juros brasileira nao esta' no kaggle e foi gentilimente fornecida por amigos da Secretaria do Tesouro Nacional.

A base de dados da curva americana foi coletada no site abaixo:

base de dados da curva americana

Para fins de replicacao do codigo, salvei as bases de dados bem como o codigo no link do github abaixo:

github com as bases para download e o codigo

Comments

Os links estao com problema. Use a ferramenta de colar links que resolve.

---

Fiz a correcao ja! Muito obrigado, Professor.

Answer 1

O objetivo desse exercicio e' prever um dado termo da curva de juros (real ou nominal) brasileira ou americana usando tecnicas de machine learning.

Como coletamos tanto os dados da curva nominal e real, podemos facilmente calcular a inflacao implicita para o Brasil e o Estados Unidos. Dados de curva de juros tendem a presentar problemas de multicolinearidade dado que normalmente um terno (eg: 5 anos) pode ser aproximado por uma cominacao de outros termos (eg: 3 e 7 anos). Desse modo, na parte de previsao do exercicio criamos a seguinte regra de bolso para os modelos benchmark (lasso e ridge) : se a curva ser estimada e' a brasileira entao as variaveis independentes serao da curva americana e vice-e-versa. Em termos de dados faltantes, optou-se por preencher tais valores usando interpolacoes lineares. O legal de trabalhar com dados da curva de juros elas capturam a precificacao de variaveis economica importantes tais como juros nominal, juros real e inflacao.

importando as bibliotecas e declarando as funcoes

from sklearn import datasets
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
plt.style.use('fivethirtyeight')

import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
import matplotlib.dates as mdates
import seaborn as sns
import scipy.stats as stats
from scipy.stats import norm
from scipy.special import boxcox1p


import statsmodels
import statsmodels.api as sm
from sklearn import linear_model
from sklearn.preprocessing import StandardScaler
from sklearn.preprocessing import scale
from sklearn.linear_model import Lasso, LassoCV, Ridge, RidgeCV
from sklearn.model_selection import cross_validate
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error

import warnings
warnings.simplefilter(action='ignore', category=FutureWarning)
from keras.models import Sequential
from keras.layers import LSTM, Dense, Dropout
from keras import optimizers
from sklearn.preprocessing import MinMaxScaler


def processData(data, lb):
    X, Y = [], []
    for i in range(len(data)-lb-1):
        X.append(data[i:(i+lb)])
        Y.append(data[i+lb])
    return np.array(X), np.array(Y)


def RMSLE(predictions, realizations):
        predictions_use = predictions.clip(0)
        rmsle = np.sqrt(np.mean(np.array(np.log(predictions_use + 1) - 
                                     np.log(realizations + 1))**2))
        return rmsle  

Construindo e limpado a base dados

# import Brazil nominal and real fixed-rate  data
bra_nominal = pd.read_csv('Downloads/Curva_zero_prefixados.csv', parse_dates=['Date'], index_col = 'Date')
bra_nominal.tail()

bra_real = pd.read_csv('Downloads/Curva_zero_indice_de_precos.csv', parse_dates=['Date'], index_col = 'Date')
bra_real.tail()

for c in bra_nominal.columns:
    bra_nominal = bra_nominal.rename(columns = {c: 'BR_'+c})
bra_nominal.tail()

for c in bra_real.columns:
    bra_real = bra_real.rename(columns = {c: 'BR_'+c})
bra_real.tail()

bra_nominal = bra_nominal.drop(bra_nominal.columns[[0,1,2,3,5,7,9,16,17,18,19,20,21,22,23,24,25,26,27,28,29,30,31,32,33,34,35,36]], axis=1) 

df_bra =  pd.merge(bra_nominal, bra_real, on='Date')
df_bra['BR_inf_6M'] = (((1+df_bra['BR_6M_N']/100)/(1+df_bra['BR_6M_R']/100))-1)*100
df_bra['BR_inf_1Y'] = (((1+df_bra['BR_1Y_N']/100)/(1+df_bra['BR_1Y_R']/100))-1)*100
df_bra['BR_inf_2Y'] = (((1+df_bra['BR_2Y_N']/100)/(1+df_bra['BR_2Y_R']/100))-1)*100
df_bra['BR_inf_3Y'] = (((1+df_bra['BR_3Y_N']/100)/(1+df_bra['BR_3Y_R']/100))-1)*100
df_bra['BR_inf_4Y'] = (((1+df_bra['BR_4Y_N']/100)/(1+df_bra['BR_4Y_R']/100))-1)*100
df_bra['BR_inf_5Y'] = (((1+df_bra['BR_5Y_N']/100)/(1+df_bra['BR_5Y_R']/100))-1)*100
df_bra['BR_inf_10Y'] = (((1+df_bra['BR_10Y_N']/100)/(1+df_bra['BR_10Y_R']/100))-1)*100

inf_bra = df_bra[['BR_inf_6M','BR_inf_1Y','BR_inf_2Y','BR_inf_3Y','BR_inf_4Y','BR_inf_5Y','BR_inf_10Y']]


bra_nominal_summ = bra_nominal.describe()
bra_nominal_summ = bra_nominal_summ.transpose()
bra_nominal_summ

bra_real_summ = bra_real.describe()
bra_real_summ = bra_real_summ.transpose()
bra_real_summ

inf_bra_summ = inf_bra.describe()
inf_bra_summ = inf_bra_summ.transpose()
inf_bra_summ



# import US nominal and real fixed-rate  data

us_nominal = pd.read_csv('Downloads/USTREASURY-YIELD.csv')
us_nominal = us_nominal.rename(columns = {
        '1 MO':'1M_N',
       '2 MO':'2M_N',
       '3 MO':'3M_N',
       '6 MO':'6M_N',
       '1 YR':'1Y_N',
       '2 YR':'2Y_N',
       '3 YR':'3Y_N',
       '5 YR':'5Y_N',
       '7 YR':'7Y_N',
       '10 YR':'10Y_N',
       '20 YR':'20Y_N',
       '30 YR':'30Y_N'})    
us_nominal['Date']=pd.to_datetime(us_nominal['Date'])
us_nominal = us_nominal.set_index('Date')

for c in us_nominal.columns:
    us_nominal[c] = pd.to_numeric(us_nominal[c], errors='coerce')

us_nominal= us_nominal['2020-10-09':'2007-01-02']
us_nominal.tail()

us_real = pd.read_csv('Downloads/USTREASURY-REALYIELD.csv')
us_real = us_real.rename(columns = {
       '5 YR':'5Y_R',
       '7 YR':'7Y_R',
       '10 YR':'10Y_R',
       '20 YR':'20Y_R',
       '30 YR':'30Y_R'})    
us_real['Date']=pd.to_datetime(us_real['Date'])
us_real = us_real.set_index('Date')

for c in us_real.columns:
    us_real[c] = pd.to_numeric(us_real[c], errors='coerce')

us_real= us_real['2020-10-09':'2007-01-02']
us_real.tail()

for c in us_nominal.columns:
    us_nominal = us_nominal.rename(columns = {c: 'US_'+c})
us_nominal.tail()

for c in us_real.columns:
    us_real = us_real.rename(columns = {c: 'US_'+c})
us_real.tail()

us_nominal = us_nominal.drop(us_nominal.columns[[0,1]], axis=1) 


df_us =  pd.merge(us_nominal, us_real, on='Date')
df_us['US_inf_5Y'] = (((1+df_us['US_5Y_N']/100)/(1+df_us['US_5Y_R']/100))-1)*100
df_us['US_inf_7Y'] = (((1+df_us['US_7Y_N']/100)/(1+df_us['US_7Y_R']/100))-1)*100
df_us['US_inf_10Y'] = (((1+df_us['US_10Y_N']/100)/(1+df_us['US_10Y_R']/100))-1)*100
df_us['US_inf_20Y'] = (((1+df_us['US_20Y_N']/100)/(1+df_us['US_20Y_R']/100))-1)*100
df_us['US_inf_30Y'] = (((1+df_us['US_30Y_N']/100)/(1+df_us['US_30Y_R']/100))-1)*100

inf_us = df_us[['US_inf_5Y','US_inf_7Y','US_inf_10Y','US_inf_20Y','US_inf_30Y']]


us_nominal_summ = us_nominal.describe()
us_nominal_summ = us_nominal_summ.transpose()
us_nominal_summ

us_real_summ = us_real.describe()
us_real_summ = us_real_summ.transpose()
us_real_summ

inf_us_summ = inf_us.describe()
inf_us_summ = inf_us_summ.transpose()
inf_us_summ

df_all =  pd.merge(df_bra, df_us, on='Date')
df_all = df_all.interpolate(limit_direction='both')

df_all['BR_inv_steepness'] = df_all['BR_10Y_N'] - df_all['BR_2Y_N']
df_all['BR_inv_steepness_D'] =np.where(df_all['BR_inv_steepness'] < 0, 1, 0)

Charts para analise (code)

#Charts
plt.plot(bra_nominal_summ['max'])

fig, ax = plt.subplots(3, 1, figsize=(7,21))

ax[0].plot(bra_nominal_summ['max'], color='red', linestyle = '-.', linewidth=1, label = 'max')
ax[0].plot( bra_nominal_summ['mean'] + bra_nominal_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean + std')
ax[0].plot( bra_nominal_summ['mean'], color='blue')
ax[0].plot(bra_nominal_summ['mean'] - bra_nominal_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean - std')
ax[0].plot( bra_nominal_summ['min'], color='red', linestyle = '-.', linewidth=1, label = 'min')
ax[0].fill_between(bra_nominal_summ.index, bra_nominal_summ['min'], bra_nominal_summ['mean'] - bra_nominal_summ['std'], color='r', alpha=.1)
ax[0].fill_between(bra_nominal_summ.index, bra_nominal_summ['mean'] + bra_nominal_summ['std'], bra_nominal_summ['max'], color='r', alpha=.1)
ax[0].fill_between(bra_nominal_summ.index, bra_nominal_summ['mean'] - bra_nominal_summ['std'], bra_nominal_summ['mean'] + bra_nominal_summ['std'], color='b', alpha=.1)
ax[0].legend(fontsize=12, loc = 'upper left');
ax[0].set_ylim(-2.5, 20)
ax[0].set_xlabel('Bond duration')
ax[0].set_ylabel('Yield [%]')
ax[0].set_title('Brazil Nominal Yield Curve')

ax[1].plot(bra_real_summ['max'], color='red', linestyle = '-.', linewidth=1, label = 'max')
ax[1].plot( bra_real_summ['mean'] - bra_real_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean - std')
ax[1].plot( bra_real_summ['mean'], color='blue')
ax[1].plot(bra_real_summ['mean'] + bra_real_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean + std')
ax[1].plot( bra_real_summ['min'], color='red', linestyle = '-.', linewidth=1, label = 'min')
ax[1].fill_between(bra_real_summ.index, bra_real_summ['min'], bra_real_summ['mean'] - bra_real_summ['std'], color='r', alpha=.1)
ax[1].fill_between(bra_real_summ.index, bra_real_summ['mean'] + bra_real_summ['std'], bra_real_summ['max'], color='r', alpha=.1)
ax[1].fill_between(bra_real_summ.index, bra_real_summ['mean'] -  bra_real_summ['std'],  bra_real_summ['mean'] +  bra_real_summ['std'], color='b', alpha=.1)
ax[1].legend(fontsize=12);
ax[1].set_ylim(-2.5, 20)
ax[1].set_xlabel('Bond duration')
ax[1].set_ylabel('Yield [%]')
ax[1].set_title('Brazil Real Yield Curve')

ax[2].plot(inf_bra_summ['max'], color='red', linestyle = '-.', linewidth=1, label = 'max')
ax[2].plot( inf_bra_summ['mean'] - inf_bra_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean - std')
ax[2].plot( inf_bra_summ['mean'], color='blue')
ax[2].plot(inf_bra_summ['mean'] + inf_bra_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean + std')
ax[2].plot( inf_bra_summ['min'], color='red', linestyle = '-.', linewidth=1, label = 'min')
ax[2].fill_between(inf_bra_summ.index, inf_bra_summ['min'], inf_bra_summ['mean'] - inf_bra_summ['std'], color='r', alpha=.1)
ax[2].fill_between(inf_bra_summ.index, inf_bra_summ['mean'] + inf_bra_summ['std'], inf_bra_summ['max'], color='r', alpha=.1)
ax[2].fill_between(inf_bra_summ.index, inf_bra_summ['mean'] -  inf_bra_summ['std'],  inf_bra_summ['mean'] +  inf_bra_summ['std'], color='b', alpha=.1)
ax[2].legend(fontsize=12);
ax[2].set_ylim(-2.5, 20)
ax[2].set_xlabel('Inflation Horizon')
ax[2].set_ylabel('Inflation [%]')
ax[2].set_title('Brazil Implicit inflation')

plt.plot(us_nominal_summ['max'])

fig, ax = plt.subplots(3, 1, figsize=(7,21))

ax[0].plot(us_nominal_summ['max'], color='red', linestyle = '-.', linewidth=1, label = 'max')
ax[0].plot( us_nominal_summ['mean'] + us_nominal_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean + std')
ax[0].plot( us_nominal_summ['mean'], color='blue')
ax[0].plot(us_nominal_summ['mean'] - us_nominal_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean - std')
ax[0].plot( us_nominal_summ['min'], color='red', linestyle = '-.', linewidth=1, label = 'min')
ax[0].fill_between(us_nominal_summ.index, us_nominal_summ['min'], us_nominal_summ['mean'] - us_nominal_summ['std'], color='r', alpha=.1)
ax[0].fill_between(us_nominal_summ.index, us_nominal_summ['mean'] + us_nominal_summ['std'], us_nominal_summ['max'], color='r', alpha=.1)
ax[0].fill_between(us_nominal_summ.index, us_nominal_summ['mean'] - us_nominal_summ['std'], us_nominal_summ['mean'] + us_nominal_summ['std'], color='b', alpha=.1)
ax[0].legend(fontsize=12, loc = 'upper left');
ax[0].set_ylim(-2.5, 7.5)
ax[0].set_xlabel('Bond duration')
ax[0].set_ylabel('Yield [%]')
ax[0].set_title('US Nominal Yield Curve')

ax[1].plot(us_real_summ['max'], color='red', linestyle = '-.', linewidth=1, label = 'max')
ax[1].plot( us_real_summ['mean'] - us_real_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean - std')
ax[1].plot( us_real_summ['mean'], color='blue')
ax[1].plot(us_real_summ['mean'] + us_real_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean + std')
ax[1].plot( us_real_summ['min'], color='red', linestyle = '-.', linewidth=1, label = 'min')
ax[1].fill_between(us_real_summ.index, us_real_summ['min'], us_real_summ['mean'] - us_real_summ['std'], color='r', alpha=.1)
ax[1].fill_between(us_real_summ.index, us_real_summ['mean'] + us_real_summ['std'], us_real_summ['max'], color='r', alpha=.1)
ax[1].fill_between(us_real_summ.index, us_real_summ['mean'] -  us_real_summ['std'],  us_real_summ['mean'] +  us_real_summ['std'], color='b', alpha=.1)
ax[1].legend(fontsize=12);
ax[1].set_ylim(-2.5, 7.5)
ax[1].set_xlabel('Bond duration')
ax[1].set_ylabel('Yield [%]')
ax[1].set_title('US Real Yield Curve')

ax[2].plot(inf_us_summ['max'], color='red', linestyle = '-.', linewidth=1, label = 'max')
ax[2].plot( inf_us_summ['mean'] - inf_us_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean - std')
ax[2].plot( inf_us_summ['mean'], color='blue')
ax[2].plot(inf_us_summ['mean'] + inf_us_summ['std'], color='blue', linestyle = '-.', linewidth=1, label = 'mean + std')
ax[2].plot( inf_us_summ['min'], color='red', linestyle = '-.', linewidth=1, label = 'min')
ax[2].fill_between(inf_us_summ.index, inf_us_summ['min'], inf_us_summ['mean'] - inf_us_summ['std'], color='r', alpha=.1)
ax[2].fill_between(inf_us_summ.index, inf_us_summ['mean'] + inf_us_summ['std'], inf_us_summ['max'], color='r', alpha=.1)
ax[2].fill_between(inf_us_summ.index, inf_us_summ['mean'] -  inf_us_summ['std'],  inf_us_summ['mean'] +  inf_us_summ['std'], color='b', alpha=.1)
ax[2].legend(fontsize=12);
ax[2].set_ylim(-2.5, 7.5)
ax[2].set_xlabel('Inflation Horizon')
ax[2].set_ylabel('Inflation [%]')
ax[2].set_title('US Implicit inflation')

fig, ax = plt.subplots(4,3,figsize=(24, 12))
for i, c in enumerate(df_us.columns):
    ax[i//3][i%3].plot(df_us[c], label = c )
    ax[i//3][i%3].legend()
    ax[i//3][i%3].set_ylabel('Yield [%]')
fig.delaxes(ax[3][1])


fig, ax = plt.subplots(5,3,figsize=(24, 12))
for i, c in enumerate(df_bra.columns):
    ax[i//3][i%3].plot(df_bra[c], label = c )
    ax[i//3][i%3].set_ylabel('Yield [%]')
    ax[i//3][i%3].legend()
fig.delaxes(ax[3][2])

f, ax = plt.subplots(1,2, figsize=(18,6))
corrmatrixus = df_us.corr()
sns.heatmap(corrmatrixus, square=True, ax=ax[0])
ax[0].set_title('US')
corrmatrixj = df_bra.corr()
sns.heatmap(corrmatrixj, square=True, ax=ax[1])
ax[1].set_title('Brazil')

f, ax = plt.subplots(figsize=(22,12))
corrmatrix = df_all.corr()
sns.heatmap(corrmatrix, square=True, ax=ax)
ax.set_title('Brazil and US')


fig, ax = plt.subplots(figsize=(18,10.5))
diff = df_bra['BR_10Y_N'] - df_bra['BR_2Y_N']
plt.plot(df_bra.index, diff)
plt.fill_between(df_bra.index, 0, diff, where= diff < 0, color='red')
plt.fill_between(df_bra.index, 0, diff, where= diff > 0, color='blue')
plt.ylabel('BR_Nominal_Yield$_{10Y}$ - BR_Nominal_Yield$_{2Y}$')

fig, ax = plt.subplots(figsize=(18,10.5))
diff = df_us['US_5Y_N'] - df_us['US_1Y_N']
plt.plot(df_us.index, diff)
plt.fill_between(df_us.index, 0, diff, where= diff < 0, color='red')
plt.fill_between(df_us.index, 0, diff, where= diff > 0, color='blue')
plt.ylabel('US_Nominal_Yield$_{10Y}$ - US_Nominal_Yield$_{2Y}$')

Charts

Answer 2

Forecasting

Como modelos benchmark, foram escolhido os modelos industriais lasso e ridge. Todavia, optou-se, para fins de balanceamento, usar um cross-validation em ambos. Ademais, o codigo foi montado com vista a selicionar o \(\alpha\) em cada modelo que minimize o erro de previsao. Tendo em vista, o problema de multicolinearidade dos termos adjacentes, optei usar como variaveis independentes os dados do outro pais.

Para fins de comporacao com os modelos, industriais usou-se o modelo Long Short-term Memory network (LSTM) que e' um tipo de rede neural recorrente (RNN). Um benefício desse tipo de rede é que ela pode aprender e lembrar em longas sequências e não depende de uma observação defasada de janela pré-especificada como entrada.Por padrão, uma camada LSTM no Keras mantém o estado entre os dados em um batch. Um batch de dados é um número de linhas de tamanho fixo do conjunto de dados de treinamento que define quantos padrões processar antes de atualizar os pesos da rede.
A camada LSTM espera que a entrada esteja em uma matriz com as dimensões: [amostras, intervalos de tempo, features].

A camada LSTM espera que a entrada esteja em uma matriz com as dimensões: [amostras, intervalos de tempo, recursos]. O intervalo de tempo (lstmtimestep) indica o numero de dias N que temos que prever um valor no tempo T, com base nos dados dos dias T-N. Ja os features indicam as caracteristicas/variaveis que levaremos em conta na previsao. No nosso caso, levaremos em conta somente as taxas historicas do termo da curva a ser previsoto com logica analoga a previsao via um processo auto regressivo.

A normalizacao dos dados foi feita via funcao MinMaxScale.

Deixo aqui algumas referencia de cunho pratico que ajudaram a aprender o pouco o que sei sobre o LSTM:

link 1
link 2

No codigo abaixo, o target escolhido foi o termo de 5 anos da curva de juros nominal brasileira. O codigo foi montado de modo que usario possa alterar os principais inputs dos modelos, como, por exemplo, o target.

Para fins de comparacao com o kaggle, usamos o rmsle como metrica de erro.

Forecasting (code)

## forecasting

if __name__ == "__main__":

    target= 'BR_5Y_N'
    cv = 10 
    test_size = 0.8
    lstm_neuro = 16
    lstm_time_step = 7
    lstm_dropout = 0.25
    lstm_learning_rate = 0.001
    lstm_epochs = 200
    y = df_all[target]
    model_list = ['cv_lasso', 'cv_ridge', 'lstm']
    rmsle=[]
    df_chart_rmsle = {}
    df_all2 = df_all
    for c in list(df_all2.columns):
       st = c
       if st[0:1]==target[0:1]:    
            df_all2.drop([c], axis = 1, inplace = True) 
    X = df_all2
    #X= df_all.drop([target], axis = 1)
    X_train, X_test , y_train, y_test = train_test_split(X, y, test_size= test_size, random_state=1)

## LASSO

    lassocv = LassoCV(alphas=None, cv=cv  , max_iter=10000, normalize=True)
    lassocv.fit(X_train, y_train)
    alphas = 10**np.linspace(6,-2,50)*0.5
    lasso = Lasso(max_iter=10000, normalize=True)
    coefs = []

    for a in alphas:
        lasso.set_params(alpha=a)
        lasso.fit(X, y)
        coefs.append(lasso.coef_)

    np.shape(coefs)

    lassocv = LassoCV(alphas=None, cv=cv , max_iter=10000, normalize=True)
    lassocv.fit(X_train, y_train)
    lasso.set_params(alpha=lassocv.alpha_)
    print("Alpha=", lassocv.alpha_)
    lasso.fit(X_train, y_train)
    pred_lasso= lasso.predict(X_test)
    print("mse = ",mean_squared_error(y_test, pred_lasso))
    print("best model coefficients:")
    print(pd.Series(lasso.coef_, index=X.columns))

    plt.figure(figsize=(6,6))
    plt.scatter(y_test, pred_lasso, s = 0.8)
    plt.xlim(5,15)
    plt.ylim(5, 15)
    plt.plot([5, 15], [5, 15], color='red', linestyle='-', linewidth=3)
    plt.suptitle('Forecast vs historical', fontsize=18)
    plt.xlabel('LASSO', fontsize = 16)
    plt.ylabel('y_test', fontsize=16)
    plt.show()

    print(RMSLE(y_test,pred_lasso))
    rmsle.append(RMSLE(y_test,pred_lasso))


## RIDGE

    ridgecv = RidgeCV(alphas=alphas, normalize=True)
    ridgecv.fit(X_train, y_train)
    print("Alpha=", ridgecv.alpha_)
    ridge_best = Ridge(alpha=ridgecv.alpha_, normalize=True)
    ridge_best.fit(X_train, y_train)
    pred_ridge= ridge_best.predict(X_test)
    print("mse = ",mean_squared_error(y_test,pred_ridge))
    print("best model coefficients:")
    print(pd.Series(ridge_best.coef_, index=X.columns))

    plt.figure(figsize=(6,6))
    plt.scatter(y_test, pred_ridge, s = 0.8)
    plt.xlim(5,15)
    plt.ylim(5, 15)
    plt.plot([5, 15], [5, 15], color='red', linestyle='-', linewidth=3)
    plt.suptitle('Forecast vs historical', fontsize=18)
    plt.xlabel('Ridge', fontsize = 16)
    plt.ylabel('y_test', fontsize=16)
    plt.show()

    print(RMSLE(y_test,pred_ridge))
    rmsle.append(RMSLE(y_test,pred_ridge))


## LSTM


    ts1 = y.dropna()
    ts1_idx = ts1.index
    ts2 = ts1.values
    scl = MinMaxScaler()
    ts2 = ts2.reshape(ts2.shape[0],1)
    ts2 = scl.fit_transform(ts2)


    X, y = processData(ts2, lstm_time_step)
    X_train, X_test = X[:int(X.shape[0]* test_size)], X[int(X.shape[0]* test_size):]
    y_train, y_test = y[:int(y.shape[0]* test_size)], y[int(y.shape[0]* test_size):]
    print('Checar se a divisao foi feita da forma correta')
    print(X_train.shape)
    print(X_test.shape)
    print(y_train.shape)
    print(y_test.shape)


    model = Sequential()
    model.add(LSTM(lstm_neuro, input_shape=(lstm_time_step,1)))
    model.add(Dropout(lstm_dropout))
    model.add(Dense(1))
    adam = optimizers.adam(lr=lstm_learning_rate, clipnorm=1.)
    model.compile(optimizer=Adam, loss='mse')
    X_train = X_train.reshape((X_train.shape[0], X_train.shape[1], 1))
    X_test = X_test.reshape((X_test.shape[0], X_test.shape[1], 1))
    history = model.fit(X_train, y_train, epochs=lstm_epochs, validation_data=(X_test, y_test), shuffle=False)


    plt.plot(history.history['loss'], label = 'Training')
    plt.plot(history.history['val_loss'], label = 'Validation')
    plt.legend()


    Xt16 = model.predict(X_test)
    fig, ax = plt.subplots(1,1,figsize=(24, 12))
    ax.plot( scl.inverse_transform(y_test.reshape(-1,1)), label='Test data')
    ax.plot( scl.inverse_transform(Xt16), label ='16L')
    ax.set_ylabel('Yield [%]')
    ax.legend()

    plt.figure(figsize=(6,6))
    plt.scatter(scl.inverse_transform(Xt16),  scl.inverse_transform(y_test.reshape(-1,1)), s = 0.8)
    plt.xlim(5,15)
    plt.ylim(5, 15)
    plt.plot([5, 15], [5, 15], color='red', linestyle='-', linewidth=3)
    plt.suptitle('Forecast vs historical', fontsize=18)
    plt.xlabel('LSTM 16', fontsize = 16)
    plt.ylabel('y_valid', fontsize=16)
    plt.show()


    print(RMSLE(scl.inverse_transform(y_test.reshape(-1,1)),scl.inverse_transform(Xt16)))
    rmsle.append(RMSLE(scl.inverse_transform(y_test.reshape(-1,1)),scl.inverse_transform(Xt16)))

    #Chart rmsle

    df_chart_rmsle={'Used model':model_list,'RMSLE':rmsle}
    print(df_chart_rmsle)
    DF_chart_rmsle = pd.DataFrame(df_chart_rmsle)
    sns.factorplot(y='Used model',x='RMSLE',data=DF_chart_rmsle,kind='bar',size=5,aspect=2)

Forecasting (Charts)
O lstm tem um erro de previsao menor do que os demais. Contudo, caso a resticao de uso da curvas domestica nao esteja ativa, o lasso e ridge poussem um desempenho melhor do que lstm.

A titulo de curiosidade, deixei tambem um grafico comparando o desempenho dos 3 modelos em quesito quando o target e' o termo de 1 ano da curva nominal brasileira. E' interessante ver que o LSTM possui um erro menor do que os demais, de mod que o gap do erro de previsao entre rede neurais e modelos industriais aumentou consideravelmente.

Comments

Análise interessante Daniel. Ótima maneira de trabalhar com dados desse tipo, não é algo fácil e os gráficos e comparações ajudam a entender o o que foi proposto. Só tenho uma dúvida: ao usar LSTM, a série se torna estacionária com o uso da RNN?

OBS: os links não estão indo direto, mas basta colar eles no navegador que funciona. Em relação ao código, ao replicar tive que ajustar o optimzer.Adam, com o \(\it A \) maiúsculo. No mais, todas imagens e resultados foram replicadas. Inclusive vou guardar o modelo para futuro estudo.

---

Oi Gustavo, fico feliz que vc tenha gostado da resolucao.  Quanto a sua duvida, a normalizacao e' feita via MinMaxScaler(). Mudei o texto da resposta para deixar isso claro.  Sobre o typo do "A" na minha versao do spyder esta rodando com o minusculo. De qualquer forma, eu modifiquei a codigo conforme a sua sugestao. Muito obrigado pelos apontamentos.


## LSTM


    ts1 = y.dropna()
    ts1_idx = ts1.index
    ts2 = ts1.values
    scl = MinMaxScaler()
    ts2 = ts2.reshape(ts2.shape[0],1)
    ts2 = scl.fit_transform(ts2)

---

Pefeito Daniel! Obrigado pelo esclarecimento.