Model-Specific to Model-General Uncertainty for Physical Properties

uncertainty-quantification
gaussian-process-regression
bayesian-inference
equation-of-state
machine-learning
Industrial & Engineering Chemistry Research (2022)
Published

2022

Doi

Ni Zhan and John R. Kitchin, “Model-Specific to Model-General Uncertainty for Physical Properties”, Industrial & Engineering Chemistry Research 61(24), 8368-8377 (2022). https://doi.org/10.1021/acs.iecr.1c04706

Read the paper

This summary was written by a language model from the paper and has not been checked line by line. Rely on the published paper for details.

Question

Physical properties such as rate constants, diffusivities or elastic moduli are usually extracted by fitting a model to data, and the reported value should carry an uncertainty. Uncertainty can be split into aleatoric (noise in observations) and epistemic (ignorance about the best model) parts, and it can be computed by frequentist confidence intervals or by Bayesian posterior distributions. The paper asks what a given uncertainty estimate actually captures: does it describe only the spread in one model’s parameters, or also the uncertainty in choosing the model form? (Zhan 2022)

Approach

The test problem is the energy-volume equation of state (EOS) of a solid, where the equilibrium volume and energy come from the function and its minimum, and the bulk modulus comes from the second derivative (B = V d2E/dV2), so derivative-based properties have higher variance. DFT energy-volume data for fcc Pd and Au were taken from earlier group work (Boes and co-workers, 2016 and 2017), restricted to a window near the minimum (24 points for Pd, 29 for Au) to keep residuals homoscedastic across models. Seven analytical EOS forms were fitted: stabilized jellium, Anton-Schmidt, a third-order polynomial (Taylor), Murnaghan, Birch-Murnaghan, Poirier-Tarantola and Vinet.

Three uncertainty methods were compared. First, nonlinear least squares with the delta method, following the authors’ earlier AIChE Journal work, using the inverse Fisher information for properties that are direct parameters and the property’s derivative with respect to parameters for derived quantities. Second, Bayesian nonlinear regression framed as a probabilistic graphical model in pyro, with normal priors on parameters (centered near the regression solution) and a normal residual whose standard deviation has a Uniform(0, 1) prior; the posterior was approximated by stochastic variational inference with a multivariate normal guide and by Hamiltonian Monte Carlo with the No-U-Turn sampler. Third, a Gaussian process (GP) with an RBF kernel whose joint covariance spans the function and its first and second derivatives, built from kernel derivatives computed with jax automatic differentiation and vectorized mapping, with hyperparameters fitted in gpytorch by maximizing the marginal likelihood. Property distributions from the GP came from 1000 joint samples of function, slope and curvature on a 1000-point volume grid, locating the minimum with a Brent root finder on the first derivative.

Key findings

All seven EOS fits had small errors, yet the point estimates differed: minimum volumes spread over about 0.01 Å3 (Pd) and 0.03 Å3 (Au), energies were nearly identical, and bulk moduli spread over about 15 GPa (Pd, roughly 165 to 180 GPa) and 20 GPa (Au, roughly 137 to 160 GPa), with the polynomial the outlier. The training window mattered: using the full data with stabilized jellium gave bulk moduli about 6 GPa larger for both metals. Delta method standard errors grew with model RMSE, for example 0.069 GPa for Pd stabilized jellium versus 4.3 GPa for the polynomial. For bulk modulus, the approximately 95% intervals of different models did not overlap, so these uncertainties are tied to the chosen model, which the authors call model-specific.

Bayesian regression gave posteriors quantitatively close to the delta method intervals (the delta intervals were somewhat wider for Pd Poirier-Tarantola), with HMC tighter and closer to the regression optimum than SVI. SVI took about 5 minutes per model and was sensitive to initialization, HMC about 20 minutes but more robust. Bayesian regression is therefore also model-specific. The derivative-aware GP gave the widest distributions for all three properties, and they covered the predictions of essentially every analytical EOS except the polynomial, whose error is one to two orders of magnitude larger. The authors describe this as model-general uncertainty, a smooth, residual-weighted combination over plausible model forms rather than a multimodal sum of individual model distributions. GP variance grows with derivative order and near the edges of the data, and the RBF kernel encodes no EOS physics, so it should not be trusted for extrapolation.

On practicality, nonlinear regression is fastest and easiest and is recommended when parameters are roughly normal; Bayesian regression may help for skewed or discrete parameter distributions; GPs are easy to fit (minutes) but scale as O(n3) in data size.

Significance

The paper gives a clear operational distinction between parameter uncertainty within one model and uncertainty that includes model selection, and shows that common frequentist and Bayesian regression tools capture only the former. The joint function-derivative GP is a practical way to report physical properties that depend on derivatives with uncertainty that spans model form, which the authors argue is essential when the true model is unknown, as in Bayesian optimization or sparse data. Suggested extensions include DFT functionals and properties derived from molecular dynamics.