.Rhistory

Results3[[i]] = list(Beta = betas[i], Communities = temp$Communities)
Mean.Score[i] = temp$Mean.Score
Number.Communities[i] = length(temp$Communities)
}
if(length(Results2) > 0){
betas = seq(0.01, 1, by = 0.01)
Results3 = list()
Number.Communities = rep(0,length(betas))
Mean.Score = rep(0,length(betas))
for(i in 1:length(betas)){
temp = cleanup(Results2, betas[i])
Results3[[i]] = list(Beta = betas[i], Communities = temp$Communities)
Mean.Score[i] = temp$Mean.Score
Number.Communities[i] = length(temp$Communities)
}
}
if(length(Results2) > 0){
betas = seq(0.01, 1, by = 0.01)
Results3 = list()
Number.Communities = rep(0,length(betas))
Mean.Score = rep(0,length(betas))
for(i in 1:length(betas)){
temp = cleanup(Results2, betas[i])
Results3[[i]] = list(Beta = betas[i], Communities = temp$Communities)
Mean.Score[i] = temp$Mean.Score
Number.Communities[i] = length(temp$Communities)
}
}
temp
length(Results2)
betas = seq(0.01, 1, by = 0.01)
Results3 = list()
Number.Communities = rep(0,length(betas))
Number.Communities = rep(0, length(betas))
Number.Communities
Mean.Score = rep(0, length(betas))
temp
for(i in 1:length(betas)){
temp = cleanup(Results2, betas[i])
Results3[[i]] = list(Beta = betas[i], Communities = temp$Communities)
Mean.Score[i] = temp$Score
Number.Communities[i] = length(temp$Communities)
}
i
temp$Score
temp
length(betas)
temp = cleanup(Results2, betas[i])
temp
cleanup = function(Results, beta){
k = length(Results)
if(k < 2){
Results <- Results[[1]]
Results$Mean.Score <- Results$Score
return(Results)
}
indx.rm = numeric()
if(k > 1){
for(i in 1:(k-1)){
for(j in (i+1):k){
match = length(intersect(Results[[i]]$B, Results[[j]]$B))*length(intersect(Results[[i]]$I, Results[[j]]$I))/(min(length(Results[[i]]$B)*length(Results[[i]]$I),length(Results[[j]]$B)*length(Results[[j]]$I)))
pot = c(i,j)
indx = pot[which.min(c(Results[[i]]$Score, Results[[j]]$Score))]
if(match > beta-0.00001){indx.rm = c(indx.rm,indx)}
}
}
}
cleanup = function(Results, beta){
k = length(Results)
if(k < 2){
Results <- Results[[1]]
Results$Mean.Score <- Results$Score
return(Results)
}
indx.rm = numeric()
if(k > 1){
for(i in 1:(k-1)){
for(j in (i+1):k){
match = length(intersect(Results[[i]]$B, Results[[j]]$B))*length(intersect(Results[[i]]$I, Results[[j]]$I))/(min(length(Results[[i]]$B)*length(Results[[i]]$I),length(Results[[j]]$B)*length(Results[[j]]$I)))
pot = c(i,j)
indx = pot[which.min(c(Results[[i]]$Score, Results[[j]]$Score))]
if(match > beta-0.00001){indx.rm = c(indx.rm,indx)}
}
}
}
Results[indx.rm] = NULL
temp = 0
for(i in 1:length(Results)){
temp = temp + Results[[i]]$Score
}
Mean.score = temp/length(Results)
return(list(Communities = Results, Mean.Score = Mean.score))
}
if(length(Results.temp) < 1){return("No Community Found")}
Scores = rep(0, length(Results.temp))
for(i in 1:length(Results.temp)){
if(length(Results.temp[[i]]$B) == 0){Scores[i] = -1000}
if(length(Results.temp[[i]]$B) > 0){
Scores[i] = Results.temp[[i]]$Score
}
}
Scores = round(Scores, 2)
#keep only unique communities with score greater than threshold
indx = which(!duplicated(Scores) == TRUE)
indx.2 = which(Scores > min.score)
Results2 = Results.temp[intersect(indx, indx.2)]
if(length(Results2) == 0){
Results = NULL
return(Object = NULL)
}
if(length(Results2) > 0){
betas = seq(0.01, 1, by = 0.01)
Results3 = list()
Number.Communities = rep(0, length(betas))
Mean.Score = rep(0, length(betas))
for(i in 1:length(betas)){
temp = cleanup(Results2, betas[i])
Results3[[i]] = list(Beta = betas[i], Communities = temp$Communities)
Mean.Score[i] = temp$Mean.Score
Number.Communities[i] = length(temp$Communities)
}
}
Z = data.frame(Beta = betas, Mean.Score = Mean.Score, Number.Communities = Number.Communities)
Object = list(Community.List = Results3, Diagnostics = Z)
Z
community.object <- multilayer.extraction(adjacency = network, seed = 123, min.score = 0, prop.sample = .10)
plot(community.object, main = "Diagnostic Plot AU_CS")
cleanup = function(Results, beta){
k = length(Results)
if(k < 2){
Results <- Results[[1]]
Results$Mean.Score <- Results$Score
Results$Number.Communities <- 1
return(Results)
}
indx.rm = numeric()
if(k > 1){
for(i in 1:(k-1)){
for(j in (i+1):k){
match = length(intersect(Results[[i]]$B, Results[[j]]$B))*length(intersect(Results[[i]]$I, Results[[j]]$I))/(min(length(Results[[i]]$B)*length(Results[[i]]$I),length(Results[[j]]$B)*length(Results[[j]]$I)))
pot = c(i,j)
indx = pot[which.min(c(Results[[i]]$Score, Results[[j]]$Score))]
if(match > beta-0.00001){indx.rm = c(indx.rm,indx)}
}
}
}
Results[indx.rm] = NULL
temp = 0
for(i in 1:length(Results)){
temp = temp + Results[[i]]$Score
}
Mean.score = temp/length(Results)
return(list(Communities = Results, Mean.Score = Mean.score))
}
object <- refine(community.object, k = 6, m = m, n = n)
object <- refine(community.object, k = 1, m = m, n = n)
vertex.set <- 1:10
layer.set <- 1:2
n = 61
if(length(layer.set) < 1 || length(vertex.set) < 1){
return(obs.score = 0)
}
if(length(layer.set) == 1){
super.mod <- mod.matrix[[layer.set]] #just a single igraph object
}
if(length(layer.set) > 1){
#merge the modularity graphs
super.mod <- graph.empty(n = n, directed = FALSE)
for(j in layer.set){
super.mod <- union(super.mod, mod.matrix[[j]]) #take union of all networks in the layer.set
}
}
super.mod
super.mod.subgraph <- induced_subgraph(super.mod, v = vertex.set)
#
edge.weights <- as.data.frame(get.edge.attribute(super.mod.subgraph))
edge.weights[is.na(edge.weights)] <- 0
edge.weights
modularity.score <- colSums(edge.weights)
modularity.score
modularity.score <- rowSums(edge.weights)
modularity.score
modularity.score[which(modularity.score < 0)] <- 0
modularity.score
tot.mod <- sum(modularity.score)
obs.score <- (tot.mod)^2 / (n^2*choose(length(vertex.set), 2)*(length(layer.set)))
obs.score
choose(length(vertex.set), 2)
source("score_temp.R")
community.object <- multilayer.extraction(adjacency = network, seed = 123, min.score = 0, prop.sample = .10)
score = function(mod.matrix, vertex.set, layer.set, n){
if(length(layer.set) < 1 || length(vertex.set) < 1){
return(obs.score = 0)
}
if(length(layer.set) == 1){
super.mod <- mod.matrix[[layer.set]] #just a single igraph object
}
if(length(layer.set) > 1){
#merge the modularity graphs
super.mod <- graph.empty(n = n, directed = FALSE)
for(j in layer.set){
super.mod <- union(super.mod, mod.matrix[[j]]) #take union of all networks in the layer.set
}
}
#take sub-graph of the modularity matrix
super.mod.subgraph <- induced_subgraph(super.mod, v = vertex.set)
#get the edge weights all together
edge.weights <- as.data.frame(get.edge.attribute(super.mod.subgraph))
edge.weights[is.na(edge.weights)] <- 0
modularity.score <- rowSums(edge.weights) #sum across vertices first
modularity.score[which(modularity.score < 0)] <- 0 #only keep positive values
tot.mod <- sum(modularity.score)
obs.score <- (tot.mod)^2 / (n^2*choose(length(vertex.set), 2)*(length(layer.set)))
return(obs.score)
}
community.object <- multilayer.extraction(adjacency = network, seed = 123, min.score = 0, prop.sample = .10)
plot(community.object)
dev.new()
plot(community.object)
mod.matrix <- modularity.matrix(adjacency)
mod.matrix
for(i in 1:m){
if(i == 1){
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- initialization(graph, prop.sample, m, n)
}else{
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- Map(c, initial.set, initialization(graph, prop.sample, m, n))
}
}
print(paste("Search Stage"))
cat(paste("Searching over", length(initial.set[[1]]), "seed sets \n"))
Results.temp <- list()
K <- length(initial.set[[1]])
#detectCores detects the number of cores available on your instance
registerDoParallel(detectCores())
Results.temp <- foreach(i=1:K,.packages="MultilayerExtraction") %dopar% {
starter <- list()
starter$vertex.set <- as.numeric(initial.set$vertex.set[[i]])
starter$layer.set <- as.numeric(initial.set$layer.set[[i]])
single.swap(starter, adjacency, mod.matrix, m, n)
}
Results.temp
if(length(Results.temp) < 1){return("No Community Found")}
Scores = rep(0, length(Results.temp))
Scores
for(i in 1:length(Results.temp)){
if(length(Results.temp[[i]]$B) == 0){Scores[i] = -1000}
if(length(Results.temp[[i]]$B) > 0){
Scores[i] = Results.temp[[i]]$Score
}
}
Scores
Scores = round(Scores, 2)
Scores
indx = which(!duplicated(Scores) == TRUE)
indx
for(i in 1:length(Results.temp)){
if(length(Results.temp[[i]]$B) == 0){Scores[i] = -1000}
if(length(Results.temp[[i]]$B) > 0){
Scores[i] = Results.temp[[i]]$Score
}
}
Scores = round(Scores, 5)
Scores
indx = which(!duplicated(Scores) == TRUE)
indx
indx.2 = which(Scores > min.score)
indx.2
Results2 = Results.temp[intersect(indx, indx.2)]
Results2
if(length(Results2) == 0){
Results = NULL
return(Object = NULL)
}
if(length(Results2) > 0){
betas = seq(0.01, 1, by = 0.01)
Results3 = list()
Number.Communities = rep(0, length(betas))
Mean.Score = rep(0, length(betas))
for(i in 1:length(betas)){
temp = cleanup(Results2, betas[i])
Results3[[i]] = list(Beta = betas[i], Communities = temp$Communities)
Mean.Score[i] = temp$Mean.Score
Number.Communities[i] = length(temp$Communities)
}
}
Mean.Score
Z = data.frame(Beta = betas, Mean.Score = Mean.Score, Number.Communities = Number.Communities)
Object = list(Community.List = Results3, Diagnostics = Z)
Z
multilayer.extraction = function(adjacency, seed = 123, min.score = 0, prop.sample = 0.05, directed = c(FALSE, TRUE)){
#adjacency should be an edgelist with three columns - node1, node2, layer
#layer should be numbered with integers
#each network has the same number of nodes
#nodes are indexed by the integers starting from 1
m <- max(adjacency[, 3]) #max of the layer index
n <- length(unique(c(adjacency[, 1], adjacency[, 2])))
directed <- directed[1]
#Calculate the modularity matrix
print(paste("Estimation Stage"))
mod.matrix <- modularity.matrix(adjacency)
#Initialize the communities #TODO
print(paste("Initialization Stage"))
for(i in 1:m){
if(i == 1){
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- initialization(graph, prop.sample, m, n)
}else{
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- Map(c, initial.set, initialization(graph, prop.sample, m, n))
}
}
#Search Across Initial sets
print(paste("Search Stage"))
cat(paste("Searching over", length(initial.set[[1]]), "seed sets \n"))
Results.temp <- list()
K <- length(initial.set[[1]])
#detectCores detects the number of cores available on your instance
registerDoParallel(detectCores())
Results.temp <- foreach(i=1:K,.packages="MultilayerExtraction") %dopar% {
starter <- list()
starter$vertex.set <- as.numeric(initial.set$vertex.set[[i]])
starter$layer.set <- as.numeric(initial.set$layer.set[[i]])
single.swap(starter, adjacency, mod.matrix, m, n)
}
#Cleanup the results: Keep the unique communities
print(paste("Cleaning Stage"))
if(length(Results.temp) < 1){return("No Community Found")}
Scores = rep(0, length(Results.temp))
for(i in 1:length(Results.temp)){
if(length(Results.temp[[i]]$B) == 0){Scores[i] = -1000}
if(length(Results.temp[[i]]$B) > 0){
Scores[i] = Results.temp[[i]]$Score
}
}
Scores = round(Scores, 5)
#keep only unique communities with score greater than threshold
indx = which(!duplicated(Scores) == TRUE)
indx.2 = which(Scores > min.score)
Results2 = Results.temp[intersect(indx, indx.2)]
if(length(Results2) == 0){
Results = NULL
return(Object = NULL)
}
if(length(Results2) > 0){
betas = seq(0.01, 1, by = 0.01)
Results3 = list()
Number.Communities = rep(0, length(betas))
Mean.Score = rep(0, length(betas))
for(i in 1:length(betas)){
temp = cleanup(Results2, betas[i])
Results3[[i]] = list(Beta = betas[i], Communities = temp$Communities)
Mean.Score[i] = temp$Mean.Score
Number.Communities[i] = length(temp$Communities)
}
}
Z = data.frame(Beta = betas, Mean.Score = Mean.Score, Number.Communities = Number.Communities)
Object = list(Community.List = Results3, Diagnostics = Z)
class(Object) = "MultilayerCommunity"
return(Object)
}
community.object <- multilayer.extraction(adjacency = network, seed = 123, min.score = 0, prop.sample = .10)
dev.new()
plot(community.object)
object <- refine(community.object, k = 6, m = m, n = n)
num.layers <- colSums(object$Layers > 0)
num.vertices <- colSums(object$Vertices)
print(data.frame(Layers = num.layers, Vertices = num.vertices))
indx1 <- which(object$Vertices[,1] == 1)
indx2 <- setdiff(which(object$Vertices[,2] == 1), indx1)
indx3 <- setdiff(setdiff(which(object$Vertices[,3] == 1), indx1), indx2)
indx4 <- setdiff(setdiff(setdiff(which(object$Vertices[,4] == 1), indx3),       indx2),indx1)
indx5 <- setdiff(setdiff(setdiff(setdiff(which(object$Vertices[,5] == 1), indx4), indx3),indx2), indx1)
indx6 <- setdiff(setdiff(setdiff(setdiff(setdiff(which(object$Vertices[,6] == 1),
indx1), indx2), indx3), indx4), indx5)
none2 <- setdiff(1:61, union(union(union(indx1, indx2), union(indx3, indx4)), union(indx5, indx6)))
image(object$Layers, main = "Layers in Communities")
image(object$Layers, main = "Layers in Communities")
re.order <- c(indx1, indx2, indx3, indx4, indx5, indx6, none2)
image(AU_CS$coauthor[re.order,re.order], main = "Co-Author")
image(AU_CS$leisure[re.order,re.order], main = "Leisure")
image(AU_CS$work[re.order,re.order], main = "Work")
image(AU_CS$lunch[re.order,re.order], main = "Lunch")
image(AU_CS$facebook[re.order,re.order], main = "Facebook")
dev.new()
image(AU_CS$coauthor[re.order,re.order], main = "Co-Author")
image(AU_CS$leisure[re.order,re.order], main = "Leisure")
image(AU_CS$work[re.order,re.order], main = "Work")
image(AU_CS$lunch[re.order,re.order], main = "Lunch")
image(AU_CS$facebook[re.order,re.order], main = "Facebook")
rm(list = ls())
install.packages("devtools")
library(devtools, quietly = TRUE)
#install and load MultilayerExtraction
devtools::install_github('jdwilson4/multilayer_extraction')
library(MultilayerExtraction, quietly = TRUE)
library(Matrix, quietly = TRUE)
library(igraph, quietly = TRUE)
library(foreach, quietly = TRUE)
library(doParallel, quietly = TRUE)
library(parallel, quietly = TRUE)
data("AU_CS")
#number of nodes and layers
n <- dim(AU_CS[[1]])[1]
m <- length(AU_CS)
#types of relationships among attendees
relationship.names <- names(AU_CS)
print(relationship.names)
image(AU_CS$coauthor, main = "Co-author")
image(AU_CS$facebook, main = "Leisure")
image(AU_CS$leisure, main = "Facebook")
Mean.Degree <- sapply(AU_CS, mean) * n
print(Mean.Degree)
network <- adjacency.to.edgelist(AU_CS)
community.object <- multilayer.extraction(adjacency = network, seed = 123, min.score = 0, prop.sample = .10)
#plot the number of communities across overlap parameter beta
plot(community.object, main = "Diagnostic Plot AU_CS")
rm(list = ls())
devtools::install_github('jdwilson4/multilayer_extraction')
library(MultilayerExtraction, quietly = TRUE)
?modularity_matrix
?adjacency_to_edgelist
data("AU_CS")
#number of nodes and layers
n <- dim(AU_CS[[1]])[1]
m <- length(AU_CS)
network <- adjacency.to.edgelist(AU_CS)
setwd("/Users/jdwilson4/Dropbox/Github/MultilayerExtraction")
document()
devtools::install_github('jdwilson4/multilayer_extraction')
library(MultilayerExtraction, quietly = TRUE)
devtools::install_github('jdwilson4/multilayer_extraction')
devtools::install_github('jdwilson4/multilayer_extraction')
document()
document()
devtools::install_github('jdwilson4/multilayer_extraction')
devtools::install_github('jdwilson4/multilayer_extraction')
document
document()
document()
devtools::install_github('jdwilson4/multilayer_extraction')
library(MultilayerExtraction, quietly = TRUE)
network <- adjacency.to.edgelist(AU_CS)
community.object <- multilayer.extraction(adjacency = network, seed = 123, min.score = 0, prop.sample = .10)
network
multilayer.modularity.matrix()
?multilayer.modularity.matrix
rm(list = ls())
devtools::install_github('jdwilson4/multilayer_extraction')
library(MultilayerExtraction, quietly = TRUE)
library(Matrix, quietly = TRUE)
library(igraph, quietly = TRUE)
library(foreach, quietly = TRUE)
library(doParallel, quietly = TRUE)
library(parallel, quietly = TRUE)
data("AU_CS")
#number of nodes and layers
n <- dim(AU_CS[[1]])[1]
m <- length(AU_CS)
network <- adjacency.to.edgelist(AU_CS)
?multilayer.extraction
community.object <- multilayer.extraction(adjacency = network, seed = 123, min.score = 0, prop.sample = .10)
class(network)
adjacency = network
m
n
for(i in 1:m){
if(i == 1){
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- initialization(graph, prop.sample, m, n)
}else{
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- Map(c, initial.set, initialization(graph, prop.sample, m, n))
}
}
directed = FALSE
for(i in 1:m){
if(i == 1){
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- initialization(graph, prop.sample, m, n)
}else{
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- Map(c, initial.set, initialization(graph, prop.sample, m, n))
}
}
prop.sample = 0.01
for(i in 1:m){
if(i == 1){
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- initialization(graph, prop.sample, m, n)
}else{
graph <- graph_from_edgelist(as.matrix(subset(as.data.frame(adjacency), layer == i))[, 1:2],
directed = directed)
initial.set <- Map(c, initial.set, initialization(graph, prop.sample, m, n))
}
}
mod.matrix <- multilayer.modularity.matrix(adjacency)
class(mod.matrix)
mod.matrix
community.object <- multilayer.extraction(adjacency = network, seed = 123, min.score = 0, prop.sample = .10)
library(MultilayerExtraction)
?multilayer.extraction