<jupyter_start><jupyter_text>Table of Contents1&nbsp;&nbsp;最近N天的走势放大图2&nbsp;&nbsp;最近N天的走势预测  使用RNN预测股票下一日的收盘价这基本是一个不可能完成的任务，不过作为RNN的练习，还是一个不错的题目：有数据，有场景，有吸引力。本案例主要参考：https://github.com/DarkKnight1991/Stock-Price-Prediction ，这是一个`many-to-one`的RNN案例，即通过前60日股价数据（open,close,high,low,volume）预测下一日的收盘价，其中feature size=5。先上运行结果图：![prediction result](https://raw.githubusercontent.com/subaochen/subaochen.github.io/master/images/stock-pred-result.png)But，这里仅仅是一个练习，不作为投资参考，请参见：https://www.blueskycapitalmanagement.com/machine-learning-in-finance-why-you-should-not-use-lstms-to-predict-the-stock-market/  导入必要的包这里使用了sklearn包的MinMaxScaler进行数据的预处理，主要是归一化处理。<jupyter_code>import tensorflow as tf
import tensorflow.keras as keras
import numpy as np
import os
import sys
import time
import pandas as pd
from tqdm._tqdm_notebook import tqdm_notebook
import pickle
#from keras.models import load_model
from tensorflow.keras.layers import Dense, Dropout
#from keras.layers import LSTM
from tensorflow.keras.callbacks import ModelCheckpoint, EarlyStopping, ReduceLROnPlateau, CSVLogger
from tensorflow.keras import optimizers
# from keras.wrappers.scikit_learn import KerasClassifier
from sklearn.preprocessing import MinMaxScaler
from sklearn.model_selection import train_test_split
from sklearn.metrics import mean_squared_error
import logging
from matplotlib import pyplot as plt
os.environ['TF_CPP_MIN_LOG_LEVEL'] = '2'
logging.getLogger("tensorflow").setLevel(logging.ERROR)
os.environ['TZ'] = 'Asia/Shanghai'  # to set timezone; needed when running on cloud
time.tzset()<jupyter_output><empty_output><jupyter_text>集中设置参数<jupyter_code>params = {
    "batch_size": 20,  # 20<16<10, 25 was a bust
    "epochs": 30, # 由于启用了earlyStopping机制，通常会提前终止
    "lr": 0.00010000, # 学习率
    "time_steps": 60, # RNN的(最大？)滑动窗口大小，这是使用60日的数据预测下一日的某个特征
    "feature_size": 5 # 5个特征
}

iter_changes = "dropout_layers_0.4_0.4"
DATA_FILE="ge.us.txt"
PATH_TO_DRIVE_ML_DATA="./"
INPUT_PATH = PATH_TO_DRIVE_ML_DATA+"inputs"
OUTPUT_PATH = PATH_TO_DRIVE_ML_DATA+"outputs/"+time.strftime("%Y-%m-%d")+"/"+iter_changes
TIME_STEPS = params["time_steps"]
FEATURE_SIZE = params["feature_size"]
BATCH_SIZE = params["batch_size"]
stime = time.time()
# check if directory already exists
if not os.path.exists(OUTPUT_PATH):
    os.makedirs(OUTPUT_PATH)
    print("Directory created", OUTPUT_PATH)
