text:
\begin{document}
\IEEEoverridecommandlockouts
\title{A Primal--Dual SGD Algorithm  for\\ Distributed Nonconvex Optimization}
\author{Xinlei Yi, Shengjun Zhang, Tao Yang, Tianyou Chai, and Karl H. Johansson
\thanks{This work was supported by the
Knut and Alice Wallenberg Foundation, the  Swedish Foundation for Strategic Research, the Swedish Research Council, and the National Natural Science Foundation of China under grants 61991403, 61991404, and 61991400.}
\thanks{X. Yi and K. H. Johansson are with the Division of Decision and Control Systems, School of Electrical Engineering and Computer Science, KTH Royal Institute of Technology, 100 44, Stockholm, Sweden. {\tt\small \{xinleiy, kallej\}@kth.se}.}
\thanks{S. Zhang is with the Department of Electrical Engineering, University of North Texas, Denton, TX 76203 USA. {\tt\small  ShengjunZhang@my.unt.edu}.}
\thanks{T. Yang and T. Chai are with the State Key Laboratory of Synthetical Automation for Process Industries, Northeastern University, 110819, Shenyang, China. {\tt\small \{yangtao,tychai\}@mail.neu.edu.cn}.}
}

\maketitle

\begin{abstract}
The distributed nonconvex optimization problem of minimizing a global cost function formed by a sum of $n$ local cost functions by using local information exchange is considered.
This problem is an important component of many machine learning techniques with data parallelism, such as deep learning and federated learning.
We propose a distributed primal--dual stochastic gradient descent (SGD) algorithm,  suitable for arbitrarily connected communication networks and any smooth (possibly nonconvex) cost functions.
We show that the proposed algorithm achieves the linear speedup convergence rate $\mathcal{O}(1/\sqrt{nT})$ for general nonconvex cost functions and the linear speedup convergence rate $\mathcal{O}(1/(nT))$ when the global cost function satisfies the Polyak--{\L}ojasiewicz (P--{\L}) condition, where $T$ is the total number of iterations.
We also show that the output of the proposed algorithm with constant parameters linearly converges to a neighborhood of a global optimum.
We demonstrate through numerical experiments the efficiency of our algorithm in comparison with the baseline centralized SGD and recently proposed distributed SGD algorithms.

\emph{Index Terms}---Distributed nonconvex optimization, linear speedup, Polyak--{\L}ojasiewicz condition, primal--dual algorithm, stochastic gradient descent
\end{abstract}

\section{Introduction}
Consider a network of $n$ agents, each of which has a local smooth (possibly nonconvex) cost function $f_i: \mathbb{R}^{p}\rightarrow \mathbb{R}$.
All agents collaboratively solve the following optimization problem
\begin{align}\label{sgd:eqn:xopt}
 \min_{x\in \mathbb{R}^p} f(x)=\frac{1}{n}\sum_{i=1}^nf_i(x).
\end{align}
Each agent $i$ only has information about its local cost function $f_i$ and can communicate with its neighbors through the underlying communication network. The communication network is modeled by an undirected graph $\mathcal G=(\mathcal V,\mathcal E)$, where $\mathcal V =\{1,\dots,n\}$ is the agent set, $\mathcal E
\subseteq \mathcal V \times \mathcal V$ is the edge set, and $(i,j)\in \mathcal E$ if agents $i$ and $j$ can communicate with each other. The set $\mathcal{N}_i=\{j\in \mathcal V :~ (i,j)\in \mathcal E\}$ is the neighboring set of agent $i$.
The optimization problem \eqref{sgd:eqn:xopt} incorporates many popular machine learning approaches with data parallelism, such as deep learning \cite{dean2012large} and federated learning \cite{McMahan2017communication}.
A star graph is a special undirected graph, in which there is one and only one agent  (hub agent) that connects to all of the other agents (leaf agents) and each leaf agent only connects to the hub agent. Such a graph corresponds to the master/worker architecture adopted by many parallel learning algorithms.

In this paper, we consider the case where each agent is able to collect stochastic gradients of its local cost function and propose a distributed stochastic gradient descent (SGD) algorithm to solve  \eqref{sgd:eqn:xopt}.
In general, SGD algorithms are suitable for scenarios where explicit expressions of the gradients are unavailable or difficult to obtain.
For example, in some big data applications, such as empirical risk minimization, the actual gradient is to be calculated from the entire data set, which results in a heavy computational burden.
A stochastic gradient can be calculated from a randomly selected subset of the data and is often an efficient way to replace the actual gradient.
Other examples when SGD algorithms are suitable include scenarios where data are arriving sequentially such as in online learning \cite{langford2009sparse}.

\subsection{Literature Review}
When the communication network is a star graph, various parallel SGD algorithms have been proposed to solve \eqref{sgd:eqn:xopt}.
A potential performance bottleneck of such algorithms lies on the communication burden of the master.
To overcome this issue, a promising strand of research is combining parallel SGD algorithms with communication reduction approaches, e.g., asynchronous parallel SGD algorithms \cite{recht2011hogwild,de2015taming,lian2015asynchronous,lian2016Comprehensive,Zhou2018distributedas},
gradient compression based parallel SGD algorithms \cite{de2015taming,pmlr-v80-bernstein18a,
jiang2018linear,reisizadeh2019fedpaq,basu2019qsparse},
periodic averaging based parallel SGD algorithms \cite{jiang2018linear,reisizadeh2019fedpaq,wang2018adaptive,yu2019parallel,
haddadpour2019trading,Yu2019on,haddadpour2019local}, and parallel SGD algorithm with dynamic batch sizes \cite{Yu2019Computation}.
Convergence properties of these algorithms have been analyzed in detail.
In particular, in \cite{jiang2018linear,yu2019parallel,Yu2019on,Yu2019Computation}, an $\mathcal{O}(1/\sqrt{nT})$ convergence rate has been established for general nonconvex cost functions, where $T$ is the total number of iterations.
This rate is $n$ times faster than the well known $\mathcal{O}(1/\sqrt{T})$ convergence rate established by SGD over a single agent, and thus a linear speedup in the number of agents is achieved.
In \cite{haddadpour2019local,Yu2019Computation}, the  convergence rate has been improved to $\mathcal{O}(1/(nT))$ when the global cost function satisfies the P--{\L} condition, which also achieves a linear speedup.
In addition to the star architecture restriction, aforementioned parallel SGD algorithms require certain restrictions on the cost functions, such as bounded gradients of the local cost functions or bounded difference between the gradients of the local and global cost functions.

Distributed algorithms executed over arbitrarily connected communication networks have been suggested to overcome communication bottlenecks for parallel SGD algorithms.
Various distributed SGD algorithms have been proposed to solve \eqref{sgd:eqn:xopt}, e.g., synchronous distributed SGD algorithms \cite{Yu2019on,jiang2017collaborative,lian2017can,george2019distributed},
asynchronous distributed SGD algorithms \cite{pmlr-v80-lian18a,Assran2019Stochastic},
compression based distributed SGD algorithms \cite{tang2018communication,reisizadeh2019robust,taheri2020quantized,Singh2020Communication},
and periodic averaging based distributed SGD algorithm \cite{wang2018cooperative}.
Convergence properties of these algorithms have been analyzed and the linear speedup convergence rate $\mathcal{O}(1/\sqrt{nT})$ has been established for general nonconvex cost functions \cite{Yu2019on,lian2017can,Assran2019Stochastic,tang2018communication,taheri2020quantized,
Singh2020Communication,wang2018cooperative}.
However, similar to aforementioned parallel SGD algorithms, these distributed algorithms require  restrictive assumptions on the cost functions.
In order to remove these restrictions, the authors of \cite{Tang2018Decentralized} proposed a variant of the distributed SGD algorithm proposed in \cite{lian2017can}, named $\mathrm{D}^2$, in which each agent stores the stochastic gradient and its local model in last iteration and linearly combines them with the current stochastic gradient and local model. For this algorithm the authors established the linear speedup convergence rate $\mathcal{O}(1/\sqrt{nT})$, but they required that the eigenvalues of the mixing matrix associated with the communication network  are strictly greater than $-1/3$.
The authors of \cite{lu2019gnsd,zhang2019decentralized} proposed  distributed stochastic gradient tracking algorithms suitable for arbitrarily connected communication networks.
However, these algorithms only achieve $\mathcal{O}(1/\sqrt{T})$ convergence rate, which is not a speedup.
Moreover, gradient tracking algorithms have the common potential drawback that in order to track the global gradient,
at each iteration each agent needs to communicate one additional $p$-dimensional variable with its neighbors.
This results in heavy communication  burden when $p$ is large.

Note that all aforementioned distributed SGD algorithms converge to stationary points, which may be local or global optima, or saddle points.
None of existing studies on distributed SGD algorithms consider finding the global optimum when the global cost function satisfies some additional property, such as the P--{\L} condition studied for the parallel algorithms in  \cite{haddadpour2019local,Yu2019Computation}. Noting above, two core theoretical questions arise.

(Q1) Are there any distributed SGD algorithms such that they not only are suitable for arbitrarily connected communication networks and any smooth cost functions but also find stations points with the linear speedup convergence rate $\mathcal{O}(1/\sqrt{nT})$?

(Q2) If the P--{\L} condition holds in addition, can the above SGD algorithms find the global optimum with the linear speedup convergence rate $\mathcal{O}(1/(nT))$?

