Alignment of dynamic networks

TL;DR

This paper introduces DynaMAGNA++, the first dynamic network alignment method, which outperforms static methods by leveraging evolving network information in biological and social systems.

Contribution

It presents the first dynamic network alignment algorithm, extending MAGNA++ to incorporate temporal information through novel dynamic conservation measures.

Findings

Dynamic NA outperforms static NA in experiments.

DynaMAGNA++ effectively aligns evolving networks.

The method is applicable to biological and social networks.

Abstract

Networks can model real-world systems in a variety of domains. Network alignment (NA) aims to find a node mapping that conserves similar regions between compared networks. NA is applicable to many fields, including computational biology, where NA can guide the transfer of biological knowledge from well- to poorly-studied species across aligned network regions. Existing NA methods can only align static networks. However, most complex real-world systems evolve over time and should thus be modeled as dynamic networks. We hypothesize that aligning dynamic network representations of evolving systems will produce superior alignments compared to aligning the systems' static network representations, as is currently done. For this purpose, we introduce the first ever dynamic NA method, DynaMAGNA++. This proof-of-concept dynamic NA method is an extension of a state-of-the-art static NA method, MAGNA++. Even though both MAGNA++ and DynaMAGNA++ optimize edge as well as node conservation across the aligned networks, MAGNA++ conserves static edges and similarity between static node neighborhoods, while DynaMAGNA++ conserves dynamic edges (events) and similarity between evolving node neighborhoods. For this purpose, we introduce the first ever measure of dynamic edge conservation and rely on our recent measure of dynamic node conservation. Importantly, the two dynamic conservation measures can be optimized using any state-of-the-art NA method and not just MAGNA++. We confirm our hypothesis that dynamic NA is superior to static NA, under fair comparison conditions, on synthetic and real-world networks, in computational biology and social network domains. DynaMAGNA++ is parallelized and it includes a user-friendly graphical interface.

Tables

Table 1. Table 1 : Network discrimination performance of DynaMAGNA++, while optimizing both dynamic edge and node conservation (E+N), dynamic edge conservation alone (E), and dynamic node conservation alone (N), for biological synthetic networks, with respect to the area under the precision-recall curve (AUPR), F-score at which precision and recall cross and are thus equal (F-score cross cross {}_{\mbox{cross}} ), maximum F-score (F-score max max {}_{\mbox{max}} ), and the area under the ROC curve (AUROC).

NA method	AUPR	F-score${}_{\text{cross}}$	F-score${}_{\text{max}}$	AUROC
DynaMAGNA++ (E+N)	0.836	0.675	0.788	0.928
DynaMAGNA++ (E)	0.551	0.400	0.750	0.878
DynaMAGNA++ (N)	0.770	0.625	0.808	0.934

Equations

$$N_c=\sum _{(u,v)\in V_1\otimes V_1}\text{mathds}1\left[(u,v)\in E_1\wedge (f(u),f(v))\in E_2\right],$$

$$\begin{array}{rl}{N}_n& =\sum _{(u,v)\in V_1\otimes V_1}\text{mathds}1[((u,v)\in E_1\wedge (f(u),f(v))\notin E_2)\vee \\ & ((u,v)\notin E_1\wedge (f(u),f(v))\in E_2)]\\ & =|E_1|+|E_2^{\prime }|-2N_c,\end{array}$$

$$\text{mbox}{S}^3=\frac{N_c}{N_c+N_n}.$$

$$\begin{array}{r}\text{mbox}CET((u_1,u_2),(v_1,v_2))=\sum _{e\in T_{u_1{u}_2}}\sum _{e^{\prime }\in T_{v_1{v}_2}}ct(e,e^{\prime }),\end{array}$$

$$\text{mbox}NCET((u_1,u_2),(v_1,v_2))=\sum _{e\in T_{u_1{u}_2}}d(e)+\sum _{e^{\prime }\in T_{v_1{v}_2}}d(e^{\prime })-2\sum _{e\in T_{u_1{u}_2}}\sum _{e^{\prime }\in T_{v_1{v}_2}}ct(e,e^{\prime }),$$

$$T_c=\sum _{(u,v)\in V_1\otimes V_1}\text{mbox}CET((u,v),(f(u),f(v)))$$

$$T_n=\sum _{(u,v)\in V_1\otimes V_1}\text{mbox}NCET((u,v),(f(u),f(v))).$$

