Python SciPy optimize.isotonic_regression用法及代码示例

本文简要介绍 python 语言中 scipy.optimize.isotonic_regression 的用法。

用法:

scipy.optimize.isotonic_regression(y, *, weights=None, increasing=True)#

非参数等渗回归。

一个(非严格)单调递增数组x与相同的长度y由池相邻违规者算法(PAVA)计算，参见[1]。有关更多详细信息，请参阅注释部分。

参数 ：：

y： (N,) 数组

响应变量。

weights： (N,) 数组 或 None

箱重。

increasing： bool

如果为 True，则拟合单调递增，即等渗回归。如果为 False，则拟合单调递减(即反调回归)。默认为 True。

返回 ：：

res： OptimizeResult

优化结果表示为 OptimizeResult 对象。重要的属性是：

• x ：等渗回归解，即长度与 y 相同的递增(或递减)数组，元素范围从 min(y) 到 max(y)。

• weights：包含每个块(或池)B 的案例权重总和的数组。

• blocks：长度为 B+1 的数组，其中包含每个块(或池)B 的起始位置的索引。j-th 块由 x[blocks[j]:blocks[j+1]] 给出，其所有值都相同。

注意：

给定数据 \(y\) 和案例权重 \(w\) ，等渗回归解决以下优化问题：

\[\operatorname{argmin}_{x_i} \sum_i w_i (y_i - x_i)^2 \quad \text{subject to } x_i \leq x_j \text{ whenever } i \leq j \,.\]

对于每个输入值 \(y_i\) ，它会生成一个值 \(x_i\) ，使得 \(x\) 增加(但不严格)，即 \(x_i \leq x_{i+1}\) 。这是由 PAVA 完成的。该解决方案由池或块组成，即 \(x\) 的相邻元素，例如\(x_i\) 和 \(x_{i+1}\) 都具有相同的值。

最有趣的是，如果平方损失被广泛的 Bregman 函数替换，那么解决方案保持不变，Bregman 函数是严格一致的均值评分函数的独特类别，请参阅 [2] 和其中的参考文献。

根据[1]，PAVA 的实现版本的计算复杂度为 O(N)，输入大小为 N。

参考：

[1] (1,2)

布辛，F. M. T. A. (2022)。单调回归：简单快速的 O(n) PAVA 实现。统计软件杂志，代码片段，102(1), 1-25。 DOI:10.18637/jss.v102.c01

[2]

Jordan, A.I.、Mühlemann, A. 和 Ziegel, J.F. 说明可识别泛函的等渗回归问题的最优解。安研究所统计数学 74, 489-514 (2022)。 DOI:10.1007/s10463-021-00808-0

例子：

此示例表明isotonic_regression 确实解决了约束优化问题。

>>> import numpy as np
>>> from scipy.optimize import isotonic_regression, minimize
>>> y = [1.5, 1.0, 4.0, 6.0, 5.7, 5.0, 7.8, 9.0, 7.5, 9.5, 9.0]
>>> def objective(yhat, y):
...     return np.sum((yhat - y)**2)
>>> def constraint(yhat, y):
...     # This is for a monotonically increasing regression.
...     return np.diff(yhat)
>>> result = minimize(objective, x0=y, args=(y,),
...                   constraints=[{'type': 'ineq',
...                                 'fun': lambda x: constraint(x, y)}])
>>> result.x
array([1.25      , 1.25      , 4.        , 5.56666667, 5.56666667,
       5.56666667, 7.8       , 8.25      , 8.25      , 9.25      ,
       9.25      ])
>>> result = isotonic_regression(y)
>>> result.x
array([1.25      , 1.25      , 4.        , 5.56666667, 5.56666667,
       5.56666667, 7.8       , 8.25      , 8.25      , 9.25      ,
       9.25      ])

与调用 minimize 相比，isotonic_regression 的一大优点是它更加用户友好，即不需要定义目标函数和约束函数，并且速度快了几个数量级。在商用硬件上(2023 年)，对于长度为 1000 的正态分布输入 y，最小化器大约需要 4 秒，而 isotonic_regression 大约需要 200 微秒。