PCL 点云中的数学

函数求导

方差&协方差矩阵
基本概念
方差（Variance）衡量的是单个随机变量的变化(比如一个人在群体中的身高)，概率论中方差用来度量随机变量和其数学期望（即均值）之间的偏离程度。

标准差（Standard Deviation）是方差的算术平方根，用σ表示。标准差能反映一个数据集的离散程度。

协方差（Covariance）衡量的是两个随机变量同时变化的离散程度(比如一个人在群体中的身高和体重)。方差是协方差的一种特殊情况。

协方差矩阵
协方差矩阵表示一组随机变量之间的两两线性相关性，其中的每个元素代表了两个随机变量之间的协方差。

(3)

对角线上的元素为各个随机变量的方差，非对角线上的元素为两两随机变量之间的协方差，根据协方差的定义

(4)

我们一起通过二维点分布，探究协方差矩阵与数据分布特点的关系。

为了简化图像，我们令均值向量

的方差均为1，我们可以生成如下随机分布，共200个点：

import numpy as np
import matplotlib.pyplot as plt
# %matplotlib inline

plt.style.use('ggplot')
plt.rcParams['figure.figsize'] = (12, 8)

# Normal distributed x and y vector with mean 0 and standard deviation 1
x = np.random.normal(0, 1, 200)
y = np.random.normal(0, 1, 200)
X = np.vstack((x, y))

plt.scatter(X[:, 0], X[:, 1])
plt.title('Generated Data')
plt.axis('equal')
plt.show()

我们可以通过如下代码计算随机生成的点的真实协方差矩阵

# Covariance
def cov(x, y):
    xbar, ybar = x.mean(), y.mean()
    return np.sum((x - xbar)*(y - ybar))/(len(x) - 1)

# Covariance matrix
def cov_mat(X):
    return np.array([[cov(X[0], X[0]), cov(X[0], X[1])],
                     [cov(X[1], X[0]), cov(X[1], X[1])]])

# Calculate covariance matrix
print(cov_mat(X)) # (or with np.cov(X))

[[ 1.08970816 -0.03960157]
 [-0.03960157  0.99871813]]

即协方差矩阵主轴上的数据描述了其各个维度的离散程度（分布情况），值越大，分散的越开。

线性变换
接下来我们为所有点进行线性变换，来观察数据分布的变化以及协方差矩阵的变化。

只进行缩放操作：

代码实现：

import numpy as np
import matplotlib.pyplot as plt

plt.style.use('ggplot')
plt.rcParams['figure.figsize'] = (12, 8)

# Normal distributed x and y vector with mean 0 and standard deviation 1
x = np.random.normal(0, 1, 200)
y = np.random.normal(0, 1, 200)
X = np.vstack((x, y)) # 2xn

sx, sy = 0.5, 2.0
Scale = np.array([[sx, 0], [0, sy]])

Y = Scale.dot(X)

原始点集

plt.scatter(X[0, :], X[1, :])

缩放后点

plt.scatter(Y[0, :], Y[1, :])
plt.title('Generated Data')
plt.axis('equal')
plt.show()
先缩放，后旋转：

代码实现：

import numpy as np
import matplotlib.pyplot as plt

plt.style.use('ggplot')
plt.rcParams['figure.figsize'] = (12, 8)

# Normal distributed x and y vector with mean 0 and standard deviation 1
x = np.random.normal(0, 1, 200)
y = np.random.normal(0, 1, 200)
X = np.vstack((x, y)) # 2xn

sx, sy = 0.5, 2.0
Scale = np.array([[sx, 0], [0, sy]])

# Rotation matrix
theta = np.pi / 6
c, s = np.cos(theta), np.sin(theta)
Rot = np.array([[c, -s], [s, c]])

# Transformation matrix
T = Rot.dot(Scale)

# Apply transformation matrix to X
Y = T.dot(X)

