*Over the next few days, I am posting portions of a paper I am currently writing about making greedy algorithms as efficient as possible. I post this material before submission because: 1) much of it has been discussed before in a variety of places, and I am merely assembling some of these methods in one place; 2) I think it can find immediate use, rather than wait to release the entire paper after the paper deadline in February. Today, we begin with some notation, and a look at the naïve implementation of orthogonal matching pursuit (OMP). This provides us a baseline complexity.*

We are interested in efficiently modeling a signal \(\vu\)

by a linear combination of the atoms defined in a dictionary

\(\mathcal{D} = \{ \vphi_\omega \in C^M \}_{\omega \in \Omega= \{1, 2, \ldots, N\}}\):

$$

\vu = \MPhi \vx + \vn

$$

where \(\MPhi := [\vphi_1 | \vphi_2 | \cdots | \vphi_N]\).

By \(\Omega_k \subset \Omega\) we denote

an ordered subset of cardinality \(k\) of indices into the dictionary.

We denote the \(i\)th element of a set by \(\mathcal{I}(i)\).

The matrix \(\MPhi_{\Omega_k}\) is thus composed of

the ordered atoms indexed by \(\Omega_k\).

We denote \(\MP_k\) the orthogonal projection matrix onto

the range space of \(\MPhi_{\Omega_k}\), or equivalently,

onto the span of the subdictionary \(\mathcal{D}_k \subset \mathcal{D}\).

The projection matrix \(\MP_k^\perp = \MI – \MP_k\) is thus the

orthogonal projection onto the subspace orthogonal to

the span of the subdictionary,

or equivalently, the left null space of \(\MPhi_{\Omega_k}\).

We denote the inner product between two vectors in \(C^M\)

as \(\langle \vu_1, \vu_2 \rangle = \vu_2^H \vu_1\),

where \(^H\) denotes the complex conjugate transpose.

The vector \(\ve_k\) is the \(k\)th element of the standard basis of \(R^N\),

i.e., all zeros with a 1 in its \(k\)th row.

Finally, all norms are implicitly the \(\ell_2\)-norm,

\(\|\vx\|^2 = \langle \vx, \vx \rangle\).

In its \(k\)th iteration, a naïve implementation of OMP

creates and searches through the set

\(\mathcal{I}_{k-1} = \{\langle \vr_{k-1}, \vphi_n \rangle/\|\vphi_n\| \}_{n \in \Omega}\),

which has a complexity of \(\mathcal{O}(NM)\).

After selecting a new atom,

OMP orthogonalizes the residual by evaluating

\(\vr_{k} = \MP_k^\perp\vu\),

which has complexity of at most \(\mathcal{O}(Mk + Mk^2 + k^3)\)

(because of evaluating the pseudoinverse).

Thus, the \(k\)th iteration of naïve OMP has complexity

\(\mathcal{O}(NM + Mk + Mk^2 + k^3)\),

which shows that as the model becomes larger,

the orthgonalization procedure dominates the complexity.

Here is my code: ompdnaive.m

How the orthogonalization can be done, when the input signal and the atoms are 2D?

LikeLike

Hello. Why not just vectorize them?

LikeLike