text:
\begin{document}
\maketitle

\begin{abstract}

Although stochastic optimization is central to modern machine learning, the precise mechanisms underlying its success, and in particular, the precise role of the stochasticity, still remain unclear. Modelling stochastic optimization algorithms as discrete random recurrence relations, we show that multiplicative noise, as it commonly arises due to variance in local rates of convergence, results in heavy-tailed stationary behaviour in the parameters. A detailed analysis is conducted for SGD applied to a simple linear regression problem, followed by theoretical results for a much larger class of models (including non-linear and non-convex) and optimizers (including momentum, Adam, and stochastic Newton), demonstrating that our qualitative results hold much more generally. In each case, we describe dependence on key factors, including step size, batch size, and data variability, all of which exhibit similar qualitative behavior to recent empirical results on state-of-the-art neural network models from computer vision and natural language processing. Furthermore, we empirically demonstrate how multiplicative noise and heavy-tailed structure improve capacity for basin hopping and exploration of non-convex loss surfaces, over commonly-considered stochastic dynamics with only additive noise and light-tailed structure.

\end{abstract}

\section{Introduction}

Stochastic optimization is the process of minimizing a deterministic objective function via the simulation of random elements, and it is one of the most successful methods for optimizing complex or unknown objectives.
Relatively simple stochastic optimization procedures---in particular, stochastic gradient descent (SGD)---have become the backbone of modern machine learning (ML) \cite{ma2017power}. To improve understanding of stochastic optimization in ML, and particularly why SGD and its extensions work so well, recent theoretical work has sought to study its properties and dynamics~\cite{liu2019deep}.

Such analyses typically approach the problem through one of two perspectives.
The first perspective, an \emph{optimization (or quenching) perspective}, examines convergence either in expectation \cite{chen2018convergence,zhou2018convergence,gower2019sgd,jain2019sgd,fontaine2020continuous} or with some positive (high) probability \cite{roosta2016sub,du2017gradient,kleinberg2018alternative,ward2018adagrad} through the lens of a deterministic counterpart.
This perspective inherits some limitations of deterministic optimizers, including assumptions (e.g., convexity, Polyak-\L{}ojasiewicz criterion, etc.) that are either not satisfied by state-of-the-art problems, or not strong enough to imply convergence to a quality (e.g., global) optimum.
More concerning, however, is the inability to explain what has come to be known as the ``generalization gap'' phenomenon:
increasing stochasticity by reducing batch size appears to improve generalization performance \cite{keskar2016large,martin2018implicit}.
Empirically, existing strategies do tend to break down for inference tasks when using large batch sizes \cite{golmant2018computational}.
The second perspective, a \emph{probabilistic (annealing) perspective}, examines algorithms through the lens of Markov process theory \cite{freidlin1998random,henderson2003theory,mandt2016variational,dieuleveut2017bridging}.
Here, stochastic optimizers are interpreted as samplers from probability distributions concentrated around optima \cite{mandt2017stochastic}, and annealing the optimizer (by reducing step size) increasingly concentrates probability mass around global optima \cite{henderson2003theory}.
Traditional analyses trade restrictions on the objective for precise annealing schedules that guarantee adequate mixing and ensure convergence.
However, it is uncommon in practice to consider step size schedules that decrease sufficiently slowly as to guarantee convergence to global optima with probability one \cite{li2019budgeted}.
In fact, SGD based methods with poor initialization can easily get stuck near poor local minima using a typical step decay schedule \cite{liu2019bad}.

More recent efforts conduct a \emph{distributional analysis}, directly examining the probability distribution that a stochastic optimizer targets for each fixed set of hyperparameters \cite{mandt2016variational,babichev2018constant,dieuleveut2017bridging}.
Here, one can assess a stochastic optimizer according to its capacity to reach and then occupy neighbourhoods of high-quality optima in the initial stages, where the step size is large and constant.
As the step size is then rapidly reduced, tighter neighbourhoods with higher probability mass surrounding nearby minima are achievable. This is most easily accomplished using a variational approach by appealing to continuous-time Langevin approximations \cite{mandt2016variational, chaudhari2018stochastic}, whose stationary distributions are known explicitly \cite{ma2015complete}; but these approaches also require restrictive assumptions, such as constant or bounded volatility \cite{mandt2017stochastic,orvieto2019continuous}.
Interestingly, these assumptions parallel the common belief that the predominant part of the stochastic component of an optimizer is an additive perturbation \cite{kleinberg2018alternative,zhang2019algorithmic}.
Such analyses have been questioned with recent discoveries of non-Gaussian noise \cite{simsekli2019tail} that leads to \emph{heavy-tailed} stationary behaviour (whose distributions have infinite Laplace transform).
This behavior implies stronger exploratory properties and an increased tendency to rapidly reach faraway basins than earlier Langevin-centric analyses suggest; but even these recent non-Gaussian analyses have relied on strong assumptions, additive noise, and continuous-time approximations of their own \cite{simsekli2020fractional}.

\paragraph{Main Contributions.}  We model stochastic optimizers as Markov random recurrence relations, thereby avoiding continuous-time approximations, and we examine their stationary distributions.
The formulation of this model is described in \S\ref{sec:Model}.
We show that \emph{multiplicative noise}, frequently assumed away in favour of more convenient \emph{additive noise} in continuous analyses, can often lead to heavy-tailed stationary behaviour, playing a critical role in the dynamics of a stochastic optimizer, and influencing its capacity to hop between basins in the loss landscape.
In this paper, we consider heavy-tailed behavior that assumes a power law functional form.
We say that the stationary distributions of the parameters/weights $W$ exhibit \emph{power laws}, with tail probabilities $\mathbb{P}(\|W\| > t) = \Omega(t^{-\alpha})$ as $t \to \infty$, for some $\alpha > 0$ called the \emph{tail exponent} (where smaller tail exponents correspond to heavier tails) --- further details, including precise definitions are in Appendix \appref{A}{\ref{sec:LinearRecurrence}}.
To inform our analysis, in \S\ref{sec:LinearOptim}, we examine the special case of constant step-size SGD applied to linear least squares, and we find it obeys a \emph{random linear recurrence relation} displaying both multiplicative and additive noise.
Using well-known results \cite{buraczewski2016stochastic}, we isolate three regimes determining the tail behaviour of SGD (shown in Table~\ref{tab:Regimes}, discussed in detail in \S\ref{sec:LinearOptim}), finding stationary behaviour \emph{always} exhibits a precise power law in an infinite data regime.
\begin{table}
\centering
\begin{tabular}{|l|c|c|c|}
\hline
\rowcolor{Gray}\textbf{Regime} & \textbf{Condition on $A$} & \textbf{Tails for $B$} & \textbf{Tails for $W_{\infty}$} \\
\hline
Light-tailed (\S\appref{A.1}{\ref{sec:GoldieGrubel}}) & $\mathbb{P}(\|A\| \leq 1) = 1$ & $\mathcal{O}(e^{-\lambda t})$ & $\mathcal{O}(e^{-\rho t})$ \\
Heavy-tailed (multiplicative) (\S\appref{A.2}{\ref{sec:KestenGoldie}}) & $\|\sigma_{\min}(A)\|_{\alpha} = 1$ & $o(t^{-\alpha})$ & $\Theta(t^{-\alpha})$ \\
Heavy-tailed (additive) (\S\appref{A.3}{\ref{sec:GrinceviciusGrey}}) & --- & $\Omega(t^{-\beta})$ & $\Omega(t^{-\beta})$ \\
\hline
\end{tabular}
\caption{\label{tab:Regimes} Summary of the three primary tail behaviour regimes for the stationary distribution of (\ref{eq:LinearRecurrence}).}
\end{table}
\begin{table}
    \centering
    \begin{tabular}{|l|c|c|c|}
    \hline
    \rowcolor{Gray}\textbf{Regime} & \textbf{Noise} & \textbf{Condition on }$r$ & \textbf{Quadratic Example} \\
