Typesetting math: 100%
Skip to main content Skip to article

### Reading Assistant

#### Questions you could ask:

#### Actions you could take:

Ask about this article
Articles left: 5
Renews on 28 February 2026

- View PDF

## Physics Reports

Volume 1018, 23 May 2023, Pages 1-64

# The structure and dynamics of networks with higher order interactions

Author links open overlay panel S. Boccaletti a b c d, P. De Lellis e,
,
, R. Criado a, S. Jalan d, M. Romance a
https://doi.org/10.1016/j.physrep.2023.04.002
Get rights and content
Full text access

## Highlights

- •
  Networks with higher-order interactions.
- •
  Processes on networked systems beyond pairwise interactions.
- •
  Emergent collective dynamics in hypernetworks.
- •
  Control in hypernetworks.

## Abstract

All beauty, richness and harmony in the emergent dynamics of a complex system largely depend on the specific way in which its elementary components interact. The last twenty-five years have seen the birth and development of the multidisciplinary field of Network Science, wherein a variety of distributed systems in physics, biology, social sciences and engineering have been modeled as networks of coupled units, in the attempt to unveil the mechanisms underneath their observed functionality. There is, however, a fundamental limit to such a representation: networks capture only pairwise interactions, whereas the functioning of many real-world systems not only involves dyadic connections, but rather is the outcome of collective actions at the level of groups of nodes. For instance, in ecological systems, three or more species may compete for food or territory, and similar multi-component interactions appear in functional and structural brain networks, protein interaction networks, semantic networks, multi-authors scientific collaborations, offline and online social networks, gene regulatory networks and spreading of consensus or contagious diseases due to multiple, simultaneous, contacts. Such multi-component interactions can only be grasped through either hypergraphs or simplicial complexes, which indeed have recently found a huge number of applications. In this report, we cover the extensive literature of the past years on this subject, and we focus on the structure and dynamics of hypergraphs and simplicial complexes. These are indeed becoming increasingly relevant, thanks to the enhanced resolution of data sets and the recent advances in data analysis techniques, which (concurrently and definitely) have shown that such structures play a pivotal role in the complex organization and functioning of real-world distributed systems.

## Keywords

Complex networks
Higher-order interactions
Hypernetworks
Simplicial complexes
Structural properties
Emergent collective dynamics
Processes on networked systems

## 1. Introduction

Extended physical systems can be categorized in two main classes: continuous and distributed.
By continuous ones we mean those systems whose extension is in a physical space (or, in the case of delayed dynamical systems, in time), and which are generally modeled by means of partial differential equations. In these cases, instabilities in space and/or in time lead to the spontaneous (or forced) emergence of “patterns”, which may coexist or compete in the physical domain, giving rise to a very complex spatiotemporal dynamics. One such example is Rayleigh–Bénard convection in fluids under the effect of gravity and heated from below [1], [2]. There, patterns spontaneously arise as those observed on the surface of the Sun, or those verified in the mantle of the Earth, or (when furthermore rotation is acting as a driving force) like the Great Red Spot seen in Jupiter or like the clouds that are forming in the Earth atmosphere. Other examples are nonlinear optical systems confined within cavities with high Fresnel numbers which give rise to pattern formation and competition [3], or solidification of materials which may result in intricate patterns such as those observed in snowflakes and dendritic crystals.
At variance, distributed systems not necessarily have extension in space, but can be seen as formed by a myriad of elementary (or unit) components, each interacting with all the others via a (in general complex and intricate) wiring of connections. Therefore, each system’s component may be identified by a node, or vertex, of a network (sometimes associated to a vector state which evolves following a set of ordinary differential equations), and the interplay with the rest of the units may be represented by the edges, or links, of that network. Thanks to their interactions, the unitary components organize into an emergent collective dynamics, whose complexity, beauty, richness and harmony largely depend on the topological properties of the underneath graph’s structure of connections.
Examples of distributed systems which have been modeled as networks are the brain (where nodes are neurons and links are synaptic connections), the WWW (where nodes are web-pages and links are hyperlinks from a page to another), the air transportation system (where nodes are airports and links are flights), food webs and ecological systems (where nodes are species and links are ecological interactions), scientific coauthorship (where nodes are scientists and links are co-Authored papers), just to mention a few cases. In fact, we ourselves, as individuals, can be seen as the units of a network of social relationships of different kinds and, as biological systems, as the delicate result of a very complicated network of biochemical reactions.
The last twenty-five years have seen the birth and development of the multidisciplinary field of Network Science, wherein a variety of distributed systems in physics, biology, social sciences and engineering have been mapped into networks of coupled units, in the attempt to discover the mechanisms underneath their observed functionality [4], [5], [6].
Network Science is nowadays one of the hottest and most successful research fields, with interdisciplinary applications in areas as diverse as systems engineering, systems biology, nonlinear and statistical physics, applied mathematics, genetics, meteorology, econophysics, neuroscience and personalized medicine.
There exists, however, a fundamental limit in the network representation of a distributed system: networks capture only pairwise interactions, whereas, in general, the way many real-world systems interplay at their component level not only involves dyadic connections, but rather entails “group interactions” among the network’s nodes. In other words, when one adopts a network representation of a system, the overall action of the entire system on each unitary component is assumed to be always exhaustively described by (or factorizable into) a combination of pairwise interactions.
Such a hypothesis may indeed find justification when one studies certain types of processes. For instance, when the nature of the interaction is linear, in virtue of the superposition principle one has that the net effect caused by the simultaneous interplay with more than one unit is always equivalent to the sum of the effects that would have been caused by each dyadic interaction individually.
However, as soon as the assumption of linear coupling is abandoned, that hypothesis is very short in representing faithfully many circumstances. Indeed, from functional [7], [8], [9] and structural [10] brain networks to protein interaction networks [11], to semantic networks [12] and co-Authorship graphs in science [13] there are a lot of practical situations which simply cannot be factorized in terms of pairwise interactions [14] and where, instead, higher-order interactions have to be taken into account for a proper description of the system.
Let us pay few words to illustrate a couple of examples.
The first is borrowed from social science, and refers to collaboration networks. Imagine to map the entire scientific literature into a giant network where nodes are fellow researchers of different disciplines and all published co-Authored papers define the interactions among scientists, and let us focus on three specific scientists, say
i,j, and
k. Let us furthermore conjecture that the literature contains at least a paper co-Authored by
i and
j but not co-Authored by
k, and at least another two papers: one co-Authored by
j and
k but not by
i, and one co-Authored by
i and
k but not by
j. The case is mapped into a triangle which connects the three nodes
i,j, and
k, and which is indeed factorizable as the sum of the three dyadic interactions. A completely different situation occurs when the literature does not contain such three papers, but contains a single published Manuscript which is co-Authored by
i,
j, and
k together. Also in this latter case a triangle would represent the interaction of the three scientists, but now such a triangle would stand for a genuine triadic interaction.
As a second example, one can think of species interactions within ecological webs. Ecologists have identified several types of dyadic (species-species) interactions as, for instance, mutualism, commensalism, competition, predation, and parasitism. Besides such two species interactions one can have, however, other more complicated situations, wherein a given species acts as a catalyst (or as a inhibitor) of the interaction between two other species, even without directly interacting with them. Imagine two species of herbivores which compete for the same food, and another species of small animals (insects, for example) which, when they die, are acting as a fertilizer for the grass on which the two species feed. Once again, this three-species interaction cannot be factorized into the sum of three dyadic interplays.
Appearance of such multi-component actions is much more generic than the two described examples, as higher-order interactions abundantly take place in functional and structural brain networks, protein interaction networks, semantic networks, offline and online social networks, gene regulatory networks and spreading of consensus or contagious diseases due to multiple, simultaneous, contacts.
To properly describe these, and many other cases, one has then to map many-body interactions into proper representations by the use of hypergraphs and/or simplicial complexes. While hypergraphs and simplicial complexes are known in classical graph theory since longtime [15], the availability of new data sets and the recent advances in topological data analysis techniques [16] renewed the interest of the scientific community [17]. On the side of modeling, progresses were made in the last years in extending to hypergraphs standard graph models, such as random graphs models [18], the configuration model [19], models of network growth [20], [21] and activity driven models [22]. Furthermore, it was shown that these kinds of group interactions have a very important role in the collective dynamics and processes taking place in networked systems, such as spreading of ideas and social contagion [23], [24], evolutionary games [25], and synchronization of networked dynamical systems [26], [27], in some cases even determining the emergence of novel states, or states which would be inherently prohibited with only dyadic interactions.
Even though some nice reports on this subject are already available in the literature under the form of monographs [28] or books [29], [30], [31], [32], [33], higher-order networks are one of the most rapidly emerging field of modern research with interdisciplinary applications in several areas. A systematic study of structural properties in network systems with higher-order interactions is still missing, as well as a complete investigation on these systems’ functionality, dynamics and control is still at its infancy.
That is why we decided to embark ourselves in the present journey, with the belief that the modern advancements regarding the intimate relationship between the structure and dynamics of network beyond pairwise interactions, as well as the many novel dynamical behaviors which have been recently unveiled as entirely due to the presence of such group interactions, are justifying the additional effort of recapitulating them within a comprehensive review paper which has, furthermore, the pretention of indicating points that are still open and delineating some promising directions which are, in our opinion, still unexplored and merit instead to be pursued in future years.
Besides the present introduction, the report is, therefore, organized in other 4 Chapters.
In Chapter 2, we start by discussing the mathematical formalism of hypergraphs and simplicial complexes, with a special attention to review the various measures which are used to characterize the connectivity structure of these objects, from generalized degrees and Laplacians, to distances, efficiency, vulnerability, resilience, centrality and ranking. The Chapter ends with an account on preferential an non preferential models that are used to generate hypernetworks endowed with specific scaling properties.
Chapter 3 is devoted to discuss the effects of such higher order interactions when different processes are studied in hypernetworks, from evolutionary games to spreading processes and social contagion.
In Chapter 4, we review the main features characterizing the emerging collective dynamics of networks with higher order interactions, with a specific attention to the setting of synchronized states among identical and non identical dynamical systems, and to consensus dynamics. The Chapter ends with an account on methods for the control of hypernetworks, and on hypergraph controllability.
Finally, in Chapter 5 we report our conclusions, and discuss few open points which constitute (in our opinion) real challenges for the future years.

## 2. The mathematical formalism of networks with higher order interactions

The discipline of Network Science is strongly intertwined with the study of heterogeneity in real-world systems, which led scientists to realize the multifold roles that nodes play in complex networks, establishing a vision that surpassed that of the classical Erdős–Rényi random-graph models [4]. The next breakthrough happened when the scientific community explored networks in which not only nodes, but also links, could have a different nature, thus starting the development of multilayer network models [34]. We are currently progressing towards a third milestone, which is the consistent study of networks whose links are not limited to joining only two nodes, but rather can represent simultaneous interactions between any finite number thereof [28]. These new models are called higher-order networks, and their underlying structure is that of hypergraphs and simplicial complexes [28], [30], [35], [36].
Such structures arise naturally when one describes relationships between groups of more than two nodes, and thus they are the natural generalization of the concept of graph. The extension to higher-order networks of the numerous concepts, measures and observables defined in classical Network Science (from degree to directionality, from distances to mesoscales, just to mention a few) is far from being an obvious task and, in fact, is one of the currently most active areas of research for the community of applied mathematicians. In some cases, the current state of the art is still in its infancy, and there are many problems that still need solutions. Therefore, in this chapter we decided to focus on the few well-established concepts, and to summarize the results that have been used in studies of dynamics and processes taking place on higher-order networks, which are reviewed in the following chapters. At the same time, here we establish the notation we use throughout the report and include a final section on models and generating algorithms for hypergraphs and simplicial complexes. A table with the list of major symbols used is enclosed. Note that, although we believe that the definitions used here are the most common and most generally accepted, differences are sometimes encountered in the literature. For example, some Authors have defined the degree of a node in a higher-order network not as the overall number of hyperedges that contain it, but as a vector whose components are the number of hyperedges of a given cardinality that contain it (see Ref. [28] and references therein).

### 2.1. Basic and preliminary definitions

A hypergraph
H is a pair
H=(X,E), where
X is a finite set of nodes (also called vertices)
X=1,…,N and
E=h 1,h 2,…,h n is a family of nonempty subsets of elements of
X [28], [29], [30], [35], [36], [37], [38], [39]. These subsets are called hyperedges or hyperlinks, and they represent an interaction taking place between elements of
X (see panel (a) of Fig. 1). Hypergraphs are the natural structures that underlie higher-order networks because any higher-order interaction can be readily made to correspond to a hyperedge
h whose cardinality
|h| is greater than 2 and equals the number of elements that participate in it.
A hyperedge consisting of a single node is called a loop. In the following, except if otherwise explicitly stated, we will only consider loopless hypergraphs.

1. Download: Download high-res image (240KB)
2. Download: Download full-size image

Fig. 1. Hypergraphs and simplicial complexes. Schematic illustration of a hypergraph [panel (a)] and of a simplicial complex [panel (b)]. The hypergraph is formed by 5 nodes, one red hyperedge of size 2, two green hyperedges of size 3 and one yellow hyperedge of size 4. The simplicial complex is made of 19 nodes, two simplices of cardinality 1 (links), 8 simplices of cardinality 2 (green triangles) and 2 simplices of cardinality 3 (yellow tetrahedra).

A loopless hypergraph that has no repeated hyperedges is called a simple hypergraph. In simple hypergraphs,
E is not just a family of hyperedges but a set thereof.
If
|h|=k the dimension of
h is said to be
k−1.
The order of
H is the cardinality of
X, and its size is
|E|. In other words, the order of a hypergraph is the number of nodes it contains, and its size is the number of hyperedges that form it.
Two hyperedges
h i and
h j are said to be incident if
h i∩h j≠0̸, i.e., if they have at least one node in common.
Two nodes are said to be adjacent if there is at least one hyperedge that contains both.
The collection of hyperedges incident on a node
i,
E(i)=h|i∈h is called the star of
i.
Two nodes
i and
j are said to be similar if their stars are the same, i.e., if
E(i)=E(j). Equivalently, two nodes are similar if they participate in all the very same hyperedges.
The degree of node
i in the hypergraph
H is denoted by
deg H(i), or simply
deg(i) if there is no ambiguity, and is the number of hyperedges that contain it, which equals the size of its star. In formulae,
deg(i)=|E(i)|. The largest degree amongst the nodes in
H is denoted by
Δ(H).
The rank of
H is the maximum cardinality of its hyperedges, i.e.,
rank(H)=max h∈E|h|. Similarly, the co-rank of
H is the minimum cardinality of its hyperedges, i.e.,
co-rank(H)=min h∈E|h|. If
rank(H)=co-rank(H)=k, then all the hyperedges have the same cardinality, as they all involve the same number of nodes
k, and the hypergraph is called
k-uniform, or simply uniform.
If all the nodes in the hypergraph have the same degree
k, the hypergraph is said to be
k-regular, or simply regular.
Given two hypergraphs
H=(X,E) and
H′=(X′,E′), if
X′⊆X and
∀h′∈E′ there is only one
h∈E, such that
h′⊂h, then
H′ is called a sub-hypergraph of
H. In other words, a sub-hypergraph of
H is a hypergraph on a subset of the nodes of
H such that for each hyperedge of the sub-hypergraph there is only one hyperedge of
H containing it. Note that
H=(X,E) is a sub-hypergraph of itself. Additionally, if
Y⊂X, the induced sub-hypergraph by
H on
Y is the hypergraph
H Y=(Y,E Y) where
E Y is defined by
h′∈E Y⇔∃h∈E such that h′=h∩Y.
Thus, the difference between a sub-hypergraph and an induced sub-hypergraph of a same hypergraph
H is that former may exclude hyperedges formed by some of its nodes that are present in
H, whereas the latter may not.

#### 2.1.1. Paths, distances and connectedness of hypergraphs

Given a hypergraph
H, a walk
W of length
r from node
x to node
y is an alternating node-hyperedge sequence
(x=k 1,h 1,k 2,h 2,…,k r,h r,k r+1=y) such that
k i≠k i+1 and
k i,k i+1∈h i. If the sequence of 2-tuples of consecutive elements of
W,
((k 1,h 1),(k 2,h 1),(k 2,h 2),…,(k r,h r),(k r+1,h r)), has no repeated elements, then
W is called a trail. Note that this definition effectively generalizes that of a trail in a graph, by requiring that the induced dyadic sub-interactions used within each hyperedge to pass from each node in
W to the next be different. A walk in which all the nodes are distinct is called a path. Note that in simple graphs the requirement of having distinct nodes implies that all the traversed edges are also different; however, this implication does not hold in hypergraphs, and repeated hyperedges are allowed in paths. Finally, if
x=k 1=k r+1=y the path is called a cycle.
A hypergraph is said to be connected if any pair of nodes is connected by at least a walk.
The distance between two nodes
i and
j in the hypergraph
H is the minimum length of a path connecting them, and it is indicated by
d H(i,j), or simply
d(i,j) if there is no ambiguity. If there is no path between
i and
j then their distance is defined to be infinite. The diameter of a hypergraph is defined as the largest distance between any two nodes.
If a hypergraph
H is not connected, its connected components, or simply components, are its maximal-order connected induced sub-hypergraphs. Informally, these are the individual connected sub-hypergraphs that partition
H and that are such that no walk exists between a node in one component and a node in another. More formally, they are the equivalence classes induced by the equivalence relation
R defined on the set of nodes by
(i,j)∈R⇔i=j or there is a walk between i and j.

#### 2.1.2. Simplices and simplicial complexes

A special type of hypergraphs that have attracted considerable attention in the study of complex systems are simplicial complexes. By definition, a simplicial complex is a hypergraph
H=(X,E) such that, for every
h∈E, if
h′⊂h and
|h′|⩾2, then
h′∈E. In other words, simplicial complexes are hypergraphs in which the presence of a hyperedge of any cardinality implies the existence of all hyperedges of any smaller cardinality on the same nodes (see panel (b) of Fig. 1). The origin of their name lies in the traditional terminology that defines a
p-simplex as the collection of a set of size
p+1 and all its proper subsets of size greater than 1. Note that this is in analogy to hyperedges, so that the largest hyperedge in a
p-simplex has dimension
p although it should be pointed out that in Algebraic Topology the hyperedges such that
|h|=1 are also included in the classic definition of simplicial complexes. Finally, the dimension of a simplicial complex
S is the maximum dimension of its simplices, which equals
rank(S)−1.
Although a
p-simplex can be considered as an abstract object determined by its nodes, traditionally simplices have been viewed through a geometric representation as polyhedra in a multidimensional space. Thus, for example, a 1-simplex is a line between two nodes (an edge), a 2-simplex is a triangle and a 3-simplex is a tetrahedron. More generally, a
p-simplex is the
p-dimensional analogue of a triangle, i.e., the convex envelope of a set of
p+1 independent nodes in a general position. It is also worth noting that simplicial complexes are similar to those complexes considered in algebraic topology as the working objects of the classic theory of homology and homotopy [40], [41], [42].
Finally, the simplicial network closure of a hypergraph
H=(X,E) is the simplicial complex
H̄=(X,Ē) on the same nodes obtained by replacing each hyperedge of
H with its induced simplex, without repetitions:
Ē=⋃h∈E P(h)∖h′∈P(h)||h′|<2,where
P(h) denotes the power set of
h. For instance, if we take
X=i,j,k,l and
E=i,j,k,j,k, then
Ē=i,j,k,i,j,i,k,j,k.

#### 2.1.3. Matrix and hypermatrix representations of hypernetworks

One of the most natural ways to represent group interactions is by using bipartite graphs: given a simple hypergraph
H=(X,E), build the bipartite graph
B(H)=(X∪E,E) defined by the two node sets
X and
E, and containing all and only the edges
(i,h) such that
i∈X,
h∈E and
i∈h. In other words, the nodes of
B represent the nodes and the hyperedges of
H, and the edges of
B link the nodes representing nodes of
H with those representing the hyperedges of
H that contain them (see Fig. 2).
An alternative representation of a hypergraph can be obtained by means of its incidence matrix. Originally introduced by Kirchhoff in 1847 for the description of electrical circuits and their applications, incidence matrices are an established way to describe relationships between two different classes of objects. The incidence matrix
I(H) of a hypergraph
H is a matrix whose rows correspond to the nodes of
H and whose columns correspond to its hyperedges. As an example, consider the hypergraph
H with nodes
X=1,2,3,4,5,6 and hyperedges
h 1=1,2,4,h 2=2,3,h 3=1,2,3,h 4=3,5,6 (illustrated in the left panel of Fig. 2). Its incidence matrix is
I(H)=1 0 1 0 1 1 1 0 0 1 1 1 1 0 0 0 0 0 0 1 0 0 0 1.It is straightforward to verify that
I(H) can be obtained directly from the adjacency matrix
A(B(H)) of the bipartite graph
B(H). In fact, one can write
A(B(H)) as the block matrix
A(B(H))=0 I(H)I T(H)0,where the superscript T indicates transposition.

1. Download: Download high-res image (285KB)
2. Download: Download full-size image

Fig. 2. Bipartite representation of hypergraphs. The hypergraph
H=({1,2,3,4,5,6},{h 1,h 2,h 3,h 4}) (left panel) is mapped into a bipartite graph (right panel) in which the top nodes correspond to the nodes of
H, the bottom nodes correspond to the hyperedges of
H and the edges link the nodes to all the hyperedges to which they belong.

A further useful matrix is the frequency matrix of relations of
H,
A(H), which is an
|X|×|X| matrix defined by (1)
A i,j(H)=deg(i)if i=j|h|i,j∈h|if i≠j.From its definition, it is clear that its diagonal elements contain the degrees of the nodes and its off-diagonal elements contain the number of times the two corresponding nodes are found in the same hyperedge. So, for the same example hypergraph as before, the frequency matrix of relations is
A(H)=2 2 1 1 0 0 2 3 2 1 0 0 1 2 3 0 1 1 1 1 0 1 0 0 0 0 1 0 1 1 0 0 1 0 1 1.
With this formalism, it is possible to verify that
I(H)I T(H)=A(H),and
I T(H)I(H)=Ã(H),where
Ã(H) is an
|E|×|E| matrix defined by
Ãi j(H)=|h i|if i=j|h i∩h j|if i≠j,so that its diagonal elements contain the cardinalities of the hyperedges and its off-diagonal elements contain the number of nodes that the two corresponding hyperedges have in common. Finally, using the previous definitions, one can write
A 2(B(H))=A(H)0 0 Ã(H).
Taking the matrix
A(H) and replacing the diagonal elements with zeros one obtains the adjacency matrix of the hypergraph
A(H). Formally, this can be written as
A(H)=A(H)−D,where
D is the diagonal matrix containing the degrees of the nodes. If the hyperedges are weighted, one can obtain a weighted adjacency matrix by exploiting the definition of
A(H) in terms of the incidence matrix. Then, introducing the diagonal matrix
W containing the weights of the hyperedges, and treating it as a metric, one can write [43]
A(H)=I(H)W I T(H).
Finally, higher-order matrix representations of higher-order networks have also been introduced in Refs. [28], [30], [35], [36], [44], [45], [46]. A common formalism is to define one
N×⋯×N︸k times adjacency hyper-matrix
A(k)(H) for each cardinality
k of the hyperedges, so that
A i 1,…,i k(k)=1 if i 1,…,i k∈E,0 otherwise.
In the following, the adjacency hyper-matrix of the hypergraph is also called adjacency tensor.

#### 2.1.4. Graph representations

One of the most frequently pursued objectives when working with higher-order networks is to try to approximate them by conventional graphs or by regular hypergraphs and compare the analysis results obtained for the original higher-order graph with those obtained on the approximating structures. Several approximations can be defined using common techniques of geometric convex analysis. In particular, approximations of a higher-dimension convex body in
R n can be obtained by using the two basic tools of sections and projections in lower dimensional spaces
R k with
k<n.
Although the relevance of this type of approximations is indisputable, the best approximation one may adopt often depends on the specificities of the problem under study. Therefore, the correct question is not which conventional graph best approximates a hypernetwork, but rather which standard graph best approximates the hypernetwork under the specific conditions of the process or dynamics under study.
Sectioning effectively consists in partitioning a hypergraph
H of rank
R into a set of maximal
k-uniform hypergraphs
s k(H)=(X,E k), which are called the
k-skeletons of
H, so that for each hyperedge cardinality
k,
E k=h∈E||h|=k, and
H=⋃k=2 R s k(H). The
k-projection of
H, instead, is the
k-uniform hypergraph
π k(H)=(X,J k) where
i 1,…,i k∈J k⇔∃h∈E|i 1,…,i k∈h.In other words, to form the
k-projection of a hypergraph, first remove all hyperedges of cardinality smaller than
k, and then replace each hyperedges
h of cardinality
x>k with hyperedges between all possible choices of
k nodes out of the
x nodes in
h. A weighted version of such projections is also possible, by providing each hyperedge in the projection with a weight equal to the number of hyperedges of
H that contain all its nodes. Note that
π 2(H) is the classic projection of
H into a simple graph.
Skeletons and projections have some basic properties, which we list below:

1. Download: Download high-res image (242KB)
2. Download: Download full-size image

Fig. 3. Walks on hypergraphs. Illustration of 1-walks [panel (a)] and of 2-walks [panel (b)].

- 1.
  s k(H)⊆π k(H)∀k∈N,2⩽k⩽R, with the equality holding for
  k=R.
- 2.
  ⋃k=2 R π k(H) is a simplicial complex
  ∀k∈N,2⩽k⩽R.
- 3.
  ⋃k=2 R π k(H) is the smallest simplicial complex that contains
  H.
- 4.
  H is a simplicial complex if and only if
  s k(H)=π k(H)∀k∈N,2⩽k⩽R.
- 5.
  H⊆⋃k=2 R π k(H), with the equality holding if and only if
  H is a simplicial complex.
- 6.
  H is a simple graph if and only if
  s 2(H)=π 2(H)=H.
- 7.
  deg π k(H)(i)⩾deg H(i)⩾deg s k(H)(i) for all nodes
  i and
  ∀k∈N,2⩽k⩽R.
- 8.
  d H(i,j)⩽d π k((H))(i,j)⩽d s k(H)(i,j) for all pairs of nodes
  i and
  j and
  ∀k∈N,2⩽k⩽R.
- 9.
  d H(i,j)=d π 2(H)(i,j) for all pairs of nodes
  i and
  j.

A consequence of the last property is that the metric structure of
H and that of
π 2((H)) coincide, and therefore all the structural metric parameters are the same for
H and for its projection to a simple graph. Also note that the property holds for the specific definition of distance we have given above. However, other definitions of walks in hypergraphs exist [28]. For example, one can define a walk as an ordered list of hyperedges such that the intersections of consecutive hyperedges contain either only one or at least two nodes, yielding so-called 1-walks and 2-walks, respectively (see Fig. 3). In this case, a different metric structure is obtained with consequent differences in the metric properties.
A very important concept is that of the line graph of a hypergraph. Line graphs were originally introduced by H. Whitney in 1932 [47], and they have been studied since then under different names [48], [49], [50]. Notably, the study of these objects and their applications has been increasingly gaining significance, especially in recent years [39], [51], [52], [53], [54], [55], [56], [57]. Given a hypergraph
H=(X,E), its line graph is the graph
L(H)=(E,E′), where
h i,h j∈E′⇔h i∩h j≠0̸. In other words, each node of the line graph corresponds to a hyperedge, and a link in the line graph exists if and only if the two corresponding edges in the original hypergraph are incident. The line graph is therefore a particular type of intersection graph [57].
Finally, the dual of a hypergraph
H(X,E) is the hypergraph
H∗=(E,E′) where
E′=h 1′,h 2′,…,h|X|′ and
h i′=h|i∈h. Thus, each of the nodes of the dual corresponds to a hyperedge of the original graph, and a hyperedge exists between all the nodes whose corresponding hyperedges in the original hypergraph share a node. Note that this marks a substantial difference with respect to duals in simple graphs. In fact, in classic graph theory, duals refer to drawings of planar graphs, whereas here the dual is effectively a different type of line graph, made possible by the higher dimension of the hyperedges. Also, it is possible to verify that
(H∗)∗=H up to isomorphism. Moreover, if
I(H) is the incidence matrix of
H, then
I T(H) is the incidence matrix of
H∗. Last, it can be shown that the line graph of a hypergraph is isomorph to the 2-projection of its dual i.e.,
L((H))≅π 2(H∗) (see Fig. 4).

1. Download: Download high-res image (195KB)
2. Download: Download full-size image

Fig. 4. Line graph and projection graphs. Panel (a) is an illustration of a hypergraph with 6 nodes and 4 hyperedges. Panel (b) is a sketch of the corresponding line graph, which is isomorph to the 2-projection of its dual. Panel (c) shows the 2-projection.

### 2.2. Centrality of nodes in a hypergraph

Several types of centrality measures of nodes in higher-order networks have been introduced in the literature. A first one consists in the extension of the eigenvector centralities used for simple graphs to the case of hypergraphs [52], [58], [59]. Specifically, Ref. [59] defines the eigenvector centrality
c(i) of a node in a
k-uniform hypergraph
H=(X,E) as a scalar field such that the values it takes on each node must be proportional to the sum of the values that an arbitrary function
F:R k−1⟶R takes on all the sets of its neighbors. In formulae, one has
c(i)=1 λ∑i,j 1,…,j k−1∈E F(c(j 1),…,c(j k−1)).For the particular choice of
F(x 1,…,x k−1)=x 1+⋯+x k−1, one recovers the classic eigenvector centrality of the weighted projection
π 2(H), while if
F(x 1,…,x k−1)=x 1 x 2⋯x k−1 then one obtains the
Z-eigenvector centrality introduced in Ref. [52]. Other choices are also possible, such as the
p-dimensional metrics
F(x 1,…,x k−1)=‖(x 1,…,x k−1)‖p. In all these cases, the existence and uniqueness of the centrality values can be proved by means of the non-linear Perron–Frobenius theory [60]. Note that the above formalism can be extended to general hypergraphs of rank
R by considering a family of functions
F k:R k−1⟶R and using
c(i)=1 λ∑k=2 R∑i,j 1,…,j k−1∈E F k(c(j 1),…,c(j k−1)).Also in this case, the same considerations about existence and uniqueness apply.
More recently, a new measure of node centrality in hypergraphs has been introduced, which is no longer a scalar value, but a vector with dimension equal to the rank of the hypergraph minus one [46]. The procedure to calculate this vector centrality starts with constructing the corresponding line graph and calculating the classic eigenvector centrality
c(h) for each of its nodes
h, which correspond to the hyperedges of the original hypergraph. Then, one defines the vector centrality of each node
i as the
R−1-dimensional vector
c i whose components are given by (2)
c i k=1 k∑h||h|=k,i∈h c(h).In other words, the
kth component of the vector centrality of node
i is the sum of the eigenvector centralities of all the nodes of the line graph of the hypergraph that correspond to hyperedges of size
k that contain
i. Ref. [46] shows rigorously that this quantity is related to the classic eigenvector centrality in graphs, but it introduces a significant added value in higher-order networks. In particular, it is able to discriminate between the roles that the same node may play at different orders of interactions, providing information that cannot be inferred from scalar measures.

### 2.3. Derivative graph of a hypergraph

One of the classical methods for comparing the degree of coincidence or similarity between two sets
A and
B is the Jaccard index [61], which is defined by (3)
J(A,B)=|A∩B||A∪B|.Different generalizations of the Jaccard index have been introduced in the literature [62], [63], [64], [65], [66], [67], [68], including the overlap index
C(A,B), defined as
C(A,B)=J(A,B)|A∩B|min|A|,|B|.Note that the values of both the Jaccard index and the overlap index are always between 0 and 1.
Recently, Ref. [37] has introduced the derivative of a hypergraph with respect to a pair of nodes
i and
j, a concept that allows one to quantify the similarity of pairs of nodes in the hypergraph structure. Specifically, if
H is a hypergraph and
A is its frequency matrix of relations, the derivative of
H with respect to the nodes
i and
j is the numerical value (4)
∂H∂i,j=A i,i−A i,j+A j,j−A i,j A i,j=A i,i−2 A i,j+A j,j A i,j.If two nodes have disjoint stars, the derivative of the hypergraph with respect to them is defined to be
∞. The derivatives of a hypergraph provide a direct comparison of the role played by two of its nodes: the smaller the derivative, the larger the similarity between the nodes with respect to which it was computed. So, if two nodes
i and
j are similar, i.e., if they participate in all the same hyperedges, then
∂H∂i,j=0. As an example, one can conveniently write the derivatives of the hypergraph in panel (a) of Fig. 5 in matrix form:
0 1 2 3 1∞∞1 2 0 1 2∞∞3 1 0∞2 2 1 2∞0∞∞∞∞2∞0 0∞∞2∞0 0.Note that a relationship exists between the Jaccard index of the stars of two nodes and the derivative of a hypergraph with respect to them [38]: (5)
J(E(i),E(j))=1 1+∂H∂i,j.Finally, if a hypergraph has rank
R one can show that if its derivative with respect to two nodes is finite, then its upper bound is
2 R−1, which implies that
1 2 R is a lower bound for the Jaccard index of the stars of the two nodes.
Using the derivatives, one can build two more graphs associated to a hypergraph, namely the derivative graph
∂H and the homogeneity graph
H G(H) [37]. To compute
∂H, first created a fully connected graph with the same nodes as
H. Next, if the derivative of
H with respect to two nodes is infinite, delete the edge between the two nodes. Then, if the derivative with respect to two nodes is 0 (i.e., if the two nodes are similar in
H), contract the edge between them. Finally, assign to each remaining edge a weight equal to the derivative of
H with respect to them. To build the homogeneity graph, follow the same procedure, but give each edge the inverse of the weight it has in the derivative graph.

1. Download: Download high-res image (209KB)
2. Download: Download full-size image

Fig. 5. Derivative of a hypergraph. By computing the derivatives of the hypergraph in panel (a) with respect to any two pairs of nodes, one can build its derivative graph
∂H, and its homogeneity graph, which are shown in panels (b) and (c), respectively. Note that the edge between the similar nodes 5 and 6 has been contracted, “collapsing” the two nodes into one, and that no edge is present between nodes that yield an infinite derivative.

### 2.4. Generalized Laplacians and entropy

