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.

> ## Sep 2024: Stopped using output because of volative behaviour
> ## of interER, power, geometric, barabasi
> 
> library(pcalg)
> suppressWarnings(RNGversion("3.5.0"))
> ## setwd("/sfs/u/kalischm/research/packages/unifDAGs/")
> ## source("aux_general.R")
> ## source("randDAG.R")
> 
> ### Check all methods: ----------------------------------------------
> 
> ## MM hack: extract them from the randDAG() function definition
> body. <- body(randDAG)
> is.switch <- function(P) !is.symbol(P) && identical(as.symbol("switch"), P[[1]])
> switchCall <- body.[vapply(body., is.switch, NA)][[1]]
> stopifnot(identical(as.symbol("switch"), switchCall[[1]]))
> (rDAGmeths <- names(switchCall)[-c(1:2, length(switchCall))])
[1] "er"        "regular"   "watts"     "bipartite" "barabasi"  "geometric"
[7] "power"     "interEr"  
> rDAGall <- function(n, d, ...)
+   sapply(rDAGmeths, function(meth) randDAG(n,d, method=meth, ...),
+          simplify=FALSE)
> set.seed(37)
> rD.10.4 <- rDAGall(10, 4) # 2024-02: no "low-level warning" anymore
> ##     , warning = function(w) {
> ##         rDAG.warn <<- conditionMessage(w); invokeRestart("muffleWarning") })
> ## ## with a low-level warning:
> ## ## IGNORE_RDIFF_BEGIN
> ## rDAG.warn
> ## ## IGNORE_RDIFF_END
> ## stopifnot(grepl("graph_molloy_.*Cannot shuffle graph", rDAG.warn))
> ##
> rD.10.4 # looks ok
$er
A graphNEL graph with directed edges
Number of Nodes = 10 
Number of Edges = 21 

$regular
A graphNEL graph with directed edges
Number of Nodes = 10 
Number of Edges = 20 

$watts
A graphNEL graph with directed edges
Number of Nodes = 10 
Number of Edges = 20 

$bipartite
A graphNEL graph with directed edges
Number of Nodes = 10 
Number of Edges = 15 

$barabasi
A graphNEL graph with directed edges
Number of Nodes = 10 
Number of Edges = 21 

$geometric
A graphNEL graph with directed edges
Number of Nodes = 10 
Number of Edges = 13 

$power
A graphNEL graph with directed edges
Number of Nodes = 10 
Number of Edges = 16 

$interEr
A graphNEL graph with directed edges
Number of Nodes = 10 
Number of Edges = 19 

> ## Show, but ignore the package startup messages:
> ## IGNORE_RDIFF_BEGIN
> stopifnot( require("graph") )
Loading required package: graph
Loading required package: BiocGenerics

Attaching package: 'BiocGenerics'

The following objects are masked from 'package:stats':

    IQR, mad, sd, var, xtabs

The following objects are masked from 'package:base':

    Filter, Find, Map, Position, Reduce, anyDuplicated, aperm, append,
    as.data.frame, basename, cbind, colnames, dirname, do.call,
    duplicated, eval, evalq, get, grep, grepl, intersect, is.unsorted,
    lapply, mapply, match, mget, order, paste, pmax, pmax.int, pmin,
    pmin.int, rank, rbind, rownames, sapply, setdiff, table, tapply,
    union, unique, unsplit, which.max, which.min

> ## IGNORE_RDIFF_END
> 
> stopifnot(vapply(rD.10.4, isDirected, NA),
+           vapply(rD.10.4, inherits, NA, what="graph"))
> ## nice plot of all 8 :
> # op <- par(mfrow=c(4,2))
> # invisible(lapply(names(rD.10.4), function(nm) plot(rD.10.4[[nm]], main=nm)))
> # par(op)
> 
> 
> str(outs <- lapply(rD.10.4, leaves, "out"))
List of 8
 $ er       : chr "3"
 $ regular  : chr [1:2] "3" "8"
 $ watts    : chr [1:2] "4" "5"
 $ bipartite: chr [1:3] "5" "7" "8"
 $ barabasi : chr "4"
 $ geometric: chr [1:3] "3" "4" "7"
 $ power    : chr [1:4] "2" "7" "8" "9"
 $ interEr  : chr [1:2] "3" "7"
> if(packageVersion("igraph") < "2.0") ## currently fails
+ stopifnot(identical( outs,
+     list(er = "3", regular = c("1", "5", "6"), watts = c("3", "4", "6"),
+          bipartite = c("1", "2", "5"), barabasi = c("4", "8"),
+          geometric = c("4", "7"), power = c("4", "5", "9"),
+          interEr = c("3", "7"))
+ ))
> 
> str(ins <- lapply(rD.10.4, leaves, "in"))
List of 8
 $ er       : chr [1:3] "1" "4" "7"
 $ regular  : chr [1:2] "1" "10"
 $ watts    : chr [1:2] "7" "10"
 $ bipartite: chr "3"
 $ barabasi : chr [1:3] "1" "6" "7"
 $ geometric: chr [1:3] "5" "6" "9"
 $ power    : chr [1:3] "4" "5" "10"
 $ interEr  : chr [1:2] "2" "4"
