Integrating GNN and Neural ODEs for Estimating Non-Reciprocal Two-Body Interactions in Mixed-Species Collective Motion
Authors: Masahito Uwamichi
Simon K. Schnyder
Tetsuya J. Kobayashi
Satoshi Sawai
Venue: NeurIPS 2024 poster (2024)
Keywords: Deep Learning
Neural Differential Equations
Graph Neural Networks
System Identification
Active Matter
Collective Motion
Non-Reciprocal
Primary area: machine_learning_for_physical_sciences
TL;DR: We present a framework to estimate non-reciprocal two-body interactions from trajectories in mixed-species collective motion. Our method accurately replicates interactions and collective behaviors, demonstrated through numerical experiments.

Abstract: Analyzing the motion of multiple biological agents, be it cells or individual animals, is pivotal for the understanding of complex collective behaviors. 
With the advent of advanced microscopy, detailed images of complex tissue formations involving multiple cell types have become more accessible in recent years. However, deciphering the underlying rules that govern cell movements is far from trivial. Here, we present a novel deep learning framework for estimating the underlying equations of motion from observed trajectories, a pivotal step in decoding such complex dynamics. 
Our framework integrates graph neural networks with neural differential equations, enabling effective prediction of two-body interactions based on the states of the interacting entities. We demonstrate the efficacy of our approach through two numerical experiments. First, we used simulated data from a toy model to tune the hyperparameters. Based on the obtained hyperparameters, we then applied this approach to a more complex model with non-reciprocal forces that mimic the collective dynamics of the cells of slime molds. Our results show that the proposed method can accurately estimate the functional forms of two-body interactions -- even when they are nonreciprocal -- thereby precisely replicating both individual and collective behaviors within these systems.

