数学狂想曲（三）——统计杂谈, PID算法, 20世纪10大算法, 矩阵&向量的积

http://antkillerfarm.github.io/

统计杂谈

统计模拟

统计模拟是数理统计中非常有用的工具之一， 它是利用计算机产生某概率模型的随机数，再通过这些随机数来模拟真实模型。

这方面的书籍首推Sheldon M.Ross所著的《Simulation》。

注：Sheldon M.Ross，斯坦福大学统计学博士（1968），UCB教授（1976～2004），南加州大学Industrial and Systems Engineering系主任。

以下仅记录一些简单的结论。

均匀分布

通过线性同余发生器可以生成伪随机数，且该随机数满足均匀分布。

xn=axn−1modulom

正态分布

Box-Muller变换：如果随机变量U1,U2独立且U1,U2∼Uniform(0,1)，

Z0Z1=−2lnU1−−−−−−−√cos(2πU2)=−2lnU1−−−−−−−√sin(2πU2)

则Z0,Z1独立且服从标准正态分布。

注：George Edward Pelham Box，1919~2013，英国统计学家。伦敦学院大学博士，先后供职于普林斯顿大学和威斯康辛-麦迪逊大学。Ronald Aylmer Fisher的女婿。英国皇家学会会员，美国统计协会主席，数理统计学会（这是一个国际组织）主席。

Mervin Edgar Muller，俄亥俄州立大学教授。

统计力学与组合优化

MCMC和Gibbs Sampling最早都是统计力学的概念，后来才被用于机器学习领域。现将统计力学与组合优化的对应关系罗列如下：

统计力学	组合优化
样本	问题实例
状态（构形）	构形
能量	代价函数
温度	控制参数
基态能量	最小代价
基态构形	最小构形

频率统计学派 vs. 贝叶斯学派

对数学史感兴趣的朋友，可以看看陈希孺院士的《数理统计学简史》一书。rickjin文章的内容有相当部分取自该书。

注：陈希孺，1934～2005，数理统计学家。1956年毕业于武汉大学数学系，1997年当选为中国科学院院士。

该书中关于频率统计学派和贝叶斯学派的争议，引起了我的注意。

频率统计学派是所谓的正统派，由于其简单且便于理解的特点，多数入门级的数理统计学教程，一般都是按照该学派的思路写的。

而贝叶斯学派可谓另辟蹊径，它和频率统计学派的差异，参见《机器学习（七）》。由于该派系的思想比较新颖，我一度以为它和频率统计学派的关系，就犹如相对论之于经典力学。

然而，陈希孺院士告诉我们，两者各有优劣，尚未到一方决出胜负的阶段。比如，贝叶斯学派的先验估计，既是其成功的奥秘，也是其不成功的软肋。比如，对于“无信息先验分布”，目前尚处于“信则灵，不信则无”的境地。

陈院士的观点是：各取所长，为我所用。

PID算法

PID算法是自动控制领域最基础应用也最广泛的算法。它是三个单词proportional（P，比例）、integral（I，积分）和differential（D，微分）的缩写。

PID算法的流程如下图所示：

图1

以下我们以水箱水龙头为例，解释一下PID算法的各个要素。

场景1：假设我们面对的系统是一个简单的水箱的液位，要从空箱开始注水直到达到某个高度，而你能控制的变量是注水龙头的开关大小。

图1中的r(t)表示期待的设定值，在这个场景中，就是水箱的液位。y(t)表示当前的液位。e(t)=r(t)−y(t)表示误差值。

针对这个场景，我们可以设计如下算法：

水箱液位离预定高度远的时候开关开大点，离的近的时候开关就开小点，随着液位逐步接近预定高度逐渐关掉水龙头。用数学表示就是：

u(t)=Kpe(t)(1)

上式中的u(t)表示需要控制的量，在这里就是水龙头的开合度。Kp被称为比例系数。

场景2：假设这个水箱不仅仅是装水的容器了，还需要持续稳定的给用户供水。

以下用c表示给用户供水的量(c≥0)。显然如果使用公式1，则系统稳定时，u(t)−c=0，即Kpe(t)=c。

由于c和Kp都不为0，因此e(t)也不为0，这就导致始终无法加注到指定水位。这种稳态误差被称为静差。

为了平衡c，我们修改算法为：

u(t)=Kpe(t)+Ki∫t0e(τ)dτ(2)

Ki被称为积分系数。积分环节的意义就相当于你增加了一个水龙头，这个水龙头的开关规则是水位比预定高度低就一直往大了拧，比预定高度高就往小了拧。如果漏水速度不变，那么总有一天这个水龙头出水的速度恰好跟漏水的速度相等了，系统就和第一小节的那个一样了。那时，静差就没有了。这就是所谓的积分环节可以消除系统静差。

场景3：假设用户的用水量是变化的，即c(t)。

这时，如果仍采用公式2，则由于c(t)的变化，导致e(t)不恒为0。为了减小c(t)的变化，对e(t)的影响，我们继续修改算法：

