A Robust Diffusion Estimation Algorithm with Self-Adjusting Step-Size in WSNs

Abstract

In wireless sensor networks (WSNs), each sensor node can estimate the global parameter from the local data in a distributed manner. This paper proposed a robust diffusion estimation algorithm based on a minimum error entropy criterion with a self-adjusting step-size, which are referred to as the diffusion MEE-SAS (DMEE-SAS) algorithm. The DMEE-SAS algorithm has a fast speed of convergence and is robust against non-Gaussian noise in the measurements. The detailed performance analysis of the DMEE-SAS algorithm is performed. By combining the DMEE-SAS algorithm with the diffusion minimum error entropy (DMEE) algorithm, an Improving DMEE-SAS algorithm is proposed for a non-stationary environment where tracking is very important. The Improving DMEE-SAS algorithm can avoid insensitivity of the DMEE-SAS algorithm due to the small effective step-size near the optimal estimator and obtain a fast convergence speed. Numerical simulations are given to verify the effectiveness and advantages of these proposed algorithms.

1. Introduction

The problem of parameter estimation, which is the indirect determination of the unknown parameters from measurements of other quantities [1,2,3,4,5,6], is a key issue in the signal processing field. Distributed estimation has become very popular for parameter estimation in wireless sensor networks. The objective is to enable the nodes to estimate a vector of parameters of interest in a distributed manner from the observed data. Distributed estimation schemes over adaptive networks can be mainly classified into incremental strategies [7,8,9], consensus strategies [10,11], and diffusion strategies [12,13,14,15,16,17,18,19,20,21,22]. In the incremental strategies, data is processed in a cyclic fashion through the network. The consensus strategies rely on the fusion of intermediate estimates of multiple neighboring nodes. In the Diffusion strategies, information is processed at all nodes while the nodes communicate with all their neighbors to share their intermediate estimates. The diffusion strategies are particularly attractive because they are robust, flexible and fully-distributed, such as the diffusion least mean squares (DLMS) algorithm [12]. In this paper, we focus on the diffusion estimation strategies.

The performance of distributed estimation degrades severely when the signals are perturbed by non-Gaussian noise. Non-Gaussian noise may be natural, due to atmospheric phenomena, or man-made, due to either electric machinery present in the operation environment, or multipath telecommunications signals [23,24,25]. Recently, some researchers focus on improving robustness for non-Gaussian noise of distributed estimation methods. The efforts are mainly directed at searching for a more robust cost function to replace the MSE criterion, which is optimal only when the measurement noise is Gaussian. To address this problem, the diffusion least mean p-power (DLMP) based on p-norm error criterion was proposed to estimate the parameters of the wireless sensor networks [26]. The correntropy as a nonlinear similarity measure has been successfully used as a robust and efficient cost function for non-Gaussian signal processing [27,28,29,30]. In [27], two robust MCC based diffusion algorithms, namely the Adapt-then-Combine (ATC) and Combine-then-Adapt (CTA) diffusion maximum correntropy criterion (DMCC) algorithms, are developed to improve the performance of the distributed estimation over network in impulsive noise environments.

The error entropy criterion based on the minimum error entropy (MEE) method also has shown its ability to achieve more accurate estimates than mean-square error (MSE) under non-Gaussian noise [31,32,33,34,35,36,37]. In [31], the diffusion minimum error entropy (DMEE) was proposed. The DMEE algorithm achieved improved performance for non-Gaussian noise with the fixed step-size, but it still suffers from conflicting requirements between convergence rate and the steady-state mean square error. A large step-size leads to a fast convergence rate but a large mean-square error at the steady state. For this problem, variable step-size techniques have been widely used to improve the convergence of diffusion LMS algorithms remarkably by adjusting the step-size appropriately [38,39,40,41]. Lee et al. [38] proposed a novel variable step-size diffusion LMS algorithm which controls the step-size suboptimally to attain the minimum mean square error at each time instant. In [41], Abdolee investigated the effect of adaptation step-sizes on the tracking performance of DLMS algorithms in networks under non-stationary signal conditions. However, to the best of our knowledge, the variable step-size technique has not been extended to the field of distributed minimum error entropy estimation for non-Gaussian noise yet.

In this paper, we incorporate the minimum error entropy criterion with self-adjusting step-size (MEE-SAS) [42] into the cost function in diffusion distributed estimation. Then, we figure out the diffusion-strategy solutions, which are referred to as the diffusion MEE-SAS (DMEE-SAS) algorithm. Numerical simulation results show that the DMEE-SAS algorithm outperforms DLMS, DLMP and DMEE algorithms when the noise is modeled to be non-Gaussian noise. We also design an Improving DMEE-SAS algorithm by using a switching scheme between DMEE-SAS and DMEE algorithms for a non-stationary environment, which tracks the changing estimator very effectively. The Improving DMEE-SAS algorithm can avoid the small effective step-size of the DMEE-SAS algorithm when it is close to the optimal estimator.

