Rで決定木・回帰木モデルをつくるとき、標準のrpart::rpart()を使ってもいいし、他にもパッケージが山ほどある。なのに、巷の解説をみると、なぜかmvpartパッケージを使って説明していることが多い。すでに公開停止されたパッケージなのに、なぜか今でも広く使われているように思う。なぜだろう？ そんなに使いやすいわけでもないんですけど。
かくいう私も、ふだんの仕事ではmvpartパッケージを使うことが多い。最大の理由は慣れているからなんですが(嗚呼、老害)、もうひとつ、ちょっとした理由がある。ふつう回帰木の目的変数はひとつだけど、複数の目的変数について同時にやりたくなることがありませんか？ 私の場合でいうと、消費者の行動を表す量的変数が複数あって、ソシオグラフィカルな属性が複数あって、行動が似通っている消費者群を属性の組み合わせで定義したい、というような場合だ。こういう多変量回帰木には、標準のrpartパッケージでは対応できないし、他のパッケージもふつうは対応できてないんじゃないかと思う。たぶん。知らんけど。
でもmvpartパッケージってろくな解説もないし、細かい点についてちょっと疑問が生じても解決する術がない。どうしたものかなあ… と思っていて、ふと気が付いたんだけど、このパッケージの作者はGlen De’athさん。多変量回帰木を最初に提唱したといわれる論文De’ath(2002)の著者その人なのである。なあんだ、そうだったの？ じゃあ元の論文をよめばいいじゃん。
De’ath, G. (2002) Multivariate regression trees: A new technique for modeling species-environment relationships. Ecology, 83(4), 1105-1117.
というわけで、読んでみました。馴染みのない分野の論文で、よくわからん点も多いのだが、そこはあまり気にしないことにして…
1. イントロダクション
環境と種の関係をどうやってモデル化するか。
ひとつの種についてなら、GLMとかGAMとかCARTとかが使われており、CARTの予測性能が高いと云われている。アンバンスとか欠損とか非線形関係とか高次交互作用とかがあるからだろう。
コミュニティの分析、つまり複数の種についての分析の場合は、オーディネーション手法が広く使われている[よくわかんないけど、観測地点を行、種と環境変数を列、セルに観測値をいれたデータ行列を主成分分析とか多重対応分析とかで視覚化するというようなことかな？]。そういう分析はふつう、環境変数の線形結合で表される環境成分に対して、個々の種のabundances[なんて訳すの?]が線形な関係を持っている、ないし単峰の関係を持っている、という仮定をおいている。ほんとに同時分析することもあれば、まず種のほうの構造を決めて、次に環境変数との関連を調べることもある。
クラスタ分析も使われてますね。その場合は、まず種によって測定地点のクラスタを作っておいて、つぎに環境変数と関連付ける。同時分析の試みは少ない。
本研究は、種と環境の関連を調べるための新手法として多変量回帰木 MRT を提案する。複数の種について分析するが、環境との関係の形式についてはなにも仮定していない、というのがポイント。クラスタ分析としてみると、種についてなんらかの形で類似したクラスタを作るんだけど、クラスタの定義は環境変数で行うわけで、一種の制約付きクラスタリングである。
なお、以下では反応変数のことを種、事例のことを観測地点、説明変数のことを環境特性と呼びますが、別にそうでなくてもいいです。
2. 単変量回帰木
[単変量回帰木 URT についての説明。中略…]
URTの場合、ノードのimpurityはふつう、反応変数のノード平均に関する二乗和として定義する[偏差二乗和ということね]。分割についてはノード内二乗和の最小化、言い換えるとノード間二乗和の最大化が目指される。
[…中略…]
木の適合度として、葉のimpurityの合計をルート・ノードのimpurityで割った値を用いる。これを相対誤差 RE という。[え？ それって\(R^2\)と同じことだという理解であってる？] 予測という観点からは RE は楽観的になっているので、交差妥当化相対誤差 CVRE を使う。
3. MRTへの拡張
多変量への拡張にあたっては、単に、ノードのimpurityを、単変量のimpurityの和として定義しなおすだけでよい。幾何的に言うと、impurityはノード重心に対する観測地点のユークリッド距離の二乗和である。
以下ではこれをSS-MRTと呼ぼう。
さて、単変量から多変量への拡張にあたって、次の難題が生じる。
- どうやって解釈するか。4節で述べる。
- MRTは制約付きクラスタリングの一形式である。これを、同じクラスタ数の無制約なクラスタリングと比べたい。比較に際してはMRTでいうimpurityと同じメトリックが使えると助かる。で、もしも無制約なクラスタリングのほうが種の分散をよりよく説明できているならば、未観察の説明変数があるということになる。[←なるほどね]
- 上記のimpurity指標を拡張し、現在オーディネーションとかクラスタリングとかで使われている種の非類似性の指標と同じものできるとありがたい。5節で述べる。
4.MRTの初歩的事例
12種類のhunting spidersと6個の環境変数との関係について調べよう。[確認してないけど、これRのmvpartパッケージについてるデータじゃないかな]
まず、種を標準化する。次に観測地点を標準化する。[え？ データ行列を最初に縦に標準化し、次に横に標準化するってこと？] 環境成分の分析ではふつうそうするし、あとで冗長性分析 RDA や正準対応分析 CCAと比べることができる。
MRTにより葉4枚の木が得られた。[…中略…]
4.1 分割についての検討
木における種の分布は各ノードの棒グラフをみてもわかるけど、種を表側、分岐を表頭にとって、各分岐で各種の分散がどれだけ減ったかを表にしてもよい。[あああ、なるほど。これmvpart()はぱっと出してくれたりするの？]
4.2 ツリー・バイプロット
[…] MRTでは、葉の平均反応は多変量であり、葉を比較するのがURTより難しい。多変量データの構造について検討する方法のひとつは、主成分バイプロットによる低次元空間上での視覚化である。距離バイプロットなら、群の平均の最良の最小二乗表現を与えてくれるので、MRT(群間の二乗和を最大化している)と同じメトリックを使っていることになる。[…]
[距離バイプロットってなんのこと？ いわゆる共分散バイプロットのことだろうか。まあいいや、とにかく、葉を行、説明変数を列にとった集計表をバイプロット表現しておいて、そこに個々の観測地をsupplimentary pointとして重ねる、というアイデアだと思う]
4.3 indicator種の同定
[その葉に特有な種を見つけるという話。これはちょっとdomain-specificな話だと思うのでパス]
4.4 木の葉と無制約なクラスタの比較
同じデータについて、完全結合法の階層クラスタ分析をやって、k-means法をやりなおしてクラスタ内平方和を最小化した。比べてみるとMRTと似た結果になった。つまり環境変数で種の分散をうまく説明できている。[うーん。MRTと全然違う結果になったとして、ただそれだけで、環境変数で種の分散をうまく説明できていない、という言い方ができるだろうか？ 二乗和を比べるってんならわかるけど]
5. SS-MRTを超えて
回帰木というのは、木の成長と剪定においてimpurityと分割基準と予測誤差がそれぞれ独立だという特徴を持っている。つまり、impurityだけ別のに切り替えてもいいわけだ。
5.1 加算多変量回帰木
平均に対する二乗和じゃなくて、中央値からの差の絶対値の合計を使ってもいい。こういうimpurityを加法距離という。加法距離をつかったMRTをA-MRTと呼ぼう。外れ値に頑健だという点が魅力である。
5.2 距離ベース多変量回帰木
[面白そうな話なんだけど、いまちょっと時間がないのでパス…]
5.3 A-MRTとdb-MRTの比較
[パス]
6. MRT, 冗長性分析、正準対応分析の比較
[パス]
7. 考察
[パス]
———
残念ながら別の用事で忙しくなって、途中までしか読んでないんだけど、まあ目的は達したから良しとしよう。
Appendixをみると、SS-MRTでは、観測地点\(i\)における種\(j\)の値を\(x_{ij}\)、平均を\(\bar{x}_j\)として、impurityは \(\sum_{ij} (x_{ij} – \bar{x}_j)^2\)だと書いてある。要するにMRTにおける最小化の対象というのは、各反応変数のクラスタ内偏差平方和の、反応変数を通した単純合計なわけね。それはそれでクリアな定義である。反応変数ごとに事前にスケーリングしとくかどうか、ユーザが熟慮のうえ決めなさい、ということであろう。この点、mvpart()の実装はどうなっているのかしらん、あとで調べてみよう。
統計手法の提案であるにも関わらず本文には数式が一切出てこないという、ちょっと変わった論文であった。心がなごみますね。