BN3 for QC Notebook. Third approach
FIRST LEVEL BN SETUP.
library(bnlearn)
library(gRain)## Loading required package: gRbase##
## Attaching package: 'gRbase'## The following objects are masked from 'package:bnlearn':
##
## ancestors, children, parentslibrary(ggplot2)
require(scales)## Loading required package: scales##
## Attaching package: 'scales'## The following object is masked from 'package:gRbase':
##
## ordinalDefining the level “l” BN
imagename
bnl_l = function(a,b,c) {
#
# a pr(A1 = 'F')
# b pr(B1 = 'F' | A1 = 'C')
# c pr(B1 = 'F' | A1 = 'F')
#
fc <- c("Failure","Correct")
# bnl = model2network('[A1][B1|A1]')
Astp = cptable(~A1, values=c(a,1-a), levels=fc)
Bstp = cptable(~B1+A1, values=c(c,1-c,b,1-b), levels=fc)
cpts = compileCPT(list(Astp,Bstp))
bnl.stp = grain(cpts)
return(list('net'=bnl.stp, 'parms'=c(a,b,c)))
}Defining the level “l+1” BN
The model now looks like A->B->C + B0->B + C0->C
bnl_lp2 = function(a,b,c,d,e,f,g,h,i,b0,c0) {
#
# a pr(A2 = 'F')
# b0 pr(B0 = 'F')
# c0 pr(C0 = 'F')
# b pr(B2 = 'F' | A2 = 'F', B0 = 'F')
# c pr(B2 = 'F' | A2 = 'F', B0 = 'C')
# d pr(B2 = 'F' | A2 = 'C', B0 = 'F')
# e pr(B2 = 'F' | A2 = 'C', B0 = 'C')
# f pr(C2 = 'F' | B2 = 'F', C0 = 'F')
# g pr(C2 = 'F' | B2 = 'F', C0 = 'C')
# h pr(C2 = 'F' | B2 = 'C', C0 = 'F')
# i pr(C2 = 'F' | B2 = 'C', C0 = 'C')
#
fc <- c("Failure","Correct")
Astp = cptable(~A2, values=c(a,1-a),levels=fc)
B0stp= cptable(~B0, values=c(b0,1-b0),levels=fc)
C0stp= cptable(~C0, values=c(c0,1-c0),levels=fc)
Bstp = cptable(~B2+A2+B0, values=c(b,1-b,c,1-c,d,1-d,e,1-e),levels=fc)
Cstp = cptable(~C2+B2+C0, values=c(f,1-f,g,1-g,h,1-h,i,1-i),levels=fc)
cpts = compileCPT(list(Astp,B0stp,Bstp,C0stp,Cstp))
bnl.stp = grain(cpts)
return(list('net'=bnl.stp, 'parms'=c(a,b,c,d,e,f,g,h,i,b0,c0)))
}Defining the level “l” BN
#
layerl = function(t,t1) {
bnl0 = bnl_l(t,t,t)[['net']]
bnl0c = compile(bnl0)
def = c('Failure','')
nodes = c('A1','B1')
sev1 = setEvidence(object=bnl0c,nodes=nodes,states=def)
p1 = pEvidence(sev1)
#
bnl1 = bnl_l(t1,t,t)[['net']]
bnl1c = compile(bnl1)
sev2 = setEvidence(object=bnl1c,nodes=nodes,states=def)
p2 = pEvidence(sev2)
return(c(p1,p2))
}
#Failure in C | to Failures in A, B0 and C0
#
t0 =0.01
t =0.5
lowp =c(0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,0.99)
prl =c()
prlp1=c()
prlp2=c()
prlp3=c()
prlp4=c()
prp1 =c()
prp2 =c()
for (i in 1:length(lowp)){
resl = layerl(t0,lowp[i])
p1 = resl[1]
p2 = resl[2]
prp1[i]=p1
prp2[i]=p2
# a b c d e f g h i b0 c0
bnlp1 = bnl_lp2(t,t,t,t,t,t,t,t,t,p1,p1)[['net']]
bnlp1c = compile(bnlp1)
nodes = c('A2','B0','B2','C0','C2')
def = c('Failure','Failure','','Failure','Failure')
sev1 = setEvidence(object=bnlp1c,nodes=nodes,states=def)
prlp1[i]= pEvidence(sev1)
#
# a b c d e f g h i b0 c0
bnlp2 = bnl_lp2(t,t,t,t,t,t,t,t,t,p2,p1)[['net']]
bnlp2c = compile(bnlp2)
sev2 = setEvidence(object=bnlp2c,nodes=nodes,states=def)
prlp2[i]= pEvidence(sev2)
#
# a b c d e f g h i b0 c0
bnlp3 = bnl_lp2(t,t,t,t,t,t,t,t,t,p1,p2)[['net']]
bnlp3c = compile(bnlp3)
sev3 = setEvidence(object=bnlp3c,nodes=nodes,states=def)
prlp3[i]= pEvidence(sev3)
#
# a b c d e f g h i b0 c0
bnlp4 = bnl_lp2(t,t,t,t,t,t,t,t,t,p2,p2)[['net']]
bnlp4c = compile(bnlp4)
sev4 = setEvidence(object=bnlp4c,nodes=nodes,states=def)
prlp4[i]= pEvidence(sev4)
}
#
df = data.frame(x=lowp,base=prlp1,ex_b0=prlp2,ex_c0=prlp3,
ex_all=prlp4,pr1=prp1,pr2=prp2)
##
p = ggplot(df, aes(x=x)) +
geom_line(aes(y = pr1, color = "cyan"), linetype="twodash") +
geom_line(aes(y = pr2, color = "orange"), linetype="twodash") +
geom_line(aes(y = base, color = "black")) +
geom_line(aes(y = ex_b0, color="green")) +
# geom_line(aes(y = ex_c0), color="darkblue") +
geom_line(aes(y = ex_all, color="red")) +
scale_y_continuous(trans='log10',
labels = trans_format("log10", math_format(10^.x))) +
theme(legend.position="top") +
scale_color_identity(name = "P(C=F|A=B0=C0=F)",
breaks = c("cyan", "orange", "black",
"green","red"),
labels = c(expression("Level L"[ref]),
expression(L[exc]),
expression("Level L+1"[ref]),
expression("L+1"["exc<B0>"]),
expression("L+1"["exc<B0+C0>"])),
guide = "legend") +
xlab("Excitation at Level L [0-1]") +
ylab("Prob. at Level L+1 for specific event") + theme_bw()
print(p)
Integration of information
#
defp = t(data.frame(
"01"=c('Failure','Failure','Failure','Failure','Failure'),
"02"=c('Failure','Failure','Failure','Correct','Failure'),
"03"=c('Failure','Failure','Correct','Failure','Failure'),
"04"=c('Failure','Failure','Correct','Correct','Failure'),
"05"=c('Failure','Correct','Failure','Failure','Failure'),
"06"=c('Failure','Correct','Failure','Correct','Failure'),
"07"=c('Failure','Correct','Correct','Failure','Failure'),
"08"=c('Failure','Correct','Correct','Correct','Failure'),
"09"=c('Correct','Failure','Failure','Failure','Failure'),
"10"=c('Correct','Failure','Failure','Correct','Failure'),
"11"=c('Correct','Failure','Correct','Failure','Failure'),
"12"=c('Correct','Failure','Correct','Correct','Failure'),
"13"=c('Correct','Correct','Failure','Failure','Failure'),
"14"=c('Correct','Correct','Failure','Correct','Failure'),
"15"=c('Correct','Correct','Correct','Failure','Failure'),
"16"=c('Correct','Correct','Correct','Correct','Failure'),
stringsAsFactors = FALSE
))
colnames(defp) = c('A2','B0','B2','C0','C2')#
sim = function(lowp,t){
prl =c()
prlp1=c()
prlp2=c()
prlp3=c()
prlp4=c()
prlp5=c()
prp1 =c()
prp2 =c()
#
for (i in 1:length(lowp)){
resl = layerl(t0,lowp[i])
p1 = resl[1]
p2 = resl[2]
prp1[i]=p1
prp2[i]=p2
# a b c d e f g h i b0 c0
bnlp1 = bnl_lp2(p1, t, t, t, t, t, t, t, t,p1,p1)[['net']]
bnlp2 = bnl_lp2(t , t, t, t, t, t, t, t, t,p2,p1)[['net']]
bnlp3 = bnl_lp2(p2,p2, t, t, t, t, t, t, t,p1,p2)[['net']]
bnlp4 = bnl_lp2(p2, t, t, t, t, t, t, t, t,p2,p2)[['net']]
bnlp5 = bnl_lp2(p2,p2, t,p2, t,p2, t,p2, t,p2,p2)[['net']]
bnlp1c = compile(bnlp1)
bnlp2c = compile(bnlp2)
bnlp3c = compile(bnlp3)
bnlp4c = compile(bnlp4)
bnlp5c = compile(bnlp5)
#
nodes = c('A2','B0','B2','C0','C2')
prlp1[i]= 0.
prlp2[i]= 0.
prlp3[i]= 0.
prlp4[i]= 0.
prlp5[i]= 0.
for (j in 1:nrow(defp)){
pp = 0.
def = defp[j,]
sev1 = setEvidence(object=bnlp1c,nodes=nodes,states=def)
pp = pEvidence(sev1)
prlp1[i] = prlp1[i] + pp
pp=0
sev2 = setEvidence(object=bnlp2c,nodes=nodes,states=def)
pp = pEvidence(sev2)
prlp2[i]= prlp2[i] + pp
pp=0
sev3 = setEvidence(object=bnlp3c,nodes=nodes,states=def)
pp = pEvidence(sev3)
prlp3[i]= prlp3[i] + pp
pp=0
sev4 = setEvidence(object=bnlp4c,nodes=nodes,states=def)
pp = pEvidence(sev4)
prlp4[i]= prlp4[i] + pp
pp=0
sev5 = setEvidence(object=bnlp5c,nodes=nodes,states=def)
pp = pEvidence(sev5)
prlp5[i]= prlp5[i] + pp
}
}
#
#
df = data.frame(x=lowp,base=prlp1,ex_b0=prlp2,ex_c0=prlp3,
ex_two=prlp4,ex_all=prlp5,pr1=prp1,pr2=prp2)
return(df)
}
#
drawpr = function(lowp,t){
res=sim(lowp,t)
p = ggplot(res, aes(x=x)) +
geom_line(aes(y = pr1, color = "cyan"), linetype="twodash") +
geom_line(aes(y = pr2, color = "orange"), linetype="twodash") +
geom_line(aes(y = base, color = "black")) +
geom_line(aes(y = ex_b0, color="green")) +
geom_line(aes(y = ex_two), color="darkblue") +
geom_line(aes(y = ex_all, color="red")) +
scale_y_continuous(trans='log10',
labels = trans_format("log10", math_format(10^.x))) +
theme(legend.position="top") +
scale_color_identity(name = paste0("P(C=F| t=",t,")"),
breaks = c("cyan", "orange", "black",
"green","darkblue","red"),
labels = c(expression("Level L"[ref]),
expression(L[exc]),
expression("Level L+1"[ref]),
expression("L+1"["exc<B0>"]),
expression("L+1"["exc<A+B0+C0>"]),
expression("L+1"["exc<AllExt>"])),
guide = "legend") +
xlab("Excitation at Level L [0-1]") +
ylab("Prob. at Level L+1 for specific event") + theme_bw()
return(p)
}
##
t0 =0.01
lt =c(0.05,0.1,0.15,0.2,0.25,0.3,0.4,0.5)
lowp =c(0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,0.99)
res=sim(lowp,t0)
for (t in lt){
p=drawpr(lowp,t)
plot(p)
}
