Abstract

Compressed sensing (CS) is an important research area of signal sampling and compression, and the essence of signal recovery in CS is an optimization problem of solving the underdetermined system of equations. Greedy pursuit algorithms are widely used to solve this problem. They have low computational complexity; however, their recovery performance is limited. In this paper, an intelligence recovery algorithm is proposed by combining the Bat Algorithm (BA) and the pruning technique in subspace pursuit. Experimental results illustrate that the proposed algorithm has better recovery performance than greedy pursuit algorithms. Moreover, applied to the microseismic monitoring system, the BA can recover the signal well.

1. Introduction

In recent years, the collected data quantity by the Internet of things (IoT) has increased dramatically. To facilitate storage and transmission, signal sampling based on the traditional Shannon–Nyquist sampling theory works as follows: massive data are collected at the sampling stage and most are discarded at the compression stage. This process is extremely wasteful. CS, which was proposed by D. Donoho, E. Candès, and T. Tao, brings a revolutionary breakthrough in signal sampling [13]. If only nonzero elements exist in a signal, the signal is called -sparse signal with sparsity . CS theory indicates that sparse signals can be recovered more accurately using less measurements than the Nyquist sampling principle.

The signal recovery problem in CS is a nonconvex combinatorial optimization problem. It can be solved by three types of algorithms. The first type is greedy pursuit algorithms, including Matching Pursuit (MP), Orthogonal Matching Pursuit (OMP) [4], Generalized Orthogonal Matching Pursuit (GOMP) [5], and Subspace Pursuit (SP) [6]. They are two-stage algorithms. In the first stage, they seek support iteratively. In the second stage, the signal can be recovered by using the least square method. This kind of algorithm is effective; however, they have weak recovery performance. The second kind of algorithm is a convex optimization algorithm, which recovers signals by converting nonconvex problems to convex problems. The most common convex optimization algorithms are Basic Pursuit (BP) and Gradient Projection for Sparse Reconstruction (GPSR) [7]. Although they have high accuracy, their computational complexity is high. The third kind of recovery algorithm is Bayesian algorithms [8], which can be divided into the following two types: Maximum A Posteriori (MAP) Estimation and Hierarchical Bayesian. The former underlines a signal distribution, and the latter introduces one or two variables, which control the sparse signal. Bayesian algorithms achieve a balance of high accuracy and short recovery time [9, 10]. In recent years, some scholars have applied Swarm Intelligence Algorithms to signal recovery in CS, such as Particle Swarm Optimization (PSO) [11] and the Grey Wolf Optimizer Algorithm [12]. These methods have good global search ability. However, they have a slow convergence velocity and are easy to be trapped into local optimum.

Recently, Yang et al. simulated bat echolocation and proposed a Bat Algorithm (BA) [13] based on stochastic optimization, which has simple structure, fast convergence, strong robustness, and parallel computing. The BA combines the advantages of PSO, Harmony Search (HS), Simulated Annealing (SA), and other algorithms. Many scholars applied it to optimization problems, such as numerical optimization, economic dispatch with random wind power [14], multirobot formation reconfiguration [15], microgrid operation management [16], and image processing [17]. In these areas, the BA shows better performance than do other intelligence optimizer algorithms.

In this paper, inspired by the BA and the greedy pursuit algorithm, an intelligence recovery algorithm for CS is proposed. OMP is used to initialize the swarm. The BA and the pruning technique are combined to update populations. The proposed algorithm outperforms the traditional greedy pursuit algorithms. Moreover, it runs faster than other swarm intelligence algorithms such as PSO.

The remainder of this paper is organized as follows: Section 2 briefly introduces the basic theory of CS and the BA. In Section 3, the proposed algorithm is introduced in detail. In Section 4, experiments are done to evaluate the algorithm performance. Finally, conclusions are drawn in Section 5.

The notations used in the rest of this paper include the following: denotes the absolute value of a real number or the cardinality of a set. represents . For the vector , is the estimate of . In a matrix , is the -th column of . , which denotes the positions of nonzero elements in x, and it is called support set. is the universal set of support. The matrix is submatrix of . The pseudoinverse of a matrix is defined as , where denotes the transpose operator.

2. Basic Theory

2.1. Compressed Sensing

Supposing an n-dimensional signal can be represented on a basis as , where , if there are only nonzero elements in , is called signal with sparsity , and the positions of nonzero elements in are called the support set, which can be denoted as . The signal can be sampled by projection onto the measurement matrix aswhere is the measurement signal, is the original signal, is the sparse signal, is the measurement matrix, and is the sparse basis. is the sensing matrix.