The classic graph Laplacian is a matrix that provides an algebraic description of a graph, in that its spectral properties reveal a number of important topological features [69], [70], [71], [72]. The extension and generalization of Laplace operators to hypergraphs and simplicial complexes has been the object of several recent studies [27], [28], [30], [45], [72], [73], [74], [75].
A first generalization can be obtained by means of the adjacency matrix
A(H). Then, introducing the Laplacian degree of node
i as
D(i)=∑j=1 N A i,j (see Ref. [75]), one can define the Laplacian matrix of
H as in Refs. [74], [76]
L(H)=D−A(H),where
D is the diagonal matrix of the Laplacian degrees. Since
L is a zero-row-sum matrix, its smallest eigenvalue
λ 0=0, which has a corresponding eigenvector whose components are all equal. Moreover, in analogy with the classic graph Laplacian, the multiplicity of
λ 0 is equal to the number of connected components of
H.
A different approach, suitable for simplicial complexes, starts from the adjacency hyper-matrix formalism described in Section 2.1.3, see Ref. [44]. Given a simplicial complex
S of dimension
D with
N nodes, consider its set of
D adjacency hyper-matrices
A(r), with
2⩽r⩽D+1. Then, for each node
i and for each
1⩽d⩽D, define the generalized
d-degree of node
i,
k i(d), as the number of
d-simplices it participates in:
k i(d)=1 d!∑i 2=1 N∑i 3=1 N⋯∑i d+1=1 N A i,i 2,i 3,…,i d+1(d+1).Next, for each edge
i,j, define its generalized
d-degree
k i,j(d) as the number of
d-simplices to which the edge belongs:
k i,j(d)=A i,j(2)for d=1,1(d−1)!∑i 3=1 N∑i 4=1 N⋯∑i d+1=1 N A i,j,i 3,i 4,…,i d+1(d+1)for 2⩽d⩽D.Finally, define the generalized Laplacian of order
d as the matrix
L(d) defined by
L i,j(d)=d!k i(d)if i=j,−(d−1)!k i,j(d)if i≠j and A i,j(2)=1,0 otherwise.Note that if
D=2, i.e., if
S is a simple graph, then
L(1) is the traditional graph Laplacian.
A related concept is that of graph entropy, first introduced by J. Körner in 1973 [77] and later extended in Ref. [30] to quantify the amount of information carried by a hypergraph when modeling systems in different areas of science and engineering [28], [30]. Given a hypergraph
H on
N nodes, consider the matrix
L′(H)=L(H)∑i=1 N D(i).Its eigenvalues
μ 0,μ 1,…,μ N−1 are related to those of
L,
λ 0,λ 1,…,λ N−1, by the relation
μ i=λ i∑i=1 N D(i).But
∑i=1 N D(i) is also the trace of
L, which means that it is equal to the sum of its eigenvalues. Thus, it is
μ i=λ i∑i=1 N λ i,which means that
μ 0,μ 1,…,μ N−1 is a discrete probability distribution. This allows one to use the eigenvalues of
L′ to define the algebraic hypergraph entropy
S(H) [30] as
S(H)=−∑i=1 N−1 μ i log 2(μ i),where the sum starts from 1 because the first eigenvalue is always 0. An expression of the spectral entropy specific for simplicial complexes has been introduced in Ref. [78].
Table 1 contains a list of the major symbols that will be used along the entire report.

Table 1. List of the major symbols used in this report.

H | Hypergraph
--- | ---
X | Set of nodes in H
E | Set of hyperedges in H
N | Number of nodes in X
h | Hyperedge
|h| | Cardinality of a hyperedge
deg H(i) | Degree of node i in the hypergraph H
Δ H | Largest degree amongst the nodes in H
d H(i,j) | Distance between two nodes i and j in H
S | Simplicial complex
d−simplex | Simplex of order d, 1≤d≤D
B(H) | Bipartite graph
I(H) | Incidence matrix
A(H) | Adjacency matrix of the hypergraph
R, rank (H) | Rank
H∗ | Dual of a hypergraph
c(i) | Eigenvector centrality of a node i in a hypergraph H
L | Laplacian matrix
k i(d) | Generalized d-degree of a node i
k i,j(d) | Generalized d-degree of an edge {i,j}
L(d) | Generalized Laplacian of order d
S(H) | Algebraic hypergraph entropy
P(k) | Degree distribution

### 2.5. Models and generating algorithms

During the last years, significant progress has been made in modeling hypergraph structures. In particular, several generating models have been introduced, among which we recall the extensions of some established static graph models (such as random graphs [18], [23], the configuration model [44] and activity-driven models [22]), as well as growth algorithm where random-graph models have been developed by aggregating hypergraphs or simplicial complexes as fundamental building blocks, inspired by quantum networks [19], [79] and complex materials [80].

#### 2.5.1. Preferential and non preferential growing mechanisms

Recently, a simple technique has also been introduced to grow hypergraphs and simplicial complexes by means of a preferential rule, a non preferential one, or a combination thereof [21]. The method generates hypergraphs characterized by a power-law scaling of their degree distribution
P(k), offering full control over the distributions
P(k(d)) of the generalized degrees. For ease of illustration, we consider here simplicial complexes of order 2, and aim to grow a complex of
N nodes endowed with a transitivity coefficient
T=1 [81], [82]. To do so, start at time
t=0 with a seed consisting of a clique of size
N 0<N. At each subsequent time step
m tri edges are selected randomly, with the condition that no pair of chosen edges have a node in common. Then, a new node is added to the graph, and new edges are placed between the new node and the
2 m tri nodes located at the ends of the selected edges, thus generating
m tri new triangles. If the edges are selected uniformly at random, then the probability of choosing any edge at time
t
P(t)=1 N 0(N 0−1)2+2 m tri(t−1).If, instead, the edges are chosen according to a preferential-attachment rule, then the probability of selecting edge
i,j at time
t is proportional to the number of triangles that the edge is part of in that moment:
P i,j(t)=2 k i,j(2)∑i=1 N∑j=1 N k i,j(2),where the factor of 2 accounts for the double-counting in the denominator.
Figs. 6a and 6c show the degree distribution for complexes of size
N=1 0 5 generated in the two ways. In both cases, one can observe a clear power-law scaling, as also found in the majority of real-world networks [4], [5], [6], [84], [85], [86], [87], [88]. Figs. 6b and 6d are visualizations of typical complexes generated using the two mechanisms with
N=200 and
m tri=1. Fig. 7 shows the distribution
P(k(2)) of the generalized degree
k(2) for the two cases. Comparing the two panels of the figure one can see that
P(k(2)) decays exponentially if the complex is grown without preferential attachment, whereas it scales as a power-law for preferential attachment. This implies that the two mechanisms yield substantially different topologies of 2-simplices to the resulting complexes.

1. Download: Download high-res image (1MB)
2. Download: Download full-size image

Fig. 6. Preferential and non-preferential models for the growth of hypernetworks. The non-preferential model (panels a and b) and the preferential-attachment model (panels c and d). Panels a and c show logarithmic plots of the degree distribution
P(k) for different values of
m tri. The data are averaged over 100 realizations of complexes with
N=1 0 5. The dashed are to the analytical predictions given by Eqs. (6), (8). Panels b and d show visualizations of two complexes generated with the described mechanisms
N=200 and
m tri=1. The size of the nodes is proportional to 50 times the square root of their eigenvector centrality [83], the width of each link is proportional to
k i,j 2 and the color of the links relates to the number of triangles.

Figure reprinted from Ref. [21].
Ref. [21] also provides rigorous expressions for
P(k) and
P(k(2)). In the case of no preferential attachment, by indicating the number of nodes of degree
k at time
t with
N(k,t), one can write down the rate equation
Δ N(k,t)=2 m tri∑k k N(k,t)[−k N(k,t)+(k−1)N(k−1,t)]+δ k,2 m tri,where
Δ N(k,t)=N(k,t+1)−N(k,t) and
δ is the Kronecker symbol. Asymptotically, the number of nodes for large
t is
N(t)≈t, and one seeks a solution of the form
N(k,t)=t P(k), where
P(k) is the steady-state distribution. Also, since the total number of edges is approximately
2 m tri t, one has
∑k k N(k,t)≈4 m tri t. Then, the rate equation becomes
P(k)=k−1 k+2 P(k−1)+2 k+2 δ k,2 m tri.The solution of this equation for
k⩾2 m tri is (6)
P(k)=4 m tri(2 m tri+1)k(k+1)(k+2)∼k−3,which fits the data remarkably well (see Fig. 6).

1. Download: Download high-res image (227KB)
2. Download: Download full-size image

Fig. 7. The difference in scaling between non-preferential and preferential-attachment growth. Panel a shows the distribution of
P(k(2)) vs.
k(2) obtained by growing complexes with
1 0 4 nodes with the non-preferential model. Panel b shows the same distribution for complexes grown with preferential attachment. Results are averaged over 100 realizations. The dashed lines are the solutions of Eq. (7) (panel a) and of Eq. (9) (panel b). For the sake of consistency note that
k l in the axis of this figure corresponds to
k(2) in our notation.

Figure reprinted from Ref. [21].
A similar procedure can be followed to obtain the asymptotic distribution
P(k(2)). One first considers the number of edges participating in
k(2) triangles at time
t,
N e(k(2),t), and then writes down its rate equation as
Δ N e(k(2),t)=−m tri N e(k(2),t)N e(t)+m tri N e(k(2)−1,t)N e(t)+2 m tri δ k(2),1±…,where
N e(t)≈2 m tri t and the unwritten terms account for the formation of triangles from two or more linked edges. Moreover, one has that
P(k(2))=N e(k(2),t)/N e(t), and the recursive equation is
P(k(2))=1 3 P(k(2)−1,t)+2 3 δ k(2),1,which admits the solution (7)
P(k(2))=2 3 k(2).The exponential scaling predicted by Eq. (7) is illustrated in Fig. 7a.
If using preferential attachment to select the edges, the analytical treatment is more involved, but it is possible to prove that the closed form for the degree distribution is (8)
P(k)=3 3+2 m tri Γ(k 2)Γ(m tri+2.5)Γ(m tri)Γ(k 2+2.5)∼k−2.5,where
Γ is the gamma function [21]. Furthermore, the closed form for
P(k(2)) is (9)
P(k(2))=3 4∏l=2 k(2)l−1 l+3=3 4 4!k(2)(k(2)+1)(k(2)+2)(k(2)+3)∼(k(2))−4.
The predictions of Eq. (8) are verified very well by the results of numerical simulations (see Fig. 6c), as are those of Eq. (9) (see Fig. 7b).

#### 2.5.2. Generating hypernetworks with desired power-law scaling properties

The generating algorithm can be generalized to a mixed model, through which it is possible to impose any desired power-law scaling in the triangular structure of the network. This is done by fixing the edge-selection probability at time
t to be
P i,j(t)=A 4∑i=1 N deg(i)+B 3 k i,j(2)∑i=1 N∑j=1 N k i,j(2),for some constants
A and
B. The requirements of normalization and well-definedness of the distribution
P i,j(t) lead to the bounds
−1⩽A⩽1 2 and
0⩽B⩽2, and to the condition
A=1 2−3 4 B. The cases
A=0,B=2 3 and
A=1 2,B=0 recover the preferential and non-preferential models, respectively. With this formalism, complexes generated with strictly positive values of
B display scale-free distributions
P(k)∼k−γ and
P(k(2))∼k(2)−γ(2), with exponents given by (10)
γ=1+1 A+B=1+4 2+B
γ(2)=1+2 B.For
B=0, instead, one has
P(k(2))=2/3 k(2). Eq. (10) implies that the values of
γ lie between 2 (
B=2) and 3 (
B=0), whereas
γ(2) is equal to 2 for
B=2, and diverges as
B goes to 0. In turn, this means that by choosing
B to be between 1 and 2, one can arbitrarily obtain values of
γ and
γ(2) between 2 and 3. Fig. 8 shows the (generalized) degree distributions for three distinct values of
B, illustrating that the analytic predictions well approximate the numerically generated data. Finally, Ref. [21] also discusses how this generating algorithm can be extended to general hypergraphs and simplicial complexes of any dimension.

1. Download: Download high-res image (323KB)
2. Download: Download full-size image

Fig. 8. The mixed model. Distributions of degrees (panels a, b and c) and generalized degrees (panels d, e and f) in complexes grown with the mixed model for different values of
B. The data (blue lines) are averages over 100 different realizations of 50000 nodes each. The dotted lines show the scaling with the exponents predicted by Eq. (10). For the sake of consistency note that
k l in the axis of this figure corresponds to
k(2) in our notation.

Figure reprinted from Ref. [21].

## 3. Dynamical processes in networks with higher-order interactions

The use of networks to model structures that support dynamical processes can be traced back to the very birth of graph theory. In fact, the problem of the existence of Eulerian trails in a graph inherently includes the idea of traversing its structure, and of doing so in a very specific way. As we bring the abstraction to the next level, a network no longer has to support an actual physical motion, but it can be taken to simply represent relations of some given nature between elements whose state is a dynamical variable. In this case, the nodes can influence each other if they share an edge, whose strength represents the local intensity of a mutual or directed interaction. It has to be furthermore remarked that the state space of the nodes is not constrained to be continuous. Thus, it is easy to understand how the network paradigm has obtained numerous successes in modeling situations in which the discrete condition of a node can propagate to its neighbors under many possible sets of rules. A typical example of such a process is the spreading of a disease from an infected individual to a susceptible neighbor. Indeed, epidemic spreading is one of the topics that Network Science has contributed most to shed light upon. More generally, this is an instance of a contact process, in which the change in the state of a node is only enabled by the direct influence of its immediate neighbors. Below, we review how our understanding of such processes has benefited from considering interactions beyond the pairwise level, touching upon epidemic spreading, as well as models of opinion formation, social contagion, and games.

### 3.1. Spreading processes

The initial realization that higher-order interactions play a significant role in numerous complex systems did not occur, strictly speaking, in the very near past. However, explicitly accounting for them within models and their related mathematical treatment is a fairly recent development. Some of the first steps in this direction were taken by a number of mathematical epidemiologists, who studied the transmission of infectious diseases within household environments [89], [90], [91], [92], [93], [94], [95], [96], [97], [98], [99]. This setting consists of networks that include local close-contact neighborhoods, which, despite their name, can represent not just actual households, but any other similar environment, such as workplaces, gyms, clubs or, in general, locations where people spend a non-negligible amount of time in close proximity to each other. These studies were quite successful in introducing methods that could both replicate real-world data and provide suggestions as to the strategies best suited for the containment of epidemics, for example via targeted vaccination campaigns and lockdowns [100], [101].

#### 3.1.1. Spreading of epidemics

Classic epidemic studies can be traced back to the work of Kermack and McKendrick [102]. In their basic incarnation, epidemic models group individual in compartments, assuming that individuals can be susceptible (they are healthy and can be infected), infectious (they are infected and can transmit the disease), or recovered (they are healed and immunized). In Susceptible–Infected–Susceptible (SIS) models, a susceptible agent can get the disease and become infectious, to then possibly heal and become susceptible again, whereas in Susceptible–Infected–Removed (SIR) models, when an infectious individual recovers, it becomes immune to the disease. The simplest epidemic models assume homogeneous mixing, whereby any individual can interact with any other, and contagion takes place with a given probability. Under this simplifying assumption, it was possible to show the existence of a threshold below which the epidemics dies out.
Network models of epidemic spreading dropped this assumption of homogeneous mixing, and considered that the interaction of the individuals takes place over a graph. In this setting, contagion was considered as the result of pairwise interactions between the individuals. In the context of epidemics, network science was able to show how the epidemic threshold may change, and possibly disappear, depending on the topological properties of the network, and brought the use of contagion models beyond the classic field of epidemic spreading, see e.g., Ref. [103], [104] and references therein.
However, it is to be noted that, whilst still being foundational works in their own right, such first, seminal articles included the contributions of the household structure in a way that did not directly account for the synergistic effects of multiple exposures, thus substantially reducing the problem to a traditional pairwise network approach.
The first work that, instead, explicitly considered higher-order interactions is that of Bodó and collaborators [105], who studied a SIS model on hypergraphs. Their use of generic hypergraphs allowed them to naturally consider the occurrence of pairwise interactions together with that of higher-order ones, providing at the same time a straightforward way to account for the household structure of a network and to include a degree of control over it. In their work, they assumed that the probabilities of changing state, i.e., of contracting an infection or of recovering from one and becoming again susceptible, are very small. More specifically, given a time interval
Δ t, they expressed the probabilities of node
i changing state within that interval as
1−e−γ i Δ t and
1−e−r i Δ t where
γ i and
r i are the recovery and infection rates, respectively. The peculiarity of the model is that, whilst
γ i is a fixed parameter,
r i is given by (11)
r i=τ∑h∋i f(k h),where
τ is a proportionality constant, the sum is over all the hyperedges
h that include node
i,
k h is the number of infected nodes in the hyperedge
h, and
f is an arbitrarily chosen function of
k h. At variance with previous works, the Authors’ choice for
f is that of a piecewise linear function that is the identity up to a threshold
c, from where it then stays constant. Thus, rather than considering the existence of a minimum value of infected neighbors that allows an infection to occur, they modeled the infection rate as having no threshold, but saturating at a critical amount of exposure. A key result from the simulations of the model is that the size of the infection is generally smaller than or equal to the one observed in the pairwise network obtained by replacing all edges of order greater than 2 with cliques. This is particularly true for household models, where the equality is reached when the saturation threshold equals the greatest household size. Moreover, the effect of higher-order edges is more pronounced in the early stages of the infection, as it is clearly illustrated by the simulation results on the preferential attachment hypergraphs reported in Fig. 3 of Ref. [105].

#### 3.1.2. Spreading of information

As already mentioned, the dynamical features of contact spreading processes are, in fact, general. Therefore, the same framework used to study disease spreading can actually be employed to model the spread of information. With this in mind, Ma and Guo considered a SIR model on hypergraphs [106]. In their vision, infected nodes are “spreaders” of some rumor or piece of information, whereas susceptible nodes are “ignorant” with respect to it and recovered nodes are “stiflers”, as they are unwilling or unable to further its propagation. Then, hyperedges represent local communities in a network, such as social circles that form naturally within the hierarchical structure of a commercial company. On this kind of network, they considered four different spreading mechanisms, which they called the probability transmission, the one-way transmission, the group transmission and the gossip transmission. Of these, the probability transmission is the classic spreading mechanism, in that all ignorant (susceptible) neighbors of each spreader (infected) node can become knowledgeable of the rumor (infected) with a given probability. The other three, instead, assume certainty of spreading but interpolate the number of neighbors of a spreader that become infected at every step, from just 1 in the one-way transmission to all neighbors in the gossip one. Analysis of these four mechanisms shows that, within the expected qualitative differences, high-degree nodes have a large importance for the spreading process, suggesting that they be the main target of control measures if the propagation has to be stopped.
In a somewhat similar manner, Suo et al. considered a SIS model for information spreading on linear
k-uniform hypergraphs [107]. Note that these are a particular type of hypergraph in which all edges connect exactly
k nodes, and any two edges have at most 1 node in common. On these structures, they considered two transmission mechanisms, which they called the reactive process (RP) and the contact process (CP). The former is effectively identical to the probability transmission of Ref. [106], whereas in the latter the infected nodes spread the information to all susceptible nodes within an edge randomly selected among those incident to them. With these governing rules, they found that structural differences only affect the initial stages of the spreading process, and, similarly to the results of Ref. [106], high-degree nodes are the fastest information spreaders. Their main finding is that the two transmission mechanisms result in different steady-state behaviors. More specifically, the CP features a continuous phase transition in the transmission probability
β, so that the number of informed individuals vanishes in the steady state if
β is smaller than a critical value (see Fig. 9).
Note that, unlike the model of Ref. [105], the dynamics in those of Refs. [106], [107] are not affected by the higher order of the edges in the hypergraph, which only acts as a supporting structure. Thus, the results found in these cases should be the same as those obtained by replacing hyperedges with cliques. Also note that, when operating this substitution on a linear hypergraph, one effectively obtains a locally treelike household network.

1. Download: Download high-res image (129KB)
2. Download: Download full-size image

Fig. 9. Presence of a phase transition in the CP transmission mechanism of information on
k-uniform linear hypergraphs. If the infection probability
β is smaller than a critical value, the fraction
ρ of infected nodes at steady state vanishes. The structural parameter
m 1 is equal to
k−1. Note that the structures corresponding to
m 1=1 are simply random trees.

Figure reprinted from Ref. [107], © 2018 with permission from Elsevier.

#### 3.1.3. Phase transitions

A related model, which also features interesting phase transition behavior, is that proposed in Ref. [23], where the Authors considered a SIS process on simplicial complexes. As discussed in the previous Chapter, simplicial complexes are a special type of hypergraphs in which the presence of a hyperedge of size
k implies the existence of all possible hyperedges of all sizes smaller than
k on the same nodes. So, for example, if the hyperedge of size 3
(v 1,v 2,v 3) is present, then also the edges
(v 1,v 2),
(v 1,v 3) and
(v 2,v 3) are part of the hypergraph. From a functional point of view, this means that if a certain
k-body interaction exists between
k elements of the system, then all other possible lower-order interactions between the same elements also exist, all the way down to pairwise ones. Note that, structurally, this implies that simple simplicial complexes are linear hypergraphs. The reason is that two edges of the same size can have at most 1 node in common, because if they had more than one, this would result in multiple edges at the lower orders of interaction. For instance, if the edges
(v 1,v 2,v 3) and
(v 2,v 3,v 4) existed, the edge
(v 2,v 3) would occur twice, once induced by
(v 1,v 2,v 3) and once induced by
(v 2,v 3,v 4). Thus, simple simplicial complexes form a subset of linear hypergraphs. On these particular structures, the Authors considered an infection process in which a susceptible node can become infected if it participates in a hyperedge where all other nodes are already infected. The infectivity rates
β k are in general different between hyperedge sizes, whereas the rate of recovery
μ is constant for all nodes. Indicating with
d k the average degree of the nodes computed on hyperedges of size
k, an endemic steady state appears for high enough values of the 2-body rescaled infectivity
λ 2=d 2 β 2/μ. The main result is that higher-order infectivities act as control parameters for the nature of this phase transition. In fact, even if one only considers up to 3-body interactions, a larger
β 3 first increases the density of infected nodes in the endemic state, and then it changes the order of the transition itself, making it switch from continuous to first-order, as it can be seen in Fig. 10, where the steady state fraction of infected nodes
ρ∗ as a function of
λ 2 is plotted for different values of
λ 3=〈d 3〉β 3/μ.
To explain these results, the Authors developed a mean-field theory, writing the evolution equation for the fraction
ρ of infected nodes as (12)
ρ̇(t)=−μ ρ(t)+∑k=1 K β k d k ρ k(t)[1−ρ(t)],where
K is the largest size of hyperedges present in the simplicial complex. Numerical solutions of this equation reproduce the model simulations rather closely, providing a way to quickly estimate the critical behavior of the model and its steady state.

1. Download: Download high-res image (128KB)
2. Download: Download full-size image

Fig. 10. Switch in phase transition order in the SIS model on simplicial complexes. The results for
λ 3≡d 3 β 3 μ=0 (blue circles) are equivalent to a traditional SIS model on graphs. As
λ 3 is increased to 0.8 (black triangles), the fraction of infected nodes in the endemic state becomes larger and the curve becomes steeper, but with an unchanged critical point. Finally, as
λ 3 is brought to 2 (orange squares), the transition becomes first-order-like and a hysteresis cycle appears.

Figure adapted from Ref. [23].
The requirement of the structure being a simplicial complex was removed by the Authors of Ref. [108], who effectively extended the previous model to the case of scale-free
k-uniform hypergraphs, in which the probability that a random node has degree
d is proportional to
d−γ. To describe the behavior of the SIS model on such structures, they employed a heterogeneous mean-field theory, starting by expressing the evolution of the fraction of infected nodes with degree
d,
ρ d, as (13)
ρ̇d(t)=−μ ρ d(t)+β d Θ k−1[1−ρ(t)],where
β=β k, and
Θ is the joint probability that a node in a randomly chosen edge has degree
d and is infected. The similarity between Eqs. (12), (13) is evident, with the main difference being the absence of the sum over hyperedge sizes in Eq. (13) which, in turn, is due to the fact that once the requirement of using simplicial complexes is removed one does not have to consider all lower-order interactions together with a
k-body one. Solving Eq. (13) by a self-consistent approach, one finds the existence of a critical value of the power-law exponent of the degree distribution
γ c=2+(k−1)−1, which yields three distinct regimes. For
γ<γ c, the epidemic threshold in the control parameter
η=β μ disappears,
η=0 is effectively a critical point, and an endemic steady state always exists for any
η>0. When
γ=γ c, the theory predicts the presence of an epidemic threshold at
η=1 min i d(i)(k−1 π sin(π k−1))k−1.Finally, for
γ>γ c, one finds a hybrid transition with increasing values of
η c.
A further model sharing similarities with that of Ref. [23] is presented in Ref. [109]. The main difference is that for a node to become infected in the model of Ref. [109], one does not need it to participate in a hyperedge where all other nodes are already infected. Rather, contagion can occur if the number of neighbors already infected is greater than an arbitrary threshold. Qualitatively, the behavior of the model is fairly similar to that of Ref. [108], in that it presents both continuous and discontinuous phase transitions. However, the specific values of quantities of interest can vary significantly with the details of the structure considered. To illustrate this point, the Authors studied two extreme cases, namely that of a random regular graph with the addition of a hyperedge containing all the nodes, which they term a “hyperblob”, and that of a star graph with the addition of a hyperedge containing all the nodes, which they call a “hyperstar”. These cases are analytically solvable, and show that the dynamics has important differences between them, such as the hyperstar admitting a vanishing critical infection rate for the lower branch of the solution in the first-order transition regime. Also, the model provides a possible explanation of the large variability observed in real-world studies that reported substantial differences in the size of the critical mass needed to move a system out of an equilibrium state. In fact, these differences could be explained as being simply an expression of another key phenomenological feature of the model, which is the presence of bistability. This coexistence of states that appears for particular choices of the parameters remains also when the model is extended to include a recovered state from where one can become susceptible again. In fact, Wang and collaborators studied a type of Susceptible–Infected–Recovered–Susceptible (SIRS) model on simplicial complexes, and found that the phenomenology is effectively the same as that of SIS models [110]. The main difference, due to the presence of the extra state, is the appearance of limit cycles, which reproduce periodic epidemic waves as often observed in actual outbreaks (see Fig. 11).
A theoretical approach to study the rich phenomenology of this latter model is offered in Ref. [111]. There, the Authors extended the microscopic Markov Chain approach and the epidemic link equations to simplicial complexes with up to 3-body interactions. The resulting equations, made treatable via moment-closure approximations, reproduce some of the already mentioned phenomenology, including the appearance of first-order transitions (as seen in Fig. 12).
However, the nature itself of the approach makes it unsuitable to precisely compute the transition points. To circumvent this problem, the Authors of Ref. [112] introduced a new formalism, based on decomposing a simplicial complex into a set of edge-disjoint cliques, on which the microscopic epidemic equations were expressed. The use of this approach shows that the dependence of the critical point from the strength of the higher-order interactions is mediated by the higher-order dynamical correlations, making it not amenable to be estimated via mean-field methods.
An alternative approach is offered in Ref. [113], where the Authors focused on studying the transitions that simplicial epidemic models have between disease-free, bistable and endemic states. Using linear dynamical systems theory, they established rigorous conditions for the asymptotic behavior of the models and showed that bistability is a direct consequence of the presence of higher-order interactions.
A more general theoretical treatment is presented in Ref. [114], in terms of a mean-field theory applicable to actual hypergraphs, and not only simplicial complexes. The Authors illustrated the applicability of their method considering the two specific limiting cases of collective contagion, in which a node can become infected only if all the other participants in a hyperedge are infected (as in Ref. [23]), and individual contagion, in which a node only needs one other infected node in a shared hyperedge to have the possibility of becoming infected (as in Ref. [109] with a threshold of 1). Their approach can also work in intermediate cases, thereby fully accounting for the model of Ref. [109]. The main result is an estimate of the competition between contagions happening on hyperedges of different order. In fact, even when considering only pairwise and 3-body interactions, they found that having a heterogeneous distribution of the pairwise degrees suppresses the occurrence of first-order transitions, whereas a heterogeneous distribution of higher-order degrees promotes them. Moreover, the method lends itself well to explaining particularly complex cases such as the “hipster effect”, which causes trends to be rejected via group dynamics but accepted via pairwise peer pressure.
Further generalizations were studied by Higham and de Kergorlay, who considered both the deterministic and the stochastic SIS model on hypergraphs with nonlinear infection rates [115]. In this work, they established a spectral condition for the disease-free state to be asymptotically stable. Specifically, the condition is (14)
Λ(W)c β μ<1,where
β is the strength of the infectious interactions,
μ is the recovery rate and
c is a constant that is determined by the infection function used. Also,
W=I I T is the frequency matrix of relations, where
I is the incidence matrix of the hypergraph, and
Λ(W) is its largest eigenvalue. Note that each element
W i,j is just the number of hyperedges in which nodes
i and
j appear together. This means that the condition in Eq. (14) does not depend explicitly on the specific size of the edges that link node pairs. Moreover, since
W is a matrix, the interactions that are accounted for are only pairwise ones, even though they are implicitly weighted by the possible belonging of pairs of nodes to edges of different size. Also, if one considers a simple graph with adjacency matrix
A as the structure that supports the epidemic, the condition reduces to the well-known expression
Λ(A)β μ<1.Thus, it can be said that Ref. [115] provides the first mathematically rigorous proof of the extinction condition for SIS epidemics on simple graphs.

1. Download: Download high-res image (144KB)
2. Download: Download full-size image

Fig. 11. Appearance of epidemic waves in the simplicial SIRS model. A quenched mean-field treatment predicts the recurrence in time of a substantial number of infections in a simplicial complex with up to 3-body interactions. For this figure, the nodes have average pairwise degree 6 and each form on average 6.645 3-body interactions, with infection probability 0.5, recovery probability 0.9, probability of becoming susceptible again of 0.01 and 3-body enhancement factor of 88.

Figure reprinted from Ref. [110], © 2021 with permission from AIP Publishing.

1. Download: Download high-res image (230KB)
2. Download: Download full-size image

Fig. 12. First-order transitions in simplicial contagion. In a simplicial complex where the nodes have average pairwise degree 12 and each form on average 5 3-body interactions, with recovery probability 0.2 and 3-body infection probability 0.1, the fraction of infected nodes at steady state features a hysteresis cycle controlled by the pairwise infection probability
β.

Figure reprinted from Ref. [111].
Along the same lines, St-Onge and collaborators explored different means to explicitly account for nonlinearity in the transmission of a disease. In Ref. [116], they introduced an explicit nonlinear integral kernel to express the probability that a node becomes infected if it belongs to an environment of a given size in which a specific fraction of neighbors are already infected. Their application of the nonlinear kernel to hypergraph structures exploits a degree-based mean-field approach, and its results reproduce some of the features already observed in other models. In particular, first-order transitions and bistable regimes appear when the dependence of the infection probability on the local infected fraction is superlinear. Also, in the same regimes, the epidemic can grow super-exponentially if the global fraction of infected nodes becomes large enough.
In Ref. [117], the Authors instead adopted a different strategy, developing approximate master equations based on mean-field infection rates. With this formalism, they showed that the third moment of the hyperdegree distribution controls the possibility of having a discontinuous transition. In fact, they found that a first-order transition is only possible, in the thermodynamic limit, if the third moment of the distribution is finite. So, for example, if the degrees are power-law-distributed, first-order transitions can only occur when the power-law exponent is greater than 4. The behavior of the epidemic is also affected by the heterogeneity of edge sizes. In this case, heterogeneous edge-size distributions cause the disease to become endemic within edges of the largest sizes, which then act as reservoirs of infection and, in turn, inhibit bistability.
Whenever one studies epidemics, a question that arises naturally is how best to stop or slow down the spread of the disease. When the structure considered is that of a hypergraph, Jhun introduced a new strategy and extended existing ones based on simple graphs [118]. More specifically, they introduced a measure of the importance of each edge for the propagation of the disease, which they called “simultaneous infection probability” (SIP). They estimated the SIP via an individual-based mean-field theory, and showed that vaccinating the nodes within an edge in decreasing order of SIP is a highly effective strategy that can achieve low thresholds for herd immunity. Also, they generalized the edge eigenscore of graphs to
k-uniform hypergraphs, starting from the H-eigenvectors, which are defined as the vectors
v whose elements
v i 1 obey the equation
∑i 2,i 3,…,i k=1 N A i 1,i 2,…,i k v i 2 v i 3⋯v i k=λ v i 1 k−1,where
A is the adjacency tensor of the hypergraph. Then, given a hyperedge
(i 1,i 2,…,i k), they defined its H-eigenscore as the product
v i 1∗v i 2∗⋯v i k∗, where
v∗ is the H-eigenvector corresponding to the largest eigenvalue
λ. An edge vaccination strategy based on the H-eigenscore is also very effective, even though, of course, it can only be applied to uniform hypergraphs.

#### 3.1.4. Time varying and adaptive hypernetworks

All the models discussed so far consider networks whose structures are fixed in time. However, in real-world systems, links can often be created or destroyed, leading to networks that, over time, become significantly different from their past states [119]. To capture the effect of temporal variability on epidemic spreading, the Authors of Ref. [120] extended the model of Ref. [23] studying simplicial complexes with up to 3-body interactions, as in Ref. [111]. The main finding is that if the structure of the complexes is randomized at every time step, so that both the pairwise edges and the 3-body ones are placed according to an Erdős–Rényi-like process, then the early onset of epidemics disappears, and the only control parameter remains the pairwise infectivity, regardless of the initial infection seed and of the 3-body infectivity. Similarly, tuning the temporal correlations between explored structures, they found that the simplicial effects are stronger for lower correlations and weaker for higher ones.
Whilst the spreading of diseases is a typical phenomenon that is modeled on networks with pairwise and higher-order interactions, it remains a specific instance of a possibly more general group of processes, in which the state of a node diffuses to its neighbors. A first study of this kind is due to Wang and collaborators, who introduced a model of information diffusion on an adaptive hypergraph, with several parameters designed to mimic group interactions in a collaborative environment [121]. Later, Peng and co-Authors developed a statistical framework to explicitly model the same kind of environment, using a hypergraph representation of multilayer networks [122]. In general, such diffusive processes depend on some difference between the states of the nodes involved. As such, their mathematical formulation involves Laplacian operators. It becomes then very important to find generalizations of the traditional graph Laplacian to higher-order structures, and to study how they affect diffusive dynamics on hypergraphs. For the case of simplicial complexes, one possible generalization is based on their incidence matrices. Specifically, if
I[k] is the
kth-order incidence matrix, whose rows correspond to simplices of order
k−1 and whose columns correspond to simplices of order
k, then the
kth-order Laplacian is [123] (15)
L[k]=I[k]T I[k]+I[k+1]I[k+1]T.With this formalism, Torres and Bianconi [124] investigated higher-order diffusion on the so-called “network geometry with flavor” (NGF) model of simplicial complexes [125], [126]. Studying the properties of the simplicial Laplacians, they found that the dynamics is strongly affected by the spectral dimension of the second term in the right-hand side of Eq. (15). In fact, the estimate of the return-time probability for the dynamics changes significantly with the order of the interactions considered, converging to a definite value for dyadic interactions, but decaying as a power-law for higher-order ones (see Fig. 13).
Generalization of such treatments to hypergraphs are still in their infancy, mainly due to the inherent difficulty of defining a Laplacian operator on hypergraphs. The principal works in this direction are due to Chen and collaborators [127], [128], [129] and to Jost and Mulas [130], who offered the first rigorous and methodic study of a hypergraph Laplacian. Nonetheless, particular features of dynamical processes, such as their linear stability, can be assessed in the vicinity of fixed points of the dynamics without necessarily having to resort to a full tensorial formalism, as shown by Arruda and collaborators, who developed a framework to carry out this type of analyses and demonstrated its use in studying social contagion and diffusion on hypergraphs [131].

