Year: 2018

Abstract:
While it has not yet been proven, empirical evidence suggests that model generalization is related to local properties of the optima which can be described via the Hessian. We connect model generalization with the local property of a solution under the PAC-Bayes paradigm. In particular, we prove that model generalization ability is related to the Hessian, the higher-order “smoothness” terms characterized by the Lipschitz constant of the Hessian, and the scales of the parameters. Guided by the proof, we propose a metric to score the generalization capability of the model, as well as an algorithm that optimizes the perturbed model accordingly.

heading: 1 INTRODUCTION
text: Deep models have proven to work well in applications such as computer vision (Krizhevsky et al., 2012) (He et al., 2014) (Karpathy et al., 2014), speech recognition (Mohamed et al., 2012) (Hinton et al., 2012), and natural language processing (Socher et al., 2013) (Graves, 2013) (McCann et al., 2018). Many deep models have millions of parameters, which is more than the number of training samples, but the models still generalize well (Huang et al., 2017).
On the other hand, classical learning theory suggests the model generalization capability is closely related to the “complexity” of the hypothesis space, usually measured in terms of number of parameters, Rademacher complexity or VC-dimension. This seems to be a contradiction to the empirical observations that over-parameterized models generalize well on the test data1. Indeed, even if the hypothesis space is complex, the final solution learned from a given training set may still be simple. This suggests the generalization capability of the model is also related to the property of the solution.
Keskar et al. (2017) and Chaudhari et al. (2017) empirically observe that the generalization ability of a model is related to the spectrum of the Hessian matrix∇2L(w∗) evaluated at the solution, and large eigenvalues of the ∇2L(w∗) often leads to poor model generalization. Also, (Keskar et al., 2017), (Chaudhari et al., 2017) and (Novak et al., 2018b) introduce several different metrics to measure the “sharpness” of the solution, and demonstrate the connection between the sharpness metric and the generalization empirically. Dinh et al. (2017) later points out that most of the Hessian-based sharpness measures are problematic and cannot be applied directly to explain generalization. In particular, they show that the geometry of the parameters in RELU-MLP can be modified drastically by re-parameterization.
Another line of work originates from Bayesian analysis. Mackay (1995) first introduced Taylor expansion to approximate the (log) posterior, and considered the second-order term, characterized by the Hessian of the loss function, as a way of evaluating the model simplicity, or “Occam factor”. Recently Smith & Le (2018) use this factor to penalize sharp minima, and determine the optimal batch size. Germain et al. (2016) connect the PAC-Bayes bound and the Bayesian marginal likelihood when the loss is (bounded) negative log-likelihood, which leads to an alternative perspective on Occam’s razor. (Langford & Caruana, 2001), and more recently, (Harvey et al., 2017) (Neyshabur et al., 2017) (Neyshabur et al., 2018) use PAC-Bayes bound to analyze the generalization behavior of the deep models.
Since the PAC-Bayes bound holds uniformly for all “posteriors”, it also holds for some particular “posterior”, for example, the solution parameter perturbed with noise. This provides a natural
1For example over-parameterized neural network can fit any function of sample size n, making the Rademacher complexity large, but empirically those neural networks generalizes well. (Zhang et al., 2016)
way to incorporate the local property of the solution into the generalization analysis. In particular, Neyshabur et al. (2017) suggests to use the difference between the perturbed loss and the empirical loss as the sharpness metric. Dziugaite & Roy (2017) tries to optimize the PAC-Bayes bound instead for a better model generalization. Still some fundamental questions remain unanswered. In particular we are interested in the following question:
How is model generalization related to local “smoothness” of a solution?
In this paper we try to answer the question from the PAC-Bayes perspective. Under mild assumptions on the Hessian of the loss function, we prove the generalization error of the model is related to this Hessian, the Lipschitz constant of the Hessian, the scales of the parameters, as well as the number of training samples. The analysis also gives rise to a new metric for generalization. Based on this, we can approximately select an optimal perturbation level to aid generalization which interestingly turns out to be related to Hessian as well. Inspired by this observation, we propose a perturbation based algorithm that makes use of the estimation of the Hessian to improve model generalization.

