# Sparse regularization of inverse problems by biorthogonal frame thresholding

Jürgen Frikel  Department of Computer Science and Mathematics, OTH Regensburg
Galgenbergstraße 32, D-93053 Regensburg, Germany
<EMAIL_ADDRESS>Markus Haltmeier  Department of Mathematics,
University of Innsbruck
Technikerstrasse 13, 6020 Innsbruck, Austria
<EMAIL_ADDRESS>
###### Abstract

We analyze sparse frame based regularization of inverse problems by means of a
biorthogonal frame decomposition (BFD) for the forward operator, which
generalizes the SVD. The BFD allows to define a non-iterative (direct) frame
thresholding approach which we show to provide a convergent regularization
method with linear convergence rates. These results will be compared to the
well-known analysis and synthesis variants of sparse $\ell^{1}$-regularization
which are usually implemented thorough iterative schemes. If the frame is a
basis (non-redundant case), the three versions of sparse regularization,
namely synthesis and analysis variants of $\ell^{1}$-regularization as well as
the BFD thresholding are equivalent. However, in the redundant case, those
three approaches are pairwise different.

##  1 Introduction

This paper is concerned with inverse problems of the form

$\bm{y}^{\delta}=\mathbf{A}\bm{x}+z\,,$ (1.1)

where $\mathbf{A}\colon\mathbb{D}(\mathbf{A})\subseteq\mathbb{X}\to\mathbb{Y}$
is a linear operator between Hilbert spaces, and $z$ denotes the data
distortion (noise). We allow unbounded operators and assume that
$\mathbb{D}(\mathbf{A})$ is dense. Moreover, we assume that the unknown object
$\bm{x}$ is an element of a closed subspace space
$\mathbb{X}_{0}\subseteq\mathbb{X}$ on which $\mathbf{A}$ is bounded. We are
particularly interested in problems, where (1.1) is ill-posed in which case
the solution of (1.1) (if existent) is either not unique or the solution
operator is not continuous (hence, the solution process is unstable with
respect to data perturbations). In order to stabilize the inversion of (1.1)
one has to apply regularization methods, cf. [11, 23]. The basic idea of
regularization is to include a-priori information about the unknown object
into the solution process.

In this paper, we use sparsity based regularization, where the a-priori
assumption on the unknown object is sparsity of $\bm{x}$ with respect to a
frame $(u_{\lambda})_{\lambda\in\Lambda}$ of $\mathbb{X}$, cf. [23, 14, 7, 17,
21, 3]. That is, we regularize the recovery of $\bm{x}$ from measurements
(1.1) by enforcing sparsity of $\bm{x}$ with respect to a suitably chosen
frame of $\mathbb{X}$. Sparse regularization is well investigated and has been
applied to many different imaging problems, and by now there are many
algorithms available that implement sparse regularization. However, when
dealing with frames, there are at least two fundamentally different concepts
implementing sparsity, namely the synthesis and the analysis variant. The
reason for this lies in the fact that expansions of $\bm{x}\in\mathbb{X}$ with
respect to frames are not unique (which is in contrast to basis expansions).
In the synthesis variant, it is assumed that the unknown is a sparse linear
combination of frame elements, whereas, in the analysis variant, it is
required that a the inner products $\left\langle
u_{\lambda},\bm{x}\right\rangle$ with respect to a given frame are sparse. The
difference between these approaches has been pointed out clearly in [10].

Sparse regularization is widely used in inverse problems as it provides good
regularization results and is able to preserve or emphasize features (e.g.,
edges) in the reconstruction. However, this often comes at the price of speed,
since most of the algorithms implementing sparse regularization are based on
variational formulations that are solved by iterative schemes.

In the present paper, we investigate a third variant of sparse regularization
that is based on an operator adapted biorthogonal frame decomposition (BFD) of
the unknown object, cf. [9, 4, 6] and which allows to define a direct (non-
iterative) sparse regularization method. In the noise-free case ($z=0$),
explicit reproducing formulas for the unknown object can be derived from the
BFD, where the frame coefficients of $\bm{x}$ are calculated directly from the
data $\bm{y}=\mathbf{A}\bm{x}$. In the presence of noise ($z\neq 0$),
regularized versions of those formulas are obtained by applying component-wise
soft-thresholding to the calculated coefficients, where the soft-thresholding
operator is defined as follows:

######  Definition 1.1 (Soft-thresholding).

Let $\Lambda$ be some index set.

  * $\blacksquare$

For $\eta,d\in\mathbb{K}$ let
$\operatorname{soft}(\eta,d)\coloneqq\operatorname{sign}(\eta)\max{\left\\{0,\left|\eta\right|-d\right\\}}$.

  * $\blacksquare$

For $\bm{\eta},\bm{d}\in\mathbb{K}^{\Lambda}$ we define the component-wise
soft-thresholding by

$\mathbb{S}_{\bm{d}}(\bm{\eta})\coloneqq\left(\operatorname{soft}\left(\eta_{\lambda},d_{\lambda}\right)\right)_{\lambda\in\Lambda}\,.$
(1.2)

Here and below we define
$\operatorname{sign}(\eta)\coloneqq\eta/\left|\eta\right|$ for
$\eta\in\mathbb{K}\setminus\\{0\\}$ and $\operatorname{sign}(0)\coloneqq 0$.
The advantage of the BFD-variant of sparse regularization lies in the fact
that it admits a non-iterative (direct) and fast implementation of sparse
regularization that can be easily implemented for several inverse problems.

We point out, that the three variants of sparse regularization (mentioned
above) are equivalent if orthonormal bases are used instead of frames, but
they are fundamentally different in the redundant case.

As the main theoretical results in this paper we show that the thrid variant
of sparse regularization, which we call _BFD-thresholding_ , defines a
convergent regularization method and we derive linear convergence rates for
sparse solutions. For the basis case, the same results follow from existing
results of $\ell^{1}$-regularization [7, 14, 15]. In the redundant case, the
results follow from [15] for the synthesis approach and from [16] for the
analysis approach. In case of BFD-thresholding, we are not aware of any
results concerning convergence analysis or convergence rates.

### Outline

This paper is organized as follows. In Section 2 we define the biorthogonal
frame decompositions of operators and give several examples of biorthogonal
frame expansions for various operators using wavelet, curvelet, and shearlet
frames. In Section 3 we review the convergence theory of
$\ell^{1}$-regularization and the convergence rates. In Section 4, we show
that BFD-thresholding is a convergent regularization method and derive its
convergence rates.

##  2 Biorthogonal frame decomposition

In this section, we introduce the concept of operator adapted biorthogonal
frame decompositions (BFD) and discuss some classical examples of such BFDs in
the case of the classical 2D Radon transform and the forward operator of
photoacoustic tomography with a flat observation surface.

###  2.1 Formal definition

We define the operator adapted biorthogonal frame decomposition as a
generalization of the wavelet vaguelette decomposition and the biorthogonal
curvelet or shearlet decompositions to general frames, cf. [9, 4, 6].

######  Definition 2.1 (Biorthogonal frame decomposition (BFD)).

Let $(u_{\lambda})_{\lambda\in\Lambda}\in\mathbb{X}^{\Lambda}$,
$(v_{\lambda})_{\lambda\in\Lambda}\in\mathbb{Y}^{\Lambda}$ and let
$(\kappa_{\lambda})_{\lambda\in\Lambda}$ be a family of positive numbers. For
a linear operator
$\mathbf{A}\colon\mathbb{D}(\mathbf{A})\subseteq\mathbb{X}\to\mathbb{Y}$, we
call $(u_{\lambda},v_{\lambda},\kappa_{\lambda})_{\lambda\in\Lambda}$ a
_biorthogonal frame decomposition (BFD) for $\mathbf{A}$_, if the following
conditions hold:

  1. (D1)

$(u_{\lambda})_{\lambda}$ is a frame of $\mathbb{X}$,

  2. (D2)

$(v_{\lambda})_{\lambda}$ is a frame of
$\overline{\operatorname{ran}(\mathbf{A})}=\overline{\mathbf{A}(\mathbb{X})}$,

  3. (D3)

