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, parents

Defining 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 + A->C

bnl_lp2 = function(a,b,c,d,e,f,g) {
  #
  # a pr(A2 = 'F')
  # b pr(B2 = 'F' | A2 = 'C')
  # c pr(B2 = 'F' | A2 = 'F')
  # d pr(C2 = 'F' | B2='F', A2='F')
  # e pr(C2 = 'F' | B2='F', A2='C')
  # f pr(C2 = 'F' | B2='C', A2='F')
  # g pr(C2 = 'F' | B2='C', A2='C')
  #
  fc <- c("Failure","Correct")
  Astp = cptable(~A2, values=c(a,1-a),levels=fc)
  Cstp = cptable(~C2, values=c(d,1-d),levels=fc)  
  Bstp = cptable(~B2+A2, values=c(c,1-c,b,1-b),levels=fc)
  Cstp = cptable(~C2+A2+B2, values=c(d,1-d,e,1-e,f,1-f,g,1-g),levels=fc)
  cpts = compileCPT(list(Astp,Bstp,Cstp))
  bnl.stp = grain(cpts)
  return(list('net'=bnl.stp, 'parms'=c(a,b,c,d,e,f,g)))
}

Defining the level “l” BN

t =0.01
t1=0.10
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)
def= c('Failure','')
nodes=c('A1','B1')
sev2 = setEvidence(object=bnl1c,nodes=nodes,states=def)
p2=pEvidence(sev2)
p2
## [1] 0.1
#

Failure in A1 => C2

#
lowp=c(0.01,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,0.99)
prb=c()
for (i in 1:length(lowp)){
  bnlp1 = bnl_lp2(lowp[i],p1,p1,lowp[i],p1,p1,p1)[['net']]
  bnlp1c= compile(bnlp1)
  nodes=c('A2','B2','D2')
  def  =c('Failure','','Failure')
  sev3 = setEvidence(object=bnlp1c,nodes=nodes,states=def)
  prb[i]=pEvidence(sev3)
}
plot(lowp,prb,type="b")

#
#