[机器学习] XGBoost 自定义损失函数-FocalLoss

一.XGBoost损失函数

损失函数：损失函数描述了预测值和真实标签的差异，通过对损失函数的优化来获得对学习任务的一个近似求解方法
boosting类算法的损失函数的作用： Boosting的框架, 无论是GBDT还是Adaboost, 其在每一轮迭代中, 根本没有理会损失函数具体是什么, 仅仅用到了损失函数的一阶导数通过随机梯度下降来参数更新
二阶导数：GBDT在模型训练时只使用了代价函数的一阶导数信息，XGBoost对代价函数进行二阶泰勒展开，可以同时使用一阶和二阶导数
牛顿法梯度更新：XGBoost是用了牛顿法进行的梯度更新。通过对损失进行分解得到一阶导数和二阶导数并通过牛顿法来迭代更新梯度。

XGBOOST是GBDT的一种高效实现方式，XGBOOST在每一次迭代时，对损失函数进行了二阶泰勒展开，在下一次迭代时同时考虑了一阶导和二阶导。这一点与GBDT仅使用一阶导有所区别。

XGBOOST是一个非常灵活的模型，允许使用者根据实际使用场景调整损失函数，对于常见的二分类问题一般使用的binary：logistic损失函数，其形式为：

这个损失函数对于正类错分和负类错分给予的惩罚时相同的，但是对于不平衡数据集，或者某些特殊情况（两类错分代价不一样）的时候这样的损失函数就不再合理了。

XGBOOST很灵活，允许使用者自定义损失函数的形式，只需要输入损失函数的一阶导和二阶导（grad，hess）即可。

基于XGBoost的损失函数的分解求导，可以知道XGBoost的除正则项以外的核心影响因子是损失函数的1阶导和2阶导，所以对于任意的学习任务的损失函数，可以对其求一阶导数和二阶导数带入到XGBoost的自定义损失函数范式里面进行处理。

def custom_obj(pred, dtrain):# pred 和 dtrain 的顺序不能弄反# STEP1 获得labellabel = dtrain.get_label()# STEP2 如果是二分类任务，需要让预测值通过sigmoid函数获得0～1之间的预测值# 如果是回归任务则下述任务不需要通过sigmoid# 分类任务sigmoid化def sigmoid(x):return 1/(1+np.exp(-x))sigmoid_pred = sigmoid(原始预测值)#回归任务pred = 原始预测值# STEP3 一阶导和二阶导grad = 一阶导hess = 二阶导return grad, hess

二.常见的损失函数

分类问题

非平衡分类学习任务，例如首笔首期30+的风险建模任务，首期30+的逾期率比例相对ever30+的逾期率为1/3左右，通过修正占比少的正样本权重来对影响正样本对损失函数的贡献度，可以进一步提升模型的效果

def weighted_binary_cross_entropy(pred, dtrain,imbalance_alpha=10):# retrieve data from dtrain matrixlabel = dtrain.get_label()# compute the prediction with sigmoidsigmoid_pred = 1.0 / (1.0 + np.exp(-pred))# gradientgrad = -(imbalance_alpha ** label) * (label - sigmoid_pred)hess = (imbalance_alpha ** label) * sigmoid_pred * (1.0 - sigmoid_pred)return grad, hess

Focal Loss

Focal Loss for Dense Object Detection 是ICCV2017的Best student paper,文章思路很简单但非常具有开拓性意义，效果也非常令人称赞。

大家还可以看知乎的讨论：如何评价 Kaiming 的 Focal Loss for Dense Object Detection？

Focal Loss的引入主要是为了解决难易样本数量不平衡（注意，有区别于正负样本数量不平衡）的问题，实际可以使用的范围非常广泛，为了方便解释，拿目标检测的应用场景来说明

Focal Loss的主要思想就是改变损失函数.Focal loss是在交叉熵损失函数基础上进行的修改

单阶段的目标检测器通常会产生高达100k的候选目标，只有极少数是正样本，正负样本数量非常不平衡。我们在计算分类的时候常用的损失——交叉熵的公式如下：

