# Fluctuation Results for Multi-species Sherrington-Kirkpatrick model in the
Replica Symmetric Regime

Partha S. Dey⋆ and  Qiang Wu† Department of Mathematics, University of
Illinois at Urbana-Champaign, 1409 W Green Street, Urbana, Illinois 61801
<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract.

We study the Replica Symmetric region of general multi-species Sherrington-
Kirkpatrick (MSK) Model and answer some of the questions raised in Ann.
Probab. 43 (2015), no. 6, 3494–3513, where the author proved the Parisi
formula under _positive-definite_ assumption on the disorder covariance matrix
$\Delta^{2}$. First, we prove exponential overlap concentration at high
temperature for both _indefinite_ and _positive-definite_ $\Delta^{2}$ MSK
model. We also prove a central limit theorem for the free energy using overlap
concentration. Furthermore, in the zero external field case, we use a
quadratic coupling argument to prove overlap concentration up to $\beta_{c}$,
which is expected to be the critical inverse temperature. The argument holds
for both _positive-definite_ and emphindefinite $\Delta^{2}$, and $\beta_{c}$
has the same expression in two different cases. Second, we develop a species-
wise cavity approach to study the overlap fluctuation, and the asymptotic
variance-covariance matrix of overlap is obtained as the solution to a matrix-
valued linear system. The asymptotic variance also suggests the de
Almeida–Thouless (AT) line condition from the Replica Symmetry (RS) side. Our
species-wise cavity approach does not require the positive-definiteness of
$\Delta^{2}$. However, it seems that the AT line conditions in _positive-
definite_ and _indefinite_ cases are different. Finally, in the case of
_positive-definite_ $\Delta^{2}$, we prove that above the AT line, the MSK
model is in Replica Symmetry Breaking phase under some natural assumption.
This generalizes the results of J. Stat. Phys. 174 (2019), no. 2, 333–350,
from 2-species to general species.

###### Key words and phrases:

Spin glass, Phase diagram, Central limit theorem, Cavity method.

###### 2010 Mathematics Subject Classification:

Primary: 82B26, 82B44, 60F05.

##  1\. Introduction and Main Results

Multi-species Sherrington-Kirkpatrick (MSK) model, as a generalization of the
classical SK model [Tal11a, Tal11b, Pan13], was introduced by Barra  _et. al._
in [BCMT15]. In the MSK model, the spins are divided into finitely many
species, and the density of each species is asymptotically fixed with
increasing system volume. The interaction parameter between different spins
depends only on the species structure. Thus, one can consider the MSK model as
an inhomogeneous generalization of the classical SK model. The detailed
definition is given in Section 1.1.

One important question in spin glasses is to understand the behavior of
limiting free energy. In the classical SK model, Parisi [Par83, Par80]
proposed his celebrated variational formula, known as the Parisi formula,
which gives the limiting free energy at all temperatures. Talagrand [Tal06]
proved this in the setting of mixed $p$-spin models for even $p$. Later
Panchenko [Pan14], using a different approach, proved the Parisi formula of
mixed $p$-spin models for all $p$.

The authors in [BCMT15] proposed a Parisi type formula for the limiting free
energy in the MSK model. They also proved the upper bound by adapting Guerra’s
interpolation method under _positive-definite_ condition on the interaction
matrix $\Delta^{2}$. Later, Panchenko [Pan15] completed the proof by proving a
multi-species version of the Ghirlanda-Guerra identities to match the lower
bound. The _positive-definite_ condition on $\Delta^{2}$ was used only in the
proof of the upper bound. In [Pan15], Panchenko also proposed several open
questions, such as understanding behaviors of the model when the interaction
matrix $\Delta^{2}$ is _indefinite_ , deriving an analog of the AT line
condition in MSK model, among others. The authors in [BSS19] derived the AT
line condition of two-species _positive-definite_ case and proved that above
the AT line, the 2-species SK model is RSB. The essential tool of their
argument is a perturbation of the Parisi formula, which is only known for
_positive-definite_ $\Delta^{2}$. The idea fails in the general species case
because of algebraic difficulties.

For the interesting question about what happens when $\Delta^{2}$ is
_indefinite_ , to the best of our knowledge, only some partial results on a
few particular models have been known very recently. In [BGG11, BGGPT14],
there are some conjectured min-max type formulas by physicists for the
bipartite SK model, which is a special case of _indefinite_ MSK model. Another
particular case of _indefinite_ MSK model known as deep Boltzmann machine
(DBM) was investigated in [ABCM20, ACCM21] in the replica symmetry regime, and
a complete solution was obtained in [ACCM20] under assumptions on the
Nishimori line. Besides that, a min-max formula for the replica symmetric
solution in the DBM model is proved in [Gen20]. Recently, Mourrat [MD20,
Mou20A, Mou20B, Mou20C, Mou18, Mou19] has reinterpreted the Parisi formula as
the solution to an infinite-dimensional Hamilton-Jacobi equation in the
Wasserstein space of probability measures on the positive half-line.
Particularly in [Mou20B], he studied the bipartite SK model and questioned the
possibility to represent the limiting free energy as a variational formula. In
the setting of the Hamilton-Jacobi equation, non-convex $\Delta^{2}$ breaks
down the application of the Hopf-Lax formula while solving the PDE. For the
spherical case, Auffinger and Chen [AC14] studied the bipartite spherical SK
(BSSK) model and proved a variational formula for the limiting free energy at
sufficiently high temperature. Later, Baik and Lee [BL17] computed the
limiting free energy and its fluctuations at all non-critical temperatures in
the BSSK model using tools from random matrix theory. However, computing
limiting free energy for general _indefinite_ model is still unknown.
Specifically for Ising spins, we do not even know the limiting free energy at
any temperature, other than the lower bound by Panchenko [Pan15].

In this paper, we do a high-temperature analysis of the Ising MSK model, and
our main contributions are:

  1. (1)

For external field $h\geqslant 0$, by extending Latala’s argument, we prove an
RS regime in the MSK model, where the overlap has exponential concentration.
Note that our approach is unified for both _positive-definite_ and
_indefinite_ MSK model, but the proved RS regimes have a different form; see
Theorem 1.7. For $h=0$, by using a quadratic coupling argument, we prove the
concentration of overlap up to
$\beta_{c}=\rho(\Delta^{2}\Lambda)^{-\nicefrac{{1}}{{2}}}$, where $\Delta^{2}$
is the interaction matrix and $\Lambda$ is the diagonal matrix containing the
species size ratios. This result is true for _indefinite_ and _positive-
definite_ $\Delta^{2}$, _i.e.,_ we prove the whole RS regime for the general
MSK model, which is given as $\beta<\beta_{c}$, see Theorem 1.8. The above
overlap concentrations also enable us to prove the limit and central limit
theorem of free energy in the corresponding RS regime; see Theorem 1.9.

  2. (2)

By developing a different species-wise cavity method, we derive a linear
system of Gibbs average of a quadratic form of overlap vectors. The system is
solved using linear algebraic methods, enabling us to compute the variance-
covariance structure of overlap vectors. The computation also suggests the AT-
line condition in the MSK model from the replica symmetry side. Note that our
species-wise cavity method does not require _positive-definite_ $\Delta^{2}$,
but the AT line condition for _indefinite_ $\Delta^{2}$ seems to be more
complicated; see the discussions below Theorem 1.13.

  3. (3)

In the case of _positive-definite_ $\Delta^{2}$, we prove the AT line
condition from the RSB side under some natural assumption. The key is still
the perturbation idea of the Parisi formula. This generalizes the result for
$2$-species SK model in [BSS19]. For the _indefinite_ $\Delta^{2}$ case, we
use our species-wise cavity approach to give a conjectured form of the AT line
condition. To get a rigorous proof of the AT line in this case is challenging.
Because, first, the Parisi formula for _indefinite_ MSK model is still in
mystery, so the classical perturbation technique can not be used. Second,
proving the uniqueness of stationary point is hard. We prove a uniqueness
result in 2-species case using an elementary approach, see Proposition 1.6.

The definition of the MSK model is given in Section 1.1, then we review the
main results of the $2$-species case in [BSS19], and conclude Section 1 by the
statement of our results. Readers can always check notations back in Section
1.4. A road map is given in Section 1.5.

###  1.1. MSK model

Fix $N\geqslant 1$. For a spin configuration on the $N$-dimensional hyper-
cube,
$\boldsymbol{\sigma}=(\sigma_{1},\sigma_{2},\ldots,\sigma_{N})\in\Sigma_{N}:=\\{-1,+1\\}^{N}$,
we consider the Hamiltonian given by

$\displaystyle
H_{N}(\boldsymbol{\sigma}):=\frac{\beta}{\sqrt{N}}\sum_{1\leqslant
i<j\leqslant N}g_{ij}\sigma_{i}\sigma_{j}+h\sum_{i=1}^{N}\sigma_{i}$ (1)

where $g_{ij}$, the disorder interaction parameters are independent centered
Gaussian random variables, $\beta>0$ is the inverse temperature and
$h\geqslant 0$ is the external field.

In the classical SK model, the variance structure of disorder $g_{ij}$ is
homogeneous, usually taken as $g_{ij}\sim\operatorname{N}(0,1)$. While in the
MSK model, the variance of $g_{ij}$ depends on the structure of species among
$N$ spins.

Assume that, there are $m\geqslant 2$ species. When $m=1$ the MSK model
reduces to the classical SK model. We partition the set of spins into $m$
disjoint sets, namely,

$I=\bigcup_{s=1}^{m}I_{s}=\\{1,2,\ldots,N\\},\qquad I_{s}\cap
I_{t}=\emptyset\text{ for }s\neq t.$

For $i\in I_{s},j\in I_{t}$, we assume

$\operatorname{\mathds{E}}g_{ij}^{2}=\Delta_{st}^{2},$

_i.e.,_ the inhomogeneity of MSK model comes from the interaction among
different species. While proving MSK Parisi formula in [Pan15, BCMT15], the
assumption that $\Delta^{2}=(\\!(\Delta_{st}^{2})\\!)_{s,t=1}^{m}$ is
symmetric and _positive-definite_ was used. But in this paper, we do not
require the positive-definiteness condition for most of the results. Besides
that, we assume that the ratio of spins in each species is fixed
asymptotically, _i.e.,_ for $s=1,2,\ldots,m$ and

$\lambda_{s,N}:={|I_{s}|}/{N},$

we have

$\lim_{N\to\infty}\lambda_{s,N}=\lambda_{s}\in(0,1).$

Since $\lambda_{s,N}$ and $\lambda_{s}$ are asymptotically the same, for the
rest of the article we will use $\lambda_{s}$ instead of $\lambda_{s,N}$ for
convenience. We denote
$\Lambda:=\operatorname{diag}(\lambda_{1},\lambda_{2},\ldots,\lambda_{m})$.

The overlap vector between two replicas
$\boldsymbol{\sigma}^{1},\boldsymbol{\sigma}^{2}\in\Sigma_{N}$ is given by

$\mathbf{R}_{12}=(R_{12}^{(1)},R_{12}^{(2)},\ldots,R_{12}^{(m)})^{\intercal},$

where

$R_{12}^{(s)}:=\frac{1}{|I_{s}|}\sum_{i\in
I_{s}}\sigma_{i}^{1}\sigma_{i}^{2}.$

is the overlap restricted to species $s\in\\{1,\ldots,m\\}$. In some cases we
write it as $R^{(s)}$ for short if there are only 2 replicas involved. All
vectors will be considered as a column vector in the rest of the article.

A central question in spin glasses is to understand the free energy

$\displaystyle F_{N}:=\frac{1}{N}\log Z_{N},\quad\text{where}\quad
Z_{N}:=\sum_{\boldsymbol{\sigma}\in\Sigma_{N}}\exp(H_{N}(\boldsymbol{\sigma}))$
(2)

is the partition function, and the Gibbs measure is given by

$\displaystyle
G_{N}(\boldsymbol{\sigma})={\exp(H_{N}(\boldsymbol{\sigma}))}/{Z_{N}},\qquad\boldsymbol{\sigma}\in\Sigma_{N}.$
(3)

Later we will also use $\langle\cdot\rangle$ to denote the Gibbs average just
for convenience; check Section 1.4. When $\Delta^{2}$ is _positive-definite_ ,
Parisi formula gives a variational representation of the limiting free energy
for all $\beta>0,h\in\mathds{R}$, which was set up in [BCMT15, Pan15] for MSK
model. This paper is mainly about the phase transitions and fluctuation of
free energy in the RS regime. Before stating the main results, let us recall
some related results in MSK model in [BSS19, Pan15, BCMT15]. We will only use
the Parisi formula in the proof of Theorem 1.15.

###  1.2. Related results in MSK model

Consider a sequence of real numbers

$\displaystyle 0=\zeta_{0}<\zeta_{1}<\cdots<\zeta_{k}<\zeta_{k+1}=1,$ (4)

and for each $s\in\\{1,2,\ldots,m\\}$, the sequences

$\displaystyle 0=q_{0}^{s}\leqslant q_{1}^{s}\leqslant\cdots\leqslant
q_{k+1}^{s}\leqslant q_{k+2}^{s}=1.$ (5)

For $0\leqslant l\leqslant k+2$, define
$\mathbf{q}_{\ell}=(q_{\ell}^{s})_{s=1}^{m}$,

$\displaystyle
Q_{\ell}:=\frac{1}{2}{\mathbf{q}_{\ell}}^{\intercal}\Lambda\Delta^{2}\Lambda\mathbf{q}_{\ell},\quad
Q_{\ell}^{s}:=(\Delta^{2}\Lambda\mathbf{q}_{\ell})_{s}\quad\text{for}\
1\leqslant s\leqslant m.$ (6)

Given these sequences, consider i.i.d. standard Gaussian random variables
$(\eta_{\ell})_{1\leqslant\ell\leqslant k+2}$. Define

$X_{k+2}^{s}:=\log\cosh\bigl{(}h+\beta\sum_{0\leqslant\ell\leqslant
k+1}\eta_{\ell+1}(Q_{\ell+1}^{s}-Q_{\ell}^{s})^{\nicefrac{{1}}{{2}}}\bigr{)},$

then recursively define for $0\leqslant\ell\leqslant k+1$,

$X_{\ell}^{s}:=\frac{1}{\zeta_{\ell}}\log\operatorname{\mathds{E}}_{\ell+1}\exp(\zeta_{\ell}X_{\ell+1}^{s}).$

where $\operatorname{\mathds{E}}_{\ell+1}$ denotes expectation w.r.t.
$\eta_{\ell+1}$. The following theorem gives the Parisi formula in MSK model.

######  Theorem 1.1 ([Pan15]*Theorem 1).

For the MSK model with positive-definite $\Delta^{2}$, the limiting free
energy is given by

$\displaystyle\lim_{N\to\infty}F_{N}=\inf_{\boldsymbol{\zeta},\mathbf{q}}\mathscr{P}(\boldsymbol{\zeta},\mathbf{q}),$
(7)

where

$\mathscr{P}(\boldsymbol{\zeta},\mathbf{q})=\log
2+\sum_{s=1}^{m}\lambda_{s}X_{0}^{s}-\frac{\beta^{2}}{2}\sum_{\ell=1}^{k+1}\zeta_{\ell}(Q_{\ell+1}-Q_{\ell}).$

is the Parisi functional and the $\inf$ in (7) is taken over all the sequences
in (4) and (5).

In the above variational formula (7), let $k=0$ and $q_{1}^{s}=q_{s}\in[0,1]$,
then the Parisi functional on the RS regime simplifies to

$\displaystyle\mathscr{P}_{\operatorname{RS}}(\mathbf{q})$
$\displaystyle:=\log
2+\sum_{s=1}^{m}\lambda_{s}\operatorname{\mathds{E}}\log\cosh(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h)+\frac{\beta^{2}}{4}(\mathbf{1}-\mathbf{q})^{\intercal}\Lambda\Delta^{2}\Lambda(\mathbf{1}-\mathbf{q})$

where $\mathbf{1}$ is the vector of all $1$’s and
$\mathbf{q}=(q_{1},q_{2},\ldots,q_{m})^{\intercal}$. Taking derivative w.r.t.
$q_{t}$ for $t=1,\ldots,m$

$\displaystyle\frac{\partial\mathscr{P}_{\operatorname{RS}}}{\partial
q_{t}}=\frac{\beta^{2}}{2}\lambda_{t}\sum_{s=1}^{m}\Delta_{s,t}^{2}\lambda_{s}\left[q_{s}-\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h)\right].$
(8)

and setting it 0, we get the set of critical points (when $\Delta^{2}$ is
invertible) to be

$\displaystyle\mathcal{C}(\beta,h)=\left\\{\mathbf{q}\in[0,1]^{m}\ \bigl{|}\
\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h)=q_{s}\
\text{for}\ s=1,2,\ldots m\right\\}.$ (9)

For the MSK model with positive-definite $\Delta^{2}$, we define the Replica
Symmetric solution as

$\displaystyle\operatorname{RS}(\beta,h):=\inf_{\mathbf{q}\in[0,1]^{m}}\mathscr{P}_{\operatorname{RS}}(\mathbf{q}).$
(10)

Note that $\mathscr{P}_{\operatorname{RS}}$ is well-defined for _indefinite_
$\Delta^{2}$. However, one expects the replica symmetric solution to be
achieved at a saddle point of $\mathscr{P}_{\operatorname{RS}}$ instead of a
minimizer, see [BGG11]*Theorem 4. The $\mathscr{P}_{\operatorname{RS}}$
functional in the _indefinite_ case seems to be not convex. We only use the
$\operatorname{RS}(\beta,h)$ expression for proving the replica symmetry
breaking in Theorem 1.15 for _positive-definite_ case. Actually, in the RS
region for _positive-definite_ $\Delta^{2}$, the infimum should be achieved at
the critical point in $\mathcal{C}(\beta,h)$ instead of some boundary point.
This fact can easily be partially checked by comparing
$\operatorname{RS}(\beta,h)$ with the replica symmetric solution we obtained
from the overlap concentration results in Theorem 1.9. However, to rigorously
prove this for the whole region can be a nasty calculus problem. The other way
is to prove the convexity of the $\mathscr{P}_{\operatorname{RS}}$ functional,
which is beyond the scope of the current paper. Therefore, we will assume that
the infimum in (10) is not achieved at the boundary for Theorem 1.15. The
following results were proved in [BSS19] for $m=2$.

######  Theorem 1.2 ([BSS19]*Theorem 1.1–1.2).

Restricted to $m=2$ species model, under the assumption

$\displaystyle\Delta_{12}^{2}=1,\Delta_{11}^{2}\Delta_{22}^{2}>1\ \text{and}\
\lambda_{1}\Delta_{11}^{2}\geqslant\lambda_{2}\Delta_{22}^{2},$ (11)

  1. (i)

If either $h>0$ or

$\displaystyle\beta^{2}<\frac{2}{\lambda_{1}\Delta_{11}^{2}+\lambda_{2}\Delta_{22}^{2}+\sqrt{(\lambda_{1}\Delta_{11}^{2}-\lambda_{2}\Delta_{22}^{2})^{2}+4\lambda_{1}\lambda_{2}}}$
(12)

then $\mathcal{C}(\beta,h)=\\{\mathbf{q}_{\star}\\}$ is a singleton.

  2. (ii)

Assume $h>0$, let $\mathbf{q}\in\mathcal{C}(\beta,h)$, and
$\gamma_{s}:=\operatorname{\mathds{E}}\operatorname{sech}^{4}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h)$
for $s=1,2$. If

$\displaystyle\beta^{2}>\frac{2}{\lambda_{1}\gamma_{1}\Delta_{11}^{2}+\lambda_{2}\gamma_{2}\Delta_{22}^{2}+\sqrt{(\lambda_{1}\gamma_{1}\Delta_{11}^{2}-\lambda_{2}\gamma_{2}\Delta_{22}^{2})^{2}+4\lambda_{1}\lambda_{2}\gamma_{1}\gamma_{2}}}$
(13)

then

$\lim_{N\to\infty}F_{N}<\operatorname{RS}(\beta,h).$

Comparing with the classical SK model, since the RHS of (12) is reduced to the
critical $\beta$ for $h=0$ in SK model when $m=1$ and the similar form of RHS
of (12) and (13), one can reasonably guess that (12) gives an RS regime of MSK
model, and (13) is the AT line condition. The authors in [BSS19] proved that
(13) indeed gives RSB phase for $m=2$ using the idea in [Tal11b]*Chapter 13 by
a 1-RSB perturbation of the Parisi formula, and the uniqueness of
$\mathbf{q}_{\star}$ is essential in the proof of Theorem 1.2 part (ii),
without that, it is hard to give an accurate description of the AT line.

In this paper, we prove the uniqueness of $\mathbf{q}_{\star}$ in the general
species case under an analogous condition of (12). To fully generalize Theorem
1.2 in the $m>2$ case (_i.e.,_ incorporating the case $h>0$), one needs to
analyze the uniqueness of the solution to a nonlinear system. Unfortunately,
we are not able to prove that. However, assuming $\Delta^{2}$ is _positive-
definite_ , $\mathbf{q}_{\star}$ is unique when $h>0$ and the infimum in (10)
is not achieved at the boundary, we can prove the AT line condition for
general $m$. The idea is still based on the perturbation technique of the
Parisi formula. On the other hand, our variance analysis in Section 5 and
Section 6 suggests the AT line condition should be true from the RS side when
$\Delta^{2}$ is _positive-definite_.

###  1.3. Statement of the main results