\subsection{Main Contributions}
This paper provides positive answers for the above two questions. More specifically, the contributions of this paper are summarized as follows.

(i) We propose a distributed primal--dual SGD algorithm to solve the optimization problem~\eqref{sgd:eqn:xopt}.
In the proposed algorithm, each agent maintains the primal and dual variable sequences and only communicates the primal variable with its neighbors.
This algorithm is suitable for arbitrarily connected communication networks and any smooth (possibly nonconvex) cost functions.

(ii) We show in Corollary~\ref{sguT:coro-sg-smT} that our algorithm finds a stationary point with the linear speedup convergence rate $\mathcal{O}(1/\sqrt{nT})$ for general nonconvex cost functions, thus answers (Q1). Compared with \cite{jiang2018linear,yu2019parallel,Yu2019on,Yu2019Computation,lian2017can,
Assran2019Stochastic,tang2018communication,taheri2020quantized,Singh2020Communication,
wang2018cooperative,Tang2018Decentralized}, we achieve the same convergence rate but under weaker assumptions related to network architectures and/or cost functions, and compared with \cite{lu2019gnsd,zhang2019decentralized}, we not only establish linear speedup but also just use half communication in each iteration.

(iii) We show in Theorem~\ref{sgu:thm-sg-diminishingt} that our algorithm finds a global optimum with the linear speedup convergence rate $\mathcal{O}(1/(nT))$   when the global cost function satisfies the P--{\L} condition,  thus answers (Q2).  Compared with \cite{haddadpour2019local,Yu2019Computation,stich2018local,Koloskova2019decentralized,
Singh2020Communication,olshevsky2019non} , we achieve the same convergence rate but under weaker assumptions related to network architectures and/or cost functions, and compared with \cite{reisizadeh2019fedpaq,jiang2017collaborative,rabbat2015multi,lan2018communication,
yuan2018optimal,jakovetic2018convergence,fallah2019robust}, we not only establish linear speedup but also relax the strong convexity by the P--{\L} condition.

(iv) We show in Theorems~\ref{sg:thm-sg-fixed-a5} and \ref{sg:thm-sg-fixed} that the output of our algorithm with constant parameters linearly converges to a neighborhood of a global optimum when the global cost function satisfies the P--{\L} condition.
Compared with \cite{pu2018swarming,jiang2017collaborative,fallah2019robust,pu2018distributed,xin2019distributed}, which used the strong convexity assumption, we achieve the similar convergence result under weaker assumptions on the cost function.

The comparison of this paper with other related studies in the literature is summarized in Table~\ref{sgd:table}.
\begin{table*}[htbp]
\caption{Comparison of this paper to some related works. }
\label{sgd:table}
\vskip -1in
\begin{center}
\begin{tiny}

\scalebox{1}{
\begin{tabular}{M{0.9cm}|M{1.0cm}|M{1.9cm}|M{1.9cm}|M{2.2cm}|M{1.8cm}|M{4.8cm}N}
\hline
Reference&Problem type&Extra assumption&Communication network&Communicated variable&Communication rounds&Convergence rate&\\[7pt]

\hline
\cite{jiang2018linear}&Nonconvex&Bounded $\|\nabla f_i-\nabla f\|$&Star graph&One quantized variable &$\mathcal{O}(n^{5/4}T^{3/4})$&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]

\hline
\multirow{2}{*}[-4.5pt]{\cite{reisizadeh2019fedpaq}}& Nonconvex& \multirow{2}{*}[-4.5pt]{Identical $\nabla f_i$}
& \multirow{2}{*}[-4.5pt]{Star graph}& \multirow{2}{*}[-4.5pt]{\parbox{2.2cm}{\centering One quantized variable}}& \multirow{2}{*}[-4.5pt]{$\mathcal{O}(T)$}&$\mathcal{O}(1/\sqrt{T})$&\\[9pt]
\cline{2-2}\cline{7-7}
&Strongly convex& & & & &$\mathcal{O}(1/T)$&\\[9pt]

\hline
\cite{yu2019parallel}&Nonconvex&Bounded $\|\nabla f_i\|$&Star graph&One full-information variable &$\mathcal{O}(n^{3/4}T^{3/4})$&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]

\hline
\multirow{2}{*}[-4.5pt]{\cite{Yu2019on}}&  \multirow{2}{*}[-4.5pt]{Nonconvex}& \multirow{2}{*}[-4.5pt]{\parbox{1.9cm}{\centering Bounded $\|\nabla f_i-\nabla f\|$}}
&Star graph& \multirow{2}{*}[-4.5pt]{\parbox{2.2cm}{\centering Two full-information variables}}& $\mathcal{O}(n^{3/4}T^{3/4})$&\multirow{2}{*}[-4.5pt]{$\mathcal{O}(1/\sqrt{nT})$}&\\[9pt]
\cline{4-4}\cline{6-6}
& & &Connected graph & & $\mathcal{O}(T)$&&\\[9pt]

\hline
\cite{haddadpour2019local}&P--{\L} condition&Identical $\nabla f_i$&Star graph&One full-information variable &$\mathcal{O}((nT)^{1/3})$&$\mathcal{O}(1/(nT))$&\\[9pt]

\hline
\multirow{2}{*}[-4.5pt]{\cite{Yu2019Computation}}&Nonconvex&\multirow{2}{*}[1.5pt]{\parbox{1.9cm}{\centering Identical $\nabla f_i$, exponentially increasing batch size}}&\multirow{2}{*}[-4.5pt]{Star graph}&\multirow{2}{*}[-4.5pt]{\parbox{2.2cm}{\centering One full-information variable}} &$\mathcal{O}(\sqrt{nT}\log(\frac{T}{n}))$&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]
\cline{2-2}\cline{6-7}
&P--{\L} condition& & & &$\mathcal{O}(\log(T))$&$\mathcal{O}(1/(nT))$&\\[9pt]

\hline
\multirow{2}{*}[-4.5pt]{\cite{jiang2017collaborative}}&Nonconvex&\multirow{2}{*}[-4.5pt]{Bounded $\|\nabla f_i\|$}&\multirow{2}{*}[-4.5pt]{Connected graph}&\multirow{2}{*}[-4.5pt]{\parbox{2.2cm}{\centering One full-information variable}} &\multirow{2}{*}[-4.5pt]{$\mathcal{O}(T)$}&$\mathcal{O}(1/T^\theta),~\forall \theta\in(0,0.5)$&\\[9pt]
\cline{2-2}\cline{7-7}
&Strongly convex& & & &&$\mathcal{O}(1/T)$; linearly to a neighbor of the global optimum (constant stepsize)&\\[9pt]

\hline
\cite{lian2017can}&Nonconvex&Bounded $\|\nabla f_i-\nabla f\|$&Connected graph&One full-information variable &$\mathcal{O}(T)$&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]

\hline
\cite{Assran2019Stochastic}&Nonconvex&Bounded $\|\nabla f_i-\nabla f\|$&Uniformly jointly
strongly connected digraph&One full-information variable &$\mathcal{O}(T)$&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]

\hline
\cite{tang2018communication}&Nonconvex&Bounded $\|\nabla f_i-\nabla f\|$&Connected graph&One compressed variable &$\mathcal{O}(T)$&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]

\hline
\cite{taheri2020quantized}&Nonconvex&Bounded $\|\nabla f_i\|$&Strongly connected digraph&One quantized variable &$\mathcal{O}(T)$&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]

\hline
\multirow{2}{*}[-4.5pt]{\cite{Singh2020Communication}}&Nonconvex&\multirow{2}{*}[-4.5pt]{Bounded $\|\nabla f_i\|$}&\multirow{2}{*}[-4.5pt]{Connected graph}&\multirow{2}{*}[-4.5pt]{\parbox{2.2cm}{\centering One compressed variable}} &\multirow{2}{*}[-4.5pt]{Event-triggered}&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]
\cline{2-2}\cline{7-7}
&Strongly convex& & & &&$\mathcal{O}(1/(nT))$&\\[9pt]

\hline
\cite{wang2018cooperative}&Nonconvex&Identical $\nabla f_i$&Connected graph&One full-information variable &$\mathcal{O}(n^{3/2}\sqrt{T})$&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]

\hline
\cite{Tang2018Decentralized}&Nonconvex&The eigenvalues of the mixing matrix are strictly greater than $-1/3$&Connected graph&One full-information variable &$\mathcal{O}(T)$&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]

\hline
\cite{lu2019gnsd,zhang2019decentralized}&Nonconvex&No&Connected graph&Two full-information variables &$\mathcal{O}(T)$&$\mathcal{O}(1/\sqrt{T})$&\\[9pt]

\hline
\cite{stich2018local}&Strongly convex&Bounded $\|\nabla f_i\|$&Star graph&One full-information variable &$\mathcal{O}(\sqrt{T/n})$&$\mathcal{O}(1/(nT))$&\\[9pt]

\hline
\cite{Koloskova2019decentralized}&Strongly convex&Bounded $\|\nabla f_i\|$&Connected graph&One compressed variable &$\mathcal{O}(T)$&$\mathcal{O}(1/(nT))$&\\[9pt]

\hline
\cite{olshevsky2019non}&Strongly convex&No&Connected graph&Two full-information variables &$\mathcal{O}(T)$&$\mathcal{O}(1/(nT))$&\\[9pt]