$\forall\lambda\in\Lambda\colon\mathbf{A}^{*}v_{\lambda}=\kappa_{\lambda}u_{\lambda}$.

######  Remark 2.2.

The BFD generalizes the singular value decomposition (SVD) and the wavelet-
vaguelette decomposition (WVD) (cf. [9]) as it allows the systems
$(u_{\lambda})_{\lambda}$ and $(v_{\lambda})_{\lambda}$ to be non-orthogonal
and redundant. Note that by (D2) and (D3), the frame $(u_{\lambda})_{\lambda}$
satisfies
$u_{\lambda}\in\mathbf{A}^{*}(\overline{\operatorname{ran}(\mathbf{A})})=\operatorname{ran}(\mathbf{A}^{*})$,
where we have made use of the identity
$\operatorname{ran}(\mathbf{A})^{\bot}=\ker(\mathbf{A}^{*})$. For typical
inverse problems this yields a notable smoothness assumption on the involved
elements of the frame $(u_{\lambda})_{\lambda}$.

Although the SVD has proven itself to be a useful tool for analyzing and
solving inverse problems it has the following drawbacks: First, ONBs that are
provided by the SVD (though optimally adapted to the operator in
consideration), in many cases, don’t provide sparse representations of signals
of interest and, hence, are not suitable for the use in sparse regularization.
In particular, frames that provide sparse representation of signals (such as
wavelets or wavelet-like systems) are often not part of the SVD. Second, SVD
is often very hard to compute and not known analytically for many practical
applications.

To overcome some of those difficulties, wavelet-vaguelette decompositions were
introduced. Nevertheless, this concept builds upon expansions of signals with
respect to orthogonal wavelet-systems, which may not provide an optimal sparse
representation of signals of interest, e.g., signals with sharp edges. Thus,
by allowing general frames, the BCD offers great flexibility in the choice of
a suitable function system for sparse regularization while retaining the
advantages.

######  Definition 2.3.

For a frame $(u_{\lambda})_{\lambda\in\Lambda}$ of $\mathbb{X}$, the
_synthesis operator_ is defined as
$\mathbf{U}\colon\ell^{2}(\Lambda)\to\mathbb{Y}\colon\bm{\xi}\mapsto\sum_{\lambda\in\Lambda}\xi_{\lambda}u_{\lambda}$
and the corresponding _analysis operator_ is defined as its adjoint,
$\mathbf{U}^{\ast}\colon\mathbb{X}\to\ell^{2}(\Lambda)\colon\bm{x}\mapsto(\left\langle\bm{x},u_{\lambda}\right\rangle)_{\lambda\in\Lambda}$.

In what follows, the synthesis operator of a frame will be always denoted with
the corresponding upper case letter, e.g., if
$(v_{\lambda})_{\lambda\in\Lambda}$ is a frame, then $\mathbf{V}$ denotes the
corresponding synthesis and $\mathbf{V}^{\ast}$ the analysis operator.

In order to simplify the notation, we will also refer to
$(\mathbf{U},\mathbf{V},\bm{\kappa})$ the as BFD instead of using the full
notation $(u_{\lambda},v_{\lambda},\kappa_{\lambda})_{\lambda\in\Lambda}$.

If a BFD exists for an operator $\mathbf{A}$, it immediately gives rise to a
_reproducing formula_

$\bm{x}=\sum_{\lambda\in\Lambda}\left\langle\bm{x},u_{\lambda}\right\rangle\bar{u}_{\lambda}=\sum_{\lambda\in\Lambda}\kappa_{\lambda}^{-1}\left\langle\mathbf{A}\bm{x},v_{\lambda}\right\rangle\bar{u}_{\lambda}=\bar{\mathbf{U}}\circ\mathbf{M}^{\bm{\texttt{+}}}_{\bm{\kappa}}\circ\mathbf{V}^{\ast}(\mathbf{A}\bm{x}),$
(2.1)

where $(\bar{u}_{\lambda})_{\lambda\in\Lambda}$ is the dual frame to
$(u_{\lambda})_{\lambda\in\Lambda}$ (cf. [5]) and $\bar{\mathbf{U}}$ the
corresponding synthesis operator. Moreover,
$\mathbf{M}^{\bm{\texttt{+}}}_{\bm{\kappa}}$ denotes the Moore-Penrose inverse
of $\mathbf{M}_{\bm{\kappa}}$ and performs point-wise division with
$\bm{\kappa}$, i.e.

$(\mathbf{M}^{\bm{\texttt{+}}}_{\bm{\kappa}}(\bm{\eta}))_{\lambda}\coloneqq\begin{cases}(\eta_{\lambda}/\kappa_{\lambda})_{\lambda\in\Lambda}&\text{
if }\kappa_{\lambda}\neq 0\\\ 0&\text{ otherwise }\,.\end{cases}$ (2.2)

Hence, from given (clean) data $\bm{y}$, one can calculate the frame
coefficient of $\bm{x}$ and obtain a reconstruction via (2.1). The key to the
practical use of this reproducing formulas is the efficient implementation of
the analysis and synthesis operators $\mathbf{V}^{\ast}$ and
$\bar{\mathbf{U}}$, respectively. For particular cases, we will provide
efficient and easy to implement algorithms for the evaluation of
$\mathbf{V}^{\ast}$ and $\bar{\mathbf{U}}$.

Note that, the reproducing formula (2.1) (similarly to the SVD) reveals the
ill-posedness of the operator equation through the decay of the _quasi-
singular values_ $\kappa_{\lambda}$. A regularized version of the reproducing
formula (2.1) can be obtained by incorporating soft-thresholding of the frame
coefficients. In section 4, we present a complete analysis of this approach as
we are not aware of any results for the general BFD-thresholding in the
context of regularization theory.

We now provide several examples of BFDs, including the wavelet vaguelette
decomposition and the biorthogonal curvelet decomposition [4] for the Radon
transform as well as a BFD for the forward operator of photoacoustic
tomography with flat observation surface.

###  2.2 Radon transform

######  Definition 2.4 (Radon transform).

The Radon transform $\mathbf{R}\colon L^{2}(B_{1}(0))\to
L^{2}(\mathbb{S}^{1}\times\mathbb{R})$ is defined by

$\forall(\theta,s)\in\mathbb{S}^{1}\times\mathbb{R}\colon\quad\mathbf{R}f(\theta,s)=\int_{\mathbb{R}}f(s\theta+t\theta^{\bot})\,\mathrm{d}s\,.$
(2.3)

It is well known that the Radon transform is bounded on $L^{2}(B_{1}(0))$, see
[20]. Let $\mathbf{F}g(\theta,\omega)=\int_{\mathbb{R}}g(\theta,s)e^{-i\omega
s}$ be the Fourier transform with respect to the first component and consider
the Riesz potential [20]

$(\mathbf{I}^{-\alpha}g)(\theta,\omega)\coloneqq\frac{1}{2\pi}\int_{\mathbb{R}}\left|\omega\right|^{\alpha}(\mathbf{F}g)(\theta,\omega)e^{i\omega
s}\,\mathrm{d}s$ (2.4)

for $\alpha>-1$. The following hold:

  1. (R1)

The commutation relation
$(\mathbf{I}^{-\alpha}\circ\mathbf{R})f=(\mathbf{R}\circ(-\Delta)^{\alpha/2})f$.

  2. (R2)

The filtered backprojection formula
$f=(4\pi)^{-1}\mathbf{R}^{*}(\mathbf{I}^{-1}\circ\mathbf{R})f=:\mathbf{R}^{\sharp}\mathbf{R}f$.

  3. (R3)

Isometry property $4\pi\left\langle
f_{1},f_{2}\right\rangle_{L^{2}}=\langle\mathbf{I}^{-1}\circ\mathbf{R}f_{1},\mathbf{R}f_{2}\rangle_{L^{2}}$.

Using these ingredients, one can obtain a BFD for the Radon transform as
follows:

######  Example 2.5.

