Isotonic Regression（保序回归）的数学原理与拟合过程
字数 3588 2025-12-07 05:08:19

Isotonic Regression（保序回归）的数学原理与拟合过程

Isotonic Regression（保序回归）是一种非参数回归方法，其核心思想是在保持预测值单调性（非递减或非递增）约束的前提下，最小化预测值与观测值之间的误差。它不假设数据服从特定的函数形式（如线性），但要求拟合的函数是单调的。下面我将详细讲解其数学原理、问题描述、求解算法（PAVA）以及计算过程。

1. 问题描述

假设我们有 \(n\) 个观测数据点 \((x_1, y_1), (x_2, y_2), \dots, (x_n, y_n)\)，其中 \(x_i\) 是自变量（通常有序，如时间、剂量），\(y_i\) 是因变量（如响应值）。我们假设 \(x_i\) 已经按升序排列（即 \(x_1 \leq x_2 \leq \dots \leq x_n\)）。如果 \(x_i\) 相等，则需特殊处理。

Isotonic Regression 的目标是找到一组拟合值 \(\hat{y}_1, \hat{y}_2, \dots, \hat{y}_n\)，满足以下约束和优化目标：

  • 单调性约束：拟合值必须是非递减的，即 \(\hat{y}_1 \leq \hat{y}_2 \leq \dots \leq \hat{y}_n\)。如果需要非递增，只需将 \(y_i\) 取负号或调整顺序。
  • 优化目标：最小化拟合值与观测值之间的平方误差和（或绝对误差和，常用平方误差）。即求解：

\[ \min_{\hat{y}_1 \leq \hat{y}_2 \leq \dots \leq \hat{y}_n} \sum_{i=1}^n (y_i - \hat{y}_i)^2 \]

这是一个带线性不等式约束的凸优化问题，可通过Pool Adjacent Violators Algorithm (PAVA) 高效求解。

2. 数学原理与直观理解

Isotonic Regression 的拟合结果可以理解为在单调约束下对数据的“平滑”。其关键性质包括：

  • 保序性：拟合值 \(\hat{y}_i\) 的顺序与 \(x_i\) 的顺序一致。
  • 分段常数性：最优解 \(\hat{y}_i\) 是分段常数函数，即在某些连续区间上，拟合值相同。这是因为如果相邻点的拟合值不同，则它们之间的约束是严格的（\(\hat{y}_i < \hat{y}_{i+1}\)），但若允许相等，则多个点可合并为一个“块”（block）。
  • 块平均性质：每个“块”内的拟合值等于该块内所有观测值 \(y_i\) 的平均值。这是最小化平方误差的自然结果。

为什么是分段常数？
假设在某个区间上，单调约束不是紧的（即可以取不同值），但为了最小化平方误差，我们会尽量让拟合值接近观测值。如果相邻观测值的顺序违反了单调性（即 \(y_i > y_{i+1}\)），则约束会强制它们相等或调整，从而形成常数块。

3. 求解算法：PAVA（Pool Adjacent Violators Algorithm）

PAVA 是一种贪心算法，通过不断合并违反单调性的相邻块，并计算块内平均值，直到整个序列满足单调性。以下是详细步骤。

步骤1：初始化

将每个数据点视为一个独立的块（block）。每个块 \(B_k\) 包含：

  • 起始索引 \(s_k\) 和结束索引 \(e_k\)（初始时每个块只包含一个点）。
  • 块内观测值的平均值 \(\mu_k = y_{s_k}\)。
  • 块的大小 \(w_k = 1\)。

初始时，我们有 \(n\) 个块：\(B_1, B_2, \dots, B_n\)，满足 \(s_k = e_k = k\)，\(\mu_k = y_k\)。

步骤2：检查并合并违反单调性的相邻块

从左到右扫描块序列，检查相邻块 \(B_k\) 和 \(B_{k+1}\) 是否满足单调性约束：\(\mu_k \leq \mu_{k+1}\)。
如果违反（即 \(\mu_k > \mu_{k+1}\)），则合并这两个块为一个新块，并更新新块的信息：

  • 新块的索引范围：\(s_{\text{new}} = s_k\)，\(e_{\text{new}} = e_{k+1}\)。
  • 新块的大小：\(w_{\text{new}} = w_k + w_{k+1}\)。
  • 新块的平均值：\(\mu_{\text{new}} = \frac{w_k \mu_k + w_{k+1} \mu_{k+1}}{w_k + w_{k+1}}\)（加权平均）。

