Accurate multiple network alignment through context-sensitive random walk

Abstract

Background

Comparative network analysis can provide an effective means of analyzing large-scale biological networks and gaining novel insights into their structure and organization. Global network alignment aims to predict the best overall mapping between a given set of biological networks, thereby identifying important similarities as well as differences among the networks. It has been shown that network alignment methods can be used to detect pathways or network modules that are conserved across different networks. Until now, a number of network alignment algorithms have been proposed based on different formulations and approaches, many of them focusing on pairwise alignment.

Results

In this work, we propose a novel multiple network alignment algorithm based on a context-sensitive random walk model. The random walker employed in the proposed algorithm switches between two different modes, namely, an individual walk on a single network and a simultaneous walk on two networks. The switching decision is made in a context-sensitive manner by examining the current neighborhood, which is effective for quantitatively estimating the degree of correspondence between nodes that belong to different networks, in a manner that sensibly integrates node similarity and topological similarity. The resulting node correspondence scores are then used to predict the maximum expected accuracy (MEA) alignment of the given networks.

Conclusions

Performance evaluation based on synthetic networks as well as real protein-protein interaction networks shows that the proposed algorithm can construct more accurate multiple network alignments compared to other leading methods.

Background

With the availability of large-scale protein-protein interactions (PPI) networks, comparative network analysis tools have been gaining increasing interest as they provide useful means of investigating the similarities and differences between different networks. As demonstrated in [1, 2], PPI networks of different species embed various conserved functional modules - such as signaling pathways and protein complexes - which can be detected through network querying [3–5] and network alignment [6–14]. Comparative network analysis methods allow us to transfer existing knowledge on well-studied organism to less-studied ones and they have the potential to detect potential functional modules conserved across different organisms and species [1, 2, 15].

There exist several different types of comparative network analysis methods, among which global network alignment methods specifically aim to predict the best overall mapping among two or more biological networks. In order to obtain biologically meaningful results, where functionally similar biomolecules across networks are accurately mapped to each other, we should consider both the molecule-level similarity between the individual molecules as well as the similarity between their interaction patterns. The former is often called the "node similarity" while the latter is typically referred to as the "topological similarity." Examination of conserved functional modules shows that many of the molecular interactions in such modules are also well conserved, clearly showing the importance of taking the topological similarity into account when comparatively analyzing biological networks. Biological networks, such as PPI networks, are typically represented as graphs, where the nodes represent individual biomolecules (e.g., proteins) and interactions (e.g., protein binding) between biomolecules are represented by edges connecting the corresponding nodes. Given these graph representations of biological networks, the network alignment problem can be formulated as an optimization problem whose goal is to find the optimal mapping - either one-to-one or many-to-many - among a set of graphs that maximizes a scoring function that assesses the goodness of a given mapping. This is essentially a combinatorial optimization problem with a exponentially large search space, which makes finding the optimal mapping practically infeasible for large networks. As a result, existing network alignment methods employ various heuristic techniques to make the network alignment problem computationally tractable.

Several network alignment algorithms have been proposed so far [6–14], many of which focus on pairwise network alignment [16]. For example, GRAAL [9] analyzes the graphlet degree signature for two PPI networks, where it can generalize the degree of node by counting the number of graphlets for each node, and then align the two networks using a seed-and-extend approach. MI-GRAAL [10] extends GRAAL by integrating further sources of information (e.g., clustering coefficient or functional similarity) to measure the similarity between two networks. PINALOG [11] is another example of pairwise network alignment algorithm, which constructs the initial mapping for protein nodes that form dense subgraphs in the respective networks. This initial mapping is further extended by subsequently finding similar nodes in the neighborhood. Recently, a number of multiple network alignment algorithms have been proposed [12–14]. For example, SMETANA [12] tries to estimate probabilistic node correspondence scores using a semi-Markov random walk model, and then uses the estimated scores to predict the maximum expected accuracy (MEA) alignment of the given networks. Given a set of networks, NetCoffee [13] generates all possible combinations of bipartite graphs for these networks, and updates the edges in each bipartite graph based on the sequence similarity of the proteins and the topological structure of the networks. Then, the algorithm finds candidate edges (i.e., mappings) in the bipartite graphs and combines qualified edges through simulated annealing. BEAMS [14] is another recent multiple network alignment algorithm, which first extracts the so-called "backbones", or the minimal set of disjoint cliques in the filtered similarity graph, and then iteratively merges these backbones to maximize the overall alignment score.

