# Optimal-order convergence of Nesterov acceleration for linear ill-posed
problems

Stefan Kindermann Industrial Mathematics Institute, Johannes Kepler University
Linz, Linz, Austria (kindermann@indmath.uni-linz.ac.at).

###### Abstract

We show that Nesterov acceleration is an optimal-order iterative
regularization method for linear ill-posed problems provided that a parameter
is chosen accordingly to the smoothness of the solution. This result is proven
both for an a priori stopping rule and for the discrepancy principle. The
essential tool to obtain this result is a representation of the residual
polynomials via Gegenbauer polynomials.

> This article is dedicated to A. Neubauer on the occasion of his 60th
> birthday. His analytical insight shall be our unmatched benchmark.

##  1 Introduction

One option to calculate a regularized solution to an ill-posed operator
equation $Ax=y^{\delta}$, with $A:X\to Y$ and $X,Y$ being Hilbert spaces, is
to employ iterative regularization schemes, where approximate solutions
$x_{k}^{\delta}$ are calculated iteratively combined with a stopping rule as
regularization parameter choice. The simplest one being Landweber iteration
(cf., e.g., [4]), which has the downside of being rather slow. To speed up
convergence, acceleration schemes may be used such as the following Nesterov
acceleration:

$\displaystyle\begin{split}x_{k+1}^{\delta}&=z_{k}^{\delta}+A^{*}(y^{\delta}-Az_{k}^{\delta}),\qquad
k\geq 1\\\
z_{k}^{\delta}&=x_{k}^{\delta}+\alpha_{k}(x_{k}^{\delta}-x_{k-1}^{\delta}),\qquad
x_{0}=0,x_{1}=A^{*}y^{\delta},\end{split}$ (1)

where $\|A^{*}A\|\leq 1$ is assumed and where the sequence $\alpha_{k}$ is
chosen, for instance, as

$\displaystyle\alpha_{k}=\frac{k-1}{k+\beta},\qquad k\geq 1,\quad\beta>-1.$
(2)

Here $\beta$ is a parameter; common choices are, for example, $\beta=1$ or
$\beta=2$. We remark that other alternatives for the sequence $\alpha_{k}$ are
possible as well, but for the main analysis of this paper we only consider
(2).

This iteration (in a general nonlinear context) was suggested by Yurii
Nesterov for general convex optimization problems [9]. It is an instance of a
method that achieves the best rate of convergence (in the sense of objective
function decrease) that is generally possible for a first-order method. The
Nesterov acceleration can be employed to speed up convergence of gradient
methods in nonlinear or convex optimization. A particular successful instance
is the FISTA algorithm of Beck and Teboulle [2] for nondifferentiable convex
optimization.

In the realm of ill-posed problems, Hubmer and Ramlau [8] performed a
convergence analysis for nonlinear problems, and shows the efficiency of the
method.

The background and main motivation of the present article is the recent
interesting analysis of Neubauer [10] for ill-posed problems in the linear
case. He showed that (1) is an iterative regularization scheme and, more
important proved convergence rates, which are of optimal order only for a
priori parameter choices and in case of low smoothness of the solution and
suboptimal else. What is puzzling is that the method shows a quite unusual
”semi-saturation” phenomenon (we explain this term below in Section 3.1).

Our contribution in this article is twofold: At first, we prove a formula for
the residuals of the iteration (1) involving Gegenbauer polynomials. On this
basis, we can build a convergence rate analysis, which improves and extends
the results of Neubauer. In particular, we show that the method can always be
made an optimal-order method if the parameter $\beta$ is chosen accordingly to
the (Hölder-)smoothness index of the solution. This result holds for both an a
priori stopping rule and for the discrepancy principle.

Our analysis also explains the quite nebulous role that this parameter plays
in the iteration, as it is related to the index of the orthogonal polynomials
appearing in the residual formula.

Moreover, the above mentioned residual representation also clearly elucidates
the semi-saturation phenomenon because the iteration can be interpreted as a
mixture of a saturating iteration (Brakhage’s $\nu$-method) and a non-
saturating one (Landweber method).

In the following we employ some standard notation of regularization theory as
in [4]: $\delta=\|Ax^{\dagger}-y^{\delta}\|$ is the noise level and
$x^{\dagger}$ denotes the minimum-norm solution to the operator equation
$Ax=y$ with exact data $y=Ax^{\dagger}$. The index $\delta$ of $y^{\delta}$
indicates noisy data, and analogous, $x_{k}^{\delta}$ denotes the iterates of
(1) with noisy data $y^{\delta}$, while the lack of $\delta$ indicates exact
data $y$ and correspondingly the iteration $x_{k}$ with exact data $y$ in
place of $y^{\delta}$ in (1).

##  2 Residual polynomials for Nesterov acceleration

Our work follows the general theory of spectral filter-based regularization
methods as in [4], where the convergence analysis results from estimates of
the corresponding filter function. The first main result, Theorem 1 is quite
useful for this purpose as it represents the residual function in terms of
known polynomials.

The iteration (1) is a Krylov-space method, and the residual can be expressed
as

$y^{\delta}-Ax_{k}^{\delta}=:r_{k}(AA^{*})y^{\delta}$

with the residual polynomials satisfying the recurrence relation (cf. [10])

$\displaystyle\begin{split}r_{k+1}(\lambda)&=(1-\lambda)\left[r_{k}(\lambda)+\alpha_{k}(r_{k}(\lambda)-r_{k-1}(\lambda))\right],\quad
k\geq 1,\\\ \qquad r_{0}(\lambda)&=1,\qquad
r_{1}(\lambda)=(1-\lambda).\end{split}$ (3)

This is a simple consequence of the definition in (1). The $k$-th iterate can
be expressed via spectral filter functions

$x_{k}^{\delta}=g_{k}(A^{*}A)A^{*}y^{\delta},\qquad
g_{k}(\lambda):=\frac{1-r_{k}(\lambda)}{\lambda}.$

Observe that the three-term recursion (3) is not of the form to apply Favard’s
theorem [5], hence $r_{k}$ is not any orthogonal polynomial with respect to
some weight functions. (Note that Favard’s theorem fully characterizes three-
term recurrence relations that lead to orthogonal polynomials).