heading: 2 PAC-BAYES AND MODEL GENERALIZATION
text: We consider the supervised learning in PAC-Bayes scenario (McAllester, 2003) (McAllester, 1998) (McAllester, 1999) (Langford & Shawe-Taylor, 2002). Suppose we have a labeled data set S = {si = (xi, yi) | i ∈ {1, . . . , n}, xi ∈ Rd, yi ∈ {0, 1}k}, where (xi, yi) are sampled i.i.d. from a distribution xi, yi ∼ Ds.
2The variables from different layers are shared so that the model only has two free parameters w1 and w2. The bound in (a) is approximated with η = 39 using inequality (8)
The PAC-Bayes paradigm assumes probability measures over the function class F : X → Y . In particular, it assumes a “posterior” distribution Df as well as a “prior” distribution πf over the function class F. We are interested in minimizing the expected loss, in terms of both the random draw of samples as well as the random draw of functions:
L(Df ,Ds) = Ef∼DfEx,y∼Ds l(f, x, y). Correspondingly, the empirical loss in the PAC-Bayes paradigm is the expected loss over the draw of functions from the posterior:
L̂(Df ,S) = Ef∼Df 1
n n∑ i=1 l(f, xi, yi) (1)
PAC-Bayes theory suggests the gap between the expected loss and the empirical loss is bounded by a term that is related to the KL divergence between Df and πf (McAllester, 1999) (Langford & Shawe-Taylor, 2002). In particular, if the function f is parameterized as f(w) with w ∈ W , when Dw is perturbed around any w, we have the following PAC-Bayes bound (Seldin et al., 2012a) (Seldin et al., 2012b) (Neyshabur et al., 2017) (Neyshabur et al., 2018): Theorem 1 (PAC-Bayes-Hoeffding Perturbation). Let l(f, x, y) ∈ [0, 1], and π be any fixed distribution over the parametersW . For any δ > 0 and η > 0, with probability at least 1 − δ over the draw of n samples, for any w and any random perturbation u,
Eu[L(w + u)] ≤ Eu[L̂(w + u)] + KL(w + u||π) + log 1δ
η +
η
2n (2)
One may further optimize η to get a bound that scales approximately as Eu[L(w+u)] . Eu[L̂(w+
u)] + 2
√ KL(w+u||π)+log 1δ
2n (Seldin et al., 2012b). 3 A nice property of the perturbation bound (2) is
it connects the generalization with the local properties around the solution w through some perturbation u around w. In particular, suppose L̂(w∗) is a local optimum, when the perturbation level of u is small, Eu[L̂(w∗ + u)] tends to be small, but KL(w∗ + u||π) may be large since the posterior is too “focused” on a small neighboring area around w∗, and vice versa. As a consequence, we may need to search for an “optimal” perturbation level for u so that the bound is minimized.

heading: 3 MAIN RESULT
text: While some researchers have already discovered empirically the generalization ability of the models is related to the second order information around the local optima, to the best of our knowledge there is no work on how to connect the Hessian matrix∇2L̂(w) with the model generalization rigorously. In this section we introduce the local smoothness assumption, as well as our main theorem.
It may be unrealistic to assume global smoothness properties for the deep models. Usually the assumptions only hold in a small local neighborhood Neigh(w∗) around a reference point w∗. In this paper we define the neighborhood set as
Neighκ(w ∗) = {w | |wi − w∗i | ≤ κi ∀i}
where κi ∈ R+ is the “radius” of the i-th coordinate. In our draft we focus on a particular type of radius κi(w∗) = γ|w∗i |+ , but our argument holds for other types of radius, too. In order to get a control of the deviation of the optimal solution we need to assume inNeighγ, (w∗), the empirical loss function L̂ in (1) is Hessian Lipschitz, which is defined as: Definition 1 (Hessian Lipschitz). A twice differentiable function f(·) is ρ-Hessian Lipschitz if:
∀w1, w2, ‖∇2f(w1)−∇2f(w2)‖ ≤ ρ‖w1 − w2‖, (3)
where ‖ · ‖ is the operator norm. The Hessian Lipschitz condition has been used in the numeric optimization community to model the smoothness of the second-order gradients (Nesterov & Polyak, 2006) (Carmon et al., 2018) (Jin et al., 2018). In the rest of the draft we always assume the following:
3Since η cannot depend on the data, one has to build a grid and use the union bound.
Assumption 1. In Neighκ(w∗) the empirical loss L̂(w) defined in (1) is convex, and ρ-Hessian Lipschitz.
For the uniform perturbation, the following theorem holds: Theorem 2. Suppose the loss function l(f, x, y) ∈ [0, 1], and model weights are bounded |wi| + κi(w) ≤ τi ∀i. With probability at least 1 − δ over the draw of n samples, for any w̌ ∈ Rm such that assumption 1 holds
Eu[L(w̌ + u)] ≤ L̂(w̌) +O √m+∑i log τiσ̌i + log 1δ n  where ui ∼ U(−σ̌i, σ̌i) are i.i.d. uniformly distributed random variables, and
σ̌i(w̌, η, γ) = min
(√ 1
√ mn(∇2i,iL̂(w̌)/3 + ρm1/2κi(w̌)/9) , κi(w̌)
) (4)
Theorem 2 says if we choose the perturbation levels carefully, the expected loss of a uniformly perturbed model is controlled. The bound is related to the diagonal element of Hessian (logarithmic), the Lipschitz constant ρ of the Hessian (logarithmic), the neighborhood scales characterized by κ (logarithmic), the number of parameters m, and the number of samples n. Also roughly the
perturbation level is inversely related to √ ∇2i,iL̂, suggesting the model be perturbed more along the coordinates that are “flat”.4
Similar argument can be made on the truncated Gaussian perturbation, which is presented in Appendix B. In the next section we walk through some intuitions of our arguments.

