Decentralized Sporadic Federated Learning: A Unified Algorithmic Framework with Convergence Guarantees

[Weakness 2] Gradient Diversity Assumptions for Convex and Non-Convex Models

We have opted to use these data heterogeneity assumptions for convex and non-convex models based on standard assumptions for each of these models that are used in other papers, e.g., see [1-4]. In the revised version of our paper, we have added a subsection in Appendix C to further elaborate on our choice of data heterogeneity assumption for convex and non-convex models, a summary of which is given below.

We note that if we are dealing with strongly-convex models, i.e., we make Assumption 4.1-(b) for $\mu$-strongly convex models, the inequality in Assumption 4.2-(b) for data heterogeneity is a much stronger assumption than 4.1-(c). To show this, we first note that making the $\beta$-smoothness assumption (made in both 4.1-(a) and 4.2-(a)) implies that $\\| \nabla{F}(\theta) \\| \le \beta \\| \theta - \theta^\star \\|,$ where $\theta^\star$ is a point where $\nabla{F}(\theta^\star) = 0$, i.e., a globally optimal point for convex models and a stationary point for non-convex models (note that we present and prove this result in Lemma D.1-(b)). Then, since we have $\\| \nabla{F}_i(\theta) \\| \le \delta_i + \zeta_i \\| \nabla{F}(\theta) \\|$ according to Assumption 4.2-(b), we can conclude that $\\| \nabla{F}(\theta) - \nabla{F}_i(\theta) \\| \le \\| \nabla{F}(\theta) \\| + \\| \nabla{F}_i(\theta) \\| \le \delta_i + (1 + \zeta_i ) \\| \nabla{F}(\theta) \\| \le \delta_i + (1 + \zeta_i) \beta \\| \theta - \theta^\star \\|.$ We observe that defining $\delta’_i = \delta_i$ and $\zeta’_i = (1 + \zeta_i) \beta > 0$, Assumption 4.1-(c) is implied with parameters $\delta’_i$ and $\zeta’_i$.

The converse is not necessarily true though, i.e., Assumption 4.1-(c) does not imply Assumption 4.2-(b). Note that some existing works assume stricter assumptions than both 4.1-(c) and 4.2-(b), and they thus might make a common assumption among convex and non-convex models. These would be specific cases of our Assumptions 4.1-(c) and 4.2-(b). For instance, setting $\zeta_i = 0$ in our Assumption 4.1-(c) gives us the assumption made in [2, 5, 6], and setting $\zeta_i = 0$ in our Assumption 4.2-(b) provides us with the assumptions made in [7, 8].

[1] Koloskova, Anastasia, Nicolas Loizou, Sadra Boreiri, Martin Jaggi, and Sebastian Stich. "A unified theory of decentralized sgd with changing topology and local updates." In International Conference on Machine Learning, pp. 5381-5393. PMLR, 2020.

[2] Wang, Shiqiang, Tiffany Tuor, Theodoros Salonidis, Kin K. Leung, Christian Makaya, Ting He, and Kevin Chan. "Adaptive federated learning in resource constrained edge computing systems." IEEE journal on selected areas in communications 37, no. 6 (2019): 1205-1221.

[3] Bottou, Léon, Frank E. Curtis, and Jorge Nocedal. "Optimization methods for large-scale machine learning." SIAM review 60, no. 2 (2018): 223-311.

[4] Lin, Frank Po-Chen, Seyyedali Hosseinalipour, Sheikh Shams Azam, Christopher G. Brinton, and Nicolo Michelusi. "Semi-decentralized federated learning with cooperative D2D local model aggregations." IEEE Journal on Selected Areas in Communications 39, no. 12 (2021): 3851-3869.

[5] Bornstein, Marco, Tahseen Rabbani, Evan Z. Wang, Amrit Bedi, and Furong Huang. "SWIFT: Rapid Decentralized Federated Learning via Wait-Free Model Communication." In The Eleventh International Conference on Learning Representations.

[6] Chen, Yiming, Kun Yuan, Yingya Zhang, Pan Pan, Yinghui Xu, and Wotao Yin. "Accelerating gossip SGD with periodic global averaging." In International Conference on Machine Learning, pp. 1791-1802. PMLR, 2021.

[7] Koloskova, Anastasiia, Tao Lin, Sebastian Urban Stich, and Martin Jaggi. "Decentralized deep learning with arbitrary communication compression." In Proceedings of the 8th International Conference on Learning Representations. 2019.

[8] Sun, Tao, Dongsheng Li, and Bao Wang. "Decentralized federated averaging." IEEE Transactions on Pattern Analysis and Machine Intelligence 45, no. 4 (2022): 4289-4301.

Thank you for your positive assessment of our paper, including our convergence results and experiments. We also appreciate you commending the novelty of our explored idea and the clarity of our presentation. Below, we respond to each comment alongside discussing how we have addressed them in the revised manuscript.

[Weakness 1] Double stochasticity of $P^{(k)}$

We apologize for not explaining all the steps required to show that this is true. The transition matrix $P^{(k)}$ being doubly stochastic and symmetric follows from the fact that $\hat{v}\_{ij}^{(k)} = \hat{v}\_{ji}^{(k)}$ and $r_{ij} = r_{ji}$ for $i,j \in \mathcal{M}$ and $i \neq j$, and $\hat{v}_{ii}^{(k)} = 0$ for all $i \in \mathcal{M}$ and $k \geq 0$. We thank you for pointing this out, and we have added this information in Sec. 3.1 of our revised manuscript.

Having these pieces at hand, it remains to show that $p_{ij}^{(k)} = p_{ji}^{(k)}$ for all $i \neq j$, and $\sum_{j \in \mathcal{M}}{p_{ij}^{(k)}} = 1$ (which will then imply that $\sum_{i \in \mathcal{M}}{p_{ij}^{(k)}} = 1$ as well). Since $p_{ij}^{(k)} = r_{ij} \hat{v}\_{ij}^{(k)}$, symmetry of $r_{ij}$ and $\hat{v}\_{ij}^{(k)}$ immediately imply the symmetry of $p_{ij}^{(k)}$. Then, we have defined $p_{ii}^{(k)} = 1 - \sum_{j \in \mathcal{M}}{r_{ij} \hat{v}\_{ij}^{(k)}}$ in such a way that we can get $\sum_{j \in \mathcal{M}}{p_{ij}^{(k)}} = 1$. Noting that $\hat{v}\_{ii}^{(k)} = 0$, this in fact implies that $p_{ii}^{(k)} = 1 - \sum_{j \in \mathcal{M}, j \neq i}{p_{ij}^{(k)}}$.

[Weakness 3] Frequency of SGDs under diminishing learning rate policy

If we do not increase the frequency of SGDs $d_i^{(k)}$ when employing a diminishing learning rate $\alpha^{(k)}$, the final model will still converge to a neighborhood of the optimal solution similar to when a constant learning rate $\alpha$ is employed, but not to the exact solution. The intuitive reasoning behind the necessity of increasing the frequency of SGDs $d_i^{(k)}$ is that any DGD-like algorithm designed to solve the optimization problem in Eq. (1) exactly, needs all clients to participate equally. Because otherwise the dataset of a client $i$ with higher SGD frequency influences the final model much more than the dataset of a client $j$ with lower participation, i.e., $d_i^{(k)} > d_j^{(k)}$.

