The adaptive interpolation method: A simple scheme to prove replica formulas in Bayesian inference

TL;DR

This paper introduces the adaptive interpolation method, a straightforward proof scheme extending Guerra and Toninelli's approach, to rigorously establish replica formulas in various Bayesian inference problems, simplifying previous complex proofs.

Contribution

The paper presents a new, simpler proof technique called the adaptive interpolation method for validating replica formulas in Bayesian inference, applicable to multiple complex problems.

Findings

Proved replica formula for symmetric rank-one matrix estimation.

Extended the method to tensor and linear estimation problems.

Provided insights into the validity of replica formulas in Bayesian inference.

Abstract

In recent years important progress has been achieved towards proving the validity of the replica predictions for the (asymptotic) mutual information (or "free energy") in Bayesian inference problems. The proof techniques that have emerged appear to be quite general, despite they have been worked out on a case-by-case basis. Unfortunately, a common point between all these schemes is their relatively high level of technicality. We present a new proof scheme that is quite straightforward with respect to the previous ones. We call it the adaptive interpolation method because it can be seen as an extension of the interpolation method developped by Guerra and Toninelli in the context of spin glasses, with an interpolation path that is adaptive. In order to illustrate our method we show how to prove the replica formula for three non-trivial inference problems. The first one is symmetric rank-one matrix estimation (or factorisation), which is the simplest problem considered here and the one for which the method is presented in full details. Then we generalize to symmetric tensor estimation and random linear estimation. We believe that the present method has a much wider range of applicability and also sheds new insights on the reasons for the validity of replica formulas in Bayesian inference.

Equations

$$w_{ij}=\frac{s_i{s}_j}{\sqrt{n}}+z_{ij}\sqrt{\Delta }\text{for}1\le i\le j\le n,$$

$$\displaystyle P({\mathbf{w}}|{\mathbf{s}})=\frac{\exp\Big{\{}-\frac{1}{2\Delta}\sum_{i\leq j}\big{(}w_{ij}-\frac{s_{i}s_{j}}{\sqrt{n}}\big{)}^{2}\Big{\}}}{(2\pi\Delta)^{\frac{n(n+1)}{2}}}.$$

$$\displaystyle P({\mathbf{x}}|{\mathbf{w}})=\frac{\prod_{i=1}^{n}P_{0}(x_{i})P({\mathbf{w}}|{\mathbf{x}})}{\int\big{\{}\prod_{i=1}^{n}dx_{i}P_{0}(x_{i})\big{\}}P({\mathbf{w}}|{\mathbf{x}})}.$$

$$\displaystyle P\Big{(}{\mathbf{x}}\Big{|}{\mathbf{w}}=\frac{{\mathbf{s}}{\mathbf{s}}^{\intercal}}{\sqrt{n}}+{\mathbf{z}}\sqrt{\Delta}\Big{)}=\frac{\prod_{i=1}^{n}P_{0}(x_{i})e^{-{\cal H}({\mathbf{x}};{\mathbf{s}},{\mathbf{z}})}}{\int\big{\{}\prod_{i=1}^{n}dx_{i}P_{0}(x_{i})\big{\}}e^{-{\cal H}({\mathbf{x}};{\mathbf{s}},{\mathbf{z}})}}$$

$$\displaystyle{\cal H}({\mathbf{x}};{\mathbf{s}},{\mathbf{z}})\vcentcolon=\frac{1}{\Delta}\sum_{i\leq j=1}^{n}\Big{(}\frac{x_{i}^{2}x_{j}^{2}}{2n}-\frac{x_{i}x_{j}s_{i}s_{j}}{n}-\frac{x_{i}x_{j}z_{ij}\sqrt{\Delta}}{\sqrt{n}}\Big{)}$$

$$\displaystyle{\cal Z}({\mathbf{s}},{\mathbf{z}})\vcentcolon=\int\big{\{}\prod_{i=1}^{n}dx_{i}P_{0}(x_{i})\big{\}}e^{-{\cal H}({\mathbf{x}};{\mathbf{s}},{\mathbf{z}})}.$$

$$f_n:=-\frac{1}{n}{\mathbb{E}}_{\mathbf{S},\mathbf{Z}}[\ln \mathcal{Z}(\mathbf{S},\mathbf{Z})]$$

$$\displaystyle f_{\rm RS}(m;\Delta)\vcentcolon=\frac{m^{2}}{4\Delta}+f_{\rm den}\big{(}\Sigma(m;\Delta)\big{)},$$

$$\Sigma (m;\Delta ):=\sqrt{\frac{\Delta }{m}}.$$

$$\displaystyle f_{\rm den}(\Sigma)\vcentcolon=-\mathbb{E}_{S,\widetilde{Z}}\Big{[}\ln\int dxP_{0}(x)e^{-\frac{1}{\Sigma^{2}}\big{(}\frac{x^{2}}{2}-xS-x\widetilde{Z}\,\Sigma\big{)}}\Big{]}.$$

$$\underset{n\to \infty }{\lim }{f}_n=\underset{m\ge 0}{\min }{f}_{\mathrm{RS}}(m;\Delta ).$$

$$\frac{I(\mathbf{S};\mathbf{W})}{n}=f_n+\frac{\mathbb{E}[S^2{]}^2}{4\Delta }+\mathcal{O}(n^{-1}),$$

$$I(\mathbf{S};\mathbf{W})=I(\mathbf{S};\mathbf{S}{\mathbf{S}}^{⊺}/\sqrt{n}+\mathbf{Z}\sqrt{\Delta })+\mathcal{O}(\sqrt{n}),$$

$${\Delta }^{-1}:=\int dw{P}_{\mathrm{out}}(w|0)({\partial }_y\ln P_{\mathrm{out}}(w|y){|}_{y=0}{)}^2.$$

$${\mathcal{H}}_{k,t}(\mathbf{x}):=$$

$$\displaystyle\qquad+h\Big{(}{\mathbf{x}},{\mathbf{s}},{\mathbf{z}}^{(k)},\frac{K\Delta}{1-t}\Big{)}+h_{\rm mf}\Big{(}{\mathbf{x}},{\mathbf{s}},\widetilde{{\mathbf{z}}}^{(k)},\frac{K\Delta}{t\,m_{k}}\Big{)},$$

$$h(\mathbf{x},\mathbf{s},\mathbf{z},{\sigma }^2)$$

$$h_{\mathrm{mf}}(\mathbf{x},\mathbf{s},\tilde{\mathbf{z}},{\sigma }^2)$$

