# A Jacobi-Davidson Method for Large Scale Canonical Correlation Analysis
## 3 The Main Algorithm
### Subspace Extraction  
Let \(\mathcal{X}\subseteq\mathbb{R}^{m}\) and \(\mathcal{Y}\subseteq\mathbb{R}^{n}\) with \(\dim(\mathcal{X})=k\) and \(\dim(\mathcal{Y})=\ell\), respectively. As stated in [31], we call \(\{\mathcal{X},\mathcal{Y}\}\) a pair of defalting subspaces of CCGEP (4) if  
\[C\mathcal{Y}\subseteq A\mathcal{X}\quad\text{and}\quad C^{\mathrm{T}}\mathcal{ X}\subseteq B\mathcal{Y}. \tag{13}\]  
Let \(X\in\mathbb{R}^{m\times k}\) and \(Y\in\mathbb{R}^{n\times\ell}\) be \(A\)- and \(B\)-orthonormal basis matrices of the subspaces \(\mathcal{X}\) and \(\mathcal{Y}\), respectively, i.e.,  
\[X^{\mathrm{T}}AX=I_{k}\quad\text{and}\quad Y^{\mathrm{T}}BY=I_{\ell}.\]The equality (13) implies that there exist \(D_{\text{A}}\in\mathbb{R}^{k\times\ell}\) and \(D_{\text{a}}\in\mathbb{R}^{\ell\times k}\)[32] (Equation (2.11)) such that  
\[CY=AXD_{\text{A}}\quad\text{and}\quad C^{\text{T}}X=BYD_{\text{B}}. \tag{14}\]  
They imply \(D_{\text{A}}=X^{\text{T}}CY=(Y^{\text{T}}C^{\text{T}}X)^{\text{T}}=D_{\text{a }}^{\text{T}}\). So (14) is equivalent to  
\[CY=AXD_{\text{A}}\quad\text{and}\quad C^{\text{T}}X=BYD_{\text{A}}^{\text{T}}.\]  
Now if \((\lambda,\hat{x},\hat{y})\) is a singular triplet of \(D_{\text{A}}\), then \((\lambda,z)\) gives an eigenpair of (4), where \(z=[x^{\text{T}},y^{\text{T}}]^{\text{T}}\), \(x=X\hat{x}\) and \(y=Y\hat{y}\). This is because  
\[Cy=C(Y\hat{y})=(CY)\hat{y}=AXD_{\text{A}}\hat{y}=\lambda AX\hat{x}=\lambda Ax.\]  
and similarly \(C^{\text{T}}x=\lambda By\). That means  
\[ \begin{bmatrix}0&C\\ C^{\text{T}}&0\end{bmatrix}\begin{bmatrix}x\\ y\end{bmatrix} =\lambda \begin{bmatrix}A&0\\ 0&B\end{bmatrix}\begin{bmatrix}x\\ y\end{bmatrix} .\]  
Hence, we have shown that a pair of deflating subspaces \(\{\mathcal{X},\mathcal{Y}\}\) can be used to recover those eigenpairs associated with the pair of deflating subspaces of CCGEP (4). In practice, pairs of exact deflating subspaces are usually not available, and one usually use Lanczos type methods [14] or Davidson type methods [15] to generate approximate ones, such as Krylov subspaces in Lanczos method [33]. Next, we will discuss how to extract best approximate eigenpairs from a given pair of approximate deflating subspaces. In what follows, we consider the simple case: \(k=\ell\). Suppose \(\{\mathcal{U},\mathcal{V}\}\) is an approximation of a pair of deflating subspaces \(\{\mathcal{X},\mathcal{Y}\}\) with \(\dim(\mathcal{U})=\dim(\mathcal{V})=k\). Let \(U\in\mathbb{R}^{m\times k}\) and \(V\in\mathbb{R}^{n\times k}\) be the \(A\)- and \(B\)-orthonormal basis matrices of the subspaces \(\mathcal{U}\) and \(\mathcal{V}\), respectively. Denote \(\theta_{i}\), \(i=1,2,\ldots,k\) the singular values of \(U^{\text{T}}CV\) in descending order with associated left and right singular vectors \(\hat{u}_{i}\) and \(\hat{v}_{i}\), respectively, i.e.,  
\[(U^{\text{T}}CV)\hat{v}_{i}=\theta_{i}\hat{u}_{i}\quad\text{and}\quad(U^{ \text{T}}CV)^{\text{T}}\hat{u}_{i}=\theta_{i}\hat{v}_{i},\quad\text{for }1\leq i \leq k.\]  
Even though \(U\) and \(V\) as \(A\)- and \(B\)-orthonormal basis matrices are not unique, these \(\theta_{i}\) are. Motivated by the maximization principle in Theorem 1, we would seek the best approximations associated with the pair of approximate deflating subspaces \(\{\mathcal{U},\mathcal{V}\}\) to the eigenpairs \((\lambda_{i},z_{i})\) (\(1\leq i\leq j\leq k\)) in the sense of  
\[\max\sum_{i=1}^{j}\sigma_{i}(\widehat{U}_{j}^{\text{T}}C\widehat{V}_{j}). \tag{15}\]  
for any \(\widehat{U}_{j}\in\mathbb{R}^{m\times j}\) and \(\widetilde{V}_{j}\in\mathbb{R}^{n\times j}\) satisfying \(\text{span}(\widehat{U}_{j})\subseteq\mathcal{U}\), \(\text{span}(\widetilde{V}_{j})\subseteq\mathcal{V}\) and \(\widehat{U}_{j}^{\text{T}}A\widehat{U}_{j}=\widetilde{V}_{j}^{\text{T}}B \widetilde{V}_{j}=I_{j}\). We claim that the quantity in (15) is given by \(\sum_{i=1}^{j}\theta_{i}\).