1.bayesian networks的一些基本概念 

贝叶斯网bayesian networks是一种有向无环图模型(DAG),可表示为G=(V,A)。其中V是节点的集合,节点表示随机变量;A是弧(或称为边)的集合,弧的箭头表示随机变量之间的概率相依性。有向无环图DAG定义了一个因子化的V中全体节点的联合概率分布,称为全局概率分布;相对的,与每个随机变量关联的,为局部概率分布。

因子化的形式由贝叶斯网的马尔科夫性质给出,对每个随机变量,其概率只依赖于其父代:

R语言 BIC 贝叶斯评分 r语言贝叶斯网络的实现_概率分布

可以通过学习算法得到贝叶斯网的结构。学习算法首先是学习网络结构,然后在此基础上估计局部分布函数的参数。 尽管全局和局部分布的选择形式很多,最常用的还是下面的两种分布: 多项分布(离散数据) 多元正态分布(连续分布)

得到网络结构,可以利用概率的条件相依,对相关专业问题进行推断。

R中有多个包可以实现贝叶斯网的创建、学习和推断,详见下面表格。

R语言 BIC 贝叶斯评分 r语言贝叶斯网络的实现_有向无环图_02

前面我曾简单介绍过gRain包,本篇讨论最常用的bnlearn包。

2.网络创建和操作 bnlearn包自带数据集marks,88学生5门课的成绩。这个数据最早由Mardia et al(1979)研究过。

library(bnlearn)

## Warning: package 'bnlearn' was built under R version 3.0.1

data(marks)
str(marks)

## 'data.frame': 88 obs. of 5 variables:
## $ MECH: num 77 63 75 55 63 53 51 59 62 64 ...
## $ VECT: num 82 78 73 72 63 61 67 70 60 72 ...
## $ ALG : num 67 80 71 63 65 72 65 68 58 60 ...
## $ ANL : num 67 70 66 70 70 64 65 62 62 62 ...
## $ STAT: num 81 81 81 68 63 73 68 56 70 45 ...



创建一个空网络,节点对应于marks的变量。然后通过指派一个两列的矩阵来添加边。

生成一个无向图:

ug <- empty.graph(names(marks))
arcs(ug, ignore.cycles = TRUE) = matrix(c("MECH", "VECT", "MECH", "ALG", "VECT", 
    "MECH", "VECT", "ALG", "ALG", "MECH", "ALG", "VECT", "ALG", "ANL", "ALG", 
    "STAT", "ANL", "ALG", "ANL", "STAT", "STAT", "ALG", "STAT", "ANL"), ncol = 2, 
    byrow = TRUE, dimnames = list(c(), c("from", "to")))
ug

## 
## Random/Generated Bayesian network
## 
## model:
## [undirected graph]
## nodes: 5 
## arcs: 6 
## undirected arcs: 6 
## directed arcs: 0 
## average markov blanket size: 2.40 
## average neighbourhood size: 2.40 
## average branching factor: 0.00 
## 
## generation algorithm: Empty

这个ug对象属于bn类,这个类用于在bnlearn包中管理网络结构。 这个对象包括三个方面的信息:
 (1)learning:结构的学习 (2)node 节点 (3)arc 边



生成一个有向图:

dg <- empty.graph(names(marks))
arcs(dg) = matrix(c("VECT", "MECH", "ALG", "MECH", "ALG", "VECT", "ANL", "ALG", 
    "STAT", "ALG", "STAT", "ANL"), ncol = 2, byrow = TRUE, dimnames = list(c(), 
    c("from", "to")))
dg

## 
## Random/Generated Bayesian network
## 
## model:
## [STAT][ANL|STAT][ALG|ANL:STAT][VECT|ALG][MECH|VECT:ALG] 
## nodes: 5 
## arcs: 6 
## undirected arcs: 0 
## directed arcs: 6 
## average markov blanket size: 2.40 
## average neighbourhood size: 2.40 
## average branching factor: 1.20 
## 
## generation algorithm: Empty

有时候也可以从邻接矩阵(adjacency matrix)来生成有向图dg。

mat <- matrix(c(0, 1, 1, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 1, 
    0, 0, 0, 0, 0), nrow = 5, dimnames = list(nodes(dg), nodes(dg)))
mat

## MECH VECT ALG ANL STAT
## MECH 0 0 0 0 0
## VECT 1 0 0 0 0
## ALG 1 1 0 0 0
## ANL 0 0 1 0 0
## STAT 0 0 1 1 0

dg2 <- empty.graph(nodes(dg))
amat(dg2) <- mat
all.equal(dg, dg2)

## [1] TRUE

手工修改一个已经存在的网络,利用加边(set.arc)、去边(drop.arc)、颠倒(rev.arc)这几个操作,也可以得到一个需要的网络。

dg3 <- empty.graph(nodes(dg))
dg3 <- set.arc(dg3, "VECT", "MECH")
dg3 <- set.arc(dg3, "ALG", "MECH")
dg3 <- set.arc(dg3, "ALG", "VECT")
dg3 <- set.arc(dg3, "ANL", "ALG")
dg3 <- set.arc(dg3, "STAT", "ALG")
dg3 <- set.arc(dg3, "STAT", "ANL")
all.equal(dg, dg3)

## [1] TRUE

all.equal(ug, moral(dg)) #moral图

## [1] TRUE

我们全面希望了解网络的结构,可以使用存储在每个节点的信息。 
(1)节点的拓扑顺序

node.ordering(dg)

## [1] "STAT" "ANL" "ALG" "VECT" "MECH"

(2)节点的邻居(nbr)和Markov毯(mb) 
Markov毯是指对一个节点A,它的所有父节点,子节点以及和A有相同子节点的其它节点 。

nbr(dg, "STAT")

## [1] "ALG" "ANL"

mb(dg, "STAT")

## [1] "ALG" "ANL"

"ANL" %in% mb(dg, "STAT")

## [1] TRUE

"STAT" %in% mb(dg, "ANL")

## [1] TRUE

(3)某个给定节点的子代(child)和父代(parents)和子代的其它父代(o.par)

child <- children(dg, "STAT")
pa <- parents(dg, "STAT")


2.绘制网络

 bnlearn包有两种对bn对象绘制网络结构的方法,一种是利用graph包和Rgraphviz包提供的接口给出的一些绘图函数,比如graphviz.plot。再一种方法是利用plot函数。下面对于dg,来绘制网络图。 使用graphviz.plot的好处在于它返回一个graph对象,便于对网络的进一步操作。

library(Rgraphviz)

graphviz.plot(dg, layout = "fdp")plot(dg, radius = 200, arrow = 30)