Summary
The paper investigates mono-implicit Runge--Kutta (MIRK) methods for learning dynamical systems from data. In particular, MIRK methods can be made explicit by introducing the external data into the solver step itself, leading to a more efficient integrator while keeping favorable stability, symmetry, and symplecticity properties. To handle noisy data, the paper proposes the "mean inverse integrator" as an efficient way to average multiple trajectories and learn meaningful vector fields from these. The methods are demonstrated in multiple numerical experiments.
Strengths
The paper is well-written, the proposed method is presented very clearly, and the main claims made are well-supported.
Weaknesses
I find the presentation of the explicit Runge--Kutta (ERK) baseline a bit confusing, and in particular counter-intuitive of my current understanding of the usage of ERK methods for inverse problems in the context of ODEs. Concrete questions are below, in "Questions".
This criticism also extends to insufficiently detailed baselines, and potentially to insufficient baselines overall.
Questions
See also "Weaknesses". Notably:
1. How exactly is the RK4 / ISO+RK4 baseline defined?
2. Why is (6) _the_ appropriate loss function? I would rather assume something of the type $\sum || \tilde{y}_n - \hat{y}_n ||$, where $\hat{y}_n$ is a state estimate computed by the method; for instance by applying a numerical solver to simulate a whole trajectory. In particular, extrapolating exactly from the last data point seems to be a strong restriction, in particular if the data point can be very noisy. I would not consider this "the inverse problem formulation", but rather a specific choice that is made, geared towards the specific algorithm of interest.
3. Related to the last question: Could you commend on the standard least-squares approach of integrating the whole trajectory from some initial value (estimate) with RK, as opposed to the described method of just integrating from one data point to the next? The latter seems to be related to "multiple shooting"-type approaches which, while they can be very beneficial, seems to be far from the standard for neural ODEs. Though, to the best of my knowledge, even in multiple shooting it is common to extrapolate some state estimates $\hat{y}$; though part of the objective is the closeness of these estimates to the data, and they are initialized on the data.
Other remarks:
- l. 201: "evaluationg the vector filed in $\tilde{y}_n$, as is done in all explicit RK methods": To the best of my understanding this is not done in all explicit RK methods, but only in multiple shooting-type approaches (and even there this might not hold). A standard least-squares approach simulates a whole trajectory and then computes an L2-loss---and thus does not need to evaluate $f$ on the data points themselves.
- l. 129: $\tilde{y}$ is missing an index
- Related work: The write-up is very helpful; but I think it could be further improved by highlighting better how the proposed paper relates to these methods and highlight similarities and differences.
- l. 72-73: I agree with the sentiment of the statement, but I think the strict _necessity for numerical integrators_ is somewhat disputed in the literature. Notably, "gradient matching" methods compute first an interpolant, and then adjust the vector field to match the known derivative of the interpolant---this circumvents using RK or any related numerical method, which is why such methods often claim to be "simulation-free".
Rating
6: Weak Accept: Technically solid, moderate-to-high impact paper, with no major concerns with respect to evaluation, resources, reproducibility, ethical considerations.
Confidence
4: You are confident in your assessment, but not absolutely certain. It is unlikely, but not impossible, that you did not understand some parts of the submission or that you are unfamiliar with some pieces of related work.
Limitations
The authors address limitations in a dedicated section, which is much appreciated.
One additional limitation that is not explicitly mentioned in the section that comes to mind is the influence of the noise on the data: Since the data is explicitly included into the numerical solver, as opposed to being just part of some L2 loss to guide some estimated trajectory, I would expect that for very noisy data other methods might be preferrable (in particular a least-squares approach with RK from a learned initial value).
Another potential limitation: I assume that this method is not able to be a plug-in replacement in some latent neural ODE setting, where the ODE trajectory is not observed directly but only in some transformed space, e.g. when having video data and modeling an ODE in latent space? This is of course far from a trivial setting and I do not expect that such specific scenarios need to be mentioned explicitly; but the necessity of actual trajectory observations, as opposed to partial or non-linear observations, could indeed be another limitation.