Trajectory Flow Matching with Applications to Clinical Time Series Modelling

Authors: xi  zhang
yuan  pu
yuki  kawamura
andrew  loza
yoshua  bengio
dennis l shung
alexander  tong
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-14 23:08:14

Abstract: Modeling stochastic and irregularly sampled time series is a challenging problem found in a wide range of applications, especially in medicine. Neural stochastic differential equations (Neural SDEs) are an attractive modeling technique for this problem, which parameterize the drift and diffusion terms of an SDE with neural networks. However, current algorithms for training Neural SDEs require backpropagation through the SDE dynamics, greatly limiting their scalability and stability. 
To address this, we propose **Trajectory Flow Matching** (TFM), which trains a Neural SDE in a *simulation-free* manner, bypassing backpropagation through the dynamics. TFM leverages the flow matching technique from generative modeling to model time series. In this work we first establish necessary conditions for TFM to learn time series data. Next, we present a reparameterization trick which improves training stability. Finally, we adapt TFM to the clinical time series setting, demonstrating improved performance on three clinical time series datasets both in terms of absolute performance and uncertainty prediction.

Main paper: 
 introduction: 
 Real world problems often involve systems that evolve continuously over time, yet these systems are usually noisy and irregularly sampled. In addition, real-world time series often relate to other covariates, leading to complex patterns such as intersecting trajectories. For instance, in the context of clinical trajectories in healthcare, patients' vital sign evolution can follow drastically different, crossing paths even if the initial measurements are similar, due to the influence of the covariates such as medication intervention and underlying health conditions. These covariates can be time-varying or static, and often sparse. Differential equation-based dynamical models are proficient at learning continuous variables without imputations [Chen et al., 2018 , Rubanova et al., 2019 , Kidger et al., 2021b] . Nevertheless, systems governed by ordinary differential equations (ODEs) or stochastic differential equations (SDEs) are unable to accommodate intersecting trajectories, and thus requires modifications such as augmentation or modelling higher-order derivatives [Dupont et al., 2019] . While ODEs model deterministic systems, SDEs contain a diffusion term and can better represent the inherent uncertainty and fluctuations present in many real world systems. However, fitting stochastic equations to real life Figure 1 : Trajectory Flow Matching trains both an estimator of the next timepoint (x θ (t, x)) and an estimation of the uncertainty (σ θ (t, x t )). Using the conditional flow matching framework, these can be used to predict the instantaneous velocity v θ (t, x t ) and future observations. Both flows are conditioned on past data x [t-h,t-1] and conditional variables c. data is challenging because they have thus far required time-consuming backpropagation through an SDE integration. In the domain of generative models, diffusion models [Ho et al., 2020 , Nichol and Dhariwal, 2021 , Song et al., 2021] and more recently flow matching models [Lipman et al., 2023 , Albergo et al., 2023 , Li et al., 2020] have had enormous success by training dynamical models in a simulation-free framework. The simulation-free framework facilitates the training of much larger models with significantly improved speed and stability. In this work we generalize simulation-free training for fitting stochastic differential equations to time-series data, to learn population trajectories while preserving individual characteristics with conditionals. We present this method as Trajectory Flow Matching. We demonstrate that our method outperforms current state of the art time series modelling architecture including RNN, ODE based and flow matching methods. We empirically demonstrate the utility of our method in clinical applications where hemodynamic trajectories are critical for ongoing dynamic monitoring and care. We applied our method to the following longitudinal electronic health record datasets from multiple clinical settings: medical intensive care unit (MICU) data of patients with sepsis, Emergency Department (ED) data of patients with acute gastrointestinal bleeding, and MICU data of patients with acute gastrointestinal bleeding. Our main contributions are:• We prove the conditions under which continuous time dynamics can be trained simulationfree using matching techniques. • We extend the approach to irregularly sampled trajectories with a time predictive loss and to estimate uncertainty using an uncertainty prediction loss. • We empirically demonstrate that our approach reduces the error by 15-83% when applied to the real world clinical data modelling. 
