Generating Full-field Evolution of Physical Dynamics from Irregular Sparse Observations

We thank the reviewer for the careful review and will improve the clarity in the latest version. We address the comments below:(C: comment; R: response)

To aid the reviewer’s understanding, we first claim our paper’s main objective, basic concepts, and experiment setup:

C1: Confusion on more general ML relevance

R1: We propose a new generative model that directly operates on off-grid multidimensional data during training and testing phase, which is broadly applicable to unstructured domains like sensor networks, point clouds, and scientific measurements.

C2: Concerns on (i_n_m, t_m, y_i_n_m) and K-mode notation

R2:We clarify that this work aims to construct a generative model directly from off-grid data, which cannot be represented as regular matrices or tensors. Instead, we treat the data as points, following common practice in the literature [12, 31].

This work builds a generative model from off-grid data, treated as point sets rather than regular tensors, following common practice [12, 31].

The value of $K$ is for notational rigor and is will be greater than 3 in practice—for example, $K = 4$ for 3D spatial + frequency radio map/seismic data, or $K = 5$ for multispectral data with 3D spatial, frequency and observation angular modes. Our experiments focus on 2D and 3D cases for simplicity.

C3: Question on the time mode

R3:The temporal mode is not included in $K$. The reason is that the temporal mode usually exhibits much sharper or more complex fluctuations that cannot be effectively captured by a single decoupled latent factor function [31].Instead, we treat it separately using our proposed GPSD module, which leverages Gaussian processes and a temporally-augmented UNet to model temporal correlations.

C4:What does “solving” eq 6 mean?What is “alternating-direction strategy”?

R4: Solving Eq. (6) means learning the mapping between data and a low-dimensional latent space by minimizing Eq. (8) using an alternating optimization that updates $\mathcal{W}$ and $f_{\theta}$ parameters in turn.

C5:The remark about VAE is rather strange

R5: We clarify that FTM is similar to VAEs in function—mapping data to a latent space—rather than in specific techniques.

C6: Visualisation for denoisng cores.

R6: We will add visualizations in the latest version.

C7: Ablations on Tucker vs other latent representations; message-passing vs other temporal consistency guidance methods.

R7:We use FTM as the latent representation to learn a mapping from irregular multidimensional physical fields to a lower-dimensional space, as other latent representations like VAE cannot handle such irregular data in our setting.

To the best of our knowledge, we are the first to introduce temporal consistency for diffusion guidance. However, we can use common interpolation methods like linear or spline to interpolate gradients for those timesteps without observations. We conduct experiments on Observation setting 2 with three datasets. The results are:

One can see that our MPDPS performs much better. We will supplement these into the latest version.

C8: Missunderstanding on factor function

R8: We kindly refer you to line 146, where we clarify the latent function as a vector-valued function. The learned latent functions can be interpreted as positional bases functions for each mode, while the core tensor serves as the weighting coefficients that interact with these bases across the $K$ modes.

C9: Interpretion of core

R9:The core tensor serves as weighting coefficients for the $K$-mode bases, summarizing the full field at each snapshot. Pls see Sec.6 of [22] for details of its structural patterns.

C10:Concern on identifiability, rotational invariance, and sensitivity to initialisation

R10:

• FTM does not produce identifiable $\mathcal{W}$ and $f$ due to the freedom of invertible transformations, as you noted, but this is not an issue since our reconstruction task—unlike blind source separation—does not require identifiability. Like deep matrix factorization and VAEs, these transformations leave the reconstruction loss unchanged. Once $f$ is learned and fixed, the rotational ambiguity of $\mathcal{W}$ is also fixed during diffusion training and sampling.

• For the sensitivity to initialisation, we acknowledge the initialisation is crucial for model training, just like most deep-based models' training. In practice, we train with Kaiming initialisation, AdamW, and (mini‑)batch gradient descent on the convex quadratic decoder of Eq. (6).Across 5 random seeds the VRMSE standard deviation is < $10^{-3}$, showing that the model converges stably regardless of the initial values.We will add this table to the appendix in the latest version.

C11:question about meaning of “drawn from a continuous function $\mathcal W(t)$ and stochastic process interpretation”

R11:

