Least Squares Problems with Absolute Quadratic Constraints

Abstract

This paper analyzes linear least squares problems with absolute quadratic constraints. We develop a generalized theory following Bookstein's conic-fitting and Fitzgibbon's direct ellipse-specific fitting. Under simple preconditions, it can be shown that a minimum always exists and can be determined by a generalized eigenvalue problem. This problem is numerically reduced to an eigenvalue problem by multiplications of Givens' rotations. Finally, four applications of this approach are presented.

1. Introduction

The least squares methods cover a wide range of applications in signal processing and system identification [1–5]. Many technical applications need robust and fast algorithms for fitting ellipses to given points in the plane. In the past, effective methods were Bookstein's conic-fitting or Fitzgibbon's direct ellipse-specific fitting, where an algebaic distance with a quadratic constraint is minimized [6, 7]. In this paper, we develop an extended theory of minimization of least squares with a quadratic constraint based on the ideas of Bookstein and Fitzgibbon. Thereby, we show the existence of a minimal solution and present the uniqueness regarding to the smallest positive generalized eigenvalue. So, arbitrary conic fitting problems with quadratic constraints are possible.

Let 𝐴∈ℝ^𝑛×𝑚 be matrix with 𝑛≥𝑚≥2, 𝐶∈ℝ^𝑚×𝑚 be symmetric matrix, and 𝑑∈ℝ a real value. We consider the problem of finding a vector 𝑥∈ℝ^𝑚 which minimizes the function 𝐹∶ℝ^𝑚→ℝ defined by𝐹(𝑥)∶=‖𝐴𝑥‖²₂subjectto𝑥^𝑡𝐶𝑥=𝑑. (1.1) The side condition 𝑥^𝑡𝐶𝑥=𝑑 introduces an absolute quadratic constraint. The problem (1.1) is not a special case of Gander's optimization as presented in [8], because in our case 𝐶 is a real symmetric matrix in contrast to the approach of Gander, where the side condition considers real symmetric matrices 𝐶^𝑡𝐶 which are positive definite. For our considerations, we require the following three assumptions

Assumption 1.1. By replacing 𝐶 by (−𝐶), we consider 𝑑≥0. For 𝑑=0, the trivial solution 𝑥=0∈ℝ^𝑚 fulfills (1.1). Therefore, we demand 𝑑>0.

Assumption 1.2. The set 𝑁∶={𝑥∈ℝ^𝑚,𝑥^𝑡𝐶𝑥=𝑑} is not empty, that is, the matrix 𝐶 has not been less than one positive eigenvalue. If we assume that 𝐶 has only nonpositive eigenvalues, it would be negative semidefinite and 𝑥^𝑡𝐶𝑥≤0 (1.2) holds. With 𝑑>0 it follows that the set 𝑁 would be empty in this case.

Assumption 1.3. In the following, we set 𝑆=𝐴^𝑡𝐴 and assume that 𝑆 is regular. 𝑆 is sometimes called scatter matrix.

In the following two sections, we introduce the theoretical basics of this optimization. The main result is the solution of a generalized eigenvalue problem. Afterwards, we reduce this system numerically to an eigenvalue problem. In Section 5, we present four typical applications for conic fitting problems with quadratic constraints. These approximations contain the ellipse fitting of Fitzgibbon, the hyperbola fitting of O'Leary, the conic fitting of Bookstein, and an optical application of shrinked aspheres [6, 7, 9, 10].

2. Existence of a Minimal Solution

Theorem 2.1. If 𝑆∈ℝ^𝑚×𝑚 is regular, then there exists a global minimum to the problem (1.1).