notation: 
 We consider the setting of a distribution of trajectories over R d denoted X := {x 1 , x 2 , . . . , x n } where each x i consists of a set of trajectories of length T i.e. x i := {x i 1 , x i 2 , . . . , x i T } with associated times t i := {t i 1 , t i 2 , . . . , t i T }. Let x i [t-h,t-1] denote a vector of the last h observed time points. We denote a (Lipschitz smooth) time dependent vector field conditioned on arbitrary conditions c ∈ R e v(t, x t , x [t-h,t-1] , c) → dx dt : ([0, 1], R d , R h×d , R e ) → R d with flow ϕ t (v) which induces the timedependent density p t = ϕ t (v) # (p 0 ) for any density p 0 : R d → R + with R d p 0 = 1. We also consider the coupling π(x 0 , x 1 ) which operates on the product space of marginal distributions p 0 , p 1 .A stochastic differential equation (SDE) can be expressed in terms of a smooth drift f : [0, T ]×R d → R d and diffusion g : [0, T ] × R d → R d 2 in the Ito sense as:dx t = f dt + g dW twhere W t : [0, T ] → R d is the d-dimensional Wiener process. A density p 0 (x 0 ) evolved according to an SDE induces a collection of marginal distributions p t (x t ) viewed as a function p : [0, T ] × R d → R + . In a Neural SDE [Li et al., 2020 , Kidger et al., 2021a,b] the drift and diffusion terms are parameterized with neural networks f θ (t, x t ) and g θ (t, x t ).dx t = f θ (t, x t )dt + g θ (t, x t )dW t (1)where the goal is to select θ to enforce x T ∼ X true for some distributional notion of similarity such as the Wasserstein distance [Kidger et al., 2021b] or Kullback-Leibler divergence [Li et al., 2020] . However, these objectives are simulation-based, requiring a backpropagation through an SDE solver, which suffers from severe speed and stability issues. While some issues such as memory and numerical truncation can be ameliorated using the adjoint state method and advanced numerical solvers [Kidger et al., 2021b] , optimization of Neural SDEs is still a significant issue. We note that in the special case of zero-diffusion (i.e. g θ (t, x t ) = 0) this reduces to a neural ordinary differential equation (Neural ODE) [Chen et al., 2018] , which is easier to optimize than SDEs, but still presents challenges to scalability.Matching algorithms are a simulation-free class of training algorithms which are able to bypass backpropagation through the solver during training by constructing the marginal distribution as a mixture of tractable conditional probability paths. The marginal density p t induced by eq. 1 evolves according to the Fokker-Plank equation (FPE):EQUATIONwhere ∆p t = ∇ Given this construction it can be shown that the minimizer ofL match (θ) := E t,q(z),pt(x|z) ∥f θ (t, x t ) -v t (x t |z)∥ 2 + λ 2 t ∥g θ (t, x t ) -σ t (x t |z)∥ 2 (3)satisfies the FPE of the marginal p t . This is especially useful in the generative modeling setting where q 0 is samplable noise (e.g. N (0, 1)) and q 1 is the data distribution. Then we can define z := (x 0 , x 1 ) as a tuple of noise and data with q(z) := q 0 (x 0 ) ⊗ q 1 (x 1 ). This makes eq. 3 optimize a model which will draw new samples according to the data distribution q 1 (x 1 ) usingx 0 ∼ q 0 ; x 1 = 1 0 f θ (t, x t )dt + g θ (t, x t )dW t (4)with the integration computed numerically using any off-the-shelf SDE solver. While this is guaranteed to preserve the distribution over time, it is not guaranteed to preserve the coupling of q 0 and q 1 (if given).In generative modeling random pairings [Liu et al., 2023c , Albergo and Vanden-Eijnden, 2023 , Albergo et al., 2023] or optimal transport [Tong et al., 2024 , Pooladian et al., 2023] pairings are constructed for the conditional distribution q(z). However, in some problems we would like to match pairs of points as is the case in image-to-image translation [Isola et al., 2017 , Liu et al., 2023a , Somnath et al., 2023] . In this case, training data comes as pairs (x 0 , x 1 ). In this case we set q(z) := q(x 0 , x 1 ) to be samples from these known pairs, and optimize eq. 3. While empirically, these models perform well, there are no guarantees that the coupling will be preserved outside of the special case when data comes from the (entropic) optimal transport coupling π * ε (q 0 , q 1 ) and defined as:π * ε (q 0 , q 1 ) = arg min π∈U (q0,q1) d(x 0 , x 1 ) 2 dπ(x 0 , x 1 ) + ε KL(π∥q 0 ⊗ q 1 ), (5)Algorithm 1 General Trajectory Flow Matching Input: Trajectories X , noise σ, initial network v θ . while Training do z ∼ U(X ), k ∼ U{1, T -1}, t ∼ U(0, 1) µ t ← tx k + (1 -t)x k+1 x t ∼ N (µ t , σ 2 t(1 -t)I) L TFM (θ) ← v θ (k + t, x t ) -x k+1 -xt 1-t 2 L σt (θ) ← ∥σ θ (k + t, x t ) -L TFM ∥ 2 θ ← Update(θ, ∇ θ L TFM (θ), ∇ θ L σt (θ)) return v θ , σ θwhere U (q 0 , q 1 ) is the set of admissible transport plans (i.e. joint distributions over x 0 and x 1 whose marginals are equal to q 0 and q 1 ) as shown in [Shi et al., 2023] for some regularization parameter ε ∈ R ≥0 . 