Before we proceed, we may compare the residual polynomials with other well-
known cases. For classical Landweber iteration [4], which is obtained by
setting $\alpha_{k}=0$ and thus $z_{k}^{\delta}=x_{k}^{\delta}$, the
corresponding residual functions $r_{k}=:r_{k}^{(LW)}$ is

$r_{k}^{(LW)}(\lambda)=(1-\lambda)^{k}.$

On the other hand, another class of well-known iteration methods for ill-posed
problems that are based on orthogonal polynomials are two-step semiiterative
methods [6]. They have the form

$\displaystyle x_{k+1}^{\delta}$
$\displaystyle=x_{k}^{\delta}+\mu_{k+1}(x_{k}-x_{k-1})+\omega_{k+1}A^{*}(y^{\delta}-Ax_{k}),\qquad
k>1,$

where $\mu_{k}$ and $\omega_{k}$ are appropriately chosen sequences. The
corresponding residual functions satisfy the recurrence relation

$\displaystyle r_{k+1}(\lambda)$
$\displaystyle=(1-\omega_{k+1}\lambda)r_{k}(\lambda)+\mu_{k+1}(r_{k}-r_{k-1}),\qquad
k>1,$ (4)

and thus, $r_{k}(\lambda)$ form a sequence of orthogonal polynomials. Of
special interest in ill-posed problems are the $\nu$-methods of Brakhage [3,
6], defined by the sequences, for $k>1$,

$\displaystyle\mu_{k+1}$
$\displaystyle=\frac{(k-1)(2k-2)(2k+2\nu-1)}{(k+2\nu-1)(2k+4\nu-1)(2k+2\nu-3)},$
$\displaystyle\qquad\omega_{k+1}$
$\displaystyle=4\frac{(2k+2\nu-1)(k+\nu-1)}{(k+2\nu-1)(2k+4\nu-1)},$

the initial values $x_{0}=0,$ $x_{1}=\frac{4\nu+2}{4\nu+1}T^{*}y^{\delta}$,
and with $\nu>0$ a user-selected parameter. The associated residual
polynomials $r_{k}=:r_{k}^{(\nu)}$ related to (4) with $r_{0}=1$,
$r_{1}=1-\lambda\frac{4\nu+2}{4\nu+1},$ have the representation [3]

$r_{k}^{(\nu)}(\lambda)=\frac{C_{2k}^{(2\nu)}(\sqrt{1-\lambda})}{C_{2k}^{(2\nu)}(1)},$

where $C_{n}^{(\alpha)}$ denotes the Gegenbauer polynomials (aka.
ultraspherical polynomials); cf. [1].

We now obtain the corresponding representation for the Nesterov residual
polynomials, which is the basis of this article.

######  Theorem 1.

Let $\beta>-1$. The residual polynomials for the Nesterov acceleration (1)
with (2) are

$r_{k}(\lambda)=(1-\lambda)^{\frac{k+1}{2}}\frac{C_{k-1}^{(\frac{\beta+1}{2})}(\sqrt{1-\lambda})}{C_{k-1}^{(\frac{\beta+1}{2})}(1)},\qquad
k\geq 1,$ (5)

with the Gegenbauer polynomials $C_{n}^{(\alpha)}$.

###### Proof.

Defining $h_{k}(\lambda)=r_{k}(\lambda)(1-\lambda)^{-\frac{k+1}{2}}$ and
multiplying (3) by $(1-\lambda)^{-\frac{k+2}{2}}$ leads to the relation

$\displaystyle h_{k+1}(\lambda)$
$\displaystyle=(1+\alpha_{k})\sqrt{1-\lambda}h_{k}(\lambda)-\alpha_{k}h_{k-1}(\lambda),\qquad
k\geq 2,$ (6) $\displaystyle h_{1}(\lambda)$ $\displaystyle=1,\qquad
h_{2}(\lambda)=\sqrt{1-\lambda}.$

We note that $C_{n}^{(\frac{\beta+1}{2})}(x)$ satisfy the recursion relation
(cf. [1, p. 782])

$\begin{split}C_{k}^{(\frac{\beta+1}{2})}(x)&=xc_{k}C_{k-1}^{(\frac{\beta+1}{2})}(x)-d_{k}C_{k-2}^{(\frac{\beta+1}{2})}(x),\qquad
k\geq 2\\\ C_{0}^{(\frac{\beta+1}{2})}(x)&=1,\qquad
C_{1}^{(\frac{\beta+1}{2})}(x)=(\beta+1)x\end{split}$ (7)

with

$c_{k}=\frac{2k+\beta-1}{k},\qquad d_{k}=\frac{k+\beta-1}{k}.$

Using the recurrence relation with $x=1$, leads to

$\displaystyle C_{k}^{(\frac{\beta+1}{2})}(1)$
$\displaystyle=c_{k}C_{k-1}^{(\frac{\beta+1}{2})}(1)\left(1-\theta_{k}^{-1}\right)=d_{k}C_{k-2}^{(\frac{\beta+1}{2})}(1)\left(\theta_{k}-1\right)$

with

$\theta_{k}:=\frac{c_{k}C_{k-1}^{(\frac{\beta+1}{2})}(1)}{d_{k}C_{k-2}^{(\frac{\beta+1}{2})}(1)}.$

Dividing (7) by $C_{k}^{(\frac{\beta+1}{2})}(1)$ and using this relation
yields

$\displaystyle\frac{C_{k}^{(\frac{\beta+1}{2})}(x)}{C_{k}^{(\frac{\beta+1}{2})}(1)}$
$\displaystyle=\frac{C_{k-1}^{(\frac{\beta+1}{2})}(x)}{C_{k-1}^{(\frac{\beta+1}{2})}(1)}\frac{1}{1-\theta_{k}^{-1}}-\frac{C_{k-2}^{(\frac{\beta+1}{2})}(x)}{C_{k-2}^{(\frac{\beta+1}{2})}(1)}\frac{1}{\theta_{k}-1}.$
(8)

By induction (or by well-known formulae [1, 11]), it can easily be verified
that
$\frac{C_{k-1}^{(\frac{\beta+1}{2})}(1)}{C_{k-2}^{(\frac{\beta+1}{2})}(1)}=\frac{k+\beta-1}{k-1},$
from which it follows that $\theta_{k}-1=\alpha_{k}^{-1}$ as well as
$1-\theta_{k}^{-1}=(1+\alpha_{k})^{-1}$. Thus,
$\frac{C_{k}^{(\frac{\beta+1}{2})}(x)}{C_{k}^{(\frac{\beta+1}{2})}(1)}$
satisfies the same recursion as $h_{k+1}(\lambda)$, and the corresponding
initial values for $k=0,1$ agree when setting $x=\sqrt{1-\lambda}$. This
allows us to conclude that

