企业官网建设流程全解析

简介：面向需要开展全局敏感性分析的科研与工程人员，这份MATLAB代码包实现了基于快速傅里叶变换的eFAST方法，可量化各输入参数对模型输出的主效应与交互效应，支持用户自定义参数分布与模型函数。压缩包内共8个m文件，整体仅6KB，涵盖参数设置、频率选择、模型封装、敏感度指数计算与显著性检验等核心环节，尤其包含针对常微分方程（ODE）模型的调用示例，便于直接套用或二次开发。eFAST相比仅改变单个参数的局部方法，能同时评估参数间交互作用，在保证精度的同时有效降低计算成本。代码结构清晰，通过调整参数设置脚本中的采样策略与重复次数，即可适配不同复杂度的模型，帮助识别关键参数、降低调优成本。已有568人学习下载，适合从事复杂系统建模与不确定性分析的研究者参考。

1. 为什么一上来就做全局敏感性分析，而不是先调参

一个人接手带 8 个参数的 ODE 模型时，最忙乱的环节不是写方程，而是判断哪些参数值得被精确估计。之前调一个带时滞的生理模型，用lsqnonlin直接拟合，梯度总是被几个输出几乎不变的参数带偏，最后先做了一次全局敏感性分析才定位到两个主导参数。eFAST（扩展傅里叶幅度灵敏度检验）的核心思路是给每个参数分配一个互质频率，沿搜索曲线同时扰动所有参数，再用快速傅里叶变换把输出方差按频率拆开，这样不仅能看到单个参数的主效应，还能暴露参数间的交互作用。这套 eFAST MATLAB code 把采样、ODE 求解、谱分解、检出限判断串成一条流水线，入口在Parameter_settings_EFAST.m和ODE_efast.m两个文件，适合手里已有可跑的 ODE 模型、但还没搞清参数优先级的建模与仿真人员。

2. eFAST 频率分配与方差分解：先看数学再写代码

2.1 互质频率表与离散采样点数

eFAST 用连续变量s ∈ [0, 2π]作为搜索曲线的自变量。第i个参数在第s个采样点处的取值是：

x_i(s) = 0.5 + (1/π) · arcsin(sin(ω_i · s + φ_i))

这个映射把参数值固定在[0,1]，之后再用任意分布的反函数转成真实参数范围。关键是每个参数必须有一个独立的角频率ω_i，而且任意两个频率的最大公约数必须是 1，否则傅里叶谱里高阶谐波会重叠，你根本分不清哪个峰属于哪个参数。

以下是常用互质频率表的示例，实际代码中SETFREQ.m会按参数个数从类似的表里取前K个：

表中M是谐波阶数，默认取 4，表示把主频附近的 1 到 4 次谐波能量都计入对应参数的边际贡献。采样点数由最大频率决定：N = 2Mω_max + 1。这个公式来自奈奎斯特条件，采样点太少，高频成分会折叠到低频区域，敏感性指数直接被污染。

2.2 曲线采样映射与参数分布自定义

在eFASTparameterdist.m中，代码做的就是把上面的解析式转换成实际参数行向量。我一般只关心两个输入：频率数组freq和分布类型dist。

% efast_sample_curve.m —— 单条搜索曲线生成示意 s = linspace(0, 2*pi, N)'; % 等间隔采样点，N 由最大频率决定 phi = rand * 2 * pi; % 随机相位偏移，重采样时改变 X = zeros(N, K); % 每行是一组参数样本 for k = 1:K u = 0.5 + asin(sin(freq(k) * s + phi)) / pi; % 标准 eFAST 曲线，值域 [0,1] if strcmp(dist{k}, 'uniform') X(:, k) = pmin(k) + (pmax(k) - pmin(k)) * u; elseif strcmp(dist{k}, 'normal') % 把均匀分布转换成指定均值方差的正态样本 X(:, k) = norminv(u, mu(k), sigma(k)); end end

这里u不是普通随机数，它是一条周期性扫描曲线的输出，因此后续 FFT 才能把方差归因到对应频率。phi是曲线随机相位，同一频率下换phi就得到一条新曲线；ODE_efast.m会重复多次phi，用来求敏感性指数的标准差和置信区间。