First, by adapting Lalata’s argument, we prove that the MSK model is in the RS
regime when $\beta<\beta_{0}$. As noted in Section 1, our argument holds for
_indefinite_ $\Delta^{2}$ and proves the asymptotics of the free energy in the
MSK model with _indefinite_ $\Delta^{2}$. To the best of our knowledge, this
is the first result dealing with general _indefinite_ $\Delta^{2}$ in the
Ising MSK model.

For the rest of the article, we will assume the following:

######  Assumption 1.3.

Assume that $\Delta^{2}$ is a symmetric and invertible $m\times m$ nonzero
matrix with non-negative entries. However, it need not be positive-definite.

Before stating the RS phase diagram results, we need to generalize that the
set $\mathcal{C}(\beta,h)$ is singleton from $2$-species [BSS19] to general
$m$-species. We define

$\displaystyle\beta_{c}:=\rho(\Delta^{2}\Lambda)^{-\nicefrac{{1}}{{2}}}.$ (14)

where $\rho(A)$ is the _spectral radius_ or the largest absolute value of the
eigenvalues of $A$. In general, one has $\rho(A)\leqslant\|A\|$. But it is
easy to check that for symmetric $A$, $\rho(A)=\|A\|$.

######  Theorem 1.4 (Uniqueness of solution).

Assume that $\beta<\beta_{c}$. Then $\mathcal{C}(\beta,h)$ is a singleton set,
_i.e.,_ the system

$\displaystyle
q_{s}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h),\qquad
s=1,2,\ldots,m$ (15)

has a unique solution, where $\eta\sim\operatorname{N}(0,1)$.

######  Remark 1.5.

If we remove the invertiblity of $\Delta^{2}$ in Assumption 1.3, then one has
to replace the equations (15) by setting (8) to 0 and uniqueness holds for the
vector $\Delta^{2}\Lambda\mathbf{q}$.

The proof is similar to the SK model, which is because the RHS map is a
contraction for small $\beta$. For large $\beta$, the contraction is not true
anymore. In the SK model with $h>0$, Latala–Guerra lemma [Tal11a]*Proposition
1.3.8 tells us $\mathbf{q}$ is still unique. The proof is based on the
monotone property of a nonlinear function and intermediate value theorem. In
the MSK model, the analog of Latala–Guerra is not obvious since we are dealing
with a system of nonlinear equations (15). The authors in [BSS19] give proof
for the $m=2$ case using an elementary approach, but the idea is hard to
generalize for $m\geqslant 3$. Their proof holds only for _positive-definite_
$\Delta^{2}$ case. For a particular _indefinite_ MSK model, the deep Boltzmann
machine, the Latala-Guerra lemma [ACM21] has been extended to arbitrary depth
with Gaussian random field $h$. For _indefinite_ $\Delta^{2}$, the following
Proposition 1.6 gives the uniqueness of $\mathbf{q}_{\star}$ when $m=2$.

######  Proposition 1.6.

For $m=2$, assume $\Delta^{2}=\begin{pmatrix}a&b\\\ c&d\end{pmatrix}$ is
_indefinite_. If $h>0$, then the system (15) has a unique solution for all
$\beta>0$.

This proposition will be used to discuss AT line for some bipartite SK models
in Example 1.16. The proof of Proposition 1.6 and Theorem 1.4 will be given in
Section 7.

Now we define the symmetric _positive-definite_ matrix

$\displaystyle\mathcal{V}:=|\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}|.$
(16)

obtained by taking absolute values of all the eigenvalues in the spectral
decomposition of the symmetric matrix
$\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}$. We
also denote by

$\displaystyle\overline{\mathbf{R}}_{12}:=\mathbf{R}_{12}-\mathbf{q}.$ (17)

the centered overlap vector, where $\mathbf{q}$ is the unique solution to
(15). The following two theorems are about the overlap concentration. Here we
use $\nu(\cdot)=\operatorname{\mathds{E}}\langle\cdot\rangle$ to denote
expectation w.r.t. the disorder and the Gibbs measure, one can always check
the notation in Section 1.4.

######  Theorem 1.7 (Overlap concentration for general $h$).

Assume that $\beta<\beta_{0}:={\beta_{c}}/{\sqrt{4\alpha}}$, where
$\alpha=\alpha(\Delta^{2}):=1+\mathds{1}\\{\Delta^{2}\text{ is
indefinite}\\}$. For $2\eta<\beta_{c}^{2}-4\alpha\beta^{2}$, we have

$\displaystyle\nu(\exp(\eta
N\mathcal{P}(\overline{\mathbf{R}}_{12})))\leqslant\det(I-(2\eta+4\alpha\beta^{2})\mathcal{V})^{-1/2}$

where

$\displaystyle\mathcal{P}(\mathbf{x}):=\mathbf{x}^{\intercal}\Lambda^{\nicefrac{{1}}{{2}}}\mathcal{V}\Lambda^{\nicefrac{{1}}{{2}}}\mathbf{x}.$

Theorem 1.7 says that the overlap vector $\mathbf{R}_{12}$ concentrates around
$\mathbf{q}$ when $\beta<\beta_{0}$ for all $h\geqslant 0$, and its proof is
given in Section 2 by adapting Latala’s argument. Depending on the
definiteness of $\Delta^{2}$, we obtained two different RS regimes. The RS
regime for _indefinite_ $\Delta^{2}$ has an extra factor because the control
of derivative of interpolated Gibbs average in the proof is somewhat crude. To
get a sharper bound, one needs to work much harder, even unknown in SK.
However, the next theorem tells us that if $h=0$, one can prove the
concentration of overlap up to $\beta_{c}$, which is true even for
_indefinite_ $\Delta^{2}$.

######  Theorem 1.8 (Overlap Concentration for $h=0$).

If $h=0,\beta<\beta_{c}$, we have

$\displaystyle\nu\left(\exp\bigl{(}\frac{1}{8}(\beta_{c}^{2}-\beta^{2})N\mathcal{P}(\mathbf{R}_{12})\bigr{)}\right)\leqslant
K$ (18)

for some constant $K<\infty$ that does not depend on $N$.

The proof of Theorem 1.8 is given in Section 4 and is based on a quadratic
coupling argument [GT02b, Tal11b]. This theorem is expected to give the whole
RS regime when $h=0$ (by comparing with the classical SK model), even for
_indefinite_ $\Delta^{2}$.

With the control of overlap, we also prove the following.

######  Theorem 1.9 (LLN and CLT for free energy).

If $\beta<\beta_{0},h\geqslant 0$, then

$\displaystyle\left|F_{N}-\operatorname{RS}(\beta,h)\right|\leqslant\frac{K}{N},$
(19)

where $K$ is a constant that does not depend on $N$. Moreover, for $h>0$, we
have

$\displaystyle
N^{\nicefrac{{1}}{{2}}}\left(F_{N}-\operatorname{RS}(\beta,h)\right)\Rightarrow\operatorname{N}(0,b(\beta,h))\quad\text{as}\quad
N\to\infty,$ (20)

where

$\displaystyle\operatorname{RS}(\beta,h)$ $\displaystyle=\log
2+\sum_{s=1}^{m}\lambda_{s}\operatorname{\mathds{E}}\log\cosh(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h)+\frac{\beta^{2}}{4}\mathcal{Q}(\mathbf{1}-\mathbf{q}),$
$\displaystyle b(\beta,h)$
$\displaystyle:=\sum_{s=1}^{m}\lambda_{s}\operatorname{Var}(\log\cosh(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h))-\frac{\beta^{2}}{2}\mathcal{Q}(\mathbf{q}).$

and
$\mathcal{Q}(\mathbf{x}):=\mathbf{x}^{\intercal}\Lambda\Delta^{2}\Lambda\mathbf{x}$
denotes the quadratic form of $\Lambda\Delta^{2}\Lambda$.

The first part of Theorem 1.9 gives the RS solution of MSK model for general
$\Delta^{2}$. It’s easy to check that $\text{RS}(\beta,h)$ here has the same
form as $\text{RS}(\beta,h)$ in (10) derived from the Parisi formula in the
_positive-definite_ $\Delta^{2}$ case. However, in the bipartite case the
critical points in $\mathcal{C}(\beta,h)$ are saddle points, so the formula in
(7) is not true anymore. However, a modified Parisi formula is likely to be
true for $\Delta^{2}$ _indefinite_ at least in the RS regime (see Example 1.16
for further discussion). The second part is a Central Limit Theorem for the
free energy, which holds at the corresponding RS regime depending on
$\Delta^{2}$. The proof of Theorem 1.9 is in Section 3. For the CLT of free
energy in classical SK model, see [GT02, ALR87, Tin05, CN95] and references
therein.

######  Remark 1.10.

The proof of CLT for the free energy in the $\beta<\beta_{c},h=0$ case, is
similar to the classical SK model, see [ALR87] and [Tal11b]*Chapter 11. In
[ALR87], the authors proved the CLT using cluster expansion approach and later
[Tal11b] reproved it using moment method. Both the cluster expansion approach
and moment method are highly nontrivial for MSK model. However, the
characteristic function approach we used in our paper, can be generalize for
the MSK model with $\beta<\beta_{c},h=0$. Since we already have a good control
of overlap in Theorem 1.8, the next question is to determine the asymptotic
mean and variance of free energy. By the species-wise cavity approach
developed in Section 5 and Section 6, we evaluate

$\displaystyle c_{N}(\beta)$ $\displaystyle:=N\left(\log
2+\frac{1}{4}\beta^{2}\mathcal{Q}(\mathbf{1})\right)+\frac{1}{4}\log\det(I-\beta^{2}\Delta^{2}\Lambda),$
(21) $\displaystyle b(\beta)$
$\displaystyle:=\frac{1}{2}\left(-\log\det(I-\beta^{2}\Delta^{2}\Lambda)-\beta^{2}\operatorname{Tr}(\Delta^{2}\Lambda)\right).$
(22)

Then $\log Z_{N}(\beta)-c_{N}(\beta)\Rightarrow\operatorname{N}(0,b(\beta))$
as $N\to\infty$ can be carried over similarly as in [Tal11b]*Chapter 11, but
we do not pursue that here.

Next, we generalize the classical cavity approach to a species-wise cavity
method to analyze the fluctuation of the overlap vector. We choose a species
$s\in\\{1,2,\ldots,m\\}$ first with probability ${\lambda_{s}}$, then select a
spin uniformly in that species to be decoupled and finally compare this
decoupled system with the original one. The complete description is given in
Section 5. In [Pan15], Panchenko introduced another cavity idea when using the
Aizenman–Sims–Starr scheme to prove the lower bound of the Parisi formula. It
adds/decouples $k>1$ spins simultaneously, and those $k$ spins are distributed
into $m$ species according to the ratios stored in $\Lambda$. It might be
possible to do the second-moment computation along with this cavity idea, but
it will be more technically challenging since one needs to control $k>1$
spins.

For $(\beta,h)$ in the RS region proved in Theorem 1.7 and 1.8 with
$\mathbf{q}$ being the unique solution to (15), define

$\hat{q}_{s}:=\operatorname{\mathds{E}}\tanh^{4}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h),\qquad
s=1,2,\ldots,m.$

Define the $m\times m$ diagonal matrices
$\Gamma,\Gamma^{\prime},\Gamma^{\prime\prime}$, whose $s$-th diagonal entries
are respectively given by

$\displaystyle\gamma_{s}=1-2q_{s}+\hat{q}_{s}$
$\displaystyle=\operatorname{\mathds{E}}\operatorname{sech}^{4}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h),$
$\displaystyle\gamma^{\prime}_{s}=1-4q_{s}+3\hat{q}_{s}\qquad$
$\displaystyle\text{and}\qquad\gamma^{\prime\prime}_{s}=2q_{s}+q_{s}^{2}-3\hat{q}_{s}\qquad\text{
for }s=1,2,\ldots,m.$

Note that $-1/3\leqslant\gamma_{s}^{\prime}\leqslant\gamma_{s}\leqslant 1$ but
$\gamma^{\prime}_{s}$ can be negative. For large $h$, one can easily check
that $\gamma_{s}^{\prime}<-\gamma_{s}$. This fact will be used when we discuss
the AT line condition in Example 1.16. Let

$\displaystyle
U(0):=\nu(\overline{\mathbf{R}}_{12}\overline{\mathbf{R}}_{34}^{\intercal}),\quad
U(1):=\nu(\overline{\mathbf{R}}_{12}\overline{\mathbf{R}}_{13}^{\intercal})\text{
and }\quad
U(2):=\nu(\overline{\mathbf{R}}_{12}\overline{\mathbf{R}}_{12}^{\intercal})$
(23)

denote the variance-covariance matrices of the overlap vectors. Note that, all
the matrices $U(i)$’s are _positive-definite_. Using the generalized cavity
method we prove the following.

######  Theorem 1.11.

For $\beta<\beta_{0}$ or $\\{\beta<\beta_{c},h=0\\}$, the matrices

$\displaystyle\hat{U}(0)$ $\displaystyle:=-3U(0)+2U(1),$
$\displaystyle\hat{U}(1)$ $\displaystyle:=3U(0)-4U(1)+U(2),$
$\displaystyle\hat{U}(2)$ $\displaystyle:=U(0)-2U(1)+U(2)$

satisfy the following equations

$\displaystyle\begin{split}\operatorname{sym}\bigl{(}\bigl{(}I-\beta^{2}\Gamma\Delta^{2}\Lambda\bigr{)}\cdot\hat{U}(2)\bigr{)}&=\frac{1}{N}\cdot\Gamma\Lambda^{-1}+\boldsymbol{\mathfrak{R}},\\\
\operatorname{sym}\bigl{(}\bigl{(}I-\beta^{2}\Gamma^{\prime}\Delta^{2}\Lambda\bigr{)}\cdot\hat{U}(1)\bigr{)}&=\frac{1}{N}\cdot\Gamma^{\prime}\Lambda^{-1}+\boldsymbol{\mathfrak{R}},\\\
\text{ and
}\operatorname{sym}\bigl{(}\bigl{(}I-\beta^{2}\Gamma^{\prime}\Delta^{2}\Lambda\bigr{)}\cdot\hat{U}(0)\bigr{)}&=\beta^{2}\operatorname{sym}\bigl{(}\Gamma^{\prime\prime}\Delta^{2}\Lambda\cdot\hat{U}(1)\bigr{)}+\frac{1}{N}\cdot\Gamma^{\prime\prime}\Lambda^{-1}+\boldsymbol{\mathfrak{R}},\end{split}$
(24)

where $\operatorname{sym}(A)=(A+A^{\intercal})/2$ and
$\boldsymbol{\mathfrak{R}}$ is some $m\times m$ symmetric matrix with

$\max_{1\leqslant p,q\leqslant m}|\boldsymbol{\mathfrak{R}}_{p,q}|\leqslant
K\biggl{(}N^{-3/2}+\sum_{s=1}^{m}\nu(\left|R_{12}^{s}-q_{s}\right|^{3})\biggr{)}$

for some constant $K$ independent of $N$.

Unfortunately, we do not have a nice interpretation of the coefficients
appearing in the system of linear equations in Theorem 1.11. This is actually
related to an open problem by Talagrand, where he asked to identify the
underlying algebraic structure (see [Tal11a, Research Problem 1.8.3]).

Before we discuss the Theorem 1.11, let’s recall the definition of stable
matrices.

######  Definition 1.12.

Square matrix $A$ is _stable_ if all the eigenvalues of $A$ have strictly
negative real part.

Note that the linear equations in (24) are examples of continuous Lyapunov
equation of the form $AX+XA^{\intercal}+Q=0$ where $A,Q$ are symmetric
matrices. Those equations also appear in control theory (see _e.g.,_
[AM07book]). Moreover, the existence and uniqueness of the solution are
equivalent to the matrix $A$ being stable. It is an exciting question to
connect equations (24) with an appropriate control problem. In our case,
solving the system of equations given in (24), we get the asymptotic variance
of overlap.

######  Theorem 1.13 (Asymptotic variance of overlap vector).

For $\beta<\beta_{0}$ or $\\{\beta<\beta_{c},h=0\\}$, we have

$\displaystyle N\cdot U(i)\to\Sigma(i)\quad\text{as}\quad N\to\infty$

where $\Sigma(i),i=0,1,2$ satisfies the following

$\displaystyle\begin{split}\Sigma(0)-2\Sigma(1)+\Sigma(2)&=\Gamma(I-\beta^{2}\Delta^{2}\Lambda\Gamma)^{-1}\Lambda^{-1},\\\
3\Sigma(0)-4\Sigma(1)+\Sigma(2)&=\Gamma^{\prime}(I-\beta^{2}\Delta^{2}\Lambda\Gamma^{\prime})^{-1}\Lambda^{-1}.\end{split}$
(25)

######  Remark 1.14.

$\Sigma(i),i=0,1,2$ also solve a third equation (76) besides (25), but it does
not possess a simple form, see the details in Section 6.

The asymptotic variance of overlap is given as solution to the linear system
in (24). We recall the definition of _stable_ in the Definition 1.12. To get a
unique solution from the system (24), one needs the matrices
$\beta^{2}\Gamma\Delta^{2}\Lambda-I,\beta^{2}\Gamma^{\prime}\Delta^{2}\Lambda-I$
to be stable, _i.e.,_

$\displaystyle\beta^{2}\max_{\lambda\in\text{spec}(\Gamma\Delta^{2}\Lambda)}\Re(\lambda)<1\quad\text{and}\quad\beta^{2}\max_{\lambda\in\text{spec}(\Gamma^{\prime}\Delta^{2}\Lambda)}\Re(\lambda)<1.$
(26)

Note that $\gamma_{s}>0$ for all $s$ and $\Gamma\Delta^{2}\Lambda$ is similar
to a symmetric matrix. Thus all eigenvalues of $\Gamma\Delta^{2}\Lambda$ are
real, and by Perron–Frobenius theorem

$\max_{\lambda\in\text{spec}(\Gamma\Delta^{2}\Lambda)}\Re(\lambda)=\rho(\Gamma\Delta^{2}\Lambda)$

and the condition (26) is equivalent to

$\displaystyle\beta^{2}\max\\{\rho(\Gamma\Delta^{2}\Lambda),\max_{\lambda\in\text{spec}(\Gamma^{\prime}\Delta^{2}\Lambda)}\Re(\lambda)\\}<1.$
(27)

For _positive-definite_ $\Delta^{2}$, the matrix
$\Gamma^{\prime}\Delta^{2}\Lambda$ is similar to a symmetric matrix and has
real eigenvalues. Moreover, using $\gamma^{\prime}_{s}\leqslant\gamma_{s}$,
one can get

$\displaystyle\max_{\lambda\in\text{spec}(\Gamma^{\prime}\Delta^{2}\Lambda)}\Re(\lambda)\leqslant\rho(\Gamma\Delta^{2}\Lambda).$
(28)

Thus the results in Theorem 1.13 will not be true unless
$\beta^{2}<\rho(\Gamma\Delta^{2}\Lambda)^{-1}$. This indicates that
$\beta^{2}\rho(\Gamma\Delta^{2}\Lambda)=1$ is the AT line condition from the
RS side for _positive-definite_ $\Delta^{2}$. To rigorously prove that
$\beta^{2}<\rho(\Gamma\Delta^{2}\Lambda)^{-1}$ gives RS phase is still open
even in classical SK model. In [JT17], the authors proved it except for a
bounded region close to the phase boundary by using the Parisi formula.
Recently, Bolthausen [Bol14, Bol19] developed an iterative TAP approach to SK
model, where he also discussed some possible ideas to prove
$\lim_{N\to\infty}F_{N}=\text{RS}(\beta,h)$ below the AT line. For $h=0$, when
$\beta^{2}<\beta_{c}^{2}=\rho(\Delta^{2}\Lambda)^{-1}$ and $\Delta^{2}$ is
general, it is easy to check that condition (26) is satisfied, more details
can be found in Section 6.

The following Theorem 1.15 states the RSB condition when $\Delta^{2}$ is
_positive-definite_ , under the assumption that $\mathcal{C}(\beta,h)$ is a
singleton for $h>0$. Without this, it is not easy to give an accurate
description of the AT line. We assume $\Delta^{2}$ is _positive-definite_ , as
the proof is based on perturbation of the Parisi formula known only in that
case in [Pan15]. Recall that we also need to assume that the infimum in (10)
can only be achieved in $\mathcal{C}(\beta,h)$ but not on the boundary.

######  Theorem 1.15 (RSB condition).

Assume that $\Delta^{2}$ is _positive-definite_ and for some $\beta,h>0$,
$\mathcal{C}(\beta,h)$ is a singleton. We further assume that the infimum in
(10) is not achieved on the boundary. Let $\mathbf{q}$ be the unique critical
point in $\mathcal{C}(\beta,h)$. Define
$\Gamma:=\operatorname{diag}(\gamma_{1},\gamma_{2},\ldots,\gamma_{m})$ where
$\gamma_{s}:=\operatorname{\mathds{E}}\operatorname{sech}^{4}(\beta\eta\sqrt{Q_{1}^{s}}+h)\
\text{for}\ s=1,2,\ldots,m.$ If
$\beta>\rho(\Gamma\Delta^{2}\Lambda)^{-\nicefrac{{1}}{{2}}}$, then

$\lim_{N\to\infty}F_{N}<\operatorname{RS}(\beta,h).$

Note that

$\displaystyle\beta_{\textsc{AT}}(h):=\rho(\Gamma\Delta^{2}\Lambda)^{-\nicefrac{{1}}{{2}}}$
(29)

seems to be the AT line in the MSK model for _positive-definite_ $\Delta^{2}$,
and when $h=0$, $\beta_{\textsc{AT}}(0)=\beta_{c}$ gives the critical
temperature defined in (14), see the left phase diagram picture of Figure 1.
When $h=0,\beta>\beta_{c}$, there seems to be at least one non-zero solution
in $\mathcal{C}(\beta,0)$. However, to prove RSB by the perturbation argument
is still technically challenging. For the general MSK model, here is a class
of examples with a non-zero solution for $h=0,\beta>\beta_{c}$.

