从零实现深度学习框架——逻辑回归中的数值稳定
引言
本着“凡我不能创造的,我就不能理解”的思想,本系列文章会基于纯Python以及NumPy从零创建自己的深度学习框架,该框架类似PyTorch能实现自动求导。
要深入理解深度学习,从零开始创建的经验非常重要,从自己可以理解的角度出发,尽量不适用外部完备的框架前提下,实现我们想要的模型。本系列文章的宗旨就是通过这样的过程,让大家切实掌握深度学习底层实现,而不是仅做一个调包侠。
本系列文章首发于微信公众号:JavaNLP
本文我们来探讨一下逻辑回归中的数值稳定问题,所谓数值稳定即不会出现数值溢出(上溢或下溢)的问题。常见的表现是除000,或返回nan
。
上溢和下溢
我们在计算机中表示实数时,几乎总会引入一些近似误差。当许多操作复合时,即使是理论上可行的算法,在实践时也可能回导致算法失效。
一种极具毁灭性的舍入误差是下溢(underflow)。当接近零的数被四舍五入为零时发生下溢。我们通常要避免被零除或避免取零的对数。
令一种数值错误形式是上溢(overflow)。当大量级的数被近似为无穷大或负无穷大时发生上溢。
Sigmoid函数稳定性问题
我们知道Sigmoid函数公式为:
σ(x)=11+exp(−x)(1)\sigma(x) = \frac{1}{1 + \exp(-x)} \tag{1} σ(x)=1+exp(−x)1(1)
对应的图像如下:
其中包含一个exp(−x)\exp(-x)exp(−x),我们看一下exe^xex的图像:
从上图可以看出,如果xxx很大,exe^xex会非常大,而很小就没事(不会溢出),变成无限接近000。
当Sigmoid函数中的xxx负的特别多,那么exp(−x)\exp(-x)exp(−x)就会变成∞\infty∞,就出现了上溢;
那么如何解决这个问题呢?σ(x)\sigma(x)σ(x)可以表示成两种形式:
σ(x)=11+exp(−x)=exp(x)1+exp(x)(2)\sigma(x) = \frac{1}{1 + \exp(-x)} = \frac{\exp(x)}{1 + \exp(x)} \tag{2} σ(x)=1+exp(−x)1=1+exp(x)exp(x)(2)
当x≥0x \geq 0x≥0时,我们根据exe^{x}ex的图像,我们取11+exp(−x)\frac{1}{1 + \exp(-x)}1+exp(−x)1的形式;
当x<0x < 0x<0时,我们取exp(x)1+exp(x)\frac{\exp(x)}{1 + \exp(x)}1+exp(x)exp(x)
# 原来的做法
def sigmoid_naive(x):return 1 / (1 + math.exp(-x))# 优化后的做法
def sigmoid(x):if x < 0:return math.exp(x) / (1 + math.exp(x))else:return 1 / (1 + math.exp(-x))
然后用不同的数值进行测试:
> sigmoid_naive(2000)
1.0
> sigmoid(2000)
1.0
> sigmoid_naive(-2000)
OverflowError: math range error
> sigmoid(-2000)
0.0
如果传入-2000
,普通的实现会出现溢出,而优化后的版本不会。
但是这里的实现包含了if
判断,同时只判断了一个标量而不是向量。
有一种更好的方法是,我们的逻辑回归只计算出logit,然后将logit传入损失函数,这里的logit说的是逻辑回归中线性变换的输出(加权和加上偏置)。
数值稳定的BCE损失
在pytorch的github中,有一段代码 https://github.com/pytorch/pytorch/issues/751,
class StableBCELoss(nn.modules.Module):def __init__(self):super(StableBCELoss, self).__init__()def forward(self, input, target):neg_abs = - input.abs()loss = input.clamp(min=0) - input * target + (1 + neg_abs.exp()).log()return loss.mean()
笔者又在一篇文章[^1]中找到了对代码进行的解释,如下图:
我们先来分析下为什么是数值稳定的。
import numpy as np# 先定义一个数值稳定的sigmoid
def sigmoid(x):if x < 0:return np.exp(x) / (1 + np.exp(x))else:return 1 / (1 + np.exp(-x))# 逻辑回归损失的常规实现
def bce_loss_naive(y, z):return -y * np.log(sigmoid(z)) - (1-y) * np.log(1 - sigmoid(z))# 数值稳定版逻辑回归损失
def bce_loss(y, z):neg_abs = - np.abs(z)return np.clip(z, a_min=0,a_max=None) - y * z + np.log(1 + np.exp(neg_abs))
接下来我们进行测试,假设zzz是一个较大的数,比如200020002000,我们知道σ(2000)\sigma(2000)σ(2000)会输出111,那么L(1,σ(2000))L(1,\sigma(2000))L(1,σ(2000))应该为000。
> z = 2000
> y = 1
> sigmoid(z)
1.0
> bce_loss(y, z) # 数值稳定版本的
0.0
> bce_loss_naive(y, z)
/usr/local/lib/python3.7/dist-packages/ipykernel_launcher.py:2: RuntimeWarning: divide by zero encountered in log
/usr/local/lib/python3.7/dist-packages/ipykernel_launcher.py:2: RuntimeWarning: invalid value encountered in double_scalars
nan
数值稳定版本的还是稳!常规实现就会碰到nan
了。
我们在假设zzz是一个负的较大的数,比如−2000-2000−2000,那么σ(−2000)=0\sigma(-2000)=0σ(−2000)=0,即L(0,σ(−2000))L(0,\sigma(-2000))L(0,σ(−2000))也应该为000。
> z = -2000
> y = 0
> sigmoid(z)
0.0
> bce_loss(y, z)
0.0
> bce_loss_naive(y, z)
/usr/local/lib/python3.7/dist-packages/ipykernel_launcher.py:2: RuntimeWarning: divide by zero encountered in log
/usr/local/lib/python3.7/dist-packages/ipykernel_launcher.py:2: RuntimeWarning: invalid value encountered in double_scalars
nan
嗯,看起来不错。下面我们仔细分析下这段代码的正确性,没问题我们才可以放心地实现。
我们把y^=σ(z)\hat y = \sigma(z)y^=σ(z)带入逻辑回归的损失函数中得:
LCE(y,z)=−y⋅log(σ(z))−(1−y)⋅log(1−σ(z))(3)L_{CE}(y,z) = -y \cdot \log(\sigma(z))-(1-y) \cdot \log(1-\sigma(z)) \tag 3 LCE(y,z)=−y⋅log(σ(z))−(1−y)⋅log(1−σ(z))(3)
就是说如果logit,也就是z≥0z \geq 0z≥0,那么代码实现中的损失就变成了z−zy+log(1+e−z)z - zy + \log(1 + e^{-z})z−zy+log(1+e−z)。我们来验证一下。
这里分为两种情况,分别是真实标签y=0y=0y=0和y=1y=1y=1。
(1) 当z≥0z \geq 0z≥0且y=0y=0y=0时,此时代码计算的是
z+log(1+e−z)=logez+log(1+e−z)=log(1+e−z)⋅ez=log(1+ez)=−log11+ez=−log1+ez−ez1+ez=−log(1−ez1+ez)=−log(1−11+e−z)\begin{aligned} z + \log(1+e^{-z}) &= \log e^z + \log(1 + e^{-z}) \\ &= \log(1 + e^{-z}) \cdot e^z \\ &= \log(1 + e^z) \\ &= - \log \frac{1}{1 + e^z} \\ &= -\log \frac{1 + e^z - e^z}{1 + e^z} \\ &= - \log \left(1 - \frac{e^z}{1+e^z} \right) \\ &= -\log \left(1 - \frac{1}{1 + e^{-z}} \right) \end{aligned} z+log(1+e−z)=logez+log(1+e−z)=log(1+e−z)⋅ez=log(1+ez)=−log1+ez1=−log1+ez1+ez−ez=−log(1−1+ezez)=−log(1−1+e−z1)
公式(3)(3)(3)计算为−log(1−σ(z))- \log(1 - \sigma(z))−log(1−σ(z)),结果是一样的。
(2)当z≥0z \geq 0z≥0且y=1y=1y=1时,代码计算的是
z−z+log(1+e−z)=log(1+e−z)=−log(11+e−z)\begin{aligned} z - z + \log(1 + e^{-z}) &= \log (1+ e^{-z}) \\ &= - \log \left( \frac{1}{1+ e^{-z}}\right) \end{aligned} z−z+log(1+e−z)=log(1+e−z)=−log(1+e−z1)
公式(3)(3)(3)计算的是−log(σ(z))- \log(\sigma(z))−log(σ(z)),结果也是一样的。
如果z<0z < 0z<0,那么代码实现中的损失就变成了−zy+log(ez+1)-zy + \log(e^z + 1)−zy+log(ez+1),我们也来验证一下。
(3)当z<0z < 0z<0且y=1y=1y=1时,代码计算的是
−z+log(ez+1)=loge−z+log(ez+1)=log(ez+1)⋅e−z=log(1+e−z)=−log(11+e−z)\begin{aligned} -z + \log(e^z + 1) &= \log e^{-z} + \log(e^z + 1) \\ &= \log (e^z+1)\cdot e^{-z} \\ &= \log(1 + e^{-z}) \\ &= - \log \left( \frac{1}{1+ e^{-z}}\right) \end{aligned} −z+log(ez+1)=loge−z+log(ez+1)=log(ez+1)⋅e−z=log(1+e−z)=−log(1+e−z1)
公式(3)(3)(3)计算的是−log(σ(z))- \log(\sigma(z))−log(σ(z)),结果也是一样的。
(4)当z<0z<0z<0且y=0y=0y=0时,代码计算的是
log(ez+1)=−log11+ez=−log1+ez−ez1+ez=−log(1−ez1+ez)=−log(1−11+e−z)\begin{aligned} \log(e^z + 1) &= - \log \frac{1}{1 + e^z} \\ &= -\log \frac{1 + e^z - e^z}{1 + e^z} \\ &= - \log \left(1 - \frac{e^z}{1+e^z} \right) \\ &= -\log \left(1 - \frac{1}{1 + e^{-z}} \right) \end{aligned} log(ez+1)=−log1+ez1=−log1+ez1+ez−ez=−log(1−1+ezez)=−log(1−1+e−z1)
公式(3)(3)(3)计算为−log(1−σ(z))- \log(1 - \sigma(z))−log(1−σ(z)),结果也是一样的。
我们证明了这种代码实现的正确性。
下面就可以为我们的metagrad实现数值稳定版的BCE损失了。
等等,还需要先实现两个函数:abs
和clip
。
实现Clip操作
clip()
像个夹子,把Tensor中的值限制在最小值和最大值之间。
class Clip(_Function):def forward(ctx, x: ndarray, x_min=None, x_max=None) -> ndarray:if x_min is None:x_min = np.min(x)if x_max is None:x_max = np.max(x)ctx.save_for_backward(x, x_min, x_max)return np.clip(x, x_min, x_max)def backward(ctx, grad: ndarray) -> ndarray:x, x_min, x_max = ctx.saved_tensorsmask = (x >= x_min) * (x <= x_max)return grad * mask
只有在[x_min,x_max]
之间的元素才有梯度。
实现Abs操作
abs()
即求绝对值,图像如下:
我们知道,按照定义绝对值在000处是不可导的
ddx∣x∣=x∣x∣\frac{d}{dx} |x| = \frac{x}{|x|} dxd∣x∣=∣x∣x
因为除000是无意义的,但是我们和PyTorch的做法一致,当x=0x=0x=0时,令其导数也为000。
class Abs(_Function):def forward(ctx, x: ndarray) -> ndarray:ctx.save_for_backward(x)return np.abs(x)def backward(ctx, grad: ndarray) -> ndarray:x, = ctx.saved_tensors# x中元素为0的位置,返回0# 否则返回+1/-1return grad * np.where(x == 0, 0, x / np.abs(x))
实现稳定版BCE损失
现在实现起来就很顺畅:
def binary_cross_entropy(input: Tensor, target: Tensor, reduction: str = "mean") -> Tensor:''':param input: logits:param target: 真实标签 0或1:param reduction: binary cross entropy loss:return:'''neg_abs = - abs(input)errors = input.clip(x_min=0) - input * target + (1 + neg_abs.exp()).log()N = len(target)if reduction == "mean":loss = errors.sum() / Nelif reduction == "sum":loss = errors.sum()else:loss = errorsreturn loss
当然,为了稳妥起见,我们写一个测试用例:
def test_binary_cross_entropy():N = 10x = np.random.randn(N)y = np.random.randint(0, 1, (N,))mx = Tensor(x, requires_grad=True)tx = torch.tensor(x, dtype=torch.float32, requires_grad=True)ty = torch.tensor(y, dtype=torch.float32)mo = torch.binary_cross_entropy_with_logits(tx, ty).mean()to = F.binary_cross_entropy(mx, y)assert np.allclose(mo.data,to.numpy())mo.backward()to.backward()assert np.allclose(mx.grad.data,tx.grad.numpy())
确保它是通过的:
============================= test session starts =============================
collecting ... collected 1 itemtest_cross_entropy.py::test_binary_cross_entropy PASSED [100%]======================== 1 passed, 1 warning in 0.48s =========================
BCE损失的实现类实际上我们不需要改:
class BCELoss(_Loss):def __init__(self, reduction: str = "mean") -> None:super().__init__(reduction)def forward(self, input: Tensor, target: Tensor) -> Tensor:''':param input: logits:param target: 真实标签 0或1:return:'''return F.binary_cross_entropy(input, target, self.reduction)
其实之前的逻辑回归实现偶尔会遇到下面这个问题,其实就是数值稳定问题。
~/metagrad/examples/logistic_regression.py:43: FutureWarning: Support for multi-dimensional indexing (e.g. `obj[:, None]`) is deprecated and will be removed in a future version. Convert to a numpy array before indexing instead.y = y[:, np.newaxis]0%| | 0/200000 [00:00<?, ?it/s]~/metagrad/metagrad/ops.py:214: RuntimeWarning: divide by zero encountered in logreturn np.log(x)
~/metagrad/metagrad/ops.py:119: RuntimeWarning: invalid value encountered in multiplyreturn x * y
~/metagrad/metagrad/ops.py:218: RuntimeWarning: divide by zero encountered in true_dividereturn grad / x
~/metagrad/metagrad/ops.py:218: RuntimeWarning: invalid value encountered in true_divide
我们把模型的初始化权重打印出来,以便重现这个错误:
# print(model.linear.weight) Tensor([[-0.29942604 0.78735491]], requires_grad=True) 有问题的权重# 使用之前有问题的权重model.linear.weight.assign([[-0.29942604, 0.78735491]])print(f"using weight: {model.linear.weight}")
然后再次训练:
using weight: Tensor([[-0.29942605 0.7873549 ]], requires_grad=True)5%|▌ | 10026/200000 [00:05<01:45, 1796.27it/s]Train - Loss: 0.6216351389884949. Accuracy: 69.696969696969710%|█ | 20023/200000 [00:11<01:39, 1810.46it/s]Train - Loss: 0.6169218420982361. Accuracy: 71.7171717171717215%|█▍ | 29933/200000 [00:16<01:33, 1812.38it/s]Train - Loss: 0.6122889518737793. Accuracy: 74.7474747474747520%|██ | 40009/200000 [00:22<01:27, 1823.05it/s]Train - Loss: 0.6077350378036499. Accuracy: 79.7979797979797925%|██▌ | 50028/200000 [00:27<01:23, 1804.55it/s]Train - Loss: 0.6032587885856628. Accuracy: 81.8181818181818130%|██▉ | 59893/200000 [00:33<01:17, 1804.93it/s]Train - Loss: 0.5988588929176331. Accuracy: 81.8181818181818135%|███▍ | 69928/200000 [00:38<01:11, 1828.28it/s]Train - Loss: 0.5945340394973755. Accuracy: 82.8282828282828240%|████ | 80000/200000 [00:44<01:06, 1811.00it/s]Train - Loss: 0.5902827978134155. Accuracy: 83.8383838383838345%|████▌ | 90037/200000 [00:50<01:01, 1782.85it/s]Train - Loss: 0.5861039757728577. Accuracy: 85.8585858585858550%|████▉ | 99996/200000 [00:55<00:55, 1789.86it/s]Train - Loss: 0.5819962620735168. Accuracy: 86.8686868686868655%|█████▍ | 109880/200000 [01:01<00:50, 1772.11it/s]Train - Loss: 0.577958345413208. Accuracy: 86.8686868686868660%|█████▉ | 119899/200000 [01:06<00:45, 1760.05it/s]Train - Loss: 0.5739889144897461. Accuracy: 87.8787878787878865%|██████▍ | 129886/200000 [01:12<00:40, 1752.00it/s]Train - Loss: 0.5700867176055908. Accuracy: 87.8787878787878870%|██████▉ | 139954/200000 [01:18<00:33, 1791.69it/s]Train - Loss: 0.5662506222724915. Accuracy: 88.8888888888888975%|███████▍ | 149921/200000 [01:24<00:31, 1572.65it/s]Train - Loss: 0.5624792575836182. Accuracy: 89.898989898989980%|████████ | 160025/200000 [01:30<00:22, 1805.44it/s]Train - Loss: 0.5587714314460754. Accuracy: 91.9191919191919285%|████████▍ | 169954/200000 [01:35<00:16, 1778.81it/s]Train - Loss: 0.5551260113716125. Accuracy: 90.909090909090990%|████████▉ | 179929/200000 [01:41<00:11, 1777.69it/s]Train - Loss: 0.551541805267334. Accuracy: 91.9191919191919295%|█████████▍| 189974/200000 [01:47<00:05, 1825.03it/s]Train - Loss: 0.5480176210403442. Accuracy: 90.9090909090909100%|██████████| 200000/200000 [01:52<00:00, 1774.94it/s]
Train - Loss: 0.544552206993103. Accuracy: 90.9090909090909
这次没有问题了。
总结
本文我们实现了数值稳定的逻辑回归损失,下篇文章我们来实现更常用的数值稳定版Softmax回归损失。
References
- How do Tensorflow and Keras implement Binary Classification and the Binary Cross-Entropy function?
从零实现深度学习框架——逻辑回归中的数值稳定相关推荐
- python学习框架图-从零搭建深度学习框架(二)用Python实现计算图和自动微分
我们在上一篇文章<从零搭建深度学习框架(一)用NumPy实现GAN>中用Python+NumPy实现了一个简单的GAN模型,并大致设想了一下深度学习框架需要实现的主要功能.其中,不确定性最 ...
- 深度学习原理-----逻辑回归算法
系列文章目录 深度学习原理-----线性回归+梯度下降法 深度学习原理-----逻辑回归算法 深度学习原理-----全连接神经网络 深度学习原理-----卷积神经网络 深度学习原理-----循环神经网 ...
- 从零实现深度学习框架——GloVe从理论到实战
引言 本着"凡我不能创造的,我就不能理解"的思想,本系列文章会基于纯Python以及NumPy从零创建自己的深度学习框架,该框架类似PyTorch能实现自动求导.
- 从零实现深度学习框架——Seq2Seq从理论到实战【实战】
引言 本着"凡我不能创造的,我就不能理解"的思想,本系列文章会基于纯Python以及NumPy从零创建自己的深度学习框架,该框架类似PyTorch能实现自动求导.
- 从零实现深度学习框架——RNN从理论到实战【理论】
引言 本着"凡我不能创造的,我就不能理解"的思想,本系列文章会基于纯Python以及NumPy从零创建自己的深度学习框架,该框架类似PyTorch能实现自动求导.
- 从零实现深度学习框架——深入浅出Word2vec(下)
引言 本着"凡我不能创造的,我就不能理解"的思想,本系列文章会基于纯Python以及NumPy从零创建自己的深度学习框架,该框架类似PyTorch能实现自动求导. 要深入理解深度学 ...
- 从零实现深度学习框架——从共现矩阵到点互信息
引言 本着"凡我不能创造的,我就不能理解"的思想,本系列文章会基于纯Python以及NumPy从零创建自己的深度学习框架,该框架类似PyTorch能实现自动求导.
- 从零实现深度学习框架——LSTM从理论到实战【理论】
引言 本着"凡我不能创造的,我就不能理解"的思想,本系列文章会基于纯Python以及NumPy从零创建自己的深度学习框架,该框架类似PyTorch能实现自动求导.
- 【深度学习】逻辑回归及其损失函数的理解
文章目录 一.什么是二分类与逻辑回归(Logistic Regression)? 二.逻辑回归基本思路 三.定义损失函数(Loss Function) 一.什么是二分类与逻辑回归(Logistic R ...
- 快速了解深度学习框架--tensorflow(更新中)
深度学习框架(工具)简单来说即库,需要import,比如tensorflow,Caffe- 深度学习框架提供了一系列的深度学习的组件(对于通用的算法,里面会有实现),当需要使用新的算法的时候就需要用户 ...
最新文章
- 同样是搞Java,年薪15W和50W的到底差在哪里?
- ArcGIS 网络分析[8.2] 资料2 使用IDatasetContainer2接口的CreateDataset方法创建网络数据集...
- Android中Context详解
- OpenJTAG调试S3C2440裸板程序
- Python机器学习:梯度下降法003线性回归中的梯度下降法
- 为什么现在很小的孩子都会玩游戏,他们真的看得懂吗?
- Mybatis(四) 高级映射,一对一,一对多,多对多映射
- MATLAB---CAD绘制Bezier曲线算法
- 内存颗粒和闪存颗粒的区别_闪存颗粒与内存颗粒的不同
- 利用雪碧图及css自制的动态变色导航栏
- C++之友元:是朋友(friend)也是破坏者
- 【CSS+HTML】实现鼠标失去鼠标焦点动画
- Android studio成品 记账本(附带文档)
- I2C接口简介和时序
- uniapp抽奖组件-动画效果之各类抽奖(跳跃)
- 使用bat脚本运行jar程序 cmd下解决乱码问题
- IHERB上婴幼儿营养补充保健系列介绍
- 【Box3引擎摄像机扩展】Box3CameraLag Box3CameraBessel
- Mac电脑安装其他系统
- 在深圳南山科技园的两年
热门文章
- jQuery之美,第一次...
- 读取本地IP地址和子网页码
- Scrapy框架学习(二)
- 【图灵学院01】Java程序员开发效率工具IntelliJ IDEA使用
- HttpClient中post请求http、https示例
- 海思hi3518 移植live555 实现H264的RTSP播放
- C编程语言中16位整型数据的取值范围介绍
- Hdoj 1064 Financial Management
- [ZJOI2014]力
- Xcode 8.1 : Unable to read from device