Reliable Probabilistic Forecasting of Irregular Time Series through Marginalization-Consistent Flows
Authors: Vijaya Krishna Yalavarthi
Randolf Scholz
Christian Klötergens
Kiran Madhusudhanan
Stefan Born
Lars Schmidt-Thieme
Venue: ICLR 2026 Poster (2026)
Keywords: Irregular Time Series
Probabilistic Forecasting
Normalizing Flows
Primary area: probabilistic methods (Bayesian methods, variational inference, sampling, UQ, etc.)
TL;DR: Propose MOSES—a mixture of separable flows over Gaussian processes—that guarantees marginalization consistency while achieving strong predictive performance for probabilistic irregular time series forecasting.

Abstract: Probabilistic forecasting of joint distributions for irregular time series with missing values is an underexplored area in machine learning. Existing models, such as Gaussian Process Regression and ProFITi, are limited: while ProFITi is highly expressive due to its use of normalizing flows, it often produces unrealistic predictions because it lacks marginalization consistency—marginal distributions of subsets of variables may not match those predicted directly, leading to inaccurate marginal forecasts when trained on joints.
We propose MOSES (Mixtures of Separable Flows), a novel model that parametrizes a stochastic process via a mixture of normalizing flows, where each component combines a latent multivariate Gaussian with separable univariate transformations. This design allows MOSES to be analytically marginalized, enabling accurate and reliable predictions for various probabilistic queries.
Thanks to its inherent marginalization consistency, MOSES significantly outperforms all baselines—including ProFITi—on marginal predictions.
For joint predictions, it beats all other consistent models and performs close to or slightly worse than ProFITi. Implementation details:~\url{https://github.com/yalavarthivk/separable_flows}

Published as a conference paper at ICLR 2026
RELIABLE PROBABILISTIC FORECASTING OF IRREGULAR TIME SERIES
THROUGH MARGINALIZATIONN-ONIISSEEEE FLOWS
Vijaya Krishna Yalavarthi' 1.2, Randolf Scholz"
Christian Klötergens 1.2, Kiran Madhusudhanan 1.2 Lars Schmidt-Thieme 1.2
'Information Systems and Machine Learning Lab (ISMLL), University of Hildesheim, Germany 2VWFS Data Analysis Research Center (DARC), University of Hildesheim, Germany (yalavarthi, rscholz, kloetergens, madbusudhanar, schmidt-than ttimeme @ismll1 de
Stefan Born
Institute of Mathematics, TU Berlin, Germany born@math. tu-berlin.ce
ABSTRACT
Prohabilistic forecasting of joint distributions for irregular time series with missing values an underexplored area in machine learning. Existing models, such as Gaussian Process Regression and ProFITi, are limited: while ProFTTi is highly expressive due to its use normalizing flows, it often produces unrealistic pre- dictions because it lacks marginalization consistency-margae subsets of variables may not match those predicted directly, leading to inaccu- rate marginal forecasts when trained on joints. We propose MOSES (Mixtures of Separable Flows), a novel model that parametrizes a stochastic process via a mixture of normalizing flows, where each component combines latent multi- variate Gaussian with separable univariate transformations. This design allows MOSES to be analytically marginalized,, enabling accurate and reliable predie- tions for various probabilistic queries. Thanks to its inherent marginalization consistency, MOSES significantly outperforms all baselines innluudig ProF- ITi on marginal predictions. For joint predictions, it beats all other consistent models and performs close to or slightly worse than ProFITi. Implementation details:
distributions of
i wttpaattww///paaaalblll fows
INTRODUCTION
Probabilistic forecasting of irregular time series requires models that make predictions at arbitrary time points and for an arbitrary subset of variables (channels). Unlike regular time series forecasting, where predictions are confined fixed grid, irregular settings require modeling stochastic process as the number of targets can vary.
core requirement of any stochastic process model marginalization consistency: the marginal distribution of a subset of variables must agree whether computed directly or by integrating out joint. This property, a precondition for Kolmogorov's Extension Theorem, ensures that model defines a mathematically coherent stochastic process. Consistency is not only necessary for mathematical rigor but also highly desirable practice: by the data processing inequality, accurate joint predictions together with consistency ensure reliable marginals, providing both correctness and performance guarantees.
Many recent methods prioritize flexibility over consistency. For example, ProFITi Yalavarthi et al. 2025) achieves strong empirical performance using normalizing flows but violates marginalization consistency, producing contradictory marginal and joint predictions (Figure I).
"Equal contribution

---

Published as a conference paper at ICLR 2026
ground
MOSES(1)
ProFITi
GPR
ground truth
MOSES(4)
ProFITi
GPR
(ours)
(ours)
NLL = 3.472
NLL 3.379
NLL 4.390
NLL = 0.363
NLL= = U.678
N 2.149
MT 0.05
MT 1.59
MT 0.05
MT - 0.01
MT 1,33
MT 0.01
p(v2) p(u1) samples
p(y)) samples
p(y2)
ground truth marginal density
predicted marginal density
marginalized joint density
Figure 1: Demonstration of marginal consistency for MOSES (ours), ProFITi Yalavarthi et al. 2025 and Gaussian Process Regression (GPR) (Bonilla et al. 2007) on two toy datasets: blast and circle. ProFITi is inconsistent w.r.t. the marginals of the second variable /2, while MOSES( (D), where is the number of mixture components, consistent with the marginals of both y1 and 3/2- GPR marginalization consistent but predicts incorrect distributions. Marginalization inconsistency (MI, eq. 17) empirically estimates the KL divergence between the predicted marginal and the marginalized joint. See Appendix C.6 more details.
Example 1.1. In an ICU monitoring system, such a model might predict that blood pressure remains stable with 90% probability when queried individually, but only 60% when inferred from the joint distribution of vitals. These contradictions undermine reliability and could mislead clinical decisions.
Due to inconsistency, ProFITi suffers significantly on marginal distributions despite having good joints. On the other hand, existing consistent approaches, such as Gaussian Process Regression (GPR) (Dürichen et al. 2015), avoid this problem but are limited to Gaussian distributions, creating false contradiction between correctness and flexibility.
We address this problem with MOSES (Mixtures of Separable Flows), model based on mixtures of normalizing flows. Each component of MOSES uses multivariate Gaussian, induced by Gaussian process, as the base distribution, combined with univariate invertible transformations to modulate the marginals-sininn
in spirit Gaussian Copula Processes (Wilson & Ghahramani 2010). Since mixtures, Gaussian processes, and univariate transformations are each marginalization consistent, MOSES is consistent by design. At the same time, retains the expressiveness of flow-based methods for the marginals.
Our contributions as follows:
1. We formalize marginalization consistency as fundamental requirement for probabilistic forecasting models and show that violating it renders model incoherent.
2. We propose MOSES which achieves the expressiveness of modern flow approaches while maintaining strict marginalization consistency through a principled architectural design.
3. We show that consistency does not come at high performance cost: MOSES matches inconsistent models on joint prediction tasks while significantly outperforming them on marginal predictions, providing a principled solution to irregular time series forecasting.
PRELIMINARIES
Let Seq(X') denote the space of finite sequences over X. We represent an irregular time series X as a sequence of N triplets (Horn et al. 2020 Yalavarthi et al. (025):
X (L
OBS
oBs)) u
n=1:N Seq(X), X R x x R,
(I)
where 1oBs is the observation time, oBS the channel, and oBs the observed value.
2

---

Published as a conference paper at ICLR 2026
A time series query Q is sequence of K pairs:
Q ((4
(ORY (RRY
Ck () k-1-K E ()q @ 2= (1,.,,, C},
where OORY and
QRY C
specify the future time and queried channel.
forecasting answer is sequence y = (y) (.)....13K. WKK with ya R predicting the value in channel at JORY All query points occur after the observations: ming JORy > m&X, 1
ORY Ck
Requirements. A marginalization consistent probabilistic irregularly sampled time series forecast- ing model must satisfy the following requirements:
R1 Joint Multivariate Prediction. The goal of probabilistic irregular time series forecasting is to model the joint distribution ppy Q.X) over responses y, given query points Q and observed series X, allowing both context length N |X] and query length K Q to vary dynamically:
p: Seq (R x Q) Seq (2) Rz0: (v. Q, X)
. yK QI...., XI.... XN)
For given (Q. X), the mapping (yI....,yK*) plui,..... UKK X) defines probability density on BIQI.
Unlike standard multivariate time series forecasting, where future time points are fixed and typically ignored, irregular time series require conditioning on future time points and channels, as these can vary across instances.
R2 Permutation Invariance. As the time stamp and channel-ID IDD are included each sample, the order of the samples does not matter, and hence any model prediction should be independent of the order of both the query or context:
ply X) p(y" th (T) VI E Sym(IQ1), Sym([X])
R3 Marginalization Consistency/Projectin
Invariance. Predicting the joint density for the sub- query Q-x given by removing the k-th item from Q should yield the same result as marginalizing the k-th variable from the complete query Q.
/ p(y Q1 dy
p(y-k Q)::
(5)
This can be extended any subset Ks S. KI whenever we marginalize out set of variables Ks with [Ks < m, then, w.l.o.g., let Ky (1, -1,-..,7 + 1} {1,...,mm m) U (m + 1). Now, simply marginalize out the first vari- ables, and apply the property second time to marginalize out m + as well (see Appendix A. for details).
K} via induction: Assuming R3 holds
For model satisfying R1
we will only have to marginalize if we try validate the marginaliza- tion consistency. For this validation we added requirement R3 Yalavarthi et al. (2025) discussed R1
R3
and R2 but did consider R3 We argue that irregularly sampled time series are realizations of stochastic process and R3 is a tundamental property of any model that mimics it.
Theorem 2.1. Any model that satisfies RI R3 realizes an RR-vllued stochastic process over the index sel T R (1,..., C).
Proof. This is direct application of Kolmogorov's extension theorem (Oksendal 2003
WHY DO WE NEED MARGINALIZATION CONSISTENCY?
1. Without marginalization consistency, probabilistic models can be unreliable. Two test cases with the same context X but different queries Q may share overlapping targets. An inconsistent model can assign different distributions to these targets, despite identical conditioning information, which is unrealistic. Such models are inapplicable many real-world applications.
By abuse of notation, p(U X) denotes Py, (v X). For Q (f)) ply () denotes Poy Yo x -.... Vt., X). Here "I Q" index selection, not prohabilistic conditioning Q.
3