BFD for the Radon transform Let $(u_{\lambda})_{\lambda\in\Lambda}$ be either
an orthonormal basis of wavelets with compact support, a (band-limited)
curvelet or a shearlet tight frame with $\lambda=(j,k,\beta)\in\Lambda$ where
$j\geq 0$ is the scale index. Then $(\mathbf{U},\mathbf{V},\bm{\kappa})$ is a
BFD with

$\displaystyle v_{\lambda}$ $\displaystyle\coloneqq
2^{-j/2}\,(4\pi)^{-1}(\mathbf{I}^{-1}\circ\mathbf{R})u_{\lambda}$ (2.5)
$\displaystyle\kappa_{\lambda}$ $\displaystyle\coloneqq 2^{-j/2}\,.$ (2.6)

These results have been obtained in [9] for wavelet bases, in [4] for curvelet
systems and in [6] for the shearlet frame. All cases are shown in similar
manner and basically follow from (R1), (R2) and the fact that that
$2^{-j/2}(-\Delta)^{1/4}(4\pi)^{-1}u_{\lambda}\simeq u_{\lambda}$ for any of
the considered systems. The limited data case has been studied in [12].

Equation (2.5) implies

$\left\langle g,v_{\lambda}\right\rangle=2^{-j/2}(4\pi)^{-1}\left\langle
g,\mathbf{I}^{-1}\circ\mathbf{R}u_{\lambda}\right\rangle=2^{-j/2}(4\pi)^{-1}\left\langle\mathbf{R}^{*}\circ\mathbf{I}^{-1}g,u_{\lambda}\right\rangle\,.$
(2.7)

This gives an efficient numerical algorithm for the evaluation of
$\mathbf{V}^{\ast}$ provided that $\mathbf{U}^{\ast}$ is associated with an
efficient algorithm. This is in particular the case for the wavelet, shearlet
and curvelet frames as above.

######  Remark 2.6.

We would like to note that in order to define a BFD for the case of curvelets
or shearlets one needs to consider the Radon transform on subspaces of
$L^{2}(\mathbb{R}^{2})$ consisting of functions that are defined on unbounded
domains (since band-limited curvelets or shearlets have non-compact support).
However, since the Radon transform is an unbounded operator on
$L^{2}(\mathbb{R}^{2})$, the reproducing formula (2.1) will not hold any more
in general. The reproducing formulas are at least available for the case that
the object $\bm{x}$ can be represented as a finite linear combination of
curvelets or shearlets (cf. [4] and [6]).

Another possibility would be to consider projections of curvelet or shearlet
frames onto the space $L^{2}(B_{1}(0))$, which would yield a frame for this
space (cf. [5]) and then define the BFD in the same way as above. Because the
Radon transform is continuous on $L^{2}(B_{1}(0))$, the reproducing formula
(2.1) will hold for general linear combinations.

######  Algorithm 2.7.

Computing BFD coefficients for the Radon transform Let $\mathbf{U}$ be a
wavelet, shearlet of curvelet frame and define $\mathbf{V}$ by (2.5).

  1. (a)

Input: $g\in L^{2}(\mathbb{S}^{1}\times\mathbb{R})$.

  2. (b)

Compute $f_{\rm FBP}\coloneqq(4\pi)^{-1}\mathbf{R}^{*}\circ\mathbf{I}^{-1}g$.

  3. (c)

Compute $\bm{\eta}\coloneqq\mathbf{U}^{*}f_{\rm FBP}$ via wavelet, curvelet or
shearlet transform.

  4. (d)

Apply rescaling
$\bm{\eta}\leftarrow(2^{-j/2}\eta_{\lambda})_{\lambda\in\Lambda}$.

  5. (e)

Output: Coefficients $\bm{\eta}$.

###  2.3 Inversion of the wave equation

We consider a planar geometry, which has been considered in our previous work
[13]. Let $C_{0}^{\infty}(\mathbb{H}_{+})$ denote the space of compactly
supported functions $h\colon\mathbb{R}^{2}\to\mathbb{R}$ that are supported in
the half space $\mathbb{H}_{+}\coloneqq\mathbb{R}\times(0,\infty)$. For $f\in
C_{0}^{\infty}(\mathbb{H}_{+})$ consider the initial value problem

