Abstract

Finding the closest Toeplitz matrix^* * This work has been supported by UNS grant 24/L031.

Maria Gabriela Eberle; Maria Cristina Maciel

Department of Mathematics, Southern National University Av. Alem 1253, 8000 Bahía Blanca, Argentina

E-mail: geberle@criba.edu.ar, immaciel@criba.edu.ar

ABSTRACT

The constrained least-squares n × n-matrix problem where the feasibility set is the subspace of the Toeplitz matrices is analyzed. The general, the upper and lower triangular cases are solved by making use of the singular value decomposition. For the symmetric case, an algorithm based on the alternate projection method is proposed. The implementation does not require the calculation of the eigenvalue of a matrix and still guarantees convergence. Encouraging preliminary results are discussed.

Mathematical subject classification: 65F20, 65F30, 65H10.

Key words: constrained least squares, Toeplitz matrix, singular value decomposition, alternate projection method.

1 Introduction

In this work the constrained least-squares Toeplitz matrix is analyzed. Let us begin considering the general constrained least squares matrix problem

where A, B Î ^mxn, with m > n and Ì ^nxn has a particular pattern. Throughout this paper A

In the last two decades an extensive research activity has been done on this kind of problems which arise in many areas. For instance, in statistics, the problem of finding the nearest symmetric positive definite patterned matrix to a sample covariance matrix is a classical example, [11]. The symmetric positive semidefinite case has analyzed by Higham [9].

In the investigation of elastic structures, the determination of the strain matrix which relates forces f_i and displacement d_i according to Xf_i = d_i is stated as (1). The orthogonal case has been studied by Schonemann [13] and the Higham [8] and the symmetric case by Higham [10]. These problems are known in the literature as the orthogonal and symmetric Procrustes problems respectively. When is the subspace of persymmetric matrices (that is, it is symmetric about the NE-SW diagonal) the resulting problem has been studied by Eberle [4, 5].

Now in this work we are interested in solving the constrained least squares Toeplitz matrix problem

where Ì ^nxn is the subspace of Toeplitz matrices.

Let us recall that a matrix T = (t_ij) Î ^nxn is said to be a Toeplitz matrix if it is constant along its diagonals, i.e there exist scalars , r_–n+1, ¼, r₀,¼, r_n–1, such that t_ij = r_j–i, for all i,j:

The rest of the paper is organized as follows. In section 2, the general Toeplitz problem is studied. The special cases of triangular and symmetric Toeplitz problems are analyzed in section 3 and 4 respectively. The numerical results are presented in section 5 and our conclusions and final remarks in section 6.

2 The general problem

In this section the feasibility set, , of problem (2) is characterized.

Given the matrices A, B Î ^mxn, the subspace can be represented by

where a_p Î and G_p Î ^nxn are matrices defined as follows

We must observe that the matrices G_p are characterized by the following fact: for each entry ij, there exists an unique index p, –(n–1) < p < (n–1), such that (G_p)_ij = 1. Clearly they form a basis for the subspace .

2.1 The transformed problem

We will transform the problem (2) in a simpler one. To do this transformation let A be the matrix have the following singular value decomposition

where P Î ^mxm and Q Î ^nxn, are orthogonal matrices and S = diag(s₁,¼, s_n), s₁> ¼ > s_n > 0.

Since the Frobenius norm remains invariant under orthogonal transformations it is obtained

with

then the problem 2 is equivalent to

where ¢ is the following subspace:

2.2 Characterization of the solution on ¢

We introduce the notation P(A) meaning the projection of the matrix A on the set . The following theorem characterizes the projection on the subspace of matrices ¢.

Theorem 2.1. If C₁Î ^nxn, then the unique solution of the problem (4) is

for all –(n–1) < l < (n–1).

Proof. The objective function f : ^2n–1

Since f is twice continuously differentiable, we compute (a) for all p such that –(n–1) < p < (n–1)

The factor in the right hand side of (5) can be written as

and (5) becomes in

Having in mind that

is the inner product between the i^th row of (SQ^T) and the j^th column of a_lG_l, we obtain

Similarly, the element (SQ^TG_p)_ij, is the inner product between the i^th row of SQ^T and the column j of G_p. Then:

(SQ^TG_p)_ij = s_iQ_(j–p)i .

Because of the sparse structure of G_p the entries (G_p)_ij = 1 if j > p.

Therefore (6) can be written as

Defining

with

we have that the sum of all the components of v is the inner product between two columns of the orthogonal matrix Q and then it is zero when k ¹ p. Therefore for k ¹ p we have:

then v^Tq = 0 for all k ¹ p and

Then (6) becomes

for all p = 1,¼, n. Therefore, from the first order necessary condition

we obtain

We observe that the denominator of (7) is non zero, because if it were zero for any i, j, it would be

Now, we need to verify that is a minimizer of f. We compute

and

It implies that the Hessian matrix Ѳf() is positive define and therefore is the minimizer of f.

3 The triangular Toeplitz problem

In this section we consider the case when the feasibility is the subspace of triangular Toeplitz matrices. We begin considering the upper triangular case.

3.1 The upper triangular Toeplitz problem

Let the minimization problem be

where A, B Î ^mxn and is the subspace of upper triangular Toeplitz matrices. If X Î , then