In the classical SK model, one can easily check that for $h=0,\beta>1$, there
exists a non-zero solution $\hat{q}(\beta)\in(0,1)$ to the equation
$q=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{q})$. Take
$\Delta^{2},\Lambda$ such that $\mathbf{1}$ is an eigenvector of
$\Delta^{2}\Lambda$ with associated eigenvalue $\theta$. Then
$\beta_{c}=\rho(\Delta^{2}\Lambda)^{-\nicefrac{{1}}{{2}}}=\theta^{-\nicefrac{{1}}{{2}}}$
and for $\beta>\beta_{c}$, we have $\beta\sqrt{\theta}>1$. It follows that
$\mathbf{q}=\hat{q}(\beta\sqrt{\theta})\mathbf{1}$ is a non-zero solution to
the equations

$\displaystyle
q_{s}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}})\text{
for }s=1,2,\cdots,m.$

Finally, we point out that the AT line condition for $\Delta^{2}$ _indefinite_
seems to be more complicated than the _positive-definite_ case. Because for
_indefinite_ $\Delta^{2}$, the inequality (28) could fail depending on the
values of $(\beta,h)$. For fixed $\beta>\beta_{c}$ and large $h$, we observe
$\max_{\lambda\in\text{spec}(\Gamma^{\prime}\Delta^{2}\Lambda)}\Re(\lambda)\geqslant\rho(\Gamma\Delta^{2}\Lambda)$
in some examples. However, it is not obvious whether in the $\Delta^{2}$
_indefinite_ case, the AT line is given by
$\beta^{2}\rho(\Gamma\Delta^{2}\Lambda)=1$ (see from the right picture in
Figure 1).

Figure 1. Phase diagram of MSK model: the left and right one are expected for
_positive-definite_ and _indefinite_ $\Delta^{2}$ respectively. The green
shaded region corresponds to RS regime proved in Theorem 1.7 and 1.8. The
expected RSB phase is labelled in both cases. In both pictures, the red curve
and blue curve correspond to $\beta^{2}\rho(\Gamma\Delta^{2}\Lambda)=1$ and
$\beta^{2}\max_{\lambda\in\text{spec}(\Gamma^{\prime}\Delta^{2}\Lambda)}\Re(\lambda)=1$
respectively. The right picture for _indefinite_ $\Delta^{2}$ is based on the
evidence from some examples, see Example 1.16.

######  Example 1.16.

Take $\Lambda=\begin{pmatrix}1/2&0\\\ 0&1/2\end{pmatrix}$ and
$\Delta^{2}=\begin{pmatrix}0&1\\\ 1&0\end{pmatrix}$. This is the bipartite SK
model, an example of _indefinite_ MSK model. For $h>0$, by Proposition 1.6,
the following system

$\begin{cases}q_{1}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{q_{2}/2}+h)\\\
q_{2}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{q_{1}/2}+h)\end{cases}$

has a unique solution
$q_{1}=q_{2}=q=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{\frac{q}{2}}+h)$.
Similarly, $\gamma_{1}=\gamma_{2}=\gamma$ and
$\gamma_{1}^{\prime}=\gamma_{2}^{\prime}=\gamma^{\prime}$.

Then we have

$\displaystyle\rho(\Gamma\Delta^{2}\Lambda)=\frac{\gamma}{2}\quad\text{and}\quad\max_{\lambda\in\text{spec}(\Gamma^{\prime}\Delta^{2}\Lambda)}\Re(\lambda)=\rho(\Gamma^{\prime}\Delta^{2}\Lambda)=\frac{\left|\gamma^{\prime}\right|}{2}.$

We know that $\gamma\geqslant\gamma^{\prime},\gamma>0$ always holds, but
$\left|\gamma^{\prime}\right|\leqslant\gamma$ is not always true (for large
$h$), and it depends on the strength of external field $h$. The suggested AT
line condition is $\beta^{2}\max\\{\gamma,|\gamma^{\prime}|\\}<1$.

For the functional $\mathscr{P}_{\operatorname{RS}}(\mathbf{q})$, the Hessian
w.r.t. $\mathbf{q}$ is a positive multiple of

$H=\Delta^{2}-\beta^{2}\Lambda\Gamma^{\prime}=\begin{pmatrix}-\beta^{2}\gamma^{\prime}/2&1\\\
1&-\beta^{2}\gamma^{\prime}/2\end{pmatrix}.$

$H$ is an indefinite matrix in the replica symmetric region, as the
eigenvalues of $H$ are given by $-\beta^{2}\gamma^{\prime}/2\pm 1$ and the RS
region seems to be $\beta^{2}\max(\gamma,|\gamma^{\prime}|)/2<1$. So, although
the Parisi formula evaluated at the stationary point gives the limit for the
free energy (proved for small $\beta$), it is not the minimizer of the Parisi
formula directly taken from the positive-definite case.

In general, we conjecture the RSB condition for _indefinite_ $\Delta^{2}$ as
follows.

######  Conjecture 1.17.

For general indefinite $\Delta^{2}$, if

$\beta^{2}\max\\{\rho(\Gamma\Delta^{2}\Lambda),\max_{\lambda\in\text{spec}(\Gamma^{\prime}\Delta^{2}\Lambda)}\Re(\lambda)\\}>1,$

the system is in RSB phase.

The above conjecture seems far to solve since the classical idea is based on
perturbation of the Parisi formula and the Parisi formula in general
_indefinite_ case is still in mystery. Besides that, we also need to
understand the solution to (15) for _indefinite_ $\Delta^{2}$. To prove the AT
line condition from the RS side for _indefinite_ $\Delta^{2}$ is even more
challenging.

In order to illustrate the generality and correctness of our results, we look
at some models studied in other papers as particular cases in our setting.

######  Example 1.18.

Take $\Lambda=\begin{pmatrix}\lambda_{1}&0\\\ 0&\lambda_{2}\end{pmatrix}$ and
$\Delta^{2}$ as in (11), then

$\beta_{c}=\rho(\Delta^{2}\Lambda)=\frac{1}{2}\left(\lambda_{1}\Delta_{11}^{2}+\lambda_{2}\Delta_{22}^{2}+\sqrt{(\lambda_{1}\Delta_{11}^{2}-\lambda_{2}\Delta_{22}^{2})^{2}+4\lambda_{1}\lambda_{2}}\right).$

Similarly, one can check

$\beta_{AT}(h)=\rho(\Gamma\Delta^{2}\Lambda)=\frac{1}{2}\left(\lambda_{1}\gamma_{1}\Delta_{11}^{2}+\lambda_{2}\gamma_{2}\Delta_{22}^{2}+\sqrt{(\lambda_{1}\gamma_{1}\Delta_{11}^{2}-\lambda_{2}\gamma_{2}\Delta_{22}^{2})^{2}+4\lambda_{1}\lambda_{2}\gamma_{1}\gamma_{2}}\right).$

Theorem 1.4 and 1.15 extend the critical temperature (12) and AT line
condition (13) to the general MSK model.

######  Example 1.19.

Take $\Lambda=\begin{pmatrix}\lambda_{1}&0\\\ 0&\lambda_{2}\end{pmatrix}$ and
$\Delta^{2}=\begin{pmatrix}0&1\\\ \theta^{2}&0\end{pmatrix}$. For $\theta=1$,
this corresponds to the bipartite SK model. In [BL17] the authors study the
BSSK model on the sphere with $h=0$ and compute
$\beta_{c}=(\lambda_{1}\lambda_{2})^{-1/4}$. In our case, we have
$\beta_{c}=\rho(\Delta^{2}\Lambda)^{-\nicefrac{{1}}{{2}}}$ in (14) and by a
simple computation, one gets
$\rho(\Delta^{2}\Lambda)^{-\nicefrac{{1}}{{2}}}=(\lambda_{1}\lambda_{2}\theta^{2})^{-1/4}$,
which agrees with $\beta_{c}$ in [BL17].

###  1.4. Notations

  1. (1)

$\rho(A)$ is the _spectral radius_ or the largest absolute value of the
eigenvalues of $A$.

  2. (2)

$\|A\|$ is the operator norm of $A$.

  3. (3)

$\Delta^{2}=(\\!(\Delta_{s,t}^{2})\\!)_{s,t=1}^{m}$ is the covariance matrix
of Gaussian disorder.

  4. (4)

$\Lambda:=\operatorname{diag}(\lambda_{1},\lambda_{2},\ldots,\lambda_{m})$ is
the diagonal matrix, whose diagonal entries is the ratio of the spins in each
species.

  5. (5)

$\mathcal{V}:=|\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}|$
which is obtained by taking absolute values of all the eigenvalues in the
spectral decomposition of the matrix
$\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}$.

  6. (6)

$\mathbf{R}_{12}=(R^{(1)}_{12},R^{(2)}_{12},\ldots,R^{(m)}_{12})$ denotes the
overlap vector.

  7. (7)

$\overline{\mathbf{R}}_{12}:=\mathbf{R}_{12}-\mathbf{q}$ is the centered
overlap vector, and $\mathbf{q}$ is defined through (15).

  8. (8)

$Q^{s}:=(\Delta^{2}\Lambda\mathbf{q})_{s}\text{ for }s=1,2,\ldots,m$.

  9. (9)

$\Gamma:=\operatorname{diag}(\gamma_{1},\gamma_{2},\ldots,\gamma_{m})$ where
$\gamma_{s}:=\operatorname{\mathds{E}}\operatorname{sech}^{4}(\beta\eta\sqrt{Q^{s}}+h)\text{
for }s=1,2,\ldots,m$.

  10. (10)

$\langle
f(\boldsymbol{\sigma})\rangle:=\sum_{\boldsymbol{\sigma}\in\Sigma_{N}}f(\boldsymbol{\sigma})\cdot{\exp(H_{N}(\boldsymbol{\sigma}))}/{Z_{N}}$
denotes the Gibbs average of function $f$ on $\Sigma_{N}$.

  11. (11)

$\nu_{t}(\cdot)$ denotes the expectation w.r.t. the Gibbs randomness with
interpolated Hamiltonian and disorder.

  12. (12)

$\mathcal{Q}(\mathbf{x})=\mathbf{x}^{\intercal}\Lambda\Delta^{2}\Lambda\mathbf{x}$
denotes the quadratic form associated with the symmetric matrix
$\Lambda\Delta^{2}\Lambda$. For a general symmetric matrix $A$, we will use
$\mathcal{Q}(\mathbf{x},A)$ for $\mathbf{x}^{\intercal}A\mathbf{x}$. We also
write
$\mathcal{P}(\mathbf{x})=\mathcal{Q}(\mathbf{x},\Lambda^{\nicefrac{{1}}{{2}}}\mathcal{V}\Lambda^{\nicefrac{{1}}{{2}}})$.

###  1.5. Roadmap

The paper is structured as follows. Section 2 is mainly about the proof of
Theorem 1.7, where the smart path interpolation argument is applied in the MSK
model to prove the concentration of overlap at high temperature. The argument
also holds for indefinite $\Delta^{2}$, while in this case, the concentration
happens in a different regime. Section 3 is mainly about the proof of Theorem
1.9, where we prove the RS solution of the MSK model and a Central Limit
Theorem for the free energy with a non-zero external field. In Section 4, we
study the MSK model without an external field, where the proof of Theorem 1.8
is presented. A Central Limit Theorem for the free energy with zero external
field can also be proved in a similar way as in [ALR87], but we do not pursue
this direction here. In Section 5, we develop a generalized cavity method to
study the MSK model. Using second-moment computations, we derive a linear
system of overlap functions and solve it using linear algebraic methods to get
the variance-covariance structure of overlap vectors in Section 6. Basically,
Section 5 and 6 contain the proofs of Theorem 1.11 and 1.13 respectively. In
Section 7, for _positive-definite_ MSK model, we give the proof of Theorem
1.15 by using the perturbation argument of the Parisi formula. Finally Section
8 contains discussion and some further questions.

##  2\. Concentration of overlap in MSK

Originally Latala’s argument [La02] in the SK model is used to prove the
concentration of overlap in part of the RS regime. The idea is based on
Guerra’s interpolation of two different spin glass models, one of which is the
fully decoupled model, _i.e.,_ the associated Gibbs measure is a product
measure, while the other model is the standard SK model of our interest. With
the concentration of overlap in the decoupled model, one can obtain the
concentration of overlap in the SK model by controlling the derivative of the
interpolated model. For the details of this classical story, see Section 1.3
and 1.4 in [Tal11a]. This section employs a similar idea in the MSK model to
prove the concentration of overlap vectors. With this concentration result, we
also prove a part of the RS phase diagram.

Following Guerra’s interpolation, given two independent centered Gaussian
vectors
$\boldsymbol{u}=(u_{\boldsymbol{\sigma}}),\boldsymbol{v}=(v_{\boldsymbol{\sigma}})$
indexed by $\boldsymbol{\sigma}\in\Sigma_{N}$, consider the interpolation
given by $u_{\boldsymbol{\sigma}}(t)=\sqrt{t}\cdot
u_{\boldsymbol{\sigma}}+\sqrt{1-t}\cdot v_{\boldsymbol{\sigma}}$, for
$t\in[0,1]$. Suppose $F$ is a twice differentiable function allowing us to use
Gaussian integration by parts. Taking
$\phi(t):=\operatorname{\mathds{E}}F(\boldsymbol{u}(t))$ and using Gaussian
integration by parts we get the following.

######  Proposition 2.1.

For $\boldsymbol{\sigma},\boldsymbol{\tau}\in\Sigma_{N}$,

$\phi^{\prime}(t)=\sum_{\boldsymbol{\sigma},\boldsymbol{\tau}}U(\boldsymbol{\sigma},\boldsymbol{\tau})\operatorname{\mathds{E}}\left(\frac{\partial^{2}F}{\partial
x_{\boldsymbol{\sigma}}x_{\boldsymbol{\tau}}}(\boldsymbol{u}(t))\right)$

where
$U(\boldsymbol{\sigma},\boldsymbol{\tau}):=\frac{1}{2}(\operatorname{\mathds{E}}u_{\boldsymbol{\sigma}}u_{\boldsymbol{\tau}}-\operatorname{\mathds{E}}v_{\boldsymbol{\sigma}}v_{\boldsymbol{\tau}})$.

Given $(\eta_{i})_{i=1,2,\ldots,N}$ i.i.d. standard Gaussian random variables
and independent of $(g_{ij})_{1\leqslant i<j\leqslant N}$, we take

$u_{\boldsymbol{\sigma}}=\frac{\beta}{\sqrt{N}}\sum_{i<j}g_{ij}\sigma_{i}\sigma_{j},\quad
v_{\boldsymbol{\sigma}}=\beta\sum_{s=1}^{m}\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}\cdot\sum_{i\in
I_{s}}\eta_{i}\sigma_{i},\quad
w_{\boldsymbol{\sigma}}=\exp\bigl{(}h\sum_{i}\sigma_{i}\bigr{)}$

and $F(\boldsymbol{x})=\frac{1}{N}\log Z_{N}$ to be the free energy, where
$Z_{N}=\sum_{\boldsymbol{\sigma}}w_{\boldsymbol{\sigma}}e^{x_{\boldsymbol{\sigma}}}$.
In this case, the interpolated Hamiltonian is

$\displaystyle
H_{N,t}(\boldsymbol{\sigma}):=u_{t}(\boldsymbol{\sigma})=\sqrt{t}u_{\boldsymbol{\sigma}}+\sqrt{1-t}v_{\boldsymbol{\sigma}}.$
(30)

Moreover, we have

$\displaystyle\operatorname{\mathds{E}}u_{\boldsymbol{\sigma}}u_{\boldsymbol{\tau}}$
$\displaystyle=\frac{N\beta^{2}}{2}\left(\sum_{s,t}\lambda_{s}R^{s}\Delta_{s,t}^{2}\lambda_{t}R^{t}-\frac{1}{N}\sum_{s}\lambda_{s}\Delta_{s,s}^{2}\right),$
$\displaystyle\operatorname{\mathds{E}}v_{\boldsymbol{\sigma}}v_{\boldsymbol{\tau}}$
$\displaystyle=N\beta^{2}\sum_{s}\lambda_{s}R^{s}Q^{s}=N\beta^{2}\sum_{s,t}\lambda_{s}R^{s}\Delta_{s,t}^{2}\lambda_{t}q_{t}.$

In particular,

$\displaystyle U(\boldsymbol{\sigma},\boldsymbol{\tau})$
$\displaystyle=\frac{1}{2}(\operatorname{\mathds{E}}u_{\boldsymbol{\sigma}}u_{\boldsymbol{\tau}}-\operatorname{\mathds{E}}v_{\boldsymbol{\sigma}}v_{\boldsymbol{\tau}})=\frac{N\beta^{2}}{4}\left(\mathcal{Q}(\overline{\mathbf{R}})-\mathcal{Q}(\mathbf{q})-\frac{1}{N}\sum_{s}\lambda_{s}\Delta_{s,s}^{2}\right)$
$\displaystyle\text{ and }U(\boldsymbol{\sigma},\boldsymbol{\sigma})$
$\displaystyle=\frac{N\beta^{2}}{4}\left(\mathcal{Q}(\mathbf{1}-\mathbf{q})-\mathcal{Q}(\mathbf{q})-\frac{1}{N}\sum_{s}\lambda_{s}\Delta_{s,s}^{2}\right).$

By a simple computation, we have

$\frac{\partial F}{\partial x_{\boldsymbol{\sigma}}}=\frac{1}{N}\cdot
w_{\boldsymbol{\sigma}}e^{x_{\boldsymbol{\sigma}}}/Z_{N},\quad\frac{\partial^{2}F}{\partial
x_{\boldsymbol{\sigma}}\partial
x_{\boldsymbol{\tau}}}=\frac{1}{N}\cdot(\mathds{1}_{\boldsymbol{\sigma}=\boldsymbol{\tau}}-w_{\tau}e^{x_{\boldsymbol{\tau}}}/Z_{N})\cdot
w_{\boldsymbol{\sigma}}e^{x_{\boldsymbol{\sigma}}}/Z_{N}$

and we use this to rewrite the formula of $\phi^{\prime}(t)$,

$\displaystyle\phi^{\prime}(t)$
$\displaystyle=\frac{1}{N}\operatorname{\mathds{E}}(\langle
U(\boldsymbol{\sigma},\boldsymbol{\sigma})\rangle_{t}-\langle
U(\boldsymbol{\sigma},\boldsymbol{\tau})\rangle_{t})=\frac{\beta^{2}}{4}\left(\mathcal{Q}(\mathbf{1}-\mathbf{q})-\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{12}))\right).$
(31)

In the following section, we will prove that $\mathbf{R}_{12}$ is concentrated
around $\mathbf{q}$ at some high temperature regime, _i.e.,_ the term
$\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{12}))$ is very small, and this
concentration property will enable us to prove the RS solution in MSK model.

First, we will prove that the quantity
$\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{12}))$ is small. The basic idea is
to prove some exponential moments is non-increasing in $t$, then by
controlling $\nu_{0}(\mathcal{Q}(\overline{\mathbf{R}}_{12}))$, we can analyze
the other side: $\nu(\mathcal{Q}(\overline{\mathbf{R}}_{12}))$ along the
interpolation path. In order to prove the Theorem 1.7, we need to study the
property of $\nu_{t}^{\prime}$, _i.e.,_ the derivative of $\nu_{t}$ with
respect to $t$. The following lemma gives the expression of
$\nu_{t}^{\prime}$.

######  Lemma 2.2 ( [Tal11a]*Page 33).

If $f$ is a function defined on $\Sigma_{N}^{n}$, then

$\displaystyle\nu_{t}^{\prime}(f)$ $\displaystyle=\sum_{1\leqslant
l,l^{\prime}\leqslant
n}\nu_{t}(U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{l^{\prime}})f)+n(n+1)\nu_{t}(U(\boldsymbol{\sigma}^{n+1},\boldsymbol{\sigma}^{n+2})f)$
$\displaystyle\qquad\qquad\qquad-2n\sum_{l\leqslant
n}\nu_{t}(U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{n+1})f)-n\nu_{t}(U(\boldsymbol{\sigma}^{n+1},\boldsymbol{\sigma}^{n+1})f).$

Note that
$\nu_{t}(U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{l})f)=\nu_{t}(U(\boldsymbol{\sigma}^{n+1},\boldsymbol{\sigma}^{n+1})f)$,
therefore the above formula can be simplified one more step as follows:

$\displaystyle\nu_{t}^{\prime}(f)$ $\displaystyle=\sum_{1\leqslant
l<l^{\prime}\leqslant
n}2\nu_{t}(U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{l^{\prime}})f)+n(n+1)\nu_{t}(U(\boldsymbol{\sigma}^{n+1},\boldsymbol{\sigma}^{n+2})f)-2n\sum_{l\leqslant
n}\nu_{t}(U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{n+1})f).$ (32)

When $n=2$, simplifying we have

$\displaystyle\nu_{t}^{\prime}(f)$
$\displaystyle=\frac{N\beta^{2}}{2}\bigl{(}\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{12})f)-2\sum_{l\leqslant
2}\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{l,3})f)+3\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{34})f)\bigr{)}.$
(33)

Next, we prove a useful lemma comparing the quadratic forms of overlap under
$\nu_{t}^{\prime}(\cdot)$.

######  Lemma 2.3.

Consider $\rho>0$, for $l,l^{\prime}\in\mathbb{N}$, we have

$\displaystyle\left|\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{\ell\ell^{\prime}})\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))\right|$
$\displaystyle\leqslant\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12})\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))$