In this paper, we propose a novel multiple network alignment algorithm based on a context-sensitive random walk (CSRW) model. The employed CSRW model adaptively switches between different modes of random walk in a context-sensitive manner by sensing and analyzing the present neighborhood of the random walker. This context-sensitive behavior improves the quantitative estimation of the potential correspondence between nodes belonging to different networks, ultimately, improving the overall accuracy of the multiple network alignment as we will demonstrate through extensive performance evaluation based on real and synthetic biological networks.

Methods

Maximum expected accuracy (MEA) alignment of biological networks

Let us assume that we have a set of N PPI networks G=\left\{{\mathcal{G}}_1,{\mathcal{G}}_2,\dots ,{\mathcal{G}}_N\right\}. Each network {\mathcal{G}}_n=\left({\mathcal{V}}_n,{\mathcal{E}}_n\right) has a set of nodes {\mathcal{V}}_n=\left\{v_1,v_2,\dots \right\} and edges {\mathcal{E}}_n=\left\{e_{i,j}\right\}, where e _i,j represents the interaction between nodes v _i and v _j in the network {\mathcal{G}}_n. For each pair of PPI networks {\mathcal{G}}_{\mathcal{U}}=\left(\mathcal{U},\mathcal{D}\right) and {\mathcal{G}}_{\mathcal{V}}=\left(\mathcal{V},\mathcal{E}\right), we denote the pairwise node similarity score for a node pair (u _i , v _j ), where u_i\in \mathcal{U} and v_j\in \mathcal{V}, as s(u _i , v _j ). In this study, we use the BLAST bit score between proteins as their node similarity score, but other types of similarity scores based on structural or functional similarity can be also utilized if available.

Suppose {\mathcal{A}}^{*} is the true alignment of the networks in the set G, which is unknown and needs to be predicted. As in [12, 17], we can define the accuracy of a given network alignment \mathcal{A} as follows

accuracy\left(\mathcal{A},{\mathcal{A}}^{*}\right)=\frac{1}{\left|\mathcal{A}\right|}{\displaystyle \sum _{u_i~v_j\in \mathcal{A}}}1\left(u_i~v_j\in {\mathcal{A}}^{*}\right),

(1)

where 1 (·) is an indicator function, whose value is 1 if the mapping u _i ~ v _j is included in the true alignment {\mathcal{A}}^{*} and 0 otherwise. The given measure assesses the goodness of the alignment \mathcal{A} based on the relative proportion of correctly aligned nodes. Of course, since the true alignment is not known, the accuracy of a network alignment \mathcal{A} cannot be measured using (1), hence we cannot directly use this measure to compare different potential alignments to choose the best one. A reasonable alternative would be to estimate the expected accuracy as follows

E_{{\mathcal{A}}^{*}}\left[accuracy\left(\mathcal{A},{\mathcal{A}}^{*}\right)\right]=\frac{1}{\left|\mathcal{A}\right|}{\displaystyle \sum _{u_i~v_j\in \mathcal{A}}}P\left(u_i~v_j|G\right),

(2)

where P (u _i ~ v _j |G) is the posterior alignment probability between the nodes u _i and v _j given the set of networks G. Based on this measure, our objective is then to predict the maximum expected accuracy (MEA) network alignment \widetilde{{\mathcal{A}}^{*}} of the networks in G as follows

\widetilde{{\mathcal{A}}^{*}}={\displaystyle \underset{\mathcal{A}}{\mathsf{\text{max}}}}{E}_{{\mathcal{A}}^{*}}\left[accuracy\left({\mathcal{A}}^{*},\mathcal{A}\right)\right].

(3)

A similar MEA approach [18] has been formerly adopted by a number of multiple sequence alignment algorithms, including ProbCons [17], ProbAlign [19], and PicXAA [20–22]. The MEA framework has been shown to be very effective in constructing accurate alignment of multiple biological sequences, making it one of the most popular approaches for sequence alignment. Recently, the MEA approach has been also applied to comparative network analysis, where RESQUE [4] performs MEA-based network querying and SMETANA [12] performs MEA-based multiple network alignment.

