【Python讲优化】S06E05 详解多元函数的梯度应用

1.多元函数偏导数的数值解

在程序当中,利用数值方法求出各个自变量偏导数的近似解,其方法和步骤同前面讲过的导数的数值解求法并无二致:把其余的自变量固定,就将偏导数的求解方法等价为了导数的数值求解方法,我们以简单的二元函数$f(x,y)=x^2-y^2$为例,分别来看看如何利用$python$求解偏导数$f_x$和$f_y$,并实际获取点$(-1,-1)$处的数值解:

代码片段:

def f(x,y):
    return x**2-y**2

def grad_x(f, x, y):
    h = 1e-4
    return (f(x + h/2, y) - f(x - h/2, y)) / h

def grad_y(f, x, y):
    h = 1e-4
    return (f(x, y + h/2) - f(x, y - h/2)) / h

print(grad_x(f, -1, -1))
print(grad_y(f, -1, -1))

运行结果:

-2.000000000002
2.000000000002

因此,我们非常轻松的获取了函数$f(x,y)=x^2-y^2$,在点$(-1,-1)$处的偏导数$f_x$和$f_x$。

2.通过梯度场直观感受梯度

在前一讲的内容中,我们知道了,二元函数中的梯度概念可以说是一元函数中导数概念的一个延伸和拓展。

他们之间的区别是一元函数的导数$f'(x)$是一个数,而梯度$\nabla f(p)$则是一个向量,对于二元函数$f(x,y)$,梯度对应的向量就是$(f_x,f_y)$,向量的概念大家都很熟悉,他是一个有方向、有大小的量。

那么很显然的,对于一个二元函数$f(x,y)$,在他的定义域内,如果我们把每一个点的梯度都求出来,将每个点的梯度向量和各个点的位置联系起来进行集中展示,就形成了一个梯度场,场的概念是一个比较新的概念,我们还是针对二元函数$f(x,y)=x^2-y^2$来实际展示一下:

代码片段:

import numpy as np
import matplotlib.pyplot as plt

def f(x, y):
    return x**2-y**2

def grad_x(f, x, y):
    h = 1e-4
    return (f(x + h/2, y) - f(x - h/2, y)) / h

def grad_y(f, x, y):
    h = 1e-4
    return (f(x, y + h/2) - f(x, y - h/2)) / h

def numerical_gradient(f,P):
    grad = np.zeros_like(P)
    for i in range(P[0].size):
        grad[0][i] = grad_x(f, P[0][i], P[1][i])
        grad[1][i] = grad_y(f, P[0][i], P[1][i])
    return grad

x = np.arange(-2, 2, 0.25)
y = np.arange(-2, 2, 0.25)

X, Y = np.meshgrid(x, y)
X = X.flatten()
Y = Y.flatten()

grad = numerical_gradient(f, np.array([X, Y]))

plt.quiver(X, Y, grad[0], grad[1])#grad[0]是一个1*X.size的数组
plt.xlim([-2, 2])
plt.ylim([-2, 2])
plt.xlabel('x')
plt.ylabel('y')
plt.grid()
plt.show()

运行结果:
图1.梯度场示意图

图1.梯度场示意图

上面这段代码我们简要分析一下:

在$x$坐标轴上,从$-2$到$2$,按照$0.25$的间隔,我们一共生成了$16$个点,$y$坐标也是同理,一共也是有$16$个点,这样就对应表示出了平面上的$256$个点。

通过:

X, Y = np.meshgrid(x, y)

进行网格化,得到的$X$是一个二维数组,$16$行$16$列,对应表示平面上的$256$个点的横坐标,同样$Y$也是一个二维数组,$16$行$16$列,依次的对应表示上面那$256$个点的纵坐标。

X = X.flatten()
Y = Y.flatten()

上面的这两行代码将这两个二维数组展平成一维数组,$X$和$Y$都变成含有$256$个元素的一维数组,其中$X[i]$和$Y[i]$分别对应表示这$256$个点中第$i$个点的横纵坐标。

def numerical_gradient(f,P):
    grad = np.zeros_like(P)
    for i in range(P[0].size):
        grad[0][i] = grad_x(f, P[0][i], P[1][i])
        grad[1][i] = grad_y(f, P[0][i], P[1][i])
    return grad
top Created with Sketch.