This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/JBHI.2014.2312374, IEEE Journal of Biomedical and Health Informatics 1

Energy-Efficient ECG Compression on Wireless Biosensors via Minimal Coherence Sensing and Weighted ℓ1 Minimization Reconstruction Jun Zhang, Zhenghui Gu∗ , Member, IEEE, Zhu Liang Yu, Member, IEEE, and Yuanqing Li, Member, IEEE

Abstract—Low energy consumption is crucial for body area networks (BANs). In BAN-enabled electrocardiography (ECG) monitoring, the continuous monitoring entails the need of the sensor nodes to transmit a huge data to the sink node, which leads to excessive energy consumption. To reduce airtime over energy-hungry wireless links, this paper presents an energyefficient compressed sensing (CS)-based approach for on-node ECG compression. At first, an algorithm called minimal mutual coherence pursuit is proposed to construct sparse binary measurement matrices, which can be used to encode the ECG signals with superior performance and extremely low complexity. Secondly, in order to minimize the data rate required for faithful reconstruction, a weighted ℓ1 minimization model is derived by exploring the multi-source prior knowledge in wavelet domain. Experimental results on MIT-BIH arrhythmia database reveals that the proposed approach can obtain higher compression ratio than the state-of-the-art CS-based methods. Together with its low encoding complexity, our approach can achieve significant energy saving in both encoding process and wireless transmission. Index Terms—ECG telemonitoring, compressed sensing (CS), weighted ℓ1 minimization, incoherence.

I. I NTRODUCTION ITH the rapid advancement of body area network (BAN), low cost BAN-enabled electrocardiography (ECG) monitor has become important equipment for healthcare service providers and users [1]–[3]. However, in common BAN-enabled ECG monitors, the sensor nodes fall short of energy efficiency due to the huge data obtained from continuous monitoring and the energy-hungry wireless links [4]–[7]. It is desirable to improve the energy efficiency by reducing the in-network data rate and shortening the airtime over wireless links. Although conventional compression techniques, such as wavelet-based algorithms [8]–[10], have a high data compression ratio, their energy cost in encoding process will offset,

W

Manuscript received XXXX, 2013; revised xxx, 2013. This work was supported in part by the National Natural Science Foundation of China under grants 61105121, 61175114, 91120305, the Natural Science Foundation of Guangdong under grants S2012020010945, the Fundamental Research Funds for the Central Universities, SCUT under grant 2013ZZ0040, the National High-tech R&D Program of China (863 Program) under grant 2012AA011601, the High Level Talent Project of Guangdong Province and 2013KJCX0009. Asterisk indicates corresponding author. Zhenghui Gu(∗ ) is with the College of Automation Science and Engineering, South China University of Technology, Guangzhou 510640, China. Email: [email protected]. J. Zhang is with the College of Information Engineering, Guangdong University of Technology, Guangzhou 510006, China. Z. L. Yu and Y. Li are with the College of Automation Science and Engineering, South China University of Technology, Guangzhou 510640, China.

or even exceed the savings. Therefore, new techniques with low encoding complexity for on-node compression are needed. Compressed sensing (CS) [11], [12] is an emerging discipline which has attracted increasing attention thanks to its ability to guarantee data compression during the signal acquisition phase. In particular, CS has the potential to meet the demand of low encoding complexity by designing the measurement matrix carefully. Recently, the idea of CS has been applied in real-time ECG compression on resource-constrained sensors [13]–[15]. Typically, Mamaghanian et al. [13], [14] proposed a measurement matrix design technique named random sparse binary matrix (R-SBM) method to encode the ECG signals on the representative Shimmer mote. Each column of RSBM contains only zero entries besides d entries of 1’s with random locations, where the value of d is small. Using the RSBM, the encoding process of CS can be realized through a few addition operations. Thus, it can greatly reduce the energy cost. Decoding by the ℓ1 minimization algorithm [13] or model-based algorithms [14], Mamaghanian’s work shows that the CS-based method extends the sensor node lifetime through reducing the data rate. Further, Zhang et al. proposed a CS-based algorithm [15]–[18], named Block Sparse Bayesian Learning (BSBL), to reconstruct fetal ECG and EEG signals. Exploiting the intra-block correlation of the signals, the BSBL has shown its outstanding performance on both fetal ECG and EEG reconstruction. Particularly, the BSBL allows the use of an R-SBM with only two nonzero entries in each column to compress signals, which can greatly reduce energy cost in encoding process. Generally, the performance of CS relies heavily on incoherence measurement and reconstruction algorithm [19]. It is easy to confirm that when the number of entries of 1’s in each columns is fewer, the R-SBM can encode ECG signals with lower energy consumption [13], [15]. However, the incoherence property of the R-SBM will deteriorate sharply at the same time [13]. Obviously, high coherence in measurement matrix means more measurements are needed to finally reconstruct the ECG signal, which will increase the energy consumption in wireless data transmission. Hence, it is important to minimize the coherence of the SBM when the number of entries of 1’s in each column is few. On the other hand, to decode ECG signals, the ℓ1 minimization [19], [20] is the most frequently used algorithm since it is a convex problem and admits efficient solution via linear programming techniques. Essentially, the ℓ1 minimization algorithm utilizes only the sparsity prior on the property of

2168-2194 (c) 2013 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.

This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/JBHI.2014.2312374, IEEE Journal of Biomedical and Health Informatics 2

