第11章 卷积稀疏编码与卷积字典学习 (CSC/CDL)

本章关键词：平移不变性、卷积稀疏编码 (CSC)、傅里叶域 ADMM、Sherman-Morrison 公式、卷积字典学习 (CDL)、算法展开 (Unrolling)

11.1 开篇：从“切块”到“全图”的范式转变

在第10章中，我们详细讨论了基于图像块（Patch-based）的稀疏模型。尽管 K-SVD 和 OMP 等算法在去噪和修复任务中表现出色，但它们在数学建模上存在一种尴尬的“不自然感”：

  1. 平移不变性（Translation Invariance）的缺失：将图像移动一个像素，Patch 的划分网格就全变了，生成的稀疏表示可能截然不同。这与自然图像的物理特性相悖。
  2. 重叠块平均（Overlapping Patch Averaging）的理论黑洞：为了抑制块效应，我们通常对重叠 Patch 的重建结果取平均。这种操作在工程上很有效，但在变分模型中很难找到对应的能量函数解释——我们到底是在优化什么？
  3. 计算与表示的冗余：相邻 Patch 共享大部分像素，独立编码导致了巨大的计算浪费和信息的重复表示。

卷积稀疏编码（Convolutional Sparse Coding, CSC） 的提出，正是为了在数学上优雅地解决这些问题。CSC 不再将图像“切碎”，而是引入了卷积算子作为构建基石。

核心直觉：如果不把图像看作是“拼图”（Tiles），而是看作是“图章”（Stamps）的叠加呢？每一个图章（滤波器）可以在图像的任意位置盖下，且可以重叠。

11.2 CSC 变分模型构建

11.2.1 符号定义与几何直观

让我们定义 CSC 的核心组件：

  • 观测图像：$y \in \mathbb{R}^{M \times N}$ （为了简化，暂只考虑单通道灰度图）。
  • 卷积字典（滤波器组）：$\{d_k\}_{k=1}^K$，其中每个 $d_k \in \mathbb{R}^{S \times S}$ 是一个小尺寸滤波器（例如 $11 \times 11$）。$K$ 是滤波器的数量。
  • 稀疏特征图（Feature Maps）：$\{z_k\}_{k=1}^K$，其中每个 $z_k \in \mathbb{R}^{M \times N}$ 是与原图尺寸相同的矩阵。

生成模型：图像是由 $K$ 个特征图与对应滤波器卷积后求和而成的： $$ y = \sum_{k=1}^K d_k * z_k + n $$ 其中 $*$ 表示二维卷积，$n$ 是噪声。

ASCII 直观图解：

      Feature Map 1 (Sparse)          Filter 1
    +----------------------+         +---+
    | . . . . 1 . . . . .  |    *    | d |
    | . . -0.5. . . . . .  |         +---+
    +----------------------+           |

              +                        |
      Feature Map 2 (Sparse)          Filter 2         Reconstructed Image
    +----------------------+         +---+           +----------------------+
    | . . . . . . . . . .  |    *    | d |     =     |  (Edges, Textures)   |
    | . . 0.8 . . . . . .  |         +---+           +----------------------+
    +----------------------+

              + ...

注意，$z_k$ 中的非零点位置指示了“哪里”出现了某种结构，而非零点的值指示了该结构的“强度”。

11.2.2 优化目标函数

CSC 的标准形式是求解以下约束优化或无约束优化问题：

$$ \min_{\{z_k\}} \underbrace{\frac{1}{2} \left| \sum_{k=1}^K d_k * z_k - y \right|_2^2}_{\text{数据保真项}} + \underbrace{\lambda \sum_{k=1}^K | z_k |_1}_{\text{稀疏正则项}} $$

这个模型看起来很简单，但在计算上极具挑战性。

  • 变量规模：如果图像是 $512 \times 512$，且有 $K=64$ 个滤波器，我们需要优化的变量总数是 $512^2 \times 64 \approx 1600 \text{万}$ 个。
  • 耦合性：卷积操作使得 $z_k$ 中的每个像素都与周围像素纠缠在一起（由 $d_k$ 尺寸决定），且不同通道 $k$ 之间通过求和 $\sum$ 紧密耦合。