Comparing and aligning networks based on context-sensitive random walk

In order to find the alignment that maximizes the expected accuracy defined in (2), we first need an accurate method for estimating the posterior node alignment probability P (u _i ~ v _j | G). For this purpose, we adopt a context-sensitive random walk (CSRW) model, motivated by the pair hidden Markov model (pair-HMM) that has been widely used in sequence alignment [23]. The pair-HMM provides a simple, yet very effective, mathematical framework for estimating the alignment probability between symbols in different biological sequences. Unlike the traditional HMM, which generates a single symbol sequence, the pair-HMM generates a pair of aligned symbol sequences. Pair-HMM makes transitions between three different internal states M, I _X , and I _Y , where the M state emits an aligned pair of symbols, one symbol in sequence X and the other in sequence Y, while I _X and I _Y emit an unaligned symbol in sequence X and sequence Y, respectively. Given two biological sequences, the pair-HMM can be used to estimate the probability whether a given symbol pair was jointly emitted at state M, hence should be aligned to each other. This probability can be computed using the forward and backward algorithms and the resulting alignment probability provides us with a measure of confidence about the (biological) relevance between the given symbols (i.e., nucleotides, amino acids).

One of the most important features of pair-HMM is that it properly recognizes that conserved sequence patterns and motifs in different species may contain inserted and/or deleted symbols (often referred to as "indels") and therefore it specifically tries to model these indels. In a similar manner, a mathematical model that can recognize node insertions and deletions in different biological networks that contain conserved subnetwork regions and network motifs may be useful for obtaining a reliable posterior node-to-node alignment probability. Recently, random walk models have been shown to be effective for estimating the node correspondence in different networks [7, 12, 15] in a way that seam-lessly integrates both node similarity and topological similarity. However, the random walk models that were used in previous network alignment algorithms did not explicitly consider indels.

In this work, we adopt a novel context-sensitive random walk model that has been recently proposed to improve on existing models by taking such indels into account [24]. In a way that is conceptually similar to the pair-HMM, the CSRW has three different internal states M, I_{\mathcal{U}}, and I_{\mathcal{V}}, each of which corresponds to a different mode of random walk. At the M state, the random walker simultaneously moves on both networks to enter a pair of "matching" nodes. On the other hand, at the I_{\mathcal{U}} state, the random walker only moves on network {\mathcal{G}}_{\mathcal{U}} to enter a potentially "inserted" node in {\mathcal{G}}_{\mathcal{U}} that may not have a corresponding node in the network {\mathcal{G}}_{\mathcal{V}}. Similarly, at the I_{\mathcal{V}} state, the random walker only moves on {\mathcal{G}}_{\mathcal{V}} to enter a potentially inserted node in {\mathcal{G}}_{\mathcal{V}}. Transitions between states take place in a context-sensitive manner, where the random walker examines the neighboring nodes to determine the mode of random walk. For example, if there are node pairs with significant node similarity (i.e., potential orthologous nodes) in the immediate neighborhood, the CSRW switches to the M state to make a simultaneous move on both networks and randomly enter one of these node pairs. Otherwise, the CSRW switches to either I_{\mathcal{U}} or I_{\mathcal{V}} and performs an individual random walk only on one of the networks. Based on this random walk model, we compute the long-run proportion of time that a given pair of nodes will be simultaneously visited (i.e., at the M state), which can be used to compute a probabilistic correspondence score between these two nodes, as we will describe in the following section.

Estimation of node correspondence scores

Suppose we want to measure the correspondence between nodes that belong to two different networks {\mathcal{G}}_{\mathcal{U}}=\left(\mathcal{U},\mathcal{D}\right) and {\mathcal{G}}_{\mathcal{V}}=\left(\mathcal{V},\mathcal{E}\right), both of which are included in G, the set of PPI networks to be aligned. For every node pair (u _i , v _j ), where u_i\in \mathcal{U} and v_j\in \mathcal{V}, our goal is to quantify the level of confidence - which we refer to as the node correspondence score - using the CSRW model discussed earlier. For this purpose, we first construct the transition probability matrix that corresponds to the random walk. Let \mathcal{M} be the set of node pairs (u _i , v _j ) with a positive pairwise node similarity score s(u _i , v _j )