• On the function view, we regard the core trajectory as a continuous map $\mathcal W:\mathbb R^{+}\to\mathbb R^{R_1\times\cdots\times R_K}$. meaning that it can be evaluated at any time point $t$ within the range of interest. The observed cores $\mathcal W_{t_m}$ are simply samples of this function at discrete time stamps, hence the phrase “drawn”.
• On the stochastic‑process view, in the generative (Stage‑2) part we place an element‑wise Gaussian‑process (GP) prior on $\mathcal W(t)$ and inject GP noise in the forward diffusion. The reverse process denoises these stochastic samples. Thus, in training, we learn the GP prior and the denoiser to model the temporal correlations of the core sequence. This is similar to how standard diffusion models learn the distribution of images or videos in pixel space, but here we do it in a latent function space. In the sampling (Stage‑3) part, we sample the core sequence from the learned GP prior and denoise it to obtain the final core sequence.

C12:Confusion about the componants of diffusion model (s,t,i)

R12: s:diffusion step; t:time step. In Stage 1, we use FTM to map all data points at time step $t$ into a Tucker core $\mathcal{W}_t$. The spatial information is captured by latent factor functions and is decoupled from $\mathcal{W}_t$. Thus, $\mathcal{W}_t$ summarizes the full-field information at time $t$ and only the temporal mode needs to be modeled during diffusion training.

C13: Why do we want to predict target times if we already observe them?

R13: Our main task is to reconstruct time-evolving full fields from sensor observations. Although the target time points have some sparse observations, they are far from enough for downstream tasks. Hence, completing the entire field is necessary.

C14: I can’t follow eqs 10-12

R14: For detailed explanations, pls refer to R5 of Reviewer wHTh.

C15: Questions on target/test/obs/train

R15: Pls refer to the experiment setup at the begining for the train/test split. In the test phase, the target time refers to the timesteps where we aim to reconstruct the full field, while the observation time refers to the timesteps where sparse observations are available.

Given that the score is slightly below the borderline, we would like to know what other steps we can take to further improve our paper or address your concerns？

We sincerely appreciate your thoughtful feedback and recognition of our work. Below, we address the comments (C: comment, R: response):

C1:I acknowledge that I am not deeply familiar with the literature on tensor-based representations for physical field modeling, so I may have missed potential limitations in that aspect of the method.

R1: Thank you for your comment! We would like to note that using tensors to model multidimensional physical fields is still a new research direction with little prior work. Tensor-based representations have been widely used in image processing. However, to the best of our knowledge, we are the first to apply them to generative modeling of physical fields and our extensive experiments deomstrate their effectiveness.

C2: As the number of observations decreases, is there a threshold below which the model fails to generate meaningful reconstructions due to a lack of information?

R2: Thank you for the question. Diffusion models generate data based on two main sources of information: (1) the prior distribution learned during pre-training, and (2) the likelihood distribution derived from observed data. As a generative model, a well-trained diffusion model can always produce meaningful outputs even when there are few or no observations, relying mostly on the learned prior when the likelihood information is limited. To illustrated that, we plot the unconditional generation results in in Fig. 6,7,8 of Appendix F. One can see that even without observations, our method still generate meaningful trajectories.

C3:How does the model perform under non-uniform observation patterns, such as those common in weather-related applications, where sensors are densely located in developed regions but sparse over oceans or remote areas? This setting might challenge the assumption of homogeneity in sampling patterns.

R3: Thank you for the great question! Following your suggestion, we equally divide the physical field into two regions. Assuming an observation rate of 1% at each timestep, we simulate non-uniform observation patterns by allocating 60% of the observation points to one region and 40% to the other. We also test a more imbalanced case, with 80% of the observations in one region and 20% in the other. These experiments are conducted on a single sample from each dataset. We report the reconstruction VRMSE below:

We find that our method remains robust under non-uniform observation patterns in different datasets, highlighting its practical usefulness. We will supplement it in the latest version.

C4: Could you clarify the main advantages of the method compared to such approaches, particularly in terms of flexibility, generalizability, or computational efficiency?

R4: We appreciate the reviewer’s insightful comment regarding the comparison with domain-specific methods such as [1], which perform field reconstruction without pretraining. Below, we highlight the advantages of our SDIFT framework in terms of flexibility, generalizability, and computational efficiency:

• SDIFT is able to handle irregular, off-grid observations over continuous spatiotemporal domains. Leveraging the FTM with the universal approximation property, SDIFT flexibly represents high-dimensional fields into latent space without requiring explicit ODE/PDE formulations. Then, in latent space, SDIFT use GPSD to flexibly model temporally coherent core sequences at arbitrary timesteps. However, [1] relies on physics-informed architectures with predefined dynamics and structured observations (e.g., periodic sampling).

• SDIFT generalizes across diverse physical systems and observation scenarios through two key components: FTM and MPDPS. FTM constructs high-quality mappings between irregular multidimensional physical data and a lower-dimensional latent space. MPDPS effectively propagates sparse observation information across the entire core sequence, ensuring robust performance with irregular, sparse, or incomplete data. In contrast, [1] relies on gradient-based solvers that require structured observations and degrade under irregularities.

• SDIFT achieves high sampling efficiency by operating in a low-dimensional latent space. The FTM inherents a linear structure between latent factor functions and core tensors, resulting in a linear likelihood model. This design makes SDIFT scalable to large, high-dimensional systems. In contrast, [1] requires extensive iterative gradient computations over full state vectors, which becomes computationally prohibitive at scale.

• Regarding physical constraints: We agree that [1] enforces physics through explicit modeling, but this comes at the cost of flexibility and generalizability. In contrast, SDIFT implicitly captures physical regularities through FTM’s low-rank spatial structure, spatial smoothness, and GP-based temporal continuity. This design preserves flexibility while allowing the future incorporation of domain-specific priors—striking a balance between physical fidelity and adaptability to under-constrained, sparse data.

We thank the reviewer for the careful review! We address the comments below:( C: comment; R: response)

C1: The comparison with baselines is not comprehensive. Although the authors compare parameter size and runtime in Table 2, this comparison is limited to CoNFiLD only and does not include other baselines. I am particularly curious about how the proposed method performs against the other reported results.

R1: Thank you for the comment. We compared our method only with CoNFiLD because both are generative approaches. The parameter sizes and runtimes for the other methods are reported below.

The parameter sizes of different models vary because they are different types of models and we fine-tune each model individually to achieve its best performance, taking into account its unique properties.

We observe that tensor-based methods (NONFAT, DEMOTE, and LRTFR) require significantly more time to produce reconstructions. This is because they directly minimize a loss function that enforces consistency between predicted and observed values, without any pretraining phase. As a result, these methods need extensive iterations to converge and are much slower than pretrained models. Among them, LRTFR is considerably faster than NONFAT and DEMOTE, as the latter two are not scalable to large tensors.

Senseiver is the most efficient model in terms of runtime, as it feeds the indices to be queried into the network and performs the forward pass only once to produce the reconstructions. However, it does not explicitly model the temporal evolution of the physical field. As a result, its performance is inferior to ours and it fails when no observations are available.

C2: The reviewer is confused about the content in Section 3.1.

R2: Thank you for your question. In your example, $y(i\_1, i\_2) = e^{i\_1 \cdot i\_2}$, we can construct a finite-dimensional FTM solution that exactly represents $y$. That is $y(i\_1,i\_2)= \text{vec}(\mathbf{W})^{\text{T}}(f\_{\theta\_1}^1(i\_1) \otimes f\_{\theta\_2}^2(i\_2))$, where $\mathbf{W}$ being a scaler 1 and $f\_{\theta\_1}^1(\cdot)=e^{i\_1}: \mathbb{R}^1 \to \mathbb{R}^1$,$f\_{\theta\_2}^2(\cdot)=e^{i\_2}: \mathbb{R}^1 \to \mathbb{R}^1$ being scaler functions. Then the FTM model becomes $1\cdot(f\_{\theta\_1}^1(i\_1) \cdot f\_{\theta\_2}^2(i\_2))$. Since both latent functions are parameterized by neural networks, they can accurately approximate the exponential functions involved.

For illustration, let's consider a finite 2-order tensor(matrix) $\mathbf{Y}$ with size $3 \times 3$ with each entry $y(i\_1,i\_2)=e^{i\_1,i\_2}$. Then $\mathbf{Y} =$

=

