On the Local Minima of the Empirical Risk

Introduction

The optimization of nonconvex loss functions has been key to the success of modern machine learning. While classical research in optimization focused on convex functions having a unique critical point that is both locally and globally minimal, a nonconvex function can have many local maxima, local minima and saddle points, all of which pose significant challenges for optimization. A recent line of research has yielded significant progress on one aspect of this problem—it has been established that favorable rates of convergence can be obtained even in the presence of saddle points, using simple variants of stochastic gradient descent (e.g., Ge et al., 2015; Carmon et al., 2016; Agarwal et al., 2017; Jin et al., 2017a). These research results have introduced new analysis tools for nonconvex optimization, and it is of significant interest to begin to use these tools to attack the problems associated with undesirable local minima.

where the error $\nu$ typically decreases with the number of samples $n$. See Figure 1(a) for a depiction of this result, and Figure 1(b) for an illustration of the effect of sampling on the optimization landscape. We wish to exploit this nearness of $F$ and $f$ to design and analyze optimization procedures that find approximate local minima (see Definition 1) of the smooth function $F$, while avoiding the local minima that exist only in the sampled function $f$.

Although the relationship between population risk and empirical risk is our major motivation, we note that other applications of our framework include two-stage robust optimization and private learning (see Section 5.2). In these settings, the error $\nu$ can be viewed as the amount of adversarial perturbation or noise due to sources other than data sampling. As in the sampling setting, we hope to show that simple algorithms such as stochastic gradient descent are able to escape the local minima that arise as a function of $\nu$.

Much of the previous work on this problem studies relatively small values of $\nu$, leading to “shallow” local minima, and applies relatively large amounts of noise, through algorithms such as simulated annealing (Belloni et al., 2015) and stochastic gradient Langevin dynamics (SGLD) (Zhang et al., 2017). While such “large-noise algorithms” may be justified if the goal is to approach a stationary distribution, it is not clear that such large levels of noise is necessary in the optimization setting in order to escape shallow local minima. The best existing result for the setting of nonconvex $F$ requires the error $\nu$ to be smaller than $O(\epsilon^{2}/d^{8})$, where $\epsilon$ is the precision of the optimization guarantee (see Definition 1) and $d$ is the problem dimension (Zhang et al., 2017) (see Figure 2). A fundamental question is whether algorithms exist that can tolerate a larger value of $\nu$, which would imply that they can escape “deeper” local minima. In the context of empirical risk minimization, such a result would allow fewer samples to be taken while still providing a strong guarantee on avoiding local minima.

We thus focus on the two central questions: (1) Can a simple, optimization-based algorithm avoid shallow local minima despite the lack of “large noise”? (2) Can we tolerate larger error $\nu$ in the optimization setting, thus escaping “deeper” local minima? What is the largest error that the best algorithm can tolerate?

In this paper, we answer both questions in the affirmative, establishing optimal dependencies between the error $\nu$ and the precision of a solution $\epsilon$. We propose a simple algorithm based on SGD (Algorithm 1) that is guaranteed to find an approximate local minimum of $F$ efficiently if $\nu\leq O(\epsilon^{1.5}/d)$, thus escaping all saddle points of $F$ and all additional local minima introduced by $f$. Moreover, we provide a matching lower bound (up to logarithmic factors) for all algorithms making a polynomial number of queries of $f$. The lower bound shows that our algorithm achieves the optimal tradeoff between $\nu$ and $\epsilon$, as well as the optimal dependence on dimension $d$. We also consider the information-theoretic limit for identifying an approximate local minimum of $F$ regardless of the number of queries. We give a sharp information-theoretic threshold: $\nu=\Theta(\epsilon^{1.5})$ (see Figure 2).

As a concrete example of the application to minimizing population risk, we show that our results can be directly used to give sample complexities for learning a ReLU unit, whose empirical risk is nonsmooth while the population risk is smooth almost everywhere.

A number of other papers have examined the problem of optimizing a target function $F$ given only function evaluations of a function $f$ that is pointwise close to $F$. Belloni et al. (2015) proposed an algorithm based on simulated annealing. The work of Risteski and Li (2016) and Singer and Vondrak (2015) discussed lower bounds, though only for the setting in which the target function $F$ is convex. For nonconvex target functions $F$, Zhang et al. (2017) studied the problem of finding approximate local minima of $F$, and proposed an algorithm based on Stochastic Gradient Langevin Dynamics (SGLD) (Welling and Teh, 2011), with maximum tolerance for function error $\nu$ scaling as $O(\epsilon^{2}/d^{8})$The difference between the scaling for $\nu$ asserted here and the $\nu=O(\epsilon^{2})$ claimed in (Zhang et al., 2017) is due to difference in assumptions. In our paper we assume that the Hessian is Lipschitz with respect to the standard spectral norm; Zhang et al. (2017) make such an assumption with respect to nuclear norm.. Other than difference in algorithm style and $\nu$ tolerance as shown in Figure 2, we also note that we do not require regularity assumptions on top of smoothness, which are inherently required by the MCMC algorithm proposed in Zhang et al. (2017). Finally, we note that in parallel, Kleinberg et al. (2018) solved a similar problem using SGD under the assumption that $F$ is one-point convex.

Previous work has also studied the relation between the landscape of empirical risks and the landscape of population risks for nonconvex functions. Mei et al. (2016) examined a special case where the individual loss functions $L$ are also smooth, which under some assumptions implies uniform convergence of the gradient and Hessian of the empirical risk to their population versions. Loh and Wainwright (2013) showed for a restricted class of nonconvex losses that even though many local minima of the empirical risk exist, they are all close to the global minimum of population risk.

Our work builds on recent work in nonconvex optimization, in particular, results on escaping saddle points and finding approximate local minima. Beyond the classical result by Nesterov (2004) for finding first-order stationary points by gradient descent, recent work has given guarantees for escaping saddle points by gradient descent (Jin et al., 2017a) and stochastic gradient descent (Ge et al., 2015). Agarwal et al. (2017) and Carmon et al. (2016) established faster rates using algorithms that make use of Nesterov’s accelerated gradient descent in a nested-loop procedure (Nesterov, 1983), and Jin et al. (2017b) have established such rates even without the nested loop. There have also been empirical studies on various types of local minima (e.g. Keskar et al., 2016; Dinh et al., 2017).

Finally, our work is also related to the literature on zero-th order optimization or more generally, bandit convex optimization. Our algorithm uses function evaluations to construct a gradient estimate and perform SGD, which is similar to standard methods in this community (e.g., Flaxman et al., 2005; Agarwal et al., 2010; Duchi et al., 2015). Compared to first-order optimization, however, the convergence of zero-th order methods is typically much slower, depending polynomially on the underlying dimension even in the convex setting (Shamir, 2013). Other derivative-free optimization methods include simulated annealing (Kirkpatrick et al., 1983) and evolutionary algorithms (Rechenberg and Eigen, 1973), whose convergence guarantees are less clear.

Preliminaries

Our goal is to find a point that has zero gradient and positive semi-definite Hessian, thus escaping saddle points. We formalize this idea as follows.