\hline
\cite{rabbat2015multi}&Strongly convex&Identical $\nabla f_i$&Connected graph&One full-information variable &$\mathcal{O}(T)$&$\mathcal{O}(1/T)$&\\[9pt]

\hline
\cite{lan2018communication}&Strongly convex&No&Connected graph&One full-information variable &$\mathcal{O}(\sqrt{T})$&$\mathcal{O}(1/T)$&\\[9pt]

\hline
\cite{yuan2018optimal}&Strongly convex&Bounded $\|\nabla f_i\|$&Uniformly jointly
strongly connected digraph&One full-information variable &$\mathcal{O}(T)$&$\mathcal{O}(1/T)$&\\[9pt]

\hline
\cite{jakovetic2018convergence}&Strongly convex&No&Connected graph in expectation&One full-information variable &$\mathcal{O}(T)$&$\mathcal{O}(1/T)$&\\[9pt]

\hline
\cite{fallah2019robust}&Strongly convex&No&Connected graph&One full-information variable &$\mathcal{O}(T)$&$\mathcal{O}(1/T)$; linearly to a neighbor of the global optimum (constant stepsize)&\\[9pt]

\hline
\cite{pu2018swarming}&Strongly convex&No&Connected graph&One full-information variable &$\mathcal{O}(T)$&Linearly to a neighbor of the global optimum (constant stepsize)&\\[9pt]

\hline
\cite{pu2018distributed}&Strongly convex&No&Connected graph&Two full-information variables &$\mathcal{O}(T)$&Linearly to a neighbor of the global optimum (constant stepsize)&\\[9pt]

\hline
\cite{xin2019distributed}&Strongly convex&No&Strongly connected digraph&Two full-information variables &$\mathcal{O}(T)$&Linearly to a neighbor of the global optimum (constant stepsize)&\\[9pt]

\hline
\multirow{3}{*}[-12pt]{\parbox{0.9cm}{\centering This paper}}&Nonconvex&\multirow{2}{*}[-14pt]{No}&\multirow{3}{*}[-14pt]{Connected graph}&\multirow{3}{*}[-14pt]{\parbox{2.2cm}{\centering One full-information variable}} &\multirow{3}{*}[-14pt]{$\mathcal{O}(T)$}&$\mathcal{O}(1/\sqrt{nT})$&\\[9pt]
\cline{2-2}\cline{7-7}
&\multirow{2}{*}[-9pt]{\parbox{1cm}{\centering P--{\L} condition}}&& & &&$~~~~~~~~~$$\mathcal{O}(1/(T^\theta)),~\forall \theta\in(0,1)$;$~~~~~~~~~$ linearly to a neighbor of the global optimum (constant stepsize)&\\[9pt]
\cline{3-3}\cline{7-7}
&&Bounded $f^*_i$ & & &&$\mathcal{O}(1/(nT))$&\\[9pt]
\hline
\end{tabular}
}

\end{tiny}
\end{center}
\vskip -0.1in
\end{table*}

\subsection{Outline}
The rest of this paper is organized as follows. Section~\ref{sguT:sec-main} presents the new distributed primal--dual SGD algorithm. Section~\ref{sguT:sec-analysis} analyzes its convergence rate. Numerical experiments are given in Section~\ref{sgd:sec-simulation}. Finally, concluding remarks are offered in Section~\ref{sgd:sec-conclusion}.
To improve the readability, all the proofs are given in the appendix.

\noindent {\bf Notations}: $\mathbb{N}_0$ and $\mathbb{N}_+$ denote the set of nonnegative and positive integers, respectively. $[n]$ denotes the set $\{1,\dots,n\}$ for any $n\in\mathbb{N}_+$. ${\bm 1}_n$ (${\bm 0}_n$) denotes the column one (zero) vector of dimension $n$. $\col(z_1,\dots,z_k)$ is the concatenated column vector of vectors $z_i\in\mathbb{R}^{p_i},~i\in[k]$.
$\|\cdot\|$ represents the Euclidean norm for vectors or the induced 2-norm for matrices. Given a differentiable function $f$, $\nabla f$ denotes the gradient of $f$.

\section{Distributed Primal--Dual SGD Algorithm}\label{sguT:sec-main}

In this section, we propose a new distributed SGD algorithm based on the primal--dual method.

Denote $\bsx=\col(x_1,\dots,x_n)$, $\tilde{f}(\bsx)=\sum_{i=1}^{n}f_i(x_i)$,  and $\bsL=L\otimes {\bm I}_p$, where $L=(L_{ij})$ is the weighted Laplacian matrix associated with the undirected communication graph $\calG$.
Recall that the Laplacian matrix $L$ is positive semi-definite and $\nullrank(L)=\{{\bm 1}_n\}$ when $\calG$ is connected \cite{mesbahi2010graph}. The optimization problem \eqref{sgd:eqn:xopt} is equivalent to the following constrained optimization problem:
\begin{mini}
{\bsx\in \mathbb{R}^{np}}{\tilde{f}(\bsx)}{\label{sgd:eqn:xoptcon}}{}
\addConstraint{\bsL^{1/2}\bsx=}{~{\bm 0}_{np}.}{}
\end{mini}

Let $\bsu\in\mathbb{R}^{np}$ denote the dual variable. Then the augmented Lagrangian function associated with \eqref{sgd:eqn:xoptcon} is
\begin{align}\label{nonconvex:lagran}
\calA(\bsx,\bsu)=\tilde{f}(\bsx)+\frac{\alpha}{2}\bsx^\top\bsL\bsx+\beta\bsu^\top\bsL^{1/2}\bsx,
\end{align}
where $\alpha>0$ and $\beta>0$ are parameters to be designed later.

Based on the primal--dual gradient method, a distributed SGD algorithm to solve \eqref{sgd:eqn:xoptcon} is
\begin{subequations}\label{sgu:kiau-algo-compact}
\begin{align}
\bm{x}_{k+1}&=\bm{x}_k-\eta_k(\alpha_k\bsL\bm{x}_k+\beta_k\bsL^{1/2}\bm{u}_k+\bsg_k^u),\\
\bm{u}_{k+1}&=\bm{u}_{k}+\eta_k\beta_k\bsL^{1/2}\bm{x}_k,~\forall \bsx_0,~\bsu_0\in\mathbb{R}^{np},
\end{align}
\end{subequations}
where $\eta_k>0$ is the stepsize at iteration $k$, $\alpha_k>0$ and $\beta_k>0$ are the values of the parameters $\alpha$ and $\beta$ at iteration $k$, respectively, and $\bsg_k^u=\col(g^u_{1,k},\dots,g^u_{n,k})$ with $g^u_{i,k}=g_{i}(x_{i,k},\xi_{i,k})$ being the stochastic gradient of $f_i$ at $x_{i,k}$ and $\xi_{i,k}$ being a random variable.
Denote $\bsv_k=\col(v_{1,k},\dots,v_{n,k})=\bsL^{1/2}\bm{u}_k$. Then the recursion \eqref{sgu:kiau-algo-compact} can be rewritten as
\begin{subequations}\label{sgu:kia-algo-dc-compact}
\begin{align}
\bm{x}_{k+1}&=\bm{x}_k-\eta_k(\alpha_k\bsL\bm{x}_k+\beta_k\bm{v}_k+\bsg_k^u),
\label{sgu:kia-algo-dc-compact-x}\\
\bm{v}_{k+1}&=\bm{v}_k+\eta_k\beta_k\bsL\bm{x}_k,~\forall \bsx_0\in\mathbb{R}^{np},~\sum_{j=1}^nv_{j,0}={\bm 0}_p.\label{sgu:kia-algo-dc-compact-v}
\end{align}
\end{subequations}
The initialization condition $\sum_{j=1}^nv_{j,0}={\bm 0}_p$ is derived from $\bsv_0=\bsL^{1/2}\bm{u}_0$, and it is easy to be satisfied, for example, $v_{i,0}={\bm 0}_p,~\forall i\in[n]$, or $v_{i,0}=\sum_{j=1}^{n}L_{ij}x_{j,0},~\forall i\in[n]$.
Note that \eqref{sgu:kia-algo-dc-compact} can  be written agent-wise as
\begin{subequations}\label{sguT:alg:stochastic}
\begin{align}
x_{i,k+1} &= x_{i,k}-\eta_k\Big(\alpha_k\sum_{j=1}^nL_{ij}x_{j,k}+\beta_k v_{i,k}+g^u_{i,k}\Big), \label{sguT:alg:stochastic-x}\\
v_{i,k+1} &=v_{i,k}+ \eta_k\beta_k\sum_{j=1}^n L_{ij}x_{j,k},~\forall x_{i,0}\in\mathbb{R}^p, ~v_{i,0}={\bm 0}_p,~
\forall i\in[n].  \label{sguT:alg:stochastic-q}
\end{align}
\end{subequations}
This corresponds to our proposed distributed primal--dual SGD algorithm, which is presented in pseudo-code as Algorithm~\ref{sguT:algorithm-sg}.
\begin{algorithm}[tb]
\caption{Distributed Primal--Dual SGD Algorithm}
\label{sguT:algorithm-sg}
\begin{algorithmic}[1]
\STATE \textbf{Input}: parameters $\{\alpha_k\}$, $\{\beta_k\}$, $\{\eta_k\}\subseteq(0,+\infty)$.
\STATE \textbf{Initialize}: $ x_{i,0}\in\mathbb{R}^p$ and $v_{i,0}={\bm 0}_p,~
\forall i\in[n]$.
\FOR{$k=0,1,\dots$}
\FOR{$i=1,\dots,n$  in parallel}
\STATE  Broadcast $x_{i,k}$ to $\mathcal{N}_i$ and receive $x_{j,k}$ from $j\in\mathcal{N}_i$;
\STATE  Sample stochastic gradient $g_{i}(x_{i,k},\xi_{i,k})$;
\STATE  Update $x_{i,k+1}$ by \eqref{sguT:alg:stochastic-x};
\STATE  Update $v_{i,k+1}$ by \eqref{sguT:alg:stochastic-q}.
\ENDFOR
\ENDFOR
\STATE  \textbf{Output}: $\{\bsx_k\}$.
\end{algorithmic}
\end{algorithm}

