Optimized Quantization in Distributed Graph Signal Filtering
Isabela Cunha Maia Nobre Pascal Frossard









Optimized Quantization in Distributed Graph Signal Filtering
    Isabela Cunha Maia Nobre
        and Pascal Frossard

Isabela Nobre and Pascal Frossard are with the Signal Processing Laboratory (LTS4), École Polytechnique Fédérale de Lausanne (EPFL), Switzerland. (e-mails: isabela.nobre@epfl.ch, pascal.frossard@epfl.ch).
Isabela Nobre is partially supported by CAPES (grant number 88881.174577/2018-01).
    
=========================================================================================================================================================================================================================================================================================================================================================





Distributed graph signal processing algorithms require the network nodes to communicate by exchanging messages in order to achieve a common objective. These messages have a finite precision in realistic networks, which may necessitate to implement message quantization. Quantization, in turn, may generate distortion and performance penalty  in the distributed processing tasks. This paper proposes a novel method for distributed graph filtering that minimizes the error due to message quantization without compromising the communication costs. It first bounds the exchanged messages and then allocates a limited bit budget in an optimized way to the different messages and network nodes. In particular, our novel quantization algorithm adapts to both the network topology and the message importance in a distributed processing task. Our results show that the proposed method is effective in minimizing the error due to quantization and that it permits to outperform baseline distributed algorithms when the bit budget is limited.  They further show that errors produced in nodes with high eccentricity or in the first steps of the distributed algorithm contribute more to the global error. Also, sparse and irregular graphs require more irregular bit distribution. Our method provides one of the first quantization solutions for distributed graph processing, which is able to adapt to the target task, the graph properties and the communication constraints. 




Graph signal processing, quantization, distributed processing, wireless sensor networks.







§ INTRODUCTION




Many different networks, for instance wireless sensor networks, transportation networks, neural networks and social networks can be modeled as graphs where nodes support signal values or features. The field of signal processing on graphs  has been providing many powerful tools to process such signals in diverse applications, such as compression, denoising or reconstruction of sensor data <cit.>. 

Numerous applications require that the signal defined on the network is however processed distributively. Decentralized methods of graph signal processing recently emerged in order to scale to large networks as a way to deal with big data applications, privacy concerns and also bandwidth/energy constraints <cit.>. They can also be necessary when centralized topologies are not viable or suffer from a bottleneck on the central processor, with nodes that are too far to reach the central processor. Another advantage of distributed algorithms compared to centralized ones is that they can add robustness to the network in case of node failures. 

In practice, most nodes in distributed systems have limited computing power and are constrained in the amount of information they can communicate; this makes the design of efficient distributed algorithms essential. In order to enable distributed graph signal processing, linear graph signals operators can for example be approximated by shifted Chebyshev polynomials <cit.>, becoming more amenable to distributive computing for applications such as smoothing, denoising, inverse filtering and semi-supervised learning. There have been many other studies on distributed processing for graph signals or networked data <cit.>, but only few deal with the fact that, in real case scenarios, the network is subject to communication constraints that limit the precision of the messages exchanged by distributed algorithms. 




In this paper, we propose an adaptive quantization scheme for distributed graph signal processing tasks. Quantization approximates a continuous range of values by a set of discrete values and helps decreasing the communication costs. At the same time, it introduces errors that degrade the performance of distributed processing methods.   We build on our previous work <cit.> and focus on the design of a quantization scheme that minimizes the quantization error in graph signal processing tasks, by bounding the transmitted messages and by optimizing the bit allocation. We first propose a distributed processing algorithm where the messages exchanged by network nodes are bounded. We then model the error due to quantization, as a function that depends on the network topology and on the characteristics of the distributed processing task. We further cast an optimization problem to efficiently distribute bit to messages and nodes, so that the total error is minimized under communication cost constraints. We finally solve the resulting bit allocation problem with the Karush-Kuhn-Tucker (KKT) conditions. The performance of our quantizer is evaluated and  compared to the performance of an uniform quantizer. The results show that the bit allocation optimization 
improves the performance in terms of final error compared to
a uniform distribution of the bit budget across messages. They also show that a more regular graph leads towards more uniform bit distribution in the optimal allocation, which confirms the proper adaptation of our solution to the network properties. Also, since the errors propagate through the successive steps of the distributed signal processing operations, we also confirm that it is efficient to allocate a large share of the bit budget in the first steps of the iterative distributed processing algorithm.

Some works have considered different aspects of quantization in graph signal processing. The works in <cit.> and <cit.> briefly studied the effect of quantization in a linear prediction filter applied to graph signals. However the main objective of those works was to design graph filters, without considering specifically the quantization effects in the design itself. On the other hand, the work in <cit.> studied the effect of quantization in the representation of graph signals, mitigating the 
numerical effects caused by the finite-precision machines that centrally realize the filtering process. The authors specifically designed graph filters that are
robust to finite precision effects.  The above works had a centralized graph filtering approach.  There are also some works that focused on distributed optimization with quantization constraints <cit.>, but these algorithms have been developed mainly using consensus protocols, and we are  interested in solving more general processing tasks in this paper, not merely average computation. Closer to our framework, the work <cit.> derived the quantization error in distributed computations of graph signal operators and then proposed an algorithm that learns graph dictionaries to sparsely approximate graph signals while staying robust to the quantization noise. The work however does not focus on improving the quantization but rather on the design of robust graph filters.  To the best of our knowledge,our work the first one that optimizes the quantization scheme for general distributed filtering tasks.