The reason that we do not need such a requirement when a constant learning rate $\alpha$ is being employed, is that using a constant learning rate will already culminate in an inexact solution in the neighborhood of the optimal one. In other words, although an increasing SGD probability will help to reduce the size of this neighborhood, it cannot make it become zero. Thus, we have preferred to keep the algorithm more simple in the case of a constant learning rate.

Regarding the dependence the SGD probabilities $d_i^{(k)}$ on the learning rate $\alpha^{(k)}$, note that in fact they can be any increasing function and do not need to necessarily depend on the learning rate. However, in order to get the conventional convergence rate $\mathcal{O}(\ln(k) / \sqrt{k})$ known in the literature, we assume it increases with the same rate that the learning rate decreases, and we assume a diminishing policy for the learning rate as $\alpha^{(k)} = \alpha^{(0)} / \sqrt{1 + k / \gamma}$. Thus, in general if we did not care about obtaining the theoretical convergence rate, it is not necessary in practice to relate $d_i^{(k)}$ and $\alpha^{(k)}$ together in any way, and exact convergence can be achieved as long as $d_i^{(k)}$ are increasing and $\alpha^{(k)}$ is decreasing. To be precise, we would only need $\alpha^{(k)}$ to satisfy the constraints in the setup of Lemma J.4, and $d_i^{(k)}$ to satisfy $\sum_{k=0}^\infty{\alpha^{(k)} (1 - d_i^{(k)})} < \infty$.

[Question 1] Definition of $p_{ii}^{(k)}$

We note that since $\hat{v}\_{ii}^{(k)} = 0$ for all $i \in \mathcal{M}$ and $k \geq 0$ in DSpodFL, the update suggested by reviewer reduces to what we already have in the paper. We have added an explanation in Sec. 3.1 of our revised manuscript on the symmetry of $r_{ij}$ and $\hat{v}\_{ij}^{(k)}$, and the fact that $\hat{v}\_{ii}^{(k)} = 0$. Note that in our response to [Weakness 1] Double stochasticity of $P^{(k)}, we also discussed about this matter.

Again, thank you for your time and efforts in reviewing our paper. We would appreciate further opportunities to answer any remaining concerns you might have.

Thank you for your follow up. We answer this question in two parts: (i) why we use a different data heterogeneity assumption for non-convex models, and (ii) how this assumption captures data heterogeneity/gradient diversity.

(i) Different data heterogeneity assumption for non-convex models: While the conventional assumption in the literature for data heterogeneity of convex models is $\\| \nabla{F}(\theta) - \nabla{F}_i(\theta) \\| \le \delta_i$ [2, 5, 6] (as referenced in our response to [Weakness 2] previously), our Assumption 4.1-(c) for strongly-convex models includes an extra term, i.e., $\\| \nabla{F}(\theta) - \nabla{F}_i(\theta) \\| \le \delta_i + \zeta_i \\| \theta - \theta^\star \\|$. This is less strict than letting $\zeta_i = 0$, but uses the value of the globally optimal point $\theta^\star \in \mathbb{R}^n$ in its formulation, where $\nabla{F}(\theta^\star) = 0$. When dealing with non-convex models, $\nabla{F}(\theta^{‘\star}) = 0$ will in general hold over a set $\mathcal{J}$ of stationary points $\theta^{‘\star}$ with $\theta^{‘\star} \in \mathcal{J}$, and $|\mathcal{J}| > 1$. Consequently, the data heterogeneity assumption of 4.1-(c), which uses the value of the single globally optimal point $\theta^\star$ for strongly-convex models, will no longer be well-defined for non-convex models, since it is not guaranteed that there will only be one $\theta^\star$ in $\mathbb{R}^n$.

For this reason, if we wanted to use the same assumption for both convex and non-convex models, we would have to make the stricter assumption of $\\| \nabla{F}(\theta) - \nabla{F}_i(\theta) \\| \le \delta_i$, i.e., $\zeta_i = 0$ in our Assumption 4.1-(c) as well. In fact, we can show that if $\zeta_i = 0$, Assumption 4.1-(c) implies Assumption 4.2-(b). Specifically, $\\| \nabla{F}_i(\theta) \\| \le \\| \nabla{F}_i(\theta) - \nabla{F}(\theta) \\| + \\| \nabla{F}(\theta) \\| \le \delta_i + \\| \nabla{F}(\theta) \\| \le \delta_i + \zeta_i’ \\| \nabla{F}(\theta) \\|,$ where $\zeta_i’ \geq 1$. Consequently, the strict assumption of $\\| \nabla{F}(\theta) - \nabla{F}_i(\theta) \\| \le \delta_i$ for non-convex models used in [2, 5, 6] (as referenced in our response to [Weakness 2] previously), implies our Assumption 4.2-(b) for the values of $\zeta_i > 1$. Thus, by making this implication as our direct assumption in Assumption 4.2-(b), we are essentially making a less strict assumption. Furthermore, allowing $\zeta_i > 0$ to be flexible for each client $i$ and not restricting them to $\zeta_i \geq 1$ better captures heterogeneity across clients and gives us the possibility to arrive at more generalized results.

(ii) Assumption 4.2-(b) capturing data heterogeneity/gradient diversity: Let us now make Assumption 4.2-(b), which is to say $\\| \nabla{F}_i(\theta) \\| \le \delta_i + \zeta_i \\| \nabla{F}(\theta) \\|$. Then, we can calculate a measure of data heterogeneity/gradient diversity in the usual form as

$\\| \nabla{F}(\theta) - \nabla{F}\_i(\theta) \\| = \\| \frac1m \sum_{j=1, j \neq i}^m{\nabla{F}\_j(\theta)} - (1 - \frac1m) \nabla{F}\_i(\theta) \\| \le \frac1m \sum_{j=1, j \neq i}^m{\\| \nabla{F}\_j(\theta) \\|} + (1 - \frac1m) \\| \nabla{F}\_i(\theta) \\|$ $\le \frac1m \sum_{j=1, j \neq i}^m{(\delta_j + \zeta_j \\| \nabla{F}(\theta) \\|)} + (1 - \frac1m) (\delta_i + \zeta_i \\| \nabla{F}(\theta) \\|) \le 2 (1 - \frac1m) (\delta + \zeta \\| \nabla{F}(\theta) \\|).$

where we used the fact that $F(\theta) = \frac1m \sum_{i=1}^m{F_i(\theta)}$, $\delta = \max_{i \in \mathcal{M}}{\delta_i}$ and $\zeta = \max_{i \in \mathcal{M}}{\zeta_i}$. Hence, we observe how $\delta_i$ and $\zeta_i$ translate to measures of data heterogeneity/gradient diversity. We can again see that our Assumption 4.2-(b) implies a more relaxed assumption than $\\| \nabla{F}(\theta) - \nabla{F}_i(\theta) \\| \le \delta_i$, since we have $\zeta_i > 0$ in general . The reviewer is correct in saying that $\zeta_i = 0$ in our Assumption 4.2-(b) reduces to the bounded gradient assumptions of [7, 8] (as referenced in our response to [Weakness 2] previously), and in this case, we have essential made a stricter assumption which does not capture the data heterogeneity precisely. As seen above, having $\zeta_i \neq 0$ more generally allows us to capture it.

Thank you again. Please let us know if you have any remaining concerns we can address.

We appreciate your comments on our paper. Thank you for acknowledging the strength of our framework and theoretical guarantees. Please find our responses below.

[Weakness 1] Matching Existing Convergence Rates Using DSpodFL