\mathcal{M}=\left\{\left(u_i,v_j\right)|s\left(u_i,v_j\right)>0,u_i\in \mathcal{U},v_j\in \mathcal{V}\right\}.

(4)

We also define the set of non-similar node pairs as follows

\mathcal{I}=\left\{\left(u_i,v_j\right)|s\left(u_i,v_j\right)=0,u_i\in \mathcal{U},v_j\in \mathcal{V}\right\}.

(5)

Let the current position of the random walker in the product graph be (u _c , v _c ), where u_c\in \mathcal{U} and v_c\in \mathcal{V}. In each time step, the random walker examines the set of similar neighboring nodes \mathcal{N}\left(u_c,v_c\right)=\left\{\left(u_i,v_j\right)|u_i\in \mathcal{N}\left(u_c\right),v_j\in \mathcal{N}\left(v_c\right),\left(u_i,v_j\right)\in \mathcal{M}\right\} to determine its mode of random walk (corresponding to one of the three possible internal states), where \mathcal{N}\left(u_c\right) is the set of neighbors of the node u _c in the network {\mathcal{G}}_{\mathcal{U}} and \mathcal{N}\left(v_c\right) is the set of neighbors of the node v _c in the network {\mathcal{G}}_{\mathcal{V}}. If there are similar node pairs among the neighboring node pairs, hence \mathcal{N}\left(u_c,v_c\right) is not empty, the random walker switches its internal state to the M state and performs a simultaneous walk on both networks, moving from (u _c , v _c ) to one of the nodes

\left(u_i,v_j\right)\in \mathcal{N}\left(u_c,v_c\right). We define the transition probability for this simultaneous walk as follows

P\left[\left(u_i,v_j\right)|\left(u_c,v_c\right)\right]=\frac{s\left(u_i,v_j\right)}{{\displaystyle \sum _{\left(u_{i^{\prime }},v_{i^{\prime }}\right)\in \mathcal{N}\left(u_c,v_c\right)}}s\left(u_{i^{\prime }},v_{i^{\prime }}\right)}.

(6)

In case there is no similar node pair around the current position of the random walker, that is \mathcal{N}\left(u_c,v_c\right)=\varnothing, the random walker randomly changes its state to either I_{\mathcal{U}} or I_{\mathcal{V}}, and performs an individual walk on the corresponding network {\mathcal{G}}_{\mathcal{U}} or {\mathcal{G}}_{\mathcal{V}}. The probability that a given network will be chosen for an individual random walk is proportional to its size (i.e., number of nodes in the network), which ensures that both networks are equally well-traversed at the I states. The random walker randomly moves to one of the neighboring nodes with equal probability on the selected network, while staying at the same node on the other network. Based on this behavior, the transition probabilities at state I_{\mathcal{U}} are given by

P\left[\left(u_i,v_c\right)|\left(u_c,v_c\right)\right]=\frac{\left|\mathcal{U}\right|}{\left|\mathcal{U}\right|+\left|\mathcal{V}\right|}\times \frac{1}{\left|\mathcal{N}\left(u_c\right)\right|}

(7a)

for u_i\in \mathcal{N}\left(v_c\right), and the transition probabilities at state I_{\mathcal{V}} are given by

P\left[\left(u_c,v_j\right)|\left(u_c,v_c\right)\right]=\frac{\left|\mathcal{V}\right|}{\left|\mathcal{U}\right|+\left|\mathcal{V}\right|}\times \frac{1}{\left|\mathcal{N}\left(v_c\right)\right|},

(7b)

for v_j\in \mathcal{N}\left(v_c\right).

Based on the transition probabilities given by (6), (7a), and (7b), we can construct the transition probability matrix P for the random walk on the two networks {\mathcal{G}}_{\mathcal{U}} and {\mathcal{G}}_{\mathcal{V}}. Given P, we can estimate the longrun proportion of time that the random walker spends in each pair of nodes (u _i , v _j ) by computing the steady state distribution π. In practice, since real PPI networks typically have a relatively small number of interactions (therefore only few edges for most nodes), the resulting transition probability matrix for the CSRW is sparse, which makes it relatively straightforward to compute the steady state distribution using the power method.

