Network Dynamics with Higher-Order Interactions: Coupled Cell Hypernetworks for Identical Cells and Synchrony

Network interactions that are nonlinear in the state of more than two nodes - also known as higher-order interactions - can have a profound impact on the collective network dynamics. Here we develop a coupled cell hypernetwork formalism to elucidate the existence and stability of (cluster) synchronization patterns in network dynamical systems with higher-order interactions. More specifically, we define robust synchrony subspace for coupled cell hypernetworks whose coupling structure is determined by an underlying hypergraph and describe those spaces for general such hypernetworks. Since a hypergraph can be equivalently represented as a bipartite graph between its nodes and hyperedges, we relate the synchrony subspaces of a hypernetwork to balanced colorings of the corresponding incidence digraph.


Introduction
Coupled dynamical processes are ubiquitous in the world and can often be modeled by systems of ordinary differential equations (ODEs). The coupled cell network formalism developed by Golubitsky, Stewart and collaborators [1,2] and Field [3] captures the network interactions by a directed graph G to elucidate how the network structure shapes the collective dynamics. More precisely, let V = R d for some d ∈ N denote the state space of each cell i ∈ {1, . . . , n}. In a classical coupled cell system, the evolution state x i of cell i is determined by an interaction function f : V q → V . If, for example, then (j, i), (k, i), (l, i) are the edges with head i of G since, for any f , the evolution of cell i depends on the cells j, k, l. The main questions regarding coupled cell networks relate to how the network structure influences the dynamics and bifurcations of the coupled cell system without making specific assumptions on f . By contrast, in many applications the links in the networks have associated numerical values 1 called weights to represent, for example, the strength or the signal of the connection between the nodes associated with the edges. These can be realized as coupled cell networks with additive input structure; cf. [4,5,6,7]. Consider the graph G associated with (1.1) and let (w ij ) ∈ R n×n be a weight matrix. For h : V → V and g : V × V → V , cell i of the corresponding coupled cell network with additive coupling structure evolves according to where g determines the pairwise interactions between cells. In this restricted framework, adding and removing edges is natural by adjusting the corresponding weights. Networks of Kuramoto phase oscillators and pulse coupled systems are examples of coupled cell systems with additive input structure. Note that the complexity of the interactions differ in traditional coupled cell networks (1.1) and those with additive coupling structure (1.2). While the former allows for generic, nonlinear interactions between all the input nodes through f , additive coupling structure only allows for interactions between pairs of nodes. Recent research has highlighted the dynamical importance of nonpairwise interactions between nodes; cf. [8,9] for reviews. For example, in networks that describe the competitive interactions between species, one has to take into account how the interaction between two species is modulated by a third species (a triplet interaction) to explain the competition dynamics. Similarly, incorporating nonpairwise interactions in phase oscillator networks exhibits dynamics that would not arise in standard Kuramoto-type equations with pairwise interactions [10,11].
In this work, we introduce a new class of coupled cell networkscoupled cell hypernetworks-whose structure is determined by a (directed) hypergraph. A hypergraph is a generalization of a graph in which a hyperedge can join any number of nodes, that is, the directed hyperedges are from a set of k nodes (cells) to a set of l nodes (cells). This coupling structure captures that the evolution of each of the l cells depends (typically nonlinearly) on an interaction involving a set of k cells. Directed hypergraphs are used to model problems arising in, for example, operations research, computer science and discrete mathematics, to describe relationships between two sets of objects. See for example Ausiello and Laura [12] and references therein. See, also, Johnson et al. [13], Kim et al. [14] and Johnson [15]. We shall remark that in some literature, as for example in Sorrentino [16], the terminology of hypernetwork is used, not to denote a hypergraph, as in our case here, but to denote a graph that has more than one edge type, that is, with more than one adjacency matrix. We illustrate our setup in an example. Example 1.1. Consider the following system of ODEs on n = 6 state variables x i , i ∈ {1, . . . , n}: Assume that Q 2 is symmetric under permutation of the last two coordinates, that is, Q 2 (y; z, w) = Q 2 (y; w, z) for all y, z, w ∈ V . We might interpret this system as a coupled cell system with form consistent with a hypergraph H shown on the left of Figure 1: Each node of the hypergraph represents a cell, and each hyperedge represents an interaction from a cell-or a group of cells-to a cell or a group of cells. The state of cell i is determined by x i ∈ V and its evolution by the corresponding differential equation; in the following we write x = (x 1 , . . . , x n ) for the state vector. The coupling functions Q 1 and Q 2 determine the influence of one or two cells, respectively, onto another cell along the directed hyperedges. Now consider subsets of the phase space where cells are synchronized, that is, there are distinct cells whose states take the same value; sometimes this is also referred to as cluster synchronization. Some synchronization patterns are robust, that is, they are dynamically invariant subsets of the phase space for any coupling functions. In our example, consider the set { x | x 1 = x 6 = x 5 , x 2 = x 4 }, where cells 1, 6, 5 as well as 2, 4 are synchronized. Note that this set is invariant under the flow of the above equations and the dynamics restricted to this space are given byẋ This illustrates some of the main questions we will address here: Given a set of dynamical equations, such as (1.3), what is the underlying hypergraph? Given a hypergraph H and an associated coupled cell hypernetwork, how can we identify the robust synchrony subspaces? Given a robust synchrony subspace, how can we describe the dynamics on the robust synchrony subspace as a coupled cell hypernetwork and how does this relate to the original hypergraph H? ✸ The main contribution of this paper is to develop the framework of coupled cell hypernetworks and apply this framework to analyze the existence and stability of synchrony in hypernetwork dynamical systems. While the dynamical equations are similar to those in [17,18,19], we explicitly discuss the role of the interaction functions Q k . Placed within the language of coupled cell networks, our approach allows to use the general ideas developed in [7] for the analysis of network dynamical systems with higher-order interactions. Specifically, the manuscript is organized as follows. Section 2 reviews some definitions and notation on directed weighted hypergraphs. The coupled cell hypernetwork formalism for coupled differential equations is introduced in Section 3. In Section 4 we define robust synchrony subspace for hypernetworks, describe those spaces for general hypernetworks and we relate them to the balanced colorings of the corresponding incidence digraph. In Section 5 we discuss a class of hypernetworks where we can relate stability of equilibria taking into account the nonpairwise interactions. We see already for this class of examples that the nonpairwise terms cannot be disregarded. We finish with Section 6 discussing the main points presented in this work and some questions that arise naturally.