$$\text{mbox}{DS}^3=\frac{T_c}{T_c+T_n}.$$

Full text

**Results: **For this purpose, we introduce the first ever dynamic NA method, DynaMAGNA++. This proof-of-concept dynamic NA method is an extension of a state-of-the-art static NA method, MAGNA++. Even though both MAGNA++ and DynaMAGNA++ optimize edge as well as node conservation across the aligned networks, MAGNA++ conserves static edges and similarity between static node neighborhoods, while DynaMAGNA++ conserves dynamic edges (events) and similarity between evolving node neighborhoods. For this purpose, we introduce the first ever measure of dynamic edge conservation and rely on our recent measure of dynamic node conservation. Importantly, the two dynamic conservation measures can be optimized using any state-of-the-art NA method and not just MAGNA++. We confirm our hypothesis that dynamic NA is superior to static NA, under fair comparison conditions, on synthetic and real-world networks, in computational biology and social network domains. DynaMAGNA++ is parallelized and it includes a user-friendly graphical interface.

Software: Available upon request.

Supplementary information: Available upon request.

Contact: tmilenko@nd.edu, vvijayan@nd.edu

00footnotetext: ${}^{\text{1}}$Department of Computer Science and Engineering, ECK Institute for Global Health, and Interdisciplinary Center for Network Science and Applications (iCeNSA), University of Notre Dame, Notre Dame, IN 46556, USA00footnotetext: ${}^{\text{2}}$Department of Physics and Astronomy, Austin Peay State University, Clarksville, Tennessee, TN 37044, USA00footnotetext: ∗To whom correspondence should be addressed.

1 Introduction

Networks can be used to model complex real-world systems in a variety of domains (Boccaletti et al., 2006). Network alignment (NA) compares networks with the goal of finding a node mapping that conserves topologically or functionally similar regions between the networks. NA has been used in many domains and applications (Emmert-Streib et al., 2016). In computer vision, it has been used to find correspondences between sets of visual features (Duchenne et al., 2011). In online social networks, NA has been used to match identities of people who have different account types (e.g., Twitter and Facebook) (Zhang et al., 2015). In ontology matching, NA has been used to match concepts across ontological networks (Bayati et al., 2013). Computational biology is no exception. In this domain, NA has been used to predict protein function (including the role of proteins in aging), by aligning protein interaction networks (PINs) of different species, and by transferring functional knowledge from a well-studied species to a poorly-studied species between the species’ conserved (aligned) PIN regions (Faisal et al., 2015a, b; Elmsallati et al., 2016; Meng et al., 2016b; Guzzi and Milenković, 2017). Also, NA has been used to construct phylogenetic trees of species based on similarities of their PINs or metabolic networks (Kuchaiev et al., 2010; Kuchaiev and Pržulj, 2011).

NA methods can be categorized as local or global (Meng et al., 2016b; Guzzi and Milenković, 2017). Local NA typically finds highly conserved but consequently small regions among compared networks, and it results in a many-to-many node mapping. On the other hand, global NA typically finds a one-to-one node mapping between compared networks that results in large but consequently suboptimally conserved network regions. Clearly, each of local NA and global NA has its (dis)advantages (Meng et al., 2016b, a; Guzzi and Milenković, 2017). In this paper, we focus on global NA, but our ideas are applicable to local NA as well. Also, NA methods can be categorized as pairwise or multiple (Faisal et al., 2015b; Guzzi and Milenković, 2017; Vijayan and Milenković, 2016). Pairwise NA aligns two networks while multiple NA can align more than two networks at once. While multiple NA can capture conserved network regions between more networks than pairwise NA, which may lead to deeper biological insights compared to pairwise NA, multiple NA is computationally much harder than pairwise NA since the complexity of the NA problem typically increases exponentially with the number of networks. This is why in this paper we focus on pairwise NA, but our ideas can be extended to multiple NA as well. Henceforth, we refer to global and pairwise NA simply as NA.