It should be pointed out that $\{\alpha_k\}$, $\{\beta_k\}$, $\{\eta_k\}$, $\bsx_0$, $\bsv_0$, and $\bsv_1$  in Algorithm~\ref{sguT:algorithm-sg} are deterministic, while $\{\bsx_k\}_{k\ge1}$ and $\{\bsv_k\}_{k\ge2}$ are random variables generated by Algorithm~\ref{sguT:algorithm-sg}. Let $\mathfrak{F}_k$ denote the $\sigma$-algebra generated by the random variables $\xi_{1,k},\dots,\xi_{n,k}$ and let $\mathcal{F}_k=\bigcup_{s=1}^{k}\mathfrak{F}_s$. It is straightforward to see that $\bsx_k$ and $\bsv_{k+1}$ depend on $\mathcal{F}_{k-1}$ and are independent of $\mathfrak{F}_s$ for all $s\ge k$.

\section{Convergence Rate Analysis}\label{sguT:sec-analysis}
In this section, we analyze the convergence rate of Algorithm~\ref{sguT:algorithm-sg}.
The following assumptions are made.

\begin{assumption}\label{sgd:assgraph}
The undirected communication graph $\mathcal G$ is connected.
\end{assumption}

\begin{assumption}\label{sgd:assoptset}
The minimum function value of the optimization problem \eqref{sgd:eqn:xopt} is finite.
\end{assumption}

\begin{assumption}\label{sgd:assfiu}
Each  local cost function $f_i$ is smooth with constant $L_{f}>0$, i.e.,
\begin{align}\label{nonconvex:smooth}
\|\nabla f_i(x)-\nabla f_i(y)\|\le L_{f}\|x-y\|,~\forall x,y\in \mathbb{R}^p.
\end{align}
\end{assumption}

\begin{assumption}\label{sgd:ass:stochastic-grad:xi}
The random variables $\{\xi_{i,k},~i\in[n],~k\in\mathbb{N}_0\}$ are independent of each other.
\end{assumption}

\begin{assumption}\label{sgd:ass:stochastic-grad:mean}
The stochastic estimate $g_{i}(x,\xi_{i,k})$ is unbiased, i.e., for all $i\in[n]$, $k\in\mathbb{N}_0$, and $x\in\mathbb{R}^p$,
\begin{align}\label{sgd:ass:stochastic-grad:mean-equ}
\mathbf{E}_{\xi_{i,k}}[g_{i}(x,\xi_{i,k})]=\nabla f_{i}(x).
\end{align}
\end{assumption}

\begin{assumption}\label{sgd:ass:stochastic-grad:variance}
The stochastic estimate $g_{i}(x,\xi_{i,k})$ has bounded variance, i.e., there exists a constant $\sigma$ such that for all $i\in[n]$, $k\in\mathbb{N}_0$, and $x\in\mathbb{R}^p$,
\begin{align}\label{sgd:ass:stochastic-grad:variance-equ}
\mathbf{E}_{\xi_{i,k}}[\|g_{i}(x,\xi_{i,k})-\nabla f_{i}(x)\|^2]\le\sigma^2.
\end{align}
\end{assumption}

\begin{remark}
The bounded variance assumption (Assumption~\ref{sgd:ass:stochastic-grad:variance}) is weaker than the bounded second moment (or bounded gradient) assumption made in \cite{recht2011hogwild,de2015taming,lian2015asynchronous,Zhou2018distributedas,basu2019qsparse,
yu2019parallel,jiang2017collaborative,george2019distributed,taheri2020quantized,yuan2018optimal,
rakhlin2011making,stich2018local,Koloskova2019decentralized,Singh2020Communication}. Moreover, note that we make no assumption on the boundedness of the deviation between the gradients of local cost functions. In other words, we do not assume that $\frac{1}{n}\sum_{i=1}^{n}\|\nabla f_i(x)-\nabla f(x)\|^2$ is uniformly bounded, which is commonly done in studies of deep learning, e.g., \cite{jiang2018linear,yu2019parallel,Yu2019on,lian2017can,pmlr-v80-lian18a,Assran2019Stochastic,tang2018communication}.
Also, we do not assume that the mean of each local stochastic gradient is the gradient of the global cost function, i.e.,  $\mathbf{E}_{\xi}[g_{i}(x,\xi)]=\nabla f(x),~\forall x\in\mathbb{R}^p,~\forall i\in[n]$,  which is commonly assumed in studies of empirical risk minimization and stochastic optimization, e.g., \cite{lian2016Comprehensive,pmlr-v80-bernstein18a,reisizadeh2019fedpaq,wang2018adaptive,haddadpour2019trading,
haddadpour2019local,Yu2019Computation,reisizadeh2019robust,wang2018cooperative,rabbat2015multi}.
\end{remark}

\subsection{Find Stationary Points}

Let us consider the case when Algorithm~\ref{sguT:algorithm-sg} is able to find stationary points.
We have the following convergence results.

\begin{theorem}\label{sguT:thm-sg-smT}
Suppose Assumptions~\ref{sgd:assgraph}--\ref{sgd:ass:stochastic-grad:variance} hold. Let $\{\bsx_k\}$ be the sequence generated by Algorithm~\ref{sguT:algorithm-sg} with
\begin{align}\label{sguT:step:eta2-sm}
\alpha_k=\kappa_1\beta_k,~\beta_k=\beta,~ \eta_k=\frac{\kappa_2}{\beta_k},~\forall k\in\mathbb{N}_0,
\end{align}
where $\kappa_1>c_1$, $\kappa_2\in(0,c_2(\kappa_1))$, and $\beta\ge c_0(\kappa_1,\kappa_2)$ with $c_0(\kappa_1,\kappa_2),~c_1,~c_2(\kappa_1)>0$  defined in Appendix~\ref{sguT:proof-thm-sg-smT}.
Then, for any $T\in\mathbb{N}_+$,
\begin{subequations}
\begin{align}
&\frac{1}{T}\sum_{k=0}^{T-1}\mathbf{E}\Big[\frac{1}{n}\sum_{i=1}^{n}\|x_{i,k}-\bar{x}_k\|^2\Big]
=\mathcal{O}(\frac{1}{T})
+\mathcal{O}(\frac{1}{\beta^2}),\label{sguT:thm-sg-sm-equ3.1}\\
&\frac{1}{T}\sum_{k=0}^{T-1}\mathbf{E}[\|\nabla f(\bar{x}_k)\|^2]
=\mathcal{O}(\frac{\kappa_2\beta}{T})
+\mathcal{O}(\frac{4\kappa_2}{n\beta})
+\mathcal{O}(\frac{1}{T})+\mathcal{O}(\frac{1}{\beta^2}),\label{sguT:thm-sg-sm-equ3}\\
&\mathbf{E}[f(\bar{x}_{T})]-f^*=\mathcal{O}(1)
+\mathcal{O}(\frac{T}{n\beta^2})+\mathcal{O}(\frac{T}{\beta^3}),\label{sguT:thm-sg-sm-equ4}
\end{align}
\end{subequations}
where $f^*$ is the minimum function value of the optimization problem \eqref{sgd:eqn:xopt} and $\bar{x}_k=\frac{1}{n}\sum_{i=1}^{n}x_{i,k}$.
\end{theorem}
\begin{proof}
The explicit expressions of the right-hand sides of \eqref{sguT:thm-sg-sm-equ3.1}--\eqref{sguT:thm-sg-sm-equ4} and the proof are given in Appendix~\ref{sguT:proof-thm-sg-smT}. It should be highlighted that the omitted constants in the first two terms in the right-hand side of  \eqref{sguT:thm-sg-sm-equ3} do not depend on any parameters related to the communication network.
\end{proof}

Noting the right-hand side of \eqref{sguT:thm-sg-sm-equ3} and $f(\bar{x}_0)-f^*=\mathcal{O}(1)$, the linear speedup (w.r.t. number of agents) can be established if we set $\beta=\kappa_2\sqrt{T}/\sqrt{n}$, as shown in the following.