$$\displaystyle\Big{\{}$$

$${\mathbf{w}}^{(k)}=\frac{\mathbf{s}{\mathbf{s}}^{⊺}}{\sqrt{n}}+{\mathbf{z}}^{(k)}\sqrt{\frac{K\Delta }{1-t}},$$

$${\mathbf{y}}^{(k)}=\mathbf{s}+{\tilde{\mathbf{z}}}^{(k)}\sqrt{\frac{K\Delta }{t{m}_k}}.$$

$$\displaystyle{\cal H}_{k,t;\epsilon}({\mathbf{x}};\bm{\theta})\vcentcolon={\cal H}_{k,t}({\mathbf{x}})+\epsilon\sum_{i=1}^{n}\Big{(}\frac{x_{i}^{2}}{2}-x_{i}s_{i}-\frac{x_{i}\widehat{z}_{i}}{\sqrt{\epsilon}}\Big{)},$$

$$P_{k,t;ϵ}(\mathbf{x}|\mathbf{\theta })$$

$$⟨A(\mathbf{X}){⟩}_{k,t;ϵ}$$

$$f_{k,t;ϵ}$$

$$|f_{1,0;ϵ}-f_{1,0;0}|\le \frac{ϵ}{2}\mathbb{E}[S^2],\text{and}|f_{K,1;ϵ}-f_{K,1;0}|\le \frac{ϵ}{2}\mathbb{E}[S^2].$$

$${\mathcal{H}}_{1,0;0}(\mathbf{x};\mathbf{\theta })=$$

$$\displaystyle f_{1,0;0}=-\frac{1}{n}\mathbb{E}_{{\mathbf{S}},{\mathbf{Z}}}\Big{[}\ln\int\big{\{}\prod_{i=1}^{n}dx_{i}P_{0}(x_{i})\big{\}}e^{-{\cal H}_{1,0;0}({\mathbf{x}};{\mathbf{S}},{\mathbf{Z}})}\Big{]}$$

$$f_{1,0;0}=f_n.$$

$${\mathcal{H}}_{K,1;0}(\mathbf{x};\mathbf{\theta })=$$

$$m_{\mathrm{mf}}^{(K)}:=\frac{1}{K}\sum _{k=1}^K{m}_k.$$

$$\displaystyle{\cal H}_{K,1;0}({\mathbf{x}};\bm{\theta})=\frac{m_{\rm mf}^{(K)}}{\Delta}\sum_{i=1}^{n}\Big{(}\frac{x_{i}^{2}}{2}-x_{i}s_{i}-x_{i}\sqrt{\frac{\Delta}{m_{\rm mf}^{(K)}}}\sum_{k=1}^{K}\widetilde{z}_{i}^{(k)}\sqrt{\frac{m_{k}}{Km_{\rm mf}^{(K)}}}\Big{)}.$$

$$f_{K,1;0}$$

Full text

†† ​​​​​​​​​​$\dagger$ Laboratoire de Théorie des Communications, Faculté Informatique et Communications, Ecole Polytechnique Fédérale de Lausanne, CH-1015, Suisse.

$\star$ International Center for Theoretical Physics, Strada Costiera, 11 I, 34151 Trieste, Italy.

1. 1 Introduction

2. 2 The adaptive interpolation method: Main ideas

3. 2.1 Symmetric rank-one matrix estimation: Setting and main result

4. 2.2 The $(k,t)$–interpolating model

5. 2.3 The initial and final models

6. 2.4 Free energy change along the adaptive interpolation path

7. 2.5 Upper bound

8. 2.6 Lower bound

9. 2.7 Proof of the fundamental sum rule

10. 3 Application to rank-one symmetric tensor estimation

11. 3.1 Symmetric rank-one tensor estimation: Setting and main result

12. 3.2 Sketch of proof of the replica symmetric formula

13. 4 Application to Gaussian random linear estimation

14. 4.1 Gaussian random linear estimation: Setting and result

15. 4.2 Proof of the RS formula

16. 5 Concentration of overlaps

17. 6 A fluctuation identity

18. 7 Concentration of the free energy

19. A Linking the perturbed and plain free energies

20. B Proof of Lemma 3

21. C Alternative argument for the lower bound

22. D A consequence of Bayes rule

23. E A stochastic calculus interpretation

1  Introduction

A very interesting development in probability theory in recent years has been the progress on a coherent mathematical theory [1, 2, 3, 4] of the predictions of the replica and cavity methods [5] in statistical physics of spin glasses. In this respect one of the most important tools is the invention of the interpolation method by Guerra and Toninelli [6, 7] which eventually led Talagrand to a remarkable proof [8] of the Parisi formula [9] for the free energy of the Sherrington-Kirkpatrick model [10].

In more recent years the interpolation method has been fruitfully extended and adapted to problems of interest in a wide range of applications such as in coding theory, communications, signal processing and theoretical computer science, well beyond the realm of traditional statistical mechanics. Among these we highlight applications of the interpolation method to error correcting codes [11, 12, 13, 14], random linear estimation and compressive sensing [15, 16, 17, 18, 19], low-rank matrix and tensor factorization [20, 21, 22] and constraint satisfaction problems [23, 24, 25, 26]. Most of these problems are inference problems and when a Bayesian framework is adopted, they can be solved with a replica symmetric scheme (constraint satisfaction is not, as such at least, an inference problem and does not fall in this category). The replica symmetric formulas for the free energies, mutual informations and error performance measures typically predict interesting first order phase transitions, with associated “metastable states with infinite lifetime”, which pose interesting algorithmic challenges of great importance in practical applications as well as challenges from the analysis point of view. It has turned out that one can learn a great deal about the fundamental limitations for important classes of (message-passing) algorithms by studying these replica solutions (we refer to [27] for a general reference and come back to this point in the conclusion).

In spite of their complexity, for all the inference problems cited above, complete proofs of the replica symmetric formulas have been found. These proofs usually combine Guerra-Toninelli interpolation bounds with some other non-trivial idea or method, namely algorithmic approaches involving so-called spatially coupled models [28, 22, 17, 18], information theoretic methods [29, 30] or rigorous versions of the cavity method [31, 32, 33, 34] using the Aizenman-Sims-Starr principle [35]. While each of these methods has its own merit and sheds interesting light, they all lead to quite long and technically involved proofs. Besides, although each method can probably be taylored for each problem, it would clearly be more satisfactory to have a more or less unified approach.

In this paper we develop a new unified and self-contained interpolation method. We illustrate how it works for three different problems, namely rank-one symmetric matrix and tensor factorization, as well as random linear estimation and compressive sensing. Our method allows to prove at the same time matching lower and upper bounds on the free energy with much less effort than all known current proofs. All these problems are “spin systems” defined for “dense graphs” (complete graphs or hypergraphs). The ideas of this paper can also be adapted to error correcting codes that are akin to spin systems on “sparse” random graphs and we plan to come back to this aspect elsewhere111Since the first version of this manuscript, the method has been successfully applied to many other problems including non-symmetric matrix and tensor factorization [36], generalized linear models and learning [37], models of deep neural networks [38, 39], random linear estimation with structured matrices [40] and even problems defined by sparse graphical models such as the censored block model [41]..

Roughly speaking, our new scheme interpolates between the original problem and the mean-field replica solution in small steps, each step involving its own set of trials parameters and Gaussian mean-fields in the spirit of Guerra and Toninelli (this idea of interpolating in small steps originated in the sub-extensive interpolation method developed by the authors in [18, 19]). We are then able to choose the set of trial parameters in various ways so that we get both upper and lower bounds that eventually match. One can interpret the set of trial parameters as a suitable “interpolation path” that we “adapt” to obtain suitable bounds, and thus we call this method the adaptive interpolation method.222In the present formulation one can also interpret the succession of Gaussian mean-fields in each step as a Wiener process. For this reason we initially called this new approach “the stochastic interpolation method”. The interpretation in terms of a Wiener process is in fact not really needed, and here we choose a more pedestrian path, but we believe this is an aspect of the method that may be of further interest (specially for diluted systems) and briefly discuss it in Appendix E.

An important aspect of our method is the need for concentration properties of the suitable “overlap”. It was already proven long ago in [42, 43] that a concentration hypothesis for overlaps implies that the replica symmetric solution is exact (an implication that was known to physicists). However for typical spin glass systems (e.g. the Sherrington-Kirkpatrick or $p$-spin spin glass) this hypothesis can only hold in some high temperature phase, and it is also difficult to prove. We refer to [44, 43] and [1] for pioneering works on such proofs with the help of cavity-like methods. In the framework of Bayesian inference the situation is more favourable. The Bayes rule immediately implies a special set of identities obeyed by suitable “correlation functions” often known as Nishimori identities [45, 46]. These identities then allow to deduce the concentration of overlaps from the concentration of the free energy in the whole phase diagram. This is also the reason why Bayesian inference problems generally lead to replica symmetric solutions.

The paper is organized as follows. Section 2 gives a pedagogic introduction to the adaptive interpolation method for one of the simplest, yet non-trivial problems, namely rank-one symmetric matrix factorization. The replica symmetric formula for the free energy or mutual information is completely proven in a self contained and direct way (see Theorem 1). As explained in the previous paragraph, for all these problems our analysis also rests on concentration properties of the overlap parameters in the whole phase diagram (Lemma 2). This analysis is the subject of sections 5, 6 and 7, and can be read independently from the rest of the paper. We then sketch the same method for symmetric tensors (see Theorem 2). Section 4 presents the method for a more difficult problem, namely random linear estimation. In particular, we provide a much simpler and transparent proof than all other existing proofs [17, 18, 29, 30] of the replica formula (see Theorem 3).

2  The adaptive interpolation method: Main ideas

Before starting let us introduce a few notations used all along this paper: Vectorial quantities will be denoted by boldface letters, random variables by capital letters and their realizations by small letters. Expectations with respect to “quenched” variables (i.e. the variables that are fixed by the realization of the problem) are denoted $\mathbb{E}$ and those with respect to “annealed” variables (i.e. the dynamical variables) are denoted by Gibbs brackets $\langle-\rangle$ possibly with appropriate subscripts. This choice follows the standards of statistical mechanics.

2.1 Symmetric rank-one matrix estimation: Setting and main result

Consider the following probabilistic rank-one matrix estimation problem: One has access to noisy observations ${\mathbf{w}}=[w_{ij}]_{i,j=1}^{n}$ of the pair-wise product of the components of a vector ${\mathbf{s}}=[s_{1},\dots,s_{n}]^{\intercal}\in\mathbb{R}^{n}$ with i.i.d components distributed as $S_{i}\sim P_{0}$, $i=1,\dots,n$ (that we simply denote {\mathbf{S}}\overset{\text{\tiny i.i.d.}}{\mathrel{\raisebox{-2.0pt}{\sim}}}P_{0}). A standard and natural setting is the case of additive white Gaussian noise of known variance $\Delta$,

[TABLE]

where ${\mathbf{z}}=[z_{ij}]_{i,j=1}^{n}$ is a symmetric matrix with i.i.d entries $Z_{ij}\sim\mathcal{N}(0,1)$ for $1\leq i\leq j\leq n$. This is denoted {\mathbf{Z}}\overset{\text{\tiny i.i.d.}}{\mathrel{\raisebox{-2.0pt}{\sim}}}{\cal N}(0,1). The goal is to estimate the ground truth ${\mathbf{s}}$ from ${\mathbf{w}}$ assuming that both $P_{0}$ and $\Delta$ are known and independent of $n$ (the noise is symmetric so that $w_{ij}=w_{ji}$).

We consider a Bayesian setting and associate to the model (1) its posterior distribution. The likelihood of the (component-wise independent) observation matrix ${\mathbf{w}}$ given ${\mathbf{s}}$ is

[TABLE]

From the Bayes formula we then get the posterior distribution333We abusively use the notation $dxP_{0}(x)$ even though $P_{0}$ is not necessarily absolutely continuous. for ${\mathbf{x}}=[x_{1},\ldots,x_{n}]^{\intercal}\in\mathbb{R}^{n}$ given the observations (it is convenient to explicitely distinguish between the ground truth signal vector ${\mathbf{s}}$ and its estimate ${\mathbf{x}}$ sampled from the posterior)

[TABLE]

Replacing the observation ${\mathbf{w}}$ by its explicit expression (1) as a function of the signal and the noise we obtain

[TABLE]

where we call

[TABLE]

the Hamiltonian of the model. In order to obtain the last form of the posterior distribution we replaced $w_{ij}$ using (1), developed the square in $P({\mathbf{w}}|{\mathbf{x}})$, and simplified the ${\mathbf{x}}$-independent terms in the numerator and denominator. The normalization factor is by definition the partition function

[TABLE]

Our principal quantity of interest is the average free energy per component444For all other models considered in this paper we directly write the explicit expression of the free energy, but the derivation is always similar. defined by

[TABLE]

where {\mathbf{S}}\overset{\text{\tiny i.i.d.}}{\mathrel{\raisebox{-2.0pt}{\sim}}}P_{0} and {\mathbf{Z}}\overset{\text{\tiny i.i.d.}}{\mathrel{\raisebox{-2.0pt}{\sim}}}{\cal N}(0,1).

Define the replica symmetric (RS) potential $f_{\rm RS}(m;\Delta)$ as

[TABLE]

with

[TABLE]

Here $f_{\rm den}(\Sigma)$ is the free energy associated with a scalar Gaussian denoising model: $y=s+\widetilde{z}\,\Sigma$ where $S\sim P_{0}$, $\widetilde{Z}\sim{\cal N}(0,1)$. The free energy $f_{\rm den}(\Sigma)$ is minus the average logarithm of the normalization of the posterior distribution $P(x|s+\widetilde{z}\,\Sigma)\propto\exp(-\Sigma^{-2}(x^{2}/2-xs-x\widetilde{z}\,\Sigma))P_{0}(x)$:

[TABLE]

Our first theorem illustrating the adaptive interpolation method is

Theorem 1** (RS formula for symmetric rank-one matrix estimation).**

Fix $\Delta>0$. For any $P_{0}$ with bounded support, the asymptotic free energy of the symmetric rank-one matrix estimation model (1) verifies

[TABLE]

Proof.

The theorem follows from Proposition 1 in section 2.5 and Proposition 2 in section 2.6. ∎

The bounded support property hypothesis for $P_{0}$ is not really a requisite of the adaptive interpolation method, but simply makes the necessary concentration proofs for the free energy simpler. There is no condition on the size of the support, and it is presumably possible to take a support equal to the whole real line by a limiting process applied to (11), as long as the first four moments of $P_{0}$ are finite.

Formulas such as (11), where a complicated statistical model is related to a scalar (and thus analyzable) statistical model are at the root of the mean-field theory in statistical mechanics. A possible intuition behind this formula (and all formulas of the same type in this article) is as follows: The estimation problem (1) is effectively “replaced” by a decoupled estimation model ${\mathbf{y}}={\mathbf{s}}+\widetilde{{\mathbf{z}}}\,\Sigma(m;\Delta)$ where the noise variance is perfectly tuned through the minimization problem (11) in order to faithfully “summarize” the complex interactions among variables in the original model; $\Sigma(m;\Delta)$ thus plays the role of a “mean-field”. See e.g. [5, 47] for more details on the mean-field theory and its applications.

This theorem has already been obtained recently in [22, 31] (with varying hypothesis on $P_{0}$) by the more elaborate methods mentionned in the introduction. In the next paragraphs we introduce the adaptive interpolation method through a pedagogical and new proof of this theorem.

Remark 1** (Free energy, mutual information and algorithms).**

In Bayesian inference the average free energy is related to the mutual information $I({\mathbf{S}};{\mathbf{W}})$ between the observation and the unknown vector (which is formally expressed as a difference of Shannon entropies: $I({\mathbf{S}};{\mathbf{W}})=H({\mathbf{W}})-H({\mathbf{W}}|{\mathbf{S}})$). For model (1), a straightforward computation shows that when $P_{0}$ has bounded first four moments

[TABLE]

where $S\sim P_{0}$. The $n\to\infty$ limit of the mutual information (or equivalently of the average free energy) is an interesting object to compute because it allows to locate the phase transition(s) occuring in the inference problem, which corresponds to its non-analyticity point(s) as a function of $\Delta$. This phase transition threshold usually separates a low-noise regime where inference is information theoretically possible from a high-noise regime where inference is impossible. In this high-noise regime the observation simply does not carry enough information for reconstructing the signal. Furthermore, remarkably, the replica formula for the mutual information (or average free energy) also allows to determine an algorithmic noise threshold, below the phase transition threshold, which separates the information theoretic possible phase in two regions: An “easy” phase where there exist low complexity message-passing algorithms for optimal inference and a “hard” phase where message-passing algorithms yield suboptimal inference. For further information and rigorous results on these issues for model (1) we refer to [22]. A few more pointers to the literature are given in the conclusion.

Remark 2** (Channel universality).**

The Gaussian noise setting (1) is actually sufficient to completely characterize the generic model where the entries of ${\mathbf{w}}$ are observed through a noisy element-wise (possibly non-linear) output probabilistic channel $P_{\rm out}(w_{ij}|s_{i}s_{j}/\sqrt{n})$. This is made possible by a theorem of channel universality [21] (conjectured in [48] and already proven for community detection in [49]). Roughly speaking this theorem states that given an output channel $P_{\rm out}(w|y)$, such that at $y=0$ the function $y\mapsto\ln P_{\rm out}(w|y)$ is three times differentiable, with bounded second and third derivatives, then the mutual information satisfies

[TABLE]

where $\Delta$ is the inverse Fisher information (at $y=0$) of the output channel:

[TABLE]

Informally, this means that we only have to compute the mutual information for a Gaussian channel to take care of a wide range of problems, which can be expressed in terms of their Fisher information.

2.2 The $(k,t)$–interpolating model

Let ${\mathbf{z}}^{(k)}=[z_{ij}^{(k)}]_{i,j=1}^{n},\widetilde{{\mathbf{z}}}^{(k)}=[\widetilde{z_{i}}^{(k)}]_{i=1}^{n}$, $Z_{ij}^{(k)}=Z_{ji}^{(k)}\sim{\cal N}(0,1)$, $\widetilde{Z}_{i}^{(k)}\sim{\cal N}(0,1)$ for $k=1,\ldots,K$ be Gaussian noise symmetric matrices and vectors. It is important to keep in mind that these are indexed both by the vertex indices $i,j$ and the discrete global interpolation parameter $k$.

The $(k,t)$–interpolating Hamiltonian is

[TABLE]

where the trial parameters $\{m_{k}\}_{k=1}^{K}$ are to be fixed later (these will be chosen ${\cal O}(1)$ with respect to (w.r.t) $n$ and can be interpreted as signal-to-noise ratios), $t\in[0,1]$ the continuous local interpolation parameter, and

[TABLE]

Here the subscript “mf” stands for “mean-field”.

A possible interpretation of the scheme is the following. The $(k,t)$–interpolating model corresponds to the following inference model. One has access to the following sets of noisy observations about the signal ${\mathbf{s}}$ where each noise realization is independent:

[TABLE]

The first and third sets of observations correspond to similar inference channels as the original model (1) but with a much higher noise variance proportional to $K$. These correspond to the first and third terms, respectively, of the $(k,t)$–interpolating Hamiltonian (2.2). The second and fourth sets instead correspond to decoupled Gaussian denoising models, with associated “mean-field” second and fourth terms in (2.2). The noise variances are proportional to $K$ because the total number of observations is $K$ and we want the total signal-to-noise ratio to be ${\cal O}(1)$. At fixed $k$, letting $t$ increase from [math] to $1$ increases by one unit the number of decoupled observations (18) by continuously adding the observation (20): Its signal-to-noise ratio that vanishes at $t=0$ (which is equivalent to not having access to this observation) becomes finite and equal to the signal-to-noise ratio of the individual observations in the set (18) at $t=1$. Simultaneously it reduces by one the number of observations of the form (17) by “removing” the observation (19): its signal-to-noise ratio, which is finite at $t=0$, vanishes at $t=1$. From (17)–(20) it is clear that the $(k,t=1)$ and $(k+1,t=0)$–interpolating models are statistically equivalent. A complementary and more graphical illustration of the interpolation scheme is found on Figure 1.

In order to use an important concentration lemma later on, we will need a slightly more general Hamiltonian, and consider the following perturbed version of (2.2):

[TABLE]

with i.i.d $\widehat{Z}_{i}\sim{\cal N}(0,1)$ and $\bm{\Theta}\vcentcolon=\{{\mathbf{S}},\{{\mathbf{Z}}^{(k)},\widetilde{{\mathbf{Z}}}^{(k)}\}_{k=1}^{K},\widehat{\mathbf{Z}}\}$ is the collection of all quenched random variables. It should be kept in mind that the signal-to-noise ratio $\epsilon$ of this additional Gaussian “side-channel” ${\mathbf{y}}={\mathbf{s}}\sqrt{\epsilon}+\widehat{\mathbf{z}}$ will tend to [math] at the end of the proof. Therefore we always consider $\epsilon\in[0,1]$.

The $(k,t)$–interpolating model has an associated posterior distribution, Gibbs expectation $\langle-\rangle_{k,t;\epsilon}$ and $(k,t)$–interpolating free energy $f_{k,t;\epsilon}$:

[TABLE]

In the following, we simply denote $\mathbb{E}_{\bm{\Theta}}$ by $\mathbb{E}$.

Lemma 1** (Linking the perturbed and plain free energies).**

Let $P_{0}$ have finite second moment. Then for the initial and final systems

[TABLE]

A short and generic proof is found in Appendix A. This statement shows in particular that if the thermodynamic limit $n\to+\infty$ exists, then it can be exchanged with the limit $\epsilon\to 0_{+}$ (as long as $P_{0}$ has bounded second moment). We stress that the existence of the thermodynamic limit is not directly used in our subsequent analysis, but rather, follows as a consequence.

2.3 The initial and final models

Let us compute the $(k,t)$–interpolating free energy $f_{1,0;0}$ associated with the initial $(k=1,t=0)$ model. Using (2.2) and (15),

[TABLE]

As the $Z_{ij}^{(k)}$, $1\leq i\leq j\leq n$, are i.i.d ${\cal N}(0,1)$ random variables, they possess the stability property, namely $Z_{ij}\vcentcolon=\sum_{k=1}^{K}Z_{ij}^{(k)}/\sqrt{K}$ are i.i.d ${\cal N}(0,1)$ random variables as well (and symmetric). Let ${\mathbf{z}}=[z_{ij}]_{i,j=1}^{n}$. Using this we obtain

[TABLE]

which is actually the free energy (7) of the original model. We thus have:

[TABLE]

Let us now consider the free energy $f_{K,1;0}$ of the final model. Using (2.2) and (16) we get

[TABLE]

Define

[TABLE]

Simple algebra leads to

[TABLE]

We now proceed as previously using again the stability property of the Gaussian noise variables. Since $\widetilde{Z}_{i}^{(k)}$ are i.i.d ${\cal N}(0,1)$, then $\widetilde{Z}_{i}\vcentcolon=\sum_{k=1}^{K}\widetilde{Z}_{i}^{(k)}\sqrt{m_{k}/(Km_{\rm mf}^{(K)})}\sim{\cal N}(0,1)$ and are i.i.d. Let $\widetilde{\mathbf{z}}=[\widetilde{z}_{i}]_{i=1}^{n}$. Using (24) we find that $f_{K,1;0}$ can also be expressed as

[TABLE]

Expression (32) is nothing else than the free energy (10) associated with the following scalar denoising model: $y=s+\widetilde{z}\,\Sigma(m_{\rm mf}^{(K)};\Delta)$, which leads to

[TABLE]

2.4 Free energy change along the adaptive interpolation path

By construction of (2.2) we have the following coherency property (see Figure 1): The $(k,t=1)$ and $(k+1,t=0)$ models are equivalent (the Hamiltonian (2.2) is invariant under this change) and thus $f_{k,1;\epsilon}=f_{k+1,0;\epsilon}$ for any $k$. This implies that the $(k,t)$–interpolating free energy (24) verifies

[TABLE]

Let us evaluate $df_{k,t;\epsilon}/dt$. Define the overlap $q_{{\mathbf{x}},{\mathbf{s}}}\vcentcolon=n^{-1}\sum_{i=1}^{n}x_{i}s_{i}$. Starting from (24), lenghty but simple algebra (see sec. 2.7 for the details) shows that as long as $P_{0}$ has bounded first four moments,

[TABLE]

This, with (34) and (33) yields

[TABLE]

where in the last equality we used (8) and introduced the non-negative variance

[TABLE]

The fundamental sum rule (36) can now be used to prove the replica symmetric formula.

2.5 Upper bound

From (36) we recover the upper bound usually obtained by the classical method of Guerra and Toninelli [50] and applied in [21] to symmetric rank-one matrix estimation (but see also [20] which already fully proved the replica formula in the binary case). Choose $m_{k}={\rm argmin}_{m\geq 0}f_{\rm RS}(m;\Delta)$ for all $k=1,\ldots,K$. This implies $m_{\rm mf}^{(K)}={\rm argmin}_{m\geq 0}f_{\rm RS}(m;\Delta)$ as well as $V_{K}(\{m_{k}\})=0$. Thus since the integrand in (36) is non-negative we get the bound

[TABLE]

Now we apply this inequality to a sequence $\epsilon_{n}\to 0$ as $n\to+\infty$. From Lemma 1 and (28) we obtain the upper bound:

Proposition 1** (Upper bound).**

Fix $\Delta>0$. For any $P_{0}$ with bounded first four moments,

[TABLE]

2.6 Lower bound

The converse bound is generally the one requiring extra technical tools, such as the use of spatial coupling [51, 28, 22, 17, 18] or the Aizenman-Sims-Starr scheme, see [35, 33, 32, 31]. Thanks to the adaptive interpolation method the proof is quite straightforward. As in all of the existing methods, we need a concentration lemma which takes the following form in the present context (see sec. 5 for the proof).

Lemma 2** (Overlap concentration).**

Let $P_{0}$ have bounded support. For any sequences $K_{n}\to+\infty$, $0<a_{n}<b_{n}<1$, and any choice of the trial parameters $m_{k}:\epsilon\mapsto m_{k}^{(n)}(\epsilon)$, $k=1,\cdots,K$ differentiable, bounded, non-decreasing with respect to $\epsilon$, we have

[TABLE]

for any $0<\alpha<1/4$ and some constant $C>0$ independent of $n$, $K_{n}$, $b_{n}$ and the set of trial parameters ($C$ depends on the second moment and the support of $P_{0}$).

Remark 3**.**

In applications of this lemma the sequence $a_{n}$ tends to zero as slowly as we wish. In practice we will set later on $b_{n}=2a_{n}$ and take $a_{n}\to 0$ slowly enough so that $a_{n}^{-3}n^{-\alpha}\to 0$. In particular the r.h.s of (40) tends to zero.

For sequences $K_{n}$, $0<a_{n}<b_{n}<1$, and $\{m_{k}=m_{k}^{(n)}(\epsilon)\}_{k=1}^{K}$ as in Lemma 2, (36) becomes

[TABLE]

where ${\cal O}(a_{n}^{-2}n^{-\alpha})$ is uniform in the choice of $K_{n}$, $b_{n}$ and trial parameters. At this point we need another important Lemma (see Appendix B for the proof) which is made possible by construction of the adaptive interpolation method.

Lemma 3** (Weak $t$-dependence at fixed $k$).**

Fix $K$, $\epsilon$ and $\{m_{k}\}_{k=1}^{K}$. For $P_{0}$ with bounded first four moments and any $k\in\{1,\ldots,K\}$ and $t\in[0,1]$,

[TABLE]

uniformly in $\epsilon$ and $\{m_{k}\}_{k=1}^{K}$. This result also applies when fixed $K$ is replaced by $K_{n}$ and $\{m_{k}=m_{k}^{(n)}(\epsilon)\}$.

Using this lemma for a sequence $K_{n}=\Omega(n^{b})$ with $b>1$ large enough, say $b=2$, (2.6) takes the following convenient form:

[TABLE]

We now use the last crucial lemma which is a fundamental property of the adaptive interpolation.

Lemma 4** (Choice for the trial parameters).**

For a given $n$ one can freely select differentiable and non-decreasing trial parameters $\{m_{k}=m_{k}^{(n)}(\epsilon)\}_{k=1}^{K}$ as

[TABLE]

Proof.

This is authorized by construction of the adaptive interpolation method. Indeed, the $(k=1,t=0)$–interpolating model (see the Hamiltonian ${\cal H}_{1,0;\epsilon}({\mathbf{x}};\bm{\theta})$ in (2.2)) is independent of $\{m_{k}\}_{k=1}^{K_{n}}$. Thus we can freely set $m_{1}=m_{1}^{(n)}(\epsilon)=\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{1,0;\epsilon}]$. Once $m_{1}$ is fixed to this value $m_{1}^{(n)}(\epsilon)$, we go to the next step and set $m_{2}=m_{2}^{(n)}(\epsilon)=\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{2,0;\epsilon}]$, which again is possible due to the fact that the Hamiltonian ${\cal H}_{2,0;\epsilon}({\mathbf{x}};\bm{\theta})$ and the Gibbs average $\langle-\rangle_{2,0;\epsilon}$ as well depend only on $m_{1}$ which has already been fixed. And so forth: As seen from Fig. 1, the Gibbs average $\langle-\rangle_{k,0;\epsilon}$ depends only on $\{m_{k^{\prime}}\}_{k^{\prime}=1}^{k-1}$ which were already fixed in the previous steps so that the choice (44) is valid. Note that $\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{k,0;\epsilon}]\geq 0$ which is important as the $m_{k}^{(n)}(\epsilon)$’s play the role of signal-to-noise ratios, and thus must be positive. Moreover the maps $\epsilon\mapsto m_{k}^{(n)}(\epsilon)$ are differentiable and non-decreasing. They are obviously differentiable since we work with $n$ finite. To see that they are non-decreasing we look at their derivative. It is easy to see from the construction of the Gibbs bracket that $\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{k,0;\epsilon}]$ is a function $\epsilon\mapsto G_{k}^{(n)}(\epsilon+\frac{1}{K}\sum_{l=1}^{k-1}m_{l})$ so that