Existing NA methods can only align static networks. This is because in many domains and applications, static network representations are often used to model complex real-world systems, independent of whether the systems are static or dynamic. However, most real-world systems are dynamic, as they evolve over time. Static networks cannot fully capture the temporal aspect of evolving systems. Instead, such systems can be better modeled as dynamic networks (Holme, 2015). For example, a complex system such as a social network evolves over time as friendships are made and lost. Static networks cannot model the changes in interactions between nodes over time, while dynamic networks can capture the times during which the friendships begin and end. Other examples of systems that can be more accurately represented as dynamic networks include communication systems, human or animal proximity interactions, ecological systems, and many systems in biology that evolve over time, including brain or cellular functioning. In particular, regarding the latter, while cellular functioning is dynamic, current computational methods (including all existing NA methods) for analyzing systems-level molecular networks, such as PINs, deal with the networks’ static representations. This is in part due to unavailability of experimental dynamic molecular network data, owing to limitations of biotechnologies for data collection. Yet, as more dynamic molecular (and other real-world) network data are becoming available, there is a growing need for computational methods that are capable of analyzing dynamic networks (Przytycka and Kim, 2010; Przytycka et al., 2010), including methods that can align such networks.

The question is: how to align dynamic networks, when the existing NA methods can only deal with static networks? To allow for this, we generalize the notion of static NA to its dynamic counterpart. Namely, we define dynamic NA as a process of comparing dynamic networks and finding similar regions between such networks, while exploiting the temporal information explicitly (unlike static NA, which ignores this information). We hypothesize that aligning dynamic network representations of evolving real-world systems will produce superior alignments compared to aligning the systems’ static network representations, as is currently done. To test this hypothesis, we introduce the first ever method for dynamic NA.

Our proposed dynamic NA method, DynaMAGNA++, is a proof-of-concept extension of a state-of-the-art static NA method, MAGNA++ (Vijayan et al., 2015). Saraph and Milenković (2014) and Vijayan et al. (2015) compared MAGNA++ to state-of-the-art static NA methods at the time, namely IsoRank (Singh et al., 2007), MI-GRAAL (Kuchaiev and Pržulj, 2011), and GHOST (Patro and Kingsford, 2012). More recently, Meng et al. (2016b) compared MAGNA++ to additional newer static NA methods: NETAL (Neyshabur et al., 2013), GEDEVO (Ibragimov et al., 2013), WAVE (Sun et al., 2015), and L-GRAAL (Malod-Dognin and Pržulj, 2015). The comparisons were made on synthetic as well as real-world PINs, in terms of both topological and functional alignment quality. MAGNA++ was found to be superior to six of the seven existing methods and comparable to the remaining method. This is exactly why we have chosen to extend MAGNA++ rather than some other static NA method to its dynamic counterpart. However, as any future static NA methods are developed (Hayes and Mamano, 2016) that are potentially superior to MAGNA++, our ideas on dynamic NA will be applicable to such methods too. Section 2 describes the method, and Section 3 confirms our hypothesis that dynamic NA is superior to static NA, under fair comparison conditions, on both synthetic and real-world networks, and on data from both computational biology and social network domains.

2 Methods

We first summarize MAGNA++, and then we describe our proposed dynamic NA method, DynaMAGNA++, as an extension of MAGNA++.

2.1 MAGNA++

Static networks and static NA. A static network $G(V,E)$ consists of a node set $V$ and an edge set $E$. An edge $(u,v)\in E$ is an interaction between nodes $u$ and $v$. There can only be a single edge between the same pair of nodes. Given two static networks $G_{1}(V_{1},E_{1})$ and $G_{2}(V_{2},E_{2})$, assuming without loss of generality that $|V_{1}|\leq|V_{2}|$, a static NA between $G_{1}$ and $G_{2}$ is a one-to-one node mapping $f\,\colon\,V_{1}\to V_{2}$, which produces the set of aligned node pairs $\{(v,f(v))\mid v\in V_{1}\}$ (Figure 1(a)).

Static edge conservation. Given an NA between two static networks, an edge in one network is conserved if it maps to an edge in the other network, and an edge in one network is non-conserved if it maps to a non-adjacent node pair (i.e., a non-edge) in the other network (Figure 1(a)). A good static NA is a node mapping that conserves similar network regions. That is, a good static NA should have a large number of conserved edges and a small number of non-conserved edges. In this context, we measure the quality of a static NA using the popular symmetric substructure score (S3) edge conservation measure (Saraph and Milenković, 2014).

S3 is defined as follows. Formally, the number of conserved edges is

[TABLE]

and the number of non-conserved edges is

[TABLE]

