Algorithms and Analyses of Molecular Interaction Networks

Ananth Grama

Center for Science of Information, Purdue University

With significant contributions from Mehmet Koyuturk, Jayesh Pandey, Young Kim, Shahin Mohammadi, and Shankar Subramaniam

Thanks to the Indian Institute of Science and the US National Science Foundation.

Function & Topology in Molecular Networks

How does function relate to network topology?

Prior Work on Topology and Function

• Conservation [ISMB 04/Bioinf. 04]

• Alignment [RECOMB 05/JCB 06]

• Modularity [RECOMB 06/JCB 07]

• Inference [Bioinf. 06]

• Pathway Annotation [ISMB 07/Bioinf. 07, PSB 08]

• Network Abstractions/ Annotations [ECCB 08/ Bioinf. 08]

• Modularity and Domain Interactions [APBC 10/ BMC Bioinf. 10]

• Pathway Interaction Maps [PSB 12]

• Pathway Inference [ISMB 12]

Evolution of Molecular Interactions

• “Evolution thinks modular”[Vespignani, Nature Gen., 2003]

• Cooperative tasks require all participating units

– Selective pressure on preserving interactions & interacting proteins

• Proteins organized in cohesive patterns are highly conserved

[Wuchty et al., Nature Gen., 2003]

– Functional modules are likely to be consistently conserved

• Orthologs of interacting proteins are likely to interact [Wagner,

Mol. Bio. Evol., 2001]

– Conservation of interactions may provide clues on conservation of

function

• Interacting proteins follow similar evolutionary trajectories

[Pellegrini et al., PNAS, 1999]

Computational Analysis of Biological Networks

• Clustering

– Interaction network: Proteins in functional modules densely interact with each other

– Gene expression: Genes coding cooperating proteins are likely to be coregulated

– Phylogenetic profiles: Interacting proteins are likely to have co-evolved

• Graph Mining

– Common topological motifs and frequent interaction patterns reveal conserved modularity

• Graph Alignment

– Conservation/divergence of pathways, complexes, and functional modules

Frequent Interaction Patterns: Computational Problem

• Given a set of proteins V a set of interactions E, and a manyto-many mapping from V to a set of ortholog groups L = {l_{1}, l_{2}, ..., l_{n}}, the corresponding interaction network is a labeled graph G = (V, E, L).

– v ∈ V (G) is associated with a set of ortholog groups L(v) ⊆ L.

– uv ∈ E(G) represents an interaction between u and v.

• S is a sub-network of G, i.e., S ⊑ G if there is an injective mapping ϕ : V (S) → V (G) such that for all v ∈ V (S), L(v) ⊆ L(ϕ(v)) and for all uv ∈ E(S), ϕ(u)ϕ(v) ∈ E(G).

• Maximal frequent sub-network discovery

– Instance: A set of interaction networks G = {G_{1} = (V_{1}, E_{1}, L), G_{2} = (V_{2}, E_{2}, L), ..., G_{m} = (V_{m}, E_{m}, L)}, each belonging to a different organism, and a frequency threshold σ^{∗}.

– Problem: Let H(S) = {G_{i} : S ⊑ G_{i}} be the occurrence set of graph S. Find all connected subgraphs S such that |H(S)| ≥ σ^{∗}, i.e., S is a frequent subgraph in G and for all S′ ⊐ S, H(S) ≠ H(S′), i.e., S is maximal.

Ortholog Contraction

• Contract orthologous nodes into a single node

• No subgraph isomorphism

– Graphs are uniquely identified by their edge sets

• Frequent sub-networks are preserved ⇒ No information loss

– Sub-networks that are frequent in general graphs are also frequent in their ortholog-contracted representation

• Discovered frequent sub-networks are still biologically interpretable!

– Interaction between proteins becomes interaction between ortholog groups

– Ortholog-contraction may be thought of as going back in evolutionary history (to what point?)

Ortholog Contraction in PPI Networks

• Interaction between proteins → Interaction between ortholog groups or protein families

PSfrag replacements Rrp4 Rrp43 Mtr3 Ski6 Csr4

KOG3013 KOG1613 KOG1068 KOG3409

Preservation of Sub-networks

Theorem: Let G̃ be the ortholog-contracted graph obtained by contracting the orthologous nodes of network G. Then, if S is a subgraph of G, S̃ is a subgraph of G̃.

Corollary: The ortholog-contracted representation of any frequent sub-network is also frequent in the set of orthologcontracted graphs.

