Maximal Relevance and Optimal Learning Machines
O Duranthon M Marsili R Xie



1]O Duranthon


2]M Marsili
3]R Xie
[1]Department of Physics, École Normale Supérieure, 24 rue Lhomond, 75005 Paris, France
[2]The Abdus Salam International Centre for Theoretical Physics, Strada Costiera 11, 34151 Trieste, Italy
[3]Key Laboratory of Quark and Lepton Physics (MOE) and Institute of Particle Physics, Central China Normal University (CCNU), Wuhan, China


Maximal Relevance and Optimal Learning Machines
    [
    
===============================================




We explore the hypothesis that learning machines extract representations of maximal relevance, where the relevance is defined as the entropy of the energy distribution of the internal representation. 
We show that the mutual information between the internal representation of a learning machine and the features that it extracts from the data is bounded from below by the relevance. This motivates our study of models with maximal relevance – that we call Optimal Learning Machines  – as candidates of maximally informative representations.
We analyse how the maximisation of the relevance is constrained both by the architecture of the model used and by the available data, in practical cases. We find that sub-extensive features that do not affect the thermodynamics of the model, may affect significantly learning performance, and that criticality enhances learning performance, but the existence of a critical point is not a necessary condition. On specific learning tasks, we find that i) the maximal values of the likelihood are achieved by models with maximal relevance, ii) internal representations approach the maximal relevance that can be achieved in a finite dataset and iii) learning is associated with a broadening of the spectrum of energy levels of the internal representation, in agreement with the maximum relevance hypothesis.




§ INTRODUCTION


Statistical mechanics models have been used in statistical learning since the pioneering works on associative memory <cit.> and on Boltzmann Learning Machines <cit.>. Many attempts to make sense of the spectacular performance of learning machines such as deep neural networks, have focused on understanding their statistical mechanics properties, as compared e.g. to spin glasses <cit.>. 
Several attempts have been made to understand the typical properties of learning machines <cit.>. 

These properties emerge after the interaction parameters have been learned from highly structured datasets <cit.>, and they depend on what measure of error or likelihood is used in supervised or unsupervised learning, respectively. 
This makes it hard to discuss the typical properties of learned models without reference to the particular dataset on which the model is trained, the objective function assumed, the algorithms or the architecture used. This suggests that, in order to understand learning, we should focus on the ingredients which makes it possible <cit.>.



This paper takes a different perspective, and it explores the claim of Refs. <cit.> that learning maximises the relevance of the representation that the machine extracts from the data. Let us briefly recall the rationale behind this conclusion, which is articulated in Ref. <cit.>, to which we refer for more details.
In order to define the relevance, let us consider a generic unsupervised learning task that maps the data to the internal states of a learning machine, in order to capture the structure of statistical dependencies of the data.
We assume that the data  is generated from an unknown distribution p(), that we shall call the generative model. 

A learning machine maps the data  to the internal states  of the hidden layer(s) of the learning machine. Training with data drawn from () induces a distribution p() on the hidden states. 
We define the energy E_=-log p() as the coding cost associated to a microscopic state . The average energy

    ⟨ E⟩=-∑_ p()log p()≡ H[]

coincides with the entropy of , and it is a measure of the resolution of the representation. Ref. <cit.> argues that the entropy of the energy distribution

    H[E]=-∑_E p(E)log p(E),     p(E)=∑_p() δ(E-E_)

can be taken as a quantitative measure of the information that the model contains on the structure of the data, or of its generative model. In order to motivate this conclusion, we observe that H[E] is typically small for Gibbs distributions, which encode states of maximal ignorance, because they have a very narrow distribution of energies. The internal representations p() of a learning machine should be as far as possible from a maximum entropy model. Hence it should feature a broad spectrum of energies E_, that corresponds to a large value of H[E]. In this representation, two states  and ' with very different energies E_ and E_' necessarily code for significantly different patterns in the data, because they cannot occur typically, within the same maximum entropy model. 

Following Cubero et al.. <cit.>, we shall henceforth call H[s] resolution and H[E] relevance. 
The above discussion leads to the hypothesis that machines that are trained on highly structured data should approach the limit where H[E] is as large as possible. In order to explore this claim, we define the class of Optimal Learning Machines (OLM) as the ideal limit of those statistical models that maximise the relevance at a given resolution H[]

    max_{E_}: H[s]=E̅ H[E]

over the distribution of energy levels E_. The optimisation in Eq. (<ref>) is constrained, in practical cases, by the architecture of the machine. When these constrains are neglected, machines that satisfy Eq. (<ref>) have a broad distribution of energies with an exponential density of states <cit.>.
This, in turn, is equivalent to statistical criticality <cit.>, which is a widely observed statistical regularity in efficient representations <cit.>. This corroborates the conjecture that the principle of maximal relevance applies to real learning machines. 


The contributions of this paper is twofold: first we show that, in a generic learning task, the relevance lower bounds the mutual information between the internal representation and the features that the machine extracts from the data. This provides a rationale for why learning machines should approach the limit of maximal relevance, in order to extract the maximal amount of information from the data. 

Secondly, we investigate how the optimisation in Eq. (<ref>) is constrained both by the architecture of the machine used and by the amount of data on which it is trained. 
We start exploring the relevance of different architectures. Gaussian learning machines <cit.> offer an interesting counter-example, where the relevance H[E] does not vary during training. This is consistent with the fact that the learned distribution p() remains a Gaussian, irrespective of the generative model . Hence, Gaussian learning machines do not learn anything on the generative model , beyond its parameters. 

Taking H[E] as a quantitative measure of learning performance, we study a class of toy learning machines with discrete energy levels, in order to characterise the architectural features associated with a high value of H[E]. Our results suggest that features that confer superior learning performance may be associated to sub-extensive features that are not accessible to standard statistical mechanics approaches. A superior learning performance is not necessarily related to the existence of a critical point separating two different phases, but when a critical point exists, learning performance improves when the model is tuned to its critical point.

Finally, we present a series of numerical experiments, in order to understand how closely real learning machines approximate the ideal limit of OLM when the dataset is finite. Within the class of toy learning machines with discrete energy levels discussed earlier,we show that the models that have superior learning performance in two distinct (unsupervised learning) tasks are those with the highest values of H[E]. Then we move to Restricted Boltzmann Machines (RBM). There we argue that the number of samples on which the machine is trained limits the resolution of the energy scale. When the resolution and relevance are evaluated at the scale dictated by the available data, our numerical experiments show that RBM converge to states of maximal relevance, in agreement with <cit.>. 
For small systems, it is possible to compute the full spectrum of energy levels and to follow its evolution. We show that during learning the distribution of energy levels indeed broadens. We summarise and comment our results in the concluding section.



§ RELEVANCE AND HIDDEN FEATURES