It can be seen that the data quantity is greatly reduced because is an underdetermined system of equations. There are infinitely many solutions to recover from . However, if is sparse, it can be recovered from by solving the following minimization problem:where represents .

2.2. Orthogonal Matching Pursuit (OMP)

In the proposed algorithm, we use OMP to determine the initial solution.

The OMP algorithm can be stated in Algorithm 1 [4].

Input: Sensing matrix , measurement signal , sparsity .
Output: Estimated signal , estimated support set
Initialization: Residual signal , estimated index set , estimated atomic set , iteration number ;
 The -thiteration:
  Step 1: Find
  Step 2: Update the estimated support set and the estimated atomic set ;
  Step 3: Estimate the signal ;
  Step 4: Update the residual , update iteration counter ;
 Judgement: If , continue iterating; if , stop the iteration, and output the estimated support set and the estimated signal .
2.3. Bat Algorithm

In 2010, Yang et al. simulated the characteristics of bats’ predatory behaviors and proposed BA [13, 18], which is a new Swarm Intelligence Optimization Algorithm.

The first step is randomly initializing the bat positions and the velocities . Then, the BA updates the velocities and positions according to certain strategies. The velocities and positions of the -th bat at iteration are updated as follows:where is the pulse frequency of the -th bat, , is a random number. and represent the positions of the -th bat at the and the iteration, respectively; and are the velocities of the -th bat at the and the iteration, respectively; and denotes the global best position of the swarm.

For each bat, the algorithm generates a random number . If , the position of this bat is updated according to formula (6), where is the rate of pulse emission. The random walk is regarded as a procedure of the local search [13].where is a random solution in the best solution set composed of the best solution of each bat, is a random number, and is the mean of the pulse loudness in the -th iteration.

In addition, with the iterative number increasing, the pulse loudness and the rate of pulse emission need to be updated. During prey searching, the pulse loudness will gradually decline, and the rate of pulse emission will also gradually increase. The updating formulas of and are found in formulas (7) and (8), respectively:where and are constants, , and .

The BA algorithm can be described in Algorithm 2 [13].

Fitness function:
Initialization:
Randomly initialize the positions and velocities of the bat population.
Initialize the pulse frequency of each bat according to formula (3).
Initialize the rate of pulse emission and the pulse loudness .
Calculate the fitness value of each bat and find the best solution in the population.
Set the stopping criterion .
Termination condition:
Iteration counter reaches or the fitness value is smaller than .
Iteration:
Step 1: Update the bat velocity and position according to formulas (3)∼(5);
Step 2: If (), select a random solution in the best solutions set and update position of this bat by the random walk around the selected best solution as shown in formula (5);
Step 3: Calculate the fitness value of each bat’s new position;
Step 4: If (), accept the new solution and adjust , according to formulas (7) and (8);
Step 5: Find the new best solution and update the best solutions set;
Step 6: If termination condition is satisfied, stop the iteration; otherwise, continue.

3. A Bat-Inspired Sparse Recovery Algorithm

3.1. Fitness Function

If the sparsity level satisfying is known as a priori, problem (2) can be approximated as follows [11, 12]:where is the measurement signal, and is the sensing matrix.

From greedy pursuit algorithms, we can see that the signal can be recovered by using a two-step strategy. First, estimate the support set and then estimate the signal using the least square method.where denotes the estimated support set and . The vector consists of the entries of indexed by . The matrix is the submatrix of , and is the -th column of . is the pseudoinverse of a matrix , where denotes the transpose operator.

Once the support set is estimated accurately, the following equation must be satisfied.where represents .

Because , we define the fitness function as

3.2. Initialization

Suppose there are bats, the position of the first bat is initialized as the estimated support set of the OMP. For the -th () bat, its position is initialized as follows.

Randomly choose elements from the set and form the set ; here, is regarded as the velocity of the -th bat at the iteration. Then, the position of the -th () bat is a set formed by indices corresponding to the maximum absolute values in . Other related parameters are initialized as follows: the pulsing frequency of each bat is initialized according to formula (3), where . The initial emission rate . The pulse loudness , is a random number between 0 and 1, and is a random number between −1 and 1, and .

3.3. Update Strategy

The main innovation of our algorithm is the update strategy. The iterative process combines the thought of BA and greedy pursuit algorithm. The update strategy can be divided into three parts.