FTM decomposes the tensor $\mathbf{Y}$ into a Tucker core and mode-wise latent functions.
Assume $R\_1 = R\_2 = 1$, then the Tucker core is a learnable matrix $\mathbf{W} \in \mathbb{R}^{1}$, and the latent functions are defined as:

• Mode-1: $f\_{\theta\_1}^1(\cdot): \mathbb{R}^1 \to \mathbb{R}^1$
• Mode-2: $f\_{\theta\_2}^2(\cdot): \mathbb{R}^1 \to \mathbb{R}^1$

Here, $f\_{\theta\_1}^1(\cdot)$ and $f\_{\theta\_2}^2(\cdot)$ are neural networks that learn positional embeddings and are parameterized by $\theta\_1$ and $\theta\_2$.

To approximate the value at index $(1, 1)$, FTM computes: $\text{vec}(\mathbf{W})^\top \left( f\_{\theta\_1}^1(1) \otimes f\_{\theta\_2}^2(1) \right) \approx y\_{11} = e^{1}$.

Other entries in $\mathbf{Y}$ are computed similarly.

FTM minimizes the reconstruction error between $\mathbf{Y}$ and the model output, learning optimal values for $\mathbf{W}$, $f\_{\theta\_1}^1(\cdot)$, and $f\_{\theta\_2}^2(\cdot)$. It still works even when some entries in $\mathbf{Y}$ are missing. After training, $\left( \mathbf{W}, f\_{\theta\_1}^1(\cdot), f\_{\theta\_2}^2(\cdot) \right)$ can be used to represent the original tensor $\mathbf{Y}$. Therefore, FTM serves as an encoder.

C3: Confusion on the case $R\_k \to \infty$ and last step of Equation (24).

R3: Thank you for the question. We would like to clarify that the universal approximation proof is based on sufficiently large but finite $R\_k$, rather than taking $R\_k \to \infty$.

As shown in [A1], the Tucker model is a universal representer: for sufficiently large $\{R\_k\}$—up to $R\_k = I\_k$, where $I\_k$ is the size of mode $k$—any finite-dimensional tensor can be exactly represented via Tucker decomposition.

For example, consider the above case where $\mathbf{Y} \in \mathbb{R}^{3\times 3}$, $I\_1 = I\_2 = 3$. If we choose $R\_1 = R\_2 = 3$, i.e., $\mathbf{W} \in \mathbb{R}^{3 \times 3}$, then the model

$\text{vec}(\mathbf{W})^\top \left( f\_{\theta\_1}^1(\cdot) \otimes f\_{\theta\_2}^2(\cdot) \right)$

can accurately represent $\mathbf{Y}$. Therefore, $R\_k$ is a finite and bounded value, and the expression in Eq. (24) consists of finite summations.

Moreover, natural physical fields typically exhibit low-rank structure, as discussed in [20,21], implying that the effective $R_k$ values are often small in practice.

[A1] N. D. Sidiropoulos, L. De Lathauwer, X. Fu, K. Huang, E. E. Papalexakis, and C. Faloutsos, “Tensor decomposition for signal processing and machine learning,” IEEE Trans. Signal Process., vol. 65, no. 13, pp. 3551–3582, Jul. 2017.

C4:There is a typographical error in the formula on page 14, Equation (18)

R4: Thank you for pointing out this error. We will fix it in the latest version.

C5:Minor issue: references [28] and [38] refer to the same publication.

R5: Thank you for pointing out this issue. We will address it in the latest version.

R1: Thanks for the feedback!

R2: Thank you for pointing this out. We sincerely apologize for our misunderstanding. We claim that there is no finite-dimensional separable expansion solution to represent $e^{i \times j}$ with FTM.

R3: Thank you for your thoughtful comment and clear proof! We agree that $e^{xy}$ does not admit an exact finite-dimensional separable expansion in the FTM form. However, this does not contradict our proposed Universal Approximation Property (UAP) theorem.

Our UAP states that FTM can approximate any continuous function to arbitrary precision (i.e., for any $\epsilon > 0$, where $\epsilon$ is controlled by the choice of bounded values $R_1$ and $R_2$). This does not imply that FTM always yields an exact finite-dimensional separable representation for every continuous function. Instead, it focuses on characterizing the asymptotic approximation error.
Typically, exact separable representations are not required in practical applications, as approximate forms are sufficient to achieve reliable performance, as demonstrated in our experiments.

