New and not so new methods for estimating a regularization parameter - - PowerPoint PPT Presentation

New and not so new methods for estimating a regularization parameter

Giuseppe Rodriguez

in collab. with C. Brezinski, C. Fenu, P.C. Hansen, T.K. Jensen,

• M. Hochstenbach, Y. Park, L. Reichel, H. Sadok, S. Seatzu, X. Yu

Department of Mathematics and Computer Science, University of Cagliari, Italy

INdAM intensive period - CMIPI 2018

Computational Methods for Inverse Problems in Imaging Como, Italy - July 16–18, 2018

• G. Rodriguez

New and not so new estimation of a regularization parameter

Integral equations of the first kind

Linear Fredholm integral equations of the first kind b

a

k(x, y)f (y) dy = g(x), x ∈ [c, d], are the main source of ill-posed problems which are encountered in many applications, e.g., tomography, image deblurring, geophysical prospecting, mechanical engineering, etc. They can be discretized in various ways (quadrature, projection, collocation, . . . ), leading to linear discrete ill-posed problems Kf = g.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Linear discrete ill-posed problems (LDIP)

We will write Kf = g and Kf = g to denote either the continuous or the discrete problems. The matrix K is not necessarily square and in general it is severely ill-conditioned, due to the fact that the singular values of the integral operator quickly converge to zero. We will denote by gδ or gδ a noisy right-hand side. The vector f† = K †g represents the minimal norm least squares (MNLS) solution, i.e., min

f∈S f,

S = {f | min

f

Kf − g}. K † is the Moore-Penrose generalized inverse of K.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Singular Value Decomposition

Let us introduce the singular value decomposition (SVD) K = UΣV T, where UTU = Im, V TV = In, Σ = diag[σ1, σ2, . . . , σn], and σ1 ≥ σ2 ≥ · · · ≥ σr > σr+1 = · · · = σn = 0. For a LDIP cond(K) = σ1 σr ≫ 1

5 10 15 20 10

−20

10

−15

10

−10

10

−5

10

• G. Rodriguez

New and not so new estimation of a regularization parameter

LDIP analysis via SVD

The minimum norm least squares solution can be written as f† =

r

• j=1

uT

j g

σj vj. If gδ = g + e, fδ =

r

• j=1

uT

j g

σj vj +

r

• j=1

uT

j e

σj vj = f† + noise propagation For white noise uT

j e does not change much with j, while σj → 0

• quickly. As a consequence, the high frequency components (j ≫ r)
• f the solution are completely overwhelmed by noise propagation.
• G. Rodriguez

New and not so new estimation of a regularization parameter

Regularization methods

Let us consider the problem of computing the MNLS solution under the assumption gδ − g ≤ δ. Definition The family of continuous operators {Rα}, α ∈ (0, α0), is a regularization operator if for any g ∈ D(K †) there is a parameter choice rule α = α(δ, gδ) such that Rα(δ,gδ)gδ − K †g

δ→0

− − → 0, whenever gδ − g ≤ δ. Therefore, a regularization method is a pair (Rα, α(δ, gδ)).

• G. Rodriguez

New and not so new estimation of a regularization parameter

Truncated Singular Value Decomposition

The minimum norm least squares solution can be written as f† =

r

• j=1

uT

j g

σj vj. The TSVD solution is produced by the decomposition fδ =

k

• j=1

uT

j gδ

σj vj +

r

• j=k+1

uT

j gδ

σj vj = regularized solution fk + high frequency component k = 1, 2, . . . , r is a regularization parameter to be determined.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Tikhonov regularization

The Tikhonov approach uses a continuous parameter λ ∈ R+ to balance between the minimization of the residual and the control

• f the solution norm

min

f {Kf − gδ2 + λ2f2}.

In this case, the regularized solution is given by fλ = (ATA + λ2I)−1ATgδ. It can be computed by means of the SVD xλ =

r

• j=1

σ2

j

σ2

j + λ2

uT

j gδ

