1.5.5 优化及拟合：scipy.optimize

优化是寻找最小化或等式的数值解的问题。

scipy.optimize 模块提供了函数最小化（标量或多维度）、曲线拟合和求根的有用算法。

In [19]:

  1. from scipy import optimize

寻找标量函数的最小值

让我们定义下面的函数：

In [20]:

  1. def f(x):
  2.     return x**2 + 10*np.sin(x)

绘制它：

In [21]:

  1. x = np.arange(-10, 10, 0.1)
  2. plt.plot(x, f(x))
  3. plt.show()

1.5.5 优化及拟合：scipy.optimize - 图1

这个函数在-1.3附近有一个全局最小并且在3.8有一个局部最小。

找到这个函数的最小值的常用有效方式是从给定的初始点开始进行一个梯度下降。BFGS算法是这样做的较好方式：

In [22]:

  1. optimize.fmin_bfgs(f, 0)
  1. Optimization terminated successfully.
  2.          Current function value: -7.945823
  3.          Iterations: 5
  4.          Function evaluations: 24
  5.          Gradient evaluations: 8

Out[22]:

  1. array([-1.30644003])

这个方法的一个可能问题是，如果这个函数有一些局部最低点，算法可能找到这些局部最低点而不是全局最低点，这取决于初始点：

In [23]:

  1. optimize.fmin_bfgs(f, 3, disp=0)

Out[23]:

  1. array([ 3.83746663])

如果我们不知道全局最低点，并且使用其临近点来作为初始点，那么我们需要付出昂贵的代价来获得全局最优。要找到全局最优点，最简单的算法是暴力算法，算法中会评估给定网格内的每一个点：

In [24]:

  1. grid = (-10, 10, 0.1)
  2. xmin_global = optimize.brute(f, (grid,))
  3. xmin_global

Out[24]:

  1. array([-1.30641113])

对于更大的网格，scipy.optimize.brute() 变得非常慢。scipy.optimize.anneal() 提供了一个替代的算法，使用模拟退火。对于不同类型的全局优化问题存在更多的高效算法，但是这超出了scipy的范畴。OpenOpt、IPOPT、PyGMO和PyEvolve是关于全局优化的一些有用的包。

要找出局部最低点，让我们用scipy.optimize.fminbound将变量限制在(0,10)区间：

In [25]:

  1. xmin_local = optimize.fminbound(f, 0, 10)    
  2. xmin_local

Out[25]:

  1. 3.8374671194983834

注：寻找函数的最优解将在高级章节中：数学优化：寻找函数的最优解详细讨论。

寻找标量函数的根

要寻找上面函数f的根，比如f(x)=0的一个点，我们可以用比如scipy.optimize.fsolve()：

In [26]:

  1. root = optimize.fsolve(f, 1)  # 我们的最初猜想是1
  2. root

Out[26]:

  1. array([ 0.])

注意只找到一个根。检查f的图发现在-2.5左右还有应该有第二个根。通过调整我们最初的猜想，我们可以发现正确的值：

In [27]:

  1. root2 = optimize.fsolve(f, -2.5)
  2. root2

Out[27]:

  1. array([-2.47948183])

曲线拟合

假设我们有来自f的样例数据，带有一些噪音：

In [28]:

  1. xdata = np.linspace(-10, 10, num=20)
  2. ydata = f(xdata) + np.random.randn(xdata.size)

现在，如果我们知道这些sample数据来自的函数（这个案例中是$x^2 + sin(x)$）的函数形式，而不知道每个数据项的系数，那么我们可以用最小二乘曲线拟合在找到这些系数。首先，我们需要定义函数来拟合：

In [29]:

  1. def f2(x, a, b):
  2.     return a*x**2 + b*np.sin(x)

然后我们可以使用scipy.optimize.curve_fit()来找到a和b：

In [30]:

  1. guess = [2, 2]
  2. params, params_covariance = optimize.curve_fit(f2, xdata, ydata, guess)
  3. params

Out[30]:

  1. array([  0.99719019,  10.27381534])

现在我们找到了f的最优解和根，并且用曲线去拟合它，我们将这些结果整合在一个图中：

1.5.5 优化及拟合：scipy.optimize - 图2

注：在Scipy >= 0.11中，包含所有最小值和寻找根的算法的统一接口：scipy.optimize.minimize()、 scipy.optimize.minimize_scalar()和 scipy.optimize.root()。他们允许通过method关键词容易的比较多种算法。

你可以在scipy.optimize中找到对于多维度问题有相同功能的算法。

练习：温度数据的曲线拟合

下面是从1月开始阿拉斯加每个月的温度极值（摄氏度）：

最大值: 17, 19, 21, 28, 33, 38, 37, 37, 31, 23, 19, 18

最小值: -62, -59, -56, -46, -32, -18, -9, -13, -25, -46, -52, -58

  • 绘制这些温度极值。
  • 定义一个函数，可以描述温度的最大值和最小值。提示：这个函数的周期是一年。提示：包含时间偏移。
  • 用scipy.optimize.curve_fit()拟合这个函数与数据。
  • 绘制结果。这个拟合合理吗？如果不合理，为什么？
  • 最低温度和最高温度的时间偏移是否与拟合一样精确？练习：2-D 最小值

1.5.5 优化及拟合：scipy.optimize - 图3

六峰驼背函数： 1.5.5 优化及拟合：scipy.optimize - 图4

有多个全局和局部最低点。找到这个函数的全局最低点。

提示：

  • 变量可以被限定在-2 < x < 2 和 -1 < y < 1。
  • 用numpy.meshgrid() 和 pylab.imshow() 来从视觉上来寻找区域。
  • scipy.optimize.fmin_bfgs() 或者另一个多维最小化。 多几个全局最小值，那些点上的函数值十多少？如果最初的猜测是$(x, y) = (0, 0)$会怎样?

看一下非线性最小二乘曲线拟合：地形机载激光雷达数据中的点抽取练习的总结，以及更高及的例子。