三阶张量HOSVD实战：Tucker分解降维、秩选择与重构误差避坑指南

简介：这份资源围绕三阶张量的高阶奇异值分解（HOSVD）与Tucker分解展开，面向从事数据分析、机器学习及多模态数据处理的学习者与研究者，帮助理解张量降维、特征提取与奇异值计算的核心思路。压缩包共2个文件，包含1个m源码文件与1份pdf参考文献，整体约4.58MB，源码可直接运行验证算法流程，pdf则提供理论推导与公式依据，便于对照学习。目前已有1638人学习下载，说明该主题在张量分解方向具有一定关注度。通过源码与文献的配合，读者可以掌握HOSVD如何构造各模式正交基、如何收缩张量并分解出核心张量与因子矩阵，进而将其应用于图像压缩、推荐系统与特征提取等场景，为高维数据处理项目提供可复用的实现参考与排错思路。

1. 三阶张量 HOSVD：为什么你的第一版 Tucker 分解总是对不上

如果你手头有一批按「时间 × 传感器 × 工况」排布的三阶张量数据，想用 HOSVD 做 Tucker 分解降维，大概率会遇到一个很别扭的现象：明明按公式算完了，重构误差却比预期大，或者因子矩阵的符号、列顺序每次跑出来都不一样。这不是你代码写错了，而是 HOSVD 本身就不是唯一分解——它给的是「子空间」，不是「唯一解」。HOSVD 全称 Higher-Order SVD，本质是把矩阵奇异值分解推广到高阶：对三阶张量沿每个 mode 做展开（unfolding），分别做一次 SVD，得到三个因子矩阵，再用它们把原张量投影到核心张量。它和 Tucker 分解的关系是：HOSVD 是 Tucker 分解的一个「非迭代、直接算」的特例，核心张量不保证是 all-orthogonal 的最优解，但胜在快、可复现、适合做初始化。这篇文章面向的是已经拿到三阶张量、想把它压成低秩表示、又不想一上来就啃 HOOI 迭代的工程师。我会把 mode-n 展开、SVD 截断、核心张量计算、重构验证这条链路拆开，参数怎么设、坑在哪，都落到能直接抄的代码上。

2. 三阶张量 HOSVD 的数学骨架与 mode-n 展开

2.1 从矩阵 SVD 到张量：为什么不能直接套

矩阵的奇异值分解是 $A = U \Sigma V^T$，其中 $U$ 和 $V$ 是正交矩阵，$\Sigma$ 是对角奇异值矩阵。到了三阶张量 $\mathcal{X} \in \mathbb{R}^{I_1 \times I_2 \times I_3}$，你没法直接对一个三维数组做 SVD，因为「对角」这个概念在高阶下变成了「核心张量」——它不再是对角，而是一个小的三阶张量 $\mathcal{G} \in \mathbb{R}^{R_1 \times R_2 \times R_3}$，满足：

$$\mathcal{X} \approx \mathcal{G} \times_1 U^{(1)} \times_2 U^{(2)} \times_3 U^{(3)}$$

这里的 $\times_n$ 是 mode-n 乘积，$U^{(n)}$ 是第 n 个 mode 的因子矩阵，列向量张成该 mode 的子空间。HOSVD 的做法是：对每个 mode 的展开矩阵做一次 SVD，取前 $R_n$ 个左奇异向量组成 $U^{(n)}$，然后反算核心张量。注意，HOSVD 得到的核心张量不满足 all-orthogonality 的最优条件，所以它是 Tucker 分解的一个「近似解」，但计算代价极低。

2.2 mode-n 展开：三阶张量到底怎么摊平

mode-n 展开（unfolding）是把张量沿第 n 个维度「切片再拼接」成矩阵。以 $\mathcal{X} \in \mathbb{R}^{3 \times 4 \times 5}$ 为例：

