【譯】時間序列建模完整教程（R語言）
對於企業時間是最重要的因素，然而絕大多數公司很難跟上時間的腳步。但是隨著技術的發展，出現了很多有效的方法，能夠讓我們預測未來。不要擔心，本文並不會討論時間機器，討論的都是很實用的東西。
本文將要討論關於預測的方法。有一種預測是跟時間相關的，而這種處理與時間相關資料的方法叫做時間序列模型。其中一種處理基於時間的資料的方法是時間序列建模。 顧名思義，它涉及到按時(年、日、小時、分鐘)的資料，以獲得隱藏的見解，從而做出明智的決策。
當我們處理時序序列資料的時候，時間序列模型是非常有用的模型。大多數公司都是基於時間序列資料來分析第二年的銷售量，網站流量，競爭地位和更多的東西。然而很多人並不瞭解的時間序列分析這個領域。 所以，如果你不瞭解時間序列模型。這篇文章將會想你介紹時間序列模型的處理步驟以及它的相關技術。
本文包含的內容如下所示：
- 時間序列模型介紹
- 使用R語言來探索時間序列資料
- 介紹ARMA時間序列模型
- ARIMA時間序列模型的框架與應用
什麼是時間序列建模
讓我們從基礎開始。這包括平穩序列、隨機行走、 Rho 係數、 Stationarity 的 Dickey Fuller 平穩性檢驗。如果這些術語已經嚇到你了，不要擔心——它們很快就會變得清晰起來，我敢打賭你會在我解釋的時候開始享受這個主題。
平穩序列
判斷一個序列是不是平穩序列有三個評判標準：
1、均值，是與時間t 無關的常數。下圖（左）滿足平穩序列的條件，下圖（右）很明顯具有時間依賴。
2、方差，是與時間t 無關的常數。這個特性叫做方差齊性。下圖顯示了什麼是方差對齊，什麼不是方差對齊。（注意右手邊圖中的不同分佈）
3、協方差，只與時期間隔k有關，與時間t 無關的常數。如下圖（右），可以注意到隨著時間的增加，曲線變得越來越近。因此紅色序列的協方差並不是恆定的。
我們為什麼要關心時間序列的“平穩性”？
除非你的時間序列是平穩的，否則不能建立一個時間序列模型。在很多案例中時間平穩條件常常是不滿足的，所以首先要做的就是讓時間序列變得平穩，然後嘗試使用隨機模型預測這個時間序列。有很多方法來平穩資料，比如消除長期趨勢，差分化。
隨機漫步
這是時間序列中最基本的概念。 你可能很瞭解這個概念。 但是，我發現在這個行業裡，很多人把隨機漫步理解為平穩過程。 在這一節中，我會使用一些數學工具，幫助理解這個概念。我們先看一個例子。
例如: 想象一個女孩在一個巨大的棋盤上隨意移動。 在這種情況下，女孩的下一個位置只取決於最後一個位置。
這裡的 
代表這這個時間點隨機干擾項。這個就是女孩在每一個時間點帶來的隨機性。
現在我們遞迴所有x時間點，最後我們將得到下面的等式：
現在，讓我們嘗試驗證一下隨機遊走的平穩性假設：
1、是否均值為常數？
我們知道由於隨機過程的隨機干擾項的期望值為0.到目前為止：E[X(t)] = E[X(0)] = 常數
2、方差是恆定的嗎？
因此，我們推斷，隨機漫步不是一個平穩過程，因為它具有時變方差。另外，如果我們檢查協方差，我們看到協方差也依賴於時間。
Rho係數
我們已經知道一個隨機遊走是一個非平穩的過程。讓我們在方程中引入一個新的係數，看看我們是否能制定一個檢查平穩性的公式。
引入係數: Rho
現在，我們將改變Rho看看我們可不可以讓這個序列變的平穩。這裡我們將從視覺上進行確定，而不是做任何測試來檢查平穩性。
讓我們從一個Rho=0的完全平穩序列開始。這裡是時間序列的圖：
將Rho的值增加到0.5，我們將會得到如下圖：
你可能會注意到，週期變長了，但基本上似乎沒有一個嚴重的違反平穩性假設。現在讓我們採取更極端的情況下ρ= 0.9
我們仍然看到，在一定的時間間隔後，從極端值返回到零。這一系列也不違反非平穩性明顯。現在，讓我們用ρ= 1隨機遊走看看
這顯然是違反平穩條件。是什麼使rho=1變得這麼特殊的呢？，這種情況並不滿足平穩性測試？我們來找找這個數學的原因
公式 
的期望為： 
這個公式很有意義。下一個X(或者時間點t)被拉到Rho*上一個x的值。
例如，如果x(t–1)=1，E[X(T)] = 0.5（Rho=0.5）。現在，如果從零移動到任何方向下一步想要期望為0。唯一可以讓期望變得更大的就是錯誤率。當Rho變成1呢？下一步沒有任何可能下降。
Dickey Fuller平穩性檢驗
這裡學習的最後一個知識點是Dickey Fuller檢驗。。在統計學裡，Dickey-Fuller檢驗是測試一個自迴歸模型是否存在單位根。這裡根據上面Rho係數有一個調整，將公式轉換為Dickey-Fuller檢驗
我們要測試如果Rho–1=0是否差異顯著。如果零假設不成立，我們將得到一個平穩時間序列。 平穩性測試和將一個序列轉換為平穩性序列是時間序列模型中最重要的部分。因此需要記住本節提到的所有概念方便進入下一節。
使用R語言探索時間序列
本節我們將學習如何使用R處理時間序列。這裡我們只是探索時間序列，並不會建立時間序列模型。 本節使用的資料是R中的內建資料：AirPassengers。這個資料集是1949-1960年每個月國際航空的乘客數量的資料。
載入資料集
下面的程式碼將幫助我們載入資料集並且能夠看到一些資料集的巨集觀指標。
> data(AirPassengers) #載入資料 > class(AirPassengers) #檢視AirPassengers資料型別，這裡的"ts"代表的是時間序列資料 [1] "ts" > start(AirPassengers) #檢視開始時間 [1] 1949 1 > end(AirPassengers) #檢視結束時間 [1] 1960 12 > frequency(AirPassengers) #檢視時間概率的頻率，這裡是一年12個月 [1] 12 > summary(AirPassengers) #檢視資料的一些巨集觀指標 Min. 1st Qu. Median Mean 3rd Qu. Max. 104.0 180.0 265.5 280.3 360.5 622.0
檢視詳細資料
> plot(AirPassengers) #繪製出時間序列 > abline(reg=lm(AirPassengers~time(AirPassengers))) #擬合一條直線
> cycle(AirPassengers) Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec 1949 1 2 3 4 5 6 7 8 9 10 11 12 1950 1 2 3 4 5 6 7 8 9 10 11 12 1951 1 2 3 4 5 6 7 8 9 10 11 12 1952 1 2 3 4 5 6 7 8 9 10 11 12 1953 1 2 3 4 5 6 7 8 9 10 11 12 1954 1 2 3 4 5 6 7 8 9 10 11 12 1955 1 2 3 4 5 6 7 8 9 10 11 12 1956 1 2 3 4 5 6 7 8 9 10 11 12 1957 1 2 3 4 5 6 7 8 9 10 11 12 1958 1 2 3 4 5 6 7 8 9 10 11 12 1959 1 2 3 4 5 6 7 8 9 10 11 12 1960 1 2 3 4 5 6 7 8 9 10 11 12 > plot(aggregate(AirPassengers,FUN=mean)) #繪製週期內（每年）的平均值 > boxplot(AirPassengers~cycle(AirPassengers)) #繪製箱形圖，瞭解週期內的資料變化（最大值、最小值、中位數、及上下四分位數。）
重要推論
- 趨勢顯示旅客的數量每年都在增加
- 七八月的均值和方差比其他月份要高很多
- 每個月的平均值並不相同，但是方差差異很小。因此，可以看出具有很強的週期性。一個週期為12個月或更少。
在時間序列模型中，探索資料是最重要的一步——如果沒有這一步，你將不知道這個序列是不是平穩序列。就像這個例子一樣，我們已經知道了很多關於這個模型的很多細節。
ARMA 時間序列建模
ARMA也叫自迴歸移動平均混合模型。ARMA模型經常在時間序列中使用。在ARMA模型中，AR代表自迴歸，MA代表移動平均。如果這些術語聽起來很複雜，不用擔心-下面將會用幾分鐘的時間簡單介紹這些概念。
在開始之前，你首先要記住，AR或者MA並不是應用在非平穩序列上的。 在實際應用中可能會得到一個非平穩序列，你首先要做的就是將這個序列變成平穩序列（通過差分化/轉換），然後選擇可以使用的時間序列模型。
首先，本文將介紹分開介紹兩個模型（AR&MA）。接下來我們看一看這些模型的特點。
自迴歸時間序列模型（Auto-Regressive Time Series Model）
讓我們用下面的例子來理解 AR 模型:一個國家當前的GDP稱x(t)取決於去年的GDP，即x(t-1)。一個國家在一個財政年度的GDP取決於前一年建立的製造廠/服務，以及本年度新成立工廠/服務。但是GDP的主要組成部分是前一部分，因此，我們可以把 GDP 的等式寫成:
這個方程稱為AR公式。公式表示下一個點完全依賴與前面一個點。α是我們尋求的係數，以便最小化誤差函式。注意，x(t-1)同樣依賴 x(t-2)。 因此，波動在未來會逐漸消退。
例如，假設x(t)是一個城市在某一天出售的果汁瓶的數量。 在冬天，很少有供應商購買果汁瓶。突然，在一個特定的日子裡，溫度上升，果汁的需求猛增到1000個。然而，幾天後，氣候又變得寒冷起來。但是，知道人們在炎熱的日子裡已經習慣了喝果汁，在寒冷的日子裡仍有50% 的人還在喝果汁。接下來的幾天，這個比例下降到25% (50%的50%)，然後在數天之後逐漸變小。下面的圖表說明了AR 系列的慣性特性：
移動平均時間序列模型（Moving Average Time Series Model）
讓我們以另一個例子來理解移動平均時間序列模型：一個公司生產了一批型別在市場上很容易買到的包，在競爭激烈的市場。這個包的銷量在很多天都是0，於是，在某一天，這家公司做了一個實驗，生產了一個不同型別的包，這種包在市場上任何地方找不到，因此，他可以一下子賣掉了1000個。這個包的需求特別高，很快庫存快要完了。這天結束了還有100個包沒賣掉。我們把這個誤差成為時間點誤差。接下來的幾天仍有幾個客戶購買這種包。隨著時間的推移，這個袋子已經失去了它的吸引力。下面通過一個簡單的公式來描述這個場景：
如果我們嘗試繪製這個圖表，它會看起來像這樣：
注意到 MA 和 AR 模型的區別了嗎？ 在 MA 模型中，噪聲/衝擊隨著時間的推移而迅速消失。 AR模型對衝擊具有持久的影響。
AR模型與MA模型的不同
AR和MA 模型之間的主要區別是基於不同時間點的時間序列物件之間的相關性。x(t)與x(t-n)的相關性總為0，這直接源於 x(t)和 x(t-n)之間的協方差對於MA模型來說是零。然而，AM模型中x(t)與x(t-1)的相關性隨著時間的推移變得越來越小。不管是否有 AR 模型或 MA 模型，這種差異都可以被利用。
利用ACF和PACF繪圖
一旦我們得到平穩時間序列，我們必須回答兩個主要問題：
- Q1：是AR 還是MA 過程？
- Q2：我們需要使用哪種順序的AR 或MA 過程？
解決這些問題的訣竅可在上一節中找到。你沒注意到嗎？
第一個問題可以用總相關圖（也稱為自相關函式/ACF）來回答。ACF是不同滯後函式之間的總關聯圖。例如，在GDP問題上，時間點T的GDP是x(t)。我們對x(t)與x(t-1)、x(t-2)等的相關性很感興趣。現在讓我們來反思一下上面所學的內容。在移動平均序列的滯後n中，x(t)與x(t-n-1)之間不存在任何相關。因此，總相關圖在第n 次滯後處截止。因此，很容易找到一個MA系列的滯後。對於AR 序列，這種相關性將逐漸下降而沒有任何截止值。那麼，如果我們是AR 系列，我們該怎麼辦？
這是第二個竅門。如果發現每個滯後的部分相關，則在AR 級數的程度之後會被切斷。例如，如果我們有AR系列，如果我們排除第一滯後x(t-1)的影響，則我們的第二滯後x(t-2)與x(t)無關。因此，偏相關函式（PACF）在第一滯後急劇下降。下面是一些例子來澄清你對這個概念的任何懷疑：
上面的藍線顯示出明顯不同於零的值。顯然，上面的圖在第二滯後之後在PACF 曲線上有一個截止，這意味著這主要是AR過程。
顯然，上面的圖在第二滯後之後在ACF 曲線上有一個截止，這意味著這主要是MA過程。
到目前為止，我們已經討論瞭如何識別型別的固定系列使用ACF 和PACF 圖。現在，我將向大家介紹一個構建時間序列模型的全面框架。此外，我們還討論了時間序列建模的實際應用。
ARIMA時間序列模型的框架與應用
到此，本文快速介紹了時間序列模型的基礎概念、使用R探索時間序列和ARMA模型。現在我們將這些零散的東西組織起來，做一件很有趣的事情。
框架概述
下圖的框架展示瞭如何一步一步的“做一個時間序列分析”
前三步我們在前文已經討論了。儘管如此，這裡還是需要簡單說明一下：
第一步：時間序列視覺化
在構建任何型別的時間序列模型之前，分析其趨勢是至關重要的。我們感興趣的細節與系列中的任何趨勢、季節性或隨機行為有關。我們已經在本系列的第二部分中涵蓋了這一部分。
第二步：序列平穩
一旦我們知道了模式、趨勢、週期。我們就可以檢查序列是否平穩。Dicky-Fuller是一種很流行的檢驗方式。在第一部分意見介紹了這種檢驗方式。在這裡還沒有結束！如果發現序列是非平穩序列怎麼辦？
這裡有三種比較常用的技術來讓一個時間序列平穩。
1、消除趨勢：這裡我們簡單的刪除時間序列中的趨勢成分。例如，我的時間序列的方程是：
我們只需刪除括號中的部分，併為其餘部分建立模型。
2、差分：這是常用的去除非平穩性的技術。這裡我們是對序列的差分的結果建立模型而不是真正的序列。例如：
這種差異稱為 AR(I)MA 中的積分部分。 現在，我們有三個引數：
p:AR d:I q:MA
3、季節性：季節性直接被納入ARIMA模型中。下面的應用部分我們再討論這個。
第三步：找到最優引數
使用ACF 和PACF 圖可以找到引數p、d、q。這種方法的一個補充是，如果ACF 和PACF 逐漸減小，則表明我們需要使時間序列平穩並將值引入到“d”。
第四步：建立ARIMA模型
找到了這些引數，我們現在就可以嘗試建立ARIMA模型了。從上一步找到的值可能只是一個近似估計的值，我們需要探索更多(p,d,q）的組合。最小的BIC和AIC的模型引數才是我們要的。我們也可以嘗試一些具有季節性成分的模型。在這裡，在ACF/PACF圖中我們會注意到一些季節性的東西。
第五步：預測
一旦我們有了最後的ARIMA 模型，我們現在就可以預測未來的時間點。如果模型工作正常，我們也可以用視覺化趨勢來驗證。
時間序列模型的應用
這裡我們用前面的例子。使用這個時間序列做預測。我們建議你在進行下一步之前，先觀察這個資料。
下圖是這些年的乘客數的圖。在往下看之前，觀察這個圖。
這裡是我的觀察：
- 乘客有著逐年增加的趨勢。
- 這看起來有季節性，每一個週期不超過12個月。
- 資料的方差逐年增加。
在我們進行平穩性測試之前我們需要解決兩個問題。第一，我們需要消除方差不齊。這裡我們對這個序列做求對數。第二我們需要解決序列的趨勢性。我們通過對時序序列做差分。現在，我們來檢驗最終序列的平穩性。
> adf.test(diff(log(AirPassengers)), alternative="stationary", k=0) Error in adf.test(diff(log(AirPassengers)), alternative = "stationary", : could not find function "adf.test"
驗證時發現沒有安裝”tseries” 包，解決方案：
> install.packages("tseries")
> library(tseries)
測試結果：
Augmented Dickey-Fuller Test data: diff(log(AirPassengers)) Dickey-Fuller = -9.6003, Lag order = 0, p-value = 0.01 alternative hypothesis: stationary Warning message: In adf.test(diff(log(AirPassengers)), alternative = "stationary", : p-value smaller than printed p-value
下一步是找到在ARIMA 模型中使用的正確引數。我們已經知道，d分量是1，因為我們需要1的差值使級數平穩。我們使用相關圖來實現這一點。以下是ACF 系列圖：
acf(log(AirPassengers))
在上面的圖表中你看到了什麼？很顯然ACF下降的十分的慢，這就意味著乘客的數量並不是平穩的。我們在前面已經討論了，我們現狀準備在序列去對數後的差分上做迴歸，而不是直接在序列去對數後的資料熵差分。讓我們看一下差分後的ACF和PACF曲線吧。
> acf(diff(log(AirPassengers)))
pacf(diff(log(AirPassengers)))
顯然，ACF 曲線在第一次滯後之後就被切斷了。因此，我們理解P 的值應該是0，因為ACF是被切斷的曲線。Q 值為1 或2。經過幾次迭代，我們發現（0，1，1）AS（p，d，q）是最小AIC 和BIC 的組合。
讓我們建立一個 ARIMA 模型，並預測未來的10年。 此外，我們將嘗試在 ARIMA 的配方中加入季節性成分。 然後，我們將把預測和訓練資料一起顯示出來。
> (fit <- arima(log(AirPassengers), c(0, 1, 1),seasonal = list(order = c(0, 1, 1), period = 12))) Call: arima(x = log(AirPassengers), order = c(0, 1, 1), seasonal = list(order = c(0, 1, 1), period = 12)) Coefficients: ma1 sma1 -0.4018 -0.5569 s.e. 0.0896 0.0731 sigma^2 estimated as 0.001348: log likelihood = 244.7, aic = -483.4 > pred <- predict(fit, n.ahead = 10*12) > ts.plot(AirPassengers,2.718^pred$pred, log = "y", lty = c(1,3))
原文連結： ofollow,noindex" target="_blank">https://www.analyticsvidhya.com/blog/2015/12/complete-tutorial-time-series-modeling/