Let us assume that data points =(x_1…, x_d) are generated as draws from an unknown distribution[
Here and in the rest of the paper, backslashed symbols (e.g. ) refer to unknown entities.] () that is characterised by a rich structure of dependence between the components x_a. 
For example, the MNIST dataset <cit.> of hand written digits has d=28× 28 = 784 pixels, for each of which, the grey scale is codified in an integer x_a in the range [0,255]. In spite of the fact that  belongs to a very high-dimensional space, the MNIST dataset spans a manifold of only d_ int≈ 13 intrinsic dimensions <cit.>. This implies that () should be a model of strongly interacting variables, and that there should be a set of coordinates ()=(_1(),…,_d_ int()) that describes the relevant variation in the dataset. We refer to  as the hidden features. As for , the backslash indicates that  is an unknown theoretical construct that represents statistically significant patterns in the dataset[In the case of the MNIST dataset, a draw from the generative model  is a theoretical abstraction for the process of hand writing a digit by a human  and features are not necessarily related to “meaningful” patterns such as the values of digits. Deep learning machines also rely on features that look like noisy patterns, such as those used in adversarial examples. As shown by Ilyas et al. <cit.>, 
machines that only use features which are reasonable for humans perform worse, in terms of generalisation, than those that use all features. 
]. 

We focus on stochastic models, such as Restricted Boltzmann Machines (RBM) or Deep Belief Networks (DBN), with a finite number of discrete internal states. The goal of learning is to find a statistical model p() over a discrete variable , and a mapping p(|) from  to , such that the generating distribution 

    p()=∑_p(|)p()

is as close as possible to the empirical one (in unsupervised learning). The objective function employed may differ, yet in all cases, training aims at approximating the unknown generative process p()≈() as closely as possible[In supervised learning the representation aims at reproducing a functional relation x_ out=f(x_ in) between two parts of the data =(x_ in,x_ out) as well as possible. This is similar to unsupervised learning with ()=(x_ in)δ(x_ out-f(x_ in)).]. Here we focus on the outcome of this process, abstracting from algorithmic details on how the representation is learned. We remark that, unlike  which is a theoretical abstraction, p() is a proper statistical model, as well as the induced distribution p() in Eq. (<ref>). We define the energy of state  as

    E_=-log p()

which is the cost (in nats) for coding state . 


If the learning machine captures correctly the structure of the data  with which it has been trained, then it must extract features ϕ that approximate well the hidden features . Hence there must be a function ϕ() defined on the states of the hidden layers – the extracted features – such that its projection 

    ⟨ϕ|⟩≡∑_ϕ()p(|)

is close to the hidden features (). We call ϕ the extracted features, to contrast them with the hidden ones .

In an ideal situation, the distribution p() should contain no further information than that contained in ϕ(). In other words, ϕ should be the sufficient statistics for the distribution p(). This means that the distribution of  on all states with ϕ()=ϕ̅ must be a maximum entropy distribution, i.e. p(|ϕ̅)=A(ϕ̅)δ_ϕ(),ϕ̅, where A(ϕ̅) is the inverse of the number of states  with ϕ()=ϕ̅. Put differently, states with the same features ϕ()=ϕ(') should be statistically indistinguishable, i.e. p()=p(') for all  and ' such that ϕ()=ϕ('). Therefore E_=E_' for any two states  and ' such that ϕ()=ϕ('), which means that E_ must be a function of ϕ(). 
In other words, the statistical dependence of the variables can be represented as the Markov chain[The notation z→ y→ w implies that, if y is known, w and z are independent, i.e. that I(z,w|y)=0 <cit.>. This holds trivially true if w is a function of y, because then w is not a random variable, when y is known.] 

    →ϕ→ E.

The data processing inequality <cit.> then implies that 

    I(,ϕ)≥ I(,E)=H[E],

where I(X,Y) is the mutual information between variables X and Y <cit.>, and the last equality comes from the fact that I(,E)=H[E]-H[E|] and H[E|]=0. Therefore, H[E] provides a lower bound to the information I(,ϕ) that the representation contains on the extracted features. Notice that ϕ is a function of , and hence I(,ϕ)=H[ϕ]. 



Going back to the hidden features , we observe that the training process can be represented by the Markov chain

    →→→ϕ

where  is the representation that unfolds from training. If ϕ are sufficient statistics, which means that  and  are independent conditional on ϕ, then the last step of the chain can be inverted, i.e. →→ϕ→. Then the
data processing inequality implies that I(,)≤ I(,ϕ)=H[ϕ]. Under the assumption that the variation H[ϕ] of the extracted features is entirely induced by the hidden features, it is reasonable to expect that this inequality is tight[If the last step of the Markov chain (<ref>) can be reversed, then I(,)=I(ϕ,)=H[ϕ]-H[ϕ|]=I(,ϕ)-H[ϕ|]. Therefore 
I(,)≈ I(,ϕ) is equivalent to the statement that the extracted features are approximately a function of the hidden features, i.e. H[ϕ|]≈ 0.] 

The inequality (<ref>) is the main result of this Section. It connects the mutual information I(,ϕ) of the extracted features with the relevance H[E].
Apart from simple cases, such as the Gaussian learning machine discussed below, it is not clear how to identify the sufficient statistics of a distribution p(). So the extracted feature ϕ are hardly accessible in practice, 
let alone the hidden ones . The distribution p() instead can be computed for learning machines trained on structured data,  and hence H[E] can be estimated, as we shall see in the following Sections.
The  inequality (<ref>) suggests that statistical models with a high value of H[E] are natural candidates for good learning machines, even though there is no guarantee that the inequality (<ref>) is tight. 

The comparison between different learning machines on the basis of H[E] makes sense only if the representations that they generate have the same resolution H[], i.e. the same coding cost.
For this reason, we shall consider the relevance H[E] as a measure of the learning capacity, at a given level of resolution H[]. An absolute maximum value for the maximum of H[E] for a given H[] can be obtained by studying the problem

    max_{W(E)}: H[] H[E]

where W(E) is the number of energy levels with E_=E, which is constrained by the two conditions

    H[]     =     -∑_p()log p()=∑_EW(E)e^-EE 
    
    N     =    ∑_ 1=∑_EW(E)

where N is the number of available states. This problem is studied in Appendix <ref>, where we show that the solution features and exponential density of states W(E)=W_0e^(1+ν)E, where ν is a constant that depends on H[]. 
This is an ideal limit, because in practical cases the maximisation is additionally constrained by the architecture of the learning machine and by the available data. 

Notice that H[|ϕ] provides a measure of that part of H[] which is not informative on the features. By Eq. (<ref>), this is upper bounded by H[|E]. Maximising H[E] at fixed resolution H[] implies minimising H[|E]≥ H[|ϕ], i.e. squeezing noise out of the representation. 
An alternative definition of OLM, are machines that extract the most compressed representations from data, at a minimal information content H[E] on the features

    min_I(,ϕ)≥ H[E]H[].


We remark that H[E] is an unambiguous measure of information content when E_ takes values in a discrete set. When E_ can take any value on the real axis (as e.g. in RBMs), it is necessary to define the relevance with respect to a precision Δ <cit.>. In order to do this, we define the probability 

    p_Δ(E)=∑_: |E_-E|≤Δ/2 e^-E_

that a state  has energy in the interval of width Δ around E. The relevance at precision Δ is defined as

    H_Δ[E]=-∑_E p_Δ(E)log p_Δ(E)

where the sum on E is restricted to integer multiples of Δ. For a real random variable E we expect that 
H_Δ[E]≃ h[E]-logΔ for Δ→ 0, where h[E] is the differential entropy <cit.>. 
We shall take H_Δ[E] as a measure of relevance, by making sure that Δ is chosen in the range where the relation H_Δ[E]≃ h[E]-logΔ is satisfied.



§ EXPLORING DIFFERENT ARCHITECTURES


In this section, we explore classes of statistical models in order to understand how the relevance H[E] depends on the architecture. Under the hypothesis that good learning machines are those with high values of H[E], this sheds light on the architectural features that confer superior learning capacity to models. 



 §.§ The Gaussian learning machine and maximum entropy models



Consider a Gaussian learning machine with m visible units =(x_1,…, x_m) and n hidden units =(s_1,…,s_n) which are all continuous variables. The joint probability distribution is given by

    p(,)=1/Z(,ŵ)e^-1/2 (-)^2 -  1/2^2 +' ŵ

where prime denote transpose, =(b_1,…, b_m) is a vector of parameters, ŵ is an n× m matrix of couplings, and Z(,ŵ) is a normalisation constant.





The couplings  and ŵ are adjusted so as to learn data generated from  (). The marginal distribution of  is given by a Gaussian

    p()=(2π)^n/2/Z(,ŵ)e^-1/2 (-)^2 +1/2'ŵŵ ' .

As shown in <cit.>, this machine learns the first two moments of the distribution of . The dynamics of the parameters  and ŵ under different training algorithms has been discussed in Ref. <cit.> in detail. In what follows,  and ŵ are the parameters of the trained machine, so they depend on the data. The distribution in the internal layer is[The largest singular value of ŵ can be shown to be less than one, so that this is well defined.]

    p()=(2π)^m/2/Z(,ŵ)e^-1/2' (1-ŵ'ŵ) +'ŵ

