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.
 A: 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. 
A: 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.