The rest of this paper is organized as follows. In Section 2, we model the graph filter implementation in a distributed way and derive the quantization error. In Section 3, we present our new bounding scheme for processing graph signals distributively, and in Section 4 we describe our optimal bit allocation algorithm. Finally, the results and analysis of the performance of our new framework are shown in Section 5. 







§ DISTRIBUTED GRAPH SIGNAL PROCESSING
 



 §.§ Background

We first consider that the network can be represented by a weighted, undirected graph  𝒢  = (𝒱, ℰ, W), where 𝒱 represents a set of vertices, ℰ represents a set of edges and W is the weight matrix whose entries W[i,j] typically depend on the distance between nodes i and j. The number of nodes in the graph is N = |𝒱|. We define D as the diagonal degree matrix whose elements are the sum of each row of W respectively. The normalized graph Laplacian operator is frequently used in graph signal processing. It is given by ℒ = I - D^-1/2WD^-1/2, which is a real symmetric positive semi-definite matrix. We further denote as  {λ_n}_n=0..N-1  the set of eigenvalues of ℒ, which we order as { 0=λ_0 < λ_1≤λ_2 ≤ ..≤λ_N-1≤ 2}, and its eigenvectors as {𝒳_0,𝒳_1,...,𝒳_N-1}.  Finally, we denote by Λ the diagonal matrix with the eigenvalues {λ_n}_n=0..N-1 on its diagonal and 𝒳 the matrix whose columns are the eigenvectors of ℒ.

A graph signal is a function f:𝒱⟶ℝ defined on the vertices of the graph, which is represented by a vector f∈ℝ^N. The graph Fourier transform of f can be defined as an expansion of the function in terms of the eigenvectors of ℒ, that is f̂(λ_n) = ∑_i=1^Nf(i)𝒳_n(i) <cit.>. Given g(λ) the transfer function of a graph filter,  we can process the signal f_in by computing  f_out = g(ℒ)f_in, where 

    g(ℒ) :=  𝒳[   g(λ_0)                 0;                           ;        0          g(λ_N-1) ]𝒳^T.




 §.§ Distributed Filtering




As many graph signal processing tasks can be represented as filtering operations, we focus on operators of the form of (<ref>) in the rest of the paper. If the signal processing operator g(ℒ ) =   ∑_k=0^Kα_k ℒ^k is a polynomial function of order K of the Laplacian (or if it can be approximated by a polynomial, which is usual in classical settings), with {α_k}_k=0..K as the polynomial coefficients, it is amenable to distributed implementation  <cit.>. The corresponding graph filter

    f_out = ∑_k=0^Kα_k ℒ^kf_in
 
can be implemented iteratively by local processing at each node, and message exchanges between nodes. That is, each node computes the local value of the function f_out by exchanging messages with neighbors in K interactions. Along this process, all nodes together participate in computing the full function f_out. 

The distributed implementation requires the computation of the different powers of the Laplacian matrix that appear in (<ref>). We firstly define  z_k=ℒ^kf_in and begin with
z_0 = f_in. The node n sends its value z_0[n] to its one-hop neighbors in the graph, and all the other nodes do the same. After all values of z_0 are exchanged in the network, all nodes update their local status with the relation z_1=ℒ z_0. To that end, each node n will only calculate its value z_1[n] by doing ℒ_n^Tz_0, where ℒ_n is the line n of ℒ; z_0 is filled with the values of the messages exchanged from the neighbors of node n. Since  ℒ_n is zero for the nodes that are not neighbors of n, the node n does not need the values of z_0 at these nodes to calculate z_1[n]. 

The messages with the values of z_1 are then exchanged between neighbor nodes in the same way as for the values of z_0. Then z_2=ℒ z_1 is obtained in a similar fashion as z_1. This procedure repeats until K iterations, which permits to compute the full function f_out in Eq. (<ref>). Specifically, after knowing {z_0[n],z_1[n],.. z_K[n] }, the node n computes the value of the filtered signal in its own node using the relation <cit.> 

    f_out[n]=( g(ℒ)f_in) [n]  = ∑_k=0^Kα_k z_k[n].

More details are given in <cit.>. 



 §.§ Quantization error


The above computation permits to perfectly compute the response of the graph filter in distributed settings. However, in realistic cases, the messages exchanged by the network nodes have a limited precision. They are typically quantized before transmission and this modifies the outcome of the distributed filtering process. 

The quantized message at node n at step k can be written as z̃_k[n] = z_k[n]+ϵ_k[n],
where  ϵ_k[n] is the quantization error for a message at iteration step k at node n. After its transmission, the distributed update equation at node n becomes z_k+1 = ℒz̃_k[n] = ℒ(z_k[n]+ϵ_k[n] ) and integrates a quantization error term.

We define by

    ϵ_k = [ ϵ_k[0];      ⋮; ϵ_k[N] ]

