{'fraction_non_alphanumeric': 0.0588377968671046, 'fraction_numerical': 0.014815563415866598, 'mean_word_length': 3.7863016349037437, 'pattern_counts': {'":': 0, '<': 21, '<?xml version=': 0, '>': 19, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 19, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Convergence of a Multi-Agent Projected Stochastic Gradient Algorithm for Non-Convex Optimization

We introduce a new framework for the convergence analysis of a class of distributed constrained non-convex optimization algorithms in multi-agent systems. The aim is to search for local minimizers of a non-convex objective function which is supposed to be a sum of local utility functions of the agents. The algorithm under study consists of two steps: a local stochastic gradient descent at each agent and a gossip step that drives the network of agents to a consensus. Under the assumption of decreasing stepsize, it is proved that consensus is asymptotically achieved in the network and that the algorithm converges to the set of Karush-Kuhn-Tucker points. As an important feature, the algorithm does not require the double-stochasticity of the gossip matrices. It is in particular suitable for use in a natural broadcast scenario for which no feedback messages between agents are required. It is proved that our result also holds if the number of communications in the network per unit of time vanishes at moderate speed as time increases, allowing for potential savings of the network's energy. Applications to power allocation in wireless ad-hoc networks are discussed. Finally, we provide numerical results which sustain our claims.


I. INTRODUCTION
Stochastic gradient descent is a widely used procedure for finding critical points of an unknown function f [33]. Formally, it can be summarized as an iterative scheme of the form θ n+1 = θ n + γ n+1 (−∇f (θ n ) + ξ n+1 ) where ∇ is the gradient operator and where ξ n+1 represents a random perturbation. Relevant selection of the step size γ n ensures that, for a well behaved function f , sequence (θ n ) n≥0 will eventually converge to a critical point.
In this paper, we investigate a distributed optimization problem which is of practical interest in many multi-agent contexts such as parallel computing [8], statistical estimation [36], [35], [4], [29], robotics [13] or wireless networks [30]. Consider a network of N agents. To each agent i = 1, . . . , N , we associate a possibly non-convex continuously differentiable utility function f i : R d → R where d ≥ 1. Let G ⊂ R d be a nonempty compact convex subset. We address the the following optimization problem: The set G is assumed to be known by all agents. However, a given agent i does not know the utility functions f j 's of other agents j = i. Cooperation between agents is therefore needed to find minimizers of (1). Moreover, any utility function f i may be unperfectly observed by agent i itself, due to the presence of random observation noise. We thus address the framework of distributed stochastic approximation.
The literature contains at least two different cooperation approaches for solving (1). The so-called incremental approach is used by [27], [24], [28], [18], [31]: a message containing an estimate of the desired minimizer iteratively travels all over the network. At any instant, the agent which is in possession of the message updates its own estimate and adds its own contribution, based on its local observation. Incremental algorithms generally require the message to go through a Hamiltonian cycle in the network. Finding such a path is known to be a NP complete problem. Relaxations of the Hamiltonian cycle requirement have been proposed: for instance, [24] only requires that an agent communicates with another agent randomly selected in the network (not necessarily in its neighborhood) according to the uniform distribution. However, substantial routing is still needed. In [21], problem (1) is solved using a different approach, assuming that agents perfectly observe their utility functions and know also the utility functions of their neighbors. This paper focuses on another cooperation approach based on average consensus techniques, see references [15], [40]. In this context, each agent maintains its own estimate. Agents separately run local gradient algorithms and simultaneously communicate in order to eventually reach an agreement over the whole network on the value of the minimizer. Communicating agents combine their local estimates in a linear fashion: a receiver computes a weighted average between its own estimate and the ones which have been transmitted by its neighbors. Such combining techniques are often refered to as gossip methods.
The idea underlying the algorithm of interest in this paper is not new. Its roots can be found in [40], [41] where a network of processors seeks to optimize some objective function known by all agents (possibly up to some additive noise). More recently, numerous works extended this kind of algorithm to more involved multi-agent scenarios, see [25], [26], [32], [23], [19], [38] for a non exhaustive list. Multi-agent systems are indeed more difficult to deal with, because individual agents does not know the global objective function to be minimized. Reference [25] addresses the problem of unconstrained optimization, assuming convex but not necessarily differentiable utility functions. Convergence to a global minimizer is established assuming that utility functions have bounded (sub)gradients. Let us also mention [38] which focuses on the case of quadratic objective functions. Unconstrained optimization is also investigated in [9] assuming differentiable but non neces-sarily convex utility functions and relaxing boundedness conditions on the gradients. Convergence to a critical point of the objective function is proved and the asymptotic performance is evaluated under the form of a central limit theorem. In [26], the problem of constrained distributed optimization is addressed. Convergence to an optimal consensus is proved when each utility function f i is assumed convex and perfectly known by agent i. These results are extended in [32] to the stochastic descent case i.e., when the observation of utility functions is perturbed by a random noise.
In each of these works, the gossip communication scheme can be represented by a sequence of matrices (W n ) n≥1 of size N × N , where the (i, j)th component of W n is the weight given by agent i to the message received from j at time n, and is equal to zero in case agent i receives no message from j. In most works (see for instance [25], [26], [32], [9]), matrices W n are assumed doubly stochastic, meaning that W T n 1 = W n 1 = 1 where 1 the N × 1 vector whose components are all equal to one and where T denotes transposition. Although row-stochasticity (W n 1 = 1) is rather easy to ensure in practice, column-stochasticity (W T n 1 = 1) implies more stringent restrictions on the communication protocol. For instance, in [12], each one-way transmission from an agent i to another agent j requires at the same time a feedback link from j to i. Double stochasticity prevents one from using natural broadcast schemes, in which a given agent may transmit its local estimate to all its neighbors without expecting any immediate feedback [3]. Very recently, [23] made a major step forward, getting rid of the column stochasticity condition, and thus opening the road to broadcastbased constrained distributed optimization algorithms. It is worth noting however that the algorithm of [23] is such that only receiving agents update their estimates. Otherwise stated, an agent deletes its local observations as long as it is not the recipient of a message. Moreover, except perhaps in some special network topologies, the algorithm of [23] strongly relies on a specific choice of the stepsize. In particular, a necessary condition for the convergence to the desired consensus is that the stepsize vanishes at speed 1/n. However, in practice, it is often desirable to have a leeway on the choice of the stepsize to avoid slow convergence issues.

Contributions
In this paper, we address the optimization problem (1) using a distributed projected stochastic gradient algorithm involving random gossip between agents and decreasing stepsize.
• Unlike previous works, utility functions are allowed to be non-convex. We introduce a new framework for the analysis of a general class of distributed optimization algorithm, which does not rely on convexity properties of the utility functions. Instead, our approach relies on recent results of reference [7] about perturbed differential inclusions. Under a set of assumptions made clear in the next section, we establish that, almost surely, the sequence of estimates of any agent shadows the behavior of a differential variational inequality, and eventually converges to the set of Karush-Kuhn-Tucker (KKT) points of (1).
• Our assumptions encompass the case of non-doubly stochastic matrices W n and, as a particular case, the natural broadcast scheme of reference [3]. Our proofs reveal that, loosely speaking, the relaxation of column stochasticity brings a "noise-like" term in the algorithm dynamics, but which is not powerful enough to prevent convergence to the KKT points. • We show that our convergence result still holds in case the number of communications in the network per unit of time vanishes at moderate speed as time increases.
As an illustration, we apply our results to the problem of power allocation in the wireless interference channel.
The paper is organized as follows. Section II introduces the distributed algorithm and the main assumptions on the network and the observation model. The main result is stated in Section III. Section IV is devoted to its proof. We discuss applications to power allocation in Section V. Section VI describes some standard communication schemes in more details, and provides numerical results.

A. Description of the Algorithm
After an [Initialization step] where each agent starts at a given θ 0,i , each node i generates a stochastic process (θ n,i ) n≥1 in R d using a two-step iterative algorithm: [Local step] Node i generates at time n a temporary estimateθ n,i given bỹ where γ n is a deterministic positive step size, Y n,i is a random variable, and P G represents the projection operator onto the set G. Random variable Y n,i is to be interpreted as a perturbed version of the negative gradient of f i at point θ n−1,i . As will be made clear by Assumption 1(e) below, it is convenient to think of Y n,i as Y n,i = −∇f i (θ n−1,i ) + δM n,i where δM n,i is a martingale difference noise which stands for the random perturbation.
[Gossip step] Node i is able to observe the valuesθ n,j of some other j's and computes the weighted average: where for any i, N j=1 w n (i, j) = 1. In the sequel, we define the N × N matrix W n := [w n (i, j)] i,j=1···N .
Define the random vectors θ n and Y n as θ n := (θ T n,1 , . . . , θ T n,N ) T and Y n = (Y T n,1 , . . . , Y T n,N ) T . The algorithm reduces to: where ⊗ denotes the Kronecker product, I d is the d×d identity matrix and P G N is the projector onto the N th order product set G N := G × · · · × G.

B. Observation and Network Models
Random processes (Y n , W n ) n≥1 are defined on a measurable space equipped with a probability P. Notation E represents the corresponding expectation. For any n ≥ 1, we introduce the σ-field F n = σ(θ 0 , Y 1 , . . . , Y n , W 1 , . . . , W n ). The distribution of the random vector Y n+1 conditioned on F n is assumed to be such that: for any measurable set A, where (µ θ ) θ∈R dN is a given family of probability measures on R dN . We denote by |x| the Euclidean norm of a vector x. We denote by E c the complementary set of any set E. Notation x ∨ y stands for max(x, y). Assumption 1: The following conditions hold: a) (W n ) n≥1 is a sequence of matrix-valued random variables with non-negative components s.t.
• W n is row stochastic for any n: W n 1 = 1, • E(W n ) is column stochastic for any n: c) Conditionaly to F n , W n+1 and Y n+1 are independent. Moreover, the probability distribution of Y n+1 conditionaly to F n depends on θ n only. It is denoted µ θn . d) For any i = 1, . . . , N , f i is continuously differentiable. e) For any n ≥ 1, We now discuss the above Assumption. Conditions 1(a) and 1(b) summarize our assumptions on matrices W n that is, on the communication scheme used in the network. Following the work of reference [12], random gossip is assumed in this paper. Each matrix W n must be row stochastic, this means that each agent i = 1, · · · , N must compute a weighted average j w n (i, j) = 1. Note that a quite classical condition in the literature is to further assume that W n is column-stochastic for any n [25], [26], [32], [9]. Column stochasticity inevitably goes with some restrictions on the communication protocol as discussed in Section I. Here, our assumption is weaker. We only require that W n is column stochastic on average. This is for instance the case in the natural broadcast scheme of reference [3] which will be discussed in the Section II-C. Assumption 1(b) is a connectivity condition of the underlying network graph which will be discussed in more details in Section II-C.
Assumptions 1(c-e) are related to the observation model. Assumption 1(c) states three different things. First, random variables W n+1 and Y n+1 are independent conditionally on the past. Second, (W n ) n≥1 form an independent sequence (not necessarily identically distributed). Finally, the distribution of Y n+1 conditionally on the past is as in (4). Assumption 1(e) means that each Y n,i can be interpreted as a noisy version of −∇f i (θ n−1,i ). The distribution of the random additive perturbation Y n,i − (−∇f i (θ n−1,i )) is likely to depend on the past through the value of θ n−1 , but has a zero mean for any given value of θ n−1 . By Markov's inequality, Assumption 1(f) implies that (µ θ ) θ∈G N is tight i.e., for any > 0, there exists a compact set K ⊂ R dN such that sup θ∈G N µ θ (K c ) < .
Assumption 2: a) The deterministic sequence (γ n ) n≥1 is positive and such that n γ n = ∞. b) There exists α > 1/2 such that: Note that, (6) implies (but is not equivalent to) n γ 2 n < ∞, which is a rather usual assumption in the framework of decreasing step size stochastic algorithms [20]. In order to have some insights on (7), first consider the case where the matrices (W n ) n≥1 form an i.i.d. sequence i.e., the spectral radius ρ := ρ n does not depend on n. Then both conditions (5) and (7) are satisfied if and only if: Nevertheless, matrices (W n ) n≥1 do not need to be i.i.d. An interesting example is when matrix W n is likely to be equal to identity with a probability that tends to one as n → ∞. From a communication point of view, this means that the exchange of information between agents becomes rare as n → ∞. This context is especially interesting in case of wireless networks, where it is often required to limit as much as possible the communication overhead. Let us emphasize that if W n is assumed i.i.d with ρ < 1, then the usual assumptions γ n = +∞ and γ 2 n < ∞ can be substituted to Assumption 2. Assumption 2 is designed to take into account the non-stationnarity of W n .