σj vj (which contains the filter factors) or by an iterative method.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Iterative regularization

Some iterative methods (Landweber, CGLS, LSQR, . . . ) exhibit semi-convergence: the low-frequency components converge first. The iteration index k plays the role of a regularization parameter.

20 40 60 80 100 10

−3

10

−2

10

−1

10 10

1

Error iterations

Gaussian(10−4) matrix: 2000 × 1000, δ = 10−6, LSQR

• G. Rodriguez

New and not so new estimation of a regularization parameter

Regularization methods in general form

A priori information about the solution can be takein into account by including a model solution f0 and a regularization matrix L, whose kernel (approximately) contains desirable solutions. Typical choices for L are the discretization of the first and second

• derivative. The problems to be solved become:

Minimal norm solution of LDIP    min

f

L(f − f0)2

• subj. to min

f

Kf − g Tikhonov min

f {Kf − g2 + λ2L(f − f0)2}

• G. Rodriguez

New and not so new estimation of a regularization parameter

Reguarization parameter choice rules

α = α(δ, gδ) A-priori rules They depend upon the noise level δ, but not upon the data gδ. A-posteriori rules They depend upon both the noise level and the data. Morozov discrepancy principle: fα such that α = min{α : Kfα − gδ ≤ τδ}, with τ ≥ 1. Heuristic (or error-free) rules They only depend upon the data gδ.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Heuristic methods

Why heuristic? A choice rule cannot depend only upon the data [Bakushinsky, 1984], so no convergence result can be proved. Motivation In many applications the noise is not white (i.e., non Gaussian, incorrelated, equally distributed) and δ is unknown. Discrepancy sometimes does not perform well for large noise levels and inconsistent problems. Common strategies:

estimate the noise level from the data; detect (semi-)norm explosion; estimate the error in the solution.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Generalized Cross Validation (GCV)

It consists of minimizing the function [Craven, Whaba 1979] V (α) =

1 mKfα − gδ2

1

m Trace(I − K(α))

2 where the influence matrix is such that K(α)g = Kf, that is K(α) = K(K ∗K + α2I)−1K ∗. The GCV is a predictive mean-square error criterion, in the sense that it estimates the minimizer of the residual T(α) = 1 mKfα − g2. It is based on the leave-out-one lemma.

• G. Rodriguez

New and not so new estimation of a regularization parameter

GCV for large scale Tikhonov regularization

For small problems, the function V (λ) is generally computed by the SVD but there are iterative procedures for large scale matrices [Golub, von Matt 1997], [Sidje, Williams, Burrage 2008]. In [Fenu, Reichel, R 2016] the numerator and the denominator of the GCV function were expressed as Stieltjes integrals, and approximated by Gauss/Gauss-Radau quadrature rules. Exploiting the connection between quadrature rules and the (global) Lanczos method, it was possible to construct a sequence

• f upper and lower bounds

Lp,q ≤ V(λ) ≤ Up,q, where (p, q) are the (Lanczos / global Lanczos) steps.

• G. Rodriguez

New and not so new estimation of a regularization parameter

GCV via global Golub–Kahan factorization

It is possible to prove that fµ(t) := µ2 t + µ2 ⇒ fµ(KK T) = Im − K(µ), so we break the computation of the GCV denominator as trace(f (KK T)) =

• m
• j=1

trace(E T

j f (KK T)Ej).

We write Ifµ := trace(W Tfµ(KK T)W ), because one can show that Ifµ = W 2

F k

• j=1

∞ fµ(λ) d wj(λ) = W 2

F

∞ fµ(λ) d w(λ).

• G. Rodriguez

New and not so new estimation of a regularization parameter

GCV via global Golub–Kahan factorization

Computing ℓ steps of the global Golub–Kahan decomposition method with U1 = W /W F gives K[V1, V2, . . . , Vℓ] = [U1, U2, . . . , Uℓ] Cℓ + σℓ+1Uℓ+1E T

ℓ ,

K T[U1, U2, . . . , Uℓ] = [V1, V2, . . . , Vℓ] C T

