# What is an analogue of PCA in the regression context?

I'm writing code to approximate a function $y=f(\vec{x})$ where $y\in\mathbb{R}$ and $\vec{x}\in\mathbb{R}^N$ for medium-sized $N$ ($N$ between 20 and 50, approx.). I have a ton of examples, however, so my number of examples greatly exceeds $N$. But, the examples probably are somewhat redundant, and some of the $N$ variables probably aren't correlated much with $y$. While I can compute $\vec{x}$ in near zero time, I need to evaluate $f(\vec{x})$ often/quickly. So, I'd like to project $\vec{x}$ onto a lower-dimensional subspace before learning/evaluating $f(\vec{x})$.

I know the standard technique for finding a lower-dimensional projection given a lot of examples is PCA. But in the regression context, I care more about the correlation of the projection of $\vec{x}$ onto each of the basis vectors with the value of $y$ then simply maximizing the variance of the projections. So my basic question is: how can I compute such a(n) (orthogonal) basis?

I've looked into "partial least squares" and "correlated components analysis," but if I understand properly, they will only compute a single basis vector since $y\in\mathbb{R}^1$. I want fewer than $N$ vectors in my (approximate) basis but more than one. If it matters, I'll be using a kernelized regression technique, like $\nu$-SVR, LS-SVR, or kernel ridge regression.

