next previous

3. Optimization without control of robustness

Let tex2html_wrap_inline3528 be any subset of tex2html_wrap_inline3530, say that generated by an aborted matching pursuit process; tex2html_wrap_inline3532 has m  elements. Let us now consider the problem of minimizing tex2html_wrap_inline3536 on the space  E generated by the tex2html_wrap_inline3540, k  spanning

By definition, E is the range of the operator:

In the case where tex2html_wrap_inline3546 is equipped with its standard scalar product, the adjoint of S is explicitly defined by the relationship:

Indeed, for any tex2html_wrap_inline3552, we have from Eq. (36):

In what follows, S  is not necessarily a one-to-one map from tex2html_wrap_inline3556 onto  E : the vectors tex2html_wrap_inline3560 lying in tex2html_wrap_inline3562 are not necessarily linearly independent.

Let tex2html_wrap_inline3564 now be a vector minimizing on tex2html_wrap_inline3566 the quantity tex2html_wrap3748. Then, the vector tex2html_wrap_inline3570 minimizes  q  on  E . From Eq. (2), the vectors tex2html_wrap_inline3576 in question are such that

These vectors are therefore the solutions of the normal equation
(the least-squares solutions of the equation tex2html_wrap3750).

In most cases encountered in image reconstruction, the conjugate-gradients method is the best suited technique for solving Eq. (38). The version of this method presented below provides tex2html_wrap3752.


#&# Step 0: &Set tex2html_wrap_inline3582 (for example) and n = 0 ; &choose a natural starting point tex2html_wrap_inline3586 in E ; &compute tex2html_wrap_inline3590, &tex2html_wrap_inline3592; &compute tex2html_wrap_inline3594(for all tex2html_wrap_inline3596; &set tex2html_wrap_inline3598(for all tex2html_wrap_inline3600. 10ptStep 1: &Compute &tex2html_wrap_inline3602, 5pt&tex2html_wrap_inline3604, 5pt&tex2html_wrap_inline3606(for all tex2html_wrap_inline3608, 5pt&tex2html_wrap_inline3610, 5pt&tex2html_wrap_inline3612, 5pt&tex2html_wrap_inline3614, 5pt&tex2html_wrap_inline3616, 5pt&tex2html_wrap_inline3618(for all tex2html_wrap_inline3620; 6pt&if tex2html_wrap_inline3622, &termination. 6pt&Compute &tex2html_wrap_inline3624, 5pt&tex2html_wrap_inline3626(for all tex2html_wrap_inline3628; 5pt&increment n and loop to step 1. 101

Throughout this algorithm, tex2html_wrap_inline3632 is the residue of the equation tex2html_wrap3754 for tex2html_wrap3756. Likewise, tex2html_wrap_inline3638 is the value of tex2html_wrap_inline3640 at the same iterate:

The iteration in tex2html_wrap_inline3642 results from the identity:

tex2html_wrap_inline3644 The sequence of vectors tex2html_wrap_inline3646 converges towards a solution of the problem with all the remarkable properties of the tex2html_wrap3758 method (see Lannes et al. 1987b). In practice, E is chosen so that tex2html_wrap_inline3650 is a tex2html_wrap3760 map. The uniqueness of the solution can easily be verified by modifying the starting point of the algorithm. The stopping criterion is based on the fact that the final residue must be practically orthogonal to all the tex2html_wrap_inline3652's (Eq. (37)); the tex2html_wrap3762 of the angle between the vectors tex2html_wrap_inline3654 and tex2html_wrap_inline3656 is equal to tex2html_wrap3764. Here, as tex2html_wrap_inline3660 is endowed with its standard scalar product, this algorithm cannot provide the condition number of tex2html_wrap_inline3662. (The transposition of what is presented in Sect. 4.2 would give the "generalized condition number" of  AS ). We therefore recommend to use algorithm 1 only when tex2html_wrap_inline3666 is approximately known.

REMARK| 3. Let us consider the special case where  A  is the identity operator on tex2html_wrap_inline3670 (which then coincides with tex2html_wrap_inline3672). The problem is then to find tex2html_wrap_inline3674, the projection of tex2html_wrap3766 onto  E . Note that tex2html_wrap_inline3680 is then equal to unity. In this case, algorithm 1 reduces to


#&# Step 0: &Set tex2html_wrap_inline3682 (for example) and n = 0 ; &set tex2html_wrap_inline3686 and tex2html_wrap_inline3688; &compute &tex2html_wrap_inline3690, &tex2html_wrap_inline3692(for all tex2html_wrap_inline3694; &set tex2html_wrap_inline3696(for all tex2html_wrap_inline3698. 10ptStep 1: &Compute &tex2html_wrap_inline3700, 5pt&tex2html_wrap_inline3702(for all tex2html_wrap_inline3704, 5pt&tex2html_wrap_inline3706, 5pt&tex2html_wrap_inline3708, 5pt&tex2html_wrap_inline3710, 5pt&tex2html_wrap_inline3712, 5pt&tex2html_wrap_inline3714(for all tex2html_wrap_inline3716); 6pt&if tex2html_wrap_inline3718, &termination. 6pt&Compute &tex2html_wrap_inline3720, 5pt&tex2html_wrap_inline3722(for all tex2html_wrap_inline3724) ; 5pt&increment n and loop to step 1. 101

This algorithm converges towards the projection of tex2html_wrap_inline3728 onto  E with all the properties of the conjugate-gradients method.

In principle, the projection operation can also be performed by using the matching pursuit iteration (25). In this case, on setting tex2html_wrap_inline3732 equal to its optimal value, this iteration reduces to

The residues tex2html_wrap_inline3734 are then obtained via iteration (27):
and likewise for tex2html_wrap_inline3736 (cf. Eq. (28)):

At each iteration, it is then natural to choose k so that tex2html_wrap_inline3740. In the general case where the tex2html_wrap_inline3742's (tex2html_wrap_inline3744) do not form an orthogonal set, the conjugate-gradients algorithm is of course preferable.

next previous

Copyright by the European Southern Observatory (ESO)