pylops.optimization.cls_sparsity.OMP

class pylops.optimization.cls_sparsity.OMP(Op, callbacks=None)[source]

Orthogonal Matching Pursuit (OMP).

Solve an optimization problem with \(L_0\) regularization function given the operator Op and data y. The operator can be real or complex, and should ideally be either square \(N=M\) or underdetermined \(N<M\).

Parameters:
Op : pylops.LinearOperator

Operator to invert

See also

ISTA
Iterative Shrinkage-Thresholding Algorithm (ISTA).
FISTA
Fast Iterative Shrinkage-Thresholding Algorithm (FISTA).
SPGL1
Spectral Projected-Gradient for L1 norm (SPGL1).
SplitBregman
Split Bregman for mixed L2-L1 norms.

Notes

Solves the following optimization problem for the operator \(\mathbf{Op}\) and the data \(\mathbf{y}\):

\[\|\mathbf{x}\|_0 \quad \text{subject to} \quad \|\mathbf{Op}\,\mathbf{x}-\mathbf{y}\|_2^2 \leq \sigma^2,\]

using Orthogonal Matching Pursuit (OMP). This is a very simple iterative algorithm which applies the following step:

\[\begin{split}\DeclareMathOperator*{\argmin}{arg\,min} \DeclareMathOperator*{\argmax}{arg\,max} \Lambda_k = \Lambda_{k-1} \cup \left\{\argmax_j \left|\mathbf{Op}_j^H\,\mathbf{r}_k\right| \right\} \\ \mathbf{x}_k = \argmin_{\mathbf{x}} \left\|\mathbf{Op}_{\Lambda_k}\,\mathbf{x} - \mathbf{y}\right\|_2^2\end{split}\]

Note that by choosing niter_inner=0 the basic Matching Pursuit (MP) algorithm is implemented instead. In other words, instead of solving an optimization at each iteration to find the best \(\mathbf{x}\) for the currently selected basis functions, the vector \(\mathbf{x}\) is just updated at the new basis function by taking directly the value from the inner product \(\mathbf{Op}_j^H\,\mathbf{r}_k\).

In this case it is highly recommended to provide a normalized basis function. If different basis have different norms, the solver is likely to diverge. Similar observations apply to OMP, even though mild unbalancing between the basis is generally properly handled.

Methods

__init__(Op[, callbacks]) Initialize self.
callback(x, *args, **kwargs) Callback routine
finalize(x, cols[, show]) Finalize solver
run(x, cols[, show, itershow]) Run solver
setup(y[, niter_outer, niter_inner, sigma, …]) Setup solver
solve(y[, niter_outer, niter_inner, sigma, …]) Run entire solver
step(x, cols[, show]) Run one step of solver
setup(y, niter_outer=10, niter_inner=40, sigma=0.0001, normalizecols=False, show=False)[source]

Setup solver

Parameters:
y : np.ndarray

Data of size \([N \times 1]\)

niter_outer : int, optional

Number of iterations of outer loop

niter_inner : int, optional

Number of iterations of inner loop. By choosing niter_inner=0, the Matching Pursuit (MP) algorithm is implemented.

sigma : list

Maximum \(L^2\) norm of residual. When smaller stop iterations.

normalizecols : list, optional

Normalize columns (True) or not (False). Note that this can be expensive as it requires applying the forward operator \(n_{cols}\) times to unit vectors (i.e., containing 1 at position j and zero otherwise); use only when the columns of the operator are expected to have highly varying norms.

show : bool, optional

Display setup log

step(x, cols, show=False)[source]

Run one step of solver

Parameters:
x : list or np.ndarray

Current model vector to be updated by a step of OMP

cols : list

Current list of chosen elements of vector x to be updated by a step of OMP

show : bool, optional

Display iteration log

Returns:
x : np.ndarray

Updated model vector

cols : list

Current list of chosen elements

run(x, cols, show=False, itershow=[10, 10, 10])[source]

Run solver

Parameters:
x : np.ndarray

Current model vector to be updated by multiple steps of IRLS

cols : list

Current list of chosen elements of vector x to be updated by a step of OMP

show : bool, optional

Display logs

itershow : list, optional

Display set log for the first N1 steps, last N2 steps, and every N3 steps in between where N1, N2, N3 are the three element of the list.

Returns:
x : np.ndarray

Estimated model of size \([M \times 1]\)

cols : list

Current list of chosen elements

finalize(x, cols, show=False)[source]

Finalize solver

Parameters:
x : list or np.ndarray

Current model vector to be updated by a step of OMP

cols : list

Current list of chosen elements of vector x to be updated by a step of OMP

show : bool, optional

Display finalize log

Returns:
xfin : np.ndarray

Estimated model of size \([M \times 1]\)

solve(y, niter_outer=10, niter_inner=40, sigma=0.0001, normalizecols=False, show=False, itershow=[10, 10, 10])[source]

Run entire solver

Parameters:
y : np.ndarray

Data of size \([N \times 1]\)

niter_outer : int, optional

Number of iterations of outer loop

niter_inner : int, optional

Number of iterations of inner loop. By choosing niter_inner=0, the Matching Pursuit (MP) algorithm is implemented.

sigma : list

Maximum \(L^2\) norm of residual. When smaller stop iterations.

normalizecols : list, optional

Normalize columns (True) or not (False). Note that this can be expensive as it requires applying the forward operator \(n_{cols}\) times to unit vectors (i.e., containing 1 at position j and zero otherwise); use only when the columns of the operator are expected to have highly varying norms.

show : bool, optional

Display logs

itershow : list, optional

Display set log for the first N1 steps, last N2 steps, and every N3 steps in between where N1, N2, N3 are the three element of the list.

Returns:
x : np.ndarray

Estimated model of size \([M \times 1]\)