Sparse Recovery with Orthogonal Matching Pursuit under RIP

Tong Zhang

I Introduction

Here, AA is an n×dn\times d matrix. If we define an objective function

then we may estimate the parameter xˉ\bar{{\mathbf{x}}} by minimizing Q(x)Q({\mathbf{x}}), subject to appropriate constraints.

If d>nd>n, then the solution of the unconstrained optimization problem

is not unique. In order to estimate xˉ\bar{{\mathbf{x}}}, additional assumptions on xˉ\bar{{\mathbf{x}}} is necessary. We are specifically interested in the case where xˉ\bar{{\mathbf{x}}} is sparse. That is ∥xˉ∥0≪n\|\bar{{\mathbf{x}}}\|_{0}\ll n, where

It is known that under appropriate conditions, it is possible to recover xˉ\bar{{\mathbf{x}}} by solving (2) with a sparsity constraint as follows:

However, this optimization problem is generally NP-hard. Therefore one seeks computationally efficient algorithms that can approximately solve (3), with the goal of recovering sparse signal xˉ\bar{{\mathbf{x}}}. This paper considers the popular orthogonal matching pursuit algorithm (OMP), which has been widely used for this purpose (for example, see ). We are specifically interested in two issues: the performance of OMP in terms of optimizing Q(x)Q({\mathbf{x}}) and the performance of OMP in terms of recovering the sparse signal xˉ\bar{{\mathbf{x}}}.

II Main Result

Our analysis considers a more general objective function Q(x)Q({\mathbf{x}}) that does not necessarily take the quadratic form in (1). However, we assume that Q(x)Q({\mathbf{x}}) is convex. For such a general convex objective function, we consider the fully (or totally) corrective greedy algorithm in Figure 1, which was analyzed in . This paper refines the analysis to show that the algorithm works under the restricted isometry property (RIP) of (the required condition will be described later in this section). This algorithm is a direct generalization of OMP which has been traditionally considered only for the quadratic objective function in (1) with F(0)=∅F^{(0)}=\emptyset. For simplicity, we assume that the number of iterations k0k_{0} is chosen a priori. The algorithm has been known in the machine learning community as a version of boosting , and has also been proposed recently in the signal processing community .

For quadratic loss, the objective function Q(x)Q({\mathbf{x}}) is given by (1) and its derivative is ∇Q(x)=2A⊤(Ax−y)\nabla Q({\mathbf{x}})=2A^{\top}(A{\mathbf{x}}-{\mathbf{y}}). Therefore j=arg⁡max⁡i∣∇Q(x(k−1))i∣j=\arg\max_{i}|\nabla Q({\mathbf{x}}^{(k-1)})_{i}| becomes j=arg⁡max⁡i∣ai⊤(Ax−y)∣j=\arg\max_{i}|{\bf a}_{i}^{\top}(A{\mathbf{x}}-{\mathbf{y}})|, where ai{\bf a}_{i} is the ii-th column of matrix AA. This, together with F(0)=∅F^{(0)}=\emptyset, leads to the standard OMP algorithm. In order to use notation consistent with the sparse recovery literature, in the current paper, we still refer to the more general algorithm in Figure 1 as OMP even though it applies to objective functions other than (1).

The general problem of optimization under sparsity constraint is NP hard. In order to alleviate the difficulty, we consider approximate optimization under the restricted strong convexity assumption introduced below.

Given any s≥0s\geq 0, define restricted strong convexity constants ρ−(s)\rho_{-}(s) and ρ+(s)\rho_{+}(s) as follows: for all ∥x−x′∥0≤s\|{\mathbf{x}}-{\mathbf{x}}^{\prime}\|_{0}\leq s, we require

The restricted isometry constant was used to define the restricted isometry property (RIP) in the analysis of L1L_{1} regularization method . We employ the slightly more general restricted strong convexity constants in (4) because our analysis only requires the ratio ρ+(s)/ρ−(s)\rho_{+}(s)/\rho_{-}(s) to be bounded, and this is useful for general machine learning problems where ρ+(s)\rho_{+}(s) can be larger than 22.

In order to recover the target xˉ\bar{{\mathbf{x}}}, we have to assume that xˉ\bar{{\mathbf{x}}} is sparse and approximately optimizes Q(x)Q({\mathbf{x}}). If a target xˉ\bar{{\mathbf{x}}} is an exact global optimal solution, then ∇Q(xˉ)=0\nabla Q(\bar{{\mathbf{x}}})=0. However, this paper deals with approximate optimal solutions, where ∇Q(xˉ)≈0\nabla Q(\bar{{\mathbf{x}}})\approx 0. In particular, we introduce the following definition, which is convenient to apply.

