\section{Introduction}\label{sect:intro}

This paper investigates restricted Boltzmann machines (RBMs) from the point of view of Hebbian learning. RBMs are widely studied generative models of machine learning, introduced long ago in \cite{smo} and their relation with models of associative memory, such as the ones proposed by Little \cite{little} or Hopfield \cite{hop}, was noted at the very early stage of the theory, see \cite{hinton}. Simple in formulation yet very rich in structure, RBMs represent a nice prototype of energy model, on which theoretical progresses could be made and new ideas could be tested.

Here we address the following natural question from the viewpoint of associative networks: is it possible to reconstruct a pattern just looking at the typical configurations of the visible layer of a RBM? Since the hidden prior determines an energy function on the visible layer, see \cite{tubiana}, \cite{PRE}, \cite{bulso}, we can analyse the corresponding associative memory, for which the natural term of comparison is the Hopfield model, obtained by a Bernoulli-Gaussian RBM via marginalisation over the hidden layer. The question is relevant for the understanding of the RBM probability distribution and in particular of how the hidden layer affects the configurations of the visible units. Nowadays the main focus on RBMs is on learning the probability distribution underlying a given dataset. In a sense, this amounts to fitting the true law of the data by tuning the weights in the density of the model, in such a way that the datapoints fall close to deep local minima of the energy function. Therefore the analysis of learning in RBMs consists of two tasks: understanding the complex landscape of the energy in the vicinity of the datapoints at given weights and devising good optimisation algorithms to fit the data. This steps, at least from a mathematical angle, are very challenging and the problem is widely open. In these respects pattern retrieval represents a simplified setting to study at first instance, as the roles of the datapoints and of the weights is undertaken by the same objects, the patterns. So there is no optimisation, but one only has to look at the energy landscape in the vicinity of the patterns.

In fact we study the retrieval of i.i.d. binary patterns with RBMs with binary visible units, as the distribution of the hidden layer varies. We do it looking at the local minima of the energy function, in what is called in statistical physics a zero temperature limit, as an analogous rigorous analysis of the RBM density (see (\ref{eq:p}) below) is still prohibitive. Anyway for what concerns pattern reconstruction the zero-temperature regime is particularly relevant, since it is exactly where retrieval is maximised (see \cite{amit}).
We prove that the tail of the hidden prior distribution determines the retrieval capability of RBMs. More precisely, for priors with tails ranging from exponential to Gaussian we prove that deep local minima are well localised about the patterns. On the other hand, if the tails of the hidden prior decay faster than Gaussian, we show that the patterns cannot be retrieved well in any case. RBMs whose hidden priors have Gaussian tails (a class including the Hopfield model) represent special threshold cases which we treat separately in either the positive (Theorem \ref{TH1}) and the negative (Theorem \ref{TH2}) result below.

\subsection{Setting}

We consider RBMs with i.i.d. Bernoulli $\pm1$ visible units $\sigma_1,\ldots,\sigma_{N_1}$ and symmetric i.i.d. hidden units $z_1,\ldots, z_{N_2}$ and we allow a certain freedom in the choice of the hidden prior $\pi$, for which we only require that
\begin{equation}\label{eq:priorZ}
\pi(|z|\geq t)\simeq e^{-|t|^{q}}\qquad\mbox{for some $q>1$}\,.
\end{equation}

Let $\xi^{(1)},\ldots, \xi^{(N_2)}$ denote independent patterns with $N_1$ i.i.d. $\pm1$ components. The (unnormalised) density of the RBM probability distribution reads
\begin{equation}\label{eq:p}
p(\sigma,z;\xi):= \exp\left(\frac{\beta}{N_1^{\frac 1{q_-}}}\sum_{\mu\in[N_2]}(\sigma,\xi^{(\mu)}) z_\mu\right)\,,
\end{equation}
where $\beta>0$ is a parameter usually called {\em inverse temperature} and $q_-:=\min(q,2)$. We shall consider the ratio between the number of visible and hidden units as follows
\begin{equation}\label{eq:alpha}
\alpha:=\frac{N_2}{N_1^\frac{p_+}{2}}\,,
\end{equation}
where $\frac1{p_+}+\frac1{q_-}=1$.
The normalisation factors in (\ref{eq:p}) and (\ref{eq:alpha}) ensure that either the energy of the single pattern and the global maximum of (\ref{eq:p}) as $\beta\to\infty$ scale linearly with $\alpha$ (with constants depending on the prior distribution); more details on that in Section \ref{sect:T0}.

\subsection{Main results}

The main object of interest in this paper will be the following two sets. These definitions are slightly formal and a more precise characterisation of it will be given in the next section.

\begin{eqnarray}
\mathtt{LM}_{N_1}&:=&\{\mbox{local maxima of (\ref{eq:p}) as $\beta\to\infty$}\}
\label{eq:LMset}\,,\\
\mathtt{dLM}_{N_1}^{(\mu)}&:=&\{\mbox{local maxima of (\ref{eq:p}) as $\beta\to\infty$ higher than the value at $\xi^{(\mu)}$}\}\label{eq:DLMset}\,.%
\end{eqnarray}

Below $\operatorname{Hamming}(a,b)$ denotes the Hamming distance between $a,b$ (\textit{i.e. } the number of different entries) and $\widehat B^{N_1}_{\mu,R}$ is the ball in this metric centred at the $\mu$-th pattern with radius $R$. Throughout the paper we will repeatedly use that two patterns are typically separated by $N_1/2$ flips, so that $\widehat B^{N_1}_{\mu,\lfloor N_1/2\rfloor}$ and $\widehat B^{N_1}_{\mu',\lfloor N_1/2\rfloor}$ typically do not overlap.

Our first result states that the error of reconstructing a given pattern is very small in terms of the size of the network if the decay of the hidden prior (\ref{eq:priorZ}) is slower than Gaussian, while a finite fraction of bits cannot be retrieved for $q=2$.

\begin{theorem}\label{TH1}
Let $q\in(1,2)$. There exists $r_0\in(0,\frac12]$ such that with high probability
\begin{equation}\label{eq:TH1-p>2}
\max_{\mu\in[N_2]}\max_{\sigma\in \mathtt{dLM}_{N_1}^{(\mu)}\cap \widehat B^{N_1}_{\mu,\lfloor r_0N_1\rfloor}}\operatorname{Hamming}\left(\sigma, \xi^{(\mu)}\right)\lesssim (\log N_1)^{\frac{q-1}{2-q}}
\end{equation}
Let $q=2$. There exists $r_0\in(0,\frac12]$ such that for all $r\in(0,r_0)$ if $\alpha\lesssim\sqrt r$ then with high probability
\begin{equation}\label{eq:TH1-p=2}
\max_{\mu\in[N_2]}\max_{\sigma\in \mathtt{dLM}_{N_1}^{(\mu)}\cap \widehat B^{N_1}_{\mu,\lfloor r_0N_1\rfloor}}\operatorname{Hamming}\left(\sigma, \xi^{(\mu)}\right)\lesssim \lfloor rN_1\rfloor\,.
\end{equation}
\end{theorem}

Theorem \ref{TH1} characterises the models with $q<2$ as excellent in pattern reconstruction: there are high maxima located few flips away from the patterns (in fact polylog flips, see (\ref{eq:TH1-p>2})) and no higher maxima appear in an extended region. Specifying the radius $r_0$ and the basin of attraction of the patterns would require a careful optimisation of the constants in the proof that we do not attempt.

Local maxima are not directly related to the typical configurations of (\ref{eq:p}). However it is well known that any algorithmic search of typical configurations will finish to find a hopefully representative local maximum. This can be done by the usual {\sc flip} algorithm, that is greedy flipping of one units at time decreasing the energy until no more decreasing is possible. Therefore $\mathtt{dLM}_{N_1}^{(\mu)}$ has a direct interpretation in terms of retrieval. Take for instance for $q<2$. By Theorem \ref{TH1} any {\sc flip} search initialised for instance at $\xi^{(\mu)}$ will end up in a point of $\mathtt{dLM}_{N_1}^{(\mu)}$ falling only $(\log N_1)^{\frac{q-1}{2-q}}$ flips away from the pattern, which means that only few bits are mis-retrieved.

Somewhat in the opposite direction, the next result shows that for $q>2$ in (\ref{eq:priorZ}), the local maxima of the density as $\beta\to\infty$ are quite far from the patterns.

\begin{theorem}\label{TH2}
Let $q\geq2$, $r\in[0,\frac12]$ and let $\alpha_q(r)=\sqrt{S(r)}$ for $q\neq2$ and $\alpha_2(r):=S(r)/(1-2r)^2$. For every $r\in(0,\frac12)$ and for all $\alpha\gtrsim\alpha_q(r)$ we have with high probability that
\begin{equation}\label{eq:TH2}
\operatorname{Hamming}\left(\mathtt{LM}_{N_1}, \xi^{(\mu)}\right)\gtrsim \lfloor rN_1\rfloor\,.
\end{equation}
\end{theorem}

According to Theorem \ref{TH2}, if $q>2$ one could still hope for retrieval with a very small amount of hidden variables, \textit{i.e. } for $\alpha$ small enough (indeed for $\alpha=0$ reconstruction is possible, see \cite{PRE}), but for $\alpha$ larger than a given constant no recovery is allowed. For $q=2$ the situation improves a bit in the sense that pattern reconstruction becomes less and less efficient as $\alpha$ grows.

The paper \cite{PRE} showed that RBMs with hidden prior interpolating between a Gaussian and a bimodal symmetric distribution exhibit retrieval at finite $\alpha>0$, which disappears in the degenerate case when the Gaussian part is switched off. It is also argued that such a lack of retrieval should persist at least for any compactly supported hidden prior. This is demonstrated using non-rigorous replica computations and numerics. We give here the first mathematical confirm of these findings, as Theorem \ref{TH1} implies patter retrieval if in the above interpolating prior the Gaussian part is present, whatever small, and Theorem \ref{TH2} is a strong indication for lack of retrieval for hidden prior with a Bernoulli $\pm1$ distribution (for which we should read $q=\infty$).

\subsection{Related literature}

The results here presented mark a neat difference in the retrieval capabilities of RBMs with hidden priors (\ref{eq:priorZ}) with $q<2$ (very good capabilities) and $q>2$ (not so good) with a transition at the sub-Gaussian case $q=2$. As already remarked, a notable instance of the case $q=2$ is the Hopfield model, for which a similar analysis at zero temperature was done in \cite{newman} (analog of Theorem \ref{TH1}), \cite{lou} (analog of Theorem \ref{TH2}) and \cite{talahop} in the attempt of proving the picture of \cite{amit}. When comparing these papers to ours, we underline that due to the generality of the models we treat, we do not seek to characterise our estimate with the best possible constants (we leave all the numerical constants in fact undetermined), which was instead a very relevant component of all these previous papers.

We exploit and make mathematically precise the heuristics of \cite{PRE}. Namely we use that the tail of the hidden prior determines the behaviour for large argument of the energy function of the associative network (around zero it is always quadratic). It is exactly this asymptotic that governs retrieval: the more convex the better. More details on that in Section \ref{sect:T0}.
Mathematically speaking the introduction of the hidden layer is a way to linearise the energy function (over the visible units) and different prior distributions for the hidden layer correspond to different associative networks. Similar ideas have been used by \cite{tubiana}, \cite{PRE}, \cite{cocco}, \cite{bulso}, \cite{treves} to study the performance of the RBMs with varying hidden unit statistics.

The {\sc flip} algorithm is a very natural choice to explore the energy landscape of RBMs and indeed it is essentially the original network dynamics proposed by \cite{hop}. This gives a nice connection with the local max-cut problem as analysed for instance by \cite{maxcut2} and \cite{maxcut1}, even though here we exploit crucially the presence of the patterns, which constitute a special class of local minima. This is even more clear by comparing with the analysis for the Sherrington-Kirkpatrick model of \cite{SK}.

Many other dynamics have been proposed alternative to the {\sc flip} algorithm mainly for the Hopfield model and we will not give here an account on that (see the recent work \cite{zamponi} and the references therein). We just mention that the dynamics analysed by \cite{marinari}, which is a zero-temperature version of the alternate Gibbs sampling typically used to train RBMs, is in spirit very close to our zero-temperature reduction.

\subsection{Notations}

Throughout the paper $p,q\geq1$ will always be H\"older conjugate, that is $\frac 1p+\frac1q=1$ and similarly for $q_-, p_+$, with $q_-:=\min(2,q)$, $p_+:=\max(2,p)$.
$C,c$ everywhere denote positive absolute constants which may change from formula to formula. We write $X\lesssim Y$ if $X\leq CY$ and $X\simeq Y$ if $Y\lesssim X\lesssim Y$. Sometimes we write $\simeq_a$ or $\lesssim_a$ to stress the dependence of the constants $C$ above on a parameter $a$. We indicate by $(\cdot,\cdot)$ the inner product in $\mathbb{R}^{N_1}$ or $\mathbb{R}^{N_2}$ and the meaning will be always clear from the context and by $\|\cdot\|_p$ the $\ell_p$-norms. $1$ may represent the vector in $\mathbb{R}^{N_1}$ or in $\mathbb{R}^{N_2}$ with all entries equal to $1$. $B_{N}^{(q)}$ is the $\ell_q$ centred ball of radius one in $\mathbb{R}^N$. $\widehat S^{N_1-1}_{\mu,R}, \widehat B^{N_1}_{\mu,R}$ denote respectively the $N_1$-dimensional Hamming ball and sphere centred at $\xi^{(\mu)}$ of radius $R$. If $v\in\mathbb{R}^N$ and $J\subset [N]$ we denote by $v_J$ a vector in $\mathbb{R}^{|J|}$ such that $(v_J)_i=v_{j_i}$ if $J=\{j_1,\ldots, j_{|J|}\}$. For any $J\subset [N]$ we also associate a {\sc flip} operator $F_J$ defined by $(F_J v)_i=-v_i$ if $i\in J$ and $(F_J v)_i=v_i$ if $i\notin J$.
We will use the following $\psi_r$-norms:
\begin{equation}
\|Z\|_{\psi_r}:=\inf\left\{c>0\,:\,E\left[\psi_r\left(\frac{|Z|}{c}\right)\right]<2\right\}\,,\qquad r>0\,,
\end{equation}
where $\psi_r(x)=e^{x^r}$ for any $x>0$ for $r\geq1$, while for $r\in(0,1)$ there are $c(r),x(r)$ such that for $x\in(0,x(r))$ it is $\psi_r(x)=c(r)x$.
We underline that, setting $q_s:=\sup\{q'>1\,:\,\|Z\|_{\psi_{q'}}<\infty\}\,$
we have $P(|Z|\geq t)\simeq e^{-|t|^q_s}$ (we convey that bounded random variables have finite $\psi_\infty$-norm).
We will often deal with the following function:
$$
S(r):=-r\log r-(1-r)\log(1-r)\,
$$
and use the standard bound
\begin{equation}\label{eq:standard-bound}
\operatorname{Card} \widehat S^{N_1-1}_{\mu,\lfloor rN_1\rfloor}= \binom{\lfloor rN_1\rfloor}{N_1}\leq e^{N_1S(r)}\,.
\end{equation}
We denote the transpose patterns $\tilde\xi^{(i)}$ by $\tilde\xi^{(i)}_\mu:=\xi_i^{(\mu)}$, $i\in[N_1]$, $\mu\in[N_2]$. Sometimes we write $\widehat \xi:=\xi/\sqrt{N_1}$.
$A^c$ is the complement of the set $A$. We say that the event $A$ occurs with high probability if $P(A)\geq 1-N^{-x}$ for all $x>0$.

