# Analysis of Random Sequential Message Passing Algorithms for Approximate
Inference

Burak Çakmak1, Yue M. Lu2 and Manfred Opper3 1 Artificial Intelligence Group,
Technische Universität Berlin, 10587, Germany  2 John A. Paulson School of
Engineering and Applied Sciences, Harvard University, Cambridge, MA 02138, USA
3 Centre for Systems Modelling and Quantitative Biomedicine, University of
Birmingham, B15 2TT, United Kingdom

###### Abstract

We analyze the dynamics of a random sequential message passing algorithm for
approximate inference with large Gaussian latent variable models in a student-
teacher scenario. To model nontrivial dependencies between the latent
variables, we assume random covariance matrices drawn from rotation invariant
ensembles. Moreover, we consider a model mismatching setting, where the
teacher model and the one used by the student may be different. By means of
dynamical functional approach, we obtain exact dynamical mean-field equations
characterizing the dynamics of the inference algorithm. We also derive a range
of model parameters for which the sequential algorithm does not converge. The
boundary of this parameter range coincides with the _de Almeida Thouless_ (AT)
stability condition of the _replica symmetric_ ansatz for the static
probabilistic model.

Keywords: Bayesian Inference, Iterative Algorithms, Approximate message
passing, TAP Equations, Random Matrices, Dynamical Functional Theory

##  1 Introduction

The analysis of the dynamics of message passing algorithms for inference in
large probabilistic models has attracted considerable interest in the fields
of statistical physics and information sciences [1, 2, 3, 4, 5, 6, 7, 8, 9,
10, 11]. From a statistical physics point of view, the fixed points of such
algorithms correspond to solutions of TAP mean field equations for disordered
systems [12, 13, 14, 15]. The latter, under some conditions on the statistics
of the disorder, can lead to exact solutions to thermal averages in the large
system limit. Hence, message passing algorithms provide efficient computation
methods for obtaining accurate solutions to high–dimensional statistical
inference problems.

So far, most of the theoretical works on the dynamics of message passing
consider a parallel update scheme, where all dynamical nodes are updated
simultaneously at each iteration of the algorithm. For large classes of the
random interaction matrices, the exact temporal progress of the algorithm can
then be described by the so-called state-evolution equations [1, 2, 16].

In many applications, the parallel dynamics of the algorithm is often replaced
by a sequential version, where only a subset of nodes is updated per
iteration. For example, Minka’s EP (expectation propagation) algorithm [17],
which is one of the motivations behind the so-called VAMP (vector approximate
message passing) approach [18, 9], is originally formulated in terms of
sequential iterations. This type of sequential algorithms have lower
computational complexity per iteration. They can also be more memory efficient
as they only need to have access to a small batch of the available data at any
given time. Moreover, in certain situations, they were found to improve the
convergence properties [19]. Our goal in this paper is to extend the
theoretical analysis of message passing dynamics from the parallel update
setting to the sequential setting. Specifically, we address the following
issues:

  1. 1.

We analyze the dynamics of a random sequential message passing algorithm for
approximate inference with a large Gaussian latent variable model. At each
iteration, a random selection of nodes are updated by the algorithm. The
probability for a given node to be included in an update is a free parameter.
Varying this parameter allows for an interpolation between a full parallel
update of all nodes and the case where on average only a single node is
updated. Relying on the technique of the dynamical functional approach of
statistical mechanics [20], we decoupled the degrees of freedom and derive an
effective _single node_ evolution equation that characterizes the limiting
dynamics of the sequential algorithm.

  2. 2.

In practice, the probabilistic model assumed by the inference algorithm may
differ significantly from the real data generating process. We take into
account this issue by allowing for a possible mismatch between the data
generating teacher model and the model used by the student. From a technical
point of view, this more general scenario requires a larger number of time
dependent order parameters to describe the dynamics of the algorithm. In
addition, unlike the case of perfect match between the student and teacher
models [21], the message-passing algorithm is no longer guaranteed to converge
in the mismatched case. We have identified a range of model parameters for
which the convergence of the sequential algorithm is impossible.
Interestingly, the boundary of this parameter range coincides with the _de
Almeida Thouless_ (AT) stability condition of the _replica symmetric_ ansatz
for the probabilistic model [14, 15].

There have been several earlier studies of sequential dynamics for solving
various statistical physics and inference problems [22, 23, 24]. The effective
single node dynamics obtained in these studies often contain memory terms that
make it difficult to evaluate the two-time correlation functions. Remarkably,
due to the construction of our message passing algorithm, its single node
dynamics has no _memory term_. As a result, the corresponding two–time
correlation functions can be obtained by tractable recursion formulas. A
similar “memory-free” property of sequential algorithms was observed in our
previous paper [25] on solving the TAP equations for the
Sherrington–Kirkpatrick model. Finally, the issue of data-model mismatch has
also been previously considered in [26] for parallel-updating message passing
algorithms. Unlike in [26] where the analysis is focused on the “single-time”
statistics of the algorithm, we characterize the full effective single-node
dynamics. This characterization provides information about the joint
statistics of the algorithm over multiple time steps, which is crucial for
analysing the convergence properties of the message passing algorithm.

The paper is organized as follows: Section 2 presents the details of the
Bayesian probabilistic model considered in this work. We introduce in Section
3 a random sequential iterative algorithm for solving the inference problem.
Its thermodynamic properties are studied in Section 4 by using the method of
dynamical functional theory. Comparisons of the theory with simulations are
given in Section 5. We conclude the paper in Section 6 with a summary and some
discussions. The derivations of our results can be found in the Appendix.

##  2 Latent Gaussian variable models

Message passing algorithms have been successfully applied to latent Gaussian
variable models [17, 27, 28]. This class of models finds widespread
applications in statistics, machine learning and signal processing. A typical
scenario is to infer an unobserved latent vector
${\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}}\in\Re^{N\times
1}$ by using the Bayesian posterior distribution

$p({\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}}|{\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}},{\mathchoice{\mbox{\boldmath$\displaystyle
K$}}{\mbox{\boldmath$\textstyle K$}}{\mbox{\boldmath$\scriptstyle
K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}})\doteq\frac{1}{Z}\mathcal{N}({\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}}|{\mathchoice{\mbox{\boldmath$\displaystyle
0$}}{\mbox{\boldmath$\textstyle 0$}}{\mbox{\boldmath$\scriptstyle
0$}}{\mbox{\boldmath$\scriptscriptstyle
0$}}},{\mathchoice{\mbox{\boldmath$\displaystyle
K$}}{\mbox{\boldmath$\textstyle K$}}{\mbox{\boldmath$\scriptstyle
K$}}{\mbox{\boldmath$\scriptscriptstyle K$}}})\prod_{i\leq
N}p(y_{i}|\theta_{i})$ (1)

