1. 从“差不多就行”到“最优解”：为什么我们需要最小二乘

做数据分析、搞机器学习，甚至只是用Excel拟合个趋势线，你大概率都听过“最小二乘法”这个名字。它太常见了，常见到很多人把它当成一个黑盒工具，数据往里一扔，拟合直线就出来了。但你真的理解它背后在做什么吗？为什么偏偏是“二乘”（平方），而不是“一乘”（绝对值）或者“四乘”？今天，我们不谈复杂的矩阵推导，就从最朴素的直觉出发，掰开揉碎了聊聊线性最小二乘。这可能是你见过最“人话”的解读，目标是让你不仅会用，更能理解它每一步的“所以然”，下次再遇到拟合问题，你能像个老手一样，知道它在干什么，甚至能预判它可能出什么幺蛾子。

想象一个经典场景：你有一堆散点数据，想找一条直线（y = ax + b）来代表它们的趋势。眼睛一看，大概能画一条。但“大概”不行，我们需要一个明确的、可计算的准则，来判断哪条直线是“最好”的。这个准则就是：让所有数据点到这条直线的垂直距离的平方和最小。为什么是平方和？简单说，就是为了避免正负距离相互抵消（距离有正有负），同时让大的误差被显著地惩罚（平方会放大误差）。这就是“最小二乘”最核心的思想：寻找一个最优的预测模型（这里是直线），使得预测值与真实值之间的总体差异最小。

这个方法几乎无处不在。从物理学中根据实验数据确定定律参数，到经济学中分析变量间的关系，再到你手机里传感器数据的校准，背后都有它的身影。它解决的是一个最基础的“估计”问题：我们观测到的世界充满噪声，如何从这些带噪声的数据中，提炼出那个最可能存在的、简洁的规律？线性最小二乘给出了一个在数学上非常漂亮、计算上相对简单的答案。

2. 核心思想拆解：误差、惩罚与优化目标

要彻底搞懂最小二乘，我们得先把它拆成几个基本概念，看看它们是怎么组合在一起的。

2.1 模型假设：我们相信世界是线性的

线性最小二乘的第一步，是做出一个假设：我们认为要拟合的数据背后，隐藏着一个线性关系。用公式表示就是：

y = β₀ + β₁*x + ε

这里，y是我们观测到的因变量（比如销售额），x是自变量（比如广告投入），β₀和β₁是我们要求解的未知参数（截距和斜率）。ε代表误差项，它包含了所有未被模型考虑的因素（比如随机波动、测量误差等）。我们的目标就是找到一对β₀和β₁，使得根据x计算出的预测值ŷ = β₀ + β₁*x，与真实的y值尽可能接近。

这个假设是模型的基石。如果真实世界的关系根本不是线性的（比如是指数增长），那你用直线去硬套，结果肯定会跑偏。这就好比用一把直尺去量一个球体的周长，工具本身就用错了。所以，在应用最小二乘前，画个散点图看看数据的大致形态，是一个绝对不能省的好习惯。

2.2 损失函数：如何定义“好坏”

有了模型，我们需要一个量化的标准来衡量这条直线拟合得好不好。这个标准就是损失函数（也叫成本函数）。对于最小二乘法，这个函数定义为所有样本点的残差平方和：

RSS(β₀, β₁) = Σ(y_i - ŷ_i)² = Σ(y_i - (β₀ + β₁*x_i))²

RSS就是我们要最小化的那个东西。残差e_i = y_i - ŷ_i，是真实值减去预测值，可以直观理解为每个数据点垂直偏离直线的距离。

注意：这里选择“平方”作为损失，有深刻的数学和实用考量。第一，平方项处处可导，光滑连续，这为我们后续使用求导这种强有力的数学工具找到了最优解铺平了道路。如果使用绝对值（一乘），在零点处不可导，求解会麻烦得多。第二，平方项对大的误差给予更大的惩罚，这符合我们的直觉：一个偏离很远的“坏点”应该比一堆稍微偏离的“好点”更严重地影响模型。第三，在误差服从正态分布的假设下，最小二乘解恰好是模型参数的最大似然估计，这在统计上是最优的。

