Manifolds, Random Matrices and Spectral Gaps: The geometric phases of generative diffusion

## Abstract

Abstract In this paper, we investigate the latent geometry of generative diffusion models under the manifold hypothesis. To this purpose, we analyze the spectrum of eigenvalues (and singular values) of the Jacobian of the score function, whose discontinuities (gaps) reveal the presence and dimensionality of distinct sub-manifolds. Using a statistical physics approach, we derive the spectral distributions and formulas for the spectral gaps under several distributional assumptions and we compare these theoretical predictions with the spectra estimated from trained networks. Our analysis reveals the existence of three distinct qualitative phases during the generative process: a trivial phase; a manifold coverage phase where the diffusion process fits the distribution internal to the manifold; a consolidation phase where the score becomes orthogonal to the manifold and all particles are projected on the support of the data. This ‘division of labor’ between different timescales provides an elegant explanation on why generative diffusion models are not affected by the manifold overfitting phenomenon that plagues likelihood-based models, since the internal distribution and the manifold geometry are produced at different time points during generation.

## 1 Introduction

Generative diffusion models have revolutionized the fields of computer vision and generative modeling, achieving state-of-the-art performance on image generation (Ho et al.,, ; Song and Ermon,, ; Yang et al.,, ) and video generation (Ho et al.,, ; Singer et al.,, ; Blattmann et al.,, ; Tim et al.,, ). These models synthesize images through a stochastic dynamical denoising process. Experimental and theoretical arguments suggest that different features such as frequency modes and class labels are generated at different times during the process. For example, it has been shown that separation between isolated classes, like in the case of mixture of Gaussian models, happens at critical phase transition points of spontaneous symmetry breaking (Raya and Ambrogioni,, ; Ambrogioni,, ) also known as speciation events (Biroli and Mézard,, ), which appears to be hierarchically organized (Sclocchi et al.,, ). It is also well known that subspaces corresponding to different frequency models emerge at different times of diffusion (Kingma and Gao,, ). This idea has been recently refined by (Kadkhodaie et al.,, ), who showed that diffusion models give rise to a local decomposition of the image manifold into a basis of geometry-adaptive harmonic basis functions. While there are conceptual seminaries with symmetry breaking phenomena, these decomposition effects cannot be directly explained in terms of critical phase transitions as they are fundamentally linear processes. In this paper, we will provide a precise theoretical analysis of the separation of subspaces for data defined on low dimensional linear manifolds. Following (Wenliang and Moran,, ), (Horvat and Pfister,, ) and (Stanczuk et al.,, ), we will analyze the singular spectrum of the score and its relation to the manifold dimensionality.

Our main contributions are: I) an in-depth theoretical random-matrix analysis of the distribution of Jacobian spectra in diffusion models on linear manifolds and II) a detailed experimental analysis of Jacobian spectra extracted from trained networks on linear manifolds and on image datasets. The analysis of this spectra is important as it provides a detailed picture of the latent geometry that guides the generative diffusion process. We show that the linear theory predicts several phenomena that we observed in trained networks. Based on our result, we divide the generative process in three qualitatively different phases: trivial phase, manifold coverage phase and manifold consolidation phase. Using these concepts, we provide an elegant explanation of why diffusion models can avoid the manifold overfitting pathology that characterizes likelihood-based models (Loaiza-Ganem et al.,, ).

## 2 The manifold hypothesis

The manifold hypothesis states that the distribution on natural data, such as images and sound recordings, is (approximately) supported on a $m-$dimensional manifolds $\mathcal{M}$ embedded in a larger euclidean ambient space $\mathds{R}^{d}$ (Peyré,, ; Fefferman et al.,, ). While a data probability distribution supported on a $m<d$ manifold $\mathcal{M}$ cannot be expressed using a proper density function, we loosely define such density as

$$ $p_{0}(\boldsymbol{x})=\delta_{\mathcal{M}}(\boldsymbol{x})\,\rho_{\text{int}}(\boldsymbol{x})~{},$ (1) $$

where $\delta_{\mathcal{M}}$ is the Dirac function for the manifold and such that $\int_{\mathds{R}^{d}}\delta_{\mathcal{M}}(\boldsymbol{x})\bullet d\boldsymbol{x}=\int_{\mathcal{M}}\bullet\,d\boldsymbol{x}$. We call $\rho_{\text{int}}(\boldsymbol{x})$ the *internal density* instead, that is the density restricted to the manifold. The density $p_{0}(\boldsymbol{x})$ is zero outside the manifold and diverges on the manifold.

Figure: Figure 1: Visualization of the gaps in the spectrum of the (negative) Jacobian of the score for data supported on a latent manifold. Blue line: idealize spectrum of distribution with uniform internal density; Orange line: spectrum of a more realistic distribution.
Refer to caption: /html/2410.05898/assets/spectral_gaps.png

### 2.1 Manifold overfitting in likelihood-based models

Likelihood-based generative models are defined by a highly parameterized likelihood function $f(\boldsymbol{x};\boldsymbol{\theta})$, whose parameters are trained by minimizing the loss

$$ $\mathcal{L}(\boldsymbol{\theta})=-\mathbb{E}_{\boldsymbol{x}\sim p_{0}(\boldsymbol{x})}\!\left[f(\boldsymbol{x};\boldsymbol{\theta})\right]~{},$ (2) $$

which maximizes the probability of the data given the model. This maximum likelihood loss is minimized if $f(\boldsymbol{x};\boldsymbol{\theta})=p_{0}(\boldsymbol{x})$. Unfortunately, the divergence of $p_{0}(\boldsymbol{x})$ on the manifold implies that the model density $f(\boldsymbol{x};\boldsymbol{\theta})$ on the manifold becomes larger and larger during training. More problematically, the optimization problem becomes almost insensitive to the internal density $\rho_{\text{int}}(\boldsymbol{x})$. This phenomenon is called *manifold overfitting* (Loaiza-Ganem et al.,, ), since the trained model fits the manifold while ignoring its internal density, resulting in poor generation. While generative diffusion models have often been characterized in terms of a likelihood function, experimental evidence suggests that they are not affected by this pathology. In the rest of the paper, we will provide an elegant explanation of their behavior by analyzing their latent geometry.

## 3 Background on generative diffusion models

Here, we will consider a simple variance-exploding forward process where the data $\boldsymbol{x}_{0}\sim p_{0}(\boldsymbol{x})$ evolves according to the equation

$$ $\text{d}\boldsymbol{x}_{t}=\text{d}\boldsymbol{Z}_{t}$ (3) $$

