Heavy-Tailed Diffusion Models

We thank the reviewer for their comments and address each point in more detail below.

Weaknesses/Questions

Weakness 1a: The addition I would be most keen to see is a more thorough discussion of the results in Tables 2 and 3 (empirical results in the unconditional setting). For example, it would appear from Table that in the test setting for VIL, their method actually underperforms compared to EDM (+PCP) and compared to EDM (+INC) in the other setting, and the authors do not discuss this in the text as far as I can tell. (Though, Figure 3 tells a somewhat different story, and I can't work out how to reconcile the two - see question below.)

Response: We provide a detailed response for this point in our response to Question 2 below.

Weakness 1b: The authors also don't really discuss the results in Table 3 in the main body, and a short commend about what we should take from those results might be appropriate.

Response: We have added a short clarification on the results in Table 3 in the main text (lines 409-411).

Question 1: How do we tune the parameter \nu in practice? It looks from Table 3 like the choice matters a lot, so it might be nice to say something about how one should choose a suitable value of \nu or know that one has found a good value.

Response: In this work, we tune $\nu$ manually for different datasets. However, we also highlight several principled ways to learn $\nu$ directly from the data.

Firstly, it is worth noting that the parameter $\nu$ directly influences training in heavy-tailed diffusions. This is because the noising process $x_t = \mu_t x_0 + \frac{z}{\sqrt{\kappa}}$ where $\kappa \sim \chi^2(\nu)/nu$. Since the Chi-squared distribution is a special case of a Gamma distribution, we can perhaps use reparameterization [1] to estimate the parameters (i.e., $\nu$) of this distribution using a neural network and backprop end-to-end to learn this parameter. This would correspond to learning $\nu$ directly from the data.

[1] The Generalized Reparameterization Gradient, Ruiz et al.

Alternatively, one can also estimate $\nu$ from the data in a separate stage and use it as a perturbation kernel parameter. We highlight one simple way to do this in Appendix D, where we attempt to derive the optimal noise schedule for EDM and t-EDM by minimizing the mutual information between the $x_0$ and the noisy data point $x_\sigma$ w.r.t $\sigma$. In principle, for t-EDM, we can also solve this minimization problem jointly over $\nu$ and $\sigma$ to learn an initial $\nu$ from the data.

We have added these points in Appendix F.2 in our revision. Overall, our proposed framework can be readily extended to learn $\nu$ in an end-to-end or separate stage fashion. Nevertheless, we leave this exploration to future work due to space constraints.

Question 2: How confident are you that the 1-dimensional summary statistics you are looking in Table 2 (and maybe Table 3 also) at appropriate? From looking at Figure 3, t-EDM seems clearly qualitatively superior to the other methods. The fact that this is not reflected in the summary statistics in Table 2 makes me suspect the summary statistics are not capturing the most relevant aspects of the distributions.

Response: This is an interesting observation. We elaborate more upon this aspect as follows.

Firstly, our method exhibits better performance for the KS metric. This could be due to the fact that the KS metric, in addition to flattening the data, is computed over 99.9th and 0.1th percentiles for the right and left tails, respectively. This implies that the KS metric specifically examines the tail values and thus represents an appropriate metric to capture the difference between the tail estimation for the generated and the data samples. In contrast, the KR and SR metrics are based on computing higher-order moments over the entire domain of flattened empirical samples and thus might be more coarse in nature, hence the discrepancy.

However, it is worth noting that while we present EDM baselines with improved preconditioning and preprocessing as additional baselines, these can also be applied to t-EDM, potentially improving results. Therefore, a more fair comparison should be between t-EDM and the baseline EDM for which t-EDM outperforms EDM by a significant margin in all metrics. We added this observation in Section 4.1.

Secondly, due to the flattening operation, the metrics considered in this work disregard the spatial structure of the generated samples, which is undesirable. The qualitative histogram plots in Figure 3 help alleviate this to some extent. However, designing metrics that can quantify heavy-tailed behavior while considering the spatial structure of the generated samples is an interesting direction for further work but is outside the scope of this paper. We discuss these limitations in more detail in Appendix F.2.

Weakness 3c: Moreover, the claims regarding t-EDM's performance appear somewhat overstated, given the tables and the likely margin of error.

Response: We would like to highlight that t-EDM largely outperforms the baseline EDM (Row 1) model across all metrics in Table 1. Moreover, our qualitative histogram comparisons in Fig. 3 illustrate that our proposed method outperforms EDM and its improved versions in modeling heavy-tailed variables, as also acknowledged by Reviewers uzb2 and zbkB. We agree that with improved preconditioning (PCP) and preprocessing (INC), EDM outperforms t-EDM in some cases on the KR and SR metrics. However, for the KS metric, which specifically examines the left and right tails by filtering 0.1 and 99.9 percentiles of the data (and is thus arguably a better metric for comparing tail estimation), t-EDM outperforms all other baselines.

Lastly, we would like to highlight that these techniques can also be readily combined with the t-EDM model and should yield similar improvements. Therefore, the main point of comparison in Table 1 should be primarily between t-EDM and the standard baseline EDM (row 1). The current claim in lines 82-83 is based on this observation. We have added this caveat in the main text in Section 4.1. In this context, we feel that we do not overclaim our contributions in this work.

Weakness 3d: In my opinion, the experiments should include toy examples for a deeper comparison with existing diffusion models.

Response: We would like to point the reviewer to Fig. 1, where we compare the performance of the baseline EDM with t-EDM for different $\nu$ values. Fig. 1 highlights that t-EDM provides a parameter $\nu$, which can help control tail estimation and is thus better able to model heavy-tailed distributions as compared to other baselines like standard Gaussian diffusion. We also include additional quantitative results in Table 8 for the same 2-D toy dataset in Appendix C.3 (section “Quantitative Results on Neals Funnel”). From Table 8, our quantitative results support our qualitative results, implying that our proposed model can definitively help with heavy-tailed data modeling. We also include the results from Table 8 below for convenience.

The quantitative results in the above table for the toy dataset support our qualitative results in Figure 1 in the main text and imply that our proposed model can definitively help with heavy-tailed data modeling as compared to baselines like Gaussian Diffusion.

Weakness 4: Finally, the authors should provide a comprehensive literature review on extensions of diffusion models, acknowledging previous relevant works. The current literature review is too limited.

Response: We thank the reviewer for pointing this out. We have updated Appendix F.1 in the paper to add a more thorough literature review.

Overall, we thank the reviewer for their constructive suggestions and request the reviewer to reconsider their evaluation. We look forward to the reviewers’ response and would be happy to respond to additional questions/comments.

Weakness 3a: My third concern is related to the experiments and the metrics used. First, there is a large gap in performance when using $\nu = 3$ versus $\nu = 5$, raising questions about the robustness and reliability of the method.

Response: To address this question, we want to highlight a couple of observations.

Firstly, as highlighted in the main text (Proposition 1), our proposed model t-EDM converges to Gaussian diffusion (EDM) in the limit of $\nu \rightarrow \infty$. This implies that as nu increases, the performance of our method on different metrics should converge towards that of EDM. This is indeed what we observe in the results for unconditional generation in Table 2, where for larger nu (nu=5 or 7), the performance of t-EDM (let's say for the VIL channel) starts moving towards that of EDM. Therefore, this performance deviation is expected, and $\nu$ needs to be tuned accordingly.

Secondly, during training, the noising process in t-EDM can be specified as, $x_t = x_0 + \frac{\sigma}{\sqrt{\kappa}}$, where $\kappa \sim \chi^2(\nu)/\nu$. For lower values of $\nu$ (like 3), the quantity $\frac{\sigma}{\sqrt{\kappa}}$ can take large values during training (we present a pdf of the Chi-squared distribution for different nu values in Fig. 8). This effect becomes less pronounced as the value of $\nu$ increases, eventually converging to Gaussian diffusion as $\nu \rightarrow \infty$. Therefore, intuitively, the variability in the model performance for lower values of $\nu$ would be more than the variability at higher $\nu$ values. This is indeed what we observe in Table 1, where the variability in the results for nu=3 and 5 is much more than between nu=5 and 7. This also underscores the importance of tuning $\nu$ correctly for different datasets.

Overall, both these observations justify that our method is reliable for capturing heavy tails, as is also demonstrated in our empirical results in Fig. 3 and Tables 1,2, provided $\nu$ can be tuned reliably for a given dataset.

Weakness 3b: This gap suggests that the chosen metrics (KR, KS, SR) might be sensitive only to the tails of the target distributions. Additional metrics such as precision, recall, and density coverage could provide a more balanced assessment.

Response: Since we focus on comparing different diffusion models based on their ability to capture the tails of the data distribution in this work, our choice of metrics reflects this goal. We agree with the reviewer that additional metrics used in standard generative modeling for images like FID, IS, precision, and recall can provide a more balanced comparison between different models. However, a key caveat in using all these metrics is the reliance on a pre-trained feature extractor (like Inceptionv3/VGG), which is commonly trained on natural images (like ImageNet) for generative modeling of images. However, this work focuses on weather data, which has a much different structure than natural images. Due to this distribution shift, it is unclear if a feature extractor pre-trained on natural images can extract meaningful structural features from weather datasets. Therefore, computing metrics like FID, precision, or recall using a pre-trained network trained on natural images for weather modeling is dubious at best.

In this spirit, training these feature extractors on large-scale weather data using classification or self-supervised losses so that downstream metrics like FID or IS can be reported reliably could be an interesting study on its own, which is outside the scope of this paper. We have also clarified this point in Appendix F.2. Alternatively, we present qualitative histograms between the generated and the train/test set distributions (see Figs. 3,11), which also provide a nice comparison between different methods and demonstrate that t-EDM outperforms other baselines in capturing heavy tails while also modeling the distribution of non-heavy tailed structures accurately.

Weakness 2a: My second concern relates to the results and mathematical formalism of the paper. The mathematical treatment seems quite weak. For instance, the proof in the first section of the appendix contains several typos and appears incorrect in its current form. The reference to Ding (2016) would be sufficient, as it already includes the result

Response: We assume that the reviewer is referring to the proof in Appendix A.1. We have fixed some typos in the proof and made it more elaborate. We also respectfully disagree with the assessment that our mathematical treatment is weak. Some of our core mathematical contributions are:

1. Designing a student-t noising and denoising process (Sections 3.1-3.2)
2. Formulating a training objective for t-diffusions regarding the gamma power divergences (Section 3.3). We show that the training objective for t-diffusions converges to Gaussian diffusion in the limit of $\nu \rightarrow \infty$ by establishing a relation between the KL and the gamma power divergences (Proposition 1). We also show the connection of our loss functions to robust estimators in statistics (Section 5).
3. We derive a novel family of continuous-time samplers for our discrete denoising process (Proposition 2) and instantiate it in a way that provides minimal implementation overhead on top of prior diffusion models like EDM (Section 3.5)
4. We highlight the proof for all propositions and claims in Appendix A.

Overall, we believe our methodological contributions are solid (as also acknowledged by Reviewers EPoh, Nnew, and uzb2).

Weakness 2b: I also found Proposition 2 challenging to interpret, and the accompanying proof does not clarify matters. For example, it is unclear what $f$ and $g$ represent, whether they are any functions satisfying the proposition's conditions, and whether such functions exist. Similarly, the term $dS_t \sim t_d(0, dt, \nu + d)$ lacks clarity. I am also unclear on the main conclusions to be drawn from Proposition 2.

Response: As also highlighted in the main text (lines 261-264), continuous-time formulations of diffusion models can be preferred over their discrete-time analogues due to advances in fast sampling for the former. Therefore, the main idea behind Proposition 2 is to develop a family of continuous-time samplers corresponding to the discrete-time formulation of the diffusion model defined in Sections 3.1-3.3. In this context, firstly, the functions $f$ and $g$ represent design choices that a user can specify to construct continuous time samplers. We mention one plausible choice of such functions in lines 274-279 where we set $f=-\frac{\dot{\sigma}_t}{\sigma_t}$ and $g=0$. This results in our ODE solver, which we use throughout our work. Overall, there do exist functions $f$ and $g$ which can provide useful instantiations for sampler design.

Weakness 2c: Similarly, the term $dS_t \sim t_d(0, dt, \nu + d)$ lacks clarity. I am also unclear on the main conclusions to be drawn from Proposition 2.

Response: The term $dS_t \sim t_d(0, dt, \nu + d)$ intuitively means that we add Student-t noise when sampling from this continuous time SDE. More formally, it denotes the stochastic noise differential, which is distributed as a Student-t random variable and ensures that independent increments are Student-t distributed. This is intuitive since our denoising process is stochastic, with the reverse process posterior modeled as a Student-t distribution. This is also clarified in lines 279-281 in our revision.

Overall, the main idea behind Proposition 2 is to design a family of continuous-time samplers for the proposed heavy-tailed diffusion models, which can provide added advantages over the discrete-time case.

We thank the reviewer for their comments and address each point in more detail below.

Weaknesses/Questions

Weakness 1a: First, a major concern is the soundness of the proposed method. In particular, it is unclear to me that the sequence of conditional distributions q introduced by the authors genuinely defines a Markov noising process.

Response: From our construction of the forward process, the distribution $q(x_t|x_{t-\Delta t})$ can be specified as follows,

$q(x_t|x_{t-\Delta t}) = \int q(x_t|x_{t-\Delta t}, x_0)q(x_0|x_{t-\Delta t}) dx_0$. Therefore, the forward process in our design is non-markovian. However, for most practical scenarios, since diffusion model training solely depends on the perturbation kernel $q(x_t|x_0)$, it does not matter if the forward noising process is Markovian or not. Indeed, prior work in diffusion models like DDIM [1] also relies on using a non-markovian forward process design to achieve sampling speedups.

[1] Denoising Diffusion Implicit Models, Song et al.

Weakness 1b: As such, the loss function they employ does not appear to serve as a meaningful metric. More precisely, it remains uncertain why this loss function would serve as a good proxy for measuring discrepancies between the model marginal $p_0^θ$ and the data distribution. The authors acknowledge that the loss function is not a proper ELBO and claim that the limiting diffusion model is a special case of their method, but further justification of the loss function is needed

Response: This is good feedback that we address in Section 3.3 (lines 220-234) in our revised paper. Below, we argue (in more detail) that minimizing our proposed loss function in Eq. 8 indeed serves as a good proxy for measuring discrepancies between the model and the data distributions for the following reasons.

1. Firstly, the loss function in Eq. 8 is obtained by replacing the KL divergence in the DDPM objective in Eq. 1 with the gamma power divergence instead. Since the gamma power divergence is a valid divergence (See [2]), in principle, the objective in Eq. 8 is minimized when matching the learnable denoiser with the forward process posterior, which in turn implies $p_0^\theta = p_\text{data}$. This divergence-based perspective of the ELBO objective is also utilized in prior works in generative models. For instance, in the context of diffusion models, [1] replace the KL divergence in the DDPM objective using a GAN instead, which provides significant sampling speedups. Similarly, [2] utilizes a similar perspective for training VAEs for heavy-tailed data distributions. Therefore, replacing the KL divergence with alternative divergences still provides a valid proxy for measuring the discrepancies between the model and the data distributions. However, it is worth noting that different divergences would indeed result in different properties of the learned score function and are suitable for different applications. In our context of heavy-tailed diffusion models, using gamma power divergences is natural since, unlike the KL divergence, they can be computed tractably for Student-t distributions.
2. Secondly, in Proposition 1, we formally establish the connection between the KL and gamma power divergence. More specifically, the following equality holds: $D_\gamma = D_{\text{kl}} + O(\gamma)$ where $\gamma = -2/(\nu + d)$ in our work. This formally relates the objective in Eq. 8 to the DDPM objective in Eq. 1. More precisely, both objectives are equivalent in the $\nu \rightarrow \infty$ limit.
3. Lastly, the use of alternative divergences like the gamma power divergences to measure the discrepancies between the model and the data distributions in the presence of outliers/contamination is common in robust statistics [3] (referred to as M-estimators). It serves as a solid alternative to MLE-based objectives like ELBO. We also highlight this connection in Section 5 in the main text.

Therefore, using a divergence like the gamma power divergence in Eq. 8 can serve as a valid proxy for measuring the discrepancies between the model and the data distributions

1] Tackling the Generative Learning Trilemma with Denoising Diffusion GANs, Xiao et al., ICLR’22.

[2] $t^ 3$-Variational Autoencoder: Learning Heavy-tailed Data with Student's t and Power Divergence, Kim et al., ICLR’24.

[3] Variational inference based on robust divergences, Futami et al.

We thank the reviewer for insightful comments and address each point in more detail below.

Weaknesses/Questions

Weakness 1: The only thing is that, there are some duplicated references in the reference list.

Response: We have removed duplicated references in our revised version.

Question 1: Can the authors discuss potential applications in other domains where rare events or extreme values are important, such as finance, seismology, or network traffic analysis. This would help demonstrate the broader impact of the method beyond just weather forecasting.

Response: Thanks for the suggestion. We have added a note mentioning other potential applications in Appendix F.2 in lines 2097-2102.

Question 2: Do you think it possible to extend your framework into other heavy-tailed distributions beyond the Student-t, such as the Cauchy or Lévy distributions, or even unifying all these classes of perturbation kernels? Are there any theoretical or practical challenges you foresee in such extensions?

Response: This is a very interesting question. There has been some work on formulating heavy-tailed diffusion models using $\alpha$-stable Levy distributions [Yoon et al.]. However, except Gaussian and Cauchy distributions, Student-t distributions are not $\alpha$-stable in general (refer to Appendix F.1). Thus, a unification of all these classes of diffusions might be non-trivial. Regardless, this is an interesting direction for further theoretical work in this area, which we have specified in our revision in Appendix F.1.

Overall, we thank the reviewer for their positive comments and feedback. We would be happy to answer any follow-up questions.

Question 3d: Finally, could you clarify how the continuous case fits in with this? I presume, both of the approximations mentioned above are made, and then we take the continuum limit?

Response: We assume that the assumptions that the reviewer refers to are the markovianity of the reverse posterior and the parameterization in Eqn. 5. In this case, the continuum limit is taken after these two design choices.

Question 4a: Connecting the specific divergence used to a parameters in the model / denoising process (i.e. $\gamma = 2/(v+d)$) might be limiting e.g. the objective now can’t be used to learn the degree of freedom, the optimisation objective (8) cannot be used to compare between models that use different settings for $\nu$

Response: Please see our response to Weakness 1, highlighting several methods to learn $\nu$ from the data.

Question 4b: Also, the correspondence to the traditional objective only then holds as q becomes Gaussian (ideally, you’d like to recover the original ELBO without having to limit the model back to the DDPM). I’d have liked to have seen more discussion of this point.

Response: In addition to Proposition 1, which formally presents a connection between the objectives in Eq. 1 (i.e. standard ELBO) and Eq. 8 (t-Diffusion objective), prior works in generative models have also considered the divergence-based perspective of the ELBO objective in Eq. 1. For instance, in the context of diffusion models, [1] replace the KL divergence in the DDPM objective using a GAN instead, which provides significant sampling speedups. Similarly, [2] utilizes a similar perspective for training VAEs for heavy-tailed data distributions. Moreover, replacing the KL divergence with robust divergences for estimation under contamination is commonplace in robust statistics (commonly referred to as M-estimators) [3]. Therefore, replacing the KL divergence with alternative divergences is still principled. We have clarified this point in more detail in Section 3.3 (lines 220-234) of our revision.

1] Tackling the Generative Learning Trilemma with Denoising Diffusion GANs, Xiao et al., ICLR’22.

