漫步线性代数十九——快速傅里叶变换（上）

傅里叶级数是无限维度上的线性代数，它里面的向量实际上是函数$f(x)$；他们分别投影到正弦和余弦上；然后乘以傅里叶系数$a_k,b_k$。用$a_k,b_k$乘以无限的正弦和余弦序列就重新构建了函数$f(x)$，这是个经典的情况，当然也是傅里叶希望看到的，但是在实际计算中，我们用的是离散傅里叶变换，傅里叶依然存在，只不过是在有限维而已。

这是基于正交的线性代数，输入是一个数列$y_0,\ldots,y_{n-1}$而不是函数$f(x)$，输出$c_0,\ldots,c_{n-1}$的长度和输入一样，$y,c$之间是线性关系，所以肯定存在一个矩阵，而这个矩阵就是我们要介绍的傅里叶矩阵$F$，整个数字信号处理的技术依赖于它，它有许多非凡的性质。

我们获取完语音或图像或电视或声呐信号后，都需要将其数字化，这些信号用傅里叶矩阵$F$进行变换，随后还能变换回去——也就是重构，这里面的关键就是$F,F^{-1}$能够快速的进行变换：

$F^{-1}$必须简单，$F,F^{-1}$做乘法计算速度要快。

$F^{-1}$很早就有了，它的形式和$F$很像。事实上，$F$对称并且正交(不考虑因子$\sqrt{n}$)，但是它有一个缺点，也是唯一缺点：存在复数。当然，这个代价相对比较少，在我们的接受范围内。这个难度被降到了最低全是基于一个事实：$F,F^{-1}$的所有元素是某个数$w$的幂次，并且$w^n=1$。

在$4\times 4$的离散傅里叶变换中使用$w=i$(注意$i^4=1$)，整个DFT的关键点取决于$F$和它复数共轭$\bar{F}$的乘积：

$$\begin{equation} F\bar{F}=\begin{bmatrix}1&1&1&1\\1&i&i^2&i^3\\1&i^2&i^4&i^6\\1&i^3&i^6&i^9\end{bmatrix} \begin{bmatrix} 1&1&1&1\\ 1&(-i)&(-i)^2&(-i)^3\\ 1&(-i)2&(-i)^4&(-i)^6\\ 1&(-i)3&(-i)^6&(-i)^9 \end{bmatrix} =4I\tag1 \end{equation}$$

$F\bar{F}$告诉我们$F^{-1}=\bar{F}/4$，$F$的列是正交的(在$4I$中的得到零)，$n\times n$矩阵会得到$F\bar{F}=nI$，然后$F$的逆仅仅是$\bar{F}/n$。其中的奥妙在于看到复数$w=e^{2\pi i/n}$(当$n=4$时它等于$i$)。

$F$的逆非常容易获得，这使得离散变换非常重要，$F,F^{-1}$的乘法可以非常快速的完成，而不需要像往常矩阵乘法那样计算$n^2$次乘法。矩阵向量乘法$Fc,F^{-1}y$只需要$\frac{1}{2}n\log n$步即可，这就是快速傅里叶变换。

这部分我们先从$w$和它的性质开始，然后转到$F^{-1}$，最后将FFT——快速傅里叶变换。在信号处理中最大的应用就是滤波，成功的关键在于卷积运算。用矩阵术语来说就是，所有的循环矩阵(circulant matrices)都可以用$F$来对角化，所以他们可以化简为两个FFT 和一个对角矩阵。

复根

实数方程可以由复根，例如方程$x^2+1=0$的出现导致了$i$(还有$-i$)的发明，我们也将其称为一个解，并且这种情况是封闭的。如果有人问$x^2-i=0$的解，依然是有答案。复数的平方根还是复数。现在考虑一种组合$x+iy$，其中$x$表示实部，$y$表示虚部，所有$n$阶的实或复多项式有$n$个根(可能存在复数，也可能存在重复的)，这个是代数的基本定理。

我们对$x^4=1$这样的方程很感兴趣，它有四个解——1的四次根。1的两个平方根是$-1,1$，四次根就是平方根的平方根，也就是$-1,1,i,-i$，因为$i^2=-1$，所以$i^4=1$。对于八次根我们需要$i$的平方根，这会求出$w=(1+i)/\sqrt{2}$，求$w$的平方后得到$(1+2i+i^2)/2$，也就是$i$(因为$1+i^2$等于零)，由此我们得到$w^8=i^4=1$。

在傅里叶矩阵中复数$\cos\theta+i\sin\theta$是非常特别的，实部我们用$x$ 轴表示，虚部用$y$轴表示(图1)，那么数字$w$位于单位圆上；它距离原点的距离是$\cos^2\theta+\sin^2\theta=1$，它与水平线之间的夹角是$\theta$，在后面的文章我们会看到复数会以特征值的方式出现，目前为了求解$w^n=1$，我们只需要这些特殊的点，而且他们都在单位圆上。

