Source-linked AI summary
Local Gaussian process approximation for large computer experiments
Robert B. Gramacy, Daniel W. Apley
TL;DR
Large computer-experiment emulation is limited by Gaussian-process computation on large datasets and by potentially inappropriate stationarity assumptions. The paper builds local sequential GP predictors with fast updates and independently applies them across prediction locations. The resulting predictor is fast, nonstationary, highly parallelizable, and competitive with related approaches, with a trade-off between local design size and flexibility.
Problem
Gaussian-process emulation offers strong prediction but becomes computationally expensive for large datasets, while stationarity may be inappropriate.
Method
The paper constructs local GP supports through sequential design, updates quantities efficiently, and performs local inference independently across prediction locations.
Results
The resulting predictor is fast, nonstationary, highly parallelizable, and competitive with related approaches, including compactly supported covariance methods.
Takeaways & Limitations
Local design size directly controls the speed–accuracy trade-off, while independent local predictions enable near-linear multicore speedups.
Takeaways & Limitations
Larger local designs improve approximation to the full-data GP but can reduce global spatial heterogeneity and worsen out-of-sample performance.
Abstract
from arXiv · showhide
We provide a new approach to approximate emulation of large computer experiments. By focusing expressly on desirable properties of the predictive equations, we derive a family of local sequential design schemes that dynamically define the support of a Gaussian process predictor based on a local subset of the data. We further derive expressions for fast sequential updating of all needed quantities as the local designs are built-up iteratively. Then we show how independent application of our local design strategy across the elements of a vast predictive grid facilitates a trivially parallel implementation. The end result is a global predictor able to take advantage of modern multicore architectures, while at the same time allowing for a nonstationary modeling feature as a bonus. We demonstrate our method on two examples utilizing designs sized in the thousands, and tens of thousands of data points. Comparisons are made to the method of compactly supported covariances.
1 Introduction
Gaussian processes offer strong emulation performance but become computationally difficult at large data sizes, while stationarity may be inappropriate. The paper develops local, sequentially constructed predictors that support scalable, nonstationary, and parallel emulation.
- Gaussian processes provide strong out-of-sample prediction and coverage but scale poorly computationally as datasets grow.
- Existing approximations reduce stationary-model computational costs through pseudo-inputs, iterative estimators, fixed-rank kriging, or compactly supported covariances.
- The method selects influential data points near each prediction input instead of forming or using full covariance matrices.
- Localized sub-designs produce more accurate predictions than local nearest neighbors at fixed design size and leverage global nonstationarity to compensate for simple local structure.
- Unlike compactly supported covariance methods, the approach is purely local, explicitly links accuracy to computational budget, and supports independent parallel inference across prediction locations.
- The paper derives greedy local-design criteria, efficient updates, global prediction procedures, and empirical comparisons with compactly supported covariances.
2 Gaussian process predictive modeling
The paper formulates computer-model emulation with Gaussian-process means and covariances, commonly using zero mean and correlation functions based on Euclidean distance. It estimates covariance parameters through Gaussian likelihood methods and obtains Student-t predictive distributions.
- A Gaussian process is defined by its mean and covariance, with the covariance commonly separated into variance and correlation components.
- The exposition uses a zero-mean convention and an isotropic one-parameter correlation family, while the main derivations assume a separable family.
- Under data D=(X,Y), responses follow a multivariate normal likelihood with covariance Σ(X)=τ^2K, enabling analytic estimation of τ^2 and profile-likelihood inference for other parameters.
- An inverse-gamma reference prior allows analytic integration over τ^2 in the empirical-Bayes formulation.
- Newton-like methods estimate the correlation parameter effectively when analytic derivatives are available and the likelihood surface is not multimodal.
- The marginalized GP predictive distribution is Student-t, with predictive variance given by σ2(x|D,K) × N/(N−2).
3 Localized approximate emulation
The paper replaces fixed nearest-neighbor GP support with locally optimized, greedily constructed sub-designs for fast approximate prediction. Sequential updates make the approach computationally practical, while empirical results show improved accuracy and scalability relative to nearest neighbors and full-GP computation.
- Motivation: Nearest-neighbor sub-designs are computationally sensible but can be suboptimal for fixed local size because the best predictive support depends on θ.Finding the best larger sub-design is a difficult high-dimensional nonconvex optimization problem.
- Greedy local design: For each predictive location x, the method greedily adds design sites by minimizing an empirical Bayes MSPE criterion that accounts for predictive variance and uncertainty in θ.The sequential design begins from an initial sub-design and evaluates each candidate using the current local data.
- Greedy local design: The MSPE approximation combines known-parameter variance reduction, parameter-learning information, Gaussian moment matching, and the current rather than future estimate of θ.These approximations avoid analytically intractable averaging over future responses and parameter estimates.
- Fast updates: Required quantities admit sequential updates, with each iteration costing O(j^2), n iterations costing O(n^3), and one final local MLE update also costing O(n^3).Predictive quantities and likelihood components are updated using partition-inverse relationships and fixed-θ recursions.
- ALC special case: A local ALC criterion based on the variance-reduction integrand avoids the costly numerical integration required by global ALC and simplifies derivative calculations.Numerical integration can otherwise add significant computational burden.
- Empirical illustration: After 50 local points, greedy/MSPE predictions were close to the truth with lower variability than nearest-neighbor predictions, which remained biased low and more confident.Both methods acknowledged uncertainty from using a much smaller local design than the full data set.
4 Global emulation on a dense design
The method builds independent local GP approximations at each predictive location, enabling parallel global emulation while iteratively updating local designs and correlation estimates. In an illustrative example, the resulting mean closely tracks truth, but localized approximation errors and uncertainty reveal limitations under a fixed stationary covariance.
- Global prediction: Serial prediction over T locations costs O(Tn^3), whereas full GP prediction requires O(TN^2 + N^3).The local formulation avoids constructing full predictive equations at every location.
- Global prediction: Independent local designs permit trivial parallelization across predictive locations.The authors report nearly linear speedup using an OpenMP-parallelized loop over x.
- Illustrative behavior: Figure 3 shows nearly indistinguishable mean surfaces, small spatially varying errors, and non-uniform predictive standard deviations across the input space.Lower predictive variability corresponds to greater inaccuracy, suggesting difficulty with the assumed stationary covariance.
- Stagewise local inference: The iterative scheme independently builds local designs, estimates local lengthscales by MLE, optionally smooths them, and repeats before producing predictions.Local designs use MSPE or ALC heuristics, while local estimates yield a globally nonstationary predictive surface.
- Stagewise local inference: Usually one repetition suffices for joint convergence, although pathological initial lengthscales can slow convergence.A sensible initial value can be chosen from lower quantiles of squared distances in X.
- Stagewise local inference: Smoothing can mitigate multimodal likelihood mode-switching and encourage slowly varying correlation parameters, but is not essential for good predictions.The authors suggest smoothing may aid convergence of the two-stage scheme.
5 Empirical demonstration and comparison
Experiments compare local approximation variants with nearest-neighbor and compactly supported covariance methods on computer-experiment benchmarks. Greedy local methods are generally more accurate and scalable, while local inference and parallelization are important for the larger comparisons.
- Illustrative comparison: Greedy methods are generally more accurate than the worst alternatives, with greater computational effort typically improving performance.All methods show high precision and over-cover in the summarized illustrative comparison.
- Illustrative comparison: Methods without local lengthscale inference perform worst, while nearest neighbors become competitive only with an order-of-magnitude larger local design.The larger NN output also requires storing a 4x larger object.
- Borehole comparison: On the borehole benchmark, local methods without local lengthscale inference fare worst, and NN with n = 50 is dominated by greedy methods across the board.A big NN with n = 200 is competitive per second only when local MLEs are calculated.
- Uncertainty comparison: Greedy methods over-cover, whereas CSC methods achieve nominal pointwise coverage; the local methods are described as more accurate and conservative.The low standard deviation of nnbig is attributed to its larger local-design denominator rather than better global variance estimation.
- Scaling comparison: Doubling the experiment size leaves relative accuracy orders similar, while local-method runtimes scale approximately linearly and CSC runtimes increase four to eight times.The authors report that a 10,000-point CSC run was estimated to take days because of memory-related slowdowns, whereas local methods required only seconds more.
6 Discussion
The discussion presents local sequential GP designs as fast, parallelizable, and tunable approximations, while emphasizing trade-offs between computation, accuracy, and nonstationary flexibility. It also identifies resource-allocation extensions and settings where competing mean structures may be preferable.
- Contributions: The method modernizes local-neighborhood kriging into a fast, nonstationary, highly parallelizable predictor with local design size n controlling speed versus accuracy.The authors note that greedy methods share NN’s local O(n3) order but have larger constants.
- Extensions: Local designs can be augmented incrementally from stored indices and parameter estimates, with larger n(x) allocated where variance reduction or accuracy gains are greater.The proposed allocation can use active-learning-based calculations to prioritize predictive locations.
- Extensions: Stopping rules based on predictive-to-reduced variance ratios could allocate computational effort differentially across predictive locations.The ratio Vj(x)/Vj+1(x) is proposed as a response-independent metric after cancellation of shared terms.
- Limitations: Nonuniform designs and spatially varying parameters make convergence time highly unpredictable across thousands of independent local searches.This unpredictability may require a low global cap on nthresh(x).
- Limitations: Larger local designs improve approximation to the full-data GP but reduce global spatial heterogeneity and may worsen out-of-sample performance.The discussion suggests early stopping or post-process pruning to balance these effects.
- Extensions: For dense designs and predictive grids with globally shared parameters, interior greedy designs are highly regular, motivating cheaper template or restricted-candidate searches.The experiments restricted searches to the nearest N′=1000 candidates, with results except time also evaluated at N′=10000.
- Extensions: An early GPU implementation of ALC-based local design search yielded 20–70x efficiency gains.The authors identify GPU thousands-core architectures as a potential source of further parallelization benefits.
- Extensions: The local approach may also support computer-model calibration when predictions are needed only over part of the input space.The authors suggest jointly building local designs for field-data locations, while noting that MSPE details remain under development.
Supplementary Materials
The supplementary materials derive analytic and approximate updates for local GP design criteria, then document numerical stability and computational comparisons with nearest-neighbor and compactly supported covariance methods.
- Derivative and information updates: The supplement derives likelihood derivatives, Fisher information updates, and predictive-equation derivatives for sequential local GP design.The derivations allow analytically differentiable correlation functions and extend beyond the scalar-parameter exposition to separable multi-parameter families.
- Local sequential design: The Bayesian MSPE criterion selects the next local design site by balancing reduced predictive variance with uncertainty in the correlation parameters.When future observations make exact expectations intractable, the supplement uses Gaussian surrogate equations with matched moments and additional derivative approximations.
- Approximation caveat: The MSPE construction deliberately approximates future-parameter and future-prediction derivatives because their dependence on the unobserved response prevents tractable analytic expectations.For larger local designs, the current and expanded predictive means and derivatives are expected to become similar because influential points are selected early.
- Numerical stability: Greedy sub-designs produce lower covariance-matrix condition numbers than nearest-neighbor designs as local design size increases, improving numerical stability and permitting smaller nuggets.The comparison uses two nugget values; lower condition numbers imply more reliable matrix decompositions.
- Empirical comparisons: On a doubled borehole experiment, accuracy remained similar while compactly supported covariance comparators required almost 10× the runtime and greedy local runtimes scaled linearly.The tables report a Nash–Sutcliffe efficiency variant, with smaller 1−NSE values indicating better accuracy.