数据拟合和函数逼近

本文介绍数据拟合和函数逼近的相关知识。

与插值问题不同的是,拟合问题不严格要求插值点的函数值相等

Least Square Fitting


最小二乘拟合要求误差在最小二乘下尽量地小,即

min⁡∥y−f(x)∥2\min \|y-f(x)\|_2

Polynomial Fitting

首先考虑多项式拟合问题。问题可以描述为:给定 (x1,y1),…,(xn,yn)(x_1,y_1),\dots,(x_n,y_n),找一个多项式 p(x)p(x) 使得

  • deg⁡p≤m<n\deg p \le m <n
  • ∑i=1n∣yi−p(xi)∣2\sum\limits_{i=1}^n|y_i-p(x_i)|^2 最小化

运用插值中使用过的方法可以转化为最小二乘问题

min⁡∥[y1y2⋮yn]−[ϕ1(x1)ϕ2(x1)⋯ϕm(x1)ϕ1(x2)ϕ2(x2)⋯ϕm(x2)⋮⋮⋱⋮ϕ1(xn)ϕ2(xn)⋯ϕm(xn)][c1c2⋮cm]∥\min\left\| \begin{bmatrix} y_1\\y_2\\ \vdots \\y_n \end{bmatrix} - \begin{bmatrix} \phi_1(x_1) & \phi_2(x_1) & \cdots & \phi_m(x_1)\\ \phi_1(x_2) & \phi_2(x_2) & \cdots & \phi_m(x_2)\\ \vdots & \vdots & \ddots & \vdots\\ \phi_1(x_n) & \phi_2(x_n) & \cdots & \phi_m(x_n) \end{bmatrix} \begin{bmatrix} c_1\\c_2\\ \vdots \\c_{m} \end{bmatrix} \right\|

即可以把真实的模型看作 y=Xcy=Xc.

假如认为只有 yy 有噪声,则问题转化为 Ordinary least squares (OLS) 问题:

min⁡y+Δy=Xc∥Δy∥F\min\limits_{y+\Delta y=Xc}\|\Delta y\|_F

使用QR分解即可以求解:X=QRX=QR,则 c=R−1QTyc=R^{-1}Q^Ty

假如认为 XX 和 yy 都有噪声,则问题转化为 Total least squares (TLS) 问题:

min⁡y+Δy=(X+ΔX)c∥[ΔX,Δy]∥F\min\limits_{y+\Delta y=(X+\Delta X)c}\|[\Delta X, \Delta y]\|_F

可以使用 SVD 进行求解。y+Δy=(X+ΔX)cy+\Delta y=(X+\Delta X)c 可以转化为 ([ΔX,Δy]+[X,Y])[c−1]=0([\Delta X, \Delta y]+[X,Y])\begin{bmatrix}c\\-1\end{bmatrix}=0,于是可以看作给 A=[X,Y]A=[X,Y] 加上扰动后变成了有 0 奇异值的矩阵。设AA 的 SVD 为 A=∑i=1nσiuivi∗A=\sum\limits_{i=1}^n\sigma_iu_iv_i^* ,则取 [ΔX,Δy]=−σnunvn∗[\Delta X,\Delta y]=-\sigma_nu_nv_n^* 时,可以取到最优值。此时 [c−1]\begin{bmatrix}c\\-1\end{bmatrix} 为 vnv_n 的倍数

Nonlinear Fitting

再来考虑基函数为非线性函数的情况。要求解的问题描述为

min⁡∑i=1n∣yi−ϕ(xi;c)∣2\min\sum_{i=1}^n |y_i-\phi(x_i;c)|^2

其中,y=ϕ(x;c)y=\phi(x;c) 关于参数 cc 是非线性的。

如果将其转化为线性相关,则可以采用线性的方法进行处理。

对于一般的非线性函数,可以使用 Gauss-Newton 方法及其变种进行求解

Gauss-Newton Method

Gauss-Newton 法的迭代过程为

c(k+1)=c(k)−(J(k))†r(k)c^{(k+1)}=c^{(k)}-(J^{(k)})^\dagger r^{(k)}

其中

r(k)=[y1−ϕ(x1;c(k))y2−ϕ(x2;c(k))⋮yn−ϕ(xn;c(k))],   J(k)=drdc∣c=c(k)r^{(k)}= \begin{bmatrix} y_1-\phi(x_1;c^{(k)})\\ y_2-\phi(x_2;c^{(k)})\\ \vdots\\ y_n-\phi(x_n;c^{(k)}) \end{bmatrix},\ \ \ J^{(k)}=\left.\frac{dr}{dc}\right|_{c=c^{(k)}}

