Higher-order correlations reveal complex memory in temporal hypergraphs

Many real-world complex systems are characterized by interactions in groups that change in time. Current temporal network approaches, however, are unable to describe group dynamics, as they are based on pairwise interactions only. Here, we use time-varying hypergraphs to describe such systems, and we introduce a framework based on higher-order correlations to characterize their temporal organization. The analysis of human interaction data reveals the existence of coherent and interdependent mesoscopic structures, thus capturing aggregation, fragmentation and nucleation processes in social systems. We introduce a model of temporal hypergraphs with non-Markovian group interactions, which reveals complex memory as a fundamental mechanism underlying the emerging pattern in the data.

Temporal networks, where links connecting pairs of nodes are not continuously active, provide a framework to model how the interactions of a complex system evolve in time [1][2][3].They have revealed key in understanding how the time-varying interaction network of real-world social and biological systems affects the properties of dynamical processes, such as epidemic spreading [4][5][6][7], diffusion [8][9][10][11][12], synchronization [13,14], and others [15][16][17][18].Recent results have highlighted the complex way in which the activity of each link depends on the activities of all other links, showing that memory [19][20][21][22] in temporal networks is inherently a multidimensional concept with a well defined microscopic shape [23].Different approaches have aimed to describe the time evolution of a network as a trajectory in graph space, by naturally extending to the case of graphs notions such as correlations [24] or even dynamical stability [25] traditionally used for scalar or vectorial time-series.
Some early works have already started to explore the temporal dimension of higher-order interactions.For instance, group interactions in real-world social systems have been found to occur in persistent bursts of activity [28], with events of different sizes close in time also spatially correlated in the network [50].Such persistent temporal higher-order interactions have been shown to anticipate the onset of endemic states in epidemic processes [51], and to affect the convergence time of nonlinear consensus dynamics [52].Theoretical frameworks for modeling temporal group activation data [53], and for constructing simplicial complexes based on topological data analysis of multivariate time-series from brain functional activity, financial markets and disease spreading [54] have also been recently developed.However, how to analyze and characterize the temporal organization of real-world complex systems with higher-order interactions is to this day still an open problem.
In this Letter, we bridge this gap by introducing a general framework to study higher-order temporal dependencies in complex systems.We represent a complex system with interactions in groups whose size and composition can change in time as a temporal hypergraph, i.e., a hypergraph with time-varying hyperedges of different orders.We then define a set of measures to extract higher-order temporal correlations, namely to characterize how the dynamics of hyperedges of different orders are correlated.By focusing on a variety of empirical social systems, our framework reveals the existence of long-range correlations at different group sizes, and their hierarchical organization.Furthermore, we uncover the presence of temporal correlations between groups of different sizes, i.e., between hyperedges of different orders, capturing emergent interactions between coherent mesoscopic structures.Finally, to gain intuition about the underlying microscopic mechanisms, we introduce novel models of temporal hypergraphs with higher-order memory, able to explain the observed empirical patterns.
Temporal correlations in hypergraphs.-Systems with higher-order interactions can be represented by hypergraphs [55].A hypergraph is a tuple (V, H), where V is a set of N nodes, and H is a set of M hyperedges.Each hyperedge is a collection of nodes representing an interaction among multiple units.A hyperedge of order 2, or 2-hyperedge, is a set of two nodes representing a two-body interaction, a 3-hyperedge is a set of three nodes representing a group interaction among three units, and so on, up to a order D. A hypergraph can be encoded using a set of adjacency matrices A (2) , A (3) , . . ., A (D) , where the element a ii are set to zero.To describe the temporal evolution of higher-order interactions, we consider temporal hypergraphs [28].A temporal hypergraph is an ordered sequence of (static) hypergraphs (V, H(t)), with t ∈ {1, . . ., T }, where H(t) is the set of hyperedges existing at time t.For each order d, we can define a sequence of adjacency matrices (2), . . ., A (d) (T )} that represents the temporal evolution of the system interactions of order d.
The presence of higher-order interactions makes the analysis of temporal correlations a multi-faceted problem.First, to quantify temporal correlations in interactions of a given order d, we introduce the intra-order correlation matrix with d ∈ {2, . . ., D}.Here, τ is the temporal lag, A denotes the transpose of A, and we have defined the annealed adjacency matrix of order d as µ Note that, for d = 2, Eq. ( 1) recovers the correlation matrix for temporal networks [24].The diagonal terms of C (d) (τ ) capture how hyperedges of order d are temporally autocorrelated, whereas the off-diagonal terms quantify cross-correlations.When the latter are negligible, one can focus on the diagonal terms and define an intra-order correlation function that provides a scalar measure of how hyperedges of order d are autocorrelated at lag τ .Second, one can inquire whether interactions of two different orders d 1 and d 2 display a temporal interdependence, or otherwise evolve independently.To estimate such cross-order temporal correlations, we introduce the cross-order correlation matrix where d 1 , d 2 ∈ {2, . . ., D}.Note that, when d 1 = d 2 = d, we recover the intra-order correlation matrix C (d) .We can then define a scalar cross-order correlation function as All the information about intra-order and cross-order correlations can be encoded in a (D−1)×(D−1) normalized interaction matrix K d1d2 (τ ) = c (d1,d2) (τ )/2 √ σ (d1 σ (d2 , where σ (d) = c (d) (0), whose entry K d1d2 (τ ) describes how interactions of order d 1 at a given time are correlated with those of order d 2 occurring τ time steps later.Notice that matrix K d1d2 (τ ) is not symmetric, as the quantity c (d2,d1) (τ ) measures how order d 1 is correlated with order d 2 at τ time steps before, and is in general different from c (d1,d2) (τ ).The presence of a significant discrepancy between these two quantities captures asymmetries in the temporal dependencies between different orders of interaction.We quantify such an asymmetry in terms of a cross-order gap function A positive value of δ (d1,d2) (τ ) implies that the presence of groups of size d 1 correlates with the presence of groups of size d 2 after a time lag τ , more than the other way around.Temporal dependencies within and across orders emerge in social systems.
-We now characterize temporal correlations in real-world systems with group interactions.We consider different social systems, for which we have high-resolution data about their temporal evolution.We here focus on a dataset describing face-to-face interactions over a period of 32h among the N = 403 participants of a scientific conference [56,57].Three further cases, namely the social interactions occurring in an office [58], in a hospital ward [59] and in a university campus [60] are described in the Supplemental Material (SM).We encode the fine-grained temporal information of the dataset in a temporal hypergraph (V, H(t)), with |V| = 403 and t ∈ {1, . . ., T }.The set H(t) is constructed by assuming that d individuals in contact at a given time t interact together in a group of size d, thus corresponding to a hyperedge of order d at time t.We begin by studying how groups of a given size are temporally correlated.Fig. 1 reports the intra-order correlation functions c (d) (τ ) for orders d ∈ {2, . . ., 4} (circles).Significant long-range temporal autocorrelations emerge for different orders of interaction, as indicated by the slow decays of c (d) (τ ) with τ in a double logarithmic scale, up to a threshold, which typically decreases with d.This indicates that groups of larger sizes generally remain autocorrelated for shorter times.Interestingly, we also observe a saturation effect for interactions in groups of size three, with a series of peaks revealing a weak periodicity at large timescales (see SM for further details).Empirical results are compared with a null model (squares) obtained by reshuffling the sequence defining the temporal hypergraph.
We now investigate whether interactions in groups of a given size d 1 can also be correlated to interactions in groups of size d 2 = d 1 .The cross-order correlation functions c (4,5) (τ ) (cyan circles) and c (5,4) (τ ) (olive circles) for groups of sizes four and five, reported in Fig. 2(a), clearly show the presence of significant cross-order correlations (see SM for an analysis of other group sizes).We observe that both functions display a peak followed by a slow decay, which can not be reproduced by the corresponding null model (squares).For clarity of presentation, we display a binned average of the cross-correlations functions, as well as the corresponding standard deviation.Fig. 2(b) shows the normalized interaction matrix K(τ ) at time lag τ = 600s (see SM for an analysis of different τ ).We observe a banded structure around the main diagonal, meaning that cross-order correlations are higher between groups of similar sizes.This indicates that, in the interactions at a scientific conference analyzed here, groups change gradually, with the loss or the addition of one or few members (see SM for the interaction matrix of different social systems, including the social contacts in a university campus where large groups reveal a more complex correlation pattern).Finally, Fig. 2(c) shows the cross-order gap function δ (4,5) (τ ) (purple circles).Values of δ (4,5) (τ ) > 0 in almost all the range of the time lag τ considered indicate that groups of size four at a given time are correlated to those of size five occurring τ time steps later, more than the other way around.This result, which again cannot be reproduced by the null model (green squares), suggests that the formation of a group of five individuals from a group of four is more probable than the loss of one member in groups of five individuals, indicating a preferred temporal direction in the dynamics of group formation/fragmentation of this social system (see SM for an analysis of δ (d1,d2) (τ ) for other group sizes).
Temporal hypergraph models with higher-order memory.
-To investigate the mechanisms shaping the onset of intra-order and cross-order correlation profiles, we introduce two models to generate temporal hypergraphs with higher-order memory, inspired by DAR processes [21,23,61].The first model, named Discrete Auto Regressive Hypergraph (DARH) model, treats the binary states of each hyperedge h α ∈ {0, 1} (absent/present) as independent stochastic processes.Each hyperedge updates its state either drawing a state from its past, or randomly sampling a new state.With probability q (d) , where d ∈ {2, . . ., D}, a hyperedge of order d samples its state uniformly at random from its m (d) s previous states.With probability 1 − q (d) , the hyperedge state is drawn randomly following a Bernoulli process with probability y (d) .Formally, the dynamics of a hyperedge of order d, h α , is given by Our first model displays intra-order, but no cross-order temporal correlations (see SM for a characterization of the DARH model).We then introduce a second model, the cross-memory DARH (cDARH) model, a variation of the DARH model where a hyperedge of order d updates its state by drawing not only from its past but also from that of a hyperedge of a different order.We will refer to this latter property as cross-order memory.When copying from memory, each hyperedge draws from the past of other hyperedges with probability p (d) , and from its own past with probability 1 − p (d) .We assume hyperedges can copy from the memory of overlapping hyperedges only.This choice is motivated by previous studies on higher-order interactions in social networks that pointed out a tendency of groups to progressively add or remove members, one step at a time [28].For illustration, let us consider groups of size two and three: a 2-hyperedge {i, j} selects one of the (N −2) possible 3-hyperedges containing nodes i and j and draws from its previous m (3,2) c states.Similarly, a 3-hyperedge {i, j, k} selects one of the three 2-hyperedges that can be formed from it, i.e., {i, j}, {j, k}, and {i, k}, and copies a state from its previous m (2,3) c steps.Such a mechanism can be straightforwardly extended to hyperedges of other orders.Formally, we can write the dynamics of a hyperedge of order d as As for the DARH model, we have that Q t ∼ Bernoulli(q (d) ) and Y t ∼ Bernoulli(y (d) ), while the new random variable ε t (α) selects the hyperedge from which the past is copied.Consequently, the variable µ does not only depend on the order of h α but also on the value of ε t (α).For instance, in the case of 2hyperedges, µ ∼ Uniform(1, m ) if ε t (α) = α, when copying from 3-hyperedges.
We generate with the cDARH model temporal hypergraphs with N = 10 nodes, maximum hyperedge order D = 3 and a temporal range of T = 3 • 10 4 time steps.Fig. 3 reports the results obtained by setting p (2) = 0 and p (3) = 0.6, meaning that hyperedges of order three can copy from the past of hyperedges of order two, while hyperedges of order two evolve independently.We also set Intra-order correlations c (d) , with d ∈ {2, 3}, for hyperedges of order two (green circles and squares) and three (orange circles and squares).The dashed vertical lines correspond to the value of the intra-order memory of hyperedges of order two (green) and three (orange), respectively.(b) Cross-order gap function δ (2,3) between hyperedges of order two and three.We compare intra-order and cross-order correlations for the hypergraph generated using the cDARH model with a randomized null model.

m
(2,3) s = 60, for the cross-order memory of 3-hyperedges, m s = 40 and m (3) s = 10 for the intra-order memories of hyperedges of order two and three, respectively.Fig. 3(a) reveals the presence in hyperedges of order two (green circles) and three (orange circles) of significant intra-order temporal correlations.In particular, we notice that the functions c (2) and c (3) remain constant for τ ≤ m (d) , while they decay exponentially after that value, with the same rate of decay (see SM for a deeper analysis).Such a minimal model reveals that memory can be the driving mechanism for the emergence of intra-order temporal correlations, with different orders possessing different degrees of memory, explaining the hierarchical structure of correlation observed in the data.The profiles of c (d) (τ ) do not match exactly those of empirical data, which are characterized by a power-law decay.This is not surprising, as in the cDARH model a single set of values m (2) s and m (3) s is considered, while real-world social interactions can be shaped by different scales of memory [23].Fig. 3(b) shows that δ (2,3) (τ ) > 0 for different values of τ (purple circles), meaning that hyperedges of order two are correlated to hyperedges of order three occurring later in time more than the other way around.We observe a striking similarity between this trend and that observed in Fig. 2(c) for the empirical data.This result indicates that cross-order memory is a fundamental factor for the emergence of cross-order gaps in real-world social interactions.The two peaks for δ (4,5) (τ ) observed in the empirical system suggest again a more complex dependence on memory, possibly due to multiple temporal scales.
Discussion.-In this Letter, we have introduced a framework to characterize different dimensions of memory in networked systems with higher-order interactions.We have shown that real-world social systems display longrange temporal correlations at different group sizes, or-ganized in a hierarchy across multiple scales.Moreover, we found that group interactions are characterized by cross-order correlations, suggesting the emergence of nontrivial interactions between mesoscopic memories, overall leading to a complex memory structure.In the context of social systems, such cross-order interactions can be interpreted in terms of the schisming phenomenon [62][63][64][65], where group sizes in human interactions, e.g.conversations, fluctuate, nucleate, and display complex dynamics.The formalism presented here can be naturally extended to other higher-order complex systems traditionally modeled in terms of networks of interactions, such as the human brain and biological ecosystems.In conclusion, our work sheds light on the unexpectedly complex and multifaceted nature of memory that emerges in real-world interacting systems and we hope can open new avenues also in the higher-order dynamics of emergent coherent structures in a variety of physical systems, from multifragmentation in nuclear physics to vortex interaction in the atmosphere or other fluid dynamical systems.
Supplemental Material of the manuscript Higher-order correlations reveal complex memory in temporal hypergraphs

INTRA-ORDER CORRELATIONS IN SOCIAL SYSTEMS WITH HIGHER-ORDER INTERACTIONS
In the main text, we focused on a dataset describing face-to-face interactions in a scientific conference.In this section, we present the analysis of the intra-order temporal correlations for three other social systems.Two of the datasets come from the Sociopatterns project, and they describe face-to-face interactions (i) in an office, among N = 92 workers, over 11 days [58]; and (ii) in a hospital ward, among N = 75 patients, doctors, nurses, and administrative staff, over 72h [59].The third dataset contains instead proximity data of N = 692 students at the campus of the Technical University of Denmark (DTU) in Copenhagen, for a period of one month [60].In this section, we analyze temporal correlations within groups of the same size.Fig. S1 depicts the values of the intra-order correlation functions for several group sizes d, for the three empirical systems examined, namely face-toface interactions in the workplace (a), the hospital ward (b), and in the university campus (c).For each system, we consider only those group sizes for which we have a sufficient statistics.In general, for each social system, a hierarchy of long-range, slowly-decaying temporal autocorrelations emerges across different group sizes.Also, we note that intra-order correlations are lost abruptly after a typical time threshold that decreases with the group size d, i.e., larger groups remain temporally correlated for shorter time.
We observe periodicity patterns in the interactions occurring in the workplace and in the university campus.For this latter in particular, we note a series of peaks at large timescales that show the same hierarchical structure observed at small timescales.See Fig. S10 for a minimal model of a temporal hypergraph with periodicity.

CROSS-ORDER CORRELATIONS SHOW NONTRIVIAL TEMPORAL ORGANIZATION OF HIGHER-ORDER INTERACTIONS
In this section, we deepen our analysis on temporal correlations between groups of different sizes.In the main text, we showed that statistically significant cross-order correlations emerge.In particular, we focused on the cross-order correlations between groups of sizes four and five in the face-to-face interactions in the scientific conference.We here show that other pairs of group sizes in the same system are temporally correlated.
Fig. S2 displays the cross-order correlation functions c (2,3) (τ ) (cyan) and c (3,2) (τ ) (olive), describing cross-order correlations between groups of sizes two and three.We also report the cross-order correlation functions for a randomized null model where correlations have been removed by time-reshuffling (lighter squares) Both functions exhibit an exponential decay that does not match the trends of the cross-order correlation functions c (4,5) (τ ) and c (5,4) (τ ) shown in Fig. 2(a) of the main text, suggesting that different pairs of group sizes can show different kinds of interdependencies.
Fig. S3 displays the normalized interaction matrix K(τ ) at different time lags, namely τ = 60s (a), τ = 300s (b), and τ = 1500s (c).We observe that cross-order temporal correlations can emerge between different group sizes.In particular, as discussed in the main text, the matrix K(τ ) shows a banded structure around the main diagonal, meaning that correlations between groups of similar sizes is higher.The banded structure is lost at large time lags, as correlations vanish at larger timescales.The banded structure of the matrix K(τ ) suggests that, in the scientific conference dataset, groups tend to change gradually, with the loss or the addition of one member.We now want to investigate whether groups of larger size display the same property.To answer this question, we focus on the interactions among students occurring in the DTU campus in Copenhagen, in which larger groups have a sufficient statistics.Fig. S4 displays the normalized interaction matrix K(τ ) at time lag τ = 300s.Once again, the matrix K(τ ) reveals a characteristic banded structure.However, compared to the scientific conference case, we note that groups of larger size show a different pattern.Indeed, we can observe significant cross-order correlations between groups of dissimilar sizes, e.g., groups of size nine are temporally correlated with groups of size five.This result suggests that larger groups evolve differently compared to smaller ones.In particular, while small-size groups change gradually, losing or gaining one member at a time, large-size groups show more complex dynamics, as they can also split in smaller groups or emerge from them.For completeness, in Fig. S5 we report the normalized interaction matrices evaluated at time lag τ = 300s relative to the interactions in groups of size two and three occurring in the office (a) and in the hospital ward (b).In the main text, we showed that the cross-order correlation functions relative to groups of sizes d 1 and d 2 can differ significantly, as groups of a given size can anticipate those of another one, while the vice versa might not be true, or the magnitude of the two effects might be distinct.This discrepancy between the two functions is quantified in terms of the cross-order gap function δ (d1,d2) (τ ).To complement the analysis provided in the main text, here we consider the interactions in the scientific conference, studying the possible gaps in the cross-order correlation functions for other pairs of group sizes.Fig. S6 shows the cross-order gap function δ (2,3) (τ ) for the real social network (purple circles).Note that the figure axis scales are the same as those of Fig. 2(c) of the main text, for a proper comparison.We observe that δ (2,3) (τ ) does not significantly differ from the null model where temporal snapshots are reshuffled.This indicates that neither the groups of size two nor those of size three tend to anticipate the other, suggesting that groups of two people gain a member with the same probability as groups of three people lose one.d) for hyperedges of orders two and three, respectively.The dashed vertical lines correspond to the value of the intra-order memory of hyperedges of order two (green) and three (orange).The solid black line displays the best fit of the exponential decay for both c (2) and c (3) .
to better appreciate the rate comparison, here we consider a lin-log scale, whereas in Fig. 3(a) of the main text we used a log-log scale.

CORRELATION ANALYSIS OF A PERIODIC TEMPORAL HYPERGRAPH
Social systems show a variety of temporal patterns, from periodicity to decaying autocorrelations.In the main text, we concentrated on this latter feature, discussing a minimal model of temporal higher-order networks with memory.In this section, we focus instead on temporal periodic patterns in systems characterized by higher-order interactions.To do so, we build a temporal hypergraph of N nodes where each order of interaction has different periodicity.For simplicity and with no loss of generality, we limit our analysis to groups of size two and three.First, we construct a sequence of length T 2 of random 2-uniform hypergraphs, namely hypergraphs where all hyperedges have order two.Each hyperedge has a probability p 2 to be created.Each hypergraph of the sequence is nothing more than a Erdös-Rényi random network, i.e., ER(N, p 2 ).We then construct another sequence, having length T 3 , and composed by random 3-uniform hypergraphs, i.e., all hyperedges connect exactly three nodes.A hyperedge has a probability p 3 to be created.For each order of interaction, we concatenate several of the corresponding sequences one after another, thus building two temporal uniform hypergraphs.We repeat the sequences so that both temporal hypergraphs have T time-steps each.Finally, for each time-step we join the two corresponding static uniform hypergraphs.By doing this, we now have a single temporal hypergraph of length T where the 2-hyperedges are repeated every T 2 time-steps, while the 3-hyperedges have instead a period T 3 .Both orders of interactions clearly show a periodic behavior.Indeed, we observe that the intra-order correlation functions peak at τ = kT 2 and τ = kT 3 , with k ∈ N, for two-body and three-body interactions, respectively, while c (2) (τ ) and c (3) (τ ) remain close to zero for other values of τ .The periodicity of the temporal hypergraph is also reflected in the cross-order correlation functions c (2,3) (τ ) (cyan) and c (3,2) (τ ) (olive).In particular, we observe that they peak at different harmonics and that in general c (2,3) (τ ) = c (3,2) (τ ), as the periodic oscillations of the two orders of interaction are not in phase.
ij of the matrix A (d) counts the number of d-hyperedges both nodes i and j belong to, while elements a (d)

FIG. 1 .
FIG. 1. Intra-order correlations in human face-to-face interactions.Circles show the value of c (d) (τ ) for interactions in groups of different sizes, d ∈ {2, . . ., 4}.Squares refer to a randomized null model where temporal correlations have been removed by reshuffling time-steps.

FIG. 3 .
FIG. 3. Temporal correlations in the cDARH model.(a)Intra-order correlations c(d) , with d ∈ {2, 3}, for hyperedges of order two (green circles and squares) and three (orange circles and squares).The dashed vertical lines correspond to the value of the intra-order memory of hyperedges of order two (green) and three (orange), respectively.(b) Cross-order gap function δ(2,3) between hyperedges of order two and three.We compare intra-order and cross-order correlations for the hypergraph generated using the cDARH model with a randomized null model.
FIG. S1.Intra-order correlation functions in human face-to-face interactions in the office (a), in the hospital ward (b), and in the university campus (c), for different group sizes.
FIG.S2.Cross-order correlation functions for groups of sizes two and three, i.e., c(2,3) (τ ) (cyan) and c(3,2) (τ ) (olive), in the face-to-face interactions in the scientific conference.The empirical system (dark circles) is compared with a randomized null model with reshuffled time-steps (light squares).
FIG. S5.Normalized interaction matrix K(τ ) at time lag τ = 300s for the face-to-face interactions occurring in the office (a) and in the ospital ward (b).
FIG.S6.Cross-order gap function δ(2,3) (τ ), for interactions in groups of size two and three in the scientific conference (purple circles), compared with the null model (green squares).
FIG.S9.Intra-order correlations in the cDARH model.Green and orange circles represent the value of c (d) for hyperedges of orders two and three, respectively.The dashed vertical lines correspond to the value of the intra-order memory of hyperedges of order two (green) and three (orange).The solid black line displays the best fit of the exponential decay for both c(2) and c(3) .
Fig. S10 shows the intra-order correlation functions (a) and the cross-order correlation functions (b) for a periodic temporal hypergraph of N = 10 nodes.The hypergraph is composed by T = 1.95 • 10 5 snapshots.Each 2-hyperedge has a probability p 2 = 0.1 of existing and it repeats with period T 2 = 15.The 3-hyperedges have instead a probability p 3 = 0.03 of being active and they have period T 3 = 8.