1. Download: Download high-res image (140KB)
2. Download: Download full-size image

Fig. 13. The long-term behavior of the return-time probability is affected by the order of the interactions. The blue line corresponds to 0-simplices, the orange dashed line to 1-simplices, the yellow dotted line to 2-simplices and the purple dot-and-dash line to 3-simplices. The data shown are for NGF simplicial complexes of 2000 nodes, with maximum dimension 3 and flavor parameter
−1.

Figure adapted from Ref. [124].

### 3.2. Social contagion

The notion of contagion went beyond its original meaning in epidemics to describe a vast range of processes that spread among individuals, including the diffusion of political opinions, financial decision, and the adoptions of new products and technologies [132]. These phenomena of social contagion are rather important, and in several cases they cannot be described in terms of pairwise interactions. For instance, there are several ways in which individuals can be convinced to adopt a given product. Individuals may be persuaded by means of successive pairwise interactions with their neighbors. However, an alternative and conceptually different contagion mechanism is peer-pressure: the individual might be part of a social group with more than two members, where the other group members are all adopters. In this case, contagion takes place through multi-body interactions, which need to be represented by a higher-order structure like a hypergraph or a simplicial complex, in place of the standard network representation that only considers pairwise interactions.

#### 3.2.1. Contagion in hypernetworks and simplicial complexes

Iacopini et al. [23] have proposed a model that combines stochastic processes of simple contagion, which could be expressed through pairwise interactions, with the complex contagion dynamics taking place through group interactions. The simplicial complex describing both types of interactions is an extensions of the classic SIS model. Each of the
N vertices of the simplicial complex has a binary state
x i(t), with 0 and 1 corresponding to the
ith individual being susceptible or infectious, respectively. The Authors associate an order parameter
ρ(t) to the system, which corresponds to the fraction of infectious nodes at time
t.
The process of contagion happens over the simplex
σ of dimension
D, and is defined by a set of control parameters
β 1,…,β D. Note that
β 1 corresponds to the standard infection probability in pairwise interactions, which corresponds to the probability per unit time that a susceptible node
i gets the infection from an infectious neighbor
j, that is, through the link
(i,j). In the same fashion,
β 2 will be the probability per unit time that node
i receives the infection from a “full” triangle, that is, a 2-simplex
(i,j,k) in which both
j and
k are infectious, and so on for
β 3,…,β D. In general, given a susceptible node that is part of a simplex of
n+1 nodes where all other nodes are infectious, it will become infectious with probability per unit time equal to
β n. The recovery dynamics are instead driven by the (node-independent) recovery probability
μ.
The Authors simulated this higher-order SIS contagion dynamics on data sets of face-to-face interactions from the SocioPatterns collaboration [133], where a simplex of dimension 2 was reconstructed. It was observed that, for certain selections of the infectivity rates, a hysteretic behavior may appear and, at variance with the standard SIS model, an endemic state may coexist with a healthy one. In this region of parameters, the final, emergent state will depend on the initial density of infectious nodes. This emerging phenomenon, confirmed in simulations on synthetic simplicial complexes, is explained through a mean field description of the model, under the assumption of a homogeneous mixing [104]. Namely, the evolution of the density of infected nodes
ρ(t) is shown to follow this general equation: (16)
ρ̇(t)=−μ ρ(t)+∑d=1 D β d〈k d〉ρ d(t)(1−ρ(t)),where
〈k d〉 is the average number of
d-dimensional simplices incident at a node of the simplicial complex. Using this approach, the Authors were able to correctly predict the nature of the observed transitions.
When the assumption of homogeneous mixing is not tenable, alternative approaches have to be employed in order to make accurate predictions. Matamalas et al. in Ref. [111] focused their analysis on the case of 1 and 2-simplices, and have shown that both the epidemic link equations (ELE) of Ref. [134] and the microscopic Markov chain approach of Ref. [135] (MMCA) yielded an improved predictions on the prevalence and critical point, with the first method providing the best accuracy. More specifically, the MMCA allows one to define the following set of
N discrete-time equations that describe the evolution of the probability
p i that a node is infected (17)
p i(t+1)=[1−p i(t)][1−q i,1(t)q i,2(t)]+p i(t)(1−μ),where
q i,1(t) and
q i,2(t) are the probability that node
i is not infected at time
t through any pairwise or three-wise interaction, respectively, that is, (18)
q i,1(t)=∏j∈E i 1(1−β 1 p j(t)),q i,2(t)=∏j,l∈E i 2(1−β 2 p j(t)p l(t)),with
E i d being the set of
d-simplices containing node
i. Eq. (17) is based on the mean field approximation that the probability that a node is infected by one neighbor is independent of the probability of being infected by any other neighbor. However, the presence of a three-body infection mechanism makes this assumption less tenable, and induced the Authors to adopt the ELE method to appropriately adjust the computations of
q i,1(t) and
q i,2(t). The benefits for using this method are illustrated in Fig. 14, which shows how the tipping point for
β 1 to an endemic state is better estimated through the ELE method.
As already mentioned in Section 3.1.2, Ref. [112] sought for a more general method that can be adopted to study social contagion on simplicial complexes of any order. Specifically, the Authors took into account the higher-order dynamical correlation that takes place in a
d-simplex, which induces the presence of
d+1 clique in the graphs. They observed that when two or more cliques were sharing
m≥2 nodes, then care should be put in computing the contribution to the infection of one of these
m nodes of the other
m−1 nodes, to avoid overcounting it. Toward this goal, Burgio et al. sought for the minimal set of edge-disjoint cliques (so that any pair of cliques only shared at most one node) cover of the underlying graph of the simplicial complex. However, since this problem is NP-hard, they proposed a heuristic to identify such a minimal cover set. Building on this heuristic, the Authors developed the Microscopic Epidemic Clique Equation (MECLE) model that accounts for the correlations among states of the nodes belonging to the same clique, yielding a more accurate prediction of both the prevalence and of the critical point, that is, the critical value
β 1 cr of the edge infection probability
β 1 at which the epidemic-free state become unstable.

1. Download: Download high-res image (147KB)
2. Download: Download full-size image

Fig. 14. Infected nodes as a function of
β 1. Stationary density of infected nodes as a function of
β 1 in the SIS model over a random 2-dimensional simplicial complex of
N=2000 nodes, with recovery probability
μ=0.2. The hypergraph has an average number of edges
〈k 1〉 and triangles
〈k 2〉 incident at a node equal to 12 and 5, respectively. Panel (a) corresponds to
β 2=0, whereas panel (b) to
β 2=0.1. Solid dots correspond to Monte Carlo simulations obtained through the quasi stationary approach, whereas solid, dashed, and dotted lines corresponds to analytical computations using ELE, MMCA, and the mean-field approximation, respectively.

Figure adapted from Ref. [111].

1. Download: Download high-res image (524KB)
2. Download: Download full-size image

Fig. 15. The SIS model over 2-dimensional simplicial complexes. Stationary density of infected nodes as a function of
β 1 in the SIS model over different 2-dimensional simplicial complexes, with recovery probability
μ=0.2. Dots correspond to the results of Monte Carlo simulations, whereas lines to analytical computations using the methods indicated in the legends, where (i) MECLE refers to the Microscopic Epidemic Clique Equation model, which is the one that best predict the critical edge infection probability
β 1 cr at which the epidemic-free state becomes unstable, (ii) MMCA and MMCA (MECLE) refer to the Microscopic Markov Chain approximation of the simplicial Epidemic Link Equations model and the Microscopic Epidemic Clique Equations (MECLE), respectively, (iii) ELE corresponds to the simplicial Epidemic Link Equations model, and (iv) MF and MF(MECLE) refer to the homogeneous mean-field approximation of the simplicial Epidemic Link Equations model and the Microscopic Epidemic Clique Equations (MECLE), respectively. Note that MF and MF(MECLE) are almost superimposed and practically indistinguishable. Each panel corresponds to a different topology and triangle infection probability: panel (a) to a periodic triangular SC with
k̄(0,1)=k̄(0,2)=0.00 and
k̄(1,2)=3.00, with
k̄(g,r) being the mean number of
(g,n+1)-cliques1 incident on a node, and triangle infection
β 2=0.25; panel (b) to a random (see Ref. [23]) SC with
k̄(0,1)=4.10,
k̄(0,2)=0.00 and
k̄(1,2)=3.95, and
β 2=0.15; panel (c) to a Dorogovtsev–Mendes (see Ref. [136]) SC with
k̄(0,1)=1.10,
k̄(0,2)=0.00 and
k̄(1,2)=1.45, and
β 2=0.25; and panel (d) to the same as panel (c) but with
β 2=0.50.

Figure adapted from Ref. [112].
To numerically demonstrate this finding, the Authors generated different 2-dimensional simplicial complexes from graphs. Namely, starting from a given graph model, they computed its edge-disjoint clique cover and converted 3-cliques in 2-faces. A such generated simplicial complex is identified by adding a ‘SC’ to the name of the graph model used for it. Fig. 15 shows that, for each of the four simplicial complexes, the Microscopic Epidemic Clique Equation (MECLE) model best approximates the prevalence and the critical value
β 1 cr with respect to alternative analytical methods.

1. Download: Download high-res image (782KB)
2. Download: Download full-size image

Fig. 16. Simplicial SIS spreading in real world networks. Simplicial SIS spreading on dolphins network made of 62 nodes and 159 edges, and such that a certain portion of its length 3 cycles are randomly selected as 2-simplices. Panel (a) describes the traditional SIS spreading process on graphs, where blue and red nodes correspond to susceptible and infected, respectively. Gray-colored triangles denote the higher-order connections, which do not affect the SIS spreading since the normalized infectious rate
λ 2=β 2〈k 2〉/μ is set to zero. Panel (b) describes the simplicial spreading process, where blue- and orange-colored triangles denote dormant and active 2-simplices, respectively. Higher-order interactions only happen in an active simplex, with a dormant
p-simplex becoming activate only when
p of its nodes have been infected.

Figure reprinted from Ref. [24], © 2021 with permission from Elsevier.
Li et al. in Ref. [24] focused their analysis on the different time-scales at which contagion happens in the presence of higher-order interaction. This can be explained focusing for instance on the case of pairwise and three-body contagion. At the onset of the spreading, the contagion mainly happens through pairwise interaction, and therefore is mainly influenced by the pairwise infection probability
β 1. When the prevalence is sufficient, higher order contagion, whose dynamics is governed by
β 2, is triggered, thereby accelerating the spreading as illustrated in the example reported in Fig. 16. The Authors also observed that, while the ignition of the contagion process requires the pairwise interactions, it can be sustained by higher order interactions alone, as shown in Fig. 17. This observation has an interesting interpretation in social contagion, whereby the spreading of an idea needs to be ignited by individual, pairwise interactions, but can be then sustained by peer-pressure alone. The Authors showed the generality of the observed phenomena by means of mean-field methods, firstly extended to higher order contagion [23].
Contagion processes can coexist and influence each other. Several instances can be found in disease spreading, whereby the presence of a certain pathogen may favor or hinder the diffusion of another one [137], [138]. However, interaction between spreading processes is not limited to infectious disease, but can also be observed in social contagion, a striking example being the coupling between the COVID-19 epidemic and the adoption of safe behaviors [139]. Motivated by these examples, Lucas et al. in Ref. [140] focused on the case of two interacting spreading processes, a simplicial contagion
A that cooperatively and unidirectionally drives a simple (pairwise) contagion
B. Fig. 18 is a schematic illustration of the interactions, where unidirectionality is guaranteed by the fact that the coupling factors are selected as
ϵ B A=1 and
ϵ A B>1. Introducing
ρ A,
ρ B, and
ρ A B as the density of individuals infected by only
A,
B, and by both
A and
B, respectively, and considering that
ρ A tot=ρ A+ρ A B and
ρ B tot=ρ B+ρ A B are the total density of individuals infected by
A and
B, respectively, the Authors derived the following mean-field equations: (19a)
ρ̇A tot(t)=ρ A tot(t)(β 1 A〈k 1〉(1−ρ A tot(t))+β 2 A〈k 2〉ρ A tot(t)(1−ρ A tot(t))−μ)/μ,(19b)
ρ̇B tot(t)=ρ B tot(t)(β 1 B〈k 1〉(1−ρ B tot(t))+β 1 B〈k 1〉(ϵ A B−1)(ρ A tot(t)−ρ A B(t))−μ)/μ,(19c)
ρ̇A B(t)=−2 ρ A B(t)+ϵ A B β 1 B〈k 1〉(ρ A tot(t)−ρ A B(t))ρ B tot(t)/μ+β 1 A〈k 1〉(ρ B tot(t)−ρ A B(t))ρ A tot(t)/μ+β 2 A〈k 2〉(ρ B tot(t)−ρ A B(t))ρ A tot 2(t)/μ, where
β i J is the probability that a node gets the infection
J from a
i-simplex it belongs to, and
ϵ A B>1 is the coefficient that regulates the interaction between the two diseases, being the transition probability in this latter case simply the product of those probabilities coming from the single infections
A and
B (see again the schematic in Fig. 18). Note that the Authors set
β i A=0 for all
i>2 (epidemic
A only diffuses through edges and triangles), and
β i B=0 for all
i>1 (epidemic
B, if decoupled from
A, would only diffuse through edges), and we remind that
ϵ B A does not appear in the equations since it is set to 1, whereby epidemic
B does not affect the diffusion of epidemic
A. Note that (i) the last choice makes the interaction between the spreading processes asymmetric, with the spreading of
A in (19a) independent of the diffusion of
B, and driving the other two Eqs. (19b), (19c), and (ii) equation (19a) - (19c) can be completed by the following conservation law, which allows to compute the density of susceptible (non-infected) individuals
ρ S as (20)
ρ S(t)=1−ρ A tot(t)−ρ B tot(t)+ρ A B(t).An interesting numerical and analytical finding from Lucas et al. was that when the normalized pairwise infection rate
λ 1 B=β 1 B〈k 1〉/μ is smaller than 1 (meaning that, in the absence of a driving process, the epidemic would die out), the driving process
A can sustain epidemic
B if
λ 2 A=β 2 A〈k 2〉/μ>1, and one can observe a discontinuous transition when the coupling parameter
ϵ A B overcomes a critical value, as illustrated in Fig. 19. These findings suggest that, when a simplicial social contagion process, such as the adoption of safe behaviors, drives a simple (pairwise) contagion of an epidemic, acting on the social contagion might yield qualitatively different outcome on the final pathogen’s spread.

1. Download: Download high-res image (131KB)
2. Download: Download full-size image

Fig. 17. 2-simplices may sustain infection also in the absence of pairwise interactions. Infection density on a synthetic network of 1000 nodes, 4140 1-simplices, and 1401 2-simplices, generated through the extended Barabási–Albert model from Ref. [21]. Pairwise interactions are shut down in correspondence of the red dashed line, when the infection density settles. Different from panel (a), where 2-simplices are inactive (
λ 2=0), in all the other panels the higher-order interactions are able to sustain the contagion.

Figure reprinted from Ref. [24], © 2021 with permission from Elsevier.
Another bulk of work focuses instead on the case in which the two interacting spreading processes are competing with each other [141], [142], [143]. This is the case, for instance, of two concurrent epidemics such that if an individual is infected by one disease cannot be infected by the other, or in the adoption of two competing, alternative products (e.g., two competing brands of top-of-the-line smartphones). Li et al. in Ref. [141] considered then a simplicial contagion for two interacting SIS spreading processes
A and
B under absolute competition, whereby being infected by
A implies immunity to
B and vice versa. This implies that an individual can be in three states: susceptible, infected by
A, and infected by
B. In the absence of simplicial contagion, this system would have, depending on the recovery probability and infectious rates, three possible steady-state behaviors: the dominance of one of the two epidemics, or the extinction of both epidemics. The presence of simplicial interactions determines instead a novel region of the parameter space where the dominance of one epidemic or the other becomes dependent on the initial conditions. Similar outcomes have been independently observed through the MMCA by Nie et al. [142], who numerically investigated the impact of the topology of the simplicial complex, showing that an increase in the average degree of 1-simplices can trigger a transition from epidemic outbreak to extinction, while an increase in the average degree of 2-simplices favors epidemic spreading.

1. Download: Download high-res image (242KB)
2. Download: Download full-size image

Fig. 18. Schematic illustration of the interacting spreading process studied in Ref. [140]. Panel (a) depicts the transition probabilities between the four compartments: susceptible (
S, gray), infected exclusively by one disease (
A or
B, blue/red-colored, respectively) or by both (
A B, black). Panels (b)–(d) show how a susceptible node
i can get A after a contact with an infectious
k-simplex (this could also includes
A B individuals). In particular, in panel (d), as
i is part of a 2-simplex composed by two other infectious nodes, it could acquire the infection both from each of the two links (1-faces) with probability
β 1 A=β A and from the 2-face with probability
β 2 A=β A Δ. Panels (e)–(g) show that when
i is already infected with
B, the probability of getting
A is affected by the coupling term
ϵ B A. The same would happen, symmetrically, by swapping B and A.
A is selected as the driving process on B by setting
ϵ A B>1 and
ϵ B A=1.

Figure reprinted from Ref. [140].

1. Download: Download high-res image (276KB)
2. Download: Download full-size image

Fig. 19. The stationary-state obtained by numerical integration of Eq. (19). When the driver process
A follows a simple contagion (
λ 2 A=0) the transition to the epidemic state is continuous, see panels (a)–(b). On the contrary, when
A follows a simplicial contagion (
λ 2 A=2.5) a discontinuous transition may occur, depending on the values of the coupling
ϵ A B. The blue dot in panel (d) identifies the critical point corresponding to the discontinuous transition displayed in panel (c).

Figure adapted from Ref. [140].
A different contagion mechanism was considered in Ref. [143], where Nie et al. grouped the agents in three categories of adopters:
Ω A (
Ω B) is the set of agents that only adopt product
A (
B), and possibly do so only if they are connected to other adopters of
A (
B) through links or triangles; the adopters in the third category
Ω A B, when susceptible, may become adopters of either product
A or
B. In addition, when they adopted one of the two product, beyond recovering (that is, adopting neither
A or
B) they can also switch adopted product: if they are adopters of
A (
B) and are in a triangle with two adopters of
B (
A), then they may decide to adopt product
B (
A), instead. A schematic of this contagion mechanism is illustrated in Fig. 20. To identify the critical transition points in the infection rates, the Authors used the MMCA as in Refs. [111], [135] and then numerically confirmed the theoretical results in the symmetric scenario in which there is not an a priori favorite product, that is, when the recovery and infection rates are the same from both
A and
B, the fraction of adopters in
Ω A and
Ω B are the same, as well as the probability to switch from product
A to product
B and vice versa. As in Ref. [141], here the Authors show that there cannot be a steady-state where both products are adopted, but one product will dominate on the other, depending on the contagion parameters and initial conditions.

1. Download: Download high-res image (231KB)
2. Download: Download full-size image

Fig. 20. The contagion mechanism described in Ref. [143], with transition rates attached to the corresponding arrows. Panel (a) describes agents belonging to
Ω A, which may only be susceptible or adopt product
A. Agents in
Ω B behave similarly, but may only adopt product
B, as illustrated in panel (b). Finally, agents in
Ω A B may become adopters of both
A and
B. In all panels, the susceptible state is identified by a blue node, yellow (green) nodes corresponds to the compartment of agents adopting product
A (
B). A gray node correspond to a product that cannot be adopted by the agents in that category.

Figure reprinted from Ref. [143], © 2022 with permission from Elsevier.

#### 3.2.2. Contagion in time-varying structures

Depending on the time scale at which we are studying the contagion process, we may observe that, as time passes, some connections are created, whereas others are eliminated from the interaction networks. The framework of temporal graphs, that is, an ordered sequence of static graphs [144], has been typically used to describe networks whose structures evolves in time. Chowdhari et al. [120] investigated how the variability of the interaction network affects the social contagion by considering a temporal simplicial complex, where hyperedges may be added and removed during the spreading (see Fig. 21). In their study, the Authors focused on networks where contagion can either happen on one-simplices (standard, pairwise interactions, with probability
β 1) and two-simplices (three-wise interactions, with probability
β 2). To perform a numerical investigation of the impact of temporality on the spreading dynamics, the Authors leveraged the extension to simplicial complexes of the microscopic Markov Chain approach [111], with the difference standing in the sets
E i d in Eq. (18) being functions of time.
To compare social contagion through static and temporal interactions, the Authors considered random simplicial complexes (RSC) of
N=500 nodes generated with the algorithm given in Ref. [23], so to control for the average pairwise degree
〈k 1〉 and the average number of two-simplices incident on a node
〈k 2〉. Eq. (18) is then run for
T=10,000 time-steps, with the initial state corresponding to a single infected node (
ρ(0)=1/N), so to compute the asymptotic density of infected nodes as
ρ∞=∑i p i(T)/N. By comparing simulations on a static simplicial complex with those on a temporal hypergraph where at each time step they considered a new realization of the RSC with same
〈k 1〉 and
〈k 2〉 of the static case, the Authors observed a critical difference. While in the static case the epidemic onset is anticipated by increasing the rescaled 2-simplices infection parameter
λ 2=β 2〈k 2〉/μ, for temporal simplicial complexes for low values of the rescaled pairwise infection parameter
λ 1=β 1〈k 1〉/μ the endemic phase is not reached independent of the value of
λ 2, as illustrated in Fig. 22. This effect is confirmed also when multiple nodes were infected at time 0. The backward transition to the infection-free state is instead almost unchanged by the temporality of the interactions.

1. Download: Download high-res image (184KB)
2. Download: Download full-size image

Fig. 21. Contagion in time-varying structures. Schematic illustration of the temporal higher-order interactions considered in Ref. [120].

Figure reprinted from Ref. [120].

1. Download: Download high-res image (138KB)
2. Download: Download full-size image

Fig. 22. The steady-state infection rate. Color map of the steady-state infection rate as a function of the normalized infection rates
λ 1 and
λ 2 for a static (left panel) and temporal (right panel) random simplicial complex of
N=500 nodes. The recovery rate
μ=0.1, and average degrees
〈k 1〉=12 and
〈k 2〉=5 are the same for both scenarios. In the static case, the epidemic onset, identified by a solid black line, is anticipated as we increase
λ 2. The same effect is not observed on temporal simplicial complexes. The backward transition to the infection-free state (dashed black lines) is instead largely unaffected by the temporality of the interactions.

Figure reprinted from Ref. [120].
The above numerical experiment considered independent realization of the RSC at any time-step. However, when the network evolution is less abrupt, a temporal correlation between the sequence of simplicial complexes may be observed. To this aim, Chowdhari et al. proposed an algorithm where a RSC is generated only at time 0, and then its simplices are rewired with a given probability
f, which can be used to tune the temporal correlation
σ between 0 and 1. This experiment showed how the simplicial effects are more evident in the presence of high temporal correlation, with the static case that can be viewed as the limit case with unitary correlation. Nonetheless, in case of heterogeneous network structure such as the scale-free simplicial complexes [21], even for uncorrelated temporal structures (
σ≈0) the effect of simplicial interaction on the epidemic onset is not negligible, since for any value of
β 1 the endemic state can be triggered by increasing
β 2.
The effect of temporality on social contagion has also been tested in the presence of two competing spreading processes
A and
B [145], where the nodes can be in four possible states: susceptible to both infections, infected by
A, infected by
B, or infected by both
A and
B, whereas contagion can happen both via pairwise interaction and via 2-simplices. Nie et al. used the MMCA in Ref. [135] to identify the steady-state outcome of the epidemic spreading, and its effectiveness has been confirmed through numerical experiments. When the initial infection density is sufficiently high, the Authors observed that a network that evolves at a time scale that is slower than that of the spreading dynamics induces a smaller outbreak threshold.

#### 3.2.3. Other higher order contagion mechanisms

##### SIR simplicial contagion.

As in contagion on graphs, SIR models consider the case when, after healing from an infection, individuals get immune from the contagion, and therefore are classified as recovered [102]. As in SIS simplicial contagion, Palafox-Castillo and Berrones-Santos [146] considered that the infection can also arise from the interaction between a susceptible individual and its
k infected neighbors belonging to the same
k-simplex. In their preliminary work, the Authors showed how even the addition of only three-way interactions can yield notable differences in the dynamics and stationary state of the epidemic. For instance, using a mean field approach, they observed that the additional higher-order layer may prevent the possibility of containing the number of infected nodes below a certain threshold.

##### SIRS simplicial contagion.

The work of Iacopini et al. [23] has been then extended in a different direction by Wang et al. [110]. Here, the Authors considered the presence of a further compartment, the re-infectable, and study simplicial contagion in the context of SIRS models [147], widely used for studying epidemic spreading of influenza, hepatitis, and many other viruses [148], [149]. The presence of this additional compartment allowed to capture the periodicity that often characterizes epidemic outbreaks. As in Ref. [23], the Authors demonstrated the observed results through a mean field approach, and showed the existence of a discontinuous transition to the epidemic state, as well as that of limit cycles explaining the periodicity of the epidemic outbreaks.
As already mentioned several times, a simplicial complex is a class of hypergraphs with an additional property: if a hyperedge is in a simplicial complex, any non-empty subset of vertices in the hyperedge is also a hyperedge of the simplicial complex. When dealing with social contagion, this means that a node can be infected not only if all other nodes in the
d-simplex are infected, but also if any smaller (non-empty) subset of its
d neighbors in the simplex is infected. Jhun et al. [108] focused instead on
(d+1)-uniform hypergraphs, that is, on hypergraphs composed only by hyperedges of cardinality
d+1, that is, connecting
d+1 nodes. The scope of the Authors was to focus on the case in which contagion of a node
i is possible only when all
d nodes belonging to the same hyperedge as
i are infected. The contagion model from Ref. [108] could be retrieved by setting all
β i-s to 0, except of
β d. For this reason, we will drop the subscript from now on and write
β d=β. Jhun et al. focused on scale-free
(d+1)-uniform hypergraphs, which therefore have a power-law
d degree distribution, and designed it so that its exponent
λ is equal to
1+1/μ. Accordingly, they relaxed the assumption of homogeneous node degree in standard mean field analysis, and rather used the heterogeneous mean field approach first introduced in Ref. [150]. Generalizing it to the case of uniform hypergraphs, the Authors obtained the following equation for the density
ρ k d of infected nodes with a given degree
k d: (21)
ρ̇k d(t)=−μ ρ k d(t)+β(1−ρ k d(t))k d Θ d−1,where (22)
Θ=∑k=k min∞k d P d(k)ρ k(t)〈k d〉,with
k min and
P d being the smallest
d-degree and the
d-degree distribution of the hypergraph, respectively. The stationary solution for
ρ k will be (23)
ρ k d=η k d Θ d−1 1+η k d Θ d−1,where
η=β/μ. Through algebraic manipulations, the Authors were able to show that, in agreement with numerical simulations, there is a critical value
λ c=2+1/d of the exponent
λ of the scale-free degree distribution. Indeed, when
λ≤λ c, the epidemic threshold is zero, whereas it becomes finite for
λ>λ c, with the susceptibility, measuring the fluctuation of the order parameter, diverging at the transition point, as illustrated in Fig. 23.
A different contagion method, framed in the context of social communication, was proposed by Wang et al. in Ref. [151]. The Authors considered a simplicial complex built starting from a graph, where each of the minimal cliques of the original graph corresponds to a simplex [152]. In such a simplicial complex, the Authors considered a SIS propagation model, where contagion may happen both at the level of connected individual (called by the Authors “information transmission”), and between cliques sharing faces (called by the Authors “idea integration”). These two contagion mechanisms are mutually coupled, whereby a clique is infected if the ratio of its infected nodes is not smaller than a threshold
T, whereas the change of the node state depends both on the communication process intra-clique and inter-clique. An interest finding is that in this contagion scheme information does not spread when the communication between cliques vanishes, that is, in the absence of idea integration.

1. Download: Download high-res image (458KB)
2. Download: Download full-size image

Fig. 23. Heterogeneous mean-field approach. Density of infection
ρ (panels (a) and (c)) and susceptibility
χ 1 (panels (b) and (d)) as a function of
η estimated through a heterogeneous mean-field approach in Ref. [108]. Panels (a) and (b) correspond to 3-uniform scale-free simplicial complexes, whereas panels (c) and (d) correspond to 4-uniform scale-free simplicial complexes. In all panels, results are reported for different values of the exponent
λ of the scale-free degree distribution. It can be noted that, beyond a certain critical value of
λ the epidemic threshold becomes finite, and the susceptibility diverges at the transition point. Figure reprinted from Ref. [108]. © IOP Publishing Ltd and SISSA Medialab srl. Reproduced by permission of IOP Publishing. All rights reserved.

### 3.3. Games on hypernetworks and simplicial structures

The idea of the state of a node propagating to its neighbors gave rise not only to studies of epidemic spreading and diffusive processes, but also to games and to a family of models that concern the adoption of different opinions within a collection of agents. A prototypical case is simply known as the voter model, which, in its simplest instance, consists of a network of nodes whose opinion changes at every time step by copying that of a randomly chosen neighbor. The voter model and other related games on networks enjoyed a long history of interest; it is therefore only natural that they are currently being investigated on higher-order structures.

#### 3.3.1. Majority rule and voter model

A pioneering study is that of Lanchier and Neufer in Ref. [153]. Their motivation to use hypergraph was to rigorously analyze a spatial version of the majority rule model, in which groups of nodes are randomly chosen and they all adopt the opinion of the majority within their group. To do so, they mapped the dynamics, whose natural underlying structure is that of a hypergraph, into lattices, and adopted a tie-breaking advantage for even-sized edges in favor of one of the two states (see Fig. 24). Using a Markov Chain formalism, they showed that the steady state depends on the size of the edges. Specifically, if they involve an odd number of nodes, the system separates into clusters of differing opinions, whereas if the number of nodes is even, the advantaged opinion wins. Notably, this behavior differs from that of the voter model, which, in higher dimensions, always converges to a well-mixed steady state.
A different approach was taken by Chung and Tsiatas [154], who started from a fully general formulation of the voter-model dynamics on hypergraphs. In their formalism, the nodes in a hyperedge
g can adopt a given assignment of opinions
τ with a probability
p(g,σ,τ) that depends on their current opinion assignment
σ. Then, the dynamics of the model can be mapped to a random walk on the state graph
H∗ of the model. If the model has
r allowed opinions and
N nodes, then
H∗ is a directed graph with
r N nodes, each corresponding to a global state of the system, and its arcs link states that can be reached one from the other. To cope with the problem, they focused on the particular case in which the transition probabilities do not depend on the current opinion assignment of the edges, so that
p(g,σ,τ)=p(g,τ). Then, they treated the state transitions as members of a left-regular band semigroup and used spectral methods to express the probabilities of observing given global states and to estimate the time of convergence to a steady state.

1. Download: Download high-res image (248KB)
2. Download: Download full-size image

Fig. 24. Typical updating process of the voter model on a hypergraph. The dynamics is mapped on a 2-dimensional lattice, and two updates are shown, each with a block of 9 nodes.

Figure based on Ref. [153].
Gradowski and Krawiecki proposed instead a voter model sharing similarities to classic spin models [155]. In their work, they considered a two-state system on a scale-free hypergraph. Each node can have opinion
±1, and two types of dynamics were studied, namely the hyperedge update and the node update. In the former, all edges are considered, systematically and in a random order, and the nodes in each edge flip their state with a probability that depends on the majority opinion within the edge. In the latter, all nodes are chosen, in random order, together with a random edge they participate in, and their states are flipped with the same probability as before. Given a node
i with state
σ i and an edge
e to which it belongs, the flipping probability is
p(i)=1 2[1−(1−2 q)σ i∑j∈e σ j|∑j∈e σ j|],where
q is a parameter akin to a temperature that represents social noise. Simulating the system and analyzing the results, the Authors found that it features a critical transition controlled by the social noise
q (as it can be seen in Fig. 25), with critical exponents depending on the dynamics. They are compatible with mean-field Ising exponents for the hyperedge update, and they are close to those of the voter model on low-dimensional hypercubic lattices for the node update.
As already mentioned, often network structures are not static, but evolve in time with a dependence on the state of the nodes. To study one such case, Horstmeyer and Kühn considered an adaptive voter model on simplicial complexes [156]. Specifically, they simulated a version of the classic voter model on simplicial complexes including up to 3-body interactions. The rule they implemented consists in a repeated random selection of individual edges. If the chosen edge is between nodes with differing opinions (an active edge) then one of the two nodes is randomly selected. With probability
p, the chosen node rewires to another random node with the same opinion. If the edge is not part of a 2-simplex, then, with probability
1−p, the chosen node adopts the opinion of the other node. If, instead, the edge does belong to a 2-simplex, then, with probability
1−p, the chosen node adopts the majority opinion within the 2-simplex with probability
q, and it adopts the opinion of the other node in the edge with probability
1−q. The introduction of the probability
q allows one to consider the effect of peer pressure within the model. The behavior observed is qualitatively very similar to the classic version on networks, in that it features a fragmentation transition that, for values of
p above a critical point, separates a phase in which the network adopts a single opinion to one in which it polarizes into two blocks that no longer communicate with each other. The main difference with the traditional model is that the peer pressure
q speeds up the dynamics. Additionally, in the appearance of the fragmented state, the higher-order interactions vanish before the individual edges do, suggesting that such a dimensional hierarchy may be a general feature of dynamical models on higher-order networks. This is consistent with the theoretical results obtained by Neuhäuser and collaborators, who explicitly showed that true higher-order effects in models of this kind can only arise if the higher-order interaction function is nonlinear [157]. Conversely, if such a function is linear, then the dynamics can always be reduced to dyadic interactions, even if the supporting structure is that of a hypergraph.

1. Download: Download high-res image (177KB)
2. Download: Download full-size image

Fig. 25. Binder cumulants for the node-update dynamics of the voter model on hypergraphs. The curves for system sizes 1000, 2000, 5000, 10,000, 20,000, 50,000 and 100,000 cross at the critical value of the social noise, estimated to be approximately 0.152.

Figure reprinted from Ref. [155], © 2015 with permission from Acta Physica Polonica.

#### 3.3.2. Evolutionary games

