1 Introduction↩︎

In recent years, the study of distributed (or decentralized) algorithms for optimization, learning, and classification over a network of computing nodes/agents has gained significant attention due to adavances in cloud-computing and parallel data processing [1]. These networks consist of multiple agents, each with limited computational and communication capabilities, working collaboratively to solve optimization problems or learn from data in a distributed manner. However, a critical challenge in these networks is the presence of time-delays [2], [3], which can arise from communication latencies, processing times, or network congestion. Time-delays can severely impact the performance and convergence of distributed algorithms, making it essential to develop robust methods that can handle such delays. This paper explores the theoretical foundations and practical implementations of distributed optimization and learning algorithms that are resilient to time-delays, providing mathematical proofs, analysis, and potential applications.

1.1 Problem and Contributions↩︎

In distributed optimization, the idea is to optimize a cost function (or loss function) over a network of computing nodes. The objective function is the sum of some local cost functions at the nodes, and the goal is to optimize this objective using locally defined gradient-based algorithms. The common form of the optimization problem is, \[\begin{align} \min_{\mathbf{z}} F(\mathbf{z}) = \sum_{i=1}^{N} f_i(\mathbf{z}) \label{eq95prob} \end{align}\tag{1}\] with state parameter \(\mathbf{z} \in \mathbb{R}^{m}\). Functions \(f_i:\mathbb{R}^m \mapsto \mathbb{R}\) are strongly convex, differentiable with Lipschitz gradients, and denote the objective function (cost, loss, etc.) at computing node \(i\). It is assumed that the optimal point \(\underline{\mathbf{z}}^* = \min_{\mathbf{z}} F(\mathbf{z})\) for this problem exists. The primary work [4] introduces subgradient algorithms to solve this problem. ADD-OPT algorithm [5] and its recent stochastic version S-ADD-OPT [6] are popular algorithms to solve problem 1 . These algorithms work over strongly-connected directed networks with irreducible column stochastic adjacency matrices, and are granted with (i) constant step-size in contrast to existing diminishing step-size algorithms, (ii) providing accelerated convergence by tuning the step-size over a wide range, and (iii) linear convergence rate for strongly convex cost functions. Other existing distributed algorithms include: event-triggered-based second-order multi-agent systems [7], [8], double step-size solutions for nonsmooth optimization [9], reduced-complexity and flexible algorithms [10], primal-dual subgradient-based solutions [11], EXTRA algorithm for first-order consensus-based optimization [12], push-pull gradient-based methods [13], and the solutions based on alternating direction method of multipliers (ADMM) [14]–[17]. The literature also includes distributed constrained optimization with application to resource allocation under time-delay. For, example, DTAC-ADMM discusses ADMM-based distributed resource allocation under time-delay [18]. Similarly, asynchronous ADMM-based resource allocation algorithms are proposed in [19]. These works consider distributed optimization subject to a coupling resource-demand balance constraint, where the objective functions are decoupled and local. Asynchronous distributed optimization is also discussed in [20], [21], where agents perform local computations and communications without requiring global synchronization. In such methods, each node updates its local model using the most recently available information from neighbors (which may be received at irregular times). Recall that scalability is a key advantage of existing distributed optimization techniques, which follows the polynomial-order complexity of the algorithms. Polynomial-order complexity ensures computationally-efficient solutions as the number of agents or decision variables increases, making it feasible to deploy these algorithms on large-scale networks.

In this work, as the main contribution, we extend such distributed optimization algorithms to further address arbitrary and bounded time-delays over multi-agent networks. Latency is primarily addressed in consensus literature including: resilient consensus with \(l\)-hop communication [22], multi-agent consensus subject to uncertainties and time-varying delays [23], group consensus over digraphs subject to noise and latency [24], continuous-time linear average consensus with constant delays at all nodes [25], discrete-time consensus algorithms with constant communication delays [26], discrete-time consensus over digraphs under heterogeneous time-delays [27]. These works are advantageous as they provide rigorous stability/convergence analysis applicable to other distributed setups; however, they mostly assume constant homogeneous delays. For a review of consensus algorithms under time-delays and their advantages/disadvantages, refer to [28]. The concept of time-delay is not sufficiently addressed in distributed optimization literature. The inherent time-delay of information exchange among communicating nodes may lead the distributed optimization algorithm to lose convergence. The delays are typically assumed to be bounded, implying that the information sent over every link eventually reaches the destination node, i.e., no packet loss over the network. In this paper, we propose augmented consensus-based algorithms to analyze the effect of time-delays while keeping the consensus matrix on the link weights column stochastic. Our solution can tolerate heterogeneous communication delays at different links. In this regard, similar to [29], this work improves the existing algorithms over non-delayed networks [5], [6], [12], [13], [30]–[32] to more advanced delay-tolerant solutions which are not well-addressed in the literature (to our best knowledge). This work also advances the existing ADMM-based solutions [14]–[16] to withstand latency and network time-delays. Our proposed delay-tolerant augmented consensus-based DTAC-ADDOPT algorithm is in single time-scale, i.e., it performs only one step of (augmented) consensus on received information per iteration/epoch. This is computationally more efficient in contrast to the double time-scale methods [33], [34] with many steps of inner-loop consensus per iteration/epoch. Heterogeneous time-delays are considered primarily for ADMM-based [18] and gradient-descent-based [3] equality-constraint distributed optimization and resource allocation, but this current paper is our first paper addressing it over unconstrained ADD-OPT.