Proof. The real regular matrix 𝑆=𝐴^𝑡𝐴 is symmetric and positive definite. Therefore, a Cholesky decomposition 𝑆=𝑅^𝑡𝑅 exists with a regular upper triangular matrix 𝑅∈ℝ^𝑚×𝑚. In (1.1), we are looking for a solution 𝑥∈ℝ^𝑚 minimizing F(𝑥)=‖𝐴𝑥‖²=𝑥^𝑡𝑆𝑥=𝑥^𝑡𝑅^𝑡𝑅𝑥subjectto𝑥^𝑡𝐶𝑥=𝑑. (2.1) With 𝑅 regular we substitute 𝑥 by 𝑅⁻¹𝑦 for 𝑦∈ℝ^𝑚. Thus, we obtain an equivalent problem to (2.1), where we want to find a vector 𝑦∈ℝ^𝑚, minimizing 𝐹𝑦=𝑦^𝑡𝑦=𝑦²subjectto𝑦^𝑡𝑅^−1𝑡𝐶𝑅⁻¹𝑦=𝑑. (2.2) Now, we define 𝐺∶ℝ^𝑚→ℝ with 𝐺(𝑦)=𝑦^𝑡(𝑅⁻¹)^𝑡𝐶𝑅⁻¹𝑦−𝑑 and look for a solution 𝑦 on the zero-set 𝑁_𝐺 of 𝐺 with minimal distance to the point of origin. Let 𝑦₀∈𝑁_𝐺≠∅ and 𝐾^𝑚_𝑟0(0) be the closed sphere of ℝ^𝑚 in 0 with radius 𝑟₀=‖𝑦₀‖. Because of 𝑦₀∈𝑁_𝐺∩𝐾^𝑚_𝑟0(0) and 𝐺 being continuous, the set 𝑁_𝐺𝐾^𝑚_𝑟0(0)=𝑦∈𝐾^𝑚_𝑟0(0),𝐺𝑦=0 (2.3) is nonempty, closed and bounded. Therefore, for the continuous function 𝐹(𝑦)=‖𝑦‖² exists a minimal value 𝑦_𝑀 on 𝑁_𝐺∩𝐾^𝑚_𝑟0(0) with min_{𝑦∈𝑁𝐺∩𝐾𝑚𝑟0(0)}𝐹𝑦≤𝑟²₀. (2.4) For all 𝑦 from 𝑁_𝐺⧵𝐾^𝑚_𝑟0(0), it is 𝐹(𝑦)>𝑟²₀. So, 𝑦_𝑀 is a minimal value of 𝐹 in 𝑁_𝐺. By the equivalence of (2.1), and (2.2) the assumption follows.

3. Generalized Eigenvalue Problem

The minimization problem in (1.1) induces a generalized eigenvalue problem. The following theorem is already proven by Bookstein and Fitzgibbon for the special case of ellipse-fitting [6, 7].

Theorem 3.1. If 𝑥_𝑠 is an extremum of 𝐹(𝑥) subject to 𝑥^𝑡𝐶𝑥=𝑑, then a positive 𝜆₀∈ℝ exists with 𝑆−𝜆₀𝐶𝑥_𝑠=0, (3.1) that is, 𝑥_𝑠 is an eigenvector to the generalized eigenvalue𝜆₀ and 𝐹𝑥_𝑠=𝜆₀⋅𝑑 (3.2) holds.

Proof. Let 𝐺∶ℝ^𝑚→ℝ be a defined as 𝐺(𝑥)∶=𝑑−𝑥^𝑡𝐶𝑥. For 𝐺(𝑥)=0 and 𝑑>0 follows 𝑥≠0. Further, 𝐺 is continuously differentiable with 𝑑𝐺/𝑑𝑥=−2𝐶𝑥≠0 for all 𝑥 of the zero-set of 𝐺. So, if 𝑥_𝑠 is a local extremum of 𝐹(𝑥) subject to 𝐺(𝑥)=0, then it is rank(𝑑𝐺/𝑑𝑥)(𝑥_𝑠)=1. Since 𝐹 is also a continuously differentiable function in ℝ^𝑚 with 𝑚>1, it follows by using a Lagrange multiplier [11]: if 𝑥_𝑠 is a local extremum of 𝐹(𝑥) subject to 𝐺(𝑥)=0, then a 𝜆₀∈ℝ exists, such that the Lagrange function 𝜙∶ℝ^𝑚+1→ℝ given as 𝜙(𝑥,𝜆)=𝐹(𝑥)+𝜆𝐺(𝑥)=𝑥^𝑡𝐴^𝑡𝐴𝑥+𝜆𝑑−𝑥^𝑡𝐶𝑥 (3.3) has a critical point in (𝑥_𝑠,𝜆₀). Therefore, 𝑥_𝑠 fulfills necessarily the equations: grad_𝑥𝜙()()𝑥,𝜆=2𝑆𝑥−2𝜆𝐶𝑥=0,𝑑𝜙𝑥,𝜆𝑑𝜆=𝑥^𝑡𝐶𝑥−𝑑=0. (3.4)(3.5) The first equation describes a generalized eigenvalue problem with (𝑆−𝜆⋅𝐶)𝑥=0. (3.6) With 𝑑>0, 𝑥_𝑠≠0 and 𝑥_𝑠 fulfills (3.6), 𝜆 must be a generalized eigenvalue, and 𝑥_𝑠 is a corresponding eigenvector to 𝜆 of (3.6), so that (𝑆−𝜆𝐶) is a singular matrix. If 𝜆₀ is an eigenvalue and 𝑥₀≠0 a corresponding eigenvector to 𝜆₀ of (3.6), then every vector 𝛼𝑥₀ is also a solution of (3.6) for 𝜆₀. Now, we are looking for 𝛼, such that 𝑥_𝑠=𝛼𝑥₀ satisfies (3.5). For 𝜆₀≠0 and (3.4), it follows 𝑑=𝛼²𝑥^𝑡₀𝐶𝑥₀=𝛼²𝑥^𝑡₀𝑆𝑥₀𝜆₀. (3.7) Because the left side and the numerator are positive, the denominator must also be chosen positive, that is, only positive eigenvalues solve (3.4) and (3.5). By the multiplication with 𝜆₀, 𝑑⋅𝜆₀=𝛼²𝑥^𝑡₀𝑆𝑥₀=𝑥^𝑡_𝑠𝑆𝑥_𝑠=𝐹𝑥_𝑠 (3.8) follows and 𝑥_𝑠=𝛼⋅𝑥₀ fulfills the constraint 𝐺(𝑥_𝑠)=0.