[TABLE]

Now, $G_{k}^{(n)\prime}$ is an expected variance, and is therefore positive, as can be directly shown from a direct calculation (see equations (118) and (119) in section 5). This implies by induction that $\frac{d}{d\epsilon}m_{k}^{(n)}(\epsilon)\geq 0$. ∎

With this particular choice of trial parameters $\{m_{k}=m_{k}^{(n)}(\epsilon)\}_{k=1}^{K_{n}}$ the sum over $k=1,\ldots,K_{n}$ in (2.6) is set to zero: The interpolation path has been adapted (thus the name of the method). Since $V_{K_{n}}$ is non-negative, (2.6) directly implies the following lower bound:

[TABLE]

Finally, setting $b_{n}=2a_{n}$ and taking $a_{n}\to 0$ slowly enough as $n\to+\infty$ so that $a_{n}^{-3}n^{-\alpha}\to 0$, using Lemma 1 and (28) and the mean value theorem, we deduce

Proposition 2** (Lower bound).**

Fix $\Delta>0$. For any $P_{0}$ with bounded support,

[TABLE]

Remark 4** (The overlap must concentrate).**

Note that it is not obvious that one can find $\{m_{k}^{(n)}(\epsilon)\}_{k=1}^{K_{n}}$ which directly cancel the integrals in the fundamental identity (36) without using the overlap concentration of Lemma 2. Overlap concentration is a fundamental requirement of the above proof. This agrees with the statistical physics assumption that a necessary condition for the validity of the replica symmetric method is precisely the overlap concentration [5].