$\mathbf{x}$ is called a second-order stationary point (SOSP) or approximate local minimum of a function $F$ if

We note that there is a slight difference between SOSP and local minima—an SOSP as defined here does not preclude higher-order saddle points, which themselves can be NP-hard to escape from (Anandkumar and Ge, 2016).

Since an SOSP is characterized by its gradient and Hessian, and since convergence of algorithms to an SOSP will depend on these derivatives in a neighborhood of an SOSP, it is necessary to impose smoothness conditions on the gradient and Hessian. A minimal set of conditions that have become standard in the literature are the following.

A function $F$ is $\rho$-Hessian Lipschitz if $~{}\forall\mathbf{x},\mathbf{y}~{}~{}\|{\nabla^{2}F(\mathbf{x})-\nabla^{2}F(\mathbf{y})}\|\leq\rho\|{\mathbf{x}-\mathbf{y}}\|.$

Another common assumption is that the function is bounded.

A function $F$ is $B$-bounded if for any $\mathbf{x}$ that $|F(\mathbf{x})|\leq B$.

For any finite-time algorithm, we cannot hope to find an exact SOSP. Instead, we can define $\epsilon$-approximate SOSP that satisfy relaxations of the first- and second-order optimality conditions. Letting $\epsilon$ vary allows us to obtain rates of convergence.

$\mathbf{x}$ is an $\epsilon$-second-order stationary point ($\epsilon$-SOSP) of a $\rho$-Hessian Lipschitz function $F$ if

Given these definitions, we can ask whether it is possible to find an $\epsilon$-SOSP in polynomial time under the Lipchitz properties. Various authors have answered this question in the affirmative.

Main Results

In the setting we consider, there is an unknown function $F$ (the population risk) that has regularity properties (bounded, gradient and Hessian Lipschitz). However, we only have access to a function $f$ (the empirical risk) that may not even be everywhere differentiable. The only information we use is that $f$ is pointwise close to $F$. More precisely, we assume

$f,F$ are $\nu$-pointwise close; i.e., $\|F-f\|_{\infty}\leq\nu$.

As we explained in Section 2, our goal is to find second-order stationary points of $F$ given only function value access to $f$. More precisely:

Given a function pair ($F,f$) that satisfies Assumption A1, find an $\epsilon$-second-order stationary point of $F$ with only access to values of $f$.

The only way our algorithms are allowed to interact with $f$ is to query a point $\mathbf{x}$, and obtain a function value $f(\mathbf{x})$. This is usually called a zero-th order oracle in the optimization literature. In this paper we give tight upper and lower bounds for the dependencies between $\nu$, $\epsilon$ and $d$, both for algorithms with polynomially many queries and in the information-theoretic limit.

There are three main difficulties in applying stochastic gradient descent to Problem 1: (1) in order to converge to a second-order stationary point of $F$, the algorithm must avoid being stuck in saddle points; (2) the algorithm does not have access to the gradient of $f$; (3) there is a gap between the observed $f$ and the target $F$, which might introduce non-smoothness or additional local minima. The first difficulty was addressed in Jin et al. (2017a) by perturbing the iterates in a small ball; this pushes the iterates away from any potential saddle points. For the latter two difficulties, we apply Gaussian smoothing to $f$ and use $\mathbf{z}[f(\mathbf{x}+\mathbf{z})-f(\mathbf{x})]/\sigma^{2}$ ($\mathbf{z}\sim\mathcal{N}(0,\sigma^{2}\mathbf{I})$) as a stochastic gradient estimate. This estimate, which only requires function values of $f$, is well known in the zero-th order optimization literature (e.g. Duchi et al., 2015). For more details, see Section 4.1.

In short, our algorithm (Algorithm 1) is a variant of SGD, which uses $\mathbf{z}[f(\mathbf{x}+\mathbf{z})-f(\mathbf{x})]/\sigma^{2}$ as the gradient estimate (computed over mini-batches), and adds isotropic perturbations. Using this algorithm, we can achieve the following trade-off between $\nu$ and $\epsilon$.

Theorem 7 shows that assuming a small enough function error $\nu$, ZPSGD will solve Problem 1 within a number of queries that is polynomial in all the problem-dependent parameters. The tolerance on function error $\nu$ varies inversely with the number of dimensions, $d$. This rate is in fact optimal for all polynomial queries algorithms. In the following result, we show that the $\epsilon,\rho,$ and $d$ dependencies in function difference $\nu$ are tight up to a logarithmic factors in $d$.

2 Information-theoretic guarantees

If we allow an unlimited number of queries, we can show that the upper and lower bounds on the function error tolerance $\nu$ no longer depends on the problem dimension $d$. That is, Problem 1 exhibits a statistical-computational gap—polynomial-queries algorithms are unable to achieve the information-theoretic limit. We first state that an algorithm (with exponential queries) is able to find an $\epsilon$-SOSP of $F$ despite a much larger value of error $\nu$. The basic algorithmic idea is that an $\epsilon$-SOSP must exist within some compact space, such that once we have a subroutine that approximately computes the gradient and Hessian of $F$ at an arbitrary point, we can perform a grid search over this compact space (see Section D for more details):

We also show a corresponding information-theoretic lower bound that prevents any algorithm from even identifying a second-order stationary point of $F$. This completes the characterization of function error tolerance $\nu$ in terms of required accuracy $\epsilon$.

3 Extension: Gradients pointwise close

Overview of Analysis

In this section we present the key ideas underlying our theoretical results. We will focus on the results for algorithms that make a polynomial number of queries (Theorems 7 and 8).

First, we introduce Gaussian smoothing, which perturbs the current point $\mathbf{x}$ using a multivariate Gaussian and then takes an expectation over the function value.

In ZPSGD, we use the stochastic gradients suggested by Lemma 12. Perturbed Gradient Descent (PGD) (Jin et al., 2017a) was shown to converge to a second-order stationary point. Here we use a simple modification of PGD that relies on batch stochastic gradient. In order for PSGD to converge, we require that the stochastic gradients are well-behaved; that is, they are unbiased and have good concentration properties, as asserted in the following lemma. It is straightforward to verify given that we sample $\mathbf{z}$ from a zero-mean Gaussian (proof in Appendix A.2).

As it turns out, these assumptions suffice to guarantee that perturbed SGD (PSGD), a simple adaptation of PGD in Jin et al. (2017a) with stochastic gradient and large mini-batch size, converges to the second-order stationary point of the objective function.

2 Polynomial queries lower bound

The proof of Theorem 8 depends on the construction of a ‘hard’ function pair. The argument crucially depends on the concentration of measure in high dimensions. We provide a proof sketch in Appendix B and the full proof in Appendix C.

Applications

In this section, we present several applications of our algorithm. We first show a simple example of learning one rectified linear unit (ReLU), where the empirical risk is nonconvex and nonsmooth. We also briefly survey other potential applications for our model as stated in Problem 1.