In order to increase the computational efficiency of the proposed network alignment method, instead of using the original transition probability matrix P, we use a reduced matrix \tilde{P}. The reduced matrix \tilde{P} is obtained by removing the rows and columns in P that correspond to node pairs in \mathcal{I} while keeping only the rows and columns that correspond to node pairs in \mathcal{M}. After the reduction, \tilde{P} is re-normalized to make it a legitimate stochastic matrix. In practice, since the CSRW is designed to spend more time at node pairs with higher similarity, the random walker spends a relatively small amount of time at node-pairs that belong to the set \mathcal{I}, and using the reduced matrix \tilde{P} instead of P only minimally affects the estimated long-run proportion of time spent at \left(u_i,v_j\right)\in \mathcal{M}. As a result, the difference in terms of network alignment performance that results from replacing the original matrix P by this reduced matrix \tilde{P} appears to be small as shown in the supplementary material (see Section S1).

We make one further modification to the CSRW in [24] by allowing the random walker to restart at a new position at each time step with a fixed restart probability λ. Note that a similar "random walk with restart" approach was used by IsoRank [6] and IsoRankN [7], although these algorithms do not utilize the CSRW adopted in our method. We allow the random walker to select its restart position according to the pairwise node similarity, such that node pairs with higher node similarity have higher chance to be the restart position of the random walker. To this aim, we normalize the pairwise node similarity scores so that they sum up to 1. Our final node correspondence score vector c is obtained from a linear combination of the steady-state distribution of the context-sensitive random walker \tilde{\pi } (estimated using the reduced transition probability matrix \tilde{P}) and the normalized node similarity score vector s as follows

c=\lambda s+\left(1-\lambda \right)\tilde{\pi }.

(8)

The above formulation, obtained by allowing the CSRW to restart the random walk at a new position, is especially useful when comparing real PPI networks, which are often incomplete and contain many isolated nodes. Simulation results show that the incorporation of the restart scheme can make our CSRW-based alignment method more robust, especially when the available topological data are either unreliable or insufficient for detecting the similarities between networks (see Section S2).

In order to determine the restart probability λ, we first analyze the structure of the reduced product graph of {\mathcal{G}}_{\mathcal{U}} and {\mathcal{G}}_{\mathcal{V}} that contains only similar node pairs included in \mathcal{M}. Intuitively, it is desirable to increase the restart probability λ if the networks are disconnected and decrease the probability if the networks are well connected. For example, if all the nodes in the reduced product graph are completely disconnected, it is desirable to restart the random walker at every step. Additionally, when we consider the following two cases - (i) most nodes in the product graph are connected and there are only a few disconnected nodes; (ii) the product graph is equally divided into N connected subnetworks of identical size - it would be desirable to assign a higher λ to the latter case. Based on these intuitions, we set the restart probability λ as the ratio of the total number of nodes in the top K% smallest subnetworks to the total number of nodes in the reduced product graph. In this work, we used K = 99% to determine the restart probability λ.

Constructing the multiple network alignment

Once we have computed the node correspondence scores in (8) for every pair of networks in G, we take a greedy approach as in [12] to construct the multiple network alignment. The overall alignment process is as follows. First, in order to improve the reliability of the node correspondence scores, we selectively apply the probabilistic consistent transformation (PCT) defined in [12]. If λ is larger than a predefined threshold λ _t , we do not apply PCT to the node correspondence scores. A large λ implies that the product graph is ill connected (e.g., containing a large number of isolated nodes), in which case applying the PCT would not be helpful and may in fact make the scores less reliable. This is because the PCT in [12] was developed based on the assumption that the product graphs for all network pairs are relatively well connected. After the potential score refinement step through PCT, we begin with an empty alignment and greedily add aligned node pairs (u _i , v _j ) to the network alignment, starting from the pairs with the highest node correspondence scores, until there is no other node pair left that can be added without creating inconsistencies in the network alignment. Assuming that the node correspondence scores in (8) obtained by the context-sensitive random walk model with restart accurately reflect the true correspondence between nodes - such that the score is proportional to the posterior node alignment probability - the proposed network alignment scheme can be viewed as a heuristic way to find the MEA alignment of the networks in G.

Results and discussion

Datasets and experimental set-up