$\displaystyle(\partial_{t}^{2}-\Delta u(x,y,t)$ $\displaystyle=0,$
$\displaystyle(x,y,t)\in\mathbb{R}^{2}\times\mathbb{R}$ (2.8) $\displaystyle
u(x,y,0)$ $\displaystyle=y^{1/2}f(x,y)$ $\displaystyle(x,y)\in\mathbb{R}^{2}$
$\displaystyle\partial_{t}u(x,y,0)$ $\displaystyle=0$
$\displaystyle(x,y)\in\mathbb{R}^{2}\,.$

The trace map $\mathbf{A}\colon f\mapsto t^{-1/2}g$ where $g(x,t)\coloneqq
u(x,y=0,t)\chi{\left\\{t\geq 0\right\\}}$ for $(x,t)\in\mathbb{R}^{2}$ is
known to be an isometry from $L^{2}(\mathbb{H}_{+})$ to
$L^{2}(\mathbb{H}_{+})$, see [2, 13, 19]. In particular, the operator
$\mathbf{A}$ is continuous.

######  Definition 2.8 (Forward operator for the wave equation).

We define $\mathbf{A}\colon L^{2}(\mathbb{H}_{+})\to L^{2}(\mathbb{H}_{+})$ by
$\mathbf{A}f\coloneqq t^{-1/2}u$, for $f\in C_{0}^{\infty}(\mathbb{H}_{+})$,
where $u$ is the solution of (2.8), and extending it by continuity to
$L^{2}(\mathbb{H}_{+})$.

The isometry property implies that any frame gives a BFD
$(\mathbf{U},\mathbf{V},\bm{\kappa})$ by setting
$v_{\lambda}=\mathbf{A}u_{\lambda}$ and $\kappa_{\lambda}=1$. This, in
particular, yields a wavelet vaguelette decomposition and a biorthogonal
curvelet decomposition for the wave equation.

######  Example 2.9.

BFD for the wave equation Let $(u_{\lambda})_{\lambda\in\Lambda}$ be either a
wavelet frame, a curvelet frame or a shearlet frame with
$\lambda=(j,k,\beta)\in\Lambda$ where $j\geq 0$ is the scale index. Then
$(\mathbf{U},\mathbf{V},\bm{\kappa})$ is a BFD with

$\displaystyle v_{\lambda}$ $\displaystyle\coloneqq\mathbf{A}u_{\lambda}$
(2.9) $\displaystyle\kappa_{\lambda}$ $\displaystyle\coloneqq 1\,.$ (2.10)

As noted in [13] this result directly follows from the isometry property and
the associated inversion formula $f=\mathbf{A}^{*}\mathbf{A}f$.

The isometry property also gives an efficient numerical algorithm for
computing analysis coefficients with respect to the frame $\mathbf{V}$ in the
case that $\mathbf{U}$ is associated with an efficient algorithm.

######  Algorithm 2.10.

Computing BFD coefficients for the wave equation Let $\mathbf{U}$ be the
curvelet frame and define $\mathbf{V}$ by (2.9).

  1. (a)

Input: $g\in L^{2}(\mathbb{H}_{+})$.

  2. (b)

Compute $f_{\rm FBP}\coloneqq\mathbf{A}^{*}g$.

  3. (c)

Compute $\bm{\eta}\coloneqq\mathbf{U}^{*}f_{\rm FBP}$ via wavelet, curvelet or
shearlet transform.

  4. (d)

Output: Coefficients $\bm{\eta}$.

The algorithm described above can be used for any problem where the forward
operator $\mathbf{A}$ is an isometry. In the case of the wavelet transform
this simple procedure has been previously used in [13].

##  3 Sparse $\ell^{1}$-regularization

There are two fundamentally different and well studied instances of sparse
frame based regularization, namely $\ell^{1}$-analysis regularization and
$\ell^{1}$-synthesis regularization. They are defined by

$\displaystyle\mathbf{B}_{\alpha}^{\rm ANA}(\bm{y}_{\delta})$
$\displaystyle\coloneqq\operatorname*{arg\,min}_{\bm{x}\in\mathbb{X}}{\left\\{\frac{1}{2}{\left\|\mathbf{A}\bm{x}-\bm{y}_{\delta}\right\|}^{2}+\alpha\sum_{\lambda}d_{\lambda}\left|\left\langle
u_{\lambda},\bm{x}\right\rangle\right|\right\\}}$ (3.1)
$\displaystyle\mathbf{B}_{\alpha}^{\rm SYN}(\bm{y}_{\delta})$
$\displaystyle\coloneqq\mathbf{W}\left(\operatorname*{arg\,min}_{\bm{\xi}\in\ell^{2}(\Lambda)}{\left\\{\frac{1}{2}{\left\|\mathbf{A}\mathbf{W}\bm{\xi}-\bm{y}_{\delta}\right\|}^{2}+\alpha\sum_{\lambda}d_{\lambda}\left|\xi_{\lambda}\right|\right\\}}\right)\,,$
(3.2)

respectively, with weights $d_{\lambda}>0$.

######  Definition 3.1.

We call $\bm{\xi}\in\ell^{2}(\Lambda)$ sparse if the set
${\left\\{\lambda\in\Lambda\mid\xi_{\lambda}\neq 0\right\\}}$ is finite.

If $\bm{\xi}=(\xi_{\lambda})_{\lambda\in\Lambda}\in\ell^{2}(\Lambda)$ is
sparse, we write
$\operatorname{Sign}(\bm{\xi})\coloneqq\\{\bm{z}=(z_{\lambda})_{\lambda\in\Lambda}\in\ell^{2}(\Lambda)\mid
z_{\lambda}\in\operatorname{Sign}(\xi_{\lambda})\\}$ where
$\operatorname{Sign}(\,\cdot\,)\colon\mathbb{R}\to\mathbb{R}$ is the multi-
valued signum function defined by $\operatorname{Sign}(0)=[-1,1]$ and
$\operatorname{Sign}(x)={\left\\{x/\left|x\right|\right\\}}$ for $x\neq 0$. We
will use the notation

$\left\|\;\cdot\;\right\|_{\bm{d},1}\colon\ell^{2}(\Lambda)\to\mathbb{R}\cup{\left\\{\infty\right\\}}\colon\bm{\xi}\mapsto\begin{cases}\sum_{\lambda}d_{\lambda}\left|\xi_{\lambda}\right|&\text{
if }(d_{\lambda}\xi_{\lambda})_{\lambda\in\Lambda}\in\ell^{1}(\Lambda)\\\
\infty&\text{ otherwise }\,.\end{cases}$ (3.3)

Any element in the set
$\operatorname*{arg\,min}\\{{\left\|\mathbf{U}^{*}(\bm{x})\right\|}_{\bm{d},1}\mid\mathbf{A}\bm{x}=\bm{y}\\}$
is called ${\left\|\mathbf{U}^{*}(\,\cdot\,)\right\|}_{\bm{d},1}$-minimizing
solution of $\mathbf{A}(\bm{x})=\bm{y}$. Note that
${\left\|\mathbf{U}^{*}(\,\cdot\,)\right\|}_{\bm{d},1}$-minimizing solutions
exists whenever there is any solution $\bm{x}$ with
${\left\|\mathbf{U}^{*}(\bm{x})\right\|}_{\bm{d},1}<\infty$, as follows from
[23, Theorem 3.25].

Below we recall well-posedness and convergence results for both variants.
These results hold under the following quite week assumptions:

  1. (A1)

$\mathbf{A}\colon\mathbb{X}\to\mathbb{Y}$ is bounded linear;

  2. (A2)

$\mathbf{U},\mathbf{W}$ are synthesis operators of frames
$(u_{\lambda})_{\lambda\in\Lambda}$, $(w_{\lambda})_{\lambda\in\Lambda}$ of
$\mathbb{X}$;

  3. (A3)

$\bm{d}=(d_{\lambda})_{\lambda\in\Lambda}\in\mathbb{R}^{\Lambda}$ satisfies
$\inf\\{d_{\lambda}\mid\lambda\in\Lambda\\}>0$.

For certain sparse elements we will state linear error estimates which have
been derived in [14, 16]. See [7, 17, 21, 3, 23] for some further works on
sparse $\ell^{1}$-regularization, and [8, 18, 22] for wavelet regularization
methods.

###  3.1 $\ell^{1}$-Analysis regularization

Let us define the $\ell^{1}$-analysis Tikhonov functional by

$\mathcal{A}_{\alpha,\bm{y}}\colon\mathbb{X}\to\mathbb{R}\cup{\left\\{\infty\right\\}}\colon\bm{x}\mapsto\frac{1}{2}{\left\|\mathbf{A}\bm{x}-\bm{y}\right\|}^{2}+\alpha{\left\|\mathbf{U}^{*}\bm{x}\right\|}_{\bm{d},1}\,.$
(3.4)

Then we have $\mathbf{B}_{\alpha}^{\rm
ANA}(\bm{y})=\operatorname*{arg\,min}\mathcal{A}_{\alpha,\bm{y}}$.

######  Proposition 3.2 (Convergence of analysis regularization).

Let (A1)–(A3) be satisfied, suppose $\bm{y}\in\mathbb{Y}$, $\alpha>0$,
$(\bm{y}^{k})_{k\in\mathbb{N}}\in\mathbb{Y}^{\mathbb{N}}$ with
$\bm{y}^{k}\to\bm{y}$, and choose
$\bm{x}^{k}\in\operatorname*{arg\,min}\mathcal{A}_{\alpha,\bm{y}^{k}}$.

  * $\blacksquare$

Existence: The functional $\mathcal{A}_{\alpha,\bm{y}}$ has at least one
minimizer.

  * $\blacksquare$

Stability: There exists a subsequence $(\bm{x}^{k(\ell)})_{\ell\in\mathbb{N}}$
of $(\bm{x}^{k})_{k\in\mathbb{N}}$ and a minimizer
$\bm{x}_{\alpha}\in\operatorname*{arg\,min}\mathcal{A}_{\alpha,\bm{y}}$ such
that $\|\bm{x}^{k(\ell)}-\bm{x}_{\alpha}\|\to 0$. If the minimizer
$\bm{x}_{\alpha}$ of $\mathcal{A}_{\alpha,\bm{y}}$ is unique, then
$\|\bm{x}^{k}-\bm{x}_{\alpha}\|\to 0$.

  * $\blacksquare$

Convergence: Assume $\bm{y}=\mathbf{A}\bm{x}$ for $\bm{x}\in\mathbb{X}$ with
$\|\mathbf{U}^{*}\bm{x}\|_{\bm{d},1}<\infty$ and suppose
$\|\bm{y}^{k}-\bm{y}\|\leq\delta_{k}$ with $(\delta_{k})_{k\in\mathbb{N}}\to
0$. Consider a parameter choice $(\alpha_{k})_{k}\in(0,\infty)^{\mathbb{N}}$
such that
$\lim_{k\to\infty}\alpha_{k}=\lim_{k\to\infty}\delta_{k}^{2}/\alpha_{k}=0$.
Then there is an
${\left\|\mathbf{U}^{*}(\,\cdot\,)\right\|}_{\bm{d},1}$-minimizing solution
$\bm{x}^{\bm{\texttt{+}}}$ of $\mathbf{A}(\bm{x})=\bm{y}$ and a subsequence
$(\bm{x}^{k(\ell)})_{\ell\in\mathbb{N}}$ with
$\|\bm{x}^{k(\ell)}-\bm{x}^{\bm{\texttt{+}}}\|\to 0$. If the
$\|\mathbf{U}^{*}(\,\cdot\,)\|_{\bm{d},1}$-minimizing solution is unique, then
$\|\bm{x}^{k}\to\bm{x}^{\bm{\texttt{+}}}\|\to 0$.

###### Proof.

See [14, Propositions 5, 6 and 7]. ∎

In order to derive convergence rates, one has to make additional assumptions
on the exact solution $\bm{x}^{\bm{\texttt{+}}}$ to be recovered. Besides the
sparsity this requires a certain interplay between $\bm{x}^{\bm{\texttt{+}}}$
and the forward operator $\mathbf{A}$.

  1. (A4)

$\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}$ is sparse;

  2. (A5)