ℓ ,

where Ui, Uj = Vi, Vj = δij, Cℓ = Cℓ ⊗ Ik, Cℓ lower bidiagonal. It turns out that Gℓfµ = W 2

F eT 1 fµ(Tℓ)e1,

is a Gaussian quad-rule, where Tℓ := Tℓ ⊗ Ik with Tℓ := CℓC T

ℓ .

We similarly obtain the Gauss–Radau rule Rℓ+1,ζfµ = W 2

FeT 1 fµ(Tℓ+1,ζ)e1,

and for a suitable ζ the two rules bound Ifµ.

• G. Rodriguez

New and not so new estimation of a regularization parameter

GCV via partial SVD decomposition

In [Fenu, Reichel, R, Sadok 2017] bounds for the GCV were

• btained by a partial SVD decomposition of K:

uk − rk ≤ Kfλ − gδ2 ≤ uk + ck

uk =

k

• j=1

λ2 g 2

j

(σ2

j + λ)2 + ck,

ck = g2 −

k

• j=1
• g 2

j ,

rk = σ2

k(σ2 k + 2λ)

(σ2

k + λ)2 ck

and wk − sk ≤ trace(I − K(λ)) ≤ wk

wk = m −

k

• j=1

σ2

j

σ2

j + λ,

sk = (n − k) σ2

k

σ2

k + λ

to obtain Lk ≤ V(λ) ≤ Uk.

• G. Rodriguez

New and not so new estimation of a regularization parameter

GCV for large scale Tikhonov regularization

10 -6 10 -4 10 -2 10 0 10 2 10 -10 10 -8 10 -6 10 -4 10 -2 10 0 10 -20 10 -16 10 -12 10 -8 10 -4 10 0 10 -6 10 -5 10 -4 10 -3 10 -2 10 -1

k=1 k=3 k=5 k=7 k=9 k=11 k=13

Shaw example: noise 10−3, left quadrature, right gcvpsvd.

• G. Rodriguez

New and not so new estimation of a regularization parameter

L-curve

The L-curve [Hansen 1992], [Hansen, O’Leary 1993] is defined as the graph that connects adjacent points in the sequence {log Kfj − g, log fj}, j = 1, 2, . . . , ℓ. The graph is “generally” L-shaped. The L-curve criterion consist of selecting the value of the regularization parameter corresponding to the corner of the curve.

Regi´ nska method [Regi´ nska 1996] Condition L-curve [Calvetti, Lewis, Reichel 2002] L-triangle [Castellanos, G´

• mez, Guerra 2002]

Corner algorithm [Hansen, Jensen, R 2006] Residual L-curve [Reichel, Sadok 2008] restricted Regi´ nska [Reichel, R 2012]

• G. Rodriguez

New and not so new estimation of a regularization parameter

L-curve

10

−8

10

−6

10

−4

10

−2

10 10 10

2

10

4

10

6

SHAW/sin2pi, n=20, σ=10−8 − TSVD, L=I

||Axk−b||2 ||Lx k||2

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20

corner maximal regularization minimal regularization

Shaw example: 40 × 20, δ = 10−8

• G. Rodriguez

New and not so new estimation of a regularization parameter

Difficult situations: clusters of points

10

−3

10

−2

10

−1

10

1

10

2

residual norm || A x − b ||2 solution (semi)norm || L x ||2 Discrete L−curve, corner at 14

Heat(1) example: 100 × 50, δ = 10−2

• G. Rodriguez

New and not so new estimation of a regularization parameter

Difficult situations: small var. in norms/residuals

10 10

10

10

−5

10 10

5

10

10

10

15

10

20

10 10

10

10 10

2

10

4

10

6

10

8

10

10

10

12

10

14

10

16

Shaw example: 40 × 20, δ = 10−8

• G. Rodriguez

New and not so new estimation of a regularization parameter

Difficult situations: small var. in norms/residuals

10

−0.02

10

−0.01

10 10

−2

10 10

2

10

4

10

6

10

