Consensus ranking for multi-objective interventions in multiplex networks

TL;DR

This paper introduces a consensus ranking method using rank aggregation to identify influential nodes in heterogeneous multiplex networks, balancing multiple objectives for targeted interventions.

Contribution

It develops a novel approach combining rank aggregation with Pareto optimization to handle multi-objective interventions in multiplex networks with diverse layer functions.

Findings

Consensus ranking effectively balances competing layer objectives.

Weighted aggregation reveals trade-offs in multi-objective interventions.

Model networks demonstrate how multiplex structure influences intervention strategies.

Abstract

High-centrality nodes have disproportionate influence on the behavior of a network; therefore controlling such nodes can efficiently steer the system to a desired state. Existing multiplex centrality measures typically rank nodes assuming the layers are qualitatively similar. Many real systems, however, are comprised of networks heterogeneous in nature, for example, social networks may have both agnostic and affiliative layers. Here, we use rank aggregation methods to identify intervention targets in multiplex networks when the structure, the dynamics, and our intervention goals are qualitatively different for each layer. Our approach is to rank the nodes separately in each layer considering their different function and desired outcome, and then we use Borda count or Kemeny aggregation to identify a consensus ranking - top nodes in the consensus ranking are expected to effectively balance the competing goals simultaneously among all layers. To demonstrate the effectiveness of consensus ranking, we apply our method to a degree-based node removal procedure such that we aim to destroy the largest component in some layers, while maintaining large-scale connectivity in others. For any multi-objective intervention, optimal targets only exist in the Pareto-sense; we, therefore, use a weighted generalization of consensus ranking to investigate the trade-off between the competing objectives. We use a collection of model and real networks to systematically investigate how this trade-off is affected by multiplex network structure. We use the copula representation of the multiplex centrality distributions to generate model multiplex networks with given rank correlations. This allows use to separately manipulate the marginal centrality distribution of each layer and the interdependence between the layers and to independently investigate the role of the two using both analytical and numerical methods.

Tables

Table 1. Table 1: Properties of the example multiplex networks. We provide the number of nodes ( N 𝑁 N ) and layers ( L 𝐿 L ), the average degree of each layer ( ⟨ k ⟩ delimited-⟨⟩ 𝑘 \left\langle k\right\rangle ); and we quantify the degree heterogeneity of each layer with κ = ⟨ k 2 ⟩ / ⟨ k ⟩ 𝜅 delimited-⟨⟩ superscript 𝑘 2 delimited-⟨⟩ 𝑘 \kappa=\left\langle k^{2}\right\rangle/\left\langle k\right\rangle (for Poisson degree distribution κ = 1 + ⟨ k ⟩ 𝜅 1 delimited-⟨⟩ 𝑘 \kappa=1+\left\langle k\right\rangle and for scale-free networks κ → ∞ → 𝜅 \kappa\rightarrow\infty as N → ∞ → 𝑁 N\rightarrow\infty ) and with the maximum degree k max subscript 𝑘 max k_{\textrm{max}} .

	$N$	$L$	$⟨k⟩$	$\kappa =⟨k^2⟩/⟨k⟩$	$k_{\text{max}}$
Airline	575	5	4.63, 2.95, 9.02, 5.99, 4.11	28.7, 39.1, 37.2, 31.9, 28.3	109, 126, 125, 95, 66
Primate Soc.	100	4	8.62, 4.90, 4.32, 4.98	12.1, 8.9, 5.8, 6.3	35, 29, 13, 12
Human Soc.	61	5	6.33, 4.07, 0.69, 2.89, 6.36	7.8, 9.5, 2.3, 5.7, 10.8	15, 15, 5, 14, 27

Equations

$$b_{\alpha ,i}=|\{j\in \{1,2,\dots ,N\}:c_{\alpha ,i}>c_{\alpha ,j}\}|$$

$$b_i=\frac{1}{L}\sum _{\alpha =1}^L{b}_{\alpha ,i},$$

$$R={\mathrm{argmin}}_s\frac{1}{L}\sum _{\alpha =1}^{L}D(r_{\alpha },s),$$

$$D(r,s)=\frac{1}{N(N-1)}\sum _{i<j}|r(i,j)-s(i,j)|,$$