\begin{corollary}[\bf Linear Speedup]\label{sguT:coro-sg-smT}
Under the same assumptions as in Theorem~\ref{sguT:thm-sg-smT}, let $\beta=\kappa_2\sqrt{T}/\sqrt{n}$. Then, for any $T>  \max\{n(c_0(\kappa_1,\kappa_2)/\kappa_2)^2,~n^3\}$,
\begin{subequations}
\begin{align}
&\frac{1}{T}\sum_{k=0}^{T-1}\mathbf{E}\Big[\frac{1}{n}\sum_{i=1}^{n}\|x_{i,k}-\bar{x}_k\|^2\Big]
=\mathcal{O}(\frac{n}{T}),\label{sguT:coro-sg-sm-equ3.1}\\
&\frac{1}{T}\sum_{k=0}^{T-1}\mathbf{E}[\|\nabla f(\bar{x}_k)\|^2]
=\mathcal{O}(\frac{1}{\sqrt{nT}})+\mathcal{O}(\frac{n}{T}),\label{sguT:coro-sg-sm-equ3}\\
&\mathbf{E}[f(\bar{x}_{T})]-f^*=\mathcal{O}(1).\label{sguT:coro-sg-sm-equ4}
\end{align}
\end{subequations}
\end{corollary}

\begin{remark}
It should be highlighted that the omitted constants in the first term in the right-hand side of  \eqref{sguT:coro-sg-sm-equ3} do not depend on any parameters related to the communication network.
The same linear speedup result as in \eqref{sguT:coro-sg-sm-equ3} was also established by the SGD algorithms proposed in \cite{jiang2018linear,yu2019parallel,Yu2019on,Yu2019Computation,lian2017can,Assran2019Stochastic,
tang2018communication,taheri2020quantized,wang2018cooperative,Tang2018Decentralized,Singh2020Communication}.
However, in \cite{jiang2018linear,Yu2019on,lian2017can,Assran2019Stochastic,tang2018communication},
the additional assumption that the deviation between the gradients of local cost functions is bounded was made;
in \cite{yu2019parallel,taheri2020quantized,Singh2020Communication}, it was required that each local stochastic gradient has bounded second moment;
in \cite{Yu2019Computation,wang2018cooperative}, it was assumed that the mean of each local stochastic gradient is the gradient of the global cost function;
and in \cite{Tang2018Decentralized}, it was  required that the eigenvalues of the mixing matrix are strictly greater than $-1/3$. Moreover, the algorithms proposed in \cite{jiang2018linear,Yu2019Computation} are restricted to a star graph;
the distributed momentum SGD algorithm proposed in \cite{Yu2019on} requires each agent $i$ to communicate one additional $p$-dimensional variable besides the communication of $x_{i,k}$ with its neighbors at each iteration;
and the algorithm proposed in \cite{Yu2019Computation} requires an exponentially increasing batch size, which is not favorable in practice.
Under the same conditions,  the well known $\mathcal{O}(1/\sqrt{T})$ convergence rate, which is not a speedup, was achieved by the distributed stochastic gradient tracking algorithm proposed in \cite{lu2019gnsd,zhang2019decentralized}.
Moreover, similar to the distributed momentum SGD algorithm proposed in \cite{Yu2019on}, one potential drawback of the distributed stochastic gradient tracking algorithms is that at each iteration each agent needs to communicate one additional variable.
The potential drawbacks of the results stated in Corollary~\ref{sguT:coro-sg-smT} are that (i) we do not consider communication efficiency, which was considered in \cite{jiang2018linear,yu2019parallel,Yu2019Computation,tang2018communication,taheri2020quantized,
wang2018cooperative,Singh2020Communication}; and (ii) we use time-invariant undirected graphs rather than directed graphs as considered in \cite{Assran2019Stochastic,taheri2020quantized}. We leave the extension to the time-varying directed graphs with communication efficiency as future research directions.
\end{remark}

\subsection{Find Global Optimum}

Let us next consider cases when Algorithm~\ref{sguT:algorithm-sg} finds global optimum.
The following assumption is crucial.

\begin{assumption}\label{sgd:assfil} The global cost function $f(x)$ satisfies the Polyak--{\L}ojasiewicz (P--{\L}) condition with constant $\nu>0$, i.e.,
\begin{align}
\frac{1}{2}\|\nabla f(x)\|^2\ge \nu( f(x)-f^*),~\forall x\in \mathbb{R}^p.\label{nonconvex:equ:plc}
\end{align}
\end{assumption}

It is straightforward to see that every (essentially or weakly) strongly convex function satisfies the P--{\L} condition.
The P--{\L} condition implies that every stationary point is a global minimizer. But unlike (essentially or weakly) strong convexity, the P--{\L} condition  alone does not imply convexity of $f$. Moreover, it does not imply that the global minimizer is a unique either \cite{karimi2016linear,zhang2015restricted}.

Many practical applications, such as least squares and logistic regression, do not always have strongly convex cost functions.
The cost function in least squares problems has the form
\begin{align*}
f(x)=\frac{1}{2}\|Ax-b\|^2,
\end{align*}
where $A\in\mathbb{R}^{m\times p}$ and $b\in\mathbb{R}^m$. Note that if $A$ has full column rank, then $f(x)$ is strongly convex. However, if $A$ is rank deficient, then $f(x)$ is not strongly convex, but it is convex and satisfies the P--{\L} condition.
Examples of nonconvex functions which satisfy the P--{\L} condition can be found in \cite{karimi2016linear,zhang2015restricted}.

Although it is difficult to precisely characterize the general class of functions for which the P--{\L} condition is satisfied, in \cite{karimi2016linear}, one important special class was given as follows:
\begin{lemma}
Let $f(x)=g(Ax)$, where $g:\mathbb{R}^p\rightarrow\mathbb{R}$ is a strongly convex function and $A\in\mathbb{R}^{p\times p}$ is a matrix, then $f$ satisfies the P--{\L} condition.
\end{lemma}

We have the following global convergence results.

\begin{theorem}\label{sguT:thm-sg-diminishingT}
Suppose Assumptions~\ref{sgd:assgraph}--\ref{sgd:assfil} hold. For any given $T\ge(c_0(\kappa_1,\kappa_2)/\kappa_2)^{1/\theta}$, let $\{\bsx_0,\dots,\bsx_T\}$ be the sequence generated by Algorithm~\ref{sguT:algorithm-sg} with
\begin{align}\label{sguT:step:eta1T}
\alpha_k=\kappa_1\beta_k,~\beta_k=\kappa_2(T+1)^\theta,
~ \eta_k=\frac{\kappa_2}{\beta_k},~\forall k\le T,
\end{align}
where $\theta\in(0,1)$, $\kappa_1>c_1$, $\kappa_2\in(0,c_2(\kappa_1))$.
Then,
\begin{subequations}
\begin{align}
&\mathbf{E}\Big[\frac{1}{n}\sum_{i=1}^{n}\|x_{i,T}-\bar{x}_T\|^2\Big]
=\mathcal{O}(\frac{1}{T^{2\theta}}), \label{sguT:thm-sg-equ2.1bounded}\\
&\mathbf{E}[f(\bar{x}_{T})-f^*]
=\mathcal{O}(\frac{1}{nT^\theta})+\mathcal{O}(\frac{1}{T^{2\theta}}).
\label{sguT:thm-sg-equ2bounded}
\end{align}
\end{subequations}
\end{theorem}
\begin{proof}
The explicit expressions of the right-hand sides of \eqref{sguT:thm-sg-equ2.1bounded} and \eqref{sguT:thm-sg-equ2bounded}, and the proof are given in Appendix~\ref{sguT:proof-thm-sg-diminishingT}. It should be highlighted that the omitted constants in the first term in the right-hand side of \eqref{sguT:thm-sg-equ2bounded} do not depend on any parameters related to the communication network.
\end{proof}

