LS Approximation(Least-squares)

Question: 给定N个点$x_i\in\mathbb{R}^d,where\ i\in[1…N]$.希望得到$f(x)$能近似出点$x_i$处的$f_i$值。
误差函数为$J_{LS}=\sum_i||f(X_i)-f_i||^2$。那么得到最小化问题：
$$
min_{f_x\in\prod_m^d}\sum_i||f(\mathbb{x}_i)-f_i||^2
$$
$f$ 取自$\prod_m^d$,d维空间上为total degree为m的多项式空间，则可以写作：
$$
f(x)=b(x)^Tc=b(x)\cdot c
$$
其中：$\mathbb{b}(x)=[b_1(x),…,b_k(x)]^T$是多项式的基向量，$\mathbb{c}=[c_1,…,c_k]^T$是需要最小化求解的未知系数向量。

多项式的基：

$(m=2,d=2)，\mathbb{b(x)}=[1,x,y,x^2,xy,y^2]^T$
$\mathbb{R^3}(m=1,d=3)$的线性拟合，$\mathbb{b(x)}=[1,x,y,z]^T$
任意维度常数拟合，$\mathbb{b(x)}=[1]$

Solution: 误差函数$J_{LS}$的偏导为0.$\bigtriangledown J_{LS}=0,\bigtriangledown=[\partial/\partial_{c_1},…,\partial/\partial_{c_k}]$

因此，可以得到一个线性方程组(LSE)：
$$
\begin{align}
\partial J_{LS}/\partial c_1=0 &: \sum_i2b_1(\mathbb{x_i})[\mathbb{b(x_i)}^T\mathbb{c}-f_i]=0\\
\partial J_{LS}/\partial c_2=0 &: \sum_i2b_2(\mathbb{x_i})[\mathbb{b(x_i)}^T\mathbb{c}-f_i]=0\\
&\ \vdots\\
\partial J_{LS}/\partial c_k=0 &: \sum_i2b_k(\mathbb{x_i})[\mathbb{b(x_i)}^T\mathbb{c}-f_i]=0\\
\end{align}
$$
可以向量表示为：
$$
\begin{aligned}
& 2\sum_i\mathbb{[b(x_i)b(x_i)}^T\mathbb{c}-\mathbb{b(x_i)}f_i]=0 \\
& \Rightarrow \\
&c=[\sum_ib(x_i)b(x_i)^T]^{-1}\sum_ib(x_i)f_i
\end{aligned}
$$
Example.

在$\mathbb{R}^2$拟合一个二元二次多项式，i.e. $d=2,m=2$，那么有$\mathbb{b(x)}=[1,x,y,x^2,xy,y^2]^T$，那么线性方程组为：
$$
\begin{bmatrix}
1 & x_i & y_i & x_i^2 & x_iy_i & y_i^2 \\
x_i & x_i^2 & x_iy_i & x_i^3 & x_i^2y_i & x_iy_i^2 \\
y_i & x_iy_i & y_i^2 & x_i^2y_i & x_iy_i^2 & y_i^3 \\
x_i^2 & x_i^3 & x_i^2y_i & x_i^4 & x_i^3y_i & x_i^2y_i^2 \\
x_iy_i & x_i^2y_i & x_iy_i^2 & x_i^3y_i & x_i^2y_i^2 & x_iy_i^3 \\
y_i^2 & x_iy_i^2 & y_i^3 & x_i^2y_i^2 & x_iy_i^3 & y_i^4
\end{bmatrix}
\begin{bmatrix}
c_1\\c_2\\c_3\\c_4\\c_5\\c_6
\end{bmatrix}=
\sum_i
\begin{bmatrix}
1\\x_i\\y_i\\x_i^2\\x_iy_i\\y_i^2
\end{bmatrix}
f_i
$$
考虑一组二维平面上的点$P_i={(1,1),(1,-1),(-1,1),(-1,-1),(0,0),(1,0),(-1,0),(0,1),(0,-1)}$有两组对应的函数值$f_i^1={1.0, -0.5, 1.0, 1.0, -1.0, 0.0, 0.0, 0.0, 0.0}$和$f_i^2={1.0, -1.0, 0.0, 0.0, 1.0, 0.0, -1.0, -1.0, 1.0}$。

图1为该数据的拟合：

fig1

图1. 二元二次多相似局部拟合：第一行为2组9个数据点，第二行是LS拟合函数。对应的参数向量$[c_1,…,c_6]^T$分别左为$[.0.834,.0.25,0.75,0.25,0.375,0.75]^T$，右为$[0.334,0.167,0.0,.0.5,0.5,0.0]^T$.