8

10

10

10

12

10

14

10

16

10

18

10 10

0.1

10 10

2

10

4

10

6

10

8

10

10

10

12

10

14

10

16

Shaw example: 40 × 20, δ = 10−8 left: f ∈ ker(L) right: min Kf − g = 1

• G. Rodriguez

New and not so new estimation of a regularization parameter

The adaptive pruning (or corner) algorithm

The key idea in [Hansen, Jensen, R 2006]: the corner can easily be found in a L-curve if we remove the right amount of points. If too few points are removed we still maintain unwanted local features, and if too many points are removed the corner will be incorrectly located or may disappear.

1

Stage one.

1

We construct a sequence of pruned L-curves (selecting the largest line segments).

2

For each curve, we select candidate corners by two algorithms:

the first is based on the local behaviour of the curve (angles); the second is based on the global behaviour of the curve.

2

Stage two.

1

We select a candidate from the list, so that

the curve is convex at that point; the point is the last one before the residuals stagnate.

• G. Rodriguez

New and not so new estimation of a regularization parameter

The corner algorithm at work

10

−3

10

−2

10

−1

10 10

1

10

2

10

3

Heat(1) example: 100 × 50, δ = 10−2

• G. Rodriguez

New and not so new estimation of a regularization parameter

The corner algorithm at work

10

−3

10

−2

10

−1

10 10

1

10

2

10

3

Heat(1) example: 100 × 50, δ = 10−2

• G. Rodriguez

New and not so new estimation of a regularization parameter

The corner algorithm at work

10

−3

10

−2

10

−1

10 10

1

10

2

10

3

Heat(1) example: 100 × 50, δ = 10−2

• G. Rodriguez

New and not so new estimation of a regularization parameter

The corner algorithm at work

10

−3

10

−2

10

−1

10 10

1

10

2

10

3

Heat(1) example: 100 × 50, δ = 10−2

• G. Rodriguez

New and not so new estimation of a regularization parameter

The corner algorithm at work

10

−3

10

−2

10

−1

10

1

10

2

residual norm || A x − b ||2 solution (semi)norm || L x ||2 Discrete L−curve, corner at 14

Heat(1) example: 100 × 50, δ = 10−2, kbest = 12

• G. Rodriguez

New and not so new estimation of a regularization parameter

Error estimates based on interpolation/extrapolation

Let f† be the normal solution to min

f

Kf − gδ2. If f is an approximate solution and r = gδ − Kf, then f† − f2 ≃ η2

ν = dν−1

d5−2ν

1

dν−3

2

ν ∈ R where d0 = r2, d1 = K Tr2, d2 = KK Tr2. We select k∗ = arg mink ην(k). Two examples: η2 = r·K T r

KK T r , η3 = r2 K T r (Auchmuty).

[Brezinski, R, Seatzu 2008, 2009], [Reichel, R, Seatzu 2009].

• G. Rodriguez

New and not so new estimation of a regularization parameter

Error estimates based on interpolation/extrapolation

To obtain η2

ν = dν−1

d5−2ν

1

dν−3

2

, ν ∈ R, we interpolate

d0 = r2 =

m

• i=1

α2

i

≃ α2 d1 = K Tr2 =

r

• i=1

σ2

i α2 i ≃ σ2α2

d2 = KK Tr2 =

r

• i=1

σ4

i α2 i ≃ σ4α2

and extrapolate at e2 =

r

• i=1

σ−2

i

α2

i ≃ σ−2α2.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Error estimates based on interpolation/extrapolation

10

−10

10

−8

10

−6

10

−4

10

−2

10 10

−6

10

−4

10

−2

10 10

2

10

4

10

6

10

8

error η2 η3 η4

TSVD/Shaw example: 40 × 20, δ = 10−8

• G. Rodriguez

New and not so new estimation of a regularization parameter

Comparison of best performance

