らくがき入門

機械学習を始めとしたコンピュータサイエンスを主に扱っています。

PandasのDataFrameから特定の値を持つ行を削除する

言われてみたら簡単なんだけど、意外に思いつかなかった。

DataFrameの特定の行に含まれている値を指定して、それ以外を抽出するイメージ。

In [1]: import pandas as pd

In [2]: df = pd.DataFrame([[1, 2], [2, 3], [3, 4]], columns=["a", "b"])

In [3]: df
Out[3]: 
   a  b
0  1  2
1  2  3
2  3  4

In [4]: df[df.a != 2]
Out[4]: 
   a  b
0  1  2
2  3  4

Pythonによるデータ分析入門 第2版 ―NumPy、pandasを使ったデータ処理

Pythonによるデータ分析入門 第2版 ―NumPy、pandasを使ったデータ処理

pythonで地図上の2時点間の位置関係を求める;

2時点間の緯度、経度、高度が与えられたときに簡易的に2時点間のローカルな位置関係を求めます。 緯度・経度・高度といった情報ではローカルな2時点間の関係性が分かりづらく扱いづらいので、変換します。 具体的には、2時点間の距離とある時点からもう一時点を見たときの方位角、仰角を求めることを目標とします。

以下、国土地理院のサイトと独立行政法人電子航法研究所の資料を参考に作成しています。

https://www.enri.go.jp/~fks442/K_MUSEN/1st/1st060428rev2.pdf

地球上の2時点間の距離

高校などで扱った平面世界での2時点間の距離は、簡単に求めることができたと思いますが、地球上の2時点間ではそれほど簡単には求めることができません。 なぜなら、地球が球体であり、局面上にある2時点間の距離を求めることになるからです。

ジオイドとは

地球は、時点による遠心力の影響で、赤道方向が少し膨らんだ楕円体のような形をしています。 測地学では、世界の界面の平均位置に最も近い「重力の等ポテンシャル面」を「ジオイド」と定めて、これを地球の形状ということにしています。標高はジオイドを基準として定められており、標高はジオイドから測った高さになります。

GPSを用いて標高を求めるには

現在広く活用されているGPSでは、緯度・経度・高度といった幾何学的な位置を求めることができますが、重力を考慮した設計になっていないため標高を直接求めることができません。GPSを用いて標高を求めるためには、ジオイド高がj必要になります。

標高 = 楕円体高 ー ジオイド

ジオイド高は国土地理院が提供している計算サイトで計算することができます。

WGS84

WGS84とは、GPSの基準座標系で、直交座標系です。これはECEFと呼ばれ、地球の重心を原点として地球の自転軸の北極方向をZ軸、Z軸に垂直にグリニッジ子午線の方向をX軸として、これらの軸と直行するように右手系でY軸とした直交座標系です。

ECEFから東をX, 北をY、上空をZとした地平直交座標(ENU)に変換します。 ECEFからENUへの変換は、行列による回転と原点移動で実現できます。詳しくは上記で示した電子航法研究所の資料をご確認ください。

githubにサンプルコードを提示してあります。coordinate_test.pyの9, 10行目で取り扱いたい2時点の緯度・経度・ジオイド高のサンプルデータを与えています。結果表示では、2時点間の距離およびある時点(コードではRx)からみたもう1時点(コードではSat)の方位角、仰角を求めることができます。

atom上でterminalを使う

atomプラグインであるPlatformio Ide Terminalを起動して、paneとの間を移動する方法を説明します。

atomのインストール

atomは以下のサイトから自分の環境に合わせてインストールしてください。

環境

ubuntu18.04

Platformio Ide Terminal

atomでターミナルを起動するプラグインです。atom上でterminalを起動するのはTerminal-Plusが一般的なようですが、私の環境ではうまく動作しなかったのでPlatformio Ide Terminalにしました。

apm install platformio-ide-terminal

