Isotonic regression--保序回归

本文介绍了 Isotonic regression 的定义、调用方法以及在模型输出校准领域的运用。
和这篇文章有关的论文可以参考：Individual Calibration with Randomized Forecasting

1 定义

The class IsotonicRegression fits a non-decreasing real function to 1-dimensional data. It solves the following problem:
${\sum }_i{w}_i(y_i-{\hat{y}}_i{)}^2$
subject to ${\hat{y}}_i\le {\hat{y}}_j$ whenever $X_i\le X_j$ ,
where the weights $w$ are strictly positive, and both $X$ and $y$ are arbitrary real quantities(任意实数).
The increasing parameter changes the constraint to ${\hat{y}}_i\ge {\hat{y}}_j$ whenever $X_i\le X_j$. Setting it to ‘auto’ will automatically choose the constraint based on Spearman’s rank correlation coefficient.
可以这么说，IsotonicRegression是为了生成一个单调递增（monotonically increasing）的函数，使对应样本特征 $x$ 的预测值 $\hat{y}$ 和真实的样本标签 $y$ 的均方误差达到最小。操作的时候是给定 { y i , y ^ i } \{y_i,\hat{y}_i\} {yi​,y^​i​} 这一串序列的，并且这串序列是从小到大排的
IsotonicRegression produces a series of predictions ${\hat{y}}_i$
for the training data which are the closest to the targets $y$ in terms of mean squared error. These predictions are interpolated for predicting to unseen data. The predictions of IsotonicRegression thus form a function that is piecewise linear（分段线性）:

举例：

就比如在下图中，当测试集出现在 ${\hat{y}}_{i-1}$ 和 ${\hat{y}}_i$ 之间的时候，模型校准前的输出可能不是那么平稳，输出会出现在绿点，因此，我们对 ${\hat{y}}_{i-1}$ 到 ${\hat{y}}_i$ 区间的输出做一个平滑处理，例如根据 $y_{i-1}$ 和 $y_i$ 的平均值来，这样可以使得预测点变为红点，保证了当模型输出随着样本真实标签递增，有效地克服了原来输出的方差影响。

2 调用方法

代码如下，效果图如上：

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.collections import LineCollection

from sklearn.linear_model import LinearRegression
from sklearn.isotonic import IsotonicRegression
from sklearn.utils import check_random_state
# generate data points
n = 100
x = np.arange(n)
rs = check_random_state(0)
y = rs.randint(-50, 50, size=(n,)) + 50. * np.log1p(np.arange(n))

# Fit IsotonicRegression and LinearRegression models:
ir = IsotonicRegression(out_of_bounds="clip")
y_ = ir.fit_transform(x, y)

lr = LinearRegression()
lr.fit(x[:, np.newaxis], y)  # x needs to be 2d for LinearRegression

# plot results
segments = [[[i, y[i]], [i, y_[i]]] for i in range(n)]
lc = LineCollection(segments, zorder=0)
lc.set_array(np.ones(len(y)))
lc.set_linewidths(np.full(n, 0.5))

fig, (ax0, ax1) = plt.subplots(ncols=2, figsize=(12, 6))

ax0.plot(x, y, 'C0.', markersize=12)
ax0.plot(x, y_, 'C1.-', markersize=12)
ax0.plot(x, lr.predict(x[:, np.newaxis]), 'C2-')
ax0.add_collection(lc)
ax0.legend(('Training data', 'Isotonic fit', 'Linear fit'), loc='lower right')
ax0.set_title('Isotonic regression fit on noisy data (n=%d)' % n)

x_test = np.linspace(-10, 110, 1000)
ax1.plot(x_test, ir.predict(x_test), 'C1-')
ax1.plot(ir.X_thresholds_, ir.y_thresholds_, 'C1.', markersize=12)
ax1.set_title("Prediction function (%d thresholds)" % len(ir.X_thresholds_))

plt.show()

Note that we explicitly passed out_of_bounds="clip" to the constructor of IsotonicRegression to control the way the model extrapolates outside of the range of data observed in the training set. This “clipping” extrapolation can be seen on the plot of the decision function on the right-hand.

2.1参数和方法

3 model calibration

保序回归用于校准伪代码如下：

保序在这其中起到的作用就是保证有序，因为我们希望校准后的概率值应该随着函数值的增大而增大。特别要注意的是，和Platt Calibration类似，为了保证不引入偏差，用作校准的数据集应该和训练模型的数据集不同。

4 应用举例

以某种药物的使用量为例子：

假设药物使用量为数组X=0,1,2,3,4….99，病人对药物的反应量为Y=y1,y2,y3……y99 ，而由于个体的原因，Y不是一个单调函数(即：存在波动)，**如果我们按照药物反应排序，对应的X就会成为乱序，失去了研究的意义。**而我们的研究的目的是为了观察随着药物使用量的递增，病人的平均反应状况。在这种情况下，使用保序回归，即不改变X的排列顺序，又求的Y的平均值状况。

从图中可以看出，最长的蓝线x的取值约是30到60，在这个区间内，Y的平均值一样，那么从经济及病人抗药性等因素考虑，使用药量为30个单位是最理想的。
当然，这个例子侧重在回归分析某一变量对输出结果的影响，而不是模型输出校准~
当前IT行业虚拟化比较流行，使用这种方式，找到合适的判断参数，就可以使用此算法使资源得到最大程度的合理利用。

参考：
https://scikit-learn.org/stable/modules/generated/sklearn.isotonic.IsotonicRegression.html#sklearn.isotonic.IsotonicRegression
https://zhuanlan.zhihu.com/p/101766505
https://zhuanlan.zhihu.com/p/88623159