• First, we note that as outlined in Table 1 of Sec. 1, we provide convergence for last iterates of model parameters when dealing with strongly-convex models in our paper, in contrast to the majority of existing literature which only provide convergence for average iterates. Therefore when dealing with strongly-convex models, we can compare our theoretical results only for the DGD baseline because the DGD-like algorithms given in [1, 2] are among the few to show convergence for the last iterates of model parameters as well. According to Theorems 3.6, 5.5, 5.7 and D.1 in [1] and Theorems 3.5, 4.5 and B.2 in [2], the convergence rate of the algorithms ProxSkip, Decentralized Scaffnew, SplitSkip [1], GradSkip, GradSkip+ and VR-GradSkip+ [2] are all geometric, i.e., $\mathcal{O}(\rho^K)$ with $0 < \rho < 1$ (note that most of these algorithms are FL algorithms, and not decentralized FL algorithms). This rate agrees with the rate we provide in Eq. (10) of Theorem 4.11.
  
  When dealing with non-convex models, we can compare our rate with all the baselines, i.e., [3] for DGD, [4] for RG and [5] for DFedAvg. Since DGD is a special case of RG with no sporadicity in aggregations, i.e., $b_{ij} = 1$ for all $(i,j) \in \mathcal{E}$, we will compare our work with the more recent paper [4]. For DGD and RG, Lemma 17 of [4] shows the convergence upper bound before tuning the constant learning rate $\alpha$, which is $\mathcal{O}{\left( \frac{\mathbb{E}{[F{( \bar{\theta}^{(0)} )}] - F^\star}}{\alpha (K+1)} \right)} + \mathcal{O}{(\alpha)} + \mathcal{O}{(\alpha^2)}.$ For DFedAvg, Theorem 1 in [5] with a zero momentum ($\theta = 0$) obtains a less tight bound as follows $\mathcal{O}{\left( \frac{\mathbb{E}{[F{( \bar{\theta}^{(0)} )}] - F^\star}}{\alpha (K+1)} \right)} + \mathcal{O}{(\alpha)} + \mathcal{O}{(\alpha^2)} + \mathcal{O}{(\alpha^3)}.$ We observe that by setting $d_{min} = 1$ in Theorem 4.12 of our paper, these convergence rates are recovered. Recall that as discussed in Fig. 1 of our paper, all of these baseline algorithms fit within the general framework of DSpodFL with $d_i = 1$ for all clients $i \in \mathcal{M}$. That is why we can substitute $d_{min} = 1$ in Theorem 4.12 to compare our analytical results with the ones provided in [3-5].

• For a step-by-step comparison of our bounds with the one given in [4] for non-convex models, let us examine the bound derived in the proof of Lemma 16 in [4] before tuning the learning rate, i.e., $\frac1{2(T+1)} \sum_{t=0}^T{{\\| \nabla{f}(\bar{x}^{(t)}) \\|}_2^2} \le \frac{\mathbb{E}{f(\bar{x}^{(0)})} - f^\star}{(T+1) \eta} + \frac{L \hat{\sigma}^2}{n}\eta + 64 \frac{L^2 [2 \hat{\sigma}^2 + 2(\frac{6 \tau}{p} + M) \hat{\zeta}^2] \tau}{p} \eta^2.$

Translating these parameters to the setup of our paper, we have $T \to K, t \to r, f \to F, x \to \theta, \eta \to \alpha, L \to \beta, \hat{\sigma} \to \sigma, n \to m, \tau \to 1, p \to 1 - \tilde{\rho}, M \to 0.$

As a consequence, the bound in [4] using the notation and setup of our paper becomes $\frac1{2(K+1)} \sum_{r=0}^K{{\\| \nabla{F}(\bar{\theta}^{(r)}) \\|}^2} \le \frac{\mathbb{E}{F(\bar{\theta}^{(0)})} - F^\star}{(K+1) \alpha} + \frac{\beta \sigma^2}{m} \alpha + 128 \frac{\beta^2 [\sigma^2 + \frac{6}{1 - \tilde{\rho}} \delta^2]}{1-\tilde{\rho}} \alpha^2.$

Now, we compare this with the non-asymptotic bound we obtain for non-convex models in Theorem 4.12, which is $\frac{w_1}{K+1} \sum_{r=0}^K{{\\| \nabla{F}(\bar{\theta}^{(r)}) \\|}^2} \le \frac{F(\bar{\theta}^{(0)}) - F^\star}{\alpha (K+1)} + \frac{1 + \Gamma_3}{1-\frac1{\Gamma_1}} \frac{\beta^2}{m (1 - \tilde{\rho})} \frac{{\\| \mathbf{\Theta}^{(0)} - \mathbf{1}\_m \bar{\mathbf{\theta}}^{(0)} \\|}^2}{K+1} + \frac{1 + \Gamma_3}{1-\frac1{\Gamma_1}} \frac{\beta^2 (16 \frac{1+\tilde{\rho}}{1 - \tilde{\rho}} \delta^2 + \sigma^2)}{1 - \tilde{\rho}} \alpha^2 + (1+\Gamma_3) (1-d_{min}) \delta^2 + \frac{\beta \sigma^2}{2m} \alpha.$ where we used the fact that $d_{max} \le 1$.

We can see how our bound compares with the one derived in [4]. The exact value of $w_1$ in our paper depends on the arbitrarily chosen scalars $\Gamma_0, …, \Gamma_4$, but we have in general that $w_1 \le \frac12$. For common terms that our analysis has with [4], we see how their coefficients are also very similar. For example, for the term proportional to $\alpha^2$, we observe that both our bound and the one in [4] are proportional to $\frac{\beta^2}{(1-\tilde{\rho})}$ and $(\sigma^2 + \frac{q}{1-\tilde{\rho}} \delta^2)$, with slight differences in the value of $q$. Furthermore, for the term proportional to $\alpha$, both our bound and [4] are proportional to $\frac{\beta \sigma^2}{2m}$.

However, we can observe that our bound consists of two extra terms compared to [4]. The first one, $\frac{1 + \Gamma_3}{1-\frac1{\Gamma_1}} \frac{\beta^2}{m (1 - \tilde{\rho})} \frac{{\\| \mathbf{\Theta}^{(0)} - \mathbf{1}\_m \bar{\mathbf{\theta}}^{(0)} \\|}^2}{K+1}$, is capturing the effect of different model initializations for the clients. The second one, $(1+\Gamma_3) (1-d_{min}) \delta^2$ is capturing the effect of sporadic SGDs in our DSpodFL framework. We observe that if all clients conduct SGDs at every iteration, i.e., $d_i^{(k)} = 1$ for all $i \in \mathcal{M}$ and $k \geq 0$, this term becomes equal to zero, giving us the bound for non-sporadic methods outlined in [4]. Therefore, as we have claimed in Sec. 4.5 of our paper, our convergence bounds recover well-known results in the literature in the degenerate case of $d_{min} = 1$.

[1] Mishchenko, Konstantin, Grigory Malinovsky, Sebastian Stich, and Peter Richtárik. "Proxskip: Yes! local gradient steps provably lead to communication acceleration! finally!." In International Conference on Machine Learning, pp. 15750-15769. PMLR, 2022.

[2] Maranjyan, Artavazd, Mher Safaryan, and Peter Richtárik. "Gradskip: Communication-accelerated local gradient methods with better computational complexity." arXiv preprint arXiv:2210.16402 (2022).