where
$\mathcal{P}(\mathbf{x})=\mathbf{x}^{\intercal}\Lambda^{\nicefrac{{1}}{{2}}}\mathcal{V}\Lambda^{\nicefrac{{1}}{{2}}}\mathbf{x}$.

###### Proof of Lemma 2.3..

Let $i:=\left|\\{\ell,\ell^{\prime}\\}\cap\\{1,2\\}\right|$, we prove the
above inequality case by case. For $i=2$, we need to prove

$\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{12})\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))\leqslant\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12})\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))$

the inequality is obvious by the fact that
$|\mathcal{Q}(\mathbf{x})|\leqslant\mathcal{P}(\mathbf{x})$. For $i=1$ and
$i=0$, it can be proved by Hölder inequalities. We provide a proof for the
case $i=1$. First, the general form of Hölder inequality w.r.t. $\nu_{t}$ is:

$\nu_{t}(f_{1}f_{2})\leqslant\nu_{t}(f_{1}^{\tau_{1}})^{1/\tau_{1}}\nu_{t}(f_{2}^{\tau_{2}})^{1/\tau_{2}}\quad\text{for}\quad
f_{1},f_{2}\geqslant 0,\quad\frac{1}{\tau_{1}}+\frac{1}{\tau_{2}}=1$

Without loss of generality, we assume $\ell=1$ and $\ell^{\prime}=3$, then one
needs to prove

$\left|\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{13})\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))\right|\leqslant\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12})\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))$

Since

$\displaystyle\left|\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{13})\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))\right|$
$\displaystyle\leqslant\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{13})\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))$ $\displaystyle=\sum_{k}\frac{(\rho
N)^{k}}{k!}\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{13})\mathcal{P}(\overline{\mathbf{R}}_{12})^{k})$
(34)

Applying Hölder inequality with $\tau_{1}=k+1,\tau_{2}=(k+1)/k$, we get

$\displaystyle\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{13})\cdot\mathcal{P}(\overline{\mathbf{R}}_{12})^{k})$
$\displaystyle\leqslant\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{13})^{k+1})^{1/(k+1)}\cdot\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12})^{k+1})^{k/(k+1)}=\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12})^{k+1})$

the last equality is by the symmetry among replicas. Combining with the
exponential expansion (2), we proved the case $i=1$. The case $i=0$ follows in
a similar fashion. $\blacksquare$

######  Remark 2.4.

In the Lemma 2.3, if $\Delta^{2}$ is non-negative definite, then
$\mathcal{P}(\overline{\mathbf{R}}_{12})=\mathcal{Q}(\overline{\mathbf{R}}_{12})$,
thus the similar inequalities of Lemma 2.3 also hold by replacing
$\mathcal{P}(\overline{\mathbf{R}}_{12})$ with
$\mathcal{Q}(\overline{\mathbf{R}}_{12})$ and the proof is similar but
simpler.

######  Corollary 2.5.

If $\kappa>0$, then

$\displaystyle\nu_{t}^{\prime}(\exp(\kappa
N\mathcal{P}(\overline{\mathbf{R}}_{12}))\leqslant 2\alpha
N\beta^{2}\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12})\exp(\kappa
N\mathcal{P}(\overline{\mathbf{R}}_{12}))),$

where $\alpha=1+\mathds{1}\\{\Delta^{2}\text{ is indefinite}\\}$ as in Theorem
1.7.

###### Proof of Corollary 2.5..

Combined with the Lemma 2.3 and the expression of $\nu_{t}^{\prime}(f)$ with
$f=\exp(\kappa N\mathcal{P}(\overline{\mathbf{R}}_{12}))$ in (33), and the
symmetry between replicas, we can do the following analysis. When $\Delta^{2}$
is non-negative definite, as in Remark 2.4,
$\mathcal{P}(\overline{\mathbf{R}}_{12})=\mathcal{Q}(\overline{\mathbf{R}}_{12})$.
In this case, by just dropping the terms

$-2\sum_{l\leqslant 2}\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{l,3})f)$

in (33), and applying the inequalities in Lemma 2.3, we get the upper bound of
$\nu_{t}^{\prime}(f)$, which is
$2N\beta^{2}\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12})\exp(\kappa
N\mathcal{P}(\overline{\mathbf{R}}_{12})))$. For the case $\Delta^{2}$ is
indefinite, similarly we just apply all the inequalities in Lemma 2.3 to get
the bound, _i.e.,_
$4N\beta^{2}\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12})\exp(\kappa
N\mathcal{P}(\overline{\mathbf{R}}_{12})))$. $\blacksquare$

######  Corollary 2.6.

For $t<\rho/(2\alpha\beta^{2})$, the function

$t\mapsto\nu_{t}(\exp((\rho-2\alpha
t\beta^{2})N\mathcal{P}(\overline{\mathbf{R}}_{12})))$

is non-increasing.

###### Proof of Corollary 2.6..

Taking the derivative of $\nu_{t}(\exp((\rho-2\alpha
t\beta^{2})N\mathcal{P}(\overline{\mathbf{R}}_{12})))$ w.r.t. $t$, we get

$\displaystyle\frac{d}{dt}$ $\displaystyle(\nu_{t}(\exp((\rho-2\alpha
t\beta^{2})N\mathcal{P}(\overline{\mathbf{R}}_{12}))))$
$\displaystyle=\nu_{t}^{\prime}(\exp((\rho-2\alpha
t\beta^{2})N\mathcal{P}(\overline{\mathbf{R}}_{12})))-2\alpha
N\beta^{2}\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12})\exp((\rho-2\alpha
t\beta^{2})N\mathcal{P}(\overline{\mathbf{R}}_{12})))\leqslant 0.$

The last step is due to Corollary 2.5. $\blacksquare$

Our goal is to control $\nu(\exp(\kappa
N\mathcal{P}(\overline{\mathbf{R}}_{12})))$, by the idea of interpolation, we
just need to obtain a bound for $\nu_{0}(\exp(\kappa
N\mathcal{P}(\overline{\mathbf{R}}_{12})))$ since the interpolation path is
non-increasing by Corollary 2.6.

######  Lemma 2.7.

For $\kappa<1/(2\|\mathcal{V}\|)$, we have

$\displaystyle\nu_{0}(\exp(\kappa
N\mathcal{P}(\overline{\mathbf{R}}_{12})))\leqslant\det(I-2\kappa\mathcal{V})^{-1/2}.$
(35)

###### Proof of Lemma 2.7..

First, take a Gaussian vector
$\boldsymbol{u}=\sqrt{2\kappa}\Lambda^{\nicefrac{{1}}{{2}}}\mathcal{V}^{\nicefrac{{1}}{{2}}}\boldsymbol{\eta}$,
where $\boldsymbol{\eta}\sim\operatorname{N}(0,I_{m})$ and is independent with
the disorder and Gibbs randomness, then

$\displaystyle\nu_{0}(\exp(\kappa
N\mathcal{P}(\overline{\mathbf{R}}_{12}))=\operatorname{\mathds{E}}\nu_{0}(\exp(\sqrt{N}\boldsymbol{u}^{\intercal}\overline{\mathbf{R}}_{12}))$
$\displaystyle=\operatorname{\mathds{E}}\nu_{0}(\prod_{s=1}^{m}\exp(\sqrt{N}u_{s}(R^{s}-q_{s})))$
$\displaystyle=\operatorname{\mathds{E}}\nu_{0}\left(\prod_{s=1}^{m}\prod_{i\in
I_{s}}\exp\bigl{(}\sqrt{N}(\sigma_{i}^{1}\sigma_{i}^{2}-q_{s}){u_{s}}/{|I_{s}|}\bigr{)}\right)$
$\displaystyle\leqslant\operatorname{\mathds{E}}\prod_{s=1}^{m}\exp(Nu_{s}^{2}/2|I_{s}|)=\operatorname{\mathds{E}}\exp(\boldsymbol{u}^{\intercal}\Lambda^{-1}\boldsymbol{u}/2),$

where the last inequality is by the fact that (see [Tal11a]*Page 39 for a
proof)

$\nu_{0}(\exp
u(\sigma_{i}^{1}\sigma_{i}^{2}-q))=\exp(-qu)(\cosh(u)+q\sinh(u))\leqslant\exp(u^{2}/2).$

Finally we use the fact that
$\operatorname{\mathds{E}}\exp(\boldsymbol{u}^{\intercal}\Lambda^{-1}\boldsymbol{u}/2)=\operatorname{\mathds{E}}\exp(\kappa\boldsymbol{\eta}^{\intercal}\mathcal{V}\boldsymbol{\eta})=\det(I-2\kappa\mathcal{V})^{-\nicefrac{{1}}{{2}}}$
and this completes the proof. $\blacksquare$

Now we are ready to prove Theorem 1.7.

###  2.1. Proof of Theorem 1.7

Take $\rho=\eta+2\alpha\beta^{2}<1/(2\|\mathcal{V}\|)$, then

$\displaystyle\nu_{t}(\exp((\eta+2\alpha(1-t)\beta^{2})N\mathcal{P}(\overline{\mathbf{R}}_{12})))=\nu_{t}(\exp((\rho-2\alpha
t\beta^{2})N\mathcal{P}(\overline{\mathbf{R}}_{12})))$

By the non-increasing property of $\nu_{t}$ in Corollary 2.6, we have

$\nu_{t}(\exp((\rho-2\alpha
t\beta^{2})N\mathcal{P}(\overline{\mathbf{R}}_{12})))\leqslant\nu_{0}(\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))\quad\text{for}\quad t\in[0,1].$

By Lemma 2.7,

$\nu_{0}(\exp(\rho
N\mathcal{P}(\overline{\mathbf{R}}_{12})))\leqslant\det(I-2\rho\mathcal{V})^{-1/2}=\det(I-2(\eta+2\alpha\beta^{2})\mathcal{V})^{-1/2}.$

Thus,

$\nu_{1}(\exp(\eta
N\mathcal{P}(\overline{\mathbf{R}}_{12})))\leqslant\det(I-2(\eta+2\alpha\beta^{2})\mathcal{V})^{-1/2}.$

This proves Theorem 1.7. $\blacksquare$

##  3\. Asymptotic for the free energy

In this section, we give the proof of Theorem 1.9.

###  3.1. Proof of Theorem 1.9

First, let us prove the RS solution (19). Recall the formula (31) in Section
2,

$\phi^{\prime}(t)=\frac{\beta^{2}}{4}\bigl{(}\mathcal{Q}(\mathbf{1}-\mathbf{q})-\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{12}))\bigr{)}$

which gives the derivative of free energy associated to
$H_{N,t}(\boldsymbol{\sigma})$. By Theorem 1.7, when
$\beta^{2}<\beta_{0}^{2}$, we have
$\nu_{t}(\mathcal{Q}(\overline{\mathbf{R}}_{12}))\leqslant\frac{K}{N}$ by
Jensen’s inequality, which implies that
$\left|\phi^{\prime}(t)-\frac{\beta^{2}}{4}\mathcal{Q}(\bf
1-\mathbf{q})\right|\leqslant\frac{K}{N}$, combining with the fact that
$\phi(1)-\phi(0)=\int_{0}^{1}\phi^{\prime}(t)dt$, then

$\displaystyle\left|\phi(1)-\phi(0)-\frac{\beta^{2}}{4}\mathcal{Q}(\bf
1-\mathbf{q})\right|\leqslant\frac{K}{N}$

where $\phi(1)=F_{N}(\beta,h)$, and $\phi(0)$ is the free energy with
Hamiltonian $H_{N,0}(\boldsymbol{\sigma})$, a simple calculation based on the
definition of free energy can give us $\phi(0)=\log
2+\sum_{r=1}^{m}\lambda_{r}\operatorname{\mathds{E}}\log\cosh(\beta\eta\sqrt{Q^{r}}+h)$,
where $Q^{r}=(\Delta^{2}\Lambda\mathbf{q})_{r}$ and the expectation is w.r.t.
$\eta\sim\operatorname{N}(0,1)$.

Next, we prove the second part of Theorem 1.9, a quantitative Central Limit
Theorem for the free energy with $h\neq 0$. Let

$X_{t}=\frac{1}{\sqrt{N}}\bigl{(}\log Z_{N,t}-\operatorname{\mathds{E}}\log
Z_{N,0}-\frac{1}{4}N\beta^{2}t\cdot\mathcal{Q}(\bf 1-\mathbf{q})\bigr{)}$

where $Z_{N,t}$ is the partition function associated with the interpolated
Hamiltonian $H_{N,t}(\boldsymbol{\sigma})$ introduced in Section 2,
specifically

$Z_{N,t}=\sum_{\boldsymbol{\sigma}}\exp\biggl{(}\frac{\beta}{\sqrt{N}}\sum_{i<j}g_{ij}\sigma_{i}\sigma_{j}\cdot\sqrt{t}+\beta\sum_{r=1}^{m}\sqrt{Q^{r}}\sum_{i\in
I_{r}}\eta_{i}\sigma_{i}\cdot\sqrt{1-t}+h\sum_{i}\sigma_{i}\biggr{)}.$

Our goal is to prove a Central Limit Theorem for $t=1$. Let us first introduce
a lemma packing up technical details of the computation.

######  Lemma 3.1.

Let
$G_{t}(\boldsymbol{\sigma}):=\frac{\exp(H_{N,t}(\boldsymbol{\sigma}))}{Z_{N,t}},\sigma\in\\{1,-1\\}^{N}$
be the interpolated Gibbs measure. For any bounded twice differentiable
function $g\in C_{b}^{2}(\mathbb{R})$ and $t\in[0,1]$, we have

$\displaystyle\frac{d}{dt}\operatorname{\mathds{E}}g(X_{t})=\frac{\beta^{2}}{4}\operatorname{\mathds{E}}\left(g^{\prime\prime}(X_{t})\bigl{(}\langle\mathcal{Q}(\overline{\mathbf{R}}_{12})\rangle_{t}-\mathcal{Q}(\mathbf{q})-N^{-1}\operatorname{Tr}(\Delta^{2}\Lambda)\bigr{)}-\sqrt{N}g^{\prime}(X_{t})\langle\mathcal{Q}(\overline{\mathbf{R}}_{12})\rangle_{t}\right).$
(36)

In particular, for $m(t,\theta):=\operatorname{\mathds{E}}e^{\i\theta X_{t}}$
the characteristic function of $X_{t}$, we have

$\displaystyle\left|\frac{\partial}{\partial
t}m(t,\theta)-v^{2}\theta^{2}m(t,\theta)\right|\leqslant\frac{\beta^{2}C}{4}\left(\frac{|\theta|}{\sqrt{N}}+\frac{2\theta^{2}}{N}\right)\text{
for all }\theta\in\mathds{R}$ (37)

where $C$ is an upper bound of $\sup_{0\leqslant t\leqslant
1}N\nu_{t}(|\mathcal{Q}(\overline{\mathbf{R}}_{12})|)$ and
$\operatorname{Tr}(\Delta^{2}\Lambda)$; and
$v^{2}:=\beta^{2}\mathcal{Q}(\mathbf{q})/4$.

###### Proof of Lemma 3.1..

For $g\in C_{b}^{2}(\mathbb{R})$, we evaluate the LHS directly,

$\displaystyle\frac{d}{dt}\operatorname{\mathds{E}}g(X_{t})$
$\displaystyle=\frac{1}{\sqrt{N}}\operatorname{\mathds{E}}g^{\prime}(X_{t})\biggl{(}\frac{\beta}{2}\sum_{\boldsymbol{\sigma}}G_{t}(\boldsymbol{\sigma})\biggl{(}\frac{1}{\sqrt{tN}}\sum_{i<j}g_{ij}\sigma_{i}\sigma_{j}-\frac{1}{\sqrt{1-t}}\sum_{r=1}^{m}\sqrt{Q^{r}}\sum_{i\in
I_{r}}\eta_{i}\sigma_{i}\biggr{)}-\frac{N\beta^{2}}{4}\mathcal{Q}(\mathbf{1}-\mathbf{q})\biggr{)}.$

Notice that
$\frac{1}{\sqrt{tN}}\sum_{i<j}g_{ij}\sigma_{i}\sigma_{j}-\frac{1}{\sqrt{1-t}}\sum_{r=1}^{m}\sqrt{Q^{r}}\sum_{i\in
I_{r}}\eta_{i}\sigma_{i}$ and $X_{t}$ are both Gaussian random variables.
Applying Gaussian integration by parts, we get that the RHS is equal to

$\displaystyle\frac{\beta^{2}}{2N}\operatorname{\mathds{E}}\Bigg{[}g^{\prime\prime}(X_{s})\sum_{\boldsymbol{\sigma},\boldsymbol{\tau}}G_{t}(\boldsymbol{\sigma})G_{t}(\boldsymbol{\tau})\left(\frac{1}{2}(N\mathcal{Q}(\mathbf{R}_{12})-\sum_{r=1}^{m}\lambda_{r}\Delta_{r,r}^{2})-N\sum_{r,s}\lambda_{r}R_{12}^{r}\Delta_{r,s}^{2}\lambda_{s}q_{s}\right)\Bigg{]}$
$\displaystyle\qquad\qquad+\frac{\beta^{2}}{2\sqrt{N}}\operatorname{\mathds{E}}\Bigg{[}g^{\prime}(X_{s})\Bigg{(}\sum_{\boldsymbol{\sigma}}G_{t}(\boldsymbol{\sigma})\Big{(}\frac{1}{2}(N\mathcal{Q}(\mathbf{1})-\sum_{r=1}^{m}\lambda_{r}\Delta_{r,r}^{2})-N\sum_{r,s}\lambda_{r}\Delta_{r,s}^{2}\lambda_{s}q_{s}$
$\displaystyle\qquad\qquad-\sum_{\boldsymbol{\tau}}G_{t}(\boldsymbol{\tau})\Big{(}\frac{1}{2}(N\mathcal{Q}(\mathbf{R}_{12})-\sum_{r=1}^{m}\lambda_{r}\Delta_{r,r}^{2})-N\sum_{r,s}\lambda_{r}R_{12}^{r}\Delta_{r,s}^{2}\lambda_{s}q_{s}\Big{)}\Big{)}-\frac{N}{2}\mathcal{Q}(\mathbf{1}-\mathbf{q})\Bigg{)}\Bigg{]}$
$\displaystyle=\frac{\beta^{2}}{4}\operatorname{\mathds{E}}\left[g^{\prime\prime}(X_{t})\left(\langle\mathcal{Q}(\overline{\mathbf{R}}_{12})\rangle_{t}-\mathcal{Q}(\mathbf{q})-N^{-1}\operatorname{Tr}(\Delta^{2}\Lambda)\right)-\sqrt{N}g^{\prime}(X_{t})\langle\mathcal{Q}(\overline{\mathbf{R}}_{12})\rangle_{t}\right].$

For the second part, taking $g(x)=e^{\i x\theta}$ for fixed
$\theta\in\mathds{R}$, equation (36) can be rewritten as

$\displaystyle\frac{\partial}{\partial t}m(t,\theta)$
$\displaystyle=\frac{\beta^{2}\theta^{2}}{4}\mathcal{Q}(\mathbf{q})m(t,\theta)$
$\displaystyle\quad+\frac{\beta^{2}}{4}\operatorname{\mathds{E}}\left(\Bigl{(}-\frac{1}{N}\theta^{2}\left(N\langle\mathcal{Q}(\overline{\mathbf{R}}_{12})\rangle_{t}-\operatorname{Tr}(\Delta^{2}\Lambda)\right)-\i\sqrt{N}\theta\langle\mathcal{Q}(\overline{\mathbf{R}}_{12})\rangle_{t}\Bigr{)}e^{\i\theta
X_{t}}\right).$

From Section 2, we know that
$N\nu_{t}|\mathcal{Q}(\overline{\mathbf{R}}_{12})|\leqslant
N\nu_{t}(\mathcal{P}(\overline{\mathbf{R}}_{12}))$ is bounded, then

$\displaystyle\left|\frac{\partial}{\partial
t}m(t,\theta)-v^{2}\theta^{2}m(t,\theta)\right|\leqslant\frac{\beta^{2}C}{4}\left(\frac{|\theta|}{\sqrt{N}}+\frac{2\theta^{2}}{N}\right)\text{
for all }t\in[0,1],\theta\in\mathds{R}$

where
$\max\big{\\{}N\nu_{t}|\mathcal{Q}(\overline{\mathbf{R}}_{12})|,0\leqslant
t\leqslant 1,\operatorname{Tr}(\Delta^{2}\Lambda)\big{\\}}\leqslant C.$
$\blacksquare$

Now we prove the Central Limit Theorem for the free energy.

###### Proof of Theorem 1.9 part (ii).

Notice that

$\displaystyle\left|e^{-v^{2}\theta^{2}}m(1,\theta)-m(0,\theta)\right|$
$\displaystyle=\left|\int_{0}^{1}e^{-v^{2}\theta^{2}t}\left(\partial_{t}m(t,\theta)-v^{2}\theta^{2}m(t,\theta)\right)dt\right|$
$\displaystyle\leqslant\frac{\beta^{2}C}{4}\left(\frac{|\theta|}{\sqrt{N}}+\frac{2\theta^{2}}{N}\right)\int_{0}^{1}e^{-v^{2}\theta^{2}t}dt\leqslant\frac{\beta^{2}C}{4}\left(\frac{|\theta|}{\sqrt{N}}+\frac{2\theta^{2}}{N}\right).$
(38)

Thus

$\displaystyle\left|m(1,\theta)-m(0,\theta)e^{v^{2}\theta^{2}}\right|\leqslant\frac{\beta^{2}C}{4}\left(\frac{|\theta|}{\sqrt{N}}+\frac{2\theta^{2}}{N}\right)e^{v^{2}\theta^{2}}.$

