update scoring system before tests

R/MBR.R

@@ -187,9 +187,10 @@ if(FALSE){
                   }
                 }
 
-}
+            }
+            }
                 ############################## Acceptance probability
-                n.edges.tilde <- sum(dag.G.1)
+                n.edges.tilde <- sum(dag.MBR)
 
                 #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))
@@ -199,7 +200,7 @@ if(FALSE){
 
                 A <- min(exp(( s.proposed - s.current) * (n.edges/n.edges.tilde) ), 1)
 
-                # if(is.nan(score.A)){score.A <- 0}
+                if(is.nan(A)){A <- 0}
                 # if((score.A)<0){score.A <- 0}
                 # A <- min(1, score.A)
 
@@ -211,7 +212,7 @@ if(FALSE){
                   score <- score.dag(dag = dag.MBR,bsc.score = score.cache,sc = sc)
 
                 }
-            }
+
     if (is.null(A)) {
         A <- 0
     }

R/mc3.R

@@ -153,7 +153,7 @@ mc3 <- function(n.var, dag.tmp, retain, ban, max.parents, sc,sc.scaled, score.ca
         score <- score.G
 
         #if (!is.null(dag.gprime) && rbinom(n = 1, size = 1, prob = alpha) == 1) {
-        if (!is.null(dag.gprime) && runif(1)<=alpha) {
+        if (!is.null(dag.gprime) && runif(1)<alpha) {
             dag.tmp <- dag.gprime
             score <- score.Gprime
             rejection <- 0

R/mcmc.R

@@ -5,13 +5,13 @@
 ##-------------------------------------------------------------------------
 
 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, heating = 0.5,
+    1000), seed = 42, verbose = FALSE, start.dag = NULL, prior.dag = NULL, prior.lambda = NULL, prob.rev = 0.05, prob.mbr = 0.05, heating = 0,
     prior.choice = 2) {
 
     #################################################### Tests
 
     .tests.mcmcabn(score.cache, data.dists, max.parents, mcmc.scheme, seed, verbose, start.dag, prior.dag, prior.lambda,
-        prob.rev, prob.mbr, prior.choice)
+        prob.rev, prob.mbr, prior.choice,heating)
 
     ## end of tests
 
@@ -117,9 +117,10 @@ mcmcabn <- function(score.cache = NULL, score = "mlik", data.dists = NULL, max.p
 
     sc.scaled <- sc
     tmp.fact <- abs(max(sc[,ncol(sc)]))
-    sc.scaled[,ncol(sc)] <- sc[,ncol(sc)] *(tmp.fact - heating * tmp.fact)/tmp.fact
+    if(heating>0){
+      sc.scaled[,ncol(sc)] <- sc[,ncol(sc)]/(heating * tmp.fact)
+    }
 
-    #print(sc.scaled)
 
     ## scoring init
     score.init <- score.dag(dag.tmp,score.cache,sc)

R/mcmcabn-utiles.R

@@ -5,7 +5,7 @@
 ##-------------------------------------------------------------------------
 
 .tests.mcmcabn <- function(score.cache, data.dists, max.parents, mcmc.scheme, seed, verbose, start.dag, prior.dag, prior.lambda,
-    prob.rev, prob.mbr, prior.choice) {
+    prob.rev, prob.mbr, prior.choice,heating) {
     # start tests
     if (is.null(score.cache))
         stop("A cache of score should be provided. You can produce it using the R package abn.")
@@ -54,8 +54,13 @@
 
     if (is.matrix(prior.dag) && dim(prior.dag)[1] != length(data.dists))
         stop("prior.dag should be a squared matrix with dimension equal to the number of variables.")
+
+  if(heating<0)
+    stop("heating parameter shoud be between zero and one. Zero corresponds to a lower probability of accepting a lower score proposal. One corresponds to a larger probability of accepting a lowering of the score (values higher than one are possible but very unlikely to produce meaningful (or helpful) results.).")
+
 }
 
+
 ##-------------------------------------------------------------------------
 ## Internal function that call multiple times strsplit() and remove space
 ##-------------------------------------------------------------------------

R/plot-mcmcabn.R

@@ -3,7 +3,7 @@
 plot.mcmcabn <- function(x, max.score = FALSE, ...) {
 
     # utils::globalVariables(c('X','method' ,'scores')) utils::globalVariables(c('.', '%>%'))
-    dta <- data.frame(x[-1])
+    dta <- data.frame(x[2:4])
     dta$X <- 1:length(x$scores)
     max.score. <- max(x$scores)