[3] Nedic, Angelia, and Asuman Ozdaglar. "Distributed subgradient methods for multi-agent optimization." IEEE Transactions on Automatic Control 54, no. 1 (2009): 48-61.

[4] Koloskova, Anastasia, Nicolas Loizou, Sadra Boreiri, Martin Jaggi, and Sebastian Stich. "A unified theory of decentralized sgd with changing topology and local updates." In International Conference on Machine Learning, pp. 5381-5393. PMLR, 2020.

[5] Sun, Tao, Dongsheng Li, and Bao Wang. "Decentralized federated averaging." IEEE Transactions on Pattern Analysis and Machine Intelligence 45, no. 4 (2022): 4289-4301.

[Weakness 2] Experiments with higher $\tau$

We conducted new experiments during the rebuttal period, by running all the algorithms discussed in our paper for more epochs than what was initially reported on our original manuscript. As requested by the reviewer, in the tables below we present the comparison between DSpodFL and the baselines for the case of FMNIST dataset with non-IID data distribution, and also CIFAR10 with IID distribution. We ran each experiment 4 times, and the information provided in the tables below are the average of those results.

FMNIST Non-IID: Algorithm | $\tau = 10000$ | $\tau = 20000$ | $\tau = 30000$ | $\tau = 40000$

• | - | - | - | - DGD | 0.31 | 0.39 | 0.42 | 0.44 DFedAvg | 0.40 | 0.47 | 0.50 | 0.57 RG | 0.53 | 0.60 | 0.65 | 0.66 Sporadic SGDs | 0.30 | 0.36 | 0.40 | 0.42 DSpodFL | 0.60 | 0.63 | 0.67 | 0.70

CIFAR IID: Algorithm | $\tau = 10000$ | $\tau = 20000$ | $\tau = 30000$ | $\tau = 40000$

• | - | - | - | - DGD | 0.45 | 0.57 | 0.65 | 0.68 DFedAvg | 0.42 | 0.57 | 0.64 | 0.65 RG | 0.62 | 0.68 | 0.74 | 0.75 Sporadic SGDs | 0.42 | 0.55 | 0.60 | 0.65 DSpodFL | 0.68 | 0.73 | 0.75 | 0.76

We observe that even for higher values of total delay $\tau$, our DSpodFL algorithm is able to outperform the baselines.

[Weakness 3] Wall-clock time

We provide experimental results for accuracy vs. wall-clock time per the reviewer’s request. However, we believe that these wall-clock time results are not informative in our setup, since we are simulating clients as different threads on the same computing cluster, as is typical in federated learning research. Therefore, resource heterogeneity does not naturally exist in this setup, and the stochastic gradient computations of all clients take roughly the same amount of wall-clock time to finish. Moreover, physical communication delays between nodes are not present in this single-cluster setup, and they do not significantly contribute to the overall wall-clock time. Hence, all of the algorithms we consider in our paper perform more or less the same in terms of accuracy vs. wall-clock time:

FMNIST IID: Algorithm | Wall-clock time = 0.5 sec | Wall-clock time = 1.0 sec | Wall-clock time = 1.5 sec

• | - | - | - DGD | 0.71 | 0.73 | 0.76 DFedAvg | 0.57 | 0.63 | 0.64 RG | 0.67 | 0.72 | 0.74 Sporadic SGDs | 0.68 | 0.73 | 0.75 DSpodFL | 0.71 | 0.75 | 0.76

FMNIST Non-IID: Algorithm | Wall-clock time = 3 sec | Wall-clock time = 6 sec | Wall-clock time = 9 sec

• | - | - | - DGD | 0.66 | 0.69 | 0.70 DFedAvg | 0.37 | 0.41 | 0.44 RG | 0.52 | 0.55 | 0.58 Sporadic SGDs | 0.64 | 0.69 | 0.71 DSpodFL | 0.60 | 0.63 | 0.66

CIFAR10 IID: Algorithm | Wall-clock time = 10 sec | Wall-clock time = 20 sec | Wall-clock time = 30 sec

• | - | - | - DGD | 0.69 | 0.73 | 0.74 DFedAvg | 0.66 | 0.74 | 0.75 RG | 0.66 | 0.71 | 0.74 Sporadic SGDs | 0.63 | 0.69 | 0.72 DSpodFL | 0.59 | 0.66 | 0.72

CIFAR10 Non-IID: Algorithm | Wall-clock time = 30 sec | Wall-clock time = 60 sec | Wall-clock time = 90 sec

• | - | - | - DGD | 0.67 | 0.77 | 0.80 DFedAvg | 0.74 | 0.78 | 0.80 RG | 0.67 | 0.76 | 0.79 Sporadic SGDs | 0.68 | 0.77 | 0.80 DSpodFL | 0.62 | 0.73 | 0.79

In order to have more informative wall-clock time results using a single computing cluster, we can use the probability of SGDs, i.e., $d_i$’s, to model resource heterogeneity in the computational capabilities of clients. Specifically, we repeat the experiment above, this time calculating a weighted average of the wall-clock time of clients via $1 / d_i$. The following table contains some sample points from the accuracy vs. wall-clock time plots. We observe that our DSpodFL outperforms other baselines when resource heterogeneity is taken into account, i.e., when it takes heterogeneous clients different wall-clock times to finish their computations based on their resource availability. As a final note, we observe that the Sporadic SGDs algorithm is performing similarly to our DSpodFL when we measure wall-clock time scaled inversely by probability of SGDs, i.e., $1 / d_i$. This is because there are no wall-clock time delays due to communications on the single cluster setup, which is necessary in order to differentiate between DSpodFL and Sporadic SGDs. Note that the reported wall-clock times are in seconds, but scaled with $1 / d_i$ for each client.

FMNIST IID: Algorithm | Scaled wall-clock time = 15 sec | Scaled wall-clock time = 30 sec | Scaled wall-clock time = 45 sec

• | - | - | - DGD | 0.66 | 0.68 | 0.71 DFedAvg | 0.50 | 0.53 | 0.54 RG | 0.62 | 0.67 | 0.69 Sporadic SGDs | 0.71 | 0.75 | 0.77 DSpodFL | 0.70 | 0.76 | 0.75

FMNIST Non-IID: Algorithm | Scaled wall-clock time = 3 sec | Scaled wall-clock time = 6 sec | Scaled wall-clock time = 9 sec

• | - | - | - DGD | 0.40 | 0.50 | 0.54 DFedAvg | 0.21 | 0.25 | 0.32 RG | 0.35 | 0.42 | 0.47 Sporadic SGDs | 0.62 | 0.65 | 0.66 DSpodFL | 0.54 | 0.62 | 0.62

CIFAR10 IID: Algorithm | Scaled wall-clock time = 75 sec | Scaled wall-clock time = 150 sec | Scaled wall-clock time = 225 sec

• | - | - | - DGD | 0.51 | 0.61 | 0.69 DFedAvg | 0.59 | 0.63 | 0.73 RG | 0.52 | 0.64 | 0.69 Sporadic SGDs | 0.66 | 0.73 | 0.74 DSpodFL | 0.65 | 0.74 | 0.74

CIFAR10 Non-IID: Algorithm | Scaled wall-clock time = 150 sec | Scaled wall-clock time = 300 sec | Scaled wall-clock time = 450 sec