the vector containing error values for all nodes at step k, and  ϵ = [ϵ_0^T,...,ϵ_K-1^T]^T the vector containing all errors at all iteration steps and nodes. 

By taking into account the  quantization errors that accumulate through all iterations of the distributed processing task, the filtered signal can finally be written as 

    g(ℒ)f_in  = ∑_k=0^Kα_k ℒ^kz_0 +  ∑_k=0^K-1[ ∑_j=1^K-kα_k+jℒ^j] ϵ_k,

as opposed to f_out in Eq. (<ref>) for the perfect settings. More details are available in <cit.>.




§ QUANTIZED DISTRIBUTED FILTERING WITH BOUNDED MESSAGES




 §.§ Quantization error analysis


We first analyze more deeply the impact of quantization. Since the first term of the above expression (<ref>), ∑_k=0^Kα_k ℒ^kz_0, is the filtered signal in a setting without quantization (as in Eq (<ref>)), we can define the second term, 

    Q = ∑_k=0^K-1[ ∑_j=1^K-kα_k+jℒ^j] ϵ_k,

as the total error caused by the accumulated effects of the quantization errors.



We can further make the following observation on the evolution of the quantization error with the iterations of the distributed processing algorithm. At step k, the maximum value (in an absolute sense) of the messages to be transmitted is || ℒ^kf_in||_∞.
Considering that, for p>r>0,  we have ∥ v ∥_r   ≥  ∥ v ∥_p, for any v pertaining to the vector space where these norms are defined, we can write

    ∥ℒ^kf_in∥_∞  ≤  ∥ℒ^kf_in∥_2  ≤  ∥ℒ^k ∥_2·∥ f_in∥_2,

where ∥ℒ^k ∥_2 is a matrix norm induced from the 2-norm for vectors, which can be computed by  ∥ℒ^k ∥_2 = σ_max(ℒ^k ),
where σ_max(ℒ^k ) represents the largest singular value of matrix ℒ^k <cit.>. It corresponds to the square root of the largest eigenvalue of the positive-semidefinite matrix (ℒ^k )^T(ℒ^k ). Since the Laplacian matrix is diagonalizable and symmetric, we can write 
(ℒ^k )^T = (ℒ^k ), and since 

    ℒ = 𝒳Λ𝒳^T,

we have

    ℒ^2k = 𝒳Λ^2k𝒳^T,

which means that the eigenvalues of ℒ^2k  are the same as the eigenvalues of ℒ to a power of 2k. Hence, since

    ∥ℒ^k ∥_2 = σ_max(ℒ^k ) = √(λ_max(ℒ^2k)) =  √(λ_N-1^2k),

and considering that all eigenvalues of the Laplacian are real and positive values, we finally have

    ∥ℒ^kf_in∥_∞  ≤  λ_N-1^k ·∥ f_in∥_2.



This means that, as k increases, the transmitted messages can increase their ranges proportionally to the eigenvalues of ℒ, as shown in Eq. (<ref>). This also means that, at a high value of k, the value of the respective error ϵ_k will be very high, hence increasing the value of the total error. We show below how to bound these transmitted messages such that the quantization errors will remain bounded as well. 




 §.§ Distributed processing with bounded messages


Based on the previous observations, we propose a modification of the classical distributed processing algorithm that permits to bound the range of the exchanged messages. Instead of using the operator given by the normalized Laplacian ℒ at every step of the distributed algorithm, we use 

    ℒ̇=ℒ-I.

Hence, the eigenvalues of ℒ̇ will be bounded in [-1,1], instead of [0,2], with ℒ. Therefore the values of the messages being transmitted at step k will surely not surpass the range of the original signal z_0=f_in, as shown in Eq. (<ref>).

In order to integrate the modified Laplacian operator ℒ̇, we modify the distributed filtering algorithm of the previous section as follows. We start with a scenario without quantization. Firstly we set z_0 = f_in, which is then exchanged with the neighbor nodes, as before. Now, instead of multiplying the received values by ℒ, the nodes rather compute ż_1 = ℒ̇z_0. The value ż_1 is then exchanged with neighbor nodes, in a similar way as the algorithm described in the previous section. The algorithm then proceeds iteratively in the same way. 

In an ideal scenario (without quantization), we have 
ż_k = ℒ̇^kz_0 = (ℒ-I)^kz_0. We observe that, at each step, we can perfectly recover  z_k = ℒ^kz_0 from ż_k. Since the identity matrix commutes with all matrices, ℒ and I also commute. Hence, we can use the Binomial formula and derive

    ż_k = ( ∑_i=0^kki(-1)^k-iℒ^i) z_0  =  ∑_i=0^kki(-1)^k-i z_i,

or equivalently,

    z_k = ż_k-∑_i=0^k-1ki(-1)^k-iz_i  =  ∑_i=0^kkiż_i.




Therefore, with Eq. (<ref>), we can build a distributed algorithm where the values of ż_k are exchanged between neighbor nodes, but only the values z_k corresponding to the iterative solutions of the original distributed filtering algorithm need to be stored. This is useful if the sensors have memory constraints. Note that the term ∑_i=0^kkiż_i is purely combinatorial, that is, it does not depend on any specific component such as the network, data, iteration or task, etc. 



 §.§ Quantization error with bounded messages