$h_{k}(\lambda)=\frac{C_{k-1}^{(\frac{\beta+1}{2})}(\sqrt{1-\lambda})}{C_{k-1}^{(\frac{\beta+1}{2})}(1)},$

which proves the theorem. ∎

This theorem relates the residual function of Nesterov acceleration to other
known iterations. In particular, the residual $r_{k}$ is roughly the product
of that of $\frac{k}{2}$ Landweber iterations and that of $\frac{k}{2}$
iterations of a $\nu$-method with $\nu=\frac{\beta+1}{4}$.

######  Remark 1.

Gegenbauer polynomials are special cases of Jacobi polynomials and themselves
embrace several other orthogonal polynomials as special cases. Certain values
of $\beta$ in (1) lead to various specializations in (5): The choice $\beta=0$
leads to Legendre polynomials, the often encountered choice $\beta=1$ leads to
Chebyshev polynomials of the second kind [1].

We note that the result of Theorem 1 even holds for $\beta=-1$. In this case,
only $\alpha_{1}$ is not well-defined, but it is always $0$ for $\beta>-1$.
Thus, we may extend the definition of the iteration to $\beta=-1$ by setting
$x_{k}^{\delta}:=\lim_{\beta\to-1}x_{k}^{\delta}$. (This just amounts to
slightly modifying (2) by setting $\alpha_{1}=0$ for $k=1$; the remaining
iteration is well-defined by (1) and (2).) In this case, we may use [1, Eq.
22.5.28] to conclude that the resulting polynomials are Chebyshev polynomials
of the first kind.

Before we proceed with the convergence analysis, we state for generality the
corresponding theorem for Nesterov iteration with a general sequence
$\alpha_{k}$.

######  Theorem 2.

Consider the iteration (1) with a positive sequence $\alpha_{k}$. Then the
corresponding residual function can be expressed as

$r_{k}(\lambda)=(1-\lambda)^{\frac{k}{2}}\frac{P_{k}(\sqrt{1-\lambda})}{P_{k}(1)},\qquad
k\geq 1,$ (9)

where $P_{k}$ is a sequence of orthogonal polynomials obeying the recurrence
relation

$\displaystyle\begin{split}P_{k+1}(x)&=c_{k}xP_{k}(x)-d_{k}P_{k-1}(x),\qquad
k\geq 1\\\ P_{0}(x)&=1,\qquad P_{1}(x)=c_{0}x\end{split}$ (10)

with $c_{n}$ and $d_{n}$ recursively defined to satisfy

$\displaystyle\begin{split}\frac{c_{1}c_{0}}{d_{1}}&=1+\frac{1}{\alpha_{1}}\\\
\frac{c_{k}c_{k-1}}{d_{k}}&=(1+\frac{1}{\alpha_{k}})(\alpha_{k-1}+1)\qquad
k\geq 2.\end{split}$ (11)

Conversely, given a sequence of orthogonal polynomials defined by the
recurrence relation (10) with given sequences $c_{n},d_{n}$. Then there exists
a sequence $\alpha_{k}$ (defined via (11)) such that the corresponding
Nesterov iteration (1) has a residual function as in (9).

###### Proof.

The function $h_{k}(\lambda):=r_{k}(\lambda)(1-\lambda)^{\frac{k}{2}}$
satisfies the recursion (6) with $h_{0}(\lambda)=1$ and
$h_{1}(\lambda)=\sqrt{1-\lambda}$ and for $k\geq 1$. As in the proof of
Theorem 1, we may conclude that (10) leads to a similar recursion as (8):

$\frac{P_{k+1}(x)}{P_{k+1}(1)}=\frac{P_{k}(x)}{P_{k}(1)}\frac{1}{1-\theta_{k}^{-1}}-\frac{P_{k-1}(x)}{P_{k-1}(1)}\frac{1}{\theta_{k}-1},\qquad
k\geq 1.$

with

$\theta_{k}=\frac{c_{k}P_{k}(1)}{d_{k}P_{k-1}(1)},\qquad k\geq 1.$

From (10) we can conclude by some algebraic manipulations that

$\theta_{k}=\frac{c_{k}c_{k-1}}{d_{k}}\left(1-\theta_{k-1}^{-1}\right),\qquad
k\geq 2.$

If (11) holds, then from the recursion for $\theta_{k}$, it follows that we
can perform an induction step following that
$\frac{1}{\theta_{k-1}-1}=\alpha_{k-1}$ implies
$\frac{1}{\theta_{k}-1}=\alpha_{k}$. Since $\frac{1}{\theta_{1}-1}=\alpha_{1}$
by definition, we obtain that $h_{k}(\lambda)$ and $\frac{P_{k}(x)}{P_{k}(1)}$
satisfy identical recursions and have identical initial conditions with the
setting $x=\sqrt{1-\lambda}$.

Conversely, if (10) is given and the sequence $\alpha_{k}$ is recursively
defined by (11), then it follows in a similar manner that
$\frac{P_{k}(x)}{P_{k}(1)}$ has the same recursion and initial conditions as
$h_{k}(\lambda)$ and thus both functions agree. ∎

The polynomials $P_{k}(x)$ in this theorem correspond to
$xC_{k-1}^{\frac{\beta+1}{2}}(x)$ in Theorem 1.

As an illustration, we may consider the peculiar choice of $\alpha_{k}$ in
Nesterov’s original paper [9], which is also used in the well-known FISTA
iteration [2]: First, a sequence is defined recursively,

$t_{k+1}=\frac{1}{2}\left(1+\sqrt{1+4t_{k}^{2}}\right),\qquad t_{1}=1,$

and then the sequence $\alpha_{k}$ is given by

$\alpha_{k}=\frac{t_{k}-1}{t_{k+1}}.$

Note that $t_{k+1}$ is the positive root of the equation
$t_{k+1}(t_{k+1}-1)=t_{k}^{2}.$ Using this identity, we may calculate that

$(1+\frac{1}{\alpha_{k}})(\alpha_{k-1}+1)=\frac{t_{k}}{t_{k-1}}(1+\frac{t_{k}}{t_{k+1}})(1+\frac{t_{k-1}}{t_{k}}).$

