UCB-CS270-组合算法笔记-全-

UCB CS270 组合算法笔记（全）

1：梯度下降法入门 🧭

在本节课中，我们将学习优化算法中最基础且应用最广泛的方法之一：梯度下降法。我们将从理解凸函数的基本概念开始，逐步介绍梯度下降的原理、算法步骤及其在不同场景下的变体。

课程概述与期望 🗺️

首先，我想告诉大家，你们中的大多数人都已经学习过本科阶段的算法课程。本科算法课程（例如伯克利的CS 170）就像在一个规划完善的国家公园中游览。那里有铺设好的高速公路和清晰的指示牌，你可以轻松驾车游览主要景点并拍照留念。

而研究生阶段的算法课程（CS 270）则旨在带领大家接近已知领域的边缘，为未来的探索做准备。这更像是一次有向导的徒步旅行。我们将离开平坦的高速公路，深入那些尚未被完全探索的山脉。旅途可能不会像本科课程那样直接展示动态规划等核心“景点”，有时需要攀爬和跋涉。但当你登上山顶，俯瞰全景时，之前学过的知识会以一种全新的方式融入你的理解。即使你不喜欢某一段徒步路线，这个过程本身也能锻炼你的“徒步”能力，让你在面对机器学习或其他领域的新“地形”时，具备更强的探索和分析能力。

因此，请以学习知识本身为重，享受探索的过程，不必过分焦虑成绩。

课程主题 📚

我为本课程选择了一系列主题，它们代表了我个人认为对未来探索有价值的工具和领域。请注意，这只是庞大理论版图中的一条路径。

1. 连续优化：这是本科算法课程中通常缺失的重要工具。我们将学习如何求解线性规划和半定规划。
2. 线性代数工具：例如奇异值分解（SVD）、主成分分析（PCA）以及一些张量方法。
3. 嵌入与映射技术：学习如何将一个问题从一个数学空间映射到另一个空间，从而简化问题。最重要的例子是降维。
4. 谱图论：学习如何利用SVD和特征向量等线性代数工具来分析图结构，例如图划分和电网络流。
5. 半定规划：这是一种强大的算法技术。掌握了连续优化（好比汇编语言）后，半定规划就像Python这样的高级语言，能更高效地解决问题，例如最大割问题。
6. 差分隐私：这是一个新兴领域，研究如何在保护隐私的前提下进行数据分析，催生了许多有趣的算法问题。
7. 其他主题：我们还将留出时间探讨其他有趣的话题。

课程评估方式 📝

以下是本课程的评估方式：

• 双周作业：每两周提交一次。由于评分资源有限，我们将采用同伴匿名互评的方式。
• 课堂笔记整理：每位同学需要负责将一堂课的笔记用LaTeX排版（我们会提供模板），并在课后两周内提交。这对所有同学都是宝贵的参考资料。
• 期末开卷考试：计划在四月的第一周进行，为期一周。
• 课程项目：项目类似于一次“侧向徒步”，要求你选择一个研究方向或开放性问题进行调研，撰写报告，探索未知领域的边缘。

强烈建议：即使课程有录像和笔记，在听课时亲手记笔记也能帮助你更好地集中注意力、加深理解。

---

梯度下降法：从直觉到算法 🧗

现在，让我们开始今天的主要内容。梯度下降法是解决最小化问题最简单、最常用的算法之一。

问题定义

我们有一个函数 F，目标是找到使 F 值最小的点。梯度下降法的思想非常直观：从一个初始点开始，计算该点的梯度（函数上升最快的方向），然后朝梯度的反方向（即函数下降最快的方向）移动一步。

凸函数：理想的“地形” ⛰️

梯度下降法可以用于任何函数，但对于凸函数，我们可以给出漂亮的理论保证。首先，让我们定义凸性。

一个集合 K 是凸集，如果对于集合内任意两点 x 和 y，连接它们的线段上的所有点也都在集合内。用公式表示，即对于任意 λ ∈ [0, 1]，有：
λx + (1-λ)y ∈ K

一个函数 F 是凸函数，如果其图像上方的区域（称为上镜图）是一个凸集。这等价于以下更常用的定义（零阶条件）：对于任意 x, y 和 λ ∈ [0, 1]，有：
F(λx + (1-λ)y) ≤ λF(x) + (1-λ)F(y)

这意味着连接函数图像上任意两点的线段总是位于函数图像的上方。

如果函数 F 是可微的，我们还有等价的一阶条件：对于任意 x, y，有：
F(y) ≥ F(x) + ∇F(x)ᵀ(y - x)

这里，∇F(x) 是 F 在点 x 处的梯度（由所有偏导数组成的向量）。这个不等式说明，凸函数在其任意一点的切线（线性近似）之下。

如果函数二阶可微，还有二阶条件：其海森矩阵（二阶导数矩阵）在任何点都是半正定的。这意味着函数在任何点附近都像一个“碗”的形状，向上弯曲。

直观理解：海森矩阵描述了函数的局部曲率。所有二次型都可以通过拉伸和旋转，变成标准形式 ∑ w_i x_i²。半正定性保证了所有拉伸系数 w_i 都是非负的，从而函数像碗一样开口向上。

基本算法：非约束优化

在非约束优化中，我们最小化定义在整个 Rⁿ 空间上的凸函数 F。

算法步骤如下：

1. 选择一个初始点 x₀。
2. 对于每一次迭代 t = 0, 1, 2, ...：
   • 计算当前点 x_t 的梯度 ∇F(x_t)。
   • 更新：x_{t+1} = x_t - η · ∇F(x_t)
     • -∇F(x_t) 是下降方向。
     • η 是步长，是一个需要选择的参数。
3. 输出：经过 T 步迭代后，输出所有迭代点的平均值（而非最后一个点）。这样分析更简单，且在实践中常能避免在最优解附近振荡。

算法背后的直觉 🔍

关于梯度下降，有两个关键的几何直觉：

1. 函数值下降：如果步长 η 足够小，沿着负梯度方向移动应该会使函数值减小。但这依赖于函数的局部性质，对于非光滑或陡峭的函数可能不成立。
2. 靠近最优解：对于凸函数，负梯度方向总是指向全局最优解 x* 所在的大致方向。更准确地说，向量 -∇F(x) 与向量 (x - x)* 之间的夹角是锐角。这意味着沿着负梯度方向移动，即使函数值不一定立刻下降很多，你到最优解的距离也可能会减小。

这个性质可以从凸函数的一阶条件推导出来。将 y 设为最优解 x，我们得到：
*∇F(x)ᵀ (x - x) ≤ F(x) - F(x) ≤ 0**
这表明内积为非正，夹角不小于90度，所以 -∇F(x) 与 (x - x)* 的夹角为锐角。

收敛性定理 📜

在一定的假设下，我们可以证明梯度下降法的收敛速度。这里我们假设函数 F 是 L-Lipschitz 的，即其梯度（或函数值变化）不会太大。一个常用的定义是：对于任意 x, y，有 ‖∇F(x) - ∇F(y)‖ ≤ L ‖x - y‖（或 |F(x)-F(y)| ≤ L‖x-y‖，后者不要求处处可微）。

定理：假设 F 是凸函数且 L-Lipschitz，初始点 x₀ 与最优解 x* 的距离不超过 R。如果选择适当的步长 η（例如 η = R / (L√T)），那么经过 T 步迭代后，输出点的平均函数值满足：
平均[ F(x_t) ] - F(x*) ≤ (R L) / √T

这个定理表明，误差以 O(1/√T) 的速度下降。要获得高精度（例如小数点后10位），需要 O(1/ε²) 量级的迭代步数。

扩展到约束优化：投影梯度下降 🔒

很多时候，我们的变量 x 必须属于某个凸集 K。问题变为在凸集 K 上最小化凸函数 F。

算法需要稍作修改：

1. 像之前一样计算临时点：y_{t+1} = x_t - η · ∇F(x_t)。
2. 将 y_{t+1} 投影回集合 K，得到新的迭代点：x_{t+1} = Proj_K( y_{t+1} )。
   • 投影操作 Proj_K(y) 找到集合 K 中离点 y 最近的点。

为什么可行？ 关键性质是：对于凸集 K，将点投影到 K 上，不会增加该点到 K 内任何一点（特别是最优解 x*）的距离。因此，之前关于“靠近最优解”的直觉仍然成立，相同的收敛定理在约束情况下也基本适用。

投影步骤的计算复杂度取决于凸集 K 的结构，有时可能成为算法的主要瓶颈。存在一些变体算法（如 Frank-Wolfe 算法）可以避免显式的投影。

更快的收敛：强凸与平滑性 🚀

上面的定理没有对函数形状做太多假设，因此收敛较慢。如果函数性质更好，我们可以获得指数级的快速收敛。

两个关键参数是：

• 强凸性 (参数 α)：函数不仅是凸的，而且有“向下弯曲”的趋势，不会太平坦。这保证了存在唯一的最小值点，且函数在最小值点附近像二次函数一样快速上升。
• 平滑性 (参数 β)：函数不会“弯曲”得太厉害，梯度变化不会太剧烈。这保证了我们可以安全地采用较大的步长。

条件数 κ = β / α 衡量了函数最陡方向与最平坦方向的比率，决定了梯度下降的收敛难度。

定理：如果 F 是 α-强凸且 β-平滑的，并且我们知道了 α 和 β，那么可以通过选择特定的步长 η，使得梯度下降法以指数速度收敛：
F(x_t) - F(x*) ≤ exp(-O(t / κ))

这意味着只需要 O(κ log(1/ε)) 步迭代就能达到 ε 精度，远快于之前的 O(1/ε²) 步。

预条件：改善条件数 🔧

如果函数的条件数 κ 很大（即在不同方向上弯曲程度差异巨大），收敛会很慢。预条件 的思想是通过一个线性变换来改变变量的空间，使得新空间中的函数在各个方向上的弯曲程度更加均匀，从而减小条件数。

一个经典的例子是牛顿法，它在每一步都根据当前点的海森矩阵进行局部缩放，相当于动态的预条件。虽然牛顿法在接近最优解时收敛极快，但其理论分析和全局收敛性更复杂。

预条件技术是许多现代快速算法的核心。例如，在计算有向图和无向图的最大流问题时，能否进行有效的预条件是导致算法效率差异的关键之一。

---

总结 📖

本节课我们一起开启了优化算法的学习之旅。我们首先了解了研究生算法课程与本科课程的差异，明确了以探索和深度思考为目标的学习心态。接着，我们深入探讨了梯度下降法这一基础优化工具：

• 我们学习了凸函数的三种等价定义（零阶、一阶、二阶），理解了它作为理想优化“地形”的特性。
• 我们掌握了非约束梯度下降的基本算法步骤，并理解了其背后“函数值下降”和“靠近最优解”的两重直觉。
• 我们看到了一个基础的收敛定理，它给出了在 Lipschitz 连续假设下 O(1/√T) 的收敛速度。
• 我们将算法推广到约束优化问题，引入了投影梯度下降法。
• 最后，我们探讨了在函数具有更好性质（强凸且平滑）时，可以获得指数级快速收敛，并介绍了通过预条件技术来改善收敛速度的核心思想。

下节课，我们将从证明梯度下降法的收敛定理开始，通过具体的符号推导，让这些直观的理解变得更加坚实。

2：凸优化与梯度下降法

在本节课中，我们将继续学习凸优化作为算法设计工具。我们将完成梯度下降法的证明，并介绍其扩展——镜像下降法。课程内容将涵盖梯度下降的分析、投影梯度下降、随机梯度下降，以及如何通过选择不同的正则化项来推导出不同的优化算法。

梯度下降法证明回顾

上一节我们介绍了梯度下降法的基本迭代公式。本节中，我们来看看如何严格证明其收敛性。

梯度下降法用于最小化一个凸函数 ( f )。我们假设函数是 ( L )-Lipschitz 连续的，即其梯度范数有界：( |\nabla f(x)| \leq L )。算法从初始点 ( x_1 ) 开始，迭代更新规则为：
[
x_{t+1} = x_t - \eta \nabla f(x_t)
]
其中 ( \eta ) 是步长。

我们的目标是证明，经过 ( T ) 次迭代后，函数值的平均值与最优值 ( f(x^*) ) 的差距存在上界。

证明的核心步骤如下：

1. 利用凸性：对于凸函数，在任意点 ( x_t ) 处，其梯度方向提供了指向最优解 ( x^* ) 的下降方向信息。具体不等式为：
   [
   f(x_t) - f(x^) \leq \nabla f(x_t)^\top (x_t - x^)
   ]

2. 代入更新规则：将梯度 ( \nabla f(x_t) ) 用迭代公式 ( (x_t - x_{t+1}) / \eta ) 替换，得到：
   [
   f(x_t) - f(x^) \leq \frac{1}{\eta} (x_t - x_{t+1})^\top (x_t - x^)
   ]

3. 应用余弦定理（极化恒等式）：对于向量 ( a, b )，有恒等式 ( 2a^\top b = |a|^2 + |b|^2 - |a - b|^2 )。将其应用于上述不等式的右边项 ( (x_t - x_{t+1})^\top (x_t - x^) )，可得：
   [
   f(x_t) - f(x^) \leq \frac{1}{2\eta} \left( |x_t - x_{t+1}|^2 + |x_t - x^*|2 - |x_{t+1} - x^*|2 \right)
   ]

4. 求和与 telescoping（裂项相消）：将上述不等式对 ( t = 1 ) 到 ( T ) 求和。右边的 ( |x_t - x^*|2 - |x_{t+1} - x^*|2 ) 项会裂项相消，最终只剩下首项 ( |x_1 - x^*|2 ) 和末项 ( -|x_{T+1} - x^*|2 )。由于末项为负，我们可以将其舍弃（这会使不等式更宽松）。于是得到：
   [
   \sum_{t=1}^{T} [f(x_t) - f(x^*)] \leq \frac{1}{2\eta} \sum_{t=1}^{T} |x_t - x_{t+1}|^2 + \frac{1}{2\eta} |x_1 - x^*|2
   ]

5. 利用 Lipschitz 条件：根据更新规则，( |x_t - x_{t+1}| = \eta |\nabla f(x_t)| \leq \eta L )。代入上式：
   [
   \sum_{t=1}^{T} [f(x_t) - f(x^)] \leq \frac{1}{2\eta} \cdot T \cdot (\eta L)^2 + \frac{1}{2\eta} R^2 = \frac{T \eta L^2}{2} + \frac{R^2}{2\eta}
   ]
   其中 ( R = |x_1 - x^| ) 是初始点与最优点的距离。

6. 选择最优步长并求平均：为了最小化右边的上界，我们选择步长 ( \eta = \frac{R}{L\sqrt{T}} )。将其代入并除以 ( T ) 得到平均 regret 的上界：
   [
   \frac{1}{T} \sum_{t=1}^{T} [f(x_t) - f(x^*)] \leq \frac{RL}{\sqrt{T}}
   ]
   这表明，经过 ( T ) 次迭代，函数值的平均值以 ( O(1/\sqrt{T}) ) 的速率收敛到最优值附近。

这个证明的关键在于利用了凸函数梯度提供的方向信息，以及通过裂项求和简化分析。该收敛速率对于仅满足 Lipschitz 条件的凸函数是最优的。

梯度下降法的变体与评论

上一节我们分析了基础梯度下降法。本节中，我们来看看它的几个重要变体和相关性质。

投影梯度下降法

当优化问题带有约束，要求解 ( x ) 必须位于一个凸集 ( K ) 中时，我们使用投影梯度下降法。其迭代步骤如下：

1. 计算梯度步：( y_{t+1} = x_t - \eta \nabla f(x_t) )。
2. 投影回可行域：( x_{t+1} = \Pi_K(y_{t+1}) = \arg\min_{z \in K} |z - y_{t+1}|^2 )。

其分析与标准梯度下降法类似，只需在应用余弦定理时，将 ( x_{t+1} ) 替换为 ( y_{t+1} )，并利用投影算子的一个关键性质：对于凸集 ( K ) 内的任意点 ( x^* )，投影操作不会增加与 ( x^* ) 的距离，即 ( |\Pi_K(y) - x^| \leq |y - x^| )。利用此性质，分析过程几乎完全相同。

梯度下降法的性质

以下是关于梯度下降法的一些重要评论：

• 维度无关性：收敛速率 ( O(RL/\sqrt{T}) ) 与问题维度 ( n ) 无关。这意味着该方法可应用于高维甚至无限维空间。
• 收敛速率：( O(1/\sqrt{T}) ) 的收敛速率对于非光滑的 Lipschitz 凸函数是最优的。要达到精度 ( \epsilon )，需要 ( O(1/\epsilon^2) ) 次迭代。从计算精度角度看，这属于次多项式时间。若要求误差为 ( 2^{-B} )（B 比特精度），则迭代次数为 ( O(2^{2B}) )，是指数级的。
• 梯度估计：算法只需要计算（次）梯度。在实践中，即使只能访问函数值，也可以通过有限差分等方式估计梯度。
• 随机梯度下降：当只能获得梯度的无偏估计 ( \tilde{\nabla} f(x) )（即 ( \mathbb{E}[\tilde{\nabla} f(x)] = \nabla f(x) )）时，可以使用随机梯度下降。上述分析框架经过适当修改后仍然适用，因为期望运算可以与证明中的内积等操作交换。
• 更快的收敛：如果函数具有更好的性质，收敛可以更快：
  • 光滑函数：如果函数梯度是 Lipschitz 连续的（即函数是光滑的），可以获得 ( O(1/T) ) 的收敛速率。
  • 光滑且强凸函数：如果函数同时是光滑且强凸的，梯度下降法可以实现线性（指数）收敛，即误差按 ( O(c^T) )（( c < 1 )）衰减。这是真正的多项式时间算法，因为要达到 B 比特精度，只需要 ( O(\log(1/\epsilon)) = O(B) ) 次迭代。

从梯度下降到镜像下降法

前面我们看到，梯度下降法在非光滑情形下收敛较慢，且其更新规则在数学上存在一个“类型不匹配”的问题：参数 ( x ) 属于原始空间，而梯度 ( \nabla f(x) ) 本质上是对偶空间中的线性泛函。直接相加在几何上并非最自然。本节我们介绍一种更通用的视角——镜像下降法。

梯度下降的另一种形式：近端梯度视角

梯度下降法可以重新表述为以下“近端”形式：
在每一步 ( t )，寻找下一个点 ( x_{t+1} )，使其最小化线性近似加上一个二次惩罚项：
[
x_{t+1} = \arg\min_{x} \left{ \eta \nabla f(x_t)^\top x + \frac{1}{2} |x - x_t|^2 \right}
]
求解这个二次函数的最小值（令梯度为零），恰好得到标准的梯度下降更新：( x_{t+1} = x_t - \eta \nabla f(x_t) )。

这种观点的优势在于，它明确地将算法分解为两部分：

1. 模型：在当前点对目标函数进行一阶（线性）近似。
2. 正则化/信任域：添加一个惩罚项 ( \frac{1}{2} |x - x_t|^2 )，防止新点离当前点太远，确保线性近似仍然有效。

镜像下降法：推广正则化项

镜像下降法的核心思想是：将上述二次惩罚项 ( \frac{1}{2} |x - x_t|^2 ) 替换为一个更一般的“距离”函数，即 Bregman 散度。

给定一个严格凸函数 ( h )（称为 距离生成函数），其对应的 Bregman 散度定义为：
[
D_h(y, x) = h(y) - h(x) - \nabla h(x)^\top (y - x)
]
几何上，( D_h(y, x) ) 度量了函数 ( h ) 在 ( y ) 点的值与其在 ( x ) 点线性近似之间的差距。对于凸函数 ( h )，这个值总是非负的，并且在 ( y = x ) 时为零。当 ( h(x) = \frac{1}{2} |x|^2 ) 时，( D_h(y, x) = \frac{1}{2} |y - x|^2 )，就回到了梯度下降的情形。

镜像下降法的迭代步骤如下：
[
x_{t+1} = \arg\min_{x \in \mathcal{K}} \left{ \eta \nabla f(x_t)^\top x + D_h(x, x_t) \right}
]
其中 ( \mathcal{K} ) 是可行域。

镜像下降法的实例

通过选择不同的距离生成函数 ( h )，我们可以恢复多种重要的算法：

1. 梯度下降法：选择 ( h(x) = \frac{1}{2} |x|_2^2 )，则 ( D_h(x, x_t) = \frac{1}{2} |x - x_t|_2^2 )。
2. 乘法权重更新/专家算法：在概率单纯形上，选择 ( h(x) = \sum_i x_i \log x_i )（负熵），则对应的 Bregman 散度是 KL 散度：( D_h(y, x) = \sum_i y_i \log(y_i / x_i) )。由此导出的镜像下降法就是乘法权重更新算法。
3. 预条件梯度下降：选择 ( h(x) = \frac{1}{2} x^\top A x )，其中 ( A ) 是一个正定矩阵。对应的 Bregman 散度是 ( D_h(x, x_t) = \frac{1}{2} (x - x_t)^\top A (x - x_t) )。求解镜像下降更新步骤会得到 ( x_{t+1} = x_t - \eta A^{-1} \nabla f(x_t) )。这相当于在不同坐标方向上使用不同的步长，可以适应目标函数在不同方向上的曲率，从而加速收敛。
4. 牛顿法：一个更激进的想法是，在每一步使用目标函数 ( f ) 本身的局部二阶近似作为正则项。即，令 ( h_t(x) = \frac{1}{2} (x - x_t)^\top \nabla^2 f(x_t) (x - x_t) )，其中 ( \nabla^2 f(x_t) ) 是 Hessian 矩阵。这等价于在每一步最小化 ( f ) 的二阶泰勒展开，这正是牛顿法。牛顿法在初始点足够靠近最优解时，具有极快的超线性收敛速度（如 ( O(\log \log(1/\epsilon)) )）。

镜像下降法的优势在于，它允许我们根据问题的几何结构选择合适的距离生成函数 ( h )。例如，在概率单纯形上，基于 KL 散度的几何（乘法权重）比基于欧氏距离的几何（普通梯度下降）更为自然和高效。通过“镜像”映射（由 ( h ) 的梯度映射及其逆映射定义），算法将对偶空间中的梯度信息“映射”回原始空间进行更新，从而解决了类型不匹配的问题，并使得算法能自适应问题的几何。

总结

本节课中我们一起学习了凸优化中梯度下降法的完整证明及其重要变体。我们首先回顾并严格证明了梯度下降法对于 Lipschitz 凸函数的 ( O(1/\sqrt{T}) ) 收敛速率。接着，我们讨论了投影梯度下降、随机梯度下降，以及函数具备光滑性或强凸性时可能获得的更快收敛速率。

随后，我们引入了更一般的镜像下降法框架。通过将梯度下降重新解释为“线性近似加二次正则”的近端形式，我们自然地将二次正则项推广为 Bregman 散度。这一推广使我们能够根据问题域的内在几何选择不同的距离函数，从而得到如乘法权重算法、预条件梯度下降等一系列重要算法。镜像下降法统一了这些看似不同的算法，并强调了在算法设计中考虑问题几何结构的重要性。

3：镜像下降与在线凸优化

在本节课中，我们将要学习镜像下降算法的证明框架，并探讨其在线凸优化中的应用，特别是如何利用它来求解线性规划问题。

上一节我们介绍了镜像下降算法的基本概念，本节中我们来看看其证明思路以及更广泛的应用。

镜像下降算法回顾

镜像下降算法的起点是一个强凸函数 H(x)。例如，H(x) = (1/2)||x||² 就是一个强凸函数。强凸函数满足以下性质：对于任意 x, y，有 H(y) ≥ H(x) + ∇H(x)ᵀ(y-x) + (μ/2)||y-x||²，其中 μ > 0。

给定强凸函数 H，我们可以定义布雷格曼散度 D(y, x)。这是一个非对称的“距离”度量，定义为：
D(y, x) = H(y) - H(x) - ∇H(x)ᵀ(y-x)

布雷格曼散度衡量了点 y 相对于点 x 的线性近似之上的高度差。

基于此，镜像下降算法定义如下。假设当前点为 x_t，函数 F 在 x_t 处的线性近似为 F(x_t) + ∇F(x_t)ᵀ(y - x_t)。算法通过最小化该线性近似加上一个由布雷格曼散度构成的惩罚项，来找到下一个点 y_{t+1}：
y_{t+1} = argmin_y [ η ∇F(x_t)ᵀ y + D(y, x_t) ]

其中 η 是步长。通过求导并令导数为零，我们可以得到关键方程：
∇H(y_{t+1}) = ∇H(x_t) - η ∇F(x_t)

如果优化问题带有约束 x ∈ K，我们还需要将 y_{t+1} 投影回可行域 K，得到 x_{t+1}。投影使用相同的布雷格曼散度：
x_{t+1} = argmin_{x ∈ K} D(x, y_{t+1})

因此，镜像下降的一步可以形象地理解为：将当前点 x_t 通过梯度映射 ∇H 送到对偶空间，在对偶空间中执行梯度步 -η∇F(x_t)，然后通过逆映射 (∇H)⁻¹ 映射回原空间得到 y_{t+1}，最后将其投影回可行域 K 得到 x_{t+1}。

镜像下降的收敛性证明（概要）

以下是镜像下降收敛性定理的证明概要。假设 H 是强凸函数，∇F 的范数以 L 为界，则经过 T 步迭代后，有：
∑_{t=1}^T [F(x_t) - F(x)] ≤ (1/η) D(x, x_1) + (η T L²) / (2μ)

证明的核心步骤利用了凸性、布雷格曼散度的“余弦定理”以及投影的“勾股定理”。

1. 利用凸性：由 F 的凸性，对于最优解 x，有 **F(x_t) - F(x) ≤ ∇F(x_t)ᵀ (x_t - x*)**。
2. 代入梯度关系：利用从算法推导出的关系 η ∇F(x_t) = ∇H(x_t) - ∇H(y_{t+1})，将梯度项替换。
3. 应用“余弦定理”：对于布雷格曼散度，有恒等式：
   ⟨∇H(a) - ∇H(b), a - c⟩ = D(c, a) + D(a, b) - D(c, b)
   这类似于向量空间中的余弦定理。
4. 应用投影的“勾股定理”：对于投影 x_{t+1} = Π_K(y_{t+1})，有不等式：
   D(x, y_{t+1}) ≥ D(x, x_{t+1}) + D(x_{t+1}, y_{t+1})
   这类似于钝角三角形的边长关系。
5. 求和与 telescoping（裂项相消）：将上述不等式对 t 从 1 到 T 求和。利用散度项的 telescoping 性质，许多中间项会相互抵消。
6. 边界处理：剩余项可以通过梯度范数界 L 和强凸系数 μ 进行界定。

最终整理即可得到定理中的 regret 上界。选择最优步长 η ∝ 1/√T，可得平均 regret 以 O(1/√T) 的速率收敛到零。

在线凸优化

上一节我们介绍了镜像下降的收敛性，本节中我们来看看如何将其框架应用于在线凸优化。

在线凸优化问题描述如下：在每一轮 t，算法从凸集 K 中选择一个点 x_t。随后，对手（或环境）揭示一个凸成本函数 f_t，算法产生成本 f_t(x_t)。经过 T 轮后，算法的目标是最小化 regret（遗憾），即其总成本与事后看来最好的固定单点决策的总成本之差：
Regret_T = ∑_{t=1}^T f_t(x_t) - min_{x ∈ K} ∑_{t=1}^T f_t(x*)*

一个关键点是，对手在选择 f_t 时，可以基于算法之前所有的决策 x_1, ..., x_t。

