# A New Generalized Orthogonal Matching Pursuit Method

Authors: Liquan Zhao; Yulong Liu
Journal: Journal of Electrical and Computer Engineering
 (2017)

---

## Abstract

To improve the reconstruction performance of the generalized orthogonal matching pursuit, an improved method is proposed. Columns are selected from the sensing matrix by generalized orthogonal matching pursuit, and indices of the columns are added to the estimated support set to reconstruct a sparse signal. Those columns contain error columns that can reduce the reconstruction performance. Therefore, the proposed algorithm adds a backtracking process to remove the low-reliability columns from the selected column set. For anyk-sparse signal, the proposed method firstly computes the correlation between the columns of the sensing matrix and the residual vector and then selects s columns that correspond to the s largest correlation in magnitude and adds their indices to the estimated support set in each iteration. Secondly, the proposed algorithm projects the measurements onto the space that consists of those selected columns and calculates the projection coefficient vector. When the size of the support set is larger than k, the proposed method will select k high-reliability indices using a search strategy from the support set. Finally, the proposed method updates the estimated support set using the selected k high-reliability indices. The simulation results demonstrate that the proposed algorithm has a better recovery performance.


---

## Body

## 1. Introduction

In recent years, a new theory named compressive sensing (CS) [1] has surpassed the limits of the Nyquist sampling rate. Because CS can recover signals at a sampling frequency far lower than the Nyquist sampling rate, CS has aroused tremendous interests over the past few years [2, 3]. CS differs from the traditional Nyquist sampling theory and includes three procedures: sparse representation, nonrelated linear measurement, and signal reconstruction. The reconstruction algorithm aims to recover signals accurately from the measurements, and this step is one of most important parts of CS.Recently, many reconstruction algorithms have been proposed to obtain the original sparse signal from measurements. Two major classes of reconstruction algorithms arel1-minimization and greedy pursuit algorithms. Common l1-minimization approaches include basis pursuit (BP) [4], Gradient projection for sparse reconstruction (GPSR) [5], iterative thresholding (IT) [6], and other algorithms. Those algorithms possess good performance in solving a convex minimization problem, but they have a higher computational complexity.Greedy algorithms have received increasing attention for their excellent performance and small cost in recovering sparse signals from compressed measurements. A greedy algorithm proposed early on was the matching pursuit algorithm (MP) [7]. Building on the MP algorithm, the orthogonal matching pursuit algorithm (OMP) [8] was proposed to optimize the MP via orthogonalization of the estimate support set. The OMP has become a well-known greedy algorithm with wide application. The regularized orthogonal matching pursuit algorithm (ROMP) [9] was developed to refine the selected columns of the measurement matrix with a regularized rule to improve the speed of OMP. The stage wise orthogonal matching pursuit (StOMP) [10] selects multiple columns in each iteration via a presupposed threshold. The subspace pursuit (SP) [11] and compressive sampling matching pursuit (CoSaMP) [12] proposed similar improvement methods. Both of these algorithms were proposed with the idea of backtracking, and the difference is that SP selects k columns from the sensing matrix for each iteration, while CoSaMP selects 2k. The generalized orthogonal matching pursuit (GOMP) was proposed by Wang et al. [13, 14]. The algorithm selects S(S≤K) columns in each iteration. When S=1, GOMP is identical to OMP. Compared to OMP, which selects only one column in each iteration, GOMP changes the number of columns that are selected in each iteration to improve the computational efficiency and recovery performance. The generalized OMP (GOMP) has received increasing attention in recent years, because the method can enhance the recovery performance of OMP. Several papers have been published on the analysis of the theoretical performance of GOMP [13–17].


## 2. Compressive Sensing Model

Compressive sensing requires that the target signal is aK-sparse signal. It means that if we regard the signal as a dimensional vector x, there should be at most K no-zero elements in x. However, in practical applications, sparse signals may not exist in all cases. The target signal has to be transformed into a sparse signal based on a set of sparse basis Ψ=ϕ1,ϕ2,ϕ3,…,ϕN. In this case, x can be defined as (1)x=∑i=1Nαiϕi=Ψα,where α0=K; ·0 denotes the number of nonzero elements in a vector. Thus, the signal x is equivalently represented by K-sparse vector α under some linear transformation Ψ in some domains. The process of compressive sensing can be regarded as a technique that automatically selects relevant information from signals by a measurement. In the theory, x is translated into M-dimensional measurements y via a matrix multiplication with Φ. We express it as (2)y=Φx,where Φ is defined as the measurement matrix with dimensions M×N. Combining (1) with (2), we can obtain (3)y=Φx=ΦΨα=ACSα,where ACS=ΦΨ.In most scenarios,M≪N. Thus, we can interpret that x is compressed into y with a dimension ranging from N to M. Clearly, (3) is an underdetermined equation, and it is difficult to obtain an accurate solution based on the equation. That is to say, it is impossible with traditional methods that obtain the inverse of the matrix Φ to reconstruct the original signal x accurately. In this case, we can obtain x by solving the lp-minimization problem:(4)α^=minαps.t.y=ACSα.Several methods exist for solving this problem. When p=1, the problem is a l1-minimization problem, which can be solved by using a convex optimization algorithm. When p=0, the problem is a l0-minimization problem, which can be solved using a greedy algorithm. An appropriate condition for exact recovery is that the matrix ACS satisfies the condition of restricted isometry property (RIP) condition [1].Definition 1.
A sensing MatrixACS is said to satisfy the RIP condition with the smallest number of the K-restricted isometry constant δK(δK∈(0,1)), if(5)1-δKα2≤ACSα≤1+δKα2holds for any K-sparse vector α∈RN×1 with α0≤K.