so the energy is

    E()     =    1/2' (1-ŵ'ŵ) -'ŵ -m/2log (2π)+log Z(,ŵ) 
    
         =     E_0(ŵ)+1/2 (-_0)' (1-ŵ'ŵ) (-_0)

where 

    E_0(ŵ)=n/2log (2π) -log√( det(1-ŵ'ŵ))

is a constant (the ground state energy) and _0(,ŵ)=(1-ŵ'ŵ)^-1ŵ is a vector. 
Notice that the resolution is given by:

    H[]=E_0(ŵ)+n/2.

The distribution of energies is

    p(E)=1/Γ(n/2)(E-E_0)^n/2-1e^-(E-E_0).

Remarkably, p(E) is independent of  and it depends on ŵ, i.e. on the data, only through the constant E_0.
Hence the differential entropy <cit.> of this distribution 

    h(E)=n/2 +logΓ(n/2)-(n/2 -1)logψ(n/2)≃1/2log (4π e n)
    
    
    +O(1/n)

is independent of the parameters , ŵ, i.e. it is independent of the data. The scaling h(E)≃1/2log n for large n is consistent with the fact that asymptotically the distribution of energies is a Gaussian with variance proportional to n. 














Let us now comment these results in terms of the discussion in Section <ref>. In the Gaussian learning machine, the extracted features ϕ()={s_i, i=1,…,n, s_i s_j, i≤ j=1,…,n} coincide with the sufficient statistics. A Gaussian learning machine is a maximum entropy model of the form

    p(|θ)=1/Z(θ)e^θ·ϕ(),

where the parameters θ=(,ŵ) are adjusted so as to reproduce the means and covariances of the data. 
As a consequence, the generative model p() is always Gaussian, irrespective of what  is. A Gaussian learning machine does not learn anything on the shape of , because it assumes its shape from the outset.
It only learns the parameters θ. The amount of information learned on the features I(,ϕ) equals H[], because  is a part of ϕ. 
Therefore, training increases the resolution H[] but it leaves the relevance H[E] unchanged. 
The inequality (<ref>) I(,ϕ)≥ H[E] remains true, but is not informative. 

An explicit manifestation of the poor performance of Gaussian learning machines, as compared e.g. to RBMs, in reproducing a dataset of handwritten digits, is given in Fig. <ref>.
These findings suggest that the Gaussian multivariate model has a special place in statistical inference
 
 and the conclusions drawn from it should be taken with care. 


With respect to the picture in Fig. <ref>, the Gaussian learning machine can only learn one point[Conversely, the expected value ϕ̅ over the whole dataset of ϕ() is sufficient to estimate the parameters of this model.] ϕ̅, and it cannot reproduce a wide variation of features, as in Fig. <ref>.  

We expect the same to be generally true for exponential models of the type of Eq. (<ref>). Indeed, for large n, the distribution of the energy E_ is very narrow, because with very high probability

    E_≈⟨ E ⟩_θ=H[].

Schwab et al. <cit.> have shown that if the parameters θ themselves vary, as e.g. in

    p()=∫dθ p(|θ)p(θ).

then the energy acquires a broad variation and it exhibits an exponential density of state[Ref. <cit.> shows that this result is a general consequence of the Asymptotic Equipartition Property in e.g. <cit.>).
If all typical states have the same probability p(|θ)∼ e^-nH[s], then their number should be W∼ e^nH[s].] log W(E)∝ E. This is a defining feature of representations p() that maximise the relevance. 
 In the present picture, the variation of the hyper-parameters θ is induced by the distribution () of the data . In other words, the hidden features  play the role of the 
 “hidden variables” of Ref. <cit.>.
 


 §.§ Ising Learning Machines


In this Section we shall consider a learning machines whose hidden variables  are described by a fully connected model of Ising spins. These are maximum entropy models. So they are not particularly efficient as learning machines. Yet they provide a well defined playground for investigating the properties that make some architectures more efficient (i.e. with a larger H[E]) than others. This also allows us to contrast the properties of models with large H[E] with those of well known models, such as the mean field Ising ferromagnet and the J_i,j=± J spin glass. These insights make it possible to test whether models with the maximal H[E] perform better than others in this class, on a specific application to real data. This will be the subject of Section <ref>.

We consider learning machines with a distribution of internal states given by 

    p()=1/Ze^-ℋ(),     Z=∑_s e^-ℋ(),

where =(s_1,…,s_n), s_i=± 1 are n spins. The Hamiltonian takes the form

    ℋ()=-∑_i<jJ_i,js_is_j,

where the matrix Ĵ of couplings has elements J_i,j∈{± J}. 
The joint distribution  over the visible variables ∈{0,1}^m and the internal state  is given by 
p(,)=p(|)p(), where 

    p(|)=1/Q_exp[∑_j=1^m  (b_j+∑_i=1^ns_iw_i,j)x_j],
         Q_=∏_i=1^m[
    1+e^b_j+∑_i s_i w_i,j].

The parameters , ŵ and J are learned from the data, e.g. by maximum likelihood maximisation. Here we concentrate on the distribution p() of the internal variables .

We define the energy as 
the coding cost E_=-log p() =ℋ()+log Z, and address the optimisation problem

    𝒥= argmax_Ĵ: H[]=E̅H[E],

where the maximum is taken over all matrices with elements J_i,j=± J, and on J. 𝒥 is the subset of such matrices that achieve a maximal value of H[E].







In order to compute H[E] as a function of Ĵ, 

we define the distribution of the energy as

    p(E)=∑_ p()δ(E-E_)=W(E)e^-E

where W(E) is the number of states  with energy E_=E. From this and Eq. (<ref>), we compute the value of H[E]. 










The rest of this section is devoted to describe the trade-off between H[E] and H[] for models in the class of Eq. (<ref>). Since H[]≤ nlog 2 within this class, we shall discuss the rescaled resolution 

    h_s=H[]/nlog 2.

Analogously, we observe that (E_-log Z)/J takes values only on the integers [-n(n-1)/2,n(n-1)/2]. This justifies the introduction of the rescaled relevance 

    h_E=H[E]/log n

that takes values in the interval [0,2). 
The sub-extensive nature of H[E] suggests that properties that contribute to the relevance of a model may not be accessible to saddle point analysis of the partition function, that focuses only on the leading extensive terms. 
The behaviour of H[E] on the number n of spins depends on the fact that the model (<ref>) contains only one parameter J[For example, in a model with ℋ(s)=-∑_i<jJ_i,js_is_j-∑_i h_i s_i, J_i,j=± J and h_i=± h, the energy takes values on a set of point E=J n_J+h n_h with n_J∈ [-n(n-1)/2,n(n-1)/2] and n_h∈ [-n,n]. So the number of possible values of E grows at most as n^3, and we expect the h_E∈ [0,3). This suggests that H[E] grows in general with log n, with a coefficient that is bounded by a constant that grows linearly with the number of parameters.]. 

The information content of a model can be quantified in terms of the Stochastic Complexity in Minimum Description Length <cit.>. This in our case is given by 

    SC=∫dJ/J√(𝕍_J[E]),

where 𝕍_J[E] is the variance of the distribution p(E), 







which is closely related to the specific heat. 

A second measure of the information content of a model can be obtained by the area under the h_E vs h_s curve. This scores models according to their relevance across different levels of resolution. We shall call this measure Multi-Scale Relevance (MSR), following  <cit.>, that introduced a similar measure on single samples.

Both h_s and h_E only depend on Ĵ through the degeneracy W(E) of energy levels. W(E) is invariant under the permutation of the spins. Hence if Ĵ∈𝒥 then also the matrix with elements J_π(i),π(j) is a solution, where π(i) is a permutation of the indices i=1,…,n. Also, matrices Ĵ and Ĵ' that are related by a gauge symmetry J_i,j'=τ_iJ_i,jτ_j, with τ_i=± 1, have the same W(E). Hence, if Ĵ∈𝒥, then also its gauge transformed Ĵ' belongs to 𝒥. We partially exploit these symmetries by fixing the gauge with the choice J_1,j=1 for all j≠ 1. Henceforth we shall focus on the reduced set 𝒥 of matrices with this choice of the gauge. 

An interesting point that we shall discuss is how the set 𝒥 of solutions Ĵ are organised, i.e. whether they correspond to isolated maxima or whether they form a connected set of nearby solutions in terms of local moves. 
We consider local moves that correspond to “flipping” one coupling J_ı,j→ - J_i,j, or that “swap” a pair of opposite couplings J_i,j=-J_k,l, changing their sign (J_i,j,J_k,l)→ (-J_i,j,-J_k,l). Two models Ĵ and Ĵ' are connected if there is a sequence of such local moves that changes Ĵ into Ĵ'. In this way, the set 𝒥 can be divided into connected components. 
This analysis may shed light on the accessibility of solutions of Eq. (<ref>) to local learning rules  during training, and it connects to the finding of wide flat minima <cit.> in real learning machines[The latter refers to the fact that the landscape of the function that is optimised during training – e.g. the error or (minus) the likelihood – is characterised by wide flat minima in the space of parameters. If optimally trained machines correspond to models 
Ĵ of maximal relevance, then we expect that the set 𝒥 has a large connected component.]. 



  §.§.§ Exact enumeration for small systems



Each model in Eq. (<ref>) is defined by a choice of the sign of the n(n-1)/2 couplings J_i,j and by the strength J of the couplings. 
There are 2^n(n-1)/2 possible ways of choosing the signs. Yet, this number can be reduced by fixing the gauge. 


We were able to compute h_E as a function of h_s for all models up to n=9 spins. This allows us to find those models that achieve a maximal h_E for a given resolution h_s.



Representative results are shown in Fig. <ref> for n=8. There are M=219 different curves in the left panel, each corresponds to a different degeneracy W(E) of energy levels[For n=3,4,5,6,7,8 and 9, we find M=2,3,7,16,54,219 and 1625, respectively.]. 
Fig. <ref>(left) shows, in particular, the six functions h_E that achieve maximal relevance for some value of h_s. A representative choice of the corresponding matrix Ĵ∈𝒥 is shown in the right part of Fig. <ref>. We remark that there are many models  with the same W(E), that cannot be reduced to one another by a permutation of the spins (for model a there are 325 different models). 
Fig. <ref>(right) shows the degeneracy W(E) for model a) (top left in Fig. <ref> right). As a comparison, Fig. <ref> also shows the curve h_E(h_s) (left) and the degeneracy W(E) for the ferromagnet (J_i,j=J ∀ i,j). 

The architectures that achieve maximal relevance, that we identify with OLM, have the following distinctive properties: i) the large majority of the couplings J_i,j are positive, which means that OLM are close to ferromagnetic models. Indeed, the ground state of all OLM are the same as that of the ferromagnet. Yet, ii) the relevance for OLM achieves a much higher value than that of the ferromagnetic model (dotted line in Fig. <ref> left). This is consistent with the fact that W(E) for the ferromagnet is concentrated on few values of E, whereas for OLM it spreads over a larger number of values of E. For model a) W(E) is consistent with an exponential behaviour, as predicted in <cit.>. Finally iii) the architectures of OLM appear to be rather inhomogeneous, with negative J_i,j impinging on a small subset of nodes. Such a low degree of symmetry, may be the origin of the fact that the number |𝒥| of matrices Ĵ that share the same W(E) of OLM is rather large (see Table <ref>). 




















