• 沒有找到結果。

應用雙季節指數平滑模型於來店人數預測之研究 - 政大學術集成

N/A
N/A
Protected

Academic year: 2021

Share "應用雙季節指數平滑模型於來店人數預測之研究 - 政大學術集成"

Copied!
38
0
0

加載中.... (立即查看全文)

全文

(1)國立政治大學 商學院 統計系 碩士 學位論文. ‧ 國. 學. 政 治 大 應用雙季節指數平滑模型於來店人數 立 預測之研究 ‧. Applications of double seasonal exponential smoothing to store traffic forecasting n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. 研究生:許展源 撰 指導教授:翁久幸 博士 中 華 民 國 107 年 6 月 DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(2) 中文摘要. 政 治 大. Holt-Winters 法是一種可以同時考慮線性趨勢以及季節性的指數平滑法,. 立. 對於單一週期性時間序列資料有不錯的預測效果。Taylor [12] 提出的雙季節. ‧ 國. 學. 指數平滑法比 Holt-Winters 法多了一個季節影響,適合用在具有兩種週期性 的資料。另外,指數平滑法雖然簡易好用,但是並無機率模型。Hyndman et. ‧. al. [6] 將指數平滑法表示成具有單一誤差來源的狀態空間模型。有了狀態. y. Nat. 空間模型表示式, 便能寫出概似函數,進行參數估計及區間估計, 而且可以. io. sit. 很自然地新增外生變數於此模型中。然而,對於單一誤差的狀態空間模型,. n. al. er. 目前文獻上並無討論如何以類似卡爾曼濾波器的方式,來進行狀態的更新。. Ch. i n U. v. 本研究的主要貢獻有兩點。首先是關於指數平滑法在來店人數預測的比. engchi. 較。本論文使用美國服飾業的來店人數資料,我們發現該筆資料具有雙季 節性,使用 Taylor [12] 的雙季節指數平滑法,在預測上明顯優於單季節指 數平滑法。第二點貢獻是針對單一誤差的狀態空間模型,目前文獻上的更 新預測步驟,仍是延續原本的指數平滑法,本論文試著提出類似卡爾曼濾 波器的迭代更新法,來進行狀態更新。有了這種更新方法後,處理外生變 數就容易許多。. 關鍵字:指數平滑、狀態空間模型、雙季節、卡爾曼濾波. i. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(3) Abstract. The Holt-Winters method is an exponential smoothing method that can simultaneously consider linear trend and seasonality, and has a good predictive performance on single periodic time series data. The double seasonal exponential smoothing method proposed by Taylor [12] has a more seasonal effect than the Holt-Winters method and is suitable for data with double seasonality. In addition, although the exponential smoothing method is simple and easy to use, there is no probability model. Hyndman et al. [6] expresses the exponential smoothing method as a state space model with a single source of error. With the state space model representation, we can obtain likelihood function, parameter estimation and interval estimation, and naturally add exogenous variables to this model. However, for a state space model with a single source of error, there is no discussion in the literature on how to update the state in a manner similar to the Kalman filter. There are two main contributions to this study. The first is a comparison of predictive performance of the exponential smoothing method in the number of visitors. This paper uses the number of visitors from the US apparel industry. We found that the data has a double seasonality. Using the double seasonal exponential smoothing method proposed by Taylor [12], it is significantly better than the single seasonal exponential smoothing method in prediction. The second contribution is a state space model for a single source of error. The current update prediction step in the literature is still the traditional exponential smoothing method. This paper attempts to propose an iterative update method similar to the Kalman filter to update the state. With this update method, it is much easier to handle exogenous variables.. 立. 政 治 大. ‧. ‧ 國. 學. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. ii. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(4) 目錄 中文摘要.................................................................................................................................. i. Abstract .................................................................................................................................... ii. 政 治 大 表目錄..................................................................................................................................... 立 圖目錄...................................................................................................................................... iv. 目錄......................................................................................................................................... iii. v. ‧ 國. 學. 第一章 緒論 ...................................................................................................................... 1. ‧. 第二章 文獻回顧............................................................................................................. 3. y. Nat. 第三章 研究方法............................................................................................................. 5 第二節. 雙季節指數平滑法................................................................................. 第三節. 卡爾曼濾波............................................................................................. 第四節. 創新狀態空間模型................................................................................ 10. n. al. er. sit. Holt-Winters 法....................................................................................... io. 第一節. Ch. engchi. i n U. v. 5 6 7. 第四章 創新狀態空間模型的估計 ........................................................................... 14 第一節. 平滑係數選取........................................................................................ 14. 第二節. 系統狀態更新........................................................................................ 15. 第三節. 模擬........................................................................................................ 20. 第五章 實證研究............................................................................................................. 21 第一節. 資料介紹與分析步驟............................................................................ 21. 第二節. 模型分析................................................................................................ 23. 第六章 結論與建議 ........................................................................................................ 29 參考文獻................................................................................................................................. 31. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(5) 表目錄 3-1 加法誤差項模型 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. 政 治 大 5-1(b)切割後結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 立. 5-1(a)原始結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. 23 23 24. ‧. ‧ 國. 學. 5-2 實驗結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. 12. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(6) 圖目錄 3-1 卡爾曼濾波 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. 政 治 大. 4-1 模擬結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. 立. 5-1 來店人數 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. 學. ‧ 國. 9 20 21. 5-2 Holt-Winters vs. DSHW 1-250 . . . . . . . . . . . . . . . . . . . . . . . . . .. 25. 5-3 Holt-Winters vs. DSHW 900-1274 . . . . . . . . . . . . . . . . . . . . . . . .. 25. 5-5 殘差圖 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. y. 27. sit. ‧. 26. 28. Nat. 5-4 Holt-Winters vs. DSHW 1200-1274 . . . . . . . . . . . . . . . . . . . . . . .. io. n. al. er. 5-6 殘差 ACF . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. i n U. v. 6-1 多季節循環 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .. Ch. engchi. 30. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(7) 第一章 緒論 立. 政 治 大. 預測未來一直是人們感興趣的議題,對具有實體店面的零售業來說,若能事先預測. ‧ 國. 學. 當天的客人數量,便能適當配置人力; 若預測時間拉長至七天甚至一個月,便能有效控 制進貨成本、制定優惠方案。隨科技快速發展,網路已是生活不可或缺的東西,營業. ‧. 模式也從實體店面轉為電子商務,實體店家面臨極大的挑戰。儘管網路行銷擁有低成. sit. y. Nat. 本、便利性等優勢,實體店面仍保有試穿以及門市人員專業意見的優勢,若能事先預. io. er. 測來店人數,便能有效控制庫存量與人力配置。即便是非實體店家,許多也設置有線 上客服或電話客服系統,如何適當配置人力也需要良好的預測。. n. al. Ch. i n U. v. 本論文討論並比較若干指數平滑法,在來店人數預測上的表現。我們使用的資料是. engchi. 施佩吟 [15] 中的服飾業店家來店資料,該資料記錄 2007 年 1 月 1 日至 2007 年 12 月 30 日的營業日 (364 天) 中,每小時的來店人數,詳細的資料說明請見第五章。這筆資料具 有雙季節性,除了每日內有週期性,每週七天也有週期性。施佩吟 [15] 依照每日的特 定時段,將原始資料切割成多組時間序列,例如挑出所有中午 12 時的資料,得到一組 含有 364 個資料點的時間序列,接著考慮這組單季節性的時間序列資料,進行分析與 預測,其他時段如下午 1 時,下午 2 時等,也是如此處理。這樣作法的好處是不需要處 理雙季節性的問題,缺點是沒有善加利用每日內的週期性,在預測上仍有改善的空間。 雖然指數平滑法廣為人知,但指數平滑法只是一種方法而非模型,Hyndman et al. [6] 將指數平滑法表示成具有單一誤差來源的狀態空間模型 (見第三章第四節)。有了模型表 示式便能寫出概似函數,進行參數估計及區間估計等。而且可以很自然地新增外生變. 1. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(8) 數於此模型中。然而,對於單一誤差的狀態空間模型,目前文獻上並無討論如何以類 似卡爾曼濾波器的方式,來進行狀態的更新。普遍看到的預測步驟,仍是原本的指數 平滑法。 本研究的主要貢獻有兩點。首先是關於指數平滑法在來店人數預測的比較。因該筆 資料具有雙季節性,使用 Taylor [12] 的雙季節指數平滑法,在預測上明顯優於單季節指 數平滑法 (Holt-Winters)。我們也與施佩吟 [15] 擷取每日特定時段的作法進行比較,發 現擷取後的時間序列,因為沒有考慮每日內的週期性,在預測表現上確實比雙季節指 數平滑法差。第二點貢獻是針對單一誤差的狀態空間模型,本研究試著提出類似卡爾 曼濾波器的迭代更新。這裡考慮的是不具任何趨勢或季節的簡單指數平滑法加上一個. 政 治 大. 外生變數的模型,我們推導出時間更新 (time update) 與測量更新 (measurement update),. 立. 並以模擬實驗評估推導結果的正確性。. ‧ 國. 學. 本論文分成以下流程: 第一章緒論介紹研究動機與目的;第二章是文獻回顧;第三 章研究方法介紹兩種指數平滑法、卡爾曼濾波以及指數平滑的狀態空間模型;第四章. ‧. 說明如何使用卡爾曼濾波來對指數平滑法對應的狀態空間模型做狀態更新,並與傳統. sit. y. Nat. 的指數平滑模型做比較,最後以模擬的方式來驗證推導的正確性;第五章實證研究介. io. al. n. 果提出結論。. er. 紹本研究採用的資料,利用該筆資料進行預測與評估;第六章結論與建議針對研究結. Ch. engchi. i n U. v. 2. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(9) 第二章 文獻回顧 立. 政 治 大. 指數平滑是一種預測的方法,不少有效的預測方法都能看見指數平滑的身影。這些. ‧ 國. 學. 方法都具有將過去的觀測值做加權組合的特性。在設定權重的部分,新的觀測值會比 舊的觀測值來的多。因此“指數平滑”可以反映隨著觀測時間增加,權重呈指數下降. ‧. 的事實。最早的指數平滑法可追溯到 1944 年,由美國海軍的分析師 Brown [2] 所發現。. sit. y. Nat. 同一時期,Holt [4] 也在美國海軍研究辦公室 (ONR) 發展指數平滑法,但 Holt 沒有公開. io. er. 發表該方法,只記錄在 ONR 備忘錄,直到 2004 年才發表。Holt 關於加法性和乘法性 的季節指數平滑則是由他的學生 Winters [14] 提出。因此具有季節性的 Holt’s linear 法也. al. n. 被稱為 Holt-Winters 法。. Ch. engchi. i n U. v. Holt-Winters 法是針對單一週期時間序列的預測方法,這個方法在許多實際應用上 有很好的預測表現。然而,隨著科技進步,我們能夠記錄大量且多樣的資訊,單位尺 度也更小。如果以小時資料來看,每 24 小時會循環一次,又因為 7 日為一週,經過 168 小時又會出現一個循環,對電力需求來說,冬季夏季的使用量大不同,一年 8760 小時又能分出一個循環。對於這種雙週期,甚至多週期的資料,Holt-Winters 法有其 限制。Taylor [12] 為了解決短期電力預測的問題,在 Holt-Winters 法中加入第二個季節 項,提出雙季節指數平滑法,並與季節性 ARIMA 模型作比較。 雙季節指數平滑法在預測短期電力需求上有不錯的表現,如 Ismail, Zahran & ElMetaal [7],使用雙季節整合移動平均自我迴歸模型 (DSARIMA) 與雙季節指數平滑模 型,預測埃及的電力需求並比較兩種方法的效果。Souza, Barros & Miranda [10] 一樣使. 3. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(10) 用雙季節指數平滑法預測電力需求,並考慮假日以及溫度的影響。另一個例子是 Au, Ma & Yeung [1],預測基地台每小時的手機網路流量,建立 DSARIMA 的自動化預測流 程,與雙季節指數平滑模型進行比較。. 立. 政 治 大. ‧. ‧ 國. 學. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. 4. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(11) 第三章 研究方法. Holt-Winters 法. ‧ 國. 學. 第一節. 立. 政 治 大. Holt-Winters 法同時考慮了三個因子,分別為水平 (Level)、成長 (Growth) 與季節. ‧. (Seasonal),能夠更有效的預測具有季節性的時間序列資料。. ( Forecast). y. sit. al. n. (Seasonal ). io. ( Growth). ℓt = α(yt − st−m ) + (1 − α)(ℓt−1 + bt−1 ). er. ( Level ). Nat. 本論文採用具有加法趨勢項及加法季節項的 Holt-Winters 法,定義如下:. bt = β(ℓt − ℓt−1 ) + (1 − β)bt−1. v ni. Csth= γ(yt − ℓ −Ub ) + (1 − γ)st−m e n g c th−1i t−1. yˆ t+h|t = ℓt + hbt + st−m+h+m. (3.1.1) (3.1.2) (3.1.3) (3.1.4). 其中 α 、β 、γ 分別為水平影響、趨勢影響、季節影響的平滑參數,(3.1.3) 式中的 m 代表一個季節週期的長度,平滑參數反映了水平、趨勢、季節項對於資料的適應速 度。yˆ t+h|t 代表向前 h 期預測,其中 h+ m = [( h − 1) mod m ] + 1. 5. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(12) 第二節. 雙季節指數平滑法. 現實生活中的時間序列可能包含不同長度的季節週期,例如一年的用電量,在一天 24 小時的週期中,上下班時間有著固定的週期,一週七天中又有假日與非假日的週期。 若將時間拉長至兩年以上,又會有冬季夏季的週期,當資料存在著多種不同的週期時, Holt-Winters 法在這種情況下很難發揮最大功用。Taylor [12] 為了解決短期電力預測的 問題,在 Holt-Winters 法中加入第二個季節項,提出可以處理雙季節性的模型。 不同於 Holt-Winters 法,雙季節指數平滑法有兩個季節週期,m1 為短週期,m2 為 長週期,則公式可改寫成:. 政 治 大. (1). (2). ℓt = (1 − α)(ℓt−1 + bt−1 ) + α (yt − (st−m1 + st−m2 ) ). 立. (1). (2). 學. ‧ 國. bt = bt − 1 + β ( y t − ℓ t − 1 − bt − 1 − ( s t − m 1 + s t − m 2 ) ) (1). = ( 1 − γd 1 ) s t − m 1 + γd 1 ( y t − ℓ t − 1 − bt − 1 − s t − m 2 ). (2). = ( 1 − γd 2 ) s t − m 2 + γd 2 ( y t − ℓ t − 1 − bt − 1 − s t − m 1 ). st. (2). (2). (1). (1). ‧. st. (1). (2). (3.2.2) (3.2.3) (3.2.4) (3.2.5). Nat. y. yˆ t = ℓt−1 + bt−1 + st−m1 + st−m2. (3.2.1). io. n. al. (2). 及長週期季節項 st. 為水平項及趨勢項的平滑參數。. er. (1). γd2 分別為短週期季節項 st. sit. (3.2.1) (3.2.2) 中的 α 、β 為水平項及趨勢項的平滑參數,(3.2.3)(3.2.4) 中的 γd1 及. i n U. v. 在長周期中有 m2 個季節項,每 m2 單位時間更新一次。在較短的周期內有 m1 個. Ch. engchi. 季節項,每 m1 個單位時間更新一次。在這模型中並沒有要求 m2 與 m1 有倍數關係。 但是,如果 k = m2 /m1 ,則較長周期內有 k 個較短的周期。因此,對於單位時間為小 時的資料,七天 (168 小時) 為一個長週期,在經過 168 小時之後,會有 168 個季節項 (1). (1). (1). (st+1 , st+2 , ..., st+168 ) 會更新一次,一天 (24 小時) 則為短週期,經過 24 小時就會有 24 (2). (2). (2). 個季節項 (st+1 , st+2 , ..., st+24 ) 會更新一次。 在參數估計的部分,對於長週期的 168 個季節項,每一項的平滑參數皆為 γd2 ,對 於短週期的 24 個季節項,每一項的平滑參數皆為 γd1 。. 6. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(13) 第三節. 卡爾曼濾波. Kalman [8] 發表著名的論文,描述了離散數據線性濾波問題的遞歸解決方案。自那 時以來,由於數值計算技術的進步,卡爾曼濾波 (Kalman filter) 一直是廣泛研究和應用 的主題,特別是在自主或輔助導航領域。卡爾曼濾波是一組數學方程,它提供了一種 有效的計算(遞歸)方法來估計過程的狀態,以最小化平方誤差的均值。該濾波器在幾 個方面非常強大:它支持對過去,現在和甚至未來狀態的估計,並且即使未知建模系 統的確切性質,也可以這樣做。 卡爾曼濾波是一種高效率的遞歸濾波器,它能夠從一系列包含雜訊的觀測值中,估. 政 治 大. 計動態系統的狀態。卡爾曼濾波會根據每個觀測值在不同時間下的值,考慮各時間下. 立. 的聯合分布,再對未知變數進行估計,因此會比只以單一觀測值為基礎的估計方式要. ‧ 國. 學. 準。. 卡爾曼濾波的一個典型實例是針對一組有限的、包含噪聲且根據物體位置的觀察序. ‧. 列(可能有偏差),預測出物體的位置、坐標及速度。在很多工程應用(如雷達、電腦. sit. y. Nat. 視覺)中都可以找到它的身影。. io. er. 首先定義卡爾曼濾波的狀態空間模型:. n. a l yt = Hxt + ε t v i n xC t = hAx e tn−1g+cBuht−i 1 +Uwt−1. (3.3.1) (3.3.2). 隨機變數 w 與 ε 分別為程序雜訊 (process noise) 以及量測雜訊 (measurement noise) p(w) ∼ N (0, Q). (3.3.3). p(ε) ∼ N (0, R). (3.3.4). 在 這 邊 假 設 w 與 ε 互 相 獨 立 且 服 從 常 態 分 配,Q 為 程 序 雜 訊 變 異 量 (process noise covariance),R 為量測雜訊變異量 (measurement noise covariance)。H 為觀測模型,將狀 態空間映射到測量空間,在現實中,Q 、R 、H 會隨著更新而變化,此處都先假設為定 值。 首先定義 xˆ t− ∈ Rn 是已知第 t 步之前的狀態估計下,第 t 步的先驗 ( a priori ) 狀態估 計。而 xˆ t ∈ Rn 為已知觀測值 yt 下,第 t 步的後驗 ( a posteriori ) 狀態估計。接著定義. 7. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(14) 先驗以及後驗的估計誤差 (estimate errors): et− = xt − xˆ t− et = xt − xˆ t 並且定義 “先驗測量誤差共變異數 Pt− ” 與 “後驗測量誤差共變異數 Pt ” 如下: Pt− = E[et− et−T ]. (3.3.5). Pt = E[et etT ]. (3.3.6). 卡爾曼濾波是利用反饋控制 ( feedback control ) 進行估計,在時間 t 估計過程狀態, 再從含有 noise 的觀測值獲得反饋。卡爾曼濾波可分成兩個部分,第一是時間更新 ( time update ),第二是測量更新 ( measurement update )。時間更新負責估計當前狀態以及. 政 治 大. 誤差共變異數,以獲得下一時間的先驗估計,也就是更新 xˆ t− 與 Pt− ,測量更新則是負. 立. 學. ‧ 國. 責反饋,也就是更新 xˆ t 與 Pt 。時間更新方程式又稱為預測方程式,測量更新方程式則 為校正方程式。經由一些推導可得: xˆ t = xˆ t− + Kt (yt − H xˆ t− ). (3.3.7). ‧. (3.3.7) 式中的 (yk − H xˆ k− ) 稱為測量革新 (measurement innovation) 或殘差 (residual)。. y. sit. al. er. io. 完全命中。. Nat. 觀察其數學涵義可以發現,其代表的意義為觀測值與預測值的差,殘差為 0 代表預測. n. (3.3.7) 式中的 Kt 為卡爾曼增益 (Kalman gain),定義如下: Pt− H T (3.3.8) Kt = HPt− H T + R 這個 xˆ t 估計子會使後驗估計誤差共變異數最小,詳細推導過程請參考 [Maybeck79,. Ch. engchi. i n U. v. Brown92, Jacobs93] 。從 (3.3.8) 式可以發現,當觀測值誤差 R 越小,K 值就越大。當 Pt− 越小,K 值就越小。 時間更新方程式 (time update equations) xˆ t− = A xˆ t−1 + But−1. (3.3.9). Pt− = APt−1 A T + Q. (3.3.10). 由於卡爾曼濾波是依賴於前一步的估計值 ( xˆ t−1 , Pt−1 ),因此在初始階段需自定一組初 始值 ( xˆ0 , P0 ). 8. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(15) 測量更新方程式 (measurement update equations) Pt− H T Kt = HPt− H T + R. (3.3.11). xˆ t = xˆ t− + Kt (yt − H xˆ t− ). (3.3.12). Pt =( I − Kt H ) Pt−. (3.3.13). 卡爾曼增益是測量更新首先計算的項目,有了卡爾曼增益之後,便能求得最佳化的 後驗估計。測量更新方程式計算完之後,將後驗估計值 ( xˆ k , Pk ) 帶入時間更新方程式,. 政 治 大 圖 3-1 卡爾曼濾波. 求得下一期的先驗估計,形成一個循環。. 立. ‧. ‧ 國. 學. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. 此範例引用自 Greg Welch and Gary Bishop [13] 綠色直線為真實值,黑點為加上白噪音後的觀測值,藍線為卡爾曼濾波模擬的結果, 可以發現其收斂速度非常快. 9. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(16) 第四節. 創新狀態空間模型. Hyndman et al. [6] 在其 書 中, 介 紹指 數 平 滑 法 如 何 表 示 成 狀 態 空 間 模 型,記 為 (E,T,S),其中 E 為誤差項 (Error),T 為趨勢 (Trend),S 為季節性 (Seasonal)。針對誤差 項可以分成加法誤差項的模型以及具有乘法誤差項的模型。舉例來說 ETS(A, N, N) 代 表加法誤差且不含趨勢及季節的指數平滑模型,也就是具有加法誤差項的簡單指數平 滑模型。 Hyndman et al. [6] 將指數平滑法寫成狀態空間模型的型式,公式如下: 觀測方程式:. 政 x治 = Ax 大 + gε. 狀態方程式:. 立. yt = H T xt + ε t t. t −1. (3.4.1) t −1. (3.4.2). 在觀測方程式中,ε t 為誤差項,代表無法被觀測到的部分,也稱為 “innovations”,通常. 學. ‧ 國. 假設其服從變異數為 σ2 的常態分配。狀態向量 xt+1 = (t, bt , st , st−1 , ..., st−m+1 ) T 包含了 趨勢及季節。於狀態方程式中,Axt−1 描述過去的狀態 xt−1 對當前狀態 xt 的影響。向. ‧. 量 g 決定了 innovations 對狀態的附加影響程度,稱為 “persistence vector”。H 與 g 的維. sit. y. Nat. 度皆為 k x 1, A 為 k × k 矩陣。. io. er. Hyndman et al. [6] 將 (3.4.1) (3.4.2) 這種模型稱為 “創新狀態空間模型 (innovations. al. state space model)”。該模型中,觀測方程式與狀態方程式中的誤差項只差了一個乘項 g. n. v i n ,不同於 (3.3.1) 與 (3.3.2) 中 ε tC 與hwt−1 有互相獨立的性質。而由於兩誤差項有線性關 engchi U. 係,以往卡爾曼濾波的更新步驟,在此並不能適用。以下舉兩個例子說明指數平滑法 如何改寫成狀態空間模型。. 10. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(17) 加法誤差項 Holt-Winters 模型: ETS(A,A,A) 透過改寫 (3.1.1) 式至 (3.1.4) 式的 Holt-Winters 法,可以得到具有加法誤差項、加法 趨勢項及加法季節項的 Holt-Winters 模型,稱為 ETS(A,A,A),模型公式如下: 觀測方程式:. y t = ℓ t − 1 + bt − 1 + s t − m + ε t. (3.4.3). 狀態方程式:. ℓt = ℓt−1 + bt−1 + αε t. (3.4.4). bt = bt−1 + βε t. (3.4.5). st = st−m + γε t. (3.4.6). 令 yˆ t = ℓt−1 + bt−1 + st−m 為向前一期預測 (one-step forecast)。接著令 ε t = yt − yˆ t ,也. 政 治 大. 就是向前一期的預測誤差,整理後變成 yt = ℓt−1 + bt−1 + st−m + ε t ,用相同手法也可. 立. 將 (3.1.1) (3.1.2) (3.1.3) 分別整理成 (3.4.5) (3.4.6) (3.4.7),最後再將 (3.4.3) 至 (3.4.6) 用. ‧ 國. yt = H T xt + ε t. ‧. 狀態方程式: xt = Axt−1 + gε t−1    1 1 0 0 ℓt     b  0 1 0 0 t        st  0 0 0 0    x t +1 =  , A =  0 0 1 0  s t −1       s t −2  0 0 0 1  ..   ... ... ... ... . s t − m +1 0 0 0 0. sit. n. al. er. io. ··· ··· ··· ··· ··· ... ···. y. Nat. 其中.   1   1     0  H= 0 ,     0  ...  1. 觀測方程式:. 學. (3.4.1) (3.4.2) 表達如下:. Ch. engchi. i n U. v. 0 0. .  0 0   0 1 , 0 0   0 0 ... ...  1 0.   α    β     γ  g= 0     0  ...  0. 表 3-1 列出部分的加法誤差項模型. 11. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(18) 表 3-1 Seasonal. 加法誤差項模型. None (N). Additive (A). Trend. None (N). ETS(A,N,N). ETS(A,N,A). yˆ t = ℓt−1. yˆ t = ℓt−1 + st−m. ℓt = ℓt−1 + αε t. ℓt = ℓt−1 + αε t st = st−m + γε t. ETS(A,A,N). ETS(A,A,A). yˆ t = ℓt−1 + bt−1 Additive (A). ℓt. yˆ = ℓ 治 政 = ℓ + b + αε 大 ℓ = ℓ t −1. 立 b =b t. t −1. t −1. t. + βε t. t. t −1. + bt − 1 + s t − m. t. t −1. + bt−1 + αε t. bt = bt−1 + βε t. ‧ 國. 學. st = st−m + γε t. ‧. 雙季節指數平滑模型: DSHW(m1 , m2 ). y. Nat. er. io. sit. 如 同 Holt-Winter 法, 本 章 第 二 節 所 介 紹 的 雙 季 節 指 數 平 滑 法, 並 非 隨 機 模 型 (stochastic model),因此無法計算參數估計值的標準誤,也無法更進一步計算預測區間。. n. al. Ch. i n U. v. 透過改寫 (3.2.1) 式至 (3.2.5) 式,可得到雙季節指數平滑模型,記為 DSHW(m1 , m2 ),. engchi. 其中 m1 與 m2 分別為短週期與長週期的長度,公式如下: (1). (2). 觀測方程式:. y t = ℓ t − 1 + bt − 1 + s t − m 1 + s t − m 2 + ε t. (3.4.7). 狀態方程式:. ℓt = ℓt−1 + bt−1 + αε t. (3.4.8). bt = bt−1 + βε t. (3.4.9). (1). = s t − m 1 + γd 1 ε t. (2). = s t − m 2 + γd 2 ε t. st. st (1). (1). (2). (3.4.10) (3.4.11). (2). 令 yˆ t = ℓt−1 + bt−1 + st−m1 + st−m2 為向前一期預測,ε t = yt − yˆ t 為向前一期預測誤 差,(3.2.1) 式至 (3.2.5) 式也能整理成上方 (3.4.7) 式至 (3.4.11) 式。最後再以 (3.4.1) 式及. 12. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(19) (3.4.2) 式表達如下: 觀測方程式:. x t +1. 狀態方程式: xt = Axt−1 + gε t−1   ℓt     b  t   1 1 0 ···     (1)  0 1 0 · · ·  st       (1)  0 0 0 · · ·    s 0 0 1 · · ·  t −1     ...   . . . . =  (1)  , A =  .. .. .. . . s  0 0 0 · · ·  t − m1 +1       s (2)  0 0 0 · · ·  t    0 0 0 · · ·  . . . .  (2)   s t −1   .. .. .. . .   ...   0 0 0 ···   (2) s t − m2 +1. 立. 政 治 大. 0 0 0 ··· 0 0. .  0 0 0 · · · 0 0   0 1 0 · · · 0 0  0 0 0 · · · 0 0  ... ... ... . . . ... ...  , 1 0 0 · · · 0 0   0 0 0 · · · 0 1  0 0 1 · · · 0 0 ... ... ... . . . ... ...   0 0 0 ··· 1 0. . α. .    β       γd 1     0   .  .  g=  .   0     γd 2     0   .   ..  0. 學 ‧. ‧ 國. io. sit. y. Nat. n. al. er. 其中.   1   1     0   0   . H =  ..  , 1     0   0 .  ..  1. yt = H T xt + ε t. Ch. engchi. i n U. v. 13. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(20) 第四章 創新狀態空間模型的估計 政 治 大 立 式這種狀態空間模型的狀態更新,以下第一節中,我們 本章主要討論 (3.4.1)(3.4.2) ‧ 國. 學. 先介紹 Hyndman et al. [6] 如何估計參數。第二節提出類似卡爾曼濾波的估計方法,討論 系統狀態的更新。過往使用的卡爾曼濾波更新步驟,其推導過程有利用到觀測與狀態. ‧. 方程式中的誤差項互相獨立的性質,該性質在創新狀態空間模型不復存在。因此,狀. sit er. io. 第一節. y. Nat. 態更新的式子將有所不同。. 平滑係數選取 a. n. iv l C n hengchi U 關於平滑係數 (α, β, γ) 初始值的選取有幾種方法,第一種是參考過去的經驗,第. 二種是試誤法。針對創新狀態空間模型 (3.4.1)(3.4.2) 有第三種方法,利用 MLE 來求得 MSE 最小的初始值。Hyndman et al. [6] 給出一般性的指數平滑模型: y t = w ( x t −1 ) + r ( x t −1 ) ε t , x t = f ( x t −1 ) + g ( x t −1 ) ε t −1 , 有了指數平滑模型便能求得概似函數 (Likelihood function): n. n. t =1. t =1. L∗ (θ, x0 ) = nlog( ∑ ε2t ) + 2 ∑ log|r ( xt−1 )| 透過最小化 L∗ 可以求得平滑係數 θ = (α, β, γ) T 與初始狀態 x0 = (ℓ0 , b0 , s0 , s−1 , ..., s−m+1 ) T 。 這邊我們令 w( xt−1 ) = HtT xt−1 、r ( xt−1 ) = 1 、 f ( xt−1 ) = Axt−1 ,可得到指數平滑模. 14. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(21) 型的矩陣表達式: yt = HtT xt + ε t , xt = Axt−1 + gεt−1 , 有了狀態空間模型的矩陣表達式,就能使用卡爾曼濾波進行狀態更新。. 第二節. 系統狀態更新. 本節討論如何使用卡爾曼濾波對指數平滑模型進行狀態更新,針對 (3.4.1)(3.4.2) 的 一個特例。先回顧卡爾曼濾波. 政 治 大. yt = HtT xt + ε t ,. 立. xt = Axt−1 + ωt−1 ,. ‧. ‧ 國. 學. 其中 wt 與 ε t 互相獨立,這邊假設兩者的誤差項存在線性關係,兩者之間的關係為 [ ] α ωt = ε t = gε t 0 這邊我們考慮簡單指數平滑加外生變數 pt 的模型,模型如下:. y. yt = HtT xt + ε t ,. (4.2.1). sit. Nat. 觀測方程式:. n. al. er. io. 狀態方程式: xt = Axt−1 + gε t−1 , (4.2.2) [ ] [ ] [ ] ℓ t −1 α 1 0 其中 Ht = , A= , xt = , g= zt p t −1 0 0 1 其狀態向量 xt 包含水平影響 ℓt 以及外生變數 pt ,而 zt 是已知數。以下先列出推導結 [ ] 1. Ch. engchi. i n U. v. 果: 時間更新 (Time update): ˆ t −1 , xˆ − t = A xˆ t−1 + ω. 其中. ωˆ t−1 =. [ ] α 0. (yt−1 − HtT−1 xˆ t−1 ). (4.2.3a). Pt− = APt−1 A T + Q + αA[(1 − zt ) xˆ t−1 pˆ t−1 − Pt−1 Ht ]2×1 [1 0]. + αA([(1 − zt ) xˆ t−1 pˆ t−1 − Pt−1 Ht ]2×1 [1 0])T. (4.2.3b). 15. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(22) 測量更新 (Measurement update): Kt = Pt− HtT [( Ht Pt− HtT + R)]−1. (4.2.4a). − xˆ t = xˆ − t + Kt ( yt − Ht xˆ t ). (4.2.4b). Pt = ( I − Kt Ht ) Pt−. (4.2.4c). 在時間更新部分,與 (3.3.9) 比較可以看出,(4.2.3a) 中的 xˆ t− 多了前一期的預測誤差 ωˆ t−1 ,這裡沒有控制變數 ut ,故沒有 But−1 這一項。而因為觀測與狀態方程式中的誤 差項有線性關係的原故,(4.2.3b) 的 Pt− 比 (3.3.10) 多了 2 項。在測量更新的部分沒有任 何變動。 由於只有時間更新受到誤差項有線性關係的影響,以下將推導時間更新中的兩個式 子 xˆ t− 與 Pt− 。. 立. 政 治 大. ‧ 國. 0. y t −1 −. Nat. 0. HtT−1 xˆ t−1. sit. y. xˆ − t = E [ x t | y t −1 ]. [ ] α. ‧. Proof.. xˆ − t = A xˆ t−1 +. [ ] α. 學. Proposition 1: xˆ t−. n. er. io. = E[ Axt−1 + ωt−1 |yt−1 ] a l [α] v i n = A xˆ t−C E[ε |y ] 1 +h e0 n gt−c1 hti−1 U [ ] α = A xˆ t−1 + E[yt−1 − HtT−1 xt−1 |yt−1 ] 0 [ ] [ ] α α = A xˆ t−1 + y t −1 − HtT−1 E[ xt−1 |yt−1 ] 0 0 [ ] [ ] α α = A xˆ t−1 + y t −1 − HtT−1 xˆ t−1 0 0. Proposition 2: Pt− Pt− = APt−1 A T + Q + αA[(1 − zt−1 ) xˆ t−1 pˆ t−1 − Pt−1 Ht−1 ]2×1 [1 0]. + αA([(1 − zt−1 ) xˆ t−1 pˆ t−1 − Pt−1 Ht−1 ]2×1 [1 0])T 16. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(23) Proof. − T Pt− = E[( xt − xˆ − t )( xt − xˆ t ) | yt−1 ]. = E[( Axt−1 + ωt−1 − A xˆ t−1 )( Axt−1 + ωt−1 − A xˆ t−1 )T |yt−1 ] = E[ [ A( xt−1 − xˆ t−1 ) + ωt−1 ][ A( xt−1 − xˆ t−1 ) + ωt−1 ]T |yt−1 ] = APt−1 A T + Q + E[ A( xt−1 − xˆ t−1 )ωtT−1 |yt−1 ] + E[ [ A( xt−1 − xˆ t−1 )ωtT−1 ]T |yt−1 ] = APt−1 A T + Q + αA[(1 − zt−1 ) xˆ t−1 pˆ t−1 − Pt−1 Ht−1 ]2×1 [1 0] + αA([(1 − zt−1 ) xˆ t−1 pˆ t−1 − Pt−1 Ht−1 ]2×1 [1 0])T (by Lemma 4). 立. 政 治 大. 在 Pt− 部分比原先的卡爾曼濾波多出了兩個期望值 E[ A( xt−1 − xˆ t−1 )ωtT−1 |yt−1 ] 與. ‧ 國. 學. E[ [ A( xt−1 − xˆ t−1 )ωtT−1 ] T |yt−1 ], 互 為 轉 置 矩 陣, 因 此 只 要 計 算 其 中 一 個 即 可。 為 了讓 推 導 過 程簡 潔 明 瞭, 會 把 E[( xt−1 − xˆ t−1 )ωtT−1 |yt−1 ] 分成 三 個 Lemma 來推 導,. ‧. E[( xt−1 − xˆ t−1 )ωtT−1 |yt−1 ] 推導部份則為 Lemma 4。. io. sit. y. Nat. Lemma 1:. n. al. er. E[ xt−1 wtT−1 |yt−1 ] =αyt−1 xˆ t−1 − αE[( HtT−1 xt−1 ) xt−1 |yt−1 ] Proof.. Ch. e[nℓ g]c h i. E[ xt−1 wtT−1 |yt−1 ] = E[. t −1. p t −1. i n U. v. [αε t−1 0]|yt−1 ]. =αyt−1 xˆ t−1 − αE[( HtT−1 xt−1 ) xt−1 |yt−1 ]. Lemma 2: xˆ t−1 E[wtT−1 |yt−1 ] =α xˆ t−1 yt−1 − α xˆ t−1 [1 1] xˆ t−1. 17. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(24) Proof.. [. αℓˆ t−1 (yt−1 − ℓˆ t−1 − pˆ t−1 zt ) 0 xˆ t−1 E[wtT−1 |yt−1 ] = α pˆ t−1 (yt−1 − ℓˆ t−1 − pˆ t−1 zt ) 0. ]. =α xˆ t−1 yt−1 − α xˆ t−1 ℓˆ t−1 − α xˆ t−1 pˆ t−1 =α xˆ t−1 yt−1 − α xˆ t−1 [1 1] xˆ t−1. Lemma 3: E[( HtT−1 xt−1 ) xt−1 |yt−1 ] =( xˆ t−1 xˆ tT−1 + Pt−1 ) Ht−1 Proof.. 政 治 大 [ ℓ. ] ℓ t −1 t −1 E[( HtT−1 xt−1 ) xt−1 |yt−1 ] = E[[1 zt−1 ] | y t −1 ] p t −1 p t −1 [ ] ℓ2t−1 + zt−1 ℓt−1 pt−1 = E[ | y t −1 ] ℓt−1 pt−1 + zt−1 p2t−1 [ ] 1 = E[ xt−1 xtT−1 |yt−1 ] z t −1. 立. ][. ‧. ‧ 國. 學. n. al. er. io. sit. y. Nat. =( xˆ t−1 xˆ tT−1 + Pt−1 ) Ht−1. Ch. engchi. i n U. v. 完成上方 3 個 Lemma 推導之後就能求出 E[( xt−1 − xˆ t−1 )ωtT−1 |yt−1 ] Lemma 4: E[( xt−1 − xˆ t−1 )ωtT−1 |yt−1 ] =α[(1 − zt−1 ) pˆ t−1 xˆ t−1 − Pt−1 Ht−1 ]. 18. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(25) Proof. E[( xt−1 − xˆ t−1 )ωtT−1 |yt−1 ] = E[ xt−1 wtT−1 |yt−1 ] − xˆ t−1 E[wtT−1 |yt−1 ]. =αyt−1 xˆ t−1 − αE[( HtT−1 xt−1 ) xt−1 |yt−1 ] (by Lemma 1 and 2) − (α xˆ t−1 yt−1 − α xˆ t−1 [1 1] xˆ t−1 ) =α xˆ t−1 [1 1] xˆ t−1 − αE[( HtT−1 xt−1 ) xt−1 |yt−1 ] =α xˆ t−1 [1 1] xˆ t−1 − α( xˆ t−1 xˆ tT−1 + Pt−1 ) Ht−1 (by Lemma 3) [ ] =α xˆ t−1 [1 1] xˆ t−1 − α( xˆ t−1 ℓˆ t−1 xˆ t−1 pˆ t−1 Ht−1 ) − αPt−1 Ht−1. 政 治 大. =α xˆ t−1 ℓˆ t−1 + α xˆ t−1 pˆ t−1 − α xˆ t−1 ℓˆ t−1 + α xˆ t−1 pˆ t−1 zt−1. 立− αP. t − 1 Ht − 1. ‧ 國. 學. =α xˆ t−1 pˆ t−1 − α xˆ t−1 pˆ t−1 zt−1 − αPt−1 Ht−1. ‧. =α[(1 − zt−1 ) pˆ t−1 xˆ t−1 − Pt−1 Ht−1 ]. y. Nat. n. er. io. al. sit. 將 Lemma 4 帶回 Pt− 完成推導。. Ch. engchi. i n U. v. 19. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(26) 第三節. 模擬. 本 節 以 模 擬 實 驗 來 討 論 上 一 節 的 理 論 推 導 之 正 確 性。 模 擬 資 料 採 用 R 語 言 “smooth” [11] 套件中的 “sim.es” 函數生成 50 筆 ETS(A, N, A) 資料,週期設為 5 ,初始 狀態設為 (ℓ0 , b0 , s0 ) = (0.5, 0.5, 0.5) ,設定誤差項服從常態分配,初始值設為 0,再使 用 (4.2.3)-(4.2.4) 進行估計。模擬結果顯示,(4.2.3)-(4.2.4) 能適當進行狀態更新。 圖 4-1 模擬結果. 立. 政 治 大. ‧. ‧ 國. 學. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. 20. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(27) 第五章 實證研究. ‧ 國. 資料介紹與分析步驟. 學. 第一節. 立. 政 治 大. 本論文使用美國某精品服飾店資料,圖 5-1 紀錄 2007 年 1 月 1 日到 2007 年 12 月 30. ‧. 日每小時的來店人數。因為 4 月 8 日復活節、11 月 22 日感恩節及 12 月 25 日聖誕節三. sit. y. Nat. 日未營業,本研究沿用施佩吟 [15] 的資料處理做法,針對放假當日以及放假隔日,用. io. er. 前五週同一時間資料的平均值進行差補。因為營業時間非 24 小時,且剛開店以及關店 前的人潮不穩定,因此只採用每日的 12 點至 18 點共 7 小時,在雙季節指數平滑模型上. al. n. v i n Ch 設定的雙季節變為小循環為 7,大循環為 小時U *7 天),也就是說經過 7 個時間點就 i e n g49(7 h c 會出現小更新,經過 49 個時間點會出現大更新。. 圖 5-1 來店人數. 21. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(28) 施佩吟 [15] 從每日每時的來店人數資料中,擷取每日同一時段的資料進行分析與預 測。例如時間點為 12 時的資料,這樣的時間序列共有 364 筆資料,透過這樣針對不同 時間點的資料進行切割,所得到的時間序列仍有每週內的週期性,但無每日內的週期 性。本論文不依時間點將資料切割,如此可以保有每日內的週期性資訊,優點請見下 一節。 本研究使用 ETS(A,N,N)、ETS(A,A,N)、ETS(A,A,A)、DSHW(7,49) ((3.4.7)-(3.4.11)) 與包含一個外生變數 day-of-week(DOW) 的 ETS(A,N,N) 模型 ((4.2.1)-(4.2.2)) 進行預測, 並對 Holt-Winters 模型與 DSHW(7, 49) (式 (3.4.7)-式 (3.4.11)) 的殘差做 ACF 圖。 分析工具使用統計程式語言 R ,預測模型採用 “forecast” [5] 套件中的 “dshw” 函數與. 政 治 大. “ets” 函數,以向前一期預測 (one-step ahead forecast) 來進行,從 t=1274 開始預測,得. 立. 到 t=1275 到 t=2548 的預測值,再用 “forecast” 套件中的 “accuracy” 函數進行評估,選. ‧. ‧ 國. 學. 擇 RMSE、MAE、MAPE 作為預測指標。. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. 22. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(29) 第二節. 模型分析. 先針對上一節提到的切割問題作討論。表 5-1(a) 使用未經切割的資料,以向前一期 預測的方式所得到的結果,再將預測結果進行切割。表 5-1(b) 為切割時間序列後的結 果。比較雙季節指數平滑模型 (ds) 的 MAPE 、RMSE 與 MAE ,原始結果明顯優於切割 後結果,Holt-Winters 模型 (hw) 則是部分時段稍微優於切割後結果。 由於每個時段的時間序列都是以日為單位,不用考慮一日內的週期性。從表 5-1(b) 切割後結果可以看出,雙季節指數平滑模型的效果比 Holt-Winters 模型差,在表 5-1(a) 原始結果中,雙季節指數平滑模型的效果明顯優於 Holt-Winters 模型。. 政 治 原始結果 大. 表 5-1(a). hwRMSE. dsRMSE 18.3876. hr13. 17.0110. 16.5195. 22.4340. hr14. 21.2681. 18.1779. hr15. 20.9838. hr16. 19.3669. hr17. 19.7582. hr18. 23.8890. dsMAPE 23.0775. 22.7794. 19.3047. 18.4896. 28.4254. 24.6962. 21.2290. 17.6267. 17.3916. 27.9896. 23.9497. 18.6574. 14.1633. 16.9492. 26.6916. 22.8884. 15.7448. 12.9095. al. 26.9127. 24.3946. 14.2363. 12.6561. 17.4071. 14.2650. io. 75.5364. n. 17.0975. Ch. 18.7278. i n U. i e30.8040 n g c h 26.1964. 表 5-1(b). sit. 53.6892. er. 12.8239. ‧. 40.7785. Nat. hr12. hwMAPE. y. 立 dsMAE. 學. ‧ 國. hwMAE. v. 切割後結果. hwMAE. dsMAE. hwRMSE. dsRMSE. hwMAPE. dsMAPE. hr12. 12.6006. 13.7750. 17.9290. 18.5899. 21.5247. 24.5132. hr13. 18.4042. 24.8268. 25.9697. 34.4116. 19.6669. 27.0415. hr14. 21.9140. 25.2312. 31.8695. 35.5935. 19.4884. 23.7716. hr15. 23.9387. 29.9521. 34.9250. 41.6261. 18.7113. 25.0338. hr16. 25.3270. 32.2892. 37.0469. 44.0801. 18.4891. 23.9149. hr17. 27.0031. 31.2908. 39.6480. 44.3243. 18.2696. 21.8242. hr18. 23.6201. 30.7584. 32.7731. 40.7028. 18.2561. 25.0782. 23. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(30) 從表 5-2 可以看出,ETS(A,N,N) 跟 ETS(A,A,N) 的預測結果差距不大,因為 ETS(A,N,N) 與 ETS(A,A,N) 兩種模型的差別在於有沒有趨勢影響,剛好這筆來店人數資料沒有明顯 的趨勢,可以預期會有這樣的預測結果。ETS(A,N,N)+DOW 中的 DOW 代表外生變數 “day-of-week” (星期一至星期五設為 0,六日設為 1),這個模型我們使用上一節所推導 的時間更新 (time update) 與測量更新 (measurement update) 來進行迭代更新。他比前兩 個模型多了一個外生變數 DOW,在預測上有略為改進。以上這三個模型都沒有考慮季 節性。 而 ETS(A,A,A) 考慮了季節性,預測結果明顯改善,再看到雙季節指數平滑模型的 三個預測指標又明顯下降,可見該筆資料的確具有雙季節影響,適合使用雙季節指數. 政 治 大. 平滑模型來預測。. 立. MAE. MAPE. ETS(A,N,N). 47.03177. 31.91049. 38.11577. ETS(A,A,N). 47.04344. 31.90937. ETS(A,N,N)+DOW. 46.66465. 31.59696. y. ETS(A,A,A). 32.43201. 23.29364. n. Ch. i n U. v. 23.43904 16.81247. engchi. 37.26808. sit. io. al. DSHW (7, 49). 38.12456. er. Nat. RMSE. ‧. 學. ‧ 國. 表 5-2 實驗結果. 26.01653. 16.16967. 在圖5-2 雙季節指數平滑模型部分,除了少部分出現高估之外,在季節循環上有不 錯的預測結果。在 Holt-Winters 模型預測上,出現不少低估的情況。在一開始虛線包住 的區段為平日週三,卻出現假日的人潮,經查詢之後發現該日為 7 月 4 日美國獨立紀念 日。 看到圖5-3,黑色為真實來店人數,在 12 月 30 日附近,因為接近跨年的關係,店家 的人潮快速上升。Holt-Winters 模型在低點經常出現低估的情況,還出現小於 0 的預測 值。雖然雙季節指數平滑模型在高點低點都有高估低估的狀況,但與真實值相距不遠。. 24. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(31) 圖 5-2. 立. Holt-Winters vs. DSHW 1-250. 政 治 大. ‧. ‧ 國. 學. n. er. io. sit. y. Nat. al. Ch. 圖 5-3. i n U. v. eHolt-Winters n g c h i vs. DSHW 900-1274. 25. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(32) 圖 5-4 縮小區間範圍至 1200-1274 ,後期出現震幅放大的情形,Holt-Winters 模型在 這段時間的走勢中有晚一期的現象,人數下降時還保持上升的趨勢,雙季節指數平滑 模型在這方面表現穩定,幾乎能抓到上升下降的時機。 圖 5-4. 立. Holt-Winters vs. DSHW 1200-1274. 政 治 大. ‧. ‧ 國. 學. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. 26. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(33) 圖 5-5 為兩模型的殘差圖,先看到 Holt-Winters 殘差圖,每隔一段時間殘差會異常 變大,在這裡推測可能是因為第二個季節因素導致預測結果不理想,在頭尾有較大 的變動,稍後會與雙季節指數平滑模型一起討論。再看到雙季節指數平滑模型殘差 圖,其他時間的殘差明顯比 Holt-Winters 模型小且穩定,頭尾部分一樣有不規則變動。 Holt-Winters 模型與雙季節指數平滑模型兩種方法的殘差在 t ∈ (300, 800) 之間有較小變 異,但在頭 (1, 250)、尾 (900, 1274) 部分,兩者的殘差就有明顯的起伏。在雙季節指數 平滑模型殘差圖第 237 的殘差突然增加且為正值,表示在該點的預測出現高估的情況, 推測該商店可能發生什麼事。 圖 5-5 殘差圖. 立. 政 治 大. ‧. ‧ 國. 學. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. 27. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(34) 從圖 5-6 Holt-Winters 殘差的 ACF 圖可以看出預測結果還有一個季節影響,週期約 為 49,符合上面討論的長週期長度。DSHW 殘差的 ACF 則落在標準內,但仍能看出具 有些微週期性,或許是平日假日的季節影響。 圖 5-6. 立. 殘差 ACF. 政 治 大. ‧. ‧ 國. 學. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. 28. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(35) 第六章 結論與建議 立. 政 治 大. 在三種傳統的指數平滑模型比較中,Holt-Winters 模型明顯優於簡單指數平滑模型. ‧ 國. 學. 與 Holt’s linear 模型。接著與 Taylor 的雙季節指數平滑模型作比較,雙季節指數平滑更 優於 Holt-Winters 模型,對於年末突然增長的時段依舊保有不錯的預測結果。. ‧. 本論文將卡爾曼濾波器的推導邏輯,應用到觀測與狀態方程式中,誤差項具有線性. sit. y. Nat. 關係的創新狀態空間模型,推導得到一組新的更新方程式。因為創新狀態空間模型是. io. er. 指數平滑法的一種模型表示法,目前文獻中,對於創新狀態的更新步驟,基本上仍沿 用指數平滑法的操作步驟,在處理外生變數的情況較為棘手。有了類似卡爾曼濾波器. n. al. Ch. i n U. v. 的更新方法後,處理外生變數就容易許多。未來只要擴展狀態向量的維度,根據指數. engchi. 平滑法的公式做調整,將相對應的參數放入矩陣 Ht ,就能套用到其他的指數平滑法。 若能套用至雙季節指數平滑模型,加上國定假日或競爭商家的資訊作為外生變數,預 期能得到更好的預測結果。 由於研究資料的時間長度只有一年,無法比較歷年的季節變化,另外從圖 6-1 可以 看出該筆資料可能存在兩個以上季節循環。設定星期一至星期日為一個區段,圖 6-1 中 以虛線作為分段點,可以看出每個區段包含七個尖點,前四個尖點約為 100 人,也就 是星期一至星期四的人數約為 100 人,第五至第七尖點明顯上升,因此一週七天又能 再分成兩個小季節循環。星期六的人數明顯高於其他天,或許能將星期六單獨視為一 個季節循環。. 29. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(36) 立. 政 治 大 圖 6-1 多季節循環. ‧. ‧ 國. 學. n. er. io. sit. y. Nat. al. Ch. engchi. i n U. v. 30. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(37) 參考文獻 [1] S. Tom Au, Guang-Qin Ma, and Shu-Ngai Yeung. Automatic forecasting of double seasonal time series with applications on mobility network traffic prediction. 2011.. 政 治 大. [2] Robert Goodell Brown. Statistical forecasting for inventory control. McGraw-Hill, New. 立. York, 1959.. ‧ 國. 學. [3] Phillip G Gould, Anne B Koehler, J Keith Ord, Ralph D Snyder, Rob J Hyndman, and. ‧. Farshid Vahid-Araghi. Forecasting time series with multiple seasonal patterns. European Journal of Operational Research, 2008.. sit. y. Nat. io. er. [4] Charles C Holt. Forecasting seasonals and trends by exponentially weighted moving averages. Carnegie Institute of Technology, Graduate school of Industrial Administration,. n. al. 1957.. Ch. engchi. i n U. v. [5] Rob J Hyndman and Yeasmin Khandakar. Automatic time series forecasting: the forecast package for R. Journal of Statistical Software, 26(3):1–22, 2008. [6] Rob J. Hyndman, Anne B. Koehler, J. Keith Ord, and Ralph D. Snyder. Forecasting with Exponential Smoothing:The State Space Approach. Springer, 2008. [7] Mohamed A. Ismail, Alyaa R. Zahran, and Eman M. Abd El-Metaal. Forecasting hourly electricity demand in egypt using double seasonal autoregressive integrated moving average model. 2015. [8] R. E. Kalman. A new approach to linear filtering and prediction problems. Journal of Basic Engineering, 82:35–45, 1960. 31. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(38) [9] R. E. Kalman and R. S. Bucy. New results in linear filtering and new results in linear filtering and prediction theory. Journal of Basic Engineering, 83:95–108, 1961. [10] Reinaldo Castro Souza, Mônica Barros, and Cristina Vidigal C. de Miranda. Short term load forecasting using double seasonal exponential smoothing and interventions to account for holidays and temperature effects. 2007. [11] Ivan Svetunkov. smooth: Forecasting Using Smoothing Functions, 2018. [12] J. W. Taylor. Short-term electricity demand forecasting using double seasonal exponential smoothing. The Journal of the Operational Research Society, 2013.. 政 治 大 [13] Greg Welch and Gary Bishop. An introduction to the kalman filter. 立. Technical report,. ‧ 國. 學. Department of Computer Science University of North Carolina at Chapel Hill, 2006. [14] Peter R Winters. Forecasting sales by exponentially weighted moving averages. Manage-. ‧. ment Science, 1960.. y. Nat. n. al. er. io. sit. [15] 施佩吟. 指數平滑模型應用於來店人數預測之研究. 碩, 國立政治大學, 2015.. Ch. engchi. i n U. v. 32. DOI:10.6814/THE.NCCU.STAT.008.2018.B03.

(39)

參考文獻

相關文件

Then, we recast the signal recovery problem as a smoothing penalized least squares optimization problem, and apply the nonlinear conjugate gradient method to solve the smoothing

In this paper, we build a new class of neural networks based on the smoothing method for NCP introduced by Haddou and Maheux [18] using some family F of smoothing functions.

According to the authors’ earlier experience on symmetric cone optimization, we believe that spectral decomposition associated with cones, nonsmooth analysis regarding

For the proposed algorithm, we establish a global convergence estimate in terms of the objective value, and moreover present a dual application to the standard SCLP, which leads to

Abstract In this paper, we consider the smoothing Newton method for solving a type of absolute value equations associated with second order cone (SOCAVE for short), which.. 1

For the proposed algorithm, we establish its convergence properties, and also present a dual application to the SCLP, leading to an exponential multiplier method which is shown

Lin, A smoothing Newton method based on the generalized Fischer-Burmeister function for MCPs, Nonlinear Analysis: Theory, Methods and Applications, 72(2010), 3739-3758..

The superlinear convergence of Broyden’s method for Example 1 is demonstrated in the following table, and the computed solutions are less accurate than those computed by