Source-linked AI summary
A-optimal design of experiments for infinite-dimensional Bayesian linear inverse problems with regularized $\ell_0$-sparsification
Alen Alexanderian, Noemi Petra, Georg Stadler, Omar Ghattas
TL;DR
The paper addresses efficient sensor placement for infinite-dimensional Bayesian linear inverse problems governed by expensive PDEs. It combines prior-preconditioned low-rank approximations, randomized trace estimation, and continuation-based ℓ0-sparsification, obtaining binary designs with PDE-solve costs independent of parameter and sensor dimensions. In the reported advection-diffusion experiments, regularized ℓ0 designs outperform ℓ1-sparsified designs.
Problem
Computing A-optimal sensor locations is challenging for PDE-governed Bayesian inverse problems with infinite-dimensional or high-dimensional parameters.
Method
The method uses a prior-preconditioned low-rank parameter-to-observable approximation, randomized trace estimation, and continuation penalties approximating the ℓ0-norm.
Results
Regularized ℓ0-sparsified designs consistently outperform ℓ1-sparsified designs and improve over uniform and random designs.
Takeaways & Limitations
An optimal design can be computed with a number of forward PDE solves independent of the parameter and candidate sensor dimensions.
Takeaways & Limitations
The method relies on a linear parameter-to-observable map and Gaussian prior and noise distributions.
Abstract
from arXiv · showhide
We present an efficient method for computing A-optimal experimental designs for infinite-dimensional Bayesian linear inverse problems governed by partial differential equations (PDEs). Specifically, we address the problem of optimizing the location of sensors (at which observational data are collected) to minimize the uncertainty in the parameters estimated by solving the inverse problem, where the uncertainty is expressed by the trace of the posterior covariance. Computing optimal experimental designs (OEDs) is particularly challenging for inverse problems governed by computationally expensive PDE models with infinite-dimensional (or, after discretization, high-dimensional) parameters. To alleviate the computational cost, we exploit the problem structure and build a low-rank approximation of the parameter-to-observable map, preconditioned with the square root of the prior covariance operator. This relieves our method from expensive PDE solves when evaluating the optimal experimental design objective function and its derivatives. Moreover, we employ a randomized trace estimator for efficient evaluation of the OED objective function. We control the sparsity of the sensor configuration by employing a sequence of penalty functions that successively approximate the $\ell_0$-"norm"; this results in binary designs that characterize optimal sensor locations. We present numerical results for inference of the initial condition from spatio-temporal observations in a time-dependent advection-diffusion problem in two and three space dimensions. We find that an optimal design can be computed at a cost, measured in number of forward PDE solves, that is independent of the parameter and sensor dimensions. We demonstrate numerically that $\ell_0$-sparsified experimental designs obtained via a continuation method outperform $\ell_1$-sparsified designs.
1. Introduction.
The paper develops scalable A-optimal sensor-placement methods for infinite-dimensional Bayesian linear inverse problems governed by PDEs. It combines prior-preconditioned low-rank approximations, randomized linear algebra, and regularized sparsification to obtain efficient, binary sensor designs.
- Motivation: The work targets optimal sensor placement for large-scale Bayesian linear inverse problems governed by PDEs.The design minimizes uncertainty in inferred parameters through the Bayesian A-optimal criterion.
- Design criterion: A-optimal designs minimize the average posterior variance of the inversion parameters.In the infinite-dimensional setting, this criterion is expressed through the trace of the posterior covariance operator.
- Computational approach: The method uses a low-rank SVD surrogate of the prior-preconditioned parameter-to-observable map to avoid repeated PDE solves.Randomized trace estimation and randomized SVD support efficient objective and derivative evaluations.
- Sparsification: A continuation sequence of penalties approximating the ℓ0-norm produces binary sensor configurations, unlike ℓ1-sparsification.The weights represent candidate sensor inclusion, with 0 denoting absence and 1 denoting placement.
- Scalability: The number of forward PDE solves is independent of the discretized parameter and sensor dimensions.Limited independent information from nearby sensors bounds the numerical rank of the prior-preconditioned map as candidate locations increase.
2. Background.
The background formulates Bayesian inverse problems for random fields in Hilbert spaces and describes consistent finite-dimensional discretization. It also introduces randomized trace estimation and randomized SVD for large-scale covariance and parameter-to-observable computations.
- Bayesian formulation: The parameter is modeled as a random field taking values in an infinite-dimensional Hilbert space.Its law is represented as a probability measure, which observations update through the Bayesian inverse problem.
- Bayesian formulation: The paper assumes a linear parameter-to-observable map, additive Gaussian noise, and a Gaussian prior.Under these assumptions, the posterior measure is Gaussian.
- Prior modeling: The prior covariance is defined through the inverse of an elliptic differential operator to support scalable PDE-based computations.The paper chooses C0 = A^-2 with a Laplacian-like operator A.
- Discretization: Finite-dimensional inference uses a subspace discretization with a weighted inner product induced by the finite-element mass matrix.The discretized posterior remains Gaussian and supports variance and sampling computations.
- Randomized linear algebra: Randomized trace estimators approximate traces using matrix-vector products, while randomized SVD constructs low-rank surrogates with independent matrix-vector applications.These operations are suited to large-scale problems involving expensive PDE solves.
3. A-optimal design of experiments for infinite-dimensional Bayesian linear inverse problems.
The paper extends Bayesian A-optimal design to infinite-dimensional random-field inference by minimizing posterior covariance trace over weighted candidate sensor locations. It relaxes binary placement to continuous weights and uses sparsifying penalties to recover sparse designs.
- A-optimal criterion: Infinite-dimensional A-optimal design minimizes the trace of the posterior covariance operator.This extends finite-dimensional average posterior-variance minimization to Bayesian inverse problems with random-field parameters.
- Design representation: A design assigns nonnegative weights to a finite set of candidate sensor locations.The weight vector enters the Bayesian likelihood through a diagonal weighting matrix.
- Binary-design relaxation: Binary sensor placement is relaxed to weights in [0, 1] because directly optimizing binary vectors is combinatorial.Sparsifying penalties and continuation are then used to recover the desired 0–1 structure.
- Design representation: The time-dependent setting uses q = NsNτ observations and a block-diagonal weight matrix across observation times.Each block contains the candidate-sensor weights on its diagonal.
- Sparsification: The OED objective combines posterior covariance trace minimization with a penalty controlling design sparsity.The ℓ1 penalty is convex and yields sparse weights, whereas the proposed continuation family approximates the ℓ0-norm.
- Interpretation: In the Gaussian linear case, minimizing posterior covariance trace is equivalent to minimizing the prior-averaged mean square error of the posterior mean.This quantity is also called the Bayes risk of the posterior mean.
4. Numerical solution of the OED problem (3.4).
The method makes A-optimal design computation tractable by combining randomized trace estimation with a low-rank surrogate of the prior-preconditioned parameter-to-observable map. A continuation sequence of sparsity penalties produces binary sensor designs while avoiding further forward or adjoint PDE solves after surrogate construction.
- Objective and derivatives: The misfit Hessian is decomposed as a weighted sum of sensor-location contributions, enabling derivatives with respect to sensor weights.The decomposition writes H_misfit(w) as the sum of w_j F*E_jF terms.
- Objective and derivatives: A randomized trace estimator approximates the trace of the inverse Hessian used as the OED objective.The estimator requires repeated applications of H(w)^−1 to random vectors.
- Computational considerations: Repeated inverse-Hessian applications remain computationally demanding without the surrogate, and repeated applications of the prior square root are still required.These are identified as residual computational costs of the implementation.
- Low-rank approximation: The prior-preconditioned map is approximated with a low-rank SVD surrogate to accelerate applications of the inverse Hessian.Prior smoothing typically yields faster-decaying singular values for the preconditioned map than for the original map.
- Low-rank approximation: Once the low-rank surrogate is available, the method requires no further forward or adjoint PDE solves for objective and gradient evaluation.The surrogate is used to compute Θ(w) and its gradient during optimization.
- Sparsity enforcement: The continuation procedure decreases ε while reusing each solution as the next initialization, driving the penalty toward an ℓ0 approximation and binary designs.The procedure begins with an ℓ1-penalized solution and then solves successively smaller-ε problems.
5. Model problem setup.
The model problem infers an initial condition from spatio-temporal point observations of a time-dependent advection-diffusion equation in two- and three-dimensional domains. Sensors sample the evolving PDE solution, while Bayesian inversion uses a Gaussian prior and additive Gaussian noise.
- Computational domains: The two- and three-dimensional computational domains contain internal building obstacles, velocity fields, and candidate sensor locations.The two-dimensional domain is [0,1]^2, while the three-dimensional domain is [0,1]^2 × [0,0.5].
- Forward model: The state evolves according to a time-dependent advection-diffusion equation with diffusion coefficient κ and divergence-free velocity field v.The numerical experiments use κ = 0.001 in two dimensions and κ = 0.003 in three dimensions.
- Forward model: The velocity field is obtained from a steady-state Navier-Stokes problem with side walls driving the flow and Reynolds number Re = 50.The flow boundary data drive opposite directions on the left and right walls.
- Parameter-to-observable map: The observation operator extracts point values of the evolving solution at selected sensor locations and measurement times.The parameter-to-observable map first solves the PDE from initial condition m and then applies the observation operator.
- Bayesian inversion: The Bayesian inverse problem uses a Gaussian prior with covariance C0 = A^-2 and additive Gaussian noise, yielding a Gaussian posterior.The posterior mean is also characterized through a regularized deterministic inversion functional.
- Discretization and optimization: Forward and adjoint equations are discretized with finite elements in space and implicit Euler time stepping, and the OED problems use an interior-point solver with BFGS Hessian approximations.The experiments include a two-dimensional numerical study and a three-dimensional scalability test.
6. Numerical results.
Numerical experiments show that prior preconditioning enables accurate low-rank approximations, while the resulting A-optimal design method remains effective as parameter and sensor dimensions increase. Regularized ℓ0-sparsified designs outperform ℓ1, uniform, and random alternatives in the reported experiments.
- Low-rank approximation: Rapid singular-value decay, especially for the prior-preconditioned map ˜F, enables efficient low-rank approximation of the parameter-to-observable map.The singular values of ˜F decay faster than those of F, indicating that prior preconditioning improves compressibility.
- Low-rank approximation: The prior-preconditioned misfit-Hessian spectra nearly coincide across spatial and temporal discretizations, indicating resolution of the dominant problem physics.The discretizations converge toward the infinite-dimensional misfit Hessian.
- Sensor discretization: Refining candidate sensor grids yields convergent singular-value curves because neighboring sensors provide correlated information through diffusion.Increasing sensor locations can increase information and numerical rank, but the curves converge as grids are refined.
- ℓ1-sparsified designs: Rank r ≥ 40 produces very similar optimal objective values, showing that substantial compression has little influence on the computed ℓ1-sparsified design.The study varies r over 10, 15, 20, 30, 40, 60, and 80.
- Scalability: Nearly constant optimization iterations result as parameter dimension increases, while iterations and objective evaluations are largely insensitive to candidate sensor count.The observed scalability is attributed to the low numerical rank of ˜F and Newton-type optimization.
- Design comparison: Φε-sparsified designs consistently outperform ℓ1-sparsified designs and significantly improve over random designs in reducing the exact posterior-covariance trace.The experiments also show diminishing returns when more than 20 sensors are used.
- Trace estimation: Randomized trace estimation errors average 15%, 7%, 5%, 2%, and 1.5% using 1, 5, 10, 20, and 100 random vectors, respectively.Using more random vectors reduces trace-estimation error and stabilizes the resulting designs.
7. Concluding remarks.
The method is scalable but relies on linear parameter-to-observable maps, Gaussian prior and noise, and effective low-rank approximations. Its continuous relaxation makes combinatorial sensor placement tractable, while nonlinear maps and alternative criteria remain important extensions.
- The method assumes a linear parameter-to-observable map and Gaussian prior and noise distributions.
- Its computational efficiency depends on low-rank approximations of the preconditioned parameter-to-observable map and properties of the forward and observation operators.
- Continuous weights combined with sparsification indirectly control sensor counts through γ but make combinatorial optimal sensor placement computationally tractable.
- Future work includes alternative infinite-dimensional optimal experimental design criteria and nonlinear parameter-to-observable maps.
Appendix A. The mass-weighted trace estimator.
The appendix establishes the mass-weighted trace estimator used for posterior-covariance calculations. It shows unbiasedness by transforming standard Gaussian vectors and using the spectral decomposition of an M-symmetric operator.
- For an M-symmetric linear mapping A, the spectral decomposition uses M-orthogonal eigenvectors and real eigenvalues.
- With y having independent standard normal entries and z = M^-1/2y, ⟨z, Az⟩_M is an unbiased estimator of tr(A).
- The transformed coordinates q = V^T M^1/2y are standard Gaussian, enabling the trace-estimator expectation calculation through squared normal variables.
- The appendix relates posterior covariance traces to expected mean-square error for Gaussian Bayesian linear inverse problems.
- The posterior formulation assumes a Gaussian prior with additive Gaussian noise and uses the posterior mean's dependence on data.