IBM SPSS ARIMA 演算法使用案例介紹

摘要： ARIMA 演算法概述 時間序列是以規律的時間間隔採集的測量值的有序集合。時間序列分析的主要目的是根據現有的歷史資料來預測未來的資料。 IBM SPSS ARIMA（差分自迴歸移動平均值）模型是一種典型的時間序列模型，可用於預測未來值，提供比指數平滑模型更復雜的趨勢和季節性成分建模方法...

ARIMA 演算法概述

時間序列是以規律的時間間隔採集的測量值的有序集合。時間序列分析的主要目的是根據現有的歷史資料來預測未來的資料。

IBM SPSS ARIMA（差分自迴歸移動平均值）模型是一種典型的時間序列模型，可用於預測未來值，提供比指數平滑模型更復雜的趨勢和季節性成分建模方法，在模型中包含預測變數。包括自迴歸模型（AR），移動平均模型（MA）和自迴歸移動平均模型（ARMA）。

Python Notebook 環境

可以通過下面兩種 notebook 環境來進行預測分析：

本地環境

本地環境需要有 IBM SPSS 演算法 jar 包及其依賴的包（如果沒有的話可以用後面介紹的IBM Watson Studio 環境，除了載入輸入資料不同其它程式碼是一樣的）。

首先安裝 Python 到路徑 C:\Python35，這裡 Python 版本是 3.5.2，選擇自定義安裝，指定安裝路徑，並加入系統 PATH 變數。如果下面三個庫還沒有安裝的話從命令列執行 pip 安裝：

• pip install pyspark
• pip install jupyter
• pip install numpy

然後把下面 IBM SPSS 兩個 jar 包放到指定路徑 C:\Python35\Lib\site-packages\pyspark\jars，這裡 com.ibm.spss.com.ibm.json.JSON4J-1.0.jar 是 IBM SPSS 演算法包依賴的庫，IBM SPSS 演算法包這裡使用的版本號是 5.1。

• spss-assembly-spark2.0-scala2.11-5.1*.jar
• com.ibm.spss.com.ibm.json.JSON4J-1.0.jar

解壓檔案 sparkspsspy-5.1*.zip, 然後把目錄 spss 複製到下面的路徑即可。 C:\Python\Python35\Lib\site-packages

任選工作目錄，從命令列執行下面的命令就會在當前預設瀏覽器中開啟 Jupyter Notebook.

C:\temp\arima\jupyter notebook

如圖 1：

圖 1. Jupyter Home

右鍵點選 New 按鈕，下拉選單中選擇 Python 3，如圖 2：

圖 2. Jupyter Notebook

接下來就可以在 notebook 的 cell 裡面輸入 Python 程式碼了。

IBM Watson Studio

IBM Watson Studio 已經集成了 IBM SPSS 演算法包及依賴包。登入進去後，在 Projects 裡面新建一個 Project，然後在 Assets->Notebooks->New notebook，選擇 run time 的環境去建立 notebook，圖 3 顯示了核心是 Python 3.5 with Spark 2.1 的 notebook.

圖 3. IBM Watson Studio Notebook

使用案例說明

某產品分類公司想要根據過去十年的銷售資料來預測來年的男裝季度銷售，還想知道哪些預測變數（predictor）對目標被預測變數（target）影響比較大，從而做出相應的決策。預測變數包括：郵遞的產品目錄數（mail）和產品目錄的頁數（page）、開通的訂購熱線數目（phone）、印刷廣告投入額（print）以及客戶服務代表人數（service）。

這裡的使用案例的輸入資料 catalog_seasfac.csv 來自檔案 catalog_seasfac.sav ， 可在任何 IBM® SPSS® Modeler 安裝程式的 Demos 目錄中找到。

輸入原始資料的轉換

首先載入資料 catalog_seasfac.csv 為 spark dataframe ，清單 1 顯示了在本地 notebook 環境下的程式碼，csv 檔案在當前 notebook 路徑"C:\temp\arima\"下面。

清單 1. 本地環境載入輸入資料為 Spark Data Frame

from pyspark.sql import SparkSession
spark = SparkSession.builder.getOrCreate()
from pyspark.sql.types import *
from pyspark import SparkContext
from pyspark import SQLContext
sc = SparkContext.getOrCreate();
sqlContext=SQLContext.getOrCreate(sc);

df = sqlContext.read.format('com.databricks.spark.csv') \
.options(header='true', inferschema='true') \
.load('./catalog_seasfac.csv') # this is your csv file
print(df.schema)
df.show(5)

