1.累加累乘与积分

1. 将向量下标为偶数的分量累加, 写出相应表达式.
   $$\sum _{i\mathrm{m}\mathrm{o}\mathrm{d}2=0}{x}_i$$
2. 各出一道累加、累乘、积分表达式的习题, 并给出标准答案.
   ①将矩阵中大于等于1的分量平方并累加：
   $$\sum _{x_{ij}\ge 1}{x}_{ij}^2$$
   ②将矩阵中不为0的分量累乘：
   $$\prod _{x_{ij}\ne 0}{x}_{ij}$$
   ③函数$x^3-x+2$在0到1上积分：
   $${\int }_0^1{x}^3-x+2\mathrm{d}x$$
3. 你使用过三重累加吗? 描述一下其应用.
   在代码中三重循环算时间复杂度时会用到三重累加：$$\sum _{i=1}^n\sum _{j=1}^i\sum _{k=1}^{j}1$$
4. 给一个常用的定积分, 将手算结果与程序结果对比.
   定积分：${\int }_0^1{x}^2+1\mathrm{d}x$
   手算结果：$\frac{4}{3}=1.33333...$
   程序结果：$1.32835$
   

2.线性回归

1. 如何获得 $\mathbf{w}$?
   推导过程:
   $\begin{array}{rl}\Vert \mathbf{X}\mathbf{w}-\mathbf{Y}{\Vert }_2^2& =(\mathbf{X}\mathbf{w}-\mathbf{Y}{)}^{\mathrm{T}}(\mathbf{X}\mathbf{w}-\mathbf{Y})\\ & =({\mathbf{w}}^{\mathrm{T}}{\mathbf{X}}^{\mathrm{T}}-{\mathbf{Y}}^{\mathrm{T}})(\mathbf{X}\mathbf{w}-\mathbf{Y})\\ & ={\mathbf{w}}^{\mathrm{T}}{\mathbf{X}}^{\mathrm{T}}\mathbf{X}\mathbf{w}-{\mathbf{w}}^{\mathrm{T}}{\mathbf{X}}^{\mathrm{T}}\mathbf{Y}-{\mathbf{Y}}^{\mathrm{T}}\mathbf{X}\mathbf{w}+{\mathbf{Y}}^{\mathrm{T}}\mathbf{Y}\end{array}$
   将该式关于$\mathbf{w}$求导 (使用向量求导法则) 并令其为 0, 可得
   ${\mathbf{X}}^{\mathrm{T}}\mathbf{X}\mathbf{w}-{\mathbf{X}}^{\mathrm{T}}\mathbf{Y}=0$
   最后
   $\mathbf{w}=({\mathbf{X}}^{\mathrm{T}}\mathbf{X}{)}^{-1}{\mathbf{X}}^{\mathrm{T}}\mathbf{Y}$
   岭回归：
   $\underset{\mathbf{w}}{\mathrm{argmin}}=\Vert \mathbf{X}\mathbf{w}-\mathbf{Y}{\Vert }_2^2-\lambda \Vert \mathbf{w}{\Vert }_2^2$
   可推导出
   ${\mathbf{X}}^{\mathrm{T}}\mathbf{X}\mathbf{w}-{\mathbf{X}}^{\mathrm{T}}\mathbf{Y}+\lambda \mathbf{w}=0$
   最后
   $\mathbf{w}=({\mathbf{X}}^{\mathrm{T}}\mathbf{X}+\lambda \mathbf{I}{)}^{-1}{\mathbf{X}}^{\mathrm{T}}\mathbf{Y}$

2. 自己写一个小例子 $(n=3,m=1)$ 来验证最小二乘法.
   令函数 $y=2x+1$, $x$取值5, 6, 7, $\hat{y}$取值10, 13, 14.
   $\hat{y}=ax+b$, 用最小二乘法求$a,b$ :
   $$a=\frac{{\sum }_{i=1}^3{x}_i{\hat{y}}_i-3\overline{x}\overline{y}}{{\sum }_{i=1}^3{x}_i^2-3{\overline{x}}^2}=2.028$$
   $$b=\overline{y}-a\overline{x}=0.132$$

3.推导, 并描述这个方法的特点 (不少于 5 条).

$$\underset{\mathbf{w}}{\mathrm{argmax}}L(x)=\prod _{i=1}^{n}P(y_i|{\mathbf{x}}_i;\mathbf{w})$$

• 相乘计算困难, 将其求一个对数, 不改变单调性
  $$\begin{array}{rl}\log L(\mathbf{w})& =\sum _{i=1}^n\log P(y_i|\mathbf{x}i;\mathbf{w})\\ & =\sum _{i=1}^n{y}_i\log P(y_i=1|{\mathbf{x}}_i;\mathbf{w})+(1-y_i)\log (1-P(y_i=1|\mathbf{x}i;\mathbf{w}))\\ & =\sum _{i=1}^n{y}_i\log \frac{P(y_i=1|{\mathbf{x}}_i;\mathbf{w})}{1-P(y_i=1|{\mathbf{x}}_i;\mathbf{w})}+\log (1-P(y_i=1|\mathbf{x}i;\mathbf{w}))\\ & =\sum _{i=1}^n{y}_i{\mathbf{x}}_i\mathbf{w}-\log (1+e^{{\mathbf{x}}_i\mathbf{w}})\end{array}$$
• 对 $\mathbf{w}$ 求偏导
  $$\begin{array}{rl}\frac{\partial \log L(\mathbf{w})}{\partial \mathbf{w}}& =\sum _{i=1}^n{y}_i{\mathbf{x}}_i-\frac{e^{{\mathbf{x}}_i\mathbf{w}}}{1+e^{{\mathbf{x}}_i\mathbf{w}}}{\mathbf{x}}_i\\ & =\sum _{i=1}^n\left(y_i-\frac{e^{{\mathbf{x}}_i\mathbf{w}}}{1+e^{{\mathbf{x}}_i\mathbf{w}}}\right){\mathbf{x}}_i\end{array}$$
• 令该偏导为 0, 无法获得解析式, 因此用梯度下降.
  $${\mathbf{w}}^{t+1}={\mathbf{w}}^t-\alpha \frac{\partial \log L(\mathbf{w})}{\partial \mathbf{w}}$$

1. sigmoid函数更符合分类问题;
2. 将累乘变换为累加;
3. 输出为概率形式;
4. sigmoid函数任意阶可导;
5. 特征离散化，增加泛化能力.

前言：本篇系列文章参照高翔老师《视觉SLAM十四讲从理论到实践》的第四讲，讲解李群与李代数。写此篇的目的是为了补足《十四讲》中的数学类基础知识，原书内容有点像直译文献，因此笔者根据自己理解，对整体内容重新划分，更改了目录章节，改写部分公式，同时扩展了部分内容，有书的同学可以参照阅读，达到事半功倍的效果。

1. 整体误差最小化引出求导问题