PSfrag replacements

G G̃

Results: Mining PPI Networks

• PPI networks for 9 eukaryotic organisms derived from BIND and DIP

– A. thaliania, O. sativa, S. cerevisiae, C. elegans, D. melanogaster, H. sapiens, B. taurus, M. musculus, R. norvegicus

– # of proteins ranges from 288 (Arabidopsis) to 8577 (fruit fly)

– # of interactions ranges from 340 (rice) to 28829 (fruit fly)

• Ortholog contraction

– Group proteins according to existing COG ortholog clusters

– Merge Homologene groups into COG clusters

– Cluster remaining proteins via BLASTCLUST

– Ortholog-contracted fruit fly network contains 11088 interactions between 2849 ortholog groups

• Mule is available at

http://www.cs.purdue.edu/homes/koyuturk/mule/

Frequent Protein Interaction Patterns

PSfrag replacements

Pattern C. elegans S. cerevisiae D. melanogaster

Small nuclear ribonucleoprotein complex (p < 2e − 43)

Frequent Protein Interaction Patterns

PSfrag replacements

Pattern H. sapiens S. cerevisiae

B. taurus

Actin-related protein Arp2/3 complex (p < 9e − 11)

– Typeset by FoilTEX – 11

Frequent Protein Interaction Patterns

PSfrag replacements

Pattern H. sapiens S. cerevisiae D. melanogaster

Endosomal sorting (p < 1e − 78)

– Typeset by FoilTEX – 12

Modular Phylogenetics

• Top eight groups of three organisms that contain most frequent

connected sub-networks and interactions

# frequent  # frequent
Organism set  sub-networks interactions
C. elegans, D. melanogaster, H. sapiens  8  134
S. cerevisiae, D. melanogaster, H. sapiens  20  126
D. melanogaster, H. sapiens, M. musculus  17  86
S. cerevisiae, C. elegans, D. melanogaster  15  77
S. cerevisiae, C. elegans, H. sapiens  6  50
S. cerevisiae, H. sapiens, M. musculus  10  26
C. elegans, H. sapiens, M. musculus  5  23
H. sapiens, M. musculus, R. norvegicus  10  23

Runtime Characteristics

Comparison with isomorphism-based algorithms

FSG [Kuramochi & Karypis, ICDM, 2001], gSpan [Yan & Han, KDD, 2003]

FSG  Mule
Minimum  Runtime  Largest  Number of  Runtime  Largest  Number of
Dataset  Support (%)  (secs.)  pattern  patterns  (secs.)  pattern  patterns
20  0.2  9  12  0.01  9  12
16  0.7  10  14  0.01  10  14
Glutamate  12  5.1  13  39  0.10  13  39
10  22.7  16  34  0.29  15  34
8  138.9  16  56  0.99  15  56
24  0.1  8  11  0.01  8  11
20  1.5  11  15  0.02  11  15
Alanine  16  4.0  12  21  0.06  12  21
12  112.7  17  25  1.06  16  25
10  215.1  17  34  1.72  16  34

Extraction of contracted patterns

Glutamate metabolism, σ = 8%  Alanine metabolism, σ = 10%
Size of  Extraction time  Size of  Size of  Extraction time  Size of
contracted  (secs.)  extracted  contracted  (secs.)  extracted
pattern  FSG gSpan  pattern  pattern  FSG gSpan  pattern
15  10.8  1.12  16  16  54.1  10.13  17
14  12.8  2.42  16  16  24.1  3.92  16
13  1.7  0.31  13  12  0.9  0.27  12
12  0.9  0.30  12  11  0.4  0.13  11
11  0.5  0.08  11  8  0.1  0.01  8
Total number of patterns: 56  Total number of patterns: 34
Total runtime of FSG alone: 138.9 secs.  Total runtime of FSG alone :215.1 secs.
Total runtime of Mule+FSG: 0.99+100.5 secs.  Total runtime of Mule+FSG: 1.72+160.6 secs.
Total runtime of Mule+gSpan: 0.99+16.8 secs.  Total runtime of Mule+gSpan: 1.72+31.0 secs.

– Typeset by FoilTEX – 14

Pairwise Alignment of PPI Networks

• Given two PPI networks that belong to two different organisms, identify sub-networks that are similar to each other

– Biological meaning

– Mathematical modeling

• Existing algorithms

– PathBLAST aligns pathways (linear chains) to simplify the problem while maintaining biological meaning [Kelley et al., PNAS, 2004]