[2] $t^ 3$-Variational Autoencoder: Learning Heavy-tailed Data with Student's t and Power Divergence, Kim et al., ICLR’24.

[3] Variational inference based on robust divergences, Futami et al.

Question 5: In table 1, it might be useful to add p(x_T) for the two cases (by the way do you mention somewhere that q(x_T) is a product of student-ts and that you set p(x_T) = q(x_T))?

Response: We added another entry in Table 1 indicating the generative prior for Gaussian and t-diffusions. We also highlight that we set the generative prior to $p(x_T)=q(x_T|x_0)$ (this is due to the first term in the t-diffusion objective in Eqn. 8).

Question 6: Why not consider surface fields, arguably where modelling extremes have most impact?

Response: We agree that many impactful surface fields like precipitation are important heavy-tailed variables in atmospheric sciences. Our initial analysis uncovered that liquid cloud concentration was especially heavy-tailed in the dataset at hand – therefore we focused on this most challenging edge case to model. We feel it is natural to expect associated benefits to precipitation modeling. Meanwhile, it is reasonable for people to disagree about how to choose a particular field to demonstrate this method. We hope the reviewer appreciates our own rationale and very much share their interest in seeing these methods extended to benefit surface hydrological modeling.

Question 7: Have you considered a baseline that trains the EDM and then uses a base distribution for generation which is zero-mean unit-variance student-t, sweeping the degree of freedom parameter to find a suitable one?