Preliminaries on directed hypergraphs
In this section, we recall some notation and definitions on directed hypergraphs; see, for example, [20]. An hypergraph is a generalization of a graph where the graph edges are replaced by hyperedges that can join any number of nodes. In contrast to traditional directed hypergraphs, we allow for the tails to be multisets, i.e., a set that can contain an element more than once. Let #A denote the cardinality of a (multi)set A.
consists of a (finite) set of nodes C and a set of directed hyperedges E. A directed hyperedge e is a pair (T (e), H(e)), where the tail T (e) of e is a multiset of elements of C and the head H(e) of e is a subset of C; we assume that both T (e) and H(e) are nonempty. ✸ Note that, a directed hypergraph where any hyperedge e satisfies the conditions #T (e) = #H(e) = 1 is a standard directed graph.
In the above definition of directed hypergraph, we do not exclude the situation of having hyperedges e where the tail multiset T (e) has repetition of nodes. This fact is due to the association of hypergraphs with coupled cell hypernetworks and it will be clarified in Section 3.
If H = (C, E) is a hypergraph, we also write C(H) = C or E(H) = E to denote the set of nodes and hyperedges, respectively. . . , e m }. Let w : E → R be the weight function that associates a weight w j to each hyperedge e j , j = 1, . . . , m. The weight matrix W ∈ R n×m of H is the n × m matrix, where the ijth entry is the weight w j of the hyperedge e j if node i belongs to the head of the directed hyperedge e j , and 0 otherwise. A weighted directed hypergraph (H, W ) consists of H and a weight matrix W . ✸ Note that the definition of weight matrix of a directed hypergraph is distinct from that of the weighted adjacency matrix of a standard n-node directed graph, which is the n × n matrix, where the ijth entry is the weight w ij of the directed edge from node j to node i if there is a directed edge from j directed to node i, and 0 otherwise.
✸ To every directed hypergraph H can be associated a bipartite digraph D H , called the incidence digraph, Levi digraph, or König digraph of H, whose nodes are the nodes and the hyperedges of H; see, for example, [21]. Here, we generalize this concept to weighted directed hypergraphs (where the tails of the hyperedges can be multisets).  tail multiset T (e j ). The weighted incidence digraph D H of H is the weighted bipartite digraph with node set C ∪ E and edges such that: there is a directed edge from node i to the hyperedge e j with weight m j i if and only if i ∈ T (e j ); there is a directed edge with weight w j from the hyperedge e j to the node i if and only if i ∈ H(e j ). ✸ The adjacency matrix A D H of the weighted incidence digraph D H associated with a weighted directed hypergraph H has the block structure where W ∈ M n×m (R) is the weight matrix for H and the matrix T ∈ M m×n (R) describes the multiplicities of the nodes in the tail multisets of the hyperedges of H. The incidence digraph D H is represented in Figure 4 and its adjacency matrix is given by Remark 2.7. Note that, in network theory, the input set of a node in a directed network corresponds to the backward star of the node. ✸ We can define paths and connectivity in hypergraphs. A directed path of lenght q between the nodes v 1 and v q+1 is a sequence of nodes, v 1 , v 2 , . . . , v q+1 , and directed hyperedges, e 1 , e 2 , . . . , e q , where v 1 ∈ T (e 1 ), v q+1 ∈ H(e q ), and v j ∈ H (e j−1 ) ∩ T (e j ) for j = 2, . . . , q .
The nodes v 1 and v q+1 are said to be connected. An hypergraph is (weakly) connected if every pair of nodes in the hypergraph is connected by a path replacing all of its directed hyperedges with undirected hyperedges.
In the following we assume that all hypergraphs have nonempty node and hyperedge sets and are connected.

