[R program]ベイジアンネットワーク

ベイジアンネットワークを推定する方法はいくつかあるけど、個人的に好きなのはPearl先生のICアルゴリズム（Inductive Causation algorithm）での方法。

教科書はこちら。

ICアルゴリズムに似ている、GSアルゴリズムはbnlearnパッケージで出来ます。

こんな感じ。

今回の擬似データでは、Height・SBP→BMI→FBSという因果関係を作ったのだけど、BMI→FBSには矢線が付かなかった（向きが付かなかった）。

GSアルゴリズムの基準に合わなかったのだろうか。

まぁベイジアンネットはもともと探索的なものだから、線で結ばれてさえいれば、あとは解釈で因果関係を探ればいいんだけど。

ちなみにICアルゴリズム自体はグラフ理論に裏づけされていて、結構単純で以下のようなフローになっている（GSアルゴリズムはまた少し違う基準です）。

• ステップ1：(a 独立 b|Sab)となるSabをみつける（Sabは変数集合、subではなくでSのab）
  • みつからない場合（どうあってもaとbに関連がある場合）はaとbを無向辺で結ぶ

• ステップ2：隣接しないaとbが共通の隣接点cをもつとき
  • c in Sabの場合：矢印を加えない
  • c not in Sabの場合：a→c←bとする

• ステップ3：無向辺にできるだけ向きを付ける(Meek, 1995)
  • 規則1：a→b―cであればa→b→cとする
  • 規則2：a→c→b、a―bであればa→bとする
  • 規則3：a―c→b、a―d→b、a―bであればa→bとする
  • 規則4：a―c→d、c→d→b、a―bであればa→bとする

単純なアルゴリズムだけど、統計的な関連だけじゃなく因果関係をデータだけから導ける（かもしれない）という大変秀逸な理論です。

ちなみにdealパッケージでやるとこうなりました。

多分dealは単純にモデルから得た値（AICとか？）だけで矢印を引いていて、GSアルゴリズムのようにグラフなど理論の裏づけはないと思う。

どっちの方法でベイジアンネットを作るかは、数値計算から攻めるのか、グラフ理論から攻めるのかの違いですかね。

今回の結果では、GSアルゴリズムの方が真の状態に近い結果を出力してますね。

Rプログラムはこちら↓

library(bnlearn)
example(plot.bn)
example(arc.strength)

norm <- rnorm(4000)
Data <- matrix(norm, nrow=1000, ncol=4, byrow=T)
colnames(Data) <- c("Height", "BMI", "SBP", "FBS")
Data <- data.frame(Data)
Data2 <- Data
Data2$Height <- 170 + Data$Height*10
Data2$SBP <- 120 + Data$SBP*10
Data2$BMI <- 22 + Data$Height*2 + Data$BMI*2 + Data$SBP*2
Data2$FBS <- 90 + (Data2$BMI-22)*5 + Data$FBS*10
Data2 <- data.frame(Data2)
res.Data2 = gs(Data2)
plot(res.Data2)
res.arc2 <- set.arc(res.Data2, "Height", "BMI")
arc.strength(res.arc2, Data2)

#deal
library(deal)
pr_net <- network(Data2)
pr_dist <- jointprior(pr_net, 20)
dist_renew <- learn(pr_net, Data2, pr_dist)
pr_net2 <- getnetwork(dist_renew)
pos_net <- autosearch(pr_net2, Data2, pr_dist)


#相関係数と偏相関係数の確認
library(Rcmdr)
cor(Data2)
partial.cor(Data2)