where $G^{\prime}_{2}(V^{\prime}_{2},E^{\prime}_{2})$ is the subgraph of $G_{2}$ induced by $V^{\prime}_{2}=\{f(u)\mid u\in V_{1}\}$, $\mathds{1}\left[p\right]=1$ if $p$ is true and $\mathds{1}\left[p\right]=0$ if $p$ is false, and $U\otimes V$ is the Cartesian product of sets $U$ and $V$. Then,

[TABLE]

Our implementation of S3 that can compute this measure in $O(|E_{1}|+|E_{2}|)$ time complexity is described in the Supplement.

Static node conservation. A good static NA should also conserve the similarity between aligned node pairs, i.e., node conservation. Node conservation accounts for similarities between all pairs of nodes across the two networks. Node similarity can be defined in a way that depends on one’s goal or domain knowledge. In this work, we use a node similarity measure that is based on graphlets, as follows.

Graphlets (in the static setting) are small, connected, induced subgraphs of a larger static network (Milenković and Pržulj, 2008). Graphlets can be used to describe the extended network neighborhood of a node in a static network via the node’s graphlet degree vector (GDV). The GDV generalizes the degree of the node, which counts how many edges are incident to the node, i.e., how many times the node touches an edge (where an edge is the only graphlet on two nodes), into the vector of graphlet degrees (i.e., GDV), which counts how many times the node touches each of the graphlets on up to $n$ nodes, accounting in the process for different topologically unique node symmetry groups (automorphism orbits) that might exist within the given graphlet. In this work, we use all graphlets with up to four nodes, which contain 15 automorphism orbits, when calculating the GDV of a node, per recommendations of the existing studies (Hulovatyy et al., 2015, 2014). Hence, the GDV of a node has 15 dimensions containing counts for the 15 orbits.

Given GDVs of all nodes in two static networks $G_{1}(V_{1},E_{1})$ and $G_{2}(V_{2},E_{2})$, where $x_{u}$ is the GDV of node $u$, we calculate similarity $s(u,v)$ between nodes $u\in V_{1}$ and $v\in V_{2}$ by relying on an existing GDV-based measure of node similarity that was used by Hulovatyy et al. (2015). The measure works as follows. First, to extract GDV dimensions that contain the most relevant information about the extended network neighborhood of the given node, the measure first reduces dimensionality of each GDV via principal component analysis (PCA). PCA is performed on the vector set $\{x_{w}\mid w\in V_{1}\cup V_{2}\}$, where as few as needed to account for at least 99% of variance in the vector set of the first $k$ PCA components are kept. Let us denote by $y_{u}$ the dimensionality-reduced vector of $x_{u}$ that contains the $k$ PCA components. Second, we define node similarity $s(u,v)$ as the cosine similarity between $y_{u}$ and $y_{v}$. Third, given a static NA $f$, we define our node conservation measure as $\sum_{u\in V_{1}}{\frac{s(u,f(u))}{\lvert V_{1}\rvert}}$.

Objective function and optimization process (also known as search strategy). MAGNA++ is a search-based algorithm that finds a static NA by directly maximizing both edge and node conservation. Namely, MAGNA++ maximizes the objective function $\alpha S_{E}+(1-\alpha)S_{N}$, where $S_{E}$ is the S3 measure of static edge conservation described above, $S_{N}$ is the graphlet-based measure of static node conservation described above, and $\alpha$ is a parameter between 0 and 1 that controls for the two measures. In several studies, it was shown that $\alpha$ of 0.5 yields the best results (Vijayan et al., 2015; Meng et al., 2016b), which is the $\alpha$ value we use in this study, unless otherwise noted. Given an initial population of random static NAs, MAGNA++ evolves the population of alignments over a number of generations while aiming to maximize its objective function. MAGNA++ then returns the alignment from the final generation that has the highest value of the objective function.

2.2 DynaMAGNA++

Dynamic networks. A dynamic network $H(V,T)$ consists of a node set $V$ and an event set $T$, where an event is a temporal edge (Figure 1(b)). An event is represented as a 4-tuple $(u,v,t_{s},t_{e})$, where nodes $u$ and $v$ interact from time $t_{s}$ to time $t_{e}$. An event $(u,v,t_{s},t_{e})\in T$ is active at time $t$ if $t_{s}\leq t\leq t_{e}$. The duration of an event is the time during which an event is active, i.e., $t_{e}-t_{s}$. There can be multiple events between the same two nodes in the dynamic network, but no two events between the same two nodes may be active at the same time. In fact, if there are two events between the same two nodes that are active at the same time, then they must be combined into a single event.

