pylops.optimization.cls_sparsity.SplitBregman

class pylops.optimization.cls_sparsity.SplitBregman(Op, callbacks=None)

Split Bregman for mixed L2-L1 norms.

Solve an unconstrained system of equations with mixed \(L_2\) and \(L_1\) regularization terms given the operator Op, a list of \(L_1\) regularization terms RegsL1, and an optional list of \(L_2\) regularization terms RegsL2.

Parameters:

Op : pylops.LinearOperator

Operator to invert

Notes

Solve the following system of unconstrained, regularized equations given the operator \(\mathbf{Op}\) and a set of mixed norm (\(L^2\) and \(L_1\)) regularization terms \(\mathbf{R}_{2,i}\) and \(\mathbf{R}_{1,i}\), respectively:

\[J = \frac{\mu}{2} \|\textbf{y} - \textbf{Op}\,\textbf{x} \|_2^2 + \frac{1}{2}\sum_i \epsilon_{\mathbf{R}_{2,i}} \|\mathbf{y}_{\mathbf{R}_{2,i}} - \mathbf{R}_{2,i} \textbf{x} \|_2^2 + \sum_i \epsilon_{\mathbf{R}_{1,i}} \| \mathbf{R}_{1,i} \textbf{x} \|_1\]

where \(\mu\) is the reconstruction damping, \(\epsilon_{\mathbf{R}_{2,i}}\) are the damping factors used to weight the different \(L^2\) regularization terms of the cost function and \(\epsilon_{\mathbf{R}_{1,i}}\) are the damping factors used to weight the different \(L_1\) regularization terms of the cost function.

The generalized Split-Bergman algorithm [1] is used to solve such cost function: the algorithm is composed of a sequence of unconstrained inverse problems and Bregman updates.

The original system of equations is initially converted into a constrained problem:

\[J = \frac{\mu}{2} \|\textbf{y} - \textbf{Op}\,\textbf{x}\|_2^2 + \frac{1}{2}\sum_i \epsilon_{\mathbf{R}_{2,i}} \|\mathbf{y}_{\mathbf{R}_{2,i}} - \mathbf{R}_{2,i} \textbf{x}\|_2^2 + \sum_i \| \textbf{y}_i \|_1 \quad \text{subject to} \quad \textbf{y}_i = \mathbf{R}_{1,i} \textbf{x} \quad \forall i\]

and solved as follows:

\[\begin{split}\DeclareMathOperator*{\argmin}{arg\,min} \begin{align} (\textbf{x}^{k+1}, \textbf{y}_i^{k+1}) = \argmin_{\mathbf{x}, \mathbf{y}_i} \|\textbf{y} - \textbf{Op}\,\textbf{x}\|_2^2 &+ \frac{1}{2}\sum_i \epsilon_{\mathbf{R}_{2,i}} \|\mathbf{y}_{\mathbf{R}_{2,i}} - \mathbf{R}_{2,i} \textbf{x}\|_2^2 \\ &+ \frac{1}{2}\sum_i \epsilon_{\mathbf{R}_{1,i}} \|\textbf{y}_i - \mathbf{R}_{1,i} \textbf{x} - \textbf{b}_i^k\|_2^2 \\ &+ \sum_i \| \textbf{y}_i \|_1 \end{align}\end{split}\]

\[\textbf{b}_i^{k+1}=\textbf{b}_i^k + (\mathbf{R}_{1,i} \textbf{x}^{k+1} - \textbf{y}^{k+1})\]

The scipy.sparse.linalg.lsqr solver and a fast shrinkage algorithm are used within a inner loop to solve the first step. The entire procedure is repeated niter_outer times until convergence.

[1]	Goldstein T. and Osher S., “The Split Bregman Method for L1-Regularized Problems”, SIAM J. on Scientific Computing, vol. 2(2), pp. 323-343. 2008.

Methods

__init__(Op[, callbacks])	Initialize self.
callback(x, *args, **kwargs)	Callback routine
finalize([show])	Finalize solver
run(x[, show, itershow, show_inner])	Run solver
setup(y, RegsL1[, x0, niter_outer, …])	Setup solver
solve(y, RegsL1[, x0, niter_outer, …])	Run entire solver
step(x[, show, show_inner])	Run one step of solver

setup(y, RegsL1, x0=None, niter_outer=3, niter_inner=5, RegsL2=None, dataregsL2=None, mu=1.0, epsRL1s=None, epsRL2s=None, tol=1e-10, tau=1.0, restart=False, show=False)

Setup solver

Parameters:

y : np.ndarray

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

RegsL1 : list

\(L_1\) regularization operators

x0 : np.ndarray, optional