signals to guarantee the exact reconstruction. Nevertheless, other than the standard sparsity prior, incorporating additional knowledge of signals in the reconstruction algorithm is an efficient way to achieve higher decoding quality at a given number of measurements, or equivalently, a lower number of measurements for a given decoding quality [21], [22]. Besides the high sparsity nature of ECG signal in wavelet domain, we also notice that the approximation coefficients of an ECG signal are often nonzero, and the detail coefficients decay exponentially fast with increasing resolution level. If these anticipated information can be explored and exploited properly in the reconstruction algorithm, one may achieve accurate decoding with substantially lower number of measurements. In this paper, we aim to develop an energy-efficient CS method for on-node ECG compression. Our main contributions lie in the following two aspects. (1) A minimal mutualcoherence pursuit (MMCP) algorithm is proposed. In contrast to the random fashion in R-SBM, the proposed algorithm optimizes the locations of entries of 1’s to construct an SBM which achieves minimal mutual coherence. (2) We have deduced the decay factors of wavelet coefficients between neighboring resolution levels, and present a weighted ℓ1 minimization (WLM) algorithm to improve the decoding performance by embedding the multi-source knowledge within penalty weights. Finally, experimental results on MIT-BIH arrhythmia database showed that the proposed method with ultra low encoding complexity allows the number of measurements to be reduced greatly while maintaining the reconstruction accuracy. The rest of the paper is organized as follows. Section II introduces basic CS, the minimal mutual-coherence sparse binary sensing and the WLM algorithm. In Section III, the experimental results are presented and compared with the results from the state-of-the-art CS-based compression methods. Conclusions are given in Section IV. II. M ATERIAL A ND M ETHODS A. Compressed Sensing and Sparse Recovery CS, as a new compression paradigm, relies on three main requirements: sparsity representation, incoherence measurement and nonlinear reconstruction, which pertain to the signals of interest, the encoding modality and the decoding method respectively. Let x ∈ RN be a real-valued N -dimensional ECG signal vector, which we expand in an orthonormal wavelet basis Ψ = [ψ 1 ψ 2 · · · ψ N ] as follows: x = Ψθ

(1)

where θ is the coefficient vector of x. If most entries of vector θ have negligible amplitude, the ECG signal has sparse representations in the wavelet domain. The basic CS model can be expressed as y = Φx , Aθ M ×N

(2)

where Φ ∈ R (M ≪ N ) is a designed measurement (or sensing) matrix, and A = ΦΨ. y ∈ RM is the compressed signal vector. In our study x is a segment from an ECG recording, and y is the compressed data that will be transmitted via a BAN to a remote terminal.

In order to improve the energy efficiency, we need to make the data to be transmitted as less as possible. There are two key problems in CS-based data compression. First, for encoding, how should we design a measurement matrix Φ with minimal amount of measurements to ensure that data y preserves all the information of vector x? Second, for decoding, how can we reconstruct the original x from the compressed data y? For the first question, the measurement matrix Φ must obey the restricted isometry property (RIP) or the null space property (NSP) [19]. Among them, the RIP requires that all column submatrices of a certain size of the sensing matrix are well-conditioned. However, it is quite difficult to analyze these properties for deterministic matrices with minimal amount of measurements. Hence, mutual coherence (MC) [23], as an alternative way, is often used to characterize the measurement matrix, which can be defined as M C = M C(Φ) =

max

1≤k,j≤N,k̸=j

|G(k, j)|

(3)

¯TΦ ¯ and Φ ¯ is a matrix whose column is the where G , Φ normalized version of the column of Φ with unit ℓ2 -norm. Because smaller MC usually means that the measurement matrix Φ can preserve the information in the vector x with less amount of measurements, it is desirable that the MC is as small as possible. In general CS applications, random matrices such as Gaussian iid matrices or Bernoulli matrices are good choices. Introducing randomness allows to show near-optimal conditions on the number of measurements [24], [25]. In terms of the sparsity of signal vector, it allows the recovery using the below ℓ1 -minimization method (Basis Pursuit).



min θˆ s.t. y = Aθˆ (4) ˆ θ

1

In addition, Bayesian CS [26] and greedy approaches such as orthogonal matching pursuit (OMP) [23], Iterative Hard Thresholding (IHT) [27], [28] and Gradient Pursuits (GP) [29], [30] are the state-of-the-art algorithms for signal recovery. Recently, BSBL [15], [18] has been proposed as an effective method for CS recovery problem. B. Minimal Mutual Coherence Sparse Binary Sensing Although the random matrices have many theoretical benefits, they have the drawbacks of high complexity and large storage in practical implementation. It is unrealistic to use them in the BAN-based ECG monitoring. Mamaghanian et al. [13] showed that the R-SBM can greatly reduce on-node energy and space consumption in the application. However, constituted only by 0’s and 1’s, the R-SBM is inferior to the random matrices in terms of the capability of sensing the signal with minimal amount of measurements, especially when the number of entries of 1’s is extremely low. Therefore, this subsection proposes an MMCP algorithm to construct the SBM. The main difference between our method and R-SBM is that contrary to the random fashion in R-SBM, our method optimizes the locations of entries of 1’s to pursuit minimal mutual coherence.

2168-2194 (c) 2013 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.

This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/JBHI.2014.2312374, IEEE Journal of Biomedical and Health Informatics 3

1

0.9

R−SBM RGM MMC−SBM

0.8

Mutual coherence

0.7

0.6

0.5

0.4

0.3

0.2

0.1 0

2

4

6

8

10

12

14

16

18

20

22

24 25

Number of nonzero entries in each column (d)

Fig. 1.

A scene of sensing matrix construction

Fig. 2.

Suppose we construct an M × N SBM Φ in which each column contains d entries of 1’s. Denote the i-th column and the j-th row of Φ as ϕi and γ j respectively. Let Ωi be the location set of entries of 1’s in ϕi . In the beginning we initialize Φ to a zero matrix, thus Ωi (1 ≤ i ≤ N ) is an empty set. The main practices of the MMCP algorithm can be summarized in two rules: 1) increasing the entries of 1’s in Φ column by column and one by one; 2) determining the location of current 1’s with the principle of minimal MC increase. According to the definition of MC, the increase in MC results from the increase in the number of superposition of 1’s between two columns. Therefore, in order to insert a 1 into a column, we should try to choose its location that can avoid the increase in the number of superposed 1’s between the current column and the previous columns. However, except for the two columns that have the most number of superposed 1’s, the increase in the number of superposed 1’s between other two columns will not cause the increase of MC immediately. Particularly, there always exist superposed 1’s between two columns in the final constructed Φ because M < N . Thus, a single superposition of 1 between any two columns will not affect the MC of the final Φ. Moreover, if more than one superposition is inevitable, we would choose to increase the number of superposed 1’s between the current column and the column that has the least superposed 1’s with the current column at the moment. In the location set of such column, the location that has superposed 1’s with the least previous columns is chosen. Based on the above analysis, to determine the k-th location (1 ≤ k ≤ d) of 1 in ϕi , we compute the following items: Ωi ∩ Ωi′

