Surrogate models and reduction methods for UQ and inference in - - PowerPoint PPT Presentation

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks

Surrogate models and reduction methods for UQ and inference in large-scale models

Olivier Le Maître

Laboratoire d’Informatique pour la Mécanique et les Sciences de l’Ingénieur Joint work with Loïc Giraldi, Omar M. Knio and Ibrahim Hoteit (KAUST) Guoto Li (Duke)

WTDMF 2017

• O. Le Maître (LIMSI)

WTDMF 2017 1 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks

Table of contents

• 1. Motivation
• 2. Global Sensitivity Analysis via Surrogates
• 3. Bayesian Inference
• 4. Conclusions and outlooks
• O. Le Maître (LIMSI)

WTDMF 2017 2 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks

Uncertainty Quantification

Many physical systems are characterized by: the presence of broad range of scales (turbulence, porous flow, . . . ) subjected to highly varying / unknown forcing and source terms (finance, wind-load, gust, waves, earthquake, . . . ) details at all scales are often not relevant for the prediction (turbulence flows, climate modeling, . . . ) Tempering the computational complexity calls for physical models restricted at the scales of interest (upscaling, homogenization, sub-scale modeling, . . . ) modeling the -known- unknown (probabilistic approaches) These procedures introduce parameters: they usually must be calibrated (inferred, learned, identified,. . . ) from experimental observations typically parameters can not be determined exactly (limited observations, because of imperfect models and setups or becauses processes are inherently uncertain Our final goal is to infer model parameters in the most objective manner and propagate them to assess the predictive uncertainty. global sensitivity analysis & Bayesian inference

• O. Le Maître (LIMSI)

WTDMF 2017 3 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks

Table of contents

• 1. Motivation
• 2. Global Sensitivity Analysis via Surrogates
• 3. Bayesian Inference
• 4. Conclusions and outlooks
• O. Le Maître (LIMSI)

WTDMF 2017 4 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Input Uncertainty

Sources of data uncertainty Inherent variability (e.g. industrial processes). Epistemic uncertainty (e.g. model constants). May not be fully reducible, even theoretically.

• O. Le Maître (LIMSI)

WTDMF 2017 5 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Input Uncertainty

Sources of data uncertainty Inherent variability (e.g. industrial processes). Epistemic uncertainty (e.g. model constants). May not be fully reducible, even theoretically. Probabilistic framework Define an abstract probability space (Ω, A, dµ). Consider model parameter Q as random quantity: Q(ω), ω ∈ Ω. Simulation output S is random and on (Ω, A, dµ).

• O. Le Maître (LIMSI)

WTDMF 2017 5 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Input Uncertainty

Sources of data uncertainty Inherent variability (e.g. industrial processes). Epistemic uncertainty (e.g. model constants). May not be fully reducible, even theoretically. Probabilistic framework Define an abstract probability space (Ω, A, dµ). Consider model parameter Q as random quantity: Q(ω), ω ∈ Ω. Simulation output S is random and on (Ω, A, dµ). Parameter Q and simulation output S are dependent random quantities (through the mathematical model M): M(S(ω), Q(ω)) = 0, ∀ω ∈ Ω.

• O. Le Maître (LIMSI)

WTDMF 2017 5 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Uncertainty Propagation

Propagation of data uncertainty Parameter density M(S, Q) = 0 Solution density

• O. Le Maître (LIMSI)

WTDMF 2017 6 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Uncertainty Propagation

Propagation of data uncertainty Parameter density M(S, Q) = 0 Solution density Variability in model output: numerical error bars. Assessment of predictability. Support decision making process. What type of information (abstract quantities, confidence intervals, density estimations, structure of dependencies, . . . ) one needs?

• O. Le Maître (LIMSI)

WTDMF 2017 6 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

UQ Methods

Deterministic methods Sensitivity analysis (adjoint based, AD, . . . ): local. Perturbation techniques: limited to low order and simple data uncertainty. Neumann expansions: limited to low expansion order. Moments method: closure problem (non-Gaussian / non-linear problems). Simulation techniques Monte-Carlo Spectral Methods

• O. Le Maître (LIMSI)

WTDMF 2017 7 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

UQ Methods

Deterministic methods Simulation techniques Monte-Carlo Generate a sample set of parameter realizations and compute the corresponding sample set of model ouput. Use sample set based random estimates of abstract characterizations (moments, correlations, . . . ). Plus: Very robust and re-use deterministic codes: (parallelization, complex data uncertainty). Minus: slow convergence of the random estimates with the sample set dimension. Spectral Methods

• O. Le Maître (LIMSI)

WTDMF 2017 7 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

UQ Methods

Deterministic methods Simulation techniques Monte-Carlo Spectral Methods Parameterization with random variables (RVs). ⊥ projection of solution on the (L2) space spanned by the RVs. Plus: arbitrary level of uncertainty, deterministic approach, convergence rate, information contained. Minus: parameterizations (limited # of RVs), adaptation of simulation tools (legacy codes), robustness (non-linear problems, non-smooth output, . . . ). Not suited for model uncertainty

• O. Le Maître (LIMSI)

WTDMF 2017 7 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Spectral methods

Input parametrization Parametrization of Q using N < ∞ independent RVs with prescribed distribution p(ξ ξ ξ): Q(ω) ≈ Q(ξ ξ ξ(ω)), ξ ξ ξ = (ξ1, . . . , ξN) ∈ Ξ. Iso-probabilistic transformations of RVs, Karhunen-Loève expansion: Q(x x x, ω) stochastic field/process, Indentification (e.g. Bayesian). Model Solution expansion

• O. Le Maître (LIMSI)

WTDMF 2017 8 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Spectral methods

Input parametrization Model We assume that for a.e. ξ ξ ξ ∈ Ξ, the problem M(S, Q(ξ ξ ξ)) = 0

1

is well-posed,

2

has a unique solution and that the random solution S(ξ ξ ξ) ∈ L2(Ξ, pξ) : E S2 =

• Ω

S2(ξ ξ ξ(ω))dµ(ω) =

• Ξ

S2(ξ ξ ξ)p(ξ ξ ξ)dξ ξ ξ < +∞. Solution expansion

• O. Le Maître (LIMSI)

WTDMF 2017 8 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Spectral methods

Input parametrization Model Solution expansion Let {Ψ0, Ψ1, . . . } be a basis of L2(Ξ, pξ) then S(ξ ξ ξ) =

• k

skΨk(ξ ξ ξ). Knowledge of the spectral coefficients sk fully determine the random solution. Makes explicit the dependence between Q(ξ ξ ξ) and S(ξ ξ ξ). Need to select the basis and efficient procedure(s) to compute the sk.

• O. Le Maître (LIMSI)

WTDMF 2017 8 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Polynomial Chaos expansions

[Wiener, 1938]

Any well behaved RV U(ω) (e.g. 2nd-order one) defined on (Ω, A, dµ) has a convergent expansion of the form: U(ω) = u0Γ0 +

∞

• i1=1

ui1Γ1(ξi1(ω)) +

∞

• i1=1

i1

• i2=1

ui1,i2Γ2(ξi1(ω), ξi2(ω)) +

∞

• i1=1

i1

• i2=1

i2

• i3=1

ui1,i2,i3Γ3(ξi1(ω), ξi2(ω), ξi3(ω)) + . . . {ξ1, ξ2, . . . }: independent normalized Gaussian RVs. Γp polynomials with degree p, orthogonal to Γq, ∀q < p. Convergence in the mean square sense

[Cameron & Martin, 1947].

• O. Le Maître (LIMSI)

WTDMF 2017 9 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Polynomial Chaos expansions

[Wiener, 1938]

Truncated PC expansion at order No and using N RVs:

U(ω) ≈

P

• k=0

ukΨk(ξ ξ ξ(ω)), ξ ξ ξ = {ξ1, . . . , ξN}, P = (N + No)! N!No! .

{uk}k=0,...,P: deterministic expansion coefficients, {Ψk}k=0,...,P: ⊥ random polynomials wrt the inner product involving the density

• f ξ

ξ ξ:

E {ΨkΨl} = Ψk, Ψl ≡

• Ω

Ψk(ξ ξ ξ(ω))Ψl(ξ ξ ξ(ω))dµ(ω) =

• Ψk(ξ

ξ ξ)Ψl(ξ ξ ξ)p(ξ ξ ξ)dξ ξ ξ = δkl Ψk, Ψk .

• O. Le Maître (LIMSI)

WTDMF 2017 9 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Polynomial Chaos expansions

[Wiener, 1938]

Truncated PC expansion at order No and using N RVs:

U(ω) ≈

P

• k=0

ukΨk(ξ ξ ξ(ω)), ξ ξ ξ = {ξ1, . . . , ξN}, P = (N + No)! N!No! .

{uk}k=0,...,P: deterministic expansion coefficients, {Ψk}k=0,...,P: ⊥ random polynomials wrt the inner product involving the density

• f ξ

ξ ξ:

E {ΨkΨl} = Ψk, Ψl ≡

• Ω

Ψk(ξ ξ ξ(ω))Ψl(ξ ξ ξ(ω))dµ(ω) =

• Ψk(ξ

ξ ξ)Ψl(ξ ξ ξ)p(ξ ξ ξ)dξ ξ ξ = δkl Ψk, Ψk .

Gaussian RVs: p(ξ ξ ξ) = N

i=1 exp(−ξ2

i /2)

√ 2π

= ⇒ Hermite polynomials (Wiener-Hermite expansions) {Ψ0, . . . , ΨP} is an orthogonal basis of SP ⊂ L2(RN, p(ξ ξ ξ)).

• O. Le Maître (LIMSI)

WTDMF 2017 9 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Generalized PC expansion

[Xiu & Karniadakis, 2002]

Askey scheme Distribution of ξi Polynomial familly Gaussian Hermite Uniform Legendre Exponential Laguerre β-distribution Jacobi Also: discrete RVs (Poisson process). U(ω) ≈ P

k=0 ukΨk(ξ

ξ ξ(ω)) where Ψk: classical (or product of) polynomials : spectral expansions

• O. Le Maître (LIMSI)

WTDMF 2017 10 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Non-intrusive methods

Ideas and properties Compute/estimate spectral coefficients via a set of deterministic model solutions Requires a deterministic solver only Overcome issues related to non-linearities. Suffers from the curse of dimensionnality

• O. Le Maître (LIMSI)

WTDMF 2017 11 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Non-intrusive methods Basics

Use code as a black-box Compute/estimate spectral coefficients via a set of deterministic model solutions Requires a deterministic solver only

1

SΞ ≡ {ξ ξ ξ(1), . . . ,ξ ξ ξ(m)} sample set of ξ ξ ξ

2

Let s(i) be the solution of the deterministic problem M s(i), D(ξ ξ ξ(i)) = 0

3

SS ≡ {s(1), . . . , s(m)} sample set of model solutions

4

Estimate expansion coefficients sk from this sample set. Complex models, reuse of determinsitic codes, planification, . . . Error control and computational complexity (curse of dimensionality), . . .

• O. Le Maître (LIMSI)

WTDMF 2017 12 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Non intrusive spectral projection: NISP

Exploit the orthogonality of the basis: E Ψ2

k

• sk = S, Ψk =
• Ξ

S(ξ ξ ξ)Ψk(ξ ξ ξ)p(ξ ξ ξ)dξ ξ ξ. Computation of (P + 1) N-dimensional integrals S, Ψk ≈

Nq

• i=1

w(i)S ξ ξ ξ(i) Ψk(ξ ξ ξ(i)).

• O. Le Maître (LIMSI)

WTDMF 2017 13 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Non intrusive projection Random Quadratures

Approximate integrals from a (pseudo) random sample set SS: S, Ψk ≈ 1 m

m

• i=1

w(i)s(i)Ψk(ξ ξ ξ(i)).

MC LHS QMC

Convergence rate Error estimate Optimal sampling strategy

• O. Le Maître (LIMSI)

WTDMF 2017 14 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Non intrusive projection Deterministic Quadratures

Approximate integrals by N-dimensional quadratures: S, Ψk ≈

Nq

• i=1

w(i)s(i)Ψk(ξ ξ ξ(i)). Quadrature points ξ ξ ξ(i) and weights w(i) obtained by full tensorization of n points 1-D quadrature (i.e. Gauss): Nq = nN partial tensorization of nested 1-D quadrature formula (Féjer, Clenshaw-Curtis): Nq << nN

• O. Le Maître (LIMSI)

WTDMF 2017 15 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Non intrusive projection Deterministic Quadratures

l = 4 l = 5 l = 6

Important development of sparse-grid methods Anisotropic and adaptivity Extension to collocation approach (N-dimensional interpolation)

• O. Le Maître (LIMSI)

WTDMF 2017 15 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Spectral methods

Least square fit “Regression”

Best approximation is defined by minimizing a (weighted) sum of squares of residuals: R2(s0, . . . , sP) ≡

m

• i=1

wi

• s(i) −

P

• k=0

skΨk

• ξ

ξ ξ(i)

2

. Advantages/issues Convergence with number of regression points m Selection of the regression points and “regressors” Ψk Error estimate (cross validation) Compressive sensing Regularization (e.g. LARS, LASSO, BPDN)

m

• i=1

wi

• s(i) −

P

• k=0

skΨk

• ξ

ξ ξ(i)

2

+ λ

P

• k=0

|sk|. ℓ1 penalty: promoting sparsity Selection of λ

• O. Le Maître (LIMSI)

WTDMF 2017 16 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Exploiting the spectral expansion

Truncated PC expansion: S(ω) ≈ P

k=0 skΨk(ξ

ξ ξ(ω)). Convention Ψ0 ≡ 1: mean mode. Expectation of S: E {S} ≡

• Ω

S(ω)dµ(ω) ≈

P

• k=0

sk

• Ξ

Ψk(ξ ξ ξ)p(ξ ξ ξ)dξ ξ ξ = s0. Variance of S: V {S} = E S2 − E {S}2 ≈

P

• k=1

u2

k Ψk, Ψk .

Extension to random vectors & stochastic processes:

 

S1 . . . Sm

  (ω,x

x x, t) ≈

P

• k=0

 

s1 . . . sm

 

k

(x x x, t) Ψk(ξ ξ ξ(ω)).

• O. Le Maître (LIMSI)

WTDMF 2017 17 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Sobol-Hoeffding decomposition

Any S ∈ L2(Xi) has a unique hierarchical orthogonal decomposition of the form S(ξ ξ ξ) = S(ξ1, . . . , ξd) =S0 +

d

• i=1

Si(ξi) +

d

• i=1

d

• j=i+1

Si,j(ξi, ξj)+

d

• i=1

d

• j=i+1

d

• k=j+1

Si,j,k(ξi, ξj, ξk) + · · · + S1,...,d(ξ1, . . . , ξd). Hierarchical: 1st order functionals (Si) → 2nd order functionals (Si,j) → 3rd order functionals (Si,j,l) → · · · → d-th order functional (S1,...,d). Decomposition in a sum of 2k functionals Using ensemble notations: S(ξ ξ ξ) =

• i⊆Q

Si(ξ ξ ξi). Ensemble notations Let D = {1, 2, . . . , d}. Given i ⊆ Q, we denote i∼ := Q \ i its complement set in Q, such that i ∪ i∼ = Q, i ∩ i∼ = ∅. Given i ⊆ Q we define ξ ξ ξi = (ξi1, . . . , ξi|i|)

• O. Le Maître (LIMSI)

WTDMF 2017 18 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Sobol-Hoeffding decomposition

Any S ∈ L2(Ξ) has a unique hierarchical orthogonal decomposition of the form S(ξ ξ ξ) =

• i⊆Q

Si(ξ ξ ξi). Orthogonal: the functionals if the S-H decomposition verify the following

• rthogonality relations:
• Ξ

Si(ξ ξ ξi)Sj(ξ ξ ξj)p(ξ ξ ξ)dξ ξ ξ = Si, Sj = 0, ∀i, j ⊆ Q, i = j. It follows the hierarchical construction Further, the integrals of S with respect to i∼ are in this context conditional expectations, E {S|ξ ξ ξi} =

• Ξ|i∼|

S(ξ ξ ξ)dξ ξ ξi∼ = G(ξ ξ ξi) ∀i ⊆ Q, so the S-H decomposition follows the hierarchical structure S∅ = E {S} S{i} = E S|ξ ξ ξ{i}

• − E {S}

i ∈ Q Si = E {S|ξ ξ ξi} −

• ji

Sj i ⊆ Q, |i| ≥ 2.

Hoeffding, W. (1948). A class of statistics with asymptotically normal distribution. The annals of Mathematical Statistics, 19,

• pp. 293-325.
• O. Le Maître (LIMSI)

WTDMF 2017 19 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Variance decomposition

Because of the orthogonality, the variance V {S} can be decomposed as V {S} =

i=∅

• i⊆Q

V {Si} , V {Si} = Si, Si . V {Si} is interpreted as the contribution to V {S} of the interaction between ξi∈i. Sensitivity indices The partial variances V {Si} are normalized by V {S} to obtain the sensitivity indices: Si(S) = V {Si} V {S} ≤ 1,

i=∅

• i⊆Q

Si(f ) = 1. 1st order sensitivity indices. S{i}∈Q characterize the fraction of the variance due to ξi only. If d

i=1 S{i}(f ) = 1 the model is additive and the impact of the ξi can be

studied separately. Total sensitivity indices. T{i} measures the variability due to the parameter ξi, including all its interactions: T{i} :=

• i∋i

Si ≥ S{i}.

Sobol, I. M. (1993). Sensitivity estimates for nonlinear mathematical models. Wiley, 1, pp. 407-414.

• O. Le Maître (LIMSI)

WTDMF 2017 20 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

• O. Le Maître (LIMSI)

WTDMF 2017 21 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

Sensitivity indices

In many uncertainty problem, the set of uncertain parameters can be naturally grouped into subsets depending on the process each parameter accounts for. For instance, boundary conditions BC, material property ϕ, external forcing F, and Q is the union of these distinct subsets: Q = QBC ∪ Qϕ ∪ QF . The notion of first order and total sensitivity indices can be extended to characterize the influence of the subsets of parameters. For instance, SQϕ =

• i⊆Qϕ

Si, measures the fraction of variance induced by the material uncertainty alone, while TQF =

• i∩QF =∅

Si. measures the fraction of variance due to the external forcing uncertainty and all its interactions.

Sobol, I. M. (1993). Sensitivity estimates for nonlinear mathematical models. Wiley, 1, pp. 407-414. Homma T., Saltelli A., (1996). Importance measures in globla sensitivity analysis of nonlinear models. Reliab. Eng. Syst. Safety, 52:1, pp. 1-17.

• O. Le Maître (LIMSI)

WTDMF 2017 22 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

S-H decomposition from PC expansions

Recall S : ξ ξ ξ ∈ Ξ ⊂ Rd → R, where ξ ξ ξ = (ξ1, . . . , ξd) are independent real-valued r.v. with joint-probability density function p(ξ1, . . . , ξd) =

d

• i=1

pi(ξi). Let { Ψα

α α} be the set of d-variate orthogonal polynomials,

Ψα

α α(ξ

ξ ξ) =

d

• i=1

ψ(i)

αi (ξi),

with ψ(i)

l≥0 the univariate orthonornal polynomials (w.r.t. pi), such that

Ψα

α α, Ψα α α′ = δα α α,α α α′.

Then S ∈ L2(Ξ, pξ

ξ ξ) has a convergent PC expansion

S(ξ ξ ξ) = lim

No→∞

• |α

α α|≤No

Ψα

α α(ξ

ξ ξ)sα

α α,

|α α α| =

d

• i=1

|αi|.

Wiener (1938), Cameron and Martin (1947).

• O. Le Maître (LIMSI)

WTDMF 2017 23 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Uncertainty Quantification

S-H decomposition from PC expansions

Given the truncated PC expansion of S, ˆ S(ξ ξ ξ) =

• |α

α α|≤A

Ψα

α α(ξ

ξ ξ)sα

α α,

• ne can readily obtain the PC approximation of the S-H functionals through

ˆ Si(ξ ξ ξi) =

• |α

α α|≤A(i)

Ψα

α α(ξ

ξ ξi)sα

α α,

where the multi-index set A(i) is given by A(i) = {α α α ∈ A; αi > 0 for i ∈ i, αi = 0 for i / ∈ i} A. For the sensitivity indices it comes Si(ˆ S) =

• α

α α∈A(i) s2 α α α Ψα α α, Ψα α α

• α

α α∈A s2 α α α Ψα α α, Ψα α α

, T{i}(ˆ S) =

• α

α α∈T (i) s2 α α α Ψα α α, Ψα α α

• α

α α∈A s2 α α α Ψα α α, Ψα α α

, where T (i) = {α α α ∈ A; αi > 0 for i ∈ i}

Crestaux T., Le Maitre O. and Martinez J.M., (2009). Polynomial Chaos expansion for sensitivity analysis. Reliabg. Eng. Syst. Safety, 94:7, pp. 1161-1172. Sobol, I. M. (2001). Global sensitivity indices for nonlinear mathematical models and their Monte Carlo estimates. Mathematics and Computers in Simulations, 55, pp. 271-281.

• O. Le Maître (LIMSI)

WTDMF 2017 24 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Application

Application to complex models

Large scale models Number of simulations to perform (and store) Non-intrusive surrogates many model output (e.g. all nodal values over large mesh and many time-steps) Computationally intensive, although embarrassingly parallel (communications?) Calls for reduction strategies: PC expansion for the reduced coordinates in properly selected basis: POD, DMD?,. . . Large number of independent parameters A priori large spectral bases Promoting sparsity is nice but needs appropriate implementation Not necessarily effective Calls for ad-hoc reduction strategies (in the parameter space!): Tensor format, low-rank, and transformations in the input space S(ξ ξ ξ) ≈

R

• r=1

Sr(w w wT

r ξ

ξ ξ), Sr ∈ πNo(R).

• O. Le Maître (LIMSI)

WTDMF 2017 25 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Application

Circulation in the Gulf of Mexico

Objective: Assess impact of uncertainties in Wind Forcing & alertInitial Condition Simulations over a 30 days period, using HYCOM (MIT) with 4km horizontal resolution and 20 layers. Wind forcing every 3h from couples atmosphere ocean prediction (COAMP) system, plus EOF of a 60 days period. Initial Condition from assimilated simulation, plus EOF of a 14 days period. UQ has 8 (4 + 4) independent random variables with uniform distributions. PCE based on 800 LHS set, using BPDN.

• Fig. 1: Bathymetry of the Gulf of Mexico (meters), SSH

(blue box) and MLD (red box) averaging domains.

SSH(m)

• 0.2
• 0.1

0.1 0.2 0.3 0.4 0.5 0.6 PDF 1 2 3 4 5 HYCOM PCE MLD(m) 5 10 15 20 25 PDF 0.05 0.1 0.15 HYCOM PCE

• Fig. 4: Comparison of the SSH (top) and MLD (bottom)

density functions estimated by KDE method. Empirical es- timations from HYCOM realizations on ξ ξ ξ i ∈ PLHS (red curves) and PC model predictions (blue curves) obtained by evaluating PC surrogates over a fine sampling of Ξ using 105 points.

• O. Le Maître (LIMSI)

WTDMF 2017 26 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Application

Circulation in the Gulf of Mexico

Sensitivity Analysis: SSH at Deep Water Horizon and MLD at loop-current site.

SSH(m)

• 0.2
• 0.1

0.1 0.2 0.3 0.4 0.5 0.6 PDF 1 2 3 4 5 HYCOM PCE MLD(m) 5 10 15 20 25 PDF 0.05 0.1 0.15 HYCOM PCE

• Fig. 4: Comparison of the SSH (top) and MLD (bottom)

density functions estimated by KDE method. Empirical es- timations from HYCOM realizations on ξ ξ ξ i ∈ PLHS (red curves) and PC model predictions (blue curves) obtained by evaluating PC surrogates over a fine sampling of Ξ using 105 points.

Stochastic dimension index

1 2 3 4 5 6 7 8

Sensitivity Indices

0.2 0.4 0.6 0.8 1 S{i} T{i}

Stochastic dimension index

1 2 3 4 5 6 7 8

Sensitivity Indices

0.2 0.4 0.6 0.8 1 S{i} T{i}

• Fig. 6: First and Total order sensitivity indices associated to

each input variables at day 30: (Top) SSH; (Bottom) MLD

• O. Le Maître (LIMSI)

WTDMF 2017 27 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Application

Circulation in the Gulf of Mexico

POD at day 30: SSH and MLD field.

• Fig. 9: First five spatial modes uk in the expansion of the

SSH field at day 30.

• Fig. 13: First five spatial modes uk in the expansion of the

MLD field at day 30.

• O. Le Maître (LIMSI)

WTDMF 2017 28 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Application

Circulation in the Gulf of Mexico

First order sensitivity indices: SSH and MLD field.

• Fig. 16: Sensitivity analysis of the SSH (left) and MLD

(right) fields at day 30. Plotted are the first order sensitiv- ity indices related to the initial condition (top row), wind forcing (center row) and interaction between the two (bot- tom row).

• Fig. 17: 1st order sensitivity of SSH to initial condition (left

column) and MLD to wind forcing (right column) perturba- tions on selected days as indicated.

Li G., Iskandarani M., LeHenaff M., Winokur J., LeMaitre O. and Knio O. (2015), Quantifying initial and wind forcing uncertainty in the Gulf of Mexico, Computational Geosciences, (2016).

• O. Le Maître (LIMSI)

WTDMF 2017 29 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Application

Table of contents

• 1. Motivation
• 2. Global Sensitivity Analysis via Surrogates
• 3. Bayesian Inference
• 4. Conclusions and outlooks
• O. Le Maître (LIMSI)

WTDMF 2017 30 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks

Table of contents

• 1. Motivation
• 2. Global Sensitivity Analysis via Surrogates
• 3. Bayesian Inference
• 4. Conclusions and outlooks
• O. Le Maître (LIMSI)

WTDMF 2017 31 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Outline

Bayesian inference

Calibration We want to calibrate the model parameter Using observations of the system, usually corrupted by measurement noise Should be able to assess the uncertainty in our estimate (and resulting model prediction) Bayes’ formula p(q|y) = p(y|q)π(q) p(y) π(q) : prior density of the parameter p(y, q): likelihood of the noisy observations The likelihood "measures" the distance between the prediction and the observation, e.g. p(y|q) = N(y − S(q), Σǫ). Sampling of the posterior yielding posterior prediction Hyperparameter (nuisance parameter) Model error?

• O. Le Maître (LIMSI)

WTDMF 2017 32 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Outline

Challenges in Bayesian inference

Bayes’ formula p(q|y) = p(y|q)π(q) p(y) Multiple evaluations of the model p(y|q) ≈ N(y − ˆ S(q), Σǫ). Deriving a likelihood accounting for both experimental and modeling errors is very difficult More is not always better, in particular in the absence of complete information regarding protocols Calls for the selection of robust and informative observations

• O. Le Maître (LIMSI)

WTDMF 2017 33 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Reduction of observations

Table of contents

• 1. Motivation
• 2. Global Sensitivity Analysis via Surrogates
• 3. Bayesian Inference
• 4. Conclusions and outlooks
• O. Le Maître (LIMSI)

WTDMF 2017 34 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Reduction of observations

Optimal Observations Reduction

Motivation Bayesian inference in the case of overabundant data Weather forecasting Seismic wave inversion

• O. Le Maître (LIMSI)

WTDMF 2017 35 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Reduction of observations

Optimal Observations Reduction

Motivation Bayesian inference in the case of overabundant data Weather forecasting Seismic wave inversion Goal Compute an optimal approximation min

V L

P(Q | Y = y), P(Q | W = V T y) L a loss function n (random) observations Y = (Yi)n

i=1

q parameters Q = (Qi)Nq

i=1, Nq ≪ n

r dimensional reduced space V ∈ Rn×r, r ≪ n

• O. Le Maître (LIMSI)

WTDMF 2017 35 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Linear Gaussian Model

Models

Gaussian model Y = BQ + E, Observations: Y ∼ N(mY , CY ) with values in Rn Parameter of interest: Q ∼ N(mQ, CQ) with values in RNq Noise: E ∼ N(mE, CE) with values in Rn Design matrix: B ∈ Rn×Nq Forward model: A(Q) = BQ ∼ N(mA, CA), and CAQ = Cov(A(Q), Q) Reduced model W = V T BQ + V T E, Reduced observations: W ∼ N(mW , CW ) with values in Rr Reduced space: V ∈ Rn×r

• O. Le Maître (LIMSI)

WTDMF 2017 36 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Linear Gaussian Model

Posterior distributions

Given a realization y of Y Full model The posterior distribution is P(Q | Y = y) ∼ N(m⋆, C⋆) where C⋆ = CQ

• CQ + CT

AQC−1 E

CAQ

−1 CQ,

m⋆ = CT

AQC−1 Y (y − mE) + C⋆C−1 Q mQ.

Reduced model The posterior distribution is P(Q | W = V T y) ∼ N(mV , CV ) where CV = CQ

• CQ + CT

AQV

V T CEV−1 V T CAQ

−1

CQ, mV = CT

AQV (V T CY V )−1V T (y − mE) + CV C−1 Q mQ.

• O. Le Maître (LIMSI)

WTDMF 2017 37 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Linear Gaussian Model

Invariance property

Proposition (Invariance property) For all invertible matrices M ∈ Rr×r

∗

, we have mVM = mV and CVM = CV . Posterior distribution invariant under rescaling, rotation or permutation of the

• bservations

Newton method can not be directly used range(V ) is more important than V Use of a Riemannian trust region algorithm on the Grassmann manifolds Gr(r, n), the set of r-dimensional subspaces of Rn (see Absil et al. 2007, Manopt and Pymanopt libraries)

• O. Le Maître (LIMSI)

WTDMF 2017 38 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Optimality criteria

Kullback-Leibler divergence

Kullback-Leibler divergence Given two distributions P(Z0) and P(Z1) with densities fZ0 and fZ1, DKL (P(Z0) P(Z1)) = EZ0

• log fZ0

fZ1

• .

Quantify the “information lost when [P(Z1)] is used to approximate [P(Z0)]” (Burnham and Anderson, 2003) Positive and null iff P(Z0) = P(Z1) Asymmetric quantity

• O. Le Maître (LIMSI)

WTDMF 2017 39 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Optimality criteria

Kullback-Leibler divergence minimization

Kullback-Leibler divergence minimization min

[V ]∈Gr(r,n) DKL

• P(Q | Y = y) P(Q | W = V T y)

Closed form of the functional available A solution to the optimization problem exists A posteriori reduction

• O. Le Maître (LIMSI)

WTDMF 2017 40 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Optimality criteria

Kullback-Leibler divergence minimization

Kullback-Leibler divergence minimization min

[V ]∈Gr(r,n) DKL

• P(Q | Y = y) P(Q | W = V T y)

Closed form of the functional available A solution to the optimization problem exists A posteriori reduction Expected Kullback-Leibler divergence minimization min

[V ]∈Gr(r,n) EY

• DKL
• P(Q | Y ) P(Q | W = V T Y )

Closed form of the functional available A solution to the optimization problem exists A priori reduction

• O. Le Maître (LIMSI)

WTDMF 2017 40 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Optimality criteria

Entropy and mutual information

Given random variables Z, Z0, and Z1, Entropy With Z ∼ P(Z), H (Z) = EZ(− log(fZ(Z))). Amount of information contained by P(Z) Mutual information With Z0 ∼ P(Z0) and Z1 ∼ P(Z1), I(Z0, Z1) = H (Z0) + H (Z1) − H (Z0, Z1) , Amount of information that P(Z0) contains about P(Z1) Symmetric quantity

• O. Le Maître (LIMSI)

WTDMF 2017 41 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Optimality criteria

Mutual information maximization

Theorem (Mutual information maximization) We have max

V ∈Rn×r

∗

I(W , Q) = 1 2

r

• i=1

log λi, where (λi)r

i=1 are the r dominant eigenvalues of the problem

CY v = λCEv, λ ∈ R, v ∈ Rn. A solution to the optimization problem is given by the matrix V with columns being eigenvectors (vi)r

i=1 associated to the eigenvalues (λi)r i=1.

• O. Le Maître (LIMSI)

WTDMF 2017 42 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Optimality criteria

Mutual information maximization

Equivalences The mutual information maximization is equivalent to: the maximization of the expected information gain max

V ∈Rn×r

∗

EW (DKL (P(Q|W ) P(Q))) the minimization of the entropy of the posterior distribution min

V ∈Rn×r

∗

H P(Q|W = V T y) Error estimate I(Y , Q) − I(W , Q) I(Y , Q) = 1 −

r

i=1 log(λi)

n

i=1 log(λi)

= 1 −

r

i=1 log(1 + νi)

Nq

i=1 log(1 + νi)

, where (νi)Nq

i=1 are the dominant eigenvalues of the generalized eigenvalue problem

CAv = νCEv.

• O. Le Maître (LIMSI)

WTDMF 2017 43 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Linear example

Inference problem

Synthetic data For (ti)n

i=1, n = 500, a uniformly drawn sample in (−1, 1),

Yref(ti) = Aref(ti) + E(ti), ∀i ∈ {1, . . . , n}, with Aref ∼ N(mref, Cref) and E ∼ N(mE, CE). Model Yi =

Nq−1

• j=0

Tj(ti)Qj + E(ti), ∀i ∈ {1, . . . , n}, with Tj the Chebyshev polynomial of order j and Nq = 30.

• O. Le Maître (LIMSI)

WTDMF 2017 44 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Linear example

Functionals versus the dimension of the reduced space

−1 −0.5 0.5 1 −4 −2 2 4 t Observations MAP estimate 10 20 30 40 50 10−2 10−1 100 101 102 Dimension of the reduced space Kullback-Leibler divergence (nat) KLE-A KLE-Y KLE-YN KLD EKLD MI 10 20 30 40 50 10−2 10−1 100 101 102 Dimension of the reduced space Expected Kullback-Leibler divergence (nat) KLE-A KLE-Y KLE-YN KLD EKLD MI 10 20 30 40 50 10−3 10−2 10−1 100 Dimension of the reduced space Relative error on the mutual information KLE-A KLE-Y KLE-YN KLD EKLD MI Error estim.

• O. Le Maître (LIMSI)

WTDMF 2017 45 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Nonlinear problem

Inference problem: nonlinear models

Synthetic data Given two random samples (si)n

i=1 and (ti)n i=1 being independent and uniformly

distributed in (−1, 1), with n = 2000, Yref(si, ti) = exp(Fref(si, ti)) + E(si, ti), ∀i ∈ {1, . . . , n}, where Fref ∼ N(0, Cref), E ∼ N(0, CE). Model Yi = Ai(Q) + E(si, ti), ∀i ∈ {1, . . . , n}, where Ai(Q) = exp((BQ)i), Q ∼ N(0, CQ), and q = 30. Columns of B: dominant eigenvectors of Cref CQ = diag(λ1, . . . , λq): dominant eigenvalues of Cref

• O. Le Maître (LIMSI)

WTDMF 2017 46 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Nonlinear problem

Error estimation

A second-order Taylor expansion around the MAP estimate yields log fQ(q|Y = y) ≈ log fQ(qMAP|Y = y) + 1 2 (q − qMAP)T ∇2 log fQ(qMAP|Y = y)(q − qMAP) Laplace approximation Local approximation by the normal distribution N(qMAP, CMAP) where CMAP = −(∇2 log fQ(qMAP|Y = y))−1

• O. Le Maître (LIMSI)

WTDMF 2017 47 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Nonlinear problem

Error estimation

A second-order Taylor expansion around the MAP estimate yields log fQ(q|Y = y) ≈ log fQ(qMAP|Y = y) + 1 2 (q − qMAP)T ∇2 log fQ(qMAP|Y = y)(q − qMAP) Laplace approximation Local approximation by the normal distribution N(qMAP, CMAP) where CMAP = −(∇2 log fQ(qMAP|Y = y))−1 Error estimation Monitor the relative error on the MAP ǫ and the Hessian ǫH ǫ = qMAP

V

− qMAP qMAP and ǫH = (CMAP

V

)−1 − (CMAP)−1Fro (CMAP)−1Fro .

• O. Le Maître (LIMSI)

WTDMF 2017 47 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks Nonlinear problem

Errors versus the dimension of the reduced space σFref = 0.301 (top), σFref = 1.501 (bottom)

20 40 60 80 100 10−3 10−2 10−1 100 Dimension of the reduced space ǫ KLE-A KLE-Y KLE-YN KLD EKLD MI 20 40 60 80 100 10−3 10−2 10−1 100 Dimension of the reduced space ǫH KLE-A KLE-Y KLE-YN KLD EKLD MI 200 400 600 800 1,000 10−4 10−3 10−2 10−1 100 Dimension of the reduced space ǫ KLE-A KLE-Y KLE-YN KLD EKLD MI 200 400 600 800 1,000 10−4 10−3 10−2 10−1 100 Dimension of the reduced space ǫH KLE-A KLE-Y KLE-YN KLD EKLD MI

• O. Le Maître (LIMSI)

WTDMF 2017 48 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks

Table of contents

• 1. Motivation
• 2. Global Sensitivity Analysis via Surrogates
• 3. Bayesian Inference
• 4. Conclusions and outlooks
• O. Le Maître (LIMSI)

WTDMF 2017 49 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks

Conclusions and outlooks

Summary Reduction approach to enable UQ analyses: both at the deterministic and stochastic level Empirical Orthogonal Decomposition (POD, KL, . . . ) is suitable for fine sensitivity analyses of large dimensional fields Information theoretic reduction approaches are promising Should consider reduction in the input (parameter) and output (deterministic) spaces Reduction strategy(ies) should be goal-oriented Outlooks Couple the approach with filtering techniques (EnKF) Features selection for Bayesian inference Goal-oriented design of experiment

• O. Le Maître (LIMSI)

WTDMF 2017 50 / 51

Motivation Global Sensitivity Analysis via Surrogates Bayesian Inference Conclusions and outlooks

Thank you Merci

• O. Le Maître (LIMSI)

WTDMF 2017 51 / 51