$w$的平方比较简单：

$$w^2=(\cos\theta+i\sin\theta)^2=\cos^2\theta-\sin^2\theta+2i\sin\theta\cos\theta$$

实部$\cos^2\theta-\sin^2\theta$是$\cos2\theta$，虚部$2\sin\theta\cos\theta$是$\sin2\theta$(注意这里没有$i$，虚部是一个实数)，因此$w^2=\cos2\theta+i\sin2\theta$。$w$的平方依然在单位院上，但是角度变成原来的两倍$2\theta$，据此我们猜测$w^n$位于夹角为$n\theta$的线上，事实上的确如此。

有一个更好的方式来计算$w$的幂，正弦和余弦的组合还可以表示成如下方式：

$$\begin{equation} \cos\theta+i\sin\theta=e^{i\theta}\tag2 \end{equation}$$

当指数$i\theta$是虚部时乘法规则(像$(e^2)(e^3)=e^5$)依然成立，$w=e^{i\theta}$的幂依然在单位圆上：

$$\begin{equation} w^2=e^{i2\theta},\quad w^n=e^{in\theta},\quad \frac{1}{w}=e^{-i\theta}\tag3 \end{equation}$$

$n$次幂的角度对应的是$n\theta$，当$n=-1$时，倒数$1/w$ 的角度是$-\theta$，如果我们用$\cos\theta+i\sin\theta$乘以$\cos-\theta+i\sin-\theta$，将得到1：

$$e^{i\theta}e^{-i\theta}=(\cos\theta+i\sin\theta)(\cos\theta-i\sin\theta)=\cos^2\theta+\sin^2\theta=1$$

这是原书作者的小批注：我曾经收到来自纽约MIT监狱的一封信，有人问我欧拉方程(2)是正确的吗，因为数学中三个最关键的函数以如此优雅的方式组合在一起太令人震惊了。最好的解答就是观察指数的幂级数：

$$e^{i\theta}=1+i\theta+\frac{(i\theta)^2}{2!}+\frac{(i\theta)^3}{3!}+\cdots$$

实部$1-\theta^2/2+\cdots$是$\cos\theta$，虚部$\theta-\theta63/6+\cdots$是$\sin\theta$，由此得证。

图1

利用这个公式，我们就可以求解$w^n=1$，此时可将其变成$e^{in\theta}=1$，这样的话$n\theta$肯定围绕着单位圆回到起点。方程的解选为$\theta=2\pi/n$：1的$n$次根就是：

$$\begin{equation} w_n=e^{2\pi i/n}=\cos \frac{2\pi}{n}+i\sin\frac{2\pi}{n}\tag4 \end{equation}$$

它的$n$次幂是$e^{2\pi i}$，也就是1，当$n=8$时，它的根是$(1+i)/\sqrt{2}$：

$$w_4=\cos\frac{\pi}{2}+i\sin\frac{\pi}{2}=i\quad w_8=\cos\frac{\pi}{4}+i\sin\frac{\pi}{4}=\frac{1+i}{\sqrt{2}}$$

四次根在$\theta=90^{\circ}$的方向，也就是$\frac{1}{4}(360^{\circ})$，其他的四次根是$i^2=-1,i^3=-i,i^4=1$，八次根是$w_8^2,w_8^3,\ldots,w_8^8$，这些根等价于单位圆在区间$2\pi/n$上的点。这里再强调一下$w_8$的平方是$w_4$，这是快速傅里叶变换的本质。这些值相加为零，首先是考虑简单的$1+i-1-i=0$，然后是复杂的

$$\begin{equation} 1+w_8+w_8^2+\cdots+w_8^7=0\tag5 \end{equation}$$

一个证明方式是左边乘以$w_8$，得到$w_8+w_8^2+\cdots+w_8^8$，其中$w_8^8=1$，数值没有发生变换，从图上来说，这八个点都旋转了$45^{\circ}$，但是依然是这八个点。因为零是乘以$w_8$后唯一不会变化的数，所以他们的和必须是零。当$n$是偶数时这些根可以成对消掉(像$1+i^2=0,i+i^3=0$)，但是1三次根相加也是零。

傅里叶矩阵和它的逆

在连续情况下，傅里叶级数在整个区间上重现函数$f(x)$，它用了无限多个正弦和余弦函数(或者指数函数)。在离散情况下，只需要选择$n$个系数$c_0,\ldots,c_{n-1}$，所以只需要$n$个点的方程。假设输入序列是$y=2,4,6,8$，输出序列是$c_0,c_1,c_2,c_3$，那么根据$Fc=y$得：

$$\begin{equation} \begin{array}{ccccccccc} c_0&+&c_1&+&c_2&+&c_3&=&2\\ c_0&+&ic_1&+&i^2c_2&+&i^3c_3&=&4\\ c_0&+&i^2c_1&+&i^4c_2&+&i^6c_3&=&6\\ c_0&+&i^3c_1&+&i^6c_2&+&i^9c_3&=&8 \end{array}\tag6 \end{equation}$$