A related model that has recently been studied in a hypergraph version is the public good game, which consists of repeated rounds in which players interact in groups. In each round, all players are asked to contribute a fixed amount
c to a common fund. Players can choose whether to cooperate, and pay the requested amount, or defect, and offer no contribution. The total amount collected is then multiplied by a synergy factor
R and shared equally amongst the
g participants, regardless of their choice. Alvarez and collaborators studied a version of this game in which the interactions happen within hyperedges in a hypergraph [158]. Their results show that the dynamics of the simulated game reproduces the theoretical expectations of the replicator dynamics if structural correlations are absent. Conversely, in heterogeneous hypergraphs, hierarchical structures are shown to inhibit cooperation, which is an effect that is not observed in dyadic networks.
Schlager and collaborators also considered a similar game, namely the snowdrift game, on simplicial complexes [159]. In its traditional, two-player version, each player may contribute an equal amount to cover a total cost
c for a shared benefit
b, which is awarded to both players if at least one of them cooperates. They generalized the game to simplicial complexes via a direct extension of the payoff functions, so that, if the number of cooperators is
n C⩾1, each cooperator receives a payoff equal to
b−c/n C, and each defector receives simply
b. The complexes include up to 2-simplices, and the dynamics is adaptive and similar to that of Ref. [156]. At every iteration of the game, an edge is chosen and one of its two nodes is randomly selected. If the edge is not part of a 2-simplex, then with probability
φ the chosen node rewires to another random cooperator node, and with probability
1−φ it adopts the strategy opposite to that of the other node in the edge. If, instead, the edge is part of a 2-simplex, then with probability
ρ the chosen node adopts a defecting strategy if at least one of the other two nodes in the simplex cooperates, and it adopts a cooperating strategy otherwise. Finally, if the edge is part of a simplex, with probability
1−ρ the chosen node rewires to a cooperator with probability
φ, and it adopts the strategy opposite to that of the other node in the edge with probability
1−φ. To analyze the system and make it treatable, the Authors expressed the evolution of the different types of nodes and edges in terms of moment equations, and proceeded to close the system via a traditional moment-closure approximation. Even though this yields a final closed system that is still not generally amenable to analytic solution, they showed that direct simulations of the model are in good agreement with numerical solutions of the closed system of equations, and observed a remarkable stability of the best-response strategy described above.
Other evolutionary games on populations that are structured beyond pairwise interactions have been recently studied. In particular, Guo and collaborators [25] introduced a general approach that allows studying situations in which indirect interactions impact the strategy of the player. The framework enables one to investigate how the main evolutionary games (prisoner’s dilemma, snowdrift game, stag hunt game, and harmony game, which display different Nash equilibria) compete in simplicial structures and how such a competition is modulated by the unbalance of 2- and 1-simplices, which accounts for the relative prevalence of pairwise or group interactions amongst the players. The main results that are reported demonstrate that (i) higher order games allow for non-dominant strategists to emerge and coexist with dominant ones, (ii) a transition from dominant defection to dominant cooperation exists as a function of the simplicial structure of the population, and (iii) the 2-simplex interactions are the source of strategy diversity, meaning that increasing the relative prevalence of group interactions always promotes diverse strategic identities of individuals.
Stable social conventions can be overturned by a minority of individuals committed to the cause when the size of this minority overcomes a critical threshold [160], [161]. Various frameworks have been proposed to clarify the nature of this phenomenon, the simplest models being inspired by statistical physics [162], [163]. Among these models, the naming game has been successful in describing the emergence of social norms [164], explaining how a shared convention emerges from local, pairwise interactions of agents with their peers.
Recently, Iacopini et al. in Ref. [165] have proposed a novel game on hypergraphs that generalizes the standard naming game on graphs. In this setting, at each time step a hyperedge
h is randomly selected, and a speaker agent is randomly elected among the nodes in
h, whereas all the other nodes in
h act as hearers. The speaker node randomly selects a name, say A, from its vocabulary and communicates it to the other nodes in
h. In the standard naming game,
h is an edge and therefore there is only one hearer. This implies that the agreement can only be reached if the hearer has A in its vocabulary. Iacopini et al. introduced instead the possibility of group interactions, whereby multiple generalized conditions for an agreement in
h can be made. In their work, the Authors focused on what they call the unanimity condition, which is the least favorable towards achieving the agreement in
h, since it requires that all the hearers should have A in their vocabularies. If an agreement is reached in
h, then with probability
β, which can be viewed as a proxy of the social influence of the speaker, all the nodes in
h converge on the name
A and remove all the other names from their vocabulary, whereas with probability
1−β there is no convergence. If an agreement instead is not achieved, then all the nodes that did not have it already add A to their vocabulary. However, a fraction
p of the agents belong to a committed minority that does not obey to these rules, but rather sticks to the same name and does not change the vocabulary.
Iacopini et al. have shown some peculiar behaviors enabled by higher-order interactions. For instance, arbitrarily small minorities (even a single individual) may overturn the social conventions of a group composed by hundreds of peers. Furthermore, they have highlighted a non-monotonic behavior, whereby either small or large groups succeed in overturning social conventions.

## 4. Collective dynamics in networks with higher order interactions

In this Chapter, we review the main results concerning the emergence and control of collective dynamics in networks with higher order interactions. In particular, we will focus on synchronization and consensus dynamics.
Synchronization is a ubiquitous phenomenon in natural and engineered systems [166], [167], which corresponds to the emergence of a collective behavior wherein all (or some of the) system’s unitary components adjust their behavior into a common evolution.
When only dyadic interactions are considered, a broad range of studies has shed light on the intimate relationships between the topology of a network, its synchronizability, and the properties of the synchronized states. In particular, synchronous behaviors have been observed and characterized in small-world [168], weighted [169], multilayer [170], and adaptive networks [171], [172]. Beyond complete synchronization, other types of synchronization were studied, including remote synchronization [173], [174], cluster states [175], synchronization of node groups [176], as well as chimera [177], [178] and Bellerophon states [179], [180]. Furthermore, the transition to synchronization has been shown to be either smooth and reversible, or abrupt and irreversible (as in the case of explosive synchronization, resembling a first-order like phase transition [181]). When abandoning the limit of pairwise interactions, a much richer scenario emerges, with novel types of synchronized states, which are often inherently prohibited when interactions are only dyadic.
In the context of higher-order interactions, the ability of the network to recover a synchronized behavior after local perturbations has been investigated, thereby providing necessary conditions for synchronization. Obtaining global, sufficient conditions that guarantee the emergence of synchronization regardless of the extent of the perturbation is still an open problem in the literature. Similar to the approaches used in networks with pairwise interactions [182], [183], the viability of Lyapunov-based, global conditions for stability has been recently investigated also in the presence of multi-body interactions.
The Chapter is divided in five sections. The first section is devoted to discuss the synchronization of identical systems, and in particular to describe the extension of the Master Stability Function approach [184] to the case of multi-body interactions of any order. The second section reviews the main results concerning synchronization in non-identical systems, with particular reference to the higher-order Kuramoto model in all its variants. The third section reviews the major contributions in the literature concerning the emergence of consensus dynamics. The fourth section then discusses how pinning control can be applied in the context of directed hypergraphs, pinpointing an explanatory analogy with signed graphs that facilitates control design and enables the derivation of global, sufficient conditions for convergence of the node states towards a desired trajectory. Finally, the fifth section expands the controllability analysis for linear networks to the case of multi-body interactions that make the network a polynomial dynamical system.

### 4.1. Synchronization of identical systems

In this section, we refer to the main results obtained in Refs. [27], [185], where the initial analyses of complete synchronization in
p−uniform hypergraphs [186] were extended to more general higher-order structures and coupling functions. Note that, at variance with the modeling framework proposed in Ref. [187], Refs. [27], [185], [186] considered higher-order interactions that cannot be decomposed as sums of pairwise ones. Namely, Jost et al. in Ref. [185] considered chemical hypergraphs of any rank [130], which are characterized by Laplacian-like interactions. Gambuzza et al. in Ref. [27] generalized the focus to simplicial complexes where the coupling function is synchronization non-invasive, that is, it vanishes when the node states coincide. In these settings, the Authors of Ref. [27], [185] demonstrated that complete synchronization exists, and provided necessary and sufficient conditions for local stability of the synchronization solution to be stable in terms of a Master Stability Function, a well-known method originally developed in Ref. [184] for pairwise coupled systems, and later extended to static [188], time-varying [189], [190], [191], and stochastic network systems [192]. In illustrating this approach, we will follow the treatment and notation used in Ref. [27].

#### 4.1.1. Existence and stability of the synchronization state

Let us then consider an ensemble of
N dynamical systems interplaying over a hypergraph where the maximal size of its hyperedges is
D. The equations ruling the evolution of the ensemble are (24)
ẋi=f(x i)+σ 1∑j 1=1 N a i j 1(1)g(1)(x i,x j 1)+σ 2∑j 1=1 N∑j 2=1 N a i j 1 j 2(2)g(2)(x i,x j 1,x j 2)+⋯
=+σ D∑j 1=1 N…∑j D=1 N a i j 1…j D(D)g(D)(x i,x j 1,…,x j D),i=1,…,N,where
x i(t) is the
m-dimensional vector state describing the dynamics of node
i,
f:R m⟶R m describes the local (and identical for all units) dynamics,
σ 1,…,σ D are real-valued coupling strengths, while
g(d):R(d+1)×m⟶R m (
d=1,…,D) are synchronization non-invasive functions (i.e.,
g(d)(x,…,x)=0 for all
d) ruling the interaction forms at different orders. Furthermore,
a i j 1…j d(d) is the entry
i j 1…j d of the adjacency tensor
A(d), for
d=1,…,D.
For ease of illustration, we focus on the case of
D=2, which is described by the following set of coupled differential equations (25)
ẋi=f(x i)+σ 1∑j=1 N a i j(1)g(1)(x i,x j)+σ 2∑j=1 N∑k=1 N a i j k(2)g(2)(x i,x j,x k),where
σ 1 and
σ 2 are the coupling strengths associated to two- and three-body interactions.
Note that the existence and invariance of the synchronized solution
x s(t)=x 1(t)=⋯=x N(t) is guaranteed by the fact that the coupling functions are non invasive, and therefore they vanish when the node states coincide.
As for the stability of such synchronized solution, one can consider perturbations around the synchronous state, i.e.,
δ x i=x i−x s, and perform a linear stability analysis of Eq. (25). The result is (26)
δ ẋi=J f(x s)δ x i+σ 1∑j=1 N a i j(1)[∂g(1)(x i,x j)∂x i|(x s,x s)δ x i+∂g(1)(x i,x j)∂x j|(x s,x s)δ x j]
=+σ 2∑j=1 N∑k=1 N a i j k(2)[∂g(2)(x i,x j,x k)∂x i|(x s,x s,x s)δ x i+∂g(2)(x i,x j,x k)∂x j|(x s,x s,x s)δ x j
=+∂g(2)(x i,x j,x k)∂x k|(x s,x s,x s)δ x k],where
J f(x s) is the
m×m Jacobian matrix of the function
f, evaluated at
x s. It is fundamental to observe that, as the coupling functions are synchronization non invasive, their value is constant and equal to zero at the synchronization manifold, and consequently their total derivative vanishes as well, which implies (27)
∂g(1)(x i,x j)∂x i|(x s,x s)+∂g(1)(x i,x j)∂x j|(x s,x s)=0,
∂g(2)(x i,x j,x k)∂x i|(x s,x s,x s)+∂g(2)(x i,x j,x k)∂x j|(x s,x s,x s)+∂g(2)(x i,x j,x k)∂x k|(x s,x s,x s)=0.It is then possible to simplify Eq. (26) as (28)
δ ẋi=J f(x s)δ x i+σ 1∑j=1 N a i j(1)J g(1)(x s,x s)(δ x j−δ x i)
=+σ 2∑j=1 N∑k=1 N a i j k(2)[J 1 g(2)(x s,x s,x s)(δ x j−δ x i)+J 2 g(2)(x s,x s,x s)(δ x k−δ x i)],where we used the notation (29)
J g(1)(x s,x s)=∂g(1)(x i,x j)∂x j|(x s,x s),J 1 g(2)(x s,x s,x s)=∂g(2)(x i,x j,x k)∂x j|(x s,x s,x s),J 2 g(2)(x s,x s,x s)=∂g(2)(x i,x j,x k)∂x k|(x s,x s,x s).Recalling the definitions of degree
k i(1)=∑j=1 N a i j(1) and generalized degree
k i(2)=∑j=1 N∑k=1 N a i j k(2)/2 of a node
i, Eq. (28) can be then recast as (30)
δ ẋi=J f(x s)δ x i−σ 1∑j=1 N L i j(1)J g(1)(x s,x s)δ x j
=−σ 2∑j=1 N∑k=1 N τ i j k[J 1 g(2)(x s,x s,x s)δ x j+J 2 g(2)(x s,x s,x s)δ x k],where
L i j(1) is the entry
i j of the order 1 Laplacian of the simplicial complex, and
τ i j k is the element
i j k of a tensor
T∈R N×N×N, defined as
τ i j k=2 k i(2)δ i j k−a i j k(2), with
δ i j k equal to 1 if
i=j=k, and 0 otherwise.
The Authors of Ref. [27] then noted that the last term on the right hand side of Eq. (30), though not necessarily referring to diffusive three-body interactions, can be in fact rewritten in terms of the generalized Laplacian matrix. Indeed, one can notice that the Jacobian matrices
J 1 g(2)(x s,x s,x s) and
J 2 g(2)(x s,x s,x s) are both independent on
k and
j, and therefore Eq. (30) becomes (31)
δ ẋi=J f(x s)δ x i−σ 1∑j=1 N L i j(1)J g(1)(x s,x s)δ x j
=−σ 2[∑j=1 N J 1 g(2)(x s,x s,x s)δ x j∑k=1 N τ i j k+∑k=1 N J 2 g(2)(x s,x s,x s)δ x k∑j=1 N τ i j k].Then, using the fact that
T is symmetric, one has that
∑k τ i j k=∑k τ i k j, which yields (32)
δ ẋi=J f(x s)δ x i−σ 1∑j=1 N L i j(1)J g(1)(x s,x s)δ x j
=−σ 2[∑j=1 N J 1 g(2)(x s,x s,x s)δ x j L i j(2)+∑k=1 N J 2 g(2)(x s,x s,x s)δ x k L i k(2)]=J f(x s)δ x i−σ 1∑j=1 N L i j(1)J g(1)(x s,x s)δ x j
=−σ 2∑j=1 N L i j(2)[J 1 g(2)(x s,x s,x s)+J 2 g(2)(x s,x s,x s)]δ x j.
Eq. (32) can be rewritten in block form by introducing the vector
δ x=[δ x 1 T,δ x 2 T,…,δ x N T]T and by using the notation
JF=J f(x s),
JG(1)=J g(1)(x s,x s) and
JG(2)=J 1 g(2)(x s,x s,x s)+J 2 g(2)(x s,x s,x s). The result is (33)
δ ẋ=I N⊗JF−σ 1 L(1)⊗JG(1)−σ 2 L(2)⊗JG(2)δ x,where
⊗ is the direct product and
L(2) is the generalized Laplacian of order 2.
The final step to study the synchronization properties of Eq. (25) can be made by noticing that the generalized Laplacians
L(d) of any order
d are symmetric real-valued zero-row-sum matrices. Therefore, they are all diagonalizable and, for each one of them, the set of their eigenvalues is made of real non-negative numbers, and the corresponding set of eigenvectors constitutes a orthonormal basis of
R N.
Being endowed with the zero-row-sum property, all these generalized Laplacians have, as the smallest of their eigenvalues,
λ 1≡0, whose associated eigenvector
1 N(1,…,1)T is aligned along the synchronization manifold.
In general, however, these generalized Laplacians are not commuting matrices, and therefore the sets of eigenvectors corresponding to all others of their eigenvalues are different from one another. Nonetheless, any perturbation to the synchronization manifold can be expanded as linear combination of the eigenvector sets associated to any of these generalized Laplacians, which can be then used as the basis of the transverse space, with all other eigenvector sets that will be mapped to such a basis by means of unitary matrix transformations.
For instance, one can take as reference basis the one constituted by the eigenvectors of the classic Laplacian
L(1), and consider the new variable
η=(V−1⊗I N)δ x, where
V is the matrix obtained by juxtaposing column-wise the eigenvectors of
L(1). By doing so, one obtains
η̇=(V−1⊗I N)[I N⊗JF−σ 1 L(1)⊗JG(1)−σ 2 L(2)⊗JG(2)](V⊗I N)η.Furthermore, one can take into account that
V−1 L(1)V=diag(λ 1,λ 2,…,λ N)=Λ(1), where
0=λ 1<λ 2≤⋯≤λ N are the eigenvalues of
L(1), and denote by
L̃(2)=V−1 L(2)V the transformed generalized Laplacian of order 2. The result is (34)
η̇=I N⊗JF−σ 1 Λ(1)⊗JG(1)−σ 2 L̃(2)⊗JG(2)η.
Considering that also
L(2) is zero-row sum, Eq. (34) can be rewritten as (35)
η̇1=JF η 1,η i̇=(JF−σ 1 λ i JG(1))η i−σ 2∑j=2 N L̃i j(2)JG(2)η j,i=2,…,N.,where
η i=[η i(1),…,η i(m)]T, for all
i=1,…,N.
Looking at Eqs. (35), one easily sees that the dynamics of the linearized system is decoupled into two parts:
η 1 accounts for the motion along the synchronous manifold, and all other variables
η i, with
i=2,…,N, describe the transverse modes and are coupled by means of the (known) coefficients
L̃i j(2) of the transformed Laplacian. Note that the dynamics of
η 2,…,η N in general cannot be decoupled due to the non commutativity of the generalized Laplacian matrices.
Therefore, the stability of the synchronization manifold can be assessed by (i) simulating a single, uncoupled, nonlinear system; (ii) using the obtained trajectory to compute the time-varying elements of the Jacobians
JG(1) and
JG(2); (iii) simulating the dynamics of a system of
N−1 coupled linear equations; and (iv) tracking the behavior of the norm
∑i=2 N∑j=1 m(η i(j))2 to compute the maximum Lyapunov exponent
Λ max associated to all transverse modes. The necessary condition for stability of the synchronous solution is that
Λ max be negative.
The maximum conditional Lyapunov exponent
Λ max can be computed as a function of the two coupling strengths
σ 1 and
σ 2, and of the pairwise and three-body interaction topology, that is,
Λ max=Λ max(σ 1,σ 2,L(1),L(2)). Similar to the classification made for complex networks (see Ref. [4], Chapter 5), once the structure of the simplicial complex and the coupling functions
g(1) and
g(2) have been given, three classes of synchronization problem can be observed:

- 1.
  class I problems, where
  Λ max is positive in all the half plane
  (σ 1≥0,σ 2≥0), and therefore synchronization is never stable;
- 2.
  class II problems, for which
  Λ max is negative within a unbounded area of the half plane
  (σ 1≥0,σ 2≥0); and
- 3.
  class III problems, for which
  Λ max is negative in a bounded region of the half plane
  (σ 1≥0,σ 2≥0), and therefore an excessive increase of the coupling strengths may make the synchronous manifold unstable.

Examples of class II and class III problems are shown in Fig. 26 for Rössler oscillators [193] coupled on a hypernetwork. The reported numerical simulations of the system Eqs. (25) confirms the accuracy of the predictions made by solving Eq. (35).
The derivation of
Λ max can be extended to hypergraphs of any order
D, as each term on the right hand side of Eq. (24) can, indeed, be manipulated following the same steps as above. Therefore, calling
JG(d)=J 1 g(d)(x s,…,x s)+J 2 g(d)(x s,…,x s)+⋯+J d g(d)(x s,…,x s), Eq. (33) can be then generalized as (36)
δ ẋ=[I N⊗JF−σ 1 L(1)⊗JG(1)−σ 2 L(2)⊗JG(2)−⋯−σ D L(D)⊗JG(D)]δ x.

1. Download: Download high-res image (331KB)
2. Download: Download full-size image

Fig. 26. Stability of synchronization in the parameter space. Contour plots of the time averaged (over an observation time
T=500) synchronization error
E=∑i,j=1 N‖x j−x i‖2/(N(N−1))T (see the vertical bars of each panel for the color code) in the plane (
σ 1,
σ 2) for the simplicial complex sketched in the top-right inset of each panel. The performed simulations correspond to coupled Rössler oscillators, with the state of the
ith node,
i=1,…,N, being
x i=(x i,y i,z i)T and the local dynamics given by
f(x i)=[−y i−z i,x i+a y i,b+z i(x i−c)]T with parameters in their chaotic regime (
a=b=0.2,
c=9). The coupling functions are
g(1)(x i,x j)=[x j−x i,0,0]T and
g(2)(x i,x j,x k)=[x j 2 x k−x i 3,0,0]T in panels (a-c),
g(1)(x i,x j)=[x j−x i,0,0]T and
g(2)(x i,x j,x k)=[0,y j 2 y k−y i 3,0]T in panel (d), whereas
g(1)(x i,x j)=[0,y j−y i,0]T and
g(2)(x i,x j,x k)=[x j 2 x k−x i 3,0,0]T in panel (e). The blue continuous lines correspond to the theoretical predictions of the synchronization thresholds obtained from Eq. (35). Notice that panels (a-c) are examples of class III regimes, whereas panels (d-e) are examples of class II regimes.

Figure reprinted from Ref. [27].
Once again, choosing the eigenvector set which diagonalizes
L(1), and introducing the new variables
δ η=(V−1⊗I n)δ x, yields (37)
η̇1=JF η 1,
η i̇=(JF−σ 1 λ i JG(1))η i−σ 2∑j=2 N L̃i j(2)JG(2)η j−⋯−σ D∑j=2 N L̃i j(D)JG(D)η j,i=2,…,N,where the coefficients
L̃i j(d) result from transforming
L(d) with the matrix
V that diagonalizes
L(1). Conceptually, one has the same reduction of the problem to a single, uncoupled, nonlinear system, plus a system of
N−1 coupled linear equations, from which the maximum Lyapunov exponent
Λ max=Λ max(σ 1,σ 2,…,σ D,L(1),L(2),…,L(D)),can be extracted and monitored in the
D−dimensional space of the coupling strength parameters.

#### 4.1.2. Special cases

Ref. [27] also discusses some special cases in which either the topology of the connectivity structure, or the coupling functions, allow for a further manipulation of the stability equations. The first case is that of an all-to-all coupling, where the classical Laplacian matrix is
L i j(1)=−1 for i≠j,N−1 for i=j.In this case, it is easy to rewrite
L(2) (and any other generalized Laplacian), because the off-diagonal term
L i j(2) (
i≠j) is the number of triangles formed by the link
(i,j), which is
N−2 for any pair
i j, and the term of the main diagonal
L i i(2) is the number of triangles having the node
i as a vertex, which is
k i(2)=N−1 2=(N−1)(N−2)2.Therefore, one has
L(2)=(N−2)L(1), and
δ ẋi=J F δ x i−∑j=1 N L i j(1)σ 1 J G(1)+σ 2(N−2)J G(2)δ x j.By expanding
δ x on the basis of the eigenvectors of
L(1), and having that
λ 2=⋯λ N=N, one obtains (38)
η̇i=[JF−σ 1 N J G(1)−σ 2 N(N−2)J G(2)]η i,for all
i=2,…,N, which implies that, in the all-to-all case, the variables
η i are all uncoupled to each other, so that one can define a Master Stability Function (MSF)
Λ max which only depends on
σ 1,
σ 2 and
N. A similar analysis can be carried out for all-to-all simplicial complexes of any order
d, one would obtain that
Λ max=Λ max(σ 1,…,σ d,N).
Another interesting case, which we will illustrate for clarity on the case
d=2, is that of generalized diffusion interactions with natural coupling functions. Namely, such functions are diffusive as they fulfill (39)
g(1)(x i,x j)=h(1)(x j)−h(1)(x i)
g(2)(x i,x j,x k)=h(2)(x j,x k)−h(2)(x i,x i),where
h(1):R m⟶R m and
h(2):R 2 m⟶R m, and have the property of natural coupling, that is, (40)
h(2)(x,x)=h(1)(x).
The meaning of Eq. (40) is that the coupling to node
i from two-body and three-body interactions is structurally similar, in that a three-body interaction where two nodes have the same state is equivalent to a two-body interaction. In this case, the Master Stability Function
Λ max assumes a particularly convenient form, because it can be written as a function of a single parameter. The consequence of Eq. (40) is indeed that
J 1 h(2)(x s,x s)+J 2 h(2)(x s,x s)=J h(1)(x s). Therefore, one has (41)
δ ẋi=J f(x s)δ x i−σ 1∑j=1 N L i j(1)J h(1)(x s)δ x j−σ 2∑j=1 N L i j(2)J h(1)(x s)δ x j(42)
=J f(x s)δ x i−∑j=1 N σ 1 L i j(1)+σ 2 L i j(2)J h(1)(x s)δ x j. Now, one can consider the zero-row-sum, symmetric, effective matrix
M, given by
M=L(1)+r L(2),r=σ 2 σ 1.The eigenvalues of
M depend on the ratio
r of the coupling coefficients, and one has that (43)
δ ẋi=J f(x s)δ x i−σ 1∑j=1 N M i j J h(1)(x s)δ x j.
Eq. (43) allows to establish a formal full analogy between the case of a hypernetwork and that of a network on a graph with weights given by the coefficients of the effective matrix
M. In particular, a single-parameter Master Stability Function can be defined from the following
m−d i m e n s i o n a l linear parametric variational equation (44)
η̇=J f(x s)−α J h(1)(x s)η,from which the maximum Lyapunov exponent is calculated:
Λ max=Λ max(α) where
α=λ(σ 1 L(1)+σ 2 L(2)) or, equivalently,
α=σ 1 λ(L(1)+r L(2))=σ 1 λ(M).
The effective Laplacian
M therefore plays a similar role to that of the classic Laplacian in complex networks and, considering the nontrivial case where
Λ max(0)>0, and one can again define three possible classes of problems:

- 1.
  class I problems, for which the curve
  Λ max=Λ max(α) is always positive. In this case synchronization cannot be attained regardless of the higher-order structure connecting the dynamical systems;
- 2.
  class II problems, for which
  Λ max=Λ max(α) intercepts the abscissa once at
  α c. In this case, the synchronization threshold can be obtained from the equation
  σ 1 critical=α c/λ 2[M(σ 1 critical,σ 2 critical)]. When
  r is given, local stability can be observed when
  σ 1>σ 1 critical=α c/λ 2(M(r)), whereby the synchronization threshold on
  σ 1 scales with the inverse of the second smallest eigenvalue of the effective matrix;
- 3.
  class III problems, for which the curve
  Λ max=Λ max(α) has multiple crossings with the abscissa. For instance, when the crossing points are at
  α 1 and
  α 2>α 1, synchronization can be observed only if the spectrum of
  M is such that the conditions
  σ 1 λ 2(M)>α 1 and
  σ 1 λ N(M)<α 2 are simultaneously verified. In this case, the ratio
  λ 2(M)/λ N(M) can be then considered a synchronizability metric.

This implies that, for the case of natural coupling, the steps required to study synchronizability reduce to two. The first consists in computing
Λ max(α), and only pertains to node dynamics and coupling functions, whereas the second depends on the higher-order coupling configurations, and requires one to check whether
Λ max is negative for
α in
{σ 1 λ 2 M,…,σ 1 λ N(M)}. Gambuzza et al., [27] illustrated the different MSF types on networks of Rössler oscillators and Lorenz systems [193], [194], respectively, which are paradigmatic examples of chaotic systems with three-dimensional state. Namely, the state of the
i-unit in the system will be
x i=(x i,y i,z i)T∈R 3, and the individual dynamics for the Rössler is described by the vector field
f(x)=[−y i−z i,x i+a y i,b+z i(x i−x)]T, while for the Lorenz system is
f(x)=[y i−x i,x i(ρ−z i)−y i,x i y i−β z i]T. For both systems, the parameters have been selected so that the resulting dynamics is chaotic (for the Rössler oscillator
a=b=0.2,
c=9, and for the Lorenz system
σ=10,
ρ=28, and
β=8/3).
The different MSF types for both individual dynamics have been then computed with nine different kinds of natural coupling functions, and are depicted in Fig. 27 for the Rössler oscillator and in Fig. 28 for the Lorenz system. Both cases exhibit a variety of behaviors that encompass all the three possible classes of problems. In the case of Rössler oscillator, there is one class III example [Fig. 27(a)], one class II example [Fig. 27(e)], while all the remaining cases correspond to class I. In the case of the Lorenz system, there are several examples of class I problems [Fig. 28(c), (f), (g), and (h)], three of class II [Fig. 27(a), (d), and (e)], and two of class III, one characterized by a narrow parameter region for synchronization [Fig. 27(b)], and the other where the MSF assumes negative values in two different intervals of
α [Fig. 27(i)].

1. Download: Download high-res image (323KB)
2. Download: Download full-size image

Fig. 27. Synchronization in hypergraphs of Rössler oscillators, for the case of a natural coupling. The Master Stability Function is computed for each of the nine selection of the coupling function
h(2). The corresponding expression for
h(1) can be derived taking into account condition (40).

Figure reprinted from Ref. [27].

1. Download: Download high-res image (327KB)
2. Download: Download full-size image

Fig. 28. Synchronization in hypergraphs of Lorenz systems, for the case of a natural coupling. The Master Stability Function is here calculated for several coupling functions. On the top of each panel, the expression used for
h(2) is reported. The corresponding expression for
h(1) can be derived taking into account condition (40).

Figure reprinted from Ref. [27].

### 4.2. Synchronization of non-identical systems

The study of synchronization of non identical systems in the presence of higher order interactions has so far mostly concentrated on the extension of the well-known Kuramoto model [195], [196]. In the 1980s, Kuramoto demonstrated that an ensemble of phase oscillators in pairwise interaction through a sinusoidal coupling may achieve synchronization depending on the intensity of their coupling strength [197]. Namely, Kuramoto discovered that an ensemble of all-to-all coupled non-identical oscillators (each one of them rotating with an intrinsic frequency drawn from a given distribution) experienced a second-order phase transition to global frequency synchronization (where all oscillators rotate with the same frequency, equal to the average of the intrinsic frequency of each oscillator) as the pairwise coupling strength increased. The elegant insight on the onset of synchronization in coupled dynamical systems, which was explained as a pitchfork supercritical bifurcation leading to a second-order phase transition, attracted an intense research effort and inspired research on other types of coupled oscillators, thereby making the Kuramoto oscillators’ model a prototype to study the origin of synchronization.
In its standard incarnation, the Kuramoto oscillator model is described by (45)
θ i̇=ω i+K 1 N∑j=1 N sin(θ j−θ i),i=1,…,N,where
θ i is the phase and
ω i the intrinsic frequency of the
ith oscillator, respectively,
K 1>0 is the coupling strength, and
N is the total number of oscillators. The complex order parameter (46)
z 1=1 N∑j=1 N e ı θ j=r 1 e i ψ 1 corresponds to the centroid of the phases, with
ψ 1 being the average phase, and
r 1 measuring the phase coherence. In particular, a completely incoherent state corresponds to
r 1∼0, whereas phase synchronization to
r 1∼1. A detailed description of the second-order transition from an incoherent to a synchronized state as a consequence of an increase in
K 1 is contained in Ref. [198].
The Kuramoto model has been extended beyond all-to-all coupling and modified to also include multiple layers of coupled oscillators [199], [200]. Its properties have been investigated for all possible graphs describing the pairwise interactions between the oscillators, and it has been applied to capture a wide range of biological and social systems [167], [195]. Recently, the model has been generalized to also include higher-order interactions, which have been shown to induce a richer behavioral repertoire. For instance, it has been shown that the Kuramoto model may exhibit abrupt desynchronization when three-body interactions among all the oscillators are added to [201] or replace [26] the all-to-all pairwise interactions of the original model. Similar results have been obtained with a non-symmetric variation of the Kuramoto model in which the microscopic details of the interactions among the phase oscillators are described in the form of a simplicial complex [202]. A different approach has been proposed by Millán et al. [203], who have formulated a higher-order Kuramoto model in which the oscillators are placed not only on the nodes but on higher-order simplices, such as links, triangles, and so on, of a simplicial complex.
In what follows, we first focus on higher-order Kuramoto models with oscillators at the nodes, and start by analyzing the case where triadic interactions replace pairwise ones, then consider the case in which general higher-order interactions combine with standard, pairwise ones; then we briefly overview the works that consider higher-order interactions in multilayer Kuramoto oscillators, and discuss the effect of adapting the interaction strengths on the synchronization/desynchronization transitions. Finally, we account for the work where oscillators are placed not only at nodes but also at higher-order simplices.

#### 4.2.1. Higher-order interactions

To illustrate the effects of higher-order interactions, let us start by considering a first instance in which, as suggested in Ref. [26], triadic interactions replace pairwise ones: (47)
θ i̇=ω i+K 2 N 2∑j=1 N∑k=1 N sin(θ j+θ k−2 θ i),where
K 2 is the coupling strength for the three-body interactions, which is divided by
N 2 to account for the double summation in the coupling protocol.
The Authors of Ref. [26] showed both numerically and analytically how the presence of three-body interactions gives rise to a different and richer set of behaviors compared to the classic model (45). Their numerical analysis focused on a system with
N=1 0 5 oscillators, intrinsic frequencies extracted from a Lorentzian distribution with zero mean and spread equal to 1, and initial phases
θ i(0) set to 0 with probability
η, and to
π with probability
1−η, whereby
η can be viewed as an asymmetry parameter for the initial phase distribution. In the numerical experiment, the simulation begins with
K 2=16 and, after the steady state was reached, the parameter
K 2 is slowly decreased to 0, and then it is slowly reset back to 16. The experiment was repeated for different values of the asymmetry parameter
η, which was varied between 0.8 and 1, with step 0.05. As illustrated in Fig. 29, when decreasing
K 2, the solutions traverse different partial synchronization branches, with higher values of the order parameters
r 1 corresponding to higher values of
η, indicating multistability. An abrupt desynchronization transition appears at different critical values for
K 2, with a continuum of transitions with the order parameter dropping to
0. On the other hand, when
K 2 is increased back to 16, no spontaneous abrupt synchronization transition is observed.
These numerical observations are mathematically backed in Ref. [26] by an analytical description obtained through a variation of the Ott–Antonsen ansatz [204]. The presence of clustered oscillators in the desynchronized states has been used then in Ref. [205] to explore the possibility of encoding and storing sequences of bits, analytically determining which of the
2 N potential steady states are stable.
At variance with Eqs. (47), where the higher order interactions replaced the pairwise ones, another bulk of work considered the simultaneous presence of pairwise and higher-order interactions. While the pairwise interaction function is always selected as in the classic Kuramoto model, different choices have been made for higher-order interaction functions [202], [206], [207]. In what follows, we start by discussing the dynamics of an extended Kuramoto model where the coupling functions have been chosen as in Ref. [202], for which the dynamical evolution of the phases can be written as; (48)
θ i̇=ω i+K 1〈k[1]〉∑j=1 N A i j sin(θ j−θ i)+K 2 2〈k[2]〉∑j=1 N∑k=1 N B i j k sin(2 θ j−θ k−θ i)+K 3 6〈k[3]〉∑j=1 N∑k=1 N∑l=1 N C i j k l sin(θ j+θ k−θ l−θ i),where
K 1 and
K 3 are the pairwise and four-body coupling strengths,
A i j is the entry
i j of the adjacency matrix
A describing the pairwise interactions,
B i j k is the entry
i j k of the adjacency tensor
B describing the three-body interactions, and
C i j k l is the entry
i j k l of the adjacency tensor
C describing the four-body interactions, and
〈k[d]〉 is the average generalized
d degree of the hypergraph describing the coupling between the oscillators.

