函数的逼近和 DFT 组合意义

期中预习。

引入

定义函数的 \(p\)-范数 \(||f(x)||_p=\left(\int_a^b |f(x)|^p\mathrm dx\right)^{1/p}\)。

对于不同的范数，有不同的多项式逼近方法。

• 对于函数 \(P(x)\)，在某个集合 \(S\) 中寻找最小化 \(||f(x)-P(x)||_{\infty}=\max_{a\le x\le b}|f(x)-P(x)|\) 的 \(f\)：最优一致逼近
• 对于函数 \(P(x)\)，在某个集合 \(S\) 中寻找最小化 \(||f(x)-P(x)||_2\) 的 \(f\)：最优平方逼近。在离散的情况下，就是最小二乘拟合。

正交多项式

定义

对于 \(x_0,\dots,x_m\) 处的离散函数，定义带权 \(w_0,\dots,w_m\) 的内积为 \((f,g)=\sum_{i=0}^m w_if(x_i)g(x_i)\)。对于 \([a,b]\) 上的连续函数，定义带权 \(\rho(x)\ge 0\)（一般认为 \(\rho\) 只在有限个点为 0）的内积为 \(\int_a^b \rho(x)f(x)g(x)\mathrm dx\)。正交多项式就是一列两两内积为 0，自身内积（范数）大于 0 的多项式。

例：\(1,\cos x,\sin x,\cos 2x,\sin 2x,\dots\) 是 \([-\pi,\pi]\) 上带权 1 的正交函数。

构造

从 \(\{1,x,x^2,\dots\}\) 开始执行 Schmidt 正交化当然可以，但也可以递推：

\[\varphi_0(x)=1,\varphi_1(x)=x-a_0 \\ \varphi_{k+1}(x)=(x-a_k)\varphi_k(x)-b_k\varphi_{k-1}(x) \\ (a_k=\frac{(x\varphi_k,\varphi_k)}{(\varphi_k,\varphi_k)},b_k=\frac{(\varphi_k,\varphi_k)}{(\varphi_{k-1},\varphi_{k-1})}) \]

这个可以理解成一种“扰动法”：已知 \(\varphi_n(x)\) 和 \(\le n-1\) 次多项式都正交，那考虑扰动一下得到 \(x\varphi_{n}(x)\)，根据定义它会和 \(\le n-2\) 次多项式都正交，所以它可以被 \(\varphi_{n-1}(x),\varphi_{n}(x),\varphi_{n+1}(x)\) 线性表示，递推式就出来了。

正交多项式的性质（还挺有趣）：

• \(\varphi_k(x)\) 线性无关（正交必定线性无关）
• \(\varphi_k(x)\ (k\ge 1)\) 在 \([a,b]\) 内有 \(k\) 个不同的零点。
  证明：首先，\(\varphi_k(x)\) 不可能没有零点（不变号），否则内积上 \(1=\varphi_0(x)\) 不会为 0。
  假设变号 \(n\) 次，分别在 \(x_1,\dots,x_n\)，考虑内积上 \(\prod (x-x_i)\)，如果 \(n<k\)，则内积得到的函数不变号但是积分为 0，矛盾。所以 \(n=k\)。
• \(\varphi_k(x)\ (k\ge 2)\) 的零点为 \(x_1,\dots,x_k\)，则 \((a,x_1),(x_1,x_2),\dots,(x_k,b)\) 内都恰有 \(\varphi_{k+1}(x)\) 的一个零点。
  证明：不会。

常见正交多项式

• Legendre 多项式：\([-1,1]\) 上带权 \(\rho(x)=1\)，递推式 \(P_0(x)=1,P_1(x)=x,(n+1)P_{n+1}(x)=(2n+1)xP_n(x)-nP_{n-1}(x)\)。
  \[P_n(x)=\frac{1}{2^nn!}\frac{\mathrm d^n}{\mathrm dx^n}(x^2-1)^n \]
  （高代习题）
  \(P_n(-x)=(-1)^nP_n(x)\)。
  \((P_n,P_n)=2/(2n+1)\)。