Response: Thanks for the suggestion. We assume the main objective with such a baseline could be to tune the $\nu$ parameter. While we did not consider such a baseline, we highlight some computational methods that are more principled and can be used for learning $\nu$ in our response to Weakness 1. In contrast, tuning this parameter during sampling could be computationally expensive due to the iterative sampling procedure of diffusion models.

We thank the reviewer for their constructive comments and address each point in more detail below.

Weaknesses/Questions

Weakness 1: I would have liked to have seen a bit more clarity on the formulation of the denoising process (the sub-optimalities here seem significantly more severe than in the Gaussian case, see questions below for more detail). The use of the power-divergence for training was creative, but connecting the parameters of the model to the divergence parameters also seemed sub-optimal as it does not allow the framework to learn the degree of freedom parameter in the student-t distribution.

Response: Regarding learning the degree of freedom parameter $\nu$, while we don’t learn this parameter in this work, there can be several ways to learn this parameter directly from the data.

Firstly, it is worth noting that the parameter $\nu$ directly influences training in heavy-tailed diffusions. This is because the noising process $x_t = \mu_t x_0 + \frac{z}{\sqrt{\kappa}}$ where $\kappa \sim \chi^2(\nu)/nu$. Since the Chi-squared distribution is a special case of a Gamma distribution, we can perhaps use reparameterization [1] to estimate the parameters (i.e., $\nu$) of this distribution using a neural network and backprop end-to-end to learn this parameter. This would correspond to learning $\nu$ directly from the data.