如果是 IBM Watson Studio notebook, 需要把 catalog_seasfac.csv 上傳使用。登入 IBM Watson Studio 後，新建一個 project，然後在 Assets 裡面新建一個 notebook,在 notebook 右上方，點選"Find and add data"圖示上傳資料，然後右鍵選擇已上傳的資料 catalog_seasfac.csv ，選擇 Insert to Code ->Insert SparkSession DataFrame，notebook 於是自動插入程式碼來建立 Spark Dataframe,這裡也可以修改自動插入程式碼加入預先定義資料的型別 df_schema，見清單 2。

清單 2. 在 IBM Watson Studio 中載入輸入資料為 Spark Data Frame

from pyspark.sql import SparkSession
from pyspark.sql.types import *
spark = SparkSession.builder.getOrCreate()
df_schema = StructType([
StructField("date",TimestampType(),True),
StructField("men",DoubleType(),True),
StructField("women",DoubleType(),True),
StructField("jewel",DoubleType(),True),
StructField("mail",DoubleType(),True),
StructField("page",DoubleType(),True),
StructField("phone",DoubleType(),True),
StructField("print",DoubleType(),True),
StructField("service",DoubleType(),True),
StructField("YEAR_",DoubleType(),True),
StructField("MONTH_",DoubleType(),True),
StructField("DATE_",StringType(),True),
StructField("Seasonal_Err_Men",DoubleType(),True),
StructField("Seasonal_AdjSer_Men",DoubleType(),True),
StructField("Seasonal_Factors_Men",DoubleType(),True),
StructField("Seasonal_TrendCycle_Men",DoubleType(),True)])
df = spark.read\
.format('org.apache.spark.sql.execution.datasources.csv.CSVFileFormat')\
.option('header', 'true')\
.schema(df_schema)\
.load(cos.url('catalog_seasfac.csv', '***********************************'))
print(df.schema)
df.show(5)

圖 4 是打印出的資料模式和前 5 行資料，從圖中可以看到輸入資料的所有欄位及其型別，注意 date 這裡是 TimestampType 型別。因為該產品分類公司希望預測來年的男裝季度銷售，而當前原始資料的時間間隔是月，所以需要對原始資料作轉換把時間間隔轉換為季度。

圖 4. Data schema 和資料

清單 3 程式碼引入 TSDP, RTSDP 和 ARIMA 模型相關的包。

清單 3. TSDP, RTSDP 和 ARIMA 相關依賴包

from spss.ml.forecasting.timeseriesdatapreparation import TimeSeriesDataPreparation
from spss.ml.forecasting.traditional.timeseriesforecasting import TimeSeriesForecastingArima
from spss.ml.forecasting.reversetimeseriesdatapreparation import ReverseTimeSeriesDataPreparation
from spss.ml.forecasting.params.predictor import Predictor,ScorePredictor
from spss.ml.forecasting.params.targetorderlist import TargetOrderList
from spss.ml.forecasting.params.tforderlist import TFOrderList
from spss.ml.forecasting.params.temporal import ForecastEs
from spss.ml.common.wrapper import LocalContainerManager

基於 IBM SPSS 當前對演算法的設計,輸入的原始資料必須由 IBM SPSS TSDP (Time Series Data Preparation) 先處理轉換為一個合適的格式供 ARIMA 模型或其它 IBM SPSS 時間序列模型來消費,然後 IBM SPSS RTSDP (Reverse Time Series Data Preparation) 再將 ARIMA 或其它 IBM SPSS 時間序列模型的輸出轉換為同原始資料一樣的格式。

在 IBM SPSS Modeler 和 IBM SPSS Statistics 中，TSDP 和 RTSDP 被嵌入在產品中對使用者是透明的，所以不需要考慮資料轉換的問題。但是在 notebook 中這種通過 API 呼叫來使用 IBM SPSS 時間序列模型的方式必須考慮資料前後的轉換。

清單 4 程式碼即呼叫 TSDP，把輸入資料的 dataframe 轉換為供 ARIMA 消費的格式。這裡可以通過 TSDP 對輸入資料做不同的轉換。

清單 4. 呼叫 TSDP

lcm = LocalContainerManager() 
tsdp = TimeSeriesDataPreparation(lcm). \
setMetricFieldList(["men","mail","page","phone","print","service"]). \
setDateTimeField("date"). \
setInputTimeInterval("MONTH"). \
setOutTimeInterval("QUARTER") 
tsdpOutput=tsdp.transform(df)