It is interesting to analyse the properties of the set 𝒥 with respect to local dynamics, as discussed earlier. 
We find (see Table <ref> for n=8) that 𝒥 is dominated by a single largest component, which is reminiscent of the wide flat minima discussed in <cit.>. 
As h_s decreases, we see a tendency of the set to shrink and to fragment in different connected components. Yet, 
the union 𝒥_a∪𝒥_b∪…∪𝒥_f of all OLM ∈[a,b,…,f], is composed of a single connected component.

Finally, the inset of Fig. <ref> left) shows that models with a higher relevance, as measured by the MSR, typically have a larger Stochastic Complexity (SC).






  §.§.§ Some particular architectures for large n



This section analyses specific architectures for which the curve h_E(h_s) can be computed for moderately large values of n. The matrices Ĵ of the four models that we consider are schematically depicted in Fig. <ref> (left), together with the corresponding estimate of their capacity h_E. We discuss here only the main results and refer to the Appendix <ref> for more details.



  
The mean field Ising ferromagnet and spin glasses

The mean field Ising ferromagnet (MFIFM) corresponds to J_i,j=J for all i,j (model a in the left panel of Fig. <ref>). For large n, this is characterised by a high temperature disordered phase for J<J_c/n and by a low temperature ordered phase for J>J_c/n where a non-zero magnetisation spontaneously appears. The energy for the MFIFM can take at most n+1 different values, which gives an upper bound h_E≤ 1. 


Fig. <ref> (right) reports the curves h_E(h_s) for several values of n. The curve exhibits a maximum in the neighbourhood of the critical point J≈ J_c/n. The decreasing part to the right of  the maximum corresponds to the high temperature phase, whereas the one to the left of the maximum to the low temperature phase. As n increases, the maximum shifts to values of h_s closer and closer to one. At the same time, the maximum gets sharper and sharper and its value slowly decreases towards a finite limit (see inset). 

