1.bayesian networks的一些基本概念
贝叶斯网bayesian networks是一种有向无环图模型(DAG),可表示为G=(V,A)。其中V是节点的集合,节点表示随机变量;A是弧(或称为边)的集合,弧的箭头表示随机变量之间的概率相依性。有向无环图DAG定义了一个因子化的V中全体节点的联合概率分布,称为全局概率分布;相对的,与每个随机变量关联的,为局部概率分布。
因子化的形式由贝叶斯网的马尔科夫性质给出,对每个随机变量,其概率只依赖于其父代:
可以通过学习算法得到贝叶斯网的结构。学习算法首先是学习网络结构,然后在此基础上估计局部分布函数的参数。 尽管全局和局部分布的选择形式很多,最常用的还是下面的两种分布: 多项分布(离散数据) 多元正态分布(连续分布)
得到网络结构,可以利用概率的条件相依,对相关专业问题进行推断。
R中有多个包可以实现贝叶斯网的创建、学习和推断,详见下面表格。
前面我曾简单介绍过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)