heading: 4 CONNECTING GENERALIZATION AND HESSIAN
text: Suppose the empirical loss function L̂(w) satisfies the local Hessian Lipschitz condition, then by Lemma 1 in (Nesterov & Polyak, 2006), the perturbation of the function around a fixed point can be bounded by terms up to the third-order,
L̂(w + u) ≤ L̂(w) +∇L̂(w)Tu+ 1 2 uT∇2L̂(w)u+ 1 6 ρ‖u‖3 for w + u ∈ Neighκ(w) (5)
For perturbations with zero expectation, i.e., E[u] = 0, the linear term in (5), Eu[∇L̂(w)Tu] = 0. Because the perturbation ui for different parameters are independent, the second order term can also be simplified, since Eu[ 12u T∇2L̂(w)u] = 12 ∑ i∇2i,iL̂(w)E[u2i ].
Considering (2),(5) and assumption 1, it is straightforward to see the bound below holds with probability at least 1− δ
Eu[L(w∗ + u)] ≤ L̂(w∗) + 1
2 ∑ i ∇2i,iL̂(w∗)E[u2i ] + ρ 6 E[‖u‖3] + KL(w∗ + u||π) + log 1δ η + η 2n
(6)
Suppose ui ∼ U(−σi, σi), and σi ≤ κi(w) ∀i. That is, the “posterior” distribution of the model parameters are uniform distribution, and the distribution supports vary for different parameters. We also assume the perturbed parameters are bounded, i.e., |wi| + κi(w) ≤ τi ∀i.5 If we choose the prior π to be ui ∼ U(−τi, τi), and then KL(w + u||π) = ∑ i log(τi/σi).
The third order term in (6) is bounded by
ρ 6 E[‖u‖3] ≤ ρm
1/2
6 E[‖u‖33] ≤
ρm1/2
6 ∑ i κi(w)E[u2i ] = ρm1/2 18 ∑ i κi(w)σ 2 i ,
4Unfortunately the bound in theorem 2 does not explain the over-parameterization phenomenon since when m n the right hand side explodes.
5One may also assume the same τ for all parameters for a simpler argument. The proof procedure goes through in a similar way.
where we use the inequality ‖u‖2 ≤ m 1 6 ‖u‖3 and m is the number of parameters. Pluging in (6), we get
Eu[L(w + u)] ≤ L̂(w) + 1
6 ∑ i ∇2i,iL(w)σ2i + ρm1/2 18 ∑ i κi(w)σ 2 i +
∑ i log τi σi
+ log 1δ η + η 2n
(7)
Solve for σ that minimizes the right hand side, and we have the following lemma: Lemma 3. Suppose the loss function l(f, x, y) ∈ [0, 1], and model weights are bounded |wi| + κi(w) ≤ τi ∀i. Given any δ > 0 and η > 0, with probability at least 1 − δ over the draw of n samples, for any w∗ ∈ Rm such that assumption 1 holds,
Eu[L(w∗ + u)] ≤ L̂(w∗) + m/2 +
∑ i log
τi σ∗i + log 1δ
η +
η
2n (8)
where ui ∼ U(−σ∗i , σ∗i ) are i.i.d. uniformly perturbed random variables, and
σ∗i (w ∗, η, γ) = min
(√ 1
η(∇2i,iL(w∗)/3 + ρm1/2κi(w∗)/9)) , κi(w
∗) ) . (9)
In our experiment, we simply treat η as a hyper-parameter. Other other hand, one may further build a weighted grid over η and optimize for the best η (Seldin et al., 2012b). That leads to Theorem 2. Details of the proof are presented in the Appendix C and D.