We organize the paper as follows. In Section 2, we briefly revisit the minimization error entropy criterion. In Section 3, firstly, we propose the DMEE-SAS algorithm and analyze the mean, mean square and instantaneous MSD performance for the DMEE-SAS algorithm. Then, we propose the Improving DMEE-SAS algorithm for a non-stationary scenario. Simulation results are shown in Section 4. Finally, we draw conclusions in Section 5.

2. Minimization Error Entropy Criterion

Considering the limited computational capability and limited memory space for nodes in real distributed networks, this paper is based on an MEE criterion, which is simple enough and has good estimation accuracy. Important properties of MEE can be found in [32,35,37]. In many real world applications, the MEE estimator can outperform significantly the well-known MSE estimator and show strong robustness to noises, especially when data are contaminated by non-Gaussian noises. In this subsection, we introduce an MEE criterion, which could be used to derive a robust diffusion estimation algorithm with a self-adjusting step-size (DMEE-SAS) algorithm.

The aim of the adaptive signal processing problem is to minimize the difference between the desired and the system outputs, which is defined as error e. For the evaluation of the error entropy, we seek to estimate entropy directly from the error samples. Therefore, system parameters can be estimated by minimizing the Renyi’s entropy of the error e. Renyi’s entropy is given by

$$H_{\alpha }\left(e\right)=\frac{1}{1‒\alpha }log\int q^{\alpha }\left(e\right)de,$$

(1)

where $q^{\alpha }\left(e\right)$ is the probability density function of a continuous error e, and $\alpha$ is a parameter. When parameter $\alpha$ is set as 2, Equation (1) is quadratic Renyi’s entropy. Using a Gaussian kernel with kernel size $\sigma$, we can obtain a convenient evaluation of the integral operator in the formulation of quadratic Renyi’s entropy as follows:

$$\begin{array}{cc}\hfill H_2(e)& =‒log\int q^2\left(e\right)de\hfill \\ & =‒log\int {\left(\frac{1}{N}{\displaystyle \sum _{i=1}^N}{G}_{\sigma }\left(e‒e_i\right)\right)}^{2}de\hfill \\ & =‒log\left(V\right(e)),\hfill \end{array}$$

(2)

where $e=[e_1,e_2,\cdots ,e_N]$ is N independent and identically distributed samples, and the Gaussian kernel is defined as

$$G_{\sigma \sqrt{2}}(e‒e_i)=\frac{1}{\sigma \sqrt{2\pi }}exp(‒\frac{1}{2{\sigma }^2}{(e‒e_i)}^2).$$

The information $V(e)$ is quadratic information potential and is defined as the expectation of probability density function, $V\left(e\right)=E\left[q\left(e\right)\right]$. The quadratic information potential $V(e)$ can be easily estimated by using a simple and effective nonparametric estimator

$$V(e)=\frac{1}{N^2}\sum _{i=1}^N\sum _{j=1}^N{G}_{\sigma \sqrt{2}}(e_j‒e_i)\le V(0)=\frac{1}{\sigma \sqrt{2\pi }}.$$

(3)

The maximum value of the quadratic information potential $V(0)$ will be achieved when $e_1=e_2=\cdots =e_N$. The above results are obtained in the case of batch mode, where the N data points are fixed. For online training methods, in order to reduce calculation costs, the estimate of quadratic information potential can be approximated stochastically by dropping the time average in (3), leading to

$$V(e_i)\approx \frac{1}{L}\sum _{j=i‒L+1}^i{G}_{\sigma \sqrt{2}}(e_i‒e_j),$$

(4)

where L is the latest L samples at time i.

Obviously, to minimize the error entropy is equivalent to maximizing the quadratic information potential since the log is a monotonic function. Therefore, the cost function for the MEE criterion is given by

$$J_{MEE}\left(e\right)=maxV(e).$$

(5)

The selection of the kernel size $\sigma$ is an important step in estimating the information potential and is critical to the success of information theoretic criteria. In particular, increasing the kernel size leads to a stretching effect on the performance surface in the weight space, which results in increased accuracy of the quadratic approximation around the optimal point [43]. In order to ensure accuracy, in the following, a large enough kernel size can be used during the adaptation process, which is commonly used in information theoretic criteria [42,44].

3. Proposed Algorithms

As mentioned in the Introduction, the diffusion minimum error entropy algorithm achieved improved performance for non-Gaussian noise with the fixed step-size, but it still suffers from conflicting requirements between convergence rate and the steady-state mean square error. Therefore, we consider a new cost function, which can achieve fast convergence speed and strong robustness against non-Gaussian noise.

3.1. Diffusion MEE-SAS Algorithm

Consider a connected wireless sensor networks with K nodes. $k\in \{1,2,\dots ,K\}$ is the node index and i is the time index. To proceed with the analysis, we assume a liner measurement model as follows:

$$d_{k,i}=u_{{}_{k,i}}^T{w}^0+v_{k,i},$$

(6)

where $w^0$ is a $M\times 1$ deterministic but unknown vector, $d_{k,i}$ is a scalar measurement of some random process, $u_{k,i}$ is the $M\times 1$ regression vector at time i with zero mean, and $v_{k,i}$ is the random noise signal at time i with zero mean. For each node k, we have

