Householder transformation

From HandWiki
Short description: Concept in linear algebra

In linear algebra, a Householder transformation (also known as a Householder reflection or elementary reflector) is a linear transformation that describes a reflection about a plane or hyperplane containing the origin. The Householder transformation was used in a 1958 paper by Alston Scott Householder.[1]

Its analogue over general inner product spaces is the Householder operator.

Definition

Transformation

The reflection hyperplane can be defined by its normal vector, a unit vector [math]\displaystyle{ v }[/math] (a vector with length [math]\displaystyle{ 1 }[/math]) that is orthogonal to the hyperplane. The reflection of a point [math]\displaystyle{ x }[/math] about this hyperplane is the linear transformation:

[math]\displaystyle{ x - 2\langle x, v\rangle v = x - 2v\left(v^* x\right), }[/math]

where [math]\displaystyle{ v }[/math] is given as a column unit vector with conjugate transpose [math]\displaystyle{ v^\textsf{*} }[/math].

Householder matrix

The matrix constructed from this transformation can be expressed in terms of an outer product as:

[math]\displaystyle{ P = I - 2vv^* }[/math]

is known as the Householder matrix, where [math]\displaystyle{ I }[/math] is the identity matrix.

Properties

The Householder matrix has the following properties:

  • it is Hermitian: [math]\displaystyle{ P = P^* }[/math],
  • it is unitary: [math]\displaystyle{ P^{-1} = P^* }[/math],
  • hence it is involutory: [math]\displaystyle{ P = P^{-1} }[/math].
  • A Householder matrix has eigenvalues [math]\displaystyle{ \pm 1 }[/math]. To see this, notice that if [math]\displaystyle{ u }[/math] is orthogonal to the vector [math]\displaystyle{ v }[/math] which was used to create the reflector, then [math]\displaystyle{ Pu = u }[/math], i.e., [math]\displaystyle{ 1 }[/math] is an eigenvalue of multiplicity [math]\displaystyle{ n - 1 }[/math], since there are [math]\displaystyle{ n - 1 }[/math] independent vectors orthogonal to [math]\displaystyle{ v }[/math]. Also, notice [math]\displaystyle{ Pv = -v }[/math], and so [math]\displaystyle{ -1 }[/math] is an eigenvalue with multiplicity [math]\displaystyle{ 1 }[/math].
  • The determinant of a Householder reflector is [math]\displaystyle{ -1 }[/math], since the determinant of a matrix is the product of its eigenvalues, in this case one of which is [math]\displaystyle{ -1 }[/math] with the remainder being [math]\displaystyle{ 1 }[/math] (as in the previous point).

Applications

Geometric optics

In geometric optics, specular reflection can be expressed in terms of the Householder matrix (see Specular reflection § Vector formulation).

Numerical linear algebra

Householder transformations are widely used in numerical linear algebra, for example, to annihilate the entries below the main diagonal of a matrix,[2] to perform QR decompositions and in the first step of the QR algorithm. They are also widely used for transforming to a Hessenberg form. For symmetric or Hermitian matrices, the symmetry can be preserved, resulting in tridiagonalization.[3]

QR decomposition

Householder reflections can be used to calculate QR decompositions by reflecting first one column of a matrix onto a multiple of a standard basis vector, calculating the transformation matrix, multiplying it with the original matrix and then recursing down the [math]\displaystyle{ (i, i) }[/math] minors of that product. To accomplish this, a hermitian unitary matrix Q is sought which takes a complex vector x into a complex multiple of a complex vector e. For the QR decomposition, e will be a unit coordinate vector, say for the kth coordinate. A complex matrix Q having the form Q = I - u u* with u* u = 2 has the desired property of being hermitian unitary. Here * denotes the conjugate transpose. Since the only two vectors involved are x and e, the vector u must have the form u = a x + b e, where a and b are complex coefficients to be determined. Since an overall phase factor for u does not matter, the coefficient a can be chosen to be positive real. Now Q x = x (1 – a (u* x)) - e (b (u* x)). For the coefficient of the vector x to be zero, the two terms in u* x must have the same phase within a multiple of 180 degrees, so we must have arg(b) = arg(e* x) within a multiple of 180 degrees. There are two solutions according to whether an even or odd multiple of 180 degrees is chosen. So for the complex coefficient of the vector x to be zero, two linear equations in the moduli of a and b must be solved. The phases of a and b have already been determined. Suppose e to be the kth unit coordinate vector so that e* e = 1 and xk= e* x and let |x|= sqrt(x* x). Then a and b may be expressed simply either as a =1/sqrt (|x| (|x|+ |xk|)) and b = a |x| exp(i*arg(xk)) or as a =1/sqrt (|x| (|x|- |xk|)) and b = - a |x| exp(i*arg(xk)). The multiplier of e is -b/a for both solutions. The first solution is better because the denominator cannot be near zero compared to |x|.