Above is the representation of a dynamic network that our study relies on. Sometimes, dynamic data is provided in a different dynamic network representation, most often as a discrete temporal sequence of static network snapshots $G_{1}(V_{1},E_{1}),\ldots,G_{k}(V_{k},E_{k})$. We can easily convert the static snapshot-based representation of a dynamic network into our event duration-based representation (i.e., into $H(V,T)$ as defined above). We do this as follows: if there is an edge connecting two nodes in the ${t}^{\textrm{th}}$ snapshot of the snapshot-based representation, then there is an event between the two nodes that is active from time $t$ to time $t+1$ in the event duration-based representation. In other words, we combine the node sets of the snapshots into a single node set $V=V_{1}\cup\ldots\cup V_{k}$. Then, for each snapshot $G_{t}$, $t=1,\ldots,k$, we convert each edge $(u,v)\in E_{t}$ into an event between nodes $u\in V$ and $v\in V$ in the dynamic network $H(V,T)$ with start time $t$ and end time $t+1$, i.e., the event $(u,v,t,t+1)$. This allows us to use the snapshot-based representation of a dynamic network in our study.

Dynamic NA. Given two dynamic networks $H_{1}(V_{1},T_{1})$ and $H_{2}(V_{2},T_{2})$, assuming without loss of generality that $|V_{1}|\leq|V_{2}|$, a dynamic NA between $H_{1}$ and $H_{2}$ is a one-to-one node mapping $f\,\colon\,V_{1}\to V_{2}$, which produces the set of aligned node pairs $\{(v,f(v))\mid v\in V_{1}\}$ (Figure 1(b)). Note the similarity between the definitions of static NA and dynamic NA (although the process of finding the actual alignments is different). This makes static NA and dynamic NA fairly comparable.

Dynamic edge (event) conservation. First, given node pair $(u_{1},u_{2})$ in $H_{1}$ that maps to node pair $(v_{1},v_{2})$ in $H_{2}$ (Figure 1(c)), we extend the notion of a conserved or non-conserved edge from static NA to dynamic NA by accounting for the amount of time that the mapping of $(u_{1},u_{2})$ to $(v_{1},v_{2})$ is conserved or non-conserved (defined below). That is, we extend the notion of a conserved or non-conserved static edge to the amount of a conserved or non-conserved dynamic edge (event), as follows.

Intuitively, we define the amount of a conserved event as follows. Similar to how an edge $(u^{\prime}_{1},u^{\prime}_{2})$ in static network $G_{1}(V_{1},E_{1})$ is conserved if it maps to an edge $(v^{\prime}_{1},v^{\prime}_{2})$ in static network $G_{2}(V_{2},E_{2})$ (and vice versa), the mapping of $(u_{1},u_{2})$ to $(v_{1},v_{2})$ is conserved at time $t$ if both $(u_{1},u_{2})$ and $(v_{1},v_{2})$ are active at time $t$. We refer to the entire amount of time during which this mapping is conserved as the conserved event time (CET) between $(u_{1},u_{2})$ and $(v_{1},v_{2})$. In other words, it is the amount of time during which both $(u_{1},u_{2})$ and $(v_{1},v_{2})$ are active at the same time. Formally, let $T_{u_{1}u_{2}}$ be the set of events between $u_{1}$ and $u_{2}$. Similarly, let $T_{v_{1}v_{2}}$ be the set of events between $v_{1}$ and $v_{2}$. Given this, the CET between $(u_{1},u_{2})$ and $(v_{1},v_{2})$ is

[TABLE]

where the conserved time $ct(e,e^{\prime})=\max(0,\min(t_{e},t^{\prime}_{e})-\max(t_{s},t^{\prime}_{s}))$ is the amount of time during which events $e=(u_{1},u_{2},t_{s},t_{e})$ and $e^{\prime}=(v_{1},v_{2},t^{\prime}_{s},t^{\prime}_{e})$ are active at the same time, i.e., $ct(e,e^{\prime})$ is the length of the overlap of the intervals $[t_{s},t_{e}]$ and $[t^{\prime}_{s},t^{\prime}_{e}]$.