直接使用 ISTA（迭代收缩阈值算法）虽然可行，但每次迭代需要计算 $K$ 次全图卷积及其转置卷积，收敛极慢。我们需要更强大的工具：傅里叶域 ADMM。

11.3 核心求解器：频域 ADMM 与 Sherman-Morrison 技巧

这是本章的技术高潮。我们将展示如何将一个看似不可解的大规模卷积问题，分解为一系列极小的独立线性方程组。

11.3.1 变量分裂与 ADMM 框架

引入辅助变量 $\{u_k\}$ 以解耦正则项与卷积数据项：

$$ \min_{\{z_k, u_k\}} \frac{1}{2} \left| \sum_{k=1}^K d_k * z_k - y \right|_2^2 + \lambda \sum_{k=1}^K | u_k |_1 \quad \text{s.t. } \quad z_k = u_k, \forall k $$

构建增广拉格朗日函数（Scaled Form）： $$ \mathcal{L}_\rho(z, u, \eta) = \frac{1}{2} \left| \sum_k d_k * z_k - y \right|_2^2 + \lambda \sum_k | u_k |_1 + \frac{\rho}{2} \sum_k | z_k - u_k + \eta_k |_2^2 $$

ADMM 的三个标准步骤是：

  1. $z$ 步（卷积最小二乘）：$\min_z \mathcal{L}_\rho(z, u^{t}, \eta^{t})$
  2. $u$ 步（稀疏投影）：$\min_u \mathcal{L}_\rho(z^{t+1}, u, \eta^{t})$
  3. $\eta$ 步（对偶更新）：$\eta^{t+1} = \eta^t + z^{t+1} - u^{t+1}$

其中 $u$ 步非常简单，就是大家熟悉的软阈值算子（Soft Thresholding）： $$ u_k = \text{prox}_{\lambda/\rho}(z_k + \eta_k) = \text{sign}(z_k+\eta_k)\max(|z_k+\eta_k| - \frac{\lambda}{\rho}, 0) $$ 真正的难点在于 $z$ 步。

11.3.2 攻克 $z$ 子问题：频域视角

$z$ 步是一个二次优化问题。令偏导为 0，得到正规方程（Normal Equations）： $$ \sum_{j=1}^K d_j^T * \left( \sum_{k=1}^K d_k * z_k - y \right) + \rho (z_j - u_j + \eta_j) = 0, \quad \forall j=1\dots K $$ 其中 $d^T$ 表示翻转后的滤波器（相关算子）。这是一个庞大的线性方程组。

利用 Parseval 定理转换到频域： 假设采用循环边界条件，卷积等价于频域乘法。设 $\hat{x}$ 为 $x$ 的二维 DFT。 $$ \left( \hat{d}_j^* \odot \sum_{k=1}^K (\hat{d}_k \odot \hat{z}_k) \right) + \rho \hat{z}_j = \hat{d}_j^* \odot \hat{y} + \rho(\hat{u}_j - \hat{\eta}_j) $$ 这里 $\odot$ 是逐元素乘法，${}^*$ 是复共轭。

关键洞察：频率解耦（Frequency Decoupling） 注意，上述方程对于图像中的每一个频率点 $(p, q)$ 都是独立的！我们不需要解一个 $MNK \times MNK$ 的大矩阵，而是对于每个频率点 $(p, q)$，求解一个 $K \times K$ 的小线性系统。