Remark 3.2. Let 𝑥₀ be a generalized eigenvector to a positive eigenvalue 𝜆₀ of problem (3.6). Then 𝑥_1,2=±𝑑𝜆₀𝑥^𝑡₀𝑆𝑥₀𝑥₀ (3.9) are solutions of (3.8).

Lemma 3.3. If 𝑆 is regular and 𝐶 is symmetric, then all eigenvalues of (3.1) are real-valued and different to zero.

Proof. With det(𝑆)≠0,𝜆₀≠0 in (3.1). The Cholesky decomposition 𝑆=𝑅^𝑡𝑅 with a regular upper triangular matrix 𝑅 yields (3.1) 𝑅^𝑡𝑅−𝜆₀𝐶𝑥_𝑠=𝑅^𝑡𝐼−𝜆₀𝑅^𝑡−1𝐶𝑅⁻¹𝑅𝑥_𝑠=0. (3.10) With 𝑅 invertible and the substitution 𝜇₀=𝜆₀⁻¹,𝑦_𝑠=𝑅𝑥_𝑠, we obtain an eigenvalue problem to the matrix (𝑅^𝑡)⁻¹𝐶𝑅⁻¹: 𝜇₀𝐼−𝑅^𝑡−1𝐶𝑅⁻¹𝑦_𝑠=0. (3.11) Furthermore, we have 𝑅^𝑡−1𝐶𝑅^−1𝑡=𝑅^−1𝑡𝐶^𝑡𝑅^𝑡−1𝑡=𝑅^𝑡−1𝐶^𝑡𝑅⁻¹=𝑅^𝑡−1𝐶𝑅⁻¹. (3.12) Therefore, the matrix (𝑅^𝑡)⁻¹𝐶𝑅⁻¹ is symmetric and all eigenvalues 𝜇₀ are real. With 𝜆₀=𝜇₀⁻¹ for 𝜇₀≠0 follows the proposition.

Remark 3.4. Because of 𝑆 regular and 𝜆≠0 in 𝑆−𝜆₀𝐶𝑥_𝑠=0, (3.13) we can consider the equivalent problem with 𝜇₀=𝜆₀⁻¹ instead of (3.11): 𝜇₀𝐼−𝑆⁻¹𝐶𝑥_𝑠=0. (3.14) This system is called inverse eigenvalue problem to (3.1). Here, the eigenspaces to the generalized eigenvalue 𝜆₀ in (3.1) and to the eigenvalue 1/𝜆₀ in (3.14) are identical. Therefore, the generalized eigenvectors in (3.1) are perpendicular.

Definition 3.5. The set of all eigenvalues of a matrix 𝐶 is called spectrum 𝜎(𝐶). We call the set of all eigenvalues to the generalized eigenvalue problem in (3.1) also spectrum and denote 𝜎(𝑆,𝐶). 𝜎⁺(𝑆,𝐶) is defined as the set of all positive values in 𝜎(𝑆,𝐶).

Remark 3.6. In case of 𝑟𝑔(𝐶)<𝑚=𝑟𝑔(𝑆), the inverse problem in (3.14) has eigenvalue 0 with multiplicity 𝑟𝑔(𝑆)−𝑟𝑔(𝐶). Otherwise, for 𝜇≠0 and 𝜇∈𝜎(𝑆⁻¹𝐶) follows 1/𝜇∈𝜎(𝑆,𝐶). Analogously, for 𝜇∈𝜎((𝑅^𝑡)⁻¹𝐶𝑅⁻¹) with 𝜇≠0,1/𝜇∈𝜎(𝑆,𝐶).

The following lemma is a modified result of Fitzgibbon [7].

Lemma 3.7. The signs of the generalized eigenvalues of (3.1) are the same as those of 𝐶.