′

f or 1 ≤ i ≤ i − 1

(5)

Mutual coherence of sensing matrices versus d

At this moment, both Λ1 and Λ2 are nonempty set. Therefore, we should try to avoid choosing the locations 4, 8, 9 and 11. However, when all of the sets Λi′ have been checked, two situations can occur: 1)a nonempty position set, denoted as ∆, exists; 2) each location is included in a Ωi′ whose Λi′ is not an empty set. For the first situation, we choose the location j for the current 1’s, where j ∈ ∆ and the sum of γ j is minimal. Otherwise, we choose the location j, where ′ s = mini′ {|Λi′ |, 1 ≤ i ≤ i − 1}, j ∈ Ωs and the sum of γ j is minimal. Algorithm 1: Minimal MC Pursuit Input: The size of matrix Φ-(M, N ); The number of nonzero elements in each column d Initialization: Zero matrix Φ ∈ RM ×N ; Empty sets Ωi , for 1 ≤ i ≤ N Output: The sensing matrix Φ for i = 1 to N do for k = 1 to d do ∆ = {1, 2, ..., M } \ Ωi ; ′ for i = 1 to i − 1 do if Ωi ∩ Ωi′ ̸= ∅ then ∆ = ∆ \ Ωi′ ; if ∆ ̸= ∅ then Ωi ←− j, where j ∈ ∆ and the sum of γ j is minimal; Φ(j, i) = 1; else Ωi ←− j, where ′ s = min {|Λi′ |, 1 ≤ i ≤ i − 1}, j ∈ Ωs and ′ i

the sum of γ j is minimal; Φ(j, i) = 1;

′

Denote Λi′ , 1 ≤ i ≤ i − 1, as the results of (5). We check the sets Λi′ one by one. If a certain Λi′ is not an empty set, we should not choose the element of Ωi′ as the location of the current 1’s in ϕi in principle, for it will produce more than one 1’s superposition between ϕi and ϕi′ . Assume M = 12 and d = 3, a scene of determining the third location of 1’s in ϕi is depicted in Fig.1. For simplicity of the description, we omit all the columns except the first, second and i-th columns.

return Φ; The proposed algorithm can be outlined as Algorithm 1. To validate the effectiveness of the proposed MMCP algorithm, we herein compare the measurement matrices constructed by

2168-2194 (c) 2013 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.

This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/JBHI.2014.2312374, IEEE Journal of Biomedical and Health Informatics 4

300

2

200

1

100

0 0 −0.4

−0.3

−0.2

−0.1

−1

0

0.1

0.2

0.3

0.4

0.1

0.2

0.3

0.4

0.1

0.2

0.3

0.4

H2 subband 100

−2 0

500

1000

1500

2000

2500

3000 50

6 4

0 −0.4

2

40

−0.3

−0.2

−0.1

0

H1 subband 30

0 20

−2

10

−4 0

500

1000

1500

2000

2500

3000

0 −0.4

−0.3

−0.2

−0.1

0

H subband 0

(a)

(b)

Fig. 3. (a) From top to bottom, the original segment of record no. 124 MIT-BIH arrhythmia data and the corresponding wavelet coefficients of 3-level wavelet transform. (b) From top to bottom, the histograms of wavelet coefficients of the H2 , H1 , H0 high frequency subbands respectively.

the MMCP algorithm, named MMC-SBM, with the R-SBM in terms of MC. The size of the matrices is 256 × 512 in our experiment. Fig.2 plots the MC versus the number of nonzero elements d for the two classes of matrices. As a baseline, the MC corresponding to random Gaussian matrix (RGM) is also reported, where the elements of the matrix are normal distributed random numbers. As shown in Fig.2, the MMCSBM is always superior to the R-SBM, even better than the RGM in the cases of 4 ≤ d ≤ 8. It is easy to verify that the MMC-SBMs have the minimal MC for 1 ≤ d ≤ 8, because there are at most a single superposition of 1 between any two columns. Although the MC of the MMC-SBM is higher than that of the RGM when d ≥ 9 due to its binary feature and the decrease of position freedom, it is worth noting that the optimal construction of extremely sparse binary matrix is what we care about in BAN-enabled ECG telemonitoring for its ultra-low energy and space consumption. Remark 1: In the worst case, the computational complexity of Algorithm 1 scales as O(N 2 log2N ). The complexity of RSBM construction is O(M N log2M ). Although the complexity of MMC-SBM construction is higher than that of R-SBM, it can be considered a tolerable defect since the measurement matrices are always obtained offline. In the BAN-enabled ECG monitoring, the most important issues of a measurement matrix are its sensing performance and computational complexity required in signal encoding process.

C. Weighted ℓ1 Minimization Method To reconstruct the original ECG signal x, we herein propose a weighted ℓ1 minimization method (WLM). In most practical situations, we can observe only the noisy data y = Aθ + n, where n represents the additive measurement noise. We assume that the components of n are independent and identically distributed (i.i.d.) Gaussian variables with unknown variance σn2 , and the entry θˆi of θˆ is independent and has a Laplacian

distribution with standard deviation σi , i.e., √ 2 θˆi 1 p(θˆi ) = √ . exp(− ) σi 2σi

(6)