Coupled cell hypernetwork formalism
Weighted directed hypergraphs provide the backbone for the coupled cell hypernetwork formalism that we develop in this work. A hypernetwork is a weighted directed hypergraph, where each node i ∈ C comes with a phase space V = R d(i) and internal dynamics f i : V → V -we refer to a node with a phase space and internal dynamics as a cell. For simplicity, we assume that all cells are identical, i.e., V = R d and f i = f for all i. Thus, we will use the same symbol for each node/cell in a graphical representation of the network. In slight abuse of notation and terminology, we write (H, W ) for the hypernetwork, i.e., the weighted, directed hypergraph together with the data on the phase space, and use the words node/cell interchangeably.

✸
We will now define a set of dynamical equations that is compatible with the hypergraph H. For an hyperedge e ∈ E with weight w e we let k denote the cardinality #T (e). The evolution of cell i ∈ H(e) will be determined by a smooth coupling function Q k : V k+1 → V such that the evolution of cell i depends on x i and on the k variables x j with j ∈ T (e). More precisely, for a hyperedge e with tail T (e) of cardinality #T (e) = k, let x T (e) denote the k variables in the tail and write Q k (x i ; x T (e) ). We assume that Q k is invariant under permutation in the last k variables, the entries of x T (e) . Note that this implies that each hyperedge e ′ with #T (e ′ ) = k is of the same type: The interactions are governed by the same coupling function. At the same time, the strength of the interaction may be different since w e may be different from w e ′ . Definition 3.2 (Admissibility). A family Q = (Q k ), k ∈ N, of coupling functions as above is admissible for the hypernetwork (H, W ) if Q k = 0 for k ∈ B(H) and Q k = 0 otherwise. The collection of admissible family of coupling functions Q define the admissible cell vector fields for i ∈ C. ✸ Definition 3.3. Every admissible family of coupling functions Q for the hypernetwork (H, W ) and corresponding cell vector fields F i defines a coupled cell system where the state x i of cell i ∈ C evolves according toẋ . For convenience, we typically identify the dynamical system and the cell vector fields that define it. ✸ Example 3.4. Consider the hypergraph H on the right of Figure 1. For a collection of admissible family of coupling functions Q 1 , Q 2 , we have that the admissible cell vector fields are given bẏ

✸
From this perspective, a coupled cell hypernetwork characterizes a set of admissible coupling functions and admissible vector fields. However, distinct hypernetworks can have the same set of admissible coupling functions and even the same set of admissible vector fields.
Example 3.5. Consider the hypernetwork defined by the hypergraph on the left of Figure 5. For a collection of admissible family of coupling functions Q 1 , Q 2 , we have that the admissible cell vector fields are given byẋ Note that the directed hypergraph on the right of Figure 1   Example 3.6. Consider the hypernetwork defined by the hypergraph on the right of Figure 5. For a collection of admissible family of coupling functions Q 1 , Q 2 , we have that the admissible cell vector fields are given byẋ Observe that the two distinct directed hypernetworks of Figure 5 have the same set of admissible coupling functions and vector fields. ✸ Example 3.7. Consider the hypernetworks (H, W 1 ) (left) and (H, W 2 ) (right) of Figure 6. Thus, the same hypergraph H and different weighted adjacency matrices and, thus, different admissible vector fields. In fact, for an admissible coupling function Q 1 , we have that the admissible cell vector fields for (H, W 1 ) are given bẏ and the admissible cell vector fields for (H, W 2 ) are given bẏ Thus, we see that (H, W 1 ) and (H, W 2 ) have distinct set of admissible vector fields. ✸   Figure 7 are identical (equivalent) as coupled cell hypernetworks. For an admissible coupling function Q 2 , we have that for both coupled cell hypernetworks, the admissible cell vector fields are given bẏ ✸ Example 3.11. The two hypernetworks in Figure 8 are not equivalent as coupled cell hypernetworks. For an admissible coupling function Q 2 , we have that the admissible cell vector fields for the hypergraph on the left are given byẋ  Figure 7. Two identical (equivalent) coupled cell hypernetworks corresponding to two distinct weighted directed hypergraphs. Here, we are assuming all hyperedges with weight 1. where Q 2 is invariant under permutation of the last two coordinates. For an admissible coupling function Q 1 , we have that the admissible cell vector fields for the hypergraph on the right are given bẏ . That is, fixing the same cell phase spaces for the two hypergraphs, we have that the set of admissible cell vector fields for the hypergraph on the right is strictly contained in the set of admissible cell vector fields for the hypergraph on the left. ✸ Proof. Replace any hyperedge e ∈ E(H) with head set H(e) = {v 1 , . . . , v k } where k > 1, and weight w e by k hyperedges e j = (T (e), {v j }), for j = 1, . . . , k, each with weight w e j = w e . The set of admissible coupling functions and vector fields remain unchanged since they only depend on the tail of any hyperedge.