The limiting form of the curve h_E(h_s) can be computed in the limit n→∞. This shows that 

    lim_n→∞ h_E={[      1/2 J≠ J_c/n;      3/4 J=J_c/n, ].

which is consistent with the bound h_E≤ 1. The analysis leading to this result (see Appendix <ref>) confirms that h_E depends on sub-leading contributions that arise from the integration of fluctuations around the saddle point value.
At the critical point J=J_c/n, the resolution converges to h_s→ 1. Yet, depending on how the limit h_s→ 1 is taken, all values of h_E∈[1/2,3/4] can be achieved. Indeed, for J=0 we have h_s=1 and h_E=1/2, exactly. Fig. <ref> (right) is fully consistent with this result, although the convergence is very slow. 

When J_i,j=± J are chosen at random with equal probability (model b in the left panel of Fig. <ref>), the model coincides with a spin glass. We know that the relevant scale for J is 1/√(n). The ground state energy is extensive so E_ spans a number of energy levels that is at most of order n^3/2. Therefore we expect that h_E≤ 3/2 for a spin glass. In fact, for J=c/√(n) with c finite, the distribution of the energy extends over a range of order Δ E∝√(n). Since the level spacing equals 2J∼ 1/√(n), this suggests that the number of energy levels on which p(E) extends grows linearly with n. Hence we conjecture that h_E→ 1 as n→∞ for a model with randomly chosen J_i,j. This conclusion is corroborated by numerical results for systems of up to n=28 spins.




  
The star model

Let us now consider a model that differs from a MFIFM by that fact that one spin has antiferromagnetic interactions with ℓ< n other spins (model c in the left panel of Fig. <ref>). 



The degeneracy W(E) of energy levels is shown in Fig. <ref>(left) for n=256, ℓ =2 and ℓ=128. As ℓ increases, the number of different values that E takes also increases. For the MFIFM (ℓ=0) E can take only n/2+1 different values (for even n), whereas for ℓ=n/2 we find that the number of possible values of E grows almost as n^2. 
At the same time, W(E) acquires a rapid variation as a function of E, so that for large n, W(E) becomes a space filling curve (see Fig. <ref> left and inset). Interestingly, the thermodynamics of the model is dominated by the convex envelope of log W(E), which is the same as that of the MFIFM. Therefore the star model and the MFIFM have indistinguishable thermodynamic properties. This is consistent with the fact that the number of J_i,j=-1 is a fraction of order 1/n of the total number of interactions and that they impinge on one out of the n spins. It is easy to check that, for ℓ≤ n/2, the ground state of the star model is the ferromagnetic state s_i=s_j for all i≠ j.

The maximal value of h_E^*(ℓ)=max_h_sh_E for even values of n is achieved at ℓ=n/2 (see inset of Fig. <ref> right), whereas for odd values of n, h_E^* reaches a lower maximum at ℓ≈ n/3. The symmetry for ℓ→ n-ℓ of these curves is a consequence of the gauge transformation discussed above.
The curves h_E(h_s) for ℓ=n/2 and for different values of n are shown in Fig. <ref>(right). 
We find that the highest relevance H[E] for even n and ℓ=n/2, grows faster than log n. 
An extrapolation of the maximal value of h_E is consistent with h_E→ 5/4 as n→∞. At fixed values of h_s instead, we observe a slower growth of H[E] with log n, that is consistent with h_E→ 1 for large n. 
Note that the degeneracy of each energy level E_m of the MFIFM is spread over n other energy levels, √(n) of which contribute to H[E] (see inset of <ref> left). This suggests a relation h_E^ Star=h_E^MFIFM+1/2, in agreement with numerical results.

Finally we remark that the set 𝒥 of star models is connected in a single component, under a local dynamics updates of Ĵ. 




















  
The nested model

In the next architecture that we consider, the couplings in the lower right corner of the matrix Ĵ are negative, and those on the diagonal or above are positive (model d) of Fig. <ref> left), i.e. J_i,j= sign (n+ℓ-i-j+1/2).
In these models, the set of negative J_i,j is organised according to a nested network <cit.>.
Their properties are discussed in Appendix <ref> where we also derive a recursion relation for W(E) that allows us to compute h_E for moderately large values of n. Fig. <ref> shows the results for different values of n and for ℓ=1 (left) and n/2 (right). 



Although we could not access values of n as large as in the previous cases, Fig. <ref> shows that this model has qualitatively distinct features. First, we found that log W(E) approaches a continuous function of E and it has support on a number of energy levels that increases as n^2. For ℓ=1 we found that W(E)=W(-E) is a symmetric function of E (see inset of Fig. <ref> left). Second, in both cases, we found that h_E increases monotonically with n and it approaches a limiting value already for relatively small values of n. The maximal limiting values of h_E is higher than that achieved with the other architectures (h_E≈ 1.57). We note that the number of anti-ferromagnetic interactions in these models is proportional to n^2, and hence the thermodynamic properties of nested models in the n→∞ limit differ from those of the MFIFM (even though the ground state is the same for ℓ≤ n/2).

There are n!/(n-ℓ)! possible models depending on how the set of n-ℓ fully ferromagnetic spins are chosen and on the ranking of the ℓ spins. This set of models forms a single connected component under local updates of the J_i,j. This is consistent with the conjecture that models with a high learning capacity are easily accessible by local updates of the couplings.



§ NUMERICAL EXPERIMENTS


This section addresses the issue of whether the maximum relevance principle applies to learning experiments on real datasets. On the basis of this principle, we expect that i) models that achieve a higher value of H[E] should describe better complex datasets, ii) within the precision allowed by a finite dataset, training converges to models of maximal relevance and iii) that learning induces broad distributions of energy levels. We address these three issues in turn.



 §.§ Complex data is best described  by models with maximal relevance



A prediction of the maximum relevance principle is that models that better fit a dataset are those with the highest relevance H[E], provided the data contains a rich statistical structure. In order to test this prediction, we compared the ability of Ising models with different matrices Ĵ to represent a dataset, maximising the likelihood at different levels of resolution. 
 
Because of the exponential growth[To give an idea of the combinatorial complexity of the class of models, for n=8 there are 219 different classes, each with many distinct models that cannot be reduced one to the other by a gauge transformation or by a permutation of the spins. For example, there are 325 models with the same H[E] - H[] curve as model a in Fig .<ref>.] in the number of models with the number n of spins, we limit our analysis to the case n=5. There are seven different classes of models (see Table <ref>). Each class contain models with the same W(E) and hence with the same curve h_E-h_s. The class that achieves maximal value of H[E] in a wide range of H[], that we shall call class A, contains two models[One with only one negative coupling and the other with three negative couplings incident on the same node, e.g. J_1,2=J_1,3=J_1,4<0.]. Class A also has the largest values of the SC and of the MSR. The other classes – that we call B, C, D, E, FM  and AF – contain one (E, AF and FM) or two models (B, C and D). Class FM contains the ferromagnetic model and AF the anti-ferromagnetic one.  














