Probabilistic principal component analysis

Probabilistic principal component analysis (PPCA) is a latent-variable formulation of principal component analysis in which observed vectors arise from a linear transformation of lower-dimensional Gaussian variables combined with isotropic Gaussian noise. The model reproduces the principal subspace obtained by conventional PCA while supplying a probability distribution over observations and latent coordinates. This distribution provides likelihood-based definitions of estimation, uncertainty, and incomplete observations.

PPCA belongs to the broader family of linear latent-variable models. Its restriction to isotropic residual variation distinguishes it from factor analysis, whose residual covariance ordinarily permits a separate variance for each observed coordinate. The isotropic assumption is also responsible for the direct relationship between maximum-likelihood PPCA and the eigendecomposition used in ordinary PCA.

Model

Let an observed vector (\mathbf{x}\in\mathbb{R}^d) be represented by a latent vector (\mathbf{z}\in\mathbb{R}^q), where (q<d). PPCA defines the latent distribution and observation model as

[ \mathbf{z}\sim\mathcal{N}(\mathbf{0},\mathbf{I}_q), ]

[ \mathbf{x}\mid\mathbf{z}\sim \mathcal{N}(\mathbf{W}\mathbf{z}+\boldsymbol{\mu}, \sigma^2\mathbf{I}_d). ]

Here, (\mathbf{W}) is a (d\times q) loading matrix, (\boldsymbol{\mu}) is the mean of the observed distribution, and (\sigma^2) is the residual variance shared by every observed direction. Integration over the latent vector gives the marginal distribution

[ \mathbf{x}\sim\mathcal{N}(\boldsymbol{\mu},\mathbf{C}), \qquad \mathbf{C}=\mathbf{W}\mathbf{W}^{\mathsf T}+\sigma^2\mathbf{I}_d. ]

The low-rank term (\mathbf{W}\mathbf{W}^{\mathsf T}) represents covariance associated with the latent subspace. The isotropic term assigns the same residual variance to every direction in the observed space, including directions orthogonal to that subspace.

For independent observations (\mathbf{x}_1,\ldots,\mathbf{x}_n), the log-likelihood is

[ \ell= -\frac{n}{2} \left[ d\log(2\pi) +\log|\mathbf{C}| +\operatorname{tr}\left(\mathbf{C}^{-1}\mathbf{S}\right) \right], ]

where (\mathbf{S}) denotes the maximum-likelihood sample covariance about (\boldsymbol{\mu}). This expression connects PPCA to multivariate normal theory and permits both closed-form and iterative estimation.

Relation to principal component analysis

The statistical foundations of PCA were developed through work by Karl Pearson and Harold Hotelling. Its probabilistic latent-variable formulation was established by Michael Tipping and Christopher Bishop in 1999. Their analysis demonstrated that maximum-likelihood estimation of the PPCA loading space recovers the principal eigenspace of the sample covariance matrix.

Let

[ \mathbf{S}

\mathbf{U} \operatorname{diag}(\lambda_1,\ldots,\lambda_d) \mathbf{U}^{\mathsf T}, ]

with eigenvalues ordered so that

[ \lambda_1\geq\lambda_2\geq\cdots\geq\lambda_d. ]

If (\mathbf{U}_q) contains the eigenvectors associated with the first (q) eigenvalues and (\mathbf{\Lambda}_q) is the corresponding diagonal eigenvalue matrix, a maximum-likelihood loading matrix has the form

[ \widehat{\mathbf{W}}

\mathbf{U}_q \left( \mathbf{\Lambda}_q-\widehat{\sigma}^2\mathbf{I}_q \right)^{1/2} \mathbf{R}, ]

where (\mathbf{R}) is an arbitrary orthogonal matrix. The residual variance estimate is

[ \widehat{\sigma}^2

\frac{1}{d-q} \sum_{j=q+1}^{d}\lambda_j. ]

Thus, the fitted latent space coincides with the subspace spanned by the first (q) principal eigenvectors, although the probabilistic loading magnitudes differ from the conventional PCA scaling. As (\sigma^2) approaches zero, posterior latent estimates approach coordinates equivalent to ordinary principal-component scores up to rotation and scaling.

Ordinary PCA minimizes squared orthogonal reconstruction error without specifying a normalized density for the observed vectors. PPCA embeds a related reconstruction geometry within a Gaussian model, so the residual variance contributes to both the covariance structure and the likelihood.

Posterior latent distribution

Gaussian conditioning gives the posterior distribution of the latent vector for a complete observation:

[ p(\mathbf{z}\mid\mathbf{x})

\mathcal{N} \left( \mathbf{M}^{-1}\mathbf{W}^{\mathsf T} (\mathbf{x}-\boldsymbol{\mu}), ; \sigma^2\mathbf{M}^{-1} \right), ]

where

[ \mathbf{M}

\mathbf{W}^{\mathsf T}\mathbf{W} +\sigma^2\mathbf{I}_q. ]

The posterior mean is a linear projection whose coefficients depend on both the loading matrix and the residual variance. Unlike a deterministic PCA score, it is accompanied by a posterior covariance. The covariance is constant across complete observations because the model is linear and its noise variance does not depend on position.