(6)中的四个方程需要找到有四项的傅里叶级数，使得可以在区间$0$到$2\pi$上将输入和四个点$x$匹配：

$$c_0+c_1e^{ix}+c_2e^{2ix}+c_3e^{3ix}= \begin{cases} 2&x=0\\ 4&x=\pi/2\\ 6&x=\pi\\ 8&x=3\pi/2 \end{cases}$$

这些就是(6)中的四个方程，在$x=2\pi$处级数回到$y_0=2$，然后继续周期性的循环。离散傅里叶级数最好写成复数形式，也就是用指数$e^{ikx}$的组合而不是$\sin kx,\cos kx$。

对于每个$n$，联系$y,c$的矩阵存在逆，用$n$个方程来表示需要有限级数$c_0+c_1e^{ix}+\cdots$($n$项目)和$y$(在$n$个点)一致，第一项在$x=0$的地方，也就是$c_0+\cdots+c_{n-1}=y_0$，其余的点都有$w$的幂，整个问题就是$Fc=y$：

$$\begin{equation} Fc=y\quad \begin{bmatrix} 1&1&1&\cdot&1\\ 1&w&w^2&\cdot&w^{n-1}\\ 1&w^2&w^4&\cdot&w^2{n-1}\\ \cdot&\cdot&\cdot&\cdot&\cdot\\ 1&w^(n-1)&w^2(n-1)&\cdot&w^{(n-1)^2} \end{bmatrix} \begin{bmatrix} c_0\\c_1\\c_2\\\cdot\\c_{n-1} \end{bmatrix}= \begin{bmatrix} y_0\\y_1\\y_2\\\cdot\\y_{n-1} \end{bmatrix}\tag7 \end{equation}$$

傅里叶矩阵$F$中每项为$F_{jk}=w^{jk}$，这样的话我们很自然的从$0$到$n-1$计算行和列，而不是用$1$到$n$。对于第一行$j=0$，第一列$k=0$，对应的元素为$w^0=1$。

为了求出$c$，我们需要知道$F$的逆。在$4\times 4$的情况下，$F^{-1}$从$1/i=-i$开始着手，这里用到了一个通用的规则：复数$w^{-1}=\bar{w}$，它位于夹角为$-2\pi/n$的直线上，其中$w$位于夹角$+2\pi/n$的直线上。

22、根据$w^{-1=1/w=\bar{w}}$，矩阵的逆为：

$$\begin{equation} F^{-1}=\frac{1}{n} \begin{bmatrix} 1&1&1&\cdot&1\\ 1&w^{-1}&w^{-2}&\cdot&w^{-(n-1)}\\ 1&w^{-2}&w^{-4}&\cdot&\cdot\\ \cdot&\cdot&\cdot&\cdot&\cdot\\ 1&w^{-(n-1)}&w^{-2(n-1)}&\cdot&w^{-(n-1)^2} \end{bmatrix}= \frac{\bar{F}}{n}\tag8 \end{equation}$$

据此

$$\begin{bmatrix} 1&1&1\\ 1&e^{2\pi i/3}&e^{4\pi i/3}\\ 1&e^{4\pi i/3}&e^{8\pi i/3} \end{bmatrix}$$

的逆为

$$F^{-1}= \frac{1}{3}\begin{bmatrix} 1&1&1\\ 1&e^{-2\pi i/3}&e^{-4\pi i/3}\\ 1&e^{-4\pi i/3}&e^{-8\pi i/3} \end{bmatrix}$$

$F$的$j$行乘以$F^{-1}$的$j$列都是$(1+1+\cdots+1)/n=1$，比较麻烦的是非对角线部分，例如$F$的$j$行乘以$F^{-1}$的$k$列为零：

$$\begin{equation} 1\cdot1+w^jw^{-k}+w^{2j}w^{-2k}+\cdots+w^{(n-1)j}w^{-(n-1)k}=0\tag9 \end{equation}$$

这个方程的关键点是看出这些项是$W=w^jw^{-k}$的幂：

$$\begin{equation} 1+W+W^2+\cdots+W^{n-1}=0\tag{10} \end{equation}$$

$W$依然是1的一个根：$W^n=w^{nj}w^{-nk}$等于$1^j1^{-k}=1$。因为$j$和$k$不相等，所以$W$和1不相等，由此知道在单位院上你还存在其他根，这些根都满足$1+W+\cdots+W^{n-1}=0$。另一个证明有下面的方程：

$$\begin{equation} 1-W^n=(1-W)(1+W+W^2+\cdots+W^{n-1})\tag{11} \end{equation}$$

因为$W^n=1$，所以左边等于零。但是$W$又不等于1，所以最后那个因式必须为零，也就是说$F$的列是正交的。