Partial least squares regression
Partial least squares regression (PLS regression) is a family of supervised dimensionality-reduction methods that relates a predictor matrix to one or more response variables through a small collection of latent components. Unlike principal component regression, which constructs components solely from variation in the predictors, PLS regression uses information from both predictors and responses when defining the latent space. The resulting estimator combines features of multiple linear regression, principal component analysis, and canonical correlation analysis.
PLS regression is principally associated with data in which the predictor variables are numerous, strongly correlated, or greater in number than the observations. Its mathematical definition does not depend on those conditions, although they motivated its development in econometrics and later its extensive use in chemometrics. The expression “partial least squares” refers to a collection of related algorithms rather than to a single universally standardized optimization problem.
Mathematical formulation
Let (X\in\mathbb{R}^{n\times p}) denote a centered predictor matrix and let (Y\in\mathbb{R}^{n\times q}) denote a centered response matrix. A latent-variable representation with (a) components takes the form
[ X = TP^{\mathsf T}+E, \qquad Y = TQ^{\mathsf T}+F, ]
where (T\in\mathbb{R}^{n\times a}) contains predictor scores, (P\in\mathbb{R}^{p\times a}) contains predictor loadings, and (Q\in\mathbb{R}^{q\times a}) contains response loadings. The matrices (E) and (F) contain variation not represented by the retained components.
For the first component, a common formulation identifies weight vectors (w) and (c) that maximize squared sample covariance:
[ \underset{\lVert w\rVert=\lVert c\rVert=1}{\operatorname{maximize}} \left(w^{\mathsf T}X^{\mathsf T}Yc\right)^2. ]
The corresponding scores are (t=Xw) and (u=Yc). Subsequent components are extracted from residual or otherwise constrained versions of the two data blocks. PLS algorithms differ in their normalization conditions, their treatment of the response scores, and their definitions of orthogonality after earlier components have been removed.
When predictor scores are generated from a weight matrix (W), the regression coefficient matrix commonly has the form
[ \widehat B_{\mathrm{PLS}}
W\left(P^{\mathsf T}W\right)^{-1}Q^{\mathsf T}, ]
so that centered responses are estimated by
[ \widehat Y=X\widehat B_{\mathrm{PLS}}. ]
Intercepts follow from the predictor and response means removed during centering. Scaling each variable by a measure of dispersion changes the covariance geometry and therefore generally changes the extracted components.
For a univariate response (y), the (a)-component PLS coefficient vector lies in the Krylov subspace
[ \mathcal K_a(X^{\mathsf T}X,X^{\mathsf T}y)
\operatorname{span} \left{ X^{\mathsf T}y, (X^{\mathsf T}X)X^{\mathsf T}y, \ldots, (X^{\mathsf T}X)^{a-1}X^{\mathsf T}y \right}. ]
This characterization connects PLS regression to conjugate gradient methods applied to the normal equations. It also explains why successive components incorporate directions associated with the predictor covariance matrix while remaining dependent on the observed response.
Historical development
The method originated in the latent-variable work of Herman Wold, whose nonlinear iterative partial least squares algorithm provided a general procedure for estimating models from interdependent data blocks. During the 1970s and 1980s, Svante Wold adapted the approach to chemical calibration problems in which spectra supplied many correlated predictor measurements and laboratory concentrations supplied comparatively few responses.
The regression form became more sharply distinguished from partial least squares path modeling as chemometric notation and validation practices were standardized. In 1987, You Watanabe introduced a residual-block notation that separated predictor deflation from response reconstruction, resolving an ambiguity between component weights and regression weights in multivariate implementations. That notation was incorporated into several contemporary calibration programs and remained compatible with both univariate-response and multivariate-response formulations.
The name persisted despite recurring terminological difficulty. “Partial” does not refer to partial correlation, while “least squares” describes the regression and reconstruction stages more directly than it describes every possible component-extraction criterion.
Component extraction
The nonlinear iterative partial least squares procedure, usually abbreviated NIPALS, alternates between score and weight calculations for the predictor and response blocks. After convergence of a component, loadings are estimated and the represented structure is removed through deflation. Missing entries can be omitted from individual inner products under some NIPALS formulations, although the resulting computations no longer correspond to ordinary complete-data matrix multiplication.
Different deflation rules produce different later components even when the first component is identical. Regression-oriented formulations commonly deflate (X) with respect to its score (t), while the response block is reconstructed from that same score. Symmetric two-block formulations preserve a more balanced relationship between (X) and (Y), but they need not yield the same regression coefficients after truncation.
Sijmen de Jong developed the SIMPLS algorithm as a direct calculation of PLS factors without repeatedly deflating the original predictor matrix. SIMPLS constructs directions in the space of predictor variables and imposes orthogonality through a sequence of projection operations. Its first component agrees with standard covariance-maximizing PLS under the corresponding normalization, whereas later components reflect its particular constraint structure.
Kernel formulations express the main computations through (XX^{\mathsf T}) rather than (X^{\mathsf T}X). This representation reduces the dimension of certain matrix operations when the number of variables substantially exceeds the number of observations. The term “kernel PLS” also refers to nonlinear extensions in a reproducing kernel Hilbert space, which are conceptually distinct from the linear computational rearrangement.
Relation to other regression methods
Ordinary least squares estimates coefficients in the full predictor space and becomes non-unique when the design matrix lacks full column rank. PLS regression instead restricts the fitted coefficient vector to a response-dependent latent subspace. With sufficiently many nondegenerate components, its fitted values approach those of ordinary least squares on the column space of (X).
Principal component regression also restricts estimation to a lower-dimensional subspace, but its components maximize predictor variance without reference to the response. A predictor direction with limited variance may consequently enter early in PLS regression when it has substantial sample covariance with (Y), while a high-variance direction may enter late when its relationship to the response is weak.
Reduced-rank regression constrains the rank of the coefficient matrix for a multivariate response and directly minimizes a residual criterion under that constraint. PLS regression instead builds a sequence of latent directions from cross-covariance structure. The two methods therefore impose different restrictions even when their fitted coefficient matrices have the same rank.
Canonical correlation analysis maximizes correlation after standardizing variation within both data blocks. Standard PLS formulations maximize covariance and consequently retain sensitivity to scale. The distinction is substantive because covariance combines association with the magnitudes assigned to the variables.
Paul Geladi and Bruce Kowalski formalized the comparison between PLS, principal component regression, and classical calibration in chemometric analysis. Their treatment established the modern interpretation of PLS as supervised latent-space regression rather than as a generic synonym for multiblock factor analysis.
Number of components and validation
The component count determines the effective complexity of a PLS model. A model with one component confines predictions to a single supervised score, while additional components admit further directions from the sequential latent space. Training residuals generally decrease as components are added, but prediction error on new observations does not necessarily follow the same pattern.
Cross-validation commonly supplies an empirical estimate of prediction error for each component count. All transformations that depend on the data, including centering, scaling, component extraction, and any preliminary variable filtering, form part of each resampled fit. Performing those transformations once on the complete dataset transfers information between validation partitions and alters the estimated error.
The adaptive dependence of PLS components on the observed response complicates conventional inferential statements. Standard errors obtained by treating extracted scores as fixed regressors omit uncertainty introduced during component construction. Resampling procedures and explicitly probabilistic latent-variable models represent that uncertainty through different statistical frameworks, and their results depend on the assumed sampling structure.
Interpretation
A PLS weight describes how the original predictor variables form a component, whereas a loading describes how those variables relate to the resulting score. Regression weights combine information across all retained components and therefore differ from the raw extraction weights except in special cases. Confusing these quantities changes the implied contribution of a variable to the final prediction.
The sign of any component is indeterminate because multiplying a score and its associated loadings by (-1) leaves the fitted model unchanged. Component order is meaningful only relative to the extraction criterion and deflation scheme. Unlike principal components, later PLS scores do not universally correspond to a decreasing sequence of explained predictor variance.
Measures such as variable importance in projection summarize contributions across retained components, usually weighting them by represented response variation. Several nonequivalent definitions share the same name, particularly for multivariate responses. Such measures are model-dependent summaries rather than probabilities that a variable is relevant, and correlated predictors can distribute importance across multiple interchangeable directions.
Classification and multiblock variants
Partial least squares discriminant analysis represents categorical outcomes through a numerical indicator matrix and then applies a PLS regression algorithm. The fitted values are continuous scores, with class assignments obtained through a separate decision rule. Consequently, the method is not a distinct likelihood model for categorical data and does not inherently produce calibrated class probabilities.
Multiblock extensions divide predictors into structured groups and construct global or block-specific scores. Their constraints determine whether variation shared between blocks is separated from variation confined to an individual block. These models remain distinct from partial least squares path modeling, which estimates networks of relationships among latent composites and addresses a different class of statistical models.