\subsection*{Acknowledgements} The author thanks David Belius for many valuable suggestions regarding the presentation of the results.

\section{Zero temperature reduction}\label{sect:T0}

Here we show how to single out an effective energy function which captures the RBM landscape at zero temperature. We exploit a reduction argument introduced in \cite{PRE}.

Integrating out the hidden layer in (\ref{eq:p}) we get a probability distribution over the visible units. This writes as
\begin{equation}
\int p(\sigma,z;\xi) \pi(dz_1)\ldots \pi(dz_{N_2})=\exp\left(\sum_{\mu\in[N_2]}u\left(\frac\beta{N_1^{\frac 1q}}(\xi^{\mu},\sigma)\right)\right)\,,
\end{equation}
where
$$
u(x):=\log E[e^{xz_1}]\,.
$$
We are interested in studying the local minima of this distribution as $\beta$ is very large. The following simple observation starts our considerations: for any $z_1$ such that \eqref{eq:priorZ} holds for some $q>1$, it is
\begin{equation}\label{eq:start}
\lim_{x\to\infty} \frac{u(x)}{|x|^{p}}\simeq_p \|z_1\|_{\psi_q}^p\,.
\end{equation}
The proof is immediate: we write
\begin{eqnarray}
E[e^{xz_1}]&=&\int_0^\infty d\lambda P(z_1\geq x^{-1}\log\lambda)\simeq\int_0^\infty d\lambda e^{-\frac{|\log\lambda|^q}{\|z_1\|^q_{\psi_q}|x|^q}}\nonumber\\
&=&\|z_1\|_{\psi_q}|x|\int_0^\infty d\lambda e^{\lambda \|z_1\|_{\psi_q}|x|-|\lambda|^q}\nonumber\\
&=&\|z_1\|_{\psi_q}|x|e^{\frac{\|z_1\|_{\psi_q}^p|x|^p}{p}}\int_0^\infty d\lambda e^{\lambda\|z_1\|_{\psi_q} x-|\lambda|^q-\frac{\|z_1\|_{\psi_q}^p|x|^p}{p}}\nonumber\,,
\end{eqnarray}
and the last integral is finite by the Young inequality. Taking the $\log$ on both sides and passing to the limit we get (\ref{eq:start}).

Therefore by \eqref{eq:start}, as $\beta\to\infty$ we are led to consider the following energy function (reminiscent of the $p$-spin Hopfield model of \cite{bov-p}, but with a different normalisation)
\begin{equation}\label{eq:Hp}
H^{(p)}(\sigma;\xi):=-\frac{1}{N_1^{\kappa(p)}}\sum_{\mu\in[N_2]}|(\xi^{(\mu)},\sigma)|^p\,,
\end{equation}
with
$
\kappa(p):=1+p-\frac{p}{p_+}\,
$
(here we include the usual normalisation factor $1/N_1$ of the internal energy directly in the definition of $H^{(p)}$).
To fix the ideas, $\kappa(p)=p$ for $p\geq 2$ and $\kappa(p)=1+\frac{p}{2}$ for $p\leq2$.

The focus of this paper will be to study the location of the minima of \eqref{eq:Hp} on $\{-1,1\}^{N_1}$ close to the pattern configurations, in the limit $N_1,N_2\to\infty$ while $\alpha_p$ remains constant. Hence the following are more precise definitions of the sets (\ref{eq:DLMset}) and (\ref{eq:LMset}) previously introduced:

\begin{eqnarray}
\mathtt{LM}_{N_1}&:=&\{\mbox{local minima of \eqref{eq:Hp}}\}\label{eq:LMset-1}\,,\\
\mathtt{dLM}_{N_1}^{(\mu)}&:=&\{\mbox{local minima $\bar\sigma$ of (\ref{eq:Hp}), s.t. $H^{(p)}(\xi^{(\mu)};\xi)-H^{(p)}(\bar\sigma;\xi)>0$}\}\label{eq:DLMset-1}\,.
\end{eqnarray}

Next we show that the single pattern energy is close to the ground state, so providing a motivation for the normalisation factors in (\ref{eq:p}) and (\ref{eq:alpha}). First of all we give a lower bound for the ground state energy linear in $\alpha$. To do so we do not actually need binary patterns.

\begin{proposition}\label{Prop:groundstate}
Let $\xi^{(1)},\ldots, \xi^{(N_2)}$ be independent vectors in $\mathbb{R}^{N_1}$ with i.i.d. centred sub-Gaussian entries. It holds for $p\geq1$
\begin{equation}\label{eq:notfar}
\inf_{\sigma\in \{-1,1\}^{N_1}} H^{(p)}(\sigma;\xi)\gtrsim_p -(1+\alpha)
\end{equation}
with probability larger than $1-e^{-c\alpha^{\frac2p}N_1}$.
\end{proposition}
To prove Proposition \ref{Prop:groundstate} we need the following auxiliary lemma.

\begin{lemma}\label{lemma:Psup}Let $\xi^{(1)},\ldots, \xi^{(N_2)}$ be independent vectors in $\mathbb{R}^{N_1}$ with i.i.d. centred sub-Gaussian entries.
Let $p\geq1$. For all $t\gtrsim_p(1+\alpha)^{\frac1p}$
\begin{equation}\label{eq:Psup}
P\left(\frac{1}{N_1^{\max(0,\frac{q-2}{2q})}}\sup_{\sigma\in \frac{1}{\sqrt N_1}\{-1,1\}^{N_1}} \sup_{\tau\in B_{N_2}^{(q)}} \frac{(\xi^{(\mu)},\sigma)\tau_{\mu}}{\sqrt{N_1}}\geq t\right)\leq 2e^{-ct^2N_1}\,
\end{equation}
where $c>0$ depends only on the distribution of $\xi_1^{(1)}$.
\end{lemma}

\begin{proof}
We introduce the transpose patterns $\tilde\xi^{(i)}$ by $\tilde\xi^{(i)}_\mu:=\xi_i^{(\mu)}$, $i\in[N_1]$, $\mu\in[N_2]$.
First of all we note that
\begin{equation}\label{eq:notethat}
\frac{1}{N_1^{\max(0,\frac{q-2}{2q})}}\sup_{\sigma\in \frac{1}{\sqrt N}\{-1,1\}^{N_1}} \sum_{i\in[N_1]}\sum_{\mu\in[N_2]}\frac{\xi_i^{(\mu)}\sigma_i\tau_\mu}{\sqrt{N_1}}=\frac1{N_1}\sum_{i\in[N_1]}\frac{\left|(\tilde\xi^{(i)},\tau)\right|}{N_1^{\max(0,\frac{q-2}{2q})}}\,.
\end{equation}
Moreover, since for all $\tau\in B_{N_2}^{(q)}$
$$
\|\tau\|_2\leq N_2^{\max(0,\frac{q-2}{2q})}\,,
$$
we have
\begin{equation}
P\left(\frac{1}{N_1^{\max(0,\frac{q-2}{2q})}}\left|(\tilde\xi^{(i)},\tau)\right|\geq t\right)\leq 2e^{-\frac{t^2N^{\max(0,\frac{q-2}{q})}_1}{2N_2^{\max(0,\frac{q-2}{q})}\|\xi^{(1)}_1\|^2_{\psi_2}}}\,.
\end{equation}
The r.h.s. of (\ref{eq:notethat}) is the sum of independent sub-Gaussian random variables, thus for all $\tau\in B_{N_2}^{(q)}$
\begin{equation}
P\left(\frac1{N_1}\sum_{i\in[N_1]}\frac{\left|(\tilde\xi^{(i)},\tau)\right|}{N_1^{\max(0,\frac{q-2}{2q})}}\geq t\right)\leq 2e^{-\frac{t^2N^{1+\max(0,\frac{q-2}{q})}_1}{2N_2^{\max(0,\frac{q-2}{q})}\|\xi^{(1)}_1\|^2_{\psi_2}}}\,
\end{equation}
for some $c>0$ depending only on the distribution of $\xi_1^{(1)}$.

Next we cover $B^{(q)}_{N_2}$ with a number of $\ell_r$-balls in $\mathbb{R}^{N_2}$ with some small radius $\varepsilon>0$. For $p\geq2$ we can use Euclidean balls and the Sudakov inequality gives a bound on the minimal number $N(B^{N_2}_q,\varepsilon B^{N_2}_2)$ of such balls
$$
N(B^{N_2}_q,\varepsilon B^{N_2}_2)\leq e^{cN_2^{\frac2p}}\,
$$
(here we used that for a gaussian vector $g$ $E[\max_{\tau\in B^{N_2}_q}(\tau,g)]=E[\|g\|_p]\simeq N_2^{\frac1p}$).

For $p\in[1,2)$ we use $\ell_q$-balls and we have
$$
N(B^{N_2}_q,\varepsilon B^{N_2}_q)\leq e^{cN_2}
$$
(in the two estimates above the constants $c$ depends on $\varepsilon$ in a way we do not keep track of).

Assume now $(1+\alpha)^{\frac1p}\lesssim_p t$ (this is to take into account also the behaviour for small $\alpha$). By the union bound
for $p\geq2$ we get
\begin{equation}
P\left(\frac{1}{N_1^{\max(0,\frac{q-2}{2q})}}\sup_{\sigma\in \frac{1}{\sqrt N}\{-1,1\}^{N_1}} \sup_{\tau\in B_{N_2}^{(q)}} (\widehat\xi\sigma,\tau)\geq t\right)\leq 2e^{cN^{\frac2p}_2-cN_1t^2}\leq 2e^{-ct^2N_1}\,.
\end{equation}

Similarly for $p\in[1,2)$
\begin{equation}
P\left(\frac{1}{N_1^{\frac{q-2}{2q}}}\sup_{\sigma\in \frac{1}{\sqrt N}\{-1,1\}^{N_1}} \sup_{\tau\in B_{N_2}^{(q)}} (\widehat\xi\sigma,\tau)\geq t\right)\leq 2e^{cN_2-c\left(\frac{N_1}{N_2}\right)^{\frac{q-2}{q}}N_1t^2}\leq 2e^{-ct^2N_1}\,.
\end{equation}
\end{proof}

\begin{proof}[Proof of Proposition \ref{Prop:groundstate}]
The role of hidden variables at zero temperature is played by duality:
\begin{equation}\label{eq:duality}
\sum_{\mu\in[N_2]}|(\xi^{(\mu)},\sigma)|^p= \left|\sup_{\tau\in B_{N_2}^{(q)}} \sum_{\mu\in[N_2]}(\xi^{(\mu)},\sigma)\tau_\mu\right|^p\,,\qquad p\geq 1\,.
\end{equation}
Therefore (here we shorten $\wh\xi:=\xi/\sqrt {N_1}$)
\begin{eqnarray}
\inf_{\sigma\in \{-1,1\}^{N_1}} H^{(p)}(\sigma;\xi)&=&-\sup_{\sigma\in \frac{1}{\sqrt {N_1}}\{-1,1\}^{N_1}}\left|\frac{1}{N_1^{\max(0,\frac{q-2}{2q})}} \sup_{\tau\in B_{N_2}^{(q)}} \sum_{\mu\in[N_2]}(\widehat\xi^{(\mu)},\sigma)\tau_\mu\right|^p\nonumber\\
&=&- \left|\frac{1}{N_1^{\max(0,\frac{q-2}{2q})}}\sup_{\sigma\in \frac{1}{\sqrt {N_1}}\{-1,1\}^{N_1}}\sup_{\tau\in B_{N_2}^{(q)}} \sum_{\mu\in[N_2]}(\widehat\xi^{(\mu)},\sigma)\tau_\mu\right|^p\,\label{eq:bysymmetry}
\end{eqnarray}
by symmetry.
It suffices to focus on the quantity inside the modulus above, which is dealt in Lemma \ref{lemma:Psup} of Appendix \ref{Sect:inf} below. We have for all $t\geq_p(1+\alpha)^{\frac1p}$
\begin{equation}\label{eq:Psup-bis}
P\left(\frac{1}{N_1^{\max(0,\frac{q-2}{2q})}}\sup_{\sigma\in \frac{1}{\sqrt N}\{-1,1\}^{N_1}} \sup_{\tau\in B_{N_2}^{(q)}} (\widehat\xi\sigma,\tau)\geq t\right)\leq 2e^{-ct^2N_1}\,.
\end{equation}
Combining (\ref{eq:bysymmetry}) with the bound above
we obtain the statement.
\end{proof}

