简介:这份资源面向具备一定MATLAB基础、希望入门在线字典学习与在线机器学习的学习者与研究者,核心解决序列数据建模中字典更新与损失评估的实现问题。压缩包共5个文件,以4个.m脚本和1个.mat数据文件为主,整体约1.11MB,脚本分别承担主流程调度、字典更新、代价函数计算与演示运行等职责,数据文件用于存放实验所需矩阵。已有692人学习下载,说明该方向具备一定关注度。读者可借此理解代价函数如何衡量预测序列与真实序列的差异并驱动参数调整,掌握字典更新所涉及的梯度下降等优化思路,同时熟悉序列数据向在线学习格式的转换方式。代码结构紧凑,适合动手运行并交互式调整参数,在实践中加深对字典学习原理与MATLAB编程的理解。
1. 在线字典学习代码:从稀疏表示到流式更新的落地路径
在线字典学习(Online Dictionary Learning)解决的是一个很具体的问题:当样本以流式方式到达、无法一次性全部载入内存时,如何增量地维护一组过完备基(字典),让每个新样本都能被这组基稀疏表示。它和批量字典学习(如 K-SVD)最大的区别在于,字典不是一次性优化出来的,而是随着样本逐批更新,每来一批数据就修正一次字典原子。这个特性让它在信号处理、图像去噪、异常检测、特征压缩等场景里非常实用——尤其是数据量大到无法全量加载、或者数据分布随时间缓慢漂移的时候。如果你正在找一份能直接跑起来的在线字典学习代码,并且希望理解每一步在做什么、参数怎么调、哪里容易翻车,那这篇内容就是围绕这条线展开的。我会从算法核心递推讲起,然后给出可复现的 Python 实现,再拆解参数和避坑点,最后落到一个具体的调优技巧上。
2. 在线字典学习的算法骨架与递推逻辑
2.1 从批量目标函数到在线递推
字典学习的核心目标函数可以写成:
min_{D, X} (1/2) * ||Y - D X||_F^2 + lambda * ||X||_1 s.t. ||d_j||_2 <= 1, j = 1, 2, ..., K其中 Y 是观测矩阵(每列一个样本),D 是字典(每列一个原子),X 是稀疏编码矩阵,lambda 控制稀疏度,K 是字典原子个数。批量方法每次用全部 Y 来优化 D 和 X,计算量和内存占用随样本数线性增长。
在线字典学习的思路来自 Mairal 等人提出的在线优化框架:把目标函数拆成对每个样本的损失之和,然后用随机近似的方法,每来一批样本就更新一次字典。具体来说,维护两个累积量:
- A:字典原子的自相关累积矩阵,尺寸 K×K
- B:字典原子与样本的互相关累积矩阵,尺寸 n×K(n 是样本维度)
每处理一批新样本 Y_t,先固定当前字典 D_{t-1},对 Y_t 做稀疏编码得到 X_t,然后更新:
A_t = A_{t-1} + X_t X_t^T B_t = B_{t-1} + Y_t X_t^T再用 A_t 和 B_t 通过块坐标下降(Block Coordinate Descent)更新字典的每一列。这个递推结构的好处是:不需要保存历史样本,只需要维护 A 和 B 两个矩阵,内存占用与样本总量无关,只与字典大小和样本维度有关。
2.2 稀疏编码步骤的可选方案
在线字典学习里,稀疏编码是每一步都要执行的内循环,它的效率直接决定整体速度。常见做法有三种:
| 方法 | 适用场景 | 特点 |
|---|---|---|
| LARS | 小规模、需要精确解 | 精度高,但每次迭代开销大 |
| OMP | 中等规模、稀疏度已知 | 贪心策略,速度快,稀疏度需预设 |
| FISTA | 大规模、需要可控收敛 | 迭代软阈值,适合 GPU 加速 |
我一般会先用 OMP 跑通流程,因为它的稀疏度参数直观(直接设非零元素个数),调试起来不容易玄学。等流程稳定后,如果数据维度很高(比如 n > 1000),再换成 FISTA 或 LARS。
2.3 字典更新的块坐标下降
字典更新这一步,目标是在固定 A 和 B 的情况下,逐列优化 D。对第 j 列 d_j,忽略与它无关的项后,子问题变成:
min_{d_j} (1/2) d_j^T A_jj d_j - b_j^T d_j s.t. ||d_j||_2 <= 1其中 A_jj 是 A 的第 j 个对角块,b_j 是 B 的第 j 列。这个子问题有解析解:先算无约束解 u_j = (A_jj)^{-1} b_j,然后如果 ||u_j||_2 > 1,就投影到单位球上;否则直接取 u_j。实际实现时,A_jj 可能接近奇异,需要加一个小的正则项(比如 1e-6)保证数值稳定。
3. 用 Python 跑通在线字典学习的最小实现
3.1 环境准备与依赖
代码只依赖 NumPy 和 SciPy,不需要额外的深度学习框架。Python 版本建议 3.8 以上,NumPy 1.20 以上。如果你用 conda,直接:
conda create -n odl python=3.9 numpy scipy conda activate odl如果你用 pip:
pip install numpy scipy不需要 GPU,CPU 就能跑。下面所有代码都在这个环境下验证过。
3.2 核心类实现
import numpy as np from scipy.linalg import solve_triangular class OnlineDictionaryLearner: def __init__(self, n_atoms, n_features, sparsity=5, lambda_reg=1e-6, batch_size=64): """ n_atoms: 字典原子个数 K n_features: 样本维度 n sparsity: OMP 稀疏度(非零元素个数) lambda_reg: 数值稳定正则项 batch_size: 每批处理的样本数 """ self.K = n_atoms self.n = n_features self.sparsity = sparsity self.lambda_reg = lambda_reg self.batch_size = batch_size # 初始化字典:随机高斯 + 列归一化 rng = np.random.RandomState(42) self.D = rng.randn(n_features, n_atoms) self.D /= np.linalg.norm(self.D, axis=0, keepdims=True) # 累积矩阵 self.A = np.zeros((n_atoms, n_atoms)) self.B = np.zeros((n_features, n_atoms)) self.t = 0 # 已处理批次数 def _omp(self, y): """对单个样本 y 做 OMP 稀疏编码""" residual = y.copy() indices = [] coefs = np.zeros(self.K) for _ in range(self.sparsity): # 找与残差最相关的原子 proj = self.D.T @ residual idx = np.argmax(np.abs(proj)) if idx in indices: break indices.append(idx) # 最小二乘求解当前支撑集上的系数 D_sub = self.D[:, indices] coefs_sub, _, _, _ = np.linalg.lstsq(D_sub, y, rcond=None) coefs = np.zeros(self.K) coefs[indices] = coefs_sub residual = y - self.D @ coefs return coefs def _update_dict(self): """块坐标下降更新字典""" for j in range(self.K): # 提取第 j 列相关的累积量 A_jj = self.A[j, j] + self.lambda_reg if A_jj < 1e-12: continue b_j = self.B[:, j] - self.D @ self.A[:, j] + self.D[:, j] * self.A[j, j] # 无约束解 u_j = b_j / A_jj # 投影到单位球 norm_u = np.linalg.norm(u_j) if norm_u > 1.0: u_j = u_j / norm_u self.D[:, j] = u_j def partial_fit(self, Y_batch): """处理一批样本 Y_batch: (n_features, n_samples)""" n_samples = Y_batch.shape[1] X_batch = np.zeros((self.K, n_samples)) for i in range(n_samples): X_batch[:, i] = self._omp(Y_batch[:, i]) # 更新累积矩阵 self.A += X_batch @ X_batch.T self.B += Y_batch @ X_batch.T self.t += 1 # 更新字典 self._update_dict() return X_batch def transform(self, Y): """对新样本做稀疏编码,返回稀疏系数矩阵""" X = np.zeros((self.K, Y.shape[1])) for i in range(Y.shape[1]): X[:, i] = self._omp(Y[:, i]) return X这段代码的核心逻辑分三块:_omp负责稀疏编码,_update_dict负责字典更新,partial_fit把两者串起来并维护累积矩阵。_omp里每次迭代都重新做最小二乘,这是 OMP 的标准做法,虽然可以增量更新 Cholesky 分解来加速,但为了代码可读性这里用了直接求解。_update_dict里的b_j计算用了一个技巧:从 B 中减去当前字典的贡献再加上第 j 列自身的贡献,避免重复计算。partial_fit每调用一次就处理一批数据,适合流式场景。
3.3 跑一个合成数据实验
# 生成合成数据:真实字典 + 稀疏系数 np.random.seed(0) n_features, n_atoms_true = 50, 100 D_true = np.random.randn(n_features, n_atoms_true) D_true /= np.linalg.norm(D_true, axis=0, keepdims=True) n_samples = 2000 X_true = np.zeros((n_atoms_true, n_samples)) for i in range(n_samples): support = np.random.choice(n_atoms_true, 3, replace=False) X_true[support, i] = np.random.randn(3) Y = D_true @ X_true + 0.01 * np.random.randn(n_features, n_samples) # 在线字典学习 learner = OnlineDictionaryLearner(n_atoms=100, n_features=50, sparsity=3) batch_size = 64 for start in range(0, n_samples, batch_size): Y_batch = Y[:, start:start+batch_size] learner.partial_fit(Y_batch) # 评估重构误差 X_est = learner.transform(Y) Y_recon = learner.D @ X_est recon_err = np.linalg.norm(Y - Y_recon, 'fro') / np.linalg.norm(Y, 'fro') print(f"相对重构误差: {recon_err:.4f}")跑下来重构误差通常在 0.05 到 0.15 之间,取决于稀疏度和噪声水平。如果误差明显偏高,先检查稀疏度是否设得太小(比如真实稀疏度是 3,你设了 1),再检查字典原子数是否足够(原子数应该大于真实稀疏度,一般取 2 到 4 倍)。
4. 参数调优与流式场景的工程化处理
4.1 字典原子数与稀疏度的匹配关系
字典原子数 K 和稀疏度 s 之间存在一个经验关系:K 应该至少是 s 的 2 到 4 倍,否则字典表达能力不足,重构误差会居高不下。但 K 也不能太大,因为 A 矩阵是 K×K,K 超过 500 后内存和计算开销增长很快。我一般会先根据数据的本征维度估计一个下界:如果数据来自某个低维子空间,K 取子空间维度的 2 倍左右就够了。实际调的时候,从 K = 2*s 开始,逐步增加到重构误差不再明显下降为止。
4.2 批大小的选择
批大小 batch_size 影响两个东西:更新频率和梯度估计的方差。批太小(比如 1),字典更新会非常频繁,但每次更新的方向噪声很大,收敛曲线会抖;批太大(比如 512),更新次数少,收敛慢,而且失去了在线学习的意义。我一般取 32 到 128 之间,具体看数据总量:如果总样本只有几千,批大小取 32 到 64;如果总样本上万,可以取 128 到 256。另外,前几批的字典还不稳定,可以适当降低更新步长(相当于学习率),但标准在线字典学习没有显式学习率,这个效果是通过累积矩阵的递推自然实现的——早期 A 和 B 的累积量小,更新幅度大;后期累积量大,更新幅度自动减小。
4.3 处理数据分布漂移
如果数据分布随时间变化(比如传感器老化、用户行为漂移),标准在线字典学习会逐渐“忘记”旧分布,但也不会完全适应新分布,因为 A 和 B 是累积的,旧数据的影响一直在。一个实用的做法是引入遗忘因子:每次更新时把 A 和 B 乘以一个略小于 1 的系数(比如 0.99),这样旧数据的影响会指数衰减。代码改动很小:
def partial_fit_with_forgetting(self, Y_batch, forgetting=0.99): n_samples = Y_batch.shape[1] X_batch = np.zeros((self.K, n_samples)) for i in range(n_samples): X_batch[:, i] = self._omp(Y_batch[:, i]) # 遗忘因子:衰减旧累积量 self.A *= forgetting self.B *= forgetting self.A += X_batch @ X_batch.T self.B += Y_batch @ X_batch.T self.t += 1 self._update_dict() return X_batchforgetting 取 0.95 到 0.999 之间。太小会导致字典只记住最近几批,太大则漂移适应慢。这个参数没有理论最优值,得根据数据漂移的速度来试。
5. 在线字典学习代码的避坑与排查记录
5.1 字典原子出现全零或 NaN
现象:训练几批后,字典 D 的某些列变成全零,或者出现 NaN,后续所有编码结果都失效。
原因:A_jj 接近零时,无约束解 u_j = b_j / A_jj 会爆炸。A_jj 是第 j 个原子被激活次数的累积量,如果某个原子在初始化后一直没被 OMP 选中,A_jj 就一直是零。
解决:在_update_dict里加一个判断,A_jj 小于阈值(比如 1e-8)时跳过该原子的更新,或者重新随机初始化该原子。另外,初始化时确保每个原子都有机会被选中——可以用数据样本的随机子集做一次批量稀疏编码来初始化字典,而不是纯随机初始化。
5.2 重构误差不下降反而上升
现象:随着批次数增加,重构误差先降后升,或者一直不降。
原因:最常见的是稀疏度设得太小,OMP 每次只能选很少的原子,字典更新时 X_batch 太稀疏,A 和 B 的累积量不足,字典更新方向不可靠。另一个原因是批大小太小,每批的梯度噪声太大,字典在最优解附近震荡。
解决:先把稀疏度调大 1 到 2,观察误差是否下降。如果还不行,把批大小翻倍。另外检查数据是否做了归一化——如果样本的范数差异很大,OMP 会偏向选择范数大的样本对应的原子,导致字典偏向。我一般会对每个样本做 L2 归一化后再送入partial_fit。
5.3 OMP 编码速度太慢
现象:每批数据的处理时间随字典原子数 K 线性增长,K 到 200 以上时慢到无法接受。
原因:_omp里每次迭代都做一次lstsq,复杂度是 O(K * s^2),s 是稀疏度。而且每次都要重新计算self.D.T @ residual,这是 O(n * K) 的矩阵向量乘。
解决:把self.D.T @ residual改成增量更新——每次选定一个原子后,残差的变化只与那个原子有关,可以增量更新投影向量。另外,lstsq可以用 Cholesky 分解增量更新,避免每次重新求解。如果这些优化还不够,换 FISTA 做稀疏编码,它可以用矩阵运算批量处理,适合用 BLAS 加速。
5.4 字典原子之间的相关性过高
现象:字典学出来的原子之间余弦相似度很高,很多原子看起来差不多,稀疏编码时选哪个都行,导致编码结果不稳定。
原因:在线字典学习没有显式的去相关约束,如果数据本身存在高度相关的模式,字典原子会趋同。另外,如果初始化时原子之间就高度相关,后续更新很难把它们拉开。
解决:初始化时用正交化方法(比如对随机矩阵做 QR 分解)生成字典。训练过程中,可以在_update_dict里加一个去相关步骤:每次更新完一列后,把它与已更新的列做正交化。但这会增加计算量,我一般只在原子相关性超过 0.9 时才做。
5.5 流式场景下内存持续增长
现象:虽然在线字典学习理论上不需要保存历史样本,但程序运行一段时间后内存占用持续上升。
原因:检查是不是在外部维护了一个历史样本列表,每次partial_fit时把新样本 append 进去。在线字典学习本身只维护 A 和 B,内存是固定的。另一个可能是_omp里创建的临时数组没有及时释放,Python 的垃圾回收有时跟不上。
解决:确保不保存历史样本,只保留 A 和 B。如果确实需要保存一部分用于验证,用固定大小的环形缓冲区。另外,在partial_fit结束后显式删除临时变量,或者用gc.collect()强制回收。
6. 用热启动和字典剪枝把在线字典学习推到实用
在线字典学习跑通之后,下一步要解决的是“怎么让它更快、更稳、更适应真实场景”。我自己的习惯是加两个机制:热启动和字典剪枝。
热启动的思路是:不要每次从零开始训练,而是先用一小批数据做几轮批量字典学习(比如 K-SVD 跑 5 到 10 轮),得到一个粗略的字典,再用这个字典初始化在线学习的 D、A、B。这样在线阶段一开始就有不错的编码质量,收敛更快。实现上,你可以用 sklearn 的MiniBatchDictionaryLearning先跑一个partial_fit,然后把它的components_拿过来当初始字典。注意 sklearn 的字典是 (n_atoms, n_features) 的,需要转置。
字典剪枝解决的是另一个问题:训练久了之后,有些原子几乎不被任何样本选中,它们对重构没有贡献,但占着 A 和 B 的空间,还拖慢 OMP。我一般每隔 100 批做一次剪枝:统计每个原子在最近若干批中被选中的次数,如果某个原子的激活次数低于阈值(比如总批数的 1%),就把它重新初始化为当前重构误差最大的样本方向。这样既保留了字典大小不变,又让“死”原子重新参与表达。
下面是一个剪枝的实现片段:
def prune_atoms(self, activation_counts, Y_recent, threshold=0.01): """ activation_counts: 每个原子在最近若干批中被选中的次数 Y_recent: 最近一批样本,用于重新初始化死原子 threshold: 激活比例阈值 """ total = activation_counts.sum() if total == 0: return dead = np.where(activation_counts / total < threshold)[0] if len(dead) == 0: return # 用重构误差最大的样本方向重新初始化死原子 X = self.transform(Y_recent) recon = self.D @ X errors = np.linalg.norm(Y_recent - recon, axis=0) worst_idx = np.argsort(errors)[-len(dead):] for i, atom_idx in enumerate(dead): new_atom = Y_recent[:, worst_idx[i]] norm = np.linalg.norm(new_atom) if norm > 1e-8: self.D[:, atom_idx] = new_atom / norm # 重置该原子的累积量 self.A[atom_idx, :] = 0 self.A[:, atom_idx] = 0 self.B[:, atom_idx] = 0这个剪枝逻辑的核心是:死原子的累积量清零后,它会在后续批次中重新参与竞争,如果新方向确实有表达价值,它会被重新激活;如果没有,它会在下一轮剪枝中再次被替换。这样字典始终保持“活跃”状态。
验证剪枝是否有效,可以监控两个指标:一是重构误差,剪枝后应该不升或略降;二是原子激活熵,剪枝后熵应该上升(说明原子使用更均匀)。如果剪枝后误差明显上升,说明阈值设得太高,把有用的原子也剪掉了,把 threshold 调小到 0.005 再试。
最后说一个我踩过的坑:剪枝时不要直接删除原子然后缩小字典,因为 A 和 B 的尺寸变了,后续所有累积量都要重新对齐,很容易出错。保持字典大小不变、只替换原子方向,是最稳妥的做法。这个习惯让我在多个流式信号处理项目里省掉了大量调试时间。希望帮到你。
本文还有配套的精品资源,点击获取