Neural Jump Ordinary Differential Equations: Consistent Continuous-Time Prediction and Filtering
Authors: Calypso Herrera
Florian Krach
Josef Teichmann
Venue: ICLR 2021 Poster (2021)
Keywords: Neural ODE
conditional expectation
irregular-observed data modelling

Abstract: Combinations of neural ODEs with recurrent neural networks (RNN), like GRU-ODE-Bayes or ODE-RNN are well suited to model irregularly observed time series. While those models outperform existing discrete-time approaches, no theoretical guarantees for their predictive capabilities are available. Assuming that the irregularly-sampled time series data originates from a continuous stochastic process, the $L^2$-optimal online prediction is the conditional expectation given the currently available information. We introduce the Neural Jump ODE (NJ-ODE) that provides a data-driven approach to learn, continuously in time, the conditional expectation of a stochastic process. Our approach models the conditional expectation between two observations with a neural ODE and jumps whenever a new observation is made. We define a novel training framework, which allows us to prove theoretical guarantees for the first time. In particular, we show that the output of our model converges to the $L^2$-optimal prediction. This can be interpreted as solution to a special filtering problem. We provide experiments showing that the theoretical results also hold empirically. Moreover, we experimentally show that our model outperforms the baselines in more complex learning tasks and give comparisons on real-world datasets.

Published as a conference paper at ICLR 2021
NEURAL JUMP ORDINARY DIFFERENTIAL EQUATIONS: CONSISTENT CONTINUOUS- TIME PREDICTION AND FILTERING
Florian Krach Josef Teichmann Department of Mathematics, ETH Zurich, Switzerland (firstname, last name ) (maath. ethz. ch
Calypso Herrera
ABSTRACT
Combinations of neural ODEs with recurrent neural networks (RNN), like GRU- ODE-Bayes or ODE-RNN are well suited to model irregularly observed time series. While those models outperform existing discrete-time approaches, no the- oretical guarantees for their predictive capabilities are available. Assuming that the irregulall-samppled time series data originates from a continuous stochastic process, the ooPtpttimll online prediction is the conditional expectation given the currently available information. We introduce the Neural Jump ODE (NJ-ODE) that provides a data-driven approach to learn, continuously in time, the conditional expectation of stochastic process. Our approach models the conditional expec- tation between two observations with a neural ODE and jumps whenever a new observation is made. We define novel training framework, which allows us prove theoretical guarantees for the first time. In particular, we show that the output our model converges to the opptmml prediction. This can be interpreted as solution to special filtering problem. We provide experiments showing that the theoretical results also hold empirically. Moreover, we experimentally show that model outperforms the baselines in more complex learning tasks and give comparisons on real-world datasets.
INTRODUCTION
Stochastic processes widely used in many fields model time series that exhibit random behaviour. In this work, we focus on processes that can be expressed as solutions of stochastic differential equations (SDE) of the form
dx, (tt,, Xt)dt + (t,, X,),W,,
with certain assumptions on the drift and the diffusion a. With respect to the 2 the best prediction of future value of the process provided by the conditional expectation given the current value. the drift and diffusion are known or a good estimation is available, the conditional expectation can be approximated by Monte Carlo (MC) simulation. However, since f and are usually unknown, this approach strongly depends on the assumptions made on their parametric form. more flexible approach is given by neural SDEs, where the drift Jæ and diffusion are modelled by neural networks (Tzen & Raginsky, 2019; Li et al., 2020; Jia & Benson, 2019). Nevertheless, modelling the diffusion can be avoided one is only interested in forecasting the behaviour instead of sampling new paths.
alternative widely used approach to use Recurrent Neural Networks (RNN), where neural network dynamically updates latent variable with the observations of a discrete input time-series. RNNs are successfully applied to tasks for which time-series are regularly sampled, as for example speech or text recognition. However, often observations are irregularly observed in time. The standard approach of dividing the time-line into equally-sized intervals and imputing or aggregating observations might lead significant loss of information (Rubanova al., 2019). Frameworks that overcome this issue are the GRU-ODE GRU-ODE-Bayes (Brouwer et al., 2019) and the ODE-RNN (Rubanova al., 2019), which combine RNN with neural ODE (Chen al., 2018). 8). In standard RNNs, the

---