heading: 5 ON THE RE-PARAMETERIZATION OF RELU-MLP
text: Dinh et al. (2017) points out the spectrum of∇2L̂ itself is not enough to determine the generalization power. In particular, for a multi-layer perceptron with RELU as the activation function, one may re-parameterize the model and scale the Hessian spectrum arbitrarily without affecting the model prediction and generalization when cross entropy (negative log likelihood) is used as the loss and w∗ is the “true” parameter of the sample distribution.
In general our bound does not assume the loss to be the cross entropy. Also we do not assume the model is RELU-MLP. As a result we would not expect our bound stays exactly the same during the re-parameterization. On the other hand, the optimal perturbation levels in our bound scales inversely when the parameters scale, so the bound only changes approximately with a speed of logarithmic factor. According to Lemma (3), if we use the optimal σ∗ on the right hand side of the bound, ∇2L̂(w), ρ, and w∗ are all behind the logarithmic function. As a consequence, for RELU-MLP, if we do the re-parameterization trick, the change of the bound is small.
In the next two sections we introduce some heuristic-based approximations enlightened by the bound, as well as some interesting empirical observations.

heading: 6 AN APPROXIMATE GENERALIZATION METRIC
text: Assuming L̂(w) is locally convex around w∗, so that ∇2i,iL̂(w∗) ≥ 0 for all i. If we look at Lemma 3, for fixed m and n, the only relevant term is ∑ i log
τi σ∗i . Replacing the optimal σ∗, and using |wi| + κi(w) to approximate τi, we come up with PAC-Bayes based Generalization metric, called pacGen,6
Ψκ(L̂, w ∗) = ∑ i log ( (|w∗i |+ κi(w∗)) max (√ ∇2i,iL̂(w∗) + ρ(w∗) √ mκi(w∗), 1 κi(w∗) )) .
6Even though we assume the local convexity in our metric, in application we may calculate the metric on every points. When∇2i,iL̂(w∗) + ρ(w∗) √ mκi(w ∗) < 0 we simply treat it as 0.
A self-explained toy example is displayed in Figure 1. To calculate the metric on real-world data we need to estimate the diagonal elements of the Hessian∇2L̂ as well as the Lipschitz constant ρ of the Hessian. For efficiency concern we follow Adam (Kingma & Ba, 2014) and approximate ∇2i,iL̂ by (∇L̂[i])2. Also we use the exponential smoothing technique with β = 0.999 as in (Kingma & Ba, 2014).
To estimate ρ, we first estimate the Hessian of a randomly perturbed model ∇2L̂(w + u), and then approximate ρ by ρ = maxi |∇2iL(w+ui)−∇ 2 iL(w)|
|ui| . For the neighborhood radius κ we use γ = 0.1 and = 0.1 for all the experiments in this section.
We used the same model without dropout from the PyTorch example 7. Fixing the learning rate as 0.1, we vary the batch size for training. The gap between the test loss and the training loss, and the metric Ψκ(L̂, w∗) are plotted in Figure 2. We had the same observation as in (Keskar et al., 2017) that as the batch size grows, the gap between the test loss and the training loss tends to get larger. Our proposed metric Ψκ(L̂, w∗) also shows the exact same trend. Note we do not use LR annealing heuristics as in (Goyal et al., 2017) which enables large batch training.
Similarly we also carry out experiment by fixing the training batch size as 256, and varying the learning rate. Figure 4 shows generalization gap and Ψκ(L̂, w∗) as a function of epochs. It is observed that as the learning rate decreases, the gap between the test loss and the training loss increases. And the proposed metric Ψκ(L̂, w∗) shows similar trend compared to the actual generalization gap. Similar trends can be observed if we run the same model on CIFAR-10 (Krizhevsky, 2009) as shown in Figure 3 and Figure 5.