For $\theta$ fixed, let $N\to\infty$. The RHS of the last step goes to 0.
Moreover, as $N\to\infty$, we have by the classical CLT $m(0,\theta)\to
e^{-c^{2}\theta^{2}/2}$ where

$c^{2}=\lim_{N\to\infty}N^{-1}\operatorname{Var}(\log
Z_{N,0})=\sum_{s=1}^{m}\lambda_{s}\operatorname{Var}(\log\cosh(\beta\eta\sqrt{Q^{s}}+h)).$

This completes the proof. $\blacksquare$

##  4\. Concentration of overlap with zero external field

In this section, we study the MSK model without an external field. We prove
that when $\beta<\beta_{c}$, the MSK model is in the replica symmetric phase.
Note that our argument holds when $\Delta^{2}$ is indefinite. In that case,
the RS regime is still given by $\beta<\beta_{c}$ as in $\Delta^{2}$
_positive-definite_ case.

The proof of Theorem 1.8 is based on a control of the free energy. Here $q=0$
is the unique solution to (15) and the averaged Gibbs measure under the
decoupled Hamiltonian is the product of i.i.d.
Bernoulli($\nicefrac{{1}}{{2}}$) measures and thus is non-random. Before
giving the proof, we introduce several useful lemmas following the proof for
the classical case. Recall that
$\mathcal{P}(\mathbf{x})=\mathbf{x}^{\intercal}\Lambda^{\nicefrac{{1}}{{2}}}\mathcal{V}\Lambda^{\nicefrac{{1}}{{2}}}\mathbf{x}$
and
$\mathcal{V}=|\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}|$.

######  Lemma 4.1.

If $\gamma+\beta^{2}<\beta_{c}^{2}$, then

$\displaystyle\operatorname{\mathds{E}}\sum_{\boldsymbol{\sigma}_{1},\boldsymbol{\sigma}_{2}}\exp\left(H_{N}(\boldsymbol{\sigma}_{1})+H_{N}(\boldsymbol{\sigma}_{2})+\frac{\gamma
N}{2}\mathcal{P}(\mathbf{R}_{12})\right)\leqslant(\operatorname{\mathds{E}}Z_{N})^{2}\cdot\det(I-(\beta^{2}+\gamma)\mathcal{V})^{-\nicefrac{{1}}{{2}}}.$
(39)

###### Proof of Lemma 4.1.

Since the Hamiltonian $H_{N}{(\boldsymbol{\sigma})}$ is a _centered_ Gaussian
field, we know that

$\displaystyle\begin{split}\operatorname{\mathds{E}}\sum_{\boldsymbol{\sigma}_{1},\boldsymbol{\sigma}_{2}}&\exp\left(H_{N}(\boldsymbol{\sigma}_{1})+H_{N}(\boldsymbol{\sigma}_{2})+\frac{\gamma
N}{2}\mathcal{P}(\mathbf{R}_{12})\right)\\\
&=\sum_{\boldsymbol{\sigma}_{1},\boldsymbol{\sigma}_{2}}\exp\left(\frac{1}{2}\operatorname{\mathds{E}}(H_{N}(\boldsymbol{\sigma}_{1})+H_{N}(\boldsymbol{\sigma}_{2}))^{2}+\frac{\gamma
N}{2}\mathcal{P}(\mathbf{R}_{12})\right)\\\
&\leqslant\sum_{\boldsymbol{\sigma}_{1},\boldsymbol{\sigma}_{2}}\exp\left(\frac{N\beta^{2}}{2}\left(\sum_{s,t}\lambda_{s}\Delta_{s,t}^{2}\lambda_{t}-\frac{2}{N}\sum_{s}\lambda_{s}\Delta_{s,s}^{2}\right)+\frac{(\beta^{2}+\gamma)N\mathcal{P}(\mathbf{R}_{12})}{2}\right)\\\
&\leqslant(\operatorname{\mathds{E}}Z_{N})^{2}\cdot
2^{-2N}\sum_{\boldsymbol{\sigma}_{1},\boldsymbol{\sigma}_{2}}\exp\left(\frac{\delta
N}{2}\mathcal{P}(\mathbf{R}_{12})\right)\end{split}$ (40)

where $\delta=\beta^{2}+\gamma$, and the last two inequalities are based on
the fact that

$\operatorname{\mathds{E}}Z_{N}=2^{N}\exp\left(\frac{N\beta^{2}}{4}\left(\sum_{s,t}\lambda_{s}\Delta_{s,t}^{2}\lambda_{t}-\frac{1}{N}\sum_{s}\lambda_{s}\Delta_{s,s}^{2}\right)\right),$

and

$\displaystyle\operatorname{\mathds{E}}H_{N}(\boldsymbol{\sigma})^{2}$
$\displaystyle=\frac{N\beta^{2}}{2}\left(\sum_{s,t}\lambda_{s}\Delta_{s,t}^{2}\lambda_{t}-\frac{1}{N}\sum_{s}\lambda_{s}\Delta_{s,s}^{2}\right)$
$\displaystyle\operatorname{\mathds{E}}H_{N}(\boldsymbol{\sigma}_{1})H_{N}(\boldsymbol{\sigma}_{2})$
$\displaystyle=\frac{N\beta^{2}}{2}\left(\sum_{s,t}\lambda_{s}R^{s}\Delta_{s,t}^{2}\lambda_{t}R^{t}-\frac{1}{N}\sum_{s}\lambda_{s}\Delta_{s,s}^{2}\right).$

Next we prove that, for fixed $\boldsymbol{\sigma}_{2}$,

$\displaystyle 2^{-N}\sum_{\boldsymbol{\sigma}_{1}}\exp\left(\frac{\delta
N}{2}\mathcal{P}(\mathbf{R}_{12})\right)=2^{-N}\sum_{\boldsymbol{\sigma}}\exp\left(\frac{\delta
N}{2}\mathcal{P}(\mathbf{v})\right)$ (41)

where $\mathbf{v}=(\frac{1}{|I_{1}|}\sum_{i\in
I_{1}}\sigma_{i},\frac{1}{|I_{2}|}\sum_{i\in
I_{2}}\sigma_{i},\ldots,\frac{1}{|I_{m}|}\sum_{i\in I_{m}}\sigma_{i})$. Note
that for $\boldsymbol{\xi}=(\xi_{1},\xi_{2},\ldots,\xi_{m})$, where
$\xi_{s}=\frac{1}{|I_{s}|}\sum_{i\in I_{s}}\eta_{i}$, and
$\eta_{i}\sim\text{Bernouli}(\nicefrac{{1}}{{2}})$ are i.i.d, we have

$\displaystyle 2^{-N}\sum_{\boldsymbol{\sigma}}\exp\left(\frac{\delta
N}{2}\mathcal{P}(\mathbf{v})\right)$
$\displaystyle=\operatorname{\mathds{E}}\exp\left(\frac{\delta
N}{2}\mathcal{P}(\boldsymbol{\xi})\right)$
$\displaystyle=\operatorname{\mathds{E}}\exp\left(\frac{\delta
N}{2}\mathcal{Q}(\boldsymbol{\xi},\Lambda^{\nicefrac{{1}}{{2}}}\mathcal{V}\Lambda^{\nicefrac{{1}}{{2}}})\right)=\operatorname{\mathds{E}}\operatorname{\mathds{E}}_{\mathbf{u}}\exp(\sqrt{N}\mathbf{u}\boldsymbol{\xi})$

where
$\mathbf{u}\sim\operatorname{N}(0,\delta\Lambda^{\nicefrac{{1}}{{2}}}\mathcal{V}\Lambda^{\nicefrac{{1}}{{2}}})$.
By a similar way as in the proof of Lemma 2.7, we have

$\displaystyle\operatorname{\mathds{E}}\operatorname{\mathds{E}}_{\mathbf{u}}\exp(\sqrt{N}\mathbf{u}\boldsymbol{\xi})$
$\displaystyle=\operatorname{\mathds{E}}_{\mathbf{u}}\operatorname{\mathds{E}}\prod_{s=1}^{m}\prod_{i\in
I_{s}}\exp(\sqrt{N}u_{s}\eta_{i}/|I_{s}|)$
$\displaystyle\leqslant\operatorname{\mathds{E}}_{\mathbf{u}}\prod_{s=1}^{m}\exp(Nu_{s}^{2}/2|I_{s}|)=\operatorname{\mathds{E}}_{\mathbf{u}}\exp(\mathbf{u}^{\intercal}\Lambda^{-1}\mathbf{u}/2)=\det(I-\delta\mathcal{V})^{-\nicefrac{{1}}{{2}}}.$

Going back to (40) and using $\delta=\beta^{2}+\gamma$, we finish the proof.
$\blacksquare$

######  Theorem 4.2.

If $\beta<\beta_{c}$, then there exists a constant $K$ such that for each
$N>0$, and each $t>0$,

$\displaystyle\operatorname{\mathds{P}}\left(\log Z_{N}\leqslant N(\log
2+{\beta^{2}}/{4}\cdot\mathcal{Q}(\mathbf{1}))-t\right)\leqslant
K\exp(-t^{2}/K).$ (42)

The following lemma is going to tell us how to choose the $K$ appropriately.

######  Lemma 4.3.

If $\beta<\beta_{c}$, we have

$\displaystyle\operatorname{\mathds{P}}\left(Z_{N}\geqslant\operatorname{\mathds{E}}Z_{N}/2,N\langle\mathcal{Q}(\mathbf{R}_{12})\rangle\leqslant
K\right)\geqslant\frac{1}{K}.$ (43)

###### Proof of Lemma 4.3.

Notice that

$Z_{N}^{2}\langle\exp(\frac{\gamma
N}{2}\mathcal{Q}(\mathbf{R}_{12}))\rangle=\sum_{\boldsymbol{\sigma}_{1},\boldsymbol{\sigma}_{2}}\exp\left(H_{N}(\boldsymbol{\sigma}_{1})+H_{N}(\boldsymbol{\sigma}_{2})+\frac{\gamma
N}{2}\mathcal{Q}(\mathbf{R}_{12})\right).$

Thus by Lemma 4.1, we have

$\displaystyle\operatorname{\mathds{E}}\left(Z_{N}^{2}\langle\exp(\frac{\gamma
N}{2}\mathcal{Q}(\mathbf{R}_{12}))\rangle\right)\leqslant(\operatorname{\mathds{E}}Z_{N})^{2}\det(I-(\beta^{2}+\gamma)\mathcal{V})^{-\nicefrac{{1}}{{2}}}$

and

$\operatorname{\mathds{P}}(Z_{N}\geqslant\operatorname{\mathds{E}}Z_{N}/2)\geqslant(\operatorname{\mathds{E}}Z_{N})^{2}/4\operatorname{\mathds{E}}Z_{N}^{2}\geqslant\det(I-\beta^{2}\mathcal{V})^{\nicefrac{{1}}{{2}}}/4.$

By Markov’s inequality,

$\displaystyle\operatorname{\mathds{P}}\left(Z_{N}^{2}\langle\exp(\frac{\gamma
N}{2}\mathcal{Q}(\mathbf{R}_{12}))\rangle\leqslant
t(\operatorname{\mathds{E}}Z_{N})^{2}\right)\geqslant
1-\frac{1}{t}\det(I-(\beta^{2}+\gamma)\mathcal{V})^{-\nicefrac{{1}}{{2}}}$

and

$\displaystyle\operatorname{\mathds{P}}\bigl{(}Z_{N}^{2}$
$\displaystyle\langle\exp(\frac{\gamma
N}{2}\mathcal{Q}(\mathbf{R}_{12}))\rangle\leqslant
t(\operatorname{\mathds{E}}Z_{N})^{2},Z_{N}\geqslant\operatorname{\mathds{E}}Z_{N}/2\bigr{)}$
$\displaystyle\geqslant\operatorname{\mathds{P}}\left(Z_{N}^{2}\langle\exp(\frac{\gamma
N}{2}\mathcal{Q}(\mathbf{R}_{12}))\rangle\leqslant
t(\operatorname{\mathds{E}}Z_{N})^{2}\right)+\operatorname{\mathds{P}}\left(Z_{N}\geqslant\operatorname{\mathds{E}}Z_{N}/{2}\right)-1$
$\displaystyle\geqslant\operatorname{\mathds{P}}\left(Z_{N}\geqslant\operatorname{\mathds{E}}Z_{N}/{2}\right)-\frac{1}{t}\det(I-(\beta^{2}+\gamma)\mathcal{V})^{-\nicefrac{{1}}{{2}}}.$

Based on this inequality, we can choose some large $K$ such that,

$\displaystyle\operatorname{\mathds{P}}\left(Z_{N}^{2}\langle\exp(\frac{\gamma
N}{2}\mathcal{Q}(\mathbf{R}_{12}))\rangle\leqslant
K(\operatorname{\mathds{E}}Z_{N})^{2},Z_{N}\geqslant\operatorname{\mathds{E}}Z_{N}/2\right)\geqslant\frac{1}{K}.$

When $Z_{N}\geqslant\operatorname{\mathds{E}}Z_{N}/{2}$,
$Z_{N}^{2}\langle\exp(\frac{\gamma
N}{2}\mathcal{Q}(\mathbf{R}_{12}))\rangle\leqslant
K(\operatorname{\mathds{E}}Z_{N})^{2}$ implies that

$\left\langle\exp(\frac{\gamma
N}{2}\mathcal{Q}(\mathbf{R}_{12}))\right\rangle\leqslant 4K$

and hence $\langle N\mathcal{Q}(\mathbf{R}_{12}))\rangle\leqslant K$. The
constant $K$ in different expressions are different, here we abuse the
notation a bit. $\blacksquare$

######  Lemma 4.4.

The following identity:

$\displaystyle
N^{2}\langle\mathcal{Q}(\mathbf{R}_{12})\rangle=\sum_{s,t}\sum_{i\in
I_{s}}\sum_{j\in
I_{t}}\Delta_{s,t}^{2}\langle\sigma_{i}\sigma_{j}\rangle^{2}\geqslant 0$ (44)

holds.

###### Proof of Lemma 4.4.

Recall that $\mathcal{Q}(\mathbf{R}_{12})$ is the quadratic form of
$\mathbf{R}_{12}$ with matrix $\Lambda\Delta^{2}\Lambda$, then we have

$\displaystyle N^{2}\mathcal{Q}(\mathbf{R}_{12})$
$\displaystyle=N^{2}\sum_{s,t}R_{12}^{s}\lambda_{s}\Delta_{s,t}^{2}\lambda_{t}R_{12}^{t}$
$\displaystyle=N^{2}\sum_{s,t}\frac{1}{|I_{s}|}\sum_{i\in
I_{s}}\sigma_{i}^{1}\sigma_{i}^{2}\lambda_{s}\Delta_{s,t}^{2}\lambda_{t}\frac{1}{|I_{t}|}\sum_{j\in
I_{t}}\sigma_{j}^{1}\sigma_{j}^{2}=\sum_{s,t}\sum_{i\in
I_{s}}\sigma_{i}^{1}\sigma_{i}^{2}\Delta_{s,t}^{2}\sum_{j\in
I_{t}}\sigma_{j}^{1}\sigma_{j}^{2},$

After applying the Gibbs average $\langle\cdot\rangle$ on both sides, for RHS
we have
$\langle\sigma_{i}^{1}\sigma_{j}^{1}\Delta_{s,t}^{2}\sigma_{i}^{2}\sigma_{j}^{2}\rangle=\Delta_{s,t}^{2}\langle\sigma_{i}\sigma_{j}\rangle^{2}$,
then $N^{2}\ \langle\mathcal{Q}(\mathbf{R}_{12})\rangle=\sum_{s,t}\sum_{i\in
I_{s}}\sum_{j\in
I_{t}}\Delta_{s,t}^{2}\langle\sigma_{i}\sigma_{j}\rangle^{2}$. $\blacksquare$

######  Lemma 4.5 ([Tal03]*Lemma 2.2.11).

Consider a closed subset $B$ of $\mathbb{R}^{M}$, and set

$d(\mathbf{x},B)=\inf\\{d(\mathbf{x},\mathbf{y});\mathbf{y}\in B\\},$

as the Euclidean distance from $\mathbf{x}$ to $B$. Then for $t>0$, we have

$\displaystyle\operatorname{\mathds{P}}\left(d(\boldsymbol{\eta},B)\geqslant
t+2\sqrt{\log\frac{2}{\operatorname{\mathds{P}}(\boldsymbol{\eta}\in
B)}}\right)\leqslant 2\exp(-t^{2}/4)$ (45)

where $\boldsymbol{\eta}=(\eta_{1},\eta_{2},\ldots,\eta_{m})$ and
$\eta_{i}\sim\operatorname{N}(0,1)$ are i.i.d.

Now we are ready to prove the Theorem 4.2

###### Proof of Theorem 4.2.

Recall the Hamiltonian with $h=0$ is

$H_{N}(\boldsymbol{\sigma})=\frac{\beta}{\sqrt{N}}\sum_{i<j}g_{ij}\sigma_{i}\sigma_{j}$

where $\operatorname{\mathds{E}}g_{ij}^{2}=\Delta^{2}_{s,t}$ for $i\in
I_{s},j\in I_{t}$. In the following content, we let
$g_{ij}=\Delta_{s,t}\eta_{ij}$, where $\eta_{ij}\sim\operatorname{N}(0,1)$.
Let $M=N(N-1)/2$, and consider Gaussian $\boldsymbol{\eta}=(\eta_{ij})_{i<j}$,
in this case, $\boldsymbol{\eta}\in\mathbb{R}^{M}$. We understand
$H_{N}(\boldsymbol{\sigma}),Z_{N}$ as functions of $\boldsymbol{\eta}$. By
Lemma 4.3, for some suitably large $K=K_{1}$, there is a subset
$B\subset\mathbb{R}^{M}$ with

$\\{\boldsymbol{\eta}\in
B\\}=\\{Z_{N}(\boldsymbol{\eta})\geqslant\operatorname{\mathds{E}}Z_{N}/2;N\langle\mathcal{Q}(\mathbf{R}_{12})\rangle\leqslant
K_{1}\\}$

i.e. the set $B$ characterizes the event on RHS, and also
$\operatorname{\mathds{P}}(\boldsymbol{\eta}\in B)\geqslant\frac{1}{K_{1}}$.
Next we will prove

$\displaystyle\log Z_{N}(\boldsymbol{\eta})\geqslant N\left(\log
2+\frac{\beta^{2}}{4}\mathcal{Q}(\bf
1)\right)-K_{1}(1+d(\boldsymbol{\eta},B))$ (46)

For $\boldsymbol{\eta}^{\prime}\in B$, we have

$\log
Z_{N}(\boldsymbol{\eta}^{\prime})\geqslant\log\frac{1}{2}\operatorname{\mathds{E}}Z_{N}=\frac{\beta^{2}}{4}\left(N\mathcal{Q}(\mathbf{1})-\sum_{s=1}^{m}\lambda_{s}\Delta_{s,s}^{2}\right)+(N-1)\log
2.$

To prove (46), it is enough to show that

$\displaystyle\log Z_{N}(\boldsymbol{\eta})\geqslant\log
Z_{N}(\boldsymbol{\eta}^{\prime})-K_{1}d(\boldsymbol{\eta},\boldsymbol{\eta}^{\prime})$
(47)

for all $\boldsymbol{\eta}^{\prime}\in B$. Here $K_{1}$ must be bigger than
$\log 2+\frac{\beta^{2}}{4}\sum_{s=1}^{m}\lambda_{s}\Delta_{s,s}^{2}$. Notice
that

$\displaystyle
Z_{N}(\boldsymbol{\eta})=Z_{N}(\boldsymbol{\eta}^{\prime})\left\langle\exp\biggl{(}\frac{\beta}{\sqrt{N}}\sum_{s,t}\sum_{\begin{subarray}{c}i\in
I_{s},j\in I_{t},\\\
i<j\end{subarray}}\Delta_{s,t}(\boldsymbol{\eta}_{ij}-\boldsymbol{\eta}_{ij}^{\prime})\sigma_{i}\sigma_{j}\biggr{)}\right\rangle^{\prime}$

where $\langle\cdot\rangle^{\prime}$ denotes the Gibbs average with disorders
$\boldsymbol{\eta}^{\prime}$. Since there is an exponential part, it’s natural
to apply Jensen’s inequality,

$\left\langle\exp\biggl{(}\frac{\beta}{\sqrt{N}}\sum_{s,t}\sum_{\begin{subarray}{c}i\in
I_{s},j\in I_{t},\\\
i<j\end{subarray}}\Delta_{s,t}(\boldsymbol{\eta}_{ij}-\boldsymbol{\eta}_{ij}^{\prime})\sigma_{i}\sigma_{j}\biggr{)}\right\rangle^{\prime}\geqslant\exp\biggl{(}\frac{\beta}{\sqrt{N}}\sum_{s,t}\sum_{\begin{subarray}{c}i\in
I_{s},j\in I_{t},\\\
i<j\end{subarray}}\Delta_{s,t}(\boldsymbol{\eta}_{ij}-\boldsymbol{\eta}_{ij}^{\prime})\langle\sigma_{i}\sigma_{j}\rangle^{\prime}\biggr{)}$

where we use $g_{ij}=\Delta_{s,t}\eta_{ij}$ to rewrite the Hamiltonian, then
applying Cauchy-Schwarz inequality and by Lemma 4.4, we have