Now we shortly show that the patterns have energy of the same order in $\alpha$ of the global minimum, even though we can already observe a difference between the models with $p\geq 2$ and $p<2$. We deal with $\pm1$ binary patterns for simplicity, but a similar argument can be easily repeated for symmetric patterns with minor modifications. We have
$$
H^{(p)}(\xi^{(1)};\xi)=-\frac{\|\xi^{(1)}\|_2^p}{N_1^{1+p-p/p_+}}-\sum_{j=2}^{N_2}\frac{|(\xi^{(1)},\xi^{(j)})|^p}{N_1^{1+p-p/p_+}}\,.
$$
This quantity concentrates around its average as $N_1$ grows:
\begin{lemma}\label{lemma:Hconcentra}
Take $t>0$ uniformly in $N_1$, small enough. It holds
\begin{equation}\label{eq:Hconcentra}
P\left(\left|H^{(p)}(\xi^{(1)};\xi)-E[H^{(p)}(\xi^{(1)};\xi)]\right|\geq t\right)\leq  \begin{cases}
\exp\left(-cN_1 t^2\right)&p\in[1,2]\\
\exp\left(-cN_1 t^{\frac{2}{p}}\right)&p>2
\end{cases}
\end{equation}
\end{lemma}
\begin{proof}
We write
$$
H^{(p)}(\xi^{(1)};\xi)=-\frac{\|\xi^{(1)}\|_2^p}{N_1^{1+p-p/p_+}}-\sum_{\mu=2}^{N_2}\frac{|(\xi^{(1)},\xi^{(\mu)})|^p}{N_1^{1+p-p/p_+}}\,.
$$
It suffices to focus on the second summand on the r.h.s. above. We have by independence
\begin{eqnarray}
&&P\left(\left|\sum_{\mu=2}^{N_2}\frac{|(\xi^{(1)},\xi^{(\mu)})|^p}{N_1^{\frac p2}}-E\left[\frac{|(\xi^{(1)},\xi^{(\mu)})|^p}{N_1^{\frac p2}}\right]\right|\geq tN_1^{1+p/2-p/p_+}\right)\nonumber\\
&=&P\left(\left|\sum_{\mu=2}^{N_2}\frac{|(1,\xi^{(\mu)})|^p}{N_1^{\frac p2}}-E\left[\frac{|(1,\xi^{(\mu)})|^p}{N_1^{\frac p2}}\right]\right|\geq tN_1^{1+p/2-p/p_+}\right)\,,\label{eq:Mpestimato}
\end{eqnarray}
where $1$ is the constant vector with all entries equal to $1$. The random variables
$$
T_p^{(\mu)}:=\frac{|(1,\xi^{(\mu)})|^p}{N_1^{\frac p2}}-E\left[\frac{|(1,\xi^{(\mu)})|^p}{N_1^{\frac p2}}\right]\,,\qquad \mu=2,\ldots, N_2
$$
are i.i.d. with $\|T_p^{(\mu)}\|_{\psi_{\frac 2p}}\simeq 1$. Thus by Proposition \ref{prop:Hoeffging-GEN} with $\ell=2/p$ we have
\begin{equation}
\eqref{eq:Mpestimato}\leq
\begin{cases}
\exp\left(-cN_1\min(t^2, t^{\frac{2}{p}})\right)&p\in[1,2]\\
\exp\left(-cN_1\min(t^2N_1^{p-2}, t^{\frac{2}{p}})\right)&p>2
\end{cases}
\end{equation}
and the proof is complete.
\end{proof}

Moreover we have
\begin{equation}\label{eq:E-Hp-xi1}
|E[H^{(p)}(\xi^{(1)};\xi)]|\leq \frac{N_1^p+N_1^{\frac p2}N_2}{N_1^{1+p-p/p_+}}\simeq -\alpha\frac p2\Gamma\left(\frac p2\right)-\frac{1}{N_1^{1-p/p_+}} \,.
\end{equation}
In fact
\begin{eqnarray}
E[|(\xi^{(1)},\xi^{(\mu)})|^{p}]&=&E[\int d\lambda P(\xi^{(\mu)}\,:\,|(\xi^{(1)},\xi^{(\mu)})|\geq \lambda^{1/p})]\nonumber\\
&\leq& E[\|\xi^{(1)}\|_1^p]\int_0^\infty d\lambda e^{-\lambda^{\frac2p}}=\frac p2\Gamma\left(\frac p2\right) N_1^{\frac p2}\label{eq:proceeding}\,.
\end{eqnarray}
Therefore
$$
H^{(p)}(\xi^{(1)};\xi)\simeq -\frac p2\Gamma\left(\frac p2\right)\alpha-\frac{1}{N_1^{1-p/p_+}}
$$
with very high probability. We see that if $p\geq 2$ this value is really of the same order of the ground state, while if $p\in(1,2)$ for $\alpha$ small and $N_1$ large the patterns have higher energy.

\section{Retrieval for $p\geq2$}\label{sect:retrieval}

In this section we prove Theorem \ref{TH1}. We look at all configurations reachable from $\xi^{(\mu)}$ by $\lfloor rN_1\rfloor$ flips and compare the energy of the pattern with the minimal energy of such configurations. By symmetry of the patters we can reduce to look at $\mu=1$ and we may and will assume that $1$, \textit{i.e. } the vector with entries all equal to one, lies in $\widehat S^{N_1-1}_{1,\lfloor rN_1\rfloor}$.

Without further explanation, we introduce some more notations.
Given $v\in \widehat S^{N_1-1}_{1,\lfloor rN_1\rfloor}$ we set
$$
J_\mu(v)=J:=\{i\in[N_1]\,:\,v_i\neq \xi^{(\mu)}\}\,.
$$
and
\begin{equation}\label{defXY}
X^{(\mu)}(v):=\frac{1}{\sqrt {N_1}}(\xi^{(\mu)}_J,v_J)\,,\quad Y^{(\mu)}(v):=\frac{1}{\sqrt {N_1}}(\xi^{(\mu)}_{J^c},v_{J^c})\,.
\end{equation}
and $X^{(\mu)}(1)=:X^{(\mu)}$, $Y^{(\mu)}(1)=:Y^{(\mu)}$.
We conveniently let
\begin{equation}\label{eqF0}
F_p(x,y):=|x+y|^p-|x-y|^p\,,\qquad \bar F_p(r):=F_p(r,1-r)=1-(1-2r)^p>0\,.
\end{equation}
For all $p\geq1$ it is
\begin{equation}\label{eqF}
|F_p(x,y)|\leq 2^{p}|x|^{\frac p2}|y|^{\frac p2}\,.
\end{equation}
The proof of \eqref{eqF} is easy. By convexity
$$
||x-y|^p-|x+y|^p|\leq 2^p|x|^p\,,\qquad||x+y|^p-|x-y|^p|\leq 2^p|y|^p\,,
$$
therefore squaring $|F_p(x,y)|$ (\ref{eqF}) follows (in fact we have $F_2(x,y)= 4xy$).

We have the following useful representation.

\begin{lemma}\label{lemma:repr-inutile}
Let $r\in(0,\frac12)$, $J\subset[N_1]$ with $|J|=\lfloor rN_1\rfloor$. It is
\begin{equation}\label{eq:formulaF}
H^{(p)}(\xi^{(1)})-H^{(p)}(F_J\xi^{(1)})=-\frac{1}{N_1^{\frac {p_+-p}{2}}} \bar F_p(r)-\frac{1}{N_1^{\frac {p_+}2}}\sum_{\mu=2}^{N_2} F_p(X_J^{(\mu)}(\xi^{(1)}), Y_J^{(\mu)}(\xi^{(1)}))\,.
\end{equation}
\end{lemma}
\begin{proof}
Compute
\begin{eqnarray}
H^{(p)}(v)-H^{(p)}(F_Jv)&=&-\frac{1}{N_1^{\kappa(p)}}\sum_{\mu\in[N_2]}|(\xi^{(\mu)}_J,v_J)+(\xi^{(\mu)}_{J^c},v_{J^c})|^p - |-(\xi^{(\mu)}_J,v_J)+(\xi^{(\mu)}_{J^c},v_{J^c})|^p\nonumber\\
&=&-\frac{1}{N_1^{\kappa(p)-\frac p2}} \sum_{\mu\in[N_2]}F_p(X^{(\mu)}(v), Y^{(\mu)}(v))\,,\nonumber
\end{eqnarray}
by the definitions (\ref{defXY}), (\ref{eqF0}). We have
$
\kappa(p)-\frac p2=1+\frac p2-\frac{p}{p_+}=\frac{p_+}{2}\,.
$
Take now $v=\xi^{(1)}$. An easy computation gives
\begin{equation}
X_J^{(1)}(\xi^{(1)})=\frac{\|\xi^{(1)}_J\|_2^2}{\sqrt N_1}=\frac{|J|}{\sqrt N_1}=r\sqrt{N_1}\,,\quad Y_J^{(\mu)}(\xi^{(1)})=\frac{\|\xi^{(1)}_{J^c}\|_2^2}{\sqrt {N_1}}=\frac{|J^c|}{\sqrt {N_1}}=(1-r)\sqrt{N_1}\,.
\end{equation}
Thus
$$
\frac{1}{N_1^{\frac {p_+}{2}}} F_p(X^{(1)}(\xi^{(1)}), Y^{(1)}(\xi^{(1)}))= \frac{1}{N_1^{\frac {p_+-p}{2}}} F_p(r, 1-r)
$$
and \eqref{eq:formulaF} follows.
\end{proof}

The necessary tail estimates in order to prove Theorem \ref{TH1} are given in the next lemma. We shorten $P_{\xi}(\cdot)=P(\cdot\,|\,\xi)$.

\begin{lemma}\label{lemma:crux1}
Let $r\in(0,\frac12)$, $p\geq 2$, $0<y\leq 1-(1-2r)^p$, $t:= 1-(1-2r)^p-y$.
Take any $\sigma\in \widehat S^{N_1-1}_{1,\lfloor rN_1\rfloor}$.
If $\alpha\lesssim N^{\frac{p-2}{2}}r^{\frac{p-1}{2}}$ then
\begin{equation}\label{eq:stimap<2-1}
P_{\xi^{(1)}}\left(H^{(p)}(\xi^{(1)})-H^{(p)}(\sigma)\geq -t\right)\leq \exp\left(-\frac{cN_1 |y|^{\frac2p}}{\sqrt{r(1-r)}}\right)\,,
\end{equation}
and otherwise
\begin{equation}\label{eq:stimap<2-2}
P_{\xi^{(1)}}\left(H^{(p)}(\xi^{(1)})-H^{(p)}(\sigma)\geq -t\right)\leq \exp\left(-\frac{cN_1^{\frac p2} y^2}{\alpha (r(1-r))^\frac p2}\right)\,.
\end{equation}
\end{lemma}
\begin{remark}
If $r>0$ uniformly in $N_1$ we flip a number of bits proportional to $N_1$ and the tail exponent is sub-linear for all $\alpha$ if $p>2$ and for $\alpha\lesssim \sqrt r$ if $p=2$. A sub-linear number of flips corresponds to take $r\simeq N^{-x}$ for some $x>0$ (modulo log-corrections, see below). In this case we see that if $x\in(\frac{p-2}{p-1},1)$ the sub-exponential estimate still prevails, while otherwise we have a Gaussian tail.
\end{remark}

The following result is used for proving Lemma \ref{lemma:crux1}. We omit the proof, which is standard.

\begin{lemma}\label{lemma:XYsubG}
$\{X^{(\mu)}\}_{\mu\in[N_2]\setminus \{1\}}$ and $\{Y^{(\mu)}\}_{\mu\in[N_2]\setminus \{1\}}$ are independent sub-Gaussian random variables, independent one from each other, with
\begin{equation}\label{eq:XYsubG}
\|X^{(\mu)}\|_{\psi_2}\leq \sqrt{\frac{r}{\log 2}}\,\qquad \|Y^{(\mu)}\|_{\psi_2}\leq \sqrt{\frac{1-r}{\log 2}}\,.
\end{equation}
Moreover $\{F_p(X^{(\mu)}, Y^{(\mu)})\}_{\mu\in[N_2]}$ are i.i.d. $\psi_{2/p}$ r.vs with
\begin{equation}\label{eq:disF}
\|F_p(X^{(\mu)}, Y^{(\mu)})\|^{\frac 2p}_{\psi_{2/p}}\leq \frac{\sqrt{r(1-r)}}{\log 2}\,.
\end{equation}
\end{lemma}

\begin{proof}[Proof of Lemma \ref{lemma:crux1}]
By \eqref{eq:formulaF} we have
\begin{eqnarray}
&&P_{\xi^{(1)}}\left(H^{(p)}(\xi^{(1)})-H^{(p)}(\sigma)\geq -t\right)\nonumber\\&=&P_{\xi^{(1)}}\left(-\sum_{\mu=2}^{N_2} F_p(X_J^{(\mu)}(\xi^{(1)}), Y_J^{(\mu)}(\xi^{(1)}))\geq N_1^{\frac {p}{2}} (\bar F_p(r)-t)\right)\nonumber\\
&=&P_{\xi^{(1)}}\left(-\sum_{\mu=2}^{N_2} F_p(X_J^{(\mu)}, Y_J^{(\mu)})\geq N_1^{\frac {p_+}{2}} y\right)\nonumber\\
&=&P\left(-\sum_{\mu=2}^{N_2} F_p(X_J^{(\mu)}, Y_J^{(\mu)})\geq N_1^{\frac {p_+}{2}} y\right)\,,\label{eq;cont}
\end{eqnarray}
because of independence of the patterns and $F_p(r,1-r)\geq0$.

Note that $y>0$ is equivalent to $0\leq t<\bar F_p(r)$.
By Lemma \ref{lemma:XYsubG} $\{F_p(X_J^{(\mu)}, Y_J^{(\mu)})\}_{\mu\in[N_2]}$
are centred i.i.d. r.vs which fit the assumptions of Proposition \ref{prop:Hoeffging-GEN} below (with $\ell=2/p\in(0,1]$).
Therefore
\begin{equation}\label{eq:p>2}
\eqref{eq;cont}\leq \exp\left(-c\min\left(\frac{N^{p}_1y^2}{N_2 (r(1-r))^\frac p2},\frac{N_1|y|^{\frac2p} }{\sqrt{r(1-r)}}\right)\right)\,.
\end{equation}
The value of this minimum depends on $\alpha$. Bearing in mind that $y<\bar F_p(r)\lesssim r^\frac p2$ we have that the sub-linear term prevails if
$$
\alpha\lesssim N_1^{\frac{p-2}{2}}r^{\frac{p-1}{2}}\,,
$$
otherwise we have a Gaussian tail.
\end{proof}

Now we are ready for the main proof.