$\exists\bm{z}\in\operatorname{Sign}(\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}})\colon\mathbf{U}\bm{z}\in\operatorname{ran}(\mathbf{A}^{*})$;

  3. (A6)

$\exists t\in(0,1)\colon\mathbf{A}\text{ is injective on
}\operatorname{span}{\left\\{u_{\lambda}\colon{\lvert\bm{z}_{\lambda}\rvert}>t\right\\}}$.

Assumption (A5) is a so-called source condition and the main restrictive
assumption. It requires that there exists an element
$\bm{z}\in\operatorname{Sign}(\mathbf{U}^{*}\bm{x})$ that satisfies the
smoothness assumption $\mathbf{U}\bm{z}\in\operatorname{ran}(\mathbf{A}^{*})$.
Because $\bm{z}\in\ell^{2}(\Lambda)$, the space
$\operatorname{span}{\left\\{u_{\lambda}\colon\left|z_{\lambda}\right|>t\right\\}}$
is finite dimensional. Therefore, Condition (A6) requires injectivity on a
certain finite dimensional subspace.

######  Proposition 3.3 (Convergence rates for analysis regularization).

Suppose (A1)–(A6) hold. Then, for a parameter choice $\alpha=C\delta$, there
is a constant $c_{\bm{\texttt{+}}}\in(0,\infty)$ such that for all
$\bm{y}_{\delta}\in\mathbb{Y}$ with
$\|\mathbf{A}\bm{x}-\bm{y}_{\delta}\|\leq\delta$ and every minimizer
$\bm{x}^{\alpha,\delta}\in\operatorname*{arg\,min}\mathcal{A}_{\alpha,\bm{y}_{\delta}}$
we have $\|\bm{x}^{\alpha,\delta}-\bm{x}^{\bm{\texttt{+}}}\|\leq
c_{\bm{\texttt{+}}}\delta$.

###### Proof.

See [16, Theorem III.8]. ∎

###  3.2 Synthesis regularization

Let use denote the $\ell^{1}$-synthesis Tikhonov functional by

$\mathcal{S}_{\alpha,\bm{y}}\colon\ell^{2}(\Lambda)\to\mathbb{R}\cup{\left\\{\infty\right\\}}\colon\bm{x}\mapsto\frac{1}{2}{\left\|\mathbf{A}\mathbf{W}\bm{\xi}-\bm{y}\right\|}^{2}+\alpha{\left\|\bm{\xi}\right\|}_{\bm{d},1}\,.$
(3.5)

Then it holds $\mathbf{B}_{\alpha}^{\rm
SYN}(\bm{y})=\mathbf{W}(\operatorname*{arg\,min}\mathcal{S}_{\alpha,\bm{y}})$.
Synthesis regularization can be seen as analysis regularization for the
coefficient inverse problem $\mathbf{A}\mathbf{W}\bm{\xi}=\bm{y}$ and the
analysis operator $\mathbf{U}^{*}=\operatorname{Id}$. Using Proposition 3.2 we
therefore have the following result.

######  Proposition 3.4 (Convergence of synthesis regularization).

Let (A1)-(A3) be satisfied, suppose $\bm{y}\in\mathbb{Y}$, $\alpha>0$,
$(\bm{y}^{k})_{k\in\mathbb{N}}\in\mathbb{Y}^{\mathbb{N}}$ with
$\bm{y}^{k}\to\bm{y}$ and take
$\bm{\xi}^{k}\in\operatorname*{arg\,min}\mathcal{S}_{\alpha,\bm{y}^{k}}$.

  * $\blacksquare$

Existence: The functional $\mathcal{S}_{\alpha,\bm{y}}$ has at least one
minimizer.

  * $\blacksquare$

Stability: There exists a subsequence
$(\bm{\xi}^{k(\ell)})_{\ell\in\mathbb{N}}$ of
$(\bm{\xi}^{k})_{k\in\mathbb{N}}$ and
$\bm{\xi}_{\alpha}\in\operatorname*{arg\,min}\mathcal{S}_{\alpha,\bm{y}}$ such
that $(\bm{\xi}^{k(\ell)})_{\ell}\to\bm{\xi}_{\alpha}$. If the minimizer of
$\mathcal{S}_{\alpha,\bm{y}}$ is unique, then
$\|\bm{\xi}^{k}-\bm{\xi}_{\alpha}\|\to 0$.

  * $\blacksquare$

Convergence: Assume $\bm{y}=\mathbf{A}\mathbf{W}\bm{\xi}$ for
$\bm{\xi}\in\ell^{2}(\Lambda)$ with $\|\bm{\xi}\|_{\bm{d},1}<\infty$ and
$\|\bm{y}^{k}-\bm{y}\|\leq\delta_{k}$ with $(\delta_{k})_{k\in\mathbb{N}}\to
0$. Consider a parameter choice
$(\alpha_{k})_{k\in\mathbb{N}}\in(0,\infty)^{\mathbb{N}}$ with
$\lim_{k\to\infty}\alpha_{k}=\lim_{k\to\infty}\delta_{k}^{2}/\alpha_{k}=0$.
Then there exist an ${\left\|\,\cdot\,\right\|}_{\bm{d},1}$-minimizing
solution $\bm{\xi}^{\bm{\texttt{+}}}$ of
$(\mathbf{A}\mathbf{W})(\bm{\xi})=\bm{y}$ and a subsequence
$(\bm{\xi}^{k(\ell)})_{\ell\in\mathbb{N}}$ with
$\|\bm{\xi}^{k(\ell)}-\bm{\xi}^{\bm{\texttt{+}}}\|\to 0$. If the
$\|\,\cdot\,\|_{\bm{d},1}$-minimizing solution is unique, then
$\|\bm{\xi}^{k}-\bm{\xi}^{\bm{\texttt{+}}}\|\to 0$.

###### Proof.

Follows from Proposition 3.2 with $\mathbf{U}=\operatorname{Id}$ and
$\mathbf{A}\mathbf{W}$ in place of $\mathbf{A}$. ∎

We have linear convergence rates under the following additional assumptions on
the element to be recovered.

  1. (S4)

$\bm{x}^{\bm{\texttt{+}}}=\mathbf{W}\bm{\xi}^{\bm{\texttt{+}}}$ where
$\bm{\xi}^{\bm{\texttt{+}}}\in\ell^{2}(\Lambda)$ is sparse;

  2. (S5)

$\exists\bm{z}\in\operatorname{Sign}(\bm{\xi}^{\bm{\texttt{+}}})\colon\bm{z}=\operatorname{ran}(\mathbf{W}^{*}\mathbf{A}^{*})$;

  3. (S6)

$\exists t\in(0,1)\colon\mathbf{A}\mathbf{W}\text{ is injective on
}\operatorname{span}{\left\\{e_{\lambda}\colon\left|\bm{z}_{\alpha}\right|>t\right\\}}$.

######  Proposition 3.5 (Convergence rates for synthesis regularization).

Suppose that (A1)–(A3) and (S4)–(S6) hold. Then, for a parameter choice
$\alpha=C\delta$, there is a constant $c_{\bm{\texttt{+}}}\in(0,\infty)$ such
that for all $\bm{y}_{\delta}\in\mathbb{Y}$ with
$\|\mathbf{A}\bm{x}-\bm{y}_{\delta}\|\leq\delta$, every minimizer
$\bm{\xi}^{\alpha,\delta}\in\operatorname*{arg\,min}\mathcal{S}_{\alpha,\bm{y}_{\delta}}$
we have $\|\bm{\xi}^{\alpha,\delta}-\bm{\xi}^{\bm{\texttt{+}}}\|\leq
c_{\bm{\texttt{+}}}\delta$.

###### Proof.

Follows from Proposition 3.3. ∎