trajectory flow matching: 
 We now describe our simulation-free method to learn SDEs on time-series data using trajectory flow matching as summarized in Alg. 1. In the case of time series we need to ensure that trajectory couplings are preserved. We first set out a general algorithm for flow matching on vector fields in §3.1 then present a numerical reparameterization which we find stabilizes training in §3.2, a next observation prediction for irregularly sampled time series in §3.3, and finally present how to learn the noise in §3.4.In this section, we assume access to fully observed and evenly spaced trajectories X = (x 1 , x 2 , . . . , x n ) with x i := (x i 1 , x i 2 , . . . , x i T ) for clarity and notational simplicty. We note that our method is easily extensible to the more general setting of irregularly sampled trajectories. In this simplified case we letEQUATIONwhere U(X ) is the uniform empirical distribution over X , ⌈•⌉, ⌊•⌋ are the ceiling and floor functions, and N (•, •) is the multivariate normal distribution. This is a valid regression in the sense that a function minimized with Alg. 1 will return a stochastic process that will match the observed marginal distributions over time as shown in the following lemma. Lemma 3.1. The SDE dx t = u t (x|z)dt + σ 2 dW t where u t is defined in eq. 9 generates p t (x|z) in eq. 8 with initial condition p 0 := δ x1 where δ is the Dirac delta function.however, while useful, this is still insufficient for time series modeling, as it does not ensure coupling preservation. For intuition why this is an issue see figure 2 .In TFM we ensure that the couplings are preserved for history lengths h > 0. i.e.T -h , x T -h+1 , . . . , x T ) = π(x T -h , x T -h+1 , . . . , x T ).We first establish a method to ensure that these couplings are preserved allowing us to use simulation-free flow matching training for the time-series modeling task. Specifically, as long as the model takes as input (x T -h , x T -h+1 , . . . , x T ) in predicting the flow from T → T + 1, then there exists a function f θ (X T -h:T ) such that the coupling is preserved. Proposition 3.2 (Coupling Preservation). Under mild regulatory criteria on u t (•|z), p t , and q, ifE t∼U (0,T ),z∼q(z),c∼q(c|z),xt∼pt(xt|z) ∥u t (x t |z, c) -u t (x t |c)∥ 2 2 = 0and z, q(z), p t (x|z) and u t (x|z) are as defined in eqs. 6-9 then Π(u) ⋆ = Π ⋆ (x 1:T ).Where Π(u) ⋆ represents the coupling of a model which attains minimal loss according to eq. 3 and Π ⋆ (x 1:T ) is the coupling of the data distribution. Intuitively, as long as no two paths cross given conditionals c, then the coupling is preserved. In prior work c = ∅, and the coupling is only preserved in special cases such as eq. 5.We next enumerate three assumptions under which the coupling is guaranteed to be preserved at the optima. We note that these are (A1) When c = x 0 and there exists T : X → X such that T (x 0 ) = x 1 iff Π ⋆ (x 0 , x 1 ). We note that this is equivalent to asserting the existence of a Monge map T ⋆ for the coupling Π ⋆ .(A2) There exist no two trajectories x i , x j such that x i t = x j t for h consecutive observations and g = 0.(A3) Trajectories are associated with unique conditional vectors c independent of t.Even in cases when (A1)-(A3) may not hold exactly, TFM is a useful model and can often still learn useful models of the data. In some sense uniqueness up to some history length is enough as it shows TFM is as powerful as discrete-time autoregressive models. Proofs and further examples are available in §A.1.While flow matching generally predicts the flow, there is a target predicting equivalent namelygiven v θ (t, x) := x⌈t⌉ θ (t,xt)-xt ⌈t⌉-tand u t (x|z) := x ⌈t⌉ -xt ⌈t⌉-t which is equivalent to x 1 -x 0 when x t : tx 1 + (1 -t)x 0 then it is easy to show that the target predicting loss is equivalent to a timeweighted flow-matching loss. Specifically let the target predicting loss beL target (θ) = E t,q(z),pt(x|z) ∥x ⌈t⌉ θ (t, x) -x ⌈t⌉ ∥ 2 (10)then it is easy to show that Proposition 3.3. There exists a scaling function c(t) : R + → R such that L target (θ) = c(t)L match (θ).We next consider irregularly sampled time series of the form x i := (x i 1 , t i 1 ), (x i 2 , t i 2 ), . . . , (x i T , t i T ) with t i 1 < t i 2 < • • • < t i T with t next denoting the next timepoint observed after time t. In this case, when combined with the target predicting reparameterization in §3.2, we can predict the time till next observation. We therefore parameterize an auxiliary model h θ (t, x t ) : [0, T ] × R d → [0, T ] which predicts the next observation time. This is useful numerically, but also, perhaps more importantly, is useful in a clinical setting, where the spacing between measurements can be as informative as the measurements themselves [Allam et al., 2021] . h θ is trained to predict the time till the next observation with the time predictive loss:L tp (θ) = t∈T i ∥h θ (t, x t ) -(t next -t)∥ 2 2 (11)where t next is the time of the next measurement. This can be used in conjunction with the x next predictor to calculate the flow at time t asEQUATIONwhich can be used for inference on new trajectories.Finally, we consider uncertainty prediction. till now we have defined conditional probability paths using a fixed noise parameter σ. However, this does not have to be fixed. Instead, we consider a learned σ θ (t, x t ) which can be learned iteratively with the loss:L uncertainty (θ, x) = t∈T σ θ (t, x t ) -∥x θ (t, x t ) -x next ∥ 2 2 2 2(13) which learns to predict the error in the estimate of x t . This loss can be interpreted as training an epistemic uncertainty predictor which is similar to that proposed in direct epistemic uncertainty prediction (DEUP) [Lahlou et al., 2023] .  