[1] The Generalized Reparameterization Gradient, Ruiz et al.

Alternatively, one can also estimate $\nu$ from the data in a separate stage and use it as a perturbation kernel parameter. We highlight one simple way to do this in Appendix D, where we attempt to derive the optimal noise schedule for EDM and t-EDM by minimizing the mutual information between the $x_0$ and the noisy data point $x_\sigma$ w.r.t $\sigma$. In principle, for t-EDM, we can also solve this minimization problem jointly over $\nu$ and $\sigma$ to learn an initial $\nu$ from the data.

We have added these points in Appendix F.2 in our revision. Overall, while manual tuning of $\nu$ might seem to be a sub-optimality of our framework, we would like to point out that the proposed framework can be readily extended to learn $\nu$ in an end-to-end or separate stage fashion.

Weakness 2: The experimentation was a little light in that only two variables were considered in the conditional and unconditional settings, there were no experiments that considered forecasting multiple variables at one time, or autoregressive rollouts.

Response: Regarding the choice of variables, we chose variables that are most heavy-tailed in the HRRR dataset to illustrate the efficacy of our framework. Regarding multiple variables, since tuning $\nu$ per variable can be quite tedious, we leave exploration of this aspect to further work. We also highlight some methods for automatically tuning this parameter in our response to Weakness 1 above. Since we can combine techniques like per-channel-preconditioning, i.e., PCP and automatic $\nu$ tuning with heavy-tailed diffusion models, the design space for experimentation when exploring these models in the context of autoregressive rollouts in weather forecasting is quite extensive and, therefore, beyond the scope of this work. We mention this as a limitation in Appendix F.2 in our revised paper.

