Learning the effective order of a hypergraph dynamical system

Dynamical systems on hypergraphs can display a rich set of behaviors not observable for systems with pairwise interactions. Given a distributed dynamical system with a putative hypergraph structure, an interesting question is thus how much of this hypergraph structure is actually necessary to faithfully replicate the observed dynamical behavior. To answer this question, we propose a method to determine the minimum order of a hypergraph necessary to approximate the corresponding dynamics accurately. Specifically, we develop a mathematical framework that allows us to determine this order when the type of dynamics is known. We use these ideas in conjunction with a hypergraph neural network to directly learn the dynamics itself and the resulting order of the hypergraph from both synthetic and real datasets consisting of observed system trajectories.

Dynamical processes on hypergraphs 1,2 have recently received significant attention.Examples include synchronisation 3,4 , consensus dynamics [5][6][7] , epidemic spread [8][9][10] , random walks 11,12 , label propagation [13][14][15] and social contagion 16,17 .Generally, it has been found that for dynamical systems supported on hypergraphs instead of graphs, important characteristics of the dynamics can change.For example, contagion processes on hypergraphs can have different epidemic thresholds for the outbreak of an epidemic 9,10 , or group reinforcement effects can significantly alter the outcome of opinion formation on hypergraphs [17][18][19] .Hypergraphs have also gained interest as extensions of graph neural networks and been used for a variety of applications ranging from 3D shape retrieval 20 and pose estimation 21 in computer vision, to group recommendation tasks 22 , or medical applications such as cancer tissue classification 23 .
Although both aspects often appear lumped together in the mathematical equations describing a dynamical system on a hypergraph, from a modeling perspective, it is important to distinguish between (i) the topological relations constraining the possible interactions in a system (as encoded in a hypergraph) and (ii) the model of the local multi-way dynamics occurring on each hyperedge (e.g., epidemic spread, diffusion, synchronisation).It is the interplay of both aspects that leads to the (possible) emergence of higher-order effects 2,18,24 .
For instance, if we consider a linear dynamics on a hypergraph, we cannot expect any higher-order effects to emerge.As any linear map between finite-dimensional spaces can be represented by a matrix, for linear dynamics we may always find an equivalent graph-based, pairwise interaction dynamics that models the system exactly, by identifying the matrix with an effective weighted graph 5,6 .In practice, this means that it is possible to use a graph instead of a hypergraph-based dynamical system, as long as linear dynamics are considered.Similarly, it has been shown that various formulations of semisupervised and unsupervised spectral clustering problems on hypergraphs lead to an effective graph-theoretic problem [25][26][27] .
More generally, we can envision that, depending on the local dynamics, we will be able to rewrite a hypergraph dynamical system on a hypergraph of general order k as a dynamics occurring on a hypergraph with hyperedges of order at most p ≤ k.Such a simplification to lower-order relations is relevant as the use of hypergraphs presents several challenges: most prominently, since the number of hyperedges can grow combinatorially with the number of nodes, the use of hypergraph models can be computationally very expensive.This is particularly relevant for large-scale systems.
Yet, in practice, we typically know neither the exact set of relations between the entities, nor the analytical form of the local interaction dynamics, but are only given observational data, e.g., in the form of trajectories.Many different approaches have been explored to approximate such time series data.For example, it has recently been proposed to treat neural networks as models equipped with a continuum of layers 28 .This view allows a reformulation of the forward pass of a neural network as the solution of an initial value problem of an ordinary differential equation (ODE).Deep learning architectures that use this reinterpretation of the forward pass are called Neural ODEs and are useful for building continuoustime time series models.Neural ODEs have also been recently generalised to graphs 29,30 .Discovering the equations governing a dynamical system from measurements is an important problem, tackled by a large body of literature [31][32][33][34] .The main challenge here is to find a model which is complex enough to describe the existing data but not too complex to avoid overfitting.
For a dynamical system on a hypergraph, there are many possible ways to abstract an observed distributed dynamics.For instance, we may model it as emerging from a simple local dynamics on a rather complex hypergraph.Alternatively, we may consider a more complicated local dynamics interacting via a more constrained set of relations between entities.This raises the question of whether we should include the complexity of our model in the topology of the interactions or in the model of the dynamics: what multi-way relations do we need to encode in our model in practice?
In this paper, we thus introduce the concept of the dynamical order of a hypergraph dynamical system, which measures how complex the dynamics are.Based on this dynamical order, we can determine the effective order of the hypergraph dynamical system, which is given by the minimum order of a hypergraph necessary to accurately represent the corresponding dynamics.In particular, we propose an analytical framework that allows us to derive the dynamical and the effective order when the functional form of the dynamics is given.Furthermore, we propose a hypergraph neural network architecture that allows learning the hypergraph dynamics and the resulting effective order directly from data, which we test on both synthetic and real data sets.In conclusion, we present an effective method to reduce the complexity of a hypergraph dynamical system and to learn its representation from data.