• Chebyshev 多项式：\([-1,1]\) 上带权 \(\rho(x)=1/\sqrt{1-x^2}\)，递推式 \(T_0(x)=1,T_1(x)=x,T_{n+1}(x)=2xP_n(x)-P_{n-1}(x)\)。同时，\(T_n(x)=\cos (n\arccos x)\)。
  \(T_n(-x)=(-1)^nT_n(x)\)。
  \(T_n(x)\) 的零点恰好是 \(\cos \frac{2k-1}{2n}\pi,k=1,2,\dots,n\)。
  \((T_n,T_n)=\frac \pi 2,n\ge 1;\pi,n=0\)。
• Laguerre 多项式：\([0,+\infty)\) 上带权 \(\rho(x)=e^{-x}\)，递推式 \(L_0(x)=1,L_1(x)=1-x,L_{n+1}(x)=(1+2n-x)L_n(x)-n^2L_{n-1}(x)\)。
  \[L_n(x)=e^x\frac{\mathrm d^n}{\mathrm dx^n}x^ne^{-x} \]
  \((L_n,L_n)=(n!)^2\)。
• Hermite 多项式：\((-\infty,\infty)\) 上带权 \(\rho(x)=e^{-x^2}\) 的多项式。递推式 \(H_0(x)=1,H_1(x)=2x,H_{n+1}(x)=2xH_n(x)-2nH_{n-1}(x)\)。
  \[H_n(x)=(-1)^ne^{x^2}\frac{\mathrm d^n}{\mathrm dx^n}e^{-x^2} \]
  \((H_n,H_n)=2^nn!\sqrt \pi\)。

逼近函数

Taylor 多项式逼近

固定 \(x_0\)，用 \(\sum_{0\le k\le n}f^{(k)}(x_0)/k!(x-x_0)^k\) 逼近 \(f\)。需要知道导数值，估计误差要用到 \((n+1)\) 阶导数。

最佳平方逼近

Problem：有函数 \(f\)，再给定一族线性无关函数 \(\varphi_0(x),\dots,\varphi_n(x)\)，求其线性组合 \(p^*(x)\) 使得 \(||p^*(x)-f(x)||_2\) 最小。

就是求下面函数极小值：

\[I=\int_a^b \rho(x)\left(\sum_{k=0}^n a_k\varphi_k(x)-f(x)\right)^2\mathrm dx \]

取极值必要条件是 \(\partial I/\partial a_i=0,\forall 0\le i\le n\)，也即

\[\int_a^b \rho(x)2\left(\sum_{k=0}^n a_k\varphi_k(x)-f(x)\right)\varphi_i(x)\mathrm dx=0,\forall i \]

也即

\[\sum_{k=0}^n a_k(\varphi_k(x),\varphi_i(x))=(f(x),\varphi_i(x)) \]

定义 \(G_n\) 为 \((n+1)\times (n+1)\) 矩阵，\((G_n)_{ij}=(\varphi_i,\varphi_j)\)。显然 \(G\) 满秩等价于 \(\varphi\) 线性无关。上面的方程组系数矩阵为 \(G\)，所以有唯一解。

还未证明这也是充分条件，也即方程组的解的确是最小值。上式变形得

\[\left(\sum_{k=0}^na_k\varphi_k(x)-f(x),\varphi_i(x)\right)=0 \]

这意味着误差函数与 \(\varphi_i(x)\) 正交，也与 \(\varphi_i(x)\) 的线性组合正交。设最优解为 \(\varphi^*(x)\)，某个另外线性组合为 \(\varphi(x)\)，根据正交性有勾股定理 \(||f-\varphi||^2=||f-\varphi^*||^2+||\varphi^*-\varphi||^2\ge ||f-\varphi^*||^2\)。所以的确在 \(\varphi^*\) 处取到最小值。

那最优误差是多少？计算得 \(||f-\varphi^*||^2=(f,f)-(\varphi^*,f)\)。