else:
    os.rename(OUTPUT_PATH, OUTPUT_PATH+str(stime))
    os.makedirs(OUTPUT_PATH)
    print("Directory recreated", OUTPUT_PATH)<jupyter_output>Directory recreated ./outputs/2019-07-22/dropout_layers_0.4_0.4<jupyter_text>构造训练数据的方法参见下图，注意颜色相同的矩形块框出了输入数据和预测数据，图中的窗口尺寸（time_steps）是3，即使用前3天的数据预测下一日的收盘价。显然，如果样本数为`N`，则可划分的输入样本数为`N-time_steps`。![build serial time data](https://raw.githubusercontent.com/subaochen/subaochen.github.io/master/images/stock-pred-data.png)<jupyter_code>def print_time(text, stime):
    seconds = (time.time()-stime)
    print(text, seconds//60,"minutes : ",np.round(seconds%60),"seconds")


def trim_dataset(mat,batch_size):
    """
    trims dataset to a size that's divisible by BATCH_SIZE
    """
    no_of_rows_drop = mat.shape[0]%batch_size
    if no_of_rows_drop > 0:
        return mat[:-no_of_rows_drop]
    else:
        return mat


def build_timeseries(mat, y_col_index):
    """
    Converts ndarray into timeseries format and supervised data format. Takes first TIME_STEPS
    number of rows as input and sets the TIME_STEPS+1th data as corresponding output and so on.
    :param mat: ndarray which holds the dataset，shape：
    :param y_col_index: index of column which acts as output
    :return: returns two ndarrays-- input and output in format suitable to feed
    to LSTM.
    """
    # total number of time-series samples would be len(mat) - TIME_STEPS
    dim_0 = mat.shape[0] - TIME_STEPS
    dim_1 = mat.shape[1]
    # 先构造空的x，y
    x = np.zeros((dim_0, TIME_STEPS, dim_1))
    y = np.zeros((dim_0,))
    print("dim_0",dim_0)
    for i in tqdm_notebook(range(dim_0)):
        x[i] = mat[i:TIME_STEPS+i] # [i,TIME_STEPS+i)，x不包含TIME_STEPS+i这一行
        y[i] = mat[TIME_STEPS+i, y_col_index] # TIME_STEPS行x对应一个y。如果存在多个预测，修改这里
#         if i < 10:
#           print(i,"-->", x[i,-1,:], y[i])
    print("length of time-series i/o",x.shape,y.shape)
    return x, y


stime = time.time()
print(os.listdir(INPUT_PATH))<jupyter_output>['ge.us.txt']<jupyter_text>构造训练数据Again，RNN的输入数据要求的shape是(batch_size, time_steps, feature_size)数据归一化处理的公式：$$X_{scaled}=\frac{X-X_{min}(axis=0)}{X_{max}(axis=0)-X_{min}(axis=0)}\times(max-min)+min$$使用MinMaxScaler进行归一化处理时要注意，在训练数据集和测试数据集上要执行同样的缩放标准，因此fit只能执行一次，或者只能执行一次fit_transform函数。<jupyter_code>df_ge = pd.read_csv(os.path.join(INPUT_PATH, DATA_FILE))
print(df_ge.shape)
print(df_ge.tail())
tqdm_notebook.pandas('Processing...')
print(df_ge.dtypes)
# 训练需要的数据字段列表，去掉了日期和OpenInt
train_cols = ["Open","High","Low","Close","Volume"]
df_train, df_test = train_test_split(df_ge, train_size=0.8, test_size=0.2, shuffle=False)
print("Train--Test size", len(df_train), len(df_test))

# scale the feature MinMax, build array
# 只需要感兴趣的字段
x = df_train.loc[:,train_cols].values
min_max_scaler = MinMaxScaler()
# 几个方法的区别，参见：https://scikit-learn.org/stable/modules/generated/sklearn.preprocessing.MinMaxScaler.html
# fit，计算给定集合的最大值和最小值，以便归一化计算
# transform，在fit的基础上对数据进行归一化处理
# fit_transform，将fit和transform结合起来，在一个函数完成
x_train = min_max_scaler.fit_transform(x)
# 注意这里不能再次fit，在训练和测试数据集上要执行同样的缩放标准
x_test = min_max_scaler.transform(df_test.loc[:,train_cols].values)

print("Deleting unused dataframes of total size(KB)",
      (sys.getsizeof(df_ge)+sys.getsizeof(df_train)+sys.getsizeof(df_test))//1024)

del df_ge
del df_test
del df_train
del x

print("Are any NaNs present in train/test matrices?",np.isnan(x_train).any(), np.isnan(x_train).any())
x_t, y_t = build_timeseries(x_train, 3) # close price是第三个字段
x_t = trim_dataset(x_t, BATCH_SIZE)
y_t = trim_dataset(y_t, BATCH_SIZE)
print("Batch trimmed x_t size",x_t.shape)
print("Batch trimmed y_t size",y_t.shape)<jupyter_output>(14058, 7)
             Date   Open    High    Low  Close     Volume  OpenInt
14053  2017-11-06  20.52  20.530  20.08  20.13   60641787        0
14054  2017-11-07  20.17  20.250  20.12  20.21   41622851        0
14055  2017-11-08  20.21  20.320  20.07  20.12   39672190        0
14056  2017-11-09  20.04  20.071  19.85  19.99   50831779        0
14057  2017-11-10  19.98  20.680  19.90  20.49  100698474        0
Date        object
Open       float64
High       float64
Low        float64
Close      float64
Volume       int64
OpenInt      int64
dtype: object
Train--Test size 11246 2812
Deleting unused dataframes of total size(KB) 3267
Are any NaNs present in train/test matrices? False False
dim_0 11186<jupyter_text>创建训练模型batch_input_shape中的第二个参数分别为None和TIME_STEPS的区别* None: 可接收任意长度* TIME_STEPS(60):长度规定为60@TODO 使用CuDNNLSTM提高执行速度@TODO 神经元的参数可配置<jupyter_code>def loss(labels, logits):
  return tf.keras.losses.sparse_categorical_crossentropy(labels, logits, from_logits=True)
# 神经网络参数
lstm_units = [100,60]
lstm1_units = 100
lstm2_units = 60
dense1_units = 20
dense2_units = 1 # 和输出数据的数量一致

if tf.test.is_gpu_available():
  #rnn = tf.keras.layers.CuDNNGRU
  # CuDNNLSTM的执行速度大概是普通的LSTM的5倍以上
  lstm = tf.keras.layers.CuDNNLSTM
else:
  import functools
  #rnn = functools.partial(tf.keras.layers.GRU, recurrent_activation='sigmoid')
  lstm = functools.partial(tf.keras.layers.LSTM, recurrent_activation='sigmoid')

def create_model():
    lstm_model = tf.keras.Sequential()
    # (batch_size, timesteps, data_dim)
    # 没有激活函数？
    lstm_model.add(lstm(lstm1_units,batch_input_shape=(BATCH_SIZE,TIME_STEPS,FEATURE_SIZE),
        return_sequences=True, 
        recurrent_initializer='glorot_uniform',
        stateful=True))
    lstm_model.add(tf.keras.layers.Dropout(0.4))
    lstm_model.add(lstm(lstm2_units))
    lstm_model.add(tf.keras.layers.Dropout(0.4))
    lstm_model.add(tf.keras.layers.Dense(dense1_units,activation='relu'))
    # 这里只预测下一个价格，因此Dense的神经元个数为1
    lstm_model.add(tf.keras.layers.Dense(dense2_units,activation='sigmoid'))
    # 在这里SGD很难得到理想的结果，RMSprop一般可以比较好的收敛
    #optimizer = optimizers.RMSprop(lr=params["lr"])
    optimizer = tf.train.AdamOptimizer();
    #optimizer = optimizers.SGD(lr=0.000001, decay=1e-6, momentum=0.9, nesterov=True)
    #lstm_model.compile(loss='mean_squared_error', optimizer=optimizer)
    lstm_model.compile(loss='mean_squared_error', optimizer=optimizer)
    #lstm_model.build()
    return lstm_model


model = None
try:
    model = pickle.load(open("lstm_model", 'rb'))
    print("Loaded saved model:",model)
    model.summary()
except FileNotFoundError:
    print("Model not found")
    model = create_model()

model.summary()<jupyter_output>Model not found
_________________________________________________________________
Layer (type)                 Output Shape              Param #   
=================================================================
cu_dnnlstm (CuDNNLSTM)       (20, 60, 100)             42800     
_________________________________________________________________
dropout (Dropout)            (20, 60, 100)             0         
_________________________________________________________________
cu_dnnlstm_1 (CuDNNLSTM)     (20, 60)                  38880     
_________________________________________________________________
dropout_1 (Dropout)          (20, 60)                  0         
_________________________________________________________________
dense (Dense)                (20, 20)                  1220      
_________________________________________________________________
dense_1 (Dense)              (20, 1)                   21        
============================================================[...]<jupyter_text>准备测试数据和验证数据<jupyter_code>x_temp, y_temp = build_timeseries(x_test, 3)
x_val, x_test_t = np.split(trim_dataset(x_temp, BATCH_SIZE),2)
y_val, y_test_t = np.split(trim_dataset(y_temp, BATCH_SIZE),2)

print("Test size", x_test_t.shape, y_test_t.shape, x_val.shape, y_val.shape)<jupyter_output>dim_0 2752<jupyter_text>BATCH_SIZE对执行速度影响很大，以下是一些测试结果：|BATCH_SIZE|耗时（s/epoch)||------|------||10|50||20|50||512|5|但是，过大的batch_size会影响预测结果，参见：https://datascience.stackexchange.com/questions/16807/why-mini-batch-size-is-better-than-one-single-batch-with-all-training-data<jupyter_code>is_update_model = True
if model is None or is_update_model:
#     from keras import backend as K
#     print("Building model...")
#     print("checking if GPU available", K.tensorflow_backend._get_available_gpus())

    es = EarlyStopping(monitor='val_loss', mode='min', verbose=1,
                       patience=40, min_delta=0.0001)
    
    mcp = ModelCheckpoint(os.path.join(OUTPUT_PATH,
                          "best_model.h5"), monitor='val_loss', verbose=1,
                          save_best_only=True, save_weights_only=False, mode='min', period=1)

    # Not used here. But leaving it here as a reminder for future
    r_lr_plat = ReduceLROnPlateau(monitor='val_loss', factor=0.1, patience=30, 
                                  verbose=0, mode='auto', min_delta=0.0001, cooldown=0, min_lr=0)
    
    csv_logger = CSVLogger(os.path.join(OUTPUT_PATH, 'training_log_' + time.ctime().replace(" ","_") + '.log'), append=True)

    history = model.fit(x_t, y_t, epochs=params["epochs"], verbose=1, batch_size=BATCH_SIZE,
                        shuffle=False, validation_data=(trim_dataset(x_val, BATCH_SIZE),
                        trim_dataset(y_val, BATCH_SIZE)),callbacks=[es,mcp,csv_logger])
    
    #print("saving model...")
    #pickle.dump(model, open("lstm_model", "wb"))
#model.evaluate(x_test_t, y_test_t, batch_size=BATCH_SIZE)<jupyter_output><empty_output><jupyter_text>预测根据x_test_t进行预测<jupyter_code># Visualize the training data
from matplotlib import pyplot as plt
plt.figure()
plt.plot(history.history['loss'])
plt.plot(history.history['val_loss'])
plt.title('Model loss')
plt.ylabel('Loss')
plt.xlabel('Epoch')
plt.legend(['Train', 'Test'])
plt.show()
plt.savefig(os.path.join(OUTPUT_PATH, 'train_vis_BS_'+str(BATCH_SIZE)+"_"+time.ctime()+'.png'))
def plot_pred(pred, real, title):
    """绘制预测和实际的比较图"""
    plt.figure()
    plt.plot(pred)
    plt.plot(real)
    plt.title(title)
    plt.ylabel('Price')
    plt.xlabel('Days')
    plt.legend(['Prediction', 'Real'])
    plt.show()

def stock_pred(x, y, debug=True, title='Prediction vs Real Stock Price', batch_size=BATCH_SIZE):
    '''
    收盘价预测
    但是，这里没有越雷池一步，即没有往前走哪怕一小步。
    x: 输入样本列表
    y：labels，真实股价列表
    batch_size: 进行预测时的batch_size。当给出的样本数量小时，需要减少batch_size以适配之
    '''
    print('batch size:',batch_size)
    y_pred = model.predict(trim_dataset(x, batch_size), batch_size=batch_size)
    if debug==True:
        print("y_pred-shape",y_pred,y_pred.shape)
    y_pred = y_pred.flatten()
    y_test_t = trim_dataset(y, batch_size)
    error = mean_squared_error(y_test_t, y_pred)
    if debug==True:
        print("Error is", error, y_pred.shape, y_test_t.shape)
        print("y_pred:",y_pred[0:15])
        print("y_test_t:",y_test_t[0:15])
    y_pred_org = (y_pred * min_max_scaler.data_range_[3]) + min_max_scaler.data_min_[3]
    y_test_t_org = (y_test_t * min_max_scaler.data_range_[3]) + min_max_scaler.data_min_[3]
    if debug==True:
        print("y_pred_org:",y_pred_org[0:15])
        print("y_test_t_org:",y_test_t_org[0:15])

    # Visualize the prediction
    plot_pred(y_pred_org, y_test_t_org, title)
    plt.savefig(os.path.join(OUTPUT_PATH, 'pred_vs_real_BS'+str(batch_size)+"_"+str(len(x))+'.png'))
stock_pred(x_test_t, y_test_t,debug=True)
print_time("program completed ", stime)<jupyter_output>batch size: 20
y_pred-shape [[0.54787594]
 [0.5478797 ]
 [0.5478779 ]
 ...
 [0.547936  ]
 [0.54793763]
 [0.5479462 ]] (1360, 1)
Error is 0.0093615331541344 (1360,) (1360,)
y_pred: [0.54787594 0.5478797  0.5478779  0.547878   0.54787624 0.54787683
 0.54787815 0.5478748  0.54787284 0.5478715  0.5478725  0.54787284
 0.54787433 0.54787296 0.54787564]
y_test_t: [0.32378063 0.32499919 0.32800358 0.32905407 0.32905407 0.33031465
 0.32943225 0.33205846 0.32659593 0.32747834 0.31770881 0.31092267
 0.31258244 0.32378063 0.32592362]
y_pred_org: [26.53626  26.53644  26.536354 26.53636  26.536276 26.536304 26.536366
 26.536207 26.536114 26.536049 26.536097 26.536114 26.536184 26.53612
 26.536247]
y_test_t_org: [15.87  15.928 16.071 16.121 16.121 16.181 16.139 16.264 16.004 16.046
 15.581 15.258 15.337 15.87  15.972]<jupyter_text>将预测数据添加到真实数据，预测未来的价格走势<jupyter_code>print("x_test+t.shape:",x_test_t.shape)<jupyter_output>x_test+t.shape: (1376, 60, 5)<jupyter_text>最近N天的走势放大图将测试集的最近N天的走势放大了来看看，似乎预测曲线更加平顺。<jupyter_code># for days in (200,100,60,20,10):
#     y_pred = y_pred_org[-days:]
#     y_test_t = y_test_t_org[-days:]
#     plot_pred(y_pred,y_test_t,'Prediction vs Real Stock Price(latest ' + str(days) + ' days)')<jupyter_output><empty_output><jupyter_text>最近N天的走势预测在执行这个block之前，最好将测试数据重新生成一遍，以免相互影响造成测试数据的改变。 在单步预测之前，需要根据checkpoint重构模型，输入的shape修正为(1,1)，否则输入必须和训练时的batch_size相同，很受限制。<jupyter_code>x_temp, y_temp = build_timeseries(x_test, 3)
_, x_test_t = np.split(x_temp,2)
_, y_test_t = np.split(y_temp,2)

for days in (200,100,60,20):
    x = x_test_t[-days:]
    y = y_test_t[-days:]
    stock_pred(x, y, debug=False, title='Prediction vs Real Stock Price(latest ' + str(days) + ' days)')<jupyter_output>dim_0 2752
{'fraction_non_alphanumeric': 0.1330239246960387, 'fraction_numerical': 0.0685710550398745, 'mean_word_length': 4.714979454613373, 'pattern_counts': {'":': 5, '<': 40, '<?xml version=': 0, '>': 39, 'https://': 6, 'lorem ipsum': 0, 'www.': 1, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '374050', 'license': 'MIT', 'n_tokens_mistral': 7019, 'n_tokens_neox': 6547, 'n_words': 924, 'path': 'resources/stock-pred-rnn-tf1.ipynb', 'repo_name': 'subaochen/subaochen.github.io', 'size': '15298'}