experimental results: 
 In this section we empirically evaluate the performance of the trajectory flow matching objective in terms of time series modeling error, but also uncertainty quantification. We also evaluate a variety of simulation-based and simulation-free methods including both stochastic and deterministic methods. Stochastic methods are in general more difficult to fit, but can be used to better model uncertainty and variance. Further experimental details can be found in §B. Experiments were run on a computing cluster with a heterogenous cluster of NVIDIA RTX8000, V100, A40, and A100 GPUs for approximately 24,000 GPU hours. Individual training runs require approximately one gpu day.Baselines In addition to different ablations of trajectory flow matching, we also evaluate Neu-ralODE [Chen et al., 2018] , NeuralSDE [Li et al., 2020 , Kidger et al., 2021b , Kidger, 2022] , Latent NeuralODE [Rubanova et al., 2019] , and an aligned flow matching method (Aligned FM) [Liu et al., 2023a , Somnath et al., 2023] where the couplings are sampled according to the ground truth coupling during training.Metrics We primarily make use of two metrics. The average mean-squared-error (Mean MSE) over left out time series to measure the time series modeling error defined asEQUATIONwhere x and x are the predicted and true trajectories respectively. We also use the maximum mean discrepancy with a radial basis function kernel (RBF MMD) which measures how well the distribution over next observation is modelled by comparing the predicted distribution to the distribution over next states in the ground truth trajectory. Specifically we compute:EQUATIONwhere∆t = xt -x t-1 , ∆ t = x t -x t-1 , and xt := t s=t-1 f θ (s, x s )ds + g θ (t,x s )dW s is a set of samples from the model prediction at time t.We begin by evaluating how trajectory flow matching performs in a simple one dimensional setting of harmonic oscillators. We show that vanilla conditional flow and bridge matching [Liu et al., 2023c,b, Albergo and Vanden-Eijnden, 2023] , specifically aligned approaches [Somnath et al., 2023 , Liu et al., 2023a] are unable to preserve the coupling even in a simple one dimensional setting. However, augmented with our trajectory flow matching approach, and specifically using (A2), which includes information on previous observations, the model is able to fit the harmonic oscillator dataset well. The harmonic oscillator dataset consists of one-dimensional oscillatory trajectories from a damped harmonic oscillator, with each trajectory distinguished by a unique damping coefficient c. Specifically we sample trajectories x from:EQUATION)where v is the velocity of the oscillator updated byv i = v i-1 + - c m v i-1 - k m x i-1 (t i -t i-1 ); v 0 = 0 (17)with t i = 0.1 • i for i = 0, 1, 2, . . . , 99, spring constant k = 1, and mass m = 1. As c increases, the trajectories evolve from underdamped scenarios with prolonged oscillations to critically and overdamped states where the oscillator quickly stabilizes. This leads to intersecting trajectories due to frequency and phase differences, despite their shared starting point. We perform overfitting experiments on three trajectories generated by varying c. As shown in Figure 2 , models without history information are unable to distinguish between the three crossing trajectories that share the same starting point, resulting in overlapping predictions. In contrast, TFM-ODE that incorporates three previous observations is able to fit the crossing trajectories with high accuracy, with the predicted trajectories almost completely overlapping the ground truth. This is because the dataset with satisfies (A2) with h = 4 (TFM-ODE), but not h = 0 (TFM-ODE no memory and Aligned FM).Next we compared the performance of TFM and TFM-ODE with the current SDE and ODE baselines, respectively, for modeling real-world patient trajectories formed with heart rate and mean arterial blood pressure measurements within the first 24 hours of admission across three different datasets. These are clinical measurements that are taken most frequently and used to evaluate the hemodynamic status of patients, a key indicator of disease severity. Additionally, we evaluated our models against flow matching on these datasets, each with distinct characteristics, to assess their ability to generalize across different distributions. A full description of the datasets are available in Appendix B.2 with the publicly available datasets used under The PhysioNet Credentialed Health Data License Version 1.5.0 and the EHR dataset with local institutional IRB approval:• ICU Sepsis: a subset of the eICU Collaborative Research Database v2.0 [Pollard et al., 2019] of patients admitted with sepsis as the primary diagnosis • ICU Cardiac Arrest: a subset of the eICU Collaborative Research Database v2.0 [Pollard et al., 2019] of patients at risk for cardiac arrest • ICU GIB: a subset of the Medical Information Mart for Intensive Care III [Johnson et al., 2016] of patients with gastrointestinal bleeding as the primary diagnosis • ED GIB: patients presenting with signs and symptoms of acute gastrointestinal bleeding to the emergency department of a large tertiary care academic health system 4.2.1 Prediction accuracy and precision: TFM and TFM-ODE TFM-ODE yields more accurate trajectory prediction Across the three datasets TFM-ODE outperformed the baseline models by 15% to 20%, as seen in table 1 . We noticed that TFM has a similar performance as TFM-ODE. In one case TFM outperformed the non-stochastic TFM-ODE, as seen in the ICU GIB dataset. For ICU sepsis, the performance improvement from the baseline is the most significant, around 83%. This coincides with the ICU sepsis dataset having the most amount of measurement per trajectory. The improvement is seen in both TFM and TFM-ODE, possibly indicating they are able to learn better given more data, resulting in a more precise flow. Not formally measured, we noted that given the same time constraint, FM based models were significantly faster and often finished training before the time limit.TFM yields better uncertainty prediction Though TFM-ODE had lower test MSE for two out of three times, TFM yielded better uncertainty prediction overall, as seen in table 2 . Notably, TFM also had less variance in the uncertainty prediction than TFM-ODE. A plausible explanation in this case is a sacrifice in bias that subsequently decreases the variance for the stochastic implementation, reflecting the bias-variance trade off. Sampled graphs of TFM can be seen in figure 3 . It is notable that the model is able to detect the measurement uncertainty at certain timepoints, matching the increase in amplitude of oscillation in patient trajectories.TFM trajectories accurately match the noise distribution in the data TFM is able to match the noise distribution in addition to the overall trajectory shape, which is useful in settings where data has high stochasticity. We compared our models to NeuralODE and NeuralSDE in matching the variance in neighboring data points, seen in table 3 . We verify that between the baseline NeuralSDE and NeuralODE, NeuralSDE has a lower MMD and is better able to match data points. We find in ICU GIB and ED GIB datasets, TFM outperforms both in matching the variance in data. Notably, the performance pattern is reversed for the MMD metrics and mean MSE metrics with respect to TFM and TFM-ODE where better MSE leads to worse MMD and vice versa. As such, this further confirms the bias-variance trade-off for both TFM and TFM-ODE implementation. We performed ablation studies on TFM and TFM-ODE to attribute importance of various model components contributing to the performance, as seen in table 4. We examined three aspects of the model, two of which were part of our main contributions: uncertainty prediction and memory. We also ablate the model hidden dimension width to infer its potential in scaling effect.TFM and TFM-ODE performance scales with model size In contrast to Neural DE based models, TFM and TFM-ODE exhibit scaling effect, in which the model performance becomes better with a larger hidden dimension. This has been observed in previous flow matching models in image generation [Tong et al., 2024] . This may pave the way for further improvements from larger models. Trajectory memory may improve performance in high frequency measurement settings We conditioned the model based on a sliding window of trajectory history to disentangle data points that otherwise look indistinguishable to FM models. This improved the interpolation performance in the ICU Sepsis and ICU GIB dataset. Notably, this modification did not improve performance the ED GIB dataset, which could be due to shorter trajectories for patients and lower measurement frequency in the defined time period. This may also be explained by the decreased severity of disease in the ED compared to the ICU. Adding memory as a condition may be more suitable for patients whose clinical trajectories have a higher frequency of measurements. 
