Multiplex networks with heterogeneous activities of the nodes
Abstract
In multiplex networks with a large number of layers, the nodes can have different activities, indicating the total number of layers in which the nodes are present. Here we model multiplex networks with heterogeneous activity of the nodes and we study their robustness properties. We introduce a percolation model where nodes need to belong to the giant component only on the layers where they are active (i.e. their degree on that layer is larger than zero). We show that when there are enough nodes active only in one layer, the multiplex becomes more resilient and the transition becomes continuous. We find that multiplex networks with a powerlaw distribution of node activities are more fragile if the distribution of activity is broader. We also show that while positive correlations between node activity and degree can enhance the robustness of the system, the phase transition may become discontinuous, making the system highly unpredictable.
I Introduction
Multilayer networks (1); (2); (3) are formed by several interacting networks. They describe a large variety of complex systems. Examples are found in social (4), technological, communication, transportation systems (5); (6); (7) but also in biological networks of the cell, or in the brain (9); (10); (11). In the last fifteen years a lot of attention has been devoted to understand the interplay between structure and dynamics in single networks (12); (13). In recent years, it has become clear that most of the networks are not isolated and that for understanding the function of complex systems as different as complex infrastructures or the brain it is important to investigate the role of their multilayer structures (1); (2). In particular, ample debate has been devoted to characterize the robustness of multilayer networks (14); (17); (20); (21); (22); (16); (15); (18); (29); (19); (25); (26); (30); (23); (24); (10); (27); (28). It has been shown that, in the presence of interdependencies, the robustness of a multilayer network can be significantly affected. Multilayer networks can be much more vulnerable to random damage with respect to considering only their single layers taken in isolation (14); (16); (15).
In the presence of interdependencies, the notion of mutually connected component has been introduced, meaning that each pair of nodes in the mutually connected component must be connected by a path on each and every layer, internal to this component. This definition is motivated by the fact that in these interdependent systems a node is not functional if any of its interdependent nodes in the other layers is not functional. Therefore, the largest (giant) mutually connected component (MCGC) describes the robustness of the system. This component has a discontinuous phase transition as a function of the initial random damage inflicted to the nodes of the network, and close to this transition the system is affected by dramatic cascades of failure events. This transition can be studied on multilayer networks of different nature, including multiplex networks (14); (17); (16); (15) and network of networks (20); (21); (22); (23); (24). Network of networks are formed by different interacting networks, as the molecular networks in the cell, in which every node is a different type of biological molecule. Multiplex networks (6); (31); (7); (32); (33); (36); (34); (35), instead, are multilayer networks in which the same set of nodes interact through different networks. Example of multiplex networks are social networks, in which individual are connected by different types of social ties, transportation networks, like airport networks in which each airport is connected to other airports though connections of different airline companies, or collaboration networks in which scientists collaborate on different topics and eventually cite each other. Several works have focused on modeling multilayer networks (32); (33); (36); (34); (35). In particular a very useful approach employs statistical ensembles which are able to generate a large variety of multiplex network topologies with desired level of structural correlation (32); (34); (35). In networks of networks as well as in multiplex networks the possible correlations existent in this structure can strongly affect their robustness properties (1); (18); (11); (10).
An interesting result has recently shown that multiplex networks are characterized by the fact that not all the nodes are connected in every layer. In fact, many networks have been shown (7) to have heterogeneous activity of the nodes. The activity of a node is given by the number of layers in which the node is at least connected to another node. The activity of the nodes has been found to be broadly distributed in a variety of multiplex network data sets (7). Moreover, the activity of the nodes has been seen to correlate in average with the degree of single node in given layers, meaning that on average nodes that are present in many layers have also typically high degree in the single layers in which they are active (7). The heterogeneous distribution of activities of the nodes and its relation with the mean degree in single layers is frequently observed in realworld multiplex networks and may encode relevant information (7). Besides, multiplex networks with heterogeneous activities can be useful theoretical tools to address other problems like modular structures in single networks (8).
Here we characterize multiplex networks ensembles with heterogeneous activities of the nodes and we investigate their robustness properties, assuming mutual interdependencies among each layer. This ensemble of multiplex networks with heterogeneous activity of the nodes can be used to model realistic multiplex network structures, on top of which dynamical processes may occur. Moreover, ensemble of multiplex networks with heterogeneous activities of the nodes can be used to estimate the role that correlation have on their robustness properties. Here we find that heterogeneous activity of the nodes can decrease the robustness of networks in the presence of interdependencies. Nevertheless the correlation between the activity of the nodes and their degree within single layers has the opposite effect and can improve the robustness of multiplex networks.
Ii The multiplex network ensemble with heterogeneous activity of the nodes
A multiplex network of nodes and layers is completely specified when the adjacency matrices of elements if node is linked to node in layer , and otherwise are given. Every node of the multiplex network is labelled as , indicating that is the th node in layer . The replica nodes of the node are defined as all the nodes labelled as in layers (23). Interestingly, it has been observed from data (7), that in many networks not all the nodes are active (i.e. are connected to at least another node) in each layer. Let us define the activity of a node as the number of layers where node has nonzero degree. The activities of the nodes are broadly distributed (7) and they can be fitted by a powerlaw with . This implies that for some multiplex networks the bipartite network between nodes and layers described by the activity adjacency matrix can be either dense or sparse but the typical number of layers in which a node is active is always subject to unbound fluctuations. In order to characterize fully the activities of the nodes in each layer, from the adjacency matrices it is possible to construct a activity matrix of elements (Fig. 1). This matrix can be viewed as an adjacency matrix between nodes and layers indicating if node is active in layer ( or not ().
The activity of a node can be therefore expressed in terms of the matrix as in the following,
(1) 
The layer activity has been defined in (7); (1) and is given by the number of nodes active in layer , i.e.
(2) 
It is therefore possible to construct ensemble of multiplex networks with heterogeneous activity of the nodes in the subsequent manner: first we construct the network between layers and nodes described by the adjacency matrix indicating if node is active in layer , then we can construct in each layer a network between the active nodes of the layer with a degree distribution . In the case in which , we have that a node is active in a given layer if it is connected at least to another node in the same layer, i.e.
(3) 
where indicates the Kronecker delta.
In order to construct the activity matrix we can consider a microcanonical ensemble in which the activities of the nodes and the layer activities are fixed. This ensemble, called also the configuration model of a bipartite network, can be constructed by maximizing the entropy of the activity networks given by
(4) 
where is the probability of a given activity network in the ensemble under the constraints that the node activity and layer activity are kept constant. In this ensemble the probability of a matrix is given by
(5) 
where is the normalization constant also called the partition function of the ensemble. In this ensemble the probability that a node is active in layer is expressed in terms of the Lagrangian multipliers , , i.e.
(6) 
The Lagrangian multipliers are fixed by the conditions
(7) 
Finally the entropy (38); (39); (40) of the ensemble of activity networks is given by
(8) 
where and are given by
We assume here that bipartite networks between nodes and layers, characterized by the adjacency matrix is uncorrelated, i.e. we assume that the activities of the nodes are not correlated with the sizes of active nodes on the layers where the nodes are active. In this hypothesis the probabilities are proportional to the product of and which are the degrees of the mentioned bipartite network, and the following condition on the the maximal activity of the nodes and the maximal layer activity must be satisfied, i.e.
(9) 
This condition is the condition that in a bipartite network corresponds to the one imposing a structural cutoff on single uncorrelated networks (networks without degreedegree correlations) (37). In this case we have the simple expression
(10) 
with .
As we said, this argument regards lack of correlation between activities and number of active nodes on the layer. In Sec. IV, we will examine a different type of correlation, namely the one between activities and the degree distribution on the layer.
Once the activity network is constructed, in order to construct the multiplex network, we assign to each active node of layer a degree within the layer. If the activity of the nodes are uncorrelated with the degree in each layer the degrees of the nodes in layer are drawn form the degree distribution and the networks of each layer are generated by the configuration model with degree distribution . In other words, the probability of the multiplex network degree sequences is given by
(11) 
Instead, if there is a correlation between the degree of the nodes within a layer and the activity of the nodes, then we need to draw the degree of the nodes in each layer from a probability which is a function of their activity .
Iii Mutually connected component in a multiplex network with given distribution of activities of the nodes
We consider the mutually connected giant component (MCGC) in a multiplex network with given distribution of activities of the nodes, described by the ensemble of multiplex networks introduced in section II.
The layers are interdependent, meaning that each node active in a given layer is interdependent on its replica nodes in all the layers where the node is active. In particular we will assume that a node active in layer belongs to the mutually connected giant component if

at least one neighbor of node belongs to the mutually connected giant component;

in each layer where the node is active, at least one neighbor belongs to the mutually connected giant component.
On a locally treelike multiplex network, this can be easily encoded in a message passing algorithm determining if a node belongs to the mutually connected giant component (1); (30); (16); (41). The indicator function indicates if a nodes belongs () or not () to the mutually connected giant component. This indicator is determined by a set of “messages” that each node active in a layer sends to the neighboring nodes in the same layer. Each message indicates if the node belongs () or not () to the mutually connected component when the link to the node is removed.
Here our goal is to characterize the size of the mutually connected giant component as a function of the probability that random nodes are damaged in the network. In order to characterize the damage initially inflicted to the network we indicate with if a node has been damaged () or not () in the multiplex network. Given the definition of the mutually connected component, the message from a node to a node both active in layer is therefore given by
(12)  
where indicates the nodes that are neighbors of node in layer and the expression indicates all the set of all the nodes belonging to except node . Therefore the message is equal to one, if the node has not been damaged (), if at least a message coming from a neighbor node different from is positive (first factor in the square brackets), and if in all the layers where node is active (), there is at least one neighbor of node sending a positive message to node (product over ). Finally, the indicator function that node active in layer , is in the MCGC is given by
(13)  
This indicator function equals one, if the node has not been damaged (), if at least a message coming from a neighbor node of layer is positive, and if in all the layers where node is active (), there is at least one neighbor of node sending a positive message to node . In this section we only consider uncorrelated activity networks in order to allow for a theoretical treatment of the percolation properties of the multiplex networks with heterogeneous activity of the nodes. Moreover we assume that the number of nodes is large in every network of the multiplex network. Given that from Eq. we have
(14) 
we have
(15) 
where the condition is required to study the percolation transition of the MCGC. Since the value of the activity of the nodes we find the condition
(16) 
This means that, for example, for finite the number of layers is much smaller than the total number of nodes.
In order to characterize how the size of the mutually connected giant component depends on the distribution of activities , here we consider the ensemble of multiplex networks with given activities of the nodes described in section II, where we assume additionally that the damage occur on a node with probability , i.e.
(17) 
Let us observe now that the probability that a node active in layer has activity is given by
(18)  
where we have assumed that is given by Eq. . Assuming that all the layers have the same activity , and the same degree distribution , by averaging the messages given by Eq. over the ensemble of multiplex networks with given activities of the nodes, we get the probability that by following a link in layer we reach a node in the mutually connected component, obtaining and
(19) 
where the generating functions and are given by
(20) 
Which can be derived as in the following. First we define as
where indicates the average over the network, in the large network limit. Using Eq. we have
(21)  
getting Eq. in the treelike local approximation.
Using a similar derivation, under similar assumptions, the fraction of the nodes that are in the mutually connected component in layer is given by
(22) 
These equations can be studied in detail as a function of the degree distribution in each layer, and the node activities distribution . In particular, Eq. can be expressed as,
(23) 
where we have indicated by the following generating function
(24) 
Similarly can be written as
(25) 
Let us consider now Poisson networks identical on each layer: .Then we have and the equation for becomes where
(26) 
and , . A discontinuous transition occurs only if has maximum, so we need to impose the condition , where
where the function is given by
(27) 
Using this equation to isolate , and substituting into , we get for
(28) 
As a general remark, we expect in this model that the nature of the percolation phase transition will be dependent on the fraction of replica nodes with activity . In fact these nodes do not have any interdependency on replica nodes in other layers. Therefore a high density of nodes with activity favors the continuous phase transition. In the extreme case in which all the nodes have activity Eq. (26) is in fact reduced to the equation determining the giant component of a single network. In fact, the generating function defined in the last equation of is equal to one if all the nodes have activity . When the not all the nodes have activity , the non trivial functional form of is determined by the nodes that have activity . Therefore we expect that the density of such nodes is sufficiently high the transition is continuous while otherwise, we expect a discontinuous phase transition. This expected phenomenon can be related with a similar effect generated by the partial interdependence in multiplex networks where all the nodes are active in any layers (17); (16).
In the following we provide evidence for this observations by considering two different distributions of activities : a scalefree activity distribution and a Poisson activity distribution.
iii.1 Scalefree distribution
Let us consider a power law activity distribution with . In this case, the function is given by Eq. (26) where
(29) 
The position of the discontinuous transition can be calculated by solving , with . The discontinuity of the transition is due to the finite value of the solution at , as one can see that from Eq. (25), where is a simple exponential. The position of the continuous transition can be instead calculated by solving . While the tricritical point, if it exist, can be found by setting . (Note, that given the functional form of , given by Eq. , the condition is always satisfied.)
Fig. 2 illustrates the phase diagram calculated for and . For we have a continuous phase transition for . For we find a discontinuous phase transition line across all the values of . The jump in the size of the MCGC at the discontinuous transition for (shown in Fig. 3) diverges as approaches the value 2 (while at the continuous transition there is no jump). Both transition lines for and diverge as approaches the value , indicating that as decreases the multiplex network becomes more fragile. As decreases, the fraction of nodes with large activity (activity hubs) becomes important. These nodes, to be in the mutually connected giant component they must be connected in each layer, and if they are damaged they affect the connectivity of multiple layers. Therefore they are at the same time more fragile than nodes with smaller activity and are able to affect the multiplex connectivity across more layers. As a result, as decreases, random damage affects more and more layers and the multiplex network becomes less robust.
iii.2 Poisson distribution
For activities that are Poisson distributed , we have
(30) 
By setting , we can find a line of discontinuous transitions that ends at a critical point defined by characterized by
(31) 
In the multiplex network, when the percolation transition is continuous. Instead for the percolation transition is discontinuous. In Figure 4 we show the phase diagram. At , the continuous phase transition is at as expected for independent layers. As the level of interdependence becomes increasingly significant, the continuous transition becomes discontinuous (for ). This can also be seen from the size of the discontinuous jump , that approaches zero at the tricritical point and is shown in Figure 5.
Iv Correlations between the activities and the degrees of the nodes
In this section, we consider the role of correlations between the activities of the nodes and their degrees in the layer in which they are active. The higher is the activity of a node, the more fragile is the node, but the higher is the degree of the node within the layer the more robust it is. Therefore, here we want to characterize the effect that the correlations between the activity and the degrees of each node have on the robustness properties of the entire multiplex network. We start from the message passing Eqs. and . We assign to each node an activity from a distribution. Then, we assume that the matrix of activities is given by with probability given by Eq. with given by Eq. . Finally, we need to extend Eq. (11) for the degree sequence to the correlated case. As node degrees on each layer are correlated with node activities, the probability that a multiplex network has degree sequences , given the activity matrix and the activity sequence , is
In particular, we assume for simplicity that
(33) 
and we take
(34) 
where are two parameters determining the correlations between the degrees of the node and its activity. The particular choice of the functional form of in Eqs.  is dictated by the intention to model positive correlations between the activity of the degree of the layers that have been observed in real dataset (7). Real multiplex network analysis are nevertheless not sufficient to suggest the exact form of the correlations observed. Therefore this ensemble of correlated multiplex networks with heterogeneous activity of the nodes has been chosen in such way describe positive correlations between activities and degree in single layers, while keeping the model sufficiently simple to allow a number of analytical calculations.
Averaging the message passing equation over this ensemble and assuming , , and , we have only one average message determined by the equation
(35) 
Setting , the equation for reads
(36) 
with
This equation can be studied as a function of , and the parameters determining the distribution. In order to find the discontinuous phase transition, we set , . The line of these discontinuous phase transition eventually stops at a critical point , that can be calculated by setting
(37) 
The continuous phase transition can be found, instead, by imposing .
Let us now consider a powerlaw distribution of activities for and correlated degrees in every layer according to the degree distribution given by Eq. (33). In Figure 6 and Figure 7 we show the percolation transitions as a function of the parameter for and , respectively. As a general remark, for , the transition appears to be always discontinuous, while for we have both continuous and discontinuous transitions, with a region of coexistence between two percolating phases, and a critical point at the end of the line of discontinuous transitions.
First, let us observe the case of (Fig. 6). We do not plot the line as it reduces to the noncorrelated case shown in Fig. 2, where the phase transition is always continuous. For , instead, it emerges a line of discontinuous phase transitions ending in a critical point, which is determined by the Eqs. . The line of continuous phase transitions encounters the line of discontinuous transitions at a critical end point determined by the simultaneous occurrence of two minima of function (Eq. IV). Fig. 8 shows the solution of equation (IV) at the critical point. This phase diagram, therefore, is characterized by the presence of two percolating phases, separated by a short line of discontinuous transitions (Fig. 6). This coexistence region shifts towards larger values of as increases, but the critical point never joins the continuous line at a tricritical point at finite . Figure 9 shows that , without joining the continuous solution at , for .
In the case (Figure 7), we observe a discontinuous phase transition for all values of . This qualitative difference is due to the peculiar role of the nodes which are active on a single layer (). In our model, this type of nodes do not need support from the other layers (as they appear only in one of them), and therefore they drive a classical continuous percolation transition. In both cases, as increases, the percolating threshold becomes smaller, as the intensity of the correlations between the activities of the nodes and the degree of the nodes within each layer increases. Therefore the multiplex networks is more robust because it can still have a mutually connected component even if the damage is significant. On the other hand, though, the phase transition becomes discontinuous, and therefore the collapse becomes more unpredictable.
V Conclusions
In this paper, we have characterized the robustness properties of multiplex networks in a new model that encapsulates heterogeneous activities of the nodes, i.e. the possibility that each node is present only on a small fraction of layers in a multiplex, as seen in real world cases (7). In this model, we employ a notion of mutual percolation where nodes must belong to a mutually connected component only on the layers where they are active and develop an analytical approach to calculate the size of the mutually connected giant component as a function of node damage and other parameters. We show that multiplex networks with very broad activity distributions are more fragile than networks with more homogeneous distribution of activities.
We also investigate the role of correlations between the activities of the nodes and their degrees. We show that these correlations generally improve the stability of the percolating phase,and the multiplex network has a smaller percolation threshold, so the multiplex network becomes more robust. However, correlations also change the order of the phase transition, that becomes discontinuous. This provides an example on how correlations can typically reduce the fragility of multiplex networks, but at the same time they can make the system more unpredictable, as the transition becomes discontinuous.
Vi Acknowledgments
We acknowledge useful discussions with Vincenzo Nicosia and Vito Latora. This work has been partially funded by: Science Foundation Ireland, grants 14/IF/2461 and 11/PI/1026; the FETProactive project PLEXMATH (FP7ICT20118; grant 317614).
References
 S. Boccaletti, G. Bianconi, R. Criado, C. I. del Genio, J. GómezGardeñes, M. Romance, I. SendiñaNadal, Z. Wang, and M. Zanin, Phys. Rep. 544, 1 (2014).
 M. Kivelä, A. Arenas, M. Barthelemy, J. P. Gleeson, Y. Moreno, and M. A. Porter, J. Complex Netw. 2, 203 (2014).
 K.M. Lee, J. Y. Kim, S. Lee, and KI. Goh. ”Multiplex networks.” In Networks of networks: the last frontier of complexity, pp. 5372. (Springer International Publishing, 2014).
 M. Szell, R. Lambiotte, and S. Thurner, PNAS, 107, 13636 (2010).
 P. J. Mucha, T. Richardson, K. Macon, M. A Porter, and J.P. Onnela, Science, 328, 876 (2010).
 A. Cardillo, J. GómezGardeñes, M. Zanin, M. Romance, D. Papo, F. del Pozo, and S. Boccaletti, Sci. Rep. 3, 1344 (2013).
 V. Nicosia and V. Latora, Phys. Rev. E 92, 032805 (2015).
 P. Colomerde Simón and M. Boguñá, Phys. Rev. X 4, 041020 (2014).
 E. Bullmore and O. Sporns, Nat. Rev. Neurosci. 10, 186 (2009).
 S. D. S. Reis, Y. Hu, A. Babino, J. S. Andrade Jr., S. Canals, M. Sigman, and H. A. Makse, Nature Phys. 10, 762 (2014).
 G. Bianconi, Nature Phys. 10, 712 (2014).
 S. N. Dorogovtsev, A. Goltsev and J. F. F. Mendes, Rev. Mod. Phys. 80, 1275 (2008).
 A. Barrat, M. Barthélemy, A. Vespignani Dynamical Processes on complex Networks (Cambridge University Press, Cambridge, 2008).
 S. V. Buldyrev, R. Parshani, G. Paul, H. E. Stanley and S. Havlin, Nature 464, 1025 (2010).
 G. J. Baxter, S. N. Dorogovtsev, A. V. Goltsev, and J. F. F. Mendes, Phys. Rev. Lett. 109, 248701 (2012).
 S.W. Son, G. Bizhani, C. Christensen, P. Grassberger, and M. Paczuski, EPL 97, 16006 (2012).
 R. Parshani, S. V. Buldyrev, and S. Havlin, Phys. Rev. Lett. 105, 048701 (2010).
 B. Min, S. D. Yi, K.M. Lee, and K.I. Goh, Phys. Rev. E 89, 042811 (2014).
 K. Zhao and G. Bianconi, J. Stat. Mech. P05005 (2013).
 J. Gao, S. V. Buldyrev, H. E. Stanley, and S. Havlin, Nature Phys. 8, 40 (2012).
 J. Gao, S. V. Buldyrev, S. Havlin, and H. E. Stanley, Phys. Rev. E 85, 066134 (2012).
 J. Gao, S. V. Buldyrev, S. Havlin, and H. E. Stanley, Phys. Rev. Lett. 107, 195701 (2011).
 G. Bianconi, S. N. Dorogovtsev, and J. F. F. Mendes, Phys. Rev. E 91, 012804 (2014).
 G. Bianconi and S. N. Dorogovtsev, Phys. Rev. E 89, 062814 (2014).
 D. Cellai, E. Lòpez, J. Zhou, J. P. Gleeson, and G. Bianconi, Phys. Rev. E 88, 052811 (2013).
 G. J. Baxter, S. N. Dorogovtsev, J. F. F. Mendes, and D. Cellai, Phys. Rev. E 89, 042801 (2014).
 M.A. Serrano, L. Buzna, and M. Boguna, arXiv:1502.04553 (2015).
 S. Hwang, S. Choi, D. Lee, and B. Kahng, Phys. Rev. E 91, 022814 (2015).
 N. AzimiTafreshi, J. GómezGardeñes, and S. N. Dorogovtsev, Phys. Rev. E 90, 032816 (2014).
 S. Watanabe and Y. Kabashima, Phys. Rev. E. 89, 012808 (2014).
 F. Battiston, V. Nicosia and V. Latora, Phys. Rev. E 89 032804 (2014).
 G. Bianconi, Phys. Rev. E 87, 062806 (2013).
 V. Nicosia, G. Bianconi, V. Latora, and M. Barthelemy, Phys. Rev. Lett. 111, 058701 (2013).
 G. Menichetti, D. Remondini, P. Panzarasa, R. Mondragon, and G. Bianconi, PloS One 9 e 97857 (2014).
 R. Mastrandrea, T. Squartini, G. Fagiolo, and D. Garlaschelli, Phys. Rev. E, 90, 062804 (2014).
 M. De Domenico, A. SoléRibalta, E. Cozzo, M. Kivelä, Y. Moreno, M.A. Porter, S. Goméz, and A. Arenas, Phys. Rev. X, 3, 041022 (2013).
 M. Boguná,R. PastorSatorras, and A. Vespignani, The European Physical Journal B 38, 205 (2004).
 G. Bianconi, A.C.C. Coolen, C. PerezVicente, Phys. Rev. E 78, 016114 (2008).
 K. Anand and G. Bianconi, Phys. Rev. E 82, 011116 (2010).
 K. Anand and G. Bianconi, Phys. Rev. E 80, 045102 (2009).
 M. Mezard and A. Montanari, Information, Physics and Computation (Oxford University Press, Oxford, 2009).