Additionally, we would like to argue that if the dimension $M$ is not infinite (i.e., bounded), we can always choose a sufficiently large $N$ (e.g., $N = M$) such that FTM can represent the matrix $\mathbf{M}$ exactly (i.e., $\epsilon = 0$), though the resulting form is not exactly separable expansion.

We sincerely appreciate your careful checking and constructive feedback and we will make this point explicit in the main text!

Thank you for your feedback! We address the comments below:( C: comment; R: response)

R0: We sincerely thank the reviewer for the insightful comments. Regarding the relationship between our method and prior tensor-decomposition work (e.g., Chen et al. [1] and Chen et al. [2]), we agree that tensor factorization has proved powerful for spatiotemporal modeling. Our contribution, however, targets a different problem setting: generating full-field evolutions of physical dynamics from sparse, noisy, and continuously indexed observations. The differences between our method with these existing tensor-based methods are:

• (i) Existing tensor-based approaches are typically trained to deterministically reconstruct or extrapolate a single trajectory on a regular grid; They often struggle to perform reliably under extremely sparse observations. However, SDIFT learns a full generative distribution over sparse, irregular spatiotemporal coordinates and it can always produce meaningful off-grid results.

• (ii) Instead of directly factorizing the observed tensor like [1][2], we use a Functional Tucker latent space, where coordinates are embedded via continuous functions alongside a learnable, time-varying core. This framework defines a flexible manifold that GPSD can effectively model.

• (iii) The time evolutions in Chen et al. [1] are based on a linear dynamic mode decomposition assumption, while Chen et al. [2] relies on decoupled mode factors. These assumptions limit their ability to represent complex dynamics. In contrast, our GP-driven sequential diffusion framework effectively captures nonlinear dynamics and provides predictive uncertainty.

• (iv) SDIFT enables conditional generations of dynamics at arbitrary continuous time steps with the proposed MPDPS mechanism whereas Chen et al. [1] and Chen et al. [2] do not.

To better clarify our positioning, we will reference both works and include a dedicated subsection in the related work discussing differences in task formulation, sampling strategies, and uncertainty modeling.

Concerning autoregressive (AR) alternatives, we acknowledge their popularity in image and video generation (i.e., MaskGIT[B1], diffusion forcing [B2]) and have already included Senseiver—an AR attention model—as a key baseline. SDIFT reduces its VRMSE by 22 % on the same sparsity level, and and Senseiver fails entirely in the observation setting 2. The gap between our method and AR approaches arises from four key factors:

• (i) Irregular, continuous timestamps force AR models to discretise Δt, which introduces interpolation bias.
• (ii) Extremely sparse spatial observations make pixel-space AR training infeasible.
• (iii) Sampling flexibility: SDIFT solves an ODE to “jump” directly to any target time, while AR models must generate frames sequentially.
• (iv) Inference strategy: our message-passing posterior sampling integrates limited observations into the full sequence in a single pass, something AR models cannot consistently achieve. Thus, we claim it's non-trival to apply the AR directly in our problem setting.

Therefore, we argue that directly applying existing autoregressive (AR) models designed for computer vision (CV) or natural language processing (NLP) to our problem setting is non-trivial. In response to the reviewer’s request, we will note in the Limitations section that, when observations are dense and uniformly sampled, a simple AR backbone may indeed offer a more efficient solution. We will also expand the discussion to consider future work on hybrid AR–diffusion frameworks.

Thank you again for your valuable feedback!

[B1] Chang, Huiwen, et al. "Maskgit: Masked generative image transformer." Proceedings of the IEEE/CVF conference on computer vision and pattern recognition. 2022.

[B2] Chen, Boyuan, et al. "Diffusion forcing: Next-token prediction meets full-sequence diffusion." Advances in Neural Information Processing Systems 37 (2024): 24081-24125.

Thank you for your feedback. However, we beg to differ with the reviewer’s comment that our work is “incremental compared with existing tensor learning methods.”

We would like to clarify that the Functional Tucker Model (FTM) is used solely as a latent representation mapping in our framework. The directly using of functional tensor is only a fraction of main contributions of our method. Beyond that, we provide a theoretical proof of the universal approximation property of FTM in this context firstly, to support its applicability.