$\sum_{s,t}\sum_{\begin{subarray}{c}i\in I_{s},j\in I_{t},\\\
i<j\end{subarray}}\Delta_{s,t}(\boldsymbol{\eta}_{ij}-\boldsymbol{\eta}_{ij}^{\prime})\langle\sigma_{i}\sigma_{j}\rangle^{\prime}\geqslant-d(\boldsymbol{\eta},\boldsymbol{\eta}^{\prime})(N^{2}\langle\mathcal{Q}(\mathbf{R}_{12})\rangle^{\prime})^{\nicefrac{{1}}{{2}}}\geqslant-
K_{1}\sqrt{N}d(\boldsymbol{\eta},\boldsymbol{\eta}^{\prime})$

where
$d(\boldsymbol{\eta},\boldsymbol{\eta}^{\prime})=\sum_{i<j}(\eta_{ij}-\eta_{ij}^{\prime})^{2}$.
The last inequality is based on $\boldsymbol{\eta}^{\prime}\in B$, then this
proves (47), and hence (46). From that, we have

$\log Z_{N}(\boldsymbol{\eta})\leqslant
N\left(\frac{\beta^{2}}{4}\mathcal{Q}(\mathbf{1})+\log 2\right)-t\Rightarrow
d(\boldsymbol{\eta},B)\geqslant\frac{t-K_{1}}{K_{1}}\geqslant
2\sqrt{\log\frac{2}{\operatorname{\mathds{P}}(\boldsymbol{\eta}\in
B)}}+\frac{t-K_{2}}{K_{1}}$

then by (45) in Lemma 4.5, it follows

$\displaystyle\operatorname{\mathds{P}}\left(\log
Z_{N}(\boldsymbol{\eta})\leqslant
N\left(\frac{\beta^{2}}{4}\mathcal{Q}(\mathbf{1})+\log
2\right)-t\right)\leqslant
2\exp\left(-\frac{(t-K_{2})^{2}}{4K_{1}^{2}}\right)$

for $t\geqslant K_{2}$. Then for the RHS, if we take $K$ large enough, when
$t\geqslant K_{2}$,

$2\exp\left(-\frac{(t-K_{2})^{2}}{4K_{1}^{2}}\right)\leqslant
K\exp\left(-\frac{t^{2}}{K}\right).$

On the other hand, when $0\leqslant t\leqslant K_{2}$, we have
$K\exp\left(-\frac{t^{2}}{K}\right)\geqslant 1$. Therefore, in any case, we
proved

$\operatorname{\mathds{P}}\left(\log Z_{N}(\boldsymbol{\eta})\leqslant
N\left(\frac{\beta^{2}}{4}\mathcal{Q}(\mathbf{1})+\log
2\right)-t\right)\leqslant K\exp\left(-\frac{t^{2}}{K}\right).$

$\blacksquare$

Now we collect all the above to prove the main theorem.

###  4.1. Proof of Theorem 1.8

By Cauchy-schwarz w.r.t the measure $\langle\cdot\rangle$

$\displaystyle\left\langle\exp\left(\frac{\beta_{c}^{2}-\beta^{2}}{8}N\mathcal{P}(\mathbf{R}_{12})\right)\right\rangle$
$\displaystyle\leqslant\left\langle\exp\left(\frac{\beta_{c}^{2}-\beta^{2}}{4}N\mathcal{P}(\mathbf{R}_{12})\right)\right\rangle^{\nicefrac{{1}}{{2}}}$
$\displaystyle=\frac{1}{Z_{N}}\left(Z_{N}^{2}\left\langle\exp\left(\frac{\beta_{c}^{2}-\beta^{2}}{4}N\mathcal{P}(\mathbf{R}_{12})\right)\right\rangle\right)^{\nicefrac{{1}}{{2}}}$

then apply Cauchy-schwarz again for $\operatorname{\mathds{E}}$,

$\displaystyle\operatorname{\mathds{E}}\left\langle\exp\left(\frac{\beta_{c}^{2}-\beta^{2}}{8}N\mathcal{P}(\mathbf{R}_{12})\right)\right\rangle$
$\displaystyle\leqslant\left(\operatorname{\mathds{E}}\frac{1}{Z_{N}^{2}}\right)^{\nicefrac{{1}}{{2}}}\left[\operatorname{\mathds{E}}\left(Z_{N}^{2}\left\langle\exp\left(\frac{\beta_{c}^{2}-\beta^{2}}{4}N\mathcal{P}(\mathbf{R}_{12})\right)\right\rangle\right)\right]^{\nicefrac{{1}}{{2}}}$
$\displaystyle\leqslant\left(\operatorname{\mathds{E}}\frac{1}{Z_{N}^{2}}\right)^{\nicefrac{{1}}{{2}}}(K\operatorname{\mathds{E}}(Z_{N})^{2})^{\nicefrac{{1}}{{2}}}=K\left(\operatorname{\mathds{E}}\frac{(\operatorname{\mathds{E}}Z_{N})^{2}}{Z_{N}^{2}}\right)^{\nicefrac{{1}}{{2}}}$

where the second inequality is due to Lemma 4.1. By Theorem 4.2,

$\displaystyle\operatorname{\mathds{P}}\left(\frac{\operatorname{\mathds{E}}Z_{N}}{Z_{N}}>t\right)\leqslant
K\exp(-(\log t)^{2}/K)$ (48)

then we have
$\operatorname{\mathds{E}}\frac{(\operatorname{\mathds{E}}Z_{N})^{2}}{Z_{N}^{2}}<K$
by standard tail estimates. $\blacksquare$

##  5\. MSK cavity solution

For the classical SK model, the cavity method is an induction on the number of
spins. In this section, we generalize the cavity method to the MSK model and
derive a linear system for the variance-covariance matrices of the overlap
vectors. More specifically, to compare the system with $N$ spins and $(N-1)$
spins, we need to choose a spin to be decoupled from others.

Since the classical SK model is the one-species case of the MSK model, the
spin can be chosen uniformly. For the multi-species case, decoupling a spin
now depends on the structure of the species. However, we can choose a species
$s\in\\{1,2,\ldots,m\\}$ first with probability ${\lambda_{s}}$, then select a
spin uniformly inside that species to be decoupled, then compare this
decoupled system and the original one. For convenience, we denote the released
spin as $\varepsilon:=\sigma_{s_{\star}}$, where $s_{\star}$ is the index of
last spin in species $s$ for a configuration $\boldsymbol{\sigma}$. Once the
species is determined, we take the convention to decouple the last spin in the
chosen species because of symmetry inside a particular species. Note that this
is equivalent to choosing a species, u.a.r. from all spins.

Let $\delta$ represent the species to be chosen, which is a random variable
taking value in $\\{1,2,\ldots,m\\}$ with probability
$\mathbb{P}(\delta=s)=\lambda_{s}$. Given $\delta=\tilde{s}$, we decouple the
last spin $\sigma_{\tilde{s}_{\star}}$ in that species, then the interpolated
Hamiltonian between the original system and the decoupled system is

$\displaystyle
H_{t}^{\tilde{s}}(\boldsymbol{\sigma})=H_{\tilde{s}_{\star}^{-}}(\boldsymbol{\rho})+\sigma_{\tilde{s}_{\star}}\biggl{(}\sqrt{t}\cdot\frac{\beta}{\sqrt{N}}\sum_{i\neq\tilde{s}_{\star}}g_{i\tilde{s}_{\star}}\sigma_{i}+\sqrt{1-t}\cdot\beta\eta\sqrt{Q^{\tilde{s}}}+h\biggr{)}$
(49)

where
$H_{\tilde{s}_{\star}^{-}}(\boldsymbol{\rho})=\frac{1}{\sqrt{N}}\sum_{i<j}g_{ij}\rho_{i}\rho_{j}+\sum_{i\neq\tilde{s}_{\star}}\rho_{i}$,
and $\boldsymbol{\rho}$ is a ($N-1$)-dimensional vector of spins in
$\boldsymbol{\sigma}$ without $\sigma_{\tilde{s}_{\star}}$.

Recall that the overlap vector between two replicas
$\sigma^{l},\sigma^{l^{\prime}}$ in MSK model is defined as

$\displaystyle\mathbf{R}_{ll^{\prime}}=(R_{ll^{\prime}}^{(1)},R_{ll^{\prime}}^{(2)},\cdots,R_{ll^{\prime}}^{(m)})^{\intercal}$
(50)

where each coordinate in the above vector represents the marginal overlap in
the corresponding species. Similarly,

$\displaystyle\mathbf{R}_{ll^{\prime}}^{(\tilde{s}-)}=(R_{ll^{\prime}}^{(1)},R_{ll^{\prime}}^{(2)},\cdots,R_{ll^{\prime}}^{(\tilde{s}-)},\cdots
R_{ll^{\prime}}^{(m)})^{\intercal}$ (51)

where
$R_{ll^{\prime}}^{(\tilde{s}-)}=R_{ll^{\prime}}^{(\tilde{s})}-\frac{1}{|I_{\delta}|}\varepsilon_{l}\varepsilon_{l^{\prime}}$,
and $\varepsilon_{l}=\boldsymbol{\sigma}^{l}_{\tilde{s}_{\star}}$ denotes the
released spin in $\boldsymbol{\sigma}^{l}$.

In the following, $\nu_{t,\tilde{s}}(\cdot)$ denotes the expectation over the
Gibbs measure and disorder $g_{ij}$ associated with the Hamiltonian (49). In
particular, $\nu_{1,\tilde{s}}(\cdot)$ corresponds to the original system and
does not depend on $\tilde{s}$, but $\nu_{0,\tilde{s}}$ and
$\nu_{t,\tilde{s}}$ depend on $\tilde{s}$ since they both involve the
decoupling procedure. In the rest of this section, we state the results for
each fixed $\delta=\tilde{s}$, unless we state clearly otherwise.

######  Lemma 5.1.

For any $f^{-}$ on $\Sigma_{N-1}^{n}$, where $\Sigma_{N-1}$ does not contain
the spin $\sigma_{\tilde{s}_{*}}$, and $I\subset\\{1,2,\cdots,n\\}$, we have

$\displaystyle\nu_{0,\tilde{s}}\biggl{(}f^{-}\prod_{i\in
I}\varepsilon_{i}\biggr{)}=\nu_{0,\tilde{s}}(f^{-})\cdot\operatorname{\mathds{E}}(\tanh
Y^{\tilde{s}})^{|I|}$ (52)

where $Y^{\tilde{s}}=\beta\eta\sqrt{Q^{\tilde{s}}}+h$, and $|I|$ is the
cardinality of the set $I$.

The proof of this lemma is same to the classical case, see Section 1.4 in
[Tal11a]. Next we turn to computation of
$\nu_{t,\tilde{s}}^{\prime}(f):=\frac{d}{dt}\nu_{t,\tilde{s}}(f)$ in cavity
method. Let

$\displaystyle
u_{\boldsymbol{\sigma}}=\frac{\beta}{\sqrt{N}}\sigma_{\tilde{s}_{\star}}\sum_{i\neq\tilde{s}_{\star}}g_{i\tilde{s}_{\star}}\sigma_{i},\qquad
v_{\boldsymbol{\sigma}}=\beta\eta\sigma_{\tilde{s}_{\star}}\sqrt{Q^{\tilde{s}}},\qquad
w_{\boldsymbol{\sigma}}=\exp(H_{\tilde{s}_{\star}^{-}}(\boldsymbol{\rho})+h\sigma_{\tilde{s}_{\star}})$

in (49). Recall that in Section 2, the derivative $\nu_{t}^{\prime}(f)$ was
computed in (32) for some $f$ on $\Sigma_{N}^{n}$:

$\displaystyle\nu_{t}^{\prime}(f)=2\biggl{(}\sum_{1\leqslant
l<l^{\prime}\leqslant
n}\nu_{t}(U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{l^{\prime}})f)-n\sum_{l\leqslant
n}\nu_{t}(U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{n+1})f)+\frac{n(n+1)}{2}\nu_{t}(U(\boldsymbol{\sigma}^{n+1},\boldsymbol{\sigma}^{n+2})f)\biggr{)}$

where
$U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{l^{\prime}})=\frac{1}{2}(\operatorname{\mathds{E}}u_{\boldsymbol{\sigma}^{l}}u_{\boldsymbol{\sigma}^{l^{\prime}}}-\operatorname{\mathds{E}}v_{\boldsymbol{\sigma}^{l}}v_{\boldsymbol{\sigma}^{l^{\prime}}})$.
In the setting of cavity method,

$\displaystyle\operatorname{\mathds{E}}u_{\boldsymbol{\sigma}^{l}}u_{\boldsymbol{\sigma}^{l^{\prime}}}$
$\displaystyle=\beta^{2}\varepsilon_{l}\varepsilon_{l^{\prime}}\biggl{(}\sum_{s\neq\tilde{s}}\Delta_{s,\tilde{s}}^{2}\lambda_{s}R_{ll^{\prime}}^{(s)}+\Delta_{\tilde{s},\tilde{s}}^{2}\lambda_{\tilde{s}}R_{ll^{\prime}}^{(\tilde{s}-)}\biggr{)}$
$\displaystyle\operatorname{\mathds{E}}v_{\boldsymbol{\sigma}^{l}}v_{\boldsymbol{\sigma}^{l^{\prime}}}$
$\displaystyle=\beta^{2}\varepsilon_{l}\varepsilon_{l^{\prime}}Q^{\tilde{s}}$

then

$\displaystyle U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{l^{\prime}})$
$\displaystyle=\frac{1}{2}\beta^{2}\varepsilon_{l}\varepsilon_{l^{\prime}}\biggl{(}\sum_{s\neq\tilde{s}}\Delta_{s,\tilde{s}}^{2}\lambda_{s}R_{ll^{\prime}}^{(s)}+\Delta_{\tilde{s},\tilde{s}}^{2}\lambda_{\tilde{s}}R_{ll^{\prime}}^{(\tilde{s}-)}-Q^{\tilde{s}}\biggr{)}=\frac{1}{2}\beta^{2}\varepsilon_{l}\varepsilon_{l^{\prime}}\cdot\Delta_{\tilde{s}}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(\tilde{s}-)}.$

We make the inner product of vectors and the dependence of
$U(\boldsymbol{\sigma}^{l},\boldsymbol{\sigma}^{l^{\prime}})$ on $\tilde{s}$
implicit in the above expression. We use
$\Delta_{s}^{2}:=\Delta^{2}_{s,\cdot}$ to denote the $s$-th row vector of
$\Delta^{2}$, and we keep using this notation for other symmetric matrices $A$
in the rest of this section. Thus in cavity solution
$\nu^{\prime}_{t,\tilde{s}}(f)$ can be written in the following way:

######  Theorem 5.2.

For $f$ on $\Sigma_{N}^{n}$, we have

$\displaystyle\begin{split}\nu^{\prime}_{t,\tilde{s}}(f)&=\beta^{2}\biggl{(}\sum_{1\leqslant
l<l^{\prime}\leqslant
n}\nu_{t,\tilde{s}}(\varepsilon_{l}\varepsilon_{l^{\prime}}f\cdot\Delta_{\tilde{s}}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(\tilde{s}-)})-n\sum_{l\leqslant
n}\nu_{t,\tilde{s}}(\varepsilon_{l}\varepsilon_{n+1}f\cdot\Delta_{\tilde{s}}^{2}\Lambda\overline{\mathbf{R}}_{l,n+1}^{(\tilde{s}-)})\\\
&\qquad\qquad\qquad+\frac{n(n+1)}{2}\nu_{t,\tilde{s}}(\varepsilon_{n+1}\varepsilon_{n+2}f\cdot\Delta_{\tilde{s}}^{2}\Lambda\overline{\mathbf{R}}_{n+1,n+2}^{(\tilde{s}-)})\biggr{)}.\end{split}$
(53)

The next proposition is about some Hölder type inequalities.

######  Proposition 5.3.

For $f$ on $\Sigma_{N}^{n}$, and $\tau_{1}>0,\tau_{2}>0$,
$1/\tau_{1}+1/\tau_{2}=1$, we have

$\displaystyle\left|\nu(f)-\nu_{0,\tilde{s}}(f)\right|$
$\displaystyle\leqslant
K(n,\beta,\tilde{s})\nu(\left|f\right|^{\tau_{1}})^{1/\tau_{1}}\nu\left(\left|\Delta_{\tilde{s}}^{2}\Lambda(\mathbf{R}^{(\tilde{s}-)}_{12}-\mathbf{q})\right|^{\tau_{2}}\right)^{1/\tau_{2}}$
(54)
$\displaystyle\left|\nu(f)-\nu^{\prime}_{0,\tilde{s}}(f)-\nu_{0,\tilde{s}}(f)\right|$
$\displaystyle\leqslant
K(n,\beta,\tilde{s})\nu(\left|f\right|^{\tau_{1}})^{1/\tau_{1}}\nu\left(\left|\Delta_{\tilde{s}}^{2}\Lambda(\mathbf{R}^{(\tilde{s}-)}_{12}-\mathbf{q})\right|^{2\tau_{2}}\right)^{1/\tau_{2}}.$
(55)

######  Remark 5.4.

The proof of this theorem is similar to the original proof in Talagrand’s
book. For the inequality (54), we just need to control
$\nu^{\prime}_{t,\tilde{s}}(f)$, and by Theorem 5.2, it suffices to control
the terms
$\nu_{t,\tilde{s}}(\varepsilon_{l}\varepsilon_{l^{\prime}}f(\Delta_{\tilde{s}}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(\tilde{s}-)}))$.
We claim that
$\left|\nu_{t,\tilde{s}}(f)\right|\leqslant\exp(4C_{\tilde{s}}n^{2}\beta^{2})\nu(f)$,
where $C_{\tilde{s}}:=\max\\{\Delta_{\tilde{s}}^{2}\Lambda\\}$ is the maximal
entry of the vector $\Delta_{\tilde{s}}^{2}\Lambda$. Because in Theorem 5.2,
$\left|R_{12}^{s}-q_{s}\right|\leqslant 2$ implies

$\left|\nu^{\prime}_{t,\tilde{s}}(f)\right|\leqslant
4C_{\tilde{s}}n^{2}\beta^{2}\nu_{t,\tilde{s}}(f)$

then integrating the above inequality will prove the claim. Finally applying
Hölder inequality, it proves (54). For (55), we just control the second order
derivative $\nu^{\prime\prime}_{t,\tilde{s}}$ in a similar way.

Recall the notation $\overline{\mathbf{R}}_{12}=\mathbf{R}_{12}-\mathbf{q}$,
using the Proposition 5.3, one can estimate the quantities
$\nu(\overline{\mathbf{R}}_{12}^{\intercal}A\overline{\mathbf{R}}_{12})$,
$\nu(\overline{\mathbf{R}}_{13}^{\intercal}A\overline{\mathbf{R}}_{12})$ and
$\nu(\overline{\mathbf{R}}_{34}^{\intercal}A\overline{\mathbf{R}}_{12})$,
where $A=(\\!(a_{s,t})\\!)_{s,t=1}^{m}$ is a symmetric $m\times m$ matrix.

We start with
$\nu(\overline{\mathbf{R}}_{12}^{\intercal}A\overline{\mathbf{R}}_{12})$, then
use the cavity method to do some second moment computation to derive a linear
system. We incorporate the randomness of $\delta$ in the following to get

$\displaystyle\begin{split}\nu\big{(}\overline{\mathbf{R}}^{\intercal}_{12}\Lambda
A\Lambda\overline{\mathbf{R}}_{12}\big{)}&=\nu\big{(}\sum_{s=1}\lambda_{s}\overline{R}_{12}^{(s)}\sum_{t=1}\lambda_{t}a_{s,t}\overline{R}_{12}^{(t)}\big{)}\\\
&=\nu\big{(}\operatorname{\mathds{E}}^{\delta}\big{(}\overline{R}_{12}^{(\delta)}\sum_{t=1}\lambda_{t}a_{\delta,t}\overline{R}_{12}^{(t)}\big{)}\big{)}=\hat{\nu}\big{(}\overline{R}_{12}^{(\delta)}\sum_{t=1}\lambda_{t}a_{\delta,t}\overline{R}_{12}^{(t)}\big{)}\end{split}$
(56)

where $\hat{\nu}:=\nu\otimes\operatorname{\mathds{E}}^{\delta}$, and

$\overline{R}_{12}^{(\delta)}=\frac{1}{|I_{\delta}|}\sum_{i\in
I_{\delta}}(\sigma_{i}^{1}\sigma_{i}^{2}-q_{\delta})=\frac{1}{|I_{\delta}|}\sum_{\begin{subarray}{c}i\in
I_{\delta},\\\
i\neq\delta_{\star}\end{subarray}}\sigma_{i}^{1}\sigma_{i}^{2}+\frac{1}{|I_{\delta}|}\varepsilon_{1}\varepsilon_{2}-q_{\delta}.$

By using the symmetry among sites inside the species $\delta$, we continue the
above expression (56):

$\displaystyle\begin{split}\hat{\nu}\biggl{(}(\varepsilon_{1}\varepsilon_{2}-q_{\delta})\sum_{t=1}\lambda_{t}a_{\delta,t}\overline{R}_{12}^{(t)}\biggr{)}&=\hat{\nu}\biggl{(}(\varepsilon_{1}\varepsilon_{2}-q_{\delta})\bigl{(}\sum_{t\neq\delta}\lambda_{t}a_{\delta,t}\overline{R}_{12}^{(t)}+\lambda_{\delta}a_{\delta,\delta}\overline{R}_{12}^{(\delta)}\bigr{)}\bigg{)}\\\
&=\hat{\nu}((\varepsilon_{1}\varepsilon_{2}-q_{\delta})A_{\delta}\Lambda\overline{\mathbf{R}}_{12}^{(\delta-)})+\frac{1}{N}\hat{\nu}((\varepsilon_{1}\varepsilon_{2}-q_{\delta})a_{\delta,\delta}\varepsilon_{1}\varepsilon_{2}).\end{split}$
(57)