ˆ A), the maximum a By maximizing the probability p(θ|y, posteriori (MAP) estimation of θ under the prior p(θˆi ) can be given as [31] 

2 √ 

ˆ − Aθ ∑ 2 θˆi   y 2 + θˆM AP = arg min   (7) 2σn2 σi ˆ θ i Further, equation (7) can be expressed as

2

1



θˆM AP = arg min y − Aθˆ + λ Wθˆ 2 2 1 ˆ θ

(8)

where W = diag( σ11 , σ12 , ..., σ1N ), λ is a tuning parameter. Hence, Basis Pursuit in (4) can be interpreted as the MAP estimation under a Laplacian prior on θˆ and the hypothesis that the entries of θˆ have equal standard deviations. However, for most ECG signals, there are at least two facts in our application that make the above hypotheses unsatisfied: 1) since the approximation coefficients of an ECG signal in wavelet domain are often nonzero, the corresponding part of θˆ does not satisfy the sparsity hypothesis and can’t be modeled by the Laplacian distribution; 2) the detail coefficients decay exponentially fast with increasing resolution level, therefore the hypothesis that the entries of θˆ have equal standard deviations does not reflect this reality. To illustrate the two arguments, as an example, Fig.3 plots the wavelet coefficients of the 3-level wavelet transform for record no.124 MIT-BIH arrhythmia data. At the same time, the corresponding histograms of the wavelet coefficients in the three detail subbands are simultaneously reported. Clearly, the approximation coefficients (range [1, 375]) are nonzero and the standard deviations of the wavelet coefficients in the three detail subbands stably decay with increasing resolution level.

2168-2194 (c) 2013 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.

This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/JBHI.2014.2312374, IEEE Journal of Biomedical and Health Informatics 5

In the proposed WLM, we construct a N × N matrix W to incorporate the aforementioned multi-source prior. Suppose an ECG signal x is decomposed into J layers. Thus θ can be partitioned into J + 1 groups, i.e., θ = [θ TL0 , θ TH0 , ..., θ THJ−1 ]T

(9)

where θ L0 consists the approximation coefficients, and θ H0 , θ H1 , ..., θ HJ−1 represent the detail coefficients at J levels respectively. First, because almost all the entries of θ L0 are nonzero, it is not necessary to employ a sparsity constraint on θˆL0 . Therefore, the corresponding diagonal elements of W for θˆL0 are set at zero. Second, we determine the other elements of matrix W by studying the decay characteristic of detail coefficients between different levels. For special choice of wavelet ψ, we assume that the wavelet ψ has q vanishing moments [32], i.e., ∫ xn ψ(t)dt = 0 f or n = 0, 1, ..., q − 1 (10) R

The k-th entry of θ Hj is denoted as θHj,k , thus θHj,k =< x, ψj,k >

(11)

where ψj,k (t) = 2j/2 ψ(2j t − k), j, k ∈ Z. Note that ψj,k lives on an interval Ij,k of length proportional to 2−j . Hence, if signal x has s ≤ q continuous derivatives on Ij,k , it can be approximated with a Taylor polynomial [33] x(t) = Ps−1 (t) + O((2−j )s )

(12)

where Ps−1 (t) is a (s-1)-order polynomial. According to (10), we have ∫ θHj,k =< x, ψj,k >= x(t)ψj,k (t)dt R ∫ (13) = Ps−1 (t)ψj,k (t)dt + O((2−j )s ) R

= O((2−j )s ) Obviously, the magnitude of the detail coefficients is governed by the resolution level j and the local smoothness s of the underlying signal x. Imposing the prior belief that detail coefficient has a 2−s rate of decay with increasing resolution level, the reconstruction of θ can be reformulated by the convex optimization problem as follows arg min ˆ θ

J−1

2

T ∑ 1



2−(J−1−j)s θˆHj

y − Aθˆ + λ 2 2 1 j=0

(14)

where s is a tuning parameter that reflects the smoothness of signal x. It is worth noting that model (14) aims to incorporate the level-dependent prior to reconstruct vector θ, but simplify the influence of the different local smoothness of the underlying signal on the detail coefficients. In our experiments, considering that ECG signal can be regarded as a piece-wise polynomial function [34], which has 1 derivative at most of positions, we chose s = 1. Experimental results verify the validity of this selection. Remark 2: In contrast to Basis Pursuit, the main advantage of model (14) is the incorporation of the multi-source prior in

wavelet domain, including the positions of nonzero coefficients and the decay relationship of wavelet coefficients between neighboring resolution levels, which will allow the number of measurements to be significantly reduced without leading to ambiguity. Remark 3: The method (14) has similar computational complexity as the Basis Pursuit. However, it is noteworthy that the decoding process is often performed on a computer equipped with high-performance processor, the computational complexity of the reconstruction method is not a very important issue in most of the BAN-enabled ECG monitoring applications. III. R ESULTS A ND D ISCUSSION In this section, experiments were carried out to verify the performance of the proposed method. We have chosen all the 48 MIT-BIH arrhythmia data as test signals. Each data file includes two-channel ambulatory ECG recordings, and one channel data instance was chosen for the evaluation of the algorithm performance. To measure the reconstruction quality, we employ the percentage root-mean-square difference (PRD) to quantify the error percentage between the original x and the reconstructed x ˆ: ∥x − ˆx∥2 P RD = × 100 (%) (15) ∥x∥2 For ECG compression, Zigel et al [35] had classified the different values of PRD based on the signal quality perceived by specialists. In accordance with their results, PRD value below 9% is classified as ”good” reconstruction quality in this paper. A. Experimental Setup In our experiments ECG vector x is extracted from the original signal at the window size N = 512. To obtain M samples, we are interested in constructing a measurement matrix that satisfies the minimal MC, but has as less number of entries of 1’s in each column as possible. Hence, the proposed MMCP algorithm was applicable to the determination of this matrix. By increasing the number of entries of 1’s in each column from d = 1, the MMC-SBMs were constructed in a sequential way. In the initial phase, the MC of the constructed matrix will decrease sharply with the increase of d. In our experiments, we found that the matrix with minimal MC can be always obtained in a very small d (96 × 512, 128 × 512 and 160 × 512 matrices can be constructed with d = 4, 5 and 6 respectively). We chose these MMC-SBM as the measurement matrix to compress (encode) the ECG signals, and the WLM algorithm was used in the remote terminal to recover (decode) the original ECG signals. To evaluate the performance of MMC-SBM, we compared the recovery results of the WLM algorithm when the measurement matrix was the MMC-SBM, R-SBM and RGM respectively. At the same time, several state-of-the-art CS algorithms have been chosen for performance comparison, including Orthogonal Matching Pursuit (OMP) [23], Iterative hard thresholding (IHT) [27], [28], Gradient Pursuit (GP) [29], [30],