• Are you looking for the loadings of the PLS rather than the regression coefficients? They give a projection into $\mathbb R^{n \leq N}$ and you can decide to orthogonalize them. Commented Aug 30, 2011 at 7:39
• It was my understanding that because $y=f(\vec{x})\in\mathbb{R}$, I would have to choose $n=1$ for PLS. Is this not the case? I'll keep reading about PLS to see if I missed something. Commented Aug 30, 2011 at 14:15
• Actually, Matlab's PLSREGRESS function provides a full basis somehow! Is there a listing of the algorithm or what is happening from a linear algebra standpoint? Commented Aug 30, 2011 at 14:22
• You can have more than one component in normal PLS. However, I don't know whether that is true as well if you orthogonalize the loadings. Afaik, the PLS Toolbox for Matlab does not allow you to read the code, but the pls package for R is open source, so you could look into the different PLS algorithms there. Commented Sep 1, 2011 at 12:14
• [just saw Denis's response] Unfortunately this project is ongoing! We're currently looking into metric learning as a potential solution. Commented Aug 23, 2012 at 21:11

I think you want to use canonical correlation analysis (I am not sure, but I guess that is the same as your "correlated component analysis"). Canonical correlation analysis allows you to compute more than one orthogonal basis vectors in the output. You can check wikipedia for the computations but probably there are better tutorials online.

If you want to use a kernel regression algorithm afterwards I would recommend to use kernel canonical correlation analysis. The reason is that you get maximally correlated directions in the feature space right away. If you first do linear canonical correlation analysis and then stuff it into a kernel you might throw away nonlinear information the kernel could have extracted. However, the problem is that the eigenvalue problem for getting the directions scales with the number of datapoints in the kernel case (and not the number of dimensions like in the linear case). Since you have tons of data this might become a problem. If you can still fit the kernel matrix into memory you might still succeed since you only have to compute the first $n$ eigenvectors. Otherwise you might have to use a subset of the data or use a more sophisticated large-scale method.

Also note that if you use the features of kernel CCA you apply a linear regression algorithm to it afterwards and don't stuff it into another kernel. Since the coefficients of kernel CCA are basis coefficients in features space you are implicitly learning in the kernel space anyway if you use a linear regressor.

Extension after comments 1 The optimal linear predictor and the CCA component for 1d $y$ are collinear, which means that first doing CCA and the learnin a regression is equivalent. Let $C_{xx}$ and $C_{yx}$ be the covariance and the crosscovariance matrix, respectively. The optimal linear predictor $v$ is given by $v = C_{yx}C_{xx}^{-1}$. If we look at the objective for CCA in the 1d output case $$\mbox{maximize} \frac{a C_{yx} u}{a\sqrt{C_{yy}}\sqrt{u^\top C_{xx} u}}$$ when can make several adjustements to make the problem easier. First, since $a$ is 1d and will be scaled to length one later (as a basis vector), we know that $a=1$ so we don't need to optimize over. Second, maximizing the fraction is the same as keeping the denominator fixed and maximizing over the numerator. Casting the problem into a Lagrangian optimization problem, we get $$\mbox{maximize}_u C_{yx} u - \lambda u^\top C_{xx} u.$$ Substituting $w=C_{xx}^{\frac{1}{2}}u$ yields $$\mbox{maximize}_w C_{yx} C_{xx}^{-\frac{1}{2}}w - \lambda w^\top w .$$ Taking the derivative and setting it to zero yields $$C_{yx} C_{xx}^{-\frac{1}{2}} = 2\lambda w^\top$$ or $$C_{yx} C_{xx}^{-1} = 2\lambda u^\top$$ which means that $u$ is just a rescaled version of $v$. In CCA $u$ would be adjusted to have length one. The regression afterwards would scale it such that it minimizes the mean squared error.

Extension after comments 2 Note that my heuristic in the comments only makes sense if you train nonlinear regressors. If you train several linear regressors and sum their outputs you can equivalently sum the linear regressors.

• Wrt time: I just gave a try on a random $10^5 \times 100$ matrix with PLS (R, package pls, one core, no optimized BLAS): 50 PLS components take about 15.5 s user time + 3.5 s system time. 10 components: 3.3 s + 1.4 s. Obviously, optimized + multi-thread BLAS helps a lot. Commented Aug 30, 2011 at 11:38
• Thanks for your help! I'll take a look at kernel CCA but am still a bit confused about the non-kernel CCA case. In particular, the Wikipedia article says "This procedure may be continued up to $\min\{m,n\}$ times." Isn't $m=1$ in my case? If so, this means I can only get one basis vector out of it, right? Commented Aug 30, 2011 at 13:32
• Also, thanks for running a timing test! I assume by PLS you mean "partial least squares." Were you able to get 50 components to span the column space of your matrix? Does this somehow model a function from $\mathbb{R}^{50}$ to $\mathbb{R}$? I would think you would need to matrices for this. I appreciate your patience! As long as the learning step takes reasonable time, I'm happy with it -- the evaluation is what must go very quickly (will be running it millions of times). (@cbeleites for this comment; @fabee for the last -- it won't let me edit it again!) Commented Aug 30, 2011 at 13:36
• Uuh, sorry, I thought your output was also multidimensional. In that case you can still use kernel CCA or you can do something like this: 1) You run CCA and get one direction in the input. 2) You train a regressor from $\mathbb R$ to $\mathbb R$ and subtract its prediction from the output $y$ to get the residuals. 3) You project the input data onto the orthogonal complement of the CCA component. 4) Start again with 1) but now with the residuals and the deflated intput. You probably have to make sure that you don't overfit by doing proper cross-validation. Commented Aug 30, 2011 at 13:49
• I was doing something like this before but wasn't sure it was principled :-) . To check, can you collapse the first two steps you mention into just solving the least-squares problem finding $\vec{v}$ such that $\vec{v}\cdot\vec{x}\approx f(\vec{x})$ or does something else come out of CCA? Commented Aug 30, 2011 at 14:06

Iterate regressions, adding orthogonality constraints at each step:
$x_1 \perp x_0, x_2 \perp x_1 x_0 \dots$

$\begin{bmatrix} A \end{bmatrix} x_0 \approx \begin{bmatrix} b \end{bmatrix}$

$\begin{bmatrix} x_0 \\ A \end{bmatrix} x_1 \approx \begin{bmatrix} 0 \\ b \end{bmatrix}$

$\begin{bmatrix} x_0 \\ x_1 \\ A \end{bmatrix} x_2 \approx \begin{bmatrix} 0 \\ 0 \\ b \end{bmatrix}$

$\dots$

This may be related to Canonical correlation analysis as suggested above, not sure.