1.2 Applications↩︎

Distributed Training for Binary Classification: Consider a group of agents to classify \(N\) data points \({\boldsymbol{\chi}_i \in \mathbb{R}^{m-1}}\), \({i=1,\ldots,N}\), labeled by \({l_i \in \{-1,1\}}\). The problem is to find the partitioning hyperplane \({\boldsymbol{\omega}^\top \boldsymbol{\chi} - \nu =0}\), for \({\boldsymbol{\chi}\in\mathbb{R}^{m-1}}\). In the linearly non-separable case, a proper nonlinear mapping \(\phi(\cdot)\) with kernel \(K(\boldsymbol{\chi}_i,\boldsymbol{\chi}_j)=\phi(\boldsymbol{\chi}_i)^\top \phi(\boldsymbol{\chi}_j)\) can be found such that \({g(\widehat{\boldsymbol{\chi}})= \text{sgn}(\boldsymbol{\omega}^\top \phi(\widehat{\boldsymbol{\chi}}) - \nu)}\) determines the class of \(\widehat{\boldsymbol{\chi}}\). Agents collaboratively solve the problem by finding the optimal \(\boldsymbol{\omega}~\) and \(\nu\) and optimize the following convex loss [35]: \[\label{eq95svm95cent} \begin{align} \displaystyle f_i(\boldsymbol{\omega},\nu) = \boldsymbol{\omega}^\top \boldsymbol{\omega} + C \sum_{j=1}^{N} \max\{x,0\}^p \end{align}\tag{2}\] with \(p \in \mathbb{N}\) as the smoothness factor (which is typically a finite number), \(C>0\) as the margin size parameter, and \({x=1-l_j( \boldsymbol{\omega}^\top \phi(\boldsymbol{\chi}^i_j)-\nu)}\). The differentiable smooth equivalent of \(f_i\) in Eq. 2 is in the following form (assuming large enough \(\mu>0\)): \[\label{eq95svm95smooth} f_i(\boldsymbol{\omega},\nu)=\boldsymbol{\omega}^\top \boldsymbol{\omega} + C \sum_{j=1}^{N_i} \tfrac{1}{\mu}\log (1+\exp(\mu x)).\tag{3}\] This problem is also known as distributed support-vector-machine (D-SVM) [31], [32].

Distributed Least Squares: In this problem, the idea is to solve the least square problem \(H\mathbf{z}=\mathbf{b}\) in a distributed manner. Every agent/node \(i\) takes measurement \(\mathbf{b}_i \in \mathbb{R}^p\) and has a \(p\)-by-\(n\) measurement matrix \(H_i\) and collaboratively optimizes the private loss function in the following form [32]: \[\begin{align} f_i(\mathbf{x}) = \frac{1}{2} \|H_i \mathbf{z}-\mathbf{b}_i\|_2^2 \label{eq95prob95ls} \end{align}\tag{4}\] This can be addressed further in the context of distributed filtering [36]–[38].

Distributed Logistic Regression: In this problem each agent \(i\) with access to \(m_i\) training data points defined by \((c_{ij},y_{ij}) \in \mathbb{R}^p \times \{-1,1\}\), where the parameter \(c_{ij}\) has \(p\) features of the \(j\)th training data and \(y_{ij}\) denotes the binary label \(\{-1,+1\}\). Each agent, collaborating with others, solves and optimizes the private loss function in the following form [30]: \[\begin{align} f_i(\mathbf{w},b) = \sum_{j=1}^{m_i} \log(1+\exp(-(\mathbf{w}^\top c_{ij}+b)y_{ij})) + \frac{\lambda}{2} \|\mathbf{w}\|_2^2 \label{eq95prob95lr} \end{align}\tag{5}\] where the last term is for regularization to avoid overfitting.

1.3 Paper Organization↩︎

Section 2 provides the preliminary notions. Section 3 gives the main DTAC-ADDOPT algorithm with proof of convergence in Section 4. Section 5 provides the simulation results on both the academic setup and the real dataset. Section 6 provides the concluding remarks.

1.4 Notations↩︎

Table 1 summarizes the notations in this paper.

0.7

2 Preliminaries↩︎

2.1 Algebraic Graph Theory↩︎

We consider the network of agents as a digraph (directed graph) of nodes denoted by \(\mathcal{G}=\{\mathcal{V},\mathcal{E}\}\) with \(\mathcal{V}\) and \(\mathcal{E}\) respectively as the node set and link set. A link \((i,j) \in \mathcal{E}\) from node \(i\) to node \(j\) implies a communication link for message passing from agent \(i\) to agent \(j\). The adjacency matrix of \(\mathcal{G}\) is denoted by \(C\) where \(C_{ij}\) is the weight on the link \((j,i)\) (or \(j \rightarrow i\)). Define the in-neighborhood of every node \(i\) as \(\mathcal{N}_i = \{j|(j,i) \in \mathcal{E}\}\) or \(\mathcal{N}_i = \{j|C_{ij} \neq 0 \}\).

Note that, in most directed network implementations, agents already know their outgoing neighbor set for column-stochastic design of matrices, and the out-degree is locally available.

2.2 Augmented Formulation↩︎

The delay model is similar to the consensus literature [27] and is clearly defined in the following assumption.