Weakness 3: Setting the degree of freedom parameter via a parameter sweep is a bit painful. It’s a shame that it can’t be learned directly, especially if it is going to be different for different dimensions of x (locations / variables / pressure levels etc.).

Response: Kindly, refer to our response to Weakness 1 above, where we clarify some ways to tune $\nu$ directly from the data.

Question 8: For the HRRR experiments, I presume that you are training off HRRR analysis (rather than, say forecasts) this should be mentioned (I didn’t see it in the main text or in the appendix). Also, is the time interval between inputs mentioned? Is it 1hr like in StormCast?

Response: The reviewer is correct that we used data similar to what was trained in StormCast. To be clear, this is not the analysis but rather model forecast fields 1-hour post-analysis. Such a time lag is necessary to allow for some additional model spin-up; cloud condensate is not directly assimilated. We have added a sentence in Appendix C.2.1 to clarify this. We expect the methods presented to work equally well with model forecast or heavy-tailed fields available in the analysis without losing generality.

Question 9: The spatial fields shown in the appendix were useful. For the forecasts results shown in figure 8 and 9, it would be useful to show the persistence forecast.

Response: Thanks for the suggestion. We will include persistence forecasts in a subsequent revision.

Question 10: For the unconditional sampling, it looks like on the test set t-EDM outperforms the baselines 1/3 of the time (and only on the KS metric).