LocalContainerManager 用來建立容器（container）的物件, 通常包含資料轉換或建模過程中產生的 Json, PMML 或 StatXML 等與模型相關的檔案。

setMetricFieldList 指定變數欄位，這裡的變數應該是連續變數，包括被預測變數欄位和預測變數欄位。

setDateTimeField 指定時間日期欄位的名字，這裡是 date。

setInputTimeInterval 指定輸入資料的時間間隔是天（DAY），星期（WEEK），月（MONTH），季度（QUARTER）或年（YEAR），通過 setOutTimeInterval 指定轉換後的時間間隔。本文指定的輸入時間間隔是月（MONTH），輸出時間間隔是季度（QUARTER）。這裡的時間間隔根據原始資料的不同，還可以是時間間隔秒（SECOND），分鐘（MINUTE）或小時（HOUR），甚至非標準時間間隔迴圈週期（CYCLIC_PERIOD）。

在選擇時間序列的模型建模之前，通常需要分析一下時間序列曲線。清單 5 的程式碼用來把 TSDP 轉換後的格式再通過 RTSDP 轉換為和輸入原始資料一樣的格式，這樣就可以把歷史每個季度的男裝銷售曲線圖畫出來，便於進一步的分析。

清單 5. 呼叫 RTSDP

rtsdpRawdata = ReverseTimeSeriesDataPreparation(lcm). \
setInputContainerKeys([tsdp.uid])
rtsdpRawdataOutput=rtsdpRawdata.transform(tsdpOutput)
print(rtsdpRawdataOutput.schema)
rtsdpRawdataOutput.show(5)

從圖 5 中可以看到，經過 TSDP 處理之後，輸入資料的時間間隔已經從月變成了季度，這裡對預測變數作了聚集，現在每個季度的值是之前三個月的值的和。

圖 5. TSDP, RTSDP 轉換後的資料

清單 6 程式碼畫出過去十年（1989 - 1999）的男裝季度銷售額曲線圖：

清單 6. 畫出時間序列圖

import matplotlib.pyplot as plt
import matplotlib.pyplot as plt
%matplotlib inline
fig = plt.figure(figsize=(16,8))
x=rtsdpRawdataOutput.select('date').rdd.map(lambda row : row[0]).collect()
y=rtsdpRawdataOutput.select('men').rdd.map(lambda row : row[0]).collect()
plt.plot(x, y,color="blue", linewidth=1.0, linestyle="-", label="men")
plt.xlabel('date')
plt.legend(loc=0)

從圖 6 中可以看到，該時間序列整體呈現上升趨勢，即序列值趨向於隨時間變化而增加，並維持線性趨勢。此序列還有一個明顯的季節模式，即年度高點在第四季度。

圖 6. 時間序列圖

因為該產品分類公司除了想預測來年每個季度的銷售額，還想知道哪些變數對該銷售額影響最大，所以考慮使用 IBM SPSS ARIMA 模型來進行預測分析。

ARIMA 演算法建模

清單 7 的程式碼用來建立 ARIMA 模型。

清單 7. 建立 ARIMA 模型

arima = TimeSeriesForecastingArima(lcm). \
setInputContainerKeys([tsdp.uid]). \
setTargetPredictorList([Predictor(targetList = [["men"]],
predictorIncludeList=[["mail"],["page"],["phone"],["print"],["service"]])]).\
setTargetOrderList([TargetOrderList(
targetList=[["men"]],
nonSeasonal=[1,0,1],
seasonal=[1,0,1],
transType="none"
)]). \
setTFOrderList([TFOrderList(
targetList=[["men"]],
predictorList=[["mail"],["page"],["phone"],["print"],["service"]],
nonSeasonal=[1,0,0],
seasonal=[0,0,0],
transType="none"
)]). \
setOutInputData(True). \
setCalcPI(True)
arima_model=arima.fit(tsdpOutput)

setInputContainerKeys 用來設定 TSDP 轉換後的容器相關資訊。

setTargetPredictorList 用來設定被預測變數（target）和預測變數（predictor）。這裡被預測變數就是要預測的男裝銷售額 men，預測變數包括 mail, page, phone, print 和 service。

setTargetOrderList 用來設定 ARIMA 的階數，針對被預測變數(target)。對於 ARIMA 的階數，所有的值都必須是非負整數。對於自迴歸和移動平均值來說，該值表示最大階數。所有低於該值的正階數都將包括在模型中。例如，如果指定 2，那麼模型包括的階數是 2 和 1。階數依次包括下面三個值,本文中階數的設定為實驗值。