Published as a conference paper at ICLR 2021
hidden state is updated at each observation and constant in between. Conversely, in the GRU-ODE- Bayes and ODE-RDERRR IDE-RNN framework, a neural ODE is trained to model the continuous evolution of the hidden state of the RNN between two observations. While GRU-ODE-Bun -B-yes and ODE-RNN both provide convincing empirical results, they lack thorough theoretical guarantees.
Contribution. In this paper, we introduce mathematical framework to precisely describe the problem statement of online prediction and filtering of a stochastic process with temporal irregular observations. Based on this rigorous mathematical description, we introduce the Neural Jump ODE (NJ-ODE). The model architecture is very similar to the one of GRU-ODE-Bayes and ODE-RNN, however we introduce novel training framework, which in contrast to them allows us to prove convergence guarantees for the first time. Moreover, we demonstrate empirically the capabilities of our model.
Precise problem formulation. We emphasize that precise definition of all ingredients is needed, to be able to show theoretical convergence guarantees, which is the main purpose of this work. Since the objects of interest are stochastic processes, we use tools from probability theory and stochastic calculus. To make the paper more readable and comprehensible also for readers without background in these fields, the precise formulations and demonstrations of all claims are given in the appendix, while the main part of the paper focuses on giving well understandable heuristics.
2 PROBLEM STATEMENT
The problem we consider in this work, is the online forecasting of temporal data. We assume that we make observations of a Markovian stochastic process described by the stochastic differential equation (SDE)
dx, = p(t, dt + (t, X,)dW,
(1)
at tiregularlyssamples time points. Between those observation times, we want to predict the stochastic process, based only on the observations that we made previously in time, excluding the possibility to interpolate observations. Due to the Markov property, only the last observation is needed for an optimal prediction. Hence, after each observation we extrapolate the current observation into the future until we make the next observation. The time which the next observation will be made is random and assumed to be independent of the stochastic process itself.
More precisely, we suppose to have a training of N independent realisations of the Rdx- dimensional stochastic process X defined in (1). Each realisation j is observed at ny random observation times t(() ......t.. tny E [0, with values x.(3) .2) R We assume that all coordinates of the vector a!!) are observed. We are interested in forecasting how new independent realization evolves in time, such that our predictions of X minimize the expected squared distance (L'-metric) to the true unknown path. The optimal prediction, i.e. the L m'iiimiizer, is the con- ditional expectation. Given that the value of the new realization at time t is I₁, we are therefore interested in estimating the function
.(i)
[z,, E| = æz), s s2 0,
which is the ooptiil prediction until the next observation is made. To learn an approximation f of f we make use of the N realisations of the training set. After training, f is applied the new realization. Henee, this can be interpreted as special type of filtering problem. The following example illustrates the considered problem.
Example. A complicated to measure vital parameter of patients in hospital is measured multiple times during the first 48 hours of their stay. For each patient, this happens at different times de- pending on the resources, hence the observation dates are irregular and exhibit some randomness.
( = (1, 14, 27, 34) where the values (0.74, 0.65, 0.78. 0.81) are measured. Patient 2 only has n2 = 2 measure- ments at hours (t(2) (2) (3, 28) where the values (2) (2, (0.56, 0.63) are measured.
Patient 1 has n1 4 measurements at hours (tf)) .2 b3
(1) (1)
(af") .2 x3
Similarly, the j-th patient has nj measurements at times (1( Lny and has the measured values (2)).
(x(1) ....,.... xng) ). Based on this data, we want to forecast the vital parameter of new patients coming to the hospital. In particular, for patient with measured values xi at time 41, want predict what
2

---

Published as a conference paper at ICLR 2021
his values will likely be any time + > t1. Importantly, we do not only focus predicting the value at some 2o > ti, but we want to know the entire evolution of the value.
BACKGROUND
Recurrent Neural Network. The imput to RNN is discrete time series of observations (Z1,-.. , ). At each observation time 11+1, neural network, the RNNCell, updates the latent variable using the previous latent variable h2 and the input 2++1 as
h++1 == (NN Tx+1)
Neural Ordinary Differential Equation. Neural ODEs (Chen et al., 2018) are family of continuous-tiin timm models defining a latent variable hę := ht) be the solution an ODE initial-value problem
h ho [[ [(h,, 0)ds, to:
where f(.,,, fo is neural network with weights A. Therefore, the latent variables can be updated continuously by solving this ODE (3). We can emphasize the dependence h on a numerical ODE solver by rewriting (3)
he = ODESolve(Sa, hog (lo, ())
ODE-RNN. ODE-RNN (Rubanova et al., 2019) is mixture of a RNN and a neural ODE. In contrast to standard RNN, we are not only interested in an output at the observation times 4, but also in between those times. In particular, time. This is achieved by using a neural ODE to model the latent dynamics between two observation times, i.e for t < < the latent variable is defined as in (3) and (4), with ho and to replaced by h4-1 and 1 At the next observation time t,, the latent variable updated by a RNN with the new observation 2;. Fixing ho. the entire latent process can be computed by iteratively solving ODE followed by applying RNN. Rubanova et al. (2019) write this as
want to have an output stream that is generated continuously in
{ M/ ODESolve (aaSolleeJo, h4- (4-),44) 14 = (k,,
(5)
GRU-ODE-Bayes. The model architecture deseribing the latent variable in GRU-ODE-Bayes (Brouwer et al., 2019) is defined as special case of the ODE-RNN architecture. In particular, gated recurrent unit (GRU) is used for the RNN-cell and continuous version of the GRU for the neural ODE fa. Therefore, we focus on explaining the difference between our model architecture and the ODE-RNN architecture, in the following section.
PROPOSED METHOD NEURAL JUMP ODE
Markovian paths. Our assumptions the stochastic process X imply that is Markov process. In particular, the optimal prediction of future state of X only depends on the current state rather than on the full history. Hence, the previous values do not provide any additional useful information for the prediction.
JumpNN instead of RNN. Using neural ODE between two observation has the advantage that allows to continuously model the hidden state between two observations. But since the underlying process Markov, there is no need use RNN-cell to model the updates of the hidden state each new observation. Instead, whenever new observation is made, the new hidden state can solely be defined from this observation. Therefore, we replace the RNN-cell used in ODE-RNN by a standard neural network mapping the observation the hidden state. We call this the jumpNN which can be interpreted as an encoder map. Compared to ODE-RNN. this architecture is easier to train.
Last observation and time increment additional inputs for the neural ODE. The neural net- work fo used in the neural ODE takes two arguments as inputs, the hidden state hę and the current time t. However, our theoretical problem analysis suggests, that instead of t the last observation time f;-1 and the time increment t f should be used. Additionally the last observation 2;-1 should also be part of the input.
3

---

Published as a conference paper at ICLR 2021
NJ-ODE. Combining the ODE-RNN architecture (5) with the previous considerations, introduce the modified architecture of Neural Jump ODE (NJ-ODE)
{ 14 = ODFSolvee[Jo, (hieeeeooo.hii,aea h4 = jumpNN(a.).
(4--1;4)
(6)
An implementable version of this method is presented in the Algorithm 1. A neural ODE fe transforms the hidden state between observations, and the hidden state jumps according to jumpNN when a new observation is available. The outputNN, standard neural network, maps any hidden state h, to the output yt- To implement the continuous-in in- time ODE evaluation, a discretization scheme is provided by the inner loop. In the training process, the weights of all three neural networks, jumpNN, the neural ODE fo and outputNN are optimized.
Algorithm 1 The NJ-ODE. A small step size Al is fixed and we denote bn+1 := T.
Input: Data points with timestamps {(zi, (i)
for i 0 to do
h, = (umpNN(x:)
Update hidden state given next observation z; compute output
yt, outputNNN (h+:)
s <- 4
while s At 1++1 do
hs ODE ha, t -1,,1,,0 (8, At))
>> get next hidden state compute output
Ms+ OutputNN(hs+Aa)
8 f + At
end while
end for
Objective function. Our goal is to train the NJ-ODE model such that its output approximates the conditional expectation (2), which is the optimal prediction of the target process X with respect to the normn Therefore,
can prove convergence. Let yi- denote the output of the NJ-ODE at & before the jump and yı the output at & after the jump. Note that the outputs depend on parameters and the previously observed 2; which are inputs the model. Then the objective function is defined as
define
objective function, with which
(0) =
(
+
)²
N
jump part at observations
continuous part
paths
dates
between observations
We give an intuitive explanation for this definition. The "jump part" of the loss function forces jumpNN to produce good updates based on new observations, while the other part forces the jump size to be small in (the empirical) 2 Since the conditional expectation minimizes the jump size with respect to the Lo-nom, this forces the neural ODE fe to continuously transform the hidden state such that the output approximates the conditional expectation. Moreover, both parts of the loss function force the outputNN to reasonably transform the hidden state h2 to the output y/-
5 MAIN RESULT THEORETICAL CONVERGENCE GUARANTEE
In the following we informally state our main result. To formally state and prove the theorem, precise definitions of all ingredients are needed. This analysis is provided the appendix where the following results is stated in Theorem E.2 and Theorem E.
Theorem 5.1 (informal).
assume that for each number of paths N for every size of the neural networks M, their weights are chosen optimally, as to minimize & x Then, if N and M tend to infinity, the output of NJ-ODE converges in mean (L -convergence) to the conditional exception of the stochastic process X given the current information.
An intuitive explanation for this theorem given with the definition of the objective function. In this result, the focus lies on the convergence analysis under the assumption that optimal weights are found. In the Appendix discuss why this assumption restrictive.
4

---

Published as a conference paper at ICLR 2021
6
EXPERIMENTS
For further details and results for all experiments see Appendix F.
6.1 TRAINING ON SYNTHETIC DATASETS
Evaluation metric. For synthetic datasets where analytic formula for the conditional expectation
exists, we can evaluate the distance of the model output to the target process (2). We use a sampling
time grid with equidistant step size A, == K. K N, [0, 7). On this grid, compare for path
at time t, the true conditional expectation a with the predicted conditional expectation (the model
output) y/ (f), For N₂ test samples, the evaluation metric defined as
eval(2,
1 N2 K
2 LA.
()) Vid.)
t=0
Black-Scholes, Ornstein-Uhleebeeck and Heston. We test our algorithm on three scalar stochastic models, Black-Scholes, Ornstein-Uhlenbeed and Heston, with fixed parameters. For each model, generate a dataset by sampling N = 20'000 paths on the time interval [0, 1] using the Euler-scheme with 100 time steps. Independently for each path, on average 10% of the grid points are randomly chosen as observation times. The NJ-ODE is trained on 80% of the data. On the remaining 20% the model is tested, by comparing the loss function (31) computed with the NJ-ODE the loss function computed with the true conditional expectation (Figure 2). During training, the relative difference becomes very small, hence the true conditional expectation nearly replicated.
10
true path
20
true conditional expectation observed
2.0
M.
0.3
2-
Figure 2: The relative difference between the loss function computed with the predicted conditional expectation and the loss function computed with the true conditional expectation is shown with respect the number of epochs.
0.0
0.2
0.4
0.6
0.8
1.0
Figure I: Predicted and true conditional expecta- tion on a test sample of the Heston dataset.
6.2 FURTHER SYNTHETIC DATASETS
Heston model without Feller condition. We also train our model on Heston dataset for which the Feller condition is not satisfied. As explained by Andersen (2007); Jean-Françee et al. (2015) this more delicate situation. We see that our model very robust, since even this critical case learns to replicate the true conditional expectation process. In the Heston model, the variance of the stochastic process X, is a stochastic process, vr. Here, we train our model to predict both processes at the same time. The training on this 2-dimensional dataset is successful as can be seen in Figure 3. The minimal evaluation metric after 200 epochs is 0.0983.
Dataset with changing regime. In this experiment we test how well model can deal with stochastic processes, that undergo an (abrupt) change of regime at a certain point in time. Many real world time series might exhibit such change of regime. Some examples are listed below.
Longitudinal patient health recordings might experience changes depending on seasonal or longer- term influences, for example due to the seasonal flue currently the Covid-19 19 pandemic. many regions climate data has strong seasonal dependencies, that can lead to relatively abrupt changes for example when the weather changes from dry to rain season.
5

---

Published as a conference paper at ICLR 2021
2.5
2.0
1.5
1.0
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
G.4
0.6
0.8
1.6
2.5
True path
true path
our mocel
model
tuue conditional expectation observed
consitional expectation
1.0
1.0
0.5
0.0
0.2
0.9
Figure 3: Heston model without Feller condition. In both plots, the upper sub-plot corresponds the 1-dimensional path of X₁ and the lower sub-plot corresponds to the 1-dimensional path of te.
A stock market that suddenly changes from bullish to a bearish market, for example due to macro-economic event. An example for this would be the start of the Covid-19 crisis in the first quarter 2020.
We test change of regime by combining synthetic datasets. On the first half of the time interval [0, 0.5] we use the Menstinn-hllenbeck and on the second half [0.5, 1] the Black-Scholes model. In Figure we that our model correctly learns the change of regime. The minimal evaluation metric after 200 epochs is 0.0463.
true
path model conditional
conditiona expectation
10
0-
4-
0.0 0.2 0.1
0.5
0.8
1.0
0.0
0.2
0.4
0.5
0.8
1.0
Figure 4: Our model evaluated o stochastic dataset that follows an time interval [0, 0.5] and an Black-Scholes model on the time interval (0.5, 1]. We see that our model correctly learns the change of regime.
SDE on the
Dataset with explicit time dependence. Many real world datasets, have an explicit time dependence, i.e. the drift and diffusion of (9) explicitly depend f. Examples all datasets that have certain periodicity, as for example weather data (seasonal and daily periodicity), intraday periodicity of stock prices (Andersen Bollerslev, 1997) or prices of certain seasonal goods. We incorporate an explicit time dependence into the Black-Scholes dataset, by replacing the drift constant f with the time dependent constant / (sin(3t) 1), for a, B > In Figure we see that the model learns to adapt to the time-dependent coefficients. The minimal evaluation metric after 100 epochs is 0.0215 (3 2x) and 0.02805 4T) respectively.
6.3 EMPIRICAL CONVERGENCE STUDY
We confirm the theoretical results of Theorem 5.1 by an empirical convergence study for growing numbers of training samples N1 and network size M, where the performance is measured by the evaluation metric (8). For each combination of the number of training samples N, and the neural network size M. the NJ-ODE is trained times the Black-Scholes, Orstein-UUlenheck and

---

Published as a conference paper at ICLR 2021
expectation
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.4
0.6
0.8
1.0
Figure 5: NJ-ODE evaluated the time-dependent Black-Scholes dataset with 2x (left) and
3 4m (right).
10
retwork 10
10
training 200
retwork 20
relaark si
training size-ice
traiiigg 1800
retwark 16
traning_irre-2000
retwork seee
E trainingsssee sizemmaoo
trainiggiiee simngggize-2200
evalneree
evalmetre
10
traiiing size
network
Black-Scholes dataset. Mean + standard deviation (black bars) of the evaluation metric for varying training samples N₁ and network size M.
Figure 6:
Heston datasets. The means of the evaluation metric over the 5 trials with their standard deviations are shown in Figure 6 for Black-Scholes. Already 200 training samples lead to small evaluation metric if the network size is big enough, suggesting that the Monte Carlo approximation of the loss function is good, already with few samples.
6.4 COMPARISON GRU-ODE-BATE BAYES ON SYNTHETIC DATA
On the Black-Scholes, ODE- -DE-Bayes (Brouwer et al., 2019), which is, to the best of our knowledge, the neural network based method with the most similar task to ours. Results of our comparison are shown in Table 6.4. On the Black-Schole and Ornstein- Uhweenbb--lhbbkkkkke dataset, our model performs similarly to GRU-ODE-Baye is unstable and
and Heston datasets, we compare our model to GRU-
However, on the more complicated Heston dataset, the training of GRU-ODE-Bay does not converge. On the other hand, our model converges during training. Although the value of the evaluation metric is much higher than for Black-Schole and Constein-Uhlebed the resulting model output is still meaningful, which is not the case for GRU-ODE- GRU-ODE-Bayes. Hence we conclude, that GRU-ODE-Bayy DDDE-Bayes cannot be applied reasonably for the Heston dataset, while method works as the theoretical results suggest.
6.5 REAL WORLD DATASETS WITH INCOMPLETE OBSERVATIONS
Self-imputation for incomplete observations. Until now we assumed that at each observation time all coordinates of the stochastic process X are observed. However, in many real world applications, the observations are incomplete, ie. only for some of the coordinates an observation is available. To deal with such incomplete observations, propose the following self-imputation method. Whenever
7

---

Published as a conference paper at ICLR 2021
Table 1: The minimal, last and average value of the evaluation metric throughout the 100 epochs of training are shown for BRU-ODE-Bayee and our method, together with the number of trainable parameters.
Black-Scholes
Orsttii- Uhlenbeck
Heston
GRU
ours
GRU
ours
GRU
ours
min
7 10
7 10
3 10 1
5- 10
4.24
1.33
last
3. 10-3
8 10-3
1. 10-1
6. 10-4
30.9
1.35
average
3. 10-0
1.10-2
7. 10-4
8. 10-4
56.2
1.65
params
112'602
10'071
112'602
10'071
112'602
10'071
observation 2; Rx is made the mask my {0, 1 tells whether the k-th coordinate of the vector x; is observed (m = 1) or not (m = 0). Then we use the prediction yi- of our model to generate the imputed observation vector
f = my O 24
(Idx m;) ④ yi--
where 3 is the element-wise multiplication (Hadamar product) and 14x Rdx the one-vector. Instead of v, we use (ãt, m;) an input for the jump part jumpNN. The intuition behind this definition is the following. In the one dimensional case, we do not make an observation, but input Yt- if it was observation, this point of view, Yt- does not provide any additional information for the model. Similarly, we expect that imputing Vr- for unobserved coordinates does not provide any information about this coordinate to the model. However, since the model might learn to transfer the information about an observed coordinate to an unobserved one, we extend the input also include the information which coordinates were observed. For the ODE part we use yi. the prediction after processing the input (á:, m,), as input instead of x;, Here, the intuition that if the model learns how to best use the incomplete observation, i.e. if the jump is good, then this is the best approximation for Our objective function adjusted by multiplying each term the sum with the mask.
expect that this does not change the output ys for 2 t. From
Climate forecast. We compare our model to GRU-ODE-Bay Bay--Bayesss on the USHCN daily dataset (Menne al., 2016), using the same experimental setting as was used by Brouwer al. (2019). We train small (S) and large (L) version of NJ NJ-ODE with different total number of parameters. The validation set was used for early stopping after the first 100 epochs, where we trained for a total 200 epochs. We see in Table 6.5 that our small version performs slightly worse than GRU-ODE-Bayee while large version slightly outperforms it.
Table 2: Mean and standard deviation of MSE on the test sets of USHCN. Result of baselines were reported by Brouwer et al. (2019). Where known, the number of trainable parameters is reported.
USHCN MSE # params
neural ODE-VAL VAE neural ADE-VAE-MASK
0.96 + 0.11
0.83 + 0.10
sequential VAE GRU-SIMPLE
0.83 + 0.07
0
GRU-D
0.53 + 0.06
T-LSTM
0.59 0.11
GRU-ODE-Bayes
0.43 + 0.07
42'640
NJ-ODE (S)
0.45 + 0.06
10'925
NJ-ODE (L)
0.40 + 0.07
571
Physionet. We compare our model to the latent ODE their extrapolation task the PhysioNet Challenge 2012 dataset (Goldberger et al., 2000), using the same experimental setting as was used by Rubanova et al. (2019). The mean and standard deviation over 5 runs starting different random initializations is reported for our model. We see in Table 6.5 that our model outperforms the latent ODE models although having only about seventh of the trainable weights.
8

---

Published as a conference paper at ICLR 2021
Table 3: Mean and standard deviation of MSE the test set of physionet. Result of baselines were
reported by Rubanova et al. (2019). Where known, the number of trainable parameters is reported.
Physionet - MSE (x10-*
params
RNN-VAE
3.055 0.145
Latent ODE (RNN enc.)
3.162 + 0.052
Latent ODE (ODE enc)
2.231 0.029
163'972
Latent ODE + Poisson
2.208 = 0.050
181'723
NJ-ODE
1.945 0.007
24'423
RELATED WORK
Stochastic filtering theory. Our main theorem is similar to the theoretical results of neural filtering (Lo, 2009), which neural network approach to stochastic filtering theory in discrete time. In stochastic filtering, potentially incomplete and noisy observations of X are available continuously in time. This observation process is usually described by the dynamics dY, h( X, )tt + dB,, where h is measurable function and B is Brownian motion. Stochastic filtering then estimates the conditional law of X, given the noisy observations (Y=)osasr (Bain & Crisan, 2008, Def. 3.2) and therefore can provide conditional expectations. Comparably, we directly compute the conditional expectation given the last observation. Similar to our assumptions, in the neural filtering approach of only available
GRU-ODE-Bur -Byyss (Brouwer et al., 2019) and (Ryder et al.. 2018) observations irregular discrete time points. They approximate the conditional law of given the last observation by Gaussian distribution and learn its mean and variance parameters. In particular, the conditional expectation is then given by this mean parameter. In contrast, we do not make normality assumptions about the conditional distribution and we theoretically prove convergence the true conditional expectation.
Neural ODEs with jumps. Except for GRU-ODE-Bayes (Brouwer et al., 2019) and ODE-RNN (Rubanova et al., 2019) another work studying neural ODE with jumps is Neural Jump SDE (NJSDE) (Jia & Benson, 2019). Similar to the NJ-ODE framework (29), the latent process of NJSDE is described by neural ODE with jumps at random times. This model is used to describe hybrid systems which evolve continuously in time but may also interrupted by stochastic events. In contrast to that, we model the conditional expectation of continuous stochastic process.
8 CONCLUSION
We presented the Neural Jump ODE, data-driven framework for modelling the conditional expecta- tion of stochastic process given the previous observations. We introduced a rigorous mathematical description of model and more generally for the class of neural ODE based models. Moreover, for the first time we provided theoretical guarantees for model falling in this category. We evaluated our model empirically on six synthetic and two real world datasets. In comparison the baselines GRU-ODE-Bayes and latent ODE, we achieved better results especially complex datasets.
ACKNOWLEDGEMENT
The authors thank Andrew Allan, Robert A. Crowell, Anastasis Kratsios and Pierre Ruyssen for helpful discussions, providing references and insights. Moreover, the authors thank the reviewers for their thoughtful feedback that contributed to significantly improve the paper. The authors gratefully acknowledge financial support coming from the Swiss National Science Foundation (SNF) under grant 179114.
REFERENCES
Leif Andersen. Efficient simulation of the heston stochastic volatility model. Risk Management, 2007.

---

Published as a conference paper at ICLR 2021
Torben G. Andersen and Tim Bollerslev. Intraday periodicity and volatility persistence in financial markets. Journal of Empirical Finance, 1997.
Alan Bain and Dan Crisan, Fundamentals of stochastic filtering, Springer Science & Business Media, 60, 2008.
Edward De Brouwer, Jaak Simm, Adam Arany, and Yves Moreau. GRU-ODE-Bayes: Continuous modeling of sporadically-observed time series. NeurfPS, 2019.
Ricky T. Q. Chen, Yulia Rubanova, Jesse Bettencourt, and David Duvenaud. Neural ordinary
differential equations. NeurIPS, 2018.
Samuel N Cohen and Robert James Elliott. Stochastic calculus and applications. Springer, 2015.
Christa Cuchiero, Martin Larsson, and Josef Teichmann. Deep neural networks, generic universal interpolation, and controlled ODEs. 1908. 2019.
Rick Durrett. Probability: Theory and Examples. Cambridge University Press, 2010.
Soheil Feizi, Hamid Javadi, Jesse Zhang, and David Tse. Porcupine neural networks: (almost) all local optima global. arXiv:/710.02196, 2017.
Ary L. Goldberger, Luis A. N. Amaral, Leon Glass, Jeffrey M. Hausdorff, Plamen Ch. Ivanov, Roger G. Mark, Joseph E. Mietus, George B. Moody, Chung-Kang Peng, and H. Eugene Stanley, Physiobank, physiotoolkit, and physionet. Circulation, 2000.
Massimiliano Gubinelli. Stochastic analysis course note 4. Lecture Notes
https://www. iam, uni-bonn, dew.ennnmieeieeemiieiiiaaaa /wtp/uuueuuuuue httppuuuueeeeeeeeeeeeeeaeeaaaaet www.uuppllduuu..oo. gubinell//
stochastic-anal ssl.:aanssssssssllllsssl sote 4. pdf, 2016.
Calypso Herrera, Florian Krach, and Josef Teichmann. Estimating full lipschitz constants of deep neural networks. arXiv:2004.13135, 2020.
Kurt Hornik, Maxwell Stinchcombe, and Halbert White. Multilayer feedforward networks universal approximators. Neural Networks, 2, 1989.
Bégin Jean- -ran-François, Bédard Mylène, and Gaillardetz Patrice, Simulating from the heston model: A gamma approximation scheme. Monte Carlo Methods and Applications, 2015.
Junteng Jia and Austin Benson. Neural jump stochastic differential equations. NeurIPS, 2019.
Diederik Kingma and Jimmy Ba. Adam: A method for stochastic optimization. International
Conference on Learning Representations, 2014.
Bernard Lapeyre and Jérôme Lelong. Neural network regression for bermudan option pricing. arXiv:1907.06474. 2019.
Andrea Lecchini-Visintini, John Lygeros, and Jan Maciejowski. Simulated annealing: Rigorous finite-time guarantees for optimization on continuous domains. Advances in Neural Information Processing Systems, 2008.
Michel Ledoux and Michel Talagrand. Probability in banach spaces. Springer-Verlag, 62:67-69, 1991.
Xuechen Li, Ting-Kam Leonard Wong, Ricky T. Q. Chen, and David Duvenaud. Scalable gradients for stochastic differential equations. TrXiv:2001.01328 2020.
Ting-Ho Lo. Neural Filtering. Scholarpedia, 4, 2009.
Marco Locatelli. Simulated annealing algorithms for continuous global optimization: convergence conditions. Journal of Optimization Theory and applications, 2000.
M. J. Menne, N. Williams, Jr., and R. Vose. Long-term daily and monthly climate records from stations across the contiguous united states historstl climatology network). 2016.
10

---

Published as a conference paper at ICLR 2021
Philip Protter. Stochastic integration and differential equations. Springer- Verlag, 2005.
Yulia Rubanova, Ricky T. Q. Chen, and David K Duvenaud. Latent ordinary differential equations for iregularly-sampled time series. NeurIPS, 2019.
Reuven Y Rubinstein and Alexander Shapiro. Discrete event systems: Sensitivity analysis and stochastic optimization by the score function method. Wiley, 1993.
Sanae Rujivan and Song-Ping Zhu. A simplified analytical approach for pricing discretely-sampled variance swaps with stochastic volatility. Applied Mathematics Letters, 2012.
Thomas Ryder, Andrew Golightly, A. Stephen McGough, and Dennis Prangle. Black-box variational inference for stochastic differential equations. In ICML, 2018.
Belinda Tzen and Maxim Raginsky. Neural stochastic differential equations: Deep latent gaussian models in the diffusion limit. arXiv: 099833 2019.
Boxuan Yue, Junwei Fu, and Jun Liang. Residual recurrent neural networks for learning sequential representations. Information, 9, 2018.

---

Published as a conference paper at ICLR 2021
APPENDIX
A SETUP
A.1 STOCHASTIC PROCESS X
In this section rigorously describe the process X and give the assumptions which are needed in order to derive the convergence results. Let dx, dw E N and T > O be the fixed time horizon. Con- sider filtered probability space (2, F == {Fe JOSEST,,,, P), which an adapted dyy-dimensional Brownian motion (ejor is defined. We define the stochastic process' X := (Xe) the solution of the stochastic differential equation (SDE)
dx, 11l, X,)dt + (t,, X,)dW,,
for all 0 < t T, where Xo x E Rdx is the starting point and the measurable functions j [0, Rx + Rdx and [0. Rd Rxx xdw are the drift and the diffussio respectively. We impose the following assumptions:
X is continuous and square integrable, i.e. for P-a.e. w E Ω the map / + Xz(w) is continuous and E[X?] < for every t [0,
/ and both globally Lipschitz continuous in their second component, i.e. for yp E () there exists constant M 0 such that for all [0,
lp((, - ptt, v) & M z y/2 and le 2212 s (1+ |x/2)
(10)
In particular, their growth is at most linear the second component.
y is bounded and continuous in its first component (t) uniformly in its second compo- nent (I), i.e. for every E j0. T and > 0 there exists a 5 > 0 such that for all E [0. T] with - s| < 8 and x E Rdx we have (u((, uls, 2)] < e.
a is càdlåg (right-continuous with existing left-limit) in the first component and La integrable with respect to W, E L'WW), i.e
dx
= = 0) dl < co,
E
supory (x, W
(11)
where | denotes the Frobenius matrix norm. This is in particular implied if a is bounded.
A.2 RANDOM OBSERVATION DATES
In this section we describe how we model the observation dates. We want to treat irregularly observed time-series, i.e. without fixed number of observations and not necessarily equally distributed. Therefore, we suppose that we have a possibly random number of observations and that those observations are made random times ,,....;4, ln- For simplicity we make the assumption that the observation dates are independent of the stochastic process itself. In future work, this will be generalized to observation dates correlated with X, that could be modelled by a point process on the same probability space. In contrast, here we hard-code the independence assumption by considering second probability space ($2, F.P),, P), on which the random observation times of the stochastic process X are defined. More precisely, we assume that:
* n : 52 - Nzo is random variable with Ep(n] co, the random number of observations, 4 Ô - [0, T) for < S n are sorted 2 random variables, the random observation times.
We denote the joint pushforward measure of n and (1:ossinn as P. = (n, to....... In #P. The random variable n can but does not have to be unbounded. If bounded, we define K =
A stochastic process is a collection of random variables Xx + Rd w + Xi(w) for 0
For all a E f 0 < (G) (n( (w) T.
12

---

Published as a conference paper at ICLR 2021
max ( E N
NNP(( P(n 2 k) [))00] 0 to be the maximal value of n, which otherwise is infinity. We use the notation [([0, TT) for the Borel o-algebra the set [0, T. Then we define for each S S K
d [[O T)) [0, 1]:
B (B F(n2k,(1x2) (lx-))B) F(n2k)
which is probability measure on the time interval shown in the lemma below. Here, (4x-)
means the left-point of Is, for example, if to 0, = we have to E [0, 1).4 & [0, 1) but to- [0, 1), E [0, 1).
Lemma A.1. For 1 < k < K + the map A defines probability measure.
Proof. First we see that As ((0, T) - since is maps to [0. T). Furthermore, for any disjoint sets B, E [(0 T),, i N, we have {n 2 k; 4 E Ú1 B:) = Ú121 {n 2 k, 4- E B, and these sets F-measurable, since they defined through pre-images of random variables. Therefore, the additivity of I implies that N (Ú1z1 B,) Eiz1 A (B:).
Moreover, define the time of the last observation before a certain time T: [0, x 82 - [], (1,w) + T(1,)) =- max(4 < i (()), < []
Example A.2. We give two examples how the random observation dates could be defined.
Let (Nt)ue(o,T) be homogeneous Poisson point process with rate r >> 0. Hence, N(0) 0. f + (t) is constant except for jumps of size 1 at discrete random times and for any fixed E [0, T), ((f) Poisson distributed, i.e.
PtN
k- (rt) e for k E N.
Then the number of observations is defined as n= N(T) and satisfies Eg[n] rT. The observation dates b1, ..... -.n In are defined as the discontinuity times N, i.e. the times where N increases by 1. In particular, for O S k S n. t == inf(t = k} < T. Therefore, for
0 S a b we can rewrite
(a, b)) = N(M
P(N(T)22)
We define by one of the following options:
as a constant Noo-
as a Binomial random variable, n Binom(p, nmmxx) for E (0,1), amar Nyo- as a Geometric random variable, n ~ (p) ( for E 1).
as a Poisson random variable, n ~ Poi(r), 0.
and ti....,tt are defined by choosing n uniform random variables on [0, T] and sorting them.
A.3 INFORMATION O ALGEBRA
In this section, we define mathematical tool, the o-algebra, that is essential to the description of the conditional expectation. This object describes which information is available at a certain time t. In the following, we leave away Ω whenever the meaning is clear. We define the filtration of the currently available information A : (Ar)sejo.77 by
At = or (X,,44 S(),
where t, are the observation times and (--) denotes the generated o-algebra. By the definition of T we have Ar = Ar(y) for all t E [0, Tj. (Ω x Ω. F @ F.F & F.P P) is the filtered product probability space which, intuitively speaking, combines the randomness of the stochastic process with the randomness of the observations. Here, FF & F consists of the tensor pooouctt o-algebras (F @ F), F, F fort E T). As explained in the remark below, A, can be identified with sub-o-algebra of (F F),
13

---

Published as a conference paper at ICLR 2021
Remark A.3. The a r-ffilds Ar depend on i Q as well. If we look the product probability space (Ω) Q2 F @ F.F @ F. x P), where F @ F consists of the tensor-product o-algebras (F F) = Fort 10,T), then
o (X,, st)
(offl E a ( E A, 2,1, ( AEE
is a well defined sub-o-algeebn of (F 00 F)+- Furthermore, we can recover the w-wise defined version of A by intersecting each set with x { and subsequently projecting the intersection on its first component. We use the notation A. == A(( =X(x,(aa) (x((5) Σ t) to distinguish this w-wise definition from the definition sub-o- -al-eeba of the product space given above. However, Lemma B.3 implies that for our considerations, both versions of this o-algebra have the same effect, therefore will simply write As for both versions, by abuse of notation.
B OPTIMAL APPROXIMATION X OF THE STOCHASTIC PROCESS X
We are interested in the "best" approximation (or prediction) X, of the process X that one can make any time & [0, T), given the currently available information A,. For us "best" refers to the ( 22 F @ F, P x (P)-minimizer, therefore, this approximation is given by the conditional expectation. Indeed, we define :=
== ((,,)) E [0, 7 + s T), and the function
x R mx ((t, ((t,r).))= + lu(t X |Xx E)
this is proven in the following proposition.
Proposition B.1. The optimal (i.e. L'-norm minimizing) aadapted process in F2,
P) approximating (Xe)ez 0,21 is given by² X == (Xt) with X, := x this process is unique to (IP P)-null-sets. In addition we have (w-wise for is E Q2) that
[ [X|]At). Moreover,
[" S.o (((t),,-(()),,,,)) x(t),x(()) ds.
= +
(12)
implying that X is càdlåg.
The first part the result follows from the elementary Proposition below (which is proven for example in (Durrett, 2010. Thm. 5.1.8) for Ravuue random variables and can casily be extended to Rd-valued random variables when using the 2-norm).
Proposition B.2. Given a probability space (2, F.P) and sub-o-algebra A C the orthogonal projection of X E (2, F, IP) on (2, is given by X == E[X|AA,, In particular, for every 7 E 1. (1, A, P) with P + X)
> O we have
E|X- n = E
X Z X E -
X12
Proof Proposition B.I. In our case this means, that the optimal A-adapted process L2( x á IF & F.P x P) approximating X is given by (Xe)1e(o.r) with X, = Eyx [A). Here, A,, is meant as sub-o-algebra of F, F. This process unique up to (IP Pn-null-sets. Moreover, the following lemma shows that it coincides w-wise with E [x, where A. A (w)) is defined in Remark A.3.
Lemma B.3. For IP-almost-every to E Ω we have Epx (x.AA.()) surely.
Ep [x(XA|(,(())) P-almost-
'hhile we give pointwise definition, (Cohen & Elliott, 2015, Theorem 7.6.5) allows to define X directly as the optional projection. By (Cohen Elliott, 2015, Remark 7.2.2) this implies that the process X is progressively measurable, in particular, jointly measurable in and w W. However, as we show below, even from the pointwise definition, it follows that X is càdlág, hence optional (Cohen & Elliott, 2015, Theorem 7.2.7).
14

---

Published as a conference paper at ICLR 2021
Proof. Otherwise we have by Fubini's theorem, Proposition B.2 and argument similar to the one in Lemma E.4
Exx Xx E (exx(Xx
Ep
Ep
Xx - Exx Exxe [x.[4.112]]
Ep
Ep
X; 1..1
which is contradiction to Proposition B.2.
This proves the first part of Proposition B.1. The second part of this Proposition, i.e. (12), should be understood ū-wise, for w E Ω. This is justified by Lemma B.3 and derived below. In particular, for the remainder of this section, all statements are meant W-wise.
With the assumption that and ø are Lipschitz, (Protter, 2005, Thm. 7, Chap. V) implies that unique solution of (9) exists. Furthermore, this solution a Markov process as soon as the starting point v is fixed (Protter, 2005, Thm. 32, Chap. V). Hence, one can define transition function P (compare (Protter, 2005, Chap. V.6)) such that for all < t and ø bounded and measurable,
(x, o) Eelo(X))
E2 E=(o(X,))Fa)
We have that Xx()) is - measurable and therefore, since AA, Ar(s) C. Fr(s))
(r(e),((x((e)))) (=(
(13)
By our additional assumption on it follows from (Protter, 2005, Lem. before Thm. 28, Chap. IV) that
M₁ Mp= = fo'ss ( X.), 0 OSIST,
is square integrable martingale, since the Browian motion W is square integrable. In particular, for < S t < T we have [s, (r, X,)dW,, [2[M, MF 0. Moreover, the same true when conditioning on Ag" 4.
Using the martingale property of M, we have (o-wise) for every t E [0, T)
X, = Epl(X, (()) + Xx(e))
= -X((1) + Ep l u(r, Xr Ar(t) + Ee le or, Xr,,, (r(1)
(14)
= =X+() + I F(1))
Ep [u(r, Xx)) dr,
where used Fubini's Theorem (for conditional expectations) in the last step. This is justified
because Ep[se [p(r, [2dr < oo follows from μ being bounded. Let us define A = ((1rr)
[0, 10,7731t 1 +r Т} and the function
Ā: A d Rdx
((1,r),)) (xett EE +r, ) = E)
then we can use (13) to rewrite (14)
X. X + I' S(())
ds.
(15)
This proves the second part of Proposition B.I. Proposition B.4. The function ji. is (jointly) continuous.
"To see this, we choose a localizing sequence )ne such that M is bounded by n (works since M
continuous). Then the Markov property implics that B[M," M. [A.] = E[M," - -M." [F.. = 0. Since M n M, Fass and since this sequence is dominated by the integrable random variable susup,,, MM (by Doob's inequality and square integrability of M), dominated convergence implies that T.. M, M, [A,]
15

---

Published as a conference paper at ICLR 2021
Proof of Proposition B.4. For any fixed [0,T], E Rdx, define
(x) [0, x Q - Rdx (t,w) + (w)
to be the solution of the SDE
(c)(5)= /. p(r, (a S: o(r, )dW,
This solution exists and is unique by (Protter. 2005, Chap. V, Thm. 7). therefore we have that
X, (1.) PP-almost surely. Furthermore, (Gubinelli, 2016, Thm. 4) implies that for 8 < I we have X, (G)))G).(())) Hence, for 1. = 8 + we have the identity E,(, Elu(t, (E))
Furthermore, by (Gubinelli, 2016, Thm. 8) we have for any E. E <R and (s. r), (s', r') with t == s str,t
+ r' E [0, that there exists some constant C' such that
(EG))(E) - (e (5) S C [ - 2'2 + (1 + 1sl2 + 15'12)² (l1 1') +1 (16)
Therefore, we have that
(,(,,r,)) Mis' r ,2) = E [utt, (,,) - B (5)
E [utt, (,,) - E (u(t', (E))l
+ |E lu(t', (o,(())) - E lu(t', (e))
$ E [la[t, (e,((E)) ((t, (E)) 1/2
(17)
+ E [lu(t', Sao(E)) - se
S E [lu(t, Se,e(E)) (E.()11 1/2
+ M E [IS,, - (E) 1/2
where we used Jensen's inequality in the second last and (10) in the last step. Hence, for (s', r',r'' - (s, g) we have that the first term of (17) goes to zero due continuity of μ in its first component uniformly in the second component. Moreover, the second term of (17) converges to zero by (16). Together, this proves continuity of j..
C RECALL: RNN, NEURAL ODE AND ODE-RNN EQUIVALENT WAYS OF WRITING
We describe our model as the solution of the SDE (29), which a compact way of writing. In the following, we first shortly recall recurrent neural networks (RNN) and the neural ODE and then recall the ODE-RNN. Furthermore, we give a step-by-step explanation, how the way ODE-RNN was formalized, can equivalently be written in terms of an SDE similar to (29). Finally give the alternative way of writing our model and a short comparison to ODE-RNN.
Recurrent Neural Network. The input to RNN is discrete time series of observations { Lt., }. At each time 4441, the latent variable is updated using the previous latent variable hi, and the input 2t,41 as
b = (NN,
(18)
where RNNCell neural network.
Neural Ordinary Differential Equation. Neural ODEs (Chen et al., 2018) are family of continuous time models defining a latent variable hç := h(t)) be the solution an ODE initial-value problem (IVP):
h = heo + [' f(h,, s, 0)ds, t t2 to.
(19)
where f(-..,) fe is neural network with weights 0. Therefore, the latent variables can be updated continuously by solving this ODE (19). We can emphasize the dependence of h, on numerical ODE solver by rewriting (19) as
h ODFSolve(fo. hoo (to, 1))
(20)
16

---

Published as a conference paper at ICLR 2021
ODE-RNN. ODE-RNN (Rubanova et al., 2019) is mixture of a RNN and neural ODE. In contrast to a standard RNN, we are not only interested in an output at the observation times t, but also in between those times. In particular,
want to have output stream that is generated continuously in time. This is achieved by using a neural ODE to model the latent dynamics between two observation times, i.e. for 4 < < 14+1 the latent variable is defined as in (19) and (20), with hig and to replaced by hi, and Eq. At the next observation time 4 the latent variable is then updated by RNN. Rubanova al. (2019) write this as
{ N MEa == ODESolve he: (4i, (i+1)) h44s == RNNCell( 1 P(i41)
(21)
Therefore, fixing hes hg. the entire latent process can be computed by iteratively solving ODE followed by applying a RNN.
GRU-ODE-Bayes. The model architecture describing the latent variable in GRU-ODE- (Brouwer et al., 2019) is defined as the ODE-RNN but with the special choice of a gated recurrent unit (GRU) for the RNN-cell and the neural network fo also being derived from GRU.
ODE-RNN as cådlåg process. Thinking about the process := (he) defined in (21) as a (stochastic) process in time, it is defined to evolve continuously for t, < 1,+1 according to the ODE dynamics So and jumps at time 4++1 according to the RNN cell. In particular, it is defined to be càdlåg process, for which h is the standard notation for the left limit, i.e. the last point before the jump at time f441. According to this notation we have fi h1 h/a hence, we can rewrite (21)
as
{
hare1- == ODESolve( fer hes (ts, 1(41)) hi:e1
(22)
RNNCell( (hr, Itaes).
C.I RESIDUAL ODE-RNN AS A SPECIAL CASE OF CONTROLLED ODE
Residual ODE-RNN. We replace the standard RNN cell by residual RNN cell (rRNN), as it was described e.g. in Yue et al. (2018). In particular, instead of applying the RNN cell such that hr, RNNCell(hr,,,
use a residual RNN cell have he, he,- rRNNCell(b,, The residual RNN is as expressive as the standard RNN and was empirically shown to perform very similarly or even better than the standard framework (Yue al., 2018). This way, the residual RNN cell models exactly the jump of the latent variable (i.e. the differences) that occurs at the time f++1 when taking into account the next observation T4,41 Therefore, we can rewrite the ODE-RNN (22) as
{ hG+1- he
ODESolve for hi, (tz, ((11) h + (RNN
(23)
Controlled Ordinary Differential Equation. We briefly recall the definition of controlled ODEs as it was given in (Herrera al., 2020, Section 4. 4.1) and used in (Cuchiero et al., 2019; Herrera et al., 2020) to describe neural networks. We fix E, d E N and define for 1 < i < d the vector fields
V θ x Rzo x Re Re (0,t,)) V x),,
which are càglàd in /. and Lipschitz continuous in æ. Furthermore, we define the scalar càdlàg control functions
u, R2o R, ut( which have finite variation and satisfy (=) = 0. Then we define the process Z := (Z4)220 as the solution of the controlled ODE
dz, Σ
Zo =
(24)
where 2 E Re. is some starting point. (24) written in Itô's differential notation for (stochastic) integrals. The solution of (24) exists and is unique under much more general assumptions than we made here (Protter, 2005, Chap. V. Thm. 7).
"i. right-continuous with existing left limits, also denoted as RCLL
i.c. left continuous with existing right limits
17

---

Published as a conference paper at ICLR 2021
Figure 7: Schematic representation of the stochastic control Re (25).
Special Case: ODE-RNN. To write the ODE-RNN as controlled ODE, we set = 2 and define VP So. 1/2 rRNNCell and in (4) == l. As pointed by Herrera al. (2020), one can take the to be semimartingales instead of deterministic functions. In line with this, we define = as the pure jump stochastic control process
" : x [0, T] + R, (wo, (4 ( = Σ o0)
(25)
We note that a an adapted process starting at 0 with finite variation on the product probability space (2 x Ω, F @0 F, @ F,P P x PP), since the total variation of u up to time T is n and Exxa [n] co. The following result shows that the residual ODE-RNN can compactly be described as a controlled ODE.
Proposition C.1. Using the vector fields and controls defined above, the latent variable process (hx)t2to of the residual ODE-RNN can equivalent written as the solution of the controlled
ODE
dhg = fo(h:-,t-)d ddt + :RNNCell((h_..:zddttt) hto = ho.
(26)
Proof of Proposition C.I. By (Protter, 2005, Chap. II. Thm. 17), the stochastic integral is indistin- guishable from the path-wise Stbeesue-Stieltjees integral if the integrator is finite variation. Hence, we can assume that some w E fixed and that the following expressions are evaluated at this whenever applicable. First note, that is constant except at the where it increases by 1 (ef. Figure 7). In particular, the Lebesgue-Stiele Stielltes integral of some càdlàg function with respect to u is a sum, EE,ge gt:-)) where Aut."
i.e. Ju
Utt 1. Therefore, integrating (26) We:
from to to È we get
he hes + [ fo s-)ds +
100
20
hto 1" fo 8 dds +
Σ (RNNCal(().,...)))
In particular, since the first integral is continuous in f, we have for every $ k: $
pla hes = =hro + f
fo(hs-, 8 s-)ds Ids+ + Σ (r)
(h (27)
-hh. TRNNColl(haa
and for to < < tx+1
h4 = hio + 1/ d2
(2 + 1' so (o s- s-)ds
(o 8- Σ
4.5
(28)
= - + S 8o
18

---

Published as a conference paper at ICLR 2021
Together, (27) and (28) prove that (26) is equivalent (23). We also emphasize that 2; is used as input only for = ti I n.
D NEURAL JUMP ODE
We propose a model framework whose architecture can be interpreted as simplification of the ODE- RNN (Rubanova et al., 2019) and the GRU-ODE-Bayes (Brouwer al., 2019) model architecture.
D.1 THE MODEL FRAMEWORK: NEURAL ODE BETWEEN JUMPS
We define X' Rdx and H. dRd be the observation and latent space for dx. du N. Moreover, we define the following feed-forward neural networks with sigmoid activation functions:
fer Rdn x Ix x [0, T] x [0, + Rd modelling the ODE dynamics,
Rd modelling the jumps when new observations are made, and 98, Rdu dy the readout map, mapping into the target space C RR for dy E N.
Rgxx
P0g
The trainable parameters of the neural networks are 0 (01, 02, 02) θ UM21 O Here, for every M E N, ①и is defined to be the set all parameters such that fo, and pe, have hidden dimension M. We define the pure jump stochastic process (cf. Appendix C)
(
u : 12 x [0, R. (w, + (0) = Σ oo)
Then the Neural Jump ODE (NJ-ODE) is defined by the latent process H (H.))ejo.T) and the output process Y = (Y)) defined as solutions of the SDE system (cf. Appendix C)
Hu p0a (Xo)
dH, fe. X
r(t) dt + (poz (Xi) He-)
(29)
Y, goa (H))
Note that, only the values of X at the times t,, 0 i used as inputs. We will use the notation H. (X) and Y," (X) to emphasize the dependence of the latent process H and the output process Y on the model parameters 0 and the input X.
Existence and uniqueness. We note that u an adapted process starting at 0 with finite variation on the product probability space (Ω x Ω & F,F & F, x PP, since the total variation of u up to time 7' is n and Epxa[n] oc. Hence, since fe, and pe, are Lipschitz continuous (as composition of Lipschitz continuous functions) unique càdlàg solution H exists, once an initial value is fixed (Protter, 2005, Thm. 7. Chap. V). Moreover, the resulting process (Y;")eejonn (e(0,, is also càdlàg and A-adapted.
Equivalent way of writing NJ-ODE. By applying the steps that were explained in Section C backwards to (29) (but without introducing the residual version of the RNN), is easy to see that our model equivalently be written (similar to (22))
{ h h
ODESolve(fag for (he, FG, t,,t tz), (t,, t++1) Po. (xtur))
(30)
In particular, the main difference (22) is that network performing the jumps does not take h as an input, i.e. it is not RNN.
modified input for fo, and that the neural
D.2 OBJECTIVE FUNCTION
We introduce a loss function, such that the output y of the NJ-ODE can be trained to approximate X, i.. to model the best online prediction of the stochastic process X. We define the theoretical loss function and its Monte Carlo and Ergodic approximation.
19

---

Published as a conference paper at ICLR 2021
Objective function. Let define D to be the set of all R -vaueed A-adapted processes the probability space (S x Q2, F F, F,P PP. Then we define our objective functions
:D R, @ : R,
Z J V Z = px
Σ (Xx, 1=1
Zu,2 + 12,, Z 21.-1222
(31)
0 - ( = v (X)
(32)
where @ will be our (theoretical) loss function. Remark that from the definition of Y it directly follows that it is an element of ID, hence @ is well-defined.
Monte Carlo approximation of the objective function. Let us assume, that observe N E N j S N, which are themselves independent realisations of the random vector (n,4, (.)) In particular, let us assume that
independent realisations of the path X at times (E( : ((ا)
P. are i.i.d. random processes (respectively variables) for 1 S j < N and that our training data is one realisation of them. Then
X X and (n² +(5) ,,())
x Y (x()))
+ +1xg, 2
y (xxo))))))
(33)
converges (IP Pa...s to (0) N oc, by the law of large numbers (cf. Theorem E.13).
Ergodic approximation of the objective function. we only observe one realization of the path X at times (f2.... fN) we still approximate the objective function by assuming that μ and are -nneeendeenndnnt and that the stochastic process X ergodic the following sense. We fix and assume that the time increments At, := := ビ i.idd realizations of the probability distribution A1 Furthermore, we consider each observation Xi, as one sample with initial condition Xi, for which yo. is the realization of y Then we approximate the objective function by
N - (x) + YE (X) - Y -x,,(())))) (x))))
Σ (iX N 1=1
(34)
which is assumed to converge by the ergodicity assumption for N + xi to
Exxp [ Y2 + Yr Y -1-
(35)
Instead of setting the random variable n = 1, one could similarly fix the time horizon T. and take for each sample all subsequent observations that lie in the time interval [ssta Istrtt + TJ. where lann is the date of the first observation of this sample. The next sample would then start with the first observation after ttartt +
E THEORETICAL CONVERGENCE RESULTS
Our main results show, that the output of NJ-ODE converges to the conditional expectation when the size of the neural networks and the number of samples go to infinity. For each network size (and number of samples) we assume to have the weights minimizing the (Monte Carlo approxima- tion of the) loss function. In particular, weights, and therefore also do not analyse backpropagation through our model. In a similar setting, backpropagation was studied by Jia & Benson (2019).
do not consider the problem of finding those optimal
For completeness we recall the definition of L cooovveegee
Definition E.1. Let < p co. Let (R,F, F) be probability space. Then sequence of random variables (X,)eeN converges to random variable X in (R. (or simply LP), if
B|X, - P
O.
We use the notation
Xn X.
I T.-convergeece is also denoted convergence in mean.
20

---

Published as a conference paper at ICLR 2021
E. CONVERGENCE WITH RESPECT TO THEORETICAL OBJECTIVE FUNCTION
Theorem E.2. Let 0N E NTT :- argmingeey. ( for every M E N. Then, for M - the value of the loss function t (32) converges the minimal value of (31) which is uniquely achieved by X. i.e.
((0m))
min (2) - (X))
ZED
Furthermore, for every 1 S k S K we have that Y converges X as random variable in (22 [o, T), P x Aτ). In particular, the limit process -im you equals X (P x Xα)- almost surely as a random variable on Q x [0, T)-
The idea of the proof is to split up the target jump process into its continuous-in-time parts and into the jumps. Both parts are continuous functions their inputs and can therefore be approximated by neural networks (Hornik et al., 1989). More precisely, by Proposition B. and B.4 the continuous-in- time part can be written as an integral over time of a function which is jointly continuous in its inputs. This jointly continuous function is approximated by neural network, which is itself integrated over time that is the neural ODE part of the NJ-ODE. The remainder of the proof shows T, -occvveeeeeece of the output of NJ-ODE with optimal weights (with respect the loss function) to the conditional expectation process.
For completeness restate the straight forward generalization of the universal approximation result (Hornik et al., 1989, Thm. 2.4) to multidimensional output.
Theorem E.3 (Hornik), Let d E N and a be sigmoid function. Let NN v,a be the set of all I-hidden layer neural networks mapping Rr to Réd. Then for every compact subset K C RR, every 6 > and every f E C(R", Rd) there exists a neural network g E NN r. such that supzek| g(x))2 <.
((R,RR)
The following Lemmas are used in the Proof of Theorem E.2.
Lemma E.4. Let 1 < k S K and let ID be a process such that (IP {x][A{ * Z] > 0. Then
there exists an e > such that B = {1 E [0, T] Br||X, Z, 2 a [| X, []] satisfies A (B) >
Proof. First remark that since X is continuous, we have X, = X, Let define : ((w, () E
Ω 10, T] 8,-(w) (w) + Z (w)) and for each & E [0, 7] let C. {w E Ω (w,t) t) C). Then
have for B == {t E [0, ((C.) > 0} that (BB) > O, since otherwise by Fubini's 's theorem
/ 10.71
< x
PPC.))da(()) =0
which is contradiction. Now Proposition B.2 yields that for each E B there exists x,
some > 0 such that [x[[X,- 7- 2 & + Er[[X,
2 This implies the claim. Indeed, assume no such > exists, then we have for each n E that
d (1 E 10. () Bp[[X,- 2, 2 2
E [|x &,- x O. Therefore,
A (t E [0,77 Epl|X,- 2 > Ee|xX,- 1111
Σ N (t [0, Er Z 2 N + EeI1X, - x.- = =0,
MEN
which is a contradiction to (B > 0.
Lemma E.5. For any l a--adppedd process Z it holds that
2
Epxi
1
Xr,
Ze, 2 = Epxi
1
X1,
x,
+
E
pxf
Σ
X
Ze.
2
i=1
i=1
1=1
21

---

Published as a conference paper at ICLR 2021
Proof. First recall that by continuity Xt, Xt.-. Then the statement is consequence of Proposition B.2, Lemma B.3 and Fubini's ss theorem, which imply
F
Xt.-
Ep
[x,- Zt:
7
2
i=1
2 [ Ze:-
Ex
Ex
X
+ IE>
[2]
X Z: 12
pxP
Xx:-
Xx-
Px
Lemma E.6. Let 1 2 p
< xc. Let (Za)aeN be sequence of random variables, and Z and Z
random variables defined on common probability space such that Zn LP and Zn Z. Then % almost surely.
LP,
Proof. With the triangle inequality follows that BIIZ 2///2/1 = 0, which implies the claim.
Lemma E.7. The random variable ST == supusisT Xel1 is square integrable and bounded in probability, i.e. for any exist some K O such that ST K] S.
Proof. From the proof of Proposition B.1 we know that
Mc = I. (() OstST,
is square integrable martingale and that Ep|suPoSeST MM S cEpll Mr < oc, for j 1.2 and some constant c by Doob's inequality. Moreover, μ is bounded, say by B, hence
x + ( ((r, S (r),W,)
Xx1|
$ | 1//
dr + 1M
$ | + BT + |M
Therefore, Ey[Sy] S claim. Indeed, if for a fixed >> 0 such K exists, then P[ST oo] 2 which is contradiction to integrability.
|=|1 + BT + cE M lı] < oo and similar Ee [e[[]] < co, which implics the
Proof of Theorem E.2. We start by showing that X E D is the unique minimizer of v up to (P x Ax)- null-sets for any k < K. First, we recall that for every 4 we have X., Xx, and by continuity of X that X = X Therefore,
2
v(X) ExxE
X
Xx:
2
1
Ep
Ep
Xx:
8c-- 2
71
min Zeb
Exp²
Lrc -Z.-12
n
min F Ex Zeb
(|Xr, Ze.1 + 12e, 2
min VZZ), ZED
22

---

Published as a conference paper at ICLR 2021
where we used Fubini's theorem for the second line, Proposition B.2 for the third line and the triangle inequality for the fourth line. Hence, X is minimizer of V. To see that it is unique (P x Aς) )
let Z E ID be process such that (P {x + Z] > 0. By Lemma E.4, this implies that there
exists an E > such that B == ( E [0, T Es |X Z [ 2 E + Br X x)1 satisfies (B) > 0.
Now recall that by definition of X we have (BB P > fn j,tx E B)) 2 k)> > O.
This implies that there exists j E N28 such that P(n = j, to- E B3 > 0. Therefore,
1 n
E
(1,
EB
2 RB
1 (m=s} n 1=1
/o Σ
1 EB-
n
i=1
2 1Eg [1=s]] 1 1
] (e0 [B-
= ! P(n = j,tx E B) > 0.
This inequality implies now that Z is not minimizer of L because
1
(2) = Exx#
(ix Zuil2 + +lZr. 222
n
2 Epp
Bp [xtr 11
(1(1-EBB) + 11 (()-)))))) Ep [|x., 21.-121
Epp
2 Ep
$ slu +Ep |xx_ X., 22
n
= E Ep.
1 (4+ EB}
+ min (Z) ZED
71
min ZED
(Z).
Next we show that (29) can approximate X arbitrarily well. Since the dimension dæ can be chosen
freely, let us fix it to du dx. Furthermore, let us fix 87 and 0 such that pa= goz id. From Theorem E.3 it follows that for any 0 there exist M E N and Of with (0i, 02 03) that
E O, such
foz (u,,,t,)) ā (t.r, Se.
sup
(36)
(V. M. xax X
where we used that f is continuous by Proposition B.4. Since μ is bounded, also ÿı is bounded, say by B- - > O. On [-M, Mfdx we approximate j by the neural network foy and outside we continuously extend for such that it is bounded by B. By abuse of notation we call this for again and use it as our neural network. Hence, [fe; 112 S 2B. Using this, we can bound the distance between y.lrr and X. In particular, if E (4,,... }. we have
| X, 2 = -1 + (Pez (X,) H, - X, = 0.
and if l not in (h,, ,lm). then (12), (36) and the previous bound yield for ST suPoseST |Xa]i
y,0's X, < 7(0) -x-0),
+ f' |f;; (H,_.Xr(s) T(S), - -(s) - i (T((s) T(8), Xr())) ds
-(0)
< $ 71 (STSM) + +2BT (S STT + 2BT1
23

---

Published as a conference paper at ICLR 2021
Here used that Sr M implies that (36) can be used for all r(t) S S t. Moreover, by equivalence of the 1- and 2-norm, there exists a constant c >> such that
y Xe/ S CE 2cBTT (S+>M).
With Lemma E.7 we know that Er (1 [5>> = P(s > M] =: CM the convergence of d ") using these two bounds and X,, Xr. Indeed,
M M 0. Now we can show
min V (Z) v ((0M - ((%rr)
ZED
E xp
[x.
-Y
yy y-03s 2)²
< E 7
(1 [x y
- 811
+ {x+ -xc-
-x1-11 + X yelie ))²
E 2x
(xx: xx:-1 + ceeT 2-7T 1222B119220 (11
EF (|x.-- -t,, -x,-xx,---
27
133 <
+ +c(TT 122-211111 (1r>)))))))))
n 11
27 1/2
27
< E2
Ep
Xt
X 4:
+ +Bp [CET + 2cBT 1eBT11ssn (8
2
where we used the triangle-inequality for the 1.2-norm in the last step. We can bound
Ep [ce + [
2cBTT (}? < 3(cs T>² + 3(cBBT) =: e2.
which is a constant converge to as M - sc. Using that for E R, we have a S a² 1, we get
min v (2) < ( < (0tr)
ZEB
27
1/2 1/2
27
Ep
Ep
X X1
+
(MM
11 i=1 1
2
Ep
1
Ep
Xr, 8,, -
12
[xe. - xce-
1/2
+ ohi
]
71
2
2
[=2
(1
2cM
T Exx
X xte-
[2]
2cM + c21
5( 2cm) min V (Z) 2cM + ch
M
min (Z)). ZEB
ZEE
where we used (X) = minzen (Z) and that eM converges to O.
In the last step we show that the limits limms you and limm- n exist as limits in the
Banach space L L: : (22 [0, T). F @ B(O, TJ), x AR). for every k < K, and that they are both equal to X. Let us fix k < K. First we note that for every B E B(10, 7J) we have
Ep 11 [22k] ((4-EB)]
Ex, [1]] (B) F(n2k.1--BB ((nk))
P(nn 2 k)
Using "measure theoretic induction" (Durrett, 2010, Case 1-4 of Proof of Thm. 1.6.9) this yields for (I -! and B(0, 7)) meeaurrbae function Z: [0, T
R. Z₁ == Z(t) that
Ex, [7] Egg [1 (x22) Zz,
(37)
24

---

Published as a conference paper at ICLR 2021
Moreover, the triangle inequality and Lemma E.5 yield
((0%r) - v(x)= - 2 Exx
X1, - V ((()
(5- 2
(38)
xxx
Σ
Xu- x
For any gddx-value Z E L the Hölder inequality, together with the fact that n 2 1, yields
1/2
Ex [
F
n
212
S pxx [m] 1/2 TE
Txx
|2
(39)
Together, this implies that lim have
y
X L-limit. Indeed, with ē == Ex
[n] oo
EpxAe
[ yous
= ccEpxt
1 (122))
-X,
2N
2
1 8ea-
y-
27 1/2
I c0 ccēEpxp
1 [n2b]]
1/2
1/2
S cē
1 (m2x)
Xr,
y-
1/2
12
22
S Scē
- 11
8u-
- h--
Seē (0(0))) v(x))
1/2
M M-+000 0,
where used first (37) and (39) followed by two simply upper bounds and (38) in the last step.
The same argument can be applied to show that limas-+oc y X L-limit. In particular,
this proves that the limit Y := lim y exists as L-limit and by Lemma E.6 it equals X (IP x d -almost surely, for any < K.
Remark E.8. This result can be extended to any other neural network architecture for which universal approximation theorem equivalent (Hornik al., 1989, Thm. 2.4) exists. Moreover, the stochustic process X defined in (9) can be chosen more general, in particular, the diffusion part od can be replaced by any martingale, as long as the resulting process still is a Markov process and jı stays continuous.
Remark E.9. we used the modified loss function & which is identical to V except that drop the factor 1.
n. everything would work similarly and we could show I2-convergence instead of L.-- convergence " to X. However, we remark that there might exists D, such that (Z) < co while (Z) = co. In particular, if some moment of n does not exist, such a process can be constructed.
Remark E.10. It follows directly from Theorem E.2 that lim you = as random variables on Q x [0, except on sets which are rull sets with respect every product measure IP A for 1 < K.
Remark E.11. The result of Theorem E.2 does not imply that X and limm y are modifi- cations or indistinguishable. For example, if C [0,T] ] is a subset such that no left-point of the observation times (tx-) lies in B with probability greater 0, i.e. Ax(B) for 1 S k < K, then Theorem E.2 does not tell us how close (lim y ) is to f. B. In particular, it does not tell us whether they are equal Plalmost surely. Furthermore, such set B always exists, since there has to exists I. E [0, T) such that B = {1} has measure 0 for all k.
In the following corollary show, that Theorem E.2 can be extended to show convergence to the conditional expectation of (X)). for some function & E CR RR) i.c. a function that is twice continuously differentiable with bounded derivatives.
25

---

Published as a conference paper at ICLR 2021
Corollary E.12. Let up E CR.N( R), then the statement of Theorem E.2 holds equivalently, when replacing X in the loss fianctions v and + by r (x() and X by the conditional expectation r, where I == [x))
Corollary E.12 combined with the monotone convergence theorem for conditional expectation theoretically enables us to make statements about the conditional law and conditional moments of X under some a priori integrability assumptions.
Proof of Corollary E. We first remark that Proposition B.2 and Lemma E.4, E.5, B.3 hold similarly for the conditional expectation f. Hence, the same argument as in the proof of Theorem E.2 implies that I' is the unique A-adapted minimizer of V.
For simplicity of the notation we assume that dx = dy = i.. that the process X and the Brownian motion W are 1-dimensional. However, the following works well in the general case, where the correlations of the Brownian motion components have to be taken into account. By Itô's Formula (Protter, 2005, Chap. I1, Thm. 32), r (xX) is the solution of the SDE
dxp(X): of (X:)dXx + (x,)
(X)(()) X,)dt + o (1, X,) (x)) a(t, 3(t, X,)dWe.
for ax(t, X.) o' (1)))((1, Xx) 3 (X.)o x,)² and (() &f (X,)o(1, Xt). Defining à
similar to as
à A Rd + Rdx, (t,r,E) + (X,,)(((, xx x,= = Er (a(t + X = E)
one can use the boundedness of of sp' and v" to show that it is continuous and that
i, = E(x(( A.] = F=(1) + [" ā (T((), - ((t) ds.
In particular, Proposition B. and B.4 hold equivalently for I. Similar to (36), the neural network parameters can be chosen such that
sup M. xdx A
[so ou,,, - à (t,r, Se,
(40)
which then implies the statement of the Corollary similar as in the proof of Theorem E.2.
E.2 CONVERGENCE OF THE MONTE CARLO APPROXIMATION
In the following, assume the size of the neural network M is fixed and study the convergence with respect to the number of samples N. Moreover, we show that both types of convergence can be combined. To do so, we define θ == {0 E θM |0] < M), which is compact subspace of OM. It is straight forward to see, that O in Theorem E.2 can be replaced by θы Indeed, if the needed neural network weights for an e-approximation have too large norm, then one can increase until it is sufficiently big. The following convergence analysis is based on (Lapeyre & Lelong, 2019, Chapter 4.3).
Theorem E.13. Let MJ E M,N == arg inf inf s ( for every M,N E N. Then, for every M E N. x P).a.s
N-000 +
uniformly on θM-
N
Moreover, for every M E E, (P x P)-a.s
No
((min)) and (O N
((Oii))
In particular, one define increasing sequence (NM)
MMMenn MEN in N such that for every S < K
M we have that y M s converges to X for M random variable in T. (0, T). PP Xx) In particular, the limit process? := lim equals (P x Almost surely as a random variable Ω [0, T).
gmin
26

---

Published as a conference paper at ICLR 2021
The following Monte Carlo convergence analysis is based (Lapeyre & Lelong, 2019, Section 4.3). In comparison to them, we do not need the additional assumptions that were essential (Lapeyre & Lelong, 2019, Section 4.3), i.c. that all minimizing neural network weights generate the same neural network output. This assumption not needed, because we do not aim to show that " converges to y
We define the separable Banach space 8 {a (zt))N E (1Exx) | 0o) with the norm ||a]]et == 2:12
E.2.1 CONVERGENCE OF OPTIMIZATION PROBLEMS
Consider a sequence of real valued functions (fn),. defined on a compact set K C. R Define, Vn inf (x) and let on be a sequence of minimizers fn (2,,) inf K fu (z).
From (Rubinstein & Shapiro, 1993, Theorem Al and discussion thereafter) we have the following Lemma.
Lemma E.14. Assume that the sequence (fn), converges uniformly on K continuous function f. Let v* = infaek f(z and S = {- E K f(x)= = v" Then U, + v² and ((TnS**) - 0 a.s.
The following lemma is a consequence of (Ledoux & Talagrand, 1991, Corollary 7.10) and (Rubinstein & Shapiro, 1993, Lemma A1).
Lemma E.15. Let (Ei)121 be a sequence i.i.d random variables with values in S and h x - R be measurable function. Assume that a.s., the function O E R& + 1(0,E1 (0) is continuous and for all C O, E(supjo)asc [h(0, Ei)l) < +oc. Then, a.s. 0 E Rd + N E h(0,E.) converges locally uniformly to the continuous function E Rd + (E(h(0,E))
Σ M(0 E(h(0,E:)) =0 N 10
lim sup
N
E.2.2 STRONG LAW OF LARGE NUMBERS Let us define
((x, 512 by =/2
and Eg = (X)
0, O,...). ). where X2 are random variables describing the realizations of the 1
training data, as defined in Section D.2. By this definition we have n' := n (E)) = max,NN(E, * 0} P-almost-surely and we know that E are i.i.d. random variables taking values in S. Furthermore, let us write (E) make the dependence of Y on the input and the weight 0 explicit. Then define
n
(2 (6))
h(0, E2) = TÚ
Lemma E.16. The fhelowwingg properties are satisfied.
(Pi) There exists >> 0 such that for all S E and A E Oss we have (52 & <1+1)(x(() for all E [0,
(P2) Almost-surely the random function E θM + Y, uniformly continuous for every e [0, T..
(Pa) We have Exxi 1 EE 1x < oo and Exxi E Xx,- oc.
Proof. By definition of the neural networks with sigmoid activation functions (in particular having bounded outputs), all neural network outputs are bounded terms the norm of the network weights, which is assumed to be bounded, not depending on the norm of the input. Since after jump at 7 Y has the value Xx(e), we can find depending on such that the claimed bound satisfied for all proving (Pt).
27

---

Published as a conference paper at ICLR 2021
Since the activation functions continuous, also the neural networks continuous with respect to their weights 0, which implies that also 8 θ Y, is continuous. Since Oas is compact, this automatically yields uniform continuity and therefore finishes the proof of (P2).
(Pa) follows directly from the stronger result in Lemma E.7.
Proof of Theorem E. 13. We apply Lemma E.15 to the sequence of i.i.d random function h(0,5,)) From (P1) of Lemma E. 1616 we have that
F(X+) 1/ Y v (xx2 - yol
1/12 12 ve 6 1 + 1
-1
s 4 (xx,) 1
- 4 (31X,/1++ x(1
X x (3)
Hence, we obtain that
n
n2
h(0,E,) n Σ
F(X, 2 Yg, x3
12 + 4
x + 4x + 1
1=1
implying that
Exxt
sup h(0,E) #
< (12 4x) Exxt
Σ i=1
X 12
+ 4x + 10l3 < <0,
(41)
using (Ps) of Lemma E.16. By (P2) of Lemma E. 16, the function 0 + h()) is continuous. Therefore, we can apply Lemma E.15, yielding that almost-sured for N + co the function
N
h(0,E,) = (x(
(42)
N
converges uniformly on θM to
0 Exxe((((0,,)))) =
(43)
We deduce from Lemma E. 14 that d(Omin O - 0 a.s. when N + 00. Then there exists
sequence (Ómin) NEN in Owin such that 10miN 0 - 0 a.s. for N + Jc. The uniform
continuity of the random functions f + Y on θas implics that [Y, M.N. Y,
l2 + a.s. when
N + oo for all E [0, T'. By continuity of I' this yields In(Omii E1) - h(mii
M,N
E1)12 - 0 a.s. as
N + co. With (41) we can apply dominated convergence which yields
[In(emii [ E1) h(@min E1)12 0.
lim Exx +000
Since for every integrable random variable we have < [||2|2 S E[] and since O Omin we can deduce
lim ((mmn))) lim Epxē [n(OiiN Et) N
lim Exp ho Ea)] ((
(44)
-
No
Now by triangle inequality,
N (Omin)> ((0mm)) < 1& (6min.
(6
+
(
( )
(45)
(42). (43) and (44) imply that both terms on the right hand side converge to when N - co. which finishes the proof of the first part of the Theorem.
We define No = 0 and for every M E N
NM = min { N EN N > NA ( (0 < M + ( - v(X)}}
28

---

Published as a conference paper at ICLR 2021
which is possibly due to (44). Then Theorem E.2 implies that
+ (0M (X)) < मे + 2| ((X) M-000 0.
Therefore, we can apply the same arguments as in the proof of Theorem E.2 (starting from (38)) to show that
ExxA
||
y
guain
[]] cca (D(ON)NN))) -V(()) 1/2 M M- O,
for every 1 K.
Corollary E.17. In the setting of Theorem E.13, we also have that (IP x a...s
((0 M- (X) and t Nu
(Omin M.N
M-+00 v(X)
MEN where (NMs MEN is another increasing sequence in N.
Proof. The first convergence result was already shown in the proof of Theorem E. 13 and the second one can be shown similarly, when defining NM by No == 0 and for every M E N
1 NM == min { E ENN/ N NM 14s-1,101 (ONiN)) ((0ni)) s M + 1(311
which is possibly due to (45).
E.3 DISCUSSION ABOUT OPTIMAL WEIGHTS
In Theorem E.2 and E.13, the focus lies on the convergence analysis under the assumption that optimal weights are found. Below we discuss, why this assumption is not restrictive in theory.
Global versus local optima. The assumption that the optimal weights are found, is typical for convergence analysis of neural network based algorithm, since the objective function is highly complex and non-convex with respect to the weights. particular, it is well known that the standard choice of (stochastic) gradient descent optimization methods do in general only find local and not global minima. Since the difference between any local minimum and the global minimum can not generally be bounded, is unrealistic to hope for theoretical proof of convergence with respect to such optimisation schemes. On the other hand, global optimization methods as for example simulated annealing provably convergence (in probability) to a global optimum (Locatelli, 2000; Lecchini- Visintini et al., 2008). Hence, combining those with our result, convergence in probability of our model output the conditional expectation can be established, without the assumption that the optimal weights are found. However, these global optimization schemes come at the cost of much slower training compared to (stochastic) gradient descent methods when applied in practice. Moreover, several works have focused on showing that most local optima of neural networks are nearly global, see for example (Feizi et al., 2017) and the related work therein. Hence, using (stochastic) gradient descent optimization methods likely yield nearly globally optimal weights much more efficiently. In our case, this is also supported by our empirical convergence studies Section 6.3.
29

---

Published as a conference paper at ICLR 2021
F EXPERIMENTAL DETAILS
All implementations were done using PyTorch. The code is available at httpp:://gttub... Cerrerakrrechheee.cmmm/////
F.1 IMPLEMENTATION DETAILS
Dataset. For each of the SDE models (Black-Scholes, Frnstein-Uhlenbedd Heston) a dataset was generated by sampling N 20'000 paths of the SDE using the Euler-scheme. We used equidistant time grid of mesh 0.01 between time 0 and T = 1. Independently for each path, observation times were sampled from Pt. by using each of the grid points with probability 0.1 an observation time. In particular, on the time grid. Hence, 10% of the grid points were used on average. This way, n and f; are defined a discretized version of those given in Example A.2 where is binomially distributed and t; chosen uniformly on [0, T For each of these datasets, the samples were used in a 80% /20% split for training and testing. The SDEs of the dataset models and the chosen parameters are described below.
Bin(1 000, 0.1). lo 0 and the observation times {4: {4isssenn were chosen uniformly
Black-Scholes:
SDE: dX, pX,dl a x,dW, where W is I-dimensional Brownian motion conditional expectation: E(X++e|XX (x)) = Xyeus
used parameters: μ 0 = 0.3, Xo 11
Orastein-Uhlenbbekk:
SDE: dX,- - m)dt + odW,, where W is I-dimensional Brownian motion conditional expectation: FE X,) X,e m (1 ehha) used parameters: k 2,
0.3, Xu =1
Heston:
SDE: for and Z I-dimensional Brownian motions
dX; pX,dt + Vo,W, du; -k(u: m)dt + vvdde
conditional expectation': E(X++>|XX ) X,eus used parameters: μ 0.3, Xo 1 k 2, m 4, Uo 4, p= Corr(W, Z) 0.5
Architecture. In our experiments we choose the dimension of the latent variable to be dæ = 10. For fe., Ge, and مe_ we use 2-hidden-layer feed-forward neural networks, with 50 nodes in each hidden layer and tanh activation functions. Then the neural networks go, and po, are defined as residual versions of ge, and Po,- by adding a residual shorteut between the input and the output of the neural networks. Dropout was applied after each non-linearity with rate of 0.1. To scale the possibly unbounded inputs the neural networks to bounded hypercube, we applied tanh component to and h in every neural network. This was done, because neural networks sometimes become unstable when their inputs become large. To solve the ODE of the neural ODE part, the simple Euler-method was used.
Training. The neural networks were trained using the Adam optimizer (Kingma & Ba, 2014) with learning rate of 0.001 and weight decay 0.0005 for 200 epochs using batch size of 200. A random initialization was used and no hyper- parameter optimization was needed.
Further training results. Further training results test samples shown in Figure 8, 9, 10.
F.2 EXPERIMENTS ON OTHER DATASETS
We test our framework on additional synthetic datasets.
see (Rujivan& & Zhu, 2012, Equation 2.9)
30

---

Published as a conference paper at ICLR 2021
true path
true path
true conditional expertation observed
true conditional expertation
observed
2
1
0.0
0.2
0.0
0.2
Figure 8: Black-Scholes
true path
.....
true conditional expectation
3.0
observed
3.0
2.5
2,0
2.0
path
1.5
model
expectation
1.0
observed
1.0
0.2
0.4
0.6
0.2
0.4
0.6
Figure 9: Ornstein-Uhlenlee
true path
true conditional expartation
observed
path our model
expectation
2-
observed
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.4
0.6
0.8
1.0
Figure 10: Heston
31

---

Published as a conference paper at ICLR 2021
F.2.1 HESTON MODEL WITHOUT FELLER CONDITION
If the Feller condition
2 2km 2 o²
satisfied the Heston model, it is known that the variance process tt is always strictly bigger than 0 and that the Euler-scheme works well to sample from it. However, if the Feller condition not satisfied, the variance process can touch 0, where the process is reflected deterministically. Then the Euler-scheme can not longer be used to sample from the model. Moreover, this situation generally considered as more delicate. Close to the distribution of v behaves differently than further away of 0 (Andersen, 2007, Section 3). As explained by Andersen (2007), there differently well performing sampling schemes for Heston model where the Feller condition not satisfied. We the simplest implement, which is slight extension of the Fuler scheme, where values of e below O) are replaced by 0 (Andersen, 2007, Section 2.3). Although there is empirical evidence, that the resulting sampling distribution of 1 close is not correctly replicating the true distribution (Jean-François et al., 2015, Figure 2,3), this method already produces sufficiently good sample paths.
Dataset. Heston model, sampled as in Section F.I with the extension described above. Used parameters: μ 2, a = 3, Xo k = 2, m = 1, Uo = 0.5, p Corr (W, Z) 0.5. Hence the Feller condition and also the weaker condition 4 m 2 x² discussed (Jean- -raan-François et al., 2015, Section 3.2). are both not satisfied. We produce two datasets, 1-dimensional one similar to before, where only X stored and 2-dimensional one, where both X and e are stored, hence also v is target for prediction. Note that has the same conditional expectation as the SDE, the 2-dimensional dataset, and v are always observed at the same time.
Architecture & Training. Same as in Section F.1, but with batch size 100.
Results. The model learns to replicate the true conditional expectation process, which analytic, hence not effected by the sampling scheme. In particular, we see that our model is very robust, since even the delicate case where the Feller condition is not satisfied and a sampling scheme is used, that does not perfectly replicate the true distribution, our model still works well. Due to the very similar results, in Figure we only show plots on test samples of the 2-dimensional dataset, where X and are predicted.
30-
0.2
0.6
0.9
0.0
0.2
0.6
10
1.5
path truc expectation observed
true conditional expectation
observed
AA
0.5
0.0-
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.6
0.8
1.0
Figure 11: Heston model without Feller condition. In both plots, the upper sub-plot corresponds to the 1-dimensional path of X, and the lower sub-plot corresponds to the 1-dimensional path of et.
F.2.2 DATASET WITH CHANGING REGIME
Dataset. To be able to evaluate the performance of our model, we test change of regime by combining two synthetic datasets. the first half of the time interval [0, 0.5) we use the Ornstein- Uhlenbeck and on the second half [0.5. 1] the Black-Scholes model. The Black-Schole process takes as starting point, the last point the Ornstein-Uhlenecck process. We use the same hyper- parameters for the dataset generation in Section F. except that we set the parameter the Monsttin-UUheebeckk model to make the two parts act on similar scales.
10
Architecture & Training. Same in Section F.1, but with batch size 100. Moreover, used 100 neurons each hidden layer, to account for the more complicated setting, where also a time dependence has to be learnt.
32

---

Published as a conference paper at ICLR 2021
Results. In the plots test samples of the dataset shown in Figure 12 we that our model correctly learns the change of regime.
trus patt
Tue 2eth
ou model
Our motel
23
oxpeddtiinn
20
2sccondiond copeciatio
10
Figure 12: Our model evaluated on a stochastic dataset that follows an Censtein-Uhlenbee SDE on time interval [0, 0.5) and an Black-Scholes model on the time interval (0.5, 1].
F.2.3 DATASET WITH EXPLICIT TIME DEPENDENCE
Dataset. We use the Black-Scholes datasets of Section F.1, where we replace the constant y by the time dependent constant _ (sin(3t (t) + 1). for a, 3 > 0. The conditional expectation changes accordingly. For the data generation we use the same typer-parameters in Section F.1, except that we use instead of the parameter a = and 3 E {2m, 4m)
Architecture & Training. Same as in Section F. but with batch size 100. Moreover, used 400 neurons in each hidden layer, to account for the more complicated setting, where an explicit time dependence has to be learnt.
Results. We show plots on test samples of the datasets in Figure 5. We see that the model learns to adapt the time-dependent coefficients.
F.3 DETAILS ON CONVERGENCE STUDY
Evaluation metric. For the sampling time grid with equidistant step size A, = 3. EEVE N, on [0, 7] and the true and predicted conditional expectation for path j N, X and y respectively, we define the evaluation metric as
N₂
cval(X),) = N2
1
Σ
Σ (2). - -YZ,)² 1=0
(46)
where N₂ is the number of test samples, and accordingly iterates over the paths in the test set.
Increasingly big training sets. We use the following procedure to create increasingly big training sets, while keeping the exactly same test set for evaluation. Out of the initial 20'000 paths, take N2 == 4' 4'000 paths, which are fixed as the testing set. Out of the remaining 16'000 paths, we randomly choose N₁ training paths for N1 E (200, 400, 800, 1'600, 3'200, 6' 400, 12'8003.
Increasing neural network sizes. The increasingly big neural networks are defined as follows. For all involved networks, we use the feed-forwant 2-layer architecture with tanh activations (cf. Appedix F.1) where each hidden layer has the same size M for M E {10, 40, 80, 160, 320}.
Results on Black-Seli Black-Scholes dataset. accordance with the theoretical results Theorem E.2 and E.13, we see that the evaluation metric decreases when N1 and M increase (Figure 13). It is important to notice, that already a quite small number of samples can lead to a good approximation of the conditional expectation, if the network is big enough. In particular, the Monte Carlo approximation of the theoretical loss function is good already with a few samples. On the other hand, even with a large number of samples, the evaluation metric does not become small, if the network size is not big enough. From a practitioners point of view, this good news, since increasing the network size often much easier than collecting more training data.
33

---

Published as a conference paper at ICLR 2021
10"1
retwork 50c-10
10
baining sice-200
retwork
size-4co
retwork 5P2-40
haining 8078-800
retmurk pice-20
buiii y pice-1000
retwork sis=11e
2
triining size
eval metric
evall metric
training sioum12000
10
10-3
101
10'
training
netoork
Figure 13: Black-Schole dataset. Mean = standard deviation (black bars) of the evaluation metric for varying N1 and M.
Results on
and Heston dataset. We get very similar results on the Ornstein
Uhlenbeck (Figure 14) and Heston dataset (Figure 15). Similar for Black Scholes, also for
Ornstein-Uhlenbeck the training size N₁ not as important as the network size M. For network sizes, increasing the training size further than 1600 hardly changes the performance, while increasing M is crucial to get better performance. In contrast to this, for the Heston dataset we see that large
s0v-10
haining siza-100
retwork zice=20 I retwork sireman
tramningsree pire-4co
training sizemtce
sire--0
lowning traninggtieettt0 1200
T retwork
relsork_seee330 s24-320
T traning trannngsteeeaae sizeme4oo siw---iii--nnnn
mome
eval mete
10
10"4
training_size
networkk_se
Figure 14: Orastein- Uhlenbeck dataset. Mean + standard deviation (black bars) of the evaluation metric for varying N, and M.
number of training samples is more important get smaller convergence metric. This reflects the fact, that the Heston dataset is more complex and more difficult to learn.
F.4 DETAILS ON COMPARISON TO GRU-ODE -AYES
GRU-ODE-Bard (Brouwer et al., 2019). To the best of our knowledge, this the neural network based method with the most similar task to ours. In particular, this continuous time model is trained to learn the unknown temporal parameters distribution of X given the previous observations. This distribution given by the Bokker-Plance equation. Brouwer et al. (2019) outlined, that their model can exactly represent the Fokker-Planck dynamics of the -hbseeinUhlenneck process, since the corresponding distribution is Gaussian.
normal distribution, best describing the conditional
Implementation of BRU-ODE-Bayes. We use the code of the official implementation of (Brouwer al., 2019) and slightly adjust it for purpose. In particular, we do servations, hence the input mask used for this task has always only 1-entries. Furthermore, we slightly changed the scheme how the time steps are taken, be the same as in our implementation. so that comparisons can be made. Besides these minor changes, the original model is used and
incomplete ob-
hhtpps://e :/ttpa////ii thub. wtt.:/cc//ceaeeeeeeee hwwbb..cm/deebrouwe//ggggt
34

---

Published as a conference paper at ICLR 2021
4x10
4x100
2X10
3X10P
eral poatmmnne
triniin_s20-200
network o2e=10
trainnggaze 400
network 200
hriningg s0
training buiniing
retwork sice30
Detning stcw2202
network Sce-160
Beliingg 52e-6403
network 320-120
training -2930
10'
Figure 15: Heston dataset. Mean + standard deviation (black bars) of the evaluation metric for varying N1 and M.
trained on all 3 datasets (Black-Scholess Ornstein-Uhlenbeck and Heston). We tried out all com- binations of the following parameters and always chose the best performing one for comparison to our model: impute E {'True, False), logvar E {TTuu, False). mixing E {0.0001, 0.5} and hiddensize E (50, 100), whereby phidden and prephidden were chosen be equal hiddensize. The first parameter choices were the ones used in the official implementation. The model was always trained using the Euler-men mulhod for solving ODEs and with dropout 0.1, to be comparable to our implementation. Furthermore, we always used the full GRU-ODE GRU-ODE-Bayes implementation, since should the more powerful one. For the comparison, we only use the estimated mean of the normal distribution (the estimated variance is not used), which is precisely the estimated conditional expectation the model.
Implementation of NJ-ODE. Same as described in Section F.1. For each dataset only one model was trained, since we already saw the convergence properties before and wanted to have qualitative comparison GRU-ODE BooeBB. In particular, we did not try to optimize best performance (e.g. by choosing different number of layers or neurons or different activation functions) and used about 10K trainable parameters compared the best performing GRU-ODE- Bayes models of our study which used 112K trainable parameters.
hyper- prper-prrametess for
Datasets. Same described in Section F.1. The train and test sets were fixed be the same for all trained models.
Training. All models were trained using the Adam optimizer (Kingma & Ba, 2014) with learning rate of 0.001 and weight decay 0.0005 for 100 epochs using a batch size of 20. A random initialization was used.
F.4.1 FURTHER RESULTS AND DISCUSSION OF THE COMPARISON
Evolution of losses and evaluation metric during training. In Figure 16 and 17 we see, that already after few epochs the NJ-ODE model finds close to optimal weights for the given network size on the Black-Scholes and Cehtein-hhlenbeck dataset and oscillates around this optimum. The Heston dataset is considerably more difficult and the model slowly converges close-to-optimal weights for the given network size. In comparison to this, it takes the GRU-ODE-Baye model longer converge to its close-to-on tosseimal weights on the Black-Scholes and Orstein- Uhieebbeb dataset. Moreover, the model does not converge to close-to-optimes weights on the Heston dataset, but rather oscillates between bad weights. Due to some very large outliers, this is not directly visible the plot, but can be deduced from Table 6.4.
Comparison of predicted paths. For each dataset we show paths that were predicted with NJ- ODE and with GRU-ODE-Bayee Bayes, first at the optimal epoch, i.e. where the test loss was minimal during training (Figures 18, 19, 20), and then at the last epoch (Figures 21, 22, 23). For the sake of comparison, in each the performance similar for NJ-ODE and GRU-ODE GRU-ODE-Bayes Bayes on the Sacck-s and Uhieebbee dataset, with sometimes the and sometimes the other having slightly better prediction. On the Heston
the same test sample is shown. The results very
35

---

Published as a conference paper at ICLR 2021
M Mannh Mo ManMM Man
oasco
Worh Maaeh won My
m
Figure 16: NJ-ODE on Black-Scholes, Ornstein-Uhlenbeck and Heston (from left to right). Blue (1st row): training loss, orange (2nd row): evaluation loss, green (3rd row): evaluation metric.
km
M. M Mado
680
Figure 17: GRU-ODE-Bar Bayes best performing models on Black-Scholes, Omstein-Uhlenberr Uhlenbeck and Heston (from left to right). Blue (Ist row): training loss, orange (2nd row): evaluation loss, green (3rd row): evaluation metric.
dataset, the results of NJ-ODE are good predictions, being correct whenever a new observation made and not too far away even after longer periods without observations. On the other hand, GRU-ODE-BayeS is not even correct at times of new observations and learns an incorrect behaviour in between observations (e.g. making kinks where there should not be any). For the last epoch, see that this malfunction is amplified.
F.4.2 FURTHER EXPERIMENTS
Since the results of GRU-ODE-Bayee additional tests. In particular, we retrained the best combinations of GRU-ODE-Bayes with the larger batch size 100. For smaller versions of the network (M), i.e. with hiddensize prephidden phidden = 50 this stabilized the training, such that the models converged. Still, even at the best epoch, the models suffered from the same difficulties as explained in Section F.4.1. Increasing the hiddensize by factor to 100 (L), made the training unstable again, and the models did not converge any more. This gives more empirical evidence, that GRU-ODE- -Byes can not be reliably trained on more complex datasets, where the target conditional distributions differs to much from normal distribution. In contrast to this, retraining NJ-ODE with batch size 100, once for the smaller version described in Section F (S), once for larger version with 100 instead of 50 neurons in all hidden layers (M) and once with 200 neurons (L), yielded the expected results of better performance with larger networks. In particular, there is instability for larger networks. In Table F.4.2 show results of our model and of the best GRU-ODE-model for the given size.
unexpectedly bad the Heston dataset, performed
Table 4: The minimal, last and average value of the evaluation metric (smaller is better) the Heston dataset throughout the 100 epochs training, together with the number of trainable parameters.
min
last
average
params
GRU (M)
2.38
4.63
8.25
28'802
GRU (L)
5.58
17.39
1877.26
112'602
ours (S)
1.71
1.71
2.05
10'071
ours (M)
1.21
1.23
1.60
35'121
ours (L)
1.07
1.08
1.36
130/221
36

---

Published as a conference paper at ICLR 2021
F.5 REAL WORLD DATASETS WITH INCOMPLETE OBSERVATIONS
F.5.1 LOSS FUNCTION FOR INCOMPLETE OBSERVATIONS
In accordance with the self-im reff-imputation scheme, the loss function is adjusted to only use the non- imputed coordinates, i.e. if m is the random process in {0, 1 describing which coordinates are observed we have
(Z) = Epxt
Σ (mm, (X Z1,)22 Im, (Z,, Zt (122)
F.5.2 CLIMATE FORECASTING DETAILS
Dataset. We use the publicly available United State Historical Climatology Network (USHCN) daily dataset (Menne et al., 2016) together with all pre-processing steps as they were provided by Brouwer al. (2019). In particular, there are sporadically observed (i.e. incomplete observations) climate variables (daily temperatures, precipitation, and snow) measured at 1'114 stations scattered over the United States during an observation window of 4 years (between 1996 and 2000) where each station has an average of 346 observations over those 4 years. For 5 folds, the data is split into train (70%), validation (20%) and test (10%) sets. The task is to predict the next measurements after the first 3 years of observation. The mean squared error between the prediction and the correct values is computed on the validation and test set.
Baselines. We use the results reported in (Brouwer et al., 2019, Table 1) as baselines for our comparison and perform the exact same 5-fold cross validation using the same folds with the same train, validation and test sets. For completeness, give all results of the table together with our results in Table 6.5, We only show the mean squared error (MSE), since our model does not provide the negative log-likelihood.
Implementation of NJ-ODE. We once use the architecture described in Section F.1 (S) and once use the exact same architecture, but with hidden size dи 50 and 400 instead of 50 nodes in each hidden layer (L). To deal with the incomplete observations, the self-imputation scheme described in Section 6.5 is used. In particular, the data is self-imputed and passed together with the observation mask as input to the network Pay. Moreover, the network fe, uses the NJ-ODE output at the last observation time instead of the last observation as input.
Training. Both versions of NJ-ODE were trained using the Adam optimizer (Kingma & Ba, 2014) with learning rate of 0.001 and weight decay 0.0005 for 100 epochs using batch size of 100. A random initialization was used. The performance on the validation set was used for early stopping after the first 100 epochs, i.e. early stopping was possible at any epoch between {101,.... 200].
Results. The results of our model and all models reported in (Brouwer et al., 2019, Table 1) are shown in Table 6.5.
F.5.3 PHYSIONET PREDICTION DETAILS
Dataset. We use the publicly available PhysioNet Challenge 2012 dataset (Goldberger et al., 2000) together with all pre- pe-processing steps as they were described and provided by Rubanova et al. (2019). In partieular, there are 41 features of 8000 patients that are observed irregularly over 48 hour time period. The observations are put on a time grid with step size 0.016 hours leading to 3000 grid points. While in their paper Rubanova al. (2019) say that they use 2880 grid points (i.e. minute wise) in their implementation they used 3000. Moreover, in contrast to what was written in the paper, the constant features were not excluded in their implementation, hence we also keep them. Furthermore. we also rescale the time-grid to [0, 1] and normalize each feature as Rubanova et al. (2019) did it for training the model. The dataset split with the same fixed seed into 80% training and 20% test set. In particular, no cross validation but only multiple runs with new random initializations are performed to be exactly comparable to the results reported by Rubanova et al. (2019). On the test set, the observation paths are split in half. The first half (first 1500 time steps) is used as input for the model, from which the second half (second 1500 time steps) should be forecast. The masked mean squared error (MSE) with certain balancing procedure between the predictions and the true values is reported. For comparability we use the exact same implementation of the MSE as Rubanova et al. (2019).
37

---

Published as a conference paper at ICLR 2021
Baselines. We compare the performance of our model latent ODE on the extrapolation task (as described above) on physionet. As baseline for our comparison we use the results reported in (Rubanova et al., 2019, Table 5). For completeness, together with our results in Table 6.5. We shortly outline the different approach of latent ODE compared to model for the given extrapolation task. Latent ODE also splits the training samples similar to the test samples in half, using the first half as input and the second half as target. It is trained as encoder-decoder, encoding the observations in the first half and reconstructing (decoding) the second half. This falls in the standard supervised learning framework. In particular, this approach cannot straight forward be extended for online forecasting. Moreover, this approach might learn certain path dependencies. On the other hand, our model is trained always, online forecasting after each observation until the next observation made. Instead of splitting the training samples we use the entire path as input for our unsupervised training framework. Our model is based on the assumption that paths are Markov, therefore cannot learn path dependencies, i.c. dependencies on more than just the last observation. However, by training the model also on the second half' of the training samples, it learns the underlying behaviour there, which should be helpful for the extrapolation task.
show all extrapolation results of the table
Implementation of NJ-ODE. We use the architecture described in Section F.I, but with hidden size dn 41. To deal with the incomplete observations, the self-imputation scheme described in Section 6.5 is used. In particular, the data is self-imputed and passed together with the observation mask as input to the network Pey. Moreover, the network fer uses the NJ-ODE output at the last observation time instead of the last observation as input.
Training. The NJ-ODE was trained using the Adam optimizer (Kingma & Ba, 2014) with learning rate of 0.001 and weight decay 0.0005 for 175 epochs using batch size of 50. 5 runs with random initialization were performed over which the mean and standard deviations were calculated. In particular, these runs always used the same training and test set, specified by the same random seed as in thee implementation of Rubanova et al. (2019).
Results. The results of our model and all models reported in (Rubanova et al., 2019, Table 5, extrapolation) are shown in Table 6.5. We report the minimal MSE on the test set during the 175 epochs, since it was not differently specified in (Rubanova et al., 2019). However, if the MSE of the epoch is used where the training loss minimal, the results are nearly the same with 1.986 + 0.058 (x (x10-3). Moreover, of 50 (187'323 parameters in total), leading to slightly better results of 1.934 + 0.007 (x10 (at minimal MSE) and 1.982 + 0.027 (x 10 10) (at minimal training loss).
trained a larger model for 120 epochs, where 200 nodes were used instead
38

---

Published as a conference paper at ICLR 2021
True our modal
Tree path
GRU 00E Bayes
true conditoral exacctation observed
I.:: true candit ona expectation obsetved
0.0
0.2
0.4
0.5
0.8
2.0
0.0
0.2
0.4
0.5
true
Uree
model
GRU CDE Bayes
tne conditordl oxectation
tue candit ond expectation otserved
atserved
0.2
0.2
true
Srue pats Dayes true condit onal expectation absetved
our model
true conditional expectation abserved
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.4
0.6
5.8
1.0
true path
Irue path
oir model
GRU-dDr--Bayws
tne conditoral oxeetation
the candit onal expedaatio abserved
abserved
Figure 18: Comparison of predictions for Black-Scholes paths at best epoch.
39

---

Published as a conference paper at ICLR 2021
35
3,5
Suyca
expectatio
expectaiion
3.0
2.5
1.0
1.0
3.0
0.2
0.6
6.6
0.8
3.0
9.0
0.2
0.4
0.6
0.8
3.0
3.5
3.5
3.0-
3.0-
2.5-
2.5-
2.0-
2.0-
1.5-
true pats
1.5-
Irue pats
our model
GRU-DDE- Sayas
trus cond expectatin
trus conditona expectatio
1.0
observed
1.0
observed
3.5
ou mode
cond Scnal expectation
TRU-OEE- Buycs tue cond tinnal expectation observel
obso
3.0
3.0
25
2.5
2.0
2.0
0.2
0.2
expectation
3.0
3.0-
2.5
2.0
2.0
1.5
1.5
1.0
1.0
2.0
0.2
0.6
0.6
0.8
2.0
3.0
0.2
0.6
0.5
0.8
2.0
Inue 2ach
cord ticnal expectation
3.0-
25
2.5
2.0
1.5
1.5
tuu Doth
GRU-ODE Bayes
Inue condiliios expecctaiio
1.0
1.0
abserved
0.2
0.6
0.8
0.2
0.8
5.5
Figure 19: Comparison of predictions for Omsttein CaieeibnnUhlebcckkkk paths at best epoch.
40

---

Published as a conference paper at ICLR 2021
True our modal
Tree path
GRU d0E Bayes
true conditooal expectation otserved
true candit onl expectation otserved
0.0
0.2
0.4
0.5
0.8
2.0
0.0
0.2
0.4
0.5
0.8
1.0
tree
Urue
model
GRU CDE Bayes
tue conditoral oxectation
Enee candit ond expectation atserved
atserved
1.2
1.75
Free path
Irue
1.0
0.8
1.00
0.25
0.8
3.0
2rue path
208 path
model
GRU-GDE Bayea
20uc Uonal expectation
200 Uonul expediaiion
abserved
W
0.0
0.0
GaU--or-- Dayas
Irue conditiona expectatin
2.0
Irue condition expectation observed
2.0-
3.5
WWW
WWW
1.0-
1.0
0.5
0.0
0.0
0.8
Figure 20: Comparison of predictions for Heston paths at best epoch.
41

---

Published as a conference paper at ICLR 2021
True our modal
Tree path
GRU 00E Bayes
true conditoral exacctation observed
true candit ona expectation obsetved
0.0
0.2
0.4
0.5
0.8
2.0
0.0
0.2
0.4
0.5
true
Uree
model
GRU CDE Bayes
tne conditordl oxectation
tue candit ond expectation otserved
atserved
0.2
0.2
true
Srue pats Dayes true condit onal expectation absetved
our model
true conditional expectation abserved
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.4
0.6
5.8
1.0
true path
Irue path
oir model
GRU-dDr--Rayws
tne conditoral oxeetation
the candit onal expedaatio abserved
abserved
Figure 21: Comparison of predictions for Black-Scholes paths at last epoch.
42

---

Published as a conference paper at ICLR 2021
3,5
Suyca
expectatoo
expectaiion
3.0
2.5
1.0
1.0
3.0
0.2
0.6
6.6
0.8
3.0
9.0
0.2
0.4
0.6
0.8
3.0
3.5
3.5
3.0-
3.0-
2.5-
2.5-
2.0-
2.0-
1.5-
true pats
1.5-
Irue pats
our model
GRU-DDE- Sayas
trus conditooal expectatin
trus conditona expectatio
1.0
observed
1.0
observed
3.5
oU mode
cond Sonal expectation
TRu-OEE- Buycs tue cond tinnal expectation observel
obso
3.0
3.0
25
2.5
2.0
2.0
0.2
0.2
expectation
3.0-
3.0
2.5
2.0
2.0
1.5
1.5
1.0
1.0
2.0
0.2
0.6
0.6
0.8
2.0
3.0
0.2
0.6
0.5
0.8
2.0
Inue 2ach
3.5
cord tional expectation
3.0
2.5
2.5
2.0
1.5
1.5
tuu DOtS
GRU-ODE Bayes
Inue condiliion expecctaion
1.0
1.0
abserved
0.2
0.6
0.8
0.2
0.8
5.5
Figure 22: Comparison of predictions for Ornstein-Uhlenbere paths at last epoch.
43

---

Published as a conference paper at ICLR 2021
12-
tue path
GRU-005 Bayes
true conditiona expectr
observed
p2th
madtll
conditonal
osserved
0.0
0.2
0.4
0.6
0.8
2.0
0.0
0.2
0.6
0.6
3.0
tree
12-
model
tue conditonal oxectation
atserved
30
conditorna
0.0
0.2
0.6
0.8
2.0
GRU-ODE
10
condidoni
0.0
0.2
0.4
0.6
0.8
2.0
0.0
0.2
0.6
0.6
2.0
10
T7e path
Erue path
our model
GRUIODE Dayes
20u8 cand torm rspectation
true conditiita ropectation
obea ved
observed
Myp
W
0.0
0.0
0.0
20.0
true pats
true condi oxpectation
our model
17.5
oce Sayas cond tinna expectation
2.0
observed
15.0
diserved
1.5
WWW
1.0
han
0.5
2.5
0.0
0.0
0.2
0.8
0.0
0.4
0.5
1.0
Figure 23: Comparison predictions for Heston paths at last epoch.
44


Reviews:
[{"review_id": "t1bNDglOq2Z", "reviewer": "AnonReviewer3", "summary": "", "questions": "", "limitations": "", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The authors propose a method for learning the conditional expectation of stochastic process in an online fashion. The paper bears a considerable theoretical treatment, derived from the stochastic filtering literature, which is present both in the main body of the paper and the appendix. Besides the model, the paper also aims to provide a theoretical justification of the convergence of their method. \n\nI find the contribution of the paper somewhat obscure, its aims to be incremental with respect to the previous literature, and the experimental validation heavily unconvincing. I support my recommendation through the following points: \n\n- Following the well known (by now) neural ODE and neural jump SDE, the contribution of the paper seems minor. The authors state that they focus on giving theoretical guarantees, however, these are specific and loosely validated experimentally\n- There is a fair amount of space dedicate to the theoretical presentation of the background, I agree with the importance of theory, but I failed to see how that theory supports the claims of the paper. \n-the experiments are limited: only 3 synthetic examples and only one real-world one.\n-the authors states that their method focuses on approximating (directly) the conditional expectation, this seems to be a different with the previous literature. However, if that's the case, the authors should consider more benchmarks such as linear filters (adapted to non-uniformly-spaced data), Gaussian processes, or general time series models. \n\nThis paper does have a contribution. My recommendation is that  the authors show it in a clearer (to the point) manner with an improved experimental validation.  ", "extra_scores": {}, "extra_text": {}}, {"review_id": "S-eg0sehh4", "reviewer": "AnonReviewer4", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "##########################\nThe paper proposes an algorithm and an analysis of its convergence.\nThe algorithm propose to learn a model of temporal data y_1 , ... , y_T given input x_1, ..., x_T\nThe observations are assumed to arise from a deterministic latent h \ngoverned by a piecewise continuous ode (in between consecutive times t_i, t_i+1)\nwith additional deterministic jumps at transitions.\n\nIn the ODE-RNN paper, the latent h can be expressed as a single ODE for the whole time horizon (rewriting the\njump with a skip transition).\n\nThis paper appears to me as taking this expression and choosing a particular bounded form for the\ndynamics, jump and readout functions\nA statistical asymptotic analysis of the convergence of the algorithms, for random times and inputs is given. \n\n##########################\n\nMethodology:\nI find the paper quite difficult to read, I blame both its structure and my lack of ease with the mathematics used here.\nHowever from what I have understood of the algorithm proposed, I find the methodological contribution very limited.\n\nClarity:\nI come from the machine learning community and read with no difficulty \npapers cited in the related work section.\nIn comparison, I find this paper extremely difficult to read and parse despite containing the same kind of information.\n\nAsymptotic analysis.\nI leave to other reviewers the evaluation of the convergence analysis.\nMy evaluation being partial, my confidence rating is set accordingly.\n\n* For a machine learning paper presenting in the end a 3 line simple algorithm, the paper contains\na lot of superfluous mathematical notation that crowds the paper and make the reading very tedious.\nMany of the papers cited Brouwer 2019, Rubanova 2019, Li 2020, offer a much smoother read in that respect.\nAs is, this paper feels better suited to a more specialist statistics venue.\n\n* For example, many elements are introduced in the main text and are not really necessary to understand what the paper does\nThe detailed section on random inputs is used only in a theorem coming later, why have it in the main text in so much details.\nOn the other end, a description of the method this paper builds on is left into appendices.\n\n#########################\n\nAdditional comments:\n* the formatting of the references is very inconsistent, please update\n\n", "extra_scores": {}, "extra_text": {}}, {"review_id": "LANJaI_15K0", "reviewer": "AnonReviewer1", "summary": "", "questions": "", "limitations": "", "rating": 6, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "The submission studies a simplified model of ODE-RNN and GRU-ODE-Bayes theoretically, proves convergence results, and presents experimental results in companion with the theoretical results. \n\n\nThe paper does a good job in defining concepts precisely. Though, this has come at a cost of highly complex notation, which may hinder the average researcher in the ML+diffeq community, who may not have a strong background in probability theory, to understand the paper. I would therefore recommend the authors to simplify the notation by deferring the precise mathematical definition of concepts such as information sigma algebras to the appendix. \n\n\nThe section (sec. 2.4) on optimal approximation of a stochastic process in the main text is somewhat vague. Optimality certainly depends on the cost function being considered, in which case, the appendix states that the 2-norm is used here. The particular norm being used is somewhat independent of the construction of the probability space, e.g. we could consider the same prob. space and evaluate the difference between the random variable and the fixed prediction using some other function, say the metric function induced by the L1-norm). This makes terms such as “L^2(omega X omega tilde, ...)-minimizer” somewhat confusing. Note my comment here is somewhat handwavy about the precise technicalities, but it should convey the relevant idea. \n\n\nMy main concern regarding the paper is about novelty. It seems that the model considered in section 3 falls broadly in line with ODE-RNN and GRU-ODE-Bayes. On the other hand, the experiments section also doesn’t compare against latent ODE, which is a strong but relevant baseline. \n\n\nThe section (sec. 4) on theoretical convergence results mostly assume that the ERM can already be found. This rather strong assumption therefore leaves the theorems in that section not unexpected, and at the same time, less relevant for practitioners. It is also unclear whether convergence rates can be derived. \n\nThe paper does a decent job in clarifying its relationship with prior work. \n\n\nPost-rebuttal:\n- I thank the authors for improving the presentation of the paper and including additional experiments comparing to latent ODE. ", "extra_scores": {}, "extra_text": {}}, {"review_id": "9G-gmwicrfq", "reviewer": "AnonReviewer2", "summary": "", "questions": "", "limitations": "", "rating": 7, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "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": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "Summary: This paper introduces Neural Jump Ordinary Differential Equations as a method for learning models of continuous-time stochastic processes sampled at random time epochs. Specifically, the paper studies the problem of estimating the marginal conditional expectation (i.e., the L2 optimal approximation conditional on the available information) by estimating an auxiliary stochastic differential equation, parameterized by  neural networks, that approximates the conditional expectation of the process of interest at each point in time. The neural networks are trained by using a “randomized” mean squared-loss objective. The main theoretical results in the paper include asymptotic consistency of the optimal objective value in the limit of a large neural network, as well as consistency of a Monte Carlo sample average estimator of the value. The paper also establishes the L2 convergence of the estimated auxiliary solution to the marginal conditional expectation.\n\nThe technical details in the paper are mostly sound, and I believe it should be of interest to a wide community. The question of estimating stochastic models sampled at regular or irregular intervals is of broad utility. There are some technical issues however, but these can be resolved I believe. In particular, in the conclusions of Theorem 4.1 and 4.2, it seems as though the authors claim almost sure convergence, unless I am misunderstanding their statement.  What the authors establish is convergence in L2, but why does is imply almost sure convergence? Wouldn’t one require uniform integrability to conclude more? Furthermore, this is not a process level convergence result, and therefore I do not believe that they can conclude (as in Remark G.3) that the limit holds almost surely. (Also, the authors seem to suggest tin Remark G.2 that they’re not establishing L2 convergence, but this could be a problem with the writing). \n\nComing to the writing, I note that I did find the paper somewhat sloppy in its use of terminology and notation. For instance on p.1 the authors state “...while stochastic processes are continuous in time...” This is not quite true, since one can define discrete-time stochastic processes. I also found the discussion around  justifying “irregular” sampling of the stochastic process to be poorly written. In particular, it is stated that “...dividing the time-line into equally-sized intervals...is again making assumptions and information is lost...” well, any sampling will involve a loss of information, and the randomized sampling process described in this paper also involves assumptions. I don’t think this comment is appropriate. Furthermore, the authors do not make a clear case for why their irregular sampling procedure is appropriate. I’m quite certain that the sampling process introduces bias into the estimation; for instance, Theorem 1 of ref. [1] below provides sufficient conditions under which an “irregularly” sampled estimator of a functional of an SDE is unbiased. The authors must do a better job of justifying their method. I would also urge them to add an example of a randomized sampling process; for instance, a Poisson process sampler would satisfy their definition, in which case the sampling time epochs form an ordered statistic.\n\nComing to the development of the stochastic model, it is unclear to me as to why all of the random “objects” cannot be defined on the same sample space. Essentially, couldn’t one view the sampling process as a point process on the same sample space supporting the SDE? \n\nNext, in Prop. 2.1, the authors state that the optimal adapted process approximating process is \\hat{X}_t — but \\hat{X_t} is only defined pointwise (i.e., at each time ‘t’) and it is not defined as a stochastic process. Indeed, for that the authors must describe the finite dimensional distributions for all finite sets of time epochs to define the stochastic process. I believe it is inappropriate to call this a stochastic process. This doesn’t affect the main results, since the authors only establish convergence in an L2 sense, where the full distribution is not necessary. \n\nSome further minor comments:\n1. Change the term “observation dates” to “observation epochs”.\n2. Change “amount of observations” to “number of observations” (or samples).\n3. On P.3 in the definition of \\lambda_k, the set \\mathcal{B}([0,T]) is undefined.\n4. The notation defining the function \\tilde{\\mu} is very confusing, please change.\n5. P.4 “...since the variation of u...” should be “...since the total variation of u...”\n6. What do you mean by “ergodic” approximation of the objective? Isn’t it simply a sample average approximation? Which ergodic theorem is playing a role here?\n7. I would also urge you clearly define what you mean by \\mathbb{L}-convergence, for completion. \n\n[1] Unbiased Estimation with Square Root Convergence for SDE Models, Chang-Han Rhee and Peter W. Glynn.", "extra_scores": {}, "extra_text": {}}]

Author rebuttal:

Meta-review: This paper proposes a refinement, and analysis of, continuous-time inference schemes.

This paper got in-depth criticism from some very thoughtful and expert reviewers, and the authors seem to have taken it to heart.  I'm still worried about the similarity to GRU-ODE-Bayes, but I feel that the clarifications to the general theory of continuous-time belief updates is a worthy contribution, and the method proposed is a practical one.  One reviewer didn't update their score, but the other reviewers put a lot of thought into the discussion and also raised their scores.

I do think the title and name of the method is a bit misleading - I would call it something like "Consistent continuous-time filtering", because the jump ODE is really describing beliefs about an SDE.
Decision: Accept (Poster)