$$V(e_k)=\frac{1}{N^2}\sum _{i=1}^N\sum _{j=1}^N{G}_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i})\le V(0)=\frac{1}{\sigma \sqrt{2\pi }},$$

(7)

where $e_{k,i}=d_{k,i}‒u_{k,i}^{T}w$. The maximum value $V(0)$ will be achieved when $e_{k,i}=e_{k,j},j=i‒L+1,i‒L+2,\cdots ,i$.

We seek an estimate of $w^0$ by minimizing a linear combination of local information. As explained in Section 2, minimizing a linear combination of the local information is equivalent to maximizing a linear combination of the local quadratic information potential $V(e_{k,i})$. To maximize the information potential is equivalent to minimizing the following cost function:

$$\begin{array}{cc}\hfill J_k(w)& ={\displaystyle \sum _{l\in N_k}}{c}_{l,k}{[V(0)‒V(e_l)]}^2\hfill \\ & ={\displaystyle \sum _{l\in N_k}}{F}_l(w),\hfill \end{array}$$

(8)

where

$$F_l(w)={[V(0)‒V(e_l)]}^2.$$

$N_k$ denotes the one-hop neighbor set of node k, and $\left\{c_{lk}\right\}$ are some non-negative cooperative coefficients satisfying $c_{lk}=0$ if $l\notin N_k$, $1_N^T\overline{C}=1_N^T$ and $\overline{C}{1}_N=1_N$. Here, $\overline{C}$ is a $N\times N$ matrix with individual entries $\left\{c_{lk}\right\}$ and $1_N$ is a $N\times 1$ all-unity vector. The gradient of the individual local cost function is given by

$$\nabla J_k(w)=\sum _{l\in N_k}{c}_{lk}{f}_l(w),$$

(9)

where

$$f_l(w)=(\frac{2}{{\sigma }^2{N}^2})(V(0)‒V(e_k))\sum _{i=1}^N\sum _{j=1}^N{G}_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i})(e_{k,j}‒e_{k,i})(u_{k,i}‒u_{k,j}).$$

(10)

We can replace the estimate of quadratic information potential by the stochastic quadratic information potential, leading to

$$\nabla {\widehat{J}}_k(w)=\sum _{l\in N_k}{c}_{lk}{\widehat{f}}_l(w),$$

(11)

where

$${\widehat{f}}_l(w)=(\frac{2}{{\sigma }^{2}L})(V(0)‒V(e_{k,i}))\sum _{j=i‒L+1}^i{G}_{\sigma \sqrt{2}}(e_{k,i}‒e_{k,j})(e_{k,i}‒e_{k,j})(u_{k,j}‒u_{k,i}),$$

(12)

where

$$V(e_{l,i})\approx \frac{1}{L}\sum _{j=i‒L+1}^i{G}_{\sigma \sqrt{2}}(e_{l,i}‒e_{l,j}).$$

(13)

Iterative steepest-descent solution for estimating $w^0$ at each node k can thus be derived as

$$\begin{array}{cc}\hfill w_{k,i+1}& =w_{k,i}‒{\mu }_k\nabla {\widehat{J}}_k(w)\hfill \\ & =w_{k,i}‒{\mu }_k\frac{2}{{\sigma }^{2}L}{\displaystyle \sum _{l\in N_k}}{c}_{lk}[V(0)‒V(e_{l,i})]{\displaystyle \sum _{j=i‒L+1}^i}{G}_{\sigma \sqrt{2}}(e_{l,i}‒e_{l,j})(e_{l,i}‒e_{l,j})(u_{l,j}‒u_{l,i}),\hfill \end{array}$$

(14)

where ${\mu }_k$ is a positive step size. Using the general framework for diffusion-based distributed adaptive optimization [13], an adapt-then-combine (ATC) strategy for diffusion MEE-SAS algorithm can be formulated as