In a scenario with quantization, the distributed filtering algorithm with the modified Laplacian operator ℒ̇ is modified as follows. After we set z_0 = f_in, we quantize it as z̃_0 = z_0 +ϵ_0. The values z̃_0 are then exchanged with the neighbor nodes, as before. Now, instead of multiplying the received values by ℒ as in the original algorithm, the nodes rather compute ż_1 = ℒ̇z̃_0 in the bounded scheme, the resulting value is then quantized as ż̃_1 = ż_1 +ϵ_1. The quantized value ż̃_1 is eventually exchanged with neighbor nodes, in a similar way as the algorithm described above.



To recover z_k from ż_k we use the same process as in Eq. (<ref>). However, the recovery is not perfect anymore due to quantization. The quantization error now accumulates through iterations and the value of ż_k becomes 

    ż_k=(ℒ-I)^kz_0+∑_l=0^k-1(ℒ-I)^k-lϵ_l,          for       k>0,
 
while it is given by Eq. (<ref>) in the ideal settings. More specifically, we compute below the quantization error in receiving the output of the filtering process. First, Eq. (<ref>) is equivalent to 

    z_k = z_0 +  ∑_i=1^kkiż_i.

If we replace ż_i in (<ref>)  with the expression in  (<ref>), we obtain

    z_k = z_0 +  ∑_i=1^kki[ (ℒ-I)^iz_0+∑_l=0^i-1(ℒ-I)^i-lϵ_l].

Notice that, if we write ℒ^k as [(ℒ-I)+I]^k and use the Binomial formula, we obtain

    ℒ^k = ∑_i=0^kki (ℒ-I)^i.

Thus, (<ref>) can be written as

    z_k=ℒ^kz_0+∑_i=1^kki∑_l=0^i-1(ℒ-I)^i-lϵ_l,          for       k>0.
 
The distributed filtering of the graph signal f_in in the quantization regime can then be written as

    ḟ_out = g(ℒ)f_in  = ∑_k=0^Kα_k ℒ^kz_0 +∑_k=1^Kα_k ∑_i=1^kki∑_l=0^i-1(ℒ-I)^i-lϵ_l,

from which we derive the total error, which corresponds to the second term in Eq. (<ref>), 

    Q̇  = ∑_k=1^Kα_k ∑_i=1^kki∑_l=0^i-1(ℒ-I)^i-lϵ_l.


It can be rewritten as

    Q̇  = ∑_k=0^K-1[ ∑_i=1^K-kα_k+1∑_j=1^ik+ik+j(ℒ-I)^j]ϵ_k.


For the sake of clarity, we now write 

    F_k[n] =	( H_k^TH_k ) [n,n],

with 

    H_k =  ∑_i=1^K-kα_k+1∑_j=1^ik+ik+j(ℒ-I)^j.


If we combine (<ref>) and (<ref>), we get

    F_k[n] = ∑_i=1^K-kα_k+i∑_j=1^ik+ik+j∑_p=1^K-kα_k+p∑_q=1^pk+pk+q·
    
    (ℒ-I)^j+q[n,n].


Hence,

    Q̇  = ∑_k=0^K-1H_kϵ_k,

and we can finally calculate the mean squared error

    ∥Q̇∥ ^2  = ∑_k=0^K-1∑_l=0^K-1ϵ_k^TH_k^TH_lϵ_l,

and its expected value 

    E_ϵ[        ∥Q̇∥ ^2    ]   = ∑_k=0^K-1∑_l=0^K-1E_ϵ[ϵ_k^TH_k^TH_lϵ_l].






§ OPTIMIZED BIT ALLOCATION



 §.§ Rate-distortion model


In this Section, we now seek to minimize the expected value of the square of the total error in Eq. (<ref>) so that the impact of communication constraints on the distributed computations is minimum. 

The quantization error is a deterministic function of the signal  f_in. Thus, given each f_in, the quantization error is deterministic. But we would like to have an optimal quantization scheme for all input signals. In order to obtain that, we model the input signal at each node as random variables and by consequence, we model the error signal as a random variable as well  <cit.>. Thus, our aim is to minimize the expected value of the total mean square error for these random variables. 


We define x[n,k] as the number of bits used to represent the message sent from node n at step k. We consider the high rate regime and use uniform quantizers for each message. In this case, the error is directly related to the number of bits used for each message. The quantization step size is then determined by the ratio of the total quantization range over the number of quantization intervals, that is

    Δ[n,k] =  2 ·∥ f_in∥_∞/2^x[n,k],
 
the quantization error is 2 ·∥ f_in∥_∞ for al k since we filter the signal with the modified Laplacian ℒ̇ to bound the messages.

It had been shown <cit.> that even though the quantization noise and input signal are deterministically related, when the quantization step size is sufficiently small (that is, in high rate conditions), the quantization noise is uncorrelated with the quantized signal. Furthermore, the probability density function of the quantization noise will be uniform with zero mean and variance Δ[n,k]^2/12. This is called the White-Noise Quantization Error Model  <cit.>. 


