R Under development (unstable) (2024-02-05 r85863 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. > library(pcalg) > possible_p <- c(seq(10,30,by=10)) > possible_neighb_size <- c(3:10) > > num_setings <-10 > num_rep <- 2 > > ##amount of bg knowledge in percent > p_bg <- seq(0,1,by=0.1) > > counterr <- 0 > > seed1 <- 1 > for (r in seed1:num_setings) { + set.seed(r) + ##Then we can choose one setting: + p <- sample(possible_p,1) + neigh <- sample(possible_neighb_size,1) + prob <- round(neigh/(p-1),3) + + ## cat("r=",r,", p=",p,", nb=",neigh,"\n") + + ## then for every setting selected we can generate i.e. 20 graph + for (d in 1: num_rep){ + ##i left the seed inside for now + ## i am not sure which seed option makes it easier for us + ## feel free to move this around + #seed <- sample(possible_seed,1) + #set.seed(seed) + + ##get graph + g <- pcalg:::randomDAG(p, prob) ## true DAG + cpdag <- dag2cpdag(g) ## true CPDAG + + ## get adjacency matrix + cpdag.amat <- t(as(cpdag,"matrix")) + dag.amat <- t(as(g,"matrix")) + dag.amat[dag.amat != 0] <- 1 + + amat.cpdag <- cpdag.amat + if (!isValidGraph(amat.cpdag,type="cpdag")){ + counterr <- counterr +1 + } + if (!isValidGraph(dag.amat,type="dag")){ + counterr <- counterr +1 + } + } + } > counterr [1] 0 > > amat <- matrix(c(0,1,0, 0,0,1, 0,0,0), 3,3) > colnames(amat) <- rownames(amat) <- letters[1:3] > ## graph::plot(as(t(amat), "graphNEL")) > stopifnot(isValidGraph(amat = amat, type = "dag"), ## is a valid DAG + !isValidGraph(amat = amat, type = "cpdag"), ## is a valid CPDAG + isValidGraph(amat = amat, type = "pdag") ) ## is a valid PDAG > > ## a -- b -- c > amat <- matrix(c(0,1,0, 1,0,1, 0,1,0), 3,3) > colnames(amat) <- rownames(amat) <- letters[1:3] > ## plot(as(t(amat), "graphNEL")) > stopifnot(!isValidGraph(amat = amat, type = "dag"), ## not a valid DAG + isValidGraph(amat = amat, type = "cpdag"), ## is a valid CPDAG + isValidGraph(amat = amat, type = "pdag") )## is a valid PDAG > > ## a -- b -- c -- d -- a > amat <- matrix(c(0,1,0,1, 1,0,1,0, 0,1,0,1, 1,0,1,0), 4,4) > colnames(amat) <- rownames(amat) <- letters[1:4] > ## plot(as(t(amat), "graphNEL")) > stopifnot(!isValidGraph(amat = amat, type = "dag"), ## not a valid DAG + !isValidGraph(amat = amat, type = "cpdag"), ## not a valid CPDAG + !isValidGraph(amat = amat, type = "pdag") ) ## not a valid PDAG > > proc.time() user system elapsed 1.78 0.15 1.92