最近ベイズモデリングについて勉強中です。 しかし、数式ばかり追いかけてしまって分かったようなわかんないような、、、という状況です。
ということで、今回は実際にモデリングして、これまで数式で勉強してきた内容の復習をしようと思います。
参考にさせていただいたのはこちら。
実践Ｄａｔａ Ｓｃｉｅｎｃｅシリーズ ＲとＳｔａｎではじめる ベイズ統計モデリングによるデータ分析入門 (ＫＳ情報科学専門書)
- 作者:馬場真哉
- 出版社/メーカー: 講談社
- 発売日: 2019/08/23
- メディア: Kindle版
非常にわかりやすく、入門書として最適と評判ですね。 今回はこの本の中で、一般化線形モデル（GLM）について、例題を通して勉強していきたいと思います。
（統計）モデル is 何？
この手の話をすると「モデル」という単語がよく出てきます。曖昧な言葉である「モデル」という言葉について先に確認します。
モデルはこのように説明されています。
確率モデルをデータに適合させたものを、本書では統計モデルと呼んでいます．統計モデルは単にモデルと呼称することもあります． 実践Data Scienceシリーズ RとStanではじめる ベイズ統計モデリングによるデータ分析入門
与えられたデータに基づいてその周辺のしくみを数式を使って表現したものを、モデルと呼んでいるようですね。
一般化線形モデル
上で紹介した”モデル”と、一般化線形モデルについては若干ニュアンスが異なり、フレームワークの意味合いが強いです。 その説明から、モデルがどのようなものがあるのか確認します。
定義・説明
wikipediaによれば、一般化線形モデルはこのような説明になっています。
一般化線形モデル (いっぱんかせんけいモデル、英: Generalized linear model、GLM)は、残差を任意の分布とした線形モデル。似たものとして一般線形モデルがあるが、これは残差が多変量正規分布に従うモデル。一般化線形モデルには線形回帰、ポアソン回帰、ロジスティック回帰などが含まれる。 wikipedia
上で紹介した参考書籍の説明を借りれば、一般化線形モデルは一種の「モデルの型」あるいは「フレームワーク」という位置づけになります。 一般化線形モデルの中によく使われるモデルがいくつかあり、それらを総称して"一般化線形モデル"と言われるそうです。
一般化線形モデルで、モデルを決める要素として
- 確率分布
- 線形予測子
- リンク関数
の３つが挙げられ、観測データに合うようにこれら３つを変更しながら一般化線形モデルの中で適切なモデルを選択していくことが、一般化線形モデルを使用したモデリングを行うことになります。
最近勉強している理論的な話とかが出てきていないのは、MCMCが大いに関係しているようです。
数学的に解く場合には、確率変数を含んだ積分を解く必要がありますが、一般にそれを解くことは簡単でないことが多いとされています。 それをMCMCによって、求める分布の乱数を得ることができるため、擬似的にその積分計算というステップをスキップしています。 そのため、我々としては難しい積分計算を考えず、純粋にモデルの構造を予想しさえすれば、MCMCを含むライブラリで近似計算ができるように成っているんですね。
一般化線形モデルの仲間
参考にしている書籍に紹介されていたモデルをまとめるとこんな感じかと思います。
- 単回帰モデル
- 分散回帰モデル
- 正規線形モデル
- ポアソン回帰モデル
- ロジスティック回帰モデル
- これだけ今回やってないです。機会があったらそのうちやります。
モデリングしてみる
モデリングはこのように定義されています。
データからモデルを推定し、推定されたモデルを考察することで、現象の解釈や予測を行います．モデルを構築することをモデリングと呼びます．
ということで、得られているデータをもとにデータにまつわる現象の仕組み（モデル）を予想して、それについて考察を行います。
今回はモデリングするにあたって、機械学習の一般的な進め方のうち一部について見ていきます。
- どんな数式モデルに従うかを予想する
- データを使用してモデルのパラメータを設定する
- 予想したモデルについて考察する
これまでなんとなく数式を追いかけてきたベイズモデリングの理屈について、今度は実際にコードを確認していきます。
環境構築にはkaggleのRのDocker Imageを使用しました。 環境構築についてはこちらをご参考ください。
単回帰モデル
単回帰モデル（simple linear regression model）は量的データの説明変数が１つだけ使用され、正規分布に従う量的データが応答変数であるモデルです。
問題設定
Rの標準でついてくるデータセットのcarsを使ってみます。
いま、自動車の速度と停止距離についてのデータが与えられたとします。
この時、速度と停止距離の関係を表現したいとします。
考え方
上のプロットから、なんとなくですが速度と停止距離には単調増加の傾向がありそうに見えます。 そこでこの問題については単回帰モデルを適用してみます。
単回帰モデルでは、構造は下記の様になります。
が線形予測子、
が平均値、標準偏差が
を表しています。
この形式に上で与えられたデータを当てはめることが単回帰モデルを構築することに当たります。
この例では、速度が説明変数、停止距離が目的変数に当たります。
線形予測子の形式から、一次関数のグラフになっていることがわかります。 停止距離はこの線形予測子を平均値として使用した正規分布に従う、と予想していることになります。
実際にこのモデルが妥当かどうかはMCMCサンプリングをするとすぐに確認することができます。 実装の詳細は最後にgithubに譲るとして、このモデルを使って最終的に生成したものはこのようになります。
実線が回帰直線、網掛けの部分がこのモデルに基づく予測が取りうる範囲(95%の予測区間)になっています。 もともとのサンプルが、ほとんど予測区間に入っているため、ある程度妥当なモデルになっていることがわかりました。
余談ではありますが、一般的に車の速度と停止距離の関係は速度の二次関数になることが知られているので、実際には二次関数としてモデリングしたほうが妥当なのかもしれませんね。
分散分析モデル
分散分析モデルは、説明変数に質的データを取り、確率分布に正規分布を仮定したモデルです。
問題設定
いま、ミツバチの駆除のために8パターンのスプレーを使用して、その効果を測定した。
この時、スプレーの種類とミツバチの減少の関係を表現したいとします。
考え方
今度は、単回帰モデルの例とは異なり、横軸が不連続な質的データになっています。 そのため、今回は分散分析モデル(ANOVA)に適用してみます。
分散分析モデルは構造は下記の様になります。
この例では、スプレーの種類が説明変数に質的データとし、どれか一つだけ1でその他は0としています。
モデルに適用した結果がこちらになります。
今回は質的データなので連続性はありませんね。 もともとのデータが比較的収まっているので、まあまあ妥当なモデリングなんでしょうか。
正規線形モデル
正規線形モデル(GLM)では量的・質的データを問わず、複数の説明変数を持ち応答変数が正規分布に従うモデルです。 正規線形モデルの条件としては、
- 量的・質的データ問わず、複数の説明変数を線形予測子に用いることができる
- 恒等関数がリンク関数
- 正規分布を確率分布に用いている
のようになっています。 そのため、単回帰モデルや分散分析モデルも正規線形モデルに含まれることになります。
問題設定
こちらはirisデータセットを使ってみます。
各地のアヤメをサンプリングしたデータが手元にあるとします。
このときの花弁の長さと因子（種類、ガクの長さ）の関係性を表現したい状況を考えます。
考え方
上のプロットから、花弁の長さは種類とガクの大きさに相関がありそうなのは確認できます。
ということで、次のような正規線形モデルに適用させてみたいと思います。
モデルに適用させた結果がコチラになります。
一応、種類・ガクの大きさごとに傾向を掴んでそうな回帰直線が引けてます。 信頼区間の外側に大量にデータが出てしまっているので、あまり良いモデリングではありませんが。。。
ただ、普通、植物の種類って大きさとかから決まってたりするんで、そもそも種類と花弁の長さの間に相関がありそうですけどね。 完全に独立変数を取り扱えばもっときれいにモデリングできるのかもしれませんね。
ポアソン回帰モデル
ポアソン回帰モデルはこれまでの正規線形モデルと異なり、確率分布が正規分布ではなく、離散型のデータを扱うポアソン分布に従い、リンク関数に対数関数が使用されているものをポアソン回帰モデルと呼んでいます。
問題設定
今度は、生息する動物の種類のデータセットgalaを使用します。 下のようなデータが得られていたとします。
ここで、縦軸は種類数、横軸は島の面積（の対数）を表しています。 このようなときの、面積と生物の種類の関係を表現したいとします。
考え方
上の例にポアソン回帰を適合させる際の構造は下記の様になります。
確率分布にポアソン分布、リンク関数が指数関数になっているところが特徴です。
モデルに適用させた結果がコチラになります。
一定の傾向はつかめている感じはしますが、信用区間の外にデータが出てしまっています。 そのため、もうちょっとモデルを見直す必要がありますね。
参考
使ったコード
感想
今回は例題を使って、一般化線形モデルに適用させてみました。 本来はきちんと問題設定して解くのが良いと思いますが、今回はとりあえずモデリングをやってみようということで結構適当にやってみました。 そのため、モデルが妥当かどうか、またモデルが妥当でなかった際にどう修正していくのかは今後勉強しないといけないですね。
普段数式を追いかけるばかりでなにがなんだかよくわからない部分も多かったですが、具体例を使って実際に手を動かすとできてしまうのが不思議です。 今後は、理論もモデリングも並行して勉強していく必要がありそうです。
まだまだ複雑なモデリングができるわけでもなければ、データの相関関係などを見抜いたりもできませんが、そのへんはやっていくうちにできるようになるかと。 ちょっとずつ勉強していきたいものですね。