It follows that the expected value of the square of the
quantization error on the message transmitted by node n at step k is equal to its variance (since it has zero mean), and if we apply (<ref>) into the variance expression of the White-Noise Quantization Error Model, we obtain

    E_ϵ[ϵ_k[n]^2] = ∥ f_in∥_∞^2/3·2^-2x[n,k].
 
With the White-Noise Quantization Error Model, we can assume that  ϵ_k[n] and ϵ_p[m] are statistically independent for k≠ p or n≠ m. Hence the expected value of the crossed terms in Eq. (<ref>) is zero and we finally obtain the expected value of the total mean squared error as  

    E_ϵ[∥Q̇∥ ^2]  = ∑_k=0^K-1∑_n=0^N-1( H_k^TH_k ) [n,n] E_ϵ[ϵ_k[n]^2],

which is equivalent to

    E_ϵ[∥Q̇∥ ^2] = ∑_k=0^K-1∑_n=0^N-1 F_k[n]  E_ϵ[ϵ_k[n]^2],
 
with (<ref>), the MSE finally reads

    E_ϵ[∥Q̇∥ ^2] =  ∥ f_in∥_∞^2/3∑_k=0^K-1∑_n=0^N-1 F_k[n]  2^-2x[n,k].
 



 §.§ Optimal Allocation


Our objective is now to minimize the total quantization error given the constraint bit budget in the network. To that end, it is necessary to find the values of x[k,n] for every combination of k and n that obey the budget constraint and minimize E_ϵ[∥Q̇∥ ^2]  of Eq. (<ref>). This optimization can be described by the following problem:


    xminimize       E_ϵ[∥Q̇∥ ^2]  
       subject to      ∑_k=0^K-1∑_n=0^N-1x[n,k]d[n]≤ B




Here, d[n] is the degree of node n, which drives the transmission costs, and B is the total bit budget constraint over the whole network. The term x[n,k]d[n] represents the total number of bits used by node n at step k to send the respective message for all d[n] neighbors. When we sum the x[n,k]d[n] term for every k and n, we obtain the total number of bits used to process the signal, which is our constraint in the optimization problem. Since the objective and the constraint functions are both continuously differentiable and convex on the values of x[n,k], we can use the Karush-Kuhn-Tucker (KKT) conditions and they will be sufficient for finding the optimal solution. Also, we observe that the problem satisfies the KKT regularity conditions since the constraint function is affine <cit.>. 

The KKT stationarity condition states that  the solution of the optimization problem (<ref>) will be the values of x[n,k] that satisfy the equation

    ∂/∂ x[n,k][ (E_ϵ[∥Q̇∥ ^2]) +μ  (∑_k=0^K-1∑_n=0^N-1x[n,k]d[n]-B)] = 0,

where ∂/∂ x[n,k] means the partial derivative with respect to x[n,k], and μ is the KKT multiplier. The solution to Eq. (<ref>) is given by

    x[n,k] = -1/2ln(2)ln( 3 μ  d[n]/2∥ f_in∥_∞^2ln(2)F_k[n]).


The complementary slackness condition states that

    μ (∑_k=0^K-1∑_n=0^N-1x[n,k]d[n]-B)= 0,

but μ cannot be equal to zero, otherwise the expression in (<ref>) would present a ln(0) term, which is undefined. This leads to the second term of Eq. (<ref>) being equal to zero instead, and applying (<ref>) into (<ref>) and solving for μ, it results into

    μ=exp( - Bln(4)+ ∑_n=0^N
    			-1 d[n] ∑_k=0^K-1ln(    3d[n]/2∥ f_in∥_∞^2
    				ln	2· F_k[n])  /K∑ _n=0^N-1d[n]) ,

which is always non-negative for being an exponential function, thus guaranteeing that the KKT dual feasibility condition is satisfied. And we can finally replace in (<ref>) to obtain the optimized number of bits for each message. 

In practice, a conventional approach <cit.> is to further round the non-integer values in Eq. (<ref>) to become integers. Then, if some of the values of x[n,k] obtained in (<ref>) are negative or zero, they are replaced by 1, to guarantee that there is a minimal communication between neighboring nodes in every step to keep the whole system synchronized. 

We can analyze the optimal bit allocation solution as follows. We observe that the number of bits x[n,k] in Eq. (<ref>) depends on d[n] and F_k[n]. As the degree d[n] grows, the communication cost in node n grows as well, since it shares its messages with more neighbors. In order to satisfy the budget constraint, x[n,k] has thus to be smaller in a node where d[n] is larger. On the other hand, the factor F_k[n] is related to the topology and to the filter coefficients. From Eq. (<ref>), we can see that it weights the contribution of each individual error ϵ_k[n] in the global error. If we have a big F_k[n], it means that the relative contribution of ϵ_k[n] becomes big, so that its contribution needs to be reduced by increasing  x[n,k]. Also, for the same node n, there are more error terms in the global error computation (in Eq. (<ref>)) when k is low, which means that F_k[n] becomes higher in this case. It further means that the relative contribution of the error terms ϵ_k in the first iterations (small k) of the distributed processing algorithm is higher compared to the error in the later iterations. This is in agreement with the fact that the first error terms lead to higher propagation effects. These observations will be further confirmed in the upcoming result section. 





§ RESULTS




 §.§ Performance of the proposed scheme




  §.§.§ Performance of the optimal quantizer
 
