Epipolar a personal journal

最小二乘问题(六)

前面提到的最小二乘有一个共同的假设,就是残差向量遵循高斯分布。

那么,如果不满足高斯分布假设呢?

通常我们模型的噪音都是符合高斯分布的,对于模型中蕴含的非高斯分布情形,我们不在这里讨论了。我们介绍一种在高斯假设下依旧会引入非高斯分布残差的情形,并介绍在最小二乘中处理它的一种手段。

通常我们在建模时认为数据点都是由模型产生并带有一定的噪音。但在真实中,数据可能掺杂着不属于当前模型的点,叫做外点(outlier)。这些外点会污染数据,表现上便是出现了不符合高斯分布的噪音。由于我们采用平方项,当与模型偏差很远的点存在时,它会对模型的估计带来严重的不良影响。

为了减小这种影响,我们可以采用加权的最小二乘,即

理想情况下,只要令外点对应的权值为 0 即可从优化中忽略掉它们。

然而真实情况中,我们可能不知道外点都有哪些,因此我们需要一个可以自动降低外点影响的加权方案。”M-估计”便是一种常用的方法,它在优化中引入了“鲁棒函数”来根据误差的大小自动调整权重:

这里的 $\rho_i$ 便是鲁棒函数,如果 $\rho_i(x) = x^2$ ,我们就得到了一般的最小二乘问题。为了减小外点(带有巨大误差的点)的影响,一般要取 $\lim_{x\to\infty}\frac{\rho(x)}{x^2} \to 0$ 的函数。也就是当误差越大,它对整体优化的影响就越弱于普通的平方项。但同时,$\rho(x)$ 还应该满足从 0 向两个方向增长,否则优化可能会收敛到错误的结果(比如直接拒绝掉全部的数据点从而得到更小的能量)。

常见的鲁棒函数选择很多,比如 Huber 函数、Cauchy 函数等等。网上资料很丰富,这里就不列举了。

加入了鲁棒函数后,我们的问题就不是最小二乘了。但前面非线性最小二乘中,我们通过线性化得到了线性最小二乘子问题。这里我们也尝试一下转化它成为我们熟悉的形式。

我们比较带鲁棒函数的能量和一般的最小二乘能量的导数,根据极值的一阶条件我们知道:

如果我们记 $w(x) = \frac{d\rho}{dx}\frac{1}{x}$,进一步地令常系数 $w’= w(\|r\|)$ ,然后构造加权的最小二乘优化 $\min \frac{w’}{2}\|r\|^2$,它的 Jacobian 是:

可以发现,它与前面带有鲁棒函数的最小二乘优化形式是相同的。

也就是说,我们可以在每一轮迭代时,采用当前的 $\|r\|$ 计算加权值。固定加权值后,将当前迭代中的优化问题看做一个加权最小二乘问题进行优化。这时,带鲁棒函数的 M-估计与带更新的加权最小二乘是等价的。

Cholesky 分解(番外二)

为了模拟 Cholesky 分解的 fill-in 现象,写了一个矩阵消元小游戏。

它是 Web 界面,地址是 http://epipolar.tumblr.com/matrix-elimination 。这也是我第一个用 javascript 和 processing (p5.js) 编写的程序。

欢迎尝试~

最小二乘问题(五)

上回末有人搞事儿,希望加入非高斯的噪声……但我们偏不!我们要先看一些更重要的扩展。

我们记 $r(x) = Ax-b$ ,这是我们的(线性)残差函数。那么标准的最小二乘问题就是 $\min \|r(x)\|^2$ 。残差函数关于变量 $x$ 的导数(Jacobian)记作 $J_r$ ,我们知道 $J_r = A$,于是最小二乘的标准方程还可以写成 $J_r^TJ_rx = J_r^Tb$ 。

如果 $r(x)$ 并非关于 $x$ 的线性函数呢?此时我们得到了年轻人的第一个非线性最小二乘问题!

它与线性最小二乘最大的区别便在于:非线性最小二乘通常不是凸优化问题。这是什么意思?恩我们不在这里解释了。那这个区别有什么用?恩它没什么用,它是来给我们制造新的困难的,这个困难就是:对于线性最小二乘问题,你可以求解标准方程得到全局最优解;然而对于非线性最小二乘问题……

  • 它有标准方程么?是什么形式?
  • 它的解如何寻找?可以求解标准方程么?
  • 它有全局最优解么?

恩……对于一般的非线性最小二乘问题,通常我们只能找到一个局部最优解,至于局部最优是不是全局最优,这就要靠手段加运气了。它的解要如何寻找呢?答案是从一个初始的猜测出发,一步步走到更优的解,直到收敛——这也叫做迭代优化方法。至于标准方程,恩我们下面来看看非线性最小二乘的标准方程。

如果残差函数有具体的表达式,那当然是可以更好地分析,不过这里我们从它的一般形式入手。一个问题……它不是线性的,怎么办?线性化咯!

于是我们将残差函数 $r(x)$ 用它的一阶近似表达,也就是

这样,当我们从一个初始值 $x_0$ 出发时,我们可以寻找一个让一阶近似更优的变化量 $\Delta x$,这样我们便“大约”将原始问题优化了一些。

这个思路写成最优化形式,就是

它是个关于 $\Delta x$ 的线性最小二乘问题,套用我们前面的知识,我们知道它的标准方程是

前面我们已经提到了,整个问题需要一步步不停地迭代,所以说这里我们只是完成了一个小目标,接下来我们要更新变量的值并进一步优化。

所以完整的非线性最小二乘优化算法是这样的:

  1. $x \gets x_0;$
  2. $\Delta x \gets \arg\min_{\Delta x} |J_r\Delta x+r(x)|^2;$
    • 求解标准方程 $J_r^TJ_r\Delta x = -J_r^Tr(x);$
  3. $x \gets x+\gamma \Delta x;$
  4. 回到 2 ,除非 $|\Delta x|<\epsilon$ 。

这个算法叫做 Gauss-Newton 算法。而其中关于 $\Delta x$ 的标准方程也叫作 Gauss-Newton 标准方程。$\gamma > 0 $ 是一个控制步长的系数。

那么如果有多项呢?如果残差遵循不同的高斯分布呢?线性化,然后回忆前面的故事吧!

Random Thoughts

写程序无非是做两件事:解一个数学问题和把数据转换成问题要求的输入格式。

在课堂上,我们学习的是解问题。然而到了工业界,更多的人做的是转格式。

Photo

大坂 クリスタルタワー

仿佛大楼带了自发光……