推导过程与求解非线性方程组的 Newton 方法类似,Taylor 展开并忽略高阶项即可。每步迭代都求解了一次 OLS 问题。需要注意的是,GN 方法对初值十分敏感,对于一般的非线性函数,如果不经过精心设置通常很难收敛。

Variant

Damped G-N method:

Δc=−α⋅J†r\Delta c=-\alpha\cdot J^\dagger r

其中 α\alpha 满足 Wolfe 条件

Levenberg-Marquardt method:

Δc=−(JTJ+λD)−1JTr\Delta c=-(J^TJ+\lambda D)^{-1}J^Tr

其中 D=diag(JTJ)D=\text{diag}(J^TJ)

Quasi-Newton like methods: 对 Jacobi 矩阵做近似

Least Square Approximation


本节讨论对函数的最小平方逼近。对于要逼近的函数 f(x)f(x),要找 ϕ(x)∈V\phi(x)\in \mathcal V,V\mathcal V 为有限维空间,使得

min⁡ϕ∈V∥f(x)−ϕ(x)∥2\min\limits_{\phi\in\mathcal V}\|f(x)-\phi(x)\|_2

设 V\mathcal V 的正交基为 ϕ1(x),…,ϕn(x)\phi_1(x),\dots,\phi_n(x),根据最小二乘可以得到

ϕ(x)=∑i=1n<ϕi(x),f(x)>ϕi(x)\phi(x)=\sum\limits_{i=1}^n\left<\phi_i(x),f(x)\right>\phi_i(x)

正交基满足

<ϕi(x),ϕj(x)>=0, i≠j\left<\phi_i(x),\phi_j(x)\right>=0, \ i\ne j

内积定义为

<f(x),g(x)>=∫−∞∞ρ(x)f(x)g(x)‾dx,  ρ(x)≥0\left<f(x),g(x)\right>=\int_{-\infty}^\infty \rho(x)f(x)\overline{g(x)}dx,\ \ \rho(x)\ge 0

Orthogonal polynomials

正交多项式 p0(x),p1(x),…p_0(x),p_1(x),\dots 满足

  • deg⁡pn=n\deg p_n = n
  • <pm(x),pn(x)>=δmn\left<p_m(x),p_n(x)\right>=\delta_{mn}

它具有以下几个性质

  • 三项递推性:xpn(x)=γnpn−1(x)+αnpn(x)+βnpn+1(x)xp_n(x)=\gamma_np_{n-1}(x)+\alpha_np_n(x)+\beta_np_{n+1}(x)
  • pn(x)p_n(x) 有 nn 个不同的根
  • 交错根:设 pn(x)=ξn∏i=1n(x−λi),pn−1(x)=ξn−1∏i=1n−1(x−μi)p_n(x)=\xi_n\prod\limits_{i=1}^n(x-\lambda_i),p_{n-1}(x)=\xi_{n-1}\prod\limits_{i=1}^{n-1}(x-\mu_i),则有 λ1<μ1<λ2<⋯<λn−1<μn−1<λn\lambda_1<\mu_1<\lambda_2<\cdots<\lambda_{n-1}<\mu_{n-1}<\lambda_n

使用代数的方法可以证明。思路是设 A:p(x)↦xp(x)A:p(x)\mapsto xp(x),则 AA 是线性对称映射。对 AA 作 Lanczos 分解就很容易推到以上几个性质。

Best Uniform Approximation


本节讨论最佳一致逼近问题。与最小平方逼近类似,最佳一致逼近要求残差的无穷范数最小。

给定 f∈C[a,b]f\in C[a,b],找到 p∈Rn[x]p\in \R_n[x] 使得

∥f(x)−p(x)∥∞=sup⁡x∈[a,b]∣f(x)−p(x)∣\|f(x)-p(x)\|_\infty=\sup\limits_{x\in [a,b]}|f(x)-p(x)|

最小。

Properties

Weierstrass theorem

lim⁡n→∞inf⁡p∈Rn[x]∥f(x)−p(x)∥∞=0\lim\limits_{n\to \infty}\inf\limits_{p\in \R_n[x]}\|f(x)-p(x)\|_\infty=0

定理说明任何函数都能用多项式一致逼近。

Existence