Note that, in the above expression (57), there is randomness inside
$\varepsilon_{1},\varepsilon_{2},\delta$, which describes the way of choosing
the species to decouple the spin, and the computation proceeds in a quenched
way, because all the randomness only comes from $\delta$, which takes finitely
many values. In the following subsections, we deal with the two terms in (57)
separately.

###  5.1. Estimation of
$\hat{\nu}((\varepsilon_{1}\varepsilon_{2}-q_{\delta})A_{\delta}\Lambda\overline{\mathbf{R}}_{12}^{(\delta-)})$

Recall the definition of $\hat{\nu}$ below (56) and notice that

$\displaystyle\hat{\nu}((\varepsilon_{1}\varepsilon_{2}-q_{\delta})A_{\delta}\Lambda\overline{\mathbf{R}}_{12}^{(\delta-)})=\sum_{s=1}^{m}\lambda_{s}\cdot\nu((\varepsilon_{1}\varepsilon_{2}-q_{s})A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)}).$
(58)

Therefore it suffices to compute
$\nu((\varepsilon_{1}\varepsilon_{2}-q_{s})A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)})$.
First let’s take

$\tilde{f}=(\varepsilon_{1}\varepsilon_{2}-q_{s})A_{s}\Lambda\overline{\mathbf{R}}^{(s-)}_{12},$

by the second Hölder inequality (55) with $\tau_{1}=3,\tau_{2}=3/2$ in the
Proposition 5.3, we have

$\nu(\tilde{f})=\nu_{0,s}^{\prime}(\tilde{f})+\mathfrak{R}.$

Because $\nu_{0,s}(\tilde{f})=0$ by Lemma 5.1, where $\mathfrak{R}$ represents
the higher order terms with

$\displaystyle\mathfrak{R}\leqslant
K\left(N^{-3/2}+\sum_{s=1}^{m}\nu(|R_{12}^{(s-)}-q_{s}|^{3})\right)\max_{s,t}|A_{s,t}|.$
(59)

For the derivative term, by Theorem 5.2, $\nu_{0,s}^{\prime}(\tilde{f})$ is a
sum of terms in the form of

$\beta^{2}\nu_{0,s}(\varepsilon_{l}\varepsilon_{l^{\prime}}(\varepsilon_{1}\varepsilon_{2}-q_{s})\cdot
A_{s}\Lambda\overline{\mathbf{R}}_{12}^{-}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(s-)}).$

In particular, we introduce a general formula as a corollary to Theorem 5.2 to
compute $\nu_{0,s}^{\prime}(\tilde{f})$.

######  Corollary 5.5.

Consider a function $f^{-}$ on $\Sigma_{N-1}^{n}$ and two integers $x\neq
y\leqslant n$. Then

$\displaystyle\begin{split}\nu_{0,s}^{\prime}\left((\varepsilon_{x}\varepsilon_{y}-q_{s})f^{-}\right)&=\sum_{1\leqslant
l<l^{\prime}\leqslant
n}b_{s}(l,l^{\prime};x,y)\nu_{0,s}(f^{-}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(s-)})\\\
&\qquad\qquad-n\sum_{l\leqslant
n}b_{s}(l,n+1;x,y)\nu_{0,s}(f^{-}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{l,n+1}^{(s-)})\\\
&\qquad\qquad\qquad\qquad+\frac{n(n+1)}{2}b_{s}(0)\nu_{0,s}(f^{-}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{n+1,n+2}^{(s-)})\end{split}$
(60)

where
$b_{s}(l,l^{\prime};x,y)=b_{s}(\left|\\{l,l^{\prime}\\}\cap\\{x,y\\}\right|)$
and

$\displaystyle\begin{split}b_{s}(2)&=\beta^{2}\nu_{0,s}(\varepsilon_{1}\varepsilon_{2}(\varepsilon_{1}\varepsilon_{2}-q_{s}))=\beta^{2}\nu_{0,s}(1-\varepsilon_{1}\varepsilon_{2}q_{s})=\beta^{2}(1-q_{s}^{2})\\\
b_{s}(1)&=\beta^{2}\nu_{0,s}(\varepsilon_{1}\varepsilon_{3}(\varepsilon_{1}\varepsilon_{2}-q_{s}))=\beta^{2}\nu_{0,s}(\varepsilon_{2}\varepsilon_{3}-\varepsilon_{1}\varepsilon_{3}q_{s})=\beta^{2}q_{s}(1-q_{s})\\\
b_{s}(0)&=\beta^{2}\nu_{0,s}(\varepsilon_{3}\varepsilon_{4}(\varepsilon_{1}\varepsilon_{2}-q_{s}))=\beta^{2}\nu_{0,s}(\varepsilon_{1}\varepsilon_{2}\varepsilon_{3}\varepsilon_{4}-\varepsilon_{3}\varepsilon_{4}q_{s})=\beta^{2}(\hat{q}_{s}-q_{s}^{2})\end{split}$
(61)

In the above Corollary,

$q_{s}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h),\qquad\hat{q}_{s}:=\operatorname{\mathds{E}}\tanh^{4}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h),$

and we used Lemma 5.1 to get

$\displaystyle\nu_{0,s}(\varepsilon_{l}\varepsilon_{l^{\prime}}(\varepsilon_{1}\varepsilon_{2}-q_{s})\cdot
A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(s-)})$
$\displaystyle\qquad=\nu_{0,s}(\varepsilon_{l}\varepsilon_{l^{\prime}}(\varepsilon_{1}\varepsilon_{2}-q_{s}))\cdot\nu_{0,s}(A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(s-)}).$

Using Corollary 5.5 with
$n=2,x=1,y=2,f^{-}=A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)}$, we have

$\displaystyle\begin{split}&\nu^{\prime}_{0,s}((\varepsilon_{x}\varepsilon_{y}-q_{s})f^{-})\\\
&=\nu_{0,s}^{\prime}\bigl{(}(\varepsilon_{1}\varepsilon_{2}-q_{s})\cdot
A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)}\bigr{)}\\\
&=b_{s}(2)\cdot\nu_{0,s}(A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{12}^{(s-)})-2b_{s}(1)\cdot\nu_{0,s}(A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{13}^{(s-)})\\\
&\qquad-2b_{s}(1)\cdot\nu_{0,s}(A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{23}^{(s-)})+3b_{s}(0)\cdot\nu_{0,s}(A_{s}\Lambda\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{34}^{(s-)}).\end{split}$
(62)

For terms like $\nu_{0,s}(\Lambda
A_{s}\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(s-)})$
in the above expression, we apply Hölder inequality (54) with
$\tau_{1}=\frac{3}{2},\tau_{2}=3$, to obtain

$\nu_{0,s}(\Lambda
A_{s}\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(s-)})=\nu(\Lambda
A_{s}\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(s-)})+\mathfrak{R}$

and changing from $\nu(\Lambda
A_{s}\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(s-)})$
to $\nu(\Lambda
A_{s}\overline{\mathbf{R}}_{12}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}})$
will result in an error term $\mathfrak{R}$, _i.e.,_

$\nu(\Lambda
A_{s}\overline{\mathbf{R}}_{12}^{(s-)}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}}^{(s-)})=\nu(\Lambda
A_{s}\overline{\mathbf{R}}_{12}\cdot\Delta_{s}^{2}\Lambda\overline{\mathbf{R}}_{\ell\ell^{\prime}})+\mathfrak{R}.$

This is because
$R_{ll^{\prime}}^{(s-)}=R_{ll^{\prime}}^{(s)}-\frac{1}{|I_{s}|}\varepsilon_{l}\varepsilon_{l^{\prime}}$,
and the fact that $\frac{1}{N}\nu(\Lambda
A_{s}\overline{\mathbf{R}}_{12})=\mathfrak{R}$ as $xy\leqslant
x^{\frac{3}{2}}+y^{3}$ for $x,y\geqslant 0$. Combining (62) and (58), we have

$\displaystyle\begin{split}\hat{\nu}\big{(}(\varepsilon_{1}\varepsilon_{2}-q_{\delta})A_{\delta}\Lambda\overline{\mathbf{R}}_{12}^{(\delta-)}\big{)}&=\nu\big{(}\overline{\mathbf{R}}_{12}^{\intercal}S_{A}B(2)\Delta^{2}\Lambda\overline{\mathbf{R}}_{12}\big{)}-4\nu\big{(}\overline{\mathbf{R}}_{12}^{\intercal}S_{A}B(1)\Delta^{2}\Lambda\overline{\mathbf{R}}_{13}\big{)}\\\
&\qquad+3\nu\big{(}\overline{\mathbf{R}}_{12}^{\intercal}S_{A}B(0)\Delta^{2}\Lambda\overline{\mathbf{R}}_{34}\big{)}+\mathfrak{R}\end{split}$
(63)

where $S_{A}:=\Lambda A\Lambda$,

$\displaystyle
B(i)=\operatorname{diag}(b_{1}(i),b_{2}(i),\ldots,b_{m}(i)),\quad\text{for}\quad
i=0,1,2,$ (64)

with $b_{s}(i)$ as given in (61) and $\mathfrak{R}$ represent the higher order
terms as defined in (59).

###  5.2. Estimation of
$\hat{\nu}((\varepsilon_{1}\varepsilon_{2}-q_{\delta})a_{\delta,\delta}\varepsilon_{1}\varepsilon_{2})$

For the second term in the last step of (57), we get

$\displaystyle\begin{split}\hat{\nu}((\varepsilon_{1}\varepsilon_{2}-q_{\delta})a_{\delta,\delta}\varepsilon_{1}\varepsilon_{2})&=\hat{\nu}_{0,\delta}(a_{\delta,\delta}(1-\varepsilon_{1}\varepsilon_{2}q_{\delta}))+\mathfrak{R}\\\
&=\operatorname{\mathds{E}}^{\delta}(a_{\delta,\delta}(1-q_{\delta}^{2}))+\mathfrak{R}=\sum_{s=1}^{m}\lambda_{s}a_{s,s}(1-q_{s}^{2})+\mathfrak{R}.\end{split}$
(65)

In the above computation, the first equality is due to the Hölder inequality
(54) with $\tau_{1}=\infty,\tau_{2}=1$ in the Proposition 5.3, the second
equality is due to the property of $\nu_{0,s}$ as in Lemma 5.1.

By collecting all the terms, we have

$\displaystyle\begin{split}\nu(\overline{\mathbf{R}}^{\intercal}_{12}S_{A}\overline{\mathbf{R}}_{12})&=\frac{1}{N}\sum_{s=1}^{m}\lambda_{s}a_{s,s}(1-q_{s}^{2})+\nu\big{(}\overline{\mathbf{R}}_{12}^{\intercal}S_{A}B(2)\Delta^{2}\Lambda)\overline{\mathbf{R}}_{12}\big{)}\\\
&\qquad-4\nu\big{(}\overline{\mathbf{R}}_{12}^{\intercal}S_{A}B(1)\Delta^{2}\Lambda\overline{\mathbf{R}}_{13}\big{)}+3\nu\big{(}\overline{\mathbf{R}}_{12}^{\intercal}S_{A}B(0)\Delta^{2}\Lambda\overline{\mathbf{R}}_{34}\big{)}+\mathfrak{R}.\end{split}$
(66)

By a similar argument, we get

$\displaystyle\begin{split}\nu(\overline{\mathbf{R}}^{\intercal}_{12}S_{A}\overline{\mathbf{R}}_{13})&=\frac{1}{N}\sum_{s=1}^{m}\lambda_{s}a_{s,s}(q_{s}-q_{s}^{2})+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S_{A}B(1)\Delta^{2}\Lambda\overline{\mathbf{R}}_{12})\\\
&\qquad+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S_{A}(B(2)-2B(1)-3B(0))\Delta^{2}\Lambda\overline{\mathbf{R}}_{13})\\\
&\qquad\qquad+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S_{A}(6B(0)-3B(1))\Delta^{2}\Lambda\overline{\mathbf{R}}_{34})+\mathfrak{R},\end{split}$
(67)

and

$\displaystyle\begin{split}\nu(\overline{\mathbf{R}}^{\intercal}_{12}S_{A}\overline{\mathbf{R}}_{34})&=\frac{1}{N}\sum_{s=1}^{m}\lambda_{s}a_{s,s}(\hat{q}_{s}-q_{s}^{2})+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S_{A}B(0)\Delta^{2}\Lambda\overline{\mathbf{R}}_{12})\\\
&\qquad+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S_{A}(4B(1)-8B(0))\Delta^{2}\Lambda\overline{\mathbf{R}}_{13})\\\
&\qquad\qquad+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S_{A}(B(2)-8B(1)+10B(0))\Delta^{2}\Lambda\overline{\mathbf{R}}_{34})+\mathfrak{R}.\end{split}$
(68)

We define the $m\times m$ diagonal matrices

$\displaystyle\Theta(i):=\beta^{-2}B(i)\Lambda^{-1},\qquad\text{ for
}i=0,1,2,$ (69)

so that

$\displaystyle\sum_{s=1}^{m}\lambda_{s}a_{s,s}(1-q_{s}^{2})=\operatorname{Tr}(S_{A}\Theta(2)),\quad$
$\displaystyle\quad\sum_{s=1}^{m}\lambda_{s}a_{s,s}(q_{s}-q_{s}^{2})=\operatorname{Tr}(S_{A}\Theta(1)),$
$\displaystyle\text{ and
}\sum_{s=1}^{m}\lambda_{s}a_{s,s}(\hat{q}_{s}-q_{s}^{2})$
$\displaystyle=\operatorname{Tr}(S_{A}\Theta(0)).$

We also define the $m\times m$ matrices (not necessarily symmetric)

$\displaystyle\hat{B}(i):=B(i)\Delta^{2}\Lambda,\qquad\text{ for }i=0,1,2,$
(70)

Thus, from equations (66), (67) and (68), for any $m\times m$ symmetric matrix
$S$ we have

$\displaystyle\begin{split}\nu(\overline{\mathbf{R}}_{12}^{\intercal}S\overline{\mathbf{R}}_{12})&=\nu(\overline{\mathbf{R}}_{12}^{\intercal}S\hat{B}(2)\overline{\mathbf{R}}_{12})-4\nu(\overline{\mathbf{R}}_{12}^{\intercal}S\hat{B}(1)\overline{\mathbf{R}}_{13})\\\
&\qquad+3\nu(\overline{\mathbf{R}}_{12}^{\intercal}S\hat{B}(0)\overline{\mathbf{R}}_{34})+\frac{1}{N}\operatorname{Tr}(S\Theta(2))+\mathfrak{R},\\\
\nu(\overline{\mathbf{R}}_{12}^{\intercal}S\overline{\mathbf{R}}_{13})&=\nu(\overline{\mathbf{R}}_{12}^{\intercal}S\hat{B}(1)\overline{\mathbf{R}}_{12})+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S(\hat{B}(2)-2\hat{B}(1)-3\hat{B}(0))\overline{\mathbf{R}}_{13})\\\
&\qquad+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S(6\hat{B}(0)-3\hat{B}(1))\overline{\mathbf{R}}_{34})+\frac{1}{N}\operatorname{Tr}(S\Theta(1))+\mathfrak{R},\\\
\nu(\overline{\mathbf{R}}_{12}^{\intercal}S\overline{\mathbf{R}}_{34})&=\nu(\overline{\mathbf{R}}_{12}^{\intercal}S\hat{B}(0)\overline{\mathbf{R}}_{12})+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S(4\hat{B}(1)-8\hat{B}(0))\overline{\mathbf{R}}_{13})\\\
&\qquad+\nu(\overline{\mathbf{R}}_{12}^{\intercal}S(\hat{B}(2)-8\hat{B}(1)+10\hat{B}(0))\overline{\mathbf{R}}_{34})+\frac{1}{N}\operatorname{Tr}(S\Theta(0))+\mathfrak{R}.\end{split}$
(71)

The system of equations (71) can be written as a system of linear equations in
the variables
$\nu((R_{12}^{(s)}-q_{s})(R^{(t)}_{\ell\ell^{\prime}}-q_{t})),s\leqslant
t,\ell<\ell^{\prime}$ as follows.

We define the three $m\times m$ matrices

$\displaystyle
U(0):=\nu(\overline{\mathbf{R}}_{12}\overline{\mathbf{R}}_{34}^{\intercal}),\
U(1):=\nu(\overline{\mathbf{R}}_{12}\overline{\mathbf{R}}_{13}^{\intercal}),\
U(2):=\nu(\overline{\mathbf{R}}_{12}\overline{\mathbf{R}}_{12}^{\intercal}).$

Thus equation (71) can be restated in the following way: for any symmetric
$m\times m$ matrix $S$

$\displaystyle\operatorname{Tr}(SU(i))=\sum_{j=0}^{2}\sum_{k=0}^{2}C_{ij}(k)\operatorname{Tr}(S\hat{B}(k)U(j))+\frac{1}{N}\operatorname{Tr}(S\Theta(i))+\mathfrak{R},\qquad
i=0,1,2,$ (72)

where

$\displaystyle C(0):=\begin{pmatrix}10&-8&1\\\ 6&-3&0\\\
3&0&0\end{pmatrix},\quad C(1):=\begin{pmatrix}-8&4&0\\\ -3&-2&1\\\
0&-4&0\end{pmatrix}\text{ and }C(2):=\begin{pmatrix}1&0&0\\\ 0&1&0\\\
0&0&1\end{pmatrix}$

(rows and columns are indexed by $0,1,2$).

The matrices, $C(0),C(1),C(2)$, commute with each other and thus are
simultaneously “upper-triangularizable”. This fact is crucially used in the
variance computation for the overlap in Section 6. Explaining the
commutativity property of the $C(\cdot)$ matrices involves understanding the
underlying algebraic structure of the Gibbs measure and is an interesting open
question. Note that this phenomenon appears in both SK and MSK models. One can
guess that this is probably due to some high-level symmetry among replicas in
the RS regime. Also, in the MSK model, it is likely to be connected to the
synchronization property of overlap discovered by Panchenko (see [Pan15]).

Note that, the upper bound on $\mathfrak{R}$ depends on $S$ only through
$\max_{1\leqslant r,t\leqslant m}|S_{r,t}|$. We denote by $\mathbf{e}_{i}$,
the column vector with $1$ in the $i$-th coordinate and zero elsewhere. Taking
$S=(\mathbf{e}_{l}\mathbf{e}_{k}^{\intercal}+\mathbf{e}_{k}\mathbf{e}_{l}^{\intercal})/2$
for all possible choices of $1\leqslant k\leqslant l\leqslant m$, from
equation (72), we get

$\displaystyle
U(i)=\sum_{j=0}^{2}\sum_{k=0}^{2}C_{ij}(k)\cdot\operatorname{sym}(\hat{B}(k)U(j))+\frac{1}{N}\Theta(i)+\boldsymbol{\mathfrak{R}},\quad\text{for}\quad
i=0,1,2,$ (73)

where $\operatorname{sym}(A):=(A+A^{\intercal})/2$ for a square matrix $A$ and
$\boldsymbol{\mathfrak{R}}$ is a $m\times m$ matrix with

$\max_{p,q}|\boldsymbol{\mathfrak{R}}_{p,q}|\leqslant
K\left(N^{-3/2}+\sum_{s=1}^{m}\nu(|R_{12}^{(s-)}-q_{s}|^{3})\right).$

######  Remark 5.6.

The above argument does not require the positive semi-definiteness of
$\Delta^{2}$. It seems that at least in a high-temperature region, the RS
solution given by Parisi formula [Pan15] is still valid even for indefinite
$\Delta^{2}$.

##  6\. Variance of overlap

In Section 5, we introduced a species-wise cavity method for the MSK model to
derive a linear system involving the quadratic form of overlap. In this
section, by studying the linear system, we solve the overlap vectors’
variance-covariance structure. Along the way, we obtain the AT-line condition
in the MSK model.

###  6.1. Proof of Theorem 1.11

We recall the notations

$\displaystyle b_{s}(0)$ $\displaystyle=\beta^{2}(\hat{q}_{s}-q_{s}^{2}),\quad
b_{s}(1)=\beta^{2}q_{s}(1-q_{s}),\quad b_{s}(2)=\beta^{2}(1-q_{s}^{2})$

for $s=1,2,\ldots,m$; and

$\displaystyle
B(i)=\operatorname{diag}(b_{1}(i),b_{2}(i),\ldots,b_{m}(i)),\quad\Theta(i)=\beta^{-2}B(i)\Lambda^{-1},\quad\hat{B}(i)=B(i)\Delta^{2}\Lambda$

for $i=0,1,2.$ Then with

$\displaystyle C(0):=\begin{pmatrix}10&-8&1\\\ 6&-3&0\\\
3&0&0\end{pmatrix},\quad C(1):=\begin{pmatrix}-8&4&0\\\ -3&-2&1\\\
0&-4&0\end{pmatrix}\text{ and }C(2):=\begin{pmatrix}1&0&0\\\ 0&1&0\\\
0&0&1\end{pmatrix}$

(row and column indexed by $0,1,2$), we have from (73)

$\displaystyle
U(i)=\sum_{j=0}^{2}\sum_{k=0}^{2}C_{ij}(k)\cdot\operatorname{sym}(\hat{B}(k)U(j))+\frac{1}{N}\Theta(i)+\boldsymbol{\mathfrak{R}},\qquad
i=0,1,2,$

where $\operatorname{sym}(A):=(A+A^{\intercal})/2$.

The matrices, $C(0),C(1),C(2)$, commute with each other and thus are
simultaneously “upper-triangularizable”. In particular, with

