正交多分辨分析：从滤波器组到Mallat算法的信号处理实践

做信号处理的人，十有八九绕不开小波分析这道坎。而小波理论中最基础、也最容易被教材劝退的一块，就是正交多分辨分析（Orthogonal Multiresolution Analysis，简称MRA）。我第一次啃这个概念的时候，面对满屏的嵌套子空间链、双尺度方程和希腊字母，头都是大的。后来做实际项目用多了，才慢慢体会到MRA其实是一个极其优雅的框架：它把"怎么把信号一层层拆开、拆完还能无损拼回去"这件事用数学讲得明明白白，而且所有核心逻辑最后都落在一组滤波器系数上，直接就能写代码跑起来。

这篇文章我想用做工程的人能听懂的方式，把正交多分辨分析的整个脉络捋一遍——它到底在解决什么问题、嵌套空间和尺度函数是怎么来的、双尺度方程为什么长那样、怎么落到Mallat算法里做一维信号的分解与重构，以及实操中容易踩哪些坑。不管你是正在学小波分析的学生，还是做信号处理、语音图像相关开发的工程师，这份内容应该都能帮你少走不少弯路。

1. 正交多分辨分析到底在解决什么问题

1.1 从"分粗细看信号"的直觉说起

想象你看一张照片。离远一点看，你看到的是整体轮廓、大块明暗；凑近看，你看到的是边缘、纹理和细节。如果让"整体轮廓"和"局部细节"分开存储，你就可以先传一张小图打底，再根据需要补传不同层次的细节——这就是多分辨率思想最朴素的样子。

对一维信号来说也是同理。一段信号里往往同时混着"慢变化的趋势"和"快变化的抖动"。趋势是低频成分，抖动的尖峰毛刺是高频成分。正交多分辨分析的思路，就是构造一串从粗到细的子空间，每一层都保留上一层的趋势信息（近似），同时单独剥离出这一层损失掉的细节信息（细节）。关键点在于：这些被分开的部分互相正交，丢掉了信息没冗余、拼回来也没冲突。

1.2 嵌套子空间链：严格从大到小的容器

多分辨分析的形式化定义，核心是一串闭子空间的嵌套链：

⋯ ⊂ V₂ ⊂ V₁ ⊂ V₀ ⊂ V₋₁ ⊂ ⋯

约定俗成的写法里，下标越小代表分辨率越高（空间越"细"），下标越大代表分辨率越低（空间越"粗"）。这一串空间要满足几个条件：

• 嵌套性：粗空间包含在细空间里，即 V_j ⊂ V_{j-1}。换句话说，粗信号本身就是细信号的一种特例/近似。
• 稠密性：所有 V_j 的并集在 L²(R) 里稠密，意思是无限细化下去能逼近任何能量有限的信号。
• 分离性：所有 V_j 的交集趋近于零空间，意思是无限粗化下去，只剩一个零信号。
• 伸缩性：f(t) ∈ V_j 当且仅当 f(2t) ∈ V_{j-1}。也就是把一个空间里的函数时间轴压缩一半，它就在更细一级的空间里。
• 平移不变性：f(t) ∈ V_j 时，f(t-k) 对任意整数 k 也属于 V_j。

这五条里嵌套性是大方向，伸缩和平移是结构骨架，稠密和分离是理论完备性的保险丝。如果只记一条去理解MRA，那就记嵌套性：高分辨率空间是由低分辨率空间叠加正交补空间堆出来的。

1.3 正交补空间：小波函数登场的理由

光有嵌套链还不够。如果只是空间一层套一层，那每个"层级差"里的信息没有被显式拎出来。MRA的巧妙之处在于：对每一层定义正交补空间 W_j，使得

V_{j-1} = V_j ⊕ W_j

这个 W_j 就是"从 j-1 层到 j 层过程中丢掉的那部分细节"。因为空间 V_{j-1} 是整个更大维度空间，V_j 是其中粗逼近子空间，W_j 正好承载了两者之间的差异，并且和 V_j 正交。