Tridiagonalization

Main page: Tridiagonal matrix

This procedure is presented in Numerical Analysis by Burden and Faires. It uses a slightly altered [math]\displaystyle{ \operatorname{sgn} }[/math] function with [math]\displaystyle{ \operatorname{sgn}(0) = 1 }[/math].[4] In the first step, to form the Householder matrix in each step we need to determine [math]\displaystyle{ \alpha }[/math] and [math]\displaystyle{ r }[/math], which are:

[math]\displaystyle{ \begin{align} \alpha &= -\operatorname{sgn}\left(a_{21}\right)\sqrt{\sum_{j=2}^n a_{j1}^2}; \\ r &= \sqrt{\frac{1}{2}\left(\alpha^2 - a_{21}\alpha\right)}; \end{align} }[/math]

From [math]\displaystyle{ \alpha }[/math] and [math]\displaystyle{ r }[/math], construct vector [math]\displaystyle{ v }[/math]:

[math]\displaystyle{ v^{(1)} = \begin{bmatrix} v_1 \\ v_2 \\ \vdots \\ v_n \end{bmatrix}, }[/math]

where [math]\displaystyle{ v_1 = 0 }[/math], [math]\displaystyle{ v_2 = \frac{a_{21} - \alpha}{2r} }[/math], and

[math]\displaystyle{ v_k = \frac{a_{k1}}{2r} }[/math] for each [math]\displaystyle{ k = 3, 4 \ldots n }[/math]

Then compute:

[math]\displaystyle{ \begin{align} P^1 &= I - 2v^{(1)} \left(v^{(1)}\right)^\textsf{T} \\ A^{(2)} &= P^1 AP^1 \end{align} }[/math]

Having found [math]\displaystyle{ P^1 }[/math] and computed [math]\displaystyle{ A^{(2)} }[/math] the process is repeated for [math]\displaystyle{ k = 2, 3, \ldots, n - 2 }[/math] as follows:

[math]\displaystyle{ \begin{align} \alpha &= -\operatorname{sgn}\left(a^k_{k+1,k}\right)\sqrt{\sum_{j=k+1}^n \left(a^k_{jk}\right)^2} \\[2pt] r &= \sqrt{\frac{1}{2}\left(\alpha^2 - a^k_{k+1,k}\alpha\right)} \\[2pt] v^k_1 &= v^k_2 = \cdots = v^k_k = 0 \\[2pt] v^k_{k+1} &= \frac{a^k_{k+1,k} - \alpha}{2r} \\ v^k_j &= \frac{a^k_{jk}}{2r} \text{ for } j = k + 2,\ k + 3,\ \ldots,\ n \\ P^k &= I - 2v^{(k)} \left(v^{(k)}\right)^\textsf{T} \\ A^{(k+1)} &= P^k A^{(k)}P^k \end{align} }[/math]

Continuing in this manner, the tridiagonal and symmetric matrix is formed.

Examples

In this example, also from Burden and Faires,[4] the given matrix is transformed to the similar tridiagonal matrix A3 by using the Householder method.

[math]\displaystyle{ \mathbf{A} = \begin{bmatrix} 4 & 1 & -2 & 2 \\ 1 & 2 & 0 & 1 \\ -2 & 0 & 3 & -2 \\ 2 & 1 & -2 & -1 \end{bmatrix}, }[/math]

Following those steps in the Householder method, we have:

The first Householder matrix:

[math]\displaystyle{ \begin{align} Q_1 &= \begin{bmatrix} 1 & 0 & 0 & 0 \\ 0 & -\frac{1}{3} & \frac{2}{3} & -\frac{2}{3} \\ 0 & \frac{2}{3} & \frac{2}{3} & \frac{1}{3} \\ 0 & -\frac{2}{3} & \frac{1}{3} & \frac{2}{3} \end{bmatrix}, \\ A_2 = Q_1 A Q_1 &= \begin{bmatrix} 4 & -3 & 0 & 0 \\ -3 & \frac{10}{3} & 1 & \frac{4}{3} \\ 0 & 1 & \frac{5}{3} & -\frac{4}{3} \\ 0 & \frac{4}{3} & -\frac{4}{3} & -1 \end{bmatrix}, \end{align} }[/math]