• | - | - | - DGD | 0.56 | 0.67 | 0.65 DFedAvg | 0.34 | 0.45 | 0.57 RG | 0.54 | 0.64 | 0.69 Sporadic SGDs | 0.70 | 0.76 | 0.79 DSpodFL | 0.74 | 0.79 | 0.80

Again, thank you for your time and efforts in reviewing our paper. We would appreciate further opportunities to answer any remaining concerns you might have.

Thank you for the additional analysis and experiments. The reviewer decides to update the score.

Finally, regarding papers that allow multiple local steps in between aggregations in decentralized FL, like [4, 5, 10], we note that we will in fact have $d_i^{(k)} = 1$ for all $i \in \mathcal{M}$ in those iterations, resulting in $d_{min} = 1$, but instead have $b_{ij}^{(k)} = 0$. Looking at the bounds in Theorems 4.11 and 4.12, we see that substituting $d_{min} = 1$ simplifies the bounds and gives us the ones provided in [4, 5, 10], up to differences in constant scalars.

[9] Mishchenko, Konstantin, Grigory Malinovsky, Sebastian Stich, and Peter Richtárik. "Proxskip: Yes! local gradient steps provably lead to communication acceleration! finally!." In International Conference on Machine Learning, pp. 15750-15769. PMLR, 2022.

[10] Maranjyan, Artavazd, Mher Safaryan, and Peter Richtárik. "Gradskip: Communication-accelerated local gradient methods with better computational complexity." arXiv preprint arXiv:2210.16402 (2022).

[11] Nedic, Angelia, and Asuman Ozdaglar. "Distributed subgradient methods for multi-agent optimization." IEEE Transactions on Automatic Control 54, no. 1 (2009): 48-61.

[Question 6] The constant scalars $w_1, …, w_5$

The values of these scalars which show up in Theorem 4.12 are as follows $w_1 = \frac12 ( 1 - \frac1{\Gamma_0} ) ( 1 - \frac1{\Gamma_2} ) ( 1 - \frac1{\Gamma_4^2} ), w_2 = \frac{\beta^2 d_{\max} ( 1 + \Gamma_3 )}{m ( 1 - \tilde{\rho} ) \left( 1 - \frac1{\Gamma_1} \right)}, w_3 = m d_{\max} ( 16 \frac{1 + \tilde{\rho}}{1 - \tilde{\rho}} \delta^2 + \sigma^2), w_4 = ( 1 + \Gamma_3 ) \delta^2, w_5 = \frac{\beta d_{\max} \sigma^2}{2m}.$ Note that the constant scalars $\Gamma_i$ with $0 \le i \le 4$ are arbitrary which only need to satisfy the conditions $\Gamma_0, \Gamma_1, \Gamma_2, \Gamma_4 > 1$ and $\Gamma_3 > 0$.

To simplify the convergence upper bounds in Theorem 4.12 further using O-notation, we can assume a specific choice of constant learning rate $\alpha$ and minimum SGD probability $d_{min}$. If we run DSpodFL for $K$ iterations using non-convex models, let us choose $\alpha \sim 1 / \sqrt{K}$ and $d_{min} = 1 - 1/K$. Then, we would have that the convergence upper bound diminishes with the rate $\mathcal{O}{\left( \frac1K \left[ {\\| \mathbf{\Theta}^{(0)} - \mathbf{1}\_m \bar{\mathbf{\theta}}^{(0)} \\|}^2 w_2 + w_2 w_3 + w_4 \right] \right)} + \mathcal{O}{\left( \frac1{\sqrt{K}} \left[ (F{( \bar{\theta}^{(0)} )} - F^\star) + w_5 \right] \right)},$ where we have absorbed the constant scalar $w_1$ in the O-notation as it does not depend on problem-related parameters. Substituting the constant scalars $w_2, …, w_5$, we get $\mathcal{O}{\left( \frac1K \left[ \frac{\beta^2 d_{\max} ( 1 + \Gamma_3 )}{m ( 1 - \tilde{\rho} ) \left( 1 - \frac1{\Gamma_1} \right)} \left( {\\| \mathbf{\Theta}^{(0)} - \mathbf{1}\_m \bar{\mathbf{\theta}}^{(0)} \\|}^2 + \frac{\beta^2 d_{\max} ( 1 + \Gamma_3 )}{m ( 1 - \tilde{\rho} ) \left( 1 - \frac1{\Gamma_1} \right)} m d_{\max} ( 16 \frac{1 + \tilde{\rho}}{1 - \tilde{\rho}} \delta^2 + \sigma^2) \right) + ( 1 + \Gamma_3 ) \delta^2 \right] \right)} + \mathcal{O}{\left( \frac1{\sqrt{K}} \left[ (F{( \bar{\theta}^{(0)} )} - F^\star) + \frac{\beta d_{\max} \sigma^2}{2m} \right] \right)}$ for the non-convex convergence rate

[Question 7] Comparison metric $\tau$

As outlined in Sec. 5.1, the average total delay $\tau_{total}^{(k)}$ at iteration $k$ is defined as the sum of average processing delay $\tau_{proc}^{(k)}$ and average transmission delay $\tau_{trans}^{(k)}$ up to iteration $k$. At each iteration, depending on whether a client $i$ computed SGDs or not, i.e., $v_i^{(k)} \in \\{ 0,1 \\}$, there will be a processing delay incurred to finish the computation proportional to $1/d_i^{(k)}$. Similarly, depending on $\hat{v}\_{ij}^{(k)} \in \\{0,1\\}$, some of the links $(i,j)$ in the network graph will be utilized for communications, incurring a transmission delay proportional to $1/b_{ij}^{(k)}$ for that link. In order to obtain comparable result regardless of the size, connectivity and the topology of the underlying graph, we normalize (i.e., average) both processing and transmission delay to obtain $\tau_{proc}^{(k)}$, $\tau_{trans}^{(k)}$. Formally, these are defined as $\tau_{trans}^{(k)} = [ \sum_{i=1}^m{( 1 / | \mathcal{N}\_i | ) \sum_j{\hat{v}\_{ij}^{(k)} / b_{ij}}} ] / [ \sum_{i=1}^m{( 1 / | \mathcal{N}\_i | ) \sum_j{1 / b_{ij}}} ]$ and $\tau_{proc}^{(k)} = [ \sum_{i=1}^m{v_i^{(k)} / d_i}] / [\sum_{i=1}^m{1 / d_i}]$. Finally, we define $\tau_{total}^{(k)} = \tau_{trans}^{(k)} + \tau_{proc}^{(k)}$ as the average total delay as a metric to compare DSpodFL with the baselines.

[Question 8] Effect of $d_i$/$b_{ij}$ on speed of a client/link

Having explained the average total delay $\tau_{total}^{(k)} = \tau_{trans}^{(k)} + \tau_{proc}^{(k)}$ in the answer to [Question 7], we can now see what happens to different baselines under our proposed heterogeneity framework, and how $d_i$ and $b_{ij}$ translate to the speed of a client and a link, respectively.