でインストールできます。
ショートカットキーはctrl + ` で起動することができます。

paneとterminalとの間の移動はctrl + alt + fで移動することができます。

atom上でterminalもbashrcの内容が反映されるので使い心地は変わらないような気がします。

参考にしたサイト

全面的にこちらのサイトを参考にしました。色々なプラグインを比較されていてわかりやすかったです。

jump-up.info

atom関連だとこの本がまとまっていた。

テキストエディタAtom入門 (OIAX BOOKS)

テキストエディタAtom入門 (OIAX BOOKS)

時系列分析における処理フロー(定常過程と単位根過程)

前回の記事では、「時系列データの定常性を確認する」というタイトルで記事を書きました。与えられた時系列データが定常過程であれば、ARMAモデルなどの基礎的なモデルに落とし込むことができます。

今回はよく用いられる単位根検定である拡張 Dickey-Fullar(augmented Dickey-Fullar test;ADF test)を取り上げます。

使用言語はPythonです。 今回も沖本本を参照しています。

経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

単位根過程とは

前回の記事で、定常過程に関しては取り上げているのでそちらを参照してください。

leisurelab.hatenablog.com

単位根過程とは、原系列 y_tが非定常過程であり、差分系列 \Delta y_t = y_t - y_{t-1}が定常過程であるとき、過程は単位根過程であるといわれる。

単位根過程は、別名がいくつかあり、差分系列が定常となるため差分定常過程と呼ばれることがあります。

さて今回Pythonのコード中の検証データとして使用しているランダムウォークに関しても記述しておきます。 なぜ今回ランダムウォークを検証データとして使用しているかと言うと、ランダムウォークが単位根過程の代表選手だからです。 ランダムウォークの定義は、

過程 y_t
 \displaystyle
y_t = \delta + y_{t-1} + \varepsilon_t
,  \varepsilon_t \sim \text{iid}(0, \sigma ^2)
と表現されるとき、 y_tランダムウォークと呼ばれる。定数項 \deltaはドリフト率と呼ばれる。

上記の定義より、繰り返し代入すると、
 \displaystyle
y_t = \delta t + \nu_t
と表現できます。ただし、 \nu_t = \varepsilon_1 + \varepsilon_2 + \cdots + \varepsilon_tという撹乱項の線形和で表現され、これは確率的トレンドと呼ばれます。 ドリフト率 \deltaは線形トレンドの傾きを表しているので、単位根過程は線形トレンドと確率的トレンドを併せ持つ過程と言えます。

定常過程と単位根過程との性質の違い

定常過程と単位根過程の違いを以下の表にまとめました。

性質 定常過程 単位根過程
トレンド もたない 線形トレンドを持つ
長期予測 過去の観測値が減衰していき過程の期待値に近づく 過去の観測値の影響が消えない
予測誤差(MSE) 予測期間が長くなるにつれて過程の分散(有限値)に近づく 長期的に線形的に増大していく
インパルス応答関数 ショックは一時的な影響力を持つ ショックは恒久的な影響力を持つ

表のトレンドに関して詳しく説明します。
定常過程(弱定常性に従う過程)は前回の記事で説明したように、期待値と自己共分散が時間によらず一定であるという仮定を置いています。 この仮定より、明らかに定常過程はトレンドを持ちません(時間によらず一定であるから)。 定常過程は平均回帰性を持ち、過程が長期的には必ず平均の方向に戻っていくことを意味します。
ただし、定常過程であっても線形トレンドを記述できるモデルとして、トレンド定常過程と言われるものが存在します。
トレンド定常過程は、文字通り定常過程にトレンド項を加えることで線形トレンドを記述できるようにしたモデルを指します。
それでは、トレンドを持つトレンド定常過程と単位根過程は線形トレンドを記述できるという点で同じなわけですが、前述したように単位根過程は確率的トレンドも持つ点で異なります。 単位根過程は不確実性を線形的に増大させていると言えます。

インパルス応答関数に関してですが、VARモデルに関して取り上げる際に書こうと思っているので一旦保留とします。 簡単に概要だけ説明すると、ある時点の影響が将来どの程度影響をもつかを定量的に表現できるツールです。

単位根検定

単位根検定として有名なADF検定について説明します。 ADF検定は真の過程をAR(p)モデルと仮定し、仮定が単位根AR(p)過程であるという帰無仮説を、過程が定常AR(p)過程であるという対立仮説に対して検定するものです。 これが厳密な定義ですが、帰無仮説「単位根がある」、対立仮説「単位根がない(=定常である)」と言ってしまってもいいと思います。

この検定を用いる場合、ある有意水準に対して帰無仮説を棄却できれば、対立仮説が支持されて与えられた時系列データが定常過程であると言って良さそうです。 ただし、帰無仮説を棄却できなかった場合は、帰無仮説を積極的に支持したことにはならないことに注意してください。これは統計学的非対称性によるのですが長くなりそうなのでまた別の機会に取り上げます。

したがって、与えられた時系列データが単位根過程であるかは、冒頭で取り上げた単位根過程の定義である過程の差分系列が定常課程であることを用います。 したがって、手順としては以下のようになります。

  1. 与えられたデータに対してADF検定→
    帰無仮説が棄却されれば、与えられたデータは定常過程
  2. 帰無仮説を棄却できなければ差分系列を求め、その系列に対してADF検定→
    帰無仮説が棄却できれば、与えられたデータは単位根過程
  3. 帰無仮説が棄却できなければ、与えられたデータは単位根過程ではない非定常過程

このようにしてデータを分類することになります。

この手順をPythonで実装した例を示します。

gist7fb5ceb03fae95dc3f0774debf880340

経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

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

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

時系列データの定常性を考える

沖本本を基に、定常性について説明します。 定常性は時系列モデルを考える際に最も重要な概念の1つです。

定常性とは

定常性は、大雑把に言うと時系列モデルの性質です。 定常性は、同時分布だったり基本統計量に関して時間的に不変 であるという性質を表現したものです。 想像していただきたいのですが、 仮に時系列データのどの時点の期待値も、時間によらず一定であると考えられるなら かなり強力な仮定になると思いませんか。 定常性という強い仮定の下で、時系列データを分析する際に有名なARMA過程などの基礎的なモデルに落とし込むことができます。 定常性には大きく分けて2つの種類があり、

  • 弱定常性
  • 強定常性

がありますが、経済・ファイナンス分野では単に定常性というと弱定常性を指すことが多い みたいです (これは強定常性が文字通り強い仮定を置くものであり、 現実のデータを用いた場合にその仮定が正しいのかを検証することが難しいからではないかと思います。) なので、今回は弱定常性に関して詳しく見てみようと思います。

自己共分散について

弱定常性を確認する前に、導入として自己共分散について説明しておきます。

統計的アプローチを取るデータ分析の場合、まずやることとしては基本統計量を用いてデータの要約を行い、 与えられたデータがどのようなデータなのかを調査します。 時系列a分析ではない一般的なデータ分析の場合、基本統計量は

などがあります。

時系列分析においては、更に自己共分散という統計量が追加されます。 自己共分散は、同じ時系列データのある時点と別の時点との共分散 です。 なぜ、時系列分析において自己共分散を確認する必要があるのか。 これは、時系列データが過去の自分自身のデータに影響される可能性があるためです。 簡単な例として株価を考えると、ある日の株価が上がっていてその株価の値が過大評価であると市場の参加者が考えた場合は、 翌日は確実に下がるでしょう。 このように時系列データにおいて、過去の時点のデータに影響された動向変化は容易に想像できます。

自己共分散の解釈は

  • 自己共分散が正であれば、期待値を基準として同じ方向に動く傾向
  • 自己共分散が負であれば、期待値を基準として異なる方向に動く傾向

となります。

一般的にk次の自己共分散は

 \displaystyle
    \gamma_{kt} = \text{Cov}(y_t, y_{k-t}) = E[(y_t - \mu_t)(y_{t-k} - \mu_{t-k}) ]

で定義されます。

弱定常性について

やっと弱定常性の定義の話ができますね。 弱定常性の定義は

 \displaystyle E(y_t) = \mu

 \displaystyle
\text{Cov}(y_t, y_{t-k}) = E[(y_t - \mu)(y_{t-k} - \mu) ]  \displaystyle = \gamma_k

となります。

弱定常性は時系列データの期待値と自己共分散が時間によらず一定であることを要求します。 その要求を満たせば弱定常性の仮定の下で分析を進めていくことができます。 時系列モデルを構築する場合、自分がどのような仮定の下で分析しているのかを意識することはとても大事だと思います。 弱定常性の定義から自己共分散は時間差 kに依存することには注意が必要です。

与えられたデータが定常過程かどうかは次回取り上げるADF(Augmented Dickey-Fuller)検定により判別することができます。 ADF検定の概要は、帰無仮説「単位根が存在する」、対立仮説は「単位根がない(=定常である)」である検定です。 この検定において帰無仮説が棄却された場合は、対立仮説である「単位根がない(=定常である)」が支持されたと結論づけます。

まとめ

データが定常であるかどうかによって分析するプロセスが変わるので最初に与えられたデータが定常過程かどうかの確認はするべきです。 書いていて疑問に思ったのは実データに対して強定常性を仮定できるのはどういったデータかということでした。 理論的なもので、応用例が存在するかどうかは気になりますね。

経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

時系列分析における系列変換

沖本本のアウトプットを共有目的と自分用のメモ目的で書きます。

時系列分析の目的

時系列分析の目的は、複雑な観測データが持つ多様な特徴のうち 分析者にとって重要な特徴のみを簡単に表現できるモデルを構築することです。 このように作成したモデルを基に、目的に応じた分析を行います。

時系列分析における変換の目的

何も加工されていない観測データそのものの時系列データは原系列 と呼ばれます。 原系列は観測データそのものですが、 自然現象を観測データにする際にもある程度の情報の取捨選択が行われていることは一応念頭に置いておくべきです。

基本的に時系列分析の目的は、この原系列の性質を明らかにして予測などを行うことが多いですが、 原系列に対して時系列分析をするよりも何らかの変換を施して分析を進めるほうがいいことが多いです。 なので、分析におけるざっくりとした方針は以下のようになります。

  1. 原系列に対して何らかの変換を施す
  2. 原系列を変換した系列に対して分析を行う
  3. 変換を逆変換して原系列に戻すことで、原系列の性質を明らかにする

変換の種類

  • 対数変換

原系列に対して対数の変換を施します。 たとえば、経済・ファイナンスデータなどは、値が大きくなるにつれてばらつきが大きくなるデータなどがあり、 時系列分析において重要な性質である定常性の仮定を満たさない場合があります。 そのような場合は対数変換を施した対数系列 にすることで定常過程に従う時系列データにする場合があります。

  • 差分変換

原系列に対して1時点離れたデータの差を新たな系列に変換して分析することもあります。 これを差分系列と呼びます。 時系列データには、定常性を持つ定常過程に従う時系列データの他に、単位根過程に従う時系列データなどが存在します。 今回は、単位根過程の詳しい説明は省略しますが、単位根過程の差分系列は定常過程 になるということだけは抑えておいても良さそうです。 基本的に時系列分析において時系列データが定常過程に従うと仮定できると解析しやすくなるので、 その点でも単位根過程の差分系列が定常過程に従うという事実は重要です。

  • 対数差分変換

上記の対数系列の差分系列である対数差分系列に変換します。 この変換を施す動機は変化率に興味がある場合です。 離散データの場合は、対数差分系列に100を掛けて変化率とする ことができます。 対数差分系列がなぜ変化率を表現しているのかに関しては1次のテイラー展開の近似から説明できます。

  • 季節調整

時系列データには季節性を含んだデータも多いです。たとえば株価などがわかりやすい例として考えられます。 たとえば、一般的に株価は年の前半は上昇トレンドですが、年の後半は減少トレンドであることが知られています。 時系列分析において、季節変動では説明できない時系列データの動きに興味のあることが多いように思います。 季節調整自体は難しい問題で、僕自身も勉強不足なので扱いきれていませんが、Pythonのstatsmodelsというライブラリの中で季節調整のできる関数が実装されています。

from statsmodels.tsa.seasonal import seasonal_decompose
seasonal_decompose(x, freq)

xにはarray-likeな時系列データを渡して、freqで季節性の周期を渡します。たとえば、日次のデータだとfreq=7 にすることで一週間の季節性を排除します。説明のために簡易的になってしまいましたが、その他も調整するパラメータが存在するので 公式ドキュメントを確認してみてください。

まとめ

今回は時系列データの変換に関してまとめました。 機械学習におけるデータのスケール変換と同様に上記で示した変換はモデルの精度に影響が及ぶため、 重要だと思います。 統計モデリングにおいてなぜこの変換をするのか、今現在どういう仮定の下でモデリングをしているのかは 常に人に説明できるようにしておくことが重要だと日々感じています。 まだまだ勉強不足なのでアウトプットと並行しながらインプットも進めていきたいですね。

経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

経済・ファイナンスデータの計量時系列分析 (統計ライブラリー)

python+SeleniumでWebページのデータを自動取得する

SeleniumWebブラウザの自動操作ライブラリです。Seleniumを用いることで人間がWebブラウザ上で行う作業を自動化することができます。

Seleniumの一般的な用途としては、単純作業の自動化の他にWebアプリケーションの自動化などが挙げられます。

今回はSeleniumでWebページにアクセスしてcsvファイルをダウンロードするという作業を自動化します。

使用するもの

今回の主役です。様々な言語に対応して提供されていますが、今回はpythonを使用します。

  • ChromeDriver

ブラウザ(Chrome, Firefox, Safari etc...)のdriverを指定することでSeleniumを各ブラウザで実行することができます。

それぞれのブラウザではヘッドレスモード(ブラウザを実際にGUI上に立ち上げずに実行するモード)が可能です。今回はデータを取得するだけなのでヘッドレスモードで実行します(オプションで--headlessで指定するだけ)。

2018年6月2日ヘッドレスブラウザの代名詞だったPhantomJSの開発が終了アーカイブ化されたので、上記のブラウザを使用することをおすすめします。

なお、今回高頻度で利用することはないのでCircleCIは使用しません。

Selenium操作の基本部分

from selenium.webdriver import Chrome, ChromeOptions
options = ChromeOptions()
options.add_argument('--headless')
driver = Chrome('chromedriverのパス', options=options)
driver.get('操作したいブラウザのURL')

この部分がSeleniumによるブラウザ自動操作の基本部分になります。 オプションで--headlessを指定することでヘッドレスモードを選択できます。

driverでChromeを指定したあとは、driver.get(URL)でブラウザを開くことができます。

スクレイピング

Seleniumを用いたスクレイピングも一般的なスクレイピング同様HTMLの要素にアクセスして実行していきます。

# idで取得
driver.find_element_by_id('ID')
# classで取得
driver.find_element_by_class_name('CLASS_NAME')
# nameで取得
driver.find_element_by_name('NAME')
# xpathで取得
driver.find_element_by_xpath('XPATH')

ChromeのDeveloper toolsを使えばWebページを構成するHTMLがわかるので地道にみながら操作していく感じです。 elementの部分をelementsにすると該当する要素が複数取得できます。 elementでアクセスした際、要素が複数ある場合は一番はじめの要素にアクセスできます。

データの自動取得

今回は日本経済新聞の経済指標ダッシュボードのcsvデータを取得します。

https://vdata.nikkei.com/economicdashboard/macro/

経済指標ダッシュボードは日本経済の動向を把握するのに役に立つ経済指標が一覧できるサイトです。 ブラウザ上で表示されているグラフをクリックするとcsvがダウンロードできるボタンが表示されます。

基本方針としては、csvデータのURLが含まれる要素まで行ってgetAttribute()で名前付けされた属性の値を取得しにいくという方針です。

そしてそこで取得したURLに対して

import urllib
urllib.request.urlretrieve('指定のURL', '保存先のファイル名')

でURLからデータをダウンロードします。

githubにコードを上げてあります。

github.com

まとめ

Seleniumでブラウザ操作が自動化できるとかなり効率化につながりますね。まだまだ勉強できていないので色々試しながらやっていきたいです。

以下は今回勉強するに当たって参考にしたサイトです。大変勉強になりました。

tech.mercari.com

HTML付近の勉強は下の本をチラ見しながら勉強しました。図が豊富でわかりやすいです。

作りながら学ぶ HTML/CSSデザインの教科書

作りながら学ぶ HTML/CSSデザインの教科書