2.3 求解过程：从几何与微积分两个视角

怎么找到使RSS最小的β₀和β₁呢？有两种理解方式，它们殊途同归。

视角一：微积分求极值。这是一个纯粹的优化问题。我们对损失函数RSS分别关于β₀和β₁求偏导数，并令其等于零。这会得到两个方程，称为正规方程。

∂RSS/∂β₀ = -2 * Σ(y_i - β₀ - β₁x_i) = 0 ∂RSS/∂β₁ = -2 * Σ[x_i * (y_i - β₀ - β₁x_i)] = 0

解这个二元一次方程组，就能得到著名的解公式：

β₁ = Σ[(x_i - x̄)(y_i - ȳ)] / Σ[(x_i - x̄)²] β₀ = ȳ - β₁ * x̄

其中x̄和ȳ分别是x和y的样本均值。这个公式清晰展示了斜率β₁的本质：是x和y的协方差除以x的方差。当x变化时，y随之变化的比例。

视角二：几何投影。这是一个更优美的视角。我们把所有y观测值看作一个高维空间中的向量，把由x的常数项和一次项张成的空间看作一个平面（在简单线性回归中就是二维平面）。寻找最优的ŷ（即β₀ + β₁*x），本质上是在这个平面上寻找一个点，使得它到真实y向量的距离最短。根据几何知识，这个最短距离是通过向平面作垂直投影实现的。而“距离最短”正好等价于“残差平方和最小”。正规方程(XᵀX)β = Xᵀy（这是多元情况下的矩阵形式）的本质，就是确保残差向量e = y - Xβ与X张成的列空间垂直，即Xᵀe = 0。

我个人更喜欢几何视角，因为它把抽象的优化问题变成了空间中看得见的投影操作，理解起来非常直观，而且能自然推广到多元线性回归乃至更复杂的模型。

3. 手算演示与代码实现：从公式到结果

理论说得再多，不如亲手算一遍。我们用一个超级简单的数据集来演示全过程。

假设我们有三个数据点：(1, 1)，(2, 2)，(3, 2)。我们的目标是拟合一条直线y = β₀ + β₁*x。

第一步：计算必要的中间量首先计算均值： x̄ = (1+2+3)/3 = 2 ȳ = (1+2+2)/3 ≈ 1.6667

然后计算分子和分母： 分子 Σ[(x_i - x̄)(y_i - ȳ)] = (1-2)(1-1.6667) + (2-2)(2-1.6667) + (3-2)(2-1.6667) = (-1)(-0.6667) + 0*(0.3333) + 1*(0.3333) = 0.6667 + 0.3333 = 1.0 分母 Σ[(x_i - x̄)²] = (1-2)² + (2-2)² + (3-2)² = 1 + 0 + 1 = 2

第二步：套用公式求解参数β₁ = 分子 / 分母 = 1.0 / 2 = 0.5 β₀ = ȳ - β₁ * x̄ = 1.6667 - 0.5 * 2 = 0.6667

所以，我们得到的最佳拟合直线是：y = 0.6667 + 0.5*x

第三步：理解结果斜率β₁=0.5意味着，x每增加1个单位，我们预测y会增加0.5个单位。截距β₀≈0.667是当x=0时y的预测值。你可以把三个点的x值代进去，得到预测值ŷ分别为 1.1667， 1.6667， 2.1667。计算残差并验证，它们的平方和确实是最小的。

当然，现实中我们不会手算。下面用Python的numpy和sklearn分别实现一下，看看工业级工具怎么用。