The performance of our quantizer is now evaluated in details. First, we create a network with N=50 nodes that are uniformly placed at random in a unit square. The weight matrix for the network edges is generated based on a thresholded Gaussian kernel function that takes into account the physical distance between nodes. The edges weights are given by W_ij=e^-l_ij^2/θ if the distance l_ij between vertices i and j is less or equal to κ, and zero otherwise. We fix κ=0.2. 

A graph signal is defined as f_in[n]= a[n]^2+b[n]^2-1, where a[n] and b[n] are the coordinates of node n,  
and a random noise with zero-mean normal distribution is added to it. We consider denoising as the distributed graph signal processing task, via a low-pass filter 

    g(λ) = τ/τ+5λ

with τ=3. In order to  implement it distributively, a Chebyshev polynomial approximation of order K is performed, and its polynomial coefficients {α_k}_k=0..K are determined as in <cit.>. 

The distributed graph signal processing task is first performed without quantization. Then, the processing is performed with uniform quantization that is used as baseline solution. In this case, the number of bits used to represent the transmitted messages is the same for every node n and iteration step k. Finally, another processing is performed, using the optimization scheme described in this paper.  In all three cases, the bounding scheme from Section <ref> is used. We calculate the MSE between the output of the quantized and the unquantized schemes for both uniform and optimized quantization. We repeat the entire experiment 1000 times for different networks and compute average performance shown in Figure <ref>. 


		
		
		







It can be seen from Fig. <ref> that the optimization of the bit allocation proposed in this paper clearly improves the performance in terms of MSE if compared to a quantizer with the uniform bit distribution. A bigger gain appears at low bit rate where effective allocation is more important as resource are scarcer. 


To evaluate the efficiency of the bounding scheme, we process the signal using the distributed algorithm proposed in <cit.>, where the messages are not bounded. The settings are the same as the previous experiment. We obtain a MSE of 0.77 for 300 exchanged bits, a MSE of 0.55 for 400 exchanged bits and a MSE of 0.26 for 600 exchanged bits. These results show that, regardless of optimizing the bit distribution or not, our bounded range algorithm yields much lower MSE values compared to the baseline algorithm for the same bit rates.  This is due to the fact that, without bounding the transmitted messages, they grow substantially and in consequence, the quantization errors grow as well. 




We can analyze the effects of the graph structure and the number of iterations on the performance of the optimization scheme. First, the same processing (with bounding scheme) is applied into a graph with a smaller θ, which results in a bigger discrepancy between edges weights, that is, a less regular graph. The communication constraints are however unchanged, since the number of edges only depends on κ, which remains constant in our experiments. The results can be seen in Fig. <ref>, where the difference between the uniform and the optimized bit distribution is slightly bigger than in the original experiments with a more regular graph (Fig. <ref>). It means that a more regular graph tends towards a more uniform bit distribution, which seems reasonable. Finally we look at the impact of the polynomial order K on the performance. We run experiments similar to the previous ones, but with a bigger value of K. The difference between the uniform and the optimized bit allocations is shown in Fig. <ref>. We see that the gain due to optimal allocation is bigger than in the previous experiments with the same network (Fig. <ref>). When k grows, the errors propagate more across iterations and the optimized bit allocation tries to compensate it by allocating more bits in the first iteration steps. This allocation improves the performance of the optimized scheme with respect to the uniform scheme, and this improvement effect is more visible for greater K, since the errors propagate for more iterations. 




  §.§.§ Performance with different filters

We now analyze how the performance of our algorithm changes with different filters implemented in the distributed processing task. Graphs are generated in the same way as in Subsection <ref> with θ=0.001. The same input graph signal is now filtered by different operators with the polynomial approximation order fixed at K=9. Two different tasks are performed with varying filter parameters. First, we perform denoising with distributed Tikhonov regularization. Such a task can be implemented by applying the graph filter

    g(λ)=τ/τ+2λ^r

to the input signal. We fix the value of τ at 10, and we choose for r the values of 1 or 3. The second task consists in distributed smoothing, which can be performed by applying the heat kernel

    g(λ) = e^-τλ/λ_max

to the input signal. We consider τ taking the values 1 or 10. The results of the filter processing using different quantization methods are plotted in Figures  <ref> and <ref>.

Notice that, for the Tikhonov regularization, the higher r, the less smooth the filter. As for the heat kernel, the higher τ, the less smooth the filter. As we can see in Figures <ref> and <ref>, when the filter becomes smoother, the values of the MSE become smaller (both for optimized and uniform bit allocation). This happens because smooth filters tend to have smaller values of α_k, the filter coefficients derived from the Chebyshev polynomial approximation. Thus, the quantization errors are multiplied by smaller factors, resulting in smaller quantization error. 







  §.§.§ Performance with different input signals


We now analyze the performance for different input graph signals. A graph is generated in the same way as in subsection <ref> with θ=0.001, but now the input signal is a uniform random vector, whose values vary from 0 to 1. A Gaussian noise with zero-mean is added to the signal and the denoising is performed by distributively applying the same low-pass filter of Eq. (<ref>) to the signal, with its Chebyshev polynomial approximation of order K=9. The experiment is repeated 1000 times, with a different input signal at each time (while the other parameters remain fixed) and the average performance is computed and shown in Figure <ref>. We can see that the optimized bit allocation improves the filtering performance compared to the uniform bit distribution. 





  §.§.§ Performance on a real dataset