In DGD, we have that computations and communications occur at every iteration, i.e., $v_i^{(k)} = \hat{v}\_{ij}^{(k)} = 1$ for all $i \in \mathcal{M}$, $(i,j) \in \mathcal{E}$ and $k \geq 0$. Thus, we will have that $\tau_{trans}^{(k)} = [ \sum_{i=1}^m{( 1 / | \mathcal{N}\_i | ) \sum_j{1 / b_{ij}}} ] / [ \sum_{i=1}^m{( 1 / | \mathcal{N}\_i | ) \sum_j{1 / b_{ij}}} ] = 1$ and $\tau_{proc}^{(k)} = [ \sum_{i=1}^m{1 / d_i}] / [\sum_{i=1}^m{1 / d_i}] = 1$, and thus each iteration of training will incur a total average delay of $\tau_{total}^{(k)} = 2$ for this baseline. For DFedAvg, we still have $\tau_{proc}^{(k)} = 1$ as SGDs occur at every iteration. But since aggregations occur every $D$ iterations, i.e., clients conduct $D$ consecutive iterations of local SGD steps before communications, we will have $\tau_{trans}^{(k)} = 0$ for $D$ iterations and then $\tau_{trans}^{(k)} = 1$ for the next iteration. This deterministic cycle will then continue for DFedAvg.

In RG and the Sporadic SGDs baselines, we only have one of these operations occurring at every iteration, i.e., computations and communications, respectively. This means that in RG, we have $v_i^{(k)} = 1$ is deterministic but $\hat{v}\_{ij}^{(k)}$ is stochastic, and for Sporadic SGDs, we have deterministic $\hat{v}\_{ij}^{(k)} =1$ but stochastic $v_i^{(k)}$. Thus, for RG we will have $\tau_{total}^{(k)} = \tau_{trans}^{(k)} + 1 < 2$ implied by $\tau_{proc}^{(k)} = 1$, and for Sporadic SGDs, we will have that $\tau_{total}^{(k)} = 1 + \tau_{proc}^{(k)} < 2$ implied by $\tau_{trans}^{(k)} = 1$.

However, note that in DSpodFL, both computation and communication operations are carried out in a stochastic way, and thus each iteration of training requires less delay incurred on the whole decentralized system to start the next round of training. Note however, that less computation and communication might come at the cost of losing performance at each iteration, but our motivation in DSpodFL was to prove that in fact if we evaluate these algorithms based on their total delay, it can outperform existing baselines. The intuition is that due to the data distribution of clients available in the network, their processing capabilities, the graph topology and the link bandwidth capabilities, it is not necessary to over utilize all of these resources at every iteration for fast convergence. In fact, a resource-aware approach like DSpodFL, takes a step at optimally utilizing the resources to achieve the same final solutions in a shorter amount of time.

We would like to thank you again for your time and efforts in providing a comprehensive review of our paper. We have revised our manuscript accordingly, and would appreciate the opportunity to respond to any remaining concerns that you might have.

We appreciate your time in reading our paper in detail and providing us with your extensive comments. Below please find our responses to each of them, and how we have addressed them in the revised manuscript.

[Weakness 1] Correlation between computations and communications

Thank you for your detailed evaluation of our methodology. In this response, we demonstrate how in fact we allow for communications and computations to be correlated, which in turn is a solution to the reviewer’s concern. In our DSpodFL methodology, we indeed allow for slow gradient computations for clients in between aggregations, so that the latest model a client is using to calculate its gradients is not interrupted.

In our paper, we actually do not need to make the assumption that the computation and communication decisions of each client are uncorrelated in Assumption 4.3. Formally, we do not need $\mathbb{E}{[v_i^{(k)} \hat{v}\_{ij}^{(k)}]} = \mathbb{E}{[v_i^{(k)}]} \mathbb{E}{[\hat{v}\_{ij}^{(k)}]}$ for each $i \in \mathcal{M}$ and $j \in \mathcal{N}_i$. Our proofs only assume that the indicator variables $v_i^{(k)}$ and $\hat{v}\_{ij}^{(k)}$ are uncorrelated across the clients and across the network links, respectively. Formally, $\mathbb{E}{[v_i^{(k)} v_j^{(k)}]} = \mathbb{E}{[v_i^{(k)}]} \mathbb{E}{[v_j^{(k)}]}$ and $\mathbb{E}{[\hat{v}\_{ij}^{(k)} \hat{v}\_{lq}^{(k)}]} = \mathbb{E}{[\hat{v}\_{ij}^{(k)}]} \mathbb{E}{[\hat{v}\_{lq}^{(k)}]}$, for all $i \neq j$ and $(i,j) \neq (l,q)$, respectively. These assumptions essentially mean that the decision to conduct SGDs for a client $i$ at any iteration $k$ ($v_i^{(k)} \in \\{ 0,1 \\}$), is uncorrelated to the decision of a different client $j$ in doing so. Similarly, the decision for clients $(i,j)$ to utilize a link between them to communicate with each other ($\hat{v}\_{ij}^{(k)} \in \\{ 0,1 \\}$) is not correlated with the decision of a different pair of clients $(l,q)$ in doing so. Moreover, within each client, we in fact only need to assume that SGD noises $\epsilon_i^{(k)}$ are uncorrelated with the decision to conduct SGD, i.e., $v_i^{(k)}$, and no further assumptions are needed to be made on the communication decisions $\hat{v}\_{ij}^{(k)}$.

We apologize for this confusion, which was introduced by our attempt to present Assumption 4.3 as succinctly as possible. We have now rewritten the assumption to make this clear. Additionally, we have provided the precise formulation of our assumption in Appendix C.

[Weakness 2] Non-asymptotic convergence rates

We agree with the reviewer that obtaining exact convergence rates in contrast to just asymptotic convergence behavior carries a higher analytical value. In fact, a big part of our convergence analysis has gone into trying to obtain such exact convergence rates for any given iteration $K$. Eq. (9) in Theorem 4.11 shows that the convergence rate of DSpodFL under strongly-convex models is geometric (sometimes also called linear) when a constant learning rate $\alpha$ is used, i.e., $\mathcal{O}{({\rho(\Phi)}^K)}$ where $\rho(\Phi) < 1$. We have also obtained the exact convergence rate for non-convex models, as outlined in Theorem 4.12, which is sub-linear $\mathcal{O}{(1 / K)}$ when a constant learning rate $\alpha$ is used. We would like to emphasize that we only let $K \to \infty$ after obtaining the exact convergence rates, in order to also investigate the asymptotic behavior of DSpodFL.

Moreover, Theorems J.5 and L.5 in the Appendix show that when a diminishing learning rate policy of $\alpha^{(k)} = \alpha^{(0)} / \sqrt{1 + k/\gamma}$ is used, exact convergence can be achieved with a sub-linear rate of $\mathcal{O}{(\ln{K} / \sqrt{K})}$, which is standard in literature [1-3].

[1] Nedić, Angelia, and Alex Olshevsky. "Stochastic gradient-push for strongly convex functions on time-varying directed graphs." IEEE Transactions on Automatic Control 61, no. 12 (2016): 3936-3947.

[2] Makhdoumi, Ali, and Asuman Ozdaglar. "Graph balancing for distributed subgradient methods over directed graphs." In 2015 54th IEEE Conference on Decision and Control (CDC), pp. 1364-1371. IEEE, 2015.

[3] Xi, Chenguang, and Usman A. Khan. "Distributed subgradient projection algorithm over directed graphs." IEEE Transactions on Automatic Control 62, no. 8 (2016): 3986-3992.

[Weakness 3] Equivalence of Eq. (2) and (3)