Intuitively, we define the amount of a non-conserved event as follows. Similar to how an edge $(u^{\prime}_{1},u^{\prime}_{2})$ in $G_{1}$ is non-conserved if it maps to a disconnected node pair $(v^{\prime}_{1},v^{\prime}_{2})$ in $G_{2}$ (or vice versa), the mapping of $(u_{1},u_{2})$ to $(v_{1},v_{2})$ is non-conserved at time $t$ if exactly one of $(u_{1},u_{2})$ or $(v_{1},v_{2})$ is active at time $t$. We refer to the entire amount of time during which this mapping is non-conserved as the non-conserved event time (NCET) between $(u_{1},u_{2})$ and $(v_{1},v_{2})$. In other words, it is the amount of time during which $(u_{1},u_{2})$ is active, or $(v_{1},v_{2})$ is active, but not both are active at the same time. Formally, the NCET between $(u_{1},u_{2})$ and $(v_{1},v_{2})$ is

[TABLE]

where $d(e)$ is the duration of event $e$, i.e., the amount of time during which $e$ is active. We make sure to subtract twice the amount of time during which $(u_{1},u_{2})$ and $(v_{1},v_{2})$ are both active due to the above “but not both are active at the same time” constraint.

Second, given the above definitions of CET and NCET between two node pairs $(u_{1},u_{2})$ and $(v_{1},v_{2})$, we extend the S3 measure of static node conservation to a new dynamic S3 (DS3) measure of dynamic edge (event) conservation, which we propose as a contribution of this study. To define DS3, we need to introduce the notion of CET between all node pairs across the entire alignment (rather than between just two aligned node pairs), henceforth simply referred to as alignment CET, which is the sum of CET between all node pair mappings between $H_{1}$ and $H_{2}$. Analogously, we need to define the notion of alignment NCET, which is the sum of NCET between all node pair mappings between $H_{1}$ and $H_{2}$. Alignment CET measures the amount of event conservation of the entire alignment and alignment NCET measures the amount of event non-conservation of the entire alignment. A good dynamic NA is a node mapping that conserves similar evolving network regions. That is, a good dynamic NA should have high alignment CET and low alignment NCET, which is what DS3 aims to capture. Formally, alignment CET is

[TABLE]

and alignment NCET is

[TABLE]

Then,

[TABLE]

Our implementation of DS3 that can compute this measure in $O(|T_{1}|+|T_{2}|)$ time complexity is described in the Supplement.

We note that there are many real-world networks that contain events with durations that are significantly less than the entire time window of the network, called “bursty” events. Examples of networks containing bursty events are e-mail communication networks, economic networks that model transactions, and brain networks constructed from oxygen level correlations as measured by fMRI scanning, each of whose events last much less than a second while the networks’ time windows span minutes to hours (Holme, 2015). Since bursty events are so short, small perturbations in the event times can greatly affect the resulting dynamic edge (event) conservation value. Thus, in order to allow our DS3 measure to be more robust to perturbations in the event times, one may simply extend the duration of each event in the network by some time $\Delta t$. Extending the duration of each event by $\Delta t$ will account for perturbations in event times of up to $\frac{\Delta t}{2}$ due to the following. Given two events $(u_{1},u_{2},t,t)$ and $(v_{1},v_{2},t^{\prime},t^{\prime})$ with durations of 0, where $t^{\prime}=t+\frac{\Delta t}{2}$, the conserved time $ct(\cdot,\cdot)$ between the two events is 0. Thus, if we want to consider the two events as conserved, we can increase the durations of both events by $\Delta t$ to create the modified events $(u_{1},u_{2},t,t+\Delta t)$ and $(v_{1},v_{2},t^{\prime},t^{\prime}+\Delta t)$, which results in a conserved time of $\frac{\Delta t}{2}$ for the two modified events. While we do not use this technique in our work since we do not use networks with bursty events, others might in the future, and if so, this needs to be considered when performing dynamic NA.

Dynamic node conservation. Just as for static NA, a good dynamic NA method should also conserve the similarity between aligned node pairs, i.e., node conservation. To take advantage of the temporal information encoded in dynamic networks that are being aligned and also to make dynamic NA as fairly comparable as possible to static NA, in this work, we rely on a measure of node similarity based on dynamic graphlets, as follows.