C. Illustration: Some Examples of Gossip schemes
Here, we focus on two standard communication schemes and give the corresponding sequence (W n ) n≥1 for each of them. We refer the reader to reference [14] for a more complete picture and for more general gossip strategies. We introduce what we shall refer to as the pairwise and the broadcast schemes. The first one can be found in the paper of Boyd et al. [12] on average consensus while the second is inspired from the broadcast scheme depicted in [3]. The network of agents is represented as a nondirected graph (V, E) where V is the set of N nodes and E corresponds to the set edges between nodes.
1) Pairwise Gossip: At time n, a single node i wakes up (node i is chosen at random, uniformly within the set of nodes and independently from the past). Node i randomly selects a node j among its neighbors in the graph. Node i and j exchange their temporary estimatesθ n,i andθ n,j and compute the weighted average θ n,i = θ n,j = βθ n,i + (1 − β)θ n,j where 0 < β < 1. Other nodes k / ∈ {i, j} simply set θ n,k =θ n,k . Set β = 1/2 for simplicity. In this case, the corresponding matrix W n is given by W n = I N − (e i − e j )(e i − e j ) T /2 where e i denotes the ith vector of the canonical basis in R N . Note that for each n, W n forms an i.i.d. sequence of doubly stochastic matrices. Assumption 1(a) is obviously satisfied. Moreover, the spectral radius ρ of matrix E(W T n (I N − 11 T /N )W n ) satisfies (8) if and only if (V, E) is a connected graph (see [12]).
2) Broadcast Gossip: At time n, a random node i wakes up. The latter node is supposed to be chosen at random w.r.t. the uniform distribution on the set of vertices. Node i broadcasts its temporary update to all its neighbors. Any neighbor j, computes the weighted average θ n,j = βθ n,i + (1 − β)θ n,j . On the other hand, any node k which does not belong to the neighborhood N i of i (this includes i itself) simply sets θ n,k =θ n,k . Note that, as opposed to the pairwise scheme, the transmitter node i does not expect any feedback from its neighbors. It is straightforward to show that the (k, )th component of matrix W n corresponding to such a scheme writes: As a matter of fact, the above matrix W n is not doubly stochas- [3]). Thus, the sequence of matrices (W n ) n≥1 satisfies the Assumption 1(a). Once again, straightforward derivations which can be found in [3] show that the spectral radius ρ satisfies (8) if and only if (V, E) is a connected graph.