$$\left\{\begin{array}{c}{\phi }_{k,i+1}=w_{k,i}‒{\mu }_k\frac{2}{{\sigma }^{2}L}[V(0)‒V(e_{k,i})]{\displaystyle \sum _{j=i‒L+1}^i}{G}_{\sigma \sqrt{2}}(e_{k,i}‒e_{k,j})(e_{k,i}‒e_{k,j})(u_{k,j}‒u_{k,i}),\hfill \\ w_{k,i+1}={\displaystyle \sum _{l\in N_k}}{c}_{lk}{\phi }_{l,i+1}.\hfill \end{array}\right.$$

(15)

According to Equation (15), the DMEE-SAS algorithm can be seen as a diffusion estimation algorithm with variable step size ${\mu }_k(i)$, where

$${\mu }_k(i)=2{\mu }_k[V(0)‒V(e_{k,i})].$$

(16)

The DMEE-SAS algorithm is described formally in Algorithm 1.

Algorithm 1: DMEE-SAS Algorithm

Initialize: $w_{k,i}=0$

for $i=1:T$

for each node k:

Adaptation

${\mu }_k(i)=2{\mu }_k[V(0)‒V(e_{k,i})]$

${\phi }_{k,i+1}=w_{k,i}‒{\mu }_k(i)\frac{1}{{\sigma }^{2}L}{\displaystyle \sum _{j=i‒L+1}^i}{G}_{\sigma \sqrt{2}}(e_{k,i}‒e_{k,j})(e_{k,i}‒e_{k,j})(u_{k,j}‒u_{k,i})$

Combination

$w_{k,i+1}={\displaystyle \sum _{l\in N_k}}{c}_{lk}{\phi }_{l,i+1}$

end for

In the adaption step of DMEE-SAS algorithm, $V(0)‒V(e_{k,i})$ is close to $V(0)$ when the algorithm starts, and it is close to 0 when the algorithm begins to converge. $V(0)‒V(e_{k,i})$ is always a non-negative scalar quantity, which can accelerate the rate of convergence and achieve small steady-state estimation errors. The fast convergence rate and the small steady-state estimation errors of the DMEE-SAS algorithm can be established against non-Gaussian noise in the measurements.

3.2. Performance Analysis

In this section, we analyze the mean, mean-square and instantaneous MSD performance of the DMEE-SAS algorithm. For tractability of the analysis, here we focus on the case of batch mode. To briefly present the convergence property of the proposed algorithm in terms of global quantities, the following notations are introduced: $M=diag\{{\mu }_1{I}_M,\dots ,{\mu }_K{I}_M\}$, $W_i=col\{w_{1,i},\cdots w_{K,i}\}$, $w^{(0)}=col\{w^0,\cdots ,w^0\}$, ${\tilde{W}}_i=col\{{\tilde{w}}_{1,i}\cdots {\tilde{w}}_{K,i}\}$, $S=col\{s_1(w^0),\cdots ,s_K(w^0)\}$, $C={\overline{C}}^T\otimes I_M$, $I_M$ is the identity matrix.

In order to make the analysis tractable, the followings are assumed:

Assumption 1: The regressor $u_{k,i}$ is independent identically distributed (i.i.d) in time and spatially independent, and $E[u_{k,i}]=0$, $R_k=E[u_{k,i}^T{u}_{k,i}]$.

Assumption 2: The input noise $v_{k,i}$ is super-Gaussian noise. In addition, $v_{k,i}$ and the regressor $u_{k,i}$ are independent from each other. We have $E[v_{k,i}]=0$ and $E[v_{k,i}^2]={\xi }_k$.

Assumption 3: The step-sizes, ${\mu }_k$, $\forall k$, are small enough such that their squared values are negligible.

3.2.1. Mean Performance

Because the input signal and output noises are generated from stationary and ergodic processes, the double time average in Equation (10) can be replaced by the expectation, leading to

$$f_k(w)\approx \frac{2}{{\sigma }^2}(V(0)‒E[G_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i})])E[G_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i})(e_{k,j}‒e_{k,i})(u_{k,i}‒u_{k,j})].$$

(17)

We consider the gradient error caused by approximating the quadratic information potential $V(e_{k,i})$ with their instantaneous values [45]. The gradient error at iteration i and each node k is defined as follows:

$$s_k(w_{k,i})={\widehat{f}}_k(w_{k,i})‒f_k(w_{k,i}).$$

(18)

Using Equation (15), the update equation of the intermediate estimate can be rewritten as

$${\phi }_{k,i+1}=w_{k,i}‒{\mu }_k(f_k(w_{k,i})+s_k(w_{k,i})).$$

(19)

According to [44], when input signal-to-noise ratio is not too low, the error should be small on the whole. Therefore, for a relative large kernal size $\sigma$, when $w=w^0$, ($(e_{k,i}‒e_{k,j})/\sigma \approx 0$ and $G_{\sigma \sqrt{2}}(e_{k,i}‒e_{k,j})\approx \frac{1}{\sigma \sqrt{2\pi }}$. Therefore, the Hessian matrix function $H_k(w^0)$ of $F_l(w)$ is calculated as:

$$\begin{array}{cc}\hfill H_k(w^0)& =\frac{\partial f_k(w)}{\partial w}{|}_{w^0}\hfill \\ & =\frac{\partial }{\partial w}\frac{2}{{\sigma }^2}(V(0)‒E[G_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i})])E[G_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i})(e_{k,j}‒e_{k,i})(u_{k,i}‒u_{k,j})]\hfill \\ & =\frac{2}{{\sigma }^2}(V(0)‒E[G_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i})])E[G_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i}){(u_{k,i}‒u_{k,j})}^T(u_{k,i}‒u_{k,j})‒\hfill \\ & \frac{1}{{\sigma }^2}{G}_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i}){(e_{k,j}‒e_{k,i})}^2{(u_{k,i}‒u_{k,j})}^T(u_{k,i}‒u_{k,j})]\hfill \\ & +\frac{2}{{\sigma }^4}E{[G_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i})(e_{k,j}‒e_{k,i})(u_{k,i}‒u_{k,j})]}^{T}E[G_{\sigma \sqrt{2}}(e_{k,j}‒e_{k,i})(e_{k,j}‒e_{k,i})(u_{k,i}‒u_{k,j})]\hfill \\ & =\frac{1}{{\sigma }^6\pi }E[v_{{}_{k,i}}^2+v_{{}_{k,j}}^2]E[u_{k,i}^T{u}_{k,i}+u_{k,j}^T{u}_{k,j}]\hfill \\ & =\frac{4{\xi }_k{R}_k}{{\sigma }^6\pi }.\hfill \end{array}$$

