Learning Spatiotemporal Dynamical Systems from Point Process Observations
Authors: Valerii Iakovlev
Harri Lähdesmäki
Venue: ICLR 2025 Spotlight (2025)
Keywords: dynamics
spatiotemporal
neural
PDE
ODE
Primary area: generative models

Abstract: Spatiotemporal dynamics models are fundamental for various domains, from heat propagation in materials to oceanic and atmospheric flows. However, currently available neural network-based spatiotemporal modeling approaches fall short when faced with data that is collected randomly over time and space, as is often the case with sensor networks in real-world applications like crowdsourced earthquake detection or pollution monitoring. In response, we developed a new method that can effectively learn spatiotemporal dynamics from such point process observations. Our model integrates techniques from neural differential equations, neural point processes, implicit neural representations and amortized variational inference to model both the dynamics of the system and the probabilistic locations and timings of observations. It outperforms existing methods on challenging spatiotemporal datasets by offering substantial improvements in predictive accuracy and computational efficiency, making it a useful tool for modeling and understanding complex dynamical systems observed under realistic, unconstrained conditions.

Published as a conference paper at ICLR 2025
LEARNING SPATIOTEMPORAL DYNAMICAL SYSTEMS FROM POINT PROCESS OBSERVATIONS
Valerii lakovlev Harri Lähdesmäki
Department of Computer Science, Aalto University, Finland weterriiiikvvleeee harri. lalllaaaaaaaii eeaaeiiieaaa.....
ABSTRACT
Spatiotemporal dynamics models are fundamental for various domains, from heat propagation in materials to oceanic and atmospheric flows. However, currently available neural network-based spatiotemporal modeling approaches fall short when faced with data that is collected randomly over time and space, as is often the case with sensor networks in real-world applications like crowdsourced earthquake detection pollution monitoring. In response, we developed new method that can effectively learn spatiotemporal dynamics from such point process observations. Our model integrates techniques from neural differential equations. neural point processes, implicit neural representations and amortized variational inference to model both the dynamics of the system and the probabilistic locations and timings observations. It outperforms existing methods on challenging spatiotemporal datasets by offering substantial improvements predictive accuracy and computa- tional efficiency, making a useful tool for modeling and understanding complex dynamical systems observed under realistic, unconstrained conditions.
INTRODUCTION
In this work, consider the modeling of spatiotemporal dynam- ical systems whose dynamics are driven by partial differential equations. Such systems ubiquitous and range from heat propagation in microstructures to the dynamics of oceanic cur- rents. We the data-driven modeling approach, where observe a system at various time points and spatial locations and use the collected data to learn model. In practice, the data is often collected by sensor networks that make measurements random time points and random spatial locations. Such measurement approach is used, for example, in crowdsourced earthquake monitoring where smartphones are used as mea- surement devices (Minson et al., 2015; Kong et al., 201 6), in oceanographic monitoring where measurements are made by floating buoys (Albaladejo al., 2010; Xu et al., 2014; Marin- Perianu et al., 2008), and for air pollution monitoring with vehicle-moune sensors (Ma et al., 2008; Ghanem et al., 2004). This approach offers several advantages as it needs no sensor synchronization and allows the sensors to move freely. How- ever, its random nature makes modeling more challenging as it requires capturing both the system dynamics and the random observation process.
f++1
6+1 p(1++1)
p(a,)
Previous models as- sume dense observations at every time point and fixed spatiotemporal grids. Bottom: Our model works with extremely sparse observations and predicts where the next obser- vations will happen.
Figure 1: Top:
We focus on the marked spatiotemporal point process setting, where the randomness in observation times and locations defined in terms of the state of the underlying spatiotemporal dynamical system and marks are generated from the system state via an observation function. To the best of our knowledge, existing neural network-based methods cannot model such randomly observed spatiotemporal dynamical systems, mainly for two reasons. First, they do not model observation times and locations, hence they cannot predict where and when the next observation will be made.
Source code and datasets found in our github repository.

---

Published as a conference paper at ICLR 2025
Second, they assume the sensors form fixed spatial grid and record data simultaneously, which is not the case in our setup where the data could come from as little as a single sensor at each time point (Fig. 1). For example, some earlier methods assume the observations spatiotemporal grid (Long et al., 2018; Geneva & Zabaras, 2020), other methods work with irregular, but still fixed, observation locations (lakovlev al., 2021; Lienen & Günnemann, 2022), and other works go further and allow the observation locations to change over time (Pfaff' et al., 2021; Yin et al., 2023) but fix the observation times and assume dense observations. Whereas, another line of research has proposed methods that model only the spatiotemporal observation process without modeling the system dynamics (Chen et al.. 2021; Zhu et al., 2021; Zhou et al., 2022; Zhou & Yu, 2023; Du et al., 2016).
made fixed and regular
Our work fills this gap and proposes model for randomly observed spatiotemporal dynamical systems. Our model incorporates techniques from amortized variational inference (Kingma & Welling, 2013), neural differential equations (Chen et al., 2018; Rackauckas al., 2020), neural point processes (Mei & Eisner, 2017; Chen et al., 2021), and implicit neural representations (Chen et al., 2023; Yin et al., 2023) to efficiently learn both the underlying system dynamics and the random observation process. Our model uses initial observations to obtain the variational estimate of the latent initial state via a transformer encoder (Vaswani et al., 2017),, simulates the latent trajectory with neural ODEs (Chen et al., 2018), and uses implicit neural representations parameterize the point process and observation distribution. Furthermore, we identify a computational bottleneck in the latent state evaluation and propose technique to alleviate it, resulting in up to 4x faster training. Our model shows strong empirical results outperforming other models from the literature on challenging spatiotemporal datasets.
BACKGROUND
2.1 SPATIOTEMPORAL POINT PROCESSES
Spatiotemporal point processes (STPP) model sequences of events occurring in space and time. Each event has an associated event time I E Rzo and event location a; RR Given an event history H4 _ {(t,, m;) f t} with all events up to time f, we can characterize STPP by its conditional intensity function
P(t, E [t, + öt) x, E (x))
(t, 2x)
lim 5e40, 60
5t B5)(5)))
where 8t denotes infinitesimal time interval, and (ar denotes 8r-ball centered at æ. Given history H₂ with i 1 events, A (t,t,, x) describes the instantaneous probability of the next, ith event occurring time and location Given sequence of N events {(t;, x.)) bounded domain A C [0, 7 x Rd,, the fog-likelihook for the STPP is evaluated as (Daley et al., 2003)
N
S X² (t, z))ddt.
(ep((((t,2:))13)
Σ log X" (t;,
Marked STPP extends the above simple STPP by mark y; Bd that is associated to each event
(fe, T.).
2.2 ORDINARY AND PARTIAL DIFFERENTAAL EQUATIONS
Given a deterministic continuous- time dynamic system with state xt) E Rd we can describe the evolution of its state in terms of an ordinary differential equation (ODE)
dz(1) dt (t)).
For an initial state 21 at time t we can solve the ODE to obtain the system state xt)) at later times t > fj. The solution exists and is unique if is continuous in time and Lipschitz continuous in state (Coddington et al., 1956), and can be obtained either analytically or using numerical ODE solvers (Heirer al., 1987). In this work we solve ODEs numerically using ODE solvers from
2

---