related work: 
 Continuous-time neural network architectures have outperformed traditional RNN methods in modeling irregularly sampled clinical time series to optimize interpolation and extrapolation. Neural ODE with latent representations of trajectories [Rubanova et al., 2019] outperformed RNN-based approaches [Lipton et al., 2016 , Che et al., 2018 , Cao et al., 2018 , Rajkomar et al., 2018] for interpolation while providing explicit uncertainty estimates about latent states. More recently, Neural SDEs appear to outperform LSTM [Hochreiter and Schmidhuber, 1997] , Neural ODE [Chen et al., 2018 , De Brouwer et al., 2019 , Dupont et al., 2019 , Lechner and Hasani, 2020] , and attention-based [Shukla and Marlin, 2021, Lee et al., 2022] approaches in interpolation performance while natively handling uncertainty using drift and diffusion terms [Oh et al., 2024] . Discrete-time approaches offer an alternative to our continuous-time model model transformers utilize a discrete-time representation with a sequential processing [Gao et al., 2024 , Nie et al., 2023 , Woo et al., 2024 , Ansari et al., 2024 , Dong et al., 2024 , Garza and Mergenthaler-Canseco, 2023 , Das et al., 2024 , Liu et al., 2024 , Kuvshinova et al., 2024 ] models for traditional time series modeling. Adaptations to the baseline transformer includes structuring observations into text with finetuning [Zhang et al., 2023 , Zhou et al., 2023] , without finetuning [Xue and Salim, 2024, Gruver et al., 2023] , or using autoregressive model vision transformers to model unevenly spaced time series data by converting time series into images [Li et al., 2023] . Continuous-time systems are of great interest for learning causal representations using assumptions by using observations to directly modify the system state [De Brouwer et al., 2022, Jia and Benson, 2019] . Variations include intervention modeling with separate ODEs for interventions and outcome processes [Gwak et al., 2020] , using liquid time-constant networks [Hasani et al., 2021 , Vorbach et al., 2021] , or modeling treatment effects with either one [Bellot and van der Schaar, 2021] or multiple interventions [Seedat et al., 2022] Question: Does the paper provide open access to the data and code, with sufficient instructions to faithfully reproduce the main experimental results, as described in supplemental material? Answer: [Yes] Justification: Anonymised code is available as supplementary material. The ICU Sepsis dataset is available publicly online via the eICU database and the ICU GIB dataset is available via the MIMIC-III database, whereas the ED GIB datasets are datasets which include identifiable information. A deidentified dataset may be made available upon reasonable request subject to guidelines set by the institutional review board and in accordance with HIPAA policy. Guidelines:• The answer NA means that paper does not include experiments requiring code. • Providing as much information as possible in supplemental material (appended to the paper) is recommended, but including URLs to data and code is permitted. 
conclusion: 
 In this work we present Trajectory Flow Matching, a simulation-free training algorithm for neural differential equation models. We show when trajectory flow matching is valid theoretically, then demonstrate its usefulness empirically in a clinical setting. The ability to model the underlying continuous physiologic processes during critical illness using irregular, sparsely sampled, and noisy data has the potential for broad impacts in care settings such as the emergency department or ICU. These models could be used to improve clinical decision making, inform monitoring strategies, and optimize resource allocation by identifying which patients are likely to deteriorate or recover. These use cases will require thorough prospective validation and calibration for specific clinical outcomes, for example using the likelihood of a patient crossing a specific heart rate or blood pressure threshold for decisions on level of care (ICU versus inpatient floors) or specific interventions such as transfusions. In these applications, it will be important to assess and control for bias that may be present due to which patient subpopulations are present in training data.Limitations Limitations of the method includes the selective utility of integrating memory in clinical settings with high measurement frequency and no current capacity for estimating causal representations, though this will be an important future research direction. Potential harms include the following: erroneous predictions that either results in delayed care or overutilization of the health system. Accurate trajectory predictions have the potential to inform clinical decision-making regarding the appropriate level of care, leading to more timely and appropriate interventions.Question: Does the paper specify all the training and test details (e.g., data splits, hyperparameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?Answer: [Yes] Justification: We provide an overview with the significance of each clinical dataset in the body of the paper. In the Appendix we have written detailed descriptions of each dataset, including the methodology for identifying the cohorts of patients, input variables extracted, number of patient encounters, summary statistics, as well as frequency distributions. We also detail the process of training and testing for all methods presented in the paper, from baseline to our novel method. We have also provided anonymized link to all code used for the experiments.Guidelines:• The answer NA means that the paper does not include experiments.• The experimental setting should be presented in the core of the paper to a level of detail that is necessary to appreciate the results and make sense of them. • The full details can be provided either with the code, in appendix, or as supplemental material. 
broader impact: 
 Our work extends flow matching into the domain of time series modeling, demonstrating a specific instance of clinical time series prediction. In contrast to the large transformer-based models, our method has fewer in parameters and less training time needed. Notably, it scales well with parameters.As well, our parameterization on Stochastic Differential Equations (SDE) allow faster training time than traditional SDE integration. Accurate timeseries modeling in healthcare has the potential for significant benefits, but also introduces risks. Benefits that could be derived from more accurate prediction of clinical courses include improved treatment decisions, resource allocation, as well as more informative discussions of prognosis with patients or family members. Risks may come from inaccuracies in predictions which could lead to harms by biasing decision making of clinical teams. In the general case of false negative prediction (prediction of trajectories with falsely favorable outcomes) this may lead to undertreatment and in the case of false positive prediction (prediction of trajectories with incorrect detrimental outcomes) or overtreating patients. These inaccuracies may also propagate biases in training data.To move towards broad impact in the clinical domain, this work will require validation and bias estimates. Furthermore, models deployed in domains with high-stakes prediction require interpretability, which can help identify biases, miscalibration, discordance with domain knowledge, as well as build trust with teams using predictions from the model. At this time, flow-based methods have limited tools for interpretability, and we recognize this as a gap in need of future work.Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments? Answer: [Yes] Justification: The paper includes a standard error of the mean for all the results reported in a numerical format. There are no figures with error bars for which the nature of error bar calculation would be relevant. There are no statistical tests that were perform which would make discussion of the validity of statistical tests relevant.Guidelines:• The answer NA means that the paper does not include experiments.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. Justification: In the conclusion section we discuss how improved modeling of Emergency Department and ICU physiology can lead to the positive societal impact of improving disease risk prediction which can enable medical teams to make better clinical decisions, provide patients and family members with increased information on the likely course of the illness, and improve resource allocation by identifying which patients may require costly or rare resources such a blood transfusions. Guidelines:• The answer NA means that there is no societal impact of the work performed.• If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact. • Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), fairness considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations. • The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is a direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that an improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that a generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster. • The authors should consider possible harms that could arise when the technology is being used as intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional or unintentional) misuse of the technology.• If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms to monitor how a system learns from feedback over time, improving the efficiency and accessibility of ML). 