(20)

Based on the Theorem 1.2.1 of [46], we obtain

$$\begin{array}{cc}\hfill f_k(w_{k,i})& =f_k(w^0)‒({\int }_0^1{H}_k(w^0‒x{\tilde{w}}_{k,i})dx){\tilde{w}}_{k,i}\hfill \\ & =‒({\int }_0^1{H}_k(w^0‒x{\tilde{w}}_{k,i})dx){\tilde{w}}_{k,i},\hfill \end{array}$$

(21)

where ${\tilde{w}}_{k,i}=w^0‒w_{k,i}$ is the weight error vector for node k. We assume that the estimate of each node converges to the vicinity of the unknown vector $w^0$. Therefore, ${\tilde{w}}_{k,i}$ is small enough such that it is negligible, yielding

$$\begin{array}{cc}\hfill f_k(w_{k,i})& \approx ‒({\int }_0^1{H}_k(w^0)dx){\tilde{w}}_{k,i}\hfill \\ & =‒H_k(w^0){\tilde{w}}_{k,i}.\hfill \end{array}$$

(22)

We can also obtain the approximation of the gradient error at the vicinity of $w^0$, which is given by

$$\begin{array}{cc}\hfill s_k(w_{k,i})& \approx s_k(w^0)\hfill \\ & ={\widehat{f}}_k(w^0)‒f_k(w^0)\hfill \\ & =\frac{2}{{\sigma }^{2}L}(V(0)‒\frac{1}{L}{\displaystyle \sum _{j=i‒L+1}^i}{G}_{\sigma \sqrt{2}}(v_{l,i}‒v_{l,j})){\displaystyle \sum _{j=i‒L+1}^i}\hfill \\ & G_{\sigma \sqrt{2}}(v_{k,i}‒v_{k,j})(v_{k,i}‒v_{k,j})(u_{k,j}‒u_{k,i}).\hfill \end{array}$$

(23)

Substituting Equations (22) and (23) into Equation (19), an approximation of intermediate estimate can be obtained at the vicinity of

$${\phi }_{k,i+1}=w_{k,i}+{\mu }_k(H_k(w^0){\tilde{w}}_{k,i}‒s_k(w^0)).$$

(24)

By substituting Equation (24) into the second equation of Equation (15), we get the estimate of unknown parameter as follows:

$$w_{k,i+1}=\sum _{l\in N_k}{c}_{lk}[w_{l,i}+{\mu }_k(H_k(w^0){\tilde{w}}_{l,i}‒s_l(w^0))].$$

(25)

Using global quantities defined above gives the update equation for the network estimate vector as

$$W_{i+1}=C(W_i+MH{\tilde{W}}_i‒MS),$$

(26)

where H collects the Hessian matrix across the network into the global vector $H=diag(H_1(w^0),\cdots ,H_N(w^0))$. Noting that $C{w}^{(0)}=w^{(0)}$, subtraction of both sides of Equation (26) from $w^{(0)}$ gives

$${\tilde{W}}_{i+1}=C(I_{MN}‒MH){\tilde{W}}_i+CMS.$$

(27)

In view of assumptions A1 and A2, ${\tilde{W}}_i$, H and C are independent of each other. Hence, taking expectation of both sides of Equation (27) leads to

$$E[{\tilde{W}}_{i+1}]=E[C](I_{MN}‒MH)E[{\tilde{W}}_i]+CME[S].$$

(28)

We can easily find that $E[S]=col\left\{E[s_1(w^0),\cdots ,s_N(w^0)]\right\}=0$, and Equation (28) has therefore been reduced to this form

$$E[{\tilde{W}}_{i+1}]=E[C](I_{MN}‒MH)E[{\tilde{W}}_i].$$

(29)

From Equation (29), we observe that, in order to be stable for Algorithm 1 in the mean sense, the matrix $E[C](I_{MN}‒MH)$ should be stable. All the entries of $E(C)$ are non-negative and all the rows of it add up to unity. Therefore, to ensure the stability in the mean, it should hold that

$$\left|{\lambda }_{max}\{I_{MN}‒MH\}\right|<1.$$

(30)

We use the notion ${\lambda }_{max}(A)$ to denote the maximum eigenvalue of a Hermitian matrix A. Thus, we note that a sufficient condition for unbiasedness is

$$0<{\mu }_k<\frac{2}{{\lambda }_{max}\left\{H_k(w^0)\right\}}=\frac{{\sigma }^6\pi }{2{\lambda }_{max}\left\{R_k{\xi }_k\right\}}.$$

(31)

3.2.2. Mean-Square Performance

In order to make the presentation clearer, we shall introduce the following notation

$$\mathbf{\Gamma }=(I_{MN}‒MH)C^T\mathbf{\Sigma }C(I_{MN}‒MH).$$

