In this note, we present a complete analytic calculation of variance-based sensitivity indices using GP regression

Kurdyaeva, Tamara

Milias Argeitis, Andreas

NARCIS 

Analytical calculation of Sobol sensitivity indices for Gaussian Processes with a squared exponential covariance function

University of Groningen Research Database

  
 University of Groningen
Analytical calculation of Sobol sensitivity indices for Gaussian Processes with a squared
exponential covariance function
Milias Argeitis, Andreas; Kurdyaeva, Tamara
IMPORTANT NOTE: You are advised to consult the publisher's version (publisher's PDF) if you wish to cite from
it. Please check the document version below.
Publication date:
2018
Link to publication in University of Groningen/UMCG research database
Citation for published version (APA):
Milias Argeitis, A., & Kurdyaeva, T. (2018). Analytical calculation of Sobol sensitivity indices for Gaussian
Processes with a squared exponential covariance function.
Copyright
Other than for strictly personal use, it is not permitted to download or to forward/distribute the text or part of it without the consent of the
author(s) and/or copyright holder(s), unless the work is under an open content license (like Creative Commons).
Take-down policy
If you believe that this document breaches copyright please contact us providing details, and we will remove access to the work immediately
and investigate your claim.
Downloaded from the University of Groningen/UMCG research database (Pure): http://www.rug.nl/research/portal. For technical reasons the
number of authors shown on this cover page is limited to 10 maximum.
Download date: 18-03-2018
Analytical calculation of Sobol sensitivity indices for Gaussian Processes
with a squared exponential covariance function
Notation and modeling assumptions
• Column vector of random inputs: X = [X1...Xd]T ∈ Rd, X′ = [X′1...X′d]T ∈ Rd,
• Column subvector of X without Xs element: X−s = [X1 . . . Xs−1 Xs+1 . . . Xd]T ∈ Rd−1
• Column vector consists of elements X′−s and Xs: Xˆ = [X′1 . . . Xs . . . X′d]T ∈ Rd
• Random inputs Xs, s = 1, ..., d are mutually independent and each input follows a standard uniform
distribution: Xs ∼ U(0, 1), s = 1, ..., d
• Function F represents a black-box model of interest: Y = F(X), where Y is the model output
• Observed n training data points is denoted by X˜ = (˜x(1), . . . , x˜(n))
• Test n∗ data points is denoted by X∗ =
(
x(1)∗ , . . . , x
(n∗)∗
)
• Vector of function F values at n training points X˜: y˜ =
[˜
y(1), . . . , y˜(n)
]T
, y˜(i) = F
(˜
x(i)
)
, i = 1, ..., n
• Gaussian process (GP) parametrization f (x) ∼ GP(0, k(x, x′))
k(x, x′) = σ2f
d∏
j=1
exp
− (x j − x′j)22 l2j
 = ∗ ∣∣∣∣∣h j = 12 l−2j
∣∣∣∣∣
= σ2f
d∏
j=1
exp
(
−h j(x j − x′j)2
)
(1)
• Posterior GP distribution f∗|X˜, y˜, X∗ ∼ N(m(X∗), v(X∗)), where
m(X∗) =K(X∗, X˜)
[
K(X˜, X˜) + σ2I
]−1
y˜, (2)
v(X∗) =K(X∗, X∗) − K(X∗, X˜)
[
K(X˜, X˜) + σ2I
]−1
K(X˜, X∗), (3)
where the additive term σ2I reflects the fact that the function values are corrupted by additive noise
with variance σ2. Even if σ2 = 0 (as is the case when the function is exactly evaluated), a small
noise term is typically inserted to regularize the inversion of the covariance matrix K(X˜, X˜).
• K = K(X˜, X˜)
• kpq for p = 1, ..., n and q = 1, ..., n denotes the entries of the matrix
[
K + σ2 I
]−1
• wp, p = 1, ..., n are the elements of the vector
[
K + σ2 I
]−1
y˜
• Expectations and variances with respect to the GP posterior distribution are denoted as E∗ and
Var∗, respectively
• erf(·) denotes the error function.
1
Sobol indices
The estimate of the Sobol indices for the output Y with respect to the input Xs
Sˆ s =
E∗{Var{E (Y |Xs)}}
E∗{Var(Y)} , s = 1, ..., d. (4)
Inference for variance in numerator E∗{Var{E (Y |Xs)}}
Let us first calculate the numerator in the formula for Sobol indices (4).
E∗{Var{E(Y |Xs)}} = E∗{E (E(Y |Xs))2} − E∗{E(Y)}2 =
= E∗{E (E(Y |Xs))2} − Var∗{E(Y)} − (E∗{E(Y)})2 =
= I − II − III (5)
First component I from (5): E∗{E (E(Y |Xs))2}
We derive the formulas for every component of (5) separately. The first term I requires the computation
of the following integrals:
E∗{E (E(Y |Xs))2} = E
[
Var∗{E(Y |Xs)} + (E∗{E(Y |Xs)})2
]
=
=
1∫
0
∫
[0,1]d−1
∫
[0,1]d−1
v(x, xˆ) dx−sdx′−sdxs +
1∫
0
∫
[0,1]d−1
∫
[0,1]d−1
m∗(x)m∗(xˆ)dx−sdx′−sdxs =
= L + M. (6)
For the simplification we also provide the derivation of the formulas for the individual terms L and M in
(6) separately.
L =
1∫
0
∫
[0,1]d−1
∫
[0,1]d−1
v(x, xˆ) dx−sdx′−sdxs
3
=
=
1∫
0
∫
[0,1]d−1
∫
[0,1]d−1
K(x, xˆ) dx−sdx′−sdxs −
1∫
0
∫
[0,1]d−1
∫
[0,1]d−1
K(x, X˜)
(
K + σ2
)−1
K(X˜, xˆ) dx−sdx′−sdxs =
=
1∫
0
Z−s(xˆ)dxs −
1∫
0
[
zT−s(x)
(
K + σ2
)−1
z−s(xˆ)
]
dxs (7)
Now we divide the calculation of the integrals into several parts to make it easier to follow
z−s(x) =
∫
[0,1]d−1
K(x, X˜) dx−s
(1
) =σ2f exp
(
−hs(xs − x˜s)2
) ∏
j=−s
1∫
0
exp
(
−h j(x j − x˜ j)2
)
dx j =
= σ2f exp
(
−hs(xs − x˜s)2
) ∏
j=−s
I j (8)
2
Next we evaluate the integral defined as I j in (8)
I j =
1∫
0
exp
(
−h j(x j − x˜ j)2
)
dx j = ∗
∣∣∣∣∣a j = √h j(x j − x˜ j)∣∣∣∣∣ = 1√h j
√
h j(1−x˜ j)∫
−√h j x˜ j
exp
(
−a2j
)
da j =
=
√
pi
2
1√
h j
(
erf
(√
h j(1 − x˜ j)
)
+ erf
(√
h j x˜ j
))
=
√
pi
2
l j
erf  1√
2l j
(1 − x˜ j)
 + erf  1√
2l j
x˜ j

(9)
We can substitute the exact formula of I j (9) into the final formula of z−s in (8)
z−s(x) = σ2f exp
(
− 1
2l2s
(xs − x˜s)2
) ∏
j=−s
√
pi
2
l j
erf  1√
2l j
(1 − x˜ j)
 + erf  1√
2l j
x˜ j
 (10)
Finally we calculate the first double integral in (7) which was marked as Z−s
Z−s =
∫
[0,1]d−1
∫
[0,1]d−1
K(x, xˆ)dx−sdx′−s = σ2f
∏
j=−s
1∫
0

1∫
0
exp
(
−h j(x j − x′j)2
)
dx j
 dx′j = σ2f ∏
j=−s
1∫
0
I jdx′j =
= σ2f
∏
j=−s
√
pi
2
l j
1∫
0
erf  1√
2l j
(1 − x′j)
 + erf  1√
