IN recent years,the problems of sparse signal recovery have been widely investigated because of the emergence of new signal sampling theory which is known as compressed sampling or compressed sensing (CS)^{[1, 2, 3]}. Using less measurements than those required in Shannon sample principle,a sparse signal can be recovered with overwhelming probability by solving a 1-norm minimization problem. This convex optimization problem can be solved by various methods,such as basis pursuit de-noising (BPDN),least absolute shrinkage and selection operator (LASSO), Dantzig selector (DS),etc. Many researchers have attempted to discuss this problem in the classic framework of signal estimation,such as Kalman filtering (KF). In ^[4],the Dantzig selector was used to estimate the support set of a sparse signal, and then a reduced-order Kalman filter was employed to recover the signal. In ^[5],the authors presented an algorithm based on a hierarchical probabilistic model that used re-weighted $\ell_1$ minimization as its core computation and propagated second order statistics through time similar to classic Kalman filtering. The dual problem of 1-norm minimization was solved in ^[6] by introducing a pseudo-measurement equation into Kalman filtering, and the Bayesian interpretation of this method was provided in ^[7].

Distributed sensor network is an important way of data acquisition in engineering,and the distributed estimation or filtering problems have been paid much attention recently^{[8, 9, 10]}. In ^[11],three types of distributed Kalman filtering algorithms were proposed. A distributed high-pass consensus filter was used to fuse local sensor measurements,such that all nodes could track the average measurement of the whole network. These algorithms were established based on Kalman filtering in information form, and analysis of stability and performance of Kalman-consensus filter was provided in ^[12]. It is well known that robustness of Kalman filtering is not satisfactory. In ^[13],a design method of distributed $H_{\infty}$ filtering for polynomial nonlinear stochastic systems was presented,and the filter parameters were derived in terms of a set of parameter-dependent linear matrix inequalities (PDLMIs) such that a desired $H_{\infty}$ performance was achieved. A $H_{\infty}$-type performance index was established in ^[14] to measure the disagreement between adjacent nodes,and the distributed robust filtering problem was solved with a vector dissipativity method. Nevertheless,the distributed sparse signal estimation problem has not been adequately researched yet in the framework of $H_{\infty}$ filtering.

In this paper,we aim to combine the pseudo-measurement method with $H_{\infty}$ filtering,and develop a distributed $H_{\infty}$ filtering algorithm to estimate a sparse signal sequence using the measurements distributed over a sensor network. A $\ell_1$-regularized $H_{\infty}$ filter is established at first and the pseudo-measurement equation can be interpreted as a 1-norm regularization term added to the classic $H_{\infty}$ performance index. To develop the distributed algorithm,a high-pass consensus filter is embedded into $\ell_1$-regularized $H_{\infty}$ filter in information form,such that all node filters can reach a consensus on the estimates of sparse signals asymptotically and satisfy the prescribed $H_{\infty}$ performance constraint.

The remainder of this paper is organized as follows. Section Ⅱ gives a brief overview of basic problems in compressed sampling,and introduces the sparse signal recovery method using Kalman filtering with embedded pseudo-measurement. In Section Ⅲ,the centralized $\ell_1$-regularized $H_{\infty}$ filtering method is established by means of Krein space Kalman filtering,and the corresponding information filter is derived. A high-pass consensus filter is employed to develop the distributed filtering algorithm in Section Ⅳ. Simulation results are given in Section V to demonstrate effectiveness of the proposed method,and Section Ⅵ provides some concluding remarks.

Notation. The notation used here is fairly standard except where otherwise stated. The support set of ${\pmb x}\in\mathbf{R}^n$ is defined as ${\rm Supp}\{{\pmb x}\} = \{i | {\pmb x}(i)\neq 0\}$. $\|{\pmb x}\|_0$ is the cardinality of ${\rm Supp}\{{\pmb x}\}$. The 1-norm and 2-norm of ${\pmb x}$ are defined as $\parallel x{\parallel _1} = \sum\nolimits_{i = 1}^n {|x(i)|} $ and $\parallel x{\parallel _2} = \sqrt {{x^{\rm{T}}}x} $,respectively. $x_P^2$ means the product ${x^{\rm{T}}}Px$. ${\rm{sgn}}( \cdot )$ is the sign function. ${\rm{E}}[\xi ]$ defines the mathematical expectation of random variable $\xi $. ${\rm{N}}(\mu ,{\sigma ^2})$ stands for normal distribution with mean $\mu $ and variance ${\sigma ^2}$. ${{\rm{U}}_{{\rm{int}}}}[a,b]$ defines the integer uniform distribution in the interval $[a,b]$. ${\rm{vec}}\{ {x_1}, \cdots ,{x_n}\} : = {\left[ {x_1^{\rm{T}} \cdots x_n^{\rm{T}}} \right]^{\rm{T}}}$. ${\rm{diag}}\{ {A_1}, \cdots ,{A_n}\} $ denotes a block-diagonal matrix whose diagonal blocks are given by ${A_1}, \cdots ,{A_n}$. $S0$ (respectively,$S \succ 0$) means a real symmetric matrix $S$ is positive semi-definite (respectively, positive definite).