Proof. With 𝑆 being nonsingular, every generalized eigenvalue 𝜆₀ of (3.1) is not zero. Therefore, it follows for the equivalent problem (3.11) that 𝜇₀=𝜆₀⁻¹ is also a positive eigenvalue to (𝑅⁻¹)^𝑡𝐶𝑅⁻¹, where 𝑅 is an upper triangular matrix to the Cholesky decomposition of 𝑆. With Sylvester's Inertia Law, we know that the signs of eigenvalues of the matrices (𝑅⁻¹)^𝑡𝐶𝑅⁻¹ are the same as those of 𝐶.

For the following proofs, we need the lemma of Lagrange (see, e.g., [12]).

Lemma 3.8 (Lemma of Lagrange). For 𝑀⊂ℝ^𝑛, 𝑓∶𝑀→ℝ, 𝑔=(𝑔₁,…,𝑔_𝑘)∶𝑀→ℝ^𝑘, and 𝑁_𝑔={𝑥∈𝑀,𝑔(𝑥)=0∈ℝ^𝑘}, let 𝜆∈ℝ^𝑘, so that 𝑥_𝑠∈𝑁_𝑔 is a minimal value of the function Φ_𝜆∶𝑀→ℝ with Φ_𝜆(𝑥)=𝑓(𝑥)+^𝑘_𝑖=1𝜆_𝑖𝑔_𝑖(𝑥). (3.15) Then 𝑥_𝑠 is a minimal solution of 𝑓 in 𝑁_𝑔.

Definition 3.9. Let 𝜆^∗ be the smallest positive value of 𝜎⁺(𝑆,𝐶) and 𝑥^∗_𝑠 a corresponding generalized eigenvector to 𝜆^∗ to the constraint 𝑥_𝑠^∗𝑡𝐶𝑥^∗_𝑠=𝑑.

Lemma 3.10. Let 𝑆−𝜆𝐶 be a positive semidefinite matrix for 𝜆∈𝜎⁺(𝑆,𝐶). Then a generalized eigenvector 𝑥_𝑠 corresponding to 𝜆 is a local minimum of (1.1).

Proof. We consider Φ∶ℝ^𝑚→ℝ with Φ(𝑥)=𝑥^𝑡𝑆𝑥+𝜆𝑑−𝑥^𝑡𝐶𝑥. (3.16) With grad_𝑥Φ(𝑥)=2(𝑆𝑥−𝜆𝐶𝑥) holds grad_𝑥Φ(𝑥_𝑠)=0 and since the Hessian matrix (𝑆−𝜆𝐶) of Φ is positive semidefinite, the vector 𝑥_𝑠 is a minimal solution of Φ. Then, 𝑥_𝑠 minimizes also 𝐹(𝑥) subject to 𝑥^𝑡𝐶𝑥=𝑑 by the Lemma 3.8.

Remark 3.11. In fact, in Lemma 3.10, we require only a semidefinite matrix, because 𝑥_𝑠≠0 fulfills 𝑥^𝑡_𝑠(𝑆−𝜆𝐶)𝑥_𝑠=0 in (3.1).

Lemma 3.12. The matrix (𝑆−𝜆^∗𝐶) is positive semidefinite.

Proof. Let 𝜇 be an arbitrary eigenvalue of ((𝜆^∗)⁻¹𝐼−(𝑅^𝑡)⁻¹𝐶𝑅⁻¹), where 𝑅 is the upper triangular matrix of the Cholesky decomposition from 𝑆. With det𝜆^∗−1𝐼−𝑅^𝑡−1𝐶𝑅⁻¹−𝜇𝐼=det𝜆^∗−1−𝜇𝐼−𝑅^𝑡−1𝐶𝑅⁻¹=0, (3.17) it follows that ((𝜆^∗)⁻¹−𝜇) is an eigenvalue of (𝑅^𝑡)⁻¹𝐶𝑅⁻¹. This value is corresponding in (3.11) with the inverse eigenvalue of problem (3.1). Furthermore, it yields 1𝜆^∗1−𝜇≤max𝜆,𝜆∈𝜎⁺(=1𝑆,𝐶)min{𝜆,𝜆∈𝜎⁺(=1𝑆,𝐶)}𝜆^∗. (3.18) So 𝜇≥0 follows, that is, ((𝜆^∗)⁻¹𝐼−(𝑅^𝑡)⁻¹𝐶𝑅⁻¹) is positive semidefinite, and for 𝑦∈ℝ^𝑚 we obtain 𝑦^𝑡𝜆^∗−1𝐼−𝑅^𝑡−1𝐶𝑅⁻¹𝑦≥0. (3.19) By setting 𝑦=𝑅𝑥 and with regular 𝑅, we get 0≤𝑥^𝑡𝜆^∗−1𝑅^𝑡𝑅−𝐶𝑥=𝑥^𝑡𝜆^∗−1𝑆−𝐶𝑥. (3.20) With 𝜆^∗>0, 𝑥^𝑡(𝑆−𝜆^∗𝐶)𝑥≥0 for 𝑥∈ℝ^𝑚, that is, 𝑆−𝜆^∗𝐶 is positive semidefinite.