From Theorem~\ref{sguT:thm-sg-diminishingT}, we see that the convergence rate is strictly greater than $\mathcal{O}(1/(nT))$. In the following we show that the linear speedup convergence rate $\mathcal{O}(1/(nT))$ can be achieved if the P--{\L} constant $\nu$ is known in advance and each $f_i^*>-\infty$, where $f_i^*=\min_{x\in\mathbb{R}^p}f_i(x)$. The total number of iterations $T$ is not needed.
\begin{theorem}[\bf Linear Speedup]\label{sgu:thm-sg-diminishingt}
Suppose Assumptions~\ref{sgd:assgraph}--\ref{sgd:assfil} hold, and the P--{\L} constant $\nu$ is known in advance, and each $f_i^*>-\infty$.
Let $\{\bsx_{k}\}$ be the sequence generated by Algorithm~\ref{sguT:algorithm-sg} with \begin{align}\label{sgu:step:eta1t1}
\alpha_k=\kappa_1\beta_k,~\beta_k=\kappa_0(k+t_1),~ \eta_k=\frac{\kappa_2}{\beta_k},~\forall k\in\mathbb{N}_0,
\end{align}
where $\kappa_0\in[\hat{c}_0\nu\kappa_2/4,\nu\kappa_2/4)$, $\kappa_1>c_1$, $\kappa_2\in(0,\hat{c}_2(\kappa_1))$, and $t_1>\hat{c}_3(\kappa_0,\kappa_1,\kappa_2)$ with $\hat{c}_0\in(0,1)$ being a constant, $\hat{c}_2(\kappa_1)$ and $\hat{c}_3(\kappa_0,\kappa_1,\kappa_2)$ defined in Appendix~\ref{sgu:proof-thm-sg-diminishingt}.
Then, for any $T\in\mathbb{N}_+$,
\begin{subequations}
\begin{align}
&\mathbf{E}\Big[\frac{1}{n}\sum_{i=1}^{n}\|x_{i,T}-\bar{x}_T\|^2\Big]
=\mathcal{O}(\frac{1}{T^{2}}), \label{sgu:thm-sg-diminishing-equ2.1bounded}\\
&\mathbf{E}[f(\bar{x}_{T})-f^*]
=\mathcal{O}(\frac{1}{nT})+\mathcal{O}(\frac{1}{T^{2}}).
\label{sgu:thm-sg-diminishing-equ2bounded}
\end{align}
\end{subequations}
\end{theorem}
\begin{proof}
The explicit expressions of the right-hand sides of \eqref{sgu:thm-sg-diminishing-equ2.1bounded} and \eqref{sgu:thm-sg-diminishing-equ2bounded}, and the proof are given  in Appendix~\ref{sgu:proof-thm-sg-diminishingt}. It should be highlighted that the omitted constants in the first term in the right-hand side of  \eqref{sgu:thm-sg-diminishing-equ2bounded} do not depend on any parameters related to the communication network.
\end{proof}
\begin{remark}
It has been shown in \cite{rakhlin2011making} that $\mathcal{O}(1/T)$ convergence rate is optimal for centralized strongly convex optimization. This rate has been established by various distributed SGD algorithms when each local cost function is strongly convex, e.g., \cite{reisizadeh2019fedpaq,jiang2017collaborative,rabbat2015multi,lan2018communication,
yuan2018optimal,jakovetic2018convergence,fallah2019robust}. In contrast, the linear speedup convergence rate $\mathcal{O}(1/(nT))$ established in Theorem~\ref{sgu:thm-sg-diminishingt} only requires that the global cost function satisfies the P--{\L} condition, but no convexity assumption is required neither on the global cost function nor on the local cost functions. The SGD algorithms in \cite{haddadpour2019local,Yu2019Computation,stich2018local,Koloskova2019decentralized,
Singh2020Communication,olshevsky2019non} also achieve the linear speedup convergence rate.
However, the algorithms in \cite{haddadpour2019local,Yu2019Computation,stich2018local} are restricted to a star graph, while our algorithm is applicable to an arbitrarily connected graph.
Moreover, \cite{haddadpour2019local,Yu2019Computation} assumed that the mean of each local stochastic gradient is the gradient of the global cost function, and $T$ has to be known to choose the algorithm parameters. The algorithm in \cite{Yu2019Computation} furthermore requires an exponentially increasing batch size, which is not favorable in practice. In \cite{stich2018local}, it was assumed that the global cost function is strongly convex. In \cite{stich2018local,Singh2020Communication}, it was assumed that  each local stochastic gradient has bounded second moment. In \cite{Koloskova2019decentralized,Singh2020Communication,
olshevsky2019non}, it was assumed that each local cost function is strongly convex. It is one of our future research directions to achieve linear speedup with reduced communication rounds and communication efficiency for an  arbitrarily connected graph.
\end{remark}

Theorem~\ref{sgu:thm-sg-diminishingt} shows that the convergence rate to a global optimum is sublinear when we allow the algorithm parameters to be time-varying.
The following theorem establishes that the output of Algorithm~\ref{sguT:algorithm-sg} with constant algorithm parameters linearly converges to a neighborhood of a global optimum.

\begin{theorem}\label{sg:thm-sg-fixed-a5}
Suppose Assumptions~\ref{sgd:assgraph}--\ref{sgd:assfil} hold. Let $\{\bsx_{k}\}$ be the sequence generated by Algorithm~\ref{sguT:algorithm-sg} with \begin{align}\label{sgu:step:fixed-a5}
\alpha_k=\alpha=\kappa_1\beta,~\beta_k=\beta,~ \eta_k=\eta=\frac{\kappa_2}{\beta},~\forall k\in\mathbb{N}_0,
\end{align}
where $\kappa_1>c_1$, $\kappa_2\in(0,c_2(\kappa_1))$, and $\beta\ge c_0(\kappa_1,\kappa_2)$ with $c_0(\kappa_1,\kappa_2),~c_1,~c_2(\kappa_1)>0$  defined in Appendix~\ref{sguT:proof-thm-sg-smT}. Then,
\begin{align}\label{sg:thm-sg-fixed-equ1-a5}
&\mathbf{E}\Big[\frac{1}{n}\sum_{i=1}^{n}\|x_{i,k}-\bar{x}_k\|^2+f(\bar{x}_k)-f^*\Big]
\le(1-\eta\varepsilon)^{k}c_4+c_5\eta\sigma^2,~\forall k\in\mathbb{N}_+,
\end{align}
where $\varepsilon\in(0,1/\eta),~c_4,~c_5>0$ are constants defined in Appendix~\ref{sg:proof-thm-sg-fixed-a5}.

\end{theorem}
\begin{proof}
The proof is given in Appendix~\ref{sg:proof-thm-sg-fixed-a5}.
\end{proof}
\begin{remark}
It should be highlighted that we  do not need to know the P--{\L} constant $\nu$ in advance. Similar convergence result as stated in \eqref{sg:thm-sg-fixed-equ1-a5} was achieved by the distributed SGD algorithms proposed in \cite{jiang2017collaborative,fallah2019robust,pu2018swarming,pu2018distributed,xin2019distributed} when each local cost function is strongly convex, which obviously is stronger than the P--{\L} condition assumed in Theorem~\ref{sg:thm-sg-fixed-a5}. In addition to the strong convexity condition, in \cite{jiang2017collaborative}, it was also assumed that each local cost function is Lipschitz-continuous. Some information related to the Lyapunov function and global parameters, which may be difficult to get, were furthermore needed to design the stepsize. Moreover, in \cite{fallah2019robust,pu2018swarming,pu2018distributed,xin2019distributed}, the strong convexity constant was needed to design the stepsize and in \cite{pu2018distributed,xin2019distributed}, a $p$-dimensional auxiliary variable, which is used to track the global gradient, was communicated between agents. The potential drawbacks of the results stated in Theorem~\ref{sg:thm-sg-fixed-a5} are that (i) we use undirected graphs rather than directed graphs as considered in \cite{xin2019distributed}; and (ii) we do not analyze the robustness level to gradient noise as \cite{fallah2019robust} did. We leave the extension to the (time-varying) directed graphs and the robustness level analysis as future research directions.
\end{remark}

Actually, the unbiased assumption, i.e., Assumption~\ref{sgd:ass:stochastic-grad:mean}, can be removed, as shown in the following.

\begin{theorem}[\bf Biased SGD]\label{sg:thm-sg-fixed}
Suppose Assumptions~\ref{sgd:assgraph}--\ref{sgd:ass:stochastic-grad:xi}, \ref{sgd:ass:stochastic-grad:variance}, and \ref{sgd:assfil} hold. Let $\{\bsx_{k}\}$ be the sequence generated by Algorithm~\ref{sguT:algorithm-sg} with
\begin{align}\label{sgu:step:fixed}
\alpha_k=\alpha=\kappa_1\beta,~\beta_k=\beta,~ \eta_k=\eta=\frac{\kappa_2}{\beta},~\forall k\in\mathbb{N}_0,
\end{align}
where $\kappa_1>c_1$, $\kappa_2\in(0,c_2(\kappa_1))$, and $\beta\ge \breve{c}_0(\kappa_1,\kappa_2)$ with $\breve{c}_0(\kappa_1,\kappa_2)>0$ and $c_1,~c_2(\kappa_1)>0$  defined in Appendices~\ref{sg:proof-thm-sg-fixed} and \ref{sguT:proof-thm-sg-smT}, respectively. Then,
\begin{align}\label{sg:thm-sg-fixed-equ1}
&\mathbf{E}[\frac{1}{n}\sum_{i=1}^{n}\|x_{i,k}-\bar{x}_k\|^2+f(\bar{x}_k)-f^*]
\le(1-\eta\varepsilon)^{k}c_4+\breve{c}_5\sigma^2,~\forall k\in\mathbb{N}_+,
\end{align}
where $\varepsilon\in(0,1/\eta),~c_4>0$ and $\breve{c}_5>0$ are constants defined in Appendices~\ref{sg:proof-thm-sg-fixed-a5} and \ref{sg:proof-thm-sg-fixed}, respectively.

\end{theorem}
\begin{proof}
The proof is given in Appendix~\ref{sg:proof-thm-sg-fixed}.
\end{proof}

\begin{remark}
By comparing \eqref{sg:thm-sg-fixed-equ1-a5} with \eqref{sg:thm-sg-fixed-equ1}, we can see that no matter the unbiased assumption holds or not,  the output of Algorithm~\ref{sguT:algorithm-sg} with constant algorithm parameters linearly converges to a neighborhood of a global optimum, but the size of neighborhood is different. Specifically, in \eqref{sg:thm-sg-fixed-equ1-a5} the size of neighborhood is in an order of $\mathcal{O}(\eta)$, while it is $\mathcal{O}(1)$ in \eqref{sg:thm-sg-fixed-equ1}. When true gradients are available, i.e., $\sigma=0$, then from \eqref{sg:thm-sg-fixed-equ1-a5} or \eqref{sg:thm-sg-fixed-equ1} we know that a global
optimum can be linearly found. It should be highlighted that this linear convergence result is established under the P--{\L} condition and the P--{\L} constant is not used. These are two advantages since in existing studies obtaining linear convergence for distributed smooth optimization, e.g., \cite{shi2015extra,nedic2017achieving,li2020revisiting}, it is standard to assume (restricted) strong convexity, which is stronger than the P--{\L} condition, and to use the convexity parameter.
\end{remark}