safeguards: 
 Question: Does the paper describe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?Answer: [NA] Justification: The current study proposes a predictive model for time-varying datasets and as such does not pose risks that are significantly higher in comparison to conventional models that fit data trends. There is no generative component in our model, and little risk that the results could be misused in a way that misleads or proves otherwise detrimental to the broader public. We therefore believe this item does not apply for this submission.Guidelines:• The answer NA means that the paper poses no such risks.• Released models that have a high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters. • Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images. • We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make a best faith effort.12. Licenses for existing assets Question: Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?Answer: [Yes] Justification: The publicly available datasets from MIMIC-III and eICU databases are properly credited, respected, mentioned and used under the PhysioNet Credentialed Health Data License Version 1.5.0. The ED GIB EHR data is owned by the authors.Guidelines:• The answer NA means that the paper does not use existing assets.• The authors should cite the original paper that produced the code package or dataset.• The authors should state which version of the asset is used and, if possible, include a URL. • The name of the license (e.g., CC-BY 4.0) should be included for each asset.• For scraped data from a particular source (e.g., website), the copyright and terms of service of that source should be provided. • If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswithcode.com/datasets has curated licenses for some datasets. Their licensing guide can help determine the license of a dataset. • For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided. • If this information is not available online, the authors are encouraged to reach out to the asset's creators. 