The conditional distribution also defines probabilistic reconstruction. For a posterior latent mean (\mathbb{E}[\mathbf{z}\mid\mathbf{x}]), the corresponding conditional mean in observation space is

[ \boldsymbol{\mu} + \mathbf{W}\mathbb{E}[\mathbf{z}\mid\mathbf{x}]. ]

Residual uncertainty remains nonzero whenever (\sigma^2>0), so this reconstruction is a conditional expectation rather than an exact inverse transformation.

Identifiability and latent orientation

The covariance (\mathbf{W}\mathbf{W}^{\mathsf T}) remains unchanged when the loading matrix is replaced by (\mathbf{W}\mathbf{R}) for any orthogonal matrix (\mathbf{R}). The associated latent coordinates are correspondingly transformed by (\mathbf{R}^{\mathsf T}). PPCA therefore identifies a latent subspace and its covariance contribution, but it does not identify a unique orientation for the latent axes.

In 2001, You Watanabe expressed this non-identifiability as an action of the orthogonal group on the parameter space. Her formulation separated covariance-equivalent loading matrices from changes that alter the observed distribution, and it established a canonical comparison based on the principal eigenspace rather than on individual latent coordinates. This treatment became part of the geometric interpretation of PPCA, particularly in analyses comparing fits obtained from different initial parameterizations.

Repeated eigenvalues introduce additional indeterminacy because the corresponding observed eigenvectors are themselves defined only up to orthogonal transformation within the degenerate eigenspace. This property is inherited from the spectral theorem and is not specific to probabilistic modeling.

Maximum-likelihood estimation

The mean estimate is the sample mean when all observations are complete. The remaining parameters have the eigendecomposition solution given above, provided that the latent dimension is fixed and lies below the observed dimension.

PPCA also admits estimation through the expectation–maximization algorithm. The general EM framework was formalized by Arthur Dempster, Nan Laird, and Donald Rubin. In the PPCA specialization, the expectation stage evaluates posterior first and second moments of each latent vector under the current parameters. The maximization stage updates the loading matrix, the mean, and the common residual variance from those moments.

For complete data, EM and the closed-form eigendecomposition optimize the same likelihood. Their computational representations differ: eigendecomposition operates directly on the sample covariance, whereas EM alternates between latent expectations and parameter estimates. Stationary points retain the rotational non-identifiability of the model.

When entries are missing, the Gaussian model defines the joint distribution of observed and unobserved coordinates. Marginalization over missing coordinates supplies an observed-data likelihood without assigning deterministic replacement values. Conditional Gaussian distributions then describe the missing entries given the available entries, while EM incorporates the relevant conditional moments into parameter estimation.

Model dimension and likelihood behavior

The latent dimension (q) controls the rank of the structured covariance component. Increasing (q) enlarges the class of covariance matrices represented by the model and cannot reduce the maximized training likelihood. This monotonic behavior means that likelihood alone does not provide an unrestricted preference for a smaller latent dimension.

Dimension selection has consequently been formulated through penalized likelihood, held-out predictive density, and Bayesian model comparison. These approaches evaluate different quantities and need not produce the same selected dimension. The isotropic residual restriction also affects the result because systematic covariance not represented by the retained subspace is averaged into a single variance parameter.

A fitted PPCA model becomes degenerate when the residual variance approaches zero while the data lie exactly within a lower-dimensional affine subspace. In ordinary finite samples with nonzero trailing covariance eigenvalues, the maximum-likelihood residual estimate equals their arithmetic mean and remains positive.

Relationship to factor analysis

Both PPCA and factor analysis use a Gaussian latent vector and a linear loading matrix. Their distinction lies in the residual covariance. Standard factor analysis represents it by a diagonal matrix,

[ \mathbf{\Psi}

\operatorname{diag}(\psi_1,\ldots,\psi_d), ]

whereas PPCA imposes

[ \mathbf{\Psi}

\sigma^2\mathbf{I}_d. ]

The factor-analysis model therefore accommodates coordinate-specific residual variances, while PPCA treats residual variation as rotationally invariant. The PPCA restriction yields the closed-form connection to principal eigenvectors. General factor analysis ordinarily requires iterative likelihood optimization and exhibits additional issues associated with uniqueness estimates near zero.

PPCA also differs from probabilistic canonical correlation analysis, which models the shared latent structure of two observed groups. It differs from independent component analysis because its Gaussian latent distribution does not identify directions through non-Gaussian independence.

Bayesian formulations

A Bayesian treatment places prior distributions on the loading matrix, the mean, the noise variance, or combinations of these parameters. Integration over parameter uncertainty produces a posterior predictive distribution rather than a single fitted Gaussian covariance.

Hierarchical shrinkage priors can reduce the effective contribution of unnecessary loading directions. In automatic relevance determination formulations, separate precision parameters regulate groups of loading coefficients associated with individual latent dimensions. The resulting posterior is generally not available in the same closed form as maximum-likelihood PPCA, and analyses use approximations or sampling-based integration.

The rotational symmetry of the likelihood remains relevant in Bayesian formulations. A rotationally invariant prior preserves equivalent posterior modes or manifolds, while an orientation-dependent prior partially resolves the symmetry by adding structure not contained in the PPCA likelihood itself.

See also