We focus on two datasets of m=9 binary variables ∈{0,1}^m. One is that used in <cit.> on the voting patterns in the US Supreme Court. This reports the votes (in favour or against) of the m=9 judges of the US Supreme Court on N=895 cases. The second is  obtained from the MNIST dataset <cit.>, focusing on a patch of 3× 3 pixels of the N=6· 10^4 hand written digits[Each MNIST data point contains 28× 28 pixels, with an integer value in the range [0,255]. We first transform the dataset into coarse grained pixels of 2× 2 original pixels, and binarise the result applying a threshold. Then we focus on the 3× 3 patch with pixels coordinates (i,j) ranging from 5 to 7.], as in <cit.>. 
The joint distribution is built by taking p() as in Eqs. (<ref>,<ref>) and p(|) as in Eq. (<ref>). 
For every model and each value of J we maximise the log-likelihood over the parameters  and ŵ by gradient ascent. Fig. <ref> shows that over a wide range of resolutions, the maximum of the likelihood is attained for models in class A. In particular, this interval contains the point where the log-likelihood achieves its maximum. 



 §.§ RBMs achieve the maximal relevance compatible with a dataset



In order to test whether learning machines converge to models of maximal relevance during training, we focus on the architecture of Restricted Boltzmann Machines (RBM). These have n hidden nodes and n_v visible nodes, where n_v matches the dimensionality of the dataset. We study a reduced MNIST dataset with n_v=144 and N=6· 10^4 data points. Details on RBMs and on the dataset are provided in Appendix <ref>.

In general, the 2^n energy levels are learned from a sample of N≪ 2^n data points. The finiteness of the data limits the precision with which the energy spectrum can be estimated. In order to take this into account, we observe that, with a change of variables f=e^-E, the differential entropies of f and E stand in the relation

    h[E]=-∫_0^1 df p(f)log[p(f)|df/dE|]=h[f]+H[]

where H[]=⟨ E⟩ is the average energy. The same relation holds for the variables at any precision Δ, i.e. 

    H_Δ[E]=H_Δ[f]+H[]

The right hand side of Eq.(<ref>) can be estimated in a finite sample , using the empirical distribution f̂_=k_/N, where k_ is the number of times that state  occurs in the sample. 
Specifically, H_Δ[f] can be estimated using the empirical distribution of frequencies P{f_=k/N}=k m_k/N, where m_k is the number of states  that occur k times in the sample. Hence

    H_Δ[f]≃ -∑_kkm_k/Nlogkm_k/N≡Ĥ[k] .

Likewise, the second term in Eq .(<ref>) can be computed as 

    H[] ≃ -∑_k_/Nlogk_/N≡Ĥ[].

The two quantities Ĥ[] and Ĥ[k] have been introduced earlier as a measure of resolution and relevance, respectively, within a sample <cit.>. 
Taken together, Eqs. (<ref>, <ref>) and (<ref>), imply

    H_Δ[E]≃Ĥ[k]+Ĥ[] .

We remark that the precision Δ is implicitly defined by the sample. The fact that both Ĥ[k] and Ĥ[] are biased estimates of the entropy is immaterial for our purposes. We remind that our goal is that of establishing whether a learning machine approaches states of maximal relevance, for a given dataset. This can be done  comparing the values computed from a sample ŝ drawn from p(), with the theoretical maximal value of Ĥ[k]+Ĥ[] that can be achieved in a sample of N points. This has been computed in Ref, <cit.>.



Fig. <ref> reports the results of numerical experiments on the reduced MNIST dataset, where the number of hidden nodes of the RBM ranges from n=15 to 70. Each training session lasts 300 epochs with persistent contrastive divergence <cit.>. 
Sample digits generated by the RBM are displayed in Fig. <ref> (right). 
As a reference, the curve of maximal Ĥ[k]+Ĥ[] <cit.> is shown as a full black line. For each value of n, we observe trajectories that converge very fast to points that are close to the theoretical curve of maximal Ĥ[k]+Ĥ[ ], as first observed in <cit.>. As learning proceeds, the internal state of the machine evolves towards more and more compressed representations (i.e. Ĥ[] decreases).





  §.§.§ RBM develop broad distribution of energies



Here we analyse the evolution of the energy distribution to test whether the distribution of energy broadens with learning, as predicted by relevance maximisation. The setting is the same of the previous section. We focus of RBM for small values of n, for which the distribution of the energy levels E_, and hence  H_Δ[E], can be computed exactly during training. This allows us to monitor the evolution of the system as well as to appreciate the constraints imposed by the finiteness of the dataset and by the architecture. In order to do that, we take as a reference an ideal learning machine, that corresponds to the unconstrained maximisation of H_Δ [E], at a fixed H[]. The way in which this is computed parallels Ref. <cit.> and is detailed in Appendix <ref>. As discussed in <cit.>, the distinguishing feature of an ideal learning machine is that the degeneracy W(E) of energy levels follows an exponential behaviour W(E)≃ W_0 e^ν E. 

RBMs with n≤ 24 hidden nodes have been trained with (K=10) contrastive divergence <cit.> on the data. We used mini-batches of 10^3 samples for the results shown in Fig. <ref>. Although these are not ideal for learning[We confirmed that, as argued in <cit.>, best performance is obtained with mini-batches of ∼ 10 samples. Yet the trajectories obtained when the mini-batch size is small are very noisy.], they allow us to monitor the dynamics of the RBM during learning, in the plane spanned by H_Δ[E] and H[].  
As a function of Δ, the resolution is expected to behave as  H_Δ[E]=h[E]-logΔ, where h[E] is the differential entropy <cit.> of the probability density function of energy levels. Indeed we find that H_Δ[E]+logΔ is nearly constant in the range 0.1<Δ<1.




At different points in time, we compute the full spectrum of energy levels and compute H_Δ [E] and H[]. Fig. <ref> (left) suggests that learning occurs in cycles, reminiscent of those discussed in Ref. <cit.>: A compression phase where H[] decreases, followed by an expansion phase where the system comes closer to the ideal limit of a machine with maximal H_Δ[E] (dashed line). 

Fig. <ref> (right) shows the evolution of the spectrum of energy levels. As the bottom panel (b) shows, the spectrum of energy levels expands more and more during training. At times (n_e=1, n_e=10) when the system is  closer to the ideal limit (dashed line in Fig. <ref> left), the degeneracy of energy levels follows approximately an exponential behaviour W(E)≃ W_0e^ν E, whereas it departs from it in the successive compression phases. In the late stages of the dynamics, for a given input , the distribution p(|) of internal states  is sharply peaked. This makes it possible to map the dataset into a set of (clamped) internals states of the machine, as done in <cit.>, and to compute the corresponding distribution N(E) of energies. As shown in the top right panel d) of Fig. <ref>, the evolution of this distribution during learning follows the one of p(E) (panel c) and converges to a similar shape.



§ CONCLUSION


Maximal relevance has been proposed as a principle that characterises optimal learning machines <cit.>. A distinctive feature of this approach is that it allows us to discuss the properties of learning machines without any reference to the data they have been trained with, as long as this contains a sufficiently rich structure. 


The first contribution of this paper is to provide a relation between the relevance and feature extraction in learning. We show that the relevance H[E] provides a lower bound to the mutual information between the state of the machine and the features that a learning machine extracts (see Section <ref>). This does not exclude the existence of efficient learning machines with low values of H[E]. Rather it ensures that learning machines with high values of H[E] provide at least that much information on the features. 

The second contribution of this paper is to investigate architectures that achieve maximal relevance, and to confirm that, when the constraints imposed by the finiteness of the dataset and by the architecture are taken into account, the internal representations of well trained learning machines approach the limit of maximal relevance. 