$V=\begin{pmatrix}-3&2&0\\\ 3&-4&1\\\ 1&-2&1\end{pmatrix}$

and defining $T(k):=VC(k)V^{-1},k=0,1,2$ we have

$\displaystyle T(0)=\begin{pmatrix}3&-3&0\\\ 0&3&0\\\
0&0&1\end{pmatrix},T(1)=\begin{pmatrix}-4&2&0\\\ 0&-4&0\\\
0&0&-2\end{pmatrix},T(2)=\begin{pmatrix}1&0&0\\\ 0&1&0\\\ 0&0&1\end{pmatrix}.$

Define $\hat{U}(l):=\sum_{i=0}^{2}V_{li}U(i)$ for $l=0,1,2$, _i.e.,_

$\displaystyle\hat{U}(0)$ $\displaystyle:=-3U(0)+2U(1),$
$\displaystyle\hat{U}(1)$ $\displaystyle:=3U(0)-4U(1)+U(2),$
$\displaystyle\hat{U}(2)$ $\displaystyle:=U(0)-2U(1)+U(2).$

Similarly, we define

$\hat{\Theta}(l):=\sum_{i=0}^{2}V_{li}\Theta(i)\text{ for }l=0,1,2.$

We have for $l=0,1,2$,

$\displaystyle\hat{U}(l)$
$\displaystyle=\sum_{i=0}^{2}V_{li}\sum_{j=0}^{2}\sum_{k=0}^{2}C_{ij}(k)\cdot\operatorname{sym}(\hat{B}(k)U(j))+\frac{1}{N}\sum_{i=0}^{2}V_{li}\Theta(i)+\boldsymbol{\mathfrak{R}}$
$\displaystyle=\sum_{j=0}^{2}\sum_{k=0}^{2}(VC(k))_{lj}\cdot\operatorname{sym}(\hat{B}(k)U(j))+\frac{1}{N}\hat{\Theta}(l)+\boldsymbol{\mathfrak{R}}$
$\displaystyle=\sum_{k=0}^{2}\operatorname{sym}\bigl{(}\hat{B}(k)\cdot\sum_{j=0}^{2}(T(k)V)_{lj}U(j)\bigr{)}+\frac{1}{N}\hat{\Theta}(l)+\boldsymbol{\mathfrak{R}}.$

Simplifying, we get

$\displaystyle\begin{split}\operatorname{sym}\biggl{(}\bigl{(}I-\sum_{i=0}^{2}V_{li}\hat{B}(i)\bigr{)}\cdot\hat{U}(l)\biggr{)}&=\frac{1}{N}\hat{\Theta}(l)+\boldsymbol{\mathfrak{R}},\qquad
l=1,2,\\\ \text{ and
}\operatorname{sym}\biggl{(}\bigl{(}I-\sum_{i=0}^{2}V_{1i}\hat{B}(i)\bigr{)}\cdot\hat{U}(0)\biggr{)}&=\operatorname{sym}\biggl{(}\sum_{i=0}^{2}V_{0i}\hat{B}(i)\cdot\hat{U}(1)\biggr{)}+\frac{1}{N}\hat{\Theta}(0)+\boldsymbol{\mathfrak{R}}.\end{split}$
(74)

Now consider the $m\times m$ diagonal matrices
$\Gamma,\Gamma^{\prime},\Gamma^{\prime\prime}$, whose $s$-th diagonal entries
are respectively given by

$\displaystyle\gamma_{s}$
$\displaystyle=\beta^{-2}(b_{s}(0)-2b_{s}(1)+b_{s}(2))=1-2q_{s}+\hat{q}_{s},$
$\displaystyle\gamma^{\prime}_{s}$
$\displaystyle=\beta^{-2}(3b_{s}(0)-4b_{s}(1)+b_{s}(2))=1-4q_{s}+3\hat{q}_{s},$
$\displaystyle\gamma^{\prime\prime}_{s}$
$\displaystyle=\beta^{-2}(-3b_{s}(0)+2b_{s}(1))=2q_{s}+q_{s}^{2}-3\hat{q}_{s}.$

Then,

$\displaystyle\sum_{i=0}^{2}V_{2i}\hat{B}(i)$
$\displaystyle=\beta^{2}\Gamma\Delta^{2}\Lambda,\quad\sum_{i=0}^{2}V_{1i}\hat{B}(i)=\beta^{2}\Gamma^{\prime}\Delta^{2}\Lambda,\quad\sum_{i=0}^{2}V_{0i}\hat{B}(i)=\beta^{2}\Gamma^{\prime\prime}\Delta^{2}\Lambda,$
$\displaystyle\hat{\Theta}(2)$
$\displaystyle=\Gamma\Lambda^{-1},\quad\hat{\Theta}(1)=\Gamma^{\prime}\Lambda^{-1},\quad\hat{\Theta}(0)=\Gamma^{\prime\prime}\Lambda^{-1}.$

Combining with (74) we have the result. $\blacksquare$

Before going to the proof of Theorem 1.13, we will prove the following lemma,
which essentially solves the continuous Lyapunov equation. Recall that a
square matrix $A$ is _stable_ if all the eigenvalues of $A$ have a strictly
negative real part.

######  Lemma 6.1.

Let $A$ be a $m\times m$ stable matrix, and $C$ be a symmetric matrix. Suppose
that the symmetric matrix $X$ satisfies the equation,

$\operatorname{sym}(AX)=-C.$

Then we have

$X=\int_{0}^{\infty}e^{\frac{t}{2}A}Ce^{\frac{t}{2}A^{\intercal}}dt.$

###### Proof of Lemma 6.1.

First we consider the case when $A$ is stable and similar to a diagonal
matrix, _i.e.,_ $A=SDS^{-1}$ for a diagonal matrix $D$ with negative diagonal
entries and an invertible matrix $S$. We can write
$\operatorname{sym}(AX)=S\operatorname{sym}(DS^{-1}XS^{-\intercal})S^{\intercal}$
and solve for $\operatorname{sym}(DY)=-S^{-1}CS^{-\intercal}$ where
$Y=S^{-1}XS^{-\intercal}$. Furthermore,

$\int_{0}^{\infty}Se^{\frac{t}{2}D}S^{-1}CS^{-\intercal}e^{\frac{t}{2}D}S^{\intercal}dt=\int_{0}^{\infty}e^{\frac{t}{2}A}Ce^{\frac{t}{2}A^{\intercal}}dt.$

Thus, w.l.o.g. $A$ can be taken as a diagonal matrix.

We write $\operatorname{vec}(A)$ as the $m^{2}\times 1$ vector formed by
stacking the columns of $A$, _i.e.,_
$\operatorname{vec}(A)=\sum_{s=1}^{m}e_{s}\otimes
Ae_{s}=\sum_{s=1}^{m}A^{\intercal}e_{s}\otimes e_{s}$, where $e_{s}$ is the
canonical basis vector with 1 at $s$-th entry and 0 elsewhere, $\otimes$ is
the Kronecker product of matrices. We define a linear map $\phi$ from $m\times
m$ matrices to $m^{2}\times m^{2}$ matrices:

$\phi(A):=\frac{1}{2}(A\otimes I+I\otimes A).$

One can easily check that, $\operatorname{vec}(AXB^{\intercal})=(B\otimes
A)\operatorname{vec}(X)$ for any symmetric matrix $X$, and thus we have

$\operatorname{vec}(\operatorname{sym}(AX))=\phi(A)\operatorname{vec}(X).$

From $\operatorname{sym}(AX)=-C$, we get

$\displaystyle\phi(A)\operatorname{vec}(X)=-\operatorname{vec}(C).$ (75)

Since $A$ is negative definite diagonal matrix, $\phi(A)$ is also a negative
definite diagonal matrix and hence invertible with inverse given by

$(\phi(A))^{-1}=-\int_{0}^{\infty}e^{t\phi(A)}dt=-\int_{0}^{\infty}e^{\frac{t}{2}A}\otimes
e^{\frac{t}{2}A}dt.$

Thus, we have

$\operatorname{vec}(X)=-(\phi(A))^{-1}\operatorname{vec}(C)=\int_{0}^{\infty}(e^{\frac{t}{2}A}\otimes
e^{\frac{t}{2}A})\operatorname{vec}(C)dt=\int_{0}^{\infty}\operatorname{vec}(e^{\frac{t}{2}A}Ce^{\frac{t}{2}A})dt$

and

$X=\int_{0}^{\infty}e^{\frac{t}{2}A}Ce^{\frac{t}{2}A}dt.$

In the general case, by Jordan decomposition, $A$ is similar to an upper
triangular matrix with diagonal entries having negative real parts. The exact
proof goes through. $\blacksquare$

###  6.2. Proof of Theorem 1.13

First we note that $0\leqslant\gamma_{s}\leqslant
1,-1<\gamma_{s}^{\prime}\leqslant\gamma_{s}\leqslant 1$ and thus
$\rho(\Gamma\Delta^{2}\Lambda)\leqslant\rho(\Delta^{2}\Lambda),\rho(\Gamma^{\prime}\Delta^{2}\Lambda)\leqslant\rho(\Delta^{2}\Lambda)$.
If $\beta^{2}<\beta_{c}^{2}=\rho(\Delta^{2}\Lambda)^{-1}$, it is easy to see
that the matrices
$\beta^{2}\Gamma\Delta^{2}\Lambda-I,\beta^{2}\Gamma^{\prime}\Delta^{2}\Lambda-I$
are stable. Furthermore, we have

$\displaystyle\hat{U}(2)$
$\displaystyle=\frac{1}{N}\cdot\int_{0}^{\infty}\exp\left(-\frac{t}{2}(I-\beta^{2}\Gamma\Delta^{2}\Lambda)\right)\Gamma\Lambda^{-1}\exp\left(-\frac{t}{2}(I-\beta^{2}\Lambda\Delta^{2}\Gamma)\right)dt+\boldsymbol{\mathfrak{R}}$
$\displaystyle=\frac{1}{N}\cdot\int_{0}^{\infty}\Gamma^{\nicefrac{{1}}{{2}}}\Lambda^{-\nicefrac{{1}}{{2}}}\exp\left(-t(I-\beta^{2}\Gamma^{\nicefrac{{1}}{{2}}}\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}\Gamma^{\nicefrac{{1}}{{2}}})\right)\Gamma^{\nicefrac{{1}}{{2}}}\Lambda^{-\nicefrac{{1}}{{2}}}dt+\boldsymbol{\mathfrak{R}}$
$\displaystyle=\frac{1}{N}\cdot\Gamma^{\nicefrac{{1}}{{2}}}\Lambda^{-\nicefrac{{1}}{{2}}}(I-\beta^{2}\Gamma^{\nicefrac{{1}}{{2}}}\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}\Gamma^{\nicefrac{{1}}{{2}}})^{-1}\Gamma^{\nicefrac{{1}}{{2}}}\Lambda^{-\nicefrac{{1}}{{2}}}+\boldsymbol{\mathfrak{R}}$
$\displaystyle=\frac{1}{N}\cdot\Gamma(I-\beta^{2}\Delta^{2}\Lambda\Gamma)^{-1}\Lambda^{-1}+\boldsymbol{\mathfrak{R}}.$

In the second equality, we used the fact that
$ABA^{T}=B^{1/2}(B^{-1/2}AB^{1/2})(B^{-1/2}AB^{1/2})^{T}B^{1/2}$ with
$A=\exp\left(-\frac{t}{2}(I-\beta^{2}\Gamma\Delta^{2}\Lambda)\right)$ and
$B:=\Gamma\Lambda^{-1}$, a diagonal matrix with positive entries.

Similarly, we get

$\displaystyle\hat{U}(1)$
$\displaystyle=\frac{1}{N}\cdot\Gamma^{\prime}(I-\beta^{2}\Delta^{2}\Lambda\Gamma^{\prime})^{-1}\Lambda^{-1}+\boldsymbol{\mathfrak{R}}$

and

$\displaystyle\hat{U}(0)$
$\displaystyle=\int_{0}^{\infty}\exp\left(-\frac{t}{2}(I-\beta^{2}\Gamma^{\prime}\Delta^{2}\Lambda)\right)\operatorname{sym}\bigl{(}\beta^{2}\Gamma^{\prime\prime}\Delta^{2}\Lambda\cdot\hat{U}(1)\bigr{)}\exp\left(-\frac{t}{2}(I-\beta^{2}\Lambda\Delta^{2}\Gamma^{\prime})\right)dt$
(76)
$\displaystyle+\frac{1}{N}\cdot\int_{0}^{\infty}\exp\left(-\frac{t}{2}(I-\beta^{2}\Gamma^{\prime}\Delta^{2}\Lambda)\right)\Gamma^{\prime\prime}\Lambda^{-1}\exp\left(-\frac{t}{2}(I-\beta^{2}\Lambda\Delta^{2}\Gamma^{\prime})\right)dt+\boldsymbol{\mathfrak{R}}.$
(77)

Finally, we use the fact that $N\boldsymbol{\mathfrak{R}}$ converges to $0$
entrywise as $N\to\infty$. $\blacksquare$

##  7\. Uniqueness of $\mathbf{q}_{\star}$ and Replica Symmetry Breaking

In this section, for _positive-definite_ $\Delta^{2}$, we prove the AT line
condition of the MSK model when $h>0$, beyond which the MSK model is in
replica symmetry breaking phase. From classical literature of the SK model, we
know that the uniqueness of $\mathbf{q}_{\star}$ is essential to characterize
the AT line condition. In the SK model, the uniqueness was proved using
contraction argument for $\beta$ small, and the Latala-Guerra Lemma
[Tal11a]*Proposition 1.3.8 gives the uniqueness result for $h>0$. In the MSK
model, we prove the uniqueness of $\mathbf{q}_{\star}$ for small $\beta$ by
extending the contraction argument, but extending Latala-Guerra Lemma becomes
more challenging, which is about analysis of a complicated nonlinear system
(15). In [BSS19], they use an elementary approach to prove this for _positive-
definite_ $2$-species model, but the idea seems difficult to be generalized.
We first prove the uniqueness result for small $\beta$ in Theorem 1.4, then we
prove Proposition 1.6 which gives the uniqueness result for _indefinite_
2-species models when $h>0$.

###  7.1. Proof of Theorem 1.4

We need to prove that the following system

$\displaystyle
q_{s}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{(\Delta^{2}\Lambda\mathbf{q})_{s}}+h),\qquad
s=1,2,\ldots,m$ (78)

has a unique solution where
$\mathbf{q}=(q_{1},q_{2},\ldots,q_{m})^{\intercal}$. We rewrite the equations
in terms of
$\mathbf{x}=(x_{1},x_{2},\ldots,x_{m})^{\intercal}:=\Lambda^{\nicefrac{{1}}{{2}}}\mathbf{q}$.
Let $f(y):=\tanh^{2}(y)$. Then the system of equations (78) is equivalent to

$\displaystyle x_{s}={\psi}_{s},\qquad s=1,2,\ldots,m$ (79)

where

$\psi_{s}:={\lambda^{\nicefrac{{1}}{{2}}}_{s}}\operatorname{\mathds{E}}f(\beta\eta\sqrt{(\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}\mathbf{x})_{s}}+h).$

We define ${\Psi}(\mathbf{x})=(\psi_{1},\psi_{2},\ldots,\psi_{m})$. It’s easy
to compute the Jacobian matrix of the map $\Psi$,

$J({\Psi})=\left(\\!\\!\\!\left(\frac{\partial\psi_{s}}{\partial
x_{t}}\right)\\!\\!\\!\right)_{s,t=1}^{m}$

where

$\frac{\partial\psi_{s}}{\partial
x_{t}}=\frac{\beta{\lambda^{\nicefrac{{1}}{{2}}}_{s}}\Delta_{s,t}^{2}{\lambda^{\nicefrac{{1}}{{2}}}_{t}}}{2\sqrt{(\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}\mathbf{x})_{s}}}\operatorname{\mathds{E}}\eta
f^{\prime}(\beta\eta\sqrt{(\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}\mathbf{x})_{s}}+h)=\frac{\beta^{2}}{2}{\lambda^{\nicefrac{{1}}{{2}}}_{s}}\Delta_{st}^{2}{\lambda^{\nicefrac{{1}}{{2}}}_{t}}\operatorname{\mathds{E}}f^{\prime\prime}(\beta\eta\sqrt{(\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}\mathbf{x})_{s}}+h)$

and the last equality follows by Gaussian integration by parts. Now, note that

$f^{\prime}(y)=2\cdot\frac{\tanh(y)}{\cosh^{2}(y)},\
f^{\prime\prime}(y)=2\cdot\frac{1-2\sinh^{2}(y)}{\cosh^{4}(y)}\in[-2,2].$

In particular, we have

$\displaystyle
J({\Psi})=\beta^{2}A\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}$

where $A=\operatorname{diag}(a_{1},a_{2},\ldots,a_{m})$ and
$a_{s}:=\frac{1}{2}\operatorname{\mathds{E}}f^{\prime\prime}(\beta\eta\sqrt{(\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}\mathbf{x})_{s}}+h)\in[-1,1]$
for $s=1,2,\ldots,m$. Thus

$\|J(\Psi)\|\leqslant\beta^{2}\|\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}}\|=\beta^{2}\rho(\Lambda^{\nicefrac{{1}}{{2}}}\Delta^{2}\Lambda^{\nicefrac{{1}}{{2}}})=\beta^{2}\rho(\Delta^{2}\Lambda).$

In particular, $\beta^{2}<\frac{1}{\rho(\Delta^{2}\Lambda)}$ implies that
$\|J({\Psi})\|<1$, _i.e.,_ ${\Psi}$ is a contraction and the system (78) has a
unique solution. $\blacksquare$

Next, we present an elementary approach to prove Proposition 1.6, which
concerns the uniqueness of $\mathbf{q}_{\star}$ for _indefinite_ $\Delta^{2}$
with $m=2$.

###  7.2. Proof of Proposition 1.6

We use an elementary approach to prove Proposition 1.6. Recall
$\mathbf{Q}=A\mathbf{q}$, where $A=\Delta^{2}\Lambda$. For _indefinite_
$\Delta^{2}$, we have $\det(A)<0$, then

$A^{-1}=\begin{pmatrix}-a&b\\\ c&-d\end{pmatrix}\quad\text{where}\quad
a,b,c,d>0.$

Then we can rewrite the equation,

$\begin{cases}q_{1}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{Q^{1}}+h),\\\
q_{2}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{Q^{2}}+h).\end{cases}$

as

$\displaystyle\begin{cases}-aQ^{1}+bQ^{2}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{Q^{1}}+h),\\\
cQ^{1}-dQ^{2}=\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{Q^{2}}+h).\end{cases}$
(80)

Let
$g_{1}(x):=\frac{1}{b}(\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{x}+h)/x+a)$
and
$g_{2}(x):=\frac{1}{c}(\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{x}+h)/x+d)$,
by the classical Latala-Guerra lemma [Tal11a]*Appendix A.14,
$g_{1}(x),g_{2}(x)$ are both strictly decreasing on $\mathbb{R}^{+}$. The
equation (80) can be rewritten further:

$\displaystyle\begin{cases}Q^{1}=Q^{2}g_{2}(Q^{2}),\\\
Q^{2}=Q^{1}g_{1}(Q^{1}).\end{cases}$ (81)

From this, we get

$g_{1}(Q^{1})g_{2}(Q^{2})=1,$

then $Q^{2}=g_{2}^{-1}\left(\frac{1}{g_{1}(Q^{1})}\right)$ is a decreasing
function of $Q^{1}$. On the other hand, taking cross product of equation (81),
we have

$\displaystyle(Q^{2})^{2}g_{2}(Q^{2})=(Q^{1})^{2}g_{1}(Q^{1}).$ (82)

Let $h_{i}(x)=x^{2}g_{i}(x),i=1,2.$ Next we show $h_{1},h_{2}$ are both
increasing functions on $\mathbb{R}^{+}$. Since

$\displaystyle h_{1}(x)=x^{2}g_{1}(x)$
$\displaystyle=\frac{x\operatorname{\mathds{E}}\tanh^{2}(\beta\eta\sqrt{x}+h)+ax^{2}}{b}$
$\displaystyle=\int_{-\infty}^{\infty}\frac{\sqrt{x}}{b}\tanh^{2}(y)\frac{1}{\sqrt{2\pi\beta^{2}}}e^{-\frac{(y-h)^{2}}{2\beta^{2}x}}dy+\frac{ax^{2}}{b},$

where we expand the expectation part as a Gaussian integral. For this integral
form, it’s easy to check $h_{1}(x)$ is an strictly increasing function of
$x\in\mathbb{R}^{+}$. Similarly, we can prove $h_{2}$ is also strictly
increasing. From (82), we have

$h_{1}(Q^{1})=h_{2}(Q^{2}).$

Thus $Q_{2}=h_{2}^{-1}\left(h_{1}(Q^{1})\right)$ is strictly increasing as a
function of $Q^{1}$. Recall
$Q_{2}=g_{2}^{-1}\left(\frac{1}{g_{1}(Q^{1})}\right)$ is a decreasing function
of $Q^{1}$. Therefore, we conclude the uniqueness of $Q^{1}$, similarly for
$Q^{2}$. The proof is complete. $\blacksquare$

In the following part, we will prove Theorem 1.15. As we mentioned before, we
assume $\Delta^{2}$ is _positive-definite_ because the Parisi formula is only
known in that case. We further assume the uniqueness of $\mathbf{q}_{\star}$
when $h>0$, in order for accurate and correct characterization of AT line
condition.

Consider $\mathbf{p}\in[0,1]^{m}$ with $p^{s}\geqslant q_{\star}^{s}$ for all
$s=1,2,\ldots,m$, define

$\displaystyle