$$r(i,j)=\{\begin{array}{ll}1,& \text{if }{r}_i<r_j,\\ 0,& \text{if }{r}_i=r_j,\\ -1,& \text{if }{r}_i>r_j,\end{array}$$

$$D_{\text{con}}^{\text{BC/KA}}=\frac{1}{L}\sum _{\alpha =1}^{L}D(r_{\alpha },R_{\text{BC/KA}}),$$

$$\begin{array}{rl}{r}_{\alpha ,i}^{\text{D}}& =N-|\{j\in \{1,\dots ,N\}:k_{\alpha ,i}>k_{\alpha ,j}\}|,\\ r_{\alpha ,i}^{\text{K}}& =N-|\{j\in \{1,\dots ,N\}:k_{\alpha ,i}<k_{\alpha ,j}\}|,\end{array}$$

$$b_i=\frac{1}{L}\sum _{\alpha =1}^L{w}_{\alpha }{b}_{\alpha ,i},$$

$$1-s=\sum _{k=1}^{\infty }\frac{k}{c}\left(1-f(k)\right)p(k)(1-s{)}^{k-1}+\sum _{k=1}^{\infty }\frac{k}{c}f(k)p(k),$$

$$1-S=\sum _{k=0}^{\infty }\left(1-f(k)\right)p(k)(1-s{)}^k+\sum _{k=0}^{\infty }f(k)p(k).$$

$$f(k)=\{\begin{array}{ll}1,& \text{if }k>K,\\ \frac{f-{\sum }_{k>K}p(k)}{p(K)},& \text{if }k=K,\\ 0,& \text{if }k<K,\end{array}$$

$$f(k)=\{\begin{array}{ll}0,& \text{if }k>K,\\ \frac{f-{\sum }_{k<K}p(k)}{p(K)},& \text{if }k=K,\\ 1,& \text{if }k<K,\end{array}$$

$$b(k_1,\dots ,k_L)=\frac{w}{L}\sum _{\alpha \in {\mathcal{L}}_{\text{D}}}{P}_{\alpha }(k_{\alpha }-1)+\frac{1-w}{L}\sum _{\alpha \in {\mathcal{L}}_{\text{K}}}\left(1-P_{\alpha }(k_{\alpha })\right),$$

$$p_{\text{B}}(b)=\sum _{k_1,\dots ,k_L=0}^{\infty }p(k_1,\dots ,k_L)\delta \left(b-b(k_1,\dots ,k_L)\right),$$

$$f(b)=\{\begin{array}{ll}1,& \text{if }b>B,\\ \frac{f-{\sum }_{b>B}{p}_{\text{B}}(b)}{p_{\text{B}}(B)},& \text{if }b=B,\\ 0,& \text{if }b<B,\end{array}$$

$$f_{\alpha }(k_{\alpha })=\frac{{\sum }_{\{k_{\beta }:\beta \ne \alpha \}}p(k_1,\dots ,k_L)f(b(k_1,\dots ,k_L))}{{\sum }_{\{k_{\beta }:\beta \ne \alpha \}}p(k_1,\dots ,k_L)},$$

$$\begin{array}{rl}{O}_{\text{D}}& =\frac{1}{|{\mathcal{L}}_{\text{D}}|}\sum _{\alpha \in {\mathcal{L}}_{\text{D}}}\frac{S_{\alpha }-S_{\alpha }^{\text{min}}}{S_{\alpha }^{\text{max}}-S_{\alpha }^{\text{min}}},\\ O_{\text{K}}& =\frac{1}{|{\mathcal{L}}_{\text{K}}|}\sum _{\alpha \in {\mathcal{L}}_{\text{K}}}\frac{S_{\alpha }^{\text{max}}-S_{\alpha }}{S_{\alpha }^{\text{max}}-S_{\alpha }^{\text{min}}},\end{array}$$

$$P(c_1,\dots ,c_L)=C(P_1(c_1),\dots ,P_L(c_L)),$$

$$D_{\alpha \beta }=1-2\underset{0}{\overset{1}{\int }}d{u}_{\alpha }\underset{0}{\overset{1}{\int }}d{u}_{\beta }C(u_{\alpha },u_{\beta })\frac{{\partial }^{2}C(u_{\alpha },u_{\beta })}{\partial u_{\alpha }\partial u_{\beta }},$$

$$s_{\alpha ,i}=N-|\{j\in \{1,\dots ,N\}:u_{\alpha ,i}>u_{\alpha ,j}\}|.$$

$$C_{\rho }(u_1,\dots ,u_L)={\Phi }_{\rho }({\Phi }^{-1}(u_1),\dots ,{\Phi }^{-1}(u_L)),$$

$${\rho }_{\alpha \beta }=\sin \left(\frac{\pi }{4}(1-D_{\alpha \beta })\right),$$

$$p(k)=\frac{1}{N}\sum _{i=1}^N\frac{c_i^k}{k!}{e}^{-c_i},$$

$$p(k)=\frac{[c(1-\alpha ){]}^{1/\alpha }}{\alpha }\frac{\Gamma (k-1/\alpha ,c[1-\alpha ])}{\Gamma (k+1)},$$

Full text

I Introduction

In complex networks, a small subset of nodes often has disproportionate influence on the behavior of the system Freeman (1978); Page et al. (1999); Albert et al. (2000); Newman (2010), and controlling such nodes allows us to steer the network to desired states Motter and Lai (2002); Liu et al. (2011); Fiedler et al. (2013). For example, vaccinating a small fraction of carefully selected nodes can suppress large-scale disease outbreaks Pastor-Satorras and Vespignani (2002), or in a social network information is often disseminated by a small set of influencers Kitsak et al. (2010). Network centrality measures – such as degree, PageRank, or betweenness centrality – rank nodes based on their importance with respect to some process of interest. Therefore, high-centrality nodes provide effective targets to influence the behavior of a system: a node with high eigenvector centrality is an ideal target for vaccination Masuda (2009), or a node in a social network with high closeness centrality is an effective influencer Freeman (1978). Nodes, however, typically participate in multiple networks simultaneously forming a multilayer or multiplex network Kivelä et al. (2014); Boccaletti et al. (2014); Bianconi (2018). For example, a person typically participates in a number of friendship and professional social networks, and these layers are often qualitatively different in nature, including affiliative and competitive interactions Harary et al. (1953); Wasserman and Faust (1994). Perturbing a node to influence one layer inadvertently affects the other layers that the node is active in, which can have drastic unintended consequences Motter and Lai (2002); Sornette and Ouillon (2012); Yu et al. (2016); Lin et al. (2018). Therefore, it is desirable to identify targets in a multiplex network that maximize our influence in certain layers, yet minimize any unwanted impact on others. More generally, we consider the scenario where the structure, the dynamics, and our intervention goals are qualitatively different for each layer of a multiplex network, and we aim to find nodes that effectively balance the trade-offs between the multiple objectives.

Recently, a number of centrality measures were introduced to rank nodes in multilayer or multiplex networks Halu et al. (2013); Estrada and Gómez-Gardenes (2014); Solé-Ribalta et al. (2016); Rahmede et al. (2017). These centralities typically assume qualitatively similar layers with similar dynamics; and therefore they may provide useful intervention targets if the structure, the dynamics, and our intervention goals are all similar for each layer. For example, the family of “versatility centralities” extend common single-layer centralities to multilayer networks; as such, PageRank versatility is useful to identify influential scientists in a multiplex co-authorship network Perc (2010); Menichetti et al. (2014); De Domenico et al. (2015a); or airports with high betweenness versatility are prone to congestion in the multiplex air traffic network De Domenico et al. (2015b). Such multiplex centralities, however, do not address the more complex, and therefore less studied, scenario where the dynamics on the layers or our objectives are qualitatively different. For example, social and political multiplex networks are often comprised of cooperative (friendship, collaboration, alliance) and competitive (hostility, fighting, distrust) layers Harary et al. (1953); Wasserman and Faust (1994); Maoz (2010); possible interventions might aim to reduce hostility while maintaining cooperation. In fact, this is a practical issue that caretakers of captive nonhuman primate groups face: sometimes unusually aggressive animals are removed with the goal of reducing overall levels of conflict, and it is desirable that this removal does not negatively affect the cohesion of the group Judge et al. (1994); McCowan et al. (2018); Pósfai and D’Souza (2018). For another example, consider the recent work modeling brain activity as a two-layer multiplex network with synchronization and transport dynamics. A possible strategy to affect the behavior of the system is to target nodes that are influential in both layers Nicosia et al. (2017). Or another recent work studied a two-layer SIS model where one layer represents the spread of a disease and the other the spread of awareness of that disease Granell et al. (2013). Although both layers are characterized by SIS dynamics, a possible intervention would have different goals for the two layers: blocking the propagation of the disease and promoting awareness; and the best spreaders and the best blockers are known to have different properties Radicchi and Castellano (2017). In these example multiplex systems, the function of each layer and/or our corresponding intervention goals are qualitatively different and ranking nodes in such networks has so far remained an uninvestigated problem.

In this paper, we explore the use of rank aggregation methods to identify target nodes for multi-objective interventions in multiplex networks. Instead of identifying influential nodes based on a single integrated multiplex centrality measure, we rank the nodes in each layer separately based on a centrality that is relevant to the function of the particular layer and our objectives. We then use rank aggregation methods to find a consensus ranking; high-ranking nodes are expected to balance all the objectives simultaneously. An expanding number of rank aggregation methods exist, and although these were originally studied by economists in the context of social choice theory Arrow (2012), many also found applications in other disciplines, from computer science to biology Dwork et al. (2001); Lin (2010). In this paper, we focus on two widely-used methods, Borda count and Kemeny aggregation Kemeny and Snell (1962); de Borda (1781). In Sec. II, we introduce consensus ranking in multiplex networks, and we demonstrate several of its properties using model and real networks. As an example application, in Sec. III, we study degree-based node removal of multiplex networks, such that we aim to destroy the largest connected component of some layers, while maintaining connectivity in the rest. In case of multiple objectives, optimality only exists in the Pareto-sense; therefore we rely on a generalized version of the Borda count algorithm that allows us to assign varying weight to the different objectives. Using the weighted Borda count, we systematically study the trade-off between the competing goals, using analytically solvable model networks and a set of real networks.

II Consensus ranking

A multiplex network consists of $L$ layers, where each layer is a network $\mathcal{G}_{\alpha}$ ($\alpha=1,2,\ldots,L$) with the same set of nodes $\{v_{1},v_{2},\ldots,v_{N}\}$. These layers represent distinct types of interactions or relationships between the nodes, and there is no restriction on the structure of the individual layers, e.g., some layers might be directed, while others undirected, or they can be weighted or unweighted. To quantify the importance of a node $i$ in an individual layer $\alpha$, we calculate the node’s centrality $c_{\alpha,i}$ as if layer $\alpha$ is in isolation. Nodes with high centrality are important for the functioning of that given layer, and targeting such nodes is an effective strategy to influence or monitor the behavior of that particular layer.The type of centrality measure that works best depends on the function of the layer, for example, betweenness centrality is used for layers representing a transportation network, or eigenvector centrality is useful for layers governed by diffusion-like dynamics. However, nodes in a multiplex network do not perform a single task, but they simultaneously participate in all layers, prompting the question: How do we identify nodes that are important in all or most layers? This is particularly challenging if layers perform qualitatively different functions or our criteria for ranking nodes is different for each layer. To overcome these difficulties, instead of directly combining the centralities of each layer, we first determine the node ranks in each layer and then we use rank aggregation methods to identify a consensus ranking. An expanding number of rank aggregation methods exist Arrow (2012); Dwork et al. (2001); Lin (2010), here we focus on two widely-used and intuitive methods: the Borda count (BC) and the Kemeny aggregation (KA) de Borda (1781); Kemeny and Snell (1962).

We define the rank of node $i$ in layer $\alpha$ as $r_{\alpha,i}=N-b_{\alpha,i}$, where

[TABLE]

is the number of nodes that have smaller centrality than node $i$. We denote the list of node ranks in layer $\alpha$ as $r_{\alpha}=(r_{\alpha,1},r_{\alpha,2},\ldots,r_{\alpha,N})$. If there are no ties in the rankings (i.e., $c_{\alpha,i}\neq c_{\alpha,j}$ if $i\neq j$) $r_{\alpha}$ provides an ordinal ranking, the top rank being $1$. If ties exist, the ranks of tied nodes are set to be equal and assigned the worst value, e.g., if nodes $i$ and $j$ have equal centralities $c_{\alpha,i}=c_{\alpha,j}$ and are top ranked, we assign $r_{\alpha,i}=r_{\alpha,j}=2$ to both nodes.

The BC algorithm is perhaps the most straight-forward rank aggregation method de Borda (1781), it works by assigning a score $b_{i}$ to each node $i$ equal to the average number of nodes that $i$ outranks in the layers of the multiplex:

[TABLE]

where $b_{\alpha,i}$ is defined in Eq. (1). The BC ranking, $R_{\textrm{BC}}$, is simply a ranking of the nodes based on the scores $b_{i}$. The algorithm requires sorting the nodes; therefore its computational complexity is $O(N\log N)$.

The KA algorithm belongs to a larger class of rank aggregation methods that aim to identify a consensus ranking $R$ such that $R$ is at the minimum average distance from the rankings $r_{\alpha}$ Irurozki (2014); Ceberio et al. (2015). Formally, it is defined as

[TABLE]

where $D(r,s)$ is some measure of distance between two rankings $r$ and $s$, such as Hamming distance, Cayley distance, etc. The KA algorithm uses the Kemeny distance Kemeny and Snell (1962), which naturally takes into account ties:

[TABLE]

where

[TABLE]

and $s(i,j)$ is defined similarly. Kemeny distance is normalized such that if $r$ and $s$ are the same, $D(r,s)=0$; and if $r$ and $s$ rank nodes in opposite order, $D(r,s)=1$. If there are no ties in the rankings, the Kemeny distance is equivalent to Kendall’s distance 111Note that if ties are present Kemeny distance is not a distance in the mathematical sense, as it does not satisfy the triangle inequality..

Unfortunately, identifying the Kemeny consensus is an NP-complete problem for $L\geq 4$ Bartholdi et al. (1989); therefore approximate methods have to be used. Over the years, a high number of such algorithms were proposed: a recent survey compared 104 algorithms and combinations of algorithms Ali and Meilă (2012). In fact, the previously introduced BC algorithm can be considered as a simple heuristic approximation of the KA problem. Here, we implemented a local search algorithm, which was identified as providing optimal trade-off between accuracy and run-time Ali and Meilă (2012). The algorithm starts from $R_{\textrm{BC}}$ and outputs a new consensus ranking $R_{\textrm{KA}}$ by finding a local minimum of Eq. (3) using a restricted set of transformations. The computational complexity of our implementation is $O(N^{2})$; we provide the details of the algorithm in Appendix A.

II.1 Model networks

We now compare the two algorithms by applying them to a multiplex network model with tunable pairwise Kemeny distance between its layers. To generate these multiplex networks, we first independently create $L$ layers, and we determine the rank $r_{\alpha,i}$ of each node $i$ in each layer $\alpha$ with respect to some centrality. We then generate a node label sequence for each layer to ensure that the distance between each pair of layers is $D$. Here, we use the Erdős-Rényi model (ER) or the scale-free static model (SF) to generate each layer, and we use degree centrality to rank the nodes. We describe the procedure in detail in Appendix B.

By choosing parameter $D$ of the model network, we set the strength of consensus: if $D=0.5$, the layers are independent and there is no intrinsic consensus; if $D=0$, the ranking in each layer is the same and consensus is perfect. We compare the two algorithms by calculating the average distance of the layers from consensus:

[TABLE]

where $R_{\textrm{BC/KA}}$ is the consensus ranking found by the $BC$ or $KA$ algorithm. Note that $D_{\textrm{con}}$ corresponds to the cost function of the KA algorithm provided by Eq. (3), and we chose $D_{\textrm{con}}$ for comparison because BC is sometimes considered as an approximation of KA (Ali and Meilă, 2012). The results, however, have to be interpreted with care, finding a ranking that corresponds to lower $D_{\textrm{con}}$ doesn’t necessarily mean that it is better, because ultimately BC and KA define the consensus in different ways. Similarly to the fact that no single best centrality measure exists to rank the nodes, the preferred rank aggregation algorithm also depends on our particular purposes.

Figure 1 shows $D_{\textrm{con}}$ as a function of the distance between layers $D$. Since our implementation of the KA algorithm works by improving the solution of the BC algorithm, we always find that $D^{\textrm{KA}}_{\textrm{con}}\leq D^{\textrm{BC}}_{\textrm{con}}$ (the equality holds for $L=2$, where both algorithms find the same consensus, corresponding to the global minimum of $D_{\textrm{con}}$). However, this improvement is marginal in terms of $D_{\textrm{con}}$. If strong consensus exists ($D\approx 0$) both algorithms closely approximate this consensus. Even if the layers are independent ($D\approx 0.5$), we are able to find a consensus with $D_{\textrm{con}}<0.5$, especially for multiplexes with only a few layers. For the latter case, the consensus does not capture an intrinsic property of the multiplex network, rather it identifies nodes that by chance have high rank in many layers; and therefore potentially provide low-cost targets for influencing multiple layers simultaneously.

Note that Fig. 1 shows results for ER layers, we found that using the SF layers produces almost indistinguishable results. In fact, the performance of the algorithms only depends on the rankings $r_{\alpha}$, and not on the details of each layer $\alpha$. The only property of $r_{\alpha}$ that is specific to ER or SF layers is the number of ties in the ranking, and our results indicate that the performance of the BC and KA algorithms in terms of $D_{\textrm{con}}$ are insensitive to this.

II.2 Real networks

The BC and KA algorithms allow us to identify and analyze consensus ranking in real systems represented by multiplex networks. Here, we investigate three examples:

Airline network. A multiplex network with 5 layers representing the United States air traffic network in 2013. The nodes are airports and links indicate direct flights between them, and the layers correspond to the 5 largest carriers 222The airline network represents the five largest US carriers in terms of number of unique flights. It is constructed based on the publicly available database of the USA Department of Transportation (https://www.transtats.bts.gov/Fields.asp?Table_ID=259), all domestic flights labeled as “scheduled passenger service” are included.. 2. 2.

Primate social network. A multiplex network with 4 layers containing interaction data from one week of observations of a captive rhesus macaque troupe. The layers represent different interactions between animal pairs, including conflict, signaling of subordination, grooming, and huddling McCowan et al. (2011); Beisner et al. (2015). 3. 3.

Human social network. A multiplex network with 5 layers, where nodes are members of the Department of Computer Science at the Aarhus University and the links indicate various social relationships: Facebook friendship, spending leisure time together, working together, co-authorship, and regularly having lunch together Magnani et al. (2013).

Table 1 summarizes some basic properties of these networks. The airline network has a clear heavy-tailed degree distribution with hubs that have significantly more connections than average nodes. It is not possible to be as definite about the social networks due to their small size; there are, however, nodes that are connected to a significant fraction of the network, for example, in the fifth layer of the human social network the largest hub is connected to almost half of the other nodes. The primate social network is composed of competitive interactions (conflict and subordination) and cooperative interactions (grooming and huddling) – the existance of antagonistic and affiliative layers is a general property of social and political networks Harary et al. (1953); Wasserman and Faust (1994); Maoz (2010). For the primate network, we find that the adversarial layers are more heterogeneous than the affiliative layers. Interestingly, similar pattern was observed for human social networks obtained from an online game that allowed competition and alliances between players Szell et al. (2010).

When ranking nodes in a network, we chose a node centrality depending on the function of the system and our goals. Here, for illustration purposes, we calculate the node rankings for each network based on degree centrality. Figure 2(a-c) shows the pairwise Kemeny distances between layers for each network. Typically, we find $D<0.5$, indicating positively correlated rankings; the correlation is the strongest in the airline network, while rankings in the social networks are less aligned with each other. We observe an interesting pattern in the primate social networks: the distance is low between affiliative layers (grooming and huddling), and also low between the interactions involving social hierarchy (conflict and signaling); but the rankings corresponding to competitive layers are independent from rankings of affiliative layers.

We also identify the consensus ranking for each network using the BC and the KA algorithms. Figure 2d shows the average distance of the layers from consensus, $D_{\textrm{con}}$, for the original networks and their randomized counterparts, where we shuffled the node labels in each layer, eliminating any correlation between the rankings. For all multiplexes, we find stronger consensus in the real instances than in their randomized versions: the largest difference was observed in the airline network, while smallest in the primate network. As expected, KA always finds a consensus with a lower $D_{\textrm{con}}$ than BC; the difference, however, is marginal with the exception of the airline network.

Overall, we found for both model and real multiplexes that the BC and KA consensus rankings are similar in terms of $D_{\textrm{con}}$, while the BC algorithm is faster and extremely simple to implement. Both algorithms effectively identify the consensus if the layers are strongly correlated, and even if the layers of a multiplex are independent, we can identify a consensus ranking that provides low-cost targets for simultaneous intervention on multiple layers. So far, we only compared the consensus rankings in terms of $D_{\textrm{con}}$, in the next section we demonstrate how consensus ranking can identify targets for a simple degree-based node removal model with multiple objectives.

III Multi-objective degree-targeted attack

In this section, we explore an example of using consensus ranking to identify effective targets for multi-objective interventions in multiplex networks. A classic result of network science is that complex networks with heavy-tailed degree distributions are robust against random node removal; while targeted removal of high-degree hubs rapidly breaks down their large-scale connectivity Albert et al. (2000). Complex networks, however, rarely exist in isolation, nodes typically participate in multiple networks simultaneously. Therefore, removing a node from a multiplex network to reduce connectivity in one layer, might also remove connections from other layers that otherwise we would like to preserve. More specifically, here, we consider the problem of removing nodes from a multiplex network with $L$ layers such that for a set of layers $\mathcal{L}_{\textrm{D}}$ we aim to reduce the size of the largest component, while we want to keep the rest of the layers $\mathcal{L}_{\textrm{K}}$ intact. We rank the nodes in each layer based on degree centrality, but to reflect our different objectives, we reverse the rankings for the layers that we keep intact:

[TABLE]

where $r^{\textrm{D}}_{\alpha,i}$ and $r^{\textrm{K}}_{\alpha,i}$ are the ranks of node $i$ in layers that we want to destroy and keep, respectively; and $k_{\alpha,i}$ is the degree of node $i$ in layer $\alpha$. Given a multiplex network, we calculate these ranks for each layer, and we use the BC and KA algorithms to identify consensus rankings. We then iteratively remove the nodes based on this consensus ranking starting with the highest ranked nodes.

Figure 3 shows the relative size of the largest connected component $S(f)$ as a function of the fraction of nodes removed $f$ for model multiplex networks with $L=4$ ER or SF layers (for details about the model networks see Appendix B). For the ER example, we find that both the BC and KA consensus-based removal reduces $S$ faster for the layers we aim to destroy, and slower for layers we aim to preserve than random node removal, but BC preforms better at reducing $S$, while KA is better at preserving connectivity. Similarly for the SF example, BC destroys the targeted layers faster than KA, but doesn’t keep the rest of the layers intact more than random removal would. In either case, we are not as effective as if the layers would be in isolation.

Comparing the BC and KA rankings in Fig. 3, we find that neither method is objectively better than the other, instead they provide a different trade-off between the two objectives: BC preforms better at destroying layers, while KA does better at keeping layers intact. In fact, this is a pattern that we widely observed varying the parameters of the multiplex model networks. Therefore, from hereon in this section, we will focus on the BC ranking and we will study the trade-off by introducing a weighted version of BC. Further reasons for focusing on the BC algorithm are that it is much faster than the KA and the simplicity of the BC algorithm allows us to analytically solve $S(f)$ for the model networks.

The original BC algorithm provides one possible trade-off between the competing objectives of destroying certain layers, while keeping others intact. To explore other possible trade-offs, we introduce a weighted Borda count (wBC) algorithm that allows us to assign varying preference to the different objectives. The weighted Borda score of node $i$ is defined as

[TABLE]

where $w_{\alpha}=w$, if we destroy layer $\alpha$ ($\alpha\in\mathcal{L}_{\textrm{D}}$); $w_{\alpha}=1-w$, if we keep layer $\alpha$ intact ($\alpha\in\mathcal{L}_{\textrm{K}}$); and $w\in[0,1]$. The choice $w=0.5$ corresponds to the unweighted Borda score (with a multiplier of 1/2); if $w=1$, we only care about the layers we aim to destroy; and if $w=0$, we only care about the layers we aim to keep intact. In the following, we first derive an analytical solution for $S(f)$, and then we systematically investigate how multiplex network structure affects the trade-off between the different objectives using model and real networks.

III.1 Analytical solution

To analytically solve the size of the giant connected component (GCC) for consensus-based removal, we first calculate the size of the GCC for a general degree-based removal strategy (DBS) on a single layer network, and then we map the consensus-based process to a degree-based one for each layer. By DBS we mean a node removal process where the probability of removing a node only depends on its degree, i.e., the probability of removing node $i$ is $f(k_{i})$. Let $s$ be the probability that a randomly selected link leads to the GCC. Assuming local tree-like structure and uncorrelated networks, we calculate $s$ using the self-consistent equation

[TABLE]

where $p(k)$ is the degree distribution and $c$ is the average degree of the network. The first term on the right hand side is the probability that following a random link leads to a node that is not removed and not part of the GCC, and the second term is the probability that the node is removed. Once $s$ is determined, we obtain the relative size of the GCC using the equation

[TABLE]

Several node removal processes can be described as a DBS. For example, the choice $f(k)\equiv f$ leads to simple random node removal. Generally the best DBS to destroy a single-layer network removes nodes starting with the highest degree. Formally we express this as

[TABLE]

meaning that we remove all nodes with degree higher than $K$ and a certain fraction of nodes with degree $K$, where $K=\inf\{k:1-f\leq P(k)\}$ and $P(k)$ is the CDF of the degree distribution. Also, the most effective DBS to keep a single-layer network intact is to remove nodes with the lowest degree first:

[TABLE]

where $K=\inf\{k:f\leq P(k)\}$.

To map the consensus-based node removal to a DBS, we first calculate the Borda score of a node with multiplex degree $(k_{1},\ldots,k_{L})$ as

[TABLE]

where $P_{\alpha}(k)$ is the CDF of the degree distribution of layer $\alpha$. Next, we calculate the Borda score distribution, i.e., the probability that a randomly selected node of the multiplex has Borda score $b$:

[TABLE]

where $p(k_{1},\ldots,k_{L})$ is the multiplex degree distribution. When removing $f$ fraction of nodes based on consensus ranking, the probability of removing a node with Borda score $b$ is

[TABLE]

where $B=\inf\{b:1-f\leq P_{\textrm{B}}(b)\}$ and $P_{\textrm{B}}(b)$ is the CDF corresponding to $p_{\textrm{B}}(b)$ provided by Eq. (14). For each layer $\alpha$, the consensus-based removal can be formulated as an effective DBS with

[TABLE]

where $p(k_{1},\ldots,k_{L})$ is the joint multiplex degree distribution, and the summations are over the degrees in all layers, except $\alpha$.

Substituting Eq. (16) into Eqs. (9) and (10) provides $S_{\alpha}(f)$, the relative size of the GCC of layer $\alpha$ for consensus-based node removal. We numerically evaluate these equations for a class of multiplex model networks where the layers are either ER or SF networks, for SF networks we use the degree distribution provided in Eq. (23). Figure 4 compares the analytical solution of $S_{\alpha}(f)$ to simulations, showing excellent agreement. We find that there is a non-trivial connection between $S_{\alpha}(f)$ and $f$ as the node removal process switches back-and-forth between destroying and preserving layers. Note that the methods that we used to numerical evaluate the necessary equations become intractable for increasing number of layers; therefore for the numerical solutions we restrict ourselves to $L=4$ ER and $L=2$ SF layers.

III.2 Model networks.

The $w$ parameter of the wBC algorithm allows us to assign different level of importance to the two competing objectives: $w\approx 1$ favors breaking down layers, while $w\approx 0$ focuses on keeping layers intact. To systematically investigate this trade-off, we introduce the coefficients

[TABLE]

where $S_{\alpha}$ is the size of the GCC in layer $\alpha$ after removing $f$ fraction of its nodes using the wBC algorithm; $S^{\textrm{min}}_{\alpha}$ ($S^{\textrm{max}}_{\alpha}$) is the size of the GCC after removing $f$ fraction of its highest (lowest) degree nodes, i.e., it is the best we could do if layer $\alpha$ was in isolation. The coefficients $O_{\textrm{D}}$ and $O_{\textrm{K}}\in[0,1]$ measure how well the objectives are achieved relative to the case when we remove nodes from each layer independently. Analytical solution of Eq. (17) is obtained by first calulcating $S^{\textrm{min}}_{\alpha}$, $S^{\textrm{max}}_{\alpha}$, and $S_{\alpha}$ using Eqs. (11), (12), and (16).

Figure 5 shows the trade-off curve between $O_{\textrm{D}}$ and $O_{\textrm{K}}$ for model multiplex networks with $L$ layers, where we aim to destroy half of the layers, while keeping the other half intact. If there would be no trade-off, the curve would be a single point at $(O_{\textrm{K}},O_{\textrm{D}})=(1,1)$; if all layers would be identical, the trade-off curve would be the diagonal line $O_{\textrm{D}}=1-O_{\textrm{K}}$.

Figure 5a shows that increasing the number of layers affects the trade-off in two distinct ways: (i) As $L$ increases the number of conflicting objectives also increase, and the trade-off becomes more severe approaching the diagonal line. (ii) For $L\geq 4$, we aim to destroy multiple layers simultaneously, which cannot be done as efficiently as if they where in isolation; therefore even if $w=1$, $O_{\textrm{D}}$ remains less than 1. Increasing $L$ has a similar, but weaker effect on $O_{\textrm{K}}$, since destroying a network is more difficult than keeping it intact. Changing the pairwise Kemeny distance $D$ between the layers (Fig. 5b), we find that strong consensus ($D\approx 0$) leads to strong trade-off between $O_{\textrm{D}}$ and $O_{\textrm{K}}$, but allows efficient destruction ($\max_{w}O_{\textrm{D}}\approx 1$) or preservation ($\max_{w}O_{\textrm{K}}\approx 1$) of layers. Figure 5c shows the effect of degree heterogeneity, we find that destroying SF layers and preserving ER layers entails the least amount of trade-off, while destroying SF and preserving ER is the most difficult. Finally, Fig. 5d shows that the trade-off is the most pronounced at the initial stages of node removal (small $f$), and the competing objectives are less restrictive in case of large-scale removals (large $f$).

III.3 Real networks.

In the previous section, we used model multiplex networks to understand how basic network properties affect the $O_{\textrm{D}}$-$O_{\textrm{K}}$ trade-off. Real networks, however, have more complex structure. To investigate their effect, we use the three example datasets introduced in Sec. II.2: the human social network, the primate social network, and the airline network. We assign half of the layers in each multiplex to be destroyed and the rest to be kept intact, and we calculate the $O_{\textrm{D}}$-$O_{\textrm{K}}$ trade-off after removing a fraction of the nodes ($f=0.4$ for the airline network and the human social network, and $f=0.55$ for the primate social network). We then compare the trade-off curves to the following randomized null models:

Full randomization (FR): We randomly rewire all links, in effect replacing each layer with an ER network with the same number of nodes and links. 2. 2.

Node label randomization (NLR): We shuffle the node labels in each layer, removing any correlation or consensus between layers, but leaving the structure of each layer unchanged. 3. 3.

Degree preserved randomization (DPR): We rewire links such that the multiplex degree of each node is unchanged. This randomization, therefore, preserves the pairwise Kemeny distances between the layers, and removes all structure within the layers beyond their degree sequence. 4. 4.

Node label and degree preserved randomization (NLR+DPR): Combining NLR and DPR removes correlations both within layers and between layers, and only preserves the degree distributions of the individual layers.

Figure 6 provides the trade-off curves for the three example networks, each network showing a distinct behavior. For the airline network, we aim to destroy two layers and keep the other three intact. In Sec. II.2, we showed that the airline network has heterogeneous degree distributions and strong consensus between its layers. For model networks, we found that strong consensus leads to significant $O_{\textrm{D}}$-$O_{\textrm{K}}$ trade-off, and indeed, we find that for the airline network the trade-off curve is close to the diagonal. We find that the DPR trade-off curve is almost indistinguishable from the original, while the trade-off for FR, NLR, and NLR+DPR is significantly weaker, meaning that both the inter-layer correlations and the degree distributions are necessary to explain the strong trade-off. We also find that $\max_{w}O_{\textrm{D}}=1$ for the original and NLR curve, but not for FR, this indicates that the simultaneous destruction of the two layers is aided by the presence of hubs and is uneffected by the inter-layer correlations.

In the case of the primate social network, we aim to destroy the two layers that are related to competition (conflict and signaling) and to preserve the affliative layers (grooming and huddling). In fact, it is common practice for the caretakers of captive primate groups to remove unusually aggressive individuals with to goal of reducing overall levels of conflict; and it is desirable that this removal doesn’t negatively affect the cohesion of the group Judge et al. (1994); McCowan et al. (2018); Pósfai and D’Souza (2018). Figure 6b shows the $O_{\textrm{D}}$-$O_{\textrm{K}}$ curves for the primate network, and we find a very weak trade-off. This is largely explained by the inter-layer correlations, which in contrast to the airline network, reduce the trade-off. To understand this recall Fig. 2b, where we showed that Kemeny distance is small between the conflict-signaling and the grooming-huddling layer pairs, while there is no or even negative correlation between the competitive and affiliative layers. This particular structure allows us to simultaneously disrupt the competitive layers without affecting the affiliative layers. Furthermore, the competitive layers have more heterogeneous degree distributions than the affiliative layers (Table 1), a property also seen in human social networks Szell et al. (2010). Using model networks, we showed that destroying heterogenous layers and keeping homogeneous layers intact reduces the trade-off (Fig. 5c). Indeed, comparing the FR and DPR+NLR curves shows that the degree distribution of the layers contributes to the weak trade-off, albeit less than the inter-layer correlations.

Finally, for the human social network, our goal is to reduce connectivity in two layers and keep the other three intact. For both the airline network and the primate social network, we found that the inter-layer correlations significantly affect the $O_{\textrm{D}}$-$O_{\textrm{K}}$ trade-off; interestingly, in case of the human social network inter-layer correlations have little effect. The original trade-off curve is best approximated by the NLR, which removes correlations between layers, but preserves all structure within each layer; and all other randomizations show stronger trade-off than the original. Therefore, we conclude that the internal structure of the layers beyond the degree distribution, such as community structure, is what reduces trade-off.

IV Conclusion

In this paper, we explored the use of consensus rankings to identify effective targets for multi-objective interventions in multiplex networks. Our strategy is to calculate a centrality for each layer that is relevant to its specific function and our specific objectives, rank the nodes based on these centralities, and then using rank aggregation methods identify a consensus ranking. As an example process, we studied the degree-based node removal process, where we aimed to destroy the largest connected component in some layers, while keeping the other layers intact. We demonstrated that removing the nodes in order of consensus ranking effectively balances these competing objectives.

The advantage of our method is that it is agnostic to the specific properties of the layers and our goals, making it a widely applicable tool. However, the price of this flexibility is that methods designed for a specific system likely outperform our general approach – future work should explore this possible trade-off. In our work, we extensively investigated how inter-layer structural correlations affect the trade-off between different objectives. We have not, however, explored the scenario when the dynamics of the layers are also directly coupled. Centrality measures that take such coupling into account were only developed for multilayer networks where all layers follow qualitatively similar dynamics De Domenico et al. (2015b); Solé-Ribalta et al. (2016). It would be interesting to extend rank aggregation techniques to directly consider coupled dynamics.

Acknowledgements.

We thank Haochen Wu for the airline network dataset. We gratefully acknowledge support from the US Army Research Office MURI Award No. W911NF-13-1-0340, DARPA Award No. W911NF-17-1-0077, and the National Institutes of Health Award No. R24-OD011136.

Appendix A Local search algorithm

Many methods exist to approximate the Kemeny consensus, a recent survey compared the performance of 104 algorithms and combinations of algorithms Ali and Meilă (2012). It found that so-called local search methods provide an optimal trade-off between accuracy and run-time, meaning that algorithms with significantly longer runtime only marginally decreased the cost function provided by Eq. (3). Local search algorithms start from an initial ranking which can be either random or an approximation of the Kemeny consensus provided by another algorithm. Then this ranking is improved on by a series of local transformations that decrease the cost function. Here, we implement a version of local search based on the simple insert sort algorithm.

As initial ranking we use $R_{\textrm{BC}}$, the output of the Borda count algorithm. Let node $i_{1}$ be the top ranked node in $R_{\textrm{BC}}$, node $i_{2}$ the second, and so on; if there are ties in $R_{\textrm{BC}}$, we randomly break them. We then construct a new ranking $R_{\textrm{KA}}$ by starting from an empty ranking and iteratively adding nodes. First, we add node $i_{1}$. We then add $i_{2}$ such that it has new rank (i) above $i_{1}$, (ii) tied with $i_{1}$, or (iii) ranked below $i_{1}$. The new rank of $i_{2}$ is chosen to minimize the cost function. We repeat this step until all nodes are assigned a new rank.

Appendix B Centrality correlated model networks

We introduce a method to generate multiplex networks with $N$ nodes, $L$ layers and tunable Kemeny distance $D_{\alpha\beta}$ between each pair of layers $\alpha,\beta\in\{1,2,\dots,L\}$. We start by generating $L$ independent layers using any single-layer network model of choice. We then calculate the rank $r_{\alpha,i}$ of each node $i$ in each layer $\alpha$ based on some centrality $c_{\alpha,i}$. Note that we are not restricted to use the same single-layer network model or the same centrality for all layers. In the following, we describe a procedure to re-label the nodes in each layer to specify the inter-layer dependency between the ranks of nodes.

Let $P(c_{1},\dots,c_{L})$ be the CDF of the multiplex centrality distribution, i.e., $P(c_{1},\dots,c_{L})$ is the probability that a randomly selected node $i$ has $c_{\alpha,i}\leq c_{\alpha}$ for every $\alpha$; and let $P_{\alpha}(c_{\alpha})$ be the CDF of the marginal centrality distribution of layer $\alpha$. Specifying the multiplex centrality distribution would allow us to control the dependency between layers; we, however, cannot arbitrarily choose $P(c_{1},\dots,c_{L})$, since the marginals $P_{\alpha}(c_{\alpha})$ are determined by the properties of the individual layers. Yet, according to Sklar’s theorem, we can always write the multiplex centrality distribution as

[TABLE]

where $C(u_{1},\ldots,u_{L})$ is an $L$-variate copula Sklar (1959). A copula $C:[0,1]^{L}\rightarrow[0,1]$ is the CDF of a random vector $(u_{1},\ldots,u_{L})$ with uniform margins Joe (1997); Nelsen (2006). The advantage of this representation is that it separates the marginal distributions of each layer, specified by $P_{\alpha}(c_{\alpha})$, and the interdependency structure between layers, specified by $C(u_{1},\ldots,u_{L})$. If the marginals $P_{\alpha}(c_{\alpha})$ and $P_{\beta}(c_{\beta})$ are continuous, the Kemeny distance $D_{\alpha\beta}$ between the two layers is the same as the Kendall’s distance and it is completely defined by the copula

[TABLE]

where $C(u_{\alpha},u_{\beta})$ is the two-dimensional marginal of the $L$-variate copula Nelsen (2006). For example, for independent layers we have $C(u_{1},\ldots,u_{L})=\prod_{\alpha=1}^{L}u_{\alpha}$, and substituting this into Eq. (19), we get $D=0.5$ for all layer pairs.

To re-label the nodes such that the multiplex centrality distribution $P(c_{1},\dots,c_{L})$ follows Eq. (18), we draw a random vector $u_{i}=(u_{1,i},\ldots,u_{L,i})$ from the $L$-variate copula for each node $i$, and from this we obtain a vector of ranks $s_{i}=(s_{1,i},\ldots,s_{L,i})$ where

[TABLE]

We re-label each node such that $r_{\alpha,i}$, the rank of $i$ in layer $\alpha$, becomes equal to $s_{\alpha,i}$. For example, if $s_{i}=(1,3)$, we re-label nodes such that the highest ranked node in layer 1 and third highest ranked node in layer 2 are labeled $i$.

Throughout this paper, we use the Gaussian copula

[TABLE]

where $\Phi_{\rho}$ is the CDF of the L-variate standard normal distribution with correlation matrix $\rho\in[-1,1]^{L\times L}$ and $\Phi^{-1}$ is the inverse CDF of the standard normal distribution. Substituting the Gaussian copula into Eq. (19), we find the relationship

[TABLE]

where $\rho_{\alpha\beta}$ is an element of the correlation matrix $\rho$. Therefore, choosing the correlation matrix $\rho$ allows us to set the Kemeny distance between layers.

Finally, note that the above procedure does not take into account ties in the ranks, i.e., we assume that $c_{\alpha,i}\neq c_{\alpha,j}$ for any $i\neq j$. Therefore, if there are any ties in the centralities (e.g., if degree centrality is used), we randomly break them. Furthermore, Eq. (19) assumes that the marginal centrality distributions are continuous, if this is not the case, the Kemeny distance $D_{\alpha\beta}$ between layers is no longer exactly provided exactly by Eq. (19) and (22); through simulations, however, we found that $D_{\alpha\beta}$ is still well approximated by them.

B.1 Static scale-free model

To generate the scale-free layers, we use the static model Goh et al. (2001). Starting from $N$ unconnected nodes, we assign a weight $w_{i}=i^{-\alpha}$ to each node $i=1,\ldots,N$, where $\alpha\in[0,1)$. We then randomly select two nodes $i$ and $j$ with probability proportional to $w_{i}$ and $w_{j}$, respectively, and if there is no link between nodes $i$ and $j$, we connect them. We repeat this step until $L$ links are added. The resulting network has average degree $c=2L/N$ and its degree distribution can be written as sum of Poisson distributions

[TABLE]

where $c_{i}=2Lw_{i}/(\sum_{j=0}^{N}w_{j})$ is the expected degree of node $i$. For large $N$, the degree distribution is approximated as

[TABLE]

where $\Gamma(z)$ is the gamma function and $\Gamma(z,a)$ is the upper incomplete gamma function. We refer to this network as scale-free, because the tail of the distribution decays as a power-law, i.e., $p(k)\sim k^{-\gamma}$, where $\gamma=1+1/\alpha$.

When numerically evaluating the equations developed in Sec. III.1, we use Eq. (23) to represent $p(k)$.