令人惊讶的是，标准的（镜像）下降算法只需稍作修改——即在每一轮 t 使用当前揭示的函数 f_t 的梯度——就能直接应用于在线设置，并且其分析几乎不变。我们之前证明的关于 ∑[F(x_t) - F(x*)] 的界，现在直接对应于 Regret_T 的上界。因此，我们立即得到，在线（镜像）下降算法能以 O(√T) 的速率控制 regret，即平均 regret Regret_T / T → 0。

这看似矛盾：算法根据历史信息更新，而未来的成本函数可能与之无关甚至恶意针对。其深层原因在于，对手只能选择凸函数。为了在算法当前点 x_t 处造成高成本，凸的成本函数不可避免地也会对可行域 K 中的其他点（包括最优固定点 x*）造成伤害，从而限制了 regret 的增长。

应用：求解线性规划

在线凸优化框架的一个巧妙应用是求解线性规划。考虑寻找满足以下线性约束的可行点 x：
l_i(x) = w_iᵀ x - b_i ≤ 0, for i = 1,..., m，且 x ∈ [-B, B]^n

我们可以设计一个算法，它同时扮演在线凸优化算法和对手的角色：

1. 算法（玩家）：运行在线镜像下降算法。在每一轮 t，它输出一个候选点 x_t。
2. 对手：检查 x_t。如果存在某个约束 i 使得 l_i(x_t) > ε（即违反超过 ε），则选择该约束对应的线性函数 l_i 作为本轮的成本函数 f_t 发送给算法。如果所有约束都被满足到 ε 精度内，则算法终止并返回 x_t 作为近似可行解。

为什么这个算法有效？

假设存在一个严格可行解 x（满足所有约束）。对于这个点，由于所有 **l_i(x) ≤ 0**，它在整个过程中的总成本 ∑ f_t(x*) ≤ 0。另一方面，只要算法还在运行，对手每一轮都能找到一个违反超过 ε 的约束，因此算法每一轮的成本 f_t(x_t) > ε。这意味着，如果算法运行了 T 轮，则其 regret Regret_T = ∑ f_t(x_t) - ∑ f_t(x*) > εT - 0 = εT。

然而，我们知道在线下降算法能保证 Regret_T ≤ O(√T)。因此，这种高 regret（线性于 T）的情况不可能持续太久。算法必然会在 T = O(1/ε²) 轮内停止，并返回一个满足所有约束到 ε 精度的解。

这个算法简单而强大，它将线性规划的可行性问题转化为在线优化问题。其运行时间是 poly(n, m, B, 1/ε)。需要注意的是，这是关于精度 1/ε 的多项式，而非关于 log(1/ε)（即比特数的多项式）。后者是内点法等更复杂算法所能达到的“强多项式时间”目标。

总结

本节课中我们一起学习了：

1. 回顾了镜像下降算法，它通过一个强凸函数和对应的布雷格曼散度来定义迭代和投影步骤。
2. 概述了镜像下降的收敛性证明，其核心利用了布雷格曼散度的几何性质（余弦定理、投影勾股定理）。
3. 将镜像下降框架扩展到在线凸优化设置，算法通过使用每一轮揭示的成本函数的梯度，自然获得次线性 regret 界。
4. 展示了如何利用在线凸优化算法来求解线性规划可行性问题。通过让算法自身扮演对手来提供违反约束的梯度，我们得到了一个简单且理论保证的近似求解算法，其运行时间关于精度为多项式级。

这个框架展示了凸优化基础工具的通用性和强大之处。

4：质心算法与椭球算法

在本节课中，我们将学习两种用于凸优化的算法：质心算法和椭球算法。这两种算法的关键特点是，它们能以 poly(log(1/ε)) 的运行时间达到 ε 精度，这与我们之前看到的梯度下降类算法（运行时间为 poly(1/ε)）不同。这种对数级别的收敛速度使得椭球算法能够证明线性规划是多项式时间可解的。

质心算法

上一节我们介绍了算法的目标。本节中，我们来看看质心算法的基本思想。该算法旨在最小化定义在凸集 K 上的凸函数 f(x)。其核心思路类似于二分查找：我们不断切割凸集，每次切割都排除掉函数值肯定不低于当前点的一半区域。

算法步骤

以下是质心算法的基本步骤：

1. 初始化：设当前凸集为 K₁ = K。
2. 循环：对于第 t 步：
   • 计算当前凸集 Kₜ 的质心 cₜ。质心的定义是：cₜ = (∫_K x dx) / vol(K)，其中 vol(K) 表示凸集 K 的体积。对于凸集，其质心位于集合内部。
   • 在质心 cₜ 处计算函数 f 的梯度 ∇f(cₜ)。
   • 利用凸性，我们知道对于所有满足 ∇f(cₜ)ᵀ (y - cₜ) ≥ 0 的点 y，有 f(y) ≥ f(cₜ)。因此，我们可以安全地排除掉这一半空间。
   • 更新凸集：Kₜ₊₁ = Kₜ ∩ { y | ∇f(cₜ)ᵀ (y - cₜ) < 0 }。
3. 输出：返回所有质心点中函数值最小的那个，即 min{ f(c₁), f(c₂), ..., f(cₜ) }。

算法分析的关键

为什么选择质心？原因类似于二分查找选择中点：我们希望每一步都能显著减小搜索空间的“大小”。这里，我们使用体积作为度量。有一个关键定理（Grünbaum 定理）保证了这种切割的有效性：对于任何凸集 K，任何通过其质心的超平面所划分出的两个部分，其体积都满足 vol(K ∩ H) ≤ (1 - 1/e) * vol(K)。这意味着每次切割后，剩余凸集的体积至少会以一个常数因子（约 1/e）减小。

收敛性解释

假设最优解 x* 附近存在一个小的“好”区域（例如，一个以 x* 为中心、按比例 δ 缩放的凸集 K_δ），该区域内函数值接近最优。这个区域的体积是 vol(K) * δⁿ。由于算法每一步都以常数因子减小体积，在约 O(n log(1/δ)) 步后，当前凸集的体积将小于这个“好”区域的体积。这意味着算法必然在某个时刻切割到了这个区域，从而保证最终找到的质心具有接近最优的函数值。

---

上一节我们介绍了质心算法的原理，但它还不是一个完全明确的算法，主要因为计算凸集的质心通常很困难。本节中我们来看看一个更实用、理论上更强大的算法——椭球算法。

椭球算法

椭球算法的重要性在于它证明了线性规划是多项式时间可解的。更一般地说，它能够处理那些约束条件可能指数级多，但我们可以通过一个“分离预言机”高效访问的凸优化问题。

分离预言机

分离预言机是我们访问凸集 K 的方式。它是一个算法，当我们输入一个点 z 时：

• 如果 z ∈ K，预言机回答“是”。
• 如果 z ∉ K，预言机回答“否”，并且额外给出一个分离超平面 (a, b)，使得对于所有 x ∈ K 有 aᵀx ≥ b，但 aᵀz < b。这个超平面将点 z 与凸集 K 分离开。

算法基本思想

椭球算法的目标是找到一个属于凸集 K 的点（可行性问题）。它通过维护一个椭球序列来工作，每个椭球都包含整个凸集 K。椭球可以用其中心 c 和一个正定矩阵 Q 简洁地表示：E = { x | (x-c)ᵀ Q (x-c) ≤ 1 }。

1. 初始化：找到一个足够大的球（椭球 E₀），确保它包含整个凸集 K。设其中心为 c₀。
2. 迭代：对于第 t 步：
   • 查询分离预言机：当前椭球中心 cₜ 是否在 K 中？
   • 如果“是”，则算法结束，返回 cₜ。
   • 如果“否”，预言机返回一个分离超平面。这个超平面将 cₜ 与 K 分开，并且 K 完全位于超平面的另一侧。
   • 关键步骤：我们保留包含 K ∩ Eₜ 的那一半椭球（即与 K 同侧的那一半），并寻找一个新的、体积更小的椭球 Eₜ₊₁，使其完全包含这一半椭球。
3. 终止：由于每一步新椭球的体积都会以一个因子（约 exp(-1/(2n+2))）减小，而凸集 K 本身包含一个小球（故体积非零），因此在多项式步数内，椭球体积将变得比 K 的体积还小。这必然导致某个椭球的中心落入 K 内，算法成功找到可行点。

从可行性到优化

如何将最小化凸函数 f 的问题转化为可行性问题？我们可以进行二分搜索：猜测一个目标值 θ，然后利用分离预言机检查凸集 K ∩ { x | f(x) ≤ θ } 是否非空。通过不断收紧 θ，即可找到最小可能的值。

技术细节与解释

为什么是椭球？ 椭球是对称的，且在任何线性变换下仍为椭球。这允许我们将任意椭球被超平面切割的问题，通过线性变换转化为一个单位球被坐标平面切割的标准问题。对于这个标准问题，我们可以显式地构造出一个体积显著减小且包含半个球的新椭球公式，然后再变换回原空间。这种“归一化”处理使得分析变得简洁。

与质心算法的对比：质心算法需要维护和计算复杂切割后的凸集，而椭球算法始终维护一个结构简单（只需中心和矩阵）的椭球，这是其可实现的关键。

---

本节课中我们一起学习了两种用于凸优化的算法。质心算法提供了清晰的迭代切割思路，但其实现依赖于质心计算。椭球算法则是一个理论上的强大工具，它通过维护一系列体积严格递减的椭球，并利用分离预言机，证明了线性规划乃至一大类凸优化问题是多项式时间可解的。尽管其实用性可能不如内点法或针对特定问题的算法，但椭球算法在理论计算机科学中具有里程碑式的意义。

5：线性代数工具入门

在本节课中，我们将从优化问题转向一个新的主题：线性代数工具。我们将首先探讨求解线性系统的不同方法，然后介绍主成分分析（PCA）和奇异值分解（SVD）的基本概念及其计算方法。

从优化到线性代数

上一节我们介绍了多种优化方法，如梯度下降和镜像下降。本节中，我们来看看线性代数中的核心问题：求解线性系统。

线性系统的基本形式是 A x = b，其中 A 是一个矩阵，x 是未知向量，b 是已知向量。求解此系统主要有两类方法。

以下是两类主要方法的概述：

• 直接法：如高斯消元法及其变体（如LU分解）。这类方法不依赖于矩阵的条件数，总能得到精确解，但通常需要 O(n³) 的时间和 O(n²) 的空间，且难以利用矩阵的稀疏性。
• 迭代法：如梯度下降、最速下降法和共轭梯度法。这类方法从一个初始解开始，通过迭代逐步改进。其优势在于只需 O(n) 的额外空间，并且能有效利用矩阵的稀疏性。但收敛速度通常依赖于矩阵的条件数 κ。

用梯度下降求解线性系统

现在，我们具体看看如何用梯度下降法求解线性系统 A x = b。首先，我们可以通过转换，将问题转化为求解一个对称正定（PSD）矩阵的系统。

具体做法是求解 AᵀA x = Aᵀb。新矩阵 AᵀA 是实对称且半正定的。我们的目标是最小化一个凸二次函数。

我们定义目标函数 f(x) = (1/2) xᵀ A x - bᵀ x。最小化此函数等价于求解 ∇f(x) = A x - b = 0。梯度下降的更新规则如下：

x_{t+1} = x_t - η * (A x_t - b)

每一步迭代的核心计算是一次矩阵向量乘法 A x_t。如果矩阵 A 是稀疏的，这一步可以非常高效。

收敛性分析表明，迭代误差与最优解 x* 的距离满足：
| x_{t+1} - x | ≤ | I - η A | * | x_t - x |**
其中，| I - η A | 是算子范数。通过选择合适的步长 η（例如 η = 1 / λ_max），误差会以指数速度衰减，收敛速率与矩阵条件数 κ 相关。

预条件梯度下降

上一节我们分析了基础梯度下降的收敛性。本节中我们来看看如何通过预条件技术来改进它，这本质上是镜像下降思想的应用。

在预条件梯度下降中，我们引入一个预条件矩阵 L，更新规则变为：
x_{t+1} = x_t - η L⁻¹ (A x_t - b)
此时，收敛速度取决于 L⁻¹A 的条件数，而非 A 本身的条件数。

因此，问题的关键转化为：寻找一个矩阵 L，使得 L⁻¹A 的条件数很小，并且能够快速求解形如 L z = w 的线性系统。这实际上将求解 A x = b 的问题，约化到了求解更简单的 L z = w 的问题。在实践中，L 常被选为近似于 A 但结构更简单的矩阵（如树状图）。

主成分分析（PCA）与奇异值分解（SVD）

现在，我们离开线性系统求解，进入数据降维领域。一个核心问题是：给定一组高维数据点，如何找到其最佳的低维近似？

假设我们有 n 个 d 维数据点。最佳的一维近似是找到一条直线（方向向量 v），使得所有数据点到该直线的投影距离平方和最小，或者说投影长度平方和最大。

更一般地，最佳 k 维子空间可以通过一种贪心算法找到，即主成分分析（PCA）。

以下是PCA贪心算法的步骤：

1. 寻找第一个方向 v₁，最大化所有数据点在该方向上投影长度的平方和。
2. 将数据投影到与 v₁ 正交的空间中。
3. 在正交空间中重复步骤1，寻找 v₂，依此类推，直到找到 k 个方向。

这些方向 v₁, v₂, ..., vₖ 张成的子空间就是数据的最佳 k 维近似。PCA能有效揭示数据的内在结构。

幂方法计算主成分

那么，如何计算这个最大化投影平方和的首个方向 v₁ 呢？这等价于求解矩阵 B = AᵀA 的最大特征值对应的特征向量。

对于大型稀疏矩阵，我们通常使用迭代法而非直接分解。幂方法 是一种经典的迭代算法。

幂方法的过程如下：

1. 随机初始化一个向量 x₀。
2. 重复进行迭代：x_{t+1} = B x_t / | B x_t |。
3. 经过足够多次迭代后，x_t 将收敛到 B 的主特征向量方向。

其原理在于，任何初始向量都可以表示为各特征向量的线性组合。每次乘以 B，各成分会乘以对应的特征值 λᵢ。最大特征值 λ_max 对应的成分将随着迭代次数的增加而占据主导地位。

收敛速度取决于最大特征值与次大特征值之间的谱间隙。间隙越大，收敛越快。

要计算前 k 个主成分，可以使用块幂方法，即同时迭代一组正交向量，并在每次迭代后对它们进行正交化。

总结

本节课中我们一起学习了线性代数中的两大核心工具。首先，我们探讨了求解线性系统的迭代方法，特别是梯度下降及其预条件变体，它们能利用稀疏性，但收敛性依赖于条件数。接着，我们介绍了主成分分析（PCA）的概念，它用于寻找数据的最佳低维近似，并讲解了如何使用幂方法迭代计算主成分。这些线性代数工具是理解和处理高维数据、图论问题以及许多机器学习算法的基础。

6：鲁棒均值估计与高维统计

在本节课中，我们将学习如何在存在对抗性异常值的情况下，从高维数据中稳健地估计均值。我们将探讨一维和高维情况下的挑战，介绍一些核心概念，并了解一些实用的算法思想。

概述

我们从一维的鲁棒均值估计问题开始。假设我们从一个方差为1的标准高斯分布中采样，但有一个对手可以添加一小部分（比例为 ε）的任意异常值。我们的目标是估计真实均值 μ，尽管存在这些异常值。

在一维情况下，一个简单而有效的方法是使用中位数。可以证明，中位数与真实均值的偏差在 O(ε) 范围内，无论对手如何添加异常值。

公式：|median_hat - μ| ≤ O(ε)

然而，当我们转向高维（D维）时，问题变得复杂。一个自然的初步尝试是计算几何中位数，即分别计算每个坐标的中位数。

代码（概念性）：

mu_hat = []
for d in range(D):
    mu_hat_d = median([x[d] for x in samples])
    mu_hat.append(mu_hat_d)

但这种方法存在问题。在每个坐标上，估计值可能偏离真实值 O(ε)。在 D 个坐标上，总的平方误差可能达到 O(ε²D)，导致欧几里得距离误差为 O(ε√D)。这种损失源于我们随意选择了一个坐标基（例如标准基），而问题本身并没有赋予任何方向特殊性。

上一节我们介绍了几何中位数及其局限性，本节中我们来看看一个更强大的概念：Tukey 中位数。

Tukey 中位数

为了获得与维度无关的误差保证，我们需要一个在每个方向上都看起来像“中间”的点。这就是 Tukey 中位数 的思想。

对于数据集中的每个点 x，我们定义其深度为：在所有可能的方向 v 上，点 x 沿该方向投影后排序位置的最小值。

公式：depth(x) = min over all directions v of (position of projection of x along v in sorted list)

Tukey 中位数就是深度最大的那个数据点。这个点在任何方向上看起来都处于中间位置（即，两侧都有相当数量的点）。这确实是一维中位数在高维的自然推广，并且能提供维度无关的误差保证。

然而，计算 Tukey 中位数在计算上是困难的，目前没有已知的高效算法。

鲁棒均值估计的算法直觉

既然 Tukey 中位数难以计算，我们需要寻找其他途径。让我们思考一下对手如何欺骗我们。

对手的目标是使用 ε 比例的异常值，将估计的均值从真实均值 μ 推离一个距离 Δ。在一维中，为了将均值移动 Δ，对手需要将异常值放置在距离大约为 Δ/ε 的地方。这会对方差产生显著影响：异常值对方差的贡献约为 ε * (Δ/ε)² = Δ²/ε。

由于真实高斯分布的方差为1（协方差矩阵近似为单位矩阵），这种对方差的巨大扰动是可以被检测到的。在高维中，对手必须将其所有“努力”集中在某个特定方向 v 上，才能有效地移动均值。这会导致数据沿方向 v 的方差（即协方差矩阵在 v 方向上的特征值）异常增大，大约为 Δ²/ε。

因此，算法的核心直觉是：

1. 计算数据样本的协方差矩阵。
2. 识别出特征值异常大的方向（这些方向可能被异常值主导）。
3. 以某种方式“修剪”或降低这些方向上的影响，直到数据看起来像一个“好的”球体（即协方差矩阵接近单位矩阵）。
4. 用剩余数据估计均值。

这个思想不依赖于维度 D，有望获得维度无关的误差界。

均值估计的扩展讨论

在深入鲁棒估计算法之前，我们暂时离开主线，讨论一个相关的基础问题：即使在没有异常值的普通情况下，如何更好地估计均值？

经验均值的局限性

给定从某个分布 D 中抽取的 n 个一维样本，估计均值最直接的方法是计算经验均值。

公式：μ_hat = (1/n) * Σ_{i=1 to n} x_i

如果 D 是高斯分布，经验均值表现优异：估计误差超过 ε 的概率小于 δ，所需样本数 n 为 O(log(1/δ)/ε²)。误差概率随样本数呈指数衰减。

然而，对于许多现实中的重尾分布（如幂律分布），经验均值只能提供多项式收敛速度，而非指数收敛。

均值中位数算法

一个简单而强大的改进是 均值中位数 算法：

1. 将 n 个样本随机分成 k 个桶。
2. 计算每个桶内样本的均值，得到 k 个均值。
3. 输出这 k 个均值的中位数。

代码（概念性）：

def median_of_means(data, k):
    # 将数据随机分成k组
    groups = split_into_k_groups(data, k)
    # 计算每组均值
    group_means = [np.mean(group) for group in groups]
    # 返回均值的中位数
    return np.median(group_means)

这个算法非常巧妙。无论底层分布 D 是什么，它都能提供接近指数衰减的置信度。对于给定的精度 ε 和置信度 1-δ，所需样本数约为 O(log(1/δ)/ε²)，与高斯分布下的最优界仅差一个常数因子。参数 k 通常设置为 O(log(1/δ))。

高维挑战与最新进展

均值中位数算法本质是一维的。在高维中，我们需要一个“中位数”的替代品，这又回到了类似 Tukey 中位数的问题。直到大约一年半前，研究者才发现了第一个高效算法，能够在高维中实现与高斯分布类似的样本复杂度（O(log(1/δ)/ε²)），达到常数因子级别的最优性。

这个领域被称为重尾统计。

高维几何的注记

在讨论高维算法时，理解高维空间的一些反直觉性质很有帮助。例如，在高维单位球面上随机选取两个向量，它们的点积（夹角余弦）的绝对值很可能非常小，大约为 O(1/√D)。这意味着随机向量几乎总是近似正交的。此外，高维球体的大部分体积都集中在表面附近一个极薄的壳层中。这些性质是许多高维统计现象的基础。

总结

本节课中我们一起学习了：

1. 鲁棒均值估计问题：在存在对抗性异常值的情况下估计分布均值。
2. 一维解决方案：使用中位数可以获得 O(ε) 的误差。
3. 高维挑战：简单的坐标中位数（几何中位数）会引入 O(ε√D) 的误差。
4. Tukey 中位数：一个理想但难以计算的高维中位数概念。
5. 算法直觉：通过检测并修剪协方差矩阵中特征值过大的方向（可能由异常值引起）来估计均值。
6. 均值估计的基础：即使在没有异常值时，经验均值对重尾分布表现不佳，而均值中位数算法能提供更稳健的指数级收敛保证。
7. 高维特性：高维空间中的随机向量倾向于近似正交，这对理解高维统计至关重要。

这些概念为设计实际可行的、具有理论保证的鲁棒学习算法奠定了基础。

7：鲁棒均值估计与张量简介 🎯

在本节课中，我们将要学习鲁棒均值估计算法的完成部分，并初步了解张量这一数据结构。我们将首先探讨如何从包含异常值的数据集中稳定地估计真实均值，然后介绍张量的基本概念及其应用。

鲁棒均值估计算法完成 🧮

上一节我们介绍了鲁棒均值估计的问题背景。本节中我们来看看如何完成该算法。

我们有一个数据集，它是一组点 X1 到 Xn。其中 1 - ε 比例的点是内点（inliers），ε 比例的点是由对手插入的异常值（outliers）。对手可以查看我们拥有的数据点，并决定添加哪些异常值以及添加多少（少于 ε 比例）。我们的目标是估计内点分布的均值。

稳定性概念 🔒

为了实现目标，我们引入一个重要的概念：稳定性。

首先为数据集定义稳定性。一个点集是 ε-δ 稳定 的，如果移除任意 ε 比例的点，其均值的变化最多为 δ。形式化地说，对于任意至少包含 (1 - ε)n 个元素的子集，该子集的均值必须接近整个数据集的均值。

对于概率分布，我们可以类似地定义稳定性。一个分布 θ 是 ε-δ 稳定 的，如果对于它的每一个 ε-过滤（即通过“删除”至多 ε 比例概率质量得到的分布），其均值的变化最多为 δ。分布的 ε-过滤 θ' 定义为：对于所有点 x，有 P_θ'(x) ≤ (1/(1-ε)) * P_θ(x)。

我们的核心假设是：内点集本身是 ε-δ 稳定的。对于从高斯分布中抽取的样本，这很合理，因为移除少量样本不会显著改变其经验均值。

算法思路与关键引理 🧠

算法的核心思路是：找到一个同样是稳定集的子集 S。如果内点集 D_in 和找到的集合 S 都是稳定的，那么它们的均值必定接近。

证明思路：考虑两个稳定集 D_in 和 S 的交集 T。由于 D_in 稳定，μ(D_in) 接近 μ(T)。由于 S 稳定，μ(S) 也接近 μ(T)。因此，根据三角不等式，μ(D_in) 接近 μ(S)。

所以，问题转化为：如何在包含异常值的数据集中找到一个稳定子集？

为此，我们需要一个稳定性的高效充分条件。关键引理如下：

对于任何分布 θ，设其协方差矩阵的算子范数（最大特征值）为 λ。那么 θ 是 ε-δ 稳定 的，其中 δ = √(ε) * λ。

这意味着，一个分布只要其协方差矩阵的最大特征值不大，它就是稳定的。标准高斯分布的协方差矩阵是单位阵，因此非常稳定。

引理证明概要：

1. 将原始分布 θ 表示为过滤后分布 θ_in 与被过滤部分 θ_out 的凸组合：θ = (1-γ)θ_in + γθ_out。
2. 写出 θ 的协方差矩阵表达式，它包含 θ_in 和 θ_out 的协方差以及均值差项。
3. 令 v 为均值差 μ(θ_in) - μ(θ_out) 方向的单位向量，将协方差等式两边同时左乘 v^T 和右乘 v。
4. 左边 v^T Cov(θ) v ≤ λ。
5. 右边，协方差项非负，均值差项变为 2γ(1-γ) * ||μ(θ_in) - μ(θ_out)||^2。
6. 结合两边得到 ||μ(θ_in) - μ(θ_out)||^2 ≤ λ / (2γ(1-γ))。
7. 最终，原始均值与过滤后均值的差 ||μ(θ) - μ(θ_in)|| = γ * ||μ(θ_in) - μ(θ_out)|| ≤ √(ε) * λ（经过参数代入）。

寻找稳定集：凸优化表述 ⚙️

现在我们知道，要找到一个稳定集 S，等价于找到一组权重 w_i（形成一个分布），使得：

1. w_i ≥ 0，且 ∑ w_i = 1。
2. w_i ≤ (1+ε)/n（这是一个 ε-过滤约束）。
3. 该权重分布下的协方差矩阵 Cov(w) = ∑ w_i (X_i - μ(w))(X_i - μ(w))^T 的算子范数很小（例如 ≤ λ I）。

这里 μ(w) = ∑ w_i X_i 是该权重下的均值。

注意，第三个约束关于权重 w 是三次的（因为 μ(w) 中也包含 w），因此不是凸约束。然而，我们可以利用椭球法来求解。

椭球法思路：

• 可行域是满足条件1和2的凸集（一个单纯形）。
• 我们知道内点的均匀权重 w* 是一个可行解（满足所有约束）。
• 对于任意一个试探解 w，如果其协方差算子范数 λ 太大，我们可以构造一个分离超平面，将 w 和理想解 w* 分开，并告诉椭球法下一步搜索的方向。
• 具体地，设 v 是 Cov(w) 的 top 特征向量。定义线性函数 L(y) = ∑_i y_i * [ (v^T(X_i - μ(w)) )^2 ]。
• 可以验证 L(w) = λ，并且可以证明 L(w*) 大约为 O(ελ)，远小于 λ。
• 因此，超平面 L(y) ≤ λ 将 w 和 w* 分离开来。椭球法利用这个切割平面逐步逼近 w*。

这样就找到了一个稳定的子集（对应权重 w），其均值可作为对真实内点均值的鲁棒估计。

其他方法：也存在更高效的过滤算法，其核心思想是反复计算协方差矩阵，移除与 top 特征向量方向“对齐”的极端点，并确保每次移除的异常点比内点多。虽然分析更复杂，但效率更高。

张量简介 📦

上一节我们完成了鲁棒均值估计。本节中我们来看看一种新的数据结构：张量。

什么是张量？

从操作上讲，张量是一个高维数组。矩阵是二维数组，而张量可以是三维或更高维。例如，一个三阶（或三模）张量 T 可以存在于空间 R^(N×N×N) 中，就像一个数字立方体。当然，各维度的尺寸可以不同，如 R^(M×N×P)。

如何理解张量？

理解张量需要像理解矩阵一样，超越其数字阵列的表象，看到其背后与向量空间相关的本质。

• 矩阵可以表示：
  • 双线性形式：M(x, y) = x^T M y，输入两个向量，输出一个标量。
  • 线性变换：y = M x，输入一个向量，输出另一个向量。
    矩阵的具体数值依赖于所选的基，但它表示的线性变换或双线性形式是独立于基的。
• 三阶张量可以类似地表示：
  • 三线性形式：T(x, y, z) = ∑_{i,j,k} T_{ijk} x_i y_j z_k，输入三个向量，输出一个标量。
  • 输入两个向量，输出一个向量。
  • 输入一个向量，输出一个矩阵（双线性形式）。