A. Framework and Assumptions
We study the case where for any i = 1, · · · , N , the set G is determined by a finite set of p inequality constraints (1 ≤ p < ∞): for some functions q 1 , . . . , q p which satisfy the following conditions. For any θ ∈ R d , we denote by A(θ) ⊂ {1, . . . , p} the active set i.e., A(θ) = {j : q j (θ) = 0, θ ∈ G}. Denote by ∂G the boundary of G. Assumption 3: a) The set G defined by (9) is nonempty and compact. b) For any j = 1, · · · , p, q j : R d → R is a convex function, continuously differentiable in a neighborhood of ∂G. c) For any θ ∈ ∂G, (∇q j (θ) : j ∈ A(θ)) is a linearly independent collection of vectors. In other terms, some regularity is imposed on G. Moreover, (9,c) is a simple qualification assumption [34]: note that the same set G can be expressed using different constraints; for instance, one can always duplicate constraints arbitrarily. The qualification assumption says that it is up to the user to remove redundant constraints.

B. Notations
Recall that |x| represents the Euclidean norm of a vector x. Denote by ∇ the gradient operator. We denote by resp. the projector onto the consensus subspace 1 ⊗ θ : θ ∈ R d and the projector onto the orthogonal subspace. For any vector x ∈ R dN , define the vector of R d Note that We denote by x ⊥ := J ⊥ x the projection of x on the orthogonal to the consensus space. Remark that x = 1⊗ x +x ⊥ . In particular, we set θ ⊥,n := J ⊥ θ n where J ⊥ is given by (10) and refer to θ ⊥,n as the disagreement vector. Denote by: f i the average of utility functions. Define the set of KKT points of f on G (also called the set of stationary points) as: where a ∈ A} for any θ ∈ R dN and any set A. We say that a random sequence (x n ) n≥1 converges almost surely (a.s.) to a set A if d(x n , A) converges a.s. to zero as n → ∞. Please note that convergence to A does not imply convergence to some point of A.
Theorem 1 establishes two points. First, a consensus is achieved as n tends to infinity, meaning that max i,j |θ n,i −θ n,j | converges a.s. to zero. Second, the average estimate θ n converges to the set L of KKT points. As a consequence, if L contains only isolated points, sequence θ n converges a.s. to one of these points.
In particular, when f is convex, θ n converges to the set of global solutions to the minimization problem (1). However, as already remarked, our result is more general and does not rely on the convexity of f . If f is not convex, sequence θ n does not necessarily converge to a global solution. Nevertheless, it is well known that the KKT conditions are satisfied by any local minimizer [11].
The condition that f (L) has an empty interior is satisfied in most practical cases. From Sard's theorem, it holds as soon as f is d times continuously differentiable.

