)
从零实现线性回归用Python揭秘房价预测的数学之美在机器学习的世界里线性回归就像是一把打开预测大门的金钥匙。它不仅是最基础的算法之一更是理解更复杂模型的必经之路。想象一下你能够通过几个简单的数学公式就能预测一座房子的价格——这正是线性回归的魅力所在。本文将带你从数学原理出发手把手实现一个完整的线性回归模型并用它来预测波士顿房价。不同于简单的调用sklearn我们会深入算法内部用NumPy一步步构建模型让你真正理解机器学习的黑箱里发生了什么。1. 线性回归的数学基础线性回归的核心思想非常简单找到一条直线在更高维度是超平面使得这条直线能够最好地拟合数据点。用数学语言来说就是找到一组参数θ使得预测值ŷ与真实值y之间的差距最小。1.1 模型公式解析多元线性回归的预测公式可以表示为ŷ θ₀ θ₁x₁ θ₂x₂ ... θₙxₙ其中ŷ是预测值θ₀是偏置项截距θ₁到θₙ是各个特征的权重系数x₁到xₙ是输入特征为了简化计算我们通常会添加一个全为1的特征x₀这样公式可以写成向量形式ŷ hθ(x) θᵀx1.2 损失函数与优化目标为了衡量模型的好坏我们需要定义一个损失函数。在线性回归中最常用的是均方误差(MSE)def mse(y_true, y_pred): return np.mean((y_true - y_pred)**2)我们的目标就是找到一组θ使得MSE最小化。这个优化问题有解析解称为正规方程θ (XᵀX)⁻¹Xᵀy其中X是特征矩阵包含x₀1的列y是目标值向量2. 数据准备与预处理2.1 加载波士顿房价数据集波士顿房价数据集是机器学习中的经典数据集包含506个样本每个样本有13个特征和1个目标值房价中位数。from sklearn.datasets import load_boston boston load_boston() X boston.data y boston.target2.2 数据标准化为了确保不同尺度的特征对模型有相同的影响我们需要对数据进行标准化处理from sklearn.preprocessing import StandardScaler scaler StandardScaler() X_scaled scaler.fit_transform(X)2.3 添加偏置项记得我们在数学推导中添加的x₀1吗现在需要在代码中实现# 添加全1列作为偏置项 X_b np.c_[np.ones((len(X_scaled), 1)), X_scaled]3. 实现线性回归3.1 正规方程解法根据前面的数学推导我们可以直接实现正规方程def linear_regression(X, y): # 计算θ (XᵀX)⁻¹Xᵀy theta np.linalg.inv(X.T.dot(X)).dot(X.T).dot(y) return theta3.2 处理矩阵不可逆情况在实际应用中XᵀX可能不是满秩矩阵导致无法求逆。我们可以使用伪逆来避免这个问题theta np.linalg.pinv(X.T.dot(X)).dot(X.T).dot(y)3.3 完整实现将上述步骤整合成一个完整的类class LinearRegression: def __init__(self): self.theta None def fit(self, X, y): X_b np.c_[np.ones((len(X), 1)), X] self.theta np.linalg.pinv(X_b.T.dot(X_b)).dot(X_b.T).dot(y) def predict(self, X): X_b np.c_[np.ones((len(X), 1)), X] return X_b.dot(self.theta)4. 模型评估与优化4.1 常用评估指标除了MSE我们还可以计算R²分数来评估模型性能def r2_score(y_true, y_pred): ss_res np.sum((y_true - y_pred)**2) ss_tot np.sum((y_true - np.mean(y_true))**2) return 1 - (ss_res / ss_tot)4.2 交叉验证为了更可靠地评估模型我们应该使用交叉验证from sklearn.model_selection import cross_val_score lr LinearRegression() mse_scores -cross_val_score(lr, X_scaled, y, scoringneg_mean_squared_error, cv5) r2_scores cross_val_score(lr, X_scaled, y, scoringr2, cv5)4.3 特征重要性分析线性回归的一个优势是我们可以直接查看各个特征的权重feature_importance pd.DataFrame({ Feature: [Bias] list(boston.feature_names), Weight: lr.theta }).sort_values(Weight, keyabs, ascendingFalse)5. 实际应用中的注意事项5.1 特征相关性检查在应用线性回归前应该检查特征之间的相关性。高度相关的特征会导致矩阵接近奇异影响模型稳定性。corr_matrix pd.DataFrame(X_scaled).corr()5.2 异常值处理线性回归对异常值敏感可以考虑使用以下方法可视化检测箱线图、散点图使用RobustScaler代替StandardScaler考虑岭回归或Lasso回归5.3 模型假设验证线性回归有几个重要假设需要验证线性关系特征和目标之间确实存在线性关系同方差性残差的方差应该是常数正态性残差应该近似正态分布可以通过残差图来验证这些假设y_pred lr.predict(X_scaled) residuals y - y_pred plt.scatter(y_pred, residuals) plt.axhline(y0, colorr, linestyle-) plt.xlabel(Predicted values) plt.ylabel(Residuals)6. 扩展与进阶6.1 梯度下降实现除了正规方程我们还可以用梯度下降来求解线性回归def gradient_descent(X, y, learning_rate0.01, n_iters1000): m len(y) theta np.zeros(X.shape[1]) for _ in range(n_iters): gradients 2/m * X.T.dot(X.dot(theta) - y) theta - learning_rate * gradients return theta6.2 正则化方法为了防止过拟合可以引入L2正则化岭回归或L1正则化Lasso回归# 岭回归 theta_ridge np.linalg.inv(X.T.dot(X) alpha*np.eye(X.shape[1])).dot(X.T).dot(y)6.3 多项式回归当数据不是线性关系时可以通过添加多项式特征来增强模型表现from sklearn.preprocessing import PolynomialFeatures poly PolynomialFeatures(degree2) X_poly poly.fit_transform(X_scaled)7. 完整代码示例以下是整合了所有功能的完整实现import numpy as np from sklearn.datasets import load_boston from sklearn.preprocessing import StandardScaler from sklearn.model_selection import train_test_split # 加载数据 boston load_boston() X, y boston.data, boston.target # 数据预处理 scaler StandardScaler() X_scaled scaler.fit_transform(X) X_train, X_test, y_train, y_test train_test_split(X_scaled, y, test_size0.2) # 线性回归实现 class LinearRegression: def __init__(self): self.theta None def fit(self, X, y): X_b np.c_[np.ones((len(X), 1)), X] self.theta np.linalg.pinv(X_b.T.dot(X_b)).dot(X_b.T).dot(y) def predict(self, X): X_b np.c_[np.ones((len(X), 1)), X] return X_b.dot(self.theta) def mse(self, X, y): y_pred self.predict(X) return np.mean((y - y_pred)**2) def r2(self, X, y): y_pred self.predict(X) ss_res np.sum((y - y_pred)**2) ss_tot np.sum((y - np.mean(y))**2) return 1 - (ss_res / ss_tot) # 训练和评估模型 lr LinearRegression() lr.fit(X_train, y_train) print(fMSE on test set: {lr.mse(X_test, y_test):.2f}) print(fR² on test set: {lr.r2(X_test, y_test):.2f})