Performing weighted energy balance on both sides of Equation (27) and taking expectations gives

$$E[{∥{\tilde{W}}_{i+1}∥}_{\mathbf{\Sigma }}^2]=E[{∥{\tilde{W}}_i∥}_{\mathbf{\Gamma }}^2]+E[S^{T}M{C}^T\mathbf{\Sigma }CMS],$$

(32)

where $\mathbf{\Sigma }$ is an arbitrary symmetric nonnegative-definite matrix, and the notion ${∥a∥}_{\mathbf{\Sigma }}^2=a^T\mathbf{\Sigma }a$ represents a weighted vector norm for any Hermitian $\mathbf{\Sigma }>0$. By defining

$$r=vec\{E[\mathbf{\Gamma }]\},\theta =vec\{\mathbf{\Sigma }\},$$

where the vec(.) notation stacks the columns of its matrix argument on top of each other. We can modify Equation (32) to

$$E[{∥{\tilde{W}}_{i+1}∥}_{\theta }^2]=E[{∥{\tilde{W}}_i∥}_r^2]+E[S^{T}M{C}^T\mathbf{\Sigma }CMS].$$

(33)

Using the following relationship of the vectorization operator and the Kronecker product [47]:

$$vec(ABC)=(C^T\otimes A)vec\left\{B\right\}.$$

We can obtain that

$$r=\varphi \theta ,$$

(34)

where

$$\varphi =E[(I_{MN}‒MH)\otimes (I_{MN}‒MH)]\beta ,$$

(35)

$$\beta =E[C^T\otimes C^T].$$

Considering Assumption 3, we can approximate Equation (35) as

$$\begin{array}{cc}\hfill \varphi & \approx (I_{M^2{N}^2}‒I_{MN}\otimes MH‒MH\otimes I_{MN})\beta \hfill \\ & =(I_{MN}‒MH)\otimes (I_{MN}‒MH)\beta .\hfill \end{array}$$

(36)

Using the following relationship of the vectorization operator and the matrix trace [47]:

$$Tr\{A^{T}B\}=ve{c}^T(B)vec(A).$$

We find that

$$E[S^{T}M{C}^T\mathbf{\Sigma }CMS]=ve{c}^T\left\{Q\right\}\beta \theta ,$$

(37)

where

$$Q=E[MS{S}^{T}M].$$

Substituting Equations (34) and (37) into Equation (33), we can then reformulate recursion as follows:

$$E[{∥{\tilde{W}}_{i+1}∥}_{\theta }^2]=E[{∥{\tilde{W}}_i∥}_{\varphi \theta }^2]+ve{c}^T\left\{Q\right\}\beta \theta .$$

(38)

It is known that Equation (38) is stable and convergent if the matrix $\varphi$ is stable [48], form the Equation

$${\beta }^T{1}_{M^2{N}^2}=E[(C\otimes I_M)1_{MN}\otimes (C\otimes I_M)1_{MN}]=1_{M^2{N}^2},$$

We know that all the entries of $\beta$ in Equation (37) are non-negative, and all its columns sum up to unity. Using the property $\lambda (\mathbf{A}\otimes \mathbf{A})={\lambda }^2(\mathbf{A})$, the stability of $\varphi$ has the same conditions as the stability of $I_{MN}‒MH$. Therefore, we choose the step size in accordance with Equation (31), which can keep the DMEE-SAS stable in the mean-square sense.

3.2.3. Instantaneous MSD

In order to analyze instantaneous mean-square-error (MSD), we can exploit the liberty of choosing $\theta$ at time i. Then, Expression (38) gives:

$$E[{∥{\tilde{W}}_{n+1}∥}_{{\varphi }^{i‒n}\theta }^2]=E[{∥{\tilde{W}}_n∥}_{{\varphi }^{i‒n+1}\theta }^2]+ve{c}^T\left\{Q\right\}\beta {\varphi }^{i‒n}\theta .$$

(39)

The sum of both sides of Equation (39) for $n=0,1,\dots ,i‒1$ can be given by

$$E[{∥{\tilde{W}}_i∥}_{\theta }^2]=E[{∥{\tilde{W}}_0∥}_{{\varphi }^i\theta }^2]+ve{c}^T\left\{Q\right\}\sum _{n=0}^{i‒1}\beta {\varphi }^n\theta .$$

(40)

We can also adopt a similar way to describe the time instant $i+1$, given by

$$E[{∥{\tilde{W}}_{i+1}∥}_{\theta }^2]=E[{∥{\tilde{W}}_0∥}_{{\varphi }^{i+1}\theta }^2]+ve{c}^T\left\{Q\right\}\sum _{n=0}^i\beta {\varphi }^n\theta .$$

(41)

Subtraction of both sides of Equation (40) from Equation (41) gives

$$E[{∥{\tilde{W}}_{i+1}∥}_{\theta }^2]=E[{∥{\tilde{W}}_i∥}_{\theta }^2]‒E[{∥{\tilde{W}}_0∥}_{{\varphi }^i(I_{M^2{N}^2}‒\varphi )\theta }^2]+ve{c}^T\left\{Q\right\}\beta {\varphi }^n\theta .$$