## 3. GOMP Algorithm

Greedy algorithms are used widely to recover signals in CS due to their simple algorithms and low computational complexity. Two methods exist for improving the OMP algorithm. The first method is based on the idea of backtracking, as in the subspace pursuit algorithm (SP) [11]. The second method selects more than one atom in each iteration, as in GOMP. The computational complexity of the backtracking method is higher, but it yields higher accuracy in most cases. In this section, we will introduce the GOMP algorithm.GOMP firstly computes the correlation between the columns of the sensing matrixΦ and the residual vector rk-1 by Φ′rk-1 and then selects S columns that correspond to the S largest correlation in magnitude adding their indices to the estimated support set Λk in each iteration. Next the projection coefficient vector x^Λk of measurements y onto space of span(ΦΛk) is obtained using the least square method (LS). The residual rk is revised by subtracting ΦΛkx^Λk from y. These operations are repeated until either the iteration number reaches the maximum k=min{K,M/S}, or the l2-norm of the residual falls below a threshold ε. For ease of understanding, we describe the GOMP algorithm in Algorithm 1 according to [13, 14].Algorithm 1:GOMP algorithm.
Input:  measurementsy∈RM×N,
sensing matrixΦ∈RM×N,
SparsityK,
Initialize: number of indices of columns for each
selectionS(S≤K).
iteration countk=0,
residual vectorr0=y,
estimated support setΛ0=⌀.
Whilerk>ε and k<min⁡{K,M/S} do
k=k+1;
(Identification)
SelectS largest entries (in magnitude) from Φ′rk-1.
Then record the indices{φ(i)}i=1,2,3,…,S corresponding
to the entries.
(Augmentation)
Λk=Λk-1∪φ1,φ2,…,φS.
(Estimation ofxΛk)
x^Λk=argminsup⁡p(u)=Λky-Φu2.
(Residual Update)
rk=y-ΦΛkx^Λk.
End
Output The estimated supportΛ^=argminT:‖T‖=K⁡x^Λk-x^T2
and signalx^Λ^=Φ†Λ^y.In Algorithm1, we observe that the only difference between GOMP and OMP is that the GOMP selects more than one atom in each iteration. The GOMP algorithm selects S(S≥1) atoms in each iteration. When S=1, the GOMP is identical to OMP.


## 4. Improved GOMP Algorithm