IV. PROOF OF THEOREM 1 A. Sketch of Proof
We first provide some insights on the proof (all statements are made rigorous in the next subsections).
The proof is decomposed in the following steps. 1) First, establish convergence to consensus with probability 1 (see Lemma 1). In this step we show that iterating matrices W n yields to consensus whatever the behavior of sequence (θ n ) n≥1 . See subsection IV-C. 2) Show then that θ n is ruled by the following discrete time dynamical system (see Proposition 1 in subsection IV-D): In order to give some insight, assume just for a moment that g is identically 0. In that case, one could write: and view θ n as a noisy discrete approximation of the well studied Ordinary Differential Equation (ODE): whereẋ stands for the derivative of function t → x(t).
See, for instance, [17]. This line of reasoning is the so-called "ODE" method (see, for instance, [10]). If function g was regular enough, one could still use the ODE method and study an ODE of the formẋ = −∇f (x) + h(x) where h is a function derived from g (we skip the details, which are not important here). Unfortunately, in our case, function g is not regular enough. In that case, the standard ODE method fails for two reasons: i) such an ODE might have several solutions, ii) it does no longer ensure that discretized versions of the ODE stay close to suitable solutions of the original ODE.
3) The framework of Differential Inclusions (DI) addresses these two issues. The ODE is replaced by the Dİ (see (14) in what follows), where, at any time instant t, F (x(t)) is a set of vectors instead of the single vector as in an ODE. An important property is that the setvalued function F should now be upper semicontinuous. It remains to prove two assertions to finish the proof: a) Show that the noisy discretized dynamical system (13) asymptotically behaves "the same way" as a solution to the DI. To formalize this, we shall rely on the notion of perturbed differential inclusions of [7]. We prove in Section IV-E that the continuous-time process obtained from a suitable interpolation of (13) is in fact a perturbed solution to the DI. Using the results of [7] which we recall in Section IV-B, the limit set of the interpolated process can be characterized by simply studying the behavior of the solutions x(t) to the DI for large t.
b) The asymptotic behavior of the DI is addressed in Section IV-F where it is shown that f is a Lyapunov function for the set L of KKT points under the dynamics induced by the DI.

B. Preliminaries: Useful Facts about Set-Valued Dynamical Systems
Before providing the details of the proof, we recall some useful facts about perturbed differential inclusions. All definitions and statements made in this paragraph can be found in reference [7]. However, for the sake of readability and completeness, it is worth recalling some facts.
Consider an arbitrary set-valued function F which maps each point θ ∈ R d to a set F (θ) ⊂ R d . Assume that F satisfies the following conditions: Condition 1: The following holds: if it is absolutely continuous and if for any solution x to (14) and for any t > 0, and where the inequality is strict whenever x(0) / ∈ Λ. Finally, a function y : [0, ∞) → R d is called a perturbed solution to (14) if it is absolutely continuous and if there exists a locally integrable function t → U (t) such that: lim t→∞ δ(t) = 0 and such that for almost every t > 0, (15) Remark 1: Note that there is no guarantee that a perturbed solution remains close to a solution from a given time (for instance, y(t) = log(1 + t) is a perturbed solution toẋ = 0). However, loosely speaking, y is a perturbed solution if the behavior of y on the interval [t, t + T ] shadows the differential inclusion for t large enough, for each fixed width T .
The following result, found in reference [7], will be used in our proofs. Denote by S the closure of a set S. Theorem 2 ( [7]): Let V be a Lyapunov function for Λ. Assume that V (Λ) has an empty interior. Let y be a bounded pertubed solution to (14). Then, t≥0 y([t, ∞)) ⊂ Λ .
Proof: The result is a consequence of Theorem 4.2, Theorem 4.3 and Proposition 3.27 in [7].
Before going into the details of the proof of Lemma 1, it is worth noting that |Z n | 2 ≤ |Y n | 2 (just note that θ n−1 = P G N [θ n−1 ] and use the fact that G N is convex). By Assump- we expand the square Euclidean norm of the latter vector: .
Integrate both sides of the above equation w.r.t. the random variable W n : Expand the righthand side and take the expectation. Using the fact that ρ n < 1 for n large enough, . As E[|Z n | 2 ] is uniformly bounded, we obtain from Cauchy-Schwartz's inequality: for some constant C > 0. Let us denote v n := E[|J ⊥ θ n | 2 ]. Usingγ n = √ Cγ n , which also fulfills Assumption 2, we get: Let u n := n 2α v n for some α > 1/2 satisfying (6) and (7). Then, This implies in turn: u n − u n−1 ≤ −a n u n−1 + b n √ u n−1 + c n n αγ n , where b n = (1 + 1 n−1 ) 2α , a n = 1− √ bnρn n αγ n , and c n = n αγ n . A straightforward analysis of function φ n : u → −a n u + b n √ u + c n shows that u > t n implies φ n (u) < 0 where t n := (b n /a n + c n /b n ) 2 . Note that, using Assumption 1(b), a n ∼ 1−ρn n αγ n and using Assumption 2, t n is bounded above, say by a constant K > 0. Moreover, when u ≤ t n , φ n (u) ≤ φ n (b n /2a n ) = c n +b n 2 /4a n . Notice again that φ n (b n /2a n ) is bounded above, say by a constant L > 0. We have proved that if u n−1 ≤ K then u n ≤ K +L and if u n−1 > K, u n ≤ u n−1 . This implies that u n ≤ max(K + L, u 0 ). Hence v n < ∞.
Lemma 1 proves that agents asymptotically reach an agreement on their estimate. Another way to state Lemma 1 is to write that max i,j=1···N |θ n,i − θ n,j | converges a.s. to zero as n tends to infinity. Therefore, the asymptotic analysis of the whole vector θ n now reduces to the study of the average