让我们固定某个频率 $(p, q)$，令 $\mathbf{z} \in \mathbb{C}^K$ 为该频率下所有 $K$ 个特征图的值组成的列向量，$\mathbf{d} \in \mathbb{C}^K$ 为该频率下 $K$ 个滤波器的值。 上述方程可以重写为矩阵形式： $$ (\mathbf{d} \mathbf{d}^H + \rho I) \mathbf{z} = \mathbf{b} $$ 其中：

  • $\mathbf{d}$ 是 $K \times 1$ 向量。
  • $\mathbf{d} \mathbf{d}^H$ 是 $K \times K$ 的秩-1 矩阵（Rank-1 matrix）。
  • $I$ 是 $K \times K$ 单位阵。
  • $\mathbf{b} = \mathbf{d}^* \hat{y} + \rho(\mathbf{u} - \mathbf{\eta})$ 是右端项（RHS）。

11.3.3 Sherman-Morrison 公式加速

直接求解上述 $K \times K$ 系统需要 $O(K^3)$ 的复杂度。由于每个像素都要算一次，且 $K$ 可能较大（如 64 或 128），这依然太慢。

观察系数矩阵 $\mathbf{A} = \rho I + \mathbf{d}\mathbf{d}^H$。这正是 Sherman-Morrison 公式 的用武之地！ $$ (\mathbf{A} + \mathbf{u}\mathbf{v}^H)^{-1} = \mathbf{A}^{-1} - \frac{\mathbf{A}^{-1}\mathbf{u}\mathbf{v}^H\mathbf{A}^{-1}}{1 + \mathbf{v}^H\mathbf{A}^{-1}\mathbf{u}} $$ 在这里，$\mathbf{A}=\rho I$，$\mathbf{u}=\mathbf{d}$，$\mathbf{v}=\mathbf{d}$。代入公式： $$ (\rho I + \mathbf{d}\mathbf{d}^H)^{-1} = \frac{1}{\rho} I - \frac{1}{\rho^2} \frac{\mathbf{d}\mathbf{d}^H}{1 + \frac{1}{\rho}\mathbf{d}^H\mathbf{d}} $$ 计算复杂度分析：

  • 原本矩阵求逆：$O(K^3)$
  • 使用 Sherman-Morrison：主要计算量在于向量内积和乘法，复杂度降为 $O(K)$！

Rule-of-Thumb：在 CSC 的实现中，Sherman-Morrison 公式是区分“玩具代码”和“高效求解器”的分水岭。它使得计算量不再随滤波器数量 $K$ 立方增长，而是线性增长。

11.4 卷积字典学习 (CDL)

如果不知道滤波器 $d_k$，我们需要从数据中学习。这就是 卷积字典学习 (CDL)。 CDL 通常是一个双凸（Bi-convex）问题，采用交替最小化 (Alternating Minimization) 策略。

11.4.1 算法流程

  1. 稀疏编码步 (Sparse Coding Step)： 固定 $d$，求解 $z$。这就是上面讨论的 CSC 问题（使用 ADMM）。

  2. 字典更新步 (Dictionary Update Step)： 固定 $z$，求解 $d$。 $$ \min_{\{d_k\}} \frac{1}{2} \left| \sum_k z_k * d_k - y \right|_2^2 \quad \text{s.t.} \quad |d_k|_2 \le 1, \quad d_k \text{ has restricted support} $$ 注意：这里交换了卷积顺序，$z_k * d_k$ 和 $d_k * z_k$ 是一样的，但把 $d_k$ 视为变量时， $z_k$ 就充当了“滤波器”的角色。

11.4.2 字典更新的细节与陷阱

  • 空间域约束：$d_k$ 必须被限制在小的空间支撑集内（例如 $11 \times 11$），其他位置必须为 0。这意味着我们不能简单地在频域求解（频域很难施加空域的紧支集约束）。
  • 归一化约束：必须约束 $|d_k|_2 \le 1$，否则 $d$ 会趋向无穷大而 $z$ 趋向 0（所谓的尺度模糊性）。
  • 求解方法：
    • 投影梯度下降 (PGD)：计算梯度 $\nabla_d$，更新后投影到球面上并 Crop 到支撑集。简单有效。
    • ADMM：如果想更严谨，也可以对 $d$ 使用 ADMM，但这通常比 PGD 慢，因为 $d$ 变量维度很小，梯度法足够快。