Is the KS metric to be more trusted? Why weren’t the other metrics applied to the forecasting case too?

Response: This is an interesting observation. One reason for our method exhibiting better performance for the KS metric could be the fact that the KS metric, in addition to flattening the data, is computed over 99.9th and 0.1th percentiles for the right and left tails, respectively. This implies that the KS metric specifically examines the tail values and thus represents an appropriate metric to capture the difference between the tail estimation for the generated and the data samples. In contrast, the KR and SR metrics are based on computing higher-order moments over the entire domain of flattened empirical samples and thus might be more coarse in nature, hence the discrepancy.

However, it is also worth noting that while we present EDM baselines with improved preconditioning and preprocessing as additional baselines, these can also be applied to t-EDM, potentially improving results. Therefore, a more fair comparison should be between t-EDM and the baseline EDM for which t-EDM outperforms EDM by a significant margin in all metrics. We added this observation in Section 4.1. (lines 406-411).

For conditional modeling, since standard metrics like CRPS, RMSE, and SSR already capture forecasting performance while the KS metric captures tail coverage, we do not use KR and SR metrics.

Question 11: It would be interesting to see how t-Flow behaved during auto-regressive rollouts. I could imagine that the under-dispersion of EDM might actually help maintain stability in these cases, whereas heavier tails could lead to greater instability.

Response: In this work we only consider the next-frame prediction task for which we did not observe any instability issues when working with t-EDM or t-Flow models. While examining this behavior for autoregressive rollouts is indeed interesting, it feels beyond the scope of this work due to space constraints.

Overall, we thank the reviewer for their constructive feedback and for highlighting many interesting points in their review. We would also be happy to answer any subsequent questions.

Question 1: In summary, I am wondering if the performance of EDM can be increased only by preconditioning to reach that of t-EDM or if t-EDM does provide a significant difference not achievable by preconditioning.

Response: This is a valid question and, in fact, was also a motivation for us to explore the PCP and INC techniques as comparison baselines in the first place. While the possible design space for constructing more efficient preconditioning and preprocessing techniques is quite large, and it is not feasible to experiment with a lot of these different combinations, several observations justify the use of heavy-tailed diffusion models for modeling extreme events rather than relying solely on preconditioner like heuristics for improved modeling performance:

1. Theoretical arguments: [1] show that a class of Lipschitz-continuous transport maps acting on the source distribution leads to target distributions with tail properties similar to the source distribution. Therefore, a diffusion/flow model trained using a Gaussian base distribution (which is light-tailed) is unlikely to capture a target heavy-tailed distribution. Moreover, unlike a heavy-tailed prior, a light-tailed distribution like a Gaussian prior lacks the capacity to encode usual and extreme events [2]. Similarly, [3] argues that standard denoising score matching fails to capture heavy tails due to spherical Gaussian shells in high dimensions. Therefore, the inability of Gaussian diffusion models to capture tail events is likely a more fundamental modeling issue rather than one that can be resolved by updating preconditioners (PCP) or data preprocessing schemes (INC).

2. Empirical arguments: The theoretical arguments presented above are also supported by empirical results. For unconditional generation, we experimented with different PCP schemes for modulating the $\pi_\text{mean}$ parameter for controlling the noise schedule in EDM training and stuck with the RBF kernel due to better performance. However, we also observed that heuristics like PCP do not scale when training on data with multiple channels, with performance degrading rapidly (in terms of qualitative histograms) as opposed to a single-channel setting.

   As an example of this observation, we trained an EDM model on HRRR data with 10 channel inputs (as opposed to single-channel training presented in the main text). We compare qualitative histograms between the 10-channel and single-channel settings for the VIL and w20 channels in Figure 7 in Appendix C.1.7. It can be seen clearly that while PCP might prove to be effective on single-channel training (possibly due to easier modeling), on multi-channel runs, modeling tails takes a severe degradation. Therefore, under such scenarios, PCP/INC can only help so much, and t-EDM can possibly provide a more principled alternative (provided there is a principled mechanism to tune the hyperparameter nu) to model heavy tail scenarios. (We include full experimental details for our 10-channel run in Appendix C.1.7).

   Moreover, while heuristics like PCP/INC are illustrated in the main text as competing baselines, it is worth noting that they can also be used under the t-EDM framework and should be considered complementary to our proposed heavy-tailed diffusion framework for more challenging modeling tasks. Below, we also provide a more detailed response to the individual questions posed by the reviewer.

[1] Tails of Lipschitz triangular flows, Jaini et al.

[2] $t^ 3$-Variational Autoencoder: Learning Heavy-tailed Data with Student's t and Power Divergence, Kim et al.