This section briefly overviews some basic issues in compressed sampling theory,and introduces the pseudo-measurement embedded Kalman filtering method proposed in ^[6]. More details can be found in ^{[1, 2, 3, 6, 7]}.

A signal ${\pmb x} \in\mathbf{R}^n$ is sparse if $\|{\pmb x}\|_0\ll n$,and ${\pmb x}$ is called $s$-sparse if $\|{\pmb x}\|_0=s$. Assume $\{{\pmb x}_k\}_{k=0}^{N}$ is an unknown discrete-time sparse signal sequence in $\mathbf{R}^n$ and ${\pmb x}_k$ evolves according to the following dynamic model

where $A_k \in \mathbf{R}^{n\times n}$ is the state transition matrix,$\{{\pmb w}_k\}_{k=0}^{N-1}$ is a zero-mean white Gaussian sequence with covariance $Q_k \succeq 0$,and ${\pmb x}_0\sim {\rm N}({\pmb x}_{0|-1},P_{0|-1})$. The $m$-dimensional linear measurement of ${\pmb x}_k$ is

where $C_k \in \mathbf{R}^{m\times n}$ is the measurement matrix, and $\{{\pmb v}_k\}_{k=0}^{N}$ is a zero-mean white Gaussian sequence with covariance $R_k \succeq 0$. Extracting signal ${\pmb x}_k$ from measurement ${\pmb y}_k$ is an ill-posed problem in general when $m<n$.

As shown in ^{[1, 2, 3]},${\pmb x}_k$ can be accurately recovered by solving the following optimization problem

But the optimization problem (3) is NP-hard and cannot be solved effectively. Fortunately,the authors of ^[1] have shown that if the measurement matrix $C_k$ obeys the so-called restricted isometry property (RIP),the solution of (3) can be obtained with overwhelming probability by solving the following convex optimization problem

This is a fundamental result in compressive sampling,and one of the deep results is that for a $s$-sparse signal in $\mathbf{ R}^n$,only on the order of $m = s\log n$ samples are needed to reconstruct it.

For the system given in (1) and (2),Kalman filtering can provide an estimate of ${\pmb x}_k$ which is a solution of the following unconstrained $\ell_2$ minimization problem

In ^[6],the authors discussed the stochastic case of (4)

and its dual problem

By constructing a pseudo-measurement equation

where $H_k = {\rm sgn}({\pmb x}_k^\mathrm{T})$,and $\varepsilon'_k\sim{\rm N}(0,\sigma_k^2)$ serves as the fictitious measurement noise,the constrained optimization problem (7) can be solved in the framework of Kalman filtering and the specific method has been summarized as CS-embedded KF with 1-norm constraint (CSKF-1) algorithm in ^[6].

In the pseudo-measurement equation (8),the measurement matrix $H_k$ is state-dependent,and it can be approximated by $\hat{H}_k = {\rm sgn}(\hat{{\pmb x}}_{k}^\mathrm{T})$. The divergence of this approximation was discussed in ^[7]. Furthermore, $\sigma_k$ is a tuneable parameter which determines the tightness of constraint on 1-norm of the state estimate $\hat{{\pmb x}}_k$.

In this section,the pseudo-measurement equation is combined with $H_{\infty}$ filtering,and a $\ell_1$-regularized $H_{\infty}$ filtering is developed for estimating a sparse signal sequence.

Define the augmented measurement equation using (2) and (8) as

