作為一種技術手段，預測在金融、證券領域的應用非常廣泛，尤其是對股票價格的預測。我們介紹一下獲得股票數據的方法，並基於此對數據進行預處理，接著使用數據分析方法，建立基礎特徵，進一步構建預測模型，且基於新數據驗證模型效果。擬使用 VAR 及LSTM兩種算法建立預測模型。
獲取股票數據
股票數據通常可從新浪股票、雅虎股票等網頁上獲取，此外還有一些炒股軟體，如同花順、通達信等都提供了非常清楚的股票數據展示和圖表呈現。如果要獲得實時的股票數據，可以考慮使用新浪股票提供的接口獲取數據。以大秦鐵路（股票代碼：601006）為例，如果要獲取它的最新行情，只需訪問新浪的股票數據接口（具體可以百度），該接口會返回一串文本，例如：
1 var hq_str_sh601006="大秦鐵路,6.980,6.960,7.010,7.070,6.950,7.010,7.020,121033256,847861533.000,18900, 7.010,214867,7.000,66500,6.990,386166,6.980,336728,6.970,273750,7.020,836066,7.030,630800,7.040,936306,7.050,579400,7.060,2016-03-18,15:00:00,00";
這個字符串由許多數據拼接在一起，不同含義的數據用逗號隔開了，按照程式設計師的思路，順序號從0開始。
0：，股票名字1：6.980>，今日開盤價2：6.960>，昨日收盤價3：7.010>，當前價格4：7.070>，今日最高價5：6.950>，今日最低價6：7.010>，競買價，即「買一」報價7：7.020>，競賣價，即「賣一」報價8：121033256>，成交的股票數，由於股票交易以一百股為基本單位，所以在使用時，通常把該值除以一百9：847861533.000>，成交金額，單位為「元」，為了一目了然，通常以「萬元」為成交金額的單位，所以通常把該值除以一萬10：18900>，「買一」申請4695股，即47手11：7.010>，「買一」報價12：214867>，「買二」13：7.000>，「買二」14：66500>，「買三」15：6.990>，「買三」16：386166>，「買四」17：6.980>，「買四」18：336728>，「買五」19：6.970>，「買五」20：273750>，「賣一」申報3100股，即31手21：7.020>，「賣一」報價(22,23),(24,25),(26,27),(28,29)分別為「賣二」至「賣四的情況」30：2016-03-18>，日期31：15:00:00>，時間
這個接口對於JavaScript程序非常方便，如果要查看該股票的日K線圖，可訪問新浪股票的K線圖接口（具體可百度），便可得到日K線圖。
日K線圖
如果要查看該股票的分時線，可訪問連結新浪股票的分時線圖接口（具體可百度），便可得到分時線圖。
分時線圖
對於周K線和月K線的查詢，可分別訪問新浪股票的周K線圖和月K線圖的接口（具體可百度）。Python中我們可以使用pandas_datareader庫來獲取股票數據，默認是訪問yahoofinance的數據，其中包括上證和深證的股票數據，還有港股數據，該庫只能獲取股票的歷史交易記錄信息：如最高價、最低價、開盤價、收盤價以及成交量，無法獲取個股的分筆交易明細歷史記錄。上證代碼是ss，深證代碼是sz，港股代碼是hk，比如茅臺：6000519.ss，萬科000002.sz，長江實業0001.hk。這裡以貴州茅臺股票為例，說明pandas_datareader庫中股票數據的獲取方法及簡單的可視化，代碼如下：
1import pandas as pd
2import pandas_datareader.data as web
3import datetime as dt
4data = web.DataReader('600519.ss','yahoo', dt.datetime(2019,8,1),dt.datetime(2019,8,31))
5data.head
6 High Low Open Close Volume Adj Close
7# Date
8# 2019-08-01 977.000000 953.020020 976.51001 959.299988 3508952 959.299988
9# 2019-08-02 957.979980 943.000000 944.00000 954.450012 3971940 954.450012
10# 2019-08-05 954.000000 940.000000 945.00000 942.429993 3677431 942.429993
11# 2019-08-06 948.000000 923.799988 931.00000 946.299988 4399116 946.299988
12# 2019-08-07 955.530029 945.000000 949.50000 945.000000 2686998 945.000000
13
14kldata=data.values[:,[2,3,1,0]] # 分別對應開盤價、收盤價、最低價和最高價
15from pyecharts import options as opts
16from pyecharts.charts import Kline
17
18kobj = Kline.add_xaxis(data.index.strftime("%Y-%m-%d").tolist).add_yaxis("貴州茅臺-日K線圖",kldata.tolist).set_global_opts(
19 yaxis_opts=opts.AxisOpts(is_scale=True),
20 xaxis_opts=opts.AxisOpts(is_scale=True),
21 title_opts=opts.TitleOpts(title=""))
22kobj.render
貴州茅臺股票日K線圖如圖：
為給定時間序列的財務圖表，代碼中對象data包含6個屬性，依次為Open（開盤價）、High（最高價）、Low（最低價）、Close（收盤價）、Volume（成交量）、Adjusted（復權收盤價）。基於收盤價的重要性，可從收盤價的歷史數據中分割訓練集、驗證集、測試集，使用適當的特徵，建立預測模型，並實施預測。
基於VAR算法的預測
向量自回歸（VAR）模型就是非結構化的多方程模型，它的核心思想不考慮經濟理論，而直接考慮經濟變量時間時序之間的關係，避開了結構建模方法中需要對系統中每個內生變量關於所有內生變量滯後值函數建模的問題，通常用來預測相關時間序列系統和研究隨機擾動項對變量系統的動態影響。VAR模型類似聯立方程，將多個變量包含在一個統一的模型中，共同利用多個變量信息，比起僅使用單一時間序列的ARIMA等模型，其涵蓋的信息更加豐富，能更好地模擬現實經濟體，因而用於預測時能夠提供更加貼近現實的預測值。此處擬基於貴州茅臺股票數據，建立VAR的預測模型。使用後30天的數據作為驗證集，剩餘的數據用於建立預測模型。本節從VAR模型的平穩性檢驗出發，依次完成VAR模型的定階及建模預測，最終通過分析驗證集上的準確率來評估預測效果。
1、平穩性檢驗
只有平穩的時間序列才能夠直接建立VAR模型，因此在建立VAR模型之前，首先要對變量進行平穩性檢驗。通常可利用序列的自相關分析圖來判斷時間序列的平穩性，如果序列的自相關係數隨著滯後階數的增加很快趨於0，即落入隨機區間，則序列是平穩的；反之，序列是不平穩的。另外，也可以對序列進行ADF檢驗來判斷平穩性。對於不平穩的序列，需要進行差分運算，直到差分後的序列平穩後，才能建立VAR模型。此處首先提取用於建立預測模型的基礎數據，並對其進行單位根檢驗，對應的Python代碼如下：
1import statsmodels.tsa.stattools as stat
4import pandas as pd
5import numpy as np
6
7data = web.DataReader('600519.ss','yahoo', dt.datetime(2014,1,1),dt.datetime(2019,9,30))
8subdata = data.iloc[:-30,:4]
9for i in range(4):
10 pvalue = stat.adfuller(subdata.values[:,i], 1)[1]
11 print("指標 ",data.columns[i]," 單位根檢驗的p值為：",pvalue)
12# 指標 High 單位根檢驗的p值為：0.9955202280850401
13# 指標 Low 單位根檢驗的p值為：0.9942509439755689
14# 指標 Open 單位根檢驗的p值為：0.9938548193990323
15# 指標 Close 單位根檢驗的p值為：0.9950049124079876
可以看到，p值都大於0.01，因此都是不平穩序列。現對subdata進行1階差分運算，並再次進行單位根檢驗，對應的Python代碼如下：
1subdata_diff1 = subdata.iloc[1:,:].values - subdata.iloc[:-1,:].values
2for i in range(4):
3 pvalue = stat.adfuller(subdata_diff1[:,i], 1)[1]
4 print("指標 ",data.columns[i]," 單位根檢驗的p值為：",pvalue)
5# 指標 High 單位根檢驗的p值為：0.0
6# 指標 Low 單位根檢驗的p值為：0.0
7# 指標 Open 單位根檢驗的p值為：0.0
8# 指標 Close 單位根檢驗的p值為：0.0
如結果所示，對這4個指標的1階差分單獨進行單位根檢驗，其p值都不超過0.01，因此可以認為是平穩的。
2、VAR模型定階
接下來就是為VAR模型定階，可以讓階數從1逐漸增加，當AIC值儘量小時，可以確定最大滯後期。我們使用最小二乘法，求解每個方程的係數，並通過逐漸增加階數，為模型定階，Python代碼如下：
1# 模型階數從1開始逐一增加
2rows, cols = subdata_diff1.shape
3aicList =
4lmList =
5
6for p in range(1,11):
7 baseData = None
8 for i in range(p,rows):
9 tmp_list = list(subdata_diff1[i,:]) + list(subdata_diff1[i-p:i].flatten)
10 if baseData is None:
11 baseData = [tmp_list]
12 else:
13 baseData = np.r_[baseData, [tmp_list]]
14 X = np.c_[[1]*baseData.shape[0],baseData[:,cols:]]
15 Y = baseData[:,0:cols]
16 coefMatrix = np.matmul(np.matmul(np.linalg.inv(np.matmul(X.T,X)),X.T),Y)
17 aic = np.log(np.linalg.det(np.cov(Y - np.matmul(X,coefMatrix),rowvar=False))) + 2*(coefMatrix.shape[0]-1)**2*p/baseData.shape[0]
18 aicList.append(aic)
19 lmList.append(coefMatrix)
20
21#對比查看階數和AIC
22pd.DataFrame({"P":range(1,11),"AIC":aicList})
23# P AIC
24# 0 1 13.580156
25# 1 2 13.312225
26# 2 3 13.543633
27# 3 4 14.266087
28# 4 5 15.512437
29# 5 6 17.539047
30# 6 7 20.457337
31# 7 8 24.385459
32# 8 9 29.438091
33# 9 10 35.785909
如上述代碼所示，當p=2時，AIC值最小為13.312225。因此VAR模型定階為2，並可從對象lmList[1]中獲取各指標對應的線性模型。
3、預測及效果驗證
基於lmList[1]中獲取各指標對應的線性模型，對未來30期的數據進行預測，並與驗證數據集進行比較分析，Python代碼如下：
1p = np.argmin(aicList)+1
2n = rows
3preddf = None
4for i in range(30):
5 predData = list(subdata_diff1[n+i-p:n+i].flatten)
6 predVals = np.matmul([1]+predData,lmList[p-1])
7 # 使用逆差分運算，還原預測值
8 predVals=data.iloc[n+i,:].values[:4]+predVals
9 if preddf is None:
10 preddf = [predVals]
11 else:
12 preddf = np.r_[preddf, [predVals]]
13 # 為subdata_diff1增加一條新記錄
14 subdata_diff1 = np.r_[subdata_diff1, [data.iloc[n+i+1,:].values[:4] - data.iloc[n+i,:].values[:4]]]
15
16#分析預測殘差情況
17(np.abs(preddf - data.iloc[-30:data.shape[0],:4])/data.iloc[-30:data.shape[0],:4]).describe
18# High Low Open Close
19# count 30.000000 30.000000 30.000000 30.000000
20# mean 0.010060 0.009380 0.005661 0.013739
21# std 0.008562 0.009968 0.006515 0.013674
22# min 0.001458 0.000115 0.000114 0.000130
23# 25% 0.004146 0.001950 0.001653 0.002785
24# 50% 0.007166 0.007118 0.002913 0.010414
25# 75% 0.014652 0.012999 0.006933 0.022305
26# max 0.039191 0.045802 0.024576 0.052800
從上述代碼第17行可以看出這4個指標的最大百分誤差率分別為3.9191%、4.5802%、2.4576%、5.28%，最小百分誤差率分別為0.1458%、0.0115%、0.0114%、0.013%，進一步，繪製二維圖表觀察預測數據與真實數據的逼近情況，Python代碼如下：
1import matplotlib.pyplot as plt
2plt.figure(figsize=(10,7))
3for i in range(4):
4 plt.subplot(2,2,i+1)
5 plt.plot(range(30),data.iloc[-30:data.shape[0],i].values,'o-',c='black')
6 plt.plot(range(30),preddf[:,i],'o--',c='gray')
7 plt.ylim(1000,1200)
8 plt.ylabel("$"+data.columns[i]+"$")
9plt.show
10v = 100*(1 - np.sum(np.abs(preddf - data.iloc[-30:data.shape[0],:4]).values)/np.sum(data.iloc[-30:data.shape[0],:4].values))
11print("Evaluation on test data: accuracy = %0.2f%% \n" % v)
12# Evaluation on test data: accuracy = 99.03%
該預測效果如下圖，其中黑色實線為真實數據，灰色虛線為預測數據，使用VAR模型進行預測的效果總體還是不錯的，平均準確率為99.03%。針對多元時間序列的情況，VAR模型不僅考慮了其他指標的滯後影響，計算效率還比較高，從以上代碼可以看到，對於模型的擬合，直接使用的最小二乘法，這增加了該模型的適應性。
預測效果
基於LSTM算法的預測
本節主要基於LSTM算法對貴州茅臺股票數據進行預測，該算法非常擅長序列數據的建模，由於引入了遺忘門等更為複雜的內部處理單元來處理上下文信息的存儲與更新，這樣既可以消除梯度問題的困擾，也可以對存在短期或長期依賴的數據建模，該算法在文本、語音等序列數據模型中廣泛使用。本節從LSTM建模的數據要求及網絡結構設計講起，通過設置合理的參數，通過訓練得到模型，並基於該模型進行預測，最後將結果與真實數據進行比較，評估預測效果。
1、數據要求
本節使用LSTM算法對貴州茅臺股票數據進行預測，可基於前N條樣本對當前樣本進行預測，因此該模型不需要像DNN那樣，將歷史數據進行複雜轉換，將基礎數據稍加處理就能用於訓練模型。對基礎數據的處理即為對該數據進行重新封裝，將樣本前N期的集合與當前樣本對應上，分別得到訓練數據的輸入與輸出。
所示數據對應關係（具體數據為示意）
2、數據預處理
首先，需要將基礎數據重構為包含歷史3周特徵數據的基礎數據，以預測日的High（最高價）、Low（最低價）、Open（開盤價）、Close（收盤價）4個指標作為輸出數據。這裡我們使用2014年1月1日至2019年8月31日的貴州茅臺股票數據作為訓練數據，使用2019年整個9月的數據作為測試數據，來驗證模型效果。用Python將對全體數據進行標準化，並將基礎數據的特徵進行重構，代碼如下：
1SEQLEN = 21
2dim_in = 4
3dim_out = 4
4pred_len = 30
5vmean = data.iloc[:,:4].apply(lambda x:np.mean(x))
6vstd = data.iloc[:,:4].apply(lambda x:np.std(x))
7t0 = data.iloc[:,:4].apply(lambda x:(x-np.mean(x))/np.std(x)).values
8X_train = np.zeros((t0.shape[0]-SEQLEN-pred_len, SEQLEN, dim_in))
9Y_train = np.zeros((t0.shape[0]-SEQLEN-pred_len, dim_out),)
10X_test = np.zeros((pred_len, SEQLEN, dim_in))
11Y_test = np.zeros((pred_len, dim_out),)
12for i in range(SEQLEN, t0.shape[0]-pred_len):
13 Y_train[i-SEQLEN] = t0[i]
14 X_train[i-SEQLEN] = t0[(i-SEQLEN):i]
15for i in range(t0.shape[0]-pred_len,t0.shape[0]):
16 Y_test[i-t0.shape[0]+pred_len] = t0[i]
17 X_test[i-t0.shape[0]+pred_len] = t0[(i-SEQLEN):i]
如上述代碼所示，SEQLEN表示使用前期數據的長度，dim_in表示輸入數據的維度，dim_out表示輸出數據的維度，pred_len表示預測數據的長度。第5~7行代碼對數據進行zscore標準化，將數據映射到標準正態分布。第12~17行代碼對基礎數據進行重構，分別得到訓練數據X_train、Y_train以及測試數據X_test、Y_test。
3、網絡結構設計
經嘗試，我們使用近3周的歷史數據來訓練LSTM模型，同時，設置隱含層神經元的數量為64。因此，我們可以將LSTM神經網絡按下面的結構進行設計（圖中N可取21，即3周對應的天數）。
LSTM神經網絡結構
4、建立模型
現基於Keras搭建LSTM神經網絡，並基於訓練集對模型進行訓練，Python代碼如下：
1from keras.layers import LSTM, Dense
2from keras.models import Sequential
3model = Sequential
4model.add(LSTM(64, input_shape=(SEQLEN, dim_in),activation='relu',recurrent_dropout=0.01))
5model.add(Dense(dim_out,activation='linear'))
6model.compile(loss = 'mean_squared_error', optimizer = 'rmsprop')
7history = model.fit(X_train, Y_train, epochs=200, batch_size=10, validation_split=0)
8# Epoch 1/200
9# 1350/1350 [==============================] - 1s 1ms/step - loss: 0.0447
10# Epoch 2/200
11# 1350/1350 [==============================] - 1s 737us/step - loss: 0.0059
12# Epoch 3/200
13# 1350/1350 [==============================] - 1s 743us/step - loss: 0.0043
14# ......
15# Epoch 200/200
16# 1350/1350 [==============================] - 1s 821us/step - loss: 9.2794e-04
如上述代碼所示，我們使用rmsprop算法來優化模型。由於當前的建模場景是數值預測，因此使用MSE（均方誤差）來定義損失函數。算法經過200次迭代，loss從0.0447降到了9.2794e-04。我們可以基於得到的模型進行進一步預測。
5、預測實現
基於上文得到的模型，進一步編寫Python代碼，對X_test對應的輸出數據進行預測。需要注意的是，直接得到的預測結果是處於標準化的數據空間中的，需要將其還原成原始數據空間的值，結果才有意義。對應的Python代碼如下：
1preddf=model.predict(X_test)*vstd.values+vmean.values
如上述代碼所示，將模型的預測結果pred_y乘以vstd再加上vmean，即可對數據進行還原。preddf即是最終得到的預測數據，可列印其值，代碼如下：
1preddf
2# array([[1069.35781887, 1038.57915742, 1056.77147186, 1053.83827734],
3# [1070.65142282, 1039.58533719, 1057.34561875, 1054.85567074],
4# [1083.58529328, 1052.70457308, 1070.78824637, 1067.49741882],
5#
6# [1186.19297789, 1161.52758381, 1172.33666591, 1170.44623263],
7# [1181.42680223, 1155.14778501, 1166.5726204 , 1165.00336968],
8# [1186.75600881, 1160.84733425, 1172.37636963, 1170.09819923]])
9
10preddf.shape
11# (30, 4)
如上述代碼所示，preddf是一個的二維數據，包含了2019年9月整月的預測結果。
6、效果評估
對貴州茅臺股票數據預測的效果評估可以採用兩種方法。一種方法是對預測的結果與真實結果進行繪圖比較，通過直觀觀察可以知道預測效果，如果預測曲線與真實曲線完全重合或相當接近，則說明預測效果較好；反之，則說明預測模型還需要改進。另一種方法是基於貴州茅臺股票數據預測的誤差累計值來計算一個誤差率，從而得到平均精度水平，該值越大說明整體預測效果也就越好，該值越小說明預測模型還存在優化空間。編寫Python代碼，同時實現預測結果與真實數據的對比圖，以及計算累計誤差，從而全面地評估預測效果，代碼如下：
10v = 100*(1 - np.sum(np.abs(preddf - data.iloc[-30:data.shape[0],:4]).values)/np.sum (data.iloc[-30:data.shape[0],: 4].values))
12# Evaluation on test data: accuracy = 99.01%
預測評估對比圖如下。
我們可以看到，黑色實線為真實數據，灰色虛線為預測數據，橫坐標為日期下標，縱坐標為對應的股票價格。使用LSTM模型進行預測的效果總體還是不錯的，平均準確率為99.01%。對於多元時間序列數據，可嘗試使用LSTM模型，該模型能夠記憶歷史較長的重要信息，可有效識別歷史數據中存在的規律和模式，如今廣泛應用於包含大量序列數據的場景中。如果大家對Python感興趣的話，可以加一下我的薇信哦：abb436574，免費領取一套學習資料和視頻課程喲~
本文節選自《Python預測之美：數據分析與算法實戰》一書。《Python預測之美：數據分析與算法實戰（雙色）》，作者遊皓麟，以Python語言為基礎，體系化介紹預測技術工程實施的必備技能。基於Python 來做預測，不僅能夠在業務上快速落地，還讓代碼維護起來更加方便。對預測原理的深度剖析和算法的細緻解讀，是本書的一大亮點。本書共分為預測基礎、預測算法、預測案例三部分。希望讀者在看完本書後，能夠將本書的精要融會貫通，進一步在工作和學習實踐中提煉價值。
#歡迎來留言#
你學會了嗎？
對此，你怎麼看？
程序人生攜手【電子工業出版社-博文視點】送出
《Python預測之美：數據分析與算法實戰》一本
截至7月17日12:00點