合并后，用新块替换原来的两个块，块总数减1。

步骤3：回溯检查

合并后，新块的平均值可能小于其左侧相邻块的平均值（因为平均值被拉低了）。因此，需要向左回溯，检查新块与左侧块的单调性。如果违反（\(\mu_{\text{left}} > \mu_{\text{new}}\)），则继续合并左侧块与新块，并重复回溯，直到整个序列满足从左到右单调非递减。

回溯的必要性：
例如，假设有三个块，平均值分别为 5、2、4。

  • 第一轮检查块1和块2：5 > 2，合并为块A，平均值 = (5+2)/2 = 3.5。
  • 现在块A与块3比较：3.5 ≤ 4，似乎满足。但注意，原始序列是 5、2、4，合并后为 3.5、4，单调性成立。然而，如果左侧有更高值，可能需要继续合并。

实际上，更严格的PAVA实现会在每次合并后立即回溯。常见做法是：在扫描过程中，一旦合并发生，就从合并后的块开始向前（左）检查，确保与之前所有块满足单调性。

步骤4：重复扫描直到无违反

继续从前往后扫描整个块序列，重复步骤2-3，直到没有违反单调性的相邻块为止。此时，每个块内的所有点共享相同的拟合值（即块平均值）。

步骤5：输出拟合值

对于每个数据点 \(i\)，其拟合值 \(\hat{y}_i\) 等于它所在块的平均值 \(\mu_k\)。

4. 一个具体计算示例

假设观测值序列 \(y = [5, 3, 4, 2, 6]\)（已按 \(x\) 排序）。目标：找到非递减拟合值 \(\hat{y}\)。

初始化：5个块，平均值分别为 5, 3, 4, 2, 6。

  1. 扫描1：

    • 比较块1(5)和块2(3)：5 > 3，合并。新块A平均值 = (5+3)/2 = 4。块序列变为：[A(4), 4, 2, 6]。
    • 回溯：A左侧无块，继续。
    • 比较A(4)和块3(4)：4 ≤ 4，满足。
    • 比较块3(4)和块4(2)：4 > 2，合并。新块B平均值 = (4+2)/2 = 3。序列变为：[A(4), B(3), 6]。
    • 回溯：比较A(4)和B(3)：4 > 3，合并A和B。新块C平均值 = (24 + 23)/4 = (8+6)/4 = 3.5。序列变为：[C(3.5), 6]。
    • 比较C(3.5)和块5(6)：3.5 ≤ 6，满足。
      扫描结束，当前块：[C(3.5), 6]。

  2. 最终拟合值：

    • 前4个点（索引1-4）在块C中，拟合值均为 3.5。
    • 第5个点在块D（即原块5）中，拟合值为 6。
      所以，\(\hat{y} = [3.5, 3.5, 3.5, 3.5, 6]\)。

验证单调性：3.5 ≤ 3.5 ≤ 3.5 ≤ 3.5 ≤ 6，满足。平方误差和 = (5-3.5)² + (3-3.5)² + (4-3.5)² + (2-3.5)² + (6-6)² = 2.25 + 0.25 + 0.25 + 2.25 + 0 = 5.0。

5. 扩展与讨论

  • 非递增情况：若需要拟合值非递增，可将 \(y_i\) 取负号，应用PAVA后再取负回来，或直接修改比较方向。
  • 权重：PAVA可扩展为加权版本，每个点有权重 \(w_i\)，最小化加权平方和 \(\sum w_i (y_i - \hat{y}_i)^2\)，合并时计算加权平均。
  • 与线性回归关系：若数据本身单调，Isotonic Regression 可能接近线性拟合；若强烈非单调，则拟合为阶梯函数。
  • 应用场景：校准分类器概率输出（保证概率随得分单调）、剂量反应分析、任何需要单调拟合的领域。

6. 总结

Isotonic Regression 通过单调约束下的最小二乘拟合，提供了一种灵活的非参数建模工具。PAVA算法通过合并违反单调性的相邻块，并以块内平均值作为拟合值，高效求得最优分段常数解。其核心在于将全局优化问题分解为局部块的不断调整，直至整个序列满足单调性。