2168-2194 (c) 2013 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.

This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/JBHI.2014.2312374, IEEE Journal of Biomedical and Health Informatics 6

120 OMP IHT GP BCS BP BSBL WLM(RGM) WLM(R−SBM) WLM(MMC−SBM)

100

Output PRD (averaged over all records)

Bayesian CS (BCS) [26], Basis Pursuit (BP) [20] and Block Sparse Bayesian Learning (BSBL) [15], [18]. Among them, OMP, IHT, GP, BCS and BP are the representatives of the greedy algorithms, Bayesian compressive sensing or convex relaxation algorithms. The BSBL proposed recently by Zhang et al. is the most outstanding CS method for Fetal ECG(FECG) and EEG reconstruction. In our comparison experiments, these algorithms used the R-SBM with d = 12 as the measurement matrix that was suggested in Mamaghanian’s work [13]. For the OMP and BP implementations, we used the solvers SolveOM P and SolveBP , respectively, from the SparseLab toolbox [36]; for the IHT and GP implementations, we used the solvers Hard 10 M term and greed gp, respectively, from the Sparsify toolbox [37]; for the BCS implementation, we used the solver BCS f ast rvm from the BCS toolbox [38]; for the BSBL, we used the solver BSBL BO [15]. We set Ψ to be the orthonormal basis of Daubechies-6 wavelet because it can yield very sparse representation of ECG and has 6 vanishing moments.

80

60

40

20

G VG 96

128

160

192

224

256

288

320

Number of Measurments (M)

Fig. 4. Output PRD averaged over all MIT-BIH ECG records for different number of measurements M for different algorithms. WLM(RGM), WLM(RSBM) and WLM(MMC-SBM) represent the results of WLM algorithm when the measurement matrices are RGM, R-SBM and MMC-SBM respectively.

B. Average PRD and the Probability of ’Good’ Reconstruction 120

100

80

PRD

For all the 48 MIT-BIH arrhythmia records, the experimental result is shown in Fig.4, where each point is the average PRD of a specified number of measurements. At the same time, Table I reports the probability of ’good’ reconstruction quality at different situations. As it is shown in Fig. 4 and Table I, the WLM related algorithms achieved the best performance in terms of averaged output PRD and the probability of ’good’ reconstruction quality. Moreover, the performance of MMC-SBM has been validated. Among the three different measurement matrices (MMC-SBM, R-SBM and RGM), the MMC-SBM can improve the reconstruction accuracy of WLM algorithm when the number of measurements is very small. In these experiments, we can also observe that the performance improvement caused by the decoding method, WLM, is more significant than that caused by the measurement matrix MMCSBM. Considering the measurement matrix, it is noteworthy that the RGMs approach the near-optimal sensing performance in theory [24], [25]. However, they are hard to be used in applications like long-term ECG monitoring due to their high computational complexity in signal encoding. Therefore, taking into account both the superior sensing performance and low computational complexity, the MMC-SBM is a very good measurement matrix for long-term ECG monitoring. It is also shown that if the MMC-SBM is chosen as the measurement matrix, the averaged output PRD of the WLM algorithm is less than 9% for the number of measurements M = 128, and it can achieve 87.5% of ’good’ reconstruction quality with M = 160. As a comparison, most other algorithms under this condition can not recover any records in MIT-BIH arrhythmia database with ’good’ reconstruction quality due to the very small number of measurements. To observe the variance across the individual records, Fig. 5 shows the box plots for these algorithms when the number of measurements M = 128. On each box, the central mark is the median, the edges of the box are the 25th and 75th percentiles, and the whiskers extend to the most extreme data points.

60

40

20

0 OMP

IHT

GP

BCS

BP

BSBL

WLM(G)

WLM(R)

WLM(M)

Fig. 5. Box plots for all database records for OMP, IHT, GP, BCS, BP, BSBL and WLM respectively, when the number of measurements M = 128. WLM(G), WLM(R) and WLM(M) are the abbreviations of WLM(RGM), WLM(R-SBM) and WLM(MMC-SBM) respectively.

Obviously the WLM algorithm adopted the MMC-SBM as TABLE I P ROBABILITY OF R ECONSTRUCTION FOR ’G OOD ’ Q UALITY UNDER D IFFERENT N UMBER OF M EASUREMENTS M (%).

M

96 128 160 192 224 256 288 320

OMP

IHT

GP

BCS

BP

BSBL

0 0 0 72.9 79.2 87.5 91.7 95.8

0 0 0 39.6 64.6 87.5 91.7 91.7

0 0 0 16.7 85.4 91.7 91.7 93.8

0 0 0 87.5 91.7 91.7 93.8 97.9

0 0 12.5 37.5 75.0 91.7 93.8 97.9

2.1 18.8 72.9 89.6 95.8 100 100 100

WLM. RGM

R-SBM

MMC-SBM

20.8 54.2 83.3 93.8 95.8 97.9 100 100

18.8 52.1 81.3 93.8 95.8 97.9 100 100

31.3 70.8 87.5 93.8 95.8 97.9 97.9 100

2168-2194 (c) 2013 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.

This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/JBHI.2014.2312374, IEEE Journal of Biomedical and Health Informatics 7

1 0 −1 10 0 −1 20 0 −2 10 0 −1 10 0 −1 10 0 −1 0

2 0 −2

200

400

600

800

1000

200

400

600

800

1000

200

400

600

800

1000

200

400

600

800

1000

200

400

600

800

1000

200

400

600

800

1000

0