Dynamic graphlets are an extension of static graphlets (Section 2.1) from the static setting to the dynamic setting by accounting for temporal information in the dynamic network. While static graphlets can be used to capture the static extended network neighborhood of a node, dynamic graphlets can be used to capture how the extended neighborhood of a node changes over time. To describe dynamic graphlets formally, we first present the notion of a $\Delta t$-time-respecting path and a $\Delta t$-connected network. A $\Delta t$-time-respecting path is a sequence of events that connect two nodes such that for any two consecutive events in the sequence, the end time of the earlier event and the start time of the later event are within $\Delta t$ time of each other (i.e., are $\Delta t$-adjacent). A dynamic network is $\Delta t$-connected if for each pair of nodes in the network, there is a $\Delta t$-time-respecting path between the two nodes. Then, just as a static graphlet is an equivalence class of isomorphic connected subgraphs (Section 2.1), a dynamic graphlet is an equivalence class of isomorphic $\Delta t$-connected dynamic subgraphs, where two graphlets are equivalent if they both have the same relative temporal order of events. We use $\Delta t=1$ in this work, per recommendations by Hulovatyy et al. (2015). Just as the GDV of a node in a static network is a topological descriptor for the extended neighborhood the node, there exists the dynamic GDV (DGDV) of a node in a dynamic network, which describes how the extended neighborhood of a node changes over time. Specifically, just as the GDV of a node counts how many times the node touches each static graphlet at each of its automorphism orbits, the DGDV of a node counts how many times the node touches each dynamic graphlet at each of its orbits. Dynamic graphlets have a similar notion of orbits as static graphlets do, which now depend on both topological and temporal positions of a node within the dynamic graphlet. To make things comparable as fairly as possible to static NA, and per recommendations by Hulovatyy et al. (2015), in this work, we use dynamic graphlets with up to four nodes and six events, which contain 3,727 automorphism orbits, when calculating the DGDV of a node. Hence, the DGDV of a node has 3,727 dimensions containing counts for the 3,727 orbits.

Given the DGDVs of all nodes in two dynamic networks $H_{1}(V_{1},T_{1})$ and $H_{2}(V_{2},T_{2})$, we calculate similarity $s(u,v)$ between nodes $u\in V_{1}$ and $v\in V_{2}$, in the same way as in Section 2.1 (by relying on the PCA-based dimensionality reduction of all nodes’ DGVDs, computing cosine similarity between the dimensionality-reduced PCA vectors, and accounting for resulting cosine similarities between all pairs of nodes across the compared networks to obtain the total dynamic node conservation).

Objective function and optimization process (also known as search strategy). DynaMAGNA++ is a search-based algorithm that finds a dynamic NA by directly maximizing both dynamic edge (event) and node conservation. Namely, DynaMAGNA++ maximizes the objective function $\alpha S_{T}+(1-\alpha)S_{N}$, where $S_{T}$ is the DS3 measure of dynamic edge conservation described above, $S_{N}$ is the DGDV-based measure of dynamic node conservation described above, and $\alpha$ is a parameter between 0 and 1 that controls for the two measures. To make DynaMAGNA++ fairly comparable to MAGNA++, here we also use MAGNA++’s best $\alpha$ value of 0.5, unless otherwise noted. Given an initial population of random dynamic NAs, DynaMAGNA++ evolves the population of alignments over a number of generations while aiming to maximize its objective function. DynaMAGNA++ then returns the alignment from the final generation that has the highest value of the objective function.

Time complexity. To align two dynamic networks $H_{1}(V_{1},T_{1})$ and $H_{2}(V_{2},T_{2})$, DynaMAGNA++ evolves a population of $p$ alignments over $N$ generations. It does so by using its crossover function (see Saraph and Milenković (2014) for details) to combine pairs of parent alignments in the given population into child alignments, for each generation. For each generation, the dynamic edge (event) conservation, dynamic node conservation, and crossover of $O(p)$ alignments are calculated. Since dynamic edge conservation takes $O(|T_{1}|+|T_{2}|)$ to compute, dynamic node conservation takes $O(|V_{1}|)$ time to compute, crossover takes $O(|V_{2}|)$ time to compute, and $\lvert V_{1}\rvert\leq\lvert V_{2}\rvert$, the time complexity of DynaMAGNA++ is $O(Np|V_{2}|+Np(|T_{1}|+|T_{2}|))$. Note that the calculation of dynamic edge and dynamic node conservation in DynaMAGNA++ is parallelized. This allows DynaMAGNA++ to be run on multiple cores, which empirically results in close to liner speedup.