where $Z$ is a normalization constant. This model assumes that the components
of the vector $\textstyle y$ of $N$ real data values are assumed to be
generated independently from a likelihood $p(y|\theta)$ based on a vector of
unknown parameters $\textstyle\theta$. Prior statistical knowledge about
$\textstyle\theta$ is introduced by the correlated Gaussian with covariance
${\mathchoice{\mbox{\boldmath$\displaystyle K$}}{\mbox{\boldmath$\textstyle
K$}}{\mbox{\boldmath$\scriptstyle K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}}\in\Re^{N\times N}$.

We will later illustrate our theory on the well known example of Bayesian
learning of a noisy perceptron—also known as probit regression [29]. This
corresponds to a binary classification problem with class labels $y_{i}=\pm
1$. For this model, one assumes a training set given by
$\\{({\mathchoice{\mbox{\boldmath$\displaystyle
x$}}{\mbox{\boldmath$\textstyle x$}}{\mbox{\boldmath$\scriptstyle
x$}}{\mbox{\boldmath$\scriptscriptstyle x$}}}_{i},y_{i})\\}_{i\leq N}$ where
${\mathchoice{\mbox{\boldmath$\displaystyle x$}}{\mbox{\boldmath$\textstyle
x$}}{\mbox{\boldmath$\scriptstyle x$}}{\mbox{\boldmath$\scriptscriptstyle
x$}}}_{i}\in\mathbb{R}^{P\times 1}$ stands for a vector of inputs. Class
labels $y_{i}$ are generated according to the observation model

$y_{i}=\epsilon_{i}\,{\rm sign}({\mathchoice{\mbox{\boldmath$\displaystyle
x$}}{\mbox{\boldmath$\textstyle x$}}{\mbox{\boldmath$\scriptstyle
x$}}{\mbox{\boldmath$\scriptscriptstyle
x$}}}_{i}^{\top}{\mathchoice{\mbox{\boldmath$\displaystyle
w$}}{\mbox{\boldmath$\textstyle w$}}{\mbox{\boldmath$\scriptstyle
w$}}{\mbox{\boldmath$\scriptscriptstyle w$}}}+n_{i}).$ (2)

Here, we allow for additive i.i.d. Gaussian noises $n_{i}$ for all $i$ with
${\mathchoice{\mbox{\boldmath$\displaystyle n$}}{\mbox{\boldmath$\textstyle
n$}}{\mbox{\boldmath$\scriptstyle n$}}{\mbox{\boldmath$\scriptscriptstyle
n$}}}\sim\mathcal{N}({\mathchoice{\mbox{\boldmath$\displaystyle
0$}}{\mbox{\boldmath$\textstyle 0$}}{\mbox{\boldmath$\scriptstyle
0$}}{\mbox{\boldmath$\scriptscriptstyle 0$}}},\sigma^{2}{\bf I})$ as well as
i.i.d. multiplicative flip noises $\epsilon_{i}=\pm 1$ for all $i$ with
$\beta\doteq{\rm Pr}(\epsilon_{i}=-1)$. We assume a latent vector $\textstyle
w$ with Gaussian prior distribution
${\mathchoice{\mbox{\boldmath$\displaystyle w$}}{\mbox{\boldmath$\textstyle
w$}}{\mbox{\boldmath$\scriptstyle w$}}{\mbox{\boldmath$\scriptscriptstyle
w$}}}\sim\mathcal{N}({\mathchoice{\mbox{\boldmath$\displaystyle
0$}}{\mbox{\boldmath$\textstyle 0$}}{\mbox{\boldmath$\scriptstyle
0$}}{\mbox{\boldmath$\scriptscriptstyle 0$}}},{\bf I})$. To map this problem
onto the model (1), we introduce the latent vector
${\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}}={\mathchoice{\mbox{\boldmath$\displaystyle
X$}}{\mbox{\boldmath$\textstyle X$}}{\mbox{\boldmath$\scriptstyle
X$}}{\mbox{\boldmath$\scriptscriptstyle
X$}}}{\mathchoice{\mbox{\boldmath$\displaystyle\omega$}}{\mbox{\boldmath$\textstyle\omega$}}{\mbox{\boldmath$\scriptstyle\omega$}}{\mbox{\boldmath$\scriptscriptstyle\omega$}}}$
with ${\mathchoice{\mbox{\boldmath$\displaystyle
X$}}{\mbox{\boldmath$\textstyle X$}}{\mbox{\boldmath$\scriptstyle
X$}}{\mbox{\boldmath$\scriptscriptstyle
X$}}}\doteq[{\mathchoice{\mbox{\boldmath$\displaystyle
x$}}{\mbox{\boldmath$\textstyle x$}}{\mbox{\boldmath$\scriptstyle
x$}}{\mbox{\boldmath$\scriptscriptstyle
x$}}}_{1}^{\top},{\mathchoice{\mbox{\boldmath$\displaystyle
x$}}{\mbox{\boldmath$\textstyle x$}}{\mbox{\boldmath$\scriptstyle
x$}}{\mbox{\boldmath$\scriptscriptstyle
x$}}}_{2}^{\top},\cdots,{\mathchoice{\mbox{\boldmath$\displaystyle
x$}}{\mbox{\boldmath$\textstyle x$}}{\mbox{\boldmath$\scriptstyle
x$}}{\mbox{\boldmath$\scriptscriptstyle x$}}}_{N}^{\top}]$. Hence, the prior
covariance of $\textstyle\theta$ equals
${\mathchoice{\mbox{\boldmath$\displaystyle K$}}{\mbox{\boldmath$\textstyle
K$}}{\mbox{\boldmath$\scriptstyle K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}}={\mathchoice{\mbox{\boldmath$\displaystyle
X$}}{\mbox{\boldmath$\textstyle X$}}{\mbox{\boldmath$\scriptstyle
X$}}{\mbox{\boldmath$\scriptscriptstyle
X$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
X$}}{\mbox{\boldmath$\textstyle X$}}{\mbox{\boldmath$\scriptstyle
X$}}{\mbox{\boldmath$\scriptscriptstyle X$}}}^{\top}$ and we have the data
likelihood function

$p(y|\theta)=(1-\beta)\Phi\left(\frac{y\theta}{\sigma}\right)+\beta\Phi\left(-\frac{y\theta}{\sigma}\right)$
(3)

where $\Phi(\cdot)$ denotes the cumulative distribution function of the
standard normal distribution.

##  3 The random sequential VAMP algorithm

Typical prediction tasks based on observed data involve the computations of
expectations of components of $\textstyle\theta$ (or of functions of these
components) using the posterior (1). Unfortunately, except for simple Gaussian
likelihoods or simple diagonal covariance matrices, such expectations lead to
multi–dimensional integrals which cannot be computed analytically. Hence, one
has to resort to approximations. To be able to obtain reliable results in the
case of high–dimensional vectors $\textstyle\theta$, so–called message passing
algorithms have been developed which provide efficient iterative computations
of generalized mean field approximations to the desired expectations.

Given the auxiliary single-site partition function

$\displaystyle Z_{\nu}(\gamma,y)\doteq\int{\rm
d}\theta\;p(y|\theta)e^{-\frac{\nu}{2}\theta^{2}+\gamma\theta}$ (4)

the logarithmic derivatives

$\displaystyle m_{\nu}(\gamma,y)$
$\displaystyle\doteq\frac{{\partial}\ln{Z_{\nu}(\gamma,y)}}{{\partial}\gamma}\approx\mathbb{E}[\theta|{\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}},{\mathchoice{\mbox{\boldmath$\displaystyle
K$}}{\mbox{\boldmath$\textstyle K$}}{\mbox{\boldmath$\scriptstyle
K$}}{\mbox{\boldmath$\scriptscriptstyle K$}}}]$ (5) $\displaystyle
m_{\nu}^{\prime}(\gamma,y)$
$\displaystyle\doteq\frac{{\partial}m_{\nu}(\gamma,y)}{{\partial}\gamma}\approx\mathbb{E}[\theta^{2}|{\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}},{\mathchoice{\mbox{\boldmath$\displaystyle
K$}}{\mbox{\boldmath$\textstyle K$}}{\mbox{\boldmath$\scriptstyle
K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}}]-\mathbb{E}^{2}[\theta|{\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}},{\mathchoice{\mbox{\boldmath$\displaystyle
K$}}{\mbox{\boldmath$\textstyle K$}}{\mbox{\boldmath$\scriptstyle
K$}}{\mbox{\boldmath$\scriptscriptstyle K$}}}]$ (6)

provide approximations to posterior mean and variances of single components
$\theta\equiv\theta_{i}$ upon convergence of the algorithm. The mean
$\gamma_{i}$ of the cavity fields of a node [14, 30] can be computed
iteratively by the VAMP algorithm [18, 9]. In the following, we introduce a
random sequential version of the usual parallel VAMP. Before the iteration
starts, we compute the spectral decomposition

${\mathchoice{\mbox{\boldmath$\displaystyle K$}}{\mbox{\boldmath$\textstyle
K$}}{\mbox{\boldmath$\scriptstyle K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}}={\mathchoice{\mbox{\boldmath$\displaystyle
O$}}{\mbox{\boldmath$\textstyle O$}}{\mbox{\boldmath$\scriptstyle
O$}}{\mbox{\boldmath$\scriptscriptstyle
O$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
D$}}{\mbox{\boldmath$\textstyle D$}}{\mbox{\boldmath$\scriptstyle
D$}}{\mbox{\boldmath$\scriptscriptstyle
D$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
O$}}{\mbox{\boldmath$\textstyle O$}}{\mbox{\boldmath$\scriptstyle
O$}}{\mbox{\boldmath$\scriptscriptstyle O$}}}^{\top}$ (7)

where $\textstyle D$ is diagonal with the diagonal entries being the
eigenvalues of $\textstyle K$. We define the iterative algorithm in discrete
time by the vector updates for $t=1,2,\ldots T$ by

$\displaystyle{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}^{(t)}$
$\displaystyle={\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}^{(t-1)}+{\mathchoice{\mbox{\boldmath$\displaystyle
P$}}{\mbox{\boldmath$\textstyle P$}}{\mbox{\boldmath$\scriptstyle
P$}}{\mbox{\boldmath$\scriptscriptstyle
P$}}}^{(t)}[{\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}-{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}^{(t-1)}]$
(8a)
$\displaystyle{\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}$
$\displaystyle=\frac{1}{\tau^{(t)}}{\mathchoice{\mbox{\boldmath$\displaystyle
O$}}{\mbox{\boldmath$\textstyle O$}}{\mbox{\boldmath$\scriptstyle
O$}}{\mbox{\boldmath$\scriptscriptstyle
O$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
D$}}{\mbox{\boldmath$\textstyle D$}}{\mbox{\boldmath$\scriptstyle
D$}}{\mbox{\boldmath$\scriptscriptstyle
D$}}}(\lambda^{(t)}{\mathchoice{\mbox{\boldmath$\displaystyle
D$}}{\mbox{\boldmath$\textstyle D$}}{\mbox{\boldmath$\scriptstyle
D$}}{\mbox{\boldmath$\scriptscriptstyle D$}}}+{\bf
I})^{-1}{\mathchoice{\mbox{\boldmath$\displaystyle
O$}}{\mbox{\boldmath$\textstyle O$}}{\mbox{\boldmath$\scriptstyle
O$}}{\mbox{\boldmath$\scriptscriptstyle
O$}}}^{\top}{\mathchoice{\mbox{\boldmath$\displaystyle\tilde{\gamma}$}}{\mbox{\boldmath$\textstyle\tilde{\gamma}$}}{\mbox{\boldmath$\scriptstyle\tilde{\gamma}$}}{\mbox{\boldmath$\scriptscriptstyle\tilde{\gamma}$}}}^{(t)}-\tilde{{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}}^{(t)}$
(8b)
$\displaystyle{\mathchoice{\mbox{\boldmath$\displaystyle\tilde{\gamma}$}}{\mbox{\boldmath$\textstyle\tilde{\gamma}$}}{\mbox{\boldmath$\scriptstyle\tilde{\gamma}$}}{\mbox{\boldmath$\scriptscriptstyle\tilde{\gamma}$}}}^{(t)}$
$\displaystyle=\frac{m_{\nu^{(t-1)}}({\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}^{(t-1)},{\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}})}{\chi^{(t)}}-{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}^{(t-1)}.$
(8c)

The scalar quantities $\chi^{(t)}$, $\lambda^{(t)}$, $\tau^{(t)}$, and
$\nu^{(t)}$ are updated as

$\displaystyle\chi^{(t)}$ $\displaystyle=\langle
m_{\nu^{(t-1)}}^{\prime}({\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}^{(t-1)},{\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle y$}}})\rangle$ (9a)
$\displaystyle\lambda^{(t)}$ $\displaystyle=\frac{1}{\chi^{(t)}}-\nu^{(t-1)}$
(9b) $\displaystyle\tau^{(t)}$ $\displaystyle=\frac{1}{N}\sum_{i\leq
N}\frac{D_{ii}}{\lambda^{(t)}D_{ii}+1}$ (9c) $\displaystyle\nu^{(t)}$
$\displaystyle=\frac{1}{\tau^{(t)}}-\lambda^{(t)}$ (9d)

where the brackets $\langle\ldots\rangle$ denote an empirical average over the
sites. Moreover, we consider a random initialization for
${\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}^{(0)}$
from an i.i.d. normal Gaussian distribution. The diagonal matrix
${\mathchoice{\mbox{\boldmath$\displaystyle P$}}{\mbox{\boldmath$\textstyle
P$}}{\mbox{\boldmath$\scriptstyle P$}}{\mbox{\boldmath$\scriptscriptstyle
P$}}}^{(t)}$ in (8a) is composed of binary decision variables
$p_{i}^{(t)}\doteq P^{(t)}_{ii}\in\\{0,1\\}$. The original parallel version of
the VAMP algorithm is obtained when
${\mathchoice{\mbox{\boldmath$\displaystyle P$}}{\mbox{\boldmath$\textstyle
P$}}{\mbox{\boldmath$\scriptstyle P$}}{\mbox{\boldmath$\scriptscriptstyle
P$}}}^{(t)}$ is equal to the unit matrix. Random sequential updates are
introduced by making the $p_{i}^{(t)}$ random variables which decide if node
$i$ is updated ($p_{i}^{(t)}=1$) at time $k$ or not ($p_{i}^{(t)}=0$). We
assume that the $p_{i}^{(t)}$ are independent for all $i,t$ and that ${\rm
Pr}(p_{i}^{(t)}=1)=\eta$. The case $\eta=1/N$ corresponds to updating only a
single node on average.

##  4 The dynamical mean-field equations