\begin{proof}[Proof of Theorem \ref{TH1}]
We shorten
\begin{equation}\label{eq:defD}
\mathtt D_{\mu,N_1}(r_0):=\{\sigma\in \widehat B^{N_1}_{\mu,\lfloor r_0N_1 \rfloor}\,:\, H^{(p)}(\xi^{(\mu)};\xi)\geq H^{(p)}(\sigma;\xi)\}
\end{equation}
and note that since for any $\mu\in[N_2]$
$$
\mathtt{dLM}_{N_1}^{(\mu)}\cap \widehat B^{N_1}_{\mu,\lfloor r_0N_1 \rfloor}\subseteq \mathtt D_{\mu,N_1}(r_0)\,,
$$
it is
\begin{equation}\label{eq:1disth1}
P\left(\forall \mu\in[N_2]\,\,\,\,\mathtt{dLM}_{N_1}^{(\mu)}\cap \widehat B^{N_1}_{\mu,\lfloor r_0N_1 \rfloor}\subseteq \widehat B^{N_1}_{\mu,R}\right)\geq
P\left(\forall \mu\in[N_2]\,\,\,\,\mathtt D_{\mu,N_1}(r_0)\subseteq \widehat B^{N_1}_{\mu,R}\right)\,.
\end{equation}
Let us introduce the sets
\begin{equation}
\Barr_{N_1,N_2,p}(n):=\left\{\min_{\mu\in[N_2]}\min_{\sigma\in\widehat S^{N_1-1}_{\mu,n}}H^{(p)}(\sigma)-H^{(p)}(\xi^{(\mu)})\geq t(n) \right\}\,
\end{equation}
on which the minimal energy gap of $n$ flips from the patterns is a given $t(n)$. Write now for $r\in(0,\frac12)$ $n=\lfloor rN_1\rfloor$ and $t(n)=t(r)=t$ and choose, according to Lemma \ref{lemma:crux1}, $t\simeq r^{\frac p2}$. Also we take a $r_0\in(0,\frac12]$ such that $S(r)\leq c\sqrt{r}$, for a constant $c>0$ to be specified below, for all $r\in[0,r_0]$.
Then, bearing in mind (\ref{eq:1disth1}), the crux is
\begin{eqnarray}
&&P\left(\forall \mu\in[N_2]\,\,\,\,\mathtt D_{\mu,N_1}(r_0)\subseteq \widehat B^{N_1}_{\mu,R}\right)
\geq P\left(\bigcap_{n=\lfloor R\rfloor}^{\lfloor r_0N_1 \rfloor}\Barr_{N_1,N_2,p}(n)\right)\nonumber\\
&\geq&1-\sum_{n= \lfloor R\rfloor }^{\lfloor r_0N_1 \rfloor} P(\Barr^{c}_{N_1,N_2,p}(n))
\gtrsim 1-N_1\min_{\lfloor R\rfloor\leq n\leq \lfloor r_0N_1 \rfloor} P(\Barr^{c}_{N_1,N_2,p}(n))\,. \label{eq:very-impo}
\end{eqnarray}
By the standard estimate (\ref{eq:standard-bound}) and the union bound we have
\begin{eqnarray}
&&P\left(\Barr^c_{N_1,N_2,p}(n)\right)=P\left(\min_{\mu\in[N_2]}\min_{\sigma\in \widehat S^{N_1-1}_{\mu,n}}H^{(p)}(\sigma)-H^{(p)}(\xi_\mu)\leq t \right)\nonumber\\
&\leq& \exp(\log N_2+N_1
S(r))E\left[\sup_{\sigma\in \widehat S^{N_1-1}_{1,\lfloor rN_1\rfloor}}P_{\xi^{(1)}}\left(H^{(p)}(\xi^{(1)})-H^{(p)}(\sigma)\geq-t \right)\right]\,.\nonumber
\end{eqnarray}

The probabilities appearing in the last line are evaluated using Lemma \ref{lemma:crux1} with $y\simeq t$.

Let us first deal with $p>2$. We consider different regimes. If $r\gtrsim N^{-x}$ (\textit{i.e. } $n\gtrsim N_1^{1-x}$) with $x>x_p:=\frac{p-2}{p-1}$ then
Lemma \ref{lemma:crux1} yields for all $\alpha>0$
\begin{equation}\label{eq:combinalo-o-o}
P\left(\Barr^c_{N_1,N_2,p}(n)\right)\leq N_2\exp\left(N_1 \left
S(r)-c\sqrt r\right)\right)\lesssim e^{-c\sqrt{r}N_1}\simeq e^{-c\sqrt{nN_1}}\,.
\end{equation}
The above estimate must hold for all $r\in[0,r_0]$, which determines $r_0$. Thus
\begin{equation}\label{eq:pallaextp>2}
\min_{\lfloor N_1^{1-x_p}\rfloor < n\leq \lfloor r_0N_1 \rfloor} P(\Barr^{c}_{N_1,N_2,p}(n))\leq e^{-cN^{1-\frac{x_p}{2}}_1}\,.
\end{equation}
For $r< N^{-x}$ with $x\in(x_p,1)$ Lemma \ref{lemma:crux1} gives for all $\alpha>0$
\begin{equation}\label{eq:u-u-u-usala}
P\left(\Barr^c_{N_1,N_2,p}(n)\right)\leq N_2e^{N_1
S(r)-\frac{N_1^{\frac p2}r^{\frac p2}}{\alpha}}\leq N_2e^{n|\log N_1|-\frac{n^{\frac p2}}{\alpha}}\,.
\end{equation}
Thus for all $\varepsilon>0$ sufficiently small
\begin{equation}\label{eq:pallaextp>2-2}
\min_{\lfloor N_1^{2\varepsilon}\rfloor<n<\lfloor N_1^{1-x_p}\rfloor} P(\Barr^{c}_{N_1,N_2,p}(n))\leq e^{-cN_1^{\varepsilon p}}\,.
\end{equation}
Moreover by \eqref{eq:u-u-u-usala} we see that also a poly-log number of flips is allowed:
\begin{equation}\label{eq:pallaextp>2-3}
\min_{\lfloor (\log N_1)^{\frac2{p-2}}\rfloor<n\leq\lfloor N_1^{2\varepsilon}\rfloor} P(\Barr^{c}_{N_1,N_2,p}(n))\leq e^{-c(\log N_1)^{1+\frac{2}{2-p}}}\,.
\end{equation}
Finally we look at $n\simeq N_1^{1-x_p}$. In this case we have to fix some $\alpha_0>0$ and we use for $\alpha\leq \alpha_0$ the bound (\ref{eq:combinalo-o-o}) and for $\alpha>\alpha_0$ the bound (\ref{eq:u-u-u-usala}). We have
\begin{equation}\label{eq:pallaextp>2-fin}
\min_{n\simeq\lfloor N_1^{1-x_p}\rfloor} P(\Barr^{c}_{N_1,N_2,p}(n))\lesssim e^{-cN^{1-\frac{x_p}{2}}_1}\,.
\end{equation}

Combining (\ref{eq:very-impo}) with $R=(\log N_1)^{\frac2{p-2}}$ and (\ref{eq:pallaextp>2}), \eqref{eq:pallaextp>2-2}, \eqref{eq:pallaextp>2-3}, \eqref{eq:pallaextp>2-fin} we get
\begin{equation}\label{eq:final-retr-p>2}
P\left(\forall \mu\in[N_2]\,\,\,\,\mathtt{dLM}_{N_1}^{(\mu)}\cap \widehat B^{N_1}_{\mu,\lfloor r_0N_1 \rfloor}\subseteq \widehat B^{N_1}_{\mu,\lfloor (\log N_1)^{\frac2{p-2}}\rfloor}\right)\gtrsim 1-N_1e^{-c(\log N_1)^{1+\frac{2}{2-p}}}
\end{equation}
whence (\ref{eq:TH1-p>2}) follows.

Now we look at $p=2$. If $\alpha\lesssim\sqrt r$ then Lemma \ref{lemma:crux1} gives
\begin{equation}\label{eq:bound0-retr-p=2}
P\left(B^c_{N_1,N_2,p}(n)\right)\leq N_2\exp\left(N_1 \left
S(r)-c\sqrt r\right)\right)\lesssim e^{-c\sqrt{nN_1}}\qquad r\in[0,r_0]\,,
\end{equation}
where again the above estimate determines $r_0$.
Thus by (\ref{eq:very-impo}) with $R=\lfloor rN_1 \rfloor$
\begin{equation}\label{eq:caso.p=2}
P\left(\forall \mu\in[N_2]\,\,\,\,\mathtt{dLM}_{N_1}^{(\mu)}\cap \widehat B^{N_1}_{\mu,\lfloor r_0N_1 \rfloor}\subseteq \widehat B^{N_1}_{\mu,\lfloor rN_1 \rfloor}\right)\gtrsim 1-N_1e^{-c\sqrt{r}N_1}\,,
\end{equation}
whence the $p=2$ part of Theorem \ref{TH1} follows.
\end{proof}

\section{Absence of retrieval for $p\in(1,2]$}\label{sect:absence}

In this section we present the proof of Theorem \ref{TH2}.

We set for brevity for $\mu\in[N_2]\setminus \{1\}$, $p\in[1,2]$, $k\in[N_1]$, $J\subseteq[N_1]$, $v\in\{-1,1\}^{N_1}$
\begin{equation}\label{eq:W}
W_{p,k,J}^{(\mu)}(v):=\frac{2p}{\sqrt{N_1}}\xi_k^{(\mu)}v_k\operatorname{sign}(Z_{J}^{(\mu)}(v))|Z_{J}^{(\mu)}(v)|^{p-1}
\end{equation}
where
\begin{equation}\label{eq:defZ}
Z_J^{(\mu)}(v):=\frac{1}{\sqrt{N_1}}(\xi^{(\mu)},F_Jv)\,.
\end{equation}

Next we give the central technical lemma employed in the proof of Theorem \ref{TH2}. In the sequel we shorten $J+k:=J\cup\{k\}$ if $k\notin J$ and $J-k:=J\setminus\{k\}$ if $k\in J$.
\begin{lemma}\label{lemma:reprH-absence}
Let $r\in(0,\frac12]$, $J\subseteq[N_1]$ with $|J|=\lfloor rN_1\rfloor$.
For any $p\in[1,2]$ we have
\begin{eqnarray}
N_1(H^{(p)}(F_{J_{\pm k}}\xi^{(1)})-H^{(p)}(F_J\xi^{(1)}))&=&\mp \sum_{\mu=2}^{N_2} W_{p,k,J}^{(\mu)}(\xi^{(1)})-\alpha\varsigma N_1^{1-\frac p2} \nonumber\\
&\pm&\frac{2p(1-2r)^{p-1}}{N_1^{1-\frac p2}}+\OO{N_1^{2-\frac p2}}\label{eq:reprH-absence}
\end{eqnarray}
where $\varsigma$ is a positive and uniformly bounded random variable depending on $\{\xi^{(1)}_k\xi^{(\mu)}_k, Z_{J}^{(\mu)}(\xi^{(1)})\}_{\mu=2...N_2}$.
In particular
\begin{equation}\label{eq:reprH-absence-p=2}
\frac{N_1}{4}(H^{(2)}(F_{J_{\pm k}}\xi^{(1)})-H^{(2)}(F_J\xi^{(1)}))=\mp\sum_{\mu\geq2}\frac{\xi_k^{(1)}\xi_k^{(\mu)}}{\sqrt{N_1}}Z_J^{(\mu)}(\xi^{(1)})-\alpha\pm(1-2r)\mp\frac{1}{N_1}\,.
\end{equation}
\end{lemma}
\begin{proof}
By Lemma \ref{lemma:repr-inutile} we have for any $k\notin J$
\begin{eqnarray}
N_1(H^{(p)}(F_{J_{+k}}\xi^{(1)})-H^{(p)}(F_J\xi^{(1)}))&=&N_1(H(F_{J_{+ k}}\xi^{(1)})-H(\xi^{(1)})-(H(F_J\xi^{(1)})-H(\xi^{(1)})))\nonumber\\
&=&N_1^{\frac p2}\left(\bar F_p\left(r+\frac1{N_1}\right)-\bar F_p(r)\right)\nonumber\\
&+&\sum_{\mu=2}^{N_2} F_p(X_{J_{+k}}^{(\mu)}(\xi^{(1)}), Y_{J_{+k}}^{(\mu)}(\xi^{(1)}))\nonumber\\
&-&\sum_{\mu=2}^{N_2} F_p(X_{J}^{(\mu)}(\xi^{(1)}), Y_{J}^{(\mu)}(\xi^{(1)}))\nonumber\\
&=&N_1^{\frac p2}\left(|1-2r|^p-\left|1-2r-\frac2{N_1}\right|^p\right)\nonumber\\
&+&\sum_{\mu=2}^{N_2}\left( |Z_J^{(\mu)}(\xi^{(1)})|^p-\left|Z_{J}^{(\mu)}(\xi^{(1)})+2\frac{\xi_k^{(\mu)}\xi_k^{(1)}}{\sqrt{N_1}}\right|^p \right)\,.\label{eq:repr-H+}
\end{eqnarray}
Similarly for all $k\in J$
\begin{eqnarray}
N_1(H^{(p)}(F_{J_{-k}}\xi^{(1)})-H^{(p)}(F_J\xi^{(1)}))&=&N_1^{\frac p2}\left(|1-2r|^p-\left|1-2r+\frac2{N_1}\right|^p\right)\nonumber\\
\!\!\!\!\!&+&\!\!\!\!\!\sum_{\mu=2}^{N_2}\left( |Z_J^{(\mu)}(\xi^{(1)})|^p-\left|Z_{J}^{(\mu)}(\xi^{(1)})-2\frac{\xi_k^{(\mu)}\xi_k^{(1)}}{\sqrt{N_1}}\right|^p\right)\label{eq:repr-H-}
\end{eqnarray}
For $p=2$ a straightforward computation gives (\ref{eq:reprH-absence-p=2}) from (\ref{eq:repr-H+}) and (\ref{eq:repr-H-}).

In general for $p\in[1,2)$ we have to use Taylor expansion. Let $(p)_0:=1$ and $(p)_k:=\prod_{j=0}^{k-1}(p-j)$ for $k\geq1$. Assuming $r\in(0,\frac12)$, $N_1$ large enough (\textit{i.e. } $N_1(1-2r)>2$) we have
\begin{equation}
\left|1-2r\pm\frac2{N_1}\right|^p=|1-2r|^p\pm\frac{2p}{N_1}|1-2r|^{p-1}+\frac{4}{N_1^2}\sum_{k\geq2}\frac{(p)_k}{k!}\frac{2^{k-2}(1-2r)^{p-k}}{N_1^{k-2}}
\end{equation}
and using $|(p)_k|\leq k!$ we get
\begin{equation}
\left|\frac{4}{N_1^2}\sum_{k\geq2}\frac{(p)_k}{k!}\frac{2^{k-2}(1-2r)^{p-k}}{N_1^{k-2}}\right| \leq \frac{4}{N_1^2(1-2r)^{2-p}}\sum_{k\geq0}\frac{2^k}{((1-2r)N_1)^{k}}\lesssim\frac{1}{N_1^2(1-2r)^{2-p}}\,.
\end{equation}
Therefore
\begin{equation}\label{eq:bound-pezzo1}
N_1^{\frac p2}\left(|1-2r|^p-\left|1-2r\pm\frac2{N_1}\right|^p\right)=\mp\frac{2p}{N_1^{1-\frac p2}}|1-2r|^{p-1}+ \OO{N_1^{2-\frac p2}}\,.
\end{equation}
On the other hand, for $r=\frac12$ this correction term is trivially of order $N^{-\frac p2}$.