Theorem 3.13. For the smallest value 𝜆^∗ in 𝜎⁺(𝑆,𝐶) exists a corresponding generalized eigenvector 𝑥^∗_𝑠, which minimizes 𝐹(𝑥) subject to 𝑥^𝑡𝐶𝑥=𝑑.

Proof. The matrix (𝑆−𝜆^∗𝐶) is positive semidefinite by Lemma 3.12, and with Lemma 3.10 it follows that 𝑥^∗_𝑠 is a local minimum of problem (1.1). Furthermore, we know by Theorem 3.1 that if 𝑥_𝑠 is a local extremum of 𝐹(𝑥) subject to 𝑥^𝑡𝐶𝑥=𝑑, then a positive value 𝜆_𝑠 exists with 𝑆−𝜆_𝑠𝐶𝑥_𝑠=0 (3.21) and it is 𝐹(𝑥_𝑠)=𝜆_𝑠𝑑. Because of the existence of a minimum 𝑥_𝐸 in Theorem 2.1 a value 𝜆_𝐸∈𝜎⁺(𝑆,𝐶) exists in problem (3.21) regarding to 𝑥_𝐸. Otherwise, for an arbitrary local minimum 𝑥_𝑠, 𝐹𝑥^∗_𝑠=𝜆^∗𝑑=min𝜆𝑑,𝜆∈𝜎⁺(𝑆,𝐶)≤𝜆_𝑠𝑑=𝐹𝑥_𝑠. (3.22) So, 𝜆^∗=𝜆_𝐸 follows and 𝑥^∗_𝑠 is a minimum of 𝐹(𝑥) subject to 𝑥^𝑡𝐶𝑥=𝑑.

Example 3.14. We minimize 𝐹∶ℝ²→ℝ with 𝐹(𝜉₁,𝜉₂)=𝜉²₁+𝜉²₂ subject to 𝜉²₁−𝜉²₂=1. So, we have 𝑑=1, 𝑆 the identity matrix 𝐼₂∈ℝ^2×2, and 𝐶∈ℝ^2×2 a diagonal matrix with values −1 and 1. Then, we get the following generalized eigenvalue problem: det𝐼₂−𝜆𝐶=(1−𝜆)(−1−𝜆)=0 (3.23) with eigenvalues 1 and −1. Because of Theorem 3.1, we consider a generalized eigenvector (𝛼,0)^𝑡 with 𝛼∈ℝ⧵{0} only for 𝜆=1. Then, (1,0)^𝑡 and (−1,0)^𝑡 are solutions subject to 𝜉²₁−𝜉²₂=1. This result is conform to the geometric interpretation, since we are looking for 𝑥=(𝜉₁,𝜉₂)^𝑡 on the hyperbola 𝜉²₁−𝜉²₂=1 with minimal distance to the origin.

4. Reduction to an Eigenvalue Problem of Dimension 𝑟𝑔(𝐶)

In numerical applications, a generalized eigenvalue problem is mostly reduced to an eigenvalue problem, for example, by multiplication with 𝑆⁻¹. Thus, we obtain the inverse problem (3.14) from (3.1) (see, e.g., [13]). But, 𝑆 may be ill-conditioned, so that a solution of (3.14) may be numerical instable. Therefore, we present another reduction of (3.1).

