Shaping the learning landscape in neural networks around wide flat minima

The Simple Example of Discrete Weights.

In the case of binary weights it is well known that even for the single-layer network the learning problem is computationally challenging. Therefore, we begin our analysis by studying the so-called binary perceptron, which maps vectors of N inputs ξ∈−1,1N to binary outputs as σW,ξ=signW⋅ξ, where W∈−1,1N is the synaptic weights vector W=w1, w2, …, wN.

Given a training set composed of αN input patterns ξμ with μ∈1,…,αN and their corresponding desired outputs σμ∈−1,1αN, the learning problem consists of finding a solution W such that σW,ξμ = σμ for all μ. The entries ξiμ and the outputs σμ are random unbiased independent and identically distributed variables. As discussed in ref. 6 (but see also the rigorous bounds in ref. 7), perfect classification is possible with probability 1 in the limit of large N up to a critical value of α, usually denoted as αc; above this value, the probability of finding a solution drops to zero. αc is called the capacity of the device.

The standard analysis of this model is based on the study of the zero-temperature limit of the Gibbs measure with a loss (or energy) function LNE that counts the number of errors (NE) over the training set:LNE=∑μ=1αNΘ−σμσW,ξμ,[1]where Θx is the Heaviside step function, Θx = 1 if x>0 and 0 otherwise. The Gibbs measure is given byPW = 1Zβ exp−β∑μ=1αNΘ−σμσW,ξμ,[2]where β≥0 is the inverse temperature parameter. For large values of β, PW concentrates on the minima of LNE. The key analytical obstacle for the computation of PW is the evaluation of the normalization factor, the partition function Zβ:Zβ = ∑wi=±1exp−β∑μ=1αNΘ−σμσW,ξμ.[3]In the zero-temperature limit (β→∞) and below αc the partition function simply counts all solutions to the learning problem,Z∞ = limβ→∞Zβ=∑WXξW,[4]where XξW = ∏μ=1αNΘσμσW,ξμ is a characteristic function that evaluates to one if all patterns are correctly classified, and to zero otherwise.

Z∞ is an exponentially fluctuating quantity (in N), and its most probable value is obtained by exponentiating the average of logZ∞, denoted by logZ∞ξ, over the realizations of the patterns:Z∞,typical≃expNlnZ∞ξ.[5]The calculation of logZ∞ξ was done in the 1980s and 1990s by the replica and the cavity methods of statistical physics and, as mentioned above, the results predict that the learning task undergoes a threshold phenomenon at αc=0.833, where the probability of existence of a solution jumps from one to zero in the large N limit (6). This result has been put recently on rigorous grounds by ref. 7. Similar calculations predict that for any α∈0,αc, the vast majority of the exponentially numerous solutions on the hypercube W∈−1,1N are isolated, separated by a ON Hamming mutual distance (8). In the same range of α, there also exist an even larger number of local minima at nonzero loss, a result that has been corroborated by analytical and numerical findings on stochastic learning algorithms that satisfy detailed balance (9). Recently it became clear that by relaxing the detailed balance condition it was possible to design simple algorithms that can solve the problem efficiently (10⇓–12).

Local Entropy Theory.

The existence of effective learning algorithms indicates that the traditional statistical physics calculations, which focus on the set of solutions that dominate the zero-temperature Gibbs measure (i.e., the most numerous ones), are effectively blind to the solutions actually found by such algorithms. Numerical evidence suggests that in fact the solutions found by heuristics are not at all isolated; on the contrary, they appear to belong to regions with a high density of nearby other solutions. This puzzle has been solved very recently by an appropriate large deviations study (3, 4, 13, 14) in which the tools of statistical physics have been used to study the most probable value of the local entropy of the loss function, that is, a function that is able to detect the existence of regions with an ON radius containing a high density of solutions even when the number of these regions is small compared to the number of isolated solutions. For binary weights the local entropy function is the (normalized) logarithm of the number of solutions W′ at Hamming distance D N from a reference solution W:EDW = −1NlnUW,D[6]withUW,D = ∑W′XξW′δW′⋅W,N1−2D[7]and where δ is the Kronecker delta symbol. In order to derive the typical values that the local entropy can take, one needs to compute the Gibbs measure of the local entropy:PLEW = 1ZLE exp−yEDW,[8]where y has the role of an inverse temperature. For large values of y this probability measure focuses on the W surrounded by an exponential number of solutions within a distance D. The regions of HLE are then described in the regime of large y and small D. In particular, the calculation of the expected value of the optimal local entropySD ≡ EDopt=maxW−1NlnUW,Dξ[9]shows the existence of extremely dense clusters of solutions up to values of α close to αc (3, 4, 13, 14).