直观上，对张量应用矩阵乘法，相当于在其某一个“模式”（维度）上进行线性组合，类似于用矩阵左乘或右乘来组合矩阵的行或列。

张量的应用场景 🎲

张量自然出现在许多场合，一个重要的例子是描述分布的高阶矩。

• 一阶矩：均值 μ = E[x]，是一个向量。
• 二阶矩：协方差矩阵 Cov = E[(x-μ)(x-μ)^T]，描述了变量对之间的相关性。
• 三阶矩：三阶中心矩张量 M3，其元素为 M3_{ijk} = E[(x_i-μ_i)(x_j-μ_j)(x_k-μ_k)]。它自然地组织成一个三阶张量，用于描述变量三元组之间的相关性。

因此，张量可以被看作是编码随机变量高阶相关性的一种方式。当我们对随机向量进行线性变换时，其高阶矩张量也会发生相应的协调变换。

---

本节课中我们一起学习了鲁棒均值估计算法的完整推导，包括稳定性概念、关键引理以及利用凸优化（椭球法）寻找稳定子集的核心思想。随后，我们介绍了张量的基本概念，将其理解为高维数组以及作用于向量空间上的多重线性映射，并举例说明了其在表示高阶矩中的应用。下一讲我们将深入探讨张量的具体算法。

8：张量分解与独立成分分析

在本节课中，我们将学习线性代数算法的最后一部分内容：张量。我们将探讨张量的基本概念、秩的定义、一种特殊的分解形式，以及一个名为Jennrich算法的张量分解方法。最后，我们将看到该算法在独立成分分析问题中的一个应用。

张量基础与秩

上一节我们介绍了矩阵，本节我们来看看更高维度的数组——张量。

一个三维张量是一个三维的数字数组，它有三个“模式”。更高维的张量，例如四维张量，则有四个模式。就像矩阵可以视为双变量多项式、线性变换或双线性形式一样，张量可以自然地视为多变量多项式。例如，一个关于变量X, Y, Z的多线性多项式可以表示为：
[
\sum_{i,j,k} T_{ijk} X_i Y_j Z_k
]
其中 (T_{ijk}) 构成了一个三维张量。张量也常出现在数据中，例如，从分布中采样的数据，其均值是向量，协方差是矩阵，而三阶矩则对应一个三维张量。

对于矩阵，一个核心概念是秩。我们可以先定义秩为1的矩阵。

一个秩为1的矩阵可以写成一个列向量与一个行向量的外积：
[
M = u v^T
]
其第 (i, j) 项为 (u_i v_j)。基于此，矩阵 (M) 的秩 (k) 定义为能将其表示为 (k) 个秩为1矩阵之和的最小数目。

我们希望将秩的概念推广到张量。首先定义向量的张量积。

对于向量 (u, v, w \in \mathbb{R}^n)，其张量积 (u \otimes v \otimes w) 是一个三维数组，其第 (i, j, k) 项为 (u_i v_j w_k)。二维版本 (u \otimes v) 即等同于矩阵 (u v^T)。

那么，一个秩为1的张量就是三个向量的张量积 (u \otimes v \otimes w)。基于此，张量 (T) 的秩定义为能将其表示为秩为1张量之和的最小数目：
[
T = \sum_{i=1}^{k} u_i \otimes v_i \otimes w_i
]

然而，张量的秩性质不如矩阵的秩良好。

首先，一个 (n \times n \times n) 张量的最大秩可以达到 (O(n^2)) 量级，这远高于矩阵的秩（至多为 (n)）。其次，张量的秩依赖于所使用的数域（实数域或复数域）。更糟糕的是，对于矩阵，一系列低秩矩阵的极限（如果存在）仍是低秩的；但对于张量，存在高秩张量可以被一系列低秩张量任意逼近。这意味着张量的秩概念不够鲁棒。最后，计算张量的秩是NP难问题，甚至我们难以显式构造出一个秩显著高于 (n) 的张量。

尽管秩的性质不佳，但某些特殊形式的张量分解具有良好的性质且易于算法处理。

特征张量分解

在矩阵中，对称矩阵可以进行特征分解。对于张量，我们考虑一种类似的特例。

假设一个张量 (T) 可以写成若干个正交向量 (a_i) 的三重张量积之和：
[
T = \sum_{i=1}^{n} a_i \otimes a_i \otimes a_i
]
其中向量 (a_i) 彼此正交。由于正交性，这样的项数最多为 (n)，因此该张量是低秩的。我们的目标是恢复这些向量 (a_i)。

以下是一种恢复方法。

首先，选取一个随机向量 (g)，并将其应用到张量 (T) 的某一个模式上。这相当于沿着该模式对张量的切片进行线性组合。具体操作后，我们得到一个矩阵 (M_g)：
[
M_g = \sum_{i=1}^{n} (a_i \cdot g) (a_i \otimes a_i)
]
可以验证，矩阵 (M_g) 的特征向量正是这些 (a_i)，对应的特征值为 (a_i \cdot g)。由于 (g) 是随机选取的，这些特征值几乎必然互不相同，因此我们可以通过计算矩阵 (M_g) 的特征分解来唯一地恢复出所有的 (a_i)。

此外，类似于矩阵，对于由正交向量构成的张量 (T)，其对应的多项式 (T(x, x, x)) 在单位球面上的局部极大值点正是这些向量 (\pm a_i)。因此，也可以使用梯度下降等优化方法来恢复它们。

总之，对于具有这种特殊正交分解形式的张量，我们可以高效地恢复其分量。

Jennrich 算法

上一节我们看了一种特殊的分解，本节我们来看一种更一般的张量分解算法——Jennrich算法。

假设一个三阶张量 (T) 具有如下分解形式：
[
T = \sum_{i=1}^{r} u_i \otimes v_i \otimes w_i
]
其中，向量组 ({v_1, ..., v_r}) 和 ({w_1, ..., w_r}) 都是线性无关的，并且 (u_i) 是互不相同的（或满足更弱的条件）。我们的目标是恢复所有的 (u_i, v_i, w_i)。

Jennrich算法步骤如下。

1. 选取两个随机向量 (g) 和 (h)。
2. 将 (g) 和 (h) 分别应用到张量 (T) 的第一个模式上，得到两个矩阵：
   [
   M_g = \sum_{i=1}^{r} (u_i \cdot g) (v_i \otimes w_i)
   ]
   [
   M_h = \sum_{i=1}^{r} (u_i \cdot h) (v_i \otimes w_i)
   ]
3. 将这两个矩阵写成矩阵分解形式。令矩阵 (V) 的列为 (v_i)，矩阵 (W) 的行为 (w_i)，(D_g) 和 (D_h) 为对角矩阵，其对角线元素分别为 (u_i \cdot g) 和 (u_i \cdot h)。则有：
   [
   M_g = V D_g W
   ]
   [
   M_h = V D_h W
   ]
4. 计算矩阵 (M_g M_h^{-1})（假设逆存在）：
   [
   M_g M_h^{-1} = V D_g D_h^{-1} V^{-1}
   ]
   这是一个特征分解的形式。因此，矩阵 (M_g M_h^{-1}) 的特征向量就是 (V) 的列，即向量 (v_i)。
5. 类似地，计算 (M_g^T (M_h^T)) 的特征向量，可以恢复出 (W) 的行，即向量 (w_i)。
6. 一旦恢复了 (V) 和 (W)，可以通过求解线性方程组来恢复向量 (u_i)。

该算法要求 (r \leq n)，以确保 (V) 和 (W) 可逆。其核心思想是，通过在不同方向施加随机线性组合，我们得到了两个共享相同特征向量（但特征值不同）的矩阵，从而可以恢复出分解中的分量。

应用：独立成分分析

现在，我们来看Jennrich算法的一个应用：独立成分分析。

问题描述：我们获得观测样本 (y)，其生成方式为：
[
y = A x + b
]
其中：

• (x) 是一个随机向量，其各分量是独立的（例如，每个分量是随机的 ±1）。
• (A) 是一个未知的可逆方阵。
• (b) 是一个未知的偏移向量。
  目标是从观测样本 ({y}) 中恢复 (A) 和 (b)，进而可以恢复出源信号 (x)。

预处理：

1. 中心化：计算样本均值并减去，这可以消除偏移 (b) 的影响。处理后，我们可假设 (b=0)，模型简化为 (y = A x)。
2. 白化（尝试）：计算样本协方差矩阵 (\mathbb{E}[y y^T])。由于 (x) 分量独立且均值为0，有 (\mathbb{E}[x x^T] = I)。因此：
   [
   \mathbb{E}[y y^T] = A A^T
   ]
   协方差矩阵只给出了 (A A^T)，无法唯一确定 (A)，因为对于任意正交矩阵 (R)，有 ((AR)(AR)^T = A A^T)。我们遇到了旋转模糊性问题。

利用高阶统计量：
为了消除旋转模糊性，我们需要利用更高阶的统计信息。考虑四阶矩张量。
定义以下四阶张量 (M_4)，其分量计算如下：
[
M_4 = \mathbb{E}[y \otimes y \otimes y \otimes y] - T
]
其中 (T) 是一个由二阶矩（协方差）构成的张量，其 ((a,b,c,d)) 分量为所有配对乘积之和：(\mathbb{E}[y_a y_b]\mathbb{E}[y_c y_d] + \mathbb{E}[y_a y_c]\mathbb{E}[y_b y_d] + \mathbb{E}[y_a y_d]\mathbb{E}[y_b y_c])。

经过推导（此处省略详细推导），(M_4) 具有非常简洁的形式：
[
M_4 = \sum_{i=1}^{n} \kappa_i (a_i \otimes a_i \otimes a_i \otimes a_i)
]
其中 (a_i) 是矩阵 (A) 的第 (i) 列，(\kappa_i = \mathbb{E}[x_i^4] - 3) 是一个常数（对于 (x_i = \pm 1)，(\kappa_i = -2)）。

关键点：现在，我们得到了一个四阶张量，它是 (n) 个秩为1的张量（(a_i^{\otimes 4})）的加权和。这正是我们可以应用张量分解算法的形式。

应用算法：
我们可以使用类似Jennrich算法的方法来恢复向量 (a_i)。例如，可以将四阶张量通过随机投影“压缩”成一个三阶张量，然后应用Jennrich算法。一旦恢复了所有的 (a_i)（可能相差一个排列顺序），我们就成功恢复了混合矩阵 (A)。

算法失效情形：
该算法要求 (\kappa_i \neq 0)。如果 (x) 的分量是高斯分布，则 (\mathbb{E}[x_i^4] = 3)，导致 (\kappa_i = 0)，此时 (M_4) 为零张量，算法失效。这符合直觉：高斯分布是球对称的，其任意正交变换后的分布保持不变，因此从观测中无法唯一确定 (A)，旋转模糊性无法被打破。

总结

本节课中我们一起学习了张量的基本概念。我们首先看到张量的秩定义虽然自然，但性质远不如矩阵的秩良好。接着，我们学习了一种特殊的张量分解形式，当分量向量正交时，可以高效恢复。然后，我们介绍了更一般的Jennrich算法，用于分解分量向量线性无关的张量。最后，我们将该算法应用于独立成分分析问题，展示了如何利用四阶矩张量来打破二阶统计量带来的旋转模糊性，从而恢复混合矩阵。这体现了高阶张量在处理盲源分离等问题中的强大能力。

9：降维与子空间嵌入

在本节课中，我们将要学习一种强大的技术：嵌入。其核心思想是，将一个难以处理的高维或复杂领域的问题，映射到一个更容易解决的（例如更低维的）问题空间。我们将首先探讨降维，特别是随机投影，然后介绍一种更强的概念——子空间嵌入，并展示其在加速线性回归等算法中的应用。

随机投影与一维嵌入

上一节我们介绍了嵌入的广泛概念。本节中，我们来看看一个具体的例子：降维。降维是指，你有一组 D 维空间中的点，你希望将它们映射到一个更低维的子空间（例如 K 维，K 很小）。当然，在这个过程中，你允许使用随机化，并且希望在映射时保留点集的某些性质。

让我们先看看最极端的情况：一维嵌入。你想将一个从 D 维到一维的映射。这个映射是随机的。在一维空间中，你无法保留点集的太多结构，但通过一维投影仍然可以实现很多功能，并且这将成为我们构建更高维嵌入的基础模块。

最简单的尝试是随机投影。随机投影的做法是：选择一个随机方向 g，然后将所有点投影到这个方向上。定义映射 G(x) = g · x，即 x 与 g 的点积。几何上，就是选取一个随机方向，并将所有点投影到该方向上。

我们需要决定如何选取这个随机方向 g。一个简洁的方法是将其视为一个高斯向量。假设 g 是一个坐标独立的标准高斯随机向量。

这个随机映射有多好呢？它已经具有一些非常好的性质。

性质：如上构造的随机 g 是平方长度的无偏估计量。这意味着，对于 R^D 中的任意向量 x，有：

E[ (g · x)^2 ] = ||x||^2

证明：展开期望，利用 g 的分量独立且方差为 1 的性质即可得证。

由于映射 G(x) = g · x 是线性的，一旦它保留了长度，也就保留了距离。对于任意两点 x 和 y，有：

E[ (G(x) - G(y))^2 ] = E[ (g · (x - y))^2 ] = ||x - y||^2

因此，这个从 D 维到一维的随机映射，在平均意义上保持了点对之间的平方距离。当然，由于是从高维压缩到一维，很多距离信息会被破坏（例如，一个 D-1 维的子空间会被投影到零点），但平均来看，它具有上述良好性质。

此外，对于固定的 x，G(x) 本身是一个高斯随机变量，其均值为 0，方差为 ||x||^2。

提升至 K 维：Johnson-Lindenstrauss 引理

一维映射在平均意义上是正确的，但方差太大。改进估计量的自然方法是进行多次独立实验并取平均。

现在，我们定义从 D 维到 K 维的嵌入。思路是：独立进行 K 次一维投影。选取 K 个独立的高斯随机向量 g1, g2, ..., gK。定义映射 G(x) 为：

G(x) = (1/√K) * [ g1·x, g2·x, ..., gK·x ]^T

这里的归一化因子 1/√K 是为了方便。G(x) 的模长平方是 K 个独立估计量的平均值。根据期望的线性，其期望仍然是 ||x||^2。由于独立性，这个平均值的方差是单个估计量方差的 1/K，因此 ||G(x)||^2 会更集中地分布在 ||x||^2 附近。

这就是著名的 Johnson-Lindenstrauss (JL) 引理的核心思想。

JL 引理（简化表述）：固定 ε ∈ (0, 1/2) 和任意向量 x ∈ R^D。对于上述定义的随机映射 G，若选取维度 K = O( log(1/δ) / ε^2 )，则有：

Pr[ (1-ε)||x||^2 ≤ ||G(x)||^2 ≤ (1+ε)||x||^2 ] ≥ 1 - δ

这意味着，以至少 1-δ 的概率，投影后的长度在真实长度的 (1±ε) 因子范围内。

证明思路：通过归一化，可设 ||x||=1。此时，||G(x)||^2 是 K 个独立的标准高斯随机变量平方的平均值。我们关心这个平均值偏离其期望值 1 的概率。利用矩生成函数和马尔可夫不等式，可以证明其尾概率上界约为 exp(-ε^2 K / 8)。通过设置 K，即可控制失败概率 δ。

JL 引理的一个强大应用是：对于 n 个点的集合，若想以高概率（如 1-δ）同时保持所有 n^2 对点之间的距离，只需将每个点对的距离视为一个向量 xi - xj，并应用联合界。所需维度 K 仅为 O( log(n/δ) / ε^2 )。这意味着，任何 n 个点的集合都可以嵌入到 O(log n) 维的空间中，并近似保持所有欧氏距离。

oblivious 子空间嵌入

JL 映射是数据无关的，它只是随机选取一个线性变换。有时，查看数据并选取映射（如 PCA）可能更优，但 PCA 计算昂贵且是数据相关的。JL 映射可以增量计算，在流式计算中非常有用。

子空间嵌入是比 JL 引理更强的概念。假设给定一个由 n×d 矩阵 A 的列张成的 d 维子空间 S ⊆ R^n。目标是构造一个映射 G: R^n → R^k，使得对于子空间 S 中的每一个向量 x，都有：

(1-ε)||x|| ≤ ||G(x)|| ≤ (1+ε)||x||

这相当于要保护子空间内所有（无限多个）点对之间的距离。这种嵌入称为 oblivious 子空间嵌入，因为映射 G 的构造不依赖于具体的子空间 S，但对任意给定的 S 都能以高概率成功。

定理：若 G 的每一行是独立的标准高斯向量（即 G 是 k×n 的高斯随机矩阵），则只要 k = O( (d + log(1/δ)) / ε^2 )，以至少 1-δ 的概率，上述性质对任意 d 维子空间 S 成立。

这个定理的关键在于，所需维度 k 仅依赖于子空间本身的维度 d，而与原始环境维度 n 无关。

应用：加速线性回归

子空间嵌入的一个直接应用是加速大规模线性回归。线性回归问题为：给定矩阵 A ∈ R^(n×d) 和向量 b ∈ R^n，求解：

min_x ||Ax - b||^2

当 n 非常大（数据点很多）而 d 较小（特征数少）时，直接求解计算量可能很大。

利用 oblivious 子空间嵌入 G，我们可以将原问题转化为一个更小的问题：

min_x || (GA)x - (Gb) ||^2

这里，GA ∈ R^(k×d)，Gb ∈ R^k，且 k = O(d/ε^2) 远小于 n。在小矩阵上求解回归得到解 x̂。

为什么有效？ 令 S 为矩阵 [A, b] 的列张成的 (d+1) 维子空间。对于该子空间中的任意向量（例如 Ax - b），嵌入 G 能以高概率保持其长度。设 x* 是原问题的最优解。因为 x̂ 是小规模问题的最优解，所以有：

||G(Ax̂ - b)||^2 ≤ ||G(Ax* - b)||^2

利用子空间嵌入的性质，可以将两边的 G 近似“去掉”，得到：

||Ax̂ - b||^2 ≤ (1+O(ε)) ||Ax* - b||^2

因此，x̂ 是原问题的一个 (1+ε) 近似解。这样，我们仅用 O(nd + poly(d/ε)) 的时间（主要花在计算 GA 和 Gb 上）就近似解决了问题，而传统方法可能需要 O(nd^2) 时间。

定理证明思路：ε-网

如何证明 oblivious 子空间嵌入定理？核心挑战在于要处理子空间中的无限多个点。技巧是使用 ε-网。

策略：

1. 在子空间的单位球面上，构造一个有限的点集 N（称为 γ-网），使得球面上任意一点，都存在网中某点与其距离不超过 γ。
2. 证明对于网 N 中的所有点，映射 G 能以高概率保持长度。
3. 论证如果 G 对网中所有点都保持长度，那么通过线性性和连续性，它对整个单位球面（从而对整个子空间）的所有点也都保持长度。

关键点：

• 一个 d 维空间中的单位球面，可以构造一个大小为 O((1/γ)^d) 的 γ-网。这是指数于维度 d 的，但仍然是有限的。
• 对于网中这 M = O((1/γ)^d) 个点，我们可以应用联合界。为了保证每个点长度失真的失败概率小于 δ/M，根据 JL 引理，需要 k = O( log(M/δ) / ε^2 ) = O( (d + log(1/δ)) / ε^2 )。
• 通过精心选择 γ（与 ε 相关），并利用 G 的线性性和利普希茨性质，可以将网上的结论推广到整个连续空间。

这就完成了定理的证明。

快速 JL 变换

原始的 JL 映射使用稠密高斯矩阵，计算 Gx 需要 O(kd) 时间。研究者们提出了更快的结构化或稀疏矩阵构造，称为快速 JL 变换。

一个例子是以下三步构造：

1. 随机符号翻转：将 x 的每个分量独立乘以 +1 或 -1。
2. 快速傅里叶变换 (FFT)：对结果应用 FFT。
3. 均匀采样：从 FFT 输出的 n 个系数中，均匀随机选取 k 个。

这个变换可以快速计算（O(n log n) 时间），并且被证明能以高概率满足 JL 引理的性质。此外，还有基于计数草图等技术的稀疏 JL 构造。

总结

本节课中我们一起学习了：

1. 随机投影：作为一维嵌入，它是平方长度的无偏估计。
2. Johnson-Lindenstrauss 引理：通过将一维投影平均到 K 维，可以高概率地将 n 个点嵌入到 O(log n) 维空间，并近似保持所有 pairwise 距离。
3. Oblivious 子空间嵌入：一种更强的嵌入，能保持整个 d 维子空间中的所有长度，所需维度仅 O(d/ε^2)，且与原始空间维度 n 无关。
4. 应用：我们看到了子空间嵌入如何用于大幅加速大规模线性回归，仅需处理压缩后的数据。
5. 证明技术：通过构造 ε-网，将无限连续空间的问题转化为有限点集的问题，再利用联合界和 JL 引理进行证明。
6. 快速实现：介绍了通过随机符号翻转、FFT 和采样实现的快速 JL 变换。

这些降维和嵌入技术是算法设计，特别是流计算、 sketching 和数值线性代数中的基础工具，具有广泛的应用。

10：压缩感知与稀疏恢复

在本节课中，我们将完成对“无意识子空间嵌入”定理的证明，并深入探讨一种被称为“压缩感知”的特定维度约减应用。我们将学习如何利用信号的稀疏性，通过远少于信号维度的测量来精确恢复原始信号。

完成无意识子空间嵌入的证明

上一节我们介绍了无意识子空间嵌入的核心思想，并开始证明一个关键定理：对于一个随机高斯矩阵 G ∈ ℝ^{t×n}，如果行数 t 足够大，那么对于任意一个 d 维子空间，G 能以高概率近似保持该子空间内所有向量的长度。

我们通过选取子空间单位球面上的一个 γ-网（一个有限点集，使得球面上任意点都离网中某点很近）来将无限问题转化为有限问题。我们已经证明了对于网中的所有点，G 能以高概率保持其内积（从而保持长度）。现在，我们需要证明这个性质能扩展到整个连续的子空间球面。

从网扩展到整个球面

核心思路：任意球面上的点 x 都可以表示为网中点的一个线性组合，且系数呈几何级数衰减。

构造过程：

