隼時系列本：ノート3 - クッキーの日記

以下の本を読みます。

時系列分析と状態空間モデルの基礎: RとStanで学ぶ理論と実装
馬場真哉

プレアデス出版 2018-02-14
売り上げランキング : 7742

Amazonで詳しく見る by G-Tools

※ 以下、キャラクターが会話します。原作とは関係ありません。上の本からそれる話も多いです。誤りがあればご指摘ください。

前回：ノート2 ／次回：まだ

f:id:cookie-box:20180305231302p:plain:w60

第2部の4章までで色々なモデルを扱って、第2部の5章からは、じゃあどのモデルを選べばいいかって話か…。

f:id:cookie-box:20180305232608p:plain:w60

ここでモデルという言葉は ARIMA(2,0,0) と ARIMA(1,1,0) は異なるモデルということですね。なので、モデル選択というのは次数の決定のことであると書いてありますね。

あれ、AR モデルか MA モデルか ARIMA モデルか SARIMA モデルか、っていう選択を最初にするのかと思ったんだけど、その選択はないの？

まず、ARIMA モデルは AR モデルと MA モデルを包含しますよ。SARIMA モデルはさらに ARIMA モデルを包含しますが、周期を与えてつかうモデルなので、対象データに周期がある場合には SARIMA をつかってみればいいんじゃないですか。

あ、そうか。じゃあやっぱり最初に取りかかるのは次数の選択なのか。確かにARIMAモデルをつかうとしても次数をどうすればいいのか全然わからないよな…これはどういう手順で決めるのか気になるぞ…えっ、「手当たり次第」？ショックなんだけど…。

ショックを受けることないでしょう。僕たちは計算機が発達した時代に生まれたんです。考えなくても解決できる手順を中途半端に考える必要はありません。下手の考え休むに似たりです。しかし、計算機にまかせるにせよ、モデルのよさの指標については把握しておく必要があるでしょう。

う、うん…まずAICっていう指標が紹介されてるよ…それでその前に、尤度？尤度って何？この61ページのコイン投げの 2/9 っていう尤度は大きいの？小さいの？

コインの表が出る確率が 1/3 ではなく 1/2 だったらその尤度はどうなります？

それは、1/2×1/2＝1/4 になる。さっきの 2/9 と比べると 1/4 の方が大きい？

なので、あるコインについて「2回投げたら1回目が表で2回目が裏が出た」ということがわかっているとき、このコインの表が出る確率は 1/3 よりも 1/2 の方が尤もらしいということです。もっと一般に、コインの表が出る確率を $p$ とすると、「2回投げたら1回目が表で2回目が裏が出た」という結果に対する $p$ の尤度は $f(p)＝p(1-p)$ ですから、これを最大にするのは $p＝1/2$ です。なので結局、2回投げた結果しかわからない範囲では、このコインの表が出る確率として最も尤もらしいのは 1/2 です。

…2回投げて1回目が表で2回目が裏だったんだから、コインの表が出る確率が 1/2 らしいって当たり前なんじゃ？

まあ平均値については、「標本平均の期待値」が「真の平均値」に一致しますからね。でも平均値ではない分布のパラメータについては決して当たり前では…この話は脇道にそれるので置いておきましょう。…ところでハヤト、このコイン投げの結果から考えられる、もっと尤もらしいコインはあるでしょうか。

もっと尤もらしい？「2回投げたら1回目が表で2回目が裏が出た」という結果から考えられる？さっき尤度 $f(p)$ を最大にするのは $p=1/2$ だってジュン自分で言ってたじゃん。だから「表が出る確率が 1/2 のコイン」より尤もらしいコインはないだろ？

$f(p)$ はコインの表が出る確率が1回目も2回目も同じ場合の尤度ですので。2回目に投げるときは1回目に投げるときとはコインの性質が変化しているかもしれません。コインを特徴付けるパラメータを $p_1$ と $p_2$ の2個に増やして、1回目に投げるときは確率1で表が出て、2回目に投げるときは確率1で裏が出るコインを考えるとどうでしょう。こうすれば尤度は 1×1＝1 です。テキストの61ページと同様の例ですが。

えーそれはズルイだろ。

ええ、ずるいんです。ずるいというよりは、このモデリングは3回目に投げた結果を予測するのに全く有用ではないでしょう。なので、なるべくパラメータを増やさずに尤度を最大化したいんです。なるべくパラメータを増やさずに尤度を大きくできたかどうかの指標が62ページの AIC です。