import numpy as np from sklearn.linear_model import LinearRegression # 数据 x = np.array([1, 2, 3]).reshape(-1, 1) # 转成二维，sklearn要求 y = np.array([1, 2, 2]) # 方法1：使用正规方程公式 (与手算一致) X_design = np.c_[np.ones(3), x] # 添加一列1，对应截距项β₀ beta = np.linalg.inv(X_design.T @ X_design) @ X_design.T @ y print(f"正规方程求解 - 截距 β₀: {beta[0]:.4f}, 斜率 β₁: {beta[1]:.4f}") # 方法2：使用numpy的polyfit（多项式拟合，1阶就是直线） beta_np = np.polyfit(x.flatten(), y, deg=1) print(f"np.polyfit求解 - 斜率 β₁: {beta_np[0]:.4f}, 截距 β₀: {beta_np[1]:.4f}") # 方法3：使用scikit-learn model = LinearRegression(fit_intercept=True) # 默认拟合截距 model.fit(x, y) print(f"sklearn求解 - 截距 β₀: {model.intercept_:.4f}, 斜率 β₁: {model.coef_[0]:.4f}") # 计算预测值和残差平方和(RSS) y_pred = model.predict(x) rss = np.sum((y - y_pred) ** 2) print(f"残差平方和 RSS: {rss:.4f}")

运行这段代码，三种方法会得到完全相同的结果。sklearn的方式最常用也最规范，它背后默认求解的就是最小二乘问题。这里有个细节：np.polyfit返回的顺序是[斜率， 截距]，而sklearn的coef_属性存放斜率，intercept_属性存放截距，不要搞混了。

4. 多元线性回归：当世界不止一个影响因素

现实问题中，影响结果y的因素通常不止一个。比如房价，可能同时受面积、房龄、地段、楼层等多个因素影响。这时，我们就需要将简单线性回归推广到多元线性回归。

模型变成了：y = β₀ + β₁*x₁ + β₂*x₂ + ... + β_p*x_p + ε

这里x₁, x₂, ..., x_p是p个自变量（特征）。我们的目标同样是找到一组参数β₀, β₁, ..., β_p，使得残差平方和最小。此时，使用矩阵表示会异常简洁。

令：

• y是一个n×1的列向量（n个样本的观测值）。
• X是一个n×(p+1)的矩阵，第一列全是1（对应截距β₀），后面p列是各个特征的值。
• β是一个(p+1)×1的列向量，包含所有待求参数[β₀, β₁, ..., β_p]ᵀ。

那么模型可以写成：y = Xβ + ε损失函数RSS(β) = (y - Xβ)ᵀ(y - Xβ)通过求导，得到正规方程的矩阵形式：(XᵀX)β = Xᵀy如果XᵀX可逆（即X列满秩，没有完全共线性的特征），则最优解为：β = (XᵀX)⁻¹Xᵀy

这个公式是线性最小二乘理论的基石。它把复杂的多变量优化问题，转化成了一个清晰的矩阵运算。

# 多元线性回归示例 import pandas as pd from sklearn.datasets import make_regression from sklearn.model_selection import train_test_split # 生成一个模拟数据集：100个样本，3个有效特征 X, y = make_regression(n_samples=100, n_features=3, noise=10, random_state=42) X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42) model_multi = LinearRegression() model_multi.fit(X_train, y_train) print("模型参数（截距 + 三个特征的系数）：") print(f"截距 β₀: {model_multi.intercept_:.4f}") for i, coef in enumerate(model_multi.coef_): print(f"特征 x_{i+1} 的系数 β_{i+1}: {coef:.4f}") # 评估 train_score = model_multi.score(X_train, y_train) test_score = model_multi.score(X_test, y_test) print(f"\n训练集 R² 分数: {train_score:.4f}") print(f"测试集 R² 分数: {test_score:.4f}")

在多元情况下，每个系数β_j的解释需要小心：它表示在其他所有特征保持不变的情况下，特征x_j每增加一个单位，y平均变化β_j个单位。这就是“控制变量”的思想。

5. 理想与现实的差距：最小二乘的“七寸”与应对策略

最小二乘看起来完美，但它建立在几个重要的假设之上（高斯-马尔可夫假设）。一旦现实数据违背这些假设，盲目使用就会出问题。理解这些“七寸”，是你从“会用”到“懂用”的关键。

5.1 异常值：一颗老鼠屎坏了一锅粥

