ARモデルの誤差分析について

その他

2019.3.7

Topics

こんにちは。データサイエンスチームもtmtkです。
この記事では、ARモデルの誤差分析について説明します。

はじめに

前の3回の記事にわたって、時系列解析の入門としてARモデルについて紹介し、AICによるARモデルのモデル選択について解説し、ARモデルの表現能力について議論しました。
ARモデルをつかうと過去の時系列データから未来のデータの値を予測することができます。実際に予測を行う場合には、その予測の信頼区間も知りたいことが多いと思います。つまり、予測がどれくらいずれうるのかまでを知りたいということです。
この記事では、ARモデルで未来のデータを予測するとき、予測の誤差について考察します。

[latex]y_{N+1}[/latex]の誤差分析

ARモデルは、式

[latex] y_t = \sum_{i = 1}^m a_iy_{t – i} + v_t[/latex]

で表され、誤差項[latex]v_t \sim N(0, \sigma^2)[/latex]が誤差をつかさどっているのでした。
そのため、いま、過去のデータ[latex]y_1, \ldots, y_N[/latex]とARモデル[latex]y_t = \sum_{i = 1}^m a_iy_{t – i} + v_t[/latex]が与えられているとき、[latex]y_{N+1}[/latex]とそのARモデルからの推定値[latex]\hat y_{N+1} = \sum_{i = 1}^m a_iy_{N + 1 – i}[/latex]の差は、[latex]v_{N+1}\sim N(0, \sigma^2)[/latex]になります。よって、[latex]y_{N+1}[/latex]に関する信頼係数95%の信頼区間は、

[latex]\hat y_{N+1} – 1.96\sigma \leq y_{N+1} \leq \hat y_{N+1} + 1.96\sigma [/latex]
と計算できます。
データ分析を実際にする状況では、誤差項[latex]v_t[/latex]のしたがう確率分布の母分散[latex]\sigma^2[/latex]は知ることができません。そのため、最尤推定による推定値[latex]\hat\sigma ^2 =\frac{\sum_{t=m+1}^N\left(y_t – \sum_{i=1}^m a_i y_{t – i}\right)^2}{N – m}[/latex]を母分散[latex]\sigma^2[/latex]の近似値として使わざるをえません。この場合、[latex]\sigma^2[/latex]と[latex]\hat\sigma^2[/latex]の間にも誤差があるため、正規分布ではなくt分布を持ち出して議論する必要があります。
ただし、この記事では、今後の計算を簡単にするため、以下では[latex]\sigma^2 = \hat \sigma^2[/latex]を仮定します。つまり、最尤推定による母分散の推定値が真の母数に等しいとみなして計算を進めることにします。あるいは、母分散[latex]\sigma^2[/latex]が既知であると仮定していると考えてもいいです。

一般の場合([latex]y_t \quad (t > N)[/latex])の誤差分析

時刻[latex]t[/latex]での[latex]\{y_t\}_t[/latex]の予測値を表す数列を、[latex]\{\hat y_t\}_t[/latex]とおきます。いま、過去のデータ[latex]\{y_t\}_{t=1}^N[/latex]が与えられているとします。すると、[latex]1 \leq t \leq N[/latex]の範囲では観測値[latex]y_t[/latex]が与えられているため、予測値[latex]\hat y_t[/latex]は[latex]\hat y_t = y_t[/latex]とおくことができます。
未来の時刻[latex] t > N [/latex]に対しては、ARモデルとそれまでの推定値[latex]y_{t’} \, (t’ < t)[/latex]から[latex]\hat y_t = \sum_{i = 1}^m a_i \hat y_{t – i}[/latex]と計算することができます。
未来の時刻[latex] t > N[/latex]に対して、予測誤差[latex]y_t – \hat y_t[/latex]のしたがう分布を考えてみましょう。前の議論から、[latex] t = N + 1[/latex]のときは[latex] y_{N+1} – \hat y_{N+1} = v_{N+1} \sim N(0, \sigma^2)[/latex]に従うのでした。次に[latex]t = N + 2[/latex]の場合を考えます。すると、

[latex]\begin{aligned}y_{N+2} – \hat y_{N+2} &= (\sum_{i=1}^m a_iy_{t-i} + v_{N+2}) – (\sum_{i=1}^ma_i\hat y_{t-i}) \\ &= (\sum_{i=1}^m a_i y_{t – i} + v_{N+2}) – (\sum_{i=1}^m a_i y_{t-i} – a_1 v_{N+1}) \\ &= v_{N+2} + a_1 v_{N+1}\end{aligned}[/latex]
と計算されます。したがって、

[latex]y_{N+2} – \hat y_{N+2} = v_{N+2} + a_1 v_{N+1} \sim N(0, \sigma^2 + a_1^2\sigma^2)[/latex]
となります。
これまで時刻[latex]t = N + 1, N + 2[/latex]の場合について考察しました。時刻[latex]t[/latex]が一般の場合を考えましょう。誤差の確率変数の列[latex]\{\varepsilon_t\}_t[/latex]を[latex] \varepsilon_t = y_t – \hat y_{t}\, (t=1, 2, \ldots)[/latex]と定義すれば、

[latex]\begin{aligned} \varepsilon_t &= y_t – \hat y_t \\ &= \left(\sum_{i=1}^m a_i y_{t-i} + v_t\right) – \left(\sum_{i=1}^m a_i \hat y_{t-i}\right) \\ &= \left(\sum_{i=1}^m a_i y_{t-i} + v_t\right) – \left(\sum_{i=1}^m a_i y_{t-i} – \sum_{i=1}^m a_i \varepsilon_{t-i} \right) \\ &= v_t + \sum_{i=1}^m a_i \varepsilon_{t – i}\end{aligned}[/latex]
となります。定義より[latex]\varepsilon_t = 0 \, (t \leq N)[/latex]であることに注意すると、[latex]\varepsilon_t \,(t > N)[/latex]は平均[latex]0[/latex]の正規分布に従い、その分散[latex]\sigma^2_t[/latex]は

[latex]\sigma^2_t = \left\{\begin{array}{ll}\sigma^2 + \sum_{i=1}^m a_i^2 \sigma^2_{t-i} & (t = N + 1, N + 2, \ldots)\\ 0 & (t = 1, \ldots, N)\end{array}\right.[/latex]
で計算されることがわかります。
以上より、[latex]y_t[/latex]に関する95%信頼区間は、上で定義した[latex]\sigma_t[/latex]を使って、

[latex]\hat y_t – 1.96 \sigma_t \leq y_t \leq \hat y_t + 1.96\sigma_t[/latex]
と表すことができます。


(真の値、予測値、95%信頼区間をプロットしたもの)

まとめ

この記事では、ARモデルの予測の信頼区間を計算しました。ただし、計算にあたって簡単のため誤差項[latex]v_t[/latex]の母分散[latex]\sigma^2[/latex]が既知であることを仮定しました。

参考文献

  • 北川源四郎『時系列解析入門』

テックブログ新着情報のほか、AWSやGoogle Cloudに関するお役立ち情報を配信中!

tmtk

データ分析と機械学習とソフトウェア開発をしています。 アルゴリズムとデータ構造が好きです。

Recommends

こちらもおすすめ

Special Topics

注目記事はこちら