62ページ…これがAICの式か…最大化対数尤度？

61ページに説明がありましたが、あるモデル、例えば ARMA(1,1) でも何でもいいです、ARMA(1,1) だと $c, \, \phi_1, \, \theta_1$ の3つのパラメータが推定対象ですね。対象時系列データの下での尤度が最大になるようにこれらのパラメータを最適化した上で、その尤度に対数をとったものが最大化対数尤度です。最後に対数をとる必要もありませんので、対数尤度が最大になるように3つのパラメータを最適化したときに対数尤度、と言った方がいいでしょうか。次数の異なるモデル間でよさを比較するのだから、パラメータセットはベストにしておかなければならないでしょう。

…時系列データに対するARIMAモデルの尤度ってどう計算するの？時系列データって、時点ごとにデータたくさんあるよな？

さっきのコイン投げの例でも、あるパラメータの下での2回のコイン投げの結果の同時確率がそのパラメータの尤度だったでしょう。2回だろうと100回だろうと一緒ですよ。与えられた100日分の気温に対するあるパラメータセットの尤度は $p($ 1日目～100日目の気温 $|$ パラメータセット $)$ です。もっとも、ARIMAモデルの尤度関数は…インターネットを探せばありますよ。例えばこのリンクなどがそうですね。

http://econ.nsysu.edu.tw/ezfiles/124/1124/img/Chapter17_MaximumLikelihoodEstimation.pdf

…まあいいや。その最大化対数尤度と、パラメータ数で AIC が決まるのか。尤度は大きいほどよくて、パラメータ数は小さいほどよくて、AICは小さいほどよい指標なんだな。…ジュンさっき、なるべくパラメータを増やさずに尤度を大きくしたいって言ってたよな。(2-41) 式をみるとさ、パラメータを1つ増やすんだったら、対数尤度は1よりも大きく増えてくれなきゃ困るってことだよな。なんでパラメータ1個と対数尤度1単位が釣り合うんだろう。あとなんで第1項も第2項も2がかかってるの？

以下のリンクに AIC の導出がありますよ。もっとも、AIC は原論文も探せばありますが。自由度 k のカイ2乗分布の期待値が k なので AIC は (2-41) 式のようになっているんですね。全体に2がかかっているのは、こうすると第1項が漸近的にカイ2乗分布にしたがうからですね。

http://www.is.titech.ac.jp/~shimo/pub/koza2000/ohp-20000703b.pdf

ふーん…まあなんやかんやで計算機が ARMA の次数を決めてくれるんだな。でも階差についてはそういうわけにいかないのか…それが63ページからの検定…帰無仮説って何？

仮説検定ですか。例えば有意水準5%で仮説検定するのだったら、「帰無仮説が正しいとしたら起きる確率が5%以下になるくらい珍しいこと」を「偶然それが起きたと考えるには低い」とみなすんです。与えられたデータが、「帰無仮説が正しいとしたらこんなデータが得られる確率は5%以下だ」というデータだったら、帰無仮説を棄却します。つまり、帰無仮説は正しくないと考えるんです。もし帰無仮説が正しいときに与えられたデータが得られる確率が5%よりも大きいのだったら、帰無仮説を受容します。これは、帰無仮説は正しくないとまではいえない、という状態ですね。

うーん…なんか具体例ない？

今調べたんですけど、高校2年生男子の平均体重は60.4kg、標準偏差は9.93kgらしいです。この平均体重が未知だったとして、僕とハヤトが高校2年生男子の無作為標本だったとき、帰無仮説「高校2年生男子の平均体重は60kg以上」は棄却されるでしょうか。有意水準5%とします。…分散既知の場合の平均の検定では、標本平均が正規分布にしたがうことを利用します。今回は標本平均が下側5%点より小さい領域が棄却域です。正規分布の下側5%点は -1.64 で、サンプル数2なので、棄却域は $60 - 1.64 \times \sqrt{9.93^2 / 2} = 48.48$ 未満ですね。僕のハヤトの体重の平均は 48.5kg なので棄却域になく、帰無仮説は棄却されません。帰無仮説が正しい下で5%も起こりえない珍しいデータではないということですね。

でもすごいギリギリ…。

（ノート4があれば）つづく