Thus, coefficients for a recurrence formula for orthogonal polynomials that
correspond to such an iteration are

$c_{k}=1+\frac{t_{k}}{t_{k+1}},\qquad d_{k}=c_{k-1}-1.$

However, this does not seem to be related to any common polynomial family, to
the knowledge of the author.

On the other hand, we may design Nesterov iterations from the recurrence
relation of classical polynomials. For instance, the Hermite polynomials obey
a relation (10) with $c_{k}=2,$ $d_{k}=2k$. Thus, the sequence $\alpha_{k}$
has to satisfy the recursion

$a_{k}:=\frac{1+a_{k-1}}{\frac{2}{k}-a_{k-1}-1}.$

We do not know if this is of any use, though.

##  3 Convergence analysis

We consider the iteration (1) with the usual $\alpha_{k}$-sequence (2) and
show that it is an optimal-order regularization methods (of course, when
combined with a stopping rule).

###  3.1 Convergence rates and semi-saturation

In the classical analysis of regularization schemes [4], one tries to bound
the error in terms of the noiselevel $\delta$:
$\|x_{k(\delta)}^{\delta}-x^{\dagger}\|\leq f(\delta)$, where $f$ is some
function decreasing to $0$ with $\delta\to 0$. Often, Hölder-type rates are
considered with $f(\delta)=C\delta^{\xi}$. For such estimates, one has to
impose smoothness conditions in form of a source condition

$\displaystyle
x^{\dagger}=(A^{*}A)^{\mu}\omega,\qquad\|\omega\|<\infty,\quad\mu>0.$ (12)

It is also well-known [4] that the optimal rate of convergence under (12) is
of the form

$\|x_{k(\delta)}^{\delta}-x^{\dagger}\|\leq O(\delta^{\frac{2\mu}{2\mu+1}}),$

and a regularization scheme that achieves this bound is called of optimal
order.

The phenomenon of saturation is the effect that for certain regularization
method, the convergence rate $f(\delta)$ does not improve even when the
smoothness is higher, i.e., $\mu$ is larger. This happens, for instance for
Tikhonov regularization at $\mu=1$ or for the $\nu$-methods at $\mu=\nu$; see
[4].

For the Nesterov iteration (1), a detailed analysis has been performed by
Neubauer [10] with the result that, assuming a usual source condition (12) and
an appropriate a priori stopping rule, the resulting iterative regularization
scheme is of optimal order for $\mu\leq\frac{1}{2}$, and, for
$\mu>\frac{1}{2}$, the convergence rates improve with $\mu$ but in a
suboptimal way. More precisely, the convergence rates proven in [10] are

$\|x_{k(\delta)}^{\delta}-x^{\dagger}\|=\begin{cases}O(\delta^{\frac{2\mu}{2\mu+1}})&\mu\leq\frac{1}{2},\\\
O(\delta^{\frac{2\mu+1}{2\mu+3}})&\mu>\frac{1}{2}.\end{cases}$

Thus, contrary to saturating methods, the order still improves beyond the
”saturation index” $\mu=\frac{1}{2}$ but in a suboptimal way. This is what we
call ”semi-saturation”, and, to the knowledge of the author, this has not been
observed yet for a classical regularization method. A further result of [10]
is that using the discrepancy principle as stopping rule, convergence rates
are proven, which are, however, always suboptimal.

Our second main contribution is an improvement of Neubauer’s result in the
sense that we show that the Nesterov iteration is of optimal order for a
smoothness index $\mu\leq\frac{\beta+1}{4}$ with an a priori stopping rule.
Moreover, contrary to [10], we also obtain optimal-order rates with the
discrepancy principle provided that $\mu\leq\frac{\beta-1}{4}$. These findings
allows one to achieve always optimal-order convergence provided $\beta$ is
chosen sufficiently large.

Moreover, the phenomenon of semi-saturation is made transparent by referring
to the representation in Theorem 1: The residual is a product of Landweber-
type and $\nu$-type residuals, and keeping in mind that Landweber iteration
does not show saturation for Hölder indices while the $\nu$-method do, it is
clear that a product as in (5) leads to the above described semi-saturation.

###  3.2 Convergence analysis

In this section we perform a convergence analysis for the iteration (1). By
Theorem 1, we may base our investigation on the known results for Landweber
iteration and the $\nu$-methods.

We collect some useful known estimates:

$\left|\frac{C_{k-1}^{(\frac{\beta+1}{2})}(\sqrt{1-\lambda})}{C_{k-1}^{(\frac{\beta+1}{2})}(1)}\right|\leq
1,\qquad 0\leq\lambda\leq 1,\beta>-1.$ (13)

This is well-known and follows from [11, Eq. (7.33.1), (4.73)]. From this we
immediately obtain that

$|r_{k}(\lambda)|\leq 1,\qquad 0\leq\lambda<1,\quad\beta>-1,$ (14)

which has already been shown in [10]. Moreover, we may conclude from (13) and
(5) as well that

$\lim_{k\to\infty}r_{k}(\lambda)\to 0,\qquad 0<\lambda<1.$ (15)

Recall that we denote by $x_{k}$ the iteration with $y^{\delta}$ replaced by
the exact data. As usual, this allows one to split the total error into an
approximation and stability term. We estimate the stability term:

######  Proposition 1.

Let $\|A^{*}A\|\leq 1$ and define $x_{k}^{\delta}$ by (1) (2) with $\beta>-1$.
Let $x_{k}$ be the corresponding noise-free iteration with $y^{\delta}$
replaced by $y=Ax^{\dagger}$. Then we have the estimate

$\|x_{k}^{\delta}-x_{k}\|\leq\sqrt{2}\sqrt{(k-1)^{2}+\frac{k+1}{2}}\delta\leq\sqrt{2}k\delta.$
(16)

###### Proof.

Following [4], it is enough to estimate

$g_{k}(\lambda)=\frac{1-r_{k}(\lambda)}{\lambda}=r_{k}^{\prime}(\tilde{\lambda}),$

where we used the mean value theorem with $\tilde{\lambda}\in(0,\lambda)$. The
derivative may be calculated from (5) as