For the usage of greedy algorithms, it is important to generate an estimate of the correct support set. We assume the correct support isT, and the estimate of T is Λk. The goal is to determine the indices for Λk that are most similar to T. GOMP selects S columns from the sensing matrix Φ and added the indices of these columns to the estimated support set Λk to reconstruct a sparse signal. Those columns contain columns that be selected by error. Although GOMP allows for a number of error indices in Λk, these error indices will lead to the size of Λk larger than K. This increases the algorithm complexity greatly in the process of estimation of xΛk and residual update in the GOMP.In the process of estimation ofxΛk, the projection coefficient vector of the measurements onto the space of span(ΦΛk) can be expressed as(6)x^Λk=ΦΛk†y,where ΦΛk†=(ΦΛkTΦΛk)-1ΦΛkT; y are the measurements of the signal vector x. The cost of (6) is(7)4S2km+-2S2+5Nm+2S3k2+-4S3+5S2k+3S3-S2-S,where S represents the number of indices of columns for each selection, k is the number of iterations, and m is the dimensionality of measurements y [13].In the process of estimation of residual update, the residuals in GOMP can be expressed as(8)rk=y-ΦΛkx^Λk,where y represents the measurements, ΦΛk is the estimated support set, and x^Λk is the projection coefficient vector of the measurements onto the space of the estimated support set. The cost of ΦΛkx^Λk is (9)2kS-1m,where k is the number of iterations and S is the number of indices of columns for each selection.GOMP selectsS indices in each iteration and adds them to Λk. On the k-th iteration, kS indices will be selected, and the dimensionality of ΦΛk is m×kS. If all the indices in the estimated support set Λk are correct, Λk is identical to T and kS≈K. Although a selection rule exists to ensure that newly added indices belong to the correct support T in GOMP, it is unavoidable that the error indices are selected and once selected, the error indices of the selected columns will remain in the support set throughout the remainder of the reconstruction process. It means kS>K. Generally speaking, when a large S is selected, the GOMP should exert higher efficiency than a small S. However, as S increases, the probability of the selected error indices increases as well. When a large number of error indices are selected, it leads to a great increase in kS, further increasing the cost of the algorithm.To overcome this problem, we propose a method to improve the performance of GOMP. The proposed method will retain the size of the estimated support set asK and update the indices using the backtracking method to reduce the number of error indices, when the number of indices is greater than K. Even if an index is deemed reliable in some iteration, when it is considered unreliable in subsequent iterations, the index will be removed from the estimated support set. This can reduce the number of error indices in the estimated support set.The main difference between the proposed and the conventional GOMP algorithm is that the new algorithm can add or remove the index from the estimated support set at any stage of the recovery process. The proposed algorithm increases the reliability of the indices of the estimated support set and furthermore reduces the cost of the reconstruction process. The process of proposed algorithm is expressed in Algorithm2.Algorithm 2:Improved method of GOMP.
Input:  measurementsy∈RM×N,
sensing matrixΦ∈RM×N,
SparsityK,
Initialize: number of indices of columns for each selection
S(S≤K).
iteration countk=0,
residual vectorr0=y,
estimated support setΛ0=⌀.
Whilerk>ε and k<min⁡K,M/S do
k=k+1;
(Identification)
SelectS largest entries (in magnitude) from Φ′rk-1. Then
record the indices{φ(i)}i=1,2,3,…,S corresponding to the
entries.
(Augmentation)
Λk=Λk-1∪{φ1,φ2,…,φ(S)}.
(Estimation ofxΛk)
xΛk=arg⁡minsuppu=Λk⁡y-Φu2.
(Backtracking) While xΛk0≥K, select K largest
elements of|xΛk|. Then recording the indices corresponding
to the elements, and renew theΛk with those indices.
(Residual Update)  rk=y-ΦΛkx^Λk.
End
OutputThe estimated support Λ^=argminT:T=K⁡x^Λk-x^T2
and signalx^Λ^=Φ†Λ^y.In Algorithm2, we observe that the proposed method adds a backtracking step to the GOMP algorithm. When S is small, it has a high accuracy for the identification process. When S=1, GOMP is identical to OMP. At this step, GOMP has the highest accuracy in the identification step. As S increases, the probability of the selected error indices increases as well. The proposed method adds a backtracking process to remove those error indices in subsequent iterations. When S≤K, the proposed method requires more than one iteration to create a K-sized support set. When kS<K, the support set cannot contain K correct indices, and executing the backtracking step is unnecessary and will cause wasteful computation. Therefore, we designed the algorithm to execute the backtracking step when kS≥K. The dimensionality of xΛk is kS×1 before executing the backtracking step. So the improved algorithm determines whether to execute the backtracking step by the size of xΛk. Experimental evidence demonstrated that our changes improved the performance. We will describe our simulation in the next section.


## 5. Simulation and Discussion

