一步一步走向锥规划 - LP
发布网友
发布时间：2022-11-10 01:16
共1个回答
热心网友
时间：2023-11-23 04:53
引子
彼得堡学派是人类近代史的一个奇迹， 是由切比雪夫（ChebyShev）一手开创， 切比雪夫研究的领域包括， 数论、概率论、函数*近论、积分学， 这些都在这个学派极大开花， 而切比雪夫之所以能成为这个开创者， 一个重要原因就是他从小精通法语， 并且大胆积极的和法国数学家沟通交流。 另外一个重要原因就是善于培养学生，马尔可夫、李亚普诺夫都是他的学生， 都是一等一的数学大师， 这绝不是偶然。这个学派对概率和各种不等式玩的出神入化， 对近代机器学习的发展影响巨大！
帕夫努季·利沃维奇·切比雪夫本人， 出身于贵族家庭．他的祖辈中有许多人立过战功． 父亲列夫·帕夫洛维奇·切比雪夫参加过抵抗拿破仑(Napoleon)入侵的卫国战争， 母亲阿格拉费娜·伊万诺夫娜·切比雪娃也出身名门，他们共生育了五男四女，切比雪夫排行第二．切比雪夫的法语就是来自母亲精心的教育！ 题外话， 日本的一般妇女都会二门外语左右， 这对子女的教育都是极大的优势！
切比雪夫对数学的兴趣，来自阅读并思考了欧几里得(Euclid)《几何原本》(Elements)。 题外话， 现在各大中小学， 从来没有时间给学生阅读数学原著和思考， 整天重复做题， 还指望培养一流大数学家， 简直不可能。
有个年轻的犹太俄罗斯数学家Leonid Vitaliyevich Kantorovich， 受到了这个学派的极大影响， 在概率， 不等式， 数值*近方法做出了极大贡献。 并且获得诺贝尔经济学奖。
Kantorovich对一般性线性规划进行了归纳，并且最先提出了解法。 从此线性规划最先走向了经济领域，走向了世界。
什么是LP， 线性规划
标准形式
LP正如其名，目标是求线性乘积的最值， 并且受到线性不等式的约束。标准形式如下：
举个例子， 经过标准化就是上面我们形式了：
把上面的所有线性要求图形化， 根据不等式知道，定义域为*的凸区域：
图形化求解这个比较容易，因为目标公式是线性函数， 所以只要考虑凸区域的最值点， 又因为凸区域本身就是直线围起来的， 所以只要考虑角就行了
这样可以得到在x=5, y=5这个点， LP问题取到最大值。
但是我们稍微讨论一下最优值的分布， 一般情况下， 不等式围成的区域有三种情况：
空集合（null set）：这种情况天然没有最优解（optimum）。
有限可行区域： Bounded Feasible Region
有唯一解（unique optimal）
有无数解（multiple optimal）
无限可行区域： Unbounded Feasible Region
无最优解
有唯一解（unique optimal）
有无数解（multiple optimal）
顺便， 我们说明一下带权求和（weighted sum）的几何含义： 对于点积p个点和它们的系数的点积。
我们讨论三种*情况：
1. 如果常量之和为1：
这是X被称为仿射求和（affine combination）： 是两点直线上所有点
2. 如果所有常数为非负：
这是X被称为锥求和（conic combination）：是两个点向量方向之间所有可能方向可能的向量。
3. 如果同时满足affine combination 和 conic combination，
这时候被称为凸求和（convex combination）：是两个点之间线段上所有点
有了凸求和的定义， 那么我们就可以通过凸求和来判断一个点是不是角：
角只能是自己和自己的凸求和，或者角不是任何其他两个点的凸求和。
如果一个点可以表示成两个点的凸求和， 那么这个点就不是角。
这样，我们给出了LP的一般定义，以及图形化解法。
几个经典问题
在LS问题里面， 我们把误差定义为欧式空间距离（Euclidean Distance）， 如果我们换成曼哈顿距离（Manhattan Distance）， 那么我们就得到Least Absolute Deviation Regression （LADR）。
利用Subgradient求导， 我们可以得到LS类似的结论：
通过简单的变换， 我们就能把LADR问题转换成LP问题了：
除了基于Manhattan距离的回归， 还有基于Chebyshev距离的回归， 他们都是典型的LP问题：
LP的几何意义
从矩阵空间看
LP的矩阵空间比较简单， 只有线性的乘法， 因此只有简单的矩阵空间的关系，意义不是很大。
从仿射空间看（affine space）
Affine space和前面的subspace最大的区别就是， Affine只保留了相对欧式距离关系。 举个例子， 下面图里有两层子区域， 在原点的可以构成suspace， 但是上层的不能构成一个subspace， 因为不能满足线性生成（span），包括线性组合和缩放，依然属于这个区域， 具体可以看“ 矩有四子 ”中子空间的定义。 但是这个区域依然保留着相对关系，是一个仿射空间。
为什么要有仿射空间呢？
有了仿射空间的定义， 我们就可以对仿射空间进行变换。 否则只能在固定的投影空间（Projection）里面， 就需要对4个点都要计算， 因为我们要保留子空间的特性（绝对位置关系）。 但是在仿射空间里面， 我们只要保留相对位置空间， 所以只需要对3个点进行定位确定相对线性关系。
有了上面的了解， 我们再来看一个我们已经很常用的affine combination。 在下面的二维空间里面， 点D就是A, B, C三点之间的相对关系计算出来的点， 而不考虑A点自己的向量关系。
有了仿射空间， 我们在求LP问题的时候， 根据仿射空间和线性目标， 我们就知道， 极值点只能在仿射空间区域的角（Corner）取到， 这样就能构建单纯形Simplex区域来进行计算了。
Simplex就是通过空间点线性仿射连接（affine combination）构成的一个凸包区域。
LP概述
在Kantorovich提出线性规划的问题之后， 美国数学家Dantzig在1947年提出了单纯形算法（Simplex），Dantzig是来自波兰的美国数学家， 他的博士导师Jerzy Neyman相当牛， 是美国统计系（伯克利统计系）的第一个创建者。
之后，Fiacco 和 McCormick在60年代提出了早期的插值算法（interior-point methods, IPM）， 到了70年代， 人们又提出了椭面法（ellipsoid method）， 和一些subgradient methods。Ellipsoid 算法最大意义是首次提出多项式时间（polynomial-time）的算法，但是实际上这些方法都不如Simplex算法好用，Simplex算法被称为20世纪10大算法之一， 直到1984年来自印度的数学家Karmarkar提出了多项式时间（polynomial-time）的插值点算法（interior-point methods, IPM)， 一下子改变了Simplex一家独大的状态。Karmarkar就是来自传说中的IIT的。
从此interior-point methods慢慢和Simplex算法并肩成为主流， 并且Simplex属于LP专业定制，很难扩展， 而Karmarkar的算法后来扩展到了非线性（nonlinear）问题中继续使用。
因此LP算法主流主要是Dantzig的Simplex和Karmarkar的IPM算法:
1994 Simplex Method：大部分中小规模问题超快，大规模问题不如IPM，偶尔情况很差。
1984Interior-point Method, IPM：大规模情况比较快，中小规模不如Simplex， 但是不会出现最坏情况。
Simplex 算法概述
Simplex思想
算法的思想比较直观，就是要在Simplex上从任意一个顶点（角）出发， 然后找到一条路，它的终点就是最优值点。 这个算法利用前面提到的最优点必然落在角上的结论。
为什么不能遍历所有的角（Corner）来找到最优点？
在n维度的空间里面，假设有k个切面围成一个多面体（polyhedron）， 那么这个多面体最多有多少个角（Vertex）？ 有如下结论：
这个结论可以通过在基本的Simplex上面添加切面来验证：
从上面的结论，我们可以知道， 这是随着k的增长指数增长的值。 所以， 根据LP的特征， 每个不等式就是一个切面， 然后把可能区域（feasible domain）围成一个凸的多面体， 而最优点一定在角上， 但是遍历所有角（Vertex）的数目是根不等式和维度个数而指数增长的， 是不可以取的策略。
在LP中， 角的数学上如何表达的？
假设角的集合为P，并且满足如下条件：
那么如果X是角点， 对A有什么要求（什么情况下是交点是角）？
为了不失一般性， 把x向量中，零向量部分分开表示。
在这种情况下， 我们把A也进行拆分成前p列为
。 那么：
那么X是角（Extreme Point）的要求是，对应的
矩阵是満秩的，或者矩阵的所有列是线性独立的（linearly independent）。
采用反证法， 证明如下：
已经知道角点X，要求
矩阵満秩， 那么満秩是不是一定是角点呢？同样反证法， 我们只要证明非角点一定不满秩就好了。
所以， 在Simplex算法中从一个角点跳到一个相邻的角点， 都必须满足非零部分对应的列是相互独立的（矩阵是満秩的）。
为了更好表示， 我们引入basis 和 non-basis来进行划分：
这样， 我们首先得找一个基本解（basic solution， bs）， 让它满足
。 如果更进一步， 还能找到
。 那么我就找到一个基本可行解了（basic feasible solution， bfs），也就是说找到一个角了。
从下面图上可以直观的看到基本解 和基本可行解的差别了， 可行基本解使我们想要的凸多面体的角（CDEM）， 而基本解还包含其多面体区域外的角。 另外（BTW），可行解（Feasible Solutions, fs）就是凸多面体内部区域（四边形CDEM和内部）。
有了上面角的数学表达形式这个结论， 我们很容易得到：
1. 如果A不是満秩（Full Rank）的， 那么就没有角， 那么就无解了。
2. 角的个数小于从列数（n）里面选出行数（m）个列的可能性。
虽然前面我们已经给出了更为精确的结论， 但是它的由来并不清楚， 而这里我们很容易推导出一个上限（bfs是bs的一个子集， bfs数量小于bs数量）。
邻居顶点是如何表示的？
有了上面的结论， 我们在来看一下如何通过矩阵来表示邻居顶点？
既然是邻居， 那么就是要共享大部分*条件，最多只有一个*条件不一致。 那么在矩阵里面， 就是要两个基集合（basic set）只有一个列不一样。
在这样的情况下， 从一个点移到另外一个邻居点就是选择一个新的非基点（non-basic）加入到基集合里面来， 同时从基里面选择原有的一个列放回， 如下图所示：
那么如何通过矩阵做到这一点呢？ 就要通过高斯消去法（Gaussian Elimination, GE）。 举个极端的例子， 大家还记通过GE来求矩阵逆么？
我们会看到整个basis从联合矩阵右边移到左边， 然后就求解到了逆矩阵。
注意， 高斯消去法（GE）应用到计算机求解的时候经常会出现round-off error，或者divide by zero， 所以会每次选择pivot点进行操作。Pivot点选择也很简单就是把一列上最大值的行换到对应位置， 避免除以0或者小的数导致误差。
为什么需要松弛形式（slack form）？
前面我们讨论过， 顶点的形式要满足：
但是标准形式下面， 并不满足这两个条件。 因此Simplex这个算法利用的并不是LP的标准形式， 而是松弛形式（slack form ）或者说增大形式（Augmented form）， 这个形式就是通过引入松弛变量（slack）， 来把矩阵乘法的不等式变成等式。 然后所有定义域变量都放到第一象限。 从而满足我们先前设定的角的简单形式。
Simplex算法的怎么设定优化目标的？
有了松弛形式之后，就是求解目标是什么？
Simplex 的目标就是基于等式， 在一个象限里面化解到一个目标表达式：
Z + c1 * x1 + c2 * x2 + c3 * x3 = b
等价于：
Z = b - c1 * x1 - c2 * x2 - c3 * x3
使得c1, c2, c3...的符号相同。
由于 X, Y, Z >= 0。 因此：
Z取最小值： 要求c1, c2, c3 <0，Z和x1, x2, x3单调性一致，x1=x2=x3=0的时候是最小值
Z取最大值： 要求c1, c2, c3 >0，Z和x1, x2, x3单调性相反，x1=x2=x3=0的时候是最小值
有了目标了， 那么求解方式是什么？
首先我们图示化来看一下， 从一个点到最优点（optimum）理想的路径的启发。 例如下图， 最理想的情况应该吧把直线Z=4x+3y，沿着梯度上升（Gradient ascent）的方向平移， 一直移到D点， 因为D和直线仅有一个交点了， 这种情况下Z取到最大值。
但是Simplex算法是从角集合之间进行移动的。 假如， 目前在A点，A点是一个基本可行点（bfs）， 因为A点是角，所以是基本点（bs）， 另外，A点在*区域内，所以是可行点（feasible solution）。 A点又两个邻居， B和E点。 那么B点和E点哪个更好呢？
根据前面最优的梯度上升的路径， 其实我们应该选择和梯度上升更为一致的方向。 那么B和E点的方向， 哪个方向更为和梯度方向一致呢？ 我们可以把梯度分解到X轴和Y轴， 看到梯度在X轴上的分量是4，而在Y轴上的分量是3， 那么可以确定沿着X轴更为合适， 因为梯度在这个方向上投影更大， 之间的夹角更小，更为一致。
如果我们再把梯度的分量对应到等式z=4x+3y的系数， 就知道4是x的系数（coefficient）， 而3是y的系数。 所以直观上来说， 我们要选系数更大的方向。
根据前面讲述的利用高斯消去法（GE）来进行从一个bfs点， 移到邻居的bfs点。 这时候， 也要选择pivot列进行操作。 前面我们提到， 只要想办法从basic 和non-sc set 里面替换任意一列，就是变换到邻居了。 这里我们又加上， 要选择最大系数的方向操作。 因此pivot列就要选择系数最大的列优先进行GE操作。
如何选择Pivot点进行高斯消去法来交换basic列？
前面我们已经讲述了基本思想， 这里我们举个例子来进行说明：
对于标准的LP算法， 求如下最大值：
对应的图如下，所示：
首先， 我们引入Slack变量， 换成Slack形式：
这样我们得到对应的Ax = b，x > 0的矩阵形式：
（注意： 目标公式化成了-40 x1 - 30 x2 + P = 0）
根据前面的矩阵基的理论， 很容易把slack变量对应的列选出来作为基：
在这种情况下， 很容易设定non-basic变量为0， 然后basic 变量就是b值。
我们发现这时候对应的是坐标0点， 图上A点。
根据我们之前讲述的理解， 我们应该往系数最大的x1轴（x2=0）的方向移动， 因此我们找到对应列， 为pivot列。
选择x1轴对应的列为pivot列， 如下图， pivot列是P行里面除了P列之外的最小的列（取了负号的最大系数）。
那么选择哪一行进行操作了呢？ 根据高斯消去法里面Pivot行的思想， 我们应该选择pivot列里面最大值的行进行操作， 但是， 这里不是， 这里是选择正的 rhs / s 最小的行进行操作（等价于把rhs全部变成1， 然后每行除以b之后，对应的最大值的行）。如下图所示， 我们选择了24/3最小的行， 或者说选择了3/24最大列，依然符合pivot的思想。 当然这只是从pivot思想出发的阐述， 具体后面会说明为什么要选择正的值的最小pivot ratio。
通过高斯消去（GE）之后， 得到如下：
这样， 我们得到了新的基集合和非基集合， 我们看到是对应的x2=0和s3=0的直线的交点：
而这个点就是我们要找的G点：
如果继续这个过程， 我们就会再找到x2是最大系数方向， s2是pivot选择的行：
进行GE后，我们得到如下：
对应的新的基和非基集合是：
所以找到了s2=0和s3=0的交点：