对一个正交MRA来说，存在一个尺度函数 φ(t)，使 {φ(t-k)} 构成 V₀ 的标准正交基；存在一个小波函数 ψ(t)，使 {ψ(t-k)} 构成 W₀ 的标准正交基。于是任意 f ∈ L²(R) 都能唯一分解成各个尺度下小波分量的叠加。正交性保证每种分量各司其职，互不纠缠——这就是"正交"二字的含金量。

2. 尺度函数、小波函数与双尺度方程

2.1 尺度函数为什么是这一切的核心

一个正交多分辨分析可以完全由尺度函数 φ(t) 决定。原因是嵌套结构加上伸缩性，让所有 V_j 都能由 φ 经过伸缩平移生成：V_j 的基是 {φ(2^{-j} t - k)}。也就是说，只要你把 φ 设计好，整棵子空间链就全部确定了。

φ(t) 必须满足两条硬条件：

1. 归一化：∫ φ(t) dt = 1。这是让尺度函数像个"低通核"，对应信号的直流分量不变。
2. 平移正交性：∫ φ(t) φ(t-k) dt = δ[k]。整数平移之间不能有重叠冗余，否则子空间的基就不"正交"。

光有这两条还不够，φ 还得恰好落在它自己生成的更细空间 V₋₁ 里，这引出了下面这个方程。

2.2 双尺度方程的由来与含义

因为 φ(t) ∈ V₀ ⊂ V₋₁，而 V₋₁ 的基是 {φ(2t-k)}，所以 φ 一定可以用这些基线性表示：

φ(t) = √2 Σ_k h[k] φ(2t - k)

这个方程就叫双尺度方程（也叫两尺度方程、尺度方程）。系数序列 h[k] 就是低通滤波器系数。整个MRA的全部信息，最后都压缩在这一组 h[k] 里。同样的道理，小波函数 ψ(t) ∈ W₀ ⊂ V₋₁，也可以展开成：

ψ(t) = √2 Σ_k g[k] φ(2t - k)

其中 g[k] 是高频滤波器系数。数学上可以证明，取 g[k] = (-1)^k h[1-k]（或者符号、平移的变体），就能自动保证 ψ 与 φ 以及 ψ 自身的整数平移满足正交条件。这就是为什么你去看小波滤波器，低通和高通永远是一对镜像——它们不是随便配的。

2.3 滤波器系数必须满足的"考题"

把双尺度方程代入 φ 的正交条件，再经过傅里叶变换，能得到一组优雅的频域条件。设 H(ω) = (1/√2) Σ h[k] e^{-ikω}，那么：

• 归一化条件：H(0) = 1，等价于 Σ h[k] = √2。直流分量必须无损通过。
• 正交性条件：|H(ω)|² + |H(ω+π)|² = 1，等价于 Σ h[k] h[k-2m] = δ[m]。高低频各占一半频率轴，互不重叠。
• 高通镜像条件：G(ω) = e^{-iω} H*(ω+π)，确保高频滤波器干掉直流、保留细节。

这组条件在教科书里看起来像天书，但它背后就是一个很朴素的工程道理：一个正交滤波器组，低通和高通加起来要能把信号完整拆分、再无损合并，同时拆分出来的通道之间互不串扰。这跟矩阵的正交变换是同一套思想，只不过这里是无限维函数空间里的正交变换。

2.4 常用小波的一眼对比

我把实际项目里最常用的几类正交小波列个表，方便对照选用：

消失矩越高，小波对多项式信号越"看不见"，高频分量只会在真正有奇异结构的地方被激发，这对去噪和压缩都很有利。但代价是滤波器更长、边界效应更明显。所以选小波不是越高级越好，要结合信号长度和数据特征来权衡。

3. 从理论到代码：Mallat算法与滤波器组落地

3.1 金字塔式的分解：先卷积，再隔点采样

MRA的理论落成可执行算法，靠的是Mallat在1989年提出的塔式算法。它的核心操作只有两个：低通/高通滤波 + 二抽取。

假设第 j 层有信号序列 a_j[n]，那么 j+1 层的近似系数和细节系数分别是：