We show that properties that confer a superior learning capacity, as measured by H[E], may be  related to sub-extensive quantities in the statistical mechanics treatment. This suggests that the properties responsible for the spectacular performance of learning machines may not be accessible to approaches that focus only on leading thermodynamic terms. A proper understanding of these systems requires a detailed analysis of saddle point fluctuations. This result also tallies with the observation that complexity reveals itself in the sub-extensive contribution to the entropy, both in statistical models[Minimum Description Length, as well as Bayesian inference, identify the leading order term that should be used in model selection to penalise the likelihood of models for their complexity. In both cases this term is d/2log N, where d is the number of parameters of the model and N is the number of data points. This term is sub-leading with respect to the log-likelihood, which is proportional to N.] <cit.> and in time series[The complexity of a time series is defined in <cit.> as the predictive information, which is the mutual information between the recent past, in a window of T time points, and the future of the time series.] <cit.>. 

For Ising systems, we find that local update rules can explore a large part of the set of OLM, suggesting that the relevance H[E] is characterised by wide flat maxima in the space of the parameters. It is tempting to relate this property to the wide flat minima discussed by Baldassi et al. <cit.>. 

Our results also contribute to the discussion on the criticality of learned models <cit.>. We confirm that models have a superior learning performance when poised at the critical point. Yet, it is not strictly necessary for an OLM to be poised at a critical point marking the transition between two different phases. Indeed, we find recursive models with a learning performance superior to that of models that feature a phase transition, in a wide range of the resolution scales. 

We hope that these results will not only contribute to our understanding of learning machines, but that they may also pave the way to applications that may improve further the performance of learning machines. 




§ ACKNOWLEDGMENTS

We benefitted from discussions with Jean Barbier, Riccardo Zecchina, Federico Ricci-Tersenghi and Yasser Roudi. We are grateful to the authors of Re. <cit.> for sharing their data with us. We thank an anonymous referee of MSML2020 for pointing out the example in Section <ref>.





§ THE IDEAL LIMIT OF AN OPTIMAL LEARNING MACHINE



For a  learning machine with many hidden units (n≫ 1), we can approximate the distribution of energies with a continuous probability density p(E), that we assumed to be defined in an interval E∈ [0,E̅]. 
For Δ small enough, the behaviour of H_Δ[E] with Δ is given by 

    H_Δ [E]≃ h[E]-logΔ

where h[E] is the differential entropy <cit.>

    h[E]=-∫_0^E̅ dE p(E)log p(E).

For a fixed Δ, the ideal limit of a learning machine with maximal H_Δ[E] is then obtained by finding the p(E) that makes h[E] maximal.

Let W(E) be the density of energy levels, i.e. W(E)dE is the number of states s with E_s∈ [E,E+dE). The probability density function p(E) is then given by p(E)=W(E)e^-E. When no other constraint is imposed, Optimal Learning Machines are defined by the solution of the optimisation problem

    W^*(E)     =     argmax_W(E) h[E] 
           s.t.   H[]=⟨ E⟩
               ∫_0^E̅ W(E)dE=N
               ∫_0^E̅ W(E)e^-EdE=1,

where

    h[E]=-∫_0^E̅ dEp(E)log p(E)=⟨ E⟩-⟨log W(E)⟩ ,

and N is the total number of states .
The solution of this problem is found introducing Lagrange multipliers to enforce the constraints (<ref>,<ref>), and fixing E̅ to satisfy Eq.  (<ref>). The solution reads

    W(E)     =    ν/e^νE̅-1e^(1+ν)E̅, 
    
    H[]     =    ⟨ E⟩ = E̅/1-e^-νE̅-1/ν,
    
    h[E]     =    loge^νE̅-1/ν-νE̅/1-e^-νE̅+1,
    
    N     =    ν/1+νe^(1+ν)E̅-1/e^νE̅-1,

where Eq. (<ref>) needs to be solved to obtain E̅ as a function of N and 
Eq. (<ref>) yields ν as a function of H[]. The limit ν→ 0 is revealing to understand the scaling of H[] and h[E] on N. We find

    E̅    =    log (1+NE̅)≃log N+loglog N+O(loglog N/log N)     (ν=1)
    
    H[]     =    1/2E̅≃1/2log N
    
    h[E]     =    logE̅≃loglog N.

At this point, the curve h[E] vs H[] is close to its maximum. So Eq. (<ref>) provides an upper bound on the possible value of h[E]. In the case of a Boltzmann machine with n hidden binary units, there are N=2^n hidden states. Hence H[] grows linearly with n and H[E]∝log n grows with the logarithm of the number of hidden units n. This is the same scaling that we have found for discrete energy levels, with Δ=J.






§ SPIN MODELS





  
The MFIFM model

The energy for the MFIFM ranges over an interval [-Jn(n-1)/2,Jn(n-1)/2] of order n^2, but energies E_ can take at most n+1 different values 

    E_m=-J/2(m^2-n),     m=∑_i=1^ns_i=-n,-n+2,…,n-2,n.

This implies that h_E≤log(n+1)/log n. Note also that there is a single matrix Ĵ which corresponds to the degeneracy of energy levels W(E) of the MFIFM (apart from gauge transformations)..

The distribution of energy levels is given by p(E)=W(E)e^-E/Z, where Z=∑_E W(E)e^-E. Hence the relevance is given by

    H[E]= log Z -⟨log W(E)⟩+⟨ E⟩.

The degeneracy of energy level E_m in Eq. (<ref>) is 

    W(E_m)=2nn+m/2≃√(2/π n(1-μ^2))
    e^n S(μ),    μ=m/n∈ [-1,1]

