Source-linked AI summary
A Fast and Scalable Method for A-Optimal Design of Experiments for Infinite-dimensional Bayesian Nonlinear Inverse Problems
Alen Alexanderian, Noemi Petra, Georg Stadler, Omar Ghattas
TL;DR
The paper addresses sensor placement for Bayesian nonlinear inverse problems governed by PDEs, aiming to reduce uncertainty in an inferred parameter field. It develops an A-optimal, PDE-constrained method using a MAP-centered Gaussian approximation and randomized trace estimation. In porous-medium flow experiments, objective and gradient PDE-solve costs, as well as quasi-Newton iterations, are essentially independent of parameter and sensor dimensions.
Problem
The paper asks how to choose sensor locations that minimize uncertainty in parameters inferred from large-scale Bayesian nonlinear inverse problems governed by PDEs.
Method
It approximates the nonlinear posterior with a MAP-centered Gaussian covariance, estimates covariance traces randomly, and solves the resulting PDE-constrained OED problem with adjoint gradients.
Results
Objective and gradient PDE-solve costs, along with quasi-Newton iterations, are essentially insensitive to parameter dimension and only weakly dependent on sensor-candidate dimension.
Takeaways & Limitations
The method provides a scalable approach to sensor placement for elliptic-PDE coefficient inference in large-scale Bayesian inverse problems.
Takeaways & Limitations
The method relies on a Gaussian posterior approximation, which is most appropriate when the parameter-to-observable map is sufficiently well approximated by a linearization over high-posterior-probability parameters.
Abstract
from arXiv · showhide
We address the problem of optimal experimental design (OED) for Bayesian nonlinear inverse problems governed by PDEs. The goal is to find a placement of sensors, at which experimental data are collected, so as to minimize the uncertainty in the inferred parameter field. We formulate the OED objective function by generalizing the classical A-optimal experimental design criterion using the expected value of the trace of the posterior covariance. We seek a method that solves the OED problem at a cost (measured in the number of forward PDE solves) that is independent of both the parameter and sensor dimensions. To facilitate this, we construct a Gaussian approximation to the posterior at the maximum a posteriori probability (MAP) point, and use the resulting covariance operator to define the OED objective function. We use randomized trace estimation to compute the trace of this (implicitly defined) covariance operator. The resulting OED problem includes as constraints the PDEs characterizing the MAP point, and the PDEs describing the action of the covariance operator to vectors. The sparsity of the sensor configurations is controlled using sparsifying penalty functions. We elaborate our OED method for the problem of determining the sensor placement to best infer the coefficient of an elliptic PDE. Adjoint methods are used to compute the gradient of the PDE-constrained OED objective function. We provide numerical results for inference of the permeability field in a porous medium flow problem, and demonstrate that the number of PDE solves required for the evaluation of the OED objective function and its gradient is essentially independent of both the parameter and sensor dimensions. The number of quasi-Newton iterations for computing an OED also exhibits the same dimension invariance properties.
1. Introduction.
The paper develops scalable A-optimal experimental design for infinite-dimensional Bayesian nonlinear inverse problems governed by PDEs. It combines Gaussian posterior approximations, randomized trace estimation, PDE-constrained optimization, and adjoint gradients to make sensor placement computationally tractable.
- Motivation: The method targets OED for large-scale Bayesian nonlinear inverse problems, where repeated inverse and forward PDE solves make sensor-placement optimization difficult.The challenge spans state, parameter, and data dimensions.
- Contributions: The paper proposes an infinite-dimensional A-optimal formulation that minimizes expected average posterior variance for Bayesian nonlinear inverse problems.The formulation is designed to preserve the structure of the underlying function-space problem.
- Method: A Gaussian approximation centered at the MAP point replaces unavailable nonlinear posterior covariances with the inverse Hessian of the regularized data-misfit functional.The approximation is exact for linear parameter-to-observable maps and can be accurate when local linearization is adequate.
- Method: Randomized trace estimators and PDE constraints for MAP and inverse-Hessian actions turn the OED objective into a scalable PDE-constrained optimization problem.The expectation over data is approximated using samples generated from the prior and noise models.
- Application: The method is applied to elliptic-PDE coefficient inference, where sensor designs are evaluated through posterior variance and MAP-estimator quality.The application interprets sensor placement as selecting well locations for inferring a log permeability field.
2. Preliminaries.
The preliminaries define the infinite-dimensional Bayesian setting, Gaussian priors, parameter-to-observable maps, posterior measures, MAP estimation, and weighted sensor designs. They also connect covariance traces to average pointwise variance and motivate randomized trace estimation for implicit operators.
- Probability measures on Hilbert spaces: The framework models parameters as random fields in an infinite-dimensional Hilbert space, with covariance operators that are positive, self-adjoint, and trace-class.The paper uses H = L2(D) for a bounded domain D.
- Probability measures on Hilbert spaces: For an H-valued random field, the covariance trace is proportional to the spatially averaged pointwise variance, providing the basis for infinite-dimensional A-optimal design.This identity links operator uncertainty to a physically interpretable average variance over the domain.
- Bayesian inversion: The Bayesian inverse problem uses a Gaussian prior with covariance Cpr = A^-2 and a parameter-to-observable map whose evaluation typically requires a forward PDE solve.The prior covariance models correlation lengths and pointwise variance while remaining trace-class in two and three dimensions.
- Bayesian inversion: The posterior is defined from the prior and likelihood, while the MAP point minimizes a regularized data-misfit functional that depends on the observations.The MAP solution need not be unique, and its dependence on data complicates OED because data are unavailable beforehand.
- Experimental design: Sensor designs assign nonnegative weights to candidate locations, and these weights enter inference through a weighted likelihood and design-dependent posterior statistics.Weights can represent sensor activation or be combined with sparsifying penalties to encourage sparse configurations.
- Trace estimation: Randomized trace estimation can approximate traces of high-dimensional implicitly defined covariance operators using a small number of random vectors.This supports scalable evaluation of the A-optimal objective.
3. A-optimal design for Bayesian linear inverse problems.
For Bayesian linear inverse problems, A-optimal design minimizes average posterior variance, with formulations that can exploit low-rank structure to avoid further PDE solves.
- A-optimal design minimizes the average posterior variance, equivalently the trace of the posterior covariance operator.
- For a linear parameter-to-observable map with Gaussian prior, the posterior covariance has a closed-form operator expression.
- A low-rank singular value decomposition of the prior-preconditioned parameter-to-observable map can be computed once upfront.
- The resulting A-optimal objective and gradient can be evaluated without further PDE solves.
- Sparsifying penalties, including ℓ1 penalties and continuation toward the ℓ0-norm, control the number of selected sensors.
4. A-optimal design for Bayesian nonlinear inverse problems.
For nonlinear Bayesian inverse problems, the paper replaces the unavailable exact posterior covariance with a MAP-centered Gaussian approximation and randomized trace estimates to form a tractable A-optimal objective.
- The expected A-optimal objective averages posterior covariance traces over parameters and data because nonlinear posterior covariance depends on the experimental data.
- The method approximates the nonlinear posterior by a Gaussian centered at the MAP point, using the inverse Hessian or its approximation as covariance.
- Monte Carlo samples of parameters and noise generate data samples whose MAP points and Hessians enter the approximated OED objective.
- Randomized trace estimation replaces the impractical complete-basis trace computation and yields a computationally tractable OED objective.
- The resulting design problem is a Hessian-constrained formulation with sparsifying penalties and differentiability assumptions for gradient-based optimization.
5. OED for coefficient field inference in an elliptic PDE.
The paper specializes its nonlinear Bayesian OED framework to elliptic-PDE coefficient inference, expressing MAP estimation, Hessian actions, and design gradients through PDE-constrained systems.
- The elliptic-PDE application infers a log coefficient field from pointwise state observations under a PDE-constrained Bayesian inverse problem.
- The OED formulation is bilevel: inner MAP estimation is constrained by state, adjoint, and gradient equations, while the outer problem selects sensor weights.
- Hessian-vector actions are computed by coupled incremental state, incremental adjoint, and parameter equations, solved iteratively with Krylov methods.
- Objective evaluation requires MAP systems for each data sample and Hessian systems for each randomized trace vector, with PDE constraints encoding both components.
- Adjoint variables enforce the PDE constraints and provide gradients with respect to the sensor-design vector for gradient-based optimization.
- Forward-like PDE solves dominate computational cost, and the framework targets counts independent of discretized parameter and sensor dimensions.
- The paper observes dimension-independent quasi-Newton iteration counts in its example, although a general dimension-independence argument is difficult for quasi-Newton methods.
6. Example 1: Idealized subsurface flow.
The idealized subsurface-flow example constructs a prior from sparse permeability measurements, computes sparse A-optimal sensor designs, and tests their accuracy and scalability. The designs reduce posterior variance and expected MAP error, while computational costs remain largely insensitive to parameter and sensor dimensions.
- Setup of forward problem: The experiment builds a Gaussian prior for log permeability from five point measurements, using a covariance operator that imposes stronger correlation in the y-direction.The discretized state, adjoint, and parameter variables use a triangular finite-element mesh with 1,121 degrees of freedom.
- Setup of forward problem: With ℓ0 sparsification, the method produces optimal configurations with 10 sensors at γ = 0.008 and 20 sensors at γ = 0.005.The OED uses five experimental data samples and 20 random vectors in the trace estimator.
- Effectiveness of A-optimal design: For a specific truth model, optimal designs correlate lower average variance with lower L2-error of the MAP estimator, with the advantage over random designs stronger for 10 than 20 sensors.The paper notes that optimal placement matters more when sensors are scarce.
- Effectiveness of A-optimal design: The A-optimal designs minimize average posterior variance and also achieve minimal expected error between the true parameter and MAP point for 10- and 20-sensor configurations.The evaluation uses 50 independent samples and compares each optimal design with 30 randomly chosen designs.
- Scalability and performance: The OED objective, gradient, and optimization require costs insensitive to parameter dimension and only weakly dependent on candidate-sensor count, while PDE solves form the main computational building block.The study measures inner and outer CG iterations and interior-point quasi-Newton iterations.
7. Example 2: Subsurface flow based on SPE10 model.
The SPE10 porous-media experiment constructs an A-optimal sensor design and evaluates it for recovering a heterogeneous log-permeability field. The optimal placement outperforms random designs in both MAP error and average posterior variance.
- Bayesian inverse problem setup: The SPE10 test uses a vertical slice of three-dimensional permeability data with one injection well and four production wells.The injection well is centered, while production wells are placed at the domain corners.
- Bayesian inverse problem setup: The forward model represents the injection well as a mollified point source and imposes zero pressure at four corner production wells.The remainder of the boundary uses homogeneous Neumann conditions.
- Bayesian inverse problem setup: The prior uses five log-permeability estimates near the injection and production wells, with linear finite elements and 10,202 degrees of freedom.The prior covariance is constructed from a differential operator with pointwise contributions.
- A-optimal design of experiments: The OED uses 128 candidate locations, 20 randomized trace vectors, and six continuation iterations to obtain a binary design.Interior-point iterations stop at relative residual 10^-5 or after 100 iterations.
- A-optimal design of experiments: The Bayesian inverse problem is solved using pressure data collected at the A-optimal locations, with a finer forward mesh of 237,573 degrees of freedom.The finer mesh captures extreme permeability variations in the truth field.
- Results: The A-optimal placement outperforms random designs with the same sensor count in relative MAP error and average posterior variance.The comparison uses 22-sensor designs.
8. Conclusions and remarks.
The paper presents a scalable PDE-constrained method for A-optimal design whose forward-like solve cost is insensitive to parameter and sensor dimensions. It also identifies approximation, prior, sparsification, solve-count, and data-sample limitations.
- Conclusions: The method’s forward-like PDE-solve cost is insensitive to the discretized parameter-field dimension and sensor dimension.The formulation uses adjoint-based gradients for efficient gradient-based optimization.
- Conclusions: The OED formulation uses an inverse problem as an inner optimization and inverse-Hessian action equations as outer constraints.A variational formulation makes the required third derivatives tractable.
- Limitations: A Gaussian posterior approximation is a limitation when the parameter-to-observable map is not sufficiently linear over the high-posterior-probability region.Relaxing this approximation is especially challenging because the inverse problem is nested inside OED.
- Limitations: Prior samples with substantially different features can produce a design that is suboptimal for the truth parameter.The authors suggest updating the prior after initial field experiments and recomputing the OED.
- Limitations: Sparsification controls sensor count only indirectly, so multiple OED solves may be needed to choose the penalty parameter experimentally.This trade-off makes a combinatorial placement problem computationally tractable.
- Limitations: Computing an optimal design still requires many forward, adjoint, or incremental PDE solves, although shared Hessians and sample-level parallelism can reduce the burden.Low-rank Hessian approximations are identified as a possible mitigation.
- Limitations: The experiments use few prior data samples because each additional sample requires another inverse problem at every OED iteration.The authors expect diminishing returns from increasing the sample count but leave the required number for future study.
Appendix B. Gradient derivation of OED objective function
Appendix B derives the gradient of the PDE-constrained OED objective with respect to sensor weights by enforcing state and adjoint equations. It also shows that OED adjoint variables inherit structure from the state system and can be eliminated to simplify the gradient computation.
- Gradient derivation: The gradient with respect to the sensor-weight vector is obtained by differentiating the OED Lagrangian after enforcing vanishing variations in the state and adjoint variables.The weight vector enters through the weight matrix Wσ, and the resulting derivatives are assembled in vector form.
- State and adjoint equations: The OED state equations are recovered by requiring the Lagrangian variation with respect to OED adjoint variables to vanish.These equations define the state variables used in the constrained objective evaluation.
- State and adjoint equations: The OED adjoint equations have the same system structure as the OED state equations, differing in their right-hand sides up to a constant.This structural relation permits elimination of the OED adjoint variables and simplifies the right-hand sides in the gradient system.
- Discretization: Finite-element discretization represents the parameter, state, and adjoint variables with boldfaced discrete operators and variables for numerical computation.The discrete OED function and its gradient are formed from these discretized constrained systems.
Appendix C. Discretization and computational details.
Appendix C describes the finite-dimensional implementation of the OED objective and gradient using mass-weighted trace estimation and block linear solves. The computation first obtains MAP states, then covariance-action variables, and finally the objective and gradient quantities.
- Trace estimation: A mass-weighted inner product is used to discretize the infinite-dimensional Hilbert space, while Gaussian vectors provide the randomized trace estimator.The estimator uses z_k = M^-1/2ν_k with ν_k drawn from N(0, I).
- Computational workflow: The implementation solves the MAP optimization problem and evaluates the associated state and adjoint variables for each design and data sample.These quantities provide the point at which the remaining Hessian and covariance-related systems are evaluated.
- Linear systems: The covariance-action variables are computed from a block linear system whose observation block is D = B^T WσB and whose remaining blocks discretize differential operators at the MAP point.Block elimination of v_ik and q_ik reduces the system before solving H y_ik = z_k.
- Objective and gradient evaluation: Once y_ik is available, the OED objective is computed; additional analogous linear solves produce the adjoint variables required for the gradient.The discretized gradient is assembled after all required state and adjoint quantities have been obtained.