200

400

600

800

1000

0

200

400

600

800

1000

0

200

400

600

800

1000

0

200

400

600

800

1000

0

200

400

600

800

1000

0

200

400

600

800

1000

1 0 −1 1 0 −1 1 0 −1 1 0 −1 2 0 −2

(a) Recovered results by BSBL

(a) Recovered results by BSBL

1 0 −1 10

200

400

600

800

1000

200

400

600

800

1000

200

400

600

800

1000

0 −1 20 0 −2 10 0 −1 10

200

400

600

800

1000

200

400

600

800

1000

0 −1 10

2 0 −2 10 0 −1 10 0 −1 10 0 −1 10 0 −1 20

200

400

600

800

1000

200

400

600

800

1000

200

400

600

800

1000

200

400

600

800

1000

200

400

600

800

1000

200

400

600

800

1000

0 −2

0 −1 0

200

400

600

800

1000

0

(b) Recovered results by WLM

(b) Recovered results by WLM

Fig. 7. First 1000 time points of the six records in MIT-BIH arrhythmia dataset with ”worst” PRD for BSBL. (a) Recovered results of BSBL. (b) Recovered results of WLM. In each figure, the red curves represent the original signal and the blue curves are the recovered signals.

Fig. 8. First 1000 time points of the six records in MIT-BIH arrhythmia dataset with ”worst” PRD for WLM. (a) Recovered results of BSBL. (b) Recovered results of WLM. In each figure, the red curves represent the original signal and the blue curves are the recovered signals.

the measurement matrix outperforms the other algorithms. In our experiment, Record 220 and 231 are the 25th and 75th percentiles respectively for the WLM algorithm. Fig.6 displays the recovered segments for these algorithms. Comparing the results we can see all the algorithms, besides the BSBL and the WLM algorithms, can hardly recover these records at such few measurements. The BSBL and WLM algorithm recovered the segments, but the distortion of the results of WLM algorithm was smaller.

records of these ECG records by WLM algorithm. Similarly, the first 1000 time points of the six recovered records with the worst PRD values for WLM algorithm are plotted in fig.8(b), and fig.8(a) shows the recovered records of these ECG records by BSBL algorithm. In both figures, the ’red’ curves represent the original records and the ’blue’ curves are the recovered records. Obviously, the recovered records by the WLM algorithm almost exactly match the original records, thus the diagnostic features can be extracted very well. By comparison, the recovered results by the BSBL algorithm have the noticeable S-wave distortion and R amplitude distortion. Hence, in the digital CS paradigm, the WLM algorithm is strongly recommended for EGC data recovery at low measurement rate.

C. Application Perspective We further compare the performance of the WLM algorithm with the BSBL algorithm from an application perspective. In clinical applications, it is crucial for an ECG monitoring system if the diagnostic features, such as P-R interval, QRS interval, ST segment and R amplitude, are preserved in the recovered ECG signal. We used the measurement matrix of size 128 × 512 to encode the MIT-BIH arrhythmia dataset, and the WLM and BSBL algorithm are employed to recover the ECG signals because they materially outperformed the other CS algorithms in the previous experiments. Here, we study the recovered records with the worst PRD by one algorithm and see how well they are recovered by the other algorithm, and vice versa. In fig.7(a), we plot the first 1000 time points of the six recovered records with the worst PRD values for BSBL algorithm. Fig.7(b) shows the recovered

D. Energy Saving To improve the energy efficiency of a BAN-enabled ECG monitor, it is desirable to reduce the energy consumption in both encoding process and wireless transmission. Comparing with the existing CS-based methods, the proposed method enables low-complexity encoding. In our experiments, each column of the 128×512 MMC-SBM contains only 5 entries of 1’s, which saves about 58% energy consumption in encoding process compared to the R-SBMs with 12 entries of 1’s that were adopted by the other CS-based methods. In terms of the energy consumption in wireless transmission, the experimental

2168-2194 (c) 2013 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.

This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/JBHI.2014.2312374, IEEE Journal of Biomedical and Health Informatics 8

Original Segment

OMP

2

Original Segment

5

OMP

2

2

0

0

0

0

−2

−5

−2

−2

0

100

200

300

400

500

0

100

200

300

400

500

5

0

0

100

200

300

100

200

400

400

500

100

200

0

100

200

1

0

0

−1

300

400

400

500

300

400

500

300

400

500

400

500

−2 0

500

100

200

BP

300

400

500

0

100

200

BCS

2

300

GP 1

−1

0

500

BCS 2

300

IHT 2

−5

−2 0

0

GP

IHT 2

BP

2

2

0 0

0

0

−2

−2

−2

−2 −4 0

100

200

300

400

500

0

100

200

300

BSBL

400

0

500

100

WLM

2

200

300

400

500

0

2

0

0

0

0

−2

−2

−2

100

200

300

400

500

0

100

200

300

400

500

200

WLM

2

−2 0

100

BSBL

2

0

100

(a) Record 220

200

300

400

500

0

100

200

300

(b) Record 231

Fig. 6. Recovery results for (a) Record 220 and (b) Record 231. In both sub-figure, there are the original segments and the recovered results by OMP, IHT, GP, BCS, BP, BSBL and WLM, respectively.

results indicate that the proposed method can achieve satisfactory decoding quality for the MIT-BIH arrhythmia dataset with only 25% measurements, while the other CS-based methods need to increase the number of measurements to decode the ECG signals with the same quality level. For example, the BP algorithm and the BSBL algorithm require more than 43% and 32% measurements, respectively, to produce the equivalent average output PRD. Therefore, the proposed method can obtain higher compression ratio in the applications than the other CS-based methods, which result in a significant energy saving in wireless transmission.

ACKNOWLEDGMENT The authors would like to thank the associate editor and the anonymous reviewers for their valuable comments and suggestions to improve the quality of the paper.

PLACE PHOTO HERE

