人工神经网络方法解方程

OctoMiao bio photo By OctoMiao Comment

数值方法解方程大多要一个把方程或者函数参数化,例如改写成差分方程,变成步长等等的函数。例如我们要解这样一个方程:

其中初始条件 $y(0)=1$ 已知。

一种普遍的手段是,我们可以把其中的待解函数参数化为 $y(t,\{n_i\})$,然后把方程重新写成

其中 $\{n_i\}$ 是用来参数化这个函数的参数。有了这个式子,我们会联想到最小二乘法。也就是说,如果我们的参数化完全正确的话,我们应该得到方程的左边的两项之和总是零,也就是 $\frac{dy(t, \{n_i\})}{dt}$ 和 $- y(t,\{n_i\}) $ 的没有偏离,即

这个量总是零。然而数值上来讲,我们有一个更加简单的方法来让这个量在任意的 $t$ 总是零,方法就是借助最小二乘法的思路,利用一个凹函数,

倘若 $\text{cost} = 0$,我们推断上面定义的 $\text{waste}$ 函数在任意的 $t$ 都为零。

人工神经网络参数化

下面我们引入一种与人工神经网络相关的参数化方法,然后通过定义一个代价,最后最小化这个代价从而获得参数化中的各种参数,从而解出函数。

人的神经元有一个特点就是,存在激活和不激活两种状态,激活的话信号传递,不激活的话,信号被毙掉。数学上我们可以通过一个类似阶梯的函数来模拟这个行为,例如:

阶梯函数图像,可以挑选出某些参数值。

然而数值上来讲,这个函数有个非常大的缺点,就是导数不连续。尤其是对于我们想要解一个微分方程来说,用来参数化的函数里面包含导数不连续的部分,那是非常不方便的。

所以我们会选取一个导数连续但是也具有类似行为的函数,比如一个 sigmoid 函数,

这个函数的图像是

sigmoid 函数,也具有挑选出某些参数值的功能。

这个函数有个非常大的优势,就是如果我们把要求解的函数利用这个参数化,所有的导数全是可以解析求解的,对于数值计算是一个福音。

我们可以模拟神经元可以被激活这种性质,来挑选出一些我们想要的参数,这个示例中用到的就是 sigmoid 函数。我们挑选一组参数 ${v_k}$, ${w_k}$ 和 ${u_k}$,对于一个函数 $y(t)$,我们可以用下面的方法来参数化,

这样做的原因是在 $t=0$ 时,我们确保参数化的函数满足给定的初始条件($t=0$ 时 $t \sum_k v_k f(t w_k+u_k)=0$)。而 $f(x)$ 使我们的启动函数,这里我们用 sigmoid 函数。也就是说,在某个时刻 $t_i$,当参量 $t w_k+u_k$ 比较大的时候,我们输入的 $t$ 才会起作用,否则我们的输入 $t$ 被压制。

$w_k$ 是对启动函数的缩放,$u_k$ 是对启动函数的平移,而 $v_k$ 的作用是放大缩小函数值。可以设想,当我们的参数足够多的时候,我们可以用这个参数化来组成任何连续函数。

cost

那么如何得到这些参数的值呢?我们一开始定义了这样一个量,

我们想要做的是得到什么的参数可以有 $\text{cost} = 0$ 的结果。

然而我们不想做积分,所以自然的选择是离散化,选择一个序列 ${t_i}$,并且考虑

当我们选到一组参数使得这个量为零的时候,我们的参数就是函数的正确的参数化(在当前的初始条件下的)。

所以我们就把一个解方程问题转换成了一个最小化问题了,因为 $\text{cost}$ 的值越小,说明我们的参数就越接近满足

当然是在一个特定初始条件下,这里我们的初始条件是 $y(0)=1$。

剩下的部分,就无需我多说了,只需要用任何你想要的方法,来把 $\text{cost}$ 最小化,得到的参数值们 $\{v_k\}$, $\{w_k\}$ 和 $\{u_k\}$ 带回到函数的参数化

这样我们就有了这个微分方程在当前初始条件下的解。

上面我们离散化 cost 的过程中,

里面的每一个时间点的输入都是一次对神经网络的训练,如同人的学习一样,经过多次训练,我们就可以掌握方程的这种行为。而所学习到的信息,存在了神经网络的参数中(也就是网络的结构中)。

代码举例

以下是这个问题的 Python 代码,这里有一份带有全面说明的 IPython Notebook

import numpy as np
from scipy.optimize import minimize
from scipy.special import expit
import matplotlib.pyplot as plt

def cost(v,w,u,t):
v = np.array(v) # Don't know why but np.asarray(v) doesn't work here.
w = np.array(w)
u = np.array(u)

fvec = np.array(trigf(t*w + u) ) # This is a vector!!!
yt = 1 + np.sum ( t * v * fvec ) # For a given t, this calculates the value of y(t), given the parameters, v, w, u.
return ( np.sum (v*fvec + t * v* fvec * ( 1 - fvec ) * w ) + yt ) ** 2

def trigf(x):
#return 1/(1+np.exp(-x)) #
return expit(x)

def costTotal(v,w,u,t):
t = np.array(t)
costt = 0
for temp in t:
costt = costt + cost(v,w,u,temp)
return costt

costTotalF = lambda x: costTotal(np.split(x,3)[0],np.split(x,3)[1],np.split(x,3)[2],tlin)

initGuess = np.zeros(30)
result = minimize(costTotalF,initGuess,method="SLSQP")

为什么找这么多麻烦

显然,上面的例子不仅可以用任何差分方法来解,而且有解析解。我们为什么要找这么多麻烦来做这样的参数化呢?

或许我们要解决一个当前的状态依赖于全空间中的所有点的问题,这时候,我们就不好再用差分法了,因为要解出当前点我们需要知道所有的其他的点。上面这种一次性解出整个方程的方法就变得更加方便,因为我们解出方程并不是依赖于获得其他地方的信息,而是通过人工神经网络一次性猜到所有空间点的解。

不过,本文开始讲的这种方法具有非常好的普适性。我们使用 ANN 是因为 Cybenko 在 1989 年证明了一个 sigmoid 可以作为很好的 universal approximator 来处理任意的 measurable functions。1

现在我们重新思考一下这个方法的普适性。倘若我们有个系统,可以用傅里叶展开,而且所幸我们只需要展开中的前 100 项就可以很好的来估算整个系统了。那么我们就把这个展开式写出来,然后取前 100 项,这样一样定义一个 cost,一样来找到使得 cost 最小化的参数,就可以求出傅里叶展开的系数们了。

不过这个方法是不是足够有效,取决于我们所挑选的参数化形式。Kolmogorov 证明过如果选的足够好,我们可以用有限个与 y 无关的函数来精确的还原 y。

comments powered by Disqus