The probability measure Eq. 8 can be written in an equivalent form that generalizes to the nonzero errors regime, is analytically simpler to handle, and leads to novel algorithmic schemes (4):PLEW ∼ PW;β,y,λ = Zβ,y,λ−1ey ΦW,β,λ.[10]where ΦW,β,λ is a “local free entropy” potential in which the distance constraint is forced through a Lagrange multiplier λ:ΦW,β,λ = ln∑W′e−βLNEW′−λ dW,W′,[11]where d⋅,⋅ is some monotonically increasing function of the distance between configurations, defined according to the type of weights under consideration. In the limit β→∞ and by choosing λ so that a given distance is selected, this expression reduces to Eq. 8.

The crucial property of Eq. 10 comes from the observation that by choosing y to be a nonnegative integer, the partition function can be rewritten asZβ,y,λ = ∑Wey ΦW,β,λ = ∑W∑W′aa=1ye−βLRW,W′a,[12]whereLRW,W′a = ∑a=1yLNEW′a − λβ∑a=1ydW,W′a.[13]These are the partition function and the effective loss of y+1 interacting real replicas of the system, one of which acts as reference system (W) while the remaining y (W′a) are identical, each being subject to the energy constraint LNEW′a and to the interaction term with the reference system. As discussed in ref. 4, several algorithmic schemes can be derived from this framework by minimizing LR. Here we shall also use the above approach to study the existence of WFMs in continuous models and to design message-passing and greedy learning algorithms driven by the local entropy of the solutions.

Two-Layer Networks with Continuous Weights.

As for the discrete case, we are able to show that in nonconvex networks with continuous weights the WFMs exist and are rare and yet accessible to simple algorithms. In order to perform an analytic study, we consider the simplest nontrivial 2-layer neural network, the committee machine with nonoverlapping receptive fields. It consists of N input units, one hidden layer with K units and one output unit. The input units are divided into K disjoint sets of Ñ=NK units. Each set is connected to a different hidden unit. The input to the ℓ-th hidden unit is given by xℓμ=1Ñ∑i=1Ñwℓiξℓiμ, where wℓi∈R is the connection weight between the input unit i and the hidden unit ℓ and ξℓiμ is the i-th input to the ℓ-th hidden unit. As before, μ is a pattern index. We study analytically the pure classifier case in which each unit implements a threshold transfer function and the loss function is the error loss. Other types of (smooth) functions, more amenable to numerical simulation, will be also discussed in a subsequent section. The output of the ℓ-th hidden unit is given byτℓμ=signxℓμ = sign1Ñ∑i=1Ñwℓiξℓiμ.[14]In the second layer all of the weights are fixed and equal to one, and the overall output of the network is simply given by a majority vote σoutμ=sign1K∑ℓτℓμ.

As for the binary perceptron, the learning problem consists of mapping each of the random input patterns (ξℓiμ), with ℓ=1,…,K, i=1,…,Ñ, μ=1,…,αN, onto a randomly chosen output σμ. Both ξℓiμ and σμ are independent random variables that take the values ±1 with equal probability. For a given set of patterns, the volume of the subspace of the network weights that correctly classify the patterns, the so-called version space, is given byV=∫∏iℓdwℓi∏ℓδ∑iwℓi2−Ñ∏μΘσμσoutμ.[15]where we have imposed a spherical constraint on the weights via a Dirac δ in order to keep the volume finite (though exponential in N). In the case of binary weights the integral would become a sum over all of the 2N configurations and the volume would be the overall number of zero error assignments of the weights.

The committee machine was studied extensively in the 1990s (15⇓–17). The capacity of the network can be derived by computing the typical weight space volume as a function of the number of correctly classified patterns αN, in the large N limit. As for the binary case, the most probable value of V is obtained by exponentiating the average of log⁡V, Vtypical≃expNlog⁡Vξ, a difficult task which is achieved by the replica method (18, 19).