Parameter choice rules for TSVD; factors 10 and 100. Square Overdetermined systems Method systems ξ = 0 ξ = 1 ξ = 10 COSE 0% ( 0%) 0% (0%) 0% ( 0%) 0% ( 0%) L-corner 7% ( 0%) 4% (0%) 5% ( 1%) 27% (11%) Res L-curve 3% ( 1%) 1% (0%) 18% ( 5%) 34% (18%) Cond L-curve 3% ( 2%) 7% (5%) 32% (18%) 38% (23%) Regi´ nska 4% ( 0%) 2% (0%) 13% ( 0%) 24% ( 1%) ResReg 3% ( 0%) 2% (0%) 13% ( 0%) 24% ( 1%) Quasiopt 9% ( 1%) 6% (1%) 6% ( 1%) 5% ( 1%) GCV 20% (16%) 8% (4%) 8% ( 0%) 11% ( 0%) Extrapolation 5% ( 1%) 9% (0%) 35% ( 0%) 25% ( 1%) Discrepancy 0% ( 0%) 0% (0%) 31% (31%) 41% (40%)

• G. Rodriguez

New and not so new estimation of a regularization parameter

Comparison of best performance

Parameter choice rules for TSVD; factors 2 and 5. Square Overdetermined systems Method systems ξ = 0 ξ = 1 ξ = 10 COSE 6% ( 0%) 7% ( 1%) 7% ( 1%) 8% ( 1%) L-corner 24% (12%) 23% (10%) 36% (10%) 63% (40%) Res L-curve 26% ( 6%) 35% ( 1%) 61% (32%) 71% (47%) Cond L-curve 40% ( 5%) 50% (11%) 69% (46%) 67% (49%) Regi´ nska 25% (10%) 30% (11%) 58% (35%) 76% (51%) ResReg 24% ( 8%) 30% (11%) 58% (35%) 76% (51%) Quasiopt 31% (14%) 31% (14%) 31% (14%) 28% (12%) GCV 29% (22%) 17% (10%) 40% (18%) 67% (36%) Extrapolation 62% (14%) 76% (34%) 83% (62%) 76% (52%) Discrepancy 17% ( 1%) 38% ( 4%) 62% (37%) 71% (51%)

• G. Rodriguez

New and not so new estimation of a regularization parameter

The idea behind the COSE approach

“All happy families are alike; each unhappy family is unhappy in its own way.” Lev Tolstoy, Anna Karenina

5 10 15 20 25 30 35 40 −0.5 0.5 1 1.5 2 2.5 3 exact k=2 k=7 k=10 5 10 15 20 25 30 35 40 −0.5 0.5 1 1.5 2 2.5 exact µ=5e−1 µ=7e−3 µ=2e−4

Two different regularization methods agree when they are able to compute a good approximation of the solution, but produce quite different results when their outcome is inaccurate.

• G. Rodriguez

New and not so new estimation of a regularization parameter

The name of the game

COSE algorithm [Hochstenbach, Reichel, R 2015], [Park, Reichel, R, Yu 2018] Actually, it is an acronym for comparison of solutions estimator We compute the solution to a LDIP by two regularization methods, and determine the parameter in order to make the two solutions as close as possible.

• G. Rodriguez

New and not so new estimation of a regularization parameter

COSE algorithm

1 For each TSVD solution fk we compute the residuals

ρk = g − Kfk, k = 1, 2, . . .

2 We compute the Tikhonov solution fλk such that

g − Kfλk = ρk. In fact, it can be proved that Theorem ϕ(λ) := g − Kfλ2 = gT(λ−2KK T + I)−2g, with K Tg = 0, is a strictly increasing function of λ > 0; moreover ϕ(λ) = γ has a unique solution λ > 0, for any g02 < γ < g2.

• G. Rodriguez

New and not so new estimation of a regularization parameter

COSE algorithm

Let kmin correspond to the first local minimum of the difference δk = fλk − fk, k = 1, 2, . . . We choose λ = λkmin as the Tikhonov regularization parameter. We can use ρkmin = g − Kfλkmin as an estimate for the norm

• f the noise in b, and use this information to apply the