heading: 7 A PERTURBED OPTIMIZATION ALGORITHM
text: Adding noise to the model for better generalization has proven successful both empirically and theoretically (Zhu et al., 2018) (Hoffer et al., 2017) (Jastrzȩbski et al., 2017) (Dziugaite & Roy, 2017) (Novak et al., 2018a). Instead of only minimizing the empirical loss, (Langford & Caruana,
7https://github.com/pytorch/examples/tree/master/mnist
2001) and (Dziugaite & Roy, 2017) assume different perturbation levels on different parameters, and minimize the generalization bound led by PAC-Bayes for better model generalization. However how to integrate the smoothness property of the local optima is not clear.
The right hand side of (2) has Eu[L̂(w + u)]. This suggests rather than minimizing the empirical loss L̂(w), we should optimize the perturbed empirical loss Eu[L̂(w+u)] instead for a better model generalization power.
We introduce a systematic way to perturb the model weights based on the PAC-Bayes bound. Again we use the same exponential smoothing technique as in Adam (Kingma & Ba, 2014) to estimate the Hessian ∇2L̂. The details of the algorithm is presented in Algorithm 1, where we treat η as a hyper-parameter.
Even though in theoretical analysis Eu[∇L̂ ·u] = 0, in applications,∇L̂ ·u won’t be zero especially when we only implement 1 trial of perturbation. On the other hand, if the gradient ∇L̂ is close to zero, then the first order term can be ignored. As a consequence, in Algorithm 1 we only perturb the parameters that have small gradients whose absolute value is below β2. For efficiency issues we used a per-parameter ρi capturing the variation of the diagonal element of Hessian. Also we decrease the perturbation level with a log factor as the epoch increases.
We compare the perturbed algorithm against the original optimization method on CIFAR-10, CIFAR-100 (Krizhevsky, 2009), and Tiny ImageNet8. The results are shown in Figure 6. We use the Wide-ResNet (Zagoruyko & Komodakis, 2018) as the prediction model.9 The depth of the chosen model is 58, and the widen-factor is set as 3. The dropout layers are turned off. For CIFAR-10 and CIFAR-100, we use Adam with a learning rate of 10−4, and the batch size is 128. For the perturbation parameters we use η = 0.01, γ = 10, and =1e-5. For Tiny ImageNet, we use SGD with learning rate 10−2, and the batch size is 200. For the perturbed SGD we set η = 100, γ = 1,
8https://tiny-imagenet.herokuapp.com/ 9https://github.com/meliketoy/wide-resnet.pytorch/blob/master/networks/
wide_resnet.py
Algorithm 1 Perturbed OPT Require: η, γ = 0.1, β1 = 0.999, β2 = 0.1, =1e-5.
1: Initialization: σi ← 0 for all i. t← 0, h0 ← 0 2: for epoch in 1, . . . , N do 3: for minibatch in one epoch do 4: for all i do 5: if t > 0 then 6: ρ[i]← |ht[i]−ht−1[i]|‖wt−wt−1‖ 7: κ[i]← γlog(1+epoch) |wt−1[i]|+
8: σi ← min (
1 log(1+epoch) √ η(ht[i]+ρ[i]·κ[i]) , κ[i]
) · 1|gt[i]|<β2
9: ut[i] ∼ U(−σi, σi)(sample a set of perturbations) 10: gt+1 ← ∇wL̂t(wt + ut) (get stochastic gradients w.r.t. perturbed loss) 11: ht+1 ← β1ht + (1− β1)g2t+1 (update second moment estimate) 12: wt+1 ← OPT(wt) (update w using off-the-shell algorithms) 13: t← t+ 1
and =1e-5. Also we use the validation set as the test set for the Tiny ImageNet. We observe the effect with perturbation appears similar to regularization. With the perturbation, the accuracy on the training set tends to decrease, but the test on the validation set increases. The perturbedOPT also works better than dropout possibly due to the fact that the it puts different levels of perturbation on different parameters according to the local smoothness structures, while only one dropout rate is set for the all the parameters across the model for the dropout method.

heading: 8 CONCLUSION
text: We connect the smoothness of the solution with the model generalization in the PAC-Bayes framework. We prove that the generalization power of a model is related to the Hessian and the smoothness of the solution, the scales of the parameters, as well as the number of training samples. In particular, we prove that the best perturbation level scales roughly as the inverse of the square root of the Hessian, which mostly cancels out scaling effect in the re-parameterization suggested by (Dinh et al., 2017). To the best of our knowledge, this is the first work that integrate Hessian in the model generalization bound rigorously. It also roughly explains the effect of re-parameterization over the generalization. Based on our generalization bound, we propose a new metric to test the model generalization and a new perturbation algorithm that adjusts the perturbation levels according to the Hessian. Finally, we empirically demonstrate the effect of our algorithm is similar to a regularizer in its ability to attain better performance on unseen data.