Because $\mathbf{W}$ is bounded, the above convergence results can be
transferred to convergence in the signal space $\mathbb{X}$. In particular, we
have stability
$\|\mathbf{W}(\bm{\xi}^{k(\ell)})-\mathbf{W}(\bm{\xi}_{\alpha})\|\to 0$ and
convergence
$\|\mathbf{W}\bm{\xi}^{k(\ell)}-\mathbf{W}\bm{\xi}^{\bm{\texttt{+}}}\|\to 0$
under the assumptions made in Proposition 3.4, and linear convergence rates
$\|\mathbf{W}\bm{\xi}^{\alpha,\delta}-\bm{x}^{\bm{\texttt{+}}}\|\leq\tilde{c}_{\bm{\texttt{+}}}\delta$
under the assumptions made in Proposition 3.5.

###  3.3 Sparse regularization using an SVD

In the special case that $\mathbf{U}$ is part of an SVD, then analysis and
synthesis regularization are equivalent and can be computed explicitly by
soft-thresholding of the expansion coefficients.

######  Theorem 3.6 (Equivalence in the SVD case).

Let $(\mathbf{U},\mathbf{V},\bm{\kappa})$ be an SVD for $\mathbf{A}$, let
$\bm{y}_{\delta}\in\mathbb{Y}$ and consider (3.1), (3.2) with
$\mathbf{W}=\mathbf{U}$. Then

$\mathbf{B}_{\alpha}^{\rm ANA}(\bm{y}_{\delta})=\mathbf{B}_{\alpha}^{\rm
SYN}(\bm{y}_{\delta})=\\{(\mathbf{U}\circ\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}}\circ\mathbf{V}^{\ast})(\bm{y}_{\delta})\\}\,,$

equals the soft-thresholding estimator in the SVD system.

###### Proof.

Because $\mathbf{U}$ is an orthonormal basis of $\mathbb{X}$, we have
$\bm{x}=\mathbf{U}\bm{\xi}\Leftrightarrow\xi=\mathbf{U}^{*}\bm{x}$ which
implies that $\mathbf{B}_{\alpha}^{\rm
ANA}(\bm{y}_{\delta})=\mathbf{B}_{\alpha}^{\rm SYN}(\bm{y}_{\delta})$. Now let
$\bm{x}_{\alpha}\in\mathbf{B}_{\alpha}^{\rm SYN}(\bm{y}_{\delta})$ be any
minimizer of the $\ell^{1}$-analysis Tikhonov functional
$\mathcal{A}_{\alpha,\bm{y}}$. Let
$\mathbf{P}_{\operatorname{ran}(\mathbf{A})^{\bot}}$ denote the orthogonal
projection on $\operatorname{ran}(\mathbf{A})^{\bot}$. We have

$\displaystyle\mathcal{A}_{\alpha,\bm{y}}(\bm{x})$
$\displaystyle=\frac{1}{2}{\left\|\mathbf{A}\bm{x}-\bm{y}\right\|}^{2}+\alpha{\left\|\mathbf{U}^{*}\bm{x}\right\|}_{\bm{d},1}$
$\displaystyle={\left\|\mathbf{P}_{\operatorname{ran}(\mathbf{A})^{\bot}}(\bm{y})\right\|}^{2}+\sum_{\lambda\in\Lambda}\frac{1}{2}\left|\left\langle\mathbf{A}\bm{x}-\bm{y},v_{\lambda}\right\rangle\right|^{2}+\sum_{\lambda\in\Lambda}\alpha
d_{\lambda}\left|\left\langle\bm{x},u_{\lambda}\right\rangle\right|$
$\displaystyle={\left\|\mathbf{P}_{\operatorname{ran}(\mathbf{A})^{\bot}}(\bm{y})\right\|}^{2}+\sum_{\lambda\in\Lambda}\frac{1}{2}\left|\kappa_{\lambda}\left\langle\bm{x},u_{\lambda}\right\rangle-\left\langle\bm{y},v_{\lambda}\right\rangle\right|^{2}+\alpha
d_{\lambda}\left|\left\langle\bm{x},u_{\lambda}\right\rangle\right|\,.$

The latter sum is minimized by componentwise soft-thresholding. This shows
$\bm{x}_{\alpha,\delta}=(\mathbf{U}\circ\mathbb{S}_{\alpha\bm{d}/\bm{\kappa}}\circ\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbf{V}^{\ast})(\bm{y}_{\delta})=(\mathbf{U}\circ\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}}\circ\mathbf{V}^{\ast})(\bm{y}_{\delta})$
and concludes the proof. ∎

In the case that $(\mathbf{U},\mathbf{V},\bm{\kappa})$ is a redundant BFD
expansion and not an SVD, then (3.1), (3.2), and the soft-thresholding
estimator

$(\mathbf{U}\circ\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}}\circ\mathbf{V}^{\ast})(\bm{y}_{\delta})$
(3.6)

are all non-equivalent. Further, in this case, (3.1) and (3.2) have to be
computed by iterative minimization algorithms. This requires repeated
application of the forward and adjoint problem and therefore is time
consuming. In the following section, we study BFD thresholding which is the
analog of (3.6) for redundant systems. Despite the non-equivalence to
$\ell^{1}$-regularization, we are able to derive the same type of convergence
results and linear convergence rates as for the analysis and synthesis
variants of $\ell^{1}$-regularization.

##  4 Regularization via BFD thresholding

Throughout this section we fix the following assumptions:

  1. (B1)

$\mathbf{A}\colon\mathbb{X}\to\mathbb{Y}$ is bounded linear.

  2. (B2)

$(\mathbf{U},\mathbf{V},\bm{\kappa})$ is s BFD for $\mathbf{A}$.

  3. (B3)

$\bm{d}=(d_{\lambda})_{\lambda\in\Lambda}\in\mathbb{R}^{\Lambda}$ satisfies
$\inf\\{d_{\lambda}\mid\lambda\in\Lambda\\}>0$.

In this section we show well-posedness, convergence and convergence rates for
BFD soft-thresholding.

###  4.1 BFD soft-thresholding

Any BFD gives an explicit inversion formula
$\bm{x}=(\bar{\mathbf{U}}\circ\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbf{V}^{*})(\mathbf{A}\bm{x})$
where $\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}$ is defined by (2.2). For
ill-posed problems, $\kappa_{\lambda}\to 0$ and therefore the above
reproducing formula is unstable when applied to noisy data $\bm{y}_{\delta}$
instead of $\mathbf{A}\bm{x}$. Below we stabilize the inversion by including
the soft-thresholding operation.

######  Definition 4.1 (BFD soft-thresholding).

Let $(\mathbf{U},\mathbf{V},\bm{\kappa})$ be a BFD for $\mathbf{A}$. We define
the nonlinear BFD soft-thresholding estimator by

$\mathbf{B}_{\alpha}^{\rm
BFD}\colon\mathbb{Y}\to\mathbb{X}\colon\bm{y}\mapsto(\bar{\mathbf{U}}\circ\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}}\circ\mathbf{V}^{*})(\bm{y})\,.$
(4.1)

If $(u_{\lambda})_{\lambda}$, $(v_{\lambda})_{\lambda}$ are ONBs, then Theorem
3.6 shows that (3.1) and (3.2) are equivalent to (4.1). In the case of general
frames, $\mathbf{B}_{\alpha}^{\rm ANA}$, $\mathbf{B}_{\alpha}^{\rm SYN}$ and
$\mathbf{B}_{\alpha}^{\rm BFD}$ are all different.

As the main result in this paper we show that BFD soft-thresholding yields the
same theoretical results as $\ell^{1}$-regularization. Assuming efficient
implementations for $\bar{\mathbf{U}}$ and $\mathbf{V}^{*}$, the BFD estimator
has the advantage that it can be calculated non-iteratively and is therefore
much faster than $\mathbf{B}_{\alpha}^{\rm SYN}$ and $\mathbf{B}_{\alpha}^{\rm
BFD}$.

Consider the $\ell^{1}$-Tikhonov functional for the multiplication operator
$\mathbf{M}_{\bm{\kappa}}$,

$\mathcal{M}_{\alpha,\bm{\eta}}\colon\ell^{2}(\Lambda)\to\mathbb{R}\cup{\left\\{\infty\right\\}}\colon\bm{\xi}\mapsto\frac{1}{2}{\left\|\mathbf{M}_{\bm{\kappa}}\bm{\xi}-\bm{\eta}\right\|}^{2}+\alpha{\left\|\bm{\xi}\right\|}_{\bm{d},1}\,.$
(4.2)