$\displaystyle r_{k}^{\prime}(\lambda)$
$\displaystyle=\frac{k+1}{2}(1-\lambda)^{\frac{k-1}{2}}\frac{C_{k-1}^{(\frac{\beta+1}{2})}(\sqrt{1-\lambda})}{C_{k-1}^{(\frac{\beta+1}{2})}(1)}-\frac{1}{2}(1-\lambda)^{\frac{k}{2}}\left(\frac{[C_{k-1}^{(\frac{\beta+1}{2})}](\sqrt{1-\lambda})}{C_{k-1}^{(\frac{\beta+1}{2})}(1)}\right)^{\prime}.$

We use Markov’s inequality (cf. [4, Eq. (6.16)]) and (13) to conclude that

$\left|\left(\frac{[C_{k-1}^{(\frac{\beta+1}{2})}](\sqrt{1-\lambda})}{C_{k-1}^{(\frac{\beta+1}{2})}(1)}\right)^{\prime}|\right|\leq
2(k-1)^{2}\max_{0\leq\lambda\leq
1}\left|\frac{C_{k-1}^{(\frac{\beta+1}{2})}(\sqrt{1-\lambda})}{C_{k-1}^{(\frac{\beta+1}{2})}(1)}\right|\leq
2(k-1)^{2}.$

Thus,

$|g_{k}(\lambda)|\leq\frac{k+1}{2}+(k-1)^{2}.$

The result now follows with [4, Theorem 4.2] and (14). ∎

Note that this estimate is a slight improvement compared to the corresponding
estimate in [10, Equation (3.2)], which has $2k\delta$ on the right-hand side,
similar as for the $\nu$-methods.

From this we may conclude convergence:

######  Theorem 3.

Let $\|A^{*}A\|\leq 1$ and $\beta>-1$. If the iteration is stopped at a
stopping index $k(\delta)$ that satisfies $k(\delta)\delta\to 0$ and
$k(\delta)\to\infty$ as $\delta\to 0$, then we obtain convergence

$x_{k(\delta)}^{\delta}\to x^{\dagger}.$

###### Proof.

We estimate

$\displaystyle\|x_{k(\delta)}^{\delta}-x^{\dagger}\|$
$\displaystyle\leq\|x_{k(\delta)}^{\delta}-x_{k(\delta)}\|+\|x_{k(\delta)}-x^{\dagger}\|\leq\sqrt{2}k(\delta)\delta+r_{k(\delta)}(A^{*}A)x^{\dagger}.$

The first term converges to $0$ by assumption on $k(\delta)$ and the second
term does so because $k(\delta)\to\infty$ and by the dominated convergence
theorem using (14), (15) as in [4]. ∎

We now consider convergence rates, and for this, the following rather deep
estimate for orthogonal polynomials is needed. It was derived by Brakhage [3]
as well as by Hanke [4, Appendix A.2], [6] on basis of Hilb-type estimates for
Jacobi polynomials.

######  Proposition 2.

Let $\beta>-1$. Then there is a constant $c_{\beta}$ with

$\left|\lambda^{\frac{\beta+1}{4}}\frac{C_{k}^{(\frac{\beta+1}{2})}(\sqrt{1-\lambda})}{C_{k}^{(\frac{\beta+1}{2})}(1)}\right|\leq
c_{\beta}k^{-2\frac{\beta+1}{4}},\qquad 0\leq\lambda\leq 1.$ (17)

###### Proof.