当 \([a,b]=[0,1]\)，\(\varphi_k(x)=x^k\)，\(\rho(x)=1\) 时，上述方程组的稀疏矩阵就是 Hilbert 矩阵。Hilbert 矩阵是“病态矩阵”，当 \(f\) 微扰时，解的变化很大。

如果用一列正交函数进行最佳平方逼近，则方程组系数矩阵是对角阵，可以直接算得 \(a_i=\frac{(f,\varphi_i)}{(\varphi_i,\varphi_i)}\)。

最佳一致逼近

Problem：有函数 \(f\)，再给定一族线性无关函数 \(\varphi_0(x),\dots,\varphi_n(x)\)，求其线性组合 \(p^*(x)\) 使得 \(||p^*(x)-f(x)||_\infty\) 最小。

定理：对于 \([a,b]\) 上的连续函数 \(f\)，只考虑 \(\le n\) 次多项式，必存在最佳一致逼近多项式。

记 \(P(x)\) 为逼近多项式，\(R(x)=f(x)-P(x)\)，\(E=\max_{a\le x\le b}|R(x)|\)。若 \(R(x_0)=E\)，则 \(x_0\) 为正偏差；若 \(R(x_0)=-E\)，则 \(x_0\) 为负偏差。由连续性，至少存在一个偏差点。

Chebyshev 定理：若 \(n\) 次多项式 \(P(x)\) 是 \(f(x)\in C[a,b]\) 的最佳一致逼近多项式，则 \(P(x)\) 至少有 \(n+2\) 个正负交替的偏差点。这样的偏差点组称为 Chebyshev 交错点组。

（事实上，设 \(T_n\) 为 Chebyshev 多项式，则 \(2^{1-n}T_n(x)\) 是最佳逼近 0 的 \(n\) 次首一多项式）

推论：

• 若 \(n\) 次多项式 \(P(x)\) 是 \(f(x)\in C[a,b]\) 的最佳一致逼近多项式，则 \(P(x)\) 同时也是 \(f(x)\) 在 \([a,b]\) 上的一个拉格朗日插值多项式。
  证明：由 Chebyshev 定理，\(R(x)\) 在 \([a,b]\) 上至少变号 \(n+1\) 次，至少有 \(n+1\) 个根，在这些根处的 \(n\) 次拉格朗日插值多项式即为 \(P\)。