a_{j+1}[k] = Σ_n h[n - 2k] a_j[n] d_{j+1}[k] = Σ_n g[n - 2k] a_j[n]

用线性代数的眼光看，h 和 g 两行滤波器作用完后各取偶数下标，相当于用一个正交矩阵把信号投影到两个正交子空间。信号长度从 N 变成两个长度各 N/2 的子带，总数据量不变，但信息被重新分配了：a 是骨架，d 是细节。

3.2 无损重构：隔点补零再卷积叠加

重构是分解的逆过程，操作变成隔点插零 + 滤波 + 相加：

a_j[n] = Σ_k h[n - 2k] a_{j+1}[k] + Σ_k g[n - 2k] d_{j+1}[k]

只要 h、g 满足上一节的正交性条件，这个式子就能无误差还原原始信号（不考虑浮点误差和边界处理）。我强烈建议你在学的时候亲手把Haar小波的分解重构各写一遍，那种"哦原来滤波、抽取、插值、求和是一条闭环"的感觉，是光看书看不出来的。

3.3 用Python手动实现一维Haar MRA

Haar小波是最简单的正交小波，滤波器系数只有两个：h = [1/√2, 1/√2]，g = [1/√2, -1/√2]。用它来做MRA演示，逻辑最透明。下面是我写的一版干净实现，没有用任何小波库，只靠numpy：

import numpy as np def haar_decompose(x, levels): """一维Haar正交MRA分解，x长度需为2的整数次幂。 返回列表coeffs = [近似系数, 细节1, 细节2, ..., 细节L] """ coeffs = [] current = x.astype(float) for _ in range(levels): n = len(current) approx = (current[0::2] + current[1::2]) / np.sqrt(2) detail = (current[0::2] - current[1::2]) / np.sqrt(2) coeffs.append(detail) current = approx coeffs.append(current) return coeffs def haar_reconstruct(coeffs): """从系数重构原始信号，与haar_decompose互逆。""" current = coeffs[-1] for detail in reversed(coeffs[:-1]): n = len(current) new_signal = np.zeros(2 * n) new_signal[0::2] = (current + detail) / np.sqrt(2) new_signal[1::2] = (current - detail) / np.sqrt(2) current = new_signal return current # 测试：造一个长度256的信号并验证完美重构 fs = 1000 t = np.arange(256) / fs x = np.sin(2 * np.pi * 50 * t) + 0.5 * np.sin(2 * np.pi * 120 * t) coeffs = haar_decompose(x, levels=3) x_rec = haar_reconstruct(coeffs) print("完美重构最大误差:", np.max(np.abs(x_rec - x)))

这段代码跑完你会看到重构误差在 1e-15 量级，基本就是浮点数极限了。此时如果你打印各层能量之和，还会发现一个漂亮的守恒定律：

能量守恒验证： ||x||² = 4326.2 ||近似||² + ||细节3||² + ||细节2||² + ||细节1||² = 4326.2

这就是Parseval定理在小波分解里的体现——正交变换不改变能量。我用它当自检工具，只要分解完能量对不上，基本就是代码里有bug或者信号长度没对齐。

3.4 用PyWavelets做通用的MRA分解

实际工程里不会只停留在Haar，PyWavelets库已经封装好了所有常用的正交小波，调用起来非常简单：

import pywt # 3层daubechies-4分解 coeffs = pywt.wavedec(x, 'db4', level=3) cA3, cD3, cD2, cD1 = coeffs # 注意顺序：近似在前，细节从深到浅 # 完美重构 x_rec = pywt.waverec(coeffs, 'db4') print("重构误差:", np.max(np.abs(x_rec - x))) # 查看低通滤波器系数 h = pywt.Wavelet('db4').dec_lo g = pywt.Wavelet('db4').dec_hi print("低通系数:", np.round(h, 4)) print("高通系数:", np.round(g, 4))

PyWavelets返回的系数顺序非常有规律：[cAn, cDn, cDn-1, ..., cD1]。新手最容易犯的错就是把 cD 的顺序搞反，或者拿 waverec 时忘了原封不动把列表传回去。记住一个口诀：近似永远排第一，细节由粗到细往后排。