\hline
        Light-tailed \cite{ma2015complete} & White noise & $\liminf_{\|w\|\to\infty} r(w) > 0$ & $g$ bounded \\
        Heavy-tailed (multiplicative) \cite{ma2015complete} & White noise & $\limsup_{\|w\|\to\infty} r(w) = 0$ & $g(w) = c \nabla f(w)$\\
        Heavy-tailed (additive) \cite{simsekli2019tail} & L\`{e}vy noise & --- & arbitrary $g$ \\
\hline
    \end{tabular}
    \caption{Summary of three tail behaviour regimes for continuous-time models (\ref{eq:LangevinModel}).}
    \label{tab:Continuous}
\end{table}
In \S\ref{sec:General}, these observations are extended by providing sufficient conditions for the existence of power laws arising in arbitrary iterative stochastic optimization algorithms on both convex and non-convex problems, with more precise results when updates are Lipschitz. Using our results, factors influencing tail behaviour are examined, with existing empirical findings supporting the hypothesis that heavier tails coincide with improved generalization performance \cite{martin2020predicting}.
Numerical experiments are conducted in \S\ref{sec:Numerics}, illustrating how multiplicative noise and heavy-tailed stationary behaviour improve capacity for basin hopping (relative to light-tailed stationary behaviour) in the exploratory phase of learning.
We finish by discussing impact on related work in \S\ref{sec:Discussion}, including a continuous-time analogue of Table~\ref{tab:Regimes} (shown in Table~\ref{tab:Continuous}).

\paragraph{Related Work.}
There is a large body of related work, and we review only the most directly related here.
Analysis of stochastic optimizers via stationary distributions of Markov processes was recently considered in \cite{mandt2016variational, babichev2018constant, dieuleveut2017bridging}.
The latter, in particular, examined first and second moments of the stationary distribution, although these can be ineffective measures of concentration in heavy-tailed settings. Heavy tails in machine learning have been observed and empirically examined in spectral distributions of weights \cite{martin2017rethinking,martin2018implicit,martin2019traditional,martin2020heavy,martin2020predicting} and in the weights themselves \cite{simsekli2019tail,panigrahi2019non}, but (ML style) theoretical analyses on the subject seem limited to continuous-time examinations \cite{simsekli2019tail,simsekli2020fractional}\footnote{Immediately prior to uploading this article to the arXiv, we became aware of \cite{gurbuzbalaban2020heavy} (submitted to the arXiv a few days ago), which conducts a detailed analysis of heavy tails in the stationary distribution of SGD applied to least-squares linear regression. Our discussion in \S\ref{sec:LinearOptim} and Appendix \appref{A}{\ref{sec:LinearRecurrence}}, which for us serves predominantly as motivation in this special case, recovers a number of their results.  Our primary contributions in \S\ref{sec:Model} and \S\ref{sec:General} provide extensions of these results to more general stochastic optimization algorithms and more general stochastic optimization problems. We encourage readers to see \cite{gurbuzbalaban2020heavy} for details in the linear least-squares~setting. }.
Connections between multiplicative noise and heavy-tailed fluctuations can be seen throughout the wider literature \cite{deutsch1993generic, frisch1997extreme, sornette1997convergent, buraczewski2016stochastic}. From a physical point of view, multiplicative noise acts as an external environmental field, capable of exhibiting drift towards low energy states \cite{volpe2016effective}. Hysteretic optimization \cite{pal2006hysteretic} is one example of a stochastic optimization algorithm taking advantage of this property. To our knowledge, no theoretical analysis of this phenomenon has been conducted in a general optimization setting.
Indeed, while multiplicative noise in stochastic optimization has been explored in some recent empirical analyses \cite{wu2020noisy,zhang2018removing,holland2018robust}, its impact appears underappreciated, relative to the well-studied and exploited effects of additive noise  \cite{ge2015escaping,jin2017escape,du2017gradient,kleinberg2018alternative}.
Section~\ref{sec:Discussion} contains a discussion of additional related work in light of our results.

\paragraph{Notation.}  We let $\m I_d$ denote the $d \times d$ identity matrix. Unless otherwise stated, we default to the following norms: for vector $x$, we let $\|x\| = (\sum_{i=1}^n x_i^2)^{1/2}$ denote the Euclidean norm of $x$, and for matrix $\m A$, the norm $\|\m A\| = \sup_{\|x\| = 1}\|\m A x\|$ denotes the $\ell^2$-induced (spectral) norm. Also, for matrix $\m A$, $\|\m A\|_F = \|\mathrm{vec}(\m A)\|$ is the Frobenius norm, and $\sigma_{\min}(A)$, $\sigma_{\max}(A)$ are its smallest and largest singular values, respectively. For two vectors/matrices $A,B$, we let $A\otimes B$ denote their tensor (or Kronecker) product. Knuth asymptotic notation \cite{knuth1976big} is adopted. The notation $(\Omega, \mathcal{E},\mathbb{P})$ is reserved to denote an appropriate underlying probability space. For two random elements $X,Y$, $X \eqd Y$ if $X,Y$ have the same distribution. Finally, for random vector/matrix $X$, for $\alpha > 0$, $\|X\|_\alpha = (\mathbb{E}\|X\|^{\alpha})^{1/\alpha}$. All proofs are relegated to Appendix \appref{C}{\ref{app:Proofs}}.

\section{Stochastic optimization as a Markov chain}
\label{sec:Model}
In this section, we describe how to model a stochastic optimization algorithm as a Markov chain --- in particular, a \emph{random recurrence relation}.
This formulation is uncommon in the ML literature, but will be important for our analysis.
Consider a general single-objective stochastic optimization problem, where the goal is to solve problems of the form $w^{\ast} = \argmin_w \mathbb{E}_\mathcal{D} \ell(w, X)$
for some scalar loss function $\ell$, and random element $X \sim \mathcal{D}$ (the data) \cite[\S12]{kroese2013handbook}.
In the sequel, we shall assume the weights $w$ occupy a vector space $S$ with norm $\|\cdot\|$.
To minimize $\ell$ with respect to $w$, some form of fixed point iteration is typically adopted. Supposing that there exists some continuous map $\Psi$ such that any fixed point of $\mathbb{E}\Psi(\cdot,X)$ is a minimizer of $\ell$, the sequence of iterations
\begin{equation}
\label{eq:StochOptExp}
W_{k+1} = \mathbb{E}_\mathcal{D} \Psi(W_k, X)
\end{equation}
either diverges, or converges to a minimizer of $\ell$ \cite{granas2013fixed}. In practice, this expectation may not be easily computed, so one could instead consider the sequence of iterated random functions
\begin{equation}
\label{eq:StochAlgoMarkov}
W_{k+1} = \Psi(W_k, X_k),\qquad X_k \iidsim \mathcal{D},\qquad \mbox{for }k=0,1,\dots.
\end{equation}
For example, one could consider the Monte Carlo approximation of (\ref{eq:StochOptExp}):
\begin{equation}
\label{eq:StochBatchMarkov}
W_{k+1} = \frac1n \sum_{i=1}^n \Psi(W_k, X_i^{(k)}),\quad \text{with }X_i^{(k)} \iidsim \mathcal{D},\mbox{ for } k=0,1,\dots,
\end{equation}
which can be interpreted in the context of \emph{randomized subsampling with replacement}, where each $(X_i^{(k)})_{i=1}^n$ is a batch of $n$ subsamples drawn from a large dataset $\mathcal{D}$. Subsampling without replacement can also be treated in the form (\ref{eq:StochBatchMarkov}) via the Markov chain $\{W_{sk}\}_{k=0}^{\infty}$, where $s$ is the number of minibatches in each~epoch.

Since (\ref{eq:StochAlgoMarkov}) will not, in general, converge to a single point, the \emph{stochastic approximation (SA)} approach of Robbins--Munro \cite{robbins1951stochastic} considers a corresponding sequence of maps $\Psi_k$ given by $\Psi_k(w, x) \coloneqq (1 - \gamma_k)w + \gamma_k \Psi(w, x)$
where $\{\gamma_k\}_{k=1}^{\infty}$ is a decreasing sequence of \emph{step sizes}. Provided $\Psi$ is uniformly bounded, $\sum_k \gamma_k = \infty$ and $\sum_k \gamma_k^2 < \infty$, the sequence of iterates $W_{k+1} = n^{-1}\sum_{i=1}^n \Psi_k(W_k, X_i^{(k)})$ converges in $L^2$ and almost surely to a minimizer of $\ell$ \cite{blum1954approximation}.
Note that this differs from the \emph{sample average approximation (SAA)} approach \cite{kleywegt2002sample}, where a deterministic optimization algorithm is applied to a random objective (e.g., gradient descent on subsampled empirical risk).

Here are examples of how popular ML stochastic optimization algorithms fit into this framework.
\begin{example}[SGD \& SGD with momentum]
Minibatch SGD with step size $\gamma$ coincides with (\ref{eq:StochBatchMarkov}) under
$\Psi(w, X) = w - \gamma \nabla \ell(w, X).$
Incorporating momentum is possible by considering the augmented space of weights and velocities together.
In the standard setup, letting $v$ denote velocity, and $w$ the weights, $\Psi((v,w),X) = (\eta v + \nabla \ell(w,X), w - \gamma (\eta v + \nabla \ell(w,X)))$.
\end{example}
\begin{example}[Adam]
Using state-space augmentation, the popular first-order Adam optimizer \cite{kingma2015adam} can also be cast into the form of (\ref{eq:StochAlgoMarkov}).
Indeed, for $\beta_1,\beta_2,\epsilon > 0$, letting $g(w,X) = \sum_{i=1}^n \nabla \ell(w,X_i)$, for $M(m,w,X) = \beta_1 m + n^{-1}(1-\beta_1) g(w,X)$, $V(v,w,X) = \beta_2 v + (1-\beta_2)n^{-2} g(w,X)^2$, $B_1(b) = 1 - \beta_1(1-b)$, $B_2(b) = 1 - \beta_2(1-b)$, and
\[
W(b_1,b_2,m,v,w,X) = w - \eta \frac{M(m, w, X) / B_1(b_1)}{\sqrt{V(v, w, X) / B_2(b_2)} + \epsilon},
\]
iterations of the Adam optimizer satisfy $\mathcal{W}_{k+1} = \Psi(\mathcal{W}_k, X_k)$, where
$
\Psi((b_1,b_2,m,v,w),X) = (B_1(b_1), \allowbreak B_2(b_2), \allowbreak M(m,w,X), \allowbreak V(v,w,X), \allowbreak W(b_1,b_2,m,v,w,X)).
$
\end{example}
\begin{example}[Stochastic Newton]
The formulation (\ref{eq:StochAlgoMarkov}) is not limited to first-order stochastic optimization. Indeed, for $g(w,X) = \sum_{i=1}^n \nabla \ell(w,X_i)$ and $H(w,X) = \sum_{i=1}^n \nabla^2 \ell(w, X_i)$, the choice $\Psi(w, X) = w - \gamma H(w,X)^{-1}g(w,X)$ coincides with the stochastic Newton method \cite{roosta2016sub}.
\end{example}

The iterations (\ref{eq:StochAlgoMarkov}) are also sufficiently general to incorporate other random effects, e.g., dropout \cite{srivastava2014dropout}.
Instead of taking the SA approach, in our analysis we will examine the original Markov chain (\ref{eq:StochAlgoMarkov}), where hyperparameters (including step size) are fixed.
This will provide a clearer picture of ongoing dynamics at each stage within an arbitrary learning rate schedule \cite{babichev2018constant}.
Also, assuming reasonably rapid mixing, global tendencies of a stochastic optimization algorithm, such as the degree of exploration in space, can be examined through the \emph{stationary distribution} of (\ref{eq:StochAlgoMarkov}).

\section{The linear case with SGD}
\label{sec:LinearOptim}
In this section, we consider the special case of ridge regression, i.e., least squares linear regression with $L^2$ regularization, using vanilla SGD.
Let $\mathcal{D}$ denote a dataset comprised of inputs $x_i \in \mathbb{R}^d$ and corresponding labels $y_i \in \mathbb{R}^m$.
Supposing that $(X,Y)$ is a pair of random vectors uniformly drawn from $\mathcal{D}$, for $\lambda \geq 0$, we seek a solution to
\begin{equation}
\label{eq:RidgeRegression}
M^{\ast} = \argmin_{M\in\mathbb{R}^{d\times m}} \tfrac12 \mathbb{E}_\mathcal{D}\|Y - M X\|^2 + \tfrac12 \lambda \|M\|_F^2.
\end{equation}
Observing that
\[
\frac{\partial}{\partial M} \left(\tfrac12 \|Y - M X\|^2 + \tfrac12 \lambda \|M\|_F^2\right) = M (X X^\top + \lambda I_d) - Y X^\top,
\]
one can show that the solution to (\ref{eq:RidgeRegression}) is
$$
M^{\ast} = \mathbb{E}_\mathcal{D}[YX^\top]\mathbb{E}_\mathcal{D}[XX^\top + \lambda I_d]^{-1}  .
$$
(Note that $\lambda$ need not be strictly positive for the inverse to exist.)
On the other hand, applying minibatch SGD to solving (\ref{eq:RidgeRegression}) with constant (sufficiently small) step size results in a Markov chain $\{M_k\}_{k=0}^{\infty}$ in the estimated parameter matrix.
We start with the following observation, which is immediate but important.
\begin{lemma}
\label{lem:LinearRecur}
Let $n$ denote the size of each minibatch $(X_{ik},Y_{ik})_{i=1}^n$ comprised of independent and identically distributed copies of $(X,Y)$ for $k=1,2,\dots$. For $W_k$ the vectorization of $M_k$, iterations of minibatch SGD undergo the following random linear recurrence relation
\begin{equation}
\label{eq:LinearRecurrence}
W_{k+1} = \m A_k W_k + \m B_k,
\end{equation}
where $
\m A_k = I_m \otimes ((1-\lambda\gamma)\m I_d - \gamma n^{-1} \sum_{i=1}^n X_{ik}X_{ik}^\top)$, and $\m B_k = \gamma n^{-1}\sum_{i=1}^n Y_{ik} \otimes X_{ik}$. If $A_k,B_k$ are non-atomic, $\log^+\|A_k\|$, $\log^+\|B_k\|$ are integrable, and $\mathbb{E}_\mathcal{D}\log\|A_k\| < 0$, then (\ref{eq:LinearRecurrence}) is ergodic.
\end{lemma}
From here, it immediately follows that, even in this simple setting, SGD possesses \emph{multiplicative noise} in the form of the factor $A_k$, as well as \emph{additive noise} in the form of the factor $B_k$.
Under the conditions of Lemma~\ref{lem:LinearRecur}, the expectations of (\ref{eq:LinearRecurrence}) converge to $M^{\ast}$.
Although the dynamics of this process, as well as the shape of its stationary distribution, are not as straightforward, random linear recurrence relations are among the most well-studied discrete-time processes, and as multiplicative processes, are well-known to exhibit heavy-tailed behaviour.
In Appendix~\appref{A}{\ref{sec:LinearRecurrence}}, we discuss classical results on the topic.
The three primary tail regimes are summarized in Table~\ref{tab:Regimes}, where each $\alpha,\beta,\lambda,\rho$ denotes a strictly positive value. In particular, we have the following result.
\begin{lemma}
\label{lem:LinearHeavy}
Assume that the distribution of $X$ has full support on $\mathbb{R}^d$. Then any stationary distribution of (\ref{eq:LinearRecurrence}) is heavy-tailed.
\end{lemma}
Lemma~\ref{lem:LinearHeavy} suggests that heavy-tailed fluctuations could be more common than previously considered.
This is perhaps surprising, given that common Langevin approximations of SGD applied to this problem exhibit light-tailed stationary distributions \cite{mandt2016variational,orvieto2019continuous}.
In (\ref{eq:LinearRecurrence}), light-tailed stationary distributions can only occur when $X$ does not have full support on $\mathbb{R}^d$ and the step size is taken sufficiently small.
Referring to Table~\ref{tab:Regimes}, there are two possible mechanisms by which the stationary distribution of (\ref{eq:LinearRecurrence}) can be heavy-tailed.
Most discussions about SGD focus on the additive noise component $B_k$ \cite{kleinberg2018alternative} (presumably since the analysis is simpler).
In this case, if $B_k$ is heavy-tailed, then the stationary distribution of (\ref{eq:LinearRecurrence}) is also heavy-tailed.
This is the assumption considered in \cite{simsekli2019tail,simsekli2020fractional}.

However, this is not the only way heavy-tailed noise can arise.
Indeed, the multiplicative heavy-tailed regime illustrates how a power law can arise from light-tailed data (in fact, even from data with finite support).
One of the most important attributes here is that the tail exponent of the stationary distribution is entirely due to the recursion and properties of the multiplicative noise $A_k$.
In fact, this is true even if $B_k$ displays (not too significant) heavy-tailed behaviour.
Here, the heavy-tailed behaviour arises due to \emph{intrinsic factors} to the stochastic optimization, and they tend to dominate over time.
Later, we shall go further and analyze factors influencing the heaviness of these tails, or more precisely, its \emph{tail exponent}, denoted by $\alpha$ in Table~\ref{tab:Regimes}.
Evidently, $\alpha$ is dependent on the dispersal of the distribution of $A_k$, which is itself dependent (in differing ways) on the step size $\gamma$, batch size $n$, and the dispersal and condition number of the covariance matrices $X X^\top$ of the input data.
These three factors have been noted to be integral to the practical performance of SGD \cite{smith2018disciplined,xing2018walk}.

\section{Power laws for general objectives and optimizers}\label{sec:General}

In this section, we consider the general case (\ref{eq:StochAlgoMarkov}), and we examine how heavy-tailed stationary behaviour can arise for any stochastic optimizer, applied to both convex and non-convex problems.
As in Section~\ref{sec:LinearOptim}, here we are most interested in the presence of heavy-tailed fluctuations due to multiplicative factors.
(The case when heavy-tailed fluctuations arise from the additive noise case is clear: if, for every $w \in \mathbb{R}^d$, $\Psi(w,X)$ is heavy-tailed/has infinite $\alpha$th-moment, then for each $k=1,2,\dots$, $W_k$ is heavy-tailed/has infinite $\alpha$th moment also, irrespective of the dynamics of $W$.)

In our main result (Theorem~\ref{thm:Lipschitz}), we illustrate how power laws arise in general Lipschitz stochastic optimization algorithms that are contracting on average and also strongly convex near infinity with positive probability. \begin{theorem}
\label{thm:Lipschitz}
Let $(S,\|\cdot\|)$ be a separable Banach space. Assume that $\Psi: S \times \Omega \to S$ is a random function on $S$ such that $\Psi$ is almost surely Lipschitz continuous and the probability measure of $\Psi$ is non-atomic. We let $K_\Psi$ denote a random variable such that $K_\Psi(\omega)$ is the Lipschitz constant of $\Psi(\cdot,\omega)$ for each $\omega \in \Omega$. Assume that $K_\Psi$ is integrable, and $\|\Psi(w)- w\|$ is integrable for some $w^{\ast} \in S$.
Suppose there exist non-negative random variables $k_\Psi, M_\Psi$ such that, with probability one,
\begin{equation}
\label{eq:LipschitzOptim}
k_\Psi \|w - w^{\ast}\| - M_\Psi \leq \|\Psi(w) - \Psi(w^{\ast})\| \leq K_\Psi \|w - w^{\ast}\|,\qquad \mbox{for all }w \in S.
\end{equation}
If $\mathbb{E}\log K_{\Psi} < 0$, then the Markov chain given by $W_{k+1} = \Psi_k(W_k), k=0,1,\dots,$
where each $\Psi_k$ is an independent and identically distributed copy of $\Psi$, is geometrically ergodic with $W_{\infty} = \lim_{k\to\infty} W_k$ satisfying the distributional fixed point equation $W_{\infty} \eqd \Psi(W_{\infty})$. Furthermore, if $k_\Psi > 1$ with positive probability and $\|\Psi(w^\ast)\|_\alpha<\infty$ for any $\alpha > 0$, then:
\begin{enumerate}[leftmargin=*]
\item There exist $\mu,\nu,C_\mu,C_\nu > 0$ such that
$
C_\mu (1+t)^{-\mu} \leq \mathbb{P}(\|W_{\infty}\| > t) \leq C_\nu t^{-\nu},\mbox{ for all }t > 0.
$
\item There exist $\alpha,\beta > 0$ such that $\|K_{\Psi}\|_{\beta} = 1$ and $\|k_{\Psi}\|_{\alpha} = 1$ and for $\delta > 0$,
\[
0 < \limsup_{t\to\infty} t^{\alpha+\epsilon}\mathbb{P}(\|W_{\infty}\| > t),\qquad \mbox{and}\qquad \limsup_{t\to\infty} t^{\beta-\epsilon}\mathbb{P}(\|W_{\infty}\| > t) < +\infty.
\]
\item For any $p \in (0,\beta)$ and $k=1,2,\dots$, the moments of $\|W_k\|_p$ are bounded above and below by
\[
\|k_\Psi\|_p^k (\|W_0\|_p + \kappa_p^-) - \kappa_p^- \leq \|W_k\|_p \leq \|K_\Psi\|_p^k (\|W_0\|_p - \kappa_p^+)+\kappa_p^+,\]
implying $\|W_\infty\|_p \leq \kappa_p^+$, where
\[
\kappa_{p}^{-}=\frac{\| k_{\Psi}\| _{p}\| w^{\ast}\| +\| \Psi(w^{\ast})\| _{p}}{1-\| k_{\Psi}\| _{p}},\qquad\mbox{and}\qquad\kappa_{p}^{+}=\frac{\| K_{\Psi}\| _{p}\| w^{\ast}\| +\| \Psi(w^{\ast})\| _{p}}{1-\| K_{\Psi}\| _{p}}.
\]
\end{enumerate}
\end{theorem}
Geometric rates of convergence in the Prohorov and total variation metrics to stationarity are discussed in \cite{diaconis1999iterated,alsmeyer2003harris}.
From Theorem~\ref{thm:Lipschitz}, we find that the presence of expanding multiplicative noise implies the stationary distribution of (\ref{eq:StochAlgoMarkov}) is stochastically bounded between two power laws. In particular, this suggests that Lipschitz stochastic optimizers satisfying (\ref{eq:LipschitzOptim}) and $\mathbb{P}(k_\Psi > 1) > 0$ can conduct wide exploration of the loss landscape. As an example, for SGD satisfying (\ref{eq:LipschitzOptim}), we have that
\begin{equation}
\label{eq:kPsiSGD}
k_\Psi = \liminf_{\|w\|\to\infty}\sigma_{\min}(I-\gamma\nabla^2\ell(w,X)),\qquad \mbox{and} \qquad K_\Psi = \sup_w \|I - \gamma\nabla^2\ell(w,X)\|.
\end{equation}
Therefore, for the linear case (\ref{eq:LinearRecurrence}), using Theorem \ref{thm:Lipschitz}, we can recover the same tail exponent seen in Theorem \ref{thm:Kesten}. To treat other stochastic optimizers that are not Lipschitz, or do not satisfy the conditions of Theorem~\ref{thm:Lipschitz}, we present in Lemma~\ref{lem:Abstract} an abstract sufficient condition for heavy-tailed stationary distributions of ergodic Markov chains.
\begin{lemma}
\label{lem:Abstract}
For a general metric space $S$, suppose that $W$ is an ergodic Markov chain on $S$ with $W_k \cd W_{\infty}$ as $k \to \infty$, and let $f: S \to \mathbb{R}$ be some scalar-valued function. If there exists some $\epsilon > 0$ such that $\inf_{w \in S}\mathbb{P}(|f(\Psi(w))| > (1 + \epsilon)|f(w)|) > 0$, then $f(W_\infty)$ is heavy-tailed.
\end{lemma}
Under our model, Lemma~\ref{lem:Abstract} says that convergent constant step-size stochastic optimizers will exhibit heavy-tailed stationary behaviour if there is some positive probability of the optimizer moving further away from an optimum, irrespective of how near or far you are from it. Analyses concerning stochastic optimization algorithms almost always consider rates of contraction towards a nearby optimum, quantifying the \emph{exploitative} behaviour of a stochastic optimizer. To our knowledge, this is the first time that consequences of the fact that stochastic optimization algorithms could, at any stage, move away from any optimum, have been examined. In the convex setting, this appears detrimental, but this is critical in non-convex settings, where \emph{exploration} (rather than exploitation to a local optimum) is important. Indeed, we find that it is this behaviour that directly determines the tails of the stochastic optimizer's stationary distribution, and therefore, its exploratory behaviour. Altogether, our results suggest that heavy-tailed noise could be more common and more beneficial in stochastic optimization than previously acknowledged.

\paragraph{Factors influencing the tail exponent.}  Recent analyses have associated heavier-tailed structures (i.e., larger tail exponents) with improved generalization performance \cite{martin2018implicit,martin2020predicting,simsekli2019tail}. From the bounds in Theorem~\ref{thm:Lipschitz}, and established theory for the special linear case (\ref{eq:LinearRecurrence}) (see Appendix \appref{A.2}{\ref{sec:KestenGoldie}}), the following factors play a role in \emph{decreasing the tail exponent} (generally speaking, increasing $k_\Psi,K_\Psi$ in expectation or dispersion implies decreased $\alpha$) resulting in \emph{heavier-tailed~noise}.
\begin{itemize}
\item
\textbf{Decreasing batch size:}
For fixed hyperparameters, $\Psi$ becomes less deterministic as the batch size decreases, resulting in an decreased (i.e., heavier) tail exponent for ergodic stochastic optimizers. This is in line with \cite{yao2018hessian,martin2018implicit,martin2020heavy,martin2020predicting}, and it suggests a relationship to the \emph{generalization gap} phenomenon.\item
\textbf{Increasing step size:}
Also in line with \cite{martin2018implicit,martin2020heavy,martin2020predicting}, increased step sizes exacerbate fluctuations and thus result in increased tail exponents.
However, step sizes also affect stability of the algorithm.
The relationship between step and batch sizes has received attention \cite{balles2016coupling,smith2017don}; choosing these parameters to increase heavy-tailed fluctuations while keeping variance sensible could yield a valuable exploration~strategy.
\item
\textbf{More dispersed data:}
Increasing dispersion in the data implies increased dispersion for the distribution of $k_\Psi,K_\Psi$, and hence heavier tails. For the same model trained to different datasets, a smaller tail exponent may be indicative of richer data, not necessarily of higher variance, but exhibiting a larger moment of some order. Data augmentation is a strategy to achieve this \cite{perez2017effectiveness}.
\item
\textbf{Increasing regularization:}
The addition of an explicit $L^2$-regularizer to the objective function (known to help avoiding bad minima \cite{liu2019bad}) results in larger $k_\Psi$, and hence, heavier-tailed noise.
\item
\textbf{Increasing average condition number of Hessian: }
The Hessian is known to play an important role in the performance of stochastic optimizers \cite{yao2018hessian,xing2018walk}. Evidently seen in (\ref{eq:kPsiSGD}) in the case of SGD, and potentially more generally, a larger condition number of the Hessian can allow for larger tail exponents while remaining ergodic.
\end{itemize}
In each case (excluding step size, which is complicated since step size also affects stability), the literature supports the hypothesis that factors influencing heavier tails coincide with improved generalization performance; see, e.g.,~\cite{martin2018implicit} and references therein. Vanilla SGD and stochastic Newton both exhibit heavy-tailed behaviour when $\nabla^2 \ell(w,X)$ and the Jacobian of $H(w,X)^{-1}g(w,X)$ have unbounded spectral distributions, respectively (as might be the case when $X$ has full support on $\mathbb{R}^d$). However, adaptive optimizers such as momentum and Adam incorporate geometric decay that can prevent heavy-tailed fluctuations, potentially limiting exploration while excelling at exploitation to nearby optima. It has recently been suggested that these adaptive aspects should be turned off during an initial warmup phase \cite{liu2019variance}, implying that exacerbating heavy-tailed fluctuations and increased tail exponents could aid exploratory behaviour in the initial stages of training.

\section{Numerical experiments}
\label{sec:Numerics}

To illustrate the advantages that multiplicative noise and heavy tails offer in non-convex optimization, in particular for exploration (of the entire loss surface) versus exploitation (to a local optimum), we first consider stochastic optimization algorithms in the non-stationary regime.
To begin, in one dimension, we compare perturbed gradient descent (GD) \cite{jin2017escape} with additive light-tailed and heavy-tailed noise against a version with additional multiplicative noise.
That is,
\[
\mbox{(a,b)}\quad w_{k+1} = w_k - \gamma (f'(w_k) + (1+\sigma) Z_k),\quad \mbox{vs.}\quad \mbox{(c)}\quad w_{k+1} = w_k - \gamma ((1+\sigma Z_k^{(1)})f'(w_k) + Z_k^{(2)}),
\]
where $f:\,\mathbb{R}\to\mathbb{R}$ is the objective function, $\gamma,\sigma > 0$, $Z_k,Z_k^{(1)},Z_k^{(2)} \iidsim \mathcal{N}(0,1)$ in (a) and (c), and in (b), $Z_k$ are i.i.d. $t$-distributed (heavy-tailed) with 3 degrees of freedom, normalized to have unit variance. Algorithms (a)-(c) are applied to the objective $f(x) = \frac1{10}x^2 + 1 - \cos(x^2)$, which achieves a global (wide) minimum at zero. Iterations of (a)-(c) have expectation on par with classical GD. For fixed step size $\gamma = 10^{-2}$ and initial $w_0 = -4.75$, the distribution of $10^6$ successive iterates are presented in Figure~\ref{fig:BasinMulAdd} for small ($\sigma = 2$), moderate ($\sigma = 12$), and strong ($\sigma = 50$) noise. Both (b) and (c) readily enable jumps between basins. However, while (a), (b) smooth the effective loss landscape, making it easier to jump between basins, it also has the side effect of reducing resolution in the vicinity of minima. On the other hand, (c) maintains close proximity (peaks in the distribution) to critical points.

Figure~\ref{fig:BasinMulAdd} also illustrates exploration/exploitation benefits of multiplicative noise (and associated heavier tails) over additive noise (and associated lighter tails). While rapid exploration may be achieved using heavy-tailed additive noise \cite{simsekli2019tail}, since reducing the step size may not reduce the tail exponent, efficient exploitation of local minima can become challenging, even for small step sizes. On the other hand, multiplicative noise has the benefit of behaving similarly to Gaussian additive noise for small step sizes \cite{rubin2014mapping}.  We can see evidence of this behaviour in the leftmost column in Figure~\ref{fig:BasinMulAdd}, where the size of the multiplicative noise is small. As the step size is annealed to zero, multiplicative noise resembles the convolutional nature of additive noise~\cite{kleinberg2018alternative}.
\begin{figure}
\centering
\includegraphics[width=0.95\textwidth]{basin_mul_vs_add.png}
\caption{\label{fig:BasinMulAdd}Histograms (blue) of $10^6$ iterations of GD with combinations of small, moderate, and strong versus light additive, heavy additive, and multiplicative noise, applied to a non-convex objective (black). (Red) Initial starting location for the optimization.}
\end{figure}
\begin{figure}
\centering
\includegraphics[width=\textwidth]{basin_sgd_vs_gdadd.png}
\caption{\label{fig:BasinSGD} Heuristic visualization.  PCA scores for trajectories of SGD over 10 epochs versus perturbed GD with additive noise across the same number of iterations along principle components $2,3,4$.} \end{figure}

This same basin hopping behaviour can be observed heuristically in SGD for higher-dimensional models with the aid of principal component analysis (PCA), although jumps become much more frequent. To see this, we consider fitting a two-layer neural network with 16 hidden units for classification of the Musk dataset \cite{dietterich1997solving} (168 attributes; 6598 instances) with cross-entropy loss without regularization and step size $\gamma = 10^{-2}$. Two stochastic optimizers are compared: (a) SGD with a single sample per batch (without replacement), and (b) perturbed GD \cite{jin2017escape}, where covariance of iterations in (b) is chosen to approximate that of (a) on average.
PCA-projected trajectories are presented in Figure~\ref{fig:BasinSGD}.
SGD frequently exhibits jumps between basins, a feature not shared by perturbed GD with only additive noise. Further numerical analysis is conducted in Appendix~\appref{B}{\ref{sec:Numerics2}} to support claims in \S\ref{sec:General}.

\section{Discussion}
\label{sec:Discussion}

\paragraph{Lazy training and implicit renewal theory.}
Theory concerning random linear recurrence relations can be extended directly to SGD applied to objective functions whose gradients are ``almost'' linear, such as those seen in lazy training \cite{chizat2019lazy}, due to the implicit renewal theory of Goldie \cite{goldie1991implicit}. However, it seems unlikely that multilayer networks should exhibit the same tail exponent between each layer, where the conditions of Breiman's lemma (see Appendix~\appref{A}{\ref{sec:LinearRecurrence}}) break down. Significant discrepancies between tail exponents across layers were observed in other appearances of heavy-tailed noise \cite{martin2020predicting,simsekli2019tail}.
From the point of view of tail exponents, it appears that most practical finite-width deep neural networks do not exhibit lazy training. This observation agrees with the poor relative generalization performance exhibited by linearized neural networks~\cite{chizat2019lazy,oymak2019generalization}.

\paragraph{Continuous-time models.}
Probabilistic analyses of stochastic optimization often consider continuous-time approximations, e.g., for constant step size, a stochastic differential equation of the form \cite{mandt2016variational, orvieto2019continuous,fontaine2020continuous}
\begin{equation}
\label{eq:LangevinModel}
\dd W_t = -\gamma \nabla f(W_t) \dd t + g(W_t) \dd X_t,
\end{equation}
where $f:\mathbb{R}^d \to \mathbb{R}$, $g:\mathbb{R}^{d} \to\mathbb{R}^{d\times d}$. The most common of these are Langevin models where $X$ is Brownian motion (white noise), although heavy-tailed noise has also been considered \cite{simsekli2019tail,simsekli2020fractional}. In the case where $g$ is diagonal, as a corollary of \cite{rubin2014mapping,ma2015complete}, we may determine conditions for heavy-tailed stationary behaviour. Letting $r(w) = \|g(w)^{-2} \nabla f(w)\|$, in Table~\ref{tab:Continuous} we present a continuous time analogue of Table~\ref{tab:Regimes} with examples of choices of $g$ that result in each regime when $f$ is quadratic.
Langevin algorithms commonly seen in the literature \cite{mandt2016variational,orvieto2019continuous} assume constant or bounded volatility ($g$), resulting in light-tailed stationary distributions for $L^2$-regularized loss functions. However, even a minute presence of multiplicative noise (unbounded $g$) can yield a heavy-tailed stationary distribution \cite{biro2005power}. Fortunately, more recent analyses are allowing for the presence of multiplicative noise in their approximations \cite{fontaine2020continuous}.
Based on our results, it seems clear that this trend is critical for adequate investigation of stochastic optimizers using the continuous-time approach.
On the other hand, \cite{simsekli2019tail,simsekli2020fractional} invoked the generalized central limit theorem (CLT) to propose a continuous-time model with heavy-tailed (additive) L\`{e}vy noise, in response to the observation that stochastic gradient noise is frequently heavy-tailed.In their model, the heaviness of the noise does not change throughout the optimization. Furthermore, \cite{panigrahi2019non} observed that stochastic gradient noise is typically Gaussian for large batch sizes, which is incompatible with generalized CLT. We have illustrated how multiplicative noise also yields heavy-tailed fluctuations. Alternatively, heavy tails in the stochastic gradient noise for \emph{fixed} weights could be log-normal instead, which are known to arise in deep neural networks \cite{hanin2019products}. A notable consequence of the L\`{e}vy noise model is that exit times for basin hopping are essentially independent of basin height, and dependent instead on basin \emph{width}, which is commonly observed with SGD \cite{xing2018walk}. In the absence of heavy-tailed additive noise, we observed similar behaviour in our experiments for large multiplicative noise, although precise theoretical treatment remains the subject of future work.

\paragraph{Stochastic gradient MCMC.}
The same principles discussed here also apply to stochastic gradient MCMC (SG-MCMC) algorithms \cite{ma2015complete}.
For example, the presence of multiplicative noise could suggest why SG-MCMC achieves superior mixing rates to classical Metropolis--Hastings MCMC in high dimensions, since the latter is known to struggle with multimodality, although there is some skepticism in this regard \cite{wenzel2020good}. Theorem~\ref{thm:Lipschitz} implies SG-MCMC yields excellent samples for satisfying the tail growth rate conditions required in importance sampling; see \cite{hodgkinson2020reproducing} for example.

\paragraph{Geometric properties of multiplicative noise.}
It has been suggested that increased noise in SGD acts as a form of convolution, smoothing the effective landscape \cite{kleinberg2018alternative}. This appears to be partially true. As seen in Figure~\ref{fig:BasinMulAdd}, smoothing behaviour is common for additive noise, and reduces resolution in the troughs. Multiplicative noise, which SGD also exhibits, has a different effect. As seen in \cite{rubin2014mapping}, in the continuous-time setting, multiplicative noise equates to conducting additive noise on a modified (via Lamperti transform) loss landscape. There is also a natural interpretation of multiplicative noise through choice of geometry in Riemannian Langevin diffusion \cite{girolami2011riemann,ma2015complete}. Under either interpretation, it appears multiplicative noise shrinks the width of peaks and widens troughs in the effective loss landscape, potentially negating some of the undesirable side effects of additive noise.

\paragraph{Heavy-tailed mechanistic universality.}
In a recent series of papers \cite{martin2018implicit, martin2019traditional, martin2020heavy, martin2020predicting}, an empirical (or phenomenological)
theory of \emph{heavy-tailed self-regularization} is developed, proposing that sufficiently complex and well-trained deep neural networks exhibit \emph{heavy-tailed mechanistic universality}: the spectral distribution of large weight matrices display a power law whose tail exponent is negatively correlated with generalization performance. If true, examination of this tail exponent provides an indication of model quality, and factors that are positively associated with improved test performance, such as decreased batch size. From random matrix theory, these heavy-tailed spectral distributions may arise either due to strong correlations arising in the weight matrices, or heavy-tailed distributions arising in each of the weights over time. The reality may be some combination of both; here, we have illustrated that the second avenue is at least possible, and relationships between factors of the optimization and the tail exponent agree with the present findings. Proving that heavy-tailed behaviour can arise in spectral distributions of the weights (as opposed to their probability distributions, which we have treated here) is a challenging problem, and remains the subject of future work. Casting the evolution of spectral distributions into the framework of iterated random functions could prove fruitful in this respect.

\section{Conclusion}
A theoretical analysis on the relationship between multiplicative noise and heavy-tailed fluctuations in stochastic optimization algorithms has been conducted. We propose that viewing stochastic optimizers through the lens of Markov process theory and examining stationary behaviour is key to understanding exploratory behaviour on non-convex landscapes in the initial phase of training. Our results suggest that heavy-tailed fluctuations may be more common than previous analyses have suggested, and they further the hypothesis that such fluctuations are correlated to improved generalization performance. From this viewpoint, we maintain that multiplicative noise should not be overlooked in future analyses on the subject.

\paragraph{Acknowledgments.} We would like to acknowledge DARPA, NSF, and ONR for providing partial support of this work.


TITLE: Differentio Integral Series
QUESTION [0 upvotes]: How to prove the following?
$$\frac{d^a{(x}^a)}{dx^a}+\frac{d^{a-1}{(x}^a)}{dx^{a-1}}+\frac{d^{a-2}{(x}^a)}{dx^{a-2}}+\dots +\frac{d^2{(x}^a)}{dx^2}+\frac{d^1{(x}^a)}{dx^1}+\frac{d^0{(x}^a)}{dx^0}+\int^x_0{{(x}^a)dx}+\int^x_0{\int^x_0{{(x}^a)dx}dx}+\int^x_0{\int^x_0{\int^x_0{{(x}^a)dx}}dxdx}+\dots =a!e^x $$
where, $a$ is a constant.

REPLY [0 votes]: $$a!\;e^x=a!\sum\limits_{n=0}^{\infty}\frac{x^n}{n!}=\sum\limits_{n=-a}^{\infty}\frac{a!\;x^{a+n}}{(a+n)!}=\sum\limits_{n=-a}^{\infty}\frac{d^{-n}(x^a)}{dx^{-n}}$$
$$=\sum\limits_{n=-a}^{0}\frac{d^{-n}(x^a)}{dx^{-n}}+\sum\limits_{n=1}^{\infty}\frac{d^{-n}(x^a)}{dx^{-n}}=\sum\limits_{n=0}^{a}\frac{d^{n}(x^a)}{dx^{n}}+\sum\limits_{n=1}^{\infty}\frac{d^{-n}(x^a)}{dx^{-n}}$$


TITLE: Perfectly flat surface
QUESTION [1 upvotes]: Looking at microscopic images of work tops I noticed that none are perfectly flat.
Is it possible to actually create a perfectly flat surface?

REPLY [2 votes]: It depends on what you mean by "surface" I suppose. If you mean general work surface then I don't know. I do know of at least one purportedly "perfectly flat" material though, due to the geometric arrangement of it's atoms:
https://en.wikipedia.org/wiki/Mica


TITLE: operator norm minimization problem
QUESTION [0 upvotes]: Given a matrix $\mathbf Z$ and $\mathbf G$, how to solve the following spectral minimization problem?
$$\min_{\mathbf m \in  \mathbb C^n} \|\mathbf Z-\hat{\mathbf Z}\|_2 \quad \text{such that} \quad \hat{\mathbf Z} = \mathbf G \operatorname{diag}(\mathbf m){\mathbf G}^*$$
Can anyone please help me with this problem? Thanks!

REPLY [3 votes]: There are some threads, like this and this, about related subjects. The suggested method is to relate $$\min_{w\in \mathbb{C}^n}\|Z-G\,diag(w) G^*\|_2,$$ when you know $Z$ and $G$, to the minimization problem  $$\left\{\begin{array}{rrr}
minimize& x\\subject\,to &\|A(w)\|_2\leq x
\end{array}\right.,$$ with $A(w)=Z-G\,diag(w) G^*$, and $x\geq 0$.
Since $\|A(w)\|_2$ is the maximum singular value $\sigma_{max}(A(w))$ of $A(w)$, it follows that $\|A(w)\|_2\leq x$ is equivalent to $$(\sigma_{max}(A(w)) )^2=\lambda_{max}(A(w)^*A(w))\leq x^2,$$ or $$x^2-\lambda_{max}(A(w)^*A(w))=\lambda_{min}(x^2I-A(w)^*A(w))\geq 0.$$ This means that the matrix $$x^2I-A(w)^*A(w)\geq 0, $$ that is, it is a positive semidefinite matrix. This implies that the block matirx $$\left[\begin{array}{rrr}
xI& A(w)\\A(w)^*&xI
\end{array}\right],$$ is positive semidefinite.
It follows that the initial minimizing problem becomes
$$\left\{\begin{array}{rrr}
 minimize& x\\\\subject\,to &\left[\begin{array}{rrr}
  xI& A(w)\\A(w)^*&xI
 \end{array}\right]\geq 0
\end{array}\right.,$$ with $x\geq 0$, $w\in \mathbb{C}^n$, and $A(w)=Z-G\,diag(w) G^*$. The last one leads you to Semidefinite programming, and you can see this closed thread about numerical calculations.


TITLE: Foliation with leaves which are and are not dense
QUESTION [7 upvotes]: Do there exist a foliation on a closed surface (i.e. real dimension 2) which has a dense leaf and also a leaf which is not dense?

REPLY [11 votes]: No, there doesn't exist such a foliation. The existence of any foliation would mean the Euler characteristic is zero, so the surface must be either a torus or a Klein bottle. Foliations for these surfaces are understood well enough to rule out having both dense and non-dense leaves.
Any foliation will contain a "Reeb component" (for which no leaf is dense) or will be a suspension of a homeomorphism of $S^1$. If a circle homeomorphism has one dense orbit, all orbits are dense, so the situation you describe can't exist.
If you allow surfaces with boundary, you need to consider a Mobius band or annulus, but it's still not possible to have a dense and non-dense leaf: the foliation will have pieces that are either Reeb components or suspensions of interval homeomorphisms, neither of which can have a dense leaf.
For the definition of a Reeb component and complete details, see this page, especially the "Foliations of Surfaces" section:
http://www.map.mpim-bonn.mpg.de/Foliations


TITLE: Approximate $\ln(2)$ from Maclaurin expansion
QUESTION [1 upvotes]: I was given a function f(x)= ln((1-x)/(1+x)) and I had to find the 9th degree polynomial and use it to approximate ln2.  I came up with the polynomial (x- (x^2)/2 + (x^3)/3 - (x^4)/4...., but I don't know what to do next to approximate ln2.

REPLY [1 votes]: To stitch the comments together, it seems that you have already found:
$$
\ln(1+x) = x - \frac{x^2}{2} + \frac{x^3}{3}  - \frac{x^4}{4} + \cdots
\tag{$*$}
$$
You were asked to find the Maclaurin series for
$$
    f(x) = \ln\left(\frac{1-x}{1+x}\right) = \ln(1-x) - \ln(1+x)
$$
but if you use (*) twice you should be able to put that together.  Now if you solve
$$
    \frac{1-x}{1+x} = \frac{1}{2}
$$
you should be able to find the right $x$ to plug in, and get a series that converges to $\ln\left(\frac{1}{2}\right) = -\ln 2$.


TITLE: Blum Blum Shub pseudorandom numbers - clarification on prerequisites
QUESTION [1 upvotes]: By reading on Wikipedia, they say the following (these is also supported by the following work: http://www.daimi.au.dk/~mg/mamian/random-bits.pdf)
The seed x0 should be an integer that is co-prime to M (i.e. p and q are not factors of x0) and not 1 or 0.
The two primes, p and q, should both be congruent to 3 (mod 4) (this guarantees that each quadratic residue has one square root which is also a quadratic residue) and gcd(φ(p − 1), φ(q − 1)) should be small (this makes the cycle length large).
Can someone clarify why the part in bold implies that each quadratic residue has one square root which is also a QR?

REPLY [1 votes]: Because $-1$ is not a quadratic residue modulo a prime $p\equiv 3 \mod 4$ and the product of two non-residues is a residue. A (non-zero) quadratic residue $r=d^2\mod p$ has the two square roots $\pm d$ which differ by a factor of $-1$. If $d$ is a non-residue, $-d$ is a residue.


TITLE: Prove that when d(x,a) and d(y,a) is equal for all $a \in A$ and A countable, dense then x = y
QUESTION [0 upvotes]: Let A be a dense countable subset of metric space $(X,d)$ and $x, y \in X$. Then the following is true:
$$
\left( \forall_{a \in A} \text{  } d(x,a) = d(y,a) \text{ } \right) \implies x = y
$$
How can I prove this?

REPLY [2 votes]: Let $(a_n)$ be a sequence in $A$ converging to $x$. Then $$0 = d(x,x) = \lim_{n\rightarrow \infty} d(x, a_n) = \lim_{n\rightarrow \infty} d(y,a_n) = d(y,x).$$
So $x=y$. The second equality is the continuity of $d$, and the second is your assumption. I don't see why you need countability here though; just $A$ being dense is enough.


TITLE: approximating sum of i.i.d poisson using CLT
QUESTION [0 upvotes]: We have $$ S = \sum _{k=1} ^{n} X_{k}$$ each $X_k$~$ Poisson(1)$. I want to find $$P\{S_{100} \ge200\}$$
I first approximate $S_{100}$ to standard normal.
$$P\{\frac{S_{100} - 100}{\sqrt{100}} \ge \frac{200 - 100}{\sqrt{1} \cdot \sqrt{100}}\}$$
Which result in $$ P\{Z \ge 10\}$$
Then I find it difficult to process it...

REPLY [0 votes]: Then I find it difficult to process it...

This probability is actually (EDIT) $\approx 0$. There is no need to use tables or calculators; just observe that the range of a Standard gaussian is for 99.99% in $[-3.8;3.8]$


TITLE: Finding the minimum volume of a Tetrahedron.
QUESTION [0 upvotes]: Suppose you have the surface $\xi$ defined in $\mathbb{R}^3$ by the equation:
$$ \xi :\frac{x^2}{a^2} + \frac{y^2}{b^2} + \frac{z^2}{c^2} = 1 $$
For $ x \geq 0$ , $ y \geq 0$ and $ z \geq 0$. Now take any point $P \in \xi$ and consider the tangent plane ($\pi_t)$ to $\xi$ at $P$. Calculate the minimum volume of the region determined by the $xy$, $yz$, $xz$ planes and $\pi_t$.

REPLY [1 votes]: Assume for the moment $a=b=c=1$. The surface $\xi$ then is the unit sphere $S^2$.  A typical point ${\bf n}$ of $S^2$ in the first octant then has coordinates $(n_1,n_2,n_3)$ with $n_i>0$ and $\sum_i n_i^2=1$. Furthermore the tangent plane to $S^2$ at ${\bf n}$ is given by the equation ${\bf n}\cdot{\bf x}=1$, i.e., $n_1x_1+n_2x_2+n_3x_3=1$. Intersecting this plane with the three coordinate axes gives the  points $\bigl({1\over n_1},0,0\bigr)$, $\bigl(0,{1\over n_2},0\bigr)$, $\bigl(0,0,{1\over n_3}\bigr)$. It follows that the simplex $S$ in question has volume
$${\rm vol}(S)={1\over 6 n_1n_2n_3}\ .$$
This volume is minimal when its reciprocal is maximal. By the AMG inequality
$$\root3\of {n_1^2 n_2^2n_3^2}\leq{1\over3}(n_1^2+n_2^2+n_3^2)={1\over3}\ ,$$
with equality iff $n_i=3^{-1/2}$ for all $i\in[3]$. It follows that
$${\rm vol}(S)\geq{1\over6} 3^{3/2}={\sqrt{3}\over2}\ .$$
In the case of arbitrary $a$, $b$, $c>0$ we therefore have
$${\rm vol}(S)\geq{\sqrt{3}\over2}abc\ ,$$
by standard geometric principles.


TITLE: Fourier transform of piecewise function
QUESTION [4 upvotes]: I am trying to calculate the Fourier transform of the following function:
$$f(x)=\begin{cases}e^{-x}&|x|<1\\0&\text{otherwise}\end{cases}$$
What I have so far is:
$$\hat{f}(\omega)=\int_{-\infty}^\infty f(x)e^{-j\omega x}dx=\int_{-1}^1e^{-x}e^{-j\omega x}dx=\int_{-1}^1e^{-(1+j\omega)x}dx$$
which I got as:
$$(1-j\omega)\frac{e^{(1+j\omega)}-e^{-(1+j\omega)}}{1+\omega^2}$$
which I have tried to simplify to the following:
$$2(1-j\omega)\frac{\sinh(1)\cos(\omega)+j\cosh(1)\sin(\omega)}{1+\omega^2}$$
which if I separate the real and imaginary parts left me with:
$$\frac{2\left[\sinh1\cos\omega+\cosh1\,\omega\sin\omega\right]+2j\left[\cosh1\sin\omega-\sinh1\,\omega\cos\omega\right]}{1+\omega^2}$$
My question is am I missing an identity that makes this a lot easier or is this as simple as I can get it? Thanks

REPLY [1 votes]: Here's one way to calculate the Fourier transform:
The distributional derivative of $f$ satisfies the equation
$$f'(x)=-f(x)+e^{1}\delta(x+1)-e^{-1}\delta(x-1).$$
Taking the Fourier transform of both sides gives
$$
j\omega\hat f(\omega) = -\hat f(\omega)+e^{1}e^{j\omega}-e^{-1}e^{-j\omega}
$$
i.e.
$$
\hat f(\omega) = \frac{e^{1+j\omega}-e^{-(1+j\omega)}}{1+j\omega}.
$$
I don't think that the splitting into real and imaginary parts can be done in a simple way.


TITLE: Proving f is continuous using epsilon and delta
QUESTION [1 upvotes]: I need to prove that f is continuous in (1,0):
$$\frac{(x-1)y^2}{(x-1)^2+y^2}-x \; for \;x \ne (1,0) \; and \; f(1,0) = -1$$
I came up with the following until now:
$$\left|\frac{(x-1)y^2}{(x-1)^2+y^2}-x+1\right|\leq \left|\frac{(x-1)y^2}{(x-1)^2+y^2}\right|+\left|(1-x) \right|$$
However, I do not know how to continue.
Thanks in advance for your help.

REPLY [0 votes]: hint
use the fact that
$$\Bigl|\frac{(x-1)y^2}{(x-1)^2+y^2}-x+1\Bigr|\le \frac{|x-1|y^2}{y^2}+|x-1|$$
or
$$|f(x,y)-f(1,0)|\le 2|x-1|$$


text:
\section{Analysis through Approximate Dynamical Isometry}
\subsection{Assumptions}\label{ssec:assumptions}
Let us review on how to achieve the dynamical isometry.
For the sake of the prospect of the theory, let $J$ be the Jacobian of the network with ignoring the last layer $\pre^L \mapsto \post^L = W_L \pre^L$.
We say that the network achieves dynamical isometry if all eigenvalues of $JJ^T$ are one.
Consider the limit spectral distribution $\nu_\ell$ of $D_\ell^2$ given by \eqref{align:dist-of-nu}.
At the deep limit, we consider the situation such that forward and backward signal propagation are stable. Hence  we adopt the following assumption,
\begin{assumption}\label{assumption:simplify}\citep{Pennington2018emergence}
For each $L \in \N$,  sequences $q_\ell, \alpha_\ell, \gamma_\ell$ and $\sigma_{\ell+1}$ ($\ell=0,1,\dots,L-1$) are constant for every $\ell$, but depend on $L$.  In addition, each constant converges to a finite value as $L\to\infty$.
\end{assumption}

Now the limit spectral distribution $\mu_{J^TJ}$ of  $J^\top J$ as $M \to \infty$ is given by $ [(\sigma_{L-1}^{2}\,\cdot\,)_* \nu_{L-1}]^{\boxtimes (L-1)}$. By \eqref{align:S-trans-decompose},
\begin{align}
S_{\mu_{J^\top J}}(z) = (\frac{1}{ \sigma_{L-1}^{2} \gamma_{L-1}}    (1 +  \frac{1 - \alpha_{L-1}}{z +\alpha_{L-1}}))^{L-1}.
\end{align}
In order to achieve dynamical isometry, we need that $[(\sigma_{L-1}^{2}\,\cdot\,)_* \nu_{L-1}]^{\boxtimes (L-1)}$ converges to a compactly supported distribution as $L \to \infty$, and $S_{JJ^T}(z)$ converges to a non-zero function of $z$.
Since the right hand side is approximated by $\exp (  L ( 1- \alpha_{L} )(z+\alpha_{L})^{-1}  - L\log \sigma_L^2\gamma_L) )$ as $L \to \infty$, we need
\begin{align}
 \log \sigma_L^2 \gamma_L  &=  O(L^{-1}) \label{align:sigma-gamma-one}\\
 1 -\alpha_L &= O(L^{-1}) \label{align:alpha-one}
\end{align}
 and  as $L \to \infty$.

Now, the exact dynamical isometry, which means that the input-output Jacobian is an orthogonal matrix, is too strong. Thus we consider the following approximate condition.
\begin{assumption}\label{assumption:eps}
We assume that the following limits exist and $|\eps_1|, |\eps_2| < 1$:
\begin{align}
\eps_1 &= \lim_{L\to\infty}L(1-\alpha_L),\\
\eps_2 &= -\lim_{L\to\infty}L \log\sigma_{L}^2\gamma_L.
\end{align}
\end{assumption}
To achieve exact dynamical isometry, we need  $\eps_1 = \eps_2 = 0$.  \cref{assumption:eps} is implicitly assumed in \cite{Pennington2018emergence}.

\subsection{Concentration around the Maximal Eigenvalue of the Conditional FIM}

We investigate $||\mu_L||_\infty$, where we denote by $||\nu||_\infty$ the maximum of the support of $\nu$ for a compactly supported probability distribution $\nu$ on $\R_+$.
We show that $||\mu_L||_\infty$ is $O(L)$ as $L \to \infty$.
We discuss the intuitive reason for this ahead of time.
Now, $\nu_\ell$  is close to a delta measure at a point to achieve dynamical isometry.
Then the recursive equation \eqref{align:reccurence_spec}  looks like an affine transform.
Therefore, $||\mu_L||_\infty$ is the result of $L$-times affine transforms, hence  it is $O(L)$ as $L \to \infty$.

A key of the proof is to show that $||\mu_L||_\infty$ is an atom of $\mu_L$ under \cref{assumption:eps}, which is an assumption to achieve approximate dynamical isometry.
Recall that $x \in \R$ is an atom of a probability distribution $\nu$ if and only if $\nu(\{x\}) > 0$.

Here we sketch the proof of the phenomenon.
Firstly, the support of $\nu_\ell$, which is the limit spectral distribution of $D_\ell^2$, consists of two atoms $0$ and $\gamma_\ell$, and the weight $\alpha_\ell$ of $\gamma_\ell$  is close to one for sufficiently large $\ell$.
Secondly, we recall the following fact of the free multiplicative convolution: if $a$ (resp.\,$b$) is an atom of a probability  distribution $\mu$ (resp.\,$\nu$) and $\mu(\{a\}) + \nu(\{b\}) - 1 >0$, then $ab$ is an atom of $\mu \boxtimes \nu$ and
\begin{align}\label{align:free-multi-conv}
    \mu \boxtimes \nu ( \{ab\}) = \mu(\{a\}) + \nu(\{b\}) - 1.
\end{align}
Thus, if $||\mu_\ell||_\infty$ is an atom of $\mu_\ell$ with sufficiently large weight, then $||\nu_\ell||_\infty || \mu_\ell||_\infty$ is an atom of $ \nu_\ell \boxtimes \mu_\ell$.
Note that \eqref{align:free-multi-conv} is different from that in the computation rule of the classical convolution.
Additionally, by the recursive equation \eqref{align:reccurence_spec}, we have the following recursive equation:
\begin{align}\label{align:recursive_max}
    ||\mu_{\ell+1}||_\infty = q_\ell  +  \sigma_{\ell+1}^2  ||\nu_\ell ||_\infty ||\mu_\ell||_\infty.
\end{align}
This recursive equation allows us to use induction.
Then we have our theorem.
\begin{thm}\label{thm:maximum}
Consider \cref{assumption:simplify}  and  \ref{assumption:eps}.
Then for sufficiently large $L$,  it holds that $||\mu_L||_\infty$  is an atom of $\mu_L$ with weight $1- (L-1)(1-\alpha_{L-1})$,  and
\begin{align}
     \lim_{L\to\infty}L^{-1}||\mu_L||_\infty  = q  \frac{  1 - \exp\left(-\eps_2\right) }{ \eps_2}.
\end{align}
In particular, the limit has the expansion $q(1 - \eps_2/2) + O(\eps_2^2)$ as the further limit $\eps_2 \to 0$.
\end{thm}
The proof is postpone to Appendix D.

Recall that $\mu_L$ is equal to the limit distribution of eigenvalues of the conditional FIM $I(\theta|x)$ except for $LM^2 -M$ zeros.
\cref{thm:maximum} shows that the maximum of the support is $O(L)$ as $L \to \infty$.
Furthermore, we emphasize that the weight $1 - (L-1)(1-\alpha_{L-1}) \sim 1 - \eps_1$ of the maximal eigenvalue $||\mu_{L}||$ is close to $1$.
Therefore, non-zero eigenvalues of the conditional FIM  asymptotically concentrates around $qL(1-\eps_2/2)$.
In particular,  the non-degenerate-part of the conditional FIM is approximated by a scaled identity operator.

config: arxiv


TITLE: Let $A$ be a subset of $\Bbb R$ with more than one element. Let $a \in A$. If $A\setminus\{a\}$ is compact. Then prove that $A$ is disconnected.
QUESTION [0 upvotes]: Let $A$ be a subset of $\Bbb R$ with more than one element. Let $a \in A$. If $A\setminus\{a\}$ is compact. Then prove that $A$ is disconnected.
My Attempt:
We know that if $A = B \cup C$ such that $B \neq \phi$,
$C \neq \phi$ with $\overline B\cap C$ = $\phi$, $B \cap \overline C$ = $\phi$ then $A$ is called disconnected. Here $A = A\setminus\{a\} \cup \{a\}$ such that $\overline A\setminus\{a\} \cap \{a\} = \phi$
, $A\setminus\{a\} \cap \overline{\{a\}}$ = $\phi$ $\implies$ $A\setminus\{a\} \cap \{a\} = \phi$, $A\setminus\{a\} \cap \{a\}$ = $\phi$ because
$A\setminus\{a\}$ and $\{a\}$ are closed sets so $\overline A\setminus\{a\} = A\setminus\{a\}$ and $\overline {\{a\}} = \{a\}$ Am I right? Please provide me alternate proof. Thanks.

REPLY [1 votes]: As $A\backslash\{a\}$ is compact in $\mathbb{R}^n$, you know that $A\backslash \{a\}$ is closed in $\mathbb{R}^n$ and is hence closed in $A$ in the subspace topology. Also $A\backslash\{a\}$ is open in $A$ as singleton sets are closed in the subspace topology. Thus you find that $A\backslash\{a\}$ is an open and closed proper subset of $A$. Thus $A=(A\backslash\{a\})\cup \{a\}$ is a separation of $A$ into disjoint non-empty open and closed sets. Thus $A$ is disconnected.
Note: The definition of disconnectedness you gave is equivalent to the condition $A=B\cup C$ where $B,C$ are non-empty open subsets of $A$ and $B\cap C\neq \emptyset$. You can try to prove it.