11.5 扩展：多通道与多维 CSC

11.5.1 彩色图像 CSC

对于 RGB 图像，输入 $y$ 有 3 个通道。

  • 滤波器：$d_k$ 变为 3D 张量，尺寸 $3 \times S \times S$。
  • 特征图：$z_k$ 依然是 2D 矩阵（单通道）。这体现了“一种模式，多种表现”的思想——一个特征图激活，对应所有的颜色通道上产生响应。
  • 模型： $$ y^{(c)} = \sum_{k=1}^K d_k^{(c)} * z_k $$ 求解算法基本不变，只是 Sherman-Morrison 中的向量内积需要累加所有颜色通道。

11.5.2 视频 CSC (3D CSC)

对于视频，$y$ 是 $T \times H \times W$。

  • 可以使用 3D 卷积核 $d_k \in \mathbb{R}^{T' \times S \times S}$。
  • 这能捕捉动态纹理和运动模式，但也带来了巨大的内存开销。

11.6 连接深度学习：算法展开 (Unrolling)

CSC 与卷积神经网络 (CNN) 有着惊人的同构性。

  1. 单层 CSC vs. 单层 CNN：

    • 如果我们运行一次 ISTA 迭代：$z^{t+1} = \text{ReLU}(z^t - \alpha D^T(D z^t - y))$。
    • 这完全等价于一个 ResNet 块或 Recurrent CNN 层。
    • 软阈值算子 $\approx$ ReLU 激活函数。
    • 卷积 $D$ 和 $D^T$ 对应 CNN 的权重。

  2. 多层 CSC (Multi-layer CSC)： Deep CSC 模型假设 $z$ 本身也是由下一层原子生成的：$z = D_2 * z_2$。这直接对应多层 CNN 结构。Papyan 等人 (2017) 证明了多层 CSC 的前向传播就是 CNN 的一次 Pass。

  3. 可解释性网络： 通过将 CSC 的 ADMM 求解步骤展开为固定深度的网络（例如 10 层），我们可以得到一个参数可学习、但结构具有明确数学含义的“深度网络”。这被称为 Deep Unrolling 或 Algorithm Unrolling。

11.7 本章小结

  1. 范式升级：CSC 克服了 Patch 方法的平移敏感性和块效应，是更自然的图像建模方式。
  2. 计算突破：直接在空域求解 CSC 是不切实际的。频域 ADMM + Sherman-Morrison 是使 CSC 具备实用价值的关键技术组合。
  3. 秩-1 结构：在一个频率点上，多个滤波器组成的系统矩阵具有秩-1 特性，这使得求逆复杂度从立方级降为线性级。
  4. 深度桥梁：CSC 是理解 CNN 内部机制（卷积+稀疏激活）的最佳理论框架之一。

11.8 练习题

基础题 (50%)

  1. 卷积尺寸计算： 设图像大小 $256 \times 256$，滤波器 $11 \times 11$，边界填充方式为 'same'。 (a) $z_k$ 的尺寸是多少？ (b) 如果不使用 Padding（'valid' 模式），$z_k$ 的尺寸是多少？这对重建图像边缘有何影响？ (Hint: CSC 通常要求 $z_k$ 与 $y$ 同尺寸以便求和。'valid' 模式会导致边缘无法重建。)

  2. Sherman-Morrison 应用： 给定矩阵 $A = I + uv^T$，其中 $u=[1, 2]^T, v=[1, 1]^T$。 (a) 直接计算 $A^{-1}$。 (b) 使用 Sherman-Morrison 公式计算 $A^{-1}$。 (Hint: 验证两者结果一致，感受秩-1 更新的计算过程。)

  3. 计算复杂度： 假设图像总像素数为 $N$，滤波器数为 $K$。 (a) 比较空域直接卷积 $\sum d_k * z_k$ 和频域计算的复杂度。 (b) 在频域 ADMM 中，哪一步不仅与 $N$ 有关，还主要受 $K$ 影响？

挑战题 (50%)

  1. 边界伪影的数学根源： 在频域使用 FFT 求解 CSC 隐含了图像是周期延拓的假设。 (a) 这种假设会导致图像四个边缘出现什么具体的视觉伪影？（尤其是当图像上下边界亮度差异大时）。 (b) 提出一种基于“Masking”或“Padding”的修正方案，并分析这是否会破坏 Sherman-Morrison 的使用条件？ (Hint: 引入空域 Mask 算子 $W$ 后，系统变成 $W^T \mathcal{F}^{-1} \dots$，不再是对角化的。可能需要改用 Conjugate Gradient (CG) 求解 $z$ 步。)

  2. 梯度推导与反向传播： 假设我们将 CSC 的求解过程展开为一个网络层 $z = f(y; d)$。为了进行端到端的判别式训练（例如分类），我们需要计算 $\frac{\partial z}{\partial y}$ 和 $\frac{\partial z}{\partial d}$。 请基于隐函数定理（Implicit Function Theorem）或 ADMM 的定点方程，推导 $\frac{\partial z}{\partial y}$ 的解析形式。 (Hint: 当 ADMM 收敛时，定点满足 KKT 条件。对 KKT 条件两边微分。)

  3. 多层 CSC 的退化： 如果我们简单地堆叠两层 CSC：$y = D_1 * z_1, z_1 = D_2 * z_2$。 证明：这在数学上等价于单层 CSC，其等效字典为 $D_{eff} = D_1 * D_2$。那么，为什么深度学习需要非线性激活函数？ (Hint: 卷积的结合律。非线性打破了这种线性坍缩。)

11.9 常见陷阱与错误 (Gotchas)

1. 内存溢出 (OOM) 的隐形杀手

  • 现象：代码在小图上跑得通，换张 4K 图直接崩溃。
  • 原因：CSC 需要存储 $K$ 个与原图同样大小的特征图 $z$ 和拉格朗日乘子 $u, \eta$。如果是 float32，一张 1080p 图片，64 个滤波器，仅 $z$ 变量就需要 $1920 \times 1080 \times 64 \times 4$ bytes $\approx 530$ MB。算上 $u, \eta$ 和中间变量，显存轻松过几 GB。
  • 对策：必须谨慎管理内存。对于超大图像，尽管 CSC 强调全图，但在工程上可能不得不进行分块（Tile）处理并处理边界融合。

2. 滤波器初始化的“赢家通吃”

  • 现象：训练出的字典里，只有一个滤波器有意义，其他都是噪声或 0；或者所有滤波器都长得一样。
  • 原因：CDL 是非凸的，极度依赖初始化。如果初始化不好，某些滤波器可能永远得不到“激活”（$z_k$ 全为 0），从而无法获得梯度更新。
  • Rule-of-Thumb：
    • 不要用全零或常数初始化。
    • 好的初始化：从图像中随机切取 Patch，减去均值，归一化，作为初始 $d_k$。或者使用 PCA/DCT 基。

3. 频域求解中的 $\rho$ 选择

  • 现象：ADMM 收敛极慢，或者残差震荡。
  • 原因：$\rho$ 选得太小，约束 $z=u$ 极难满足；$\rho$ 选得太大，步长极小。而且在频域中，不同频率分量的能量差异巨大（低频能量极高）。
  • 技巧：使用动态 $\rho$ 调整策略（根据原始残差和对偶残差的比例自动调整）。或者对图像进行预处理（如减均值、归一化），使能量分布更可控。

4. 边界处理的偷懒代价

  • 陷阱：直接对图像做 padding 0，然后做 FFT。
  • 后果：图像边缘会出现强烈的“振铃效应”（Ringing），因为 padding 0 制造了人为的高频阶跃。
  • 对策：使用 edge (replicate) padding 或 tapered padding（使边缘平滑衰减到 0）来减少高频干扰。