Principal component analysis replaces correlated measurements by orthogonal linear coordinates ordered by their variance. For centered observations collected in an matrix , write the sample covariance matrix as . A unit loading vector gives score vector with sample variance . The Rayleigh quotient is maximized by a unit eigenvector of for its largest eigenvalue.
The subsequent principal components maximize the same quadratic form subject to orthogonality to the earlier loading vectors. By the spectral theorem, if with descending eigenvalues, the scores have diagonal covariance . Their fractions of explained variance are . Repeated eigenvalues identify an eigenspace rather than a unique axis.
Keeping the first principal components gives the rank- reconstruction . Its squared reconstruction error is
For any rank- orthogonal projection , retained variance is . The weights lie in and sum to , so this is at most the sum of the largest eigenvalues. Projecting onto any candidate reconstruction subspace is its best least-squares reconstruction; this proves the minimum-error property and agrees with the singular value decomposition. Thus PCA finds the linear subspace preserving the greatest variance, equivalently minimizing squared reconstruction error. The upper-left sketch shows the first axis aligned with the elongated cloud; the second is perpendicular. A scree plot and substantive interpretability help choose .
Figure 1.
Original sketches of principal components, classical scaling, hierarchical clustering and multivariate mean comparison
. Details.
These methods depend on measurement scale. A variable measured in large numerical units can dominate covariance-based principal component analysis; standardizing variables gives principal component analysis on a correlation matrix. Scores summarize patterns, but the method neither establishes causation nor guarantees that the largest-variance directions best predict a separate response.
The command centres the log measurements and performs principal component analysis of their covariance matrix. With cor=F, it does not first standardize the variables to variance one. For a unit principal component loading vector , the score has sample variance . By the spectral theorem for real symmetric matrices, expand in an orthonormal eigenbasis; that variance is a weighted average of the eigenvalues. Its maximum is the largest eigenvalue, reached at the corresponding eigenvector. Repeating subject to orthogonality to previous directions gives the other principal components.
The software uses the divisor- sample covariance matrix for its component standard deviations, whereas var used divisor . Diagonalizing the displayed matrix gives approximately
Taking square roots of the latter gives , agreeing with the output to the rounding of the displayed sample covariance matrix. The normalization changes these numerical standard deviations but leaves the eigenvectors and proportions of explained variance of a principal component unchanged.
The first principal component explains about of the total log-variable variance, the second another , and the last about . Keeping two preserves about . An original scree plot displays the sharp drop after the first component:
Figure 1.
Principal-component variances and cumulative explained variance for the log-tree measurements
. Details.
Consequently the centered cloud is very elongated along one direction: a one-dimensional approximation represents almost all of its squared variation. From the displayed covariance matrix, a unit first eigenvector, with its arbitrary sign chosen positive, is approximately
Its principal component score is therefore . All coefficients have the same sign, so this direction describes overall tree size, with log volume dominating the unstandardized variation. The summary alone supplies proportions, not the loading directions; the displayed sample covariance matrix supplies the latter. This is a geometric approximation, not evidence of an exact deterministic relation or independent original variables.
One further use is principal component analysis on a correlation matrix, obtained with cor=T, which equalizes marginal variances and answers a different scaling-sensitive question. Another is principal component regression: fit a response to selected principal component scores rather than highly correlated original predictors, choosing the number using prediction assessment. The fitted transformation also projects new measurements into the same low-dimensional coordinates, provided the training centering and scaling are retained.
princomp(): principal component analysis. Supply observations in rows and numeric variables in columns. After centering, let denote the resulting matrix and its sample covariance matrix, with an explicitly stated normalization. Write
The first direction maximizes the Rayleigh quotient over unit vectors; successive directions maximize it subject to orthogonality to the earlier directions. The principal component scores are , with variances and zero pairwise sample covariance. Zero covariance alone does not imply independence. Each principal component loading describes a direction in variable space; its sign is arbitrary, and repeated eigenvalues leave the basis inside the corresponding eigenspace arbitrary.
The function returns principal component loadings, component standard deviations sdev and, when requested, principal component scores. By default it uses a covariance matrix; cor=T instead performs principal component analysis on a correlation matrix, standardizing nonconstant variables to variance one. This choice matters when variables have different measurement units. These interface details are documented in the S-Plus help for princomp.
A scree plot shows the ordered eigenvalues; a plot of the first two principal component scores can reveal groups, unusual observations or a dominant gradient. The first components explain of the total variance. Reconstructing each centered observation by its orthogonal projection onto these directions minimizes average squared reconstruction error among -dimensional linear subspaces. Thus PCA gives an optimal linear compression for squared reconstruction error, with scaling and the retained dimension chosen explicitly. It is an unsupervised learning method, rather than a classification rule fitted to supplied labels.
tree(): classification trees and regression trees. A call such as tree(response ~ predictors, data=dat) fits a recursive partitioning model. At a current node, test numerical predictors against candidate thresholds, or divide the levels of a categorical predictor into two sets. Select a split that most improves the fitting criterion, then repeat within its children. Numerical splits give axis-aligned regions. A categorical response produces a classification tree; a numeric response produces a regression tree, as specified in the S-Plus tree documentation.
For a regression tree, the constant prediction in a terminal region is its sample mean, minimizing the sum of squared residuals there. For a classification tree, estimate class probabilities by the observed proportions in a leaf and predict a most probable class under equal misclassification costs. The classification-tree deviance of a node is
with ; splitting maximizes the decrease from parent to children. This fitting criterion differs from the number of wrong classifications. The middle sketches show the same three-region partition as a binary tree: first split on variable one, then split the right-hand branch on variable two.
Growing to small leaves can overfit. Use cost-complexity tree pruning, balancing fitted loss against leaf count,
and select complexity by cross-validation. Plotting the fitted tree and annotating nodes displays the rules; prediction routes each new observation through the fitted splits. A tree gives an interpretable, piecewise-constant prediction rule; training fit alone does not measure prediction accuracy.
The rpart() alternative. This also fits classification and regression trees by recursive partitioning, with squared-error splitting for ordinary numeric responses and, by default for classification, Gini impurity . Classification settings can instead use information-based splitting and specify prior probabilities and a loss matrix. It supports surrogate splits to approximate a node's primary left/right assignment when the primary predictor is missing. Its complexity-parameter path and cross-validated errors support cost-complexity tree pruning; for example, select a complexity parameter and apply prune(fit, cp=...). The maintainer's rpart documentation describes the criteria and the distinction in handling surrogate predictors. The partition and tree sketches apply to either implementation; the fitted tree need not be the same because splitting criteria and controls can differ.
cmdscale(): classical multidimensional scaling. Its input is a dissimilarity matrix or stored distance object. It outputs coordinates for the observations in a chosen dimension, usually two for visualization. Let contain the squared entries of , and put
For centered Euclidean coordinates , the distance identity gives : multiplying by removes the row and column norm terms. Diagonalize this Gram matrix and retain its largest positive eigenvalues. The coordinates are
Thus classical scaling exactly recovers Euclidean distances when all positive eigenvalues are retained. Reducing dimension approximates the centered Gram matrix; it does not in general minimize squared discrepancies in the original distances. Negative eigenvalues identify a failure of the Euclidean distance matrix criterion, so their size matters when interpreting a display. The historical S-Plus cmdscale help describes the dimension argument, eigenvalue output and optional additive adjustment to dissimilarities.
Plot the two coordinate columns on equally scaled axes and label observations. Nearby plotted points suggest small dissimilarities only insofar as the dimension reduction fits the input well. The bottom sketches give an exact Euclidean example: a configuration and its independently recovered two-dimensional coordinates have identical pairwise distances, despite different coordinate axes. Classical scaling converts distances into geometry; translations, rotations and reflections do not change the represented distances. For Euclidean distances between observations, its coordinates agree with the corresponding principal component scores up to these harmless choices of axes.
Figure 1.
Illustrative PCA directions and scree plot, equivalent recursive partition and tree, and exact classical-scaling recovery of Euclidean distances
. Details.
Factor analysis models a large set of correlated measurements using a smaller collection of shared latent factors and variable-specific variation. A classical orthogonal factor model is
Take , centered errors, and positive specific variances for the inverse-covariance and likelihood formulas below. The factor loadings form the matrix , and the model covariance is
Thus for , ; the common latent factors explain off-diagonal correlation. The communality of variable is , and . For standardized variables this is . The specific variance represents measurement noise or genuinely unshared variation; it is not another common latent factor.
This differs from principal component analysis. A principal component maximizes projected total variance and is a determined linear score after choosing its direction; a factor model attempts to explain shared covariance while reserving diagonal specific variance. Low-dimensional PCA truncation can help initialize a factor fit, but it does not justify setting all specific variances to zero. Latent factors are unobserved, not directly observed scores. When the latent factors and the errors are jointly Gaussian, their conditional mean is
which supplies one factor-score estimate, with remaining uncertainty. This follows from the block conditional normal formula in Question 1; it does not make the latent factor uniquely recoverable from one observation.
For the Gaussian likelihood, use the empirical covariance . Fit the factor model by minimizing subject to admissible specific variances. Other approaches estimate communalities and extract latent factors from the reduced correlation matrix. Choose using substantive plausibility, a scree plot, fitted residual correlations and model-fit assessment. Under local identifiability, a likelihood goodness-of-fit count for the covariance model is
where the subtracted rotational dimension must be removed from the loading parameter count. This count alone does not establish identifiability, and boundary fits with zero specific variance need special care.
Factor rotation addresses interpretability and nonuniqueness. For every orthogonal matrix , replace by ; then , so the fit is unchanged. The absolute orientation of unrestrained latent factors is therefore not identified by the covariance. Varimax rotation chooses orthogonal axes favoring a simple pattern of large and small squared loadings. Oblique rotation can permit correlated latent factors, provided their covariance is transformed consistently. Sign changes and permutations are further harmless relabelings; an interpretation must state its chosen convention.
Figure 1.
Two-factor loading directions, an equivalent rotated representation and the covariance scree plot
. Details.
The sketch uses four standardized hypothetical measurements sharing two orthogonal latent factors. Rotating the factor coordinates changes the loading descriptions while preserving every implied covariance and communality; the scree plot belongs to that same constructed covariance matrix. Factor analysis separates shared variation from specific variance, and its latent factors require an identifying or interpretive orientation.