• mode-1 展开：固定第 1 维索引，把剩下的 $(4 \times 5)$ 拉成行，得到 $3 \times 20$ 矩阵。
• mode-2 展开：固定第 2 维索引，把 $(3 \times 5)$ 拉成行，得到 $4 \times 15$ 矩阵。
• mode-3 展开：固定第 3 维索引，把 $(3 \times 4)$ 拉成行，得到 $5 \times 12$ 矩阵。

用 NumPy 实现时，np.reshape配合transpose就能完成，但顺序很容易搞反。我一般会写一个通用函数，用np.moveaxis把目标 mode 移到最前，再 reshape：

import numpy as np def unfold(tensor, mode): """ 对三阶张量做 mode-n 展开。 tensor: shape (I1, I2, I3) mode: 0, 1, 2 对应 mode-1, mode-2, mode-3 返回: shape (I_mode, 其余维度乘积) """ # 把目标 mode 移到第 0 轴 moved = np.moveaxis(tensor, mode, 0) # 第 0 轴保留，其余轴合并 return moved.reshape(moved.shape[0], -1)

逻辑说明：np.moveaxis(tensor, mode, 0)把第 mode 维换到最前面，此时数组形状变为(I_mode, ...)，再reshape成二维。参数mode从 0 开始计数，对应数学上的 mode-1 到 mode-3。这个函数是后续所有步骤的基础，写错一个轴，后面因子矩阵的维度全乱。

2.3 HOSVD 主流程：三步拿到因子矩阵和核心张量

完整流程分三步：

1. 对每个 mode 做展开，得到 $X_{(n)}$。
2. 对 $X_{(n)}$ 做 SVD，取前 $R_n$ 个左奇异向量组成 $U^{(n)}$。
3. 计算核心张量 $\mathcal{G} = \mathcal{X} \times_1 U^{(1)T} \times_2 U^{(2)T} \times_3 U^{(3)T}$。

第 3 步的 mode-n 乘积可以用np.tensordot实现。下面是一个完整的最小可运行示例：

def hosvd(tensor, ranks): """ HOSVD 分解。 tensor: shape (I1, I2, I3) ranks: 元组 (R1, R2, R3)，每个 mode 保留的秩 返回: factors 列表 [U1, U2, U3], core 张量 """ factors = [] for mode in range(3): Xn = unfold(tensor, mode) # mode-n 展开 U, S, Vt = np.linalg.svd(Xn, full_matrices=False) factors.append(U[:, :ranks[mode]]) # 截断到目标秩 # 计算核心张量：依次做 mode-n 乘积 core = tensor for mode in range(3): core = np.tensordot(core, factors[mode].T, axes=([0], [0])) # tensordot 后 mode 维跑到最后，需要移回原位 core = np.moveaxis(core, -1, 0) return factors, core

逻辑说明：np.linalg.svd返回的U是左奇异向量矩阵，列按奇异值降序排列，取前ranks[mode]列就是截断。核心张量计算时，每次tensordot会把被乘的 mode 维放到最后，所以用np.moveaxis移回第 0 轴，保证下一轮 mode 顺序正确。参数ranks直接决定压缩率，设得越小核心张量越小，但重构误差越大。

3. 参数怎么设：秩选择、截断策略与重构误差评估

3.1 秩 (R1, R2, R3) 的三种定法

HOSVD 最关键的参数就是每个 mode 保留多少秩。常见做法有三种：

• 奇异值能量占比：对每个 mode 的奇异值序列，取累计平方和达到总平方和 90%～95% 的位置。这是最稳的，适合没有先验知识的场景。
• 固定秩：如果下游任务对维度有硬约束（比如核心张量必须小于某个尺寸），直接指定。
• 肘部法：画奇异值下降曲线，找拐点。适合奇异值衰减明显的张量。

我一般先用能量占比跑一遍，看看每个 mode 需要多少秩，再根据下游任务调整。下面这段代码可以一次性输出三个 mode 的推荐秩：