The constant ϵs(xˉ)\epsilon_{s}(\bar{{\mathbf{x}}}) measures how close is ∇Q(xˉ)\nabla Q(\bar{{\mathbf{x}}}) to zero. If ∇Q(xˉ)=0\nabla Q(\bar{{\mathbf{x}}})=0, then ϵs(xˉ)=0\epsilon_{s}(\bar{{\mathbf{x}}})=0. If ∇Q(xˉ)≈0\nabla Q(\bar{{\mathbf{x}}})\approx 0, then ϵs(xˉ)\epsilon_{s}(\bar{{\mathbf{x}}}) is small. Moreover, similar to the definition of restricted strong convex constants, we are only interested in the value of ∇Q(xˉ)\nabla Q(\bar{{\mathbf{x}}}) in any subset of {1,…,d}\{1,\ldots,d\} with ss elements. The following proposition provides some estimates of ϵs(xˉ)\epsilon_{s}(\bar{{\mathbf{x}}}) using quantities that are easier to understand.

We have ϵs(xˉ)≤s∥∇Q(xˉ)∥∞\epsilon_{s}(\bar{{\mathbf{x}}})\leq\sqrt{s}\|\nabla Q(\bar{{\mathbf{x}}})\|_{\infty} and ϵs(xˉ)≤∥∇Q(xˉ)∥2\epsilon_{s}(\bar{{\mathbf{x}}})\leq\|\nabla Q(\bar{{\mathbf{x}}})\|_{2}. Moreover, if

The first two inequalities are straight-forward. For the third inequality, we note that for ∥u∥0≤s\|{\mathbf{u}}\|_{0}\leq s:

The result follows by rearranging the above inequality. ∎

The following theorem is the main result of this paper, which shows that OMP can approximately recover a sparse signal xˉ\bar{{\mathbf{x}}} in 2-norm if the condition (5) in the theorem involving strong convexity constants can be satisfied. As we shall discuss later, this condition is closely related to the RIP condition for the quadratic objective (1).

then when k=k0=s−∣Fˉ∪F(0)∣k=k_{0}=s-|\bar{F}\cup F^{(0)}|, we have

The detailed proof relies on a number of technical lemmas that are left to the appendix.

The first inequality of the theorem is a direct consequence of Lemma 9. The second inequality is a consequence of the first inequality and Lemma A.2:

Note that (5) can be satisfied as long as (ρ+(1)/ρ−(s))ln⁡(ρ+(kˉ)/ρ−(s))(\rho_{+}(1)/\rho_{-}(s))\ln(\rho_{+}(\bar{k})/\rho_{-}(s)) grows sub-linearly as a function of ss. With appropriate assumptions, this allows the ratio ρ+(s)/ρ−(s)\rho_{+}(s)/\rho_{-}(s) to be significantly larger than 11 but bounded from above (such a condition is sometimes referred to as sparse eigenvalue condition in the statistics literature). In this context, Theorem II.1 is useful for estimation problems encountered in machine learning, where ρ+(s)/ρ−(s)\rho_{+}(s)/\rho_{-}(s) may be large.

In compressed sensing, one can often control the ratio of ρ+(s)/ρ−(s)\rho_{+}(s)/\rho_{-}(s) to be not much larger than 11 using random projection. In this context, the following result gives a simpler interpretation of the above theorem, where the condition (5) of the theorem is replaced by ρ+(kˉ)≤2ρ−(31kˉ)\rho_{+}(\bar{k})\leq 2\rho_{-}(31\bar{k}).

If ρ+(kˉ)≤2ρ−(31kˉ)\rho_{+}(\bar{k})\leq 2\rho_{-}(31\bar{k}) holds, then we can let s=31kˉs=31\bar{k}, which implies that

This means that the condition (5) holds, and the corollary follows directly from Theorem II.1. ∎

For the quadratic objective (1), the condition ρ+(kˉ)≤2ρ−(31kˉ)\rho_{+}(\bar{k})\leq 2\rho_{-}(31\bar{k}) is analogous to the RIP condition in . In particular, if the matrix AA has the restricted isometry constant δ31kˉ≤1/3\delta_{31\bar{k}}\leq 1/3, then the condition ρ+(kˉ)≤4/3\rho_{+}(\bar{k})\leq 4/3 and ρ−(31kˉ)≥2/3\rho_{-}(31\bar{k})\geq 2/3 holds, with ρ+(s)\rho_{+}(s) and ρ−(s)\rho_{-}(s) defined according to (4). In this case, Corollary II.1 can be directly applied.

