Network Analysis

Input

Load the gene coordinate annotations for protein expressions and genomic CNA.

library(Biobase)
yinfo <- pData(featureData(readRDS(file = "data/prot-expression-set.rds")))
xinfo <- pData(featureData(readRDS(file = "data/cna-expression-set.rds")))

Load the Boot.Vote CNA-protein network.

library(Matrix)
scggm <- readRDS(file = "model-fits/scggm_boot.rds")
net <- scggm[c("yy","xy")]
rownames(net$xy) <- xinfo$id; colnames(net$xy) <- yinfo$id;
rownames(net$yy) <- yinfo$id; colnames(net$yy) <- yinfo$id;

Load the degree distributions for the ensemble of bootstrapped networks.

bdeg <- scggm[c("deg_xy", "deg_yy")]
names(bdeg) <- c("xy", "yy")
colnames(bdeg$xy) <- xinfo$id; colnames(bdeg$yy) <- yinfo$id;

Load the Gene Ontology mappings for enrichment analysis.

go2eg <- readRDS("data/prot-go-bp-to-entrez-gene-list.rds")
library(AnnotationDbi)
#human readable 
process_alias <- Term(names(go2eg))

Map Annotations

Convert the Boot.Vote network into an igraph object and map the annotations onto the network.

library(spacemap)
ig <- spacemap::adj2igraph(yy = net$yy, xy = net$xy, yinfo = yinfo, xinfo = xinfo)

If we query the attribute names of the nodes in the graph, we notice that the columns of xinfo and yinfo have been applied.

vertex_attr_names(graph = ig)

## [1] "name"        "alias"       "chr"         "start"       "end"        
## [6] "strand"      "description" "levels"

The igraph package has a number of ways to access the annotation information. For example, if we wish to confirm the chromosome location of GRB7, we can easily query:

vertex_attr(graph = ig, name = "chr", index = V(ig)[alias %in% "GRB7"])

## [1] "chr17"

Hub Analysis

We first prioritize the CNA- and protein- hubs. If the bdeg argument is specified, then we rank the hubs according to the average degree rank. Accordingly, the most highly ranked hubs will have the most consistently high degree across network ensemble.

To rank the protein nodes, use the rankHub command and simply specify the level = "y" argument.

ig <- rankHub(ig = ig, bdeg = bdeg$yy, level = "y")

To rank the CNA nodes, specify the level = "x" argument.

ig <- rankHub(ig = ig, bdeg = bdeg$xy, level = "x")

Alternatively, if the bdeg argument is not available, the ranking is according to degree in the final network (coded by “ig”).

tmp <- rankHub(ig = ig,  level = "x")

ig <- cisTrans(ig = ig, level = "x-y", cw = 2e6)

xhubs <- reportHubs(ig, top = 40, level = "x")

vidx <- V(ig)[match(rownames(xhubs), V(ig)$name)]
m <- vertex_attr(graph = ig, name = "mean_rank_hub", index = vidx)
s <- vertex_attr(graph = ig, name = "sd_rank_hub", index = vidx)
st <- vertex_attr(graph = ig, name = "start", index = vidx)
end <- vertex_attr(graph = ig, name = "end", index = vidx)
rnk <- paste0(round(m, 2), " (", round(s,2), ")")
xhubs$`cis genes`
tosupp <- cbind(xhubs$hub, st, end, xhubs$`# cis/ # trans`, xhubs$`Potential # cis`, xhubs$`cis genes`, rnk)
library(xtable)
print(xtable(tosupp), include.rownames = FALSE)

yhubs <- reportHubs(ig, top = 10, level = "y")

hgp <- xHubEnrich(ig = ig, go2eg = go2eg)

Module Analysis

There are many criteria to define modules of a network. This toolkit allows users to import the module membership information by themselves (see mods argument in modEnrich).

In the spaceMap publication, we use the edge-betweenness algorithm in igraph.

library(igraph)
mods <- cluster_edge_betweenness(ig)

The main goal of module analysis is identifying modules that are functionally enriched. The modEnrich function will test for significantly over-represented GO-terms (or any other valid functional grouping) within a module using hyper-geometric tests.

In the current example, only the protein nodes have functional mapping and we specify this through the levels = "y" argument. If both predictors and responses have functional mapping in the go2eg argument, then we can specify levels = c("y","x"). Other arguments are available to control the enrichment testing; see modEnrich for more details.

outmod <- modEnrich(ig = ig, mods = mods, levels = "y", go2eg = go2eg, process_alias = process_alias, prefix = "C")

The output of the module analysis is a list of 3 objects.

names(outmod)

## [1] "ig"   "etab" "eraw"

• ig is the input igraph network object updated with a “process_id” attribute for nodes affiliated with a significant GO-term. The “process_id” and “module” attributes together are useful for visualizing which nodes are enriched for a specific biological function.

• etab is the polished module enrichment table to report significant GO terms, the representation of the GO term in the module relative to the size of the GO term, and which CNA hubs belong to the module. The top ten hits as in Table S.9 of the spaceMap publication’s supplementary materials are as follows:

outmod$etab[1:20,]

• eraw contains details for each (module, GO-term) pair that was subjected to the hyper-geometric test. This output gives the user more control by reporting all tests, the relative over-representation of a GO-term in that module, the raw P-value, and the adjusted P-value.

outmod$eraw[1:5,]

Export for Visualization

There are many tools for network visualization. In the publication of spaceMap, we exported the annotated igraph network object ig to the “graphml” format, which maintained all the attributes associated with nodes when read into Cytoscape for visualization.

filename <- "~/repos/neta-bcpls/neta/spacemap-scggm-boot-vote.graphml"
#delete nodes without edges from the graph object
vis <- delete_vertices(graph = outmod$ig, v = V(outmod$ig)[igraph::degree(outmod$ig) == 0])
igraph::write_graph(graph = vis, file = filename, format = "graphml")

Here we list all the attributes associated with the nodes that can be used in tandem with Cytoscape’s filtering functions to identify specific nodes of interest.

vertex_attr_names(outmod$ig)

##  [1] "name"             "alias"            "chr"             
##  [4] "start"            "end"              "strand"          
##  [7] "description"      "levels"           "rank_hub"        
## [10] "mean_rank_hub"    "sd_rank_hub"      "ncis"            
## [13] "ntrans"           "regulates_in_cis" "potential_cis"   
## [16] "module"           "C"

We describe some of the most useful attributes for visualization:

• ‘name’: the unique node ID
• ‘alias’: the node alias (e.g. gene symbol ERBB2)
• ‘chr,start,end,strand’: the gene coordinates of nodes
• ‘description’: any note (e.g. breast cancer oncogene)
• ‘levels’: indicates whether the node belongs to predictors “x” or responses “y”
• ‘rank_hub’: the rank of the hub within its level (e.g. a value of “1” for a node of level “x” means that it is the most consistently high degree \(x\) node in the network. )
• ‘regulates_in_cis’: list of genes regulated in cis
• ‘module’: the module ID that the node belongs to.
• ‘process_id’: the significant GO-term(s) associated with the node.

Also the edge attributes are exported to ‘graphml’ format.

edge_attr_names(outmod$ig)

## [1] "levels"    "cis_trans"

• ‘levels’ indicates whether an edge is \(x-y\) or \(y-y\).
• ‘cis_trans’ indicates whether an edge is regulated in cis or in trans.

Summary

This toolkit is useful for analyzing and summarizing the output of the spaceMap model fit in the context of integrative genomics. The biological context mapped onto the network object can easily be exported to existing network visualization tools like Cytoscape.