Many times, 𝐶 is a sparse matrix with 𝑟∶=rank(𝐶)≤rank(𝑆). This symmetric matrix 𝐶 is diagonalizable in 𝐶=𝑃^𝑡𝐷𝑃 with 𝑃 orthogonal and 𝐷 diagonal. Further, we assume that the first 𝑟 diagonal entrees in 𝐷 are different to 0. For the characteristic polynomial in (3.1), it follows𝑝(𝜆)=det(𝑆−𝜆𝐶)=det𝑃𝑆𝑃^𝑡−𝜆𝐷=0. (4.1) The order of 𝑝 is 𝑟. We decompose these matrices in 𝑃𝑆𝑃^𝑡=S₁S₂S^t₂S₃D,𝐷=¹000 (4.2) with 𝑆₁,𝐷₁∈ℝ^𝑟×𝑟, 𝑆₂∈ℝ^{𝑟×(𝑚−𝑟)}, and 𝑆₃∈ℝ^{(𝑚−𝑟)×(𝑚−𝑟)}. Now, we eliminate 𝑆₂ in 𝑃𝑆𝑃^𝑡 by multiplications with Givens rotations 𝐺_𝑘∈ℝ^𝑚×𝑚,𝑘=1,…,𝑙, so that it follows: 𝐺_𝑙⋯𝐺₂𝐺₁𝑃𝑆𝑃^𝑡=Σ₁0Σ₂Σ₃,𝐺_𝑙⋯𝐺₂𝐺₁Δ𝐷=¹0Δ₂0 (4.3) with Σ₁,Δ₁∈ℝ^𝑟×𝑟, Σ₂,Δ₂∈ℝ^{(𝑚−𝑟)×𝑟}, and Σ₃∈ℝ^{(𝑚−𝑟)×(𝑚−𝑟)}. In (4.1), we achieve with orthogonal 𝐺_𝑘,𝑘=1,…,𝑙𝑝(𝜆)=det𝑃𝑆𝑃^𝑡−𝜆𝐷=det𝐺_𝑙⋯𝐺₂𝐺₁𝑃𝑆𝑃^𝑡Σ−𝜆𝐷=det¹−𝜆Δ₁0Σ₂−𝜆Δ₂Σ₃=detΣ₃detΣ₁−𝜆Δ₁=0. (4.4) Because of 𝑝(0)=det(𝑃𝑆𝑃^𝑡)=det(𝑆)≠0 and 𝑝(0)=det(Σ₃)det(Σ₁), the submatrices Σ₁,Σ₃ are regular and the generalized eigenvalues of det(Σ₁−𝜆Δ₁) are different to zero. So, with 𝑦∈ℝ^𝑟Σ₁−𝜆Δ₁𝑦=0 (4.5) can be transformed to an equivalent eigenvalue problem withΔ₁⁻¹Σ₁−𝜆𝐼𝑦=0. (4.6) This system can be solved by finding the matrix 𝑋 with Δ₁𝑋=Σ₁ using the Gaussian elimination and determining the eigenvalues of 𝑋 computing the 𝑄𝑅 algorithm [13]. Because all steps are equivalent, we have 𝜎(Δ₁⁻¹Σ₁)=𝜎(𝑆,𝐶), that is, the eigenvalues of (3.1) and (4.6) are the same.

With Theorem 3.13, we are looking for the smallest value 𝜆^∗∈𝜎⁺(𝑆,𝐶) and a corresponding generalized eigenvector 𝑥^∗_𝑠 to minimize the problem (3.1). So, 0=𝑆−𝜆^∗𝐶𝑥^∗_𝑠=𝑆𝑃^𝑡−𝜆^∗𝐶𝑃^𝑡𝑃𝑥^∗_𝑠 (4.7) yields. By substitution of 𝑦^∗_𝑠 for 𝑃𝑥^∗_𝑠, we obtain0=𝐺_𝑙⋯𝐺₂𝐺₁𝑃𝑆𝑃^𝑡−𝜆^∗𝐷𝑦^∗_𝑠=Σ₁−𝜆^∗Δ₁0Σ₂−𝜆^∗Δ₂Σ₃𝑦^∗_𝑠. (4.8) We decompose 𝑦^∗_𝑠 into the subvectors 𝑦^∗_𝑠,𝑟∈ℝ^𝑟 and 𝑦^∗_{𝑠,𝑚−𝑟}∈ℝ^𝑚−𝑟 with (𝑦^∗_𝑠)^𝑡=(𝑦^∗_𝑠,𝑟|𝑦^∗_{𝑠,𝑚−𝑟})^𝑡. Then, 𝑦^∗_𝑠,𝑟 is a generalized eigenvector for 𝜆^∗ of the problems (4.5) and (4.6).

Let 𝑦^∗_𝑠,𝑟 be an eigenvector to the smallest positive eigenvalue 𝜆^∗ of (4.6). Since Σ₃ is regular, it follows in (4.8) 𝑦^∗_{𝑠,𝑚−𝑟}=Σ₃⁻¹Σ₂−𝜆^∗Δ₂𝑦^∗_𝑠,𝑟 (4.9)

and a generalized eigenvector 𝑥^∗_𝑠 for 𝜆^∗ in (3.1) is given as: 𝑥^∗_𝑠=𝑃^𝑡𝑦^∗_𝑠,𝑟Σ₃⁻¹Σ₂−𝜆^∗Δ₂𝑦^∗_𝑠,𝑟. (4.10)

5. Applications in Conic Fitting

5.1. Fitzgibbon's Ellipse Fitting