In Appendix C we present an alternative useful, albeit not completely rigorous, argument to obtain the lower bound.

2.7 Proof of the fundamental sum rule

In this paragraph we derive the formula (35). We will need a simple but fundamental identity555This identity has been abusively called “Nishimori identity” in the statistical physics literature. One should however note that it is a simple consequence of Bayes formula (see e.g appendix B of [18]). The “true” Nishimori identity [52] concerns models with one extra feature, namely a gauge symmetry which allows to eliminate the input signal, and the expectation over ${\mathbf{S}}$ in (47) can therefore be dropped (see e.g. [20]). which is a straightforward consequence of the Bayes law. Let ${\mathbf{X}}$, ${\mathbf{X}}^{\prime}$ be two i.i.d “replicas” drawn according to the product distribution $P_{k,t;\epsilon}({\mathbf{x}}|\bm{\theta})P_{k,t;\epsilon}({\mathbf{x}}^{\prime}|\bm{\theta})$. Recall the notation $\bm{\theta}\vcentcolon=\{{\mathbf{s}},\{{\mathbf{z}}^{(k)},\widetilde{{\mathbf{z}}}^{(k)}\}_{k=1}^{K},\widehat{\mathbf{z}}\}$ for the quenched variables and $\mathbb{E}=\mathbb{E}_{\bm{\Theta}}$ for the expectation with respect to these. Then for any function $g$ which does not depend on the Gaussian noise random variables,

[TABLE]

We give a proof of this identity in Appendix D for completeness.

Let us now compute $df_{k,t;\epsilon}/dt$. Starting from (2.2), (21), (24) one obtains

[TABLE]

Now we integrate by part the Gaussian noise using the elementary formula $\mathbb{E}_{Z}[Zf(Z)]=\mathbb{E}_{Z}[f^{\prime}(Z)]$ where $f^{\prime}$ is the derivative of $f$. This leads to

[TABLE]

where ${\mathbf{X}}$, ${\mathbf{X}}^{\prime}$ are the two i.i.d replicas drawn according to (22). An application of identity (47) then leads to

[TABLE]

The Cauchy-Schwarz inequality and (47) imply that $\mathbb{E}[\langle n^{-2}\sum_{i=1}^{n}X_{i}^{2}S_{i}^{2}\rangle_{k,t;\epsilon}]={\cal O}(n^{-1})$ as long as $P_{0}$ has bounded fourth moment. Indeed, by Cauchy-Schwarz

[TABLE]

and by (47) we have $\mathbb{E}[\langle X_{i}^{4}\rangle_{k,t;\epsilon}]=\mathbb{E}[S_{i}^{4}]$ for $i=1,\dots,n$, thus we get

[TABLE]

Finally, expressing the two other terms in (49) uisng the overlap $q_{{\mathbf{x}},{\mathbf{s}}}=n^{-1}\sum_{i=1}^{n}x_{i}s_{i}$ we find (35).

3  Application to rank-one symmetric tensor estimation

The present method can be extended to cover rank-one symmetric tensor estimation, which amounts to treat the $p$-spin model on the Nishimori line. For binary spins the Guerra-Toninelli bound was proven in [20] for any value of $p$, the replica symmetric formula was proved in the whole phase diagram for $p=2$, and also in a restricted region away from the first order phase transition for $p\geq 3$. A complete proof for $p=2$ and general spins (that can thus be real) was achieved using the spatial coupling technique in [22] and in [32] by a rigorous version of the cavity method. The case $p\geq 2$ and general spins has been treated using again the cavity method in [31].

3.1 Symmetric rank-one tensor estimation: Setting and main result

The symmetric tensor problem is very close to the matrix case presented in full details in sec. 2 so we only sketch the main steps. The observed symmetric tensor ${\mathbf{w}}\in\mathbb{R}^{n_{1}\times n_{2}\times\ldots\times n_{p}}$ is obtained through the following estimation model:

[TABLE]

where ${\mathbf{s}}\in\mathbb{R}^{n}$ with i.i.d components distributed according to a known prior $P_{0}$, ${\mathbf{Z}}\in\mathbb{R}^{n_{1}\times n_{2}\times\ldots\times n_{p}}$ is a symmetric Gaussian noise tensor with i.i.d (up to the symmetry constraint) ${\cal N}(0,1)$ entries. We note that, like in the case of symmetric matrix estimation of sec. 2.1, the channel universality property (see remark 2) is valid in the present setting. This means that by covering the case of additive white Gaussian noise (52), we actually treat a wide range of (component-wise) inference channels

[TABLE]

We refer to [48, 53, 32] for more details on this point. The free energy of the model is

[TABLE]

where the Hamiltonian $\mathcal{H}({\mathbf{x}};{\mathbf{s}},{\mathbf{z}})$ is

[TABLE]

For a $P_{0}$ with bounded first four moments the free energy is related to the mutual information $I({\mathbf{S}};{\mathbf{W}})$ through

[TABLE]

We define the replica symmetric potential for symmetric tensor estimation as

[TABLE]

where $\Sigma(m;\Delta)^{2}\vcentcolon=\Delta/m^{p-1}$ and $f_{{\rm den}}(\Sigma)$ is given by (10). Next we prove the RS formula.

Theorem 2** (RS formula for symmetric rank-one tensor estimation).**

Fix $\Delta>0$. For any $P_{0}$ with bounded support, the asymptotic free energy of the symmetric tensor estimation model (52) verifies

[TABLE]

Again, we note that the bounded support property of $P_{0}$ is only needed for concentration proofs and does not impose any upper limit on the size of the support. We believe this can be removed by a limiting process as long as $P_{0}$ has bounded first four moments.

3.2 Sketch of proof of the replica symmetric formula

We prove Theorem 2. Since this proof is similar to the one of Theorem 1 for the matrix case, we only give the main ideas. The starting point is the introduction of a (perturbed) $(k,t)$–interpolating Hamiltonian:

[TABLE]

where the trial parameters $\{m_{k}\}_{k=1}^{K}$ are to be fixed later and

[TABLE]

The associated $(k,t)$–interpolating model, Gibbs expectation and $(k,t)$–interpolating free energy are defined respectively by (22), (23) and (24). Using the stability property of the Gaussian noise variables, one can check that the intial and final $(k,t)$–interpolating models are such that

[TABLE]

where

[TABLE]

By a trivial generalization of the calculations of sec. 2.7, we obtain the variation of the $(k,t)$–interpolating free energy:

[TABLE]

where the overlap is again $q_{{\mathbf{x}},{\mathbf{s}}}\vcentcolon=n^{-1}\sum_{i}^{n}x_{i}s_{i}$. This result holds as long as $P_{0}$ has finite first four moments. Proceeding similarly to section 2.4 we get the sum rule

[TABLE]

where

[TABLE]

Note that $V_{K,p}$ is non-negative by Jensen’s inequality applied to the convex function $x\in\mathbb{R}_{+}\mapsto x^{\frac{p}{p-1}}\in\mathbb{R}_{+}$ (here the trial parameters are all non-negative). For $p=2$ the identity (65) reduces to (36).

The proof of Theorem 2 proceeds from this fundamental sum rule in much the same way as in the case $p=2$ of sections 2.5 and 2.6. Here we only give a brief summary of the arguments insisting only on the essential differences. We start with the upper bound.