首先当 f(x)∈Rn[x]f(x)\in \R_n[x] 时,显然存在;那么考虑 f(x)∉Rn[x]f(x)\notin \R_n[x] 的情况。设 p(x)=∑i=0ncixip(x)=\sum\limits_{i=0}^nc_ix^i,令 c=[c0,c1,…,cn]c=[c_0,c_1,\dots,c_n],F(c)=∥f(x)−p(x)∥∞F(c)=\|f(x)-p(x)\|_\infty. 由于 FF 为连续函数,考虑 ∥c∥2=1\|c\|_2=1 的集合,为有界闭集,所以 F(c)F(c) 在 ∥c∥2=1\|c\|_2=1 上有界,设

0<m≤F(c)≤M0<m\le F(c)\le M

那么考虑任意 p(x)=αp0(x)p(x)=\alpha p_0(x),其中 p0p_0 的系数 cc 在单位球上,那么有

∥f(x)−p(x)∥∞≥∥p(x)∥∞−∥f(x)∥∞=∣α∣∥p0(x)∥∞−∥f(x)∥∞≥m∣α∣−∥f(x)∥∞,      when m∣α∣>2∥f(x)∥∞≥∥f(x)∥∞≥∥f(x)−p∗(x)∥∞\begin{aligned} \|f(x)-p(x)\|_\infty&\ge\|p(x)\|_\infty-\|f(x)\|_\infty\\ &=|\alpha|\|p_0(x)\|_\infty-\|f(x)\|_\infty\\ &\ge m|\alpha|-\|f(x)\|_\infty,\ \ \ \ \ \ when \ m|\alpha|>2\|f(x)\|_\infty\\ &\ge\|f(x)\|_\infty\\ &\ge\|f(x)-p_*(x)\|_\infty \end{aligned}

所以对于最佳一致逼近,∣α∣≤2∥f(x)∥∞m|\alpha|\le\frac{2\|f(x)\|_\infty}{m} 即 ∥c∥2≤2∥f(x)∥∞m\|c\|_2\le\frac{2\|f(x)\|_\infty}{m},所以 F(c)F(c) 在这个有界闭集上有最小值,对应的多项式即为最佳一致逼近。

ChebyShev alternation

首先介绍切比雪夫交错点组。满足以下条件的 {x0,x1,…,xm}\{x_0,x_1,\dots,x_m\} 为函数 gg 的切比雪夫交错点组:

  • a≤x0<x1<⋯<xm≤ba \le x_0 < x_1<\cdots<x_m\le b
  • g(xi)=(−1)i∥g(x)∥∞g(x_i)=(-1)^i\|g(x)\|_\infty 或 g(xi)=(−1)i+1∥g(x)∥∞g(x_i)=(-1)^{i+1}\|g(x)\|_\infty,对于任意的 ii 都成立

满足以下条件的 {x0,x1,…,xm}\{x_0,x_1,\dots,x_m\} 为函数 gg 的非一致交错点组:

  • a≤x0<x1<⋯<xm≤ba \le x_0 < x_1<\cdots<x_m\le b
  • g(xi)g(xi+1)g(x_i)g(x_{i+1}),对于任意的 ii 都成立
De La Vallée-Poussin

令 f(x)∈C[a.b],q∈Rn[x]f(x)\in C[a.b],q\in\R_n[x]. 如果 f−qf-q 有非一致交错点组 {x0,x1,…,xn+1}\{x_0,x_1,\dots,x_{n+1}\},那么有

min⁡p∈Rn[x]∥f(x)−p(x)∥∞≥min⁡i∣f(xi)−q(xi)∣\min_{p\in\R_n[x]}\|f(x)-p(x)\|_\infty\ge \min_i|f(x_i)-q(x_i)|

利用反证法可以证明

Chebyshev Therorem

令 f(x)∈C[a.b],q∈Rn[x]f(x)\in C[a.b],q\in\R_n[x]. 那么

∥f(x)−p∗(x)∥∞=min⁡p∈Rn[x]∥f(x)−p(x)∥∞\|f(x)-p_*(x)\|_\infty=\min_{p\in\R_n[x]}\|f(x)-p(x)\|_\infty

当且仅当 f−p∗f-p_* 有大小为 n+2n+2 的交错点组。此交错点组不唯一

充分性可以通过 De La Vallée-Poussin 定理得到。假设 f−qf-q 有大小为 n+2n+2 的交错点组,则

