在跨库使用线性回归时,不少工程师遇到过诡异现象:同一份数据集,在Python的scikit-learn、statsmodels以及自行实现的梯度下降脚本中训练,得到的回归系数甚至预测值都有明显偏差。这种现象通常不是算法逻辑错误,而是特征尺度差异造成的数值不稳定在背后作祟。不同数值计算库对矩阵求逆或分解的底层实现各有选择,当输入特征量纲悬殊,设计矩阵的条件数恶化,微小的浮点误差会被不成比例地放大,最终让各个库走向不同的局部数值结果。

一、为什么特征尺度会导致结果不一致
线性回归闭式解依赖对设计矩阵求逆或正交分解。若特征A的取值范围在零点几到几之间,而特征B的取值高达数万,那么设计矩阵中对应列的数值规模相差悬殊。这种列间尺度失衡会直接推高矩阵的条件数,条件数越大,矩阵越接近奇异,求解过程对舍入误差越敏感。
不同库对此的处理方式并不统一。scikit-learn的LinearRegression底层调用LAPACK的奇异值分解或QR分解,statsmodels走的是OLS基于QR的精确最小二乘,而手写梯度下降依赖学习率和迭代收敛,对特征尺度更是极度敏感。当矩阵本身已处于数值边缘,分解算法中微小的优先级差异就会让系数估计出现肉眼可见的偏移。
1.1 条件数直观理解
条件数可理解为矩阵对输入扰动的放大倍数。假设条件数为一千,那么数据或计算中的相对误差可能被放大一千倍。在特征未标准化时,混合量纲特征轻松把条件数推到九次方级别,这时任何库都很难给出完全一致的稳定解。
我们可以用numpy快速算一下条件数,观察尺度影响:
import numpy as np
# 未标准化特征:一列小数值,一列大数值
X_raw = np.array([[0.1, 20000], [0.3, 50000], [0.2, 30000]], dtype=float)
y = np.array([1.0, 2.0, 1.5])
# 计算设计矩阵条件数
X_aug = np.hstack([np.ones((X_raw.shape[0], 1)), X_raw])
cond_raw = np.linalg.cond(X_aug)
print("原始条件数:", cond_raw)
# 简单标准化
mean = X_raw.mean(axis=0)
std = X_raw.std(axis=0)
X_scaled = (X_raw - mean) / std
X_aug_s = np.hstack([np.ones((X_scaled.shape[0], 1)), X_scaled])
cond_scaled = np.linalg.cond(X_aug_s)
print("标准化后条件数:", cond_scaled)
上述代码运行后,你会发现原始条件数远大于标准化后的数值。条件数下降代表矩阵良态化,后续求解自然更稳。
二、多库结果对比与诊断方法
遇到系数不一致,第一步应固定数据并输出各库的设计矩阵条件数,以及是否开启了截距项。很多差异来自某个库默认中心化而另一个没有。确认接口一致后,再观察标准化前后的变化。
下面用scikit-learn和statsmodels做同一份数据的对照,注意我们在两个库中均先标准化:
from sklearn.linear_model import LinearRegression
import statsmodels.api as sm
from sklearn.preprocessing import StandardScaler
X = np.array([[0.1, 20000], [0.3, 50000], [0.2, 30000]], dtype=float)
y = np.array([1.0, 2.0, 1.5])
scaler = StandardScaler()
Xs = scaler.fit_transform(X)
# sklearn
sk = LinearRegression(fit_intercept=True).fit(Xs, y)
print("sklearn系数:", sk.coef_, sk.intercept_)
# statsmodels
Xs_aug = sm.add_constant(Xs)
sm_model = sm.OLS(y, Xs_aug).fit()
print("statsmodels系数:", sm_model.params)
标准化后,两者系数会高度吻合。若仍不一致,可检查statsmodels是否用了不同的求解容差,或sklearn是否因多线程BLAS产生非确定性累加顺序。
2.1 手写梯度下降的额外陷阱
梯度下降对尺度更敏感。未标准化的特征让损失函数等高线呈狭长椭圆,导致震荡或收敛极慢,最终停在偏离解析解的位置。下列代码展示尺度如何影响收敛:
def grad_descent(X, y, lr=0.01, epochs=1000):
n, d = X.shape
w = np.zeros(d)
for _ in range(epochs):
pred = X.dot(w)
grad = X.T.dot(pred - y) / n
w -= lr * grad
return w
# 未标准化直接下降
X_aug_raw = np.hstack([np.ones((X.shape[0], 1)), X])
w_raw = grad_descent(X_aug_raw, y, lr=0.0001, epochs=5000)
print("原始尺度梯度下降:", w_raw)
X_aug_s = np.hstack([np.ones((Xs.shape[0], 1)), Xs])
w_s = grad_descent(X_aug_s, y, lr=0.1, epochs=5000)
print("标准化后梯度下降:", w_s)
可以看到,原始尺度下必须极小学习率才能不发散,而标准化后正常学习率即可逼近稳定解。这解释了为什么自制脚本和成熟库差距最大。
三、修复方案与工程建议
最根本的修复是在建模前统一做特征标准化或归一化。对线性回归而言,标准化不改变系数解释的相对关系,只把量纲拉到同一区间,从根源降低条件数。若业务要求系数对应原始量纲,可在训练后反变换回原始空间。
当特征间存在共线性且尺度难调,可引入岭回归等带L2惩罚的模型,提升数值稳定性。下表列出常见处理动作:
| 问题场景 | 推荐动作 | 预期效果 |
|---|---|---|
| 多库系数偏差大 | 先标准化再训练 | 条件数下降,结果对齐 |
| 梯度下降不收敛 | 特征缩放加合适学习率 | 损失曲面变圆,快速收敛 |
| 存在共线性 | 岭回归或主成分回归 | 矩阵可逆性增强 |
3.1 反变换保留可解释性
标准化训练后,若需向业务方汇报原始量纲系数,可用公式回推。设标准化为(x-mean)/std,则原始系数等于标准化系数除以std,截距相应调整。这样既有数值稳定,又不丢解释力。
在工程落地时,建议把缩放器与模型一起持久化,避免线上输入未走同样变换而再次引发偏差。整套流水线封装后,多库不一致的问题基本可彻底规避。
linear_regressionfeature_scalingnumerical_stability修改时间:2026-08-06 16:42:35