WLS Approximation(Weighted Least Squates)

在WLS方法中使用误差函数$J_{WLS}=\sum_i\theta(||\mathbb{\bar x-x_i}||)||f(\mathbb{x_i})-f_i||^2$，$\bar x \in \mathbb{R}^d$为固定点。因此需要最小化：
$$
min_{f_x\in\prod_m^d}\sum_i\theta(||\mathbb{\bar x-x_i}||)||f(\mathbb{x}_i)-f_i||^2\tag 1
$$
其中：权重$\theta(d)$由数据点$\mathbb{x_i}$到$\mathbb{\bar x}$之间的欧氏距离表示。

未知系数也与上述权重有关，则：

$$
\begin{aligned}
f_{\mathbb{\bar x}}(\mathbb{x})&=\mathbb{b(x)^Tc(\bar x)}=\mathbb{b(x)\cdot c(\bar x)},||x-\bar x||<h \\
b(x)&=[b_1(x),…,b_k(x)]^T \\
Gaussian:\theta (d)&=e^{-\frac{d^2}{h^2}} \\
Wendland:\theta (d)&=(1-\frac{d}{h})^4(4\frac{d}{h}+1) \\
c(\bar x)&=[\sum_i\theta (d_i)b(x_i)b(x_i)^T]^{-1}\sum_i\theta (d_i)b(x_i)f_i
\end{aligned}
$$

MLS Approximation(Moving least squares)

与WLS很类似，不同点：

n move this point ($\bar x$) over the entire parameter domain, where a weighted least squares fit is computed and evaluated for each point individually.

因此：
$$
\begin{aligned}
f(x)&=f_x(x),min_{f_x\in\prod_m^d}\sum_i\theta(||x-x_i||)||f_x(x_i)-f_i||^2 \\
f_{\bar x}(x)&=b(x)\cdot c(\bar x) \\
c(\bar x)&=[\sum_i\theta (d_i)b(x_i)b(x_i)^T]^{-1}\sum_i\theta (d_i)b(x_i)f_i
\end{aligned}
$$

拟合函数的建立

在拟合区域的局部子域上，拟合函数
$$
f(x)=\sum_{i=1}^m\alpha_{i}(x)p_i(x)=p^T(x)\alpha(x)\tag{1}
$$
其中：

$\alpha(x)=[\alpha_1(x),\alpha_2(x),…,\alpha_m(x)]^T$为待求系数，它是坐标x的函数

$p(x)=[p_1(x),p_2(x),…,p_m(x)]^T$称为基函数，它是一个$k$阶完备的多项式

$m$是基函数的项数

对于二维问题：

线性基：$p(x)=[1,x,y]^T，m=3$

二次基：$p(x)=[1,x,y,x^2,xy,y^2]^T，m=6$