Reducibility of hypergraph dynamical systems
To illustrate our ideas of the dynamical and effective order of a dynamical system on a hypergraph, let us start with a Kuramoto-type dynamics 35 on a hypergraph as a concrete example.Kuramoto oscillator dynamics have been applied to various synchronization phenomena of phase oscillators 36 , ranging from power networks 37 to brain activity 38 .Several works have been working on its generalisation to simplicial complexes 3,4 and hypergraphs 39 .We will compare two different formulations of Kuramoto dynamics on a hypergraph here.
Consider a hypergraph H consisting of a set V = {1, 2, . . ., N} of N nodes, and a set E = {E 1 , E 2 , . . ., E M } of M hyperedges.Each hyperedge E α is a subset of the nodes, i.e., E α ⊆ V for all α = 1, 2, . . ., M. Each hyperedge may have a different cardinality |E α | .We define the topological order k = max α |E α | as the cardinality of the largest hyperedge.Let x(t) ∈ R N be the vector of dynamical state variables of the nodes at time t.For simplicity, we will suppress the time dependency and write x i = x i (t) for the ith component of the node state vector in the following.
There are several ways in which the well-known Kuramoto dynamics can be extended to hypergraphs.First, inspired by Adhikari et al. 39 we can write the time evolution of node i as Since the non-linear sine function acts within each entire hyperedge, for general node states x i it is not clear how to further reduce the system to a dynamical system on a lowerorder hypergraph.An alternative formulation for Kuramoto oscillator dynamics is where within each hyperedge, every pair of nodes interacts.Here, B ∈ R M×N represents the node-to-edge incidence matrix Hence, even though we are dealing here with a nonlinear dynamics on a hypergraph, the dynamics in each hyperedge is pairwise and we can reduce it to a pairwise network dynamics: the hypergraph topology simply scales the system by A = B ⊤ B. Overall, our example of Kuramoto oscillator dynamics shows that whether the hypergraph can be projected onto a lower-order system depends on the form of the dynamics supported on each hyperedge: the dynamics is reducible if the dynamics on the hyperedges can be rewritten as a linear combination of lower-order functions.More generally, for certain functional forms the non-linear dynamics can always be reduced to a lower-order hypergraph system.In the Supplementary Materials in Section 1, we show the linear-like properties that dynamics must have in order to be reducible to a network dynamical system.
In the following, we formalise this form of dynamical reduction by distinguishing between the topological order of a hypergraph and the dynamical order of a dynamics.Combining topology and dynamics to a hypergraph dynamical system, its effective order is given by the minimum of these two orders.

Dynamical order and effective order
Here we introduce a framework to separate the impact of topology and dynamics for a general hypergraph dynamical system.To this end, we group the edges of a hypergraph H by size, so that its edge set is given by E = E 2 ∪ • • • ∪ E k where E d denotes the set of all hyperedges of size d.We refer to the values in the state vector x(t) of nodes in the hyperedge E at time t excluding node i with { {x j (t) | j ∈ E, j ̸ = i} }.This is a multiset, which we denote by double braces -i.e. a set where elements may occur multiple times -since the node may have the same values.
Formally, we now define a hypergraph dynamical system as follows.
Definition 1 (Hypergraph dynamical system).Consider a family of multivariate update functions ( f d ) d∈1...K where each function f d (y 1 , . . ., y d ) acts on d variables, and f d is symmetric in its last d − 1 arguments, i.e., any permutation the last d − 1 variable does not change the function.Together with the topology of a hypergraph H with topological order k ≤ K, these update functions f d describe the local dynamics on each hyperedge.The dynamics on a node can then be written in the form of: Hypergraph dynamics model and visualisation of the topological, dynamical and effective order.When modelling complex dynamical systems, it is necessary to take into account both the system's topology and the local dynamics which are acting on the topology.We introduce a hypergraph dynamics model which separates the impact of topology and dynamics.In A, we show the update process of a single node according to the full hypergraph topology.All nodes in the hypergraph are updated synchronously.In particular, the values of the nodes that are part of hyperedges of size d are collected.Then an update function on these hyperedges which defines the dynamics is computed to obtain the update and then projected back to the node space by summing over each update.With the help of this framework, we can define and analyse different orders of the system, which we visualise in B. The topological order is given by the size of the largest hyperedge in the hypergraph.To derive the dynamical order, we look at the update functions on a hyperedge of topological order 4 in more detail.In particular, the update resulting from a hyperedge of topological order 4 can consist of a linear combination of pairwise or three-way functions which results in dynamical order 2 or 3, respectively.We can then derive the effective order of the hypergraph dynamics which is the minimum of the topological and dynamical order.The illustration highlights that the topological order of a hypergraph is only an upper bound, not a lower bound, on the effective order of the system.