IV. C ONCLUSION In this paper, an energy-efficient CS method was developed for on-node ECG compression. It consists of the improvements in both encoder and decoder. Firstly, an MMCP algorithm was proposed and applied to construct an SBM that can encode ECG recordings with high sensing performance and ultralow energy consumption. Furthermore, an WLM algorithm is derived for further reducing the measurements needed for high fidelity decoding of ECG recordings. The proposed method realizes energy efficiecy for BAN-enabled ECG monitoring from the aspects of minimizing energy consumption in encoding procedure and shortening the airtime over wireless links by efficient data compression. To evaluate the performance of the proposed method, it was applied to all of the MITBIH arrhythmia databases. Experimental results showed that the proposed method performs better compared to several other CS-based algorithms. By applying the constructed 128 × 512 SBM, in which each column contains only 5 non-zero entries, to encode the ECG recordings and decoding by the WLM algorithm, the reconstructed signals of all the MITBIH arrhythmia recordings almost exactly matched the original recordings and preserved the diagnostic features. In the future work, we will investigate the theoretic analysis and the other applications of the proposed method.

PLACE PHOTO HERE

Jun ZHANG received the bachelors and masters degrees in computer science from the Xiangtan University, Xiangtan, China, in 2002 and 2005, respectively, and the doctoral degree in pattern recognition and intelligence system from the South China University of Technological, Guangzhou, China, in 2012. He joined the College of Information Engineering, Guangdong University of Technology, in 2005 as a lecturer. His research interests include the fields of compressive sensing and biomedical signal processing.

Zhenghui GU received the Ph.D. degree from Nanyang Technological University in 2003. From 2002 to 2008, she was with Institute for Infocomm Research, Singapore. She joined the College of Automation Science and Engineering, South China University of Technology, in 2009 as an associate professor. Her research interests include the fields of biomedical signal processing and pattern recognition.

Zhu Liang YU received his BSEE in 1995 and MSEE in 1998, both in electronic engineering from the Nanjing University of Aeronautics and Astronautics, China. He received his Ph.D. in 2006 PLACE from Nanyang Technological University, Singapore. PHOTO He joined Center for Signal Processing, Nanyang HERE Technological University from 2000 as a research engineer, then as a Group Leader. In 2008, he joined the College of Automation Science and Engineering, South China University of Technology, China. He is promoted to be a full Professor in 2009. His research interests include signal processing, machine learning, computer vision and applications in biomedical engineering.

2168-2194 (c) 2013 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.

This article has been accepted for publication in a future issue of this journal, but has not been fully edited. Content may change prior to final publication. Citation information: DOI 10.1109/JBHI.2014.2312374, IEEE Journal of Biomedical and Health Informatics 9

Yuanqing Li received the B.S. degree in applied mathematics from Wuhan University, Wuhan, China, in 1988, the M.S. degree in applied mathematics from South China Normal University, Guangzhou, PLACE China, in 1994, and the Ph.D. degree in control PHOTO theory and applications from South China University HERE of Technology, Guangzhou, China, in 1997. Since 1997, he has been with South China University of Technology, where he became a full professor in 2004. In 2002-04, he worked at the Laboratory for Advanced Brain Signal Processing, RIKEN Brain Science Institute, Saitama, Japan, as a researcher. In 2004-08, he worked at the Laboratory for Neural Signal Processing, Institute for Infocomm Research, Singapore, as a research scientist. His research interests include blind signal processing, sparse representation, machine learning, brain-computer interface, EEG and fMRI data analysis.