考虑加权离散 $L_2$ 范式(向量 $\mathbf{x}=[x_1,x_2,…,x_n]$ 的 $L_2$ 范式 $L_2=(\sum_{i=1}^nx_i^2)^{1/2}$
$$
\begin{aligned}
J&=\sum_{i=1}^nw(x-x_i)[f(x)-y_i]^2 \\
&=\sum_{i=1}^nw(x-x_i)[p^T(x_i)\alpha(x)-y_i]^2
\end{aligned}\tag{3}
$$
其中：

$n$是影响区域内节点的数目

$f(x)$是拟合函数

$y_i$是$x=x_i$处的节点值，$y_i=y(x_i)$

$w(x-x_i)$是节点$x_i$的权函数

为确定系数$\alpha(x)$，$J$应取极小值。对$\alpha$求导：
$$
\frac{\partial J}{\partial \alpha}= A(x)\alpha(x)-B(x)y=0\tag{4}
$$
$$
\alpha(x)= A^{-1}(x)B(x)y\tag{5}
$$
其中：
$$
A(x)= \sum_{i=1}^nw(x-x_i)p(x_i)p^T(x_i)\tag{6} \\
$$

或者(线性二维曲线)：
$$
\begin{aligned}
A(x)&= \sum_{i=1}^np^T(x_i)w(x-x_i)p(x_i) \\
&= \begin{bmatrix}
\sum_{i=1}^nw(x-x_i) & \sum_{i=1}^nx_iw(x-x_i) \\
\sum_{i=1}^nx_iw(x-x_i) & \sum_{i=1}^nx_i^2w(x-x_i)
\end{bmatrix}
\end{aligned}
$$
$$
B(x)= [w(x-x_1)p(x_1),w(x-x_2)p(x_2),…,w(x-x_n)p(x_n)]\tag{7}
$$
$$
y^T= [y_1,y_2,…,y_n]\tag{8}
$$

把式 (5) 代入式 (1) ，就可以得到 MLS 拟合函数：
$$
f(x)=\sum_{i=1}^n\Phi_i^k(x)y_i=O^k(x)y\tag{9}
$$
其中，$O^k(x)$叫 形函数 ，$k$表示基函数的阶数：
$$
O^k(x)=[\Phi_1^k,\Phi_2^k,…,\Phi_n^k]=p^T(x)A^{-1}(x)B(x)\tag{10}
$$
如果$k=0$，则基函数$p(x)={1}$，这时的形函数为Shepard函数：
$$
O_i^{Shepard}(x)=\frac{w(x-x_i)}{\sum_{j=1}^nw(x=x_j)}\tag{11}
$$

即使基函数$p(x)$为多项式，式 (9) 中的$f(x)$也不再是多项式。

若基函数$p\in C^r$，权函数$w\in C^s$，则拟合函数$f\in C^{min(r,s)}$。

权函数

移动最小二乘中的权函数$w(x-x_i)$，其在$x$的一个子领域内不等于0，在这个子域外全为0，这个子域称为权函数的支持域(即$x$的影响区域)。

一般选择圆形作为权函数的支持域，其半径记为$s_{max}$

只有包含在支持域内的数据对$x$的取值有影响

权函数$w(x-x_i)$应该是非负的，并且随着$||x-x_i||_2$的增加单调递减。

权函数还应该具有一定的光滑性，因为拟合函数会继承权函数的连续性。

连续性：若权函数$w(x-x_i)$是$C^1$阶连续的，则拟合函数也是$C^1$阶连续的。

常用的权函数是样条函数，记$s=x-x_i，\bar s=\frac{s}{s_{max}}$则三次样条函数如式 (12) 所示：
$$
\begin{aligned}
w (\bar s) &= \frac{2}{3}-4\bar s^2 + 4\bar s^3 & (\bar s\le\frac{1}{2}) \\
w (\bar s) &= \frac{4}{3}-4\bar s+4\bar s^2-\frac{4}{3}\bar s^3 & (\frac{1}{2}<\bar s \le 1) \\
w (\bar s) &= 0 & (\bar s>1)
\end{aligned}
\tag{12}
$$

支持域应该包含足够多的节点，使得式 (6) 中$A(x)$可逆。

如果式 (1) 使用的是线性基函数，则 曲线拟合 的支持域内应该至少包含不重叠的 2个节点，曲面拟合 的支持域内应该至少包含不在同一条直线上的 3个节点

MLS拟合流程

基本思想：先讲拟合区域网格化，然后用公式 (9) 求出网格上节点值，最后连接网格节点形成拟合曲线（曲面）。

将拟合区域网格化；
对每个网格点x进行循环：
- 确定网格点x的支持域(半径)的大小；
- 确定包含在x的影响区域内的节点；
- 计算行函数$O^k(x)$；
- 计算网格点x处的函数值
结束网格点循环；
连接网格点形成拟合曲线(曲面)。

误差分析

为分析曲面拟合误差，定义误差项：
$$
error=\int_U[f^{EXACT}(x,y)-f^{MLS}(x,y)]^2dU\tag{14}
$$
其中，$f^{EXACT}(x,y)$是精确值，$f^{MLS}(x,y)$是移动最小二乘的计算值。

Applications

The MLS surface of a point-set with varying density (the
density is reduced along the vertical axis from top to bottom). The
surface is obtained by applying the projection operation described
by Alexa et. al. [2003]. Image courtesy of Marc Alexa.

Coding

References

asapmls
https://blog.csdn.net/baidu_38127162/article/details/82380914
https://blog.csdn.net/liumangmao1314/article/details/54179526
https://blog.csdn.net/hjimce/article/details/46550001
https://wenku.baidu.com/view/fe7a74976f1aff00bed51eb1.html
https://blog.csdn.net/liumangmao1314/article/details/89421806
https://en.wikipedia.org/wiki/Moving_least_squares
Final Project: Image Deformation Using Moving Least Squares

本博客所有文章除特别声明外，均采用 CC BY-SA 4.0 协议，转载请注明出处！

正则表达式上一篇

MTCNN算法与代码理解下一篇

移动最小二乘的曲线曲面拟合