[3] Heavy-tailed denoising score matching, Deasy et al.

We thank the reviewer for insightful comments and address each point in more detail below.

Weaknesses/Questions

Weakness 1*: The experiments for unconditional generation in Table 2 show that the EDM (using preconditioning) outperforms the t-EDM for the KR and SR metric in the VIL channel with test data. This hinders the claim the authors make that standard diffusion models, even with preconditioning, fail to capture the heavy-tailed distributions (line 043) and that t-EDM outperforms standard diffusion models. Qualitatively, the plots in Figure 3 do show that t-EDM outperforms EDM with preconditioning in capturing the heavy-tailed behavior. However, I believe a justification of why the KR and SR metrics do not reflect this should be included. As addressed in Appendix F.2, this metrics rely on flattening the data. This may be one reason, but the metric KS also relies on flattening the data but does show that t-EDM outperforms EDM. In any case, I believe a more throrough justification or explanation is needed here.*

Response: This is an interesting observation. We elaborate upon this aspect as follows.

Firstly, another reason for our method exhibiting better performance for the KS metric could be the fact that the KS metric, in addition to flattening the data, is computed over 99.9th and 0.1th percentiles for the right and left tails, respectively (see Appendix C.1.3). This implies that the KS metric specifically examines the tail values and thus represents an appropriate metric to capture the difference between the tail estimation for the generated and the data samples. In contrast, the KR and SR metrics are based on computing higher-order moments over the entire domain of flattened empirical samples and thus might be more coarse in nature, hence the discrepancy.

However, it is worth noting that while we present EDM baselines with improved preconditioning and preprocessing as additional baselines, these can also be applied to t-EDM, potentially improving results. Therefore, a more fair comparison should be between t-EDM and the baseline EDM for which t-EDM outperforms EDM by a significant margin in all metrics. We also clarify this in the main text in Section 4.1 (lines 406-411).

Secondly, due to the flattening operation, the metrics considered in this work disregard the spatial structure of the generated samples, which is undesirable. The qualitative histogram plots in Figure 3 help alleviate this to some extent. However, designing metrics that can quantify heavy-tailed behavior while considering the spatial structure of the generated samples is an interesting direction for further work but is outside the scope of this paper. We discuss these limitations in more detail in Appendix F.2.

Question 4: It seems unintuitive to define a diffusion model via its finite increments, only to show it converges to a continuous time diffusion and later make use of this perspective. Can the authors comment on why the did not simply define the noising process as a student-t diffusion in the first place and later discretise this?

Response: We understand different presentation styles might be possible for the methods introduced in this work. However, our presentation of the discrete-time version followed by the continuous-time formulation is based on the following rationale:

1. Firstly, presenting the discrete-time version highlights key design choices involved in the noising process design. For instance, when designing the forward noising process, we directly compute the perturbation kernel $q(x_t|x_0)$ without specifying the transition kernel $q(x_t|x_{t-\Delta t})$, which is non-markovian under our design. This is justified since only the former is utilized during diffusion model training. We also highlight this aspect in the main text in lines 149-150. Similarly, the discrete-time perspective highlights the choice of $\sigma_{12}(t)$ and $\sigma_{21}(t)$ as important for sampler design. The latter might be ignored in the continuous time presentation.
2. Secondly, presenting the discrete-time version also highlights the interpretation of our “x0-prediction” objective as minimizing the gamma-power divergence between the forward and the reverse posterior distributions, providing additional intuition on why our proposed diffusion model can capture heavy tails. On the contrary, this observation is non-trivial when directly defining the continuous-time denoising process followed by the “x0-prediction” training objective.
3. Lastly, presenting the discrete-time version first might be more accessible to readers with a background in discrete-time DDPM [Ho et al.] and not very familiar with the continuous-time theory of diffusion processes. This is because continuous treatment of t-diffusion models would most likely involve dealing with Kolmogorov equations for Levy-like processes, making the presentation very dense.

Overall, we thank the reviewer for their constructive suggestions and questions. We look forward to the reviewers’ response and would be happy to respond to additional questions/comments.

Question 2 proposition 2 to me is not well defined. The construction of the stochastic process with student t increments. Can the authors discuss how such a process is constructed in more detail?

Response: We re-highlight the overall construction for our heavy-tailed diffusion process:

1. Forward Process: We parameterize the forward process using Student-t distributions. More specifically, we define a joint distribution over $x_t$ and $x_{t-\Delta t}$ and then compute the perturbation kernel $q(x_t|x_0)$ using properties of multivariate Student-t distributions, which is utilized during training. Kindly refer to Sections 3.1-3.3 for more details on constructing the forward noising process. We would also like to refer the reviewer to Figure 2 (Algorithm 1) for an overview of the noising process and how it is utilized during training.
2. Reverse Process: We model the reverse process using a Student-t distribution based stochastic process. More specifically, we parameterize the denoising diffusion posterior $p(x_{t-\Delta t}|x_t)$ using a denoiser $D_\theta(x_t, t)$. Given a pre-trained denoiser, we can perform ancestral sampling from the posterior $p(x_{t-\Delta t}|x_t)$ during sampling. Alternatively, we can derive a continuous time formulation of our reverse process and perform sampling using off-the-shelf fast samplers for diffusion models like DDIM, EDM-ODE, etc. The latter approach is preferable since it provides better sampling speedups. Proposition 2 provides a way to construct such continuous time samplers in this context. While the proof for Proposition 2 is provided in Appendix A.5, the main idea is to “de-discretize” the posterior update in Eqn. 5 to design a family of continuous-time samplers parameterized using user-specified functions $f$ and $g$, which satisfy the constraint (as stated in Proposition 2). We also provide an instantiation of these functions and derive an ODE solver, which we use throughout this work. An overview of the sampling process using our ODE solver is given in Figure 2 (Algorithm 2).