2l j
x′j
 dx′j = σ2f ∏
j=−s
J j (11)
We present the exact analytic calculation of the J j integral in (11)
J j =
√
pi
2
l j

1∫
0
erf
 1√
2l j
(1 − x′j)
 dx′j + 1∫
0
erf
 1√
2l j
x′j
 dx′j
 =
= ∗
∣∣∣∣∣∣∣b j = 1√2l j (1 − x′j), c j = 1√2l j x′j
∣∣∣∣∣∣∣ =
=
√
pi
2
l j
−√2l j
0∫
1/
√
2l j
erf(b j)db j +
√
2l j
1/
√
2l j∫
0
erf(c j)dc j
 =
= 2
√
pil2j
1/
√
2l j∫
0
erf(t j)dt j =
= ∗
∣∣∣∣∣∣u = erf(t j), dv = dt j, du = 2√pi exp(−t2j )dt j, v = t j
∣∣∣∣∣∣ =
= 2
√
pil2j
erf(t j)t j
∣∣∣∣1/√2l j
0
− 2√
pi
1/
√
2l j∫
0
exp
(
−t2j
)
t jdt j
 =
= 2
√
pil2j
 1√
2l j
erf
 1√
2l j
 + 1√
pi
exp
(
−t2j
) ∣∣∣∣1/√2l j
0
 =
= 2
√
pi
2
l jerf
 1√
2l j
 + 2l2j
exp
− 12l2j
 − 1
 (12)
3
Then we can write down a closed form to compute Z−s in (11)
Z−s = σ2f
∏
j=−s
2
√
pi
2
l jerf
 1√
2l j
 + 2l2j
exp
− 12l2j
 − 1

 (13)
Finally, the analytic formula for L in (7) is simplified to the following form
L = σ2f
∏
j=−s
J j − σ4f
n∑
p=1
n∑
q=1
kpq ∏
j=−s
I j(xp)
∏
j=−s
I j(xq)
1∫
0
exp
(
−hs
[
(xs − x˜ps )2 + (xs − x˜qs)2
])
dxs
 =
= σ2f
∏
j=−s
J j − σ4f
n∑
p=1
n∑
q=1
kpq ∏
j=−s
I j(xp)
∏
j=−s
I j(xq) Rpqs
 (14)
The derivation of the integral Rpqs is presented separately:
Rpqs =
1∫
0
exp
(
−hs
[
(xs − x˜ps )2 + (xs − x˜qs)2
])
dxs =
= ∗
∣∣∣∣∣a2 + b2 = 12 [(a + b)2 + (a − b)2]
∣∣∣∣∣ =
= exp
[
−hs
2
(
x˜ps − x˜qs
)2] 1∫
0
exp
[
−hs
2
(
2xs − x˜ps − x˜qs
)2]
dxs =
= ∗
∣∣∣∣∣∣∣gs =
√
hs
2
(
2xs − x˜ps − x˜qs
)
, dgs =
√
2hsdxs, gs(0) = −
√
hs
2
(
x˜ps + x˜
q
s
)
, gs(1) =
√
hs
2
(
2 − x˜ps − x˜qs
)∣∣∣∣∣∣∣ =
= exp
[
−hs
2
(
x˜ps − x˜qs
)2] 1√
2hs
gs(1)∫
gs(0)
exp
(
−g2s
)
dgs =
= exp
[
−hs
2
(
x˜ps − x˜qs
)2] 1√
2hs
√
pi
2
(
erf
[
gs(1)
] − erf [gs(0)]) =
=
√
pi
2
ls exp
[
− 1
4l2s
(
x˜ps − x˜qs
)2] (
erf
[
1
2ls
(
2 − x˜ps − x˜qs
)]
+ erf
[
1
2ls
(
x˜ps + x˜
q
s
)])
(15)
This concludes the evaluation of the integral L in (14) with the exact analytic form for the components
I j(x) in (9), J j in (12) and Rpqs in (15).
On the next step to finalize the calculation of (6), we evaluate the second part M in the equation.
M =
1∫
0
∫
[0,1]d−1
∫
[0,1]d−1
m∗(x)m∗(xˆ)dx−sdx′−sdxs =
=
1∫
0
[
zT−s(x)
(
K + σ2
)−1
y˜ zT−s(xˆ)
(
K + σ2
)−1
y˜
]
dxs =
= σ4f
n∑
p=1
n∑
q=1
wp wq ∏
j=−s
I j(xp)
∏
j=−s
I j(xq)
1∫
0
exp
(
−hs
[
(xs − x˜ps )2 + (xs − x˜qs)2
])
dxs
 =
= σ4f
n∑
p=1
n∑
q=1
wp wq ∏
j=−s
I j(xp)
∏
j=−s
I j(xq)Rpqs
 (16)
4
Second component II from (5): Var∗{E(Y)}
Let us return to the equation in (5) and focus on the second term II
Var∗{E(Y)} =
∫
[0,1]d
∫
[0,1]d
k∗(x, x′) dxdx′ =
=
∫
[0,1]d
∫
[0,1]d
K(x, x′)dxdx′ −
∫
[0,1]d
∫
[0,1]d
K(x, X˜)
(
K + σ2
)−1
K(X˜, x′) dxdx′ =
= Z − zT
(
K + σ2
)−1
z (17)
Again we recast the calculation in more simple steps and use the formulas that was presented previously.
For instance, the computation of the integral z in (17) is reduced to the evaluation of I j from the equation
(9).
z(x) =
∫
[0,1]d
K(x, X˜)dx = σ2f
d∏
j=1
I j(x) (18)
Analogously, to compute Z we address to the integral J j that was calculated in (12)
Z =
∫
[0,1]d
∫
[0,1]d
K(x, x′)dxdx′ = σ2f
d∏
j=1
J j (19)
Third component III from (5): E∗{E(Y)}
We address now to the calculation of the final third component in (5).
E∗{E(Y)} =

∫
[0,1]d
K(x, X˜)dx
 [K + σ2I]−1 y˜ = zT (K(X˜, X˜) + σ2)−1 y˜ (18)= σ2f
d∏
j=1
I j
(
K(X˜, X˜) + σ2
)−1
y˜
(20)
This derivation finalize the computation of the numerator from the formula (4), all three components,I, II
and III , from (5) have been derived in analytic form.
Inference for variance in denominator E∗{Var(Y)}
The final step in the calculation of the estimate for the Sobol indices is the derivation of the analytic
formula for the denominator in (4).
E∗{var(Y)} = E∗{E(Y2)} − E∗{E(Y)}2 = E{E∗(Y2)} −
(
Var∗{E(Y)} + (E∗{E(Y)})2
)
=
= E
(
Var∗{Y} + {E∗(Y)}2
)
−
(
Var∗{E(Y)} + (E∗{E(Y)})2
)
= IV − II − III (21)
Similar to the previous routine for the numerator in (4), we divide the calculations in several steps. There
are three components II, III and IV in (21), where II is summarized by formulas (17), (18) and (19),
while III can be calculated using (20). The derivation of the forth component IV is presented in the next
section.
5
Forth component IV from (21): E∗{E(Y)}
Here we perform the last calculation of the forth component in (21).
E{Var∗{Y}} + E{E∗(Y)}2 =
∫
[0,1]d
K(x, x)dx −
∫
[0,1]d
K(x, X˜)
[
K + σ2 
]−1
K(X˜, x) dx+
+
∫
[0,1]d
[
K(x, X˜)
(
K + σ2
)−1
y˜ K(x, X˜)
[
K + σ2 
]−1
y˜ dx
]
=
= σ2f − σ4f
n∑
p=1
n∑
q=1
kpq d∏
s=1
Rpqs
 + σ4f n∑
p=1
n∑
q=1
wp wq d∏
s=1
Rpqs
 , (22)
where the derivation for Rpqs is shown in (15).
6