Used [math]\displaystyle{ A_2 }[/math] to form

[math]\displaystyle{ \begin{align} Q_2 &= \begin{bmatrix} 1 & 0 & 0 & 0 \\ 0 & 1 & 0 & 0 \\ 0 & 0 & -\frac{3}{5} & -\frac{4}{5} \\ 0 & 0 & -\frac{4}{5} & \frac{3}{5} \end{bmatrix}, \\ A_3 = Q_2 A_2 Q_2 &= \begin{bmatrix} 4 & -3 & 0 & 0 \\ -3 & \frac{10}{3} & -\frac{5}{3} & 0 \\ 0 & -\frac{5}{3} & -\frac{33}{25} & \frac{68}{75} \\ 0 & 0 & \frac{68}{75} & \frac{149}{75} \end{bmatrix}, \end{align} }[/math]

As we can see, the final result is a tridiagonal symmetric matrix which is similar to the original one. The process is finished after two steps.

Computational and theoretical relationship to other unitary transformations

The Householder transformation is a reflection about a hyperplane with unit normal vector [math]\displaystyle{ v }[/math], as stated earlier. An [math]\displaystyle{ N }[/math]-by-[math]\displaystyle{ N }[/math] unitary transformation [math]\displaystyle{ U }[/math] satisfies [math]\displaystyle{ UU^* = I }[/math]. Taking the determinant ([math]\displaystyle{ N }[/math]-th power of the geometric mean) and trace (proportional to arithmetic mean) of a unitary matrix reveals that its eigenvalues [math]\displaystyle{ \lambda_i }[/math] have unit modulus. This can be seen directly and swiftly:

[math]\displaystyle{ \begin{align} \frac{\operatorname{Trace}\left(UU^*\right)}{N} &= \frac{\sum_{j=1}^N\left|\lambda_j\right|^2}{N} = 1, & \operatorname{det}\left(UU^*\right) &= \prod_{j=1}^N \left|\lambda_j\right|^2 = 1. \end{align} }[/math]

Since arithmetic and geometric means are equal if the variables are constant (see inequality of arithmetic and geometric means), we establish the claim of unit modulus.

For the case of real valued unitary matrices we obtain orthogonal matrices, [math]\displaystyle{ UU^\textsf{T} = I }[/math]. It follows rather readily (see orthogonal matrix) that any orthogonal matrix can be decomposed into a product of 2 by 2 rotations, called Givens Rotations, and Householder reflections. This is appealing intuitively since multiplication of a vector by an orthogonal matrix preserves the length of that vector, and rotations and reflections exhaust the set of (real valued) geometric operations that render invariant a vector's length.

The Householder transformation was shown to have a one-to-one relationship with the canonical coset decomposition of unitary matrices defined in group theory, which can be used to parametrize unitary operators in a very efficient manner.[5]

Finally we note that a single Householder transform, unlike a solitary Givens transform, can act on all columns of a matrix, and as such exhibits the lowest computational cost for QR decomposition and tridiagonalization. The penalty for this "computational optimality" is, of course, that Householder operations cannot be as deeply or efficiently parallelized. As such Householder is preferred for dense matrices on sequential machines, whilst Givens is preferred on sparse matrices, and/or parallel machines.

See also

Notes

  1. "Unitary Triangularization of a Nonsymmetric Matrix". Journal of the ACM 5 (4): 339–342. 1958. doi:10.1145/320941.320947. https://hal.archives-ouvertes.fr/hal-01316095/file/p339householderb.pdf. 
  2. Taboga, Marco. "Householder matrix, Lectures on matrix algebra.". https://www.statlect.com/matrix-algebra/Householder-matrix. 
  3. Schabauer, Hannes; Pacher, Christoph; Sunderland, Andrew G.; Gansterer, Wilfried N. (2010-05-01). "Toward a parallel solver for generalized complex symmetric eigenvalue problems" (in en). Procedia Computer Science 1 (1): 437–445. doi:10.1016/j.procs.2010.04.047. 
  4. 4.0 4.1 Burden, Richard; Faires, Douglas; Burden, Annette (2016). Numerical analysis (10th ed.). Thomson Brooks/Cole. ISBN 9781305253667. 
  5. Renan Cabrera; Traci Strohecker; Herschel Rabitz (2010). "The canonical coset decomposition of unitary matrices through Householder transformations". Journal of Mathematical Physics 51 (8): 082101. doi:10.1063/1.3466798. Bibcode2010JMP....51h2101C. 

References