2/13
These update equations for i = 1 . . .N define a hypergraph dynamical system on H .
Note that Equation (4) separates the hypergraph topology from the dynamics on the edges: For each size d, a set of hyperedges defines the topological coupling of entities within the system and the update functions f d act on the d state variables associated to the nodes within each of these hyperedges.This is visualised in Figure 1 A. Details of the mathematical framework can be found in Section 5.As we are primary interested in interacting dynamics, we will for simplicy set f 1 = 0 in the following, even though all our arguments can be naturally extended to the nonzero case.
Most commonly used hypergraph dynamical processes considered in the recent literature can be written in the above form.For example, the dynamics of the Kuramoto oscillator model on hyperedges in Equation ( 1) can be written as above if we define f 1 = 0 and: for d = 2, . . ., k.
We now introduce the dynamical order of such a hypergraph dynamical system, determined by the form of the update functions f d .
Definition 2 (Dynamical order).Consider a hypergraph dynamical system.Each update function f d can be generically decomposed into a sum of functions φ p d of only p ≤ d variables: where s denotes a multi-set of variables { {y 2 , . . ., y d } } of size d −1 and v denotes a non-empty subset of s of dimension p−1.
We call the minimal p ≥ 21 for which this decomposition is possible for all update functions f d the dynamical order p dyn of the system.
Trivially, by choosing p = d we have that φ d d (y i , s) = f d (y i , s) such that the above decomposition always exists.However, for specific functional forms f d , it can be possible to decompose higher-order interactions into a combination of p-ary interactions.Let us illustrate this with a concrete example.Consider the family of update functions defined for all d via f d (y 1 , . . ., y d ) = log(y 1 . . .y d ), i.e., f 1 : y 1 → log(y 1 ), f 2 : (y 1 , y 2 ) → log(y 1 y 2 ), and so on.By the properties of the logarithm we can write f d for any d ≥ 2 as: (6) and thus the dynamical order is p dyn = 2.
From our above discussion, it should be apparent that for a dynamical order p dyn , the dynamics on the larger hyperedges always consist of linear combinations of functions of order p dyn .Hence, one may rewrite Equation (4) more compactly in terms of dynamics on (effective) hyperedges of size p dyn .This is akin to a higher-order equivalent of a clique expansion. 40.
This leads us to the definition of the effective order, which combines topological and dynamical order (see fig. 1B).
Definition 3 (Effective order).Consider a hypergraph dynamical system which consists of a hypergraph H with topological order k and a dynamical system ( f d ) d∈N with dynamical order p dyn .The effective order p min of the system is the minimum order of a hypergraph necessary to describe the corresponding dynamics accurately.It is bounded by Indeed, if k = p dyn , the dynamics on none of the hyperedges decompose and the effective order of the system is given by the topological order.Consequently, we cannot reduce the hypergraph dynamical system any further.This was the case in eq. ( 1).Thus, the effective order of this hypergraph Kuramoto model is p min = k.
In contrast, if p dyn < k, the dynamics on the hyperedges larger than the dynamical order factorise into p dyn -ary functions, and the hypergraph dynamical system can be projected onto a (new, effective) hypergraph of order p min = p dyn < k An example for an effective order of p min = 2 is given by eq. ( 2) where the update functions for d = 2, . . ., k are: This dynamics has a dynamical order of p dyn = 2 < k and the dynamics thus factorise into pairwise functions on all hyperedges so that they are effectively network dynamics.Overall, the exchange of the order of the sum and the sin function is key to make the dynamical order of eq. ( 1) and eq.( 2) differ, leading to different effective order of the systems.Though conceptually useful, the derivation of the effective order as described in this section requires knowledge of the analytical form of the dynamics.In the next section, we thus present a method to learn the corresponding functions by deriving a hypergraph neural network model.Using this computational model, we can learn the dynamics from dataand thus implicitly the effective order of the system -without prior knowledge on the functional form of the dynamics.

Learning local dynamics and effective order from data
In this section, we present a method to learn the local dynamics and thus the dynamical and effective order of the system directly from data.
We train our HyDy-GNN using datasets X = {(x (1) , ẋ(1) ), . . ., (x (S) , ẋ(S) )} consisting of pairs (x (i) , ẋ(i) ) of state vectors x (i) and gradients ẋ(i) of the dynamics.Our objective is for our model to output the state derivative vector ẋ, when the state vector x is given as input: (9)   where fd describes the (estimated) update functions of our HyDy-GNN .In particular, we use empirical risk minimization with a regularised L 1 loss function and fit the parameters θ p model of HyDy-GNN , which consist of the vector of all weights of the MLPs, by minimising the absolute difference of the prediction outputs and the provided target values of the training set x: (10)   The hyperparameter λ of the L 2 penalty term is optimised through a hyperparameter search.
A technical detail we need to take care of in this context is that we must assign an arbitrary ordering to the values in the multiset s (and thus also in the subsets v).To alleviate this problem, we train the MLPs based on all possible orderings of the values in the multiset s (or the subsets v) and take the mean to obtain the approximation of the p-ary function φ p d (see Section 5).