最小二乘的损失函数对误差进行平方，这意味着异常值（Outliers）会对模型产生巨大的、不成比例的影响。一个偏离很远的点，其残差的平方会非常大，模型会为了“讨好”这个点而严重扭曲，牺牲其他大多数正常点的拟合效果。

如何诊断？

• 可视化：绘制预测值-残差图。如果残差图中存在个别点远离零点聚集区，可能就是异常值。
• 统计量：计算标准化残差或学生化残差。绝对值大于2或3的点需要警惕。

如何应对？

1. 检查数据：首先确认异常值是否是数据录入错误或测量失误。如果是，直接修正或删除。
2. 鲁棒回归：如果异常值是真实但特殊的，可以考虑使用更鲁棒的回归方法，如Huber回归或RANSAC。Huber回归在误差较小时使用平方损失，误差较大时切换为线性损失，从而减弱异常值的影响。RANSAC则通过迭代随机采样来拟合一个受内点（正常点）支持的模型。
3. 变换变量：有时对因变量y进行变换（如取对数）可以压缩数据的尺度，减弱异常值的影响。

# 对比普通最小二乘与Huber回归对异常值的鲁棒性 import numpy as np import matplotlib.pyplot as plt from sklearn.linear_model import LinearRegression, HuberRegressor # 生成干净数据 np.random.seed(0) X = np.random.randn(100, 1) y = 2 * X.squeeze() + 1 + np.random.randn(100) * 0.5 # 添加一个异常值 X_outlier = np.append(X, [[5]]).reshape(-1, 1) y_outlier = np.append(y, [-10]) # 拟合模型 lr = LinearRegression().fit(X_outlier, y_outlier) huber = HuberRegressor(epsilon=1.35).fit(X_outlier, y_outlier) # epsilon是Huber损失从平方转为线性的阈值 # 绘图对比 plt.scatter(X_outlier, y_outlier, label='数据点') x_plot = np.linspace(X_outlier.min(), X_outlier.max(), 100).reshape(-1,1) plt.plot(x_plot, lr.predict(x_plot), 'r-', label='OLS (受异常值影响大)') plt.plot(x_plot, huber.predict(x_plot), 'g--', linewidth=2, label='Huber回归 (更鲁棒)') plt.legend() plt.xlabel('X') plt.ylabel('y') plt.title('异常值对OLS和Huber回归的影响对比') plt.show()

运行这段代码，你会清晰看到，那个在右下角的异常点，如何把OLS拟合的直线（红线）狠狠地往下拉，而Huber回归的线（绿虚线）则基本保持了原有趋势。

5.2 多重共线性：当特征们“抱团取暖”

多重共线性是指模型中的自变量之间存在高度相关关系。比如在预测房价时，同时使用了“房屋面积”和“房间数量”，这两个变量通常是相关的。这会导致什么问题？

1. 系数估计不稳定：(XᵀX)矩阵接近奇异，其逆矩阵变得数值上不稳定，微小的数据变动可能导致系数估计值发生巨大变化，甚至符号相反。
2. 系数难以解释：由于特征共享信息，很难区分每个特征对y的独立贡献。一个特征的系数可能因为另一个高度相关特征的存在而被低估或变得不显著。

如何诊断？

• 方差膨胀因子：计算每个特征的VIF。VIF = 1 / (1 - R²)，其中R²是将该特征对其他所有特征回归得到的决定系数。通常，VIF > 5 或 10 就表明存在较强的共线性。
• 条件数：计算XᵀX矩阵的条件数。条件数很大（比如 > 30）也提示存在多重共线性。

如何应对？

1. 特征选择：剔除冗余特征。可以使用领域知识，或者利用正则化方法（如Lasso）自动进行特征选择。
2. 主成分回归：使用主成分分析将相关的原始特征转换为一组不相关的主成分，然后用主成分作为新特征进行回归。
3. 岭回归：在损失函数中加入系数的L2范数惩罚项λΣβ_j²。这相当于对系数进行收缩，虽然会引入一点偏差，但能显著降低方差，稳定估计。λ是超参数，需要通过交叉验证选择。