D. Average Estimate
. Under Assumption 1(c), this means that: Proposition 1: Under Assumptions 1, 2, 3, there exists two stochastic processes (ξ n ) n≥1 , (r n ) n≥1 such that for any n ≥ 1: θ n = θ n−1 −γ n ∇f ( θ n−1 )+γ n g γn (θ n−1 )+γ n ξ n +γ n r n (19) and satisfying w.p.1: The proof is provided in Appendix I. Note that the third term in the righthand side of (19) is zero whenever θ n−1 + γ n Y n lies in G N i.e., when the projector is inoperant. In order to have some insights, assume just for a moment that this holds for any n after a certain value. In this case, equation (19) simply becomes θ n = θ n−1 − γ n ∇f ( θ n−1 ) + γ n ξ n + γ n r n . (21) In this case, by the continuity of ∇f and using the above conditions on the sequences ξ n and r n , the asymptotic behavior of sequence θ n can be directly characterized using classical stochastic approximation results (see [20], [1], [6], [16]). Indeed, a sequence θ n satisfying (21) converges to the set of critical points of f . Nevertheless, the projector P G N is generally active in practice, so that the term g γn (θ n−1 ) may be nonzero infinitely often. This additional term raises at least two problems. First, it depends on the whole vector θ n and not only on the average θ n : equation (19) looks thus nothing like a usual iteration of a stochastic approximation algorithm. Second, g γ (θ) is not a continous function of θ, whereas standard approaches often assume the continuity of the mean field of the stochastic approximation algorithm.

E. Interpolated Process
Define µ := sup θ∈G N |y|dµ θ (y). Define the following set-valued function F on R d which maps any θ to the set: Using the fact that f is continuously differentiable and that G is closed and convex, it can be shown that F satisfies Condition 1 1 . Recall notation F δ (θ) in (15). Consider stochastic processes (ξ n , r n ) n≥1 as in Proposition 1. Proposition 2: Under Assumptions 1, 2, 3, there exists a sequence of random variables (δ n ) n≥1 converging a.s. to zero and an integer n 0 such that for any n ≥ n 0 , Proof: From the definition of g γ , it is straightforward to show that |g γ (θ)| ≤ 2µ. Note that for any θ ∈ G and any y ∈ R d , the vector θ +γy −P G [θ +γy] belongs to the normal cone N G (P G [θ+γy]) at point P G [θ+γy]. Otherwise stated, P G [θ+ γy] − θ − γy can be written as a linear combination of the gradient vectors associated with the active constraints, where the coefficients of the linear combination are nonnegative (see for instance, [34], theorem 6.14). The latter linear combination is moreover unique due to the qualification constraint given by Assumption 3(c). More precisely, if A(P G [θ + γy]) represents the active set at point P G [θ + γy] for any θ, γ, y, there exists a unique collection of nonnegative coefficients (λ j (θ, γ, y) : j ∈ A(P G [θ + γy])) such that: Throughout the paper, we use the convention that λ j (θ, γ, y) = 0 in case j / ∈ A(P G [θ + γy]). The following technical lemma is proved in Appendix II.
The following function φ : [0, +∞) → [0, +∞) will be useful. Define: 1 As a purely technical point, note that the third point in Condition 1 is satified only if |∇f (θ)| increases at most at linear speed when |θ| → ∞, which has of course no reason to be true in general. This is however unimportant, as the values of θ will always be restricted to the bounded set G in the sequel. Moreover, for θ / ∈ G, one can always redefine F (θ) as in (22) but replacing f with ϕ • f where ϕ is a slowly increasing map chosen such that Condition 1 holds for ϕ • f . Since each gradient ∇q j is continuous, it is uniformly continuous on the compact set G. Thus φ(x) tends to zero as x ↓ 0. Loosely speaking, when γ is small and when all θ i 's are close to the average θ , the point P G [θ i + γy i ] is close to θ . In this case, the uniform continuity of ∇q j implies that ∇q j (P G [θ i + γy i ]) ∇q j ( θ ). Lemma 3 below states a somewhat stronger result. For any ≥ 0 and θ ∈ G, define A(θ, ) as the set of constraints which are active at least for some point in an -neighborhood of θ: Lemma 3: Under Assumption 3, there exists a constant C > 0 and a function γ → (γ) on [0, +∞) satisfying lim γ↓0 (γ) = 0 such that the following holds. For any θ ∈ G N and any γ > 0, there exists (α 1 , · · · , α p ) ∈ [0, C] p such that The proof is provided in Appendix III. The sum in the lefthand side of (27) is a (nonnegative) linear combination of the gradient vectors of the constraints at point θ . However, this does not necessarily imply that this term belongs to the normal cone N G ( θ ) because, for a fixed > 0, the set A( θ , ) is in general larger than the active set A( θ ). Nevertheless, the following lemma states that A( θ , ) is no larger than a certain active set A(θ ) for some θ in a neighborhood of θ . |θ − θ | < δ( ) and A(θ, ) ⊂ A(θ ) .
The first three terms of the righthand side of the above inequality converge to zero as γ n → 0, r n → 0 and δ n → 0. The fourth term tends to zero as well because |g γ (θ)| is uniformly bounded in (γ, θ), as remarked in the proof of Proposition 2. Finally, γ n ξ n tends to zero by (20). Thus η(t) tends to zero as t → ∞. This completes the proof of Proposition 3.