In this part, we will demonstrate the performance of the proposed algorithm with sparse signals. We used the same sparse signal sources as in GOMP withK-sparse to compare the performances of the different algorithms. The components of the sensing matrix were generated randomly with Gaussian distribution, and the size of the sensing matrix was 128×256. We used MATLAB 7.0 with a quad-core 64-bit processor in a Windows 10 environment. We executed each algorithm 1000 times and recorded the probability of the exact reconstructions. We set the threshold ε=10-6 used in [13] in both GOMP and the proposed method.In Figure1, we compare the probability of successful reconstruction of the proposed method with the ROMP [9], OMP [8], SP [11], CoSaMP [12], and StOMP [10] algorithms. It is evident from Figure 1 that the reconstruction performance of ROMP decreased rapidly with a rise in K. When K=20, the recovery probability was near zero. The probability of successful recovery of OMP decreased at a smaller rate than that of ROMP, and its critical sparsity was reached at K=55. CoSaMP and StOMP had similar performances. Both retained a high recovery performance when K≤40. The probability of successful recovery of SP was higher than for ROMP, OMP, CoSaMP, and StOMP with the same K. For our proposed method, we tested the performance with different S values, representing the number of selected columns in each iteration. We set S=10 and S=20. Figure 1 indicates that the smaller the value of S, the higher the probability of successful recovery of the proposed algorithm with the same K. Compared to the other algorithms, under the same conditions, the proposed method with S=10 had the highest probability of successful recovery, followed by the proposed method with S=20, and followed by the SP algorithm. The proposed algorithm with different S values had a higher probability than the other algorithms.Figure 1
Reconstruction performance forK-sparse Gaussian signal vector as a function of sparsity K.In Figures2 and 3, we compare the performances of the proposed method and GOMP with different values for the parameter S. In Figure 2, we selected smaller S values to compare the algorithms. All algorithms resulted in high levels of probability of recovery when K≤40. The probability of exact recovery decreased with an increase in K and reached zero when K=70. Both the proposed algorithm and GOMP had a higher probability of exact recovery with a smaller S value. When S=6, the two algorithms had nearly the same probability values. When S=9, the proposed method had a higher probability of recovery. The probability of exact reconstruction of both algorithms was 100% when K<40 and zero when K>70 for S=6 and S=9, and the difference between two algorithms mainly concentrated on the sparsity level from 40 to 70. Therefore, we compared two algorithms with K=40 to K=70 in Figure 2.Figure 2
Reconstruction performance of GOMP and the proposed method with a smallS value.Figure 3
Reconstruction performance of GOMP and the proposed method with a largeS value.In Figure3, we selected a larger value of S to compare the proposed algorithm and GOMP. The probability of exact reconstruction of both algorithms was 100% for K<10 and zero for K>70, so we compared the two algorithms with K=10 to K=70 for S=20 and S=30 in Figure 3.WhenK=40 and S=20, the probability of exact recovery of GOMP was 48.7%, while the probability value for the proposed method was 99.4%. When K=40 and S=30, the resulting probability values were 4.5% for GOMP and 99.4% for the proposed method. These results indicated that, for different S values, the proposed method had a higher probability of exact recovery than the GOMP method.Figures2 and 3 demonstrate that, for different S values, the curves of probability for exact recovery were closer to each other for the proposed method than the GOMP method. This indicated that the proposed method had a more stable performance than the GOMP method for different S values.Based on the analysis and comparison, we determined that even though the proposed algorithm had an additional backtracking process compared to GOMP, the proposed method demonstrated an excellent performance with regard to running time. In order to compare running times for all the algorithms, we ran each algorithm one thousand times to calculate the average running time. The computing environment was the same as for the determination of the probability of exact reconstruction.In Figure2, we can see that the probability of extract reconstruction of both algorithms was below 100% for K>40 with S=6 and S=9, and the effective running time should be computed for successful reconstructing. Therefore, we just compared the running time with K<40 for two methods. The difference was not obvious with K<5 in running time for both methods with S=6 and S=9, so we compared the running time with K>5 and K<40 in Figure 4.Figure 4
The average running time of GOMP and the proposed method with a smallS value.In Figure3, we can see that probability of extract reconstruction of GOMP algorithm with S=20 was below 100% for K<20, and with S=30 for K<10, and the probability of extract reconstruction of proposed algorithm with S=20 and S=30 all were below 100% for K<40. The effective running time should be computed for successful reconstructing, so we just compared the running time of both algorithms for S=20 with K<20, and S=30 with K<10. The difference was not obvious in running time for both methods with smaller K. Based on the above analysis, we just compared their running time with K>5 and K<20 for S=20, and K>5 and K<10 for S=30. In order to unify the sparsity level scope, we set K from 5 to 20 in Figure 5.Figure 5
The average running time of GOMP and the proposed method with a largeS value.In Figure4, we set S=6 and S=9 to compare the running times for GOMP and the proposed method. The running times for GOMP with S=6 and the proposed method with S=6 and S=9 were nearly identical. The running times under those conditions were less than those for GOMP with S=9. When K=40, the running time of the proposed method with S=6 and S=9 was less than that of GOMP with S=6.In Figure5, we set S=20 and S=30. From Figure 5, we can observe that the running time of the proposed algorithms was less than GOMP for S=20 with K=5 to 20, and for S=30 with K=5 to 10. This shows that the proposed algorithm has faster speed for reconstructing.These results show that the proposed method resulted not only in better performance with regard to running time, but also in a higher probability of exact reconstruction. The advantage is more obvious when largerS is taken in both algorithms.


## 6. Conclusion

In this paper, a novel method for sparse signal reconstruction has been proposed. The proposed method adds a backtracking step to maintain aK-sized estimated support set, avoiding the extra computation cost for an oversized estimated set. At the same time, the proposed method can increase the reliability of the estimated support set by removing the low-reliability columns from the estimated support set. The proposed reconstruction algorithm performs well at not only reconstructing the sparse signal (i.e., when K is small), but also the less sparse signal (i.e., when K is large). The simulation results showed that the proposed method had a better performance than GOMP with regard to both the probability of exact reconstruction and running time, especially with larger S values. In future research, we want to improve on the proposed algorithm to optimize the greedy algorithm further.



---

*Source: 3458054-2017-08-13.xml*