from statsmodels.stats.outliers_influence import variance_inflation_factor import statsmodels.api as sm import pandas as pd # 假设我们有一个包含共线性特征的数据框 df # 这里用模拟数据演示 np.random.seed(42) n = 100 x1 = np.random.randn(n) * 10 x2 = x1 + np.random.randn(n) * 0.5 # x2 与 x1 高度相关 x3 = np.random.randn(n) * 5 # x3 独立 y_sim = 2*x1 + 3*x3 + np.random.randn(n)*2 df = pd.DataFrame({'x1': x1, 'x2': x2, 'x3': x3, 'y': y_sim}) X_for_vif = df[['x1', 'x2', 'x3']] X_for_vif = sm.add_constant(X_for_vif) # 添加常数项 # 计算VIF vif_data = pd.DataFrame() vif_data['feature'] = X_for_vif.columns vif_data['VIF'] = [variance_inflation_factor(X_for_vif.values, i) for i in range(X_for_vif.shape[1])] print(vif_data)

你会看到x1和x2的VIF值非常高（可能几百甚至上千），而x3的VIF接近1。这就是多重共线性的明确信号。

5.3 异方差性：误差的“音量”不恒定

最小二乘假设误差项ε的方差是常数（同方差性）。如果误差方差随着x的增大而增大（或减小），就称为异方差性。这在金融、经济数据中很常见（如高收入群体的消费波动更大）。

后果：虽然系数估计仍然是无偏的，但标准误的估计不再有效，从而导致假设检验（如t检验、F检验）和置信区间不可靠。

如何诊断？

• 残差图：绘制预测值ŷ与残差e的散点图。如果散点呈现明显的漏斗形、扇形或任何非随机的模式，则可能存在异方差。

如何应对？

1. 变量变换：对因变量y进行变换，如取对数log(y)、平方根sqrt(y)，常能稳定方差。
2. 加权最小二乘法：如果知道方差的变化规律，可以给每个样本点赋予不同的权重（方差小的点权重大，方差大的点权重小），然后最小化加权残差平方和。
3. 使用稳健标准误：在统计推断时，使用如“White异方差稳健标准误”等方法，即使存在异方差，也能得到有效的检验结果。

5.4 模型误设：用直线拟合曲线

这是最根本的问题。如果真实关系是非线性的，那么无论你怎么优化，线性模型都无法很好地捕捉数据的结构。残差图中会显示出明显的系统性模式（如U型或倒U型）。

应对策略：

• 添加多项式项或交互项：将x²,x³或x1*x2作为新特征加入模型，使其能够拟合更复杂的曲线关系。这本质上还是在用线性模型（对特征非线性，对参数仍是线性）去逼近。
• 使用非线性模型：如决策树、支持向量机（带核函数）、神经网络等。

6. 评估与诊断：你的模型真的靠谱吗？

拟合出模型后，不能只看一个R²就万事大吉。你需要一套组合拳来评估和诊断模型。

6.1 核心评估指标

1. R²（决定系数）：最常用的指标，表示模型能解释的因变量变异比例。R² = 1 - RSS/TSS，其中TSS是总平方和。R²越接近1越好。但要注意，增加无关变量总会让R²提高，因此需要调整。
2. 调整R²：考虑了自变量个数p和样本量n的惩罚，Adj-R² = 1 - [(1-R²)(n-1)/(n-p-1)]。用于比较不同变量数的模型，比R²更可靠。
3. 均方误差及其变体：
   • MSE（均方误差）：MSE = RSS / n。数值越小越好，但量纲是y的平方，解释性稍差。
   • RMSE（均方根误差）：RMSE = sqrt(MSE)。与y同量纲，更直观。比如房价预测RMSE为5万，意味着平均预测误差在5万左右。
   • MAE（平均绝对误差）：MAE = Σ|y_i - ŷ_i| / n。对异常值不如MSE敏感，解释更直接。

6.2 统计诊断：超越R²