3.2.
Hyperedge-Maximality, hyperedge-minimality, and symmetries. In the previous section, we characterized a hypernetwork based on its set of admissible coupling functions/vector fields. In this section, we will now change perspective and focus on a specific choice of coupling function. Indeed, for a specific choice of coupling functions, we obtain a specific vector field. Conversely, we can assign a hypernetwork coupling to a dynamical system. Definition 3.14. A network dynamical system determined by x i ∈ V , i ∈ C, evolving according to is a coupled cell system for a hypernetwork coupling (H, W, Q) if X i = F i for an admissible cell vector field F i with respect to (H, W, Q) as defined in (3.6). ✸ Note that the assignment of a hypernetwork coupling to a dynamical system is not unique since the hypergraph and coupling function go hand in hand. Lemma 3.12 already indicated that even on the level of admissible vector fields, there are different hypergraphs that give rise to the same set of admissible coupling functions/vector fields. See Example 3.10 and Figure 7.
These hypernetwork couplings are identical. ✸ This implies that we can get equivalent hypernetwork couplings by splitting, or conversely combining hyperedges.
The hypernetwork couplings in Example 3.16 can be obtained by splitting/combining hyperedges.
Note that we do not require e to be distinct from e ′ 1 , . . . , e ′ k , we do not require {e ′ 1 , . . . , e ′ k } to be disjoint from E(H), nor do we necessarily have Q = Q ′ . If Q = Q ′ then the splitting/combining an hyperedge is purely structural.

✸
The following property is immediate: Then the hyperedge e is redundant. Note that redundancy here depends on the specific form of the coupling functions. ✸ To any arbitrary hypernetwork coupling we can associate a maximal and minimal dynamically equivalent hypernetwork coupling.
Definition 3.21. A hypernetwork coupling (H, W, Q) is hyperedgemaximal if no hyperedge can be split to obtain an equivalent hypergraph coupling. Conversely, a hypernetwork coupling is hyperedgeminimal if no hyperedges can be joined to obtain an equivalent hypergraph coupling structure. ✸ Note that without further assumptions, neither hyperedge-minimal nor -maximal associated hypernetwork couplings need to be unique: For example, if a hypernetwork coupling has an associated minimal hypernetwork coupling that has a single hyperedge e with weight w e then we get an infinite family of minimal hypernetwork couplings for  Note that we can always split hyperedges whose heads have cardinality greater than one. The following is an immediate consequence of Lemma 3.12: Lemma 3.24. Consider a coupled cell system with associated maximal hypernetwork coupling (H, W, Q). If (t, h) ∈ E(H) then #h = 1.
We now explore some straightforward consequences of equivalent hypernetwork couplings and how they relate to the symmetry of the coupled cell hypernetworks they define. Let S N denote the symmetric group of N elements that acts by permuting the node indices.
all cells are globally and identically coupled. These equations are S Nequivariant.
More generally we can make the following statement.
Then the coupled cell system is S k -equivariant where k = #A and S k acts by permuting the vertices in A.
Proof. By definition of a coupled cell system, Property (a) ensures that any node in A receives the same input. At the same time, Property (b) ensures that the input of any node depends in the same way on all nodes contained in A consequently, permuting nodes with indices in A does not affect the dynamical equations which proves S k -equivariance.
Of course Proposition 3.25 is a special case of the previous statement with A = C(H).

Synchrony in Coupled Cell Hypernetworks
Synchrony and synchrony patterns-where different nodes in the network evolve identically-is an essential collective phenomenon in network dynamical systems. Given a hypernetwork, what are the possible synchrony patterns for any admissible vector field? In the following we describe the synchrony patterns of a coupled cell hypernetwork and their associated balanced relations and quotient hypernetworks. 4.1. Input sets. As a first step, we generalize the concept of input equivalence relation for networks to the hypernetworks. For standard n-node directed graphs, Definition 3.2 of [1] introduces the concept of input equivalence of nodes. Roughly, two nodes c and c ′ are said to be input equivalent when besides the number of directed edges to c and c ′ is the same there is also a bijection between those sets of directed edges which preserves the edge types. In the above definition for two cells to be input equivalent, condition (iia) imposes that the sets of all cardinalities of the tail sets of the hyperedges of both cells must coincide. Moreover, condition (iib) says that, for a fixed cardinality of the tail set of a hyperedge of a cell, the summation of the weights of all the edges with the same tail set cardinality must coincide for both cells.
where Q 2 and Q 3 are invariant under permutation of the last two and three variables, respectively.
Observe that the set is flow-invariant for the above equations and the restriction of those equations to ∆ is given bẏ These equations are admissible by the hypernetwork on the right in Figure 10. This motivates the notion of a quotient hypernetwork; we make this explicit in the following section. ✸