\section{Simulations}\label{sgd:sec-simulation}

In this section, we evaluate the performance of the proposed distributed primal--dual SGD algorithm through numerical experiments.
All algorithms and agents are implemented and simulated in MATLAB R2018b, run on a desktop with \texttt{Intel Core i5-9600K processor, Nvidia RTX 2070 super, 32~GB RAM, Ubuntu 16.04}.

\subsection{Neural Networks}
We consider the training of neural networks (NN) for image classification tasks of the database MNIST \cite{lecun2010mnist}.
The same NN is adopted as in \cite{george2019distributed} for each agent and the communication graph is generated randomly. The graph is shown in Fig.~\ref{sgd:fig:graph} and the corresponding Laplacian matrix $\bm L$ is given in \eqref{sgd:lap}. The corresponding mixing matrix $\bm W$ is constructed by metropolis weight, which is given in \eqref{sgd:lapW}.

\begin{align} \label{sgd:lap}
\bm L =
\begin{bmatrix}
1 & -1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\
-1 & 3 & -1 & -1 & 0 & 0 & 0 & 0 & 0 & 0\\
0 & -1 & 3 & -1 & 0 & 0 & -1 & 0 & 0 & 0\\
0 & -1 & -1 & 4 & -1 & -1 & 0 & 0 & 0 & 0\\
0 & 0 & 0 & -1 & 2 & -1 &  0 & 0 & 0 & 0\\
0 & 0 & 0 & -1 & -1 & 2 &  0 & 0 & 0 & 0\\
0 & 0 & -1 & 0 & 0 & 0 &  2 & -1 & 0 & 0\\
0 & 0 & 0 & 0 & 0 & 0 &  -1 & 2 & -1 & 0\\
0 & 0 & 0 & 0 & 0 & 0 & 0 &  -1 & 2 & -1 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 &  0 & -1 & 1\\
\end{bmatrix}.
\end{align}
\begin{align}\label{sgd:lapW}
 \bm W =
\begin{bmatrix}
3/4 & 1/4 & 0 & 0 & 0 & 0 & 0 & 0 & 0 & 0\\
1/4 & 3/10 & 1/4 & 1/5 & 0 & 0 & 0 & 0 & 0 & 0\\
0 & 1/4 & 3/10 & 1/5 & 0 & 0 & 1/4 & 0 & 0 & 0\\
0 & 1/5 & 1/5 & 1/5 & 1/5 & 1/5 & 0 & 0 & 0 & 0\\
0 & 0 & 0 & 1/5 & 7/15 & 1/3 &  0 & 0 & 0 & 0\\
0 & 0 & 0 & 1/5 & 1/3 & 7/15 &  0 & 0 & 0 & 0\\
0 & 0 & 1/4 & 0 & 0 & 0 &  5/12 & 1/3 & 0 & 0\\
0 & 0 & 0 & 0 & 0 & 0 &  1/3 & 1/3 & 1/3 & 0\\
0 & 0 & 0 & 0 & 0 & 0 & 0 &  1/3 & 1/3 & 1/3 \\
0 & 0 & 0 & 0 & 0 & 0 & 0 &  0 & 1/3 & 2/3\\
\end{bmatrix}.
\end{align}

Each local neural network consists of a single hidden layer of 50 neurons, followed by a sigmoid activation layer, followed by the output layer of 10 neurons and another sigmoid activation layer. In this experiment, we use a subset of MNIST data set. Each agent is assigned 2500 data points randomly, and at each iteration, only one data point is picked up by the agent following a uniform distribution.

We compare our proposed distributed primal--dual SGD algorithm with time-varying and fixed  parameters (DPD-SGD-T and DPD-SGD-F)  with state-of-the-art algorithms: distributed momentum SGD algorithm (DM-SGD) \cite{Yu2019on},  distributed SGD algorithm (D-SGD-1) \cite{jiang2017collaborative,lian2017can},  distributed SGD algorithm (D-SGD-2) \cite{george2019distributed},  $\mathrm{D}^2$ \cite{Tang2018Decentralized}, distributed stochastic gradient tracking algorithm (D-SGT-1)  \cite{lu2019gnsd,xin2019distributed},  distributed stochastic gradient tracking algorithm  (D-SGT-2) \cite{zhang2019decentralized,pu2018distributed}, and  the baseline centralized SGD algorithm (C-SGD). We list all the parameters\footnote{Note: the parameter names are different in each paper.} we choose in the NN experiment for each algorithm in Table~\ref{tab:nn-par}.

We demonstrate the result in terms of the empirical risk function \cite{bottou2012stochastic}, which is given as
\begin{align*}
R(\bm z) &= - \frac{1}{n}\sum_{i = 1}^{n}\frac{1}{m_{n}} \sum _{j=1}^{m_{n}} \sum _{k=0}^{9}(t_{k} \ln y_k(\bm{x}, \bm z) + (1-t_{k})\ln (1 - y_k(\bm{x}, \bm z)))
\end{align*}
where $m_{n}$ indicates the size of data set for each agent, $t_k$ denotes the target (ground truth) of digit $k$ corresponding to a single image, $\bm x$ is a single image input, $\bm z=(z^{(1)}, z^{(2)})$ with $z^{(1)}$ and $z^{(2)}$ being the weights in the 2 layers separately, and $y_k \in [0, 1]$ is the output which expresses the probability of digit $k = 0, \dots, 9$.
The mapping from input to output is given as:
\begin{align*}
y_k(\bm{x}, \bm z) = \sigma\left( \sum_{j=0}^{50} z_{k, j}^{(2)} \sigma \left( \sum_{i = 0}^{28 \times 28} z_{j, i}^{(1)}x_{i}\right)\right),
\end{align*}
where $\sigma (s) = \frac{1}{1+\exp(-s)}$ is the sigmoid function.

\begin{figure}
\centering
\hspace{10mm}\centering{
{
\begin{tikzpicture}[-,node distance=1.4cm,
  thick,main node/.style={circle,fill=yellow!20,draw,font=\sffamily\normalsize\bfseries}]
  \node[main node] (1) {1};
  \node[main node] (2) [above right of=1] {2};
  \node[main node] (3) [right of=2] {3};
  \node[main node] (4) [below right of=2] {4};
  \node[main node] (5) [below right of=1] {5};
  \node[main node] (6) [ below of=5] {6};
  \node[main node] (7) [ right of=3] {7};
  \node[main node] (8) [ below of=7] {8};
  \node[main node] (9) [ below of=8] {9};
  \node[main node] (10) [ left of=9] {10};
\draw (1) -- (2)
(2) -- (3)
(3) -- (7)
(2) -- (4)
(3) -- (4)
(4) -- (5)
(4) -- (6)
(5) -- (6)
(7) -- (8)
(8) -- (9)
(9) -- (10);
\end{tikzpicture}}}$\qquad\qquad$
\vspace{1mm}
\caption{\label{sgd:fig:graph} Connection Topology.}
\end{figure}

\begin{table*}[htbp]
\caption{{Parameters in each algorithm in NN experiment.} }
\label{tab:nn-par}
\vskip 0.15in
\begin{center}
\begin{normalsize}
\begin{tabular}{M{3.5cm}|M{3.5cm}|M{3.5cm}|M{3.5cm}}
\hline
Algorithm & $\eta_k$ & $\alpha_k$ & $\beta_k$ \\

\hline
DPD-SGD-T & $0.08/{k^{10^{-5}}}$ & $4k^{10^{-5}}$ & $3k^{10^{-5}}$\\

\hline
DPD-SGD-F &0.03 & 5 & 20\\

\hline
DM-SGD \cite{Yu2019on} & 0.1 & \ding{55} &0.8\\

\hline
D-SGD-1 \cite{jiang2017collaborative,lian2017can}& 0.1 &\ding{55} &\ding{55}\\

\hline
D-SGD-2 \cite{george2019distributed} &\ding{55} & $0.1/(10^{-5}k + 1)$&$0.2/(10^{-5}k + 1)^{0.3}$\\

\hline
$D^{2}$ \cite{Tang2018Decentralized}& 0.01 &\ding{55} &\ding{55}\\

\hline
D-SGT-1 \cite{lu2019gnsd,xin2019distributed} & 0.01& \ding{55}&\ding{55}\\

\hline
D-SGT-2 \cite{zhang2019decentralized,pu2018distributed} & 0.01 & \ding{55}&\ding{55}\\

\hline
C-SGD & 0.1 &\ding{55} &\ding{55}\\
\hline
\end{tabular}
\end{normalsize}
\end{center}
\vskip -0.1in
\end{table*}

\begin{figure}[htbp]
\centering
  \includegraphics[width=0.85\textwidth]{nn_risk.eps}
  \caption{Empirical Risk.}
  \label{sgd:fig:nn_risk}
\end{figure}

Fig.~\ref{sgd:fig:nn_risk} shows that the proposed distributed primal--dual SGD algorithm with time-varying parameters converges almost as fast as the distributed SGD algorithm in \cite{jiang2017collaborative,lian2017can} and faster than the distributed SGD algorithms in \cite{george2019distributed,Tang2018Decentralized,lu2019gnsd,zhang2019decentralized,
pu2018distributed,xin2019distributed} and the centralized SGD algorithm. Note that our algorithm converges slower than the distributed momentum SGD algorithm \cite{Yu2019on}. This is reasonable since that algorithm is an accelerated algorithm with extra requirement on the cost functions, i.e., the deviations between the gradients of local cost functions is bounded, and it requires each agent to communicate three $p$-dimensional variables with its neighbors at each iteration. The slope of the curves are however almost the same. The accuracy of each algorithm is given in Table~\ref{tab:nn-acc}.