def recommend_ranks(tensor, energy=0.95): """ 根据奇异值能量占比推荐每个 mode 的秩。 energy: 累计能量阈值，默认 0.95 """ ranks = [] for mode in range(3): Xn = unfold(tensor, mode) _, S, _ = np.linalg.svd(Xn, full_matrices=False) cum_energy = np.cumsum(S**2) / np.sum(S**2) r = np.searchsorted(cum_energy, energy) + 1 ranks.append(r) print(f"mode-{mode+1}: 推荐秩={r}, 累计能量={cum_energy[r-1]:.4f}") return tuple(ranks)

参数energy设 0.95 还是 0.99，取决于你对重构精度的要求。0.95 通常能压掉一半以上的维度，0.99 会保留更多细节但压缩率下降。注意，三个 mode 的秩不需要相同，强行设成一样往往不是最优。

3.2 重构误差怎么算才不骗自己

HOSVD 做完必须验证重构误差。重构公式是：

$$\hat{\mathcal{X}} = \mathcal{G} \times_1 U^{(1)} \times_2 U^{(2)} \times_3 U^{(3)}$$

误差用相对 Frobenius 范数：

def reconstruct(core, factors): """从核心张量和因子矩阵重构张量""" result = core for mode in range(3): result = np.tensordot(result, factors[mode], axes=([0], [0])) result = np.moveaxis(result, -1, 0) return result def relative_error(original, reconstructed): """相对 Frobenius 误差""" return np.linalg.norm(original - reconstructed) / np.linalg.norm(original)

这里有个容易翻车的点：tensordot的轴顺序和moveaxis必须严格对应，否则重构出来的张量形状对但数值全错。我习惯在重构后打印original.shape和reconstructed.shape做一次断言，再算误差。如果误差大于 0.2，说明秩设得太小，或者展开方向搞错了。

3.3 和 HOOI 迭代的差别：什么时候该换

HOSVD 的核心张量不满足 all-orthogonality，所以它不是给定秩下的最优 Tucker 分解。HOOI（Higher-Order Orthogonal Iteration）通过交替最小二乘迭代，能进一步降低重构误差，通常比 HOSVD 低 10%～30%。但 HOOI 需要初始化，而 HOSVD 的结果正好是 HOOI 最常用的初始化。我的习惯是：先用 HOSVD 快速拿到因子矩阵和误差基线，如果误差可接受就直接用；如果不够，再把 HOSVD 结果喂给 HOOI 迭代几轮。这样既省时间，又不会陷入随机初始化的玄学。

4. 避坑与排查：HOSVD 落地时最容易翻车的 5 个点

4.1 现象：因子矩阵符号每次跑出来不一样 → 原因：SVD 符号不确定性 → 解决：固定符号或只比较子空间

np.linalg.svd返回的奇异向量符号是不确定的，同一份数据两次运行可能得到相反的列。这不是 bug，是 SVD 的固有性质。如果下游任务对符号敏感（比如要做因子解释），可以在得到U后，对每一列强制让绝对值最大的元素为正：

def fix_sign(U): """固定因子矩阵符号：每列绝对值最大元素为正""" for j in range(U.shape[1]): idx = np.argmax(np.abs(U[:, j])) if U[idx, j] < 0: U[:, j] = -U[:, j] return U

如果下游只关心子空间（比如做投影、分类），符号不影响结果，不用管。

4.2 现象：重构误差突然变大 → 原因：mode 展开顺序和 tensordot 轴不匹配 → 解决：用断言校验形状

这是最常见的翻车点。unfold时 mode 从 0 开始，tensordot时axes=([0], [0])乘的是当前第 0 轴，但乘完之后 mode 维跑到最后，必须moveaxis移回。少移一次或者移错位置，形状可能碰巧对，但数值全乱。我的做法是在reconstruct里加一行assert result.shape == original.shape，再算误差。如果误差大于 0.5，先检查展开和乘法的轴顺序，而不是急着调秩。

4.3 现象：核心张量某个维度比原张量还大 → 原因：秩设得比原始维度大 → 解决：秩上限取 min(I_n, 其余维度乘积)