1. F检验：检验模型整体是否显著。原假设是所有斜率系数都为0。如果p值很小（如<0.05），则拒绝原假设，认为至少有一个自变量与y存在线性关系。
2. t检验：检验单个系数是否显著不为0。原假设是β_j = 0。同样看p值。
3. 置信区间：给出系数估计值的一个范围，例如“我们有95%的信心认为，真实的斜率在[0.4， 0.6]之间”。区间越窄，估计越精确。

6.3 残差分析：模型的“体检报告”

这是模型诊断中最重要的一环。你需要绘制并检查以下残差图：

• 残差 vs. 拟合值图：检查同方差性和模型误设。理想情况是残差随机均匀分布在0线周围，无明显模式。
• 残差的正态概率图：检查残差是否近似正态分布。理想情况是点大致落在一条对角线上。严重的偏离会影响系数检验的有效性。
• 残差 vs. 自变量图：检查某个特定自变量与残差是否存在关系，这可能提示需要添加该自变量的高次项或交互项。

import matplotlib.pyplot as plt import statsmodels.api as sm from statsmodels.graphics.gofplots import qqplot # 使用之前的多元回归模型示例 model_sm = sm.OLS(y_train, sm.add_constant(X_train)).fit() print(model_sm.summary()) # 查看包含R²、F检验、t检验等的详细报告 # 绘制诊断图 fig = plt.figure(figsize=(12, 8)) # 1. 残差 vs. 拟合值 ax1 = fig.add_subplot(2, 2, 1) ax1.scatter(model_sm.fittedvalues, model_sm.resid) ax1.axhline(y=0, color='r', linestyle='--') ax1.set_xlabel('Fitted values') ax1.set_ylabel('Residuals') ax1.set_title('Residuals vs Fitted') # 2. 正态概率图 ax2 = fig.add_subplot(2, 2, 2) qqplot(model_sm.resid, line='s', ax=ax2) ax2.set_title('Normal Q-Q') # 3. 标准化残差平方根 vs. 拟合值 (检查同方差性) ax3 = fig.add_subplot(2, 2, 3) ax3.scatter(model_sm.fittedvalues, np.sqrt(np.abs(model_sm.get_influence().resid_studentized_internal))) ax3.set_xlabel('Fitted values') ax3.set_ylabel('Sqrt(|Standardized Residuals|)') ax3.set_title('Scale-Location') # 4. 残差 vs. 杠杆值 (检查高杠杆点) ax4 = fig.add_subplot(2, 2, 4) sm.graphics.influence_plot(model_sm, ax=ax4, criterion="cooks") ax4.set_title('Residuals vs Leverage') plt.tight_layout() plt.show()

运行这段诊断代码，你会得到四张图，它们共同构成了模型的“体检报告”。仔细阅读这份报告，比单纯看一个R²值要有用得多。

7. 从理论到实战：一个完整的端到端案例

让我们用一个接近真实的数据集，走一遍完整的线性回归建模流程。我们使用经典的波士顿房价数据集（虽然由于其伦理问题已不推荐用于实际研究，但作为教学示例仍很合适）。

目标：根据房屋的各种特征（如犯罪率、房间数、税率等）预测其价格。