使用李代数的一大动机是进行优化，而在优化过程中导数是非常必要的信息。下面来讨论一个带有李代数的函数，以及关于该李代数求导的问题，该问题有很强的实际应用背景。在SLAM中，我们使用$SO(3)$上的旋转矩阵或$SE(3)$上的变换矩阵来描述如何估计一个相机的位置和姿态（位姿）。不妨设某个时刻机器人的位姿为$\mathbf{T}$，它观察到一个世界坐标位于$p$的点，产生了一个观测数据$z$。那么，由坐标变换关系知：$$z=\mathbf{T}p+w.\text{\hspace{1em}(1.1)}$$其中，$w$为随机噪声。由于它的存在，$z$往往不可能精确地满足$z=Tp$的关系。所以，我们通常会计算理想的观测与实际数据的误差：$$e=z-\mathbf{T}p.\text{\hspace{1em}(1.2)}$$假设一共有$N$个这样的路标点和观测，于是就有$N$个上式。那么，对机器人进行位姿估计，相当于找到一个最优的$\mathbf{T}$，使得整体误差最小化：$$\underset{\mathbf{T}}{\min }J(\mathbf{T})=\sum _{i=1}^N{\Vert z_i-\mathbf{T}{p}_i\Vert }^2.\text{\hspace{1em}(1.3)}$$求解此问题，需要计算目标函数$J$关于变换矩阵$\mathbf{T}$的导数。我们把具体的算法留到后面讲（精力有限，本系列博客不再涉及，感兴趣的同学去读原书吧）。这里的重点是，我们经常会构建与位姿有关的函数，然后讨论该函数关于位姿的导数，以调整当前的估计值。然而，$SO(3),SE(3)$上并没有良好定义的加法，它们只是群，如果我们把$\mathbf{T}$当成一个普通矩阵来处理优化，就必须对它加以约束，就会比较繁琐。
而从李代数角度来说，由于李代数由向量组成，具有良好的加法运算，因此，可以使用李代数解决求导问题，求导思路分为两种：

1. 用李代数表示姿态，然后根据李代数加法对李代数求导。
2. 对李群左乘或右乘微小扰动，然后对该扰动对应的李代数求导，称为李群左扰动或右扰动模型。

第一种方式对应到李代数的求导，而第二种方式则对应到扰动模型。在讨论这两种思路的异同之前，我们先看一下BCH（Baker-Campbell-Hausdorff）线性近似。

2. BCH公式与近似形式

2.1 BCH公式

下面来考虑一个问题。虽然我们已经清楚了$SO(3)$和$SE(3)$上的李群和李代数关系，但是当在$SO(3)$中完成两个矩阵乘法时，李代数中$\mathfrak{s}\mathfrak{e}(3)$上发生了什么改变呢？反过来说，当$\mathfrak{s}\mathfrak{o}(3)$做两个李代数的加法时，$SO(3)$上是否对应着两个矩阵的乘积？如果成立，相当于：$$\exp ({\varphi }_1^{\wedge })\exp ({\varphi }_2^{\wedge })=\exp (({\varphi }_1+{\varphi }_2{)}^{\wedge })$$如果${\varphi }_1,{\varphi }_2$为标量，那么显然该式成立；但此处我们计算的是矩阵的指数函数，而非标量的指数。换言之，我们在研究下式是否成立：$$\ln (\exp (A)\exp (B))=A+B?$$很遗憾，该式在矩阵时并不成立。两个李代数指数映射乘积的完整形式，由Baker-Cambell-Hausdorff公式（BCH公式）¹给出。由于其完整形式较复杂，我们只给出其展开式的前几项：$$\begin{array}{rl}\ln (\exp (A)\exp (B))& =A+B+\frac{1}{2}[A,B]+\frac{1}{12}[A,[A,B]]-\frac{1}{12}[B,[A,B]]-\frac{1}{24}[B,[A,[A,B]]]+\cdots \\ & =A+B+\frac{1}{2}(AB-BA)+\frac{1}{12}(A^{2}B+A{B}^2-2ABA+B^{2}A+B{A}^2-2BAB)+\\ & \frac{1}{24}(A^2{B}^2-2ABAB-B^2{A}^2+2BABA)\end{array}.\text{\hspace{1em}(2.1)}$$其中[]为李括号。BCH公式告诉我们，当处理两个矩阵指数之积时，它们会产生一些由李括号组成的余项。

2.2 BCH线性近似

特别的，考虑$SO(3)$上的李代数$(\ln \exp ({\varphi }_1^{\wedge })\exp ({\varphi }_2^{\wedge }){)}^{\vee }$，当${\varphi }_1$或${\varphi }_2$为小量时，小量二次以上的项可以被忽略，而旋转或变换更新的基本都是小量，因此BCH具有实际应用意义。此时，BCH拥有线性近似表达²