ARTS repository - University of Groningen

   University of GroningenAnalytical calculation of Sobol sensitivity indices for Gaussian Processes with a squaredexponential covariance functionMilias Argeitis, Andreas; Kurdyaeva, TamaraIMPORTANT NOTE: You are advised to consult the publisher's version (publisher's PDF) if you wish to cite fromit. Please check the document version below.Publication date:2018Link to publication in University of Groningen/UMCG research databaseCitation for published version (APA):Milias Argeitis, A., & Kurdyaeva, T. (2018). Analytical calculation of Sobol sensitivity indices for GaussianProcesses with a squared exponential covariance function.CopyrightOther than for strictly personal use, it is not permitted to download or to forward/distribute the text or part of it without the consent of theauthor(s) and/or copyright holder(s), unless the work is under an open content license (like Creative Commons).The publication may also be distributed here under the terms of Article 25fa of the Dutch Copyright Act, indicated by the “Taverne” license.More information can be found on the University of Groningen website: https://www.rug.nl/library/open-access/self-archiving-pure/taverne-amendment.Take-down policyIf you believe that this document breaches copyright please contact us providing details, and we will remove access to the work immediatelyand investigate your claim.Downloaded from the University of Groningen/UMCG research database (Pure): http://www.rug.nl/research/portal. For technical reasons thenumber of authors shown on this cover page is limited to 10 maximum.Download date: 01-11-2023Analytical calculation of Sobol sensitivity indices for Gaussian Processeswith a squared exponential covariance functionNotation and modeling assumptions• Column vector of random inputs: X = [X1...Xd]T ∈ Rd, X′ = [X′1...X′d]T ∈ Rd,• Column subvector of X without Xs element: X−s = [X1 . . . Xs−1 Xs+1 . . . Xd]T ∈ Rd−1• Column vector consists of elements X′−s and Xs: X̂ = [X′1 . . . Xs . . . X′d]T ∈ Rd• Random inputs Xs, s = 1, ..., d are mutually independent and each input follows a standard uniformdistribution: Xs ∼ U(0, 1), s = 1, ..., d• Function F represents a black-box model of interest: Y = F(X), where Y is the model output• Observed n training data points are denoted by X̃ = (̃x(1), . . . , x̃(n))• Test n∗ data points are denoted by X∗ =(x(1)∗ , . . . , x(n∗)∗)• Vector of function F values at n training points X̃: ỹ =[̃y(1), . . . , ỹ(n)]T, ỹ(i) = F(̃x(i)), i = 1, ..., n• Gaussian process (GP) parametrization f (x) ∼ GP(0, k(x, x′))k(x, x′) = σ2fd∏j=1exp− (x j − x′j)22 l2j = ∗ ∣∣∣∣∣h j = 12 l−2j∣∣∣∣∣= σ2fd∏j=1exp(−h j(x j − x′j)2)(1)• Posterior GP distribution f∗|X̃, ỹ, X∗ ∼ N(m(X∗), v(X∗)), wherem(X∗) =K(X∗, X̃)[K(X̃, X̃) + σ2I]−1ỹ, (2)v(X∗) =K(X∗, X∗) − K(X∗, X̃)[K(X̃, X̃) + σ2I]−1K(X̃, X∗), (3)where the additive term σ2I reflects the fact that the function values are corrupted by additive noisewith variance σ2. Even if σ2 = 0 (as is the case when the function is exactly evaluated), a smallnoise term is typically inserted to regularize the inversion of the covariance matrix K(X̃, X̃).• K = K(X̃, X̃)• kpq for p = 1, ..., n and q = 1, ..., n denotes the entries of the matrix[K + σ2 I]−1• wp, p = 1, ..., n are the elements of the vector[K + σ2 I]−1ỹ• Expectations and variances with respect to the GP posterior distribution are denoted as E∗ andVar∗, respectively• erf(·) denotes the error function.1Sobol indicesThe estimate of the Sobol indices for the output Y with respect to the input XsŜ s =E∗{Var{E (Y |Xs)}}E∗{Var(Y)}, s = 1, ..., d. (4)Inference for variance in numerator E∗{Var{E (Y |Xs)}}Let us first calculate the numerator in the formula for estimates of Sobol indices (4).E∗{Var{E(Y |Xs)}} = E∗{E (E(Y |Xs))2} − E∗{E(Y)}2 == E∗{E (E(Y |Xs))2} − Var∗{E(Y)} − (E∗{E(Y)})2 == I − II − III (5)First component I from (5): E∗{E (E(Y |Xs))2}We derive the formulas for every component of (5) separately. The first term I requires the computationof the following integrals:E∗{E (E(Y |Xs))2} = E[Var∗{E(Y |Xs)} + (E∗{E(Y |Xs)})2]==1∫0∫[0,1]d−1∫[0,1]d−1v(x, x̂) dx−sdx′−sdxs +1∫0∫[0,1]d−1∫[0,1]d−1m∗(x)m∗(x̂)dx−sdx′−sdxs == L + M. (6)For the simplification we also provide the derivation of the formulas for the individual terms L and M in(6) separately.L =1∫0∫[0,1]d−1∫[0,1]d−1v(x, x̂) dx−sdx′−sdxs3==1∫0∫[0,1]d−1∫[0,1]d−1K(x, x̂) dx−sdx′−sdxs −1∫0∫[0,1]d−1∫[0,1]d−1K(x, X̃)(K + σ2)−1K(X̃, x̂) dx−sdx′−sdxs ==1∫0Z−s(x̂)dxs −1∫0[zT−s(x)(K + σ2)−1z−s(x̂)]dxs (7)Now we divide the calculation of the integrals into several parts to make it easier to followz−s(x) =∫[0,1]d−1K(x, X̃) dx−s(1) =σ2f exp(−hs(xs − x̃s)2) ∏j=−s1∫0exp(−h j(x j − x̃ j)2)dx j == σ2f exp(−hs(xs − x̃s)2) ∏j=−sI j (8)2Next we evaluate the integral defined as I j in (8)I j =1∫0exp(−h j(x j − x̃ j)2)dx j = ∗∣∣∣∣∣a j = √h j(x j − x̃ j)∣∣∣∣∣ = 1√h j√h j(1−x̃ j)∫−√h j x̃ jexp(−a2j)da j ==√π21√h j(erf(√h j(1 − x̃ j))+ erf(√h j x̃ j))=√π2l jerf  1√2l j(1 − x̃ j) + erf  1√2l jx̃ j(9)We can substitute the exact formula of I j (9) into the final formula of z−s in (8)z−s(x) = σ2f exp(−12l2s(xs − x̃s)2) ∏j=−s√π2l jerf  1√2l j(1 − x̃ j) + erf  1√2l jx̃ j (10)Finally we calculate the first double integral in (7) which was marked as Z−sZ−s =∫[0,1]d−1∫[0,1]d−1K(x, x̂)dx−sdx′−s = σ2f∏j=−s1∫01∫0exp(−h j(x j − x′j)2)dx j dx′j = σ2f ∏j=−s1∫0I jdx′j == σ2f∏j=−s√π2l j1∫0erf  1√2l j(1 − x′j) + erf  1√2l jx′j dx′j = σ2f ∏j=−sJ j (11)We present the exact analytic calculation of the J j integral in (11)J j =√π2l j1∫0erf 1√2l j(1 − x′j) dx′j + 1∫0erf 1√2l jx′j dx′j == ∗∣∣∣∣∣∣∣b j = 1√2l j (1 − x′j), c j = 1√2l j x′j∣∣∣∣∣∣∣ ==√π2l j−√2l j0∫1/√2l jerf(b j)db j +√2l j1/√2l j∫0erf(c j)dc j == 2√πl2j1/√2l j∫0erf(t j)dt j == ∗∣∣∣∣∣∣u = erf(t j), dv = dt j, du = 2√π exp(−t2j )dt j, v = t j∣∣∣∣∣∣ == 2√πl2jerf(t j)t j∣∣∣∣1/√2l j0−2√π1/√2l j∫0exp(−t2j)t jdt j == 2√πl2j 1√2l jerf 1√2l j + 1√πexp(−t2j) ∣∣∣∣1/√2l j0 == 2√π2l jerf 1√2l j + 2l2jexp− 12l2j − 1 (12)3Then we can write down a closed form to compute Z−s in (11)Z−s = σ2f∏j=−s2√π2l jerf 1√2l j + 2l2jexp− 12l2j − 1 (13)Finally, the analytic formula for L in (7) is simplified to the following formL = σ2f∏j=−sJ j − σ4fn∑p=1n∑q=1kpq ∏j=−sI j(xp)∏j=−sI j(xq)1∫0exp(−hs[(xs − x̃ps )2 + (xs − x̃qs)2])dxs == σ2f∏j=−sJ j − σ4fn∑p=1n∑q=1kpq ∏j=−sI j(xp)∏j=−sI j(xq) Rpqs (14)The derivation of the integral Rpqs is presented separately:Rpqs =1∫0exp(−hs[(xs − x̃ps )2 + (xs − x̃qs)2])dxs == ∗∣∣∣∣∣a2 + b2 = 12 [(a + b)2 + (a − b)2]∣∣∣∣∣ == exp[−hs2(x̃ps − x̃qs)2] 1∫0exp[−hs2(2xs − x̃ps − x̃qs)2]dxs == ∗∣∣∣∣∣∣∣gs =√hs2(2xs − x̃ps − x̃qs), dgs =√2hsdxs, gs(0) = −√hs2(x̃ps + x̃qs), gs(1) =√hs2(2 − x̃ps − x̃qs)∣∣∣∣∣∣∣ == exp[−hs2(x̃ps − x̃qs)2] 1√2hsgs(1)∫gs(0)exp(−g2s)dgs == exp[−hs2(x̃ps − x̃qs)2] 1√2hs√π2(erf[gs(1)]− erf[gs(0)])==√π2ls exp[−14l2s(x̃ps − x̃qs)2] (erf[12ls(2 − x̃ps − x̃qs)]+ erf[12ls(x̃ps + x̃qs)])(15)This concludes the evaluation of the integral L in (14) with the exact analytic form for the componentsI j(x) in (9), J j in (12) and Rpqs in (15).On the next step to finalize the calculation of (6), we evaluate the second part M in the equation.M =1∫0∫[0,1]d−1∫[0,1]d−1m∗(x)m∗(x̂)dx−sdx′−sdxs ==1∫0[zT−s(x)(K + σ2)−1ỹ zT−s(x̂)(K + σ2)−1ỹ]dxs == σ4fn∑p=1n∑q=1wp wq ∏j=−sI j(xp)∏j=−sI j(xq)1∫0exp(−hs[(xs − x̃ps )2 + (xs − x̃qs)2])dxs == σ4fn∑p=1n∑q=1wp wq ∏j=−sI j(xp)∏j=−sI j(xq)Rpqs (16)4Second component II from (5): Var∗{E(Y)}Let us return to the equation in (5) and focus on the second term IIVar∗{E(Y)} =∫[0,1]d∫[0,1]dk∗(x, x′) dxdx′ ==∫[0,1]d∫[0,1]dK(x, x′)dxdx′ −∫[0,1]d∫[0,1]dK(x, X̃)(K + σ2)−1K(X̃, x′) dxdx′ == Z − zT(K + σ2)−1z (17)Again we recast the calculation in more simple steps and use the formulas that was presented previously.For instance, the computation of the integral z in (17) is reduced to the evaluation of I j from the equation(9).z(x) =∫[0,1]dK(x, X̃)dx = σ2fd∏j=1I j(x) (18)Analogously, to compute Z we address to the integral J j that was calculated in (12)Z =∫[0,1]d∫[0,1]dK(x, x′)dxdx′ = σ2fd∏j=1J j (19)Third component III from (5): E∗{E(Y)}We address now to the calculation of the final third component in (5).E∗{E(Y)} =∫[0,1]dK(x, X̃)dx [K + σ2I]−1 ỹ = zT (K(X̃, X̃) + σ2)−1 ỹ (18)= σ2fd∏j=1I j(K(X̃, X̃) + σ2)−1ỹ(20)This derivation finalize the computation of the numerator from the formula (4), all three components,I, IIand III , from (5) have been derived in analytic form.Inference for variance in denominator E∗{Var(Y)}The final step in the calculation of the estimate for the Sobol indices is the derivation of the analyticformula for the denominator in (4).E∗{var(Y)} = E∗{E(Y2)} − E∗{E(Y)}2 = E{E∗(Y2)} −(Var∗{E(Y)} + (E∗{E(Y)})2)== E(Var∗{Y} + {E∗(Y)}2)−(Var∗{E(Y)} + (E∗{E(Y)})2)= IV − II − III (21)Similar to the previous routine for the numerator in (4), we divide the calculations in several steps. Thereare three components II, III and IV in (21), where II is summarized by formulas (17), (18) and (19),while III can be calculated using (20). The derivation of the fourth component IV is presented in thenext section.5Fourth component IV from (21): E∗{E(Y)}Here we perform the last calculation of the fourth component in (21).E{Var∗{Y}} + E{E∗(Y)}2 =∫[0,1]dK(x, x)dx −∫[0,1]dK(x, X̃)[K + σ2 ]−1K(X̃, x) dx++∫[0,1]d[K(x, X̃)(K + σ2)−1ỹ K(x, X̃)[K + σ2 ]−1ỹ dx]== σ2f − σ4fn∑p=1n∑q=1kpq d∏s=1Rpqs + σ4f n∑p=1n∑q=1wp wq d∏s=1Rpqs , (22)where the derivation for Rpqs is shown in (15).6Final formulas for Sobol indicesHere we present the final formulas for Sobol indices’ estimates Ŝ s defined as:Ŝ s =E∗{Var{E (Y |Xs)}}E∗{Var(Y)}, s = 1, ..., d. (4)We start with the analytic expression for the numerator E∗{Var{E (Y |Xs)}} in 4.E∗{Var{E(Y |Xs)}} = E∗{E (E(Y |Xs))2} − Var∗{E(Y)} − (E∗{E(Y)})2 = I − II − III. (5)For the first term I in 5, we getI : E∗{E (E(Y |Xs))2} = σ2f∏j=−sJ j − σ4fn∑p=1n∑q=1(kpq − wp wq)Rpqs ∏j=−sI j(x̃(p))I j(x̃(q)) .whereJ j = 2 l j√π2erf 1√2l j + l jexp− 12l2j − 1 ,Rpqs =√π2ls exp(−14l2s(x̃(p)s − x̃(q)s)2) (erf[12ls(2 − x̃(p)s − x̃(q)s)]+ erf[12ls(x̃(p)s + x̃(q)s)]),andI j(x̃(p))=√π2l jerf  1√2l j(1 − x̃(p)j) + erf  1√2l jx̃(p)j .The second term II in 5 can then be written asII : Var∗{E(Y)} = σ2fd∏j=1J j − σ4f d∏j=1I j(X̃)′ [K + σ2f I]−1  d∏j=1I j(X̃) .Finally, for the third III term of 5 we getIII : (E∗{E(Y)})2 =σ2f d∏j=1I j(K(X̃, X̃) + σ2)−1ỹ2The expression for the denominator E∗{Var(Y)} of 4 has the following form:E∗{Var(Y)} = E∗{E(Y2)} − Var∗{E(Y)} − (E∗{E(Y)})2 = IV − II − III. (21)As we can see, the formula 21 contains one further term compared to the expression of the numerator5, which is given byIV : E∗{E(Y2)} = σ2f − σ4fn∑p=1n∑q=1kpq d∏s=1Rpqs + σ4f n∑p=1n∑q=1wp wq d∏s=1Rpqs .7