where S(μ)=-1+μ/2log1+μ/2-1-μ/2log1-μ/2. In the regime where h_s∈ [0,1] is finite, J is of order 1/n and, for n→∞, the partition function is dominated by the saddle point 
μ^*= argmin_μ f(μ) with f(μ)=Jnμ^2/2-S(μ). 
For Jn≠ 1, the first term in Eq. (<ref>) can be computed with integration over the gaussian fluctuations 

    Z≃√(2n/π)∫_-1^1 dμ/√(1-μ^2) e^-nf(μ)≃√(2/π (1-μ^*^2) f"(μ^*))e^-nf(μ^*)

Here, the √(n) factor is canceled because f(μ)-f(μ^*)∼ (μ-μ^*)^2, and the change of variables z=√(n)(μ-μ^*) generates a 1/√(n) term.

Hence log Z≃ -nf(μ^*)+ const. The extensive term in log Z is canceled by an analogous term in

    ⟨log W(E)⟩-⟨ E⟩≃ n f(μ^*)+log√(2π n(1-μ^2))

so, to leading order, H[E]≃1/2log n+ const, which implies h_E→ 1/2 as n→∞.

For Jn=1, instead, f(μ)-f(μ^*)≃ a(μ-μ^*)^4+…. This implies that in the calculation of Z we need a change of variables z=n^1/4(μ-μ^*) that yields log Z≃1/4log n+ const. This additional term is responsible for the asymptotic behaviour h_E→ 3/4 for Jn=1.



  
The spin glass model

We computed numerically the value of max_J h_E for 100 realisations of models with random J_i,j=± 1 (with P{J_i,j=+1}=1/2), for systems of up to n=28 spins. A linear fit of H[E] versus log n yields h_E≃ 1.02 ± 0.01. 



  
The star model

In the start model we take J_i,j=-1 if i=1 and j≤ℓ+1, and J_i,j=+1 otherwise. We divide the n-1 group of spins, except spin s_1, into two groups: the ℓ ones with anti-ferromagnetic interaction with s_1 and the remaining n-ℓ-1. If q and k are the number of positive spins in the first and second group, respectively, the energy is given by

    ϵ(k,q,s_1) = J[s_1(n-1-2ℓ-2k+2q)-1/2(2k+2q-n+1)^2+n-1/2].

The degeneracy of energy levels is given by

    W(E) = ∑_k=1^n-1-ℓ∑_q=1^ℓ∑_s_1=± 1ℓ qn-1-ℓ kδ_E,ϵ(k,q,s_1).

Notice that for ℓ=0 one recovers the MFIFM. A gauge transformation J_i,j→τ_i J_i,jτ_j, with τ_1=-1 and τ_i=+1 for i>1, maps a model with  ℓ>n/2 into a model with ℓ'=n-ℓ≤ n/2. So it is sufficient to study the model for ℓ≤ n/2.

Each models 𝒥 with the star structure can be “reached” from any other model, under a dynamics where a single J_i,j↦ -J_i,j changes sign. If we consider a “rewiring” dynamics where the signs of a pair (J_i,j,J_i,k)=(-1,+1) of links are swapped (J_i,j,J_i,k)↦ (-J_i,j,-J_i,k), each of the ℓ negative links can be rewired to any of the n-ℓ-1 spins with all ferromagnetic interactions. This dynamics, spans the whole subset of models in 𝒥 with the same hub node. These subsets are related by a permutation of the nodes.



  
The nested model

The triangle model has 

    J_i,j={[                 -1 for  i>ℓ, j> n+ℓ-i;                 +1          otherwise ].

where ℓ=1,…,n-1. For ℓ=1 this corresponds to a matrix Ĵ where all the elements below the diagonal i+j=n+1 are negative, as shown in Fig. <ref>(left). 
In these models, the set of negative J_i,j is organised according to a nested network <cit.>.



These matrices can be obtained through a recursive procedure (see Fig. <ref> right): we start from n=2 with J_1,2=+1. At each iteration, from a system of n spins, we add two spins, s_0 and s_n+1.
The spin s_0 has J_0,j=-1 with all other spins j≤ n. The spin s_n+1 has J_i,n+1=+1 with all spins, including s_0. This generates a matrix with n+2 spins with the desired structure[
For each state at energy E and magnetisation M=∑_i=1^n s_i, and each choice of s_0, s_n+1, we obtain a state of the system with n+2 spins with energy

    E'=E-J[s_0s_n+1-s_0 M+s_n+1M)]

and magnetisation M'=M+s_0+s_n+1. 
Thus the energy degeneracy of a system of n spins can be computed as W_n(E)=∑_Mw_n(E,M),
where w_n satisfies the recursive equation

    w_n+2(E',M')= ∑_s_0,s_n+1=± 1∑_E,M 
    w_n(E,M)δ_E'(E)δ_M',M+s_0+s_n+1.

Here the shorthand δ_E'(E) reduces the sum to only those terms where Eq. (<ref>) is satisfied.].

It is possible to generalise the recursion relation to obtain matrices with different values of ℓ. For example, the case ℓ=n/2 can be obtained by separating the set of spins into two equal parts, ℐ_≤={i≤ n/2} and ℐ_>={i> n/2}. To each of the two sets we add two spins at each iteration, thereby obtaining a system with n+4 spins. Of the two spins added to ℐ_≤, one has J_0,j=-1 with all spins j∈ℐ_≤ and J_0,j=+1 for all j∈ℐ_>. The other three spins have ferromagnetic interactions with all other spins, as well as among themselves. It is easy to see by induction, that the ground state of this model is the ferromagnetic one s_i=s_j for all i,j. This holds true for all ℓ≤ n/2.

Numerical iteration of the recursion relations allowed us to compute the degeneracy W_n(E) and the curves h_E(h_s) for different values of n and ℓ. 

For a given value of ℓ, there are n!/(n-ℓ)! nested models depending on how the set of n-ℓ fully ferromagnetic spins are chosen and on the ranking of the ℓ remaining spins in the nested hierarchy. This set of models forms a single connected component under a dynamics where a negative link J_i,j=-1 is “rewired” to a different node k, i.e. (J_i,j,J_k,j)=(-1,+1)↦ (+1,-1). In order to show this, let ∂_i be the set of spins with J_i,j=-1. If we sort spins in such a way that δ_1=(1,…,ℓ), δ_2=(1,…,ℓ-1),…,δ_k=(1,…,ℓ-k+1),…, then the rewiring of the link (k,ℓ-k+1) to (k+1,ℓ-k+1) interchanges the order of spins k and k+1 in the hierarchy. The rewiring of the link (1,ℓ) to (1,k) with k>ℓ changes the membership of spins in the network connected by J_i,j=-1. Repeated application of these two moves in the appropriate order makes it possible to reach any model in the set of nested models, from any other model.



§ RBMS AND THE REDUCED MNIST DATASET



We focus on the architecture of Restricted Boltzmann Machines (RBM). These have n_v visible  units =(x_1,…, x_n_v), and n hidden binary units =(s_1,…,s_n) with s_i=0,1. The joint probability distribution is 

    p(,)=1/Ze^∑_i b_i s_i+∑_j c_j x_j+∑_i,js_i w_i,j x_j

where Z=Z(,,ŵ) is the partition function. The parameters , and ŵ are learned from the data, in such a way as to maximise the log-likelihood p(). The energy of hidden state  is defined as

    E_=-log p()=log Z-∑_i b_i s_i-∑_j log[1+e^c_j+∑_is_i w_i,j].

Since E_ can take any real value, we measure the relevance with respect to a precision Δ of the energy spectrum, as in Eq. (<ref>). 

In what follows, we study a reduced MNIST dataset <cit.>, where we first coarse grain the original 28× 28 pixels data in cells of 2× 2 pixels, that are transformed into binary values by applying a threshold[If the sum of the grey levels of the four pixels exceeds 400, the coarse grained pixel is assigned the value one, otherwise it is zero.]. Finally we discard the first and the last rows and column of cells, thereby leaving us with a dataset of N=60000 images of 12× 12 binary pixels, corresponding to n_v=144. The reason for working with such a reduced dataset is that it allows us to compute the exact distribution of energies, for  values of n in a  relevant range. Indeed, the evaluation of the energy requires the calculation of the partition function Z, which can be computed only for small systems.  Sect. <ref> presents results for RBMs with a number of hidden units such that the exact calculation of p_Δ(E) on all  2^n of states is within reach of our computational resources. 
For larger systems we sample N≪ 2^n states  from p() with MCMC method and estimate H_Δ(E) from the sample (see Sect. <ref>). 



In order to test whether learning machines converge to models of maximal relevance during training, we focus on the architecture of Restricted Boltzmann Machines (RBM). These have n_v visible  units =(x_1,…, x_n_v), and n hidden binary units =(s_1,…,s_n) with s_i=0,1. The joint probability distribution is 

    p(,)=1/Ze^∑_i b_i s_i+∑_j c_j x_j+∑_i,js_i w_i,j x_j

where Z is the partition function, which depends on the parameters , and ŵ that are learned from the data, in such a way as to maximise the log-likelihood p(). The energy of hidden state  is defined as

    E_=-log p()=log Z-∑_i b_i s_i-∑_j log[1+e^c_j+∑_is_i w_i,j].

Since E_ can take any real value, we measure the relevance with respect to a precision Δ of the energy spectrum, as in Eq. (<ref>). 

ieeetr