import pandas as pd import numpy as np from sklearn.datasets import fetch_california_housing # 使用更现代的加州房价数据集 from sklearn.model_selection import train_test_split, cross_val_score from sklearn.linear_model import LinearRegression, Ridge from sklearn.preprocessing import StandardScaler from sklearn.metrics import mean_squared_error, r2_score import matplotlib.pyplot as plt import seaborn as sns # 1. 加载数据 housing = fetch_california_housing() df = pd.DataFrame(housing.data, columns=housing.feature_names) df['MedHouseVal'] = housing.target # 中位数房价，单位是十万美元 print("数据形状:", df.shape) print("\n前几行数据:") print(df.head()) print("\n数据描述:") print(df.describe()) # 2. 探索性数据分析 # 查看目标变量分布 plt.figure(figsize=(12, 4)) plt.subplot(1, 2, 1) sns.histplot(df['MedHouseVal'], kde=True) plt.title('Target Distribution') # 查看特征与目标的相关性 plt.subplot(1, 2, 2) corr_matrix = df.corr() sns.heatmap(corr_matrix, annot=True, fmt='.2f', cmap='coolwarm', center=0) plt.title('Feature Correlation Matrix') plt.tight_layout() plt.show() # 3. 数据准备 X = df.drop('MedHouseVal', axis=1) y = df['MedHouseVal'] X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.2, random_state=42) # 特征标准化：对于线性模型，特别是考虑正则化时，标准化很重要 scaler = StandardScaler() X_train_scaled = scaler.fit_transform(X_train) X_test_scaled = scaler.transform(X_test) # 4. 训练基准模型（普通最小二乘） lr_baseline = LinearRegression() lr_baseline.fit(X_train_scaled, y_train) # 5. 模型评估 y_pred_train = lr_baseline.predict(X_train_scaled) y_pred_test = lr_baseline.predict(X_test_scaled) print("\n=== 基准模型 (OLS) 性能 ===") print(f"训练集 R²: {r2_score(y_train, y_pred_train):.4f}") print(f"测试集 R²: {r2_score(y_test, y_pred_test):.4f}") print(f"训练集 RMSE: {np.sqrt(mean_squared_error(y_train, y_pred_train)):.4f}") print(f"测试集 RMSE: {np.sqrt(mean_squared_error(y_test, y_pred_test)):.4f}") # 6. 使用交叉验证评估稳定性 cv_scores = cross_val_score(lr_baseline, X_train_scaled, y_train, cv=5, scoring='r2') print(f"\n5折交叉验证 R² 得分: {cv_scores.mean():.4f} (+/- {cv_scores.std()*2:.4f})") # 7. 尝试带正则化的模型（岭回归）以应对可能的共线性 ridge = Ridge(alpha=1.0) # alpha是正则化强度 ridge.fit(X_train_scaled, y_train) y_pred_test_ridge = ridge.predict(X_test_scaled) print(f"\n=== 岭回归 (alpha=1.0) 测试集 R²: {r2_score(y_test, y_pred_test_ridge):.4f}") # 8. 特征重要性分析（基于系数绝对值） coef_df = pd.DataFrame({ 'feature': housing.feature_names, 'coef_ols': lr_baseline.coef_, 'coef_ridge': ridge.coef_ }) coef_df['abs_coef_ols'] = np.abs(coef_df['coef_ols']) coef_df = coef_df.sort_values('abs_coef_ols', ascending=False) print("\n特征重要性排序 (基于OLS系数绝对值):") print(coef_df[['feature', 'coef_ols', 'coef_ridge']]) # 9. 诊断：绘制预测值与真实值对比图 plt.figure(figsize=(6, 6)) plt.scatter(y_test, y_pred_test, alpha=0.5, label='OLS Predictions') plt.scatter(y_test, y_pred_test_ridge, alpha=0.5, label='Ridge Predictions', marker='x') # 绘制理想对角线 max_val = max(y_test.max(), y_pred_test.max()) min_val = min(y_test.min(), y_pred_test.min()) plt.plot([min_val, max_val], [min_val, max_val], 'r--', label='Perfect Prediction') plt.xlabel('True Values') plt.ylabel('Predicted Values') plt.title('True vs Predicted Values') plt.legend() plt.show()

通过这个完整案例，你实践了从数据加载、探索、预处理、建模、评估到诊断的全过程。你会发现，即使是一个简单的线性模型，在真实数据上也能取得不错的效果，但同时也暴露了一些问题（比如测试集性能略低于训练集，可能存在轻微过拟合或特征共线性）。这时，岭回归的引入就是一个自然的优化尝试。

线性最小二乘是一个强大而基础的起点。它教会我们的，不仅仅是如何拟合一条直线，更是一种从数据中寻找规律、量化关系、评估不确定性的系统化思维方式。理解了它的原理、假设、优势和局限，你就掌握了数据分析中最重要的一把钥匙。下次当你面对一堆散点图时，希望你能自信地说：我知道那条“最佳”直线从何而来，也知道它何时会“说谎”。