To assess the performance of the proposed method, we tested the proposed network alignment method based on PPI networks in NAPAbench [25] and IsoBase [26]. NAPAbench is a network alignment benchmark that consists of 3 different datasets, referred to as the pairwise alignment dataset, 5-way alignment dataset, and 8-way alignment dataset. Each dataset contains three different subsets of 10 network families, each subset created using a different network growth model - CG (crystal growth), DMC (duplication-mutation-complementation), and DMR (duplication with random mutation). Each network family consists of 2, 5, or 8 PPI networks depending on the alignment dataset. For network families in the pairwise alignment dataset, each family contains one network with 3,000 nodes and the other with 4,000 nodes. In the 5-way network alignment dataset, a network family consists of 5 networks with 1,000, 1,500, 2,000, 2,500, and 2,500 nodes. Finally, in the 8-way alignment dataset, every network family consists of 8 networks, where each network contains 1,000 nodes. To evaluate the performance of the proposed method on real PPI networks, we utilized IsoBase datasets [26], which was constructed by integrating the following databases: BioGRID [27], DIP [28], HPRD [29], MINT [30], and IntAct [31]. IsoBase contains the PPI networks of five species: H. sapiens, M. musculus, D. melanogaster, C. elegans, and S. cerevisiae. Currently, the PPI network of H. sapiens in [26] has 22,369 proteins and 43,757 interactions, the PPI network of M. musculus has 24,855 proteins and 452 interactions, the PPI network of D. melanogaster has 14,098 proteins and 26,726 interactions, the PPI network of C. elegans has 19,756 proteins and 5,853 interactions, and the PPI network of S. cerevisiae has 6,659 proteins and 38,109 interactions. In our analysis, we excluded the M. musculus network as it currently contains only a small number of interactions.

Based on our simulations, we report the following performance metrics: correct nodes (CN), specificity (SPE), mean normalized entropy (MNE), conserved interaction (CI), coverage, and computation time. CN is the total number of nodes in the correct equivalence classes. Given a network alignment, an equivalence class is defined as the set of aligned nodes, and if all nodes in the equivalence class have the same functionality the given equivalence class is said to be correct. SPE is the relative number of correct equivalence classes to the total number of equivalence classes in a network alignment. For each equivalence class C, the normalized entropy can be computed by H\left(\text{C}\right)=-\frac{1}{\mathsf{\text{log}}d}{\sum }_{i=1}^d{p}_i\mathsf{\text{log}}{p}_i, where p _i is the relative proportion of nodes in C with functionality i and d is the total number of different functionalities in the given equivalence class. As a result, a network alignment that accurately maps functionally similar nodes, hence being functionally consistent, will have lower mean normalized entropy. CI is defined as the total number of edges between equivalence classes. We also count the total number of edges between correct equivalence classes, which we refer to as the conserved orthologous interactions (COI), to assess the biological relevance of the conserved interactions that have been identified by the network alignment method. Finally, for 5-way and 8-way alignment datasets, we measure the equivalence class coverage and the node coverage, where the former is the number of equivalence classes that include nodes from k different networks, and the latter is the number of nodes in an equivalence class whose equivalence class coverage is k. For the performance evaluation based on real PPI networks in IsoBase, we determined the functionality of each protein using the KEGG protein annotation [32, 33]. Note that nodes without any functional annotation in each equivalence class and equivalence classes that consist of a single node or nodes from a single network were removed before computing the performance metrics.

We compared the performance of the proposed multiple network alignment method against a number of state-of-the-art algorithms: SMETANA [12], PINALOG [11], BEAMS [14], NetCoffee [13], and IsoRankN [7]. NetCoffee was not included in pairwise network alignment experiments, since it requires at least 3 networks. For multiple network alignment experiments, PINALOG was excluded as the algorithm can only handle pairwise alignments. For IsoRankN, we set the parameter α to 0.6 as in the original paper [7]. For BEAMS, we set the filtering threshold to 0.4 for IsoBase and 0.2 for NAPAbench as in the original paper [14], and set the parameter α to 0.5. The parameter α for NetCoffee was set to 0.5. We used default parameters for SMETANA (i.e., n_max = 10, α = 0.9, and β = 0.8), and the same parameters were used in the proposed network alignment method as well. Finally, in the proposed method, we used λ _t = 0.7 to determine whether or not to apply PCT to the estimated node correspondence scores.