Due to nonsmoothness at $\bm{0}$ even for population risk, we focus on a compact region $\mathfrak{B}=\{\mathbf{w}|\mathbf{w}^{\top}\mathbf{w}^{\star}\geq\frac{1}{\sqrt{d}}\}\cap\{\mathbf{w}|\|{\mathbf{w}}\|\leq 2\}$ which excludes $\bm{0}$. This region is large enough so that a random initialization has at least constant probability of being inside it. We also show the following properties that allow us to apply Algorithm 1 directly:

These properties show that the population loss has a well-behaved landscape, while the empirical risk is pointwise close. This is exactly what we need for Algorithm 1. Using Theorem 7 we immediately get the following sample complexity, which guarantees an approximate population risk minimizer. We defer all proofs to Appendix G.

2 Other applications

Data privacy is a significant concern in machine learning as it creates a trade-off between privacy preservation and successful learning. Previous work on differentially private machine learning (e.g. Chaudhuri et al., 2011) have studied objective perturbation, that is, adding noise to the original (convex) objective and optimizing this perturbed objective, as a way to simultaneously guarantee differential privacy and learning generalization: $f=F+p(\varepsilon)$. Our results may be used to extend such guarantees to nonconvex objectives, characterizing when it is possible to optimize $F$ even if the data owner does not want to reveal the true value of $F(\mathbf{x})$ and instead only reveals $f(\mathbf{x})$ after adding a perturbation $p(\varepsilon)$, which depends on the privacy guarantee $\varepsilon$.

Motivated by the problem of adversarial examples in machine learning, there has been a lot of recent interest (e.g. Steinhardt et al., 2017; Sinha et al., 2018) in a form of robust optimization that involves a minimax problem formulation: $\min_{\mathbf{x}}\max_{\mathbf{u}}G(\mathbf{x},\mathbf{u}).$ The function $F(\mathbf{x})=\max_{\mathbf{u}}G(\mathbf{x},\mathbf{u})$ tends to be nonconvex in such problems, since $G$ can be very complicated. It can be intractable or costly to compute the solution to the inner maximization exactly, but it is often possible to get a good enough approximation $f$, such that $\sup_{\mathbf{x}}|F(\mathbf{x})-f(\mathbf{x})|=\nu$. It is then possible to solve $\min_{\mathbf{x}}f(\mathbf{x})$ by ZPSGD, with guarantees for the original optimization problem.

Acknowledgments

We thank Aditya Guntuboyina, Yuanzhi Li, Yi-An Ma, Jacob Steinhardt, and Yang Yuan for valuable discussions.

References

Appendix A Efficient algorithm for optimizing the population risk

As we described in Section 4, in order to find a second-order stationary point of the population loss $F$, we apply perturbed stochastic gradient on a smoothed version of the empirical loss $f$. Recall that the smoothed function is defined as

In this section we will also consider a smoothed version of the population loss $F$, as follows:

This function is of course not accessible by the algorithm and we only use it in the proof of convergence rates.

In the proof, we frequently require the following lemma (see e.g. Zhang et al. (2017)) that gives alternative expressions for the gradient and Hessian of a smoothed function.

Using the density function of a multivariate Gaussian, we may compute the gradient of the smoothed function as follows:

and similarly, we may compute the Hessian of the smoothed function:

For a twice-differentiable function, its gradient Lipshitz constant is also the upper bound on the spectral norm of its Hessian.

The last inequality follows from Lemma 26. ∎

A.1.2 Hessian Lipschitz

The last inequality follows from Lemma 27 and 28. ∎

A.1.3 Gradient Difference

Then the result follows from Lemma 30 and 31 ∎

A.1.4 Hessian Difference

A.2 Properties of the stochastic gradient

We prove the properties of the stochastic gradient, $\mathbf{g}(\mathbf{x};\mathbf{z})$, as stated in Lemma 24, restated as follows. Intuitively this lemma shows that the stochastic gradient is well-behaved and can be used in the standard algorithms.

Note that $\langle\frac{\mathbf{u}}{\|{\mathbf{u}}\|},\frac{\mathbf{z}}{\sigma}\rangle\sim\mathcal{N}(0,1)$. Thus, for $X\sim\mathcal{N}(0,\frac{B\|{u}\|}{\sigma})$,

This shows that $\mathbf{g}$ is sub-Gaussian with parameter $\frac{B}{\sigma}$.

implies $\mathbf{x}^{*}$ is an $O(\epsilon)$-SOSP of $F$.

By applying Lemma 13 and Weyl’s inequality, we have that the following inequalities hold up to a constant factor:

We know Eq.(2), (3) $\implies\sigma\leq\frac{\sqrt{\rho\epsilon}}{\rho\sqrt{d}}=\sqrt{\frac{\epsilon}{\rho d}}$ and $\sigma\leq\sqrt{\frac{\epsilon}{\rho d}}$.

Also Eq. (2), (3) $\implies\nu\leq\sigma\epsilon\leq\sqrt{\frac{\epsilon^{3}}{\rho d}}$ and $\nu\leq\sqrt{\rho\epsilon}\sigma^{2}\leq\sqrt{\rho\epsilon}\frac{\epsilon}{\rho d}=\sqrt{\frac{\epsilon^{3}}{\rho d^{2}}}$.

Finally Eq.(4) $\implies\nu\leq\rho\sqrt{d}\sigma^{3}\leq\frac{\epsilon^{1.5}}{\rho^{0.5}d}$.

Thus the following choice of $\sigma$ and $\nu$ ensures that $\mathbf{x}^{*}$ is an $O(\epsilon)$-SOSP of $F$:

where $c_{1},c_{2}$ are universal constants. ∎

A.4 Technical Lemmas

In this section, we collect and prove the technical lemmas used in the proofs of the above.

By the Hessian-Lipschitz property of $F$:

For brevity, denote $h=\frac{1}{(2\pi\sigma^{2})^{\frac{d}{2}}}$.

where $\omega(\Delta):=\frac{(\mathbf{z}+\Delta)(\mathbf{z}+\Delta)^{\top}-\sigma^{2}\mathbf{I}}{\sigma^{4}}e^{-\frac{\|{\mathbf{z}+\Delta}\|^{2}}{2\sigma^{2}}}$ and $\Delta=\frac{\mathbf{y}-\mathbf{x}}{2}$. Equality (5) follows from a change of variables. Now, denote $g(\mathbf{z}):=(f-F)(\mathbf{z}+\frac{\mathbf{x}+\mathbf{y}}{2})$.

Using $\omega(\Delta)=\frac{(\mathbf{z}+\Delta)(\mathbf{z}+\Delta)^{\top}-\sigma^{2}\mathbf{I}}{\sigma^{4}}e^{-\frac{\|{\Delta}\|^{2}+2\langle\Delta,\mathbf{z}\rangle}{2\sigma^{2}}}e^{-\frac{\|{\mathbf{z}}\|^{2}}{2\sigma^{2}}}$, we have the following

By a Taylor expansion up to only the first order terms in $\Delta$,

The last inequality follows from Lemma 29. ∎