For the case of even $p$ it is straightforward to derive an upper bound. One chooses all trial parameters as $m_{k}=m_{*}={\rm argmin}_{m\geq 0}f_{\rm RS}(m;\Delta)$, $k=1,\cdots,K$, which yields $V_{K,p}(\{m_{k}\}_{k=1}^{K})=0$. From there one can use the classic argument of Guerra-Toninelli: By convexity of $x\in\mathbb{R}\mapsto x^{p}\in\mathbb{R}_{+}$ for even $p$ we see that $q_{{\mathbf{X}},{\mathbf{S}}}^{p}-p\,m_{*}^{p-1}q_{{\mathbf{X}},{\mathbf{S}}}+(p-1)m_{*}^{p}\geq 0$, which implies the upper bound analogous to (38). Since Lemma 1 holds verbatim here, by taking a sequence $\epsilon_{n}\to 0_{+}$, $n\to+\infty$, we deduce as before the upper bound $\limsup_{n\to+\infty}f_{n}\leq\min_{m\geq 0}f_{\rm RS}(m;\Delta)$ (when $P_{0}$ has finite first four moments).

For the case of odd $p$ we cannot immediately apply the convexity argument. We first need to apply a concentration result. As it will become clear from its proof, Lemma 2 is generic and the same statement applies to the present tensor setting.666Here we use Lemma 2 but a weaker form of concentration is enough for this argument, namely it suffices to control the following type of “thermal” fluctuation $\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}^{2}\rangle_{k,t,\epsilon}-\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{k,t,\epsilon}^{2}]$. Moreover it is not necessary to allow for an $\epsilon$-dependence in $m_{k}$’s. Therefore for any sequence $K_{n}\to+\infty$, $0<a_{n}<b_{n}<1$, $b_{n}\to 0$, and trial parameters $\{m_{k}=m_{k}^{(n)}(\epsilon)\}_{k=1}^{K}$ which are non-decreasing functions of $\epsilon$ we have

[TABLE]

for $0<\alpha<1/4$ and ${\cal O}(a_{n}^{-2}n^{-\alpha})$ uniform in $K_{n}$, $b_{n}$, and trial parameters. Furthermore by the Nishimori identity (47) we see that $\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{k,t;\epsilon}]\geq 0$, so convexity of $x\in\mathbb{R}_{+}\mapsto x^{p}\in\mathbb{R}_{+}$ shows that the term under the integral is positive, which allows to deduce the upper bound as above (of course this argument works for any $p$ even or odd).

let us finally briefly discuss the lower bound. Lemma 3 and its proof hold for the tensor setting as well which means that in (67) we can replace $\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{k,t;\epsilon}]$ by $\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{k,0;\epsilon}]$. This then allows to choose (adapt) the sequence of trial parameters as $\{m_{k}=m_{k}^{(n)}(\epsilon)\}$ where $m_{k}=\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{k,0;\epsilon}]$, $k=1,\cdots,K_{n}$ just as in Lemma 4. Thus we obtain

[TABLE]

where we used the non-negativity of $V_{K,p}$ to get the inequality and the change of variable $\Sigma^{2}=\Delta/m^{p-1}$ to get the last line. The usual limiting argument, taking $b_{n}=2a_{n}$ with $a_{n}\to 0$ slowly enough so that $a_{n}^{-3}n^{-\alpha}\to 0$ when $n\to+\infty$, implies $\liminf_{n\to+\infty}f_{n}\geq\min_{m\geq 0}f_{\rm RS}(m;\Delta)$. We have proven this lower bound under the assumption of boundedness of the support of $P_{0}$ (used in the proof of overlap concentration)

Combining the upper and lower bounds yields Theorem 2. As a final note we remark that the alternative route to the lower bound proposed in Appendix C for the matrix case also holds essentially unchanged here.

4  Application to Gaussian random linear estimation

4.1 Gaussian random linear estimation: Setting and result

In Gaussian random linear estimation (RLE) one is interested in reconstructing a signal ${\mathbf{s}}=[s_{i}]_{i=1}^{n}\in\mathbb{R}^{n}$ from few noisy measurements ${\mathbf{y}}=[y_{\mu}]_{\mu=1}^{m}\in\mathbb{R}^{m}$ obtained from the projection of ${\mathbf{s}}$ by a random Gaussian measurement matrix $\bm{\phi}=[\phi_{\mu i}]_{\mu,i=1}^{m,n}\in\mathbb{R}^{m\times n}$ with i.i.d entries $\phi_{\mu i}\sim\mathcal{N}(0,1/n)$. The measurement rate is $\alpha\vcentcolon=m/n$. We consider i.i.d additive white Gaussian noise of known variance $\Delta$. Let the standardized noise components be $Z_{\mu}\sim\mathcal{N}(0,1)$, $\mu=1,\ldots,m$. Then the measurement model is

[TABLE]

The signal has i.i.d components distributed according to a discrete prior $P_{0}(s_{i})=\sum_{b=1}^{B}p_{b}\delta(s_{i}-a_{b})$ with a finite number $B$ of terms and $\max_{b}|a_{b}|\leq s_{\rm max}$. Note that the more general case where the signal has i.i.d vectorial components, as considered in [17, 18], can be tackled with our proof technique exactly in the same way but we consider the scalar case for the sake of notational simplicity.

The free energy of the RLE model (69) (which is also equal to the mutual information per component $I({\mathbf{S}};{\mathbf{Y}})/n$ between the noisy observation and the signal) is defined as

[TABLE]

where $[\bm{\phi}({\mathbf{x}}-{\mathbf{s}})]_{\mu}\vcentcolon=\sum_{i=1}^{n}\bm{\phi}_{\mu i}(x_{i}-s_{i})$. Let

[TABLE]

Define the following RS potential:

[TABLE]

where $i_{\rm den}(\Sigma)=I(S;S+\widetilde{Z}\,\Sigma)$ is the mutual information of a scalar Gaussian denoising model $y=s+\widetilde{z}\,\Sigma$ with $S\sim P_{0}$, $\widetilde{Z}\sim{\cal N}(0,1)$, and $\Sigma^{-2}$ an effective signal to noise ratio:

[TABLE]

We will prove the RS formula (already proven in [17, 18, 29, 30]):

Theorem 3** (RS formula for Gaussian RLE).**

Fix $\Delta>0$. For any discrete $P_{0}$, the asymptotic free energy of the RLE model (69) verifies

[TABLE]

Proof.

The result follows from Propositions 3 and 4 in sec. 4.2. ∎

4.2 Proof of the RS formula

Let ${{\mathbf{z}}}^{(k)}=[z_{\mu}^{(k)}]_{\mu=1}^{m},\widetilde{{\mathbf{z}}}^{(k)}=[\widetilde{z_{i}}^{(k)}]_{i=1}^{n}$ and $\widehat{{\mathbf{z}}}=[\widehat{z_{i}}]_{i=1}^{n}$ all with i.i.d ${\cal N}(0,1)$ entries for $k=1,\ldots,K$. Define $\Sigma_{k}\vcentcolon=\Sigma(E_{k};\Delta)$ where the trial parameters $\{E_{k}\}_{k=1}^{K}$ are fixed later on. The (perturbed) $(k,t)$–interpolating Hamiltonian for the present problem is

[TABLE]

Again, the last term is a small perturbation needed to use an important concentration result (here equation (93)). Here $\bm{\theta}\vcentcolon=\{{\mathbf{s}},\{{\mathbf{z}}^{(k)},\widetilde{{\mathbf{z}}}^{(k)}\}_{k=1}^{K},\widehat{\mathbf{z}},\bm{\phi}\}$, $k\in\{1,\ldots,K\}$, $t\in[0,1]$ and

[TABLE]

where $\bar{\mathbf{x}}={\mathbf{x}}-{\mathbf{s}}$, $\bar{x}_{i}=x_{i}-s_{i}$. Moreover the “signal-to-noise functions” $\{\gamma_{k}(t),\lambda_{k}(t)\}_{k=1}^{K}$ verify

[TABLE]

as well as the following constraint (see [18] for an interpretation of this formula)

[TABLE]

We also require $\gamma_{k}(t)$ to be strictly decreasing with $t$. The associated $(k,t)$–interpolating model, Gibbs expectation and $(k,t)$–interpolating free energy are defined respectively by (22), (23) and (24) with the Hamiltonian (4.2). Note that Lemma 1 remains valid for the present model (with the same proof).

Similarly as in sec. 2.3, and using again the stability property of the Gaussian random noise variables, it is easy to verify that the initial and final $(k,t)$–interpolating models correspond to the RLE and denoising models respectively, that is

[TABLE]

where

[TABLE]

As before we use the identity (34) and compute the free energy change along the adaptive interpolation. Straightforward differentiation leads to (with $\bar{\mathbf{X}}={\mathbf{X}}-{\mathbf{S}}$)

[TABLE]

where as before $\mathbb{E}$ denotes the average w.r.t to all quenched random variables $\bm{\theta}$ and $\langle-\rangle_{k,t;\epsilon}$ the Gibbs average with Hamiltonian (4.2). The two quantities (86) and (87) can be simplified using Gaussian integration by parts. For example, integrating by parts w.r.t $Z_{\mu}^{(k)}$,

[TABLE]

It allows to simplify $\mathcal{A}_{k,t;\epsilon}$ as follows,

[TABLE]

where we recognized the “measurement minimum mean-square-error”

[TABLE]

For $\mathcal{B}_{k,t;\epsilon}$ we proceed similarly with an integration by parts w.r.t $\widetilde{Z}_{i}^{(k)}$, and find

[TABLE]

using (81) for the last equality, and the minimum mean-square-error (MMSE) defined as

[TABLE]

The free energy can be shown to concentrate by generalizing the computations of Appendix E in [18] taking into account that the noise variables $\{Z_{\mu}^{(k)},\widetilde{Z}_{i}^{(k)}\}$ are indexed by the discrete interpolation parameter (the techniques of [18] use a discrete $P_{0}$ with bounded support for the free energy concentration). Since the free energy at fixed quenched random variables realization concentrates, both sec. VIII of [18] or sec. 5 of the present paper apply here (these are perfectly equivalent analyses and only require the identity (47) and the free energy concentration to be valid). Thus the overlap $q_{{\mathbf{x}},{\mathbf{s}}}\vcentcolon=n^{-1}\sum_{i}x_{i}s_{i}$ concentrates too. As a consequence an analog of Lemma 4.6 in [18] can be shown here: Fix a discrete $P_{0}$ with bounded support. For any sequence $K_{n}\to+\infty$, and $0<a_{n}<b_{n}<1$ (that tend to zero slowly enough in the application), and trial parameters $\{E_{k}=E_{k}^{(n)}(\epsilon)\}_{k=1}^{K_{n}}$ which are differentiable, bounded and non-increasing in $\epsilon$, we have

[TABLE]

for some $0<\alpha<1$ and $C>0$.

Now combining (34), (82), (83), (85) and (89), (91), together with (93), we obtain

[TABLE]

We need the following useful identity which can easily be checked using (72), (79), (80), (81):

[TABLE]

Let us define

[TABLE]

With the help of (95) and (96) the identity (94) becomes

[TABLE]

This is the fundamental sum rule which forms the basis for the proof of Theorem 3.

We start with the upper bound. As in sec. 2.5 we choose $E_{k}=E_{*}\vcentcolon={\rm argmin}_{E\geq 0}f_{\rm RS}(E;\Delta)$ for all $k=1,\ldots,K_{n}$ (here $\epsilon$-independent) which implies that $\Sigma_{\rm mf}(\{E_{k}=E_{*}\}_{k=1}^{K_{n}};\Delta)=\Sigma(E_{*};\Delta)$ and thus, as seen from (96), $\widetilde{f}_{\rm RS}(\{E_{k}=E_{*}\}_{k=1}^{K_{n}};\Delta)=\min_{E\geq 0}f_{\rm RS}(E;\Delta)$. Thus since the integrand in (97) is non-positive (recall that $d\gamma_{k}(t)/dt\leq 0$) and using the arguments similar to sec. 2 in order to take the $n\to+\infty$ limit, we get:

Proposition 3** (Upper bound).**

Fix $\Delta>0$. For $P_{0}$ discrete and with bounded support:

[TABLE]

Let us now prove the lower bound. This bound required the use of spatial coupling in [17, 18] or “conditional central limit theorems” in [29, 30]. Here we derive the bound in a direct and much simpler manner following the same steps as in sec. 2.6. We first need the following identity: For any discrete $P_{0}$ with bounded support, any $k\in\{1,\ldots,K\}$ and $\epsilon,t\in[0,1]$,

[TABLE]

Its proof is very similar to the one of Lemma 3. Using this identity with $K_{n}=\Omega(n^{b})$, $b>2$, in (97) and constructing $E_{k}^{(n)}(\epsilon)={\rm mmse}_{k,0;\epsilon}$ (which is indeed non-increasing with $\epsilon$ being a MMSE) for all $k=1,\ldots,K_{n}$ (by the same arguments than those in the proof of Lemma 4), we reach

[TABLE]

Recall $\Sigma^{-2}\vcentcolon=\alpha/(E+\Delta)$ and thus $E=\alpha/\Sigma^{-2}-\Delta$. For given $\Delta$ we set $\widetilde{\psi}_{\Delta}(\Sigma^{-2})\vcentcolon=\psi(\alpha/\Sigma^{-2}-\Delta;\Delta)$ and note that $\widetilde{\psi}_{\Delta}(\cdot)$ is a convex function. Thus from (96)

[TABLE]

Thus (100) becomes

[TABLE]

Taking $b_{n}=2a_{n}$, $a_{n}\to 0$, such that $a_{n}^{-3}n^{-\alpha}\to 0$ as $n\to+\infty$ we obtain (recall Lemma 1)

Proposition 4** (Lower bound).**

Fix $\Delta>0$. For any dicrete $P_{0}$ with bounded support,

[TABLE]

5  Concentration of overlaps

The main goal of this section is the proof of Lemma 2. The proof strategy outlined here is very general and it will appear to the reader that it applies to essentially any inference problem for which the identity (47) is valid and as long as the free energy can be shown to concentrate. In the framework of inference problems such proofs go back to [12, 16, 20] for binary signals (in coding, CDMA and the gauge symmetric p-spin model) and have been extended more recently in random linear estimation for arbitrary signal distributions [18]. The results and exposition given here slightly generalize and streamlines the one of the previous works.

From now on the trial parameters are chosen of the form $\{m_{k}=m_{k}^{(n)}(\epsilon)\}_{k=1}^{K}$. It will be convenient to adopt the notation $\tilde{\epsilon}=\epsilon+(K\Delta)^{-1}(\sum_{l=1}^{k-1}m_{l}^{(n)}(\epsilon)+tm_{k}^{(n)}(\epsilon))$. Here $\tilde{\epsilon}$ depends on $k,t$ but we do not write this dependence explicitly as it does not play a role (we work at fixed $k,t$ in the rest of this section). Let

[TABLE]

We will show that Lemma 2 is a direct consequence of the following:

Proposition 5** (Concentration of $\mathcal{L}$ on $\mathbb{E}\langle\mathcal{L}\rangle_{k,t;\epsilon}$ ).**