What's more, our primary contributions lie in the design of GPSD and MPDPS on the side of generative methods, which together enable to infer entire physical fields at arbitrary coordinates and timesteps in a generative manner. Such setting and paradigm are practically significant, technically underexplored, and fundamentally different compared with any existing tensor decomposition methods.

Thank you for your strong support! Here are our responses. C: comments; R: response.

C1:in section 3 124, for the diffusion prior pre-training phase, what's the definition of homogenous in "B batches of homogeneous dynamics" mean?

R1: Thank you for your question. In generative modeling, people assume that the training data is extracted from a specific distribution. Here, we also assume the training data is extracted from a specific distribution of dynamics. The homogeneous dynamics refer to the same type of dynamics, such as the same PDE or ODE, but with different initial conditions.

C2:in section 3.1, line 166, in equation 8 , do you need a coefficient in front of the TV term in practice? If so, how do you choose it?

R2: Good point! We add a coefficient in front of the TV term. We have done cross-validation experiments to choose the coefficient and we find that $10^{-6}, 10^{-7}, 10^{-8}$ all work. We will clarify it in the latest version.

C3:section 3.2, line 192, equation 9, there is a denoiser-matching objective to train the diffusion model for fixed W. Have you tried training (8) and (9) simultaneously rather than in pieces. Any +/-'s to report?

R3:Good question. Theoretically, (8) and (9) can be trained simultaneously. However, the objectives of FTM and GPSD differ significantly:

• FTM aims to learn a high‑fidelity mapping between the data domain and latent space, minimizing reconstruction error.
• GPSD aims to model the distribution of latent‑space dynamics over time.

Jointly optimizing these distinct objectives leads to slower training and reduced stability.

To address this, we adopt a two‑stage training pipeline：first training the FTM, then freezing it and training the GPSD. This strategy is widely used in the computer vision literature [A1], which has been shown to yield superior generation quality, faster sampling, and more stable convergence compared to end‑to‑end training.

[A1] Rombach, Robin, et al. "High-resolution image synthesis with latent diffusion models." Proceedings of the IEEE/CVF conference on computer vision and pattern recognition. 2022.

C4:does temporally augmented Unet just mean that the Unet gets the noised ($\mathbf{W}\_{t}^{s}$ for data time t at diffusion time s) as well the times $t=[t\_1,\cdots,t\_M]$ themselves? If not, what does it mean?

R4: The answer is yes. Apart from that, we also use Conv1D layers to capture temporal correlations among $M$ frames. The detailed archetecture is provided in Appendix B.

C5:Main question: What is proved that the specific message passing algorithm introduced, which uses $\nabla\_{W\_{\setminus t\_l}^s}\log p(O\_{t\_l}|W\_{\setminus t\_l}^s)$ . I get intuitively that we are passing some info from the observed to the rest, but abstracting a bit, how can I see more easily that updating the whole $\mathbf{W}\_t$ according to these gradients produces a sample of $\mathbf{W}\_t^{\text{all}}|\mathbf{W}\_t^{\text{obs}}$.

R5: Thank you for your question! We want to highlight that the primary goal of MPDPS is to exploit the temporal continuity of the core sequence to smoothly propagate guidance from limited observations across the entire sequence.