Other parameters. Given an initial population of dynamic NAs, DynaMAGNA++ evolves the population for up to a specified number of generations or until it reaches a stopping criterion. For each generation, DynaMAGNA++ keeps an elite fraction of alignments from the current generation’s population for the next generation’s population. In addition to the dynamic edge and node conservation measures, and the $\alpha$ parameter that controls for the contribution of the two measures, the remaining parameters of DynaMAGNA++ are (i) the initial population, (ii) the size of the population, (iii), the maximum number of generations, (iv) the elite fraction, and (v) the stopping criterion. For DynaMAGNA++, we use a population of 15,000 alignments initialized randomly, as in the original MAGNA++ paper, unless otherwise noted. We specify a maximum of 10,000 generations, since the alignments that we test all converge by 10,000 generations. The elite fraction is 0.5, as in the original MAGNA++ paper. The algorithm stops when the highest objective function value in the population has increased less than 0.0001 within the last 500 generations, since the alignments that we test do not increase by a significant amount after this point.

To fairly compare DynaMAGNA++ against MAGNA++, we aim to set the parameters of both methods to be as similar as possible. So, other than MAGNA++’s edge and node conservation measures, the remaining parameters of MAGNA++ are the same as for DynaMAGNA++. This way, any differences that we see between results of DynaMAGNA++ and results of MAGNA++ will be the consequence of the differences of the two methods’ edge and node conservation measures, i.e., of accounting for temporal information in the network with DynaMAGNA++ and ignoring this information with MAGNA++. In other words, any differences that we see between results of DynaMAGNA++ and results of MAGNA++ will fairly reflect differences between dynamic NA and static NA.

3 Results and discussion

Since there are no other dynamic NA methods to compare against, we compare DynaMAGNA++ to the next best option, namely its static NA counterpart. That is, we compare DynaMAGNA++ when it is used to align two dynamic networks, to MAGNA++ when it is used to align static versions of the two dynamic networks. By “static versions”, we mean that we “flatten” or “aggregate” a dynamic network into a static network that will have the same set of nodes as the dynamic network and a static edge will exist between two nodes in the static network if there is at least one event between the same two nodes in the dynamic network. This network aggregation simulates the common practice where network analysis of time-evolving systems is done in a static manner, by ignoring their temporal information (Hulovatyy et al., 2015; Holme, 2015).

We evaluate DynaMAGNA++ and MAGNA++ on synthetic and real-world dynamic networks, as described in the following sections. Note that there is no need to compare DynaMAGNA++ to any other static NA method besides MAGNA++, because MAGNA++ was recently shown in a systematic and comprehensive manner to be superior to seven other state-of-the-art static NA methods (Section 1). So, by transitivity, to demonstrate that dynamic NA is superior to static NA, it is sufficient to demonstrate that DynaMAGNA++ is superior to MAGNA++.

3.1 Evaluation using synthetic networks

Motivation. A good NA approach should be able to produce high-quality alignments between networks that are similar and low-quality alignments between networks that are dissimilar (Yaveroğlu et al., 2015). In this test on synthetic networks, “similar” means networks that originate from the same network model, and “dissimilar” means networks that originate from different network models. So, we refer to this test as network discrimination. Thus, in this section, we evaluate the network discrimination performance of DynaMAGNA++ and MAGNA++.

Data. We perform this evaluation on a set of biologically inspired synthetic networks. Specifically, we generate 20 dynamic networks using four biologically inspired network evolution models (or versions of the same model with different parameter values) that simulate the evolution of PINs, resulting in five networks per model (Hulovatyy et al., 2015). The four models we use are (i) GEO-GD with $p=0.3$, (ii) GEO-GD with $p=0.7$, (iii) SF-GD with $p=0.3$ and $q=0.7$, and (iv) SF-GD with $p=0.7$ and $q=0.6$, where GEO-GD is a geometric gene duplication model with probability cut-off and SF-GD is a scale-free gene duplication model (Pržulj et al., 2010). Hulovatyy et al. (2015) generalized the static versions of these models to their dynamic counterparts, and we rely on the same model networks as those used by Hulovatyy et al. (2015) (see their paper for details). The synthetic networks are represented as snapshots, with 1,000 nodes in each of the networks, an average of 24 snapshots per network, and an average of 162 edges per snapshot, where any variation in network sizes is caused by the different parameter values of the considered network models.