> if(packageVersion("igraph") < "2.0") ## currently fails
+ stopifnot(identical( ins,
+     list(er = c("1", "4", "7"), regular = c("3", "7", "10"),
+          watts = c("1", "8"), bipartite = c("4", "6"),
+          barabasi = c("6", "7"), geometric = c("5", "10"),
+          power = c("2", "7"), interEr = c("8", "10"))
+   ))
> 
> set.seed(47)
> (rD.12.2 <- rDAGall(12, 2))
$er
A graphNEL graph with directed edges
Number of Nodes = 12 
Number of Edges = 9 

$regular
A graphNEL graph with directed edges
Number of Nodes = 12 
Number of Edges = 12 

$watts
A graphNEL graph with directed edges
Number of Nodes = 12 
Number of Edges = 12 

$bipartite
A graphNEL graph with directed edges
Number of Nodes = 12 
Number of Edges = 12 

$barabasi
A graphNEL graph with directed edges
Number of Nodes = 12 
Number of Edges = 13 

$geometric
A graphNEL graph with directed edges
Number of Nodes = 12 
Number of Edges = 6 

$power
A graphNEL graph with directed edges
Number of Nodes = 12 
Number of Edges = 11 

$interEr
A graphNEL graph with directed edges
Number of Nodes = 12 
Number of Edges = 10 

> if(packageVersion("igraph") < "2.0") ## currently fails
+ stopifnot(exprs = {
+     vapply(rD.12.2, isDirected, NA)
+     vapply(rD.12.2, numNodes, 1) == 12
+     identical(vapply(rD.12.2, numEdges, 1),
+               setNames(c(9, 12, 12, 11, 11, 11, 13, 8), rDAGmeths))
+ })
> 
> ##---------------------------------------------------------------------------
> 
> ## Use the output here
> require(Matrix)
Loading required package: Matrix
> # lapply(rD.10.4, function(g) as(as(g, "Matrix"),"nMatrix"))
> # lapply(rD.12.2, function(g) as(as(g, "Matrix"),"nMatrix"))
> 
> ## Minimal checks on graphs generated via igraph
> stopifnot( require("graph") )
> set.seed(37)
> dagList <- vector(mode = "list", length = 6)
> dagList[[1]] <- randDAG(10, 4, "regular")
> dagList[[2]] <- randDAG(10, 4, "watts")
> dagList[[3]] <- randDAG(10, 4, "er")
> dagList[[4]] <- randDAG(10, 4, "bipartite")
> dagList[[5]] <- randDAG(10, 4, "geometric")
> dagList[[6]] <- randDAG(10, 4, "interEr", par2 = 0.5)
> 
> ## number of nodes
> stopifnot(all.equal(
+   sapply(dagList, numNodes),
+   rep(10,6)
+ ))
> ## number of edges
> stopifnot(all.equal(
+   sapply(dagList, numEdges),
+   c(20,20,16,15,10,18)
+ ))
> 
> ## Use the output here -- FIXME: check more
> require(Matrix)
> # lapply(rD.10.4, function(g) as(as(g, "Matrix"),"nMatrix"))
> # lapply(rD.12.2, function(g) as(as(g, "Matrix"),"nMatrix"))
> 
> ## check weights
> set.seed(123)
> n <- 100
> g <- randDAG(n=n,d=3, wFUN=list(runif,min=0,max=1))
> m <- wgtMatrix(g)
> stopifnot(sum(m != 0) == 137)
> v <- as.numeric(m)
> v <- v[v!=0]
> ## dput(as.vector(summary(v, digits=7)))
> stopifnot(all.equal(as.vector(summary(v, digits=7)),
+                     c(0.008103577, 0.2589966, 0.5287397,
+                       0.5232445, 0.8159941, 0.9915566)))
> ct <- cut(x=v, breaks=seq(0,1,by=0.1))
> stopifnot(all.equal(chisq.test(as.numeric(table(ct)), p = rep(0.1,10))$p.value,
+                     0.3101796548))
> 
> ## check generation of negative weights (fixed Bug)
> set.seed(123) ; tmp1 <- randDAG(3,2, wFUN = list(runif, min =  2, max =  2))
> set.seed(123) ; tmp2 <- randDAG(3,2, wFUN = list(runif, min = -2, max = -2))
> stopifnot(unlist(tmp1@edgeData@data) ==  2,
+ 	  unlist(tmp2@edgeData@data) == -2 )
> 
> proc.time()
   user  system elapsed 
   2.25    0.21    2.45