– NetworkBLAST compares conserved complex model with null model to identify significantly conserved subnets [Sharan et al., J. Comp. Biol., 2005]

• Our approach [Koyutürk, Kim, Topkara, Subramaniam, Szpankowski, & Grama, J. Comp. Biol., 2006]

– Guided by models of evolution

– Scores evolutionary events

– Identifies sets of proteins that induce high-scoring sub-network pairs

Evolution of PPI Networks

• Duplication/divergence models for the evolution of protein interaction networks

– Interactions of duplicated proteins are also duplicated

– Duplicated proteins rapidly lose interactions through mutations

• Allows defining and scoring evolutionary events as graphtheoretical concepts

PSfrag replacements u_{1} u′_{1} u_{2} u_{3}

Duplication Elimination Emergence

Match, Mismatch, and Duplication

• Evolutionary events as graph-theoretic concepts

– A match ∈ M corresponds to two pairs of homolog proteins from each organism such that both pairs interact in both PPIs. A match is associated with score μ.

– A mismatch ∈ N corresponds to two pairs of homolog proteins from each organism such that only one pair is interacting. A mismatch is associated with penalty ν.

– A duplication ∈ D corresponds to a pair of homolog proteins that are in the same organism. A duplication is associated with score δ.

PSfrag replacements

u_{1}
u_{2}
u_{3}
u_{4}
v_{1}
v_{2}
v_{3}
e_{1}
e_{2}
e_{3}
G:
H:

Pairwise Alignment of PPIs as an Optimization Problem

• Alignment score:

σ(A(P)) = ∑_{M∈M} μ(M) − ∑_{N∈N} ν(N) + ∑_{D∈D} δ(D)

– Matches are rewarded for conservation of interactions

– Duplications are rewarded/penalized for functional conservation/differentiation after split

– Mismatches are penalized for functional divergence (what about experimental error?)

• Scores are functions of similarity between associated proteins

• Problem: Find all protein subset pairs with significant alignment score

– High scoring protein subsets are likely to correspond to conserved modules

• A graph equivalent to BLAST

Weighted Alignment Graph

• G(V, E) : V consists of all pairs of homolog proteins v = {u ∈ U, v ∈ V }

• An edge vv′ = {uv}{u′v′} in E is a

– match edge if uu′ ∈ E and vv′ ∈ V, with weight w(vv′) = μ(uv, u′v′)

– mismatch edge if uu′ ∈ E and vv′ ∉ V or vice versa, with weight w(vv′) = −ν(uv, u′v′)

– duplication edge if S(u, u′) > 0 or S(v, v′) > 0, with weight w(vv′) = δ(u, u′) or w(vv′) = δ(v, v′)

PSfrag replacements {u_{1}, v_{1}}

{u_{2}, v_{1}}
{u_{3}, v_{3}}
{u_{4}, v_{2}}
{u_{4}, v_{4}}
μ
-ν
-δ

Maximum Weight Induced Subgraph Problem

• Definition: (MaWISh)

– Given graph G(V, E) and a constant ϵ, find Ṽ ∈ V such that ∑_{v,u∈Ṽ} w(vu) ≥ ϵ.

– NP-complete

• Theorem: (MaWISh ≡ Pairwise alignment)

– If Ṽ is a solution for the MaWISh problem on G(V, E), then P = {Ũ, Ṽ } induces an alignment A(P) with σ(A) ≥ ϵ , where Ṽ = Ũ × Ṽ.

• Solution: Local graph partitioning

– Greedy graph growing + iterative refinement

– Linear-time heuristic

• Source code available at

http://www.cs.purdue.edu/homes/koyuturk/mawish/

Alignment of Yeast and Fruit Fly PPI Networks

Rank Score z-score # Proteins # Matches  # Mismatches  # Dups.
1  15.97  6.6  18 (16, 5)  28  6  (4, 0)

protein amino acid phosphorylation (69%) JAK-STAT cascade (40%)

2 13.93 3.7 13 (8, 7) 25 7 (3, 1)

endocytosis (50%) / calcium-mediated signaling (50%)

5 8.22 13.5 9 (5, 3) 19 11 (1, 0)

invasive growth (sensu Saccharomyces) (100%) oxygen and reactive oxygen species metabolism (33%)

6 8.05 7.6 8 (5, 3) 12 2 (0, 1)

ubiquitin-dependent protein catabolism (100%) mitosis (67%)