Finding the effective model order
Using the above outlined learning paradigm, we fit a series of HyDy-GNN s with increasing orders p model ≤ k (i.e., we fit models of orders up to k).Note that, if the chosen model order is less than the true dynamical order of the system, the HyDy-GNN will not fit well, since the MLPs of order p model < p dyn cannot yield a good approximation to the nonlinear dynamics induced by the larger hyperedges.However, if the model order is equal to or greater than the dynamical order of the system, the dynamics will be learned accurately.As we take the topological order k as an upper bound, the true (a priori) dynamical order of the dynamics may be larger than p min .However, if we find a model with an effective order that is smaller than the topological order, this means we have found a model order that approximates the true dynamics sufficiently well.This may happen, e.g., if the dynamics are inherently reducible to a smaller order, i.e. their dynamical order is smaller than the topological order of the hypergraph.
Hence, to find the effective model order we need to find the model order that approximates the dynamics sufficiently well, while being as small as possible.Depending on the application considered, there are multiple possible ways to select such an order.For instance, we can select the order by manual inspection, or by setting a particular approximation threshold for the approximation quality of the dynamics, and choosing the smallest model order that achieves this threshold.
While this is not the main focus of our work here, we discuss in the following a simple scheme to automate this process of model selection.Specifically, we define a modelcorrected performance score to select the model with the lowest order p model that produces accurate results on a given dataset, as follows: where L max = max p model L(X , θ p model ) .We provide more details on the derivation of this model-corrected performance score in the Supplementary Materials in Section 2.

Results
In this section, we show that we can indeed use our framework to learn the update function on the hyperedges from data.We then show that we can use these results to infer the effective order of the hypergraph dynamical system.

Datasets
We perform our experiments on synthetic Erdős-Réyni hypergraphs and on a real hypergraph of high school students' contact patterns (SocioPatterns dataset).Both types of hypergraphs have topological order k = 4. On these hypergraphs, we simulate (higher-order variants of) four common linear and non-linear model: Kuramoto dynamics (synchronisation), SI dynamics (epidemic spread), multi-way consensus (MCM) dynamics (opinion dynamics), and linear diffusion.For each of these four dynamics we consider different variations of the dynamics, where we restrict the update functions f d to be a sum of at most p-ary update functions.Details on the generation of the hypergraphs and the dynamics can be found in Section 5.
We look into two possible scenarios which may occur in the real world: In the first scenario (I) we observe the state vectors and their derivatives for particular time-points.In the second scenario (II), we observe a set of trajectories of the dynamics.
For each of these settings, we create a corresponding dataset which we split into a training and test set.For the setting (I), each dataset consists of 500 independently sampled random initial state vectors and the corresponding derivatives of the dynamics at that point (see Section 5.3).This data can directly be used for training of the HyDy-GNN via empirical risk minimization (see eq. ( 10)).For setting (II), we generate a set of 25 trajectories, whose initial state vectors are uniformly sampled (for details see section 5.3).Starting from the random initial conditions, we simulate 100 time steps with a time delta of ∆ = 0.01 via a simple forward Euler scheme.From these trajectories we then extract a training sequence for our model, by approximating the time-derivatives ẋ via a simple forward time difference at each time point of the trajectories.We thus create a training set of 2500 tuples of state vectors and associated derivatives which we can use to train our model.We remark that in contrast to scenario (I) the training samples in scenario (II) are in fact not independent, as they are coupled over time via the governing equations of the dynamics.To differentiate between the two scenarios we will refer to the training data in the first setting (I) as point-based training data, and the training data in the second setting (II) as trajectorybased training data.

Numerical experiments: Learning the update functions
As a first step, we test whether a HyDy-GNN can accurately learn the dynamics from data.To this end, we first consider all dynamics restricted to dynamical order p dyn = 2 such that HyDy-GNN of any order should be sufficient to fit the dynamics.We then train a HyDy-GNN with order p model ∈ {2, 3, 4} to see whether we can indeed learn the update function well.
In Figure 2 we see the true and learned update functions of a hypergraph neural network trained on a point-based training set of (a) pairwise Kuramoto dynamics (p dyn = 2) and (b) pairwise MCM dynamics (p dyn = 2).We plot the true update functions as a function of the first and second input variable in the left columns.We show the learned dynamics for our HyDy-GNN s with model orders p model ∈ {2, 3, 4} (from top to bottom) in the middle columns.When p model > 2, we plot the output of the MLP 2 d with input dimension 2 for simplicity.The absolute error of the approximation, i.e., the absolute difference of the individual function values of the ground truth and learned dynamics, is shown in the right columns in Figure 2 (a) and (b).As expected, we observe a good approximation of the ground truth dynamics for all HyDy-GNN orders irrespective of the types of dynamics.In the Supplementary Materials in Figure 1, we show that for other dynamics we obtain equally good results.

