Merge branch 'issue-17' into develop

DESCRIPTION

@@ -1,7 +1,7 @@
 Package: EMJMCMC
 Type: Package
 Title: Flexible Bayesian Nonlinear Model Configuration
-Version: 1.4.4.9030
+Version: 1.4.4.9031
 Date: 2022-03-06
 Authors@R:
   c(

R/EMJMCMC2016-method-modejumping_mcmc.R

@@ -96,7 +96,8 @@ EMJMCMC2016$methods(
       j.a <- j.a + 1
       if (glob.model$print.freq > 0 && j %% glob.model$print.freq == 0) {
         cat(
-          formatC(j, width = 4L), "iterations completed up to now after",
+          formatC(j, width = 4L, format = "d"),
+          "iterations completed up to now after",
           formatC(delta.time, digits = 6L, flag = "-", format = "f"),
           "cpu minutes",
           "best MLIK found",

R/LogicRegr.R

@@ -42,19 +42,23 @@
 #' @export
 LogicRegr = function(
   formula, data, family = "Gaussian",prior = "J",report.level = 0.5, d = 20,
-  cmax = 5, kmax = 20, p.and = 0.9, p.not = 0.05, p.surv = 0.1,ncores = -1,
-  n.mods = 1000 ,advanced = list(presearch = T,locstop = F ,
-  estimator = estimate.logic.bern.tCCH,
-  estimator.args =  list(data = data.example,n = 1000, m = 50,r=1),
-  recalc_margin = 250, save.beta = FALSE, interact = TRUE,
-  relations = c("","lgx2","cos","sigmoid","tanh","atan","erf"),
-  relations.prob =c(0.4,0.0,0.0,0.0,0.0,0.0,0.0),
-  interact.param = list(
-    allow_offsprings=1,mutation_rate = 300,last.mutation = 5000,
-    max.tree.size = 1, Nvars.max = 100,p.allow.replace=0.9,p.allow.tree=0.2,
-    p.nor=0.2,p.and = 1),n.models = 10000,unique = TRUE,max.cpu = 4,
-    max.cpu.glob = 4,create.table = FALSE, create.hash = TRUE,
-    pseudo.paral = TRUE, burn.in = 50, outgraphs = FALSE, print.freq = 1000,
+  cmax = 5, kmax = 20, p.and = 0.9, p.not = 0.05, p.surv = 0.1, ncores = -1,
+  n.mods = 1000,
+  advanced = list(
+    presearch = TRUE,locstop = FALSE,
+    estimator = estimate.logic.bern.tCCH,
+    estimator.args =  list(data = data.example,n = 1000, m = 50,r=1),
+    recalc_margin = 250, save.beta = FALSE, interact = TRUE,
+    relations = c("","lgx2","cos","sigmoid","tanh","atan","erf"),
+    relations.prob =c(0.4,0.0,0.0,0.0,0.0,0.0,0.0),
+    interact.param = list(
+      allow_offsprings=1, mutation_rate = 300, last.mutation = 5000,
+      max.tree.size = 1, Nvars.max = 100, p.allow.replace=0.9, p.allow.tree=0.2,
+      p.nor=0.2, p.and = 1
+    ),
+    n.models = 10000, unique = TRUE, max.cpu = 4, max.cpu.glob = 4,
+    create.table = FALSE, create.hash = TRUE, pseudo.paral = TRUE, burn.in = 50,
+    outgraphs = FALSE, print.freq = 1000,
     advanced.param = list(
       max.N.glob=as.integer(10),
       min.N.glob=as.integer(5),
@@ -108,6 +112,10 @@ LogicRegr = function(
   if(ncores<1)
     ncores = parallel::detectCores()
 
-  return(pinferunemjmcmc(n.cores = ncores,report.level =  report.level, simplify = T, num.mod.best = n.mods, predict = F, runemjmcmc.params = advanced))
-
+  return(
+    pinferunemjmcmc(
+      n.cores = ncores, report.level = report.level, simplify = TRUE,
+      num.mod.best = n.mods, predict = FALSE, runemjmcmc.params = advanced
+    )
+  )
 }

tests/testthat/test-BLR-tutorial.R

@@ -8,16 +8,18 @@
 # omp_set_num_threads(1)
 # #***********************IMPORTANT******************************************************
 
-n_threads <- min(parallel::detectCores() - 1, 7L)
+n_threads <- 1L
 set.seed(040590)
 
 # construct a correlation matrix for M = 50 variables
 M <- 50
 m <- clusterGeneration::rcorrmatrix(M, alphad = 2.5)
 
 # simulate 1000 binary variables with this correlation matrix
-sample_size <- 100L
-X <- suppressWarnings(bindata::rmvbin(sample_size, margprob = rep(0.5, M), bincorr = m))
+sample_size <- 1000L
+X <- suppressWarnings(
+  bindata::rmvbin(sample_size, margprob = rep(0.5, M), bincorr = m)
+)
 
 # prepare the correlation matrix in the melted format
 melted_cormat <- reshape2::melt(cor(X))
@@ -41,28 +43,44 @@ formula1 <- as.formula(
 )
 
 # Bayesian logic regression with the robust-g-prior
-# FIXME: returns NULL objects even if ran at full power
 res4G <- LogicRegr(
   formula = formula1, data = df, family = "Gaussian", prior = "G",
   report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.1,
-  p.surv = 0.2,
-  ncores = n_threads,
-  n.mods = 10L
+  p.surv = 0.2, ncores = n_threads
 )
 
 # Bayesian logic regression with the Jeffreys prior
-# res4J <- LogicRegr(
-#   formula = formula1, data = df, family = "Gaussian", prior = "J",
-#   report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.1,
-#   p.surv = 0.2, ncores = n_threads
-# )
+res4J <- LogicRegr(
+  formula = formula1, data = df, family = "Gaussian", prior = "J",
+  report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.1,
+  p.surv = 0.2, ncores = n_threads
+)
+
+# NULLs are expected because predict = FALSE on LogicRegr
 
 # # print the results for the robust g-prior
-# print(base::rbind(c("expressions","probabilities"),res4G$feat.stat))
+test_that("Results for the G-prior are sensible", {
+  for (i in seq_len(nrow(res4G$allposteriors))) {
+    expect_gte(res4G$allposteriors[i, 2], 0)
+    expect_lte(res4G$allposteriors[i, 2], 1)
+  }
+  expect_gte(res4G$threads.stats[[1]]$post.populi, 0)
+  expect_gt(res4G$threads.stats[[1]]$cterm, 1000)
+  expect_equal(res4G$threads.stats[[1]]$preds, NULL)
+  expect_length(res4G, 4L)
+})
 
 # #print the results for the Jeffreys prior
-# print(base::rbind(c("expressions","probabilities"),res4J$feat.stat))
-
+test_that("Results for the Jeffrey's prior are sensible", {
+  for (i in seq_len(nrow(res4J$allposteriors))) {
+    expect_gte(res4J$allposteriors[i, 2], 0)
+    expect_lte(res4J$allposteriors[i, 2], 1)
+  }
+  expect_gte(res4J$threads.stats[[1]]$post.populi, 0)
+  expect_gt(res4J$threads.stats[[1]]$cterm, -1000)
+  expect_equal(res4J$threads.stats[[1]]$preds, NULL)
+  expect_length(res4J, 4L)
+})
 
 # # simulate Gaussian responses from a model with two-way interactions
 # # and an age effect which is an extension of S.4 of the paper