R Under development (unstable) (2024-09-06 r87103 ucrt) -- "Unsuffered Consequences" Copyright (C) 2024 The R Foundation for Statistical Computing Platform: x86_64-w64-mingw32/x64 R is free software and comes with ABSOLUTELY NO WARRANTY. You are welcome to redistribute it under certain conditions. Type 'license()' or 'licence()' for distribution details. R is a collaborative project with many contributors. Type 'contributors()' for more information and 'citation()' on how to cite R or R packages in publications. Type 'demo()' for some demos, 'help()' for on-line help, or 'help.start()' for an HTML browser interface to help. Type 'q()' to quit R. > ####' Tests the causal inference algorithms for interventional data: > ####' GIES, GES, DP > ####' > ####' @author Alain Hauser > ####' $Id: test_gies.R 454 2017-07-11 15:03:05Z mkalisch $ > > cat("Testing the causal inference algorithms for interventional data:\n") Testing the causal inference algorithms for interventional data: > > library(pcalg) > > source(system.file(package="Matrix", "test-tools-1.R", mustWork=TRUE)) Loading required package: tools > ##--> showProc.time(), assertError(), relErrV(), ... > > load("test_bicscore.rda") # in directory tests/ i.e., typically *not* installed > # str(gauss.data) > p <- ncol(gauss.data) > > (doExtras <- pcalg:::doExtras()) [1] FALSE > DBG <- if(doExtras) TRUE else FALSE # no debugging by default > ## Tolerance for numerical comparison > tol <- sqrt(.Machine$double.eps) # = default for all.equal() > > ## Define all test settings > settings <- expand.grid( + fcn = c("gies", "gds"), + cpp = c(FALSE, TRUE), + format = c("scatter", "raw"), + stringsAsFactors = FALSE) > nreps <- 5 > > for (m in seq_along(settings)) { + cat(sprintf("Algorithm: %s, C++: %s, storage format: %s\n", + settings$fcn[m], settings$cpp[m], settings$format[m])) + + for (i in 1:nreps) { + perm <- 1:nrow(gauss.data) + + ## Randomly permute data + if (i > 1) { + set.seed(i) + perm <- sample(perm) + } + + score <- new("GaussL0penIntScore", + targets = gauss.targets, + target.index = gauss.target.index[perm], + data = gauss.data[perm, ], + format = settings$format[m], + use.cpp = settings$cpp[m]) + f <- get(settings$fcn[m]) + est.graph <- f(score, verbose = DBG) + for (j in 1:p) { + if (!isTRUE(all.equal(est.graph$essgraph$.in.edges[[j]], + gauss.parents[[j]], tolerance = tol))) + stop("Parents are not estimated correctly.") + } + showProc.time() + } + cat("[Ok]\n") + } Algorithm: gies, C++: FALSE, storage format: scatter Time (user system elapsed): 0.21 0 0.2 Time (user system elapsed): 0.03 0.02 0.05 Time (user system elapsed): 0.04 0 0.04 Time (user system elapsed): 0.05 0 0.05 Time (user system elapsed): 0.05 0 0.05 [Ok] Algorithm: gds, C++: FALSE, storage format: scatter Time (user system elapsed): 0.06 0 0.06 Time (user system elapsed): 0.06 0 0.06 Time (user system elapsed): 0.07 0 0.06 Time (user system elapsed): 0.06 0 0.07 Time (user system elapsed): 0.06 0 0.06 [Ok] Algorithm: gies, C++: TRUE, storage format: scatter Time (user system elapsed): 0.02 0 0.02 Time (user system elapsed): 0 0 0 Time (user system elapsed): 0.01 0 0.01 Time (user system elapsed): 0.02 0 0.02 Time (user system elapsed): 0 0 0 [Ok] > > ## Test compatibility with deprecated calling conventions > cat("Compatibility with deprecated calling conventions... ") Compatibility with deprecated calling conventions... > score <- new("GaussL0penIntScore", + targets = gauss.targets, + target.index = gauss.target.index, + data = gauss.data) > > warningIssued <- FALSE > tryCatch(est.graph <- gies(p, gauss.targets, score), + warning = function(w) warningIssued <<- TRUE) > if (!warningIssued) { + stop("No warning issued for old calling conventions.") + } else { + for (j in 1:p) { + if (!isTRUE(all.equal(est.graph$essgraph$.in.edges[[j]], + gauss.parents[[j]], tolerance = tol))) + stop("Parents are not estimated correctly.") + } + } > warningIssued <- FALSE > tryCatch(est.graph <- gies(p = p, targets = gauss.targets, score = score), + warning = function(w) warningIssued <<- TRUE) > if (!warningIssued) { + stop("No warning issued for old calling conventions.") + } > cat("[OK]\n") [OK] > > > ## Test stepwise execution of GIES > cat(if(doExtras)"\n\n", "GIES stepwise", if(doExtras)":\n" else ": ... ", + if(doExtras) paste0(paste(rep("=", 14), collapse=""), "\n"), + sep = "") GIES stepwise: ... > for (cpp in c(FALSE, TRUE)) { + ## Randomly permute data + for (i in 1:nreps) { + perm <- 1:nrow(gauss.data) + if (i > 1) { + set.seed(i) + perm <- sample(perm) + } + score <- new("GaussL0penIntScore", + targets = gauss.targets, + target.index = gauss.target.index[perm], + data = gauss.data[perm, ], + use.cpp = cpp) + + ## Stepwise execution + essgraph <- new("EssGraph", nodes = as.character(1:p), score = score) + cont <- TRUE + while(cont) { + cont <- FALSE + while(essgraph$greedy.step("forward")) cont <- TRUE + while(essgraph$greedy.step("backward")) cont <- TRUE + while(essgraph$greedy.step("backward")) cont <- TRUE + } + for (i in 1:p) { + if(doExtras) cat(" use.cpp = ", cpp,"; i = ", i, "\n", sep="") + if (!isTRUE(all.equal(est.graph$essgraph$.in.edges[[i]], + gauss.parents[[i]], tolerance = tol))) + stop("Parents are not estimated correctly.") + } + showProc.time() + } + } Time (user system elapsed): 0.11 0 0.11 Time (user system elapsed): 0.06 0 0.06 Time (user system elapsed): 0.06 0 0.06 Time (user system elapsed): 0.06 0 0.06 Time (user system elapsed): 0.07 0 0.07 Time (user system elapsed): 0.01 0 0.01 Time (user system elapsed): 0 0 0 Time (user system elapsed): 0.02 0 0.02 Time (user system elapsed): 0 0 0 Time (user system elapsed): 0.01 0 0.01 > > cat(if(doExtras) "\n", "Done.\n") Done. > > ## add test to confirm bug-fix in 2.5-0 > new("GaussParDAG", "a") ## produced a warning prior to 2.5-0 Reference class object of class "GaussParDAG" Field ".nodes": [1] "a" Field ".in.edges": [[1]] integer(0) Field ".params": [[1]] [1] 0 0 > > proc.time() user system elapsed 1.96 0.09 2.04