优化仿真或常微分方程

什么是优化仿真或 ODE？

有时，目标函数或非线性约束函数值只能通过仿真或常微分方程 (ODE) 的数值解获得。这种优化问题有几个共同的特征和挑战，具体请参阅可能的问题和解决方法。

有关优化 ODE 的基于问题的示例，请参阅使用优化变量拟合 ODE 参数。有关基于求解器的示例，请参阅拟合常微分方程 (ODE)。

有关避免其他方法遇到的许多问题的方法，请参阅基于问题的离散化最优轨迹。该方法可以在优化过程中使用自动微分。然而，该方法的精度可能相对较低，因为它基于定步长和可能的低阶 ODE 求解算法。

要轻松优化 Simulink^® 模型，请尝试使用 Simulink Design Optimization™。

有限差分中的问题

Optimization Toolbox™ 求解器在内部使用目标函数和约束函数的导数。默认情况下，它们使用以下形式的有限差分逼近来估计这些导数

$\frac{F (x + δ) - F (x)}{δ}$

或

$\frac{F (x + δ) - F (x - δ)}{2 δ} .$

这些有限差分逼近可能不准确，因为：

大的 δ 值会对有限差分产生更大的非线性影响。
小的 δ 值会因数值精度有限而导致误差。

具体来说，对于 ODE 的仿真和数值解：

仿真通常对参数的微小变化不敏感。这意味着，如果您使用过小的扰动 δ，仿真可能返回虚假的估计导数 0。
ODE 的仿真和数值解在其函数计算中可能都有误差。这些误差在有限差分逼近中可能会被放大。
ODE 的数值解引入比计算机精度大得多的噪声值。这种噪声在有限差分逼近中可能会被放大。
如果 ODE 求解器使用变步长，则有时 F(x + δ) 计算中的 ODE 步数可能不同于 F(x) 计算中的步数。这种差别可能导致返回值中出现虚假差异，从而给出误导性的导数估计。

有限差分的建议

使用直接搜索避免有限差分
设置更大的有限差分
使用梯度计算函数
使用更严格的 ODE 容差

使用直接搜索避免有限差分. 如果您有 Global Optimization Toolbox 许可证，您可以尝试使用 patternsearch (Global Optimization Toolbox) 作为求解器。patternsearch 不尝试估计梯度，因此不会受到有限差分中的问题的限制。

如果您使用 patternsearch 进行昂贵（耗时）的函数计算，请使用 Cache 选项：

options = optimoptions('patternsearch','Cache','on');

如果您无法使用 patternsearch，并且有相对低维的无约束最小化问题，请改为尝试 fminsearch。fminsearch 不使用有限差分。但是，fminsearch 不是快速或可调的求解器。

设置更大的有限差分. 有时，您可以采取大于默认值的有限差分步长来避免有限差分中的问题中的问题。

如果您有 MATLAB^® R2011b 或更高版本，请将有限差分步长选项设置为大于默认的 sqrt(eps) 或 eps^(1/3) 的值，例如：
- 对于 R2011b–R2012b：
  options = optimset('FinDiffRelStep',1e-3);
- 对于 R2013a–R2015b 和名为 'solvername' 的求解器：
  options = optimoptions('solvername','FinDiffRelStep',1e-3);
- 对于 R2016a 及更高版本和名为 'solvername' 的求解器：
  options = optimoptions('solvername','FiniteDifferenceStepSize',1e-3);
如果不同分量中具有不同尺度，请将有限差分步长设置为与分量尺度成比例的向量。
如果您有 MATLAB R2011a 或更早版本，请将 DiffMinChange 选项设置为比默认的 1e-8 更大的值，也可以设置 DiffMaxChange 选项，例如：
```
options = optimset('DiffMinChange',1e-3,'DiffMaxChange',1);
```

注意

很难知道如何设置这些有限差分大小。

您也可以尝试设置中心有限差分：

options = optimoptions('solvername','FiniteDifferenceType','central');

使用梯度计算函数. 为了避免有限差分估计的问题，您可以在目标和非线性约束中提供逼近梯度函数。请记住使用 optimoptions 将 SpecifyObjectiveGradient 选项设置为 true；如果相关，还要将 SpecifyConstraintGradient 选项设置为 true。

