Asynchronous adaptive networks

In a recent article [1] we surveyed advances related to adaptation, learning, and optimization over synchronous networks. Various distributed strategies were discussed that enable a collection of networked agents to interact locally in response to streaming data and to continually learn and adapt to track drifts in the data and models. Under reasonable technical conditions on the data, the adaptive networks were shown to be mean-square stable in the slow adaptation regime, and their mean-square-error performance and convergence rate were characterized in terms of the network topology and data statistical moments [2]. Classical results for single-agent adaptation and learning were recovered as special cases. Following the works [3]-[5], this chapter complements the exposition from [1] and extends the results to asynchronous networks. The operation of this class of networks can be subject to various sources of uncertainties that influence their dynamic behavior, including randomly changing topologies, random link failures, random data arrival times, and agents turning on and off randomly. In an asynchronous environment, agents may stop updating their solutions or may stop sending or receiving information in a random manner and without coordination with other agents. The presentation will reveal that the mean-square-error performance of asynchronous networks remains largely unaltered compared to synchronous networks. The results justify the remarkable resilience of cooperative networks in the face of random events.


I. INTRODUCTION
Adaptive networks consist of a collection of agents with learning abilities. The agents interact with each other on a local level and diffuse information across the network to solve inference and optimization tasks in a decentralized manner. Such networks are scalable, robust to node and link failures, and are particularly suitable for learning from big data sets by tapping into the power of collaboration among distributed agents. The networks are also endowed with cognitive abilities due to the sensing abilities of their agents, their interactions with their neighbors, and the embedded feedback mechanisms for acquiring and refining information. Each agent is not only capable of sensing data and experiencing the environment directly, but it also receives information through interactions with its neighbors and processes and analyzes this information to drive its learning process.
As already indicated in [1], [5], there are many good reasons for the peaked interest in networked solutions, especially in this day and age when the word "network" has become commonplace whether one is referring to social networks, power networks, transportation networks, biological networks, or other networks. Some of these reasons have to do with the benefits of cooperation over networks in terms of improved performance and improved robustness and resilience to failure. Other reasons deal with privacy and secrecy considerations where agents may not be comfortable sharing their data with remote fusion centers. In other situations, the data may already be available in dispersed locations, as happens with cloud computing. One may also be interested in learning and extracting information through data mining from large data sets. Decentralized learning procedures offer an attractive approach to dealing with such data sets. Decentralized mechanisms can also serve as important enablers for the design of robotic swarms, which can assist in the exploration of disaster areas.

