Numerical methods for computing angles between linear subspaces

Abstract

Assume that two subspaces <italic>F</italic> and <italic>G</italic> of a unitary space are defined as the ranges (or null spaces) of given rectangular matrices <italic>A</italic> and <italic>B</italic>. Accurate numerical methods are developed for computing the principal angles <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="theta Subscript k Baseline left-parenthesis upper F comma upper G right-parenthesis"> <mml:semantics> <mml:mrow> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msub> <mml:mi>θ</mml:mi> <mml:mi>k</mml:mi> </mml:msub> </mml:mrow> <mml:mo stretchy="false">(</mml:mo> <mml:mi>F</mml:mi> <mml:mo>,</mml:mo> <mml:mi>G</mml:mi> <mml:mo stretchy="false">)</mml:mo> </mml:mrow> <mml:annotation encoding="application/x-tex">{\theta _k}(F,G)</mml:annotation> </mml:semantics> </mml:math> </inline-formula> and orthogonal sets of principal vectors <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="u Subscript k Baseline element-of upper F"> <mml:semantics> <mml:mrow> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msub> <mml:mi>u</mml:mi> <mml:mi>k</mml:mi> </mml:msub> </mml:mrow> <mml:mo>∈</mml:mo> <mml:mi>F</mml:mi> </mml:mrow> <mml:annotation encoding="application/x-tex">{u_k} \in F</mml:annotation> </mml:semantics> </mml:math> </inline-formula> and <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="v Subscript k Baseline element-of upper G comma k equals 1 comma 2 comma midline-horizontal-ellipsis comma q equals dimension left-parenthesis upper G right-parenthesis less-than-over-equals dimension left-parenthesis upper F right-parenthesis"> <mml:semantics> <mml:mrow> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msub> <mml:mi>v</mml:mi> <mml:mi>k</mml:mi> </mml:msub> </mml:mrow> <mml:mo>∈</mml:mo> <mml:mi>G</mml:mi> <mml:mo>,</mml:mo> <mml:mi>k</mml:mi> <mml:mo>=</mml:mo> <mml:mn>1</mml:mn> <mml:mo>,</mml:mo> <mml:mn>2</mml:mn> <mml:mo>,</mml:mo> <mml:mo>⋯</mml:mo> <mml:mo>,</mml:mo> <mml:mi>q</mml:mi> <mml:mo>=</mml:mo> <mml:mi>dim</mml:mi> <mml:mo>⁡</mml:mo> <mml:mo stretchy="false">(</mml:mo> <mml:mi>G</mml:mi> <mml:mo stretchy="false">)</mml:mo> <mml:mo>≦</mml:mo> <mml:mi>dim</mml:mi> <mml:mo>⁡</mml:mo> <mml:mo stretchy="false">(</mml:mo> <mml:mi>F</mml:mi> <mml:mo stretchy="false">)</mml:mo> </mml:mrow> <mml:annotation encoding="application/x-tex">{v_k} \in G,k = 1,2, \cdots ,q = \dim (G) \leqq \dim (F)</mml:annotation> </mml:semantics> </mml:math> </inline-formula>. An important application in statistics is computing the canonical correlations <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="sigma Subscript k Baseline equals cosine theta Subscript k"> <mml:semantics> <mml:mrow> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msub> <mml:mi>σ</mml:mi> <mml:mi>k</mml:mi> </mml:msub> </mml:mrow> <mml:mo>=</mml:mo> <mml:mi>cos</mml:mi> <mml:mo>⁡</mml:mo> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msub> <mml:mi>θ</mml:mi> <mml:mi>k</mml:mi> </mml:msub> </mml:mrow> </mml:mrow> <mml:annotation encoding="application/x-tex">{\sigma _k} = \cos {\theta _k}</mml:annotation> </mml:semantics> </mml:math> </inline-formula> between two sets of variates. A perturbation analysis shows that the condition number for <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="theta Subscript k"> <mml:semantics> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msub> <mml:mi>θ</mml:mi> <mml:mi>k</mml:mi> </mml:msub> </mml:mrow> <mml:annotation encoding="application/x-tex">{\theta _k}</mml:annotation> </mml:semantics> </mml:math> </inline-formula> essentially is <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="max left-parenthesis kappa left-parenthesis upper A right-parenthesis comma kappa left-parenthesis upper B right-parenthesis right-parenthesis"> <mml:semantics> <mml:mrow> <mml:mo movablelimits="true" form="prefix">max</mml:mo> <mml:mo stretchy="false">(</mml:mo> <mml:mi>κ</mml:mi> <mml:mo stretchy="false">(</mml:mo> <mml:mi>A</mml:mi> <mml:mo stretchy="false">)</mml:mo> <mml:mo>,</mml:mo> <mml:mi>κ</mml:mi> <mml:mo stretchy="false">(</mml:mo> <mml:mi>B</mml:mi> <mml:mo stretchy="false">)</mml:mo> <mml:mo stretchy="false">)</mml:mo> </mml:mrow> <mml:annotation encoding="application/x-tex">\max (\kappa (A),\kappa (B))</mml:annotation> </mml:semantics> </mml:math> </inline-formula>, where <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="kappa"> <mml:semantics> <mml:mi>κ</mml:mi> <mml:annotation encoding="application/x-tex">\kappa</mml:annotation> </mml:semantics> </mml:math> </inline-formula> denotes the condition number of a matrix. The algorithms are based on a preliminary <italic>QR</italic>-factorization of <italic>A</italic> and <italic>B</italic> (or <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="upper A Superscript upper H"> <mml:semantics> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msup> <mml:mi>A</mml:mi> <mml:mi>H</mml:mi> </mml:msup> </mml:mrow> <mml:annotation encoding="application/x-tex">{A^H}</mml:annotation> </mml:semantics> </mml:math> </inline-formula> and <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="upper B Superscript upper H"> <mml:semantics> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msup> <mml:mi>B</mml:mi> <mml:mi>H</mml:mi> </mml:msup> </mml:mrow> <mml:annotation encoding="application/x-tex">{B^H}</mml:annotation> </mml:semantics> </mml:math> </inline-formula>), for which either the method of Householder transformations (HT) or the modified Gram-Schmidt method (MGS) is used. Then <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="cosine theta Subscript k"> <mml:semantics> <mml:mrow> <mml:mi>cos</mml:mi> <mml:mspace width="thickmathspace"/> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msub> <mml:mi>θ</mml:mi> <mml:mi>k</mml:mi> </mml:msub> </mml:mrow> </mml:mrow> <mml:annotation encoding="application/x-tex">\cos \;{\theta _k}</mml:annotation> </mml:semantics> </mml:math> </inline-formula> and <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="sine theta Subscript k"> <mml:semantics> <mml:mrow> <mml:mi>sin</mml:mi> <mml:mspace width="thickmathspace"/> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msub> <mml:mi>θ</mml:mi> <mml:mi>k</mml:mi> </mml:msub> </mml:mrow> </mml:mrow> <mml:annotation encoding="application/x-tex">\sin \;{\theta _k}</mml:annotation> </mml:semantics> </mml:math> </inline-formula> are computed as the singular values of certain related matrices. Experimental results are given, which indicates that MGS gives <inline-formula content-type="math/mathml"> <mml:math xmlns:mml="http://www.w3.org/1998/Math/MathML" alttext="theta Subscript k"> <mml:semantics> <mml:mrow class="MJX-TeXAtom-ORD"> <mml:msub> <mml:mi>θ</mml:mi> <mml:mi>k</mml:mi> </mml:msub> </mml:mrow> <mml:annotation encoding="application/x-tex">{\theta _k}</mml:annotation> </mml:semantics> </mml:math> </inline-formula> with equal precision and fewer arithmetic operations than HT. However, HT gives principal vectors, which are orthogonal to working accuracy, which is not generally true for MGS. Finally, the case when <italic>A</italic> and/or <italic>B</italic> are rank deficient is discussed.