We apologize for omitting the details due to page constraints. We would like to clarify that Eq. (3) is equivalent to Eq. (2). Please find below why this is the case using the re-definition of the transition weights given in Eq. (4). We have $\theta_i^{(k+1)} = \theta_i^{(k)} + \sum_{j \in \mathcal{M}}{r_{ij} \left( \theta_j^{(k)} - \theta_i^{(k)} \right) \hat{v}\_{ij}^{(k)}} - \alpha^{(k)} g_i^{(k)} v_i^{(k)} = \theta_i^{(k)} (1 - \sum_{j \in \mathcal{M}}{r_{ij} \hat{v}\_{ij}^{(k)}}) + \sum_{j \in \mathcal{M}}{r_{ij} \hat{v}\_{ij}^{(k)} \theta_j^{(k)}} - \alpha^{(k)} g_i^{(k)} v_i^{(k)}.$ Now, using the definition of $p_{ij}^{(k)}$ and $p_{ii}^{(k)}$ as given in Eq. (4) and combining the vectors $\theta_i^{(k)}$ and $g_i^{(k)}$ into rows of the matrices $\Theta^{(k)}$ and $G^{(k)}$, respectively, and forming the doubly-stochastic symmetric matrix $P^{(k)}$ by using the elements $p_{ij}^{(k)}$ and also forming the diagonal matrix $V^{(k)}$ by elements $v_i^{(k)}$, we will arrive at Eq. (3).

[Weakness 4] Convergence rate for non-convex models

We apologize for not being able to present the constant scalars $w_1, …, w_5$ in the main text due to space constraints, and we can see how excluding them from Theorem 4.12 has made it harder to fully grasp. In the revision of our paper, we have directly written them down in Theorem 4.12.

In the original presentation of Theorem 4.12, our main goal was to directly show the influence of the learning rate $\alpha$ and minimum SGD probability $d_{min}$, to emphasize how our bound generalizes the results of existing literature under the addition of sporadic computations. As we have also discussed in the paper, the bound in Theorem 4.12 includes all terms present in the bounds of works such as [4, 5], i.e., $\mathcal{O}{\left( \frac{F{( \bar{\theta}^{(0)} )} - F^\star}{\alpha (K+1)} \right)}$, $\mathcal{O}{\left( \frac{{\\| \mathbf{\Theta}^{(0)} - \mathbf{1}\_m \bar{\mathbf{\theta}}^{(0)} \\|}^2}{K+1} \right)}$, $\mathcal{O}{( \alpha^2 )}$ and $\mathcal{O}{(\alpha)}$, with the addition of an extra term $\mathcal{O}{(1 - d_{min})}$ due to sporadic computations. We observe how this extra term in fact becomes zero if we remove the sporadic computations component of DSpodFL, i.e., choose $d_i^{(k)} = 1$ for all $i \in \mathcal{M}$ so that $d_{\min} = 1$.

Therefore, we have delegated all other parameters, including $\delta$ for data heterogeneity, $\beta$ for smoothness, $m$ for number of clients, $\sigma^2$ for the SGD noise variance and $\tilde{\rho}$ for the spectral radius of the mixing matrix to the constant scalars $w_1, …, w_5$. In our response to [Question 6], please see how each of these constant scalars are dependent on our problem-related parameters.

[4] Koloskova, Anastasia, Nicolas Loizou, Sadra Boreiri, Martin Jaggi, and Sebastian Stich. "A unified theory of decentralized sgd with changing topology and local updates." In International Conference on Machine Learning, pp. 5381-5393. PMLR, 2020.

[5] Sun, Tao, Dongsheng Li, and Bao Wang. "Decentralized federated averaging." IEEE Transactions on Pattern Analysis and Machine Intelligence 45, no. 4 (2022): 4289-4301.

[Weakness 5] Common learning rate for all clients

The choice for the learning rate $\alpha = 0.01$ that we have made for the baselines DGD [4], RG [4] and DFedAvg [5], comes from the original papers. For this reason, we believe the comparison we have made between our DSpodFL algorithm and the baselines, is a fair comparison since we have utilized the same learning rate reported in [4, 5] to properly work in training their algorithms.

On another note, our choice of using the same learning for all baselines, and not individually tuning it for each of them, follows existing research in the decentralized FL literature [4, 6, 7] (where [4] is referenced in our response to [Weakness 4]).

[6] Pu, Shi, and Angelia Nedić. "Distributed stochastic gradient tracking methods." Mathematical Programming 187, no. 1 (2021): 409-457.

[7] Bornstein, Marco, Tahseen Rabbani, Evan Z. Wang, Amrit Bedi, and Furong Huang. "SWIFT: Rapid Decentralized Federated Learning via Wait-Free Model Communication." In The Eleventh International Conference on Learning Representations.

[Weakness 6] Numerical comparison with Even et al., 2024 [8]

The reason that we have not compared our DSpodFL framework to [8] directly, is that [8] does not provide any numerical experiments using the methodology developed within. Without any specific parameter values presented in the paper, we found it hard to make a fair comparison between our method and [8]. In other words, [8] is mostly a theoretical contribution to the field of decentralized optimization, which we have referenced in our work.

[8] Even, Mathieu, Anastasia Koloskova, and Laurent Massoulié. "Asynchronous sgd on graphs: a unified framework for asynchronous decentralized and federated optimization." In International Conference on Artificial Intelligence and Statistics, pp. 64-72. PMLR, 2024.

[Question 1] Time-invariant resources in Even et al., 2024 [8]

We have stated that [8] (as referenced in our response to [Weakness 6]) does not allow time variations in resource heterogeneity, as they have assumed that resource availability of clients stay the same across the training process. Indeed, they allow for resource heterogeneity in their formulation as the reviewer has correctly mentioned, it is just that any initial resource considered for a client $i \in \mathcal{M}$ will remain fixed throughout the training. In contrast, through our formulation of DSpodFL using simple yet elegant indicator variables $v_i^{(k)}$ and $\hat{v}\_{ij}^{(k)}$ and their expected values $\mathbb{E}{[v_i^{(k)}]} = d_i^{(k)}$ and $\mathbb{E}{[\hat{v}\_{ij}^{(k)}]} = b_{ij}^{(k)}$, we allow for time-varying $d_i^{(k)}$ and $b_{ij}^{(k)}$ in DSpodFL. This in turn corresponds to dynamic resource heterogeneity in the system, which is not considered in [8]. Allowing for $d_i^{(k)}$ and $b_{ij}^{(k)}$ to be time-varying makes the theoretical analysis more challenging, e.g., see Lemma D.4.

Also note that in Appendix P.4, we have conducted several experiments under the setup where both computational and communication resources are time-varying over the training process. We observe that DSpodFL is able to outperform the baselines even in these scenarios (some more than the other).

[Question 2] Clarification on $d_{max}^{(k)} = 0$ in Proposition 4.10

We can definitely substitute $d_{max}^{(k)} = 0$ in the constraint on $\alpha^{(k)}$ given in Proposition 4.10, and since this constraint is the minimum of three values, we will get that $\alpha^{(k)} < \min\\{ 1/\mu, \infty, \infty \\} = 1/\mu$. However, we also note that $d_{max}^{(k)} = 0$ essentially means that we have $d_i^{(k)} = 0$ for all $i \in \mathcal{M}$ for a particular iteration $k$. Therefore, none of the clients will compute stochastic gradients at that iteration, and thus there is no need to choose any learning rate $\alpha^{(k)}$ to begin with.

[Question 3] Clarification on $\lim_{k \to 0}{\Phi^{(k:0)}} = 0$ in Proposition 4.10