21 4.36 6.2 9 (5, 4) 18 13 (0, 5)

cytokinesis (100%, 50%)

30 3.76 39.6 6 (3, 5) 5 1 (0, 6)

DNA replication initiation (100%, 80%)

Subnets Conserved in Yeast and Fruit Fly

Proteosome regulatory particle subnet

PSfrag replacements

S.Cerevisiae D.Melanogaster

Calcium-dependent stress-activated signaling pathway

PSfrag replacements

S.Cerevisiae D.Melanogaster

Statistical Significance of Modularity

• Existing techniques

– Mostly computational (e.g., Monte-Carlo simulations)

– Compute probability that the pattern exists rather than a pattern with the property (e.g., size, density) exists

– Overestimation of significance

• Random graph models

– PPI networks generally exhibit power-law property (or exponential, geometric, etc.)

– Analysis simplified through independence assumption

– Independence assumption may cause problems for networks with arbitrary degree distribution

– P (uv ∈ E) = d_{u}d_{v}/|E|, where d_{u} is expected degree of u, but generally d_{max}^{2} > |E| for PPI networks

• Analytical techniques based on simplified models [Koyut^{}ürk,

Grama, & Szpankowski, RECOMB, 2006]

Significance of Dense Subgraphs

• A subnet of r proteins is said to be ρ-dense if F (r) ≥ ρr^{2}, where F (r) is the number of interactions between these r proteins

• What is the expected size of the largest ρ-dense subgraph in a random graph?

– Any ρ-dense subgraph with larger size is statistically significant!

• G(n, p) model

– n proteins, each interaction occurs with probability p

– Simple enough to facilitate rigorous analysis

– If we let p = d_{max}/n, largest ρ-dense subgraph in G(n, p) stochastically dominates that in a graph with arbitrary degree distribution

Largest Dense Subgraph

• Theorem: If G is a random graph with n nodes, where every

edge exists with probability p, then

lim_{n→∞} \frac{R_{ρ}}{log n} = \frac{1}{κ(p, ρ)} (pr.), (1)

where

κ(p, ρ) = ρ log \frac{ρ}{p} + (1 − ρ) log \frac{1 − ρ}{1 − p}.  (2)

More precisely,

P (R_{ρ} ≥ r_{0}) ≤ O (\frac{log n}{n^{1/κ(p,ρ)}}) ,  (3)

where

r_{0} = \frac{log n − log log n + log κ(p, ρ)}{κ(p, ρ)}  (4)

for large n.

– Typeset by FoilTEX – 25

Generalizing Results to Complex Models

• Piecewise G(n, p) model

– Few proteins with many interacting partners, many proteins with few interacting partners

– Captures the basic characteristics of PPI networks

– The size of largest dense subgraph is still proportional to log n

• More general models

– Increasing the number of pieces, we approach models with characteristic degree distributions

– Analysis of power-law graphs in progress

• Multiple networks: Conservation

– Superpose graphs based on sequence homology

Algorithms Based on Statistical Significance

• Identification of topological modules

• Use statistical significance as a stopping criterion for graph clustering heuristics

• HCS Algorithm [Hartuv & Shamir, Inf. Proc. Let., 2000]

– Find a minimum-cut bipartitioning of the network

– If any of the parts is dense enough, record it as a dense cluster of proteins

– Else, further partition them recursively

• Use statistical significance to determine whether a subgraph is sufficiently dense

– For given number of proteins and interactions between them, we can determine whether those proteins induce a significantly dense subnet

Largest Dense Subgraph for Varying Density

PSfrag replacements

Density threshold

Size of largest dense subgraph

Size of largest conserved dense subgraph PSfrag replacements

Density threshold

Size of largest dense subgraph

Size of largest conserved dense subgraph

Yeast PPI network Yeast & Fruit Fly PPI networks

Pathway Organization: Genetic Interactome

Double mutants exhibit unexpected phenotypes, as compared

to joint single mutations.

Definition 1. • Negative Interactions: more severe phenotype than expected

– Also known as aggravating or synergistic

• Positive Interactions: Less severe phenotype than expected

– Also known as alleviating or epistatic

Most commonly used:

Phenotype : Growth rate

Model : Multiplicative null model

Organization of Genetic Interactions

Definition 2. • Between-Pathway Model – Among genes participating in redundant functions

• Within-Pathway Model – Among genes with additive effect

• Indirect Effect

– Among genes with distant functions that are not directly related

Between-Pathway Model (BPM)

