jason brownlee：時間序列預測法是一個過程，而獲得良好預測結果的唯一途徑是實踐這個過程。
在本教程中，您将了解如何利用python語言來預測波士頓每月持械搶劫案發生的數量。
本教程所述為您提供了一套處理時間序列預測問題的架構，包括方法步驟和工具，通過實踐，可以用它來解決自己遇到的相關問題。
本教程結束之後，您将了解：
如何核查python環境并準确地定義一個時間序列預測問題。
如何建構一套測試工具鍊，用于評估模型，開發預測原型。以及如何通過時間序列分析工具更好地了解你的問題。
如何開發自回歸積分滑動平均模型（arima），将其儲存到檔案，并在之後加載它對新的時間步驟進行預測。
讓我們開始吧。
波士頓
在本教程中，我們将端到端地來解析一個時間序列預測工程，從下載下傳資料集、定義問題到訓練出最終模型并進行預測。
該工程并不面面俱到，但展示了如何通過系統性地處理時間序列預測問題，來快速獲得好的結果。
我們将通過如下步驟來解析該工程 ：
環境配置.
問題描述
測試工具鍊
持續性.
資料分析
arima 模型
模型驗證
該部分提供一個解決時間序列預測問題的模闆，你可以套用自己的資料集。
本教程假設在已安裝好了scipy及其相關依賴庫的環境下進行，包括：
scipy
numpy
matplotlib
pandas
scikit-learn
statsmodels
我使用的是python2.7。該腳本将幫助你确認你安裝這些庫的版本。
我編寫該教程的所用的開發環境顯示的結果如下：
我們研究的問題是預測美國波士頓每月持械搶劫案的數量。
該資料集提供了從1966年1月到1975年10月在波士頓每月發生的的武裝搶劫案的數量， 或者說，是這十年間的資料。
這些值是一些記錄的觀察計數，總有118個觀察值。
資料集來自于mccleary＆hay（1980）。
您可以了解有關此資料集的更多資訊，并直接從 datamarket 下載下傳。
下載下傳位址：https://datamarket.com/data/set/22ob/monthly-boston-armed-robberies-jan1966-oct1975-deutsch-and-alt-1977#!ds=22ob&display=line
将資料集下載下傳為csv檔案，并将其放在目前工作目錄中，檔案名為“robberies.csv”。
我們必須開發一套測試工具鍊來審查資料和評估候選模型。
這包含兩個步驟：
1.定義好驗證資料集。
2.開發一套模型評估方法。
這個資料集不是目前時間段的資料集。 這意味着我們不能輕輕松松收集更新後的資料來驗證該模型。
是以，我們假設現在是1974年10月，并将分析和模型選擇的資料集裡扣除最後一年的資料。
最後一年的資料将用于驗證生成的最終模型。
下面的代碼會加載該資料集為pandas序列對象并将其切分成兩部分，一個用于模型開發（dataset.csv）和另一個用于驗證（validation.csv）。
運作該示例将會建立兩個檔案并列印出每個檔案中的觀察數。
1. dataset 106, validation 12
這些檔案的具體内容是：
dataset.csv：1966年1月至1974年10月的觀察（106次觀察）
validation.csv：1974年11月至1975年10月的觀察（12次觀察）
驗證資料集大小是原始資料集的10％。
請注意，由于儲存的資料集中沒有标題行，是以，稍後處理這些檔案時，我們也不需要滿足此要求。
模型評估将僅對上一節中準備的dataset.csv中的資料進行處理。
模型評估涉及兩個要素：
1. 評價名額。
2. 測試政策。
這些觀察值是搶劫案的計數。
我們将使用均方根誤差（rmse）的方式來評價預測的好壞。這将給予錯誤的預測更多的權重，并且将具有與原始資料相同的機關。
在rmse進行計算和報告之前，必須反轉對資料的所有轉換，以使不同方法的效果可直接比較。
我們可以使用scikit-learn庫的輔助函數mean_squared_error（）來執行rmse計算，該函數計算出了預期值清單（測試集）和預測清單之間的均方誤差。 然後我們可以取這個值的平方根來作為一個rmse分數。
例如：
将使用步進驗證的方式(walk-forward)來評估候選模型。
這是因為該問題定義需要滾動式預測的模型。給定所有可用資料，這需要逐漸進行預測。
步進驗證的方式将按照如下幾步執行：
資料集的前50％将被用來訓練模型。
剩餘的50％的資料集将被疊代并測試模型。
對于測試資料集中的每個步驟：
訓練該模型
利用該模型預測一次，并經預測結果儲存下來用于之後的驗證
測試資料集的資料作為下一個疊代的訓練集資料
對在測試資料集疊代期間進行的預測結果進行評估，并報告rmse得分。
由于較小的資料規模，我們将允許在每次預測之前，在給定所有可用資料的情況下重新訓練模型。
我們可以使用簡單的numpy和python來編寫測試工具的代碼。
首先，我們可以将資料集直接分成訓練和測試集。 我們總要将加載的資料轉換為float32類型，以防止仍然有一些資料為string或integer資料類型。
接下來，我們可以時間步長周遊測試資料集。訓練資料集存儲在一個python的list對象中，因為我們需要在每次疊代中很容易的附加一些觀察值，numpy數組連接配接有些overkill。
模型進行的預測，按慣例稱為yhat。這是由于結果或觀察被稱為y，而yhat（有上标“^”的y）是用于預測y變量的數學符号。
列印預測和觀察結果，對每個觀察值進行完整性檢查，以防模型存在問題。
在資料分析和模組化陷入困境之前的第一步是建立一個評估原型。
這裡将為利用測試工具進行模型評估和性能測量，分别提供一套模版。通過該性能測量，可以将目前模型和其他更複雜的預測模型進行比較。
時間序列預測的預測原型被稱為天真預測或persistence。
在這裡，先前時間步驟的觀察結果将被用作于下一時間步長的預測。
我們可以将其直接插入上一節中定義的測試工具中。
下面羅列了完整的代碼：
運作測試工具，列印測試資料集中每次疊代的預測值和實際觀察值。
示例到進行列印模型的rmse這一步結束。
在該示例中，我們可以看到持久性模型實作了51.844的rmse。 這意味着，對于每個預測，該模型中平均有51次搶劫值預測是錯誤的。
我們現在有了一個預測方法和性能評估的原型; 現在我們可以開始深入資料。
我們可以使用統計彙總和資料的繪圖，來快速了解該預測問題的資料結構。
在本節中，我們将從四個角度來看資料：
1.摘要統計。
2.線圖。
3.密度圖。
4.盒型-虛線圖。
利用文本編輯器中打開資料dataset.csv檔案或原始的robberies.csv檔案，并檢視資料。
快速檢查表明，資料集裡沒有明顯丢失的觀察資料。 我們可能已經更早的發現了像nan或'？'等值數，如果我們試圖強制序列化資料為浮點值。
摘要統計提供了觀察值的快速預覽。 它可以幫助我們快速了解我們正在使用的東西。
以下示例為計算并列印時間系列的摘要統計資訊。
運作該示例會提供大量的摘要資訊供回顧。
這些統計資料包含的一些資訊為：
觀察次數（計數）與我們的期望相符，這意味着我們将資料正确處理了。
平均值約為173，我們可以想象出我們在這個序列中的級别。
在112次搶劫時，标準偏差（平均值的平均分布）相對較大。
百分位數與标準差一緻表明資料有很大的差異。
count 106.000000 mean 173.103774 std 112.231133 min 29.000000 25% 74.750000 50% 144.500000 75% 271.750000 max 487.000000
對于序列中比較大的變化将可能模型很難進行高度準确的預測，如果它是由随機波動（例如，非系統的）引起，
畫一個時間序列的線圖可以對該問題提供很多深入的資訊。
下面的示例建立并顯示出資料集的線圖。
運作示例并檢視繪圖。 注意系列中的任何明顯的時間結構。
圖中的一些觀察結果包括：
随着時間的推移，劫案有增加的趨勢。
似乎沒有任何明顯的異常值。
每年都有相對較大的波動，上下波動。
後幾年的波動似乎大于前幾年的波動。
整體趨勢表明該資料集幾乎肯定不是穩定的，并且波形的明顯變化也可能對此有所提示。
每月波士頓搶劫案資料線圖
這些簡單的觀察結果表明，我們可以看到對趨勢進行模組化和從時間序列中将它提取出來的好處。
或者，我們可以使用差分法進而使序列的模組化平穩。如果後期波動有增長趨勢，我們甚至可能需要兩個差分到不同級别。
提取出密度圖，可以提供對資料結構的進一步了解。
下面的示例建立了沒有任何時間結構的觀測值的直方圖和密度圖。
運作示例并檢視繪圖。
圖中展示的一些分析結果包括：
分布不是高斯分布。
分布是左移的，可以是指數或雙高斯分布
每月波士頓搶劫密案數量密度圖
我們可以看到在模組化之前進行一些資料的次元變換的一些益處。
我們可以按年份對月度資料進行分組，并了解每年觀察值的分布情況，以及觀察值可能如何變化。
我們确實希望從資料中看到一些趨勢（增加平均值或中位數），但看到其餘的分布是如何變化的可能會更有趣。
下面的示例将觀察值按照年份劃分，并為每個觀察年建立一個盒型-虛線圖。 最後一年（1974年）隻包含10個月，可能與其他有12個月的觀察值的年份相比會無效。 是以，隻有1966年和1973年之間的資料被繪制。
from pandas import series from pandas import dataframe from pandas import timegrouper from matplotlib import pyplot series = series.from_csv('dataset.csv') groups = series['1966':'1973'].groupby(timegrouper('a')) years = dataframe() for name, group in groups: years[name.year] = group.values years.boxplot() pyplot.show()
圖中的一些觀察結果包括:
每年的中間值（紅線）顯示可能不是線性的趨勢
散列值和中間50％的資料（藍色框）不同，但時間長了可能會有變化。
較早的年份，也許是前2年，與資料集的其餘部分有很大的不同。
波士頓每月搶劫案數目的盒形-虛線圖
觀察結果表明，年度波動可能不是系統性的，難以模拟。 他們還建議将模組化的前兩年的資料從模型分離出來可能有一些好處，如果确實不同。
這種對資料進行年度觀察是一個有趣的方法，而且可以通過檢視年度彙總統計資料和其中的變化，來進一步跟進未來的變化。
接下來，我們可以開始檢視該序列的預測模型。
在本節中，我們将會針對該問題開發自回歸積分滑動平均模型，即arima模型。
我們将其分四個步驟：
1.開發一套手動配置的arima模型。
2.使用arima的網格搜尋來找到優化的模型。
3.預測殘差的分析，以評估模型中的所有偏差。
4.使用能量變換探索模型的改進。
nasonasonal arima（p，d，q）需要三個參數，并且通常是傳統地進行手動配置。
時間序列資料的分析。假設我們使用的是靜态時間序列。
時間序列基本上肯定是非靜态的。我們可以通過首先區分序列，并使用統計測試，來確定結果是靜态的。
下面的示例建立一個靜态版本的序列，并将其儲存到檔案stationary.csv。
from statsmodels.tsa.stattools import adfuller # create a differe def difference(dataset): diff = list() for i in range(1, len(dataset)): value = dataset[i] - dataset[i - 1] diff.append(value) return series(diff) x = series.values # difference data stationary = difference(x) stationary.index = series.index[1:] # check if stationary result = adfuller(stationary) print('adf statistic: %f' % result[0]) print('p-value: %f' % result[1]) print('critical values:') for key, value in result[4].items(): print('\t%s: %.3f' % (key, value)) # save stationary.to_csv('stationary.csv')
運作示例輸出統計顯着性測試的結果，以表明序列是否靜止。 具體來說，強化迪基-福勒檢驗。
結果表明，測試的統計值-3.980946小于-2.893的5％時的臨界值。 這表明我們可以否決具有小于5％的顯着性水準的零假設（即，結果為統計結果的機率很低）。
拒絕零假設意味着過程沒有機關根，反過來說就是，該時間序列是固定的或者不具有時間依賴性結構。
adf statistic: -3.980946 p-value: 0.001514 critical values: 5%: -2.893 1%: -3.503 10%: -2.584
這表明至少需要一個級别的差分。 我們的arima模型中的d參數應至少為1的值。
下一步是分别選擇自回歸（ar）和移動平均（ma）參數的滞後值，p和q。
我們可以通過審查自相關函數（acf）和部分自相關函數（pacf）圖來做到這一點。
以下示例為該序列建立acf和pacf圖。
from statsmodels.graphics.tsaplots import plot_acf from statsmodels.graphics.tsaplots import plot_pacf pyplot.figure() pyplot.subplot(211) plot_acf(series, ax=pyplot.gca()) pyplot.subplot(212) plot_pacf(series, ax=pyplot.gca())
運作示例并檢視繪圖，詳細了解一下如何設定arima模型的p和q變量。
下面是從圖中得到的一些觀察結論。
acf顯示1個月的顯著滞後。
pacf顯示了大約2個月的顯著滞後，滞後至大約12個月以後。
acf和pacf在同一點顯示出下降，可能是ar和ma參數的混合。
p和q值的正确的起點是1或2。
每月波士頓搶劫案數量的acf和pacf圖
這個快速分析結果表明arima（1,1,2）對原始資料可能是一個很好的切入點。
雷鋒網了解，實驗表明，這套arima的配置不會收斂，并引起底層庫的錯誤。一些實驗表明，該模型表現不太穩定，同時定義了非零的ar和ma階數。
該模型可以簡化為arima（0,1,2）。 下面的示例示範了此arima模型在測試工具上的性能。
from sklearn.metrics import mean_squared_error from statsmodels.tsa.arima_model import arima from math import sqrt # load data # prepare data x = x.astype('float32') train_size = int(len(x) * 0.50) train, test = x[0:train_size], x[train_size:] # walk-forward validation history = [x for x in train] predictions = list() for i in range(len(test)): # predict model = arima(history, order=(0,1,2)) model_fit = model.fit(disp=0) yhat = model_fit.forecast()[0] predictions.append(yhat) # observation obs = test[i] history.append(obs) print('>predicted=%.3f, expected=%3.f' % (yhat, obs)) # report performance mse = mean_squared_error(test, predictions) rmse = sqrt(mse) print('rmse: %.3f' % rmse)
運作此示例結果的rmse值為49.821，低于持續性模型（persistence）。
... >predicted=280.614, expected=287 >predicted=302.079, expected=355 >predicted=340.210, expected=460 >predicted=405.172, expected=364 >predicted=333.755, expected=487 rmse: 49.821
這是一個好的開始，但是我們可以達到改進結果通過一個配置更好的arima模型。
許多arima配置在此資料集上不穩定，但可能有其他超參數導向一個表現良好的模型。
在本節中，我們将搜尋p，d和q的值，以找出不會導緻錯誤的組合，并找到可獲得最佳性能的組合。 我們将使用網格搜尋來探索整數值子集中的所有組合。
具體來說，我們将搜尋以下參數的所有組合：
p: 0 - 12.
d: 0 - 3.
q: 0 -12.
這是（13 * 4 * 13），或676，測試工具的運作結果，并且将花費需要一些時間用來執行。
下面是利用網格搜尋的測試工具的完整示例：
import warnings # evaluate an arima model for a given order (p,d,q) and return rmse def evaluate_arima_model(x, arima_order): # prepare training dataset # make predictions for t in range(len(test)): model = arima(history, order=arima_order) history.append(test[t]) # calculate out of sample error return rmse # evaluate combinations of p, d and q values for an arima model def evaluate_models(dataset, p_values, d_values, q_values): dataset = dataset.astype('float32') best_score, best_cfg = float("inf"), none for p in p_values: for d in d_values: for q in q_values: order = (p,d,q) try: mse = evaluate_arima_model(dataset, order) if mse < best_score: best_score, best_cfg = mse, order print('arima%s mse=%.3f' % (order,mse)) except: continue print('best arima%s mse=%.3f' % (best_cfg, best_score)) # load dataset # evaluate parameters p_values = range(0,13) d_values = range(0, 4) q_values = range(0, 13) warnings.filterwarnings("ignore") evaluate_models(series.values, p_values, d_values, q_values)
運作該示例并運作所有組合，并将無錯誤收斂的結果報出來。
該示例目前硬體上運作至少需要2小時。
結果表明，最佳配置為arima（0，1，2），巧合的是，這在前面的部分中已經被證明。
arima(6, 1, 0) mse=52.437
arima(6, 2, 0) mse=58.307
arima(7, 1, 0) mse=51.104
arima(7, 1, 1) mse=52.340
arima(8, 1, 0) mse=51.759
best arima(0, 1, 2) mse=49.821
對一個模型來說，一個不錯的最終檢查，是檢查其殘差預測誤差。
理想情況下，殘差的分布應該是符合零均值的高斯分布。
我們可以通過用直方圖和密度圖繪制殘差來檢查這一點。
以下示例為計算測試集上預測的殘差，并建立該密度圖。
# walk-foward validation # errors residuals = [test[i]-predictions[i] for i in range(len(test))] residuals = dataframe(residuals) residuals.hist(ax=pyplot.gca()) residuals.plot(kind='kde', ax=pyplot.gca())
運作該示例會建立兩個圖表。
這些圖表征了一個有長右尾的高斯樣分布。
這可能表明該預測是有偏差的，并且在這種情況下，在模組化之前對原始資料進行基于能量的變換可能是有用的。
殘差密度圖
通過檢查所有時間序列的殘差中是否有某種類型的自相關性也是一個好辦法。 如果存在，它将表明模型有更多的機會來模組化資料中的時間結構。
下面的示例是重新計算殘差，并建立acf和pacf圖以檢查出比較明顯的的自相關性。
plot_acf(residuals, ax=pyplot.gca()) plot_pacf(residuals, ax=pyplot.gca())
結果表明，該時間序列中存在的少量的自相關性已被模型捕獲。
殘差acf和pacf圖
box-cox變換是可以對一組能量變換進行評估的方法，包括但不限于資料的對數，平方根和互逆變換。
下面的示例對資料的執行一個對數變換，并生成一些曲線以檢視對時間序列的影響。
from scipy.stats import boxcox
from statsmodels.graphics.gofplots import qqplot
transformed, lam = boxcox(x)
print('lambda: %f' % lam)
pyplot.figure(1)
# line plot
pyplot.subplot(311)
pyplot.plot(transformed)
# histogram
pyplot.subplot(312)
pyplot.hist(transformed)
# q-q plot
pyplot.subplot(313)
qqplot(transformed, line='r', ax=pyplot.gca())
運作示例建立三個圖：經變換的時間序列的線圖，展示變換後值的分布的直方圖，以及展示變換後的值得分布與理想化高斯分布相比的q-q圖。
這些圖中的一些觀察結果如下：
大的波動已在時間序列的線圖中被移除。
直方圖顯示出了數值更平坦更均勻（表現良好）的分布。
q-q圖是合理的，但仍然不适合高斯分布。
每月波士頓搶劫案資料的boxcox變換圖
毫無疑問，box-cox變換對時間序列做了一些事情，而且是有效的。
在使用轉換資料測試arima模型之前，我們必須有一種方法來進行逆轉轉換，以便可以将對在轉換後的尺度上進行訓練的模型所做的預測轉換回原始尺度。
示例中使用的boxcox（）函數通過優化損失函數找到理想的lambda值。
lambda在以下函數中用于資料轉換：
transform = log(x), if lambda == 0 transform = (x^lambda - 1) / lambda, if lambda != 0
這個變換函數可以直接逆過來，如下：
x = exp(transform) if lambda == 0 x = exp(log(lambda * transform + 1) / lambda)
這個逆box-cox變換函數可以在python中如下實作：
# invert box-cox transform from math import log from math import exp def boxcox_inverse(value, lam): if lam == 0: return exp(value) return exp(log(lam * value + 1) / lam)
我們可以用box-cox變換重新評估arima（0,1,2）模型。
這包括在拟合arima模型之前的初次變換，然後在存儲之前對預測進行反變換，以便稍後與期望值進行比較。
boxcox（）函數可能失敗。 在實踐中，我已經認識到這一點，它似乎是傳回了一個小于-5的lambda值。 按照慣例，lambda值在-5和5之間計算。
對于小于-5的lambda值添加檢查，如果是這種情況，則假定lambda值為1，并且使用原始資料來拟合模型。 lambda值為1與未轉換相同，是以逆變換沒有效果。
下面列出了完整的示例。
# transform transformed, lam = boxcox(history) if lam < -5: transformed, lam = history, 1 model = arima(transformed, order=(0,1,2)) # invert transformed prediction yhat = boxcox_inverse(yhat, lam)
運作此示例将列印每次疊代的預測值和預期值
注意，當使用boxcox（）轉換函數時，可能會看到一些警告; 例如：
runtimewarning: overflow encountered in square
llf -= n / 2.0 * np.log(np.sum((y - y_mean)**2. / n, axis=0))
不過現在可以忽略這些。
在轉換資料上生成的模型最終的rmse為49.443。 對于未轉換的資料，這是一個比arima模型更小的誤差，但隻是很細微的，它可能是，也可能不是統計意義上的不同。
>predicted=276.253, expected=287 >predicted=299.811, expected=355 >predicted=338.997, expected=460 >predicted=404.509, expected=364 >predicted=349.336, expected=487 rmse: 49.443
我們将使用該box-cox變換的模型作為最終模型。
在開發完模型并選擇最終模型後，必須對其進行驗證和确定。
驗證是一個可選的過程，為我們提供了“最後的檢查手段”，以確定我們沒有被自己欺騙。
此部分包括以下步驟：
1. 确定模型：訓練并儲存最終模型。
2. 進行預測: 加載最終模型并進行預測。
3. 驗證模型: 加載和驗證最終模型。
模型的确定涉及在整個資料集上拟合arima模型，值得注意的一點是，該步驟中的資料集已經經過了變換。
一旦拟合好之後，模型可以儲存到檔案以供後續使用。這裡由于我們對資料執行了box-cox變換，是以我們需要明确所選擇的lamda值，以便使得來自模型的任何預測都可以被轉換回原始的未轉換的尺度。
下面的示例展示了在box-cox變換資料集上使用arima（0,1,2）模型，并将整個拟合對象和lambda值儲存到檔案中的過程。
import numpy
# transform data
# fit model
# save model
model_fit.save('model.pkl')
numpy.save('model_lambda.npy', [lam])
運作該示例會建立兩個本地檔案：
model.pkl這是通過arima.fit() 調用生成的arimaresult對象。 包括了各種系數和在拟合模型時傳回的所有其他内部資料。
model_lambda.npy這是作為一行一列numpy數組對象存儲的lambda值。
這些檔案或許會儲存的過于詳細。其實在操作中真正需要的是：來自模型的ar和ma系數、用于差異數量的d參數、以及滞後的觀察值、模型殘差和用于變換的lamda值。
這裡，最直接的方法是加載模型并進行單一預測。
這是相對直接的，涉及恢複儲存的模型、lambda值和調用forecast() 方法。
然而，在目前穩定版本的statsmodels庫（v0.6.1）中存在一個bug，它的錯誤提示如下：
typeerror: __new__() takes at least 3 arguments (1 given)
據我測試發現，這個bug也似乎存在于statsmodels的0.8版本候選1版本中。
更多詳細資訊，請參閱zae myung kim的讨論和在github上的解決辦法。
連結：https://github.com/statsmodels/statsmodels/pull/3217
我們可以使用一個猴子更新檔（monkey patch）來解決這個問題，在儲存之前将一個 __getnewargs __() 執行個體函數添加到arima類中。
下面的示例與上一小節的變化相同。
# monkey patch around bug in arima class def __getnewargs__(self): return ((self.endog),(self.k_lags, self.k_diff, self.k_ma)) arima.__getnewargs__ = __getnewargs__
我們現在可以加載模型并進行單一預測。
下面這個示例展示了加載模型，對下一時間步進進行了預測，反轉box-cox變換，并列印預測結果的過程。
from statsmodels.tsa.arima_model import arimaresults
model_fit = arimaresults.load('model.pkl')
lam = numpy.load('model_lambda.npy')
print('predicted: %.3f' % yhat)
最後得到的預測值為452。
如果我們在validation.csv中檢視的話，我們可以看到下一個時間段的第一行的數值是452。模型達到了100% 的正确率，這是非常驚人的（也是幸運的）。
predicted: 452.043
在這一步驟中，我們可以加載模型并以模拟操作的方式使用它。
在測試工具部分中，我們将原始資料集的最後12個月儲存在單獨的檔案中，以驗證最終模型。
我們現在可以加載此validation.csv檔案，并用它來檢驗我們的模型對真正“未知”的資料究竟表現如何。
這裡可以采用兩種方式：
加載模型并使用它預測未來的12個月。 超過前一兩個月的預測将很快開始降低技能（start to degrade in skill）。
加載模型并以步進預測的方式使用該模型，更新每個時間步長上的變換和模型。 這是首選的方法，在實踐中使用這個模型，它将表現出最佳的性能。
與前面章節中的模型評估一樣，我們将以步進預測的方式進行預測。 這意味着我們将在驗證資料集中逐漸引導時間點，并将觀察結果作為曆史記錄的更新。
# load and prepare datasets dataset = series.from_csv('dataset.csv') x = dataset.values.astype('float32') history = [x for x in x] validation = series.from_csv('validation.csv') y = validation.values.astype('float32') # load model # make first prediction history.append(y[0]) print('>predicted=%.3f, expected=%3.f' % (yhat, y[0])) # rolling forecasts for i in range(1, len(y)): obs = y[i] mse = mean_squared_error(y, predictions) pyplot.plot(y) pyplot.plot(predictions, color='red')
運作以上示例将列印出在驗證資料集中所有時間步長上的每個預測值和預期值。
修正階段最終的rmse預測結果為53次搶劫。 這與預期的49次并沒有太大的出入，我希望它在簡單的持久性模型上也能有類似的表現。
>predicted=452.043, expected=452 >predicted=423.088, expected=391 >predicted=408.378, expected=500 >predicted=482.454, expected=451 >predicted=445.944, expected=375 >predicted=413.881, expected=372 >predicted=413.209, expected=302 >predicted=355.159, expected=316 >predicted=363.515, expected=398 >predicted=406.365, expected=394 >predicted=394.186, expected=431 >predicted=428.174, expected=431 rmse: 53.078
我們在這裡同時提供了預測結果與驗證資料集之間互相比較的圖表。
該預測具有持續性預測的特性。這表明雖然這個時間序列有明顯的趨勢，但它仍然是一個相當困難的問題。
驗證時間步驟的預測
需要說明的是，本教程并不完美，你可以有很多改善結果的方法。
本節列出了一些或許可行的改進措施。
統計顯着性檢驗。使用統計學測試來檢查不同模型之間的結果差異是否具有統計學的意義。t檢驗将是一個很好的入手點。
使用資料變換的網格進行搜尋。通過box-cox變換在arima超參數中重複網格搜尋過程，并檢視是否可以實作一個與此前不同的，并且更好的參數集。
檢查殘差。通過box-cox變換檢測最終模型的預測誤差的殘差，以檢視是否存在可以進一步優化的偏差或自相關。
精簡模型存儲。嘗試簡化模型儲存，僅僅存儲那些必要的系數，而不是整個arimaresults對象。
手動處理趨勢。使用線性或非線性模型直接對趨勢進行模組化，并明确地從序列中删除它。如果該趨勢是非線性的，并且更适合非線性模組化，那這種方式可能會帶來更好的性能表現。
置信區間。顯示驗證資料集中預測的置信區間。
資料選擇。考慮在沒有前兩年資料的情況下對問題進行模組化，并檢視這是否對預測技能（forecast skill）有影響。
通過本教程的學習，你可以初步了解到基于python的時間序列預測項目的具體實作步驟和工具。
我們在本教程中讨論了很多内容，具體來說包括以下三個方面：
如何通過性能名額和評估方法開發測試工具，以及如何快速開發基準預測原型和方法。
如何通過時間序列分析來提出最好的模拟預測問題的方法。
如何開發arima模型，并儲存它，然後加載它用來對新資料進行預測。
本文作者：ai研習社