在 Python 中实现线性回归模型有多种方式，可分为**手动推导（基于最小二乘法的解析解）**、**使用 `scikit-learn`（工业级封装）** 和**使用 `statsmodels`（侧重统计推断）**。下面我将从原理出发，逐步给出三种主流实现方式，并附详细解释与对比。 --- ### ✅ 1. 手动实现（解析解：正规方程 Normal Equation） 线性回归的目标是最小化残差平方和（RSS）： \[ \min_{\mathbf{w}} \|\mathbf{X}\mathbf{w} - \mathbf{y}\|^2 \] 当设计矩阵 $\mathbf{X} \in \mathbb{R}^{n \times (d+1)}$（含偏置列）满秩时，闭式解为： \[ \mathbf{w} = (\mathbf{X}^\top \mathbf{X})^{-1} \mathbf{X}^\top \mathbf{y} \] > ⚠️ 注意：需添加全1列（截距项），或用 `np.column_stack([np.ones(n), X])` ```python import numpy as np import matplotlib.pyplot as plt # 生成模拟数据：y = 2x + 1 + noise np.random.seed(42) X = np.random.randn(100, 1) * 2 y = 2 * X.flatten() + 1 + np.random.randn(100) * 0.5 # 添加偏置项（截距列） X_with_bias = np.column_stack([np.ones(X.shape[0]), X]) # 正规方程求解：w = (X^T X)^{-1} X^T y w = np.linalg.inv(X_with_bias.T @ X_with_bias) @ X_with_bias.T @ y # 提取参数 intercept, slope = w[0], w[1] print(f"手动实现结果 — 截距: {intercept:.4f}, 斜率: {slope:.4f}") # 可视化 plt.scatter(X, y, alpha=0.6, label='Data') plt.plot(X, slope * X + intercept, 'r-', label=f'Fit: y = {slope:.2f}x + {intercept:.2f}') plt.legend() plt.title("Linear Regression via Normal Equation") plt.show() ``` ✅ **优点**：原理清晰、无迭代、适合中小规模数据 ❌ **缺点**：矩阵求逆计算复杂度 $O(d^3)$，对病态矩阵敏感（可用伪逆 `np.linalg.pinv` 改进） --- ### ✅ 2. 使用 `scikit-learn`（推荐用于工程实践） `sklearn.linear_model.LinearRegression` 默认使用 **SVD 分解**（比直接求逆更稳定），支持正则化（如 `Ridge`）、多输出、标准化等。 ```python from sklearn.linear_model import LinearRegression from sklearn.metrics import mean_squared_error, r2_score # 注意：sklearn 的 LinearRegression 默认包含截距，无需手动加 bias 列 model = LinearRegression(fit_intercept=True) model.fit(X, y) # 自动处理 X 为 2D array（n_samples × n_features） # 获取参数 print(f"sklearn 结果 — 截距: {model.intercept_:.4f}, 系数: {model.coef_[0]:.4f}") # 预测与评估 y_pred = model.predict(X) print(f"MSE: {mean_squared_error(y, y_pred):.4f}") print(f"R²: {r2_score(y, y_pred):.4f}") ``` ✅ **优点**：API 统一、高效稳定、支持 pipeline、交叉验证、扩展性强（Lasso/Ridge/ElasticNet） ❌ **缺点**：黑盒程度略高（但源码开源可查） --- ### ✅ 3. 使用 `statsmodels`（适合统计分析与推断） 提供详细的统计摘要（p 值、置信区间、F 检验、残差诊断等），类似 R 的 `lm()`。 ```python import statsmodels.api as sm # 必须显式添加常数项（截距） X_sm = sm.add_constant(X) # 等价于 np.column_stack([np.ones(n), X]) model_sm = sm.OLS(y, X_sm).fit() print(model_sm.summary()) # 输出完整统计报告（含 t-stat, p-value, R²_adj, etc.） ``` 📌 输出中关键字段： - `coef`：系数估计值 - `P>|t|`：系数是否显著非零（<0.05 通常认为显著） - `R-squared`：决定系数 - `Prob(F-statistic)`：整体模型显著性检验 --- ### 🔁 补充：梯度下降实现（理解优化过程） 适用于大规模数据或后续扩展（如逻辑回归、神经网络）： ```python def linear_regression_gd(X, y, lr=0.01, epochs=1000): n, d = X.shape w = np.random.randn(d + 1) * 0.01 X_bias = np.column_stack([np.ones(n), X]) for _ in range(epochs): y_pred = X_bias @ w error = y_pred - y grad = (2 / n) * X_bias.T @ error w -= lr * grad return w[0], w[1:] # intercept, coefficients intercept_gd, coef_gd = linear_regression_gd(X, y) print(f"梯度下降结果 — 截距: {intercept_gd:.4f}, 系数: {coef_gd[0]:.4f}") ``` --- ### 📌 总结对比表 | 方法 | 适用场景 | 是否需归一化 | 是否提供统计推断 | 计算稳定性 | 学习价值 | |------------------|------------------------------|--------------|------------------|------------|----------| | 正规方程（手动） | 教学/小数据/理解原理 | 否 | ❌ | 中（病态易失败） | ⭐⭐⭐⭐⭐ | | `sklearn` | 工程部署/快速建模/生产环境 | 推荐（尤其含量纲差异） | ❌（但可配合 `sklearn` 生态做检验） | ⭐⭐⭐⭐⭐（SVD） | ⭐⭐⭐ | | `statsmodels` | 学术研究/假设检验/报告撰写 | 否 | ✅（完整统计量） | ⭐⭐⭐⭐ | ⭐⭐⭐⭐ | | 梯度下降 | 大数据/自定义损失/深度学习基础 | ✅强烈推荐 | ❌ | ⭐⭐（需调参） | ⭐⭐⭐⭐⭐ | ---