MATLABで『時系列解析入門』方針・1章
最終更新:2021/07/14
非専門家が書いています. 十分に批判的に読んで頂くようお願いいたします. 間違い・疑問点などあれば, ぜひコンタクトフォームへ連絡いただけると幸いです.
方針
僕は音波・気候変動データ・動画などと様々な時系列データを扱っているが, Rなど時系列解析には便利なツールが多々あるので自分でアルゴリズムを実装する機会は少なかった. しかし, 将来を見据えると, 一度しっかり自分の手でも実装し理解を深める必要があると感じた. そこで, かつて読んだことがあり, お気に入りである北川源四郎先生の『時系列解析入門』をMATLABで実装してみることにした.
この本のいいところはたくさんあるが, 一つ例が多い点がある. なので, 実装してそれがしっかり動いているのか確認しやすい. また, この本にのっとった講義動画が以下のサイトにあり, 非常にわかりやすい ( 東京大学OCWx数理手法Ⅶ ) . 日本 ( いや世界? ) を代表する時系列解析の大家の講義動画をこのように見ることができるのは本当にすごいことだと思う.
とはいえ, すべて ( 例えば分散の計算など ) を一から書くのもやりすぎなきがするし, ここで作ったコードは実際の解析に使いたいと思っているのであまりに複雑にもしたくない. そこで, ”基本的な方針として”, MATLABの初等的な関数, 具体的にはどのTOOLBOXも必要としない関数, は積極的に使うこととした.
利用データ
データは, 北川先生のウェブサイトや講義サイトで公開されている. これらをかき集めてきた.
例:教科書p.2にある太陽黒点数 ( 実はこのデータはMATLABにも入っており, load sunspot.mat
で得られる ) .
あと, せっかくなので新しいデータも手持ちから一つとってきた ( 一部拡大 ) .
1章の内容
簡単な例で練習.
p.9 対数差分による前処理
d = nikkei; subplot(211) plot((nikkei)) title('日経225平均株価') subplot(212) plot(diff(log(nikkei))) title('対数差分')
図1.2 ( p.9 ) が再現されている.
p.10 前期比・前年同期比
d = whard; subplot(311); plot(whard) title('WHARD original') xticks(0:12:length(d)); subplot(312); dd = d(2:end)./d(1:end-1); % main Eq.(1.6) plot(dd) title('前期比') xticks(0:12:length(d)); subplot(313); dd = d(13:end)./d(1:end-12); % main Eq.(1.6) plot([ones(1, 12) dd]) % 頭に1をパディング title('前年同期比') xticks(0:12:length(d)); yline(1.0) ylim([0.7 1.3])
図1.3 ( p.11 ) が再現されている.
p.11 移動平均
移動平均には, movmean
という関数がある.
d = tokyo; subplot(221) plot(d); ylim([0 40]) legend('original', 'Box', 'off') subplot(222) plot(movmean(d, 5)); ylim([0 40]) legend('k = 5', 'Box', 'off') subplot(223) plot(movmean(d, 17)); ylim([0 40]) legend('k = 17', 'Box', 'off') subplot(224) plot(movmean(d, 29)); ylim([0 40]) legend('k = 29', 'Box', 'off')
図1.4 ( p.12 ) が再現されている.
movmean
は, はしっこの処理がどうなっているかしっかり認識して使うこと.
例:
>> movmean([1 2 3 4 5 6], 3) ans = 1.5000 2.0000 3.0000 4.0000 5.0000 5.5000