We consider an average case analysis of the algorithm in the limit
$N\to\infty$ assuming that the data $\textstyle y$ is generated from a given
likelihood model, where for generality, we consider a data-model mismatching
scenario. In general, we assume that the components of the vector $\textstyle
y$ are generated independently from a likelihood $p_{0}(y|\theta)$ which is
not necessarily equal to $p(y|\theta)$. For the noisy perception model this
would correspond to different sets of hyperparameters $\beta$ and $\sigma$ in
(3). But we also assume that the prior covariance matrix $\textstyle K$ is the
same for both true parameter and the parameter in the inference (student)
model. We choose $\textstyle K$ to be a random matrix with a rotational
invariant distribution. This means that the matrix $\textstyle O$ in (8) is
assumed to be a Haar matrix, i.e. a random rotation. In this way, it is
possible to model matrices $\textstyle K$ with nontrivial (weak) dependencies
between entries.

Following previous studies [5, 21], we derive an effective dynamics of a
single node. This is obtained by averaging the generating functional of the
dynamics over the randomness of $\textstyle y$, $\textstyle O$ and
$\\{{\mathchoice{\mbox{\boldmath$\displaystyle P$}}{\mbox{\boldmath$\textstyle
P$}}{\mbox{\boldmath$\scriptstyle P$}}{\mbox{\boldmath$\scriptscriptstyle
P$}}}^{(t)}\\}$ and a subsequent decoupling of the degrees of freedom. This
involves order parameter functions which are self-averaging for $N\to\infty$.
Generating functionals are partition functions for the computation of
expectations of dynamical variables where the dynamics is included in terms of
Dirac $\delta$ functions. To avoid cluttered notation, and with
${\mathchoice{\mbox{\boldmath$\displaystyle J$}}{\mbox{\boldmath$\textstyle
J$}}{\mbox{\boldmath$\scriptstyle J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}}\doteq{\mathchoice{\mbox{\boldmath$\displaystyle
K$}}{\mbox{\boldmath$\textstyle K$}}{\mbox{\boldmath$\scriptstyle
K$}}{\mbox{\boldmath$\scriptscriptstyle K$}}}^{-1}$ 111 Unless the covariance
matrix $\textstyle K$ has an inverse, one can consider the substitution
${\mathchoice{\mbox{\boldmath$\displaystyle K$}}{\mbox{\boldmath$\textstyle
K$}}{\mbox{\boldmath$\scriptstyle K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}}\to{\mathchoice{\mbox{\boldmath$\displaystyle
K$}}{\mbox{\boldmath$\textstyle K$}}{\mbox{\boldmath$\scriptstyle
K$}}{\mbox{\boldmath$\scriptscriptstyle K$}}}+\epsilon{\bf I}$ for
$\epsilon>0$ and perform the limit $\epsilon\to 0$ at the end of the analysis.
In any case, the need for the inverse
${\mathchoice{\mbox{\boldmath$\displaystyle K$}}{\mbox{\boldmath$\textstyle
K$}}{\mbox{\boldmath$\scriptstyle K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}}^{-1}$ will be bypassed in the analysis. Hence, without loss of
generality, we can assume that $\textstyle K$ has an inverse., the dynamical
functionals corresponding to (8) and (9) for $T$ discrete time steps can be
written in the form

$\displaystyle Z_{i}(\\{l^{(t)}\\})=$
$\displaystyle\int\prod_{t=1}^{T}\left\\{{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(t)}{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle m$}}{\mbox{\boldmath$\textstyle
m$}}{\mbox{\boldmath$\scriptstyle m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}\;\delta\left[{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}-f_{t}\left(\\{{\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(t)},{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(l)},{\mathchoice{\mbox{\boldmath$\displaystyle
P$}}{\mbox{\boldmath$\textstyle P$}}{\mbox{\boldmath$\scriptstyle
P$}}{\mbox{\boldmath$\scriptscriptstyle
P$}}}^{(l)}\\}_{l=1}^{t};{\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle y$}}}\right)\right]\right.$
$\displaystyle\left.\qquad\qquad\qquad\qquad\times\delta({\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(t)}-{\mathchoice{\mbox{\boldmath$\displaystyle
J$}}{\mbox{\boldmath$\textstyle J$}}{\mbox{\boldmath$\scriptstyle
J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle m$}}}^{(t)}){\rm
e}^{i\psi_{i}^{(t)}l^{(t)}}\right\\}$ (10)

where $\\{f_{t}\\}$ is an appropriate sequence of non-linear scalar functions.
Using the Fourier representation of the Dirac measures the averaged generating
functional is of the form

$\displaystyle\mathbb{E}[Z_{i}(\\{l^{(t)}\\})]=\int{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}}{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle y$}}{\mbox{\boldmath$\textstyle
y$}}{\mbox{\boldmath$\scriptstyle y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}}{\rm dP}({\mathchoice{\mbox{\boldmath$\displaystyle
O$}}{\mbox{\boldmath$\textstyle O$}}{\mbox{\boldmath$\scriptstyle
O$}}{\mbox{\boldmath$\scriptscriptstyle
O$}}})\;p_{0}({\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}}|{\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}})\mathcal{N}({\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}}|{\mathchoice{\mbox{\boldmath$\displaystyle
0$}}{\mbox{\boldmath$\textstyle 0$}}{\mbox{\boldmath$\scriptstyle
0$}}{\mbox{\boldmath$\scriptscriptstyle
0$}}},{\mathchoice{\mbox{\boldmath$\displaystyle
K$}}{\mbox{\boldmath$\textstyle K$}}{\mbox{\boldmath$\scriptstyle
K$}}{\mbox{\boldmath$\scriptscriptstyle K$}}})\prod_{t\leq T}{\rm
dP}({\mathchoice{\mbox{\boldmath$\displaystyle P$}}{\mbox{\boldmath$\textstyle
P$}}{\mbox{\boldmath$\scriptstyle P$}}{\mbox{\boldmath$\scriptscriptstyle
P$}}}^{(t)})\;Z_{i}(\\{l^{(t)}\\})$ (11) $\displaystyle=c\int{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}}{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle y$}}{\mbox{\boldmath$\textstyle
y$}}{\mbox{\boldmath$\scriptstyle y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}}\;p_{0}({\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}}|{\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}})\prod_{t\leq
T}{\rm dP}({\mathchoice{\mbox{\boldmath$\displaystyle
P$}}{\mbox{\boldmath$\textstyle P$}}{\mbox{\boldmath$\scriptstyle
P$}}{\mbox{\boldmath$\scriptscriptstyle P$}}}^{(t)}){\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(t)}{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle m$}}{\mbox{\boldmath$\textstyle
m$}}{\mbox{\boldmath$\scriptstyle m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle\hat{\psi}$}}{\mbox{\boldmath$\textstyle\hat{\psi}$}}{\mbox{\boldmath$\scriptstyle\hat{\psi}$}}{\mbox{\boldmath$\scriptscriptstyle\hat{\psi}$}}}^{(t)}\delta\left[{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}-f_{t}\left(\\{{\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(l)},{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(l)},{\mathchoice{\mbox{\boldmath$\displaystyle
P$}}{\mbox{\boldmath$\textstyle P$}}{\mbox{\boldmath$\scriptstyle
P$}}{\mbox{\boldmath$\scriptscriptstyle
P$}}}^{(l)}\\}_{l=1}^{t};{\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle y$}}}\right)\right]$
$\displaystyle\quad\times{\rm e}^{{\rm
i}\sum_{k}({\mathchoice{\mbox{\boldmath$\displaystyle\hat{\psi}$}}{\mbox{\boldmath$\textstyle\hat{\psi}$}}{\mbox{\boldmath$\scriptstyle\hat{\psi}$}}{\mbox{\boldmath$\scriptscriptstyle\hat{\psi}$}}}^{(t)})^{\top}{\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(t)}}{\rm
e}^{{\rm
i}\sum_{t}\psi_{i}^{(t)}l^{(t)}}\mathbb{E}_{\mathchoice{\mbox{\boldmath$\displaystyle
O$}}{\mbox{\boldmath$\textstyle O$}}{\mbox{\boldmath$\scriptstyle
O$}}{\mbox{\boldmath$\scriptscriptstyle O$}}}\left[{\rm
e}^{-\frac{1}{2}{\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}}^{\top}{\mathchoice{\mbox{\boldmath$\displaystyle
J$}}{\mbox{\boldmath$\textstyle J$}}{\mbox{\boldmath$\scriptstyle
J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}}{\mathchoice{\mbox{\boldmath$\displaystyle\theta$}}{\mbox{\boldmath$\textstyle\theta$}}{\mbox{\boldmath$\scriptstyle\theta$}}{\mbox{\boldmath$\scriptscriptstyle\theta$}}}-{\rm
i}\sum_{t\leq
T}({\mathchoice{\mbox{\boldmath$\displaystyle\hat{\psi}$}}{\mbox{\boldmath$\textstyle\hat{\psi}$}}{\mbox{\boldmath$\scriptstyle\hat{\psi}$}}{\mbox{\boldmath$\scriptscriptstyle\hat{\psi}$}}}^{(t)})^{\top}{\mathchoice{\mbox{\boldmath$\displaystyle
J$}}{\mbox{\boldmath$\textstyle J$}}{\mbox{\boldmath$\scriptstyle
J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}}{\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(t)}}\right]$
(12)

where ${\rm dP}({\mathchoice{\mbox{\boldmath$\displaystyle
O$}}{\mbox{\boldmath$\textstyle O$}}{\mbox{\boldmath$\scriptstyle
O$}}{\mbox{\boldmath$\scriptscriptstyle O$}}})$ stands for the Haar invariant
measure of the orthogonal group $O(N)$ and $c$ stands for a nonrandom term to
ensure the normalization property $\mathbb{E}[Z(\\{l^{(t)}=0\\})]=1$.

Appendix A gives a short summary of details and references needed for the
computations of the expectations and the subsequent decoupling of the degrees
of freedom. We find that the effective statistics of an arbitrary single node
$\gamma^{(t)}\equiv\gamma_{i}^{(t)}$ (with similar definitions for other
variables) of the algorithm (8) and (9) is given by the stochastic process

$\displaystyle\left(\theta,y,\phi^{(1:T)}\right)$
$\displaystyle\sim\mathcal{N}(\theta|0,q)p_{0}(y|\theta)\mathcal{N}(\phi^{(1:T)}|\theta\mathcal{\hat{B}},\mathcal{C})$
(13a) $\displaystyle\gamma^{(t)}$
$\displaystyle=\gamma^{(t-1)}+p^{(t)}(\phi^{(t)}+\gamma^{(t-1)})$ (13b)

where for short $\phi^{(1:T)}\doteq(\phi^{(1)},\phi^{(2)},\ldots,\phi^{(T)})$
and $q\doteq\lim_{N\to\infty}\frac{1}{N}{\rm
tr}({\mathchoice{\mbox{\boldmath$\displaystyle K$}}{\mbox{\boldmath$\textstyle
K$}}{\mbox{\boldmath$\scriptstyle K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}})$.

Luckily, similar to previous results [21] obtained for the simpler scenario of
parallel dynamics and matching teacher–student models, the effective dynamics
does not contain memory terms. These terms are often encountered for the
stochastic dynamics of disordered systems [23, 24] and would render the
driving process $\phi^{(1:T)}$ non Gaussian. This would preclude the
computation of explicit analytical results for averages at finite time $t$ and
one would have to resort to Monte–Carlo simulations [31] of the effective
process (13).

The entries of the $T\times 1$ vector $\mathcal{\hat{B}}$ and the $T\times T$
covariance matrix $\mathcal{C}$ are recursively computed according to

$\displaystyle\mathcal{\hat{B}}^{(t)}$
$\displaystyle=\frac{\tau^{(t)}\zeta^{(t)}\mathbb{E}[\theta{\tilde{\gamma}}^{(t)}]}{1-q\tau^{(t)}\zeta^{(t)}}$
(14) $\displaystyle\mathcal{C}^{(t,t^{\prime})}$
$\displaystyle=\frac{\mathcal{D}^{(t,t^{\prime})}+\mathcal{Q}^{(t,t^{\prime})}\left(q\mathcal{\hat{B}}^{(t)}\mathcal{\hat{B}}^{(t^{\prime})}+\mathcal{\hat{B}}^{(t)}\mathbb{E}[\theta\tilde{\gamma}^{(t^{\prime})}]+\mathcal{\hat{B}}^{(t^{\prime})}\mathbb{E}[\theta\tilde{\gamma}^{(t)}]+\mathbb{E}[\tilde{\gamma}^{(t)}\tilde{\gamma}^{(t^{\prime})}]\right)}{1-\mathcal{Q}^{(t,t^{\prime})}}$
(15)

where we have introduced the auxiliary dynamical order parameters

$\displaystyle\zeta^{(t)}$ $\displaystyle=\frac{{\rm
R}(-\tau^{(t)})-1/q}{q-\tau^{(t)}}$ (16a)
$\displaystyle{\mathcal{Q}}^{(t,t^{\prime})}$
$\displaystyle=\tau^{(t)}\tau^{(t^{\prime})}\left\\{\begin{array}[]{cc}\frac{{\rm
R}(-\tau^{(t^{\prime})})-{\rm
R}(-\tau^{(t)})}{\tau^{(t)}-\tau^{(t^{\prime})}}&t\neq t^{\prime}\\\ {\rm
R}^{\prime}(-\tau^{(t)})&\text{else}\end{array}\right.$ (16d)
$\displaystyle\mathcal{D}^{(t,t^{\prime})}$
$\displaystyle=\frac{\hat{\mathcal{B}}^{(t)}\hat{\mathcal{B}}^{(t^{\prime})}}{\zeta^{(t)}\zeta^{(t^{\prime})}}\left\\{\begin{array}[]{cc}\frac{\zeta^{(t)}-\zeta^{(t^{\prime})}}{\tau^{(t)}-\tau^{(t^{\prime})}}&t\neq
t^{\prime}\\\ \frac{\zeta^{(t)}-{\rm
R}^{\prime}(-\tau^{(t)})}{q-\tau^{(t)}}&\text{else}.\end{array}\right.$ (16g)

Here, the function $\rm R(\omega)$ is defined as

${\rm R}(\omega)\doteq{\rm G}^{-1}(\omega)-\frac{1}{\omega}$ (17)

where ${\rm G}^{-1}$ denotes the functional inverse (w.r.t. decomposition) of
the function

${\rm G}(z)\doteq\lim_{N\to\infty}\frac{1}{N}{\rm
tr}({\mathchoice{\mbox{\boldmath$\displaystyle K$}}{\mbox{\boldmath$\textstyle
K$}}{\mbox{\boldmath$\scriptstyle K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}}(z{\mathchoice{\mbox{\boldmath$\displaystyle
K$}}{\mbox{\boldmath$\textstyle K$}}{\mbox{\boldmath$\scriptstyle
K$}}{\mbox{\boldmath$\scriptscriptstyle K$}}}-\bf I)^{-1}).$ (18)

Note that, when $\textstyle K$ has an inverse, the function ${\rm R}$ stands
for the R-transform [32] of the limiting spectral distribution of
${\mathchoice{\mbox{\boldmath$\displaystyle K$}}{\mbox{\boldmath$\textstyle
K$}}{\mbox{\boldmath$\scriptstyle K$}}{\mbox{\boldmath$\scriptscriptstyle
K$}}}^{-1}$. In general, ${\rm R}(\omega)$ is well-defined and it is related
to the limiting distribution of the non-zero eigenvalues of $\textstyle K$.
Finally, the random field $\tilde{\gamma}^{(t)}$ stands for the effective
stochastic process of an arbitrary component of
${\mathchoice{\mbox{\boldmath$\displaystyle\tilde{\gamma}$}}{\mbox{\boldmath$\textstyle\tilde{\gamma}$}}{\mbox{\boldmath$\scriptstyle\tilde{\gamma}$}}{\mbox{\boldmath$\scriptscriptstyle\tilde{\gamma}$}}}^{(t)}$
in (8). Specifically, we have

$\tilde{\gamma}^{(t)}=\frac{m_{\nu^{(t-1)}}(\gamma^{(t-1)},y)}{\chi^{(t)}}-\gamma^{(t-1)}$
(19)

where for convenience we have replaced the empirical averages in representing
the dynamical order parameter in the algorithm, such as $\chi^{(t)}$,
$\lambda^{(t)}$ and etc, by the averages w.r.t. the effective stochastic
process (13), e.g.
$\chi^{(t)}=\mathbb{E}[m^{\prime}_{\nu^{(t-1)}}(\gamma^{(t-1)},y)]$. We can
see that the explicit computations of order parameter functions require
expectations of nonlinear functions of pairs of correlated Gaussian random
variables. These can be performed easily by numerical quadrature. To obtain a
recursion for such order parameters, we note that the first line of (13)
implies

$\mathcal{C}_{\phi}^{(t,t^{\prime})}=\mathcal{C}^{(t,t^{\prime})}+q\hat{\mathcal{B}}^{(t)}\hat{\mathcal{B}}^{(t^{\prime})}$
(20)

where
$\mathcal{C}_{\phi}^{(t,t^{\prime})}\doteq\mathbb{E}[\phi^{(t)}\phi^{(t^{\prime})}]$.
Finally, the covariance of the $\gamma^{(t)}$ variables can be obtained from
the second line of (13) by averaging over the decision variables $p^{(t)}$

$\displaystyle\mathcal{C}_{\gamma}^{(t,t^{\prime})}=(1-\eta)^{2}\mathcal{C}_{\gamma}^{(t-1,t^{\prime}-1)}+\eta^{2}\left[\mathcal{C}_{\phi}^{(t,t^{\prime})}\right.$
$\displaystyle\left.+\sum_{l^{\prime}=1}^{t^{\prime}-1}(1-\eta)^{t^{\prime}-l^{\prime}}\mathcal{C}_{\phi}^{(t,l^{\prime})}+\sum_{l=1}^{t-1}(1-\eta)^{t-l}\mathcal{C}_{\phi}^{(t^{\prime},l)}\right]\;.$
(21)

Combined with (15) and (19), we obtain a closed set of equations for the
iterative computation of two time correlation functions. We give explicit
results of such computations together with comparisons to simulations of the
algorithm for the perceptron model (3) in section 5. In the following section,
we will analyse the local convergence properties of the algorithm based on a
recursion for the necessary single time order parameters.

###  4.1 The fixed point solution

We assume in the following that parameters of the probabilistic model and
initial conditions are chosen in such a way that asymptotically for large
times, the algorithm will converge to a fixed point. We will then analyze the
consistency of this assumption and establish a necessary criterion for
convergence and show its relation to the AT line of the static learning model.
Translated to the case of the single node dynamics of the dynamical mean field
approach, we assume that $\gamma^{(t)}$ converges to a static random variable
$\gamma^{\star}$ for $t\to\infty$. It follows from (13), that
$\gamma^{\star}=\phi^{\star}$ which is the limit of the Gaussian random
variable $\phi^{(t)}$ which drives the dynamics. Specifically, we have

$\left(\theta,y,\gamma^{\star}\right)\sim\mathcal{N}(\theta|0,q)p_{0}(y|\theta)\mathcal{N}(\gamma^{\star}|\theta\mathcal{\hat{B}}^{\star},\mathcal{C}^{\star})$
(22)

where $\mathcal{\hat{B}}^{\star}$ stands for the stationary solution of
$\hat{\mathcal{B}}^{(t)}$ etc. It is easy to see that

$\displaystyle\hat{\mathcal{B}}^{\star}=\zeta^{\star}\mathbb{E}[\theta
m_{\nu^{\star}}(\gamma^{\star},y)].$ (23)

Then, it follows from (15) that

$\mathcal{C}^{\star}=\frac{(\hat{\mathcal{B}}^{\star})^{2}}{\zeta^{\star}(q-\chi^{\star})}+\left(\mathbb{E}[m_{\nu^{\star}}(\gamma^{\star},y)^{2}]-\frac{(\hat{\mathcal{B}}^{\star})^{2}}{(\zeta^{\star})^{2}(q-\chi^{\star})}\right){\rm
R}^{\prime}(-\chi^{\star}).$ (24)

For example, in the teacher–student matching case, i.e. when
$p_{0}(y|\theta)=p(y|\theta)$, we have

$\mathbb{E}[m_{\nu^{\star}}(\gamma^{\star},y)^{2}]=\mathbb{E}[\theta
m_{\nu^{\star}}(\gamma^{\star},y)]=q-\chi^{\star}$ (25)

and the general solutions (23) and (24) simplify to
$\mathcal{\hat{B}}^{\star}=\mathcal{C}^{\star}={\rm R}(-\chi^{\star})-1/q$.
which agrees with our previous results [21].

###  4.2 Single-time recursion of dynamics and the AT instability criteria

In order to study the convergence towards the fixed point over time, we need
the covariances between the random variables $\gamma^{(t)}$, $\phi^{(t)}$ and
their asymptotic limits. This can be obtained from recursions of the single
time order parameters defined from the limits
$\mathcal{C}_{\phi,\gamma}^{(t)}\doteq\lim_{t^{\prime}\to\infty}\mathcal{C}_{\gamma,\phi}^{(t,t^{\prime})}$
(assuming the limits exist). Using (20), (21) and (15), we get

$\displaystyle\mathcal{C}_{\gamma}^{(t)}$
$\displaystyle=(1-\eta)\mathcal{C}_{\gamma}^{(t-1)}+\eta\mathcal{C}_{\phi}^{(t)}$
(26) $\displaystyle\mathcal{C}_{\phi}^{(t)}$
$\displaystyle=\mathcal{C}^{(t)}+q\hat{\mathcal{B}}^{(t)}\hat{\mathcal{B}}^{\star}$
(27) $\displaystyle\mathcal{C}^{(t)}$
$\displaystyle=\frac{\mathcal{D}^{(t)}+\mathcal{Q}^{(t)}\left(q\mathcal{\hat{B}}^{(t)}\mathcal{\hat{B}}^{\star}+\mathcal{\hat{B}}^{(t)}\mathbb{E}[\theta\tilde{\gamma}^{\star}]+\mathcal{\hat{B}}^{\star}\mathbb{E}[\theta\tilde{\gamma}^{(t)}]+\mathbb{E}[\tilde{\gamma}^{(t)}\tilde{\gamma}^{\star}]\right)}{1-\mathcal{Q}^{(t)}}$
(28)

where e.g.
$\mathcal{Q}^{(t)}\doteq\lim_{t^{\prime}\to\infty}\mathcal{Q}^{(t,t^{\prime})}$.
In a similar way, we can show that

$\displaystyle\mathbb{E}[\theta\tilde{\gamma}^{(t)}]$
$\displaystyle=(1-\eta)\mathbb{E}[\theta\tilde{\gamma}^{(t-1)}]+\eta\mathbb{E}[\theta
f_{t}(\phi^{(t-1)})]$ (29)
$\displaystyle\mathbb{E}[\tilde{\gamma}^{(t)}\tilde{\gamma}^{\star}]$
$\displaystyle=(1-\eta)\mathbb{E}[\tilde{\gamma}^{(t-1)}\tilde{\gamma}^{\star}]+\eta\mathbb{E}[f_{t}(\phi^{(t-1)})f_{\star}(\phi^{\star})]$
(30)

where we have introduced the function

$f_{t,\star}(x,y)\doteq\frac{m_{\nu^{(t-1)},\nu^{\star}}(x,y)}{\mathbb{E}[m^{\prime}_{\nu^{(t-1)},\nu^{\star}}(x,y)]}-x.$
(31)

Based on these recursions, we will derive a condition on the parameters of the
model for which the assumption of convergence leads to a contradiction. To
this end, we study the asymptotic speed of convergence
$\gamma^{(t)}\to\gamma^{\star}$ which we define as

$\mu_{\gamma}\doteq\lim_{t\to\infty}\frac{~{}\mathcal{C}_{\gamma}^{\star}-\mathcal{C}_{\gamma}^{(t)}~{}}{\mathcal{C}_{\gamma}^{\star}-\mathcal{C}_{\gamma}^{(t-1)}}\;.$
(32)

The condition $\mu_{\gamma}\geq 1$ implies that the algorithm no longer
converges. Note, that a divergence of the algorithm was not observed for the
teacher–student matching scenario discussed in [21]. In B, we derive the
explicit formula

$\mu_{\gamma}=1-\eta\frac{1-\mathbb{E}[(m^{\prime}_{\nu^{\star}}(\gamma^{\star},y))^{2}]{\rm
R}^{\prime}(-\chi^{\star})}{1-(\chi^{\star})^{2}{\rm
R}^{\prime}(-\chi^{\star})}.$ (33)

Here, from definition of the function ${\rm R}(\omega)$ (17) it follows that
the term $1-(\chi^{\star})^{2}{\rm R}^{\prime}(-\chi^{\star})$ is always
positive. Hence, $\mu_{\gamma}\geq 1$ if and only if

$\mathbb{E}[(m^{\prime}_{\nu^{\star}}(\gamma^{\star},y))^{2}]{\rm
R}^{\prime}(-\chi^{\star})\geq 1.$ (34)

Following the arguments in [14, 15] we conclude that equation (34) coincides
with the stability condition of the replica symmetric ansatz for the static
probabilistic model- known as the de Almeida Thouless (AT) criterion.
Remarkably, the stability criterion is independent of the update parameter
$\eta$ . This indicates that (at least within our theoretical setting), a
diverging parallel iterated algorithm cannot be made convergent by reverting
to a random sequential version.

##  5 Simulation results

In the following, we compare our analytical results to numerical simulations
of the algorithm for the perceptron model. We assume that the teacher model
from which data are generated is of the general form (3) with teacher
parameters denoted by $\beta$ and $\sigma$. For the student likelihood used in
the inference algorithm we restrict ourselves to the simple noise free
likelihood model

$p(y|\theta)=\Theta(y\theta)$ (35)

where $\Theta$ stands for the unit-step function. This is a special case of
(3) corresponding to the limits $\beta,\sigma\to 0$. We specialise on the
following random matrix models for $\textstyle X$: (i) the entries of
$\textstyle X$ are independent Gaussian with zero mean variance $1/N$; (ii)
${\mathchoice{\mbox{\boldmath$\displaystyle X$}}{\mbox{\boldmath$\textstyle
X$}}{\mbox{\boldmath$\scriptstyle X$}}{\mbox{\boldmath$\scriptscriptstyle
X$}}}={{\mathchoice{\mbox{\boldmath$\displaystyle
O$}}{\mbox{\boldmath$\textstyle O$}}{\mbox{\boldmath$\scriptstyle
O$}}{\mbox{\boldmath$\scriptscriptstyle
O$}}}}{\mathchoice{\mbox{\boldmath$\displaystyle
S$}}{\mbox{\boldmath$\textstyle S$}}{\mbox{\boldmath$\scriptstyle
S$}}{\mbox{\boldmath$\scriptscriptstyle S$}}}$ where $\textstyle S$ is the
$N\times P$ projection matrix with $N\geq P$ and
${\mathchoice{\mbox{\boldmath$\displaystyle S$}}{\mbox{\boldmath$\textstyle
S$}}{\mbox{\boldmath$\scriptstyle S$}}{\mbox{\boldmath$\scriptscriptstyle
S$}}}_{ij}=\delta_{ij}$ for all $i,j$, and $\textstyle O$ is an $N\times N$
Haar random matrix. The function ${\rm R}(\omega)$ in (17) for these models
reads as

${\rm
R}(\omega)=\left\\{\begin{array}[]{cc}\frac{q-1-\sqrt{(q-1)^{2}-4\omega}}{2\omega}&\text{model}~{}(i)\\\
1+\frac{q-1}{\omega}&~{}\text{model}~{}(ii)\end{array}\right.\;.$ (36)

We will first present non asymptotic (finite times) results. In order to
demonstrate that our analytical approach also applies to non convergent
dynamics of the algorithm, we consider model parameter settings from the
_unstable_ region. Specifically, we set $\beta_{0}=\frac{1}{2}$ and $\alpha=2$
and $\eta=0.8$. In this case, we obtain the following results for two–time
correlations over 5 time steps with

$\displaystyle\mathbb{E}[\phi^{(1:5)}(\phi^{(1:5)})^{\top}]$
$\displaystyle=\tiny{\left[\begin{array}[]{ccccc}6.47&11.73&17.45&23.55&29.98\\\
11.73&24.03&36.60&50.15&64.46\\\ 17.45&36.60&58.6041&80.29&103.80\\\
23.55&50.15&80.2957&114.26&146.86\\\
29.98&64.46&103.80&146.86&194.73\end{array}\right]}$ (42)
$\displaystyle\frac{1}{N}{\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(1:5)}({\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(1:5)})^{\top}$
$\displaystyle={\tiny\left[\begin{array}[]{ccccc}6.48&11.74&17.38&23.36&29.62\\\
11.74&24.09&36.50&49.82&63.80\\\ 17.38&36.50&58.15&79.39&102.27\\\
23.36&49.82&79.39&112.72&144.46\\\
29.62&63.80&102.27&144.46&190.87\end{array}\right]}.$ (48)

We clearly see a strong increase in autocorrelations over time, indicating the
divergence of the algorithm. Here, the simulation result is based on a single
instance of the model with $N=2^{14}$. The results were obtained from the
random matrix (i). For the random matrix model (ii) we have similar theory-
experiment agreement.

Secondly, we present results on the error of estimating the true teacher
parameter $\textstyle\theta$ at each iteration step in Figure 1. The
parameters correspond to the region of convergence. In contrast to typical
results for cases of teacher–student model matching (with optimally chosen
variance of initial conditions), the prediction error turns out to be non-
monotonic.

Figure 1: Predicting of the estimation error for a given iteration time-step
$t$. E.g. the figure with label (i) is for the random matrix model (i). The
model parameters are chosen as $\sigma_{0}^{2}=10^{-2}$, $\beta_{0}=0.2$,
$\eta=0.5$, $N=3P/2$ and $P=2^{12}$.

Finally, we illustrate the asymptotic speed of convergence predicted by the
theory compared to a single simulation of the algorithm. To show the
robustness of our results, we have chosen the parameters yielding large values
of static order parameters. Nevertheless, we find a remarkably good prediction
of the exponential convergence.

Figure 2: Asymptotic of the algorithm. The model parameters are chosen as
$\sigma_{0}^{2}=0.01$, $\beta_{0}=0.3$, $\eta=0.5$, $N=3P$ and $P=2^{12}$. In
this case, we have $\mathcal{C}_{\phi}^{\star}\approx 8\times 10^{15}$ and
$\hat{\mathcal{B}}^{\star}\approx 10^{7}$.

##  6 Summary and discussion

We have analysed the dynamics of a message passing algorithm for inference in
large latent Gaussian variable models. Our analysis is based on a teacher-
student scenario together with random matrix assumptions for data. We have
focused on the problem of student-teacher mismatch and random sequential
updates. Using a dynamical functional approach we have decoupled the degrees
of freedom and have derived an effective stochastic dynamics for single nodes.
The absence of memory terms in the single node dynamics leads to tractable
recursions for two–time correlation functions. Comparison between our theory
and simulations on single instances of large systems show excellent agreement.

We have shown that a teacher-student mismatch opens the possibility of a
divergence of the algorithm. We have identified the range of model parameters
for which convergence to a fixed point is impossible. Our main result is that
the critical set of parameters is identified as the AT line of the static
replica symmetric solution. It would be interesting to see if one could prove
global convergence in the stable region. For this one would have to go beyond
the local stability analysis presented in this paper and study the full
temporal development of a set of coupled order parameters. A possible
simplification could be the construction of a Lyapunov function for the single
node dynamics.

Since the static AT stability criterion is independent of the update schedule
of the dynamics, we were not able to show that a divergent (parallel)
algorithm can be made convergent using random sequential updates. One might
argue that this negative result could be related to the random matrix
distributions used in the modeling of the data. A second possibility is the
simplicity of the node update used in our model. As an alternative one could
define random updates of the
${\mathchoice{\mbox{\boldmath$\displaystyle\tilde{\gamma}$}}{\mbox{\boldmath$\textstyle\tilde{\gamma}$}}{\mbox{\boldmath$\scriptstyle\tilde{\gamma}$}}{\mbox{\boldmath$\scriptscriptstyle\tilde{\gamma}$}}}^{(t)}$
in last line of (8) variables instead. The dynamical mean field analysis of
the corresponding model will be given elsewhere.

## Acknowledgment

This work was supported by the German Research Foundation, Deutsche
Forschungsgemeinschaft (DFG), under Grant “RAMABIM” with No. OP 45/9-1, by the
US National Science Foundation under grant CCF-1910410, and by the Harvard FAS
Dean’s Competitive Fund for Promising Scholarship.

##  Appendix A The dynamical functional analysis

The disorder average in (12) can be computed using the saddle-point method.
Specifically, we can follow the steps [21, Eq. (B.6)–(B.34)], by essentially
replacing all averages over the matrix $\textstyle A$ by averages over
$\textstyle J$ and read off the result. In our case, the variables
$\textstyle\theta$, ${\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle m$}}}^{(k)}$ and
${\mathchoice{\mbox{\boldmath$\displaystyle\hat{\psi}$}}{\mbox{\boldmath$\textstyle\hat{\psi}$}}{\mbox{\boldmath$\scriptstyle\hat{\psi}$}}{\mbox{\boldmath$\scriptscriptstyle\hat{\psi}$}}}^{(k)}$
play the roles of $\textstyle u$,
${\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}(k)$
and
${\mathchoice{\mbox{\boldmath$\displaystyle\hat{\rho}$}}{\mbox{\boldmath$\textstyle\hat{\rho}$}}{\mbox{\boldmath$\scriptstyle\hat{\rho}$}}{\mbox{\boldmath$\scriptscriptstyle\hat{\rho}$}}}(k)$
in [21], respectively. Doing so leads to the large $N$ limit approximation of
the the averaged generating functional as

$\displaystyle\mathbb{E}[Z_{i}(\\{l(t)\\})]\simeq\int{\rm d}\theta{\rm d}y{\rm
d}\phi^{(1:T)}\;p_{0}(y|\theta)\mathcal{N}(\theta|0,q)\mathcal{N}(\phi^{(1:T)}|\theta\mathcal{\hat{B}},\mathcal{C})\prod_{t\leq
K}{\rm dP}(p^{(t)}){\rm d}\psi^{(t)}{\rm d}m^{(t)}$
$\displaystyle\quad\times\delta\left[m^{(t)}-f_{t}\left(\\{\psi^{(l)},m^{(l)},p^{(l)}\\}_{l=1}^{t};y\right)\right]\delta\left[\psi^{(t)}-\phi^{(t)}-\sum_{t\leq
T}\mathcal{\hat{G}}^{(t,t^{\prime})}m^{(t)}\right]{\rm e}^{{\rm
i}\psi^{(t)}l^{(t)}}.$ (49)

Here, $\hat{\mathcal{G}}^{(t,l)}$ denotes the $(t,l)$th indexed entries of the
$T\times T$ memory matrix $\hat{\mathcal{G}}$ which is defined in terms of the
R-transform and its power series expansion as

$\hat{\mathcal{G}}={\rm
R}(\mathcal{G})=\sum_{n=1}^{\infty}c_{n}\mathcal{G}^{n-1}.$ (50)

Te entries of the $T\times T$ response matrix $\mathcal{G}$ are given by

$\mathcal{G}^{(t,t^{\prime})}\doteq\mathbb{E}\left[\frac{\partial
m^{(t)}}{\partial\phi^{(t^{\prime})}}\right].$ (51)

Moreover, the Gaussian process $\\{\phi^{(t)}\\}$ has the $T\times 1$ mean
vector $\theta\mathcal{\hat{B}}$ and $T\times T$ covariance matrix
$\mathcal{C}$ which are computed by

$\displaystyle\mathcal{\hat{B}}$
$\displaystyle=\left(\sum_{n=2}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}(-q)^{n^{\prime}}\mathcal{G}^{n-n^{\prime}-2}\right)\mathcal{B}$
(52) $\displaystyle\mathcal{C}$
$\displaystyle=\sum_{n=2}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}\mathcal{G}^{n^{\prime}}\mathcal{C}_{m}(\mathcal{G}^{\top})^{n-2-n^{\prime}}$
$\displaystyle-\sum_{n=3}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}(-q)^{n^{\prime}}\sum_{l=0}^{n-n^{\prime}-3}(\mathcal{G}^{\top})^{l}\mathcal{B}\mathcal{B}^{\top}\mathcal{G}^{n-n^{\prime}-l-3}$
(53)

where we have defined

$\displaystyle\mathcal{C}_{m}^{(t,t^{\prime})}\doteq\mathbb{E}[m^{(t)}m^{(t^{\prime})}]~{}~{}\text{and}~{}~{}\mathcal{B}^{(t)}\doteq-\mathbb{E}[\theta
m^{(t)}].$ (54)

###  A.1 The analysis of sequential dynamics

By using the property of Dirac-delta function
$\delta({\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle
y$}}})=|{\mathchoice{\mbox{\boldmath$\displaystyle
X$}}{\mbox{\boldmath$\textstyle X$}}{\mbox{\boldmath$\scriptstyle
X$}}{\mbox{\boldmath$\scriptscriptstyle
X$}}}|\delta({\mathchoice{\mbox{\boldmath$\displaystyle
X$}}{\mbox{\boldmath$\textstyle X$}}{\mbox{\boldmath$\scriptstyle
X$}}{\mbox{\boldmath$\scriptscriptstyle
X$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
y$}}{\mbox{\boldmath$\textstyle y$}}{\mbox{\boldmath$\scriptstyle
y$}}{\mbox{\boldmath$\scriptscriptstyle y$}}})$ we note that