HOSVD 的秩 $R_n$ 不能超过展开矩阵的秩，而展开矩阵的秩上限是 $\min(I_n, \prod_{i \neq n} I_i)$。如果你设的 $R_n$ 超过这个值，U[:, :R_n]会取不满，NumPy 不报错但核心张量维度异常。建议在hosvd函数开头加一行校验：

for mode in range(3): max_rank = min(tensor.shape[mode], np.prod(tensor.shape) // tensor.shape[mode]) assert ranks[mode] <= max_rank, f"mode-{mode+1} 秩超限: {ranks[mode]} > {max_rank}"

4.4 现象：数据量一大内存爆掉 → 原因：full SVD 计算全量奇异向量 → 解决：用截断 SVD 或随机 SVD

np.linalg.svd默认算全量奇异向量，对于大张量，展开矩阵可能非常大。如果只需要前 $R_n$ 个奇异向量，用scipy.sparse.linalg.svds或sklearn.utils.extmath.randomized_svd可以只算前 k 个，内存和时间都省一个量级。参数k设成目标秩加 10 左右，留一点余量再截断。

4.5 现象：重构误差和理论值对不上 → 原因：忘了 HOSVD 不是最优分解 → 解决：用 HOOI 做基线对比

HOSVD 的误差天然比 HOOI 高，如果你拿 HOSVD 的误差去对标「Tucker 分解的最优误差」，会觉得怎么调都差一截。正确的做法是：先跑 HOSVD 记录误差，再用它的因子矩阵初始化 HOOI，跑 10～20 轮，看误差能降多少。如果降幅很小，说明 HOSVD 已经接近最优；如果降幅很大，说明你的秩设得偏大，HOSVD 没抓住主要结构。

5. 进阶技巧：用 HOSVD 做张量补全与秩自适应

5.1 用 HOSVD 初始化张量补全

HOSVD 除了降维，还能用来做缺失值补全。思路是：对含缺失的张量，先用均值填充，做一次 HOSVD，用重构结果替换缺失位置，再重新做 HOSVD，迭代几轮。这就是所谓的「HOSVD 补全」：

def hosvd_completion(tensor, mask, ranks, iters=20): """ tensor: 含缺失值的张量（缺失位置填 0） mask: 布尔张量，True 表示观测到，False 表示缺失 ranks: HOSVD 秩 """ filled = tensor.copy() for it in range(iters): factors, core = hosvd(filled, ranks) recon = reconstruct(core, factors) # 只在缺失位置更新 filled[~mask] = recon[~mask] err = np.linalg.norm((filled - recon)[mask]) / np.linalg.norm(filled[mask]) if it % 5 == 0: print(f"iter {it}, 观测位置误差={err:.6f}") return filled

参数iters一般 20～50 轮就收敛，ranks可以比纯降维时略大，因为补全需要更多自由度。注意，观测位置的误差才是真正要监控的指标，缺失位置的「误差」没有意义。

5.2 秩自适应：别一次定死

固定秩的 HOSVD 需要你提前知道每个 mode 保留多少。如果不知道，可以用「秩递增」策略：从较小的秩开始，每次增加 1，监控观测位置误差，直到误差不再明显下降。下面这个表格是我在几个典型三阶张量上总结的经验值：

压缩率 = 核心张量元素数 / 原张量元素数。可以看到，秩不需要设得很大就能压到 10% 以下，但误差会随压缩率上升。我的习惯是先把误差目标定在 0.15 以内，再反推秩。

5.3 一个我踩过的坑：别在归一化之前做 HOSVD

如果三个 mode 的量纲差异很大（比如一个 mode 是时间秒级，一个是温度摄氏度，一个是压力帕斯卡），直接做 HOSVD 会让量纲大的 mode 主导奇异值分布，秩分配严重偏斜。正确做法是先对每个 mode 做标准化（减均值除标准差），再做 HOSVD，重构后再反标准化。这一步不做，后面调秩怎么调都不对。我现在养成的习惯是：拿到三阶张量先检查每个 mode 的数值范围，差两个数量级以上就先标准化，再进 HOSVD 流程。希望帮到你。

本文还有配套的精品资源，点击获取