データサイエンティスト(仮)

元素粒子論博士。今はデータサイエンティスト(仮)。

Rでスパースモデリング:Group Lasso

はじめに

スパース推定の代表例として、Lasso(Least Absolute Shrinkage and Selection Operator)があります。
これは、正則化項として L_1 ノルムを加えるようなパラメータ推定法です。

 \displaystyle {
  \Omega(\bf{w}) = \| \bf{w} \|_1 = \sum _{i = 1}^p |w_i| 
  }

このような正則化項を導入すると、大部分のパラメータは0に潰れてスパースな解が得られます。

それに加えて、近年、Lassoは様々な拡張モデルが提案されています。
その中の一つが、グループ化された変数がある場合に適用されるGroup Lassoというモデルです。
各グループを指定するラベルを gとすると、Group Lassoの罰則項は以下のように書かれます。

 \displaystyle {
  \Omega_{\text{GLasso}}(\bf{w}) = \sum _{g} \|w_g \| _2 
  }

ここで、各グループに関するノルムは L_1ではなく L_2ノルムです。
このように設定することにより、同じグループ内の変数がすべて残る/一斉に0になるという傾向をもつようになります。このあたりの直感的(図形的)説明は、Group Lassoでグループごと重みが0に潰れる理由 | Preferred Researchを参照してください。

Group Lassoの使いどころ

カテゴリカル変数そのものの効果をモデルに反映させたい

例えば、都道府県を説明変数に入れる場合があります。素朴には、各都道府県をフラグ化(and 適切にreferenceを設定)して、変数選択やLassoなどでどのような都道府県がモデルの変数として選択させる、などを調べることがあると思います。しかしながら、この方法だと、「都道府県じたい」が効果のある変数かどうかを判定したい場合にはうまく働きません。このようなあるカテゴリカル変数じたいの効果の有無を調べる場合にGroup Lassoは効果を発揮します。

先見的に似た変数の効果がわかっていて、それらを同時にモデルに反映させたい

遺伝子情報などの中には、あらかじめ効果がわかっているものがあることが少なくありません。特に、似たような効果をもつ遺伝子を同時にモデルに反映させたい場合に、Group Lassoが有効です。

Group Lassoの種類

Group Lassoには大雑把に重複なしGroup Lasso、重複ありGroup Lassoの2種類があります。前者は説明変数間で重複がないようにグループ化されている場合で、後者は説明変数間で重複を許すようにグループ化されている場合になります。この2種類は、どういった目的のモデルを作成したいかで使い分ける形になるでしょう。

RでGroup Lasso

Rに限っただけでもLassoの拡張モデルを扱えるパッケージはたくさんあります。ここでは、{grpreg}というパッケージを使ってGroup Lassoの挙動を見てみましょう。

{grpreg}のreferenceは以下になります。
https://cran.r-project.org/web/packages/grpreg/grpreg.pdf

パッケージをインストールし、サンプルデータ`Birthwt`を使用します。

if(!require(grpreg)){
    install.packages("grpreg", dependencies = TRUE)
}; require(grpreg)

data(Birthwt)

ここで、Birthwtの中身はざっくり以下になっています。

  • X : matrix 説明変数(妊婦の属性) 標準化の操作は済
  • bwt : numeric 乳児の体重[kg] 線形回帰とかを想定した変数
  • low : integer 低体重(2.5kgがthreshold)なら1, そうでなければ0 ロジスティック回帰とかを想定した変数
  • group : factor あらかじめ変数を手でグループ分けしておいた結果を格納

説明変数の雰囲気を見ておきます。

f:id:tekenuko:20160829225417p:plain
それぞれの変数は以下です。

  • age1~3 : 妊婦の年齢に関連する変数
  • lwt1~3 : 妊婦の体重に関連する変数
  • white, black : 妊婦の人種 その他(黄色人種とか)がreference
  • smoke : 妊婦が喫煙歴あるか
  • ptl1, ptl2m : 早産の回数フラグ 1回、2回以上
  • ht : 高血圧の有無
  • ut : 子宮過敏性の有無
  • ftv1, 2, 3m : 妊娠の3ヶ月目までで医者を呼んだ回数フラグ 1回、2回、3回以上


Group Lassoで乳児の体重を妊婦の状態から説明するモデルを構築してみます。モデル構築には、grpreg関数を利用します。

model.grplasso <- grpreg(X = Exp, y = Obj, group, penalty="grLasso")

grpreg関数は、penaltyで罰則項の選択が可能です。grLassoはいわゆる普通のGroup Lassoで、他にはGroup SCAD, MCP, Bridgeなども選べる余地があります。さらに、{glnmet}のようにalphaをいじってElastic net的なモデリングも可能になっています。残差分布はGausian, Binomial, Poissonから選択できます。罰則の強さは何点かでグリッドサーチして求めており、基準はBIC (free energyの鞍点近似)がデフォルトになっています。ただし、自分で細かい調整がしたい場合はcv.grpregという関数を使うとより詳細なグリッドサーチが可能です。

実際に選択された変数を見てみましょう。

coef(model.grplasso)
(Intercept)
3.02538345542349
age1
0.12939932166821
age2
0.514120194706666
age3
0.307216803283638
lwt1
0.628666306376423
lwt2
-0.151307766326747
lwt3
0.49392658258674
white
0.167709378969912
black
-0.0534652473427624
smoke
-0.174393102956969
ptl1
-0.157879134519634
ptl2m
0.0444118454553915
ht
-0.26668574037556
ui
-0.36973322736348
ftv1
0
ftv2
0
ftv3m
0

確かに、グループ単位で0になる/ならないが表現されているモデルになっています。

おわりに

カテゴリカル変数の効果を見たい、先見的に似たような変数がわかっていてそれらの効果を取り込みたい、といった状況は実社会でもたびたびあると思われます。そのような場合は、Lassoだけでなく、Group Lassoを試してみると、新たな価値が生まれるかもしれません。