The first part is the update of velocity and position: a set is formed with elements selected randomly from the best solution ; we define . If , a set consists of elements selected randomly from the set , update . If , the set consists of elements selected randomly from the set , update . The position of the -th bat is a set composed of the indices corresponding to the maximum absolute values in .

The second part is the local search: generate a random number ; if , select a random solution from the best solutions set and update the position of this bat by randomly replacing elements of this solution with other elements in . After updating all positions, calculate the new fitness of each bat. Specifically, in order to accelerate the convergence speed of this algorithm, all the new solutions will be accepted in the algorithm.

The third part is the adjustment of the emission rate and the pulse loudness as formulas (7) and (8). If , update the best solutions set . If , update the global best solution .

3.4. Stopping Criterion

If , the iteration is terminated; here, is the termination threshold. Moreover, in order to avoid too many iterations, set a maximum allowed iteration number . If the iteration number is larger than , the iteration is also terminated.

The algorithmic process in Algorithm 3.

Input: Sensing matrix , measurement vector , sparsity level , termination threshold , maximum allowed iterative number .
Output: The estimated signal
Fitness function:
Initialization: Initialize the positions and velocities of bats (Section 3.2).
Termination condition: If or , terminated the iteration.
Iteration:
Step 1: A set is formed with elements selected randomly from the best solution , define .
Step 2: If , a set consists of elements selected randomly from the set , update . If , the set consists of elements selected randomly from the set , update .
Step 3: The position of the -th bat is a set composed of the indices corresponding to the maximum absolute values in .
Step 4: If , select an optimal solution from the best solutions set and then replace elements of this solution with other elements in .
Step 5: Calculate the new fitness .
Step 6: If , adjust and .
Step 7: If , update the best solutions set . If , update the best solution .
Step 8: If termination condition is satisfied, stop iteration; otherwise, continue.
3.5. Efficiency Analysis

The computational complexity of the algorithm is mainly dependent on the initialization and iteration.

The algorithm initializes solutions using OMP, which has the complexity [19].

In each iteration, the complexity is mainly dependent on three operations. The first is the multiplication of matrix, which has the computational complexity . The second is sorting of an n-dimensional vector, which has the computational complexity . The third is the least square, which has the computational complexity . In summary, the computational complexity in the iteration phase is , where is the population size.

In general, the computational complexity of the proposed algorithm is similar to OMP when the signal sparsity is small. The complexity depends on the number of iterations when the sparsity is larger or when the measurement number is small. Moreover, this algorithm uses the idea of the BA, which combines the advantages of PSO, HS, SA, and other intelligence algorithms to optimize the convergence speed and search accuracy. It usually runs faster than other intelligence algorithms. OMP can generate an initial solution closed to the target, which accelerates the convergence of the BA. In addition, as swarm intelligence algorithm, it can run in parallel to reduce the running time [12, 20].

4. Experiments

In this section, we compare this algorithm with OMP, GOMP, SP, Fast Laplace, and PSO, where Fast Laplace is a Bayesian-based algorithm, and PSO is a swarm intelligence algorithm for CS. Then, we apply this algorithm to recover the real microseismic signals.

Candes and Tao proved that the Gaussian random matrix [2, 21, 22] with independent identically distribution can be a universal measurement matrix. Therefore, the measurement matrices in these experiments are the Gaussian random matrix with the size of . The number of iterations of the OMP, GOMP, and SP algorithms is . The maximum number of iterations of PSO and the proposed algorithm are . The algorithm will stop when or the number of iterations reaches .

The software and hardware environment of the experiment is as follows:Processor: Intel (R) Pentium (R) CPU G3220 @ 3.00 GHz 3.0Memory: 4.00 GOperating system: 64-bit Windows7Simulation software: MATLAB R2014a

4.1. Comparison with Other Algorithms

In the first experiment, we compare the recovery performance of different algorithms against the change of sparsity. We use a 128 × 256 random Gaussian matrix to sample a 256-length sparse signal. The sparsity level is set as . The parameters are set as , , , ; the rate of pulse emission and the pulse loudness are set between 0 and 1; , , and . Here, we compare the algorithm performance with the following two metrics: recovery error and exact recovery rate.

Recovery error is a metric to evaluate the error between the original sparse signal and the recovered sparse signal. The recovery error is defined as follows (13):where is the original sparse signal, is the recovered sparse signal, and is the cycle number to obtain the average performance. The simulation result of recovery error is shown in Figure 1(a). The exact recovery rate is obtained by calculating the number of times that the estimated values are exact. If , this experiment is considered successful. In 100 experiments, if there are successful experiments, the exact recovery rate is . The exact recovery rates of different algorithms are shown in Figure 1(b).

