[計量経済学] [機械学習] Generalized Random Forest (GRF) について
概要
Athey, Tibshirani, & Wager (2016, Generalized Random Forests) で提案されている Generalized Random Forest (GRF) について解説してみる.
[1610.01271] Generalized Random Forests 2019/7/4 追記: この論文は Annals of Statistics にアクセプトされたようだ.
projecteuclid.org
計量経済学と機械学習の両方の文脈を追う必要が出てくるので, 機械学習を学んできた人, (計量) 経済学を学んできた人, それぞれに対して伝わりやすいように説明を試みる.
先日の Tokyo R #71 で以下のような発表があった.
ここで公開されているソースを確認していると, grf
というパッケージが読み込まれていた. しかし, ソースコードではこのパッケージの関数は使われていない. なんのパッケージなのか気になって調べたところ, 上記で提案されたアルゴリズムの実装らしい. ということで調べてみると, 機械学習を経済学で応用することを研究している経済学者らが提案した方法である. よって, (計量) 経済学の知識と機械学習の知識の両方がないと何をやろうとしているのかわかりづらいと思われる. たとえば, 一般化ランダムフォレスト (GRF) という名前は, ランダムフォレストの一般化なのか? 何をどう一般化したのか? と思うかもしれないが, これは計量経済学の文脈を知っているとその名前の意図に気づける. というわけで, 予備知識を補いながら解説していく.
この論文で提案されている GRF をできる限り短くまとめると, 「個別効果の異質性 (heterogeneity) のある因果推論モデルのパラメータ (異質的処置効果: Heterogeneous Treatment Effect) を推定するために, モデルの誤差を基準にした Ramdom Forest で個体をクラスタリングし, 得られた推定クラスタをもとに一般化モーメント法 (GMM) を 局所回帰してパラメータを識別するアルゴリズム」である. これで何をやってるかわかった人はこれ以降読まなくても大丈夫だが, たぶん多くの人は分からないと思うので, 説明をしていく.
予備知識のセットアップ
経済学寄りの人, 機械学習寄りの人のどちらに対しても説明になるように, GRFの理解のために特に重要な概念を順に説明していく*1. 以下では, 「目的意識」「一般化モーメント法」「カーネル回帰」「ランダムフォレスト」について簡単に説明しておく.
目的は因果推論
単に因果推論, と言うだけでは, 文脈を共有しない人には少しわかりにくいと思われる. 特に機械学習寄りの人向けに経済学の問題意識を解説する. というデータにたいして, パラメータ (重みベクトル)
をもつモデル
これは典型的な因果推論である. 以下, 因果推論を形式的に説明しておく. この応用は, 独立変数 の増減に対して目的変数
がどれだけ変化するか, という含意を得たい, ということである*4.
というバイナリ変数とすると, これはある政策 (施策) を打ったかどうかを表す.
のうち,
以外の残りの独立変数
は
としてまとめておく. そこで
ただし, このような解釈は本質的に反実仮想的 (counterfactual) であり, が本当にこのような因果効果の大きさとみなせるかどうかは, さらに追加の仮定が必要になる. これについては, モデルをどう定式化するかだけでなく, サンプルの取り方の問題などが絡んでくる. この仮定を無視した因果推論の「濫用」はしばしば問題となり, 慎重に分析しなければならないが, ここでそこまで詳細に説明するのは面倒*5なので, ここではこの問題をクリアしていると仮定し, 単なるパラメータ同定と同じ問題としてみなすことにする.
一般化モーメント法 (GMM)
一般化モーメント法 (GMM) は計量経済学ではとてもメジャーなパラメータの推定方法であり, 様々な派生バージョンが存在する. しかし, 経済学以外の分野では私の知る範囲ではほとんど見られない.
典型的な回帰モデル,
ここで, や, それ以外の変数が
と無相関ならば, それらを
とすると,
GMMの枠組みではモーメント条件からどうやってパラメータを推定するのか, 簡単に説明する.
モーメント条件は期待値だから、適切な仮定のもとでは標本平均 (標本モーメントともいう) で推定することができる. よって, 標本平均がゼロという条件から逆算してパラメータを推定することが考えられる. もちろん, 実際のデータにはノイズがあるから, 標本平均が正確にゼロとなるわけではない. そこで, GMM によるパラメータの推定は, 距離最小化問題として定式化される. 独立変数が複数ある場合に一般化すると、モーメント条件は, 以下のようなベクトルで表される.
これ以上の GMM の詳しい説明は, 以前書いた記事や, 記事で挙げられている参考文献, あるいはぐぐったら出てくる各大学の講義ノートなどで確認してほしい.
カーネル回帰
カーネル回帰*9あるいはカーネル平滑化法という名称は分野によってばらつきがありそうだ. 統計学/計量経済学では局所回帰 (local regression) やノンパラメトリック回帰と呼ばれることが多い*10ようで, 本文中でも局所回帰という表現が使われている. カーネル回帰は, 機械学習の教科書でもよく取り上げられているが, 機械学習の典型的な問題である予測 (外挿) には使いにくいため, 実際に使う人は少ないかもしれない. 計量経済学では, 部分的に使われることもあるが, 標準的な教科書にはあまり取り上げられることのないトピックのように思われる.
カーネル回帰はナダラヤ・ワトソン (Nadaraya-Watson) 推定量が代表例で, 矩形3次カーネルやエパネチニコフ (Epanechnikov) ・カーネル (Епанечников, В. А., 1969) などのカーネル関数を使った方法のことを指す*11.
局所回帰の基本的なアイディアは, ある点の周囲に分布するデータだけでフィットさせる方法のことである. 例えば, 線形回帰で
その際に, 点が離れるほど重みが減衰するようにしたい. この重み付けをカーネル関数と呼ばれる関数で行うのがカーネル回帰である. 例えばエパネチニコフ・カーネルならば
これは, ある点の周囲で最も近い観測点から順に個を取り出して計算に使うk-最近傍法 (k-NN) と似ている. 再近傍法ベースのアルゴリズムには B-スプライン近似や
ggplot2
の geom_smooth()
でデフォルトで使われる LOESS がある*12.
カーネル回帰は入力された点 ごとにパラメータを計算し, 予測値
を返すという処理をする, 怠惰学習*13と呼ばれるタイプのアルゴリズムである. よって, 一旦パラメータを学習させれば後は任意の
に対してただちに予測値を計算, ということができない. よって, 活用できる場面が限られるが, 上記のように線形回帰モデルを定義していながら滑らかな回帰曲線を得られるという特徴がある. たとえば, 3次曲線をもとに生成させたデータ に, 上記の単純なカーネル線形回帰を適用した結果が以下になる.
同様に, GMM (より一般化するなら, 尤度関数にも) カーネル関数でウエイトをつけて計算することもできる. 本文中では, カーネル関数で重み付けした GMM は local GMM と呼ばれている.
カーネル回帰の詳しい話も, 文献が充実しているので適当に確認してほしい.
例えば, 日本学術会議 経済学委員会 数量的経済・政策分析分科会の発表スライド
西山慶彦 『ノンパラメトリック,セミパラメトリック計量経済分析:概観』
や, 機械学習の分厚い教科書などでは説明がなされていることが多い. 例えばHastie, Tibshriani, & Friedman (2009, The elements of statistical learning: Data mining, inference, and prediction)など.
ランダムフォレスト
経済学寄りの人はおそらくランダムフォレストを知らないことが多いので, 今度は主に経済学寄りの人向けに説明する. 本文にも簡単な解説があるが, 一応.
ランダムフォレストは, パラメータの識別ではなく, 変数 に対して目的変数
の予測値
を求めるためのアルゴリズムである. よって, 予測値
がどれだけ
を正確に予測するかが課題であり, パラメータをバイアスなく推定するといったことは必ずしも要求されない. これはランダムフォレストに限らず, 機械学習全般に当てはまる*14. このランダムフォレストというアルゴリズムも, 線形回帰のような形でパラメータを明示的に求めることなく, 予測値を返すアルゴリズムである.
ランダムフォレストの説明の前に, まず基礎となる決定木について説明する. 決定木は, 条件分岐の組み合わせ構造で記述されるモデルを作成するアルゴリズムである. 決定木と呼ばれるアルゴリズムは実は複数あるが, 一番広く使われているのは CART アルゴリズム*15で, 今回提案される方法でも CART を使っているため, これを説明する. CART は2択分岐をネストさせて作るアルゴリズムで, 分類にも回帰にも適用できる. 決定木という名称は, モデルがノードと枝で構成された木構造で表せることに由来している. この木構造を特定するには, 数学的には組合せ最適問題になり, 全体最適化の計算が難しい. そこで, 局所最適化を繰り返すヒューリスティックな方法 (貪欲法) で計算するしかない. 具体的には, 分岐のたびに, 分岐の誤分類を最小限に抑えるように分割方法を決める, ということである.
例えば以下の図が決定木の例である.
これは という3通りの値を取る目的変数に対して決定木を適用したものである*16. この図を例に, 決定木のアルゴリズムの流れを説明する. (1) は, 木の出発点である根ノードである. このノードにはデータ
が存在する. これを, 2つのデータに分割する. この分類によって,
の誤分類が最小になるように*17分割の条件を決めなければならない. ここでは,
ならば
, そうでなければ
, という条件が選ばれた. この条件式は, 必ず1つの変数をつかった不等号で表される. これによって, データはノード(2) と ノード (3) に分割された. ノード (2) に振り分けられたデータは, 全て
なので, これ以上分割する必要はない. 一方, ノード (3) に分けられたデータには B, C の例が五分五分で含まれている. 今度は
という条件で分割している. これによって
で誤分類率 0.9 % のノード (6) と
で誤分類率 0.2 % のノード (7) が出来た. このようにして, 誤分類がなくなるまで分割する*18.
このようにしてできた木構造に対して, データの例を入力し, 終端のノード (リーフノード) のどれにたどり着くかでの予測値を出力できる.
貪欲法で解かれる決定木は, 計算量が比較的少ないが, 結果のばらつきが大きくなってしまいがちである. そこで, 結果のばらつきの大きい決定木を複数作って, その平均を取る*19ことでばらつきを小さくする, というのがランダムフォレストのアイディアである*20. 実際に, ランダムフォレストと, これを構成するそれぞれの木の決定領域を図示してみる.