• 自迴歸（Autoregressive, p）：模型中的自迴歸階數。自迴歸階數指定序列中哪些以前的值用於預測當前值。例如，自迴歸階數 1 指定序列中過去一個時限的值用於預測當前值。
• 差分（Difference, d）：指定在估計模型之前應用於序列的差分的階。當趨勢出現時（具有趨勢的序列通常是不穩定的，而 ARIMA 建模時假定是穩定的），差分是必需的並可用於去除這些趨勢的影響。差分階數與序列趨勢度一致，1 即一階差分表示線性趨勢，2 即二階差分表示二次趨勢，0 表示沒有差分。
• 移動平均值（Moving Average, q）：模型中移動平均值階數的值。移動平均值階數指定如何使用與序列以前值均值之間的偏差來預測當前值。例如，移動平均值階數 1 指定在預測序列的當前值時，可考慮與序列（來自過去一個時限中的每一個）均值之間的偏差。

setTFOrderList 用來設定轉換函式的階數,針對預測變數（predictor）。對於轉換函式的階數，所有的值都必須是非負整數。對於分子和分母來說，該值表示最大階數。所有低於該值的正階數都將包括在模型中。階數依次包括下面三個值,本文中階數的設定為實驗值。

• 分子（Numerator）: 轉換函式的分子階數指定選定的獨立（預測變數）序列中哪些先前的值用於預測相依序列的當前值。例如，分子階數 1 指定獨立序列過去一個時限的值以及獨立序列的當前值可用於預測每個相關序列的當前值。
• 差分（Difference）:指定在估計模型之前應用於所選獨立（預測變數）序列的差分的階數。當趨勢出現時，差分是必需的並可用於去除這些趨勢的影響。
• 分母（Denominator）: 轉換函式的分母階數指定如何使用與選定獨立（預測變數）序列的先前值均值之間的偏差來預測相依序列的當前值。例如，分母階數 1 指定當預測每個相關序列的當前值時，可考慮與獨立序列過去一個時限的均值之間的偏差。

nonSeasonal 指非季節性的階數，取值可以參考上面階數的定義。

seasonal 指季節性的階數，取值可以參考上面階數的定義。

注意：季節性自迴歸成分、移動平均值成分和差分成分與其非季節性對應成分起著相同的作用（季節性分子、分母和差分成分與其非季節性對應成分起著相同的作用）。但是對於季節階數，當前的序列值會受到由一個或多個季節週期分隔的以前序列值的影響。例如，對於以月為時間單位的資料（季節週期為 12），季節階數 1 表示當前序列值會受到當前序列之前的 12 個週期內的序列值的影響。因此，對於以月為時間單位資料，將季節階數指定為 1 相當於將非季節階數指定為 12。

setOutInputData 用來將輸入資料包含在預測的輸出中，從而可以比較輸入值和輸出預測值。

setCalcPI 用來計算預測變數的重要性,只有在存在預測變數的時候才有效。

ARIMA 模型解析及預測變數重要性

建模後還需要對模型相關檔案進行解析，清單 8 的程式碼用來輸出容器中所有模型相關的檔案的列表：

清單 8. 輸出模型相關檔案列表

arima_container = lcm.importContainerSeq(arima.uid)
arima_container.entriesNames()

圖 7 為輸出檔案列表的例子，這裡 TSDPOutput.json 是 TSDP 的輸出，datamodel.xml 是資料模型檔案，containers[0]中的 *.es.xml 和 0.xml 是 PMML 檔案，containers[1]中的 StatXML.xml 和 0.xml 都是 StatXML。StatXML.xml 中可以看到構建模型的各種設定及其值，0.xml 可以檢視模型的各個統計量的資訊進而判斷模型的好壞。

圖 7. 模型相關檔案列表

清單 9 中的程式碼輸出 PMML 檔案，資料模型檔案和 StatXML 檔案，其它檔案可仿照例子輸出。

清單 9. 輸出模型相關檔案

containers = arima_container.containers();
print("containers[0] 0.xml:\n"+containers[0].containerEntry('0.xml').stringContent()+"\n")
print("containers[0] datamodel.xml:\n"+containers[0].containerEntry('datamodel.xml').stringContent()+"\n")
print("containers[1] StatXML.xml:\n"+containers[1].containerEntry('StatXML.xml').stringContent()+"\n")
print("containers[1] 0.xml:\n"+containers[1].containerEntry('0.xml').stringContent()+"\n")