For $k$ even, this is [4, Eq. (6.22)] (with $k$ there meaning $2k$ here), or
[6, Th. 4.1]. However, the result there is based on the Hilb-type formula
([11, Theorems 8.21.12, 8.21.13] which holds for all $k$ as in [3, p. 170].
Thus, by following the steps in [4, Appendix A.2], the result is obtained. ∎

Note that in case $-1<\beta<1$, the constant $c_{\beta}$ may be explicitly
calculated from [11, Eq. (7.33.5)].

The corresponding estimates for the residuals of Landweber iteration are
standard; cf. [4, Eq. (6.8)]:

$|\lambda^{\mu}(1-\lambda)^{k}|\leq c_{\mu}(k+1)^{-\mu}.$ (18)

As a consequence, we may state our main convergence rate result for an a
priori stopping rule:

######  Theorem 4.

Let $\|A^{*}A\|\leq 1$ and $\beta>-1$, and suppose that a source condition
(12) is satisfied with some $\mu>0$.

  1. 1.

If $\mu\leq\frac{\beta+1}{4}$ and the stopping index is chosen as

$k(\delta)=O(\delta^{-\frac{1}{2\mu+1}}),$

then we obtain optimal order convergence

$\displaystyle\|x_{k(\delta)}^{\delta}-x^{\dagger}\|\leq
O(\delta^{\frac{2\mu}{2\mu+1}}).$

  2. 2.

If $\mu>\frac{\beta+1}{4}$ and the stopping index is chosen as

$k(\delta)=O(\delta^{-\frac{1}{\mu+\frac{\beta+1}{4}+1}}),$ (19)

then we obtain suboptimal order convergence

$\displaystyle\|x_{k(\delta)}^{\delta}-x^{\dagger}\|\leq
O(\delta^{\frac{\mu+\frac{\beta+1}{4}}{\mu+\frac{\beta+1}{4}+1}}).$

###### Proof.

For $\lambda\leq 1$ the estimate (17) yields (by interpolation) and
$(1-\lambda)^{\frac{k+1}{2}}\leq 1$ that

$|r_{k}(\lambda)\lambda^{\mu}|\leq Ck^{-2\mu},\qquad\mu\leq\frac{\beta+1}{4}.$
(20)

In case of $\mu>\frac{\beta+1}{4}$, we have with additionally using (18)

$\displaystyle|r_{k}(\lambda)\lambda^{\mu}|$
$\displaystyle\leq|(1-\lambda)^{\frac{k+1}{2}}\lambda^{\mu-\frac{\beta+1}{4}}|c_{\beta}k^{-2\frac{\beta+1}{4}}$
$\displaystyle\leq
c_{\mu,\beta}(\frac{k+1}{2})^{-(\mu-\frac{\beta+1}{4})}|c_{\beta}k^{-2\frac{\beta+1}{4}}\leq
Ck^{-\left(\mu+\frac{\beta+1}{4}\right)}.$

The result now follows by standard means:

$\displaystyle\|x_{k}^{\delta}-x^{\dagger}\|$
$\displaystyle\leq\|x_{k}^{\delta}-x_{k}^{\delta}\|+\|x_{k}^{\delta}-x^{\dagger}\|\leq\sqrt{2}k\delta+\|r_{k}(A*A)(A^{*}A)^{\mu}\omega\|$
$\displaystyle\leq\begin{cases}\sqrt{2}k\delta+Ck^{-2\mu}&\mu\leq\frac{\beta+1}{4},\\\
k\delta+Ck^{-(\mu+\frac{\beta+1}{4})}&\mu>\frac{\beta+1}{4}.\end{cases}$

Solving for $k$ by equating the two terms in the last bounds yields the a
priori parameter choice and the corresponding rates. ∎

These results correspond to those of Neubauer when $\beta=1$. However, for
$\beta>1$ this is an improvement as we obtain optimal-order convergence if
$\beta$ is chosen larger than $4\mu-1$. We note that in the optimal order
case, the number of iteration needed is of order
$O(\delta^{-\frac{1}{2\mu+1}})$, which is the same order as for semiiterative
methods and for the conjugate gradient method. Thus, the Nesterov acceleration
certainly qualifies being called a fast method.

###  3.3 Discrepancy principle

With the improved estimates, we can as well strengthen the result of [10] when
the iteration is combined with the well-known discrepancy principle. Recall
that it defines a stopping index $k(\delta)$ a posteriori by the first
(smallest) $k$ that fulfils the inequality

$\|Ax_{k}^{\delta}-y^{\delta}\|\leq\tau\delta,$ (21)

where $\tau>1$ is fixed. The corresponding convergence rates can be obtained
by a slight modification of the proof in [10] and the general theory in [4].

######  Theorem 5.

Let $\|A^{*}A\|<1$, $\beta>-1$, and assume a source condition (12) satisfied.
If the iteration (1) is stopped by the discrepancy principle (21), then we
obtain the following convergence rates:

  1. 1.

If $\mu+\frac{1}{2}\leq\frac{\beta+1}{4}$, then we achieve optimal order
convergence rates

$\displaystyle\|x_{k(\delta)}^{\delta}-x^{\dagger}\|\leq
O(\delta^{\frac{2\mu}{2\mu+1}})$

with stopping index being of the same order as in (19).

  2. 2.

$\mu+\frac{1}{2}\geq\frac{\beta+1}{4}$, then we obtain that

$k(\delta)=O(\delta^{-\frac{1}{\frac{1}{2}+\mu+\frac{\beta+1}{4}}})$

and a rate of

$\displaystyle\|x_{k(\delta)}^{\delta}-x^{\dagger}\|\leq
O(\delta^{\frac{\mu+\frac{\beta+1}{4}-\frac{1}{2}}{\mu+\frac{\beta+1}{4}+\frac{1}{2}}}).$

###### Proof.

The proof [10, Theorem 4.1] only needs minor modifications. The estimate [10,
Eq. (4.3)]

$\|x_{k(\delta)}-x^{\dagger}\|\leq\|r_{k(\delta)}(T^{*}T)w\|^{\frac{1}{2\mu+1}}\left((\tau+1)\delta\right)^{\frac{2\mu}{2\mu+1}}$

is valid independent of our new rate results, hence it follows as in [10, Eq.
(4.4)] that $\|x_{k(\delta)}-x^{\dagger}\|\leq
o(\delta^{\frac{2\mu}{2\mu+1}})$. It remains to estimate
$\|x_{k(\delta)}^{\delta}-x^{\dagger}\|$ by (16) combined with an upper bound
for $k(\delta)$. Estimate [10, Eq. (4.2)] and the discrepancy principle yields

$\tau\delta\leq\delta+\|(T^{*}T)^{\frac{1}{2}+\mu}r_{k}(T^{*}T)w\|$

for $k=k(\delta)$. Using (20) in case that
$\mu+\frac{1}{2}\leq\frac{\beta+1}{4}$, we obtain

$(\tau-1)\delta\leq Ck(\delta)^{-2(\frac{1}{2}+\mu)},$

which yields (19), and with (16) we obtain
$\|x_{k(\delta)}^{\delta}-x^{\dagger}\|=O(\delta^{\frac{2\mu}{2\mu+1}})$,
which proves the result in the optimal case.

In case that $\mu+\frac{1}{2}>\frac{\beta+1}{4}$, the corresponding estimate
is

$(\tau-1)\delta\leq Ck(\delta)^{-(\frac{1}{2}+\mu+\frac{\beta+1}{4})},$

from which the result in the second case follow. ∎

These rates agree with those of [10] when setting $\beta=1$. There, however,
only the suboptimal case 2. was possible. Our improvement is to show that we
may achieve optimal order results even with the discrepancy principle provided
$\beta$ is sufficiently large.

######  Remark 2.

It is clear that in practice $\beta$ should be selected in the regime of
optimal rates, i.e. $\beta>4\mu-1$ for a prior choices and $\beta>4\mu+1$ for
the discrepancy principle. However, it is a rule of thumb to choose such
parameter also as small as possible, or more precisely, in such a way to come
close to the saturation point, i.e. $\beta\sim 4\mu-1$, respectively
$\beta\sim 4\mu+1$.

######  Remark 3.

For semiiterative methods, a modified discrepancy principle [6, 4] has been
defined, where the residual in (21) is replaced by an expression of the form
$(y^{\delta},s_{k}(AA^{*})y^{\delta})$ with a constructed function $s_{k}$.
This yields an order-optimal method as for the a priori stopping rule. An
adaption of this strategy for Nesterov iteration is certainly possible and
this should yield order-optimal rates for all $\mu\leq\frac{\beta+1}{4}$.
However, the strategy is quite involved and it is not completely clear to us
how to include this into the iteration efficiently. We thus do not intend to
investigate such modifications in this article.

##  4 Numerical results

In this section we present some small numerical experiments to illustrate the
semi-saturation phenomenon and to investigate the performance of Nesterov’s
iteration, in particular, with respect to the optimal-order results.

In a first example we consider a simple diagonal operator $A={\rm
diag}(\frac{1}{n^{2}})$, for $n=1,\ldots 1000$, as well as an exact solution
$x^{\dagger}=(\frac{1}{n^{4}}(-1)^{n})_{n=1}^{1000}$, which amounts to a
source condition being satisfied with index $\mu=0.75$. Thus, we are in a case
of higher smoothness, where the results of the present article really improve
those of [10]. We add standard normally distributed Gaussian noise to the
exact data and performed various iterative regularization schemes: Landweber
iteration, the $\nu$-method, and the Nesterov iteration, the latter two with
various settings of the parameters $\nu$ and $\beta$, respectively.

We calculated the stopping index either by the discrepancy principle (21) with
$\tau=1.01$ or, since we have the luxury of an available exact solution in
this synthetic example, we also calculate the oracle stopping index, which is
defined as

$k_{opt}={\rm argmin}_{k}\|x_{k}^{\delta}-x^{\dagger}\|.$

In other words, $k_{opt}$ is the theoretically optimal possible stopping
index.

Figure 1: Log-log plot of the error $\|x_{k(\delta)}^{\delta}-x^{\dagger}\|$
versus the noiselevel $\delta$ for Nesterov iteration (full line, blue),
Landweber iteration (dotted line, black), and the $\nu$-method (dashed dotted
line, red). Left: optimal stopping rule. Right: stopping by discrepancy
principle. The parameters $\beta,\nu$ are in an optimal-order regime.

Figure 2: The similar plot as in Figure 1(left) for various iteration
parameter in a suboptimal-order regime. Left: $\nu=0.4$ and $\beta=0$. Right:
$\nu=0.3$ and $\beta=-0.5$. Stopping by optimal stopping rule $k_{opt}$.

In Figure 1, we display the error $\|x_{k(\delta)}^{\delta}-x^{\dagger}\|$
against various noiselevels on a log-log scale. The curves correspond to
convergence rates for Nesterov iteration (full line, blue), Landweber
iteration (dotted line, black), and the $\nu$-method (dashed dotted line,
red). The parameter were chosen as $\beta=4$ and $\nu=1$, i.e., we are in the
optimal-order case covered by item 1 in Theorem 4 and Theorem 5. On the left-
hand side we employ the oracle stopping rule using $k_{opt}$ and on the right-
hand side we use the discrepancy principle.

As can be observed, all three methods show a similar (optimal-order) rate, as
stated in Theorems 4 and 5. In particular, this verifies one of our findings
that the discrepancy principle for Nesterov’s iteration leads to an optimal-
order method provided $\beta$ is chosen appropriately.

In Figure 2, we illustrate the semi-saturation phenomenon: Here $\beta$ and
$\nu$ are deliberately chosen as too small ($\beta=0$, $\nu=0.4$ on the left-
hand side and $\beta=-0.5$, $\nu=0.3$ on the right-hand side). We observe that
for small $\nu$, the convergence rate of the $\nu$-method is slow as a result
of its saturation. On the other hand, the Nesterov iteration also has a slower
rate than the non-saturating Landweber iteration, but, as can be expected from
our residual polynomial representation, it is in between the other two.

We remark that the $\nu$-methods show some unpleasant behaviour when $\nu$ is
chosen small. The residual is highly oscillating and for small noiselevel we
could not even reach the prescribed discrepancy, and if we did, then the
number of iteration was quite high, even higher than for Landweber iteration.
This might be attributed to our quite aggressive setting of the discrepancy
principle with $\tau=1.01$. In that respect, the Nesterov iteration was very
well-behaved, and we had no problem with a small $\beta$, which is probably
due to the robust Landweber-component in the representation (5).

The optimal-order convergence only partly illustrates the effective
performance of the methods. In Tables 1 we therefore provide the ratio of
errors values, i.e., the numbers in the table are
$\frac{\|x_{method,k}^{\delta}-x^{\dagger}\|}{\|x_{Nesterov,k_{opt}}^{\delta}-x^{\dagger}\|}$,
where $x_{Nesterov,k_{opt}}^{\delta}$ denotes Nesterov iteration with the
optimal stopping rule and $x_{method,k}^{\delta}$ the iteration of the
respective method with the respective stopping rule. All results correspond to
an optimal-order regime of parameters (those of Figure 1). The number of
iterations (both for the oracle stopping rule and the discrepancy principle)
are given in Table 2. In these tables, we also include the corresponding
results for the conjugate gradient iteration CGNE [7].

Table 1: Errors compared to Nesterov iteration: $\frac{\|x_{method,k}^{\delta}-x^{\dagger}\|}{\|x_{Nesterov,k_{opt}}^{\delta}-x^{\dagger}\|}$.  |  | $\delta$
---|---|---
Method | Stopping | $10^{-5}$ | $10^{-4}$ | $10^{-3}$ | $10^{-2}$ | $10^{-1}$
Nesterov | $k_{opt}$ | 1 | 1 | 1 | 1 | 1
Landweber | $k_{opt}$ | 1.15 | 0.83 | 0.96 | 1.05 | 1.06
$\nu$-Method | $k_{opt}$ | 1.02 | 1.06 | 1.01 | 1.26 | 0.97
CGNE | $k_{opt}$ | 1.02 | 0.82 | 1.05 | 1.02 | 0.84
Nesterov | Discrepancy | 1.58 | 1.10 | 1.41 | 2.84 | 1.90
Landweber | Discrepancy | 2.23 | 1.17 | 1.41 | 2.80 | 1.98
$\nu$-Method | Discrepancy | 1.02 | 1.13 | 1.00 | 1.56 | 1.88
CGNE | Discrepancy | 1.81 | 1.19 | 1.05 | 2.51 | 1.97
Table 2: Number of iterations for various methods; setting as in Table 1. |  | $\delta$
---|---|---
Method | Stopping | $10^{-5}$ | $10^{-4}$ | $10^{-3}$ | $10^{-2}$ | $10^{-1}$
Nesterov | $k_{opt}$ | 371 | 163 | 65 | 26 | 15
Landweber | $k_{opt}$ | 11000 | 2193 | 512 | 145 | 36
$\nu$-Method | $k_{opt}$ | 190 | 82 | 33 | 22 | 9
CGNE | $k_{opt}$ | 10 | 6 | 4 | 3 | 2
Nesterov | Discrepancy | 260 | 111 | 39 | 13 | 1
Landweber | Discrepancy | 5106 | 1080 | 220 | 37 | 1
$\nu$-Method | Discrepancy | 190 | 96 | 33 | 10 | 1
CGNE | Discrepancy | 8 | 5 | 4 | 2 | 1

In terms of the number of iteration, the Nesterov iteration is slightly slower
than the $\nu$-methods (approximately by a constant factor of 1.5) but both
have a similar modest increase of iterations when $\delta$ is decreased. Both
need more iteration than the CGNE-method, which, of course, is the fastest one
by design. The slightly higher number of iterations might be attributed to the
better error estimate in (16). (Note that the $\nu$-methods have a $2$ in
place of $\sqrt{2}$ there). It might appear a little bit paradoxical that a
better estimate leads to slower convergence, but this is clear from the theory
as the number of iteration is a decreasing function of $\delta$ and thus also
of any factor in front of $\delta$. This factor, however, pays off when
considering the total error of the method, and we observe that Nesterov
iteration with the optimal choice $k_{opt}$ indeed has almost always a
slightly smaller error than the $\nu$-method. Surprisingly, it is in several
instances also better than the CGNE-method. However, the Nesterov method
sometimes loses some of its advantages against the $\nu$-method, when using
the discrepancy principle, but the performance is still acceptable.

Some further experiments indicate that the results are rather insensitive to
overestimating $\beta$. As stated in Remark 2, the best choice is usually
related to the smoothness index, but there was arose no serious problems when
$\beta$ was larger.

Further numerical experiments have been performed in [10]: Even though the
value of $\beta$ was not reported there, the results are consistent with our
theory with the choice $\beta=1$. The forward operator there was the Green’s
function for the solution of the 1D boundary value problem
$-u^{\prime\prime}=f$ with homogeneous boundary conditions. Exact solutions
with various smoothness are stated there: Example 5.1 with $\mu=\frac{1}{8}$,
Example 5.2 with $\mu=\frac{5}{8}$, and Example 5.3 with $\mu=\frac{17}{8}$.
We used the same problem and the same examples, but we calculated $A$ by using
a FEM-discretization of the boundary value problem and $A$ as the
corresponding solution operator. For simplicity we ignored discretization
errors and took the discretized (projected) solution as $x^{\dagger}$.

Figure 3: Convergence rates for the examples in [10]. Left: Ex. 1, smoothness
index $\mu=\frac{1}{8}$. Center: Ex. 2, smoothness index $\mu=\frac{5}{8}$.
Right: Ex. 3, smoothness index $\mu=\frac{17}{8}$. Displayed are the errors
versus the noiselevel on a logarithmic scale. A marker ’x’ indicates optimal
choice of $\beta$, and ’+’ indicates suboptimal choice $\beta=1$. The full
line indicates the optimal order rate.

The main purpose of this experiment is to verify that the discrepancy
principle ($\tau=1.1$) can be made an optimal-order method. We choose
$\beta=3.5$ for the first two examples and $\beta=9.5$ for the third, which
should in any case lead to an optimal-order situation. In Figure 3, we plotted
the error versus the relative noiselevel on a logarithmic scale for the three
examples with this choice of $\beta$, indicated by the marker ’x’. As a
comparison, we also indicated the predicted optimal rate by a solid line.
Furthermore, also shown and marked with ’+’ are the corresponding results for
$\beta=1$, i.e., in the suboptimal case.

These results clearly illustrate that for the discrepancy principle we may
achieve the optimal order rates with the correct choice of $\beta$ and for a
wrong choice of $\beta$ the rate deteriorates. For low-smoothness as in
Example 1 (left picture in Figure 3, however, there seems to occur almost no
deterioration contrary to expectation.

##  5 Conclusion

We have provided a representation of the residual polynomials for Nesterov’s
acceleration method for linear ill-posed problems as a product of Gegenbauer
polynomials and Landweber-type residuals. This allowed us to prove optimal-
order rates for an a priori stopping rule and the discrepancy principle as
long as $\beta$ in (2) is sufficiently large. The number of iteration is shown
to be of the same order as for other fast methods such as the $\nu$-method or
the conjugate gradients methods. Moreover, our representation clearly explains
the observed semi-saturation phenomenon.

Within the class of linear iterative methods, the Nesterov acceleration is an
excellent choice, as it is a fast method as well as a quite robust one.
Although, it must be conceded, that it cannot compete with the conjugate
gradient method in terms of number of iterations. However, this is compensated
by its flexibility and simplicity of use, which also allows one to easily
integrate it into existing gradient methods and also to apply it in nonlinear
cases.

## References

  * [1] M. Abramowitz and I. A. Stegun, Handbook of mathematical functions with formulas, graphs, and mathematical tables, vol. 55 of National Bureau of Standards Applied Mathematics Series, U.S. Government Printing Office, Washington, D.C., 1964.
  * [2] A. Beck and M. Teboulle, A fast iterative shrinkage-thresholding algorithm for linear inverse problems, SIAM J. Imaging Sci., 2 (2009), pp. 183–202.
  * [3] H. Brakhage, On ill-posed problems and the method of conjugate gradients, in Inverse and ill-posed problems (Sankt Wolfgang, 1986), vol. 4 of Notes Rep. Math. Sci. Engrg., Academic Press, Boston, MA, 1987, pp. 165–175.
  * [4] H. W. Engl, M. Hanke, and A. Neubauer, Regularization of inverse problems, vol. 375 of Mathematics and its Applications, Kluwer Academic Publishers Group, Dordrecht, 1996.
  * [5] J. Favard, Sur les polynomes de tchebicheff, C.R. Acad. Sci. Paris, 200 (1935), pp. 2052–2053.
  * [6] M. Hanke, Accelerated Landweber iterations for the solution of ill-posed equations, Numer. Math., 60 (1991), pp. 341–373.
  * [7] , Conjugate gradient type methods for ill-posed problems, vol. 327 of Pitman Research Notes in Mathematics Series, Longman Scientific & Technical, Harlow, 1995.
  * [8] S. Hubmer and R. Ramlau, Nesterov’s accelerated gradient method for nonlinear ill-posed problems with a locally convex residual functional, Inverse Problems, 34 (2018), pp. 095003, 30 pages.
  * [9] Y. E. Nesterov, A method for solving the convex programming problem with convergence rate $O(1/k^{2})$, Dokl. Akad. Nauk SSSR, 269 (1983), pp. 543–547.
  * [10] A. Neubauer, On Nesterov acceleration for Landweber iteration of linear ill-posed problems, J. Inverse Ill-Posed Probl., 25 (2017), pp. 381–390.
  * [11] G. Szegő, Orthogonal polynomials, American Mathematical Society, Providence, R.I., fourth ed., 1975.  American Mathematical Society, Colloquium Publications, Vol. XXIII.