時系列解析と疫学

Transcript

1. NIID Winter seminer 2023 時系列解析と疫学 国⽴感染症研究所 第12室⻑ ⽶岡 ⼤輔
2. ⽬次 • 時系列の回帰分析 ADF検定、DW検定、GLS、潜在成⻑曲線モデル，ITS • 状態空間モデル ローカル線形トレンドモデル等，隠れマルコフモデル，変化点モデル • 時系列データ同⼠の関係 ラグ相互相関、DTW，コヒーレンス， ウェーブレットクロススペクト ル解析，時系列クラスタリング，グレンジャー因果性検定，直交化イ ンパルス応答 • 多変量時系列データの要約 時系列因⼦分析，関数主成分分析 2/60 かるく概観するだけならとても良い本です！ 本書の内容をカバーします。
3. 時系列分析の⽬的 【記述】 時系列変化を図⽰し、特徴を表現する 【説明】 時系列変化をモデル化し特徴を明らかにする 【予測】 現在までのデータから将来の値を推計する 【制御】 操作可能な変数を変化させ、アウトカムを 望ましい状態にする 3/48 疫学ではこの辺がメインか？ ⾦融や地震などは明確に 「予測」に重点 = black-box 的予測でも当てれればOK (deepでポン)
4. 最も基礎的な指標: ⾃⼰相関関数 (ACF)とコレログラム Auto-correlation function (ACF): ⾃⼰相関係数 時系列の元データとラグk だけずらしたデータの相関 ⾒えること： 隣り合った時点がどれくらい類似している? 周期性がある？その周期性は？ 4/48 条件：yt が定常過程(任意の時点tで確率分布が同じ) 標本分散 <latexit sha1_base64="/micsJPtjP4hMk9QkmEHHKknM2Y=">AAACqHicfVDLbhMxFHWGVwmvFJZsLCKktNBopuK1oFKhCLFBBKlpijLpyOM4E2v8GNl3ECPLP8afsGMLX4EnDYi2iCtZOj7nXl+fk1eCW4jjb53o0uUrV69tXO/euHnr9p3e5t0jq2tD2Zhqoc1xTiwTXLExcBDsuDKMyFywSV4etPrkMzOWa3UITcVmkhSKLzglEKisd/jq4G1W4j2czheGUJd4N1A75Va6JOBSywtJ/MmuT20tMwd7iT9xQfaDJoOdNCfGNX4rXBw8Kv0fIuv142G8KnwRJGvQR+saZZudl+lc01oyBVQQa6dJXMHMEQOcCua7aW1ZRWhJCjYNUBHJ7Myt7Hv8MDBzvNAmHAV4xf494Yi0tpF56JQElva81pL/0qY1LF7MHFdVDUzR00WLWmDQuM0Sz7lhFEQTAKGGh79iuiQhRgiJn9nSvg1aCxusvGHBomHvA/WhYoaANtsuJaaQXPlguUgft+h/jeTL78aAut2Qd3I+3YvgaHeYPBs+/fikv/96nfwGuo8eoAFK0HO0j96hERojir6i7+gH+hltR6NoEn06bY0665l76ExF+S+TLdSQ</latexit> ACFk = 1 (n k)ˆ2 n k X t=1 (yt ¯ y)(yt+k ¯ y) 標本平均 <latexit sha1_base64="nwku9awieEJziCGY3eUmvNHxGCk=">AAACenicfVDbahsxEJW3t9S9Oc1jX0RMStI6Zrekl4cEQtNAX0pSqJOAd1lm5bEtopUWaTZ0WfwH/Zq+Jj/Sf+lDtY4DTVI6QnA4c0ajc7JCSUdh+KsV3Ll77/6DpYftR4+fPH3WWX5+5ExpBQ6EUcaeZOBQSY0DkqTwpLAIeabwODvda/rHZ2idNPobVQUmOUy0HEsB5Km083J/WKWU8B0e9njsz545a5ger9KaNqezZEesTzfSTjfsh/Pit0G0AF22qMN0ubUdj4woc9QkFDg3jMKCkhosSaFw1o5LhwWIU5jg0EMNObqknhua8TXPjPjYWH818Tn790QNuXNVnnllDjR1N3sN+a/esKTxh6SWuigJtbhcNC4VJ8ObdPhIWhSkKg9AWOn/ysUULAjyGV7b0rxNxijnrXxCb9HiF08dFGiBjH1Vx2AnudQzb3kS9xr0PyF8vxJ61G77vKOb6d4GR2/60bv+269b3d2Pi+SX2Au2ytZZxN6zXfaZHbIBE+wH+8nO2UXrd7AabASvL6VBazGzwq5VsPUHLDDAyA==</latexit> E[yt] = 0, Cov[yt, yt h] = c(h) 重要な仮定 弱定常性 (Weak stationarity) or ⼆次定常性 期待値と分散共分散が時刻t によらず同⼀の関数に従う なんか周期性が⾒え ますねぇ ↓ 周期性を考慮したモ デリング なんか⻑期的な⾃⼰ 相関が⾒えますねぇ ↓ 階差{yt -yt-1 }をとって みよう
5. 定常性、⼤事 (詳しいことは後ほど) とってもざっくり⾔うと「薄⽬で⾒ると平べったい」 定常じゃない = ⾮定常じゃないとどうなる？ • ⾒せかけの回帰が発⽣しうる：単位根過程であった場合、回帰分析にお いて（本当は違うのに）有意であるかのような結果が出てしまう。 • 回帰による誤差がトレンドをもつ可能性があり、その場合、残差は定常 という(線形)回帰の仮定を満たすことができず、結果⾃体の信頼性が⼤ きく揺らぐことになる。 5/48 ⾮定常性の種類 【トレンド⾮定常性】 平均が時間とともに変化する⾮定常性(i.e., なんかトレンドがある)。 → 差分を取る 【季節性⾮定常性】 ⼀定の周期を持つ季節変動に存在する⾮定常性(i.e., なんか季節性 がある) → 季節差分を取る 【分散⾮定常性】分散が時間とともに変化する⾮定常性(i.e., なんか段々広がっていく) 。 →対数変換などを⽤いる
6. 時系列の分解 時系列データyt を「トレンドmt 、季節成分st 、誤差et 」に分解 STL: Seasonal and Trend decomposition using LOESS が有名 局所多項式回帰（LOESS）を使⽤してデータのトレンドと季節成分を抽出 それ以外の成分として誤差を推定 成分抽出の精度は、LOESS内の近傍サイズや 多項式次数などのパラメータに依存します。 6/48 <latexit sha1_base64="AHIHSXErya3gicrdXcx17vY9MTs=">AAACaXicfVBdS1tBEN3c2mqjrUlLQezLYiiUVsK9otaHFkR98KXUglEhCZe5m0lc3I/L7txiuOTH9LX9Rf4G/4R7Ywp+FAd2OJw5szNzslxJT3F8VYuezT1/Mb/wsr649Or1cqP55sTbwgnsCKusO8vAo5IGOyRJ4VnuEHSm8DS72K/qp7/QeWnNMY1z7GsYGTmUAihQaePdOCX+jeuQP3M/zZhS2mjF7Xga/DFIZqDFZnGUNmtfewMrCo2GhALvu0mcU78ER1IonNR7hcccxAWMsBugAY2+X073n/APgRnwoXXhGeJT9m5HCdr7sc6CUgOd+4e1ivxfrVvQcKdfSpMXhEbcDhoWipPllRl8IB0KUuMAQDgZduXiHBwICpbdm1L9TdYqH045wHCiw++B+pGjA7LuU9kDN9LSTMLJo956hZ4SwuU/YUD1evA7eejuY3Cy0U6221s/N1u7ezPnF9h7tsY+soR9YbvskB2xDhOsZL/ZH/a3dh01o5Vo9VYa1WY9b9m9iFo3Zgq6lw==</latexit> yt = mt + st + et Rだと簡単、stl関数に⼊れるだけ
7. ARモデル ⾃⼰回帰モデル (Autoregression model: AR model) 気持ち：今の値は1期前の値に依存しているよ! (今⽇が⾼けりゃ明⽇も⾼い) 重要な条件: |φ|<1 (定常過程になる) (ちなみにφ=1でランダムウォーク(単位根仮定)) 7/48 <latexit sha1_base64="A9owpiGwcWp/Z26eLujI1d+PUN0=">AAACmHicfVBdb9MwFHUzPkb56rY3eLGokAaUKplg2wNIFUwCHoCB6DapKdGNe5tas+PIvpmIov4hfg2v49fgdEWwDXEtS0fnnuvrc9JCSUdheNoKVq5cvXZ99Ub75q3bd+521tYPnCmtwKEwytijFBwqmeOQJCk8KiyCThUepsevm/7hCVonTf6FqgLHGrJcTqUA8lTS2asS4i95XMwkr5KankbzJ/EJWCycVF5APR435w/FYyc1/7AZ9jzINHzdepR0umE/XBS/DKIl6LJl7SdrrRfxxIhSY05CgXOjKCxoXIMlKRTO23HpsABxDBmOPMxBoxvXC7tz/tAzEz411t+c+IL9e6IG7VylU6/UQDN3sdeQ/+qNSprujmuZFyVhLs4WTUvFyfAmOz6RFgWpygMQVvq/cjEDC4J8wue2NG+TMcp5K3voLVp876mPBVogYx/XMdhMy3zuLWdxr0H/E8K330KP2m2fd3Qx3cvgYKsfbfeff3rWHbxaJr/K7rMHbJNFbIcN2Fu2z4ZMsO/sBztlP4N7wSB4E7w7kwat5cwGO1fB5188jMyW</latexit> yt = yt 1 + "t, " ⇠ N(0, 2) AR(1) <latexit sha1_base64="59i+NiJaXShHFAEs5Dx6Yh6GS9k=">AAACwXicfVFLbxMxEPZueZTwSuHIxSJCKjREuxGvA0gV5cAFKBJpK2XDyutMNlb9kj1bsVrlV3Lip3DDmw2CtoixLH3zzWeP53NhpfCYJD+ieOvK1WvXt2/0bt66feduf+fekTeV4zDhRhp3UjAPUmiYoEAJJ9YBU4WE4+L0oK0fn4HzwugvWFuYKVZqsRCcYaDyvqpzpG9oZpciT2mdN/g0Xe2t03GXjkM6N+jpXqeyHW0DfcYcWC9kuAiHNGvXH4pmXij6cTcZBlAq9nX8uJf3B8koWQe9DNINGJBNHOY70evQnFcKNHLJvJ+micVZwxwKLmHVyyoPlvFTVsI0QM0U+Fmz9mVFHwVmThfGha2Rrtm/TzRMeV+rIigVw6W/WGvJf9WmFS5ezRqhbYWgeddoUUmKhrYm07lwwFHWATDuRHgr5UvmGMfwFee6tHejMdKHUd5BGNHBh0B9suAYGvekyZgrldCrMHKZDVv0PyH79lsYUK/1O73o7mVwNB6lL0bPPz8b7L/dOL9NHpCHZJek5CXZJ+/JIZkQTr6Tn1EcbcUHsYht7DppHG3O3CfnIm5+AYaE2gE=</latexit> yt = 1yt 1 + 2yt 2 + · · · + pyt p + "t, " ⇠ N(0, 2) AR(p) pは⾃分で与 えてね
8. MAモデルとARMAモデル 移動平均モデル (Moving-average model: MA model) 気持ち：今の値は1期前のノイズに依存しているよ! 昨⽇、偶然インスタでバズって売上急上昇（= 急に来た誤差、ノイズ）は数⽇間(q⽇間)影響があるよ ⾃⼰回帰移動平均モデル (ARMA model) ARとMAを合体させたら強い 8/48 MA(q) qは⾃分で 与えてね <latexit sha1_base64="jb9r55NNH4JeKhnps+7AhdAghX8=">AAAC5HicfZFNb9QwEIad8FXC17Yc4WCxKiq0rJIVBQ4gVcCBC1Aktq20WSLHO8latePUnlSsovwDDkiIK/+LP8IZZ7srdbeIkSy9ep+xxzOTllJYDMPfnn/p8pWr19auBzdu3rp9p7O+cWB1ZTgMuJbaHKXMghQFDFCghKPSAFOphMP0+E3LD0/BWKGLzzgtYaRYXohMcIbOSjrfpwnSVzRW1TaNT5mB0grpAG7HOAFkSXTerfFJ1CxIf4X0HRlrtAt+ssJPGhrHwVIRGluh6IetcMeJXLEv/UdJpxv2wlnQiyKaiy6Zx36y7r10ZXmloEAumbXDKCxxVDODgktogriyUDJ+zHIYOlkwBXZUz0bX0E3njGmmjTsF0pl7/kbNlLVTlbpMxXBiV1lr/osNK8xejGpRlBVCwc8KZZWkqGm7BzoWBjjKqROMG+H+SvmEGcbRbWupSvs2ai2ta+UtuBYNvHfWxxIMQ20e1zEzuRJF41rO451W/S+RfV0kOhUEbt7R6nQvioN+L3rW2/30tLv3ej75NXKPPCBbJCLPyR55R/bJgHDyx7vvbXoP/cz/5v/wf56l+t78zl2yFP6vvz906qM=</latexit> yt = µ + "t + ✓1"t 1 + ✓2"t 2 + · · · + ✓q"t q "t ⇠ N(0, 2) <latexit sha1_base64="+LcBmQG6r3ifQP6i1fdr1U70Jwo=">AAADFnicfVFNj9MwEHXC1xK+unDkYlEhFbpUSQS7HEBaAQcuLItEt5WaEjmum1prJ157siKK8j/4NRyQEFeu/BucppW2XcRIlt689zzjGSdKcAO+/8dxr1y9dv3Gzk3v1u07d+91du+fmLzQlA1pLnI9TohhgmdsCBwEGyvNiEwEGyWnbxt9dM604Xn2GUrFppKkGZ9zSsBScedHGQN+jSO14HGAy7iCZ0HdX6Zhm4Y2neVgcL91qZZWdUOcE82U4cKWgn4ECwYkDi6yq3qtEm4p69Jr/WxLP6txFHkbTXBkuMRHPX/PglSSL+GTuNP1B/4y8GUQrEAXreI43nVe2ba0kCwDKogxk8BXMK2IBk4Fq72oMEwRekpSNrEwI5KZabVcdo0fW2aG57m2JwO8ZC/eqIg0ppSJdUoCC7OtNeS/tEkB85fTimeqAJbRttG8EBhy3PwcnnHNKIjSAkI1t2/FdEE0oWD/d6NLUxvyXBg7yjtmR9Tsg6U+KqYJ5PppFRGdSp7VduQ02mvQ/4zk69pokefZfQfb270MTsJBsD948el59/DNavM76CF6hHooQAfoEL1Hx2iIqNNzjpyRM3a/ud/dn+6v1uo6qzsP0Ea4v/8CvMb8ng==</latexit> yt = 1yt 1 + 2yt 2 + · · · + pyt p + "t + ✓1"t 1 + ✓2"t 2 + · · · + ✓q"t q "t ⇠ N(0, 2) AR(p) MA(q) ARMA(p,q) 常に定常性 を満たす
9. ARIMAモデル ⾃⼰回帰和分移動平均モデル (Auto-regressive Integrated Moving average model: ARIMA model) AR、MA、ARMAは定常なモデルには有効だが、⾮定常なモデルには無⼒ 階差を取ると定常(後述のd 次和分過程 I (d )) になることが多々ある 階差をとったあとのデータに対してARMA(p,q)を考える 9/48 <latexit sha1_base64="wYd8u2KGTQobcdUofwknQlQKz5Q=">AAADGnicfVFNj9MwEHXC1xK+unDkYlEhFtqtkoivA0gr4MAFWBDdXakpkeO6qbV27LUnK6oo/4RfwwEJceXCv8FpU2nbRYxk6c17zzOecaYFtxCGfzz/wsVLl69sXQ2uXb9x81Zn+/aBVaWhbEiVUOYoI5YJXrAhcBDsSBtGZCbYYXb8utEPT5mxXBWfYa7ZWJK84FNOCTgq7fyYP0gBv8SJnvE0wi6rYDeqe4s8bvPY5RMFFveWPt3yum6YU2KYtly4ctBLYMaApNFZtq24VOINZVV7pZ9s6Cc1TpJgrQlOLJf4/U7YdyCX5Ev8MO10w0G4CHweRC3oojb2023vhWtLS8kKoIJYO4pCDeOKGOBUsDpISss0occkZyMHCyKZHVeLhdf4vmMmeKqMOwXgBXv2RkWktXOZOackMLObWkP+SxuVMH0+rnihS2AFXTaalgKDws3v4Qk3jIKYO0Co4e6tmM6IIRTcH691aWqDUsK6Ud4wN6Jh7xz1QTNDQJlHVUJMLnlRu5HzpN+g/xnJ15XRoSBw+442t3seHMSD6OngycfH3b1X7ea30F10D+2gCD1De+gt2kdDRL2+98kbeYn/zf/u//R/La2+1965g9bC//0XDxz9Yg==</latexit> y0 t = 1y0 t 1 + 2y0 t 2 + · · · + py0 t p + "t + ✓1"t 1 + ✓2"t 2 + · · · + ✓q"t q "t ⇠ N(0, 2) <latexit sha1_base64="l+KCMW9MB/0Fx6iypyW6iAl35q4=">AAACZnicfVDRShtBFJ1stWqsGi3igy9DQ1FKDbuibR8siPrQF1HBqJCE5e7kJg7Oziwzd4vLkl/xtf6Sf9DP6GyMULX0wMDh3HPn3nuSTElHYfhQC95MTb+dmZ2rz79bWFxqLK9cOJNbgW1hlLFXCThUUmObJCm8yixCmii8TG4Oq/rlT7ROGn1ORYa9FIZaDqQA8lLcWCk2YuLfeRHTVhGXtBWN4kYzbIVj8NckmpAmm+A0Xq7tdftG5ClqEgqc60RhRr0SLEmhcFTv5g4zEDcwxI6nGlJ0vXK8/Ih/9EqfD4z1TxMfq393lJA6V6SJd6ZA1+5lrRL/VevkNPjWK6XOckItHgcNcsXJ8CoJ3pcWBanCExBW+l25uAYLgnxez6ZUf5MxyvlTjtCfaPHYSycZWiBjP5VdsMNU6pE/edj9XLH/GeH2yehZve7zjl6m+5pcbLeiL63ds53m/sEk+Vm2zj6wTRaxr2yf/WCnrM0Eu2V37Be7r/0OFoPVYO3RGtQmPe/ZMwT8D0fHuic=</latexit> y0 t = yt yt 1 1階の階差 <latexit sha1_base64="kmP7xlM50fDzWLAhTNWxWvGCUhU=">AAACbHicfVBdSxtBFJ2s2mq0bdQ+FEQYGtpIqWG39OuhBak+9EW00KiQpMvdyU0cnJ1ZZu5Kl2X9NX1t/0//RH9DZ2MK9QMvzHA499y5c06SKekoDH83grn5hXv3F5eayysPHj5qra4dOZNbgT1hlLEnCThUUmOPJCk8ySxCmig8Ts526/7xOVonjf5KRYbDFCZajqUA8lTcelJ8KzudKib+kRedmLb9VdJ2VMWtdtgNp8VvgmgG2mxWh/Fq48NgZESeoiahwLl+FGY0LMGSFAqr5iB3mIE4gwn2PdSQohuWUwsVf+aZER8b648mPmX/nyghda5IE69MgU7d9V5N3tbr5zR+PyylznJCLS4XjXPFyfA6Dz6SFgWpwgMQVvq/cnEKFgT51K5sqd8mY5TzVvbQW7S476mDDC2QsS/KAdhJKnXlLU8GL2t0lxC+/xN61Gz6vKPr6d4ER6+60dvumy+v2zufZskvsg32lG2xiL1jO+wzO2Q9JtgF+8F+sl+NP8HjYCPYvJQGjdnMOrtSwfO/yfS8Lg==</latexit> y00 t = y0 t y0 t 1 2階の階差 ARIMA(p,d,q)と表現する dは階差の回数 Q&A Q1. p,qはどうやって選べばいいの？ A1. Rだとauto.arima関数でAICベースでやってくれるよ Q2.予測できる？ A2. forecast関数で簡単にできるよ
10. ここまでは、yt (⾃分の値)だけを 使ってきたけど、共変量も⼊れ たいよ。だって疫学だもの。 10/48 時系列回帰〜
11. 単位根過程と⾒せかけの回帰 単位根回帰 AR(p)モデルで特性⽅程式のd個の根が1となる場合にその確率過程はd次の和分過程I (d )とする 簡単に⾔えば： yt は⾮定常だけど、d回の階差をとれば定常になる時系列のこと ⼤問題：⾒せかけの回帰 (Spurious regression) y とx が共に単位根のとき、ほぼ確実に有意なβが得られる 11 <latexit sha1_base64="m3DRxUPTUOWFkf2jz7UlAXZwjJc=">AAACc3icfVBdSxtBFJ2sVm1sa2wf89DBKJS2hl3R2gcLoj74IlVoVMiGcHdykwzOziwzdyVhyaO/pq/tj+kP8b2zMQU/ihcGDuecmTvnJJmSjsLwTyWYm3+xsLj0srr86vWbldrq23NnciuwJYwy9jIBh0pqbJEkhZeZRUgThRfJ1WGpX1yjddLoHzTOsJPCQMu+FECe6tbej7vEv/E4QQI+8vgTj6/BYuak8jp1a42wGU6HPwXRDDTYbE67q5W9uGdEnqImocC5dhRm1CnAkhQKJ9U4d5iBuIIBtj3UkKLrFNMkE77hmR7vG+uPJj5l798oIHVunCbemQIN3WOtJP+ntXPqf+0UUmc5oRZ3i/q54mR4WQvvSYuC1NgDEFb6v3IxBAuCfHkPtpRvkzHK+ShH6CNaPPHU9wwtkLEfixjsIJV64iMP4s8les4Io39Gj6pV33f0uN2n4HyrGX1p7pxtN/YPZs0vsTpbYx9YxHbZPjtmp6zFBLthP9kv9rtyG9SDtWD9zhpUZnfesQcTbP4FQvi/pw==</latexit> yt = xt + "t 例えば、2つの⾮定常な時系列があり、それぞれが時間ととも に増加する傾向にあるとする。これらの系列をそのまま回帰分 析すると、系列間に強い正の相関があるように⾒えるかもしれ ない。しかし、これは⾒せかけの回帰。 なんにも考えずに、時系列デー タ同⼠を回帰してみたらなんか 有意な結果得られた！やったね ！！！！！ <latexit sha1_base64="A9owpiGwcWp/Z26eLujI1d+PUN0=">AAACmHicfVBdb9MwFHUzPkb56rY3eLGokAaUKplg2wNIFUwCHoCB6DapKdGNe5tas+PIvpmIov4hfg2v49fgdEWwDXEtS0fnnuvrc9JCSUdheNoKVq5cvXZ99Ub75q3bd+521tYPnCmtwKEwytijFBwqmeOQJCk8KiyCThUepsevm/7hCVonTf6FqgLHGrJcTqUA8lTS2asS4i95XMwkr5KankbzJ/EJWCycVF5APR435w/FYyc1/7AZ9jzINHzdepR0umE/XBS/DKIl6LJl7SdrrRfxxIhSY05CgXOjKCxoXIMlKRTO23HpsABxDBmOPMxBoxvXC7tz/tAzEz411t+c+IL9e6IG7VylU6/UQDN3sdeQ/+qNSprujmuZFyVhLs4WTUvFyfAmOz6RFgWpygMQVvq/cjEDC4J8wue2NG+TMcp5K3voLVp876mPBVogYx/XMdhMy3zuLWdxr0H/E8K330KP2m2fd3Qx3cvgYKsfbfeff3rWHbxaJr/K7rMHbJNFbIcN2Fu2z4ZMsO/sBztlP4N7wSB4E7w7kwat5cwGO1fB5188jMyW</latexit> yt = yt 1 + "t, " AR(1)でρ=1 全く関係ないランダムウォークxとy を回帰すると係数が何故か、有意に なっちゃう!!
12. じゃあ、どうすんの？→ 検定しようぜ 単位根過程かどうかの検定は沢⼭ある Dickey-Fuller (DF)検定、Augmented DF検定、Phillips‒Perron検定、KPSS検定 ADF検定が最もよく使われる AR(p)を仮定し、H0: γ=0 (単位根あり) vs H1: γ<0 (単位根なし)を検定する （基本的）⼿順： 1. 元データyt でADF検定をする 2. 階差列をつくってADF検定する 3. １へ戻る (有意差ありになるまでコレを繰り返す) 12/48 有意差あり (= 単位根なし) → OK 有意差なし (= 単位根あるかも) → Step 2に移⾏
13. 系列相関 (Serial correlation) 定義 回帰分析において残差(実測と理論値の差）に⾃⼰相関があること 普通の分析⼿法はサンプルのi.i.d.(独⽴同分布)を仮定しているが、時系列 データはそうではないよね。。。(前のデータに今期のデータが依存) 統計的には、OLSがBLUE（最良線形不偏推定量）ではなくなる Durbin-Watson検定 【使い⽅】回帰分析した後に、残差に系列相関があるかどうかを検定 13/48 <latexit sha1_base64="kaXuX88KXdbGiKG5X7UfxyANKt4=">AAAClnicfVBdixMxFE3Hr3X86uqL4EuwuIhdy4z49eDKsir6Ilawu4WmDHfS2zZsJhmSO4tl6P/x1/iq/hszbQV3V7wQcjj33Juck5daeUqSX63owsVLl69sXY2vXb9x81Z7+/aht5WTOJBWWzfMwaNWBgekSOOwdAhFrvEoP37T9I9O0HllzRdalDguYGbUVEmgQGXtg51FRnyPixwJsqS7vtNhRl1eZbTLhYh3AtgTbm4DU9PjdNkVJ+Cw9EqHFZS1O0kvWRU/D9IN6LBN9bPt1isxsbIq0JDU4P0oTUoa1+BISY3LWFQeS5DHMMNRgAYK9ON6ZXbJHwRmwqfWhWOIr9i/J2oovF8UeVAWQHN/tteQ/+qNKpq+HNfKlBWhkeuHppXmZHmTHJ8oh5L0IgCQToW/cjkHB5JCvqdeaXaTtdoHK28xWHT4MVCfSnRA1j2qBbhZocwyWJ6J3Qb9Twhf/wgDiuOQd3o23fPg8Ekvfd579vlpZ/9gk/wWu8fus4csZS/YPvvA+mzAJPvGvrMf7Gd0N3odvYver6VRazNzh52qqP8bz/vLnQ==</latexit> yt = 0 + 1Xt + ut, ut = ⇢ut 1 + "t H0: ρ=0 (系列相関なし) vs H1: ρ≠0 (系列相関あり) 注意: 有意差あれば系列相関あるよってシグナル <latexit sha1_base64="yDLZ9zeJvEKYK/Q9FRfznr0zJu8=">AAACkHicfVBdaxNBFJ1s/WjjV1offRkMQhUbdoO2fWixah9EECs0TSGbLrOT2XTofCwzd4ph2T/TX+OrvvlvvJtGsa14YeDMOffOnXPyUkkPcfyzFS3dun3n7vJK+979Bw8fdVbXjrwNjosBt8q645x5oaQRA5CgxHHpBNO5EsP87H2jD8+F89KaQ5iVYqzZ1MhCcgZIZZ2d/SHdpemkcIxXqQ86q2C3X58c0vVZBhsBrxtJ/fykX/9RE1RDBkhlnW7ci+dFb4JkAbpkUQfZamsnnVgetDDAFfN+lMQljCvmQHIl6nYavCgZP2NTMUJomBZ+XM1t1vQZMhNaWIfHAJ2zf09UTHs/0zl2agan/rrWkP/SRgGK7XElTRlAGH65qAiKgqVNZnQineCgZggYdxL/Svkpw8QAk72ypXkbrFUerewLtOjEJ6Q+l8IxsO5FlTI31dLUaHmavmzQ/xrZ19+NiNptzDu5nu5NcNTvJZu9119edffeLZJfJk/IU7JOErJF9sgHckAGhJML8o18Jz+itWg7ehO9vWyNWouZx+RKRR9/AYiOyk4=</latexit> DW = PT t=2 (yt ut 1)2 PT t=1 u2 t Durbin-Watson⽐ 0-4までの値をとり0に近いほど 系列相関が強い
14. 中断時系列デザイン (Interuppeted time serise; ITS) 介⼊のアウトカムへの影響を時系列のレベル（切⽚）とトレンド（傾き）の 変化で定量化する⽅法 ITSの前提条件 (Kontopantelis et al. (2015, BMJ)) 1. Pre-interventionとpost-interventionで集団の分布が同⼀ 2. Pre-interventionの時系列変化が線形 3. 介⼊以外に時系列の変化を説明できるものがない 14/48 トレンド の変化 レベルの 変化 <latexit sha1_base64="4FdMjk6IY1LHzFlGnWuGSviQTt0=">AAACqHicfVBdaxNBFJ2sVWv8aFoffRkMgrQSdutXHxSK9sEXMULSVLJhuTu5SYfOxzJztxiW/DH/iW++6q9wNtmibcULwxzOOXfu3JMXSnqK4++t6MbGzVu3N++07967/2Crs71z7G3pBA6FVdad5OBRSYNDkqTwpHAIOlc4ys/e1/roHJ2X1gxoUeBEw9zImRRAgco6g0VG/C1PcyTI4r31nQwasM/71tNRqYNprzE95wOektTo/4hBOweHhZcqPEpZpxv34lXx6yBpQJc11c+2W2/SqRWlRkNCgffjJC5oUoEjKRQu22npsQBxBnMcB2ggjJ9Uq/WX/ElgpnxmXTiG+Ir9u6MC7f1C58GpgU79Va0m/6WNS5odTCppipLQiPWgWak4WV5nyafSoSC1CACEk+GvXJyCA0Eh8UtT6rfJWuXDKkcYVnT4MVCfCnRA1u1WKbi5lmYZVp6nz2r0PyN8vTAG1G6HvJOr6V4Hx/u95FXv5ecX3cN3TfKb7BF7zJ6yhL1mh+wD67MhE+wb+8F+sl/RbtSPRtGXtTVqNT0P2aWK8t9BitL+</latexit> yt = 0 + 1T + 2PostDumt + 3T ⇥ PostDumt + "t Bernal et al. (2016, IJE) Preのトレンド レベルの変化量 トレンドの変化量 レベルだけ変化 β3 = 0 トレンドだけ変化 β2 = 0 レベルとトレンド 両⽅変化 β2 ≠ 0 β3 ≠ 0
15. ITS (Cont.) 介⼊の効果は様々な形で表現可能 関数系は⾃分で決めてOK(是⾮統計家に相談してください) 15/48 Lagを付きで介⼊ 効果が出る ⼀時的なレベル 変化 ⼀時的なトレンド の変化とそれに続 くレベルの変化 Q&A Q1: 前後、何時点とればいいん です？ A1: 多ければ多いほどパワーは 増⼤。それぞれ最低8時点づつ？ (Penfold and Zhang, 2013) Q2: じゃあ、多ければ多いほどい いんですね！？ A2: ITSの前提である、期間中の 分布が同じという仮定が成り⽴た なくなっちゃうので注意！
16. ITSの実際の解析⼿順 Jandoc et al. (2015, JCE) のガイドラインに沿って報告する 実際の⼿順（簡略版） 1. 研究期間と解析前のタイムポイントを定義する 2. 介⼊がどういったモデルで表現できるか考える 3. 介⼊効果が漸進的か遅延している場合はラグ効果も考慮する 4. 残差にACFやDW検定などで系列相関がないか確かめる 5. もしあれば、季節成分などを含めて、もう⼀度、系列相関check 6. Sub-analysisとして、介⼊を受けないアウトカムや集団で同じモデルをチェ ックしてみる。 季節性調整について シンプルな⼿法：⽉や曜⽇ダミーを⼊れる 頻出⼿法：フーリエ項を⼊れる 16/48 <latexit sha1_base64="MscMu9VhxcTlCqCcH2SJXnVkbXY=">AAADDnicjVFNa9tAEF2pH0ndjzjNsZelppCWYCSTtDn0ENIQekniQp0ELGNW65G8ZLUrdkehRug/9Nf0UAi55i/k33TluNDYpu3AwuPNm5mdN3EuhcUguPX8Bw8fPV5ZfdJ4+uz5i7Xm+stTqwvDoce11OY8ZhakUNBDgRLOcwMsiyWcxRef6vzZJRgrtPqKkxwGGUuVSARn6Khh86cVKpKQ4GY0SgzjZSfKBcWqDDtVZEQ6xrdblGv7T818n+3/6LNME4002iUjjxe0w2YraAfToIsgnIEWmUV3uO59dN15kYFCLpm1/TDIcVAyg4JLqBpRYSFn/IKl0HdQsQzsoJyaXNE3jhnRRBv3FNIp+2dFyTJrJ1nslBnDsZ3P1eSyXL/AZHdQCpUXCIrfDUoKSVHT+mJ0JAxwlBMHGDfC/ZXyMXOmoLvrvSl1b9RaWrfKAbgVDRw56iQHw1Cbd2XETJoJVbmV02irRn8Tsm+/hQ41Gs7vcN7dRXDaaYfv2ztftlt7+zPnV8kr8ppskpB8IHvkM+mSHuEe9Q69E6/rf/d/+Ff+9Z3U92Y1G+Re+De/AC1l+WQ=</latexit> sin ✓ 2⇡t 12 ◆ , cos ✓ 2⇡t 12 ◆ , sin ✓ 4⇡t 12 ◆ , cos ✓ 4⇡t 12 ◆ , . . . , cos ✓ 2N⇡t 12 ◆ ⽉ごとなら12。週ごとなら52など ⽉次データなら1⽉、２⽉とか。週次データなら第23週とか
17. 状態空間モデル 観測値のモデルと、直接観測できない潜在的な「状態」のモデルを同時 に使う 何が嬉しいの？→とても柔軟にモデリング可能 (⽋損値OK，ARIMAなども表現可能、定常性の仮定がいらない) But, 推定がしばしば⼤変。。。 線形ガウスとかだとKalman-Filterみたいな⾼速な⼿法がある ⾊々複雑にするとMCMC 17/48 <latexit sha1_base64="/OChslWUmHqP33uS95+if07Ws0k=">AAACgHicfVBdaxQxFM2OH63j19Y++hJcLMXWdaa0WkShqGhfxApuW9hZwp3s3W1oJhmSO3WXYf+Ev8ZX/Rf+GzPbFWwrXggczjnJzTl5qZWnJPnViq5dv3FzaflWfPvO3Xv32ysPDr2tnMSetNq64xw8amWwR4o0HpcOocg1HuWnbxv96AydV9Z8oWmJgwLGRo2UBAqUaG+uTQTx1/y9oImo6Wk62/gg6ExQlsVr07m030i08VWQaHeSbjIffhWkC9BhizkQK61X2dDKqkBDUoP3/TQpaVCDIyU1zuKs8liCPIUx9gM0UKAf1PNYM/44MEM+si4cQ3zO/n2jhsL7aZEHZwF04i9rDfkvrV/RaHdQK1NWhEaeLxpVmpPlTUd8qBxK0tMAQDoV/srlCTiQFJq8sKV5m6zVPkR5hyGiw4+B+lSiA7LuSZ2BGxfKzELkcbbZoP8ZYfLHGFAch77Ty+1eBYdb3fR5d+fzdmfvzaL5ZfaQPWLrLGUv2B7bZwesxyT7xr6zH+xnFEXr0bMoPbdGrcWdVXZhope/AdOjwuM=</latexit> xt = Ftxt 1 + Gtvt yt = Htxt + wt システムモデル (観測されない) 観測モデル (観測されない) 細かいことを話そうと思 ったけど、⻑くなりそう なんで今回は勘弁
18. 時系列データ同⼠の関係 時系列同⼠の距離(類似度)をどう測るか？ 1. CCF (cross-correlation function): ⼀番シンプル Window を制限して移動させていくとCCFの時間変化も⾒れる 18/48 <latexit sha1_base64="C/hCAEH+5l8I9M6Fk/Ax7smQ4DQ=">AAACdnicfZBdaxNBFIYn61cbv1K9LMjQUKxawq74ddFCsV54I1Zo2kg2LmcnJ+mQ+VhmzkrDsvf9Nd7qX/GfeOlsGsG24oGBh/e8M2fOmxdKeorjn63o2vUbN2+trLZv37l7735n7cGRt6UT2BdWWTfIwaOSBvskSeGgcAg6V3icz/ab/vFXdF5ac0jzAkcapkZOpAAKUtbZ2M+qwed6a/aE7/LUlzqraDepvxzOMzoN/GxWZ51u3IsXxa9CsoQuW9ZBttbaScdWlBoNCQXeD5O4oFEFjqRQWLfT0mMBYgZTHAY0oNGPqsUyNd8MyphPrAvHEF+of9+oQHs/13lwaqATf7nXiP/qDUuavBlV0hQloRHngyal4mR5kwwfS4eC1DwACCfDX7k4AQeCQn4XpjRvk7XKh1XeYVjR4YcgfSzQAVn3tErBTbU0dVh5mm439D8jnP4xBmq3Q97J5XSvwtHzXvKq9/LTi+7e22XyK2ydbbAtlrDXbI+9ZweszwQ7Y9/Yd/aj9St6FG1Gj8+tUWt55yG7UFH8GzqhwQE=</latexit> CXY (k) = T X t=1 ytxt+k 時系列y の時刻t のときの値と時系列x の時刻t+k のときの値の相関
19. Dynamic time warping: DTW CCFの問題（の⼀つ）は、⻑さの違う時系列では計算できない DTW: ようは距離 (⼤きい値だと遠い) 時系列データx とy の各要素を距離が⼩さくなるようにアライメント 19/48
20. 時系列のクラスタリング 時系列間の距離(類似度)が定義できるとそれ使⽤してクラスタリング可能 特徴ベースアプローチ ⽣の時系列データをそのまま従来のクラスタリング⼿法に適⽤ Ex. DTWを使ったk-means 形状ベースアプローチ ⽣の時系列を低次元の特徴ベクトルに変換し，抽出した特徴ベクトルに対して従来のク ラスタリング⼿法を適⽤ Ex. K-shapelet モデルベースアプローチ クラスターごとにモデルを仮定し， それらのモデルをデータに適合させる ことでクラスタリングを⾏う Ex. 混合ARIMA過程 20/48
21. コヒーレンス (Coherence) 周波数領域での相関 特に周期性同⼠の類似性に興味があるとき 21/48 <latexit sha1_base64="QD0lYv6vzdh4XXI655UsCIZxtU4=">AAAClnicfVDbahsxEJW3t9S9Oe1LoS+ippCGYnZNm/ShDSFJSV5KXagTB69ttPKsLaKVFkkbYpT9n3xNXpv+TWcdB3IpHRBzdOaMRnOSXArrwvBPLbh3/8HDR0uP60+ePnv+orH8ct/qwnDoci216SXMghQKuk44Cb3cAMsSCQfJ0XZVPzgGY4VWv9wsh0HGJkqkgjOH1Kixta2nw/bIn8zKlfQ9/UrjcWoY96c9vA79anmI+XTYLv0VUeXDa8Vy1GiGrXAe9C6IFqBJFtEZLde+xGPNiwyU45JZ24/C3A08M05wCWU9LizkjB+xCfQRKpaBHfj5siV9h8yYptrgUY7O2esdnmXWzrIElRlzU3u7VpH/qvULl34eeKHywoHil4PSQlKnaeUcHQsD3MkZAsaNwL9SPmVolkN/b0yp3nZaS4ur7ACuaOA7Uj9yMMxps+pjZiaZUCWuPIk/VOh/QnZyJURUr6Pf0W1374L9ditaa336+bG5ubVwfom8IW/JConIOtkke6RDuoSTM3JOfpOL4HWwEXwLdi+lQW3R84rciKDzF2ZSytQ=</latexit> Coh2 xy (f) = |X(f)⇤Y (f)|2 X(f)⇤X(f)Y (f)⇤Y (f) パワースペクトル 周波数ごとの波の強さ クロススペクトル xt のフーリエ変換をX (f ) X* (f )は複素共役 <latexit sha1_base64="aMJni/vgnnfFk/WajU+CHwZKRdI=">AAAChHicfVDbbhMxEHW2UEq4NC2PvFhESClqo93S2wOgCnjgBVEk0kbKptGsM5ta9dorexY1Wu1v8DV9bf+hf4M3CRJtESNZPjpzxuNzklxJR2F40wiWHjxcfrTyuPnk6bPnq6219WNnCiuwJ4wytp+AQyU19kiSwn5uEbJE4Uly/qnun/xE66TRP2ia4zCDiZapFECeGrXCfifd4O95LDWNyi1/pTStTuc3x9NyS27HueQpVRcd2hjTqNUOu+Gs+H0QLUCbLepotNZ4F4+NKDLUJBQ4N4jCnIYlWJJCYdWMC4c5iHOY4MBDDRm6YTmzVvHXnhnz1Fh/NPEZ+/dECZlz0yzxygzozN3t1eS/eoOC0oNhKXVeEGoxX5QWipPhdU58LC0KUlMPQFjp/8rFGVgQ5NO8taV+m4xRzlv5jN6ixa+e+pajBTL2TRmDnWRSV97yJN6s0f+EcPFH6FGz6fOO7qZ7Hxxvd6O97u73nfbhx0XyK+wle8U6LGL77JB9YUesxwT7xS7ZFbsOloPN4G2wO5cGjcXMC3argg+/AR2RxVE=</latexit> X(f) = Z 1 1 e i2⇡ftx(t)dt 80HzのCos波と 120HzのSin波の合成 80HzのCos波と 200HzのSin波の合成 うんうん、80と 120のところにき ちんとピーク⽴っ てるね 80Hzが被っている うんうん、（カブ ってる）80のとこ ろにきちんとピー ク⽴ってるね Fourier変換なのでSin波(とCos 波)の重ね合わせで表現できる 波を考える
22. ウェーブレットクロススペクトル解析 （ざっくり⾔うと）時間により周波数成分が変化する場合の解析⼿法 Fourieのsin波とcos波じゃなくもっと⼩さい波(wavelet)の重ね合わせで表現すればいいじゃない パラメータ: a が周波数帯域(スケール)に関連、b が時間(シフト)に関連 スカログラム: a、b、W をplotしたもの 22/48 <latexit sha1_base64="G44uK+mP+1/j6eZHjrK2mBTUdRw=">AAACrnicfVDdbtMwGHXD3yh/HVxyY1EhtdOoEsQGFyBNwAU3iCHRdqLOIsdxUmuOHewvqJWVd+M1eAFu4RFw2iCxDfFJto/Od44/+6SVFBbC8HsvuHL12vUbOzf7t27fuXtvsHt/ZnVtGJ8yLbU5SanlUig+BQGSn1SG0zKVfJ6evWn786/cWKHVJ1hXPC5poUQuGAVPJYPP82Q1ovvpGL/CJMsNZS5qHLFfDDjaNEQoSJzfc1g3p9sTr0YwJpUVp3tE8hxGnQ+epI33ECOKJYwzSAbDcBJuCl8GUQeGqKvjZLf3kmSa1SVXwCS1dhGFFcSOGhBM8qZPassrys5owRceKlpyG7tNCA1+7JkM59r4pQBv2L8djpbWrsvUK0sKS3ux15L/6i1qyF/ETqiqBq7YdlBeSwwat4niTBjOQK49oMwI/1bMltQnAj73c1Pau0Fraf1X3nL/RcPfe+pDxQ0FbfYcoaYohfIxmoLst+h/Qrr6I/So3/d5RxfTvQxmTyfR4eTg47Ph0esu+R30ED1CIxSh5+gIvUPHaIoY+oZ+oJ/oVxAGsyAOkq006HWeB+hcBcvf3zLYIw==</latexit> Wx(a, b) = 1 p a Z 1 1 x(t) ⇤ ✓ t b a ◆ dt マザーウェーブレット: 試⾏錯誤して選んでく ださい トランプ⽀持者のインスタ投稿数のスカログラム a b <latexit sha1_base64="PLeKBUYMkH4TTR/V1X29I6e4F5g=">AAACeHicfVDLThsxFHWGPiD0EWDJxmqKeAhFM4jXAiREu+imKpUIQUrS6I5zEyw89si+g4hG8wN8DVv4k/5KV/WEqVSg6pUsnXvOvbbPiVMlHYXhz1ow8+Llq9ezc/X5N2/fvW8sLJ45k1mBbWGUsecxOFRSY5skKTxPLUISK+zEl59KvXOF1kmjT2mSYj+BsZYjKYA8NWh87Azy60mxBpvxOj/kZVf82Ji2HlfCoNEMW+G0+HMQVaDJqjoZLNQOekMjsgQ1CQXOdaMwpX4OlqRQWNR7mcMUxCWMseuhhgRdP5/aKfiKZ4Z8ZKw/mviU/Xsjh8S5SRL7yQTowj3VSvJfWjej0X4/lzrNCLV4eGiUKU6Gl9nwobQoSE08AGGl/ysXF2BBkE/w0Svl3WSMct7KZ/QWLX711LcULZCxG3kP7DiRuvCWx73NEv1vEK7/DHpUr/u8o6fpPgdnW61ot7Xzfbt5dFwlP8uW2Qe2xiK2x47YF3bC2kywG3bL7th97VfAg9Vg/WE0qFU7S+xRBVu/AXrCwIE=</latexit> Wxy(a, b) = W⇤ x (a, b)Wy(a, b) ウェーブレットク ロススペクトル クリントン vs トランプ⽀持者のインスタ投稿数のクロススペクトル
23. ヒルベルト変換 Sin波とcos波は相関=1 でも位相はπ/2だけずれている 時系列データから位相を計算するにはヒルベルト変換を⽤いればOK 23/48 どれくらいズレてるん ですか？
24. ⽅向データ 普通の平均などが⾃然な定義とならない！ 23時に寝た時と1時に寝た時、平均⼊眠時刻は(23+1)/2= 12時？ ⾓度データ[-π, π)に変換して統計解析(することが多い) 24/48 1987~1998 年に、ピッツバーグで銃犯罪が発⽣した時刻 (Gill & Hangartner, 2010) ⼀様に分布しているか？→レイリー検定 24 時間を[−π, π) の⾓度 として表している。 例えば0は0時、πは正午に 対応 1/1 6/15 12/31 例えばカレンダーに対して 回帰する。1/1の次は12/31 なんで、実は両端は⼀致し てほしい。。。 円周上S1の分布: フォン・ミーゼス分布、wrapped正規分布、wrappedコーシー分布 球⾯上Sdの分布 : フォン・ミーゼス・フィッシャー分布、フィッシャー‒ビングハム分布、ケント分布 <latexit sha1_base64="CsrNQjMwwwBwzG4S97YyYoMYzHs=">AAACl3icfVBdSxtBFJ1s1dq09aM+FV8GQ6FoCbtitQ8VpYr4IlVo/CAbwt3JTTI4O7PM3BXDkh/UX9PX+m+cjStELb0wcObcc+/MOUmmpKMwvKsFr2Zm517Pv6m/ffd+YXFp+cO5M7kV2BJGGXuZgEMlNbZIksLLzCKkicKL5Pqg7F/coHXS6F80yrCTwkDLvhRAnuouHVzxXR4nSNAN+UaFIi6M4zEN/WWKdFJPkTdgMXNSlVsaYTOcFH8Jogo0WFWn3eXa97hnRJ6iJqHAuXYUZtQpwJIUCsf1OHeYgbiGAbY91JCi6xQTt2P+yTM93jfWH018wk5PFJA6N0oTr0yBhu55ryT/1Wvn1P/WKaTOckItHh7q54qT4WV0vCctClIjD0BY6f/KxRAsCPIBP3ml3E3GKOetHKK3aPHEUz8ztEDGrhcx2EEq9dhbHsRfSvQ/Idw+Cj2q133e0fN0X4LzzWa03fx6ttXY/1ElP89W2Rr7zCK2w/bZMTtlLSbYb/aH/WV3wcdgLzgKjh+kQa2aWWFPKji7Bzj9y0c=</latexit> Y = 0 + 1cos✓ + 1sin✓ + " Mardia (1976)
25. VAR Vector Auto-regressive model (VAR): 複数時系列のAR(p )モデル 多くの時系列を同時に考慮したほうが精度良さそうでしょ？ 過去の⾃分だけでなく、他の時系列データの過去も考慮されている 推定：個別の回帰式をOLS(⼀致性アリ、これがVARが流⾏した理由の⼀つ by 沖本本) ⼀般のSUR(Seemingly unrelated regression)ではダメ (VARは全⽅程式が同じ共変量の組を持つから特別にOK) ワクチンと超過死亡と移動量と… など幾つでもOK But 多すぎると推定が不安定になる (n変量VAR(p)はn(np+1)+n(n+1)/2個のパラメータ) 定常条件などはAR過程のときと同じ様なもの（Yule-Walker⽅程式） 25/48 <latexit sha1_base64="1M02/Se5cuOa50QBOQCkrnI1VOQ=">AAAC93icfVHbitRAEO3E2268zeqjL42Di+g6poOXFRQW9cEXcQVnd2EyhEpPzUyznXTo7iyGkG/xQRBf/Rz/xs4ki3vTgobDqVOnuqrSQgpjw/C351+6fOXqtbX14PqNm7duDzbu7BlVao5jrqTSBykYlCLHsRVW4kGhEbJU4n56+K7N7x+hNkLlX2xV4DSDRS7mgoN1VDL4vlklNduyDX1DY5DFEpLwcZyihZZmTZd9wpq/ZNSSUU8egcbCCOm8Ops4psFmJzjhyY7Lo4s8o397tjbJYBiOwlXQ84D1YEj62E02vNfxTPEyw9xyCcZMWFjYaQ3aCi6xCeLSYAH8EBY4cTCHDM20Xm2zoQ8cM6Nzpd3LLV2xJytqyIypstQpM7BLczbXkhflJqWdb09rkRelxZx3jealpFbR9jR0JjRyKysHgGvh/kr5EjRw6w54qkvrbZWSxo3yHt2IGj866lOBGqzSj+oY9CITeeNGXsRbLfqfEL4eCx0KArdvdna758FeNGIvRs8/PxvuvO03v0bukfvkIWHkJdkhH8guGRPurXtPvW3vlV/53/wf/s9O6nt9zV1yKvxffwC34uxA</latexit> y1,t = ↵0 + 1,1y1,t 1 + 1,2y2,t 1 + "1,t y2,t = ↵1 + 2,1y1,t 1 + 2,2y2,t 1 + "2,t 現在の値も含めたSVARという モデルもある (Sims, 1980) ⾒せかけの回帰のときと同じで、単位根同⼠を回帰しちゃうとダメ！ Q. じゃあどうすんの？ A. (前説明した通り)検定して差分をとってればOK。でも共和分関係にあるとVARダメってい うトラップも存在
26. グレンジャーの因果“性” 普通の因果ではない！！ (グレンジャーの意味での因果、予測的因果) ある時系列y の未来の予測性能(MSE)がx を考慮したときの⽅がうまくい くこと グレンジャーの因果性検定y2 → y1 : かどうか？(残差のF検定) 26/48 『⽇本でのCOVID-19禍におけるPCR陽性者数と感染 死亡者数のGranger因果性の検討』濱⽥悦⽣,2021 <latexit sha1_base64="1M02/Se5cuOa50QBOQCkrnI1VOQ=">AAAC93icfVHbitRAEO3E2268zeqjL42Di+g6poOXFRQW9cEXcQVnd2EyhEpPzUyznXTo7iyGkG/xQRBf/Rz/xs4ki3vTgobDqVOnuqrSQgpjw/C351+6fOXqtbX14PqNm7duDzbu7BlVao5jrqTSBykYlCLHsRVW4kGhEbJU4n56+K7N7x+hNkLlX2xV4DSDRS7mgoN1VDL4vlklNduyDX1DY5DFEpLwcZyihZZmTZd9wpq/ZNSSUU8egcbCCOm8Ops4psFmJzjhyY7Lo4s8o397tjbJYBiOwlXQ84D1YEj62E02vNfxTPEyw9xyCcZMWFjYaQ3aCi6xCeLSYAH8EBY4cTCHDM20Xm2zoQ8cM6Nzpd3LLV2xJytqyIypstQpM7BLczbXkhflJqWdb09rkRelxZx3jealpFbR9jR0JjRyKysHgGvh/kr5EjRw6w54qkvrbZWSxo3yHt2IGj866lOBGqzSj+oY9CITeeNGXsRbLfqfEL4eCx0KArdvdna758FeNGIvRs8/PxvuvO03v0bukfvkIWHkJdkhH8guGRPurXtPvW3vlV/53/wf/s9O6nt9zV1yKvxffwC34uxA</latexit> y1,t = ↵0 + 1,1y1,t 1 + 1,2y2,t 1 + "1,t y2,t = ↵1 + 2,1y1,t 1 + 2,2y2,t 1 + "2,t <latexit sha1_base64="jGQPIENCAp577xIhCzmekfpARhY=">AAACZHicfVBdaxNBFJ1sta0x1rTFJ0EGgyASwm6xHw8tFPXBF7GCaQPZEO5ObpIhszPLzN3SsOSf+Gr/U/+Av8PZNIJJSg8MHM49d+69J8mUdBSGd5Vg48nTza3tZ9XntRc7L+u7e5fO5FZgWxhlbCcBh0pqbJMkhZ3MIqSJwqtk8rmsX12jddLonzTNsJfCSMuhFEBe6tfrcYIE/SJqHsz4GQ95v94IW+EcfJ1EC9JgC1z0dyun8cCIPEVNQoFz3SjMqFeAJSkUzqpx7jADMYERdj3VkKLrFfPVZ/ydVwZ8aKx/mvhc/b+jgNS5aZp4Zwo0dqu1Unyo1s1peNIrpM5yQi3uBw1zxcnwMgc+kBYFqaknIKz0u3IxBguCfFpLU8q/yRjl/Clf0J9o8ZuXvmdogYz9UMRgR6nUM3/yKG6W7DEj3Pwzelat+ryj1XTXyeVBKzpqHf742Dj/tEh+m71mb9l7FrFjds6+sgvWZoJds1/sN7ut/AlqwX7w6t4aVBY9+2wJwZu/S7e4sw==</latexit> 1,2 = 0 <latexit sha1_base64="Hp3H0To/9FR5B1otMHdCjolmtnI=">AAACf3icfVBdTxNBFJ2uqFi/ijz6MqExFrVl16hgIglBY3wxolAg6TbN3endMmF2ZjNzl9hs9kf4a3iFn8G/cbbURMB4k5mcnHvu3DknyZV0FIYXjeDWwu07dxfvNe8/ePjocWvpyb4zhRXYF0YZe5iAQyU19kmSwsPcImSJwoPk+GPdPzhB66TRezTNcZjBRMtUCiBPjVovP/NNHo9TC6Ls7O7+GIXd+o5W12xVztBaZ6+r8260Wo1a7bAXzorfBNEctNm8dkZLjQ/x2IgiQ01CgXODKMxpWIIlKRRWzbhwmIM4hgkOPNSQoRuWM1cVf+aZMU+N9UcTn7F/T5SQOTfNEq/MgI7c9V5N/qs3KCjdGJZS5wWhFpeL0kJxMryOiI+lRUFq6gEIK/1fuTgCnxD5IK9sqd8mY5TzVj6ht2jxq6e+5WiBjH1RxmAnmdSVtzyJX9Xof0L4+UfoUbPp846up3sT7L/uRe96b7+/aW9tz5NfZE/ZCuuwiK2zLfaF7bA+E+wXO2Vn7DxoBM+DXhBeSoPGfGaZXang/W8oR8EZ</latexit> F = (SSR0 SSR1)/r SSR1/(T np 1) 制約(β=0)を⼊れた ときの残差平⽅和 制約(β=0)を⼊れな いときの残差平⽅和 n変量VAR(p)の場合のF値 制約の数 (今は1)
27. インパルス応答関数: IRF (Non-orthgonalized) Impulse Response Function (ラグk の関数) (Orthgonalized) Impulse Response Function 27/48 <latexit sha1_base64="6MmCyJ+rPMbb7A6dn2x+aXWhFDU=">AAACmHicfVBdaxNBFJ2sXzV+NNU3fRkMQtUQdsWvB4WgRe2DWItpC9kQ7k7uptPMzgwzd4th2T/kr/G1/hpn04i2FS8sHM45d+6ek1klPcXxSSu6dPnK1Wtr19s3bt66vd7ZuLPnTekEDoVRxh1k4FFJjUOSpPDAOoQiU7ifzd81+v4xOi+N/koLi+MCZlrmUgAFatLZ2t59P6nkUb05f8Tf8HSaOxBVasGRBMUXQevRk3ld/+HSY3BovVRhvzrqUV1POt24Hy+HXwTJCnTZanYmG63X6dSIskBNQoH3oyS2NK6aC0Jh3U5LjxbEHGY4ClBDgX5cLePW/GFgpjw3Lnya+JL9e6OCwvtFkQVnAXToz2sN+S9tVFL+alxJbUtCLU4P5aXiZHjTHZ9Kh4LUIgAQToZ/5eIQQmEUGj5zpXmbjFE+RNnCENHhp0B9tuiAjHtcpeBmhdR1iDxLew36nxG+/TYG1G6HvpPz7V4Ee0/7yYv+8y/PuoO3q+bX2H32gG2yhL1kA/aR7bAhE+w7+8FO2M/oXjSIPkTbp9aotdq5y85MtPsLyKfO4g==</latexit> IRFij(k) = @yi,t+k @"j,t j 番⽬の変数における時刻t のショックがi 番⽬の変数 における時刻t +k をどれだけ変化させるか でも、変数i とj が相関してたらj における ショックって測りにくくないですか？ <latexit sha1_base64="OvK2sjNuYNRMiMFp3hgfRxbe4o4=">AAADEHicfVFda9RAFJ3Erxq/tvroy9DFIlqXJGjrg0JRH/oituK2hc0SbmbvZodOMmFmUlxC/oS/xgdBfPUf+G+cbLLY7RYvDBzOOffM3DtJIbg2vv/Hca9dv3Hz1sZt787de/cf9DYfHmtZKoZDJoVUpwloFDzHoeFG4GmhELJE4Ely9r7RT85RaS7zL2Ze4DiDNOdTzsBYKu792J7HVbBjavqWRiCKGcT+8yhBAw0d1K36Iqj/kWFDhh15DgoLzYXNamOiiHrbreFCZrBsD6/KDFczU8gy27KWvVRCWrb5ca/vD/xF0XUQdKBPujqMN5030USyMsPcMAFajwK/MOMKlOFMYO1FpcYC2BmkOLIwhwz1uFqsuaZPLDOhU6nsyQ1dsBc7Ksi0nmeJdWZgZvqy1pBXaaPSTF+PK54XpcGctRdNS0GNpM2f0QlXyIyYWwBMcftWymaggBn7syu3NNlGSqHtKB/Qjqjwo6U+FajASPWsikClGc9rO3Ia7TTof0b4ujRa5Hl238Hl7a6D43AQ7A5eHb3s77/rNr9BHpMt8pQEZI/skwNySIaEOVvOgXPkfHa/ud/dn+6v1uo6Xc8jslLu779MpvYB</latexit> y1,t = ↵0 + 1,1y1,t 1 + 1,2y2,t 1 + "1,t y2,t = ↵1 + 2,1y1,t 1 + 2,2y2,t 1 + 1"1,t + 2u2,t 誤差項をy1 と相関する部分と そうでないu とに分解 = 直交化 <latexit sha1_base64="RERpRETg5ogHdyoQIJT4pMEcnaw=">AAACjnicfVBdaxNBFJ2sX238SvXRl8EgVA1hV7SKWCytSH0Qq5i2kA3L3cnddJrZmWXmbmlY9r/01/RVX/tvnE0j2la8MHA4537MOWmhpKMwPGsF167fuHlrabl9+87de/c7Kw92nSmtwIEwytj9FBwqqXFAkhTuFxYhTxXupdOtRt87Quuk0d9pVuAoh4mWmRRAnko6bz99+5hU8rBenT7l6zweZxZEFRdgSYLiM6/16Pm0rv9wZVId9qiuk0437Ifz4ldBtABdtqidZKX1Lh4bUeaoSShwbhiFBY2qZq9QWLfj0mEBYgoTHHqoIUc3quYma/7EM2OeGeufJj5n/56oIHdulqe+Mwc6cJe1hvyXNiwpezOqpC5KQi3OD2Wl4mR4kxgfS4uC1MwDEFb6v3JxAD4m8rleuNLsJmOU81Y+oLdo8bOnvhRogYx9VsVgJ7nUtbc8iXsN+l8jHP9u9Kjd9nlHl9O9CnZf9KO1/quvL7sbm4vkl9gj9pitsoi9Zhtsm+2wARPshJ2yH+xn0AnWgvXg/Xlr0FrMPGQXKtj+BVSTylQ=</latexit> IRFij(k) = @yi,t+k @uj,t 無相関な部分へのショックがどれだけ他に影響するか？を定量化
28. 次元削減 多次元のデータを低い次元に落とすこと 有名なところにPCA、t-SNE、UMAPなど 多変量時系列データも低い次元に落とすが可能 時系列的な前後関係(や時間依存関係)を保持しながら、低次元に落とすのはなかなか困難 普通にPCAとかやるとだめ! (離散Fourier変換により漸近的に独⽴な周波数領域のデータに変換してからとかはOK) 有名なところにDynamic factor analysis、Functional PCAなど 28/48 PCAの復習 Kernel法的考え⽅ 逆に⾼次元に射影すると⾒通しが良くなる場合があるよね
29. 関数データ解析 ちょっと研究中 要は、点データじゃなくって、関数をデータとしてみる 均⼀で計測数が多い場合はぶっちゃけ古典的な経時データ解析⼿法使えばOK 観測時点、観測地点数などが個体ごとに異なっていたりすると効果を発揮 29/48 元データ 関数データ化 基底関数の係数に関してPCA <latexit sha1_base64="2Ll+s0wMPKhMJA0xdkopwWbWw7o=">AAACi3icfZBdSxtBFIYn29aPtGq0l94MDQVbS9j1sxQL0g/ojahgVMiG5ezkJA7OziwzZ8Ww5Kf013irP8B/09kkharFAwMP7zlnZt43zZV0FIb3teDFy1czs3Pz9ddvFhaXGssrp84UVmBbGGXseQoOldTYJkkKz3OLkKUKz9LL71X/7Aqtk0af0DDHbgYDLftSAHkpaewOE8m/8jhFgiTk6zyWmpITfp3INfowkSNPPVqPr8Bi7qTye5Q0mmErHBd/CtEUmmxaR8lybS/uGVFkqEkocK4ThTl1S7AkhcJRPS4c5iAuYYAdjxoydN1y7HDE33ulx/vG+qOJj9V/N0rInBtmqZ/MgC7c414l/q/XKaj/uVtKnReEWkwe6heKk+FVXLwnLQpSQw8grPR/5eICLAjyoT54pbqbjFHOW/mB3qLFAy8d5miBjP1YxmAHmdQjb3kQf6rouUG4/jvoqV73eUeP030KpxutaKe1fbzV3P82TX6OrbJ3bI1FbJfts1/siLWZYL/ZDbtld8FCsBl8CfYmo0FtuvOWPajg5x84zseA</latexit> yi = 0 + Z T xi(t) 1(t)dt + "t Ramsey & Silverman (2005) y がスカラーで x が関数 <latexit sha1_base64="c9EVRVpW/adeqqWrbRWJXPCit+k=">AAAClHicfVBdaxQxFM2OX3X86FbBF1+Ci6WKLDPFr4cWFqvgi1rBbQubdbiTvbsbmkyG5I52Gfbv+Gt8VfDfmJmuYFvxQuDk3HNzc05eauUpSX51okuXr1y9tnY9vnHz1u317sadA28rJ3EorbbuKAePWhU4JEUaj0qHYHKNh/nxXtM//ILOK1t8okWJYwOzQk2VBApU1h1snmRqix7t8q+fBdkyU6KcqywNlBDxpsiRoL3tcjEDY6BVtZrtVpN1e0k/aYtfBOkK9Niq9rONzo6YWFkZLEhq8H6UJiWNa3CkpMZlLCqPJchjmOEowAIM+nHdWl3yh4GZ8Kl14RTEW/bviRqM9wuTB6UBmvvzvYb8V29U0fTluFZFWREW8nTRtNKcLG9y4xPlUJJeBADSqfBXLufgQFJI98yW5m2yVvtg5TUGiw7fBepDiQ7Iuse1ADczqlgGyzPxpEH/E8LJH2FAcRzyTs+nexEcbPfT5/1nH5/2Bq9Wya+x++wB22Ipe8EG7C3bZ0Mm2Tf2nf1gP6N70U60F705lUad1cxddqai978BGGzKWA==</latexit> xi(t) = w> i 1(t) 1(t) = > 2(t) x とβ の基底関数展開