In [1]:
from __future__ import print_function
import pandas as pd
import numpy as np
import tensorflow as tf
import matplotlib.pyplot as plt
%matplotlib inline
plt.rcParams['figure.figsize'] = (10, 10)
path='/tmp/'
In [2]:
import os
os.environ["CUDA_DEVICE_ORDER"]="PCI_BUS_ID"
os.environ["CUDA_VISIBLE_DEVICES"]="0"
from tensorflow.python.client import device_lib
print(device_lib.list_local_devices())
In [4]:
# Generate seasonal data dataframe
n_cases = 365*10
data=[]
for i in range(n_cases):
t = (i/100.)
s = abs((i%7)-3.5)*2.
data += list(t + s + np.random.randn(1))
index = pd.date_range('1/1/2000', periods=n_cases)
serie = pd.Series(data=data, index=index, name='value')
serie.head
Out[4]:
In [5]:
plt.plot(serie)
Out[5]:
In [6]:
plt.plot(serie['1/1/2002':'2/1/2002'])
Out[6]:
In [7]:
def extract_sequences(serie_np, maxlen=14, step=1):
'''
Cut the serie in redundant sequences of maxlen data
One sequence of length 14 for each day
'''
sequences = []
nex_data = []
for i in range(0, serie_np.shape[0] - maxlen, step):
sequences.append(serie_np[i: i + maxlen])
nex_data.append(serie_np[i + maxlen])
return np.array(sequences), np.array(nex_data)
#X_trn , y_trn = extract_sequences(seqs, maxlen=14*24, step=1)
In [8]:
# Extract sequences
n_weeks = 20
maxlen=7*n_weeks
X_trn , y_trn = extract_sequences(serie.as_matrix()[:n_cases-300], maxlen=maxlen, step=1)
X_tst , y_tst = extract_sequences(serie.as_matrix()[n_cases-300:], maxlen=maxlen, step=1)
print(X_trn.shape, y_trn.shape)
print(X_tst.shape, y_tst.shape)
print(y_tst[0])
In [9]:
# Normalize
range_trn = np.max(X_trn) - np.min(X_trn)
print(range_trn)
X_trn = (X_trn / range_trn) -0.5
y_trn = (y_trn / range_trn) -0.5
X_tst = (X_tst / range_trn) -0.5
y_tst = (y_tst / range_trn) -0.5
In [10]:
def generate_batch(X, y, batch_size=4, limit=-1, shuffle_data=True):
'''
Generate batches for one epoch
Randomize order for each epoch
'''
shuffle_index = [i for i in range(0, X.shape[0], batch_size)]
if shuffle_data:
from random import shuffle
shuffle(shuffle_index)
for i in shuffle_index[:limit]:
yield X[i:i+batch_size], y[i:i+batch_size]
gb = generate_batch(X_trn, y_trn)
In [11]:
gpu_options = tf.GPUOptions(allow_growth = True)
sess = tf.InteractiveSession(config=tf.ConfigProto(gpu_options=gpu_options, log_device_placement=True))
In [12]:
def dense(x, input_size=10, otput_size=1):
W = tf.Variable(tf.truncated_normal([input_size, otput_size], stddev=0.1))
b = tf.Variable(tf.constant(0.1, shape=[otput_size]))
return tf.matmul(x,W) + b
In [14]:
# Parameters
#features = 1
lstm_feat = 512
#Inputs
x_input = tf.placeholder(tf.float32, shape=[None, maxlen], name='x')
x_input_lstm = tf.reshape(x_input, [-1,maxlen, 1])
y_input = tf.placeholder(tf.float32, shape=[None], name='y')
# Recurrent layer
lstm1 = tf.contrib.rnn.BasicLSTMCell(lstm_feat)
lstm_out, _ = tf.nn.dynamic_rnn(lstm1, x_input_lstm, dtype=tf.float32, scope='lstm04')
print(lstm_out)
#Final dense layer
y_pred = tf.reshape(dense(lstm_out[:,-1,:], input_size=lstm_feat, otput_size=1),[-1])
print(y_pred)
# Loss function
cost = tf.reduce_sum(tf.square(y_pred - y_input))
# Trainer
learning_rate = tf.placeholder(tf.float32, name='learning_rate')
train_step = tf.train.AdamOptimizer(learning_rate=learning_rate).minimize(cost)
In [15]:
sess.run(tf.global_variables_initializer())
In [16]:
# Train graph
num_epoch=200
batchSize=256
i=0
lr=0.001
for epoch in range(num_epoch):
c_trn = []
gb = generate_batch(X_trn, y_trn, batch_size=batchSize)
for x_b, y_b in gb:
feed_dict={x_input: x_b, y_input: y_b, learning_rate: lr}
_, c, prediction, l_out = sess.run([train_step, cost, y_pred, lstm_out], feed_dict=feed_dict)
c_trn += [c]
i += 1
if i%10==0:
c_tst = cost.eval(feed_dict={x_input: X_tst, y_input: y_tst})
print('Epoch: ', epoch, ' - LR: ',lr, ' - Cost: ',np.mean(c_trn, axis=0), ' - Cost test: ',c_tst )
lr *= 0.99
In [17]:
#Score test data
predict_tst = []
real_tst = []
gb_test = generate_batch(X_tst, y_tst, batch_size=X_tst.shape[0]-1 , shuffle_data=False)
for x_b, y_b in gb_test:
p = y_pred.eval(feed_dict={x_input: x_b})
predict_tst = np.concatenate((predict_tst, p))
real_tst = np.concatenate((real_tst, y_b))
# Plot for 1 step forecast
plt.plot((real_tst[:50]+0.5)*range_trn)
plt.plot((predict_tst[:50]+0.5)*range_trn)
plt.show()
In [ ]:
In [18]:
#Recursive aplication of forecast
n_forecast = 100
x = [X_tst[0]]
y = y_tst[0]
y_fore = []
y_real = np.empty([0])
for i in range(n_forecast):
pred = y_pred.eval(feed_dict={x_input: x})
y_real = np.concatenate((y_real, np.array([y_tst[i]])))
y_fore = np.concatenate((y_fore, pred))
x = [np.concatenate((x[0][1:],pred))]
# Plot results
plt.plot((y_real+0.5)*range_trn)
plt.plot((y_fore+0.5)*range_trn)
plt.show()
# RMSE out of sample
from sklearn.metrics import mean_squared_error
from math import sqrt
rmse = sqrt(mean_squared_error((y_real+0.5)*range_trn, (y_fore+0.5)*range_trn))
print('RMSE: ', rmse)
In [ ]:
In [ ]: