Zero-shot Imputation with Foundation Inference Models for Dynamical Systems
Authors: Patrick Seifner
Kostadin Cvejoski
Antonia Körner
Ramses J Sanchez
Venue: ICLR 2025 Poster (2025)
Keywords: Zero-shot imputation
foundation models
time series imputation
dynamical systems
amortized inference
zero-shot interpolation
foundation models for time series
Primary area: learning on time series and dynamical systems
TL;DR: We introduce a framework for zero-shot imputation of missing time series data, whose underlying dynamics are assumed to be determined by ODEs. Our foundation models often outperform state-of-the-art methods, which are trained on the target datasets.

Abstract: Dynamical systems governed by ordinary differential equations (ODEs) serve as models for a vast number of natural and social phenomena. In this work, we offer a fresh perspective on the classical problem of imputing missing time series data, whose underlying dynamics are assumed to be determined by ODEs. Specifically, we revisit ideas from amortized inference and neural operators, and propose a novel supervised learning framework for *zero-shot time series imputation*, through parametric functions satisfying some (hidden) ODEs. Our proposal consists of two components. First, a broad probability distribution over the space of ODE solutions, observation times and noise mechanisms, with which we generate a large, synthetic dataset of (hidden) ODE solutions, along with their noisy and sparse observations. Second, a neural recognition model that is trained *offline*, to map the generated time series onto the spaces of initial conditions and time derivatives of the (hidden) ODE solutions, which we then integrate to impute the missing data. We empirically demonstrate that *one and the same* (pretrained) recognition model can perform zero-shot imputation across 63 distinct time series with missing values, each sampled from widely different dynamical systems. Likewise, we demonstrate that it can perform zero-shot imputation of missing high-dimensional data in 10 vastly different settings, spanning human motion, air quality, traffic and electricity studies, as well as Navier-Stokes simulations — *without requiring any fine-tuning*. What is more, our proposal often outperforms state-of-the-art methods, which are trained on the target datasets.

Our pretrained model, repository and tutorials are available online.

Published as a conference paper at ICLR 2025
ZERO-SHOT IMPUTATION WITH FOUNDATION INFERENCE MODELS FOR DYNAMICAL SYSTEMS
Patrick Seifner'.2 Kostadin Cvejoski Antonia Körner & Ramsés J. Sánchez! Lamarr Institute', University of Bonn² & Fraunhofer IAIS seifnerêcs. uni-bonn.de, sanchezéca, uni-bonn.de
ABSTRACT
Dynamical systems governed by ordinary differential equations (ODEs) serve as models for vast number of natural and social phenomena. In this work, we offer a fresh perspective on the classical problem of imputing missing time se- ries data, whose underlying dynamics are assumed be determined by ODEs. Specifically, we revisit ideas from amortized inference and neural operators, and propose a novel supervised learning framework for zero-shot time series imputa- tion, through parametric functions satisfying some (hidden) ODEs. Our proposal consists of two components. First, a broad probability distribution over the space of ODE solutions, observation times and noise mechanisms, with which we gener- ate large, synthetic dataset of ODE solutions, along with their noisy and sparse observations. Second, a neural recognition model that is trained offline, to map the generated time series onto the spaces of initial conditions and time derivatives of the (hidden) ODE solutions, which we then integrate to impute the missing data. We empirically demonstrate that one and the same (pretrained) recognition model can perform zero-shot imputation across 63 distinct time series with miss- ing values, each sampled from widely different dynamical systems. Likewise, we demonstrate that it can perform zero-shot imputation of missing high-dimensional data in 10 vastly different settings, spanning human motion, air quality, traffic and electricity studies, as well Navier-Stokes simulations without requiring any fine-tuning. What is more, our proposal often outperforms statc-of-the-an art meth- ods which are trained on the target datasets.
Our pretrained model, repository and tutorials are available online!
INTRODUCTION
Dynamical systems are mathematical systems that change with time according to fixed evolution rule, and serve representational and analytical tools for phenomena which generate patterns that change over time. Very often, the recorded changes of these empirical patterns are such that they can viewed as occurring continuously time, and thus can be represented mathematically by systems whose evolution rule is defined through differential equations. Dynamical systems governed by ordinary differential equations (ODEs) correspond an important subset of these models, and describe the rate of change of single parametric function x R - R which represents the state of the (D-dimensional) system, as time evolves, by means vector field f. TR+ RR. - RRD. In equations, we write
dx(t) dt
xt)) = ft,xx(())) where x(()
(1)
These deceptively simple systems have had fundamental role in our understanding of many natural processes across nearly every scientific discipline - from their very introduction and application to celestial mechanics in the late seventeenth century (Newton, 1687; Bernoulli, 1712), to their function as models concentration changes in molecular reaction networks (Hoff, 1986); models of population oscillations in biology (Lotka, 1925; Volterra, 1927); of atmospheric convection and its chaotic features (Lorenz, 1963); and of the coherent, high energy modes within turbulent flows
htttps://iimdsciencce github..hhb... 100..1///paaa o/pp/ IM/intro.htmi

---

Published as a conference paper at ICLR 2025
Fample
Emepie
80
f(t)
0,0
X
3
V3
$ w1
Un
u
N
Foundation Inference Model (FIM). Left: Graphical model of the data generation model (Eq. 2). Filled (empty) circles represent observed (unobserved) random variables. Cen- ter: Schematic representation of imputation mechanism for temporal-wise missing patterns (Eq. 7). Right: Neural operator module processing the context vector ug and query time t (Eq. 4 and 8).
Figure 1:
(Noack et al., 2003), just to name a few - and continue to be the go-to mathematical objects for the representation of dynamic phenomena today.
In this work, we consider the general problem of imputing missing values in time series data, recorded from some empirical process (y" IR+ - RD) whose dynamics are assumed be gov- emedd by some unknown ODE. In other words, we assume that both available and missing values in the series y" (F1),... ............. y" (Tt) correspond to the values taken by the solution x(t) of some hidden ODE, at the observation times Ti,..., 7. potentially corrupted by some noise signal of which only a few statistics are known. Therefore, the goal to infer the ODE solution x(t) that best interpolates the noisy time series y"
and hence imputes its missing values.
The current machine learning paradigm tackles this problem by (implicitly) constraining models to handle single process only. That is. practitioners typically encode their inductive biases into either the model architectures or the training objectives, and optimize the model parameters to fit a single empirical distribution (see e.g. Section 2). One disadvantage of this approach is that models trained to fit a single process tend to be overly specific to its distribution, and thus can rarely be reused to impute the missing values of second one, even when both processes are assumed to be governed by e.g. similar ODEs. Another disadvantage is that, to succeed, the paradigm requires practitioners to have access to enough observations the process they study, to train and test their models from scratch, and that they also have the experience and expertise to face the trials and tribulations of their intricate training procedures.
In this paper, we instead frame the imputation problem as an instance of amortized inference, in the sense introduced by Stuhlmüller et al. (2013). Indeed, in lieu of training one complex model single empirical process, we train neural recognition model offline to infer a large and varied set of ODE solutions x(t), from a synthetic dataset that composed of noisy series of observation on those solutions, displaying different missing value patterns. Somewhat more precisely, we train our model supervised fashion, to infer both the (latent) initial conditions x(0) and (latent) time derivatives (l) that determine the target set of ODE solutions. However, opposite to Stuhlmüller et al. (2013) and other follow-up works, like that by Paige & Wood (2016), who treat their recognition models as auxiliary Monte Carlo methods, we employ our pretrained models to directly impute different synthetic, simulation and experimental datasets, without any parameter fine-tuning. We therefore adopt the "zero-shot" terminology introduced by Larochelle et al. (2008) Let us briefly note here that we have recently used this general amortized inference framework to train models that perform zero-shot inference of both Markov jump processes (Berghaus et al., 2024) and stochastic differential equations (Seifner et al., 2025).
In what follows, first review previous work on time series imputation in Section 2. Section 3 introduces our main ideas, the synthetic dataset encoding our assumptions, and our recognition
We invite the reader to check Appendix A, where we also comment on the differences hetween our method- ology and other meta-learning methods.
2

---

Published as a conference paper at ICLR 2025
model, which we name Foundation Inference Model (FIM) for dynamical systems. In Section 4, we report our experimental findings and empirically demonstrate that: (i) the hierarchical structure underlying FIM - which treats x(0) and x(t) as latent variables allows to reconstruct, in a zero-shot fashion, the phase portrait of complex dynamical systems; (ii) FIM is able to impute, in zero-shot mode, missing values in set of 63 noisy time series, each of which is sampled from dynamical systems of different dimensionalities; and (iii) the same (pretrained) FIM can perform zero-shot imputation of vastly different, high-dimensional, experimental and simulation data, while often out- out-performing state-of- of-the-art models which are trained on the target datasets. Finally, Section 5 comments on the main limitations of our methodology, and closes the paper with some concluding remarks about future work.
RELATED WORK
In its most general form, the problem of imputing missing time series data with an ODE model involves inferring the vector field (f Rt x RD + PR that defines the (hidden) ODE in question. In practice, however, one does not need explicitly infer the vector field in order to impute the missing data- - and we will take this perspective in the present work. Prominent examples are the variants by Che al. (2018) and Cao al. (2018) that model the dynamics interpolating the data in some latent, high-dimensional space via linear ODEs, or the neural ODE model of Chen al. (2018) that learns latent albeit nonlinear ODEs. The latter has in fact been modified suit very different imputation and "smoothing" scenarios (Rubanova et al., 2019; Yildiz et al., 2019; Norcliffe al., 2021; Seifner Sanchez, 2023). Other works depart from ODEs and assume the dynamics are stochastic. Examples thereof include the works by Fortuin et al. (2020) and Fang al. (2024), which leverage Gaussian processes to represent the time evolution, and that by Tashiro et al. (2021), which instead utilizes conditional score-based diffusion models. Researcher have also recently dispensed with continuous-timg models altogether, and deployed self-attention mechanisms to impute missing data (Du et al., 2023). Wang al. (2024) provides a recent and comprehensive review on these and many other imputation methods, and we refer the reader to it for completeness.
Regardless of whether they rely on ODEs or not, all the models above find themselves under the umbrella of the classical paradigm, insofar as they are all optimized with respect to a single em- pirical process. There are, nevertheless, two recent exceptions. Similar us, Becker al. (2023) and d' Ascoli et al. (2024) generate large datasets of ODE systems and their (noisy) observations. Opposite us, they attempt explicitly infer the vector fields, and do so in symbolic form. The work of Becker et al. (2023) is however limited to one dimensional ODEs, whereas that of d' Ascoli al. (2024) is limited to six -i-dimensioonll ones. To the best of our knowledge, we present the first zero-shot solution that is applicable to real-world empirical processes of any dimensionaliry.
FOUNDATION INFERENCE MODELS FOR DYNAMICAL SYSTEMS
In this section, we introduce novel methodology for zero-shot imputation of missing time series data. Let the series dimensional empirical process y" (t), where each observation is represented by vector y (T;) e RR, with j 1,...,I. Suppose that some of the components of these observation vectors are missing, and the goal is to impute them back. Our proposal frames the problem of estimating these missing values as an inference task, in which one seeks to infer the ODE solution x(t) that best interpolates the series (T))....... ym")., (n), and thus imputes its missing values.
correspond to a sequence of 1 observations on some D-
The classical formulation of the imputation problem typically involves different missing patterns and in this work we focus on two of them. The first one is the so-called point-wise missing pattern, where individual vectors in the series randomly lack some of their components. The second one is the temporal missing pattern, where certain components of the vectors in the series are missing over consecutive observation times. To handle them, make the following two simple assumptions.
'e choose to name our model foundation model because it goes in line with the definition proposed by Bommasani et al. (2021). To wit: foundation model is any model that is trained on broad data (generally using self-supervision scale) that can he adapted wide range of downstream tasks.
3

---

Published as a conference paper at ICLR 2025
First, we assume that for every time series of observations featuring point-wise missing patterns, one can always find a certain time scale Tample or some (sequential) subset of observations, for which the best interpolating ODE solution is "simple" Furthermore, we assume that the set of all such simple parametric functions can be well-represented by a heuristically constructed synthetic distribution. Second, we assume that time series featuring temporal missing patterns involve more complex interpolating functions, meaning that no such Tample is to be found in this case. Although more complex in nature, we assume that these functions are locally "simple", and that they of- exhibit generic secular and seasonal structures, which encode important information about the missing values and can be well-represented by a second, synthetic distribution over parametric func- tions. Should these two general assumptions hold true, a model trained infer both our synthetic set of "simple" parametric functions, and that of functions exhibiting generic, secular and seasonal structures, from noisy and sparse observations on them, will automatically interpolate any unseen sequence of empirical observations and, consequently, impute all of its missing values. In the ex- perimental section below, we empirically demonstrate that this is indeed the case in a variety of scenarios.
Our methodology thus consists of two components. The first comprises synthetic data generation model, and encodes our beliefs about both the "simple" ODE solutions that interpolate the data locally, and the more complex ones that feature global, generic structures. The second corresponds to a neural recognition model of minimal inductive biases, that maps sets of noisy observations onto the space of parametric functions. In what follows, we delve into the details of these two components and name recognition model Foundation Inference Model (HIM) for dynamical systems. Figure illustrates the FIM framework.
3.1 SYNTHETIC DATA GENERATION MODEL
In this subsection, we describe the synthetic data generation model we use to sample large and varied set of ODE solutions, together with their noisy and sparse observations. Given that every ODE solution is parametric function of time, and that each component of any such D- D-dimensional function is itself a -ieddimenioonal, parametric function of time, we focus only on the space of 1D time series. In other words, we opt for a channel independent strategy (Nie et al., 2023; Han et al., 2024). Let us then define the probability of observing the 17 noisy time series /1, U2. ..... I1 at the observation times 0 S T1 < T2 < any of the components of some D-dimensional process- - as
T. S 1 - which might correspond to the values taken by
II P (Unlar O)(()) [a:- xn f."' [(s)ds)) F=1
(2)
where 8(.) represents the Dirac delta function, which identifies the ODE solution xt)), evaluated at time The with zo J {sdds. That is, we understand our ODE solutions being determined by some initial condition Zo and the parametric function f(t), which represents the time derivative of (t). Note that we use the notation z; to denote 2(T.), and that we denote both random variables and their values with the same symbol. Let us now specify each term in Eq. 2, starting from the right.
DISTRIBUTION OVER INITIAL CONDITIONS. We define the prior p(zo) over initial conditions as standard Gaussian distribution. A standard Gaussian suffices, because the values of every time series we process are first normalized to lie on the unit interval (see Section 3.2.1).
DISTRIBUTION OVER PARAMETRIC FUNCTIONS OF TIME. The prior distribution p(f.ny) factor- izes as p(ffl )P(s)) with p(ny) the prior over the hyperparameter set 7). The conditional dis- tribution p(f is a distribution over the space of parametric functions of time, defined on the unit interval. It encodes our beliefs about the class of interpolating functions we expect to find in practice. Indeed, we briefly motivated earlier, we design two such distributions. One represents "simple" parametric functions that are assumed to be typical interpolating functions imputing point- wise missing patterns the characteristic time scale Tsimple The other represents functions that are locally simple, but that exhibit secular and seasonal structures at longer time scales (i.e. 7 > Tample)- Structures that, in turn, are assumed to carry crucial information about temporal missing patterns.
By simple we loosely mean functions that can be approximated by low-degree polynomials, and functions whose Fourier transform has support at low frequencies only.
4

---

Published as a conference paper at ICLR 2025
We define these distributions by means of random Chebyshev expansions and Gaussian processes with different kernels, and refer the reader to Appendix 1.1 for details.
DISTRIBUTION OVER OBSERVATION TIMES. We define the prior distribution Pani (Ti......... Tr, 7g with hyperparameter set 7g. to represent missing data patterns that we expect to be relevant in real- world imputation tasks. Again, we consider two such distributions. The first one represents point- wise missing patterns and allows for regular and irregular observation grid instances with different observation count. The second one represents temporal missing patterns and combines point-wise patterns with randomly located observation gaps. We allow the latter to be as large as one-third of the unit interval. Appendix B.2 provides details regarding our implementations.
DISTRIBUTION OVER NOISE PROCESSES. When recording any empirical process, one typically only has access to the mean square error of those measurements. According to the maximum entropy principle, a Gaussian distribution is the best guess one can make about the noise distribution actually, it the most likely distribution - given the available information (that is, given those first two moments) (Jaynes, 2003). We choose our noise model Puuise (ya [z;; o)) accordingly, and set p(o) to also be a Gaussian distribution zero mean and variance 10
We use the generative model, Eq. 2 above, to sample large and varied set of noisy and sparse ODE solutions. We refer the reader to Appendix 4 for details on the specifics of the sampling procedure.
3.2 FOUNDATION INFERENCE MODEL
In this subsection, we introduce recognition model that exploits ideas from neural operators (Lu et al., 2021; Kovachki et al., 2023) to map time series data onto parametric functions. Indeed, given set of noisy observations (y1, ..,, (yi, 72) on some ODE solution s=(t) sampled from the data generation model, Eq. 2 above - our goal is to infer both the 1 D parametric function f(t) and initial condition zo that specified ztt) in the first place. In other words, we want to reverse the data generation process, Below, we first introduce a neural interpolation model that is trained to infer the distribution p(f(n)) over "simple" functions, from time series exhibiting point-wise missing patterns. Note that in this setting Tsimple is, by construction, of order one. Later, we introduce a second interpolation model that trained to infer p(f(m)) over functions that are locally "simple" but display global structures at time scales >>> Tsimple from time series characterized by temporal missing patterns. In what follows, we denote the first interpolation model with FIM- É, for it handles local dynamic features, and use F IM to refer to the the second one.
3.2.1 1
FIM FOR INTERPOLATING POINT- WISE MISSING PATTERNS
In order to interpolate time series featuring point-wiss missing patterns, and values of every scale, we first need to normalize every input sequence and rescale their target parametric functions accordingly (see Appendix C.1 for details). Let us label the set of normalized, noisy observations with Y, the space of rescaled 1D parametric functions with F and that of rescaled initial conditions with Xn- Our interpolation problem can then be understood as the problem of mapping ) onto both F and Xg. We begin with the first of these two maps.
Let us use d² and 1h to denote feedforward (FFN) and sequence processing neural networks, re- spectively. Let us also denote the trainable network parameters with 0. We now define the function
h (u", (i)), with u² 2g (v-l (y,, (3,,, T,,-.., (t,
(3)
where o and the composition ob wif can be interpreted as the trunk and branch nets of DeepONets (Lu et al., 2021), while uz is the context vector encoding the (context) points (y1, Ti...., U1, 7)) Given the vector-valued, parametric function (t). we now define the mean and variance of a Gaussian random variable over the possible values of the estimated time derivative as follows
Î(1) (1)= (1)), log Var( og (h²(1)),
(4)
where, similar to the works of Lakshminarayanan et al. (2017) and Valdenegro- & Mori (2022), the variance Var(f) is used to represent the model's uncertainty in the estimation of f(t (see e.g. Figure 3 for an illustration). Next, to perform the map between y and Xy, we also model the initial condition f()) as Gaussian random variable, whose mean and variance are given by
so =v((1)) log Var(fo) =b((
(5)
5

---

Published as a conference paper at ICLR 2025
GT trajectory
GT Initicial values
FN trajectory
Velocity
Valocity
Accala Accalaration
Position
Position
Valocity
Figure 2: Zero-shot phase portrait reconstruction of two dynamical systems. Left: Van der Pol oscillator. Center and Right: Rössler attractor in the position-velocity and veleeity-acceleration planes, respectively.
Equations 4 and allow us to express the function interpolating the time series (y1, ...... (yt, m)
at any desired time E (0, 1) - as &. + J fs))ds. We train F FIM- in a supervised fashion, to maximize the log-likelihood of Eqs. 4 and 5, and minimize the one-step reconstruction error of the integrated solution, both with respect to the distribution over "simple" functions. Thus, FTM-f can be understood an estimator of the posterior distribution over the space of interpolating functions, given the context points (3/, (1,,,,, 7t). We refer the reader to Appendices C.2 and C.3 for details regarding the model architecture and training objective.
3.2.2 BEYOND SIMPLE FUNCTIONS: PROCESSING DATA ANY LENGTH AND DIMENSIONALITY WITH FIM-E
Let us briefly comment how use our (pretrained) FIM-& process time series of any length and dimensionality. We provide further details (and limitations) of our strategies' implementation in Appendix C.6. During training. F WI-E processes only 1D time series with at least L-min and at most Lim observations (see Appendix B.2). We define L as the minimum number of context points FIM generally needs to function, shorter time series are considered out-of-distribution Similarly, time series with more than Lmux observations also fall outside the distribution, not just because of their length, but because their (hidden) interpolating function may not be well-represented by our distribution of "simple" functions. To address this limitation, we first split any target time series longer' than Lmux into successive and overlapping time windows, ensuring that both their observation count and dynamic features remain within distribution. That is, we assume that each time window spans time scale of order (TTameee. Then, we combine the local FIM- estimates obtained from each window into global one for the entire target time series (see Appendix C.6). In an analogous manner, we adopt a channel independent strategy and process each component of any target, D D-dimensional process independently with F FIM-(.
We empirically demonstrate the efficacy of these approaches in Section 4 below. We also investigate the effect of changing the length (and number) of the overlapping windows in Appendix F.4.4.
3.2.3 FIM FOR INTERPOLATING TEMPORAL MISSING PATTERNS
In this subsection, tackle the interpolation of time series displaying temporal missing patterns by decoupling trends and seasonality from local fluctuations. Suppose we are given a noisy time series with observations, whose values and observation times lie within the unit interval. Suppose that this time series has been split into K sequential (i.e. ordered) sets
V1.
yar +1...., ..... wwa
y1
V1,
where wi is the number of observations within the ktth set. Suppose now that the gth set is missing, and the goal is to impute it back. See e.g. the center image Figure 1.
Note that this approach can also be applied to time series of length shorter than Tamax but whose dynamic features believed to he out-of-diitrriution.

---

Published as a conference paper at ICLR 2025
Table 1: MAE of the inferred time derivative xtt) and ODE solution x(t) on ODEBench. The
standard deviation is calculated across 10 samplings of the corruption schemes.
Inferred time derivative x(t)
Reconstructed ODE solution x(t)
Model
= 0
7 = 0.05
= 0
? == 0.05
ODEForer
8.00 + 0.40
7.90
0.60
1.18 + 0.05
1.16 1 0.05
FIM-
2.44 1 0.05
3.79 = 0.05
0.17 1 0.01
0.34 1 0.01
We assume that locally, within every (even the missing one), the functions underlying the data are well represented by our synthetic distribution of simple functions. That is, we assume that each set spans time scale of the order (Tampl) and thus, that its underlying function can be modelled well with our pretrained FIM-f. We also assume that beyond those time scales (T 2 Tsimple), there exist inter-set, global structures and correlations that carry information about the missing (i.e. gth) set. Our task is to define a model that encodes precisely this information. In other words, we require a model that extends the context of FIM- to longer time scales.
Since F FIM- is trained to deal with normalized data, we first normalize each of the available sets, and denote with 5j the statistics containing information about the local scale (i.e. the norms) of the jthh set (see Appendix D for details). Let us now process each (normalized and available) set with our pretrained encoding networks ol o v/f to obtain the sequence
(uf, si)..... (u 1), (u 1q+1), (u&: sk) with u² o8 ( Ti,..., Muy, Twrs
where, for simplicity, we relabelled the sub-indices (of the jth set) on the right hand side as m < E' W; + m, with mn an integer between 1 and wy.. Our second interpolation model F FIM consists of sequence processing network v/f with trainable parameter set yo, that computes
u 2 (Luf: ..... (. (u
$4+1). ................................ (u sk)),
i.e. the context vector for the missing (i.e. qth) set. At this point, we can use the pretrained Φ1 og dg and og networks of F IM-I to estimate the function ftt and its variance along the gap. That
f(t) = of (t))))), (t)), log Var(f)(t) og (t), with h (t) = of og of (t)
(8)
We optimize yp in a supervised manner - while keeping fixed - to maximize the likelihood of ftt) along the missing set, with respect to our synthetic dataset of complex functions that feature global, generic structures and temporal missing patterns. We refer the reader to Appendix D, where we provide additional details and discuss about how integrate f(() to infer 2()) along the gap.
EXPERIMENTS
In this section, we methodology widely different imputation tasks, which involve datasets of varying complexity, different dimensionalities and noise signals of very varied nature. We use FIM-f and FIM to impute- in zero-shot mode- missing data featuring point-wise and temporal missing patterns, respectively. To be precise, we apply our pretrained models directly to the test sets of the target datasets, without any parameter fine-tuning. F IM- was pretrained to infer 2M ODE solutions from time series with (Imin- Lmax) set to (4, 128). F FIM was pretrained to infer 500K local ODE solutions from time series with observation gaps that amounted to one-third of the data. Addi- tional information regarding model architecture and hyperparameters, training details and ablation studies can all be found in Appendices C and D.
METRICS. Below we evaluate the performance of our models wrt. the mean-absolute error (MAE). In the Appendix, we also report our results wrt. the root-mean-square error (RMSE) and, in some cases, the R² coefficient of determination and the mean-relative error (MRE). Formulas for these metrics can be found in Appendix E. BASELINES. Depending the task, compare findings against the (symbolic) ODE- Former (d' (d'Ascoli et al., 2024) model; the LatentODE (Chen et al., 2018), NeuralODEPro- cesses (Norcliffe al., 2021) and BRITS (Cao et al., 2018) models; the GP-VAC -P-VAE (Fortuin et al.,
7

---

Published as a conference paper at ICLR 2025
Observations
Imputation
FIR
Std(/)
20
2000
0.0
0.2
0-4
0.6
0.8
1.0
0.0
0.1
0.2
0.4
0.0
0.2
0.6
Time
Time
Time
Figure 3: Zero-shot imputation with FIM. Left and Center: Point-wise missing imputation in (one dimension of) the Beijing and GuangZhou datasets, respectively. (Right): Temporal missing impu- tation in single PCA dimension of the Motion Capture dataset.
2020) and the ecore-based CSDI (Tashiro al., 2021) models; the (self-ateention) SAITS (Du al., 2023) model; and the (Gaussian process) nPODE (Heinonen et al., 2018) model, among others. Besides ODEForeer, all other baselines trained on their larget datasets.
4.1 PHASE PORTRAIT RECONSTRUCTION
Before testing our methodology in real-world imputation scenarios proper, we explore its ability to accurately reconstruct the phase portrait of complex dynamical systems, using noisy and sparse observations on only one of their coordinates.
Phase portraits are geometric representations of the trajectories of dynamical systems in phase space. They help unveil attractor sets or limit circles and allow to determine, inter alia, how chaotic the sys- tems in question might be (Benettin et al., 1976). In the seminal paper by Packard al. (1980), the authors demonstrated that one can obtain a faithful phase-portrait representation of any D- dimensional dynamical system, through one of its coordinates, say x (t), and all its time derivatives (ir (t), 31 (t)..... up to order D - 1. Carrying out such a reconstruction from noisy and sparse observations on x1 (t) alone, not only entails interpolating the data, but also numerically computing the time derivatives of the interpolating function. In this section, we empirically demonstrate that one can use the hierarchical structure underlying FIM which treats the time derivatives á (t) (f()) in our notation) as a latent variable to reconstruct, in zero-shot fashion, the phase portrait of complex dynamical systems.
Suppose we simulate dynamical system to obtain the solution x(t) over the interval J0, T'. rep- resented on a fine-grid of points. We then introduce two types of data corruption. The first is where (N,
multiplicative noise Vij = (1 + e)x.y. with second is random subsampling, where a fraction p of the fine-grid is removed (and thus corresponds point-wise missing data pattern). For concreteness, set T 10 and 0.05. The task is to reconstruct the phase portrait of the dynamical system from the 1D time series /11,-- V/1 alone.
1.......... D and j -1..L. L. The
Let start with the (nonlinear) Van der Pol oscillator, whose dynamics are given by the second- order ODE ë distorted as increases the strength of the nonlinear term p. We μ 0.5, simulate the system with 12 different initial conditions and record only the noisy position of the oscillator 128 times per trajectory. The left panel of Figure 2 shows the 12 trajectories inferred by FIM- on (&1. 21) phase space. together with the ground-truth, on plotting grid of 2048 points, which amounts to an imputation of 1920 missing points. The agreement good and the visible deviations are due to small discrepancies in the estimation of the oscillator's velocity, for some initial conditions. We suspect that such rapidly changing functions are not well-represented by our synthetic distribution of "simple" functions. Next we consider the Rössler system, which is 3D dynamical system that features a chaotic attractor, and record up to 2048 noisy observations on 21(t) in order to discem its main features. The center panel of Figure 2 displays the trajectory inferred by F FIM-( l in (ø,, I1) space, on a plotting grid of 8192 points, which is in very good agreement with the ground- truth. Similarly, the right panel of the same figure displays the inferred trajectory in (Ø,, i1 space,
= O. This system features a limit circle around [x|| = 1, which gets
8

---

Published as a conference paper at ICLR 2025
Table 2: MAE (at missing values) on 8 datasets featuring 50% point-wise missing patterns. Base- lines scores for GuangZhou and Solar are extracted from Fang et al. (2024). The rest are extracted from Du al. (2024). See Appendix H for standard deviations and additional baselines.
Air quality
Traffic PeMS
Electricity ETT
Method
Beijing
Italy
GuangZhou
Pedestrian
Solar
Electricity
BRITS
0.169
0.321
3.335
0.287
0.259
1.985
0.238
1.124
SAITS
0.194
0.285
3.391
0.302
0.205
1.827
0.223
1.399
GP-VAE
0.258
0,453
3,419
0.346
0,451
1.810
0.414
1.099
CSDI
0.144
0,958
3,202
0.288
0,351
0.804
0.318
0.798
BayOTIDE
2.687
0.734
FIK-(
0.166
0.215
2.427
0.365
0.273
0.595
0.279
0.083
obtained by applying FIM- twice to the noisy observations on x (1). Overall, the agreement is still strong and F FIM- only struggles to fit some very rapid changes in the acceleration function, which likely lie -ustof-distributtin. We provide details of these computations and explore the systems further in Appendices F.1 and F.2, respectively.
The discussion above yielded a qualitative picture of the inference capabilities of FIM- To obtain a more quantitative perspective, we now compare F IM- against the ODEFormer model of d' A'cooliii al. (2024). ODEFormer is trained offline infer the vector field f(x) of dynamical systems in symbolic form. We estimate x() with ODEFormer by first solving their inferred ODE, to obtain an estimate of the ODE solution x(),, and then evaluating their inferred vector field along their esti- mated (t). As target dataset we analyse ODEBench, which was also introduced by d' d'Ascoli et al. (2024). ODEBench consists of 63 autonomous ODEs of different dimensionalities (specifically. 23 fine-grid
1D, 28 2D, 10 3D and 4D equations) and their solutions. The latter are represented of 512 points. We set = 0.5 (which defines point-wise missing patterns) and let y be either 0 or 0.05. We then compute the MAE of ODEFormer and FIM on the estimation of both xt() and x(t) with respect to the ground-truth trajectories across all 63 target ODEs. Table I reports our results averaged over all ODEs and shows that FIM-C outperforms ODEFormer in every case.
Some final remarks are in order. First, d' Ascoli et al. (2024) originally report n² scores in their paper. We also report these per dimension in Table 9 of Appendix F.3. Our conclusions remain un- changed. Second, these results also reveal that, despite being trained on additive noise only, FIM-& can handle multiplicative noise well. Third, we empirically demonstrate that F IM- is also superior to LatentODE trained on the iconic Lorenz systems in Appendix F.4. Fourth, we also demonstrate that FIM-& outperforms the very recent NeuralODEProcesee model on low-data regimes in Ap- pendix F.5. All together, these results reflect the capabilities of FIM- to not only impute point- wise missing data in dynamical systems, but also accurately reconstruct their phase portraits, both in zero-shot mode. Finally, Appendix I contains preliminary results on fine-tuning FIM.
4.2 IMPUTATION OF POINT-WISE MISSING PATTERNS
In this subsection we evaluate F FIM- on 8 real-world, high-dimensional datasets featuring point- wise missing patterns. Specifically, we study the Guangzhou dataset, which contains traffic speed records with 214 channels and 500 observations, and the Solar dataset, which consists of solar-power generation records with 137 channels and 52560 observations. We obtained the (preprocessed) datasets from Fang et al. (2024). We also study two popular air quality datasets from Beijing and Italy, which have 132 channels with 1458 observations, and 13 channels with 774 observations, respectively; the Electricity and ETT-h/ datasets of electricity consumption, common in forecasting studies, with 370 (1 (1457) and (358) channels (observations) each; and two additional traffic-related datasets, namely the PeMS dataset of road occupancy with 862 channels and 727 observations, and the single-channel Pedestrian dataset, which reports pedestrian activity in Australia and consists of 3633 observations. We obtained this second set of (preprocessee) datasets from Du et al. (2024).
After being split into train, validation and test sets, fifty percent of these subsets is randomly re- moved, defined missing and set aside for evaluation. We only make use of the (available 50% of the) test subsets with F FIM-f. Figure illustrates the type of zero-shot ODE solutions inferred by FTM-( on the Beijing and Guangzhou datasets (left and center). The bottom panels portrait instead

---

Published as a conference paper at ICLR 2025
Table 3: MAE (at missing values) on Motion Capture (MC) and Navier-Stokes datasets featuring temporal missing patterns of 20%. The large error bars in the MC dataset have been reported before (Heinonen et al., 2018). (Cubic) (Fl includes Savgol filter.
Motion Capture
Navier Stokes
Model
PCA
No PCA
PCA
No PCA
LatentODE
1.658 1 0.989
0,076 + 0.030
(Cubic) spline
3.362 1.175
4.209 = 1.436
0.085 i 0.003
0.083 0.003
(Cubic) spline(F)
2.897 1 0.871
2.998 - 0.881
0.084 1 0.000
0.075 1 0.002
FIM
1.765 + 0.627
1.611 = 0.453
0.063 = 0.003
0.051 + 0.002
the inferred (time) derivatives of the interpolating functions, together with the confidence of the model, which increases in regions with high, local information. Table 2 reports the average MAE the missing values for all models. Remarkably, FIM- outperforms all baselines in 4 out of datasets, comes second in one, and third in two, which suggests that there is indeed enough local information to perform the data imputatiio with "simple" functions in these use cases. We close this subsection by referring the reader to Appendix G, where we report the scores of additional, albeit less common baselines, as well as our results on the data splits investigated by Du et al. (2023).
4.3 IMPUTATION OF TEMPORAL MISSING PATTERNS
In this subsection, we look into the harder problem of imputing time series data featuring temporal missing patterns. Indeed, we explore the problem setup proposed by Heinonen et al. (2018), 8), in which (about) 20% of the data from the middle of the time series is missing completely. More precisely. we consider their human Motion Capture dataset, which consists of 43 trajectories, each with 50 channels and 100 observations. We also apply their setup the Navier-Stokes simulation of Course & Nair (2023), which instead contains 596602 channels. To be able to handle the (high) dimensionality of the datasets, Heinonen et al. (2018) projects the data to latent space using PCA, and trains their models to perform the imputation there (see Appendix H for details). We follow their methodology and train a LatentODE model, to outperform the results reported by Heinonen et al. (2018). Additionally, compare against naive spline interpolation.
Having set the stage. we leverage our pretrained FIM to infer the ODE solution that best imputes the missing data in PCA space, and report the MAE at the missing values, after projecting back to data space, in Table 3. Since FIM (and spline) can be applied to each channel independently, also report the MAE we obtain by imputing the data directly in high-dimensional (data) space. Again, one and the same F FIM performs comparably to (or even better than) LatentODE in the Motion Capture dataset, and outperforms it on the Navier-Stokes case. Note that imputing the data in high-dimensional space helps F FM avoid the errors introduced by the PCA projections. The right panel of Figure illustrates the class of zero-shot ODE solutions inferred by F FIM in the Motion Capture dataset, and demonstrates how FIM leverages the global, seasonal structures outside the gap to impute the (locally simpler) missing data. The agreement is strong.
CONCLUSIONS
In this work, we introduced a novel methodology for zero-shot imputation of time series data, whose underlying dynamics are assumed to be governed by ordinary differential equations (ODEs). We em- pirically demonstrated that one and the same Foundation Inference Model (FIM) can impute datasets of any dimensionality, featuring noise signals of very different nature, even in cases which do not naturally admit a description in terms of ODEs. In fact, we showed that FIM often outperforms SOTA models that are trained to the target distributions.
The main limitation of our methodology is clearly imposed by our synthetic distributions. Evaluating FIM on empirical datasets whose distribution significantly deviates from our synthetic distribution will inevitably yield poor estimates. The left panel of Figure 2 provides such an example. Indeed, for some initial conditions, the velocity of the Van der Pol oscillator features rapid changes that are not well represented by our pretraining distribution. Future work shall extend our decoupling of local and global features of the imputation model to the problem of zero-shot forecasting.
10

---

Published as a conference paper at ICLR 2025
6 REPRODUCIBILITY STATEMENT
The Appendix includes all information required to reproduce our presented methods and results. Let us therefore provide short overview of its contents, focused on the parts related reproducibility.
An important part of our methodology is our synthetically generated training dataset, which Ap- pendix B is dedicated to. The generation and hyperparameter choices are described extensively in Appendix B.1 (distributions over functions), Appendix B.2 (observation grids) and Appendix B.3 (noise processes). Appendix B.4 schematizes the whole generation algorithm.
Appendix C includes all details related to FIM-f, our model for point-wise missing pattern impu- tation. particular, Appendix C.2 describes the complete model architecture, including hyperpa- rameters of our single trained F FM-& model, and model inputs and outputs. Its training procedure is discussed in Appendix C.4, including hyperparameters for the optimizer and computing resources, whereas Appendix C.3 states the training objective.
The details for FIM, our model for temporal missing pattern imputation, are given Appendix D. Appendix D.2 describes the complete model architecture, where also include the hyperparame- ters used for our single trained IM model. The training details are given in Appendix D.4, including hyperparameters of the optimizer and computing resources, and Appendix D.3 states the associated training objective.
Equations of the evaluation metrics used for comparisons in experiments, and remarks about their usage in the different applications, are provided in Appendix E. Training details and com- puting resources for Laten(ODE, of our baselines, stated (for each dataset individually) in Appendices F.4.2, H.1.2 and H.2.2.
The datasets in our experiments are addressed individually in their own subsection in Appen- dices F, G and H. Each subsection contains either generation hyperparameters or links to their sources, and the applied pre-processing steps.
Finally, our pretrained model, repository and tutorials are available online
ACKNOWLEDGMENTS
This research has been funded by the Federal Ministry of Education and Research of Germany and the state of North-Rhine Westphalia as part of the Lamarr- Institut for Machine Learning and Artificial Intelligence.
REFERENCES
Sören Becker, Michal Klein, Alexander Neitz, Giambattista Parascandolo, and Niki Kilbertus. Pre- dicting ordinary differential equations with transformers. In International Conference on Machine Learning, pp. 1978-2002. PMLR, 2023.
Giancarlo Benettin, Luigi Galgani, and Jean-Marie Strelcyn. Kolmogorov entropy and numerical experiments. Physical Review A, 14(6):2338, 1976.
David Berghaus, Kostadin Cvejoski, Patrick Seifner, Cesar Ojeda, and Ramses Sanchez. Founda- 2024.
tion inference models for markov jump processes. arXiv preprint
Johann Bernoulli. Extrait de la Réponse de M. Bernoulli à M. Herman, datée de Basle le 7. Octobre 1710. Mémoires de Ac. Royale des Sciences, Boudot, Paris, 1710, 1712.
Rishi Bommasani, Drew A Hudson, Ehsan Adeli, Russ Altman, Simran Arora, Sydney von Arx, Michael S Bernstein, Jeannette Bohg, Antoine Bosselut, Emma Brunskill, et al. On the opportu- nities and risks of foundation models. arXiv preprint
2021.
Wei Cao. Dong Wang, Jian Li, Hao Zhou, Lei Li, and Yitan Li. Brits: Bidirectional recurrent
imputation for time series. Advances in neural information processing systems, 31, 2018.
Phttps://finee ienceeeneeeeeo www.eeeeiiieeeeeeeeeeeeo.oeeeeee .iteeiittti thtb... /0 IM/intro.html

---

Published as a conference paper at ICLR 2025
Stephanie Chan, Adam Santoro, Andrew Lampinen, Jane Wang. Aaditya Singh, Pierre Richemond, James McClelland, and Felix Hill. Data distributional properties drive emergent in-context learn- ing in transformers. Advances in Neural Information Processing Systems, 35:18878 11888888-1899) 2022.
Zhengping Che, Sanjay Purushotham, Kyunghyun Cho, David Sontag, and Yan Liu. Recurrent neural networks for multivariate time series with missing values. Scientific reports, 8(1):6085, 2018.
Ricky TQ Chen, Yulia Rubanova, Jesse Bettencourt, and David K. Duvenaud. Neural ordinary differential equations. Advances in neural information processing systems, 31, 2018.
Kevin Course and Prasanth Nair. State estimation of a physical system with unknown governing equations. Nature, 622(7982):261-26)) 2023.
Stéphane d' Ascoli, Sören Becker, Philippe Schwaller, Alexander Mathis, and Niki Kilbertus, ODE- Former: Symbolic regression of dynamical systems with transformers. In The Twelfth Interna- tional Conference on Learning Representations, 2024.
Wenjie Du, David Côté, and Yan Liu. Saits: self-attention bassed imputation for time series. Expert Systems with Applications, 219:119619, 2023.
Wenjie Du, Jun Wang, Linglong Qian, Yiyuan Yang, Fanxing Liu, Zepu Wang, Zina Ibrahim, Haoxin Liu, Zhiyuan Zhao, Yingjie Zhou, Wenjia Wang, Kaize Ding, Yuxuan Liang, B. Aditya Prakash, and Qingsong Wen. Tsi-bench: Benchmarking time series imputation. arXiv preprint arXiv:2406. 12747, 2024.
Emilien Dupont, Arnaud Doucet, and Yee Whye Teh. Augmented neural odes. Advances in neural information processing systems, 32, 2019.
Harrison Edwards and Amos Storkey.
Towards neural statistician.
arXiv
preprint
arXiv:1606.02185, 2016.
Shikai Fang, Qingsong Wen, Yingtao Luo, Shandian Zhe, and Liang Sun. BayOTIDE: Bayesian on- line multivariate time series imputation with functional decomposition. In Proceedings of the 41st International Conference on Machine Learning, volume 235, pp. 12993-13009 PMLR, 2024.
Vincent Fortuin, Dmitry Baranchuk, Gunnar Rätsch, and Stephan Mandt. Gp-vae: Deep probabilis- tic multivariate time series imputation. In Proc. AISTATS, pp. 1651-1661, 2020.
Marta Garnelo, Dan Rosenbaum, Christopher Maddison, Tiago Ramalho, David Saxton, Murray Shanahan, Yee Whye Teh, Danilo Rezende, and SM Ali Eslami. Conditional neural processes. In International conference on machine learning, pp. 1704-1713, PMLR, 2018a
Marta Garnelo, Jonathan Schwarz, Dan Rosenbaum, Fabio Viola, Danilo J Rezende, SM Eslami, and Yee Whye Teh. Neural processes. arXiv preprint 1rX.v:180..01622, 2018b.
Lu Han, Han-Jia Ye, and De-Chuan Zhan. The capacity and robustness trade-off: Revisiting the channel independent strategy for multivariate time series forecasting. IEEE Transactions Knowledge and Data Engineering, 2024.
Markus Heinonen, Cagatay Yildiz, Henrik Mannerström, Jukka Intosalmi, and Harri Lähdesmäki. Learning unknown ode models with gaussian processes. International conference on machine learning, pp. 1959-1968 PMLR, 2018.
Luke B Hewitt, Maxwell Nye, Andreea Gane, Tommi Jaakkola, and Joshua B Tenenbaum. The variational homoencoder: Learning to learn high capacity generative models from few examples. arXiv preprint arXiv: 1807.0881 2018.
Jacobus Henricus van 't Hoff. Studies in Chemical Dynamics. Amsterdam: F. Mueller and Co. London: Williams and Norgate, 1986.
Edwin Thompson Jaynes. Probability Theory: the Logic of Science. Cambridge University Press, 2003.
12

---

Published as a conference paper at ICLR 2025
Xiajun Jiang, Ryan Missel, Zhiyuan Li, and Linwei Wang. Sequential latent variable models for few-shot high- high-dimensional time-series forecasting, In The Eleventh International Conference on Learning Representations, 2023.
Hyunjik Kim, Andriy Mnih, Jonathan Schwarz, Marta Garnelo, Ali Eslami, Dan Rosenbaum, Oriol Vinyals, and Yee Whye Teh. Attentive neural processes. arXiv preprint arXiv::901.05700 2019.
Günter Klambauer, Thomas Unterthiner, Andreas Mayr, and Sepp Hochreiter. Self-normalizing neural networks. Advances in neural information processing systems, 30, 2017.
Nikola Kovachki, Zongyi Li, Burigede Liu, Kamyar Arizzadenesheli, Kaushik Bhattacharya, An- drew Stuart, and Anima Anandkumar. Neural operator: Learning maps between function spaces with applications to pdes, Journal of Machine Learning Research, 24(89):1-97, 2023.
Balaji Lakshminarayanan. Alexander Pritzel, and Charles Blundell. Simple and scalable predictive uncertainty estimation using deep ensembles. Advances neural information processing systems, 30, 2017.
Hugo Larochelle, Dumitru Erhan, and Yoshua Bengio. Zero-data learning of new tasks. In AAAI, volume 1, pp. 3, 2008.
Edward N Lorenz. Deterministic nonperiodic flow. Journal of atmospheric sciences, 2002);300-14,, 1963.
Ilya Loshchilov and Frank Hutter. Sgdr: Stochastic gradient descent with warm restarts. arXiv preprint arXiv:1608.00 03993, 2016.
llya Loshchilov and Frank Hutter. Decoupled weight decay regularization. arXiv preprint
arXiv:1711.05101, 2017.
Alfred James Lotka. Elements of physical biology. Williams & Wilkins, 1925.
Lu Lu, Pengzhan Jin, Guofei Pang, Zhongqiang Zhang, and George Em Karniadakis. Learning nonlinear operators via deeponet based on universal approximation theorem of operators. Nature machine intelligence, 3(3):218-229 2021.
Isaac Newton. Philosophiae Naturalis Principia Mathematica. ef. especially Liber de Motu Sects. II, III, VII. London, 1687.
Yuqi Nie, Nam H Nguyen, Phanwadee Sinthong, and Jayant Kalagnanam. A time series is worth 64 words: Long-term forecasting with transformers. In The Eleventh International Confer- ence on Learning Representations, 2023. URL /ttp:://///eeeeeeeee/e ew. net htt./formmmiie JbdcOvTOcol.
Bernd R Noack, Konstantin Afanasiev, Marek Morzyúski, Gilead Tadmor, and Frank Thiele. A hierarchy of low-dimensional models for the transient and post-transient cylinder wake. Journal of Fluid Mechanics, 497:335-363, 2003.
Alexander Norcliffe, Cristian Bodnar, Ben Day, Jacob Moss, and Pietro Liò. Neural {ode} processes. In International Conference on Learning Representations, 2021. URL https:
//openzeview. wttp://wwwwwee.tttneeeeeeoee wtwee/foomm/iieeeeee att.e/oooommm1d222//yylal
Norman H Packard, James P Crutchfield, J Doyne Farmer, and Robert S Shaw. Geometry from time series. Physical review letters, 45(9):712, 1980.
Brooks Paige and Frank Wood. Inference networks for sequential monte carlo in graphical models. In International Conference on Machine Learning, pp. 3040-3049, PMLR, 2016.
Yulia Rubanova, Tian Qi Chen, and David Duvenaud. Latent Ordinary Differential Equations for Irregularly-Sample Time Series. In Advances in Neural Information Processing Systems 32, pp- 5321-5331, 2019.
Abraham Savitzky and Marcel JE Golay, Smoothing and differentiation of data by simplified least squares procedures. Analytical chemistry, 36(8):1627-1639, 1964.
13

---

Published as a conference paper at ICLR 2025
Patrick Seifner and Ramses Sanchez. Neural Markov jump processes. In Proceedings of the 40th International Conference on Machine Learning, volume 202. PMLR, 2023.
Patrick Seifner, Kostadin Cvejoski, David Berghaus, Cesar Ojeda, and Ramses Sanchez. Foun- dation inference models for stochastic differential equations: A transformer-based approach for zero-shot function estimation. arXiv preprint arXiv:2502. 19049, 2025.
Satya Narayan Shukla and Benjamin Marlin. Multi-time attention networks for irregularly sampled time series. In ICML Workshop on the Art of Learning with Missing Values (Artemiss), 2020.
Gautam Singh, Jaesik Yoon, Youngsung Son, and Sungjin Ahn. Sequential neural processes. Ad- vances in Neural Information Processing Systems, 32, 2019.
Andreas Stuhlmüller, Jacob Taylor, and Noah Goodman. Learning stochastic inverses. Advances in neural information processing systems, 26, 2013.
Yusuke Tashiro, Jiaming Song, Yang Song, and Stefano Ermon. Csdi: Conditional score-based diffusion models for probabilistic time series imputation. Advances in Neural Information Pro- cessing Systems, 34:24804-2453 2021.
Matias Valdenegro-Towe Toro and Daniel Saromo Mori. A deeper look into aleatoric and epistemic uncer- tainty disentanglement. In 2022 IEEE/CVF Conference on Computer Vision and Pattern Recog- nition Workshops (CVPRW), pp. 1508-1516. IEEE, 2022.
Jean-Philippe Vert, Koji Tsuda, and Bernhard olkoo.. primer on kernel methods. Kernel methods in computational biology, 47:35-70, 2004.
Vito Volterra. Variazioni fluttuazioni del numero d' d'individui in specie animali conviventi, vol- ume 2. Societá anonima tipografica" Leonardo da Vinci", 1927.
Jack M Wang, David Fleet, and Aaron Hertzmann. Gaussian process dynamical models for human motion. IEEE transactions pattern analysis and machine intelligence, 30(2):283-298, 2007.
Jun Wang, Wenjie Du, Wei Cao, Keli Zhang, Wenjia Wang, Yuxuan Liang, and Qingsong Wen. Deep learning for multivariate time series imputation: A survey. arXiv preprint arXiv:2402.04059, 2024.
Rui Wang, Robin Walters, and Rose Meta-learning dynamics forecasting using task inference.
Advances in Neural Information Processing Systems, 35:21640-21653. 2022.
Christopher Williams and Carl Rasmussen. Gaussian processes for regression. Advances in neural information processing systems, 8, 1995.
Cagatay Yildiz, Markus Heinonen, and Harri Lahdesmaki. Ode2vae: Deep generative second order odes with bayesian neural networks. Advances in Neural Information Processing Systems, 32, 2019.
14

---

Published as a conference paper at ICLR 2025
A COMPARISON WITH OTHER META-LEARNING APPROACHES
(The discussion that follows took place during the rebuttal process].
Conditional neural network models such as the Neural Statistician (Edwards & Storkey, 2016; He- witt et al., 2018), or members of the Neural Process family (Garnelo al., 2018b;a; Kim et al,, 2019) and their later extensions to dynamical process (Singh al., 2019; Wang et al., 2022; Jiang et al., 2023), consider the problem of (meta)training a single model on a set of different, albeit related datasets Di...., Dm. each of which is characterized by corresponding (that is, shared) context latent variable C1,...,Crm Cm-
Main difference with FIM:
The first and most important difference between all these approaches and ours is that they need to be trained on datasets from their target domains, which makes both their inferred representa- tions and optimized weights "problem specific". Indeed, every one of these works focuses, by construction, meta-learning among similar systems (or datasets) only. To illustrate, let con- sider the work of Jiang et al. (2023). They first forecast bouncing ball simulations under different gravity conditions. For them, each gravity corresponds to different dataset, and their inferred representations must encode their corresponding gravities (or at least some aspects of them). Later, they consider forecasting turbulent flow dynamics under different buoyant forces and train their model anew, in order to encode the forces characterizing each dataset. Thus, the parameters of the neural network building of their model are not the same for the bouncing ball experiment and the turbulent flow experiment. Similarly, the semantic content of the shared representations inferred by their model not the same for the bouncing ball experiment and the turbulent flow experiment. In the former case it corresponds to gravity, whereas in the latter it corresponds to buoyant force.
In sharp contrast, proposal is not problem specific, because it entirely relies on the gen- eral assumptions (inductive biases, or prior knowledge) of Section 3. These simple assumptions pertain the nature of interpolating functions only, and are therefore independent of the actual data generation mechanisms underlying our target datasets. Given that our model is only trained on synthetic dataset encoding these two assumptions, can only recognize patterns that help it reconstruct the best interpolating function given the available (context) data, regardless of the underlying data generation mechanism. In other words, our method maintains the same network parameters and representation semantics throughout all experiments. The representations consis- tently correspond to the time derivative and initial conditions of the hidden interpolating function, regardless of the target dataset.
2. A second important difference lies in how these works leverage their prior knowledge for transfer. Specifically, they assume that their target domains allow for the collection their different yet related datasets D,,... Dm, on which to train their (meta)models.
We instead fully rely on our two assumptions of Section equation 3, and therefore on our synthetic dataset which encodes them.
3. A third difference is that all these works assumed there is shared latent variable em encoding the main features of the mith dataset D,, in their collection. Most of these works rely on neural variational inference optimize their encoder-decoder pairs, and infer their representations in an unsupervised manner (see however Wang et al. (2022)).
In our proposal our latent variables are not shared, for each time series has associated to it a single interpolating function. Transfer learning takes place because of the general assumptions we encode into our synthetic dataset.
B SYNTHETIC DATA GENERATION MODEL: SPECIFICS
B. ON THE DISTRIBUTION OVER PARAMETRIC FUNCTIONS OF TIME
In this section we define two distributions over the space of parametric functions. We use the first one to train FIM for interpolation. We use the second one to train extended FIM for imputation.
15

---

Published as a conference paper at ICLR 2025
B.1.1 DISTRIBUTION OVER "SIMPLE" INTERPOLATION FUNCTIONS
To train our model for point-wise missing patterns, define p(f(n)) either via Gaussian processes (GP) with Radial Basis Function (RBF) kernels (Williams & Rasmussen, 1995; Vert et al., 2004), or as truncated Chebyshev expansions, whose coefficients and degree are both randomly sampled. Each alternative has different Ty hyperparameters.
In the case of RBF kernels, there is single free hyperparmeter, the lengthscale Nif, which controls the scale at which variations take place. A small Mi results in parametric functions with short-range fluctuations, whereas larger T/f produce smoother functions that capture broader trends. We define p(ns) as mixture of two Beta distributions, with equal mixing coefficients. To wit
1 p(r)) (eta 12, 10) (eta
Here, the first component returns functions with faster change (i.e. higher frequency), whereas the second component returns smoother functions.
In the case of (the M M-order) truncaled Chebyshev expansions, write
M
f(t)
= Σ m=1
am Tm
(10)
where Tm (l) is the mth Chebyshev polynomial with (real) coefficient Om. Now, in order to generate a random parametric function with Eq. 10, we sample both the degree M and the set coefficients fai...., av) } from the prior p(ny). ). Intuition says would like high-order polynomials to occur rarely. It also says one would like the scale of their coefficients to be small (see also the work of e.g. Chan et al. (2022)). We therefore (implicitly) define the distribution p(ny) over hyperparameters as
aM) (aI....,am) ~N N (1, 1 )T with M Zipf(2).
(a), (01.-..,aM))
(11)
In practice, we generate a synthetic dataset of 2M parametric functions (200K for the test set), each of which is evaluated on a fine grid of Lmaz = 128 points. This fine grid is defined regularly on the unit interval [0, 1]. Half of this dataset consists of random Chebyshev expansions. The other half consists of parametric functions sampled via GPs.
B.1.2 DISTRIBUTION OVER IMPUTATION FUNCTIONS WITH GLOBAL PATTERNS
To train our model for temporal missing patterns, we require data with trends and seasonality, such that model learns to these (global) patterns to impute the missing data. We hence define p(f[ny) as again via GPs, but opt for Periodic kernels, which exhibit the required seasonality.
Periodic kernels are specified by two free hyperparameters: the lengthscale, denoted by n/y, and the period, denoted by rif- The lengthscale determines the local fluctuations of the function, as described in Appendix 1.1.1. The period determines the frequency of repetitions of the function, i.e. the seasonality.
We sample both hyperparameters independently. In other words, our distribution over Periodic kernel GPs p(ns) factorizes as p(nf) p( rff) p(n)))((f))) We define uniform distributions
for both hyperparameters:
p(n))) = u(r) (0.)5. p(nf)) = 21 [[033, 0.5))
We generate dataset of 500k parametric functions (50K for the test set). Each function is evaluated on a regular fine grid of Lmmx = 256 points on the unit interval [0, 1].
Note that we do not need to explicitly address trends in the parametric functions over time. A sample from the defined GPs will (a.s.) have a non-zero mean. During data generation (see Ap- time, such that the corresponding ODE solutions
pendix B.4), these deviations will accumulate exhibit trends naturally. In combination with instance normalization during model input processing (see Appendix D.1), 1), this setup covers variety of trends.
16

---

Published as a conference paper at ICLR 2025
B.2 ON THE DISTRIBUTION OVER OBSERVATION GRIDS
In this section define distributions observation grids. We the first train FIM for interpolation. We use the second one to train our extended FIM for imputation.
B.2.1 OBSERVATION GRIDS ENCODING POINT-WISE MISSING PATTERNS
We define the prior distribution over observation grids (Tna(T1...... Tb ny) in a way that it allows for both regular and irregular observation grid instances, each with at most Lmm and at least Luin observations. Note that the latter number defines the minimum number of context points needed for FIM to function.
Let us recall that our parametric functions f(t) are evaluated on fine grid of 128 points (see Ap- pendix B.1 above) over the unit interval [0, 1]. This fine grid is subsequently subsampled randomly, to define the random observation times ....... Tr.
We employ two subsampling schemes regular and irregular which occur with equal probability within our entire dataset.
Regular scheme. In this scheme the observation times are obtained from the fine grid by strides of regular length, where the stride length sampled from the (ni,orm[[12,,, 16)) distribution. Note that this distribution defines p(Ng)-
Irregular scheme. In this scheme the observation times are instead defined via Bernoulli masks, whose survival probabilities (here denoted by ngy) are sampled from categorical distribution p(ng): The latter is defined over the (0.0625, 0.25, 0.5} with class probabilities f0.5, 0.25, 0.25], re- spectively.
Note that ensure that there are at least Lmin observations in each time series, by rejecting observation grid instances with less than eight points. In contrast, there are time series with at most T 128 observations.
B.2.2 OBSERVATION GRIDS ENCODING TEMPORAL MISSING PATTERNS
To train our model for temporal missing patterns, the observation grids of training data must exhibit such patterns. We sample such grids in a two-step process.
In the point-wise step, for the observation grids sampled in Appendix B.2.1, only employ posnt-wise subsampling schemes. These are required for our model can handle irregular grid time series data of variable lengths.
In the temporal step we introduce the temporal missing patterns, by dropping consecutive range of observations. We train our model to impute these dropped values.
More formally, let т, denote times in the point-wise observation grid and 7, denote times in the temporal missing pattern observation grid. Then their joint distribution depends on hyperpa- rameters, sppi and r/mp and factorizes as
Pana ( 3k. T1
T1,
11
jamp
11g 1ppii
Premp (
TI. ( y(1 P Ppain (TI...., " W(
where ( S (Ti)L
Let us now report the specifics for both steps:
Point-wise step. One half of the observation grids are generated with regular subsampling scheme, sampling stride lengths from Uniform([L, 2, 3, 4)). The other half is generated with irregular subsampling scheme, sampling from a Bernoulli distribution with survival probability 0.5.
Temporal step. To generate missing pattern in the time series, we first sample the position of the missing pattern from a uniform distribution, specifically Uniform(1, 3), which determines the window that will be masked. The length of the missing pattern is randomly sampled from uniform distribution in the range [10, 30). The remaining part the time series then split into four equal windows.
17

---

Published as a conference paper at ICLR 2025
B.3 ON THE DISTRIBUTION OVER NOISE PROCESSES
Let x R4 + R be time-dependent trajectory. We define the distribution over noise processes on observations s(() R of such trajectory to factorize, such that the noise process is fixed for each trajectory. More formally, denoting the noise process by a and the noisy trajectory by y, we assume
Paoise (y(t), rtt), (t)) = (y(() ært), p)(o)) for all t E R.
As outlined in Section 3.1, we use additive Gaussian noise
Pasise (y(() xxt), a) (x((t),
for our synthetic datasets. Here, the noise process is defined only by single value: the standard deviation of Gaussian distribution. Hence, our choice of distribution over noise processes reduces to a choice of suitable distribution for σ. We use p(o) N (0, X) in our synthetic datasets, interpreting negative samples for as their absolute value
The value of À differs between our two synthetic datasets. In our dataset for point-wise missing patterns, we use = 0.1 and our dataset for temporal missing patterns, we use 0.05. We reduce À for the temporal-missing patterns retain more structure of the trajectory, such that the model can learn to impute based on the structure of the time series.
B.4 ON THE GENERATION OF THE SYNTHETIC DATASET
In this subsection schematize the data generation process. Underlying this process is sampling procedure of noisy time series data, based on samples random functions, observation grids and noise processes. The corresponding distributions are described in Appendices 1-I...
To generate the jth instance of our synthetic datasets, we utilise these distributions in the following generation steps:
1. Sample a function
J3 p(flnss) where Mfj ~ p(ns): and record its values (f)(4))E on the fine grid (t:) the regular grid on the unit interval [0, 1].
2. Sample initial value roj NN0 1) that defines the ODE solution
(t) = Toj + S. fz(s)ds
and record its values {t)(t) on the fine grid.
3. Sample observation grid
Tljs Tú Posa(TI, (os), where Va. p(1g))
that defines a time series with (clean) observations ((rij = a(T;); Top) of the ODE
solution.
4. Sample noisy observations
Vij ~ Pussee Fijs Gz); where of -j~ p(o), of the ODE solution that define the time series {wy Tyy
Note that each step consists of i.i.d. samples of their respective (hierarchical) distributions. Thus, each instance of our synthetic datasets is also i.i.d.
C FIM FOR LOCAL INTERPOLATION: ADDITIONAL DETAILS
In this section, we begin with instance normalization in section C.I as pre-processing step for handling times series of varying scales. We then outline the training objective section C.3, focusing on the loss function design. Next, we discuss methods for processing time series of any length and dimensionality (section C.6).
This distrihution also called folded Gaussian distribution.
18

---

Published as a conference paper at ICLR 2025
C.1 INPUT PRE- PREEPRRCESSING AND OUTPUT POST-PROCESSIN
To accommodate time series of all scales, we employ min-maax-normaliztion for the observation
values and times per time series before processing them with FIM-& and renormalize the model outputs accordingly. To express this instance normalization in equations, let (y1, (80,,),,,,)) Ti) be
noisy observations, following the notation of Section 3.2.1 and denote
min
Tiax
max i=
Ts
Vmin
min s=
Mir Umax
max yi
(12)
v=,....
Before applying FIM-( l described in Section 3.2.1, we replace the inputs by their normalized values:
yi
Mmin
Tmin
y/
(13)
Ymax
Y/min
Tmax
7min
Let f(t), log Var (t). Zo and Var(&o) be the outputs of FIM-L, following Equations and 5 of Section 3.2.1, given the normalized inputs. Under the renormalization maps
y (ymax Ymin ) + Ymin+ f
(Tmax
Tmin + Tmin
(14)
the model outputs are transformed to the original time and value scale as follows:
Î(1)= < f(1)
Mmax Mmin
(15)
Tiax
Tmin
log (ar(f)(t) + log Var + 2 log
Imax
Mmin
(16)
Tmax
Tmin
to To(ymax - Y/min) + ymin
(17)
log Var($o) < log Var( ) + log(Umax ymin)
(18)
While the transformations of i.g and log Var( follow immediately from the value renormalization map, the transformations of f()) and (ar( follow from the linearity of the derivative (in case
of the value renormalization map) and the chain rule (in case of the time renormalization map).
This described instance normalization approach enables FIM- to be used a range of (real-world) applications, we show in experiments (Section 4). But there are also limitations, which we will now address.
Firstly, some signal-to-noise ratio patterns are not accurately resolved because of this normalization. For example, observations (almost) constant functions with additive noise are not interpolated very accurately, because the normalization parameters mainly determined by the noise, not the underlying signal.
Secondly, the internal embeddings of FIM-( do not contain any information about the scale of the imput. Downstream tasks, like imputation of temporal missing patterns could benefit from this information. In such scenarios, opted for separate embeddings of statisties about the input data scales, to reintroduce them alongside the F FIM-& embeddings (see Appendix D.2).
C.2 MODEL ARCHITECTURE
Let us provide more details about the interpolation model architecture, which was already outlined in Section 1, and hyperparameters of the single FIM- model all experimental results were derived with. Note that we use SeLU (Klambauer et al., 2017) as activation for feed-forward neural networks and employ dropout rate of 0.1.
In total, FIM- has roughly 20 million learnable parameters. The results of an ablation study on architecture and hyperparameter choices are described in Appendix C.5.
Model inputs:
noisy time series (y1, " (31, T.) with observation values yi E R and ordered, but potentially irregular, observation times E R4 with Tı ... < T/.
iii) A time t R4 at which evaluate the functions f and log Var(f)
Model architecture:
19

---

Published as a conference paper at ICLR 2025
(i) Temporal embedding. We the learnable time embedding from Shukla & Marlin (2020) with output in 12512 Its ith dimension is defined as
(0(() wot bo if 0 sin(w,t b,) otherwise
where w; and b; are learnable parameters.
(ii) Trunk net equivalent: We use the time embedding of in combination with feed-forward neural network 0% with 4x1024 hidden layers and output in R to encode the evaluation time t. The composition of % can thus be understood as the Trunk DeepOnet. Let us denote
t. (0f (1)
(iii) Individual observations embedding. us denote the ith element of the time series in our input with
y Concat(ys, (T)))
(iv) Branch net equivalent. The embedded time series processing network, which can be under- stood as the branch net of DeepONets, is a composition of sequence processing network of and a feed-forward neural network 02. The sequence processing network vi is defined as bi-directional LSTM with hidden states in RR512, The feed-forward neural network @b with 4x 1024 hidden layers and outputs in R512 Let us denote
u pf(vf(y?......)))
(v) Time derivative projection. We combine the output of the trunk net and branch net equivalent networks via seed-forwane neural network oy with 4x 1024 hidden layers and outputs in RR522 Let us denoted its output by
(1) og (Concat ,,²)).
The Gaussian distribution over the values of the estimated time derivative is parameterized by two linear projections from (t) to R, denoted by
f(1) o (h (t)) and log Var(f) (h"
(vi) Intial condition projection: Parameters of the Gaussian distribution for the initial condition
sg = o (1 and log Var(ão) = (f("1)
are modeled with feed-forward neural networks og and df each with 4x 1024 hidden layers and outputs in R.
Model outputs: To summarize, our model returns
a time derivative f(t), with uncertainty estimate log Var(f) for time t E R, and iii) initial value žo, with uncertainty estimate log Var($o).
Combined, these outputs define a function 2(r) interpolating the noisy (input) time series (I..T)...... T)::
#(T) to + S f(s)ds
C.3 TRAINING OBJECTIVE
The parameter set of FIM-E is trained in a supervised fashion, utilising the synthetic dataset for point-wise missing patterns introduced in Appendix B.
Following the notation of Appendix B.4, let f the sampled function from one instance in synthetic dataset, { its values on the fine grid, zo the initial value of the ODE solution, recorded on the fine grid {x(t+)) and {Vi, ri)!_ the associated noisy time series. Using ((,, T inputs for FIM-f, denote its outputs by f. log Var(f), to and log Var($o) as in Appendix C.2.
The training objective of FTM-f consists of three parts:
20

---

Published as a conference paper at ICLR 2025
Maximizing the Gaussian log- ligllikeiihood log C(f(t.) f(t.), Var(f)(t:)) of the time derivative at all L points t; of the fine grid.
(ii) Minimizing the one-step-ahead reconstruction error (.r(4;+1) - (x((+) + +(())))) ((4)((-11 1:))) of the ODE solution at - points on the fine grid.
(iii) Maximizing the Gaussian oogllikelihood log L(ro Lo, Var(fo)) of the initial value.
Expressed as equation, train F PIM- E to minimize
(f(t)) f(t.)))
C
+ log loog(Vai (aa(Var() (va(())ttt))
(
2Var(
E (~p(snt)) memp(zt)
|x(t++1) - (x((s) + f(t))ttt) +.mi
1=1
(I) 20)²
E F0mpp(50)
+ log (og(Var((o)))
2 2Vr (t)
Note that the one-step-ahead reconstruction error is simply the absolute error of single step of the Fuler method, when using the inferred time derivative f(44), but the ground-truth starting point at4)). Such term is only viable in our supervised training regime on synthetically generated data, where we have access the ground-truth ODE solution
C.4 TRAINING PROCEDURE AND IMPLEMENTATION
The parameters 0 of FIM- - were optimized with AdamW (Loshchilov & Hutter, 2017), using learning rate of 16 and weight decay 1e-1, Using a batch size of 1024, the loss on a validation set converged after approximately 500 epochs.
We used four A100 80GB GPUs to train FIM-(. The implementation is done in Jax². Its code and the trained model weights are provided in the supplementary material.
C.5 ABLATION STUDIES: ARCHITECTURE AND DATASET DESIGN
In this section, we present ablation studies of the design of the synthetic training dataset dataset for and the architecture of FIM-f. We evaluate the ablations by the Ra-accuracy the ODEBench dataset, corrupted with = 0.5, o = 0.05 (see Appendix F.3 for more details). It is important to that these models trained only on synthetic data, and not on the ODEBench dataset, in accordance with our zero-shot approach.
First, experiment with the size of model (2M, 20M, 50M parameters) and the number of trajectories in the the training dataset (2M, 8M trajectories) and summarize the results in Table For 2M trajectories, increasing the number of parameters from 2M to 50M does increase the performance by around 4%. For 8M trajectories, the trend persists with the 20M parameter model, reverts with the 50M parameter model, resulting no significant improvement. The architecture specifics for the models of different sizes follows:
(i) 2M parameters: MLP: 4x256, hidden-dim: 256;
(ii) 20M parameters: MLP: 4x1024, hidden-dim: 512:
(iii) 50M parameters: MLP: 3x2048, hidden-dim: 1024.
Next, we experiment with the architecture of the sequence processing network vf. considering BiL.TMM and transformer networks. The architecture specifics in this case are:
Transformer: MLP: 4x256, N-of-layers: 4, QKV-dim: 256, N-of-heads 8, output dim: 256, total-parammcte--coutt:: 20M:
21

---

Published as a conference paper at ICLR 2025
iii BiLSTM: MLP: 4x1024, hidden-dim: 512 (i.e. 256 per each directions), total Ootal-parameter count: 20M.
The results are presented in Table 5 and show that there no significant difference in performance when evaluated on the ODEBench dataset.
Table 4: I8 accuracy with standard devia- tion, ablating parameter count and number of time series train set.
Table 5: R² Accuracy with standard devia-
tion, ablating the architecture of the sequence processing network v9.
Parameters
Time Series in Train Set 2M
Architecture
R² Accuracy (%)
8M
BiLSTM
87.30 1 1.35
2M
82.6 = 0.7
82.3 + 0.9
Transformer
87.38 + 1.42
20M
87.3 1.3
86.2 1.0
50M
87.7 = 2.0
84.1 + 1.3
An additional study we conducted involves varying the size of the maximum and minimum training contexts points in time series of the synthetic train set, which denote by Lmax and Lmi- The results for ablating Lma, presented in Table 6, indicate that, at least on the ODEBench dataset, the model performs better with more context points. However, the results are within error bars. Similarly, the results for ablating Luais shown in Table 7, reveal that the model performs better when only trained on time series with more observations. Again, this result is only valid for the ODEBench dataset, as other applications require small context windows, e.g. because there are too few observations in the time series.
Table 6: R² Accuracy with standard devia- tion, ablating Lixx. the maximum number of observations in a time series of the synthetic train dataset.
Table 7: R² accuracy with standard devia- tion, ablating Lmins the minimum number of observations in a time series of the synthetic train dataset.
Lmnx
R² Accuracy (%)
Luan
R Accuracy (%)
128
87.3 + 1.3
87.3 + 1.3
64
86.7 + 1.2
32
89.0 = 1.0
C.6 PROCESSING TIME SERIES OF ANY LENGTH AND DIMENSIONALITY
As a zero-shot inference model, FIM- should be able to cope with dynamic phenomena of very diverse nature. In particular, it should handle observations with different dimensionalities, lengths and scales. Our approach handle different scales with instance normalization was discussed in Appendix C.I. Now want to address our approach to the other two phenomena and offer some limitation of these approaches.
C. COMPOSITION ACROSS DIMENSIONS
To handle different (in fact, arbitrary) dimensional observations, we process each feature of a time series separately with our pretrained model FIM-(. The synthetic training dataset only contains ID time series (see Appendix B. and the model architecture only returns ID time derivatives (see Appendix C.2). Such channel independent strategy has been used previously, e.g. by Nie et al. (2023) and Han et al. (2024).
Regarding the inference problem of the time derivative on ODE solution path, this reduction to ID systems is exact, not an approximation. Indeed, consider a vector field f R+ x RD - IR" and the solution vector field along the solution path, f(t, x(t)) is a purely demendennn function, which naturally splits into D independent coordinate functions (f (t)) ... fn (t. x(t)).
R - RD to some initial value problem (rp, f). Then its time derivative, i.e. the
Coordinates of datasets generated from real-world scenarios, which we consider in our experiments, might not be independent or uncorrelated. Still, we show that approach can perform well in
22

---

Published as a conference paper at ICLR 2025
riety of such scenarios. However, extreme cases like the PHYSIONET2012 dataset (Du al., 2024), with very sparse observations in some coordinates, are a natural limitation channel independent approaches like ours.
C.6.2 COMPOSITION ALONG TIME
Time series in the training set of FIM-& are bounded in length by fixed upper bound Lmax. as described in Section 3.1. Therefore, FIM-& can only (reasonably) process time series which lengths do not exceed this upper bound. Similarly, FIM-f can only (reasonably) infer time derivatives approximately contained in the broad distributions over parametric functions of time, specified in Appendix B.1.
Datasets from real-world scenarios will likely not adhere to these limitations derived from our syn- thetic training dataset.
To adapt the pretrained FIM-& to handle such scenarios, we employ a windowing scheme. Con- cretely, split the time series into successive, overlapping windows that are still processable by the FIM-&. f. Applying (instances of) FIM-& to each window individually yields local estimates of the solution ((t).
As windows are processed separately, instance normalization (see Appendix C.1) and corresponding renormalization occurs individually at their respective scales.
To get a global estimate out of the local estimates, one can combine them on the overlaps with two approaches.
Interpolating the ODE solutions. Consider two successive, overlapping time windows A and B, defined on the intervals [t ft'l and [t t respectively, that
10 B < <
Applying instance of F to each window yields local solutions at( and t (t) in the respective windows. By defining
(1),
i 18 if st
a(t)
4 (t(
+ (
1- (t),,
(19)
(t).).
if t
for E [16 (f]. we combine (1) and (() to a global solution x(1) with linear inter- polation on the overlap. The time derivative corresponding to â can be found (a.e.) by taking the derivative of equation 19.
2. Interpolating the time derivatives and solving the ODE. Alternatively we could interpolate the estimated time derivatives f*(() and (!) of the time windows A and B, just as done in Eq. 19 but with â replaced with f. Once such interpolated time derivative is available one can integrate the ODE over the complete interval [t tt.. However, this second approach accumulates the errors of the ODE solutions of each interval.
In practice, and for all experiments described in Section 4 and the Appendix, we used the first approach.
Notation. There are two methods for specifying the overlapping windows:
1. Windows specified by the number of windows that cover the observations. For windows we write FIM- FIM-E(w.n. = m). identifying this method.
2. Windows specified by the number of observations they contain. For observations per windows, we write F FM FIM {(o.n. =m), identifying this method.
Limitations. Our experiments in Section 4 show that compositional approach along time can be quite effective. Still, there are natural limitations, which we want address.
parametric functions of time are broad, they do not cover all local patterns exhibited by real-world time series. For example, FIM-& struggles with high- frequency time series with sparse observations. Many windows required to accurately interpolate
Firstly, although distributions
23

---

Published as a conference paper at ICLR 2025
such patterns. Future work will explore the design of the synthetic training data distributions, to cover more of these patterns with F FIM-
Secondly, (subsequent) temporal missing patterns can not be accurately imputed with the window- ing approach. As windows are processed individually, and their outputs merged manually on the overlaps, information is passed between them. Therefore, imputation that requires global, long range patterns can not be resolved. We handle the temporal missing pattern imputation with an additional imputation module FIM (see Appendix D).
D FIM FOR INTERPOLATING TEMPORAL MISSING PATTERNS: ADDITIONAL DETAILS
D.1 INPUT
AND OUTPUT POST-PROCESSING
Analogous the local interpolation model discussed in Appendix C, our zero-shot model for inter- polating temporal missing patterns FIM should accommodate time series of all scales. Hence, we again employ min-max-nowe normaaliliiniion for the observation values per time series, as described for FIM-6 in Appendix C.I.
Our approach to the temporal missing pattern imputation task involves one more step of instance normalization, which we want to address now. For this sake, let us briefly recall the notation in- troduced in Section 3.2.3. (ordered) time series with observations {(), T.)) split into K sequentially ordered sets
Vie V ywr+1. :
V1 wl
where the gth set empty, meaning it contains no observations, but defines a consecutive time range, which we call the time gap. Our goal is to impute values in this gap.
Recall that one part of our model (described in Section 3.2.3 and Appendix D.2) applies a pretrained FIM-f to all sets, except the 9th. Internalyy, FIM-E normalizes each set locally, processes it and renormalizes its outputs accordingly. the crucial part of our approach, the model predicts an embedding u, for the qth set, which is processe with F FIM E, to return f. log Var(f) and Po inside the gap.
It is only possible to locally renormalize these outputs wrt. the time range of the gap. because there are no observation values. The global renormalization is not effected by this deficit, because on this scale, the predicted output at the gap is renormalized with the global normalization parameters, equivalent to the outputs at all other sets.
Because the model is subjected to this deficit during training, and inputs are instance normalized globally, the model learns to adjust these conditions.
D.2 MODEL ARCHITECTURE
Let us provide more details about the architecture of FIM, the imputation model for temporal miss- ing patterns, which was already outlined Section 3.2.3 and hyperparameters the single FIM model all experimental results were derived with. Following the model architecture of FIM-( (see Appendix C.2), we use SeLU (Klambauer et al., 2017) as activation for all feed-forward neural networks and employ dropout rate of 0.1.
FTM contains a pretrained FTM-E model, which parameters 0 are frozen. Additionally, FTM has roughly 6M trainable pararameters, which denote by ip.
Model inputs:
A noisy time series (y1, (1).T.)...... R)) with observation values y/ E R and ordered, but potentially irregular, observation times R4 with 7ı <... T7
(ii) An ordered splitting of the interval [m, n]] into K = 5 sequentially ordered sets, which induces splitting of the observations into K sets
F1,-.- Don U yen +1... Bent twg ... yt- uy +1..... :
24

---

Published as a conference paper at ICLR 2025
where the ęth set only defines consecutive time range (or imputation gap) and does not contain any actual observations.
(iii) A time i E R2 in the imputation gap at which to evaluate the functions f and log Var( at.
Model architecture:
Local scale statistics. Defining wo 0 and w
Σk= j-1
Wk for convenience, set
min
Fur +i
V wy
yj a
i= 1...,013
Vj lrst Wuger4
V) Vu +w0
U/
y
Brst
nist
rj
lust
T Fuggee +00g
T
and denote by
$j 1t J 1/
range
2//
160 R2
the statistics about the position and local scale for the jth set, where j * q. (ii) Local scale embedding. A learnable linear layer og embedds the local scale statistics for all
j # q. Let us denote
= =f E 12512
As mentioned in Appendix D.I, we also have access to the boundary times of the imputation gap: 7a . and -iil We can thus locally time instance normalize as I. <- (L - fisst ) //Tditt
such that on the correct scale to apply the Trunk net of the pretrained F FTM- Let us again
denote the output of the Trunk net t" (0f (t).
(iii) FTM-L embedding. We apply the temporal encoder oi and the Branch net equivalent network og v of the underlying, pretrained IM-! model all sets j * q. Concretely, ob yields the individual observation embeddings
Concat Mii (T)) 513 for all i 1, I....,
and ob o vf yields
un = 22 iy
y
twy E R122 for all q.
(iv) Gap embedding estimation. A sequence processing network wig transformer with 4 layers, heads and (self) attention dimension 512, processes the sequence
(uf sf).... ... (u.- - (u sit s4+1)..................... (uk sk)
and returns an embedding for the gth set, that is
up 2/ (u ... (
(ug-1 stt2). (kk,sk)) E RR512
(v) FIM-l time derivative projection. We apply feed-forward neural networks Φg. of and of of the underlying, pretrained FIM-( model to the estimated embedding up and the time in the imputation gap. We write
ft) (h) log (ar(f.)(t) (h)(()), with ht(t) da dt uff).
Finally, we revert the local time instance normalization from (ii) inside the imputation gap, with the formulas presented in Appendix C.
Model outputs: In conclusion, model returns time derivative fT)), with uncertainty estimate log Var f. (t), for a time E R4 in the evaluation gap.
Let us note passing that the intermediate embeddings u², for j + q can be processed as usual with the underlying F FIM-E. Therefore, the temporal missing pattern imputation model can also return (t), log Var( (f) and initial conditions for sets # q, via applications of the underlying FIM-L. In Appendix D.6 we discuss how to combine these outputs from all sets to infer continuous interpolating path (t) defined on 1).
25

---

Published as a conference paper at ICLR 2025
D.3 TRAINING OBJECTIVE
The parameters & of FIM are trained in a supervised fashion, utilising the synthetic dataset for temporal missing patterns introduced in Appendix B.
We briefty recall the notation of Appendix B.4 for the data generated by the synthetic dataset. Let f be the sampled function from one instance in the synthetic dataset, {f((4)) its values on the fine grid, To the initial value of the ODE solution, recorded on the fine grid {a(t+)) and {(tt, T the associated noisy time series.
The noisy time series is split into K sequentially ordered sets, outlined in Appendix D.2. Follow- ing the notation therein, we can identify the subset
{ts first s (ti)
of points on the fine grid contained in the imputation window.
Using {(V,, n=) as inputs for FIM, we denote its outputs evaluated on fine grid points contained in the imputation window by f(t)) and log Var( (tt) as in Appendix D.2.
The training objective of FIM is then maximize the Gaussian log- likelihood of the time derivative log (() f(t)), (()())(()) at all M points E of the fine grid contained in the imputation gap. Expressed as equation, train FIM to minimize
M
(1(() - ((9))) log Var ((t))
C =
E f
(V)r(f)(t))
2
D.4 TRAINING PROCEDURE
The parameters y of FIM were optimized with AdamW (Loshchilov & Hutter, 2017), using a weight decay of 1e-3 We where the learning rate decays from le to le over 400 epochs. Note that we do not deploy warm restarts. With batch size of 1024, we trained FIM single A A100 80GP GPU.
cosine annealing schedule as introduced by Loshchilov & Hutter (2016),
D.S ABLATION STUDIES: DATASET DESIGN
The dataset is generated using GP model with periodic kernels. The grid spans from 0 to 1 with resolution 512 points. Each GP created by combining kernels, and the parameters are sampled from uniform distributions. The lengthscale É is sampled from (4(0.55, 1), while the period is sam- pled from (0.3,, 0.5). The training dataset contains 500,000 samples and the test dataset contains 50,000 samples.
Our experimental datasets for the temporal missing pattern imputation (see Appendix H) contain significant trends. Moreover, the temporal size of the imputation gap is 20% for all samples in the time series Yet our synthetic generated training data includes both: trends and variable window sizes.
To assess the impact the trend in our synthetic generated training data on the evaluation perfor- mance, train two sets of models: trend. For both datasets, we train two models, varying the size of the imputation window. Specifi- cally, the window size is sampled from either (10, 30) or (d,, 50). In Table we present the results of the four different FIM models trained on the four different dataset. One can see that the model trained on the dataset with smaller imputation window performs the best. Increasing the imputation window size leads to worse performance on the Motion Capture dataset for both case (PCA and no PCA). The performance of the FTM model also decreases when trend is included in the training dataset.
with trend in the training data and the other without
D.6 CONTINUITY ALONG THE GAP
Similar to the windowing scheme of F IM-f f (see Appendix C.6.2), FIM processes a time series split into K sets of sequentially ordered observations. Continuing the notation from Appendix D.2,
26

---

Published as a conference paper at ICLR 2025
Table 8: Ablation study on the size of the imputation window size and the occurrence of trends in the synthetic training data of FIM. We report RMSE in the imputation gap of the Motion Capture dataset.
Motion Capture
Dataset
PCA
No PCA
No Trend (10,0, 30)
3.27 + 1.12
2.97 1 0.96
Trend (10,, 30)
3.41 +
1 1.07
No Trend (5, 50)
3.59 1 1.14
3.52
1 1.18
Trend - (15 50)
3.67 = 1.13
3.59
+ 1.19
now detail our approach to connecting these outputs, such that the combined interpolating path ((()) continuous [0, 1].
We can group sets containing observations, i.e. j outside of the imputation gap, into two disjoint groups: sets S < q to the lefi of the imputation window and sets 4 < j to right of the imputation window. Model outputs for each group be combined with the ideas from FIM-(, presented in Appendix C.6.2. Indeed, these sets are just windows that can be combined on some predefined overlap.
Let f and jr denote the inferred trajectories to the left and right of the imputation gap [r s/m), which are designed to be continuous by the efforts described in Appendix C.6.2. To connect the two trajectories to continuous path inside the imputation gap, we first extend both of them individually by integrating the time derivative fy inferred from F IM with initial values glen and fright respectively.
Expressed in equations, define for & fin [lls]]
the following extensions:
a d (t) t (T [ frrtt
fs(dds
f." (s)dss
d (t) sright
(T
Jost
These extensions are then combined via weighted temporal interpolation, where define the com- bined solution at time t E [0, by
(t) Jnst
Best for < < fisst
a(1)
plon
+
{ 1---he lust
rright (t) for E [T [Io
(20)
lust
grighe
for > lnst
Note that & is continuous by design and its corresponding combined time derivative f can be found (a.e.) by taking the derivative of equation 20.
D.7 PROCESSING TIME SERIES OF ANY DIMENSIONALITY
FIM contains a pretrained FIM- model. In fact, FIM just extends that (frozen) model by two networks: a linear layer Og to embed local scale statistics, and sequence processing network 1/2. processing a sequence of FIM-f embeddings of sets (consecutive) observations. Therefore, also only returns ID time derivatives. To impute arbitrary dimensional time series data with FIM, we process each component of a time series separately. The same channel independent strategy is used for FIM-6. In Appendix C.6.1 we discuss limitations of this approach, which are also valid for the temporal missing pattern imputation model FIM.
E EVALUATION METRICS
Depending on the experiment and the baselines, we use one or more of the following evaluation metrics for multi-dimensional time series.
27

---

Published as a conference paper at ICLR 2025
Following the notation of Du al. (2024), consider multi-dimensional target time series I1...., LL, with 2, E RD. and corresponding estimated time series E1,..., -L. Some baselines only compute the metries below certain elements and components of each time series. Let m......, m be an associated sequence of masks my E {0, 1) indicating which elements and components are part of the metric calculation. Finally, denote by aq a4 and my the dth component of each vector.
We compute per trajectory the Mean Absolute Error (MAE), (Root) Mean Squared Error (R)MSE))
Mean Relative Error (MRE) and the R22 score, defined as follows:
MAE(2, m)) Di=l= d=1
D
1
2
m
Σ 2d=1 m
D
MSE(2, 2t, m)) E d=1 (1
) m
2
Σ/= Ed=1 m
RMSE 2, m) VMSE (z, ft, m)
E Ed=1
1
m
MRE(Z, 2)
Ei= Ex=1
[9 m
1 (x D
E/e
(x
1
El (2) xd
where nd = 1 EL 2. Note that we only require the R² score for the comparison on ODEBench, where no targets are masked, so we do not add the mask in the formula above.
For the imputation task of point-wise missing patterns, usually compute the metrics only at the missing points, following the methodology of our baselines, e.g. (Du et al., 2024). One exception is the ODEBench evaluation, where the baseline values computed against the complete target trajectory, including observed and missing time points. We adapt our computation to the baseline methodology for valid comparison.
For the imputation task of temporal missing patterns, metrics are only computed in the imputation gap, in accordance with our baselines.
Usually, we report mean and standard deviation of the metrics calculated over all trajectories in a dataset. The only exception is the -ccuuurayy defined by d' Ascoli et al. (2024) the percentage of time series with 12 score larger than 0.9. Note that Du et al. (2024) and Fang et al. (2024) only report metrics averaged over all missing values over all time series, not standard deviation time series.
F PHASE PORTRAIT RECONSTRUCTION: ADDITIONAL RESULTS
F.1
VAN DER POLI OSCILLATOR
F.1.1 DATA DESCRIPTION AND PRE- PRECESSCESSIG
One particular (second-order) differential equation from the ODEBench dataset is the Van der Pol oscillator:
dx =
du = p(1 - (1-22)vu- x
We simulate the system in the time interval [0, 10] for set of initial positions z(0), initial velocities (0) and parameters μ::
-3, x00) -3,,-3,--11 1,,1,
(0) E E -455,4455)
# (0.1, 1.5)
28

---

Published as a conference paper at ICLR 2025
Time
Time
1ime
Figure 4: Phase portraits of Van der Pol systems with parameters fr 0.1, 0.5, 1.5. 128 noisy observation (black circles) get interpolated by FIM-f(w.n. 4) (blue line). The inferred time derivative matches (magenta line) matches the ground truth vector field the solution path (black dashed lines).
Thme
Time
10me
Figure 5: Phase portraits of Van der Pol systems with parameters 0.1, 0.5, 1.5. 32 noisy observation (black circles) get interpolated byF IM-((w... = 4) (blue line). The inferred time derivative matches (magenta line) matches the ground truth vector field the solution path (black dashed lines).
We generate 1 / time series by observing the position s=(t) only at a regular grid of either 128 or 32 points the interval [0, 10] and corrupt them with multiplicative noise of level = 0.05 inspired by the ODEBench setup (see Appendix F.3).
F.1.2 MODELLING AND RESULTS
We apply FIM- ((w.n. 4) and visualise the inferred path and the inferred time derivative by means of phase portraits. Figure contains the visualisation for 128 observations and Figure 5 for 32 observations. In each figure, the parameter f ranges from μ 0.1 in the left, to = 0.5 in the middle and μ = 1.5 in the right subplot.
FTM-6 interpolates the observed path and the inferred time derivative f(() matches the ground truth vector field v(t) at the path closely. The accuracy of the (implied) estimation f(0) of the initial velocity v(0) is not perfect and degrades with increasing μ.. In other words. based purely the (limited) observations, the model can not infer the sharp change in velocity at the beginning of trajectory. The bottom, 1D subplots reveal that FIM-6 adjusts to quickly.
accurate estimation
Even with only 32 observations (considerably less than the roughy 256 observations in the ODEBench), FIM-E recovers the ground truth time derivative quite well, showcasing that FIM-E does not require many observations approximate simpler regions of the dynamics.
29

---

Published as a conference paper at ICLR 2025
F.2 RÖSSLER ATTRACTOR
F.2.1 DATA DESCRIPTION AND PRE-PROCESSING
In Figure 2 we consider the Rössler attractor differential equation from the ODEBench in its chaotic variation:
dr -5(y + 2) dy = 5(0.2y + x) d= 5(0.2 + 2 5.7 2)))
Also following the ODEBench, simulate the system in the time interval [0, 10] with initial value (2.3, 1.1, 0.8). From 4096 observations on a regular time grid in [0, 10), we subsample irregular observation grid by dropping each observation (independently) with probability 0.5. The subsam- pled data is then corrupted with multipiccative noise of level 0.05. This corruption scheme follows one corruption scheme of ODEBench, which described see Appendix F.3.
F.2.2 MODELLING AND RESULTS
We apply FIM-E(w.n.. 32) recover the trajectory from the corrupted data. Our model also returns the time derivative along the recovered trajectory, which is the velocity of the particle in the Rössler attractor system. The centered plot in Figure 2 shows the ground-truth and model inference by means of phase portrait. It contains position and velocity for the first component, i.e. x and dr, along the trajectory.
Our model can also approximate the acceleration of the particle. First, we discretize the inferred velocity via a regular grid of size 8192 in the time interval [0. 10). To this discretized function, we apply FIM-I(w.n. = 64) and recover an approximation of the acceleration of the particle. The rightmost plot in Figure shows the velocity and acceleration the first component along the observed trajectory, again by means of phase portrait.
F.3
ODE BENCH
F.3.
DATA DESCRIPTION AND PRE-PROCESSING
The ODEBench dataset (d' Ascoli et al., 2024) contains the solution of 63 ODEs with 2 initial con- ditions each. The solutions are available on regular grid length 512 the time interval [0, 10).
We follow the pre-processing of d' Ascoli et al. (2024), who corrupt both the observations values and the observation grid. Let us recall their corruption scheme here for completeness.
The obvervation values corrupted by multiplicative noise, yielding noisy observations of the ground truth solution x; = xt,,) at some time via yi = (1 + e)2;, where e ~ (), The observation grid is corrupted by dropping each observation (independently) with probability
All combinations of y E {0, 0.05} and p E {0, 0.5) are considered.
F.3.2 BASELINES
We interested in two tasks on this dataset:
reconstruction of the underlying ODE solution (which is the task considered by d' A'Acolli et al.
(2024)) and
(ii)
inference of the vector field on the solution path, i.e. the time derivative.
Note that in the of p 0.5, this setup is very similar to the missing point imputation task considered other works.
While d' Ascoli et al. (2024) consider the first task (and we compare against their results below), they do not report any results related to the second task. We therefore reevaluate ODEFormer on the ODEBench, gain access the time derivative information needed for comparison. We use
30

---

Published as a conference paper at ICLR 2025
the implementation and weights provided by the authors", apply it to 10 samplings of the corruption schemes (as we do for our own model) and average the results.
On the solution reconstruction task, re-evaluation yields different results than reported by the authors (see e.g. Table 9). We report both results, if they are available, and name the results from re-evaluation ODEFormer Re-ev..
Results of other well performing models on the ODEBench have been extracted from Figure 4 of d' d'Ascoli al. (2024).
F.3.3 MODELLING AND RESULTS
We apply F FIM- to all equations of ODEBench and compare our results to the best performing models from d' d'Ascoli et al. (2024), as well as our re-evaluation of ODEFormer outlined in Ap- pendix F.3.2.
The authors of ODEBench compare several methods based on the reconstruction of the ground- truth ODE solution all 512 available points of each trajectory. This comparison is based on the R² accuracy the percentage of predictions of which the R score exceeds 0.9. We compare our zero-shot model FIM- with different numbers of windows to the best performing models from d' A'Accoli et al. (2024). We report our result wrt. the Re-accuracy (see Table 9), the RMSE (see Table 10) and the MAE (see Table 11. For our models, and if available for ODEFormer, we provide the metrics for equations grouped by their dimensionality, and average our results over 10 samplings of the corruption schemes detailed above.
Averaged over all equations, FIM- improves the baselines on all corruption schemes, in particular with larger number of windows. The optimal number of windows depends on the difficulty of the underlying solution. MMore complex dynamics require higher local resulution, thus more windows, than easier dynamics, which benefit from the smooting effect of low number of windows. Gener- ally speaking, equations of the same dimensionality in the ODEBench are of similar complexity. We therefore denote by FIM-A( (Weighted Sum) the average the best performing number of windows per equation, weighted by the number of equations of that dimension.
Noticeably, the performance drops for equations dimension Figure shows the path recon- struction for one initial value problem of each equation. Trajectories of 3D systems included in ODEBench are inherently more complex than for all other dimensions. Still, FIM - É(w.n. = 8), used for Figure 6, infers suitable initial conditions and vector fields that approximate the ground- truth solution well.
Let now consider the accuracy of the inferred time derivative along the solution path. We provide the results for our model and for our re-evaluation of ODEFormer on ODEBench. Similar to the original ODEBench task, compute the metrics against the time derivative along the solution path of the ground-truth solution. To extract an estimation of the this time derivative from ODEFormer. we can evaluate the inferred equation along two possible paths:
the ground-truth solution of the ground truth equation
(ii) the inferred trajectory from ODEFormer
While we report both evaluation approaches, we find that the second performs consistently better than the first. We mark the results with ODEFormer (g.1. traj) and ODEFormer (inf. traj.) respec- tively.
We report the results wrt. the RMSE (see Table 12) and the MAE (see Table 13). Our model performs better than ODEFormer on both evaluation approaches. The inferred time derivative is much closer the ground-truth vector field along the solution path, especially when we use more than 2 windows.
Pttpps: thub, wttp:/coo/cccco///oeeco/ former
31

---

Published as a conference paper at ICLR 2025
Table 9: ODEBench Rº accuracy of FTM-& l for different numbers of window, split by dimensions. The standard deviation is calculated across 10 samplings of the corruption schemes.
0.5
0.5
Dim
Model
Y
Y
ODERsmmer
69.9
60.3
ODERocmer-gg
ODERomer Re-ev.
64.3
64.7 =
62.9 =
61.6 6 3.0
PySR
82.5
63.2
77.0
38.2
AlL
69.0
61.0
71.5
76.1 + 1.3
f(wr.n.
86.6
FI3- E(w... 8)
91.3
88.3 0.5
91.6
87.5 1.0
FMM F[mn... 16)
160.0 - 0,0
96.2 0.5
97.2 1.2
91.6 2 1.9
(Woightdd
100,0
96.9
1.2
ODEFvmmer Re
89.1 0.0
87.8 2.0
87.5 2.5
84.3 3.8
1.7
FIM-E(w.n. =
100.0 0.0
95.4
99.8 0.7
95.2 0.9
FIM-E[w.n. =
100.0 - 0.0
96.4 0.7
99.8 =: 0.7
93.5 1:5-::1 1.0
2IM-E(w.n.
100,0 0.0
93.7
100.0
+: 1.4
ODERamer Re-ev.
62.5 + 0.0
63.8
3.5
61.4 = 3.5
60.7
E(w.n.
625:3:00
61.4
80.0
11.3-22
946
93.8
95.0
9364
FIM f(.....
100.0 0.0
98.6
0.8
99.3 =
97.5 12.5--00 0.9
FIM ((.. 1G)
160.0 + 0.0
98.4
0.6
99.5 +: 0.9
97.0 1.5
ODE Re
10.0 4 0.0
13.1 = 3.7
9.4 8: 3.2
12.5-4:477
FIM-E(
15.0 0.0
15.0
17.5 1: 2.5
180:11
FIM-F(w.n.
3000::00
1:::::2 2.4
FIH-E(w.n. F12-E(w.n. 16)
45.0
40.5
1.6
49.5 6.0
43.0 6.3
100.0 -t: 0.0
95.0 0.0
84.0 +: 6.1
74.5 1:: 6.9
ODEicmmer Re
75.0 +
08.8
11.6
68.8 + 11.6
50.3 + 17:7
FI9 Flm....
T00.0 = 0.3)
TinLI
11.00
11IL 100
100. 2: 0J)
f(w...
160,0
100.0
0.0
100.0
100.0
Fim--iw...
100.0
100.0
FI9 f(n.n.. 16)
100.0 - 0.0
100.01 0.0
100.0 0.0
100.0 0.0
Table 10: ODEBench RMSE of FIM- for different numbers windows, split by dimensions. The standard deviation is calculated across 10 samplings of the corruption schemes.
p
0.5
Dim.
Medal
7
7 = 0.05
7 0.05
149561
FIM E(mron. FI3 f(w.n.
149276 2
L.3191 = 0.01377
1.1169
1.020%
117046 012288
AlL
F(an.n.
0/47144
0389672 + 0.0
057397 0.00558
0.747 + 0,01644
0.451
086858 + 0.01
0.1079
f(mn.n. 16)
0.14572
0.0059
0,24992
FI44 (Weighted Sun)
0.14572 + 0.0
0.25858 0.0059
0,28492 0.000
0,41192 + 0.0205
ODLivcmer Re ev.
0.50954 0.0
0.52900
0,13294
0,58829 0.04373
0.66961 4 0.09289
F14-[(w.n.
0.47952 0.0
0.54023 1 0.01949
O 012225 1 0.00566
0.26812 1 6.01333
FIM-E[w.n.
0.09113 - 0.0
0.22905 0.01549
0.04558
0.00269
0,27359 1: 0.06622
FI34-((w... FI4-f(w.n. 16)
0.01472
0.01135- 0.0
0.36661 0.01391
0.02188 0.00081
0.5:251
0.01067
Revev.
0.0
0.55307
0
1194-2um-n =
0.4882 = DAD0S
1L15043 01010
0.420 0.011
FIM f(w.n.
0.14509 = 0.0
0.15664 0.00234
13077 = 0.01999
0.5587 + 0.01388
0.0
0L060O8 0.00199
0.48811
0,9019 4 000698
FTM- f(m.n. 16)
0.0857
0.003327
0
0.01066
+ 0
6.34588
6,65072
0.38931
6344523 1
700731
0,07836
79348 0.14398
3
(..-- FIM-f(W.N.
5.15268 :: 0.0
5.16363 4 0.02147
4.11386 à: 0.0952
4.4104
0.9907
2,8077 1:
2,79512 0.01428
2.5567
0.05278
2,57153 1: 0.02113
f(w.n. 16)
0.0U5S
0.13815
ODEKemer Re-ev.
B/3579
0.31700 0.00393
0.3305
0,00303
0,0530 + 0.97477
f[w.n.
FIM-Fiu... =40
0/0213 - 0.0
0,09696 2: 0.00964
0.00447 2 0.0009
0,00786 1 0.08079
FT9 f(m.n. =8)
000079 : 0.0
0.00715
0.00041
0.0144 - 0.00065
0.0961 + 0.0046
f(w.n.
3,0046
0.0011 0.00083
32

---

Published as a conference paper at ICLR 2025
Table 11: ODFBench MAH of FTM-f f for different numbers of windows, split by dimensions. The standard deviation is calculated across 10 samplings of the corruption schemes.
Dia.
Model
p 0
p=
p p= 0.5
P = 0.5
0
0.05
0.06
ODEiommer Re-ev.
1,0675 = 0.0
1X8829 # 0.15297
1.17828 # 0,03353
L,1533 + 0.04755
FI8-E( = 2)
1,13992 = 0.0
116/99 + 0.0090
0,83511 + 0.01963
0,85S8 + 0.01983
FIM-E(w
0.65753 ::
0.71483 0.00482
0,57207
0.01368
0.64132 0.01129
FIM- ffw.
0,33755
0,90221
0.41244
FIH- 3IN-E(w.N. n...... 16)
0.06228 0.0
021283 0.00345
0.17416
0.01029
0,34229
0.01006
ODERamer Revev.
0.3223
0.34722 +
0.99660
047208 0.4834
FIM-Z(w-n, n. = 2)
0.3878 = 0.0
0A3727 1 0.01N0
010355 Σ 0,0006
0,3000 I 00105
77M--[mm...
0.7499 0.0
0.16663 + 0.01033
0,1213 0.00204
0.20268
0.01535
f(w-n.
0,040
0.20011
0.26112 0.01003
FIM-E[w.n. 16;
0.08011 0.0
0,26645 1 0.00826
0.0692 - 0.00025
0.38727
0.00971
ODEFomer
0,40300
0,07357
F1H-E[w.n.
0.3346 0.0
0.4077 0.00229
022652 1 0.01166
0,23367
000281
2 FIH-f(w
0,08653 = 0.0
0.99166 + 0.00135
0 0,00575
0,09457 + 0L006UT
FIM-E(w.N) FIM f(w.n. 16)
0,02190
0,00982 0.00094
0.02557
0.00285
0.05379 0.00299
0,0485
0.03821 0.00141
0.01448
0.06105 0.00215
ODEFamer
4.13377
4404779
4,88917 0.11551
411050
4,42416
3
FIH-f( n. = 4)
3.84243 : 0.0
3,84711 : 0.01762
3,31 1 0.08397
3.30732 1: 0.08331
FIM-F(w f.--f(w...
2,01825 = 0.0
2,01172 + 0.01012
183143 0.03501
1,84830 0.04100
f(w.n.. 16)
0.57452
1,09200 0.06818
ODEKomer FiH-t(w.n.
0,0299
010361 0JIE2SK
0,02508 0.00251
0,02356 OLORD4S
0.0255
0,00900 000129
4 FIN-E(w.n.
0.00132 :: 0.0
0.00331 0.00054
0.00131
0.00037
0,005S8
0.00061
FIM f(w.n.
0.00066 = 0.0
0,00111
0.00X00
0.00IK6
+ 0.00020
0,00685
+ 0.00E3T
=IH-E(w.n. 166)
0,0027
0.00500
0.80017
0.00053
0.00117
0,00833
0.00441
Table 12: ODEBench RMSE of the inferred time derivative of F FI-E for different number of win- dow, split by dimensions. The standard deviation is calculated across 10 samplings of the corruption schemes.
Model
p
P 0
0 = 0.5
0.05
0.UG
ODEForeer (g.1. tmai.
30,89650 0.0
71,30478 +: 131,29822
69.30647 :: 50.34435
93.26899
105 60553
ODEFumer (iní. (ra..)
16,50882 0.0
13.10672 +: 2.30755
19,82354 1: 2,94507
17.77896 +: 4.61166
FIN-7[w.n.
10,09933 1.00
10,07:11 1: 0.01023
980062 1: 0,008
933593 Σ 0.00R33
FIN-C(w.o.
9,01272
FIN-C(w... 8)
684600 1:
7.177 0.03336
691581 1: 0.03668
7A4845 1 0.04305
2I..-( 16)
3.49679 + 0.0
4,79088 + 0.06319
4.46558 + 0.13783
6.48983 + 0.17314
PIN-E (Weighted Sum)
3.49679 1:
3,78833 1 0,0204
4,46358 0,1321
4,78349
0.136708
ODEPormer (p.t. tnj))
8,58175 + 0.0
6,98587 + 29939
23 2306 + 18.09459
20A331 + 28J9377
ODfFormer (inf Imj
2.4357
12.3632
PIR Z(w.n. 2]
0,35337 0102785
INSALE =E
PIN f(w.n. 4)
0.08764 + 0.0
0.34632 + 0.03219
008866 d: 0100263
0.43992 + 0.04078
PIN-f(w...
0.06123 1 0.0
0,87099 1 0.08969
0.66833 2 0.00991
1,34817 = 10076
fow... 16)
0,04306
2.85851 0.14102
0.5609 0.00669
4.54185 0.25996
ODEFomer (g.. tnj)
9,48347 1: 0.0
115,23801 1: 298.68335
9.06334 E: 3445004
100.5754 E: 239 8902
ODEtormer m()))
199699
10577 2.00177
5.8334 3.27421
2
FIN-((w.n = 2]
FIN-f(w.n. = 4)
0,70673 1: 0.0
L43282 1:
0,73517 1 0.00655
1,41955 1 0000186
0.69129 1: 0.03759
121969 1 0.01069
1,22766 1 0,11812
0.7396 1 0.03952
PIN-f(w.n.
0,21854 4
0,32851 0.01379
0,38857
0,65065 + 0.05582
FIN-C(w.o. 16)
0.51666
0.18537
0,8219
ODEFFomer (g.t. tmj.
145, 845883444 +
110457446 + 67.51671
35% 12564 284-26566
51259 199.26689
ODFFormer (inf tra))
86,05997 0.00
84346 + 4-6072
96053
58068
72.54148 + 74.72853
3 PIR Mu..n. = 2]
5827993 2.00
58.28155 4 00003S
S7.7226
2: 010192
ST.TITY = 01U1773
PIN-f(w.n.
54,04507 1.00
53,92284 1 0.01339
94947 0.05699
52 52.58163 1 0.0522
PIN-ffw.n.
42.37861 1.00
50.597566 0.22360
PTN f(w.n. 16)
21.67166 0.0
22.147%
0.09837
48430 0.83561
28.11654 0,89331
ODEFormer (g.t tnji
0,12387
0.16521 : 0,04465
0.43155 + D.4TST
ODEifrmer (inf. tral.)
0,0639 0.0
0,08313 0.01145
0.06529 + 0.01458
0.92522 + 2.12509
4
PIN PIN-f(w..
FIN-t(w.n. = 21
0.03748 1 0.0
0.06423 4
0.03646 1 0.00119
0,0116
0.00122
0.0766 1: 0.00126
0.00546
0,00844
0.01566
0,02066 1. 0.00155
0.00131
FIN-f(w.n.
0.0245
0.00529
0.00078
0.00243
FIN-C[w.n. = 18)
0.00283 1 0.0
0.0707
1 0.00888
0.00366
1: 0.00053
0.10999
0.00559
33

---

Published as a conference paper at ICLR 2025
Nusy dbservation
FIM
Ground-irah
80
-15
-25
0.25
www
050
10
Figure 6: Reconstruction from noisy observations (dots) of all equations in the ODEBench dataset, using the first set of initial conditions, with corruptions 0.5 and y 0.05. FTM-f (line) recovers the ground truth path (dashed line) closely for (almost) equations.
34

---

Published as a conference paper at ICLR 2025
Table 13: ODEBench MAE of the inferred time derivative of FIM-& for different numbers of win- dows, split by dimensions. The standard deviation is calculated across 10 samplings of the corrup- tion schemes.
0
U.D
0.5
0.05
0.05
ODEFamer (g Irj])
9.08246 = 0.0
10,3509 =: 5.12400
11:43992, 1: 3.80992
14,56679 1: 5162873
ODE (inf. mra)))
7.47623
7.59557 :: 0.80806
8,03686 0.38837
7,87898 0.61724
FIH-E[w.m
6,71538
6.1500
FIM-f[w.n.
5.01699
5,64520
5,74910 0.01754
FIH-E(w.N. n. 8)
4.21719 0.0
4.45152 1: 0.01683
4,21716 1: 0.01543
4,58770 1: 0.02192
FIH-f(w.n. 16)
1,9056
2.75833
0.03090
2.44375
0,05500
3,79191 0.A4225
ODEFomer (g.. tra.)
1002333 = 0.0
0,85939 + 0.28399
193692 0.79555
L69K55 1.27994
ODERamer (inf. 3ma
1.22521
1.79378
1,74344
118--[[....
0DG9ST
FIM-F(m.n.
0.3440
0.0
0.21994 = 0.01723
0,02922 + 0.00100
0,29445 + 0.07833
FIM f(w...
0.08447
0.0
(L57124 0.04161
0.02071 0,00083
085952
0.04000
FIM f(w... 16)
000120
1,3990
0.12402
ODEiomeer taj.)
1.14244
4,60410 :: 9.66250
1,29322
0.11240
4,67499 7.50841
mm)))
0.55929
2
F1M-2[w.n. =
0.78924
0178135 2: 0.00215
0.39850
0.00569
0.M0813
000050
FIH-E(w.n.
0.33064 =: 0.0
0.15679 :: 0.00296
0,29782 + 0.00601
0,33816 4: 0.00718
FIH-E(w.n.
0.08201 : 0.0
0,15581 1: 0.00496
0.10640
0.00910
0.21836 0.01454
FIM-E(w.n. 16;
0.00003
ODEFamer (gt. tnja
$7,4099
50,32888 + 20. 20.88564
64,04977 244457
74.75737
ODEFomer mma
42:71222
0.0
43 86978
44.45747 2,04799
43,39972 +
3
F1n-E[w.n. = 2)
39125040 3: 0.0
39 25000
U.INSOL
3M92776 E ULNIPSS
F1M f(m...
35.2175
= 0.0
35 50043 0.014422
3466156 + 0.09962
34,59348 + 0.8666
f(w.n.
26,29693
0.0
26,29145
0.10493
FIM-f(w.n. = 16)
11 0.0
12.40848 0.03408
15 15,20400 0.34251
15,83900 1 0.32635
(gt. trj.)
ODE (inf. mi))
0.01773
0.01.333 1: 0.00195
0.01558
0.00135
0,06136
0.0929]
4. FIN-t(w.n.
0.01224 = 0.0
0.01785 + 0.00347
0,0E562 + 0.0006
0,00K22
0.00007
FIH-E(w.n.
0,00178
0,00667 0,00026
0,00179 0.0001
0,00892
0.00061
FIM-E(w.n
0.00135
0.00134
FIM-E(w.n) 16)
0.00117 0.0
0.04355 0.00189
0.00143
0,00005
0,06961
0.00279
F.4 LORENZ SYSTEM AND LATENTÔDE
F.4.1 DATA DESCRIPTION AND PRE-PROCESSING
To compare our zero-shot model to a LatentODE, a specialised black-box model requiring training, we consider the chaotic Lorenz system
dx ay x) dy = pz F2
dz ry Az
with parameters (a, 3) (10, 28, 28. We simulate the system in the time interval [0, 10), starting with initial values sampled from N ([2.3, 8.1, 12.4], 1). Note that the Lorenz system with these parameters and the mean initial value part of the ODEBench dataset.
To generate time series, we sample an initial value, simulate the system and observe it at the regular grid on (o, 1] with 512 points. The observations are subsampled with survival probability of 0.5 and corrupted with additive gaussian noise sampled from (0, 0.05).
We generate 1024 time series for the training of LatentODE and additional 128 time series each for validation and test.
F.4.2 TRAINING LATENTODE
We train LatentODE (Rubanova et al., 2019) on the standardised train set, selecting the model based on the lowest loss the validation set. The encoder is LSTM, the emission model diagonal gaussian, fixing its standard deviation to 0.01. The target objective the likelihood the observa- tions.
35

---

Published as a conference paper at ICLR 2025
Cican
Observations
Ground Truth Vectorficld
Vectorfiald
Xmmmy
www
-15
-20
-20
200
Figure 7: Inference of FIM f(w.n. 16) on a Lorenz system time series, split int their individual dimensions. (Top) Noisy observations (black dots) of the system (black line) are interpolated by FIM-f (blue line). (Bottom) Inferred values of the time derivative (blue dashed line) match the ground truth vector field values along the solution path (black dashed line) closely.
We train for 6000 epochs, with minibatches of size 32, using AdamW with learning rate 1e-3 and weight decay le 1e-2. To help the model learn, we slowly anneal the input time series length over the initial 2000 epochs, starting at 25 observations.
We selected the model architecture by grid search. See Table 14 for the search grid and final param- eters. The model trained roughly 9 hours A100 40GB GPU.
Table 14: Hyperparameters, including grid search range over some, for the LatentODE baseline on all datasets. Selected hyperparameters, determined by the loss on the validation set, are underlined. For simplicity, we set the dimension of the forward neural ODE to the hidden size of the encoder.
Hyperparameter
Lorenz System
Motion Capture
Navier Stokes
Hidden size
64, 128, 256 [64, 64], [128, 128]
32, 64, 128
64, 128, 256
NeuralODE hidden layers
[64, 64
[128, 128]
Initial condition hidden layers
64, 64 [128, 128
64, 64
128, 128
Emission model hidden layers
64, 64 [128, 128]
64, 64
[28, 128]
F.4.3 MODELLING AND RESULTS
Note that Latent ODE is notoriously difficult to train and we only managed to fit the data on 64-dimensional hidden space (Dupont et al., 2019). The time derivative inference is therefore not directly feasible with our trained Latent ODE, We nevertheless compare FIM-É(w.n. = 16) against said Latent ODE model the reconstruction task, well against the Scipy implementation of the cubic spline.
Table 15 displays results in terms of RMSE calculated against the ground-truth. clean solution path on all 512 of the regular grid. They averaged over the 128 trajectories in the test set for described data corruption scenario. FIM-( not only outperforms both baselines, but also perfectly infers the hidden vector field along the solution path, as can be seen in Figure 7.
Table 15: Reconstruction task the Lorenz system data. RMSE is calculated the whole trajec- tory.
Model
Lorenz system
LatentODE
3.25 0.99
Cubic spline
3.97 5.8
FIM-F
2.01 + 0.33
"https://scipy"
36

---

Published as a conference paper at ICLR 2025
F.4.4 ABLATION STUDY: NUMBER OF WINDOWS
Additionally, we perform an ablation study on the numbers of windows processed by FIM-é on the Lorenz test set.
Figure shows the RMSE as the number of windows increases from 1 to 32.
Initially, the RMSE considerably higher than the LatentODE baseline, demonstrating that FIM-& on it own can not handle arbitrarily complex systems, because it has not been trained on such com- plex data. With the addition of more windows, the RMSE decreases sharply. However, beyond approximately 14 windows, the RMSE reduction tapers off, suggesting diminishing returns with further increases in the number of windows.
12
LD
FIM-t
10
8-
RMSE
4-
2
of Windows
Figure 8: RMSE of LatentODE of FIM-((w.n. m) for numbers of windows m between 1 and 32 on the Lorenz test set. The standard deviation across 128 time series is represented by the shaded area for the LatentODE and by the whiskers for each F FIM-((w.n. m). The performance of FTM-f increases, as the number of windows increase.
F.5 COMPARISON AGAINST NEURAL ODE PROCESSES
F.5.1 DATA DESCRIPTION AND PRE-PROCESSING
Norcliffe al. (2021) introduce four ID synthetic datasets of different functional type, including sines, exponentials, straight lines, and harmonic oscillators. Each task is defined by a parameterized function, where the parameters trajectory is generated by sampling from these parameter distributions and then sampling from the function evenly spaced timestamps, i, over fixed range produce 100 data points (t. v). The equations for these tasks, along with the ranges for their defining parameters, are provided in Ta- ble 16.
sampled from predefined uniform distributions. An example
10 random context points are selected out of all 100 points. They serve as the model inputs, where the task is reconstruct the original function all 100 points.
Table 16: Description functions considered by Norcliffe al. (2021), including the mathematical form, ranges for parameters and b, range for and the number of test samples.
Task
Form
b
test
Sines
= a sin(t - b)
(-1,1)
(-1/2 1-1/2,1/2
(-m, m
10
Exponentials
/60 x exp(t b)
(-1.1
(- 1-/2,11/2
(-1,4)
10
Straight lines
at
-1,1)
1-2,,
(0.5)
10
Oscillators
y = a sin(t b) exp(-t/22
(-1,1
(-1/2,1/2)
(0,5)
10
37

---

Published as a conference paper at ICLR 2025
FIM
Ground Truth
Observations
0.5
1.75
0.
1.0
9.73
0.0
-61
3.00
Figure 9: Samples from the one dimensional regression dataset of Norcliffe et al. (2021). The ground-truth path (black dashed lines) observed at 10 irregular observation times (black dots). FIM-8 8 interpolates (blue line) well and recovers the ground-truth dynamics.
F.5.2 MODELLING AND RESULTS
In Table 17, we present the performance results our pre-trained model FIM IIM-((w.n. 1) on ID Regression Task proposed by Norcliffe et al. (2021). MSE, MAE Rta are calculated over all the points (not just a subset of the target points) for different samples from each function type. The metrics from 10 samples per function type averaged, and we report the mean and standard deviation.
It is important to reiterate that model has not been trained this dataset. Still, our model outperforms the baseline from Norcliffe et al. (2021) three out the four tasks. The Oscillator set is the only one where model underperforms. Its trajectories exhibit rapid changes the edges of the considered interval, which F FIM-& l can not identify based solely on the sparse observations.
Figure shows four samples from the test set and the corresponding inferences of model.
Table 17: Results the 1D) regression dataset from Norcliffe al. (2021). We report the and standard deviation over 10 samples of each function class. The Neural ODE Processes baselines are selected as the top performing model from Table of Norcliffe al. (2021).
Function
Model
R²t
MAEX 10 24
MSEX 10 -21 21
Neural ODE Processes
0.25 + 0.04
Exponentials
F 1]
0.9656
0.0079
1.11
0.11
0.07 0.02
Neural ODE Processes
3.16 0.76
Linear
F FIM
W.n.. (....
1)
0,9992
0.0002
3.85
1 0.47
0.34 0.09
Neural ODE Processes
0.55 1 0.03
Oscillators
F FIM (w...
0.8822
= 0.028
20.43
= 2.54
22.62 = 7.22
Neural ODE Processes
2.09 0.12
Sine
FIM- w.n.
1)
0.9679 0.0042
4.84 + 0.53
0.67 + 0.16
G POINT-WISE MISSING PATTERN IMPUTATION: ADDITIONAL RESULTS
G.1 COMPARING AGAINST BAYOTIDE
G.1.1 DATA DESCRIPTION AND PRE- PRO-ESSSNSINNG
For the point-wise missing pattern imputation model, first consider two available datasets from Fang et al. (2024), namely the Praffi-GGangzzhe and the Solar dataset.
is single time series with 500 observations and 214 components. Solar is single time series with 52560 observations and 137 components.
For each of these datasets, the authors drop either 30% 50% of the observed values. Importantly, components are dropped independent from each other, such that at any given observation time, components (likely) only partially observed.
38

---

Published as a conference paper at ICLR 2025
For experiments, the rre-peeeesed datasets provided by the authors' including the sampled mask indicating dropped values.
G.1.2 MODELLING AND RESULTS
We apply FIM-((w.n. m) with different number windows to both datasets and both obser- vation ratios and compare our results to the baselines extracted from Tables and from Fang et al. (2024).
Following the experimental setup of Fang et al. (2024), we compute the metrics RMSE and MAE only the missing values. Table 18 contains the results with an observation ratio of 50% and Table 19 the results for an observation ratio of 70%.
The two datasets are of vastly different lengths. Still, with a suitable number of windows, FIM-6 can impute the missing values well in both datasets. It even outperforms all available baselines, even though they have been trained on the datasets.
The center plot of Figure 3 displays partial time series of one component of the Cuaffc-Guangzhou dataset, including the missing values and imputation by FI-EC. We can see that FIM-& captures the local patterns well and uses them effectively to impute the missing values.
Table 18:
datasets from Fang et al. (2024) with observed ratio 50%. Baselines have be extracted from Table 2 of Fang et al. (2024). We highlight the overall best model per dataset in bold and underline our best-performing model.
Performance
Traffic- Praffi-GuangZZhor
Solar-Power
Method
RMSE
MAE
RMSE
MAE
Simple.ccn
9.852
7.791
3.213
2.212
BRITS
4.874
3.335
2.842
1.985
NAOMI
5.986
4.543
2.918
2.112
SAITS
4.839
3.391
2.791
1.827
TIDER
4.708
3.469
1.679
0.838
Multi-Task GP
4.887
3.530
2.847
1.706
GP-VAE
4.844
3.419
3.720
1,810
CSDI
4.813
3,202
2.276
0.804
CSBI
4.790
3.182
2.097
1.033
HayOTIDE-fix weight
11.032
9,294
5.245
2.153
BayOTIDE-trend only
4.188
2,875
1.789
0.791
BayOTIDE
3,820
2.687
1.699
0.734
F FIM-f(w.n.. 1)
10.614
8.113
FIM-E(w.n. 2)
7.690
5,648
FIM-E(w n. = 4)
5.377
3.861
FIM-L(w.n. = 8)
4.440
3,090
FIM-f(w.n. ... 16)
3.741
2.562
FIM-E(w.n. 32)
3.600
2.427
FIM-E(w.n. 800)
2.044
1.188
FIM-É(W. n. 1600)
1.667
0.755
F FIM-E(w.n. = 3200)
1.550
0.595
G.2 COMPARING AGAINST SAITS
G.2.1 DATA DESCRIPTION AND PRE- PROCESSING
Du et al. (2023) collect four time series datasets (PhysioNer-200 20011 Air- AiirQQuaity, Electricity and ETT) and benchmark a range of imputation models, including their own, on several imputation tasks. For detailed information about each of these datasets, including the sequence lengths and number of samples, we refer the reader to Table of Du et al. (2023).
Phhttps://o thub. com /uamm/uuaauu--aaaa fang/ BayOTi
39

---

Published as a conference paper at ICLR 2025
Table 19: Performance on two datasets from Fang et al. (2024) with observed ratio 70%. Baselines have be extracted from Table 3 of Fang et al. (2024). We highlight the overall best model per dataset in bold and underline our best-performing model.
Traftic-Gannghhou
Solar-Power
Method
RMSE
MAE
RMSE
MAE
SimpleMean
10.141
8.132
3.156
2.319
BRITS
4.416
3.003
2.617
1.861
NAOMI
5.173
4.013
2.702
2.003
SAITS
4,407
3,025
2.359
1,575
TIDER
4.168
3.098
1.676
0.874
Multi Task GP
4.471
3.223
2.618
1.418
GP-VAE
4.373
3.156
3,561
1.723
CSDI
4.301
2,991
2.132
1.045
CSBI
4.201
2,955
1.987
0.926
BayOTIDE-fix weight
13.319
9.290
5.238
2.026
BayOTIDE-tree only
4.002
2.759
1.651
0.712
BayOTIDE
3.724
2.611
1.621
0.709
FIM-L(w.n.
12.533
9.680
FIM-&(w.r) = 2)
9.020
6,810
FIM-E(w.n. = 4)
5.584
4,024
FIM-E(w.r) m. 8)
4.401
3,084
FIM-((w..r 16)
3.387
2.360
PIIM-((W.n. 32)
3.087
2.148
FIM-E(w n. = 800)
2.016
1.181
F FIM-E(w.n. = 1600)
1,470
0.655
FIM- 8 (w... 3200)
1.282
0.474
Here, we only consider one of their imputation tasks: point-wise imputation of 10% missing ob- servations. The authors provide pre-processed datasets for this pattern which we use for the following evaluation. Note that their pre-processed data includes pre-sampled masks for the missing values, to enable fair comparison to their results.
G.2.2 MODELLING AND RESULTS
Du at. al. (2023) report MAE, RMSE and MRE for all datasets and wide range imputation models and methods in their Table 2, which we include here for completeness. We apply FIM-((w.n. = m) with different numbers windows m and FIM-((o.n. n) for different window lengths n to all datasets and report our results, together with the baselines, in Table 20.
With the right number or size of windows, FIM-( competitive on three out of four datasets. The PhysioNet-2012 dataset is naturally sparse, where some components of time series only contain few observations. As a zero-shot approach that employs a channel independent strategy, FIM does not perform well in such situations. In other words, PhysioNet-2012 demonstrates the limitations of FIM-l outlined in Appendix C.6.
G.3 TSI-BENCH
G.3.1 DATA DESCRIPTION AND PRE- PROCESSEGS
TSI-Bench is a collection of eight time series imputation datasets from different application do- mains: air-quality, traffic, energy and healthcare, assembled by Du et al. (2024). Its authors provide pre- pro-psecedeed splits of each 1 for several imputation patterns.
Here, consider their 10% and 50% point missing observation pattern. Note that the pre- rre-procesed datasets include sampled masks for missing values, enabling a fair comparison to their results.
B phttpps:/ggithbb..mm/////// Phttps: thub. htttt//////gtaab..oo////////////
40

---

Published as a conference paper at ICLR 2025
Table 20: Performance on datasets from Du et al. (2023) with 10% imputation data. Baselines have be extracted from Table of Du et al. (2023). We highlight the overall best model per dataset in bold and underline our best-performing model.
PhysioNet-2012
Air-Quality RMSE
Electricity RYSSE:
ETT RMSE
Method
MAR
RYSE:
MRF
MAE
MRE
MAR
MRR
MAE
MRE
Median
0.726
0.968
103.5%
0.763
1.175
107.4%
2,056
2,732
110.1%
1.145
1.847
139.1%
Last
0.862
1.207
123.0%
0.967
1.408
136.38
1,006
1,535
53.9%
1.007
1.365
96.45
1.179
EPGAN
0.702
0.964
100.1%
0,750
1.126
105.6%
0,584
0.703
89.0%%
ANN
0,776
0.294
0.643
06.01%
0.376
0.428
BRITS
0.256
0.767
36.39%
0.155
01.225
21,698
0.347
1.322
45.39%
0.130
0.259
12.59
Transfommer
0.190
0.445
26.91
0.158
0.521
22.3%6
0.823
1,301
44.0%
0,114
0.173
10.9%
SAfTS-base
0.521
44.0%
SAITS
0.186
8.431
26.6%
0.137
0.518
19.3%
0.735
1,162
39.45
0.092
0.139
8.8%
FI8- (....
0.443
0.899
63.3%
0.167
0.433
23.0%
0.118
0.251
6.3%
0.108
0.178
10.5%
f(w...
0.146
4.9%
f(w.n.
58.0%
0:145
4.2%
FIN-t(w.n. 5)
0.409
0,751
58.4%
0.185
0.541
26.0%
0.011
0.166
3.8%
0.134
0.252
13.0%
FIN f(w... 16)
0.499
0.842
68.5%
0.394
0.747
55.7%
0,700
0.162
3.7%
0.489
0.912
47.5%
3.9%
FIN-t(o.n. =4)
0364
0.784
51.0%
0.137
0.435
19.3%
0,009
0.100
3.7%
0.106
0.134
10.3%
FIN-f(o.n. =6)
04483
1.047
69.0%
0.163
0,460
23.1%
0.126
0.329
6.6%
0.117
0.194
11.0%
FIN n..-.(.... n.m.=
0.478
1.083
68 756
0.165
0.443
23.39
0.118
0.311
6.3%
0.114
0.195
11.19
FIN-f(o.n. 12]
FIN-f(o.n. m. 12)
0.477
1.016
63.49
0.165
0.447
23.35
0.112
0.285
6.0%
0.109
0.182
10.85
Io... 16)
0.453
0.987
64.7%
0.167
0.433
23.6%6
0.111
0.271
6.0%
0,108
0.178
10.5%
For details about each individual time series dataset, including the pre-processing, we refer the reader Appendix A of (Du al., 2024).
G.3.2 MODELLING AND RESULTS
We apply F FIM-((w.n. m) with different numbers of windows and FIM-((o.n. n) for different window lengths n to all available datasets both the 10% and 50% missingness patterns. In accordance with Du et al. (2024), compute the performance metrics only the missing values.
Table 21 contains the MAE of all baseline models (extracted from Table of (Du et al., 2024)) and F FIM-f in the 10% missingness pattern. Table 22 contains the MAE and MSE of all baseline models (extracted from Table 11 of (Du et al., 2024)) and FTM- in the 50% missingness pattern.
Because of the short, and sometimes very sparse data (in particular in the 50% missingness pattern), we experimented with different application strategies for F FIM- and report the corresponding results in both tables.
Let us now compare the two window specifying methods, FIM- f(w.n. = m) and FIM-((o... = n) on the TSI-Bench datasets. Here, specifying the windows by their number of observations performs better overall, in particular in the 50% missingness pattern of Table 22. In the 10% missingness pattern, the (relative) difference between the two approaches is smaller, although still present.
With its windowing scheme, FIM- is not limited by the length of time series it can process. Some datasets in the TSI-Bench based on a single long time series, that is split up into smaller chunks during pre-processing, to accommodate methods with such limitations. For these datasets, we ex- perimented with reassembling the original time series via concatenation, before applying FIM-(. We denote these models by Nong-FIM-((w..r. =m).
The reassembled time series could provide more context for our zero-shot method to better extract the (local) patterns. Table 22 reveals that such strategy can indeed improve the performance, in particular in terms of the MSE. However, it does not improve the performance all datasets and for all metrics.
41

---

Table 21: MAE on datasets from Du et al. (2024) with 10% point missingness. Baselines have been extracted from Table 2 of Du et al. (2024). Parenthesis indicate the standard deviation of five training rounds of neural models. We highlight the overall best model per dataset in bold and underline our best-performing model.
Published as a conference paper at ICLR 2025
Method
BeljingAir
ItalyAir
PeMS
Pedestrian
ET1.h1
Electricity
PhystoNet2012
PhysdeNet2019
ITransfommrr
0.133
0.233 (0.014)
0.226 (0.001)
0.148 (0 000))
0.263 (0.004)
0.571 (0.178)
0.379 (0.002)
0:462 (0) (8.006)
SAITS
0.155 (0.004
0.185 (8.010)
0.287 (0.001)
0.131 (0.006)
0.144 (0.000)
1,377 (1.026)
0.257 (0.019)
0.352 (0.0K15)
Nonstatioary
0.209 (0.)
0.266 (0.007)
0.331 (0.017)
0.453 (0
0.359 (0.03))
0.213 (0.010
0,410 (0.002)
0.458 (0.001)
ETSformer
0.187 (0.02)
0:250
0.947 (0.006)
0.237 (0.00)
0:412 (0.005
0.373
0.451
PatchTS
0.358 (0.011)
0.274 (0:036)
0.330 (0:013)
0.126 (0 003)
0.340
0.550 (0.039)
0.301 (0.011)
(0.00))
Grossfummer
0.184 (0.004)
0.246 (0.011)
0.337 (0.007)
0.119 (0.005)
0.233 (0 (0.008)
0.540 (0.034)
0.525 (0.302)
0.378 (0.007)
Informer
0.148 (0.002)
0.205 (0.008)
0.302 (0:005)
154 (0.00))
0.167 (0.006
1,291 (0.031)
0.297 (0.001)
0.403 (0.007)
Autofommer
0.257 (0.012)
0.598 (0.074)
0.197
0.748
(0.902)
Pyraformer
0.178 (0.004)
0.217 (0.00)
0.285 (0.003)
0.153 (0.012)
0.182 (0.00)
1.096 (0.033)
0.294 002)
0.387 (0.004)
Transfemmer
0.142 (0.00))
0.191 (0.010)
0.294 (0.002)
0.136 (0.009)
0.178 (0.01.5
1,316 (1.036)
0,250 (0.005
0.341 (0.002)
BRITS
0.127 (0.01))
0.235 (0.007)
0.149
(0.002)
0.971 (0.016)
0.297
0,355 (0.0011)
MRNN
0,568 (0.002)
0.638 (0.003)
0.624 (0,000)
0,735 (0.001)
0.789 (3.019)
1,824 (0.005)
0.708 (0.029)
(0.015)
GRUD
0.233 (0.002)
0.368 (0.012)
0.355 (0.002)
0.204 (0.008)
0.325 (0.004)
0,976 (0.015)
0,450 (0.004)
0.471 (0.001)
TimrsNei
0,230 (0,00)
0,260 (0.004
0.312 (0:001)
0.157 (0.00)
0,251000083
1.011 (0.016)
0,353 (0.001)
0,394 (0.0033
MICN
(0,001)
(0.004)
(0.003)
(0.010)
0.392 (0.006)
SCINei
0.191 (0.011)
0,268 (0.010)
0.487 (0.101)
0.149 (0.01))
(0.099
0,581 (0.015)
0.341 (0.005)
0.427 (0.002)
StemKiNN
0.361 (0.0022
0.200 (0.008)
0.495 (0:079)
0.127 (0.000)
0.248 (0.012)
1,300 (0.078)
0,351 (0.001)
0.416 (0.0K2)
0.211 (0008)
0.273 (0,005)
0.396 (0:02))
0.138 (0005)
0.262 (0.099)
0.718 (0.023)
0.315 (0.008)
0:406 (0.0173
Kuupa
0.363 (0.108)
0.307 (0.041)
0.532 (0.122)
0.173 (0 020)
0.435 (0. (3.)
1.309 (0.531)
0.413 (0000)
0.451 (0. (0.019)
TM incar
0.215 (0.016)
0.242 (0 (0.009)
0. 0.362 (0.009)
0.179 (0 004)
0.227 (0.006)
0.519 (0.006)
0.370 (0.00 (0.000)
0.432 (0.001)
FILM
0.318 (0.010)
(0.011)
0.764 (0.064)
0.413 (0.010)
0,583 (0.008)
0.834 (0.031)
0.458 (0.001)
0,494 (0.003)
(0.01))
0.539 (0.418)
0.33% (0:047)
0:231 (0.064)
0.151 (0.000)
1,483 (0,439)
0.252 (0.002))
0:408 (0.019)
US--AN
0.137 (0.002)
0.266 (0.012)
0.296 (0:001)
0.151 (0.016)
0.458 (0 (..000
0.958 (0.009)
0,310 (0.007)
0. (0.002)
GI-VAE
0.240 (0.066)
0.369 (0.012)
0.341 (0:007)
0.319 (0.010)
0.329 (0.017)
1,152 (0.0774
0.445 (0.005)
0,562 (0.0044)
0.721
0.798
0.728
0,737
0.422
0,708
0.762
Median
0.581
0,518
0.778
0.667
0.710
0.408
0.690
0.747
LOCF
0.188
0,233
0.375
0.257
0,315
0,104
0,469
0.478
Lineir
0,135
0.111
0.167
0.197
0.065
0,366
0,387
FIM flm.n.
0.162
0.317
0.234
0.105
0,489
0.448
FIM f[w.n. 2)
0.131
0,150
0.247
1.194
0.265
0.088
0.460
0.434
FIM-F(w.nt 4)
0.127
0.196
02237
0. 170
0,234
0.078
0.444
0.476
f[w.n.
0.165
0311
0.245
0.234
0.072
0.445
FIM f[w.n. 16]
0,401
0.444
0.449
0.071
0.521
0.538
FIM- FIM-f(w.n. fm 32)
0.501
1.505
0.580
0,612
Fin-f(o.n. =4
0,122
62088
0.171
0,230
01069
0.402
0.418
FIN U.N-t[o.n.
0277
0.177
0.295
0.119
0.536
0,499
FIN-t(o.n. = 8)
0.145
0.162
0.301
0.184
0,302
0,117
0.532
0.487
fIn-f(o.n. = =8]
0,144
0.162
0.307
0.186
0,308
0,114
0,523
0.481
ZIN-f(o.n.
0.147
0.162
0272
0.295
0.469
fl.. 12)
0.162
0274
0.110
0.507
0.463
FIM FIM-f [o.n. = 16)
0.148
0.162
0.317
0.234
0,518
0,110
0.502
0.461
Long-FTH S-N--F---fopp FIN-f(w... 16)
0.425
Long f(w.n. 322)
0.455
0.733
Long-FIN-t[w.r. = 664)
0.397
0,235
0.705
0.495
long-FiN-f(w.n.- = =128)
0.258
0.174
0.477
0.291
Long-FIN-t(w.er -256]
0.167
0,145
0.273
0,240
0.172
Long-FIN-Rlw.n. -512)
0217
0.225
0,101
long-FTH-f[w.r. -1024)
0.120
0.081
Long-FIN-([w.r. = =2048)
0.072
Long-FIN-t[w.r. 4096)
0.070

---

Table 22: MAE and MSE on datasets from Du et al. (2024) with 50% point missingness. Baselines have been extracted from Table 11 of Du et al. (2024). Parenthesis indicate the standard deviation of five training rounds of neural models. We highlight the overall best model per dataset in bold and underline our best- dest-performing model.
Published as a conference paper at ICLR 2025
BeijiugAir
ItalyAir
PeMS
ET
Electricity
Pedestrian
Method
MSE
MAE
MAE
MSE
MSE
MAE
MSE
Transformer
0.163 (0.003)
0.233 (0.00))
0.321 (0.007)
0.327 (0.011)
(0:007)
0.539 (0.016)
(0.002)
0.233 (0
0.893
1.584 (0, (60)))
(0,65))
SAITS
0.94 (0 003)
0.193 (0.007)
0.285 (0.010)
0.236 (0.014)
0.302 (0.001)
0.595 (0.003)
0.223 (0.007)
0.107 (0 (0.005)
1.399 (0.069)
3.837 (0.316)
0.208 (0.01)
0.392 (0.927)
Monelarisetar
0:231 (0.001)
0.271 (0.007)
0.314 (0.005)
0.361 (0.010)
0.394 (0.013)
0.688 (0.016)
0.382 (0/004)
0.292 (0.00)
02177 (0.031)
0.191 (0.048)
0.487 (0.033)
0,859 (0.099)
BTSfommer
0,249 (0.004)
0.261 (0.009)
0.401 (0.007)
0:421 (0.011)
0.366 (0,007)
0.586 (0.009)
0.364 (0.01.3
0.269 (0.022)
0.78 (0.008)
1.687 (0.02)
0.330 (0.00)
(..19(0.010)
PatchTS
(0.009)
0.206 (0.007)
(0.011)
(0.110)
(0.006)
0.609
(0.023)
0.149 (0.017)
0.856
141)
(0
0.35) (0.005)
Crossformer
0.215 (0.007)
0.221 (0.004)
0.325 (0.009)
0.293 (0.009)
0.357 (0.003)
0.607 (0.008)
(0.021)
0.146 (0.01))
0.980 (0.344)
2:255 (1.656)
0.191 (0.008)
0,356 (0.014)
Informer
0.184 (0.005)
0.213 (0.003)
0.304 (0.007)
0,247 (0.015)
0.330 (0.005)
0.600 (0.009)
0.279 (0.008)
0.162 (0.007)
1.277 (0.026)
3,239 (0.080)
0.210 (0.006)
0.37% (0.02))
Antoformer
0,898 (0.001)
1.554 (0.003)
0.833 (0.017)
1,800 (0.044)
0.602 (0.068)
1.242 173)
0.964 (0.008)
1,553 (0.025)
2.164 (0.001)
5,092 (0.00)
1.033 (0.015)
2.273 (0.067)
Pyraformer
0.198 (0.005)
0.223 (0.011)
(0.012)
(0.015)
(0.002)
0.580 (0.004)
0.29) (0.026)
(0.021)
1.131
2.711 (0.079)
(0.006)
(0.907)
Tranformer
0.85 (0.007)
0.192 (0.005)
(0.01)
0,230 (0.017)
0,316 (0.004)
0,588 (0.005)
0.274 (0.012)
0.162 (0.017
1.365 (0.034)
3.554 (0.085)
0.194 (0.0144
0.342 (0,033)
BRIIS
0.169 (0001)
0.194 (00023)
0,321 (0.008)
0.283 (0,005)
0.287 (0.001)
0,561 (0.02)
0.238 (0000)
0.127 (0.004)
1.124 (0.10))
2828 (0.023)
0.23 (0.077)
0.433 (0.021
MRNN
0.633 (0.006)
(0.00))
0.734
(0.006)
0.645 (0:001)
(0
0.816
1,219 (0.004)
1810
(0.01))
(0.001)
1.258 (0.000)
GRUD
0.29 (0.001)
0.303 (0.002)
0,476 (0.009)
(0:011)
0.372 (0.0022
0619 (0.002)
0.477 (0.008)
0.357 (0.005)
1 (0.011)
2.45% (0.034)
0 (0.005)
0.507 (0.007)
TimesSNt
0.265 (00015)
0:233 (0007))
0,570 (0.010)
0,323 (0.012)
0,348 (0.002)
0.367 (0001)
0.339 (8004)
0.210
1.131 (0017)
2.44 (0.077)
0.209 (0006)
0,392 (0.017)
MICN
0:456 (0.006)
0.553 (0.013)
0,548 (0.001
0.852 (0.012)
0.392 (0.006)
0.608 (0.010)
0,606 (0.073)
0.688 (0.152)
0.965 (0.008)
2.018 (0.032)
SCTNet
(0:012)
0.230 (0.036)
0,357 (0.008)
(0.006)
(0.093)
0,849
(0.013)
0.778 (0.023)
(0 (15)
0 0.251 (0.005)
0.391 (0.015)
StamGNNN
0.186 (0.004)
0.363 (0.005)
0.307 (0.014)
0.280 (0.019)
0.446 (0.021)
0.562 (0.064)
0.325 (0.019)
0.200 (0.025)
1.362 (0.187)
5,803 (0 (0.920)
0.200 (0.099
0.343 (0.014)
FrrTS
02335 (0:015)
0.246 (0.010)
0,349 (0.015)
0,345 (0.055)
0,422 (0:019)
0.686 (0.027)
(0019)
0,195 (0.030)
0871 (0.084)
1.326 (025)
0.234 (0.004)
0.314 (.314(0.0066
Koopu
0.373 (0.079)
0.445 (0.119)
0.345 (0.032)
(0.056)
0.506 (0:114)
0.855 10.
0,515 (0.159)
0,577 (0.351)
1.255 (0.250)
7,290 11.677)
0.246 (0.017)
0.330 (0:030)
DLinear
0.245 (0:005)
0.242 (0.006)
0.340 (0.004)
0.337 (0.005)
0.389 (0:013)
0.601 (0.020)
0.311 (0.000)
0.186 (0. (0.000)
0.734 (0.011)
0.988 (0.038))
0.310 (0.002)
0:455 (0.906)
HILM
0,331 10.009)
0.409 (0.009)
0,402 (0.018)
0.468 (0.052)
0.781 (0.0599
1.499 124)
0.589 (30005)
0.793 (0.00B)
09907 (0.024)
1.434 (0.078)
0.453 (0.007)
0.664 (0.000
CSDI
(0.007)
(0.
0.958
183)
(0:042)
0.651
(0.06))
(0.011)
21,850 140)
1.117
US-GAN
0.192 (0.001)
0.187 (0.005)
0.357 (0.009)
0,278 (0:011)
0.530 (D:001)
0.566 (0.001)
0,755 (0.973)
(3.955)
1.119 (0.007)
2.600 (0.018)
(0.005)
0.528 (0.007)
GP-VAE
0,258 (0.004)
0.234 (0.008)
0,453 (0.014)
0,495 (0.022)
0,346 (0.0155
0.617 (0.015)
0,414 (0.013)
0,301 (0.011)
1,099 (0.032)
2.973 (0 (0.040)
0.451 (0.022))
0,677 (0.031)
Meun
0,708
1,073
0,368
1,096
0.799
1,416
0.735
0.971
0.123
0.581
0,763
1.25%
Median
0.677
1.143
0.583
1.116
0.777
1.476
0:708
1,022
0.408
0/627
1.386
LOCF
0.264
0.429
0.345
0.511
0,547
1.094
0.425
0.491
0.140
0.181
0.365
0.635
Linear
0.165
0,231
0,214
0.252
0.343
0.539
0,267
0.178
0
0.035
0.247
0.279
FIM- (wr..
0.325
0.254
0
0.496
0.370
0.110
0.057
0.366
F10- m..... 2)
0.185
0,350
0.236
0.371
0.434
0.840
0.337
0.275
0.006
0.047
0.305
FIM-4(w. a. = 47
0.990
0
0.352
0.609
0.384
0.714
0.317
0,251
0.093
0.045
0,273
0
FIM-M W..n. 81
0
0.520
0,563
1 0.49
0.473
0.866
0,326
0,305
0.091
0.046
0.368
0.450
Mw.n.
0.886
0.626
1.1%6
0.445
0.508
0.122
0,140
0.763
FIM-Lww.m. 32)
0.650
0.91
0.33
1.499
FIM- [oon. = 4)
0.100
01.445
0:225
0.249
0.375
0.70.
0.279
0.19K
0103
0.274
0,272
Mo...
0.340
0,506
0.250
1,205
0.469
0.310
0.336
0.287
FIM Mo....
0,422
0.254
0.287
0.583
1.327
0.467
0.661
0.199
0.294
0.314
21.--(o.n. = 8)
0:233
0
0.254
0.287
0.568
1.235
0.467
0.649
0.191
0.257
0.502
0.293
PIM Ko.n. 12)
0.308
0.370
0.254
0.287
0.533
1,058
0.420
0.453
0.187
0.227
0.327
0.334
foo.. 122)
0.212
0.410
0.254
0.277
0.546
1.136
0.448
0.182
0.211
0.524
0.325
FINs-6(o.n. 16)
0.198
0.325
0.254
0.287
0.496
0.960
0.417
0.500
0.198
0.240
0.343
0,366
Long-F 101-2 w.m. = 16]
0.414
0.158
Iong-F2M-t(w.nn. 32)
0.358
0.38 0.383
0.905
1,692
0.693
1.084
Long-E28-E(w.n. 644
0.512
0.283
0.709
1.193
0.576
Long-FI8-E(w.m. = 128)
0.266
0.397
0.238
0.245
0,516
0.848
0.362
0,303
1ong-Fi8-t(w.n. = 256)
0.201
0.279
0.215
0.231
0.417
0.721
0.328
0.281
0.157
0.103
1ong-F28-((w.n. 512)
0.177
0.265
0.365
9.651
0,310
0.265
0.106
0.052
Long-F38-E(w.n. 1024)
0.279
0.011
0.044
Long-FI8-I(w... 2048)
0.089
0.043
Long-FI8-t(w:) = 4098)
0.097
0.066

---

Published as a conference paper at ICLR 2025
H TEMPORAL MISSING PATTERN IMPUTATION: ADDITIONAL RESULTS
H.1 MOTION CAPTURE
H.1.1 DATA DESCRIPTION AND PRO-PROCESSING
Let us now consider the imputation problem setup proposed by Heinonen et al. (2018) on a human motion capture dataset, consisting of 50-dimensional pose measurements of walking subjects. We take the data provided by Yildiz et al. (2019), which was pre- rre-processed according to previous work of Wang al. (2007).
The dataset contains 43 trajectories of maximal length of 125. Following Heinonen et al. (2018), we remove 20% out of the center of each trajectory. To compare to Heinonen al. (2018), we also apply the provided PCA projection and consider the first 3 PCA components. In the following, we call this the PCA setup.
Because F IM can be applied to arbitrary dimensional data, due to our channel independent strategy (see Appendix D.7), we also consider the dimensional data directly. In the following, we call this the No PCA setup.
H.1.2 TRAINING LATENTODE
To obtain another baseline model, train LatentODE (Rubanova et al., 2019) the first PCA components on the likelihood outside the interpolation window, i.c. it has never seen and is not trained on the (missing) data inside the imputation window. We train and test all 43 trajectories, selecting the model based on the performance on all trajectories. This approach similar to the approach of models trained in Heinonen et al. (2018).
The training data is standardised for training. We use a LSTM as as the encoder and diagonal gaussian emission model, fixing its standard deviation to 0.01. We train for 120.000 gradient descent iterations over the whole dataset, using Adam W with learning rate le-3 and weight decay 1e-2, To 2 help the model learn, we slowly anneal the input time series length over the initial 60.000 epochs, starting at 25 observations.
We did ablation over the hidden size of the model. See Table 14 for the final hyperparameters. The model trained roughly 3.5 hours on A100 40GB GPU.
H.1.3 MODELLING AND RESULTS
In accordance with Heinonen et al. (2018), we compute the performance metrics in 50 dimensional space and only on the missing points inside the imputation gap. We evaluate FIM, LatentODE and also cubic spline composed with Savitzky-Golay filters (Savitzky & Golay, 1964) in the PCA setup. Table reports the RMSE, including baselines extracted from Heinonen et al. (2018), while Table 24 reports the MAE.
FIM performs almost as well as the specialized LatentODE model, especially when considering the large standard deviation for both approaches. These results indicate that a zero-shot imputation approach is indeed viable.
PCA dimensionality reduction induces some level of intrinsic emrrr In contrast to baseline models, FIM is not restricted by dimensionality of the data. By applying our model to the 50 dimensional data without PCA projection, can, in principle, avoid this emror.
We report the results of this No PCA setup in Tables 23 and 24. Indeed, our model is performing better without the PCA projection, showcasing the strengths of our flexible, zero-shot methodology
44

---

Published as a conference paper at ICLR 2025
Table 23: RMSE of imputation in the Motion Capture and Navier Stokes datasets. The RMSE is calculated only in the imputation window, in accordance with Heinonen et al. (2018).
Motion Capture
Navier Stokes
Model
Filter
PCA
No PCA
PCA
No PCA
LatentODE
3.066 1: 1.767
0.133 0.053
Cubic spline
7.333 1 2.49
8.95 + 2.78
0.174 1 0.003
0.186 1 0.004
Cubic spline
Savgol(15, 3)
6.317 = 1.813
6.752 = 1.676
0.151 + 0.003
0.143 + 0.002
Cubic spline
Savgol( 8, 3)
7.078 = 2.147
8.271 2.282
0.148 0.003
0.147 0.003
Cubic spline
Savgol( 4,
7.251 = 2.184
8.783 = 2.508
0.171 + 0.003
0.18 + 0.004
npODE
3.94 + 3.50
GPDM
5.31 + 3.39
VGPLVM
3.91 + 1.89
FIM
3.271 = 1.22
2.977 0.96
0.103 0.004 0.971 = 0.005
Table 24: MAE imputation in the Motion Capture and Navier Stokes datasets. The MAE is calculated only in the imputation window, in accordance with Heinonen et al. (2018).
Motion Capture
Navier Stokes
Model
Filter
PCA
No PCA
PCA
No PCA
LatentODE
1,658 + 0.989
0.076 + 0.03
Cubic spline
3.362 + 1.175
4.209 + 1.436
0.085 = 0.003
0.083 + 0.003
Cubic spline
Savgol 15, 3)
2.897 0.871
2.998 + 0.881
0.084 0.00
0.075 0.002
Cubic spline
Savgol( 3)
3.229 1.029
3.695 + 1.22
0.081 + 0.003
0.072 0.002
Cubic spline
Savgol( 4, 3)
3.298 + 1.035
4.061 + 1.354
0.085 + 0.003
0.082 + 0.003
FIM
1.765 + 0.627
1.611 + 0.453
0.062 + 0.003
0.051 + 0.002
H.2
NAVIER STOKES
H.2.1
DATA DESCRIPTION AND PRE- PROCESSING
As an application of FIM to high-dimensional data, we consider the simulation of two-dimensional, incompressible Davier-Stokes equation from (Course Nair, 2023) The equation is simulated on wwo-dimensional grid of size 199 x 1499 for total of 596, 602 states. Following (Course & Nair, 2023), we remove the first 20% of the trajectory for warmup and are left with 2441 simulation steps.
We use this preprocessed simulation create another imputation dataset, with similar setup in the motion capture dataset. We drop the last observation and cut the remaining trajectory into 61 time series of length 40. Then remove the central 20% of each time series, creating temporal missing pattern imputation task.
While model can handle this high-dimensional data, as will show below, we need to apply PCA dimensionality reduction to train a (specialized) baseline model. Following (Course & Nair, 2023), we project the data to 38 dimensional space with randomised PCA, which already captures the high dimensional dynamics well.
In the following, we will again refer to these two setups as PCA and No PCA respectively.
H.2.2 TRAINING LATENTODE
As a baseline, we train LatentODE (Rubanova et al., 2019) in the PCA setup on the likelihood outside of the imputation window, with the same reasoning for the motion capture imputation problem in Section H.1.2.
LatentODE is trained all 61 (standardised) time series, uses LSTM as the encoder and diag- onal gaussian emission model with fixed standard deviation of 0.01. We train for 300.000 gradient descent iterations the whole dataset, using AdamW with learning rate 1 and weight decay
Phttps://gg wwtphhc.w.w.ouseeeeeeeeeee h/siiiiiiivviiie
45

---

Published as a conference paper at ICLR 2025
Table 25: Performance comparison of different models across window counts.
Window count
Pre-trained
Fine-tuned noisy obs
Fine-tuned on clean obs
10.32 + 1.478
8.983 0.932
9.001 + 0.945
4
7.136 + 0.697
5.936 + 0.432
5.947 + 0.433
2.860 + 0.237
2.708 0.176
2.704 + 0.175
12
1.956 0.161
1.913 + 0.144
1.911 = 0.143
16
1.652 + 0.126
1.702 0.127
1.703 + 0.126
20
1.526 = 0.125
1,586 + 0.112
1.588 = 0.112
24
1.504 + 0.138
1.532 0.108
0.100
28
1.516 + 0.133
1.505 + 0.098
1.506 + 0.098
32
1.561 + 0.142
1.500 + 0.118
1.501 + 0.118
1e-2, We slowly anneal the input time series length over the initial 60.000 epochs, starting at just observations.
We studied some ablation over the hidden size. The search grid, including the final hyperparameters, shown in Table Table 14. The model trained roughly hours on a A100 40GB GPU.
H.2.3 MODELLING AND RESULTS
Let us first consider the PCA setup where we trained the LatentODE baseline. For both FIM and LatentODE we report the RMSE, in Table 23, and MAE, in Table 24, inside the imputation gap. In this particular task, F IM was again able to match the performance of LatentODE.
As FIM can be applied to data of any dimensionality, we also experimented with the No PCA setup, imputing the missing data 596, 602 dimensional space directly. As in the Motion Capture dataset, FIM performs even better without the errors induced by the PCA projection.
Finally, we report the computational load of inference in dataset, as it is considerably larger than in all other experiments, because of its high dimensionality. The application of F FMM took roughly 9 hours on A 100 40GB GPU. As a comparison, the application cubic spline on the same data took roughly 1.1 hours on 32 CPU cores.
I ON FINE-TUNING FIM
Although outside the scope of the present work, fine-tuning FIM to specific applications would be of great interest and value.
In practice, one does not have access to any ground-truth initial condition, nor to any time derivative (of the hidden interpolating function). However, the model can still be optimized to reconstruct available data. To illustrate the plausibility of this approach in controlled setting, have fine- tuned F FIM- to reconstruct the Lorentz system data we studied in Appendix F.4. The details this dataset can be found in that Appendix.
Setting up the fine-tuning:
1. FIM-( l takes as input the noisy observations along one path and returns (its best guess of) the interpolation function, which can be evaluated at any time. We fine-tune FIM- in two settings, namely
(a) to reconstruct the noisy input data, and
(b) to reconstruct the ground-truth, clean values the Lorentz system's solution, the input observation times.
Clearly, case (a) is the realistic setting. We consider case (b) as a control case, because FIM-& was pretrained reconstruct the clean path.
2. As we described in Subsection 3.2.2 of the main text, FIM-E can be used process either sets of up to Laux observations simultancously, or subsets (or windows) with fewer observations, fol- lowed combination of the local window estimates. This trick is what allows interpolate
46

---

Published as a conference paper at ICLR 2025
time series of any length and handle interpolation functions which lie outside our distribution of "simple" functions.
The number of windows into which we split the target dataset is a hyperparameter of the model, it clearly problem dependent. Indeed, Figure 8 of Appendix F.4.4 demonstrates the effect of modifying it on the Lorentz system data. To fine-tune FIM-C, we therefore need to specify the number windows hyperparameter. That is, we need to specify the temporal scale at which we want to fine-tune the model. In our simple experiment below we chose to fine-tune F FIM- on four (4) windows.
We fine-tuned FIM-& with learning rate of le-6, which was the learning rate at which we stopped pretraining. We fine-tuned it for two epochs only, on all components of the Lorenz systems (i.e. channel independence strategy).
Other details:
Table 25 (which expands on Fig. 8) contains our results.
The first thing we notice is that the performance of FIM-& clearly improves at the time scale (i.e. at the window count) on which it was fine-tuned. Indeed, we fine-tuned it on 4 windows. Very importantly, the fine-tuning process does not negatively perturb the interpolation performance at different time scales (i.c. at fewer or larger window counts). Therefore, there is no catastrophic forgetting. This observation alone demonstrates the plausibility of fine-tuning with FIM-C.
The second thing we notice is that fine-tuning on the noisy and clean data had a similar effect. The continuous nature of the interpolation function appears to make the fine-tuning robust to noise in the target data. However, we have also observed that if one fine-tunes FIM- on the noisy data for many more epochs, simply with larger learning rate, its performance different time scales does deteriorate. This occurs because the model starts fitting the noise. Therefore, to properly address fine-tuning within our framework, we would need define an emission model that handles the noise.
One could picture an emission model, (pre-)trained on our synthetic dataset, which would take as input both the interpolation function xtt) represented on some sensor grid, like in DeepONets and the noisy data .1..... Wt- The emission model would then output the probability of observing the noisy value y at any desire time т,, that is p(y: (xr(:)).
Coupling F FIM- together with this emission model would allow to finetune FIM-& (or the emission model, or both) to any target dataset. We plan to pursue this general direction in the near future.
47


Reviews:
[{"review_id": "IcDlGQVjjW", "reviewer": "Reviewer_6Udz", "summary": "This paper proposes a supervised learning framework for zero-shot time series imputation in dynamical systems of any dimensionality. The framework uses a synthetic data generation model based on two key assumptions: \n(1) time series with point-wise missing data have simple ODE-based interpolation solutions, and \n(2) time series with temporal missing patterns are locally simple. \nThe authors introduce two neural interpolation models, FIM and FIM-l, designed for point-wise and temporal missing patterns, respectively.\nThe approach is evaulated over 8 datasets", "questions": "— The FIM-l model only works for time series that follow a ‘simple’ distribution. Since existing imputation models already perform well on real-world and synthetic data, the practical value of the proposed model might be limited.\n— The authors did not explain clearly how the FIM could handle processes of any dimensionality, and how the performance of the model varies with the dimensionality in practice.", "limitations": "", "rating": 5, "confidence": 3, "soundness": 2, "presentation": 3, "contribution": 2, "strengths": "— The proposed FIM has great zero-shot imputation performance and achieves SOTA results on multiple datasets.\n— Unlike prior work, this paper proposes a zero-shot approach that can be used on processes of any dimensionality.\n— The synthetic data generation model offers a general method for creating meaningful synthetic time series, which could support future research.\n- The authors are honest about the limitation of the approach, especially those related to the synthetic distributions used in the analysis and setting (ODEs) as opposed to real-world data which can exhibit any (or no ) distrubution.", "weaknesses": "- Fig 2 is difficult to read, the colors are too  similar\n- The paper is too \"buzzwordy\", I would have  preferred more of a proper technical discussion\n- The results of the proposed method are not always close to SOTA or better than SOTA", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "qZd6YfIe3J", "reviewer": "Reviewer_nBHN", "summary": "This paper introduces foundation model for time-series imputation, using several parametric functions of ODEs for training the model. The authors generate many synthetic noisy, irregular observations to train recognition model in offline manner.  They use combination of objective function with respect to initial value and data points to train this model. This trained model can conduct zero-shot imputation on unseen time series dataset with different dimensions. They validate their method on various different benchmarks in several domains.", "questions": "- I think some ablation studies for measuring the generalization capability of this model will be beneficial.", "limitations": "", "rating": 8, "confidence": 3, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "- This paper is well written and organized. I enjoy reading this paper.\n- To the best of my knowledge, this paper is novel in that this firstly introduces the zero-shot imputation of time-series with dimension free inference. \n- I think the training objective and the way they give supervision to the model is also novel.\n- The experimental results are impressive, covering a wide range of datasets and scenarios enough to prove the zero-shot capability of this method. Additionally, the paper faithfully includes details for reproducing the experiments which is helpful", "weaknesses": "- Although this model impute the data in zero-shot manner, experimental results are not that impressive in Table 2. Also, the authors use reported value from the other papers which can question whether comparison is really fair.\n-  I think it would be better to include ODE based methods for the baseline considering concept of this paper.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "zkR3UIu4US", "reviewer": "Reviewer_8Y1g", "summary": "The paper develops a ODE-based foundation model for time series imputation. The idea is to simulate a large collection of incomplete time series including both point-wise and temporal missing patterns. A recognition model is trained to map such time series to the initial condition, and time derivatives at each time step, with both mean and log-variance prediction. The design of the recognition model is inspired by deep ONet. There is no need to fine tune the model on application-specific data. The experiments show the improvement over many sota methods specifically trained on application data.", "questions": "see above.", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1. very interesting work on time series imputation, potentially pointing out a new direction\n2. design of neural recognition model is novel and interesting, albeit inspired by deepOnet. Such \"borrowing\" is still interesting. \n3. zero-shot performance is surprisingly good", "weaknesses": "1. as admitted by the authors, one limitation is that when the actual application does not well match the assumption of the simulation data, the performance can deteriorate. Though the paper focuses on zero-shot prediction, it will be great to discuss  possible methods for fine tuning or adaptation to specific applications. Since your training assumes the ground-truth of the initial and time derivatives of each sampled trajectories are known, and the neural recognition models directly learn to fit them, such training cannot directly apply to real dataset where these ground-truth is unknown. \n2.  The design of the neural recognition model seems not well justified, especially for the branch-net component. why do you use RNN? Why not using attention mechanism instead? Given RNN/LSTM's are replaced by attention  nearly every where, the authors should explain their rationale of such choices, either theoretically or empirically. \n3. How are the hyperparameters selected? Why using 1024 dimensional embeddings? Are this done by a rigorous validation process? If so, how?  It will be good for the authors to specify the details of hyperparmeter selection and validation process.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "SOrfu8yXU9", "reviewer": "Reviewer_6Fr3", "summary": "This paper addresses the challenges of imputing missing time-series data. A framework for zero-shot time-series imputation is proposed with 1) a synthetic data generation model for sampling a set of ODE solutions and 2) a neural recognition model mapping the time series data onto parametric functions. Empirical results present improved performance of imputation compared to baseline models.", "questions": "Please find the questions in the Weaknesses section above.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "1. The problem of imputing time-series data is important and the idea of linking amortized inference and neural operator is interesting.\n2. The methodology of the foundation inference model is well-structured.\n3. Experiments were performed on the imputation of point-wise and temporal missing patterns to present the effectiveness of the proposed method.", "weaknesses": "1. I think the idea of zero-shot amortized inference framework is somewhat incremental. There have been works on amortized inference for few-shot time-series forecasting (e.g. [1] [2] [3][4]) in terms of lacking sufficient observation and learning single dynamics, where the idea of amortized inference learning the prior knowledge of dynamics is similar to the proposed method. Could the authors add a discussion about these works and the benefits of the proposed method? Also, it would be good if the authors could make a comparison of some of these methods in experiments.\n2. In Equation 3, 4, and 5, do $\\phi^{\\theta}_i$ and $\\psi^{\\theta}$ share the same parameters $\\theta$? Similar problem to the parameter $\\psi$ in FIM model in Equation 7 and 8.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:

Decision: Accept (Poster)