1 2context("Graphlets") 3 4sortgl <- function(x) { 5 cl <- lapply(x$cliques, sort) 6 n <- sapply(cl, length) 7 list(cliques=cl[order(n)], thresholds=x$thresholds[order(n)]) 8} 9 10test_that("Graphlets work for some simple graphs", { 11 library(igraph) 12 13 g <- make_full_graph(5) 14 E(g)$weight <- 1 15 gl <- graphlet_basis(g) 16 17 expect_that(names(gl), equals(c("cliques", "thresholds"))) 18 expect_that(length(gl$cliques), equals(1)) 19 expect_that(sort(gl$cliques[[1]]), equals(1:vcount(g))) 20 expect_that(gl$thresholds, equals(1)) 21 22 g2 <- make_full_graph(5) 23 E(g2)$weight <- 1 24 E(g2)[1%--%2]$weight <- 2 25 gl2 <- sortgl(graphlet_basis(g2)) 26 27 expect_that(gl2, equals(list(cliques=list(1:2, 1:5), thresholds=c(2,1)))) 28}) 29 30test_that("Graphlets filtering works", { 31 library(igraph) 32 gt <- data.frame(from =c("A", "A", "B", "B", "B", "C", "C", "D"), 33 to =c("B", "C", "C", "D", "E", "D", "E", "E"), 34 weight=c( 8 , 8 , 8 , 5 , 5 , 5 , 5 , 5 )) 35 36 g <- graph_from_data_frame(gt, directed=FALSE, vertices=data.frame(LETTERS[1:5])) 37 gl <- sortgl(graphlet_basis(g)) 38 39 expect_that(gl$cliques, equals(list(1:3, 2:5))) 40 expect_that(gl$thresholds, equals(c(8, 5))) 41}) 42 43## Naive version of graphlets 44 45unvs <- function(x) lapply(x, as.vector) 46 47threshold.net <- function(graph, level) { 48 N <- vcount(graph) 49 graph.t <- delete_edges(graph, which(E(graph)$weight < level)) 50 51 clqt <- unvs(max_cliques(graph.t)) 52 clqt <- lapply(clqt, sort) 53 clqt[order(sapply(clqt, length), decreasing=TRUE)] 54} 55 56graphlets.old <- function(graph) { 57 58 if (!is_weighted(graph)) { stop("Graph not weighted") } 59 if (min(E(graph)$weight) <= 0 || any(!is.finite(E(graph)$weight))) { 60 stop("Edge weights must be non-negative and finite") 61 } 62 63 ## Do all thresholds 64 cl <- lapply(sort(unique(E(graph)$weight)), function(w) { 65 threshold.net(graph, w) 66 }) 67 68 ## Put the cliques in one long list 69 clv <- unlist(cl, recursive=FALSE) 70 71 ## Sort the vertices within the cliques 72 cls <- lapply(clv, sort) 73 74 ## Delete duplicate cliques 75 clu <- unique(cls) 76 77 ## Delete cliques that consist of single vertices 78 clf <- clu[sapply(clu, length) != 1] 79 80 clf 81} 82 83test_that("Graphlets work for a bigger graph", { 84 library(igraph) 85 set.seed(42) 86 g <- make_graph("zachary") 87 E(g)$weight <- sample(1:5, ecount(g), replace=TRUE) 88 89 gl <- graphlet_basis(g) 90 gl2 <- graphlets.old(g) 91 92 glo <- sort(sapply(unvs(gl$cliques), paste, collapse="-")) 93 gl2o <- sort(sapply(gl2, paste, collapse="-")) 94 95 expect_that(glo, equals(gl2o)) 96}) 97 98graphlets.project.old <- function(graph, cliques, iter, Mu=NULL) { 99 100 if (!is_weighted(graph)) { stop("Graph not weighted") } 101 if (min(E(graph)$weight) <= 0 || any(!is.finite(E(graph)$weight))) { 102 stop("Edge weights must be non-negative and finite") 103 } 104 if (length(iter) != 1 || !is.numeric(iter) || 105 !is.finite(iter) || iter != as.integer(iter)) { 106 stop("`iter' must be a non-negative finite integer scalar") 107 } 108 109 clf <- cliques 110 111 ## Create vertex-clique list first 112 vcl <- vector(length=vcount(graph), mode="list") 113 for (i in 1:length(clf)) { 114 for (j in clf[[i]]) { 115 vcl[[j]] <- c(vcl[[j]], i) 116 } 117 } 118 119 ## Create edge-clique list from this, it is useful to have the edge list 120 ## of the graph at hand 121 el <- as_edgelist(graph, names=FALSE) 122 ecl <- vector(length=ecount(graph), mode="list") 123 for (i in 1:ecount(graph)) { 124 edge <- el[i,] 125 ecl[[i]] <- intersect(vcl[[edge[1]]], vcl[[edge[2]]]) 126 } 127 128 ## We will also need a clique-edge list, the edges in the cliques 129 system.time({ 130 cel <- vector(length=length(clf), mode="list") 131 for (i in 1:length(ecl)) { 132 for (j in ecl[[i]]) { 133 cel[[j]] <- c(cel[[j]], i) 134 } 135 } 136 }) 137 138 ## OK, we are ready to do the projection now 139 if (is.null(Mu)) { Mu <- rep(1, length(clf)) } 140 origw <- E(graph)$weight 141 w <- numeric(length(ecl)) 142 a <- sapply(clf, function(x) length(x) * (length(x) + 1) / 2) 143 for (i in 1:iter) { 144 for (j in 1:length(ecl)) { 145 w[j] <- sum(Mu[ ecl[[j]] ]) 146 } 147 for (j in 1:length(clf)) { 148 Mu[j] <- Mu[j] * sum(origw[cel[[j]]] / (w[cel[[j]]] + .0001)) / a[j] 149 } 150 } 151 152 ## Sort the cliques according to their weights 153 Smb <- sort(Mu, decreasing=TRUE, index=TRUE) 154 list(cliques=clf[Smb$ix], Mu=Mu[Smb$ix]) 155} 156 157test_that("Graphlet projection works", { 158 library(igraph) 159 160 D1 <- matrix(0, 5, 5) 161 D2 <- matrix(0, 5, 5) 162 D3 <- matrix(0, 5, 5) 163 D1[1:3, 1:3] <- 2 164 D2[3:5, 3:5] <- 3 165 D3[2:5, 2:5] <- 1 166 g <- graph_from_adjacency_matrix(D1 + D2 + D3, mode="undirected", weighted=TRUE) 167 g <- simplify(g) 168 169 gl <- graphlet_basis(g) 170 glp <- graphlets(g) 171 glp2 <- graphlets.project.old(g, cliques=gl$cliques, iter=1000) 172 173 glp$cliques <- unvs(glp$cliques) 174 expect_that(glp, equals(glp2)) 175}) 176