For clarity, consider the setup in Fig. 2, where we aim to generate the core tensors at three target timesteps: $t\_l$, $t'$, and $t\_{l+1}$, denoted as $\mathcal{W}\_{t\_l}$, $\mathcal{W}\_{t'}$, and $\mathcal{W}\_{t\_{l+1}}$, respectively. In this example, observations are only available at $t\_l$ and $t\_{l+1}$ (observation timesteps),denoted as $\mathcal{O}\_{t\_l}$ and $\mathcal{O}\_{t\_{l+1}}$. Standard DPS methods cannot provide gradient guidance for generating the intermediate core $\mathcal{W}\_{t'}$, due to the lack of direct observations at $t'$.

MPDPS addresses this limitation by propagate the information of $\mathcal{O}\_{t\_l}$ to {$\mathcal{W}\_{t\_l}, \mathcal{W}\_{t'},\mathcal{W}\_{t\_{l+1}}$} and the information of $\mathcal{O}\_{t\_{l+1}}$ to {$\mathcal{W}\_{t\_l}, \mathcal{W}\_{t'},\mathcal{W}\_{t\_{l+1}}$}, which means each observation set will contribute to the generation of all target sequences.

Without loss of generality, we take $\mathcal{O}\_{t\_l}$ as an example. At diffusion step $s$, we use $\mathcal{O}\_{t\_l}$ twice to obtain gradient guidance for the three cores:

• First use to compute gradient for $\mathcal{W}\_{t\_l}^s$: This follows standard DPS. We directly compute the likelihood gradient with respect to $\mathcal{W}\_{t\_l}^s$ as $\nabla\_{\mathcal{W}\_{t\_l}^s} \log p(\mathcal{O}\_{t\_l} | \mathcal{W}\_{t\_l}^s)$, as illustrated by the green& red arrows at $t\_l$ in Fig.2(a).

• Second use to compute gradients for $\mathcal{W}\_{t'}^s$ and $\mathcal{W}\_{t\_{l+1}}^s$:
  It includes three steps:

  1. Compute denoised cores at $t'$ and $t\_{l+1}$: we first use $\mathcal{W}\_{t'}^s$ and $\mathcal{W}\_{t\_{l+1}}^s$ to estimate the corrsponding denoised cores $\hat{\mathcal{W}}\_{t'}^0 = D\_{\theta}(\mathcal{W}\_{t'}^s)$ and $\hat{\mathcal{W}}\_{t\_{l+1}}^0 = D\_{\theta}(\mathcal{W}\_{t+1}^s)$ via the pretrained denoiser $D\_{\theta}(\cdot)$.

  2. Estimate pseudo cores at $t\_l$ via GPR: With the temporal continuity assumption, we use Gaussian Process Regression (GPR) to estimate the pseudo denoised core $\hat{\mathcal{W}}\_{t\_l}^0 = GPR (\hat{\mathcal{W}}\_{t'}^0, \hat{\mathcal{W}}\_{t\_{l+1}}^0) = GPR(D\_{\theta}(\mathcal{W}\_{t'}^s), D\_{\theta}(\mathcal{W}\_{t+1}^s))$ , which can be viewed as a regression-based prediction using two denoised cores.

  3. Gradient computation: With pseudo core $\hat{\mathcal{W}}\_{t\_l}^0$ $\quad$, we compute an additional likelihood of $\mathcal{O}\_{t\_l}$. Then we can then compute the gradient of the likelihood respect to ${\mathcal{W}}\_{t'}^s$ and ${\mathcal{W}}\_{t\_{l+1}}^s$, as they directly determine $\hat{\mathcal{W}}\_{t\_l}^0$ via $D\_{\theta}$ and GPR. We further derive a closed-form solution in Appendix C for gradients with respect to $\mathcal{W}\_{t'}^s$ and $\mathcal{W}\_{t\_{l+1}}^s$ (i.e., $\nabla\_{W\_{\setminus t\_l}^s}\log p(O\_{t\_l}|W\_{\setminus t\_l}^s)$ where $W\_{\setminus t\_l}^s=${$\mathcal{W}\_{t'}^s,\mathcal{W}\_{t\_{l+1}}^s$}), which can be efficiently computed. The results are shown in Eq.(10)(11). The process is illustrated by the green arrows pointing to $\mathcal{O}\_{t\_l}$ in Fig. 2(b).

Similarly, information from $\mathcal{O}\_{t\_{l+1}}$ is propagated to $\mathcal{W}\_{t\_l}$ 、$\mathcal{W}\_{t'}$, and $\mathcal{W}\_{t\_{l+1}}$ in the same manner, and we sum the gradients from all observations to obtain the final gradient guidance for each core, as shown in Eq. (12).

In this way, each observation contributes to the generation of all target timesteps. In other words, every target core—regardless of whether it has a corresponding observation—receives information from all observed timesteps. To sample the entire core sequence, we simply add the gradients from all observations to the posterior, as shown in Eq. (12). This mechanism encourages smooth and coherent generation across the entire sequence.

We hope this explanation clarifies how MPDPS updates the entire core sequence. We will include this in the latest version.