2.3 FFT 谱分解与一阶敏感性指数

对一条曲线的输出y，去掉均值后做 FFT，再把单侧功率谱的每个频点能量除以总能量，得到该频率对输出方差的贡献比例。

% fft_variance.m —— 一条曲线上的谱分解 Y = fft(y - mean(y)); % 去掉直流分量 power = abs(Y(2:floor(N/2)+1)).^2; % 单侧功率谱 S_total = sum(power); for i = 1:K idx = freq(i) * (1:M); % 主频的 1 到 M 次谐波对应位置 S_i(i) = sum(power(idx)) / S_total; % 一阶敏感性指数 end

代码中idx计算的是第i个参数的主频带上所有谐波位置。S_i是该参数的贡献占总方差的百分比。值得注意的是，总效应指数不能只在单条曲线上用频带法硬拆，因为交互效应会分布在其他参数的频率组合上；常见实现是把phi重新随机化，做多组独立曲线后由efast_sd.m和efast_ttest.m给出稳定估计。M 值不宜超过 6，否则高频混叠会把S_i抬得很虚。

3. ODE_efast.m 的接入点：把你的微分方程装进全局敏感性分析的入口

3.1 优先修改Model_efast.m，主程序不要动

ODE_efast.m是调度器，它只认一个接口：输入一行1×K参数向量，输出一个标量结果。要接你自己的模型，改Model_efast.m即可。以 SIR 传染模型为例：

% Model_efast.m —— 用户自定义 ODE 模型 function out = Model_efast(p) % p 是 1×3 参数行向量，顺序由 Parameter_settings_EFAST.m 决定 beta = p(1); % 传染率 gamma = p(2); % 恢复率 I0 = p(3); % 初始感染比例 tspan = [0 100]; y0 = [1 - I0; I0; 0]; % S, I, R 初值 [~, Y] = ode45(@(t,y) sir_rhs(t,y,beta,gamma), tspan, y0); out = Y(end, 2); % 取第 100 时刻的感染比例作为敏感性分析指标 end function dydt = sir_rhs(t, y, beta, gamma) dydt = [-beta * y(1) * y(2); beta * y(1) * y(2) - gamma * y(2); gamma * y(2)]; end

这段代码里p是唯一的输入，内部把向量里的元素解析成模型参数，最后只返回一个标量。eFAST 的全部统计学计算依赖这个标量，所以不要在这里输出整个状态向量。你要评估峰值时刻、累计病例数还是末端状态，都先算成数再返回。

3.2ODE_efast.m内部的循环骨架

主程序的循环结构并不复杂，每个参数轮流成为“被扫描参数”，其他参数在同一曲线上同步波动：

% ODE_efast.m 内部的核心循环示意 for rep = 1:R % R 次随机相位重采样 phase = rand(1, K) * 2 * pi; for i = 1:K % 生成第 i 条搜索曲线，矩阵 X 的每行是一组参数 X = gen_search_curve(freq, phase, N, K); for n = 1:N y(n) = Model_efast(X(n, :)); % 每次调用都求解一次 ODE end spec = abs(fft(y - mean(y))).^2; % 只取单侧谱，索引要 +1 抵消直流 Si(rep, i) = sum(spec(freq(i) * (1:M) + 1)) / sum(spec(2:end)); end end

这段代码的关键是gen_search_curve生成的X里，第i列在整条曲线上使用freq(i)扫描，其余列的频率不同，但共用同一个phase。因此输出序列y中每个参数的影响都被调制到不同频带上，FFT 之后才能按频率分离。每次n循环都是一次完整的 ODE 积分，N一般在两三百到一千之间，重采样次数R决定了总求解次数。

3.3 多个输出指标与时间序列的处理方式

如果你的 ODE 输出是一整段时间序列，不要直接把向量交给 eFAST。敏感性指数定义在单个标量输出上，常见做法是取几个关键指标分别跑一轮分析，比如峰值、达到某个阈值的时间、末端值、曲线下面积。参数设置文件里可以加一个out_select字段，在Model_efast.m中根据它切换返回指标。每多一个输出指标，计算量都会成倍增加。

4. 全局敏感性分析完整运行：从参数设置到结果解释

4.1 在Parameter_settings_EFAST.m里定义参数空间