從 containers[0]中的 PMML 檔案 0.xml 可以得到各個預測變數的重要性，見圖 8：

圖 8. 預測變數重要性

從資料模型 datamode.xml 中可以得到對應的變數的編碼，１對應的是 mail，3 對應的是 phone，見圖 9。因為這兩個預測變數的 importance 值最大，所以這兩個預測變數對被預測變數 men 的影響最大。

圖 9. 資料模型中的變數

從 containers[1]中的 StatXML 檔案 0.xml 中可以看到模型相關的統計量，可以對模型的好壞做出評估，見圖 10：

圖 10. 模型相關統計量

ARIM 演算法預測

清單 10 的程式碼用來作預測及設定預測引數:

清單 10. 預測及引數設定

transformer = arima_model. \
setTargets(ScorePredictor(targetIncludedList=[["men"]])). \
setForecast(ForecastEs(outForecast = True, forecastSpan = ４,outCI = True))
predictions = transformer.transform(tsdpOutput)

setTargets 用來指定被預測變數 men。

setForecast 用來指定預測相關設定，outForecast 指定是否作預測，forecastSpan 指定預測的時間間隔，這裡時間間隔是季度，所以 4 指定預測來年的 4 個季度，outCI 指定是否計算輸出預測置信區間 LCI (Lower Confidence Interval) 和 UCI (Upper Confidence Interval)。

清單 11 的程式碼用來將預測結果通過 RTSDP 轉換為和輸入原始資料一樣的格式。對於輸出的 dataframe, 為了便於顯示沒有選取所有欄位而通過 select 選取了部分欄位。

清單 11. RTSDP 轉換

rtsdp = ReverseTimeSeriesDataPreparation(lcm). \
setInputContainerKeys([tsdp.uid]). \
setDeriveFutureIndicatorField(True)

score = rtsdp.transform(predictions).select("date","$FutureFlag","men","$TS-men","$TSLCI-men","$TSUCI-men")
score.show(score.count())

setInputContainerKeys 用來設定 TSDP 轉換後的容器相關資訊。

setDeriveFutureIndicatorField 是否輸出未來指示欄位。

圖 11 為預測的輸出，這裡$FutureFlag 列為 0 顯示歷史時間，為 1 顯示未來的時間。$TS-men 列為男裝銷量的預測值，$TSLCI-men 列和$TSUCI-men 列表明瞭預測值的置信區間。

圖 11. 預測輸出

清單 12 的程式碼用來畫出男裝銷量的歷史值和預測值的比較圖：

清單 12. 畫出預測圖

from datetime import datetime,time
x=score.select('date').rdd.map(lambda row : row[0]).collect()
yObserve=score.select('men').rdd.map(lambda row : row[0]).collect()
yPredict=score.select('$TS-men').rdd.map(lambda row : row[0]).collect()

fig = plt.figure(figsize=(16,8))
plt.plot(x, yObserve,color="blue", linewidth=1.0, linestyle="-", label="men")
plt.plot(x, yPredict,color="red", linewidth=1.0, linestyle="-", label="$TS-men")
plt.xlabel('date')
plt.legend(loc=0)

圖 12 展示了男裝季度銷售額 men 和其預測值的比較圖$TS-men，從圖中可以看到 1999 年的整體銷售情況預測較好，1998 年第四季度高峰期後 1999 年上半年再次回覆正常，下半年一直保持平穩的上升趨勢。

圖 12. 預測圖

結束語

本文通過在 notebook 中轉換原始資料，設定模型引數，建立模型，並根據預測變數的重要性（Predictor Importance）找到對預測影響最大的預測變數，做出預測，完整地展示瞭如何通過 Python 呼叫 IBM SPSS ARIMA API 進行建模、預測及分析。

參考資料

• 檢視 IBM SPSS Modeler 關於時間序列演算法的幫助文件 ，瞭解 IBM SPSS 時間序列演算法。
• 參考 IBM SPSS Modeler 中應用 ARIMA 來進行預測分析的教程 ，瞭解如何使用 IBM SPSS Modeler 來進行 ARIMA 模型的預測分析。
• 檢視 IBM SPSS 預測分析演算法在 IBM Watson Studio Notebook 中的幫助文件 ，這裡 Forecasting 包括多種時間序列相關的演算法。
• 檢視IBM Watson Studio首頁，檢視相關內容及登入後使用 Notebook。