We now consider a dataset containing the rain gauge time series components of a curated set of historical daily rainfall data for the Amazonian rainforest region in Brazil <cit.>. The data is provided by the Brazillian water management agency Agencia Nacional de Aguas (ANA). Daily rainfall intensities in [mm/day] were recorded by 850 rain gauges spread in a large region spanning the southern Amazonian rainforest to the Cerrado biomes of Brazil from 1926 to 2013. We define a geographical graph, where the nodes of the graph consist of the rain gauge stations. The weight matrix is generated based on a thresholded Gaussian kernel function that takes into account the geographical distance between stations. The year of 2010 was chosen for containing the most complete data, and the average daily rain of that year was computed for each rain site to form the signal value at each node. This results in a graph signal of dimension 850. This graph can be visualized in Fig. <ref>.

A Gaussian noise with zero-mean is added to the rain signal and denoising is performed by distributively applying the same low-pass filter as in  Eq. (<ref>)  (actually its Chebyshev polynomial approximation of degree K) to the signal and results are shown in Fig <ref>. When the signal is  processed with the algorithm proposed in <cit.>, where the messages are not bounded,  we obtain a MSE of 93.7 for 60000 exchanged bits and a MSE of 87.5 for 80000 exchanged bits. Comparing these results with Fig. <ref>, we can see that our proposed modified algorithm yields much lower MSE values compared to the baseline algorithm. Also in Fig. <ref> we see that the optimized bit allocation further improves the MSE in the bounded scheme, when compared to the uniform bit distribution. 




		



The rain signal is smooth, and the filter is the same as the one used in Subsection <ref>, so the differences in amplitude observed between Fig. <ref>  and the results in Subsection <ref> stem from the range of the input signal and the topology. The graph of the rain topology dataset has more nodes, which increases the number of computations necessary to filter the signal. With the increase in the number of computations, more errors are committed and accumulated, hence why the MSE values tend to be higher. 




 §.§ Analysis of the bit allocation
 



  §.§.§ Evolution of F_k[n]
 

Since F_k[n] influences the optimal bit allocation x[n,k], it is important to understand which factors influence it. A new graph is generated in the same way as in Subsection <ref>, also the same graph signal is considered to be filtered by the same filter of Eq. (<ref>) with an approximation of order K = 9 and edge weight scale factor θ=2.  We compute F_k[n] with Eq. (<ref>). Now we observe the relationship of F_k[n] with n and k. In order to analyze the relation with n first, we plot F_0[n] over the graph in Fig. <ref>, namely F_k[n] at step k=0. At this initial step we have the highest variance of F_k[n], so we can visualize its relationship with n better. The colors represent the values of F_0[n] at different nodes. It can be noted that groups of nodes that are more isolated tend to have higher values of F_k[n]. As a result, a higher number of bits needs to be allocated to those nodes.  This means that the contribution to the global error of these nodes becomes higher. That happens because isolated nodes tend to amplify their errors, whereas more connected nodes tend to dissolve their errors into the network.

We now plot in Fig. <ref> the values of F_k[n]  versus k. For a given k, we show different values of F_k[n], one value for each n. We observe that the general trend is to F_k[n] decrease with k, as expected since initial errors indeed propagate more than the later quantization errors as discusses earlier.


















  §.§.§ Variance of x[n,k]


We now have a deeper look at the actual optimal bit allocation in order to understand for which topologies the variance of x[n,k] is high or low. A low variance means that the optimal solution approaches the one of a uniform solution. In this case, it might be simpler to perform a uniform bit allocation even with slightly worse MSE results.

In another set of experiments, 300 graphs generated from a weight matrix based on a thresholded Gaussian kernel function were generated with the edge weight scale factor θ fixed at 2. The edge threshold value κ varies from 0.15 to 0.4, so that we can have graphs with different values of mean degree and variance. The unconnected graphs were removed from the experiment. The graphs were chosen with N=50 nodes. The same graph signal as in Subsection <ref> is generated and filtered with the same filter of Eq. (<ref>) with polynomial approximation order K=9.  We use Eq. (<ref>) for computing x[n,k]. 


In Fig. <ref>, we can see the relationship between the degrees' variance and mean with the variance of x[n,k] with respect to n (x[n,k] is averaged with respect to k). Each dot corresponds to the experiment resulted with the use of a different graph. It can be noted from Figs. <ref>  that if the graph has low degree variance, that is, it is a more regular graph, it tends to have a more uniform bit distribution. At the same time, graphs with higher degree mean (that is, denser graphs) also tend to have more uniform bit distribution, which is expected since more connected graphs tend to have the quantization errors smoothed out more uniformly in the network. These facts mean that when the graph is both sparser and irregular it will require a more irregular bit distribution hence optimal bit allocation is important. 