1. Download: Download high-res image (264KB)
2. Download: Download full-size image

Fig. 29. Multistability in the de-synchronization transition. The order parameter
r 1 as a function of the coupling strength
K 2=K for different values of the asymmetry parameter
η. Blue, red, green, orange, and purple circles are obtained from simulation simulations of Eq. (47) with
N=1 0 5 oscillators for
η equal to 1, 0.95, 0.90, 0.85, and 0.8, respectively.

Figure reprinted from Ref. [26], © 2019 by the American Physical Society.
To obtain an analytic insight on the transitions taking place in Eq. (48), the Authors of Ref. [202] performed a mean-field (all-to-all) approximation, and rewrote the oscillator dynamics as (49)
θ i̇=ω i+1 2 ı(H e−ı θ i−H⋆e ı θ i),where
H=K 1 z 1+K 2 z 2 z 1⋆+K 3 z 1 2 z 1⋆, with
⋆ denoting the complex conjugate, and (50)
z 2=1 N∑j=1 N e 2 ı θ j=r 2 e ı ψ 2,being the 2-cluster order parameter.
In the thermodynamic limit, one can introduce the density function
f(θ,ω,t), so that
f(θ,ω,t)d θ d ω is the fraction of oscillator with phase in
[θ,θ+d θ) and intrinsic frequency in
[ω,ω+d ω) at time
t. As the oscillators are conserved,
f needs to fulfill the continuity condition, whereby (51)
∂f∂t+∂∂θ f(θ̇i)=0.Expanding the density function in terms of its Fourier series, one obtains
f(θ,ω,t)=g(ω)2 π 1+∑n=1∞f n(ω,t)e ι n θ+c.c,with
c.c. being the complex conjugate of the previous term, and
g(ω) an intrinsic frequency distribution.
Next, one can substitute the above Fourier series expansion of the density function and the expression of
θ̇i in Eq. (49) in the continuity equation, and then equate the coefficients of
e ı n θ to calculate the coefficients
f n, which, in principle, are infinite. To tackle this problem, Ott and Antonsen in their landmark paper Ref. [208] proposed the Ansatz
f n(ω,t)=α n(ω,t) for some function
α which is analytic in the complex
ω plane, with
|α|<1 to avoid divergence of the series. Upon inserting this Ansatz into the continuity equation, equating any power of
e ı n θ collapses to a single dynamical equation for
α, given by (52)
α̇=−ı ω α+1 2(H⋆−H α 2).
Additionally, in the thermodynamic limit, the complex order parameter is
z 1=∬f(θ,ω,t)e ı θ d θ d ω=∫α(ω,t)g(ω)d ω.For Lorentzian intrinsic frequency distribution with mean
ω 0 and spread
Δ, one has
z 1⋆=α(ω 0−ı Δ t) and
z 2⋆=z 1⋆2. Substituting this in Eq. (52) evaluated at
ω=ω 0−ı Δ, and taking
Δ=1, one finally obtains the following dynamics for the order parameter
r 1 and for the average phase
ψ: (53)
ṙ1=−r 1+K 1 2 r 1(1−r 1 2)+K 2+3 2 r 1 3(1−r 1 2),
ψ̇=ω 0,where
K 2+3=K 2+K 3. Note that the dynamics of the order parameter can be studied independently of that of
ψ, and that
r 1=0 is always a fixed point of the above equation, which is stable until a critical coupling strength and becomes unstable through a subcritical or supercritical pitchfork bifurcation. The other fixed points are (54)
r 1=K 2+3−K 1±(K 1+K 2+3)2−8 K 2+3 2 K 2+3.Note that the specific form of Eq. (53), with the fixed points depending on the sum of
K 2 and
K 3, depends on the choice of the coupling functions for the triadic and four-body interaction terms in Eq. (48), and alternative choices of these functions may not yield a closed-form expression for the dynamics and fixed-points of
r 1.
In an all-to-all setting with
N=1 0 4 oscillators, Ref. [202] showed an excellent matching between the numerical simulations and the theoretical predictions from Eq. (53). This is illustrated in Fig. 30 by plotting the order parameter
r 1 as a function of the pairwise coupling strength
K 1 for higher-order couplings such that
K 2+3 is equal to 0, 2, 5, and 10. This numerical experiment illustrates that, when the higher order coupling is sufficiently small, a second order transition to synchronization occurs through a supercritical pitchfork bifurcation, taking place at a critical value of
K 1 sync=2 for the pairwise coupling; however, when
K 2+3 overcomes a critical threshold (see, for instance, the cases
K 2+3=5,8,10 in Fig. 30) the synchronization branch folds over itself, thus giving rise to hysteresis and abrupt transitions between an incoherent and a synchronized behavior. In this regime, at
K 1 sync a subcritical (rather than supercritical) pitchfork bifurcation occurs. Furthermore, at a lower value
K 1 desync<K 1 sync the synchronized branch first emerges through a saddle node bifurcation. Interestingly, it can be noted that when
K 2+3 is larger than 8, the synchronized branch emerges at negative values of
K 1, this meaning that higher-order interactions may stabilize synchronized states even in the presence of repulsive pairwise interactions.
In a later work Ref. [207], Kovalenko et al. have been able to numerically illustrate and analytically explain that a certain degree of synchronization may be observed also when both pairwise and triadic interactions are repulsive. This intriguing result has been obtained in an ensemble of
N globally coupled oscillators, described by (55)
θ i=ω i+K 1 N∑j=1 N sin(θ j−θ i)+K 2 2 N 2∑j=1 N∑k=1 N sin(θ j+θ k−2 θ i),for
i=1,…,N. Note that the triadic coupling function considered here is the same as model (47), but different from the one selected for model (48). System (55) was studied in Ref. [207] by rewriting its dynamics as (56)
σ̇i=W i σ i+K 1(ρ−(ρ⋅σ i)σ i)+K 2(ρ⋅σ i)(ρ−(ρ⋅σ i)σ i),where
σ i is the unit vector
[cos(θ i),sin(θ i)]T,
ρ is the mean vector
∑i=1 N σ i/N, and
W i is the antisymmetric matrix
0 ω i−ω i 0.Note that the use of Eq. (56) can be extended to capture the dynamics not only of
2-dimensional unit vectors, but also of general
D-dimensional vectors whose dynamics lie on the unit sphere
S D−1. Kovalenko et al. have been able to analytically study Eq. (56) for both
D=2 and
D=3, obtaining a self-contained equation for the order parameter
r 1=|ρ|=|z 1|.

1. Download: Download high-res image (234KB)
2. Download: Download full-size image

Fig. 30. Abrupt synchronization in simplicial complexes: All-to-all case. Synchronization profiles describing the macroscopic system state: the order parameter
r=r 1 as a function of 1-simplex coupling
K 1 for higher order coupling
K 2+3=0,2,5,8, and 10 (blue to red). Solid and dashed curves represent stable and unstable fixed points given by Eq. (54), respectively, and circles denote results taken from direct simulations of Eq. (48) specified for the all-to-all case, and with
N=1 0 4 oscillators. For
K 2+3<2 and
K 2+3>2 the pitchfork bifurcation is supercritical and subcritical, respectively.

Figure reprinted from Ref. [202].
In Fig. 31 both the theoretical prediction and the numerical evaluation of the backward and forward transition to synchronization are described. The theoretical predictions have been obtained by computing the fixed point from the analytic expression for the order parameter
r 1, and then performing a parametric stability analysis that leverages the findings reported in the Supplementary Information of Ref. [207]. From a numerical standpoint, the simulations have been performed on an ensemble of
N=5,000 oscillators, varying the coupling strengths
K 1 and
K 2 between
−30 and
30. The backward transition has been inspected by setting the initial condition
σ i(0) to be all equal to a random unitary vector
e D in
R D, whereas the forward transition has been sought by setting
σ i(0) to
e D or
−e D with equal probability 0.5.
The theoretical and numerical analyses offer further backing to previously observed phenomena, but also support the existence of novel behaviours induced by higher-order interactions. As noted in Ref. [26], in the absence of pairwise interactions, the forward transition is inhibited, see Fig. 31, panels (b) and (f), and this is also true for
D=3, see panels (d) and (h). However, the most intriguing evidence from Fig. 31 is that in the backward transition some synchronization features are retained, and this is illustrated both analytically and numerically for
D=2, see panels (a) and (e), respectively, whereas for the 3-dimensional case the phenomenon was only captured analytically, see panels (c) and (e).
We remark that alternative functional forms for the sinusoidal coupling terms have been investigated, see for instance the four-body interaction function selected in Ref. [206], which differs from the one proposed in model (48). With this choice of the interaction function, and under the assumption of identical intrinsic frequencies, the Authors were able to detect the emergence of chimera states without introducing phase-lag, and used the Ott–Antonsen manifold to obtain an insight on such phenomenon. For the case of identical Kuramoto oscillators, the local stability of the synchronized state has been studied on simplicial complexes of any order [204]. In an effort to generalize the treatment for non-identical oscillators, Bick et al. [209] proposed a framework based on a mean-field formulation involving probability measures, which allowed to consider virtually any type of higher-order coupling functions, and the presence of multiple interacting populations of oscillators. Using this formalism, the Authors identified invariant subspaces of synchrony patterns, and provided a suitable definition of stability that allowed to show how the all-synchronized state can be shown to be locally stable.

1. Download: Download high-res image (428KB)
2. Download: Download full-size image

Fig. 31. The synchronization scenario in parameter space. Theoretical predictions (panels (a) to (d)) and numerical simulations (panels (e) to (h)) for
D=2 and
D=3 of the order parameter as a function of the coupling strengths (
K 1=λ 1,K 2=λ 2). Panels (a), (c), (e), (g) refer to the backward transition, whereas panels (b),(d),(f),(h) to the forward transition. In panels (a)–(d), the values around the line
λ 1=−λ 2 of the quadrant defined by
λ 1>0 and
λ 2<0 are interpolated. In each panel, the synchronized (incoherent) the state is represented by the red (blue) color, and the values of the order parameter
r 1=R are coded according to the vertical color bar reported at the rightmost of the figure. For better representing the results, in the semi-plane
λ 2<0 a vertical logarithmic scale is adopted for
|λ 2|.

Figure reprinted from Ref. [207], © 2021 by the American Physical Society.

#### 4.2.2. Multilayer and adaptive wirings

The first work discussing the simultaneous presence of two layers and of higher order interactions was Ref. [26], where the Authors investigated the impact that the higher-order system Eq. (47) has when acting as a driver layer on a standard, pairwise layer of Kuramoto oscillators. Ref. [26] numerically illustrated how the higher-order layer effectively drives the pairwise one only when the coupling between the layers dominates the pairwise coupling term.
A later work then focused on more general multilayer structures, where all layers may exhibit higher-order interactions described by all-to-all simplicial complexes [210]. Namely, the Authors focused on two interconnected layers, and introduced an order parameter for each of the layers. To reveal the multiple routes to first-order transition to synchronization associated with different basins of initial conditions, Jalan and Suman extended the Ott–Antonsen approach to two-layer networks, thereby analytically obtaining the entire phase diagram for the order parameters as a function of the coupling strengths.
Adaptation in pairwise interactions has been an intensive area of research grounded on an empirical observation that connection strengths between pairs of dynamical units in real-world systems change with time. In networks whose structure and weights evolve with time, the dynamics of evolving units affects their connections and vice versa [211], [212], [213], [214], [215], [216]. One of the adaptation schemes widely considered in Kuramoto oscillators with pairwise coupling is based on the Hebbian learning rule, first postulated by Hebb [217] and later supported by experimental evidence [218]. The Hebbian learning rule relies on the principle that “cells that fire together, wire together”. Indeed, it has been observed that the synaptic coupling between a pair of neurons is strengthened or weakened if the two neurons are or are not simultaneously firing as a consequence of presynaptic and postsynaptic spikes. If the relative spike timing is coded in terms of the phases of the oscillators, a neural network can be described as a network of phase oscillators. The adaptation implemented on connection-weights ensures that such Hebbian learning mechanism on pairwise coupled Kuramoto oscillators drives e.g., the emergence of cluster states [219], mesoscale structures [220], and explosive synchronization [221] in complex networks. In the presence of higher-order interactions, a natural research question is to understand the impact that the adaptation of the hyperedge weights has on the evolution of the oscillator dynamics.
A natural extension of the Hebbian learning rule was proposed in Ref. [222], [223], where the phase-evolution of
N non-identical Kuramoto oscillators coupled on a simplicial complex is described by (57)
θ̇i=ω i+K 1〈k[1]〉∑j=1 N A i j sin(θ j−θ i)+K 2 2〈k[2]〉∑j,k=1 N B i j k sin(2 θ j−θ k−θ i),and, at variance with the static case, the elements
A i j and
B i j k of the adjacency matrices and tensors, respectively, that describe the simplicial complex are adapted as (58)
Ȧi j=α cos(θ j−θ i)−μ A i j,
Ḃi j k=β cos(2 θ j−θ k−θ i)−ν B i j k,where
α,
β,
μ, and
ν are positive scalars.
In the absence of adaptation, Eq. (57) would exhibit a first-order phase transition to global synchronization for various synthetic and real-world higher-order coupling [80], [224], [225], [226]. This transition is illustrated in Fig. 32, panel (a), which reports simulations of model (57) on a random 2-simplicial complex of
N=1 0 3 oscillator, with average degrees
〈k[1]〉=14 and
〈k[2]〉=10, and intrinsic frequency uniformly drawn from in the interval
[−1,1]; and where
λ≔K 1+K 2 is first adiabatically increased from 0, and then adiabatically decreased until 0. The simulations were carried out for different value of the parameter
0≤p=λ 2/λ≤1, which allows to balance the relevance of pairwise interactions, whereby large (small) values of
p corresponds to the dominance of triadic (pairwise) interactions. In the static setting, where all
A i js and
B i j ks are equal to 1, a first-order forward and backward transition is observed for all
p. Increasing
p, the forward transition to synchronization is observed at larger values, and a prominent hysteretic behavior is observed.
On the other hand, when the strength of the interactions is adapted according to Eq. (58), the phase coordination in a unique cluster is prevented, whereby one observes the order parameter
r 1 (denoted
R 1 in Fig. 32) to be approximately zero for all
λ and
p, as illustrated in panel (b). Nonetheless, monitoring the
2-cluster order parameter
r 2 (
R 2 in Fig. 32), one sees that the adaptive model admits a 2-cluster state as a stable solution instead of global synchronized state of the static case. The abrupt transition to a 2-cluster states happens at larger values of
λ, which increases with
p; larger values of
p also correspond to a more prominent hysteretic behavior. It has to be remarked that, as noted in Ref. [222], in the absence of pairwise interactions (that is, when
p=1), there is no forward transition to synchronization for random initial conditions.
To capture the onset of the 2-cluster state, Kachhvah and Jalan [223] used the Ott–Antonsen approach that, for Lorentzian intrinsic frequency distribution with mean
ω 0 and spread
Δ, allowed them to obtain the following differential equation for the 2-cluster order parameter
z 2 [for its definition, see Eq. (50)]: (59)
z 2̇=2 ι ω 0 z 2−2 Δ z 2−λ z 2 2 2 α μ(1−p)z 2⋆+β ν p z 2 z 4⋆+λ 2 α μ(1−p)z 2+β ν p z 2⋆z 4.Looking at Eq. (50), one can then obtain the dynamics of
r 2 and
ψ 2, which can be shown to be decoupled. Studying the fixed points of
r 2 and their stability, a bifurcation analysis can be performed, showing that the forward transition to a 2-cluster state happens through a sub-critical pitchfork bifurcation, with the fixed point
r 2=0 becoming unstable. In Fig. 32, panels (c) and (d) show that the theoretical predictions provide meaningful insight on the numerical results.

1. Download: Download high-res image (411KB)
2. Download: Download full-size image

Fig. 32. Comparing the
r 1 and
r 2 order parameters’ behavior. Panels (a) and (b) report the order parameters
r 1=R 1 and
r 2=R 2 as a function of
λ=K 1+K 2 computed from numerical simulations of
N=1 0 3 oscillators intrinsic frequency uniformly drawn from in the interval
[−1,1], whose dynamics are described by model (57), and coupled on a random 2-simplicial complex with average degrees
〈k[1]〉=14 and
〈k[2]〉=10. In the legend, a superscript
f (
b) corresponds to points obtained when
λ is adiabatically increased (decreased). Panel (a) corresponds to the case of static interactions, where all
A i j and
B i j are 1, whereas panel (b) corresponds to interaction strengths adapted through Eq. (58). In panels (c) and (d), numerical and analytical predictions of the onset of a stable 2-cluster state is illustrated in the
(r 2=R 2,λ) plane for an all-to-all connected 2-simplicial complex of
N=1 0 4 oscillators, with Lorentzian distribution of the intrinsic frequencies with spread
Δ=0.1, and parameters
μ and
ν of the adaptive law (58). In panel (c), the parameters
α and
β of Eq. (58) are set to 1, whereas
p is varied between 0.2 and 0.9; in panel (d),
p is set to 0.9, whereas
α is set equal to
β and varied between 0.25 and 1. The solid and dashed lines are the respective theoretical traces of stable and unstable solutions.

Figures reprinted from Ref. [223], © 2022 by the American Physical Society.

#### 4.2.3. Dynamics on hyperedges

The hypergraph formalism allows one to model dynamical processes where the signals are associated not only to the nodes, but also to hyperedges of different dimensions, which include links, triangles, and so on. Such signals, defined on higher-order structures, are called topological signals, and are of interest in the context of higher-order synchronization [227] and in signal processing [123], [124]. The interaction between topological signals defined on hyperedges of different dimensions in non-trivial. For instance, the synaptic activity and cell body activity of a neuron can be affected by the gliomas in brain tumors [228].
Based on the formalism introduced in Ref. [203], Ghorbanchian et al. in Ref. [229] studied a topological Kuramoto model where the topological signals are defined on the nodes and links, as illustrated in Fig. 33. Focusing on all-to-all simplicial complexes, they leveraged algebraic topology to gather an understanding of the emerging dynamics, whereby phase dynamics was projected on nodes and on triangles by means of suitably defined incidence matrices. Namely, the Authors of Ref. [229] showed that for such a model, the coupling between the signals associated to nodes and links produces an explosive synchronization phenomenon, in which phases defined on nodes synchronize simultaneously to phases defined on the hyperlinks. The proposed model was numerically tested both on synthetic simplicial complexes and on real connectomes [230], [231].
In a later work, Arnaudon et al. [232] further generalized the model to consider weights on any simplices, and both linear and nonlinear frustration terms, yielding a natural simplicial generalization of the Sakaguchi–Kuramoto model [233]. The Authors discovered novel dynamical phenomena, which include the emergence of simplicial phase re-locking in regimes of high frustration.

1. Download: Download high-res image (155KB)
2. Download: Download full-size image

Fig. 33. Schematic comparison of the standard Kuramoto model on a graph with the topological Kuramoto model on a simplicial complex. Panel (a) shows a graph where to each node, say
i, is associated to a dynamical variable, namely its phase
θ i. The synchronization between these variables is captured by the standard Kuramoto model. Panel (b) shows a simplicial complex formed by nodes, links, and triangles (here shaded in orange) in which a dynamical variable is associated not only to nodes, but also to the links, e.g., a variable
ϕ i j will be associated to a link
(i,j). The coupled synchronization dynamics on the simplicial complex are described by the higher-order topological Kuramoto model formalized in Ref. [203], [229].

Figure reprinted from Ref. [229].

### 4.3. Consensus dynamics

Consensus dynamics arises in many fields of science and engineering when a group of agents needs to agree upon a certain quantity of interest. Notable examples include opinion dynamics [234] (where groups of individual interact to find an agreement on a certain topic of discussion), wireless sensor networks [235] (where a given quantity, e.g., the temperature of a room, needs to be estimated in a distributed fashion), or rendez-vous problems in robotics [236] (where an ensemble of robots needs to meet in a certain point of the state space).
Starting from various seminal works of the seventies and eighties [237], [238], [239], [240], researchers have pointed out the topological properties of the interaction graph that were required to achieve consensus in the network. One of the crucial findings was that, for static graphs, the existence of a directed spanning tree in the graph would suffice for consensus in networks of simple integrators, and the analysis was later extended for time-varying [241], state-dependent [242], and dynamic topologies [214]. These models consider pairwise, independent interactions, and do not take into account the multibody effects on agents that belong to the same group, albeit, for instance, the dynamics in a social clique is also shaped by mechanisms of peer influence and reinforcement [22]. Recently, several studies have been conducted to take into account these higher-order interactions, by leveraging the formalism of simplicial complexes and hypergraphs [157], [227], [243], [244], [245], [246], [247].
As we will show in the following, there are two directions in which consensus-like diffusion dynamics can be extended on hypergraphs. One is to maintain the association between states and nodes, whose evolution will depend on the combined states of all the nodes inside each multibody interacting unit. Another alternative is to associate a state variable and a dynamic flow to each simplex.
A first extension of standard, pairwise dynamics is to account for the interactions between triangles of users. With this goal, Neuhäuser et al. [157] considered an ensemble of
N nodes connected through a hypergraph
H composed by a set of undirected triangles
T, which can be described by the adjacency tensor
A∈R N×N×N whose element
i j k is (60)
A i j k=1, if{i,j,k}∈T,0,otherwise.Notice that, since the triangles in
T are undirected,
A is symmetric in all dimensions, and that the Authors of Ref. [157] assumed the absence of pairwise interactions. Each node in the system, say
i, is associated to a scalar state variable
x i and its nonlinear consensus dynamics can be written as (61)
ẋi=∑j,k=1 N A i j k f i(j,k)(x i,x j,x k),i=1,…,N,where
f i(j k):R 3→R is the interaction function associated to the triangle
(i,j,k) and its symmetries are aligned with those of
A, whereby
f i(j k)=f i(k j)=f i{j k} for all
i. The Authors observed that, when the interaction functions are linear, the three-body interactions can be expressed as linear combinations of pairwise interactions, and therefore (61) could be studied as a standard consensus problem on an equivalent undirected graph. With the goal of mimicking multibody phenomena that are relevant in sociology, such as peer pressure [248] and homophily [249], Neuhäuser et al. then focused on the following form for the interaction functions: (62)
f i{j,k}=1 2 s(|x j−x k|)(x j−x i)+(x k−x i),thereby assuming the interaction functions being the same for each triangle. In analogy with the classic consensus model on graphs, model (61)–(62) was called three-way consensus model (or 3CM) by the Authors. If the function
s is monotonically decreasing, then the effect of nodes
j and
k on node
i is reinforced if their states are similar, whereas it is reduced if their states are very diverse.
Note that the interaction function (62) is symmetric with respect to
i and
j, and therefore it is consistent with the choice of an underlying undirected topology. It is possible to then describe the network dynamics in terms of a time-varying, weighted graph. To this aim, the Authors introduced a state-dependent adjacency matrix
A(t), whose element
i j can be defined as (63)
A i j(t)=∑k A i j k s(|x j−x k|)=∑k∈I i j s(|x j−x k|),where
I i j is the index set of the nodes that form a triangle with
i and
j. Denoting
ℒ the state-dependent Laplacian matrix associated to
A, it is then possible to rewrite the dynamics (61)- (62) as (64)
ẋi=−∑j ℒ i j(t)x j,i=1,…,N,which, in compact matrix form, becomes (65)
ẋ=−ℒ(t)x,where
x=[x 1;…;x N]. Interestingly, the Authors showed that, albeit the hypergraph is undirected, the time-varying adjacency matrix
A is not symmetric, and therefore the average state is not constant as in consensus over undirected graphs. The way the average state shifts will depend on the interplay between the influence function
s in Eq. (62), the initial conditions of the node states, and the hypergraph topology. In particular,

- 1.
  In the mean field, that is, when the hypergraph is fully connected, the shift only depends on the initial state of the node. In particular, when the initial state is uniformly distributed in a finite interval, the average state does not change, whereas an asymmetric distribution yields a variation in the average state. In the latter case, the shift will depend on the selection of the interaction function
  s, as illustrated in Fig. 34. In particular, the Authors showed that, assigning binary initial conditions, with 80% of the nodes initialized with 0 and 20% with 1, a monotone decreasing interaction function yields a shift of the average value towards the majority, whereas a monotone increasing interaction function tends to balance opinions.
- 2.
  When the mean field assumption is dropped, the topology of the hypergraphs will affect the final consensus value and the rate of convergence. For instance, Neuhäuser et al. considered a network of 100 nodes, with two fully-connected clusters of 10 nodes, named A and B where the nodes of cluster A had initial condition 0, whereas the nodes of B were initialized with 1. The two clusters were connected by 80 randomly placed triangles, so that a fraction
  p∈[0,1] of them is directed towards
  A, in the sense that two of its nodes belong to
  B and one to
  A. The Authors varied
  p from 0 to 1 with step 0.25, and computed the average and standard deviation of the final consensus value. Despite the fact that the hyperedges are undirected,
  p≠0.5 determines an asymmetric flow of influence between the two clusters, with the average consensus value shifting from 0 to 1 as
  p varies from 0 to 1, and with a higher standard deviation the closer
  p is to 0.5, as it can be seen in the left panel of Fig. 35.

The study of consensus dynamics over hypergraphs was later extended to topologies with hyperedges of any dimension [244], [246], and to the case in which the hypergraph topology changes over time [243], [247]. Sahasrabuddhe et al. [246] introduced a general model of consensus dynamics over hypergraphs which do not pose any restriction on the cardinality of the hyperedge. Let
H={V,E} a generic hypergraph, with
V and
E being the sets of its nodes and hyperedges, respectively. The dynamics of the 3CM in Eq. (61) on
H then becomes (66)
ẋi=∑h∈E f i h(x i,x j,x k,…),i=1,…,N,where (67)
f i h(x i,x j,x k,…)=∑j∈h s i j(x i,x j,x k,…)(x j−x i),i∈h,0, otherwise,where
s i j(h), invariant with respect to any permutation of the indices
k∈h that are neither
i nor
j, captures the overall influence of the hyperedge
h. In the special case where
s i j=s i for all
i,j, that is, when the scaling functions are independent of acting node, model (66) can be rewritten as (68)
ẋi=∑n=2 N 1(n−2)!∑j k…A i j k…(n)s i(x i,x j,x k,…)(x j−x i),i=1,…,N,where
A(i) is the adjacency tensor describing all interactions through hyperedges of cardinality
i. As in the case of three-way interactions, from Eq. (68) it is possible to introduce an equivalent state-dependent Laplacian matrix
ℒ(t) so that (69)
ẋ=−ℒ(t)x.Alternative choices of the nonlinear scaling function model different types of multi-way consensus dynamics. In particular, the Authors of Ref. [246] considered two conceptually different submodels of Eq. (68). The first, called multi-way consensus model I (or MCM-I) captures homophily, whereby given a node
i of a hyperedge
h, the interaction function is chosen as (70)
s i(x i,x j,x k,…)=g i(|〈x〉h−x i|),where
〈x〉h is the average state of the nodes belonging to
h, and
g i is a nonlinear function. In the context of opinion dynamics, Eq. (70) can be interpreted as a function of the distance between the opinion of
i and the average opinion of the group
h.

1. Download: Download high-res image (190KB)
2. Download: Download full-size image

Fig. 34. The numerical simulation of model (61)–(62) in the mean field approximation. The function
s(z) has been selected as
2 exp(λ z), and the initial configuration such that 80% of the nodes are initialized with 0 and 20% with 1. The average state (red solid line) is conserved only for
λ=0 (center panel). On the contrary, or
λ<0 (left panel), the average state shifts towards the majority following a reinforcement effect, whereas for
λ>0 (right panel), the opposite happens, with the dynamics drifting towards balance. In all panels, the dotted red line corresponds to the initial average state, whereas the dashed blue line corresponds the final state approximation proposed in Ref. [157].

Figure reprinted from Ref. [157], © 2020 by the American Physical Society.

1. Download: Download high-res image (122KB)
2. Download: Download full-size image

Fig. 35. Numerical simulation of model (61)–(62). The network is here is composed by two fully connected clusters,
A and
B, interconnected by 80 randomly placed triangles, so that a fraction
p∈[0,1] of them is directed towards
A. The function
s(z) has been selected as
2 exp(λ z), with
λ=−100 so that pairs of nodes with similar state exert a much stronger influence on other nodes. The final consensus value (left panel) and the rate of convergence (right panel), averaged over 20 simulations, are reported. The error bars correspond to one standard deviation.

Figure reprinted from Ref. [157], © 2020 by the American Physical Society.
The second selection of
s i captures instead conformity and constitutes what the Authors called multi-way consensus model II (or MCM-II): (71)
s i(x i,x j,x k,…)=g i(|〈x〉h/i−x j|),which can be interpreted as a function of the distance between the opinion of an individual
j with respect to the average opinion of the group except
i, whereby the influence of node
j within hyperedge
h depends on how close its opinion is with respect to that of the rest of the group.
Sahasrabuddhe et al. [246] illustrated the fundamental differences between these two submodels of (68) on a fully connected hypergraph of
N=10 nodes where they chose
g i(z)=exp(λ z) for all
i, with
λ being a real number, and set the initial state of
n 0 nodes to 0, and of the remaining
N−n 0 nodes to 1. As illustrated in Fig. 36, the two submodels evolve in opposite ways: the left panel shows that, for MCM-I, when
g i is monotonically increasing (
λ>0), the average node state drifts towards the opinion of the initial majority, whereas when
g i is monotonically decreasing (
λ<0) the opinions drift toward those of the initial minority. The exact opposite happens when MCM-II is considered, see the right panel of Fig. 36.

1. Download: Download high-res image (336KB)
2. Download: Download full-size image

Fig. 36. Numerical simulations of models MCM-I (left panel) and MCM-II (right panel) on a fully connected hypergraph.
N=10 nodes, and the function
g i(z) in Eqs. (70), (71) has been selected as
exp(λ z). The final consensus value as a function of the initial average node state is reported.

Figure reprinted from Ref. [246].

#### 4.3.1. Consensus on time-varying hypergraphs

Most of the works on higher-order interaction, with the exception of Chowdhari et al. [120] in the context of social contagion, have focused on static hypergraphs. However, in several applications, the topology of the interactions between the agents of the network evolves on a time scale that is comparable to that of the dynamics taking place on the nodes, and therefore cannot be neglected. Two alternative approaches, which we will review in the following, have been considered when the hypergraphs evolve in time. In the first, explored in Ref. [243], the Authors resort to temporal hypergraphs to describe a time-varying network topology. A second approach was followed by Schawe and Hernández [247], who considered instead a state-dependent evolution of the hypergraph in a generalization of the Deffuant–Weisbuch opinion dynamics model [250].

##### Dynamics on temporal hypergraphs.

A temporal hypergraph is a natural extension of the concept of temporal network over a standard graph (see, e.g., Ref. [144] and references therein). Neuhäuser et al. [243] extended the three-way consensus model (61) - (62) over a temporal hypergraph described by a sequence of adjacency tensors
A[1],A[2],…, each describing the topology of a hypergraph in a time-period of length
τ. As in the static case, these adjacency tensors can be associated to the corresponding state-dependent Laplacians, namely
ℒ[1],ℒ[2],…, therefore yielding to the following temporal 3CM: (72)
ẋ(t)=−ℒ[ℓ](t)x(t),(ℓ−1)τ≤t≤ℓ τ,ℓ=1,2,…Note that when the interaction function
s in Eq. (62) is linear,
ℒ[ℓ](t) becomes a constant Laplacian, model (72) can be solved analytically, and conditions for consensus can be given, see Ref. [251]. For general nonlinear interaction functions, an analytical examination of the properties of Eq. (72) becomes difficult, thus preventing the comparison of consensus dynamics over temporal hypergraphs with those over the corresponding time-aggregated static counterpart, obtained by averaging the adjacency tensors over time. For this reason, the Authors of Ref. [243] decided to rather compare the state dynamics on the time-aggregated hypergraphs with the average state dynamics on what they called a random temporal hypergraph, that is, a hypergraph whose topology, in each time-window, is drawn at random from a multiset of
r, possibly repeating Laplacian matrices. Neuhäuser et al. showed that, in the case of a linear interaction function, the temporal dynamics are slower than the time-aggregated dynamics. Additionally, when the interaction function is nonlinear, they also observed a first-mover advantage, whereby groups that are in the local majority in early-active hyperedges have a higher influence on the final consensus value.

##### Deffuant–weisbuch (DW) on hypergraphs.

Schawe and Hernández in Ref. [247] focused on the higher-order interaction that may take place on the classic Deffuant–Weisbuch model. In the DW mode, at every time step an edge from a backbone, static graph
G is selected, and the corresponding pair of agents update their opinion only if their are sufficiently close to each other. Therefore, the effective network topology in the DW model is a state-dependent graph
G e(t), since the interaction between any pair of agents
(i,j) can take place only if
(i,j) belongs to
G and their opinions are not too far from each other.
The DW interaction mechanism has been then generalized to hypergraphs in Ref. [247], which considered as backbone topology a hypergraph
H, and selected at every time step a hyperedge
h in
H, whose nodes would update their opinion only if the maximum and minimum opinion in the hyperedge
h were not too far away. This means that the effective network topology in the DW model is a state-dependent hypergraph
H e(t), whereby the interactions between a group of nodes belonging to a hyperedge
h∈H take place only if the extremal opinions in
h are sufficiently close. In this setting, at every time-step
k, a hyperedge
h∈E is selected, and the opinion of each member of
h is updated as (73)
ẋi(t+1)=1|h|∑j∈h x j(t),if h∈H e(t),x i(t),otherwise,for all
i∈h, where (74)
h∈H e(t)⟺max j∈h x j(t)−min j∈h x j(t)<ɛ,with
ɛ>0 being the confidence interval for the interaction to take place. This interaction rule models situations where an individual with a very different opinion in a subgroup may halt the discussion and prevent the achievement of an otherwise possible agreement between the other members of the subgroup. This implies that large groups have a low probability to achieve an agreement.
To explore the properties of such a modified DW model, the Authors of Ref. [247] performed numerical simulations on different types of
k-uniform hypergraphs. In all simulations, they considered 1000 independent realizations of the system for each value of
ɛ, which is varied between 0 and 0.6 with step 0.002.
Fig. 37 illustrates the effect of higher-order interactions on
k-uniform Erdős–Rényi (ER) hypergraphs. Compared with the classic DW model, which corresponds to the case
k=2, it can be noted that the average size of the largest opinion cluster is approximately naught for a larger range of confidence
ɛ, since, as discussed above, the interaction is more difficult the larger is the cardinality of the hyperedge. Furthermore, the phase transition from polarization to consensus, which is sharper for larger system size, disappears for
k>4. Schwaze and Hernández then showed that these numerical findings are confirmed on homogeneous Barabási–Albert hypergraph, with the difference that the polarization plateau disappears already for
k=3. The Authors suggest that this could be due to the presence of hubs that, belonging to a large number of hyperedges, are more likely to interact and form a large cluster compared to the other nodes. On the contrary, in regular, spatial hypergraphs, where the hyperedges are composed by the
k-closest neighbors, the behavior is qualitatively different, and presence of higher-order interactions induces a sharp, albeit continuous, transition from complete fragmentation to consensus.