4.2.
Robust synchrony subspaces. Consider a hypernetwork (H, W ) with n cells that take their state in V . Let ∆ ⊂ V n be a subspace of the hypernetwork total phase space defined by equality of cell states-a polydiagonal subspace. Define an equivalence relation ⊲⊳ on the cells of the hypernetwork in the following way: If x i = x j is an equality defining ∆ then i ⊲⊳ j. To highlight the underlying equivalence relation, we write ∆ = ∆ ⊲⊳ . We say that ∆ ⊲⊳ is a hypernetwork synchrony subspace when it is left invariant under the flow of every coupled cell system with form consistent with the hypernetwork, as defined above, that is for any admissible vector field. In slight abuse of notation and terminology, we will forget about the phase space and call ∆ a synchrony subspace of the weighted hypergraph (H, W ) if it is a hypernetwork synchrony subspace for any hypernetwork on (H, W ). Finally, if ∆ ⊆ R n is a polydiagonal subspace and K ∈ M n×n (R) leaves ∆ invariant, we also say that ∆ is a synchrony space of K. By Lemmas 3.12 and 3.19, we have the following result. Recall that for traditional coupled cell networks there is the notion of a balanced equivalence relation ⊲⊳ on the set of cells [1,2]. The balanced equivalence relations ⊲⊳ are in one-to-one correspondence with synchrony patterns: ∆ ⊲⊳ is a synchrony space for the network (that is, it is left invariant under the flow of every coupled cell system with form consistent with the network) if and only if ⊲⊳ is balanced. Motivated by the definition of balanced relation of a network introduced in [1,2] and generalized to the weighted network setup in [6,7], we now define balanced equivalence relation in the hypernetwork setup.
Consider a hypernetwork (H, W ) with set of cells C and set of hyperedges E. The hypernetwork is the union of constituent hypernetworks (H k , W k ) with identical set of cells C and hyperedges E k that contain the hyperedges whose tail sets have cardinality k 1 ; note that E k = ∅ if and only if k ∈ B(H) with B(H) as in (3.5). For simplicity, we will just write H k for (H k , W k ) (and H for (H, W )) in the following. Trivially, the input equivalence relation of H is a refinement of the input equivalence relation of every H k . (iii) We say that ⊲⊳ is balanced for the constituent hypernetwork H k if for every two distinct cells c, c ′ ∈ C such that c ⊲⊳ c ′ , the set of patterns determined by the hyperedges of the sets BS(c) and BS(c ′ ) coincide and each pattern has the same pattern weight on both cells. ✸ Definition 4.6. Consider a hypernetwork H with cells C, hyperedges E, and constituent hypernetworks H k as defined above. Let ⊲⊳ be an equivalence relation on C refining ∼ I . We say that ⊲⊳ is balanced if it is balanced for every constituent hypernetwork H k . ✸ Note that input equivalence is not always a balanced relation; this was already noted by Stewart [22,Section 6] for standard n-node directed graphs. That is, the coarsest balanced equivalence relation refines ∼ I but does not need not to coincide with ∼ I . See also Aldis [23] for the description of a polynomial-time algorithm to compute the coarsest balanced equivalence relation of a graph. Since it is a necessary condition for an equivalence relation on the nodes to be balanced is to refine ∼ I , we include that assumption at the above definition. The coarsest partition corresponds to the most synchrony that is possible.  In Figure 11, cells in the class 1 have white color, cells in the class 2 have blue color, and those in the class 4 have pink color. Consider the equivalence classes ordered as 1, 2, 4 . We have that ⊲⊳ determines two types of patterns, (2, 1, 0) and (1, 2, 0), for the hyperedges in both BS (4) and BS (14). The pattern (2, 1, 0) corresponds to a hyperedge with tail  Figure 11. The equivalence relation with three classes represented by the three colours is balanced for the hypernetwork.
set consisting of two white cells and one blue cell; the pattern (1, 2, 0) corresponds to a hyperedge whose tail set has two blue cells and one white cell. For cell 4, the incoming hyperedge with pattern (2, 1, 0) has weight 1 and the hyperedge with pattern (1, 2, 0) has weight 2. Proof. Let H be a hypernetwork which is a network, that is, the tail sets of all the hyperedges have cardinality 1. Thus H 1 = H. Given an  where each edge has weight w e = 1. The resulting hypernetwork is identical as a coupled cell system to the hypernetwork with underlying hypergraph H ′ = (C, E ′ ) such that the head H(e) of any hyperedge e ∈ E ′ has cardinality 1. Specifically, by splitting the head sets of hyperedges e 2 and e 4 we have By assumption in the beginning of this section, we will identify H = (C, E) with H ′ = (C, E ′ ) and drop the ′ . For the balanced coloring indicated by the shading of the nodes in Figure 1, the cells of the quotient are given by the equivalence classes The sets BS(1), BS(2), BS(3) are obtained from BS(1), BS(2) and BS (3), respectively, and thus all with weight equal to 1. Note that H/⊲⊳ is identical as coupled cell hypernetwork to the hypernetwork shown in Figure 1 to the right. ✸ Remark 4.14. The (somewhat nonstandard) convention to allow multisets as tails of directed hyperedges becomes essential in the coupled cell hypernetwork formalism presented in this work that considers generic features for all admissible vector fields simultaneously. By contrast, if one considers a specific hypernetwork coupling (H, W, Q), then one may be able to identify edges whose tail sets have cardinality k with edges with lower tail set cardinalities. For example, consider cells whose phase space is R and hypergraph coupling with  Figure 13. The quotient hypernetwork of the hypernetwork in Figure 11 by the balanced equivalence relation on the set of nodes whose classes are represented by the three colours.
Consider the decomposition of H into its constituent hypernetworks H k , for k = j 1 , . . . , j r according to the (positive and integer) cardinalities k of the tail sets of its hyperedges. Given two distinct cells c, c ′ such that   Figure 11 and the balanced equivalence relation ⊲⊳ presented in Example 4.8(i). Consider coupled cell systems consistent with H, where the cell phase space is V , the internal dynamics is given by f : V → V and the coupling by Since the equivalence relation ⊲⊳ is balanced, then the polydiagonal space is a synchrony space of H, that is, equations for H leave ∆ ⊲⊳ invariant. The restriction of those equations to ∆ ⊲⊳ gives rise to coupled cell systems consistent with the quotient hypernetwork H/⊲⊳ in Figure 13 with cells evolving according tȯ x 2 = f (x 2 ), ✸ Remark 4.19. Due to Lemma 3.19, the results in Theorem 4.17, concerning the restriction of the dynamcs to the synchrony subspace ∆ ⊲⊳ , apply to every hypernetwork obtained from the quotient hypernetwork Q by one (or more) purely structural combining hyperedge operations, since they are identical as coupled cell systems.  (i) Given an equivalence relation ⊲⊳ on C for H, we define the equivalence relation ⊲⊳ D on C ∪ E for D H in the following way: , for e i , e j ∈ E. with − → m(e) the pattern determined by ⊲⊳ on the hyperedge e.
(ii) Given an equivalence relation⊲⊳ on C ∪ E for D H , we define the equivalence relation⊲⊳ H on C for H through (a) c⊲⊳ H c ′ iff c⊲⊳ c ′ , for c, c ′ ∈ C. We say that the relation⊲⊳ H is the projection of the relation⊲⊳. ✸ Given the definition above, we have then the following result. Proof. Let (H, W ) be an hypernetwork with cells C and hyperedges E, and let D H be the associated incidence digraph with nodes C ∪ E.
(i) Let ⊲⊳ be a balanced equivalence relation on the set of cells C of the hypernetwork (H, W ) and consider the corresponding equivalence relation ⊲⊳ D on the set of nodes C ∪ E of the bipartite network D H , as in Definition 4.20. By definition, two nodes e i , e j ∈ E of D H such that e i ⊲⊳ D e j correspond to two hyperedges e i and e j of H that have the same pattern determined by ⊲⊳. Also, note that the input set of a node e i ∈ E of D H corresponds to the tail T (e i ) of the hyperedge e i in H. Thus: (a) for every two nodes e i , e j ∈ E of D H such that e i ⊲⊳ D e j there is a bijection between their input sets in D H preserving the ⊲⊳ D -classes.
Consider now two nodes c, d ∈ C of D H such that c ⊲⊳ D d, and thus with c ⊲⊳ d. Then, since ⊲⊳ is balanced, the pattern sets determined by the hyperedges of the sets BS(c) and BS(d) coincide and each pattern has the same weight on both cells. Note that the input set I B (i) of a node i ∈ C of D H is given by the backward stars BS(i) of i in H. We have then: (b) for any two nodes c, c ′ ∈ C such that c ⊲⊳ D c ′ , for every ⊲⊳ D -class, the sum of the weights of the edges in D H directed to nodes c and c ′ , from the nodes in that ⊲⊳ D -class, is the same. From (a) and (b), it follows that the equivalence relation ⊲⊳ D , as defined in Definition 4.20, is balanced. Thus, we have shown that, for every balanced equivalence relation ⊲⊳ for the hypernetwork (H, W ), we can associate a balanced equivalence relation ⊲⊳ D for the incidence digraph D H .
(ii) Let⊲⊳ be a balanced equivalence relation on the set of nodes C ∪E for the incidence digraph D H and consider the equivalence relation⊲⊳ H that is a projection on the set of cells For the projected equivalence relation⊲⊳ H on C for H obtained from⊲⊳ consider the usual polydiagonal subspace In terms of synchrony subspaces for the hypernetwork (H, W ) we have then the following result.
In terms of synchrony subspaces for the hypernetwork (H, W ) we have then that they can be obtained via the 'projection' of the synchrony subspaces of the adjacency matrix of the incidence digraph D H . Remark 4.24. A relevant consequence of the results in this section is that the existing results regarding balanced relations and synchrony spaces for networks can be used to obtain analogous results for hypernetworks. For example, the work of Aldis [23] with the description of a polynomial-time algorithm to compute the coarsest balanced equivalence relation of a graph and the work of Aguiar and Dias [24] describing an algorithm to compute the lattice of synchrony subspaces for the adjacency matrix of a network. ✸ where f : V → V , Q 1 : V 2 → V , Q 2 : V 3 → V are smooth functions and Q 2 is symmetric under permutation of the last two coordinates. Looking at the equations, we can conclude that the set of nontrivial synchrony subspaces for the hypernetwork H is given by Now, let us see how we can get this set of synchrony subspaces using the incidence digraph D H associated with H. The digraph D H is represented in Figure 4 and its adjacency matrix given by For an eigenvalue λ of a matrix let W λ denote the associated (generalized) eigenspace. Moreover, write v 1 , . . . , v k for the span of vectors v 1 , . . . , v k . The eigenvalues of the matrix A D H are λ ∈ {0, ±1, ±0.5 ± i0.866}; the algebraic multiplicity of λ = 0 is three and that of λ = ±1 is two. The corresponding (generalized) eigenspaces are where v 1 = (0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0) v 2 = (0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 1) and The polydiagonal subspaces given by equalities of cell coordinates and equalities of edge coordinates that are invariant by the adjacency matrix A D H arẽ These now relate to the synchrony spaces of H: We have that∆ 1 'projects into' the synchrony subspace ∆ 1 of H and∆ 2 and∆ 3 'project into' the synchrony subspace ∆ 2 of H. ✸ We stress that our results are valid for both unweighted and weigthed hypernetworks; the previous example can be seen as a hypernetwork where all weights are equal to one.
Remark 4.26. Note that there is no need to consider more than one adjacency matrix for the incidence digraph D H in order to separate the hyperedges with tails with different multiplicites since those hyperedges as nodes in D H cannot synchronize given that the row sum of the corresponding rows in the submatrix T of adjacency matrix A D H is different. ✸