heading: A DETAILS OF THE TOY EXAMPLE
text: This section discusses the details of the toy example shown in Figure (1). We construct a small 2- dimensional sample set from a mixture of 3 Gaussians, and then binarize the labels by thresholding them from the median value. The sample distribution is shown in Figure 1b. For the model we use a 5-layer MLP with sigmoid as the activation and cross entropy as the loss. There are no bias terms in the linear layers, and the weights are shared. For the shared 2-by-2 linear coefficient matrix, we treat two entries as constants and optimize the other 2 entries. In this way the whole model has only two free parameters w1 and w2.
The model is trained using 100 samples. Fixing the samples, we plot the loss function with respect to the model variables L̂(w1, w2), as shown in Figure 1a. Many local optima are observed even in this simple two-dimensional toy example. In particular: a sharp one, marked by the vertical green line, and a flat one, marked by the vertical red line. The colors on the loss surface display the values of the generalization metric scores (pacGen) defined in Section 6. Smaller metric value indicates better generalization power.
As displayed in the figure, the metric score around the global optimum, indicated by the vertical green bar, is high, suggesting possible poor generalization capability as compared to the local optimum indicated by the red bar. We also plot a plane on the bottom of the figure. The color projected on the bottom plane indicates an approximated generalization bound, which considers both the loss and the generalization metric.10 The local optimum indicated by the red bar, though has a slightly higher loss, has a similar overall bound compared to the “sharp” global optimum.
On the other hand, fixing the parameter w1 and w2, we may also plot the labels predicted by the model given the samples. Here we plot the prediction from both the sharp minimum (Figure 1c) and the flat minimum (Figure 1d). The sharp minimum, even though it approximates the true label better, has some complex structures in its predicted labels, while the flat minimum seems to produce a simpler classification boundary.