discrepancy principle to choose the regularization parameter. It turns out that the solution fλkmin is an excellent approximation of the solution for many test problems. The computations of the algorithm are simple and inexpensive when the SVD of K, or GSVD of (K, L), is available.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Notation for the problem in general form

Assuming N(K) ∩ N(L) = {0}, we consider the two problems    min

f

Lf2

• subj. to min

f

Kf − g min

f {Kf − g2 + λ2Lf2}

Let the GSVD of (K, L) be (m ≥ n ≥ p) K = U Σ In−p

• Z −1,

L = V

• M
• Z −1,

where U and V have orthonormal columns, Σ = diag[σ1, σ2, . . . , σp], M = diag[µ1, µ2, . . . , µp] and σ2

i + µ2 i = 1 for i = 1, 2, . . . , p.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Expression of the two solutions

TGSVD: fk =

p

• i=p−k+1

uT

i g

σi zi +

n

• i=p+1

(uT

i g) zi

Tikhonov: fλ =

p

• i=1

σi uT

i g

σ2

i + λ2µ2 i

zi +

n

• i=p+1

(uT

i g) zi.

We compute ρk = Kfk − UUTg, k = 1, 2, . . . , ℓ, where UUT is the orthogonal projection in the range of K, and seek for the solution fλk which produces a residual identical to ρk.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Determining the Tikhonov parameter λ

The zero λk of Kfλ − UUTg2 = ρ2

k

must be found for k = 1, 2, . . . It can be determined by a Newton method, i.e., by computing the limit τ = 1/λ2 of the sequence τq+1 = τq + 1 2  

p

• j=1

µ4

j (uT j g)2

(σ2

j τq + µ2 j )2 − ρ2 k

  ·  

p

• j=1

σ2

j µ4 j (uT j g)2

(σ2

j τq + µ2 j )3

 

−1

. The value of the regularization parameter is given by λk = τ −1/2. At this point we search for the minimum of δk = fk − fλk.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Numerical results

10 20 30 40 0.5 1 k=1,

1=7.19e-01

10 20 30 40 0.5 1 k=3,

3=2.71e-01

10 20 30 40 0.5 1 k=5,

5=1.08e-01

10 20 30 40 0.5 1 k=7,

7=2.06e-01

10 20 30 40 0.5 1 k=9,

9=4.24e-01

10 20 30 40 0.5 1 k=11,

11 =3.42e+00

Gravity, 40 × 40, L = D1, δ = 10−2: TGSVD, Tikhonov

• G. Rodriguez

New and not so new estimation of a regularization parameter

Numerical results

5 10 15 20 10 -1 10 0 10 1 10 2 10 3 10 4

Error

k

Gravity, same example: error and δk’s, both minimal at k = 5.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Large scale computations

Application of r steps of Golub–Kahan bidiagonalization to K with initial vector b yields the decompositions KWr = Wr+1Br+1,r, K T Wr = WrBT

r,r,

where Wr+1 ∈ Rm×(r+1) and Wr ∈ Rn×r have orthonormal columns, w1 = g/g, Br+1,r :=        ρ1 σ2 ρ2 ... ... σr ρr σr+1        ∈ R(r+1)×r. and Br,r is the leading r × r submatrix of Br+1,r.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Large scale computations

The columns of Wr span the Krylov subspace Kr(K TK, K Tb). Following [Hochstenbach, Reichel 2010], we project the Tikhonov problem in Kr, by min

y∈Rr{KWry − g2 + λ2LWry2},

and introduce the (compact) QR factorization LWr = QrRr,

• btaining the small-scale problem

min

y∈Rr{Br+1,ry − ge1 2 + λ2Rry2}.

The GSVD approach is then applied to min

y Br+1,ry − ge1 2.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Large scale computations: an alternative approach

In case the structure of the null space of L is known, an alternative, more effective, technique can be adopted [Morigi, Reichel, Sgallari 2006], [Baglama, Reichel 2007]. Let the orthonormal columns of the matrix ˘ Ws ∈ Rn×s span N(L) and introduce the QR factorization, K ˘ Ws = ˘ Qs ˘ Rs and the orthogonal projectors P ˘