The subspace , can be written as:

where a_pÎ and the matrices G_p Î ^nxn are defined for all 1 < i, j, p < n as follow

We remark that the matrices G_p are characterized by the property that for each entry ij, there exists an unique p, 1 < p < n such that (G_p)_ij = 1.

The set , is a basis for the subspace of the upper triangular Toeplitz matrices of ^nxn. Similarly to the general Toeplitz case, the singular value decomposition is used in order to obtain a solution of an equivalent problem

with

P Î ^mxm y Q Î ^nxn, are orthogonal matrices, S = diag(s₁,¼, s_n) with s₁> ¼ > s_n> 0, and the subspace

3.1.1 Characterization of the solution on ¢

The following theorem, similar to Theorem 2.1 for the dense case characterizes the projection on ¢.

Theorem 3.1. Let be C₁Î ^nxn, then the unique solution of the problem

is given by

for all 1< l < n.

Proof. Similar to the proof of theorem 2.1 having in mind that the first derivative is

The element (SQ^T((a_lG_l))_ij is the inner product between the i^th row of (SQ^T) and the j^th column of (a_lG_l), then it results

Analogously, the element ( SQ^TG_p)_ij, is the inner between the row i of SQ^T and the column j of G_p, then

3.2 The lower triangular Toeplitz problem

Let us consider the problem

where A, B Î ^mxn and is the subspace of the lower triangular Toeplitz matrices. If X Î , then

The subspace , can be expressed as

where b_pÎ and the matrices G_p Î ^nxn are defined as

for all 1 < i, j, p < n. The set , is a basis for the subspace . Again, applying the singular value decomposition of A the problem can be transformed in

with

where P Î ^mxm y Q Î ^nxn, are the orthogonal matrices, S = diag(s₁,¼, s_n), s₁> ¼ > s_n > 0, and the subspace

As in the general and upper triangular Toeplitz cases the projection is characterized on ¢.

Theorem 3.2. If C₁Î ^nxn, then the unique solution of the problem

is given by

for all 1 < l < n.

Proof. Similar to the proof of Teorema 2.1.

4 The symmetric Toeplitz problem

This section is concerned with the problem

where A, B Î ^mxn , is the subspace of the Toeplitz matrices and is the subspace of the symmetric matrices of ^nxn. Clearly, the solution is a matrix which is symmetric not only with respect to the main diagonal but also to the SW-NE diagonal. Since the feasibility region is the intersection of two subspace and because we know how to solve the persymmetric Procrustes problem [5] and the general Toeplitz problem, see section 2, we develop an algorithm based on the alternating projection method [3, 1, 6, 7]. Therefore the solution of (16) is obtained projecting alternating on the subspaces y .

Given the problem

according with Escalante and Raydán [7] it can be transformed, via the singular value decomposition of A in the following equivalent problem

where

and the feasibility region

¢ = {Z Î ^nxn : Z = SY, Y = Y^T}.

Hence if the solution of the last problem is

where is obtained by the methods proposed by Higham for the symmetric case [10], the solution of the first problem is

On the other hand, we know how to solve the general Toeplitz problem, see section 2. Combining both results, to solve the problem 16 is equivalent to solve

where

Remark. The set ¢¢, comes from the fact that the symmetric and Toeplitz problems have the same objective function.

The alternating projection algorithm for the symmetric Toeplitz case (19) is as follows:

Algorithm 4.1. Given A, B Î ^mxn. X₀ = C₁Î ^nxn, where C₁ is obtained via the singular value decomposition of A, for i = 0, 1, 2, ¼, until convergence repeat

The algorithm terminates when two consecutive projection on the same subspace are close enough. The convergence is guaranteed by the general theory established by von Neumann [12] for the alternating projection algorithm.

5 Numerical experience

In this section we show some examples. The algorithms have been implemented in Fortran, by using FORTRAN POWER STATION (1993) in an environment PC with a processor Pentium(R2) Intelmmx(TM) Technology, 31.0 MB de RAM. and 32 bits of virtual memory. The LINPACK library dsvdc subroutine have been used to obtain the singular value decomposition.

Example 5.1. The general problem (2) has been solved with the following matrices

The solution is

We project on the sets y .

a) The projection on is

b) The projection on is

Example 5.2. We illustrate the behavior of the algorithm 4.1 for the symmetric case carrying out the following example

The solution is

The algorithm terminates when two consecutive projections on the subspace of the Toeplitz matrices is less or equal to a fixed tolerance. We set = 10^-6 for this experiment.

6 Concluding remarks

We have studied the Procrustes Toeplitz problem, analyzing the general, the triangular and the symmetric cases. In the last case an algorithm based on the alternating projection method is proposed. We illustrate our results with some examples.

The Toeplitz systems Tx = b where T Î ^nxn is a Toeplitz matrix, are present in many disciplines (signal and image processing, times series analyzes, queuing networks, partial differential equations problems and control problems). When these systems are solved via preconditioned conjugate gradient methods, the search of adequate preconditioners is the main issue. Among the possible preconditioners for a Toeplitz matrix, Chan proposes the preconditioner c(T) defined to be the minimizer of C – T

7 Acknowledgments

The authors wish to thank two anonymous referees for helpful suggestions.

#537/01.

Publication Dates