We thank you again for your thorough examination of our paper. It is correct that the values that the consensus error and average model error converge to, depend on $\Psi^{(k)}$. However, as we have stated in the paper, $\lim_{k \to 0}{\Phi^{(k:0)}} = 0$ is a sufficient condition just for convergence, and the exact values that $\nu^{(k)}$ will converge to will depend on $\Psi^{(k)}$. Since the elements of $\Psi^{(k)}$, i.e., $\psi_1^{(k)}$ and $\psi_2^{(k)}$, are bounded values as seen in their definitions in Lemmas 4.7 and 4.8, we do not need to worry about divergence of $\nu^{(k)}$ as long as $\lim_{k \to 0}{\Phi^{(k:0)}} = 0$. However, convergence of $\nu^{(k)}$ to non-zero values is indeed possible, and occurs in Theorems 4.11 and 4.12 due to using a constant learning rate. Achieving zero optimality gap requires DSpodFL to use diminishing learning rate, which is standard in literature, and is presented in detail in Appendices J and L.

[Question 4] Double stochasticity of matrices $R$ and $P^{(k)}$

The adjustment of transition weights $p_{ii}^{(k)}$ to make the matrix $P^{(k)}$ doubly stochastic, is primarily to account for the sporadicity of communications at each iteration $k$ of DSpodFL. In other words, since we have $p_{ij}^{(k)} = r_{ij} \hat{v}\_{ij}^{(k)}$ and $\hat{v}\_{ij}^{(k)} \in \\{ 0,1 \\}$, we essentially zero out some of the elements at each row of the matrix $R$ to form the off-diagonal elements of $P^{(k)}$, and thus need to adjust the diagonal elements to retain double stochasticity. However, note that the guarantee that $0 < p_{ii}^{(k)} \le 1$ relies on double stochasticity of $R$ and Metropolis-Hastings weights that are used to design its entries $r_{ij}$. Therefore, although the reviewer is correct in that we do not use the diagonal values $r_{ii}$ of $R$ directly and therefore they can have any arbitrary value, but we still need off-diagonal elements of each row in matrix $R$ to satisfy $0 \le \sum_{j \in \mathcal{M}, j \neq i}{r_{ij}} < 1$.

For coherence in moving from $R$ to $P^{(k)}$, we have opted to complete the constraint $0 \le \sum_{j \in \mathcal{M}, j \neq i}{r_{ij}} < 1$ by defining $r_{ii}$ such that the matrix $R$ would become doubly-stochastic. This is also because if for a client $i \in \mathcal{M}$, we have that it is connected to all other clients in the network, i.e., $\mathcal{N}\_i = \mathcal{M} \setminus \\{ i \\}$, then if at some iteration it communicates with all of those clients, i.e., $\hat{v}\_{ij}^{(k)} = 1$ for all $j \in \mathcal{N}\_i$, then $p_{ii}^{(k)} = r_{ii}$.

[Question 5] Theoretical comparison with Koloskova et al., 2020 [4]

First, we note that as outlined in Table 1 of Sec. 1, we provide convergence for last iterates of model parameters when dealing with strongly-convex models in our paper, in contrast to the majority of existing literature which only provide convergence for average iterates. Therefore when dealing with strongly-convex models, we can compare our theoretical results only for the DGD baseline because the DGD-like algorithms given in [9, 10] are among the few to show convergence for the last iterates of model parameters as well. We apologize for any confusion that we might have caused by citing the work by Koloskova et al., 2020 [4] (referenced in our response to [Weakness 4]) in lines 386 and 387 of our original manuscript, as our goal was to cite this work as a paper which presents a unified framework for DGD methods. However, we will be able to compare our results in the non-convex case with [4], which we will talk about later in this response. In the revised manuscript when presenting our results for convex models, we have changed the paper we cite on DGD-like methods to [10], which also shows convergence for last iterates and we are able to better compare our analytical bounds with it when dealing with strongly-convex models.

• For strongly-convex models, according to Theorems 3.6, 5.5, 5.7 and D.1 in [9] and Theorems 3.5, 4.5 and B.2 in [10], the convergence rate of the algorithms ProxSkip, Decentralized Scaffnew, SplitSkip [9], GradSkip, GradSkip+ and VR-GradSkip+ [10] are all geometric, i.e., $\mathcal{O}(\rho^K)$ with $0 < \rho < 1$ (note that most of these algorithms are FL algorithms, and not decentralized FL algorithms). This rate agrees with the rate we provide in Eq. (10) of Theorem 4.11.

• When dealing with non-convex models, we can compare our rate with all the baselines, i.e., [11] for DGD, [4] for RG and [5] for DFedAvg (referenced in our response to [Weakness 4]). Since DGD is a special case of RG with no sporadicity in aggregations, i.e., $b_{ij} = 1$ for all $(i,j) \in \mathcal{E}$, we will compare our work with the more recent paper [4]. For DGD and RG, Lemma 17 of [4] shows the convergence upper bound before tuning the constant learning rate $\alpha$, which is $\mathcal{O}{\left( \frac{\mathbb{E}{[F{( \bar{\theta}^{(0)} )}] - F^\star}}{\alpha (K+1)} \right)} + \mathcal{O}{(\alpha)} + \mathcal{O}{(\alpha^2)}.$ For DFedAvg, Theorem 1 in [5] with a zero momentum ($\theta = 0$) obtains a less tight bound as follows $\mathcal{O}{\left( \frac{\mathbb{E}{[F{( \bar{\theta}^{(0)} )}] - F^\star}}{\alpha (K+1)} \right)} + \mathcal{O}{(\alpha)} + \mathcal{O}{(\alpha^2)} + \mathcal{O}{(\alpha^3)}.$ We observe that by setting $d_{min} = 1$ in Theorem 4.12 of our paper, these convergence rates are recovered. Recall that as discussed in Fig. 1 of our paper, all of these baseline algorithms fit within the general framework of DSpodFL with $d_i = 1$ for all clients $i \in \mathcal{M}$. That is why we can substitute $d_{min} = 1$ in Theorem 4.12 to compare our analytical results with the ones provided in [3-5].

We are grateful for your follow up questions regarding our paper. Below, please find our responses to each comment.

Follow-up to [Weakness 1]:

It is correct that our formulation allows for $v_i^{(k)}$ and $\hat{v}_{ij}^{(k)}$ to be uncorrelated only at iteration $k$, and our expectation calculations at any iteration $k$ in fact depend on conditioning over past realizations of random variables involved in our paper. However, we will demonstrate that in our formulation of DSpodFL, it is sufficient for those indicator random variables to be uncorrelated for only one iteration. This is mainly because of the way we define an iteration of training in our work, which is elaborated in several steps as follows:

(i) Start iteration $k=0$ at time $t=0$.

(ii) Let $T_i(t)$ be the processing time it takes for client $i$ to finish its gradient computation at time $t$. Then, start all inter-client communications at time $t’ = t + \max_{i \in \mathcal{M}}{T_i(t)}$ after all gradient computations have finished. This essentially means having correlated $v_i$ and $\hat{v}_{ij}$ for a single iteration, because we are delaying the communications until all gradients have been computed.

(iii) Define $T_{(i,j)}'(t')$ be the total transmission time it takes for clients $i, j$ to use the link $(i, j)$ to transmit their models to each other at time $t'$. Then wait until $t'' = t' + \max_{(i,j) \in \mathcal{E}}{T_{(i,j)}'(t')}$ before starting any other gradient computation for the clients. This again is possible due to $v_i$ and $\hat{v}_{ij}$ being allowed to be correlated for a single iteration.