1. 令 x₁ = x。
2. 在网中找到离 x₁ 最近的点 y₁。根据 γ-网的定义，有 ||x₁ - y₁|| ≤ γ。
3. 将 x₁ 分解为：x₁ = y₁ + ||x₁ - y₁|| * ( (x₁ - y₁) / ||x₁ - y₁|| )。令 x₂ = (x₁ - y₁) / ||x₁ - y₁||（这是一个单位向量）。
4. 对 x₂ 重复此过程：找到网中最近点 y₂，分解 x₂ = y₂ + ||x₂ - y₂|| * x₃。
5. 持续迭代，我们可以将 x 写成一个无限级数：x = Σ_{i=1}^{∞} α_i y_i，其中系数 α_i 满足 |α_i| ≤ γ^{i-1。

现在，考虑变换后的向量 Gx 的长度：
||Gx||² = <Gx, Gx> = Σ_{i,j} α_i α_j <Gy_i, Gy_j>

由于 y_i 和 y_j 都是网中的点，根据我们的假设，它们的内积被很好地保持：<Gy_i, Gy_j> = <y_i, y_j> ± ε。
因此：
||Gx||² = Σ_{i,j} α_i α_j <y_i, y_j> ± ε * Σ_{i,j} |α_i α_j|

第一项正是 ||x||² = <x, x>。第二项是误差项。由于系数 α_i 几何衰减，其绝对值和 Σ_i |α_i| 是一个有界常数（例如，小于 2）。因此，误差项被控制在 O(ε) 级别。

这就证明了对于子空间球面上的任意点 x，其长度在经过 G 变换后都能在 (1 ± O(ε)) 因子内得到保持。定理得证。

总结：通过使用 γ-网，我们将对无限连续集合的联合界问题，转化为对有限点集（其大小与子空间维度 d 呈指数关系）的联合界问题。所需的测量数 t 与 d + log(1/δ) 成正比，这解释了定理中的参数。

引入压缩感知 🎯

无意识子空间嵌入处理的是子空间，而压缩感知处理的是另一类结构：稀疏向量。许多自然信号（如图像、音频）在某个合适的基（如傅里叶基、小波基）下是稀疏的，即只有少数坐标值显著非零。

核心问题：能否利用这种稀疏性，仅通过 m << n 次线性测量（远少于信号维度 n）来完整恢复一个 n 维稀疏信号 x？

线性测量模型

我们通过一个测量矩阵 M ∈ ℝ^{m×n} 进行测量。我们不是直接观测 x，而是观测线性测量值 b = Mx ∈ ℝ^m。
恢复任务：在已知 M 和测量结果 b 的前提下，恢复出原始的稀疏信号 x。

这里存在两个关键问题：

1. 可识别性：矩阵 M 需要具备什么性质，才能保证稀疏解 x 是唯一确定的？
2. 算法：如何从 b 中有效地恢复出 x？

ℓ₁ 最小化：一种启发式方法

最直接的稀疏恢复优化问题是最小化 x 的 ℓ₀ “范数”（即非零元素个数）：

minimize ||x||₀ subject to Mx = b

然而，这是一个非凸、NP难的问题。

一个广泛使用的启发式方法是将其松弛为凸优化问题，即最小化 ℓ₁ 范数：

minimize ||x||₁ subject to Mx = b

其中 ||x||₁ = Σ_i |x_i|。

直观解释：ℓ₁ 范数球（例如三维中的八面体）是“尖”的，其顶点和低维面对应于稀疏向量（顶点是只有一个非零元素的向量，棱边是有两个非零元素的向量等）。最小化 ℓ₁ 范数相当于在约束超平面 Mx = b 上寻找与 ℓ₁ 球面首次接触的点。对于“大多数”超平面，首次接触点很可能落在这些“尖”的顶点或棱边上，从而产生一个稀疏解。

受限等距性质与精确恢复

为了理论证明 ℓ₁ 最小化的有效性，我们引入一个关键概念。

定义（受限等距性质，RIP）：一个矩阵 M ∈ ℝ^{m×n} 满足 (k, ε)-RIP，如果对于所有至多 k 个非零元素的向量 x（k-稀疏向量），都有：
(1 - ε) ||x||₂² ≤ ||Mx||₂² ≤ (1 + ε) ||x||₂²

这意味着 M 近似保持所有 k-稀疏向量的长度。

为什么随机矩阵满足 RIP？
所有 k-稀疏向量的集合可以表示为 C(n, k) 个 k 维子空间的并集。我们已经知道，对于一个固定的 k 维子空间，随机高斯矩阵能以高概率（通过调整行数 m）保持其上的长度。通过对这 C(n, k) 个子空间进行联合界，我们可以证明：如果 m = O( k + log(C(n, k)) ) = O( k log(n/k) )，那么随机高斯矩阵 M 将以高概率满足 (k, ε)-RIP。这从信息论角度看也是直观的：我们需要约 log(C(n, k)) 比特来编码非零值的位置，以及 O(k) 次测量来编码这些位置上的值。

ℓ₁ 最小化精确恢复定理

定理：如果测量矩阵 M 满足 (2k, ε)-RIP，且 ε 足够小（例如小于某个常数如 0.4），那么对于任意 k-稀疏信号 x，通过测量 b = Mx 并求解 ℓ₁ 最小化问题：

minimize ||x'||₁ subject to Mx' = b

得到的解 x' 将精确等于原始信号 x。

证明思路（反证法）：

1. 假设存在另一个解 x' ≠ x，且 ||x'||₁ ≤ ||x||₁。令 Δ = x' - x，则 MΔ = 0。
2. 设 S 是 x 的非零坐标索引集，|S| = k。由于 x' 的 ℓ₁ 范数更小，可以推导出 Δ 在 S 上的 ℓ₁ 质量至少占其总 ℓ₁ 质量的一半：||Δ_S||₁ ≥ 0.5 ||Δ||₁。
3. 将 Δ 在 **S补集上的坐标按绝对值降序排列，并分成若干块，每块大小约为k`。
4. 关键步骤：
   • 向量 [Δ_S; Δ_{第一块}] 是 3k-稀疏的。由于 M 满足 (2k, ε)-RIP（实际上需要 3k），其范数被大致保持，因此 ||M(Δ_S + Δ_{第一块})||₂ 相对较大。
   • 利用坐标递减的性质，可以证明剩余各块的 ℓ₂ 范数之和很小。
5. 然而，由于 MΔ = 0，我们有 0 = MΔ_S + Σ_{所有块} MΔ_{块}。这意味着 MΔ_S 必须能被后续各块的像抵消。但根据 RIP，MΔ_S 的“能量”较大，而后续各块像的“能量”总和很小，无法将其抵消至零，从而产生矛盾。因此，假设不成立，x' 必须等于 x。

总结

本节课中我们一起学习了：

1. 完成了无意识子空间嵌入的证明：通过构造 γ-网，我们将对无限子空间的保证转化为对有限点集的联合界，证明了随机高斯矩阵能以接近最优的参数保持任意子空间的结构。
2. 引入了压缩感知框架：利用信号在特定基下的稀疏性，我们可以用远少于信号维度的线性测量来捕获信号的全部信息。
3. 分析了 ℓ₁ 最小化算法：我们将难解的 ℓ₀ 最小化问题松弛为凸的 ℓ₁ 最小化问题。通过几何直观和受限等距性质，我们证明了对于满足 RIP 的随机测量矩阵，ℓ₁ 最小化能精确恢复原始稀疏信号。

压缩感知将维度约减的思想应用于稀疏向量集，在医学成像、单像素相机等领域有重要应用。其核心在于，通过巧妙的线性测量和高效的凸优化算法，我们可以突破传统采样定理的限制。

11：近似最近邻搜索与局部敏感哈希

在本节课中，我们将学习度量空间和嵌入的基本概念，并探讨一个核心的计算问题：近似最近邻搜索。我们将看到，通过一种称为“局部敏感哈希”的巧妙技术，可以高效地解决这个问题。

度量空间与嵌入

上一节我们介绍了度量空间和嵌入的概念。本节中，我们来看看它们的正式定义。

度量空间 是一个集合 X 及其上的一个距离函数 d。距离函数 d: X × X → ℝ⁺ 满足以下性质：

• d(x, x) = 0
• d(x, y) = d(y, x) （对称性）
• d(x, y) + d(y, z) ≥ d(x, z) （三角不等式）

常见的例子包括欧几里得空间（L2范数）、曼哈顿空间（L1范数）、字符串的编辑距离以及图上的最短路径距离。

嵌入 是从一个度量空间 (X, d) 到另一个度量空间 (X‘, d’) 的映射函数 φ。如果该映射能以因子 D 保留距离，则称为失真度为 D 的嵌入。即，对于所有 x, y ∈ X，满足：

d(x, y) / D ≤ d‘(φ(x), φ(y)) ≤ d(x, y)

失真度为1的嵌入称为等距嵌入。

最近邻搜索问题

现在，我们来看一个度量空间上的基础计算问题：最近邻搜索。

给定一个度量空间（例如 ℝᵈ 和 L2 范数）中的一个点集 P（数据点），以及一个查询点 q，目标是找到 P 中距离 q 最近的点。

一个朴素的解决方案是存储所有点，查询时计算 q 到每个点的距离。这需要 O(nd) 的存储空间和 O(nd) 的查询时间（d 为维度）。在低维度下，可以使用 KD 树等数据结构，但其空间复杂度通常随维度 d 指数增长，对于高维数据（如图像、文本向量）不适用。即使先使用约翰逊-林登斯特劳斯引理进行降维，再构建数据结构，空间需求仍然很大。

因此，我们将关注一个放松版本的问题。

近似最近邻问题

为了解决高维下的效率问题，我们考虑一个近似且带尺度约束的版本：C-近似 R-近邻 问题。

输入：

• 点集 P（数据点）。
• 距离尺度 R。
• 近似因子 C > 1。

查询：

• 一个查询点 q。

目标：

• 如果存在一个数据点 p ∈ P 满足 d(p, q) ≤ R，则算法必须返回一个点 p’ ∈ P，满足 d(p’, q) ≤ C * R。
• 如果不存在这样的点（即所有点距离都大于 R），算法可以不返回任何结果或返回任意点。

一个更强的版本是 C-近似最近邻 问题，它要求返回的点距离 q 至多是最近点距离的 C 倍。通常，可以通过在尺度 R 上二分搜索，利用解决前一个问题的算法来解决后一个问题。

局部敏感哈希

解决近似近邻搜索的一个核心思想是使用哈希，但不是普通的哈希。我们需要一种能反映数据点之间“距离”的哈希函数，即局部敏感哈希。

一个 (R, CR, P_close, P_far) 局部敏感哈希族 是一族哈希函数 H，满足对于任意两点 x, y：

• 若 d(x, y) ≤ R，则 Pr[h(x) = h(y)] ≥ P_close。
• 若 d(x, y) ≥ CR，则 Pr[h(x) = h(y)] ≤ P_far。

其中 P_close > P_far。哈希函数 h 是从该族中随机选取的。

哈希族构造示例

以下是几个LSH族的例子：

1. 汉明距离（L1 on {0,1}ᵈ）

• 哈希函数：随机选取一个坐标 i，h(x) = x_i。
• 性质：Pr[h(x) = h(y)] = 1 - (汉明距离(x, y) / d)。

2. 单位球面上的余弦距离（角度）

• 哈希函数：随机选取一个超平面（通过随机高斯向量 v 定义），h(x) = sign(v·x)。
• 性质：两点被超平面分开的概率正比于它们之间的角度。

3. 杰卡德相似度（用于集合/文档）

• 距离：d(A, B) = 1 - |A ∩ B| / |A ∪ B|。
• 哈希函数（最小哈希）：随机排列整个字典（词表），h(A) 定义为排列后 A 中第一个元素。
• 关键性质：Pr[h(A) = h(B)] = |A ∩ B| / |A ∪ B| = 1 - d(A, B)。这个相关的随机选择（同一排列用于所有文档）至关重要。

4. L2 距离（欧几里得距离）

• 构造分两步：
  1. 随机投影：选取随机高斯向量 v，将点 x 映射到一维实数 v·x。期望距离得以保留。
  2. 离散化：在一维线上随机设置一个阈值 θ，h(x) = 1 if (v·x ≥ θ) else 0。
• 最终，两点哈希值不同的概率正比于它们的 L2 距离。

从LSH到数据结构

仅仅有LSH族还不够，我们需要利用它来构建高效的数据结构。思路是放大 P_close 和 P_far 之间的差距。

1. 放大差距
构造新的哈希函数 g(x) = (h₁(x), h₂(x), ..., h_t(x))，即串联 t 个独立的LSH函数。那么：

• 若 d(x, y) ≤ R，则 Pr[g(x) = g(y)] ≥ (P_close)^t。
• 若 d(x, y) ≥ CR，则 Pr[g(x) = g(y)] ≤ (P_far)^t。
  通过增大 t，可以使相近点碰撞概率与相远点碰撞概率的差距指数级扩大。

2. 设置参数并重复
我们希望对于任意查询点 q，没有远点（距离 ≥ CR）与它碰撞。这要求 (P_far)^t ≈ 1/n，解得 t ≈ log(1/n) / log(P_far) = O(log n)。
此时，一个近点（距离 ≤ R）与 q 碰撞的概率至少为 (P_close)^t ≈ n^{-ρ}，其中 ρ = log(P_close) / log(P_far) < 1。

由于单次成功概率低（n^{-ρ}），我们构建 L = n^ρ 个独立的哈希表，每个哈希表使用不同的 g 函数。这样，对于一个近邻查询点，它在至少一个哈希表中与近点碰撞的概率就变成了常数。

3. 算法流程

• 预处理：构建 L = n^ρ 个哈希表。对于每个数据点 p，计算它在所有 L 个哈希表中的索引 g_i(p)，并将其存入对应的桶中。
• 查询：对于查询点 q，计算它在 L 个哈希表中的索引 g_i(q)。检查每个对应桶中是否存在数据点。如果找到，则计算其与 q 的实际距离，若满足 ≤ C*R 则返回。

4. 复杂度分析

• 空间：O(n * L) = O(n^{1+ρ})，用于存储 L 个哈希表。
• 查询时间：O(L) = O(n^ρ) 次哈希计算和桶访问，通常还需进行常数次距离计算验证。

扩展到复杂度量

对于编辑距离、地球移动距离等更复杂的度量，直接构造LSH可能困难。一个常见的方法是分两步：

1. 嵌入：将原度量空间嵌入到更简单的空间（如 L1）中，虽然会引入失真 D（可能是 log n 等增长缓慢的函数）。
2. 应用LSH：在目标空间（如 L1）上使用已知的高效LSH方案。

总结

本节课中我们一起学习了：

1. 度量空间与嵌入 的形式化定义。
2. 最近邻搜索 问题及其在高维下面临的挑战。
3. 近似最近邻 问题的放松定义（C-近似 R-近邻）。
4. 局部敏感哈希 的核心思想：哈希碰撞的概率反映数据点的距离。
5. 如何通过 串联哈希函数 和 构建多个哈希表 来放大概率差距，从而构造出解决近似近邻搜索的数据结构，其查询时间为 O(n^ρ)，其中 ρ 是LSH族的一个参数。
6. 几个经典度量（汉明距离、余弦距离、杰卡德相似度、L2距离）的LSH构造示例。
7. 处理复杂度量的一般策略：先嵌入到简单空间（如L1），再应用LSH。

局部敏感哈希是处理高维相似性搜索的强大工具，在文档去重、图像检索、推荐系统等领域有广泛应用。

12：度量空间嵌入与低直径分解 🧭

在本节课中，我们将深入学习度量空间，特别是如何将复杂的度量空间嵌入到更简单的树度量中。我们将介绍低直径分解这一核心工具，并展示如何利用它来构建高效的随机树嵌入。

---

度量空间与树嵌入 🌳

上一节我们介绍了度量空间的基本概念。本节中，我们来看看如何简化复杂的度量空间。

一个度量空间由一组点集 X 和一个距离函数 d(x, y) 定义。一个自然的例子是图上的最短路径度量：给定一个带权图，点之间的距离定义为两点间最短路径的权重之和。

处理一般的度量空间可能相当复杂。因此，我们希望将其简化。具体来说，我们今天要寻找的是嵌入，它能将一般的度量空间映射到树度量。树度量，顾名思义，就是定义在树结构上的最短路径度量，它显然要简单得多。

进行这种嵌入的优势在于，许多NP完全问题在树度量上可以高效求解。因此，解决近似完全问题的一个自然方法是：获取一个一般度量空间（或从图中得到的最短路径度量），将其嵌入到一个树度量中，然后在树度量上解决问题。这为你最初关心的问题提供了一个近似算法。事实上，许多近似算法和在线算法都使用了这种方法。

理想与现实

理想情况下，对于一个度量空间 (X, d)，我们希望找到一个嵌入到树度量空间 (T, d_T)，使得对于所有点 x, y，有 d(x, y) = d_T(x, y)。但正如你所料，这通常是不可能的。

例如，考虑一个 n 个点的环。如果你试图将其嵌入到任何一棵树中，你必然会引入 Ω(n) 的失真。这意味着某些距离会被极大地拉伸。问题的关键在于，如果只使用一棵确定的树，可能会产生如此大的失真。

随机树嵌入

因此，我们实际感兴趣的是概率树嵌入，即寻找一个随机的或概率性的树嵌入。

给定度量空间 (X, d)，一个随机树嵌入定义为一个树度量上的分布。我们约定，对于从这个分布中采样的每一棵树，距离都不会减少（即 d_T(x, y) ≥ d(x, y)），但也不会增加太多。具体来说，我们希望对于所有点对 (x, y)，树距离的期望值最多是原始距离的某个因子 α 倍：

E[d_T(x, y)] ≤ α * d(x, y)

有时，如果原始度量空间来自一个图 G=(V, E)，我们可能希望嵌入树是原图的生成树，这被称为低拉伸生成树。今天我们主要关注更一般的随机树嵌入。

一个重要结论是：对于任意包含 n 个点的度量空间 (X, d)，都存在一个树度量分布，其失真 α 为 O(log n)。这意味着我们可以在对数因子内用树来近似任何度量空间。

---

低直径分解 🔍

构建树嵌入的一个关键工具是低直径分解。这个概念本身在属性测试、分布式算法和许多近似算法中都非常有用。

问题定义

给定一个度量空间 (X, d)，我们的目标是将它分割成若干小块（即一个划分 P = {X₁, X₂, ...}），使得每一块的直径（块内任意两点间的最大距离）至多为某个小值 Δ。

这很容易做到，例如将每个点单独作为一块。但我们真正关心的是一个权衡：在将空间分割成小块的同时，我们不希望“切碎”得太过分。具体来说，如果两个点 x 和 y 原本很接近，我们希望它们有很高的概率落在同一块中。

我们希望分离两个点的概率与它们之间的距离成正比。用 P(x) 表示点 x 所属的划分块。我们希望存在一个参数 β，使得对于任意两点 x, y：

Pr[P(x) ≠ P(y)] ≤ β * (d(x, y) / Δ)

分母中的 Δ 是自然的：如果你要求更小的直径（更小的 Δ），那么分离邻近点的概率就会相应增加。

直观例子：实数线

考虑实数线 R，其度量为通常的绝对值距离。假设我们希望将实数线分割成直径为1的块。

一个自然的确定性方法是：在整数坐标处进行划分，例如 [0,1), [1,2), ...。但这样会有一个问题：如果两个点 x 和 y 恰好位于某个划分边界的两侧（例如 0.999 和 1.001），即使它们距离非常近（0.002），也总是会被分离。这违反了我们的概率上界。

一个简单的改进方法是：随机平移划分的起点。例如，随机均匀地选择一个偏移量 s ∈ [0, 1)，然后在 s, s+1, s+2, ... 处进行划分。这样，两个距离为 r 的点被分离的概率恰好是 r（当 r ≤ 1 时）。这里我们得到了 β = 1 的理想情况。

通用算法

对于一般的度量空间，我们有一个非常简洁的算法。目标是生成一个直径至多为 Δ 的划分。

以下是算法步骤：

1. 随机选择一个半径 R，其值在 [Δ/4, Δ/2] 区间内均匀分布。
2. 将点集 X 随机排列（生成一个随机排列）。
3. 按照随机排列的顺序处理每个点 x_i：
   • 考虑以 x_i 为中心、半径为 R 的“球” B(x_i, R)，即所有与 x_i 距离不超过 R 的点。
   • 将 B(x_i, R) 中尚未被分配到任何块的所有点，分配到一个新的块中（这个块以 x_i 为“中心”）。

这个算法直观上类似于在空间中随机撒点并“吞噬”周围区域的过程。随机半径 R 的引入是为了避免度量空间在某个固定半径附近出现不良行为。

算法分析（核心结论）

该算法可以保证，对于任意两点 x, y，它们被分离的概率满足：

Pr[P(x) ≠ P(y)] ≤ O(log n) * (d(x, y) / Δ)

其中 n 是总点数。这意味着我们得到了一个 β = O(log n) 的低直径分解。这个对数因子是不可避免的，并且与我们最终树嵌入的失真因子相匹配。

---

层次化树分解与嵌入 🌲

现在，我们回到最初的目标：如何构建一个随机树嵌入？我们将使用一种递归的策略，称为层次化树分解。

构建过程

1. 递归分解：从整个度量空间开始，我们应用低直径分解算法，将其划分为直径约为 Δ 的块。然后，将每个块视为一个新的子度量空间，递归地对每个子块应用低直径分解，但将目标直径减半（例如设为 Δ/2）。继续这个过程，每次递归都将目标直径除以一个常数因子（如2或8），直到每个块只包含一个点。
2. 构建树结构：这个递归分解过程自然定义了一棵树。
   • 叶子节点：是原始的点。
   • 内部节点：对应递归过程中产生的每一个块（划分）。
   • 边与权重：连接一个块节点与其子块节点（即该块在下一层递归中被划分成的更小块）。我们为这条边赋予一个权重，其值等于该父块节点的直径上界（例如，在递归层 i，直径为 Δ/2^i，则权重可设为 Δ/2^i）。

为什么这是好的嵌入？

考虑任意两个点 x 和 y，设它们的原始距离 d(x, y) 介于 8^{j-1} 和 8^j 之间。

• 下界（距离不被压缩）：在递归分解中，直到高度约为 j 的层级，块的直径都小于 d(x, y)。因此，x 和 y 在高度 j 或更早的层级就必定被分离到不同的块中。这意味着在树中，x 和 y 的最近公共祖先至少位于高度 j。树中 x 到 y 的路径必须经过至少一条权重约为 8^j 的边，因此树距离 d_T(x, y) ≥ 8^j ≥ d(x, y)。所以嵌入不会减少距离。
• 上界（距离不被过度拉伸）：我们需要证明期望的树距离不会比原始距离大太多。关键观察是：x 和 y 在高度 h 的层级被分离的概率，由我们的低直径分解引理控制，大约为 O(log n) * (d(x, y) / 8^h)。树距离可以表示为，对所有可能分离它们的高度 h，权重 8^h 乘以在该高度分离的概率之和。通过计算这个和式，并利用低直径分解的概率上界，我们可以证明：
  E[d_T(x, y)] ≤ O(log n) * d(x, y)

结合上下界，我们得到了一个失真为 O(log n) 的随机树嵌入。

---

应用示例：批量购买网络设计 💻

让我们看一个随机树嵌入的应用实例：批量购买网络设计问题。

问题描述

假设你有一个图（代表通信网络），以及一系列流量需求 (s_i, t_i, f_i)，表示需要在源点 s_i 和汇点 t_i 之间支持 f_i 单位的流量。你可以在图的边上购买容量，每条边 e 上每单位容量的成本为 c_e。容量可以被所有流量需求共享（分时复用）。目标是购买总成本最低的容量配置，以满足所有流量需求。

树度量的优势

在一般的图上，这是一个复杂的优化问题。然而，在树结构上，这个问题变得非常简单。因为在树中，任意两点之间只有唯一路径。为了满足需求 (s_i, t_i, f_i)，你只需在从 s_i 到 t_i 的唯一路径上的每条边，增加 f_i 单位的容量。所有需求独立处理，然后将每条边上所有需求要求的容量相加，再乘以边成本，即得总成本。这是一个确定性的、容易计算的解。

近似算法

利用我们学到的随机树嵌入，可以设计一个近似算法：

1. 将原始图（或其最短路径度量）通过随机树嵌入，失真为 O(log n) 地映射到一棵随机树 T 上。
2. 在树 T 上求解上述批量购买网络设计问题（这很容易）。
3. 将树 T 上的解（即各边应购买的容量）解释回原图。由于嵌入的失真，这个解可能无法在原图上直接支持所有流量，但可以证明，通过适当缩放，它能给出原问题的一个 O(log n) 近似解。

这个例子展示了如何将复杂问题简化到树结构上求解，是随机树嵌入威力的一个典型体现。

---

总结 📚

本节课中我们一起学习了：

1. 随机树嵌入的概念：通过一个树度量上的分布来近似复杂的度量空间，并保证期望失真为 O(log n)。
2. 低直径分解这一核心工具：一种将度量空间划分为小块的方法，同时保证邻近点被分离的概率与其距离成正比（乘以一个 O(log n) 因子）。我们介绍了一个简洁的随机算法来实现它。
3. 层次化树分解的构造：通过递归应用低直径分解，并赋予边权重，可以系统地构建出满足要求的随机树嵌入。
4. 一个应用实例：批量购买网络设计问题，展示了如何利用树嵌入将复杂图上的难题转化为树上的易解问题，从而获得高效的近似算法。

这些工具在组合优化、在线算法和分布式计算等领域有着广泛的应用，是处理度量空间相关问题的强大基础。

13：稀疏割问题与度量嵌入

在本节课中，我们将学习一个名为“稀疏割”的问题。我们将看到如何通过将问题嵌入到 L1 空间或度量嵌入中，来紧密地联系并解决这个问题，并利用这种方法构建一个算法。

图扩展概述

图扩展是衡量图连通性好坏的一种度量。

假设我们有一个图 G，它由一组顶点和边组成。我们假设图是无权且 d-正则的，即所有顶点的度均为 D。这使得定义更简单。

考虑图中的一个顶点子集 S。子集 S 的扩展定义为穿过该子集的边数，即从 S 到其补集 S^c 的边数。

我们通常将这个数除以与 S 关联的边的总数，即 S 的体积。在 d-正则图中，S 的体积就是 d 乘以 S 的基数 |S|。因此，在 d-正则图中，扩展可以表示为：
Φ(S) = (从 S 到 S^c 的边数) / (d * |S|)

另一种理解方式是：从 S 中随机选择一个顶点，然后沿着该顶点的一条随机边移动一步，离开集合 S 的概率就是 Φ(S)。这是一个介于 0 和 1 之间的数。

图的扩展（也称为图的电导）定义为所有大小不超过 n/2 的子集 S 中，Φ(S) 的最小值。

为什么这个度量能反映连通性？让我们看几个例子。

在完全图中，任意大小小于 n/2 的子集 S，其扩展至少为 1/2。这意味着完全图是高度连通的。

在一个 n×n 的网格图中，如果我们从中间切开，穿过割的边数仅为 n，而左侧子集的体积约为 n^2，因此扩展约为 1/n。这表明网格图的连通性较差。

与最小割相比，扩展是一个更好的连通性度量。最小割可能只分离一个顶点，不能很好地反映整体连通性。

扩展图

扩展图是一种具有常数度（例如度为 3）但扩展至少为某个常数的图。这意味着尽管边数稀疏（O(n) 条边），但其连通性却像完全图一样好。

随机 d-正则图以高概率是扩展图。此外，也存在明确的扩展图构造。

以下是一个明确的扩展图构造示例：
取顶点为 {0, 1, ..., p-1}，其中 p 为素数。对于每个顶点 x，连接以下三条边：

1. x 连接到 (x+1) mod p
2. x 连接到 (x-1) mod p
3. x 连接到 (1/x) mod p（当 x≠0 时，0 可特殊处理）
   这样就得到了一个 3-正则图，并且它是一个扩展图。

扩展图在算法中有广泛应用，如路由、最短路径计算、拉普拉斯矩阵求解等。

稀疏割问题

我们关注的问题是：给定一个图 G，计算或近似其扩展值 Φ(G)。计算精确值是 NP 难的，因此我们寻求近似算法。

为了便于处理，我们引入一个与 Φ(G) 相近的量 Ψ(G)：
Ψ(G) = min_{S: |S|≤n/2} (从 S 到 S^c 的边数) / (|S| * |S^c|)

可以证明，近似 Ψ(G) 与近似 Φ(G) 在常数因子内是等价的，因为 |S^c| 总是在 n/2 和 n 之间。

我们可以将 Ψ(G) 用割度量来表示。对于子集 S，定义其指示函数 1_S(v)，当顶点 v 在 S 中时为 1，否则为 0。这定义了一个割度量 d_S(u, v) = |1_S(u) - 1_S(v)|。

利用这个度量，Ψ(G) 可以重写为：
Ψ(G) = min_S ( Σ_{(u,v)∈E} d_S(u, v) ) / ( Σ_{所有 u,v} d_S(u, v) )

我们的目标是找到一个割度量，使得这个比值最小化。

线性规划松弛

直接处理割度量很困难，因此我们考虑一个线性规划松弛。我们寻找一个满足度量条件（非负性、对称性、三角不等式）的任意距离函数 d(u, v)，来最小化相同的比值。

由于比值是齐次的，我们可以固定分母至少为 1，然后最小化分子。这给出了以下线性规划：

最小化： Σ_{(u,v)∈E} d(u, v)
约束条件：

1. d 是一个度量（满足所有度量公理）。
2. Σ_{所有 u,v} d(u, v) ≥ 1

求解这个线性规划，我们得到一个距离函数 d，它最小化了目标比值。

通过度量嵌入进行舍入

现在我们得到了一个任意的度量 d，需要将其“舍入”回一个割度量。我们的策略是先将这个度量嵌入到 L1 空间（即具有 L1 距离的欧几里得空间），然后再从 L1 嵌入中得到一个割。

关键定理是：任何 n 点度量空间都可以嵌入到 L1 空间中，且失真度仅为 O(log n)。这意味着嵌入后的距离与原距离之比最多为 O(log n) 倍。

Frechet 嵌入：基本构件

我们使用一种称为 Frechet 嵌入的技术作为构建块。给定度量空间和其一个子集 A，Frechet 嵌入将每个点 x 映射到它到集合 A 的距离：f_A(x) = d(x, A) = min_{y∈A} d(x, y)。

Frechet 嵌入是收缩的：对于任意两点 x, y，有 |f_A(x) - f_A(y)| ≤ d(x, y)。

Bourgain 嵌入构造

我们的嵌入是随机化的。我们构造 log n 个随机子集 A_1, A_2, ..., A_log n。子集 A_t 通过以概率 1/(2^t) 独立包含每个顶点来构建。

然后，我们将一个点 x 嵌入到一个 log n 维向量中：
F(x) = ( d(x, A_1), d(x, A_2), ..., d(x, A_log n) )

这个嵌入的 L1 距离满足：

• 上界： 由于每个 Frechet 分量都是收缩的，所以 ||F(x) - F(y)||_1 ≤ log n * d(x, y)。
• 下界（期望）： 我们可以证明 E[ ||F(x) - F(y)||_1 ] ≥ Ω(1) * d(x, y)。

下界的证明思路是：对于每一对点 x, y，考虑围绕它们的一系列半径递增的球。通过精心选择半径，使得对于每个尺度 t，以常数概率，A_t 包含 x 附近小球中的点，但不包含 y 附近大球中的点。当这个事件发生时，第 t 个分量的贡献至少是某个半径差。对所有尺度求和，这些半径差会叠加为 d(x, y) 的一部分。

从高概率到确定嵌入

上述构造给出了一个期望意义上的好嵌入。为了得到一个对所有点对都同时保持低失真的嵌入，我们可以独立重复此构造 O(log n) 次，并将所有嵌入连接起来。通过浓度不等式和并集界限，可以高概率保证最终嵌入的失真为 O(log n)。

从 L1 嵌入到割

最后一步是将 L1 嵌入转换回一个割度量。假设我们有一个将顶点映射到实数线（一维 L1）的嵌入。我们可以通过随机阈值来产生一个割：随机选择一个阈值 θ，定义 S = { v | f(v) ≤ θ }。可以证明，在这个随机割下，两点 u, v 被分开的概率正好等于 |f(u) - f(v)|。

对于更高维的 L1 嵌入（即 R^k 中的点），我们可以随机选择一个维度，然后在该维度上应用上述随机阈值法。这样得到的随机割，其期望割度量距离与 L1 距离成正比。

算法总结

综上所述，我们得到了一个近似稀疏割的算法：

1. 求解线性规划松弛，得到一个度量 d。
2. 应用 Bourgain 嵌入，将度量 d 嵌入到 L1 空间，失真为 O(log n)。
3. 通过随机阈值法，从 L1 嵌入中抽取一个割 S。
4. 这个割 S 的 Ψ(S) 值（以及 Φ(S) 值）是图最优稀疏割的 O(log n) 因子近似。

课程总结

在本节课中，我们一起学习了稀疏割问题及其与度量嵌入的深刻联系。我们首先定义了图扩展，并讨论了扩展图。然后，我们将稀疏割问题表述为一个优化问题，并通过线性规划进行松弛。为了解决舍入问题，我们深入探讨了 Bourgain 嵌入定理，该定理表明任何度量都可以低失真地嵌入到 L1 空间。最后，我们展示了如何从 L1 嵌入中恢复出一个实际的割，从而得到一个高效的 O(log n) 近似算法。这套将组合优化问题与几何嵌入相联系的方法，是算法设计中的一个强大工具。

14：引言与切格不等式

在本节课中，我们将开始学习一个新的主题：谱图理论。我们将了解如何通过研究与图相关的矩阵的特征值和特征向量来理解图的结构。本节课将介绍拉普拉斯矩阵、其物理意义，并证明一个核心定理——切格不等式。

概述：什么是谱图理论？

谱图理论的核心思想是，通过研究与图相关的矩阵（如拉普拉斯矩阵或邻接矩阵）的特征值和特征向量，来理解图的性质和结构。

给定一个图 G = (V, E)，我们可以定义一个与之相关的二次型，即图的拉普拉斯算子。

图的拉普拉斯矩阵

图的拉普拉斯算子是一个二次型，我们称之为 L_G(x)。其定义如下：

公式： L_G(x) = Σ_{(i,j)∈E} (x_i - x_j)^2

这个二次型总是大于等于零，因此它是一个半正定矩阵。我们可以将其写成矩阵形式：

公式： L_G = D - A

其中，D 是度矩阵（一个对角矩阵，对角线元素为每个顶点的度数），A 是邻接矩阵（如果顶点 i 和 j 之间有边，则 A_{ij} = 1）。

对于 d-正则图（每个顶点的度数都是 d），度矩阵 D 简化为 d 乘以单位矩阵。

由于 L_G 是一个实对称半正定矩阵，它有 n 个特征值和 n 个正交的特征向量。我们记特征值为 λ_1 ≤ λ_2 ≤ ... ≤ λ_n，对应的特征向量为 v_1, v_2, ..., v_n。

我们总是有 λ_1 = 0，其对应的特征向量是常数向量 v_1 = (1, 1, ..., 1)^T。

特征值的变分刻画

理解特征值的一个非常有用的方式是通过瑞利商和 Courant-Fischer 定理。对于一个实对称矩阵 M，其特征值 λ_k 可以刻画为：

公式： λ_k = min_{dim(V)=k} max_{x∈V, x≠0} (x^T M x) / (x^T x)

这意味着 λ_k 是你能找到的 k 维子空间中，所有向量瑞利商的最大值的最小可能值。特别地，λ_1 是使得瑞利商最小的最佳一维子空间。

将这个定理应用于拉普拉斯矩阵，λ_k 可以理解为在某种意义下，图在 k 维空间中“最佳嵌入”的度量。

物理直观：弹簧网络

理解拉普拉斯特征向量的一个直观方式是想象一个弹簧网络。将图的每个顶点看作一个点，每条边看作一根弹簧。这个弹簧系统的总能量与拉普拉斯二次型成正比：

公式： 能量 ∝ Σ_{(i,j)∈E} (x_i - x_j)^2

因此，寻找拉普拉斯矩阵的特征向量和特征值，本质上是在寻找这个弹簧网络在满足正交性约束下的低能量状态（平衡态）。第一特征向量对应所有顶点重合的零能量状态。后续的特征向量则给出了将图“展开”到一维、二维或更高维空间的方式，同时最小化弹簧系统的总能量。

通过计算图的第二、第三特征向量，我们可以将图的顶点映射到实数线上或平面上，从而得到图的一种自然“嵌入”。这种嵌入仅依赖于图的连接信息，却能神奇地反映出图的一些几何或聚类结构。

上一节我们介绍了拉普拉斯矩阵和特征值的直观意义，本节中我们来看看特征值与图的一个基本性质——连通性之间的关系。

特征值与连通分量

图 G 至少有 k 个连通分量，当且仅当其拉普拉斯矩阵的第 k 个特征值 λ_k = 0。

证明思路：

• 如果图有 k 个连通分量：设这 k 个分量为 S_1, S_2, ..., S_k。我们可以构造 k 个向量 w_1, w_2, ..., w_k，其中 w_i 是集合 S_i 的指示向量（在 S_i 中为1，否则为0）。这些向量是正交的，并且每个向量的瑞利商都为 0。根据 Courant-Fischer 定理，这意味着 λ_k ≤ 0。又因为拉普拉斯矩阵半正定，所以 λ_k = 0。
• 如果 λ_k = 0：这意味着存在一个 k 维子空间，其中所有向量的瑞利商为 0。瑞利商为 0 意味着对于该子空间中的任何向量 x，所有相连的顶点 i, j 都有 x_i = x_j。由此可以推导出图至少被分成了 k 个连通分量。

这个结论虽然简单，但它引出了一个更深刻的问题：如果 λ_2 很小但不为零，那意味着什么？这就引出了图“扩展性”的概念。

图的扩展性

为了量化图的“连通程度”，我们引入扩展性的概念。对于一个图 G 中的顶点子集 S（其大小不超过 |V|/2），其扩展性 φ(S) 定义为：

公式： φ(S) = |E(S, V\S)| / (d * |S|) （对于 d-正则图）

其中，|E(S, V\S)| 是横跨集合 S 和其补集之间的边数。d * |S| 是集合 S 的“体积”。

图 G 的扩展性 φ(G) 定义为所有大小不超过 |V|/2 的集合 S 中，φ(S) 的最小值。

扩展性 φ(G) 衡量了图的最薄弱环节。φ(G) 越大，图整体连通性越好；φ(G) 越小，说明图容易被切成两个联系稀疏的部分。

上一节我们定义了图的扩展性，本节的核心是建立扩展性与第二特征值之间的定量关系，即切格不等式。

切格不等式

我们通常使用归一化的拉普拉斯矩阵来简化表述，对于 d-正则图，定义为：

公式： L̃_G = (1/d) * L_G

其第二特征值记为 λ₂(L̃_G)。切格不等式指出：

公式： φ(G)² / 2 ≤ λ₂(L̃_G) ≤ 2 φ(G)

这个不等式建立了图的最优切割（扩展性 φ(G)）与谱信息（第二特征值 λ₂）之间的桥梁。它意味着：

• 如果 λ₂ 很小，那么图存在一个扩展性很差的切割（即容易分成两部分）。
• 反之，如果图的扩展性很好（φ(G) 大），那么 λ₂ 也一定不会太小。

我们通常将 λ₂ ≤ 2 φ(G) 称为不等式的“简单方向”，而将 φ(G)² / 2 ≤ λ₂ 称为“困难方向”。简单方向证明了一个弱上界，而困难方向证明了一个更强的下界。

证明简单方向 (λ₂ ≤ 2 φ(G))

这个方向的证明是构造性的。

1. 找到最小扩展集：设集合 S 是达到最小扩展性 φ(G) 的集合，即 φ(G) = φ(S)。
2. 构造测试向量：考虑集合 S 的指示向量（在 S 中为1，否则为0）。但这个向量不与全1向量正交。
3. 正交化：我们构造向量 w = 指示向量(S) - (|S|/|V|) * 全1向量。这个向量 w 与全1向量正交。
4. 计算瑞利商：可以计算出向量 w 的瑞利商恰好为 2 φ(S)。因为 w 与全1向量正交，根据 Courant-Fischer 定理，第二特征值 λ₂ 不会超过任何此类向量的瑞利商，因此 λ₂ ≤ 2 φ(S) = 2 φ(G)。

证明困难方向 (φ(G)² / 2 ≤ λ₂) 与谱划分算法

这个方向的证明不仅是一个存在性证明，更给出了一个高效的随机算法（谱划分算法），可以从一个具有小瑞利商的向量（如近似第二特征向量）中，构造出一个扩展性较小的切割。

算法步骤如下：

1. 获取向量：计算（或近似计算）拉普拉斯矩阵的第二特征向量 x。每个顶点 i 被赋予一个实数值 x_i。
2. 中心化：将向量 x 平移一个常数，得到新向量 y = x + α * 1，使得恰好有一半的 y_i ≤ 0，一半的 y_i ≥ 0。平移不影响瑞利商。
3. 缩放：将向量 y 缩放为 z = β * y，使得 Σ_i z_i² = 1。缩放也不影响瑞利商。此时，z 的瑞利商仍然是 λ₂。
4. 构造平方嵌入：对于每个顶点 i，考虑其带符号的平方值：sign(z_i) * z_i²。这将所有顶点映射到区间 [-1, 1] 上，且区间总长度为 1。
5. 随机阈值切割：在区间 [-1, 1] 上均匀随机地选取一个阈值 τ。定义集合 S 为所有满足 sign(z_i) * z_i² < τ 的顶点 i 的集合（如果 |S| > |V|/2，则取补集）。

算法分析：
通过计算顶点落入集合 S 的概率以及边被切割的概率，并应用柯西-施瓦茨不等式等技巧，可以证明该算法产生的切割的期望扩展性不超过 √(2 λ₂)。由概率方法可知，至少存在一个阈值 τ 能产生一个扩展性 ≤ √(2 λ₂) 的切割。因此，φ(G) ≤ √(2 λ₂)，平方后即得 φ(G)² / 2 ≤ λ₂。

这个算法非常实用，是谱聚类和图划分的基础。它表明，即使不精确计算特征向量，只要找到一个瑞利商较小的向量，就能高效地找到图的一个稀疏切割。

总结

本节课中我们一起学习了谱图理论的入门知识：

1. 我们引入了图的拉普拉斯矩阵及其特征值和特征向量。
2. 我们了解了特征值的变分刻画及其与弹簧网络模型的物理直观联系。
3. 我们证明了图的连通分量数量与其拉普拉斯矩阵零特征值重数之间的关系。
4. 我们定义了衡量图连通性强弱的扩展性 φ(G)。
5. 我们重点学习并证明了谱图理论的核心定理之一——切格不等式，它建立了图的扩展性 φ(G) 与第二特征值 λ₂ 之间的定量关系。
6. 在证明切格不等式的困难方向时，我们实际上得到了一个非常实用的谱划分算法，该算法能够利用特征向量信息来寻找图中的稀疏切割。

切格不等式是连接图的结构性质与代数性质的桥梁，而谱划分算法则是理论应用于实际的一个典范。在接下来的课程中，我们将继续探索谱图理论的其他有趣主题和应用。

15：谱图理论与随机游走 🧮

在本节课中，我们将继续探索谱图理论，并了解它与图的结构性质（如连通性和扩展性）之间的深刻联系。我们还将学习谱图理论的一个重要应用：分析随机游走在图上的混合时间。

课程概述与回顾 📚

上一节我们介绍了图的拉普拉斯矩阵及其二次型，并学习了Cheeger不等式。该不等式将图的第二小特征值（谱隙）与图的扩展性紧密联系起来。本节中，我们将看看Cheeger不等式的更多推论，并探讨谱图理论如何帮助我们理解随机游走的性质。

Cheeger不等式的推论与应用 🔍

Cheeger不等式告诉我们，图的第二小特征值（λ₂）与图的扩展性（φ(G)）满足以下关系：

λ₂ / 2 ≤ φ(G) ≤ √(2λ₂)

这个不等式在路径图（Path）和环图（Cycle）上是紧的。对于路径图，λ₂ 约为 O(1/n²)，而其扩展性 φ(G) 约为 O(1/n)，这体现了不等式右侧的平方根间隙。

一个自然的问题是，对于更高的特征值（如 λ₃, λ₄...）是否有类似的结论？答案是肯定的。存在一个更广义的定理：如果图的第 k 个特征值 λₖ 很小（例如为 ε），那么图中存在 k 个互不相交的集合 S₁, S₂, ..., Sₖ，使得每个集合的扩展性都很小（至多为 √(λₖ) 乘以某个多项式因子）。这是Cheeger不等式在高维上的推广。

从算法角度看，我们之前学过通过求解线性规划（如Bourgain嵌入）来寻找近似最小割的算法，它能获得 O(log n) 的近似比。而基于Cheeger不等式的谱方法虽然可能给出较差的绝对近似比（例如 √φ(G)），但在扩展性 φ(G) 本身是常数的情况下，它能给出常数近似，并且计算效率极高。

扩展图与谱验证 🌳

扩展图是一类具有高度连通性的稀疏图。从组合角度定义，一个 d-正则图 G 如果其扩展性 φ(G) 至少是一个常数，则它是一个组合扩展图。然而，验证一个图是否具有常数扩展性是NP难的。

Cheeger不等式为我们提供了一种可验证的替代定义：谱扩展图。一个 d-正则图 G 如果其归一化拉普拉斯矩阵的第二小特征值 λ₂ 至少是一个常数，则它是一个谱扩展图。根据Cheeger不等式，组合扩展性和谱扩展性在常数因子内是等价的。因此，在实践中，我们通常通过构造和验证谱扩展图来获得组合扩展图。

平面图的谱性质 ✈️

谱方法不仅能分析抽象图，也能分析具有几何结构的图。例如，对于平面图，我们可以证明其第二小特征值必然很小。

定理：对于任何平面图 G，其归一化拉普拉斯矩阵的第二小特征值 λ₂ 满足 λ₂ ≤ 8/(d n)，其中 d 是图的度（假设为正则或最大度），n 是顶点数。

这个定理符合直觉：平面图（如网格）很容易被稀疏地切割开。证明巧妙地运用了平面图的圆填充表示和球极平面投影。

1. 圆填充定理：任何平面图都可以用一系列互不相交的圆来表示，使得两个顶点之间有边当且仅当对应的圆相切。
2. 球极投影：将平面上的圆填充投影到球面上。该投影保持“圆”的性质，并将整个图形约束在有限大小的球面上。
3. 构造嵌入：以球面上各圆的圆心作为对应顶点的嵌入向量（在三维空间中）。
4. 计算瑞利商：通过计算嵌入向量的瑞利商来界定 λ₂。利用圆的半径和球面面积有限的事实，可以证明瑞利商的上界为 O(1/(d n))，从而 λ₂ 也很小。

这个证明中一个微妙之处在于确保嵌入向量的和为零（即正交于全1向量），这需要用到布劳威尔不动点定理等拓扑学工具。此结论可以推广到更高亏格的曲面上的图。

随机游走与混合时间 🚶‍♂️

现在，让我们看看扩展性/谱隙的一个核心应用：分析随机游走的混合时间。考虑一个 d-正则图 G 上的简单随机游走：从某个顶点 u 开始，每一步均匀随机地走向当前顶点的一个邻居。

设 p_t 是一个向量，其第 v 个分量表示游走在 t 步后位于顶点 v 的概率。那么，概率分布的演化满足以下线性关系：

p_t = (A / d) * p_{t-1} = M * p_{t-1}

其中 A 是邻接矩阵，M = A/d 是随机游走矩阵。进而，我们有 p_t = M^t * p_0。

我们的目标是了解 p_t 何时会接近均匀分布 π = (1/n, ..., 1/n)。我们考察差值向量 p_t - π。由于 π 是 M 的特征值为1的特征向量（Mπ = π），因此差值向量的演化也满足：

p_t - π = M^t * (p_0 - π)

将初始差值向量在 M 的特征向量基上展开。设 M 的特征值为 1 = λ₁ > λ₂ ≥ ... ≥ λ_n ≥ -1，对应的特征向量为 v₁ (=π), v₂, ..., v_n。那么，经过 t 步后，差值向量在每个特征方向上的分量会乘以 λ_i^t。由于 λ₁=1 对应均匀分布方向，而其他特征值 |λ_i| < 1，因此差值向量会以指数速度衰减。

具体地，我们可以得到混合时间的上界。定义 ε-混合时间 T_mix(ε) 为满足 ||p_t - π||₂ ≤ ε 的最小 t。通过分析可得：

T_mix(ε) = O( log(1/ε) / log(1 / |λ|) )

其中 |λ| = max{ |λ₂|, |λ_n| } 是仅次于1的特征值的最大绝对值（对于正则图，常称为谱隙）。对于扩展图（λ₂ 是常数），混合时间是 O(log n) 的，非常快。对于任意连通图，根据Cheeger不等式，λ₂ ≥ φ(G)²/2 ≥ Ω(1/(d² n²))，因此混合时间至多是多项式级别的（如 O(n³ log n)）。

随机游走算法应用示例 💡

上述分析引出了一个简单而重要的对数空间算法：用随机游走判定 s-t 连通性。

算法：给定图 G 和顶点 s, t，从 s 开始进行随机游走，步数为 T = poly(n)（例如 n⁵）。如果在游走过程中访问到 t，则输出“连通”，否则输出“不连通”。

原理：如果 s 和 t 连通，则图是连通的，随机游走将在多项式步内以高概率接近均匀分布，从而以至少 1/n 的概率访问到 t。重复多项式次尝试，访问到 t 的概率极高。如果 s 和 t 不连通，则游走永远无法从 s 所在的连通分量到达 t。

空间复杂度：算法只需记录当前顶点位置和一个步数计数器，每个都需要 O(log n) 比特，因此总空间复杂度为 O(log n)。这比传统的 BFS/DFS 需要 O(n) 空间要节省得多。是否存在确定性的对数空间算法解决 s-t 连通性曾是一个重要问题，已在2005年得到解决。

总结与回顾 🎯

本节课我们一起深入学习了谱图理论的几个关键应用：

1. 我们回顾了Cheeger不等式，它建立了图的谱隙与扩展性之间的桥梁。
2. 我们看到了该不等式如何导致对扩展图的谱验证，并讨论了平面图必然具有小特征值（即容易分割）的定理及其巧妙证明。
3. 我们探讨了谱隙如何控制随机游走的混合时间。扩展图上的随机游走混合极快（O(log n)步），而任何连通图上的随机游走都会在多项式步内混合。
4. 最后，我们以此为基础，介绍了一个经典的对数空间随机算法：使用随机游走判定图的 s-t 连通性。

谱图理论将线性代数的工具与图的结构和随机过程深刻联系起来，是算法设计与分析中一个强大而优美的框架。

16：有效电阻与随机游走

在本节课中，我们将深入探讨有效电阻的概念，并揭示它与随机游走之间优美而深刻的联系。我们将学习如何计算有效电阻，理解其基本性质，并最终证明一个将有效电阻与随机游走的“击中时间”联系起来的核心定理。

概述

有效电阻是电路网络中的一个基本概念。给定一个由电阻（边）和节点组成的网络，以及两个特定的节点A和B，有效电阻是指，如果将整个网络替换为连接A和B的单个电阻，该电阻的值应为多少。我们将看到，这个概念不仅可以通过串联和并联的物理公式计算，还可以通过分析网络中的能量耗散或电势差来理解。更重要的是，有效电阻与图上随机游走的性质有着紧密的联系，这为我们分析图算法提供了强大的工具。

有效电阻的定义与性质

首先，我们明确有效电阻的定义。考虑一个电阻网络，其中每条边的电阻值为1（对应于图的边）。对于任意两个节点 u 和 v，有效电阻 R_eff(u, v) 可以通过以下两种等价方式计算：

1. 能量视角：在 u 点注入1安培的电流，在 v 点引出1安培的电流。此时，整个网络消耗的能量即为 R_eff(u, v)。
2. 电势视角：在上述电流配置下，u 点和 v 点之间的电势差即为 R_eff(u, v)。

基于能量最小化的特性（即电流总是选择能量消耗最小的路径流动），我们可以得到一个直观但重要的性质：向图中添加边只会降低（或保持不变）任意两点间的有效电阻。这是因为添加边扩大了电流可选路径的集合，最小化问题在一个更大的集合上进行，其最优值（最小能耗）不会增加。

随机游走与击中时间

为了建立与随机游走的联系，我们需要引入“击中时间”的概念。

定义（击中时间 H_{uv}）：H_{uv} 表示从一个随机游走从节点 u 出发，首次到达节点 v 所需的期望步数。

需要注意的是，击中时间通常不是对称的，即 H_{uv} 不一定等于 H_{vu}。例如，在一个星形图中，从中心节点到叶节点的击中时间很短，但从叶节点返回中心节点的击中时间则可能很长。

核心定理：有效电阻与击中时间

现在，我们来看连接这两个概念的核心定理。

定理：对于任何无向图 G=(V, E)，设其边数为 m。对于任意两个节点 u 和 v，其击中时间满足：
H_{uv} + H_{vu} = 2m * R_eff(u, v)

等式的左边 H_{uv} + H_{vu} 被称为 u 和 v 之间的“通勤时间”，即从 u 出发走到 v 再回到 u 的期望时间。这个定理表明，通勤时间与两点间的有效电阻成正比。

定理证明思路

证明巧妙地构造了三种不同的电流配置，并利用了电势差与击中时间满足相同线性方程的事实。

1. 情景A：在每个节点 i 注入 d_i（节点 i 的度数）单位的电流，并在特定节点 v 引出总计 2m 单位的电流。可以证明，在此情景下，任意节点 u 与 v 之间的电势差 P_{uv}^{(A)} 恰好等于击中时间 H_{uv}。
2. 情景B：在每个节点 i 注入 d_i 单位电流，但在节点 u（而非 v）引出 2m 单位电流。类似地，此时 v 与 u 之间的电势差 P_{vu}^{(B)} 等于 H_{vu}。
3. 情景C：将情景B的所有电流方向反转。此时，电流变为从 u 点注入 2m 单位，并从每个节点 i 引出 d_i 单位。根据电路理论，反转电流方向会反转电势差符号，因此有 P_{uv}^{(C)} = -P_{vu}^{(B)} = H_{vu}（经过适当调整）。

由于电流分布构成一个线性空间，我们可以将情景A和情景C的电流叠加。

• 叠加后，除了在 u 点有净注入 2m 单位电流、在 v 点有净引出 2m 单位电流外，其他节点的注入和引出相互抵消。这正是测量 u 和 v 之间有效电阻的标准设置：注入 2m 单位电流，测量电势差。
• 因此，叠加情景下的电势差 P_{uv}^{(A+C)} = 2m * R_eff(u, v)。
• 另一方面，根据线性叠加，P_{uv}^{(A+C)} = P_{uv}^{(A)} + P_{uv}^{(C)} = H_{uv} + H_{vu}。

联立上述等式，即得 H_{uv} + H_{vu} = 2m * R_eff(u, v)。证明完毕。

定理的推论与应用

这个优美的定理带来了几个重要的推论：

1. 有效电阻是度量：由于通勤时间满足三角不等式（绕路总比直路耗时更长），因此有效电阻 R_eff(·,·) 也满足三角不等式，可以作为图节点间的一种距离度量。
2. 边有效电阻的上界：对于图中的任何边 (u, v)，其有效电阻 R_eff(u, v) ≤ 1。因为如果只有这一条边，电阻为1；添加其他并联路径只会降低总电阻。这直接推出通勤时间 H_{uv} + H_{vu} ≤ 2m。
3. 图的覆盖时间上界：从任一节点出发，随机游走访问图中所有节点所需的期望时间（覆盖时间）最多为 2m(n-1)。思路是固定一个节点访问顺序，利用通勤时间上界进行求和。

有效电阻与生成树

最后，我们简要介绍有效电阻另一个迷人的表征：它与均匀随机生成树密切相关。

表征定理：考虑在图中发送1单位电流从 u 到 v 的电网络流。以下方法可以得到这个电流：

1. 从图的所有生成树中均匀随机地选取一棵生成树 T。
2. 在树 T 中，找到从 u 到 v 的唯一路径，并沿该路径发送1单位流。
3. 对所有可能的生成树 T 重复此过程，并将得到的流取平均（期望）。

这个平均流就是所求的电网络流。

基于此，一个惊人的结论是：对于任意边 (u, v)，其有效电阻等于该边被包含在一棵均匀随机生成树中的概率。
R_eff(u, v) = Pr[(u, v) 属于随机生成树 T]

直观理解：如果边 (u, v) 是连接 u 和 v 的唯一捷径，那么它几乎总会出现在生成树中，有效电阻接近1。如果存在很多并联路径，那么有些生成树会选择其他路径而不包含 (u, v)，该边出现在树中的概率就会降低，有效电阻也随之减小。

总结

本节课我们一起学习了有效电阻的核心概念及其在图论中的深刻应用。

• 我们首先从物理和能量角度定义了有效电阻，并了解了其单调性等基本性质。
• 接着，我们引入了随机游走的击中时间概念，并证明了连接二者的关键定理：H_{uv} + H_{vu} = 2m * R_eff(u, v)。
• 利用这个定理，我们推导出了有效电阻作为度量、边电阻上界以及图覆盖时间上界等重要推论。
• 最后，我们揭示了有效电阻与均匀随机生成树之间的优美联系：边的有效电阻等于它出现在随机生成树中的概率。

这些内容展示了来自不同领域（电路理论、概率论、图论）的概念如何交织在一起，为解决算法问题提供强大而统一的视角。在下一节课中，我们将回到计算电网络流本身的问题，并由此引向求解拉普拉斯系统这一更广泛的课题。

17：求解拉普拉斯系统 🧮

在本节课中，我们将要学习如何高效求解拉普拉斯线性系统。拉普拉斯系统在科学计算、物理模拟和组合优化中非常普遍。我们将探讨一种近乎线性时间的求解方法，并理解其背后的核心思想。

概述

拉普拉斯矩阵 L 定义为 L = D - A，其中 D 是度矩阵，A 是邻接矩阵。求解拉普拉斯系统 Lx = b 是一个基础问题。传统的求解方法如高斯消元法或梯度下降法，对于大规模图来说效率不高。我们将介绍一种基于随机 Kaczmarz 方法的算法，它利用图的组合结构，能够以近乎线性时间求解此类系统。

从线性系统到随机 Kaczmarz 方法

上一节我们介绍了拉普拉斯系统。本节中我们来看看一种用于求解一般线性系统的迭代算法——随机 Kaczmarz 方法。

该算法非常直观。假设我们有一个线性方程组 A_i x = b_i，我们从一个初始点 x_0 开始。在每一步，我们随机选择一个约束方程 i，然后将当前点 x_t 投影到该约束所定义的超平面上，得到 x_{t+1}。

投影公式如下：
x_{t+1} = x_t + (b_i - A_i \cdot x_t) * A_i / ||A_i||^2

分析表明，如果我们以概率 p_i = ||A_i||^2 / (sum_j ||A_j||^2) 选择第 i 个约束，那么算法会指数收敛。收敛速度取决于矩阵 A 的一个条件数参数 κ(A)。

拉普拉斯系统与电流流

上一节我们介绍了一个通用的线性系统求解器。本节中我们来看看如何将其应用于求解拉普拉斯系统。

首先，我们需要将拉普拉斯系统 Lv = b 重新解释。可以将图视为一个电阻网络，每条边的电阻为 1。向量 b 表示注入每个顶点的边界电流（需满足 sum(b_i) = 0）。我们的目标是求解每个顶点上的电压 v，或者等价地，求解每条边上的电流 f。

根据欧姆定律和基尔霍夫电流定律，边 (i, j) 上的电流 f_{ij} = v_i - v_j，并且对于每个顶点 i，流入的边界电流 b_i 等于从该顶点流出的所有边电流之和。这正是 Lv = b 所描述的。

因此，求解拉普拉斯系统等价于在给定边界电流 b 的情况下，寻找满足基尔霍夫电压定律（即电势差定义一致）的电流流 f。

构建约束系统

上一节我们将问题转化为寻找一个特定的电流流。本节中我们来看看如何用线性约束来描述这个电流流。

所有满足基尔霍夫电流定律（流量守恒）的流构成一个线性空间。而电学流（即可由电势差导出的流）是这个空间的一个子空间，其附加条件是满足基尔霍夫电压定律：对于图中任意一个环，环上各边的电流代数和为零。

为了应用 Kaczmarz 方法，我们需要一个有限的约束集。我们通过以下步骤构建：

1. 选择图的一个生成树 T。
2. 对于每条不在树中的边 e = (u, v)，将其加入树中会形成一个唯一的环 C_e。
3. 对于每个这样的环 C_e，我们施加约束：环上所有边的电流代数和为零。

可以证明，这组环约束足以保证整个流场满足基尔霍夫电压定律。因此，我们的线性系统变为 A f = 0，其中矩阵 A 的每一行对应一个非树边构成的环，每一列对应一条边。

算法效率与低拉伸生成树

上一节我们得到了约束矩阵 A。本节中我们来看看算法的运行时间如何依赖于图的结构。

回顾随机 Kaczmarz 方法的收敛时间与参数 κ(A) 有关。对于我们的矩阵 A，可以证明其最小奇异值至少为 1。而矩阵 A 的 Frobenius 范数 ||A||_F 等于所有约束环的长度之和。

对于一个由非树边 (u, v) 定义的环，其长度为 1 + dist_T(u, v)，其中 dist_T(u, v) 是 u 和 v 在生成树 T 中的距离。我们定义边 (u, v) 的拉伸为 stretch_T(u, v) = dist_T(u, v)（因为原图中边长为1）。那么总 Frobenius 范数约为 m + sum_{e not in T} stretch_T(e)。

因此，算法的迭代次数（正比于 ||A||_F^2）取决于所有非树边的总拉伸。为了获得近乎线性的运行时间，我们需要一个低拉伸生成树，使得总拉伸为 O(m polylog(n))。

幸运的是，对于任何图，都存在这样的低拉伸生成树，并且可以在线性时间内构建。其构造思想类似于我们之前学过的度量空间分解和层次树嵌入。

算法实现与总结

上一节我们分析了理论效率的关键在于低拉伸生成树。本节中我们简要讨论实现并总结全课。

在获得低拉伸生成树后，随机 Kaczmarz 方法的每一步需要随机选择一个环（即非树边），概率与其长度（拉伸+1）成正比，然后更新流 f。直接操作维度为 m（边数）的流向量成本很高。但由于所有更新都围绕着树上的路径进行，我们可以利用高效的数据结构（如树状路径更新）来在近常数时间内完成每次迭代，从而使得总运行时间达到近乎线性。

本节课中我们一起学习了求解拉普拉斯系统的近乎线性时间算法。其核心步骤是：

1. 将拉普拉斯系统 Lv = b 转化为在给定边界电流下寻找电学流 f 的问题。
2. 利用生成树将电学流需满足的基尔霍夫电压定律转化为一系列环约束 A f = 0。
3. 对约束系统应用随机 Kaczmarz 方法。
4. 通过为图构造一个低拉伸生成树，确保约束矩阵 A 的条件良好，从而使迭代次数近乎线性。
5. 利用树的数据结构高效实现每次迭代更新。

这种方法巧妙地将连续优化（Kaczmarz/梯度下降）与离散图论（生成树、低拉伸）结合起来，是理论计算机科学中一个非常优美的成果，并催生了许多其他图算法上的突破。

18：无向图最大流算法

在本节课中，我们将学习如何利用拉普拉斯线性系统求解器（即“电网络流”计算）来设计一个求解无向图最大流的算法。我们将从基础概念出发，逐步构建算法，并最终分析其性能。

概述

最大流问题是一个经典的组合优化问题。我们将在无向、无权（所有边容量为1）的图中寻找最大流。核心思想是使用上一讲介绍的拉普拉斯线性系统求解器，它能高效计算满足流守恒且能量最小的“电网络流”。我们将通过乘性权重更新算法，将满足所有边容量约束的困难问题，转化为一系列只需满足“平均”约束的较简单问题，并用电网络流求解器来实现。

基础概念与问题形式化

首先，我们形式化最大流问题。给定源点 s 和汇点 t，目标是找到从 s 到 t 的最大流量 F。

一个流 f 是一个定义在边上的向量。它必须满足两类约束：

1. 流守恒约束（简单约束）：对于除 s 和 t 外的所有顶点，流入等于流出。对于 s，净流出为 F；对于 t，净流入为 F。所有满足此约束的流构成集合 K。
2. 容量约束（困难约束）：对于每条边 e，流量 f_e 不得超过其容量（这里假设为1）。

因此，最大流问题等价于在集合 K 中找到一个点（流），使其同时满足所有容量约束。直接求解这个线性规划是困难的。

乘性权重算法框架

我们将使用乘性权重算法来协调“寻找可行流”和“检查约束违反”这两个任务。

上一节我们介绍了最大流问题的形式。本节中，我们来看看如何用乘性权重算法框架将其转化为更易处理的问题。

算法的基本设定如下：

• 专家：每条边对应一个“专家”，代表一个需要检查的容量约束。
• 算法（乘性权重）：在每一轮 t，算法维护一个在边（专家）上的概率分布 p^t。
• 预言机（Oracle）：观察到分布 p^t 后，预言机需要给出一个流 f^t。这个流必须满足流守恒约束（即 f^t ∈ K），同时还要满足一个关于 p^t 的“平均”容量约束。
• 收益（Gain）：如果预言机给出的流 f^t 在边 e 上的流量为 f_e^t，则该边在此轮的收益为 g_e^t = f_e^t - 1（因为容量是1，f_e^t - 1 度量了违反程度）。
• 算法的收益：算法在第 t 轮的总收益是期望值 ⟨p^t, g^t⟩，即各边违反程度的加权平均。

乘性权重算法保证，在经过足够多轮后，算法的平均收益不低于任何单个“专家”（即固定检查某条边）在整个过程中的平均收益。

预言机的任务是：给定概率分布 p，找到一个流 f ∈ K，使得平均违反程度 ∑_e p_e * f_e 不超过 1 + ε。这比要求所有 f_e ≤ 1 要容易得多，因为它只是一个线性约束。

如果这样的预言机存在，乘性权重算法经过 T ≈ (ρ log m) / ε^2 轮后（其中 ρ 是“宽度”，即单轮中可能出现的最大违反值），输出的平均流 f* = (1/T) ∑ f^t 将是一个 2ε-近似可行解（即每条边上的流量最多为 1 + 2ε）。

因此，问题的核心转化为：如何构建一个高效的预言机？ 这个预言机需要能快速找到满足流守恒且平均违反程度小的流。

构建预言机：从最短路到电网络流

我们首先尝试一个简单的预言机：最短路预言机。

给定边上的概率分布 p（可视为边长），计算从 s 到 t 的最短路径，并将全部流量 F 沿这条路径发送。

• 优点：如果最短路径长度 ≤ 1/F，则总成本 ∑ p_e f_e ≤ 1，满足要求。
• 缺点：宽度 ρ 可能高达 F（当所有流量挤在一条边上时，f_e - 1 可达 F-1）。在最坏情况下 F ≈ m，导致迭代轮数 T ≈ O(m^2/ε^2)，效率不高。

为了降低宽度，我们需要让流更均匀地分散开。这引出了更聪明的预言机：电网络流预言机。

以下是电网络流预言机的步骤：

1. 将每条边 e 的电阻 r_e 设置为 r_e = p_e + ε/m。
2. 使用拉普拉斯求解器，计算在此电阻网络中，从 s 注入 F 单位电流、从 t 流出 F 单位电流时，各边上的电网络流 f̂。

电网络流的关键性质是：在所有满足流守恒的流中，它最小化总能量消耗 ∑_e f_e^2 * r_e。

接下来我们分析这个预言机的性能，即它的平均违反程度和宽度。

首先估计宽度，即单边上可能的最大流量。设存在一个满足所有容量约束（即 f_e* ≤ 1）的最优流 f*。其能量为：
∑_e (f_e*)^2 * r_e ≤ ∑_e 1 * (p_e + ε/m) = 1 + ε
由于电网络流 f̂ 最小化能量，因此其能量也满足 ∑_e (f̂_e)^2 * r_e ≤ 1 + ε。
对于任何边 e，有 (f̂_e)^2 * (ε/m) ≤ 1 + ε，由此可得 f̂_e ≤ O(√(m/ε))。因此，宽度 ρ = O(√(m/ε))，相比最短路预言机的 O(m) 有了显著改进。

其次，检查平均违反程度。我们的目标是 ∑_e p_e * f̂_e。利用柯西-施瓦茨不等式：
∑_e p_e * f̂_e ≤ (∑_e p_e)^{1/2} * (∑_e p_e (f̂_e)^2)^{1/2} ≤ 1 * (1 + ε)^{1/2} ≤ 1 + ε
因此，电网络流 f̂ 确实满足了平均违反程度 ≤ 1 + ε 的要求。

算法总结与复杂度分析

现在我们将所有部分组合起来，形成完整的最大流算法。

1. 外层框架：使用乘性权重算法。专家是图的边。
2. 内层预言机：在每一轮，根据乘性权重算法产生的边概率分布 p，构造电阻网络并求解电网络流 f̂，将其返回给乘性权重算法。
3. 更新与迭代：乘性权重算法根据返回的流 f̂ 计算各边收益 g_e = f̂_e - 1，并据此更新边的概率分布。
4. 输出：经过 T 轮后，输出所有轮次流的平均值 f*。

复杂度分析：

• 宽度：ρ = O(√(m/ε))
• 乘性权重迭代轮数：T = O( (ρ log m) / ε^2 ) = O( (√m log m) / ε^{2.5} )
• 每轮计算成本：使用近线性时间的拉普拉斯求解器，计算一次电网络流需要 Õ(m) 时间。
• 总时间复杂度：Õ(m * T) = Õ(m^{3/2} / ε^{2.5})。

这是一个针对无向图的近似最大流算法。通过更精巧的设计（例如，当某边流量过大时将其删除并递归处理），可以将复杂度进一步改进到约 Õ(m^{4/3})，最终甚至能发展出近线性时间的算法。这些改进的核心思想是利用组合预条件子来加速优化过程。

总结

本节课我们一起学习了如何利用电网络流求解器和乘性权重算法来求解无向图的最大流问题。关键步骤包括：

1. 用乘性权重框架将复杂的容量约束转化为一系列平均约束。
2. 构建电网络流预言机来高效解决每个平均约束子问题，其优势在于能自动分散流量，从而降低算法的宽度参数。
3. 最终得到一个时间复杂度为 Õ(m^{3/2}) 的近似算法。

这种方法展示了将连续优化工具（拉普拉斯求解器）与组合优化经典问题相结合的强大能力，并为后续更快的算法奠定了基础。

19：半正定规划入门 🧮

在本节课中，我们将要学习一种强大的优化工具——半正定规划。它是线性规划的推广，结合了线性规划、凸优化和谱方法（如特征值计算）的思想，为我们提供了一种表达和解决复杂问题的高级“编程语言”。

什么是半正定规划？🤔

上一节我们介绍了课程概述，本节中我们来看看半正定规划的具体定义。

半正定规划本质上是一个定义在半正定矩阵集合上的线性规划。具体来说，我们有一个由变量构成的 n x n 矩阵 X，我们要求这个矩阵是半正定的。

一个矩阵 M 是半正定的，当且仅当满足以下等价条件之一：

• 所有特征值非负。
• 对于所有实向量 x，有 xᵀMx ≥ 0。
• 存在矩阵 V，使得 M = VVᵀ。
• 存在系数 cᵢⱼ，使得 xᵀMx 可以表示为平方和：∑ (∑ cᵢⱼ xⱼ)²。

所有半正定矩阵构成 ℝ^(n²) 空间（或考虑对称性后的 ℝ^(n + n choose 2)）中的一个凸锥。

一个标准的半正定规划问题形式如下：

• 变量：一个 n x n 的对称矩阵 X。
• 约束：X ≽ 0（表示 X 是半正定矩阵），以及 m 个关于矩阵元素的线性约束：〈Aₗ, X〉 ≤ bₗ，其中 〈A, X〉 = trace(AᵀX) = ∑ᵢⱼ Aᵢⱼ Xᵢⱼ。
• 目标：最小化（或最大化）一个线性目标函数：〈C, X〉。

向量视角的解释 📐

上一节我们介绍了矩阵形式的半正定规划，本节中我们来看看等价的向量形式解释。

根据半正定矩阵的性质，存在向量 v₁, …, vₙ ∈ ℝⁿ，使得矩阵 X 的元素满足 Xᵢⱼ = vᵢ · vⱼ。因此，半正定规划可以重新解释为寻找一组向量 {vᵢ}，并满足以下条件：

• 目标：最小化 ∑ᵢⱼ Cᵢⱼ (vᵢ · vⱼ)。
• 约束：对于 l = 1 … m，有 ∑ᵢⱼ Aᵢⱼ⁽ˡ⁾ (vᵢ · vⱼ) ≤ bₗ。

在这种视角下，我们只能对向量之间的内积（即旋转不变的量）施加线性约束，而不能直接约束向量的具体坐标。

经典应用：最大割问题 ✂️

上一节我们理解了半正定规划的基本形式，本节中我们来看一个经典应用案例——最大割问题。

最大割问题描述如下：给定一个无向图 G = (V, E)，目标是找到一个划分 (S, V\S)，使得被切割的边数（即连接 S 和 V\S 的边）最多。这是一个NP难问题。

Goemans-Williamson SDP松弛

我们为最大割问题设计一个半正定规划松弛。思路是为每个顶点 i 分配一个单位向量 vᵢ（即 ‖vᵢ‖² = 1）。对于最优割 (S, V\S)，我们期望的“理想”解是：如果 i ∈ S，则 vᵢ = u（某个单位向量）；如果 i ∉ S，则 vᵢ = -u。这样，对于边 (i, j)，‖vᵢ - vⱼ‖² 在边被切割时为4，否则为0。

因此，切割的边数等于 (1/4) ∑_{(i,j)∈E} ‖vᵢ - vⱼ‖²。我们构建如下SDP：

• 最大化：(1/4) ∑_{(i,j)∈E} ‖vᵢ - vⱼ‖²
• 约束：对于所有 i，‖vᵢ‖² = 1。

这个SDP是原问题的一个松弛，因为任何实际的割（对应 ±1 赋值）都是该SDP的可行解，因此SDP的最优值至少是最大割的值。

随机超平面舍入法

求解SDP后，我们得到一组单位向量 {vᵢ}。为了得到一个实际的割，我们采用一个简单而优雅的舍入方法：

1. 随机选取一个通过原点的超平面（等价于随机选取一个法向量 g）。
2. 对于每个顶点 i，如果 vᵢ · g ≥ 0，则将其放入 S；否则放入 V\S。

以下是算法性能的分析：
对于一条边 (i, j)，其被切割的概率等于向量 vᵢ 和 vⱼ 之间的夹角 θ 与 π 的比值，即 Pr[边 (i,j) 被切割] = θ / π，其中 cos θ = vᵢ · vⱼ。

因此，算法切割边数的期望值为 ∑_{(i,j)∈E} (θᵢⱼ / π)。我们需要证明这个值至少是SDP最优值（进而也是最大割最优值）的 α 倍。

通过逐项比较，我们需要最小化比值：(θ/π) / ( (1/4)‖vᵢ - vⱼ‖² )。由于 ‖vᵢ - vⱼ‖² = 2 - 2 cos θ = 4 sin²(θ/2)，该比值化为 (θ/π) / sin²(θ/2)。通过数值计算，这个比值的最小值 α_GW ≈ 0.878。

因此，Goemans-Williamson 算法是一个 0.878-近似算法。在唯一游戏猜想下，这个近似比是紧的。

另一个应用：介数约束问题 🔢

上一节我们看到了SDP在组合优化中的成功应用，本节中我们再看一个简洁的例子——介数约束问题。

问题描述：我们需要寻找 {1, …, n} 的一个排列 π，满足一系列介数约束，形式为“元素 i 在元素 j 和 k 之间”。即使存在满足所有约束的排列，找到它也是NP难的。我们的目标是满足尽可能多的约束。

SDP建模与舍入

我们为每个元素 i 分配一个向量 vᵢ。对于一个约束“i 在 j 和 k 之间”，其几何意义是向量 (vⱼ - vᵢ) 和 (vₖ - vᵢ) 方向相反（夹角为钝角）。因此，我们添加约束：(vⱼ - vᵢ) · (vₖ - vᵢ) ≤ 0。

我们求解这个可行性SDP得到向量 {vᵢ}。为了得到一个排列，我们随机选取一个方向 g，并将每个向量投影到该方向上，得到标量 xᵢ = vᵢ · g。然后按照 xᵢ 的升序对元素进行排序，得到最终的排列。

性能分析

对于单个介数约束，由于 (vⱼ - vᵢ) 和 (vₖ - vᵢ) 的夹角 θ ≥ 90°，在随机投影下，元素 i 的投影值落在 j 和 k 的投影值之间的概率至少为 θ/π ≥ 1/2。因此，期望至少有一半的约束被满足。这是一个简单的常数因子近似算法。

总结 📝

本节课中我们一起学习了半正定规划的基础知识及其两个经典应用。

• 我们首先定义了半正定规划，它是在半正定矩阵锥上进行的线性优化，并可以等价地从向量内积的角度理解。
• 接着，我们深入探讨了最大割问题的Goemans-Williamson算法，该算法通过SDP松弛和随机超平面舍入，获得了突破性的0.878近似比，展示了SDP将谱方法（特征值）与线性约束相结合的强大能力。
• 最后，我们看了介数约束问题，通过将顺序关系转化为向量间的几何约束，并用随机投影舍入，得到了一个简单的常数因子近似算法。

这些例子表明，半正定规划为我们提供了一套系统且强大的框架，可以将复杂的组合优化问题松弛为可高效求解的凸优化问题，并通过巧妙的舍入技术获得高质量的近似解。在接下来的课程中，我们将学习如何更系统地构建这类SDP松弛。

20：半正定规划与多项式系统

在本节课中，我们将继续学习半正定规划。我们将看到如何将组合优化问题（如顶点覆盖）表达为多项式不等式系统，并介绍一种强大的工具——平方和半正定规划层次结构，它能系统地生成这些问题的凸松弛。

上一节我们介绍了半正定规划及其在最大割问题中的应用。本节中，我们来看看如何将更一般的问题表述为多项式系统，并利用半正定规划进行松弛求解。

从组合问题到多项式系统

许多组合优化问题可以自然地写成一个多项式（不）等式系统。我们以最小顶点覆盖问题为例。

输入：一个图 G = (V, E)。
目标：找到一个最小的顶点子集 S ⊆ V，使得每条边至少有一个端点在 S 中。

为了将其表述为多项式系统，我们为每个顶点 i 引入一个变量 x_i，其意图是：

• x_i = 1 表示顶点 i 在覆盖集 S 中。
• x_i = 0 表示顶点 i 不在覆盖集 S 中。

于是，我们可以写出以下约束：

1. 布尔约束：对于所有 i ∈ V，有 x_i^2 - x_i = 0。这确保了 x_i ∈ {0, 1}。
2. 覆盖约束：对于所有边 (i, j) ∈ E，至少有一个端点在 S 中。这等价于 (1 - x_i)(1 - x_j) = 0。
3. 目标函数：我们希望最小化 ∑ x_i。我们可以先考虑判定版本：是否存在一个大小不超过 C 的顶点覆盖？即 ∑ x_i ≤ C。

这样，我们就得到了一个关于变量 x_1, ..., x_n 的多项式系统。类似地，最大割问题也可以表述为多项式系统（例如，约束 x_i^2 = 1 和优化 ∑_{(i,j)∈E} (x_i - x_j)^2）。

然而，直接求解多项式系统是 NP 难的。因此，我们需要寻找一种可计算的松弛方法。

凸化与矩方法

多项式系统求解困难的原因在于其解空间是非凸的。一个自然的想法是将其“凸化”。给定解集 S，一个明显的凸化是考虑 S 上的所有概率分布。概率分布的集合是凸的。

但直接处理整个概率分布是不可行的，因为我们需要为指数多个赋值（x ∈ S）指定概率。一个更可行的思路是只处理概率分布的低阶矩，即变量及其乘积的期望值。

对于顶点覆盖问题，我们关注一阶矩和二阶矩：

• 一阶矩：X_i = E[x_i]
• 二阶矩：X_{ij} = E[x_i x_j]

我们的目标是找到满足以下条件的矩 {X_i, X_{ij}}：

1. 原始约束的线性化：原始多项式约束会对矩产生线性约束。例如：
   • 由 x_i^2 - x_i = 0 可得 X_{ii} - X_i = 0。
   • 由 (1 - x_i)(1 - x_j) = 0 可得 1 - X_i - X_j + X_{ij} = 0。
   • 由 ∑ x_i ≤ C 可得 ∑ X_i ≤ C。
2. 矩矩阵的半正定性：我们要求由这些矩构成的矩矩阵 M_2 是半正定的。M_2 的行和列对应单项式 {1, x_1, ..., x_n}，其 (A, B) 位置的元素是 E[Mon_A * Mon_B]。例如：

   M_2 = [ 1,    X_1,    X_2,    ..., X_n;
           X_1,  X_{11}, X_{12}, ..., X_{1n};
           X_2,  X_{21}, X_{22}, ..., X_{2n};
           ...   ...     ...     ...  ...;
           X_n,  X_{n1}, X_{n2}, ..., X_{nn} ]
   

   要求 M_2 ≽ 0 等价于要求对于所有系数向量 p，有 E[(p_0 + ∑ p_i x_i)^2] ≥ 0，这显然是任何真实概率分布都必须满足的性质。

这样，我们就得到了一个关于矩变量 {X_i, X_{ij}} 的半正定规划：在满足一系列线性约束和矩矩阵 M_2 半正定的条件下进行优化。这就是二阶平方和半正定规划松弛。

平方和半正定规划层次结构

我们不必停留在二阶矩。可以固定一个偶数 D = 2d，考虑所有阶数不超过 D 的矩 {E[x^α]}，其中 α 是一个多重指数且 |α| ≤ D。相应的松弛称为 d 阶平方和半正定规划松弛。

以下是其更简洁的表述方式。我们寻找一个伪期望算子 Ẽ，它是一个线性泛函，作用于阶数不超过 D 的多项式，并满足：

1. 对于原系统中的每个约束 P_i(x) ≥ 0，有 Ẽ[P_i(x)] ≥ 0。
2. 对于任何阶数不超过 d 的多项式 Q(x)，有 Ẽ[Q(x)^2] ≥ 0。这等价于相应的 D 阶矩矩阵是半正定的。

对于每个 d = 1, 2, ...，我们都可以写出这样一个 SDP。这就形成了一个平方和 SDP 层次结构。

• d=1 通常对应基本的线性规划松弛。
• 随着 d 增大，SDP 的规模（变量数和约束数）快速增长（约 n^O(d)），计算代价也增大，但松弛也变得更紧。
• 当 d 足够大（例如，对于 n 个布尔变量的问题，d = n）时，该 SDP 松弛是精确的，但求解它需要指数时间。

这个层次结构的力量在于它能自动地、系统地推导出许多有效的约束，包括一些非局部的全局约束，而不仅仅是手工添加的局部约束（如“三角形不等式”）。

层次结构的能力与意义

平方和 SDP 层次结构比早期的线性规划层次结构更强大。一个关键的性质是：对于定义在布尔域（x_i^2 = x_i）上的问题，d 阶 SOS SDP 的解，在任意 2d 个变量的子集上，都对应着一个真实的概率分布。这意味着该 SDP 能够“捕获”所有关于至多 2d 个变量的、在所有真解上都成立的约束。

对于许多经典的组合优化问题（如顶点覆盖、最大割），其近似比在层次结构的低阶（如二阶）就达到了已知的最好结果，并且没有证据表明更高阶能带来改进。然而，在平均情况问题、统计推断和某些随机模型中的算法设计方面，平方和层次结构展现了巨大的潜力，它仿佛一个“编译器”，能将某种形式的证明转化为高效的算法。我们将在下一讲中探讨对偶性和证明的概念。

本节课中我们一起学习了如何将组合优化问题表述为多项式系统，并引入了强大的平方和半正定规划层次结构作为系统的松弛框架。我们看到了它如何通过矩和伪期望算子将问题凸化，以及其层次结构如何提供一系列从宽松到紧密的松弛。下一节，我们将深入探讨如何分析这些基于平方和 SDP 的算法。

21：SOS证明系统与伪期望

在本节课中，我们将学习如何分析从平方和（SOS）半定规划（SDP）层次结构中得到的“伪期望”或“伪矩”。我们将探讨如何从这些伪矩中提取原始问题的解，并理解其背后的核心思想——对偶性与证明系统。我们将看到，SOS SDP的对偶对象正是“平方和证明”，而“弱对偶性”为我们分析SDP解提供了一个强大的工具。

---

回顾：平方和（SOS）SDP层次结构

上一节我们定义了度数为 D 的平方和（SOS）SDP层次结构。它计算一个多项式系统解的分布（或伪分布）的矩（或伪矩）。

给定一个多项式系统：
[
P_i(x) \geq 0, \quad i = 1, \ldots, m
]
度数为 D 的 SOS SDP 解本质上会给你一组伪矩。我们称之为“伪矩”，因为你永远无法确定是否存在一个真实的分布以这些数为矩。这类似于顶点覆盖线性规划（LP）的分数解：可能不存在一个真实的顶点覆盖解分布恰好具有这些分数值，但这就是SDP给出的结果。

我们尚未分析如何证明关于这些伪矩的陈述，或者如何从伪矩中提取原始问题的解。

---

对偶性与证明系统

在深入SOS SDP之前，让我们暂时忘记SDP，重新思考对偶性的概念。对偶性通常在线性规划（LP）的优化背景下讨论，但它背后有一个更广泛的思想：对偶性与证明。

线性规划中的证明系统

假设我们有一个线性规划，其可行域由以下约束定义：
[
A_i \cdot y \leq b_i, \quad i = 1, \ldots, m
]
这些约束共同定义了一个多面体 (P)，即LP的可行域。

如果我们知道一个点 (y) 在这个可行域内，那么关于 (y) 我们所知的唯一事实就是这些约束。我们可以问：是否存在其他关于 (y) 的线性不等式陈述 (c \cdot y \leq d) 也成立？

如何证明这样一个陈述？我们可以从公理（即原始约束）出发，通过非负线性组合来推导新的陈述。例如：
[
\lambda_1 (A_1 \cdot y) + \lambda_2 (A_2 \cdot y) \leq \lambda_1 b_1 + \lambda_2 b_2, \quad \lambda_1, \lambda_2 \geq 0
]
这为我们提供了一个证明系统：规则是取公理的非负线性组合。

• 健全性： 每个能被证明的陈述都是真的。这对应于弱对偶性，在这个系统中是显然成立的。
• 完备性： 每个真的陈述都能被证明。这对应于强对偶性。对于线性规划，这意味着如果对于所有可行解 (y) 都有 (c \cdot y \leq d) 成立，那么总存在一组非负系数 (\lambda_i)，使得：
  [
  c = \sum_{i} \lambda_i A_i \quad \text{且} \quad d \geq \sum_{i} \lambda_i b_i
  ]
  这个恒等式就是证明。对偶线性规划本质上就是在寻找这样的证明（即系数 (\lambda_i)）。

关键要点是：在线性规划的可行域上，每个真的线性不等式都可以通过取约束的非负线性组合来证明。

---

多项式等式系统的证明

这个概念在数学中非常基础。考虑一个多项式等式系统：
[
P_1(x) = 0, \quad P_2(x) = 0, \ldots, P_m(x) = 0
]
假设对于所有满足该系统的 (x)，都有 (Q(x) = 0) 成立。我们如何证明这一点？

一个自然的证明（证书）是展示一个恒等式：
[
Q(x) = \sum_{i=1}^{m} R_i(x) P_i(x)
]
其中 (R_i(x)) 是多项式。如果我们看到这样的恒等式，就知道只要 (P_i(x)=0)，就必有 (Q(x)=0)。希尔伯特零点定理指出，只要陈述为真，这样的证明总是存在的（尽管可能需要对 (Q(x)) 取幂）。

---

多项式不等式系统的证明：平方和（SOS）证明

现在考虑包含不等式的多项式系统：
[
P_i(x) \geq 0, \quad i = 1, \ldots, m
]
假设该系统蕴含另一个不等式 (Q(x) \geq 0)。我们如何证明它？

一个多项式非负的最自然证明是展示它是平方和：
[
Q(x) = \sum_{j} S_j(x)^2
]
但这还没有使用我们的公理 (P_i(x) \geq 0)。为了结合公理，我们考虑以下形式的证明：
[
Q(x) = \underbrace{\sum_{j} S_{0j}(x)^2}{\text{平方和项}} + \underbrace{\sum^{m} P_i(x) \left( \sum_{k} S_{ik}(x)^2 \right)}{\text{带系数的公理}}
]
这里，(S(x)) 和 (S_{ik}(x)) 都是多项式。由于：

1. 平方和项显然非负。
2. 每个 (P_i(x) \geq 0)（公理）。
3. 每个 (\sum_k S_{ik}(x)^2 \geq 0)。
   因此整个表达式非负，从而证明了 (Q(x) \geq 0)。

定义（平方和证明）：
给定多项式系统 (P_i(x) \geq 0)，一个证明 (Q(x) \geq 0) 的度数为 D 的平方和证明是一个如下形式的恒等式：
[
Q(x) = S_0(x) + \sum_{i=1}^{m} P_i(x) S_i(x)
]
其中 (S_0(x)) 和每个 (S_i(x)) 都是平方和多项式，且等式中出现的所有多项式的最高次数不超过 D。

Positivstellensatz 定理指出，只要 (Q(x) \geq 0) 在系统下为真，就总存在这样的平方和证明（可能需要乘以一个多项式因子）。同样，一般没有实用的度数上界。

---

SOS SDP 的对偶与平方和证明

现在，让我们将SOS证明系统与SOS SDP的对偶性联系起来。

考虑最大割（Max Cut）的SDP松弛（用矩阵形式表示）：
[
\begin{aligned}
& \max \quad L \cdot X \
& \text{s.t.} \quad X_{ii} \leq 1, \quad \forall i \
& \qquad \quad X \succeq 0 \quad \text{(X是半正定矩阵)}
\end{aligned}
]
半正定约束 (X \succeq 0) 等价于无穷多个线性约束：(v^T X v \geq 0, \ \forall v)。

如果我们为这个（无穷维）LP写对偶，引入对偶变量 (\alpha_i)（对应 (X_{ii} \leq 1)）和 (\beta_v)（对应 (v^T X v \geq 0)），对偶问题会寻找一个如下形式的恒等式来证明目标函数的上界 (B)：
[
B - L \cdot X = \sum_i \alpha_i (1 - X_{ii}) + \sum_v \beta_v (v^T X v)
]
其中 (\alpha_i \geq 0, \beta_v \geq 0)。项 (\sum_v \beta_v (v^T X v)) 是半正定矩阵的非负组合，其本身也是一个半正定矩阵，记作 (Y)。

现在，从伪期望的角度解释这个恒等式。在度数为2的SOS SDP中，矩阵 (X) 的条目对应于伪期望：(X_{ij} = \tilde{\mathbb{E}}[x_i x_j])。目标 (L \cdot X) 对应伪期望 (\tilde{\mathbb{E}}[\sum_{(i,j)\in E} (x_i - x_j)^{2])。类似地，(v}T X v = \tilde{\mathbb{E}}[(\sum_i v_i x_i)^2])。

因此，上述对偶恒等式在多项式世界中转化为：
[
B - \sum_{(i,j)\in E} (x_i - x_j)^2 = \sum_i \alpha_i (1 - x_i^2) + \sum_v \beta_v (\sum_i v_i x_i)^2
]
这正是我们之前定义的度数为2的平方和证明！ 右边第一项是平方和（常数），第二项是公理 ((1 - x_i^2) \geq 0) 乘以平方和系数 ((\sum_i v_i x_i)^2)。

核心关联：

• 度数为 D 的 SOS SDP 的对偶正是在寻找度数为 D 的平方和证明。
• 随着度数 D 增加，证明系统的能力增强。

---

伪期望与弱对偶性

现在我们可以形式化地定义伪期望，并陈述连接伪期望与平方和证明的关键性质。

定义（伪期望）：
对于一个多项式系统 (P_i(x) \geq 0)，一个度数为 D 的伪期望 (\tilde{\mathbb{E}}) 是一个从次数不超过 D 的多项式到实数的线性泛函，满足：

1. (\tilde{\mathbb{E}}[1] = 1)。
2. 对于所有 (i) 和所有次数足够低的多项式 (S(x))，有 (\tilde{\mathbb{E}}[P_i(x) \cdot S(x)^2] \geq 0)。
3. 对于所有次数足够低的多项式 (S(x))，有 (\tilde{\mathbb{E}}[S(x)^2] \geq 0)。

度数为 D 的 SOS SDP 会输出这样一个伪期望。

弱对偶性（关键分析工具）：
如果不等式 (Q(x) \geq 0) 拥有一个度数为 D 的平方和证明（即可以写成前述形式），那么对于该多项式系统的任何度数为 D 的伪期望 (\tilde{\mathbb{E}})，都有：
[
\tilde{\mathbb{E}}[Q(x)] \geq 0
]
这个陈述非常基础但强大。它意味着：如果你使用一个“简单的”（低度数）平方和证明来证明某个事实 (Q(x) \geq 0)，那么这个事实在SDP返回的“伪世界”解中也必然成立。

强对偶性（大致上说，每个对伪期望成立的事实也能被低度平方和证明）在大多数关心的情况下也成立，但需要更仔细的论证。

---

总结与应用展望

本节课我们一起学习了：

1. 证明系统的概念： 在线性规划中，证明通过对约束取非负线性组合完成。
2. 平方和证明： 对于多项式不等式系统，证明通过将目标写成平方和与带平方和系数的公理组合来完成。
3. SOS SDP的对偶： 度数为 D 的 SOS SDP 的对偶问题正是在寻找度数为 D 的平方和证明。
4. 伪期望与弱对偶性： 这是分析SOS SDP算法的核心工具。如果我们能用低度平方和证明某个性质对于真实解成立，那么该性质对于SDP输出的伪期望也成立。

在下一节课中，我们将把这个理论应用到具体算法中。策略如下：

• 对于一个难题（例如鲁棒线性回归），我们将其表述为一个多项式系统。
• 我们求解该系统的度数为 D 的 SOS SDP，获得伪期望 (\tilde{\mathbb{E}})。
• 为了证明 (\tilde{\mathbb{E}}) 给出的解接近真实答案，我们将在真实解的世界里，使用低度数的平方和证明来证明“任何真实解都接近真实答案”。
• 根据弱对偶性，这个证明意味着伪期望 (\tilde{\mathbb{E}}) 也满足同样的性质，从而保证SDP返回的解是高质量的。

这种将代数证明“编译”为SDP解性能保证的框架，是分析基于SOS的算法的强大语言。

22：平方和（SOS）SDP的应用

在本节课中，我们将学习平方和（SOS）半正定规划（SDP）的两个具体应用。我们已经介绍了SOS SDP及其对偶性，并构建了相关理论框架。今天，我们将看到如何利用这些工具解决实际问题。

应用一：鲁棒线性回归

首先，我们来看第一个应用：鲁棒线性回归问题。

问题定义

假设存在一个未知的线性函数，由一个单位向量 L̂ 定义，即 L̂ · x。我们从这个函数中获取样本。真实的样本形式为 (x̂, L̂ · x̂ + γ)，其中 x̂ 是随机向量（例如随机的±1向量），γ 是高斯噪声。

在标准线性回归中，我们的目标是通过最小化平方误差来找到这个线性函数 L。然而，在鲁棒线性回归中，我们面临一个更复杂的情况：我们接收到的样本中，有 (1-ε) 比例是来自真实分布的真实样本，但有 ε 比例是任意且对抗性的异常值。我们的目标仍然是恢复出接近原始 L̂ 的线性函数。

多项式系统建模

为了应用SOS SDP，我们首先将问题表述为一个多项式系统。

决策变量：

• W： 对于每个样本 i，变量 Wi 表示该样本是否为真实样本（意图是 Wi = 1 表示真实，Wi = 0 表示异常）。
• L： 我们想要寻找的线性函数向量。

多项式系统 P：

1. 二元约束： Wi² - Wi = 0，对于所有 i。这强制 Wi 为 0 或 1。
2. 真实样本数量约束： Σᵢ Wi = (1 - ε)n。这表示至少 (1-ε)n 个样本被标记为真实。
3. 归一化约束： ||L||₂² ≤ 1。这是一个简单的归一化。
4. 最小化目标： 最小化在标记为真实的样本上的平方误差：Minimize Σᵢ Wi (yᵢ - L · xᵢ)²。

这里，xᵢ 和 yᵢ 是已知的样本数据，是多项式中的常数。多项式变量是 W 和 L。

SOS SDP 松弛与解提取

接下来，我们求解这个多项式系统的一个低阶（例如，10阶）SOS SDP松弛。求解后，我们会得到一组“伪期望”值，例如 Ẽ[Wi], Ẽ[Lⱼ], Ẽ[WiWj] 等，它们对应于系统 P 的某个解分布的矩。

如何从中提取解呢？方法非常简单：我们直接输出伪期望值 Ẽ[L] 作为我们找到的线性函数 L。即，L** = Ẽ[L]。

可识别性证明与平方和论证

现在，我们需要证明这个输出的 L* 确实是好的。这涉及到证明一个“可识别性”声明：多项式系统 P 的任何解（或伪期望解）都对应一个在真实样本上误差较小的线性函数。

具体来说，我们希望证明以下声明：对于系统 P 的任何解（或伪期望），在真实样本集 S 上的误差 Σ_{i∈S} (yᵢ - L · xᵢ)² 很小。

我们将通过构造一个平方和（SOS）证明来证明这个声明。SOS证明是一系列遵循特定代数规则（例如，平方项非负、不等式可加、可乘以平方多项式等）的推导。其美妙之处在于，一旦我们为多项式系统 P 构造了一个SOS证明来推导出某个不等式，那么同样的不等式对于SDP松弛得到的伪期望值 Ẽ[·] 也必然成立。

定理：假设 Ẽ[·] 是系统 P 的一个阶数为 O(1) 的伪期望函数子，且 L* = Ẽ[L]。那么，L* 在真实样本上的误差最多比真实线性函数 L̂ 的误差大 O(√ε)。

证明思路（通过SOS论证）：

1. 我们从目标误差表达式开始，并利用恒等式将其分解为两部分：一部分对应被SDP解正确标记的真实样本，另一部分对应被错误标记的真实样本。
2. 对于第一部分（正确标记的部分），我们利用SDP是最小化问题的松弛这一事实，将其误差与真实函数 L̂ 的误差联系起来。
3. 对于第二部分（错误标记的部分），我们应用伪期望版本的柯西-施瓦茨不等式。这是一个可以通过SOS证明的代数不等式。
4. 应用柯西-施瓦茨后，我们得到两个因子的乘积。其中一个因子与错误标记的样本比例 ε 有关，可以证明其大小为 O(ε n)。
5. 另一个因子涉及真实样本上的四阶矩。通过利用样本分布的性质（如 x 是随机±1向量）和SDP中的约束（如 ||L||₂² ≤ 1），我们可以证明这个因子的大小为 O(n)。
6. 综合以上步骤，最终得出误差上界为 L̂ 的误差加上 O(√ε) 项。

通过这个SOS证明，我们确保了从SDP解中提取的 L* 满足所需的误差界限。

---

关于平方和（SOS）能力的讨论

上一节我们看到了如何用SOS SDP解决鲁棒线性回归。本节中，我们来探讨SOS方法的能力与局限。

SOS证明系统非常强大，可以证明许多事实：

• 单变量多项式：一个单变量多项式 p(x) 对所有实数 x 非负（p(x) ≥ 0），当且仅当它可以写成平方和的形式。这意味着对于单变量多项式优化，SOS SDP 是精确的。
• 许多经典不等式：如柯西-施瓦茨不等式、霍尔德不等式、AM-GM不等式等，都有低阶的SOS证明。

然而，也存在一些看似简单的事实，但需要非常高阶的SOS才能证明，甚至可能没有低阶SOS证明：

1. 奇偶性示例：考虑一个包含 (2n-1) 个顶点的完全图。多项式系统要求为每个顶点分配一个关联边变量（意图表示完美匹配），并满足每个顶点恰好关联一条边、变量为0/1等约束。由于顶点数为奇数，完美匹配不存在，系统不可行。但证明此系统的不可行性需要 SOS 的阶数为 Ω(n)。

2. 网格路径不相交示例：在一个 n×n 网格上，要求找到从左下到右上的路径和从左上到右下的路径，且两条路径不相交。直观上，这样的两条路径必然相交。但将这个事实形式化为多项式系统的不可行性，并给出SOS证明，可能需要非常高的阶数（如 n 的量级）。

这些例子表明，尽管SOS框架很强大，能够形式化地推导许多复杂的不等式并设计算法，但它并非万能。对于一些组合结构固有的、基于奇偶性或拓扑的简单矛盾，SOS证明可能需要指数级（高阶）的规模才能捕捉到。

---

总结

本节课中，我们一起学习了平方和（SOS）半正定规划（SDP）的两个核心应用方面。

首先，我们深入探讨了鲁棒线性回归问题。通过将问题建模为多项式系统，应用SOS SDP松弛，并利用SOS证明来论证解的可识别性，我们得到了一个能容忍一定比例异常值的有效算法。关键步骤包括将问题转化为代数形式、利用伪期望的柯西-施瓦茨不等式等SOS可证明的引理进行推导。

其次，我们讨论了SOS方法的能力与局限性。我们了解到SOS可以精确处理单变量多项式非负性问题，并能捕获许多经典不等式。然而，对于一些基于组合或拓扑原理的简单不可行问题，可能需要非常高阶的SOS证明，这揭示了该方法的计算边界。

这种基于SOS SDP和可证明代数推理的范式非常通用，可应用于其他鲁棒估计问题，如鲁棒均值估计或高斯混合模型学习。其核心思想是一致的：将学习问题表述为多项式系统，证明任何（伪）解都接近真实答案，并通过SOS论证使该证明适用于SDP松弛解。

23：线性规划在组合优化中的应用 🧮

在本节课中，我们将要学习线性规划在组合优化中的一个重要应用。我们将从一个经典问题——二分图匹配开始，介绍其线性规划松弛，并证明该松弛是“整的”。接着，我们会将这种方法推广到更复杂的问题，如有度限制的生成树问题，并介绍一种称为“迭代舍入”的强大技术。

二分图匹配的线性规划松弛

上一节我们介绍了课程概述，本节中我们来看看第一个经典例子：二分图匹配问题。

给定一个二分图 G = (V1 ∪ V2, E)，其中 |V1| = |V2| = n，以及边上的权重函数 w(e)。目标是找到一个最小权重的完美匹配 M。

我们可以为这个问题建立一个线性规划松弛。为每条边 e ∈ E 引入一个变量 x_e，其意图是：如果边 e 在匹配中，则 x_e = 1；否则 x_e = 0。约束条件是每个顶点的关联边之和不超过1（即度数约束）。

公式：

• 变量：对于所有 e ∈ E，有 x_e ≥ 0。
• 约束：对于所有顶点 v ∈ V1 ∪ V2，有 ∑_{e incident to v} x_e ≤ 1。
• 目标：最小化 ∑_{e ∈ E} w_e * x_e。

我们称这个线性规划的可行域为二分图匹配多面体。一个关键性质是，这个多面体是整的，这意味着它的所有极点（角点）坐标都是整数。因此，当我们求解这个线性规划时，如果算法返回一个极点解，那么这个解自动对应一个合法的匹配，无需额外的舍入操作。

为了证明这个多面体是整的，我们采用一种算法化的思路：对于任意权重函数 w，我们都可以将线性规划的最优解“舍入”为一个具有相同成本的实际匹配。如果存在一个分数极点，我们总能构造一组权重使其成为唯一最优解，而此时无法找到同等成本的整数匹配，从而产生矛盾。因此，证明的核心在于设计这样一个舍入算法。

以下是该算法的步骤描述：

算法从求解线性规划开始，得到一个极点解 x。

1. 处理零边：如果存在边 e 使得 x_e = 0，则从图中删除边 e，然后重新求解线性规划。
2. 处理整边：如果存在边 e 使得 x_e = 1，则将边 e 加入匹配 M，并从图中删除这条边及其两个端点，然后重新求解线性规划（相应地调整其他顶点的度数约束）。
3. 矛盾情况：如果既没有 x_e = 0 的边，也没有 x_e = 1 的边，即所有 x_e 都是分数（严格介于0和1之间）。我们将证明，对于一个极点解，这种情况不可能发生。

现在，我们来分析为什么“所有边都是分数”的情况会导致矛盾。因为 x 是一个极点解，它必须满足：紧约束（取等号的约束）的数量至少等于变量的数量（即边数 |E|）。由于所有 x_e > 0，非负约束 x_e ≥ 0 都不是紧的。因此，所有紧约束都来自度数约束 ∑ x_e ≤ 1。

设 W 是度数约束取等号（即和为1）的顶点集合。那么有 |W| = |E|（因为紧约束数至少为 |E|）。对于 W 中的任意顶点 v，由于其关联边的 x_e 之和为1，且每个 x_e 都是分数，那么顶点 v 的度数至少为2。因此，所有 W 中顶点的度数 ≥ 2。

考虑边数 |E| 与度数之和的关系：2|E| = ∑_{v} degree(v) ≥ 2|W|。结合 |E| = |W|，可知所有不等式都必须取等号。这意味着：

• 所有 W 中顶点的度数恰好等于2。
• 所有不在 W 中的顶点度数为0。

因此，由 W 中的顶点和分数边构成的子图是一个每个顶点度数均为2的图，即一系列环的并集。然而，在二分图中，环的长度必须是偶数。如果我们分别对左侧顶点 V1 和右侧顶点 V2 的度数约束求和，会得到：
∑_{v in V1} (∑_{e incident to v} x_e) = ∑_{v in V2} (∑_{e incident to v} x_e) = |E|。
这表明这两组约束线性相关，独立紧约束的数量至多为 |W| - 1。这与极点解需要至少 |E| = |W| 个线性独立紧约束相矛盾。因此，情况3不可能发生，算法总能在前两步取得进展，最终得到一个整数匹配。

生成树与有度限制的生成树 🌳

上一节我们证明了二分图匹配多面体的整性，本节中我们来看看另一个经典问题：生成树。我们为其定义一个称为“子图消除”的线性规划。

对于图 G=(V,E)，为每条边 e 引入变量 x_e，意图是 x_e=1 表示边在生成树中。约束条件有两个部分：

1. 对于任意顶点子集 S ⊆ V, S ≠ ∅，有 ∑_{e ∈ E(S)} x_e ≤ |S| - 1，其中 E(S) 是 S 内部的边集。这保证了无环。
2. ∑_{e ∈ E} x_e = |V| - 1。这保证了边数正确。

公式：

• 对于所有 S ⊆ V, S ≠ ∅：∑_{e ∈ E(S)} x_e ≤ |S| - 1
• ∑_{e ∈ E} x_e = |V| - 1
• 对于所有 e ∈ E：x_e ≥ 0

这个多面体也是整的，即其极点对应着图的生成树。现在，我们考虑一个更难的问题：有度限制的生成树。除了生成树的基本约束外，我们还对某些顶点 v 的度数施加限制：deg_T(v) ≤ B_v。

在线性规划中，这转化为添加约束：∑_{e incident to v} x_e ≤ B_v（对于指定的顶点集合 W 中的 v）。这个问题是NP难的，因为当 B_s = B_t = 1 且其他顶点 B_v = 2 时，它就等价于哈密顿路径问题。

然而，我们可以设计一个近似算法。我们证明存在一个多项式时间算法，可以找到一棵生成树 T，其成本不超过线性规划最优值，但每个顶点的度数至多为 B_v + 2（后续改进可达到 B_v + 1）。算法采用了迭代舍入的技术。

算法的框架如下，我们不断求解当前图的线性规划，得到一个极点解 x，然后根据 x 的结构做出决定，修改问题实例，并迭代进行。

1. 删除零边：若存在边 e 使 x_e = 0，则从图 G 中删除 e。
2. 处理叶子顶点：若存在顶点 v 只有一条关联边 e 满足 x_e > 0（即在该解中度数为1），则将边 e 加入树 T，从图中删除顶点 v，并将边 e 另一端点 u 的度限制 B_u 减1。
3. 放松度限制：若存在顶点 v（其度限制 B_v 原本存在）在当前解中的度数（即关联边 x_e 之和）至少为3，则直接删除该顶点 v 的度限制约束 ∑ x_e ≤ B_v。
4. 重复上述过程，直到构建出完整的生成树。

算法的关键在于证明，我们永远不会陷入一种“僵局”，即当前极点解 x 不满足以上任何一条可操作的条件。假设陷入僵局，则意味着：没有边为0，没有顶点度数为1，且所有有度限制的顶点 v ∈ W 的度数都至少为4（否则会被条件3处理），其他顶点度数至少为2。通过计算边数并分析紧约束的数量，可以导出矛盾。

紧约束主要来自两部分：子图消除约束和度限制约束。对于子图消除约束，一个关键引理是：如果两个集合 A 和 B 对应的约束是紧的，那么 A ∪ B 和 A ∩ B 对应的约束也是紧的。通过这种“非交叉化”操作，可以证明所有紧的子图消除约束张成的空间，等价于一个层叠族（Laminar Family）的约束所张成的空间。一个不含单元素集的层叠族最多包含 n-1 个集合。因此，线性独立的子图消除紧约束最多有 n-1 个。

在僵局假设下，设 |W|=k，总边数 |E| ≥ n + k。作为一个极点解，至少需要 |E| 个线性独立紧约束。然而，我们最多只有 (n-1) 个来自层叠族的独立紧约束，加上最多 k 个度限制紧约束，总数最多为 n - 1 + k，这小于 n + k。矛盾产生。因此，算法永远不会陷入僵局，总能不断推进直至得到生成树。

迭代舍入与顶点覆盖 🛡️

上一节我们看到了迭代舍入技术在有度限制生成树问题上的应用，本节中我们简要探讨该方法在另一个经典问题——顶点覆盖上的体现，以展示其通用性。

顶点覆盖问题的标准线性规划松弛如下：为每个顶点 v 引入变量 x_v。
公式：

• 对于所有边 (u,v) ∈ E：x_u + x_v ≥ 1
• 对于所有顶点 v ∈ V：0 ≤ x_v ≤ 1
• 目标：最小化 ∑_{v} w_v * x_v

一个简单的迭代舍入算法是：

1. 求解线性规划，得到极点解 x。
2. 若存在顶点 v 满足 x_v ≥ 1/2，则将 v 加入顶点覆盖，并从图中删除 v 及其所有关联边。
3. 重复上述过程。

可以证明，在极点解中，总存在某个变量的值至少为 1/2（否则所有约束无法线性独立到足以定义极点）。因此，算法总能进行。每次我们将一个 x_v ≥ 1/2 的顶点加入覆盖，其成本最多是线性规划中该顶点所贡献成本的两倍，从而最终得到一个成本不超过线性规划最优值两倍的整数解。这实际上就是经典的2-近似算法。迭代舍入的观点为这个经典算法提供了一个新的证明角度。

总结 📝

本节课中我们一起学习了线性规划在组合优化中的核心应用。我们从二分图匹配的线性规划松弛出发，证明了其多面体的整性，并引入了通过分析极点解紧约束数量来进行“舍入”的证明思想。接着，我们将问题推广到生成树及其带度限制的变体，介绍了迭代舍入这一强大技术。该技术通过反复求解线性规划、固定那些接近整数的变量，并利用极点解的紧约束数量必须足够多这一维度论证，来保证算法进度并得到近似解。最后，我们看到了迭代舍入观点如何自然地解释顶点覆盖的2-近似算法。这些内容展示了线性规划作为工具，在精确求解和近似算法设计中的深远影响。

24：匹配问题的高级算法

在本节课中，我们将要学习匹配问题的两种高级算法：基于线性规划的多面体方法和基于多项式检验的代数方法。我们将看到如何通过扩展约束来精确描述一般图的完美匹配，以及如何利用行列式和随机赋值来高效地检测匹配的存在性。

多面体方法：完美匹配的线性规划

上一节我们介绍了二分图匹配的线性规划。本节中，我们来看看如何将线性规划推广到一般图的完美匹配问题。

对于一个图 G=(V, E)，我们希望找到一组边，使得每个顶点恰好与一条边关联，且边之间不共享顶点，这称为完美匹配。我们定义变量 x_e，意图是当边 e 在匹配中时 x_e = 1，否则为 0。

首先，我们写出基本的线性规划约束：

• x_e ≥ 0，对于所有边 e ∈ E。
• 对于每个顶点 v ∈ V，其关联边的变量之和为 1：∑_{e ∈ δ(v)} x_e = 1。

然而，对于一般图，仅这些约束是不够的。例如，在一个三角形（三个顶点的环）中，给每条边赋值 1/2 满足上述约束，但三角形显然没有完美匹配。因此，我们需要引入额外的约束，称为奇集不等式。

以下是奇集不等式的定义：

• 对于任意顶点子集 U ⊆ V，如果 |U| 是奇数，则至少有一条匹配边离开该集合。用公式表示为：∑_{e ∈ δ(U)} x_e ≥ 1。

这个不等式是直观的：如果一个奇数大小的顶点集 U 中的所有顶点都在内部相互匹配，由于奇数的限制，至少会有一个顶点无法在内部找到配对，因此必须有一条边连接到集合外部。

现在，我们得到了定义完美匹配多面体 P 的三组约束：

1. 非负性：x_e ≥ 0， ∀ e ∈ E。
2. 度约束：∑_{e ∈ δ(v)} x_e = 1， ∀ v ∈ V。
3. 奇集不等式：∑_{e ∈ δ(U)} x_e ≥ 1， ∀ U ⊆ V， |U| 为奇数。

一个关键定理是：这个多面体 P 的所有极值点恰好对应图的完美匹配。这意味着，如果我们在这个多面体上优化一个线性目标函数（例如，最大权重匹配），得到的最优解将直接是一个完美的（或最大权重的）匹配，而不会是分数解。

定理证明思路（迭代舍入法）

我们可以使用“最小反例法”和“迭代舍入”的思路来证明这个定理。

1. 假设存在最小反例：假设图 G 是满足以下条件的最小图（按顶点数加边数衡量）：其多面体 P 存在一个极值点 x，而 x 不是一个匹配的指示向量（即不是 0-1 解）。
2. 分析极值点性质：
   • 首先，x_e 不可能为 0。如果某条边 x_e = 0，我们可以删除这条边，得到一个更小的反例图。
   • 其次，每个顶点的度至少为 2。如果某个顶点度为 1，那么关联的唯一那条边 e 的 x_e 必须为 1（以满足度约束）。这样，我们可以将这条边及其两个顶点从图中移除，得到一个更小的反例。
3. 应用紧约束：由于 x 是极值点，它必须使至少 |E| 个约束取等号（紧）。度约束有 |V| 个。因此，除了所有度约束可能都是紧的之外，至少还有一个奇集不等式是紧的。即存在一个奇数大小的顶点集 U，使得 ∑_{e ∈ δ(U)} x_e = 1。
4. 收缩与归纳：考虑将集合 U 收缩为一个超级顶点，得到图 G1；同时，将集合 V\U 收缩为一个超级顶点，得到图 G2。这两个图都比原图 G 小。
   • 将解 x 分别限制到 G1 和 G2 上，得到 x¹ 和 x²。可以验证它们是 G1 和 G2 对应多面体的可行解。
   • 由于 G1 和 G2 更小，根据最小反例假设，它们对应的多面体是“好的”（即极值点都是匹配）。因此，x¹ 和 x² 都可以表示为各自图中一组完美匹配的凸组合。
5. 合并匹配：关键点在于，因为原图中跨越割 (U, V\U) 的边变量之和为 1，所以在 G1 和 G2 的匹配凸组合中，每一个匹配都恰好包含一条跨越收缩后割的边。这样，我们可以将 G1 中选取某条跨越边 e 的匹配，与 G2 中选取同一条边 e 的匹配“粘贴”起来，从而将 x¹ 和 x² 的凸组合合并为原图 G 中完美匹配的一个凸组合。
6. 得出结论：这意味着原图的极值点 x 本身也是完美匹配的凸组合，这与 x 是极值点（不能是凸组合除非是单个点）相矛盾。因此，假设不成立，多面体 P 的所有极值点都是完美匹配。

这个证明展示了组合优化中一个强大的技巧：利用线性规划的极值点结构和组合问题的特性，通过归纳和收缩来证明规划的精确性。

算法意义与备注

• 这个线性规划虽然约束数量是指数级的（所有奇子集），但存在多项式时间的分离预言机。这意味着我们可以使用椭球法在多项式时间内求解它。
• 一个有趣的事实是：对于完美匹配多面体，不存在变量和约束都是多项式规模的线性规划刻画（除非某些复杂性假设不成立）。这与生成树多面体不同，后者可以通过增加流变量来构造多项式规模的线性规划。这使得匹配多面体在某种程度上更为复杂和有趣。

---

代数方法：通过行列式检测匹配

上一节我们探讨了使用线性规划精确求解匹配问题。本节中，我们转向一种完全不同的、基于代数和随机化的方法，用于检测匹配的存在性。

二分图匹配的代数检测

给定一个二分图 G=(V1 ∪ V2, E)，我们想知道它是否存在完美匹配。

我们构造一个 |V1| × |V2| 的矩阵 B，其元素定义如下：

• 如果 (i, j) ∈ E，则 B_{ij} = x_{ij}（一个变量）。
• 否则，B_{ij} = 0。

现在，考虑这个矩阵的行列式 det(B)。根据行列式的定义：
det(B) = ∑_{π ∈ S_n} sign(π) ∏_{i=1}^{n} B_{i, π(i)}
其中 S_n 是所有排列的集合。

观察这个多项式：

• 如果一个排列 π 对应的边 (i, π(i)) 不全在图中（即某个 B_{i, π(i)} = 0），那么该项为 0。
• 因此，只有那些使得所有 (i, π(i)) 都是图中边的排列 π 才会贡献非零项。这样的排列 π 恰好定义了从 V1 到 V2 的一个完美匹配。

所以，det(B) 作为一个多项式，是非零多项式当且仅当图 G 存在完美匹配。

算法思路：我们不需要显式展开这个可能包含指数多项式的行列式。相反，我们可以使用随机赋值的方法。

1. 为每个变量 x_{ij} 随机选择一个值（例如，从一个足够大的整数集合中随机选取）。
2. 计算赋值后数值矩阵 B 的行列式（这是一个数字，可以在 O(n³) 时间内用高斯消元法计算）。
3. 如果行列式结果非零，则输出“存在完美匹配”；如果为零，则输出“不存在完美匹配”。

这个算法的正确性依赖于一个重要的引理——施瓦茨-齐佩尔引理：

对于一个非零的 n 元多项式 P，如果每个变量从一个有限集合 S 中随机独立选取，那么 P 求值为 0 的概率不超过 deg(P) / |S|，其中 deg(P) 是多项式的总次数。

在我们的例子中，det(B) 的次数最多为 n。因此，如果我们从大小为 poly(n) 的集合中随机赋值，算法出错的概率可以控制在很小的范围内（例如 1/poly(n)）。这是一个单边错误的蒙特卡洛随机算法：如果图有完美匹配，它几乎总能发现；如果它说“没有”，那一定没有。

从检测到搜索：上述算法只能检测是否存在完美匹配。要实际找到一个完美匹配，我们可以使用标准的“搜索归约到检测”框架：依次考虑每条边 e，从图中移除 e 后检查是否仍存在完美匹配。如果不存在，则 e 必在某个完美匹配中，将其加入匹配集并对剩余图递归操作；如果存在，则丢弃 e。这需要多次调用检测算法，但通过适当的错误概率控制，我们可以在多项式时间内以高概率找到一个完美匹配。

一般图匹配的代数检测

对于一般图 G=(V, E)，我们也可以构造一个矩阵 T 来检测完美匹配。T 是一个 |V| × |V| 的斜对称矩阵（或称为 Tutte 矩阵），定义如下：

• 如果 (i, j) ∉ E，则 T_{ij} = 0。
• 如果 (i, j) ∈ E，则：
  • 若 i < j，令 T_{ij} = x_{ij}。
  • 若 i > j，令 T_{ij} = -x_{ji}。

关键定理是：图 G 存在完美匹配当且仅当 det(T) 作为多项式是非零的。

证明思路：

1. （必要性） 如果存在完美匹配 M，考虑一个特殊的排列 π：对于匹配边 (i, j)，令 π(i)=j 且 π(j)=i。这个排列由若干个对换（长度为 2 的环）构成。计算 det(T) 中对应 π 的项，会发现它包含形如 ∏_{(a,b)∈M} (x_{ab})² 的单项式。由于矩阵的斜对称性，任何其他排列都无法产生完全相同的这个单项式，因此它不会被抵消，det(T) 非零。
2. （充分性） 如果 det(T) 非零，则存在某个排列 π 使得对应项非零。将排列 π 分解为环。可以证明，如果 π 包含一个奇数长度的环，那么通过“反转”这个奇环可以得到另一个排列 π‘，且 det(T) 中 π 和 π’ 对应的项相互抵消。因此，在非零项对应的排列 π 中，所有环的长度必须是偶数。而一个所有环都是偶数的排列，其边集自然包含了若干个偶环，从每个偶环中可以很容易地选出一个完美匹配（间隔选边即可）。

因此，我们同样可以通过随机赋值给变量 x_{ij}，然后计算数值矩阵 T 的行列式，来以高概率检测一般图是否存在完美匹配。

---

应用：红黑匹配问题

上一节的方法展示了代数技术的威力。本节中，我们来看一个更复杂的变体问题，它目前只有基于代数方法的随机算法。

问题描述：给定一个二分图，每条边被染成红色或黑色。再给定一个整数 k。目标是判断是否存在一个完美匹配，其中恰好包含 k 条红边。

这是一个非常自然的约束匹配问题，但令人惊讶的是，我们不知道任何解决该问题的确定性多项式算法。

随机算法思路：
我们修改之前二分图的矩阵 B 的定义：

• 如果 (i, j) 是黑边，则 B_{ij} = x_{ij}。
• 如果 (i, j) 是红边，则 B_{ij} = y * x_{ij}。这里我们引入了一个新变量 y。
• 如果 (i, j) 不是边，则 B_{ij} = 0。

现在，计算行列式 det(B)。将其视为变量 x_{ij} 和 y 的多项式。我们可以将其按 y 的幂次展开：
det(B) = ∑_{j=0}^{n} y^j * Q_j(x)
其中 Q_j(x) 是一个只关于 x_{ij} 的多项式。

关键观察：多项式 Q_k(x) 非零当且仅当图中存在一个恰好包含 k 条红边的完美匹配。因为 y^j 的系数 Q_j(x) 本质上汇总了所有恰好包含 j 条红边的完美匹配的贡献。

算法步骤：

1. 为所有 x_{ij} 变量随机赋值（例如，从一个大集合中随机选数），得到数值矩阵 B(α)。
2. 现在 det(B(α)) 变成了一个单变量 y 的多项式。我们计算这个多项式在多个 y 值（例如 y = 1, 2, ..., n+1）上的值。
3. 通过多项式插值，我们可以恢复出 det(B(α)) 作为 y 的多项式的各个系数，特别是系数 Q_k(α)。
4. 如果 Q_k(α) ≠ 0，则输出“存在恰好 k 条红边的完美匹配”；否则输出“不存在”。

根据施瓦茨-齐佩尔引理，如果原问题的答案是“存在”（即 Q_k(x) 非零），那么随机赋值后 Q_k(α) 为零的概率很小。因此，这是一个正确的随机算法。

这个例子凸显了代数方法的灵活性：通过引入额外的变量（如 y），我们可以让行列式编码更复杂的组合信息（如红边的数量）。

---

核心挑战：多项式恒等检验

上述所有代数算法的核心都归结为一个更基础的问题：多项式恒等检验。

问题定义：给定一个由某个算术电路（或公式，如行列式）计算出的多元多项式 P(x1, ..., xn)，判断这个多项式是否恒等于零（即展开后所有项是否完全抵消）。

我们有一个非常简单的随机算法：为所有变量随机赋值，然后计算电路输出。如果结果非零，则 P 肯定非零；如果结果为零，则猜测 P 恒为零。根据施瓦茨-齐佩尔引理，当 P 非零时，随机赋值得到零的概率很低。

然而，为 PIT 问题寻找一个确定性多项式时间算法，是理论计算机科学中一个长期悬而未决的重大开放问题。如果这个问题被解决，它将直接给出红黑匹配等问题的确定性算法，并且会带来复杂性理论领域的许多突破性进展（例如电路下界）。目前，我们相信 PIT 应该有确定性算法（因为随机算法通常可以去随机化），但尚未找到。

---

总结

本节课中我们一起学习了匹配问题的两种高级视角：

1. 多面体方法：通过添加奇集不等式，我们得到了完美匹配的精确线性规划刻画。虽然规划规模是指数级的，但存在多项式时间的分离预言机。我们使用迭代舍入和最小反例法证明了该多面体的整数性。
2. 代数方法：通过为图的邻接矩阵（或 Tutte 矩阵）的边引入变量，其行列式多项式编码了完美匹配的信息。利用施瓦茨-齐佩尔引理和随机赋值，我们可以在多项式时间内以高概率检测匹配的存在性。这种方法非常灵活，可以用于解决像红黑匹配这样的约束匹配问题。
3. 基础问题：这些代数算法的核心依赖于多项式恒等检验问题，该问题的确定性求解目前仍是开放的挑战。

这些内容展示了组合优化与代数、概率方法之间深刻而美妙的联系。

25：属性测试入门 🧪

在本节课中，我们将学习一个名为“属性测试”的算法与复杂性领域。属性测试的核心思想是，我们希望通过随机检查输入的一小部分，来快速判断一个大型对象（如列表、函数或图）是否满足某个特定属性，或者是否“远离”满足该属性。这对于处理大数据或需要快速预检查的场景非常有用。

---

什么是属性测试？ 🤔

属性测试的基本设定如下：我们有一个对象，例如一个长度为 n 的列表或字符串。我们想知道这个对象是否满足某个属性 P（例如，列表是否已排序）。通常，完全验证属性需要读取整个输入，时间复杂度为 O(n)。而属性测试的目标是设计一个随机化算法，它只进行 o(n) 次查询（例如 O(log n) 次），就能以高概率做出正确判断。

具体来说，我们期望的算法保证是：

• 如果对象满足属性 P，则测试算法以至少 0.99 的概率“接受”。
• 如果对象“远离”满足属性 P（例如，需要修改很多元素才能使其满足属性），则测试算法以至少 0.99 的概率“拒绝”。

这里的关键在于“远离”的定义。对于几乎满足属性的对象（例如只有一个元素错位的“几乎排序”列表），我们允许算法犯错，因为我们只检查了输入的一小部分。我们只对明显违反属性的情况做出强力保证。

---

示例一：测试列表是否“接近排序” 📊

上一节我们介绍了属性测试的概念，本节中我们来看看一个具体例子：测试一个数字列表是否“接近”已排序状态。

问题定义

输入是一个包含 n 个数字的列表 A[1..n]。我们希望测试该列表是否“接近”已排序。首先，我们需要定义“接近”的含义。

一个递增子序列是指从原列表中按顺序取出的一些元素，它们满足递增关系。列表的最长递增子序列长度是衡量其“有序程度”的好方法。

我们定义：一个序列是 l-远离排序 的，如果其最长递增子序列的长度为 n - l。也就是说，我们需要修改至少 l 个元素才能使其完全排序。我们的目标是设计一个算法，能够区分“列表已排序”（l=0）和“列表是 εn-远离排序”（l ≥ εn）这两种情况。

初步尝试与失败

首先，我们可能会尝试一些简单的检查方法。

以下是两种直观但会失败的测试方法：

1. 随机检查连续对：随机选择一个索引 i，检查是否 A[i] < A[i+1]。如果列表是 1,2,...,n/2,1,2,...,n/2，它距离排序很远（最长递增子序列长度仅为 n/2），但此测试只有 1/n 的概率能发现那个唯一的逆序对。
2. 随机检查任意对：随机选择两个索引 i < j，检查是否 A[i] < A[j]。可以构造一个序列（例如分成 √n 块，每块内递减，但块之间递增），使得该序列远离排序，但随机检查两个相距较远的元素时，它们看起来是递增的，因此测试失败。

有效算法：基于二分搜索的测试

现在，我们来看一个真正有效的属性测试算法。该算法非常优雅。

算法描述：
重复以下步骤 O(1/ε) 次：

1. 在 {1, ..., n} 中均匀随机选择一个索引 i。
2. 尝试在列表 A 中执行标准的二分搜索算法来查找值 A[i]。
3. 在执行二分搜索的过程中，算法会查询一系列位置并与 A[i] 比较。如果在此过程中发现任何违反“搜索路径上值应有序”的情况（例如，在应该更小的位置遇到了更大的值），则立即拒绝（判定列表不接近排序）。
   如果所有重复测试都未触发拒绝，则接受列表。

该算法总共进行 O((log n)/ε) 次查询。

算法正确性证明

为什么这个算法有效？我们来证明其核心思想。

首先，定义一个“好元素”：元素 A[i] 是好的，如果针对 A[i] 执行二分搜索时，没有发现任何违规。否则，它是坏的。

引理 1：如果算法以高概率（如 0.99）接受，则好元素的数量至少为 (1 - ε)n。

• 证明：如果坏元素超过 εn 个，那么每次随机选择一个索引 i 时，有至少 ε 的概率选到坏元素。一旦选到坏元素，二分搜索必定会发现违规并导致拒绝。重复 O(1/ε) 次后，算法仍然接受的概率极低（小于 (1-ε)^(O(1/ε)) ≈ e^{-O(1)}），与高概率接受的假设矛盾。

引理 2：所有好元素构成了原列表的一个递增子序列。

• 证明：任取两个好元素 A[i] 和 A[j]，且 i < j。考虑对它们分别进行二分搜索。这两个搜索路径在二分搜索树中会有一个最低公共祖先节点，设其索引为 k。在 A[i] 的搜索路径上，由于没有违规，且路径在 k 之后转向左子树，可知 A[i] < A[k]。在 A[j] 的搜索路径上，路径在 k 之后转向右子树，可知 A[k] < A[j]。因此，A[i] < A[j]。所以所有好元素是递增的。

结合两个引理：如果算法高概率接受，则存在一个长度至少为 (1 - ε)n 的递增子序列（即所有好元素）。根据定义，这意味着列表至多是 εn-远离排序的。反之，如果列表是 εn-远离排序的（即最长递增子序列长度小于 (1 - ε)n），那么好元素不可能太多，算法有很大概率在重复测试中选到坏元素并拒绝。

---

示例二：估计最大匹配大小 🕸️

上一节我们看到了一个关于列表属性的测试算法，本节中我们来看看一个关于图属性的例子：估计一个极大匹配的大小。这个问题可以等价于为顶点覆盖问题提供一个常数时间的 2-近似算法。

问题定义

输入是一个最大度不超过 d 的图 G（d 被视为常数）。一个匹配是边集的一个子集，其中任意两条边没有公共顶点。一个极大匹配是指不能再添加任何边而不破坏匹配性质的匹配。注意，极大匹配的大小不是唯一的，它取决于构造匹配的顺序。我们的目标是输出一个数字 t，使得存在某个极大匹配 M，其大小 |M| 满足 |t - |M|| ≤ εn。

基础构造算法

如果我们不关心时间复杂度，一个标准的随机化贪心算法可以构造一个极大匹配：

1. 为所有边生成一个均匀随机排列 π（即给每条边分配一个 1 到 |E| 的唯一随机优先级）。
2. 按照 π 中优先级从高到低（或从低到高）的顺序扫描边。
3. 对于每条边，如果它能被加入当前匹配（即其两个端点都未被匹配），则加入；否则跳过。

这个算法必然产生一个极大匹配，记作 M。我们的目标是不显式构造 M，而是快速估计 |M|。

估计思路与核心挑战

一个自然的估计思路是：

1. 随机采样 O(1/ε^2) 条边。
2. 对于每条采样边 e，判断它是否属于上述算法产生的极大匹配 M。
3. 计算属于 M 的采样边的比例，乘以总边数，即可估计 |M|。

核心挑战在于：如何在不运行完整贪心算法的情况下，快速判断一条给定边 e 是否在 M 中？我们需要一个期望常数时间（关于 d 和 1/ε）的查询算法。

常数时间成员查询算法

以下是判断边 e 是否在 M 中的算法 Test(e)：

算法 Test(e)：

1. 检查边 e 的所有邻边（即共享端点的边）。如果存在一条邻边 e' 满足 π(e') < π(e)（即 e' 的优先级更高，在贪心算法中更早被考虑），那么 e 能否加入 M 取决于 e' 是否被加入了 M。
2. 因此，我们需要递归地调用 Test(e') 来判断 e' 的状态。
3. 如果所有优先级更高的邻边 e' 都未被加入 M（即 Test(e') 返回 False），那么当贪心算法扫描到 e 时，其端点都空闲，e 会被加入，故 Test(e) 返回 True。
4. 如果存在某个优先级更高的邻边 e' 已加入 M（即 Test(e') 返回 True），那么 e 的某个端点已被占用，e 无法加入，故 Test(e) 返回 False。
5. 递归的基准情况：如果某条边 e* 在其局部邻域中具有最小的 π 值（即没有优先级更高的邻边），那么 Test(e*) 直接返回 True。

运行时间分析

关键在于分析 Test(e) 的期望查询次数。当 π 是随机排列时：

• 为了查询边 e，算法可能需要沿着一条路径 (e, e1, e2, ..., ek) 进行递归，其中每条边的 π 值严格递减。
• 对于一条长度为 k 的路径，π 值在该路径上严格递减的概率是 1/k!（所有 k! 种顺序等可能）。
• 在最大度为 d 的图中，从 e 出发长度为 k 的路径（按边计）最多有约 2 * d^k 条。
• 因此，Test(e) 访问的总边数的期望值最多为：
  ∑_{k≥1} (2 * d^k) * (1/k!) = 2 * (e^d - 1) = O(e^d)
  这是一个只依赖于常数 d 的常数。

因此，判断单条边是否在 M 中的期望时间是常数。我们总共需要判断 O(1/ε^2) 条采样边，故总期望运行时间为 O(e^d / ε^2)，对于常数 d 和 ε，这也是常数。

---

总结与拓展 🎓

本节课中我们一起学习了属性测试的基本概念和两个经典示例。

1. 我们首先了解了属性测试的目标：用次线性时间的随机查询，判断对象是否满足属性或远离该属性。
2. 然后，我们深入探讨了测试列表接近排序性的算法，该算法巧妙地利用二分搜索路径的一致性进行检测。
3. 接着，我们研究了估计图中极大匹配大小的问题，并给出了一个常数时间的近似算法。其核心是设计了一个高效的“成员查询”过程，通过分析随机排列下递归查询的期望代价证明了其高效性。

属性测试与许多其他领域紧密相关，例如：

• 图属性测试：测试图是否接近二分图、是否无三角形等。基于西泽雷吉引理，人们对哪些图属性可以常数时间测试有了深入的理解，尽管所需常数可能极大。
• 概率可检查证明：这是属性测试的重要起源之一，研究如何用极少量的随机查询来验证数学证明的正确性。
• 委托计算验证：如何让计算能力弱的验证者高效检查强大云服务器的计算结果是否正确。

属性测试展示了在面对海量数据时，随机化和近似计算所能带来的巨大效率提升。

26：Lovász局部引理及其构造性算法

在本节课中，我们将学习Lovász局部引理，并重点探讨一个用于K-SAT问题的构造性算法。我们将从一个具体的例子入手，逐步理解算法的核心思想及其巧妙的分析过程。

概述

Lovász局部引理是一个概率论中的重要工具，它用于证明在多个“坏事件”并非完全独立，但每个事件仅依赖于少数其他事件的情况下，所有坏事件都不发生的概率大于零。本节课我们将通过K-SAT问题，学习该引理的一个构造性证明版本，即不仅证明解的存在性，还能通过一个高效的算法找到这个解。

K-SAT问题定义

首先，我们回顾K-SAT问题。给定n个布尔变量 X1, X2, ..., Xn 和m个子句 C1, C2, ..., Cm。每个子句是k个文字的析取（逻辑或）。例如，一个子句可能形如 Xi1 ∨ ¬Xi2 ∨ Xi3 ∨ ... ∨ Xik。目标是找到一个变量的赋值，使得所有子句都被满足（即为真）。

对于一个随机赋值，任意一个特定的K-SAT子句被满足的概率是 1 - 1/(2^k)，因为只有唯一的一种赋值（所有文字为假）会使其不满足。

核心定理

我们旨在证明以下定理：

对于任意一个K-SAT公式，如果每个变量参与的子句数量小于 2^k / k（减去一个常数，如20），那么该公式必然是可满足的，并且可以在多项式时间内找到一个满足的赋值。

这里的“变量度”指的是任意一个变量出现在不同子句中的最大次数。

算法描述

上一节我们介绍了定理的目标，本节中我们来看看证明该定理的算法。这是一个非常简洁的算法。

算法伪代码如下：

输入：一个K-SAT公式 φ
1. 随机生成一个初始赋值 X。
2. While (存在被违反的子句 C)：
      调用 Fix(C) 函数。

Fix(C) 函数定义如下：

Fix(子句 C):
   1. 对子句C中的所有变量，重新随机赋值。
   2. While (存在与C共享变量的、被违反的子句 D)：
          调用 Fix(D)。

算法从随机赋值开始，然后递归地“修复”那些不满足的子句。修复一个子句意味着将其包含的所有变量重新随机化。这个过程可能会破坏其邻近子句的满足性，因此算法会递归地修复那些新产生的违反子句。

算法正确性分析：第一部分

算法的分析分为两部分。首先，我们证明算法的外层循环（第2步）执行次数有限。

关键观察： 在每次调用 Fix(C) 并返回后，整个公式中被满足的子句数量至少增加一个。

原因如下：
考虑任意一个在 Fix(C) 调用前就被满足的子句D。在 Fix(C) 的执行过程中，D的状态可能从满足变为不满足。但是，算法规定，一旦一个子句因邻居子句的修改而变为不满足，算法会立即递归调用 Fix(D) 来修复它。因此，在 Fix(C) 最终返回时，D必然是满足的。同时，最初触发修复的子句C在修复后也必然被满足（尽管中间可能又被破坏，但最终会被修复）。因此，每次成功执行完一次外层循环（即调用并返回一次Fix），至少有一个子句从永久不满足变为永久满足。

由此可得，外层循环最多执行m次（m为子句总数）。因此，只要内层的递归过程 Fix 能够终止，整个算法就会在多项式步数内结束。

算法正确性分析：第二部分

接下来，我们需要证明递归函数 Fix 本身不会无限运行下去。我们将采用反证法，并借助一个精巧的“通信游戏”思想。

证明思路：
假设算法运行了非常多的步骤S仍未终止。我们将模拟算法与一个“调试器”之间的单向通信。

1. 算法发送给调试器的信息包括：
   • 每次调用 Fix(C) 时，发送子句C的标识。
   • 每次 Fix(C) 返回时，发送一个“返回”信号。
   • 在算法运行了S步后（或假设的终止点），发送当前的完整赋值Y。
2. 调试器的能力： 在收到以上信息后，调试器可以逆向推导出算法的全部随机比特，包括初始赋值X和每次重随机化子句时所用的随机数。

逆向推导如何工作？
调试器从最终赋值Y开始，按照调用 Fix 的逆序进行处理。关键点在于：算法只在子句被违反时才重随机化它。对于一个K-SAT子句，有且仅有一种赋值会违反它。因此，当调试器知道一个子句C被重随机化，并且知道重随机化后的状态，它就能唯一确定重随机化之前该子句的变量赋值（即那个唯一的违反赋值）。通过这种方式，调试器可以一步步倒退回初始状态X，并恢复出所有中间步骤使用的随机比特。

通信比特数计算：
现在我们来计算算法发送了多少比特，以及调试器恢复出了多少比特。

• 算法发送的比特数：
  • 最终赋值Y：n 比特。
  • S次调用/返回的信号：O(S) 比特。
  • 子句标识：
    • 对于外层循环调用的子句（最多m个），每个需要 log m 比特来标识。
    • 对于递归调用的子句，由于它一定是当前正在处理的子句的邻居，所以只需用 log d 比特来指明是哪个邻居即可，其中 d 是一个子句的最大邻居数（与变量度相关）。
  • 总发送比特数约为：n + O(S) + m*log m + S*log d
• 调试器恢复的比特数：
  • 初始赋值X：n 比特。
  • S次子句重随机化：每次需要 k 比特随机数来给k个变量赋值。
  • 总恢复比特数为：n + S*k

导出矛盾：
如果算法运行了非常多的步骤S，使得 S*k（恢复的比特数主要增长项）显著大于 S*log d（发送比特数的主要增长项），我们就得到了一个矛盾：调试器用较少的通信量恢复出了更多的信息。这违反了信息论的基本原理。

具体来说，当 k > log d + O(1) 时，对于足够大的S，就会产生矛盾。在我们的设定中，变量度小于 2^k / k，因此子句的邻居数 d 也小于约 2^k，从而 log d < k。条件满足，故假设不成立，算法必然在多项式步数内终止。

与Lovász局部引理的联系

上一节我们分析了一个具体的构造性算法，本节中我们来看看它背后的普遍原理——Lovász局部引理（LLL）。

LLL的非构造性表述：
考虑一个概率空间和一组“坏事件” E1, E2, ..., Em。假设每个坏事件发生的概率至多为 p，且每个坏事件最多与d个其他坏事件相关（即独立于其他事件）。如果满足：
p * e * (d + 1) < 1 （其中 e 是自然常数）
那么，所有坏事件都不发生的概率大于零。

在K-SAT的例子中，坏事件是“子句i被违反”，其概率 p = 1/(2^k)，依赖数 d 与变量度相关。LLL的非构造性证明直接保证了满足条件的公式必有解，但没有给出找解的方法。

构造性与非构造性证明：
我们之前学习的算法，为K-SAT这一特殊情况提供了LLL的一个构造性证明。这意味着它不仅证明了解的存在，还给出了一个高效的算法来实际找到这个解。这是算法研究中的一个美妙成果，因为很多存在性定理（如纳什均衡的存在性）并不容易转化为构造性算法。

总结

本节课我们一起学习了Lovász局部引理及其在K-SAT问题上的一个精彩应用。我们首先描述了一个简单的递归随机算法，然后通过一个基于信息论思想的巧妙分析，证明了该算法能够在多项式时间内找到低度K-SAT公式的解。最后，我们将此结果与更一般的、非构造性的Lovász局部引理联系起来，理解了构造性证明的价值所在。这个算法展示了概率方法与算法设计结合所能产生的强大力量。