Background

Condition-specific network activation is characteristic for cellular systems and other real-world interaction networks. If measurements of network states are available across a versatile set of conditions or time points, it becomes possible to construct a global view of network activation patterns. Different parts of the network respond to different conditions, and in different ways. Systematic, data-driven identification of these responses will help to obtain a holistic view of network activity [1-2]. This package provides robust probabilistic algorithms for functional network analysis [1, 3].

The methods are based on nonparametric probabilistic modeling and variational learning, and provide general exploratory tools to investigate the structure

NetResponse detects and characterizes subnetworks that exhibit context-specific activation patterns across versatile collections of functional measurements, such as gene expression data. The implementations are partially based on the agglomerative independent variable group analysis (AIVGA) and variational Dirichlet process Gaussian mixture models (Kurihara et al. 2007). The tools are particularly useful for global exploratory analysis of genome-wide interaction networks and versatile collections of gene expression data, and in general to discover subnetworks with alternative states.

Usage examples

Examples on running NetResponse algorithm and visualizing the results. The algorithm combines network and functional information to detect coherent subnetworks that reveal distinct activation modes across conditions.

Generate simulated data.

library(netresponse)
res <- generate.toydata(Dim = 3, Nc = 3, Ns = 200, sd0 = 3, rgam.shape = 1, rgam.scale = 1)

D <- res$data
component.means <- res$means
component.sds   <- res$sds
sample2comp     <- res$sample2comp

# Use fully connected network
network <- matrix(rep(1, 9), nrow = 3) 

Fit NetResponse model.

# Various network formats are supported, see help(detect.responses) for
# details. With large data sets, consider the 'speedup' option.
set.seed(4243)
res <- detect.responses(D, network, mixture.method = "vdp", pca.basis = TRUE)

# List subnets (each is a list of nodes)
subnet.id <- names(get.subnets(res))[[1]]

library(ggplot2)
vis <- plot_responses(res, subnet.id, plot_mode = "pca")

# Modify the resulting ggplot2 object to enhance visualization
p <- vis$p # Pick the ggplot2 object from results
p <- p + geom_point(size = 3) # Modify point size
print(p) 

vis <- plot_responses(res, subnet.id, plot_mode = "network")

vis <- plot_responses(res, subnet.id, plot_mode = "heatmap")

vis <- plot_responses(res, subnet.id, plot_mode = "boxplot_data")

plot_scale(vis$breaks, vis$palette, two.sided = TRUE)

subnet.id <- 'Subnet-1'

# Sample - response probabilities (soft cluster assignment)
response.probs <- sample2response(res, subnet.id)
tail(round(response.probs, 6))

# Sample - response hard assignments
hard.clusters <- response2sample(res, subnet.id)
print(hard.clusters)

params <- get.model.parameters(res, subnet.id) 
names(params)

Nonparametric Gaussian mixture models

Nonparametric Gaussian mixtures with variational Dirichlet processes based on implementations by Kurihara et al. (2007) and Honkela et al..

# Generate 2-dimensional simulated data with 3 clusters
res <- generate.toydata(Dim = 2, Nc = 3, Ns = 200, sd0 = 3, rgam.shape = 1, rgam.scale = 1)

D <- res$data
real.means <- res$means
real.sds   <- res$sds
real.sample2comp     <- res$sample2comp

# Infinite Gaussian mixture model with       
# Variational Dirichlet Process approximation       
mixt <- vdp.mixt( D )
            
# Centroids of the detected Gaussian components       
estimated.means <- mixt$posterior$centroids

# The colors denote the known clusters
# The blue ball denotes the original (known) cluster centroids and
# the triangle denotes the estimated centroids
plot(D, col = real.sample2comp, pch = 1)
points(real.means, col = "blue", pch = 16, cex = 2)
points(estimated.means, col = "blue", pch = 17, cex = 2)
            
# Hard mixture component assignment for each sample
estimated.sample2comp <- apply(mixt$posterior$qOFz, 1, which.max)

# Compare known and estimated mixture components
# (note that cluster indices may have switched due to unidentifiability)
# nearly all samples have one-to-one match between the real and estimated 
# clusters
head(table(estimated.sample2comp, real.sample2comp))

citation("netresponse")

sessionInfo()