Linearization and stability-a case study
In the previous sections, we considered the question what type of synchrony patterns can robustly exist for coupled cell hypernetworks and how they depend on the properties of the underlying hypergraph. We now consider linear stability of solutions on synchrony subspaces; asymptotic stability is crucial to actually observe synchrony patterns in real-world systems. We show that in a class of examples that linear stability may or may not depend on higher-order interactions.
Here we consider weighted directed hypernetworks (H, W ) with n nodes and directed hyperedges of the two types shown in Figure 14: There is an edge between nodes i, j with weight K ij and for each pair of nodes k, l in {1, . . . , n} there is a hyperedge ({k, l}, {i}) for i = 1, . . . , n with weight H kl . Note that we do not assume any relationship between the weights K ij of the pairwise interactions and the weights H kl between the nonpairwise interactions. For the remainder of this section, we fix a hypernetwork coupling through the coupling functions  [26]), that is, equations (5.9) with K = H and K is skewsymmetric. Here, K ij represents the effect of species j on the growth rate of species i. The dynamics of species i is determined by the fitness of species i given by n j=1 K ij p j and the average fitness for the system n k=1 n l=1 K kl p k p l ; this ensures that no species can increase in density without other species decreasing. The condition n i=1 p i = 1 ensures that total abundance conservation is maintained for all time. In this model terms of the form K ij p i p j represent pairwise interactions between the species i and j and n k=1 n l=1 K kl p k p l represents an average of nonpairwise interactions between all the species.
In [27], it is shown that for a skew-symmetric n × n matrix K is skew symmetric the system has a unique equilibrium solution p, which is linearly neutrally stable. For a skew-symmetric matrix K, the quadratic form w → w T Kw is null and with K = H the system (5.9) reduces tȯ p i = (Kp) i p i (5.11) for i = 1, . . . , n. Chawanya and Tokita [27] reports that the condition of skew symmetry of K (on the interactions between the species) can be used to yield and stabilize a large complex ecosystem. The antisymmetry model assumption is based on the fact that many species interact with each other in prey-predator or parasitic relationships. ✸ We can make the following two observations. We show two examples of system (5.10), one with no nonpairwise interactions and one with nonpairwise interactions, admitting an equilibrium whose stability does depend on the nonpairwise interactions terms.
Examples 5.4. Consider the system (5.10) where n = 4 and Note that K is a skew symmetric matrix. The eigenvalues of K are λ = 0 (double) and a pair of nonzero imaginary eigenvalues λ = ±i √ 11/2. Moreover, W 0 = (1, 1, 1, 1), (0, 2, 1, 0) . (a) Assume that in (5.10) there are no nonpairwise interactions, that is, H = 0. We have that p * = 1 4 (1, 1, 1, 1) is an equilibrium of the system (5.10) with stability determined by K (by Lemma 5.3), that is, the equilibrium p * has neutral linear stability in the sense that all eigenvalues have zero real part.  Figure 15. The equivalence relation with three classes represented by the three colours is balanced for the network.
is a synchrony space for K and thus, by Lemma 5.2, also for the system (5.13). Moreover, which has eigenvalues 0, − 1 2 , and 1 4 (1 ± i √ 3). That is, it has the same stability as for the system without nonpairwise interactions, H = 0. ✸

