Polynomial Chaos Expansion (PCE)

Theoretical Foundations of PCE

PCE as a Spectral Method in Probability Space

Polynomial chaos expansion (PCE) represents a response with random inputs as a series in orthogonal polynomials. Wiener formulated it in 1938 for Gaussian variables and Hermite polynomials, and the generalisation to arbitrary distributions (the Wiener-Askey scheme: uniform ↔ Legendre, gamma ↔ Laguerre, and so on) opened it up to engineering. At heart it is a spectral method in probability space — where a Fourier series expands the spatial direction, orthogonal polynomials expand the stochastic direction. If the response is a smooth function of the inputs, the coefficients decay quickly and a handful of terms already gives high accuracy (spectral convergence). That character — fast when smooth, painful when discontinuous — carries over intact, so ordinary spectral-method intuition applies directly.

What Orthogonality Buys You

Because the basis \( \Psi_{\boldsymbol{\alpha}} \) is orthogonal under the weight of the input distribution (\( \mathbb{E}[\Psi_i \Psi_j] = \delta_{ij}\gamma_i \)), statistics turn into algebra the moment the coefficients are known.

$$ f \approx \sum_{\boldsymbol{\alpha}} c_{\boldsymbol{\alpha}} \Psi_{\boldsymbol{\alpha}}(\mathbf{X}) \;\Rightarrow\; \mathbb{E}[f] = c_{\mathbf{0}}, \quad \mathrm{Var}[f] = \sum_{\boldsymbol{\alpha}\ne\mathbf{0}} c_{\boldsymbol{\alpha}}^2 \gamma_{\boldsymbol{\alpha}} $$

Sobol' sensitivity indices likewise follow exactly, as partial sums of the coefficients. Statistics and sensitivities are read off rather than waited for through Monte Carlo convergence (\( 1/\sqrt{N} \)) — this is the strength that sets PCE apart from other surrogates (for the practical workflow see sparse PCE and LAR).

Intrusive and Non-Intrusive — Two Lineages

This article deals mainly with the design of non-intrusive schemes. Intrusive PCE is an option in a research context — in-house codes and dedicated UQ frameworks.

Numerical Routes to the Coefficients

Projection — Quadrature Design Is Everything

Orthogonality gives each coefficient as the inner product \( c_{\boldsymbol{\alpha}} = \mathbb{E}[f\,\Psi_{\boldsymbol{\alpha}}]/\gamma_{\boldsymbol{\alpha}} \), and evaluating it by numerical quadrature is the projection method. In one dimension, Gaussian quadrature (with points matched to the distribution) delivers the highest accuracy for the fewest points. In several dimensions a full tensor product explodes as \( n^d \), so a Smolyak sparse grid thins out the quadrature of the cross terms. Projection has a distinctive character: you do not choose the sample points (the quadrature rule does), and their number follows from order and dimension. For \( d \lesssim 10 \) with a smooth response it is still the most dependable route. Two cautions apply: (1) sparse grids carry negative weights and are sensitive to non-smooth responses, and (2) quadrature points often land at the edges of the range, so designing input ranges whose edges do not break the analysis matters even more here than with regression.

When to Use Regression Instead

Regression (least squares, sparse regression) places samples freely, which is what lets it combine with discarding failed runs, adding samples sequentially, and sparsification. Rough guidance:

• Low dimension (d≤5), smooth response, stable solver → projection (Gaussian quadrature or sparse grid). Deterministic and highly reproducible
• Moderate dimension, some failed runs, sparsity worth exploiting → regression (LAR with corrected LOO)
• High dimension (d>20) → cut the dimension down first with screening

Verifying Convergence — Never Skip the Order Study

Verification of a PCE rests on an order study: raise the order \( p \) and check that the statistics settle. The mean is fixed at low order, the variance lags a little, and tail probabilities (reliability analysis) demand higher order still — every output statistic has its own order requirement, so confirm convergence directly on the quantity you actually need. Report pointwise error at hold-out points (or LOO) alongside it and the result takes exactly the form of a mesh convergence study.

Guidance for Practical Use

Where PCE Wins and Where It Loses

Input Modelling Dominates the Result

Everything PCE outputs — mean, variance, sensitivities, probabilities — is an answer conditional on the assumed input distributions. Distribution family (normal or lognormal), tail weight, and the presence of correlation change the result substantially, and that uncertainty does not disappear no matter how accurate the expansion becomes. The discipline in practice is threefold: (1) record the evidence behind each distribution (measurement, standard, literature), (2) run doubtful distribution assumptions as several cases and look at the sensitivity, and (3) treat correlation explicitly through a Nataf transform or similar, since ignoring it distorts variance and sensitivities systematically. Validating the surrogate and validating the input model are separate exercises, and the latter drives the result more — it is the first thing to look at when reviewing a UQ report.

From Scalars to Vector and Field Outputs

When the response is a field (a stress distribution, a temperature field), the standard approach is not to fit an independent PCE at every node but to compress the field into a few POD mode coefficients first and build a PCE per coefficient — the point where this meets ROM. Statistics of the whole field (mean field, variance field, confidence intervals at any point) then come out of a single consistent framework.

Tools and Implementation

Tool Options

A Minimal Implementation Sketch (ChaosPy, Projection)

import chaospy as cp
import numpy as np

dist = cp.J(cp.Normal(210e9, 6e9), cp.Uniform(0.28, 0.32))  # E, Poisson's ratio
expansion = cp.generate_expansion(4, dist)                   # order-4 orthogonal basis
nodes, weights = cp.generate_quadrature(5, dist, rule="gaussian")
evals = [run_fem(*n) for n in nodes.T]                       # run CAE at the quadrature points
model = cp.fit_quadrature(expansion, nodes, weights, evals)  # coefficients by projection
print(cp.E(model, dist), cp.Std(model, dist))                # mean, standard deviation
print(cp.Sens_m(model, dist))                                # first-order Sobol' indices

The regression version swaps in fit_regression and reuses the same post-processing. Being able to step through projection, then regression, then sparse regression and compare them is the educational value of a script-based implementation.

Research Frontiers

Stochastic Galerkin Reconsidered

Long dismissed as "beautiful theory that loses to non-intrusive methods on implementation cost", the intrusive route is being reconsidered as differentiable programming and automatic code transformation mature. If automatic solver rewriting lowers the implementation barrier, intrusive PCE could take the lead on long-time accuracy in time-dependent problems — precisely the weak point of non-intrusive methods — and it continues as one strand of UQ infrastructure research.

Arbitrary Distributions, Data-Driven Bases (aPC), and Bayesian Integration

Arbitrary PC builds orthogonal bases numerically from the moments of measured data; Bayesian coefficient estimation carries the uncertainty of the PCE itself (arising from small samples) as a credible interval; and evidence-based model comparison handles the choice of order and basis. The field is deepening in a statistically honest direction. Combining Bayesian calibration with a PCE surrogate is where this trend touches practice.

Long-Time Accuracy in Time-Dependent Problems

Applying PCE to unsteady responses (vibration, transients) runs into the classic "long-time problem": nonlinearity in the stochastic space grows with time and the order requirement diverges. Time-interval splitting (time-adaptive PCE), separation of phase and amplitude, and flow-map-based methods are all under study, and application to UQ of dynamic systems (seismic response, flutter) is spreading. In practice, the single most effective way to dodge the difficulty is to narrow the question down: which quantity, at which instant, is the UQ actually about?

Troubleshooting

Symptoms, Causes, and Fixes

Getting the Hang of It

Related: sparse PCE and LAR (the practical workflow), Kriging surrogates, reduced-order models (ROM).