We justify the above assumption. Note that, in practice, delays may change on a time-scale much slower than the algorithm step-size (or are upper-bounded by the same constant); therefore, the derived bounds using the maximum delay remain valid in practical cases. Also, in many networks, communication paths and routing remain stable for long periods. Communication latencies in these settings are dominated by propagation and queuing delays that are (on the algorithm time-scale) nearly constant and hence well modeled as time-invariant. Moreover, treating delays as fixed (but heterogeneous) provides a conservative worst-case analysis useful for algorithm design and safety guarantees.

For every set of connected nodes \((i,j)\) and \((i,k)\), the communication delay implies the heterogeneous scenario. Define the augmented state vectors \(\widehat{\mathbf{z}}_k = (\mathbf{z}_k; \mathbf{z}_{k-1}; \dots; \mathbf{z}_{k-\overline{\tau}})\) as the column-concatenation of delayed state vectors (";" denotes column concatenation). Given the column stochastic consensus matrix \(C\) and maximum delay \(\overline{\tau}\), its augmented matrix is defined as, \[\begin{align} \label{eq95aug95C} \small \overline{C} = \left( \begin{array}{cccccc} C_0 & I_n & 0_n & \ldots & 0_n & 0_n \\ C_1 & 0_n & I_n &\ldots & 0_n& 0_n\\ \vdots & \vdots & \vdots & \ddots & \vdots & \vdots \\ C_{\overline{\tau}-1} & 0_n & 0_n & \ldots & I_n & 0_n \\ C_{\overline{\tau}} & 0_n & 0_n & \ldots & 0_n & I_n \end{array} \right), \normalsize \end{align}\tag{6}\] with \(I_{n}\) and \(0_{n}\) respectively as \(n\)-by-\(n\) identity and zero matrices. The non-negative matrices \(C_r\) with \(r \in \{0,\dots,\overline{\tau}\}\) are defined based on delay \(0 \leq r \leq \overline{\tau}\) as, \[\begin{align} C_{r,ij} = \left\{ \begin{array}{ll} C_{ij}, & \text{If}~ \tau_{ij}=r \\ 0, & \text{Otherwise}. \end{array}\right. \end{align}\] Assuming time-invariant delays, for every \((i,j) \in \mathcal{E}\), only one of the entries \(C_{0,ij},C_{1,ij}, \ldots, C_{\overline{\tau},ij}\) is equal to \(C_{ij}\) and the rest are zero. This implies that the column-sum of the first \(n\) columns of \(\overline{C}\) and \(C\) are equal. Note that, given a column-stochastic \(C\) matrix, the augmented matrix \(\overline{C}\) is also column stochastic from the definition. It should be noted that this large augmented matrix is only used in proof analysis of the proposed algorithm, and it is not practically used in the agents’ dynamics (see the iterative dynamics 11 14 in the next section).

3 The Main Algorithm↩︎

The main algorithm in compact matrix form (subject to no time-delay) is given as follows: \[\begin{align} \tag{7} \mathbf{x}_{k+1} &= {C}_k \mathbf{x}_{k} - \alpha \mathbf{g}_k \\ \tag{8} \mathbf{y}_{k+1} &= {C}_k \mathbf{y}_{k} \\ \tag{9} \mathbf{z}_{k+1} &= \frac{\mathbf{x}_{k+1}}{\mathbf{y}_{k+1}} \\ \tag{10} \mathbf{g}_{k+1} &= {C}_k \mathbf{g}_{k} + \nabla \mathbf{f}_{k+1} - \nabla \mathbf{f}_{k} \end{align}\] For delayed case, define the augmented vectors \(\widehat{\mathbf{x}}_k,\widehat{\mathbf{y}}_k,\widehat{\mathbf{g}}_k\) of size \(n(\overline{\tau}+1)\). Let, \(\mathbf{Y}_k = diag(\widehat{\mathbf{y}}_k).\) Further, define the auxiliary matrix \(\Xi^{n}_{i,\overline{\tau}}\) is an \(n \times (\overline{\tau}+1)n\) matrix defined as \(\Xi^{n}_{i,\overline{\tau}}= (\mathbf{b}^{\overline{\tau}+1}_i \otimes I_{n})^\top\) with \(\mathbf{b}^{\overline{\tau}+1}_i\) as the unit column-vector of the \(i\)’th coordinate (\(1\leq i \leq {\overline{\tau}+1}\)), i.e., \(\mathbf{b}^{\overline{\tau}+1}_i =( \underbrace{\overbrace{0;\dots;0}^{i-1};1;0;\dots;0}_{\overline{\tau}+1})\) In case \(\mathbf{x} \in \mathbb{R}^{np}\) then \(\Xi^{np}_{i,\overline{\tau}} = \Xi^{n}_{i,\overline{\tau}} \otimes I_p\). Then, putting \(i=1\), we have \(\mathbf{x}_k = \Xi^{np}_{1,\overline{\tau}} \widehat{\mathbf{x}}_k, \mathbf{y}_k = \Xi^{np}_{1,\overline{\tau}} \widehat{\mathbf{y}}_k, \mathbf{g}_k = \Xi^{np}_{1,\overline{\tau}} \widehat{\mathbf{g}}_k\). In fact, \(\Xi^{np}_{1,\overline{\tau}}\) returns the first \(np\) rows of the augmented vector of size \(np(\overline{\tau}+1)\).

The main distributed optimization dynamics under communication time-delays are in the following vector form, \[\begin{align} \tag{11} \mathbf{x}_{k+1,i} &= \sum_{j \in \mathcal{N}_i} \sum_{r=0}^{\overline{\tau} } C_{k,ij} \mathcal{I}_{k-r,ij}(r) \mathbf{x}_{k-r,j} - \alpha \mathbf{g}_{k,i} \\ \tag{12} \mathbf{y}_{k+1,i} &= \sum_{j \in \mathcal{N}_i} \sum_{r=0}^{\overline{\tau} } C_{k,ij} \mathcal{I}_{k-r,ij}(r) \mathbf{y}_{k-r,j}\\ \tag{13} \mathbf{z}_{k+1,i} &= \frac{\mathbf{x}_{k+1,i}}{\mathbf{y}_{k+1,j}} \\ \tag{14} \mathbf{g}_{k+1,i} &= \sum_{j \in \mathcal{N}_i} \sum_{r=0}^{\overline{\tau} } C_{k,ij} \mathcal{I}_{k-r,ij}(r) \mathbf{g}_{k-r,j} + (\nabla \mathbf{f}_{k+1,i} - \nabla \mathbf{f}_{k,i}) \end{align}\] where \(\nabla \mathbf{f}_{k+1,i}\) denotes \(\nabla \mathbf{f}_{i}(\mathbf{z}_{k+1,i})\) and \(\mathcal{I}\) is the indicator function, \[\begin{align} \label{eq95mcI} \mathcal{I}_{k,ij}(\tau) = \left\{ \begin{array}{ll} 1, & \text{if}~ \tau_{ij}(k) = \tau,\\ 0, & \text{otherwise}. \end{array}\right. \end{align}\tag{15}\] In practice, agents use Eqs. 11 14 to update their states, i.e., the recipient agents use the neighboring data as they arrive with some possible delays. The proposed solution is summarized in Algorithm 1.

Equivalently in the matrix form, \[\begin{align} \tag{16} \widehat{\mathbf{x}}_{k+1} &= \overline{C}_k \widehat{\mathbf{x}}_{k} - \alpha \widehat{\mathbf{g}}_k \\ \tag{17} \widehat{\mathbf{y}}_{k+1} &= \overline{C}_k \widehat{\mathbf{y}}_{k} \\ \tag{18} \mathbf{z}_{k+1} &= \Xi^{np}_{i,\overline{\tau}} \widehat{\mathbf{z}}_{k+1},~ \widehat{\mathbf{z}}_{k+1} = \frac{\widehat{\mathbf{x}}_{k+1}}{\widehat{\mathbf{y}}_{k+1}} \\ \tag{19} \widehat{\mathbf{g}}_{k+1} &= \overline{C}_k \widehat{\mathbf{g}}_{k} + (\overline{\nabla \mathbf{f}}_{k+1} - \overline{\nabla \mathbf{f}}_{k}) \end{align}\] where, \(\overline{\nabla \mathbf{f}}_{k} := (\nabla \mathbf{f}_{k}; \nabla \mathbf{f}_{k-1}; \dots ; \nabla \mathbf{f}_{k-\overline{\tau}})\). Augmented matrix \(\overline{C}_k\) (defined in Section 2) is column-stochastic [27]. For the proposed solution, one can substitute the strong-connectivity of \(\mathcal{G}\) (irreducibility of \(C_k\)) in [5], [13] by the irreducibility of \(\overline{C}_k\). Note that given a column-stochastic consensus matrix \(C\), its augmented consensus version \(\overline{C}_k\) is also column-stochastic.

The following lemma from our previous work [41] relates the spectral property of the delayed and deay-free system matrices.

Define, \[\begin{align} y &:= \sup_k \| \mathbf{Y}_k\| \\ y_- &:= \sup_k \| \mathbf{Y}^{-1}_k\| \end{align}\] and recall the column-stochasticity of \(\overline{C}_k\). Then, the following lemma holds.

Finding an exact relation between \(\sigma\) and \(\sigma_1\) as a function of \(\overline{\tau}\) could be one direction of future research. Here, the relation between \(\sigma=\left\| \overline{C} - \overline{C}_{\infty} \right\|\) (the augmented case) and \(\sigma_1=\left\| C - C_{\infty} \right\|\) (the delay-free case) is approximated in the following lemma.

This lemma implies that \(\sigma \leq \sigma_1^{\frac{1}{1+\overline{\tau}}} < 1\), which says that by increasing \(\overline{\tau}\), \(\sigma\) gets closer to \(1\).

In the absence of delays, from [5] we have, \[\begin{align} \underline{\overline{\mathbf{x}}}_k &:= \frac{1}{n} (\mathbf{1}_n \otimes I_{p})(\mathbf{1}^\top_n \otimes I_{p})\mathbf{x}_k \\ \underline{\overline{\mathbf{g}}}_k &:= \frac{1}{n} (\mathbf{1}_n \otimes I_{p})(\mathbf{1}^\top_n \otimes I_{p})\mathbf{g}_k \\ \mathbf{z}^* &:= \mathbf{1}_n \otimes \underline{\mathbf{z}}^* \\ \underline{\mathbf{h}}_k &:= \frac{1}{n} (\mathbf{1}_n \otimes I_{p})(\mathbf{1}^\top_n \otimes I_{p}) \nabla \mathbf{f}_{k} \\ \underline{\mathbf{q}}_k &:= \frac{1}{n} (\mathbf{1}_n \otimes I_{p})(\mathbf{1}^\top_n \otimes I_{p}) \nabla \mathbf{f}(\underline{\overline{\mathbf{x}}}_k). \end{align}\] and, in the presence of communication time-delays (with max delay \(\overline{\tau}\)), the variables change to the augmented version as \[\begin{align} \overline{\mathbf{x}}_k &:= (\underline{\overline{\mathbf{x}}}_k;\underline{\overline{\mathbf{x}}}_{k-1}; \dots; \underline{\overline{\mathbf{x}}}_{k-\overline{\tau}}) \\ \overline{\mathbf{g}}_k &:= (\underline{\overline{\mathbf{g}}}_k;\underline{\overline{\mathbf{g}}}_{k-1}; \dots; \underline{\overline{\mathbf{g}}}_{k-\overline{\tau}}) \\ \mathbf{z}^* &:= \mathbf{1}_{n(\overline{\tau}+1)} \otimes \underline{\mathbf{z}}^* \\ \mathbf{h}_k &:= (\underline{\mathbf{h}}_k;\underline{\mathbf{h}}_{k-1}; \dots; \underline{\mathbf{h}}_{k-\overline{\tau}}) \\ \mathbf{q}_k &:= (\underline{\mathbf{q}}_k;\underline{\mathbf{q}}_{k-1}; \dots; \underline{\mathbf{q}}_{k-\overline{\tau}}). \end{align}\] Let’s define the following variables for the proof analysis: \[\begin{align} \mathbf{t}_k &:= \left( \begin{array}{c} \|\widehat{\mathbf{x}}_k - \mathbf{Y}_\infty \overline{\mathbf{x}}_k\| \\ \| \overline{\mathbf{x}}_k - \mathbf{z}^*\|_2 \\ \|\widehat{\mathbf{g}}_k - \mathbf{Y}_\infty \mathbf{h}_k\| \end{array} \right), \\ \mathbf{s}_k &:= \left( \begin{array}{c} \|\mathbf{x}_k \|_2 \\ 0 \\ 0 \end{array} \right),\\ \tag{20} G &:= \left( \begin{array}{ccc} \sigma & 0 & \alpha\\ \alpha c l y_- & \eta & 0 \\ c d \epsilon l y_-(\kappa + \alpha l y y_-) & \alpha d \epsilon l^2 y y_- & \sigma + \alpha c d \epsilon l y \end{array} \right),\\ \tag{21} H_k &:= \left( \begin{array}{ccc} 0 & 0 & 0\\ \alpha l y_- T \gamma_1^{k-1} & 0 & 0 \\ (\alpha l y +2)d\epsilon l y_-^2 T \gamma_1^{k-1} & 0 & 0 \end{array} \right), \end{align}\] where \(\kappa := \|\overline{C}_k - I_{np\overline{\tau}}\|_2 = \|C - I_{np}\|_2\), \(\epsilon := \|I_{np\overline{\tau}} - \overline{C}_\infty \|_2 := \|I_{np} - C_\infty \|_2\), \(\eta := \max\{1-n(\overline{\tau}+1)l,1-n(\overline{\tau}+1)s\}\). \(y=\sup_k\|\mathbf{Y}_k\|_2\), \(y_-=\sup_k\|\mathbf{Y}^{-1}_k\|_2\), and \(c,d\) are positive constant from the equivalence of \(\|\cdot\|\) and \(\|\cdot\|_2\). These variables are used in the following lemmas,

It should be clarified that Eq. ?? provides a linear iterative relation between \(t_k\) and \(t_{k+1}\) via matrices, \(G\) and \(H_{k-1}\). Therefore, the convergence of \(t_k\) follows from spectral analysis of matrices \(G\) and \(H\). In other words, to prove linear convergence of \(\|t_k\|_2\) toward zero, the sufficient condition is to prove \(\rho(G)<1\), as well as the linear decay of matrix \(H_{k-1}\), which is straightforward from Eq. 21 since \(0<\gamma_1<1\).

So, we need to prove that \(\rho(G)<1\) as a sufficient condition to bound \(\alpha\) (the spectral radius of \(G\) defined in Eq. 20 being less than \(1\)). This is discussed in the following lemma and proved by matrix perturbation theory.

It should be noted that the bound on \(\alpha\) holds for both heterogeneous delays \(\tau \leq \overline{\tau}\) and homogeneous maximum delays \(\tau = \overline{\tau}\). For both cases, the gradient-tracking rate \(\alpha\) satisfying Eq. ?? ensures the convergence.

The following lemma provides the main results on the linear convergence of Algorithm 1.

It should be noted that large delays may cause considerable computational overhead as the dimension of the augmented matrices scales with the time-delay bound \(\overline{\tau}\). However, this trade-off is inherent to worst-case delay handling in this paper; handling delayed messages explicitly enables delay-tolerant convergence and explicit stability margins for gradient-tracking rate \(\alpha\) (as shown in Eq. ?? ) while, on the other hand, results in higher computational costs for large delays. In this paper, considering reasonable and sufficiently small delay bounds (to avoid packet loss), the convergence analysis and computational complexity are justified.

4 Convergence and Proof of Lemma 5↩︎

This section presents the proof of Lemma 5 and the convergence analysis in three separate steps.

Step-I:

First, from Eq. 16 , Lemma 3 and Lemma 7 we bound \(\|\widehat{\mathbf{x}}_k - \mathbf{Y}_\infty \overline{\mathbf{x}}_k\|\) as, \[\begin{align} \nonumber \|\widehat{\mathbf{x}}_k - \mathbf{Y}_\infty \overline{\mathbf{x}}_k\| \leq& \|\overline{C}_k \widehat{\mathbf{x}}_{k-1} - \mathbf{Y}_\infty \overline{\mathbf{x}}_{k-1}\| \\&+ \alpha \| \widehat{\mathbf{g}}_{k-1} - \mathbf{Y}_\infty \mathbf{h}_{k-1}\| \label{eq95step195final} \\ \nonumber \leq& \sigma \|\widehat{\mathbf{x}}_{k-1} - \mathbf{Y}_\infty \overline{\mathbf{x}}_{k-1}\| \\&+ \alpha \| \widehat{\mathbf{g}}_{k-1} - \mathbf{Y}_\infty \mathbf{h}_{k-1}\| \end{align}\tag{23}\]

Step-II:

Next we bound \(\|\overline{\mathbf{x}}_{k} - \mathbf{z}^*\|_2\). From Lemma 7, \[\begin{align} \label{eq95qh} \overline{\mathbf{x}}_{k} = (\overline{\mathbf{x}}_{k-1} - \alpha \mathbf{q}_{k-1}) - \alpha (\mathbf{h}_{k-1} - \mathbf{q}_{k-1}) \end{align}\tag{24}\] Let’s define \(\mathbf{x}_+ = \overline{\mathbf{x}}_{k-1} - \alpha \mathbf{q}_{k-1}\) as the augmented version of centralized GD step. Redefining Lemma 9 and Eq. ?? for the augmented variables, we get \[\begin{align} \|\mathbf{x}_+ - {\mathbf{z}}^*\|_2 \leq \eta \| \widehat{\mathbf{x}}_{k-1} - {\mathbf{z}}^* \|_2 \end{align}\] For the second term in 24 , from Lipschitz condition we obtain, \[\begin{align} \|\mathbf{h}_{k-1} - \mathbf{q}_{k-1}\|_2 \leq \|\frac{1}{n(\overline{\tau}+1)} (\mathbf{1}_{n(\overline{\tau}+1)} \mathbf{1}^\top_{n(\overline{\tau}+1)}) \otimes I_{p})\|_2 l \|\widehat{\mathbf{z}}_{k-1} - \overline{\mathbf{x}}_{k-1}\|_2 \end{align}\] Then, \[\begin{align} \nonumber \|\overline{\mathbf{x}}_{k} - \mathbf{z}^*\|_2 &\leq \|\mathbf{x}_+ - {\mathbf{z}}^*\|_2 + \alpha l \|\mathbf{h}_{k-1} - {\mathbf{q}}_{k-1}\|_2 \\ &\leq \eta \| \widehat{\mathbf{x}}_{k-1} - {\mathbf{z}}^* \|_2 + \alpha l \| \widehat{\mathbf{z}}_{k-1} - \overline{\mathbf{x}}_{k-1} \|_2 \label{eq:proofxz} \end{align}\tag{25}\] From Eq. 13 (or Eq. 18 ) and recalling Lemma 8 we get, \[\begin{align} \nonumber \|\widehat{\mathbf{z}}_{k-1} - \overline{\mathbf{x}}_{k-1}\|_2 \leq& \|\mathbf{Y}^{-1}_{k-1}(\widehat{\mathbf{x}}_{k-1} - \mathbf{Y}_\infty \overline{\mathbf{x}}_{k-1})\|_2 \\\nonumber &+ \|\mathbf{Y}^{-1}_{k-1}\mathbf{Y}_\infty - I_{np(\overline{\tau}+1)})\overline{\mathbf{x}}_{k-1}\|_2 \\ \nonumber \leq& y_-\|\widehat{\mathbf{x}}_{k-1} - \mathbf{Y}_\infty \overline{\mathbf{x}}_{k-1}\|_2 \\&+ y_-T \gamma_1^{k-1} \|\widehat{\mathbf{x}}_{k-1}\|_2 \label{eq:step2semifinal} \end{align}\tag{26}\] where we also used the fact that \(\|\overline{x}_{k-1}\|_2 \leq \|\widehat{x}_{k-1}\|_2\). Then, by substituting the above in Eq. 25 we get, \[\begin{align} \nonumber \|\overline{\mathbf{x}}_{k} - \mathbf{z}^*\|_2 \leq& \alpha c l \mathbf{y}_- \|\widehat{\mathbf{x}}_{k-1} \mathbf{Y}_\infty \overline{\mathbf{x}}_{k-1} \|_2 \\ &+ \eta \|\overline{\mathbf{x}}_{k-1} - \mathbf{z}^*\|_2 + \alpha l \mathbf{y}_- T \gamma_1^{k-1} \| \widehat{\mathbf{x}}_{k-1} \|_2 \label{eq95step295final} \end{align}\tag{27}\]