All experiments were performed on a personal computer with a 2.4 GHz Intel i7 processor and 8 GB memory.

Performance assessment based on NAPAbench network alignment benchmark

We first evaluated the performance of the proposed algorithm using the NAPAbench network alignment benchmark and compared it to other leading algorithms. The evaluation results are summarized in Table 1, 2, and 3, which show the average CN, SPE, and MNE of various network alignment algorithms.

Table 1 Performance comparison for pairwise network alignment.

Table 2 Performance comparison for 5-way network alignment.

Table 3 Performance comparison for 8-way network alignment.

As we can see in Table 1 in most cases, the proposed algorithm yields a significantly higher CN and SPE compared to other algorithms, which shows that the algorithm is capable of finding conserved nodes with both high sensitivity and specificity. Furthermore, the mean normalized entropy (MNE) is also much lower, indicating that the proposed algorithm yields network alignment results that are more functionally coherent. This table shows that BEAMS yields higher CN for the CG dataset, although its SPE is lower and its MNE is higher than the proposed method. Both SMETANA and the proposed algorithm shows similar performance on the CG dataset, but we can also see that the proposed algorithm consistently outperforms SMETANA on the DMC/DMR datasets. Multiple network alignment results obtained using the 5-way alignment dataset and the 8-way alignment dataset show similar trends. Tables 2 and 3 show that, in most cases, our proposed algorithm outperforms other algorithms with higher CN, higher SPE, and lower MNE. For multiple network alignment, we further compared different network alignment algorithms based on their capability of predicting equivalence classes that span all networks, since one of the main goals of multiple network alignment is to find functionally homologous proteins that are conserved in the networks of all target species. Simulation results show that, in most cases, our proposed method also yields much higher CN and SPE as well as lower MNE for equivalence classes that span all networks.

Next, we compare the number of conserved (orthologous) interactions identified by different network alignment algorithms. As Figure 1 shows, the proposed method was able to identify the largest number of conserved interactions as well as conserved orthologous interactions in most cases, resulting in higher CI and COI. The performance of SMETANA was comparable to the proposed method, while other algorithms typically resulted in lower CI and COI. It is worth noting that more than 95% of the conserved interactions that were detected by our proposed network alignment algorithm were between correct equivalence classes (i.e., conserved orthologous interactions). This certainly shows that our method can effectively detect biologically meaningful conserved interactions through network alignment.

Figure 1

The total number of conserved orthologous interactions (COI) and conserved interactions (CI): (a) pairwise network alignment; (b) 5-way network alignment; (c) 8-way network alignment.

We also analyzed the overall coverage of the predicted alignment results for the 5-way and 8-way network alignments. The results are shown in Figure 2 for the 5-way alignment and in Figure 3 for the 8-way alignment. For the 5-way network alignment, we can see that around 40% of the equivalence classes predicted by the proposed method contained nodes from all 5 networks. SMETANA shows a similar level of coverage, while for the remaining algorithms, only about 30% of the predicted equivalence classes included nodes from all 5 networks. The overall node coverage also shows similar trends. The 8-way alignment results summarized in Figure 3 show that the proposed algorithm can effectively find equivalence classes with good coverage, which include nodes from a large number of networks. For example, we can see that around 40% of the equivalence classes predicted by the proposed method contained nodes from all 8 networks.

Figure 2

Equivalence class coverage and node coverage for 5-way network alignment: (a) equivalence class coverage for each network growth model; (b) node coverage for each network growth model.

Figure 3

Equivalence class coverage and node coverage for 8-way network alignment: (a) equivalence class coverage for each network growth model; (b) node coverage for each network growth model.

Table 4 shows the mean computation time of the respective algorithms for aligning the network families in the NAPAbench datasets. As we can see in Table 4 SMETANA requires the least amount of time for aligning the networks in NAPAbench, while IsoRankN needs the most computation time. In our simulations, we observed that NetCoffee runs relatively fast, although its computation time varies significantly depending on the network structure. For example, it took much longer to align networks in the DMR dataset using NetCoffe, compared to networks in the DMC or CG datasets.

Table 4 Mean computation time for aligning PPI networks in the NAPAbench datasets (in seconds).