运行的第一件事是打开Parameter_settings_EFAST.m，把参数数量、取值范围、分布类型和采样规模写清楚。名字可能随代码版本略有差异，但核心字段是这几个：

% Parameter_settings_EFAST.m 常用配置项 K = 3; % 参数个数，必须和 Model_efast.m 里 p 的长度一致 N = 377; % 每条搜索曲线的采样点数，参照 2M*max(freq)+1 M = 4; % 谐波阶数，默认 4 R = 5; % 独立重采样次数，用于计算标准差和 t 检验 freq = [11 21 37 47 61 79]; % 互质频率表，实际只取前 K 个 range = [0.01 0.5; 0.05 0.8; 0.1 0.9]; % 每个参数的上下界 dist = {'uniform', 'uniform', 'uniform'}; % 每个参数的分布类型

N不是越大越好。它直接决定 ODE 求解次数，K=5、N=489、R=5时一共要跑5×489×5 = 12225次ode45，如果模型本身要好几秒算完，总时间会很难看。我通常先用N=129或N=257快速跑一轮，确认量级后再加大N和R做正式分析。

4.2 标准命令序列

按顺序在 MATLAB 命令行执行：

Parameter_settings_EFAST; % 把上一步配置写入工作区 ODE_efast; % 生成搜索曲线并求解 ODE，得到敏感性指数 CVmethod; % 交叉验证，检查指数是否随样本数收敛 efast_sd; % 计算多组重采样下 Si 和 STi 的标准差 efast_ttest; % 对重采样结果做 t 检验，标记不显著参数

CVmethod.m的作用是计算敏感性指标的变异系数，如果CV > 0.1，说明当前N或R不足以得到稳定结果。efast_ttest输出的是像是“这个参数是否得到了显著非零的敏感度”的结论，而不是参数本身是否显著；它帮你把结果从“数值排序”变成“统计判断”。

4.3 输出结果怎么看

运行完成后，工作区里常见的输出数组是S、ST和它们各自的标准差。下面是一张典型的汇总结果：

当STi明显大于Si，代表该参数主要通过与其他参数耦合影响输出，局部敏感性分析看不到这种关系。最后一个参数虽然主效应只有 0.05，但总效应达到 0.47，所以也不能简单固定，否则模型在参数组合变化时的行为会被低估。

5. 用敏感性排序做降维：固定参数与再拟合的技巧

5.1 把低敏感参数固定为中位数，重新标定少参数模型

全局敏感性分析的直接收益不是一张排名表，而是帮你把高维 ODE 模型降到可标定的低维模型。常见做法是：先用一轮 eFAST 拿到STi，把总效应低于阈值（比如 0.05）的参数固定为当前范围中位数，只对剩余高敏感参数做拟合：

% reduce_model.m —— 根据 STi 固定低敏感参数 threshold = 0.05; low_idx = find(STi < threshold); % 找出可固定参数索引 fixed_value = pmedian(low_idx); % 取当前范围中位数 % 重新组装 p_fixed，只留下高敏感参数作为估计目标 p_est = p_high_only(fixed_value);

固定参数之前要确认它不会导致 ODE 刚性变化。有些参数单独看对末端输出不敏感，但会影响数值稳定性，比如积分步长下限或反应速率常数。如果固定后ode45的计算时间变长或结果发散，就把该参数重新放回拟合集合。

5.2 用efast_ttest结果区分“不敏感”和“没算出来”

efast_ttest.m的输出能帮你减少一种常见错误：某些参数的真实效应很小，但由于N太小，FFT 频谱里出现假峰，导致排序虚高。我一般把p<0.01且STi排名前 40% 的参数作为待估计参数；如果 t 检验不显著，即使STi数值不小，也先当作采样波动处理，适当增加R再跑一轮。

另一个容易踩的坑：Si和STi都是相对方差的占比，它们之和可以大于 1，因为交互效应会在多个参数的总效应里重复计算。所以不要拿两个加和去校正模型误差，直接按STi排序即可。最终完成一轮筛选后，把高敏感参数单独提出来重跑 ODE 参数估计，会明显快于在全部参数上盲目搜索，而且参数间相关性带来的局部极小值也会少很多。

本文还有配套的精品资源，点击获取