Let $P_{0}$ with finite second moment and bounded support in $[-M,M]$. For any choice of trial parameters $\{m_{k}=m_{k}^{(n)}(\epsilon)\}_{k=1}^{K}$ that are non-decreasing bounded and differentiable functions of $\epsilon\in]0,1[$, and any sequences $0<a_{n}<b_{n}<1$, we have

[TABLE]

for any $0<\alpha<1/4$ with $C>0$ a constant uniform in $k,t$ and the trial parameters and depending only on the second moment of $P_{0}$ and $M$.

The proof of this proposition is broken in two parts. Notice that

[TABLE]

Thus it suffices to prove the two following lemmas. The first lemma expresses concentration w.r.t the posterior distribution (or “thermal fluctuations”) and is an elementary consequence of concavity properties of the free energy.

Lemma 5** (Concentration of $\mathcal{L}$ on $\langle\mathcal{L}\rangle_{k,t;\epsilon}$ ).**

Let $P_{0}$ with finite second moment. For any choice of trial parameters $\{m_{k}=m_{k}^{(n)}(\epsilon)\}_{k=1}^{K}$ that are non-decreasing bounded and differentiable functions of $\epsilon\in]0,1[$, and any sequences $0<a_{n}<b_{n}<1$,

[TABLE]

The second lemma expresses the concentration of the Gibbs average w.r.t the realizations of quenched disorder variables.

Lemma 6** (Concentration of $\langle\mathcal{L}\rangle_{k,t;\epsilon}$ on $\mathbb{E}\langle\mathcal{L}\rangle_{k,t;\epsilon}$ ).**

Let $P_{0}$ with finite second moment and bounded support in $[-M,M]$. For any choice of trial parameters $\{m_{k}=m_{k}^{(n)}(\epsilon)\}_{k=1}^{K}$ that are non-decreasing bounded and differentiable functions of $\epsilon\in]0,1[$, and any sequences $0<a_{n}<b_{n}<1$,

[TABLE]

for any $0<\eta<1/2$ and where $C>0$ depends only on the second moment of $P_{0}$ and $M$. In particular $C$ is independent of $k,t$ and the trial parameters.

Remark 5**.**

Thanks to the identity (111) below, that we will show in section 6, the statements of Proposition 5 and Lemmas 5 and 6 hold if we replace $\mathcal{L}$ by the overlap $q_{{\mathbf{x}},{\mathbf{s}}}$.

The proof of this last lemma is based on an important but generic result concerning the concentration of the $(k,t)$–interpolating free energy for a single realization of quenched variables. Let

[TABLE]

Recall that $f_{k,t;\epsilon}=\mathbb{E}[F_{k,t;\epsilon}(\bm{\Theta})]$.

Proposition 6** (Concentration of the $(k,t)$–interpolating free energy).**

Let $P_{0}$ with bounded support in $[-M,M]$. One can find $c>0$ which depends only on $M$ and $\Delta$ such that for all $k=1,\ldots,K$, $t\in[0,1]$ and $\epsilon\in[0,1]$,

[TABLE]

where $u>0$. Explicit expressions for $c$ can be derived from (176) in sec. 7.

This proposition is proved in sec. 7. In the rest of this section we prove Lemmas 2, 5 and 6. The parameters $k$ and $t$ stay fixed and do not play any role, but it is important to be careful about the $\epsilon$ dependence.

Proof of Lemma 2

The proof is based on the remarkable identity (here $S\sim P_{0}$)

[TABLE]

Its derivation is found in sec. 6 and involves lengthy algebra using identity (47) and integrations by parts w.r.t the Gaussian noise. This formula implies

[TABLE]

and using Fubini’s theorem

[TABLE]

Then applying Proposition 5 we obtain (since the bounds are uniform in $k,t$)

[TABLE]

so that (40) is verified for any $0<\alpha<1/4$. $\blacksquare$

We now turn to the proof of Lemmas 5 and 6. The main ingredient is a set of formulas for the first two derivatives of the free energy w.r.t $\tilde{\epsilon}$. For any given realisation of the quenched disorder we have the equalities (here $Z_{i}\sim\mathcal{N}(0,1)$ i.i.d)

[TABLE]

Averaging (115) and (116) and using a Gaussian integration by parts w.r.t $z_{i}$ and the identity $\mathbb{E}[\langle X_{i}\rangle_{k,t;\epsilon}S_{i}]=\mathbb{E}[\langle X_{i}\rangle_{k,t;\epsilon}^{2}]$ (again a special case of (47)), we find (see Appendix A)

[TABLE]

There is another useful formula for $d^{2}f_{k,t;\epsilon}/d{\tilde{\epsilon}}^{2}$ that can be worked out directly (see sec. 6) by differentiating the second expression in (117) instead of the first:

[TABLE]

This formula clearly shows that $f_{k,t;\epsilon}$ is a concave function of $\tilde{\epsilon}$.

Proof of Lemma 5

From (118) we have

[TABLE]

where we used $\tilde{\epsilon}\geq\epsilon$ and $\mathbb{E}[\langle X_{i}^{2}\rangle_{k,t;\epsilon}]=\mathbb{E}[S^{2}]$ (an application of (47)). We perform an integration of this inequality over $\epsilon\in[a_{n},b_{n}]$. Note that the map $\epsilon\in[a_{n},b_{n}]\mapsto\tilde{\epsilon}\in[\tilde{\epsilon}(a_{n}),\tilde{\epsilon}(b_{n})]$ is differentiable and the inverse map is well defined and also differentiable since we have assumed that the trial parameters are differentiable and non decreasing. Obviously the Jacobian $J=d\tilde{\epsilon}/d\epsilon\geq 1$ since the trial parameters are non-decreasing. Integrating over $\epsilon\in[a_{n},b_{n}]$ and performing the change of variables $\epsilon\mapsto\tilde{\epsilon}$, and using $J\geq 1$, we obtain

[TABLE]

From (117) combined with the convexity of the square and an application of the Nishimori identity, we see that the first term is certainly smaller in absolute value than $\frac{1}{n}\mathbb{E}[S^{2}]$. The second term is smaller than $\frac{\mathbb{E}[S^{2}]}{4n}|\ln a_{n}|$. This concludes the proof of Lemma 5. $\blacksquare$

Proof of Lemma 6

In what follows we view $F_{k,t;\epsilon}(\bm{\theta})$ and $f_{k,t;\epsilon}$ as functions of $\tilde{\epsilon}$. Recall that $P_{0}$ has bounded support in $[-M,M]$. Define the two functions of $\tilde{\epsilon}$

[TABLE]

Because of (116) we see that the second derivative of $\widetilde{F}(\tilde{\epsilon})$ is negative, so this is a concave function of $\tilde{\epsilon}$ (without this extra term $F_{k,t;\epsilon}$ is not necessarily concave, although $f_{k,t;\epsilon}$ is concave). Note also that $\widetilde{f}(\tilde{\epsilon})$ is concave. Concavity implies for any $\delta>0$

[TABLE]

The difference between the derivatives appearing on the r.h.s of these inequalities cannot be considered small because at a first order transition point the derivatives have jump discontinuities. Set

[TABLE]

where the signs of these quantities follow from concavity of $\widetilde{f}(\tilde{\epsilon})$. From (123), (124) and (125) we get

[TABLE]

Now we will cast this inequality in a more usable form. From (122)

[TABLE]

with

[TABLE]

and from (115), (117),

[TABLE]

From (127), (129) it is easy to show that (5) implies

[TABLE]

At this point we use Proposition 6. A standard argument given at the end of this proof shows that this proposition implies

[TABLE]

for any $0<\eta<1$. Squaring, then taking the expectation of (5) and using $\mathbb{E}[A^{2}]=\mathcal{O}(n^{-1})$, $\tilde{\epsilon}\geq\epsilon$, and $(\sum_{i=1}^{p}v_{i})^{2}\leq p\sum_{i=1}^{p}v_{i}^{2}$,

[TABLE]

We now take $\epsilon\in[a_{n},b_{n}]$ and $0<\delta<a_{n}$. Using the change of variables $\epsilon\mapsto\tilde{\epsilon}(\epsilon)$, that the Jacobian $J=d\tilde{\epsilon}/d\epsilon\geq 1$, $|d\widetilde{f}(\tilde{\epsilon})/d\tilde{\epsilon}|\leq(\mathbb{E}[S^{2}]+M/\sqrt{\tilde{\epsilon}})/2$ from (117) and (122), $C^{\pm}(\tilde{\epsilon})\geq 0$ from (125), and the mean value theorem

[TABLE]

Thus, integrating (5) over $\epsilon\in[a_{n},b_{n}]$ yields with $0<\delta<a_{n}$

[TABLE]

where $B\geq\tilde{\epsilon}$, because $\tilde{\epsilon}$ is bounded by assumption of the boundedness of the $m_{k}$’s and $\epsilon\leq 1$. Finally we choose $\delta=a_{n}n^{-\frac{1}{4}+\frac{\eta}{2}}$, $0<\eta<1/2$, and obtain for $n$ large enough (and $a$ fixed positive small)

[TABLE]

for some constant $C>0$ depending only on $M$ and $\mathbb{E}[S^{2}]$.

It remains to justify (131). By the Cauchy-Schwarz inequality and Proposition 6 we have

[TABLE]

If we can show that the moments of the (random) free energy $F_{k,t;\epsilon}(\bm{\Theta})$ are bounded uniformly in $n$, then the choice $u=n^{-1/2+\eta}$ for any $0<\eta<1/2$ allows to conclude the proof. Let us briefly show how the moments are estimated. By the Jensen’s inequality

[TABLE]

The expectation over ${\mathbf{X}}$ is computed from (2.2) and one finds a polynomial in $\{s_{i},\{z_{ij}^{(k)},\widetilde{z}_{i}^{(k)}\}_{k=1}^{K},\widehat{z}_{i}\}_{i=1}^{n}$ which all have bounded moments. On the other hand from (15), (16) by completing the squares we have

[TABLE]

and find that $\mathcal{H}_{k,t;\epsilon}({\mathbf{x}};\bm{\theta})$ is lower bounded by a polynomial in $\{s_{i},\{z_{ij}^{(k)},\widetilde{z}_{i}^{(k)}\}_{k=1}^{K},\widehat{z}_{i}\}_{i=1}^{n}$. This is also the case for $F_{k,t;\epsilon}(\bm{\theta})$. With these upper and lower bounds on $F_{k,t;\epsilon}(\bm{\theta})$ it is easy to show that for any integer $p$

[TABLE]

where $C_{p}$ is independent of $n$ and depends only on $\Delta$ and moments of $P_{0}$. $\blacksquare$

6  A fluctuation identity

The purpose of this appendix is to prove the identity (111) relating the various fluctuations. This identity is quite powerful and holds in quite some generality and in particular for the three applications presented in this paper. To alleviate the notation we denote $\langle-\rangle_{k,t;\epsilon}$ simply by $\langle-\rangle$. It actually follows from the exact formula

[TABLE]

that we derive next. But before doing so, let us show how (140) implies (111). First we note that by (47) the last subdominant sum equals $\mathbb{E}[S^{2}]/4n\tilde{\epsilon}={\cal O}(1/n)$. We then express the first two terms in terms of the overlap $q_{{\mathbf{x}},{\mathbf{s}}}$. From (47) we have $\mathbb{E}[\langle X_{i}X_{j}\rangle^{2}]=\mathbb{E}[S_{i}S_{j}\langle X_{i}X_{j}\rangle]$ and therefore

[TABLE]

Similarly $\mathbb{E}[\langle X_{i}\rangle^{2}]=\mathbb{E}[S_{i}\langle X_{i}\rangle]$, so

[TABLE]

and $\mathbb{E}[\langle X_{i}X_{j}\rangle\langle X_{i}\rangle\langle X_{j}\rangle]=\mathbb{E}[S_{i}S_{j}\langle X_{i}\rangle\langle X_{j}\rangle]$ which implies

[TABLE]

Replacing the three last identities in (140) leads to (111).

We now summarise the main steps leading to the formula (140), using the identity (47) and integrations by parts w.r.t the Gaussian noise. This formula follows by summing the two following identities

[TABLE]

We first derive the second identity which requires somewhat longer calculations.

Derivation of (145)

First we compute $\mathbb{E}[\langle\mathcal{L}\rangle]^{2}$. From (104) we have

[TABLE]

From (47) we have $\mathbb{E}[\langle X_{i}\rangle S_{i}]=\mathbb{E}[\langle X_{i}\rangle^{2}]$ and by an integration by parts

[TABLE]

Thus we find

[TABLE]

which is formula (117). Squaring, we have

[TABLE]

Now we compute $\mathbb{E}[\langle\mathcal{L}\rangle^{2}]$. From (104) we have

[TABLE]

Taking the expectation and using (47) (for the terms that do not contain explicit $z$-factors) we find

[TABLE]

In order to simplify this expression we now integrate by parts all terms that contain explicit $Z$-factors:

[TABLE]

Then

[TABLE]

and finally

[TABLE]

where in the last two identities we used (47) after the integration by parts. Replacing the last five identities (152)–(156) into (151) we get

[TABLE]

Subtracting (157) and (149) we finally find (145). $\blacksquare$

Derivation of (144)

Acting with $n^{-1}d/d{\tilde{\epsilon}}$ on both sides of (148) we find

[TABLE]

Computing the derivative of $\mathcal{L}$ and using (147) we find that (158) is equivalent to

[TABLE]

Now we compute the terms in the first sum. We have

[TABLE]

Then from (47),

[TABLE]

It remains to integrate by parts the two terms involving the explicit $\widehat{Z}_{j}$ dependence (these can be found in the previous integrations by parts). This leads to

[TABLE]

The formula (144) then follows from (159) and (6).

7  Concentration of the free energy

In this section we prove Proposition 6. We will call $\mathbb{E}_{{\mathbf{Z}}}$, $\mathbb{P}_{{\mathbf{Z}}}$ the expectation and probability law over all Gaussian variables, $\mathbb{E}_{{\mathbf{S}}}$, $\mathbb{P}_{{\mathbf{S}}}$ the ones over the input signal variables, and $\mathbb{E}$, $\mathbb{P}$ the ones over the joint law. The proof is broken up in two lemmas. We first show a lemma which expresses concentration w.r.t all Gaussian sources of disorder uniformly in the input signal.

Lemma 7** (Concentration w.r.t the Gaussian quenched disorder).**

Take $P_{0}$ with bounded support in $[-M,M]$. For any signal realisation ${\mathbf{s}}$ and all $k=1,\ldots,K$, $t\in[0,1]$ and $\epsilon>0$ we have

[TABLE]

where $u>0$.

Proof.

The proof method is again based on an interpolation (of a different kind) that goes back to a beautiful work of Guerra and Toninelli [54]. We fix the input signal realisation ${\mathbf{s}}$ and consider two i.i.d copies for the Gaussian quenched variables ${\mathbf{z}}^{(k,1)}=[z_{ij}^{(k,1)}]_{i,j=1}^{n},\widetilde{{\mathbf{z}}}^{(k,1)}=[\widetilde{z_{i}}^{(k,1)}]_{i=1}^{n}$ and ${\mathbf{z}}^{(k,2)}=[z_{ij}^{(k,2)}]_{i,j=1}^{n},\widetilde{{\mathbf{z}}}^{(k,2)}=[\widetilde{z_{i}}^{(k,2)}]_{i=1}^{n}$. We also need two copies of the extra Gaussian noise introduced in the perturbation term (21), namely $\widehat{{\mathbf{z}}}^{(1)}=[\widehat{z}_{i}^{(1)}]_{i=1}^{n}$ and $\widehat{{\mathbf{z}}}^{(2)}=[\widehat{z}_{i}^{(2)}]_{i=1}^{n}$. We define an Hamiltonian interpolating between the two realizations of the Gaussian disorder, with new interpolating parameter $\tau\in[0,1]$:

[TABLE]

Let ${\cal Z}_{k,t;\epsilon}(\tau)\vcentcolon=\int\{\prod_{i=1}^{n}dx_{i}P_{0}(x_{i})\}\exp(-{\cal H}_{k,t,\tau;\epsilon})$ the partition function associated to ${\cal H}_{k,t,\tau;\epsilon}$. Let $s>0$ be a trial parameter to be fixed later on and let

[TABLE]

where $\mathbb{E}_{1}$ and $\mathbb{E}_{2}$ are the expectations w.r.t the two independent sets of Gaussian variables (note that $\varphi_{k,t;\epsilon}(\tau)$ depends on the fixed signal instance ${\mathbf{s}}$). Using the union bound for the first inequality and Markov’s inequality together with $\exp(\varphi_{k,t;\epsilon}(1))=\mathbb{E}_{\mathbf{Z}}[\exp(-snF_{k,t;\epsilon}(\bm{\Theta}))]$ and $\exp(\varphi_{k,t;\epsilon}(0))=\exp(-sn\mathbb{E}_{\mathbf{Z}}[F_{k,t;\epsilon}(\bm{\Theta})])$ for the second one, one deduces that

[TABLE]

Our essential task is now to prove an upper bound on $|\varphi_{k,t;\epsilon}^{\prime}(\tau)|$. We have

[TABLE]

where

[TABLE]

We then replace this expression in the numerator of (166) and integrate by parts over all standard Gaussian variables of type $z^{(1)}$ and $z^{(2)}$. Doing so generates partial derivatives of the form $\mathbb{E}_{2}[\frac{\partial}{\partial z^{(1)}}\langle-\rangle]$ and $\mathbb{E}_{2}[\frac{\partial}{\partial z^{(2)}}\langle-\rangle]$ as well as derivatives of the form $\frac{\partial}{\partial z^{(1)}}\exp(s\,\mathbb{E}_{2}[\ln{\cal Z}_{k,t;\epsilon}(\tau)])$. A lengthy but straightforward calculation shows that only the later survive. The numerator of (166) becomes

[TABLE]

Working out the partial derivatives yields

[TABLE]

For bounded signals we have $|x_{i}|<M$ as well as $m_{k}\leq M^{2}$. Thus the sum of these four terms is bounded by

[TABLE]

for all $k=1,\dots,K$. This is an upper bound for the numerator of (166), which implies $|\varphi_{k,t;\epsilon}^{\prime}(\tau)|\leq s^{2}n(2M^{4}/\Delta+\epsilon M^{2}/2)$. From (165)

[TABLE]

and the best possible value $s=u(M^{4}/\Delta+\epsilon M^{2}/2)^{-1}$ yields (163) and ends the proof. ∎

The second lemma expresses concentration w.r.t the input signal of the free energy averaged over the Gaussian disorder. Recall that $\mathbb{P}_{{\mathbf{S}}}$ is the probability law w.r.t the signal realisation.

Lemma 8** (Concentration w.r.t the signal realisation).**

Take $P_{0}$ with bounded support in $[-M,M]$. For all $k=1,\ldots,K$, $t\in[0,1]$, and $\epsilon\in[0,1]$ we have

[TABLE]

where $u>0$.

Proof.

We first prove a bounded difference property on $\mathbb{E}_{{\mathbf{Z}}}[F_{k,t;\epsilon}(\bm{\Theta})]$ and then apply the McDiarmid inequality [55, 56]. Let ${\mathbf{s}}$ and ${\mathbf{s}}^{\prime}$ two signal realisations that differ at the component $i$ only, i.e. $s_{j}=s_{j}^{\prime}$ for $j\neq i$. We first consider the difference of Hamiltonians corresponding to these two realisations. From (2.2)–(21) we have

[TABLE]

For a signal distribution with bounded support $[-M,M]$ we get (recall $|m_{k}|\leq M^{2}$)

[TABLE]

Now set $g(s_{1},\ldots,s_{n})\vcentcolon=\mathbb{E}_{{\mathbf{Z}}}[F_{k,t;\epsilon}(\bm{\Theta})]$. We have (here ${\mathbf{X}}\sim P_{0}$)

[TABLE]

and since from (171)

[TABLE]

we readily obtain

[TABLE]

with $c_{i}=2(M^{4}/\Delta+\epsilon M^{2})/n$, $i=1,\ldots,n$. McDiarmid’s inequality states that

[TABLE]

which here reads (169) and ends the proof of the lemma. ∎

Proof of Proposition 6

From the triangle inequality and the union bound

[TABLE]

where the last inequality comes from Lemmas 7 and 8. $\blacksquare$

Appendix A Linking the perturbed and plain free energies

The purpose of this appendix is to prove Lemma 1. We first note that differentiating the function $\epsilon\mapsto f_{k=1,t=0;\epsilon}$ in (24)

[TABLE]

By a Gaussian integration by parts the last term becomes

[TABLE]

By an application of the identity (47) we have $\mathbb{E}[\langle X_{i}\rangle_{1,0;\epsilon}S_{i}]=\mathbb{E}[\langle X_{i}\rangle_{1,0;\epsilon}^{2}]$. Therefore we find

[TABLE]

Now by convexity and (47) we have $\mathbb{E}[\langle X_{i}\rangle_{1,0;\epsilon}^{2}]\leq\mathbb{E}[\langle X_{i}^{2}\rangle_{1,0;\epsilon}]=\mathbb{E}[S^{2}]$. Therefore

[TABLE]

and the first inequality of the Lemma follows from an application of the mean value theorem.

The second inequality follows from the Lipschitz continuity of the free energy $f_{k=K,t=1;\epsilon}$ of the decoupled scalar system. We refer to [57] for the proof of this standard fact.

Appendix B Proof of Lemma 3

The proof of this lemma uses another interpolation:

[TABLE]

where ${\mathbf{X}},{\mathbf{X}}^{\prime},{\mathbf{X}}^{\prime\prime}$ etc are i.i.d replicas distributed according to (22). Computations similar to those in sec. 2.7 lead to

[TABLE]

where we define

[TABLE]

Finally from (182) and Cauchy-Schwarz, one obtains

[TABLE]

The last equality is true as long as the prior $P_{0}$ has bounded first four moments. We prove this claim now. Let us start by studying $\mathbb{E}[\langle q_{{\mathbf{X}},{\mathbf{S}}}^{2}\rangle_{k,s;\epsilon}]$. Using Cauchy-Schwarz for the inequality and (47) for the subsequent equality,

[TABLE]

where the last equality is valid for $P_{0}$ with bounded second and fourth moments. For $\mathbb{E}[\langle g({\mathbf{X}},{\mathbf{X}}^{\prime};{\mathbf{S}})^{2}\rangle_{k,s;\epsilon}]$ we proceed similarly by decoupling the expectations using Cauchy-Schwarz and then using (47) to make appear only terms depending on the signal ${\mathbf{s}}$. One finds that under the same conditions on the moments of $P_{0}$ we have $\mathbb{E}[\langle g({\mathbf{X}},{\mathbf{X}}^{\prime};{\mathbf{S}})^{2}\rangle_{k,s;\epsilon}]={\cal O}(n^{2})$. Combined with (185) leads to the last equality of (184) and ends the proof.

Appendix C Alternative argument for the lower bound

We present an alternative useful, albeit not completely rigorous, argument to obtain the lower bound (46). With enough work the argument can be made rigorous. Note that defining

[TABLE]

the identity (2.6) is equivalent to

[TABLE]

Setting $b_{n}=2a_{n}$, taking a sequence $a_{n}\to 0$ slowly enough as $n\to+\infty$, using Lemma 1 and (28), we obtain

[TABLE]

Simple algebra starting from $\partial_{m_{k}}\widetilde{f}_{\rm RS}(\{m_{k}\}_{k=1}^{K_{n}};\Delta)=0$ implies, under the assumption that the extrema are attained at interior points of $\mathbb{R}^{K_{n}}_{+}$ (the point to work out to make the argument rigorous), that the minimizer of $\widetilde{f}_{\rm RS}(\{m_{k}\}_{k=1}^{K_{n}};\Delta)$ satisfies

[TABLE]

The right hand side is independent of $k$, thus the minimizer is $m_{k}=m_{*}$ for $k=1,\cdots,K_{n}$ where

[TABLE]

Thus

[TABLE]

From (188) we get

[TABLE]

which is the inequality (46).

Appendix D A consequence of Bayes rule

The purpose of this appendix is to prove the identity (47). Recall that the Gibbs bracket $\langle-\rangle_{k,t;\epsilon}$ is the average with respect to the posterior $P_{k,t;\epsilon}({\mathbf{x}}|\bm{\theta})$ where $\bm{\theta}\vcentcolon=\{{\mathbf{s}},\{{\mathbf{z}}^{(k)},\widetilde{{\mathbf{z}}}^{(k)}\}_{k=1}^{K},\widehat{\mathbf{z}}\}$. Using Bayes law we have:

[TABLE]

It remains to notice that

[TABLE]

where the Gibbs bracket on the right hand side is an average with respect to the product measure of two posteriors $P_{k,t;\epsilon}({\mathbf{x}}|\bm{\theta})P_{k,t;\epsilon}({\mathbf{x}}^{\prime}|\bm{\theta})$.

Appendix E A stochastic calculus interpretation

We note that the proofs do not require any upper limit on $K$. This suggests that it is possible to formulate the adaptive interpolation method entirely in a continuum language. Here we informally show this for the simplest problem, namely symmetric rank-one matrix factorisation, and plan to come back to a rigorous formulation of the continuum formulation in future work.

It is helpful to first write down explicitly the $(k,t)$–interpolating Hamiltonian (2.2) (leaving out the perturbation in (21) which is irrelevant for the argument here)

[TABLE]

and to define the step-wise function $m(u)=m_{k^{\prime}}$ for $k^{\prime}/K\leq u<(k^{\prime}+1)/K$, $k^{\prime}=1,\dots,K$.

Let us first look at the terms that do not involve Gaussian noise and become simple Riemann integrals. We have for the contribution coming from (195) and (197),

[TABLE]

Similarly, we have for the terms coming from (196) and (198),

[TABLE]

Now we treat the more interesting contributions involving the Gaussian noise. Let $B(u)$ be the Wiener process defined by $B(0)=0$, $\mathbb{E}[B(u)]=0$, $\mathbb{E}[B(u)B(v)]=\min(u,v)$ for $u,v\in\mathbb{R}_{+}$. We introduce independent copies $B_{ij}(u)$, $i,j=1,\dots,n$ and consider the sum of increments (also written as an Ito integral)

[TABLE]

Since the increments are independent and $\mathbb{E}[(B(u)-B(v))^{2}]=|u-v|$, this is a Gaussian random variable with zero mean and variance $(K+1-k-t)/K$. It is therefore equal in distribution to

[TABLE]

and the contribution of the (random) Gaussian noise in (195) and (197) becomes

[TABLE]

To represent the contributions of (196), (198) we introduce independent copies of the Wiener process $\widetilde{B}_{i}(u)$, $i=1,\dots,n$ and form the Ito integral

[TABLE]

which has the same variance than

[TABLE]

Indeed

[TABLE]

Therefore the contribution of (196) and (198) can be represented as

[TABLE]

Finally, collecting (199), (200), (203), (207), setting $\tau\vcentcolon=(t+k)/K$ and $K\to\infty$, we obtain a continuous form of the random $(k,t)$–interpolating Hamiltonian,

[TABLE]

where $m(u)$ is an arbitrary trial function and $\mathbf{B}$ denotes the collection of all Wiener processes. Note that $\int_{\tau}^{1}du\,B_{ij}(u)=B_{ij}(1)-B_{ij}(\tau)$ which is distributed as $\sqrt{1-\tau}Z_{ij}$ for $Z_{ij}\sim\mathcal{N}(0,1)$, and $\int_{0}^{\tau}\sqrt{m(u)}d\widetilde{B}_{i}(u)$ is distributed as $\sqrt{\int_{0}^{\tau}m(u)}\widetilde{Z}_{i}$ for $\widetilde{Z}_{i}\sim\mathcal{N}(0,1)$. Therefore (208) is equal in distribution to

[TABLE]

Clearly, the usual Guerra-Toninelli interpolation appears as a special case where one chooses a constant trial function $m(u)=m$ constant. When we go from (208) to (209) we eliminate completely the Wiener process, however we believe it is useful to keep in mind the point of view expressed by (208) which may turn out to be important for more complicated problems.

Starting from (208) or (209) it is possible to evaluate the free energy change along the interpolation path. We define the free energy

[TABLE]

For $\tau=0$ using we recover the original Hamiltonian $\mathcal{H}_{k=1,t=0}$ (see (26)) and $f(0)=f$ given in (7). For $\tau=1$ setting $\int_{0}^{1}du\,m(u)=m_{\rm mf}$ we recover the mean-field Hamiltonian $\mathcal{H}_{k=K,t=1}$ (see (31)) and $f(1)=f_{\rm den}(\Sigma(\int_{0}^{1}du\,m(u));\Delta)$. Then proceeding similarly to sec. 2.7 one finds the identity

[TABLE]

where $\langle-\rangle_{\tau}$ is the Gibbs average w.r.t (208).

Of course this immediately gives the upper bound in Proposition 1. The matching lower bound is obtained by the same ideas used in the discrete version. We briefly review them informally in the continuous language. One first introduces the $\epsilon$-perturbation term (21) and proves a concentration property for the overlap analogous to Lemma 2. Starting with the continuous version of the interpolating Hamiltonian the proof of the free energy concentration is essentially identical (even simpler) than in sec. 7, which implies the overlap concentration through sec. 5 that is unchanged. Then, the square in the remainder term is approximately equal to $(\mathbb{E}_{{\mathbf{S}},\mathbf{B}}[\langle q_{{\mathbf{X}},{\mathbf{S}}}\rangle_{\tau,\epsilon}]-m(\tau))^{2}$ and we make it vanish by choosing

[TABLE]

This continuous setting thus allows to avoid proving Lemma 3. This then easily yields the lower bound in Proposition 1. One must still check that (212) has a solution. The right hand side is a function $G_{n,\epsilon}(\tau;\int_{0}^{\tau}du\,m(u))$ so setting $x(\tau)=\int_{0}^{\tau}du\,m(u)$, $dx/d\tau=m(\tau)$, we recognize that (212) is a first order differential equation with initial condition $x(0)=0$. The existence of a unique global solution on $\tau\in[0,1]$ is then proved using the Cauchy-Lipschitz theorem. Moreover this solution is differentiable and monotone increasing with respect to $\epsilon$. This last step of the analysis replaces Lemma 4.

Acknowledgments

Jean Barbier acknowledges funding by the Swiss National Science Foundation grant no. 200021-156672. We thank Thibault Lesieur for providing us the expression of the RS potential for tensor estimation. We also acknowledge helpful discussions with Olivier Lévêque and Léo Miolane on the stochastic calculus interpretation and continuous version of Appendix E.