First, we would like to find an ellipse for a given set of points in ℝ². Generally, a conic in ℝ² is implicitly defined as the zero set of 𝑓∶ℝ⁶×ℝ²→ℝ to a constant parameter 𝑎=(𝛼₁,…,𝛼₆)^𝑡∈ℝ⁶: 𝑓𝑎,𝜉,𝜂=𝛼₁𝜉²+𝛼₂𝜉𝜂+𝛼₃𝜂²+𝛼₄𝜉+𝛼₅𝜂+𝛼₆. (5.1) The equation 𝑓(𝑎,𝜉,𝜂)=0 can also be written with 𝑥=(𝜉,𝜂)^𝑡 as 0=𝑥^𝑡𝐴𝑥+𝑏^𝑡𝛼𝑥+𝑐with𝐴=¹𝛼₂𝛼/2₂/2𝛼₃𝛼,𝑏=⁴𝛼₅,𝑐=𝛼₆. (5.2) The eigenvalues 𝜆₁,𝜆₂ of 𝐴 characterize a conic uniquely [14]. Thus, we need 4𝜆₁𝜆₂=4𝛼₁𝛼₃−𝛼²₂>0 for ellipses in 𝑓(𝑎,𝑥)=0. Furthermore, every scaled vector 𝜇𝑎 with 𝜇∈ℝ⧵{0} describes the same zero-set of 𝑓. So, we can impose the constraint for ellipses with 4𝛼₁𝛼₃−𝛼²₂=1. For 𝑛(𝑛≥6) given points (𝜉_𝑖,𝜂_𝑖)^𝑡∈ℝ², we want to find a parameter 𝑎∈ℝ⁶, which minimizes 𝐹∶ℝ⁶→ℝ with𝐹(𝑎)=^𝑛_𝑖=1𝑓𝑎,𝜉_𝑖,𝜂_𝑖²subjectto4𝛼₁𝛼₃−𝛼²₂=1. (5.3) This ellipse fitting problem is established and solved by Fitzgibbon [7]. With the following matrices 𝐷∈ℝ^𝑛×6, 𝐶∈ℝ^6×6𝜉𝐷=²¹𝜉₁𝜂₁𝜂²₁𝜉₁𝜂₁1𝜉²₂𝜉₂𝜂₂𝜂²₂𝜉₁𝜂₂1𝜉⋮⋮⋮⋮⋮⋮_2𝑛𝜉_𝑛𝜂_𝑛𝜂²_𝑛𝜉_𝑛𝜂_𝑛1,𝐶=0020000−10000200000000000000000000000, (5.4) and 𝐹(𝑎)=∑^𝑛_𝑖=1𝑓(𝑎,𝜉_𝑖,𝜂_𝑖)²=‖𝐷𝑎‖²₂, we achieve the equivalent problem:min_𝑎∈ℝ6‖𝐷𝑎‖²₂subjectto𝑎^𝑡𝐶𝑎=1. (5.5) For 𝑆=𝐷^𝑡𝐷, we have a special case of (1.1). Assuming 𝑆 is a regular matrix and since the eigenvalues of 𝐶 are −2, −1, 0, and 2, by lemma 3.3 we know that the generalized eigenvalue problem (𝑆−𝜆𝐶)𝑎=0 (5.6) has exactly one positive solution 𝜆^∗∈ℝ. Because of Theorem 3.13 a corresponding generalized eigenvector 𝑎^∗ to 𝜆^∗ minimizes the problem (5.5) and 𝑎^∗ consists of the coefficients of an implicit given ellipse.

A numerically stable noniterative algorithm to solve this optimization problem is presented by Halir and Flusser [15]. In comparison with Section 4, their method uses a special block decomposition of the matrices 𝐷 and 𝐶.

5.2. Hyperbola Fitting

Instead of ellipses, O'Leary and Zsombor-Murray want to find a hyperbola to a set of scattered data 𝑥_𝑖∈ℝ² [9]. A hyperbola is a conic, which can uniquely be characterized by 4𝜆₁𝜆₂=4𝛼₁𝛼₃−𝛼²₂<0 [14]. So, we consider the constraint 𝛼²₂−4𝛼₁𝛼₃=1 and obtain the optimization problem: min_𝑎∈ℝ6‖𝐷𝑎‖²₂subjectto𝑎^𝑡(−𝐶)𝑎=1 (5.7) with 𝐷 and 𝐶 being chosen in 5.1. The matrix (−𝐶) has two positive eigenvalues. In this case, a solution is given by a generalized eigenvector to the smallest value in 𝜎⁺(𝑆,−𝐶). But O'Leary and Zsombor-Murray determine the best hyperbolic fit by evaluation of 𝜅_𝑖=𝛼²_2,𝑖−4𝛼_1,𝑖𝛼_3,𝑖, where the eigenvector 𝑎_𝑖=(𝑎_𝑖,1,…,𝑎_𝑖,6)^𝑡 is associated to a positive value of 𝜎⁺(𝑆,−𝐶).

