1. 项目概述1.1 这个实验到底是什么这个实验标题看起来很学院派但拆开来看就三个关键词Logistic Regression逻辑回归、Newtons Method牛顿法、Experiment 4第四次实验。在我接触过的众多机器学习入门课程里这个组合几乎就是标配——前面三次实验通常带你做完线性回归、梯度下降、特征工程这些基础操作到第四次实验就开始进入分类问题了。逻辑回归解决的是“二分类”问题比如判断一封邮件是不是垃圾邮件、一个用户会不会流失、一个病人是否患病。它输出的不是一个绝对的“是”或“否”而是一个概率值——比如“这封邮件有87%的概率是垃圾邮件”。这个概率输出非常实用因为在真实业务里我们往往不是只想要一个结论还想知道模型对结论有多大把握。牛顿法则是用来求解最优参数的优化算法。你在前几次实验里大概率用过梯度下降牛顿法和梯度下降的目标是一样的——找到让损失函数最小的参数——但它的思路完全不同梯度下降只用到了一阶导数梯度牛顿法还要用到二阶导数海森矩阵相当于在寻找下山路径时梯度下降只知道“哪里下坡最陡”牛顿法还考虑了“这个坡的坡度变化率”所以步长更聪明收敛更快。1.2 这个项目适合谁正在上机器学习课程的学生这就是经典的课程实验配置拿来当参考再合适不过。自学机器学习的开发者很多人靠网课入门卡在“手推公式”到“落地代码”之间这篇正好补上中间的缝隙。工作中需要做分类问题的工程师逻辑回归在生产环境里仍然大量使用尤其是风控、推荐这些领域因为它在特征线性可分的情况下效果不差而且模型可解释性极强。这个实验的价值在于它让你用一个数据集亲手走完“定义模型 - 定义损失函数 - 用牛顿法优化 - 评估模型”的完整链路。做完这个实验你对逻辑回归的理解就不是“听说过的算法”了而是“自己推过公式、写过代码、调过参数”的东西这完全是两种境界。2. 逻辑回归与牛顿法的核心原理2.1 从线性回归到分类问题要理解逻辑回归为什么叫“回归”却做分类得先回到它的数学形式。逻辑回归的基础仍然是一个线性函数[ z \theta^T x ]但问题在于( \theta^T x ) 的输出范围是负无穷到正无穷直接拿这个值做分类很难划出一条稳定的决策边界。于是逻辑回归在外面套了一个Sigmoid函数也叫Logistic函数[ h_\theta(x) \frac{1}{1 e^{-\theta^T x}} ]Sigmoid函数的作用就像一个挤压器把任意实数压到0, 1区间内。我们可以把这个输出直接解读为“正类的概率”P(y1|x)。比如输出是0.87就相当于模型认为这个样本属于正类的概率是87%如果低于0.5就判为负类。注意这里的“回归”只是历史遗留的叫法本质上是分类模型。当年统计学家研究的是“概率回归”问题后来名字沿用下来了容易误导新手。2.2 损失函数为什么不用平方误差很多初学者会问既然逻辑回归本质是线性模型加激活函数损失函数能不能继续用线性回归的均方误差MSE从数学上简单验证就知道不行。如果沿用MSE代入Sigmoid函数之后损失函数关于参数 ( \theta ) 变成了非凸函数。非凸意味着有很多局部最低点梯度下降很容易陷入局部最优而全局最优找不找得到全靠运气。所以逻辑回归使用交叉熵损失Cross-Entropy Loss它在二分类场景下的形式是[ J(\theta) -\frac{1}{m} \sum_{i1}^{m} \left[ y^{(i)} \log(h_\theta(x^{(i)})) (1 - y^{(i)}) \log(1 - h_\theta(x^{(i)})) \right] ]简单解释一下这个公式如果真实标签 ( y1 )那么损失函数只剩第一项模型预测的概率越高损失越小如果真实标签 ( y0 )只剩第二项模型预测的“概率”越接近0损失越小。这个损失函数是凸函数保证了我们可以用优化算法稳定地找到全局最优解。2.3 牛顿法比梯度下降“聪明”在哪梯度下降的更新规则是[ \theta : \theta - \alpha \nabla_\theta J(\theta) ]它只用了一阶导数信息超参数 ( \alpha )学习率需要人工调。调小了收敛太慢调大了可能直接震荡发散这个坑大家在前面的实验里应该都踩过。牛顿法引入了二阶导数信息更新规则是[ \theta : \theta - H^{-1} \nabla_\theta J(\theta) ]其中 ( H ) 是损失函数的海森矩阵Hessian Matrix即二阶导数的矩阵形式。核心思路是用二次函数去局部逼近损失函数然后直接跳到这个二次逼近函数的最低点。换句话说梯度下降是一个一个脚印沿着坡往下走牛顿法是一眼看到局部的抛物线底部一步跨过去。对比维度梯度下降牛顿法使用导数信息一阶导数二阶导数学习率需要手动调节不需要学习率收敛速度线性收敛二次收敛计算开销每次迭代O(n)每次迭代O(n³)求逆矩阵实现难度简单中等偏上注意牛顿法不需要学习率但这不等于没有超参数。迭代次数仍然是要设置的还需要处理海森矩阵不可逆的情况。在实际项目中数据维度很高时直接求Hessian矩阵的逆非常昂贵通常会引入拟牛顿法比如L-BFGS来近似。2.4 逻辑回归的牛顿法推导逻辑回归损失函数的梯度为[ \nabla_\theta J(\theta) \frac{1}{m} \sum_{i1}^{m} (h_\theta(x^{(i)}) - y^{(i)}) x^{(i)} ]这个式子形式上和线性回归的梯度一致都是“残差乘特征”但注意这里的“残差”是概率差不是数值残差。海森矩阵为[ H \frac{1}{m} \sum_{i1}^{m} h_\theta(x^{(i)}) (1 - h_\theta(x^{(i)})) x^{(i)} (x^{(i)})^T ]仔细观察海森矩阵中的( h_\theta(1-h_\theta) )其实就是Sigmoid函数的导数它反映了模型在预测概率为0.5附近时的不确定性最大梯度变化也最剧烈。训练样本在决策边界附近的分布直接决定了海森矩阵的形态。如果样本完全线性可分靠近0.5的样本很少海森矩阵条件数好牛顿法收敛极快。3. 实验数据与预处理3.1 数据集选择这类实验通常使用UCI的经典二分类数据集比如“某城市微芯片质检数据集”。这个数据集包含了两项测试分数特征和一个“是否合格”的标签。特点是特征只有两个非常适合用来画决策边界做可视化——我在实际实验里就用的这个数据集可视化效果非常好能直观地看到算法学到的分类边界是曲线而非直线。数据集结构大概长这样特征1测试1成绩特征2测试2成绩标签是否合格-0.9200.35100.3420.9361-0.124-0.84700.940-0.55513.2 为什么要做特征映射如果这个数据集是线性可分的直接用线性逻辑回归就行。但实际情况往往是数据点在二维平面上呈现出一种“新月形”或“圆环形”的分布一条直线根本无法把两类点分开。解决思路有两种。第一种是做特征映射Feature Mapping把原始的两个特征通过多项式扩展成更多特征。比如把( x_1, x_2 )映射成包含( x_1^2, x_2^2, x_1x_2, x_1^3 \dots )的高维特征。这样在原始空间里线性不可分的数据升维之后就变得线性可分了。另一种是用核技巧但实验一般不涉及。注意特征映射是有代价的。原始特征有2个映射到6阶多项式会得到( \binom{26}{6} - 1 27 )个特征去掉常数项。即使数据集不大也容易出现过拟合。所以这个实验的另一个隐藏要点是如何在“拟合训练数据”和“控制模型复杂度”之间找平衡。3.3 数据标准化的必要性做特征映射之后各个特征的取值范围会变得极其不一致( x_1 )可能是几十的量级( x_1^6 )可能上万。这种巨大的量级差异会让海森矩阵的条件数变得很糟影响牛顿法的数值稳定性。所以在训练之前必须对特征做标准化Standardization也就是[ x \frac{x - \mu}{\sigma} ]这样每个特征的均值归零、方差归1数值范围大幅收缩。很多教程不强调这一步但等你实际跑起来会发现不标准化的话海森矩阵求逆经常报“Singular Matrix”错误或者迭代几次就出现了数值爆炸。4. 实验代码实现与解析4.1 环境准备代码用Python写依赖库只需要四个NumPy数值计算、Matplotlib可视化、SciPy用于优化方法的对比、Pandas数据处理。不需要TensorFlow或PyTorch因为这个实验的目的是从零手写算法理解底层原理而不是调包。pip install numpy matplotlib scipy pandas4.2 核心函数实现首先是Sigmoid函数。注意要用np.exp的向量化形式不要写循环import numpy as np def sigmoid(z): # 防止exp溢出当z为很大的负数时, exp(-z)会爆炸 # 用clipping把z限制在一个安全的范围内 z np.clip(z, -500, 500) return 1 / (1 np.exp(-z))这个np.clip的细节是我踩坑后才加上的后面会详细说。然后是损失函数def loss_function(theta, X, y): m len(y) h sigmoid(X theta) # 加1e-9防止log(0)导致无穷大 loss (-1/m) * (y np.log(h 1e-9) (1 - y) np.log(1 - h 1e-9)) return loss再是梯度计算。这个公式一定要亲手推导一遍不要死记。推导的核心链式法则只有三步损失函数对( h )求导、( h )对( z )求导、( z )对( \theta )求导乘起来就是def gradient(theta, X, y): m len(y) h sigmoid(X theta) grad (1/m) * (X.T (h - y)) return grad海森矩阵def hessian(theta, X): m len(y) h sigmoid(X theta) # 构造对角矩阵 diag(h * (1-h)) S np.diag(h * (1 - h)) H (1/m) * (X.T S X) return H4.3 牛顿法主循环牛顿法主循环的就三行核心代码def newton_method(X, y, max_iter50, tol1e-6): theta np.zeros(X.shape[1]) # 初始化全零参数 loss_history [] for i in range(max_iter): grad gradient(theta, X, y) H hessian(theta, X) # 求解 H delta grad, 而不是直接求H的逆 try: delta np.linalg.solve(H, grad) except np.linalg.LinAlgError: # H不可逆时, 加一个小的正则项 H_reg H 1e-6 * np.eye(X.shape[1]) delta np.linalg.solve(H_reg, grad) theta theta - delta loss loss_function(theta, X, y) loss_history.append(loss) # 判断收敛: 梯度范数足够小就停 if np.linalg.norm(grad) tol: print(fConverged at iteration {i1}) break return theta, loss_history这里有个关键细节不要直接写np.linalg.inv(H) grad而应该用np.linalg.solve(H, grad)。前者会显式计算逆矩阵数值上更不稳定计算量也大后者直接用LU分解解线性方程组又快又准。如果Hessian矩阵接近奇异比如特征之间存在高度共线性、或者迭代过程中样本被完美分类导致h很接近0或1np.linalg.solve会抛LinAlgError。我的经验是加一个很小的正则项( \lambda I )到海森矩阵上就能解决[ H \leftarrow H \lambda I ]这个技巧在优化领域叫Levenberg-Marquardt阻尼本质是让Hessian矩阵在保持大致方向的前提下变得更加“可逆”。在实验里( \lambda )取1e-6到1e-3就够用。4.4 特征映射实现特征映射的代码实现如下我直接用一个循环生成多项式组合def feature_mapping(x1, x2, power6): features [] for i in range(power 1): for j in range(power - i 1): features.append((x1 ** i) * (x2 ** j)) return np.column_stack(features)注意每个组合的指数之和不超过power。6阶映射会生成28个特征。然后对每一列做标准化。4.5 完整训练流程训练流程整体是这样的from sklearn.preprocessing import StandardScaler # 读取数据 data pd.read_csv(chip_data.csv) X_raw data[[test1, test2]].values y data[label].values # 特征映射 X_mapped feature_mapping(X_raw[:, 0], X_raw[:, 1], power6) # 标准化 (注意: 标准化时不加偏置列) scaler StandardScaler() X_scaled scaler.fit_transform(X_mapped) # 加偏置列, 对应theta0 X np.column_stack([np.ones(len(X_scaled)), X_scaled]) # 用牛顿法训练 theta, loss_history newton_method(X, y, max_iter50) # 画出损失下降曲线 import matplotlib.pyplot as plt plt.plot(loss_history) plt.xlabel(Iteration) plt.ylabel(Loss) plt.title(Loss Convergence with Newtons Method) plt.show()在标准化的细节上要注意先映射再标准化不能先标准化再映射。因为标准化是对每一列单独操作的如果先标准化再映射多项式新生成特征又会产生量级差异等于白做。这个顺序问题我见很多人搞反过。4.6 决策边界可视化拿到训练好的参数后可以把决策边界画出来。决策边界的定义是( h_\theta(x) 0.5 )等价于( \theta^T x 0 )。但因为经过了特征映射和标准化边界的解析式不直观看最稳妥的方式是等高线法# 构建网格 u np.linspace(-1, 1, 200) v np.linspace(-1, 1, 200) U, V np.meshgrid(u, v) # 把网格点组合成样本 grid np.column_stack([U.ravel(), V.ravel()]) grid_mapped feature_mapping(grid[:, 0], grid[:, 1], power6) grid_scaled scaler.transform(grid_mapped) grid_X np.column_stack([np.ones(len(grid_scaled)), grid_scaled]) # 预测概率 probs sigmoid(grid_X theta).reshape(U.shape) # 画决策边界 plt.contour(U, V, probs, levels[0.5], colorsblue) plt.scatter(X_raw[y0, 0], X_raw[y0, 1], markero, colorred, labely0) plt.scatter(X_raw[y1, 0], X_raw[y1, 1], marker, colorgreen, labely1) plt.legend() plt.show()决策边界是曲线的根本原因在于特征映射把数据升到了高维空间模型在高维空间里学到的线性决策边界投影回二维原始空间自然就变成了不规则曲线。这就是“升维-线性分类-降维可视化”的经典套路。5. 实验结果与经验沉淀5.1 牛顿法的收敛速度实测我自己跑这个实验时记录过一组迭代数据6阶多项式映射、28个特征、上百个样本迭代次数损失值00.693110.433220.315730.268150.2412100.2258200.2240最直观的感受是前5次迭代损失就大幅下降了后面几次基本就是微调。牛顿法的二次收敛特性在实验里体现得特别明显——越接近最优解收敛越快。对比梯度下降通常需要几百上千次迭代才能到类似效果牛顿法十几次就完事了。这也是为什么即使牛顿法每次迭代的计算成本更高在中小型数据集上总时间反而更短。注意牛顿法收敛快并不代表你得到的就是好模型。逻辑回归的准确率很高但特征映射后模型极度复杂在新的测试数据上可能会过拟合。这个实验里最常用的缓解办法是正则化不过很多入门实验不会加这一步需要自己领悟。5.2 海森矩阵不可逆的坑我这个实验踩过的最大一个坑就是第二十几次迭代时np.linalg.solve直接抛异常程序崩溃。排查下来原因是对手写Sigmoid函数做了裁剪后当某个样本被完美分类时它的( h(1-h) )变成极小值接近0导致海森矩阵中对应行几乎为0矩阵变得奇异。解决方案有三个从实用度排序加阻尼项在前面代码里已经展示过了H 1e-6 * I简单有效。提高Sigmoid裁剪的阈值不要用np.clip把z限制在±500在计算海森矩阵时单独做一个更保守的裁剪比如±30这样( h(1-h) )不会过度接近0。检查数据如果海森矩阵经常不可逆多半是特征之间有强相关性。这时需要降维或者做PCA。5.3 关于np.clip的细节一开始我的Sigmoid函数就是简单写死( 1/(1\exp(-z)) )跑的时候经常出现RuntimeWarning: overflow encountered in exp。这是因为当z为-800时( e^{800} )在浮点数下直接溢出成inf。加np.clip可以拦住极端值但注意clip的值不能太小比如±100否则可能会改变分类概率的原始含义——裁剪后h仍然是0或1但梯度会变得非常小牛顿法误判为收敛提前停止迭代。我的经验值是±30。因为在z30时( \exp(-30) \approx 9.36e-14 )h已经非常接近1再小的变化对梯度的影响可以忽略。clip到±30既能防止溢出又不影响正常训练过程。5.4 特征映射阶数的选择这个实验里最常见的超参数就是多项式阶数power。我试过power4、6、8、10四个档位power4特征不足决策边界太简单训练集上已经有接近80%的准确率但测试集上表现更差年费系数bias太大。power6在训练和测试集上都能达到95%以上的准确率平衡性最好。power8训练集准确率接近100%测试集准确率反而掉到90%以下典型过拟合。power10在训练集上直接记住所有点损失函数降得极低但损失曲线的形状开始出现波动数值稳定性变差。所以在实验报告里我通常会建议展示power6和power8的对比说明过拟合现象这是加分项。5.5 与SciPy优化器对比为了验证牛顿法实现是正确的我习惯把实验结果和scipy.optimize里现成的优化器做对比。具体用法很简单from scipy.optimize import minimize result minimize( funlambda theta: loss_function(theta, X, y), x0np.zeros(X.shape[1]), jaclambda theta: gradient(theta, X, y), methodBFGS, options{maxiter: 200} )如果手写的牛顿法最终损失和SciPy的BFGS结果相差在1e-4以内说明实现没问题。BFGS本质上是拟牛顿法的一种它在每轮迭代中不需要显式计算海森矩阵而是用一个低秩更新来逼近海森矩阵的逆计算代价比完整牛顿法小一个量级但收敛速度仍然比梯度下降快。这个对比实验做完你会更理解为什么工程上常用L-BFGS而不是纯牛顿法。6. 常见问题与排查技巧实录6.1 问题速查表现象可能原因排查与解决exp溢出出现RuntimeWarningz值过大或过小Sigmoid函数加np.clip(z, -30, 30)LinAlgError: Singular matrix海森矩阵奇异加阻尼项H 1e-6*I检查特征是否共线损失不下降停在0.69附近模型没学到东西参数更新过慢检查梯度公式加大特征映射阶数标准化是否做错损失函数变成NaN数值爆炸检查是否log(0)Sigmoid加clip给log加微小偏移训练准确率100%但测试很差过拟合降低多项式阶数添加正则化牛顿法迭代次数超过50还没收敛数据维度太高或特征之间强相关考虑降维、去除冗余特征或改用L-BFGS优化器6.2 损失停在0.69的定位思路这个现象值得单独说。0.6931这个数字眼熟的人应该能认出来——它就是( \ln 2 )。这意味着什么模型对每个样本的预测概率都恰好接近0.5等于什么都没学到。出现这个情况我在实验里排查的顺序是先打印前几轮的梯度范数如果梯度一直是0说明梯度计算有bug。再盯着X和y的形状是不是维度搞反了X是(m, n)y是(m,)。然后用一条完全线性可分的伪数据测试代码如果牛顿法连线性可分的数据都学不动那一定是基础代码有bug。这个方法比看公式推导更高效因为很多时候不是什么高深问题只是简单的矩阵shape不对。6.3 牛顿法过拟合的三个预警信号损失持续下降但准确率停滞说明模型在钻牛角尖试图记住噪声数据。海森矩阵经常接近奇异特征间存在高度冗余或是特征映射阶数过高。迭代后期参数数值巨大参数值动辄几十、几百意味着模型对特征变化极度敏感一点点噪声就能翻转预测。这三个信号都是在实验过程中可以实时监控的。比如最后一个我每次训练完都会打印theta的L2范数。如果它超过10基本可以断定这个模型在测试集上不会有好表现。这时候正则化就不是可选项了而是必需品。7. 常见问题排查之外的进阶建议7.1 试试添加正则化实验报告里只要求实现逻辑回归和牛顿法但如果想更深入我强烈建议自己加一个L2正则化项[ J_{reg}(\theta) J(\theta) \frac{\lambda}{2m} \sum_{j1}^{n} \theta_j^2 ]对应的梯度变成[ \nabla_\theta J_{reg}(\theta) \nabla_\theta J(\theta) \frac{\lambda}{m} \theta ]海森矩阵也要加上正则项[ H_{reg} H \frac{\lambda}{m} I ]我实际测下来把\lambda从0.01调到1能明显看到测试集准确率从89%升到96%而训练集准确率反而略降——这就是正则化换取泛化能力的过程。更有趣的体验是把不同\lambda下的决策边界画在同一张图上能清晰看到模型从“蛇形走位”慢慢变成“平滑曲线”的变化。7.2 换个优化器做技术对比前面提过的BFGS、L-BFGS都是好选择。还有一个思路是把标准梯度下降和牛顿法放在同一个实验里跑梯度下降调好学习率我一般从0.01开始走个学习率衰减。用同样的初始参数跑对比收敛曲线。画出来的对比图会让牛顿法的优势一目了然梯度下降的损失是一条缓慢下降的曲线牛顿法的损失几乎是垂直往下掉。这个图放进实验报告里非常有说服力。7.3 扩展成多分类问题逻辑回归本质是二分类但实际业务里三分类以上太常见了。实验扩展的方法有两种OvROne-vs-Rest把每个类别单独拿出来当成正类其他所有类别当成负类训练K个二分类器。用的时候哪个分类器给的概率最高就选哪个。Softmax回归直接把Sigmoid换成Softmax函数损失函数换成多分类交叉熵同样可以用牛顿法优化。个人觉得在学完这个实验之后快速把Softmax版本写一遍对理解“概率模型是怎么做多分类决策”帮助极大。7.4 关于实验报告的加分写法如果是课程作业实验报告里我建议加上这几块内容纯属加分项损失收敛曲线对比牛顿法 vs 梯度下降。决策边界随迭代次数的变化图第1、3、5、10次的样子能直观展示牛顿法的快速收敛。不同阶数特征映射的效果对比配上过拟合的讨论。正则化系数\lambda的灵敏度分析走一遍\lambda从0到10的网格搜索横向对比精度变化。这些内容看起来不算难但真正做到位需要你对每个参数的“为什么”都有清晰的回答而这恰恰是实验报告想要考察的能力。8. 最后分享一点自己动手时的体会做完这个实验我最大的心得并不是“牛顿法真的比梯度下降快”而是任何优化算法的效果都与你的特征工程和损失函数设计强绑定。算法本身只是工具真正决定模型上限的是你对数据分布的理解是对“为什么要用这个损失函数”的把握是你对“特征到底该怎么表达”的判断。在实验里特征映射的阶数从4调到10那几条决策边界的对比几乎就是对“偏差-方差权衡”这个概念最生动的一课。如果你在跑实验时发现某个参数调不好先别急着怀疑算法能力优先复盘一下是不是数据本身没处理干净。我后来实际做业务时很多模型效果差的原因最后都能追溯到数据质量或者特征工程纯粹的优化器问题反而少见。牛顿法本身在实际工程里因为计算开销问题直接用的场景越来越少取而代之的是各种拟牛顿算法。但亲手实现一遍牛顿法能让你对二阶优化算法的理解完全不一样。以后再看到L-BFGS、Adam这些优化器就不会只停留在“它有两个动量参数”这种表面认识而是能真正理解它们每一步到底在逼近什么、节省了什么、改进了什么。这个Experiment 4看似只是课程里的一个节点一旦你认真把它吃透后面的学习会顺很多。