where $\bar{{\pmb y}}_k=\left[ \begin{array}{c} {\pmb y}_k \\ 0 \\ \end{array} \right]$, $\bar{C}_k = \left[ \begin{array}{c} C_k \\ H_k \\ \end{array} \right]$, $\bar{{\pmb v}}_k = \left[ \begin{array}{c} {\pmb v}_k \\ \varepsilon'_k \\ \end{array} \right]$, and denote $\bar{R}_k = \left[ \begin{array}{cc} R_k & 0 \\ 0 & \sigma_k^2 \\ \end{array} \right]$. Consider the system described by (1) and (9) with Gramian matrix

where $\delta_{kj}$ is Kronecker delta function. We aim to design a full-order filter in the form of

such that for all non-zero ${\pmb w}_k$ and ${\pmb v}_k$,the filtering error $\tilde{{\pmb x}}_k={\pmb x}_k-\hat{{\pmb x}}_k$ satisfies the following $\ell_1$-regularized $H_{\infty}$ performance constraint

where $\gamma>0$ is a prescribed scalar.

There exists a filter in the form of (11) achieving the performance (12),if and only if the filtering error covariance matrix $P_{k|k}$ satisfies

for $0\leq k < N$,where the initial value is $P_{0|-1}$,and the predicted error covariance matrix $P_{k|k-1}$ satisfies the Riccati recursion

The filtered estimates $\hat{{\pmb x}}_{k|k}$ are recursively computed as

where

and the predicted estimates

Proof. Define the quadratic performance function

where $\breve{{\pmb y}}_k=\left[ \begin{array}{c} {\pmb y}_k \\ 0 \\ \hat{{\pmb x}}_{k|k} \\ \end{array} \right]$, $\breve{C}_k=\left[ \begin{array}{c} C_k \\ H_k \\ I \\ \end{array} \right]$ and $\breve{R}_k=\left[ \begin{array}{ccc} R_k & 0 & 0 \\ 0 & \sigma_k^2 & 0 \\ 0 & 0 &-\gamma^2I \\ \end{array} \right]$. Then the proposed filter is not hard to obtain according to Krein space Kalman filtering in ^[15]. Details of the proof are omitted.

In ^[7],the pseudo-measurement equation (8) was interpreted in Bayesian filtering framework,and a semi-Gaussian prior distribution was discussed. Here,according to Theorem 1,the pseudo-measurement equation in $H_{\infty}$ filtering can be interpreted as a 1-norm regularization term added to the classic $H_{\infty}$ performance index.

To establish a distributed filtering algorithm,the result in Theorem 1 will be rebuilt in information form,so by denoting

we can obtain the $\ell_1$-regularized $H_{\infty}$ information filter.

The filter established in Theorem 1 is equivalent to the following information form.

Measurement update:

Pseudo-measurement update:

Time update:

Proof. According to matrix inversion lemma,we have

then (15) can be rewritten as

On the other hand,by substituting (22) into (24),we have

It is easy to obtain the following equation

so

which means (15) is equivalent to (22) and (24). Equation (13) is obtained immediately by combining (21) and (23).

In this section,we will develop a distributed $H_{\infty}$-consensus filtering based on the results presented in Section Ⅲ.

Consider a sensor network whose topology is represented by an undirected graph $\mathcal{G}=(\mathcal{V,E,A})$ of order $r$ with the set of nodes $\mathcal{V}=\{1,2,\cdots,r\}$,the set of edges $\mathcal{E} \subseteq \mathcal{V}\times\mathcal{V}$,and the adjacency matrix $\mathcal{A}=[a_{ij}]_{r\times r}$ with nonnegative adjacency elements $a_{ij}$. An edge of $\mathcal{G}$ is denoted by an unordered pair $(i,j)$. The adjacency elements associated with the edges are positive,i.e., $a_{ij}>0\Leftrightarrow (i,j)\in \mathcal{E}$. Node $j$ is called a neighbor of node $i$ if $(i,j)\in \mathcal{E}$ and $j\neq i$. The neighbor set of node $i$ is denoted by $\mathcal{N}_i$. Assume $\mathcal{G}$ is strongly connected.

Assume the measurement of sensor node $i$ is described by the linear model (2). Denote $\bar{{\pmb y}}_{i,k}=\left[ \begin{array}{c} {\pmb y}_{i,k} \\ 0 \\ \end{array} \right]$, $\bar{C}_{i,k}=\left[ \begin{array}{c} C_{i,k} \\ H_{k} \\ \end{array} \right]$, $\bar{{\pmb v}}_{i,k}=\left[ \begin{array}{c} {\pmb v}_{i,k} \\ \varepsilon'_{i,k} \\ \end{array} \right]$, we can get the following augmented measurement equation of sensor $i$

where $C_{i,k}$ is the measurement matrix of sensor $i$,$H_{k} = {\rm sgn}({\pmb x}_k^\mathrm{T})$,${\pmb v}_{i,k}\sim{\rm N}(0,R_{i,k})$,$\varepsilon'_{i,k}\sim{\rm N}(0,\sigma_{i,k}^2)$. And let $\bar{R}_{i,k}={\rm diag}\{R_{i,k},\sigma_{i,k}^2\}$. Define $\mathcal{Y}_k = {\rm vec}\{\bar{{\pmb y}}_{1,k},\cdots,\bar{{\pmb y}}_{r,k}\}$,$\mathcal{V}_k = {\rm vec}\{\bar{{\pmb v}}_{1,k},\cdots,\bar{{\pmb v}}_{r,k}\}$, $\mathcal{C}_k = {\rm vec}\{\bar{C}_{1,k},\cdots,\bar{C}_{r,k}\}$ and ${\mathcal{R}_k} = {\rm diag}\{\bar{R}_{1,k},\cdots,\bar{R}_{r,k}\}$. Then we have the augmented global measurement

Denoting

where $\bar{S}^[1]_{k}=\frac{1}{r}\sum_{i=1}^rC_{i,k}^\mathrm{T}R_{i,k}^{-1}C_{i,k}$ and $\bar{S}^[2]_{k}=\frac{1}{r}\sum_{i=1}^r\sigma_{i,k}^{-2}H_k^\mathrm{T}H_k$, we can obtain a local filter for each sensor node,which has the same performance as the centralized filter presented in Theorem 1 by using the augmented global measurement (28).

Suppose every node of the network applies the following filter

where $P_{i,0|-1}=rP_{0|-1}$,then the local and central state estimates for all nodes are the same,i.e.,$\hat{{\pmb x}}_{i,k|k}=\hat{{\pmb x}}_{k|k}$.

Proof. The proof is a direct combination of Corollary 1 and the result in ^[11].

In Corollary 2,the average measurement ${\pmb z}_k$ and the average inverse covariance matrix $\bar S_k^{^{\left[ 1 \right]}}$ can be computed approximately on each node using the following high-pass consensus filter proposed in ^[11]

It has been verified that,if the network is strongly connected, ${\pmb \xi}_i\rightarrow \frac{1}{r}\sum_{i=1}^r{\pmb u}_i$ as $t\rightarrow\infty$ for any $i \in \mathcal{V}$. $\bar S_k^{^{\left[ 2 \right]}}$ can be approximated on each node by $\bar S_{i,k}^{^{\left[ 2 \right]}}=\sigma_{i,k}^{-2}\hat{H}_{i,k}^\mathrm{T}\hat{H}_{i,k}$, where $\hat{H}_{i,k} = {\rm sgn}(\hat{{\pmb x}}_{i,k|k}^-)^\mathrm{T}.$p>

The distributed $H_{\infty}$ filtering can be summarized as the following Algorithm 1.

Download: JPG larger image

In this section,we will illustrate effectiveness of the method proposed in Section Ⅳ.

Without loss of generality,consider a sensor network with 6 nodes as shown in Fig. 1,whose topology is represented by an undirected graph $\mathcal{G}=(\mathcal {V},\mathcal {E}, \mathcal {A})$ with the set of nodes $\mathcal{V}=\{1,2,3,4,5,6\}$,the set of edges $\mathcal{E}=\{(1,1),(1,2),(1,3),(2,2),(2,3),(2,4),(2,5),(3,3),(3,5),(3,6),$ $(4,4),(4,5),(5,5),(5,6),(6,6)\}$,and the adjacency matrix

Fig. 1

	Download: JPG larger image
Fig. 1. Topological structure of the sensor network.

Here,we attempt to estimate a 10-sparse signal sequence $\{{\pmb x}_k\}$ in $\mathbf{R}^{256}$ and assume the sequence has a constant support set. The dynamics of $\{{\pmb x}_k\}$ can be described by the following linear system

where $\{{\pmb w}_{k}(i)\}$ is a zero-mean Gaussian white noise sequence with standard deviation 0.01. For system (1),we have $A = I_{256}$ and $Q_k = 10^{-4} \times I_{256}$. Both the index $i \in {\rm Supp}\{{\pmb x}_k\}$ and the value of ${\pmb x}_k(i)$ are unknown. The initial sparse signal ${\pmb x}_0$ is generated by creating the support set as $i\sim {\rm U}_\mathrm{int}[1, 256]$ and setting the value as ${\pmb x}_0(i)\sim {\rm N}(0,1)$. For sensor node $i$,the sensing matrix $C_{i,k} \in \mathbf{R}^{12\times 256}$ consists of entries sampled according to ${\rm N}(0,\frac{1}{72})$. Then,we have the global measurement ${\pmb y}_{\mathrm{cent}} = {\rm vec}\{{\pmb y}_1,\cdots,{\pmb y}_6\} \in \mathbf{R}^{72}$

Table I

Table I NMSEs WHEN $k=80$ Algorithms KF CSKF-1 Node 1 Node 2 Node 3 Node 4 Node 5 Node 6 Errors 0.7800 0.0293 0.0109 0.0074 0.0078 0.0087 0.0075 0.0086

Table I NMSEs WHEN $k=80$

and the global measurement matrix $C_{\mathrm{cent}} = {\rm vec}\{C_1,\cdots,\\C_6\} \in \mathbf{R}^{72\times 256}$. This type of matrix like $C_{\mathrm{cent}}$ has been shown to satisfy the restricted isometry property with overwhelming probability. The measurement noise ${\pmb v}_{i,k} \sim {\rm N}(0,R_{i,k})$, where $R_{i,k}=10^{-4} \times I_{12}$,and denote $R_{\mathrm{cent}} = {\rm diag}\{R_1,\cdots,R_6\}\in \mathbf{R}^{72\times 72}$. Set filter initial states as $\hat{{\pmb x}}_{i,0|-1}=0~(i = 1,\cdots,6)$,and filter parameters as $\beta=0.2$,$\lambda=1$,$\gamma=10$,$\sigma_{i,k} = 10+0.5k$,and $N=80$.

Now we are ready to design the distributed $H_{\infty}$ filters for system (38). The centralized Kalman filtering and CSKF-1 algorithm proposed in ^[6] are implemented respectively using ${\pmb y}_{\mathrm{cent}}$,$C_{\mathrm{cent}}$ and $R_{\mathrm{cent}}$,and we will compare their performances with that of Algorithm 1. The following normalized mean squared error (NMSE) is employed to evaluate the performance,

The results are presented in Figs. 2$\sim$4. In Fig. 2,the first figure shows the actual signal ${\pmb x}_{80}$,the second figure gives $\hat{{\pmb x}}_{80}$ from centralized Kalman filtering,the third figure shows $\hat{{\pmb x}}_{80}$ from centralized CSKF-1 algorithm. Fig. 3 gives the estimates of ${\pmb x}_{80}$ from sensor nodes using Algorithm 1. It is obvious that all nodes give satisfactory estimates of the actual sparse signal. Fig. 4 presents the normalized mean squared errors,which indicates that all node filters are stable and reach a consensus on the estimates of sparse signals. Moreover,there exist smaller steady errors by using Algorithm 1 than centralized Kalman filtering and CSKF-1 algorithm, and the errors when $k=80$ are presented in Table I. All of these results demonstrate effectiveness of the distributed filter presented in this paper.

Fig. 2

	Download: JPG larger image
Fig. 2. Actual signal ${\pmb x}_{80}$ and its estimates $\hat{{\pmb x}}_{80}$ using centralized KF and CSKF-1.

In this paper,the problem of distributed sparse signal estimation in sensor networks has been considered. The pseudo-measurement equation was interpreted as a 1-norm regularization term in the classic $H_{\infty}$ performance index,and the $\ell_1$-regularized $H_{\infty}$ filtering method was proposed. By means of a high-pass consensus filter,the distributed $H_{\infty}$ filtering algorithm was established. A numerical example verified effectiveness of the proposed method. However, the algorithm presented in this paper is only suitable to deal with the time-varying sparse signal with an invariant or slowly changing support set,and the more general methods need to combine a support set estimator with the filter. The decomposition of sensing matrix needs further research as well. These problems will be discussed in future work.

Fig. 3

	Download: JPG larger image
Fig. 3. Estimates of ${\pmb x}_{80}$ using Algorithm 1.

Fig. 4

	Download: JPG larger image
Fig. 4. MSEs of centralized KF,CSKF-1 and Algorithm 1.