It is interesting to observe that except for constants, the result of this paper for OMP is as strong as those for more sophisticated greedy algorithms such as ROMP or CoSaMP . For example, Corollary II.1 can be applied when δs≤1/3\delta_{s}\leq 1/3 with s=31kˉs=31\bar{k}, while a similar result for CoSaMP in applies when δs≤0.1\delta_{s}\leq 0.1 with s=4kˉs=4\bar{k}. Nevertheless, the difference in the constants may still suggest possible advantages for more complex algorithms such as CoSaMP under suitable conditions.

For quadratic objective function, a simple instantiation of ϵs(xˉ)\epsilon_{s}(\bar{{\mathbf{x}}}) using Proposition II.1 leads to the following sparse recovery result that is relatively simple to interpret.

III Discussion

Our result for signal recovery is stronger than previous results for OMP that relied on different conditions. For example, considered the problem of recovering the support set of a sparse signal under a stronger condition (also see for recovery properties under stochastic noise). A similar analysis was employed in , where it was shown that for any fixed sparse signal xˉ\bar{{\mathbf{x}}} with kˉ=∥xˉ∥0\bar{k}=\|\bar{{\mathbf{x}}}\|_{0}, OMP can recover the signal with large probability using O(kˉln⁡d)O(\bar{k}\ln d) measurements. A more refined analysis in shows that a lower bound of n=2kˉln⁡(d−kˉ)n=2\bar{k}\ln(d-\bar{k}) measurements is enough for recovery. However, the above results are not uniform with respect to all kˉ\bar{k}-sparse signals xˉ\bar{{\mathbf{x}}} (that is, for any set of random projections, there exist kˉ\bar{k}-sparsity signals that fail the analysis). In comparison, the RIP condition holds uniformly by definition, and hence our result applies uniformly to all kˉ\bar{k}-sparse signals. Although our result is stronger than previous results in terms of signal recovery in 2-norm, the result requires running the OMP algorithm for more than kˉ\bar{k} iterations, and hence doesn’t recover the true support set of the ideal signal. In comparison, results such as also imply exact recovery of the correct support set (but under stronger assumptions) using only kˉ\bar{k} OMP iterations. It is also known that it is impossible to uniformly recover the support set (in kˉ\bar{k} iterations) with the OMP algorithm with O(kˉln⁡d)O(\bar{k}\ln d) measurements . This means that it is necessary to run OMP for more than kˉ\bar{k} iterations in order to achieve the best 2-norm recovery performance with as few meausrements as possible.

It is worth mentioning that some previous results apply uniformly to all kˉ\bar{k}-sparse signals. For example, results in depend on the stronger mutual incoherence condition. Unfortunately the mutual incoherence condition can only be satisfied with Ω(kˉ2ln⁡d)\Omega(\bar{k}^{2}\ln d) random projections. Therefore in recent years there have been significant interests in studying OMP under the RIP. In addition to the current paper, a number of recent papers investigated this issue, reaching varying conclusions . For example, the RIP-based analysis for sparse signals (but without noise) was considered in , with the conclusion that under a sufficiently strong assumption on the RIP constant (in fact, the resulting condition is similar to the mutual incoherence condition), exact recovery is possible in kˉ\bar{k} iterations. The condition required for the RIP constant was weakened in , where the author showed that by running the OMP algorithm more than kˉ\bar{k} iterations, it is possible to achieve exact recovery (again assuming no noise). The condition in can be satisfied with only O(kˉ1.6ln⁡d)O(\bar{k}^{1.6}\ln d) measurements, which is a significant improvement over the traditional Ω(kˉ2ln⁡d)\Omega(\bar{k}^{2}\ln d) measurements. The result obtained in the current paper is along the same line as , but reduced the required number of measurements to the optimal order of O(kˉln⁡d)O(\bar{k}\ln d).

It is also interesting to compare the new OMP result in this paper to that of Lasso, which is also known to work under the RIP. However, a more refined comparison illustrates differences between the known theoretical results for these two methods. For OMP, the result in Theorem II.1 can be applied as long as the condition

