Learning Spatiotemporal Dynamical Systems from Point Process Observations

**Q1:** 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? **A1:** In Section 4.1, we describe the structure and function of each component of our model, which makes it easier to understand and interpret. However, interpreting the meaning and structure of the learned latent space and dynamics is more difficult. Indeed, deep generative models, despite their flexibility and strong predictive performance, are generally hard to interpret. However, this is a limitation of deep generative models in general, and not something specific to our approach. Interpretability of deep generative models would require a separate study, which is outside the scope of our present work. ------------------------------- **Q2:** A little more discussion around the limitation of the Poisson process, and potential solution would have been welcome. **A2:** In line 497 we discussed what we believe are the main limitations of using the Poisson process. These limitations could be alleviated by using different point process models and adapting the architecture to efficiently account for potential interactions between the observation events and system dynamics. We will extend the discussion on limitations but leave the extensions to be considered in future work where such modeling assumptions are needed. We will extend the discussion of these two topics and add it to the revised manuscript.

**Q1:** Observation function is constrained to a normal distribution with fixed variance... **A1:** Our model supports and can be easily extended to other observation likelihood models. This would only require a reparameterization such that the observation function $g$ (in Eq. 11) returns parameters of the new distribution. -------------------- **Q2:** It would be valuable to understand the limits of the method, its computational cost, and the time this architecture needs to train... **A2:** As we mention in line 318, our method requires at most 1.5 hours of training. We also discuss limitations of our model at the end of Section 7. Our method, as many others, requires hyperparameter tuning to achieve best predictive performance, but it is rather robust to many hyperparameters, with the latent state dimension $d_z$, and complexities of the dynamics function $f$ and latent state decoder $\phi$ being most important hyperparameters to tune. -------------------- **Q3:** No runtime comparison of the different models provided ... Please provide information about the runtime of the entire model when generating a rollout of a spatiotemporal sequence of frames. **A3:** In our experiments, system dynamics simulation is the largest contributor to the entire model runtime. For example, on the Navier-Stokes dataset, which has the largest number of positions of interest, the runtime split across the model components is: encoder - 10\%, dynamics - 63\%, decoder - 27\%. Note that these numbers are with the latent space interpolation. Without it, dynamics simulation effectively takes all of the model runtime. -------------------- **Q4:** More details about the differences between the introduced method and AutoSTPP would be valuable... **A4:** The major advantage of our model over AutoSTPP is that our model leverages the observed system states $y$ to better model the observation process of $t$ and $x$, while AutoSTPP doesn't. -------------------- **Q5:** 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? **A5:** Modeling observation times and locations is fundamental for spatio-temporal point processes, which our method combines with simultaneously modeling the observations of the underlying spatio-temporal dynamics. Table 2 indeed shows that if the sensor locations are known, modeling them does not improve the model's predictive accuracy. However, systems we consider in our work produce observations at random spatiotemporal locations, so the sensor locations are unknown at test time. This means we need to predict where the observations will be made. Modeling sensor locations takes ≈15\% of the total runtime. -------------------- **Q6:** Have you tried higher-order interpolations? What is the error that incurs from the interpolation compared to modeling the full temporal grid? **A6:** We didn't consider higher order methods because the nearest neighbor and linear interpolation approach the full temporal grid errors for sufficient number of interpolation points ($n=50$ in our case). -------------------- **Q7:** Have you explored other solvers beyond dopri5, such as Euler, which is much cheaper? ... 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? **A7:** Yes, the time steps $\tau_1, \dots, \tau_n$ are separated by a constant time, which allows to use e.g., the Euler solver. However, our method is not constrained by this assumption and works with both regular and irregular time grids as well as fixed and adaptive step sizes. We tried both Euler and dopri5 solvers with 20 interpolation points ($n=20$), and the training time difference was small (41 min for Euler, and 49 min for dopri5); the predictive performance for both solvers was similar as well (within 3\% of each other). -------------------- **Q8:** How does the latent space dimensionality affect the runtime? Might be interesting to report along with its effect on the parameter count around line 375. **A8:** To answer this question we measure the runtime using on our Navier-Stokes dataset as this is the largest dataset in our work and we use the model with the largest number of parameters for it. We measure the total training time and the total number of model parameters as a function of the latent space dimension $d_z$: | $d_z$ | Time | Parameters | |-----|--------|---------| | 368 (default) | 45 min | 3022210 | | 184 | 43 min | 2721002 | | 92 | 43 min | 2570398 | | 46 | 43 min | 2495096 | We see that the both the training time and number of parameters are rather weakly affected by the latent space dimension. We thank the reviewer for careful reading and comments. We will incorporate the above comments and clarifications into the revised manuscript.

**Q1:** 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)? **A1:** As we describe on line 184, we use a low-dimensional latent state because it allows to simulate the dynamics considerably faster than a full-grid spatiotemporal discretization. -------------------------- **Q2:** 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. **A2:** Indeed, deep generative models, despite their flexibility and strong predictive performance, are generally hard to interpret. However, this is a limitation of deep generative models in general, and not something specific to our approach. We also note that, if needed, various model properties might be enforced, for example with penalty losses, architectural choices or even expert-defined parametric model components, thus giving a degree of control over the model and its interpretability. -------------------------- **Q3:** This encoder-based latent space modeling approach might not be able to model general PDE systems. **A3:** Following previous studies, we addressed this concern by testing our model on a wide range of challenging PDE systems, where it showed strong predictive performance. While this is not a proof that our model works for all conceivable PDE systems, this a strong evidence that it works for a wide range of realistic and practical scenarios. -------------------------- **Q4:** Why is the model claimed to be continuous? **A4:** Our model is continuous because it defines the system state at any arbitrary spatial and temporal location, rather than restricting it to predefined discrete grids. We note that all differential equation systems (except selected toy examples that allow closed-form solutions) need some numerical solution methods, such as the adaptive dopri5 solver that we used in our study.

**Q1:** Can you give me a possible explanation of why it works [on sparse data]? **A1:** Our model was designed to work on sparse data by 1) modeling the observation times t and locations x via a point process (Eq. 10), and 2) using the observed states y to improve the modeling of observation events (Eq.11). These features allow us to handle sparse observations made at arbitrary spatiotemporal locations, which is in contrary to all previous methods. To the best of our knowledge, previous methods can model only observation times and locations using a point process, or they can model spatio-temporal dynamics assuming that the observations are available at pre-defined time points or dense spatial locations, but but not both. ---------------- **Q2:** 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? **A2:** Transformer was used as a flexible general-purpose architecture supporting training-time parallelization. From our experience, RNN-based encoders could give results on-par with Transformer, but require longer training times. ---------------- **Q3:** 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? **A3:** To answer this question we compared the training times with neural ODE and a map-based dynamics (where the next state $z_{i+1}$ is evaluated as $z_{i+1} = f(z_i)$). We used the Navier-Stokes dataset as this is the largest of our datasets. As the ODE solver we used dopri5 with high absolute and relative tolerances of $10^{-5}$, which ensured accurate solutions. As the result, training with neural ODE dynamics took 45 minutes, while using map-based dynamics took 40 minutes. We also measured training times with lower absolute and relative tolerances set to $10^{-3}$, which resulted in 42 minutes-long training. So, we see that while using neural ODEs as the dynamics model increases the training time, the increase is quite modest.