The proof strategy used in this paper is based on the following Lemma.

######  Lemma 4.2 ($\ell^{1}$-minimization for multiplication operators).

  1. (a)

$\forall\alpha\in\mathbb{R}_{>0}\forall\bm{\eta}\in\ell^{2}(\Lambda)\colon$
$\mathcal{M}_{\alpha,\bm{\eta}}$ has a unique minimizer.

  2. (b)

$\forall\alpha\in\mathbb{R}_{>0}\forall\bm{\eta}\in\ell^{2}(\Lambda)\colon$
$(\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}})(\bm{\eta})=\operatorname*{arg\,min}\mathcal{M}_{\alpha,\bm{\eta}}$.

  3. (c)

$\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}}\colon\ell^{2}(\Lambda)\to\ell^{2}(\Lambda)$
is continuous.

  4. (d)

$\mathbf{V}^{*}\mathbf{A}=\mathbf{M}_{\bm{\kappa}}\mathbf{U}^{*}$.

###### Proof.

Because $(\operatorname{Id},\operatorname{Id},\bm{\kappa})$ is an SVD for
$\mathbf{M}_{\bm{\kappa}}$, Items (a), (b) follow from Theorem 3.6, the
equivalence of $\ell^{1}$-regularization and soft-thresholding in the SVD
case. Item (c) follows from Proposition3.4. Moreover, the equality
$(\mathbf{V}^{*}\mathbf{A}\bm{x})_{\lambda}=\left\langle
v_{\lambda},\mathbf{A}\bm{x}\right\rangle=\left\langle\mathbf{A}^{*}v_{\lambda},\bm{x}\right\rangle=\kappa_{\lambda}\left\langle
u_{\lambda},\bm{x}\right\rangle=(\mathbf{M}_{\bm{\kappa}}\mathbf{U}^{*}\bm{x})_{\lambda}$
shows Item (d). ∎

Note that the continuity of
$\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}}$
(see Item (b) in the above lemma) is not obvious as it is the composition of
the soft thresholding $\mathbb{S}_{\alpha\bm{d}}$ with the discontinuous
operator in $\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}$. The characterization
in Item (b) of
$\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}}$ as
minimizer of the $\ell^{1}$-Tikhonov functional
$\mathcal{M}_{\alpha,\bm{\eta}}$ and the existing stability results for
$\ell^{1}$-Tikhonov regularization yields an elegant way to obtain the
continuity of
$\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}}$.
Verifying the continuity directly would also be possible but seems to be a
harder task. A similar comment applies to the proof of Theorem 4.3 where we
use the convergence of $\mathcal{M}_{\alpha,\bm{\eta}}$ to show convergence
the BFD soft-thresholding estimator $\mathbf{B}_{\alpha}^{\rm BFD}$.

###  4.2 Convergence analysis

In this section, we show that $(\mathbf{B}_{\alpha}^{\rm BFD})_{\alpha>0}$ is
well-posed and convergent.

######  Theorem 4.3 (Well-posedness and convergence).

Let (B1)-(B3) be satisfied, suppose $\bm{y}\in\mathbb{Y}$ and let
$(\bm{y}^{k})_{k\in\mathbb{N}}\in\mathbb{Y}^{\mathbb{N}}$ satisfy
$\bm{y}^{k}\to\bm{y}$.

  1. (a)

Existence: $\mathbf{B}_{\alpha}^{\rm BFD}\colon\mathbb{Y}\to\mathbb{X}$ is
well-defined for all $\alpha>0$.

  2. (b)

Stability: $\mathbf{B}_{\alpha}^{\rm BFD}\colon\mathbb{Y}\to\mathbb{X}$ is
continuous for all $\alpha>0$.

  3. (c)

Convergence: Assume $\bm{y}=\mathbf{A}\bm{x}$ for some $\bm{x}\in\mathbb{X}$
with $\|\mathbf{U}^{*}\bm{x}\|_{\bm{d},1}<\infty$, suppose
$\|\bm{y}^{k}-\bm{y}\|\leq\delta_{k}$ with $(\delta_{k})_{k\in\mathbb{N}}\to
0$ and consider a parameter choice
$(\alpha_{k})_{k\in\mathbb{N}}\in(0,\infty)^{\mathbb{N}}$ with
$\lim_{k\to\infty}\alpha_{k}=\lim_{k\to\infty}\delta_{k}^{2}/\alpha_{k}=0$.
Then $\|\mathbf{B}_{\alpha_{k}}^{\rm
BFD}(\bm{x}^{k})-\bm{x}^{\bm{\texttt{+}}}\|\to 0$.

###### Proof.

(a), (b): According to Lemma 4.2, the mapping
$\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}}\colon\ell^{2}(\Lambda)\to\ell^{2}(\Lambda)$
is well defined and continuous. Moreover, by definition we have
$\mathbf{B}_{\alpha}^{\rm
BFD}=\bar{\mathbf{U}}\circ(\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}})\circ\mathbf{V}^{*}$
which implies existence and stability of BFD thresholding.
(c): We have

$\|\mathbf{B}_{\alpha_{k}}^{\rm
BFD}(\bm{y}^{k})-\bm{x}^{\bm{\texttt{+}}}\|=\|\bar{\mathbf{U}}\circ(\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}})\circ\mathbf{V}^{*}(\bm{y}^{k})-\bar{\mathbf{U}}\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}\|\\\
\leq{\left\|\mathbf{U}\right\|}\;\|(\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}})(\mathbf{V}^{*}\bm{y}^{k})-\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}\|\,.$
(4.3)

Moreover,
$\|\mathbf{V}^{*}\bm{y}^{k}-\mathbf{M}_{\bm{\kappa}}\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}\|=\|\mathbf{V}^{*}\bm{y}^{k}-\mathbf{V}^{*}\mathbf{A}\bm{x}^{\bm{\texttt{+}}}\|\leq{\left\|\mathbf{V}\right\|}\|\bm{y}^{k}-\mathbf{A}\bm{x}^{\bm{\texttt{+}}}\|\leq{\left\|\mathbf{V}\right\|}\delta_{k}$.
Therefore, Proposition 3.4 and the equality
$\operatorname*{arg\,min}\mathcal{M}_{\alpha,\mathbf{V}^{*}\bm{y}^{k}}=(\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}})(\mathbf{V}^{*}\bm{y}^{k})$
shown in Lemma 4.2 imply
$\|(\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}})(\mathbf{V}^{*}\bm{y}^{k})-\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}\|\to
0$ for $k\to\infty$. Together with (4.3) this yields (c) and completes the
proof. ∎

###  4.3 Convergence rates

Next we derive linear convergence rates for sparse solutions. Let use denote
by
$\operatorname{supp}(\xi)\coloneqq{\left\\{\lambda\in\Lambda\colon\xi_{\lambda}\neq
0\right\\}}$ the support of $\xi\in\ell^{2}(\Lambda)$. To derive the
convergence rates, we assume the following for the exact solution
$\bm{x}^{\bm{\texttt{+}}}$ to be recovered.

  1. (B4)

$\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}$ is sparse.

  2. (B5)

$\forall\lambda\in\operatorname{supp}(\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}})\colon\kappa_{\lambda}\neq
0$.

Note that assumptions (B4), (B5) imply the source condition

$\bm{z}\in\partial\|\,\cdot\,\|_{\bm{d},1}\cap\operatorname{ran}(\mathbf{M}_{\bm{\kappa}}^{*})\neq\emptyset\,.$

is satisfied for some element $\bm{z}\in\ell^{2}(\Lambda)$ that can be chosen
such that $\left|z_{\lambda}\right|<1$ for
$\lambda\not\in\operatorname{supp}(\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}})$.
Moreover, it follows that $\mathbf{M}_{\bm{\kappa}}$ is injective on
$\operatorname{span}{\left\\{e_{\lambda}\mid\left|\bm{z}_{\lambda}\right|>t\right\\}}$
with
$t\coloneqq\max{\left\\{\left|\kappa_{\lambda}\right|\mid\lambda\not\in\operatorname{supp}(\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}})\right\\}}$.
Because $\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}\in\ell^{2}(\Lambda)$, we have
$t<1$.

