【MATLAB第107期】基于MATLAB的Morris全局敏感性分析模型(无目标函数)
一、原理介绍
1.基本原理:
Morris方法采用概率均匀抽样的方式估计每个模型输入因子在输出结果中的重要性,通过比较系统在不同输入参数值上的输出结果变化来度量每个参数的敏感性 。
该方法的原理与局部敏感性分析中的“一次变化一个因子”方法相似,通过计算所谓的基本效应(elementary effect, EE)来衡量每个参数对模型输出的全序影响。
如式(1)所示, 对初始状态进行多次反复采样, 求各次采样的基效应,
并用式(2)求其均值 μi, 衡量每个参数对模型输出的全序影响。
在具体实施过程中,Morris方法通过对每个输入变量进行随机抽样(本文使用LHC抽样,并提供多种抽样方法如sobol /halton),并在这些抽样点上加上随机扰动来分析输入变量对输出的影响。计算每个输入变量的平均值和标准差,然后计算变量在不同取值下的输出变化率。通过比较各个变量的输出变化率来确定其对输出的敏感性。
Morris方法将每个输入变量分为几个等级,通过随机选取初始值、加上随机扰动得到新值、计算新值对应的输出值及其变化量,重复这一过程以获得一定数量的变化量,再计算这些变化量的平均值和标准差。最终,根据变量的平均变化量和标准差来计算变量的元素敏感度和总体敏感度,从而确定不同输入变量的相对重要性。
2.算法步骤:
进行参数采样,通常采用拉丁超立方抽样方法以保证参数的均匀分布 。
对每个候选参数向量进行扰动计算,即改变其中一个参数的值按照一定的步长。
计算扰动后的输出值,并使用所谓的“元距离”(Elementary Effect)来度量输出结果的变化 。
重复计算以提高结果的准确性,多次重复上述步骤以获取一组灵敏度度量值 。
灵敏度度量:Morris方法通过计算每个参数的平均元距离和标准差来评估参数的敏感性,这些度量值反映了参数变化对输出结果的影响程度 。
二、morris方法与sobol方法对比
Morris方法和Sobol方法都是全局敏感性分析的常用技术,它们各自有不同的特点和应用场景:
Morris方法:上文已介绍
Sobol方法:Sobol方法基于方差分解,通过计算每个参数或参数组合对输出方差的份额来确定其对输出的敏感性。这种方法可以识别参数之间的相互作用,并量化它们对输出总方差的贡献。Sobol方法通常采用蒙特卡洛随机采样法生成样本 。
两种方法的主要差异包括:
抽样方法:Morris方法通常使用分层抽样,如拉丁超立方抽样,而Sobol方法则多采用随机采样。
计算成本:Morris方法通常计算成本较低,适合初步筛选重要的输入参数,而Sobol方法计算成本较高,但可以提供更全面和深入的敏感性分析。
结果解释:Morris方法的结果较容易解释,通常使用均值和标准差来表示参数的敏感性;Sobol方法则提供了一阶和总效应指数,可以区分参数的独立影响和相互作用。
适用性:Morris方法适用于模型输入参数数量较多的情况,而Sobol方法适用于需要深入理解参数影响和参数间相互作用的情况。
三、无目标函数实现思路
借鉴往期文章, 如用BP神经网络对输入与输出数据进行非线性函数拟合,通过LHC函数抽样参数,调用BP训练模型进行计算。
nfac = 7; % 输入参数数量
npoint = 1000; % 每个参数的采样点数
VarMin=[137 0 0 160 4.4 708 650];%各个参数下限
VarMax=[374 193 260 240 19 1049.90 902];%各个参数上限
四、代码展示
clear all
% 调用Morris敏感性分析函数
nfac = 7; % 输入参数数量
npoint = 1000; % 每个参数的采样点数
VarMin=[137 0 0 160 4.4 708 650];%各个参数下限
VarMax=[374 193 260 240 19 1049.90 902];%各个参数上限
[mu, order] = morris(@fun, VarMin, VarMax,nfac, npoint);
% 输出结果
disp('敏感度指数(降序排列):');
disp(mu);
disp('参数顺序(按敏感度):');
disp(order);
figure()
bar([order]);
xlabel('输入变量')
ylabel('评价指标')
title('morris全局敏感性分析(无目标函数)')
function [mu, order] = morris(model_function,VarMin, VarMax,nfac, npoint)
% 1.设置参数变化量
% 2.初始化输出矩阵
% 3.LHC生成随机采样点
% 4.计算基础输出
% 5.计算每个参数的基本效应
% 5.1正向变化
% 5.2负向变化
% 6.计算敏感度指数
五、代码获取
1.阅读首页置顶文章
2.关注CSDN
3.根据自动回复消息,回复“107期”以及相应指令,即可获取对应下载方式。