R EFERENCES [1] R. L. Ashok and D. P. Agrawal, “Next-generation wearable networks,” IEEE Computer Mag., vol. 36, no. 11, pp. 31–39, 2003. [2] M. A. Hanson, H. C. Powell, A. T. Barth, K. Ringgenberg, B. H. Calhoun, J. H. Aylor, and J. Lach, “Body area sensor networks: Challenges and opportunities,” IEEE Computer Mag., vol. 42, no. 1, pp. 58–65, 2009. [3] F. Rinc´on, J. Recas, N. Khaled, and D. Atienza, “Development and evaluation of multilead wavelet-based ecg delineation algorithms for embedded wireless sensor nodes,” IEEE Trans. Inf. Technol. Biomed., vol. 15, no. 6, pp. 854–863, 2011. [4] T. Ma, P. L. Shrestha, M. Hempel, D. Peng, H. Sharif, and H.-H. Chen, “Assurance of energy efficiency and data security for ECG transmission in basns,” IEEE Trans.Biomed. Eng., vol. 59, no. 4, pp. 1041–1048, 2012. [5] S. Lee, J. Kim, and M. Lee, “A real-time ECG data compression and transmission algorithm for an e-health device,” IEEE Trans.Biomed. Eng., vol. 58, no. 9, pp. 2448–2455, 2011. [6] A. M. Dixon, E. G. Allstot, D. Gangopadhyay, and D. J. Allstot, “Compressed sensing system considerations for ECG and EMG wireless biosensors,” IEEE Trans. Biomed. Circ. Sys., vol. 6, no. 2, pp. 156–166, 2012. [7] C. Lin, Y. Yang, and J. Wang, “A wearable sensor module with a neural-network-based activity classification algorithm for daily energy expenditure estimation.” IEEE Trans. Inf. Technol. Biomed., vol. 16, no. 5, pp. 991–998, 2012. [8] B. A. Rajoub, “An efficient coding algorithm for the compression of ECG signals using the wavelet transform,” IEEE Trans. Biomed. Eng., vol. 49, no. 4, pp. 355–362, 2002. [9] S.-G. Miaou and S.-N. Chao, “Wavelet-based lossy-to-lossless ecg compression in a unified vector quantization framework,” IEEE Trans. Biomed. Eng., vol. 52, no. 3, pp. 539–543, 2005. [10] C.-T. Ku, K.-C. Hung, T.-C. Wu, and H.-S. Wang, “Wavelet-based ECG data compression system with linear quality control scheme,” IEEE Trans. Biomed. Eng., vol. 57, no. 6, pp. 1399–1409, 2010. [11] D. L. Donoho, “Compressed sensing,” IEEE Trans. Inf. Theory, vol. 52, no. 4, pp. 1289–1306, 2006. [12] E. J. Cand`es, J. Romberg, and T. Tao, “Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information,” IEEE Trans. Inf. Theory, vol. 52, no. 2, pp. 489–509, 2006. [13] H. Mamaghanian, N. Khaled, D. Atienza, and P. Vandergheynst, “Compressed sensing for real-time energy-efficient ECG compression on wireless body sensor nodes,” IEEE Trans.Biomed. Eng., vol. 58, no. 9, pp. 2456–2466, 2011. [14] H. Mamaghanian, N. Khaled, D. Atienza, and P. Vandergheynst, “Structured sparsity models for compressively sensed electrocardiogram signals: A comparative study,” in Biomedical Circuits and Systems Conference (BioCAS), 2011 IEEE. IEEE, 2011, pp. 125–128. [15] Z. Zhang, T.-P. Jung, S. Makeig, and B. D. Rao, “Compressed sensing for energy-efficient wireless telemonitoring of noninvasive fetal ECG via block sparse bayesian learning,” IEEE Trans.Biomed. Eng., vol. 60, no. 2, pp. 300–309, 2013. [16] Z. Zhang and B. D. Rao, “Recovery of block sparse signals using the framework of block sparse bayesian learning,” in Acoustics, Speech and Signal Processing (ICASSP), 2012 IEEE International Conference on. IEEE, 2012, pp. 3345–3348. [17] Z. Zhang and B. D. Rao, “Sparse signal recovery with temporally correlated source vectors using sparse bayesian learning,” IEEE J. Sel. Topics Siginal Proc., vol. 5, no. 5, pp. 912–926, 2011.

[18] Z. Zhang, T. Jung, S. Makeig, and B. Rao, “Compressed sensing of EEG for wireless telemonitoring with low energy consumption and inexpensive hardware,” IEEE Trans.Biomed. Eng., vol. 60, no. 1, pp. 221–224, 2013. [19] E. J. Cand`es and M. B. Wakin, “An introduction to compressive sampling,” IEEE Signal Proc. Mag., vol. 25, no. 2, pp. 21–30, 2008. [20] S. S. Chen, D. L. Donoho, and M. A. Saunders, “Atomic decomposition by basis pursuit,” SIAM REV., vol. 43, no. 1, pp. 129–159, 2001. [21] C. Miosso, R. von Borries, M. Arg`aez, L. Velazquez, C. Quintero, and C. Potes, “Compressive sensing reconstruction with prior information by iteratively reweighted least-squares,” IEEE Trans. Signal Process., vol. 57, no. 6, pp. 2424–2431, 2009. [22] Y. Wang and W. Yin, “sparse signal reconstruction via iterative support detection,” SIAM J Imag. Sci., vol. 3, no. 3, pp. 462–491, 2010. [23] J. Tropp and A. Gilbert, “Signal recovery from random measurements via orthogonal matching pursuit,” IEEE Trans. Inf. Theory, vol. 53, no. 12, pp. 4655–4666, 2007. [24] H. Rauhut, “Compressive sensing and structured random matrices,” Theoretical foundations and numerical methods for sparse recovery, vol. 9, pp. 1–92, 2010. [25] E. J. Candes and Y. Plan, “A probabilistic and ripless theory of compressed sensing,” IEEE Trans. Inf. Theory, vol. 57, no. 11, pp. 7235– 7254, 2011. [26] S. Ji, Y. Xue, and L. Carin, “Bayesian compressive sensing,” IEEE Trans. Signal Proc., vol. 56, no. 6, pp. 2346–2356, 2008. [27] T. Blumensath and M. E. Davies, “Iterative hard thresholding for compressed sensing,” Applied and Computational Harmonic Analysis, vol. 27, no. 3, pp. 265–274, 2009. [28] T. Blumensath and M. E. Davies, “Normalized iterative hard thresholding: Guaranteed stability and performance,” IEEE J. Sel. Topics Siginal Proc., vol. 4, no. 2, pp. 298–309, 2010. [29] T. Blumensath and M. E. Davies, “Gradient pursuits,” IEEE Trans. Signal Proc., vol. 56, no. 6, pp. 2370–2382, 2008. [30] T. Blumensath and M. E. Davies, “Stagewise weak gradient pursuits,” IEEE Trans. Signal Proc., vol. 57, no. 11, pp. 4333–4346, 2009. [31] J. Zhang, Z. L. Yu, Y. Li, G. N. Liu, and Z. Gu, “Weighted regularized sparse recovery method for optical power monitoring,” IEEE Photonics Tech. Let., vol. 24, no. 1, pp. 55–57, 2012. [32] D. Ingrid, Ten lectures on wavelets. Society for Industrial and Applied Mathematics, 1992. [33] J. Bergh, F. Ekstedt, and M. Lindberg, Wavelets: Mit Anwendungen in Signal-und Bildverarbeitung. Springer London, Limited, 2007. [34] R. Nygaard and D. Haugland, “Compressing ECG signals by piecewise polynomial approximation,” in Acoustics, Speech and Signal Processing, 1998. Proceedings of the 1998 IEEE International Conference on, vol. 3. IEEE, 1998, pp. 1809–1812. [35] Y. Zigel, A. Cohen, and A. Katz, “The weighted diagnostic distortion (WDD) measure for ECG signal compression,” IEEE Trans. Biomed. Eng., vol. 47, no. 11, pp. 1422–1430, 2000. [36] SparseLab Toolbox [Online]. Available:http://sparselab.stanford.edu/. [37] Sparsify Toolbox [Online]. Available:http://www.personal.soton.ac.uk/ tb1m08/. [38] BCS Toolbox [Online]. Available:http://people.ee.duke.edu/ ∼lcarin/ BCS.html.

2168-2194 (c) 2013 IEEE. Personal use is permitted, but republication/redistribution requires IEEE permission. See http://www.ieee.org/publications_standards/publications/rights/index.html for more information.