The customary perspective to reason about epidemic mitigation in temporal networks hinges on the identification of nodes with specific features or network roles. The ensuing individual-based control strategies, however, are difficult to carry out in practice and ignore important correlations between topological and temporal patterns. Here we adopt a mesoscopic perspective and present a principled framework to identify collective features at multiple scales and rank their importance for epidemic spread. We use tensor decomposition techniques to build an additive representation of a temporal network in terms of mesostructures, such as cohesive clusters and temporally-localized mixing patterns. This representation allows to determine the impact of individual mesostructures on epidemic spread and to assess the effect of targeted interventions that remove chosen structures. We illustrate this approach using high-resolution social network data on face-to-face interactions in a school and show that our method affords the design of effective mesoscale interventions. Keywords: complex networks, temporal networks, tensor decomposition, epidemic spread, high-resolution social networks, human contact networks, targeted interventions Many natural and artificial systems are adequately represented in terms of interaction networks between their components 1 . The network representations has led to key insights on the relation between the structure of these systems and the dynamics of diverse processes that take place over the network, such as information diffusion, epidemic spread, and much more 2;3 . In recent years, the availability of time-resolved network data has pushed network science beyond the static graph representation and has prompted new research on understanding and modeling time-varying graphs, commonly referred to as temporal networks 4 . This has spurred an intense activity on studying how dynamical processes are affected by the temporal evolution of the network on which they unfold, focusing on heterogeneous distributions of interevent times ("burstiness"), heterogeneous activity distributions, causality constraints, etc. [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] . Highresolution social network data [17] [18] [19] [20] , in particular, have opened the way to richer individual-based models of epidemic spread 18;21-23 and have raised new questions on designing control strategies for epidemic processes in temporal networks 22;24-26 . However, micro-interventions that operate on individual nodes, such as targeted vaccination, are difficult to implement, both for lack of high-resolution data in general cases and because action at the level of individuals is subjected to many informational and operational constraints. Therefore, it is natural to consider interventions that target important collective patterns and structural units 32 to represent a temporal network as an additive superposition of mesostructures, each of which can encode complex structural and temporal correlations. To assess the importance of individual mesostructures, we take advantage of the additivity of the tensor decomposition, which allows us to dissect the original network and reassemble altered temporal networks where chosen mesostructures are selectively removed. By comparing the dynamics of a given process on the original network and on the altered ones we can then determine the specific role of individual mesostructures and we can rank the effectiveness of interventions that target those structures for removal. We illustrate our approach using high-resolution social network data on face-to-face interactions in a school setting and investigating the importance of individual mesostructures for simple epidemic processes. Schools are actually an interesting context for models of epidemic spread, as they are thought to play an important role in the community spread of infectious diseases [33] [34] [35] [36] [37] . Macro-scale interventions such as school closure are considered a viable strategy for epidemic mitigation 38 but come with steep socio-economic costs 39;40 , thus calling for more targeted approaches. We show that our approach identifies effective strategies for epidemic mitigation that involve the targeted removal of individual mesostructures. Remarkably, we find that the most important mesostructures for epidemic spread correspond neither to the strongest structures nor to cohesive structures that can be easily detected by community detection methods. Mesoscale structure of temporal networks We first summarize a general framework that we introduced in ref. 31 to expose the mesoscale structure of a temporal network. The approach combines a tensor representation of a temporal network 41 and a dimensionality reduction technique based on tensor decomposition 32 . We consider an undirected and unweighted temporal network, hence the state of the network at time t can be represented by a binary-valued adjacency matrix M(t) ∈ R N ×N , where the matrix entry M ij (t) indicates the status of link i-j at time t, and N is the number of nodes of the network. The matrices describing the state of the network at different times can be combined into a 3-mode tensor, T ∈ R N ×N ×L , where the first two dimensions are the customary node indices of the adjacency matrix and the third dimension is a temporal index. L is the number of network snapshots, each for a different point in time. To detect structures, i.e., correlated activity/connectivity patterns of the temporal network, we set up a latent factor analysis based on non-negative tensor factorization 31 : The central idea is to approximate the tensor T with a sumT of rank-1 non-negative tensors (see Methods): with R < min {N, S}. This decomposition allows to describe the original temporal network as a purely additive superposition of a chosen number of activity/connectivity patterns, encoded by the individual tensor components. Each component S r is regarded as a mesoscale structure of the temporal network and can encode cohesive clus-2 ters (i.e., communities), temporally-localized mixing patterns, and more. Such a parts-based representation is a general and powerful feature of tensor decompositions with non-negativity constraints and greatly enhances the interpretability of individual components. A similar approach based on nonnegative matrix factorization has been proposed to detect overlapping communities in static networks 43 . The additive decomposition of Eq. 1 also suggests a natural way to investigate the role played by individual mesoscale structures, both on the overall network architecture and on the dynamics of processes occurring over the network, e.g., epidemic spread. In particular, we can assess the importance of individual mesostructures for dynamical processes unfolding over the temporal network by selectively removing those structures: We can start from the decomposition of Eq. 1, remove a chosen term (i.e., a chosen structure), and sum all the other terms, obtaining a modified temporal network where the chosen structure has been selectively erased. We can then simulate a given dynamical process over the original temporal network and over the modified temporal network, and compare the dynamics or the outcome of the process, hence learning about the impact of the removed structure. As a case study we consider an empirical temporal network of human face-to-face interactions in a primary school 44 , measured by using wearable sensors 19 . The dataset describes the interactions (or "contacts") between 231 children aged 6 to 12 and 10 teachers, organized in 10 classes. It comprises the contacts that occurred at the school premises during two consecutive days in October 2009, from 8:30am to 5:15pm of each day. We chose this dataset because it is publicly available 23;44 and provides not only a high-resolution temporal network, but also ground truth information on communities, namely, class attributes for all nodes. Moreover, the school schedule is known and can be used to understand the activity patterns that mix multiple classes, as done in Ref. 31. In the following this will allow us to relate the effect of removing an individual component to specific knowledge about the social behavior associated with that component. The school data are naturally represented as a temporal network, where nodes are individuals and links are interactions between pairs of nodes. The networks snapshots of our dataset are recorded every 20 seconds, however the school activity is scheduled on a coarser temporal scale. Hence, we aggregate the data over consecutive time intervals of 15 minutes to obtain unweighted network snapshots: for every 15-minute interval, a link is drawn between two nodes if the corresponding individuals had at least one contact during that interval. This yields a temporal network represented by binaryvalued valued tensor T with N = 241 and L = 70. Following ref. 31, we carry out the non-negative tensor factorization of Eq. 1 and approximate the tensor T with the sumT of R = 14 components (see Methods and Supplementary Information): the results of the decomposition are reported in Fig. 1 . Components S 1 through S 9 and component S 11 correspond to the 10 school classes of the school. They can be validated using the know class attributes for students and teachers, and they can also be detected to some extent using community detection algorithms on the temporally-aggregated network 31 . The components S 10 , S 12 , S 13 and S 14 , however, mix individuals from different classes and -contrary to the above class components -they are only active during lunch breaks. Such components cannot be found by community detection of the temporally-aggregated networks as they correspond to temporally-localized mixing patterns of the temporal network, due to mixing of students at lunch time. Tensor decomposition, however, can naturally detect these patterns as well as the customary cohesive network communities. We remark that whereas the original temporal network T is unweighted, the approximating networkT is in general weighted (see Methods). For any given network snapshot, corresponding to a 15-minute time interval, we will interpret the weight of the edge between two nodes as the cumulated time that those nodes spent in contact over that 15-minute interval. In Table 1 we report the fraction of total tensor weight corresponding to each of the R = 14 components, S r / T . Left-hand color-coded matrix: Node weights for each component. Each row corresponds to a vector a r of the decomposition, normalized to unity in Euclidean norm for clarity of visualization. Nodes are ordered to expose the block structure of the matrix due to school classes. Right-hand plots: Component activity as a function of time over the 2-day span of the dataset. Each row correspond to a vector c r , normalized to unity in Euclidean norm. Components S 1 to S 9 and component S 11 are school classes, while components S 10 , S 12 , S 13 and S 14 mix individuals from different classes and are active during lunch breaks. Mesoscale targeted interventions for epidemic spreading The decomposition of a temporal network into a superposition of components, interpreted as mesostructures, makes it possible to pinpoint and assess the contribution of individual mesoscopic features to the overall structural and functional properties of the temporal network. Indeed, the additive decomposition of Eq. 1 allows to selectively remove a chosen component s by excluding the corresponding term, yielding an altered temporal net- which can be compared to the original network to elucidate the specific role of component s. When component s can be interpreted in terms of a specific behavioral pattern (e.g., the lunch breaks of Fig. 1 ) its removal can be regarded as the effect of an intervention strategy that selectively targets that behavior (e.g., removing lunch breaks from the school schedule). As a first step in exploring mesoscale interventions, here we focus on the case study of the previous section and investigate the effect of removing individual mesostructures on the dynamics of simple epidemic processes unfolding over the school temporal network. We study how the timing and size of the epidemic are affected by the removal of individual network components, that is, by different targeted interventions aimed at removing them. We start with the simple case of a susceptibleinfected (SI) process. Each network node can be in either of two states: susceptible (S) or infected (I). A susceptible node in contact with an infected one becomes infected with a fixed probability λ per unit time. Once infected, a node stays indefinitely in that state. The system is initialized with all 4 r 1 2 3 4 5 6 7 weight fraction 11.3% 8.6% 8.8% 7.1% 8.3% 6.9% 5.7% r 8 9 10 11 12 13 14 weight fraction 6.8% 7.5% 7.5% 4.8% 6.8% 5.7% 6% Table 1 : Component weight. Fraction of total tensor weight corresponding to each term S r of the decompositionT . nodes in the S state, except for a single infected node (seed node). The timing and duration of contacts between nodes is described by the temporal networkT , which is a sequence of weighted networks: for each network snapshot, the weight of an edge between two nodes gives the total duration of the contacts between those nodes during the corresponding time interval. Although the SI process does not describe any realistic infectious disease, it represents a paradigmatic dynamical process and it is frequently used to investigate the structural properties of temporal networks. To investigate the effect on epidemic spread of removing a chosen mesostructure S r , we simulate the epidemic process on the altered networkT r , and compare its dynamics to that observed on the full temporal networkT . We quantify the effect of removing component r by measuring the epidemic delay ratio where t j is the half-infection time onT when the seed is node j, and t r j is the same quantity for the altered temporal networkT r . The half-infection time is measured from the first time the seed node infects another node, and the average · is computed over all possible seed nodes j and different starting times for the SI process (see Supplementary Information for details). We also study the impact of mesoscale features on the dynamics of a more realistic epidemic process, the Susceptible-Infected-Recovered (SIR) model. In this model, S → I transitions occur as in the SI case when a susceptible node is in contact with an infected one, with probability λ per unit time. Infected nodes recover with a constant rate µ, and Recovered (R) nodes no longer take part in the epidemic process. To assess the effect of removing a given mesostructure r we compare the size of the epidemics at the end of the process, measured as the final number of nodes that are in states I or R. We compute the ratio ρ between the average size Ω r of the epidemic on the altered networkT r and the average size Ω of the epidemic on the full networkT : The average is computed over all possible seed nodes and over different seeding times near the beginning of the data set (see SI for details). Note that, due to the finite temporal span of our network data, the epidemic process may not be over by the end of the simulation (i.e., some nodes may still be in the I state). In this case Eq. 4 only yields an estimate of how much the epidemic has been mitigated during the time interval covered by the data. We start from the decomposition of the school temporal network described in Fig. 1 and study the effect of removing individual mesostructures on the dynamics of an SI process. We consider several values of the transmission rate λ: here we present results for λ = 0.7 and in the Supplementary Information we show that qualitatively similar results hold for other values of the parameter. Figure 2 shows the delay ratio τ r of Eq. 3 obtained by selectively removing each of the components S r , with r = 1, . . . , 14: that is, comparing the dynamics of the SI process over each of the altered networks T r with the dynamics observed for the full temporal networkT . Since the altered networks have by definition a smaller total weight T r than the full network, we expect epidemic spread to be mitigated for the altered networks. To assess whether a given component r plays a structural role for epidemic spread that goes beyond its weight, it is thus important to compare the delay ratio τ r observed 5 on removing that component with a suitable null model. Such a null model is obtained by measuring the delay ratio for several stochastic realizations of altered temporal networks obtained fromT by removing, at random, a fraction of weight equal to that of component r (as given by Table 1 ). The distributions of τ r for the null models are shown by the boxplots of Fig. 2 . Figure 2 : Epidemic delay ratio as a function of removed component. The epidemic delay ratio τ r of Eq. 3 is plotted (solid circles) as a function of the removed network component and compared to the null model described in the main text (boxplots). In each boxplot, the central band is the median, the bottom and top of the box are the first and third quartiles. Whiskers extent to the most extreme data point within 1.5 IQR of the inner quartiles. Points beyond whiskers are outliers. We observe that the delay ratio obtained by selectively removing a component (solid circles) is almost always larger than the one obtained by removing links at random (boxplots): Such targeted interventions have a stronger effect than random removals of links. For the first 9 structures, the delay ratio is positively correlated with the component weight of Table 1 , as we could naively expect. Strikingly, the removal of components 12 or 14 slows down the epidemic spread considerably. These two mesostructures, despite carrying a smaller total weight than most other components (see Table 1 ), thus appear to play a more important role for SI spreading. As shown in Fig. 1 both components 12 and 14 mix nodes from different classes and have activity concentrated during the lunch break of the first day. The other components correspond either to individual classes (components 1 through 9 and 11) or to components, such as 10 and 13, with little activity on the first day. Overall, the SI process is slowed down the most by removing comparatively weak structures (in terms of weight) that mix classes and occur early in the data. We remark that the most effective intervention for mitigating epidemic spread thus involves the removal of seemingly minor mesostructures that correspond to complex correlations between link activity and network structure, and cannot be uncovered by standard community detection approaches. We now turn to the case of an SIR process and investigate the effect of mesostructure removal on the dynamics of an SIR epidemic over the temporal network. We compute the epidemic size ratio ρ λ,µ (r) of Eq. 4 for a wide range of values of the parameters λ and µ, and for all components r. The results are shown in Fig 3 for selected components and in the Supplementary Information for all the other components. We notice that for large λ values and small µ values (bottom right of the heat maps) the epidemic size ratio ρ λ,µ (r) 1: in this region the epidemics spreads and finishes fast, hence mesostructure removal cannot mitigate its size. The parameter region where mesostructure removal can have a strong effect, therefore, is the arc-shaped region visible, e.g., in the heat map for r = 1 (lighter blue/green). The removal of class components r = 1 and r = 2 mildly mitigates the epidemic for a broad range of parameter within that region (almost no effect is obtained for the other class components, as shown in Fig. S4 of the Supplementary Information). Removing the class component 1 has a comparatively higher effect because of the higher overall weight of that component (Table 1) . On taking a closer look at Fig. 3 we notice that, depending on the removed component, significant epidemic mitigation can be achieved in two main parameter regions within the arc. The first region falls within the black contour (high values of λ and µ) and corresponds to parameter values yielding epidemics that finish by the end of the second day. In this region the epidemic can be strongly mitigated (ρ λ,µ 0.5) by removing mesostructures that involve mixing of classes on the first school day, namely r = 12 and r = 14. This is con-6 Figure 3 : Epidemic size ratio as a function of SIR model parameters. Each heat map corresponds to a targeted intervention that selectively removes one component. For each removed component r the heat map shows the epidemic size ratio ρ λ,µ (r) as a function of the SIR parameters λ and µ. ρ λ,µ (r) = 1 indicates that the intervention does not affect epidemic size. The white area is the region where the epidemic dies out, i.e., it fails to affect more than 1% of the network nodes. The region inside the black contour line corresponds to parameter values such that the SIR epidemic finishes within the finite span of the school dataset (2 days), both for the full and the altered temporal networks. That is, for those parameter values the epidemic size ratio is not affected by the finite temporal span of the dataset. sistent with the results of Fig. 2 for the SI case. Conversely, removing the class-mixing components r = 10 and r = 13, which are mostly active on the second day (Fig. 1) has a negligible effect in the same parameter region. For a sample parameter choice within this region, in Fig. 4 we illustrate the temporal evolution of the epidemics and the effect of removing individual components: the epidemic on the unmodified network peaks on the first day and, consistently with the above remarks, removing r = 12 or r = 14 strongly mitigates the fraction of affected nodes, whereas the removal of r = 10 and r = 13 has a negligible effect. We remark that removing the early class-mixing mesostructures also mitigates the SIR epidemic for a wide range of parameter values outside the black contour, as visible in Fig. 3 , and that removing these mesostructures is more effective, in general, than removing class components, even though the latter account for a larger network weight. The second parameter region of interest lies in the lower part of the arc (low values of λ and µ) and is visible in Fig. 3 for r = 10 and r = 13. This region corresponds to parameter values for which the SIR process is slower and the epidemic does not finish over the 2-day temporal span of the data. In this case, removing either of the class-mixing mesostructures r = 10, 12, 13, 14 as a similar and limited effect on the epidemic size reached at the end of the second day. We have put forward a methodology to systematically assess the relevance of mesostructures for dynamical processes such as epidemic spread in temporal networks. The proposed method only uses time-resolved connectivity information without metadata and allows to uncover complex patterns of correlated link activity whose removal has a strong impact on simulated epidemic spreading for Fraction of infected nodes as a function of time for a sample parameter choice. Temporal evolution of the average fraction of infected nodes (top) and of the cumulative fraction of infected nodes (bottom) as a function of time for a sample parameter choice of λ = 0.11 and µ = 0.16. The curve marked as none corresponds to the spread on the unmodified temporal networkT . The other curves show the effect of selectively removing components r = 1, 10, 11, 12, 14. For the chosen parameter values, epidemic spread occurs on the first day, hence removing the early class-mixing components r = 12 or r = 14 has a strong mitigation effect. a broad range of model parameters. The removal of an individual mesostructure can be regarded as a targeted intervention at the mesoscale, such as the removal of a specific group of nodes or the suppression of a specific activity pattern (e.g., a lunch break). These targeted interventions have a potentially high relevance, as they do not involve the system as a whole, nor they need actions at microscopic level of individual nodes. Hence, they might strike a balance between cost and effectiveness that might make them actionable in a variety of applications. As an illustration of the method, we have considered a case study where simple but paradigmatic spreading processes occur on top of an empirical temporal network describing face-to-face interactions in a school. In this case, the most important mesostructures for epidemic spread are not those involving most of the link activity, but instead consist of weaker, temporally-localized mixing patterns corresponding to scheduled social activities. It is important to remark that these mesostructures cannot by found by means of community detection techniques on static networks, and that our methods reveals them in an unsupervised fashion. We also notice that the actionability of an intervention that suppresses a given mesostructure hinges on the recognizability of that mesostructure in terms of relevant contextual information, such as a known or understandable grouping of nodes, temporal activity localized at specific times of a known schedule, or a specific relation to other metadata for the context at hand. For example, the selective removal of component S 12 of Fig. 1 only becomes an actionable targeted intervention when its temporal activity profile is identified as a mixing pattern associated with the lunch break on the first school day. Our work is a first step towards the systematic design of mesoscale interventions in temporal networks. Natural future steps include the use of our methodology on different temporal networks, both empirical and synthetic, the investigation of the impact of incomplete or noisy data, and the study of other dynamical processes. 8 A 3-mode tensor T , with entries T ijk , can be approximated by a sumT of R rank-1 tensors 32;45;46 , where each component tensor S r is the outer products of three vectors, S r = a r • b r • c r . The Frobenius norm of the difference between T and T , = T −T , quantifies the error in approximating the original tensor with the sumT . The number of components R is interpreted as the number of sought mesostructures in the decomposition of the temporal network. Its value needs to be set by balancing the conflicting goals of minimizing , i.e., recovering as much as possible of the original tensor, and avoiding overfitting 31;32;47 . The sets of vectors a {1,2,...,R} , b {1,2,...,R} , c {1,2,...,R} can be combined into matrices A ∈ R N ×R , B ∈ R N ×R and C ∈ R S×R . These matrices all have R columns, one for each component of the decomposition, i.e., one for each sought mesoscale structure. For an undirected temporal network (such as the one of our case study) the adjacency matrix represented on each tensor slice is symmetric and the decomposition yields A B. In practice it is possible to directly seek a decomposition with A ≡ B, i.e., with component tensors of the form S r = a r •a r •c r . Hence in the main text we only refer to the matrices A and C: the matrix elements a ir of A associate each component r = 1, . . . , R with individual nodes of the network, while the matrix elements c kr of C encode the temporal activity pattern of each component, with time indexed by k. The decomposition of Eq. 5 can be written as a ir a jr c kr . The matrix elements a ir indicate notion of membership of node i to component r, so this representation is suitable to describe overlapping mesostructures and complex correlations between connectivity patterns and activity patterns over time. Notice that even when the original temporal network is unweighted, that is, the tensor T is binary-valued, the component tensors S r correspond in general to weighted networks, and the approximating tensor T is in general real-valued. To compute the decomposition of Eq. 5 we need to solve an optimization problem that minimizes the residual . A standard way to do so is to convert the 3-mode problem yielded by Eq. 5 into three coupled 2-mode sub-problems: This is done by unfolding the original tensor T along its three modes, a technique also known as matricization 45 . Using this representation, the original tensor approximation problem is cast into 3 coupled matrix approximation problems that involve the factor matrices A, B, and C. Since the corresponding optimization problems are convex with respect to either of the matrices -but not with respect to them allit is possible to solve the optimization through a method known as Alternative Non-negative Least Squares 48 . The minimization of is usually carried out with non-negativity and/or sparsity constraints on the factor matrices 32;45 , as this is known to yield decompositions that can be interpreted as parts-based representations of the original data 42 . Here we carry out the decomposition with non-negativity constraints, using Alternating Nonnegative Least Squares and a block-coordinate descent method to improve convergence speed 49;50 . Our implementation builds on the MATLAB code of ref. 51. The number of components R is set by using the core consistency diagnostic 31;47 optimized over 5 stochastic realizations of the decomposition that differ for the random initial conditions of the factor matrices. The components S r are ordered so that the Euclidean norm of S r decreases for increasing r. That is, the first column of matrices A, B and C correspond to the component with highest Euclidean norm, and successive columns encode weaker and weaker components.
Based on time-evolving networks, the work in REF identifies and ranks collective features for epidemic spread.
Revealing latent factors of temporal networks for mesoscale intervention in epidemic spread
ArXiv
ArXiv
Mathematics Physics Computer Science