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 微秒。