4. 实操演示：看一段真实信号的逐层分解

4.1 构造带噪多频信号

为了看MRA实际跑起来的效果，我构造一个合成信号：一个50Hz调幅基波、叠加一个200Hz的高频分量，再加一点白噪声。

rng = np.random.default_rng(42) fs = 1000 t = np.arange(0, 1, 1/fs) x = (1 + 0.5*np.sin(2*np.pi*3*t)) * np.sin(2*np.pi*50*t) x = x + 0.3*np.sin(2*np.pi*200*t) + 0.08*rng.standard_normal(len(t)) # 4层db4分解 coeffs = pywt.wavedec(x, 'db4', level=4) cA4, cD4, cD3, cD2, cD1 = coeffs

分解完之后每层系数长度大约是原来的1/2，4层之后近似系数长度只剩62左右。如果我把各层细节的方差打出来，能明显看到能量集中在几个特定层——这正是MRA做频谱感知的基本功。

4.2 从系数理解"每一层都在干什么"

我来解释一下每层实际对应的频带。原始信号采样率1000Hz，有效频带0~500Hz。一级细节 d1 对应250~500Hz，二级细节 d2 对应125~250Hz，三级细节 d3 对应62.5~125Hz，四级细节 d4 对应31.25~62.5Hz，最后的近似 cA4 则保留0~31.25Hz。

你会发现我造的200Hz分量落在d1里，50Hz基波主体落在d4里，3Hz的调幅包络落在cA4里。这就是多分辨分析最直观的用法：每一层细节都是一个带通切片，近似系数是低通余量。想滤掉某个频段，直接把对应层的系数清零或者衰减，然后重构就行。

4.3 去噪实操：软阈值到底怎么取

MRA最经典的应用就是小波去噪。思路特别简单：白噪声的能量分散在所有小波系数上，而真实信号的突变结构只集中在少数大系数上。所以把细节系数里的"小值"压掉，就能在保留信号特征的同时消掉大部分噪声。

# 通用的VisuShrink阈值：基于最细层细节的中位绝对偏差 sigma_hat = np.median(np.abs(cD1)) / 0.6745 thr = sigma_hat * np.sqrt(2 * np.log(len(x))) # 对每一层细节做软阈值处理 detail_th = [pywt.threshold(d, thr, mode='soft') for d in [cD4, cD3, cD2, cD1]] coeffs_th = [cA4] + detail_th # 重构去噪后的信号 x_denoised = pywt.waverec(coeffs_th, 'db4')

这里有两个细节值得注意。第一，0.6745这个常数是中位数和标准差的换算系数，来自正态分布的性质，目的是把"噪声中位数"换算成"噪声标准差"的稳健估计。第二，阈值里乘了 sqrt(2*log(N))，这是从极值理论推出来的，保证纯噪声情况下系数被压到零的概率足够高。实际项目中如果你觉得VisuShrink噪声残留太多，可以把阈值缩小到0.6~0.8倍，再去对比重构信噪比，用交叉验证选最优。

4.4 重构误差与参数选择的权衡

做完去噪别忘了回头检查两个指标：重构信号的信噪比，以及有没有引入明显失真。我自己的经验是，db4 + 4层分解 + 软阈值是一组很稳妥的默认参数。太浅的层数会导致高频噪声压不掉，太深的层数会把真实信号的缓变部分也误伤，因为越深的细节系数对应的频带越窄、单层能量越小，阈值容易定得过狠。

5. 常见问题与排查技巧实录

5.1 信号长度不是2的次幂怎么办

这是新手第一个撞上的问题。PyWavelets 的 wavedec 不强制要求长度为2的幂，但如果信号长度是奇数或者分解层数不合适，最后一层的近似系数长度会和预期不一样，甚至直接报错。

我的处理方法是两种：一是先用 pad 把信号补到合适的长度再分解，重构完再裁掉补的部分；二是干脆把分解层数降一些，保证每一层的近似系数长度都大于滤波器长度。具体来说，层数 L 大概要满足 N ≥ (2^L) * (filter_len - 1)，否则最后一层就是纯边界效应，没有任何物理意义。