\begin{table*}[htbp]
\caption{{Accuracy on each algorithm in NN experiment.} }
\label{tab:nn-acc}
\vskip 0.15in
\begin{center}
\begin{normalsize}
\begin{tabular}{M{4.1cm}|M{3.2cm}}
\hline
Algorithm & Accuracy  \\

\hline
DPD-SGD-T & $93.04\%$\\

\hline
DPD-SGD-F & $92.76\%$\\

\hline
DM-SGD \cite{Yu2019on} & $93.44\%$\\

\hline
D-SGD-1 \cite{jiang2017collaborative,lian2017can}& $92.96\%$\\

\hline
D-SGD-2 \cite{george2019distributed} & $92.88\%$\\

\hline
$D^{2}$ \cite{Tang2018Decentralized}& $90.44\%$\\

\hline
D-SGT-1 \cite{lu2019gnsd,xin2019distributed} & $92.88\%$\\

\hline
D-SGT-2 \cite{zhang2019decentralized,pu2018distributed} & $92.96\%$\\

\hline
C-SGD & $93\%$\\
\hline
\end{tabular}
\end{normalsize}
\end{center}
\vskip -0.1in
\end{table*}

\subsection{Convolutional Neural Networks}
Let us consider the training of a convolutional neural networks (CNN) model.
We build a CNN model for each agent  with five 3$\times$3 convolutional layers using ReLU as activation function, one average pooling layer with filters of size 2$\times$2, one sigmoid layer with dimension 360, another sigmoid layer with dimension 60, one softmax layer with dimension 10. In this experiment, we use the whole MNIST data set. We use the same  communication graph as in above NN  experiment. Each agent is assigned 6000 data points randomly. We set the batch size as 20, which means at each iteration, 20 data points are chosen by the agent to update the gradient, which is also following a uniform distribution. For each algorithm, we do 10 epochs to train the CNN model.

We compare our algorithms DPD-SGD-T and DPD-SGD-F  with the fastest one above: DM-SGD, D-SGD-1, and C-SGD. We list all the parameters we choose in the CNN experiment for each algorithm in Table~\ref{tab:cnn-par}.

\begin{table*}[ht!]
\caption{{Parameters in each algorithm in CNN experiment.} }
\label{tab:cnn-par}
\vskip 0.15in
\begin{center}
\begin{normalsize}
\begin{tabular}{M{3.5cm}|M{2.5cm}|M{2.5cm}|M{2.5cm}}
\hline
Algorithm & $\eta_k$ & $\alpha_k$ & $\beta_k$ \\

\hline
DPD-SGD-T & $0.5/{k^{10^{-5}}}$ & $0.5k^{10^{-5}}$ & $0.1k^{10^{-5}}$\\

\hline
DPD-SGD-F &0.5 & 0.5 & 0.1\\

\hline
DM-SGD \cite{Yu2019on} & 0.1 & \ding{55} &0.8\\

\hline
D-SGD \cite{jiang2017collaborative,lian2017can}& 0.1 &\ding{55} &\ding{55}\\

\hline
C-SGD & 0.1 &\ding{55} &\ding{55}\\
\hline
\end{tabular}
\end{normalsize}
\end{center}
\vskip -0.1in
\end{table*}

We demonstrate the training loss and the test accuracy of each algorithm in Fig.~\ref{sgd:fig:cnn_loss} and Fig.~\ref{sgd:fig:cnn_acc}. Here we use Categorical Cross-Entropy loss, which is a softmax activation plus a Cross-Entropy loss. We can see that our algorithms perform almost the same as the DM-SGD and better than the D-SGD-1 and  the centralized C-SGD. The accuracy of each algorithm is given in Table~\ref{tab:cnn-acc}.

\begin{figure}[!ht]
\centering
        \includegraphics[width=0.85\textwidth]{cnn_loss.eps}
        \caption{CNN training loss.}
        \label{sgd:fig:cnn_loss}
\end{figure}

\begin{figure}[!ht]
\centering
        \includegraphics[width=0.85\textwidth]{cnn_acc.eps}
        \caption{CNN accuracy.}
        \label{sgd:fig:cnn_acc}
\end{figure}

\begin{table*}[ht!]
\caption{{Accuracy on each algorithm in CNN experiment.} }
\label{tab:cnn-acc}
\vskip 0.15in
\begin{center}
\begin{normalsize}
\begin{tabular}{M{3.5cm}|M{2.5cm}}
\hline
Algorithm & Accuracy  \\

\hline
DPD-SGD-T & $94.75\%$\\

\hline
DPD-SGD-F & $93.17\%$\\

\hline
DM-SGD \cite{Yu2019on} & $94.29\%$\\

\hline
D-SGD \cite{jiang2017collaborative,lian2017can}& $92.96\%$\\

\hline
C-SGD & $89.91\%$\\
\hline
\end{tabular}
\end{normalsize}
\end{center}
\vskip -0.1in
\end{table*}

\section{Conclusions}\label{sgd:sec-conclusion}
In this paper, we studied distributed nonconvex optimization. We proposed a distributed primal--dual SGD algorithm and derived its convergence rate. More specifically, the linear speedup convergence rate $\mathcal{O}(1/\sqrt{nT})$ was established for smooth nonconvex cost functions under  arbitrarily connected communication networks.  The convergence rate was improved to the linear speedup convergence rate $\mathcal{O}(1/(nT))$ when the global cost function additionally satisfies the P--{\L} condition. It was also shown that the output of the proposed algorithm with constant parameters linearly converges to a neighborhood of a global optimum. Interesting directions for future work include achieving linear speedup under the P--{\L} condition while considering communication reduction.

\bibliographystyle{IEEEtran}
\bibliography{refextra}


TITLE: How do eigenvalues of a matrix X change if we linear transform the matrix X?
QUESTION [2 upvotes]: I have a matrix $X$ which has eigenvalues $U$.
Now create a new matrix $Y = AX$ where $A$ is a nonsingular matrix.
How do the eigenvectors and eigenvalues of $Y$ change in relation to the eigenvectors and eigenvalues of $X$ and the matrix $A$?

REPLY [0 votes]: If the eigenvector equation for $X$ is
$$
Xv=\lambda v
$$
then
$$
Yv=AXv=A\lambda v=\lambda Av
$$
so it comes down to what effect $A$ has on $v$ which will vary according to $A$.
In special cases, such as $A$ being a constant diagonal matrix with value $k$ on the diagonal then the above equation continues
$$
\lambda Av=\lambda k Iv = \lambda kv
$$
In such a special case the eigenvectors of $X$ and $A$ are the same and the eigenvalues are multiplied by $k$.


TITLE: Prove the existence or the non-existence of a couple of numbers ($n$,$m$) such that $n^2=m!$
QUESTION [2 upvotes]: In recent days, while I was doing exercises on combinatorics, I thought if a number $m!$ could be a perfect square. I proved to demonstrate it through the prime factorization. My attempt: $$k=m!=p_1^{a_1}\cdot p_2^{a_2}\cdot p_3^{a_3}....\cdot p_n^{a_n}$$
with $a_1,a_2,..,a_n$ that have to be even. Now if there exists a prime number $p_j$ between $\frac {k}{2}<p_j<k$ this prime number in the prime factorization the exponent ($a_j$) relative to the prime number $p_j$ is equal to $1$ ($a_j=1$) because the number $k$ isn't a multiple of $p_j$; for example $12!$ isn't a perfect square because the prime numbers $7$ and $11$ $6<7,11<12$ aren't divisors of $12$.
My question
There exist a prime number $p_j$ between $k$ and $\frac {k}{2}$ and how can I prove if there exist or not exist a number $m!=n^2$?

REPLY [0 votes]: Subtracting $m^2$ from both sides of your equation we get $$n^2-m^2=m!-m^2 \tag{0} $$ $$ (n-m)(n+m)=m\left((m-1)!-m\right), $$ that is, assuming $m\ne0$, $$ \frac{n}{m}+1=\frac{(m-1)!-m}{n-m}.\tag{1}$$ We do know that $n$ is a multiple of $m$, therefore both sides of $(1)$ are integer. Now, $(1)$ is equivalent to $$\frac{n}{m}=\frac{(m-1)!-n}{n-m} $$ $$\frac{n}{m^2}=\frac{\frac{(m-1)!}{m}-\frac{n}{m}}{n-m} \tag{2}$$ as well as $$\frac{n}{m^2}+\frac{1}{m}=\frac{\frac{(m-1)!}{m}-1}{n-m} \tag{3}.$$ But while both the right-hand sides of  $(2)$ and $(3)$ are integer, since $(0)$ being rewritable as $$m^2\left(\left(\frac{n}{m}\right)^2 -1\right) =m( (m-1)! - m)$$ yields that $m$ divides $(m-1)!-m$, and so $(m-1)!$, their left-hand sides can be simultaneously integer if and only if $m=1$. Therefore, the only couples of solutions are $(m,n)=(0,\pm1),(1,\pm1)$.