We hope that this high-level overview clarifies the construction of the process with Student-t increments and the significance of the result in Proposition 2.

Question 3: In the proof of proposition 2 the authors derive the last term converging to Brownian motion, not a student process, can they please comment on this discrepancy.

Response: There might be some misunderstanding here. The noise increments in Proposition 2 are indeed distributed as Student-t random variables, as also stated in the statement for Proposition 2 on line 276. We do notice that there is a final step missing in the proof in Appendix A.5 where it can be stated that the random variable $dS_t = dW_t / \sqrt{\kappa_t}$ is Student-t distributed as $\kappa_t \sim \chi^2$ distributed. We have clarified this line in more detail in the proof in Appendix A.5 to make this more explicit for the readers.

We thank the reviewer for insightful comments and address each point in more detail below.

Weaknesses/Questions

Weakness 1: The paper fails to discuss previous work, for example "Heavy-tailed denoising score matching" covering the subject matter.

Response: We thank the reviewer for mentioning this work. However, several key differences exist between our proposed methods and the Generalized Normal distribution-based score matching introduced in the mentioned paper, which we have clarified in Section F.1 in our revision with a more detailed explanation (also see our response to Question 1 below).

Weakness 2: The paper lacks experiments on synthetic data which is generated as heavy tailed. While not needed, I believe it would help make the case that this method definitively helps with modeling heavy tails by demonstrating it on definitively heavy tailed data.

Response: We would like to point the reviewer to Fig. 1, where we compare the performance of the baseline EDM with t-EDM for different $\nu$ values on a 2-D toy dataset (also referred to as the Neal’s Funnel dataset) with heavy tails along the x-axis dimension. Fig. 1 highlights that t-EDM provides a parameter $\nu$, which can help control tail estimation and is thus better able to model heavy-tailed distributions as compared to other baselines like standard Gaussian diffusion. We also include additional quantitative results in Table 8 for the same 2-D toy dataset in Appendix C.3 (section “Quantitative Results on Neals Funnel”). We also include the results from Table 8 below for convenience.

The quantitative results in the above table for the toy dataset support our qualitative results in Figure 1 in the main text, and imply that our proposed model can definitively help with heavy-tailed data modeling as compared to baselines like Gaussian Diffusion.

Question 1: Can the authors comment on the drawbacks of using student t noise discussed in "Heavy-tailed denoising score matching" and why this makes sense in the context proposed?

Response: We assume the reviewer refers to the paper [1] "Heavy-tailed denoising score matching" by Deasy et al. The motivation to use heavy-tailed noising distributions is similar between [1] and our work, i.e., using base Gaussian distribution for noising concentrates on a narrow shell and neglects the tails (lines 045-046 in the main text). However, in this context, it is worth noting a few key differences from [1]:

1. Relevance of using Student-t distributions: The choice of a Generalized Normal Distribution in [1] as a noising process is mostly inspired by the convenience of working with finite moments of the underlying noisy state vector. The choice of the Student-t distribution in [1] is discouraged due to non-finite moments. However, we argue that this is not representative of whether a Student-t kernel can/cannot be used for modeling heavy-tailed distributions. In fact, as we show in this work, the Student-t based noising process outperforms traditional Gaussian diffusion models for modeling heavy-tailed data distributions (for instance, see Fig. 1). Moreover, as we discuss in Appendix F.1, several prior works have utilized Student-t distributions for other model families like VAEs [2] and likelihood-based flows [3], justifying our choice of the heavy-tailed distribution.
2. Easy Sampler Design: In [1], the authors rely on annealed Langevin dynamics for sampling and briefly discuss the existence of a reverse SDE for sampling in Appendix B.3. However, an exact formulation of such a reverse process is missing from [1] and would involve levy formulations of Kolmogorov’s forward and reverse equations which can be quite complex. In contrast, due to the conditional properties of multivariate Student-t distributions, our method does not rely on a reverse process formulation based on Levy processes (see Proposition 2 in the main text). Therefore, designing samplers using our proposed Student-t diffusion formulation is easier than the Generalized Normal distribution-based noising formulation discussed in [1].

Therefore, using a Student-t formulation makes sense in the proposed context. We have included these details in Appendix F.1.

[2] Tails of Lipschitz triangular flows, Jaini et al.

[3] $t^ 3$-Variational Autoencoder: Learning Heavy-tailed Data with Student's t and Power Divergence, Kim et al.