ページ

2016年8月10日水曜日

【R】 地域志向と就職希望地 : 連鎖グラフによる分析(その2)


前回は、連鎖グラフを得るところまでできた。しかし、本当は予測まで行いたいところ。しかしながら、前回使ったパッケージでは、予測に手間がかかることは分かった(出来ないわけではない)。

○方針
  1. 同時分布の共分散行列を取得し、それからシュール補行列の方法を用いて期待値を計算する。(以前やったことがある。すべて量的変数として扱うので可能。)
  2. 新しいパッケージbnlearnを使う。これは試したことがないが、質的変数も量的変数も扱える強力なパッケージ。

今回は2番目の、bnlearnパッケージを使う方法でアプローチしてみる。まずは、インストールして使い方に慣れる必要がある。



bnlearnパッケージのインストール

gRimパッケージはインストール方法が特殊であり、バイオコンダクターのbiolite関数を使う必要があった。bnlearnは普通にCRANからインストールすれば良い。


bnlearnを試してみる

> library(bnlearn)
 警告メッセージ: 
 パッケージ ‘bnlearn’ はバージョン 3.2.5 の R の下で造られました  
> data(learning.test)
> str(learning.test)
'data.frame':   5000 obs. of  6 variables:
 $ A: Factor w/ 3 levels "a","b","c": 2 2 1 1 1 3 3 2 2 2 ...
 $ B: Factor w/ 3 levels "a","b","c": 3 1 1 1 1 3 3 2 2 1 ...
 $ C: Factor w/ 3 levels "a","b","c": 2 3 1 1 2 1 2 1 2 2 ...
 $ D: Factor w/ 3 levels "a","b","c": 1 1 1 1 3 3 3 2 1 1 ...
 $ E: Factor w/ 3 levels "a","b","c": 2 2 1 2 1 3 3 2 3 1 ...
 $ F: Factor w/ 2 levels "a","b": 2 2 1 2 1 1 1 2 1 1 ...
> 

まずはライブラリを呼び出して、サンプルデータを読み込んでみる。strで内容を表示させてみると、すべてカテゴリー変数で構成されているようだ。

res = gs(learning.test)
plot(res)

bnlearnパッケージではデータに基づいてグラフを推定することを「Structure learning」と呼んでいる。得られるグラフはUGではなく、一挙にDAGまで行くようだ。複数のアルゴリズムがサポートされており、それぞれ別種の関数になっている。簡単なところでは、

  • gs:Grow-Shrink
  • iamb:Incremental Association Markov Blanket
  • fast.iamb:Fast Incremental Association
  • inter.iamb:Interleaved Incremental Association
  • hc:Hill Climbing
  • tabu:Tabu Search

(なぜ、簡単なところでは・・・と言っているかというと、他にもハイブリッド方式とか、追加オプションとかいろいろ入っているため。)

Grow-Shrinkアルゴリズムで計算して、グラフを表示させると冒頭の図が得られる。またresというオブジェクトの内容は以下のようにbnlearn特有のオブジェクトになっている。

> res

  Bayesian network learned via Constraint-based methods

  model:
    [partially directed graph]
  nodes:                                 6 
  arcs:                                  5 
    undirected arcs:                     1 
    directed arcs:                       4 
  average markov blanket size:           2.33 
  average neighbourhood size:            1.67 
  average branching factor:              0.67 

  learning algorithm:                    Grow-Shrink 
  conditional independence test:         Mutual Information (disc.) 
  alpha threshold:                       0.05 
  tests used in the learning procedure:  60 
  optimized:                             TRUE 

bnlearnオブジェクトの構築

さて、私のやりたいことは、上のサンプルのようにデータから一挙にStructure learningすることではなく、空のオブジェクトの中に手でStrucyureを構築して、これとデータとを合わせてParameter learningすることなのであった。Structure learningで期待する結果が得られれば、それに越したことはないのだが。それは後で比較してみれば良い。

まずは空のオブジェクトの作り方。

> vnames <- c(c1lab,c2lab,c3lab)
> e = empty.graph(vnames)
> e

  Random/Generated Bayesian network

  model:
   [sex][from][grd][a1][a2][a3][b1][b2][q5] 
  nodes:                                 9 
  arcs:                                  0 
    undirected arcs:                     0 
    directed arcs:                       0 
  average markov blanket size:           0.00 
  average neighbourhood size:            0.00 
  average branching factor:              0.00 

  generation algorithm:                  Empty 

変数名だけ指定してempty.graphという関数で空のオブジェクトが生成される。枝が設定されていないので、図示すると9つの変数がばらばらになっているグラフが表示される。

ここに枝を入れ込む。igD3は見かけだけ連鎖グラフなので、改めて反対向きの枝を追加することを忘れない。

> arc.set <- matrix(vnames[t(get.edgelist(igD3))],ncol = 2, byrow = TRUE,dimnames = list(NULL, c("from", "to")))
> intraarc1 <- arc.set[arc.set[,1] %in% c1lab & arc.set[,2] %in% c1lab,]
> intraarc2 <- arc.set[arc.set[,1] %in% c2lab & arc.set[,2] %in% c2lab,]
> arc.set <- rbind(arc.set,intraarc1[,c(2,1)],intraarc2[,c(2,1)])
> arcs(e) = arc.set
> e

  Random/Generated Bayesian network

  model:
    [partially directed graph]
  nodes:                                 9 
  arcs:                                  10 
    undirected arcs:                     5 
    directed arcs:                       5 
  average markov blanket size:           2.22 
  average neighbourhood size:            2.22 
  average branching factor:              0.56 

  generation algorithm:                  Empty  

これを図示すると、以下のようになる。前回の結果と同じ(当然だが)。


ここまでうまくいったように思えたのだが・・・bnlernのパラメータ・ラーニングはDAGのみを受け付ける(PDAGは、したがって連鎖グラフは受け付けない)ということが判明した。

やはり、第一のアプローチによるしかないようである。次回。