heading: B TRUNCATED GAUSSIAN
text: Because the Gaussian distribution is not bounded but the inequality (5) requires bounded perturbation, we first truncate the distribution. The procedure of truncation is similar to the proof in (Neyshabur et al., 2018) and (McAllester, 2003).
Let ui ∼ N(0, σ2i ). Denote the truncated Gaussian as Nκi(0, σ2i ). If ũi ∼ Nκi(0, σ2i ) then
Pκi(ũ) = 1
Zi
{ p(ui) if |ui| < κi(w)
0 o.w. (10)
Now let’s look at the event
E = {u | |ui| < κi(w) ∀ i} (11)
If ∀i σi < κi(w)√2erf−1( 12m ) , by union bound P(E) ≥ 1/2. Here erf−1 is the inverse Gaussian error function defined as erf(x) = 2√ π ∫ x 0 e−t 2
dt, andm is the number of parameters. Following a similar procedure as in the proof of Lemma 1 in (Neyshabur et al., 2018),
KL(w + ũ||π) ≤ 2(KL(w + u||π) + 1) (12)
Suppose the coefficients are bounded such that ∑ i w 2 i ≤ τ , where τ is a constant. Choose the prior π as N(0, τI), and we have
KL(w + u||π) ≤ 1 2 (m log τ − ∑ i log σ2i −m+ 1 τ ∑ i σ2i + 1) (13)
10the bound was approximated with η = 39 using inequality (8)
Notice that after the truncation the variance only becomes smaller, so the bound of (6) for the truncated Gaussian becomes
Eu[L(w + ũ)] ≤ L̂(w)+ 1
2 ∑ i ∇2i,iL(w)σ2i + ρm1/2 6 ∑ i κi(w)σ 2 i
+ m log τ −
∑ i log σ 2 i −m+ 1τ ∑ i σ 2 i + 1 + 2 log 1 δ
2η +
η
2n (14)
Again when L̂(w) is convex around w∗ such that ∇2L̂(w∗) ≥ 0, solve for the best σi and we get the following lemma: Lemma 4. Suppose the loss function l(f, x, y) ∈ [0, 1], and model weights are bounded ∑ i w 2 i ≤ τ . For any δ > 0 and η, with probability at least 1 − δ over the draw of n samples, for any w∗ ∈ Rm such that assumption 1 holds,
Eu[L(w∗ + ũ)] ≤ L̂(w∗) + m log τ −
∑ i log σ 2 i + 1 + 2 log 1 δ
2η +
η
2n (15)
where ũi ∼ Nκi(0, σ∗i ) are i.i.d. random variables distributed as truncated Gaussian,
σ∗i = min
(√ 1
η∇2i,iL̂(w∗) + ρηm1/2 3 κi(w ∗) + 1τ
, κi(w ∗)√ 2erf−1( 12m )
) (16)
and σ∗2i is the i-th diagonal element in Σ ∗.
Again we have an extra term η, which may be further optimized over a grid to get a tighter bound. In our algorithm we treat η as a hyper-parameter instead.

heading: C PROOF OF LEMMA 3
text: Proof. We rewrite the inequality (7) below
Eu[L(w + u)] ≤ L̂(w) + 1
6 ∑ i ∇2i,iL(w)σ2i + ρm1/2 18 ∑ i κi(w)σ 2 i +
∑ i log τi σi
+ log 1δ η + η 2n
(17)
The terms related to σi on the right hand side of (17) are
1 6 ∇2i,iL(w)σ2i +
ρm1/2
18 κi(w)σ
2 i − log σi η
(18)
Since the assumption is ∇2i,iL̂(w∗) ≥ 0 for all i, ∇2i,iL̂(w) + ρm1/2κi(w)/3 > 0. Solving for σ that minimizes the right hand side of (17), and we have
σ∗i (w, η, γ) = min
(√ 1
η(∇2i,iL̂(w)/3 + ρm1/2κi(w)/9) , κi(w)
) (19)
The term 16 ∑ i∇2i,iL(w)σ2i + ρm1/2 18 ∑ i κi(w)σ 2 i on the right hand side of (7) is monotonically increasing w.r.t. σ2, so
1
6 ∑ i ∇2i,iL(w)σ∗2i + ρm1/2 18 ∑ i κi(w)σ ∗2 i
≤ ∑ i ( 1 6 ∇2i,iL(w) + ρm1/2 18 κi(w) ) 1 η(∇2i,iL̂(w)/3 + ρm1/2κi(w)/9) = m
2η (20)
Combine the inequality (20), and the equation (19) with (17), and we complete the proof.

heading: D PROOF OF THEOREM 2
text: Proof. Combining (4) and (7), we get
Eu[L(w̌ + u)] ≤ L̂(w̌) + 1
2
√ m
n +
∑ i log τi σ̌i
+ log 1δ η + η 2n
The following proof is similar to the proof of Theorem 6 in (Seldin et al., 2012b). Note the η in Lemma (3) cannot depend on the data. In order to optimize η we need to build a grid of the form
ηj = e j √ 2n log 1
δj
for j ≥ 0. For a given value of ∑ i log τi σ̌i , we pick ηj , such that
j =
⌊ 1
2 log
(∑ i log τi σ̌i
log 1δj + 1 )⌋ where bxc is the largest integer value smaller than x. Set δj = δ2−(j+1), and take a weighted union bound over ηj-s with weights 2−(j+1), and we have with probability at least 1− δ,
Eu[L(w̌ + u)] ≤ L̂(w̌) + 1
2
√ m
n + (1 + 1/e)
√√√√√∑i log τiσ̌i + log 1δ + log 22 ( 2 + log (∑ i log τi σ̌i log 1δ + 1 )) 2n
Simplify the right hand side and we complete the proof.

heading: E PROOF OF LEMMA 4
text: Proof. We first rewrite the inequality (14) below:
Eu[L(w∗ + ũ)] ≤ L̂(w∗)+ 1
2 ∑ i ∇2i,iL(w∗)σ2i + ρm1/2 6 ∑ i κi(w ∗)σ2i
+ m log τ −
∑ i log σ 2 i −m+ 1τ ∑ i σ 2 i + 1 + 2 log 1 δ
2η +
η
2n
The terms related to σi on the right hand side of (14) is( 1
2 ∇2i,iL(w∗) +
ρm1/2
6 κi(w
∗) + 1
2τη
) σ2i −
log σ2i 2η
(21)
Take gradients w.r.t. σi, when ∇2i L̂ ≥ 0, we get the optimal σ∗i ,
σ∗i = min
(√ 1
η∇2i,iL̂(w∗) + ρηm1/2 3 κi(w ∗) + 1τ
, κi(w ∗)√ 2erf−1( 12m )
)
Note the first term in (21) is monotonously increasing w.r.t. σi, so( 1
2 ∇2i,iL(w∗) +
ρm1/2
6 κi(w
∗) + 1
2τη
) σ∗2i
≤ ( 1
2 ∇2iL(w∗) +
ρm1/2
6 κi(w
∗) + 1
2τη
) 1
η∇2i L̂(w∗) + ρηm1/2 3 κi(w ∗) + 1τ
= 1
2η (22)
Summing over m parameters and combine (14), we complete the proof.

heading: F A LEMMA ABOUT EIGENVALUES OF HESSIAN AND GENERALIZATION
text: By extrema of the Rayleigh quotient, the quadratic term on the right hand side of inequality (5) is further bounded by
uT∇2L̂(w)u ≤ λmax(∇2L̂(w))‖u‖2. (23)
This is consistent with the empirical observations of Keskar et al. (2017) that the generalization ability of the model is related to the eigenvalues of ∇2L̂(w). The inequality (23) still holds even if the perturbations ui and uj are correlated. We add another lemma about correlated perturbations below.
Lemma 5. Suppose the loss function l(f, x, y) ∈ [0, 1]. Let π be any distribution on the parameters that is independent from the data. Given δ > 0 η > 0, with probability at least 1− δ over the draw of n samples, for any local optimal w∗ such that ∇L̂(w∗) = 0, L̂(w) satisfies the local ρ-Hessian Lipschitz condition in Neighκ(w∗), and any random perturbation u, s.t., |ui| ≤ κi(w∗) ∀i, we have
Eu[L(w∗ + u)] ≤ L̂(w∗) + 1
2 λmax
( ∇2L̂(w∗) )∑ i E[u2i ] + ρ 6 E[‖u‖3]
+ KL(w∗ + u||π) + log 1δ
η +
η
2n . (24)
Proof. The proof of the Lemma 5 is straightforward. Since∇L̂(w∗) = 0, the first order term is zero at the local optimal point even if E[u] 6= 0. By extrema of the Rayleigh quotient, the quadratic term on the right hand side of inequality (5) is further bounded by
uT∇2L̂(w)u ≤ λmax ( ∇2L̂(w) ) ‖u‖2. (25)
Due to the linearity of the expected value, E[uT∇2L̂(w)u] ≤ λmax ( ∇2L̂(w) )∑ i E[u2i ], (26)
which does not assume independence among the perturbations ui and uj for i 6= j.

heading: G PERTOPT V.S. DROPOUT: AN EMPIRICAL COMPARISON
text: This section contains several figures comparing dropout and the proposed perturbation algorithm. Dropout can be viewed as multiplicative perturbation using Bernoulli distribution. It has already been widely used in almost every deep models. For comparison we present results using the exact same wide resnet architectures except the dropout layers are turned on or off. We report the accuracy with dropout rate of 0.0, 0.1, 0.3, and 0.5 on CIFAR-10 and CIFAR-100. For Tiny ImageNet we report the result with dropout rate being 0.0, 0.1, and 0.3. Again for the pertOPT algorithm all the dropout layers are turned off.
The depth of the chosen wide resnet model (Zagoruyko & Komodakis, 2018) is 58, and the widenfactor is set as 3. For CIFAR-10 and CIFAR-100, we use Adam with a learning rate of 10−4, and the batch size is 128. For the perturbation parameters we use η = 0.01, γ = 10, and =1e-5. For Tiny ImageNet, we use SGD with learning rate 10−2, and the batch size is 200. For the perturbed SGD we set η = 100, γ = 1, and =1e-5. Also we use the validation set as the test set for the Tiny ImageNet.
Figure (7), (8), and (9) show the accuracy versus epochs for training and validation in CIFAR10, CIFAR-100, and Tiny ImageNet respectively. It is pretty clear that with added dropout the validation/test accuracy got boosted compared to the original method. For CIFAR-10, dropout rate 0.3 seems to work best compared to all the other dropout configurations. For CIFAR-100 and Tiny
ImageNet, dropout 0.1 seems to work better. This may be due to the fact that CIFAR-10 has less training samples so more regularization is needed to prevent overfit.
Although both perturbedOPT and dropout can be viewed as certain kind of regularization, in all experiments the perturbed algorithm shows better performance on the validation/test data sets compared to the dropout methods. One possible explanation is maybe the perturbed algorithm puts different levels of perturbation on different parameters according to the local smoothness structures, while only one dropout rate is set for all the parameters across the model.

Review summaries:
The authors study generalization capabilities of neural networks local minimums thanks to a PAC-Bayesian analysis that grasps the local smoothness properties. Even if some assumptions are made along the way, their analysis provides a metric that gives insight on the accuracy of a solution, as well as an optimization algorithm. Both of these result show good empirical behavior.