where $\text{d}\boldsymbol{Z}_{t}$ is a standard Brownian motion. The formal solution of Eq. ([3](#S3.E3)) can be given in term of the heat kernel:
$p(\boldsymbol{x}_{t},t)=\mathbb{E}_{\boldsymbol{x}_{0}\sim p_{0}}\!\left[\frac{1}{\sqrt{2\pi t}}e^{-\frac{\left\lVert\boldsymbol{x}_{t}-\boldsymbol{x}_{0}\right\rVert_{2}^{2}}{2t}}\right]~{}.$
The *target distribution* $p_{0}(\boldsymbol{x})$ is then recovered by reversing the diffusion process (Anderson,, ). We initialize this reverse process from $\boldsymbol{x}_{t_{f}}\sim\mathcal{N}\left(0,t_{f}I_{d}\right)$ at some large time $t_{f}$, which evolves backward according to

$$ $\text{d}\boldsymbol{x}_{t}=-\nabla_{\boldsymbol{x}}\log p_{t}(\boldsymbol{x}_{t})\text{d}t+d\boldsymbol{Z}_{t}$ (4) $$

The function $s(\boldsymbol{x},t)=\nabla_{\boldsymbol{x}}\log p_{t}(\boldsymbol{x})$ is is the so-called score function. From a set of training points $\{\boldsymbol{y}^{1},\dots,\boldsymbol{y}^{N}\}\stackrel{{\scriptstyle\text{iid}}}{{\sim}}p_{0}$, we can train a neural approximation of $s(\boldsymbol{x},t)$ by learning a denoising autoencoder $\hat{\boldsymbol{\epsilon}}_{\boldsymbol{\theta}}(\mathbf{x},t)$, which is trained to recover the standardized noise $\epsilon_{t}$ from the noisy state $x_{t}=x_{0}+\sqrt{t}\epsilon_{t}$ (Hyvärinen and Dayan,, ; Vincent,, ; Ho et al.,, ). In order to avoid singularities in neural network output, the learned score is parametrized as
$\hat{s}_{\boldsymbol{\theta}}(\boldsymbol{x},t)=-\frac{\hat{\boldsymbol{\epsilon}}_{\boldsymbol{\theta}}(\boldsymbol{x},t)}{\sqrt{t}}~{}.$

It is also convenient to define the support score $\tilde{s}_{\mathcal{M}}(\boldsymbol{x},t)$, defined as the score function obtained from the uniform data distribution $\tilde{p}_{0}(\boldsymbol{x})=\frac{1}{|\mathcal{M}|}\delta_{\mathcal{M}}(\boldsymbol{x})$.

## 4 Dynamic latent manifolds and spectral gaps

Consider a generative diffusion model with $p_{0}(\boldsymbol{x})$ defined on a $d$-dimensional manifold $\mathcal{M}_{0}$ according to Eq. ([2](#S2.E2)). In the course of the diffusion process, we can define a time-dependent family of stable latent manifolds

$$ $\mathcal{M}_{t}=\{\boldsymbol{x}^{*}\mid\tilde{s}_{\mathcal{M}}(\boldsymbol{x}^{*},t)=0,\text{ with }J_{\mathcal{M}}(\boldsymbol{x}^{*},t)~{}\text{ n.s.d.}\}~{},$ (5) $$

where the negative semi-definiteness (n.s.d.) is a stability condition on the Jacobian matrix $J_{\mathcal{M}}(\boldsymbol{x},t)$ of the support score $\tilde{s}_{\mathcal{M}}(\boldsymbol{x}^{*},t)$. Due to the noise, the diffusing particles are likely to be found in shells of radius $\sqrt{t}$ around each point of the latent manifold.

For a small perturbation $\boldsymbol{p}$ around a point $\boldsymbol{x^{*}}$ on the latent manifold at time $t$, the score function is well approximated by its linearization:

$$ $s(\boldsymbol{p},t)\approx J(\boldsymbol{x}^{*},t)~{}\boldsymbol{p}=-\sum_{j}\left(\boldsymbol{v}_{j}\cdot\boldsymbol{p}\right)\lambda_{j}(\boldsymbol{x}^{*},t)\boldsymbol{v}_{j}~{}~{},$ (6) $$

where $J(\boldsymbol{x}^{*},t)$ is the Jacobian of the score and the $\boldsymbol{v}_{j}$ and $\lambda_{j}(\boldsymbol{x}^{*},t)$ are respectively the $j$-th eigenvector and the associated eigenvalue of $-J(\boldsymbol{x}^{*},t)$. The spectrum of eigenvalues provides detailed information concerning the local geometry of the manifold around. Perturbations aligned to the tangent space of $\mathcal{M}_{t}$ correspond to small eigenvalues, while orthogonal perturbations correspond to high eigenvalues, as the score tends to push the stochastic dynamics towards its fixed-points. Therefore, we can estimate the dimensionality of the manifold from the location of a drop (i.e. a sharp change) in the sorted spectrum of eigenvalues (Stanczuk et al.,, ). This is visualized in Fig. [1](#S2.F1). This drop corresponds exactly to a gap (i.e. a separation) in the eigenvalues spectrum; in the following we will refer to both as gaps.

### 4.1 subspaces and intermediate gaps

Consider the situation where the internal density $\rho_{\text{int}}(\boldsymbol{x})$ is not locally flat around a point $\boldsymbol{x}^{*}\in\mathcal{M}_{0}$. In this case, at a finite time $t$ the actual score function does not vanish on the latent manifold $\mathcal{M}_{t}$ as there is a gradient of the log-density along the tangent directions. This implies that the spectrum of tangent eigenvalues can have a series of sub-gaps with separate different tangent subspaces with different ‘local variance’. In image generation tasks, these subspaces are often associated with different frequency modes, as noted in (Kingma and Gao,, ).

Consequently, we can quantify the sensitivity to the internal density at time $t$ by studying the statistics and temporal evolution of intermediate gaps

$$ $\Delta^{\text{GAP}}_{k}(\boldsymbol{x}^{*},t)=\lambda_{k+1}(\boldsymbol{x}^{*},t)-\lambda_{k}(\boldsymbol{x}^{*},t)~{},$ (7) $$

where the indices $k$ depend on the dimensionality of the subspaces. Note however that under realistic data distributions it is unlikely to find sharp intermediate discontinuities since each subspace will have a different eigenvalue, resulting in a smooth gradient.

## 5 Phenomenology of generative diffusion on manifolds

This section contains an intuitive picture that follows from our theoretical results on linear models, which we will fully outline in the next section. The theory considers the case of a linear manifold with Gaussian internal distribution. While only linear models are theoretically tractable, we conjecture that their phenomenology captures the main features of subspace separation in the tangent space of curved manifolds. We validated the theory using networks trained on both linear data and highly non-linear data such as natural images (see section [8](#S8)). Based on the dynamics of the spectral gaps, we found that the generative dynamics of $\boldsymbol{x}_{t}$ according to Eq. ([4](#S3.E4)) can be separated into three distinct phases. The phase separation do not correspond to singularities as there are cross-over events, not genuine phase transitions.
During all our analysis we will exclusively work with eigenvalues since, in the linear manifold model, the Jacobian of the true score is symmetric. The same phenomenology is nevertheless fully appreaciable when using the singular values in our experiments with neural approximations of the score.

### 5.1 Phase I: The trivial phase

In the trivial phase, the diffusing particle moves according to the noise distribution without strong biases towards the manifold directions. In this dynamic regime, the latent manifold $\mathcal{M}_{t}$ is a single point surrounded by an isotropic quadratic well of potential. The spectral gaps are not visible and all eigenvalues have approximately the same value due to the isotropy of the noise distribution. This trivial phase is analogous to the initial phases described in (Raya and Ambrogioni,, ) and (Biroli et al.,, ).

### 5.2 Phase II: Manifold coverage

The manifold coverage phase begins with the opening of the first of a series of spectral gaps corresponding to local subspaces. In this phase, different subspaces with different variances can therefore be identified by intermediate gaps in the spectra, as sketched in Fig. [2](#S5.F2). When the intermediate gaps are opened, the diffusing particles spread across the relevant manifold directions according to their relative variances. In other words, during this regime of generative diffusion, the process fits the distribution of the data internal to the manifold.

We assume low-rank covariance $\Sigma=FF^{\top}$ for the data distribution.
In terms of random matrix theory, the gap-forming phenomenology has two distinct processes: the emergence of intermediate gaps (i.e. steps in the dimensionality plot) between separated bulks of the spectrum, the opening of a final gap that allows to infer the dimensionality of the full manifold.
Our analysis gives us the time scale at which such intermediate gaps are maximally opened, i.e.

$$ $t_{\text{max}}^{(k)}=\sqrt{\gamma_{+}(\sigma_{k})\,\gamma_{-}(\sigma_{k+1})},$ (8) $$

where $\gamma_{-}(\sigma_{k+1})$ and $\gamma_{+}(\sigma_{k})$ are specific eigenvalues of $\Sigma$ (see Fig. [9](#A1.F9)) associated to two hierarchically consecutive variances (see Appendix [A.2](#A1.SS2) for an exhaustive analysis).
In most cases, when $\sigma_{k+1}^{2}\ll\sigma_{k}^{2}$, the dependence on the two variances is $\mathcal{O}(\sigma_{k}\cdot\sigma_{k+1})$. This is the timescale where the score is maximally sensitive to the relative variance of the two subspaces, which guides the particles toward the correct internal distribution.

### 5.3 Phase III: Manifold consolidation

Finally, the manifold consolidation phase is characterized by the asymptotic closure of the intermediate gaps and the sharpening of the total manifold gap, indicating the full dimensionality of $\mathcal{M}$. In this final regime, the score assumes the form

$$ $\nabla_{\mathbf{x}}\log p_{t}(\mathbf{x})\simeq\frac{1}{t}\Big{[}\Pi-I_{d}\Big{]}\mathbf{x}.$ (9) $$

where $\Pi=F(F^{\top}F)^{-1}F^{\top}$ is the projection matrix over the manifold.

Figure: Figure 2: Visualization of the trivial (I), coverage (II) and consolidation phase (III).
Refer to caption: /html/2410.05898/assets/x1.png

The component of the score orthogonal to the manifold diverges proportionally to $t^{-1}$, while the tangent components converge to a constant and become therefore negligible in this regime.
This results in the consolidation of the gap corresponding to the manifold dimensionality $m$ and to the (relative) closure of the intermediate gaps. Therefore, in this final phase the dynamics of the model simply projects the particles into the manifold $\mathcal{M}_{t}\rightarrow\mathcal{M}_{0}$. In the generative modeling literature, this phenomenon is also known as *manifold overfitting* as the terms corresponding the the internal distribution is negligible (Loaiza-Ganem et al.,, ). Fortunately, this internal coverage occurs in the model throughout intermediate phases across the previous diffusion process, resulting in a balanced coverage of the internal distribution.

## 6 Theoretical analysis of the spectral of linear diffusion models

In this section, we provide our main theoretical results concerning the spectral distribution for random linear subspaces and the relative spectral gaps formulas. We start by reviewing diffusion with data supported on linear manifolds, where the exact score function can be computed.

### 6.1 Linear manifolds

Normally the distribution $p_{0}(\boldsymbol{x})$ is unknown. It is however interesting to investigate tractable special case where the distribution is a multivariate Gaussian defined on a linear manifold:

$$ $\boldsymbol{y}^{\mu}=F\boldsymbol{z}^{\mu}$ (10) $$

where $F\in\mathbb{R}^{d\times m}$ is an arbitrary projection matrix that implicitly define the structure of the latent manifold. and $\boldsymbol{z}^{\mu}\sim\mathcal{N}(0,I_{m})$ the latent space vector. In this setting, the distribution can be explicitly written as

$$ $p_{0}(\mathbf{x})=\mathbb{E}_{\mathbf{z}\sim\mathcal{N}(0,I_{m})}\delta(\mathbf{x}-F\mathbf{z})=\mathcal{N}(\mathbf{x};0,FF^{\top})$ (11) $$

Therefore, the density of the process at a given time $t$ is again Gaussian and can be computed from

$$ $p_{t}(\boldsymbol{x})=\mathbb{E}_{\mathbf{z}}\frac{1}{\sqrt{(2\pi t)^{d}}}\ e^{-\frac{1}{2t}\|\mathbf{x}-F\mathbf{z}\|^{2}}.$ (12) $$

While linear manifolds are very simple when compared with real data, they still exhibit a rich and non-trivial phenomenology that elucidate several universal phenomena of diffusion under the manifold hypothesis. In fact, these linear models capture the structure of tangent spaces of smooth manifolds.

The score function of the linear model is solvable analytically since we only have to perform Gaussian integrals, from which we obtain a quadratic form in $\mathbf{x}$ that we can rewrite as

$$ $\log p_{t}(\mathbf{x})=\frac{1}{2t}\mathbf{x}^{\top}J_{t}\mathbf{x}+const.$ (13) $$

where the constant does not depend on $\mathbf{x}$ and

$$ $J_{t}=\frac{1}{t}F\left[I_{m}+\frac{1}{t}F^{\top}F\right]^{-1}F^{\top}-I_{d}.$ (14) $$

The score function is thus derived as
$\nabla_{\mathbf{x}}\log p_{t}(\mathbf{x})=\frac{1}{t}J_{t}\mathbf{x}.$
It is then useful to analyze the spectrum of the matrix $J_{t}$, since $J_{t}$ is proportional to the Jacobian of the score function. In fact, since the gradient of the score is orthogonal to the manifold sufficiently close to it, the number of null eigenvalues of $J_{t}$ will correspond to the manifold dimension and we should expect to see a drop in the spectrum.

### 6.2 The linear random-manifold score in large dimensions: A random matrix analysis

In the following, we provide an outline of our theoretical results on the distribution of spectral gaps in the matrix $J_{t}$ under random linear manifolds. This choice reflects the fact that the distribution and support of the data is usually not known in advance, and it is therefore important to quantify the statistical variability induced by this uncertainty. Specifically, we will consider random projection matrices $F$ distributed according to different distributions. To ensure tractability, we perform the analysis in the limit of large $d$ (visible) and $m$ (latent) dimensions while keeping the ratio $\alpha_{m}=m/d$ constant.

### 6.3 The isotropic case

If the elements of the projection matrix $F$ are sampled as $F_{ij}\sim\mathcal{N}(0,\sigma^{2}/m)$, we are able to derive analytically the full expression of the distribution of the eigenvalues of $J_{t}$. It is given by a simple transformation of the distribution of the eigenvalues of $F^{\top}F$, which is known to be the Marchenko-Pastur distribution reported in Appendix [A.1](#A1.SS1).

Figure: (a) t = 10
Refer to caption: /html/2410.05898/assets/x2.png

In Fig. [3](#S6.F3) we show the shape of the spectrum at different times. The bulk of the distribution, inherited from the density of the eigenvalues of $F^{\top}F$, gradually shifts from left to right in the support.
By measuring the cumulative function of the spectrum, one can isolate a drop in the effective dimensionality of the manifold, as also plotted in Fig. [3](#S6.F3). The step is present at any time in the process and it is implied by the gap between the left bound of the bulk and the spike in $-1$.
The width of this gap evolves in time according to

$$ $\Delta^{\text{GAP}}_{\text{fin}}(t;\sigma)=\frac{\sigma^{2}(1+\alpha_{m}^{-1/2})^{2}}{t+\sigma^{2}(1+\alpha_{m}^{-1/2})^{2}}.$ (15) $$

If we name $\gamma_{+}(\sigma)$ the left bound eigenvalue of the bulk in the spectrum of $F^{\top}F$ (see Fig. [8](#A1.F8)), one can recover a more general expression for the gap, being

$$ $\frac{\gamma_{+}(\sigma)}{t+\gamma_{+}(\sigma)}=\Delta.$ (16) $$

Hence we can resolve the gap at a scale $\Delta$ at the time

$$ $t_{\text{in}}=\gamma_{+}(\sigma)\left(\frac{1-\Delta}{\Delta}\right).$ (17) $$

### 6.4 Intermediate gaps and subspaces with different variances

Another relevant case for our study is the one where we consider a manifold having multiple subspaces with different variances. Here we will focus on the instance of two distinct variances. This scenario is reproduced by considering a number $f\cdot m$ of columns of $F$ to have element Gaussian distributed with zero mean and variance $\sigma_{1}^{2}/m$, and the remaining $(1-f)\cdot m$ columns with elements from a Gaussian with variance $\sigma_{2}^{2}/m$. The spectrum of $J_{t}$ can be computed also in this case, as explained in Appendix [A.2](#A1.SS2).
The density function of the eigenvalues shows a transient behavior of the spectrum in the form of an intermediate drop in the estimated dimensionality of the hidden data manifold. This behavior is reported in Fig. [4](#S6.F4).
Even though the expression of the spectral density doesn’t have an explicit analytical form and has to be computed numerically, one can adopt a special assumption on the behavior of the density of the eigenvalues of $F^{\top}F$ to estimate the typical times at which the intermediate drop occurs.
Generally speaking, the spectrum of $F^{\top}F$ can be composed of two separated bulks, as observable in Fig. ([9](#A1.F9)). This happens when $\sigma_{2}^{2}$ and $\sigma_{1}^{2}$ are significantly different. In analogy with the single variance scenario, we name $\gamma_{+}$ the left bound of the bulk associated to higher eigenvalues, i.e. with the higher variance, and $\gamma_{-}$ the right bound of the bulk associated with smaller eigenvalues, i.e. smaller variance. Most commonly, $\gamma_{-}=\gamma_{-}(\sigma_{2})$ and $\gamma_{+}=\gamma_{+}(\sigma_{1})$.
In this case the gap-forming phenomenology provides for two distinct processes: the emergence of intermediate gaps (i.e. steps in the dimensionality plot) between spearated bulks of the spectrum, the opening of a final gap that allows to infer the dimensionality of the full manifold.
The width of the intermediate gap between two bulks can be obtained from Eq. ([27](#A1.E27)) as

$$ $\Delta^{\text{GAP}}_{\text{inter}}(t;\sigma_{1},\sigma_{2})=\frac{t}{t+\gamma_{-}(\sigma_{2})}-\frac{t}{t+\gamma_{+}(\sigma_{1})}.$ (18) $$

By imposing $\Delta^{\text{GAP}}_{\text{inter}}(t;\sigma_{1},\sigma_{2})=\Delta$ one finds the following quadratic form

$$ $\Delta t^{2}+\big{[}(\Delta-1)\gamma_{-}+(\Delta+1)\gamma_{+}\big{]}t+\Delta\gamma_{-}\gamma_{+}=0.$ (19) $$

Considering $\Delta\ll 1$ and $\gamma_{+}\ll\gamma_{-}$ the opening time for the intermediate gap can be found by

$$ $t_{\text{in}}(\Delta)\simeq\Delta^{-1}\gamma_{-}(\sigma_{1}),$ (20) $$

that is a reference time at which the gap becomes visible. On the other hand, by assuming the closure time to be close to zero, it can be obtained as

$$ $t_{\text{fin}}(\Delta)\simeq\Delta\gamma_{+}(\sigma_{2}).$ (21) $$

Furthermore, the time at which the gap is maximum in width, and so maximally visible, is located in between $t_{\text{in}}$ and $t_{\text{fin}}$. This is the most important time scale for the problem, it is obtained by imposing $\partial\Delta^{\text{GAP}}/\partial t=0$ and it measures

$$ $t_{\text{max}}=\sqrt{\gamma_{-}(\sigma_{1})\gamma_{+}(\sigma_{2})}.$ (22) $$

Indeed, when $\sigma_{1}^{2}\gg\sigma_{2}^{2}$ the total spectrum can be approximated by a mixture of two separated Marchenko-Pastur distributions, with variances $\sigma_{1}^{2}$ and $\sigma_{2}^{2}$, and parameters $\alpha_{m}$ and $\gamma$ to be rescaled with respect to $f$ and $(1-f)$.
This approximation becomes exact under a slight modification of $F$ which does not imply any loss of the quality of the description. Now the relevant quantities for the gap become

$$ $\displaystyle\Delta^{\text{GAP}}_{\text{inter}}(t;\sigma_{1},\sigma_{2})$ $\displaystyle=\frac{t\big{[}f\sigma_{1}^{2}(1-\sqrt{\frac{1}{f\alpha_{m}}})^{2}-(1-f)\sigma_{2}^{2}(1+\sqrt{\frac{1}{(1-f)\alpha_{m}}})^{2}\big{]}}{\Big{[}t+(1-f)\sigma_{2}^{2}(1+\sqrt{\frac{1}{(1-f)\alpha_{m}}})^{2}\Big{]}\Big{[}t+f\sigma_{1}^{2}(1-\sqrt{\frac{1}{f\alpha_{m}}})^{2}\Big{]}}$ (23) $\displaystyle t_{\text{in}}(\Delta)$ $\displaystyle\simeq\Delta^{-1}f\left(1-\sqrt{\frac{1}{f\alpha_{m}}}\right)^{2}\sigma_{1}^{2},$ (24) $\displaystyle t_{\text{fin}}(\Delta)$ $\displaystyle\simeq\Delta(1-f)\left(1+\sqrt{\frac{1}{(1-f)\alpha_{m}}}\right)^{2}\sigma_{2}^{2},$ (25) $\displaystyle t_{\text{max}}$ $\displaystyle=\sqrt{f(1-f)}\left(1-\sqrt{\frac{1}{f\alpha_{m}}}\right)\left(1+\sqrt{\frac{1}{(1-f)\alpha_{m}}}\right)\sigma_{1}\cdot\sigma_{2}.$ (26) $$

Figure: (a) t = 10
Refer to caption: /html/2410.05898/assets/x8.png

This same analysis can be extended to the more general case where the spectral density is known to be formed by different detached bulks, associated to hierarchically smaller variances of the data. The evolution of the intermediate gaps in a double-variance diffusion model is reported in Fig. [4](#S6.F4): notice that $t_{\text{max}}=\mathcal{O}(\sigma_{2})$ is consistent with Fig. [4](#S6.F4)b and [4](#S6.F4)f, where the gap was found to be maximum in width.
It is worth noting that subspaces with higher variances are the first ones to be explored by diffusion, and to be learned by the model.
This point suggests that the model is sensitive to the parameters of the probability distribution on the manifold as recently suggested by other works in the literature (see section [9](#S9) for further details).

## 7 Experiments with synthetic linear datasets

We first measure the spectrum of the singular values of the Jacobian of a score function trained through a neural network on a linear manifold data-model generated by two variances $\sigma_{1}^{2}$, $\sigma_{2}^{2}$ as described in Section [6.4](#S6.SS4). Results are reported in Figure [5](#S7.F5) (left and central panels). The opening of the gaps is consistent with the theory for the exact score: an intermediate gap associated to the subspace with higher variance first opens; subsequently, the gap relative to the lower variance subspace, which here corresponds to the final gap, opens. We can infer the dimensions of the subspaces by subtracting the location of the dashed vertical lines to $d$. We underline the fact that higher-variance subspaces are learned first by repeating the experiment after swapping the values of the variances. Eventually, the right panel in Figure [5](#S7.F5) reports the same experiments where variances are uniformly generated in the interval $[10^{-3},1]$: it is evident that the $d$ intermediate gaps is now a continuous line, as it is expected to be in more realistic natural data-sets.

Figure: Figure 5: Ordered singular values obtained with the trained score model, for difference variances on the subspaces. Data are generated according to the linear-manifold model with $d=100$ and $m=40$. Left: $\sigma_{1}^{2}=1$, $\sigma_{2}^{2}=0.01$, $f=0.75$; Center: $\sigma_{1}^{2}=0.01$, $\sigma_{2}^{2}=1$, $f=0.75$; Right: variances sampled uniformly between $10^{-3}$ and $1$.
Refer to caption: /html/2410.05898/assets/x16.png

Figure: Figure 6: Comparison between spectra obtained with the trained score model and with the numerical analysis, for difference variances on the subspaces. Data are generated according to the linear-manifold model with $d=100$ and $m=40$, $\sigma_{1}^{2}=1$, $\sigma_{2}^{2}=0.01$, $f=0.75$; from left to right, the spectra are evaluated a time $t\approx 0.45$, $t\approx 0.3$, $t\approx 0.11$.
Refer to caption: /html/2410.05898/assets/x19.png

We will now compare the gaps computed analytically with ones obtained from real neural networks trained on the same linear data-model.
The results of such comparison are presented in Figure [6](#S7.F6), and they show a good agreement between the ordered distribution of the singular values obtained through empirical methods, and the relative analytical counterpart, computed through the replica method. The opening of the predicted intermediate gaps signal the right dimension of the linear subspaces as verified from the experiments. One can notice from the figure that the analytical profile shows the shape of a sharp step between the zero value along the $x$-axis and the first appearing gap: this shape is related to the Dirac-delta spike that the spectrum of the eigenvalues presents at $r=-1$ (see Appendix [A.2](#A1.SS2) for details about the spectrum); on the other hand, the numerical profile looks different in the same region, and this behaviour is associated to the absence of the spike in the distribution of the singular values, that leaves room to a separated bulk from the other ones.
This evident discrepancy between theory and experiment is probably due to the final configuration of the trained neural network and leaves space for further investigations.

## 8 Experiments with natural image datasets

While our theoretical analysis is limited to linear random-manifold models, several qualitative figures of its phenomenology can be observed in networks trained in natural images. Fig. [7](#S8.F7) shows the temputal evolution of the spectrum estimated numerically from the Jacobian of models trained on MNIST, Cifar10 and CelebA and averaged over $100$ images. The details of this experiment are given in Supp. [C](#A3) and Supp. [B](#A2).
We see that for low times (at the end of the reverse process), the network trained on MNIST shows a sharp total manifold gap at around the $700$-th singular value. Consistently with our theory, for larger time the spectrum becomes smoother. The non averaged result obtained around single images also reveal the presence of several discernible small intermediate gaps (see left panel in Fig. [7](#S8.F7)). The results on Cifar10 and CelebA are more difficult to interpret, probably due to a less well-defined manifold and subspace structure. However, the inflection points of the spectra reveal a more suble gap structure, which recedes for small values of $t$ as more local subspaces are revealed.

Figure: Figure 7: Jacobian spectra of diffusion models trained on MNIST, Cifar10 and CelebA.
Refer to caption: /html/2410.05898/assets/x22.png

## 9 Related work & Discussion

Generative diffusion on manifold was first studies mathematically in (De Bortoli,, ) and in (Pidstrigach,, ), which provided rigorous error bounds. The use of spectral gaps to quantify the dimensionality of the manifold was introduced in (Stanczuk et al.,, ), where the total manifold gap is analyzed. Similar results were presented experimentally in (Stanczuk et al.,, ), where the authors described the phenomenology of on-manifold and on-manifold singular values in trained models and studied the deviations from the conservative vector field predicted by the exact score. These concepts were further explored theoretically in (Horvat and Pfister,, ), which offers an analysis of how violations of conservativity affect generation. Several recent studies investigated the local linear structure of trained diffusion models. For example, (Kadkhodaie et al.,, ) studied the expansion of the Jacobian of trained models and described it as an optimal geometry-adaptive Harmonic representation. Similarly, () characterized the linear expansion of the Jacobian of trained networks and characterized the resulting components in terms of their frequency content.
Our work can be seen as a theoretical complement to this more applied line of research, as we provide a comprehensive random-matrix analysis of the phenomenon in tractable models. The dynamic geometry of diffusion manifolds were also investigated in () using techniques inspired by research on latent generators such as GANs and autoencoders. Another recent work (Peng et al.,, ) uses parameterized low-rank score functions to show an equivalence between training diffusion model and subspace clustering. In our work, we prove results for the exact score for high-dimensional datasets.
Finally, (Sakamoto et al.,, ) studies the dynamic geometry of tubular neighborhoods of the latent manifold and connect their singularities with spontaneous symmetry breaking events.

In the last few years several authors employed the use of statistical physics techniques for the analysis of generative diffusion. Montanari, (), Benton et al., () and Alaoui et al., () studied the reverse generative process as a form of stochastic localization (El Alaoui and Montanari,, ; Brown et al.,, ). In parallel, Biroli and Mézard, () investigated the behavior of generative diffusion processes in the thermodynamic limit. The evolution of the fixed-points of the exact score was first studied in (Raya and Ambrogioni,, ) for the analysis spontaneous symmetry breaking and further analyzed in (Biroli and Mézard,, ) for the analysis of memorization and glassy phase transitions. The theory of critical phase transitions in the fixed-point dynamics of generative diffusion was further studied in Ambrogioni, (), where it was shown that these transitions are always of the mean-field type.

Generative diffusion models exhibit rich geometric structures that have the potential to explain their impressive generative capabilities. Our work introduces the use of random-matrix theory techniques for the analysis of their dynamic local geometry and gives room to the use of advanced statistical physics techniques, which may in the future unveil more global, topological and non-linear aspects of the dynamic geometry of diffusion generative models.

## Appendix A Analytical derivation of the Spectrum of J t subscript 𝐽 𝑡 J_{t}

### A.1 Single variance scenario

We want to compute the spectrum of the matrix in ([14](#S6.E14)). Let us first consider the case in which $F$ is a $d\times m$ matrix with Gaussian entries, and call $\gamma$ the eigenvalues of $FF^{\top}$.

The function that gives the eigenvalues $r$ of $J_{t}$ as function of $\gamma$ is

$$ $r_{j}=\frac{1}{t}\frac{\gamma_{j}}{1+\frac{1}{t}\gamma_{j}}-1=-\frac{t}{t+\gamma_{j}}$ (27) $$

Thus, knowing that the distribution of $\gamma$ is Marchenko-Pastur, we can obtain the distribution of $r$

$$ $\rho_{t}(r)=-\frac{\alpha_{m}}{2\pi}\frac{1}{r(1+r)}\sqrt{\left(r_{+}-r\right)\left(r-r_{-}\right)}+\left(1-\alpha_{m}\right)\delta\left(r+1\right)\theta\left(\alpha_{m}^{-1}-1\right)$ (28) $$

for $r\in[r_{-}(t),r_{+}(t)]$, with $r_{\pm}(t)=-\frac{t}{\left(1\pm\frac{1}{\sqrt{\alpha_{m}}}\right)^{2}+t}$.

One could ask whether the bulk of $J_{t}$ separates from $r=-1$ at a discrete time. This separation corresponds to a drop in the histogram of eigenvalues. According to Eq. ([28](#A1.E28)), the bulk is always separated from the spike at finite time $t$, because $(1+\alpha_{m}^{-1/2})^{2}+t$ for every $t$, and the width of the gap is given by $\Delta^{\text{GAP}}(t)=r_{-}(t)+1$

$$ $\Delta^{\text{GAP}}(t)=\frac{(1+\alpha_{m}^{-1/2})^{2}}{t+(1+\alpha_{m}^{-1/2})^{2}}.$ (29) $$

With respect to the starting spectrum of $F^{\top}F$, this condition reads

$$ $\frac{\gamma_{+}}{t+\gamma_{+}}=\Delta$ (30) $$

so the time when we see the drop at a scale $\Delta$ is $t=\gamma_{+}^{2}\frac{(1-\Delta)}{\Delta}$.

Figure: Figure 8: Spectrum of the eigenvalues of $F^{\top}F$ as obtained from random matrix theory with eigenvalue $\gamma_{+}$ indicated by red arrow. $\gamma_{+}$ is provided by the Marchenko-Pastur density function. Control parameters are chosen to be $\alpha_{m}=0.5$, $\sigma^{2}=1$.
Refer to caption: /html/2410.05898/assets/x25.png

### A.2 Double variance scenario

We want to compute the spectrum of $J_{t}$ when $F_{i\mu}\sim\mathcal{N}(0,\sigma_{1}^{2})$ for $\mu<fm/2$ and $F_{i\mu}\sim\mathcal{N}(0,\sigma_{2}^{2})$ for $\mu>(1-f)m/2$, with $f\in[0,1]$. We use the replica method to compute the spectrum of $A=\frac{1}{m}FF^{\top}$, then with a transform we obtain the spectrum of $J_{t}$. In order to obtain the spectrum we need to compute the expectation of the resolvent of $A$ in the $d\to+\infty$ limit, and to do this we will rely on the replica method

$$ $\displaystyle\mathbb{E}\left[g_{A}(z)\right]$ $\displaystyle=-\frac{2}{d}\frac{\partial}{\partial z}\mathbb{E}\left[\log\frac{1}{\sqrt{\det\left(zI_{d}-A\right)}}\right]$ (31) $\displaystyle=-\frac{2}{d}\frac{\partial}{\partial z}\lim_{n\to 0}\mathbb{E}\left[\frac{Z^{n}-1}{n}\right]$ (32) $$

with

$$ $\displaystyle Z^{n}$ $\displaystyle=\det\left(zI_{d}-A\right)^{-n/2}$ (33) $\displaystyle=\int\prod_{a=1}^{n}\prod_{i=1}^{d}\frac{d\phi_{i}^{a}}{\sqrt{2\pi}}\;e^{-\frac{1}{2}\sum_{a=1}^{n}\sum_{i,j=1}^{d}\phi_{i}^{a}\left(z\delta_{ij}-\frac{1}{m}\sum_{\mu}F_{i\mu}F_{j\mu}\right)\phi_{j}^{a}}$ (34) $$

and taking the expectation

$$ $\displaystyle\mathbb{E}\left[Z^{n}\right]$ $\displaystyle=\int\prod_{a,i}\frac{d\phi_{i}^{a}}{\sqrt{2\pi}}\;e^{-\frac{z}{2}\sum_{a}\sum_{i}\text{($\phi_{i}^{a})^{2}$}}\mathbb{E}\left[e^{\frac{1}{2m}\sum_{a}\sum_{\mu}(\sum_{i}\phi_{i}^{a}F_{i\mu})^{2}}\right]$ (35) $$

$$ $=\int\prod_{a,\mu}\frac{d\eta_{\mu}^{a}}{\sqrt{2\pi}}\;e^{-\frac{1}{2}\sum_{a}\sum_{\mu}\text{($\eta_{\mu}^{a})^{2}$}}\int\prod_{a,i}\frac{d\phi_{i}^{a}}{\sqrt{2\pi}}\;e^{-\frac{z}{2}\sum_{a}\sum_{i}\text{($\phi_{i}^{a})^{2}$}}\prod_{\mu}\mathbb{E}\left[e^{\frac{1}{\sqrt{m}}\sum_{a}\left(\sum_{i}\phi_{i}^{a}F_{i\mu}\right)\eta_{\mu}^{a}}\right]$ (36) $$

where in the last step we have used the independence of the rows of $F$ and applied a Hubbard-Stratonovic transform.

We can separate the product over $\mu$ and integrate over the distribution of $F$

$$ $\displaystyle\mathbb{E}\left[Z^{n}\right]$ $\displaystyle=\int\prod_{a,\mu}\frac{d\eta_{\mu}^{a}}{\sqrt{2\pi}}e^{-\frac{1}{2}\sum_{a,\mu}(\eta_{\mu}^{a})^{2}}\int\prod_{a,i}\frac{d\phi_{i}^{a}}{\sqrt{2\pi}}e^{-\frac{z}{2}\sum_{a,i}(\phi_{i}^{a})^{2}}$ (37) $\displaystyle\times\prod_{\mu=1}^{fm-1}\mathbb{E}\left[e^{\frac{1}{2\sqrt{m}}\sum_{a}\left(\sum_{i}\phi_{i}^{a}F_{i\mu}\right)\eta_{\mu}^{a}}\right]\prod_{\mu=fm}^{m}\mathbb{E}\left[e^{\frac{1}{2\sqrt{m}}\sum_{a}\left(\sum_{i}\phi_{i}^{a}F_{i\mu}\right)\eta_{\mu}^{a}}\right]$ (38) $\displaystyle=\int\prod_{a,\mu}\frac{d\eta_{\mu}^{a}}{\sqrt{2\pi}}\;e^{-\frac{1}{2}\sum_{a,\mu}(\eta_{\mu}^{a})^{2}}\int\prod_{a,i}\frac{d\phi_{i}^{a}}{\sqrt{2\pi}}e^{-\frac{z}{2}\sum_{a,i}(\phi_{i}^{a})^{2}}$ (39) $\displaystyle\times e^{\frac{\sigma_{1}^{2}}{2m}\sum_{i}\sum_{\mu<fm}\left(\sum_{a}\phi_{i}^{a}\eta_{\mu}^{a}\right)^{2}+\frac{\sigma_{2}^{2}}{2m}\sum_{i}\sum_{\mu\geq fm}\left(\sum_{a}\phi_{i}^{a}\eta_{\mu}^{a}\right)^{2}}$ (40) $\displaystyle=\int\prod_{a,\mu}\frac{d\eta_{\mu}^{a}}{\sqrt{2\pi}}\;e^{-\frac{1}{2}\sum_{a,\mu}(\eta_{\mu}^{a})^{2}}\;\int\prod_{a,i}\frac{d\phi_{i}^{a}}{\sqrt{2\pi}}\;e^{-\frac{z}{2}\sum_{a,i}(\phi_{i}^{a})^{2}})$ (41) $\displaystyle\times e^{\frac{\sigma_{1}^{2}}{2m}\sum_{ab}\left(\sum_{i}\phi_{i}^{a}\phi_{i}^{b}\right)\left(\sum_{\mu<fm}\eta_{\mu}^{a}\eta_{\mu}^{b}\right)+\frac{\sigma_{2}^{2}}{2m}\sum_{ab}\left(\sum_{i}\phi_{i}^{a}\phi_{i}^{b}\right)\left(\sum_{\mu>fm}\eta_{\mu}^{a}\eta_{\mu}^{b}\right)}.$ (42) $$

Introducing $q_{ab}=\frac{1}{d}\sum_{i}\phi_{i}^{a}\phi_{i}^{b}$

$$ $\displaystyle\mathbb{E}\left[Z^{n}\right]$ $\displaystyle=\int\prod_{a,b}\frac{dq_{ab}d\hat{q}_{ab}}{2\pi}\int\prod_{ai}\frac{d\phi_{i}^{a}}{\sqrt{2\pi}}\;e^{-\sum_{ab}\frac{1}{2}\hat{q}_{ab}\left(dq_{ab}-\sum_{i}\phi_{i}^{a}\phi_{i}^{b}\right)-\frac{z}{2}\sum_{a}\sum_{i}\text{($\phi_{i}^{a})^{2}$}}$ (43) $\displaystyle\times\left[\int\prod_{a=1}^{n}\frac{d\eta^{a}}{\sqrt{2\pi}}\;e^{-\frac{1}{2}\sum_{a}\text{($\eta^{a})^{2}$}+\frac{\sigma_{1}^{2}}{2\alpha_{m}}\sum_{ab}q_{ab}\eta^{a}\eta^{b}}\right]^{fm}$ (44) $\displaystyle\times\left[\int\prod_{a=1}^{n}\frac{d\eta^{a}}{\sqrt{2\pi}}\;e^{-\frac{1}{2}\sum_{a}\text{($\eta^{a})^{2}$}+\frac{\sigma_{2}^{2}}{2\alpha_{m}}\sum_{ab}q_{ab}\eta^{a}\eta^{b}}\right]^{(1-f)m}$ (45) $\displaystyle=\int\prod_{a,b}\frac{dq_{ab}d\hat{q}_{ab}}{2\pi}\;e^{nd\Phi(q,\hat{q})}$ (46) $$

with

$$ $\Phi(q,\hat{q})=-\frac{1}{2n}\sum_{a,b}q_{ab}\hat{q}_{ab}+G_{S}(\hat{q})+f\alpha_{m}G_{E}(q,\sigma_{1})+(1-f)\alpha_{m}G_{E}(q,\sigma_{2})$ (47) $$

where

$$ $\displaystyle G_{S}(\hat{q})$ $\displaystyle=\frac{1}{n}\log\int\prod_{a=1}^{n}\frac{d\phi^{a}}{\sqrt{2\pi}}\;e^{-\frac{z}{2}\sum_{a}(\phi^{a})^{2}+\frac{1}{2}\sum_{ab}\hat{q}_{ab}\phi^{a}\phi^{b}}$ (48) $\displaystyle G_{E}(q,\sigma)$ $\displaystyle=\frac{1}{n}\log\int\prod_{a=1}^{n}\frac{d\eta^{a}}{\sqrt{2\pi}}\;e^{-\frac{1}{2}\sum_{a}\text{($\eta^{a})^{2}$}+\frac{\sigma^{2}}{2\alpha_{m}}\sum_{ab}q_{ab}\eta^{a}\eta^{b}}$ (49) $$

Using the replica symmetric ansatz $q_{ab}=\delta_{ab}q$, $\hat{q}_{ab}=-\delta_{ab}\hat{q}$

$$ $G_{S}(\hat{q})=-\frac{1}{2}\log(z+\hat{q})$ (50) $$

$$ $G_{E}(q,\sigma)=-\frac{1}{2}\log(1-\frac{\sigma^{2}q}{\alpha_{m}}).$ (51) $$

Putting all together we have

$$ $\Phi(z)=\frac{1}{2}\hat{q}q-\frac{1}{2}\log(z+\hat{q})-f\frac{\alpha_{m}}{2}\log\left(1-\frac{\sigma_{1}^{2}q}{\alpha_{m}}\right)-(1-f)\frac{\alpha_{m}}{2}\log\left(1-\frac{\sigma_{2}^{2}q}{\alpha_{m}}\right).$ (52) $$

The integral can be evaluated by the saddle point method

$$ $\displaystyle q$ $\displaystyle=\frac{1}{z+\hat{q}}$ (53) $\displaystyle\hat{q}$ $\displaystyle=-f\frac{\alpha_{m}\sigma_{1}^{2}}{\alpha_{m}-\sigma_{1}^{2}q}-(1-f)\frac{\alpha_{m}\sigma_{2}^{2}}{\alpha_{m}-\sigma_{2}^{2}q}.$ (54) $$

We can find the Stieltjes transform

$$ $\displaystyle\mathbb{E}[g_{A}(z)]$ $\displaystyle=-2\alpha_{m}\frac{\partial\Phi(z)}{\partial z}$ (55) $\displaystyle=\alpha_{m}q^{*}(z)$ (56) $$

where $q^{*}$ is found by solving the saddle point equation

$$ $zq^{3}+q^{2}\left(\alpha_{m}-1-\frac{z\alpha_{m}}{\sigma_{1}^{2}}-\frac{z\alpha_{m}}{\sigma_{2}^{2}}\right)+\\ +q\left(\frac{\alpha_{m}^{2}}{\sigma_{1}^{2}\sigma_{2}^{2}}(z-f\sigma_{1}^{2}-(1-f)\sigma_{2}^{2})+\frac{\alpha_{m}}{\sigma_{1}^{2}}+\frac{\alpha_{m}}{\sigma_{2}^{2}}\right)-\frac{\alpha_{m}^{2}}{\sigma_{1}^{2}\sigma_{2}^{2}}=0.$ (57) $$

The asymptotic distribution of eigenvalues can be obtained from the Stieltjes transform as

$$ $\displaystyle\rho(\gamma)$ $\displaystyle=\frac{1}{\pi}\lim_{\epsilon\to 0^{+}}\mathrm{Im}\left(g_{A}(\gamma-i\epsilon)\right)$ (58) $\displaystyle=\frac{1}{\pi}\lim_{\epsilon\to 0^{+}}\mathrm{Im}\Phi^{\prime}(\gamma-i\epsilon)$ (59) $\displaystyle=\frac{1}{\pi}\alpha_{m}\lim_{\epsilon\to 0^{+}}\mathrm{Im}[q^{*}]$ (60) $$

Figure: Figure 9: Spectrum of the eigenvalues of $F^{\top}F$ as obtained from random matrix theory with eigenvalues $\gamma_{-}$ and $\gamma_{+}$ indicated by red arrows. Control parameters are chosen to be $\alpha_{m}=0.5$, $f=0.5$, $\sigma_{1}^{2}=1,\sigma_{2}^{2}=0.1$.

Once the density of the eigenvalues is computed one can perform the same change of variables described in Appendix [A.1](#A1.SS1) for the case of single variance, and obtain the density $\rho_{t}(r)$ for the eigenvalues of $J_{t}$.
Figures ([4](#S6.F4)) and ([4](#S6.F4)) report the evolution in time of the spectral density, as well as its cumulative function, when $f=0.25$ and $f=0.75$. The cumulative function has been used to estimate the formation of the intermediate gaps to compare with the experiments for the estimation of the data-manifold dimension.

## Appendix B Network training and Model Architecture Details

**Table 1: Table displaying both model and training configurations for each dataset.**
| Dataset | Image Size | Latent Dim. | Channel Mult. | Param. Count | Batch size | Iterations |
| --- | --- | --- | --- | --- | --- | --- |
| Cifar10 | 32 | 128 | (1, 2, 2, 2) | 35.7M | 128 | 500,000 |
| MNIST | 28 | 128 | (1, 2, 2) | 24.5M | 128 | 400,000 |
| CelebA | 64 | 64 | (1, 1, 2, 2, 4, 4) | 27.4M | 64 | 800,000 |

For the image datasets, we used the diffusion setting in (Ho et al.,, ). We use the variance scheduler with $\beta_{\text{min}}=10^{-4}$ and $\beta_{\text{min}}=2\times 10^{-2}$, $T=1000$ time steps, and score model backbone (PixelCNN++ (Salimans et al.,, )). Furthermore, for each of the datasets, we adjusted the partameters to account for the different complexity see [Table 1](#A2.T1).

For the linear models, we used a Variance Exploding continuous score model trained with 2M steps (batch size 128). The model had a Residual architecture with size 128 hidden channels in each layer, two residual blocks comprised by two linear layers with SiLu.

We primarily utilized NVIDIA Tesla V100 GPUs with 32 GB of memory.

## Appendix C Experimental methodology:
Computing the Singular Values of the Jacobian of the Score Function

For computing the singular values, we use the procedure from (Stanczuk et al.,, ) reported in algorithm [1](#alg1). For the linear models and MNIST models we used a symmetrized version which we empirically found to be more stable, reported in algorithm [2](#alg2).

Figure: Algorithm 1 Estimate singular values at $x_{0}$

Figure: Algorithm 2 Estimate singular values at $x_{0}$ with central difference

## Peer Reviews

confidence: 3
contribution: 3
initial_rating: 6
presentation: 3
questions: - I find the definition of the "support score", and later of the "stable latent manifolds" in Eq. (5), very confusing. Why are the manifolds defined in this way? What is the role of this definition for the rest of the paper? It seems to me that from Eq. (6) onward, this support score is not mentioned anymore and we come back to the normal score.
- Line 131: "radius $\sqrt{t}$" should be "radius proportional to $\sqrt{t}$"?
- I don't know how to reconcile the phenomenology described in Section 5 with the results presented in the later sections. Indeed, from the plots shown in Fig. 3 and 4, it seems that at all times we are in the so-called "Manifold Coverage" phase, with the intermediate gaps open. Does this mean that the "trivial phase" is only at $t=\infty$ and the "consolidation phase" is only at $t=0$? In general, I think the point of separation between these three phases is not well defined.
soundness: 3
strengths: - The paper is fairly-well written.
- The authors adopt an interesting new perspective on the study of generative diffusion processes.
summary: The paper studies the generative process of diffusion models, meaning how samples of the target distribution are generated starting from random initial values, by looking at spectral gaps in the Jacobian matrix of the score function. 

The manifold hypothesis is adopted and, in order to conduct a mathematic analysis, the authors consider the simplified case of linear manifolds and Gaussian target measures. In this setting, the temporal evolution of the spectra of the Jacobian is computed either analytically or through a Replica computation, and the resulting spectral gaps are used to estimate the dimensionality of the data-manifold.

The main result of the paper is that, during the generative process, there is a series of spectral gaps opening, each corresponding to a different local subspace (this phase is denoted Manifold Coverage). Then, the more the process gets close to t=0, the more these gaps shrink, except from the principal one, which is the only one sharpening and remaining asymptotically.

The authors verify this theoretical prediction on synthetic datasets, by learning the score with a Neural Network and looking at the empirical distribution of the respective Jacobian.

Finally, some experiments on image datasets are presented, to look at whether the phenomenology described for linear manifolds is valid, in some way, also for more realistic settings.
weaknesses: - Section 7: The justification of the discrepancy in Fig. 6 between the analytical (replica) and the empirical results is weak, in my opinion. If the Neural Network is expressive enough to learn the exact score, shouldn't we expect the two curves to overlap? If one does the same experiments for the isotropic case for which an analytic solution is available, does one get the same discrepancy?
- Section 8: The results for Cifar10 and CelebA do not show the appearance of any spectral gap. I didn't understand the justification at the end of Section 8. Also for MNIST, while a spectral gap seems to appear at small values of $t$, it seems to indicate a dimensionality of the manifold of around 700. How is this compatible to the fact that the dataset considers only 10 classes?
- [minor]: In the Introduction, the Jacobian is mentioned, but it should be said it is the Jacobian of the score.

The paper presents many typos, here is a list of the ones I found (there could be many more):
- Line 40: models should be modes
- Line 66: idealize should be idealized
- Line 113: is is
- Line 126: Eq. (2) should be Eq. (1)
- Line 229: the the
- Line 357: spearated 
- Line 501: temputal
- Line 510: suble

I suggest a careful reading to correct any additional typo present.

---

confidence: 4
contribution: 2
initial_rating: 5
presentation: 2
questions: - The experimental detail of section 7 is very light and lacks much detail needed for recreation. For example the number of samples from the distribution from training is missing.
- Section 7 leaves a poignant question unexplored - how does the number of samples from the data distribution used to train the score function affect the empirical results in estimating the dimension of the data manifold. The distribution of eigenvalues obtained in the results appear to not mach perfectly the theoretical results. Is this due to modeling issues, the amount of data used, or some other issues? Exploring a range of dataset sizes may shed some light on why the image experiments fail to show the definitive eigenspectra cliff needed to detect the intrinsic data dimensionality. Perhaps these datasets simply don't contain enough examples to learn the score well.
- If the conclusion of section 8 is that for real datasets the proposed method of looking at the distribution of eignevalues of the score function fails to detect the intrinsic dimensionality of the data, can the authors comment on how this menthod might be of practical utility in training or understanding diffusion models?
soundness: 3
strengths: - The results are of some interest and validate observational results of the change in distribution of the eigen spectra of the score function of diffusion models through time.
summary: The paper provides a random-matrix theory derived result computing the exact eigenspectra of the score of a diffusion modeling a distribution supported on a linear submanifold of a sample space. This distribution is a Gaussian supported on a hyperplane in Euclidean space. A three phase process on the eigenvalues is observed in the diffusion placed on the modeled density. This theoretical computation matches the observation of a number of previous works studying the eigenspectra of the score function of diffusion models.
weaknesses: - The setting studied is a rather limited setting, and given the linear nature of the submanifold the density is supported on, and the Gaussian nature of the density.
- The experimental section of the paper is very light and leaves many questions unexplored.

---

confidence: 4
contribution: 3
initial_rating: 6
presentation: 3
questions: Following the above weakness, I have the following questions. 

1. Can the authors either show some results for nonlinear manifold? If not, at least certain discussion is expected.
2. Can the authors further discuss the results on CIFAR-10 and CelebA? These results are not sufficient to support the claim in the paper, nor the conjecture. 
3. This is not really a question. At a first glance I thought this is a theoretical paper, but after reading it carefully I feel it's more like a summary and conjecture based on empirical observation, as the theoretical part is only restricted to linear models.
soundness: 3
strengths: 1. The perspective on generative diffusion models, focusing on the spectrum of the score function's Jacobian, appears novel. 

2. The three-phase framework in generative diffusion models provides insight into model behavior across different generative stages.

3. The approach is supported by mathematical derivations under a linear assumption and is further validated with illustrative examples on toy datasets like MNIST.
summary: This paper provides a theoretical perspective on generative diffusion models under the manifold hypothesis, focusing on the spectrum of the score function's Jacobian. The authors propose that spectral gaps indicate the presence and dimensionality of distinct sub-manifolds. They derive analytical forms of the spectrum under the assumption that the manifold is a linear subspace, revealing three phases in the generative process: a trivial phase, a manifold coverage phase, and a consolidation phase. THey argue that this division explains why generative models avoid the manifold overfitting issues seen in likelihood-based models.
weaknesses: 1. The linear assumption is overly restrictive and unlikely to hold in practical scenarios, limiting the generalizability of the theoretical results.
2. The real data results on CIFAR-10 and CelebA (Figure 7) are difficult to interpret and lack sufficient clarity to be fully convincing.
3. In Section 5, the conjecture that "their phenomenology captures the main features of subspace separation in the tangent space of curved manifolds" lacks a clear intuition or justification, making it unconvincing based on the current presentation.

A minor issue:
There are some citation formatting issues. For example. the third reference in the first paragraph of the intro, should be Song et al, instead S et al. See also the B. et all in the same paragraph, and so on.

---

confidence: 3
contribution: 2
initial_rating: 6
presentation: 2
questions: * In section "The linear random-manifold score in large dimensions: A random matrix analysis", is the approach aimed to study the behaviour of diffusion models when encountering random datasets?
* In Equation (5), $M_t$ appears to represent the manifold where noised data reside. Is the idea that particles are found in shells of radius $\sqrt{t}$ around the latent manifold consistent with this? If $M_t$ is not merely a union of spheres centered on data points (I don't think I believe it is the case), what precisely is its mathematical function? I understand you want to consider a specific space on which the noised data lies, in order to use the intepretations deriving from (6), but I don't think introducing $M_t$ is a good mathematical strategy
soundness: 3
strengths: * Clear exposition, good fit in the existing literature
* The analysis of the Jacobian spectrum and its phase-wise evolution provides a meaningful interpretation of the manifold structure and how diffusion models capture it
* Using linear manifolds and random matrices enables tractable and insightful theoretical analysis supporting the findings.
* Experiments on complex datasets, such as MNIST and CIFAR-10, align with theoretical predictions, showing the three observable phases.
summary: This paper investigates the spectral properties of generative diffusion models within the framework of the manifold hypothesis. By analyzing the spectrum of the Jacobian of the score function, the authors identify and interpret three distinct phases in the generative process: trivial phase, manifold coverage phase, and manifold consolidation phase. These phases explain the internal distribution alignment and resilience against manifold overfitting, which is often seen in likelihood-based models. Through a mix of theoretical analysis using linear manifolds and random matrix theory, alongside high-dimensional experiments on image datasets, the paper demonstrates how these phases manifest and validate the theory.
weaknesses: * Introducing the mathematics for manifolds earlier could clarify the manifold overfitting problem and improve the logical flow.

See questions

---

confidence: 3
contribution: 2
initial_rating: 5
presentation: 1
questions: 1. (044-045) These decomposition phenomena cannot be directly explained in terms of critical phase transitions as they are fundamentally linear processes. 
   - This reasoning is not clear to me. 
   - I assume authros wanted to explain the separation of subspaces in diffusion models using the phase transition idea, much like how B. et al. 2024a explained the separation of data classes. 
   - what's the main difference of this work from B. et al. 2024a?
2. 091-095 meantioned that the divergence of data distribution implies the model density becomes larger during training. Could authros elaborate on this? Moreover, when claiming diffusion models having three phases are unaffected by manifold overfitting, could authors provide the analysis accordingly in comparison to the likelihood based models? 
3. 096-097: it's pointed out in this work, diffusion models are not affected by the manifold overfitting. Could authors provide some evidence on this, e.g., literatures? 
4.  140-146: it's not straightforward to see how perturbations aligned to the tangent space of the manifold correcpond to small eigenvalues, while the orthogonal ones correspond to high eigenvalues. Could authors provide more details on this?
5. It's hard to parse Fig 1. Is this at a fixed time t?  I checked the cited paper Stanczuk et al. 2022, but it's still unclear to me. Why the x-axis is latent dimensionality? I do not directly see the link that perturbations more aligned to the tangent space of the manifold are located at the regim of smaller latent dimentionality, not the orthogonal ones have larger dimensionality. How do you find the one-to-one correspondence between the eigenvalues and the latent dimensionality? In Stanczuk et al. 2022, the authors pointed out that the vanishing singular values correspond to the tangent space, while the large ones correspond to the normal space. But it seems there is some nontrivial steps needed from there to here. 
7. 172: The phase separation do not correspond to singularities as there are cross-over events, not genuine phase transitions. 
   1. what are singularities in this context? 
   2. what are cross-over events?
8. It seems $\sigma_k$ is not defined. I assume it is the variance of the data distribution over manifolds with different dimensions. When would $\sigma_{k+1}^2\ll \sigma_k^2$ happen? 
9.  What is asymptotic closure? Do you mean the end time of the phase II can be infinite and one considers an approximate of when this closure happens? 
10. Please also define $r$ in line 481.

Other suggestions:
- I think the template of the paper is different from the provided one.
- Please correct and unify the bibliography.
soundness: 2
strengths: I think this paper is the first work using random matrix to analyze the jacobian in diffusion models, though the idea of analyzing high-dim jacobians using random matrix theory is not new, assuming I read enough literature on understanding diffusion models. Understanding diffusion models, especially from the latent space geometry, is increasingly researched and this paper provides a new way of approaching this analysis.
summary: This paper aims to understand the diffusion models from the geometry of the latent space. The authors did analysis in terms of the spectrum of the score functions' Jacobian, from which the spectral gaps indicates three learning phases of the models. Further, the authors reasoned that the fitting of the distribution and the manifold geometry are disentangled so to explain the indifference of diffusion models to the manifold overfitting.
weaknesses: 1. The paper is delivered in a rather vague way and there are multiple important points which should be clarified in terms of either equations, or citations on other literatures. When carefully reading the paper, I couldn't really parse the ideas at several important locations, though some of these unclear points are only explained later on in a more concrete way in Section 6. I have detailed a few of them in the questions section. 
2. The paper strives to better understand the generative diffusion models. While I appreciate the effort using RMT (random matrix theory), the analysis does not really offer deep insights to diffusion models.  In the practical diffusion models, the analysis in Figure 7 merely demonstrated there is a eigenspectrum transition in the Jacobian of the trained model on MINST, which does not have an apparent link to the three phases proposed in the paper. Similar observations have already been shown by _Diffusion Models Encode the Intrinsic Dimension of Data Manifolds_. It would be great if authors could point out the difference or novel insights compared to this work. 
3. The paper _Dynamical Regimes of Diffusion Models_ showed the three phases of the generative process in diffusion models using phase transition analysis where these key transitions are also analyzed via the distributional solution of the forward dynamics. These transitions are determined by not only the eigenvalues of the covariance matrix of the data, but also the entropy of the noised data distribution. It really interests me to see the difference between the analysis in this work and the one in _Dynamical Regimes of Diffusion Models_. So it's necessary to point out the new results from this paper compared to the previous analyses.
4. I try to locate in the paper where these different phases explains why diffusion models are prone to manifold overfitting. Unfortunately, I could not find it. It's mentioned in section 2 that likelihood based models lead to increasing model density during training (091-092), thus, resulting in the manifold overfitting. Could authors explain how the three phases of diffusion models help them to be unaffected by the manifold overfitting from this perspective? For example, add a section or paragraph that explicitly connects their three-phase model to the issue of manifold overfitting, explaining step-by-step how each phase contributes to avoiding this problem

Meta-review: Summary:
This paper presents a novel theoretical analysis of generative diffusion models through the lens of differential geometry and random matrix theory. By examining the spectral gaps in the Jacobian of the score function, the authors identify three distinct phases of generative diffusion: a trivial phase, a manifold coverage phase, and a consolidation phase. This work provides valuable insights into the generative process, particularly in understanding how diffusion models avoid manifold overfitting—a common issue in likelihood-based models. The theoretical analysis is supported by experiments on synthetic and real-world datasets.

Strengths:

Novel Perspective: The use of random matrix theory to analyze the Jacobian in generative diffusion models is novel, contributing a fresh theoretical angle to the field.

Insightful Contributions: The identification of the three geometric phases offers an explanation for why diffusion models do not suffer from manifold overfitting, which has been a critical issue in generative modeling.

Rigorous Theoretical Foundation: The paper develops analytical tools for the linear manifold case, providing a clear and elegant connection between theory and empirical observations.

Improved Clarity and Experiments: Following reviewer feedback, the authors significantly revised the paper, clarifying key points, addressing methodological concerns, and expanding the discussion of real-world experiments.

Weaknesses:

Simplified Assumptions: The linear manifold assumption limits the generality of the theoretical results. While the authors provide justification and some experiments to extend the findings to nonlinear manifolds, these remain conjectural.

Experimental Limitations: Some reviewers found the experiments on real-world datasets to be insufficiently linked to the theoretical claims, with weaker evidence of the three-phase structure in CIFAR-10 and CelebA datasets.

Presentation Concerns: The initial presentation of the paper was noted to be fragmented and difficult to follow. While the revised version addresses many of these concerns, there is room for further improvement in connecting sections and streamlining definitions.

Discussion:

The reviewers largely agree on the novelty and value of the theoretical contributions, particularly the identification of the three phases of generative diffusion and the insights into manifold overfitting. While there are concerns about the linearity assumption and the real-world experimental evidence, the authors have made substantial improvements to the paper during the discussion phase, addressing many of these criticisms. The analysis has significant potential to inspire further research on the geometry of diffusion models and their generative properties.

Suggestions for Improvement Before Camera-Ready Submission:

Further refine the connections between sections to improve the paper's readability and logical flow.

Expand the discussion of how the insights from this paper could be used to improve training or evaluation of diffusion models in practice.

Include additional commentary on the limitations of the linear manifold assumption and future directions for extending the theory to nonlinear manifolds.

Conclusion:
This paper provides an innovative and rigorous theoretical framework for understanding the geometric properties of generative diffusion models. Despite some limitations, its contributions are significant and warrant publication. The constructive feedback from reviewers has been well-addressed, and the revised manuscript is a valuable addition to the field.

Decision: Accept (Poster)