In Fig. <ref>, we can further see the relationship with the variance of the nodes' eccentricities. The eccentricity of a graph node n is the maximum distance between n and any other node of the graph.  In  Fig. <ref> we can also notice the tendency that graphs with more irregular eccentricities require more irregular bit distribution, because more eccentric nodes require more bits than the less eccentric ones, as discussed previously. Which means that, for graphs with more irregular eccentricities, the benefit of solving the optimal bit allocation problem is higher. 












§ CONCLUSION


In this paper, we have shown how the quantization error in distributed graph signal processing tasks can be minimized by bounding the transmitted messages and by optimizing the bit allocation in the network. Our method is important in cases where we have a bit budget constraint. The optimal bit allocation varies with the network nodes and the steps of the iterative filtering process; more specifically, nodes at the border of the network and the initial steps of the iterative algorithm tend to require more bits. Its variance varies with the topology, where irregular and sparse topologies lead to high variance of the optimal bit allocation. Experimental results show that our distributed processing algorithm substantially decreases the quantization error with regard to previous solutions and to a uniform bit allocation scheme. Future work will investigate non-uniform quantization for each message (varying quantization step size), which will improve the performance in low-rate regime. 













IEEEbib






§ BIT ALLOCATION IN THE NETWORK 


We now observe the relationship between F_k[n] and x[n,k]. In the same experiment of Subsection <ref>, we show in Fig. <ref> the values of x[n,k] plotted as the colors of the graph. Comparing Figs. <ref> and <ref>, we notice that, when the values of x[n,k] are high,  they are also high in F_k[n]; but the opposite does not necessarily hold true. This happens because the values of F_k[n]  have the tendency of being higher on the more isolated nodes, but the values of x[n,k] take into consideration communication costs too. Thus they tend to be higher on those nodes that are not only at the border of the graph but have also small degrees. 






§ ERROR PROPAGATION






We do an analysis of the error propagation for illustration purpose. Three graphs are generated in the same way as in subsection <ref>. The same input graph signal is filtered by distributively applying the same low-pass filter of Eq. (<ref>) to the signal, with its Chebyshev polynomial approximation of order K = 17. We fix θ=2 and vary κ. When κ is small, the resulting graph is sparse. Thus, the three generated graphs have the same set of nodes but have different sparsity values, according to the chosen value of κ.  Figs. <ref>, <ref> and <ref> show the three graphs for the values of κ = 0.18, 0.25 and 0.3 respectively.

We are interested in observing how the individual error ϵ_0[n], at a specific value of n and generated in the first iteration step k=0, is propagated through the network.  To achieve this, we perform an experiment where all the other nodes send their messages without quantization  (that is, as the perfect representation of their true value) and node n only quantizes at step k=0, while for the next values of k it sends the messages unquantized. By doing this we can assure that all observed errors in the distributed signal processing stem from  ϵ_0[n] and its propagation. The bounded scheme proposed in this paper is used for all cases. 




The chosen node n  is highlighted within an orange circle. The colors represent the absolute difference between z_17[n] in an unquantized processing (true value) and its value in the experiment mentioned where only  z_0[n] is quantized but all the other values (different k's or different n's) are not. We choose to represent this difference at k=17 since it is the filtering final step and we can thus see the accumulated propagation of ϵ_0[n].  The difference is plotted in log scale. 

We can see in Fig. <ref>  that the node n is only connected to another node and disconnected from the rest. Thus, ϵ_0[n] does not propagate to most of the nodes resulting in these unconnected nodes having no errors on z_17[n]. On the other hand, node n and its neighbor present very high error, which suggests that the error was amplified among them. When we compare to Figs. <ref> and <ref>, where the node n is connected to the rest of the network, ϵ_0[n] propagates to the other nodes and thus the value of the error at n drops considerably. When we consider the last Figure in particular, we notice that there is barely any error in any node. The global MSE for each case is 6.97e-07 for κ=0.18, 1.10e-07 for κ=0.25 and 0.47e-7 for κ =0.3, respectively. Notice that Fig. <ref> has the highest eccentricity for node n, which seems to suggest that nodes with high eccentricity, that is, more isolated nodes, tend to have their errors amplified, whereas the ones that are more central tend to have their errors dissolved into the network. This is coherent to what we previously observed in Section <ref>.












§ VARIANCE OF  THE BIT ALLOCATION FOR BINOMIAL DISTRIBUTION



Another 300 graphs were generated considering node degrees that follow a binomial distribution. This distribution was chosen was because the mean and variance of this distribution can be chosen independently (as long as the variance is below the mean). The graphs were chosen with N=50 nodes. The same graph signal as in Subsection <ref> is generated and filtered with the same filter of Eq. (<ref>) with polynomial approximation order K=9.  We use Eq. (<ref>) for computing x[n,k]. 

In Fig. <ref>, we can see the relationship between the degrees' variance and mean with the variance of x[n,k] with respect to n (x[n,k] is averaged with respect to k). Each dot corresponds to the experiment resulted with the use of a different graph. We can notice that the results of Fig. <ref> are similar to those of Fig. <ref>, which means our observations are consistent for different types of graphs. 

In Fig. <ref>, we can further see the relationship with the variance of the nodes' eccentricities. In orange, we have graphs generated from binomial distribution node degrees and in blue we have Gaussian kernel generated graphs. Here, the same tendency is again noticed for both types of graphs.