$$(\ln \exp ({\varphi }_1^{\wedge })\exp ({\varphi }_2^{\wedge }){)}^{\vee }\approx \{\begin{array}{c}{\mathbf{J}}_l({\varphi }_2{)}^{-1}{\varphi }_1+{\varphi }_2当{\varphi }_1为小量,\\ {\mathbf{J}}_r({\varphi }_1{)}^{-1}{\varphi }_2+{\varphi }_1当{\varphi }_2为小量.\end{array}\text{\hspace{1em}(2.2)}$$
以第一个近似为例。该式告诉我们，当对一个旋转矩阵$R_2$（李代数为${\varphi }_2$）左乘一个微小旋转矩阵$R_1$（李代数为${\varphi }_1$）时，可以近似地看作，在原有的李代数${\varphi }_2$上加了一项${\mathbf{J}}_l({\varphi }_2{)}^{-1}{\varphi }_1$。同理，第二个近似描述了右乘一个微小位移的情况。于是，李代数在BCH近似下，分成了左乘近似和右乘近似两种，在使用时我们须注意使用的是左乘模型还是右乘模型。
本文以左乘为例。左乘BCH近似雅克比${\mathbf{J}}_l$事实上就是《李群与李代数1》中的式$(4.6)$的内容：$${\mathbf{J}}_l=\mathbf{J}=\frac{\sin \theta }{\theta }\mathbf{I}+(1-\frac{\sin \theta }{\theta })\mathbf{a}{\mathbf{a}}^T+\frac{1-\cos \theta }{\theta }{\mathbf{a}}^{\wedge }.\text{\hspace{1em}(2.3)}$$它的逆为：$${\mathbf{J}}_l^{-1}=\frac{\theta }{2}\cot \frac{\theta }{2}\mathbf{I}+(1-\frac{\theta }{2}\cot \frac{\theta }{2})\mathbf{a}{\mathbf{a}}^T+\frac{\theta }{2}{\mathbf{a}}^{\wedge }.\text{\hspace{1em}(2.4)}$$而右乘雅克比仅需要对自变量取负号即可：$${\mathbf{J}}_r(\varphi )={\mathbf{J}}_l(-\varphi ).\text{\hspace{1em}(2.5)}$$这样，我们就可以谈论李群乘法与李代数加法的关系了。

2.3 BCH近似的意义

为了方便读者理解，我们重新叙述BCH近似的意义。假定对某个旋转$R$，对应的李代数为$\varphi$，我们给它左乘一个微小旋转，记作$\Delta R$，对应的李代数为$\Delta \varphi$。那么，在李群上，得到的结果就是$\Delta R\cdot R$，而在李代数上，根据BCH近似表达公式$(2.2)$，为：${\mathbf{J}}_l(\varphi {)}^{-1}\Delta \varphi +\varphi$。合并起来，可以简单的写成：$$\exp (\Delta {\varphi }^{\wedge })\exp ({\varphi }^{\wedge })=\exp \left(({\mathbf{J}}_l(\varphi {)}^{-1}\Delta \varphi {)}^{\wedge }+\varphi \right).\text{\hspace{1em}(2.6)}$$同理可得右乘扰动$R\cdot \Delta R$对应的表达式为：$$\exp ({\varphi }^{\wedge })\exp (\Delta {\varphi }^{\wedge })=\exp \left(({\mathbf{J}}_r(\varphi {)}^{-1}\Delta \varphi {)}^{\wedge }+\varphi \right).\text{\hspace{1em}(2.7)}$$反之，如果我们在李代数上进行加法，让一个$\varphi$加上$\Delta \varphi$，那么可以近似为李群上带左右雅克比的乘法：$$\exp \left((\varphi +\Delta \varphi {)}^{\wedge }\right)=\exp \left(({\mathbf{J}}_l(\varphi )\Delta \varphi {)}^{\wedge }\right)\exp ({\varphi }^{\wedge })=\exp ({\varphi }^{\wedge })\exp \left(({\mathbf{J}}_r(\varphi )\Delta \varphi {)}^{\wedge }\right).\text{\hspace{1em}(2.8)}$$这就为之后李代数上做微分提供了理论基础。

同样地，对于$SE(3)$，也有类似的BCH近似：$$\exp (\Delta {\xi }^{\wedge })\exp ({\xi }^{\wedge })\approx \exp \left(({ȷ}_l^{-1}\Delta \xi +\xi {)}^{\wedge }\right).\text{\hspace{1em}(2.9)}$$$$\exp ({\xi }^{\wedge })\exp (\Delta {\xi }^{\wedge })\approx \exp \left(({ȷ}_r^{-1}\Delta \xi +\xi {)}^{\wedge }\right).\text{\hspace{1em}(2.10)}$$这里$ȷ$形势比较复杂，它是一个$6\times 6$的矩阵，读者可以参考文献[2]中式$(7.82)$和$(7.83)$的内容。
反之，如果我们在李代数上进行加法，让一个$\xi$加上$\Delta \xi$，那么可以近似为李群上带左右雅克比的乘法：$$\exp \left((\xi +\Delta \xi {)}^{\wedge }\right)\approx \exp \left(({\mathbf{ȷ}}_l(\xi )\Delta \xi {)}^{\wedge }\right)\exp ({\xi }^{\wedge })=\exp ({\xi }^{\wedge })\exp \left(({\mathbf{ȷ}}_r(\xi )\Delta \xi {)}^{\wedge }\right).\text{\hspace{1em}(2.11)}$$
下面我们正式开始讨论李代数的求导、李群扰动对应李代数的扰动模型求导两种思路的异同。

3. 微分模型——李代数求导

首先，考虑SO(3)上的情况。假设我们对一个空间点$p$进行了旋转，得到$Rp$。现在，计算旋转之后点的坐标相对于旋转的导数，我们非正式地记为：$$\frac{\partial (Rp)}{\partial R}.\text{\hspace{1em}(3.1)}$$由于$SO(3)$没有加法，所以该导数无法按照导数的定义进行计算。设$R$对应的李代数为$\varphi$，我们转而计算³：

$$\frac{\partial (\exp ({\varphi }^{\wedge })p)}{\partial \varphi }$$按照导数的定义有：$$\begin{array}{rl}\frac{\partial (\exp ({\varphi }^{\wedge })p)}{\partial \varphi }& =\underset{\delta \varphi \to 0}{\lim }\frac{\exp ((\varphi +\delta \varphi {)}^{\wedge })p-\exp ({\varphi }^{\wedge })p}{\delta \varphi }\\ & =\underset{\delta \varphi \to 0}{\lim }\frac{\exp (({\mathbf{J}}_{\mathbf{l}}\delta \varphi {)}^{\wedge })\exp ({\varphi }^{\wedge })p-\exp ({\varphi }^{\wedge })p}{\delta \varphi }\\ & =\underset{\delta \varphi \to 0}{\lim }\frac{(\mathbf{I}+({\mathbf{J}}_{\mathbf{l}}\delta \varphi {)}^{\wedge })\exp ({\varphi }^{\wedge })p-\exp ({\varphi }^{\wedge })p}{\delta \varphi }\\ & =\underset{\delta \varphi \to 0}{\lim }\frac{({\mathbf{J}}_{\mathbf{l}}\delta \varphi {)}^{\wedge }(\exp ({\varphi }^{\wedge })p)}{\delta \varphi }\\ & =\underset{\delta \varphi \to 0}{\lim }\frac{-(\exp ({\varphi }^{\wedge })p{)}^{\wedge }{\mathbf{J}}_{\mathbf{l}}\delta \varphi }{\delta \varphi }\\ & =-(Rp{)}^{\wedge }{\mathbf{J}}_{\mathbf{l}}.\end{array}$$推导说明：第二行的变换为BCH线性近似，第三行为泰勒展开舍去高阶项后的近似（由于取了极限，可以写等号），第四行至第五行将反对称符号看做叉积，交换之后变号，即$a^{\wedge }b=-b^{\wedge }a$。于是，我们推导出了旋转后的点相对于李代数的导数：$$\frac{\partial Rp}{\partial \varphi }=(-Rp{)}^{\wedge }{\mathbf{J}}_{\mathbf{l}}.\text{\hspace{1em}(3.2)}$$不过，由于这里仍然含有形势比较复杂的${\mathbf{J}}_{\mathbf{l}}$，我们不希望计算它。而下面要讲的扰动模型则提供了更简单的导数计算公式。

4. 扰动模型求导（左乘）

4.1 SO(3)上的扰动模型求导

另一种求导方式是对R进行一次扰动$\Delta R$，然后求结果相对于扰动的变化率。这个扰动可以乘在左边也可以乘在右边，最后结果会有一点儿微小的差异，我们以左扰动为例。设$\Delta R$对应的李代数为$\phi$，$R$左乘微小旋转$\Delta R$后就是$\Delta R\cdot R$。
这里需要着重强调的是，$\Delta R\cdot R$非常重要，理解了它就等于理解了扰动，本章的扰动模型求导就是根据$\Delta R\cdot R$进行的。
根据BCH近似，其对应的李代数就是$\phi +\varphi$，对$\phi$求导得：$$\begin{array}{rl}\frac{\partial (Rp)}{\partial \phi }& =\underset{\phi \to 0}{\lim }\frac{(\Delta R\cdot R)p-Rp}{\phi }\\ & =\underset{\phi \to 0}{\lim }\frac{\exp ({\phi }^{\wedge })\exp ({\varphi }^{\wedge })p-\exp ({\varphi }^{\wedge })p}{\phi }\\ & =\underset{\phi \to 0}{\lim }\frac{(\mathbf{I}+{\phi }^{\wedge })\exp ({\varphi }^{\wedge })p-\exp ({\varphi }^{\wedge })p}{\phi }\\ & =\underset{\phi \to 0}{\lim }\frac{{\phi }^{\wedge }\exp ({\varphi }^{\wedge })p}{\phi }\\ & =\underset{\phi \to 0}{\lim }\frac{-(Rp{)}^{\wedge }\phi }{\phi }=-(Rp{)}^{\wedge }.\end{array}\text{\hspace{1em}(3.3)}$$推导说明：第一行是求导公式的直接展开，第二行直接将李群替换为对应的李代数，第三行为泰勒展开舍去高阶小项后的近似（取极限后写等号），第四到第五行即$a^{\wedge }b=-b^{\wedge }a$。可见，相比于直接对李代数求导，省去了一个雅克比${\mathbf{J}}_{\mathbf{l}}$的计算。这使得扰动模型更为实用。请读者务必理解这里的求导运算，这在位姿估计中具有重要的意义。

4.2 SE(3)上的扰动模型求导

下边，我们给出$SE(3)$上的扰动模型。假设某空间点p经过一次变换$\mathbf{T}$（对应李代数为$\xi$），得到$\mathbf{T}p$⁴。现在给$\mathbf{T}$左乘一个扰动$\Delta T=\exp ({\chi }^{\wedge })$，我们设扰动项的李代数为$\chi =[\Delta \rho ,\Delta \varphi {]}^T$，那么：

$$\begin{array}{rl}\frac{\partial \mathbf{T}p}{\partial \chi }& =\underset{\chi \to 0}{\lim }\frac{(\Delta T\cdot T)p-Tp}{\chi }\\ & =\underset{\chi \to 0}{\lim }\frac{\exp ({\chi }^{\wedge })\exp ({\xi }^{\wedge })p-\exp ({\xi }^{\wedge })p}{\chi }\\ & =\underset{\chi \to 0}{\lim }\frac{(I+{\chi }^{\wedge })\exp ({\xi }^{\wedge })p-\exp ({\xi }^{\wedge })p}{\chi }\\ & =\underset{\chi \to 0}{\lim }\frac{{\chi }^{\wedge }\exp ({\xi }^{\wedge })p}{\chi }\\ & =\underset{\chi \to 0}{\lim }\frac{\left[\begin{array}{cc}\Delta {\varphi }^{\wedge }& \Delta \rho \\ 0^T& 0\end{array}\right]\left[\begin{array}{c}Rp+t\\ 1\end{array}\right]}{\chi }\\ & =\underset{\chi \to 0}{\lim }\frac{\left[\begin{array}{c}\Delta {\varphi }^{\wedge }(Rp+t)+\Delta p\\ 0^T\end{array}\right]}{[\Delta \rho ,\Delta \varphi {]}^T}\\ & =\left[\begin{array}{cc}\mathbf{I}& -(Rp+t{)}^{\wedge }\\ 0^T& 0^T\end{array}\right]\stackrel{\underset{\mathrm{d}\mathrm{e}\mathrm{f}}{}}{=}(Tp{)}^{\odot }.\end{array}\text{\hspace{1em}(3.4)}$$我们把最后的结果定义成一个算符$\odot$⁵，它把一个齐次坐标的空间点变换成一个$4\times 6$的矩阵。
上式推导中，需要稍微解释的是第五行的矩阵求导顺序，假设$\mathbf{a},\mathbf{b},\mathbf{x},\mathbf{y}$都是列向量，那么在我们的符号写法下，有如下规则：$$\frac{d\left[\begin{array}{c}a\\ b\end{array}\right]}{d\left[\begin{array}{c}x\\ y\end{array}\right]}={\left(\frac{d[a,b{]}^T}{d\left[\begin{array}{c}x\\ y\end{array}\right]}\right)}^T={\left[\begin{array}{cc}\frac{da}{dx}& \frac{db}{dx}\\ \frac{da}{dy}& \frac{db}{dy}\end{array}\right]}^T=\left[\begin{array}{cc}\frac{da}{dx}& \frac{da}{dy}\\ \frac{db}{dx}& \frac{db}{dy}.\end{array}\right]\text{\hspace{1em}(3.5)}$$

至此，我们已经介绍了李群与李代数上的微分运算。下面再补充一点其它性质。

4.3 伴随性质

对于公式：$$R{p}^{\wedge }{R}^T=(Rp{)}^{\wedge }.\text{\hspace{1em}(3.6)}$$有类似的公式：$$R\exp (p^{\wedge })R^T=\exp ((Rp{)}^{\wedge }).\text{\hspace{1em}(3.7)}$$该式称为$SO(3)$上的伴随性质。同样地，在$SE(3)$上也有伴随性质：$$T\exp ({\xi }^{\wedge })T^{-1}=\exp ((Ad(T)\xi {)}^{\wedge }),\text{\hspace{1em}(3.8)}$$其中：$$Ad(T)=\left[\begin{array}{cc}R& t^{\wedge }R\\ 0& R\end{array}\right].$$本性质是书中习题，感兴趣的同学可以研究下。

5. 相似变换群相关

最后，我们要提一下在单目视觉中使用的相似变换群$Sim(3)$，它对应的李代数，以及他们之间的指数/对数映射关系，最后看一下它对应的扰动模型。如果你只对双目SLAM或RGB-D SLAM感兴趣，可以跳过本节。因仿射变换和射影变换不再具备指数/对数映射性质和扰动模型，所以不再考虑它们。

5.1 相似变换群$Sim(3)$

我们已经介绍过单目的尺度不确定性。如果在单目SLAM中使用$SE(3)$表示位姿，那么由于尺度不确定性与尺度漂移，整个SLAM过程中的尺度会发生变化，这在$SE(3)$中未能体现出来。因此，在单目情况下，我们一般会显示地把尺度因子表达出来。用数学语言来说，对于位于空间的点$p$，在相机坐标系下要经过一个相似变换，而非欧氏变换：$$p^{\prime }=\left[\begin{array}{cc}s\mathbf{R}& \mathbf{t}\\ 0^T& 1\end{array}\right]p=s\mathbf{R}\mathbf{p}+\mathbf{t}.\text{\hspace{1em}(5.1)}$$在相似变换中，我们把尺度$s$表达出来了，它同时作用在$p$的3个坐标之上，对$p$进行了一次缩放。
与$SO(3)、SE(3)$类似，相似变换也对矩阵乘法构成群，称为相似变换群$Sim(3)$：$$Sim(3)=\left\{\mathbf{S}=\left[\begin{array}{cc}s\mathbf{R}& \mathbf{t}\\ 0^T& 1\end{array}\right]\in {\mathbb{R}}^{4\times 4}\right\}.\text{\hspace{1em}(5.2)}$$同样地，$Sim(3)$也有对应的李代数、指数映射、对数映射、扰动模型等，下边分别介绍。

5.2 李代数$\mathfrak{s}\mathfrak{i}\mathfrak{m}(3)$及指数/对数映射

5.2.1 李代数$\mathfrak{s}\mathfrak{i}\mathfrak{m}(3)$

李代数$\mathfrak{s}\mathfrak{i}\mathfrak{m}(3)$是一个7维向量$\zeta$，它的前6维与$\mathfrak{s}\mathfrak{e}(3)$相同，最后多了一项相似变换系数$\sigma$：$$\mathfrak{s}\mathfrak{i}\mathfrak{m}(3)=\left\{\mathbf{\zeta }=\left[\begin{array}{c}\mathbf{\rho }\\ \mathbf{\varphi }\\ \sigma \end{array}\right]\in {\mathbb{R}}^7,\mathbf{\rho }\in {\mathbb{R}}^3,\mathbf{\varphi }\in \mathfrak{s}\mathfrak{o}(3),{\mathbf{\zeta }}^{\wedge }=\left[\begin{array}{cc}\sigma \mathbf{I}+{\mathbf{\varphi }}^{\wedge }& \mathbf{\rho }\\ 0^T& 0\end{array}\right]\in {\mathbb{R}}^{4\times 4}\right\}.\text{\hspace{1em}(5.3)}$$

5.2.2 李代数到李群的指数映射

关联$Sim(3)$和$\mathfrak{s}\mathfrak{i}\mathfrak{m}(3)$的仍是指数映射和对数映射。先看李代数$\mathfrak{s}\mathfrak{i}\mathfrak{m}(3)$到李群$Sim(3)$的指数映射：$$\exp ({\mathbf{\zeta }}^{\wedge })=\left[\begin{array}{cc}{e}^{\sigma }\exp ({\mathbf{\varphi }}^{\wedge })& {\mathbf{J}}_s\mathbf{\rho }\\ 0^T& 1\end{array}\right].\text{\hspace{1em}(5.4)}$$其中：$$\begin{array}{rl}{\mathbf{J}}_s& =\frac{e^{\sigma }-1}{\sigma }\mathbf{I}+\frac{\sigma e^{\sigma }\sin \theta +(1-e^{\sigma }\cos \theta )\theta }{{\sigma }^2+{\theta }^2}{\mathbf{a}}^{\wedge }+\\ & \left(\frac{e^{\sigma }-1}{\sigma }-\frac{(e^{\sigma }\cos \theta -1)\sigma +(e^{\sigma }\sin \theta )\theta }{{\sigma }^2+{\theta }^2}\right){\mathbf{a}}^{\wedge }{\mathbf{a}}^{\wedge }\end{array}$$通过上式可知，已知李代数$\zeta$，它与李群中元素的对应关系为：$$s=e^{\sigma },\mathbf{R}=\exp ({\mathbf{\varphi }}^{\wedge }),\mathbf{t}={\mathbf{J}}_s\mathbf{\rho }.\text{\hspace{1em}(5.5)}$$说明：李群中的尺度因子$s$即为李代数中$\sigma$的指数函数；旋转部分与$SO(3)$中一致；平移部分，在$\mathfrak{s}\mathfrak{e}(3)$中需要乘一个雅克比$\mathbf{J}$，而相似变换中是更为复杂的${\mathbf{J}}_s$。

5.2.3 李群到李代数的对数映射

已知李群$Sim(3)$时，先将矩阵简化为如下形式：左下角最后一行前6个元素为0，最后一个元素为1。求解各元素的步骤叙述如下：

1. $\sigma$。对于矩阵左上角$3\times 3$的行列式，由于$\Vert \mathbf{R}\Vert =1$，故行列式的值就是尺度因子$s$，此时$\sigma =ln(s)$；
2. $\mathbf{\varphi }$。 简化后的李群$Sim(3)$的左上角$3\times 3$矩阵除以尺度因子$s$，得到旋转矩阵$\mathbf{R}$，根据$\mathbf{R}$和罗德里格斯公式可求得$\mathbf{a}$和$\theta$，进而由$\mathbf{\varphi }=\mathbf{a}\theta$求得$\mathbf{\varphi }$；
3. $\mathbf{\rho }$。由$\sigma 、\mathbf{a}、\theta$可求得${\mathbf{J}}_s$，简化后的李群$Sim(3)$右上角一列为$\mathbf{t}$，此时$\mathbf{\rho }=\mathbf{t}/{\mathbf{J}}_s$。

最后由$\sigma$、$\mathbf{\varphi }$和$\mathbf{\rho }$组合，得到对应的李代数$\mathbf{\zeta }$。

5.3 扰动模型

$Sim(3)$的BCH近似与$SE(3)$是类似的，我们可以讨论一个点$p$经过相似变换$Sp$后，相对于$S$的导数。同样地，存在微分模型和扰动模型两种方式，而扰动模型不用计算雅克比，较为简单。因此我们省略推导过程，直接给出扰动模型的结果。
设给予$Sp$左侧一个小扰动$\exp ({\eta }^{\wedge })$，并求$Sp$对于扰动的导数。因为$Sp$是4维的齐次坐标，$\eta$是7维向量，所以该导数应该是$4\times 7$的雅克比。方便起见，记$Sp$的前3维组成向量为$q$，那么：$$\frac{\partial \mathbf{S}\mathbf{p}}{\partial \mathbf{\eta }}=\left[\begin{array}{ccc}\mathbf{I}& -{\mathbf{q}}^{\wedge }& \mathbf{q}\\ 0^T& 0^T& 0\end{array}\right].\text{\hspace{1em}(5.6)}$$此式可对比公式$(3.4)$理解记忆。关于$Sim(3)$，我们就介绍到这里。关于其更详细的资料，建议读者参见文献[3]。

5.4 总结

本讲《李群与李代数》引入了李群$SO(3)$和$SE(3)$，以及它们对应的李代数$\mathfrak{s}\mathfrak{o}(3)$和$\mathfrak{s}\mathfrak{e}(3)$。我们介绍了位姿在他们上面的表达和转换，然后通过BCH的线性近似，就可以对位子进行扰动并求导。这给之后讲解位姿的优化打下了理论基础，因为我们需要经常对某一个位姿的估计值进行调整，使它对应的误差减小，因此只有弄清楚如何对位姿进行调整和更新之后，才能继续下一步的内容。

本讲内容比较偏理论，讲的过程也非常精简，速度也相对快，但本讲内容是解决后续许多问题的基础，请读者务必理解，特别是位姿估计部分。值得一提的是，除了李代数，同样也可以用四元数、欧拉角等方式表示旋转，只是后续处理要麻烦一些。在实际应用中，也可以使用$SO(3)$加上平移的方式来代替$SE(3)$，从而回避一些雅克比的计算。

为方便理解，下面来看两个具体的例子。第一个例子是李群和李代数的基本操作，第二个例子是关于最小误差的计算。

6. 编码实战

6.1 李群与李代数的运算库：Sophus

我们已经介绍了李代数的入门知识，现在是通过实践演练并巩固所学知识的时候了。第3讲中的Eigen提供了几何模块，但没有提供李代数的支持。一个较好的李代数库Strasdat维护的Sophu库⁶，它支持本章讨论的$SO(3)$和$SE(3)$，此外，还有二维运动$SO(2)$和$SE(2)$及相似变换$Sim(3)$的内容。它是直接在Eigen基础上开发的，因此我们不需要安装额外的依赖库。读者可以直接从GitHub上获取Sophus，原书目录slambook/3rdparty也提供Sophus源代码。

由于历史原因，Sophus早期版本只提供双精度的李群/李代数类，后续版本改写成了模板类。模板类的Sophus可以使用不同精度的李群/李代数，但同时也增加了使用难度。本书中使用带模板类的Sophus库，Sophus本身也是一个cmake工程，它只需编译即可，无需安装。

下面来演示Sophus库中的$SO(3)$和$SE(3)$运算：

#include<iostream>
#include<cmath>
#include<eigen3/Eigen/Core>
#include<eigen3/Eigen/Geometry>
#include "sophus/se3.hpp"

using namespace std;
using namespace Eigen;
using namespace Sophus;

int main(int argc, char **argv) {
    Matrix3d R = AngleAxisd(M_PI/2, Vector3d(0, 0, 1)).toRotationMatrix();
    Quaterniond q(R);
    Sophus::SO3d SO3_R(R);
    Sophus::SO3d SO3_q(q);
    cout<<"SO(3) from matrix:\n"<<SO3_R.matrix()<<endl;
    cout<<"SO(3) from quaternion:\n"<<SO3_q.matrix()<<endl;
    cout<<"they are equal"<<endl;
    
    Vector3d so3 = SO3_R.log();
    cout<<"so3 = "<<so3.transpose()<<endl;
    cout<<"so3 hat = \n"<<Sophus::SO3d::hat(so3)<<endl;
    cout<<"so3 hat vee = "<<Sophus::SO3d::vee(Sophus::SO3d::hat(so3)).transpose()<<endl;
    
    Vector3d update_so3(1e-4, 0, 0);
    Sophus::SO3d SO3_updated = Sophus::SO3d::exp(update_so3) * SO3_R;
    cout<<"SO3 updated = \n" <<endl<<SO3_updated.matrix()<<endl;
    
    cout<<"******************************"<<endl;
    
    Vector3d t(1, 0, 0);
    Sophus::SE3d SE3_Rt(R, t);
    Sophus::SE3d SE3_qt(q, t);
    cout<<"SE(3) from R, t:\n"<<SE3_Rt.matrix()<<endl;
    cout<<"SE(3) from q, t:\n"<<SE3_qt.matrix()<<endl;
    cout<<"they are equal"<<endl;
    
    typedef Eigen::Matrix<double, 6, 1> Vector6d;
    Vector6d se3 = SE3_Rt.log();
    cout<<"se3 = "<<se3.transpose()<<endl;
    cout<<"se3 hat = \n"<<Sophus::SE3d::hat(se3)<<endl;
    cout<<"se3 hat vee = "<<Sophus::SE3d::vee(Sophus::SE3d::hat(se3)).transpose()<<endl;
    
    Vector6d update_se3;
    update_se3.setZero();
    update_se3(0, 0) = 1e-4;
    cout<<"update_se3 = \n" <<endl<<update_se3.matrix()<<endl;
    Sophus::SE3d SE3_updated = Sophus::SE3d::exp(update_se3) * SE3_Rt;
    cout<<"SE3 updated = \n" <<endl<<SE3_updated.matrix()<<endl;
    
    return 0;
}

该演示程序分为两部分：前半部分介绍$SO(3)$上的操作，后半部分为$SE(3)$。我们演示如何构造$SO(3),SE(3)$对象，对它们进行指数、对数映射，以及当知道更新量后，如何对李群元素进行更新。为了编译它，请在CMakeLists.txt添加以下几行代码：

cmake_minimum_required(VERSION 3.0)

project(useSophus)

# 为使用 sophus，需要使用find_package命令找到它
find_package(Sophus REQUIRED)
find_package(Pangolin REQUIRED)
include_directories(${Sophus_INCLUDE_DIRS})
include_directories(${Pangolin_INCLUDE_DIRS})

add_executable(trajectoryError trajectoryError.cpp)
target_link_libraries(trajectoryError ${Pangolin_LIBRARIES})

include_directories("/usr/include/eigen3")
include_directories("/home/frank/slam14/slambook2/3rdparty/Sophus")

add_executable(useSophus useSophus.cpp)

# 把工程调为Debug编译模式
set(CMAKE_BUILD_TYPE "Debug")

find_package命令是cmake提供的寻找某个库的头文件和库文件的指令。如果cmake能找到它，就会提供头文件和库文件所在目录的变量，在这个例子中，就是Sophus_INCLUDE_DIRS。基于模板的Sophu库和Eigen库一样，是仅含头文件而没有源文件的，根据它们，我们就能将Sophus库引入自己的cmake工程。程序运行结果如下：

SO(3) from matrix:
2.22045e-16          -1           0
          1 2.22045e-16           0
          0           0           1
SO(3) from quaternion:
2.22045e-16          -1           0
          1 2.22045e-16           0
          0           0           1
they are equal
so3 =      0      0 1.5708
so3 hat = 
      0 -1.5708       0
 1.5708       0      -0
     -0       0       0
so3 hat vee =      0      0 1.5708
SO3 updated = 

          0          -1           0
          1           0     -0.0001
     0.0001 2.03288e-20           1
******************************
SE(3) from R, t:
2.22045e-16          -1           0           1
          1 2.22045e-16           0           0
          0           0           1           0
          0           0           0           1
SE(3) from q, t:
2.22045e-16          -1           0           1
          1 2.22045e-16           0           0
          0           0           1           0
          0           0           0           1
they are equal
se3 =  0.785398 -0.785398         0         0         0    1.5708
se3 hat = 
        0   -1.5708         0  0.785398
   1.5708         0        -0 -0.785398
       -0         0         0         0
        0         0         0         0
se3 hat vee =  0.785398 -0.785398         0         0         0    1.5708
update_se3 = 

0.0001
     0
     0
     0
     0
     0
SE3 updated = 

2.22045e-16          -1           0      1.0001
          1 2.22045e-16           0           0
          0           0           1           0
          0           0           0           1

6.2 评估轨迹的误差

在实际工程中，我们经常需要评估一个算法的估计轨迹与真实轨迹的差异来评价算法的精度。真实轨迹往往通过某些更高精度的系统获得，而估计轨迹则是由待评价的算法计算得到的。在三维刚体运动中，我们演示了如何显示存储在文件中的某条轨迹，本节我们考虑如何计算两条轨迹的误差。

考虑一条估计轨迹${\mathbf{T}}_{esti,i}$和真实轨迹（ground truth）${\mathbf{T}}_{gt,i}$，其中$i=1,\cdots ,N$，那么我们可以定义一些误差指标来描述它们之间的差别。常见的误差指标有很多种，介绍如下：

1. 绝对轨迹误差（Absolute Trajectory Error, ATE），计算公式如下：$$AT{E}_{all}=\sqrt{\frac{1}{N}\sum _{i=1}^N{\Vert \log ({\mathbf{T}}_{gt,i}^{-1}{\mathbf{T}}_{esti,i}{)}^{\vee }\Vert }^2}.\text{\hspace{1em}(6.1)}$$这实际上是每个位姿李代数的均方根误差（Root-Mean-Squared Error，RMSE），它可以刻画两条轨迹的旋转和平移误差。
2. 平均平移误差（Average Translational Error），也有的文献仅考虑平移误差，计算公式如下：$$AT{E}_{trans}=\sqrt{\frac{1}{N}\sum _{i=1}^N{\Vert trans({\mathbf{T}}_{gt,i}^{-1}{\mathbf{T}}_{esti,i})\Vert }^2}.\text{\hspace{1em}(6.2)}$$其中trans表示取括号内部变量的平移部分。因为从整条轨迹来看，旋转出现误差后，随后也会在平移上出现误差，所以以上两种指标在实际中都适用。
3. 相对位姿误差（Relative Pose Error，RPE），除此之外，还可以定义相对误差。例如，考虑$i$时刻到$i+\Delta t$时刻的运动，那么RPE可定义为：$$RP{E}_{all}=\sqrt{\frac{1}{N}\sum _{i=1}^N{\Vert \log {\left(({\mathbf{T}}_{gt,i}^{-1}{\mathbf{T}}_{gt,i+\Delta t}{)}^{-1}({\mathbf{T}}_{esti,i}^{-1}{\mathbf{T}}_{esti,i+\Delta t})\right)}^{\vee }\Vert }^2}.\text{\hspace{1em}(6.3)}$$同样地，也可只取平移部分：$$RP{E}_{trans}=\sqrt{\frac{1}{N}\sum _{i=1}^N{\Vert trans{\left(({\mathbf{T}}_{gt,i}^{-1}{\mathbf{T}}_{gt,i+\Delta t}{)}^{-1}({\mathbf{T}}_{esti,i}^{-1}{\mathbf{T}}_{esti,i+\Delta t})\right)}^{\vee }\Vert }^2}$$利用Sophus库，很容易实现这部分计算。下面我们演示绝对轨迹误差的计算。在这个例子中，我们有groundtruth.txt和estimated.txt两条轨迹，下面的代码江都区这两条轨迹，计算RMSE然后显示到3D窗口中：

#include<iostream>
#include<fstream>
#include<unistd.h>
#include<pangolin/pangolin.h>
#include<sophus/se3.hpp>

using namespace Sophus;
using namespace std;

string groundtruth_file = "../groundtruth.txt";
string estimated_file = "../estimated.txt";

typedef vector<Sophus::SE3d, Eigen::aligned_allocator<Sophus::SE3d>> TrajectoryType;

void DrawTrajectory(const TrajectoryType &gt, const TrajectoryType &esti);

TrajectoryType ReadTrajectory(const string &path);

int main(int argc, char **argv) {
    TrajectoryType groundtruth = ReadTrajectory(groundtruth_file);
    TrajectoryType estimated = ReadTrajectory(estimated_file);
    cout<<groundtruth_file<<endl;
    assert(!groundtruth.empty() && !estimated.empty());
    assert(groundtruth.size() == estimated.size());
    
    double rmse = 0;
    for (size_t i=0; i<estimated.size(); i++) {
        Sophus::SE3d p1 = estimated[1], p2 = groundtruth[i];
        double error = (p2.inverse() * p1).log().norm();
        rmse += error * error;
    }
    rmse = rmse/double(estimated.size());
    rmse = sqrt(rmse);
    cout<<"RMSE = "<<rmse<<endl;
    
    DrawTrajectory(groundtruth, estimated);
    return 0;
}

TrajectoryType ReadTrajectory(const string& path) {
    ifstream fin(path);
    TrajectoryType trajectory;
    if (!fin) {
        cerr<<"trajectory: "<<path<<" not found."<<endl;
        return trajectory;
    }
    
    while (!fin.eof()) {
        double time, tx, ty, tz, qx, qy, qz, qw;
        fin>>time>>tx>>ty>>tz>>qx>>qy>>qz>>qw;
        Sophus::SE3d p1(Eigen::Quaterniond(qx, qy, qz, qw), Eigen::Vector3d(tx, ty, tz));
        trajectory.push_back(p1);
    }
    
    return trajectory;
}

void DrawTrajectory(const TrajectoryType &gt, const TrajectoryType &esti) {
  // create pangolin window and plot the trajectory
  pangolin::CreateWindowAndBind("Trajectory Viewer", 1024, 768);
  glEnable(GL_DEPTH_TEST);
  glEnable(GL_BLEND);
  glBlendFunc(GL_SRC_ALPHA, GL_ONE_MINUS_SRC_ALPHA);

  pangolin::OpenGlRenderState s_cam(
      pangolin::ProjectionMatrix(1024, 768, 500, 500, 512, 389, 0.1, 1000),
      pangolin::ModelViewLookAt(0, -0.1, -1.8, 0, 0, 0, 0.0, -1.0, 0.0)
  );

  pangolin::View &d_cam = pangolin::CreateDisplay()
      .SetBounds(0.0, 1.0, pangolin::Attach::Pix(175), 1.0, -1024.0f / 768.0f)
      .SetHandler(new pangolin::Handler3D(s_cam));


  while (pangolin::ShouldQuit() == false) {
    glClear(GL_COLOR_BUFFER_BIT | GL_DEPTH_BUFFER_BIT);

    d_cam.Activate(s_cam);
    glClearColor(1.0f, 1.0f, 1.0f, 1.0f);

    glLineWidth(2);
    for (size_t i = 0; i < gt.size() - 1; i++) {
      glColor3f(0.0f, 0.0f, 1.0f);  // blue for ground truth
      glBegin(GL_LINES);
      auto p1 = gt[i], p2 = gt[i + 1];
      glVertex3d(p1.translation()[0], p1.translation()[1], p1.translation()[2]);
      glVertex3d(p2.translation()[0], p2.translation()[1], p2.translation()[2]);
      glEnd();
    }

    for (size_t i = 0; i < esti.size() - 1; i++) {
      glColor3f(1.0f, 0.0f, 0.0f);  // red for estimated
      glBegin(GL_LINES);
      auto p1 = esti[i], p2 = esti[i + 1];
      glVertex3d(p1.translation()[0], p1.translation()[1], p1.translation()[2]);
      glVertex3d(p2.translation()[0], p2.translation()[1], p2.translation()[2]);
      glEnd();
    }
    pangolin::FinishFrame();
    usleep(5000);   // sleep 5 ms
  }

}

代码中文本文件groundtruth.txt和estimated.txt请从博客资源部分下载。
该程序输出结果为RMSE = 3.09706。读者也可以尝试将旋转部分去掉，仅计算平移部分的误差。这个例子中的数据已经经过预处理，包括轨迹的时间对齐、外参预估等，请读者知悉。程序运行后的绘制图如下：

后记：

SLAM系列到此完结。笔者认为数学部分难以理解的部分已经补充完毕，剩下的部分，高翔老师的书讲解的没那么混乱了，仔细阅读可以理解。后边十讲的部分，笔者会以思维导图的形式继续发博客总结，喜欢的童鞋记得关注我吧。

参考文献：

1. 《视觉SLAM十四讲：从理论到实践》，高翔、张涛等著，中国工信出版社。
2. T.Barfoot, State estimation for robotics: A matrix lie group approach, 2016.
3. H. Strasdat, Local accuracy and global consistency for efficient visual slam. PhD thesis, Citeseer, 2012.

---

1. 参见https://groupprops.subwiki.org/wiki/Baker-Campbell-Hausdorff_formula。 ↩︎

2. 对BCH具体形式和近似表达的具体推导，本文不作具体讨论，请参考文献[2]。 ↩︎

3. 严格来说，在矩阵微分中，只能求行向量关于列向量的导数，所得结果是一个矩阵。但本文写成列向量对列向量的导数，读者可以认为先对分子进行转置，再对最后结果进行转置，这使得式子变得简洁。在这种意义下，可以认为$\frac{d(Ax)}{dx}=A$ ↩︎

4. 请注意，为了使乘法成立，$p$必须使用齐次坐标。 ↩︎

5. 作者将其读作“咚”，像一个石子掉在井里的声音。 ↩︎

6. 最早提出李代数的是Sophus Lie，因此库以他的名字命名。 ↩︎

本文来源于个人的凸优化学习笔记参考cs229 cvxoptnote,写成笔记的原因仅仅是想通过个人的笔记自己讲述与推导一下这些数学公式，内容可能会很简单，强力建议想得到一手资料的人好好学习文末参考资料

凸集合

定义就直接跳过了，这里简单写一些常见的凸集

• 凸集的交, 设 $C_i,i=1,2,3,...,n$ 是凸集，那么我们有${\cap }_{i=0}^n{C}_i$ 也是凸集
• affine subspace 即 $Ax+b=0,x\in R^n$, 这既是 convex,也是concave
• 凸函数的sublevel 即 i f   f ( x )   i s   c o n v e x   f u n c t i o n , t h e n { f ( x ) &lt; = α } i s   c o n v e x   s e t if \ f(x)\ is\ convex\ function,then \{f(x) &lt;=\alpha\} is\ convex\ set if f(x) is convex function,then{f(x)<=α}is convex set

凸函数

Define

性质，两条重点

1. first-order-approximation
   
   $$f(y)\ge f(x)+(\nabla f(x){)}^T(y-x)$$
2. 二阶导数大于0
   
   $${\nabla }_x^{2}f(x)\ge 0$$

一些简单的凸函数

1. Affine funciton
   
   $f(x)=A^{T}x+b,x\in R^n$
2. quadratic function
   
   $f:R^n->R,f(x)=2^{-1}{x}^{T}Ax+b^{T}x+c,A\in S^n,b\in R^n,c\in R$
   
   这里简单说一下 ${\nabla }^{2}f(x)=A$, 因此如果 $A$ 是半正定(positive-semidefinite)的,那么$f(x)$ is convex
3. Norm,e.g $f(x)=||x|{|}_2^2$
4. 凸函数的非负权重和
   
   $f(x)={\sum }_{i=0}^n{c}_i{f}_i(x),f_i(x)\in convexfunctionclass$
5. 凸函数集合的最大值函数
   
   $f(x)=\max f_i(x),f_i(x)\in convexfunctionclass$
   
   证明:
   $$\begin{array}{rl}f(\theta x+(1-\theta )y)& =\max f_i(\theta x+(1-\theta )y)\\ & \le \max \theta f_i(x)+(1-\theta )f_i(y)\\ & \le \theta f(x)+(1-\theta )f(y)\end{array}$$

凸优化问题

define

凸优化问题的局部最优解一定是全局最优解

接下来是凸优化问题最重要，也是最有趣的一个性质

证明

我们用反证法，
设 $x\in R^n,f(x)\in convexfunctionclass$,且， $x$ 是 凸优化问题 $f(x)$ 的一个局部最优解，那么我们有，在$x$ 的邻域中,i.e.$||x-z|{|}_2<R$ ,不存在一个点 $z$,使得$f(z)<f(x)$, 同时，假设 $x$ 不是全局最优解，那么我们有,在定义域中，存在一点 $y,s.t.f(y)<f(x)$

设 $z=\theta y+(1-\theta )x,\theta =\frac{R}{2||x-y|{|}_2}$,可以证明$||x-z|{|}_2=2^{-1}R<R$

同时,

$$\begin{array}{rl}f(z)& =f(\theta y+(1-\theta )x)\\ & \le \theta f(y)+(1-\theta )f(x)\\ & \le f(x)\end{array}$$

矛盾

一些机器学习的example

1. 线性回归
   
   $J(\theta )=||W^T\Phi (X)-Y|{|}_2^2=\sum (w^t\varphi (x_i)-y_i{)}^2$
   因为 $f_i(w)=(w^t\varphi (x_i)-y_i{)}^2$ 是一个复合affine 凸函数(i.e. x^2 是凸函数，而 $w^t\varphi (x_i)-y_i$ 是对 $w$ 的affine function)
2. logistic reg
   
   $J(\theta )=-\sum y_i\log g(\theta x_i)+(1-y_i)\log 1-g(\theta x_i)$
3. SVM
   
   SVM 的目标函数是一个 QP(quadratic opt) 问题，$O(n^3)$
   这个问题还引出一些凸问题有趣的性质，比如lagrange duality 和 KKT,我们后面分析

lagrange & KKT

lagrange mutiplier

对于上面的凸优化问题来说，我们有它的lagrange function

$$\begin{array}{rl}\mathcal{L}(x,\alpha ,\beta )& =f(x)+\sum {\alpha }_i{g}_i(x)+\sum {\beta }_i{h}_i(x)\\ & s.t.{\alpha }_i>=0\end{array}$$

其中 $x$ 称为 primal variable, $\alpha ,\beta$ 是dual variable

NOTE 其实可以将 $\alpha ,\beta$ 理解成打破原来约束的代价

primal & dual 原问题与对偶问题

中间那段 english 不用看，直接看形式就行了

我们下面试图证明这样一个弱对偶性(weak duality)

即 $d^{\ast }={\theta }_D$