$\displaystyle\delta\left[{\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}-\left(\frac{1}{\tau^{(t)}}(\lambda^{(t)}{\bf
I}+{\mathchoice{\mbox{\boldmath$\displaystyle J$}}{\mbox{\boldmath$\textstyle
J$}}{\mbox{\boldmath$\scriptstyle J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}})^{-1}-{\bf
I}\right)\tilde{{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}}^{(t)}\right]=\frac{1}{c^{(t)}}\delta\left[\tilde{{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}}^{(t)}-\tau^{(t)}(\lambda^{(t)}{\bf
I}+{\mathchoice{\mbox{\boldmath$\displaystyle J$}}{\mbox{\boldmath$\textstyle
J$}}{\mbox{\boldmath$\scriptstyle J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}})({\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}+\tilde{{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}}^{(t)})\right]$
(55) $\displaystyle=\frac{1}{c^{(t)}}\int{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle m$}}{\mbox{\boldmath$\textstyle
m$}}{\mbox{\boldmath$\scriptstyle m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}{\rm
d}{\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(t)}\delta[\tilde{{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}}^{(t)}+\lambda^{(t)}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}+{\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(t)}]\delta[{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}+\tau^{(t)}({\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}+\tilde{{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}}^{(t)})]\delta[{\mathchoice{\mbox{\boldmath$\displaystyle\psi$}}{\mbox{\boldmath$\textstyle\psi$}}{\mbox{\boldmath$\scriptstyle\psi$}}{\mbox{\boldmath$\scriptscriptstyle\psi$}}}^{(t)}-{\mathchoice{\mbox{\boldmath$\displaystyle
J$}}{\mbox{\boldmath$\textstyle J$}}{\mbox{\boldmath$\scriptstyle
J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle m$}}}^{(t)}]$ (56)

where for short we define the dynamical determinant
$c^{(t)}\doteq\left|\tau^{(t)}(\lambda^{(t)}{\bf
I}+{\mathchoice{\mbox{\boldmath$\displaystyle J$}}{\mbox{\boldmath$\textstyle
J$}}{\mbox{\boldmath$\scriptstyle J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}})\right|$ which do not depend on the disorder variables and thereby they
solely play the role of appropriate constant terms in the disorder average.
Indeed, one can express the dynamics
${\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}$
in (8b) in terms of the system of equations

$\displaystyle{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle m$}}}^{(t)}$
$\displaystyle=-\tau^{(t)}({\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}+\tilde{{\mathchoice{\mbox{\boldmath$\displaystyle\gamma$}}{\mbox{\boldmath$\textstyle\gamma$}}{\mbox{\boldmath$\scriptstyle\gamma$}}{\mbox{\boldmath$\scriptscriptstyle\gamma$}}}}^{(t)})$
(57a)
$\displaystyle{\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}$
$\displaystyle={\mathchoice{\mbox{\boldmath$\displaystyle
J$}}{\mbox{\boldmath$\textstyle J$}}{\mbox{\boldmath$\scriptstyle
J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}-\nu^{(t)}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle m$}}}^{(t)}.$ (57b)