Step-III:

Next, we bound \(\|\widehat{\mathbf{g}}_{k} - \mathbf{Y}_\infty \mathbf{h}_k\|_2\). From Eq. 19 \[\begin{align} \nonumber \|\widehat{\mathbf{g}}_{k} - \mathbf{Y}_\infty \mathbf{h}_k\|_2 &\leq \|\overline{C}_k \widehat{\mathbf{g}}_{k-1} - \mathbf{Y}_\infty \mathbf{h}_{k-1}\|_2 \\ \label{eq95gh} &+ \|(\overline{\nabla \mathbf{f}}_k - \overline{\nabla \mathbf{f}}_{k-1} - \mathbf{Y}_\infty( \mathbf{h}_{k} - \mathbf{h}_{k-1})\|_2 \end{align}\tag{28}\] From Lemma 3 and Lemma 7, \[\begin{align} \|\overline{C}_k \widehat{\mathbf{g}}_{k-1} - \mathbf{Y}_\infty \mathbf{h}_{k-1}\|_2 \leq \sigma \|( \widehat{\mathbf{g}}_{k-1} - \mathbf{Y}_\infty \overline{\mathbf{g}}_{k-1})\|_2 \end{align}\] Further, the second term in 28 can be recalculated as, \[\begin{align} \nonumber &\|(\overline{\nabla \mathbf{f}}_k - \overline{\nabla \mathbf{f}}_{k-1}) - \mathbf{Y}_\infty( \mathbf{h}_{k} - \mathbf{h}_{k-1})\|_2 =\\ \nonumber &\|(I_{np(\overline{\tau}+1)} - \frac{\mathbf{Y}_\infty}{n}(\mathbf{1}_{n(\overline{\tau}+1)} \otimes I_{p})(\mathbf{1}^\top_{n(\overline{\tau}+1)} \otimes I_{p})) (\overline{\nabla \mathbf{f}}_k - \overline{\nabla \mathbf{f}}_{k-1}) \|_2 \\ &\leq \epsilon l \|\widehat{\mathbf{z}}_k - \widehat{\mathbf{z}}_{k-1}\|_2 \end{align}\] which follows from the Lipschitz condition. Therefore, \[\begin{align} \label{eq95gh2} \|\widehat{\mathbf{g}}_{k} - \mathbf{Y}_\infty \mathbf{h}_k\|_2 &\leq \sigma \|\widehat{\mathbf{g}}_{k-1} - \mathbf{Y}_\infty \mathbf{h}_{k-1}\|_2 + d \epsilon l \|\widehat{\mathbf{z}}_k - \widehat{\mathbf{z}}_{k-1}\|_2 \end{align}\tag{29}\] To bound \(\|\widehat{\mathbf{z}}_k - \widehat{\mathbf{z}}_{k-1}\|_2\) we have, \[\begin{align} \nonumber \|{\mathbf{h}}_{k-1} \|_2 &= \|\frac{1}{n(\overline{\tau}+1)}(\mathbf{1}_{n(\overline{\tau}+1)} \otimes I_{p})(\mathbf{1}^\top_{n(\overline{\tau}+1)} \otimes I_{p}\nabla \mathbf{f}(\overline{\mathbf{x}}_{k-1}) \|_2 \\ &\leq l \|\overline{\mathbf{x}}_{k-1} - \mathbf{z}^*\|_2 \end{align}\] Therefore, using Eq. 26 , we obtain, \[\begin{align} \nonumber \|\mathbf{Y}_k^{-1} \widehat{\mathbf{g}}_{k-1}\|_2 \leq& y_-\| \widehat{\mathbf{g}}_{k-1} - \mathbf{Y}_\infty \mathbf{h}_{k-1} \|_2 \\\nonumber &+ y_- y l \| \overline{\mathbf{x}}_{k-1} - \mathbf{z}^* \|_2 \\\nonumber &+ y_-^2 y l \| \widehat{\mathbf{x}}_{k-1} - \mathbf{Y}_\infty \overline{\mathbf{x}}_{k-1} \|_2 \\ &+ y_-^2 y l T \gamma_1^{k-1}\| \widehat{\mathbf{x}}_{k-1} \|_2 \end{align}\] Recall that \((\overline{C}-I_{n(\overline{\tau}+1)p})\mathbf{Y}_\infty \overline{\mathbf{x}}_{k-1} = \mathbf{0}\). Then, \[\begin{align} \nonumber \|\widehat{\mathbf{z}}_k - \widehat{\mathbf{z}}_{k-1}\|_2 \leq& (y_- \kappa + \alpha y_-^2 y l) \|\widehat{\mathbf{x}}_{k-1} - \mathbf{Y}_\infty \overline{\mathbf{x}}_{k-1}\|_2 \\\nonumber &+ \alpha y_- \|\widehat{\mathbf{g}}_{k-1} - \mathbf{Y}_\infty {\mathbf{h}}_{k-1}\|_2 \\\nonumber &+ \alpha y_- y l \|\overline{\mathbf{x}}_{k-1} - {\mathbf{z}}^*\|_2 \\ &+ (\alpha y l +2) y_-^2 T \gamma_1^{k-1} \|\widehat{\mathbf{x}}_{k-1}\|_2 \end{align}\] Substitute the above in Eq. 29 , \[\begin{align} \nonumber \|{\mathbf{g}}_{k} - \mathbf{Y}_\infty \mathbf{h}_k\|_2 \leq& (c d \epsilon l \kappa y_- + \alpha c d \epsilon l^2 y y_-^2) \|\widehat{\mathbf{x}}_{k-1} \\ \nonumber &- \mathbf{Y}_\infty \overline{\mathbf{x}}_{k-1}\|_2 + \alpha d \epsilon l^2 y y_- \|\overline{\mathbf{x}}_k - \mathbf{z}^*\|_2 \\ \nonumber &+ (\sigma + \alpha c d \epsilon l y_-) \|\widehat{\mathbf{g}}_{k-1} - \mathbf{Y}_\infty \mathbf{h}_{k-1}\|_2 \\ &+ (\alpha y l +2)d \epsilon l y_-^2 T \gamma_1^{k-1} \|\widehat{\mathbf{x}}_{k-1}\|_2 \label{eq95step395final} \end{align}\tag{30}\] Finally, combining Eqs. 23 , 27 , and 30 results in Lemma 5 and proves the convergence.

5 Simulations↩︎

5.1 Academic Example↩︎

For the experimental simulation, we consider a quadratic cost function as Eq. 5 similar to [44] with randomly set parameters. The number of agents is set as \(n=10\) nodes. The bound on the time-delay is set \(\overline{\tau}=5\) and \(\alpha=0.005\). We consider convergence over two cases of random Erdos-Renyi (ER) networks: (i) static networks where the structure of the multi-agent network is time-invariant, and (ii) dynamic (switching) networks where the network structure randomly changes every \(2\) iterations. The simulations are shown in Fig. 2. For the switching case, there exist some oscillations in the residual decay due to changes in the network topology.

Next, we redo the simulations over an ER network to check the convergence for different values of bound on the time-delays, i.e., \(\overline{\tau}= 5,10,15,20\). We set gradient-tracking rate \(\alpha = 0.001\) and \(\alpha = 0.005\) for this simulation. The mean-square-error (MSE) residuals at agents for different bounds on the time-delay are shown in Fig. 3. As it can be seen from the figure, for large value of \(\overline{\tau}\), the residual decay becomes unstable and the optimization diverges (see the residual for \(\overline{\tau}= 15,20\) in the right figure for \(\alpha = 0.005\)).

5.2 Real Data-Set Example↩︎

We use the MNIST dataset for distributed optimization, which is a well-known dataset in the field of machine learning and image classification. It consists of handwritten digits from \(0\) to \(9\) and is commonly used to train and test various classification algorithms. The dataset includes \(70000\) images of handwritten digits. Each image is a \(28 \times 28\) grayscale image, resulting in \(784\) pixels per image, and is associated with a label from \(0\) to \(9\), indicating the digit it represents. A set of sampled images is shown in Fig. 4.

The data set and image processing algorithms are taken from [45]. We randomly select \(N = 12000\) labelled images from the MNIST data set to be classified using logistic regression with a convex regularizer. The data are distributed among \(n=16\) agents to be cooperatively classified. In our cost optimization setup, define \[\begin{align} \label{eq95F95mnist} \min_{\mathbf{b},c} & F(\mathbf{b},c) = \frac{1}{n}\sum_{i=1}^{n} f_i(\mathbf{x}) \end{align}\tag{31}\] with every node \(i\) taking a batch of \(m_i=750\) sample images. Each node \(i\), then, locally minimizes the following classification cost: \[\begin{align} \label{eq95fij95regression} f_i(\mathbf{x}) = \frac{1}{m_i}\sum_{j=1}^{m_i} \ln(1+\exp(-(\mathbf{b}^\top x_{i,j}+c)y_{i,j}))+\frac{\lambda}{2}\|\mathbf{b}\|_2^2. \end{align}\tag{32}\] with \(\mathbf{b},c\) as the classifier parameters. The residual is defined as \(F(\overline{\mathbf{x}}^k)-F(\mathbf{x}^*)\) with \(\overline{\mathbf{x}}^k = \frac{1}{n} \sum_{i=1}^{n}\mathbf{x}^k_i\). We run and compare the residual of distributed training for different existing distributed optimization techniques in the literature over an exponential network. The following optimization algorithms are used for comparison: GP [13], SGP [46], [47], S-ADDOPT [6], and PushSAGA [48]. The simulation results are given in Fig. 5 for an exponential graph of \(n=16\) nodes (the graph is shown in the figure). It should be mentioned that GP, SGP, S-ADDOPT, and Push-SAGA are not delay-tolerant and, thus, are simulated for delay-free case. Therefore, as expected, they show better performance in the absence of time-delays, while practically they do not converge in the presence of time-delays. On the other hand, our DTAC-ADDOPT algorithm converges in the presence of heterogeneous time-delays. For this simulation, we set \(\overline{\tau}=3\). The slower rate of convergence for DTAC-ADDOPT is due to time-delays in the data sharing as compared to the other delay-free optimization techniques.

6 Conclusions and Future Works↩︎

Delay-tolerant distributed optimization over digraphs is proposed in this work. We present a distributed algorithm over a multi-agent network that is robust to time-delayed information-exchange among the agents. The delay-tolerance is shown both by mathematical proofs and experimental simulations. Future research direction includes finding a tighter bound between \(\sigma_1\) and \(\sigma\) based on \(\overline{\tau}\) in Lemma 3. One can extend the convergence analysis to find maximum delay \(\tau_{\max}\) for which the algorithm fails to converge when \(\tau > \tau_{\max}\). Our analysis is based on bounded delays, considering no packet loss over the network, where the extension to certain classes of packet losses via standard buffering/retransmission or stochastic-analysis variants are left for future research. Moreover, distributed optimization subject to asynchronous and event-triggered operation, privacy-preserving distributed optimization [49], and adding nonlinearities and momentum terms to reach faster convergence (similar to coupling-constrained optimization and resource allocation in [50]) are other open problems and directions of future research. Different applications in machine learning setups can also be considered for future research.

References↩︎