English

   University of GroningenAnalytical calculation of Sobol sensitivity indices for Gaussian Processes with a squaredexponential covariance functionMilias Argeitis, Andreas; Kurdyaeva, TamaraIMPORTANT NOTE: You are advised to consult the publisher's version (publisher's PDF) if you wish to cite fromit. Please check the document version below.Publication date:2018Link to publication in University of Groningen/UMCG research databaseCitation for published version (APA):Milias Argeitis, A., & Kurdyaeva, T. (2018). Analytical calculation of Sobol sensitivity indices for GaussianProcesses with a squared exponential covariance function.CopyrightOther than for strictly personal use, it is not permitted to download or to forward/distribute the text or part of it without the consent of theauthor(s) and/or copyright holder(s), unless the work is under an open content license (like Creative Commons).The publication may also be distributed here under the terms of Article 25fa of the Dutch Copyright Act, indicated by the “Taverne” license.More information can be found on the University of Groningen website: https://www.rug.nl/library/open-access/self-archiving-pure/taverne-amendment.Take-down policyIf you believe that this document breaches copyright please contact us providing details, and we will remove access to the work immediatelyand investigate your claim.Downloaded from the University of Groningen/UMCG research database (Pure): http://www.rug.nl/research/portal. For technical reasons thenumber of authors shown on this cover page is limited to 10 maximum.Download date: 21-06-2022Analytical calculation of Sobol sensitivity indices for Gaussian Processeswith a squared exponential covariance functionNotation and modeling assumptions• Column vector of random inputs: X = [X1...Xd]T ∈ Rd, X′ = [X′1...X′d]T ∈ Rd,• Column subvector of X without Xs element: X−s = [X1 . . . Xs−1 Xs+1 . . . Xd]T ∈ Rd−1• Column vector consists of elements X′−s and Xs: X̂ = [X′1 . . . Xs . . . X′d]T ∈ Rd• Random inputs Xs, s = 1, ..., d are mutually independent and each input follows a standard uniformdistribution: Xs ∼ U(0, 1), s = 1, ..., d• Function F represents a black-box model of interest: Y = F(X), where Y is the model output• Observed n training data points are denoted by X̃ = (̃x(1), . . . , x̃(n))• Test n∗ data points are denoted by X∗ =(x(1)∗ , . . . , x(n∗)∗)• Vector of function F values at n training points X̃: ỹ =[̃y(1), . . . , ỹ(n)]T, ỹ(i) = F(̃x(i)), i = 1, ..., n• Gaussian process (GP) parametrization f (x) ∼ GP(0, k(x, x′))k(x, x′) = σ2fd∏j=1exp− (x j − x′j)22 l2j = ∗ ∣∣∣∣∣h j = 12 l−2j∣∣∣∣∣= σ2fd∏j=1exp(−h j(x j − x′j)2)(1)• Posterior GP distribution f∗|X̃, ỹ, X∗ ∼ N(m(X∗), v(X∗)), wherem(X∗) =K(X∗, X̃)[K(X̃, X̃) + σ2I]−1ỹ, (2)v(X∗) =K(X∗, X∗) − K(X∗, X̃)[K(X̃, X̃) + σ2I]−1K(X̃, X∗), (3)where the additive term σ2I reflects the fact that the function values are corrupted by additive noisewith variance σ2. Even if σ2 = 0 (as is the case when the function is exactly evaluated), a smallnoise term is typically inserted to regularize the inversion of the covariance matrix K(X̃, X̃).• K = K(X̃, X̃)• kpq for p = 1, ..., n and q = 1, ..., n denotes the entries of the matrix[K + σ2 I]−1• wp, p = 1, ..., n are the elements of the vector[K + σ2 I]−1ỹ• Expectations and variances with respect to the GP posterior distribution are denoted as E∗ andVar∗, respectively• erf(·) denotes the error function.1Sobol indicesThe estimate of the Sobol indices for the output Y with respect to the input XsŜ s =E∗{Var{E (Y |Xs)}}E∗{Var(Y)}, s = 1, ..., d. (4)Inference for variance in numerator E∗{Var{E (Y |Xs)}}Let us first calculate the numerator in the formula for estimates of Sobol indices (4).E∗{Var{E(Y |Xs)}} = E∗{E (E(Y |Xs))2} − E∗{E(Y)}2 == E∗{E (E(Y |Xs))2} − Var∗{E(Y)} − (E∗{E(Y)})2 == I − II − III (5)First component I from (5): E∗{E (E(Y |Xs))2}We derive the formulas for every component of (5) separately. The first term I requires the computationof the following integrals:E∗{E (E(Y |Xs))2} = E[Var∗{E(Y |Xs)} + (E∗{E(Y |Xs)})2]==1∫0∫[0,1]d−1∫[0,1]d−1v(x, x̂) dx−sdx′−sdxs +1∫0∫[0,1]d−1∫[0,1]d−1m∗(x)m∗(x̂)dx−sdx′−sdxs == L + M. (6)For the simplification we also provide the derivation of the formulas for the individual terms L and M in(6) separately.L =1∫0∫[0,1]d−1∫[0,1]d−1v(x, x̂) dx−sdx′−sdxs3==1∫0∫[0,1]d−1∫[0,1]d−1K(x, x̂) dx−sdx′−sdxs −1∫0∫[0,1]d−1∫[0,1]d−1K(x, X̃)(K + σ2)−1K(X̃, x̂) dx−sdx′−sdxs ==1∫0Z−s(x̂)dxs −1∫0[zT−s(x)(K + σ2)−1z−s(x̂)]dxs (7)Now we divide the calculation of the integrals into several parts to make it easier to followz−s(x) =∫[0,1]d−1K(x, X̃) dx−s(1) =σ2f exp(−hs(xs − x̃s)2) ∏j=−s1∫0exp(−h j(x j − x̃ j)2)dx j == σ2f exp(−hs(xs − x̃s)2) ∏j=−sI j (8)2Next we evaluate the integral defined as I j in (8)I j =1∫0exp(−h j(x j − x̃ j)2)dx j = ∗∣∣∣∣∣a j = √h j(x j − x̃ j)∣∣∣∣∣ = 1√h j√h j(1−x̃ j)∫−√h j x̃ jexp(−a2j)da j ==√π21√h j(erf(√h j(1 − x̃ j))+ erf(√h j x̃ j))=√π2l jerf  1√2l j(1 − x̃ j) + erf  1√2l jx̃ j(9)We can substitute the exact formula of I j (9) into the final formula of z−s in (8)z−s(x) = σ2f exp(−12l2s(xs − x̃s)2) ∏j=−s√π2l jerf  1√2l j(1 − x̃ j) + erf  1√2l jx̃ j (10)Finally we calculate the first double integral in (7) which was marked as Z−sZ−s =∫[0,1]d−1∫[0,1]d−1K(x, x̂)dx−sdx′−s = σ2f∏j=−s1∫01∫0exp(−h j(x j − x′j)2)dx j dx′j = σ2f ∏j=−s1∫0I jdx′j == σ2f∏j=−s√π2l j1∫0erf  1√2l j(1 − x′j) + erf  1√2l jx′j dx′j = σ2f ∏j=−sJ j (11)We present the exact analytic calculation of the J j integral in (11)J j =√π2l j1∫0erf 1√2l j(1 − x′j) dx′j + 1∫0erf 1√2l jx′j dx′j == ∗∣∣∣∣∣∣∣b j = 1√2l j (1 − x′j), c j = 1√2l j x′j∣∣∣∣∣∣∣ ==√π2l j−√2l j0∫1/√2l jerf(b j)db j +√2l j1/√2l j∫0erf(c j)dc j == 2√πl2j1/√2l j∫0erf(t j)dt j == ∗∣∣∣∣∣∣u = erf(t j), dv = dt j, du = 2√π exp(−t2j )dt j, v = t j∣∣∣∣∣∣ == 2√πl2jerf(t j)t j∣∣∣∣1/√2l j0−2√π1/√2l j∫0exp(−t2j)t jdt j == 2√πl2j 1√2l jerf 1√2l j + 1√πexp(−t2j) ∣∣∣∣1/√2l j0 == 2√π2l jerf 1√2l j + 2l2jexp− 12l2j − 1 (12)3Then we can write down a closed form to compute Z−s in (11)Z−s = σ2f∏j=−s2√π2l jerf 1√2l j + 2l2jexp− 12l2j − 1 (13)Finally, the analytic formula for L in (7) is simplified to the following formL = σ2f∏j=−sJ j − σ4fn∑p=1n∑q=1kpq ∏j=−sI j(xp)∏j=−sI j(xq)1∫0exp(−hs[(xs − x̃ps )2 + (xs − x̃qs)2])dxs == σ2f∏j=−sJ j − σ4fn∑p=1n∑q=1kpq ∏j=−sI j(xp)∏j=−sI j(xq) Rpqs (14)The derivation of the integral Rpqs is presented separately:Rpqs =1∫0exp(−hs[(xs − x̃ps )2 + (xs − x̃qs)2])dxs == ∗∣∣∣∣∣a2 + b2 = 12 [(a + b)2 + (a − b)2]∣∣∣∣∣ == exp[−hs2(x̃ps − x̃qs)2] 1∫0exp[−hs2(2xs − x̃ps − x̃qs)2]dxs == ∗∣∣∣∣∣∣∣gs =√hs2(2xs − x̃ps − x̃qs), dgs =√2hsdxs, gs(0) = −√hs2(x̃ps + x̃qs), gs(1) =√hs2(2 − x̃ps − x̃qs)∣∣∣∣∣∣∣ == exp[−hs2(x̃ps − x̃qs)2] 1√2hsgs(1)∫gs(0)exp(−g2s)dgs == exp[−hs2(x̃ps − x̃qs)2] 1√2hs√π2(erf[gs(1)]− erf[gs(0)])==√π2ls exp[−14l2s(x̃ps − x̃qs)2] (erf[12ls(2 − x̃ps − x̃qs)]+ erf[12ls(x̃ps + x̃qs)])(15)This concludes the evaluation of the integral L in (14) with the exact analytic form for the componentsI j(x) in (9), J j in (12) and Rpqs in (15).On the next step to finalize the calculation of (6), we evaluate the second part M in the equation.M =1∫0∫[0,1]d−1∫[0,1]d−1m∗(x)m∗(x̂)dx−sdx′−sdxs ==1∫0[zT−s(x)(K + σ2)−1ỹ zT−s(x̂)(K + σ2)−1ỹ]dxs == σ4fn∑p=1n∑q=1wp wq ∏j=−sI j(xp)∏j=−sI j(xq)1∫0exp(−hs[(xs − x̃ps )2 + (xs − x̃qs)2])dxs == σ4fn∑p=1n∑q=1wp wq ∏j=−sI j(xp)∏j=−sI j(xq)Rpqs (16)4Second component II from (5): Var∗{E(Y)}Let us return to the equation in (5) and focus on the second term IIVar∗{E(Y)} =∫[0,1]d∫[0,1]dk∗(x, x′) dxdx′ ==∫[0,1]d∫[0,1]dK(x, x′)dxdx′ −∫[0,1]d∫[0,1]dK(x, X̃)(K + σ2)−1K(X̃, x′) dxdx′ == Z − zT(K + σ2)−1z (17)Again we recast the calculation in more simple steps and use the formulas that was presented previously.For instance, the computation of the integral z in (17) is reduced to the evaluation of I j from the equation(9).z(x) =∫[0,1]dK(x, X̃)dx = σ2fd∏j=1I j(x) (18)Analogously, to compute Z we address to the integral J j that was calculated in (12)Z =∫[0,1]d∫[0,1]dK(x, x′)dxdx′ = σ2fd∏j=1J j (19)Third component III from (5): E∗{E(Y)}We address now to the calculation of the final third component in (5).E∗{E(Y)} =∫[0,1]dK(x, X̃)dx [K + σ2I]−1 ỹ = zT (K(X̃, X̃) + σ2)−1 ỹ (18)= σ2fd∏j=1I j(K(X̃, X̃) + σ2)−1ỹ(20)This derivation finalize the computation of the numerator from the formula (4), all three components,I, IIand III , from (5) have been derived in analytic form.Inference for variance in denominator E∗{Var(Y)}The final step in the calculation of the estimate for the Sobol indices is the derivation of the analyticformula for the denominator in (4).E∗{var(Y)} = E∗{E(Y2)} − E∗{E(Y)}2 = E{E∗(Y2)} −(Var∗{E(Y)} + (E∗{E(Y)})2)== E(Var∗{Y} + {E∗(Y)}2)−(Var∗{E(Y)} + (E∗{E(Y)})2)= IV − II − III (21)Similar to the previous routine for the numerator in (4), we divide the calculations in several steps. Thereare three components II, III and IV in (21), where II is summarized by formulas (17), (18) and (19),while III can be calculated using (20). The derivation of the fourth component IV is presented in thenext section.5Fourth component IV from (21): E∗{E(Y)}Here we perform the last calculation of the fourth component in (21).E{Var∗{Y}} + E{E∗(Y)}2 =∫[0,1]dK(x, x)dx −∫[0,1]dK(x, X̃)[K + σ2 ]−1K(X̃, x) dx++∫[0,1]d[K(x, X̃)(K + σ2)−1ỹ K(x, X̃)[K + σ2 ]−1ỹ dx]== σ2f − σ4fn∑p=1n∑q=1kpq d∏s=1Rpqs + σ4f n∑p=1n∑q=1wp wq d∏s=1Rpqs , (22)where the derivation for Rpqs is shown in (15).6Final formulas for Sobol indicesHere we present the final formulas for Sobol indices’ estimates Ŝ s defined as:Ŝ s =E∗{Var{E (Y |Xs)}}E∗{Var(Y)}, s = 1, ..., d. (4)We start with the analytic expression for the numerator E∗{Var{E (Y |Xs)}} in 4.E∗{Var{E(Y |Xs)}} = E∗{E (E(Y |Xs))2} − Var∗{E(Y)} − (E∗{E(Y)})2 = I − II − III. (5)For the first term I in 5, we getI : E∗{E (E(Y |Xs))2} = σ2f∏j=−sJ j − σ4fn∑p=1n∑q=1(kpq − wp wq)Rpqs ∏j=−sI j(x̃(p))I j(x̃(q)) .whereJ j = 2 l j√π2erf 1√2l j + l jexp− 12l2j − 1 ,Rpqs =√π2ls exp(−14l2s(x̃(p)s − x̃(q)s)2) (erf[12ls(2 − x̃(p)s − x̃(q)s)]+ erf[12ls(x̃(p)s + x̃(q)s)]),andI j(x̃(p))=√π2l jerf  1√2l j(1 − x̃(p)j) + erf  1√2l jx̃(p)j .The second term II in 5 can then be written asII : Var∗{E(Y)} = σ2fd∏j=1J j − σ4f d∏j=1I j(X̃)′ [K + σ2f I]−1  d∏j=1I j(X̃) .Finally, for the third III term of 5 we getIII : (E∗{E(Y)})2 =σ2f d∏j=1I j(K(X̃, X̃) + σ2)−1ỹ2The expression for the denominator E∗{Var(Y)} of 4 has the following form:E∗{Var(Y)} = E∗{E(Y2)} − Var∗{E(Y)} − (E∗{E(Y)})2 = IV − II − III. (21)As we can see, the formula 21 contains one further term compared to the expression of the numerator5, which is given byIV : E∗{E(Y2)} = σ2f − σ4fn∑p=1n∑q=1kpq d∏s=1Rpqs + σ4f n∑p=1n∑q=1wp wq d∏s=1Rpqs .7