For the smallest nontrivial value of K, K=3, it has been found that above α0≃1.76 the space of solutions changes abruptly, becoming clustered into multiple components.* Below α0 the geometrical structure is not clustered and can be described by the simplest version of the replica method, known as replica symmetric solution. Above α0 the analytical computation of the typical volume requires a more sophisticated analysis that properly describes a clustered geometrical structure. This analysis can be performed by a variational technique which is known in statistical physics as the replica-symmetry-breaking (RSB) scheme, and the clustered geometrical structure of the solution space is known as RSB phase.

The capacity of the network, above which perfect classification becomes impossible, is found to be αc≃ 3.02. In the limit of large K (but still with Ñ≫1), the clustering transition occurs at a finite number of patterns per weight, α0≃2.95 (15), whereas the critical capacity grows with K as αc∝ln⁡K (20).

The Existence of WFM.

In order to detect the existence of WFM we use a large deviation measure which is the continuous version of the measure used in the discrete case: Each configuration of the weights is reweighted by a local volume term, analogously to the analysis in Local Entropy Theory. For the continuous case, however, we adopt a slightly different formalism which simplifies the analysis. Instead of constraining the set of y real replicas^† to be at distance D from a reference weight vector, we can identify the same WFM regions by constraining them directly to be at a given mutual overlap: For a given value q1∈−1,1 , we impose that Wa⋅Wb=Nq1 for all pairs of distinct replicas a,b. The overlap q1 is bijectively related to the mutual distance among replicas (which tends to 0 as q1→1). That, in turn, determines the distance between each replica and the normalized barycenter of the group N∑aWa/∑aWa, which takes the role that the reference vector had in the previous treatment. Thus, the regime of small D corresponds to the regime of q1 close to 1, and apart from this reparametrization the interpretation of the WFM volumes is the same. As explained in Materials and Methods, the advantage of this technique is that it allows to use directly the first-step formalism of the RSB scheme (1-RSB). Similarly to the discrete case, the computation of the maximal WFM volumes leads to the following results: For K=3 and in the large y limit, we find‡ Vq1 = maxWaa=1y1Ny⟨ln⁡VWaa=1y,q1⟩ξ = GSq1 + αGEq1withGSq1 = 121+ln2π+ln1−q1GEq1 = ∫∏ℓ=13Dvℓmaxu1,u2,u3−∑ℓ=13uℓ22+ + lnH∼1H∼2+H∼1H∼3+H∼2H∼3−2H∼1H∼2H∼3,where H∼ℓ≡Hd01−q1uℓ+q01−q1vℓ, Hx ≡ ∫x∞Dv, Dv≡dv12πe−12v2, d0≡yq1−q0 and q0 satisfies a saddle point equation that needs to be solved numerically. Vq1 is the logarithm of the volume of the solutions, normalized by Ny, under the spherical constraints on the weights, and with the real replicas forced to be at a mutual overlap q1. In analogy with the discrete case, we still refer to Vq1 as to the local entropy. It is composed by the sum of 2 terms; the first one, GSq1, corresponds to the log-volume at α=0, where all configurations are solutions and only the geometric constraints are present. This is an upper bound for the local entropy. The second term, V1q1 ≡ αGEq1, is in general negative and it represents the log of the fraction of solutions at overlap q1 (among all configurations at that overlap), and we call it normalized local entropy. Following the interpretation given above, we expect (in analogy to the discrete case at small D; see Fig. 3) that in a WFM this fraction is close to 1 (i.e., V1 is close to 0) in an extended region where the distance between the replicas is small, that is, where q1 is close to 1; otherwise, WFMs do not exist in the model. In Fig. 1 we report the values of V1q1 vs. the overlap q1 for different values of α. Indeed, one may observe that the behavior is qualitatively similar to that of the binary perceptron: Besides the solutions and all of the related local minima and saddles predicted by the standard statistical physics analysis (15⇓–17, 20) there exist absolute minima that are flat at relatively large distances. Indeed, reaching such wide minima efficiently is nontrivial, and different algorithms can have drastically different behaviors, as we will discuss in detail in Numerical Studies.

• Download figure
• Open in new tab
• Download powerpoint

Fig. 1.