Synthetic Data
We now demonstrate how we can use our framework to learn the effective order of an observed dynamics.Specifically, we examine hypergraph dynamical systems with Kuramoto, SI, MCM and Diffusion dynamics defined in terms of p-ary update functions with p ∈ {2, 3, 4} (see Section 5.3).Note that, in general, p does not have to be the dynamical order as we only ensure that there exists some decomposition of the dynamics into p-variate functions, but this decomposition does not have to be the minimal one.In fact, for the diffusion dynamics we always have p dyn = 2 due to linearity.However, for the other update functions chosen here, the dynamical order is indeed always equal to the chosen decomposition, i.e., p dyn = p in case of the Kuramoto, SI and MCM dynamics, as we show in the Supplementary Materials in Section 1.As the topological order of all considered hypergraphs is k = 4, the effective order is thus equal to the dynamical order for all hypergraph dynamical systems.
The results of the experiments on synthetic datasets are shown in Figure 3, where we compare the learning results for the derivative datatset in Figure 3 (a)-(b) and the trajectory dataset in Figure 3 (c)-(d).While we always trained the model with the pointwise L 1 -loss, as discussed in section 2.1, we evaluate the performance with two measures: (a.1) pointwise mean absolute error (MAE) (c.1) trajectory MAE (mean MAE over a trajectory) (see Supplementary Materials for details).Both MAEs are calculated on a cross-validation with 10 folds.Our results show that we are indeed able to learn the effective order from data.In Figure 3 (a.1),we show the performance of the different models with p model ∈ {2, 3, 4}, which is given by the lowest mean absolute error (MAE) of a cross-validation with 10 folds.
We see in Figure 3 that the models, whose order is equal or larger than the effective order of the hypergraph dynamical system, clearly outperform the models with an order smaller than the effective order.As the dynamical order of the diffusion dynamics is p dyn = 2 for all hypergraph diffusion systems, all models can adequately approximate the dynamics in this case.We plot the model corrected performance in Figure 3 (a.2).We see that our performance score almost always selects the model with p model = p min , which demonstrates that our method is able to infer the effective order for the synthetic training data set.In the case of SI dynamics, we find that the observed dynamics is also approximated well by an order 3 model, which results in the selection of p min = 3 as effective order in this particular case.
In Figure 3(c), we evaluate the performance of our model   We observe a very good approximation to the ground truth dynamics for all model orders and both types of dynamics.This implies that our framework allows us to learn the update function of a dynamical system from data.As this is true for all model orders, we can additionally conclude that the framework can separate the effects of topology and dynamics, as it accurately captures the pairwise update functions even for higher model orders.
for long-term trajectory predictions based on the trajectorybased data set.In Figure 3 (c.1),we thus plot the trajectory MAE of ground truth and predicted trajectories.We see that the results are similar to (a).However, in this case the model-corrected performance score is calculated by substituting L MAE into eq.( 11), as described in the Supplementary Materials in Section 3. The results also suggests the selection of an effective order of p min = 3 for MCM and SI.An evaluation using the pointwise MAE which leads to similar results can be found in the Supplementary Materials in Figure 2.
In Figure 3 (b) and (d), we show an example of the longterm prediction of trajectories with 200 time steps for the two datasets.Long-term prediction here means that we provide one initial condition and integrate over 200 time steps with a time delta of ∆ = 0.01 via a forward Euler scheme.In particular, we plot ground truth Kuramoto trajectories (blue) with p dyn ∈ {2, 3, 4} (top to bottom) and the predicted trajectories given by a HyDy-GNN with p model ∈ {2, 3, 4} (left to right).
We find that only a HyDy-GNN with large enough order approximates the long-term behaviour of the dynamics well, as indicated also by the results in (a) and (c).Although the trajectory training set derived from a smaller amount of different hypergraph topologies (only 25 different hypergraphs are considered instead of 500 in the derivative dataset), it provides a better approximation of the trajectories (e.g., compare (b) and (d) for p model = 4, p dyn = 3).

Real-word datasets
Similar to the results on the synthetic datasets, the results in Figure 4 show that using our method we can learn the effective order of dynamics unfolding on also on a real-world hypergraph of high-school student contacts.Here, we consider the three interaction dynamics which would be expected on a contact dataset, epidemic spread (SI), opinion dynamics (MCM) and diffusion (e.g.diffusion of ideas or information).For diffusion, a linear dynamics, the inferred effective order is correctly inferred as p min = 2.For nonlinear dynamics such as epidemic spread and opinion formation, the best-performing model corresponds to an effective order p min = {2, 3, 4} rather than the order of the hypergraph.This highlights the importance of carefully considering the type of dynamics when working with real-world hypergraph structures.Although this will not always hold, a lower-order hypergraph will be sufficient to capture the relevant dynamics in many cases.Our method enables researchers to derive how much they can simplify their models from their specific datasets of interest.