##### Majority rule (MR) models.

In opinion dynamics, a large bulk of work has considered agents that have to focus on a finite number of options, and often the focus was restricted to the paradigmatic two-option case, see for instance the already mentioned voter model of Ref. [252], [253], [254], where at every step an agent adopts the opinion of one of its neighbors in a graph. In the majority rule (MR) model [255], instead, at each update event a group of
m≥3 agent is selected in the population, and the group simultaneously adopts the opinion of the majority (if
m is even, then the opinion is either selected randomly, or a bias for one opinion is introduced). Since the MR model introduces group, rather than pairwise interactions, the use of graphs to encode such interactions is questionable, and therefore Noonan and Lambiotte in Ref. [256] proposed to model and study the model on hypergraphs, and restricted their investigation on the case where groups of
m=3 agents interact at every time-step.

1. Download: Download high-res image (469KB)
2. Download: Download full-size image

Fig. 37. Size of the largest cluster. Mean size
〈S〉 (averaged over 1000 realizations) of the largest cluster as a function of the confidence
ɛ for model (73)–(74) on
k-uniform ER hypergraphs with expected mean degree
10 and for different number of nodes
N. Panels (a), (b), (c), (d), and (e) correspond to
k equal to 2, 3, 4, 5, and 6, respectively. The insets show the variances
V ar(S) of the largest cluster across the realizations.

Figure adapted from Ref. [247].
The Authors considered therefore a hypergraph
H with node set
V(H)={1,…,N} and whose hyperedges are only triangles. Denoting with
T(H) the set of its triangles, at every time-step a hyperedge
h∈T is selected, and the majority rule is applied to the nodes in
h. The MR dynamics were then derived in terms of Fokker–Planck equations, and the Authors showed that, for large number of agents, in the stochastic dynamics the drift term dominates over the diffusion term, thus making a deterministic approximation accurate.
The effectiveness of this approach was illustrated over different hypergraph topologies. In the mean-field case where
T contains all possible triangles, the Authors were able to accurately predict the density
ρ of agents having state 1 at time
t. In particular, defining
κ=(2 ρ(0)−1)2/(ρ(0)(1−ρ(0))) they obtained the following approximation (75)
ρ(t)=1 2 1−1−4/(4+4 κ exp(3 t)),ρ(0)<0.5,1 2 1+1−4/(4+4 κ exp(3 t)),ρ(0)>0.5,which, for
N=1 0 4, showed an excellent agreement with numerical simulations, with consensus achieved on 0 for
ρ(0)<0.5, or on 1 if
ρ(0)>0.5.
Then, Noonan and Lambiotte analyzed the dynamics over a tripartite hypergraph, where the agents were grouped in three sets
V a,
V b, and
V c, and the set of triangle
T(H) consisted of all possible triplets with one agent in each set, that is,
T(H)={{i,j,k}:i∈V a,j∈V b,k∈V c}. For large
N, they have been able to approximate the dynamics of the densities
ρ a,
ρ b and
ρ c of the agents in
V a,
V b, and
V c, respectively, that have opinion 1, and found two asymptotically stable fixed point (
[0,0,0] and
[1,1,1]) and a saddle point
[1/2,1/2,1/2]. Furthermore, they were able to accurately predict the exit probability
E(ρ 0), that is, the probability that all agents reach consensus on 1 if the initial conditions are given by
ρ 0, as (76)
E(ρ 0)=0,if ρ a(0)+ρ b(0)+ρ c(0)<3/2,1,if ρ a(0)+ρ b(0)+ρ c(0)<3/2,which is in excellent agreement with the numerical simulations illustrated in Fig. 38.
Furthermore, the Authors also extended the analysis to the case of modular hypergraphs and to hypergraphs with different node weighting. Interestingly, when considering a modular hypergraph with two fully connected communities linked by a number of intra-community triangle
c, they were able to predict the threshold on
c so that the MR model transitions from the coexistence of different opinion in different communities to global consensus.

1. Download: Download high-res image (246KB)
2. Download: Download full-size image

Fig. 38. Exit probabilities. Exit probability on a tripartite hypergraph of
N=900 as a function of
ρ a(0) and
ρ b(0), which are varied on a grid with uniform mesh size 0.01. Left, center, and right panels correspond to
ρ c(0) equal to 0.25, 0.5, and 0.75, respectively. The line
ρ a(0)+ρ b(0)=3/2−ρ c(0) is depicted in dashed red, and blue (below the dashed line) and yellow (above the dashed) lines correspond to consensus achieved onto 0 or 1, respectively. The blue dots represents the outlying simulations that are not in agreement with approximation Eq. (76).

Figure reprinted from Ref. [256], © 2021 by the American Physical Society.

#### 4.3.2. Consensus in simplicial complexes

In all consensus models presented so far, it is assumed that the dynamics takes place on the nodes of a hypergraph (or of a simplicial complex), whose topology only affects the interactions taking place between the nodes. An alternative take on higher-order consensus dynamics is to consider that a state variable is associated to each simplex (or hyperedge) and formulate flows on simplices of any dimension, and this is the path followed by DeVille in Ref. [227].
Given a simplicial complex
X of maximal dimension
D, and denoting
X d the set of its simplices of dimension
d, it can be defined, for all
d=0,1,…,D, the
R-vector space with basis given by the elements of
X d. For each
d, the following boundary map can be defined: (77)
∂d:C d(X)⟶C d−1(X),[v 0,…,v d]⟼∑l=0 d(−1)l[v 0,…,v ˆ i,…,v d],where
v i (a 0-simplex) is a vertex of
X, square brackets indicate that a simplex is ordered, and
[v 0,…,v ˆ i,…,v d] denotes
[v 0,…,v d]/v i. Furthermore, the Authors introduced a function
w:C d(X)→C d(X) for each
d, which associates the weight
w(F) to each simplex
F∈X. Picking a basis that corresponds to (any) ordering of the simplices in
X, one can then obtain
B d and
W d as the matrix representations of
∂d and
w, respectively.
Let
f:R→R be a function such that
f(0)=0 and
f′(0)>0, which induces a vector function
f:R n→R n such that the
ith component of
(x) is defined as
f(x i). Given some
ω d∈C d(X),
d=0,…,D, DeVille defined the consensus dynamics on the simplicial complex
X as (78a)
d d t θ 0=ω 0−W 0−1 B 1 W 1 f(B 1 T θ 0),(78b)
d d t θ d=ω d−W d−1 B d+1 W d+1 f(B d+1 T θ d)+B d T W d−1−1 f(B d W d θ d),1≤d<D,(78c)
d d t θ D=ω D−B D T W D−1−1 f(B D W D θ D), and called the
dth equation as the
d-simplex flow (the 0-simplex will be the “vertex flow” and the 1-simplex the “edge flow”). This choice for the simplicial flow Eqs. (78) ensures that (i) it reduces to standard consensus models when the simplicial complex is a graph, and (ii) the dynamics of each simplex are coupled only to those of the simplices to which it is related by boundary operations. Interestingly, with respect to point (i), DeVille showed that an appropriate choice of its parameters makes the vertex flow Eq. (78a) equivalent to the flow of the network dynamical system (79)
d d t x i=ζ i−∑j=1 N a i j f(x j−x i),on a graph
G={V,E}, where
a i j is the element
(i,j) of the symmetric adjacency matrix associated to
G. The vertex flow Eq. (78a) is equivalent to Eq. (79) when
B 1 is the incidence matrix of
G,
W 1 is a
|E|×|E| matrix with
(W 1)e e=a e corresponding to the weight of edge
e∈E,
W 0=I n, and
ω 0=[ζ 1;…;ζ n].
DeVille then characterized the solutions of the simplicial flow Eq. (78), and their stability properties. In particular, noting that the composition map
∂d∘∂d+1=0, it is natural to define the
dth homology group
H d(X)≔ker∂d/Im∂d+1 for
d=0,…,D−1. Solutions
x d∈ker B d+1 T∩ker B d W d are then called homological solutions since they form a basis for
H d(X). The Author of Ref. [227] showed that, linearizing around a homological solution, a negative semidefinite operator is obtained, but asymptotic stability cannot be concluded in general, since they might have multiple zero eigenvalues. However, when
f(0)=0 and
f′(0)>0, it can be shown that the set of homological solutions is asymptotically stable in the sense that any sufficiently small perturbation would relax to the subspace of homological solutions exponentially fast.
Nijholt and DeVille [245] then further explored the dynamics of systems defined over simplicial flows. In particular, they provide a method to build simplicial flows that are equivalent to any given dynamical systems, and provide examples of the Guckenheimer–Holmes cycle [257], and of a combination of the Lorenz system [258] and Sel’kov model [259]. Additionally, they discussed the impact of changing the orientations in simplicial complexes, and how the algebraic symmetries in the simplicial complex reflect onto the symmetries of the dynamical system defined on it.

### 4.4. Pinning control of hypergraphs

The main contributions of control theory in the study of network systems coupled through pairwise interactions have been the development of

- 1.
  Mathematical tools for steering the network toward a desired state (see, e.g., Ref. [260] and references therein). Among the wide range of control schemes, one of the most popular has been the so-called pinning control, which prescribes that a virtual node, the pinner, injects the control signal only to a (small) fraction of the network nodes [261], [262], [263], [264], [265], [266], [267].
- 2.
  Analytical tools to study controllability of arbitrary directed networks, and identify the set of driver nodes where a time-dependent control can be injected to guarantee full (or partial) controllability [268], [269], [270], [271], possibly minimizing the spent control energy [272], [273], [274], [275], [276].

Recent works have expanded the findings on pinning control and controllability in the presence of multi-body interactions, as we will discuss in the following. More specifically, we will review how directed hypergraphs have been used to model the case in which the pinner node cannot measure the individual state of the controlled nodes, and how the controllability analysis on linear networks on graphs has been expanded to the case of multi-body interactions that can be translated into polynomial dynamical systems [277].
In standard pinning control on networks, a virtual node, the pinner, exerts a control action on a limited fraction of the network nodes, denoted “pinned”. This action is proportional to the difference between the state of the pinner and that of the pinned node [261], [267]. The underlying assumption in this control strategy is that the pinner is capable of measuring the state of a given network node, and therefore this control signal can be suitably represented by a directed edge from the pinner to the controlled node. However, such a scheme is not always applicable as it may not be possible to feed back the state of the pinned node. For example, in microbial consortia the control designer only has available an aggregated measurement of the fluorescence of a cell group, whereby the fluorescence of a single cell (i.e., of a single node) cannot always be measured [278]. Similar challenges might be posed by constraints on the actuation that force the controller to inject the same control signal to more nodes (i.e., more cells). To account for these constraints, De Lellis et al. [279] modeled the interaction between coupled units through a directed hypergraph [280], using directed hyperedges from the pinner (the tail of the hyperedge) to a subset of nodes (the heads of the hyperedge) to describe the control signals. Namely, the Authors of Ref. [279] considered the dynamics of the nodes of a controlled hypergraph
ℋ c={V c,E c}, where
V c={ν 1,…,ν N} is the set of controlled nodes, and
E c={h 1,…,h M} the set of directed hyperedges of
ℋ c, with the
ith directed hyperedge
h i being an ordered pair
(T(h i),H(h i)) of (possibly empty) disjoint subsets of the hypergraph nodes. The ordered subsets
T(h i) and
H(h i) of
V are the set of tails and heads of the hyperedge
h i, with cardinality
|T(h i)| and
|H(h i)|, respectively, and such that
T(h i)∩H(h i)=0̸. De Lellis et al. associated a node
ν i∈V c to a state variable
x i∈R n, and to a hyperedge
h∈E c a tail state matrix
x h τ=[x t(h,1),…,x t(h,|T(h)|)] and a head state matrix
x h h=[x h(h,1),…,x h(h,|H(h)|)], with
x t(h,i) and
x h(h,j) being the state of the
ith tail and of the
jth head of
h, respectively. The node dynamics can be then written as (80)
ẋi=f(x i,t)+∑h∈E c⋅,i σ h(x h τ α h−x h h β h)+u i,where
f:R n×R≥0⟶R n is the continuous and differentiable vector field describing the individual dynamics,
σ h is the coupling gain associated to the hyperedge
h;
α h=[(α h)t(h,1),…,(α h)t(h,|T(h)|)]T and
β h=[(β h)h(h,1),…,(β h)h(h,|H(h)|)]T are the (ordered) vectors stacking the weights associated to the tails and heads of
h, respectively, such that
α h T 1|T(h)|=β h T 1|H(h)|=1; the control input
u i tries to steer the dynamics of the controlled nodes onto the trajectory of the pinner, defined by the following Cauchy problem: (81)
ẋs=f(x s,t),x s(0)=x s 0,where
x s∈R n is the state of the pinner.
Note that the coupling protocol in Eq. (80), called hyperdiffusive in Ref. [279], is such that each head of a hyperedge receives a signal which is proportional to the difference between a convex combination of the states of the tails and a convex combination of the states of the heads of the hyperedge. Such a protocol is synchronization noninvasive according to the definition given in Ref. [27]. To describe the higher-order links connecting the pinner with the controlled network Eq. (80), it is possible to consider an augmented hypergraph
ℋ={V,E}, of which
ℋ c is a proper sub-hypergraph. More specifically, the node set
V also includes the pinner, that is,
V=V c∪{v s}, and
E=E c∪⋃i=1 N E s,i,where
E s,i is the set of hyperedges having the pinner
v s as a tail, and node
v i as a head. Note that the state evolution of the pinner, by definition Eq. (81), is not influenced by those of the controlled nodes, and therefore it is a node that cannot be a head of any hyperedge, that is,
E is such that
s∉⋃h∈E H(h). In the context of hypergraphs, the standard definition of the set of pinned nodes
P as the nodes who are directly influenced by the pinner extends as
P={i∈V:E s,i≠0̸}, allowing to define the control action as (82)
u i=∑h∈E s,i k h(x h τ α h−x h h β h),i∈P,0,otherwise,where
k h is the control gain associated to the hyperedge
h, and where the pinner’s state
x s is one of the columns of the tail state matrix
x h τ. Examples of pinning control actions associated to hyperedges of size 3 are reported in Fig. 39.
As in the standard pinning control of graphs, the goal of the control action Eq. (82) is to steer the dynamics of the nodes in the controlled network toward the pinner’s trajectory. Defining the pinning error of node
i as
e i=x i−s, and the network pinning error as
e=[e 1;…;e N], the control objective is to drive
e(t) asymptotically to zero, that is, (83)
lim t→+∞‖e(t)‖=0.

1. Download: Download high-res image (53KB)
2. Download: Download full-size image

Fig. 39. Sample pinning hyperedges. Denoting with
k h the hyperedge weight, in panel (a), the inputs to nodes 1 and 2 would be
u 1=u 2=k h(x s−(β 1 x 1+β 2 x 2)), whereas in panel (b) the input to node 2 would be
u 2=k h((α 1 x 1+α s x s)−x 2).

#### 4.4.1. Equivalence with signed graphs

It can be noted that, given a hyperedge
h∈E, the following equality holds (84)
x h τ α h−x h h β h=∑j∈T(h)((α h)j(x j−x i))−∑j∈H(h)((β h)j(x j−x i)).This observation allows to rewrite network (80) - (81) over the hypergraph
ℋ as an equivalent network over a weighted signed graph
S=S(ℋ). Namely, Eqs. (80) - (81) are equivalent to (85)
żi=f(z i,t)+∑j=1 N+1 a i j(z j−z i),where
z=[z 1;…;z N+1], with
z i=x i for
i=1,…,N, and
z N+1=x s, stacks the states of the nodes of the controlled network and that of the pinner, and
a i j is the (off-diagonal) entry
i j of the adjacency matrix associated to
S, defined as (86)
a i j=∑h∈E c j,i(α h)j σ h−∑h∈E c⋅,{i,j}(β h)j σ h+I P(i)(∑h∈E{s,j},i(α h)j k h−∑h∈E c s,{i,j}(β h)j k h),where
I P(i)=1 if
i∈P, and zero otherwise, and where, given two node subsets
V 1,V 2⊆V,
E V 1,V 2 is the set of hyperedges having a tail in
V 1 and a head in
V 2, that is,
{h∈E:V 1⊆T(h)∧V 2⊆H(h)}; finally,
E⋅,V 1={h i∈E:V 1⊆H(h i)} is the set of hyperedges having all the nodes in
V 1 as heads. A sample 5-node hypergraph
ℋ and its equivalent signed graph
S(ℋ) is depicted in Fig. 40.

1. Download: Download high-res image (151KB)
2. Download: Download full-size image

Fig. 40. The signed graph. A sample hypergraph
ℋ is depicted in panel (a), and its equivalent associated signed graph
S(ℋ) in panel (b).

Figure reprinted from Ref. [279].

#### 4.4.2. Pinning control design

The equivalent representation (85) - (86) of the higher-order network (80) - (81) allowed the Authors of Ref. [279] to derive sufficient conditions for driving the state of the controlled node onto the trajectory of the pinner. Sorting the eigenvalues of the Laplacian matrix associated with the signed graph
S in ascending order of their real part, such that
ℜ(λ 1)≤…≤ℜ(λ N) (where
ℜ stands for the real part), a sufficient condition for global convergence of
e(t) to zero can be obtained when
ℜ(λ 2)>0 and the vector field
f fulfills the so-called QUAD assumption [279]. Interestingly, the condition
ℜ(λ 2)>0 on the interaction topology is not necessarily fulfilled if there is a hyperpath on
ℋ from the pinner to the rest of the nodes, different from what happens on directed graphs [264]. Therefore, the standard guideline for pinning control of digraphs, which prescribes to pin at least one node in each root strongly connected component of the controlled graph, cannot be trivially extended to graphs.
De Lellis et al. then built a heuristic algorithm to choose the pinning hyperedges, which starts by extracting from
S a (positively-signed) graph
S+ that only contains the positive edges of
S. The design procedure first ensures that at least one node in each root strongly connected component in
S+ is pinned through edges or hyperedges, which is a necessary condition for having
ℜ(λ 2)>0. The heuristic continues with the identification a critical negative edge in
S, say
h−, that, when added to
S+, would make the second smallest eigenvalue associated to the Laplacian matrix cross the imaginary axes. Then, the algorithm proceeds by adding a hyperedge with heads in one of the endpoints, say
i, of
h− and only in the strongly connected components to which
i belongs. The algorithm terminates when there are no more critical edges in
S, and therefore the final selection of pinning hyperedges guarantees that
ℜ(λ 2)>0. Fig. 41 illustrates the effectiveness of the heuristic on a sample hypergraph.

1. Download: Download high-res image (311KB)
2. Download: Download full-size image

Fig. 41. Pinning control. Application of the pinning control design heuristic from Ref. [279] to a network of Chua’s circuit coupled through the hypergraph
ℋ depicted in panel (a), where the pinner is node 9, and the we report the Laplacian matrix
L of the associated signed graph
S(ℋ). Panel (b) and (c) shows the dynamics of the error norm for each of the controlled node before and after the application of the pinning selection algorithm.

Figure adapted from Ref. [279].

### 4.5. Controllability of hypergraphs

In linear networks over digraphs, the classic controllability problem tackled in the literature studies the following dynamics (87)
ẋ=A x+B u,where
x=[x 1;…;x n]∈R n is the network state, with
x i∈R being the (scalar) state of node
i,
A∈R n×n is the adjacency matrix of a graph,
B∈R n×m is the input matrix, and
u∈R m is the control input. In this context, starting from the classic Kalman criterion [281], the problem of selecting the driver nodes where the input should be injected to achieve complete or partial controllability of the network has been widely investigated in the literature (see, e.g., Refs. [268], [269], [270], [271], [274] and references therein).
Chen et al. in Ref. [277] generalized the controllability problem over network (87) to the case of interactions taking place over an undirected hypergraph
ℋ=(V,E). Specifically, they considered a
k-uniform hypergraph, whose adjacency tensor is a
kth order
n-dimensional supersymmetric tensor, whereby its element
j 1 j 2…j k, defined as (88)
A j 1 j 2…j k=1(k−1)!,if {j 1,j 2,…,j k}∈E,0,otherwise,is invariant under any permutation of the indices. The tensor multiplication of
A along mode
1≤p≤k for a vector
v∈R n is a
k−1 order
n-dimensional sensor, denoted
A×p v, defined as
(A×p v)j 1 j…j p−1 j p+1…j k=∑j p=1 n T j 1…j p…j k v j p.Given
k vectors
v 1,…,v k, the tensor multiplication can be extended as
A v 1 v 2…v k=A×1 v 1×2 v 2×3…×k v k.Using tensor algebra, the Authors of Ref. [277] generalized system (87) as (89)
ẋ=A x k−1+∑j=1 m b j u j,where
b j is the
jth column of
B, and (90)
A x k−1=A×1 x×2 x…×k−1 x.We note that, for
u=0, Eq. (89) would describe a homogeneous polynomial system, whereby
A x k−1∈R n is a homogeneous polynomial of degree
k−1. Therefore, the higher order interactions in the hypergraph translate into multiplicative terms in the interaction model proposed by Chen et al., which were motivated by the nonlinear dynamical systems that model protein–protein interactions (see, e.g., Ref. [282]).
Leveraging classic results from controllability of polynomial systems [283], Chen et al. derived a condition for controllability of system (89). Letting
C 0 being the linear span of
{b 1,…,b m}, for a generic integer
q,
C q can be recursively defined as the linear span of
C q−1∪{A v 1 v 2…v k−1|v l∈C q−1}.The Authors of Ref. [277] proved that, for an even
k, the polynomial system (89) is strongly controllable, in the sense that it can be driven from any initial state to any target state in any finite time, if the
C(A,B)=∪q≥0 C q=R n. Note that, for linear networks, that is, when
ℋ is a graph, this condition reduces to the well-known Kalman condition. Using this result for assessing controllability can be computationally prohibitive, and therefore the Authors proposed an efficient approach to compute a reduced controllability matrix, denoted
C r such that, if
rank(C r)=n, system (89) is strongly controllable.
When
k is odd, that is, when we consider odd uniform hypergraphs, conditions for controllability have not been derived yet. Indeed, such a hypergraph would give rise to a polynomial system of even degree, whose controllability is an open problem in the literature [284]. Nonetheless, using the results of Ref. [285], Chen et al. showed that the same condition used for controllability for even
k, that is,
C(A,B)=R n, can be used for an odd
k to guarantee that system (89) is accessible, which is a weaker form of controllability that implies that, for any initial state
x 0, the reachable set in a positive time
T contains a non-empty open set.
Based on their results, Chen et al. then derived a heuristic algorithm to determine the minimal selection of driver nodes to enforce controllability (accessibility in the case of
k odd): the driver nodes were chosen based on the maximum change in the ranking of the reduced controllability matrix
C r. Fig. 42 reports the driver node selection for select 4-uniform
r-hyperchains,
r-hyperrings and
r-hyperstars, where
r≤k−1 is the number of nodes common between two intersecting hyperedges.2 Interestingly, the Authors point out that the selection of the driver nodes is linked to the degree distribution, whereby numerical evidence shows that nodes with the highest degree are always picked as drivers.

1. Download: Download high-res image (494KB)
2. Download: Download full-size image

Fig. 42. Driver node selection for different types of hypergraphs. Panels A, B, and C correspond to 4-uniform 1-hyperchains, 2-hyperchains, and 3-hyperchains, respectively. Panels D, E, and F correspond to 4-uniform 1-hyperrings, 2-hyperrings, and 3-hyperrings, respectively. Panels G, H, and I correspond to 4-uniform 1-hyperstars, 2-hyperstars, and 3-hyperstars, respectively. The selected driver nodes are identified by arrows, with cyan arrows corresponding to the nodes with the highest degree, which are always selected as drivers.

Figure based on Ref. [277].

## 5. Conclusions and perspectives

At the end of this long journey, we dare to make a few concluding remarks with the aim to identify questions that remain still open to future progresses, and to point this way the reader towards some issues and problems that we consider to be of key relevance, and that we believe that would most likely attract soon the attention of scientists in the area.
Our Chapter 2 has concentrated on giving the essential mathematical concepts, and in describing some measures and properties that have been used in the rest of the report. In fact, a full extension to hypergraphs of the measures and definitions which are currently used for describing the structural properties of classical networks is a real challenge for applied mathematicians. The extent of this difficulty can be easily understood by considering that, while classic networks are described in terms of matrices, hypergraphs requires a tensor representation, and the transition from matrix to tensor algebra is a passage which, roughly speaking, corresponds to the opening of a Pandora’s box.
It is unclear, for instance, whether the proper approach to describe networks with interactions at different orders is still that of defining scalar measures (duly weighted to account for the various orders of interplay), or instead to adopt vectorial quantities (having different components, one per each order of interaction) that may carry a more refined information on the overall structure of the hypergraph. This is the case of Ref. [46], where it was proved that a vector centrality offers a way to differentiate nodes in different rankings, and therefore to distinguish the role that a given node may play in the system at a given order of interaction.
A second absolutely relevant question is describing group interactions of heterogeneous nature. A clear issue is the problem of directionality: in classical networks directionality simply means an arrow indicating the direction where the pairwise interaction occurs. However, a definition of asymmetric hyperlinks may be not unique. Let us consider, for example, the case of (scientific) collaboration networks. We already saw that this case can be mapped into a higher order network, where nodes represent scientists and each scientific paper is an hyperlink among the co-Authors of the Manuscript [13], [286]. However, the role played by each Author in the development of a Manuscript is not the same: in many disciplines, the first and last Authors are the scientists who mostly contributed to the Manuscript, as well as in multidisciplinary collaborations different kind of scientists (such as theoretical researchers and experimental scientists) may play a different role in the development of the collaboration. This would require a completely new formalism which would introduce “structured” hyperlinks, i.e.. new mathematical objects which are not simply the entries of a tensor, but are endowed with internal properties that would describe the complexity of the interactions among the nodes contained in the hyperlink. And a second Pandora’s box opens up.
Chapter 3 of our report is an account on different processes that may take place among units interacting beyond a pairwise level, from epidemic spreading, to models of opinion formation, social contagion, and games.
Here, several open questions exist that, we are sure, will soon attract the interest of fellow researchers. A fundamental point is the study of spreading models that mix contact spreading and spatial diffusion, as well as the study of models of multiple diseases spreading, with nontrivial interactions between them. In turn, this would lead to a full understanding on the effect of evolving higher-order structures on the spreading processes. As a particularly important case, the study of how reactive structural adaptivity affects epidemic spreading.
Another missing point is a full equilibrium thermodynamic theory for contagion models, as well as the determination of the rigorous mathematical conditions for the applicability of the numerous approximations that have been made in the literature (in other words, what are the actual underlying or implicit assumptions used, and in what systems would they fail?). At the same time, a pressing issue is the estimation of the errors caused by the presence of structural and dynamical correlations when using mean-field methods to study epidemics and related strategies to stop or slow down the spreading.
Finally, we mention the need of developing (structural) renormalization techniques for contact and dynamical processes on hypergraphs, which, in turn, would help with the inclusion of correlation effects, which is currently numerically difficult for networks with hyperedges of high dimension, as well as the study of the effects of stochasticity on spreading processes.
Chapter 4 is the longest of this report, and accounts for the emergence of collective dynamics in hypernetworks and control of hypergraphs.
In particular, the first part of the Chapter is dedicated to synchronization of networked dynamical units, i.e., to the fact that an ensemble of dynamical systems may give rise to a macroscopic synchronous motion as a result of the specific network structure of interactions among them.
Synchronization of networked units is a behavior observed far and wide in natural and man made systems, and correspond to the fact that, as a control parameter (typically the coupling strength in each link or hyperlink) is increased, the ensemble passes from a fully disordered state (at small coupling values, in which each unit evolves independently on the others) to a fully ordered one characterized by the fact that all units follow the same trajectory.
The transition to synchronization of identical units can be either abrupt and discontinuous (and, in this case, the term used in the literature is explosive synchronization) or second-order like i.e., continuous and progressive. This latter state is called cluster synchronization and corresponds, for classical networks, to the emergence of structured states where the ensemble splits into different subsets (the clusters) each one evolving in unison. Moreover, it is known that the underlying symmetries of a network are responsible for the way nodes split in functional clusters during cluster synchronization. A clear challenge is, therefore, describing the transition to synchronization in the case of hypernetworks, where the concept of clusters and symmetries need still to be properly defined.
When the units are non identical, as in the case of the Kuramoto model, we have reviewed a bulk of studies which have pointed out that higher order interactions actually introduce a very rich synchronization scenario, including states that would be inherently prevented for pairwise interactions. Almost all analytical studies, however, are so far confined to the case of global coupling, and the extension of analytical approaches to structured hypernetworks is certainly a problem to be faced within the next years.
Another crucial problem which needs thorough investigation is how inclusion of higher-order interaction facilitates first-order transition to synchronization in non identical systems. Specifically, how much higher-order interactions are sufficient to cause a first-order-like transition in the system, and which kind of distribution of higher-order interactions is supportive or detrimental for causing abrupt synchronization in hypernetworks.
Chapter 4 contains moreover a thorough review on control and controllability of hypernetworks.
The full exploitation of the hypergraph formalism in control theory is yet far to come. In the presence of diffusive-like higher order coupling, for instance, the analogy with signed graphs could be further explored to optimally select the units to be pinned in pinning control problems, minimize the control energy required to achieve the desired coordination task in the network, and to identify the driver nodes guaranteeing controllability of the overall system.
On the other hand, a theoretical framework that uses directed hypergraphs to translate the actuation and sensing constraints on networked dynamical systems is far from having been rigorously developed. In particular, control strategies should be developed to tackle coordination problem for more general higher-order interaction of the type described in Ref. [27], where the protocol may not be hyperdiffusive as that considered in Ref. [279].
Another challenging problem in control of networked systems is to steer the trajectory of a co-evolving network, when dynamics take place both at the nodes and at the edges. In the presence of higher-order interactions, dynamics may also take place on higher-order structures such as triangles and tetrahedra, see e.g., Ref. [227]. In this framework, therefore, an open challenge is that of designing control inputs that would guarantee the emergence of a desired behavior in the overall higher-order system.
As for possible relevant applications, we here briefly quote some circumstances where we believe that the formalism of higher-order networks can be used successfully in future years, especially in power grids, neuroscience, social science and linguistics.
Coupled Kuramoto oscillators with inertia [287] provide a model that successfully describes the dynamics of power-grids [288], [289], [290]. As we have seen in the previous Chapter, a large body of literature exists on Kuramoto oscillators interplaying (without inertia) on hypergraphs, and recently an analytical framework has been developed which incorporates inertia in Kuramoto oscillators coupled interacting on simplicial complex [291]
The brain is another example where the importance of higher-order interactions has been revealed for information processing [10], and the existence of its various functional responses [292].
Modern linguistics is characterized by its innovative approach of considering a language as a system or a complex network, this way developing fresh and efficient quantitative methods for its study. In particular, the higher order network formalism entitles one to build new models of linguistic that consider sentences or paragraphs as groups or collections of certain words, and makes it possible to capture and analyze information on the mesoscopic structure of texts (relationships between words, sentences, paragraphs, chapters and texts). Additionally, new tools for detecting plagiarism and obtaining elements that allow characterizing the style of an author can be introduced.
On the other hand, considering each hyper-arista as a property or a characteristic that a node may or may not have, or even as an event or an affair to which a node may or may not participate, applicability of high-order networks in the field of social network models is beyond any doubt.
Finally, let us conclude this report by remarking, once again, that Network Science is one of the hottest and most successful fields of modern research. Given the enormous success and the gigantic number of applications that complex networks have already encountered in almost all realms of nonlinear science, it is easy therefore to see that studies on hypergraphs and other higher order networks constitute the next step that network scientists will make within the next decade. At the same time, we believe that the years to come will possibly bring to light other phenomena whose description will instead require to mix continuous and distributed approaches, and for which therefore the only representation of the system as a single network could be too short a limitation.

## Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Acknowledgments

The Authors would like to gratefully acknowledge D. Aleja, R. E. Amritkar, A. Arenas, B. Barzel, I. Belykh, V. Belykh, G. Bianconi, E. M. Bollt, J. Bragard, J. M. Buldú, F. Bullo, J. Burguete, A. Buscarino, A. Cardillo, T. Carroll, H. Cerdeira, M. Clerc, G. Contreras-Aso, S. K. Dana, F. Della Rossa, M. Di Bernardo, A. Diaz Guilera, T. Di Matteo, R. D’Souza, E. Estrada, S. Focardi, U. Feudel, L. Fortuna, M. Frasca, N. Frolov, L.V. Gambuzza, J. Garcia-Ojalvo, A. Garcimartin, D. Ghosh, G. Giacomelli, J. Gómez-Gardeñes, N. Gupta, S. Havlin, C. Hens, J. Hizanidis, P. Hövel, T. House, A. E. Hramov, I. Iacopini, S. Jafari, M. Jusup, B. Kahng, T. Kapitaniak, P. Khanra, K. Kovalenko, C. Kühn, S. Kurkin, J. Kurths, M. Lakshmanan, S. Lepri, D. Liuzza, V. Latora, F. Lo Iudice, S. Majhi, C. Masoller, H. Mancini, D. Maza, R. Meucci, L. Minati, G. Mindlin, Y. Moreno, D. Musatov, S. Olmi, I. Omelchenko, O. Omel’chenko, G. V. Osipov, P. Pal, P. Parmananda, L. M. Pecora, M. Perc, N. Perra, A. N. Pisarchik, A. Politi, M. Porfiri, A. Prasad, E. Primo, I. Procaccia, S. Rakshit, A. Raigorovskii, R. Ramaswamy, P.L. Ramazza, A. Rizzo, R.G. Rojas, I. Samoylenko, E. Schöll, M. D. Shrimali, P. Simon, S. Sinha, F. Sorrentino, B. Tadić, A. Torcini, J. Tredicce, E. Vasilyeva, S. Yanchuk, A. Zakharova, I. Zuriguel, for the many discussions on the various subjects covered in this review article, which greatly inspired our writing.
S.B. acknowledges support from the Italian Ministry of Foreign Affairs and International Cooperation (grant number ZA23GR02). P.D. acknowledges support from the Research Project PRIN 2017 “Advanced Network Control of Future Smart Grids” funded by the Italian Ministry of University and Research (2020–2023). C.I.D.G. acknowledges support from UKRI under Future Leaders Fellowship grant number MR/T020652/1. K.A.B. acknowledges support from 2023/00004/001M2978, 2023/00005/016M3033 (URJC Grants). R.C. and M.R. acknowledge support from projects n.PGC2018-101625-B-I00 (Spanish Science Ministry, AEI/FEDER, UE), 2019/00034/001M1967, 2023/00004/001M2978, 2023/00005/016M3033 (URJC Grants). S.J. acknowledges support from SERB Power grant SPF/2021/000136 and the computational facilities received from the Department of Science and Technology (DST), Government of India , under the FIST scheme (Grant No. SR/FST/PSI-225/2016).