F. Study of the Differential Inclusion Dynamics
Proposition 4: [2] Any solution x of the differential inequality dx dt ∈ F (x) with F defined by eq. (22) such that ∀t ∈ R, x(t) ∈ G satisfies for almost all t ∈ R: where P T G (x) stands for the projection onto the tangent cone T G (x) at point x.
Proof: For the sake of completeness, we reproduce here the proof of reference [2, pp. 266]. Let x(t) be a solution of eq. (14); it is differentiable for almost every t by definition. At such t one has dx since x(t) ∈ G for all t. For the same reason (using < 0), one has − dx dt (t) ∈ T G (x(t)). By convexity of G, T G (θ) is the dual cone of N G (θ) for every θ. Hence ∀v ∈ N G (x(t)) one has at the same time dx dt (t), v ≤ 0 and − dx dt (t), v ≤ 0. As a consequence, ∀v ∈ N G (x(t)), dx dt (t), v = 0. Now, considering that x is a solution of eq. (14) with F defined by eq. (22): dx , v = 0 and it is a classical fact from convex analysis (see, for instance, reference [34]) that, N G (x(t)) and T G (x(t)) are dual cones: v = P T (−∇f (x(t))).
The following proposition is straightforward to prove but has an important role.
Proposition 5: Any solution x of eq. (30) admits f as a Lyapunov function for the set L of KKT points defined by eq. (12).
Proof: One has d dt f (x(t)) = ∇f (x(t)), dx dt (t) . From Proposition 4, we deduce: When T and N are dual cones, the decomposition v = This gives the sought result.
We are now in a position to apply Theorem 2 with Lyapunov function f itself and Λ the set of KKT points of our optimization program (see Proposition 5).

V. APPLICATION: POWER ALLOCATION IN AD-HOC WIRELESS NETWORKS
The context of power allocation for wireless networks has recently raised a great deal of attention in the field of distributed optimization and game theory [37], [5], [22]. Application of distributed optimization to power allocation has been previously investigated in [30].
Consider an ad-hoc network composed of N sourcedestination pairs. We focus on the so-called interference channel i.e., the signal received by the destination of a given pair i = 1, . . . , N is corrupted both by an additive white Gaussian noise of variance σ 2 i and by the interference produced by other sources j = i. Denote by p i ≥ 0 the transmission power of source i. The power of the useful received signal at the destination i is proportional to p i A i,i where A i,i represents the channel gain between the ith source and the corresponding destination. On the other hand, the level of interference endured by the destination i is proportional to j =i p j A j,i where A j,i is the (positive) channel gain between source j and destination i. As will be explained below, the reliability of the link between the ith source and its destination is usually expressed as an increasing function of the signal to interferenceplus-noise ratio, defined as A i,i p i /(σ 2 i + j =i A j,i p j ). We assume that there is no Channel State Information at the Transmitter (no CSIT) i.e., all channel gains are unknown at all transmitters. However, we assume that the destination associated with the ith source-destination pair • knows the set of channel gains A i := (A 1,i , · · · , A N,i ) T , • ignores all other channel gains A j for j = i. Figure 1 below illustrates the interference channel with N = 2 transmit-destination pairs. We assume that 0 ≤ p i ≤ P i .

D1
No CSIT D2 Destination 1 observes A 1,1 A 2,1 where P i is the maximum allowed power for user i. Define θ = (p 1 , · · · , p N ) T as the vector of all powers of all users. The aim is to select a relevant value for parameter θ. We assume that destinations are able to communicate according to an underlying connected graph. The proposed algorithm works as follows.
1) In a first step, the set of destination nodes cooperate and jointly search for a relevant global power allocation θ.
The desired vector θ corresponds to a local minimizer of an optimization problem which will be made clear below. 2) Once an agreement is found on the power allocation vector θ, each destination i provides its own source with the corresponding power p i using a dedicated channel. Consider fixed deterministic channels. As a performance metric, consider the error probability observed at each destination. Assuming for instance that each transmitter uses a 4-QAM modulation, the error probability at the ith destination is given by [39,Section 3.1]: where Q(x) = ( √ 2π) −1 ∞ x e −t 2 /2 dt. We investigate the following minimization problem: where β i is an arbitrary positive deterministic weight known only by agent i and where G := {(p 1 , · · · , p N ) ∈ R N : ∀i = 1, · · · , N, 0 ≤ p i ≤ P i }. The above optimization problem is non-convex. Note that, utility functions (31) can of course be replaced by any other continuously differentiable functions of the signal-to-interference-plus-noise ratio without changing the results of this section. Section II suggests the following deterministic distributed gradient algorithm. Each user i has an estimate θ n,i of the whole vector θ at the nth iteration. Here, we stress the fact that a given user has not only an estimate of its own power allocation p i , but has also an estimate of what should be the power allocation of other users j = i. Denote by θ n = (θ T n,1 , · · · , θ T n,N ) T the vector of size N 2 which gathers all local estimates. Denote by A := ((A 1 ) T , · · · , (A N ) T ) T the vector which gathers all N 2 channel gains. The distributed algorithm writes: where for any θ = (θ T 1 , · · · , θ T N ) T in R N 2 , we set and where ∇ θ is the gradient operator with respect to the first argument θ of P e,i (θ, A i ). Corollary 1: Under the stated assumptions on the sequences (γ n ) n≥1 and (W n ) n≥1 , the algorithm (33) is such that sequence ( θ n ) n≥1 converges to the set of KKT points of (32).
Remark 2: In many situations, the channel gains are random and rapidly time-varying. In this case, it is more realistic to assume that each destination i observes a sequence of random i.i.d. channel gains (A i n ) n≥1 . The algorithm (33) can be extended to this context without difficulty. This yields an algorithm which for solving the following optimization problem:

A. Scenario #1
As a benchmark, we address the convex optimization scenario formulated in [23]. Define G ⊂ R 2 as the unit disk in R 2 centered at the origin. Consider the minimization of Here, (R 1 , · · · , R N ) is a collection of i.i.d. real Gaussian distributed random variables with mean 0.5 and unit variance, and (s 1 , · · · , s N ) is a collection of deterministic elements of R 2 . The number of agents is set as N = 10 or N = 50. We used different graphs: the complete graph where any agent is connected to all other agent, and the cycle. We evaluate the performance of both pairwise and broadcast algorithms described in Section II-C. The weighting coefficient β used to compute the average is set to β = 0.5. As for comparison, we also evaluate the performance of the broadcast-based algorithm of [23]. The common point between the algorithm of [23] and the broadcast algorithm described in Section II-C is that they both rely on the broadcast gossip scheme of [3] but the core of the algorithms is rather different as explained in Section I. In order to distinguish both broadcast algorithms, we will designate the algorithm of [23] as the broadcast algorithm with sleeping phases, refering to the fact that each agent does not update its estimates as long as it is not the recipient of a message. On the other hand, we refer to the broadcast algorithm of Section II-C as the broadcast algorithm without sleeping phases.
It is worth remarking that a fair comparison between different stochastic approximation algorithms is generally a delicate task, because the behavior of each particular algorithm is sensitive to the choice of the stepsize. In this paragraph, we simply set γ n = γ 0 /n ξ for all n, where γ 0 > 0 and 0.5 < ξ ≤ 1 are parameters chosen in an ad-hoc fashion. More degrees of freedom are of course possible when choosing γ n , but a complete discussion would be out the scope of this paper. Recall that the algorithm of [23] requires a more specific choice of the stepsize which solely depends on the initial step. We shall denote by γ s 0 the latter initial stepsize used with the algorithm [23], where the upperscript s stands for sleeping phases.
For each algorithm, we evaluate the deviation of the estimates from the global minimizer θ : Note that ∆ n depends on the parameters (s 1 , · · · , s N ). We consider 50 Monte-Carlo runs, each of them consisting of 10000 iterations of each algorithm. For each run, we randomly select the parameters (s 1 , · · · , s N ) according to the uniform distribution on the unit disk G. The kth Monte-Carlo run yields a sequence (∆ (k) n : 1 ≤ n ≤ 10000) for each algorithm. Figure 2 represents the average deviation (1/50) 50 k=1 ∆ (k) n as a function of the number n of iterations. In Figure 2(a), we set N = 50 and the graph is a cycle. In Figure 2(b), we set N = 10 and the graph is a complete graph. It is worth noting that the pairwise gossip algorithm behaved at least as well as both broadcast based algorithms in our experiments. This fact might seem surprising at first glance. Indeed, in the framework of average consensus i.e., when the aim is not to optimize an objective function but simply to compute an average in a distributed fashion [12], the broadcast gossip algorithm of [3] is known to i) reach a consensus faster than the pairwise algorithm of [12] and ii) fail to converge to the desired value. In the context of distributed optimization, a different phenomenon happens: our theoretical analysis showed that broadcast based optimizers do converge to the desired value. However, in this example, there is no clear gain in using a broadcast-based algorithm. Convergence has been established in Theorem 1.Appendix I reveals that part of the perturbation (denoted ξ n ) is due to the fact that 1 T W n = 1 T (see the term ξ (1) n at equation (35)). This part of the perturbation is clearly zero when W n is doubly stochastic. This is the case for the pairwise algorithm, but not for the broadcast algorithm.
As a conclusion, there should be interesting comparisons to make between pairwise optimizer and the broadcast ones.

B. Scenario #2
Consider the distributed power allocation algorithm of Section V. In order to validate the proposed algorithm, we study the 2×2 interference channels shown in Figure 1. As a toy but revealing example, first assume fixed channel gains chosen as must not exceed a maximum power of P 1 = P 2 = 10. The aim is to minimize the weighted sum of the error probabilities as in (32) where β 1 = 2/3, β 2 = 1/3. Strictly speaking, we actually implement a distributed gradient descent w.r.t. to the parameter vector θ in log-scale in order to avoid slow convergence. Figure 3(a) represents the objective function (32) w.r.t. (p 1 , p 2 ) in dB (the x-axis and y-axis are 10 log 10 p 1 and 10 log 10 p 2 respectively). On this example, there exists a unique minimum achieved at point (p 1 , p 2 ) = (10, 5.4). Figure 3(b) represents, on a single run, the trajectory of the estimates θ n,1 = (p 1 n,1 , p 2 n,2 ) of the first agent as a function of the number of iterations. We compare the pairwise and the broadcast gossip schemes. Note that we only plot the result for the broadcast scheme without sleeping phase, as we observed slow convergence of the algorithm of [23] on this particular example. The two upper curves represent the estimate of power p 1 (using a pairwise and a broadcast scheme respectively) while the two lower curves represent the estimate of power p 2 . Each algorithm converges to the desired value First agent's estimates of p 1 and p 2 as a function of the number of iterations γn = 200/n 0.7 for n ≤ 3000 -γn = 30/n 0.7 for n > 3000. (10, 5.4). However, the convergence curve is rather smooth in the pairwise case, and is more erratic in the broadcast case. Indeed, matrices W n are non doubly stochastic in the broadcast scheme. As already explained above, non doubly stochastic matrices introduce an artificial noise term which is the main cause of the erratic shape of the trajectory.
We finally provide numerical results in the case where channel gains are random and time-varying. We assume Rician fading [39,Section 2.4]. For any n, we set EA 1,1 n = EA 2,2 n = 2, EA 1,2 n = EA 2,1 n = 1. The variance of each channel gain is 0.5. Components of A n are assumed independent. Figure 4 represents the average trajectory of the estimates θ n,1 = (p 1 n,1 , p 2 n,2 ) of the first agent as a function of the number of iterations. Trajectories have been averaged based on 50 Monte-Carlo runs. Once again, we observe convergence of the distributed algorithms. The convergence is faster in the pairwise case.

