自己回帰移動平均モデル
自己回帰移動平均モデル（じこかいきいどうへいきんモデル、英: Autoregressive moving average model、ARMAモデル）は、統計学において時系列データに適用されるモデルである。George Box と G. M. Jenkins の名をとって "ボックス・ジェンキンスモデル" とも呼ばれる。
時系列データ Xt について、ARMAモデルはその将来の値を予測するためのツールとして機能する。モデルは自己回帰(AR)部分と移動平均(MA)部分からなる。一般に ARMA(p,q)モデルと表記され、p は自己回帰部分の次数、q は移動平均部分の次数を表す（定義は後述）。
Contents
自己回帰モデル
AR(p) という表記は次数 p の自己回帰モデルを表す。AR(p)モデルは次の式で表される。
[math] X_t = c + \sum_{i=1}^p \varphi_i X_{t-i}+ \varepsilon_t .\,[/math]
ここで[math]\varphi_1, \ldots, \varphi_p[/math]はモデルのパラメータ、[math]c[/math] は定数項、[math]\varepsilon_t[/math] は誤差項（後述）である。定数項は単純化するために省かれることが多い。
自己回帰モデルは基本的に無限インパルス応答フィルタに一種の変形を加えたものである。
モデルとして定常的であるために、パラメータの値には何らかの制約が必要である。例えば、|φ1| > 1 となる AR(1)モデルは定常的ではない。
例: AR(1)過程
AR(1)過程は次の式で表される。
[math]X_t = c + \varphi X_{t-1}+\varepsilon_t,\,[/math]
ここで、[math]\varepsilon_t[/math] は、[math]\sigma^2[/math]の分散に従うホワイトノイズである（[math]\varphi_1[/math] のような添え字は省いてある）。この過程は [math]|\varphi|\lt 1[/math] であれば、共分散定常性を有する。[math]\varphi=1[/math] であれば、[math]X_t[/math] は単位根を表し、ランダムウォークと見なされ、共分散定常性を有しない。そうでない場合、[math]X_t[/math] の期待値の計算は単純である。ここで共分散定常性を以下のように定式化する。
[math]\mbox{E}(X_t)=\mbox{E}(c)+\varphi\mbox{E}(X_{t-1})+\mbox{E}(\varepsilon_t)\Rightarrow \mu=c+\varphi\mu+0.[/math]
従って、次のようになる。
[math]\mu=\frac{c}{1-\varphi},[/math]
ここで [math]\mu[/math] は平均である。c = 0 なら、平均も 0 になり、分散は次のようになる。
[math]\textrm{var}(X_t)=E(X_t^2)-\mu^2=\frac{\sigma^2}{1-\varphi^2}.[/math]
自己共分散は次の式で表される。
[math]B_n=E(X_{t+n}X_t)-\mu^2=\frac{\sigma^2}{1-\varphi^2}\,\,\varphi^{|n|}.[/math]
この自己共分散関数は減衰時間 [math]\tau=-1/\ln(\varphi)[/math] で減衰する（これを確かめるには、[math]B_n=K\phi^{|n|}[/math] で [math]K[/math] が [math]n[/math] に独立な場合を考えればよい。[math]\phi^{|n|}=e^{|n|\ln\phi}[/math] であり、指数関数的減衰の法則 [math]e^{-n/\tau}[/math] に適合することに注意されたい）。スペクトル密度関数は自己共分散関数の逆フーリエ変換である。離散系では、離散時間逆フーリエ変換が適用される。
[math]\Phi(\omega)= \frac{1}{\sqrt{2\pi}}\,\sum_{n=-\infty}^\infty B_n e^{-i\omega n} =\frac{1}{\sqrt{2\pi}}\,\left(\frac{\sigma^2}{1+\varphi^2-2\varphi\cos(\omega)}\right). [/math]
[math]X_j[/math] が離散的であるため、この式の分母にあるコサインの項が折り返し雑音（エイリアス）を表している。標本化間隔([math]\Delta t=1[/math])が減衰時間([math]\tau[/math])より十分に小さいと仮定すると、[math]B_n[/math] に連続体近似を適用できる。
[math]B(t)\approx \frac{\sigma^2}{1-\varphi^2}\,\,\varphi^{|t|}[/math]
この場合、スペクトル密度はローレンツ分布に従う。
[math]\Phi(\omega)= =\frac{1}{\sqrt{2\pi}}\,\frac{\sigma^2}{1-\varphi^2}\,\frac{\gamma}{\pi(\gamma^2+\omega^2)}[/math]
ここで [math]\gamma=1/\tau[/math] は減衰時間 [math]\tau[/math] に関する角周波数である。
[math]X_t[/math] の別の表現方法として、最初の式で [math]X_{t-1}[/math] を [math]c+\varphi X_{t-2}+\varepsilon_{t-1}[/math] に置き換える方法がある。これを再帰的に N回繰り返すと次の式になる。
[math]X_t=c\sum_{k=0}^{N-1}\varphi^k+\varphi^NX_{\varphi-N}+\sum_{k=0}^{N-1}\varphi^k\varepsilon_{t-k}.[/math]
N が無限大に近づくと、[math]\varphi^N[/math] はゼロに近づき、最終的に次の式が得られる。
[math]X_t=\frac{c}{1-\varphi}+\sum_{k=0}^\infty\varphi^k\varepsilon_{t-k}[/math]
ARパラメータの計算
AR(p)モデルは次の方程式で与えられる。
[math] X_t = \sum_{i=1}^p \varphi_i X_{t-i}+ \varepsilon_t.\,[/math]
これはパラメータ [math]\varphi_i[/math]（i = 1, ..., p）に基づいている。これらパラメータは以下の Yule-Walker方程式で計算できる可能性がある。
[math] \gamma_m = \sum_{k=1}^p \varphi_k \gamma_{m-k} + \sigma_\varepsilon^2\delta_m [/math]
ここで m = 0, ... , p であり、p + 1 個の方程式となる。[math]\gamma_m[/math] は X の自己共分散関数、[math]\sigma_\varepsilon[/math] は入力ノイズ過程の標準偏差、δm はクロネッカーのデルタである。
この式の最後の部分は m = 0 のときだけ 0 でない値となるので、この方程式は一般に m > 0 のときの行列式で表すことで解ける。
[math]\begin{bmatrix} \gamma_1 \\ \gamma_2 \\ \gamma_3 \\ \vdots \\ \end{bmatrix} = \begin{bmatrix} \gamma_0 & \gamma_{-1} & \gamma_{-2} & \dots \\ \gamma_1 & \gamma_0 & \gamma_{-1} & \dots \\ \gamma_2 & \gamma_{1} & \gamma_{0} & \dots \\ \dots & \dots & \dots & \dots \\ \end{bmatrix} \begin{bmatrix} \varphi_{1} \\ \varphi_{2} \\ \varphi_{3} \\ \vdots \\ \end{bmatrix} [/math]
これにより [math]\varphi[/math] が全て求められる。また、m = 0 のときは次のようになる。
[math] \gamma_0 = \sum_{k=1}^p \varphi_k \gamma_{-k} + \sigma_\varepsilon^2 [/math]
これにより [math]\sigma_\varepsilon^2[/math] が求められる。
導出
AR過程を定義する方程式は次の通りである。
[math] X_t = \sum_{i=1}^p \varphi_i\,X_{t-i}+ \varepsilon_t.\,[/math]
両辺に Xt-m をかけて、期待値を求めるとしたとき、次のようになる。
[math]E[X_t X_{t-m}] = E\left[\sum_{i=1}^p \varphi_i\,X_{t-i} X_{t-m}\right]+ E[\varepsilon_t X_{t-m}].[/math]
自己共分散関数の定義から、[math]E[X_t X_{t-m}] =\gamma_m[/math] である。ノイズ関数の値は互いに独立であり、ゼロより大きい m について Xt − m は εt に独立である。m ≠ 0 の場合、[math]E[\varepsilon_t X_{t-m}] = 0[/math] となる。m = 0 の場合、次のようになる。
[math]E[\varepsilon_t X_{t}] = E\left[\varepsilon_t (\sum_{i=1}^p \varphi_i\,X_{t-i}+ \varepsilon_t)\right] = \sum_{i=1}^p \varphi_i\, E[\varepsilon_t\,X_{t-i}] + E[\varepsilon_t^2] = 0 + \sigma_\varepsilon^2, [/math]
従って、次が得られる。
[math]\gamma_m = E\left[\sum_{i=1}^p \varphi_i\,X_{t-i} X_{t-m}\right] + \sigma_\varepsilon^2 \delta_m.[/math]
さらに
[math]E\left[\sum_{i=1}^p \varphi_i\,X_{t-i} X_{t-m}\right] = \sum_{i=1}^p \varphi_i\,E[X_{t} X_{t-m+i}] = \sum_{i=1}^p \varphi_i\,\gamma_{m-i}, [/math]
これにより次の Yule-Walker方程式が導かれる。
[math]\gamma_m = \sum_{i=1}^p \varphi_i \gamma_{m-i} + \sigma_\varepsilon^2 \delta_m.[/math]
移動平均モデル
MA(q)という表記は、次数 q の移動平均モデルを表す。以下の数式で表される。
[math] X_t = \varepsilon_t + \sum_{i=1}^q \theta_i \varepsilon_{t-i}\,[/math]
ここで θ1, ..., θq はモデルのパラメータ、εt, εt-1,... は誤差項である。移動平均モデルも無限インパルス応答フィルタに一種の変形を加えたものである。
自己回帰移動平均モデル
ARMA(p, q)という表記は、p次の自己回帰とq次の移動平均を組合わせたモデルを指す。以下の数式で表される。
[math] X_t = \varepsilon_t + \sum_{i=1}^p \varphi_i X_{t-i} + \sum_{i=1}^q \theta_i \varepsilon_{t-i}.\,[/math]
誤差項
誤差項 εt は一般に「独立かつ同一の分布に従う」(i.i.d.)無作為変数であり、ゼロを平均値とする正規分布に従う。すなわち εt ~ N(0,σ2) で、σ2 は分散である。このような仮定を弱めることもあるが、そうするとモデルとしての性質が変化する。特に、i.i.d. という仮定を変更すると根本的な性質が変化する。
ラグ（遅れ）作用素を使った記法
ARMAモデルをラグ作用素（遅れ作用素） L を使って表す場合もある。この場合、AR(p)モデルは次のように表される。
[math] \varepsilon_t = \left(1 - \sum_{i=1}^p \varphi_i L^i\right) X_t = \varphi X_t\,[/math]
ここで、φ は次の多項式で表される。
[math] \varphi = 1 - \sum_{i=1}^p \varphi_i L^i.\,[/math]
また、MA(q)モデルは次のように表される。
[math] X_t = \left(1 + \sum_{i=1}^q \theta_i L^i\right) \varepsilon_t = \theta \varepsilon_t\,[/math]
ここで θ は次の多項式で表される。
[math] \theta= 1 + \sum_{i=1}^q \theta_i L^i.\,[/math]
以上から、ARMA(p, q)モデルは次のように表される。
[math] \left(1 - \sum_{i=1}^p \varphi_i L^i\right) X_t = \left(1 + \sum_{i=1}^q \theta_i L^i\right) \varepsilon_t\,[/math]
あるいは、もっと簡潔に記せば、次のようになる。
[math] \varphi X_t = \theta \varepsilon_t.\,[/math]
ラグ作用素とは、時系列データのある時点のデータで他の時点のデータを表すように係数化したもの。上記の式はいずれも Xt しか出現しない（他の時点のデータが出てこない）ことに注意されたい。他の時点のデータは全てラグ作用素によって表されている。
実データへの適用
実データに適用する場合、ARMAモデルの p と q を選択後、誤差項を最小化するパラメータを探るため最小二乗法を使うのが普通である。また、実データに適合する最小の p および q を見つけることでよい結果が得られることが知られている。純粋なARモデルでは、これに Yule-Walker 方程式を利用することができる。
一般化
Xt の過去の値や誤差項 εt との依存関係は特に何かが判明するまで線形な関係であると見なすのが普通である。依存関係が非線形である場合、このモデルを「非線形移動平均」(NMA)、「非線形自己回帰」(NAR)、「非線型自己回帰移動平均」(NARMA)モデルと呼ぶ。
自己回帰移動平均モデルは他の方法でも一般化される。例えば、ARCH（自己回帰条件付き分散変動）モデルやARIMA（自己回帰和分移動平均）モデルがある。複数の時系列を扱うベクトルARIMAモデルもある。データに季節変動効果がある場合、SARIMA（季節ARIMA）モデルでモデル化される。
その他の一般化として、「多変量自己回帰」(MAR)モデルがある。
関連項目
参考文献
- George Box and Gwilym M. Jenkins. Time Series Analysis: Forecasting and Control, second edition. Oakland, CA: Holden-Day, 1976.
- Mills, Terence C. Time Series Techniques for Economists. Cambridge University Press, 1990.
- Percival, Donald B. and Andrew T. Walden. Spectral Analysis for Physical Applications. Cambridge University Press, 1993.
- Yoshitsugu Hayashi,Hiroshi Ohkama,Yoshitaka Fujiwara. An Estimation Method of Auto-Regressive Parameters with Time-varying Cost. Faculty of Enginnering, Kitami Institute of Technology, 1997.