(42)

By setting

$$\theta =vec\left\{I_{MN}\right\}$$

in Equation (44) and dividing both sides of it by N, the instantaneous MSD for the whole network are computed by:

$$\frac{1}{N}\sum _{k=1}^{N}E\left[{∥{\tilde{w}}_{k,i}∥}^2\right]=\frac{1}{N}E\left[{∥{\tilde{W}}_i∥}^2\right],$$

(43)

where ${\tilde{W}}_i$ can be obtained by the following iteration:

$$\frac{1}{N}E[{∥{\tilde{W}}_{i+1}∥}^2]=\frac{1}{N}E[{∥{\tilde{W}}_i∥}^2]‒\frac{1}{N}E[{∥{\tilde{W}}_0∥}_{{\varphi }^i(I_{M^2{N}^2}‒\varphi )vec\left\{I_{MN}\right\}}^2]+\frac{1}{N}ve{c}^T\left\{Q\right\}\beta {\varphi }^{i}vec\left\{I_{MN}\right\}.$$

(44)

3.3. An Improving Scheme for the DMEE-SAS Algorithm

The too small effective step size near the optimal estimator will hinder the tracking ability of the DMEE-SAS algorithm in a non-stationary environment. In a non-stationary environment, the optimal estimator has small changes. A random-walk model is commonly used in the literature to describe the non-stationarity of the weight vector [48].

Therefore, we try to combine the DMEE-SAS algorithm with the DMEE algorithm [31] in a non-stationary environment where tracking is important. The DMEE-SAS algorithm should be used due to the faster convergence when the algorithm starts, and the DMEE algorithm should be used when the algorithm begins to converge. We use the Lyapunov stability theory [49] to analyze the switching time for each node.

The Lyapunov energy function is a method for analyzing the convergence characteristics of dynamic systems. The cost function can be viewed as a Lyapunov energy function. For the DMEE-SAS algorithm, the continuous-time learning rule is

$$\dot{w}=‒{\mu }_{DMEE‒SAS}\frac{\partial J_k{(w)}_{DMEE‒SAS}}{\partial w}.$$

(45)

The temporal dynamics for the Lyapunov energy that describes the DMEE-SAS algorithm can be obtained as follows:

$$\begin{array}{cc}\hfill {\dot{J}}_k{(w)}_{DMEE‒SAS}& ={\displaystyle \sum _{l\in N_k}}{c}_{lk}(‒2)[V(0)‒V(e_{l,i})]\frac{\partial V{(e_{l,i})}^T}{\partial w}\dot{w}\hfill \\ & ={\displaystyle \sum _{l\in N_k}}{c_{lk}(‒4){\mu }_{k,DMEE‒SAS}{[V(0)‒V(e_{l,i})]}^2∥\frac{\partial V(e_{l,i})}{\partial w}∥}^2.\hfill \end{array}$$

(46)

The individual local energy function for DMEE algorithm can be written as

$$J_k{(w)}_{DMEE}=‒\sum _{l\in N_k}{c}_{lk}V(e_{l,i}).$$

(47)

For the DMEE algorithm, the continuous-time learning rule is

$$\dot{w}=‒{\mu }_{DMEE}\frac{\partial J_k{(w)}_{DMEE}}{\partial w}.$$

(48)

In a similar way, the temporal dynamics for the Lyapunov energy that describes the DMEE algorithm can be obtained as follows:

$$\begin{array}{cc}\hfill {\dot{J}}_k{(w)}_{DMEE}& ={\displaystyle \sum _{l\in N_k}}{c}_{lk}\frac{\partial V{(e_{l,i})}^T}{\partial w}\dot{w}\hfill \\ & ={\displaystyle \sum _{l\in N_k}}{c_{lk}(‒{\mu }_{l,DMEE})∥\frac{\partial V(e_{l,i})}{\partial w}∥}^2.\hfill \end{array}$$

(49)

The switching time is determined as

$$\left|{\dot{J}}_k{(w)}_{DMEE‒SAS}\right|>\left|{\dot{J}}_k{(w)}_{DMEE}\right|\iff V(e_{l,i})<V(0)‒\frac{1}{2}\sqrt{\frac{{\mu }_{l,DMEE}}{{\mu }_{l,DMEE‒SAS}}}(l\in N_k).$$

(50)

When the condition of Equation (50) is met, we should switch from the DMEE-SAS algorithm to the DMEE-SAS algorithm. We introduce the following auxiliary variable:

$$s_{k,i}=\left\{\begin{array}{c}1,V(e_{k,i})<V(0)‒\frac{1}{2}\sqrt{\frac{{\mu }_{k,DMEE}}{{\mu }_{k,DMEE‒SAS}}},\hfill \\ 0,\mathrm{otherwise}.\hfill \end{array}\right.$$

This yields the following algorithm, which we refer to as the improving DMEE-SAS algorithm:

$$\left\{\begin{array}{c}{\phi }_{k,i+1}=w_{k,i}‒s_{k,i}{\mu }_{{}_{k,DMEE‒SAS}}\frac{2}{{\sigma }^{2}L}[V(0)‒V(e_{k,i})]{\displaystyle \sum _{j=i‒L+1}^i}{G}_{\sigma \sqrt{2}}(e_{k,i}‒e_{k,j})(e_{k,i}‒e_{k,j})\hfill \\ (u_{k,j}‒u_{k,i})‒(1‒s_{k,i}){\mu }_{{}_{k,DMEE}}\frac{1}{{\sigma }^{2}L}{\displaystyle \sum _{j=i‒L+1}^i}{G}_{\sigma \sqrt{2}}(e_{k,i}‒e_{k,j})(e_{k,i}‒e_{k,j})(u_{k,j}‒u_{k,i})\hfill \\ w_{k,i+1}={\displaystyle \sum _{l\in N_k}}{c}_{lk}{\phi }_{l,i+1}.\hfill \end{array}\right.$$

(51)

For the purpose of clarity, we summarize the procedure of the Improving DMEE-SAS algorithm in Algorithm 2.

Algorithm 2: Improving DMEE-SAS Algorithm

Initialize:

$w_{k,i}=0$

for $i=1:T$

for each node k:

Adaptation

each node calculates the switching time using Equation (50).

each node updates intermediate estimate ${\phi }_{k,i}$ according to the first equation of Equation (51).

Combination

$w_{k,i+1}={\displaystyle \sum _{l\in N_k}}{c}_{lk}{\phi }_{l,i+1}$

end for

4. Simulation Results

Twenty sensors are randomly placed in a square $100\times 100$ shown in Figure 1. The communication distance is set as 50. In this paper, the performance of the steady-state network MSD [12] is adopted for performance comparison. All of the performance measures are averaged over 100 trials.

We employ the super-Gaussian distribution as the noise model in our simulations. We generate the noise from the zero-mean generalized Gaussian distribution of probability density function $q_V(v)=\propto exp(‒{\left|v\right|}^p)$, where p is a positive shape parameter of probability density function [50]. We set $p=0.6$ to make the noise distribution be super-Gaussian.

(a) In Stationary Environment

Here, the proposed DMEE-SAS algorithm performance is compared with that of some existing algorithms in the literature. We assume the communication link is the ideal link. The unknown parameter vector $w^0$ is set to ${[\frac{1}{\sqrt{6}},\frac{1}{\sqrt{6}},\frac{1}{\sqrt{6}},\frac{1}{\sqrt{6}},\frac{1}{\sqrt{6}},\frac{1}{\sqrt{6}}]}^T$.

We set the window length $L=8$ and kernel size $\sigma =1.5$ for both DMEE and DMEE-SAS algorithms. Furthermore, the p is 1.2 for the DLMP algorithm. The steady state MSD curves are plotted in Figure 2. It is found that the DMEE-SAS algorithm is robust to the non-Gaussian noises and performs better than the DLMP algorithm [26] and DLMS [12]. The DMEE-SAS algorithm achieves a better convergence performance than the DMEE [31] algorithm when the DMEE-SAS and DMEE algorithms achieve comparable performance.

(b) In Non-stationary Environment

Here, the simulations are carried out in the same environments as those shown in stationary environment, except for the optimal estimator $w^0$. We compare the proposed Improving DMEE-SAS algorithm with other algorithms.

Motivated by [51], we assume a time-varying $w^0$ of length 6 as follows:

$$w_i^0=\frac{1}{2}{[a_{1,i},a_{2,i},a_{3,i},a_{4,i},a_{5,i},a_{6,i}]}^T,$$

where $a_{k,i}=[cos(wi+\frac{(k‒1)}{2}\pi )]$ for $k=1,2,3,4,5,6$ and $w=\frac{\pi }{3000}$.

The unknown link is assumed to change at time 6000. In Figure 3, the Improving DMEE-SAS algorithm can detect the weight vector change and the performance of it is better than the DLMS algorithm. We observe that Improving DMEE-SAS and DMEE algorithms achieve comparable performance and Improving DMEE-SAS achieves better convergence performance than the DMEE algorithm. When compared with the DMEE-SAS algorithm, the Improving DMEE-SAS algorithm exhibits a significant improvement in performance when the estimate is close to the optimal estimator. The Improving DMEE-SAS algorithm achieves a low MSD and fast rate of convergence in the non-stationary environment.

5. Conclusions

In this paper, a robust diffusion estimation algorithm with self-adjusting step-size is developed which is called the DMEE-SAS algorithm. The mean and mean square convergence analysis of this new algorithm are carried out, and a sufficient condition for ensuring the stability is obtained. Simulation results illustrate that the DMEE-SAS algorithm can achieve better performance than the DLMS, robust DLMP, and DMEE algorithms in non-Gaussian noise scenario. In addition, we propose the Improving DMEE-SAS algorithm, using it in the non-stationary scenario where the unknown parameter is changing over time. The Improving DMEE-SAS algorithm combined the DMEE-SAS algorithm with the DMEE algorithm, and it can avoid the small effective step-size of the DMEE-SAS algorithm when close to the optimal estimator.