Consequently, by the general results of the dynamical functional theory (49),
$\\{\phi_{i}^{(t)}\\}_{t=1}^{T}$ (for an arbitrary component $i$) can be
transformed into a Gaussian random sequence by appropriate subtractions. The
subtractions define an auxiliary dynamical system which is obtained by
replacing the variable
${\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}$
by

${\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}_{aux}^{(t)}={\mathchoice{\mbox{\boldmath$\displaystyle
J$}}{\mbox{\boldmath$\textstyle J$}}{\mbox{\boldmath$\scriptstyle
J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle m$}}}^{(t)}-\sum_{l\leq
t}\hat{\mathcal{G}}^{(t,l)}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle m$}}}^{(t)}$ (58)

for $t=1,2,\ldots T$. The entries of the response matrix $\mathcal{G}$ read

$\mathcal{G}^{(t,t^{\prime})}\doteq\mathbb{E}\left[\frac{\partial
m^{(t)}}{\partial\phi_{aux}^{(t^{\prime})}}\right].$ (59)

Moreover, by construction we have

$\frac{\partial{\gamma}^{(t-1)}}{\partial{\phi}_{aux}^{(t^{\prime})}}=\underbrace{p^{(t^{\prime})}\prod_{l=t^{\prime}+1}^{t-1}(1-p^{(l)})}_{\doteq
p^{(t,t^{\prime})}}\quad~{}t^{\prime}<t.$ (60)