ランダムフォレストの決定領域
この図では, tree1
から tree8
までの8つの木と, それらを合成した forest
の決定領域と, 実際のデータの散布図を重ねたものである. 色分けされた領域は, それぞれの木の予測値を表し, 破線は真の境界を表している. 個々の木の決定領域はかなり異なるが, 的中率 (accuracy) はそれぞれ大差ない. つまり, それぞれの木にはばらつきがあるが, それぞれ同じくらい正しい決定領域の近傍をさまよっている. このような時, 複数のモデルの予測値を合成したランダムフォレストの予測する決定領域は、ばらつきのあるそれぞれの木の決定領域を平均するため, より正しい決定領域に近づき, 予測精度が向上する*21. 実際, この図では 8つの木の的中率がそれぞれ 0.63 から 0.87 程度だったものが, ランダムフォレストに合成することで 0.94 にまで上昇している.
かなり大雑把にアイディアだけを述べたが, GRF を説明するのに最低限必要なランダムフォレストへの理解は「2分割を繰り返す」「ブートストラップ標本で複数回計算し、結果の平均を取る」「複雑な分類ができるわりに計算量が少なく, 実装も比較的簡単なアルゴリズム」程度で良いと思う.
本題
以上, 「目的意識」「一般化モーメント法」「カーネル回帰」「ランダムフォレスト」の4つを踏まえて, 論文の主題について説明する.
この論文でも, メインはモーメント条件
主な目的はパラメータ の同定を通して因果効果の大きさ知ることで, それ以外のパラメータ
は, あまり関心のない局外母数 (nuisance parameter) ということになる. パラメータの同定が目的と言ったが, ここでは最低限
さえ同定できればいい*22.
は, そのような2種類のパラメータを含むモデルを表している.
は, 問題に直接関係のあるデータを表し, 例えば,
として, そこに共変量
を加えた,
さらに, より正確には, モーメント条件は,
さて, GMMとは, モーメント条件の期待値を標本平均に置き換えた標本モーメントの距離最小化問題を解くことだった. ここでも同様だが, 実際のデータには異質効果が存在するため, 通常と同じように推定するとパラメータの推定にバイアスがかかる. よって, パラメータを推定するには, 新たに重み を加えた, 以下のような最適化問題を解くことになる (GMM の重みは, 観測例
ごとに異なる値であるという仮定を置いていなかった).
という状況でどうすればいいか, というのがこの論文の問題設定である.
この式は, をカーネル関数とすれば, カーネル回帰ともみなせる. もちろん,
が任意の
で同一であると仮定すれば, 通常の回帰モデルにもなる. さらに,
の部分については特に仮定してないため, 理論上は GMM や最小二乗法だけでなく, 分位点回帰, サポートベクター回帰など色々な回帰アルゴリズムを適用できる.
既存の研究でも, こういった問題に対処する方法が全く提案されなかったわけではない. 例えば経済学分野の生産関数の推定では, Olley & Pakes (1996), Levinsohn & Petrin (2003), Ackerberg, Caves, & Frazer (2006) という一連の研究で, 企業ごとに異なる生産性を表現するために, カーネル多項式回帰が使われている. これらは生産関数の推定という特定のジャンルの問題を対象とした方法なので問題にならなかったが, 一般にカーネル回帰は の次元が増えると近傍の点が減るため (次元の呪い), 計算がうまくいかないという問題がある*23.
そこで, この論文では, をカーネル関数ではなく, ランダムフォレストを使って推定する, という方法を提案している. 決定木で, 分岐を作る際に, 分割による局所回帰の当てはまりを最適化するように決定する. そのようにして作成したランダムフォレストについて, 「分類結果」ではなく, 最終的にどのリーフノードに行き着くか, という情報をクラスタリングに使う*24. 本来のランダムフォレスト (決定木) は, 誤分類を基準に分割するものだが GRF では回帰モデルの平均2乗損失関数を代わりに当てはめている. この点では一般化かもしれないが, 個人的には GMM の派生アルゴリズムという意味ということにしたほうがしっくりくる気がする.
この GRF のアイディアをそのままアルゴリズムに落とし込むと, 次のようになる.
本の決定木について, それぞれ以下を繰り返す.
データの初期値として, ブートストラップ標本を生成する.
親ノード
のデータを,
の2つに分割する. 分割は決定木と同様に, いずれかの軸に平行になるようにする.
分割後の2集合それぞれに対応する部分データで上の最適化問題
を解き, パラメータを得る.
(2) から誤分類率を計算する. 誤分類率は以下のように, パラメータの平均2乗誤差の期待値として定義されるので, これを推定する.
はノード
に属する部分データ,
は親ノード
の時点でのデータ, つまり
を使って (b) の式から推定したパラメータである.
(2), (3) の方法で計算した誤分類率が最小になるような分割
を選択する.
以上を
本の木について行ったら, 結果を合成する.
を,
とする.は, 木
において
と同じリーフノードに行き着く点の集合を意味する. この
をもとに, (A) を解く.
ところで, 本来の決定木やランダムフォレストでは, 誤分類率の計算はこのような複雑ではなかった. ノードの分割のたびに, 最適な分割を探索するためにモデルのパラメータをいちいち推定する必要があるので, 計算量が膨大になり実用的ではない. そこで, この論文では, 計算を早めるために誤分類率を近似計算で置き換える方法を提案している*25.
さらに, この方法で得られた推定量の理論的性質について調べ,
がレギュラーな条件を満たすなら, GRFで推定された
は一致性を持つことや, 漸近分布をを示している. レギュラーな条件は最小二乗法をはじめ, 分位点回帰*26などロバストな回帰など, 多くのアルゴリズムと矛盾しないため, 応用できる範囲は広い.
2019/7/4 追記: 私はこの論文の研究の背景をあまり確認してなかったのだが, 中村知繁氏の以下のスライドの44ページ目以降から読むと, この研究のモチベーションがわかりやすい (スライドの前半部分も因果推論フレームワークの理解に有益である).
speakerdeck.com
参考文献
- Ackerberg, D. a, Caves, K., & Frazer, G. (2006). Structural Identification of Production Functions. Ariel, 12(Iv), 1–35. Retrieved from http://folk.uio.no/rnymoen/Ackerberg_Caves_Frazer.pdf
- Athey, S., Tibshirani, J., & Wager, S. (2016). Generalized Random Forests, doi:https://doi.org/1610.01271
- Hastie, T., Tibshriani, R., & Friedman, J. (2009). The elements of statistical learning: Data mining, inference, and prediction (2nd ed.). Springer. Retrieved from https://web.stanford.edu/~hastie/ElemStatLearn/
(杉山守・井出剛・神嶌敏弘・栗田多喜夫・前田英作他監訳『統計的学習の基礎—データマイニング・推論・予測—』, 共立出版, 2014 年)- 作者:Trevor Hastie,Robert Tibshirani,Jerome Friedman
- 発売日: 2014/06/25
- メディア: 単行本
- Levinsohn, J., & Petrin, A. (2003). Estimating Production Functions Using Inputs to Control for Unobservables. Review of Economic Studies, 70(2), 317–341. doi:https://doi.org/10.1111/1467-937X.00246
- Olley, G. S., & Pakes, A. (1996). The Dynamics of Productivity in the Telecommunications Equipment Industry. Econometrica, 64(6), 1263–1297. doi:https://doi.org/10.2307/2171831
- Zhou, Z.-H. (2012). Ensemble methods: Foundations and algorithms. CRC Press. (宮岡悦良・下川朝有訳,『アンサンブル法による機械学習 — 基礎とアルゴリズム —』,近代科学社,2017 年)
- 作者:Zhou,Zhi‐Hua
- 発売日: 2017/07/04
- メディア: 単行本
- 星野崇宏 (2009) 『調査観察データの統計科学 --因果推論・選択バイアス・データ融合』岩波書店.
調査観察データの統計科学―因果推論・選択バイアス・データ融合 (シリーズ確率と情報の科学)
- 作者:星野 崇宏
- 発売日: 2009/07/29
- メディア: 単行本
- Епанечников, В. А. (1969). Непараметрическая оценка многомерной плотности вероятности. Теория Вероятностей и Ее Применения, 14(1), 156–161. English version: doi:https://doi.org/10.1137/1114019
*1:元論文も経済学者を対象とした論文であるため, 馴染みの薄い機械学習に関しては, 基本的なことがらも順序だてて説明している. しかしそれだけだと論文のリンク貼って「読め」で終わってしまうのでもう少し丁寧に説明しておく.
*2:予測値を得ることが目的なのか, パラメータの同定が目的なのか, という前提を共有しておらずに噛み合わない議論が時々インターネットで見られる (と書いたが自分も以前、論点のはっきりしていない文章を書いていた).
*3:なので, ニューラルネットであっても, その理論的挙動や感度分析の方法論の整備が進めば, もしかしたら使われるようになるかもしれない. あまり自信はないが.
*4:つまり, 一般化して言うと,
を知りたいということになる.
*5:星野 (2009) などが日本語の因果推論の教科書として支持されている.
*6:なお, 操作変数アプローチが純粋な統計学的に見ると間違っている, などという主張がまれに現れるが, 実際のところ, SEM のように複数の方程式を集約したモデルを推定しているというだけであり, 多変量分布で表されるモデルの最尤推定と実質的に変わらない. ではなぜ最尤法を使わないのかと言うと, 理由の1つとして, モデルを線形近似した GMM のほうが, モデルの仮定の誤りに対して最尤法より頑健である, という性質があることが挙げられる. それ以外にも, どういう問題に応用するかで GMM と最尤法の優位点はいろいろと変わってくるため, 一概に論じることは難しい. 具体的なアルゴリズムどうしを比較するのではなく, 「ベイズ統計」や「ディープラーニング」といった広いカテゴリでまとめて, どちらが優位か比較するのがあまり意味をなさないのと同じように.
*7:SEM の定義にコンセンサスがあるのか, 私はよく分かっていない. 非線形なモデルも SEM に含むという人もいるかもしれない.
*8:パラメータの同定という観点から, 変数のスケール調整は好まれない
*9:カーネル密度推定 (KDE) と名前が似ているが, カーネル回帰が教師あり学習とすれば KDE は教師なし学習と位置づけられる.
*10:こういうのもノンパラメトリックと言えるのか昔から疑問に思っている. そもそもノンパラメトリックを厳密にどう定義するかという問題かもしれないが. 頑健性 (robustness) と同じように明確に定義しないものなのだろうか.
*11:名前にカーネルとあるが, サポートベクターマシン (SVM) などに利用されるカーネルトリックとは直接の関係はない.
*12:LOESS などの最近傍法で重み付けするアルゴリズムは, カーネルによる重み付けとは理論上の性質が異なる. よって、狭義には最近傍法ベースのアルゴリズムだけを局所回帰というほうが正確であるように思えるが、カーネル回帰も含めて局所回帰と呼ばれることも多い.
*13:lazy learning.「遅延学習」というほうが適切か?
*14:よって, 機械学習で「バイアス」と言ったら, 通常は予測値と目的変数の差の大きさについて言及しているはずである.
*15:R の rpart や Python の scikit-laern でも CART が採用されている.
*16:実はこれは iris データの名前を変えただけである.
*17:この誤分類の計算には何通りかある. 単純に誤分類の発生割合とすることもできるが, CART では, ジニ多様性指数を使うことが多い.
*18:分割の終了条件にもいろいろなやり方があり, 誤分類がなくなるまで続けるという方法だけが正解というわけではない. 分割を細かくしすぎると過剰適合になる可能性もある.
*19:今回は予測対象がカテゴリカル変数なので, 平均ではなく多数決で最終的な予測値を決める.
*20:より正確には, ランダムフォレストとはブートストラップ標本を生成し, それぞれの標本で決定木を作り, 結果を合成するアルゴリズムのことを言う.
*21:このような手法は集団学習とかアンサンブル学習と総称される. ここでは大雑把な説明をしたが, より厳密な理論を知りたい場合は, Hastie et al. (2009, The elements of statistical learning: Data mining, inference, and prediction) や Zhou (2012, Ensemble methods: Foundations and algorithms) などを読むと良い.
*22:あらゆる変数が制御可能ということはあまりないため, 一部の変数だけでもよいから特定したい, という場合はよくある.
*23:さらに言えば, OP, LP, ACF の論文では, モデルに追加の仮定を設けることで式をシンプルにすることで推定しやすくしている, という面もある
*24:このようなアイディアはこの論文が最初ではない. よく似たアイディアとして, リーフノードを特徴量に使うことで複雑な構造のデータに対しても線形モデルの予測性能を向上できるという方法が https://research.fb.com/publications/practical-lessons-from-predicting-clicks-on-ads-at-facebook/ で紹介されている.
*25:具体的な計算手順まで細かく解説するのはしんどいので, 興味があれば各自で調べてほしい.
*26:quantile regression forest というアルゴリズムが過去に提案されているようだが, これは, 予測値の期待値を平均ではなく中央値に寄せるようにしたランダムフォレストのバリエーションの1つであり, GRF とは直接の関係はない.