• Bi-cliquish structure

• Have been used to:

1. Predict co-pathway membership of gene pairs

2. Extract redundant pathways

The Genetic Landscape of a Cell

Adopted from Costanzo et al., 2010

• Baker’s yeast, Saccharomyces cerevisiae

• Synthetic Genetic Array (SGA)

• 1712 query genes

1. 1378 null alleles of nonessential genes

2. 334 hypomorphic or conditional alleles of essential genes

• 3885 array strains

Functional Annotations

• KEGG Pathway Database

• Annotations for 1026 genes in the experiment

• 96 Pathways

– 80 pathways after filtering pathways with less than 10 genes.

Local Neighborhood Similarity: A Predictor of Co-Pathway Membership

Similarity prediction methods

• Number of Shared Neighbors

• Congruence Score

• Pearson Correlation of Interaction Profiles

Both vi and vj have three shared neighbors. However, in the first case their congruence score is almost 0.6, while in the second case it is approximately 2 (assuming a graph of size 10).

Evaluating Ranking Methods

Given a pathway P_{A} and a cut size (target set) l.

Definition 3.

P − value(X = k) = P rob(k ≤ X)
= HGT (k|N, N_{A}, l)
= ∑_{x=k}^{min(N_{A},l)} \frac{C(l, x)C(N − l, N_{A} − x)}{C(N, N_{A})}

X: Random variable denoting the number of true positives in a random sample, N: Total number of gene pairs, N_{A}: Number of

gene pairs in pathway A, l: Size of target set

Minimum HyperGeometric (mHG) Score

Target size unknown:

Definition 4. The Minimum HyperGeometric (mHG) score mHG(λ) = min_{1≤l≤N}HGT (b_{l}(λ); N, N_{A}, l),

where b_{l}(λ) = ∑_{i=1}^{l} λ_{i}

λ_{i} is 1 if both of the genes in the i^{th} ranked gene pair are members of P_{A}, and 0 otherwise.

• mHG Adjusted for Multiple Comparison

– Typeset by FoilTEX – 37

Highlights

Basic Idea

Heterogeneous performance of copathway membership predictions

⟺

Existence of specific structure around enriched pathways

• Decomposing neighborhood of each pathway

• Inferring lethal crosstalk among pathways

Modified Congruence Score (MCS)

Evaluating Neighborhood Overlap of Gene Pairs With Respect to

a Given Pathway

Definition 5.

P − value(X = k_{ij}^{B}) = P rob(k_{ij}^{B} ≤ X)
= HGT (k_{ij}^{B}|n_{B}, d_{i}^{B}, d_{j}^{B})

= ∑_{x=k_{ij}^{B}}^{min(d_{i}^{B},d_{j}^{B})} \frac{C(d_{j}^{B}, x)C(n_{B} − d_{j}^{B}, d_{i}^{B} − x)}{C(n_{B}, d_{i}^{B})}

MCS is defined as −log_{10} of the P-value.

Modified Congruence Score (MCS): continued

Example 1

A sample neighborhood configuration for v_{i} and v_{j}. Here n = 15, D_{i} = 6, D_{j} = 5, n_{B} = 6, d_{i} = 3, d_{j} = 4 and k = 2.

Constructing Neighborhood Overlap Graph For a Given Pathway Pair

Definition 6. The neighborhood overlap graph (NOG) of a given pathway P_{A} with respect to pathway P_{B}, denoted by H_{A→B} = (V_{H}, E_{H}), is an unweighted, undirected graph defined over same vertices as P_{A}. In this graph, there is a link between vertices v_{i} and v_{j} if the network structure around them with respect to P_{B} is statistically significant .

Pruning neighborhood overlap graph, finding cohesive subgraphs, and identifying interaction ports

Adopted from Batagelj and Zaversnik, 2002

1. Iterative peeling of K-shells Pruning hairy components

2. Connected components in each core

3. Evaluating the significance of components

• Evaluating significance using ER random graph model

KEGG Crosstalk Map

– Typeset by FoilTEX – 43

Interaction Port Case Study

Crosstalk Between Protein Processing in ER and Proteasome

Ongoing Work: The Interaction Map of Aging

– Typeset by FoilTEX – 45

Ongoing Work

• The surprisingness of choice in networks.

• Tissue-specific alignments.

Thanks

• To the US National Science Foundation for their funding.

• To our hosts at IISc for all of their support and funding.

• To all of my collaborators.

• To you for your attention!