## References

1. [1]
   Koschmieder E.L.
   Bénard Cells and Taylor Vortices
   Cambridge University Press (1993)
   Google Scholar
2. [2]
   Getling A.V.
   Rayleigh-Bénard Convection: Structures and Dynamics, Vol. 11
   World Scientific (1998)
   Google Scholar
3. [3]
   Arecchi F.T., Boccaletti S., Ramazza P.
   Pattern formation and competition in nonlinear optics
   Phys. Rep., 318 (1–2) (1999), pp. 1-83
   View PDF View article View in Scopus
   Google Scholar
4. [4]
   Boccaletti S., Latora V., Moreno Y., Chavez M., Hwang D.-U.
   Complex networks: Structure and dynamics
   Phys. Rep., 424 (4–5) (2006), pp. 175-308
   View PDF View article View in Scopus
   Google Scholar
5. [5]
   Boccaletti S., Bianconi G., Criado R., Del Genio C.I., Gómez-Gardenes J., Romance M., Sendina-Nadal I., Wang Z., Zanin M.
   The structure and dynamics of multilayer networks
   Phys. Rep., 544 (1) (2014), pp. 1-122
   View PDF View article View in Scopus
   Google Scholar
6. [6]
   Latora V., Nicosia V., Russo G.
   Complex Networks: Principles, Methods and Applications
   Cambridge University Press (2017)
   Google Scholar
7. [7]
   Petri G., Expert P., Turkheimer F., Carhart-Harris R., Nutt D., Hellyer P.J., Vaccarino F.
   Homological Scaffolds of brain functional networks
   J. R. Soc. Interface, 11 (101) (2014), Article 20140873
   Crossref
   View in Scopus
   Google Scholar
8. [8]
   Lord L.-D., Expert P., Fernandes H.M., Petri G., Van Hartevelt T.J., Vaccarino F., Deco G., Turkheimer F., Kringelbach M.L.
   Insights into brain architectures from the homological Scaffolds of functional connectivity networks
   Front. Syst. Neurosci., 10 (2016), p. 85
   View in Scopus
   Google Scholar
9. [9]
   Lee H., Kang H., Chung M.K., Kim B.-N., Lee D.S.
   Persistent brain network homology from the perspective of dendrogram
   IEEE Trans. Med. Imaging., 31 (12) (2012), pp. 2267-2277
   View in Scopus
   Google Scholar
10. [10]
    Sizemore A.E., Giusti C., Kahn A., Vettel J.M., Betzel R.F., Bassett D.S.
    Cliques and cavities in the human connectome
    J. Comp. Neurosci., 44 (1) (2018), pp. 115-145
    Crossref
    View in Scopus
    Google Scholar
11. [11]
    Estrada E., Ross G.J.
    Centralities in simplicial complexes. Applications to protein interaction networks
    J. Theoret. Biol., 438 (2018), pp. 46-60
    View PDF View article View in Scopus
    Google Scholar
12. [12]
    Sizemore A.E., Karuza E.A., Giusti C., Bassett D.S.
    Knowledge gaps in the early growth of semantic networks
    (2017)
    arXiv preprint arXiv:1709.00133
    Google Scholar
13. [13]
    Patania A., Petri G., Vaccarino F.
    The shape of collaborations
    EPJ Data Sci., 6 (1) (2017), p. 18
    View in Scopus
    Google Scholar
14. [14]
    Petri G., Scolamiero M., Donato I., Vaccarino F.
    Topological strata of weighted complex networks
    PLoS One, 8 (6) (2013), Article e66506
    Crossref
    View in Scopus
    Google Scholar
15. [15]
    Aleksandrov P.S.
    Combinatorial Topology, Vol. 1
    Courier Corporation (1998)
    Google Scholar
16. [16]
    Carlsson G.
    Topology and data
    Bull. Amer. Math. Soc., 46 (2) (2009), pp. 255-308
    View in Scopus
    Google Scholar
17. [17]
    Salnikov V., Cassese D., Lambiotte R.
    Simplicial complexes and complex systems
    Eur. J. Phys., 40 (1) (2018), Article 014001
    Google Scholar
18. [18]
    Costa A., Farber M.
    Random simplicial complexes
    Configuration Spaces, Springer, Berlin (2016), pp. 129-153
    Crossref
    View in Scopus
    Google Scholar
19. [19]
    Courtney O.T., Bianconi G.
    Generalized network structures: The configuration model and the canonical ensemble of simplicial complexes
    Phys. Rev. E, 93 (6) (2016), Article 062311
    Crossref
    View in Scopus
    Google Scholar
20. [20]
    Bianconi G., Rahmede C.
    Complex quantum network manifolds in dimension d
    > 2 are scale-free
    Sci. Rep., 5 (2015), p. 13979
    View in Scopus
    Google Scholar
21. [21]
    Kovalenko K., Sendiña-Nadal I., Khalil N., Dainiak A., Musatov D., Raigorodskii A.M., Alfaro-Bittner K., Barzel B., Boccaletti S.
    Growing scale-free simplices
    Commun. Phys., 4 (1) (2021), pp. 1-9
    Crossref
    Google Scholar
22. [22]
    Petri G., Barrat A.
    Simplicial activity driven model
    Phys. Rev. Lett., 121 (22) (2018), Article 228301
    View in Scopus
    Google Scholar
23. [23]
    Iacopini I., Petri G., Barrat A., Latora V.
    Simplicial models of social contagion
    Nature Commun., 10 (1) (2019), pp. 1-9
    Crossref
    Google Scholar
24. [24]
    Li Z., Deng Z., Han Z., Alfaro-Bittner K., Barzel B., Boccaletti S.
    Contagion in simplicial complexes
    Chaos Solitons Fractals, 152 (2021), Article 111307
    View PDF View article View in Scopus
    Google Scholar
25. [25]
    Guo H., Jia D., Sendiña-Nadal I., Zhang M., Wang Z., Li X., Alfaro-Bittner K., Moreno Y., Boccaletti S.
    Evolutionary games on simplicial complexes
    Chaos Solitons Fractals, 150 (2021), Article 111103
    View PDF View article View in Scopus
    Google Scholar
26. [26]
    Skardal P.S., Arenas A.
    Abrupt desynchronization and extensive multistability in globally coupled oscillator simplexes
    Phys. Rev. Lett., 122 (2019), Article 248301
    View in Scopus
    Google Scholar
27. [27]
    Gambuzza L.V., Di Patti F., Gallo L., Lepri S., Romance M., Criado R., Frasca M., Latora V., Boccaletti S.
    Stability of synchronization in simplicial complexes
    Nature Commun., 12 (1) (2021), pp. 1-13
    Google Scholar
28. [28]
    Battiston F., Cencetti G., Iacopini I., Latora V., Lucas M., Patania A., Young J.-G., Petri G.
    Networks beyond pairwise interactions: Structure and dynamics
    Phys. Rep., 874 (2020), pp. 1-92
    View PDF View article View in Scopus
    Google Scholar
29. [29]
    Berge C.
    Hypergraphs: Combinatorics of Finite Sets, Vol. 45
    Elsevier (1984)
    Google Scholar
30. [30]
    Bretto A.
    Hypergraph Theory. An Introduction. Mathematical Engineering
    Springer, Cham (2013)
    Google Scholar
31. [31]
    Voloshin V.I.
    Introduction to Graph and Hypergraph Theory
    Nova Kroshka Books, UK (2013)
    Google Scholar
32. [32]
    Bianconi G.
    Higher-Order Networks
    Cambridge University Press (2021)
    Google Scholar
33. [33]
    Battiston F., Petri G.
    Higher-order systems
    (2022)
    Google Scholar
34. [34]
    Zhang X., Boccaletti S., Guan S., Liu Z.
    Explosive synchronization in adaptive and multilayer networks
    Phys. Rev. Lett., 114 (3) (2015), Article 038701
    Crossref
    View in Scopus
    Google Scholar
35. [35]
    Estrada E., Rodriguez-Velazquez J.A.
    Complex networks as hypergraphs
    (2005)
    arxiv preprint physics/0505137
    Google Scholar
36. [36]
    Johnson J.
    Hypernetworks in the Science of Complex Systems, Vol. 3
    World Scientific (2013)
    Google Scholar
37. [37]
    Criado-Alonso A., Aleja D., Romance M., Criado R.
    Derivative of a hypergraph as a tool for linguistic pattern analysis
    Chaos Solitons Fractals, 163 (2022), Article 112604
    View PDF View article View in Scopus
    Google Scholar
38. [38]
    Criado-Alonso Á., Aleja D., Romance M., Criado R.
    A new insight into linguistic pattern analysis based on multilayer hypergraphs for the automatic extraction of text summaries
    Math. Methods Appl. Sci. (2023)
    Google Scholar
39. [39]
    Tyshkevich R.I., Zverovich V.E.
    Line hypergraphs: A survey
    Acta Appl. Math., 52 (1) (1998), pp. 209-222
    Crossref
    View in Scopus
    Google Scholar
40. [40]
    Eilenberg S., Zilber J.A.
    Semi-simplicial complexes and singular homology
    Ann. of Math. (1950), pp. 499-513
    Crossref
    Google Scholar
41. [41]
    Massey W.S.
    A Basic Course in Algebraic Topology, Vol. 127
    Springer (2019)
    Google Scholar
42. [42]
    Rotman J.J.
    An Introduction to Algebraic Topology, Vol. 119
    Springer Science & Business Media (2013)
    Google Scholar
43. [43]
    Zhou D., Huang J., Schölkopf B.
    Learning with hypergraphs: Clustering, classification, and embedding
    Adv. Neural Inf. Process. Syst., 19 (2006)
    Google Scholar
44. [44]
    Courtney O.T., Bianconi G.
    Generalized network structures: The configuration model and the canonical ensemble of simplicial complexes
    Phys. Rev. E, 93 (2016), Article 062311
    Crossref
    View in Scopus
    Google Scholar
45. [45]
    Gambuzza L., Di Patti F., Gallo L., Lepri S., Romance M., Herrero R., Frasca M., Latora V., Boccaletti S.
    The Master Stability Function for Synchronization in Simplicial Complexes
    Springer International Publishing (2022), pp. 249-267
    Crossref
    View in Scopus
    Google Scholar
46. [46]
    Kovalenko K., Romance M., Vasilyeva E., Aleja D., Criado R., Musatov D., Raigorodskii A., Flores J., Samoylenko I., Alfaro-Bittner K., et al.
    Vector centrality in hypergraphs
    Chaos Solitons Fractals, 162 (2022), Article 112397
    View PDF View article View in Scopus
    Google Scholar
47. [47]
    Whitney H.
    Congruent graphs and the connectivity of graphs
    Amer. J. Math., 54 (1) (1932), pp. 150-168
    Crossref
    Google Scholar
48. [48]
    Hemminger R., Beineke L.
    Lowell W.B., Wilson R.J. (Eds.), Line Graphs and Line Digraphs, Selected Topics in Graph Theory, Academic Press, New York (1978)
    Google Scholar
49. [49]
    Bermond J.-C., Heydemann M.-C., Sotteau D.
    Line graphs of hypergraphs I
    Discrete Math., 18 (3) (1977), pp. 235-241
    View PDF View article View in Scopus
    Google Scholar
50. [50]
    Heydemann M.-C., Sotteau D.
    Line graphs of hypergraphs II
    Colloq. Math. Soc. J. Bolyai, vol. 18 (1976), pp. 567-582
    Google Scholar
51. [51]
    Bagga J.
    Old and new generalizations of line graphs
    Int. J. Math. Math. Sci., 2004 (29) (2004), pp. 1509-1521
    View in Scopus
    Google Scholar
52. [52]
    Benson A.R.
    Three hypergraph eigenvector centralities
    SIAM J. Math. Data Sci., 1 (2) (2019), pp. 293-312
    Crossref
    View in Scopus
    Google Scholar
53. [53]
    Criado R., Flores J., García del Amo A., Romance M.
    Centralities of a network and its line graph: An analytical comparison by means of their irregularity
    Int. J. Comput. Math., 91 (2) (2014), pp. 304-314
    Crossref
    View in Scopus
    Google Scholar
54. [54]
    Criado R., Flores J., Garcia del Amo A., Romance M., Barrena E., Mesa J.A.
    Line graphs for a multiplex network
    Chaos, 26 (6) (2016), Article 065309
    View in Scopus
    Google Scholar
55. [55]
    Evans T.S., Lambiotte R.
    Line graphs, link partitions, and overlapping communities
    Phys. Rev. E, 80 (1) (2009), Article 016105
    Crossref
    View in Scopus
    Google Scholar
56. [56]
    Evans T.S., Lambiotte R.
    Line graphs of weighted networks for overlapping communities
    Eur. Phys. J. B, 77 (2) (2010), pp. 265-272
    Crossref
    View in Scopus
    Google Scholar
57. [57]
    Naik R.N.
    On intersection graphs of graphs and hypergraphs: A survey
    (2018)
    arXiv preprint arXiv:1809.08472
    Google Scholar
58. [58]
    Arrigo F., Higham D.J., Tudisco F.
    A framework for second-order eigenvector centralities and clustering coefficients
    Proc. R. Soc. Lond. Ser. A Math. Phys. Eng. Sci., 476 (2236) (2020), Article 20190724
    View in Scopus
    Google Scholar
59. [59]
    Tudisco F., Higham D.J.
    Node and edge nonlinear eigenvector centrality for hypergraphs
    Commun. Phys., 4 (1) (2021), pp. 1-10
    View in Scopus
    Google Scholar
60. [60]
    Lemmens B., Nussbaum R.
    Nonlinear Perron-Frobenius Theory, Vol. 189
    Cambridge University Press (2012)
    Google Scholar
61. [61]
    Jaccard P.
    Distribution de la flore alpine dans le bassin des dranses et dans quelques régions voisines
    Bull. Soc. Vaudoise Sci. Nat., 37 (1901), pp. 241-272
    Google Scholar
62. [62]
    Brusco M., Cradit J.D., Steinley D.
    A comparison of 71 binary similarity coefficients: The effect of base rates
    PLoS One, 16 (4) (2021), Article e0247751
    Crossref
    View in Scopus
    Google Scholar
63. [63]
    Costa L.d.F.
    Further generalizations of the Jaccard I ndex
    (2021)
    arXiv preprint arXiv:2110.09619
    Google Scholar
64. [64]
    Costa L.d.F.
    On the effects of text preprocessing on paragraph similarity networks
    (2022)
    Google Scholar
65. [65]
    Costa L.d.F.
    On similarity
    Physica A, 599 (2022), Article 127456
    View PDF View article View in Scopus
    Google Scholar
66. [66]
    Costa L.d.F.
    Coincidence complex networks
    J. Phys.: Complexity, 3 (1) (2022), Article 015012
    Google Scholar
67. [67]
    Talaga S., Nowak A.
    Structural measures of similarity and complementarity in complex networks
    Sci. Rep., 12 (1) (2022), pp. 1-18
    Google Scholar
68. [68]
    Vijaymeena M., Kavitha K.
    A survey on similarity measures in text mining
    Mach. Learn. Appl.: Int. J., 3 (2) (2016), pp. 19-28
    Google Scholar
69. [69]
    Bollobás B.
    Random graphs
    Modern Graph Theory, Springer (1998), pp. 215-252
    Crossref
    Google Scholar
70. [70]
    Chung F.R.
    Spectral Graph Theory, Vol. 92
    American Mathematical Soc. (1997)
    Google Scholar
71. [71]
    Mohar B.
    The Laplacian Spectrum of Graphs. Graph Theory, Combinatorics and Applications, Vol. 2 (1988)
    Wiley-Intersci. Publ., Wiley, New York (1991), pp. 871-898
    Google Scholar
72. [72]
    Schaub M.T., Benson A.R., Horn P., Lippner G., Jadbabaie A.
    Random walks on simplicial complexes and the normalized Hodge 1-Laplacian
    SIAM Rev., 62 (2) (2020), pp. 353-391
    Crossref
    View in Scopus
    Google Scholar
73. [73]
    Mulas R., Horak D., Jost J.
    Graphs, simplicial complexes and hypergraphs: Spectral theory and topology
    Higher-Order Systems, Springer (2022), pp. 1-58
    Crossref
    View in Scopus
    Google Scholar
74. [74]
    Rodriguez J.A.
    On the Laplacian eigenvalues and metric parameters of hypergraphs
    Linear Multilinear Algebra, 50 (1) (2002), pp. 1-14
    View in Scopus
    Google Scholar
75. [75]
    Rodriguez J.A.
    Laplacian eigenvalues and partition problems in hypergraphs
    Appl. Math. Lett., 22 (6) (2009), pp. 916-921
    View PDF View article View in Scopus
    Google Scholar
76. [76]
    Rodriguez J.A.
    On the Laplacian spectrum and walk-regular hypergraphs
    Linear Multilinear Algebra, 51 (3) (2003), pp. 285-297
    View in Scopus
    Google Scholar
77. [77]
    J. Körner, Coding of an information source having ambiguous alphabet and the entropy of graphs, in: 6th Prague Conference on Information Theory, 1973, pp. 411–425.
    Google Scholar
78. [78]
    Maletić S., Rajković M.
    Combinatorial Laplacian and entropy of simplicial complexes associated with complex networks
    Eur. Phys. J. Spec. Top., 212 (1) (2012), pp. 77-97
    Crossref
    View in Scopus
    Google Scholar
79. [79]
    Courtney O.T., Bianconi G.
    Weighted growing simplicial complexes
    Phys. Rev. E, 95 (6) (2017), Article 062301
    Crossref
    View in Scopus
    Google Scholar
80. [80]
    Šuvakov M., Andjelković M., Tadić B.
    Hidden geometries in networks arising from cooperative self-assembly
    Sci. Rep., 8 (1) (2018), pp. 1-10
    Google Scholar
81. [81]
    Barrat A., Weigt M.
    On the properties of small-world network models
    Eur. Phys. J. B, 13 (3) (2000), pp. 547-560
    View in Scopus
    Google Scholar
82. [82]
    Newman M.E.
    Scientific collaboration networks. I. Network construction and fundamental results
    Phys. Rev. E, 64 (1) (2001), Article 016131
    Google Scholar
83. [83]
    Bonacich P.
    Factoring and weighting approaches to status scores and clique identification
    J. Math. Sociol., 2 (1) (1972), pp. 113-120
    Crossref
    View in Scopus
    Google Scholar
84. [84]
    Dorogovtsev S.N., Mendes J.F.
    Evolution of networks
    Adv. Phys., 51 (4) (2002), pp. 1079-1187
    View in Scopus
    Google Scholar
85. [85]
    Newman M.E.
    The structure and function of complex networks
    SIAM Rev., 45 (2) (2003), pp. 167-256
    View in Scopus
    Google Scholar
86. [86]
    Estrada E.
    The Structure of Complex Networks: Theory and Applications
    Oxford University Press (2012)
    Google Scholar
87. [87]
    Albert R., Barabási A.-L.
    Statistical mechanics of complex networks
    Rev. Modern Phys., 74 (2002), pp. 47-97
    View in Scopus
    Google Scholar
88. [88]
    Barabási A.-L., Albert R.
    Emergence of scaling in random networks
    Science, 286 (5439) (1999), pp. 509-512
    View in Scopus
    Google Scholar
89. [89]
    Becker N.G., Dietz K.
    The effect of household distribution on transmission and control of highly infectious diseases
    Math. Biosci., 127 (2) (1995), pp. 207-219
    View PDF View article View in Scopus
    Google Scholar
90. [90]
    Ball F., Mollison D., Scalia-Tomba G.
    Epidemics with two levels of mixing
    Ann. Appl. Probab. (1997), pp. 46-89
    View in Scopus
    Google Scholar
91. [91]
    Ball F.
    Stochastic and deterministic models for SIS epidemics among a population partitioned into households
    Math. Biosci., 156 (1–2) (1999), pp. 41-67
    View PDF View article View in Scopus
    Google Scholar
92. [92]
    Ghoshal G., Sander L., Sokolov I.
    SIS epidemics with household structure: The self-consistent field method
    Math. Biosci., 190 (1) (2004), pp. 71-85
    View PDF View article View in Scopus
    Google Scholar
93. [93]
    Ferguson N.M., Cummings D.A., Cauchemez S., Fraser C., Riley S., Meeyai A., Iamsirithaworn S., Burke D.S.
    Strategies for containing an emerging influenza pandemic in Southeast Asia
    Nature, 437 (7056) (2005), pp. 209-214
    Crossref
    View in Scopus
    Google Scholar
94. [94]
    Neal P.
    Stochastic and deterministic analysis of SIS household epidemics
    Adv. Appl. Probab., 38 (4) (2006), pp. 943-968
    View in Scopus
    Google Scholar
95. [95]
    Wu J.T., Riley S., Fraser C., Leung G.M.
    Reducing the impact of the next influenza pandemic using household-based public health interventions
    PLoS Med., 3 (9) (2006), Article e361
    Crossref
    Google Scholar
96. [96]
    Ball F., O’Neill P.D., Pike J.
    Stochastic epidemic models in structured populations featuring dynamic vaccination and isolation
    J. Appl. Probab., 44 (3) (2007), pp. 571-585
    View in Scopus
    Google Scholar
97. [97]
    Dodd P., Ferguson N.
    Approximate disease dynamics in household-structured populations
    J. R. Soc. Interface, 4 (17) (2007), pp. 1103-1106
    Crossref
    View in Scopus
    Google Scholar
98. [98]
    Fraser C.
    Estimating individual and household reproduction numbers in an emerging epidemic
    PLoS One, 2 (8) (2007), Article e758
    Crossref
    View in Scopus
    Google Scholar
99. [99]
    House T., Keeling M.J.
    Deterministic epidemic models with explicit household structure
    Math. Biosci., 213 (1) (2008), pp. 29-39
    View PDF View article View in Scopus
    Google Scholar
100. [100]
     St-Onge G., Thibeault V., Allard A., Dubé L.J., Hébert-Dufresne L.
     Social confinement and mesoscopic localization of epidemics on networks
     Phys. Rev. Lett., 126 (9) (2021), Article 098301
     View in Scopus
     Google Scholar
101. [101]
     St-Onge G., Thibeault V., Allard A., Dubé L.J., Hébert-Dufresne L.
     Master equation analysis of mesoscopic localization in contagion dynamics on higher-order networks
     Phys. Rev. E, 103 (3) (2021), Article 032301
     View in Scopus
     Google Scholar
102. [102]
     Kermack W.O., McKendrick A.G.
     Contributions to the mathematical theory of epidemics – I
     Proc. R. Soc., 115A (1927), pp. 700-721
     Crossref
     View in Scopus
     Google Scholar
103. [103]
     Barrat A., Barthelemy M., Vespignani A.
     Dynamical Processes on Complex Networks
     Cambridge University Press (2008)
     Google Scholar
104. [104]
     Kiss I.Z., Miller J.C., Simon P.L., et al.
     Mathematics of Epidemics on Networks, Vol. 598
     Springer, Cham (2017)
     Google Scholar
105. [105]
     Bodó Á., Katona G.Y., Simon P.L.
     SIS epidemic propagation on hypergraphs
     Bull. Math. Biol., 78 (4) (2016), pp. 713-735
     Crossref
     View in Scopus
     Google Scholar
106. [106]
     Ma T., Guo J.
     Study on information transmission model of enterprise informal organizations based on the hypernetwork
     Chinese J. Phys., 56 (5) (2018), pp. 2424-2438
     View PDF View article View in Scopus
     Google Scholar
107. [107]
     Suo Q., Guo J.-L., Shen A.-Z.
     Information spreading dynamics in hypernetworks
     Physica A, 495 (2018), pp. 475-487
     View PDF View article View in Scopus
     Google Scholar
108. [108]
     Jhun B., Jo M., Kahng B.
     Simplicial SIS model in scale-free uniform hypergraph
     J. Stat. Mech. Theory Exp., 2019 (12) (2019), Article 123207
     Crossref
     View in Scopus
     Google Scholar
109. [109]
     de Arruda G.F., Petri G., Moreno Y.
     Social contagion models on hypergraphs
     Phys. Rev. Res., 2 (2) (2020), Article 023032
     View in Scopus
     Google Scholar
110. [110]
     Wang D., Zhao Y., Luo J., Leng H.
     Simplicial SIRS epidemic models with nonlinear incidence rates
     Chaos, 31 (5) (2021), Article 053112
     View in Scopus
     Google Scholar
111. [111]
     Matamalas J.T., Gómez S., Arenas A.
     Abrupt phase transition of epidemic spreading in simplicial complexes
     Phys. Rev. Res., 2 (1) (2020), Article 012049
     View in Scopus
     Google Scholar
112. [112]
     Burgio G., Arenas A., Gómez S., Matamalas J.T.
     Network clique cover approximation to analyze complex contagions through group interactions
     Commun. Phys., 4 (1) (2021), pp. 1-10
     View in Scopus
     Google Scholar
113. [113]
     Cisneros-Velarde P., Bullo F.
     Multi-group SIS epidemics with simplicial and higher-order interactions
     IEEE Trans. Control Netw. Syst. (2021)
     Google Scholar
114. [114]
     Landry N.W., Restrepo J.G.
     The effect of heterogeneity on hypergraph contagion models
     Chaos, 30 (10) (2020), Article 103117
     View in Scopus
     Google Scholar
115. [115]
     Higham D.J., De Kergorlay H.-L.
     Epidemics on hypergraphs: Spectral thresholds for extinction
     Proc. R. Soc. Lond. Ser. A Math. Phys. Eng. Sci., 477 (2252) (2021), Article 20210232
     View in Scopus
     Google Scholar
116. [116]
     St-Onge G., Sun H., Allard A., Hébert-Dufresne L., Bianconi G.
     Universal nonlinear infection kernel from heterogeneous exposure on higher-order networks
     Phys. Rev. Lett., 127 (15) (2021), Article 158301
     View in Scopus
     Google Scholar
117. [117]
     St-Onge G., Iacopini I., Latora V., Barrat A., Petri G., Allard A., Hébert-Dufresne L.
     Influential groups for seeding and sustaining nonlinear contagion in heterogeneous hypergraphs
     Commun. Phys., 5 (1) (2022), pp. 1-16
     Google Scholar
118. [118]
     Jhun B.
     Effective epidemic containment strategy in hypergraphs
     Phys. Rev. Res., 3 (3) (2021), Article 033282
     View in Scopus
     Google Scholar
119. [119]
     Ghosh D., Frasca M., Rizzo A., Majhi S., Rakshit S., Alfaro-Bittner K., Boccaletti S.
     The synchronized dynamics of time-varying networks
     Phys. Rep., 949 (2022), pp. 1-63
     View PDF View article Crossref
     View in Scopus
     Google Scholar
120. [120]
     Chowdhary S., Kumar A., Cencetti G., Iacopini I., Battiston F.
     Simplicial contagion in temporal higher-order networks
     J. Phys.: Complexity, 2 (3) (2021), Article 035019
     Crossref
     View in Scopus
     Google Scholar
121. [121]
     Wang J.-P., Guo Q., Yang G.-Y., Liu J.-G.
     Improved knowledge diffusion model based on the collaboration hypernetwork
     Physica A, 428 (2015), pp. 250-256
     View PDF View article View in Scopus
     Google Scholar
122. [122]
     Peng G., Wang H., Zhang H., Huang K.
     A hypernetwork-based approach to collaborative retrieval and reasoning of engineering design knowledge
     Adv. Eng. Inform., 42 (2019), Article 100956
     View PDF View article View in Scopus
     Google Scholar
123. [123]
     Barbarossa S., Sardellitti S.
     Topological signal processing over simplicial complexes
     IEEE Trans. Signal Process., 68 (2020), pp. 2992-3007
     Crossref
     View in Scopus
     Google Scholar
124. [124]
     Torres J.J., Bianconi G.
     Simplicial complexes: Higher-order spectral dimension and dynamics
     J. Phys.: Complexity, 1 (1) (2020), Article 015002
     Crossref
     View in Scopus
     Google Scholar
125. [125]
     Bianconi G., Rahmede C.
     Network geometry with flavor: From complexity to quantum geometry
     Phys. Rev. E, 93 (3) (2016), Article 032315
     Crossref
     View in Scopus
     Google Scholar
126. [126]
     Bianconi G., Rahmede C.
     Emergent hyperbolic network geometry
     Sci. Rep., 7 (1) (2017), pp. 1-9
     Google Scholar
127. [127]
     Chan T.-H.H., Louis A., Tang Z.G., Zhang C.
     Spectral properties of hypergraph Laplacian and approximation algorithms
     J. ACM, 65 (3) (2018), pp. 1-48
     Crossref
     Google Scholar
128. [128]
     Chan T.-H.H., Tang Z.G., Wu X., Zhang C.
     Diffusion operator and spectral analysis for directed hypergraph Laplacian
     Theoret. Comput. Sci., 784 (2019), pp. 46-64
     View PDF View article View in Scopus
     Google Scholar
129. [129]
     Chan T.-H.H., Liang Z.
     Generalizing the hypergraph laplacian via a diffusion process with mediators
     Theoret. Comput. Sci., 806 (2020), pp. 416-428
     View PDF View article View in Scopus
     Google Scholar
130. [130]
     Jost J., Mulas R.
     Hypergraph Laplace operators for chemical reaction networks
     Adv. Math., 351 (2019), pp. 870-896
     View PDF View article View in Scopus
     Google Scholar
131. [131]
     Ferraz de Arruda G., Tizzani M., Moreno Y.
     Phase transitions and stability of dynamical processes on hypergraphs
     Commun. Physics, 4 (1) (2021), pp. 1-9
     Crossref
     Google Scholar
132. [132]
     Ugander J., Backstrom L., Marlow C., Kleinberg J.
     Structural diversity in social contagion
     Proc. Natl. Acad. Sci., 109 (16) (2012), pp. 5962-5966
     Crossref
     View in Scopus
     Google Scholar
133. [133]
     SocioPatterns Collaboration
     (2022)
     http://www.sociopatterns.org/
     . (Accessed June 2022)
     Google Scholar
134. [134]
     Matamalas J.T., Arenas A., Gómez S.
     Effective approach to epidemic containment using link equations in complex networks
     Sci. Adv., 4 (12) (2018), p. eaau4212
     View in Scopus
     Google Scholar
135. [135]
     Gómez S., Arenas A., Borge-Holthoefer J., Meloni S., Moreno Y.
     Discrete-time Markov chain approach to contact-based disease spreading in complex networks
     Europhys. Lett., 89 (3) (2010), p. 38009
     Crossref
     View in Scopus
     Google Scholar
136. [136]
     Dorogovtsev S.N., Mendes J.F., Samukhin A.N.
     Size-dependent degree distribution of a scale-free growing network
     Phys. Rev. E, 63 (6) (2001), Article 062101
     Google Scholar
137. [137]
     Cai W., Chen L., Ghanbarnejad F., Grassberger P.
     Avalanche outbreaks emerging in cooperative contagions
     Nat. Phys., 11 (11) (2015), pp. 936-940
     Crossref
     View in Scopus
     Google Scholar
138. [138]
     Poletto C., Meloni S., Van Metre A., Colizza V., Moreno Y., Vespignani A.
     Characterising two-pathogen competition in spatially structured environments
     Sci. Rep., 5 (1) (2015), pp. 1-9
     Google Scholar
139. [139]
     Perra N.
     Non-pharmaceutical interventions during the COVID-19 pandemic: A review
     Phys. Rep., 913 (2021), pp. 1-52
     View PDF View article View in Scopus
     Google Scholar
140. [140]
     Lucas M., Iacopini I., Robiglio T., Barrat A., Petri G.
     Simplicially driven simple contagion
     Phys. Rev. Res., 5 (1) (2023), Article 013201
     View in Scopus
     Google Scholar
141. [141]
     Li W., Xue X., Pan L., Lin T., Wang W.
     Competing spreading dynamics in simplicial complex
     Appl. Math. Comput., 412 (2022), Article 126595
     View PDF View article View in Scopus
     Google Scholar
142. [142]
     Nie Y., Li W., Pan L., Lin T., Wang W.
     Markovian approach to tackle competing pathogens in simplicial complex
     Appl. Math. Comput., 417 (2022), Article 126773
     View PDF View article View in Scopus
     Google Scholar
143. [143]
     Nie Y., Zhong X., Lin T., Wang W.
     Homophily in competing behavior spreading among the heterogeneous population with higher-order interactions
     Appl. Math. Comput., 432 (2022), Article 127380
     View PDF View article View in Scopus
     Google Scholar
144. [144]
     Holme P., Saramäki J.
     Temporal networks
     Phys. Rep., 519 (3) (2012), pp. 97-125
     View PDF View article View in Scopus
     Google Scholar
145. [145]
     Nie Y., Zhong X., Wu T., Liu Y., Lin T., Wang W.
     Effects of network temporality on coevolution spread epidemics in higher-order network
     J. King Saud Univ.-Comput. Inf. Sci. (2022)
     Google Scholar
146. [146]
     Palafox-Castillo G., Berrones-Santos A.
     Stochastic epidemic model on a simplicial complex
     Physica A, 606 (2022), Article 128053
     View PDF View article View in Scopus
     Google Scholar
147. [147]
     Hethcote H.W.
     The mathematics of infectious diseases
     SIAM Rev., 42 (4) (2000), pp. 599-653
     View in Scopus
     Google Scholar
148. [148]
     Jódar L., Villanueva R.J., Arenas A.
     Modeling the spread of seasonal epidemiological diseases: Theory and applications
     Math. Comput. Modelling, 48 (3–4) (2008), pp. 548-557
     View PDF View article View in Scopus
     Google Scholar