7/13
(1) (2)  and show that we are indeed able to learn the effective order of a hypergraph dynamical system from data.In (a.1), we display the mean absolute error (MAE) of cross-validation with 10 folds for all three models with p model ∈ {2, 3, 4} which were trained on the respective dynamics datasets (top row).We can see that only if the model order is larger or equal to the effective order of the system (which is p min = p for Kuramoto, SI and MCM and p min = 2 for Diffusion), the dynamics can be learned by the HyDy-GNN .In (a.2) we plot the model-corrected performance which is maximised for a model that optimises for both accuracy and model complexity.We see that in almost all cases our performance score selects the model with p model = p min .This shows that our method is indeed able to infer the effective order for the synthetic derivative training set.Only in the case of SI dynamics, we can see that dynamics with effective order p min = 4 can already be roughly captured by an order 3 model, which results in the selection of p model = 3 as the optimal order.In (c), we evaluate the performance of our model for long-term trajectory predictions based on the trajectory training set.In (c.1), we thus plot the MAE of ground truth and a predicted trajectory.This evaluation thus includes accumulated errors.We see that the results are similar to (a), whereas in this case the model-corrected performance score also suggests p min = 3 for a dynamics with p dyn = 4 for MCM and SI, even though the model of order 4 leads to a lower MAE.If a more accurate model is preferred, one could reweight the model-corrected performance score accordingly.In the right column (b) and (d), we show an example of the long-term prediction of trajectories with 200 time steps for the two datasets.In particular, we plot ground truth Kuramoto trajectories (blue) with p ∈ {2, 3, 4} (top to bottom) and the predicted trajectories given by a HyDy-GNN with p model ∈ {2, 3, 4} (left to right).We can see that only models which are equal to or larger than the effective order of the system (which is given by p min = p for Kuramoto dynamics) can capture their long-term behaviour.Moreover, the trajectory training set provides slightly better results.This is because the trajectory dataset is biased towards the long-term behaviour of the dynamics.As a result, it is more suitable for this type of long-term trajectory prediction.Generally, the results confirm that our method is capable of learning the effective order, both for the derivative and trajectory training set and the resulting models can capture the long-term behaviour of the system.In particular, we simulate epidemic spread (SI), opinion dynamics (MCM) and diffusion with order p ∈ {2, 3, 4} as all of these dynamics can be expected to appear on a contact pattern hypergraph.In (1), we again show the MAE of cross-validation with 10 folds for models with p model ∈ {2, 3, 4} which are trained on derivative datasets.As in the case of the synthetic datasets, we can see that only if the model order is larger or equal to the effective order, the trajectories can be learned by the hypergraph neural networks.This is undermined by the results of the model-corrected performance in (2).We see that the model order of the best model corresponds to the effective order of the system.Again, only in the case of SI dynamics, we can see that p dyn = 4 can already be roughly captured by a model with p model = 3.Overall, the results highlight the importance of carefully considering the type of dynamics when working with real-world hypergraph structures.

Discussion
We presented a framework to infer a hypergraph dynamical system, trading off topological and dynamical complexity.We used this framework to derive an effective hypergraph dynamical system for a given dynamics both analytically and in a data-driven way, based on observational data.In particular, using a neural network as a flexible way to approximate the local interaction dynamics, we were able to accurately learn the hypergraph dynamics and reduce the hypergraph order, while respecting the observed dynamics.In this context, finding an effective dynamical order that is smaller than the topological order of the hypergraph, indicates that we can "prune down" some hyperedges to be of smaller order as they do not play a role for the specific dynamics of interest.This is interesting from a point of view of model complexity reduction.Moreover, the learned models are capable of predicting the long-term behaviour of dynamical systems.We believe that our methodological approach has a wide range of potential applications and provides a starting point for a research question that is currently somewhat understud-ied: namely, the trade-off between topological and dynamical complexity on higher-order domains such as hypergraphs, simplicial and cellular complexes.Future research is needed to further investigate this connection between (hypergraph) topology, dynamics and effective order.

Methods
In this section, we formalise our hypergraph dynamical system model in more detail.In particular, we consider a dynamics in terms of the topological structure and the local dynamics.The topology of the system determines which components of the system interact locally.The dynamics and specifically the update functions on the edges determine how the components interact (update), and how these local interactions are globally combined (projection) is again determined by the topology.

Hypergraph dynamical systems: Analytical
framework Let V = {1, . . ., N} be a set of nodes in the system, represented by vertices.Let x(t) ∈ R n be the vector of dynamical state variables of the nodes and t the time.The hypergraph topology is given by a set of hyperedges E = {E 1 , . . ., E M }.For simplicity we group the hyperedges according to their size d to obtain the whole hyperedge set Topology For each d-hyperedge E d α = {i, j, . . ., k} we define its indicator matrix such that S d α x maps the node state vector x to the node state vector x E d α which only contains the state variables of the nodes included in hyperedge E d α .We can then construct a (lifted) state vector by concatenating all the individual indicator vectors of the hyperedges by multiplying the original state vector x with a matrix L d ∈ R d|E d |×n defined as: The resulting vector corresponds to lifting the initial state vector into a higher-dimensional state space, hence we call the matrix a lifting matrix.The lifting matrix determines which components of the system interact, represented by its non-zero entries (see Figure 5).