5.3. Bookstein's Conic Fitting

In Bookstein's method, the conic constraint is restricted to2trace²𝐴−2det𝐴=2𝜆²₁+𝜆²₂=2𝛼²₁+𝛼²₂+2𝛼²₃=1, (5.8) where 𝜆₁,𝜆₂ are the eigenvalues of 𝐴 in 𝑓 [6]. There, it is 𝜆_1,2∈[−√2/2,√2/2] and at least one of them different to 0. But the constraint (5.8) is not a restriction to a class of conics. Here, we determine an arbitrary conic, which minimizes𝐹(𝑎)=^𝑛_𝑖=1𝑓𝑎,𝜉_𝑖,𝜂_𝑖²subjectto2𝛼²₁+𝛼²₂+2𝛼²₃=1. (5.9) The resulting data matrix 𝐷∈ℝ^𝑛×6 is the same as for Fitzgibbon's problem. The constraint matrix 𝐶∈ℝ^6×6 has a diagonal shape with the entrees (2,1,2,0,0,0), that is, all eigenvalues of 𝐶 are nonnegative. In the case of a regular matrix 𝑆, the problem (5.9) is solved for a generalized eigenvector to the smallest value in 𝜎⁺(𝑆,𝐶).

5.4. Approximation of Shrinked Aspheres

After the molding process in optical applications, the shrinkage of rotation-symmetric aspheres is implicitly defined for 𝑥=(𝜉,𝜁)^𝑡 in 𝛼₁𝜁²+𝛼₂𝜁+𝛼₃+𝛼₄𝜉²=0subjectto𝛼²₂−4𝛼₁𝛼₃=4𝑟², (5.10) where 𝑟∈ℝ⧵{0} and 𝑎=(𝑎₁,…,𝑎₄)^𝑡 are aspheric-specific constants [10]. For 𝑖=1,…,𝑛 with 𝑛≥4, the scattered data 𝑥_𝑖=(𝜉_𝑖,𝜁_𝑖)^𝑡∈ℝ² of a shrinked asphere are given in this approximation problem. Here, we are looking for the conic parameter 𝑎=(𝛼₁,…,𝛼₄)^𝑡 for a fixed value 𝑟_ref, which minimizes 𝐹(𝑎)=^𝑛_𝑖=1𝛼₁𝜁²_𝑖+𝛼₂𝜁_𝑖+𝛼₃+𝛼₄𝜉²_𝑖²subjectto𝛼²₂−4𝛼₁𝛼₃=4𝑟²_ref. (5.11) Analogously to Fitzgibbon, we have the matrices 𝐷∈ℝ^𝑛×4 and 𝐶∈ℝ^4×4 with 𝜁𝐷=²¹𝜁₁1𝜉²₁𝜁²₂𝜁₂1𝜉²₂𝜁⋮⋮⋮⋮_2𝑛𝜁_𝑛1𝜉²_𝑛,𝐶=00−200100−20000000 (5.12) and with 𝐹(𝑎)=‖𝐷𝑎‖²₂ we get the following optimization problem:min_𝑎∈ℝ4‖𝐷𝑎‖²₂subjectto𝑎^𝑡𝐶𝑎=4𝑟²_ref. (5.13) This is also an application of (1.1). The matrix 𝐶 has the eigenvalues −2,0,1, and 2. So, the generalized eigenvalue problem in (3.1) with regular 𝑆=𝐷^𝑡𝐷∈ℝ^4×4 has two positive values in 𝜎⁺(𝑆,𝐶). With Theorem 3.13, a generalized eigenvector 𝑎^∗∈ℝ⁴ to the smaller of both values solves (5.13).

The coefficients 𝛼_𝑖 in the problems (5.5) and (5.13) correspond not to the same monomials𝜉^𝑘𝜁^𝑙. Hence, we have different matrices 𝐷 and 𝐶.

6. Conclusion

In this paper, we present a minimization problem of least squares subject to absolute quadratic constraints. We develop a closed theory with the main result that a minimum is a solution of a generalized eigenvalue problem corresponding to the smallest positive eigenvalue. Further, we show a reduction to an eigensystem for numerical calculations. Finally, we study four applications about conic approximations. We analyze Fitzgibbon's method for direct ellipse-specific fitting, O'Leary's direct hyperbola approximation, Bookstein's conic fitting, and an optical application of shrinked aspheres. All these systems are attribute to the general optimization problem.