---

Published as a conference paper at ICLR 2026
()SES(1)
M2SES(2)
MOSES(3)
MOSES(4)
MOSES(5)
MOSES(
Dataset
GMM(1)
GMM(2)
GMM(3)
GMM(4)
GMM(5)
GMM(6)
Figure (Top) Importance of multiple mixture components: M1SES( cannot represent the correct distribution, but MOSES(4) can. (Bottom) Limitation of Gaussian Mixture Models: more components are needed single ellipsoids cannot model curvature. Color represents mixture component; see Appendix or details.
2. Marginalization consistency provides performance guarantees. For consistent model, if the joint prediction over K points is accurate, predictions for any subset of size < K are also guaranteed to be accurate. This follows from the data processing inequality (DPI; Murphy 2022), which ensures marginalizing cannot increase KL divergence:
DEL (plvi...., VK Q1 ................................ QK x) plyi...., QI....,Q...... xx.
.., DEL (ply...... ... -K-1 Q1, ... Qx j(n..................... (.u.....,KK- 11,--., (K-,,-,,K-I,XX
Intuitively, marginalization "averages out" errors, so accuracy cannot worsen on subsets. Consis- tent models thus remain reliable for marginal predictions, as confirmed in Tables [ and 2 where inconsistent model, ProFITí, performs poorly on marginals despite good joint predictions.
CONSTRUCTING MARGINALIZATION CONSISTENT CONDITIONAL DISTRIBUTIONS
Our goal is to build model for the conditional joint distribution p(,,,...,, Qi,.... Q Kx, X), as in (3). Since the model should satisfy R3 it follows that the marginal distribution of yn must only depend on Qs and X.
Separably Parametrized Gaussians. A simple way model a permotation-invariaat con- ditional distribution for variable number of targets is multivariate Normal distribution Nyy (X, Q). Σ(X, Q)), where the mean and covariance are separably parametrized:
μ4 X), Σn.e E(,x,Qe,X)
(6)
with functions ju and Σ as in Gaussian Processes. Such a setup is marginalization consistent by design, since marginalizing Gaussian simply involves selecting the relevant rows and columns of Σ and elements of μ.. However, Gaussian Processes are limited because they can only model joint Gaussian distributions. For more flexible distributions, normalizing flows are common alternative (Rezende & Mohamed 2015).
Separable Normalizing Flows. transformation f R - K to a source distribution pz on IR K The density of the resulting target distribution is given by the change-of-variable formula:
Normalizing flows model distribution by applying an invertible
0s (05"(v:0))
(2) = (y;
det
By
Most existing flows use simple source distributions, typically standard multivariate normal, and model interactions between variables through the transformation (Rezende & Mohamed 2015
4

---

Published as a conference paper at ICLR 2026
Pr Q.X)
ovalustion
hul
Ph
Prity bu)
ha.1
b..K
H 10 hots is MMXM w/ w(h""
pa ha) (z. hu)
R1 C
Fac C
P2,(2 lt]
Figure 3: Illustration of MOSES. D flows (fixed). K variables (variable at inference time). The Encoder (enc) takes the observed series X and query Q as input, and outputs an embedding (depends on both X and Q) and w (depends only on X). ha is used to parametrize both μaa Ea of Pza and the flow transformation fa- Each flow is separable and applies a spline transform Φ each component 2 of (z,(z] ha). yielding yv of FLLow (y ha).
Pd
Papamakarios et al. 2021 Conditional flows for variable number of targets follow the same Idea, designing expressive transformations for vectors of arbitrary size (Liu et al.) 2019 Bilos & Günnemann 2021 Yalavarthi al., 2025). However, these models generally do not guarantee marginalization consistency, and no simple condition on the transform alone can provide it.
We propose different approach: simple, separable transformations combined with richer source distribution, namely a Gaussian Process with full covariance. In this setup, dependencies between variables are captured by the source distribution rather than the transformation, enabling marginalization consistency. This is similar in spirit to copula- based architectures, where dependen- cies between variables can be separated from their univariate marginal distributions.
Lemma 3.1 (separable flow condition). Assume each variable's 's transformation depends only on its own query and the shared context, not on other queries. A conditional flow model over TR.K or Seq(R) is separable if is expressed in the form
12, = (0(21 Q1. QK,
for some univariate function Φ: R x Q x Seq(X)) - R that is invertible in the first argument. Any model that consists of such separable flow transformation, combined with a marginalization consistent model for the source distribution, is itself marginalization consistent. (Proof: Appendix A.2
Conditional Mixtures of Flows. Using Gaussians as base distributions, we can model only linear dependencies between variables. To increase expressiveness, we combine multiple separable flows into a mixture. Even few components can significantly improve the model, reaching performance comparable to simple GMM without flow transformations (Figure 2 Appendix C.6).
Lemma 3.2. Given probabilistic models (Pa)d=I:D that satisfy RI R3
then a mixture model
ply Q. = Σ (x Pa(3) Q.X)
with permutation invariant weight function w: (eq + A where A denotes probability
simplex in D variables A (w R wa 0, Ea Wa 1}, also satisfies R1 R3 (Proof: Ap- pendix A.3)
4 MIXTURES OF SEPARABLE FLOWS (MOSES)
Based on the constructions from the previous section, we build marghnalization-coniisett model for forecasting irregular time series using four components (Figure 3):
Separable encoder: shared encoding of the observations, hous oBs (OB)) and /) vectorized, query-specific encodings baąe == enc (Qk, X; X000RY)
2. Separable Gaussians:
D multivariate Gaussians, (2, (2 Has Ex), each parametrized separately by the corresponding query encoding ha-
5

---

Published as a conference paper at ICLR 2026
3. Separable transformations: D separable normalizing flows fa are parametrized using the encoded queries hą. Each fa is applied on top of the Gaussian base distribution pz., providing a separable normalizing flow p.
4. Mixture of flows: The D flows are combined using mixing weights w :- w(h which only depend on the shared observation encoding h not on the queries.
1. Separable Encoder. To encode both the observations = ((t2us Co (D ) I= and queries Q (tt 0 ))==1:K, we apply positional embedding with learnable parameters (a,, f=1: to the time component (Kazemi et al. 2019).
f=1:F
pos_embed()) = Sayt+ sin(aat by)
if =1 else
(10)
We use one-hot encodings for the channel component. The value component is passed through.
[pos_embed (t, one-hot(c..") 1 1
(11a)
1==1:W
q [pos_emml bos_embed (( one-hot(ck (ou 1=
(11b)
The observations R. (F+C+)) are further encoded via self-attention, and the queries
E (x((++C) via cross-attention w.r.t. the encoded observations:
h MHA(x, X, x; 00RS)
(E (ERMXM)
(12a)
== MHA(q. h OBS h ORS, (ORY)
(E RER.D.
(12b)
reshape(h)
(E (DDDxXXM)
(12c)
where MHA denotes multi-head attention. For the encoding of the queries, we use an encoding dimension D M and reshape each hş into D encodings hd, dimension M.
2. Separable Gaussians pz, (z pa, Ea). We model the means and covariances of the base multivariate Gaussian using simple linear and quadratic functions of the encoded queries ha:
p(ba)= MaOMEAN
u(ba)x k = ha. MMEAN
(13a)
E(ha) T + (haBcov) (hooTory
E(haea = Ens (Cd, Qhoy
(13b)
VM'
VM
where JMEAN Mxx11 and gcov MM'' are trainable weights shared across all D mixture components, 1x is the K x K identity matrix, and 6ạ1 is the Kronecker delta.
The scaling by /M in (13b) ensures stable learning, following (Vaswani et al. 2017). Since E(ha) the sum of a positive semidefinite and positive definite matrix, it remains positive definite. The encoding h, incorporates both the context X and the queries Q, consistent with the separable Gaussian setup in (6). Note that although the base Gaussian defined via Gaussian processes (GPs), we do not perform GP-style inference.
Separable transformations fa. To obtain separable invertible transformations, we apply univariate bijective function to each variable independently. ppline-based functions are particularly popular due to their expressiveness and good generalization (Durkan al., 2019 (Dolatabadi et al. 2020). We use computationally efficient Linear Rational Spline (LRS) transformations (Dolatabadi ct al. 2020). For a conditional LRS Φ(zz ha,s; (FLOM), parameters such as bin widths and heights, derivatives at the knots, and À are computed from the conditioning input ha and shared model parameters gFLow, The same gFLOw is used for all variables 21:K, allowing the transformation to handle varying numbers of variables, and is also shared across all D mixture components. For details, see Appendiix A.5 In summary, the conditional flow model is separable across the query size JK with
fa(2 Q,X) = f(z ba) = (Ø(21,
((zK, ha,k))
(14)
4. Mixture Model. We model the mixture weights via cross-attention, using trainable parameters 3 E DDM as attention queries, and a softmax to ensure the weights sum to 1:
softmax(MHAl) h obs (OBB,8MIX)
(15)
Theorem 4.1. Our model. MOSES, satisfies RI R3 and hence realizes a stochastic process via Kolmogorov's 's Extension Theorem (see Theorem 2. Proof: Appendix A.4
R3

---

Published as a conference paper at ICLR 2026
Computational Complexities. The D separable flows computationally efficient: since they are separable, their Jacobian matrix is diagonal, and computing its determinant requires O(K) operations. The main computational cost lies in evaluating 2 and det 2 for the base distribution, which typically requires ((K*) operations. However, for large K, our low-rank modification Ea = Ix + UU (see (13b)) reduces their computation to (M using the Woodbury and
Weinstein Weinstein-Aronzzajnn identities. This approach scales well for large values of > M', as M' is independent of K.
Training. Given batch B of training instances (Q. X)y), we minimize the normalized joint negative log-likelihood (njNLL) Yallavrthii et al. 2025
CONNLL
log ply x,
(16)
(Q
where / includes all the model parameters. njNLL generalizes NLL to varying target sizes.
Relationship to Copulas. Each component of MOSES is conceptually similar to Gaussian copula process Wilson & Ghahramani (2010), in that separates the linear dependency structure from the marginal distributions. The base Gaussian covariance matrix can be decomposed into correlation matrix and diagonal variances. The diagonal variances, together with the base mean and the univariate spline transformations, determine the contribution of each component to the marginal distribution of each variable. On the other hand, since MOSES is a mixture of D such components, each with its own transformations, the overall model does not have a simple separation between dependency structure and marginal distributions.
Many copula-based architectures train marginals and dependencies separately avoid identifiability issues. However, MOSES is trained end-to-end, and because of the restricted domain of the spline transformations, identifiability issues are largely mitigated. shift in the mean of the base distribution cannot be compensated by the spline transformation.
To the best of our knowledge, copula-based architectures have not been applied to probabilistic forecasting of irregular time series. Models such TACTIS and TACTiS-2 (Drouin et al. 2022 Ashok et al., (202) are designed for regularly sampled, fully observed data. While their principles could extend to irregular time series, their implementations are not readily adaptable, and they lack marginalization consistency due to their non-separable encoder and copula structure.
5 RELATED WORK
There have been multiple works that deal with point forecasting of irregular time series Ansari
jet al. 2023 Che et al.] 2018 Scholz et al. 2023: Chen et al.) 2024, Yalavarthi et al., 2024 Zhang
eet al. 2024 Klötergens
al. 2224 2025a a b) In this work, deal with probabilistic forecasting
of irregular time series. Models such as NeuralFlows (Biloš et al. 20.T). GRU-ODE (De Brouwer
eet al., 2019), and CRU (Schirmer et al., 2022 predict only the marginal distribution for single time stamp. Additionally, interpolation models like HetVAE (Shukla & Marlin, (022) and Triplet- former (Yalavarthi et al. 2023) can also be applied for probabilistic forecasting. However, they also produce only marginal distributions. All the above models assume the underlying distribution is Gaussian which is not the case for many real-world datasets. On the other hand, Gaussian Process Regression (GPR Dürichen et al. 2015). and ProFITi (Yalavarthi et al.) 2025) can predict proper joint distributions. ProFITi is not marginalization consistent because of its separabbe encoder and probabilistic component.
There have been works on models for tractable and consistent marginals for fixed size variables (e.g. tabular data). Probabilistic Circuits (Choi et al. 2020) create sum-prod network on the marginal distributions. Later, Sidheekh et al. (2023) added univariate normalizing flows to the leaf nodes of the circuit for better expressiveness. However, it is not trivial to extend such circuits to deal with sequential data of variable size. Gaussian Mixture Models (GMMs) (Duda & Hari 1974) are often used only for unconditional density estimation, but be extended to conditional density estimation. They can provide tractable and consistent marginal distributions. However, GMMs are not expressive enough and often require very large number of components to approximate even
7

---

Published as a conference paper at ICLR 2026
Table I: Comparing NLL. Lower is better, best results in bold, second best in italics.
Model
USHCN
PhysioNet
MIMIC-III
MIMIC-IV
inconsistent
ProFITi
-3.226 0.225
-0.647 0.078
-0.377 0.032
-1.777 0.066
GRU-ODE
0.766 = 0.159
0.501 0.001
0.961 = 0.064
0.823 0.318
consistent
NeuralFlows
0.152
0.496 = 0.001
0,998 + 0.177
0.689 0.087
univariate
CRU
0,761 0.191
1.057 = 0.007
1.234 + 0.076
OOM
Triplctformer+
4.632 + 8.179
0.519 # 0.112
1,051 # 0.141
0.686 0.115
GPR
2.011 + 1.376
1.367 = 0.074
3.146 + 0.359
2.789 0.057
consistent
GMM
1.050
1,063 = 0.002
1,160 + 0.020
1.076 0.003
multivariate
MOSES (ours)
-3.357 0.176
-0.491 0.041
-0.305 # 0.027
-1.668 = 0.097
simple distributions, see Figure
Note that normalizing flow models such as Dinh et al. (2017) 7):
Paaamaaarrios et (2.1 2021 neither provide tractable marginals nor are applicable to varying
number variables.
Prior work on mixtures of normalizing flows has focused on fixed-length sequences. Pires & Figueiredo (2020) and Ciobanu (202T) used affine coupling and masked autoregressive flows for density estimation, Postels al. (202I) applied them reconstruction tasks. These models face challenges with variable-length sequences and intractable marginals, [Sendera et al. (202T) introduced non-Gauee non-Gaussian Gaussian Processes for few-shot learning but only support vinglevvrrabll prediction. In contrast, MOSES handles multiple variables and missing values.
6 EXPERIMENTS
6.1 TOY EXPERIMENT
We illustrate marginalization consistency using two synthetic bivariate distributions (Blast and Circle; see Figure Appendix B). The task is to estimate the unconditional joint distribution. MOSES accurately models both joint and marginal distributions while preserving consistency. In contrast, ProFITi captures the joint well - especially for Blast - but fails on marginals due to its triangular attention mechanism, which enforces a fixed dependency order. Gaussian Process Regression maintains consistency but lacks predictive accuracy. To predict the marginalization inconsistency use 2- Wasserstein distance (WD). For each variable yk. we compare:
1. p(yk Qk X): the predicted marginal,
2. fMA (yk Qk, X): the marginal obtained by integrating the joint p(y Q X).
Since sampling directly from fMak is difficult, we sample from the joint and extract the k-th compo- nent. The marginalization inconsistency is defined as the average WD across variables:
K WD M(n Qk, x), MAR (yk
MI
Qu, X)).
(17)
K
We use 1000 samples to compute MI. Currently, we compute only univariate marginals; multivariate marginals are possible in principle but computationally prohibitive.
6.2 MAIN EXPERIMENT
We evaluate our model on four real-world datasets: one climate dataset (USHCN) and three medical datasets (PhysioNet" 12, MIMIC-III, and MIMIC-IV). Following prior work (Yalavarthi et al. 20255] Biloš et al. 2021). observe the first 36h and predict the next 3 time steps for medical datasets, and observe years and predict time steps for USHCN. Both the number of observations (N) and queries (K) vary across samples (see Table 5). We split each dataset into training, validation, and test sets using a 70: 10:10:20 ratio. We train MOSES using the Adam optimizer with a learning rate of 0.001 and batch size of 64. Hyperparameter search is over mixture components D E {1,3,5, 5, 7, 10).
8

---

Published as a conference paper at ICLR 2026
Table 2: Trained for njNLL and evaluate for mNLL, lower is better.
Model
USHCN
PhysioNet' 12
MIMIC-I11
MIMIC-IV
inconsistent
ProFITi
-3.324 # 0.206
-0.016 0.085
0.408 0.030
0.500 0.322
GRU-ODE
0.776 = 0.172
0.504 0.061
0.839 = 0.030
0.876 = 0.589
consistent
Neural-Flows
0.775 0.180
0.492 0.029
0.866 + 0.097
0.796 0.053
univariate
CRU
0.762 0.180
0.931 0.019
1.209 1 0.044
OOM
Tripletformer+
0,411 = 7.506
0.110
0.894 : 0.083
0.751 = 0.063
GPR
1.235 0.096
1.161 0.065
1.341 = 0.009
1.161 = 0,010
consistent multivariate
GMM
1.042 # 0.021
1.069 = 0.002
1.124 # 0.007
1.075 = 0.007
MOSES (ours)
-3.355 0.156
-0,271 0.028
0.163 = 0.026
-0.634 0.017
MOSES
ProFITi
CRU
NerualFlow
GRU-ODE
GMM
33
-1
10 10"h
10
10
10
10*
10-A
10
10-
10
10
10*
MIL USHCN
MTJ
MT MIMIC-ITI
MIJ MIMIC-IV
PhysioNet' 12
Figure 4: njNLL vs. MI. MOSES is marginalization consistent within sampling error.
attention heads E (1, 2), and latent sizes M, F E (16,32, 64, 128]. All models are implemented in Poroh and trained on NVIDIA RTX 3090, A40 and GTX 1080 Ti GPUs.
Baselines. We use NeuralFlows Bilos et al., 2021 GRU-ODE (De Brouwer et al.
2019
CRU (Schirmer al. 2022), GPR (Durichen et al. 2013), and ProFITi (Yalavarthi et al. 2025
Our encoder is similar to Tripletformer (Yalavarthi et al. (223) that predict marginal distributions for interpolation. We used it for forecasting and called the model Tripletformer+. NeuralFlows, GRU-ODE, CRU, and Tripletformer+ predict only marginals and are marginalization consistent, as their joint distribution the product marginals. GPR also marginalization consistent. We also compare with Gaussian Mixture Model (GMM) which is MOSES without flows attached to highlight the advantage of flows in MOSES.
Results. To highlight the importance of Marginalization Consistency in probabilistic forecasting models train the model for njNLL (Equation (16)) and evaluate for two metries: 1. Normalized Joint Negative Log-Likelihood (njNLL: Table 1) and 2. Marginal Negative Log-Likelihood (mNLL: Table 2). While njNLL measures joint density, mNLLL assesses univariate marginal density (Biloš
et al. 2021 Schirmer et al. 2022). A good model must perform well on both the metrics.
MOSES outperforms all narginalization-coniiin similarly or slightly worse than ProFITi on njNLL, it significantly surpasses ProFITi on mNLL. For USHCN, both models perform comparably within standard deviation. From Figure our model achieves similar likelihoods to ProFITi, with MI values close to zero, while MI values of ProFITi are up to an order of magnitude larger. This difference arises because ProFITi prioritizes joint distributions but neglects marginalization consistency, leading to performance degradation when queried on single instances (Figure T). In contrast, our model maintains consistency. Similar results observed when comparing Energy Score (for multivariate distributions) and CRPS (for univariate distributions) as seen Table
models across both metrics. Although performs
We noticed that ProFITi's gains primarily stem from its encoder. When using the same encoder (ProFITi-TF), MOSES achieves superior accuracy in MIMIC-III and MIMIC-IV TTall 4). Further- more, when we replaced our encoder with the encoder of ProFITi which is GraFITi, MOSES-GraFIT has significantly better results than ProFITi in both MIMIC-II and MIMIC-IV IV. However, since GraFITi is an inconsistent model, GraFT has worse mNLL for all the datasets.

---

Published as a conference paper at ICLR 2026
Table 3: Comparing models w.r.t. Energy Score and CRPS. Lower is better, best results in bold.
Energy Score
CRPS
Dataset
ProFITi
MOSES (ours)
ProFITi
MOSES (ours)
USHCN
0.452 + 0.044
0,552 # 0.044
0.182 # 0.007
0.220 + 0.019
PhysioNet'12 12
0.879 + 0.303
1.599 0.013
0.271 0.003
0.260 = 0.002
MIMIC-III
1.606 0.168
1.353 0.033
0,319 0.003
0.296 + 0.005
MIMIC-VV
0.808 + 0.003
0.906 = 0.029
0.279 = 0.012
0.245 + 0.010
Table 4: Comparing njNLL and mNLL across datasets verify the contribution of probabilistic component of ProFITi. "ProFITi-TF" denotes ProFITiiTTrasformn using same encoder as MOSES. "MOSES-Gratee GraFITi" is MOSES with GraFITi encoder.
Dataset
USHCN
PhysioNet' 12
MIMIC-I1
MIMIC-IV
njNLL
ProFITi
-3.226++0222 2.25
-0.647=00078
-0.377=0.032
-1.777=00066
ProFITi-TF
-3.415+0.271
-0.657=0.034
0.51620.111
-1.405=0:220
MOSES-Gra GraFITi
-3,282+0.22 260
-0.615=0.030
-0.396=0.025
-2.197=0.027
MOSES
-3.35740.176
00.491400041
-0.305=0.027
-1.66820.097
mNLL
ProFITi
-3.324+0.206 206
-0.01620.085
0,408=0.030
0,500=0.322
ProFITi-TF
+3.440+0.243
0.017=0.042
1.279=0.057
0.345=0.325
MOSES-GraFIl MOSES
-3.252+0.25
-0.227=0.113
0.520=0.125
-0.299=0.1 1114
-3.355+0.156
-0.271+0.028
0.163=0.026
-0.634=0.017
Discussion. While MOSES provides strong theoretical guarantees through marginalization consis- tency, this comes with an inherent trade-off in expressivity. The separability constraint necessary for tractable marginalization limits the model's ability to capture complex dependencies directly through the flow transformations. In MOSES, dependencies are primarily captured through the low-rank covariance structure latent Gaussian distribution and mixtures, while the expressive power of normalizing flows focused modeling flexible marginal distributions.
This architectural choice explains why MOSES achieves superior marginal performance and perfect consistency, while ProFITi-wiit its fully non-separable flows-can sometimes achieve better joint distribution modeling at the cost of consistency. In applications like healthcare forecasting, where reliability and interpretability are crucial, the consistency guarantees of MOSES provide great value despite this trade-off.
In future work, we aim explore more flexible architectures that maintain marginalization consis- tency, such as through copula models or probabilistic circuits. A key challenge will be adapting these methods handle the irregular time series setting.
CONCLUSIONS
We introduced MOSES, marginalization connssstttttnnnnniittetnttt mixture of separable flows for probabilistic forecasting of irregular time series with missing values. By carefully parametrizing its components, we ensured both decomposability and marginalization consistency. Experiments on four real-world irregularly sampled datasets show that MOSES achieves substantial gains in marginal predictions while preserving competitive accuracy in joint predictions compared to more flexible but inconsistent models. This work provides initial step toward addressing marginalization inconsistency in probabilistic forecasting. Future research will likely focus on enhancing performance while preserving consistency. A key limitation of MOSES is that mixture weights independent of the query Q. Although adapting weights to query time points could improve joint predictions, requirements RI enforce muery-nndependence
necessary condition for marginalization consistency.
10

---

Published as a conference paper at ICLR 2026
USAGE OF LARGE LANGUAGE MODELS
We would like to note that all scientific contributions in this manuscript such as conceptual ideas, problem identifications, model development, theorems, proofs, experimental results, discussion and so on are entirely the work of the authors. The usage of Large Language Models is limited to minor linguistic refinement such as correcting grammar and typos.
REFERENCES
Abdul Fatir Ansari, Alvin Heng, Andre Lim, and Harold Soh. Neural continuous-discrete state space models for irregularl-samppl time series. In Andreas Krause, Emma Brunskill, Kyunghyun Cho, Barbara Engelhardt, Sivan Sabato, and Jonathan Scarlett (eds.), Proceedings of the 40th International Conference on Machine Learning, volume 202 of Proceedings of Machine Learning Research, pp. 926-951. PMLR, 2023.
Arjun Ashok, Étienne Marcotte, Valentina Zantedeschi, Nicolas Chapados, and Alexandre Drouin. TACTiS-2: Better, faster, simpler attentional copulas for multivariate time series. In The Twelfth International Conference on Learning Representations, 2024.
Marin Biloš and Stephan Günnemann. Normalizing flows for permutation invariant densities. In Marina Meila and Tong Zhang (eds.), Proceedings of the 38th International Conference on Machine Learning, volume 139 of Proceedings of Machine Learning Research, pp. 957-967. PMLR, 2021.
Marin Biloš, Johanna Sommer, Syama Sundar Rangapuram, Tim Januschowski, and Stephan Gün- nemann. Neural flows: Efficient alternative to neural ODEs. Advances in Neural Information Processing Systems, 34:21325-21333 2021.
Edwin V Bonilla, Kian Chai, and Christopher Williams. Multi-task gaussian process prediction. In Advances in Neural Information Processing Systems, volume 20, 2007.
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):1-12, 2018. doi: 1111038/41159--18--2444
Yuqi Chen, Kan Ren, Yansen Wang, Yuchen Fang, Weiwei Sun, and Dongsheng Li. ContiFormer: Continuous-time transformer for irregular time series modeling. Advances in Neural Information Processing Systems, 36, 2024.
YooJung Choi, Antonio Vergari, and Guy Van den Broeck. Probabilistic Circuits: A Unifying Framework for Tractable Probabilistic Models+. 2020.
Sebastian Ciobanu. Mixtures of normalizing flows. In Proceedings of ISCA 34th International Conference on Computer Applications in Industry and Engineering, EPiC Series in Computing, volume 79, pp. 82-90, 2021.
Edward De Brouwer, Jaak Simm, Adam Arany, and Yves Moreau. GRU-ODE- Continuous modeling of sporadically-obseree time series. Advances in Neural Information Processing Systems, 32, 2019.
Laurent Dinh, Jascha Sohk-Dickstein, and Samy Bengio. Density estimation using real NVP. In International Conference on Learning Representations, 2017.
Hadi Mohaghegh Dolatabadi, Sarah Erfani, and Christopher Leckie. Invertible generative modeling using linear rational splines. International Conference on Artificial Intelligence and Statistics, PMLR, 2020.
pp.
Alexandre Drouin, faienne Marcotte, and Nicolas Chapados. Tactis: Transformmer-attentional copulas for time series. In International Conference Machine Learning, pp. 5447-5493. PMLR, 2022.
Richard O. Duda and Peter E. Hart. Pattern classification and scene analysis. In A Wiley- Interssintercceeece Publication, 1974.
11

---

Published as a conference paper at ICLR 2026
Robert Dürichen, Marco A. F. Pimentel, Lei Clifton, Achim Schweikard, and David A. Clifton. Multitask gaussian processes for multivariate physiological time-series analysis. IEEE Transactions on Biomedical Engineering,
2015. doi: 1LL10////ME.221225553333
Conor Durkan, Artur Bekasov, lain Murray, and George Papamakarios. Neural spline flows. Advances in neural information processing systems, 32, 2019.
Max Horn, Michael Moor, Christian Bock, Bastian Rieck, and Karsten Borgwardt. Set functions for time series. In International Conference on Machine Learning, pp. 4353 1433-433.. PMLR, 2020.
Johnson, Bulgarelli, T Pollard, S Horng, and LA Celi. Mark. MIMIC-IV (version 1.0).
PhysioNer, 2021.
Alistair EW Johnson, Tom Pollard, Lu Shen, Li-wei Lehman, Mengling Feng, Mohammad Ghassemi, Benjamin Moody, Peter Szolovits, Leo Anthony Celi, and Roger G Mark. MIMIC-III, a freely accessible critical care database. Scientific data, 3(1):1-9, 2016...
Seyed Mehran Kazemi, Rishab Goel, Sepehr Eghbali, Janahan Ramanan, Jaspreet Sahota, Sanjay Thakur, Stella Wu, Cathal Smyth, Paseal Poupart, and Marcus Brubaker. Time2 Learning Vector Representation of Time, July 2019. URL http:::/aaxxiv.orr///aaaaa// 11Xi01107..0321 [cs].
Christian Klötergens, Vijaya Krishna Yalavarthi, Maximilian Stubbemann, and Lars Schmidt-Thiees Functional latent dynamics for irregularly sampled time series forecasting. In Joint European Conference Machine Learning and Knowledge Discovery in Databases, pp. 421-436. Springer, 2024.
Christian Klötergens, Vijaya Krishna Yalavarthi, Tim Dernedde, and Lars Schmidt-Thieme. Mixing it up: Exploring mixer networks for irregular multivariate time series forecasting, 2025a. URL https://arxiv. wttp://wwww..rr/p/////////////// 11816
Christian Klötergens, Vijaya Krishna Yalavarthi, Randolf Scholz, Maximilian Stubbemann, Stefan Born, and Lars Sehmidt-Thiemme Physiome-ODE A benchmark for irregularly sampled multivari- ate time-series forecasting based on biological ODEs. The Thirteenth International Conference Learning Representations, 2025b.
Balaji Lakshminarayanan, Alexander Pritzel, and Charles Blundell. Simple and scalable predictive uncertainty estimation using deep ensembles. Advances in Neural Information Processing Systems (NIPS), 30, 2017.
Jenny Liu, Aviral Kumar, Jimmy Ba, Jamie Kiros, and Kevin Swersky. Graph normalizing flows. Advances in Neural Information Processing Systems, 32, 2019.
Matthew J Menne, CN Williams Jr, and Russell S Vose. United States historical climatology network daily temperature, precipitation, and snow data. Carbon Dioxide Information Analysis Center, Oak Ridge National Laboratory, Oak Ridge, Tennessee, 2015.
Kevin P. Murphy.
Probabilistic Machine Learning: An Introduction. Adaptive Computation and Machine Learning Series. MIT Press, Cambridge, MA, USA, March 2022. ISBN 978-0-262- 04682-4.
Bernt Øksendal. Stochastic Differential Equations. Universitext. Springer, Berlin, Heidelberg, 2003. ISBN 978-3-540-04758-2 1---3-642-14394--- doi: 10.1007/978-33-42--444 14394-6.
George Papamakarios, Theo Pavlakou, and lain Murray. Masked autoregressive flow for density estimation. In I. Guyon, U. Von Luxburg. S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett (eds.), Advances in Neural Information Processing Systems, volume 30. Curran Associates, Inc., 2017.
George Papamakarios, Eric Nalisnick, Danilo Jimenez Rezende, Shakir Mohamed, and Balaji Lakshminarayanan. Normalizing flows for probabilistic modeling and inference. Journal of Machine Learning Research, 22(1), January 2021. ISSN 1532-4435.
12

---

Published as a conference paper at ICLR 2026
Guilherme G. P. Freitas Pires and Mário A. T. Figueiredo. Variational mixture of normalizing flows. In 28th European Symposium on Artificial Neural Networks, Computational Intelligence and Machine Learning, ESANN 2020, Bruges, Belgium, October 2-4, 2020, pp. 205-210, 2020.
Janis Postels, Mengya Liu, Riccardo Spezialetti, Lue Van Gool, and Federico Tombari. Go with the flows: Mixtures of normalizing flows for point cloud generation and reconstruction. In 2021 International Conference on 3D Vision (3DV), pp. 1258. IEEE, 2021. doi: 10.1109/ 3DVS3792.2021.000
Kashif Rasul, Abdul-Saboor Sheikh, Ingmar Schuster, Urs Bergmann, and Roland Vollgraf. Multi- variate probabilistic time series forecasting via conditioned normalizing flows. In International Conference on Learning Representations, 2021.
Danilo Rezende and Shakir Mohamed. Variational Inference with Normalizing Flows. In International Conference on Machine Learning, pp. 1530-1538. PMLR, June 2015.
Mona Schirmer. Mazin Eltayeb, Stefan Lessmann, and Maja Rudolph. Modeling irregular time series with continuous recurrent units. In Proceedings of the 39th International Conference on Machine Learning, volume 162, pp. 19388-19406 PMLR, 2022.
Randolf Scholz, Stefan Born, Nghia Duong-Trung, Mariano Nicolas Stuz-Bouraazou, and Lars Schmidt-Thiieme Latent linear ODEs with neural kalman filtering for irregular time series forecasting, 2023.
Maximilian Seitzer, Arash Tavakoli, Dimitrije Antic, and Georg Martius. On the Pitfalls of Het- eroscedastic Uncertainty Estimation with Probabilistic Neural Networks. In International Confer- ence on Learning Representations, 2021.
Marcin Sendera. Jacek Tabor, Aleksandra Nowak, Andrzej Bedychaj. Massimiliano Patacchiola, Tomasz Trzcinski, Przemystaw Spurek, and Maciej Zieba. Non-gaussian gaussian processes for few-shot regression. Advances in Neural Information Processing Systems, 34: 34:10285-10298, 2021.
Satya Narayan Shukla and Benjamin Marlin, Heteroscedastic temporal variational autoencoder for irregularly sampled time series. In International Conference on Learning Representations, 2022.
Sahil Sidheekh, Kristian Kersting, and Sriraam Natarajan. Probabilistic Flow Circuits: Towards Unified Deep Models for Tractable Probabilistic Inference. In Proceedings of the Thirty-Ninth Conference on Uncertainty in Artificial Intelligence, pp. 1964-1973. PMLR, July 2023.
Ikaro Silva, George Moody, Daniel J Scott, Leo A Celi, and Roger G Mark. Predicting in in-hospital mortality of icu patients: The physionet/computing in cardiology challenge 2012. In 2012 Computing in Cardiology, pp. 245 245-248. IEEE, 2012.
Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N Gomez, Lukasz Kaiser, and Illia Polosukhin. Attention is All you Need. In I. Guyon, U. V. Luxburg, S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett (eds.), Advances in Neural Information Processing Systems 30, pp- 5998-6008. Curran Associates, Inc., 2017.
Andrew G Wilson and Zoubin Ghahramani. Copula processes. Advances in Neural Information
Processing Systems, 23, 2010.
Vijaya Krishna Yalavarthi, Johannes Burchert, and Lars Schmidt- Tripletformer for Proba- bilistic Interpolation of Irregularly sampled Time Series. In 2023 IEEE International Conference Big Data (BigData), pp. 986-995, December 2023. doi: 10. D0.1109/BigDaaa550444.0000 20230003
Vijaya Krishna Yalavarthi, Kiran Madhusudhanan, Randolf Scholz, Nourhan Ahmed, Johannes Burchert, Shayan Jawed, Stefan Born, and Lars Schmidt-Thiee GraFITi: Graphs for Forecasting Irregularly Sampled Time Series. In Michael J. Wooldridge. Jennifer G. Dy, and Sriraam Natarajan (eds.). Thirty-Eightt AAAI Conference on Artificial Intelligence, AAAI 2024. Thirty-Sixth Con- ference Innovative Applications of Artificial Intelligence, IAAI 2024, Fourteenth Symposium on Educational Advances in Artificial Intelligence, EAAI 2014, February 20-27, 2024, Vancouver, Canada, pp. (6255-16263 AAAI Press, 2024. doi: 10. 10.1609/AAA1 13811555555
13

---

Published as a conference paper at ICLR 2026
Vijaya Krishna Yalavarthi, Randolf Scholz, Stefan Born, and Lars Schmidt-Thieme. Probabilistic Forecasting of Irregularly Sampled Time Series with Missing Values via Conditional Normalizing Flows. Proceedings of the AAAI Conference Artificial Intelligence, 39(20): April doi: 10. 10.9/aa/aaii.... 1.09///aii.339920....... v39120.35494.
2025. ISSN 2374-3468.
Weijia Zhang, Chenlong Yin, Hao Liu, Xiaofang Zhou, and Hui Xiong. Irregular Multivariate Time Series Forecasting: Transformable Patching Graph Neural Networks Approach. In Firtt-First International Conference on Machine Learning. June 2024.
14

---

Published as a conference paper at ICLR 2026
A THEORY
A.1 PROOF OF EXTENSION OF EQUATION (5) TO ARBITRARY FINITE SETS
We show that Equation (5) is equivalent to the more general statement that for any finite subset Ks C {l,..., K holds:
/ |2|kgl p(y Q) dyks
p(y-Ks Ks X)
(18)
Proof. Obviously. (8) is equivalent to (5) when Kgl 1 via Ks [k}. The other di-
rection can shown by induction on size of Ks. Assume that (8) holds for any sub- set of size at most and let Ks be subset of size m + 1. w.l.o.g., we may assume that Kg = {L...., . 1] = [L...., m] ) {m + 1}. Let ÿ = y [l....mm (fm-1..... yk) and
Q = Q- (1 = (Qm+1, Qx), then y- K2 セ- and Q K5 Q-1 1. It follows that
ply-Ks 1QQQks. Ks, X) = p 1 QQ
plym+1...... yK Qmtl..... ..K X) dy
(f Qi:K X) dyn dyn..... dym+1
/ Munx Qu.K. X) dytmmt1
A.2 PROOF OF LEMMA 3.1
Proof. Since X is a common conditional to all the marginals, we can ignore it. So, assume that f is separable transformation:
f(x Q) (d(21 Q1)....,. 6(zK Qk))
(19)
and that (z(2) Q) is marginalization consistent model. Then, the predictive distribution is
ply Q) P2 (y Q) Q) det df df-1(y Q)
(20)
dy
Since is separable, it follows that the Jacobian is diagonal:
df (y Q)
a a (yn Q1). ... Φ-'(xx (yK (Qx))
dy
d(yi.... ..... MK)
de d Q1)
doz Qk))
diag
(21)
dy1
dyk
Hence, the determinant of the Jacobian is the product of the diagonal elements:
df d (y Q)
dø
(yk Qx)
do (yl Qk)
det
det
II
(22)
dy
A=I:K dyj
dys
k=I:
15

---

Published as a conference paper at ICLR 2026
Using this fact, integrate the joint density over y/ to get the marginal density:
f p(y Q)dyy
/ oz Q) 1
a
IQ)
det
dyk
dy
-f 1z (y) 11)))
dΦ
(yk
II
dy
22
dy
k=1:k
II 1
de
(un Qi)
f pz
(f- (1)
dd dy
d$
[Qx)
dyk
dyn
dØ
(yn
Qi)
-f pz(( (Q) dzk
II 1
transf.-thm
dyr
dØ
(Mi
Qi)
II 1
(2(# 1Q-k)
dyı
df
(V- 1-x)
p2(2=k Q-+k
det
dy-k
p(y-k k (Q-k)
A.3 PROOF OF LEMMA 3.2
Proof. Consider a mixture model of the form
ply Q,
Σ (a) Q)
(23)
satisfying the conditions from Lemma 3.2 i.e. the component models pa satisfy the requirements R1
R3 and the weight function w: Seq(X)) A is permutation invariant with respect to X.
1. p satisfies R1
By construction of the mixture model, it has the same domain and codomain
as the component models.
2. P satisfies R2
Let E Sym([Q)) and E (y)([X)), then
p(y Q x²)
wa(X pa(s"" 1Q,x*)
D Σ
wa(X))a(v Q,
permutation invariance w and Pa
d=1
ply Q)
3. P satisfies R3
/ p(y Q, x) dyn / Σ wa((X)ala 1)
Σ wa(X)) / Paly 11X)) dyk
d=1
Σ waa(X (X))x(X) X) Pa is marginalization consistent
d=1
ply- & QQ-k,
16

---

Published as a conference paper at ICLR 2026
A.4
PROOF OF HEOREM 4.
Proof. Due Lemma 3.2) it sufficient to show that all the component models satisfy the require-
ments RI R3 and that the mixture weights are permutation invariant with respect to X.
The latter immediately follows from (aa and (3): since hoas is computed by self-attention over x,
the encoding of X, it is permutation equivariant with respect the sequence elements, i.c.
OBS (x') = (ho*s(x))
for any T E Sym([XX)
(24)
Then, w is computed by cross-attention over hoss, and hence is permutation invariant with respect to X, as the sequence length is summed over in Equation (15)
joint multivariate prediction, by construction.
R3
marginalization consistency, is satisfied by the component models due to Lemma 3.1 since separable Gaussians as source distributions and separable flows as transformations.
use
R2
permutation invariance of the component model Pa can be seen follows:
As already mentioned, hsss permutation equivariant with respect to X. By Equation (12b)
follows that b2 is permutation invariant with respect to X, since the observation length summed over, and permutation equivariant with respect to Q:
ha("", x*) ha(Q, x)" for any и Sym(IQl), E Sym(IX))
(25)
Now, consider the d-th component model Py, (y Q. Let : fo f [v Q, By Equation (14)
and the gguivatiaace//ivariancc of ha w.r.t. Q2 and X. our flow model satisfies
'v" 1=, = fy" h = (f- ha))" =
(26)
I.e. z is permutation equivariant with respect to Q and invariant with respect to X. Therefore:
Pva Q",
det d (v 1Q)
of (v" (T) Q" (7)
dy
df (yy h7)
=N((f- h7) ((h)), E(hT))
dot
by remark above
dy
df
(y/ hg)
= =N((" p((T), E(h))) det
dy
NN(= det df (y/T h
Di permutation invariance of GP
dy
det d(f
(y) hx)) 8
-N E)
by 26
dy
det df (y ha)
-N(2 (2,)).
by 22 (product is permutation invariant)
dy
PY, (v) (2
A.5 LINEAR RATIONAL SPLINES
Linear Rational Splines (LRS) are computationally efficient spline functions Dolatabadi et al. (2020).
Formally, given set of monotonically increasing points {wum Um hmmm)mmm... 1:M called knots, that is
Um < um+11 and Unp < Um+11 along with their corresponding derivatives (Am > 0) mm then the
LRS transformation (uu) within bin u E [um, Um+1] is:
F'm (Am -6] -6)+mmm Bm
00 Am
Φ(u)=
aon (Am -8) amm
where á
Um
E |0,
Om Pm 400m
Am S u u
Um Um
&m (1 +m Am))
(27)
17

---

Published as a conference paper at ICLR 2026
1) signifies the location of automatically inserted virtual knot between Um and Um+11 with value Unc- The values of Am. O', am and Em are all automatically derived from the original knots and their derivatives Dolatabadi et al. (2220). For conditional LRS o(2g; h 0), the function parameters such as width and height of each bin, the derivatives at the knots, and X are computed from the conditioning input a and some model parameters 0. 0 helps project ha, the function parameters, and is common to all the variables 21:K so that the transformation & can be applied for varying number of variables K. Additionally, we 0 common all the components as well. Since, each component has separate embedding for variable z (ha,)), we achieve different transformations in different components for same variable.
Here, Am 10
In summary, the conditional flow model is separable across the query size f = f1 x ... x fK with
fa(z) == f(= ha) = (0(21, (4,),...,((2K,, ha,k))
(28)
18

---

Published as a conference paper at ICLR 2026
Table Statistics of the datasets used in our experiments. Sparsity means the percentage of missing observations in the time series. N is the total number of observations and K is the number of queries in our experiments Section 6
Name
#Samples
#Channels
Sparsity
K
USHCN
1100
77.9%
322
PhysioNet'12
12,000
37
85.7%
519
53
MIMIC-III
21,000
96
94.2%
4-709
MIMIC-IV
18,000
102
97.8%
1-1382
1 79
B DATASETS
REAL-WORLD DATASETS.
real-world datasets are used in our experiments, based on prior works in irregular time series forecasting. Their statistics are summarized in Table 5
USHCN Menne et al. 2015).. This is a climate dataset consisting of 5 climate variables such as daily temperatures, precipitation and measured 150 years 1218 meteorological stations the USA. Following De Brouwer et al. (2019); Yalavarthi et al. 2025), we selected 1114 stations and observation window of years from 1996 until 2000.
PhysioNet2012 Silva et al. (2012). This physiological dataset consists of the medical records of 12,000 patients who are admitted into ICU. 37 vitals are recorded for 48 hrs. Following the protocol of Yalavarthi et al. (2024) Che et al. (2018), dataset consists of hourly observations in each series.
MIMIC-III Johnson et al. (2016). This also a physiological dataset. it a collection of readings of the vitals of the patients admitted to ICU at Beth Israeli Hospital. Dataset consists of 18,000 instances and 96 variables are measured for 48 hours. Following De Brouwer et al. (2019): Bilos
al. 2021): Yalavarthi et al. (2025) observations are rounded 30 minute intervals.
MIMIC-IVJJohnonn et al. (2021). The successor of the MIMIC-II dataset. Here, 102 variables from patients admitted to ICU tertiary academic medical center in Boston are measured for 48 hours. Following De Brouwer et al. (2019) Biloš et al. (2021); Yalavarthi et al. (2025) observations,
MIMIC-IV Johnson
are rounded to minute intervals.
SYNTHETIC DATASETS.
We also synthetic datasets evaluate the marginalization consistency of the models.
Circle (toy dataset). The circle dataset consists of samples from a unit circle in R² with added Gaussian noise. Formally, the data is generated as follows:
y
0.05 (N I2) where z ~~ N0,
(29)
|=|2
Blast distribution (toy dataset). This is bivariate, non-isotropic distribution defined as follows:
y sign(=) ② ④ where 2~~ (80
[ 1
(30)
2
19

---

Published as a conference paper at ICLR 2026
Table We compare models with respect MSE, where lower values indicate better performance. The perf-peroogning model is shown in bold, the second-best in italics, and models with performance within one standard deviation of the best are underlined
Model
USHCN
PhysioNet'12
MIMIC-III
MIMIC-IV
inconsistent
ProFITi
0.308 0.061
0.305 0.007
0.548 0.063
0.389 + 0.015 0.365 + 0.012
GRU-ODE
0.410 0.106
0.329 0.004
0.479 0.044
consistent
Neural-Flows
0.424 0.110
0.331 0.006
0.479 0.045
0.374 + 0.017
univariate
CRU
0.290 0.060
0.725 = 0.037
OOM 0.369 + 0.030
Tripletformer+
0.349 = 0.131
0.293 0.018
0.547 # 0.068
GPR
0.597 0.1
0.575
0.059
0.862 + 0.016
0.609 + 0.014
consistent
GMM
0.294 0.083
0.293 0.005
0.535 # 0.064
0.332 + 0.015
multivariate
MOSES (ours)
0.411 0.099
0.307 0.006
0.517 + 0.057
0.342 + 0.028
C ADDITIONAL EXPERIMENTS
C.I COMPARING FOR POINT FORECASTING
While point forecasting important task in time series analysis, the goal of probabilistic forecasting fundamentally different. Probabilistic forecasting aims to capture the full predictive distribution rather than just single-point estimate. Nonetheless, one might intuitively expect that the best probabilistic model would also yield the most accurate point estimates, However, this is not always
the case in practice, as noted in prior works (Lakshminarayanan al. 2017 Seitzer et al. 2021 Rasul et al. 2021 Yalavarthi et al. 2025).
We compare probabilistic models in terms of point prediction accuracy using Mean Squared Error (MSE), as reported in Tabll 6) Our results show that no single model consistently outperforms the others across all datasets. GMM is the overall best ranked model with MOSES and GRU-ODE being the next best. We believe there are two primary reasons for this phenomenon:
(1.) MSE is related to the Negative Log-Likelihood (NLL) of Gaussian distribution with a fixed standard deviation. Therefore, models explicitly trained by minimizing Gaussian Negative Log- Likelihood (even if they predict more than just the mean) are naturally optimized for this metric.
(2.) Probabilistic models trained to predict the underlying data distribution, not solely the optimal point estimate (c.g., the conditional mean). Their objective is to accurately capture the uncertainty and dependencies in the data, which involves learning the (co)variance structure. This focus on the full distribution can sometimes lead to point estimates that are not strictly optimized for minimizing the squared error, even if the overall probabilistic forecast is superior.
Except for MOSES, GMM and ProFITi, all the other probabilistic models are designed to predict Gaussian distributions. MOSES and GMM have competing accuracy in all the datasets other than USHCN. In USHCN, MOSES performs worse because of outliers in the data. This is the reason for large standard deviations for all the models in USHCN dataset.
C.2 EXPERIMENT ON VARYING OBSERVATION AND FORECAST HORIZONS
We would like to see if MOSES scalable to long observations and forecast horizons. For this, we performed an experiment on varying length observation and forecasting horizons on PhysioNet') 12 dataset and compared against the published results from (Yalavarthi et al., 2025) in Table The observation and forecasting horizons are: [(36h, 12h), (24h, 36h), (12h, 36h)).
Tables and present the njNLL and mNLL results for ProFITi and MOSES. The results follow the trends observed in Tables [ and 2 ProFITi performs best when predicting joint distributions. However, its lack of marginalization consistency leads to a severe performance drop when predicting marginal distributions. In contrast. MOSES maintains stable performance from njNLL to mNLL. While it performs slightly worse than ProFITi on njNLL, it significantly outperforms ProFITi on MNLL.
20

---

Published as a conference paper at ICLR 2026
Table Experiment on varying observation and forecast horizons. Evaluation metric-njNLL, Lower the better. Best results are in bold.
36/12
24/24
12/36
NeuralFlows
0.708 0.048
1.097:0.044
1.436+0. 187
ProFITi
-0.768=0.041
-0.355=0.243
-0,291=0.415
MOSES (ours)
-0.315#0.016
-0.298+0.027
-0.06330.04p
Table 8: Experiment on varying observation and forecast horizons. Evaluation metric-mNLL, Lower the better. Best results are in bold.
36/12
24/24
12/36
NeuralFlows
0.709=0.047
1.06==0.022
1.439=0.198
ProFITi
1.376=1.764
0.705=0.177
2.977=2.978
MOSES (ours)
-0.083+0.025
-0.020+0.060
0.040=0.131
C.3 ABLATION STUDY.
Using PhysioNet 2012, we show the importance of different model components. As summarized in Table 9 the performance reduced by removing the flows (MOSES - f) which is same GMM. is expected that normalizing flows are more expressive compared to simple mixture of Gaussians. On the other hand, by using only isotropic Gaussian the base distribution (MOSES COV) model performance worsened. Similarly, parameterizing the components weights have a slight advantage over fixing them to 1/D with being the number of components. One interesting observation is even using single component (MOSES(1)) gives similar results compared to mixture of such components. This could be because the dataset we have may not require multiple components. We note that we have D = 1 in our hyperparameter space, and select the best D) based on validation dataset. Note that (M)SES(1)) can be seen as representative of Copula Processes (Wilson & Ghahramani 2010) and Non-Gaussian Gaussian Process (Sendera et al. 2021).
When we conduct this study on other datasets Table 10 we realized that mixer model indeed helpful if there is a requirement.
C.4 COMPARING THE NUMBER OF PARAMETERS AND RUNTIME FOR GMM AND MOSES
Since MOSES is built upon GMM, Table 11 presents the number of parameters and runtime for both MOSES and GMM. For reference, we also include ProFITi.
The results show that GMM has relatively low number of parameters for MIMIC-III and MIMIC-IV, whereas for USHCN and PhysioNet' 12, the number of parameters is significantly higher. The primary difference between GMM and MOSES is the inclusion of flows. Given that all other factors remain the same, MOSES is expected to have a slightly higher number of parameters than GMM due to these additional flows. Also, the parameters for the flows are shared among all the variables and the components, their number does not grow with increase in components or variables. However, differences in the chosen hyperparameters for GMM and MOSES lead to some discrepancies from this expectation. Moreover, the inclusion of flows in MOSES results in slightly higher runtime compared to GMM.
C.5 SENSITIVITY ANALYSIS
We conducted a sensitivity analysis on key hyperparameters: the number mixture components (D), the latent embedding size (M), and the number of flow layers. The results shown in Figure 5) For each dataset, we started from the best hyperparameter set and varied one parameter at a time.
Our analysis reveals that the model is most sensitive to the number of flow layers, followed by the latent embedding size. This aligns with our architectural design: increasing the number flow layers directly enhances the expressiveness of the univariate transformations, which crucial for modeling
21

---

Published as a conference paper at ICLR 2026
Table 9: Ablation study on PhysioNet' 12
Table 10: Importance of mixer components
Model
njNLL (J)
Dataset
MOSES(()
(OSES(D)
MOSES
-0.491 0.041
USHCN
-3.255=0.200
-3.357=0.176
MOSES-/
1.063 0.002
PhysioNet' 12 MIMIC-III
-0.493=0.029
-0.491=0.041
MOSES-COM
-0.308 = 0.024
0.173=0.332
-0.305=0.027
MOSES-W
-0.451 = 0.038
MIMIC-IV
-1.68620.072
-1.668=0.097
MOSES (I)
-0.493 # 0.029
Table 11:
Comparing the number of parameters and run-time per epoch for results in Table 1 GMM and MOSES, and ProFITi for reference
for
USHCN
PhysioNet'12
MIMIC-II III
MIMIC-IV IV
Parameters
Time
Parameters
Run Time
Parameters
Run
Parameters
Time
ProFITi GMM
1.093.0K
3.8s
75.8K
42.14s
59.7K
66.8s
285.9K
70.2s
416.0K
0.9s
390,9K
5.9s
33.0K
18.5s
101.1K K
21.3s
MOSES (ours)
167.6K
2.46
134.6K
14.1s
112.6K
25.46
398.6K
33.3s
complex marginal distributions. The latent embedding size is also important, as it parameterizes both the marginal distributions (via the flow conditioning) and the covariance structure of the latent Gaussian, which captures dependencies between variables.
In contrast, the model exhibits robustness to the number of mixture components. A single component (D 1) is often sufficient to capture the primary dependency structure in the data, and increasing D does not always improve the accuracy. This indicates that the core dependency structure can often be captured by a single, well-structured Gaussian component for many datasets. However, the mixture framework remains vital, as allows the model to scale its capacity to represent more complex, multi-modal joint distributions when necessary, demonstrated for MIMIC-III and in Figure
C.6 UNCONDITIONAL DENSITY ESTIMATION: FIGURES AND
In both Figures T and 2, all models are trained on unconditional density estimation tasks using either the circle or blast toy datasets. For the time series models, the history set to be empty or constant, and the query points are fixed. While MOSES and ProFITi can model arbitrary distributions, GPR is restricted to predicting multivariate Gaussians.
Figure illustrates the marginalization consistency property for two-variable distributions on both datasets. Both MOSES and GPR are consistent models, meaning their predicted marginals align closely with the marginals obtained from the joint distribution (resulting in very low Marginalization Inconsistency value). In contrast, ProFITi is inconsistent, leading to incoherence between the numerically integrated joint distribution and the predicted marginal of the second variable. This inconsistency arises in the second variable because ProFITi parameterizes the joint distribution using lower triangular matrix, which enforces dependency of the second variable on the first while treating the first variable as independent. Consequently, when predicting the marginal of the second variable, its dependency on the first variable lost, producing a marginal that differs from the true distribution.
Figure 2 shows the advantage of MOSES compared to GMM. GMM requires more number of components compared to MOSES for predicting similar distributions. In this experiment, we use the unconditional version the model meaning the flow architecture independent. Hence, mean and covariance of the Gaussian distributions which are same for both MOSES and GMM are trainable parameters. addition to them, we have the parameters for the separable flows which are not shared among the variables and components. Due to this, MOSES have significantly higher number of parameters compared to GMM. Although one could condition the models upon constant vectors and reduce the number of parameters, did not use this approach and stick the general practice.
22

---

Published as a conference paper at ICLR 2026
-3.222
-2.49
1.474
3.321
2.990
njNLL
3.400
33.20
-2.005
-.479
33.471
-3.744
10
64
128
# #Mxxture Components (D)
Latent Embedding (M)
#How Layers
(a) USHCN
-0.474
00.622
0.608
NINLL
-0.485
00547
-0.616
64
128
Components (D)
Embedding (M)
Layers
(b) PhysioNet2012
-0.217
-0056
1.288
-0.251
-0.18
0.713
n/NUL
-0.265
-0.221
0.138
-0319
0.304
10
16
32
64
128
2 Flow Layers
# Mixture Components (D)
Latent Embedding (M)
MIMIC-III
-1.495
0.547
1.337
-1.576
00.777
N/NLL
11657
11.407
.08222
11.79
11.836
-1.92
64
128
#Mixture Components (D)
Latent Embedding (M)
#Flow Layers
(d) MIMIC-IV
Figure 5: Sensitivity Analysis for MOSES. Model performance with varying number of mixture components, latent embedding size, and flow layers. Evaluation metric: njNLL, lower the better. For each dataset, we fixed optimal hyperparameters and varied one parameter at time.
23

---

Published as a conference paper at ICLR 2026
4
0
0
0
(a) Ground truth Gaussian (b) Distribution generated (c) Distribution generated (d) Distribution generated
Mixture Distribution
by ProFITi
by M1SES(1)
by MOSES(2)
Figure 6: Demonstrating the advantage mixture components
C.7 DEMONSTRATING THE ADVANTAGES OF MIXTURE COMPONENTS IN MOSES
M1SES(1) is restricted to modeling linear dependencies in the joint distribution, as it uses single Gaussian base. While the spline transformations can indeed create complex, multi-modal univariate marginal distributions, they cannot create model complex dependencies between variables.
When we use multiple components (D 1), MOSES gains expressivity as each component can capture different dependency structures. To demonstrate this, we evaluated ture distribution (Figure 6a). As shown in Figure 6b. ProFITi fits this distribution easily due to its non-separable architecture. In contrast, MOSES(1) fails to capture the multi-modalii of the joint distribution Figuue 6c). as its unimodal base cannot represent multiple correlation patterns. However, MOSES with multiple components successfully models the distribution by having different components capture the different modes and their respective dependencies (Figure 6d).
Gaussian mix-
C.8 DEMONSTRATING THE COMPUTATIONAL COSTS OF MOSES
In Figure we compare the computational costs of MOSES and ProFITi across varying forecast horizons on medical datasets. Both models were configured with similar parameter counts (~550K) for fair comparison. While MOSES incurs computational overhead from covariance matrix operations (as analyzed in [Section 4), it demonstrates superior runtime efficiency relative to ProFITi. This efficiency stems from MOSES's s separable flow architecture, where mixture components process in parallel, contrasting with ProFITi's sequential, cascaded flow computations.
Figure 8 further analyzes computational scaling with respect to the number of Gaussian components. As expected, increasing the mixture count D raises computational requirements, though the paral- lelizable nature of these components allows only slight increase in run time. Note that we use different GPU (A40) order run this experiment and hence the results on efficiency is different from Table 11
24

---

Published as a conference paper at ICLR 2026
Peak Memory (GB)
7.67
Peak Merory (GB)
20.82
Peak Memory (GB)
5.00
ProFITi
ProFITi
ProFITi
13.24
11.98
3.17
3.56
2.19
7.69
5.97
6.49
1.80
1,20
0.89
5708
#512
=264
S1336
#906
#454
12-36
24-24
36:12
12:36
24-24
36-12
12-36
24:74
36:12
Configuration (observed-foeecast)
Configuration (observed-forecat))
Configuration (observed-foeccst)
MIMIC-I1I
MIMIC-IV
PhysioNet 2012
(a) Peak memory usage with varying observed and forecast lengths. Max queries (K) indicated below cach bar. Number of components of MOSES set to 10.
Epoch Time [sec)
40.9s
39.65
38.15
Epoch Time [sac)
77.6s
Epoch Time [sec)
18.1s
32.7
64.6s
59.9s
14.5
19.55
9.8s
16.3
12:55
28.5s
23.0c
6.5s
s1336
4459
#297
s147
12-36
24-24
36-12
0.0
12-36
24-24
36-12
12:36
24-24
36-12
Configuration (obses)ed-forecast))
Configuration (observed-forecstt)
Configuration (observed-forecast)
MIMIC-III
MIMIC-IV
PhysioNet 2012
(b) Run time per epoch with varying observed and forecast lengths. Max queries (K) indicated below each bar. Number of components of MOSES set to 10.
Figure 7:
Demonstrating (a) peak memory usage and (b) run time per epoch across the medical datasets with varying observation and forecast lengths.
Peak Memary (G8)
7.78
Paak Memary (GB)
20.92
Peak Peak Memory (GB)
4.90
14.00
11.82
7.54
0.99
0.60
3
10
3
10
3
10
Number of Gaussians
Number of Gaussians
Number of Gaussians
MIMIC-III
MIMIC-IV
PhysioNet 2012
(a) Peak memory usage with varying number of Gaussian components, Observed 12 hours and forecast 36 hours.
23.2
12.0
Epoch Time (sec)
19.3s
Epoch Time (cac)
28.9s
Epoch Time (sac)
10.0s
15.4
14.0s
22.4s
7.45
8.2s
10.5s
11.9s
14.5s
18.65
5.5s
11.6
10.7s
0.0
3
10
10
3
5
10
Number of Gaussians
Number of Gaussians
Number of Gaussians
MIMIC-III
MIMIC-IV
PhysioNet 2012
(b) Run time epoch with varying number of Gaussian components. Observed 12 hours and forecast 36 hours.
Figure Demonstrating (a) peak memory usage and (b) run time per epoch across the medical datasets with varying observation number of components. Observed 12 hours and forecast 36 hours.
25


Reviews:
[{"review_id": "kTP9LcwZoW", "reviewer": "Reviewer_b5pF", "summary": "This paper provide a time series prediction method applied in irregular (missing) dataset. Requirements 1~3 satisfy the problems that irregular dataset facing. The key contribution of R3 is crucial for irregular forecasting because the consistency of probability of one variable if some observations are missed. For implementing abovementioned proposition, MOSES use Separable Normalizing Flows through dealing each variable individually. For forecasting, MOSES train a weight matrix to combine variables. The experiments results suggest MOSES take a large margin than baselines.", "questions": "1.\tAblation experiments may be conducted in ENC to demonstrate its effects.\n2.     As Weaknesses.2 says, is there some method to make probability consistency not only in hidden variable z but also in observation y, or to bridge them to some degree?\n3.\tThough the frequency-time embedding method is employed in MOSES, I'm still concerned about the accuracy degeneration when query-t has a big gap between observe-t, because the query-t must be greater than observe-t which is used in the training stage.\n4.\tThe time series generation model showed impressive performance in irregular datasets, the comparison with them may be interesting. \nLi Y, Lu X, Wang Y, et al. Generative time series forecasting with diffusion, denoise, and disentanglement[J]. Advances in Neural Information Processing Systems, 2022, 35: 23009-23022.\nNarasimhan S S, Agarwal S, Akcin O, et al. Time weaver: A conditional time series generation model[J]. arXiv preprint arXiv:2403.02682, 2024.", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1.\tThe probability consistency is crucial for irregular forecasting.\n2.\tThe pave to accomplish probability consistency is efficient.\n3.\tVerify the importance of probability consistency fully in toy dataset, and MOSESE make a big gap to baseline in four real world dataset.", "weaknesses": "1.\tThough Kolmogorov’s extension theorem make sure the exchangeable of R3, it still need handle different size of subset for ensuring probability consistency.\n2.\tMOSES only ensure the probability consistency of hidden variable z, then bridge this property to y through function f. So the representability of enc might be vital.\n3.\tIn line 257, author claim w depend only on hobs, not on queries. Whereas hobs come from X, and Q is a subset of X. So queries might effect w.", "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": "ynyeqXnc8Q", "reviewer": "Reviewer_1gzJ", "summary": "The paper proposes MOSES, a mixture of separable flows for probabilistic forecasting on irregular multivariate time series. The core idea is to guarantee marginalization consistency by (i) using a consistent Gaussian base with low-rank-plus-identity covariance, (ii) applying per-target invertible spline transforms whose parameters depend only on that target’s query and the shared context, and (iii) mixing such components with query-independent weights. Experiments on several datasets show competitive joint NLL and markedly better marginal NLL and near-zero marginal-inconsistency.", "questions": "- Can you show that your R3 implies full projective consistency for any two finite index sets, and formalize the “conditional on $X$” version of Kolmogorov?\n- Can you provide scaling plots (time/memory) vs ProFITi and others as $K$ and component count $D$ grow.\n- Add other metrics (e.g., CRPS, Energy Score) in the main text to complement NLL.\n- Include ablations: sensitivity to mixture size $D$, latent sizes, and spline bins; discuss the marginal-vs-joint performance trade-off and when query-independent mixing may hurt.", "limitations": "", "rating": 4, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "- Clear formalization of requirements (e.g., marginalization consistency) and an architecture that enforces them by construction.\n- Clean, scalable base: low-rank plus diagonal covariance and diagonal Jacobian flows give tractability while retaining flexible marginal shapes.", "weaknesses": "- Architectural assumption is critical. Consistency hinges on each target’s transform depending only on its own query + context (not other queries). This should be stated as a formal assumption in the lemmas.\n- The low-rank plus covariance shared via a global parameter limits cross-dimensional dependence. The “flow-level expressiveness” claim should acknowledge this trade-off.\n- Complexity claims are theoretical; wall-clock/memory vs. strong baselines at larger $K$ are missing.\n- “Fully identifiable” is too strong.", "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": "gc65bN35GB", "reviewer": "Reviewer_C64z", "summary": "This paper treats the problem of probabilistic forecasting of multivariate irregular time series with missing values. The emphasis is put specifically on obtaining the joint distribution over all covariates at new unobserved query timestamps. The authors identify two main issues with the existing methods: Inconsistent marginals in the case of baselines that directly learn the joint distribution with no consistency guarantess (e.g. ProFitI), and Lack of expressivity of baselines that provide consistent marginals by design (e.g. GPs). In order to fill this gap, the authors propose MOSES, a method that learns a mixture of separable conditional normalizing flows, enabling marginals consistency through separability, and expressivity through the flows' invertible transformations. Experiments in toy datasets and 4 real world datasets show that the proposed method indeed achieves better marginals consistency, trading-off errors in the joint distribution estimation compared to baselines (ProFITi).", "questions": "- Regarding the illustrative examples in Figure 1, although they do a good job in showing the marginal consistency feature of MOSES, I believe that they fail in providing a comprehensive comparison with ProFITi. Specifically, I'm referring to the fact that in these examples, MOSES appear to have a better fit even for the joint density (first line of plots showing the samples). Dont the authors think that it would be better to illustrate a case where ProFITi achieves a better joint density, all while having inconsistent marginals (similar to the results on Physionet'12, MIMIC-III, and MIMIVC-IV)? In this case, the tradeoffs between fitting the joint density and having consistent marginals would be clearer in my opinion. \n- I have a question regarding how the transformation parameters $\\theta_{\\text{flow}}$ are learned, it's mentioned in Appendix A.4 that these parameters are set based on the linear parameters of the LRS model, the knots, their derivatives, etc, and that they are common between variables and mixture components. Doesn't this violates the separability characteristic of the NF in the sense that the transformation $\\phi$ now depends on other query embeddings through the knots?\n- Correct me if I'm wrong but you learn $D$ query embeddings (so $D$ different sets of parameters $\\theta_{\\text{d}}^{\\text{query}}$) then you train everything end-to-end, so I guess in the final likelihood term on $y$, the gradients flowing to a given set of parameters $\\theta_{\\text{d}}^{\\text{query}}$ are coming only from the corresponding mixture terms $p_{\\mathrm{Y}_d}(y|\\mathbf{h}_d)$. If this is correct, how do you prevent these query representations from collapsing into the same representation?", "limitations": "", "rating": 6, "confidence": 4, "soundness": 3, "presentation": 4, "contribution": 3, "strengths": "- The paper is clearly written, highlighting all the important aspects of the method, its comparison with existing baselines, and the experimental results.\n- The approach is theoretically well-motivated, and innovative in imposing a dependency structure in the source distribution of the NF, conserving separability afterwards through component-wise transformations $\\phi$.\n- The experimental results are consistent with the paper's narrative in the sense that they show the proposed approach indeed improves the marginal distributions negative log-likelihood.", "weaknesses": "- My main concern is that unlike the claim in the abstract \"it matches ProFITi in joint prediction performance\", the proposed method achieves statistically-significant worse joint negative log-likelihood than ProFITi in 3 out of 4 real-world tasks considered in the paper (referring to Table 1). Furthermore, in the only dataset where it outperforms ProFITi, it's shown in the ablation study in the appendix (Table 11) that the mixture components ($D>1$) are important for this dataset, suggesting the improvement could be related to the mixture rather than the marginals consistent design of MOSES. The last matter make me wonder whether ProFITi can be augmented to fit a mixture joint density, in which case it could maybe win also in the 4th dataset (and perhaps in univariate marginals as well).", "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": "WVnEI82zsx", "reviewer": "Reviewer_8T7r", "summary": "This paper proposes a model named MOSEs, designed to achieve marginalization consistency in probabilistic irregular time series forecasting. Formally, given an input $X$ and a full query set $Q$, marginalization consistency requires that for any subset $S \\subseteq Q$, $p(y_{S} \\mid S, X) = \\int p(y \\mid Q, X), dy_{Q \\setminus S}$ i.e., the marginal distribution obtained by direct querying should be consistent with that obtained through marginalization. MOSEs enforces this property using a mixture of separable normalizing flows, where each flow is constructed from a parameterized source Gaussian distribution followed by separable univariate transformations. Multiple flows are then combined with softmax-weighted mixing coefficients. Experimental results demonstrate that MOSEs achieves superior marginalization consistency compared with baseline models.", "questions": "1. To further evaluate the expressiveness of separable flows, please compare **MOSES(1)** with ProFITi on the following toy example: a 2D Gaussian mixture distribution with two components—one centered at $(1, 1)$ and the other at $(-1, -1)$. In this case, $p(x_1, x_2) \\ne p(x_1)p(x_2)$, making it suitable to test whether the separability assumption limits modeling capacity.\n\n2. Since ProFITi performs better in terms of njNLL, is it possible to marginalize its joint distribution to obtain the corresponding marginal distributions? If so, how does this method perform on njNLL compared with MOSES?", "limitations": "", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "1. The issue studied in this paper—marginalization consistency—is crucial for real-world probabilistic forecasting applications.\n\n2. The paper is well-structured and easy to follow.\n\n3. The proposed model theoretically guarantees marginalization consistency in closed form.", "weaknesses": "1. The separability constraint (Equation 8) is overly restrictive, substantially limiting the flexibility of the normalizing flow. A flow defined as $f(x_1, x_2) = (f(x_1), f(x_2))$ effectively models only the marginal distributions and implicitly assumes $p(x_1, x_2) \\ne p(x_1)p(x_2)$. Although the use of parameterized source Gaussians and mixture modeling can partially mitigate this issue, the restriction still reduces the expressive power that typically makes flow-based models advantageous.\n\n2. Consequently, the prediction accuracy (as measured by njNLL) is inferior to that of non-separable flow models. This limitation is further supported by the results shown in the Figure 2, where the GMM achieves better approximations, despite MOSEs employing more complex neural architectures.", "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)