Initial guess of size \([M \times 1]\). If None, initialize internally as zero vector

niter_outer : int, optional

Number of iterations of outer loop

niter_inner : int, optional

Number of iterations of inner loop of first step of the Split Bregman algorithm. A small number of iterations is generally sufficient and for many applications optimal efficiency is obtained when only one iteration is performed.

RegsL2 : list, optional

Additional \(L^2\) regularization operators (if None, \(L^2\) regularization is not added to the problem)

dataregsL2 : list, optional

\(L^2\) Regularization data (must have the same number of elements of RegsL2 or equal to None to use a zero data for every regularization operator in RegsL2)

mu : float, optional

Data term damping

epsRL1s : list

\(L_1\) Regularization dampings (must have the same number of elements as RegsL1)

epsRL2s : list

\(L^2\) Regularization dampings (must have the same number of elements as RegsL2)

tol : float, optional

Tolerance. Stop outer iterations if difference between inverted model at subsequent iterations is smaller than tol

tau : float, optional

Scaling factor in the Bregman update (must be close to 1)

restart : bool, optional

Initialize the unconstrained inverse problem in inner loop with the initial guess (True) or with the last estimate (False). Note that when this is set to True, the x0 provided in the setup will be used in all iterations.

show : bool, optional

Display setup log

Returns:

x : np.ndarray

Initial guess of size \([N \times 1]\)

step(x, show=False, show_inner=False, **kwargs_lsqr)

Run one step of solver

Parameters:

x : list or np.ndarray

Current model vector to be updated by a step of OMP

show_inner : bool, optional

Display inner iteration logs of lsqr

show : bool, optional

Display iteration log

**kwargs_lsqr

Arbitrary keyword arguments for scipy.sparse.linalg.lsqr solver used to solve the first subproblem in the first step of the Split Bregman algorithm.

Returns:

x : np.ndarray

Updated model vector

run(x, show=False, itershow=[10, 10, 10], show_inner=False, **kwargs_lsqr)

Run solver

Parameters:

x : np.ndarray

Current model vector to be updated by multiple steps of IRLS

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.

show_inner : bool, optional

Display inner iteration logs of lsqr

**kwargs_lsqr

Arbitrary keyword arguments for scipy.sparse.linalg.lsqr solver used to solve the first subproblem in the first step of the Split Bregman algorithm.

Returns:

x : np.ndarray

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

finalize(show=False)

Finalize solver

Parameters:

show : bool, optional

Display finalize log

Returns:

xfin : np.ndarray

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

solve(y, RegsL1, x0=None, niter_outer=3, niter_inner=5, RegsL2=None, dataregsL2=None, mu=1.0, epsRL1s=None, epsRL2s=None, tol=1e-10, tau=1.0, restart=False, show=False, itershow=[10, 10, 10], show_inner=False, **kwargs_lsqr)

Run entire solver

Parameters:

y : np.ndarray

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

RegsL1 : list

\(L_1\) regularization operators

x0 : np.ndarray, optional

Initial guess of size \([M \times 1]\). If None, initialize internally as zero vector

niter_outer : int, optional

Number of iterations of outer loop

niter_inner : int, optional

Number of iterations of inner loop of first step of the Split Bregman algorithm. A small number of iterations is generally sufficient and for many applications optimal efficiency is obtained when only one iteration is performed.

RegsL2 : list, optional

Additional \(L^2\) regularization operators (if None, \(L^2\) regularization is not added to the problem)

dataregsL2 : list, optional

\(L^2\) Regularization data (must have the same number of elements of RegsL2 or equal to None to use a zero data for every regularization operator in RegsL2)

mu : float, optional

Data term damping

epsRL1s : list

\(L_1\) Regularization dampings (must have the same number of elements as RegsL1)

epsRL2s : list

\(L^2\) Regularization dampings (must have the same number of elements as RegsL2)

tol : float, optional

Tolerance. Stop outer iterations if difference between inverted model at subsequent iterations is smaller than tol

tau : float, optional

Scaling factor in the Bregman update (must be close to 1)

restart : bool, optional

Initialize the unconstrained inverse problem in inner loop with the initial guess (True) or with the last estimate (False). Note that when this is set to True, the x0 provided in the setup will be used in all iterations.

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.

show_inner : bool, optional

Display inner iteration logs of lsqr

**kwargs_lsqr

Arbitrary keyword arguments for scipy.sparse.linalg.lsqr solver used to solve the first subproblem in the first step of the Split Bregman algorithm.

Returns:

x : np.ndarray

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

iiter : int

Iteration number of outer loop upon termination

cost : numpy.ndarray

History of cost function through iterations