Hence, the response terms read

$\displaystyle\mathcal{G}^{(t,t^{\prime})}$
$\displaystyle=-\tau^{(t)}\delta_{tt^{\prime}}-\tau^{(t)}\mathbb{E}\left[f_{t}(\gamma^{(t-1)};y)p^{(t,t^{\prime})}\right]$
$\displaystyle=-\tau^{(t)}\delta_{tt^{\prime}}-\tau^{(t)}{\rm
Pr}(p^{(t,t^{\prime})}=1)\mathbb{E}\left[f_{t}^{\prime}(\gamma^{(t-1)})p^{(t,t^{\prime})}|p^{(t,t^{\prime})}=1\right]$
$\displaystyle=-\tau^{(t)}\delta_{tt^{\prime}}-\tau^{(t)}\eta(1-\eta)^{t-1-t^{\prime}}\mathbb{E}\left[f_{t}^{\prime}(\phi_{aux}^{(t^{\prime})})\right]$
$\displaystyle=-\tau^{(t)}\delta_{tt^{\prime}}$ (61)

where for convenience we have introduced the function

$f_{t}(x,y)\doteq\frac{m_{\nu^{(t-1)}}(x,y)}{\mathbb{E}[m^{\prime}_{\nu^{(t-1)}}(x,y)]}-x$
(62)