u(t)=Kpe(t)+Ki∫t0e(τ)dτ+Kdde(t)dt(3)

Kd被称为微分系数。由于c(t)的变化不可能事先得知，因此，微分环节只能减小c(t)的变化所造成的影响，而不能消除。

公式3在Laplace域可写做：

L(s)=Kp+Ki/s+Kds(4)

从公式4可以看出，PID controller的数学原理和锁相环（Phase-locked loop）非常类似，它们实际上都是Feedback Control系统。



参考：

https://en.wikipedia.org/wiki/PID_controller

《Feedback Control of Dynamic Systems》，Gene F. Franklin，J David Powell，Abbas Emami-Naeini著。

注：Gene F. Franklin，1927~2012，美国控制论学家。哥伦比亚大学博士，斯坦福大学教授。

J David Powell，美国航空航天学家。斯坦福大学博士和教授。

Abbas Emami-Naeini，斯坦福大学博士和讲师，SC Solutions公司总监。

https://www.zhihu.com/question/23088613/answer/23942834

http://blog.163.com/suyanqiang0613@126/blog/static/1100531332012111475650867/

20世纪10大算法

2000年，IEEE评选出20世纪10大算法。名单如下：

1.Metropolis Algorithm for Monte Carlo

2.Simplex Method for Linear Programming

3.Krylov Subspace Iteration Methods

4.The Decompositional Approach to Matrix Computations

5.The Fortran Optimizing Compiler

6.QR Algorithm for Computing Eigenvalues

7.Quicksort Algorithm for Sorting

8.Fast Fourier Transform

9.Integer Relation Detection

10.Fast Multipole Method

详细内容参见：

http://www.uta.edu/faculty/rcli/TopTen/topten.pdf

中文版本：

http://blog.csdn.net/v_JULY_v/article/details/6127953

矩阵&向量的积

矩阵&向量有很多种积的定义，特罗列如下：

向量的点积

Dot product，又称inner product。

代数定义：

a⋅b=∑i=1naibi=a1b1+a2b2+⋯+anbn

几何定义：

a⋅b=∥a∥ ∥b∥cos(θ)

复变积分定义：

⟨ψ,χ⟩=∫baψ(x)χ(x)¯¯¯¯¯¯¯dx

矩阵的积

matrix product的定义如下（以3阶方阵为例）：

AB=⎛⎝⎜apubqvcrw⎞⎠⎟⎛⎝⎜αλρβμσγντ⎞⎠⎟=⎛⎝⎜aα+bλ+cρpα+qλ+rρuα+vλ+wρaβ+bμ+cσpβ+qμ+rσuβ+vμ+wσaγ+bν+cτpγ+qν+rτuγ+vν+wτ⎞⎠⎟

可以看出，积矩阵的每个元素是矩阵A、B相应行列向量的内积。

向量的向量积

Cross product是一个向量，其定义如下：

a×b=∥a∥∥b∥sin(θ) n

它还有个更出名的定义：

u×v=∣∣∣∣iu1v1ju2v2ku3v3∣∣∣∣

向量的混合积

triple product，又称mixed product或box product。

a⋅(b×c)≡det⎡⎣⎢a1b1c1a2b2c2a3b3c3⎤⎦⎥=det(a,b,c)

向量的笛卡尔积

Cartesian product实际上是集合论中的概念。

A×B={1,2}×{3,4}={(1,3),(1,4),(2,3),(2,4)}

向量的张量积

Tensor product，又称outer product。

u⊗v=uvT=⎡⎣⎢⎢⎢u1u2u3u4⎤⎦⎥⎥⎥[v1v2v3]=⎡⎣⎢⎢⎢u1v1u2v1u3v1u4v1u1v2u2v2u3v2u4v2u1v3u2v3u3v3u4v3⎤⎦⎥⎥⎥

可以看出，Tensor product和Cartesian product，虽然形式上都是各向量的组合，然而前者是2维的，而后者是1维的。

矩阵的张量积

张量积推广到矩阵，即所谓Kronecker product。

注：Leopold Kronecker，1823～1891，德国数学家。柏林大学博士和教授。导师Dirichlet。他最牛逼的地方是，当Riemann去世的时候，拒绝了哥廷根大学的offer。而这个位置的前任分别是：Carl Gauss、Dirichlet、Riemann。

A⊗B=⎡⎣⎢⎢a11B⋮am1B⋯⋱⋯a1nB⋮amnB⎤⎦⎥⎥

Hadamard product

又叫Schur product或entrywise product。

注：Jacques Salomon Hadamard，1865～1963，法国数学家。巴黎高等师范学校博士，波尔多大学教授。他曾于1936年访华，执教于清华大学。中国偏微分方程研究事业的主要创始人之一——吴新谋教授，就是他的学生。

⎛⎝⎜a11a21a31a12a22a32a13a23a33⎞⎠⎟∘⎛⎝⎜b11b21b31b12b22b32b13b23b33⎞⎠⎟=⎛⎝⎜a11b11a21b21a31b31a12b12a22b22a32b32a13b13a23b23a33b33⎞⎠⎟