1

I am starting to use the pcalg package, and I am having some doubts regarding the functions pdag2allDags and addBgKnowledge:

I am using the sample data gmG provided by the package

library(pcalg)
library(Rgraphviz)
data("gmG")
df<-gmG$x
suffStat <- list(C = cor(df), n = nrow(df))
fci.fit<-fci(suffStat, indepTest = gaussCItest, p = ncol(df),alpha = 0.01)
plot(fci.fit)

I would like to get the all equivalent DAGs. From documentation, it should work using the function pdag2allDags (from here). We should need just to get the amat (adjacent matrix) data.

From the specifications on the documentation, I thought the following should work...

plotAllDags <- function(res) {
    require(graph)
    p <- sqrt(ncol(res$dags))
    nDags <- ceiling(sqrt(nrow(res$dags)))
    par(mfrow = c(nDags, nDags))
    for (i in 1:nrow(res$dags)) {
        tmp <- matrix(res$dags[i,],p,p)
        colnames(tmp) <- rownames(tmp) <- res$nodeNms
        plot(as(tmp, "graphNEL"))
    }
}

res1<-pdag2allDags(as(fci.fit,"amat"))
plotAllDags(res1)

But, instead, it returns:

Error in sqrt(ncol(res$dags)) : non-numeric argument to mathematical function

We see also the amat in the fci's object. So, I tried:

res2<-pdag2allDags(fci.fit@amat)
plotAllDags(res2)

It also returns the same:

Error in sqrt(ncol(res$dags)) : non-numeric argument to mathematical function

But if I use the pc algorithm, it works:

pc.fit<-pc(suffStat, indepTest = gaussCItest, p = ncol(df),alpha = 0.01)
plot(pc.fit)
res0<-pdag2allDags(as(pc.fit,"amat"))
plotAllDags(res0)

What is going on? Is not pdag2allDags intended to deal with all amat objects (pc, fci, rfci, etc.)? I was not able to find any other ...allDags function in the documentation. How to get the all equivalent DAGs from the output of the fci function?

The same happens for the function addBgKnowledge. It works for pc:

pc.amat2<-addBgKnowledge(gInput = pc.fit@graph,x=1,y=2)
plot(pc.amat2)

but not for fci, even the documentation saying it uses an amat

fci.amat2<-addBgKnowledge(gInput = as(fci.fit,"amat"),x=1,y=2)
plot(as(t(fci.amat2),"graphNEL"))

It delivers:

Error in h(simpleError(msg, call)) : error in evaluating the argument 'x' in selecting a method for function 'plot': argument is not a matrix

hamagust
  • 728
  • 2
  • 10
  • 28

1 Answers1

1

I believe it should not be an ideal answer, it should have a better way to do it, but as I was not getting any answer, I have made some attempts, and I am sharing how I had circumvented the issue for now...

If we assess the structure, both are indeed dealing with adjacency matrix (amat)... but they are different... the pc output comprises a 'cpdag' type amat, and fci output, a 'pag' type amat...

as(pc.fit,"amat")

we get:

Adjacency Matrix 'amat' (8 x 8) of type ‘cpdag’:
  1 2 3 4 5 6 7 8
1 . 1 . . . . . .
2 1 . 1 . 1 . . .
3 . 1 . . . . . .
4 . . . . . . . .
5 . 1 . . . . . .
6 1 . . . 1 . . .
7 . . . . . 1 . .
8 1 . . . 1 . . .

and with

as(fci.fit,"amat")

we get:

Adjacency Matrix 'amat' (8 x 8) of type ‘pag’:
  1 2 3 4 5 6 7 8
1 . 1 . . . 2 . 2
2 1 . 1 . 1 . . .
3 . 1 . . . . . .
4 . . . . . . . .
5 . 1 . . . 2 . 2
6 3 . . . 3 . 2 .
7 . . . . . 3 . .
8 3 . . . 3 . . .

or, with

fci.fit@amat

we get

  1 2 3 4 5 6 7 8
1 0 1 0 0 0 2 0 2
2 1 0 1 0 1 0 0 0
3 0 1 0 0 0 0 0 0
4 0 0 0 0 0 0 0 0
5 0 1 0 0 0 2 0 2
6 3 0 0 0 3 0 2 0
7 0 0 0 0 0 3 0 0
8 3 0 0 0 3 0 0 0

The pcalg documentation on pdag2allDags does not explicitly mention that...

  • function definition says 'a Given Partially Directed Acyclic Graph (PDAG)'.
  • the example shows an 'amat' defined as matrix of zeros and ones, without any connection to fci or pc functions

But the solution I found starts with taking the 'pag' type amat and converting it to 'cpdag'...

Well... the 'pag' type amat indeed is not comprised just by zeroes and ones... But we now from here:

  • 3 means the column is the cause and line the consequence...
  • and 2 means the column is the consequence and line is the cause...

so, we replace them creating a new function

pag2cpdag<-function(fciAlgo){
  res<-as(fciAlgo,"amat")#a amat type pag
  res[res==3]<-1
  res[res==2]<-0
  return(res)
}

So we get with

pag2cpdag(fci.fit)

the following output:

Adjacency Matrix 'amat' (8 x 8) of type ‘pag’:
  1 2 3 4 5 6 7 8
1 . 1 . . . . . .
2 1 . 1 . 1 . . .
3 . 1 . . . . . .
4 . . . . . . . .
5 . 1 . . . . . .
6 1 . . . 1 . . .
7 . . . . . 1 . .
8 1 . . . 1 . . .

I was not able to change the type, but it delivers a matrix in the cpdag's amat format. And can be applied to pdag2allDags...

res1<-pdag2allDags(pag2cpdag(fci.fit))
plotAllDags(res1)

to return adequatelly the all equivalent DAGs...

But it still does not work for addBgKnowledge with objects generated by fci even the documentation saying it uses an amat

fci.amat2<-addBgKnowledge(gInput = pag2cpdag(fci.fit),x=1,y=2)
plot(as(t(fci.amat2),"graphNEL"))

it returns:

Error in h(simpleError(msg, call)) : error in evaluating the argument 'x' in selecting a method for function 'plot': no method or default for coercing “amat” to “graphNEL”

The output of

fci.amat2

appears to be ok...

Adjacency Matrix 'amat' (8 x 8) of type ‘pag’:
  1 2 3 4 5 6 7 8
1 . . . . . . . .
2 1 . . . . . . .
3 . 1 . . . . . .
4 . . . . . . . .
5 . 1 . . . . . .
6 1 . . . 1 . . .
7 . . . . . 1 . .
8 1 . . . 1 . . .

Looking at example, a simple matrix work for plotting it as graphNEL... so, just converting to a simple matrix solves the issue:

pc.amat2<-addBgKnowledge(pag2cpdag(fci.fit),x=1,y=2)
class(pc.amat2)<-"matrix"
plot(as(t(pc.amat2),"graphNEL"))
hamagust
  • 728
  • 2
  • 10
  • 28