Commit 5ea6ba2e authored by Gilles Kratzer's avatar Gilles Kratzer
Browse files

update scoring system

parent 33dbbeb5
Pipeline #2437 passed with stage
in 2 seconds
Package: mcmcabn
Title: Flexible Implementation of a Structural MCMC Sampler for DAGs
Version: 0.2
Version: 0.3
Authors@R: c(person("Gilles", "Kratzer", role = c("aut", "cre"),
email = "gilles.kratzer@math.uzh.ch",
comment = c(ORCID = "0000-0002-5929-8935")),
......
......@@ -14,9 +14,8 @@ print.mcmcabn
importFrom(cowplot, ggdraw, axis_canvas, insert_yaxis_grob)
importFrom(graphics, plot)
importFrom(stats, as.formula, dist, quantile, sd)
importFrom(matrixStats, logSumExp)
#methods for class mcmcabn
S3method(plot, mcmcabn)
S3method(summary, mcmcabn)
S3method(print, mcmcabn)
\ No newline at end of file
S3method(print, mcmcabn)
......@@ -5,3 +5,7 @@
## mcmcabn 0.2:
* patch to make compatible with abn 2.0
## mcmcabn 0.3:
* update mcmcabn to make it compatible with constraints imported from the cache of scores
MBR <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score, verbose) {
MBR <- function(n.var, dag.tmp, retain, ban, max.parents, sc,sc.scaled, score.cache, score, verbose) {
## scores
score.G <- 0
......@@ -47,7 +47,7 @@ MBR <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
new.parent.i <- new.parent.i[,-length(new.parent.i)]
z.i.G.0 <- (logSumExp(sc.tmp[, ncol(sc.tmp)]))
#z.i.G.0 <- (logSumExp(sc.tmp[, ncol(sc.tmp)]))
dag.G.1 <- dag.G.0
......@@ -81,8 +81,8 @@ MBR <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
}
new.parent.j <- new.parent.j[-length(new.parent.j)]
score.G[j] <- logSumExp(sc.tmp[, ncol(sc.tmp)])
score.G[j] <- sum(sc.tmp[, ncol(sc.tmp)])
#score.G[j] <- logSumExp(sc.tmp[, ncol(sc.tmp)])
}
}
......@@ -90,7 +90,7 @@ MBR <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
dag.MBR <- dag.G.1
############################ complementary inverse move
if(FALSE){
# parent set
parent.set.G.1 <- dag.G.1[selected.node, ]
......@@ -187,17 +187,24 @@ MBR <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
}
}
}
############################## Acceptance probability
n.edges.tilde <- sum(dag.G.1)
score.A <- ( z.i.G.0 / z.star.i.G.prime.0 * sum(is.finite(score.G)) / sum(is.finite(score.G.prime)))
#score.A <- n.edges/n.edges.tilde * ((z.star.x.i.M.dot) / (z.star.x.j.M.dot) * (z.x.j.M.cross) / (z.x.i.M.tilde.cross))
if(is.nan(score.A)){score.A <- 0}
if((score.A)<0){score.A <- 0}
A <- min(1, score.A)
#score.A <- ( z.i.G.0 / z.star.i.G.prime.0 * sum(is.finite(score.G)) / sum(is.finite(score.G.prime)))
#score.A <- n.edges/n.edges.tilde * ((z.star.x.i.M.dot) / (z.star.x.j.M.dot) * (z.x.j.M.cross) / (z.x.i.M.tilde.cross))
s.proposed <- score.dag(dag.MBR,score.cache,sc.scaled)
s.current <- score.dag(dag.tmp,score.cache,sc.scaled)
A <- min(exp(( s.proposed - s.current) * (n.edges/n.edges.tilde) ), 1)
# if(is.nan(score.A)){score.A <- 0}
# if((score.A)<0){score.A <- 0}
# A <- min(1, score.A)
if (rbinom(n = 1, size = 1, prob = A) == 1) {
if (runif(1)<A) {
rejection <- 0
dag.tmp <- dag.MBR
......@@ -208,6 +215,6 @@ MBR <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
if (is.null(A)) {
A <- 0
}
############################## Return ##score.MBR <- score
############################## Return
return(list(dag.tmp = dag.tmp, score = score, alpha = A, rejection = rejection))
} #EOF
REV <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score, verbose) {
REV <- function(n.var, dag.tmp, retain, ban, max.parents, sc,sc.scaled, score.cache, score, verbose) {
rejection <- 1
A <- 0
......@@ -9,10 +9,10 @@ REV <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
# randomly select one
if (sum(dag.tmp) != 0){
if (sum(dag.tmp*(1-retain)*(1-t(ban))) != 0){
#i->j
selected.edge <- which(x = ((dag.tmp*(1-retain)*(1-t(ban))) == 1), arr.ind = TRUE)[sample(x = 1:n.edges, size = 1), ,drop=FALSE]
selected.edge <- which(x = ((dag.tmp*(1-retain)*(1-t(ban))) == 1), arr.ind = TRUE)[sample(x = 1:sum(dag.tmp*(1-retain)*(1-t(ban))), size = 1), ,drop=FALSE]
# store current parent set (j)
parent.set <- dag.tmp[selected.edge[1], ,drop=FALSE]
......@@ -42,7 +42,8 @@ REV <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
sum(sc.tmp[, ncol(sc.tmp)]))), ,drop=FALSE]
new.parent.i <- new.parent.i[-length(new.parent.i)]
# store partition function
z.star.x.i.M.dot <- (logSumExp(sc.tmp[, ncol(sc.tmp)]))
#z.star.x.i.M.dot <- (logSumExp(sc.tmp[, ncol(sc.tmp)]))
#z.star.x.i.M.dot <- (sum(sc.tmp[, ncol(sc.tmp)]))
# M direct sum
dag.M.cross <- dag.M.dot
......@@ -67,7 +68,8 @@ REV <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
new.parent.j <- new.parent.j[-length(new.parent.j)]
# store partition function
z.x.j.M.cross <- (logSumExp(sc.tmp[, ncol(sc.tmp)]))
#z.x.j.M.cross <- (logSumExp(sc.tmp[, ncol(sc.tmp)]))
#z.x.j.M.cross <- (sum(sc.tmp[, ncol(sc.tmp)]))
# M tilde
......@@ -78,12 +80,13 @@ REV <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
############################## computing acceptance probability ##############################################
if(FALSE){
# score node j that do not include i
sc.tmp <- sc[score.cache$children == selected.edge[1], ,drop=FALSE]
sc.tmp <- sc.tmp[sc.tmp[, selected.edge[2]] == 1, ,drop=FALSE]
sc.tmp <- sc.tmp[rowSums(sc.tmp[, descendents.M.dot.j,drop=FALSE] == 0) == length(descendents.M.dot.j), ,drop=FALSE]
z.star.x.j.M.dot <- (logSumExp(sc.tmp[, ncol(sc.tmp)]))
#z.star.x.j.M.dot <- (logSumExp(sc.tmp[, ncol(sc.tmp)]))
dag.M.tilde.cross <- dag.M.dot
dag.M.tilde.cross[selected.edge[1],] <- parent.set
......@@ -97,26 +100,27 @@ REV <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
,drop=FALSE]
}
z.x.i.M.tilde.cross <- (logSumExp(sc.tmp[, ncol(sc.tmp)]))
}
############################## Acceptance probability
score.A <- n.edges/n.edges.tilde * ((z.star.x.i.M.dot) / (z.star.x.j.M.dot) * (z.x.j.M.cross) / (z.x.i.M.tilde.cross))
if(is.nan(score.A)){score.A <- 0}
if((score.A)<0){score.A <- 0}
A <- min(1, score.A)
#score.A <- min(exp((z.star.x.i.M.dot / (z.star.x.j.M.dot) * (z.x.j.M.cross) / (z.x.i.M.tilde.cross))*n.edges/n.edges.tilde),1)
s.proposed <- score.dag(dag.M.tilde,score.cache,sc.scaled)
s.current <- score.dag(dag.tmp,score.cache,sc.scaled)
A <- min(exp(( s.proposed - s.current) * (n.edges/n.edges.tilde) ), 1)
#if(is.nan(score.A)){score.A <- 0}
if (rbinom(n = 1, size = 1, prob = A) == 1) {
#if((score.A)<0){score.A <- 0}
#A <- min(1, score.A)
#if (rbinom(n = 1, size = 1, prob = A) == 1) {
if (runif(1)<A){
rejection <- 0
dag.tmp <- dag.M.tilde
score.M.tilde <- 0
for (a in 1:n.var) {
sc.tmp <- sc[score.cache$children == a, , drop=FALSE]
score.M.tilde <- sum(min(sc.tmp[(apply(sc.tmp, 1, function(x) identical(unname(x[1:n.var]), unname(dag.tmp[a,
])))), n.var + 1]), score.M.tilde)
}
score <- score.M.tilde
score <- score.dag(dag.M.tilde,score.cache,sc)
}
}#eoif
......
mc3 <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score, prior.choice, prior.lambda, prior.dag, verbose) {
mc3 <- function(n.var, dag.tmp, retain, ban, max.parents, sc,sc.scaled, score.cache, score, prior.choice, prior.lambda, prior.dag, verbose, scaling.factor) {
## construction of neighbours list
neighbours.list <- NULL
......@@ -140,15 +140,20 @@ mc3 <- function(n.var, dag.tmp, retain, ban, max.parents, sc, score.cache, score
])))), n.var + 1]), score.Gprime)
}
# alpha <- min(exp((score.Gprime + log(prior.Gprime) + log(n.Gprime)) - ((score.G + log(prior.G)) + log(n.G))),1)
alpha <- min(exp(score.Gprime - score.G) * (n.G/n.Gprime) * (prior.Gprime/prior.G), 1)
if(!is.numeric(alpha) | is.nan(alpha)) alpha <- 0
score.Gprime.scaled <- score.dag(dag.gprime,score.cache,sc.scaled)
score.G.scaled <- score.dag(dag.tmp,score.cache,sc.scaled)
alpha <- min(exp(( score.Gprime.scaled - score.G.scaled) * (n.G/n.Gprime) * (prior.Gprime/prior.G)), 1)
#alpha <- min(exp(( score.Gprime.scaled - score.G.scaled) * (n.G/n.Gprime) * (prior.Gprime/prior.G)), 1)
#alpha <- min(exp(( score.Gprime - score.G) * (n.G/n.Gprime) * (prior.Gprime/prior.G)), 1)
if(!is.numeric(alpha) | is.nan(alpha)) alpha <- 0
score <- score.G
if (!is.null(dag.gprime) && rbinom(n = 1, size = 1, prob = alpha) == 1) {
#if (!is.null(dag.gprime) && rbinom(n = 1, size = 1, prob = alpha) == 1) {
if (!is.null(dag.gprime) && runif(1)<=alpha) {
dag.tmp <- dag.gprime
score <- score.Gprime
rejection <- 0
......
......@@ -5,7 +5,7 @@
##-------------------------------------------------------------------------
mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.parents = 1, mcmc.scheme = c(100, 1000,
1000), seed = 42, verbose = FALSE, start.dag = NULL, prior.dag = NULL, prior.lambda = NULL, prob.rev = 0.05, prob.mbr = 0.05,
1000), seed = 42, verbose = FALSE, start.dag = NULL, prior.dag = NULL, prior.lambda = NULL, prob.rev = 0.05, prob.mbr = 0.05, heating = 0.5,
prior.choice = 2) {
#################################################### Tests
......@@ -32,6 +32,8 @@ mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.p
prob.mc3 <- 1 - (prob.rev + prob.mbr)
## output
out.dags <- array(data = NA, dim = c(n.var, n.var, mcmc.scheme[1] + 1))
out.scores <- rep(NA, length = mcmc.scheme[1] + 1)
......@@ -96,10 +98,12 @@ mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.p
colnames(dag.tmp) <- rownames(dag.tmp) <- 1:n.var
}
## scores
if (score %in% c("bic", "aic", "mdl")) {
tmp <- score.cache$node.defn[, as.numeric(colnames(dag.tmp))]
colnames(tmp) <- as.numeric(colnames(dag.tmp))
sc <- cbind(tmp, -score.cache[[score]])
}
......@@ -107,17 +111,24 @@ mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.p
tmp <- score.cache$node.defn[, as.numeric(colnames(dag.tmp))]
colnames(tmp) <- as.numeric(colnames(dag.tmp))
sc <- cbind(tmp, score.cache[[score]])
}
## scoring init
score.init <- 0
for (a in 1:n.var) {
sc.tmp <- sc[score.cache$children == as.numeric(colnames(dag.tmp)[a]), ,drop = FALSE]
score.init <- sum(min(sc.tmp[which(apply(sc.tmp, 1, function(x) identical(unname(x[1:n.var]), unname(dag.tmp[a,
])))), n.var + 1]), score.init)
sc.scaled <- sc
tmp.fact <- abs(max(sc[,ncol(sc)]))
sc.scaled[,ncol(sc)] <- sc[,ncol(sc)] *(tmp.fact - heating * tmp.fact)/tmp.fact
}
#print(sc.scaled)
## scoring init
score.init <- score.dag(dag.tmp,score.cache,sc)
# for (a in 1:n.var) {
# sc.tmp <- sc[score.cache$children == as.numeric(colnames(dag.tmp)[a]), ,drop = FALSE]
# score.init <- sum(min(sc.tmp[which(apply(sc.tmp, 1, function(x) identical(unname(x[1:n.var]), unname(dag.tmp[a,
# ])))), n.var + 1]), score.init)
#
# }
## out
......@@ -144,7 +155,7 @@ mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.p
method.choice <- sample(x = c("MC3", "REV", "MBR"), size = 1, prob = c(prob.mc3, prob.rev, prob.mbr))
switch(method.choice, MC3 = {
out <- mc3(n.var, (dag.tmp), retain, ban, max.parents, sc, score.cache, score, prior.choice, prior.lambda, prior.dag,
out <- mc3(n.var, (dag.tmp), retain, ban, max.parents, sc,sc.scaled, score.cache, score, prior.choice, prior.lambda, prior.dag,
verbose)
dag.tmp <- out$dag.tmp
score <- out$score
......@@ -152,14 +163,14 @@ mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.p
if (verbose) {
print("REV move")
}
out <- REV(n.var, (dag.tmp),retain, ban, max.parents, sc, score.cache, score, verbose)
out <- REV(n.var, (dag.tmp),retain, ban, max.parents, sc,sc.scaled, score.cache, score, verbose)
dag.tmp <- out$dag.tmp
score <- out$score
}, MBR = {
if (verbose) {
print("MBR move")
}
out <- MBR(n.var, dag.tmp,retain, ban, max.parents, sc, score.cache, score, verbose)
out <- MBR(n.var, dag.tmp,retain, ban, max.parents, sc,sc.scaled, score.cache, score, verbose)
dag.tmp <- out$dag.tmp
score <- out$score
})
......@@ -185,7 +196,7 @@ mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.p
method.choice <- sample(x = c("MC3", "REV", "MBR"), size = 1, prob = c(prob.mc3, prob.rev, prob.mbr))
switch(method.choice, MC3 = {
out <- mc3(n.var, (dag.tmp), retain, ban, max.parents, sc, score.cache, score, prior.choice, prior.lambda, prior.dag,
out <- mc3(n.var, (dag.tmp), retain, ban, max.parents, sc,sc.scaled, score.cache, score, prior.choice, prior.lambda, prior.dag,
verbose)
dag.tmp <- out$dag.tmp
score <- out$score
......@@ -195,7 +206,7 @@ mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.p
if (verbose) {
print("REV move")
}
out <- REV(n.var, (dag.tmp),retain, ban, max.parents, sc, score.cache, score, verbose)
out <- REV(n.var, (dag.tmp),retain, ban, max.parents, sc,sc.scaled, score.cache, score, verbose)
dag.tmp <- out$dag.tmp
score <- out$score
alpha <- out$alpha
......@@ -204,7 +215,7 @@ mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.p
if (verbose) {
print("MBR move")
}
out <- MBR(n.var, (dag.tmp), retain, ban, max.parents, sc, score.cache, score, verbose)
out <- MBR(n.var, (dag.tmp), retain, ban, max.parents, sc,sc.scaled, score.cache, score, verbose)
dag.tmp <- out$dag.tmp
score <- out$score
alpha <- out$alpha
......
......@@ -53,6 +53,8 @@ ___
- 18/02/2019 - new pre-print [Is a single unique Bayesian network enough to accurately represent your data?](https://arxiv.org/pdf/1902.06641.pdf) on arXiv
- 01/07/2019 - mcmcabn 0.2 available on CRAN
____
**`mcmcabn` is developed and maintained by [Gilles Kratzer](https://gilleskratzer.netlify.com/) and [Prof. Dr. Reinhard Furrer](https://user.math.uzh.ch/furrer/) from
......
......@@ -60,7 +60,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="../index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......
......@@ -30,7 +30,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="../index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......@@ -150,7 +150,7 @@
<a class="sourceLine" id="cb3-19" data-line-number="19"> <span class="dt">Dyspnea =</span> <span class="st">"binomial"</span>)</a>
<a class="sourceLine" id="cb3-20" data-line-number="20"></a>
<a class="sourceLine" id="cb3-21" data-line-number="21"><span class="co">#plot the BN as described in the paper</span></a>
<a class="sourceLine" id="cb3-22" data-line-number="22">dag &lt;-<span class="st"> </span><span class="kw"><a href="https://malco.io/ggdag//reference/dagify.html">dagify</a></span>(Tuberculosis<span class="op">~</span>Asia,</a>
<a class="sourceLine" id="cb3-22" data-line-number="22">dag &lt;-<span class="st"> </span><span class="kw"><a href="https://www.rdocumentation.org/packages/ggdag/topics/dagify">dagify</a></span>(Tuberculosis<span class="op">~</span>Asia,</a>
<a class="sourceLine" id="cb3-23" data-line-number="23"> Either<span class="op">~</span>Tuberculosis,</a>
<a class="sourceLine" id="cb3-24" data-line-number="24"> XRay<span class="op">~</span>Either,</a>
<a class="sourceLine" id="cb3-25" data-line-number="25"> Dyspnea<span class="op">~</span>Either,</a>
......@@ -159,7 +159,7 @@
<a class="sourceLine" id="cb3-28" data-line-number="28"> Either<span class="op">~</span>LungCancer,</a>
<a class="sourceLine" id="cb3-29" data-line-number="29"> Dyspnea<span class="op">~</span>Bronchitis)</a>
<a class="sourceLine" id="cb3-30" data-line-number="30"></a>
<a class="sourceLine" id="cb3-31" data-line-number="31"><span class="kw"><a href="https://malco.io/ggdag//reference/ggdag_classic.html">ggdag_classic</a></span>(dag, <span class="dt">size =</span> <span class="dv">6</span>) <span class="op">+</span><span class="st"> </span><span class="kw"><a href="https://malco.io/ggdag//reference/theme_dag_blank.html">theme_dag_blank</a></span>()</a></code></pre></div>
<a class="sourceLine" id="cb3-31" data-line-number="31"><span class="kw"><a href="https://www.rdocumentation.org/packages/ggdag/topics/ggdag_classic">ggdag_classic</a></span>(dag, <span class="dt">size =</span> <span class="dv">6</span>) <span class="op">+</span><span class="st"> </span><span class="kw"><a href="https://www.rdocumentation.org/packages/ggdag/topics/theme_dag_blank">theme_dag_blank</a></span>()</a></code></pre></div>
<p><img src="mcmcabn_files/figure-html/unnamed-chunk-3-1.png" width="672" style="display: block; margin: auto;"></p>
<p>For the first search, it is advised to limit the number of parents per node to 2, which is the maximum number of parent of the network.</p>
<div class="sourceCode" id="cb4"><pre class="sourceCode r"><code class="sourceCode r"><a class="sourceLine" id="cb4-1" data-line-number="1"><span class="co">#loglikelihood scores</span></a>
......
......@@ -60,7 +60,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......
......@@ -30,7 +30,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......@@ -125,6 +125,7 @@
<ul>
<li><p>08/03/2019 - mcmcabn is available on CRAN (v 0.1)</p></li>
<li><p>18/02/2019 - new pre-print <a href="https://arxiv.org/pdf/1902.06641.pdf">Is a single unique Bayesian network enough to accurately represent your data?</a> on arXiv</p></li>
<li><p>01/07/2019 - mcmcabn 0.2 available on CRAN</p></li>
</ul>
<hr>
<p><strong><code>mcmcabn</code> is developed and maintained by <a href="https://gilleskratzer.netlify.com/">Gilles Kratzer</a> and <a href="https://user.math.uzh.ch/furrer/">Prof. Dr. Reinhard Furrer</a> from <a href="https://www.math.uzh.ch/as/index.php?id=as">Applied Statistics Group</a> from the University of Zurich.</strong></p>
......
......@@ -60,7 +60,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="../index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......@@ -135,6 +135,13 @@
<ul>
<li>patch to make compatible with abn 2.0</li>
</ul>
</div>
<div id="mcmcabn-0-3" class="section level2">
<h2 class="hasAnchor">
<a href="#mcmcabn-0-3" class="anchor"></a>mcmcabn 0.3:</h2>
<ul>
<li>update mcmcabn to make it compatible with constraints imported from the cache of scores</li>
</ul>
</div>
</div>
......@@ -144,6 +151,7 @@
<ul class="nav nav-pills nav-stacked">
<li><a href="#mcmcabn-0-1">0.1</a></li>
<li><a href="#mcmcabn-0-2">0.2</a></li>
<li><a href="#mcmcabn-0-3">0.3</a></li>
</ul>
</div>
</div>
......
......@@ -63,7 +63,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="../index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......
......@@ -63,7 +63,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="../index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......
......@@ -60,7 +60,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="../index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......
......@@ -63,7 +63,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="../index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......
......@@ -63,7 +63,7 @@
</button>
<span class="navbar-brand">
<a class="navbar-link" href="../index.html">mcmcabn</a>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.2</span>
<span class="version label label-default" data-toggle="tooltip" data-placement="bottom" title="Released version">0.3</span>
</span>
</div>
......@@ -263,28 +263,28 @@
#&gt; Number of MCMC steps: 100
#&gt; Thinning: 0
#&gt;
#&gt; Maximum score: -11252.53
#&gt; Empirical mean: -12577.33
#&gt; Empirical standard deviation: 1420.952
#&gt; Maximum score: -11511.4
#&gt; Empirical mean: -13093.47
#&gt; Empirical standard deviation: 1155.325
#&gt; Quantiles of the posterior network score:
#&gt; 0.025 0.25 0.5 0.75 0.975
#&gt; BN score -15127.73 -13458.25 -12044.84 -11338.42 -11254.25
#&gt; 0.025 0.25 0.5 0.75 0.975
#&gt; BN score -15129.49 -13795.44 -12702.89 -12576.32 -11511.4
#&gt;
#&gt;
#&gt; Global acceptance rate: 0.3366337
#&gt; Global acceptance rate: 0.1980198
#&gt; Accepted Rejected
#&gt; MBR 2 2
#&gt; MC3 32 63
#&gt; REV 0 2
#&gt; MBR 1 0
#&gt; MC3 18 79
#&gt; REV 1 2
#&gt;
#&gt;
#&gt; Sample size adjusted for autocorrelation: 1.429643
#&gt; Sample size adjusted for autocorrelation: 1.873662
#&gt;
#&gt; Autocorrelations by lag:
#&gt; 0 1 2 3 4 5 6 7
#&gt; acf 1 0.9718074 0.9436429 0.9153732 0.8849037 0.8544765 0.823888 0.7931788
#&gt; 0 1 2 3 4 5 6 7
#&gt; acf 1 0.9632124 0.926465 0.8894539 0.8524517 0.8154436 0.7784273 0.7414108
#&gt; 8 9 10
#&gt; acf 0.7625209 0.7275857 0.6933869</div><div class='input'>
#&gt; acf 0.7044402 0.6674692 0.6278747</div><div class='input'>
<span class='co'># Uniquelly with MC3 moves</span>
<span class='no'>mcmc.out.asia.small</span> <span class='kw'>&lt;-</span> <span class='fu'>mcmcabn</span>(<span class='kw'>score.cache</span> <span class='kw'>=</span> <span class='no'>bsc.compute.asia</span>,
......@@ -304,26 +304,26 @@
#&gt; Number of MCMC steps: 100
#&gt; Thinning: 0
#&gt;
#&gt; Maximum score: -13046.48
#&gt; Empirical mean: -13667.21
#&gt; Empirical standard deviation: 777.3341
#&gt; Maximum score: -11643.48
#&gt; Empirical mean: -14129.55
#&gt; Empirical standard deviation: 947.4905
#&gt; Quantiles of the posterior network score:
#&gt; 0.025 0.25 0.5 0.75 0.975
#&gt; BN score -15147.68 -14640.66 -13136.8 -13136.28 -13046.48
#&gt; 0.025 0.25 0.5 0.75 0.975
#&gt; BN score -15147.68 -14744 -14689.33 -12829.51 -12828.45
#&gt;
#&gt;
#&gt; Global acceptance rate: 0.3960396
#&gt; Global acceptance rate: 0.1287129
#&gt; Accepted Rejected
#&gt; MC3 40 61
#&gt; MC3 13 88
#&gt;
#&gt;
#&gt; Sample size adjusted for autocorrelation: 1.732075
#&gt; Sample size adjusted for autocorrelation: 2.050773
#&gt;
#&gt; Autocorrelations by lag:
#&gt; 0 1 2 3 4 5 6 7
#&gt; acf 1 0.965945 0.9316219 0.8967652 0.8616907 0.8313675 0.8032237 0.7761743
#&gt; 0 1 2 3 4 5 6 7
#&gt; acf 1 0.9428316 0.8772444 0.8444724 0.8117004 0.7789152 0.748184 0.7182643
#&gt; 8 9 10
#&gt; acf 0.7491203 0.7197729 0.6902472</div><div class='input'>
#&gt; acf 0.6883445 0.6573608 0.6254991</div><div class='input'>
<span class='co'>#let us define a starting DAG</span>
<span class='no'>startDag</span> <span class='kw'>&lt;-</span> <span class='fu'><a href='https://www.rdocumentation.org/packages/base/topics/matrix'>matrix</a></span>(<span class='kw'>data</span> <span class='kw'>=</span> <span class='fu'><a href='https://www.rdocumentation.org/packages/base/topics/c'>c</a></span>(<span class='fl'>0</span>, <span class='fl'>0</span>, <span class='fl'>0</span>, <span class='fl'>1</span>, <span class='fl'>0</span>, <span class='fl'>0</span>, <span class='fl'>0</span>, <span class='fl'>0</span>,
<span class='fl'>0</span>, <span class='fl'>0</span>, <span class='fl'>1</span>, <span class='fl'>0</span>, <span class='fl'>0</span>, <span class='fl'>0</span>, <span class='fl'>0</span>, <span class='fl'>0</span>,
......@@ -354,26 +354,26 @@
#&gt; Number of MCMC steps: 100
#&gt; Thinning: 0
#&gt;
#&gt; Maximum score: -12351.24
#&gt; Empirical mean: -12726.16
#&gt; Empirical standard deviation: 574.7013
#&gt; Maximum score: -11654.16
#&gt; Empirical mean: -12049.47
#&gt; Empirical standard deviation: 706.0496
#&gt; Quantiles of the posterior network score:
#&gt; 0.025 0.25 0.5 0.75 0.975
#&gt; BN score -14160.15 -12492.22 -12461.74 -12449.72 -12351.25
#&gt; 0.025 0.25 0.5 0.75 0.975
#&gt; BN score -14114.03 -11980 -11654.16 -11654.16 -11654.16
#&gt;
#&gt;
#&gt; Global acceptance rate: 0.3960396
#&gt; Global acceptance rate: 0.1287129
#&gt; Accepted Rejected