VII. CONCLUSION
We introduced a new framework for the analysis of a class of constrained optimization algorithms for multi-agent systems. The methodology uses recent powerful results about dynamical systems which do not rely on the convexity of the objective function, thus addressing a wider range of practical distributed optimization problems. Also, the proposed framework allows to alleviate the common assumption of doublestochasticity of the gossip matrices, and therefore encompasses the natural broadcast gossip scheme. The algorithm has been proved to converge to a consensus. The interpolated process of average estimates is proved to be a perturbed solution to a differential variational inequality, w.p.1. As a consequence, the average estimate converges almost surely to the set of KKT points.

ACKNOWLEDGMENTS
We are grateful to the anonymous Reviewers for their useful comments. We thank Philippe Ciblat for relevant remarks regarding Section V, as well as Eric Moulines and Filippo Santambrogio for fruitful discussions.

APPENDIX I PROOF OF PROPOSITION 1
From (3) and Assumption 1(e), it is straightforward to show that the decomposition (19) holds if one sets: and ξ n = ξ n where: where Z n is given by (16). We first prove that r n tends to zero. Remark that: Each gradient ∇f i is continuous, and thus uniformly continuous on the compact set G. By Lemma 1, | θ n−1 − θ n−1,i | converges to zero a.s. for any i. Therefore, r n converges a.s. to zero. To prove Proposition 1, it is thus sufficient to show that sup k≥n k =n γ ξ (j) tends to 0 when n → ∞ for j = 1, 2.
First, consider ξ (1) n . Recalling that W n is row-stochastic, it follows that {(1 T W n −1 T ) ⊗ I d }J = 0. Thus, one may write: where the random vector Z n is given by (16). Define M n := n k=1 γ k ξ (1) k . It is straightforward to show that M n is a martingale adapted to (F n ) n≥1 . Indeed, by Assumption 1(c), W n and Z n are independent conditioned on F n−1 . Therefore: We derive: ⊗ I d and where the last equality is obtained by integrating the inner terms w.r.t.
As the spectral radius of matrix E[W T k 11 T W k ] is uniformly bounded, there exists a constant C > 0 such that: By Lemma 1, the first term in the righthand side of the above inequality is finite. Recalling that |Z k | ≤ |Y k |, we deduce from Assumption 1(f) that E|Z k | 2 is uniformly bounded. As k γ 2 k < ∞, we conclude that sup n E|M n | 2 < ∞. This implies that the martingale converges a.s. to a finite random variable M ∞ . Thus, for any k ≥ n, Thus, sup k≥n k =n γ ξ (1) tends a.s. to zero as n → ∞.
We now study ξ (2) n . Clearly, ξ n is a martingale difference noise sequence. Therefore, where we used the fact that |Y k | ≤ |Z k | for any k. Note that E[|Y n | 2 |F n−1 ] is uniformly bounded with sup θ |y| 2 dµ θ by Assumption 1(c). The latter constant is finite by Assumption 1(f). Therefore, sup n≥1 E |Y n | 2 < ∞. Since k γ 2 k < ∞, we conclude that sup k≥n k =n γ ξ (2) tends to zero using the same arguments. This completes the proof of Proposition 1.

APPENDIX IV PROOF OF LEMMA 4
For any ≥ 0, j = 1, · · · , p, define ∂G j := {θ ∈ G : ∃θ ∈ q −1 j ({0}), |θ − θ| ≤ }. It is useful to remark that ∂G 0 j = q −1 j ({0}) ∩ G is the set of points in G for which the jth constraint is active. In particular, that ∂G 0 j ⊂ ∂G j for any ≥ 0. Denote by d H the Hausdorff distance between sets. Define: The key point is to show that lim ↓0 δ( ) = 0. By contradiction, assume that this is not the case. Then there exists a constant c > 0 and a sequence n ↓ 0 such that δ( n ) > c for each n. As there is a finite number of subsets of {1, · · · , p}, it is straightforward to show that there exists a certain subset E ⊂ {1, · · · , p} such that for any n ≥ 1, First note that j∈E ∂G n j is nonempty. Indeed, if it was empty, j∈E ∂G 0 j would be empty as well, so that the Hausdorff distance in the lefthand side of (37) would be d H (∅, ∅) = 0 < c. Thus, for any n ≥ 1, there exists θ n ∈ j∈E ∂G n j such that The sequence (θ n ) n≥1 lies in the compact set G. Thus, there exists a subsequence which converges to some point θ ∈ G. Without loss of generality, we shall still denote this subsequence by (θ n ) n≥1 in order to simplify the notations. We thus consider that lim n→∞ θ n = θ . We shall now prove that θ ∈ j∈E ∂G 0 j . For any n ≥ 1, θ n ∈ ∂G n j . Thus, there exists θ (j) n ∈ G such that q j (θ (j) n ) = 0 and |θ n − θ (j) n | ≤ n . As q j is convex, it is also Lipschitz on the compact set G. Denote by K j its Lipschitz constant on G: |q j (θ n )| = |q j (θ n ) − q j (θ (j) n )| ≤ K j n .
It is useful to remark that, as a by-product of the above proof, we also obtained the following result. Consider any set E ⊂ {1, · · · , p} and assume that there exists a sequence n ↓ 0 s.t. for any n ≥ 1 there exists θ n ∈ j∈E ∂G n j . Due to the above arguments, any limit point of such a sequence (θ n ) n≥1 belongs to the set j∈E ∂G 0 j which is thus nonempty. Let us state this result the other way around: for any E such that j∈E ∂G 0 j = ∅, there exists E 0 > 0 such that for any < E 0 , j∈E ∂G j = ∅. We set 0 = min{ E 0 : j∈E ∂G 0 j = ∅}. It remains to prove (28). Let 0 < < 0 and θ ∈ G. Trivially, θ ∈ j∈A(θ, ) ∂G j .