∥f−q∥∞≥∥f−q∗∥∞≥min⁡i∣f(xi)−q(xi)∣=∥f−q∥∞\|f-q\|_\infty\ge\|f-q_*\|_\infty\ge\min_i|f(x_i)-q(x_i)|=\|f-q\|_\infty

所以 qq 为最佳一致逼近

必要性利用反证法可以证明。假设 g=f−p∗g=f-p_* 有大小为 m+1m+1 的交错点组,m≤nm\le n. 对于每个区间 [xi,xi+1][x_i,x_{i+1}],取最右端的零点 tit_i,这样就可以得到一系列点

a≤x0<t0<x1<⋯<xm<tm<xm+1≤ba\le x_0<t_0<x_1<\cdots<x_m<t_m<x_{m+1}\le b

令 q(x)=α∏i=0m(x−ti)q(x)=\alpha\prod\limits_{i=0}^m(x-t_i). 令 gg 在每个区间 [ti−1,ti][t_{i-1},t_i] 上的次大值为 mim_i,即

mi=max⁡{∣g(x)∣:g(x)g(xi)<0,x∈[ti−1,ti]}m_i=\max\{|g(x)|:g(x)g(x_i)<0,x\in[t_{i-1},t_i]\}

令 ξi=M+mi2,M=∥g∥∞\xi_i=\frac{M+m_i}{2},M=\|g\|_\infty. 控制 α\alpha 的值能使得 ∥q∥∞≤min⁡ξi\|q\|_\infty\le \min\xi_i,那么此时有

∥f−(p∗+q)∥∞=∥g−q∥∞≤M−12min⁡ξi\|f-(p_*+q)\|_\infty=\|g-q\|_\infty\le M-\frac12\min\xi_i

于是就找到了更优的逼近,产生矛盾。

Uniqueness

有了 Chebyshev 定理,就可以来证明唯一性了

假设存在两个最佳一致逼近 p1,p2p_1,p_2,设 ∥f(x)−p1(x)∥∞=∥f(x)−p2(x)∥∞=M, p0=p1+p22\|f(x)-p_1(x)\|_\infty=\|f(x)-p_2(x)\|_\infty=M,\ p_0=\frac{p_1+p_2}{2}

那么 M≤∥f(x)−p0(x)∥∞=∥12(f(x)−p1(x))+12(f(x)−p2(x))∥∞≤M2+M2=MM\le\|f(x)-p_0(x)\|_\infty=\|\frac12(f(x)-p_1(x))+\frac12(f(x)-p_2(x))\|_\infty\le\frac M2+\frac M2=M,所以 p0p_0 也是最佳一致逼近多项式。

将 n+2n+2 个交错点代入 f−p0f-p_0 有

f(xi)−p0(xi)=12(f(xi)−p1(xi))+12(f(xi)−p2(xi))f(x_i)-p_0(x_i)=\frac12(f(x_i)-p_1(x_i))+\frac12(f(x_i)-p_2(x_i))

若 f(xi)−p0(xi)=Mf(x_i)-p_0(x_i)=M,则有 f(xi)−p1(xi)=f(xi)−p2(xi)=Mf(x_i)-p_1(x_i)=f(x_i)-p_2(x_i)=M;f(xi)−p0(xi)=−Mf(x_i)-p_0(x_i)=-M 时同理. 所以 nn 次多项式 p1−p2p_1-p_2 至少有 n+2n+2 个根,p1−p2p_1-p_2 只能为零多项式。

Remez Algorithm

Remez 算法通过迭代法对最佳一致逼近问题进行求解,基于 Chebyshev 定理

算法流程:

  • 初始化点 {x0,…,xn+1}⊂[a,b]\{x_0,\dots,x_{n+1}\}\sub [a,b],常用 Chebyshev 结点
  • 重复以下步骤
  • 求解 p(x)∈Rn[x]p(x)\in \R_n[x] 和 EE 使得 p(xi)+(−1)iE=f(xi),∀ip(x_i)+(-1)^iE=f(x_i),\forall i
  • 使用 f−pf-p 的一个或部分极值点替代 {x0,…,xn+1}⊂[a,b]\{x_0,\dots,x_{n+1}\}\sub [a,b] 中的点,替代的过程要保证取代的点和极值点符号相同
  • ∥f(x)−p(x)∥∞\|f(x)-p(x)\|_\infty 和 ∣E∣|E| 很接近,则终止得到结果

Remez 算法对初值不敏感。


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