Performance assessment based on protein-protein interaction networks in IsoBase

For further evaluation, we performed additional experiments using real PPI networks in IsoBase. Table 5 shows the pairwise network alignment performance of the tested algorithms for several PPI network pairs. As we can see in this table, the proposed algorithm consistently performs fairly well in all cases, outperforming the other algorithms. We can make similar observations in Table 6 which summarizes the performance evaluation results for aligning 3 PPI networks. The proposed algorithm attains high CN, high SPE, and low MNE across all cases, showing that it can effectively compare and accurately align real PPI networks. BEAMS shows good performance on multiple alignment of real networks that is comparable to the proposed method, with a slightly lower SPE and a slightly higher MNE. Additionally, although BEAMS and IsoRankN achieve higher CN in some cases, the proposed method consistently yields higher CN than these methods with comparable SPE and MNE when we consider multiple network alignment results for regions that are conserved across all networks. Another observation we can make in Table 5 is that IsoRankN performs very well on real PPI networks compared to the other more recent algorithms. This is especially interesting, if we consider the fact that the performance of IsoRankN lagged behind the other algorithms according to the large-scale evaluations using NAPAbench. One possible explanation is that, for constructing the network alignment, IsoRankN relies on node similarity (i.e., sequence similarity in this case) more strongly compared to the other algorithms. In order to find out whether this is indeed a plausible explanation, we performed network alignment experiments solely using node similarity scores (i.e., without considering network topology), where we constructed the network alignment in a greedy manner by iteratively adding protein pairs with the highest node similarity scores. The alignment results are shown in Tables 5 and 6 right below the results for IsoRankN (labeled as "Node Similarity"). Surprisingly, these results show that this simple greedy network alignment approach that uses node similarity alone outperforms IsoRankN in most cases and surpasses all the other algorithms in all cases. In fact, currently available PPI networks are known to be very incomplete and these network typically contain a large number of isolated nodes. They are suspected to include a large number of spurious interactions while still missing many potential protein-protein interactions [34, 35]. Furthermore, only a small proportion of proteins in these PPI networks have reliable functional annotations (e.g., according to KEGG orthology), making it difficult to reliably assess the quality of a predicted network alignment. As a result, for current PPI networks, utilization of topological similarity between networks may not be necessarily helpful for improving the overall quality of the network alignment across the entire network. Moreover, since only a few large real PPI networks are available at the moment, we risk overtraining network alignment algorithms if they are mainly evaluated solely based on real PPI networks.

Table 5 Pairwise network alignment results for real PPI networks.

Table 6 Multiple network alignment results for real PPI networks (for 3 species).

Figure 4 shows the computation time for aligning the PPI networks in IsoBase. SMETANA required the least computation time for pairwise network alignment and NetCoffee was the fastest among all for aligning the PPI networks of 3 species. Although IsoRankN yielded accurate alignment results for real PPI networks in IsoBase, it also required the largest amount of computation time in most cases. Figure 4 shows that our proposed network alignment algorithm requires relatively longer running time compared to other algorithms, in exchange for the improved alignment accuracy. Currently, the main bottleneck is the time required to construct the transition probability matrix \tilde{P} of the context-sensitive random walker, and we are currently optimizing the code for our algorithm to make it computationally more efficient.

Figure 4

Computation time for aligning real PPI networks (in seconds): (a) pairwise network alignment; (b) multiple network alignment (for 3 species).

Conclusions

In this paper, we proposed a novel network alignment algorithm based on a context-sensitive random walk model that has been recently introduced. The CSRW provides an effective mathematical framework for comparing different biological networks and quantifying the node-to-node correspondence between nodes that belong to different networks. In our proposed method, we combined the CSRW model with a restart scheme, where the restart probability is automatically adjusted based on the characteristics of the networks under comparison. Furthermore, the proposed network alignment algorithm employs adaptive probabilistic consistency transformation, where the PCT is adaptively activated or deactivated based on the overall structure of the given networks. As we have shown through extensive performance evaluations based on biologically realistic PPI networks in NAPAbench as well as real PPI networks in IsoBase, the novel network alignment algorithm proposed in this paper can significantly improve the overall accuracy of pairwise as well as multiple network alignment.