149. [149]
     Thornley S., Bullen C., Roberts M.
     Hepatitis B in a high prevalence New Zealand population: A mathematical model applied to infection control policy
     J. Theoret. Biol., 254 (3) (2008), pp. 599-603
     View PDF View article View in Scopus
     Google Scholar
150. [150]
     Pastor-Satorras R., Vespignani A.
     Epidemic spreading in scale-free networks
     Phys. Rev. Lett., 86 (14) (2001), p. 3200
     View in Scopus
     Google Scholar
151. [151]
     Wang D., Zhao Y., Leng H., Small M.
     A social communication model based on simplicial complexes
     Phys. Lett. A, 384 (35) (2020), Article 126895
     View PDF View article View in Scopus
     Google Scholar
152. [152]
     Kozlov D.
     Combinatorial Algebraic Topology, Vol. 21
     Springer Science & Business Media (2008)
     Google Scholar
153. [153]
     Lanchier N., Neufer J.
     Stochastic dynamics on hypergraphs and the spatial majority rule model
     J. Stat. Phys., 151 (1) (2013), pp. 21-45
     Crossref
     View in Scopus
     Google Scholar
154. [154]
     Tsiatas A., Chung F.
     Hypergraph coloring games and voter models
     Internet Math., 10 (1) (2014), p. 1558
     Google Scholar
155. [155]
     Gradowski T., Krawiecki A.
     Majority-vote model on scale-free hypergraphs
     Acta Phys. Polon. A, 127 (3A) (2015)
     Google Scholar
156. [156]
     Horstmeyer L., Kuehn C.
     Adaptive voter model on simplicial complexes
     Phys. Rev. E, 101 (2) (2020), Article 022305
     View in Scopus
     Google Scholar
157. [157]
     Neuhäuser L., Mellor A., Lambiotte R.
     Multibody interactions and nonlinear consensus dynamics on networked systems
     Phys. Rev. E, 101 (3) (2020), Article 032310
     View in Scopus
     Google Scholar
158. [158]
     Alvarez-Rodriguez U., Battiston F., de Arruda G.F., Moreno Y., Perc M., Latora V.
     Evolutionary dynamics of higher-order interactions in social networks
     Nat. Hum. Behav., 5 (5) (2021), pp. 586-595
     Crossref
     View in Scopus
     Google Scholar
159. [159]
     Schlager D., Clauß K., Kuehn C.
     Stability analysis of multiplayer games on adaptive simplicial complexes
     Chaos, 32 (5) (2022), Article 053128
     View in Scopus
     Google Scholar
160. [160]
     Granovetter M.
     Threshold models of collective behavior
     Am. J. Sociol., 83 (6) (1978), pp. 1420-1443
     Google Scholar
161. [161]
     Xie J., Sreenivasan S., Korniss G., Zhang W., Lim C., Szymanski B.K.
     Social consensus through the influence of committed minorities
     Phys. Rev. E, 84 (1) (2011), Article 011130
     Crossref
     View in Scopus
     Google Scholar
162. [162]
     Castellano C., Fortunato S., Loreto V.
     Statistical physics of social dynamics
     Rev. Modern Phys., 81 (2) (2009), p. 591
     View in Scopus
     Google Scholar
163. [163]
     Vespignani A.
     Modelling dynamical processes in complex socio-technical systems
     Nat. Phys., 8 (1) (2012), pp. 32-39
     Crossref
     View in Scopus
     Google Scholar
164. [164]
     Dall’Asta L., Baronchelli A., Barrat A., Loreto V.
     Nonequilibrium dynamics of language games on complex networks
     Phys. Rev. E, 74 (3) (2006), Article 036105
     Crossref
     Google Scholar
165. [165]
     Iacopini I., Petri G., Baronchelli A., Barrat A.
     Group interactions modulate critical mass dynamics in social convention
     Commun. Phys., 5 (1) (2022), pp. 1-10
     Crossref
     Google Scholar
166. [166]
     Pikovsky A., Rosenblum M., Kurths J.
     Synchronization: A Universal Concept in Nonlinear Sciences, Vol. 12
     Cambridge University Press (2003)
     Google Scholar
167. [167]
     Boccaletti S., Pisarchik A.N., Del Genio C.I., Amann A.
     Synchronization: From Coupled Systems to Complex Networks
     Cambridge University Press (2018)
     Google Scholar
168. [168]
     Barahona M., Pecora L.M.
     Synchronization in small-world systems
     Phys. Rev. Lett., 89 (5) (2002), Article 054101
     Google Scholar
169. [169]
     Chavez M., Hwang D.-U., Amann A., Hentschel H.G.E., Boccaletti S.
     Synchronization is enhanced in weighted complex networks
     Phys. Rev. Lett., 94 (2005), Article 218701
     Crossref
     View in Scopus
     Google Scholar
170. [170]
     del Genio C.I., Gómez-Gardeñes J., Bonamassa I., Boccaletti S.
     Synchronization in networks with multiple interaction layers
     Sci. Adv., 2 (11) (2016)
     Google Scholar
171. [171]
     Gutiérrez R., Amann A., Assenza S., Gómez-Gardenes J., Latora V., Boccaletti S.
     Emerging meso-and macroscales from synchronization of adaptive networks
     Phys. Rev. Lett., 107 (23) (2011), Article 234103
     Crossref
     View in Scopus
     Google Scholar
172. [172]
     Avalos-Gaytán V., Almendral J.A., Leyva I., Battiston F., Nicosia V., Latora V., Boccaletti S.
     Emergent explosive synchronization in adaptive complex networks
     Phys. Rev. E, 97 (2018), Article 042301
     View in Scopus
     Google Scholar
173. [173]
     Gambuzza L.V., Cardillo A., Fiasconaro A., Fortuna L., Gómez-Gardenes J., Frasca M.
     Analysis of remote synchronization in complex networks
     Chaos, 23 (4) (2013), Article 043103
     View in Scopus
     Google Scholar
174. [174]
     Nicosia V., Valencia M., Chavez M., Díaz-Guilera A., Latora V.
     Remote synchronization reveals network symmetries and functional modules
     Phys. Rev. Lett., 110 (2013), Article 174102
     Crossref
     View in Scopus
     Google Scholar
175. [175]
     Pecora L.M., Sorrentino F., Hagerstrom A.M., Murphy T.E., Roy R.
     Cluster synchronization and isolated desynchronization in complex networks with symmetries
     Nature Commun., 5 (1) (2014), pp. 1-8
     Crossref
     Google Scholar
176. [176]
     Gambuzza L.V., Frasca M., Latora V.
     Distributed control of synchronization of a group of network nodes
     IEEE Trans. Automat. Control, 64 (1) (2018), pp. 365-372
     Google Scholar
177. [177]
     Abrams D.M., Strogatz S.H.
     Chimera states for coupled oscillators
     Phys. Rev. Lett., 93 (17) (2004), Article 174102
     Google Scholar
178. [178]
     Panaggio M.J., Abrams D.M.
     Chimera states: Coexistence of coherence and incoherence in networks of coupled oscillators
     Nonlinearity, 28 (3) (2015), p. R67
     View in Scopus
     Google Scholar
179. [179]
     Bi H., Hu X., Boccaletti S., Wang X., Zou Y., Liu Z., Guan S.
     Coexistence of quantized, time dependent, clusters in globally coupled oscillators
     Phys. Rev. Lett., 117 (20) (2016), Article 204101
     Crossref
     View in Scopus
     Google Scholar
180. [180]
     Xu C., Boccaletti S., Guan S., Zheng Z.
     Origin of Bellerophon states in globally coupled phase oscillators
     Phys. Rev. E, 98 (5) (2018), Article 050202
     View in Scopus
     Google Scholar
181. [181]
     Boccaletti S., Almendral J., Guan S., Leyva I., Liu Z., Sendiña-Nadal I., Wang Z., Zou Y.
     Explosive transitions in complex networks’ structure and dynamics: Percolation and synchronization
     Phys. Rep., 660 (2016), pp. 1-94
     View PDF View article View in Scopus
     Google Scholar
182. [182]
     Belykh V.N., Belykh I.V., Hasler M.
     Connection graph stability method for synchronized coupled chaotic systems
     Physica D, 195 (1–2) (2004), pp. 159-187
     View PDF View article View in Scopus
     Google Scholar
183. [183]
     Li Z., Chen G.
     Global synchronization and asymptotic stability of complex dynamical networks
     IEEE Trans. Circuits Syst. II, 53 (1) (2006), pp. 28-33
     View in Scopus
     Google Scholar
184. [184]
     Pecora L.M., Carroll T.L.
     Master stability functions for synchronized coupled systems
     Phys. Rev. Lett., 80 (10) (1998), p. 2109
     View in Scopus
     Google Scholar
185. [185]
     Mulas R., Kuehn C., Jost J.
     Coupled dynamics on hypergraphs: Master stability of steady states and synchronization
     Phys. Rev. E, 101 (6) (2020), Article 062313
     View in Scopus
     Google Scholar
186. [186]
     Krawiecki A.
     Chaotic synchronization on complex hypergraphs
     Chaos Solitons Fractals, 65 (2014), pp. 44-50
     View PDF View article View in Scopus
     Google Scholar
187. [187]
     Carletti T., Fanelli D., Nicoletti S.
     Dynamical systems on hypergraphs
     J. Phys.: Complexity, 1 (3) (2020), Article 035006
     Crossref
     View in Scopus
     Google Scholar
188. [188]
     Sun J., Bollt E.M., Nishikawa T.
     Master stability functions for coupled nearly identical dynamical systems
     Europhys. Lett., 85 (6) (2009), p. 60011
     Crossref
     View in Scopus
     Google Scholar
189. [189]
     Stilwell D.J., Bollt E.M., Roberson D.G.
     Sufficient conditions for fast switching synchronization in time-varying network topologies
     SIAM J. Appl. Dyn. Syst., 5 (1) (2006), pp. 140-156
     Crossref
     View in Scopus
     Google Scholar
190. [190]
     Frasca M., Buscarino A., Rizzo A., Fortuna L., Boccaletti S.
     Synchronization of moving chaotic agents
     Phys. Rev. Lett., 100 (4) (2008), Article 044102
     Crossref
     View in Scopus
     Google Scholar
191. [191]
     Zhou J., Zou Y., Guan S., Liu Z., Boccaletti S.
     Synchronization in slowly switching networks of coupled oscillators
     Sci. Rep., 6 (1) (2016), pp. 1-8
     Google Scholar
192. [192]
     Della Rossa F., DeLellis P.
     Stochastic master stability function for noisy complex networks
     Phys. Rev. E, 101 (5) (2020), Article 052211
     View in Scopus
     Google Scholar
193. [193]
     Rössler O.E.
     An equation for continuous chaos
     Phys. Lett. A, 57 (5) (1976), pp. 397-398
     View PDF View article View in Scopus
     Google Scholar
194. [194]
     Strogatz S.H.
     Nonlinear Dynamics and Chaos with Student Solutions Manual: With Applications to Physics, Biology, Chemistry, and Engineering
     CRC Press (2018)
     Google Scholar
195. [195]
     Acebrón J.A., Bonilla L.L., Vicente C.J.P., Ritort F., Spigler R.
     The Kuramoto model: A simple paradigm for synchronization phenomena
     Rev. Modern Phys., 77 (1) (2005), p. 137
     View in Scopus
     Google Scholar
196. [196]
     Rodrigues F.A., Peron T.K.D., Ji P., Kurths J.
     The Kuramoto model in complex networks
     Phys. Rep., 610 (2016), pp. 1-98
     View PDF View article View in Scopus
     Google Scholar
197. [197]
     Kuramoto Y.
     Chemical turbulence
     Chemical Oscillations, Waves, and Turbulence, Springer (1984), pp. 111-140
     Crossref
     Google Scholar
198. [198]
     Strogatz S.H.
     From Kuramoto to Crawford: Exploring the onset of synchronization in populations of coupled oscillators
     Physica D, 143 (1–4) (2000), pp. 1-20
     View PDF View article View in Scopus
     Google Scholar
199. [199]
     Kumar A., Jalan S.
     Explosive synchronization in interlayer phase-shifted Kuramoto oscillators on multiplex networks
     Chaos: Interdiscipl. J. Nonlinear Sci. (Fast Track), 3 (2021), Article 041103
     View in Scopus
     Google Scholar
200. [200]
     Kachhvah A.D., Jalan S.
     Explosive synchronization and chimera in interpinned multilayer networks
     Phys. Rev. E, 104 (2021)
     Google Scholar
201. [201]
     Tanaka T., Aoyagi T.
     Multistable attractors in a network of phase oscillators with three-body interactions
     Phys. Rev. Lett., 106 (22) (2011), Article 224101
     Crossref
     View in Scopus
     Google Scholar
202. [202]
     Skardal P.S., Arenas A.
     Higher order interactions in complex networks of phase oscillators promote abrupt synchronization switching
     Commun. Phys., 3 (2020), p. 218
     View in Scopus
     Google Scholar
203. [203]
     Millán A.P., Torres J.J., Bianconi G.
     Explosive higher-order Kuramoto dynamics on simplicial complexes
     Phys. Rev. Lett., 124 (218301) (2020)
     Google Scholar
204. [204]
     Lucas M., Cencetti G., Battiston F.
     Multiorder Laplacian for synchronization in higher-order networks
     Phys. Rev. Res., 2 (3) (2020), Article 033410
     View in Scopus
     Google Scholar
205. [205]
     Skardal P.S., Arenas A.
     Memory selection and information switching in oscillator networks with higher-order interactions
     J. Phys.: Complexity, 2 (1) (2020), Article 015003
     Google Scholar
206. [206]
     Kundu S., Ghosh D.
     Higher-order interactions promote chimera states
     Phys. Rev. E, 105 (L042202) (2022)
     Google Scholar
207. [207]
     Kovalenko K., Dai X., Alfaro-Bittner K., Raigorodskii A.M., Perc M., Boccaletti S.
     Contrarians synchronize beyond the limit of pairwise interactions
     Phys. Rev. Lett., 127 (2021), Article 258301
     View in Scopus
     Google Scholar
208. [208]
     Ott E., Antonsen T.M.
     Low dimensional behavior of large systems of globally coupled oscillators
     Chaos, 18 (3) (2008), Article 037113
     View in Scopus
     Google Scholar
209. [209]
     Bick C., Böhle T., Kuehn C.
     Multi-population phase oscillator networks with higher-order interactions
     NoDEA Nonlinear Differential Equations Appl., 29 (64) (2022)
     Google Scholar
210. [210]
     Jalan S., Suman A.
     Multiple first-order transitions in simplicial complexes on multilayer systems
     Phys. Rev. E, 106 (4) (2022), Article 044304
     View in Scopus
     Google Scholar
211. [211]
     Philip Seliger S.C.Y., Tsimring L.S.
     Plasticity and learning in a network of coupled phase oscillators
     Phys. Rev. E, 65 (041906) (2002)
     Google Scholar
212. [212]
     Zhou C., Kurths J.
     Dynamical weights and enhanced synchronization in adaptive complex networks
     Phys. Rev. Lett., 96 (16) (2006), Article 164102
     Crossref
     View in Scopus
     Google Scholar
213. [213]
     De Lellis P., di Bernardo M., Garofalo F.
     Novel decentralized adaptive strategies for the synchronization of complex networks
     Automatica, 45 (5) (2009), pp. 1312-1318
     View PDF View article View in Scopus
     Google Scholar
214. [214]
     Šiljak D.
     Dynamic graphs
     Nonlinear Anal. Hybrid Syst., 2 (2) (2008), pp. 544-567
     View PDF View article View in Scopus
     Google Scholar
215. [215]
     De Lellis P., Garofalo F., Porfiri M., et al.
     Evolution of complex networks via edge snapping
     IEEE Trans. Circuits Syst. I. Regul. Pap., 57 (8) (2010), pp. 2132-2143
     View in Scopus
     Google Scholar
216. [216]
     Zhang X., Boccaletti S., Guan S., Liu Z.
     Explosive synchronization in adaptive and multilayer networks
     Phys. Rev. Lett., 114 (038701) (2015)
     Google Scholar
217. [217]
     Hebb D.O.
     The Organization of Behavior: A Neuropsychological Theory
     Psychology Press (2005)
     Google Scholar
218. [218]
     H. Markram M.F., Sakmann B.
     Regulation of synaptic efficacy by coincidence of postsynaptic APs and EPSPs
     Science, 275 (213) (1997)
     Google Scholar
219. [219]
     Berner R., Sawicki J., Schöll E.
     Birth and stabilization of phase clusters by multiplexing of adaptive networks
     Phys. Rev. Lett., 124 (8) (2020), Article 088301
     View in Scopus
     Google Scholar
220. [220]
     Pitsik E., Makarov V., Kirsanov D., Frolov N., Goremyko M., Li X., Wang Z., Hramov A., Boccaletti S.
     Inter-layer competition in adaptive multiplex network
     New J. Phys., 20 (7) (2018), Article 075004
     Crossref
     View in Scopus
     Google Scholar
221. [221]
     Kachhvah A.D., Dai X., Boccaletti S., Jalan S.
     Interlayer hebbian plasticity induces first-order transition in multiplex networks
     New J. Phys., 22 (12) (2020), Article 122001
     Crossref
     View in Scopus
     Google Scholar
222. [222]
     Kachhvah A.D., Jalan S.
     Hebbian plasticity rules abrupt desynchronization in pure simplicial complexes
     New J. Phys. (Fast Track), 24 (2022), Article 052002
     Crossref
     View in Scopus
     Google Scholar
223. [223]
     Kachhvah A.D., Jalan S.
     First-order route to antiphase clustering in adaptive simplicial complexes
     Phys. Rev. E, 105 (6) (2022), p. L062203
     View in Scopus
     Google Scholar
224. [224]
     Malayaja Chutani B.T., Gupte N.
     Hysteresis and synchronization processes of Kuramoto oscillators on high-dimensional simplicial complexes with competing simplex-encoded couplings
     Phys. Rev. E, 104 (2021), Article 034206
     Google Scholar
225. [225]
     Motter A.E., Myers S.A., Anghel M., Nishikawa T.
     Spontaneous synchrony in power-grid networks
     Nat. Phys., 9 (3) (2013), pp. 191-197
     Crossref
     View in Scopus
     Google Scholar
226. [226]
     Schnitzler A., Gross J.
     Normal and pathological oscillatory communication in the brain
     Nat. Rev. Neurosci., 6 (4) (2005), pp. 285-296
     Crossref
     View in Scopus
     Google Scholar
227. [227]
     DeVille L.
     Consensus on simplicial complexes: Results on stability and synchronization
     Chaos, 31 (2) (2021), Article 023137
     View in Scopus
     Google Scholar
228. [228]
     Araque A., Carmignoto G., Haydon P.G., Oliet S.H., Robitaille R., Volterra A.
     Gliotransmitters travel in time and space
     Neuron, 81 (4) (2014), pp. 728-739
     View PDF View article View in Scopus
     Google Scholar
229. [229]
     Ghorbanchian R., Restrepo J.G., Torres P.G., Bianconi G.
     Higher-order simplicial synchronization of coupled topological signals
     Commun. Phys., 4 (1) (2021), pp. 1-13
     View in Scopus
     Google Scholar
230. [230]
     Hagmann P., Cammoun L., Gigandet X., Meuli R., Honey C.J., Wedeen V.J., Sporns O.
     Mapping the structural core of human cerebral cortex
     PLoS Biol., 6 (7) (2008), Article e159
     Crossref
     Google Scholar
231. [231]
     Varshney L.R., Chen B.L., Paniagua E., Hall D.H., Chklovskii D.B.
     Structural properties of the Caenorhabditis elegans neuronal network
     PLoS Comput. Biol., 7 (2) (2011), Article e1001066
     Crossref
     View in Scopus
     Google Scholar
232. [232]
     Arnaudon A., Peach R.L., Petri G., Expert P.
     Connecting hodge and Sakaguchi-Kuramoto through a mathematical framework for coupled oscillators on simplicial complexes
     Commun. Phys., 5 (1) (2022), pp. 1-12
     Google Scholar
233. [233]
     Sakaguchi H., Kuramoto Y.
     A soluble active rotater model showing phase transitions via mutual entertainment
     Progr. Theoret. Phys., 76 (3) (1986), pp. 576-581
     Google Scholar
234. [234]
     Anderson B., Ye M.
     Recent advances in the modelling and analysis of opinion dynamics on influence networks
     Int. J. Autom. Comput., 16 (2) (2019), pp. 129-149
     Crossref
     View in Scopus
     Google Scholar
235. [235]
     Lewis F.L.
     Wireless sensor networks
     Smart Environ.: Technol., Protoc., Appl. (2004), pp. 11-46
     Crossref
     Google Scholar
236. [236]
     Lin J., Morse A.S., Anderson B.D.
     The multi-agent rendezvous problem
     42nd IEEE International Conference on Decision and Control, Vol. 2, IEEE (2003), pp. 1508-1513
     View in Scopus
     Google Scholar
237. [237]
     Bertsekas D., Tsitsiklis J.
     Parallel and Distributed Computation: Numerical Methods
     Prentice-Hall, Englewood Cliff, NJ (1989)
     Google Scholar
238. [238]
     Chatterjee S., Seneta E.
     Towards consensus: Some convergence theorems on repeated averaging
     J. Appl. Probab., 14 (1) (1977), pp. 89-97
     Google Scholar
239. [239]
     DeGroot M.H.
     Reaching a consensus
     J. Amer. Statist. Assoc., 69 (345) (1974), pp. 118-121
     View in Scopus
     Google Scholar
240. [240]
     Siljak D.D.
     Large-Scale Dynamic Systems: Stability and Structure
     North-Holland New York (1978)
     Google Scholar
241. [241]
     Olfati-Saber R., Murray R.M.
     Consensus problems in networks of agents with switching topology and time-delays
     IEEE Trans. Automat. Control, 49 (9) (2004), pp. 1520-1533
     View in Scopus
     Google Scholar
242. [242]
     Blondel V.D., Hendrickx J.M., Tsitsiklis J.N.
     On Krause’s multi-agent consensus model with state-dependent connectivity
     IEEE Trans. Automat. Control, 54 (11) (2009), pp. 2586-2597
     View in Scopus
     Google Scholar
243. [243]
     Neuhäuser L., Lambiotte R., Schaub M.T.
     Consensus dynamics on temporal hypergraphs
     Phys. Rev. E, 104 (6) (2021), Article 064305
     View in Scopus
     Google Scholar
244. [244]
     Neuhäuser L., Lambiotte R., Schaub M.T.
     Consensus dynamics and opinion formation on hypergraphs
     Higher-Order Systems, Springer (2022), pp. 347-376
     Crossref
     View in Scopus
     Google Scholar
245. [245]
     Nijholt E., DeVille L.
     Dynamical systems defined on simplicial complexes: Symmetries, conjugacies, and invariant subspaces
     Chaos, 32 (9) (2022), Article 093131
     View in Scopus
     Google Scholar
246. [246]
     Sahasrabuddhe R., Neuhäuser L., Lambiotte R.
     Modelling non-linear consensus dynamics on hypergraphs
     J. Phys.: Complexity, 2 (2) (2021), Article 025006
     Crossref
     View in Scopus
     Google Scholar
247. [247]
     Schawe H., Hernández L.
     Higher order interactions destroy phase transitions in Deffuant opinion dynamics model
     Commun. Phys., 5 (1) (2022), pp. 1-9
     Google Scholar
248. [248]
     Asch S.E.
     Effects of group pressure upon the modification and distortion of judgments
     Guetzkow H. (Ed.), Groups, Leadership and Men, Carnegie Press (1951), pp. 177-190
     Google Scholar
249. [249]
     McPherson M., Smith-Lovin L., Cook J.M.
     Birds of a feather: Homophily in social networks
     Annu. Rev. Sociol. (2001), pp. 415-444
     View in Scopus
     Google Scholar
250. [250]
     Deffuant G., Neau D., Amblard F., Weisbuch G.
     Mixing beliefs among interacting agents
     Adv. Complex Syst., 3 (01n04) (2000), pp. 87-98
     Google Scholar
251. [251]
     Olfati-Saber R., Fax J.A., Murray R.M.
     Consensus and cooperation in networked multi-agent systems
     Proc. IEEE, 95 (1) (2007), pp. 215-233
     View in Scopus
     Google Scholar
252. [252]
     Clifford P., Sudbury A.
     A model for spatial conflict
     Biometrika, 60 (3) (1973), pp. 581-588
     Google Scholar
253. [253]
     Holley R.A., Liggett T.M.
     Ergodic theorems for weakly interacting infinite systems and the voter model
     Ann. Probab. (1975), pp. 643-663
     Google Scholar
254. [254]
     Redner S.
     A Guide to First-Passage Processes
     Cambridge University Press (2001)
     Google Scholar
255. [255]
     Galam S.
     Minority opinion spreading in random geometry
     Eur. Phys. J. B, 25 (4) (2002), pp. 403-406
     View in Scopus
     Google Scholar
256. [256]
     Noonan J., Lambiotte R.
     Dynamics of majority rule on hypergraphs
     Phys. Rev. E, 104 (2) (2021), Article 024316
     View in Scopus
     Google Scholar
257. [257]
     Guckenheimer J., Holmes P.
     Structurally stable heteroclinic cycles
     Mathematical Proceedings of the Cambridge Philosophical Society, vol. 103, Cambridge University Press (1988), pp. 189-192
     View in Scopus
     Google Scholar
258. [258]
     Lorenz E.N.
     Deterministic nonperiodic flow
     J. Atmos. Sci., 20 (2) (1963), pp. 130-141
     View in Scopus
     Google Scholar
259. [259]
     Sel’Kov E.
     Self-oscillations in glycolysis 1. A simple kinetic model
     Eur. J. Biochem., 4 (1) (1968), pp. 79-86
     Crossref
     View in Scopus
     Google Scholar
260. [260]
     Liu Y.-Y., Barabási A.-L.
     Control principles of complex systems
     Rev. Modern Phys., 88 (3) (2016), Article 035006
     Crossref
     View in Scopus
     Google Scholar
261. [261]
     Chen T., Liu X., Lu W.
     Pinning complex networks by a single controller
     IEEE Trans. Circuits Syst. I. Regul. Pap., 54 (6) (2007), pp. 1317-1326
     View in Scopus
     Google Scholar
262. [262]
     Della Rossa F., De Lellis P.
     Synchronization and pinning control of stochastic coevolving networks
     Annu. Rev. Control, 53 (2022), pp. 147-160
     View PDF View article View in Scopus
     Google Scholar
263. [263]
     Grigoriev R., Cross M., Schuster H.
     Pinning control of spatiotemporal chaos
     Phys. Rev. Lett., 79 (15) (1997), p. 2795
     View in Scopus
     Google Scholar
264. [264]
     Lu W., Li X., Rong Z.
     Global stabilization of complex networks with digraph topologies via a local pinning algorithm
     Automatica, 46 (1) (2010), pp. 116-121
     View PDF View article View in Scopus
     Google Scholar
265. [265]
     Porfiri M., di Bernardo M.
     Criteria for global pinning-controllability of complex networks
     Automatica, 44 (12) (2008), pp. 3100-3106
     View PDF View article View in Scopus
     Google Scholar
266. [266]
     Sorrentino F., di Bernardo M., Garofalo F., Chen G.
     Controllability of complex networks via pinning
     Phys. Rev. E, 75 (4) (2007), Article 046103
     Crossref
     View in Scopus
     Google Scholar
267. [267]
     Wang X.F., Chen G.
     Pinning control of scale-free dynamical networks
     Physica A, 310 (3–4) (2002), pp. 521-531
     View PDF View article View in Scopus
     Google Scholar
268. [268]
     Gao J., Liu Y.-Y., D’souza R.M., Barabási A.-L.
     Target control of complex networks
     Nature Commun., 5 (1) (2014), pp. 1-8
     Google Scholar
269. [269]
     Liu Y.-Y., Slotine J.-J., Barabási A.-L.
     Controllability of complex networks
     Nature, 473 (7346) (2011), pp. 167-173
     Crossref
     View in Scopus
     Google Scholar
270. [270]
     Lo Iudice F., Garofalo F., Sorrentino F.
     Structural permeability of complex networks to control signals
     Nature Commun., 6 (1) (2015), pp. 1-6
     Google Scholar
271. [271]
     Yuan Z., Zhao C., Di Z., Wang W.-X., Lai Y.-C.
     Exact controllability of complex networks
     Nature Commun., 4 (1) (2013), pp. 1-9
     View in Scopus
     Google Scholar
272. [272]
     Yan G., Ren J., Lai Y.-C., Lai C.-H., Li B.
     Controlling complex networks: How much energy is needed?
     Phys. Rev. Lett., 108 (21) (2012), Article 218703
     Crossref
     View in Scopus
     Google Scholar
273. [273]
     Sun J., Motter A.E.
     Controllability transition and nonlocality in network control
     Phys. Rev. Lett., 110 (20) (2013), Article 208701
     Crossref
     View in Scopus
     Google Scholar
274. [274]
     Pasqualetti F., Zampieri S., Bullo F.
     Controllability metrics, limitations and algorithms for complex networks
     IEEE Trans. Control Netw. Syst., 1 (1) (2014), pp. 40-52
     View in Scopus
     Google Scholar
275. [275]
     Li A., Cornelius S.P., Liu Y.-Y., Wang L., Barabási A.-L.
     The fundamental advantages of temporal networks
     Science, 358 (6366) (2017), pp. 1042-1046
     Crossref
     View in Scopus
     Google Scholar
276. [276]
     De Lellis P., Di Meglio A., Garofalo F., Lo Iudice F.
     The inherent uncertainty of temporal networks is a true challenge for control
     Sci. Rep., 11 (1) (2021), pp. 1-7
     Crossref
     View in Scopus
     Google Scholar
277. [277]
     Chen C., Surana A., Bloch A.M., Rajapakse I.
     Controllability of hypergraphs
     IEEE Trans. Netw. Sci. Eng., 8 (2) (2021), pp. 1646-1657
     Crossref
     View in Scopus
     Google Scholar
278. [278]
     Salzano D., Fiore D., di Bernardo M.
     Ratiometric control of cell phenotypes in monostrain microbial consortia
     J. R. Soc. Interface, 19 (192) (2022), Article 20220335
     View in Scopus
     Google Scholar
279. [279]
     De Lellis P., Della Rossa F., Iudice F.L., Liuzza D.
     Pinning control of hypergraphs
     IEEE Control Syst. Lett., 7 (2022), pp. 691-696
     Google Scholar
280. [280]
     Gallo G., Longo G., Pallottino S., Nguyen S.
     Directed hypergraphs and applications
     Discrete Appl. Math., 42 (2–3) (1993), pp. 177-201
     View PDF View article View in Scopus
     Google Scholar
281. [281]
     Kalman R.E.
     Mathematical description of linear dynamical systems
     J. Soc. Ind. Appl. Math., Ser. A: Control, 1 (2) (1963), pp. 152-192
     Crossref
     Google Scholar
282. [282]
     Fujarewicz K., Kimmel M., Swierniak A.
     On fitting of mathematical models of cell signaling pathways using adjoint systems
     Math. Biosci. Eng., 2 (3) (2005), p. 527
     Crossref
     Google Scholar
283. [283]
     Jurdjevic V., Kupka I.
     Polynomial control systems
     Math. Ann., 272 (3) (1985), pp. 361-368
     View in Scopus
     Google Scholar
284. [284]
     Melody J., Basar T., Bullo F.
     On nonlinear controllability of homogeneous systems linear in control
     IEEE Trans. Automat. Control, 48 (1) (2003), pp. 139-143
     View in Scopus
     Google Scholar
285. [285]
     Aeyels D.
     Local and global controllability for nonlinear systems
     Systems Control Lett., 5 (1) (1984), pp. 19-26
     View PDF View article View in Scopus
     Google Scholar
286. [286]
     Zeng A., Shen Z., Zhou J., Wu J., Fan Y., Wang Y., Stanley H.E.
     The science of science: From the perspective of complex systems
     Phys. Rep., 714 (2017), pp. 1-73
     View PDF View article View in Scopus
     Google Scholar
287. [287]
     Tanaka H.-A., Lichtenberg A.J., Oishi S.
     First order phase transition resulting from finite inertia in coupled oscillator systems
     Phys. Rev. Lett., 78 (11) (1997), p. 2104
     View in Scopus
     Google Scholar
288. [288]
     Rohden M., Sorge A., Timme M., Witthaut D.
     Self-organized synchronization in decentralized power grids
     Phys. Rev. Lett., 109 (6) (2012), Article 064101
     Crossref
     View in Scopus
     Google Scholar
289. [289]
     Dörfler F., Chertkov M., Bullo F.
     Synchronization in complex oscillator networks and smart grids
     Proc. Natl. Acad. Sci., 110 (6) (2013), pp. 2005-2010
     Crossref
     View in Scopus
     Google Scholar
290. [290]
     Filatrella G., Nielsen A.H., Pedersen N.F.
     Analysis of a power grid using a Kuramoto-like model
     Eur. Phys. J. B, 61 (2008), pp. 485-491
     Crossref
     View in Scopus
     Google Scholar
291. [291]
     Sabhahit N.G., Khurd A.S., Jalan S.
     Self-consistent method for Kuramoto oscillators with inertia having higher-order interactions
     (2023)
     arXiv preprint arXiv:2303.08363
     Google Scholar
292. [292]
     Guo T., Zhang Y., Xue Y., Qiao L., Shen D.
     Brain function network: Higher order vs. more discrimination
     Front. Neurosci., 15 (2021), Article 696639
     View in Scopus
     Google Scholar

## Cited by (217)

- ### Patterns of neuronal synchrony in higher-order networks

  2025, Physics of Life Reviews

- ### Reinforcement learning and collective cooperation on higher-order networks

  2024, Knowledge Based Systems

- ### Synchronization in multiplex networks

  2024, Physics Reports

- ### Epidemic spreading on higher-order networks

  2024, Physics Reports

- ### Contagion dynamics on higher-order networks

  2024, Nature Reviews Physics

- ### Multi-chimera states in a higher order network of FitzHugh–Nagumo oscillators

  2024, European Physical Journal Special Topics

View all citing articles on Scopus

1
  According to [112], given a simplicial complex
  K, a clique of
  n+1 nodes is a
  (1,n+1)-clique if its nodes form a
  n-simplex in
  K, whereas it is a
  (0,n+1)-clique. otherwise. Since a 1-simplex is equivalent to a 2-clique, the Authors choose to assign
  g=0 to any 2-clique.

2
  In standard
  k-uniform hyperchains, hyperrings and hyperstars,
  r=k−1.

View Abstract
© 2023 Elsevier B.V. All rights reserved.
The structure and dynamics of networks with higher order interactions

### Citations

- Citation Indexes 214

### Captures

- Mendeley Readers 74

### Social Media

- Shares, Likes & Comments 121

View details