9/13
Dynamics The (possibly non-linear) dynamics are then defined by the update operator F : R d|E d | → R d|E d | .This operator consists of components F : R d → R d which defines the update for each hyperedge: Even though the nodes are ordered, the update function is effectively a function on the interaction set and structured as: where we recall that the individual interaction functions f d are invariant to a permutation of their last d − 1 arguments.In other words, the function f d (x i , { {x j (t) | j ∈ E d α , j ̸ = i} }) computes the update for node i resulting from the hyperedge E d α .
Following this generally non-linear update step, the updates are projected back onto the nodes via the transposed lifting operator, which simply sums up the updates on all hyperedges that a node is a part of.This procedure is displayed in Figure 5.This results in the following global dynamics: For each node i the dynamics can thus be written as Note that, since F d is equivariant under permutation of the ordering of the nodes in each hyperedge, the dynamics is globally equivariant under permutation of the nodes.This restricts the dynamics that are possible.For example, is the only possible dynamics in this framework that can be called linear.Still, many commonly used dynamical systems can be rewritten in our framework.

Technical details of HyDy-GNN
In our analytical framework for hypergraph dynamical systems (eqs.(13) to (15)) the values of nodes of each hyperedge are collected by the lifting operator L d , then transformed by a possibly non-linear function f d and this update is projected back down to the node space by summing over each update.show the update process of a single node according to the hyperedges of order 3 in detail.Using the lifting operator L 3 , the values of nodes that are part of hyperedges of size 3 are collected.Then, a possibly non-linear function on these hyperedges is computed to obtain the update for a specific node, which is centred in these update functions.Subsequently, the obtained update values are projected back down to the node space by summing over each update.
As already stated in the main part of this work, we adapt this framework to introduce a neural network-based learning approach for dynamics on hypergraphs.We call our hypergraph dynamics learning model (HyDy-GNN ).
A technical issue in this context is that the input for the update functions in our analytical framework are sets of state variable, which is a type of input Neural Networks cannot deal with.We thus derive equivalent formulations for eqs.( 14 and the global update function in eq. ( 14) is then redefined as F ′ d by substituting in F′ p d instead of Fd accordingly.This substitution results in the same dynamics as long as f d is decomposable in p-dimensional functions as in definition 2. The model then still computes the global updates, i.e. the

Figure 1 .
Figure 1.Hypergraph dynamics model and visualisation of the topological, dynamical and effective order.When modelling complex dynamical systems, it is necessary to take into account both the system's topology and the local dynamics which are acting on the topology.We introduce a hypergraph dynamics model which separates the impact of topology and dynamics.In A, we show the update process of a single node according to the full hypergraph topology.All nodes in the hypergraph are updated synchronously.In particular, the values of the nodes that are part of hyperedges of size d are collected.Then an update function on these hyperedges which defines the dynamics is computed to obtain the update and then projected back to the node space by summing over each update.With the help of this framework, we can define and analyse different orders of the system, which we visualise in B. The topological order is given by the size of the largest hyperedge in the hypergraph.To derive the dynamical order, we look at the update functions on a hyperedge of topological order 4 in more detail.In particular, the update resulting from a hyperedge of topological order 4 can consist of a linear combination of pairwise or three-way functions which results in dynamical order 2 or 3, respectively.We can then derive the effective order of the hypergraph dynamics which is the minimum of the topological and dynamical order.The illustration highlights that the topological order of a hypergraph is only an upper bound, not a lower bound, on the effective order of the system.
Learned MCM update functions.

Figure 2 .
Figure 2. Learned update functions.Here we see the ground truth (left columns) and learned update functions (middle columns) resulting from a HyDy-GNN of order p model ∈ {2, 3, 4} (top to bottom) trained on a synthetic derivative dataset of pairwise Kuramoto dynamics (a) and pairwise MCM dynamics (b) of order p = p dyn = 2.In particular, the x-and y values correspond to a range of two-dimensional inputs of the functions and we plot the corresponding function value.In the case when p model > 2, we only plot the output of the MLP 2d with input dimension 2. The ground truth update functions are shown in the left and the absolute errors of the approximations, which are given by the absolute difference of the individual function values of the ground truth and learned dynamics, are displayed in the right columns.We observe a very good approximation to the ground truth dynamics for all model orders and both types of dynamics.This implies that our framework allows us to learn the update function of a dynamical system from data.As this is true for all model orders, we can additionally conclude that the framework can separate the effects of topology and dynamics, as it accurately captures the pairwise update functions even for higher model orders.

( c )Figure 3 .
Figure 3. Experimental results on synthetic hypergraphs.In this figure, we compare the learning results for the derivative ((a) and (b)) and the trajectory dataset ((c) and (d))and show that we are indeed able to learn the effective order of a hypergraph dynamical system from data.In (a.1), we display the mean absolute error (MAE) of cross-validation with 10 folds for all three models with p model ∈ {2, 3, 4} which were trained on the respective dynamics datasets (top row).We can see that only if the model order is larger or equal to the effective order of the system (which is p min = p for Kuramoto, SI and MCM and p min = 2 for Diffusion), the dynamics can be learned by the HyDy-GNN .In (a.2) we plot the model-corrected performance which is maximised for a model that optimises for both accuracy and model complexity.We see that in almost all cases our performance score selects the model with p model = p min .This shows that our method is indeed able to infer the effective order for the synthetic derivative training set.Only in the case of SI dynamics, we can see that dynamics with effective order p min = 4 can already be roughly captured by an order 3 model, which results in the selection of p model = 3 as the optimal order.In (c), we evaluate the performance of our model for long-term trajectory predictions based on the trajectory training set.In (c.1), we thus plot the MAE of ground truth and a predicted trajectory.This evaluation thus includes accumulated errors.We see that the results are similar to (a), whereas in this case the model-corrected performance score also suggests p min = 3 for a dynamics with p dyn = 4 for MCM and SI, even though the model of order 4 leads to a lower MAE.If a more accurate model is preferred, one could reweight the model-corrected performance score accordingly.In the right column (b) and (d), we show an example of the long-term prediction of trajectories with 200 time steps for the two datasets.In particular, we plot ground truth Kuramoto trajectories (blue) with p ∈ {2, 3, 4} (top to bottom) and the predicted trajectories given by a HyDy-GNN with p model ∈ {2, 3, 4} (left to right).We can see that only models which are equal to or larger than the effective order of the system (which is given by p min = p for Kuramoto dynamics) can capture their long-term behaviour.Moreover, the trajectory training set provides slightly better results.This is because the trajectory dataset is biased towards the long-term behaviour of the dynamics.As a result, it is more suitable for this type of long-term trajectory prediction.Generally, the results confirm that our method is capable of learning the effective order, both for the derivative and trajectory training set and the resulting models can capture the long-term behaviour of the system.

4 Figure 4 .
Figure 4. Experimental results on a real-world contact pattern hypergraph.In this figure we display the results of experiments on three real-world datasets.In particular, we simulate epidemic spread (SI), opinion dynamics (MCM) and diffusion with order p ∈ {2, 3, 4} as all of these dynamics can be expected to appear on a contact pattern hypergraph.In (1), we again show the MAE of cross-validation with 10 folds for models with p model ∈ {2, 3, 4} which are trained on derivative datasets.As in the case of the synthetic datasets, we can see that only if the model order is larger or equal to the effective order, the trajectories can be learned by the hypergraph neural networks.This is undermined by the results of the model-corrected performance in(2).We see that the model order of the best model corresponds to the effective order of the system.Again, only in the case of SI dynamics, we can see that p dyn = 4 can already be roughly captured by a model with p model = 3.Overall, the results highlight the importance of carefully considering the type of dynamics when working with real-world hypergraph structures.
and E d α is an arbitrary hyperedge with index α and size d i.e. |E d α | = d.The topological order of the hypergraph is defined by its maximal hyperedge size k.In the case of networks, we have that d = 2 as the sets represent edges and thus only consist of two nodes E 2 α = {i, j}.For d > 2 the sets E d α = {i, j, . . ., m} represent hyperedges in which more than two nodes interact locally.

Figure 5 .
Figure 5. Hypergraph dynamics model.In this figure weshow the update process of a single node according to the hyperedges of order 3 in detail.Using the lifting operator L 3 , the values of nodes that are part of hyperedges of size 3 are collected.Then, a possibly non-linear function on these hyperedges is computed to obtain the update for a specific node, which is centred in these update functions.Subsequently, the obtained update values are projected back down to the node space by summing over each update.

4/13 2.1 Learning a hypergraph dynamical system
In our formulation of hypergraph dynamical systems, the values of the nodes of each hyperedge are transformed by a possibly non-linear function f d , and this update is projected back onto the nodes by summing over each update.In the spirit of this framework, we introduce a neural network-based learning approach for dynamics on hypergraphs.We call our framework Hypergraph Dynamics Graph Neural Network (HyDy-GNN ).Apart from minor technical adaptations, which are described in more detail in section 5, the main difference between HyDy-GNN and the analytical formulation of our dynamics is that in HyDy-GNN , we approximate the update function f d using multi-layer perceptrons (MLPs).We thus do not need to specify the functional form of the dynamics a priori, but can learn the dynamics entirely from observational data, as MLPs are universal function approximators that can approximate any function to arbitrary accuracy.In order to estimate the effective (minimal) order of the dynamics, we consider HyDy-GNN s operating on MLPs with different input dimensions.For a HyDy-GNN of order p model only MLPs with up to p model input variables are used to model the p-dimensional functions φ p d in definition 2. For a general multi-set s of size d − 1 we then have y 1 , { {y 1 , ..., y d } }\y 1 ) f d (y 2 , { {y 1 , ..., y d } }\y 2 ) . . .f d (y d , { {y 1 , ..., y d } }\y d )