对于某些 ODE，您可以在求解 ODE 的同时以数值方式计算梯度。例如，假设目标函数 z(t,x) 的微分方程为
$\frac{d}{d t} z (t, x) = G (z, t, x),$
其中，x 是您要最小化的参数向量。假设 x 是标量。则其导数 y 的微分方程，
$y (t, x) = \frac{d}{d x} z (t, x)$
为
$\frac{d}{d t} y (t, x) = \frac{\partial G (z, t, x)}{\partial z} y (t, x) + \frac{\partial G (z, t, x)}{\partial x},$
其中，z(t,x) 是目标函数 ODE 的解。您可以在与 z(t,x) 相同的微分方程组中求解 y(t,x)。此解为您给出逼近的导数，而不使用有限差分。对于非标量 x，对每个分量求解一个 ODE。
有关此方法的理论和计算介绍，请参阅 Leis 和 Kramer 合著的·[2]。有关此方法和有限差分方法的计算经验，请参阅 Raue et al. 所著 [3] 中的图 7。
对于某些仿真，您可以在仿真内估计导数。例如，在 Reiman 和 Weiss 合著的 [4] 中所述的似然比方法或在 Heidelberger、Cao、Zazanis 和 Suri 合著的 [1] 中分析的无穷小扰动分析方法，就是在估计目标或约束函数的同一仿真中估计导数。

使用更严格的 ODE 容差. 您可以使用 odeset 将 AbsTol 或 RelTol ODE 求解器容差设置为低于默认值的值。不过，选择太小的容差可能会导致求解缓慢、无法收敛或其他问题。对于 RelTol，切勿选择小于 1e-9 的容差。AbsTol 的每个分量的下限取决于您的问题的规模，此处无法提供建议。

随机函数中的问题

如果仿真使用随机数，则对目标函数或约束函数进行两次求值会返回不同结果。这会影响函数估计和有限差分估计。有限差分的值可能主要源于随机性带来的变化，而不是源于不同计算点 x 和 x + δ 带来的变化。

对随机函数的建议

如果您的仿真使用您控制的流中的随机数，请在每次计算目标或约束函数之前重置随机流。这种做法可以减少结果的变异性。例如，在以下目标函数中：

function f = mysimulation(x)
rng default % or any other resetting method
...
end

有关详细信息，请参阅生成可重复的随机数。

目标和约束的常用计算方法

通常，在同一仿真运行期间，仿真会同时计算目标函数和约束。或者，目标和非线性约束函数都使用相同的成本高昂的计算。fmincon 等求解器会分别计算目标函数和非线性约束函数。这可能导致极大的效率损失，因为求解器会调用两次成本高昂的计算。要避免此问题，请使用同一函数中的目标和非线性约束中的方法，或当使用基于问题的方法时使用具有串行或并行共同函数的目标和约束，基于问题中的方法。

目标或约束函数计算失败

对于某些参数值，您的仿真或 ODE 可能会失败。

针对计算失败的建议

设置适当的边界. 虽然您可能不知道对参数空间的所有限制，但请尝试对所有参数设置适当的上界和下界。这可以加快优化速度，并有助于求解器避免出现有问题的参数值。

使用遵守边界的求解器. 如迭代可能违反约束中所述，一些算法可能在中间迭代中违反边界约束。为了优化仿真和 ODE，请使用始终遵守边界约束的算法。请参阅满足边界约束的算法。

返回 NaN. 如果您的仿真或 ODE 求解器无法在 x 点处成功计算目标或非线性约束函数，请让您的函数返回 NaN。大多数 Optimization Toolbox 和 Global Optimization Toolbox 求解器在遇到 NaN 值时，都具有尝试不同迭代步的稳健性。这些稳健的求解器包括：

fmincon interior-point、sqp 和 trust-region-reflective 算法：
fminunc
lsqcurvefit
lsqnonlin
patternsearch

有些人倾向于在不成功、不可行或其他不理想的点上返回任意大的目标函数值。但是，这种做法可能会使求解器困惑，因为求解器没有意识到返回的值是任意的。当您返回 NaN 时，求解器可能会尝试在不同的点进行计算。

参考书目

[1] Heidelberger, P., X.-R. Cao, M. A. Zazanis, and R. Suri. Convergence properties of infinitesimal perturbation analysis estimates. Management Science 34, No. 11, pp. 1281–1302, 1988.

[2] Leis, J. R. and Kramer, M.A. The Simultaneous Solution and Sensitivity Analysis of Systems Described by Ordinary Differential Equations. ACM Trans. Mathematical Software, Vol. 14, No. 1, pp. 45–60, 1988.

[3] Raue, A. et al. Lessons Learned from Quantitative Dynamical Modeling in Systems Biology. Available at https://journals.plos.org/plosone/article?id=10.1371/journal.pone.0074335, 2013.

[4] Reiman, M. I. and A. Weiss. Sensitivity analysis via likelihood ratios. Proc. 18th Winter Simulation Conference, ACM, New York, pp. 285–289, 1986.