University of Groningen

   University of GroningenAnalytical calculation of Sobol sensitivity indices for Gaussian Processes with a squaredexponential covariance functionMilias Argeitis, Andreas; Kurdyaeva, TamaraIMPORTANT NOTE: You are advised to consult the publisher's version (publisher's PDF) if you wish to cite fromit. Please check the document version below.Publication date:2018Link to publication in University of Groningen/UMCG research databaseCitation for published version (APA):Milias Argeitis, A., & Kurdyaeva, T. (2018). Analytical calculation of Sobol sensitivity indices for GaussianProcesses with a squared exponential covariance function.CopyrightOther than for strictly personal use, it is not permitted to download or to forward/distribute the text or part of it without the consent of theauthor(s) and/or copyright holder(s), unless the work is under an open content license (like Creative Commons).Take-down policyIf you believe that this document breaches copyright please contact us providing details, and we will remove access to the work immediatelyand investigate your claim.Downloaded from the University of Groningen/UMCG research database (Pure): http://www.rug.nl/research/portal. For technical reasons thenumber of authors shown on this cover page is limited to 10 maximum.Download date: 12-11-2019Analytical calculation of Sobol sensitivity indices for Gaussian Processeswith a squared exponential covariance functionNotation and modeling assumptions• Column vector of random inputs: X = [X1...Xd]T ∈ Rd, X′ = [X′1...X′d]T ∈ Rd,• Column subvector of X without Xs element: X−s = [X1 . . . Xs−1 Xs+1 . . . Xd]T ∈ Rd−1• Column vector consists of elements X′−s and Xs: Xˆ = [X′1 . . . Xs . . . X′d]T ∈ Rd• Random inputs Xs, s = 1, ..., d are mutually independent and each input follows a standard uniformdistribution: Xs ∼ U(0, 1), s = 1, ..., d• Function F represents a black-box model of interest: Y = F(X), where Y is the model output• Observed n training data points are denoted by X˜ = (˜x(1), . . . , x˜(n))• Test n∗ data points are denoted by X∗ =(x(1)∗ , . . . , x(n∗)∗)• Vector of function F values at n training points X˜: y˜ =[˜y(1), . . . , y˜(n)]T, y˜(i) = F(˜x(i)), i = 1, ..., n• Gaussian process (GP) parametrization f (x) ∼ GP(0, k(x, x′))k(x, x′) = σ2fd∏j=1exp− (x j − x′j)22 l2j = ∗ ∣∣∣∣∣h j = 12 l−2j∣∣∣∣∣= σ2fd∏j=1exp(−h j(x j − x′j)2)(1)• Posterior GP distribution f∗|X˜, y˜, X∗ ∼ N(m(X∗), v(X∗)), wherem(X∗) =K(X∗, X˜)[K(X˜, X˜) + σ2I]−1y˜, (2)v(X∗) =K(X∗, X∗) − K(X∗, X˜)[K(X˜, X˜) + σ2I]−1K(X˜, X∗), (3)where the additive term σ2I reflects the fact that the function values are corrupted by additive noisewith variance σ2. Even if σ2 = 0 (as is the case when the function is exactly evaluated), a smallnoise term is typically inserted to regularize the inversion of the covariance matrix K(X˜, X˜).• K = K(X˜, X˜)• kpq for p = 1, ..., n and q = 1, ..., n denotes the entries of the matrix[K + σ2 I]−1• wp, p = 1, ..., n are the elements of the vector[K + σ2 I]−1y˜• Expectations and variances with respect to the GP posterior distribution are denoted as E∗ andVar∗, respectively• erf(·) denotes the error function.1Sobol indicesThe estimate of the Sobol indices for the output Y with respect to the input XsSˆ s =E∗{Var{E (Y |Xs)}}E∗{Var(Y)} , s = 1, ..., d. (4)Inference for variance in numerator E∗{Var{E (Y |Xs)}}Let us first calculate the numerator in the formula for estimates of Sobol indices (4).E∗{Var{E(Y |Xs)}} = E∗{E (E(Y |Xs))2} − E∗{E(Y)}2 == E∗{E (E(Y |Xs))2} − Var∗{E(Y)} − (E∗{E(Y)})2 == I − II − III (5)First component I from (5): E∗{E (E(Y |Xs))2}We derive the formulas for every component of (5) separately. The first term I requires the computationof the following integrals:E∗{E (E(Y |Xs))2} = E[Var∗{E(Y |Xs)} + (E∗{E(Y |Xs)})2]==1∫0∫[0,1]d−1∫[0,1]d−1v(x, xˆ) dx−sdx′−sdxs +1∫0∫[0,1]d−1∫[0,1]d−1m∗(x)m∗(xˆ)dx−sdx′−sdxs == L + M. (6)For the simplification we also provide the derivation of the formulas for the individual terms L and M in(6) separately.L =1∫0∫[0,1]d−1∫[0,1]d−1v(x, xˆ) dx−sdx′−sdxs3==1∫0∫[0,1]d−1∫[0,1]d−1K(x, xˆ) dx−sdx′−sdxs −1∫0∫[0,1]d−1∫[0,1]d−1K(x, X˜)(K + σ2)−1K(X˜, xˆ) dx−sdx′−sdxs ==1∫0Z−s(xˆ)dxs −1∫0[zT−s(x)(K + σ2)−1z−s(xˆ)]dxs (7)Now we divide the calculation of the integrals into several parts to make it easier to followz−s(x) =∫[0,1]d−1K(x, X˜) dx−s(1) =σ2f exp(−hs(xs − x˜s)2) ∏j=−s1∫0exp(−h j(x j − x˜ j)2)dx j == σ2f exp(−hs(xs − x˜s)2) ∏j=−sI j (8)2Next we evaluate the integral defined as I j in (8)I j =1∫0exp(−h j(x j − x˜ j)2)dx j = ∗∣∣∣∣∣a j = √h j(x j − x˜ j)∣∣∣∣∣ = 1√h j√h j(1−x˜ j)∫−√h j x˜ jexp(−a2j)da j ==√pi21√h j(erf(√h j(1 − x˜ j))+ erf(√h j x˜ j))=√pi2l jerf  1√2l j(1 − x˜ j) + erf  1√2l jx˜ j(9)We can substitute the exact formula of I j (9) into the final formula of z−s in (8)z−s(x) = σ2f exp(− 12l2s(xs − x˜s)2) ∏j=−s√pi2l jerf  1√2l j(1 − x˜ j) + erf  1√2l jx˜ j (10)Finally we calculate the first double integral in (7) which was marked as Z−sZ−s =∫[0,1]d−1∫[0,1]d−1K(x, xˆ)dx−sdx′−s = σ2f∏j=−s1∫01∫0exp(−h j(x j − x′j)2)dx j dx′j = σ2f ∏j=−s1∫0I jdx′j == σ2f∏j=−s√pi2l j1∫0erf  1√2l j(1 − x′j) + erf  1√2l jx′j dx′j = σ2f ∏j=−sJ j (11)We present the exact analytic calculation of the J j integral in (11)J j =√pi2l j1∫0erf 1√2l j(1 − x′j) dx′j + 1∫0erf 1√2l jx′j dx′j == ∗∣∣∣∣∣∣∣b j = 1√2l j (1 − x′j), c j = 1√2l j x′j∣∣∣∣∣∣∣ ==√pi2l j−√2l j0∫1/√2l jerf(b j)db j +√2l j1/√2l j∫0erf(c j)dc j == 2√pil2j1/√2l j∫0erf(t j)dt j == ∗∣∣∣∣∣∣u = erf(t j), dv = dt j, du = 2√pi exp(−t2j )dt j, v = t j∣∣∣∣∣∣ == 2√pil2jerf(t j)t j∣∣∣∣1/√2l j0− 2√pi1/√2l j∫0exp(−t2j)t jdt j == 2√pil2j 1√2l jerf 1√2l j + 1√piexp(−t2j) ∣∣∣∣1/√2l j0 == 2√pi2l jerf 1√2l j + 2l2jexp− 12l2j − 1 (12)3Then we can write down a closed form to compute Z−s in (11)Z−s = σ2f∏j=−s2√pi2l jerf 1√2l j + 2l2jexp− 12l2j − 1 (13)Finally, the analytic formula for L in (7) is simplified to the following formL = σ2f∏j=−sJ j − σ4fn∑p=1n∑q=1kpq ∏j=−sI j(xp)∏j=−sI j(xq)1∫0exp(−hs[(xs − x˜ps )2 + (xs − x˜qs)2])dxs == σ2f∏j=−sJ j − σ4fn∑p=1n∑q=1kpq ∏j=−sI j(xp)∏j=−sI j(xq) Rpqs (14)The derivation of the integral Rpqs is presented separately:Rpqs =1∫0exp(−hs[(xs − x˜ps )2 + (xs − x˜qs)2])dxs == ∗∣∣∣∣∣a2 + b2 = 12 [(a + b)2 + (a − b)2]∣∣∣∣∣ == exp[−hs2(x˜ps − x˜qs)2] 1∫0exp[−hs2(2xs − x˜ps − x˜qs)2]dxs == ∗∣∣∣∣∣∣∣gs =√hs2(2xs − x˜ps − x˜qs), dgs =√2hsdxs, gs(0) = −√hs2(x˜ps + x˜qs), gs(1) =√hs2(2 − x˜ps − x˜qs)∣∣∣∣∣∣∣ == exp[−hs2(x˜ps − x˜qs)2] 1√2hsgs(1)∫gs(0)exp(−g2s)dgs == exp[−hs2(x˜ps − x˜qs)2] 1√2hs√pi2(erf[gs(1)] − erf [gs(0)]) ==√pi2ls exp[− 14l2s(x˜ps − x˜qs)2] (erf[12ls(2 − x˜ps − x˜qs)]+ erf[12ls(x˜ps + x˜qs)])(15)This concludes the evaluation of the integral L in (14) with the exact analytic form for the componentsI j(x) in (9), J j in (12) and Rpqs in (15).On the next step to finalize the calculation of (6), we evaluate the second part M in the equation.M =1∫0∫[0,1]d−1∫[0,1]d−1m∗(x)m∗(xˆ)dx−sdx′−sdxs ==1∫0[zT−s(x)(K + σ2)−1y˜ zT−s(xˆ)(K + σ2)−1y˜]dxs == σ4fn∑p=1n∑q=1wp wq ∏j=−sI j(xp)∏j=−sI j(xq)1∫0exp(−hs[(xs − x˜ps )2 + (xs − x˜qs)2])dxs == σ4fn∑p=1n∑q=1wp wq ∏j=−sI j(xp)∏j=−sI j(xq)Rpqs (16)4Second component II from (5): Var∗{E(Y)}Let us return to the equation in (5) and focus on the second term IIVar∗{E(Y)} =∫[0,1]d∫[0,1]dk∗(x, x′) dxdx′ ==∫[0,1]d∫[0,1]dK(x, x′)dxdx′ −∫[0,1]d∫[0,1]dK(x, X˜)(K + σ2)−1K(X˜, x′) dxdx′ == Z − zT(K + σ2)−1z (17)Again we recast the calculation in more simple steps and use the formulas that was presented previously.For instance, the computation of the integral z in (17) is reduced to the evaluation of I j from the equation(9).z(x) =∫[0,1]dK(x, X˜)dx = σ2fd∏j=1I j(x) (18)Analogously, to compute Z we address to the integral J j that was calculated in (12)Z =∫[0,1]d∫[0,1]dK(x, x′)dxdx′ = σ2fd∏j=1J j (19)Third component III from (5): E∗{E(Y)}We address now to the calculation of the final third component in (5).E∗{E(Y)} =∫[0,1]dK(x, X˜)dx [K + σ2I]−1 y˜ = zT (K(X˜, X˜) + σ2)−1 y˜ (18)= σ2fd∏j=1I j(K(X˜, X˜) + σ2)−1y˜(20)This derivation finalize the computation of the numerator from the formula (4), all three components,I, IIand III , from (5) have been derived in analytic form.Inference for variance in denominator E∗{Var(Y)}The final step in the calculation of the estimate for the Sobol indices is the derivation of the analyticformula for the denominator in (4).E∗{var(Y)} = E∗{E(Y2)} − E∗{E(Y)}2 = E{E∗(Y2)} −(Var∗{E(Y)} + (E∗{E(Y)})2)== E(Var∗{Y} + {E∗(Y)}2)−(Var∗{E(Y)} + (E∗{E(Y)})2)= IV − II − III (21)Similar to the previous routine for the numerator in (4), we divide the calculations in several steps. Thereare three components II, III and IV in (21), where II is summarized by formulas (17), (18) and (19),while III can be calculated using (20). The derivation of the fourth component IV is presented in thenext section.5Fourth component IV from (21): E∗{E(Y)}Here we perform the last calculation of the fourth component in (21).E{Var∗{Y}} + E{E∗(Y)}2 =∫[0,1]dK(x, x)dx −∫[0,1]dK(x, X˜)[K + σ2 ]−1K(X˜, x) dx++∫[0,1]d[K(x, X˜)(K + σ2)−1y˜ K(x, X˜)[K + σ2 ]−1y˜ dx]== σ2f − σ4fn∑p=1n∑q=1kpq d∏s=1Rpqs + σ4f n∑p=1n∑q=1wp wq d∏s=1Rpqs , (22)where the derivation for Rpqs is shown in (15).6Final formulas for Sobol indicesHere we present the final formulas for Sobol indices’ estimates Sˆ s defined as:Sˆ s =E∗{Var{E (Y |Xs)}}E∗{Var(Y)} , s = 1, ..., d. (4)We start with the analytic expression for the numerator E∗{Var{E (Y |Xs)}} in 4.E∗{Var{E(Y |Xs)}} = E∗{E (E(Y |Xs))2} − Var∗{E(Y)} − (E∗{E(Y)})2 = I − II − III. (5)For the first term I in 5, we getI : E∗{E (E(Y |Xs))2} = σ2f∏j=−sJ j − σ4fn∑p=1n∑q=1(kpq − wp wq)Rpqs ∏j=−sI j(x˜(p))I j(x˜(q)) .whereJ j = 2 l j√pi2erf 1√2l j + l jexp− 12l2j − 1 ,Rpqs =√pi2ls exp(− 14l2s(x˜(p)s − x˜(q)s)2) (erf[12ls(2 − x˜(p)s − x˜(q)s)]+ erf[12ls(x˜(p)s + x˜(q)s)]),andI j(x˜(p))=√pi2l jerf  1√2l j(1 − x˜(p)j) + erf  1√2l jx˜(p)j .The second term II in 5 can then be written asII : Var∗{E(Y)} = σ2fd∏j=1J j − σ4f d∏j=1I j(X˜)′ [K + σ2f I]−1  d∏j=1I j(X˜) .Finally, for the third III term of 5 we getIII : (E∗{E(Y)})2 =σ2f d∏j=1I j(K(X˜, X˜) + σ2)−1y˜2The expression for the denominator E∗{Var(Y)} of 4 has the following form:E∗{Var(Y)} = E∗{E(Y2)} − Var∗{E(Y)} − (E∗{E(Y)})2 = IV − II − III. (21)As we can see, the formula 21 contains one further term compared to the expression of the numerator5, which is given byIV : E∗{E(Y2)} = σ2f − σ4fn∑p=1n∑q=1kpq d∏s=1Rpqs + σ4f n∑p=1n∑q=1wp wq d∏s=1Rpqs .7

https://pure.rug.nl/ws/files/65061085/sobols_indices_v2.pdf

Analytical calculation of Sobol sensitivity indices for Gaussian Processes with a squared exponential covariance function

Abstract

Similar works

Full text

Available Versions

NARCIS

University of Groningen Research Database

ARTS repository - University of Groningen

ARTS repository - University of Groningen

University of Groningen