Normalized local entropy V1 vs. q1, that is, logarithm of the fraction of configurations of y real replicas at mutual overlap q1 in which all replicas have zero error loss LNE. The curves are for a tree-like committee machine with K=3 hidden units trained on αN random patterns, for various values of α, obtained from a replica calculation in the limit of large N and large number of real replicas y. When the curves approach 0 as q1→1 it means that nearly all configurations are minima, and thus that the replicas are collectively exploring a wide minimum (any q1<1 implies distances of O(N) between replicas). The analogous figure for the limiting case of a large number of hidden units K can be found in SI Appendix.

The case K=3 is still relatively close to the simple perceptron, although the geometrical structure of its minima is already dominated by nonconvex features for α>1.76. A case that is closer to more realistic ANNs is K≫1 (but still N≫K), which, luckily enough, is easier to study analytically. We findGSq1 = 121+ln2π+ln1−q1GEq1 =  ∫  Dv maxu−u22+lnHΔq1e1−q1eu+q0e1−q1ev,where Δq1e=q1e−q0e with q1e≡1−2πarccosq1, q0e≡1−2πarccosq0, and q0 is fixed by a saddle point equation. The numerical results are qualitatively similar to those for K=3: We observe that indeed WFM still exist for all finite values of α. The analogue of Fig. 1 for this case is reported in the SI Appendix.

The results of the above WFM computation may require small corrections due to RSB effects, which, however, are expected to be very tiny due to the compact nature of the space of solutions at small distances.

A more informative aspect is to study the volumes around the solutions found by different algorithms. This can be done by the BP method, similarly to the computation of the weight enumerator function in error correcting codes (21).

Not All Devices Are Appropriate: The Parity Machine Does Not Display HLE/WFM Regions.

The extent by which a given model exhibits the presence of WFM can vary (see, e.g., Fig. 1). A direct comparison of the local entropy curves on different models in general does not yet have a well-defined interpretation, although at least for similar architectures it can still be informative (22). On the other hand, the existence of WFM itself is a structural property. For neural networks, its origin relies on the threshold sum form of the nonlinearity characterizing the formal neurons. As a check of this claim, we can analyze a model that is in some sense complementary, namely the so-called parity machine. We take its network structure to be identical to the committee machine, except for the output unit, which performs the product of the K hidden units instead of taking a majority vote. While the outputs of the hidden units are still given by sign activations, Eq. 14, the overall output of the network reads σoutμ=∏ℓ=1Kτℓ. The volume of the weights that correctly classifies a set of patterns is still given by Eq. 15.

Parity machines are closely related to error-correcting codes based on parity checks. The geometrical structure of the absolute minima of the error loss function is known (20) to be composed of multiple regions, each in one to one correspondence with the internal representations of the patterns. For random patterns such regions are typically tiny and we expect the WFM to be absent. Indeed, the computation of the volume proceeds analogously to the previous case^§ , and it shows that in this case for any distance the volumes of the minima are always bounded away from the maximal possible volume, that is, the volume one would find for the same distance when no patterns are stored. The log-ratio of the 2 volumes is constant and equal to −α⁡log2. In other words, the minima never become flat, at any distance scale.

The Connection between Local Entropy and CE.

Given that dense regions of optimal solutions exist in nonconvex ANN, at least in the case of independent random patterns, it remains to be seen which role they play in current models. Starting with the case of binary weights, and then generalizing the result to more complex architectures and to continuous weights, we can show that the most widely used loss function, the so-called CE loss, focuses precisely on such rare regions (see ref. 23 for the case of stochastic weights).

For the sake of simplicity, we consider a binary classification task with one output unit. The CE cost function for each input pattern readsLCEW = ∑μ=1M fγσμN∑i=1Nwiξiμ,[16]where fγx = −x2+12γlog2⁡coshγx. The parameter γ allows one to control the degree of “robustness” of the training (Fig. 2, Inset). In standard machine learning practice γ is simply set to 1, but a global rescaling of the weights Wi can lead to a basically equivalent effect. That setting can thus be interpreted as leaving γ as implicit, letting its effective value, and hence the norm of the weights, to be determined by the initial conditions and the training algorithm. As we shall see, controlling γ explicitly along the learning process plays a crucial role in finding HLE/WFM regions.