$y{}'$ 是经过激活函数的输出，所以在0-1之间。可见普通的交叉熵对于正样本而言，输出概率越大损失越小。对于负样本而言，输出概率越小则损失越小。此时的损失函数在大量简单样本的迭代过程中比较缓慢且可能无法优化至最优。

为了解决正负样本不平衡的问题，我们通常会在交叉熵损失的前面加上一个参数平衡因子alpha，用来平衡正负样本本身的比例不均. 文中alpha取0.25，即正样本要比负样本占比小，这是因为负例易分。

只添加alpha虽然可以平衡正负样本的重要性，, 对难易样本的不平衡没有任何帮助

根据正、负、难、易，样本一共可以分为以下四类：

但这并不能解决全部问题, 实际上目标检测中大量的候选目标都是像下图一样的易分样本。

这些样本的损失很低，但是由于数量极不平衡，易分样本的数量相对来讲太多，最终主导了总的损失。而本文的作者认为，易分样本（即置信度高的样本）对模型的提升效果非常小，模型应该主要关注与那些难分样本

那么Focal loss是怎么改进的呢？

一个简单的思想：把高置信度(p)样本的损失再降低一些不就好了吗！

首先在原有的基础上加了一个因子，其中gamma>0使得减少易分类样本的损失。使得更关注于困难的、错分的样本。

举个例，gamma为2，对于正类样本而言，预测结果为0.95肯定是简单样本，所以（1-0.95)的gamma次方就会很小，这时损失函数值就变得更小。而预测概率为0.3的样本其损失相对很大。

对于负类样本而言同样，预测0.1的结果应当远比预测0.7的样本损失值要小得多。对于预测概率为0.5时，损失只减少了0.25倍，所以更加关注于这种难以区分的样本。这样减少了简单样本的影响，大量预测概率很小的样本叠加起来后的效应才可能比较有效。

Focal Loss的最终形式结合了上面的公式解决了难易样本的不平衡，解决了正负样本的不平衡，结合使用，同时解决正负难易2个问题！

#focal 损失
# LGB+Focal Loss 其中alpha：为不能让容易分类类别的损失函数太小, 默认值0.25；
#gamma：更加关注困难样本 即关注y=1的样本 默认值2
from scipy.misc import derivative
def focal_loss_lgb_sk( y_pred,dtrain, alpha=0.25, gamma=2):label = dtrain.get_label()a,g = alpha, gammadef fl(x,t):p = 1/(1+np.exp(-x))return -( a*t + (1-a)*(1-t) ) * (( 1 - ( t*p + (1-t)*(1-p)) )**g) * ( t*np.log(p)+(1-t)*np.log(1-p) )partial_fl = lambda x: fl(x, label)grad = derivative(partial_fl, y_pred, n=1, dx=1e-6)hess = derivative(partial_fl, y_pred, n=2, dx=1e-6)return grad, hess

最终的Focal Loss形式如下：

gamma调节简单样本权重降低的速率，当gamma为0时即为交叉熵损失函数，当gamma增加时，调整因子的影响也在增加。实验发现gamma为2是最优。

这样以来，训练过程关注对象的排序为正难>负难>正易>负易。

1. alpha $\in$ (0,1) 二分类场景下，类似于正例的sample_weight概念，可以按照样本占比，适度加权, 假如有5条正例、95条负例，则建议 alpha=0.95，取0.5 相关于关掉此功能

2. gamma $\in$ [0， + $\infty$ ) 对于难度的区分程度, 取 gamma = 0 相当于关掉该功能, gamma越大则越重视难度，即专注于比较困难的样本。建议在 (0.5, 10.0) 范围尝试

自定义损失函数 XGBoost+Focal Loss

xgboost如果想应用focal loss可以使用自定义损失函数来实现。

首先介绍一个比较好用的求导的python包 sympy。用这个包来求导，可以比较轻松求出一阶导和二阶导。

下面的alpha 与上一节中平衡因子alpha 不太一样，相当于上一节中只在y=1处加入了alpha，然后y=0不加入。 # alpha指的是用来处理非平衡的参数, label为1的样本的损失是label为0的多少倍（1-N）

