Skip to content

Commit

Permalink
Rewritten output as unit tests (#17)
Browse files Browse the repository at this point in the history
  • Loading branch information
wleoncio committed Feb 8, 2024
1 parent 64cdf7e commit 6de6170
Showing 1 changed file with 98 additions and 28 deletions.
126 changes: 98 additions & 28 deletions tests/testthat/test-inference-help.R
Original file line number Diff line number Diff line change
@@ -1,48 +1,118 @@
# make sure that you are using Mac Os or Linux (mclapply is currently not supported for Windows unless some mclapply hack function for Windows is preloaded in your R session)

#*********************** IMPORTANT******************************************************
# if a multithreaded backend openBLAS for matrix multiplications
# is installed on your machine, please force it to use 1 thread explicitly
# unless ncores in LogrRegr is reasonably small, in the latter case
# you might want to experiment with the combinations of blas_set_num_threads and ncores
library(RhpcBLASctl)
blas_set_num_threads(1)
omp_set_num_threads(1)
#*********************** IMPORTANT******************************************************

# simulate Gaussian responses

threads <- 1L
set.seed(040590)
X1 <- as.data.frame(array(data = rbinom(n = 50 * 1000, size = 1, prob = runif(n = 50 * 1000, 0, 1)), dim = c(1000, 50)))
Y1 <- rnorm(n = 1000, mean = 1 + 0.7 * (X1$V1 * X1$V4) + 0.8896846 * (X1$V8 * X1$V11) + 1.434573 * (X1$V5 * X1$V9), sd = 1)
X1 <- as.data.frame(
array(
data = rbinom(n = 50 * 1000, size = 1, prob = runif(n = 50 * 1000, 0, 1)),
dim = c(1000, 50)
)
)
Y1 <- rnorm(
n = 1000,
mean = 1 + 0.7 * (X1$V1 * X1$V4) + 0.8896846 * (X1$V8 * X1$V11) + 1.434573 * (X1$V5 * X1$V9),
sd = 1
)
X1$Y1 <- Y1

# specify the initial formula
formula1 <- as.formula(paste(colnames(X1)[51], "~ 1 +", paste0(colnames(X1)[-c(51)], collapse = "+")))
formula1 <- as.formula(
paste(colnames(X1)[51], "~ 1 +", paste0(colnames(X1)[-c(51)], collapse = "+"))
)
data.example <- as.data.frame(X1)

# run the inference with robust g prior
res4G <- LogicRegr(formula = formula1, data = data.example, family = "Gaussian", prior = "G", report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.01, p.surv = 0.2, ncores = 32)
print(res4G$feat.stat)
res4G <- suppressMessages(
LogicRegr(
formula = formula1, data = data.example, family = "Gaussian", prior = "G",
report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.01,
p.surv = 0.2, ncores = threads, print.freq = 0L
)
)
# run the inference with Jeffrey's prior
res4J <- LogicRegr(formula = formula1, data = data.example, family = "Gaussian", prior = "J", report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.01, p.surv = 0.2, ncores = 32)
print(res4J$feat.stat)

res4J <- suppressMessages(
LogicRegr(
formula = formula1, data = data.example, family = "Gaussian", prior = "J",
report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.01,
p.surv = 0.2, ncores = threads, print.freq = 0L
)
)

# change to Bernoulli responses
X1 <- as.data.frame(array(data = rbinom(n = 50 * 1000, size = 1, prob = 0.3), dim = c(1000, 50)))
X1 <- as.data.frame(
array(data = rbinom(n = 50 * 1000, size = 1, prob = 0.3), dim = c(1000, 50))
)
Y1 <- -0.7 + 1 * ((1 - X1$V1) * (X1$V4)) + 1 * (X1$V8 * X1$V11) + 1 * (X1$V5 * X1$V9)
X1$Y1 <- round(1.0 / (1.0 + exp(-Y1)))

# specify the initial formula
formula1 <- as.formula(paste(colnames(X1)[51], "~ 1 +", paste0(colnames(X1)[-c(51)], collapse = "+")))
formula1 <- as.formula(
paste(colnames(X1)[51], "~ 1 +", paste0(colnames(X1)[-c(51)], collapse = "+"))
)
data.example <- as.data.frame(X1)


# run the inference with robust g prior
res1G <- LogicRegr(formula = formula1, data = data.example, family = "Bernoulli", prior = "G", report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.2, p.surv = 0.2, ncores = 32)
print(res1G$feat.stat)
res1G <- suppressWarnings(
suppressMessages(
LogicRegr(
formula = formula1, data = data.example, family = "Bernoulli", prior = "G",
report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.2,
p.surv = 0.2, ncores = threads, print.freq = 0L
)
)
)

# run the inference with Jeffrey's prior
res1J <- LogicRegr(formula = formula1, data = data.example, family = "Bernoulli", prior = "J", report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.2, p.surv = 0.2, ncores = 32)
print(res1J$feat.stat)
res1J <- suppressWarnings(
suppressMessages(
LogicRegr(
formula = formula1, data = data.example, family = "Bernoulli", prior = "J",
report.level = 0.5, d = 15, cmax = 2, kmax = 15, p.and = 0.9, p.not = 0.2,
p.surv = 0.2, ncores = threads, print.freq = 0L
)
)
)
test_that("outputs are correct", {
expect_equal(
res4G$feat.stat,
matrix(
c(
"I(V9)", "1", "I(V5)", "1", "I(((V5))|((V9)))", "1", "I(V8)",
"0.999798166079656", "I(V11)", "0.998980768894478", "I(V1)",
"0.996827860563252", "I(V4)", "0.989626671838795"
),
nrow = 7, ncol = 2, byrow = TRUE
)
)
expect_equal(
res4J$feat.stat,
matrix(
c(
"I(((V9))&((V5)))", "0.999999999999948", "I(((V8))&((V11)))",
"0.9999999984096", "I(((V4))&((V1)))", "0.99896404225898"
),
nrow = 3, ncol = 2, byrow = TRUE
)
)
expect_equal(
res1G$feat.stat,
matrix(
c(
"I(((V9))&((V5)))", "1", "I(V8)", "0.999999999999935",
"I(((V1))&(1-(V4)))", "0.999986154466532",
"I(V11)", "0.999445448885581", "I(V1)", "0.997930154511171",
"I(V4)", "0.997375287784621"
),
nrow = 6, ncol = 2, byrow = TRUE
)
)
expect_equal(
res1J$feat.stat,
matrix(
c(
"I(V1)", "1", "I(V4)", "1", "I(V11)", "1", "I(((V4))|((V8)))", "1",
"I(((V5))&((V9)))", "1"
),
nrow = 5, ncol = 2, byrow = TRUE
)
)
})

0 comments on commit 6de6170

Please sign in to comment.