net.1 = custom.fit(net, dist = list(x1=cpt1.x1, x2=cpt1.x2, x3=cpt1.x3, x4=cpt1.x4)) query=rep(NA,length(net.1)) names(query) <- names(net.1) query[2]='false' my.bn.inference(net.1,query) ## $logp.evi ## [1] 0 ## $logp.query ## [1] -0.1743534 ## $p ## [1] 0.84 ## $evidence ## [1] "" ## $query ## [1] "x2=false"