where the last step is correct due to the second inequality we proved. The proof of the fourth inequality directly follows from the third inequality. ∎

The last inequality follows from Lemma 32. ∎

Inequality (6) follows by applying a generalization of mean-value theorem to vector-valued functions. ∎

Appendix B Overview for polynomial queries lower bound

The first property is due to the concentration of measure in high dimensions. The latter two properties are intuitively shown in Figure 5. These properties suggest a natural construction for $f$:

To see why this construction gives a hard instance of Problem 1, recall that the direction $\mathbf{v}$ is uniformly random. Since the direction $\mathbf{v}$ is unknown to the algorithm at initialization, the algorithm’s first query is independent of $\mathbf{v}$ and thus is likely to be in region $S_{\mathbf{v}}$, due to property 1. The queries inside $S_{\mathbf{v}}$ give no information about $\mathbf{v}$, so any polynomial-time algorithm is likely to continue to make queries in $S_{\mathbf{v}}$ and eventually fail to find $\mathbf{v}$. On the other hand, by property 2 above, finding an $\epsilon$-SOSP of $F$ requires approximately identifying the direction of $\mathbf{v}$, so any polynomial-time algorithm will fail with high probability.

We deal with first problem by constructing a function $\bar{F}(\mathbf{y})$ on $^{d}$, ignoring the boundary condition, and then composing it with a smooth periodic function. For the second problem, we carefully construct a smooth function $h$, as shown in Figure 5, to have zero function value, gradient and Hessian at the boundary of the ball and outside the ball, so that no algorithm can make use of the padding region to identify an SOSP of $F$. Details are deferred to section C in the appendix.

Appendix C Constructing Hard Functions

In this section, we prove Theorem 8, the lower bound for algorithms making a polynomial number of queries. We start by describing the hard function construction that is key to the lower bound.

We will first present a “scale-free” version of the hard function, where we assume $\rho=1$ and $\epsilon=1$. In section C.2, we will show how to scale this hard function to prove Theorem 8.

where $h(\mathbf{y})=h_{1}(\mathbf{v}^{\top}\mathbf{y})\cdot h_{2}(\sqrt{\|{\mathbf{y}}\|^{2}-(\mathbf{v}^{\top}\mathbf{y})^{2}})$, and

and the vector $\mathbf{v}$ is uniformly distributed on the $d$-dimensional unit sphere.

We will state the properties of the hard instance by breaking the space into different regions:

“hypercube” $H=[-\frac{\pi}{2},\frac{\pi}{2}]^{d}$ be the $d$-dimensional hypercube with side length $\pi$.

“band” $S_{\mathbf{v}}=\{\mathbf{x}\in S:\langle\sin\mathbf{x},\mathbf{v}\rangle\leq\frac{\log d}{\sqrt{d}}\}$

We also call the union of $S_{2}$ and $S_{\mathbf{v}}$ the “non-informative” region.

Our construction happens within the ball. However it is hard to fill the space using balls, so we pad the ball into a hypercube. Our construction will guarantee that any queries to the non-informative region do not reveal any information about $\mathbf{v}$. Intuitively the non-informative region is very large so that it is hard for the algorithm to find any point outside of the non-informative region (and learn any information about $\mathbf{v}$).

Let $F,f,\mathbf{v}$ be as defined in equations (7), (8). Then $F,f$ satisfies:

$f$ in the non-informative region $S_{2}\cup S_{\mathbf{v}}$ is independent of $\mathbf{v}$.

$F$ has no SOSP in the non-informative region $S_{2}\cup S_{\mathbf{v}}$.

$F$ is $O(d)$-bounded, $O(1)$-Hessian Lipschitz, and $O(1)$-gradient Lipschitz.

These properties will be proved based on the properties of $h(\mathbf{y})$, which we defined in (7) to be the product of two functions.

Property 1. On $S_{\mathbf{v}}$, $f(\mathbf{x})=F(\mathbf{x})=\|{\sin\mathbf{x}}\|^{2}$, which is independent of $\mathbf{v}$. On $S_{2}$, we argue that $h(\sin\mathbf{x})=0$ and therefore $f(\mathbf{x})=\|{\sin\mathbf{x}}\|^{2}~{}\forall\mathbf{x}\in S_{2}$. Note that on $S_{2}$, $\|{\mathbf{x}}\|>3/\mu$ and $(\sin x)^{2}>(\frac{2x}{\pi})^{2}~{}\forall|x|<\frac{\pi}{2}$, so

Therefore, $h(\sin\mathbf{x})=h_{1}(\mathbf{v}^{\top}\sin\mathbf{x})\cdot h_{2}(\sqrt{\|{\sin\mathbf{x}}\|^{2}-(\mathbf{v}^{\top}\sin\mathbf{x})^{2}})=0$.

For $\mathbf{x}\in S_{\mathbf{v}}$, we have $\mathbf{v}^{\top}\sin\mathbf{x}\in[-\frac{\log d}{\sqrt{d}},\frac{\log d}{\sqrt{d}}]$. By symmetry, we may just consider the case where $\mathbf{v}^{\top}\sin\mathbf{x}>0$.

Here $C>0$ is a large enough universal constant.

Property 3. (Part I.) We show that there are no SOSP in $S_{2}$. For $\mathbf{x}:\|{\mathbf{x}}\|>3/\mu$, we argue that either the gradient is large, due to contribution from $\|{\sin\mathbf{x}}\|^{2}$, or the Hessian has large negative eigenvalue (points close to the boundary of $H$). Denote $G(\mathbf{x})=\|{\sin\mathbf{x}}\|^{2}$. We may compute the gradient of $G$ as follows:

where $\xi\approx 0.95$ is the positive root of the equation $\sin 2x=x$. On $S_{2}$, $\nabla F(\mathbf{x})=\nabla G(\mathbf{x})$, so $\|{\nabla F(\mathbf{x})}\|>\frac{3}{\mu}=1\times 10^{-2}$ for $\mathbf{x}\in S_{2}\cap[-\xi,\xi]^{d}$. We may also compute the Hessian of $G$:

Since $(\frac{\pi}{4}-\xi)<-0.15$, $\lambda_{\min}(\nabla^{2}F(\mathbf{x}))<-0.3$.

(Part II.) We argue that $F$ has no SOSP in $S_{\mathbf{v}}$. For $\mathbf{y}=\sin(\mathbf{x})$, we consider two cases: (i) $z=\sqrt{\|{\mathbf{y}}\|^{2}-(\mathbf{v}^{\top}\mathbf{y})^{2}}$ large and (ii) $z$ small.

Write $g(\mathbf{x})=h(\sin\mathbf{x})$, and denote $\left.\nabla h(\mathbf{x})\right|_{\sin\mathbf{x}},\left.\nabla^{2}h(\mathbf{x})\right|_{\sin\mathbf{x}}$ with $\nabla h(\mathbf{y}),\nabla^{2}h(\mathbf{y})$. Let $\mathbf{u}\circ\mathbf{v}$ denote the Schur product of $\mathbf{u}$ and $\mathbf{v}$. We may compute the gradient and Hessian of $g$:

Now we change the coordinate system such that $\mathbf{v}=(1,0,\cdots,0)$. $\|{\nabla h(\mathbf{y})}\|$ and $\lambda_{\min}(\nabla^{2}h(\mathbf{y}))$ are invariant to such a transform. Under this coordinate system, $h(\mathbf{y})=h_{1}(y_{1})\cdot h_{2}(\sqrt{\|{\mathbf{y}}\|^{2}-(y_{1})^{2}})$

(i): $z\geq\frac{1}{2\mu}$. We show that $\|{\nabla F}\|$ is large.

Let $\mathcal{P}_{-1}(\mathbf{u})$ denote the projection of $\mathbf{u}$ onto the orthogonal component of the first standard basis vector.

Since $\forall~{}i\neq 1,\frac{\partial}{\partial y_{i}}h(\mathbf{y})=h_{1}(y_{1})h_{2}^{\prime}(z)\frac{y_{i}}{z}$, we have

(ii): $z<\frac{1}{2\mu}$. We show that $\nabla^{2}F(\mathbf{x})$ has large negative eigenvalue. First we compute the second derivative of $h$ in the direction of the first coordinate:

Now we use this to upper bound the smallest eigenvalue of $\nabla^{2}F(\mathbf{x})$.

Property 4. $O(1)$-bounded: Lemma 34 shows that $|h(\mathbf{y})|\leq 1$. $\|{\sin\mathbf{x}}\|^{2}\leq d$. Therefore $|F|\leq 1+d$.

$O(1)$-gradient Lipschitz: $\|{\nabla^{2}F(\mathbf{x})}\|\leq\|{\nabla^{2}G(\mathbf{x})}\|+\|{\nabla^{2}g(\mathbf{x})}\|$. We know $\|{\nabla^{2}G(\mathbf{x})}\|\leq 2$.

$O(1)$-Hessian Lipschitz: First bound the Hessian Lipschitz constant of $G(\mathbf{x})$.

Now we bound the Hessian Lipschitz constant of $g(\mathbf{x})$. Denote $\mathbf{A}(\mathbf{x})=diag(\cos\mathbf{x})$ and $\mathbf{B}(\mathbf{x})=diag(\sin\mathbf{x})$.

Therefore $F(\mathbf{x})$ is $(2.8\times 10^{10})$-Hessian Lipschitz.

Now we need to prove smoothness properties of $h(\mathbf{y})$ that are used in the previous proof. In the following lemma, we prove that $h(\mathbf{y})$ as defined in equation (7) is bounded, Lipschitz, gradient-Lipschitz, and Hessian-Lipschitz.

$h(\mathbf{y})$ as given in Definition 7 is O(1)-bounded, O(1)-Lipschitz, O(1)-gradient Lipschitz, and O(1)-Hessian Lipschitz.

WLOG assume $\mathbf{v}=(1,0,\cdots,0)^{\top}$. Denote $u=y_{1},\mathbf{w}=(y_{2},\cdots,y_{d})^{\top}$. Let $\otimes$ denote tensor product.

Note that $|h_{1}^{\prime}|\leq 3\mu,|h_{2}^{\prime}|\leq 2\mu,|h_{1}^{\prime\prime}|\leq 32\mu^{2},|h_{2}^{\prime\prime}|\leq 12\mu^{2},|h_{1}^{\prime\prime\prime}|\leq 300\mu^{3},|h_{2}^{\prime\prime\prime}|\leq 48\mu^{3}$. Assume $\mu>2$.

O(1)-bounded: $|h|\leq|h_{1}|\cdot|h_{2}|\leq 1$.

O(1)-Lipschitz: $\|{\nabla h(\mathbf{y})}\|=\sqrt{h_{2}(\|{\mathbf{w}}\|)h_{1}^{\prime}(u)+h_{1}(u)h_{2}^{\prime}(\|{\mathbf{w}}\|)}\leq 3\mu\leq O(1)$.

$\|{\nabla h_{1}(u)}\|\leq 3\mu$. Notice that the following are also O(1):

Therefore, $\|{\nabla^{2}h(\mathbf{y})}\|\leq 24\mu^{2}+32\mu^{2}+2\cdot 2\mu\cdot 3\mu\leq 68\mu^{2}$.

O(1)-Hessian Lipschitz: We first argue that $\nabla^{2}h_{2}(\mathbf{w})$ is Lipschitz. For $\|{\mathbf{w}}\|\geq 1/\mu$, $\nabla^{2}h_{2}(\mathbf{w})=0$. So we consider $\|{\mathbf{w}}\|<\mu$. We obtain the following by direct computation.

We may easily check that indeed $\lim_{\|{\mathbf{w}}\|\to 1}\|{\nabla^{2}h_{2}(\mathbf{w})}\|=0$.

Therefore $\nabla^{2}h_{2}(\mathbf{w})$ is $144\mu^{3}$-Lipschitz.

By triangle inequality, using the above, we obtain

This proves that $h(\mathbf{y})$ is $1000\mu^{3}$-Hessian Lipschitz.

C.2 Scaling the Hard Instance

Now we show how to scale the function we described in order to achieve the final lower bound with correct dependencies on $\epsilon$ and $\rho$.

where $r=\sqrt{\epsilon/\rho}$ and $F,f$ are defined as in Equation 7. Define the ‘scaled’ regions:

This is implied by Lemma 33. To see this, notice

We have simply scaled each coordinate axis by $r$.

C.3 Proof of the Theorem

We are now ready to state the two main lemmas used to prove Theorem 8.

Denote unit vector $\hat{\mathbf{y}}=\sin\frac{\mathbf{x}}{r}/\|{\sin\frac{\mathbf{x}}{r}}\|$. This gives:

In last inequality, we used Lemma 36, which finishes the proof. ∎

Now we have all the ingredients to prove Theorem 8, restated below more formally.

Note that because $\|{\mathbf{x}_{t}}\|\leq r/100$, $\|{\sin\frac{1}{r}\mathbf{x}_{t}}\|\leq\|{\mathbf{x}_{t}}\|/r\leq 1/100$.

For completeness, we now state the classical result showing that most of the surface area of a sphere lies close to the equator; it was used in the proof of Lemma 36.

Appendix D Information-theoretic Limits

In this section, we prove upper and lower bounds for algorithms that may run in exponential time. This establishes the information-theoretic limit for problem 1. Compared to the previous (polynomial time) setting, now the dependency on dimension $d$ is removed.

We first restate our upper bound, first stated in Theorem 9, below.

The algorithm is based on a procedure to estimate the gradient and Hessian at point $x$. This procedure will be applied to a exponential-sized covering of a compact space to find an SOSP.

where $r$ is scalar in the order of $O(\sqrt{\epsilon/\rho})$.

We will first show that any solution of this problem will give good estimates of the gradient and Hessian of $F$.

Any solution $(\mathbf{g},\mathcal{H})$to the above feasibility problem, Eq.(10), gives

When we have $\|{f-F}\|_{\infty}\leq\nu$, above feasibility problem is equivalent to solve following:

Due to the Hessian-Lipschitz property, we have $|F(\mathbf{y})-F(\mathbf{x})-\nabla F(\mathbf{x})^{\top}(\mathbf{y}-\mathbf{x})-(\mathbf{y}-\mathbf{x})^{\top}\nabla^{2}F(\mathbf{x})(\mathbf{y}-\mathbf{x})|\leq\frac{1}{6}\rho r^{3}$, this means above feasibility problem is also equivalent to:

Given $\nu\leq\frac{1}{c}\sqrt{\frac{\epsilon^{3}}{\rho}}$ for large enough constant $c$, and picking $r=c^{\prime}\sqrt{\frac{\epsilon}{\rho}}$ with proper constant $c^{\prime}$, we prove the lemma. ∎

We then argue that (10) always has a solution.

$(\nabla F(\mathbf{x}),\nabla^{2}F(\mathbf{x}))$ is clearly one solution to the feasibility problem Eq.(10). Then, this lemma is true due to Hessian Lipschitz and gradient Lipschitz properties of $F$. ∎

Now, since the algorithm can do an exhaustive search over a compact space, we just need to prove that there is an $\epsilon$-SOSP within a bounded distance.

Combining all these lemmas we are now ready to prove the main theorem of this section:

We show that Algorithm 2 is guaranteed to succeed within a number of function value queries of $f$ that is exponential in all problem parameters. First, by Lemma 43, we know that at least one of $\{\mathbf{x}_{t}\}_{t=1}^{N}$ must be an $O(\epsilon)$-SOSP of $F$. It suffices to show that for any $\mathbf{x}$ that is an $O(\epsilon)$-SOSP, Algorithm 2’s subroutine will successfully return $\mathbf{x}$, that is, it must find a solution $\mathbf{g},\mathcal{H}$, to the feasibility problem 10 that satisfies $\|{\mathbf{g}}\|\leq O(\epsilon)\text{ and }\lambda_{\min}(\mathcal{H})\geq-O(\sqrt{\rho\epsilon})$.

Next, notice that because all the covers in Algorithm 2 have size at most $O((d/\epsilon)^{d^{2}})$ must terminate in $O(e^{d^{2}\log\frac{d}{\epsilon}})$ steps. ∎

In the following two lemmas, we provide simple methods for constructing an $\epsilon$-cover for a ball (as well as a sphere), and for matrices with bounded spectral norm.

For any point $\mathbf{y}$ in the ball, we can find $\mathbf{x}\in C$ such that $|\mathbf{y}_{i}-\mathbf{x}_{i}|\leq\epsilon/\sqrt{d}$ for each $i\in[d]$. By the Pythagorean theorem, this implies $\|{\mathbf{y}-\mathbf{x}}\|\leq\epsilon$. ∎

For any matrix $M$ in $\mathcal{M}$, we can find $N\in C$ such that $|N_{i,k}-M_{i,k}|\leq\epsilon/d$ for each $i,k\in[d]$. Since the Frobenius norm dominates the spectral norm, we have $\|{N-M}\|\leq\|{N-M}\|_{F}\leq\epsilon$. ∎

D.2 Information-theoretic Lower bound

Appendix E Extension: Gradients pointwise close

Given function pair ($F,f$) that satisfies Assumption A1, find an $\epsilon$-second-order stationary point of $F$ with only access to function values of $\mathbf{g}=\nabla f$.

Appendix F Proof of Extension: Gradients pointwise close

Recall the definition of the gradient smoothing of a function given in Definition 11. In this section we will consider a smoothed version of the (possibly erroneous) gradient oracle, defined as follows.

We proceed by exchanging the order of differentiation. The last equality follows from applying lemma 19 to the function $\frac{\partial}{\partial x_{j}}f(\mathbf{x}+\mathbf{z})$

We will prove the 4 claims of the lemma one by one, in the following 4 sub-subsections.

The last inequality follows from Lemma 55. ∎

F.1.2 Hessian Lipschitz

The last inequality follows from Lemmas 27 and 56. ∎

F.1.3 Gradient Difference

The inequality at (13) follows from Lemma 31. ∎

F.1.4 Hessian Difference

The last inequality follows from Lemma 23 and 49. ∎

F.2 Properties of the stochastic gradient

For the second claim, since function $f$ is L-Lipschitz, we know $\|{\mathbf{g}(\mathbf{x},\mathbf{z})}\|=\|{\nabla f(\mathbf{x}+\mathbf{z})}\|\leq L$. This implies that $\mathbf{g}(\mathbf{x};\mathbf{z})$ is sub-Gaussian with parameter $L$. ∎

F.3 Proof of Theorem 47

implies $\mathbf{x}^{*}$ is an $O(\epsilon)$-SOSP of $F$.

By Lemma 48 and Weyl’s inequality, we have that the following inequalities hold up to a constant factor:

We know Eq.(14), (15) $\implies\sigma\leq\frac{\sqrt{\rho\epsilon}}{\rho\sqrt{d}}=\sqrt{\frac{\epsilon}{\rho d}}$ and $\sigma\leq\sqrt{\frac{\epsilon}{\rho d}}$.

Thus the following choices ensures $\mathbf{x}^{*}$ is an $O(\epsilon)$-SOSP of $F$:

F.4 Technical lemmas

In this section, we collect and prove the technical lemmas used in section F.

For brevity, denote $h=\frac{1}{(2\pi\sigma^{2})^{\frac{d}{2}}}$. We have:

where $\Delta=\frac{\mathbf{y}-\mathbf{x}}{2}$. The last equality follows from a change of variables. Now denote $\mathbf{g}(\mathbf{z}):=(\nabla f-\nabla F)(\mathbf{z}+\frac{\mathbf{x}+\mathbf{y}}{2})$. By a Taylor expansion up to only the first order terms in $\Delta$, we have

The last inequality follows from Lemma 55.

Appendix G Proof of Learning ReLU Unit

In this section we analyze the population loss of the simple example of a single ReLU unit.

Recall our assumption that $\|{\mathbf{w}^{\star}}\|=1$ and that the data distribution is $\mathbf{x}\sim\mathcal{N}(0,\mathbf{I})$; thus,

We use the squared loss as the loss function, hence writing the empirical loss as:

The main tool we use is a closed-form formula for the kernel function defined by ReLU gates.

(Cho and Saul, 2009) For fixed $\mathbf{u},\mathbf{v}$, if $\mathbf{x}\sim\mathcal{N}(0,\mathbf{I})$, then

where $\theta$ is the angel between $\mathbf{u}$ and $\mathbf{v}$ satisfying $\cos\theta=\mathbf{u}^{\top}\mathbf{v}/(\|{\mathbf{u}}\|\|{\mathbf{v}}\|)$.

Then, the population loss has the following analytical form:

and so does the gradient ($\hat{\mathbf{w}}$ is the unit vector along $\mathbf{w}$ direction):