new assets: 
 Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?Answer: [No] Justification: We will not be releasing new datasets. Guidelines:• The answer NA means that the paper does not release new assets.• Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc. • The paper should discuss whether and how consent was obtained from people whose asset is used. • We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurIPS Code of Ethics and the guidelines for their institution. • For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review. 
 
 Appendix: We hope to extend our method to cover other types of time series that have periodicity in the components, potentially incorporating Fourier transform [Li et al., 2021] and Physics-Inspired Neural Networks (PINN). Since interpretability is an important factor for clinical reliability, we are developing methods to further elucidate key components affecting the prediction. As well, we hope to incorporate functional flow matching for fully continuous setting [Kerrigan et al., 2024] .We first prove a Lemma which shows TFM learns valid flows between distributions with the target prediction reparameterization trick. Lemma A.1. If p t (x) > 0, δ data is Lipschitz continuous for all x ∈ R d and t ∈ [0, 1], L F M and L T F M are equal,∇ θ L F M (θ) = ∇ θ L T F M (θ)Proof. This proof is a simple extension of Lipman et al. [2023] , Tong et al. [2024] which proved L CF M and L F M are equal under similar constraint.Given δ data = t 1 -t 0 , we have u t (x) = x1-x0 δ data where t 0 is the previous time in the time series, and t 1 is the current time for inference. For the time series data, we are assuming it to be Lipschitz continuous there exist L ≥ 0 such that for all x, y ∈ R n , |f (x) -f (y)| ≤ L∥x -y∥.EQUATIONBy bilinearity of the 2-norm and since x 1 is independent of θ. Next,E pt(x) 1 (1 -t) 2 ∥x 1 θ (t, x)∥ 2 = ∥x 1 θ (t, x)∥ 2 p t (x)dx = ∥x 1 θ (t, x)∥ 2 p t (x|z)q(z)dzdx = E q(z),pt(x|z) ∥x 1 θ (t, x)∥ 2 Finally, E pt(x) x1 θ (t, x), x 1 = x1 θ (t, x), x 1 p t (x|z)q(z)dz p t (x) p t (x)dx = x1 θ (t, x), x 1 p t (x|z)q(z)dz dx = x1 θ (t, x), x 1 p t (x|z)q(z)dzdx = E q(z),pt(x|z) x1 θ (t, x), x 1Where we first substitute then change the order of integration for the final equality. Since at all times t the gradients of L FM andL T F M are equal, ∇ θ L FM (θ) = ∇ θ L T F M by substitution. LE t,q(z),pt(x|z) ∥v θ (t, x) -u t (x|z)∥ 2 = E t,q(z),pt(x|z) 1 (⌈t⌉ -t) 2 ∥x ⌈t⌉ θ (t, x) -x ⌈t⌉ ∥ 2 (21) E t,q(z),pt(x|z) ∥v θ (t, x) -u t (x|z)∥ 2 = E t,q(z),pt(x|z) x⌈t⌉ θ (t, x) -x ⌈t⌉ -t - x ⌈t⌉ -x ⌈t⌉ -t 2 (22) = E t,q(z),pt(x|z) 1 (⌈t⌉ -t) 2 ∥x ⌈t⌉ θ (t, x) -x ⌈t⌉ ∥ 2 (23)Lemma 3.1 The SDE dx t = u t (x|z)dt + σ 2 dW t where u t is defined in eq. 9 generates p t (x|z) in eq. 8 with initial condition p 0 := δ x1 where δ is the Dirac delta function.Proof. For simplicity of notation we first show the case where ⌈t⌉ = 1.dx t = u t (x|z)dt + σ 2 dW t = 1 -x t 1 -t dt + σ 2 dW t (24)which is equivalent to the d dimensional Brownian bridge which has marginalEQUATIONwhich completes the proof for ⌈t⌉ = 1.Proposition 3.2 (Coupling Preservation) Under mild regulatory criteria on u t (•|z), p t , and q, if E t∼U (0,T ),z∼q(z),c∼q(c|z),xt∼pt(xt|z) ∥u t (x t |z, c) -u t (x t |c)∥ 2 2 = 0 and z, q(z), p t (x|z), and u t (x|z) are as defined in eqs. 6-9 then Π(u) ⋆ = Π ⋆ (x 1:T ).Proof. We prove the deterministic case with T = 1. The extensions to stochastic and T > 1 are evident. The couplings are equal if the marginal vector field u t (x t |c) = u t (x t |z, c) everywhere as the coupling is governed by the push forward flows ϕ(x 0 , c) = 1 0 u t (x t |c)dt, and ϕ(x 0 , c, z)= 1 0 (u t (x t |z, c). If E t∼U (0,T ),z∼q(z),c∼q(c|z),xt∼pt(xt|z) ∥u t (x t |z, c) -u t (x t |c)∥ 2 2 = 0 then ϕ(x 0 , c, z) = ϕ(x 0 , c) for all x 0 and therefore the couplings of the optimal map are equivalent. We note that this requires exchange of integrals under the same conditions as Lemma A.1.Next we show how (A1)-(A3) satisfy Prop. 3.2.(A1) When c = x 0 and there exists T : X → X such that T (x 0 ) = x 1 if and only if Π ⋆ (x 0 , x 1 ).We note that this is equivalent to asserting the existence of a Monge map T ⋆ for the coupling Π ⋆ .In the two timepoint case, c = x 0 is sufficient as long as there aren't two trajectories that have the same x 0 but different x 1 s. Conditioning on this way ensures the conditions of of Prop. 3.2 as the uniqueness property ensures the uniqueness of u t (x t |c).(A2) There exist no two trajectories x i , x j such that x i t = x j t for h + 1 consecutive observations. In this case notice that this is simply a multi-timepoint extension of A1 to c = x t-h-1:t-1 , i.e. conditioned on a history of length h. If this is the case then the same reasoning as A1 applies. (A3) Trajectories are associated with unique conditional vectors c independent of t.This satisfies Prop 3.2 by definition.Proposition 3.3 There exists a scaling function c(t) : R + → R such that L target (θ) = c(t)L match (θ).Proof. We start with the matching loss.E t,q(z),pt(x|z) ∥v θ (t, x) -u t (x|z)∥ 2 = E t,q(z),pt(x|z) 1 (⌈t⌉ -t) 2 ∥x ⌈t⌉ θ (t, x) -x ⌈t⌉ ∥ 2 (26)by substitution,E t,q(z),pt(x|z) ∥v θ (t, x) -u t (x|z)∥ 2 = E t,q(z),pt(x|z) x⌈t⌉ θ (t, x) -x ⌈t⌉ -t - x ⌈t⌉ -x ⌈t⌉ -t 2 (27) = E t,q(z),pt(x|z) 1 (⌈t⌉ -t) 2 ∥x ⌈t⌉ θ (t, x) -x ⌈t⌉ ∥ 2 (28)completing the proof. The three oscillation trajectories correspond to c = 0.25 (the red trajectory in Figure 2 ), c = 2 (blue), and c = 3.75 (green). Before used as an input, t was scaled to between 0 and 1 by dividing by 10.In order to accurately model the perturbations in the physiologic signals (mean arterial pressure and heart rate) of the underlying patient states, we need to learn beyond the general trend of the data. While the physiologic measurements themselves reflect patient status and drive clinical decision making, the degree of variation holds information that goes beyond the snapshot at a single time point. Our approach models the data distribution and stochasticity rather than just fitting the average trajectory. Other time-varying such as treatment conditions and non-time-varying covariates such as underlying disease states may also hold information that may impact the underlying state generating the physiologic signals. Our approach also incorporates this information to inform the trajectory modeling.The data distribution in the ICU datasets reflect its status as the most resource-intensive clinical setting with increased measurement frequency and data distribution shift towards more abnormal physiologic values (Figure 4 ). The ED dataset reflects its status as the clinical setting focused on triaging patients, with sparser and physiologic measurements that fall in the normal range.For each clinical dataset, we modeled patient trajectories formed with heart rate and blood pressure measurements during the first 24 hours following admission. The timeline for each trajectory, originally in minutes, was scaled to a range between 0 and 1 by dividing by 1440. Additionally, heart rate and blood pressure values were z-score normalized to standardize the data.Intensive Care Unit Sepsis (ICU Sepsis) Dataset The eICU Collaborative Research Database v2.0 [Pollard et al., 2019] Norepinephrine infusion rates were calculated by converting drug doses or infusion rates to µg/kg/min, and where drug doses were not explicitly available, the dose was inferred from the free text given in the drug name. Start and end times for norepinephrine infusion were calculated by dividing the dose by the infusion rate. Where there appeared to be multiple infusions at the same time, the maximum infusion rate was taken as the infusion rate. As a conditional input to the models, the norepinephrine infusion doses are then scaled to between 0 and 1 by dividing by the maximum norepinephrine value in the dataset.The APACHE-IV score, a validated critical care risk score, predicts individual patient mortality risk [Zimmerman et al., 2006] . In data preprocessing, we uses logistic regression of the score against binary hospital mortality data to generate a probability for each patient, serving as an additional input condition for models.Intensive Care Unit Cardiac Arrest (ICU Cardiac Arrest) Dataset This dataset was extracted from the to eICU Collaborative Research Database v2.0 [Pollard et al., 2019] Age, initial lab test measurements (three labs omitted due to missing data), and pre-existing medical conditions were used to train an XGBoost model [Chen and Guestrin, 2016] to predict the binary outcome variable indicating the need for hospital-based care. The resulting probabilities of requiring hospital-based care (outcome of 1) for each patient were then calculated using the trained model and used as conditional input to conditional models in experiments on this dataset. Of note, the outcome variable was defined as 1 if a patient (1) requires red blood cell transfusion, (2) requires urgent intervention (endoscopic, interventional radiologic, or surgical) to stop bleeding or (3) all-cause 30-day mortality. Labs and medical conditions included in this dataset are listed below. Labs in bold were excluded from the XGBoost risk score calculation due to missing data. Since the trajectories in this dataset are deterministic and regularly sampled, we deployed only TFM-ODE and applied solely the L match loss (i.e, no uncertainty or time predictive loss), as these methods sufficiently address the structured nature of the data to generate proof-of-concept results.The three models presented in Figure 2 all have hidden size of 256, σ of 0.1, trained under seed=0 with Adam optimizer with learning rate 1 × 10 -3 for a maximum of 1000 epochs with early stopping (patience=3) validation loss.All the models for clinical data experiments are trained with Adam optimizer. A maximum training time and epochs are set to 48 hours and 300, with early stopping (patience=3) monitoring validation loss. All metrics reported were ran with 5 seeds (0,1,2,3,4) to ensure it is reproducible.TFM, TFM-ODE, and ablations The TFM models were trained with learning rate 1 × 10 -6 and had σ of 0.1. The complete models have hidden size of 256 and memory of 3, while ablation study with a hidden size of 64 and/or no memory was performed (Table 4 ). The noise parameter for the SDE implementation was set to 0.1 for ablations without L uncertainty . The hyperparameters σ = 0.1 and memory=3 for full models were selected through experiments with different values of σ and memory (Figure 5 and 6 ).FM The FM baseline models were trained with learning rate 1 × 10 -6 . All models had a hidden size of 64 with σ of 0.1.Latent Neural ODE The latent Neural ODE models were trained with a learning rate of 1 × 10 -3 . 100 GRU units were used for the encoder model and the number of latent dimensions was 2.Baseline Neural SDE and Neural ODE Both baseline models were trained with learning rate 1 × 10 -5 and had a hidden size of 64.• If the contribution is a dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable. • Depending on the contribution, reproducibility can be accomplished in various ways.For example, if the contribution is a novel architecture, describing the architecture fully might suffice, or if the contribution is a specific model and empirical evaluation, it may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed. • While NeurIPS does not require releasing code, the conference does require all submissions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. (d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome to describe the particular way they provide for reproducibility.In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers to have some path to reproducing or verifying the results. 


Decision: The paper presents a significant advancement by bridging flow matching and differential equation-based models, a novel and rigorous approach that is well-supported by theoretical proofs and extensive experiments across multiple real-world datasets. The method’s ability to provide uncertainty prediction and its potential to address the scalability issues of stochastic differential equations (SDEs) while delivering state-of-the-art performance highlight its originality and practical significance. The clarity of the paper, with its logically flowing ideas, and its application to a domain with notable societal benefits further underscore its strengths.

All reviewers are highly supportive of this paper. It stands out by combining substantial algorithmic and theoretical contributions with practical significance.

Reviews: review: paper_summary: paper_summary: This paper presents Trajectory Flow Matching (TFM), a simulation-free training algorithm for neural differential equation models. This enables modeling of continuous physiologic processes using irregular, sparsely sampled, and noisy data, all with better scalability.

The authors provide theoretical proof that matching techniques can allow for simulation-free training, and empirical proof of the effectiveness of TFM in clinical settings. An ablation study is also provided to validate components of TFM.
strength_weakness: strength_weakness: - The motivation for the paper is that SDEs are not scalable. While this is a known fact, the paper lacks a proper discussion of the impact of simulation-free training on the scalability of TFM.
- The provided code lacks instructions to reproduce the results (empty README). 
- Explanability studies would have been welcome to help clinical usability.
questions: questions: - Did you compare TFM-ODE with other models that are able to process irregular clinical time series?
- What is your reasoning behind modeling specifically the heart rate + MAP pair ?
limitations: limitations: The authors have adequately addressed the technical limitations (memory not always suited depending on the dataset, no causal representations estimations) as well as societal impacts (false predictions).
paper_summary: paper_summary: This paper presents Trajectory Flow Matching, an extension of flow matching to time series. It can model irregular, sparsely sampled, and noisy time series. It trains in a simulation-free manner, bypassing backpropagation through the dynamics. The method is tested on ICU physiological time series, demonstration SoTA performance, and uncertainty prediction.
strength_weakness: strength_weakness: None of particular note.
questions: questions: 1) It would be interesting to understand if this could be applied to periodic time series like ECG. Do you see ways to introduce this behaviour, and if so, could you call this out as a future research direction?
limitations: limitations: The authors clearly describe the limitations.
paper_summary: paper_summary: This paper proposes a trajectory flow matching method for time series data, by using flow matching at each time points. To preserve the coupling of time series, the vector fields are conditioned on history lengths (or even more general $c$). The method also provides ways for model stability, irregularly sampled trajectories and uncertainty prediction. After evaluating the model performance on simple harmonic oscillators, they further apply their method to three clinical time series datasets.
strength_weakness: strength_weakness: The TFM is the model for discrete time series. It would potentially be more useful to model things in a continuous way.
questions: questions: The paper is clearly written. No major questions.
limitations: limitations: The same as weakness