which fulfills the divergence-free property
$\mathbb{E}[f^{\prime}_{t}(x,y)]=0$. Thereby, we get

$\displaystyle{\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}_{aux}^{(t)}$
$\displaystyle={\mathchoice{\mbox{\boldmath$\displaystyle
J$}}{\mbox{\boldmath$\textstyle J$}}{\mbox{\boldmath$\scriptstyle
J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}-\mathcal{\hat{G}}^{(t,t)}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle m$}}}^{(t)}$ (63)
$\displaystyle={\mathchoice{\mbox{\boldmath$\displaystyle
J$}}{\mbox{\boldmath$\textstyle J$}}{\mbox{\boldmath$\scriptstyle
J$}}{\mbox{\boldmath$\scriptscriptstyle
J$}}}{\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle m$}}}^{(t)}-{\rm
R}(-\tau^{(t)}){\mathchoice{\mbox{\boldmath$\displaystyle
m$}}{\mbox{\boldmath$\textstyle m$}}{\mbox{\boldmath$\scriptstyle
m$}}{\mbox{\boldmath$\scriptscriptstyle
m$}}}^{(t)}={\mathchoice{\mbox{\boldmath$\displaystyle\phi$}}{\mbox{\boldmath$\textstyle\phi$}}{\mbox{\boldmath$\scriptstyle\phi$}}{\mbox{\boldmath$\scriptscriptstyle\phi$}}}^{(t)}.$
(64)

The effective stochastic process (w.r.t. dynamical functional analysis) of
$\phi_{i}^{(1:T)}$ becomes then a Gaussian process as

$\phi^{(1:K)}\sim\mathcal{N}(\theta\mathcal{\hat{B}},\mathcal{C})$ (65)

We next use the result
$\mathcal{G}^{(t,t^{\prime})}=-\tau^{(t)}\delta_{tt^{\prime}}$ to compute the
the necessary order parameters $\mathcal{\hat{B}}$ and $\mathcal{C}$ in (52)
and (53), respectively.

####  A.1.1 Computation of $\mathcal{\hat{B}}$

We have

$\displaystyle\mathcal{\hat{B}}^{(t)}$
$\displaystyle=\left(\sum_{n=2}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}(-q)^{n^{\prime}}(-\tau^{(t)})^{n-n^{\prime}-2}\right)\mathcal{B}^{(t)}$
(66)
$\displaystyle=\frac{\mathcal{B}^{(t)}}{\tau^{(t)}-q}\sum_{n=2}^{\infty}c_{n}[(-q)^{n-1}-(-\tau^{(t)})^{n-1}]$
(67)
$\displaystyle=\frac{\mathcal{B}^{(t)}}{\tau^{(t)}-q}\sum_{n=1}^{\infty}c_{n}[(-q)^{n-1}-(-\tau^{(t)})^{n-1}]$
(68) $\displaystyle=\frac{\mathcal{B}^{(t)}}{\tau^{(t)}-q}({\rm R}(-q)-{\rm
R}(-\tau^{(t)}))$ (69)
$\displaystyle=\frac{\mathcal{B}^{(t)}}{\tau^{(t)}-q}(1/q-{\rm
R}(-\tau^{(t)}))$ (70)

On the other hand, we have

$\displaystyle\mathcal{B}^{(t)}$
$\displaystyle=\tau^{(t)}\mathbb{E}[\theta(\phi^{(t)}+\tilde{\gamma}^{(t)})]$
(71)
$\displaystyle=\tau^{(t)}(q\mathcal{\hat{B}}^{(t)}+\mathbb{E}[\theta\tilde{\gamma}^{(t)}]).$
(72)

Combining both result we easily obtain that

$\mathcal{\hat{B}}^{(t)}=\frac{\tau^{(t)}\zeta^{(t)}\mathbb{E}[\theta\tilde{\gamma}^{(t)}]}{1-q\tau^{(t)}\zeta^{(t)}}\quad\text{with}\quad\zeta^{(t)}\doteq\frac{{\rm
R}(-\tau^{(t)})-1/q}{q-\tau^{(t)}}.$ (73)

####  A.1.2 Computation of $\mathcal{C}$

Recall that

$\displaystyle\mathcal{C}$
$\displaystyle=\sum_{n=2}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}\mathcal{G}^{n^{\prime}}\mathcal{C}_{m}(\mathcal{G}^{\top})^{n-2-n^{\prime}}$
$\displaystyle\underbrace{-\sum_{n=3}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}(-q)^{n^{\prime}}\sum_{l=0}^{n-n^{\prime}-3}(\mathcal{G}^{\top})^{l}\mathcal{B}\mathcal{B}^{\top}\mathcal{G}^{n-n^{\prime}-l-3}}_{\mathcal{D}}.$
(74)

We then write

$\displaystyle\mathcal{C}^{(t,t^{\prime})}-\mathcal{D}^{(t,t^{\prime})}$
$\displaystyle=\mathcal{C}_{m}^{(t,t^{\prime})}\sum_{n=2}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}(-\tau^{(t)})^{n^{\prime}}(-\tau^{(t^{\prime})})^{n-2-n^{\prime}}$
(75)
$\displaystyle=\mathcal{C}_{m}^{(t,t^{\prime})}\sum_{n=2}^{\infty}c_{n}\frac{(-\tau^{(t^{\prime})})^{n-1}-(-\tau^{(t)})^{n-1}}{\tau^{(t)}-\tau^{(t^{\prime})}}$
(76) $\displaystyle=\mathcal{C}_{m}^{(t,t^{\prime})}\frac{{\rm
R}(-\tau^{(t^{\prime})})-{\rm
R}(-\tau^{(t)})}{\tau^{(t)}-\tau^{(t^{\prime})}}$ (77)

On the other hand, by construction we have

$\displaystyle\mathcal{C}_{m}^{(t,t^{\prime})}$
$\displaystyle=\tau^{(t)}\tau^{(t^{\prime})}(\mathbb{E}[\phi^{(t)}\phi^{(t^{\prime})}]+\mathbb{E}[\phi^{(t)}\tilde{\gamma}^{(t^{\prime})}]+\mathbb{E}[\phi^{(t^{\prime})}\tilde{\gamma}^{(t)}]+\mathbb{E}[\tilde{\gamma}^{(t)}\tilde{\gamma}^{(t^{\prime})}])$
(78)
$\displaystyle=\tau^{(t)}\tau^{(t^{\prime})}(\mathcal{C}^{(t,t^{\prime})}+q\mathcal{\hat{B}}^{(t)}\mathcal{\hat{B}}^{(t^{\prime})}+\mathbb{E}[\phi^{(t)}\tilde{\gamma}^{(t^{\prime})}]+\mathbb{E}[\phi^{(t^{\prime})}\tilde{\gamma}^{(t)}]+\mathbb{E}[\tilde{\gamma}^{(t)}\tilde{\gamma}^{(t^{\prime})}])$
(79)
$\displaystyle=\tau^{(t)}\tau^{(t^{\prime})}(\mathcal{C}^{(t,t^{\prime})}+q\mathcal{\hat{B}}^{(t)}\mathcal{\hat{B}}^{(t^{\prime})}+\mathcal{\hat{B}}^{(t)}\mathbb{E}[\theta\tilde{\gamma}^{(t^{\prime})}]+\mathcal{\hat{B}}^{(t^{\prime})}\mathbb{E}[\theta\tilde{\gamma}^{(t)}]+\mathbb{E}[\tilde{\gamma}^{(t)}\tilde{\gamma}^{(t^{\prime})}])$
(80)

where in the last line we have invoked the results

$\displaystyle\mathbb{E}[\phi^{(t)}\tilde{\gamma}^{(t^{\prime})}]$
$\displaystyle=\mathbb{E}[\phi^{(t)}f_{t^{\prime}}(\gamma^{(t^{\prime}-1)})]$
(81)
$\displaystyle=(1-\eta)\mathbb{E}[\phi^{(t)}f_{t^{\prime}}(\gamma^{(t^{\prime}-2)})]+\eta\mathbb{E}[\phi^{(t)}f_{t^{\prime}}(\phi^{(t^{\prime}-1)})]$
(82)
$\displaystyle=(1-\eta)\mathbb{E}[\phi^{(t)}f_{t^{\prime}}(\gamma^{(t^{\prime}-2)})]+\eta\mathcal{\hat{B}}^{(t)}\mathbb{E}[\theta
f_{t^{\prime}}(\phi^{(t^{\prime}-1)})]$ (83)
$\displaystyle=\mathcal{\hat{B}}^{(t)}\mathbb{E}[\theta\tilde{\gamma}^{(t^{\prime})}].$
(84)

The equation (83) follows from the Stein’ lemma. Then, by invoking (80) in
(77) we get

$\mathcal{C}^{(t,t^{\prime})}=\frac{\mathcal{D}^{(t,t^{\prime})}+\mathcal{Q}^{(t,t^{\prime})}\left(q\mathcal{\hat{B}}^{(t)}\mathcal{\hat{B}}^{(t^{\prime})}+\mathcal{\hat{B}}^{(t)}\mathbb{E}[\theta\tilde{\gamma}^{(t^{\prime})}]+\mathcal{\hat{B}}^{(t^{\prime})}\mathbb{E}[\theta\tilde{\gamma}^{(t)}]+\mathbb{E}[\tilde{\gamma}^{(t)}\tilde{\gamma}^{(t^{\prime})}]\right)}{1-\mathcal{Q}^{(t,t^{\prime})}}$
(85)

We complete the derivation by simplifying $\mathcal{D}^{(t,t^{\prime})}$:
Firstly, for $\tau^{(t)}\neq\tau^{(t^{\prime})}$ we have

