I have an adjacency matrix from the package 'bnlearn' using the function amat (the matrix is acyclical). For example:
+---+-------------------------------+
| | 1 2 3 4 5 |
+---+-------------------------------+
| 1 | 0 1 0 0 0 |
| 2 | 0 0 1 0 0 |
| 3 | 0 0 0 1 0 |
| 4 | 0 0 0 0 1 |
| 5 | 0 0 0 0 0 |
+---+-------------------------------+
I need to find the complete dependency matrix from this. For one lag dependency matrix I can use:
New_matrix<- if(old_matrix+old_matrix*old_matrix)>0 then 1 else 0
For two lag dependency matrix I can use:
New_matrix_2<- if(new_matrix+new_matrix*old_matrix)>0 then 1 else 0
The problem is I don't know where the adjacency is complete, that is for how many iterations do I run this to get to the final matrix with all interdependencies incorporated?
+---+-------------------------------+
| | 1 2 3 4 5 |
+---+-------------------------------+
| 1 | 0 1 1 1 1 |
| 2 | 0 0 1 1 1 |
| 3 | 0 0 0 1 1 |
| 4 | 0 0 0 0 1 |
| 5 | 0 0 0 0 0 |
+---+-------------------------------+
For this, the answer is 3 iterations. But the matrix that I need to solve this for is 500x500. Is there a direct way to arrive at the complete adjacency matrix?
To find the paths from all nodes, it is probably easier to use the igraph package.
Using your example,
library(bnlearn)
library(igraph)
# Create BN in your example
g <- empty.graph(LETTERS[1:5])
amat(g) <- rbind(cbind(0, diag(4)),0)
amat(g)
# A B C D E
# A 0 1 0 0 0
# B 0 0 1 0 0
# C 0 0 0 1 0
# D 0 0 0 0 1
# E 0 0 0 0 0
# Convert to igraph object using BN adj. matrix as input
g1 <- graph_from_adjacency_matrix(amat(g))
# You can find all ancestors for each node by using
# the mode="in" argument, and order to specify the depth of the search
neighborhood(g1, order=nrow(amat(g)), mode="in")
# Similarly, you can get the full connected graph
# using the same options
ances <- connect(g1, order=nrow(amat(g)), mode="in" )
get.adjacency(ances, sparse=FALSE)
# A B C D E
# A 0 1 1 1 1
# B 0 0 1 1 1
# C 0 0 0 1 1
# D 0 0 0 0 1
# E 0 0 0 0 0
Alternatively, you can use matrix exponential
m <- amat(g)
1* as.matrix((Matrix::expm(m) - diag(ncol(m))) > 0)