# A Jacobi-Davidson Method for Large Scale Canonical Correlation Analysis
## 3 The Main Algorithm
### Correction Equation
Suppose \(\lambda_{i}\) for \(1\leq i\leq j\), \(X_{j}=[x_{1},\ldots,x_{j}]\) and \(Y_{j}=[y_{1},\ldots,y_{j}]\) have been computed where \(j\leq k\). We can consider the generalized eigenvalue problem \[\widetilde{K}z=\lambda Mz,\] (24) where \[\widetilde{K}= \begin{bmatrix}I_{m}-AX_{j}X_{j}^{\mathrm{T}}&0\\ 0&I_{n}-BY_{j}Y_{j}^{\mathrm{T}}\end{bmatrix}\begin{bmatrix}0&C\\ C^{\mathrm{T}}&0\end{bmatrix}\begin{bmatrix}I_{m}-X_{j}X_{j}^{\mathrm{T}}A&0 \\ 0&I_{n}-Y_{j}Y_{j}^{\mathrm{T}}B\end{bmatrix} .\] (25) By (11), it is clear that the eigenvalues of (24) consist of two groups. Those eigenvalues associated with the eigenvectors \([x_{1}^{\mathrm{T}},y_{1}^{\mathrm{T}}]^{\mathrm{T}},\ldots,[x_{j}^{\mathrm{T }},y_{j}^{\mathrm{T}}]^{\mathrm{T}},[x_{1}^{\mathrm{T}},-y_{1}^{\mathrm{T}}]^{ \mathrm{T}},\ldots,[x_{j}^{\mathrm{T}},-y_{j}^{\mathrm{T}}]^{\mathrm{T}}\) are shifted to zero and the others remain unchanged. Furthermore, for the correction equation, we find \(s\) and \(t\) subject to additional \(A\)- and \(B\)-orthogonality constrains for \(s\) and \(t\) against \(X_{j}\) and \(Y_{j}\), respectively. By a similar derivation of (22), the correction equation after deflation becomes \[ \begin{bmatrix}I_{m}-A\tilde{x}\tilde{x}^{\mathrm{T}}&0\\ 0&I_{n}-B\tilde{y}\tilde{y}^{\mathrm{T}}\end{bmatrix}\begin{bmatrix}-\theta_{1 }A&\widehat{C}\\ \widehat{C}^{\mathrm{T}}&-\theta_{1}B\end{bmatrix}\begin{bmatrix}I_{m}-\tilde {x}\tilde{x}^{\mathrm{T}}A&0\\ 0&I_{n}-\tilde{y}\tilde{y}^{\mathrm{T}}B\end{bmatrix}\begin{bmatrix}s\\ t\end{bmatrix} =- \begin{bmatrix}r_{s}^{0}\\ r_{s}^{0}\end{bmatrix} ,\] (26) where \(\widehat{C}=(I_{m}-AX_{j}X_{j}^{\mathrm{T}})C(I_{n}-Y_{j}Y_{j}^{\mathrm{T}}B)\). Notice that \(s\perp_{A}X_{j}\) and \(t\perp_{B}Y_{j}\) mean \(U\perp_{A}X_{j}\) and \(V\perp_{s}Y_{j}\) in Algorithm 1, respectively. It follows that \(U^{\mathrm{T}}\widehat{C}V=U^{\mathrm{T}}CV\). 4. At step 5, LAPACK's routine xGESVD can be used to solve the singular value problem of \(U^{\mathrm{T}}CV\) because of its small size, where \(U^{\mathrm{T}}CV\) takes the following form:This form is preserved in the algorithm during refining the basis \(U\) and \(V\) at step 8. The new basis matrices \(U\widehat{U}\) and \(V\widehat{V}\) are reassigned to \(U\) and \(V\), respectively. Although a few extra costs are incurred, this refinement is necessary in order to have faster convergence for eigenvectors as stated in [36; 37]. Furthermore, the restart is easily executed by keeping the first \(s_{\min}\) columns of \(U\) and \(V\) when the dimension of the subspaces \(\operatorname{span}\{U\}\) and \(\operatorname{span}\{V\}\) exceeds \(s_{\max}\). The restart technique appears at step 8 to keep the size of \(U\), \(V\) and \(U^{\mathsf{T}}CV\) small. There are many ways to specify \(s_{\max}\) and \(s_{\min}\). In our numerical examples, we just simply take \(s_{\max}=3k\) and \(s_{\min}=k\). ```
0: Initial vectors \(u_{0}\), \(v_{0}\), \(s=u_{0}\), \(t=v_{0}\) and \(V=U=[\,]\). 0: Converged canonical weight vectors \(\tilde{x}_{i}\) and \(\tilde{y}_{i}\) for \(i=1,\ldots,k\). 1:for\(iter=1,2,\ldots,\)until convergence do
2:\(A\)- and \(B\)-orthogonal \(s\) and \(t\) against \(U\) and \(V\), respectively, to obtain \(\tilde{s}\) and \(\tilde{t}\). **Theorem 3**. : _Assume the condition of Lemma 1, \(\sin\bigtriangleup_{\iota}(x,\tilde{x})=\mathcal{O}(\epsilon)\) and \(\sin\bigtriangleup_{\iota}(y,\tilde{y})=\mathcal{O}(\epsilon)\). Let \([s^{\mathrm{T}},t^{\mathrm{T}}]^{\mathrm{T}}\) be the exact solution of the correction Equation (22). Then,_

\[|\sin\bigtriangleup_{\iota}(x,\tilde{x}+s)|=\mathcal{O}(\epsilon^{3})\quad \text{and}\quad|\sin\bigtriangleup_{\iota}(y,\tilde{y}+t)|=\mathcal{O}( \epsilon^{3}). \tag{27}\]

## 4 Numerical Examples

In this section, we present some numerical examples to illustrate the effectiveness of Algorithm 1. Our goal is to compute the first few canonical weight vectors.