$\displaystyle\mathcal{D}^{(t,t^{\prime})}$
$\displaystyle=-\mathcal{B}^{(t)}\mathcal{B}^{(t^{\prime})}\sum_{n=3}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}(-q)^{n^{\prime}}\sum_{l=0}^{n-n^{\prime}-3}(-\tau^{(t)})^{l}(-\tau^{(t^{\prime})})^{n-n^{\prime}-l-3}$
(86)
$\displaystyle=\frac{\mathcal{B}^{(t)}\mathcal{B}^{(t^{\prime})}}{\tau^{(t)}-\tau^{(t^{\prime})}}\sum_{n=3}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}(-q)^{n^{\prime}}\\{(-\tau^{(t)})^{n-n^{\prime}-2}-(-\tau^{(t^{\prime})})^{n-n^{\prime}-2}\\}$
(87)
$\displaystyle=\frac{\mathcal{B}^{(t)}\mathcal{B}^{(t^{\prime})}}{\tau^{(t)}-\tau^{(t^{\prime})}}\left[\frac{{\rm
R}(-\tau^{(t)})-1/q}{q-\tau^{(t)}}-\frac{{\rm
R}(-\tau^{(t^{\prime})})-1/q}{q-\tau^{(t^{\prime})}}\right]$ (88)

where we used the fact that

$\displaystyle\sum_{n=3}^{\infty}c_{n}\sum_{n^{\prime}=0}^{n-2}(-q)^{n^{\prime}}(-\tau^{(t)})^{n-n^{\prime}-2}$
$\displaystyle=\frac{1}{\tau^{(t)}-q}\sum_{n=3}^{\infty}c_{n}[(-q)^{n-1}-(-\tau^{(t)})^{n-1}]$
(89)
$\displaystyle=\frac{1}{\tau^{(t)}-q}\sum_{n=1}^{\infty}c_{n}[(-q)^{n-1}-(-\tau^{(t)})^{n-1}]-1$
(90) $\displaystyle=\frac{{\rm R}(-\tau^{(t)})-{\rm
R}(-q)}{q-\tau^{(t)}}-c_{2}$ (91) $\displaystyle=\frac{{\rm
R}(-\tau^{(t)})-1/q}{q-\tau^{(t)}}-c_{2}$ (92)

Moreover, in the equal-time case, we have

$\mathcal{D}^{(t,t)}=\frac{\mathcal{B}^{(t)}\mathcal{B}^{(t)}}{q-\tau^{(t)}}\left[\frac{{\rm
R}(-\tau^{(t)})-1/q}{q-\tau^{(t)}}-{\rm R}^{\prime}(-\tau^{(t)})\right].$ (93)

##  Appendix B Derivation of (33)

For convenience, we define the single time deviations
$\Delta_{\gamma,\phi}^{(t)}\doteq\mathcal{C}_{\gamma}^{\star}-\mathcal{C}_{\gamma,\phi}^{(t)}$.
Furthermore, from (26) we write the recursion

$\Delta_{\gamma}^{(t)}=(1-\eta)\Delta_{\gamma}^{(t-1)}+\eta\Delta_{\phi}^{(t)}.$
(94)

Moreover, we introduce the rate
$\mu_{\phi}\doteq\lim_{t\to\infty}\frac{\Delta_{\phi}^{(t)}}{\Delta_{\phi}^{(t-1)}}$.
From (94) it follows that $\mu_{\gamma}=\mu_{\phi}$. We next compute
$\mu_{\phi}$. To this end, from (27) and (28) we firstly write
$\Delta_{\phi}^{(t)}$ in the form

$\displaystyle\Delta_{\phi}^{(t)}=\mathcal{C}_{\gamma}^{\star}-c^{(t)}\hat{\mathcal{B}}^{\star}\hat{\mathcal{B}}^{(t)}-\frac{\mathcal{Q}^{(t)}}{1-\mathcal{Q}^{(t)}}\mathbb{E}[\tilde{\gamma}^{(t)}\tilde{\gamma}^{\star}]$
(95)

for an appropriately computed constant sequence $c^{(t)}$ (which does not
depend on $\Delta_{\phi}^{(t-1)}$). Then, from (30) we further write

$\displaystyle\Delta_{\phi}^{(t)}$
$\displaystyle=\mathcal{C}_{\gamma}^{\star}-c^{(t)}\hat{\mathcal{B}}^{\star}\hat{\mathcal{B}}^{(t)}-(1-\eta)\frac{\mathcal{Q}^{(t)}(1-\mathcal{Q}^{(t-1)})}{(1-\mathcal{Q}^{(t)})\mathcal{Q}^{(t-1)}}(\mathcal{C}_{\gamma}^{\star}-c^{(t-1)}\hat{\mathcal{B}}^{\star}\hat{\mathcal{B}}^{(t-1)})$
$\displaystyle+(1-\eta)\frac{\mathcal{Q}^{(t)}(1-\mathcal{Q}^{(t-1)})}{(1-\mathcal{Q}^{(t)})\mathcal{Q}^{(t-1)}}\Delta_{\phi}^{(t-1)}-\eta\frac{\mathcal{Q}^{(t)}}{1-\mathcal{Q}^{(t)}}\mathbb{E}[f^{\prime}_{t}(\phi^{(t-1)})f^{\prime}_{\star}(\phi^{\star})].$
(96)

Thereby, we get the derivative

$g_{t}(\Delta_{\phi}^{(t-1)})\doteq\frac{\partial\Delta_{\phi}^{(t)}}{\partial\Delta_{\phi}^{(t-1)}}=(1-\eta)\frac{\mathcal{Q}^{(t)}(1-\mathcal{Q}^{(t-1)})}{(1-\mathcal{Q}^{(t)})\mathcal{Q}^{(t-1)}}+\eta\frac{\mathcal{Q}^{(t)}}{1-\mathcal{Q}^{(t)}}\mathbb{E}[f^{\prime}_{t}(\phi^{(t-1)})f^{\prime}_{\star}(\phi^{\star})].$
(97)

We then obtain the rate as

$\displaystyle\mu_{\phi}$
$\displaystyle=\lim_{t\to\infty}g_{t}(0)=(1-\eta)+\eta\frac{(\chi^{\star})^{2}{\rm
R}^{\prime}(-\chi^{\star})}{1-(\chi^{\star})^{2}{\rm
R}^{\prime}(-\chi^{\star})}\mathbb{E}[(f^{\prime}_{\star}(\phi^{\star}))^{2}]$
(98) $\displaystyle=(1-\eta)+\eta\frac{{\rm
R}^{\prime}(-\chi^{\star})}{1-(\chi^{\star})^{2}{\rm
R}^{\prime}(-\chi^{\star})}(\mathbb{E}[(m_{\nu^{\star}}^{\prime}(\phi^{\star}))^{2}]-(\chi^{\star})^{2})$
(99)
$\displaystyle=1-\eta\frac{1-\mathbb{E}[(m_{\nu^{\star}}^{\prime}(\phi^{\star}))^{2}]{\rm
R}^{\prime}(-\chi^{\star})}{1-(\chi^{\star})^{2}{\rm
R}^{\prime}(-\chi^{\star})}.$ (100)

## References

  * [1] Kabashima Y 2003 Journal of Physics A: Mathematical and General 36 11111
  * [2] Bolthausen E 2014 Communications in Mathematical Physics 325 333–366.
  * [3] Bayati M and Montanari A 2011 IEEE Transactions on Information Theory 57 764–785
  * [4] Mimura K and Okada M 2014 IEEE Transactions on Information Theory 60 3645–3670
  * [5] Opper M, Çakmak B and Winther O 2016 Journal of Physics A: Mathematical and Theoretical 49 114002
  * [6] Çakmak B and Opper M 2019 Phys. Rev. E 99(6) 062140
  * [7] Çakmak B and Opper M 2020 Journal of Statistical Mechanics: Theory and Experiment 2020 103303
  * [8] Takeuchi K 2020 IEEE Transactions on Information Theory 66 368–386
  * [9] Rangan S, Schniter P and Fletcher A K 2019 IEEE Transactions on Information Theory 65 6664–6684
  * [10] Fletcher A K, Rangan S and Schniter P 2018 Inference in deep networks in high dimensions 2018 IEEE International Symposium on Information Theory (ISIT) (IEEE) pp 1884–1888
  * [11] Fan Z 2020 arXiv preprint arXiv:2008.11892
  * [12] Thousless D J, Andersen P W and Palmer R G 1977 Philosophical Magazine 35 593–601
  * [13] Parisi G and Potters M 1995 Journal of Physics A: Mathematical and General 28 5267
  * [14] Opper M and Winther O 2001 Physical Review E 64 056131–(1–14)
  * [15] Kabashima Y 2008 Journal of Physics: Conference Series 95
  * [16] Donoho D L, Maleki A and Montanari A 2009 Proceedings of the National Academy of Sciences 106 18914–18919
  * [17] Minka T P 2001 Expectation propagation for approximate Bayesian inference Proc. 17th Conference on Uncertainty in Artificial Intelligence (UAI)
  * [18] Ma J and Ping L 2017 IEEE Access 5 2020–2033
  * [19] Manoel A, Krzakala F, Tramel E and Zdeborová L 2015 Swept approximate message passing for sparse estimation Proceedings of the 32nd International Conference on Machine Learning (Proceedings of Machine Learning Research vol 37) ed Bach F and Blei D (Lille, France: PMLR) pp 1123–1132 URL https://proceedings.mlr.press/v37/manoel15.html
  * [20] Martin P C, Siggia E D and Rose H A 1973 Physical Review A 8 423
  * [21] Çakmak B and Opper M 2020 Journal of Physics A: Mathematical and Theoretical 53 274001
  * [22] Sommers H J 1987 Physical review letters 58 1268
  * [23] Sollich P and Barber D 1997 EPL (Europhysics Letters) 38 477
  * [24] Mignacco F, Krzakala F, Urbani P and Zdeborová L 2020 Dynamical mean-field theory for stochastic gradient descent in gaussian mixture classification (Preprint 2006.06098)
  * [25] Çakmak B and Opper M 2021 Phys. Rev. E 103(3) L030101 URL https://link.aps.org/doi/10.1103/PhysRevE.103.L030101
  * [26] Takahashi T and Kabashima Y 2020 Macroscopic analysis of vector approximate message passing in a model mismatch setting 2020 IEEE International Symposium on Information Theory (ISIT) (IEEE) pp 1403–1408
  * [27] Opper M and Winther O 2000 Neural computation 12 2655–2684
  * [28] Rasmussen C E and Williams C K 2006 Gaussian processes for machine learning vol 1 (Cambridge: MIT press)
  * [29] Neal R M 1997 arXiv preprint physics/9701026
  * [30] Çakmak B and Opper M 2018 Expectation propagation for approximate inference: Free probability framework 2018 IEEE International Symposium on Information Theory (ISIT) (Piscataway, NJ, USA: IEEE) pp 1276–1280 ISSN 2157-8117
  * [31] Eisfeller H and Opper M 1992 Physical Review Letters 68 2094
  * [32] Akemann G, Baik J and Di Francesco P (eds) 2011 The Oxford Handbook of Random Matrix Theory (Oxford University Press)