# 原始点集
plt.scatter(X[0, :], X[1, :])
# 缩放、旋转后
plt.scatter(Y[0, :], Y[1, :])
plt.title('Generated Data')
plt.axis('equal')
plt.show()

# Calculate covariance matrix
print(cov_mat(Y)) # (or with np.cov(Y))
其协方差矩阵输出结果为：


[[ 1.1830068  -1.66028421]
 [-1.66028421  3.10071868]]

可验证与之前代码输出的协方差结果近似。现在我们已经得到协方差矩阵与缩放矩阵、旋转矩阵的关系，由于协方差矩阵是正定对称矩阵即可对角化矩阵，则接下来我们对协方差矩阵进行特征值分解，观察其特征向量矩阵、特征值矩阵与缩放矩阵、旋转矩阵的关系

特征分解
根据线性代数相关定义，对于任意可对角化矩阵，都存在一个特征分解（eigen decomposition）

我们通过以下代码对Σ
进行特征分解、对结果进行验证：

import numpy as np

sigma = np.array([
    [19 / 16, -15 * np.sqrt(3) / 16,],
    [-15 * np.sqrt(3) / 16, 49 / 16]
])

eVa, eVe = np.linalg.eig(sigma)
print("eigen value:\n", eVa)
print("eigen vector:\n", eVe)

输出

eigen value - SS^T:
 [0.25 4.  ]
eigen vector - U:
 [[-0.8660254  0.5      ]
 [-0.5       -0.8660254]]

可验证该输出结果与上述一致。

我们可以将特征值和特征向量分别相乘，得到变换后椭圆区域的两个基向量：

# Calculate covariance matrix
C = cov_mat(Y)
print(C)  # (or with np.cov(Y))

eVa, eVe = np.linalg.eig(C)

plt.scatter(Y[:, 0], Y[:, 1])
for value, eigen in zip(eVa, eVe.T):
    plt.plot(
        [0, 3 * np.sqrt(value) * eigen[0]],
        [0, 3 * np.sqrt(value) * eigen[1]],
         lw=5)

plt.title('Transformed Data')
plt.axis('equal')
plt.show()

输出：

总结：

一个多元正态分布的概率密度函数中，其协方差矩阵
包含了分布的旋转角度和缩放信息，其均值向量
包含了密度中心的位置信息。

雅克比矩阵
理解雅克比向量

将多个维度作为整体，来理解函数与变量的关系。

雅克比矩阵即通过线性代数的向量来描述高纬度的数据在空间中的表现。

通过向量将每一个函数的变量变化（偏导数）合并为一个整体，在高维空间去观察变化。

(7)

雅克比寻找极值

当函数进入极大值或极小值时，其梯度会比较平缓，接近或等于0。当极值只有一个时，找到对应的x，y位置比较容易。

当函数比较复杂时，极值可能有多个，会有全局最高的山顶、全局最低的山谷、局部山顶、局部的山谷。这些地方的梯度也都是0。如果将x，y等于0，则无法确定是哪个极值。

如果我们任意选择一点，无法知道最高峰和最低谷的位置。但是可以知道其最近的较高点和较低点。雅克比Vector的方向指向了最近的较高点，其模表述了幅度。即可以帮我们高效地找到局部的较高和较低点，但不能帮我们找到全局的最高点和最低点。

黑塞矩阵
理解黑塞矩阵

数学上，黑塞矩阵（Hessian matrix又称海森矩阵）是一个自变量为向量的实值函数的二阶偏导数组成的方阵。即可以理解为对Jacobian Vector再次进行求导，把每一项偏导项再做一次求导。其矩阵是对称矩阵。

那么如何在J处于极值情况下，判断这个极值是最大值、最小值、还是鞍值？这就需要使用黑塞矩阵了。

对于二次函数，黑塞矩阵为常数矩阵，与参数无关。

(8)

黑塞矩阵常通过牛顿法来处理大规模优化问题。