is satisfied. With F(0)=∅F^{(0)}=\emptyset, this roughly requires (ρ+(1)/ρ−(s))ln⁡(ρ+(kˉ)/ρ−(s))(\rho_{+}(1)/\rho_{-}(s))\ln(\rho_{+}(\bar{k})/\rho_{-}(s)) to grow sub-linearly as a function of ss in order to apply the theory. In comparison, the known condition for Lasso (e.g., this has been made explicit in ) requires ρ+(s)/ρ−(s)\rho_{+}(s)/\rho_{-}(s) to grow sub-linearly as a function of ss. To compare the two conditions, we note that the condition for OMP is weaker in terms of of the upper convexity constant as there is no explicit dependency on ρ+(s)\rho_{+}(s); however, the dependency on ρ−(s)\rho_{-}(s) is stronger in OMP than Lasso due to the logarithmic term. Although it is unclear how tight these conditions are, the comparison nevertheless indicates that even though both algorithms work under the RIP, there are still finer differences in their theoretical analysis: Lasso is slightly more favorable in terms of its dependency on the lower strong convexity constant, while OMP is more favorable in terms of its dependency on the upper strong convexity constant. We further conjecture that the extra logarithmic dependency ln⁡(ρ+(kˉ)/ρ−(s))\ln(\rho_{+}(\bar{k})/\rho_{-}(s)) in OMP is necessary. In practice, some times Lasso performs better while other times OMP performs better (for example, see experimental results in ). Therefore some discrepancy in their theoretical analysis is expected. More specifically, for sparse recovery, one often observes that Lasso is superior when the nonzero coefficients have a similar magnitude (which happens to be the case that the extra ln⁡(ρ+(kˉ)/ρ−(s))\ln(\rho_{+}(\bar{k})/\rho_{-}(s)) factor is required in our OMP analysis) while OMP performs better when the nonzero coefficients exhibit rapid decay (which happens to be the case that the extra ln⁡(ρ+(kˉ)/ρ−(s))\ln(\rho_{+}(\bar{k})/\rho_{-}(s)) factor can be removed from our analysis). The theory in this paper significantly narrows the previous theoretical gap between these two sparse recovery methods by positively answering the open question of whether OMP can recover sparse signals under the RIP. Therefore our result allows practitioners to apply OMP with more confidence than previously expected.

Acknowledgements

The author would like to thank the anonymous referees for pointing out many relevant references and for suggestions to improve the presentation.

Appendix A Technical Lemmas

Let x′=xˉFˉ∩F{\mathbf{x}}^{\prime}=\bar{{\mathbf{x}}}_{\bar{F}\cap F}, then by the definition of x{\mathbf{x}}, we know that Q(x)≤Q(x′)Q({\mathbf{x}})\leq Q({\mathbf{x}}^{\prime}). Therefore

which implies the lemma. The first inequality is by the definitions of ρ+(s)\rho_{+}(s) and ϵs(xˉ)\epsilon_{s}(\bar{{\mathbf{x}}}). The last inequality follows from the fact that ab≤0.5a2+0.5b2ab\leq 0.5a^{2}+0.5b^{2} with a=ϵs(xˉ)/ρ+(s)a=\epsilon_{s}(\bar{{\mathbf{x}}})/\sqrt{\rho_{+}(s)} and b=ρ+(s)∥xˉFˉ∖F∥2b=\sqrt{\rho_{+}(s)}\|\bar{{\mathbf{x}}}_{\bar{F}\setminus F}\|_{2}. ∎

we obtain the desired inequality. The first inequality is by the definitions of ρ−(s)\rho_{-}(s) and ϵs(xˉ)\epsilon_{s}(\bar{{\mathbf{x}}}). The last inequality again follows from the fact that ab≤0.5a2+0.5b2ab\leq 0.5a^{2}+0.5b^{2} with a=ϵs(xˉ)/ρ+(s)a=\epsilon_{s}(\bar{{\mathbf{x}}})/\sqrt{\rho_{+}(s)} and b=ρ+(s)∥xˉFˉ∖F∥2b=\sqrt{\rho_{+}(s)}\|\bar{{\mathbf{x}}}_{\bar{F}\setminus F}\|_{2}. ∎

The next lemma shows that each greedy search makes reasonable progress. This proof is essentially identical to a similar result in but with refined notations used in the current paper. We thus include the proof for completeness. It allows the readers to verify more easily that the proof in remains unchanged with our new definitions.