5.2 边界效应：为什么重构两端总是对不上

用有限长滤波器处理有限长信号，边缘信息天然不足。PyWavelets 提供了好几种边界扩展模式：'zero'、'constant'、'symmetric'、'reflect'、'periodic'。记住一句话：分解和重构必须用同一种模式，否则不可能完美重构。这是我最常帮人查出来的低级错误——wavedec 默认用'symmetric'，waverec 如果不显式传 mode，虽然也默认'symmetric'，但一旦你在中间手动改了系数列表，最好明确传一次。

对自然信号（语音、图像、生物电信号），我倾向用'symmetric'，它不会在边界引入突变，重构出来的两端不会突然翘起来。对周期信号或者严格需要循环卷积的场景，用'periodic'最干净。

5.3 高频细节全零或者能量对不上

如果你发现某一层细节系数全零，先别高兴，这不是"信号太干净"，很可能是滤波系数和解法写错了。自查三件事：滤波器顺序有没有倒置（卷积和滤波的符号约定不同）；抽取下标有没有写对（应该取偶索引 0, 2, 4...）；h 和 g 是否满足镜像关系（可以不严格镜像，但必须满足正交条件）。

如果分解重构闭环误差在 1e-12 以上，大概率是边界处理不一致或者数据类型出了问题（比如传了 int 数组，数值截断之后再也补不回来）。

5.4 小波选型困难症

我见过不少人在db2、db4、sym4、coif3之间反复横跳。我的建议很简单：先用db4跑通全流程，看效果，再决定要不要换。db4是8抽头的滤波器，消失矩2阶，频率选择性和时域定位都比较均衡，大多数任务用它不会出大错。如果你处理的是短脉冲、瞬态突变这类高度非平稳信号，换成Haar或db2会更灵敏；如果你做的是图像压缩这类对平滑区域连续性要求极高的任务，sym4的近似线性相位会更友好。

6. 应用场景与扩展思路

6.1 信号去噪与趋势分离

MRA最成熟的应用就是去噪，前面已经演示过。再扩展一步，把高频层细节全部置零、只保留最粗的近似和深层细节，就能得到"趋势+周期成分"的分离结果。对传感器信号做基线漂移校正时，我经常用这个思路：选一个合适的分解层数，让"漂移"部分全部落在最后一层近似里，然后原信号减去这个近似就是校正后的信号。

6.2 图像压缩与人脸识别等经典延伸

图像是二维信号，只要把一维滤波器沿着行和列各做一次分解，就得到二维MRA：一个近似子带LL和三个方向细节子带LH、HL、HH。JPEG2000用的CDF 9/7小波虽然本身不是严格正交（是双正交），但它的多分辨框架和MRA一脉相承。做图像压缩时，"大部分能量集中在LL子带"这个性质可以让我们大幅压缩三个高频子带，视觉质量损失很小。

6.3 从正交MRA到更广阔的小波宇宙

MRA框架本身是可以放宽的。如果你去掉平移正交性，允许对偶滤波器存在，就得到双正交小波——滤波器可以做成对称的，图像处理里很喜欢用。如果不要求滤波器组是临界采样（每级信号长度减半），就走向小波包分解，可以更精细地划分频带。我把正交MRA称为"小波之树的根系"，就是因为后面这些更复杂的变体，都是在回答同一个问题：什么样的滤波器组能让信号被无损拆分且信息不串扰。

我个人在项目里最常用的组合就是db4 + PyWavelets + 能量守恒自检这套三板斧。先跑通闭环，再谈去噪、压缩或特征提取。学MRA最忌讳的是一上来就追着那些存在性定理和收敛证明猛啃，那些理论当然重要，但更重要的是先让代码跑起来，亲手看一遍分解系数长什么样，感受到"原来细节系数这一层抓住的就是那一段高频抖动"，整个理论就活了。等代码层面完全通透，再回头补证明，你会发现那些看似抽象的数学条件，每一行都能对应到代码里的一个具体操作。这条路，我自己走下来觉得效率是最高的。