Discussion
Here we developed a framework for coupled cell systems with higherorder interactions. In contrast to other approaches to dynamics on hypergraphs-including [17,19]-our framework allows for directionality of the interactions and coupling weights. The framework is restricted by the assumption of homogeneity in the kth order coupling: The interaction is mediated by a single coupling function Q k for any edge of tail size k. These assumptions do shape the set of admissible vector fields. Recall the hypernetwork of Example 4.8(ii), which is depicted in Figure 12. As an example, the admissible evolution equations for nodes 4 and 12 take the shapė x 4 = f (x 4 ) + Q 3 (x 4 ; x 1 , x 2 , x 8 ) + Q 3 (x 4 ; x 5 , x 6 , x 7 ), x 12 = f (x 12 ) + Q 3 (x 12 ; x 1 , x 2 , x 3 ) + Q 3 (x 12 ; x 9 , x 10 , x 11 ) .
By contrast, if we forget the hyperedge structure and consider the related network shown in Figure 15 then the equations for cells 4 and 12 in the formalism of Golubitsky, Stewart and collaborators [1,2] have the formẋ 4 = g(x 4 ; x 1 , x 2 , x 5 , x 6 , x 7 , x 8 ), x 12 = g(x 12 ; x 1 , x 2 , x 3 , x 9 , x 10 , x 11 ), where g is invariant under permutations of the last six arguments. Even though the combinatorial representation of the equations is a network (a directed graph), the admissible vector fields that are determined by the interaction function g can have nonlinear dependencies between the cell coordinates x k . By contrast, in the additive input setup [4,5,6,7] no nonlinear interactions beyond pairs of cells are possible and the admissible equations for cells 4 and 14 have the forṁ x 12 = f (x 12 ) + h(x 12 , x 1 ) + h(x 12 , x 2 ) + h(x 12 , x 3 ) + h(x 12 , x 9 ) + h(x 12 , x 10 ) + h(x 12 , x 11 ).
The admissible vector fields of our framework are richer than the additive setup. Moreover, they explicitly capture higher-order interaction structure, which is only implicit in the classical formalism of Golubitsky, Stewart, and collaborators but important from a dynamical point of view; cf. Section 5.
What is an appropriate combinatorial structure to encode higherorder interactions in network dynamical systems (cf. [9])? The framework developed above is phrased in terms of (directed) hypergraphs. First, the hypergraphs employed are nonstandard: The tails of each hyperedge is a multiset rather than a set. This is crucial to define a quotient of a hypernetwork without making further assumptions on the coupling functions as arguments on the synchrony subspace can appear multiple times. Second, different hypergraphs can represent the same coupled cell hypernetwork. This is due to the fact that hyperedgeheads can contain more than one element which may allow to easily identify symmetries (cf. Proposition 3.25).
It is worth pointing out that in the formalism developed above we typically consider all admissible vector fields at the same time. More specifically, we ask: What are the dynamical features of all ordinary differential equations (ODE) that are compatible with the hypernetwork structure? This elucidates the constraints network structure imposes.
For example, Theorem 4.16 allows to translate structural properties (balanced relations on a hypergraph) into dynamical properties (any ODE consistent with the hypernetwork will have a particular synchrony subspace). Consequently, these properties are not specific to any choice of coupling function. While this is the same approach as in traditional coupled cell systems, the approach is in contrast to some applications where a fixed coupling function is considered: A specific coupling function may be imposed by a particular physical system. But a nongeneric choice of coupling function can lead to nongeneric dynamical behavior and nonproper hypernetwork couplings (Definition 3.22).
The importance of higher-order interactions in network dynamical systems has repeatedly been highlighted. The framework presented here bridges coupled cell systems and higher-order interaction networks. Specifically, it allows to characterize synchrony patterns (whether global or localized/clustered). While other approaches are possible, our framework strikes a balance between generality and results that can elucidate synchronization phenomena in real-world systems.