A. Asynchronous Behavior
The survey article [1] and monograph [5] focused on the case of synchronous networks where data arrive at all agents in a synchronous manner and updates by the agents are also performed in a synchronous manner. The network topology was assumed to remain largely static during the adaptation process. Under these conditions, the limits of performance and stability of these networks were identified in some detail for two main classes of distributed strategies: consensus and diffusion constructions. In this chapter, we extend the overview from [1] to cover asynchronous environments. In such environments, the operation of the network can suffer from the occurrence of various random events including randomly changing topologies, random link failures, random data arrival times, and agents turning on and off randomly. Agents may also stop updating their solutions or may stop sending or receiving information in a random manner and without coordination with other agents. Results in [2]- [4] examined the implications of such asynchronous events on network performance in some detail and under a fairly To appear as book chapter in Cooperative and Graph Signal Processing, P. Djuric and C. Richard, editors, Elsevier, 2018. The work in this book chapter was supported in part by NSF grants ECCS-1407712 and CCF-1524250. The authors are grateful to IEEE for allowing reproduction of substantial material from [1]- [4] in this book chapter. Co-author X. Zhao was a PhD student in Electrical Engineering at UCLA. The first author is with the Ecole Polytechnique Federale de Lausanne, EPFL, CH-1015 Lausanne. Contact email: ali.sayed@epfl.ch. general model for the random events. The purpose of this chapter is to summarize the key conclusions from these works in a manner that complements the presentation from [1] for the benefit of the reader. While the works [2]- [4] consider a broader formulation involving complex-valued variables, we limit the discussion here to real-valued variables in order not to overload the notation and to convey the key insights more directly. Proofs and derivations are often omitted and can be found in the above references; the emphasis is on presenting the results in a motivated manner and on commenting on the insights they provide into the operation of asynchronous networks.
We indicated in [2]- [4] that there already exist many useful studies in the literature on the performance of consensus strategies in the presence of asynchronous events (see, e.g., [11]- [21]. There are also some studies in the context of diffusion strategies [22], [23]. However, with the exception of the latter two works, the earlier references assumed conditions that are not generally favorable for applications involving continuous adaptation and learning. For example, some of the works assumed a decaying step-size, which turns off adaptation after sufficient iterations have passed. Some other works assumed noise free data, which is a hindrance when learning from data perturbed by interferences and distortions. A third class of works focused on studying pure averaging algorithms, which are not required to respond to continuous data streaming. In the works [2]- [4], we adopted a more general asynchronous model that removes these limitations by allowing for various sources of random events and, moreover, the events are allowed to occur simultaneously. We also examined learning algorithms that respond to streaming data to enable adaptation. The main conclusion from the analysis in these works, and which will be summarized in future sections, is that asynchronous networks can still behave in a mean-square-error stable manner for sufficiently small step-sizes and, interestingly, their steady-state performance level is only slightly affected in comparison to synchronous behavior. The iterates computed by the various agents are still able to converge and hover around an agreement state with a small mean-square-error. These are reassuring results that support the intrinsic robustness and resilience of network-based cooperative solutions.

B. Organization of the Chapter
Readers will benefit more from this chapter if they review first the earlier article [1]. We continue to follow a similar structure here, as well as a similar notation, since the material in both this chapter and the earlier reference [1] are meant to complement each other. We organize the presentation into three main components. The first part (Sec. II) reviews fundamental results on adaptation and learning by single stand-alone agents. The second part (Sec. III) covers asynchronous centralized solutions. The objective is to explain the gain in performance that results from aggregating the data from the agents and processing it centrally at a fusion center. The centralized performance is used as a frame of reference for assessing various implementations. While centralized solutions can be powerful, they nevertheless suffer from a number of limitations. First, in real-time applications where agents collect data continuously, the repeated exchange of information back and forth between the agents and the fusion center can be costly especially when these exchanges occur over wireless links or require nontrivial routing resources. Second, in some sensitive applications, agents may be reluctant to share their data with remote centers for various reasons including privacy and secrecy considerations. More importantly perhaps, centralized solutions have a critical point of failure: if the central processor fails, then this solution method collapses altogether.
For these reasons, we cover in the remaining sections of the chapter (Secs. IV and V) distributed asynchronous strategies of the consensus and diffusion types, and examine their dynamics, stability, and performance metrics. In the distributed mode of operation, agents are connected by a topology and they are permitted to share information only with their immediate neighbors. The study of the behavior of such networked agents is more challenging than in the single-agent and centralized modes of operation due to the coupling among interacting agents and due to the fact that the networks are generally sparsely connected.

II. SINGLE-AGENT ADAPTATION AND LEARNING
We begin our treatment by reviewing stochastic gradient algorithms, with emphasis on their application to the problems of adaptation and learning by stand-alone agents.

A. Risk and Loss Functions
Thus, let J(w) : IR M×1 → IR denote a twice-differentiable real-valued (cost or utility or risk) function of a real-valued vector argument, w ∈ IR M×1 . When the variable w is complex-valued, some important technical differences arise that are beyond the scope of this chapter; they are addressed in [2]- [5]. Likewise, some adjustments to the arguments are needed when the risk function is non-smooth (non-differentiable); as explained in the works [24], [25]. It is sufficient for our purposes in this chapter to limit the presentation to real arguments and smooth risk functions without much loss in generality.
We denote the gradient vectors of J(w) relative to w and w T by the following row and column vectors, respectively: These definitions are in terms of the partial derivatives of J(w) relative to the individual entries of w = col{w 1 , w 2 , . . . , w M }, where the notation col{·} refers to a column vector that is formed by stacking its arguments on top of each other. Likewise, the Hessian matrix of J(w) with respect to w is defined as the following M × M symmetric matrix: which is constructed from two successive gradient operations. It is common in adaptation and learning applications for the risk function J(w) to be constructed as the expectation of some loss function, Q(w; x), where the boldface variable x is used to denote some random data, say, J(w) = E Q(w; x) (2) and the expectation is evaluated over the distribution of x.
Example II.1 (Mean-square-error (MSE) costs). Let d denote a zero-mean scalar random variable with variance σ 2 d = E d 2 , and let u denote a zero-mean 1 × M random vector with covariance matrix R u = E u T u > 0. The combined quantities {d, u} represent the random variable x referred to in (2). The cross-covariance vector is denoted by r du = E du T . We formulate the problem of estimating d from u in the linear least-mean-squares-error sense or, equivalently, the problem of seeking the vector w o that minimizes the quadratic cost function: This cost corresponds to the following choice for the loss function: Such quadratic costs are widely used in estimation and adaptation problems [26]- [30]. They are also widely used as quadratic risk functions in machine learning applications [31], [32]. The gradient vector and Hessian matrix of J(w) are easily seen to be: Example II.2 (Logistic or log-loss risks). Let γ denote a binary random variable that assumes the values ±1, and let h denote an M × 1 random (feature) vector with R h = E hh T . The combined quantities {γ, h} represent the random variable x referred to in (2). In the context of machine learning and pattern classification problems [31]- [33], the variable γ designates the class that feature vector h belongs to. In these problems, one seeks the vector w o that minimizes the regularized logistic risk function: where ρ > 0 is some regularization parameter, ln(·) is the natural logarithm function, and w 2 = w T w. The risk (4a) corresponds to the following choice for the loss function: Once w o is recovered, its value can be used to classify new feature vectors, say, {h ℓ }, into classes +1 or −1. This can be achieved by assigning feature vectors with h T ℓ w o ≥ 0 to one class and feature vectors with h T ℓ w o < 0 to another class. It can be easily verified that: where I M denotes the identity matrix of size M × M .

B. Conditions on Cost Function
Stochastic gradient algorithms are powerful iterative procedures for solving optimization problems of the form While the analysis that follows can be pursued under more relaxed conditions (see, e.g., the treatments in [34]- [37]), it is sufficient for our purposes to require J(w) to be strongly-convex and twice-differentiable with respect to w. The cost function J(w) is said to be ν-strongly convex if, and only if, its Hessian matrix is sufficiently bounded away from zero [35], [38]- [40]: for all w and for some scalar ν > 0, where the notation A > 0 signifies that matrix A is positive-definite. Strong convexity is a useful condition in the context of adaptation and learning from streaming data because it helps guard against ill-conditioning in the algorithms; it also helps ensure that J(w) has a unique global minimum, say, at location w o ; there will be no other minima, maxima, or saddle points. In addition, it is well-known that strong convexity helps endow stochastic-gradient algorithms with geometric convergence rates in the order of O(α i ), for some 0 ≤ α < 1 and where i is the iteration index [35], [36]. In many problems of interest in adaptation and learning, the cost function J(w) is either already strongly-convex or can be made strongly-convex by means of regularization. For example, it is common in machine learning problems [31], [32] and in adaptation and estimation problems [28], [30] to incorporate regularization factors into the cost functions; these factors help ensure strong convexity. For instance, the mean-square-error cost (3a) is strongly convex whenever R u > 0. If R u happens to be singular, then the following regularized cost will be strongly convex: where ρ > 0 is a regularization parameter similar to (4a). Besides strong convexity, we shall also assume that the gradient vector of J(w) is δ-Lipschitz, namely, there exists δ > 0 such that for all w 1 and w 2 . It can be verified that for twice-differentiable costs, conditions (6) and (8) combined are equivalent to For example, it is clear that the Hessian matrices in (3c) and (4d) satisfy this property since in the first case and in the second case, where the notation λ min (R) and λ max (R) refers to the smallest and largest eigenvalues of the symmetric matrix argument, R, respectively. In summary, we will be assuming the following conditions [2], [5], [41], [42].
Assumption II.1 (Conditions on cost function). The cost function J(w) is twice-differentiable and satisfies (9) for some positive parameters ν ≤ δ. Condition (9) is equivalent to requiring J(w) to be ν-strongly convex and for its gradient vector to be δ-Lipschitz as in (6) and (8), respectively.

C. Stochastic-Gradient Approximation
The traditional gradient-descent algorithm for solving (5) takes the form: where i ≥ 0 is an iteration index and µ > 0 is a small step-size parameter. Starting from some initial condition, w −1 , the iterates {w i } correspond to successive estimates for the minimizer w o . In order to run recursion (11), we need to have access to the true gradient vector. This information is generally unavailable in most instances involving learning from data. For example, when cost functions are defined as the expectations of certain loss functions as in (2), the statistical distribution of the data x may not be known beforehand. In that case, the exact form of J(w) will not be known since the expectation of Q(w; x) cannot be computed. In such situations, it is customary to replace the true gradient vector, ∇ w T J(w i−1 ), by an instantaneous approximation for it, and which we shall denote by ∇ w T J(w i−1 ). Doing so leads to the following stochastic-gradient recursion in lieu of (11): We use the boldface notation, w i , for the iterates in (12) to highlight the fact that these iterates are now randomly perturbed versions of the values {w i } generated by the original recursion (11). The random perturbations arise from the use of the approximate gradient vector. The boldface notation is therefore meant to emphasize the random nature of the iterates in (12).
We refer to recursion (12) as a synchronous implementation since updates occur continuously over the iteration index i. This terminology is meant to distinguish the above recursion from its asynchronous counterpart, which is introduced later in Sec. II-E. We illustrate construction (12) by considering a scenario from classical adaptive filter theory [26]- [28], where the gradient vector is approximated directly from data realizations. The construction will reveal why stochastic-gradient implementations of the form (12), using approximate rather than exact gradient information, become naturally endowed with the ability to respond to streaming data.
Example II.3 (LMS adaptation). Let d(i) denote a streaming sequence of zero-mean random variables with variance σ 2 d = E d 2 (i). Let u i denote a streaming sequence of 1 × M independent zero-mean random vectors with covariance matrix R u = E u T i u i > 0. Both processes {d(i), u i } are assumed to be jointly wide-sense stationary. The cross-covariance vector between d(i) and u i is denoted by r du = E d(i)u T i . The data {d(i), u i } are assumed to be related via a linear regression model of the form: for some unknown parameter vector w o , and where v(i) is a zero-mean white-noise process with power σ 2 v = E v 2 (i) and assumed independent of u j for all i, j. Observe that we are using parentheses to represent the time-dependency of a scalar variable, such as writing d(i), and subscripts to represent the time-dependency of a vector variable, such as writing u i . This convention will be used throughout the chapter. In a manner similar to Example II.1, we again pose the problem of estimating w o by minimizing the mean-square error cost where now the quantities {d(i), u i } represent the random data x i in the definition of the loss function, Q(w; x i ). Using (11), the gradient-descent recursion in this case will take the form: The main difficulty in running this recursion is that it requires knowledge of the moments {r du , R u }. This information is rarely available beforehand; the adaptive agent senses instead realizations {d(i), u i } whose statistical distributions have moments {r du , R u }. The agent can use these realizations to approximate the moments and the true gradient vector. There are many constructions that can be used for this purpose, with different constructions leading to different adaptive algorithms [26]- [29].
It is sufficient to focus on one of the most popular adaptive algorithms, which results from using the data {d(i), u i } to compute instantaneous approximations for the unavailable moments as follows: By doing so, the true gradient vector is approximated by: Observe that this construction amounts to replacing the true gradient vector, ∇ w T J(w), by the gradient vector of the loss function itself (which, equivalently, amounts to dropping the expectation operator). Substituting (13e) into (13c) leads to the well-known (synchronous) least-mean-squares (LMS, for short) algorithm [26]- [28], [43]: The LMS algorithm is therefore a stochastic-gradient algorithm. By relying directly on the instantaneous data {d(i), u i }, the algorithm is infused with useful tracking abilities. This is because drifts in the model w o from (13a) will be reflected in the data {d(i), u i }, which are used directly in the update (13f).
If desired, it is also possible to employ iteration-dependent step-size sequences, µ(i), in (12) instead of the constant step-size µ, and to require µ(i) to satisfy Under some technical conditions, it is well-known that such step-size sequences ensure the convergence of w i towards w o almost surely as i → ∞ [5], [35]- [37]. However, conditions (14) force the step-size sequence to decay to zero, which is problematic for applications requiring continuous adaptation and learning from streaming data. This is because, in such applications, it is not unusual for the location of the minimizer, w o , to drift with time. With µ(i) decaying towards zero, the stochastic-gradient algorithm (12) will stop updating and will not be able to track drifts in the solution. For this reason, we shall focus on constant step-sizes from this point onwards since we are interested in solutions with tracking abilities. Now, the use of an approximate gradient vector in (12) introduces perturbations relative to the operation of the original recursion (11). We refer to the perturbation as gradient noise and define it as the difference: The presence of this perturbation prevents the stochastic iterate, w i , from converging almost-surely to the minimizer w o when constant step-sizes are used. Some deterioration in performance will occur and the iterate w i will instead fluctuate close to w o . We will assess the size of these fluctuations by measuring their steady-state mean-square value (also called mean-squaredeviation or MSD). It will turn out that the MSD is small and in the order of O(µ) -see (28c) further ahead. It will also turn out that stochastic-gradient algorithms converge towards their MSD levels at a geometric rate. In this way, we will be able to conclude that adaptation with small constant step-sizes can still lead to reliable performance in the presence of gradient noise, which is a reassuring result. We will also be able to conclude that adaptation with constant step-sizes is useful even for stationary environments. This is because it is generally sufficient in practice to reach an iterate w i within some fidelity level from w o in a finite number of iterations. As long as the MSD level is satisfactory, a stochastic-gradient algorithm will be able to attain satisfactory fidelity within a reasonable time frame. In comparison, although diminishing step-sizes ensure almost-sure convergence of w i to w o , they nevertheless disable tracking and can only guarantee slower than geometric rates of convergence (see, e.g., [5], [5], [35], [36]). The next example from [41] illustrates the nature of the gradient noise process (15) in the context of mean-square-error adaptation.
Example II.4 (Gradient noise). It is clear from the expressions in Example II.3 that the corresponding gradient noise process is where we introduced the error vector: Let the symbol F i−1 represent the collection of all possible random events generated by the past iterates {w j } up to time i − 1 (more formally, F i−1 is the filtration generated by the random process w j for j ≤ i − 1): It follows from the conditions on the random processes for some constant c ≥ 0. If we take expectations of both sides of (16e), we further conclude that the variance of the gradient noise, E s i (w i−1 ) 2 , is bounded by the combination of two factors. The first factor depends on the quality of the iterate, E w i−1 2 , while the second factor depends on σ 2 v . Therefore, even if the adaptive agent is able to approach w o with great fidelity so that E w i−1 2 is small, the size of the gradient noise will still depend on σ 2 v .

D. Conditions on Gradient Noise Process
In order to examine the convergence and performance properties of the stochastic-gradient recursion (12), it is necessary to introduce some assumptions on the stochastic nature of the gradient noise process, s i (·). The conditions that we introduce in the sequel are similar to conditions used earlier in the optimization literature, e.g., in [35, pp. 95-102] and [44, p. 635]; they are also motivated by the conditions we observed in the mean-square-error case in Example II.4. Following the developments in [2], [41], [42], we let denote the conditional second-order moment of the gradient noise process, which generally depends on i. We assume that, in the limit, the covariance matrix tends to a constant value when evaluated at w o and is denoted by For example, comparing with expression (16a) for mean-square-error costs, we have Assumption II.2 (Conditions on gradient noise). It is assumed that the first and second-order conditional moments of the gradient noise process satisfy (17b) and almost surely, for some nonnegative scalars β 2 and σ 2 s .
Condition (19a) ensures that the approximate gradient vector is unbiased. It follows from conditions (19a)-(19b) that the gradient noise process itself satisfies: It is straightforward to verify that the gradient noise process (16a) in the mean-square-error case satisfies conditions (19a)-(19b). Note in particular from (16e) that we can make the identifications σ 2 s → 4σ 2 v Tr(R u ) and β 2 → 4c.

E. Random Updates
We examined the performance of synchronous updates of the form (12) in some detail in [1], [5]. As indicated earlier, the focus of the current chapter is on extending the treatment from [1] to asynchronous implementations. Accordingly, the first main digression in the exposition relative to [1] occurs at this stage.
Thus, note that the stochastic-gradient recursion (12) employs a constant step-size parameter, µ > 0. This means that this implementation expects the approximate gradient vector, ∇ w T J(w i−1 ), to be available at every iteration. Nevertheless, there are situations where data may arrive at the agent at random times, in which case the updates will also be occurring at random times. One way to capture this behavior is to model the step-size parameter as a random process, which we shall denote by the boldface notation µ(i) (since boldface letters in our notation refer to random quantities). Doing so, allows us to replace the synchronous implementation (12) by the following asynchronous recursion: Observe that we are attaching a time index to the step-size parameter in order to highlight that its value will now be changing randomly from one iteration to another. For example, one particular instance discussed further ahead in Example II.5 is when µ(i) is a Bernoulli random variable assuming one of two possible values, say, µ(i) = µ > 0 with probability p µ and µ(i) = 0 with probability 1 − p µ . In this case, recursion (21) will be updating p µ −fraction of the time. We will not limit our presentation to the Bernoulli model but will allow the random step-size to assume broader probability distributions; we denote its first and second-order moments as follows.
Assumption II.3 (Conditions on step-size process). It is assumed that the stochastic process {µ(i), i ≥ 0} consists of a sequence of independent and bounded random variables, µ(i) ∈ [0, µ ub ], where µ ub > 0 is a constant upper bound. The mean and variance of µ(i) are fixed over i, and they are denoted by: withμ > 0 and σ 2 µ ≥ 0. Moreover, it is assumed that, for any i, the random variable µ(i) is independent of any other random variable in the learning algorithm.
The following variable will play an important role in characterizing the mean-square-error performance of asynchronous updates and, hence, we introduce a unique symbol for it: where the subscript "x" is meant to refer to the "asynchronous" mode of operation. This expression captures the first and second-order moments of the random variations in the step-size parameter into a single variable, µ x . While the constant step-size µ determines the performance of the synchronous implementation (12), it turns out that the constant variable µ x defined above will play an equivalent role for the asynchronous implementation (21). Note further that the synchronous stochastic-gradient iteration (12) can be viewed as a special case of recursion (21) when the variance of µ(i) is set to zero, i.e., σ 2 µ = 0, and the mean valueμ is set to µ. Therefore, by using these substitutions, we will able to deduce performance metrics for (12) from the performance metrics that we shall present for (21). The following two examples illustrate situations involving random updates.
Example II.5 (Random updates under a Bernoulli model). Assume that at every iteration i, the agent adopts a random "on-off" policy to reduce energy consumption. It tosses a coin to decide whether to enter an active learning mode or a sleeping mode. Let 0 < p µ < 1 denote the probability of entering the active mode. During the active mode, the agent employs a step-size value µ. This model is useful for situations in which data arrives randomly at the agent: at every iteration i, new data is available with probability p µ . The random step-size process is therefore of the following type: In this case, the mean and variance of µ(i) are given by: Example II.6 (Random updates under a Beta model). Since the random step-size µ(i) is limited to a finite-length interval [0, µ ub ], we may extend the Bernoulli model from the previous example by adopting a more general continuous Beta distribution for µ(i). The Beta distribution is an extension of the Bernoulli distribution. While the Bernoulli distribution assumes two discrete possibilities for the random variable, say, {0, µ}, the Beta distribution allows for any value in the continuum [0, µ]. Thus, let x denote a generic scalar random variable that assumes values in the interval [0, 1] according to a Beta distribution. Then, according to this distribution, the pdf of x, denoted by f x (x; ξ, ζ), is determined by two shape parameters {ξ, ζ} as follows [45], [46]: where Γ(·) denotes the Gamma function [47], [48]. Figure 1 plots f x (x; ξ, ζ) for two values of ζ. The mean and variance of the Beta distribution (25a) are given by: We note that the classical uniform distribution over the interval [0, 1] is a special case of the Beta distribution for ξ = ζ = 1 -see Fig. 1. Likewise, the Bernoulli distribution with p µ = 1/2 is recovered from the Beta distribution by letting ξ = ζ → 0. In the Beta model for asynchronous adaptation, we assume that the ratio µ(i)/µ ub follows a Beta distribution with parameters {ξ, ζ}. Under this model, the mean and variance of the random step-size become:

F. Mean-Square-Error Stability
We now examine the convergence of the asynchronous stochastic-gradient recursion (21). In the statement below, the notation a = O(µ) means a ≤ bµ for some constant b that is independent of µ.   the gradient noise process, and the random step-size process hold. Let µ o = 2ν/(δ 2 + β 2 ). For any µ x satisfying it holds that E w i 2 converges exponentially (i.e., at a geometric rate) according to the recursion where the scalar α satisfies 0 ≤ α < 1 and is given by It follows from (28a) that, for sufficiently small step-sizes: Proof: We subtract w o from both sides of (21) to get We now appeal to the mean-value theorem [5], [35], [49] to write: where we are introducing the symmetric and random time-variant matrix H i−1 to represent the integral expression. Substituting into (29a), we get so that from Assumption II.2: It follows from (9) that since ν ≤ δ. In the first line above, the notation ρ(A) denotes the spectral radius of its matrix argument (i.e., ρ(A) = max k |λ k (A)| in terms of the largest magnitude eigenvalue of A). From (31a) we obtain Taking expectations of both sides of (30), we arrive at (28a) from (20b) and (31b) with α given by (28b). The bound in (27) on the moments of the random step-size ensures that 0 ≤ α < 1. For the O(µ 2 x ) approximation in expression (28b) note from (23) that Iterating recursion (28a) gives Since 0 ≤ α < 1, there exists an iteration value I o large enough such that It follows that the variance E w i 2 converges exponentially to a region that is upper bounded by 2(

G. Mean-Square-Error Performance
We conclude from (28c) that the mean-square error (MSE) can be made as small as desired by using small step-sizes, µ x . In this section we derive a closed-form expression for the asymptotic mean-square error, which is more frequently called the mean-square deviation (MSD) and is defined as: Strictly speaking, the limit on the right-hand side of the above expression may not exist. A more accurate definition for the MSD appears in Eq. (4.86) of [5], namely, However, it was explained in [5][Sec. 4.5] that derivations that assume the validity of (34a) still lead to the same expression for the MSD to first-order in µ x as derivations that rely on the more formal definition (34b). We therefore continue with (34a) for simplicity of presentation. We explain below how an expression for the MSD can be obtained by following the energy conservation technique of [5], [28], [29], [50], [51]. For that purpose, we need to introduce two smoothness conditions.

Assumption II.4 (Smoothness conditions).
In addition to Assumptions II.1 and II.2, we assume that the Hessian matrix of the cost function and the noise covariance matrix defined by (17a) are locally Lipschitz continuous in a small neighborhood around w = w o : for small perturbations δw ≤ r and for some τ, τ 2 ≥ 0 and 1 ≤ κ ≤ 2.
The range of values for κ can be enlarged, e.g., to κ ∈ (0, 4]. The only change in allowing a wider range for κ is that the exponent of the higher-order term, O(µ 3/2 x ), that will appear in several performance expressions, as is the case with (41a)-(41b), will need to be adjusted from 3 2 to min{ 3 2 , 1 + κ 2 }, without affecting the first-order term that determines the MSD [3], [5], [52]. Therefore, it is sufficient to continue with κ ∈ [1, 2] to illustrate the key concepts though the MSD expressions will still be valid to first-order in µ x .
Using (9), it can be verified that condition (35a) translates into a global Lipschitz property relative to the minimizer w o , i.e., it will also hold that [3], [5]: for all w and for some τ ′ ≥ 0. For example, both conditions (35a)-(35b) are readily satisfied by mean-square-error costs. Using property (35c), we can now motivate a useful long-term model for the evolution of the error vector w i after sufficient iterations, i.e., for i ≫ 1. Indeed, let us reconsider recursion (29c) and introduce the deviation matrix: where the constant (symmetric and positive-definite) matrix H is defined as: Using (35c) and the fact that (E a) 2 ≤ E a 2 for any real-valued random variable, we can bound the conditional expectation of the norm of the perturbation term as follows: so that using (28c), we conclude that: We can deduce from this result that µ(i) c i−1 = O(µ 2 x ) asymptotically with high probability [3], [5]. To see this, let r c = mµ 2 x , for any constant integer m ≥ 1. Now, calling upon Markov's inequality [53]- [55], we conclude from (38b) that for i ≫ 1: This result shows that the probability of having µ(i) c i−1 bounded by r c can be made arbitrarily close to one by selecting a large enough value for m. Once the value for m has been fixed to meet a desired confidence level, then r c = O(µ 2 x ). This analysis, along with recursion (37a), motivate us to assess the mean-square performance of the error recursion (29c) by considering instead the following long-term model, which holds with high probability after sufficient iterations i ≫ 1: Working with iteration (39) is helpful because its dynamics is driven by the constant matrix H as opposed to the random matrix H i−1 in the original error recursion (29c). If desired, it can be shown that, under some technical conditions on the fourth-order moment of the gradient noise process, the MSD expression that will result from using (39) is x ) of the actual MSD expression for the original recursion (29c) -see [3], [5], [52] for a formal proof of this fact. Therefore, it is sufficient to rely on the long-term model (39) to obtain performance expressions that are accurate to first-order in µ x . Figure 2 provides a block-diagram representation for (39). Before explaining how model (39) can be used to assess the MSD, we remark that there is a second useful metric for evaluating the performance of stochastic gradient algorithms. This metric relates to the mean excess-cost; which is also called the excess-risk (ER) in the machine learning literature [31], [32] and the excess-mean-square-error (EMSE) in the adaptive filtering literature [26]- [28]. We denote it by the letters ER and define it as the average fluctuation of the cost function around its minimum value: Using the smoothness condition (35a), and the mean-value theorem [35], [49] again, it can be verified that [3], [5], [52]: ER Lemma II.2 (Mean-square-error performance). Assume the conditions under Assumptions II.1, II.2, II.3, and II.4 on the cost function, the gradient noise process, and the random step-size process hold. Assume further that the asynchronous step-size parameter µ x is sufficiently small to ensure mean-square stability as required by (27). Then, the MSD and ER metrics for the asynchronous stochastic-gradient algorithm (21) are well-approximated to first-order in µ x by the expressions: where R s and H are defined by (17b) and (36b), and where we are adding the superscript "asyn" for clarity in order to distinguish these measures from the corresponding measures in the synchronous case (mentioned below in (43a)-(43c)). Moreover, we derived earlier in (28b) the following expression for the convergence rate: Proof: We introduce the eigen-decomposition H = U ΛU T , where U is orthonormal and Λ is diagonal with positive entries, and rewrite (39) in terms of transformed quantities: . Let Σ denote an arbitrary M × M diagonal matrix with positive entries that we are free to choose. Then, equating the weighted squared norms of both sides of (42a) and taking expectations gives for i ≫ 1: From (17b), (19b), (28c), and (35b) we obtain: Therefore, substituting into (42b) gives for i → ∞: Since we are free to choose Σ, we let Σ = 1 2 Λ −1 and arrive at (41a) since w i 2 = w i 2 and U ΣU T = 1 2 H −1 . On the other hand, selecting Σ = 1 4 I M leads to (41b).
We recall our earlier remark that the synchronous stochastic-gradient recursion (12) can be viewed as a special case of the asynchronous update (21) by setting σ 2 µ = 0 and µ x =μ ≡ µ. Substituting these values into (41a)-(41c), we obtain for the synchronous implementation (12): which are the same expressions presented in [1] and which agree with classical results for LMS adaptation [56]- [61]. The matrices R s that appear in these expressions, and in (41a)-(41c), were defined earlier in (17b) and they correspond to the covariance matrices of the gradient noise processes in their respective (synchronous or asynchronous) implementations. The examples that follow show how expressions (41a) and (41b) can be used to recover performance metrics for mean-square-error adaptation and learning under random updates.
Example II.7 (Performance of asynchronous LMS adaptation). We reconsider the LMS recursion (13f) albeit in an asynchronous mode of operation, namely, We know from Example II.3 and (18a) that this situation corresponds to H = 2R u and R s = 4σ 2 v R u . Substituting into (41a) and (41b) leads to the following expressions for the MSD and EMSE of the asynchronous LMS filter: where here, and elsewhere, we will be using the notation ≈ to indicate that we are ignoring higher-order terms in µ x . For example, let us assume a Bernoulli update model for µ(i) where the filter updates with probability p µ using a step-size value µ or stays inactive otherwise. In this case, we conclude from (24b) that µ x = µ so that the above performance expressions for the MSD and EMSE metrics will coincide with the values obtained in the synchronous case as well. In other words, the steady-state performance levels are not affected whether the algorithm learns in a synchronous or asynchronous manner. However, the convergence rate is affected sinceμ = µp µ and, therefore, It follows that asynchronous LMS adaptation attains the same performance levels as synchronous LMS adaptation albeit at a slower convergence rate.
We may alternatively compare the performance of the synchronous and asynchronous implementations by fixing their convergence rates to the same value. Thus, consider now a second random update scheme with meanμ and let us set µ =μ. That is, the step-size used by the synchronous implementation is set equal to the mean step-size used by the asynchronous implementation. Then, in this case, we will get α sync LMS = α asyn LMS so that the convergence rates coincide to first-order. However, it now holds that MSD asyn LMS > MSD sync LMS sinceμ x > µ so that some deterioration in MSD performance occurs.
Example II.8 (Performance of asynchronous online learners). Consider a stand-alone learner receiving a streaming sequence of independent data vectors {x i , i ≥ 0} that arise from some fixed probability distribution X . The goal is to learn the vector w o that optimizes some ν-strongly convex risk function J(w) defined in terms of a loss function [62], [63]: In an asynchronous environment, the learner seeks w o by running the stochastic-gradient algorithm with random step-sizes: The gradient noise vector is still given by Since ∇ w J(w o ) = 0, and since the distribution of x i is stationary, it follows that the covariance matrix of s i (w o ) is constant and given by The excess-risk measure that will result from this stochastic implementation is then given by (41b) so that ER = O(µ x ).

III. CENTRALIZED ADAPTATION AND LEARNING
The discussion in the previous section establishes the mean-square stability of stand-alone adaptive agents for small step-sizes (Lemma II.1), and provides expressions for their MSD and ER metrics (Lemma II.2). We now examine two situations involving a multitude of similar agents. In the first scenario, each agent senses data and analyzes it independently of the other agents. We refer to this mode of operation as non-cooperative processing. In the second scenario, the agents transmit the collected data for processing at a fusion center. We refer to this mode of operation as centralized or batch processing. We motivate the discussion by considering first the case of mean-square-error costs. Subsequently, we extend the results to more general costs.

A. Non-Cooperative MSE Processing
Thus, consider separate agents, labeled k = 1, 2, . . . , N . Each agent, k, receives streaming data where we are using the subscript k to index the data at agent k. We assume that the data at each agent satisfies the same statistical properties as in Example II.3, and the same linear regression model (13a) with a common w o albeit with noise v k (i). We denote the statistical moments of the data at agent k by We further assume in this motivating example that the R u,k are uniform across the agents so that In this way, the cost , which is associated with agent k, will satisfy a condition similar to (9) with the corresponding parameters {ν, δ} given by (cf. (10a)): Now, if each agent runs the asynchronous LMS learning rule (44) to estimate w o on its own then, according to (45a), each agent k will attain an individual MSD level that is given by where we are further assuming that the parameter µ x is uniform across the agents to enable a meaningful comparison. Moreover, according to (28b), agent k will converge towards this level at a rate dictated by: The subscript "ncop" is used in (49b) and (49c) to indicate that these expressions are for the non-cooperative mode of operation. It is seen from (49b) that agents with noisier data (i.e., larger σ 2 v,k ) will perform worse and have larger MSD levels than agents with cleaner data. We are going to show in later sections that cooperation among the agents, whereby agents share information with their neighbors, can help enhance their individual performance levels.

B. Centralized MSE Processing
Let us now contrast the above non-cooperative solution with a centralized implementation whereby, at every iteration i, the N agents transmit their raw data {d k (i), u k,i } to a fusion center for processing. In a synchronous environment, once the fusion center receives the raw data, it can run a standard stochastic-gradient update of the form: where µ is the constant step-size, and the term multiplying µ can be seen to correspond to a sample average of several approximate gradient vectors. The analysis in [1], [5] showed that the MSD performance that results from this implementation is given by (using expression (64a) with H k = 2R u and R s, Moreover, using expression (64b) given further ahead, this centralized solution will converge towards the above MSD level at the same rate as the non-cooperative solution: In an asynchronous environment, there are now several random events that can interfere with the operation of the fusion center. Let us consider initially one particular random event that corresponds to the situation in which the fusion center may or may not update at any particular iteration (e.g., due to some power-saving strategy). In a manner similar to (21), we may represent this scenario by writing: with a random step-size process, µ(i); this process is again assumed to satisfy the conditions under Assumption II.3. The analysis in the sequel will show that the MSD performance that results from this implementation is given by: where we are adding the superscript "1" to indicate that this is a preliminary result pertaining to the particular asynchronous implementation (51a). We will be generalizing this result soon, at which point we will drop the superscript "1". Likewise, using expression (64b) given further ahead, this version of the asynchronous centralized solution converges towards the above MSD level at the same rate as the non-cooperative solution (49c) and the synchronous version (50c): Observe from (50b) and (51b) that the MSD level attained by the centralized solution is proportional to 1/N times the average noise power across all agents. This scaled average noise power can be larger than some of the individual noise variances and smaller than the remaining noise variances. This example shows that it does not generally hold that centralized stochastic-gradient implementations outperform all individual non-cooperative agents [64]. Now, more generally, observe that in the synchronous batch processing case (50a), as well as in the asynchronous implementation (51a), the data collected from the various agents are equally aggregated with a weighting factor equal to 1/N . We can incorporate a second type of random events besides the random variations in the stepsize parameter. This second source of uncertainty involves the possibility of failure (or weakening) in the links connecting the individual agents to the fusion center (e.g., due to fading or outage, or perhaps due to the agents themselves deciding to enter into a sleep mode according to some power-saving policy). We can capture these possibilities by extending formulation (51a) in the following manner: where µ(i) continues to be a random step-size process, but now the coefficients {π k (i); k = 1, 2, . . . , N } are new random fusion coefficients that satisfy: Assumption III.1 (Conditions on random fusion coefficients). It is assumed that π k (i) is independent of µ(i) and of other random variables in the learning algorithm. The fusion coefficients are also independent over time, namely, π k (i) and π ℓ (j) are independent for any i = j. For the same time i, the coefficients {π k (i)} are correlated over space in view of the first requirement in (52b). The mean and co-variance(s) of each π k (i) are denoted by: for all k, ℓ = 1, 2, . . . , N and all i ≥ 0. When k = ℓ, the scalar c π,kk corresponds to the variance of π k (i) and, therefore, we shall also use the alternative notation σ 2 π,k for this case: It is straightforward to verify that the first and second-order moments of the coefficients {π k (i)} satisfy: c π,kℓ = 0, for any k (54c) Note that the earlier asynchronous implementation (51a) is a special case of the more general formulation (52a) when the mean of the random fusion coefficients are chosen asπ k = 1/N and their variances are set to zero, σ 2 π,k = 0. In order to enable a fair and meaningful comparison among the three centralized implementations (50a), (51a), and (52a), we shall assume that the means of the fusion coefficients in the latest implementation are set tō π k = 1/N (results for arbitrary mean values are listed in Example III.1 further ahead and appear in [4]). We will show in the sequel that for this choice ofπ k , the MSD performance of implementation (52a) is given by: Moreover, using expression (62b) given further ahead, this asynchronous centralized solution will converge towards the above MSD level at the same rate as the non-cooperative solution (49c) and the synchronous version (50c):

C. Stochastic-Gradient Centralized Solution
The previous two sections focused on mean-square-error costs. We now extend the conclusions to more general costs. Thus, consider a collection of N agents, each with an individual convex cost function, J k (w). The objective is to determine the unique minimizer w o of the aggregate cost: It is now the above aggregate cost, J glob (w), that will be required to satisfy the conditions of Assumptions II.1 and II.4 relative to some parameters {ν c , δ c , τ c }, with the subscript "c" used to indicate that these factors correspond to the centralized implementation. Under these conditions, the cost J glob (w) will have a unique minimizer, which we continue to denote by w o . We will not be requiring each individual cost, J k (w), to be strongly convex. It is sufficient for at least one of these costs to be strongly convex while the remaining costs can be convex; this condition ensures the strong convexity of J glob (w). Moreover, although some individual costs may not have a unique minimizer, we require in this exposition that w o is one of their minima so that all individual costs share a minimum at location w o ; the treatment in [1], [5], [65] considers the more general case in which the minimizers of the individual costs {J k (w)} can be different and need not contain a common location, w o . There are many centralized solutions that can be used to determine the unique minimizer w o of (56), with some solution techniques being more powerful than other techniques. Nevertheless, we shall focus on centralized implementations of the stochastic-gradient type. The reason we consider the same class of stochastic gradient algorithms for non-cooperative, centralized, and distributed solutions in this chapter is to enable a meaningful comparison among the various implementations. Thus, we first consider a synchronous centralized strategy of the following form: with a constant step-size, µ. When the fusion center employs random step-sizes, the above solution is replaced by: where the process µ(i) now satisfies Assumption II.3. More generally, the asynchronous implementation can employ random fusion coefficients as well such as: where the coefficients {π k (i)} satisfy Assumption III.1 with means

D. Performance of Centralized Solution
To examine the performance of the asynchronous implementation (57c)-(57d), we proceed in two steps. First, we identify the gradient noise that is present in the recursion; it is equal to the difference between the true gradient vector for the global cost, J glob (w), defined by (56) and its approximation. Second, we argue that (57c) has a form similar to the single-agent stochastic-gradient algorithm (21) and, therefore, invoke earlier results to write down performance metrics for the centralized solution (57c)-(57d).
We start by introducing the individual gradient noise processes: for k = 1, 2, . . . , N . We assume that these noises satisfy conditions similar to Assumption II.2 with parameters and Additionally, we assume that the gradient noise components across the agents are uncorrelated with each other: Using these gradient noise terms, it is straightforward to verify that recursion (57c) can be rewritten as: where s i (w i−1 ) denotes the overall gradient noise; its expression is given by Since iteration (59) has the form of a stochastic gradient recursion with random update similar to (21), we can infer its mean-square-error behavior from Lemmas II.1 and II.2 if the noise process s i (w i−1 ) can be shown to satisfy conditions similar to Assumption II.2 with some parameters {β 2 c , σ 2 s }. Indeed, starting from (60), some algebra will show that The following result now follows from Lemmas II.1 and II.2 [4], [5].
. For any µ x /N < µ o , the iterates generated by the asynchronous centralized solution (57c)-(57d) satisfy: where the scalar α satisfies 0 ≤ α < 1 and is given by It follows from (62a) that for sufficiently small step-size parameter µ x ≪ 1: Moreover, under smoothness conditions similar to (35a) for J glob (w) for some parameter τ c ≥ 0, and similar to (35b) for the individual gradient noise covariance matrices, it holds for small µ x that: where We can recover from the expressions in the lemma, performance results for the particular asynchronous implementation described earlier by (57b) by setting σ 2 π,k = 0 so that It is seen from (63) and (64a) that when the fusion center operates under the more general asynchronous policy (57c), the additional randomness in the fusion coefficients {π k (i)} degrades the MSD performance relative to (64a) due to the presence of the factor 1 + N 2 σ 2 π,k > 1, i.e., to first-order in µ x we have: MSD asyn,1 < MSD asyn (65a) In comparison, the added randomness due to the {π k (i)} does not have significant impact on the convergence rate since it is straightforward to see that, to first order in the step-size parameter µ x : Likewise, settingμ k = µ and σ 2 π,k = 0, the performance of the synchronous centralized solution (57a) is obtained as a special case of the results in the lemma: These expressions agree with the performance results presented in [1], [5].
Example III.1 (Case of general mean-values). Although not used in this chapter, we remark in passing that if the mean values {π k } are not fixed at 1/N , as was required by (57d), then the MSD performance of the asynchronous centralized solution (57c) will instead be given by

E. Comparison with Non-Cooperative Processing
We can now compare the performance of the asynchronous centralized solution (57c) against the performance of non-cooperative processing when agents act independently of each other and run the recursion: This comparison is meaningful when all agents share the same unique minimizer so that we can compare how well the individual agents are able to recover the same w o as the centralized solution. For this reason, we re-introduce the requirement that all individual costs {J k (w)} are ν-strongly convex with a uniform parameter ν. Since J glob (w) is the aggregate sum of the individual costs, then we can set the lower bound ν c for the Hessian of J glob (w) at ν c = N ν. From expressions (28b) and (64b) we then conclude that, for a sufficiently small µ x , the convergence rates of the asynchronous non-cooperative solution (68a) and the asynchronous centralized solution with random update (57b) will be similar: Moreover, we observe from (41a) that the average MSD level across N non-cooperative asynchronous agents is given by so that comparing with (64a), some simple algebra allows us to conclude that, for small step-sizes and to first order in µ x : MSD asyn,1 cent That is, while the asynchronous centralized solution (57b) with random updates need not outperform every individual non-cooperative agent in general, its performance outperforms the average performance across all non-cooperative agents. The next example illustrates this result by considering the scenario where all agents have the same Hessian matrices at w = w o , namely, H k ≡ H for k = 1, 2, . . . , N . This situation occurs, for example, when the individual costs are identical across the agents, say, J k (w) ≡ J(w), as is common in machine learning applications. This situation also occurs for the mean-square-error costs we considered earlier in this section when the regression covariance matrices, {R u,k }, are uniform across all agents, i.e., R u,k ≡ R u for k = 1, 2, . . . , N . In these cases with uniform Hessian matrices H k , the example below establishes that the asynchronous centralized solution (57b) with random updates improves over the average MSD performance of the non-cooperative solution (68a) by a factor of N .
Example III.2 (N -fold improvement in MSD performance). Consider a collection of N agents whose individual cost functions, J k (w), are ν-strongly convex and are minimized at the same location w = w o . The costs are also assumed to have identical Hessian matrices at w = w o , i.e., H k ≡ H. Then, using (64a), the MSD of the asynchronous centralized implementation (57b) with random updates is given by: Example III.3 (Random fusion can degrade MSD performance). Although the convergence rates of the asynchronous centralized solution (57c) and the non-cooperative solution (68a) agree to first-order in µ x , the relation between their MSD values is indefinite (contrary to (68d)), as illustrated by the following example.
Consider the same setting of Example III.2 and assume further that the variances of the random fusion coefficients are uniform, i.e., σ 2 π,k ≡ σ 2 π . Then, using (63), the MSD of the asynchronous centralized implementation (57c) is given by Therefore, to first-order in µ x , we find that In other words, if the variance of the random fusion coefficients is large enough, then the centralized solution will generally have degraded performance relative to the non-cooperative solution (which is an expected result).
Example III.4 (Fully-connected networks). In preparation for the discussion on networked agents, it is useful to describe one extreme situation where a collection of N agents are fully connected to each other -see Fig. 3. In this case, each agent is able to access the data from all other agents and, therefore, each individual agent can run a synchronous or asynchronous centralized implementation, say, one of the same form as (57b): When this happens, each agent will attain the same performance level as that of the asynchronous centralized solution (57b). Two observations are in place. First, note from (72) that the information that agent k is receiving from all other agents is their gradient vector approximations. Obviously, other pieces of information could be shared among the agents, such as their iterates {w ℓ,i−1 }. Second, note that the right-most term multiplying µ(i) in (72) corresponds to a convex combination of the approximate gradients from the various agents, with the combination coefficients being uniform and all equal to 1/N . In general, there is no need for these combination weights to be identical. Even more importantly, agents do not need to have access to information from all other agents in the network. We are going to see in the sequel that interactions with a limited number of neighbors is sufficient for the agents to attain performance that is comparable to that of the centralized solution. Figure 4 shows a simple selection of connected topologies for five agents. The leftmost panel corresponds to the noncooperative case and the rightmost panel corresponds to the fully-connected case. The panels in between illustrate some other topologies. In the coming sections, we are going to present results that allow us to answer useful questions about such networked agents such as: (a) Which topology has best performance in terms of mean-square error and convergence rate? (b) Given a connected topology, can it be made to approach the performance of the centralized solution? (c) Which aspects of the topology influence performance? (d) Which aspects of the combination weights (policy) influence performance? (e) Can different topologies deliver similar performance levels? (f) Is cooperation always beneficial? and (g) If the individual agents are able to solve the inference task individually in a stable manner, does it follow that the connected network will remain stable regardless of the topology and regardless of the cooperation strategy?

IV. SYNCHRONOUS MULTI-AGENT ADAPTATION AND LEARNING
In this section, we describe distributed strategies of the consensus (e.g., [7], [11], [15] and [68]- [77]) and diffusion (e.g., [1], [5]- [7], [41], [78]- [80]) types. These strategies rely solely on localized interactions among neighboring agents, and they can be used to seek the minimizer of (56). We first describe the network model. Figure 5 shows an example of a network consisting of N connected agents, labeled k = 1, 2, . . . , N . Following the presentation from [5], [7], the network is represented by a graph consisting of N vertices (representing the agents) and a set of edges connecting the agents to each other. An edge that connects an agent to itself is called a self-loop. The neighborhood of an agent k is denoted by N k and it consists of all agents that are connected to k by an edge. Any two neighboring agents k and ℓ have the ability to share information over the edge connecting them.

A. Strongly-Connected Networks
We assume an undirected graph so that if agent k is a neighbor of agent ℓ, then agent ℓ is also a neighbor of agent k. We assign a pair of nonnegative scaling weights, {a kℓ , a ℓk }, to the edge connecting k and ℓ. The scalar a ℓk is used by agent k to scale the data it receives from agent ℓ; this scaling can be interpreted as a measure of the confidence level that agent k assigns to its interaction with agent ℓ. Likewise, a kℓ is used by agent ℓ to scale the data it receives from agent k. The weights {a kℓ , a ℓk } can be different so that the exchange of information between the neighboring agents {k, ℓ} need not be symmetrical. One or both weights can also be zero.
A network is said to be connected if paths with nonzero scaling weights can be found linking any two distinct agents in both directions, either directly when they are neighbors or by passing through intermediate agents when they are not neighbors. In this way, information can flow in both directions between any two agents in the network, although the forward path from an agent k to some other agent ℓ need not be the same as the backward path from ℓ to k. A strongly-connected network is a connected network with at least one non-trivial self-loop, meaning that a kk > 0 for some agent k.
The strong connectivity of a network translates into a useful property on the combination weights. Assume we collect the coefficients {a ℓk } into an N × N matrix A = [a ℓk ], such that the entries on the k-th column of A contain the coefficients used by agent k to scale data arriving from its neighbors ℓ ∈ N k ; we set a ℓk = 0 if ℓ / ∈ N k . We refer to A as the combination matrix or policy. It turns out that combination matrices that correspond to strongly-connected networks are primitive -an N × N matrix A with nonnegative entries is said to be primitive if there exists some finite integer n o > 0 such that all entries of A no are strictly positive [5], [7], [81].

B. Distributed Optimization
Network cooperation can be exploited to solve adaptation, learning, and optimization problems in a decentralized manner in response to streaming data. To explain how cooperation can be achieved, we start by associating with each agent k a twice-differentiable cost function J k (w) : IR M ×1 → IR. The objective of the network of agents is to seek the unique minimizer of the aggregate cost function, J glob (w), defined by (56). Now, however, we seek a distributed (as opposed to a centralized) solution. In a distributed implementation, each agent k can only rely on its own data and on data from its neighbors.
We continue to assume that J glob (w) satisfies the conditions of Assumptions II.1 and II.4 with parameters {ν d , δ d , τ d }, with the subscript "d" now used to indicate that these parameters are related to the distributed implementation. Under these conditions, the cost J glob (w) will have a unique minimizer, which we continue to denote by w o . For simplicity of presentation, we will also assume in the remainder of this chapter that the individual costs J k (w) are strongly convex as well. These costs can be distinct across the agents or they can all be identical, i.e., J k (w) ≡ J(w) for k = 1, 2, . . . , N ; in the latter situation, the problem of minimizing (56) would correspond to the case in which the agents work together to optimize the same cost function. If we let w o k denote the minimizer of J k (w), we continue to assume for this exposition that each J k (w) is also minimized at w o : The case where the individual costs are only convex and need not be strongly convex is discussed in [1], [5], [42]; most of the results and conclusions continue to hold but the derivations become more technical. Likewise, the case in which the individual costs need not share minimizers is discussed in these references. In that case, as was already shown in [65], the iterates w k,i by the individual agents will not approach the minimizer w o of (56), but rather the minimizer w ⋆ of a weighted aggregate function with some positive scaling weights {q k }. This function was defined in Eq. (57b) of [1]. It was also explained in [1] how this convergence property adds a useful degree of freedom to the operation of the network, and how it can be exploited advantageously to steer the network to converge to desirable limit points, including to w o , through the selection of the combination policy A (which determines the scaling weights {q k } in the weighted aggregate cost). It is, nevertheless, sufficient for the purposes of this chapter to continue with the case (73). There are many important situations in practice where the minimizers of individual costs coincide with each other. For instance, examples abound where agents need to work cooperatively to attain a common objective such as tracking a target, locating a food source, or evading a predator (see, e.g., [6], [82], [83]). This scenario is also common in machine learning problems [31], [32], [84]- [87] when data samples at the various agents are generated by a common distribution parameterized by some vector, w o . One such situation is illustrated in the next example.
Example IV.1 (Mean-square-error (MSE) networks). Consider the same setting of Example II.3 except that we now have N agents observing streaming data {d k (i), u k,i } that satisfy the regression model (13a) with regression covariance matrices R u,k = E u T k,i u k,i > 0 and with the same unknown w o , i.e., The individual mean-square-error costs are defined by and are strongly convex in this case, with the minimizer of each J k (w) occurring at If we multiply both sides of (74a) by u T k,i from the left, and take expectations, we find that w o satisfies This relation shows that the unknown w o from (74a) satisfies the same expression as w o k in (74c), for any k = 1, 2, . . . , N , so that we must have w o = w o k . Therefore, this example amounts to a situation where all costs {J k (w)} attain their minima at the same location, w o .
We shall use the network model of this example to illustrate other results in the chapter. For ease of reference, we shall refer to strongly-connected networks with agents receiving data according to model (74a) and seeking to estimate w o by adopting the mean-square-error costs J k (w) defined above, as mean-square-error (MSE) networks. We assume for these networks that the measurement noise process v k (i) is temporally white and independent over space so that in terms of the Kronecker delta δ k,ℓ . Likewise, we assume that the regression data u k,i is temporally white and independent over space so that E u T k,i u ℓ,j = R u,k δ k,ℓ δ i,j Moreover, the measurement noise v k (i) and the regression data u ℓ,j are independent of each other for all k, ℓ, i, j. These statistical conditions help facilitate the analysis of such networks.
In the next subsections we list synchronous distributed algorithms of the consensus and diffusion types for the optimization of J glob (w). We only list the algorithms here; for motivation and justifications, the reader may refer to the treatments in [1], [5]. Moreover, Sec. V.D of [1] provides commentary on several other related works in the literature, in addition to the history and evolution of the consensus and diffusion strategies.

C. Synchronous Consensus Strategy
Let w k,i denote the iterate that is available at agent k at iteration i; this iterate serves as the estimate for w o . The consensus iteration at each agent k is described by the following construction (see, e.g., [11], [15], [36], [69], [71], [77]): where the {µ k } are individual step-size parameters, and where the combination coefficients {a ℓk } that appear in (78) are nonnegative scalars that are required to satisfy the following conditions for each agent k = 1, 2, . . . , N : Condition (79a) implies that the combination matrix A = [a ℓk ] satisfies where 1 denotes the vector with all entries equal to one. We say that A is left-stochastic. One useful property of left-stochastic matrices is that their spectral radius is equal to one [7], [81], [88]- [90]: An equivalent representation that is useful for later analysis is to rewrite the consensus iteration (78) as shown in the following listing, where the intermediate iterate that results from the neighborhood combination is denoted by ψ k,i−1 . Observe that the gradient vector in the consensus implementation (80) is evaluated at w k,i−1 and not ψ k,i−1 .
Consensus strategy for distributed adaptation for each time instant i ≥ 0: each agent k = 1, 2, . . . , N performs the update: We remark that one way to motivate the consensus update (78) is to start from the non-cooperative step (68a) and replace the first iterate w k,i−1 by the convex combination used in (78).
Example IV.2 (Consensus LMS networks). For the MSE network of Example IV.1, the consensus strategy reduces to:

D. Synchronous Diffusion Strategies
There is an inherent asymmetry in the consensus construction. Observe from the computation of w k,i in (80) that the update starts from ψ k,i−1 and corrects it by the approximate gradient vector evaluated at w k,i−1 (and not at ψ k,i−1 ). This asymmetry will be shown later, e.g., in Example VII.4, to be problematic when the consensus strategy is used for adaptation and learning over networks. This is because the asymmetry can cause an unstable growth in the state of the network [80] -see also the explanations in [1] and [5][Sec. 10.6]. Diffusion strategies remove the asymmetry problem.
Combine-then-Adapt (CTA) Diffusion. In the CTA formulation of the diffusion strategy, the same iterate ψ k,i−1 is used to compute w k,i , thus leading to description (82a) where the gradient vector is evaluated at ψ k,i−1 as well. The reason for the name "Combine-then-Adapt" is that the first step in (82a) involves a combination step, while the second step involves an adaptation step. The reason for the qualification "diffusion" is that the use of ψ k,i−1 to evaluate the gradient vector allows information to diffuse more thoroughly through the network. This is because information is not only being diffused through the aggregation of the neighborhood iterates, but also through the evaluation of the gradient vector at the aggregate state value.
Diffusion strategy for distributed adaptation (CTA) for each time instant i ≥ 0: each agent k = 1, 2, . . . , N performs the update: Adapt-then-Combine (ATC) Diffusion. A similar implementation can be obtained by switching the order of the combination and adaptation steps in (82a), as shown in the listing (82b). The structure of the CTA and ATC strategies are fundamentally identical: the difference lies in which variable we choose to correspond to the updated iterate w k,i . In ATC, we choose the result of the combination step to be w k,i , whereas in CTA we choose the result of the adaptation step to be w k,i .
Diffusion strategy for distributed adaptation (ATC) for each time instant i ≥ 0: each agent k = 1, 2, . . . , N performs the update: One main motivation for the introduction of the diffusion strategies (82a) and (82b) is the fact that they enable single time-scale distributed learning from streaming data under constant step-size adaptation and in a stable manner [2], [65], [78], [79] and [91]- [95] -see also [5][Chs. [9][10][11]; the diffusion strategies further allow A to be leftstochastic, which permit larger modes of cooperation than doubly-stochastic policies. The CTA diffusion strategy (83a) was first introduced for mean-square-error estimation problems in [78], [91]- [93]. The ATC diffusion structure (83b), with adaptation preceding combination, appeared in the work [96] on adaptive distributed least-squares schemes and also in the works [79], [95], [97], [98] on distributed mean-square-error and state-space estimation methods. The CTA structure (82a) with an iteration dependent step-size that decays to zero, µ(i) → 0, was employed in [13], [99], [100] to solve distributed optimization problems that require all agents to reach agreement. The ATC form (82b), also with an iteration dependent sequence µ(i) that decays to zero, was employed in [101], [102] to ensure almost-sure convergence and agreement among agents.

Example IV.3 (Diffusion LMS networks).
For the MSE network of Example IV.1, the ATC and CTA diffusion strategies reduce to: Example IV.4 (Diffusion logistic regression). We revisit the pattern classification problem from Example II.2, where we consider a collection of N networked agents cooperating with each other to solve the logistic regression problem. Each agent receives streaming data {γ k (i), h k,i }, where the variable γ k (i) assumes the values ±1 and designates the class that the feature vector h k,i belongs to. The objective is to use the training data to determine the vector w o that minimizes the cost under the assumption of joint wide-sense stationarity over the random data. It is straightforward to verify that the ATC diffusion strategy (82b) reduces to the following form in this case: V. ASYNCHRONOUS MULTI-AGENT ADAPTATION AND LEARNING There are various ways by which asynchronous events can be introduced into the operation of a distributed strategy. Without loss in generality, we illustrate the model for asynchronous operation by describing it for the ATC diffusion strategy (82b); similar constructions apply to CTA diffusion (82a) and consensus (80).

A. Asynchronous Model
In a first instance, we model the step-size parameters as random variables and replace (82b) by: In this model, the neighborhoods and the network topology remain fixed and only the µ k (i) assume random values. The step-sizes can vary across the agents and, therefore, their means and variances become agent-dependent. Moreover, the step-sizes across agents can be correlated with each other. We therefore denote the first and secondorder moments of the step-size parameters by: When ℓ = k, the scalar c µ,kk coincides with the variance of µ k (i), i.e., c µ,kk = σ 2 µ,k ≥ 0. On the other hand, if the step-sizes across the agents happen to be uncorrelated, then c µ,kℓ = 0 for k = ℓ.
More broadly, we can allow for random variations in the neighborhoods (and, hence, in the network structure), and random variations in the combination coefficients as well. We capture this more general asynchronous implementation by writing: where the combination coefficients {a ℓ k (i)} are now random and, moreover, the symbol N k,i denotes the randomlychanging neighborhood of agent k at time i. These neighborhoods become random because the random variations in the combination coefficients can turn links on and off depending on the values of the {a ℓ,k (i)}. We continue to require the combination coefficients {a ℓk (i)} to satisfy the same structural constraint as given before by (79a), i.e., ℓ∈N k,i a ℓk (i) = 1, and a ℓk (i) > 0, if ℓ ∈ N k,i a ℓk (i) = 0, otherwise Since these coefficients are now random, we denote their first and second-order moments by: When ℓ = n and k = m, the scalar c a,ℓk,nm coincides with the variance of a ℓk (i), i.e., c a,ℓk,nm = σ 2 a,ℓk ≥ 0.
whereas the asynchronous ATC diffusion strategy (88) reduces to: We can view implementation (91a) as a special case of the asynchronous update (91b) when the variances of the random combination coefficients {a ℓk (i)} are set to zero.
Example V.2 (Asynchronous consensus LMS networks). Similarly, for the same MSE network of Example IV.1, the asynchronous consensus strategy is given by

B. Mean Graph
We refer to the topology that corresponds to the average combination coefficients {ā ℓk } as the mean graph, which is fixed over time. For each agent k, the neighborhood defined by the mean graph is denoted by N k . It is straightforward to verify that the mean combination coefficientsā ℓk satisfy the following constraints over the mean graph (compare with (79a) and (89)): One example of a random network with two equally probable realizations and its mean graph is shown in Fig. 6 [2]. The letter ω is used to index the sample space of the random matrix A i . There is one useful result that relates the random neighborhoods {N k,i } from (88) to the neighborhoods {N k } from the mean graph. It is not difficult to verify that N k is equal to the union of all possible realizations for the random neighborhoods N k,i ; this property is already illustrated by the example of Fig. 6: for any k, where Ω denotes the sample space for N k,i .

C. Random Combination Policy
The first and second-order moments of the combination coefficients will play an important role in characterizing the stability and mean-square-error performance of the asynchronous network (88). We collect these moments in matrix form as follows. We first group the combination coefficients into a matrix: The sequence {A i , i ≥ 0} represents a stochastic process consisting of left-stochastic random matrices whose entries satisfy the conditions in (89) at every time i. We subsequently introduce the mean and Kronecker-covariance matrix of A i and assume these quantities are constant over time; we denote them by the N × N matrixĀ and the N 2 × N 2 matrix C A , respectively:Ā The matrix C A is not a conventional covariance matrix and is not necessarily Hermitian. The reason for its introduction is because it captures the correlations of each entry of A i with all other entries in A i . For example, for a network with N = 2 agents, the entries ofĀ and C A will be given by: We thus see that the (ℓ, k)−th block of C A contains the covariance coefficients of a ℓk with all other entries of A i . One useful property of the matrices {Ā, C A } so defined is that their elements are nonnegative and the following matrices are left-stochastic: Assumption V.1 (Asynchronous network model). It is assumed that the random processes {µ k (i), a ℓm (j)} are independent of each other for all k, ℓ, m, i, and j. They are also independent of any other random variable in the learning algorithm. The asynchronous network model described in this section covers many situations of practical interest. For example, we can choose the sample space for each step-size µ k (i) to be the binary choice {0, µ} to model random "on-off" behavior at each agent k for the purpose of saving power, waiting for data, or even due to random agent failures. Similarly, we can choose the sample space for each combination coefficient a ℓ k (i), ℓ ∈ N k \{k}, to be {0, a ℓk } to model a random "on-off" status for the link from agent ℓ to agent k at time i for the purpose of either saving communication cost or due to random link failures. Note that the convex constraint (89) can always be satisfied by adjusting the value of a kk (i) according to the realizations of {a ℓk (i); ℓ ∈ N k,i \{k}}.
Example V.3 (The spatially-uncorrelated model). A useful special case of the asynchronous network model of this section is the spatially-uncorrelated model. In this case, at each iteration i, the random step-sizes {µ k (i); k = 1, 2, . . . , N } are uncorrelated with each other across the network, and the random combination coefficients {a ℓk (i); ℓ = k, k = 1, 2, . . . , N } are also uncorrelated with each other across the network. Then, it can be verified that the covariances {c µ,kℓ } in (87c) and {c a,ℓk,nm } in (90c) will be fully determined by the variances {σ 2 µ , σ 2 a,ℓk }: and c a,ℓk,nm =

D. Perron Vectors
Now that we introduced the network model, we can move on to examine the effect of network cooperation on performance. Some interesting patterns of behavior arise when agents cooperate to solve a global optimization problem in a distributed manner from streaming data [1], [5]. For example, one interesting result established in [2], [3], [42], [65], [103] is that the effect of the network topology on performance is captured by the Perron vector of the combination policy. This vector turns out to summarize the influence of the topology on performance so much so that different topologies with similar Perron vectors will end up delivering similar performance. We explained the role of Perron vectors in the context of synchronous adaptation and learning in [1], [5]. Here we focus on asynchronous networks. In this case, two Perron vectors will be needed since the randomness in the combination policy is now represented by two moment matrices,Ā and C A . In the synchronous case, only one Perron vector was necessary since the combination policy was fixed and described by a matrix A. Although we are focusing on the asynchronous case in the sequel, we will be able to recover results for synchronous networks as special cases.
Let us first recall the definition of Perron vectors for synchronous networks, say, of the form described by (82b) with combination policy A. We assume the network is strongly-connected. In this case, the left-stochastic matrix A will be primitive. For such primitive matrices, it follows from the Perron-Frobenius Theorem [81] that: (a) the matrix A will have a single eigenvalue at one; (b) all other eigenvalues of A will be strictly inside the unit circle so that ρ(A) = 1; and (c) with proper sign scaling, all entries of the right-eigenvector of A corresponding to the single eigenvalue at one will be positive. Let p denote this right-eigenvector with its entries {p k } normalized to add up to one, i.e., Ap = p, 1 T p = 1, p k > 0, k = 1, 2, . . . , N We refer to p as the Perron eigenvector of A. It was explained in [1], [5] how the entries of this vector determine the mean-square-error performance and convergence rate of the network; these results will be revisited further ahead when we recover them as special cases of the asynchronous results.
On the other hand, for an asynchronous implementation, the individual realizations of the random combination matrix A i in (88) need not be primitive. In this context, we will require a form of primitiveness to hold on average as follows.
Definition V.1 (Strongly-connected asynchronous model). We say that an asynchronous model with random combination coefficients {a ℓk (i)} is strongly-connected if the Kronecker-covariance matrix given byĀ ⊗Ā + C A is primitive.
Observe that if we setĀ = A and C A = 0, then we recover the condition for strong-connectedness in the synchronous case, namely, that A should be primitive.
Definition V.1 means that the directed graph (digraph) associated with the matrixĀ⊗Ā+C A is strongly-connected (e.g., [89, pp. 30,34] and [5]). It is straightforward to check from the definition of C A in (96b) that so that the digraph associated withĀ ⊗Ā + C A is the union of all possible digraphs associated with the realizations of A i ⊗ A i [104, p. 29]. Therefore, as explained in [2]- [4], definition V.1 amounts to an assumption that the union of all possible digraphs associated with the realizations of A i ⊗ A i is strongly-connected. As illustrated in Fig. 7, this condition still allows individual digraphs associated with realizations of A i to be weakly-connected with or without self-loops or even to be disconnected [3].
is strongly-connected due to the existence of the source and sink nodes. However, the digraph associated with E (A i ⊗ A i ), which is the union of the first two digraphs, is strongly-connected, where information can flow in any direction through the network. This is a modified version of Fig. 7 from [3]; extracted with permission.
It follows from property (98) and Definition V.1, the matrixĀ ⊗Ā + C A is left-stochastic and primitive. It can be verified that mean matrixĀ is also primitive ifĀ ⊗Ā + C A is primitive (although the converse is not true). Therefore, the matrixĀ is both left-stochastic and primitive. We denote the Perron eigenvector ofĀ ⊗Ā + C A by p c ∈ IR N 2 ×1 , which satisfies: Likewise, we denote the Perron eigenvector ofĀ byp ∈ IR N ×1 , which satisfies: Observe that if we set C A = 0 andĀ = A, then we recover the synchronous case information, namely,p = p and p c = p ⊗ p.
The entries of the Perron vectors {p c ,p} are related to each other, as illustrated by the following explanation (proofs appear in [3][App. VII]). Since the vector p c is of dimension N 2 × 1, we partition it into N sub-vectors of dimension N × 1 each: where p k denotes the k-th sub-vector. We construct an N × N matrix P c from these sub-vectors: We use p ℓk to denote the (ℓ, k)-th element of matrix P c , which is equal to the ℓ-th element of p k . It can be verified that the matrix P c in (103b) is symmetric positive semi-definite and it satisfies wherep is the Perron eigenvector in (102b). We can further establish the following useful relations: We can also verify that the matrix difference is symmetric, positive semi-definite, and satisfies C c 1 = 0. Moreover, it is straightforward to verify that where c c,kk denotes the (k, k)-th entry in C c .

VI. ASYNCHRONOUS NETWORK PERFORMANCE
We now comment on the performance of asynchronous networks and compare their metrics against both noncooperative and centralized strategies.

A. MSD Performance
We denote the MSD performance of the individual agents and the average MSD performance across the network by: where the error vectors are measured relative to the global optimizer, w o . We further denote the gradient noise process at the individual agents by and define where F i−1 now represents the collection of all random events generated by the iterates from across all agents, {w k,j , k = 1, 2, . . . , N }, up to time i − 1: It is stated later in (149a), under some assumptions on the gradient noise processes (106a), that for stronglyconnected asynchronous diffusion networks of the form (86), and for sufficiently small step-sizes µ x [3], [5]: The same result holds for asynchronous consensus and CTA diffusion strategies. Observe from (107) the interesting conclusion that the distributed strategy is able to equalize the MSD performance across all agents for sufficiently small step-sizes. It is also instructive to compare expression (149a) with (63) and (67) in the centralized case.
Observe how cooperation among agents leads to the appearance of the scaling coefficients {p k , p c,kk }; these factors are determined byĀ and C A . Note further that if we setμ k = µ, σ 2 µ,k = 0,p k = p k and p c,kk = p 2 k , then we recover the MSD expression for synchronous distributed strategies: This result agrees with expression (62) from [1].
Example VI.1 (MSE networks with random updates). We continue with the setting of Example IV.1, which deals with MSE networks. We assume the first and second-order moments of the random step-sizes are uniform, i.e.,μ k ≡μ and c µ,kk ≡ σ 2 µ , and also assume uniform regression covariance matrices, i.e., R u,k ≡ R u for k = 1, 2, . . . , N . It follows that H k = 2R u ≡ H and R s,k = 4σ 2 v,k R u . Substituting into (149a), and assuming a fixed topology with fixed combination coefficients set to a ℓk , we conclude that the MSD performance of the diffusion strategy (91a) with random updates is well approximated by: where and p k is the k−th entry of the Perron vector p defined by (100). If the combination matrix A happens to be doubly stochastic, then its Perron eigenvector becomes p = 1/N . Substituting which agrees with the centralized performance (51b). In other words, the asynchronous diffusion strategy is able to match the performance of the centralized solution for doubly stochastic combination policies, when both implementations employ random updates. Since the centralized solution can improve the average MSD performance over non-cooperative networks, we further conclude that the diffusion strategy can also exceed the average performance of non-cooperative networks.
Example VI.2 (Asynchronous MSE networks). We continue with the setting of Example VI.1 except that we now employ the asynchronous LMS diffusion network (91b). Its MSD performance is well approximated by: If the mean combination matrixĀ happens to be doubly stochastic, then its Perron eigenvector becomesp = 1/N . Substitutinḡ p k = 1/N into (109a), and using p c,kk =p 2 k + c c,kk , where c c,kk is from (103e), gives It is clear that if c c,kk = σ 2 π,k , then the MSD performance in (110b) will agree with the centralized performance (55a). In other words, the distributed diffusion strategy is able to match the performance of the centralized solution.

VII. NETWORK STABILITY AND PERFORMANCE
In this section, we examine more closely the performance and stability results that were alluded to in the earlier sections. We first examine the consensus and diffusion strategies in a unified manner, and subsequently focus on diffusion strategies due to their enhanced stability properties, as the ensuing discussion will reveal.

A. MSE Networks
We motivate the discussion by presenting first some illustrative examples with MSE networks, which involve quadratic costs. Following the examples, we extend the framework to more general costs.
Example VII.1 (Error dynamics over MSE networks). We consider the MSE network of Example IV.1, which involves quadratic costs with a common minimizer, w o . The update equations for the non-cooperative, consensus, and diffusion strategies are given by (68a), (81), and (83a)-(83b). We can group these strategies into a single unifying description by considering the following structure in terms of three sets of combination coefficients {a o,ℓk (i), a 1,ℓk (i), a 2,ℓk (i)}: In (111) CTA diffusion: where A i denotes some generic combination policy satisfying (89). We associate with each agent k the following three errors: which measure the deviations from the global minimizer, w o . Subtracting w o from both sides of the equations in (111) we get In a manner similar to (16a), the gradient noise process at each agent k is given by In order to examine the evolution of the error dynamics across the network, we collect the error vectors from all agents into N × 1 block error vectors (whose individual entries are of size M × 1 each): Motivated by the last term in the second equation in (114a), and by the gradient noise terms (114b), we also introduce the following N × 1 column vectors whose entries are of size M × 1 each: We further introduce the Kronecker products and the following N × N block diagonal matrices, whose individual entries are of size M × M each: From (114a) we can then easily conclude that the block network variables satisfy the recursions: so that the network weight error vector, w i , evolves according to: For comparison purposes, if each agent operates individually and uses the non-cooperative strategy (68a), then the weight error vector would instead evolve according to the following recursion: where the matrices {A o,i , A 1,i , A 2,i } do not appear any longer, and with a block diagonal coefficient matrix It is also straightforward to verify that recursion (117b) can be equivalently rewritten in the following form in terms of the gradient noise vector, s i , defined by (115b): Example VII.2 (Mean-error behavior). We continue with the setting of Example VII.1. In mean-square-error analysis, we are interested in examining how the quantities E w i and E w i 2 evolve over time. If we refer back to the data model described in Example IV.1, where the regression data {u k,i } were assumed to be temporally white and independent over space, then the stochastic matrix R i appearing in (117b)-(118) is seen to be statistically independent of w i−1 . We further assume that, in the unified formulation, the entries of the combination policies {A o,i , A 1,i , A 2,i } are independent of each other (as well as over time) and of any other variable in the learning algorithm. Therefore, taking expectations of both sides of these recursions, and invoking the fact that u k,i and v k (i) are also independent of each other and have zero means (so that E z i = 0), we conclude that the mean-error vectors evolve according to the following recursions: The matrixB controls the dynamics of the mean weight-error vector for the distributed strategies. Observe, in particular, that B reduces to the following forms for the various strategies (non-cooperative (68a), consensus (81), and diffusion (83a)-(83b)): Example VII.3 (MSE networks with uniform agents). The results of Example VII.2 simplify when all agents employ step-sizes with the same mean value,μ k ≡μ, and observe regression data with the same covariance matrix, R u,k ≡ R u [7], [80]. In this case, we can expressM and R from (152b) and (119c) in Kronecker product form as follows: so that expressions (122a)-(122d) reduce toB Observe thatB atc =B cta , so we denote these matrices byB diff . Using properties of the eigenvalues of Kronecker products of matrices, it can be easily verified that the M N eigenvalues of the aboveB matrices are given by the following expressions in terms of the eigenvalues of the component matrices {Ā, R u } for k = 1, 2, . . . N and m = 1, 2, . . . , M : Example VII.4 (Potential instability in consensus networks). Consensus strategies can become unstable when used for adaptation purposes [5], [80]. This undesirable effect is already reflected in expressions (125a)-(125c). In particular, observe that the eigenvalues ofĀ appear multiplying (1 − 2µλ m (R u )) in expression (125c) for diffusion. As such, and since ρ(Ā) = 1 for any left-stochastic matrix, we conclude for this case of uniform agents that ρ(B diff ) = ρ(B ncop ). It follows that, regardless of the choice of the mean combination policyĀ, the diffusion strategies will be stable in the mean (i.e., E w i will converge asymptotically to zero) whenever the individual non-cooperative agents are stable in the mean: individual agents stable =⇒ diffusion networks stable (126a) The same conclusion is not true for consensus networks; the individual agents can be stable and yet the consensus network can become unstable. This is because λ k (Ā) appears as an additive (rather than multiplicative) term in (125b) (see [5], [6], [80] for examples): individual agents stable consensus networks stable (126b) The fact that the combination matrixĀ T appears in an additive form in (122b) is the result of the asymmetry that was mentioned earlier in the update equation for the consensus strategy. In contrast, the update equations for the diffusion strategies lead tō A T appearing in a multiplicative form in (122c)-(122d).
Example VII.5 (Useful stability result). It is observed from expressions (122c)-(122d) in the asynchronous case, as well as from the corresponding expressions (81c)-(81d) in the synchronous case studied in [1], that the mean stability of diffusion strategies usually involves examining the stability of a matrix product of the form: where D is a block diagonal symmetric matrix with blocks of size M × M , while A 1 and A 2 are Kronecker product matrices defined in terms of N × N left-stochastic matrices A 1 and A 2 as Matrix products of the form (127) are induced by the cooperation mechanism that is inherent to diffusion learning. They have a useful property: it turns out that these matrix products are stable, regardless of A 1 and A 2 , as long as D is stable (i.e., has all its eigenvalues strictly inside the unit disc). This useful result is easy to establish for symmetric left-stochastic matrices A 1 and A 2 , as already noted in [78]. This is because for symmetric matrices, their spectral radii coincide with their 2−induced norms and, hence, Consequently, since we already know that ρ(A 1 ) = ρ(A 2 ) = 1, it follows that which confirms that a stable D guarantees a stable B, regardless of A 1 and A 2 . The conclusion that B in (127) is stable whenever D is stable continues to hold even when the matrices A 1 and A 2 are not necessarily symmetric. However, the argument leading to (129) will need to be adjusted because property (128) need not hold anymore. This more general result was established in [7, App. D] and also in [105, App. A, pp. 3471-3473], where it was shown that multiplication of a symmetric block diagonal matrix D by any (not necessarily symmetric) left-stochastic Kronecker-product transformations from left and right generally reduces the spectral radius, i.e., Accordingly, a stable D again ensures a stable B. This conclusion was established in the above references by replacing the 2−induced norm used to arrive at (129) by a more convenient block-maximum norm, denoted by · b,∞ and defined as follows. Let x = col{x 1 , x 2 , . . . , x N } denote an N × 1 block column vector whose individual entries are themselves vectors of size M × 1 each. Following [7], [36], [106], the block maximum norm of x is denoted by x b,∞ and is defined as That is, x b,∞ is equal to the largest Euclidean norm of its block components (this definition extends the regular notion of the ∞−norm of a vector to block vectors). The vector norm (131) induces a block-maximum matrix norm. Let A denote an arbitrary N × N block matrix with individual block entries of size M × M each. Then, the block-maximum norm of A is defined as The block-maximum norm has several useful properties -see [7]. In particular, when A is N × N left-stochastic and A = A ⊗ I M , then it can be verified that A T b,∞ = 1. Likewise, when D is block diagonal and symmetric, then D b,∞ = ρ(D). Consequently, repeating the argument leading to (129) and replacing the 2−induced norm used there by the block-maximum norm we have and we again conclude that a stable symmetric D guarantees a stable B for general left-stochastic matrices A 1 and A 2 since, for symmetric D, it holds that ρ(D) = D b,∞ .
Remark 1. The same argument (133) can be used to relax the requirement of symmetry and stability on the block diagonal matrix D. Actually, as long as D b,∞ < 1, which is guaranteed by requiring the block diagonal entries of D to have their 2−induced norms bounded by one, we can again conclude that ρ(B) < 1 so that B is stable.
Remark 2. The validity of property (130) for general left-stochastic matrices was already noted and exploited in earlier works, e.g., in Lemma 1 of [79] and in Lemma 2 of [125]. However, the statement of Lemma 1 in [79] left out the qualification "diagonal" for the center matrix, and the norm · ρ that was used in the proof of the lemma in Appendix I of [79] should be replaced by the · ∞ −norm. Although these corrections were already noted in the references [7], [105], we restate below the correct form of Lemma 1 from [79] for accuracy, and provide its adjusted proof using the current notation. That lemma deals with left-stochastic matrices A 1 and A 2 prior to extension by Kronecker products. Its correct statement should read as follows (the word "diagonal" is missing from the statement in [79]).
Restatement of Lemma 1 [79]: Let A 1 , A 2 , and D denote arbitrary N × N matrices, where A 1 and A 2 have real non-negative entries, with columns adding up to one, i.e., 1 T A 1 = 1 T , 1 T A 2 = 1 T . Then, the matrix B = A T 2 DA T 1 is stable for any choice of A 1 and A 2 if, and only if, the diagonal matrix D is stable.
Proof: One immediate derivation is to employ the ∞−norm and to note that and, hence, It follows that a stable D guarantees a stable B. A second derivation that fixes the argument from Appendix I of [79] relies on using the ∞−norm instead of the ρ−norm used there. Thus, note that we can alternatively argue that when D is stable. This completes the proof of Lemma 1 from [79].
Remark 3. These same arguments establish the validity of Lemma 2 from [125], which deals with a special case involving a matrix of the form B = A T 2 D using the current notation (the matrix in [125] is denoted by F = C T M with B ← F , A 2 ← C, and D ← M; moreover, the matrix M is now generally non-symmetric but still block diagonal with stable entries of the form M kk = (I − P k S k )F for some matrices {P k , S k , F } defined in [125]. When M kk 2 < 1, these block diagonal entries will have 2−induced norms smaller than one so that D b,∞ < 1. Again, the ρ−norm used in the proof of Lemma 2 in [125] should be replaced by the · b,∞ norm and the argument adjusted as in (133) or, alternatively, we note that networks, which deal with mean-square-error cost functions. We now consider networks with more general costs, {J k (w)}, and apply diffusion strategies to seek the global minimizer, w o , of the aggregate cost function, J glob (w), defined by (56). Without loss in generality, we consider ATC diffusion implementations of the form (88): Similar conclusions will apply to CTA diffusion implementations.
Assumption VII.1 (Conditions on cost functions). The aggregate cost J glob (w) in (56) is twice differentiable and satisfies a condition similar to (9) in Assumption II.1 for some positive parameters ν d ≤ δ d . Moreover, all individual costs {J k (w)} are assumed to be strongly-convex with their global minimizers located at w o , as indicated earlier by (73).
As explained before following (73), references [1], [5], [42] present results on the case in which the individual costs are only convex and need not be strongly convex. These references also discuss the case in which the individual costs need not share minimizers. With each agent k in (88), we again associate a gradient noise vector: Assumption VII.2 (Conditions on gradient noise). It is assumed that the first and second-order conditional moments of the gradient noise components satisfy: almost surely for some nonnegative scalars β 2 k and σ 2 s,k , and where F i−1 represents the collection of all random events generated by the iterates from across all agents, {w ℓ,j , ℓ = 1, 2, . . . , N }, up to time i − 1. Moreover, it is assumed that the limiting covariance matrix of s k,i (w o ) exists: We collect the error vectors and gradient noises from across all agents into N × 1 block vectors, whose individual entries are of size M × 1 each: and where we are dropping the argument w k,i−1 from the s k,i (·) for compactness of notation. Likewise, we introduce the following N × N block diagonal matrices, whose individual entries are of size M × M each: where Now, in a manner similar to (29b), we can appeal to the mean-value theorem [5], [35], [49] to note that (143) so that the approximate gradient vector can be expressed as: Subtracting w o from both sides of (138a)-(138b), and using (144), we find that the network error vector evolves according to the following stochastic recursion: Recursion (145a) describes the evolution of the network error vector for general convex costs, J k (w), in a manner similar to recursion (119a) in the mean-square-error case. However, recursion (145a) is more challenging to deal with because of the presence of the random matrix H i−1 ; this matrix is replaced by the constant term R in the earlier recursion (119a) because that example deals with MSE networks where the individual costs, J k (w), are quadratic in w and, therefore, their Hessian matrices are constant and independent of w. In that case, each matrix H k,i−1 in (142c) will evaluate to 2R u,k and the matrix H i−1 in (145b) will coincide with the matrix R defined by (119c). The next statement ascertains that sufficiently small step-sizes exist that guarantee the MSE stability of the asynchronous diffusion strategy (138a)-(138b) [2], [5].
Result (147) shows that the MSD of the network is in the order of µ x,max . Therefore, sufficiently small stepsizes lead to sufficiently small MSDs. As was the case with the discussion in subsection II-G, we can also seek a closed-form expression for the MSD performance of the asynchronous diffusion network and its agents. To do that, we first introduce the analog of Assumption II.4 for the network case.

Assumption VII.3 (Smoothness Conditions).
The Hessian matrix of the individual cost functions {J k (w)}, and the noise covariance matrices defined for each agent in a manner similar to (17a) and denoted by R s,k,i (w), are assumed to be locally Lipschitz continuous in a small neighborhood around w o : for small perturbation δw ≤ r d and for some τ k,d , τ k,s ≥ 0 and 1 ≤ κ ≤ 2.
Lemma VII. where H k = ∇ 2 w J k (w o ). Moreover, for large enough i, the convergence rate towards the above steady-state value is well approximated by the scalar: Proof: See App. XII in [3].
Example VII.6 (Gaussian regression data). MSE performance expressions of the form (149a) are accurate to first-order in the step-size parameters, i.e., they are on the order of O(μ k ). The same is true of expression (112) from [1] in the synchronous case. There are situations, however, where exact expressions for the MSE performance can be derived for multi-agent networks. We illustrate this possibility here for the case of MSE networks of the type described earlier in Example IV.1. We consider first synchronous networks and comment later on how the results should be adjusted to handle asynchronous behavior. We refer to the same setting of Example 9 from [1] where we have N agents observing streaming data {d k (i), u k,i } that satisfy the regression model: We assume the regression vectors are zero-mean Gaussian-distributed with diagonal covariance matrices denoted by Λ k , say, We further assume that the regression data is temporally white and independent over space so that We also assume that the measurement noise process v k (i) is temporally white and independent over space so that E v k (i)v ℓ (j) = σ 2 v,k δ k,ℓ δ i,j in terms of the Kronecker delta sequence δ m,n . Likewise, the measurement noise v k (i) and the regression data u ℓ,j are assumed to be independent of each other for all k, ℓ, i, j. The Gaussian assumption on the regression data is useful because, for this case, fourth-order moments of the regression vectors can be evaluated in closed-form (and such fourth-order moment calculations will arise in the process of determining closed-form expressions for the MSE). Specifically, for any matrix Q ≥ 0 of size M × M , it holds for independent real-valued Gaussian regressors that [28] (there is a typo in reproducing this expression in the second equation of App. II in [79]): E u T k,i u k,i Qu T ℓ,i u ℓ,i = R u,k QR u,ℓ + δ k,ℓ {R u,k Tr (QR u,k ) + R u,k QR u,k } (150d) We assume the network is running a diffusion strategy of the form: which includes both the ATC and CTA LMS diffusion forms as special cases. We know from (77b) in [1] that the error network vector evolves according to the dynamics which is defined in terms of the following N × 1 column vectors whose entries are of size M × 1 each: are N × N block diagonal matrices, whose individual entries are of size M × M each. Let Σ be any N × N non-negative definite block matrix that we are free to choose, with blocks of size M × M . Computing the Σ−weigthed squared norm of the error vector in (151b) under expectation gives (see the derivation leading to (269) in [7] or (38)- (39) in [79]): where the deterministic weighting matrix Σ ′ is given by: We can evaluate the last expectations in (153a)-(153b) in closed-form. But first we need to introduce a convenient blockvector notation, denoted by bvec(·). Thus, given an N × N block matrix, with blocks of size M × M each, say for N = 3, its block vectorization is obtained as follows. We first vectorize each of the block entries and define the column vector where the notation A ⊗ b B denotes the block Kronecker product of two block matrices A and B (assumed here to be both of size N × N with M × M blocks); the kℓ−th block of A ⊗ b B has size N M 2 × N M 2 and is given by [153]: in terms of the traditional Kronecker product operation. Using the block Kronecker properties (154c)-(154d) we now find that the last expectation in (153a) is given by: This is the same expression (56) from [78] for the case of CTA diffusion, and the same expression (42) from [79]) for CTA and ATC diffusion (except that this latter reference used the traditional vec(.) and Kronecker notation ⊗ instead of bvec(.) and the block Kronecker notation ⊗ b ).
Let us now evaluate the last expectation in (153a). Let Q = MA 2 ΣA T 2 M so that we can rewrite more compactly: The kℓ−th block entry in this matrix is given by in terms of the kℓ−th block of Q. Using property (150d) for Gaussian regressors, we get K kℓ = 4R u,k Q kℓ R u,ℓ + 4δ kℓ {R u,k Tr [Q kk R u,k ] + R u,k Q kk R u,k } (156c) It is clear from the above expression that the matrix K has the following general form involving two block diagonal matrices: Introduce the block diagonal matrices (written for N = 3): where λ k = vec(Λ k ). Then, it can be verified that bvec(Z 1 ) = L 1 (R ⊗ b R)bvec(Q) (158a) bvec(Z 2 ) = L 2 bvec(Q) (158b) Noting that bvec(Q) we conclude that Note that the matrix X has the following block-diagonal structure: k , . . . , X which is the same structure derived through equations (57)- (67) in [78]. Substituting (160a) into (153b) and using again (154c)-(154d) we find that the block vectorized forms of the weighting matrices {Σ, Σ ′ } are related via (the expression for F below fixes the typo in Eq. (44) from [79]; in the derivation in this example we considered the special case in which S m = I in (44) from [79]): This is the same expression for F in Eq. (69) from [78] for the case of CTA diffusion (where A 2 = I N ). This is also the same expression for F in Eq. (41) from [79] for CTA and ATC diffusion (except that this latter reference wrote ⊗ instead ⊗ b ). Substituting (155a) and the above expression for σ ′ into (153a), and using the compact notation x 2 σ for x 2 Σ , we rewrite (153a) in the form In seady-state, as i → ∞, the mean-square-error approaches As explained in Sec. 6.6 of [7], the network MSD can be assessed by selecting σ to satisfy which is expression (105a) from [78].
To illustrate how these results are adjusted for asynchronous behavior, we consider the case in which the step-size parameters {µ k } in (151a) are replaced by random values {µ k (i)}. We denote the mean and variances of these random variables as follows: the constant step-size case in order to examine solutions that are able to adapt and learn continuously from streaming data. The presentation complements the results from [1], [5]. We introduced a fairly general model for asynchronous behavior that allows for random step-sizes, link failures, random topology variations, and random combination coefficients. We examined the mean-square-error performance and stability properties under asynchronous events and recovered results for synchronous operation as a special case. The results indicate that asynchronous networks are robust, resilient to failure, and remain mean-square stable for sufficiently small step-sizes.