Ws = ˘

Ws ˘ W T

s , P⊥ ˘ Ws = I− ˘

Ws ˘ W T

s , P ˘ Qs = ˘

Qs ˘ QT

s , P⊥ ˘ Qs = I− ˘

Qs ˘ QT

s .

• G. Rodriguez

New and not so new estimation of a regularization parameter

Large scale computations: an alternative approach

Then, we can write the residual as Kf − g2 = P ˘

QsKP ˘ Wsf − (P ˘ Qsg − P ˘ QsKP⊥ ˘ Wsf)2

+ P⊥

˘ QsKP⊥ ˘ Wsf − P⊥ ˘ Qsg2.

The two components of the solution can be computed by solving min

f∈Rn{P⊥ ˘ QsAP⊥ ˘ Wsf − P⊥ ˘ Qsg2 + λ2LP⊥ ˘ Wsf2}

and ˘ Rsy − ( ˘ QT

s g − ˘

QT

s AP⊥ ˘ Wsf) = 0,

where y = ˘ W T

s f. Then,

f = P ˘

Wsf + P⊥ ˘ Wsf.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Comparison of large scale computations

100 200 300 400 500

• 0.05

0.05 0.1 0.15 0.2 0.25 0.3 0.35

solution Algorithm 2 Algorithm 3

100 200 300 400 500

• 1
• 0.5

0.5 1 1.5 2 2.5

solution Algorithm 2 Algorithm 3

Phillips and Shaw examples: 500 × 500, δ = 10−2, L = D2.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Numerical results: image deblurring

• bserved image

best error, k=22 COSE, k=34 L−corner, k=52

Figure: barbara test image, size 256 (n2 = 65536), noise 10−2.

• G. Rodriguez

New and not so new estimation of a regularization parameter

Numerical results: image deblurring

• bserved image

best error, k=22 COSE, k=33 L−corner, k=52

Figure: spacestat test image, size 256 (n2 = 65536), noise 10−2.

• G. Rodriguez

New and not so new estimation of a regularization parameter

1

P.C. Hansen, T.K. Jensen, and G. Rodriguez. An adaptive pruning algorithm for the discrete L-curve criterion.

• J. Comput. Appl. Math., 198(2):483-492, 2007.

2

• C. Brezinski, G. Rodriguez, and S. Seatzu. Error estimates for linear systems

with applications to regularization.

• Numer. Algorithms, 49(1-4):85-104, 2008.

3

• C. Brezinski, G. Rodriguez, and S. Seatzu. Error estimates for the regularization
• f least squares problems.
• Numer. Algorithms, 51(1):61-76, 2009.

4

• L. Reichel, G. Rodriguez, and S. Seatzu. Error estimates for large-scale ill-posed

problems.

• Numer. Algorithms, 51(3):341-361, 2009.

5

• L. Reichel and G. Rodriguez. Old and new parameter choice rules for discrete

ill-posed problems.

• Numer. Algorithms, 63(1):65-87, 2013.

6

• M. E. Hochstenbach, L. Reichel, and G. Rodriguez. Regularization parameter

determination for discrete ill-posed problems. J. Comput. Appl. Math., 273:132-149, 2015.

7

• C. Fenu, L. Reichel, and G. Rodriguez. GCV for Tikhonov regularization via

global Golub-Kahan decomposition.

• Numer. Linear Algebra Appl., 23:467-484,

2016.

8

• C. Fenu, L. Reichel, G. Rodriguez, and H. Sadok. GCV for Tikhonov

regularization by partial SVD. BIT, 57:1019-1039, 2017.

9

• Y. Park, L. Reichel, G. Rodriguez, and X. Yu. Parameter determination for

Tikhonov regularization problems in general form. J. Comput. Appl. Math., 343:12-25, 2018.

• G. Rodriguez

New and not so new estimation of a regularization parameter