Isotonic Regression（保序回归）的数学原理与拟合过程 Isotonic Regression（保序回归）是一种非参数回归方法，其核心思想是在保持预测值单调性（非递减或非递增）约束的前提下，最小化预测值与观测值之间的误差。它不假设数据服从特定的函数形式（如线性），但要求拟合的函数是单调的。下面我将详细讲解其数学原理、问题描述、求解算法（PAVA）以及计算过程。 1. 问题描述 假设我们有 \( n \) 个观测数据点 \((x_ 1, y_ 1), (x_ 2, y_ 2), \dots, (x_ n, y_ n)\)，其中 \(x_ i\) 是自变量（通常有序，如时间、剂量），\(y_ i\) 是因变量（如响应值）。我们假设 \(x_ i\) 已经按升序排列（即 \(x_ 1 \leq x_ 2 \leq \dots \leq x_ n\)）。如果 \(x_ i\) 相等，则需特殊处理。 Isotonic Regression 的目标是找到一组拟合值 \(\hat{y}_ 1, \hat{y}_ 2, \dots, \hat{y}_ n\)，满足以下约束和优化目标： 单调性约束 ：拟合值必须是非递减的，即 \(\hat{y}_ 1 \leq \hat{y}_ 2 \leq \dots \leq \hat{y}_ n\)。如果需要非递增，只需将 \(y_ i\) 取负号或调整顺序。 优化目标 ：最小化拟合值与观测值之间的平方误差和（或绝对误差和，常用平方误差）。即求解： \[ \min_ {\hat{y}_ 1 \leq \hat{y}_ 2 \leq \dots \leq \hat{y} n} \sum {i=1}^n (y_ i - \hat{y}_ i)^2 \] 这是一个带线性不等式约束的凸优化问题，可通过 Pool Adjacent Violators Algorithm (PAVA) 高效求解。 2. 数学原理与直观理解 Isotonic Regression 的拟合结果可以理解为在单调约束下对数据的“平滑”。其关键性质包括： 保序性 ：拟合值 \(\hat{y}_ i\) 的顺序与 \(x_ i\) 的顺序一致。 分段常数性 ：最优解 \(\hat{y}_ i\) 是分段常数函数，即在某些连续区间上，拟合值相同。这是因为如果相邻点的拟合值不同，则它们之间的约束是严格的（\(\hat{y} i < \hat{y} {i+1}\)），但若允许相等，则多个点可合并为一个“块”（block）。 块平均性质 ：每个“块”内的拟合值等于该块内所有观测值 \(y_ i\) 的平均值。这是最小化平方误差的自然结果。 为什么是分段常数？ 假设在某个区间上，单调约束不是紧的（即可以取不同值），但为了最小化平方误差，我们会尽量让拟合值接近观测值。如果相邻观测值的顺序违反了单调性（即 \(y_ i > y_ {i+1}\)），则约束会强制它们相等或调整，从而形成常数块。 3. 求解算法：PAVA（Pool Adjacent Violators Algorithm） PAVA 是一种贪心算法，通过不断合并违反单调性的相邻块，并计算块内平均值，直到整个序列满足单调性。以下是详细步骤。 步骤1：初始化 将每个数据点视为一个独立的块（block）。每个块 \(B_ k\) 包含： 起始索引 \(s_ k\) 和结束索引 \(e_ k\)（初始时每个块只包含一个点）。 块内观测值的平均值 \(\mu_ k = y_ {s_ k}\)。 块的大小 \(w_ k = 1\)。 初始时，我们有 \(n\) 个块：\(B_ 1, B_ 2, \dots, B_ n\)，满足 \(s_ k = e_ k = k\)，\(\mu_ k = y_ k\)。 步骤2：检查并合并违反单调性的相邻块 从左到右扫描块序列，检查相邻块 \(B_ k\) 和 \(B_ {k+1}\) 是否满足单调性约束：\(\mu_ k \leq \mu_ {k+1}\)。 如果违反（即 \(\mu_ k > \mu_ {k+1}\)），则合并这两个块为一个新块，并更新新块的信息： 新块的索引范围：\(s_ {\text{new}} = s_ k\)，\(e_ {\text{new}} = e_ {k+1}\)。 新块的大小：\(w_ {\text{new}} = w_ k + w_ {k+1}\)。 新块的平均值：\(\mu_ {\text{new}} = \frac{w_ k \mu_ k + w_ {k+1} \mu_ {k+1}}{w_ k + w_ {k+1}}\)（加权平均）。 合并后，用新块替换原来的两个块，块总数减1。 步骤3：回溯检查 合并后，新块的平均值可能小于其左侧相邻块的平均值（因为平均值被拉低了）。因此，需要向左回溯，检查新块与左侧块的单调性。如果违反（\(\mu_ {\text{left}} > \mu_ {\text{new}}\)），则继续合并左侧块与新块，并重复回溯，直到整个序列满足从左到右单调非递减。 回溯的必要性 ： 例如，假设有三个块，平均值分别为 5、2、4。 第一轮检查块1和块2：5 > 2，合并为块A，平均值 = (5+2)/2 = 3.5。 现在块A与块3比较：3.5 ≤ 4，似乎满足。但注意，原始序列是 5、2、4，合并后为 3.5、4，单调性成立。然而，如果左侧有更高值，可能需要继续合并。 实际上，更严格的PAVA实现会在每次合并后立即回溯。常见做法是：在扫描过程中，一旦合并发生，就从合并后的块开始向前（左）检查，确保与之前所有块满足单调性。 步骤4：重复扫描直到无违反 继续从前往后扫描整个块序列，重复步骤2-3，直到没有违反单调性的相邻块为止。此时，每个块内的所有点共享相同的拟合值（即块平均值）。 步骤5：输出拟合值 对于每个数据点 \(i\)，其拟合值 \(\hat{y}_ i\) 等于它所在块的平均值 \(\mu_ k\)。 4. 一个具体计算示例 假设观测值序列 \(y = [ 5, 3, 4, 2, 6 ]\)（已按 \(x\) 排序）。目标：找到非递减拟合值 \(\hat{y}\)。 初始化 ：5个块，平均值分别为 5, 3, 4, 2, 6。 扫描1 ： 比较块1(5)和块2(3)：5 > 3，合并。新块A平均值 = (5+3)/2 = 4。块序列变为：[ A(4), 4, 2, 6 ]。 回溯：A左侧无块，继续。 比较A(4)和块3(4)：4 ≤ 4，满足。 比较块3(4)和块4(2)：4 > 2，合并。新块B平均值 = (4+2)/2 = 3。序列变为：[ A(4), B(3), 6 ]。 回溯：比较A(4)和B(3)：4 > 3，合并A和B。新块C平均值 = (2 4 + 2 3)/4 = (8+6)/4 = 3.5。序列变为：[ C(3.5), 6 ]。 比较C(3.5)和块5(6)：3.5 ≤ 6，满足。 扫描结束，当前块：[ C(3.5), 6 ]。 最终拟合值 ： 前4个点（索引1-4）在块C中，拟合值均为 3.5。 第5个点在块D（即原块5）中，拟合值为 6。 所以，\(\hat{y} = [ 3.5, 3.5, 3.5, 3.5, 6 ]\)。 验证单调性：3.5 ≤ 3.5 ≤ 3.5 ≤ 3.5 ≤ 6，满足。平方误差和 = (5-3.5)² + (3-3.5)² + (4-3.5)² + (2-3.5)² + (6-6)² = 2.25 + 0.25 + 0.25 + 2.25 + 0 = 5.0。 5. 扩展与讨论 非递增情况 ：若需要拟合值非递增，可将 \(y_ i\) 取负号，应用PAVA后再取负回来，或直接修改比较方向。 权重 ：PAVA可扩展为加权版本，每个点有权重 \(w_ i\)，最小化加权平方和 \(\sum w_ i (y_ i - \hat{y}_ i)^2\)，合并时计算加权平均。 与线性回归关系 ：若数据本身单调，Isotonic Regression 可能接近线性拟合；若强烈非单调，则拟合为阶梯函数。 应用场景 ：校准分类器概率输出（保证概率随得分单调）、剂量反应分析、任何需要单调拟合的领域。 6. 总结 Isotonic Regression 通过单调约束下的最小二乘拟合，提供了一种灵活的非参数建模工具。PAVA算法通过合并违反单调性的相邻块，并以块内平均值作为拟合值，高效求得最优分段常数解。其核心在于将全局优化问题分解为局部块的不断调整，直至整个序列满足单调性。