Recall now \eqref{eq:W} and compute
\begin{eqnarray}
&&\left|Z_{J}^{(\mu)}(\xi^{(1)})\pm2\frac{\xi_k^{(\mu)}\xi_k^{(1)}}{\sqrt{N_1}}\right|^p\nonumber\\
&=&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|\geq 2N_1^{-\frac12}\}}\Big( |Z_{J}^{(\mu)}(\xi^{(1)})|^p\pm W_{p,k,J}^{(\mu)}(\xi^{(1)}) \nonumber\\
&+& \sum_{\ell\geq 2}\frac{(p)_\ell}{\ell!}\frac{|Z_{J}^{(\mu)}(\xi^{(1)})|^p}{(Z_{J}^{(\mu)}(\xi^{(1)}))^\ell} \frac{(\pm 2\xi^{(\mu)}_k\xi^{(1)}_k)^\ell}{N_1^{\frac \ell2}} \Big)\nonumber\\
&+&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|< 2N_1^{-\frac12}\}}\frac{2^p}{N_1^{\frac p2}}\Big(1\pm \frac p2\xi_k^{(\mu)}\xi_k^{(1)}Z_{J}^{(\mu)}(\xi^{(1)})\sqrt{N_1}\nonumber\\
&+&\sum_{\ell\geq 2}\frac{(p)_\ell}{\ell!} \left(\pm\frac{\xi_k^{(\mu)}\xi_k^{(1)}}{2}Z_{J}^{(\mu)}(\xi^{(1)})\sqrt{N_1}\right)^\ell\Big)\nonumber\,.
\end{eqnarray}
Thus
\begin{eqnarray}
&&\left( |Z_J^{(\mu)}(\xi^{(1)})|^p-\left|Z_{J}^{(\mu)}(\xi^{(1)})\pm2\frac{\xi_k^{(\mu)}\xi_k^{(1)}}{\sqrt{N_1}}\right|^p\right)=\mp W_{p,k,J}^{(\mu)}(\xi^{(1)})\nonumber\\
&-&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|\geq 2N_1^{-\frac12}\}}\sum_{\ell\geq 2}\frac{(p)_\ell}{\ell!}\frac{|Z_{J}^{(\mu)}(\xi^{(1)})|^p}{(Z_{J}^{(\mu)}(\xi^{(1)}))^\ell} \frac{(\pm2\xi^{(\mu)}_k\xi^{(1)}_k)^\ell}{N_1^{\frac \ell2}}\label{eq:riga2}\\
&-&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|< 2N_1^{-\frac12}\}}\left(\frac{2^p}{N_1^p}\pm\frac{2^p}{N_1^p}\frac{\xi^{(\mu)}_k\xi^{(1)}_k}{2}Z_{J}^{(\mu)}(\xi^{(1)})\sqrt{N_1}- |Z_{J}^{(\mu)}(\xi^{(1)})|^p\mp W_{p,k,J}^{(\mu)}(\xi^{(1)})\right)\nonumber\\\label{eq:riga1}\\
&-&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|< 2N_1^{-\frac12}\}}\frac{2^p}{N_1^{\frac p2}}\sum_{\ell\geq 2}\frac{(p)_\ell}{\ell!} 2^{-\ell}(\pm Z_{J}^{(\mu)}(\xi^{(1)})\xi_k^{(\mu)}\xi_k^{(1)}\sqrt{N_1})^\ell\,.\label{eq:riga3}
\end{eqnarray}

The contributions (\ref{eq:riga2}) and (\ref{eq:riga3}) are very similar and will be dealt together.
Using that for $\ell\geq2$ it is $\ell-p>0$ and $\ell^2(p)_\ell\leq 2p^2\ell!$ we have
\begin{equation}\label{eq:bound-riga2}
|\eqref{eq:riga2}|\leq \frac{2p^2}{N_1^{\frac p2}}\sum_{\ell\geq2} \ell^{-2}\,.
\end{equation}
Similarly
\begin{equation}\label{eq:bound-riga3}
|\eqref{eq:riga3}|\leq \frac{2^{p+1}p^2}{N_1^{\frac p2}}\sum_{\ell\geq2} \ell^{-2}\,.
\end{equation}

Furthermore, depending on the value of $\operatorname{sign}(Z_{J}^{(\mu)}(\xi^{(1)}))\xi^{(1)}_k\xi^{(\mu)}_k$ we can have either
\begin{eqnarray}
&&\eqref{eq:riga2}=1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|< 2N_1^{-\frac12}\}}\frac{1}{N_1^{\frac p2}}\sum_{\ell\geq 2}\frac{(p)_\ell}{\ell!}\frac{|Z_{J}^{(\mu)}(\xi^{(1)})|^p}{|Z_{J}^{(\mu)}(\xi^{(1)})|^\ell} \frac{1}{N_1^{\frac \ell2}}\nonumber\\
&=&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|< 2N_1^{-\frac12}\}}\frac{|Z_{J}^{(\mu)}(\xi^{(1)})|^p}{N_1^{\frac p2}}\left[(1+|Z_{J}^{(\mu)}(\xi^{(1)})\sqrt N_1|^{-1})^p-(1+|Z_{J}^{(\mu)}(\xi^{(1)})\sqrt N_1|^{-1})\right]\geq0\nonumber\\\,
\end{eqnarray}
and
\begin{eqnarray}
&&\eqref{eq:riga3}=1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|\geq 2N_1^{-\frac12}\}}\frac{2^p}{N_1^{\frac p2}}\sum_{\ell\geq 2}\sum_{\ell\geq 2}\frac{(p)_\ell}{\ell!} 2^{-\ell}| Z_{J}^{(\mu)}(\xi^{(1)})\sqrt{N_1}|^\ell\nonumber\\
&=&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|\geq 2N_1^{-\frac12}\}}\frac{|Z_{J}^{(\mu)}(\xi^{(1)})|^p}{N_1^{\frac p2}}\left[(1+|2^{-1}Z_{J}^{(\mu)}(\xi^{(1)})\sqrt N_1|)^p-(1+|2^{-1}Z_{J}^{(\mu)}(\xi^{(1)})\sqrt N_1|)\right]\geq0\nonumber\\\,
\end{eqnarray}
(equality is achieved only if $p=1$) or
\begin{eqnarray}
&&\eqref{eq:riga2}=1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|< 2N_1^{-\frac12}\}}\frac{1}{N_1^{\frac p2}}\sum_{\ell\geq 2}(-1)^\ell\frac{(p)_\ell}{\ell!}\frac{|Z_{J}^{(\mu)}(\xi^{(1)})|^p}{|Z_{J}^{(\mu)}(\xi^{(1)})|^\ell} \frac{2^\ell}{N_1^{\frac \ell2}}\nonumber\\
&=&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|< 2N_1^{-\frac12}\}}\frac{|Z_{J}^{(\mu)}(\xi^{(1)})|^p}{N_1^{\frac p2}}\sum_{\ell\geq 1}\left(\frac{\sqrt{N_1}Z_{J}^{(\mu)}(\xi^{(1)})}{2}\right)^{-2\ell}\left(\frac{(p)_{2\ell}}{{2\ell}!}-\frac{(p)_{2\ell+1}}{{2\ell+1}!}\left|\frac{\sqrt{N_1}Z_{J}^{(\mu)}(\xi^{(1)})}{2}\right|^{-1}\right)\nonumber\\
\end{eqnarray}
(split the sum over even and odd $\ell\geq2$ and rename the indices) and
\begin{eqnarray}
&&\eqref{eq:riga3}=1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|\geq 2N_1^{-\frac12}\}}\frac{2^p}{N_1^{\frac p2}}\sum_{\ell\geq 2}(-1)^\ell\frac{(p)_\ell}{\ell!}2^{-\ell}| Z_{J}^{(\mu)}(\xi^{(1)})\sqrt{N_1}|^\ell\nonumber\\
&=&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|\geq 2N_1^{-\frac12}\}}\frac{2^p}{N_1^{\frac p2}}\sum_{\ell\geq 1}\left(\frac{\sqrt{N_1}Z_{J}^{(\mu)}(\xi^{(1)})}{2}\right)^{2\ell}\left(\frac{(p)_{2\ell}}{{2\ell}!}-\frac{(p)_{2\ell+1}}{{2\ell+1}!}\frac{|\sqrt{N_1}Z_{J}^{(\mu)}(\xi^{(1)})|}{2}\right)\,.\nonumber\\
\end{eqnarray}

The quantities in the two displays above are positive, since  for $\ell\geq1$ $(p)_{2\ell}\geq0$ and $(p)_{2\ell+1}\leq0$, which can be shown by observing that $(p)_{\ell\geq3}=(-1)^{\ell}p(p-1)\prod_{j=2}^{\ell-1}(j-p)$ (again the equality is achieved only for $p=1$). We conclude
$$
\eqref{eq:riga2}<0,\quad \eqref{eq:riga3}<0\quad\mbox{for $p\in(1,2)$ and}\quad \eqref{eq:riga2}=\eqref{eq:riga3}=0\quad\mbox{for $p=1$.}
$$

With a bit of algebra we rewrite the term in (\ref{eq:riga1}) as
\begin{eqnarray}
&&1_{\{|Z_{J}^{(\mu)}(\xi^{(1)})|< 2N_1^{-\frac12}\}}\frac{2^p}{N_1^\frac{p}{2}}\left(\left(1-\frac{|2Z_{J}^{(\mu)}(\xi^{(1)})|^p}{N_1^{\frac p2}}\right)\right.\nonumber\\
&\mp&\left. p\xi^{(\mu)}_k\xi^{(1)}_k\operatorname{sign}(Z_{J}^{(\mu)}(\xi^{(1)}))\left(\frac{|2Z_{J}^{(\mu)}(\xi^{(1)})|^{p-1}}{N_1^{\frac {p-1}2}}-\frac{|2Z_{J}^{(\mu)}(\xi^{(1)})|}{N_1^{\frac 12}}\right)\right)\label{eq:bitofalgebra}
\end{eqnarray}
According to the value of $\xi^{(\mu)}_k\xi^{(1)}_k\operatorname{sign}(Z_{J}^{(\mu)}(\xi^{(1)}))$ the term inside the parenthesis can be either
$$
1-|x|^p+p(|x|^{p-1}-|x|)\quad\mbox{ or }\quad 1-|x|^p-p(|x|^{p-1}-|x|)\,,
$$
where we shortened $|x|:=\frac{|2Z_{J}^{(\mu)}(\xi^{(1)})|}{N_1^{\frac 12}}<1$. The first expression above is clearly positive, while the second one is positive thanks to Lemma \ref{lemma:tecnico}. From the representation (\ref{eq:bitofalgebra}) we also get the bound
\begin{equation}\label{eq:bound-algebra}
\eqref{eq:riga1}\leq \frac{C}{N_1^\frac p2}\,.
\end{equation}

Therefore there is a positive and uniformly bounded random variable $$\varsigma:=\varsigma(\{\xi^{(1)}_k\xi^{(\mu)}_k, Z_{J}^{(\mu)}(\xi^{(1)})\}_{\mu=2...N_2})$$ such that
\begin{equation}\label{eq:bound-pezzo2}
\sum_{\mu=2}^{N_2}\left( |Z_J^{(\mu)}(\xi^{(1)})|^p-\left|Z_{J}^{(\mu)}(\xi^{(1)})\pm2\frac{\xi_k^{(\mu)}\xi_k^{(1)}}{\sqrt{N_1}}\right|^p\right)=\mp \sum_{\mu=2}^{N_2} W_{p,k,J}^{(\mu)}(\xi^{(1)})+ \alpha\varsigma N_1^{1-\frac p2}\,.
\end{equation}
This and (\ref{eq:bound-pezzo1}) give (\ref{eq:reprH-absence}).
\end{proof}

Now we turn to the proof of Theorem \ref{TH2}, which we conveniently split in severals steps.

{\em Step 1: reduction.}
Due to the exchangeability of the patterns and their entries we have
\begin{eqnarray}
&&P\left(\mathtt{LM}^{(\mu)}_{N_1}\cap \widehat B^{N_1}_{\mu, \lfloor r N_1\rfloor}\neq \emptyset\right)\nonumber\\
\!\!\!\!&\leq& \!\!\!\!N_2 \sum_{\ell=1}^{\lfloor rN_1\rfloor} \sum_{\substack{J\subset [N_1]\\|J|=\ell}} P\left(\bigcap_{k\notin J}\{H(F_{J_{+k}}\xi^{(1)})-H(F_{J}\xi^{(1)})>0\}\,,\,\,\bigcap_{k\in J}\{H(F_{J-k}\xi^{(1)})-H(F_{J}\xi^{(1)})>0\}\right)\,\nonumber\\
\!\!\!\!&=&\!\!\!\!\!\!N_2\sum_{\ell=1}^{\lfloor rN_1\rfloor}\binom{N_1}{\ell}P\left(\bigcap_{k>\ell}\{H(F_{[\ell]_{+k}}\xi^{(1)})-H(F_{[\ell]}\xi^{(1)})>0\}\,,\,\,\bigcap_{k\leq\ell}\{H(F_{[\ell]-k}\xi^{(1)})-H(F_{[\ell]}\xi^{(1)})>0\}\right)\nonumber\,.\\\label{eq:proof1}
\end{eqnarray}
We set for brevity
\begin{eqnarray}
Z_{\ell}^{(\mu)}(\xi^{(1)})&:=&\frac{1}{\sqrt{N_1}}\left(\sum_{k=1}^{\ell}\xi_k^{(1)}\xi_k^{(\mu)}-\sum_{k\geq\ell+1}\xi_k^{(1)}\xi_k^{(\mu)}\right)\,,\\
W_{p,k,\ell}^{(\mu)}(\xi^{(1)})&:=&\frac{2p}{\sqrt{N_1}}\xi_k^{(\mu)}\xi_k^{(1)}\operatorname{sign}(Z_{[\ell]}^{(\mu)}(\xi^{(1)}))|Z_{[\ell]}^{(\mu)}(\xi^{(1)})|^{p-1}\,,\quad W_{p,k,\ell}^{(\mu)}:=W_{p,k,\ell}^{(\mu)}(1)\,,\\
M^{(\mu)}&:=&\frac1{N_1}\sum_{i\in[N_1]}\xi_i^{(\mu)}\,,\quad Q^{(\mu)}:=\operatorname{sign}(M^{(\mu)})|M^{(\mu)}|^{p-1}\,.
\end{eqnarray}
By Lemma \ref{lemma:tecnico} in Appendix \ref{Sect:AppTH2} we have for $N_1$ large enough ($r':=\ell/N_1$)
\begin{eqnarray}
\eqref{eq:proof1}&=&N_2\sum_{\ell=1}^{\lfloor rN_1\rfloor}\binom{N_1}{\ell}P\left(\forall k>\ell\,\,\,-\sum_{\mu=2}^{N_2} W_{p,k,\ell}^{(\mu)}(\xi^{(1)})\gtrsim\alpha N_1^{1-\frac p2} -\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}}\,,\,\right.\nonumber\\
&&\left.\forall k\in[\ell] \,\,\,\sum_{\mu=2}^{N_2} W_{p,k,\ell}^{(\mu)}(\xi^{(1)})\gtrsim\alpha N_1^{1-\frac p2}+\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}} \right)\nonumber\\
&=&N_2\sum_{\ell=1}^{\lfloor rN_1\rfloor}\binom{N_1}{\ell}P\left(\forall k>\ell\,\,\,-\sum_{\mu=2}^{N_2} W_{p,k,\ell}^{(\mu)}\gtrsim\alpha N_1^{1-\frac p2} -\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}}\,,\,\right.\nonumber\\
&&\left.\forall k\in[\ell] \,\,\,\sum_{\mu=2}^{N_2} W_{p,k,\ell}^{(\mu)}\gtrsim\alpha N_1^{1-\frac p2}+\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}} \right)\nonumber\\
&=&N_2\sum_{\ell=1}^{\lfloor rN_1\rfloor}\binom{N_1}{\ell}P\left(\forall k>\ell\,\,\,(Q,\tilde\xi^{(k)})\gtrsim\alpha N_1^{1-\frac p2} -\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}}\,,\,\right.\nonumber\\
&&\left.\forall k\in[\ell] \,\,\,(Q,\tilde\xi^{(k)})\gtrsim\alpha N_1^{1-\frac p2}+\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}} \right)
\end{eqnarray}
(recall that $\tilde\xi^{(k)}$ denotes the $k$-th transposed pattern). In the second identity above we have exploited independence of the pattern $\xi^{(1)}$ to replace $W_{p,k,\ell}^{(\mu)}(\xi^{(1)})$ by $W_{p,k,\ell}^{(\mu)}$ and in the third one the independence of the first $\ell$ entries from all the others and the flip-symmetry to replace $Z_{\ell}^{(\mu)}$ by $M^{(\mu)}$. We also used that the independent random variables $Q^{(2)}\,\ldots Q^{(N_2)}$ are symmetric.