where j=arg⁡max⁡i∣∇Q(x)i∣j=\arg\max_{i}|\nabla Q({\mathbf{x}})_{i}|.

For all i∈{1,…,d}i\in\{1,\ldots,d\} and η>0\eta>0, we define

It follows from the definition of ρ+(1)\rho_{+}(1) that min_α Q(x+ αe_j) ≤ Q(x+ η sgn( ¯ x _j) e_j) ≤ Q_j(η) . Since the choice of j=arg⁡max⁡i∣∇Q(x)i∣j=\arg\max_{i}|\nabla Q({\mathbf{x}})_{i}| achieves the minimum of min⁡imin⁡ηQi(η)\min_{i}\min_{\eta}Q_{i}(\eta), the lemma is a direct consequence of the following stronger statement:

with an appropriate choice of η\eta; this is because

Therefore, we now turn to prove that (6) holds. Denoting u=∑i∈Fˉ∖F∣xˉi∣u=\sum_{i\in\bar{F}\setminus F}|\bar{{\mathbf{x}}}_{i}|, we obtain that

Since we assume that x{\mathbf{x}} is optimal over FF, we get that ∇Q(x)i=0\nabla Q({\mathbf{x}})_{i}=0 for all i∈Fi\in F. Additionally, xi=0{\mathbf{x}}_{i}=0 for i∉Fi\not\in F and xˉi=0\bar{{\mathbf{x}}}_{i}=0 for i∉Fˉi\not\in\bar{F}. Therefore,

Combining the above with the definition of ρ−(s)\rho_{-}(s), we obtain that

and rearranging the terms, we conclude our proof of (6). ∎

The direct consequence of the previous lemma is the following result, which is critical in our analysis. The idea of using a nesting approximating sequence has appeared in , but the current version is improved. The change is necessary for the purpose of this paper. In the following μ\mu can be chosen as any positive number if L=1L=1.

we have from (8) that if ∣Fˉ1∖Fˉ(k)∣≠0|\bar{F}_{1}\setminus\bar{F}^{(k)}|\neq 0, then

Now assume that the lemma holds at L=m−1L=m-1 for some m>1m>1. That is, with

We thus obtain from (8) that if ∣FˉL∖Fˉ(k)∣≠0|\bar{F}_{L}\setminus\bar{F}^{(k)}|\neq 0, then

The following lemma is a slightly stronger version of the theorem, which we can prove more easily by induction.

Consider the OMP algorithm. If there exist kk and ss such that ∣Fˉ∪F(k)∣≤s|\bar{F}\cup F^{(k)}|\leq s and

We prove this result by induction on ∣Fˉ∖F(0)∣|\bar{F}\setminus F^{(0)}|. If ∣Fˉ∖F(0)∣=0|\bar{F}\setminus F^{(0)}|=0, then the bound in (9) holds trivially because Q(x(k))≤Q(x(0))≤Q(xˉ)Q({\mathbf{x}}^{(k)})\leq Q({\mathbf{x}}^{(0)})\leq Q(\bar{{\mathbf{x}}}).

where μ=10ρ+(m)/ρ−(s)\mu=10\rho_{+}(m)/\rho_{-}(s). We have L≤⌊log⁡2m⌋+1L\leq\lfloor\log_{2}m\rfloor+1 because the second inequality is automatically satisfied when L=⌊log⁡2m⌋+1L=\lfloor\log_{2}m\rfloor+1 (the right hand side is zero in this case). Moreover, if the second inequality is always satisfied for all L≥1L\geq 1, then we can simply take L=1L=1 (and ignore the first inequality).

where (10) is used to derive the second inequality.

then (12) implies that (9) holds automatically (since μ≥10\mu\geq 10), which finishes the induction. Therefore in the following, we only consider the case (13) does not hold, which implies that

Therefore m−∣Fˉ∖F(k)∣+1>2L−1m-|\bar{F}\setminus F^{(k)}|+1>2^{L-1}. That is, ∣Fˉ∖F(k)∣≤m−2L−1|\bar{F}\setminus F^{(k)}|\leq m-2^{L-1}. It follows from the induction hypothesis that after another

OMP iterations, (9) holds. Therefore by combining this estimate with (11), we know that the total number of OMP iterations for (9) to hold (starting with F(0)F^{(0)}) is no more than

This finishes the induction step for the case ∣Fˉ∖F(0)∣=m|\bar{F}\setminus F^{(0)}|=m. ∎

References