Comments: comment: Thanks for the authors for their response and I would keep my initial rating.
Thank you for your response, particularly expanding future work and adding the README.md detailed in another review.
Thank you for your response and for clarifying the methods you had in mind. At this time, we have not compared TFM to these specific methods, and we will look into adding comparisons for the revised version of our paper.
Thank you for your response and actions to improve your paper even further. I appreciate your clarification about the choice of clinical variables.


When I mentioned "other models that are able to process irregular clinical time series", I should have been clearer that I had in mind non-differential equations-based models such as STraTS [1],  RAINDROP [2], or Warpformer [3] for example, or even older models such as SAnD [4] or InterpNet [5]. 


Whether you compared your model to these or not, your work still impacts the ODE domain, and my review thus remains positive (Accept), as you have addressed my other concerns.

---

[1] S. Tipirneni and C. K. Reddy, “Self-Supervised Transformer for Sparse and Irregularly Sampled Multivariate Clinical Time-Series,” ACM Trans. Knowl. Discov. Data, vol. 16, no. 6, p. 105:1-105:17, Jul. 2022, doi: 10.1145/3516367.

[2] X. Zhang, M. Zeman, T. Tsiligkaridis, and M. Zitnik, “Graph-Guided Network for Irregularly Sampled Multivariate Time Series,” presented at the International Conference on Learning Representations, Oct. 2021.

 [3] J. Zhang, S. Zheng, W. Cao, J. Bian, and J. Li, “Warpformer: A Multi-scale Modeling Approach for Irregular Clinical Time Series,” in Proceedings of the 29th ACM SIGKDD Conference on Knowledge Discovery and Data Mining, Aug. 2023, pp. 3273–3285. doi: 10.1145/3580305.3599543.

[4] H. Song, D. Rajan, J. Thiagarajan, and A. Spanias, “Attend and Diagnose: Clinical Time Series Analysis Using Attention Models,” Proceedings of the AAAI Conference on Artificial Intelligence, vol. 32, no. 1, Art. no. 1, Apr. 2018, doi: 10.1609/aaai.v32i1.11635.

[5] S. N. Shukla and B. M. Marlin, “Interpolation-Prediction Networks for Irregularly Sampled Time Series,” ArXiv, Sep. 2019.