Published as a conference paper at ICLR 2025
torchdi (Chen, 2018) package. Similarly, the dynamics of spatiotemporal systems with state *(t, x) defined over both space and time is described in terms of a partial differential equation (PDE)
(z x) F(1, =(, (l)
0 which incorporates both temporal and spatial derivatives, indicated by ae and V, respectively.
PROBLEM SETUP
In this work, we model spatiotemporal dynamical systems from data. The data consists of multiple trajectories collected by observing a system over a period of time. trajectory consists of N triplets {(t,, :. y:) ) where E R' is a system observation, with & Ry0 and I; E R being the corresponding observation time and location. The number of observations N can change across different trajectories in the dataset. Due to randomness of the observation process we assume the observation times and locations do not overlap (neither within a trajectory nor across trajectories), resulting in single observation y: per time point and location. For brevity, we describe our method for a single observed trajectory, but extension to multiple trajectories is straightforward.
We assume the data generating process (DGP) consists of a latent spatiotemporal state uul1, a)) E Rd with E 2g0 and æ E Hd*, whose dynamics are governed by a PDE. The latent state directly affects the times and locations of the observations, and is only partially observed:
@u(1, F(u(t, Vu(t, x))),
-(t,,)) ~ p(u).
(5)
Öt
t, x ~ nhpp()(u(t,x)))) y: (ut(
1
N
(6)
N.
(7)
According to this process, each trajectory is generated as follows. We first sample the latent field .lt1,)) initial time point t and specify its dynamics by PDE (Eq. 5). Next, sample the observation times t; and locations æ; from a non-honooeeee non-homogeneous Poisson process (nhpp) with intensity sample y; from the observation distribution
that is function of the latent state (Eq. 6). Finally, parameterized by a function g (Eq. 7). All datasets in this work were generated according to this process (see Appendix for details). We assume the data generating process fully unknown, and our goal is to construct and learn its model from the data.
4 METHODS
In this section we describe our proposed model, its components, and the parameter inference method.
4.1 MODEL
Our goal is to model the true DGP (Eqs. 5-7). To this end, we define our generative model as:
dz(r) f(z(r)),
2(t1) p(z(tz)).
(8)
dt
u(t, x) o(z(t), x),
(9)
t. x, ~ (bpp(>(u((t())))
= 1,..., N.
(10)
y: ~ ply,,g(u(tt, #))),
(11)
and the corresponding joint distribution is (details in App. B):
N II
p(z1, (ti, 2;, M = ( 12))
a
p(y 21, 14, m.).
(12)
i=1
Below, we descrihe each component in detail, and diagram of our model is shown in Figure 2.
3

---

Published as a conference paper at ICLR 2025
Latent dynamics (Eq. 8). Our goal is to model the latent PDE dynamics of the DGP (Eq. 5). To do that, we introduce cow-dimensional temporal latent state x(T) E Rd with 7 E R>0 that encodes the spatiotemppral field u(T, )), and introduce neural network- parameterizee ODE governing its dynamics (Eq. 8). We use a low-dimensional latent state because it allows to simulate the dynamics considerably faster than a full-grid spatiotemporal discretization (Wu et al., 2022; Yin et al., 2023). During our experiments, we observed that it takes ODE solvers long time to solve the ODE and that this time increases with the number of observations. We found this happens because of bottleneck in existing ODE solvers caused by extremely dense time grids (such as in Fig. 3). Such time grids foree the solver to choose small steps which prevents it from adaptively selecting the optimal step size thus reducing its efficiency. To alleviate this bottleneck, we move away from the original dense time grid l1....,4y iN to an auxiliary sparse grid TI Ta. with n. E N, where the first and last time points coincide with the original grid: T1 = [1 Tn Ly (see Fig. 2). Then, we solve for the latent state (t) only at the sparse grid, which allows the ODE solver to choose the optimal step size and results in approximately an order of magnitude (up to 9x) faster simulations (see Sec. 5). Finally, we use the computed latent states z(n))....,z) ,((T)) to approximate the original latent state x() at intermediate time points via interpolation: Et) interpolate(t; ..r)),....,.......... (,,7),...,,..(..
Latent state decoding (Eq. 9). Next, we need to recover the latent spatiotemporal state u(t, m) from the low-dimensional representation (()) We define the latent spatiotemporal state utt, æ;) as a function of (t)(t), x) of the latent state (t) and evaluation location æ (Eq. 9). In particular, we the input and maps it directly to IRd a trainable linear mapping of x to R This formulation results in space-time continuous 2utt, æ) which can be evaluated at any time point and spatial location, which is
parameterize of (t)( x) by an MLP which takes (t) + proj(z)
utt, x), where proj(z)
required to define the intensity function (tu(t, x)) and observation function g g(u(t, x)).
Intensity function (Eq. 10). We parameterize the intensity function X(u(t, æ:) by an MLP which takes utt, æ)) as the input and maps it to intensity of the Poisson process (Eq. 10). To ensure non- negative values of the intensity function and improve its numerical stability we further exponentiate the output of the MLP and add small constant to it.
Observation function (Eq. 11). The observation function (tt(t, x)) maps the latent spatiotemporal state utt, x) parameters of the observation model (Eq. 11). Therefore, its exact structure depends on the observation model. In this work the observation model is a normal distribution, where the variance is fixed and the observation function g(t(t, returns the mean.
Details of the model specification and parameterizaaion can be found in Appendix B and C.
4.2 PARAMETER AND LATENT STATE INFERENCE
To infer the model parameters and posterior of the latent initial state 21 we use variational inference (Blei et al., 2017). We define an approximate posterior ((11; Φ) with variational parameters e to approximate the true posterior p(zi)(4, n; Ni and then minimize the Kullback-Leible divergence
KL
(4i, Fz (.)
(13)
over the model and variational parameters obtain an estimate of the model parameters and approximation of the posterior.
To avoid optimizing the local variational parameters di for each trajectory in the dataset, we use amortization (Kingma & Welling, 2013) and define ub as the output of an encoder:
3 EEnoder 24, y y.
(14)
which maps observations (t,, 3e, y ) to the local variational parameters th This allows to optimize a fixed set of encoder parameters instead of a dataset size-dependent set of local variational parameters. We discuss the structure of the encoder in the next section.
4

---

Published as a conference paper at ICLR 2025
Latent Space
*(n) 21
a(2))
23
Encoder
24
2(4)
9
Context
72
LN-T1
Figure 2: Model diagram. Black crosses on the time axis are observation times, while white crosses on field images the corresponding observation locations. The initial observations (context) mapped by the encoder to the latent initial state distribution q(21). ). The initial state 2 is then sampled from that distribution and evolved through time using the dynamics model. The latent trajectory is evaluated only at the sparse time grid Tl..... Fn and the latent state x(t) at other time points is evaluated via interpolation. The latent state (() then mapped by é to the spatiotemporal state (t, æ), which is used to parameterize the point process and observation distribution (via mappings À and g. respectively) for predicting subsequent observation times and locations.
In practice, instead of minimizing the KL divergence we, equivalently, maximize the evidence lower bound (ELBO) which for our model is defined as:
c
Σ (=isV) [In (ly t4 24, #1) + Bg(
Σ In ((u(t,, - - (u(t,x)))d
(15)
(8) Expected observation log-lik.
Expected STPP log-lik.
KT [[][((21)) (1 +
(16)
(86) KL between prior and posterior
The ELBO is maximized wrt. the model and encoder parameters. Appendix B contains detailed derivation of the ELBO, and fully specifies the model and the approximate posterior. While the term (iii) can be computed analytically, computation of terms (i) and (ii) involves approximations: Monte Carlo integration for the expectations and intensity integral, and numerical ODE solvers for the solution of the initial value problems. Appendix B details the computation of ELBO.
4.3 ENCODER
Our encoder maps the observations {4i, 2;, ) to the local variational parameters & of ()) using Eq. 14. This process begins by embedding each observation into a high-dimensional feature space. The embedded sequence is then processed by a stack of transformer encoder layers (Vaswani al., 2017), after which the output is mapped the variational parameters v. Below, we describe this process in more detail.
Fi y: into sequence of vectors (b:) where First, we convert the observation sequence (4,
b: proj((.) + proj(æ;)) + proj(y:),
(17)
and proj are separate trainable linear projections. This additive decomposition is simple yet shows strong empirically performance in experiments. Since t; and æ; already provide positional information, we omit positional encodings. Additionally, we introduce an aggregation token, [AGG], with learnable representation to indicate where the encoder should aggregate information about the observations. This gives us the following encoder input sequence: {(b,) J [AGG]]. Next, stack of transformer encoder layers processes the input sequence, producing the output sequence [AGG]].
{1b, [AGGT] Finally, map the output aggregation token [AGG] to the variational parameters V.. This mapping depends on the form of the approximate posterior (see Appendix B for specification of the mapping).
5

---

Published as a conference paper at ICLR 2025
K01
Figure 3: Examples of trajectories from our datasets. White crosses indicate the observation locations marks on the horizontal time grid indicate the corresponding observation times f. and color-coded fields denote the system state with observations y at white crosses (t, 2:). Row Burger's (1D) dataset. The full trajectory is shown, with 1-D measurement location æ along the vertical and time along the horizontal directions. Rows 2, 3 and 4: Shallow Water, Navier-Stokes, and Scalar Flow datasets, respectively. Due to the 2D nature of the systems, only snapshots of the system state are plotted.
4.4 FORECASTING
Given set of initial observations (ti,0; y/ Ji=l of a test trajectory (we call it "context")), want to predict how the system behaves at some time point > fNex and spatial location æ. Given the context, predict the observation y at and æ the expectation of the following approximate posterior predictive distribution:
p(ylt, æ. (1i 27 N = p(y]t, 21) (z1))=(2))dz1
(18)
where 3/* = Encoder((t), a vE Nr), and the expectation is approximated via Monte Carlo integration. In Eq. 18, we omit explicitly conditioning on training data for brevity. Similarly, the distribution over the observation times and locations is evaluated as:
p(t. x [t], (Nrr
p(t(t pz)) (2))))))-())))2)) (21)dz1.
(19)
In our experiments, the context for each trajectory defined as all observations withing fixed initial time period (see Appendix for details).
EXPERIMENTS
In this section we demonstrate properties our method and compare against other methods from the literature. Our datasets are generated by three commonly-used PDE systems: Burgers' (models nonlinear ID wave propagation). Shallow Water (models 2D wave propagation under the gravity), and Navier-Stoke with transport (models the spread of a pollutant in a liquid over 2D domain). In addition to the synthetic data, we include a real-world dataset Scalar Flow (Eckert al., 2019), which contains observations of smoke plumes raising in warm air. See Appendix A for details about the dataset generation. Figure shows examples of trajectories from the datasets.
In all cases model has most 3 million parameters, and training takes most 1.5 hours single GeForce RTX 3080 GPU. The training is done for 25k iterations with learning rate 3e-4 and batch size 32. We use the adaptive ODE solver (dopri5) from torchdi f ffeq package with relative

---

Published as a conference paper at ICLR 2025
and absolute tolerance set to le-5. As performance metrics we mean absolute error (MAE) for the observations yi, and event-averaged log-liielhoood of the point process for t, and 2;, both evaluated on the test set. See Appendix C for details about training setup.
Context Size. Our encoder maps the initial observations (con-
Norter Stotees
Stalow
Scelzz 7low
text) {tz xf, yi /evx to parameters & of the approximate pos-
0.16
0.12
0.12
0.11
terior ((11 v). Here we look at how accuracy of the state
0.10
predictions and process likelihood is affected by the context size. We train and test model with different context sizes:
0.05
0.09
0.01
0.08
6.40
8.91
full (using all Nexx points), half of the context (using points from Natx/ 22 Netx), etc., and show the results in Figure 4.
5.58
8.55
4.75
0.20
3.93
7.85
We see that both MAE and process likelihood improve as we increase the context size, but the improvements tend saturate for larger contexts. This indicates that the encoder mostly uses observations that are close to the time point at which the latent state is inferred, and does not utilize observations that are too far away from it. This effect is especially visible with synthetic data, but the plots suggest the real-world dataset Scalar Flow could still benefit from larger context.
3.10
7.49
Nesuf2
Context Size
Figure 4: Test MAE (1) and log-
likelihood (t) vs. context size. For Scalar Flow we use the right axis for better visibility.
Strrpes Burgecs'
(nn)
Scalor Ctni
Stallow Water flia)
Scalar Cial
Dataset
Interp.
Seq.
0.12
0.10
6.00
8.91
Burgers
0.007
0.018
MAB
0.09
0.09
5.25
8.55
Shallow Water
0.011
0.045
0.07
0.08
4.50
8.20
Navier-Stokes
0.012
0.108
0.04
0.07
3.75
7.85
Scalar Flow
0.011
0.089
0.01
0.06
3.00
7.49
20
Grid Resolution (n)
Grid Resolution
Table 1: Latent state evalua- tion time (sec), interpolation sequential method. Aver- aged time over all trajectories.
Test MAE (L) and process log-likelihood (1) for dif- ferent temporal grid resolutions (n) and interpolation methods (nn = nearest neighbor, lin = linear). For Scalar Flow we use the right axis for better visibility.
Figure 5:
Latent State Interpolation. As discussed in Section 4. we do not evaluate the latent state at the full time grid f1....,t,, iN using the ODE solver (we call it the sequential method). Instead, we use adaptive ODE solver evaluate the latent state only at a sparser time grid 7i.-..y Tu, and then interpolate the evaluations to the full grid. Here we look at how this approach affects training times, and investigate the effects of the grid resolution n and different interpolation methods. We vary the resolution from coarse n = 2 to fine n 50, and test nearest neighbor and linear interpolation methods. Figure 5 shows that both MAE and process log- likelihhoood improve with the grid resolution, but only when we use nearest neighbor interpolation. With linear interpolation, the model achieves its best performance with just n 2, meaning that interpolation is done between only two points in the latent space, and further increasing the resolution does not improve the results. To study this observation in more detail, we provide additional experiments in Appendix D. Both interpolation methods achieve the same optimal performance when n > 20. This effect is visible for both synthetic and, to a smaller extent, real data. Furthermore, Table shows that our method results in up to 9x faster latent state evaluation than the sequential method, which in case translates in up to 4x faster training.
Latent Space Dimension. Our model assumes the system dynamics can be accurately modeled in a low-dimensional latent space. In this experiment we show that this assumption is valid at least for the four systems that we consider and best predictive accuracy is achieved for relatively oow-dimensional latent space dimension dz. In Figure we show how MAE and process log- likelihoooo improve increase the latent state dimension. For most datasets the predictive accuracy converges for as low as 16-dimensional latent space, except for the Navier-Stoker dataset, where d 64 is required, likely
7

---

Published as a conference paper at ICLR 2025
due to its more intricate state (as can be qualitatively observed in Fig. 3). We note that d2 affects the number of model parameters, but the difference between 1- and 64-dimensional latent space less than 3%.
Burgurs
Shallew
Navier Stokes
Scalar
0.15
8.50
6.00
STPP log-lik.
0.11
5.38
MAE
0.08
5.25
4.75
0.04
3.62
4.12
0.01
2.00
3.50
1
4
16
64
368
16
368
Latent Space Dim.
Latent Space Dim.
Figure 6: Test MAE (L) and process log-likelihood (f) for different latent space dimensions (dz). For Burgers' and Navier-Stokes we use the right axis for better visibility.
Interaction Between Observation and Process Models. In our model the observation times and locations are modeled by a point process model (Eq. 10), while observations are modeled by an observation model (Eq. 11). Intuitively, one could expect some interaction between them and ask the question of how important is one model for the performance of the other. In this experiment we look at the effect of removing the point process model on the performance of the observation model and vice versa. We remove the point process model by removing Eq. 10 from our model, so that we do not model the observation times and locations. Similarly, we remove the observation model by removing Eq. 11 from our model and also removing observations y; from our data, so that there is no information about the observation values. Table 2 we show that removing the point process model seems to have practically no effect on MAE, while removing the observation model decreases the process likelihood. This shows that if we know the system states we can model the observation times and locations more accurately.
MAE (-STPP)
Log-lik
Dataset
MAE
Log-lik.
(-Obs.)
Burgers'
0.042
0.043
3.77
3.70
Shallow Water
0.012
0.012
3.67
3.56
Navier-Stoker
0.071
0.071
5.96
5.92
Scalar Flow
0.090
0.089
8.43
6.87
Table 2: Test MAE (L) and log-likelihood (t) of the full model, and after removing the point process (-STPP) or the observation model (- Obs.).
Comparison to Other Methods. We compare our method against two groups of methods from the literature. The first group consists of neural spatiotemporal dynamic models that can model only the observations y: but not the time points and spatial locations. This group includes FEN (Lienen Günneeman, 2022) (graph neural network- based model, simulates the system dynamics every point on the spatial grid), NSPDE (Salvi et al., 2022) (space- and time-continuous model simulating the dynamics in the spectral domain), and CNN-ODE (our simple baseline that uses a CNN encoder/decoder to map observations to/from the latent space, and uses neural ODEs (Chen et al., 2018) to model the latent space dynamics). Note that due to their assumption about multiple observations per time point, dynamic models require time binning (see Appendix C.4). The second group consists of neural spatiotemporal point processes that model only the observation times and locations. This group consists of DSTPP (Yuan al., 2023), NSTPP (Chen al., 2021), and AutoSTPP (Zhou & Yu, 2023). We also provide simple baselines for both groups: median predictor for dynamic models, and constant (learnable) intensity for the point process group. Hyperparameters of each method were tuned for the best performance. See Appendix C.4 for details about the models and hyperparameters used.
Table 3 shows the comparison results. We see that most methods from the first group perform rather poorly and fail to beat even the MAE of the median predictor, with the CNN-ODE model showing the best results. However, CNN-ODE still performs considerably worse than model synthetic data, likely because of the loss of information caused by time binning. On the real-world Scalar Flow dataset CNN-ODE and model have similar performance. In the second group most methods
8

---

Published as a conference paper at ICLR 2025
Model
Burgers'
Shallow
Navier-Sioles
Scalar
MAE (1)
Log-lk. (1)
MAE (Q)
Log (1)
MAB (1)
Log. (1)
MAB (4)
Log-lis. (
Modian Prod.
0.161
0.148
0.155
0,140
FEN
NIA
0.536 0,008
0,158 1. 0.006
0.197 0.008
NSPDE
0.166 8,002
0.098 0.001
0.142 0.008
CNN- ODE
0,078 0.000
0.932 & 0.001
0.077 4: 0.002
0.091 0.001
Const, Insensity
3.59
2.02
5.86
6.87
DSTPP
0.67 1 0.13
0.08 0.01
0:41 0.01
3.02 = 0.01
NSTPP
2.35 N/A
0.13
0.03
AUDOSTPP
0.91
0.01
8.49 0.05
Ours
0.043 0.002
3.77 + 0.02 0.012 + 0.001
3.62 0.05
0.071 + 0.001
5.96 0.01
0.092
0.002
8:41 0.04
Table 3: Model comparisons. MAE (L) and Log-lik (per event) (1) on test data. The first group of methods (FEN, NSPDE, CNN-ODE) contains neural spatiotemporal dynamical models and are evaluated using MAE. The second group (DSTPP, NSTPP, AutoSTPP) contains neural spatiotemporal point process models and are evaluated using Log-lik. Error bars represent one standard error calculated over realizations of a random seed controlling the model initialization.
showed poor results and failed to beat the simple constant intensity baseline, however AutoSTPP shows very strong performance, although still not outperforming method on the synthetic data, and being close to on the Scalar Flow dataset.
These experiments show the challenging nature of modeling randomly observed dynamical systems, and demonstrate our method's properties and strong performance.
6 RELATED WORK
Traditionally, intensity functions simple parametric functions incor- porating domain knowledge about the system being modeled. While simple and interpretable, this approach might require strong domain knowledge and might have limited expressivity due overly simplistic form of the intensity function. These limitations lead to the development of neural point processes which parameterize the intensity function by a neural network, leading to improved flexi- bility and expressivity. For example Mei & Eisner (2017); Omi et al. (2019); Jia & Benson (2019); Zuo al. (2020) use neural networks to model dime-dependent intensity functions in a flexible and efficient manner, with neural architectures ranging from recurrent neural networks to transformers and neural ODEs. These techniques consistently outperform classical intensity parameterizations and do not require strong domain expertise to choose the right form of the intensity function. Other works extend this idea to marked point processes where intensity is a function of time and also of a discrete or continuous mark. For example, Chen et al. (2021); Zhu et al. (2021); Zhou et al. (2022); Zhou & Yu (2023); Yuan et al. (2023) use spatial coordinates as a mark and model the intensity using a wide range of methods aimed at improving flexibility and efficiency. Works such Du et al. (2016); Xiao et al. (2019); Boyd et al. (2020) use discrete observations as marks, and Sharma et al. (2018) use continuous marks with underlying dynamic model, which makes their work most similar ours, with major differences related to the model, training process, and them working only with temporal processes, whereas we are dealing with more general spatiotemporal processes.
Neural Point Processes.
Neural Spatiotemporal Dynamic Models. To model the spatiotemporal dynamics our method uses the "encode-process -recod"" approach employed in many other works. The main idea of the approach consists of taking a set of initial observations and using it to estimate the latent initial state. A dynamics function is then used to map the latent state forward in time to other time points either discretely (Long al., 2018: HAN et al., 2022; Wu et al., 2022) or continuously (Yildiz et al., 2019; Salvi et al., 2022). Finally, the latent state is mapped to the observations using discrete (Yildiz et al., 2019; Wu et al., 2022) or continuous (Yin et al., 2023; Chen et al., 2023) decoder. Such methods assume multiple observations per time point and do not model the observation times and locations. However, in contrary to many previous works that use the encode- encode-process decode approach to learn a discriminative model, we use the auto-encoding variational Bayes to learn generative model of the underlying dynamical system.
Implicit Neural Representations. Our model represents the continuous latent spatiotemporal state ult, a)) in terms of function ((ā(1), a:). This approach representing continuous field is called

---

Published as a conference paper at ICLR 2025
implicit neural representations (Park et al., 2019; Chen & Zhang, 2019; Mescheder et al., 2019; Chibane et al., 2020), and it has found applications in many other works related to modeling of spatiotemporal systems (Chen et al., 2023; Yin et al., 2023).
CONCLUSION
In this work, developed method for modeling spatiotemporal dynamical systems from marked point process observations. Along with our method, we proposed an interpolation-basee technique to greatly speed up the training time. In the experiments we showed that our method effectively utilizes the context of various lengths to improve accuracy of the initial state inference, uses the system state observations to better model the observation times and locations, and that low-resolution linear interpolation is sufficient to accurately represent the latent trajectories. We further demonstrated that our method achieves strong performance on challenging spatiotemporal datasets and outperforms other methods from the literature.
Limitations. Using a Poisson process for sampling event times and locations can be limiting in certain applications. First, the assumption of non-overlapping time points may hold in scenarios where effectively simultaneous observations occur, such as in high-frequency sensor networks. Second, the lack of interaction between the event occurrences and system dynamics may not hold if the act of observation influences the system (e.g., triggering human actions or affecting sensor availability). Finally, during data collection faithfully capturing the event intensity requires sufficiently dense sensor placement, especially in high-intensity regions. Future work could relax these assumptions to better model such scenarios, for example, by using different point process models and adapting the architecture to efficiently account for potential interactions between the observation events and system dynamics.
REFERENCES
Cristina Albaladejo, Pedro Sánchez, Andrés Iborra, Fulgencio Soto, Juan López, and Roque Torres. Wireless sensor networks for oceanographic monitoring: A systematic review. Sensors, 10(7): 6948-6968, 2010.
David M Blei, Alp Kucukelbir, and Jon D McAuliffe. Variational inference: A review for statisticians. Journal of the American statistical Association. 112(518):859 877, 2017.
Alex Boyd, Robert Bamler, Stephan Mandt, and Padhraic Smyth. User-dependent neural sequence models for continuous-time event data. Advances in Neural Information Processing Systems, 33: 21488-21499 2020.
Peter Yichen Chen, Jinxu Xiang, Dong Heon Cho, Yue Chang, G A Pershing, Henrique Teles Maia, Maurizio M Chiaramonte, Kevin Thomas Carlberg, and Eitan Grinspun. CROM: Continuous reduced-order modeling of PDEs using implicit neural representations. In The Eleventh Interna- tional Conference on Learning Representations, 2023. URL https: hutpss///oennsereeew.eeee
forum? 1.d-FUOR1tGG00s
Ricky T. Q Chen. torchdiffeq. 2018. torchdi ffeq.
URL https://ggithub..cm//ttt//
Ricky T. Q. Chen, Brandon Amos, and Maximilian Nickel. Neural spatio-temporal point processes. In International Conference on Learning Representations, 2021.
Ricky TQ Chen, Yulia Rubanova, Jesse Bettencourt, and David Duvenaud. Neural ordinary
differential equations. Advances in neural information processing systems, 31, 2018.
Zhiqin Chen and Hao Zhang. Learning implicit fields for generative shape modeling. Proceedings of the IEFF/CVF Conference on Computer Vision and Pattern Recognition, pp. 5939-5949 2019.
Julian Chibane, Thiemo Alldieck, and Gerard Pons-Moll. Implicit functions in feature space for 3d shape reconstruction and completion. In Proceedings of the IEEE/CVF conference on computer vision and pattern recognition, pp. 6970-6981, 2020.
10

---

Published as a conference paper at ICLR 2025
Earl A Coddington, Norman Levinson, and T Teichmann. Theory of ordinary differential equations, 1956.
Daryl Daley, David Vere-Jones, et al. An introduction to the theory of point processes: volume I: elementary theory and methods. Springer, 2003.
Nan Du, Hanjun Dai, Rakshit Trivedi, Utkarsh Upadhyay, Manuel Gomez- Rodriguez. and Le Song- Recurrent marked temporal point processes: Embedding event history to vector. In Proceedings of the 22nd ACM SIGKDD international conference on knowledge discovery and data mining, pp- 1555 1564, 2016.
Marie-Lena Eckert, Kiwon Um, and Nils Thuerey. Scalarflow: A large-scale volumetric data set of real-world scalar transport flows for computer animation and machine learning. ACM Transactions on Graphics, 38(6):239, 2019.
Nicholas Geneva and Nicholas Zabaras. Modeling the dynamics of pde systems with physics- constrained deep auto-regressive networks. Journal of Computational Physics,
2020.
Moustafa Ghanem, Yike Guo, John Hassard, Michelle Osmond, and Mark Richards. Sensor grids for air pollution monitoring. In Proc. 3rd UK e-Science All Hands Meeting, Nottingham, UK, 2004.
Pablo Gómez, Gabriele Meoni, and Håvard Hem Toftevaag. torchquad, 2024. URL https:
//github. toom/eesattet essatorchquudd
XU HAN, Han Gao, Tobias Pfaff, Jian-Xun Wang, and Liping Liu. Predicting physics in mesh- reduced space with temporal attention. In International Conference on Learning Representations, 2022. URL https:tpp//oo openooveeevveeeee ew net/ fo rum? ? LdNNffCC
E Heirer, SP Nørsett, and G Wanner. Solving ordinary differential equations i: Nonstiff problems, 1987.
Dan Hendrycks and Kevin Gimpel. Gaussian error linear units (gelus), 2016.
Philipp Holl, Vladlen Koltun, Kiwon Um, and Nils Thuerey. phiflow: A differentiable pde solving framework for deep learning via physical simulations. In NeurIPS Workshop on Dif- ferentiable vision, graphics, and physics applied to machine learning, 2020. URL https: //montreal robot ics. .ww.c.cca.c..oo................. / di wfcvge/ a /ttpssssstsseepsseeeeeesssssssss 3.pdf.
Valerii lakovlev, Markus Heinonen, and Harri Lähdesmäki. Learning continuous-tima {pde)s from sparse data with graph neural networks. In International Conference on Learning Representations, 2021. URL https: www.reveeevveeeeeeeee
net www.m/orrmmmmmo fow.fommmmmmuuaaa 1a
Junteng Jia and Austin R Benson. Neural jump stochastic differential equations. Advances in Neural Information Processing Systems, 32, 2019.
Diederik P Kingma and Max Welling. Auto-encoding variational bayes. arXiv preprint
arXiv:1312.61 14,6111111010 2013.
Qingkai Kong, Richard M Allen, Louis Schreier, and Young- Woo Kwon. Myshake: A smartphone
Woong-Woo seismic network for earthquake early warning and beyond. Science advances, 2(2)se1501055, 2016.
PA W Lewis and Gerald S Shedler. Simulation of nonhomogeneous poisson processes by thinning. Naval research logistics quarterly, 96(3):403-413 1979.
Marten Lienen and Stephan Günnemann. Learning the dynamics of physical systems from sparse ob- servations with finite element networks. In International Conference on Learning Representations, 2022. URL https: : / opw.reveewe..eee wwwpereveww.eee/forrrr.. net/forum wtt. -d-HHmuuk----222
Zichao Long, Yiping Lu, Xianzhong Ma, and Bin Dong. Pde-net: Learning pdes from data. In International conference machine learning, pp. 3208-3216 PMLR, 2018.
Ilya Loshchilov and Frank Hutter. Decoupled weight decay regularization. In International Confer- ence on Learning Representations, 2019. URL https:/ http///pppe//eiieee//e/ee net forum?id= Bkg6Ricqq) Cqy7.
11

---

Published as a conference paper at ICLR 2025
Yajie Ma. Mark Richards, Moustafa Ghanem, Yike Guo, and John Hassard. Air pollution monitoring and mining based on sensor grid in london. Sensors, 8(6):3601-300 -3623, 2008.
Mihai Marin-Perianu, Supriyo Chatterjea, Raluca Marin- Stephan Bosch, Stefan Dulman, Stuart Kininmonth, and Paul Havinga. Wave monitoring with wireless sensor networks. In 2008 International Conference on Intelligent Sensors, Sensor Networks and Information Processing, pp. 611-616. IEFE, 2008
Hongyuan Mei and Jason M Eisner. multivariate point process. Advances in neural information processing systems, 30, 2017.
The neural hawkes process: A neurally self-modulating
Lars Mescheder, Michael Oechsle, Michael Niemeyer, Sebastian Nowozin, and Andreas Geiger. Occupancy networks: Learning 3d reconstruction function space. In Proceedings of the IEEE/CVF conference computer vision and pattern recognition, pp. 4460-4470, 2019.
Sarah Minson, Benjamin A Brooks, Craig L Glennie, Jessica R Murray, John O Langbein, Susan Owen, Thomas H Heaton, Robert A lannucci, and Darren L Hauser. Crowdsourced earthquake early warning. Science advances, 13) 1500336, 2015.
Takahiro Omi, Kazuyuki Aihara, al. Fully neural network based model for general temporal point processes. Advances in neural information processing systems, 32, 2019.
Jeong Joon Park, Peter Florence, Julian Straub, Richard Newcombe, and Steven Lovegrove. Deepsdf: Learning continuous signed distance functions for shape representation. In Proceedings of the IEEE/CVF conference computer vision and pattern recognition, pp. 165-174, 2019.
Tobias Pfaff, Meire Fortunato, Alvaro Sonzaez-Gonzalez, and Peter Battaglia. Learning mesh-based simulation with graph networks. In International Conference Learning Representations, 2021.
URL https:://///mm /tstt:::/oppeee wwpp/ppeeneeeeeeeeeeeeeeeeeeeeee net forum? iediiidooo redd=ooNNYOL_. XP.
Christopher Rackauckas, Yingbo Ma. Julius Martensen, Collin Warner, Kirill Zubov, Rohit Supekar, Dominic Skinner, Ali Ramadhan, and Alan Edelman. Universal differential equations for scientific machine learning. August 2020. doi: 1.11203//....s55555555... URL http: /tx. doi org/ 10.21203/RS3 RSS 3. R--RRS55155 55125/V1.
Cristopher Salvi, Maud Lemercier, and Andris Gerasimovics. Neural stochastic pdes: Resolution- invariant learning of continuous spatiotemporal dynamics. Advances in Neural Information Processing Systems, 15:1333-1344 2022.
Anuj Sharma, Robert Johnson, Florian Engert, and Scott Linderman. Point process latent variable models of larval zebrafish behavior. Advances in Neural Information Processing Systems, 31, 2018.
Makoto Takamoto, Timothy Praditia, Raphael Leiteritz, Dan MacKinlay, Francesco Alesiani, Dirk Pflüger, and Mathias Niepert. PDEBENCH: AN EXTENSIVE BENCHMARK FOR SCI- EN- TIFIC MACHINE LEARNING. ICLR 2023 Workshop on Physics for Machine Learning, 2023.
URL htttt::///peereeee//
net forum? 1d=b85wox202k1,
Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan Gomez, Lukasz Kaiser, and Illia Polosukhin. Attention is all you need. Advances in neural information processing systems, 30, 2017.
Tailin Wu, Takashi Maruyama, and Jure Leskovee. Learning to accelerate partial differential equations via latent global evolution. In Alice H. Oh, Alekh Agarwal, Danielle Belgrave, and Kyunghyun Cho (eds.), Advances in Neural Information Processing Systems, 2022. URL https: /epenreviivveeeeeeeeeeeeeeeeevee ew net forrmmid==x7ttgpee
Shuai Xiao, Junchi Yan, Mehrdad Farajtabar, Le Song, Xiaokang Yang, and Hongyuan Zha. Learning time series associated event sequences with recurrent point process networks. IEEE transactions on neural networks and learning systems, 30( 10):312433336, 2019.
Guobao Xu, Weiming Shen, and Xianbin Wang. Applications of wireless sensor networks in marine environment monitoring: A survey. Sensors, 14(9):1693--1655. 2014.
12

---

Published as a conference paper at ICLR 2025
Cagatay Yildiz, Markus Heinonen, and Harri Lahdesmaki. Ode2vae: Deep generative second order odes with bayesian neural networks. Advances in Neural Information Processing Systems, 32, 2019.
Yuan Yin, Matthieu Kirchmeyer, lean-Yves Franceschi, Alain Rakotomamoniy, and patrick gallinari. Continuous PDE dynamics forecasting with implicit neural representations. In The Eleventh International Conference on Learning Representations, 2023, URL httpps//openzeeien.
Yuan Yuan, Jingtao Ding, Chenyang Shao, Depeng Jin, and Yong Li. Spatio-temporal diffusion point processes. In Proceedings of the 29th ACM SIGKDD Conference on Knowledge Discovery and Data Mining, KDD '23. ACM, August 2023, doi: 10.1 (0.11455358030053399555 URL http:
//dx. doi. org/ 10. 11 145535800055.... 3599511.
Zihao Zhou and Rose Yu. Automatic integration for spatiotemporal neural point processes, In Thirty-seventh Conference Neural Information Processing Systems, 2023. URL https:
/ //pppe/eeeew.neeeeeeeeee net. prrpen/wiiee.nttttoooooeeeeee-pe
Zihao Zhou, Xingyi Yang, Ryan Rossi, Handong Zhao, and Rose Yu. Neural point process for learning spatiotemporal event dynamics. Learning for Dynamics and Control Conference, pp. 777-789. PMLR, 2022.
Shixiang Zhu, Shuang Li, Zhigang Peng, and Yao Xie. Imitation learning of neural spatio-temporal point processes. IEEE Transactions on Knowledge and Data Engineering, 34(11):5391-544 2021.
Simiao Zuo, Haoming Jiang, Zichong Li, Tuo Zhao, and Hongyuan Zha. Transformer hawkes process. In International conference on machine learning, pp. 1 11692-1 11092-11700 11702. PMLR, 2020.
13

---

0, Published as a conference paper at ICLR 2025
A DATA
As described in Section 3, we assume the following data generating process
@u(t,
u(t1,-) ~ p(u), 4, 0; ~ ulpp(((u(t;x)))) y: ((ull((((tt,
(20)
ôt
1 N,
(21)
N,
(22)
We selected three commonly used PDE systems: Burgers', Shallow Water, and Navier-Stokes with transport. Next we discuss the data generating process for each dataset.
A.1 BURGERS'
The Burgers' system in ID is characterized by I-dimensional state u(l,2) a) and the following PDE dynamics:
Ou Ou
8uu 01."
1400
(23)
01
dr
where v is the diffusion coefficient. We obtain data for this system from the PDEBench dataset (Takamoto al., 2023). The system
600
400
simulated on a time interval [0 and on a spatial domain [0,
with periodic boundary conditions. The dataset contains 10000 trajectories, each trajectory is evaluated on a uniform temporal grid with 101 time points, and uniform spatial grid with 256 nodes.
100
2010
4041
Figure 7: Trajectory lengths
for Burgers' dataset.
Since the data is discrete, we define continuous field u(t, x) by linearly interpolating the data across space and time. Next, we use the interpolant define the intensity function
Aut(,, .)) (t,))u(t,
(24)
and use it to simulate the non-homogeneous Poisson process using the thinning algorithm (Lewis & Shedler, 1979). the result, for each the 10000 trajectories, we obtain sequence of time points and spatial locations tt.,2;44 J where N can be different for each trajectory. Finally, compute the observations as
y: = attt,
(25)
We use the first 0.5 seconds as the context for the initial state inference, and further filter the dataset to ensure that each trajectory has at least 10 points in the context, resulting in approximately 7000 trajectories. We use 80%\/1%%/10% split for traniingvalidatiion//singg
A.2 SHALLOW WATER
The Shallow Water system in 2D is characterized by 3-dimensional state
Tojj. leng bs (Shallow Welser
u(t, x) [h(t, x), att, x), (t)
(26)
and the following PDE dynamics:
0h @hhu)
Öhvv)
0t +
+
= 0,
(27)
dx
8y
du
Ou
Ou
ah 0,
Ot +
u
dr
Dy
Or
(28) Figure 8:
600 500
Trajectory lengths
800
for Shallow Water dataset.
dv dv
U du
@h
9 Oy
(29)
0t
Ox
Oy
where is the gravitational acceleration, hit, x) is wave height, and u(t, x) and ett, x) are horizontal and vertical velocities, respectively. We obtain data for this system from the PDEBench dataset
14

---

Published as a conference paper at ICLR 2025
(Takamoto al., 2023). The system is simulated time interval [o. 1] and spatial domain [-2.5, 2.5]². The dataset contains 1000 trajectories, each trajectory is evaluated on a uniform temporal grid with 101 time points, and uniform 128-by-125 spatial grid.
Since the data is discrete, we define continuous field u(t, z) by linearly interpolating the data across space and time. Next, we use the interpolant define the intensity function
((u(1, æ)) = (t, x) 1/,
(30)
indicating that measurements are made only when the wave height below above the baseline of 1. We use the intensity function to simulate the non-homogeneous Poisson process using the thinning algorithm (Lewis Shedler, 1979). As the result, for each of the 1000 trajectories, we obtain sequence of time points and spatial locations {&,, x; ";) where N can be different for each trajectory. Finally, we compute the observations as
y = h(4,
(31)
observing only the wave height.
We use the first 0.1 seconds as the context for the initial state inference, and further filter the dataset to ensure that each trajectory has at least 10 points in the context, resulting in approximately 800 trajectories. We use 80%/10%/100% split for rraining/validation//ssinng
A.3 NAVIER-STOKES
The Savier-Stokes system with a transport equation is characterized
by a 4-dimensional state
Tra. lesgths (Navter Staues)
u(l, x) [o(t, æ), u(t, æ), o(1,,),,
(32)
and the following PDE dynamics:
153
Du
dv
0,
(33)
Fx
Oy
50
Ou + u ax
du
Ou
Op + v
82u
v
1250 1500 1750 2000 2250 2000 2750 (34) Figure 9: Trajectory lengths for Navier-Stokes dataset.
Ot
dy
dr
2x
oy²
du
dv
dv
Op
02
10 + u 8x
+ U
+ C,
(35)
Du
Oy
Ou2
01²2
Oc = dr
dc
dc
0-c
0c
σt
v Dy
+
d.2
0y²
(36)
where the diffusion constant v = 0.002, and c(t, x) is concentration of the transported species, p(t, 2) is pressure, and u(t, a) and v(t, w) are horizontal and vertical velocities, respectively. For each trajectory, we start with zero initial velocities and pressure, and the initial scalar field eo (I. y) is generated as:
do(x,
Σ
Aut (occos22TT20) Iy)) Th sin Iy)).
(37)
k.1=-N
(0 y) - min(Co)
co(x, v) max ((o) min(Cn)
(38)
where = 2 and λaa 7kk ~ (0.
We use PhiFlow (Holl et al., 2020) to solve the PDEs. The system is simulated on a time interval [0, 1] and on a spatial domain [0, 1]2 with periodic boundary conditions. The solution is evaluated at randomly selected spatial locations and time points. We use 1089 spatial locations and 25 time points. The spatial and temporal grids are the same for all trajectories. The dataset contains 1000 trajectories.
Since the data is discrete, define continuous field u(t, 2) by linearly interpolating the data across space and time. Next, we use the interpolant to define the intensity function
A(u(t, 2000)c(t, c(t,
(39)
15

---

Published as a conference paper at ICLR 2025
indicating that measurement intensity is proportional to the species concentration. We the intensity function to simulate the non-homogeneous Poisson process using the thinning algorithm (Lewis & Shedler. 1979). As the result, for each of the trajectories, spatial locations {4, 24 ) where N can be different for each trajectory. Finally, we compute the observations
obtain sequence of time points and
y: = c(t;,
(40)
observing only the species concentrations.
We use the first 0.5 seconds as the context for the initial state inference, and further filter the dataset to ensure that each trajectory has at least 10 points in the context, resulting in all 1000 trajectories satisfying this condition. We use 80%/10%/10% split for training/validatioonttsnnngg
A.4 SCALAR FLOW
To demonstrate capabilities of our model on a real-world process, we use the Scalar Flow Dataset (Eckert et al., 2019). This dataset consists of videos of smoke plumes rising in hot air. Each video has 150 frames recorded over 2.5 seconds, and each frame has the resolution of 1080x 1920 pixels. The observations are postprocessed camera images of the smoke plumes taken from multiple views. For simplicity, we use only the front view. We downsample the original frames from 1080x1922 to 120x213 pixels and smooth out the high-frequency noise by applying Gaussian smoothing with the standard deviation of 1.5.
Then, we emulate real-world measurement process where each "sensor" (pixel) makes mea- surement with the probability directly proportional to the observed smoke density (higher den- sity implies more frequent measurements). In particular we set the observation probability to 0.03 x smoke density, resulting in 2000 observations per trajectory on average. As a result, we have a sequence of smoke density observations made at different time points and spatial locations. Note that we do not make any Poisson process assumptions here. Instead, we use realistic measurement process where sensors are utilized only when some action in the measured process stars to happen, which is the application
targeting in our work.
We use the first half (1.25 sec.) of each trajectory as the context for the initial state inference. In total, the dataset contains 100 trajectories. We train our and other models 80 trajectories, and use 10 trajectories for validation, and 10 for testing.
B
MODEL, POSTERIOR, AND ELBO
B. MODEL
As described in Section 4.1 our model defined as
Zi p(21), dz(r)
[(z(T)),
(41)
dt
(t) interpolate (t)e(t);(()).....)) (l),x) a))
(42)
(43)
4, .. wwwppp((p(lltlla
1...., N,
(44)
~ plyllllguu æ:))),
N.
(45)
The corresponding joint distribution is
p(21, (4i, 2,, yi = 2
II p(y P xi),
(46)
i=1
where
p(z1) = =N(O
(47)
p((1,, az 1211 II X ((u(t, m:) exp (-/ x(u(t, (t)(tt
(48)
p (,, = (g xx)), o2
(49)
where A is the normal distribution, 0 is vector, and / is the identity matrix. We set by 10-3
16

---

Published as a conference paper at ICLR 2025
B.2 POSTERIOR
We define the approximate posterior with vi - [W [oa] as
q( v/) = (Z 1/ diag(No=)),
(50)
where N is the normal distribution, and d is matrix with vector ب) on the diagonal.
As discussed in Section 4.3, the encoder maps the context {4,, c, ye Nexx to the local variational parameters & (which we break up into [4/,. 2oa)). The encoder first converts the context sequence {t,, v.. yt the following input sequence: { Je [AGG]]. We pass it though the transformer encoder and read the value of the aggregation token at the last layer, which we denote by [AGG]. Then, we map [AGG] [e) th as follows:
N.a ) to a sequence of vectors {b, Noxx and then adds an aggregation token, which gives
1/ Linear([AGG]),
(51)
2 exp [Linear[[AGG]]
(52)
where Linear separate linear mappings.
B.3 ELBO
Given definitions in the previous sections, can write the ELBO
/ ((21) 1(21, (4,, (0, Vila dz
(53)
((21)
-f ((2))) In I1 p(y: hi ) x;
/ 1 ((21))
dz1
(54)
(2)1)
(())) II p(y 21, t, xz, dzz
(55)
/ ((21) In In x; ) 1 (2))
(56)
+ / ((z1) In p(z1) dz1
(57)
((21)
N
By(x.) [ln |21 m7))
(58)
In X(u((,, 2x)) f X(ult, dxodt
E Eg(ss)
(59)
C2
KL(g(z1) [|
(60)
La
C1 + Lg -c- Lz.
(61)
Computing ELBO. We compute the ELBO as follows:
1. Compute local variational parameter from the context as
Vu 2 Encoder((t, x,, y: Nax)
2. Sample the latent initial state as 21 ~ N(21 W/, diag(1ba)) using reparameterization
3. Evaluate 2(T)...... (TT)) by solving dz(r) solve the ODE using torchdi f feq package (Chen et al., 2018).
= ff (x((T) with ø; as the initial condition. We
dt
17

---

Published as a conference paper at ICLR 2025
4. Compute the latent states ((t)) ...f),................ z(tx) via interpolation
(t) = interpolate x (TT) (T))) (Tn))
5. Compute u(f1, æ1). .... ultN NN) as ((t;, x.) (tt(t)..
6. Compute L1 using Monte Carlo integration
7. Compute L2 using Monte Carlo integration, with the intensity integral computed using torchquadi (Gómez al., 2024) package.
8. Compute L3 analytically
Sampling is done using reparametrization (Kingma & Welling, 2013) to allow for exact gradient evaluation via backpropagation. Monte Carlo integration done using sample size of one.
C EXPERIMENTS SETUP
C.I DATASETS
See Appendix A for all details about dataset generation. For all datasets we use 80%/10% 10% vainn/vlidattonntestt splits. The context size Nex is determined for each trajectory separately based on the time interval that we consider to be the context, for Burgers' and Navier-Stokes it is the first 0.5 seconds (out of 2 seconds), while for Shallow Water it is 0. seconds (out of second). Since the number of events occurring within this fixed time interval is different for each trajectory, the context size consequently differs across trajectories. See encoder architecture description below for details on how different context lengths are processed.
C.2 TRAINING, VALIDATION, AND TESTING
We use Adam W (Loshchilov & Hutter, 2019) optimizer with constant learning rate 3e-4 (we use linear warmup for first 250 iterations). We train for 25000 iterations, sampling a random minibatch from the training dataset at every iteration. We do not any kind of scaling for the ELBO terms.
We compute validation error (MAE on observations) on single random minibatch from the validation set every 125 iterations. We save the current model weights if the validation error averaged over last 10 validation runs is the smallest so far.
To simulate the model's 's dynamics (Chen et al., 2018). In particular, we use the dopri5 solver with rtol = atol = 10 without the adjoint method. To evaluate the intensity integral f >(u(t,, ax) dxdd
differentiable ODE solvers from torchdi ffeq package
Monte Carlo integration from the torchquad (Gómez et al., 2024) package with 32 randomly sampled points for training, and 256 randomly sampled points for testing, found these sample sizes be sufficient for producing small variance of the integral estimation in our case.
C.3 ENCODER AND MODEL ARCIITECTURES
Encoder: We use padding to efficiently handle input contexts of varying length. The time points, coordinates, and observations are linearly projected to the embedding space (128-dim for Burgers' and Shallow Water, and 192-dim for Navier-Stokes). Then, the projections are summed and passed through stack of transformer encoder layers (4 layers for Burgers' and Shallow Water, and 5 layers for Favier-Stokes. 4 attention heads were used all cases). We introduce a learnable aggregation token which we append at the end of each input sequence. The final layer's output is read from the last token of the output sequence (corresponds to the aggregation token) and is mapped to the local variational parameters as discussed in Appendix B.
Model. The dynamics function fx((t)) is an MLP with 3 layers 368 neurons each, and GeLU (Hendrycks & Gimpel, 2016) nonlincarities. We set resolution of the uniform temporal grid .i...., Ty to n = 50. The dynamics function maps the current latent state xt) to its time derivative dt dy-dimensional vectors.
dz(t) both
The function ((t(t,, x)) is an MLP with hidden layers with width 368 for Burgers', 256 for Shallow Water, and 512 for Navier-Stokes. We GeLU (Hendrycks & Gimpel, 2016) nonlinearities. As the
18

---

Published as a conference paper at ICLR 2025
input the MLP takes the sum of (t) E Rd with proj(z) E R where proj is a linear projection, and maps the sum to the output u(t, x) E R
The intensity function A(u(t, æ)) is an MLP with 3 layers, 256 neurons each, and GeLU (Hendrycks & Gimpel, 2016) nonlincarities. We further exponentiate the output of the MLP and add a small constant (1e-4) to it to ensure the intensity is positive and to avoid numerical instabilities.
We set the observation function g(u(t,z)) to the identity function, thus g x)) = =u(t x)).
C.4 MODEL. COMPARISON
C.4.1 DYNAMIC SPATIOTEMPORA MODELS
In this section we discuss neural spatiotemporal models we used for the comparisons. But first, we discuss data preprocessing that we do in order to make these models applicable.
Data Preprocessing. All models in this section assume that dense observations are available at every time point, and some further assume the observations are located on a uniform spatial grid. To make our data satisfy these requirements we use time binning and spatial interpolation. First, we use time binning and divide the original time grid into 10 bins and group all time points into 10 groups based on which bin they belong to. We used 10 time bins as it ensured that there was a sufficient number of time points in each bin while also being sufficient to capture changes of the system state. We denote by 1 the time point that represents the bin Next, we use all observations that belong to bin and interpolate them to uniform 32x32 spatial grid. In many cases some points of the spatial grid were outside of the convex hull of the interpolation points, so we set values of such points -1 (values of the state range from 0 to 1). We used linear interpolation. We denote the resulting interpolated values at time tj as u, which is the interpolated system state at the 32x32 spatial grid. When we compute the loss between uu and the corresponding model prediction use only those spatial locations that had proper interpolation values (i.e. were inside the convex hull of the interpolation points).
This is a simple baseline where we compute the median of the observations y; on the training data and use it the prediction on the test data.
Median predictor.
CNN-ODE. This is another relatively simple baseline based on latent neural ODEs (Chen et al., 2018). It operates in some sense similarly our model it takes sequence of initial observations (context) .i,..., UN... and uses a CNN encoder to map the context to the latent initial state 21- Then maps 21 to 22..... EN via a latent ODE, and finally decodes z; to obtain u, via CNN decoder. In our case we used the first two observations u and U2 as the context.
FEN. This model discretizes the spatial domain into mesh and uses graph neural network to model the dynamics of the state at each node of the mesh. The model starts at u, and maps it directly u2, , uN via an ODE solver without using any context. We simulated the dynamics as stationary and autonomous, with the free-form term active. The dynamics function is a 2-layer MLP with the width of 512 neurons and ReLU nonlinearities. All other hyper-parameters were left as defaults as we found that changing them did not improve the results. We used the official code from Lienen & Günnemann (2022).
This model applies point- wiie transformation to map the system state at every spatial node to a sigh--imensional latent representation, then uses the fixed-point method to simulate the dynamics at every node, and finally maps the nigh-dimensional latent states back to the observation space to obtain the predictions. We used the fxed-point method with maximum number of available modes and single iteration as we found that larger number of iterations does not improve the results. All other hyper-parameters were left as defaults as we found that changing them did not improve the results. We used the official code from Salvi al. (2022).
NSPDE.
C.4.2 NEURAL SPATIOTEMPORAN PROCESSES
Constant Intensity Baseline. As a simple baseline we use the following intensity function X(t, x) = where is a learnable parameter that we fit on the training data.
19

---

Published as a conference paper at ICLR 2025
hy
dz
his
23
92
368
512
40
128
0,090
0.091
0.084
0.081
40
512
0.075
0.073
0.074
0.072
368
128
0.084
0.083
0.084
0.082
368
512
0.072
0.073
0.073
0.072
Table 4: Dependence of the test MAE on latent space dimension dz, dynamics complexity b and decoder complexity hg.
DSTPP.
This model uses a transformer-based encoder to condition a diffusion-based density model on the event history. We use 500 time steps and 500 sampling steps, and train for at most 12 hours. All other hyper- pprammerrseters were left as defaults as we found that changing them did not improve the results or made the training too slow. We used the official code from Yuan et al. (2023).
NSTPP. This model represents the vent history in terms of latent state governed by neural ODE, with the latent state being discontinuously updated after each event. The latent state is further used to represent the distribution of spatial locations via normalizing flow model. We used the attentive CNF version of the model. All other hyper-parameters were left as defaults as we found that changing them did not improve the results. We used the official code from Chen et al. (2021).
AutoSTPP.
This model uses dual network approach for flexible and efficient modeling of spatiotem- poral processes. We used 5 product networks with 2 layers 128 neurons each and hyperbolic tangent nonlincarities. The step size was set to 20 and number of steps was set to 101 in each direction. All other cryper-parameter were left as defaults as we found that changing them did not improve the results. We used the official code from Zhou & Yu (2023).
D EFFECTS OF DYNAMICS FUNCTION AND DECODER COMPLEXITY
Our model (similarly to previous competing models) consists of a composition of mappings that are parameterized by deep neural networks. Therefore, the role and importance of different model components depends on a number of hyperparameters. This, in turn, is related to the result reported in the main text (Figure 5), where we observed that, in the case of linear interpolation, decreasing the number of interpolation points did not decrease predictive performance of our model. Decreasing the number of interpolation points to n = 2 corresponds to replacing non-linear dynamics function with linear model, but this decrease in complexity of the dynamics model can, in some cases, be compensated by the capacity of decoder function. In this section, we conduct a series of experiments to better understand the role of the dynamics and decoder functions in our model. In particular, look at the model's predictive accuracy as a function of the latent space dimension dz. and complexity of the dynamics function ʃ and decoder cb. Since both / and Φb are 3-layer MIPs, we define their complexities as the width of the hidden layers hy and hg. respectively. In particular, we look at low- and dibenniinniii latent spaces (dz 40 and d2 = 368); we consider medium and large decoder function o (with hg = 128 and h 512), and dynamics functions of different complexities ranging from h 23 to hy 512. We use linear interpolation with 50 interpolation points to allow the model represent complex latent trajectories if it helps to fit the data better. We use the -tokes-Stokes data (since this is the most complex synthetic dataset in our work) and report only the mean absolute error (as log- likelihoo follows the same pattern). The results are in Table 4.
As can be seen, for the low-dimensional latent space (dz 40) and redium--sized decoder (with h 128) increasing complexity of the dynamics function h, Ieads to a meaningful reduction of MAE. However, after increasing complexity of the decoder from It = 128 to he = 512, increasing h results in significant improvements. Then, if we switch to high-dimensional latent space (d2 = 368), we see that complexity of the decoder he becomes the main determinant of the model's s predictive accuracy, with dynamics function complexity h, having observable effect.
In summary, we see that increasing complexity of the dynamics function of the neural ODE can be useful in some cases dow-diineli latent space and medium-sized decoder), but with a sufficiently
20

---

Published as a conference paper at ICLR 2025
large decoder and large latent spaces its benefits become less apparent, meaning that in such cases the dynamics can be modeled using simple function approximators. While the experimental results show that in some cases accurate predictions can be achieved with simple dynamies, enabling the model to represent complex dynamics via Neural ODEs can be useful also for cases where complex latent trajectories expected or enforced.
21


Reviews:
[{"review_id": "LJp7m8PvU0", "reviewer": "Reviewer_UEak", "summary": "This is an engineering oriented work that model STPP with intensity driven by a continuous latent states governed by a Neural-ODE, with initial states generated by a transformer encoder. The formulation sounds valid and the proposed benefits are for sparsely distributed data. The main contributions are the new formulation and the interpolation-based speedup technique.", "questions": "* Can you give me a possible explanation of why it works?\n* There is no ablation study on why transformer are used for generating the initial states. Or do you have evidence the initial state part is robust to architecture choice?\n* Despite the proposed speedup method, I believe neural-ODE is still untolerably slow and does not scale well. Do you have actual training/ inference time comparison?", "limitations": "", "rating": 8, "confidence": 4, "soundness": 2, "presentation": 4, "contribution": 3, "strengths": "- The challenge seems well grounded as the sparse data over a large spatial domain is common for many types of problems, e.g., few agent trajectories over a large geographical domain.\n- The method looks (possibly) scalable with low-resolution linear interpolation.\n- The math formulation is clear and the empirical results are fair. \n- A lot of ablation study, accounting for context size, resolution, removal of components", "weaknesses": "* I don't think the paper really answer the question of why it work on sparse data. There is no theoretical analysis / visualization of how the low-dimensional latent space captures the full dynamics from sparse observations. No discussion of information-theoretic bounds on what can be learned from sparse observations. It is reasonable to expect normalizing-flow based method (like Neural-STPP) not working well because the distribution is too localized, but I don't see why your method have an theoretical advantage over SOTA with kernel-based or closed spatial distribution.", "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": {}}, {"review_id": "FPi2ZQmSnB", "reviewer": "Reviewer_f2Uy", "summary": "This article is devoted to the problem of learning spatiotemporal dynamics from randomly sampled points in space and time. This problem is particularly well suited for the situation where we have sensors that record a system, and we have to predict also the behavior of the sensors during the dynamics (e.g. meteorological sensors that are carried by currents). The method proposed in this article is based on the use of neural ODEs in a learned latent space.", "questions": "My main questions relate the points raised in the weaknesses above.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 3, "strengths": "-- The problem of learning spatiotemporal point processes is rather important, and any contribution to this problem should be well welcomed by the scientific community.\n-- The overall idea of the article is meaningful. \n-- The numerical results are rather good.", "weaknesses": "-- Some explanations are not properly given. For instance, I assume that they are using an ODE solver in a latent space because a direct approach would immediately incur into stiffness problems. Why not using a neural PDE solver? Why is it better to learn a latent space and use an ODE solver for a problem that is formulated as a PDE (as in Eqn 5 of the paper)? This is unclear.\n-- The latent approach makes the approach less clear, and more out of the control of the user. I suppose the authors have no idea why the encoder creates a certain latent space rather than another. A theoretical approach seems very complicated, in fact the authors limit themselves mostly to empirical results. \n-- It is unclear if a general system can be learned in this way. In a sense, we might think of the encoded latent space as a low-degree approximation of the system, but it might be that certain PDE models stemming from Eqn 5 might not be suitably tackled by such approach. \n-- One of the main claims is that the model is continuous. An interpolation task should be performed in this case to show that they can handle continuity well. They use interpolation in the method, but it is unclear if in an experiment where portions of the trajectories are completely hidden during training, could be recovered during evaluation.", "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": {}}, {"review_id": "L9F8XMRSM9", "reviewer": "Reviewer_BCH2", "summary": "A composition of different ML methods is presented to simulate spatiotemporal processes from point observations without access to the whole spatial field. The proposed approach encodes sparse context point observations into a latent state. In this latent state, the dynamic process evolution is integrated over large time steps, while a fine temporal resolution is obtained via interpolation in the latent state. A decoder projects the predicted latent states back to the observation space.", "questions": "1. For what reason is the distribution of the next sensor signal's location predicted? What is the benefit of such a prediction and what computational cost does it impose? If I understand correctly, Table 2 suggests removing the point process model (which simulates the next sensor signal position and time, if I'm correct). At least according to a minimal model when following Occams Razor.\n2. The interpolation ablation is very illustrative. Have you tried higher-order interpolations to infer $\\hat{z}(t_i)$, i.e., quadratic, cubic? What is the error that incurs from the interpolation compared to modeling the full temporal grid $t_1, \\dots, t_N$? Table 1 demonstrates the time improvement when using interpolations; it would be great to see the error associated with the two techniques (Interp. vs Seq.).\n3. Have you explored other solvers beyond dopri5, such as Euler, which is much cheaper? Or does the method depend on an adaptive solver to account for potentially different deltas between time steps? Figure 2 somehow suggest that the effectively processed time steps $\\tau_m$ are separated by a constant time delta. Is this a requirement of the architecture?\n4. How does the latent space dimensionality $d_z$ affect the runtime? Might be interesting to report along with its effect on the parameter count around line 375.", "limitations": "", "rating": 8, "confidence": 4, "soundness": 4, "presentation": 3, "contribution": 3, "strengths": "_Originality:_ The combination of various different branches from machine learning is original. They are composed in an elegant and versatile way to solve spatiotemporal problems efficiently. The use of latent states enforces abstractions that push for generalizability. Intriguingly, the introduced method does not rely on complex architectures, but mostly emerges out of various MLPs that are well placed and wired.\n\n_Quality:_ Claims are well supported with experimental results, which are contrasted against several recent  and competitive ML architectures. Figures are well designed and support the message conveyed in the manuscript.\n\n_Clarity:_ The manuscript is well organized, structured, and written. A rich appendix provides details about the model design, yet a supplementary material to validate the results is missing.\n\n_Significance:_ Results appear significant in terms of how sparse the data is sampled. Three synthetic problems of varying difficulty, as well as a real-world problem demonstrate the applicability of the method. Results are reported in two metrics along with standard deviations, which helps assessing the quality of the forecasts.", "weaknesses": "1. Observation function is constrained to a normal distribution with fixed variance. It would be helpful to add arguments of this design choice, to what extend it limits the expressivity of the model, as well as for what problems this formulation is sufficient.\n2. Ablations showing the performance under different spatial and temporal sparsities would be highly informative to understand the quality and limitations of the model at different tasks. Presumably, e.g., Navier-Stokes likely depends on more dense samples compared to Shallow Water. Extending this ablation to the other benchmarked methods would also provide valuable insights about the models' data efficacy.\n3. Limitations are not reported. It would be valuable to understand the limits of the method, its computational cost, and the time this architecture needs to train. Also, it is unclear to which extend the method can directly be applied to a task at hand or how much fine tuning is involved.\n4. No runtime comparison of the different models provided. If I'm not mistaken, the model must be called for each spatial position of interest in each time step, which amounts to a large number of model calls. Thus, to extend on Table 1, please provide information about the runtime of the entire model when generating a rollout of a spatiotemporal sequence of frames.\n5. More details about the differences between the introduced method and AutoSTPP would be valuable, given that these two approaches perform almost equally well. For what reason is your method superior to AutoSTPP?\n\n_Minor Comments_\n- Typo in line 306, \"withing\"\n- $N_{\\text{ctx}}$ is unclear in Figure 4. What value does the variable take? Would be good to have the actual value. EDIT: C.1 provides this information; I thus suggest to refer to C.1 in the Caption of Figure 4.", "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": {}}, {"review_id": "yczfWOnso9", "reviewer": "Reviewer_Nyvc", "summary": "This paper presents a novel method for modeling spatiotemporal dynamical systems from point process observations. The model integrates techniques from neural differential equations, neural point processes, implicit neural representations, and amortized variational inference. The authors also introduce a technique to speed training by addressing a computational bottleneck in latent state evaluation. The experimental results demonstrate the effectiveness of the model on challenging spatiotemporal datasets.", "questions": "Questions are related to the weaknesses:\nCould you address the issue of interpretability and the Poisson process a bit more", "limitations": "", "rating": 8, "confidence": 2, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- The paper is well-written and technically sound. The methodology is clearly presented, and the experimental setup is detailed\n- The proposed model is technically sound. It effectively combines techniques from various fields, including neural differential equations, neural point processes and amortized variational inference\n- experiments and \"ablations studies\" are comprehensive, showing the impact of many parameters of the model", "weaknesses": "- While focusing on predictive capability and computational efficiency, discussing the interpretability of the model would enhance its value. Can something be said about the dynamical system?\n- A little more discussion around the limitation of the Poisson process, and potential solution would have been welcome.", "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:

Decision: Accept (Spotlight)