We first prove the properties of the population loss, which were stated in Lemma 16 and we also restate the lemma below. Let $\mathfrak{B}=\{\mathbf{w}|\mathbf{w}^{\top}\mathbf{w}^{\star}\geq\frac{1}{\sqrt{d}}\}\cap\{\mathbf{w}|\|{\mathbf{w}}\|\leq 2\}$.

The population and empirical risk $R,\hat{R}_{n}$ of learning a ReLU unit problem satisfies:

If $\mathbf{w}_{0}\in\mathfrak{B}$, then runing ZPSGD (Algorithm 1) gives $\mathbf{w}_{t}\in\mathfrak{B}$ for all $t$ with high probability.

Inside $\mathfrak{B}$, $R$ is $O(1)$-bounded, $O(\sqrt{d})$-gradient Lipschitz, and $O(d)$-Hessian Lipschitz.

Inside $\mathfrak{B}$, $R$ is nonconvex function, $\mathbf{w}^{\star}$ is the only SOSP of $R(\mathbf{w})$.

To prove these four claims, we require following lemmas.

The first important property we use is that the gradient of population loss $l$ has the one-point convex property inside $\mathfrak{B}$, stated as follows:

Note that inside $\mathfrak{B}$, we have the angle $\theta\in[0,\pi/2)$. Also, let $\mathfrak{W}_{\theta}=\{\mathbf{w}|\angle(\mathbf{w},\mathbf{w}^{\star})=\theta\}$, then for $\theta\in[0,\pi/2)$:

On the other hand, note that $\theta\leq 2\sin\theta$ holds true for $\theta\in[0,\pi/2)$; thus we have:

where the second last inequality used the fact that $\sin\theta\leq\|{\mathbf{w}-\mathbf{w}^{\star}}\|$ for all $\mathbf{w}\in\mathfrak{B}$. ∎

One-point convexity guarantees that ZPSGD stays in the region $\mathfrak{B}$ with high probability.

ZPSGD (Algorithm 1) with proper hyperparameters will stay in $\mathfrak{B}$ with high probability.

The algorithm always moves towards $\mathbf{x}^{\star}$ in the region $\mathfrak{B}-\{\|{\mathbf{w}-\mathbf{w}^{\star}}\|\leq 1/10\}$.

The algorithm will not jump from $\{\|{\mathbf{w}}\|\leq 1/10\}$ to $\mathfrak{B}^{c}$ in one step.

Let $\mathbf{w}(t)=\frac{1}{5}(\mathbf{w}^{\star}+t\mathbf{e})$ where $\mathbf{e}$ is any direction so that $\mathbf{e}^{\top}\mathbf{w}^{\star}=0$

which is nonconvex in domain $t\in$. Therefore $f(\mathbf{w}(t))$ is nonconvex along this line segment inside $\mathfrak{B}$.

Note that in above setup, $\tan\theta=t$, so the population loss can be calculated as:

It’s easy to show $\mathbf{w}(t)\in\mathfrak{B}$ for all $t\in$ and if $g(t)=R(\mathbf{w}(t))$, then $g^{\prime\prime}(0.6)<0$ and thus the function is nonconvex. ∎

Next, we show that the empirical risk and the population risk are close by a covering argument.

For sample size $n\geq d$, with high probability, we have:

Let $\{\mathbf{w}^{j}\}_{j=1}^{J}$ be a $\epsilon$-covering of $\mathfrak{B}$. By triangular inequality:

where $\mathbf{w}^{j}$ is the closest point in the cover to $\mathbf{w}$. Clearly, the $\epsilon$-net of $\mathfrak{B}$ requires fewer points than the $\epsilon$-net of $\{\mathbf{w}|\|{\mathbf{w}}\|\leq 2\}$. By the standard covering number argument, we have $\log N_{\epsilon}=O(d\log\frac{1}{\epsilon})$. We proceed to bound each term individually.

Term $T_{2}$: For a fixed $j$, we know $\hat{R}_{n}(\mathbf{w}^{j})=\frac{1}{n}\sum_{i=1}^{n}(y_{i}-\text{ReLU}(\mathbf{x}_{i}^{\top}\mathbf{w}^{j}))^{2}$, where $y_{i}-\text{ReLU}(\mathbf{x}_{i}^{\top}\mathbf{w}^{j})$ is sub-Gaussian with parameter $O(1)$, thus $(y_{i}-\text{ReLU}(\mathbf{x}_{i}^{\top}\mathbf{w}^{j}))^{2}$ is sub-Exponential with parameter $O(1)$. We have the concentration inequality:

That is, with $n\geq d$, and probability $1-\delta$, we have:

Term $T_{3}$: Since the population loss is $O(1)$-Lipschitz in $\mathfrak{B}$, we have:

Term $T_{1}$: Note that for a fixed pair $(\mathbf{x}_{i},\mathbf{y}_{i})$, the function $g_{i}(\mathbf{w})=(y_{i}-\text{ReLU}(\mathbf{x}_{i}^{\top}\mathbf{w}))^{2}$ is $O(\|{\zeta_{i}}\|\|{\mathbf{x}_{i}}\|+\|{\mathbf{x}_{i}}\|^{2})$-Lipschitz. Therefore,

With high probability, $\frac{1}{n}\sum_{i}\left[\|{\zeta_{i}}\|\|{\mathbf{x}_{i}}\|+\|{\mathbf{x}_{i}}\|^{2}\right]$ concentrates around its mean, $O(d)$.

By picking $\epsilon$ (for the $\epsilon$-covering) small enough, we finish the proof. ∎

Finally we prove the smoothness of population risk in $\mathfrak{B}$, we have $1/\sqrt{d}\leq\|{\mathbf{w}}\|\leq 2$.

For population loss $R(\mathbf{w})=\frac{1}{4}\|{\mathbf{w}}\|^{2}+\frac{5}{4}-\frac{1}{2\pi}\|{\mathbf{w}}\|[\sin\theta+(\pi-\theta)\cos\theta]$, its gradient and Hessian are equal to:

where $\hat{\mathbf{w}}$ is the unit vector along the $\mathbf{w}$ direction, and $\hat{\mathbf{u}}$ is the unit vector along the $\mathbf{w}^{\star}-\hat{\mathbf{w}}\cos\theta$ direction.

Note $\|{\mathbf{w}^{\star}}\|=1$. Let $z(\mathbf{w},\mathbf{w}^{\star})=\frac{\mathbf{w}^{\top}\mathbf{w}^{\star}}{\|{\mathbf{w}}\|}$, we have:

Since $\cos\theta=z(\mathbf{w},\mathbf{w}^{\star})$, we obtain:

Therefore, the Hessian (when $\theta\neq 0$):

where $\hat{\mathbf{u}}$ is the unit vector along $\mathbf{w}^{\star}-\hat{\mathbf{w}}\cos\theta$ direction.

And for $\theta=0$, Hessian $\nabla^{2}R(\mathbf{w})=\frac{1}{2}\mathbf{I}$. We prove this by taking the limit. For $\hat{\mathbf{v}}=\mathbf{w}^{\star}$