from sympy import diff, symbols, log# y指的是true label
# p指的是预测出来的概率值（是经过sigmoid转换之后的）
# gamma指的是local foss的参数
# alpha指的是用来处理非平衡的参数
# 同时考虑了focal loss参数和处理非平衡参数的损失函数
y, p, gamma, alpha = symbols('y p gamma alpha')
loss = alpha * (-y * log(p) * (1 - p) ** gamma) - (1 - alpha) * (1 - y) * log(1 - p) * p ** gamma
grad = diff(loss, p) * p * (1 - p)  # 求一阶导
print('grad:')
print(grad)
hess = diff(grad, p) * p * (1 - p)  # 求二阶导
print('hess:')
print(hess)'''
grad:
p*(1 - p)*(alpha*gamma*y*(1 - p)**gamma*log(p)/(1 - p) - alpha*y*(1 - p)**gamma/p - gamma*p**gamma*(1 - alpha)*(1 - y)*log(1 - p)/p + p**gamma*(1 - alpha)*(1 - y)/(1 - p))
hess:
p*(1 - p)*(p*(1 - p)*(-alpha*gamma**2*y*(1 - p)**gamma*log(p)/(1 - p)**2 + alpha*gamma*y*(1 - p)**gamma*log(p)/(1 - p)**2 + 2*alpha*gamma*y*(1 - p)**gamma/(p*(1 - p)) + alpha*y*(1 - p)**gamma/p**2 - gamma**2*p**gamma*(1 - alpha)*(1 - y)*log(1 - p)/p**2 + 2*gamma*p**gamma*(1 - alpha)*(1 - y)/(p*(1 - p)) + gamma*p**gamma*(1 - alpha)*(1 - y)*log(1 - p)/p**2 + p**gamma*(1 - alpha)*(1 - y)/(1 - p)**2) - p*(alpha*gamma*y*(1 - p)**gamma*log(p)/(1 - p) - alpha*y*(1 - p)**gamma/p - gamma*p**gamma*(1 - alpha)*(1 - y)*log(1 - p)/p + p**gamma*(1 - alpha)*(1 - y)/(1 - p)) + (1 - p)*(alpha*gamma*y*(1 - p)**gamma*log(p)/(1 - p) - alpha*y*(1 - p)**gamma/p - gamma*p**gamma*(1 - alpha)*(1 - y)*log(1 - p)/p + p**gamma*(1 - alpha)*(1 - y)/(1 - p)))
'''

注意: xgboost版本为 0.90

在xgboost里边设置自定义损失函数，将上边生成的一阶导和二阶导贴过来。

当设置alpha=0.5, gamma=0时，就是原始的损失函数, 相当于logloss*0.5，对比输出的结果可以判断求导是否正确。

