sobol全局灵敏性分析
合集下载
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
Sobol 全局灵敏性分析
最近在研究全局敏感性分析方法中的 Sobol 方法，看了一些国内的论文，发现一个通病，就是 公式一挂就可以得出结果了，真心觉得这种论文很“恶心”，主要原因是自己看不太懂。直到在维 基百科上面找到了这种方法的详细解释，今天我们就根据网页上的步骤用一个例子来走一遍。
1. 假设现在有一个函数 : Y = sin(x1) + 7 (sin(x2 ))2 + 0.1 x34 siN*D 的矩阵 ABi (i = 1,2, …,D) ，即用矩阵 B 中的第 i 列替换矩阵 A 的第 i 列，以本体为例：
0.5 0.5 0.5
0.5 0.5 0.5
0.5 0.5 0.5
0.25 0.25 0.25
0.75 0.75 0.25
0.75 0.25 0.75
AB1
AB2
AB3
1.11036605855
3.50765176932
0.67596163481
YA
YB 3.50765176932
1.11036605855 YAB1
3.95562603058
1.31095036329
1.70665120003
1.71834424557
2.09136387768
2.09136387768
YA 和 Y B 即
2.09136387768 1.11036605855 3.50765176932 1.31095036329 Y 2.09136387768 3.50765176932 1.11036605855 1.70665120003
这样我们就可以通过上述的公式求解 x1 的一介影响指数。下面是求解过程。 8.Y 的均值等于 mean(Y)=(2.09136387768+1.11036605855+3.50765176932+ 1.31095036329+2.09136387768+3.50765
原理在这就不说了，大家可以自行谷歌。为了方便讲解例子我们设置采样的样本数为
4(N = 4) ，自
变量数目为 3(D = 3) 。我们按照上述网页的步骤。
4. 生成 N * 2D （即 4 行 6 列）的样本矩阵。这个就是我们 Sobol sequence 做的事情。这边我们生
成的矩阵为：
0.5 0.5 0.5 0.5 0.5 0.5 0.75 0.25 0.25 0.25 0.75 0.75 m 0.25 0.75 0.75 0.75 0.25 0.25 0.375 0.375 0.625 0.875 0.375 0.125
0.06941973968
0.06941973968，那么 ST1
0.0831043122 0.835332581542
至此 x1 的一介影响指数和全局影响指数都可以求出来了，然后 x2、x3 的一介影响指数与全局影 响指数的求解过程和上面一样。
当然在实际的求解过程中肯定不会用这么少的样本量的，但是求解过程是一样的。希望这篇文 章会对做敏感性分析的小伙伴们有些帮助。
其中 :
1 VarXi EX : i (Y Xi ) N
N
j 1 f (B) j
( f (AB i ) j
f (A) j )
EX: i (VarXi (Y X : i ))
1 2N
N
j
(
1
f (A)
j
2
f (AB i ) j )
Var (Y) Var (YA YB ) 这里的“加”不是普通意义上的相加而是构造了一个新的矩阵包括
0.75 0.75 0.75
0.25 0.25 0.75
0.25 0.75 0.25
0.875 0.375 0.625
0.375 0.375 0.625
0.375 0.375 0.125
6. 经过这三步我们构造了 A、B、AB1、AB2、AB3 这五个矩阵，这样我们就有 (D + 2) * N ( 即 20) 组
其视作一个黑盒子，只有输入和输出，这时我们对其进行敏感性分析就很有必要了。经过敏感性分
析我们就能找出对结果影响较大的参数。这样对于调整结果是很有帮助的。
3. 接着上面的例子，首先我们得根据三个自变量的范围进行采样，这边采样的方法一般都是蒙特卡
洛采样以及一系列基于蒙特卡洛采样的变种，这个例子中我们采用了 Sobol sequence ，具体的采样
176932+1.11036605855+1.70665120003)/8 = 2.0545456218 Var(Y) = 0.835332581542
1 VarX 1
N
N
Y (Y j 1 Bj
AB1 j
YAj )
本例中 VarX1就等于 (2.09136387768 * (2.09136387768 - 2.09136387768) + 3.50765176932* (0.67596163481 - 1.11036605855) + 1.11036605855 * (3.95562603058 - 3.50765176932) +
1.70665120003 * (1.71834424557 - 1.31095036329) )/4 = -0.0827611931562 那么
0.0827611931562
S1
0.09907573 0.835332581542
1 EX: i 2N
N
2
j 1 (YAj YAB1 j )
本例中 EX~1 就等于 ((2.09136387768 - 2.09136387768)^2 + (1.11036605855 - 0.67596163481)^2 + (3.50765176932 - 3.95562603058)^2 + (1.31095036329 - 1.71834424557)^2 )/(2 * 4)=
x1、 x2、x3输入数据，因此我们将有 20 组 Y 值。将上述的数据带入函
Y = sin(x1) + 7 (sin(x2 ))2 + 0.1 x34 sin(x1) ，这里详细的计算过程就不描述了。根据输入我们得出
对应的 Y 值矩阵。
2.09136387768
2.09136387768
2.09136387768
x1、x2、 x3 三个自变量对应变量 Y 有影响。
2. 然后一般会给这三个参数一个 取值范围 ，这里假设三个自变量的取值范围都设为 [0,1] 。
敏感性分析的目的就是求取这三个参数对于 Y 值得贡献。当然我们这边可能有人一下子就可以
分析出那个参数对于 Y 值影响最大，但是在解决实际问题时，这个函数一般都是未知，我们只能将
5、 将矩阵的前 D 列设置为矩阵 A，后 D列设置为 B 列，在我们的例子中就是矩阵 m的前 3 列设置 为矩阵 A，后 3 列设置为矩阵 B。
0.5 0.5 0.5 0.75 0.25 0.25 A 0.25 0.75 0.75 0.375 0.375 0.625
0.5 0.5 0.5 0.25 0.75 0.75 B 0.75 0.25 0.25 0.875 0.375 0.125
3.93432481933
1.1316672698
YAB2 0.683693008537 YAB3 3.49992039559
1.31095036329
1.30537043023
7. 根据一介影响指数公式：
Si VarXi EX : i (Y Xi ) Var (Y)
总效应指数 : STi
EX : i (VarXi (Y X i )) Var (Y)
最近在研究全局敏感性分析方法中的 Sobol 方法，看了一些国内的论文，发现一个通病，就是 公式一挂就可以得出结果了，真心觉得这种论文很“恶心”，主要原因是自己看不太懂。直到在维 基百科上面找到了这种方法的详细解释，今天我们就根据网页上的步骤用一个例子来走一遍。
1. 假设现在有一个函数 : Y = sin(x1) + 7 (sin(x2 ))2 + 0.1 x34 siN*D 的矩阵 ABi (i = 1,2, …,D) ，即用矩阵 B 中的第 i 列替换矩阵 A 的第 i 列，以本体为例：
0.5 0.5 0.5
0.5 0.5 0.5
0.5 0.5 0.5
0.25 0.25 0.25
0.75 0.75 0.25
0.75 0.25 0.75
AB1
AB2
AB3
1.11036605855
3.50765176932
0.67596163481
YA
YB 3.50765176932
1.11036605855 YAB1
3.95562603058
1.31095036329
1.70665120003
1.71834424557
2.09136387768
2.09136387768
YA 和 Y B 即
2.09136387768 1.11036605855 3.50765176932 1.31095036329 Y 2.09136387768 3.50765176932 1.11036605855 1.70665120003
这样我们就可以通过上述的公式求解 x1 的一介影响指数。下面是求解过程。 8.Y 的均值等于 mean(Y)=(2.09136387768+1.11036605855+3.50765176932+ 1.31095036329+2.09136387768+3.50765
原理在这就不说了，大家可以自行谷歌。为了方便讲解例子我们设置采样的样本数为
4(N = 4) ，自
变量数目为 3(D = 3) 。我们按照上述网页的步骤。
4. 生成 N * 2D （即 4 行 6 列）的样本矩阵。这个就是我们 Sobol sequence 做的事情。这边我们生
成的矩阵为：
0.5 0.5 0.5 0.5 0.5 0.5 0.75 0.25 0.25 0.25 0.75 0.75 m 0.25 0.75 0.75 0.75 0.25 0.25 0.375 0.375 0.625 0.875 0.375 0.125
0.06941973968
0.06941973968，那么 ST1
0.0831043122 0.835332581542
至此 x1 的一介影响指数和全局影响指数都可以求出来了，然后 x2、x3 的一介影响指数与全局影 响指数的求解过程和上面一样。
当然在实际的求解过程中肯定不会用这么少的样本量的，但是求解过程是一样的。希望这篇文 章会对做敏感性分析的小伙伴们有些帮助。
其中 :
1 VarXi EX : i (Y Xi ) N
N
j 1 f (B) j
( f (AB i ) j
f (A) j )
EX: i (VarXi (Y X : i ))
1 2N
N
j
(
1
f (A)
j
2
f (AB i ) j )
Var (Y) Var (YA YB ) 这里的“加”不是普通意义上的相加而是构造了一个新的矩阵包括
0.75 0.75 0.75
0.25 0.25 0.75
0.25 0.75 0.25
0.875 0.375 0.625
0.375 0.375 0.625
0.375 0.375 0.125
6. 经过这三步我们构造了 A、B、AB1、AB2、AB3 这五个矩阵，这样我们就有 (D + 2) * N ( 即 20) 组
其视作一个黑盒子，只有输入和输出，这时我们对其进行敏感性分析就很有必要了。经过敏感性分
析我们就能找出对结果影响较大的参数。这样对于调整结果是很有帮助的。
3. 接着上面的例子，首先我们得根据三个自变量的范围进行采样，这边采样的方法一般都是蒙特卡
洛采样以及一系列基于蒙特卡洛采样的变种，这个例子中我们采用了 Sobol sequence ，具体的采样
176932+1.11036605855+1.70665120003)/8 = 2.0545456218 Var(Y) = 0.835332581542
1 VarX 1
N
N
Y (Y j 1 Bj
AB1 j
YAj )
本例中 VarX1就等于 (2.09136387768 * (2.09136387768 - 2.09136387768) + 3.50765176932* (0.67596163481 - 1.11036605855) + 1.11036605855 * (3.95562603058 - 3.50765176932) +
1.70665120003 * (1.71834424557 - 1.31095036329) )/4 = -0.0827611931562 那么
0.0827611931562
S1
0.09907573 0.835332581542
1 EX: i 2N
N
2
j 1 (YAj YAB1 j )
本例中 EX~1 就等于 ((2.09136387768 - 2.09136387768)^2 + (1.11036605855 - 0.67596163481)^2 + (3.50765176932 - 3.95562603058)^2 + (1.31095036329 - 1.71834424557)^2 )/(2 * 4)=
x1、 x2、x3输入数据，因此我们将有 20 组 Y 值。将上述的数据带入函
Y = sin(x1) + 7 (sin(x2 ))2 + 0.1 x34 sin(x1) ，这里详细的计算过程就不描述了。根据输入我们得出
对应的 Y 值矩阵。
2.09136387768
2.09136387768
2.09136387768
x1、x2、 x3 三个自变量对应变量 Y 有影响。
2. 然后一般会给这三个参数一个 取值范围 ，这里假设三个自变量的取值范围都设为 [0,1] 。
敏感性分析的目的就是求取这三个参数对于 Y 值得贡献。当然我们这边可能有人一下子就可以
分析出那个参数对于 Y 值影响最大，但是在解决实际问题时，这个函数一般都是未知，我们只能将
5、 将矩阵的前 D 列设置为矩阵 A，后 D列设置为 B 列，在我们的例子中就是矩阵 m的前 3 列设置 为矩阵 A，后 3 列设置为矩阵 B。
0.5 0.5 0.5 0.75 0.25 0.25 A 0.25 0.75 0.75 0.375 0.375 0.625
0.5 0.5 0.5 0.25 0.75 0.75 B 0.75 0.25 0.25 0.875 0.375 0.125
3.93432481933
1.1316672698
YAB2 0.683693008537 YAB3 3.49992039559
1.31095036329
1.30537043023
7. 根据一介影响指数公式：
Si VarXi EX : i (Y Xi ) Var (Y)
总效应指数 : STi
EX : i (VarXi (Y X i )) Var (Y)