形式陈述
给定固定实矩阵 公理库 矩阵 Matrix 以有限行列集合为索引、取值于半环,并以中间指标求和定义乘法的函数。 A ∈ R m × n 、目标秩 k ≥ 1 和过采样数 p ≥ 0 ,令 ℓ = k + p ≤ min ( m , n ) 。随机化 SVD 先找一个至多 ℓ 维的输出子空间,再在该空间内求秩至多 k 的最佳近似。基本的两遍算法如下,以下等式先在精确算术中理解。
抽取 Ω ∈ R n × ℓ ,各元素相互独立 公理库 独立性 Statistical independence 从概率表理解独立性,区分两两、相互和条件独立,并用可计算反例澄清零协方差与条件均值的限度。 ,并分别服从标准正态分布 公理库 正态分布 Normal distribution · Gaussian distribution · 高斯分布 具有指数平方密度、在仿射变换与独立求和下封闭的概率分布族。 ,形成 Y = A Ω ∈ R m × ℓ 。
令 r = rank Y ,求 col ( Y ) 的正交规范基 Q ∈ R m × r ,满足 Q T Q = I r 。满列秩时可用经济型QR 分解 公理库 QR 分解 QR factorization · QR decomposition · Economy-size QR 把长方矩阵分解为正交列与上三角因子,并区分经济型表示、数值算法和秩亏边界。 ;秩亏时只保留真实列空间的基,不把补齐 QR 的额外方向计入其中。
再访问 A ,形成 B = Q T A ∈ R r × n ,对这个行数较少的矩阵做SVD 公理库 奇异值分解 Singular value decomposition · SVD 任意有限维线性映射都可在正交规范基下表示为非负对角伸缩。 :B = U ~ Σ V T 。
取 k ′ = min ( k , r ) ,保留前 k ′ 个奇异项,输出
A ^ k = Q B k ′ = ( Q U ~ : , 1 : k ′ ) Σ 1 : k ′ , 1 : k ′ V : , 1 : k ′ T . 若 r = 0 ,输出零矩阵即可。输出通常保存左因子、奇异值和右因子,不必重新形成整个 m × n 矩阵。还应记录 k , p 、随机种子、采用的秩容差与误差诊断,便于复现并解释数值结果。为了统一公式,下文用 B k 表示“至多保留 k 项”的截断,故 k ≥ r 时 B k = B 。
P = Q Q T 是到采样列空间的正交投影 公理库 正交投影 Orthogonal projection 把向量映到子空间上最近点并使误差与子空间正交的线性算子。 。算法产生的两个近似要分清:Q B = P A 的秩至多为 ℓ ;最终的 Q B k 才保证秩至多为 k 。对任意这样选取的 Q ,都有精确的 Frobenius 范数 公理库 矩阵范数与诱导算子范数 Matrix norm · Induced matrix norm · Operator norm of a matrix 用诱导范数和常用可计算矩阵范数度量线性映射的放大能力,并区分算子范数、Frobenius 范数与谱半径。 分解
‖ A − Q B k ‖ F 2 = ‖ ( I − P ) A ‖ F 2 + ∑ j > k σ j ( B ) 2 . 第一项是采样空间没有捕获的能量,第二项是进入这个空间后为满足目标秩而再丢弃的能量。这一恒等式不要求 Ω 随机;随机性用于分析第一项通常有多大。
直觉
A Ω 的每一列都是 A 各列的一个随机线性组合。若某些左奇异方向的伸缩量较大,它们往往也在这些组合中较明显。收集 k + p 个组合而非仅 k 个,为目标子空间提供额外余量;正交化随后去掉尺度和重复方向,让 Q 只表达子空间本身。
一旦固定 Q ,大矩阵中能由该空间解释的全部信息就是 B = Q T A 。SVD 在较小的坐标系统中挑选最有价值的 k 个方向,再用 Q 搬回原来的输出空间。这个截断步骤对固定空间是最优的;若采样空间本身遗漏了重要方向,后续 SVD 无法补回。
图片加载失败 随机 SVD 的两层误差 图中的数值来自下面的固定算例。它把“找空间”与“空间内截断”分开,同时保留两层误差的实际大小;箭头不表示每次随机抽样都会取得这些数值。
例子与边界
一个可完整复算的五阶例子
取标准基 e 1 , … , e 5 ,设
A = diag ( 4 , 3 , 2 , 1 , 1 2 ) , k = 2 , p = 2 , ℓ = 4 , 并指定 以下便于算术的采样矩阵:
Ω = ( e 1 + 2 e 5 , e 2 , e 3 , e 4 ) = ( 1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1 2 0 0 0 ) . 这是确定性演示,不是用于验证高斯抽样表现的试验。直接相乘得到
Y = ( 4 e 1 + e 5 , 3 e 2 , 2 e 3 , e 4 ) , Q = ( q , e 2 , e 3 , e 4 ) , q = 4 e 1 + e 5 17 . 四列两两正交,故 r = 4 。第二遍压缩得到
B = Q T A = ( 16 / 17 0 0 0 1 / ( 2 17 ) 0 3 0 0 0 0 0 2 0 0 0 0 0 1 0 ) . B 的行也两两正交,所以 B B T = diag ( 1025 / 68 , 9 , 4 , 1 ) 。其降序奇异值为 1025 / 68 , 3 , 2 , 1 ;可取 U ~ = I 4 ,对应右奇异向量依次为
v 1 = 32 e 1 + e 5 1025 , e 2 , e 3 , e 4 . 保留前两项便得到完整的秩二近似
A ^ 2 = q 16 e 1 T + 1 2 e 5 T 17 + 3 e 2 e 2 T = ( 64 / 17 0 0 0 2 / 17 0 3 0 0 0 0 0 0 0 0 0 0 0 0 0 16 / 17 0 0 0 1 / 34 ) . 采样空间的正交补由单位向量 w = ( − e 1 + 4 e 5 ) / 17 张成。因此 I − P = w w T ,并且
‖ ( I − P ) A ‖ F 2 = ‖ w T A ‖ 2 2 = ( − 4 ) 2 + 2 2 17 = 20 17 . 空间内部再丢掉奇异值 2 , 1 ,损失 2 2 + 1 2 = 5 。两项相加给出
‖ A − A ^ 2 ‖ F 2 = 20 17 + 5 = 105 17 ≈ 6.17647 . 直接从矩阵 A − A ^ 2 求平方和也得到同一结果。全空间内的最佳秩二近似则是 A 2 = diag ( 4 , 3 , 0 , 0 , 0 ) ,误差平方为 2 2 + 1 2 + ( 1 / 2 ) 2 = 21 / 4 = 5.25 。差距来自把 e 1 与 e 5 混成了 q :Q 的列空间并不包含 e 1 。而四秩近似 Q B 的误差平方仅为 20 / 17 ,不能把这个较小数字当作最终秩二误差。
采样秩不足意味着什么
对于任意给定的 Ω ,r < k 只说明样本不足以张成 k 维空间,不说明 A 的秩小于 k 。例如 A = I 3 、k = 2 、p = 0 ,指定 Ω = ( e 1 , e 1 ) ,则 r = 1 ,输出 e 1 e 1 T ,仍遗漏两个方向。准确重构至少要求 col ( A ) ⊆ col ( Q ) ,最终秩截断还要求 rank A ≤ k 。
独立高斯抽样额外保证 rank ( A Ω ) = min ( rank A , ℓ ) 几乎处处成立:在 A 的非零右奇异子空间中,投影后的高斯矩阵仍为独立标准高斯矩阵,而其适当最大阶子式是非零多项式,取零的概率为零。因此若 rank A ≤ k ,基本算法在精确算术下几乎必然重构 A 。浮点中的数值秩由容差决定,这一精确秩结论不会自动决定该容差。
推论与应用
为什么两种误差平方可以相加
对任意 C ∈ R r × n ,写成
A − Q C = ( I − P ) A + Q ( B − C ) . 左边第一项的每一列都在 col ( Q ) 的正交补中,第二项每一列都在该空间中。更明确地,其 Frobenius 内积为
tr ( A T ( I − P ) Q ( B − C ) ) = 0 , 因为 ( I − P ) Q = 0 。又由 Q T Q = I r ,有 ‖ Q ( B − C ) ‖ F = ‖ B − C ‖ F ,故
‖ A − Q C ‖ F 2 = ‖ ( I − P ) A ‖ F 2 + ‖ B − C ‖ F 2 . 第一项与 C 无关。截断 SVD 的最佳低秩逼近定理说明,秩至多为 k 的 C 中,B k 使第二项最小,值为 ∑ j > k σ j ( B ) 2 。因此既得到了前述恒等式,也证明了 Q B k 在所有列空间包含于 col ( Q ) 、秩至多为 k 的矩阵中最优:每个这样的矩阵 X 都有 X = Q C 、C = Q T X ,且 rank C = rank X 。这一限制空间内的最优性,并不声称它等于全空间最佳近似 A k 。
高斯抽样的期望保证
令全空间最优误差为 τ k = ( ∑ j > k σ j ( A ) 2 ) 1 / 2 。Halko–Martinsson–Tropp 的 Theorem 10.5 对固定实矩阵、k ≥ 2 、p ≥ 2 、k + p ≤ min ( m , n ) 及独立标准高斯 Ω 给出投影误差的期望 公理库 期望 Expectation · Expected value 实值或复值随机变量关于概率测度的 Lebesgue 积分,概括加权平均与总体质量平衡。 界
E ‖ ( I − P ) A ‖ F ≤ 1 + k p − 1 τ k . 该定理证明中的平方误差估计还给出
E ‖ ( I − P ) A ‖ F 2 ≤ ( 1 + k p − 1 ) τ k 2 . 后一个式子不能由前一个式子简单平方推出;它来自论文 Theorem 9.1 的确定性界和 Propositions 10.1–10.2 对高斯矩阵及其伪逆的二阶矩估计。这里引用这一概率分析,下面单独推导最终秩 k 输出的保证。
由于 Q T A k 的秩至多为 k ,可将它作为 B 的低秩候选。利用截断最优性及正交投影不增范数,得到
‖ B − B k ‖ F ≤ ‖ B − Q T A k ‖ F = ‖ Q T ( A − A k ) ‖ F ≤ τ k . 将此代入确定性的平方和分解,再对随机样本取期望,便得到本页推导的保守估计
E ‖ A − Q B k ‖ F 2 ≤ ( 2 + k p − 1 ) τ k 2 , E ‖ A − Q B k ‖ F ≤ 2 + k p − 1 τ k . 最后一步使用 E Z ≤ E Z 2 。这两个式子针对最终截断,与 Theorem 10.5 针对投影的常数不同。期望描述重复抽样的平均表现,不保证某一次输出满足同一个上界;需要逐次判断时,应计算或估计该次残差,并明确诊断的精度。
计算成本与使用条件
对稠密矩阵,两次矩阵乘法分别形成 A Ω 和 Q T A ,每次成本为 O ( m n ℓ ) ;正交化、小矩阵 SVD 和提升左因子合计为 O ( ( m + n ) ℓ 2 ) 。总成本为 O ( m n ℓ + ( m + n ) ℓ 2 ) ,当 ℓ ≪ min ( m , n ) 时才体现低秩计算的优势。Q 和 B 需要 O ( ( m + n ) ℓ ) 存储;显式形成输出矩阵另需 O ( m n ) 存储。
这一版本必须能第二次访问 A 。只允许单遍读取时,无法直接沿用 B = Q T A 这一步;结构化随机变换、单遍算法各有自己的误差分析。奇异值衰减慢时,可考虑交替乘 A , A T 并在中间正交化的子空间迭代,但这里的成本与定理均针对上述基本算法。
有限精度实现应使用稳定的正交化,检查 ‖ Q T Q − I ‖ ,并核对重构残差。精确算术中可由 ‖ A ‖ F 2 − ∑ j ≤ k σ j ( B ) 2 计算最终误差平方;当误差很小时,两个接近的数相减可能损失有效位数,应改为直接计算残差或采用另有误差控制的估计。小残差和重复抽样的一致性是有用诊断,严格的单次证书则还需相应数值误差界。
参考资料