import numpy as np
import xgboost as xgbgamma = 0
alpha = 0.5# y指的是true label
# p指的是预测出来的概率值（是经过sigmoid转换之后  p = 1.0 / (1.0 + np.exp(-p))
# gamma指的是local foss的参数
# alpha指的是用来处理非平衡的参数, label为1的样本的损失是label为0的多少倍
def logistic_obj(p, dtrain):y = dtrain.get_label()p = 1.0 / (1.0 + np.exp(-p))grad = p * (1 - p) * (alpha * gamma * y * (1 - p) ** gamma * np.log(p) / (1 - p) - alpha * y * (1 - p) ** gamma / p - gamma * p ** gamma * (1 - alpha) * (1 - y) * np.log(1 - p) / p + p ** gamma * (1 - alpha) * (1 - y) / (1 - p))hess = p * (1 - p) * (p * (1 - p) * (-alpha * gamma ** 2 * y * (1 - p) ** gamma * np.log(p) / (1 - p) ** 2 + alpha * gamma * y * (1 - p) ** gamma * np.log(p) / (1 - p) ** 2 + 2 * alpha * gamma * y * (1 - p) ** gamma / (p * (1 - p)) + alpha * y * (1 - p) ** gamma / p ** 2 - gamma ** 2 * p ** gamma * (1 - alpha) * (1 - y) * np.log(1 - p) / p ** 2 + 2 * gamma * p ** gamma * (1 - alpha) * (1 - y) / (p * (1 - p)) + gamma * p ** gamma * (1 - alpha) * (1 - y) * np.log(1 - p) / p ** 2 + p ** gamma * (1 - alpha) * (1 - y) / (1 - p) ** 2) - p * (alpha * gamma * y * (1 - p) ** gamma * np.log(p) / (1 - p) - alpha * y * (1 - p) ** gamma / p - gamma * p ** gamma * (1 - alpha) * (1 - y) * np.log(1 - p) / p + p ** gamma * (1 - alpha) * (1 - y) / (1 - p)) + (1 - p) * (alpha * gamma * y * (1 - p) ** gamma * np.log(p) / (1 - p) - alpha * y * (1 - p) ** gamma / p - gamma * p ** gamma * (1 - alpha) * (1 - y) * np.log(1 - p) / p + p ** gamma * (1 - alpha) * (1 - y) / (1 - p)))return grad, hessdtrain = xgb.DMatrix('../demo_data/agaricus.txt.train')
dtest = xgb.DMatrix('../demo_data/agaricus.txt.test')print("---------no focal loss-----------")watchlist = [(dtrain, 'train'), (dtest, 'eval')]params = {'max_depth': 2, 'eta': 0.01, 'silent': 1, 'eval_metric': 'auc', "objective": "binary:logistic"}
xgb.train(params=params, dtrain=dtrain, num_boost_round=3, early_stopping_rounds=50,evals=[(dtrain, 'train'), (dtest, 'test')], verbose_eval=1)print("---------focal loss-----------")
params = {'max_depth': 2, 'eta': 0.01, 'silent': 1, 'eval_metric': 'auc', "objective": "binary:logitraw"}
xgb.train(params=params, dtrain=dtrain, num_boost_round=3, early_stopping_rounds=50,evals=[(dtrain, 'train'), (dtest, 'test')], verbose_eval=1, obj=logistic_obj)

回归问题

美团的ETA预估

深度学习在美团配送ETA预估中的探索与实践 - 美团技术团队

通过对损失函数进行修正大幅度提高了业务预测效果（使得结果倾向于提前到达）

1.带惩罚的最小二乘损失函数

def custom_normal_train( y_pred,dtrain):label = dtrain.get_label()residual = (label - y_pred).astype("float")grad = np.where(residual<0, -2*(residual)/(label+1), -10*2*(residual)/(label+1))#对预估里程低于实际里程的情况加大惩罚hess = np.where(residual<0, 2/(label+1), 10*2/(label+1))#对预估里程低于实际里程的情况加大惩罚return grad, hess

https://www.kaggle.com/c/allstate-claims-severity/discussion/24520#140255

kaggle 的保险定价竞赛中评价函数为MAE（存在不可导区域），通过使用下述的损失函数来获得MAE的近似损失函数来求解该问题

1.fair loss

def fair_obj( preds,dtrain):labels = dtrain.get_label() con = 2residual = preds-labelsgrad = con*residual / (abs(residual)+con)hess = con**2 / (abs(residual)+con)**2return grad,hess

https://www.kaggle.com/c/allstate-claims-severity/discussion/24520#140255 提供了非常高质量的近似MAE损失函数优化的讨论

2.mae的近似损失函数- Pseudo Huber loss

def huber_approx_obj( preds,dtrain):  d = preds -dtrain.get_label()  h = 1 #h is delta in the formula  scale = 1 + (d / h) ** 2  scale_sqrt = np.sqrt(scale)  grad = d / scale_sqrt  hess = 1 / scale / scale_sqrt  return grad, hess

3.mae的近似损失函数 ln(cosh) loss

def logcoshobj( preds,dtrain): label =  dtrain.get_label()d = preds - labelgrad = np.tanh(d)/label  hess = (1.0 - grad*grad)/label  return grad, hess