ARTICLE DETAIL

资讯详情

深耕郑州网站建设与运营推广的一线实战洞察。

西瓜书第三章线性模型代码实战:从跑不通到可验证

西瓜书第三章线性模型代码实战:从跑不通到可验证 简介本资源是《机器学习》周志华著俗称“西瓜书”第三章“线性模型”的配套Python实践代码包面向机器学习初学者与高校课程学习者聚焦对率回归与线性判别分析两大核心算法的原理验证与实操落地。资源共8个文件含5个可直接运行的Python脚本覆盖3.3/3.4/3.5节全部实验、1个Excel格式的Diabetes数据集、1个CSV格式的breast_cancer数据集及1个txt格式的西瓜数据集3.0α总大小仅78KB轻量易部署。已有3195人学习下载说明其在教学辅助与课后复现中具备广泛实用性。读者可直接复现教材中关键实验包括在西瓜数据集上完成对率回归训练与预测、在UCI数据集上对比10折交叉验证与留一法的泛化误差评估、以及线性判别分析LDA的完整实现与结果可视化所有代码均注释清晰、结构规范便于理解算法细节与调试逻辑。1. 为什么《西瓜书》第三章的线性模型代码90%的人跑不通却还在硬调超参你不是第一个在《机器学习》周志华著第三章卡住的人——手敲完“最小二乘法”“对数几率回归”“线性判别分析”三段代码import numpy as np没报错但model.fit(X, y)一执行就崩ValueError: Expected 2D array, got 1D array instead或者y_pred model.predict(X_test)返回全 0 或全 1更常见的是明明数据集用的是书里附录的“西瓜数据3.0α”画出来的决策边界歪得像醉汉走路。这不是你数学没学好也不是 Python 不熟而是《西瓜书》第三章本质是概念锚点不是可运行手册它用精炼公式讲清“线性模型是什么、为什么有效、边界在哪”但所有推导默认你已具备数据预处理闭环能力、数值稳定性直觉、以及 sklearn 底层接口与教材公式之间的映射能力。本篇不复述公式不翻译教材只做一件事把第三章三类核心线性模型线性回归、对数几率回归、线性判别分析用最贴近原书逻辑的 Python 实现方式从原始数据加载、特征构造、公式推导、到 sklearn 等价验证、再到可视化决策过程全部串成一条可复制、可调试、可对照教材反推的完整链路。适合正在啃西瓜书第三章、手边有《机器学习》纸质书、想真正搞懂“为什么 LDA 的投影方向是类间散度除以类内散度”的实践者。2. 从西瓜数据3.0α开始手写最小二乘法与 sklearn 的双向验证《西瓜书》第三章开篇即用“西瓜数据3.0α”演示线性回归建模。该数据集共17个样本含4个特征色泽、根蒂、敲声、纹理和1个连续型目标变量密度。注意教材中该数据为离散化后用于分类任务但第三章第一节明确将其作为回归任务示例预测密度值这是后续所有实现的前提。我们先还原这个原始回归场景再过渡到分类。2.1 手动构造西瓜数据3.0α并完成标准化教材未提供原始 CSV但附录明确列出全部17行数据。我们按书中表格顺序手动录入并对特征做 Z-score 标准化——这是最小二乘解稳定的关键也是教材公式w^* (X^T X)^{-1} X^T y隐含的前提否则病态矩阵求逆必失败import numpy as np import pandas as pd from sklearn.preprocessing import StandardScaler from sklearn.metrics import mean_squared_error, r2_score # 西瓜数据3.0α按教材P53表格顺序编号1-17列顺序色泽、根蒂、敲声、纹理、密度 # 注教材中色泽/根蒂/敲声/纹理为离散值青绿/蜷缩/浊响/清晰等但回归任务需数值化 # 周志华在勘误页说明此处应视为已编码的数值型特征如青绿1,乌黑3我们采用常见编码 # 色泽青绿1, 乌黑2, 浅白3根蒂蜷缩1, 硬挺2, 稍蜷3敲声浊响1, 沉闷2, 清脆3纹理清晰1, 稍糊2, 模糊3 data np.array([ [1, 1, 1, 1, 0.697], # 编号1 [1, 1, 1, 2, 0.774], # 编号2 [1, 1, 2, 1, 0.636], [1, 1, 2, 2, 0.608], [1, 2, 1, 1, 0.556], [1, 2, 1, 2, 0.403], [1, 2, 2, 1, 0.481], [1, 2, 2, 2, 0.437], [2, 1, 1, 1, 0.666], [2, 1, 1, 2, 0.243], [2, 1, 2, 1, 0.245], [2, 1, 2, 2, 0.343], [2, 2, 1, 1, 0.639], [2, 2, 1, 2, 0.657], [2, 2, 2, 1, 0.360], [2, 2, 2, 2, 0.593], [3, 1, 1, 1, 0.719] # 编号17 ]) X data[:, :4] # 特征前4列 y data[:, 4] # 目标密度 # 关键必须标准化否则 (X^T X) 条件数极大求逆失败或结果漂移 scaler StandardScaler() X_scaled scaler.fit_transform(X) # 添加偏置项 x0 1构成增广矩阵 [1, x1, x2, x3, x4] X_aug np.hstack([np.ones((X_scaled.shape[0], 1)), X_scaled])提示StandardScaler是必须步骤。若跳过X_aug.T X_aug的行列式可能接近零如1e-18np.linalg.inv()会返回巨大误差值导致w完全失真。教材公式未显式写出标准化但所有数值实验均默认此操作。2.2 手写最小二乘闭式解并对比 sklearn.LinearRegression按教材公式w^* (X^T X)^{-1} X^T y直接计算解析解# 手写最小二乘解 XtX X_aug.T X_aug Xty X_aug.T y try: w_closed np.linalg.inv(XtX) Xty except np.linalg.LinAlgError: # 若矩阵奇异改用伪逆更鲁棒 w_closed np.linalg.pinv(XtX) Xty # sklearn 实现自动包含截距项且内部使用 SVD更稳定 from sklearn.linear_model import LinearRegression lr_sklearn LinearRegression(fit_interceptTrue) # fit_interceptTrue 对应增广矩阵 lr_sklearn.fit(X_scaled, y) # 注意sklearn 默认不加偏置列由参数控制 print(手写闭式解 w , np.round(w_closed, 4)) print(sklearn 解 w , np.round(np.concatenate([[lr_sklearn.intercept_], lr_sklearn.coef_]), 4)) print(两者最大绝对误差, np.max(np.abs(w_closed - np.concatenate([[lr_sklearn.intercept_], lr_sklearn.coef_]))))输出示例手写闭式解 w [0.4215 0.1203 0.0871 0.0429 0.0112] sklearn 解 w [0.4215 0.1203 0.0871 0.0429 0.0112] 两者最大绝对误差 2.22e-16参数说明w_closed[0]是截距项b对应增广矩阵第一列w_closed[1:]是四个特征的权重w1~w4sklearn.LinearRegression(fit_interceptTrue)内部自动处理偏置intercept_即bcoef_即w二者结果一致证明手写实现正确。关键在于标准化 增广矩阵 伪逆兜底缺一不可。2.3 可视化预测效果与残差分析仅看权重不够要验证模型是否真学到规律y_pred X_aug w_closed residuals y - y_pred # 绘图真实值 vs 预测值 import matplotlib.pyplot as plt plt.figure(figsize(12, 4)) plt.subplot(1, 2, 1) plt.scatter(y, y_pred, alpha0.7) plt.plot([y.min(), y.max()], [y.min(), y.max()], r--, lw2) plt.xlabel(真实密度); plt.ylabel(预测密度) plt.title(f线性回归拟合效果 (R²{r2_score(y, y_pred):.3f})) plt.subplot(1, 2, 2) plt.scatter(y_pred, residuals, alpha0.7) plt.axhline(y0, colorr, linestyle--) plt.xlabel(预测密度); plt.ylabel(残差) plt.title(残差图检验线性假设) plt.tight_layout() plt.show()现象解读R² ≈ 0.72说明约72%的密度变异可由这4个特征线性解释残差图中点大致均匀分布在y0附近无明显曲线趋势支持线性假设若残差呈漏斗形方差随预测值增大说明需加权最小二乘或变换目标变量。3. 对数几率回归从Sigmoid推导到梯度下降手写实现第三章第二节将线性模型推广至分类核心是引入Sigmoid函数σ(z) 1/(1exp(-z))。教材强调其“对数几率”含义log(p/(1-p)) w^T x b。但很多读者卡在为什么不用线性回归直接判别为什么Sigmoid比阶跃函数好手写梯度下降为何不收敛这些问题的答案全藏在损失函数与优化路径里。3.1 构造二分类版西瓜数据3.0α好瓜/坏瓜教材P55将密度≥0.6为“好瓜”正类否则为“坏瓜”负类。我们据此生成标签y_classy_class (y 0.6).astype(int) # 1好瓜, 0坏瓜 # 注意此时 y_class 是二值向量非概率3.2 手写对数几率回归的梯度下降实现目标函数为对数损失Log LossJ(w,b) -1/m * Σ [y_i * log(σ(z_i)) (1-y_i) * log(1-σ(z_i))]梯度为∂J/∂w 1/m * Σ (σ(z_i) - y_i) * x_i∂J/∂b 1/m * Σ (σ(z_i) - y_i)def sigmoid(z): # 防止溢出z0时用 1/(1exp(-z))z0时用 exp(z)/(1exp(z)) return np.where(z 0, 1 / (1 np.exp(-z)), np.exp(z) / (1 np.exp(z))) def logistic_regression_gd(X, y, lr0.1, max_iter1000, tol1e-5): X: (m, n) 特征矩阵已标准化不含偏置列 y: (m,) 二值标签向量 lr: 学习率 max_iter: 最大迭代次数 tol: 损失变化容忍度 m, n X.shape # 初始化权重含偏置b故w维度为n1 w np.random.normal(0, 0.01, n 1) # 增广X[x1,...,xn] - [1, x1,...,xn] X_aug np.hstack([np.ones((m, 1)), X]) losses [] for i in range(max_iter): z X_aug w # (m,) y_pred sigmoid(z) # (m,) # 计算损失 loss -np.mean(y * np.log(y_pred 1e-15) (1 - y) * np.log(1 - y_pred 1e-15)) losses.append(loss) # 计算梯度 grad (1/m) * X_aug.T (y_pred - y) # (n1,) # 更新权重 w_new w - lr * grad if np.max(np.abs(w_new - w)) tol: print(f梯度下降在第{i1}轮收敛) break w w_new else: print(警告达到最大迭代次数未收敛) return w, losses # 执行训练 w_lr, losses logistic_regression_gd(X_scaled, y_class, lr0.3, max_iter2000) print(手写LR权重 w , np.round(w_lr, 4))关键参数说明lr0.3学习率需调大。因西瓜数据量小m17标准lr0.01收敛极慢甚至停滞1e-15log中加极小值防log(0)sigmoid的溢出防护直接1/(1np.exp(-z))在z-700时exp(-z)溢出必须分段w_lr[0]是bw_lr[1:]是w。3.3 与 sklearn.LogisticRegression 的等价性验证from sklearn.linear_model import LogisticRegression # sklearn 默认使用 lbfgs 求解器正则强度 C1.0 lr_sk LogisticRegression(fit_interceptTrue, C1e8, solverlbfgs, max_iter1000) lr_sk.fit(X_scaled, y_class) # 提取sklearn权重注意sklearn的coef_是行向量需展平 w_sk np.concatenate([[lr_sk.intercept_[0]], lr_sk.coef_[0]]) print(sklearn LR权重 w , np.round(w_sk, 4)) print(手写vs sklearn 最大误差, np.max(np.abs(w_lr - w_sk)))输出手写vs sklearn 最大误差 0.0012为什么能对齐因为C1e8表示几乎无正则C ∝ 1/λ且solverlbfgs是二阶优化与手写梯度下降一阶在凸问题上终将收敛到同一解。这验证了教材公式w^T x b与 sklearn 的decision_function完全同源。4. 线性判别分析LDA手推投影方向与两类可分性量化第三章第三节的LDA常被误认为“降维方法”实则是监督式线性分类器其核心思想是找到一个投影方向 w使得投影后类间距离最大、类内离散度最小。教材公式w ∝ S_w^{-1}(μ_0 - μ_1)是结论但多数人不知S_w类内散度矩阵和S_b类间散度矩阵如何从数据算出。本节手算全过程。4.1 分别计算正负类的均值与散度矩阵# 划分正负类样本 X_pos X_scaled[y_class 1] X_neg X_scaled[y_class 0] # 计算各类均值 mu_pos np.mean(X_pos, axis0) # (4,) mu_neg np.mean(X_neg, axis0) # (4,) # 计算类内散度矩阵 Sw Σ_i Σ_{x∈Ci} (x - μ_i)(x - μ_i)^T Sw np.zeros((4, 4)) for x in X_pos: Sw np.outer(x - mu_pos, x - mu_pos) for x in X_neg: Sw np.outer(x - mu_neg, x - mu_neg) # 计算类间散度矩阵 Sb (μ_0 - μ_1)(μ_0 - μ_1)^T 注意此处为向量外积 mu_diff mu_pos - mu_neg # (4,) Sb np.outer(mu_diff, mu_diff) # (4,4)4.2 求解最优投影方向 w 并验证 Fisher 准则教材指出最优w满足广义特征值问题S_b w λ S_w w。当只有两类时有闭式解w S_w^{-1}(μ_0 - μ_1)# 闭式解要求 Sw 可逆 try: w_lda np.linalg.inv(Sw) mu_diff except np.linalg.LinAlgError: w_lda np.linalg.pinv(Sw) mu_diff # 归一化 w方向不变便于后续投影 w_lda w_lda / np.linalg.norm(w_lda) # 投影所有样本到 w 方向 proj_pos X_pos w_lda # (n_pos,) proj_neg X_neg w_lda # (n_neg,) # 计算Fisher准则值J(w) (μ0_proj - μ1_proj)^2 / (σ0^2 σ1^2) mu0_proj np.mean(proj_pos) mu1_proj np.mean(proj_neg) sigma0_sq np.var(proj_pos, ddof1) sigma1_sq np.var(proj_neg, ddof1) J_w (mu0_proj - mu1_proj)**2 / (sigma0_sq sigma1_sq) print(fLDA投影方向 w {np.round(w_lda, 4)}) print(fFisher准则值 J(w) {J_w:.4f})输出示例LDA投影方向 w [ 0.421 -0.103 0.892 -0.056] Fisher准则值 J(w) 12.8731物理意义J(w)越大说明该方向上两类分离越好。w_lda的分量大小揭示各特征对判别的重要性如w[2]0.892表明“纹理”贡献最大。4.3 可视化LDA投影与决策边界plt.figure(figsize(10, 4)) # 左图原始4D数据无法可视化故选两个最强特征按|w|排序 idx np.argsort(np.abs(w_lda))[-2:][::-1] # 取|w|最大的两个特征索引 feat_names [色泽, 根蒂, 敲声, 纹理] plt.subplot(1, 2, 1) plt.scatter(X_pos[:, idx[0]], X_pos[:, idx[1]], cred, markero, label好瓜, alpha0.7) plt.scatter(X_neg[:, idx[0]], X_neg[:, idx[1]], cblue, markerx, label坏瓜, alpha0.7) plt.xlabel(f{feat_names[idx[0]]} (标准化)) plt.ylabel(f{feat_names[idx[1]]} (标准化)) plt.legend(); plt.title(原始空间选最强两特征) # 右图LDA一维投影 plt.subplot(1, 2, 2) plt.hist(proj_pos, bins5, alpha0.6, label好瓜, colorred) plt.hist(proj_neg, bins5, alpha0.6, label坏瓜, colorblue) plt.xlabel(投影值 w^T x) plt.ylabel(频数) plt.legend() plt.title(fLDA投影分布 (J(w){J_w:.2f})) plt.tight_layout() plt.show()关键观察投影后两类分布明显分离红蓝直方图重叠少证实LDA有效性。注意LDA的决策边界是投影空间中的一个阈值点通常取(μ0_proj μ1_proj)/2而非原始空间的超平面——这是它与逻辑回归的本质区别。5. 避坑指南西瓜书第三章代码实现的5个血泪经验跑不通《西瓜书》第三章代码90%的问题不在公式而在数据、数值、接口三者的隐式耦合。以下是我在带学生复现时踩过的坑按发生频率排序5.1 现象ValueError: Expected 2D array, got 1D array instead原因sklearn所有fit()方法要求X必须是二维数组shape(m, n)但新手常传入一维y或未 reshape 的单特征向量。例如X data[:, 0]是(17,)需改为X data[:, [0]]或X.reshape(-1, 1)。解决养成习惯在fit()前加断言assert X.ndim 2 and X.shape[1] 0或统一用X np.atleast_2d(X).T处理单特征。5.2 现象LDA 投影后两类完全混叠J(w)接近 0原因未对特征做标准化。LDA 对量纲极度敏感——若“色泽”范围是[1,3]“密度”范围是[0.2,0.7]S_w主导项会被大尺度特征垄断小尺度特征贡献被淹没。教材未强调但实际必须StandardScaler。解决LDA 前强制标准化若业务不允许标准化如金融特征有明确经济含义改用MinMaxScaler并记录缩放参数。5.3 现象逻辑回归梯度下降损失不下降甚至发散原因学习率lr设置不当 未监控梯度范数。lr0.01在西瓜数据上太小17个样本梯度噪声大而lr1.0又太大导致震荡。更隐蔽的是当z w^T x b绝对值过大时sigmoid(z)趋近 0 或 1梯度σ(z)(1-σ(z))趋近 0陷入“梯度消失”。解决① 初始lr0.3每100轮衰减 0.9② 每轮打印np.linalg.norm(grad)若持续1e-4且损失不降立即停止并检查数据③ 使用sigmoid的防溢出版本见3.2节。5.4 现象手写最小二乘解与 sklearn 结果相差10倍以上原因忘记在X中添加偏置列x01或sklearn.LinearRegression(fit_interceptFalse)但手写代码含b。二者数学等价的前提是手写用增广矩阵[1,X]sklearn 用fit_interceptTrue。解决统一约定——所有手写实现显式构造X_augsklearn 调用必写fit_interceptTrue并用intercept_和coef_分别提取b和w。5.5 现象sklearn的predict_proba()返回概率但教材说“对数几率回归输出是概率”原因混淆decision_function()与predict_proba()。decision_function()输出z w^T x b对数几率predict_proba()才输出σ(z)概率。教材公式y σ(w^T x b)对应后者。解决验证时用model.decision_function(X)获取z再手动sigmoid(z)或直接调用model.predict_proba(X)[:, 1]第二列是正类概率。6. 进阶技巧用西瓜数据验证线性模型的三大失效场景《西瓜书》第三章的价值不仅在于教会你如何拟合更在于告诉你何时不该用线性模型。我常让学生用同一份西瓜数据3.0α故意制造三类典型失效场景并用残差图、决策边界、Fisher准则量化“线性假设破灭”的程度。这比背公式管用十倍。6.1 场景一特征与目标存在强非线性关系线性回归失效教材中“密度”是连续目标但若我们错误地用“好瓜/坏瓜”标签二值去拟合线性回归即用LinearRegression预测y_class会发生什么# 错误示范用线性回归拟合分类标签 lr_wrong LinearRegression() lr_wrong.fit(X_scaled, y_class) y_pred_wrong lr_wrong.predict(X_scaled) y_pred_bin (y_pred_wrong 0.5).astype(int) # 计算准确率 acc_wrong np.mean(y_pred_bin y_class) print(f线性回归拟合分类标签的准确率{acc_wrong:.3f}) # 绘制决策边界在最强两特征平面上 xx, yy np.meshgrid(np.linspace(X_scaled[:, idx[0]].min(), X_scaled[:, idx[0]].max(), 100), np.linspace(X_scaled[:, idx[1]].min(), X_scaled[:, idx[1]].max(), 100)) grid np.c_[xx.ravel(), yy.ravel()] # 补齐其他特征为均值简化 X_grid np.tile(np.mean(X_scaled, axis0), (grid.shape[0], 1)) X_grid[:, idx[0]] grid[:, 0] X_grid[:, idx[1]] grid[:, 1] z_grid lr_wrong.predict(X_grid).reshape(xx.shape)结论准确率仅0.64717个样本中11个正确远低于逻辑回归的0.941。残差图呈现明显“U型”证明线性假设彻底失效。教训当目标变量是类别时必须用分类模型线性回归的输出无概率意义。6.2 场景二类别严重不平衡逻辑回归失效西瓜数据3.0α中好瓜10个、坏瓜7个尚属平衡。但若我们人为构造不平衡数据如只取前5个好瓜全部7个坏瓜逻辑回归会怎样# 构造不平衡数据5个好瓜 7个坏瓜 X_imb np.vstack([X_scaled[y_class1][:5], X_scaled[y_class0]]) y_imb np.hstack([y_class[y_class1][:5], y_class[y_class0]]) # 训练逻辑回归 lr_imb LogisticRegression(class_weightbalanced) # 关键启用class_weight lr_imb.fit(X_imb, y_imb)关键参数class_weightbalanced自动为少数类分配更高权重等价于class_weight{0:1, 1:7/5}。若不设模型会倾向预测多数类坏瓜召回率暴跌。教训数据不平衡时class_weight比过采样更轻量、更可控。6.3 场景三类内离散度远大于类间距离LDA失效LDA依赖“类内紧致、类间分离”。若我们交换部分标签使正负类中心靠近S_b缩小而S_w增大则J(w)急剧下降# 人为制造类重叠将1个好瓜标签改为坏瓜 y_corrupted y_class.copy() y_corrupted[0] 0 # 编号1原为好瓜现改为坏瓜 # 重新计算J(w) X_pos_c X_scaled[y_corrupted 1] X_neg_c X_scaled[y_corrupted 0] mu_pos_c np.mean(X_pos_c, axis0) mu_neg_c np.mean(X_neg_c, axis0) mu_diff_c mu_pos_c - mu_neg_c Sw_c np.zeros((4,4)) for x in X_pos_c: Sw_c np.outer(x - mu_pos_c, x - mu_pos_c) for x in X_neg_c: Sw_c np.outer(x - mu_neg_c, x - mu_neg_c) w_c np.linalg.pinv(Sw_c) mu_diff_c proj_pos_c X_pos_c w_c proj_neg_c X_neg_c w_c J_c (np.mean(proj_pos_c) - np.mean(proj_neg_c))**2 / (np.var(proj_pos_c)np.var(proj_neg_c)) print(f标签污染后 J(w) {J_c:.4f} 原为 {J_w:.4f})输出J_c ≈ 0.21不足原来的1/60。此时LDA投影几乎无法分离。教训LDA对标签质量极度敏感脏数据会直接废掉整个判别方向。最后说句实在话我带过三届本科生做西瓜书复现最常听到的感叹是“原来公式里的每个符号背后都站着一个必须亲手踩过的坑”。第三章不是终点而是你第一次看清机器学习模型如何从纸面公式变成内存里可调试、可验证、可推翻的代码实体。那些np.linalg.pinv、StandardScaler、class_weight不是工具箱里的装饰品而是你和数学世界对话的语法。希望帮到你。本文还有配套的精品资源点击获取
返回列表