######  Theorem 4.4 (Convergence rates).

Suppose that (B1)–(B5) hold. Then, for the parameter choice $\alpha=C\delta$
with $C\in(0,\infty)$, there is a constant $c_{\bm{\texttt{+}}}\in(0,\infty)$
such that for all $\bm{y}_{\delta}\in\mathbb{Y}$ with
$\|\mathbf{A}\bm{x}-\bm{y}_{\delta}\|\leq\delta$ we have
$\|\mathbf{B}_{\alpha}^{\rm
BFD}(\bm{y}_{\delta})-\bm{x}^{\bm{\texttt{+}}}\|\leq
c_{\bm{\texttt{+}}}\delta$.

###### Proof.

As in the proof of Theorem 4.3 one obtains

$\displaystyle\|\mathbf{B}_{\alpha}^{\rm
BFD}(\bm{y}^{k})-\bm{x}^{\bm{\texttt{+}}}\|$
$\displaystyle\leq{\left\|\mathbf{U}\right\|}\;\|(\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}})(\mathbf{V}^{*}\bm{y}^{k})-\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}\|$
(4.4)
$\displaystyle\|\mathbf{V}^{*}\bm{y}^{k}-\mathbf{M}_{\lambda}\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}\|$
$\displaystyle\leq{\left\|\mathbf{V}\right\|}\delta\,.$ (4.5)

According to the considerations below (B4), (B5) the conditions (S4)–(S6) are
satisfied for the operator $\mathbf{M}_{\bm{\kappa}}$ in place of $\mathbf{A}$
and with $\mathbf{W}=\operatorname{Id}$. The convergence rates result in
Proposition 3.3, estimate (4.5), and the identity
$\operatorname*{arg\,min}\mathcal{M}_{\alpha,\mathbf{V}^{*}\bm{y}}=(\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}})(\mathbf{V}^{*}\bm{y})$
shown in Lemma 4.2 imply
$\|(\mathbf{M}_{\bm{\kappa}}^{\bm{\texttt{+}}}\circ\mathbb{S}_{\alpha\bm{d}})(\mathbf{V}^{*}\bm{y})-\mathbf{U}^{*}\bm{x}^{\bm{\texttt{+}}}\|\leq
c\delta$. Together with (4.4) this implies $\|\mathbf{B}_{\alpha}^{\rm
BFD}(\bm{y}_{\delta})-\bm{x}^{\bm{\texttt{+}}}\|\leq
c{\left\|\mathbf{U}\right\|}\delta$ and concludes the proof. ∎

##  5 Conclusion

To overcome the inherent ill-posedness of inverse problems, regularization
methods incorporate available prior information about the unknowns to be
reconstructed. In this context, a useful prior is sparsity with respect to a
certain frame. There are at least two different regularization strategies
implementing sparsity with respect to a frame, namely $\ell^{1}$-analysis
regularization and $\ell^{1}$-synthesis regularization. In this paper, we
analyzed BFD-thresholding as a third variant of sparse regularization. One
advantage of BFD-thresholding compared to other sparse regularization methods
is its non-iterative nature leading to fast algorithms. Besides having a BFD,
actually computing the BFD soft-thresholding estimator (4.1) requires the dual
frame $\bar{\mathbf{U}}$. While in the general situation, the BFD and the dual
frame have to be computed numerically, we have shown that for many practical
examples (see Section 2) they are known explicitly and efficient algorithms
are available for its numerical evaluation.

The BFD-approach presented in this paper is well studied in the context of
statistical estimating using certain multi-scale systems. However, its
analysis in the context of regularization theory has not been given so far. In
this paper we closed this gap and presented a complete convergence analysis of
BFD-thresholding as regularization method.

## Acknowledgement

The work of M.H has been supported by the Austrian Science Fund (FWF), project
P 30747-N32.

## References

  * [1]
  * [2] A. L. Buhgeim and V. B. Kardakov.  Solution of an inverse problem for an elastic wave equation by the method of spherical means.  Sibirsk. Mat. Z., 19(4):749–758, 1978.
  * [3] S. Bürger, J. Flemming, and B. Hofmann.  On complex-valued deautoconvolution of compactly supported functions with sparse fourier representation.  Inverse Probl., 32(10):104006, 2016.
  * [4] E. J. Candès and D. Donoho.  Recovering edges in ill-posed inverse problems: Optimality of curvelet frames.  Ann. Statist., 30(3):784–842, 2002.
  * [5] O. Christensen.  An introduction to frames and Riesz bases.  Applied and Numerical Harmonic Analysis. Birkhäuser Boston Inc., Boston, MA, 2003.
  * [6] F. Colonna, G. Easley, K. Guo, and D. Labate.  Radon transform inversion using the shearlet representation.  Appl. Comput. Harmon. Anal., 29(2):232–250, 2010.
  * [7] I. Daubechies, M. Defrise, and C. De Mol.  An iterative thresholding algorithm for linear inverse problems with a sparsity constraint.  Comm. Pure Appl. Math., 57(11):1413–1457, 2004.
  * [8] V. Dicken and P. Maass.  Wavelet-Galerkin methods for ill-posed problems.  J. Inverse Ill-Posed Probl., 4:203–221, 1996.
  * [9] D. L. Donoho.  Nonlinear solution of linear inverse problems by wavelet-vaguelette decomposition.  Appl. Comput. Harmon. Anal., 2(2):101–126, 1995.
  * [10] M. Elad, P. Milanfar, and R. Rubinstein.  Analysis versus synthesis in signal priors.  Inveerse Problems, 23(3):947, 2007.
  * [11] H. W. Engl, M. Hanke, and A. Neubauer.  Regularization of Inverse Problems.  Kluwer Academic Publishers, Dordrecht, 1996.
  * [12] J. Frikel.  Sparse regularization in limited angle tomography.  Appl. Comput. Harmon. Anal., 34(1):117–141, 2013.
  * [13] J. Frikel and M. Haltmeier.  Efficient regularization with wavelet sparsity constraints in photoacoustic tomography.  Inverse Problems, 34(2):024006, 2018.
  * [14] M. Grasmair, M. Haltmeier, and O. Scherzer.  Sparse regularization with $l^{q}$ penalty term.  Inverse Problems, 24(5):055020, 13, 2008.
  * [15] M. Grasmair, M. Haltmeier, and O. Scherzer.  Necessary and sufficient conditions for linear convergence of $\ell^{1}$-regularization.  Comm. Pure Appl. Math., 64(2):161–182, 2011.
  * [16] M. Haltmeier.  Stable signal reconstruction via $\ell^{1}$-Minimization in redundant, non-tight frames.  IEEE Trans. Sig. Processing, 61(2):420–426, 2013.
  * [17] D. A. Lorenz.  Convergence rates and source conditions for tikhonov regularization with sparsity constraints.  Journal of Inverse and Ill-Posed Problems, 16(5):463–478, 2008\.
  * [18] A.K. Louis, P. Maass, and A. Rieder.  Wavelets. Theorie und Anwendungen.  Teubner, Stuttgart, 1998.
  * [19] E. K. Narayanan and Rakesh.  Spherical means with centers on a hyperplane in even dimensions.  Inverse Probl., 26(3):035014, March 2010.
  * [20] F. Natterer.  The Mathematics of Computerized Tomography, volume 32 of Classics in Applied Mathematics.  SIAM, Philadelphia, 2001.
  * [21] R. Ramlau and G. Teschke.  A tikhonov-based projection iteration for nonlinear ill-posed problems with sparsity constraints.  Numerische Mathematik, 104(2):177–203, 2006.
  * [22] A. Rieder.  A wavelet multilevel method for ill-posed problems stabilized by Tikhonov regularization.  Numer. Math., 75:501–522, 1997.
  * [23] O. Scherzer, M. Grasmair, H. Grossauer, M. Haltmeier, and F. Lenzen.  Variational methods in imaging.  Applied Mathematical Sciences, 167, 2009.