It can be seen that when the sparsity is large, the recovery errors of three commonly used greedy pursuit algorithms such as OMP, GOMP, and SP increase rapidly, and the exact recovery rate decreases quickly. The recovery error of the Fast Laplace algorithm is lower than that of traditional greedy pursuit algorithms. However, experiments demonstrate that the exact recovery rates of the Bayesian algorithm are 0 because there exists no experiment satisfying . The swarm intelligence algorithms also perform well. The recovery error of the proposed algorithm is lower than 0.4, and the exact recovery rate is higher than other five algorithms when the sparsity is very large.

In the second experiment, we compare the recovery performances against the change of measurement number. We use a Gaussian matrix to sample a 256-length signal with sparsity 30. The measurement number is set as , and other parameters are the same to the first experiment.

Figures 2(a) and 2(b) illustrate the recovery error and the exact recovery rate, respectively, of the six algorithms under different measurement numbers. As shown in the figure, the recovery of the proposed algorithm performs better than those commonly used greedy pursuit algorithms such as OMP, GOMP, and SP, when the number of measurement vectors is small. It also has more advantages compared with the PSO algorithm and the Fast Laplace algorithm.

In the third experiment, we compare the recovery efficiency of different algorithms. Figure 3(a) is the running time of the six algorithms with different sparsity. The sparsity level is set as . Figure 3(b) is the running time of the six algorithms with different measurement number. The measurement number is taken as . Other parameters are the same as those of the first experiment, and the average time of each algorithm in recovering the signal is counted. As shown in Figure 3, when the sparsity is small or the measurement number is large, the recovery efficiency of this algorithm is similar to the commonly used greedy pursuit algorithms such as OMP. Although the running time of the algorithm increases when the sparsity is large or the measurement number is small, the algorithm is more efficient than other swarm intelligence algorithms such as PSO. However, as a swarm intelligence algorithm, it can run in parallel to reduce the running time [12, 20].

4.2. Application in Microseismic Signals

Microseismic monitoring is an effective technology for accident prediction in coal mines. Under the background of coal mine IoT, the collected data quantity by numerous microseismic nodes in the monitoring area are massive [23]. The effective compression of microseismic data can reduce the storage space and energy consumption. It plays an important role in improving the data transmission rate. Participants of the project collected three microseismic signals from a coal mine in Anhui Province, China. In this section, we use these signals as the original signal and apply the proposed algorithm to recover signal.

We perform experiments on three different microseismic signals. The length of the signal is , as shown in Figure 4. It can be seen that the microseismic signals are one-dimensional time-frequency signals. Microseismic signals are not sparse in the time domain; however, their representations under a certain basis are sparse or compressible.

In this experiment, the 2-layer db2 wavelet transform [24] is used as sparse basis. Wavelet transform uses a cluster of functions to approximate a signal which can be seen as a decomposition of the signal. It divides the signal into a low-frequency signal and multiple high-frequency signals. The low-frequency signal is the stationary part of the signal. The high-frequency signal indicates the details of the signal. Figure 5 shows the wavelet decomposition coefficients of three groups of microseismic signals.

The measurement matrix used in this experiment is a 256 × 512 random Gaussian matrix. The sparsity level is and . Figure 6 shows the three groups of recovered signals. The recovery errors in the three groups of experiments are 0.1862, 0.1250, and 0.1512.

5. Conclusions

In this paper, an Intelligence Sparse Recovery Algorithm is proposed for CS inspired by the BA and the pruning technique of SP. The algorithm inherits the global search ability and the local search ability of the BA, and we use OMP to obtain a better initial solution to accelerate the convergence. Compared with the traditional pursuit algorithm, the algorithm in this paper has better recovery performance under large sparsity and small measurement. At the same time, it runs faster than the intelligent optimization algorithm such as PSO. The proposed algorithm can be applied to the microseismic monitoring system, which can recover the signal well. In future research, performance and efficiency of the algorithm will be more rigorously analyzed.

Data Availability

The data of microseismic signals used to support the finding of this study have not been made available for the time being because the data are obtained through cooperation between us and the Coal Mine Group, and it must be approved by both parties before it can be disclosed.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

This work was financially supported by the National Science Foundation of China (grant no. 51574232).