• 若 \(n\) 次多项式 \(P(x)\) 是 \(f(x)\in C[a,b]\) 的最佳一致逼近多项式，\(f\) 有 \(n+1\) 阶导数且 \(f\) 的 \(n+1\) 阶导不变号，则 \(|R(a)|=|R(b)|=E\)。
  证明：反证法，假设 \(a\) 不是偏差点，则在 \((a,b)\) 内至少有 \(n+1\) 个交错偏差点。在这些点处，\(R(x)\) 取极值，故 \(R'(x)=0\)。求 \(n+1\) 阶导，用 Rolle 定理，得 \(f^{(n+1)}(x)\) 在 \((a,b)\) 上存在零点，矛盾。

一般情况下最佳一致逼近多项式的计算比较困难。

当 \(n=1\) 且 \(f(x)\) 在 \((a,b)\) 上有二阶导且二阶导不变号时，有简单的算法：从上述推论知道存在三个点 \(a,x_0,b\) 满足 \(P(a)-f(a)=-(P(x_0)-f(x_0))=P(b)-f(b)\)。问题就是求出 \(x_0\)：设 \(P(x)=a_1x+a_0\)，则 \(R'(x)=f'(x)-a_1\) 在 \((a,b)\) 上单调，所以其在 \((a,b)\) 上至多有一个零点，而偏差点就是 \(R'\) 零点，所以 \(x_0\) 满足 \(f'(x_0)=a_1\)。综上列出方程：

\[\begin{cases} a_1a+a_0-f(a)=a_1b+a_0-f(b)\\ a_1a+a_0-f(a)=f(x_0)-(a_1x_0+a_0)\\ f'(x_0)=a_1 \end{cases} \]

解得

\[\begin{cases} a_1=\frac{f(b)-f(a)}{b-a}\\ a_0=\frac{f(a)+f(x_0)}2-\frac{f(b)-f(a)}{b-a}\frac{a+x_0}{2} \end{cases} \]

对于一般情况，有一种迭代算法叫 Remes 算法：任选 \(n+2\) 个点钦定为初始 Chebyshev 交错点组，将 \(P(x)\) 的 \(n+1\) 个系数和最大误差 \(E\) 作为未知数，共 \(n+2\) 个未知数 \(n+2\) 个方程，解出 \(P(x)\)，并利用某种算法计算 \(R(x)\) 的真实极值点，在这些极值点处进行下一轮迭代。

最小二乘拟合

Problem：有 \(m+1\) 个离散的点 \((x_i,y_i)\ (0\le i\le m)\)，每个点分别带权 \(w_i\)，再给定一族线性无关函数 \(\varphi_0(x),\dots,\varphi_n(x)\)，求其线性组合 \(p^*(x)\) 使得 \(\sum_{i=0}^m w_i(y_i-p^*(x_i))\) 最小。通常 \(m\) 远大于 \(n\)。

此时，仍有类似最佳一致逼近的方程组。但是现在是权值给定的离散内积，\(\varphi_i\) 无关不能推出离散内积矩阵可逆。为了使离散内积矩阵 \(G\) 可逆，需要 \(n+1\) 个长为 \(m+1\) 的向量 \(\phi_k=(\varphi_k(x_0),\varphi_k(x_1),\dots,\varphi_k(x_m))\) 线性无关。在这些向量线性无关前提下，方程组有唯一解 \(\varphi^*\)，且类似地有以下结论：

• 该唯一解的确取到最小平方误差。
• 误差与基函数正交。
• 误差的值为 \((f,f)-(f,\varphi^*)\)。

定理：设 \(\varphi_0(x),\dots,\varphi_n(x)\) 的任意非零线性组合在给定点集 \((x_0,y_0),\dots,(x_m,y_m)\) 上都至多只有 \(n\) 个零点（\(m\ge n\)），则称这组函数满足 Haar 条件。若满足 Haar 条件，则上面提到的向量组 \(\phi_k\) 线性无关，\(G\) 可逆。

证明：设 \((m+1)\times (n+1)\) 矩阵 \(H\) 是依次把 \(\phi_k\) 按列排列得到的矩阵，\(b\) 为任一 \((n+1)\times 1\) 列向量，则 \(Hb\) 的位置 \((i,1)\) 就是 \(\sum b_j\varphi_j(x)\) 这个线性组合在 \(x_i\) 处的点值。\(Hb\) 有 \(m+1\) 行，根据 Haar 条件，\(Hb=0\) 没有非零解，所以 \(H\) 列满秩，\(\phi_k\) 线性无关。记 \(D\) 为将权重 \(w_0,\dots,w_m\) 排在对角线上得到的对角阵，则根据内积定义 \(G=H'DH\)。若 \(b'Gb=0\)，则 \((Hb)'D(Hb)=0\)，\(Hb=0\)，必有 \(b=0\)，所以 \(G\) 满秩（而且还正定）。

如果 \(\varphi_i(x)\) 正交，同样可以直接算得系数 \(a_i==\frac{(f,\varphi_i)}{(\varphi_i,\varphi_i)}\)。

实际运用时，可以先对数据变形，例如 \(y'_i=\ln y_i\)，再用多项式拟合。

最佳平方三角逼近

设对于 \(a=t_0<t_1<\dots<t_k=b\)，\(f(x)\) 在 \((t_{i-1},t_i)\) 内连续，在端点处左右极限均存在，则称 \(f\) 在 \([a,b]\) 内分段连续。

对于周期为 \(2\pi\) 且在 \([-\pi,\pi]\) 上分段连续的函数，设

\[a_j=\frac 1\pi \int_{-\pi}^\pi f(x)\cos (jx)\mathrm dx\\ b_j=\frac 1\pi \int_{-\pi}^\pi f(x)\sin (jx)\mathrm dx\\ \]

定义 \(S(x)=\frac {a_0}2+\sum_{j\ge 1} a_j\cos jx+b_j\sin jx\) 为 \(f(x)\) 的 Fourier 级数。

定理：在上述条件下，

• 若 \(x\) 是 \(f(x)\) 的连续点，则 \(S(x)\) 收敛于 \(f(x)\)
• 若 \(x\) 是 \(f(x)\) 的间断点，则 \(S(x)\) 收敛于 \((f(x+0)+f(x-0))/2\)

定理：对于 \([-\pi,\pi]\) 上的连续函数 \(f(x)\)，若要求用 \(1,\sin x,\cos x,\dots,\sin nx,\cos nx\) 的线性组合逼近 \(f\)，则存在唯一的最佳平方逼近函数 \(S_n(x)=\frac {a_0}2+\sum_{1\le j\le n} a_j\cos jx+b_j\sin jx\)。换句话说，最佳平方三角逼近就是 Fourier 级数的部分和。

定理：对于 \([-\pi,\pi]\) 上的连续函数 \(f(x)\)，设

\[a_n=\frac 1{2\pi}\int_{-\pi}^{\pi} f(x)e^{-inx}\mathrm dx \]

则 \(f(x)=\sum_{n=-\infty}^{+\infty}a_ne^{inx}\)。

这里，\(e^{inx}\ (n\in \Z)\) 是 \([-\pi,\pi]\) 上的正交函数，内积 \((f,g)\) 定义为 \(\int_{-\pi}^\pi f(x)\overline {g(x)}\mathrm dx\)。

这个定理其实就是把 Fourier 级数展开换了一种形式写出来。

最小二乘三角逼近

给定复值函数 \(f(x)\) 在 \([0,2\pi]\) 的 \(n\) 等分点上的值 \(f_k=f(\frac{2\pi k}{n})\ (0\le k<n)\)，在这些点上定义 \((f,g)=\sum_{k=0}^{n-1}f_k\overline{g_k}\)，同时可以定义正交性和范数 \(||f||_2=\sqrt{(f,f)}\)。

注意到在这样的定义下，\(1,e^{ix},e^{2ix},\dots,e^{(n-1)ix}\) 是一族正交函数。记 \(\phi_j=(1,e^{ij2\pi/n},\dots,e^{ij2\pi(n-1)/n})\)，则 \((\phi_i,\phi_j)=n[i=j]\)。这就是所谓“单位根反演”。

此时套用最小二乘逼近中的结论，得到

定理：设有给定的离散点 \((x_i,y_i)\ (0\le i<n)\)，\(x_i=\frac{2\pi i}{n}\)，\(y_i\) 为复数。若要用 \(1,e^{ix},\dots,e^{i(n-1)x}\) 的线性组合逼近 \(f\)，则最小二乘逼近函数为

\[\begin{aligned} &\sum_{k=0}^{n-1}\frac{(f,e^{ikx})}{(e^{ikx},e^{ikx})}e^{ikx}\\ ={}&\frac 1n\sum_{k=0}^{n-1}e^{ikx}\sum_{j=0}^{n-1}f_je^{-ij\frac{2\pi k}n} \end{aligned} \]

可以发现，从 \(f\) 的点值到最小二乘逼近线性组合的系数的变换是线性变换，将其称为离散 Fourier 变换（DFT）。具体地，\((f_0,f_1,\dots,f_{n-1})\) 的离散 Fourier 变换 \((c_0,c_1,\dots,c_{n-1})\) 满足

\[c_k=\frac 1n\sum_{j=0}^{n-1}f_je^{-ij\frac{2\pi k}n} \]

该变换的矩阵称为 Fourier 矩阵（DFT 矩阵）。称 \(f\) 为 \(c\) 的离散 Fourier 逆变换，满足

\[f_k=\sum_{j=0}^{n-1}c_je^{ij\frac{2\pi k}n} \]

这样，我们为 DFT 找到了一个组合意义：用 \(e^{ikx}\) 去逼近等距节点时进行的系数变换。