Integrating GNN and Neural ODEs for Estimating Non-Reciprocal Two-Body Interactions in Mixed-Species Collective Motion
Masahito Uwamichi' Simon K. Schnyder", Tetsuya J. Kobayashi Satoshi Sawai
*Graduate School of Arts and Sciences, The University of Tokyo
IInstitute of Industrial Science, The University of Tokyo
Universal Biology Institute, The University of Tokyo twamichi@g ecc. w-c... ac. jp. simon@sat. t. u- jp. tetsuyadsat. .t. u-tokyo ac jp, casawai aaaaeiimaallee..ee ecc. u-tokyo.... ac. ac.
Abstract
Analyzing the motion of multiple biological agents, be it cells or individual animals, is pivotal for the understanding of complex collective behaviors. With the advent advanced microscopy, detailed images of complex tissue formations involving mul- tiple cell types have become more accessible in recent years. However, deciphering the underlying rules that govern cell movements is far from trivial. Here, we present a novel deep learning framework for estimating the underlying equations of motion from observed trajectories, pivotal step in decoding such complex dynamics. Our framework integrates graph neural networks with neural differential equations, enabling effective prediction of two-body interactions based the states of the interacting entities. We demonstrate the efficacy of our approach through two numerical experiments. First, we used simulated data from toy model to tune the hyperparameters. Based on the obtained hyperparameters, we then applied this approach a more complex model with non-reciprocal forces that mimic the collective dynamics of the cells of slime molds. Our results show that the proposed method can accurately estimate the functional forms of two-body interactions even when they are nonreciprocal thereby precisely replicating both individual and collective behaviors within these systems.
Introduction
Collective motion, phenomenon observed in various biological systems, is characterized by coordinated movement of multiple entities. This behavior is prevalent in wide range of self- propelled systems, collectively referred to as active matter ( from flocks of birds 2
3 4 and schools of fish (5) to cellular slime molds (6), microswimmers (7), swarming leukocytes (8), and human crowds (.) Understanding the underlying mechanisms collective motion is crucial for elucidating the principles governing the dynamics of these systems. In particular, the interactions between individual entities play a pivotal role in shaping the collective behavior of the system.
Recent advances in imaging technologies have enabled detailed observations at the cellular level, providing insights into the dynamics of complex tissue formation involving multiple cell types. For example, cellular slime molds, model organism for studying collective motion, exhibit intricate behaviors such as aggregation, migration, and differentiation (6). These behaviors are driven by the interactions between different cell types, which are mediated by chemical signals and physical forces. Decoding the underlying equations of motion that govern these interactions essential for understanding the emergent properties of these systems.
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

In this work, we present a novel deep leamning framework for estimating two-body interactions in a mixed-species collective motion. Our framework integrates graph neural networks (GNNs) with neural differential equations (neural DEs) to predict interactions between pairs of entities based on their states. GNNs are well-suited for modeling complex interactions graph-structured data, while neural DEs provide a flexible framework for learning the dynamics of the system. By combining these two approaches, we can effectively capture the interactions between individual entities and predict their collective behavior. We demonstrate the efficacy of our framework through two numerical experiments. The first experiment uses a toy model designed to generate data for refining the hyperparameters of our framework. The second experiment explores complex scenario that partially mimics the collective motion of cellular slime molds, where two different cell types interact with each other. Our results show that our method can accurately estimate the two-body interactions, thereby replicating both individual and collective behaviors within these systems.
The rest of this paper is organized as follows. In section 2 we introduce the study of collective motion. In section [] we provide an overview of related work on collective motion and deep learning for dynamical systems. In section[4 we describe our deep learning framework for estimating two-body interactions in mixed-species collective motion. In section 5
we successfully apply this framework to two numerical experiments. Finally, in section 6 we discuss the implications of our work and outline potential future research directions.
2 Background
Let us first describe the formulation of collective motion in active matter. Starting with the Vicsek model (I0), the collective behavior in active matter is described based on the centroids, velocities, and orientations of each individual. A unique aspect of active matter is that it allows for spontaneous generation of forces and torques, which is justified by the ability to extract and utilize energy from the external environment (L.. The Viesek model itself is a multi-particle model in discrete time, where each individual possesses a velocity along its orientation and adjusts its direction based on interactions with nearby individuals. This model assumes that interactions among individuals are local and can be represented using dynamical graph structure (12). Stemming from this model, numerous other models have been proposed, differing in the nature of interactions between individuals and the forms of their equations of motion. Particularly models with continuous time motion equations often use Langevin equations where the motion of individuals is described by the summation of pairwise interactions (13). Furthermore, since active matter can utilize energy from the external environment, focusing solely on moving agents categorizes as an open system. This implies that the conservation of the total system energy is not guaranteed, making these systems inherently unsuitable for Hamiltonian descriptions. Given this background, the formulation of collective motion in active matter typically involves direct descriptions of forces rather than using free energy.
3 Related Work
System Identification methods: This section discusses methods for estimating governing laws from data, commonly referred to as system identification. System identification aims to estimate the equations of motion of a system from data. One well-known method for estimating differential equations is Sparse Identification of Nonlinear Dynamics (SINDy) (14). SINDy estimates nonlinear differential equations from data using LASSO regression to sparsely estimate the terms of the equations. However, when applying this method many-body systems with pairwise interactions, the number of parameters to estimate increases exponentially with the number of individuals, which raises computational costs and leads to potential instabilities. Other system identification methods use Bayesian optimization (15) 16 17), c.g. representing the functions in the motion equations through basis function expansion, and estimating their coefficients via Bayesian optimization. This method has been successfully used in Vicsek model, reproducing its formation of orientational order. However, the effectiveness of this approach for more complex systems, such as mixed- species systems, remains unclear. Deep learning approaches to analyzing collective motion have also been proposed (18) 19 20). These incorporate Attention mechanisms to analyze behaviors limited to turning right or left (18), systems that limit the scope to Vicsek-type models to estimate orientational order parameters (19), or to estimate two-body interactions by expressing forces as neural networks (20). Another possibility is to express Hamiltonians (21 or Lagrangians (22) as
2

---

neural networks, in order to preserve known symmetries of the studied system. Some proposed including knowledge or assumptions about the physical context into the loss function (23). So far, none of these methods have been applied to systems in which multiple species interact with each other in unknown ways.
Integration of graph neural network with neural differential equations: In this paper, we integrate GNN with neural ordinary differential equations (neural ODEs) to estimate interactions through evaluating equations of motion with dynamic edge structure. Most of the previous studies on integration of the two methods solve and train ODEs on graph with fixed edges and consider This approach requires fully or almost fully connected graph for collective motion where the adjacency changes drastically over time, and causes memory inefficiency. On the other hand, some introduced given time series of edge structure (29), making the system difficult to be extrapolated. Thus we resolved these issues by setting a rule to define edges at each time point instead of explicitly giving graph structure time series.
changes in the edge weight (24 25 26 27 28).
4 Method
We aim to estimate the rules of motion for individual entities within collective motion data. This section presents the general framework that is used to generate the training data from physical model for two numerical experiments, as well as the learning algorithm and neural network architecture.
4.1 Physical model
We represent the state of each entity at time t as Z denotes the variables described by the motion equations in the state space X' C R", and y', a non-temporal variable, represents auxiliary attributes such the type of each entity within the feature space V We define a distance function d time if (t) d()) 2 (t)) < dg, where do R is predefined threshold. The motion of each entity governed by the following Langevin equation:
2 3 = (z)))) E X x V, where x' X
Z + R, and assume that entities and j interact at
dx' (2)()
Σ (=10,.)))) + (t).
(1)
jst dd)())/d <donj
with given self-driven forces F(!)
R", forces due to interactions between pairs F
2 x Z + R", intensity of noise R, and Wiener process (' R". Given the initial state (to) E Z at time to E R. the state Att E after time interval At E IR can be determined by solving the motion equation using the operator S(F(1), F(2), a. do) (Z(t)), At) + Z(to + At). Here, (to + At) (to + At), v²) defined by the following integral:
tot At
x' (to + At) = x' 2'()) +
do' dx'
(2)
which is computed numerically (section 4.3 For o we evaluated eq with the neural ODE method (30 3T), and for a * 0 the neural stochastic differential equation (neural SDE)
method 32 33).
4.2
Interpretation of the system as graph
The system at any time / can be represented as a directed graph G(I) = (V, E(), Z()), where the set of entities represent the vertices V = 2,, ... N), pairs of interacting entities constitute edges F(()= = j(, E V x Vide(() < do A j]. and Z())= = (='()) E V} the set of states of all entities. Importantly, since the entities change states over time, the graph in general evolves over time. The forces ((I) and F(2) are evaluated on each vertex and each edge, respectively.
4.3 Solving the equations of motion of the graph through neural DEs
To obtain a time series of collective entities (eq. 2). we have to evaluate the equation of motion (eq. We achieved this with made wrapper that connects neural ODE (SDE) to GNN (Figure
3

---

Initial State (z())
Neural ODE
State at time t+ (z(++dtt)
State at time (2(())
GNN
Define edges (E(t)
Caleulate self-generated force (F (z(())
Calculate interactions (FA (())
Total Force
Total Interaction CΣ fs uor P (z((0,(()))
Figure 1: Schematic figure of the integration of GNN and neuralODE. For details, see section 4.3
Dp. The wrapper is equipped with the distance function d and the threshold do to define the edges Ett each time t, the functions to calculate self- propppooplliion F() and interaction F=22, and the noise intensity The computation S(F(), p(2), d, do) is then run through the following steps:
0. Create a graph object only with nodes V and set the static states {s') } of the entities. Store this graph in the wrapper object. Pass the initial state (x)) into the neural ODE (SDE) object (Figure T bluu).
1. The neural ODE (SDE) object passes the current state (x) to the forward function of the wrapper object.
2. The forward function places the state values (z) at the nodes of the graph. Subsequently, the graph is passed GNN (Figure rreg).
3. GNN first defines and updates the edges Ett) using the distance function d and the threshold do stored in the wrapper. Next, it calculates the interactions F(2) the edges and adds them to the self-pro propupropulsion F(1) on the nodes to return the forces to the ODE (SDE) object.
If SDE, GNN also returns the noise intensity in the wrapper.
4. The neural ODE (SDE) object updates the state (ar(t + dt
5. Repeat steps 4 until the end of the time series.
When inferring the forces FO), F(2) we model them using neural networks (section 4.5)
4.4 Generation of training data
To construct S.m =((F()) ((2)) Osi, dhim, d. distance functions dain provided, suitable choices for the noise magnitude asim, threshold do.sim and forms for the forces f:(1) f=(2) were specified (see section 5). Thus defined S was used to conduct simulations according to the motion equations, from initial states assigned uniformly at random z0) E Xo with appropriate values for y¹ y. The simulation was repeated M times, and the m-th result denoted as Zm :t + Sxim ( (0), f) with the m-th initial condition Z (0). The combined data set Zuim
{7mh im 11...,MM was then used to apply the subsequent learning algorithm, aiming to estimate [1)
sim and (2)
4.5 Learning algorithm and neural network architecture for the forces
Subsequently, developed a learning algorithm to estimate the simulations conducted.
(1) sim and p(()) using the results from
The self propulsoolion and interaction functions, F() and NN were modeled using neural networks.
(2)
Unless otherwise specified, these functions consisted of a three-layer fully connected network
followed by scaling layer. Each fully connected layer consisted of 128 nodes, with the Exponential Linear Unit (ELU) serving as the activation function. The scaling layer involved a scalar R and vector B IR"" as learnable parameters, and it transformed an input vector E E to e a B. This

---

Training data; (A)
1-10
1-60
1-110
(B)
Estimated Forces and Dynamics:
(C)
1=60
11.0
N.
V=1x10"
A00
(8)
1-10
1-60
11.0
-5
X10
(K)
X30
400
X10
Figure 2: Snapshots of the simulation results in the harmonic interaction model with friction constants.
(A) The case with a friction constant o = x 10-2 and (B) with 7 = x 10-1, Positions of individuals are indicated by blue circles, and velocities by black arrows. (C-H) The functions estimated from the data with harmonic interaction plotted against the true values. The upper panel displays F(1) and the lower panel (22) for different cases: (C-D) with ? x 10-2 and Ni.. = 270, 1 x 10-1 and N = 3. Blue and
(E-F) with y 1 x 10-2 and N 3, and (G-H) with orange indicates x- y- element F. A black line serves as guide indicating where the estimated values equal the true values. (I-K) Snapshots of the simulation results in the model estimated from the data with harmonic interaction. The panels represent different scenarios: (I) with 7 = 1 x 10-2 and Nu. 270, (J) with x 10 and Nu, 3, and (K) with 7 = 1 10- and Nu.. 3.
configuration allowed the output of the fully connected network to be scaled appropriately, reflecting physical scales.
The parameters of fully connected networks and the scaling layer, collectively denoted were optimized to approximate 1 (0) and 1x (0) to f:(1) and fa(2) respectively. To evaluate the sim
deviation of F (0) and F(2 (f), we first solved the motion equations using these functions. In practice, we constructed wrapper object SNx(O) S(F()) (0), r(2 (0), 0, daim do, and used perform simulations. Here, dp,NN was a suitably set threshold, and for simplicity, ONN = In the simulation, the initial state was snapshot in training data zm (to), and we used this
5

---

to compute the state T time units later as (0: to + +) = SNN (Z) (fo), r), which would be compared with a truth 2im (tg + T).
was used to evaluate the discrepancy between the solution by 0 and the simulation data through L (0: (o) L _( (to + r), Z (0; to 7)). To optimize 0, the average [()) of this loss function over various m and to was minimized:
For each numerical experiment, a loss function L: 2 x
arg min (()
The M simulation data sets were split into two parts, with Mr sets used as training data and Myal sets used as validation data. This optimization was conducted using the LAMB optimization algorithm (34), leveraging gradient information with respect to 8. The parameter 8 was evaluated using the minimum of the loss function calculated on the validation data set.
5
Experiments
Here, we demonstrate the efficacy of the proposed method two numerical experiments.
5.1 Underdamped Brownian Motion with Harmonic Interaction
Training data: First, to demonstrate that the self-propulsion and interaction forces could be estimated using our method, we performed simulations using a toy model. We used N 100 particles. each with position r' (t) E IR² and velocity (f) E R?, such that 2 (t) = (t) = (r' (t), (t)) and y' is empty. The particles were connected by harmonic oscillators. and each velocity was subjected to damping force, so the equations of motion read:
dr v'dt,
du'
(-70²-
Σ JEL
V di + odw"(()
(5)
U(r) = 2 -r)):
(6)
where & the strength of the interaction, re the natural spring length, and y R the friction coefficient. The distance function was dsim(2 22) = 1r' - r²). with threshold set at do.sim 5.
In order to fit these terms into the aforementioned framework, the equations were reformulated as:
F!1) (2')
(F((1) (x*) (11)) ( (2*)) = (v', -r0"),
F(2) sim
(2 2d) = (Fr(2) sim, (2 2²) (22)
sim (2 2²)) = (O. -))).
(8)
Here, h = 1 10 = =, 1 x 10 and was set to 1 10 2 or 1 10-1 The initial states x² (0) were uniformly sampled from Xo = [0, 5)2 x (0, x 10 31,, and simulations were performed to obtain Zzim-
For these simulations, the Buler-Maruyama method with time step of 0.1 was used, and data were collected at = 0, 1..... 50. The data showed that all particles cyclically gathered towards the center before dispersing (Figure 2 A-B)). When & = 1 10 2, the amplitude of uggregation-dispersee
was
nearly constant (Figure 2 A), Supplemental Movie S1), whereas for
x 10 the amplitude of
segregation-disperiionn decreased (Figure 2 B), Supplemental Movie S2).
Learning: With these data, we applied our proposed method to estimate F(2) were modeled follows:
and F1) and sim
F(2) sim
F
UL (=':0) = (v', (u', NN (2 20:0 = (0. (1- - - 0))
(9)
(10)
. and P were each independent neural networks as described in section 4.5 We evaluated the predictive accuracy for position and velocity follows, using Z (0; to +T) predicted for 7 10

---

and the simulation data Z (to):
(2 (to + T), LNN 10;00 to T))
= N
Σ
IA (0: +7) - sim (to + r)1
(11)
IEV
= 1
(Z (to + +), N to 7))
(0NN(0:00+7)) EVV
r'si (to + r)1?.
(12)
We sampled these metrics for 60 randomly selected pairs (m, In)- The total loss function L()) for one batch was defined as the dimensionless sum of these metrics, normalized by the variance in the simulation data:
(,)
Σ (m,ta)
(Z (to T), Z (0; to + T)) Var,,, tm. Ir (to (to + T]],
(13)
(o)
E
(Z (to +T), ZNN (0; F))/ M (to + T)],
(14)
(m,ta)
((0) = (.) + (0)
(15)
To minimize the loss functtion, we manually searched for optimal hyperparameters for the LAMB
optimizer. We determined that 31 = 0.5, B2 = 0.4, = x 10-6, and a learning rate of 1 x 10-3 effectively minimized the loss function, although convergence was notably slow. important to note that we were unable to identify hyperparameters that would allow for both faster convergence and adequate estimation of the functional forms of F(1) and F(2), As a result of minimizing this loss function, we optimized 0 to estimate F and F (Figure (C-H)). When
10-2 and the (() training set size was Mua 270 and validation set size Mvai 30, after 300 epochs (15 days), NNN, and
(2( NNoo were observed to approximate simme and p respectively (Figure 2 C-D)). In contrast,
F[1)
with Mun = 3 and Mval = 3, especially F(1) did not approximate F(1) sim, even after 6000 epochs (4 days; Figure F)).). However, under the same dataset size but with y x 10-1, (()) and (=)2) NN,,v
were confirmed to approximate sim, and p(2) sim, (Figure 2 G-H).
The accuracy of the estimation results is quantified in Supplemental Table S For each particle or pair in the dataset Zsi, we calculated the forces F. and F!??) and computed the Mean Squared Error (MSE) and Mean Absolute Error (MAE). These errors were then normalized by the I.1 or L2 norm of F(1) sim.e and F(2) sim,e respectively, to provide dimensionless measures of accuracy. Due to the extensive computation time required for estimation, comprehensive statistics could not be gathered. Instead, we present the results of all trials to illustrate the trends in estimation accuracy. As observed, when Mya = 3 with a low friction constant = x 10-2, significaat estimation errors
occurred. Conversely, the estimation errors for the interaction force F(?) were largely unaffected by y, suggesting that the estimation of interaction and friction are somewhat independent. However, when Mra 270, slight improvement in the accuracy of interaction estimates was noted compared to the Mf = 3 scenario. This observation indicates that the accuracy of estimation depends on the number of data points, for the interaction well as the friction.
Furthermore, to verify how well these estimates fit the training data, we conducted simulations from random initial values Z(0) sampled in a similar manner to training data creation, and visualized (Figure IIKK); Supplemental Movies S3-5). In all cases, the estimates
the results (Z confirmed to adequately replicate the training data for all cases, even with y x 10-2 and Mμa = The ability to replicate training data suggests that friction had almost no effect in this case. However, this implies that the proposed method may not adequately estimate very weak effects.
5.2 Mixed Species Collective Motion with Overdamped Self- eeff-propulsion
Next, to test the proposed method in complex systems involving interactions among multiple species, we conducted simulations using more complex model that emulates real collective movements.
Training data: this model, each of the N 400 individuals has position r(() E [0, I.J., orientation 5'(t) E [0, 2T), and species type cf e {0, 1}. Both r' and o² are subject to periodic
7

---

Training data: (A)
20000
(B)
Pammeter
Pammrntter set
Estimated Forces and Dynamics:
(K)
2001
70-
Paraneter set
5
6
15
X10
(H))
(I)
(L)
70
2000
Parameter set (l)
-2
15
100
Figure 3: Snapshots of the simulation results in the mixed-species model. Panel (A) displays a representative snapshot of the training data for the case with parameter (i), and panel (B) with parameter set (ii) (see Supplemental Table s2). Each individual is represented by an arrow located their position and directed toward their polarity. Black arrows indicate individuals with cf
while red arrows represent those with e = 1. (C-J) The functions estimated from the data with the mixed-species model plotted against the parameter set (i). and (G-J) with parameter set (ii). A gray line serves as guide indicating where the estimated values equal the true values. (K-L) Snapshots of the simulation results in the model estimated from data with the mixed sixee-species model. (K) The case with parameter set (i), and (L) with parameter set (ii).
values. The rows display different cases: (C-F) with
boundary conditions, ensuring continuous and consistent movement dynamics across the defined space. This defines x² (1) = (t" o' (()) and y² = (c*), constituting =²()) (x* v1), The motion equations for individual are described as follows:
dr = =Vopp'
Σ
BJ
(16)
ELL)E))))
dø
Σ
(acF(c') + oculd')()) (r" .pi) dt + od)
(17)
E(
JN
(r. - 1 pis
(18)
1 r pd
5(1 |rid
(19)
pàù p²
Jy
(x)
(20)
|rés|
p (cos Φ², sin ) pi (- sin Φi cos Φ¹),
(21)
3 = r - E |- 1-L22,
(22)
with self ceff-propulsion speed vn 2 0, the strength of the excluded volume interaction / 2 0, the strengths of contact following and chemotaxis acr, ach 0, respectively, the diffusion length of the chemoattractant > 0. the noise strength 2 0, and the modified Bessel function of the second kind Kj. The model was obtained by introducing chemotaxis term into a preexisting model (35). to make it more appropriate for cellular slime molds. The chemotaxis term assumes rapid diffusion
8

---

of chemotactic substances secreted by each individual (36). We varied the sensitivity of contact following and chemotaxis depending on the species type, thus modeling a system where different species interact non-reciprocally. The interaction terms were adapted into our framework as follows:
(1)1) sim
(x*)
(F F (x')
(23)
(vop', 0),
(24)
F sim
(2' x)) = (m (2', (=',,") F o (2')=)))
(25)
(BJEV (ocr(c")
aci(d") (ea)))))) (r" .11)
(26)
The distance function was dsmm(2 20) = |r with periodic boundary conditions of width L = 20
and a threshold do,sim = 1. The parameters were vn = L, 3 = L, A = 0.5, and a = v0.2, with 200 individuals each species cf 0, 1. Parameter dependencies based on cf are shown in Supplemental Table $2
The initial states x' (0) were uniformly sampled from [0, L] x [0, 2- and simulations were carried out to generate Zim using the Heun method with timestep of 0.1, collecting data for t = 0,1,... 300. Simulation results for parameter set (i) showed the formation of mixed-species circular clusters ex- hibiting rotational behavior (Figure ([A) Supplemental Movie S6). Conversely, results for parameter set (ii) formed centipede-like clusters of mixed species that moved translationally local alignments (Figure B), Supplemental Movie S7).
In order to estimate p(()) and Ell)) we modeled ENN and E as follows:
Learning:
(2 -R(00) F FN,, (c: (c: 0), Fl (ct
(27)
(2) NN (=", 22;0) ( (R(O)) (2)
F(2) NN. (R (R(- . - 2 0,C,0,0)
(28)
F12)
( (R(-)R)())))))) pr'r Φi - ,.c.,,,, e,
F N are independent neural networks with their respective parameters, and R() is the rotation matrix for ø', which converts all angles into relative angles with respect to o¹, thereby maintaining the model's invariance to coordinate rotation. Additionally, when inputting the categorical variable c' into the neural network, is transformed into two-dimensional vector via an embedding layer. The predictions Z (0: to + 7) made by this neural network for 1 and the simulation data z (to) defined the prediction errors for position and orientation as follows:
Here, ,N NN. NN,
F( (1)
C, (to + T), ZNN (0; to +T)
1
Σ
d(rá + ram + 7))
(29)
N
Lo (to + T) 7m (0; +
N (1
1
Σ NeV
cos (0; +7) D'im (to r))) (30)
Predictions were computed using the Buler-Maruyama method with timestep of 0.1. These metrics were sampled for 60 randomly selected pairs (m, 4)) and each was normalized by the variance in the simulation data to compose the loss function for one batch:
(0) = Σ
L. 7 (to + +7, NN (0; 40)) /arr,, ( ("o + +r)].
(31)
(m.
(0) = Σ
(Z (1o + T), Z (0; (0))
Σ (1-cos((lae cos( (40 + +)))
(32)
(m.)oo
((,,to)
((0) (0) + (0).
(33)
To minimizing this loss function, we tried the same hyperparameters as those used for the harmonic in- teraction model, and found we optimized 0 to estimate F and (22)) (Figure C-J). All experiments
were conducted with Mm = 3, Myal = 3. After 20,000 epochs (25 days), I'N,, P P ( were found to approximate F(())_. r. di sim.r of our estimation results with the training data, sampled random initial values ((0) similarly to the
(1) (2)
sim, closely (Figure 3 C-J)). To further verify the fit

---

training data creation process and conducted simulations to visualize the outcomes (NN(0)) t) (Figure B-L; Supplemental Movies S8-9). The results confirmed that all estimates adequately reproduced the training data.
To further test the predictive power of the model trained on the mixed state, performed simulations where all individuals were of the same type. For initial conditions Zo( where all _N = 400 individuals shared the same type C¹ c, either c = 0 or c 1, we performed simulations using the aforementioned model Ssim and the neural network SNN ( (0*) trained above (Supplemental Figure $D For parameter set (i) with c O. both Sim (0), t) and SNN(#" (Z) t) showed the formation of two clusters that exhibited rotational movements (Supplemental Figure AT])-)). Conversely, for parameter set (i) with c 1 and parameter set (ii) with c 0, 1. both (zi (Zo()), t) and SNN (0), t) formed centipede-like clusters that moved translationally in alignment with their local direction (Supplemental Figure C-]C These results suggest that the interaction of each type was approximated accurately enough to explain the dynamics of training data.
As a negative control, we conducted simulations where both types of individuals shared the same parameters, as in parameter sets (iii)-(v), and trained the neural networks for each case (Supplemental
(1) p(2) Figure $2.. In all instances, it was confirmed that (1) NN, NN.o: NN. NN. closely approximated
p(2)
sim.es F F(22 ml()
sim.o- Notably, these estimation results showed minimal dependency on ef in
(2 (*), (2 (2', 0*)) This outcome suggests that the species dependency of interactions estimated by the proposed method in mixed assemblies is based on actual data, rather than being an artifact of the methodology.
Finally, to summarize the results above, we quantified the estimation accuracy for each trial, similar to the methods described in section 5.1 (Supplemental Table $3 To note, since F 0, the MSE and MAE for F(] not normalized. Our results demonstrate that Fy!) is estimated with high accuracy across all conditions. Also, the estimation of F(1) is maintained sufficiently small magnitudes relative the order of magnitude (10") typically seen with 1x(2) (Figure FJJ).. Regarding the interactions, both F(2) and F exhibited somewhat higher errors, approximately of the order of 10-1 These errors likely stem from the high dimensionality of the inputs to the F(2 which may prevent the neural network from achieving sufficient learning, or from the training process neglecting the noise.
Conclusion
In this study, we proposed novel method for estimating interactions among individuals within models of collective motion. This method combines dynamic GNNs with neural ODEs to estimate interactions among individuals. We demonstrated that this method was able to estimate interactions in both simple models and complex mixed-species collective motion models. In models with simple harmonic interactions, the proposed method identified the relevant interactions effectively. Importantly, the method successfully inferred the non-reciprocal interactions between the different species in the complex mixed-specier model for the collective dynamics of slime molds. For the latter case, the continuous updating of the edge structures in our approach substantially reduces memory needs. In our simulations of 400 bodies with an edge density of 2%. the required memory of 200 GB for a fully connected graph is reduced to 30GB thus making them feasible on off-the-shelf GPUs. Nevertheless, we should note that this approach is limited due to the long time required for estimation. Moreover, the method currently estimates deterministic motion equations and only considers pairwise interactions, thereby not accommodating interactions among three or more bodies or the effects of noise. Future research is expected to extend this approach to estimate more general interactions and develop methods for stochastic motion equations. Applying the present methods to real data in systems such as immune cells should help clarify the complex rules behind their migration stategies.
Acknowledgements
This research was supported by ISPS KAKENHI Grants JP22H05673 (awarded to MU), JP24H01464 (MU), JP22H04841 (SKS), JP22K14012 (SKS), JP JP19H05799 (TJK), JP19 (SS) and
10

---

JP19KKO282 (SS), as well as JST CREST Grants JPMJCR2011 (TJK) and JPMICRI923 (SS). We acknowledge the use of ChatGPT and GitHub Copilot for their assistance in swiftly articulating ideas in clear English. Passages written with the help of these tools were reviewed and rewritten by the authors.
References
[1] T. Viesek and A. Zafeiris, "Collective motion," Physics Reports, vol. 517, no. 3-4, pp. 71-140, 2012.
[2] A. Cavagna, A. Cimarelli, I. Giardina, G. Parisi, R. Santagati, F. Stefanini, and M. Viale, "Scale-free correlations in starling flocks," Proceedings of the National Academy of Sciences of the United States of America, vol. 107, pp. 11865-11870, 6 2010.
[3] A. Flack, M. Nagy, W. Fiedler, I. D. Couzin, and M. Wikelski, "From local collective behavior to global migratory patterns in white storks," Science, vol. 360, no. 6391, pp. 911-914 2018.
[4] H. J. Charlesworth and M. S. Turner, "Intrinsically motivated collective motion," Proceedings of the National Academy of Sciences, vol. 1 116, pp. 15362-15367, 2019.
[5] Y. Katz, K. Tunstrom, C. C. Ioannou, C. Huepe, I. D. Couzin, S. A. Levin, I. D. C. Designed, and I. D. C. Performed, "Inferring the structure and dynamics of interactions in schooling fish." Proceedings of the National Academy of Sciences, vol. 108, pp. 18720-18725, 2011.
[6] T. Fujimori, A. Nakajima, N. Shimada, and S. Sawai, "Tissue self-organization based on collective cell migration by contact activation locomotion and chemotaxis," Proceedings of the National Academy of Sciences of the United States of America, vol. 116, pp. 4291-4296, Mar. 2019.
[7] J. Elgeti, R. G. Winkler, and G. Gompper, "Physics of microswimmers single particle motion and collective behavior: A review," Reports on Progress in Physics, vol. 78, 2015.
[8] E. Strickland, D. Pan, C. Godfrey, J. S. Kim, A. Hopke, W. Ji, M. Degrange, B. Villavicencio, M. K. Mansour, C. S. Zerbe, D. Irimia, A. Amir, and O. D. Weiner, "Self--xtinguishing relay waves enable homeostatic control of human neutrophil swarming," Developmental Cell, July 2024.
[9] J. Kok, M. K. Lim, and C. S. Chan, "Crowd behavior analysis: A review where physics meets biology," Neurocomputing, vol. 177, pp. 342-362, 2016.
[10] T. Vicsek, A. Czirók, E. Ben-Jacob, I. Cohen, and O. Shochet, "Novel Type of Phase Transition in System of Self-Driven Particles," Physical Review Letters, vol. 75, pp. 1226-1229, Aug. 1995.
[11] G. I. Menon, Active Matter, pp. 193-218. New York, NY: Springer New York, 2010.
[12] S. Boccaletti, V. Latora, Y. Moreno, M. Chavez, and D. Hwang. "Complex networks: Structure and dynamics," Physics Reports, vol. 424, pp. 175-308. Feb. 2006.
[13] H. Chaté, F. Ginelli, G. Grégoire, and F. Raynaud, "Collective motion of self-propelled particles interacting without cohesion," Physical Review E, vol. 77, no. 4, p. 046113, 2008.
[14] K. Champion, B. Lusch, J. N. Kutz, and S. L. Brunton, "Data-driven discovery of coordinates and governing equations," Proceedings of the National Academy of Sciences of the United States of America, vol. 116, pp. 22445-2212 22451, Nov. 2019.
[15] D. B. Brückner, A. Fink, C. Schreiber, P. J. F. Röttgermann. J. O. Rädler, and C. P. Broedersz, "Stochastic nonlinear dynamics of confined cell migration in two-state systems," Nature Physics, vol. 15, pp. 595-601, June 2019.
[16] D. B. Brückner, P. Ronceray, and C. Broedersz, "Inferring the Dynamics of Underdamped Stochastic Systems," Physical Review Letters, vol. 125, p. 058103, July 2020.

---

[17] D. B. Brückner, N. Arlt, A. Fink, P. Ronceray, J. O. Rädler, and C. P. Broedersz, "Learning the dynamies of cell-cell interactions confined cell migration." Proceedings of the National Academy of Sciences of the United States of America, vol. 118, p. e2016602118 18, Feb. 2021.
[18) J. H. Heras, F. Nomero-Ferrero, R. C. Hinz, and G. G. de Polavieja, "Deep attention networks reveal the rules of collective motion in zebrafish," PLoS Computational Biology, vol. 15,
p. el007354, Sept. 2019.
[19] M. Ruiz-Garcia, C. M. Barriuso G., L. C. Alexander, D. G. A. L. Aarts, L. M. Ghiringhelli, and C. Valeriani, "Discovering dynamic laws from observations: The case of self--ropede propelled, interacting colloids," Physical Review E, vol. 109, p. 064611, 1, June 2024.
[20] H. Koyama, H. Okumura, A. M. Ito, K. Nakamura, T. Otani, K. Kato, and T. Fujimori, "Effective mechanical potential of cell-cell interaction explains three-dimensional morphologies during early embryogenesis." PLoS Computational Biology, vol. 19, p. e1011306, Aug. 2023.
[21] S. Greydanus, M. Dzamba, and J. Yosinski, "Hamiltonian Neural Networks." Advances in Neural Information Processing Systems, vol. 32, pp. 1-16, June 2019.
[22] M. Cranmer, S. Greydanus, S. Hoyer, P. Battaglia, D. Spergel, and S. Ho, "Lagrangian neural networks," ICLR 2020 Workshop on Integration of Deep Neural Models and Differential Equations, Apr. 2019.
[23] K. Linka, A. Schäfer, X. Meng, Z. Zou, G. E. Karniadakis, and E. Kuhl, "Bayesian Physics Informed Neural Networks for real-world nonlinear dynamical systems," Computer Methods in Applied Mechanics and Engineering, vol. 402, p. 115346, Dec. 2022.
[24] L.-P. Xhonneux, M. Qu, and J. Tang. "Continuous graph neural networks," Proceedings of the 37th International Conference on Machine Learning, vol. 119, pp. 10432-10441, July 2020.
[25] M. Poli, S. Massaroli, J. Park, A. Yamashita, H. Asama, and J. Park, "Graph Neural Ordinary Differential Equations," arXiv: (EIi:0991..7532 June 2021.
[26] L. Chen, K. Wu, J. Lou, and J. Liu, "Signed Graph Neural Ordinary Differential Equation for Modeling Continuous-Time Dynamics," AAAI, vol. 38, pp. 8292 8301, Mar. 2024.
[27] A. Han, D. Shi, L. Lin, and J. Gao, "From Continuous Dynamics to Graph Neural Networks: Neural Diffusion and Beyond," arXiv:2310.1012.. Oct 2023.
[28] Y. Zhang, S. Gao, J. Pei, and H. Huang, "Improving Social Network Embedding via New Second-Order Continuous Graph Neural Networks," Proceedings of the 28th ACM SIGKDD Conference on Knowledge Discovery and Data Mining, pp. 2515-2522 Aug. 2022.
[29] L. Pan, J. Lu, and X. Tang, "Spatial-temporal graph neural ODE networks for skeleton- action recognition." Scientific Reports, vol. 14, p. 7629, Apr. 2024.
[30] R. T. Chen, Y. Rubanova, J. Bettencourt, and D. K. Duvenaud, "Neural ordinary differential equations." Advances in neural information processing systems, vol. 31, Dec. 2018.
[31] M. Poli, S. Massaroli, A. Yamashita, H. Asama, J. Park, and S. Ermon, "Torchdyn: implicit models and neural numerical methods in pytorch," " Neural Information Processing Systems, Workshop on Physical Reasoning and Inductive Biases for the Real World, vol. 2, Dec. 2021.
[32] X. Li, T.-K. L. Wong, T. Q. Chen, and D. Duvenaud, "Scalable gradients for stochastic differential equations," " International Conference on Artificial Intelligence and Statistics, Aug. 2020.
[33] P. Kidger, J. Foster, X. Li, H. Oberhauser, and T. Lyons, "Neural SDEs as GANs," International Conference on Machine Learning, July 2021.
[34] Y. You, J. Li, S. Reddi, J. Hseu, S. Kumar, S. Bhojanapalli, X. Song, J. Demmel, K. Keutzer, and C.-J. Hsieh, "Large batch optimization for deep learning: Training BERT 76 minutes." International Conference on Learning Representations, Apr. 2020.
12

---

[35] T. Hiraiwa, "Dynamic Deff-Organization of Idealized Migrating Cells by Contact Communica- tion," Physical Review Letters, vol. 125, no. 26, p. 268104, 2020.
[36] B. Liebchen and H. Löwen, "Modeling chemotaxis of microswimmers: From individual to collective behavior," in Chemical Kinetics: Beyond the Textbook, (K. Lindenberg, R. Metzler, and G. Oshanin, eds.), p. 493-516, London, United Kingdom: World Scientific Publishing Europe, 2019.
[37] M. Wang, D. Zheng, Z. Ye, Q. Gan, M. Li, X. Song, J. Zhou, C. Ma, L. Yu, Y. Gai, T. Xiao, T. He, G. Karypis, J. Li, and Z. Zhang, "Deep Graph Library: A Graph-Centric, Highly-Performant Package for Graph Neural Networks," arXiv:1909.013/5 Aug. 2019.
[38] M. Novik, "torch-optimizer collection of optimization algorithms for PyTorch.." " Jan. 2020.
13

---

A Supplemental Material
Computational Resources
Detailed specifications of the computational setup are provided here to ensure reproducibility. We employed the Dell Precision 7920 Tower, which comprises 64 GB RAM, 2 Intel Xeon CPUs, and NVIDIA RTX A6000 GPUs. The system operates under Linux Ubuntu 22.04, with all neural network training processes managed through pyenv with python 3.8.13 8.13 maintain environment consistency. To use the GPUs in the computation, we used CUDA 11.5 and PyTorch 1.13.1.
Existing Assets
To construct our framework, utilized the Deep Graph Library 1.1.1 (https://www, dg1 ai/) (37), TorchSDE 0.2.5
http://hhbb.coooooggeeeeeeeeeeee Tor 1.0.4 hutpp:://github..cm/////ffff///
(32) (3),
(31), and torch-optimizer 0.3.0 luttps://githen com/ cowp:j/tiittyppoooroooooioiiiiii (38), all of which are publicly available under the MIT License and Apache License 2.0.
Supplemental Tables
Trial ID
Y
Mira
MSE(F!))
MAE(F.9)
MSE(F, (PE
MAE(F??
1
1 10
270
3.6 x 10
5.9 X 10-2
5.1 x 10
4.8 x 10-2
2
1 10 -
1.7
1.4
2.2 x 10
1.4 x
3
1 10-2
3
9.4
3.7
2.2 x 10
1.3 10-1
4
X 10-1
3
1.2 x 10-4
1.7 10-
2.2 x 10-2
1.0 x 10-1
5
x 10-1
3
5.4
x 10
3.8 x 10-2
2.3 x 10-2
1.1 x 10-1
Supplemental Table S1: Normalized MSE and MAE of the estimation of F()) and F(2) for the harmonic interaction model. Each row indicates independent experiment.
Parameter set
aCF(()
acF(l)
aci
acC (1)
0.1
0.9
2.0
0.2
0.9
0.5
0.5
0.5
(ii1)
0.9
0.9
0.2
0.2
(iv)
0.9
0.9
0.5
0.5
(v)
0.1
0.1
2.0
2.0
Supplemental Table S2:
Parameter used the simulations of the mixed-species model.
Trial ID
Parameter set
MSE(F!))
M)
MSE
MAE(F 9
MSE
MAE(F!
MSE(F
MAE(F
1.5 x 10
3.7x 10
4.5 10
6.7 10
5.6 10
277 10
24
28 x 10
2
(i)
1.7 x 10
3.5x 10
3.0 x 10
3,9 10
9.9 10
3.0 x 10-1
16 10
28 x 10
(i)
20x 10-3
42x 10
3.8 10
1 10
19x 10
16 10
5.4
20 x 10
5.8 10
23 10
7.1x
2.7 10
2.7
17 10
5.8
21 10
5
(18)
1.0 10-3
3.1 10
233 10-6
1.1 10
4.2 10
1. 10
6.8 10
2,5 x 10
6
(V)
1.6 10
13x 10
3.2 10
5.0 10
37 10
17 10
73 10
24 10-1
2 10
48x 10
1.9 10
28
19
8.8
1.6
8
6.1 10
8.0
3.1
10
3.9 10
2.3 10
1.1% 10
2.0 x 10-1
Supplemental Table S3: Normalized MSE and MAE of the estimation of F()) and F(2) for the
mixed-species model. Each row indicates an independent experiment. To note, MSE (F(1) and
MAE(F()) are not normalized since the true value is zero in any case.
14

---

Supplemental figures
Original Model
Prodiction
Original Model
Prodiction
(A)
(B)
200.0
200.0
200.0
== 200.0
15
Punnnder adt
6 19
(F)
(G)
(HI
200..0
200.0
200.0
Paramete set (1)
1e
10
20
10 20
Supplemental Figure S1: Snapshots of the simulation results in the mixed-species model and the
estimated model. Panels (A, C, E, G) depict results from the mixed-species model, and panels (B, D, F, H) from the estimated model. The upper rows (A-D) represent simulations with parameter set (i), and the lower rows (E-H) with parameter (ii). Panels (A, B, E, F) include individuals with ct
and panels (C, D, G. H) with cf 1. Colors correspond to those used in Figure
15

---

(B)
d/w0
eww1
Puramets
(F)
(G)
0000
Parameter
830
(K)
Parameter
2411
x10
Supplemental Figure S2: The functions estimated from data with the mixed-species model, without dependency on species type, plotted against the true values. The rows display different cases: (A-D) with parameter set (iii), (E-H) with parameter set (iv), and (I-L) with parameter set (v). A gray line serves as guide to indicate where the estimated values equal the true values.
Supplemental Movies
Supplemental Movie S1:
A visualization of the simulation based on the harmonic interaction model with friction constant 7 = x 10
Supplemental Movie S2:
A visualization of the simulation based on the harmonic interaction model with friction constant
7 = 11x x 10-1
Supplemental Movie S3:
The results of simulations using the estimated functions (1) and F (2) Ni = 270 data from the harmonic interaction model with ? 10
NN. which were trained using
Supplemental Movie S4:
The results of simulations using the estimated functions P and F NN which were trained using
Nm = data from the harmonic interaction model with 7 x
Supplemental Movie S5:
The results of simulations using the estimated functions (1) and (2)
NN which were trained using
Nm data from the harmonic interaction model with y 10
Supplemental Movie S6:
A visualization of the simulation based on the mixed-species model with parameter set (i).
16

---

Supplemental Movie S7:
visualization of the simulation based on the mixed-species model with parameter (ii).
Supplemental Movie S8:
The results of simulations using the estimated functions F and FNN which were trained using Nua = 3 data from the mixed-species model with parameter set (i).
Supplemental Movie $9:
The results of simulations using the estimated functions P and NN which were trained using Na 3 data from the mixed-species model with parameter set (ii).
17

---

NeurIPS Paper Checklist
Claims
Question: Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification:
See Section 1.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made in the paper and important assumptions and limitations. A No or
NA answer to this question will not perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results be expected to generalize other settings.
is fine include aspirational goals as motivation as long is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification: See Section
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those not discussed in the paper.
The authors are encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification. asymptotic approximations only holding locally). The authors should reflect how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low images are taken in low lighting. Or speech-to-text system might be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach address problems of privacy and faimess.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play an impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory Assumptions and Proofs
Question:
For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [NA]
18

---

Justification:
The paper do not include theoretical results.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated referenced in the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but they appear in the supplemental material, the authors are encouraged to provide short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental Result Reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that affects the main claims and/or conclusions of the paper (regardless whether the code and data are provided or not)?
Answer: [Yes]
Justification: See Section 4-5.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g.. in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make clear how to reproduce that algorithm. (b) the contribution is primarily a new model architecture, the paper should describe the architccture clearly and fully.
(c) If the contribution is new model (e.g., a large language model), then there should either way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reprodueibility. the case of closed-source models, may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers have some path to reproducing or verifying the results.
5. Open access to data and code
Question: Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
19

---

Answer: [Yes]
Justification:
See the Supplemental Material.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurIPS code and data submission guidelines (https://nips. cc/
for more details.
While we encourage the release of code and data, we understand that this might not be possible, so "No" an acceptable answer. Papers cannot be rejected simply for not including code, unless this is central to the contribution (e.g., for new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https
/t.pss.ccppullc//uuii / Swweleee//ooSSeeeeeeeeeeee
www.SSubbissioooiiiooioooooo for more details.
The authors should provide instructions on data access and preparation, including how access the raw data, preprocessed data, intermediate data, and generated data,
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing as much information as possible in supplemental material (appended the paper) is recommended, but including URLs to data and code is permitted.
6. Experimental Setting/Details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: See Section 4-5.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core the paper to a level detail that is necessary appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, supplemental material.
7. Experiment Statistical Significance
Does the paper report error bars suitably and correctly defined other appropriate information about the statistical significance of the experiments?
Question:
Answer: [No]
Justification: Due the extensive computation time required for estimation, comprehensive statistics could not be gathered. Instead, we present the results of all trials to illustrate the trends in estimation accuracy (see Supplemental Tables 1,3).
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper.
The factors of variability that the error bars capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.)
The assumptions made should be given (e.g., Normally distributed errors).
20

---

It should be clear whether the error bar is the standard deviation or the standard error of the mean.
is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not to show in tables figures symmetric error bars that would yield results that are out range (e.g. negative error rates).
If error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures tables in the text.
8. Experiments Compute Resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?
Question:
Answer: [Yes]
Justification:
See the Supplemental Material.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper).
9. Code Of Ethics
Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics http:///wara.cccca//a/////aacaaa
Answer: [NA]
Justification:
Guidelines:
The answer NA means that the authors have reviewed the NeurIPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- eration due to laws or regulations in their jurisdiction).
10. Broader Impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [No]
Justification:
Guidelines:
The answer NA means that there no societal impact of the work performed.
If the authors answer NA or No. they should explain why their work has no socictal impact why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
21

---

The conference expects that many papers will be foundational research and not tied particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology.
If there negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms monitor how system learns from feedback over time, improving the efficiency and accessibility of MI..
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [No]
Justification:
Guidelines:
The answer NA means that the paper poses such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Lieenses for existing assets
Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper. properly credited and are the license and terms of use explicitly mentioned and properly respected?
Answer: [Yes]
Justification: See the Supplemental Material
Guidelines:
The answer NA means that the paper does existing assets.
The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset is used and, if possible, include URL.
The name of the license (c.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (c.g., website), the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswi thcode. com com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided.
22

---

If this information is not available online, the authors are encouraged to reach out to the asset's 's creators.
13. New Assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [Yes]
Justification: See the Supplemental Materials.
Guidelines:
The answer NA means that the paper does not release new assets,
Researchers should communicate the details of the dataset/code/moodl as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsourcing and Research with Human Subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well details about compensation (if any)?
Answer: [NA]
Justification:
The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsoureing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According the NeurIPS Code of Ethies, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional Review Board (IRB) Approvals or Equivalent for Research with Human Subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or an equivalent approval/review based on the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: The paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurlPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
23


Reviews:
[{"review_id": "HsDIs21vDr", "reviewer": "Reviewer_KAZi", "summary": "This paper seeks to model dynamic interactions of biological agents (primarily cellular slimes) by adapting the (continuous time) Langevin model, which models the collective motion of the system from considering individual force terms and the summation of pairwise interactions. Such pairwise interactions can be modelled as a graph; thus the others model the interaction terms with message passing on a graph and the solution to the differential equations are modelled with a neural ODE. The main contribution seems to be the learning of the force terms with fully connected neural networks. The authors claim that the main problem with previous works is that no models have been proposed that can model systems that interact in unknown ways; however, the system is trained in a supervised way based on known simulations. The results suggest the learning model can replicate the simulation,  I don’t believe there was any attempt to extrapolate beyond the training regime.", "questions": "1. The main contribution seems to be to use a fully connected network to learn the force terms through a data term that forces the model to match the simulation results. Why not use a Physics informed neural network to embed the physics into the model, and thus hopefully extrapolate beyond the training regime?\n2. Wouldn’t it make more sense to model this as a graph ODE? http://arxiv.org/pdf/1911.07532\n3. It would help future readers if the authors explain in more detail how the whole network is trained - including where the gnn and neural ODE are called (it took me a long time to work this out as they are only briefly mentioned in the text) \n4. This paper seems like it might be relevant and seems to be going beyond what is proposed here through learning models of interactions https://arxiv.org/abs/2303.09906 with PiNNs - can the authors explain the advs/disadvs of their work relative to an approach such as this", "limitations": "Limitations are not really discussed beyond runtimes", "rating": 3, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "1. The paper describes a very interesting problem which (to my knowledge) hasn’t been extensively explored through learning based approaches\n2. It uses state-of-the-art message passing gnns and Neural ODEs as part of the \n3. The model is accurately replicating simulation behaviour\n4. Code is available", "weaknesses": "1. The paper is really hard to read. It took me a very long time to understand how it was implemented and although the code was generously submitted, this isn’t commented making it difficult to follow in the time available. Gnns and ODEs are in the title of the paper but are barely mentioned in the text. \n2. The methodological contribution is limited\n3. The model is trained in a supervised way meaning that it cannot learn new behaviours or extrapolate beyond the training regime. As this was stated as the main limitation of previous works it suggests this isn’t adding anything beyond these?\n4. Results are based on validating relative to simple simulations - and run time is apparently anyway slow, so what is the value of the proposed approach", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for the feedback and comments. Your careful reading was very helpful in identifying places that we missed to explain in the draft.\n\nRe – Question 1 and Weakness 2 (Why not use a Physics informed neural network):\n\nThis was briefly described in the Background section, however it was not fully conveyed in relation to Physics-Informed-Neural-Networks (PINN). While there is no strict definition of PINN, we believe our model absolutely falls into that category, since we build our architecture on equations of motion for active matter (a class of coupled equations of motion important in non-equilibrium and biological systems) with forces being represented by neural networks. Most works on PINN centers on the estimation of relatively simple close to equilibrium systems such as diffusion and fluid dynamics where semi or almost non-supervised learning is partially possible based on first principles. Instead, actively propelled particles as such studied here are far from equilibrium and Hamiltonian is not available. It is that we are exploring this difficult class of systems which is nevertheless fundamental in biological physics that makes the study novel and important.\n\nRe – Question 2 & Weakness 2 (Relationship with Graph ODEs):\n\nOur model is in fact a graph ODE except that the graph is time evolving due to change in particle positions. Valid edges were detected at every timepoint based on particle-to-particle distance. Existing graph ODEs instead have explicitly defined (i.e. fixed or semi-fixed) graphs. We realized that this wasn’t fully described and would be happy to clarify in the final edition.  \n\nRe – Question 3 & Weakness 1 (Clarify the details of training):\n\nWe recognize now that the section that describes this (p.4 L. 130-136) should have been clearer. Apart from adding more comments in the codes, we propose to reorganize and extend the Method section, and add a flow chart in the final version. The reorganized section would make the following clear:\n- Both for generating the training data with predetermined forces and the learning approach with neural networks, we solve a differential equation, eq. (1) in the original draft, over time. Depending on whether there is noise in the system, we employed either neural ODE or neural SDE. By neural ODE we simply mean that we solve an ODE in which the forces are given by neural networks. In either case, the r.h.s. of eq. (1) has to be evaluated at least once per time step.\n- To generate the training data, predetermined force functions were simply evaluated. During the learning process, the forces are represented by neural networks. Graph neural networks (GNN) techniques were employed to make this evaluation efficient. We note that the state of the system at any given time can be expressed by a graph, with each entity representing a vertex, and interactions between pairs of entities representing edges. As the positions evolve over time, the graph is time-dependent, too.\n- Whenever the r.h.s. of eq. (1) needs to be evaluated by the solver, we pass the system state to the GNN wrapper. The GNN wrapper defines the graph structure and performs a message passing process, which evaluates the neural network (or predetermined function to generate training data) of force F^{(1)} for each vertex, evaluates the neural network for the force F^{(2)} for each edge, and finally adds the forces on the edges to each connecting vertex.\n- Loss function: In order to infer the forces F^{(1)} and F^{(2)}, we split training data into pairs of states which are a time tau apart. Then, for each pair of states, we integrate the equation of motion from the first state for a time tau using the neural network forces and then compare the resulting state to the second state. This exploits that the noise in the training data has the tendency to partially cancel out over time.\n\nRe – Question 4 (advantages and disadvantages of our method in comparison with arxiv:2303.09906)\n\nThe preprint addresses learning the dynamics of an order parameter of alignment of agents instead of the underlying full trajectory of the interacting agents. The order parameter encodes some of the information about the collective dynamics of the agents, since it is an average over the agent’s states at a given time, but it does not come close to capturing their full dynamics. Their approach infers the dynamics of a one-dimensional quantity, whereas we infer the full dynamics of all the agents in the system. Their approach is useful in situations where one is only interested in understanding the dynamics of a mesoscopic or macroscopic quantity such as an order parameter and where having to keep track of the complete dynamics of the system is a hindrance, e.g. in terms of data storage or computational cost. Their approach however cannot be applied to inferring the microscopic interactions between the agents (because they are intentionally averaged out), which is the aim of this study.\n\nRe – Weaknesses 3 & 4 (Unable to extrapolate as other previous works due to supervised approach. The value of the proposed approach):\n\nThis is related to Question1. We should note that in complex systems where first principle is not easily available, an unsupervised approach is not feasible. To the best of our knowledge, the only studies of learning the interaction of actively propelled particles is on the Viscek model (cited refs [4, 24]) which still has an order parameter that can be exploited. The aim of our work instead is to go beyond such system where no clear order parameter is given. From application point of view, too, the interaction term in Viscek model leaves out crucial aspects of interaction between biological cells and thus we turned to the current model equation. Whether learning such models is even possible was not at all trivial and this work demonstrates its first realization and thus novel. We realized that this aspect was not sufficiently explained and will be amended in the corrected edition.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "yiZvY3aakk", "reviewer": "Reviewer_Shpk", "summary": "The paper introduces a novel deep learning framework for estimating two-body interactions in mixed-species collective motion. The authors combine graph neural networks (GNNs) with neural ordinary differential equations (Neural ODEs) to predict interactions between pairs of entities based on their states. This approach represents the system as a dynamic graph, using GNNs to efficiently calculate interactions and Neural ODEs to learn system dynamics. Authors demonstrate their method via two experiments: a simple harmonic interaction model and a more complex model simulating cellular slime molds. The framework successfully estimates interaction functions and replicates both individual and collective behaviors in these systems. The authors provide detailed quantitative analyses of estimation accuracy and demonstrate the method’s ability to generalize to single-species scenarios after training on mixed-species data.", "questions": "Questions for the authors:\n1. How does the performance of your method compare to existing approaches for estimating collective motion dynamics, such as SINDy or Bayesian optimization methods?\n2. Have you explored the method’s sensitivity to noise in the input data? How robust is the estimation process to measurement errors or stochastic fluctuations in the trajectories?\n3. I’m intrigued by the use of LAMB optimizer. What specific characteristics of your problem or empirical results led you to choose LAMB over AdamW?\n4. What are the primary factors limiting the scalability of your approach to systems with stochastic motion equations and higher-order interactions?", "limitations": "yes", "rating": 7, "confidence": 4, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "Strengths:\n- Novel integration of GNNs and Neural ODEs for collective motion analysis\n- Ability to handle mixed-species systems and estimate species-dependent interactions (and demonstrating its effectiveness across simple and complex models)\n- Successful generalization to single-species scenarios after training on mixed-species data\n- Potential applicability to a wide range of biological systems exhibiting collective motion", "weaknesses": "Weaknesses:\n- High computational cost and long training times .\n- Current limitation to deterministic motion equations and pairwise interactions.\n- Lack of ablative studies on hyper parameter optimizations and method’s sensitivity to noise.\n- Lack comparative analysis with existing methods for estimating collective motion dynamics.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for the feedback and comments.    We are encouraged by your overall evaluation and many of the technical suggestions which are very helpful in future extension of the work.\n\nRe – Question 1 and Weakness 4 (Comparative analysis with SINDy and Bayesian optimization methods):\n\nThank you for pointing out the lack of comparative analysis. Although both SINDy and Bayesian optimization are in principle applicable to our target system, they require drastic change of codes. Thus we have not tested them. \n\nRe – Question 2 and Weakness 3 (The method’s sensitivity to noise in the input data):\n\nNo, we have not systematically studied dependency on the noise parameter, nor introduced additional noise, and this is on-going and part of our future plans.  Due to the very fact that stochastic differential equations are employed to generate the dataset and to train ODEs, we expect the method itself to be robust against stochasticity. In the mixed-species case, linear relationship between the true vs estimation was robust (fig. 2F and 2J).  We do however see consistent under-estimation, albeit slight, of the interaction term of \\phi (F^{(2)}_\\phi) which presumably is caused by the stochasticity that is added to d\\phi/dt. \n\nRe – Question 3 (Advantage to use LAMB optimizer in comparison with AdamW):\n\nOur module has a layer where the output of the fully connected network is multiplied by a trainable scalar. When using Adam, the scalar behaved erratically compared with weights and biases in fully connected layers. The ability of LAMB to automatically adjust the learning rate between layers helped us stably train the module. We are aware that there is also an option to set a different learning rate for the scalar using AdamW, which we did not test.\n\nRe – Question 4 and Weakness 2 (Primary factors which limits the scalability of our approach to systems with stochastic motion equations and higher-order interactions):\n\nComputational load is the primary factor when extending the present approach to stochastic systems.  With the current stand-alone off-the shelf GPU in the market (NVIDIA A6000), it takes weeks to train ODEs.  Since the bottleneck seems to be at neural ODE calculation, performing bayesian inference to fit SDEs would take an unrealistically long time.  Some form of parallelization and use of super computers may help resolve this issue in the future. As for higher-order interactions, the limitation is the lack of a good library to deal with hypergraphs, at least to the best of our knowledge. The other option is to use message passing twice, which we have not tried, because of the caveat of making the training and interpretation of neural networks more difficult.\n\nRe – Weakness 1 (High computational cost and long training times):\n\nSince our target of estimation is pairwise force which would work as a total force to drive the entities, we are essentially dealing with an inverse problem. The costs are necessary for stable and accurate estimation of each pairwise force against noise.  Given the nature of the problem, we believe costs of memory and time for estimation is a difficult issue to resolve.  Still, as we mentioned in the last section, we agree that it is an important avenue to explore in future research.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "9NGnam15N1", "reviewer": "Reviewer_T7k9", "summary": "This paper presents a model combining GNN and Neural ODEs for the task of two-body interaction prediction. Two numerical experiments were conducted to evaluate the model's performance.", "questions": "The author should discuss in their response regarding the missing related works that involve integrating GNN+Neural ODEs, and highlight their differences/contribution from previous work.", "limitations": "Limitations were discussed in the paper.", "rating": 5, "confidence": 4, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "- Using a GNN and Neural ODE to describe the dynamics of two-body interaction state changes is well-motivated and reasonable.\n- Experiments from interaction models are implemented from simple to complex setups, showing the efficiency and capacity of the proposed method.\n- The paper structure and writing is clear and easy to follow.", "weaknesses": "One of the main contribution claimed in this paper (as well as in its title) is the combination of GNN and Neural ODE, but this approach has been studied in multiple works back since 2020, the authors seem to overlook the important and related body of works [1,2,3,4]. [2] is a direct competing method (GDE) where the authors proposed a model integrating graph NN and Neural ODEs. GDE-related applications exist as well, for example, social network embedding [5], action recognition [6]. \n\nGiven the missing discussion of related work, and my concern about the major technical contribution claim in this paper, I do not see much merits of this study to be novel and interesting to the community.\n\n[1] L. Xhonneux et al.: Continuous Graph Neural Networks. PMLR 2020.\n[2] M. Poli et al.: Graph Neural Ordinary Differential Equations. AAAI 2020.\n[3] L. Chen et al.: Signed Graph Neural Ordinary Differential Equation for Modeling Continuous-Time Dynamics. AAAI 2024.\n[4] A. Han et al.: From Continuous Dynamics to Graph Neural Networks: Neural Diffusion and Beyond. TMLR 2024.\n[5] Y. Zhang et al.: Improving Social Network Embedding via New Second-Order Continuous Graph Neural Networks. KDD 2022.\n[6] L. Pan et al.: Spatial-temporal graph neural ODE networks for skeleton-based action recognition. Scientific Reports 2024.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "Thank you for bringing to our attention these papers.  We will cite these papers in the final edition.\n\nRe – Question and Weakness (Relationship to Graph ODE):\n\nThe pre-existing methods [1-5] solve and train ODEs on a graph with fixed edges and consider changes in the edge weight. These approaches require a fully or almost fully connected graph in solving ODE for an application to collective motion where the adjacency drastically changes in time, such as collective cell migration, flocks, and pedestrian dynamics. In such systems addressing the adjacency at each time is critical because the behavior of individuals is determined by sensible information, that is, the state of neighboring others. This is computationally expensive making it difficult to be extrapolated to situations with unknown graph structures.   In contrast, since edge structures are continuously being updated in our approach, the pruning effect substantially reduces memory requirement - in our case of 400-body simulations with the edge density of 2%,  the required memory of 200 GB for fully connected graph is reduced to 30GB making it feasible on an off-the-shelf GPU board. As for reference [6], although the method introduces a dynamic edge structure, the time series of the graph (edge) structure needs to be estimated first using an entire time series and thus the trained model cannot be further extrapolated to situations outside of the particular graph structure. Our method resolves this issue by setting a rule to define edges at each time point instead of explicitly giving graph structure time series.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "jj0kjsaBN0", "reviewer": "Reviewer_Sm49", "summary": "This work proposes a graph neural network-based deep learning model for two-body interactions, which has better efficiency based on two numerical experiments.", "questions": "see strengths", "limitations": "see strengths", "rating": 7, "confidence": 1, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "Two-body interactions are not my research field. Thus, I don't have enough expertise to evaluate this work. I sent an email to AC to re-assign this paper when I received the assignment. Unfortunately, I didn't receive any feedback. My justification is just based on the writing without any methodology consideration.", "weaknesses": "see strengths", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:
We thank reviewers for careful consideration of our work and are thankful for their astute criticism, which was very helpful in vastly improving the quality of the manuscript.  The manuscript was reviewed by 4 reviewers, with two of them (Shpk and Sm49) recommending acceptance. The other two, KAZi and T7k9, expressed disapproval out of concerns of novelty of the work which we believe are mainly due to the relatively short description provided in our original draft on the research background and the methodology.  

We have provided a point-by-point response to the reviewers’ comments which we incorporate in an updated draft.  We feel strongly that we have been able to address the concerns of the referees. In particular, we would like to briefly point out the following modifications:
- Referees KAZi and T7k9 both asked about the relationship to graph ODE.  We clarify the relationship of our method with pre-existing approaches to combine GNNs and neural ODE. Our new approach does not require any given graph structure, but defines it using a given rule at every time step of neural ODE calculation. This enabled us to extrapolate the trained model to unknown initial conditions, even in case the adjacency drastically changes in time.
- Describe Methods more clearly. In particular, how Graph Neural Networks (GNNs)  are used to evaluate the equations of motion effectively. We added a flow chart which describes how our model evaluates the equations of motion (attached PDF file). 
- Describe Physics Informed Neural Network and where the present work sits relative to related literatures.
- Referee KAZi expressed difficulty in readability while referee T7k9 found the writing to be rather strong.  Grammar and expression will be checked and corrected.

Meta-review: The paper studies the problem of learning models for collective behavior from data. While the initial reviews were quite mixed, the rebuttal and discussion phase lead to some convergence in terms of understanding, clarity and the overall score (including a score increase that is clearly reflected in the comments, but not entered correctly). From my point of view, this is a quite interesting paper, and while there are certainly some edges to be ironed out in terms of presentation and clarity (see reviewer comments), I am recommending acceptance at this point. I am also strongly encouraging the authors to take all comments seriously and adjust the manuscript accordingly.
Decision: Accept (poster)