For any $\hat{\mathbf{v}}\perp\mathbf{w}^{\star}$, the angle $\theta$ between $\mathbf{w}+\epsilon\hat{v}$ and $\mathbf{w}^{\star}$ is $\Theta(\frac{\epsilon}{\|{\mathbf{w}}\|})$ up to first order in $\epsilon$, we have:

The population loss function $R$ is $O(1)$-bounded, $O(1)$-Lipschitz, $O(\sqrt{d})$-gradient Lipschitz, and $O(d)$-Hessian Lipschitz.

The bounded, Lipschitz, and gradient Lipschitz are all very straightforward given the formula of gradient and Hessian. We will focus on proving Hessian Lipschitz. Equivalently, we show upper bounds on following quantity:

Note that the change in $\theta$ is at most $O(\frac{\epsilon}{\|{\mathbf{w}}\|})$, we have:

For four claims in Lemma 16, claim 1 follows from Lemma 60; claim 2 follows from Lemma 64; claim 3 follows from Lemma 62; claim 4 follows from Lemma 61 and Lemma 59. ∎

G.2 Proof of Theorem 17

Appendix H Proof of Stochastic gradient descent

Here for completeness we give the result for perturbed stochastic gradient descent, which is a adaptation of results in Jin et al. (2017a) and will be formally presented in Jin et al. (2018).

function $f$ satisfies following property:

for any $\lambda>0,\delta>0$, if minibatch size $m\geq\frac{2\lambda^{2}\sigma^{2}}{\epsilon^{2}}\log\frac{d}{\delta}$, then for a fixed $\mathbf{x}$, with probability $1-\delta$, we have:

This lemma means, when mini-batch size is large enough, we can make noise in the stochastic gradient descent polynomially small.

Consider the setting of Theorem 65, if $\|{\nabla f(\mathbf{x}_{t})}\|\geq\epsilon$, then by running Algorithm 1, with probability $1-\delta$, we have $f(\mathbf{x}_{t+1})-f(\mathbf{x}_{t})\leq-\eta\epsilon^{2}/4$.

By gradient Lipschitz, and the fact $\|{\xi_{t}}\|\leq\epsilon/20$ and with minibatch size $m$ large enough, with high probability we have $\|{\nabla f(\mathbf{x}_{t})-\mathbf{g}_{t}}\|\leq\epsilon/20$. Let $\zeta_{t}=\mathbf{g}_{t}-\nabla f(\mathbf{x}_{t})+\xi_{t}$, by triangle inequality, we have $\|{\zeta_{t}}\|\leq\epsilon/10$ and update equation $\mathbf{x}_{t+1}=\mathbf{x}_{t}-\eta(\nabla f(\mathbf{x}_{t})+\zeta_{t})$:

Consider the setting of Theorem 65, if $\|{\nabla f(\mathbf{x}_{t})}\|\leq\epsilon$ and $\lambda_{\min}(\nabla^{2}f(\mathbf{x}_{t}))\leq-\sqrt{\rho\epsilon}$, then by running Algorithm 1, with probability $1-\delta$, we have $f(\mathbf{x}_{t+\mathscr{T}})-f(\mathbf{x}_{t})\leq-\mathscr{F}$.

Combining lemma 67 and 68, we know with probability $1-\frac{\Delta_{f}}{\mathscr{F}}\delta$, algorithm will find $\epsilon$-second order stationary point in following iterations:

where $\zeta_{t}=\mathbf{g}_{t}-\nabla f(\mathbf{x}_{t})+\xi_{t}$.

(Improve or Localize) Suppose $\{\mathbf{x}_{t}\}_{t=0}^{T}$ is a SGD sequence, then for all $t\leq T$:

where $\zeta_{t}=\mathbf{g}_{t}-\nabla f(\mathbf{x}_{t})+\xi_{t}$

For any $t\leq T-1$, by Lemma 69, we have:

Finally, by Cauchy-Schwarz, we have for all $t\leq T$:

To study escaping saddle points, we need a notion of coupling. Recall the PSGD update has two source of randomness: $\mathbf{g}_{t}-\nabla f(\mathbf{x}_{t})$ which is the stochasticity inside the gradient oracle and $\xi_{t}$ which is the perturbation we deliberately added into the algorithm to help escape saddle points. Let $\text{SGD}_{\xi}^{(t)}(\cdot)$ denote the update via SGD $t$ times with perturbation $\xi=\{\xi_{2},\cdots\}$ fixed. Define Stuck region:

Intuitively, the later perturbations of coupling sequence are the same, while the very first perturbation is used to escape saddle points.

Since $\mathbf{x}_{\mathscr{T}}$ and $\mathbf{x}^{\prime}_{0}$ are independent, we have

In the remaining proof, we will proceed proving Eq.(22) by showing two steps:

$\min\{f(\mathbf{x}_{\mathscr{T}})-f(\mathbf{x}_{0}),f(\mathbf{x}^{\prime}_{\mathscr{T}})-f(\mathbf{x}^{\prime}_{0})\}\leq-2\mathscr{F}$ with probability $1-\delta$

The final result immediately follow from triangle inequality.

Part 2. Assume the contradiction $\min\{f(\mathbf{x}_{\mathscr{T}})-f(\mathbf{x}_{0}),f(\mathbf{x}^{\prime}_{\mathscr{T}})-f(\mathbf{x}^{\prime}_{0})\}\geq-2\mathscr{F}$, by Lemma 70 (note $\|{\zeta_{t}}\|\leq\|{\mathbf{g}_{t}-\nabla f(\mathbf{x}_{t})}\|+\|{\xi_{t}}\|\leq 2r$ with high probability when $m$ is large enough), with $1-\delta/2$ probability, this implies localization:

That is, the first term is always the dominating term. It is easy to check for base case $t=0$; we have $0\leq\|{\mathbf{w}_{0}}\|/2$. Suppose for all $t^{\prime}\leq t$ the induction holds, this gives:

where the third last inequality use the fact $\mathbf{w}_{0}$ is along minimum eigenvector direction of $\mathcal{H}$, the last inequality uses the fact $\eta\rho\mathscr{S}T=c^{-1}\leq 1/4$ for $c$ large enough.

On the other hand, with $1-\delta/2$ probability, we also have:

where the last inequality requires $\max_{\tau}\|{\bm{h}_{\tau}}\|\leq\gamma\|{\mathbf{w}_{0}}\|$ which can be achieved by making minibatch size $m$ large enough. Now, by triangular inequality, we finishes the induction.

Let $r_{0}=\delta r\sqrt{\frac{2\pi}{d}}$ and applying Lemma 71, we know $\mathcal{X}^{\xi}_{\text{stuck}}(\mathbf{x}_{t})$ has at most width $\eta r_{0}$ in the minimum eigenvector direction of $\nabla^{2}f(\mathbf{x}_{t})$ and thus,

Therefore the probabilty of escaping saddle point is $(1-\delta)(1-\sqrt{\delta})\geq 1-2\sqrt{\delta}$. Reparametrizing $\delta^{\prime}=2\sqrt{\delta}$ only affects constant factors in $\chi$, hence we finish the proof. ∎