• Download figure
• Open in new tab
• Download powerpoint

Fig. 2.

Mean error rate achieved when optimizing the CE loss in the binary single-layer network, as predicted by the replica analysis, at various values of γ (increasing from top to bottom). The figure also shows the points αc≈0.83 (up to which solutions exist) and αU≈0.76 (up to which nonisolated solutions exist). (Inset) Binary CE function fγ(x) for various values of γ (increasing from top to bottom). For low values of γ, the loss is nonzero even for small positive values of the input, and thus the minimization procedure tends to favor more robust solutions. For large values of γ the function tends to max(−x,0). The dotted line shows the corresponding NE function, which is just 1 in case of an error and 0 otherwise (cf. Eq. 1).

For the binary case, however, the norm is fixed and thus we must keep γ as an explicit parameter. Note that since limγ→∞fγx = max−x,0 the minima of LCE below αc at large γ are the solutions to the training problem, that is, they coincide with those of LNE.

We proceed by first showing that the minimizers of this loss correspond to near-zero errors for a wide range of values of α and then by showing that these minimizers are surrounded by an exponential number of zero error solutions.

In order to study the probability distribution of the minima of LCE in the large N limit, we need to compute its Gibbs distribution (in particular, the average of the log of the partition function; see Eq. 5) as it has been done for the error loss LNE. The procedure follows standard steps and it is detailed in SI Appendix. The method requires one to solve 2 coupled integral equations as functions of the control parameters α, β, and γ. In Fig. 2 we show the behavior of the fraction of errors vs. the loading α for various values of γ. Up to relatively large values of α the optimum of LCE corresponds to extremely small values of LNE, virtually equal to zero for any accessible size N.

Having established that by minimizing the CE one ends up in regions of perfect classification where the error loss function is essentially zero, it remains to be understood which type of configurations of weights are found. Does the CE converge to an isolated point-like solution in the weight space (such as the typical zero energy configurations of the error function)^¶ or does it converge to the rare regions of HLE?

In order to establish the geometrical nature of the typical minima of the CE loss, we need to compute the average value of EDW (which tells us how many zero energy configurations of the error loss function can be found within a given radius D from a given W; see Eq. 6) when W is sampled from the minima of LCE. This can be accomplished by a well-known analytical technique (24) which was developed for the study of the energy landscape in disordered physical systems. The computation is relatively involved, and here we report only the final outcome. For the dedicated reader, all of the details of the calculation, which relies on the replica method and includes a double analytic continuation, can be found in SI Appendix. As reported in Fig. 3, we find that the typical minima of the CE loss for small finite γ are indeed surrounded by an exponential number of zero error solutions. In other words, the CE focuses on HLE regions. The range of γ for which this happens is generally rather wide, but neither too-small nor too-large values work well (additional details on this point are provided in SI Appendix). On the other hand, this analysis does not fully capture algorithmic dynamic effects, and in practice using a procedure in which γ is initialized to a small value and gradually increased should be effective in most circumstances (4).

• Download figure
• Open in new tab
• Download powerpoint

Fig. 3.

Average local entropy around a typical minimum of the LCE loss for various values of α and γ. The gray upper curve, corresponding to α=0, is an upper bound since in that case all configurations are solutions. For α=0.4, the 2 curves with γ=1 and 2 nearly saturate the upper bound at small distances, revealing the presence of dense regions of solutions (HLE regions). There is a slight improvement for γ=2, but the curve at γ→∞ shows that the improvement cannot be monotonic: In that limit, the measure is dominated by isolated solutions. This is reflected by the gap at small D in which the entropy becomes negative, signifying the absence of solutions in that range. For α=0.6 we see that the curves are lower, as expected. We also see that for γ=1 there is a gap at small D, and that we need to get to γ=3 in order to find HLE regions.

As an algorithmic check we have verified that while a simulated annealing approach gets stuck at very high energies when trying to minimize the error loss function, the very same algorithm with the CE loss is indeed successful up to relatively high values of α, with just slightly worse performance compared to an analogous procedure based on local entropy (13). In other words, the CE loss on single-layer networks is a computationally cheap and reasonably good proxy for the LE loss. These analytical results extend straightforwardly to 2-layer networks with binary weights. The study of continuous weight models can be performed resorting to the BP method.