{\em Step 2: disentangling by the FKG inequality.}
Let us now call $\mathcal M$ the $\sigma$-field generated by $M^{(2)},\ldots, M^{(N_2)}$ and notice that $Q^{(2)}\,\ldots Q^{(N_2)}$ are $\mathcal M$-measurable. We shorten $P_{\mathcal M}(\cdot):=P(\cdot\,|\,\mathcal M)$.
The crucial observation here (first done for the Hopfield model in \cite{lou}) is that for each $\mu\in[N_2]$ the law of $\xi^{(\mu)}$ conditionally on  $\mathcal M$ is a permutation distribution (as the increments of a simple random walk given the position).
Therefore the FKG inequality applies (see for instance \cite{brics}) and we have
\begin{eqnarray}
&&\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!\!P_{\mathcal M}\left(\forall k>\ell\,\,\,(Q,\tilde\xi^{(k)})\gtrsim\alpha N_1^{1-\frac p2} -\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}}\right.\,,\nonumber\\
&&\left.\,\forall k\in[\ell] \,\,\,(Q,\tilde\xi^{(k)})\gtrsim\alpha N_1^{1-\frac p2}+\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}} \right)\nonumber\\
&\leq&\prod_{k>\ell}P_{\mathcal M}\left((Q,\tilde\xi^{(k)})\gtrsim\alpha N_1^{1-\frac p2} -\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}}\right)\nonumber\\
&&\prod_{k\in[\ell]}P_{\mathcal M}\left((Q,\tilde\xi^{(k)})\gtrsim\alpha N_1^{1-\frac p2} +\frac{2p(1-2r')^{p-1}}{N_1^{1-\frac p2}}\right)\label{eq:continua}
\end{eqnarray}
Note that $E[\xi_k^{(\mu)}\,|\,\mathcal M]=M^{(\mu)}$ and $\Var[\xi_k^{(\mu)}\,|\,\mathcal M]=1-(M^{(\mu)})^2$.
Therefore introducing $\overline Q\in\mathbb{R}^{N_2}$ with components $\overline Q^{(\mu)}:=Q^{(\mu)}\sqrt{1-(M^{(\mu)})^2}$, we get
\begin{equation}\label{eq:Hoeff}
P_{\mathcal M}\left((Q,\tilde\xi^{(k)})\geq t\right)\leq \exp\left( -c\frac{(t-\|M\|_p^p)^2}{\|\overline Q\|_2^2} \right)\,,
\end{equation}
by the Hoeffding inequality. Note that this quantity is independent on $k$ and also we have the simple bound $\|\overline Q\|_2^2\leq \| Q\|_2^2=\|M\|^{2p-2}_{2p-2}$. Then
\begin{eqnarray}
&&\eqref{eq:continua}_{p\in[1,2)}\leq \exp\left( -cN_1\frac{(\alpha N_1^{1-\frac p2}-\|M\|_p^p)^2}{\|M\|^{2p-2}_{2p-2}} \right)\,,\label{eq:rhs1}\\
&&\eqref{eq:continua}_{p=2}\leq \exp\left( -\frac{cN_1}{\|M\|^{2}_{2}}\left(r'(\alpha-(1-2r')-\|M\|_2^2)^2+(1-r')(\alpha+(1-2r')-\|M\|_2^2)^2\right) \right)\,.\nonumber\\\label{eq:rhs2}
\end{eqnarray}
{\em Step 3: concentration.}
Now we have to take the global expectations of the r.h.s. above. A similar computation as in \eqref{eq:proceeding} gives $E[\|M\|_p^p]\simeq\alpha N_1^{1-\frac p2}$ (this is an identity for $p=2$). Moreover $\|(M^{(\mu)})^p-E[(M^{(\mu)})^p]\|_{\psi_{\frac2p}}\simeq N^{-\frac p2}$. Hence by Proposition \ref{prop:Hoeffging-GEN} in Appendix \ref{Sect:Tail} (with $\ell=2/p\geq1$)
\begin{equation}\label{eq:conc-Mp}
P\left(\left|\sum_{\mu\geq2} (M^{(\mu)})^p-E[(M^{(\mu)})^p]\right|\geq t \right)\leq \exp\left(-c\min\left(\frac{t^2N_1^{p-1}}{\alpha},\frac{t^{\frac2p}N^{2-\frac2p}}{\alpha^{\frac 2p-1}}\right)\right)\,.
\end{equation}
{\em Step 4: finalising the argument for $p\in(1,2)$.}
Using (\ref{eq:conc-Mp}) with $t= \frac12(\alpha N^{1-\frac p2}-E[\|M\|_p^p]$ in the first line below and likewise but for $2p-2$ in the second line we have
\begin{eqnarray}
E[\mbox{r.h.s. of \eqref{eq:rhs1}}]&\leq &E\left[\exp\left( -c\frac{\alpha^2N_1^{3-p}}{\|M\|^{2p-2}_{2p-2}} \right)\right]+e^{-c\min(\alpha,\alpha^2)N_1}\nonumber\\
&\leq &e^{-c\alpha N_1}+2e^{-c\min(\alpha,\alpha^2)N_1}\label{eq:final-bound1}
\end{eqnarray}
Combining (\ref{eq:proof1}), (\ref{eq:continua}), (\ref{eq:rhs1}), (\ref{eq:final-bound1}) and using the standard bound (\ref{eq:standard-bound}) we have for $p\in[1,2)$
\begin{eqnarray}
P\left(\mathtt{LM}^{(\mu)}_{N_1}\cap \widehat B^{N_1}_{\mu, \lfloor r N_1\rfloor}\neq\emptyset\right)&\leq& N_2\sum_{\ell=1}^{\lfloor rN_1\rfloor}\binom{N_1}{\ell}\left(e^{-c\alpha N_1}+2e^{-c\min(\alpha,\alpha^2)N_1}\right)\nonumber\\
&\lesssim& r N_2N_1 e^{N_1(S(r)-c\min(\alpha,\alpha^2))}\,,
\end{eqnarray}
which yields the assertion for $p\in[1,2)$.

{\em Step 5: finalising the argument for $p=2$.} Using (\ref{eq:conc-Mp})$_{p=2}$ with $t=\tau:=(1+\alpha)(1-2r')$ in the first line in the following display and $t=C\alpha$ in the second line, we have
\begin{eqnarray}
E[\mbox{r.h.s. of \eqref{eq:rhs2}}]&\leq &E\left[\exp\left( -cN_1\frac{((1-2r')-\tau)^2}{\|M\|^{2}_{2}} \right)\right]+e^{-cN_1(1+\alpha)(1-2r')\min\left(1,\frac{(1+\alpha)(1-2r')}{\alpha}\right)}\nonumber\\
&\leq &e^{-c\alpha N_1(1-2r')^2}+e^{-C\alpha N_1}+e^{-cN_1(1+\alpha)(1-2r')\min\left(1,\frac{(1+\alpha)(1-2r')}{\alpha}\right)}\nonumber\\
&\leq &e^{-c\alpha N_1(1-2r)^2}+e^{-C\alpha N_1}+e^{-cN_1(1+\alpha)(1-2r)\min\left(1,\frac{(1+\alpha)(1-2r)}{\alpha}\right)}\label{eq:final-bound2}\,,
\end{eqnarray}
as $r'$ ranges from $1/N_1$ to $r$.
So combining (\ref{eq:proof1}), (\ref{eq:continua}), (\ref{eq:rhs1}), (\ref{eq:final-bound2}) and using again (\ref{eq:standard-bound}) we obtain
\begin{eqnarray}
P\left(\mathtt{LM}^{(\mu)}_{N_1}\cap \widehat B^{N_1}_{\mu, \lfloor r N_1\rfloor}\neq\emptyset\right)&\leq& N_2\sum_{\ell=1}^{\lfloor rN_1\rfloor}\binom{N_1}{\ell}\left(e^{-c\alpha N_1(1-2r)^2}+e^{-c\alpha N_1}+e^{-cN_1(1+\alpha)(1-2r)\min\left(1,\frac{(1+\alpha)(1-2r)}{\alpha}\right)}\right)\nonumber\\
&\lesssim& r N_2N_1 \left(e^{N_1(S(r)-c\alpha(1-2r)^2)}+e^{N_1(S(r)-c\alpha)}\right.\nonumber\\
&+&\left.e^{N_1\left(S(r)-c(1+\alpha)(1-2r)\min\left(1,\frac{(1+\alpha)(1-2r)}{\alpha}\right)\right)}\right)\,.\nonumber
\end{eqnarray}
The first two summands are negative if $\alpha>S(r)/(1-2r)^2$. The fact that also the third one is so is verified in Lemma \ref{lemma:verify}, Appendix \ref{Sect:AppTH2}. The proof is complete.


\section{\bf Introduction}
\label{section:1}

The extremely rapid increase in viscosity or relaxation time when a liquid cools is a universal hallmark of glassy liquids \cite{RevModPhys.83.587,doi:10.1146/annurev-conmatphys-031113-133848,doi:10.1126/science.1120714,Karmakar_2015,acsomega.0c04831}. The temperature dependence of relaxation time, which is characterized by its ``fragility,'' \cite{Angell1924} is, however, non-universal, differing widely among glassy liquids. Glassy liquids which show an Arrhenius temperature dependence (exponential in inverse temperature) in relaxation time are called strong glassformers and those with stronger-than-Arrhenius (super-Arrhenius) dependence on temperature are known as fragile glassformers. An infamous property of the glass transition is that the increase of relaxation time upon cooling is unaccompanied by any qualitative change in the liquid structure, and very little quantitative change in it. As a result, structural insight into the differences between strong and fragile liquids has been lacking.

Many experimental \cite{Mauronatcomm,wei_linking_2015,doi:10.1063/1.4964362} and numerical \cite{sun_structural_2022,sastry_nature,Kaori_nature,PhysRevLett.126.208001} studies have been conducted to better understand the correlation of fragility to structure.
Recent efforts to characterize structure in terms of locally preferred structures \cite{doi:10.1063/1.2773716,Hallett,PhysRevX.8.011041}, medium-range crystalline order \cite{PhysRevLett.99.215701,tanaka_2010,Kawasaki_2011,PhysRevLett.121.085703}, and complex topological clusters~\cite{Mauro201131} have been successful for some systems but have not been applied across the strong-to-fragile spectrum of glassformers. However, a machine learning approach \cite{cubuk_identifying_2015,Schoenholz_Nat_Physics} shows promise. This approach identifies a microscopic structural variable called ``softness'' \cite{Schoenholz_Nat_Physics} that describes each particle's local structural environment in a way that is designed to be strongly correlated with its dynamics. It has been applied successfully across a spectrum of systems \cite{cubuk_identifying_2015,Schoenholz_Nat_Physics,Schoenholz263,doi:10.1126/science.aai8830,Sharp10943,ma_heterogeneous_2019,harrington_machine_2019,PhysRevE.101.010602,cubuk_unifying_2020,D0SM01575J,PhysRevX.11.041019,https://doi.org/10.48550/arxiv.2011.13049}.
Most important for our purpose here, softness appears to be associated with the typical energy barrier for rearrangement, with particles of higher softness having lower barriers. This association between structure and typical energy barrier holds across systems spanning a broad range of disorder spanning from glassforming liquids~\cite{Schoenholz_Nat_Physics} to crystals~\cite{Sharp10943} as well as a wide range of fragilities ranging from fragile (super-Arrhenius) behavior~\cite{Schoenholz_Nat_Physics} to strong (Arrhenius)~\cite{cubuk_unifying_2020} and even beyond (sub-Arrhenius behavior)~\cite{D0SM01575J}. Thus the softness of a particle appears to provide a useful structural measure of the typical energy barrier that must be surmounted for the particle to rearrange.

Here we extend the interpretation that softness provides a measure of such energy barriers to study the implied entropic barriers as well for a family of systems that span a wide range of fragilities. In order to study behavior across the family, we use the same definition of softness for all the systems.
We show that the average softness $\langle S \rangle$ and free energy barrier $\Delta F(\langle S \rangle)$ are stronger functions of temperature for fragile liquids than strong ones. Moreover, the deviation of the relaxation time from simple proportionality to $1/P_R(\langle S\rangle)$, where $P_R(S)$ is the probability that a particle of softness $S$ will rearrange, is compressed into a smaller range of temperatures relative to $T_g$ as fragility increases. Together, these three factors lead to a much stronger increase of the relaxation time with decreasing temperature for fragile liquids than for strong ones.

Sec.~\ref{section:2} describes the model and simulation details, Sec.~\ref{section:3} describes the machine learning analysis methodology, Sec.~\ref{section:4} contains our results and Section~\ref{section:5} summarizes our conclusions.

\section{\bf Model and Simulations}
\label{section:2}
We have conducted equilibrium molecular dynamics (MD) simulations in three dimensions for model glass forming liquids with harmonic repulsive interactions. Specifically, we study $50:50$ binary mixtures of $N = 10976$ particles that interact via the harmonic pair potential \cite{Durian_PRL1995,PhysRevLett.88.075507,doi:10.1146/annurev-conmatphys-070909-104045,Berthier_2009}:
\begin{equation}
V_{\alpha\beta}(r) = \epsilon_{\alpha\beta} \left ( 1-\frac{r}{\sigma_{\alpha\beta}}\right)^2
\end{equation}
for $r < \sigma_{\alpha\beta}$ and $V_{\alpha\beta}(r) = 0$ otherwise.
Here $(\sigma_{\alpha\alpha}, \sigma_{\alpha\beta}, \sigma_{\beta\beta})$ = $(1.00,\ 1.20,\ 1.40)$, and   ($\epsilon_{\alpha\alpha}, \epsilon_{\alpha\beta}, \epsilon_{\beta\beta}) = (1.0,\ 1.0,\ 1.0) .$ Each simulation is run at fixed number density ($\rho$) and temperature ($T$).  All data was collected after the system reaches thermal equilibrium. This model shows a strong increase of fragility with increasing density \cite{Berthier_2009,kinetic_fra}.  In particular, at low densities near the jamming transition the relaxation time exhibits Arrhenious temperature dependence (strong behavior) while at high density the relaxation time is super-Arrhenious (fragile behavior)~\cite{Berthier_2009}. We performed our simulations in the  density range $\rho \in [0.65, 0.82]$ to span a wide range of fragilities. In the following, length is measured in units of $\sigma_{\alpha\alpha}$, energy in units of $\epsilon_{\alpha\alpha}$, temperature in units of $\epsilon/\kappa_B$, time scale in units of $(\frac{m\sigma_{\alpha\alpha}^2}{\epsilon_{\alpha\alpha}})^{1/2}$,  masses and $\kappa_B$ are set to unity.

\section{\bf Softness Analysis}
\label{section:3}
Following earlier work, we use machine learning to identify a structural variable, softness, that characterizes the local structural environment surrounding each particle and correlates strongly with rearrangements. To do this, we must first characterize the local structural environment, then identify rearranging and non-rearranging particles for the training set. We then use a Support Vector Machine (SVM) \cite{SVM1,SVM2} to construct a hyperplane that best separates the rearranging particles in the training set from the non-rearranging ones.

In contrast to previous work~\cite{Schoenholz_Nat_Physics,Schoenholz263,Sussman10601} we do not apply the softness analysis to inherent structures (configurations obtained by quenching snapshots to zero temperature). Instead, we apply it directly to configurations of the thermal system. We find that this diminishes the classification accuracy slightly (to around $82\%$) compared to typical results for inherent structures in Lennard-Jones systems (around 88-90\%) but that softness still captures enough to be useful.

The local structural environment is quantified by three sets of structure functions. The first and second sets of structure functions are  closely related and are inspired by the symmetry functions that were proposed by Behler and Parrinello \cite{PhysRevLett.98.146401}. Both these structural descriptors are widely used  to quantify the particle local structure environment in previous studies \cite{cubuk_identifying_2015,Schoenholz_Nat_Physics,doi:10.1126/science.aai8830,Schoenholz263,Sussman10601,Sharp10943,PhysRevX.11.041019}. The first set of structural function measures the local density of a reference particle $\kappa$ at distance $\mu \pm L$.
\begin{equation}
G(\kappa;\mu) = \sum_i \exp^{-(r_{ik}-\mu)^2/L^2}
\label{density}
\end{equation}
$\kappa$ is the reference particle whose local environment we want to probe.
$\mu$ is the distance at which we want to measure the local density of particle $\kappa$. $r_{ik}$ is the distance between particle $\kappa$ and $i$. $G(\kappa;\mu)$ provides the radial density  of reference particle $\kappa$ at distance $\mu$. Here we keep $L = 0.1 \sigma_{\alpha\alpha}$ fixed. We take $\mu$ in between
$0$ to $3.0 \sigma_{\alpha\alpha}$ with an increment of $0.1 \sigma_{\alpha\alpha}$. The second set of structural descriptor captures bond angle information.
\begin{equation}
\Psi(\kappa;\xi,\lambda,\zeta)=\sum_{i,j} \exp^{-(r_{i\kappa}^2+r_{j\kappa}^2+r_{ij}^2)/\xi^2} (1+\lambda\cos{\theta_{\kappa ij}})^\zeta
\label{angular}
\end{equation}
where $\theta_{\kappa ij}$ is the angle between particles $\kappa$, $i$ and $j$ and $\lambda = \pm 1$. The angular resolution of the structure functions is controlled by the parameter $\zeta$.  The third set of structural descriptors includes the number of Delauney neighbors for the reference particle $\kappa$, the number of Delauney neighbors that are in contact (within interaction range) and the number of neighbors that are not in contact (beyond the interaction range) with the reference particle $\kappa$~\cite{doi:10.1063/5.0035395,https://doi.org/10.48550/arxiv.2011.13049}. Altogether we retain 73 structural descriptors. Note that this is overkill; our goal is not to obtain the most parsimonious description of local structure that yields useful information, but to have a definition expansive enough to yield high prediction accuracy across the entire range of system densities studied.

Once each particle's structural descriptors ($\mathbf{F_i}$ = $\{F_i^1,..,F_i^M\}$, here M = 73) have been chosen, the next step is to choose the training set, which consists of two equal-sized subsets, one containing rearranging particles and the other subset containing non-rearranging particles. Stochasticity plays a role in whether a particle will rearrange or not; to reduce the effects of stochasticity on the training set, we choose particles from the extremes of the distribution: to be labeled as a rearranging particle it must rearrange within the very short  cage-breaking time scale, $\tau_\beta$ \cite{doi:10.1063/1.5033555} in the future; to be a non-rearranging particle it must not have rearranged for a long time, $10\tau_\alpha$, in the past. We then use support-vector machine (SVM) method  with a linear kernel to construct the hyperplane $w \cdot \mathbf{F} - b =0 $ that best separates the rearranging from the non-rearranging particles.

The first step is to define what we mean by a rearrangement. We calculate $D_{min}^2(i)$ \cite{PhysRevE.57.7192}, the non-affine displacement of a particle $i$
\begin{multline}
D^2_{\mathrm{min}}(i,t,\Delta t )= \underset{ \epsilon_{\alpha \beta}}{\text{min}}\Big\{ \sum_\alpha \sum_j \big[ (r^\alpha_j(t)-r^\alpha_i(t))-\\ \sum_{\beta}(\delta_{\alpha \beta} -\epsilon_{\alpha \beta}) (r^\beta_j(0)-r^\beta_i(0)) \big ]^2 \Big\}
\end{multline}
Here $\alpha$ and $\beta$ denote spatial coordinates and $j$ indexes the neighbors of particle $i$ out to a distance of $R_D = 2.5$.  The first term represents the actual displacements of nearby particles relative to particle $i$, while the second term represents the relative displacements of neighboring particles if they were in a region of uniform strain. We then find the $\epsilon_{\alpha \beta}$ that minimizes $D^2 (t,\Delta t )$. The local deviation from affine deformation during the time interval $t$ and $\Delta t$, referred to as $D_{min}^2$, measures the non-affine component of local deformation. The  time interval $\Delta t$ is set to be longer than the ballistic time scale and close to the cage breaking time scale; we use $\Delta t = 12$.

We can now define a rearrangement. A particle is rearranging if $D_{min}^2$ is above than a threshold value $D^2_{\mathrm{min},0}$ \cite{cubuk_identifying_2015}.
We choose $D^2_{\mathrm{min},0}$ so that roughly 1 \% of small particles undergo a rearrangement at its low temperature state point (the low-temperature state point for each density has been chosen so that the systems' relaxation times are relatively similar). Following that, we use the same threshold for all subsequent temperatures for a specific density system. Table ~\ref{table:1} lists the values of $D^2_{\mathrm{min},0}$ for each density.

\begin{table}
\caption{}
 \begin{center}
  \begin{tabular}{ | c | c | c | }

 \hline
 Density & $D^2_{\mathrm{min},0}$ \\
 \hline
 $\rho = 0.82$  & 0.28     \\
\hline
$\rho = 0.78$  &  0.26    \\
\hline
$\rho = 0.72$  &  0.24    \\
\hline
$\rho = 0.69$  &  0.235    \\
\hline
$\rho = 0.65$  &  0.21    \\
\hline
\end{tabular}
\label{table:1}
\end{center}
\end{table}

\vskip +0.05in

We show the results for small particles in our binary mixture here, but our results have been verified for large particles (see Appendix ~\ref{Appendix_sce_1}).

Now that we have defined rearrangements, we must next assemble a training set and construct an SVM hyperplane. Our study requires obtaining a hyperplane that yields high prediction accuracy across the range of temperatures and densities of interest. To do this, we start by assembling a training set from the highest density system studied ($\rho=0.82$) at a low temperature (T = 0.0045). We then test the resulting hyperplane on systems at lower densities. We find that the classification accuracy is highest for the same density system and that it drops with decreasing density (see Table~\ref{table:highrho}) for the cross-validation accuracies at different densities for this hyperplane.

\begin{table}
\caption{}
\begin{center}
  \begin{tabular}{ | c | c | c | }
 \hline
 Density & Classification Accuracy \\
 \hline
 $\rho = 0.82$  & 85 \%    \\
\hline
$\rho = 0.78$  & 60 \%    \\
\hline
$\rho = 0.72$  & 51 \%    \\
\hline
$\rho = 0.69$  & 50.8 \%    \\
\hline
$\rho = 0.65$  & 50.6 \%    \\
\hline
\end{tabular}
\label{table:highrho}
\end{center}
\end{table}

We then repeat this procedure, assembling the training set from the lowest density system ($\rho=0.65$) at a low temperature (T = 0.00065) where the relaxation time $\tau_\alpha$ is fairly similar to that of the high-density system at $T=0.0045$. The cross-validation accuracy of the resulting hyperplane is shown in Table~\ref{table:lowrho}. In this case, the cross-validation accuracy is significantly lower at the higher densities.

\begin{table}
\caption{}
\begin{center}
  \begin{tabular}{ | c | c | c | }
 \hline
 Density & Classification Accuracy \\
 \hline
 $\rho = 0.82$  & 50 \%    \\
\hline
$\rho = 0.78$  & 50.1 \%    \\
\hline
$\rho = 0.72$  & 51 \%    \\
\hline
$\rho = 0.69$  & 56.5 \%    \\
\hline
$\rho = 0.65$  & 84 \%  \\
\hline
\end{tabular}
\label{table:lowrho}
\end{center}
\end{table}

Finally, we assemble a training set with half of the particles in the rearranging subset as well as half from the non-rearranging subset taken from the first training set, and the remaining half for each subset taken from the second training set. The resulting hyperplane performs very well in classifying rearranging and non-rearranging particles in a variety of fragile liquids at all the density state points studied (Table~\ref{table:allrho}).

\begin{table}
\caption{}
\begin{center}
  \begin{tabular}{ | c | c | c | }
 \hline
 Density & Classification Accuracy \\
 \hline
 $\rho = 0.82$  & 83 \%    \\
\hline
$\rho = 0.78$  &  82.5\%    \\
\hline
$\rho = 0.72$  &  82.4\%    \\
\hline
$\rho = 0.69$  &  82.1\%    \\
\hline
$\rho = 0.65$  &  82\%    \\
\hline
\end{tabular}
\label{table:allrho}
\end{center}
\end{table}
This hyperplane successfully classifies rearranging/non-rearranging particles based on their local structural descriptors over a wide spectrum of fragile systems. We use this hyperplane throughout our studies to obtain the results reported here.
The continuous variable softness of particle $i$ ($S_i = w \cdot \mathbf{F_i} - b $) is defined as the signed distance between the hyperplane and the position of local structural environment, as specified by a point in M-dimensional structure-function space.

\begin{figure*}
\centering
\includegraphics[scale=0.18]{Fig/Fig_1/Fig_1.pdf}
\caption{The probability of that a particle of softness $S$ will rearrange, $P_R(S)$ vs.~$S$ for systems ranging from strong to  fragile liquids (Fig.1(a) to Fig.1(e)); lines are to guide the eye and temperatures $T$ are indicated in the legend. Arrhenius plots of $P_R(S)$ as a function of $1/T$ for different softnesses ranging from $S \sim -3.5$ (blue) to $S \sim 4.5$ (red) for systems ranging from strong to fragile liquids (Fig.1(f) to Fig.1(j)). Lines are fits to the Arrhenius form (Eq.~\ref{eq:Arrhenius}). Vertical line in Fig.1(f)  indicates the $1/T_0 = 0$ location.}
  \label{fig:1}
\end{figure*}
\section{\bf Results}
\label{section:4}

It has been shown previously that fragility increases strongly with density in this family of systems~\cite{Berthier_2009,kinetic_fra}. We first recapitulate earlier analyses~\cite{Schoenholz_Nat_Physics} to calculate the probability that a particle is rearranging ($P_R$) as a function of its softness at each density (Fig.~\ref{fig:1} (a) to Fig.~\ref{fig:1} (e)). The probability of rearrangement is the fraction of particles that undergo rearrangement in the interval $\Delta t$ for a given softness $S$. At the lowest temperature studied for each system, $P_R(S)$ increases by several orders of magnitude as softness increases, showing that the softness is highly correlated with the probability of rearranging. Next,
in Fig.~\ref{fig:1} (f) to Fig.~\ref{fig:1} (j) we plot $P_R(S)$ as a function of $1/T$ for a range of softness values. As in previous studies~\cite{Schoenholz263,Sussman10601,Sharp10943,PhysRevE.101.010602,D0SM01575J} we see that $P_R(S)$ exhibits Arrhenius behavior for a given softness $S$:

\begin{equation}
P_R(S) = \exp \Big(\Sigma (S) -  \frac{\Delta E(S)}{T} \Big)
\label{eq:Arrhenius}
\end{equation}
where the effective configurational entropy barrier $(\Sigma (S))$ and energy barrier $(\Delta E(S))$ terms do not depend on temperature. The values of $(\Delta E(S))$ and $(\Sigma (S))$ are plotted in Fig.~\ref{fig:2} (a) and Fig.~\ref{fig:2} (b). Note that $\Delta E (S)$ and $\Sigma (S)$ become stronger functions of softness with increasing density (increasing fragility). Note also that $\Delta E(S)$ decreases with decreasing density; this is expected since particle overlaps decrease as the density is lowered below the unjamming transition. Furthermore, although $T_g$ decreases at low density, rescaling $\Delta E(S)$ by $T_g$ does not change these trends.

A consistent feature at every density studied is that the Arrhenius curves for $P_R(S)$ at different $S$ intersect at a reasonably well-defined temperature. In previous work, it was shown that this temperature corresponds to the onset temperature, $T_0$, below which liquids begin to show characteristic glassy dynamics, including dynamic heterogeneities, effective barrier heights that increase with decreasing temperature (for fragile systems) and non-exponential relaxation in time~\cite{Schoenholz_Nat_Physics}. The plots for $P_R(S)$ show that $T_0$ corresponds to the temperature at which structure no longer affects dynamics. As $T$ decreases below $T_0$, structure (softness) affects dynamics increasingly strongly.

Figs.~\ref{fig:1}(f)-(j) show that $T_0$ is a function of density, $\rho$. The density dependence is plotted in Fig.~\ref{fig:2}(c). Interestingly, $T_0$ is independent of $\rho$ for sufficiently high $\rho$, corresponding to fragile systems. Below $\rho \approx 0.69$,
however, $1/T_0$ decreases rapidly and even goes negative for the strongest systems studied.

Below $T_0$ (above $1/T_0$) the variation with $S$ of $P_R(S)$ or equivalently the effective energy barrier $\Delta E(S)$, implies that particles of higher softness face lower energy barriers and relax more rapidly, leading to heterogeneity of the relaxation. This heterogeneity is typically characterized using $\chi_4(t)$, the variance of the overlap function ($Q(t)$)
~\cite{doi:10.1146/annurev.physchem.51.1.99,PhysRevResearch.2.022067} as
\begin{equation}
\chi_4(t) = N \big [\langle Q^2(t) \rangle - \langle Q(t) \rangle ^2 \big].
\label{EQ:1}
\end{equation}
Where the self-overlap correlation function or two-point overlap correlation function is defined as
\begin{equation}
Q(t) = \frac{1}{N} \sum_{i=1}^{N} w\big( |\vec{r}_i(t) - \vec{r}_i(0)| \big).
\label{EQ:2}
\end{equation}
It calculates the amount of overlap between the initial ($t=0$) and subsequent configurations at time $t$. The window function $w(x)$ is defined as $w(x) = 1$, if $x<a$ and $0$ otherwise. The window function $w(x)$ was created to eliminate the decorrelation caused by particles' vibrational motions inside the local cages that their neighbors create. We chose $a=0.3$, which corresponds to the mean square displacement plateau value. The reported results are not sensitive to the precise choice of $a$.

 \vskip +0.05in
\begin{figure*
\includegraphics[scale=0.13]{Fig/Fig_2/Fig_2.pdf}
\caption{(a) Energetic and (b) entropic barriers implied by the Arrhenius behavior of $P_R(S)$ (Eq.~\ref{eq:Arrhenius}) as functions of softness $S$. Densities $\rho$ are indicated in the legend. (c) Inverse onset temperature, $1/T_{0}$ as a function of system density, with fragility increasing with density.
}
  \label{fig:2}
\end{figure*}

Fig.~\ref{fig:3}(a) shows representative examples of $\chi_4(t)$ for the highest- and lowest-density systems studied, showing that there is a peak in this quantity. The peak position, $\tau_\alpha$, is a measure of the relaxation time, while peak magnitude, $\chi_4^P$ quantifies the degree of dynamical heterogeneity.

The temperature dependences of $\chi_4^P$ and $\tau_\alpha$ are shown in Fig.~\ref{fig:3}(b)-(c), where we have rescaled $T$ by the corresponding glass transition temperature ($T_g$) for the system at each density. Here we define the glass transition temperature, $T_g$ as the temperature at which $\tau_\alpha = 10^7$. Values of $T_g$ are listed in Table \ref{table:Tg}.
\begin{table}
\caption{}
 \begin{center}
  \begin{tabular}{ | c | c | c | }
 \hline
 Density & $T_g$ \\
\hline
$\rho = 0.65$  &  0.000357   \\
\hline
$\rho = 0.69$  &  0.001195   \\
\hline
$\rho = 0.72$  &  0.001812    \\
\hline
$\rho = 0.78$  &  0.003163   \\
\hline
 $\rho = 0.82$  & 0.004031  \\
 \hline
\end{tabular}
\label{table:Tg}
\end{center}
\end{table}

\begin{figure}[htp]
\includegraphics[scale=0.245]{Fig/Fig_3/Fig_3.pdf}
\caption{(a) The dynamical heterogeneity function $\chi_4(t)$ as a function of time $t$ for strong ($\rho=0.65$, $T = 0.00065$ dotted line) and fragile ($\rho = 0.82$, $T = 0.0045$, solid line) liquids with the same relaxation times. The inset shows the peak value of $\chi_4$, $\chi_4^P$, as a function of the relaxation time $\tau_\alpha$ (time at the peak) for the strongest (solid line) and most fragile (dotted line) liquids studied. (b) The peak value  $\chi_4^P$ is plotted as a function of $T_g/T$, the inverse temperature made dimensionless by the glass transition temperature $T_g$. (c) Angell plot of $\tau_\alpha$ as a function of $T_g/T$. Densities are indicated in the legend. The behavior is Arrhenius at $\rho=0.65$ and becomes increasingly super-Arrhenius with increasing $\rho$. Curves are fits to the Vogel-Fulcher-Tamann form, Eq.~\ref{VFT}. (d) Kinetic fragility  (left y-axis) and $T_{VFT}$ (right y-axis), as defined in Eq.~\ref{VFT}, as a function of density.}
  \label{fig:3}
\end{figure}

Fig.~\ref{fig:3}(c) shows that fragility increases with density, as shown earlier in Ref.~\cite{Berthier_2009,kinetic_fra}. We fit the relaxation time $\tau_\alpha(T)$ using the Vogel-Fulcher-Tamann (VFT) formula \cite{10004192038,https://doi.org/10.1111/j.1151-2916.1925.tb16731.x,https://doi.org/10.1002/zaac.19261560121}
\begin{equation}
\tau_\alpha(T)=\tau_0 \exp\left [\frac{1}{K_{VFT} \ (T/T_{VFT}-1)}\right],
\label{VFT}
\end{equation}
where $\tau_0$ is the relaxation time at the onset temperature $T_0$. Here $T_{VFT}$ is the apparent divergence temperature for relaxation time and $K_{VFT}$ determines the degree of kinetic fragility, known as ``fragility index''. In Fig.~\ref{fig:3}(c) the curves are the VFT fits. Fig.~\ref{fig:3}(d) shows the kinetic fragility ($K_{VFT}$) (left axis, black points) and $T_{VFT}$ (right axis, blue points) as a function of density ($\rho$). We see that kinetic fragility ($K_{VFT}$) and $T_{VFT}$ provide essentially the same information and increase with increasing fragility.

Returning to Fig.~\ref{fig:3}(a), we note that it shows $\chi_4(t)$ for the strongest (solid line) and most fragile (dotted line) systems (the lowest- and highest-density systems), with temperatures chosen for each so that the two systems have the same value of $\tau_\alpha$. In the inset we plot $(\chi_4^P)$  as a function of $\tau_\alpha$ for the two systems. Evidently the more fragile system has a higher value of $\chi_4^P$ than the stronger system at the same relaxation time, as expected from the stronger dependence of $\Delta E$ and $\Sigma$ on $S$ for the more fragile system.

For fragile liquids the relaxation time increases more rapidly with decreasing temperature $T$ than for strong liquids. Because softness corresponds to an energy barrier, one might expect that the mean softness, $\langle S \rangle$ of a fragile system is a stronger function of temperature than that of a strong system. In Fig.~\ref{fig:4} (a), we plot $\langle S \rangle$ as a function of $T_g/T$. For fragile liquids, $\langle S \rangle$ decreases at a significantly faster rate with decreasing temperature than for strong liquids, as expected. Note that the definition of softness is identical for all systems studied, so this result indicates that structure changes more rapidly with temperature, leading to a more rapid change of relaxation time in fragile liquids than in strong ones.

\begin{figure}[htp]
\includegraphics[scale=0.126]{Fig/Fig_4/Fig_4.pdf}
\caption{(a) Average softness $\langle S \rangle$ as a function of $T_g/T$ for systems at densities indicated in the legend.  The temperature dependence of $\langle S \rangle$ increases with increasing fragility (density).
(b) The relation between $\Sigma(\langle S \rangle)$ and $\Delta E(\langle S \rangle)/T_0$ (described after Eq.~\ref{eq:sigmavsdeltae}) for $\rho=0.82$.}
  \label{fig:4}
\end{figure}

To better understand the connection between structure, energy barrier and dynamics, we study the free energy barrier:
\begin{equation}
    \Delta F(S) = \Delta E(S) - T \Sigma (S)
    \label{eq:free_energy}
\end{equation}

Note that at the onset temperature, $T_0$, $P_R(S)$ is independent of softness so from Eq.~\ref{eq:Arrhenius} we have
\begin{equation}
\Sigma(S) = \Delta E(S)/T_0 + C^\prime.
    \label{eq:sigmavsdeltae}
\end{equation}
where $C^\prime$ is a constant. In particular, at any temperature $T$, Eq.~\ref{eq:sigmavsdeltae} implies that $\Sigma(\langle S \rangle)= C^\prime + \Delta E(\langle S \rangle)/T_0$. This relation is verified in Fig.\ref{fig:4}(b) for the most fragile system, $\rho=0.82$; we find excellent agreement with $C^\prime=-0.35$.
Combining this with Eq.~\ref{eq:free_energy} we see that
\begin{equation}
    \Delta F(S)=(1-T/T_0) \Delta E(S) - T C^\prime.
      \label{eq:free_energy_new}
\end{equation}
If we were to assume that $\tau \propto 1/P_R(\langle S \rangle)$ as suggested by earlier work~\cite{Schoenholz263},  this would imply
$\tau \propto \exp({\Delta F(\langle S \rangle)}/T)$. We check the relation $\tau \propto 1/P_R(\langle S \rangle)$ by plotting the product $\tau_\alpha \times P_R(\langle S \rangle)$ as a function of inverse temperature in Fig.~\ref{fig:5} (b). We see that the product is not constant with temperature, implying that deviations from the simple relation $\tau \propto 1/P_R(\langle S \rangle)$. These deviations are pushed to lower temperatures with increasing fragility (density), partly due to the decrease in $T_0$, indicated by the vertical lines in Fig.~\ref{fig:5}(b). Thus, deviations are compressed into a narrower range, leading to a stronger $T$-dependence of the relaxation time, with increasing fragility. In all cases, however, such deviations are observed, showing that the dynamics are not simple, likely due to spatial and temporal correlations such as those that lead to kinetic heterogeneities (\textit{e.g.}~facilitation), which are stronger for more fragile systems (Fig.~\ref{fig:3}(a)).

\begin{figure}[htp]
\includegraphics[scale=0.13]{Fig/Fig_5/Fig_5.pdf}
\caption{(a) Data for the dimensionless free energy barrier $(\Delta F(\langle S \rangle)/T_g)$ for systems of different densities against inverse temperature ($T_g/T$). (b) Product of $\tau_\alpha$ and $P_R(\langle S \rangle )$ against inverse temperature ($T_g/T$). Vertical lines represent the re-scaled onset temperature ($T_g/T_0$) for systems at densities $\rho=0.69$ to $\rho=0.82$.}
  \label{fig:5}
\end{figure}

\section{\bf Discussion}
\label{section:5}

In summary, our results show that structure makes a strong contribution to the wide variation in fragility observed in glassy liquids. Our analysis suggests that three factors control fragility. Compared to strong liquids, in fragile ones (1) local structure ($\langle S \rangle$) decreases more strongly with decreasing temperature, (2) the free energy barrier $\Delta F(\langle S \rangle)$ rises more rapidly with decreasing $\langle S \rangle$, and (3) the upwards deviation of $\tau$ from $1/P_R(\langle S \rangle)$ is compressed into a smaller temperature range relative to $T_g$.  All of these factors contribute to the greater sensitivity of the relaxation time $\tau_\alpha$ to temperature in fragile liquids.

Our interpretation of $\Sigma(S)$ and $\Delta E(S)$ as entropic and energetic barriers stems from the Arrhenius behavior of $P_R(S)$. Direct confirmation that these quantities can be interpreted as barriers is difficult because calculations of energy barriers are extremely costly in a system with $O(e^{O(N)})$ minima once the number of particles, $N$, is of reasonable size. Direct calculations of entropic barriers are even more challenging but could be worthwhile in small systems to compare against our results

The result of Ref.~\cite{Schoenholz_Nat_Physics} that $T_0$ marks the onset of the dependence of $P_R(S)$ on softness, or local structure, is surprisingly powerful. It allows direct extraction of the entropic barrier (Eq.~\ref{eq:sigmavsdeltae}) and free energy barrier (Eq.~\ref{eq:free_energy_new}) to rearrangement that only include an temperature independent constant, given the energy barrier. Note that $T_0$ marks a crossover, not a transition, in low dimensions, so it is not sharply defined. Nonetheless, our results show that it is sufficiently well defined to be potentially useful as a way of extracting barriers.

Remarkably, we find that $T_0$ is actually negative for the strongest glassformers. Eq.~\ref{eq:sigmavsdeltae} shows that when $T_0$ is positive, a decreasing energy barrier with increasing $S$ implies an increasing entropic barrier with $S$. Thus, softer particles have lower energetic barriers and higher entropic barriers. However, when $T_0<0$, the opposite is true--softer particles have lower energy barriers and lower entropic barriers as well.
Note that for the strongest system at $\rho=0.65$, the energy barrier is essentially independent of softness (Fig.~\ref{fig:2}(a) ) and $T_0<0$ with a very small magnitude (Fig.~\ref{fig:2}(c)). As a result, the heterogeneity in particle mobilities comes from variation of the entropic barriers with $S$ (Fig.~\ref{fig:2}(b)), not of the energy barriers. In other words, at high densities/fragilities, heterogeneous mobility stems from the softness dependence of the energy barriers while at low densities/fragilities, it stems from the softness dependence of the entropic barriers. It is notable that a single hyperplane is predictive at all these densities and captures this change of behavior.

\subsection*{Conflicts of interest}
There are no conflicts of interest to declare.

\subsection*{Acknowledgements}
This research was supported by NSF through the University of Pennsylvania Materials Research Science and Engineering Center (MRSEC) (DMR-1720530) (IT) and through the Extreme Science and Engineering Discovery Environment (XSEDE) \cite{xsede} via TG-BIO200091. Additional support was provided by the Simons Foundation through the Simons Collaboration on ``Cracking the glass problem'' (Award 454945 to SAR and AJL) and through Investigator Award 327939 to AJL.

\subsection*{Data Availability Statement}
Raw data used to make plots and other data will be published at the following url [???] after the revision based on reviewer comments.