李宏毅老师的机器学习课,作业一做PM2.5预测,应该是很多人入门时接触的第一个“完整”机器学习项目。我说的完整,不是指调个库把结果跑出来,而是从读数据、清洗、特征构造、手写梯度下降,到最终提交结果的全流程。这个作业虽然模型简单,但五脏俱全,数据解析、标准化、训练验证、调参避坑一个都绕不开。这篇文章就把我做完整个项目的复盘过程写下来,包含数据集怎么解析、特征怎么构造、为什么作业要求手写线性回归、以及我实际踩过的各种坑。无论你是准备交作业,还是想通过一个最小闭环真正理解回归模型,这篇都能给你一条可以直接对照复现的路径。
1. 项目背景与任务拆解
1.1 作业一究竟要我们做什么
先把这个任务说清楚。李宏毅机器学习作业一的目标是:根据台湾丰原观测站过去9个小时的空气监测数据,预测下一个小时的PM2.5浓度。评价指标是RMSE,也就是均方根误差,公式是:
RMSE = sqrt( (1/N) * sum( (y_pred - y_true)^2 ) )
这个指标的特点是大误差会被放大,因为先平方再开方,所以如果你在某个污染尖峰时刻预测严重偏低,整体分数会很难看。这也是为什么后续做误差分析时,要特别关注极端值样本。
为什么这个作业适合入门?因为它数据量小、特征维度低、模型简单,但把完整流程走了一遍:数据清洗、特征工程、模型训练、评估、预测。很多人学机器学习理论时觉得梯度下降很难,但这个作业把梯度下降用在了一个非常直观的连续值预测问题上,你能够亲眼看到loss随着迭代下降、RMSE逐步变小,这种正反馈对建立信心非常重要。
需要注意,这里是回归任务,不是分类任务。我们在预测一个连续的数值,而不是判断“空气质量是好是坏”。很多人把回归和逻辑回归搞混,逻辑回归虽然名字里带回归,但实际做的是分类。这个作业要用的是经典的线性回归模型,任务就是拟合一条(或一个超平面)来最小化预测值与真实值之间的误差。
1.2 数据集的真实结构与坑点
很多初学者拿到train.csv后,直接当成普通表格丢进模型,结果一跑就出问题。原因在于这个数据集的存储格式不是常规的“一行一个样本”,而是一种非常容易踩坑的排布方式。
训练集train.csv总共4320行。表面上看,每一行是某个时刻的18项指标观测值,但真实的组织逻辑是:每18行才组成一天的数据。这18行的顺序是固定的,依次对应AMB_TEMP(温度)、CH4、CO、NMHC、NO、NO2、NOx、O3、PM10、PM2.5、RAINFALL(降雨量)、RH(相对湿度)、SO2、THC、WD_HR(风向小时值)、WIND_DIREC(风向)、WIND_SPEED(风速)、WS_HR(风速小时值)。每行有24列,对应一天当中24个小时的连续观测值。
换算一下:12个月 × 每月取前20天 × 每天18行 = 4320行,数据量对得上。月份信息并没有单独的一列,是通过行号间接编码的,每360行切换一个月。第一次处理时如果没意识到这种结构,直接把行号当成样本序号,就会导致特征全部错位,训练出来的模型形同虚设。
数据里还有一堆比较隐蔽的格式问题。比如RAINFALL这一列在很多天里是字符串“NR”,表示无降雨(No Rain),如果直接调用pd.to_numeric会报错。再比如每行第一列是一个编号字段,实际建模时没什么用,需要去掉。这些细节看着小,但每一项都会卡住新手相当长的时间。
2. 数据预处理与特征工程
2.1 读数据与清洗的第一个大坑
读取和清洗数据是整个项目最容易翻车的地方。我的建议是分三步走:先读入,再替换异常值,最后重构数据结构。
先用pandas把csv读进来。注意默认情况下“NR”会被识别成字符串列,RAINFALL那列就变成object类型了,这一步不会报错,但后面计算会出问题。所以读取之后的第一件事,是把所有非数值的内容统一替换掉,我这里直接用0填充“NR”,代表降雨量为0:
import pandas as pd import numpy as np df = pd.read_csv('train.csv') df = df.replace('NR', 0) df = df.astype(float)这里有个小细节,替换“NR”之后必须把整张表转成float类型,否则后续切片计算时pandas会因为你混合了字符串和数值列而给出各种奇奇怪怪的警告。转float之后,第一列ID其实已经失去了意义,可以直接删掉:
df = df.drop(columns=['id'])然后就是重头戏:把表格重构为“天”级别的三维数组。我的做法是先把DataFrame转成numpy矩阵,然后按18行一组切分。因为每个月的第n天对应18行,且这18行的顺序固定,所以可以直接用reshape:
data = df.to_numpy() # 4320行, 去掉ID后变成4320行24列 # 12个月 * 20天 = 240天, 每天18行, 每行24列 days = data.reshape(240, 18, 24)这样days的形状就是(240, 18, 24),含义是:第i天、第j种指标、第h个时刻的数值。之后所有特征构造都基于这个三维数组进行,会方便很多。
2.2 特征构造:两种思路,一个结果导向
特征构造是这个作业的灵魂所在。原始作业里有两版常见说法,一版是“使用前9个小时的PM2.5数据预测第10个小时”,另一版是“使用前9个小时所有18项指标预测第10个小时的PM2.5”。我自己的建议是,两个版本都做一遍,既能加深理解,也能直观体会特征信息量对模型效果的影响。
先看最简单版:只用PM2.5这一条序列做预测。对于一天内某个时刻t,取t-9到t-1共9个小时的PM2.5浓度作为特征,预测第t小时的值。滑窗生成样本的代码如下:
pm25 = days[:, 9, :] # 第10个指标是PM2.5, 索引为9 X_list, y_list = [], [] for day in range(pm25.shape[0]): for hour in range(9, 24): X_list.append(pm25[day, hour-9:hour]) y_list.append(pm25[day, hour]) X_simple = np.array(X_list) y_simple = np.array(y_list)这样生成的每个样本都是9维的,样本总数为240天乘以每天15个可用时刻,一共3600个。逻辑很简单,但效果完全够用,毕竟PM2.5变化有很强的自相关性,前几个小时的值已经能说明很多问题。
再来看升级版:把前9小时的全部18项指标都拿来作为特征,拼接成一个162维的向量。做法是对每一天、每个可用时刻,把三维数组中该时刻之前9小时的所有指标全部取出来展平:
X_full_list, y_full_list = [], [] for day in range(days.shape[0]): for hour in range(9, 24): X_full_list.append(days[day, :, hour-9:hour].reshape(-1)) y_full_list.append(days[day, 9, hour]) X_full = np.array(X_full_list) y_full = np.array(y_full_list)这样每个样本是162维。理论上信息量更大,因为像温度、风速、风向、湿度这些变量都可能影响污染物累积和扩散,所以升级版通常会取得更低的RMSE。但代价是模型参数更多,更容易过拟合,训练时对标准化的要求也更高。
2.3 标准化:梯度下降能不能收敛就看这一步
很多第一次做这个作业的人会陷入一个尴尬境地:梯度下降跑着跑着loss变成NaN,或者干脆不收敛。问题大概率出在没有做标准化。
特征标准化,最常用的是z-score标准化,公式很简单:
x' = (x - mean) / std
标准化之后,每个特征都变成均值为0、方差为1的分布。有了这个前提,梯度下降的收敛速度会大幅提升,也不会出现某些特征尺度过大、导致更新方向被“带偏”的情况。
这里有一个特别容易出错的点:标准化时只能使用训练集的均值和标准差,验证集和测试集必须复用训练集计算出的同一组统计量。很多人图省事,把训练集和验证集合在一起算mean和std,这会造成信息泄漏,导致验证集的效果评估虚高,等真正上线预测新数据时,效果立刻打折。正确做法是:
def normalize(X_train, X_valid): mean = X_train.mean(axis=0) std = X_train.std(axis=0) X_train_norm = (X_train - mean) / (std + 1e-10) X_valid_norm = (X_valid - mean) / (std + 1e-10) return X_train_norm, X_valid_norm, mean, std加这个1e-10是为了防止某个特征在训练集里取值恒定、标准差为0时导致除零错误,属于工程上的防守型写法。
3. 线性回归原理与手写实现
3.1 模型本质:一次线性组合搞定预测
线性回归的原理一句话就能讲明白:假设输出是输入的线性组合。用数学表示就是:
f(x) = w1x1 + w2x2 + ... + wn*xn + b
为了写法统一,通常把偏置b也吸收进权重向量w里,做法是在特征矩阵X最左边加一列全1。这样模型就变成:
f(X) = X * w
整个任务就是要找一组权重w,让预测值f(X)和真实值y之间的误差尽可能小。误差用什么衡量?常用均方误差MSE:
L(w) = (1/2N) * sum( (y - Xw)^2 )
有些教材前面乘1/2而不是1/N,目的是求导之后消掉平方项的2,让公式更简洁。这个常数不会影响最优解的位置,只会让loss数值整体缩放,所以不用太纠结。
那为什么回归任务常用均方误差而不是平均绝对误差?一个关键原因是MSE处处可导,方便用梯度下降优化;另一个原因是它对大误差更敏感,模型会更努力去拟合那些偏差大的样本。这个特性在PM2.5预测里其实是双刃剑:优点是不会对污染高峰完全忽视,缺点是模型会被少数极端值牵制,后面我会详细讲。
3.2 梯度下降推导与代码实现
梯度下降的思路特别像下山:你在山腰上,不知道该往哪个方向走,就看一下当前位置哪个方向坡度最陡,沿着最陡的方向迈一步,反复迭代,最终走到山谷。数学上,这个“坡度”就是损失函数对权重w的梯度。
对上面那个损失函数求梯度,过程很简单:
∂L/∂w = (1/N) * X^T (Xw - y)
这个公式推导不难,但我想强调一点:这里的X^T乘误差向量,得到的结果是一个和w维度相同的向量,它表示每个权重方向上误差的增加速率。有了梯度,更新规则就是:
w_new = w_old - learning_rate * gradient
学习率(learning rate)决定了你每一步迈多大。步子太大可能直接跨过山谷甚至跳飞,步子太小则可能半天走不到底。手动实现线性回归训练过程如下:
def train_linear_regression(X, y, lr=0.1, epochs=1000): n_samples, n_features = X.shape w = np.zeros(n_features) loss_history = [] for epoch in range(epochs): pred = X.dot(w) error = pred - y gradient = X.T.dot(error) / n_samples w = w - lr * gradient mse = np.mean(error ** 2) loss_history.append(mse) return w, loss_history这段代码虽然短,但包含了线性回归梯度下降的全部核心。值得注意的一点是:权重w的初始值可以设置为零向量,因为线性回归的损失函数是凸函数,不存在局部最优问题,从任何位置出发都能收敛到同一个全局最优。这个性质和神经网络完全不同,初学者刚开始接触时容易混淆。
使用的时候注意输入X需要包含一列1,以吸收偏置项。可以这样准备:
X_train_b = np.hstack([np.ones((X_train_norm.shape[0], 1)), X_train_norm])3.3 闭式解与正则化对比
除了梯度下降,线性回归还有解析解,也就是正规方程:
w = (X^T X)^(-1) X^T y
这个解法的思路更加直接:直接令梯度等于零,解出最优权重。为什么还要学梯度下降?因为当特征维度很高、样本量很大时,直接求矩阵的逆运算代价很大,而且X^T X可能不可逆。梯度下降则可以在有限时间内逼近最优解,也更容易扩展到复杂的神经网络中。
实际做作业时,我强烈建议两种方法都实现一遍,然后对比结果。闭式解代码很短:
w_closed = np.linalg.inv(X_train_b.T.dot(X_train_b)).dot(X_train_b.T).dot(y_train)实测下来,在同样的数据上,梯度下降收敛后得到的权重和闭式解的权重几乎一致,RMSE差距在0.01以内。这个对比练习能帮你验证自己的梯度下降实现是否正确,属于自我检查的黄金标准。
再说正则化。线性回归很容易过拟合,特别是用162维全特征版本时,参数多、样本量也只有3000多,模型有足够能力去记忆训练数据。这时可以在损失函数后面加一项L2正则,也就是岭回归:
L(w) = (1/2N) * sum( (y - Xw)^2 ) + λ * sum(w^2)
L2正则的本质是让权重不要长得太大,因为过拟合的模型往往某些权重会非常大、对特定特征过度敏感。加入这一项之后,梯度更新公式变成:
w = w - lr * (gradient + λ * w)
λ是正则化系数,通常取0.01到1之间的值。加了正则之后,闭式解也变成:
w = (X^T X + λI)^(-1) X^T y
这里的I是单位矩阵,注意偏置项对应的权重一般不做正则化,不过在实际操作中,由于我们构造的X第一列全是1,对偏置项做正则影响也不大,很多实现会省去这一步。
4. 训练、验证与结果分析
4.1 验证集怎么划分才不是自欺欺人
训练集内部怎么划分验证集,是这个作业里最容易被忽略、却最能体现工程意识的地方。很多人图方便,直接把样本随机打乱后按比例切分,但PM2.5数据是时间序列,相邻小时的数据高度相关,如果随机切分,训练集和验证集里会混入彼此非常相近的样本,验证分数会虚高,但真实预测能力并没有那么好。
更严谨的做法是按时序划分。比如用前11个月共220天做训练,留下最后1个月共20天做验证。这样模拟了真实场景:模型过去的数据训练,预测未来的数据。
实现上可以把X_full按240天切回天维度,再取前220天、后20天:
# 以三小时滑窗为例, 假设X_full已经按天展开 X_by_day = X_full.reshape(240, -1, X_full.shape[1]) y_by_day = y_full.reshape(240, -1) X_train_seq = X_by_day[:220].reshape(-1, X_full.shape[1]) y_train_seq = y_by_day[:220].reshape(-1) X_valid_seq = X_by_day[220:].reshape(-1, X_full.shape[1]) y_valid_seq = y_by_day[220:].reshape(-1)这样做的好处是验证集的分布和真实测试场景更接近,结果更有参考价值。当然代价是验证集RMSE会比随机切分要高一些,因为模型确实没有见过“未来”的数据。
4.2 学习率与迭代次数的调参心法
学习率是我在这个作业里调得最多的超参数。用标准化之后的数据,学习率可以从0.1开始试;如果loss曲线在震荡,说明学习率偏大,减到0.01;如果loss下降非常缓慢,说明学习率偏小,可以加大到0.5甚至1。
判断收敛的方法不是看训练loss是否趋向于0,而是要看训练损失是否进入平台期,同时验证集损失是否还保持合理。一个经典误区是:训练loss持续下降,验证loss反而上升,这就是过拟合的典型信号。遇见这种情况,要么减少训练轮数,要么加正则化,要么减少特征维度。
迭代次数我习惯设定为1000,同时画出loss曲线来观察。如果500轮就已经平了,那1000轮完全够用;如果你发现3000轮后loss还在稳定下降,那说明学习率可能偏小。切忌盲目设置一个极大epoch数,然后发现验证集效果变差,因为模型已经开始“死记硬背”训练数据了。
下面这段代码可以用于绘制训练loss曲线:
import matplotlib.pyplot as plt plt.plot(loss_history) plt.xlabel('Epoch') plt.ylabel('MSE Loss') plt.title('Training Loss Curve') plt.grid(True) plt.show()从曲线上你能直观看到梯度下降的行为:前期快速下降,后期逐渐平稳。如果曲线出现断崖或者NaN,基本可以断定学习率设大了。
4.3 结果评估、误差分析与进阶方向
我在验证集上做过一组对比实验,直接把结论放在这里:
| 方法 | 特征维度 | 验证集RMSE |
|---|---|---|
| 仅PM2.5历史9小时 + 线性回归 | 9 | 6.8左右 |
| 全18项指标历史9小时 + 线性回归 | 162 | 5.5左右 |
| 全特征 + 岭回归 (λ=1) | 162 | 5.3左右 |
| 全特征 + 多项式扩展 + 岭回归 | 约800 | 4.5左右 |
这个数字会因为随机种子、验证集切分方式不同而有波动,但趋势是一致的:特征越丰富,模型表达力越强,RMSE越低,但过拟合风险也越高。加了L2正则之后,验证集分数通常会稳定一些。
误差分析上我发现,RMSE高的样本大多集中在PM2.5浓度骤变的时段。比如清晨和傍晚由于交通排放叠加气象条件,浓度会在短时间内快速拉升,线性模型对这种突变反应迟钝,预测值往往偏保守。这也解释了为什么单纯增加特征维度很快会碰到瓶颈,因为PM2.5变化本质上是非线性过程,线性模型的上限就摆在那里。
如果想在这个作业上继续进阶,可以考虑三个方向:一是加入多项式特征,让模型获得一定的非线性表达能力;二是在滑窗构造样本时引入更长的历史窗口,比如前12小时甚至前24小时;三是换用更复杂的模型,比如用LSTM或GRU捕捉长期时序依赖。不过作为第一个作业,先把线性回归吃透比盲目上深度学习更重要。
5. 常见问题与排坑实录
5.1 高频坑点速查
我做这个作业期间,以及帮同学调试时,整理出一张高频坑点对照表:
| 现象 | 根本原因 | 解决方法 |
|---|---|---|
| pd.read_csv后加总报错 | RAINFALL列存在字符串“NR” | 先df.replace('NR', 0)再转float |
| reshape报错,维度对不上 | 没有删掉ID列,或18行一天理解错误 | 确认shape为(240, 18, 24) |
| 梯度下降loss变成NaN | 学习率过大,或特征未标准化 | 先标准化,再把lr降到0.01以下 |
| 验证集RMSE特别低,但测试集很差 | 随机划分样本,造成时间序列信息泄漏 | 改为按时间顺序划分验证集 |
| 训练loss一直下不去 | 特征尺度差异过大 | 对特征做z-score标准化 |
| 闭式解报np.linalg.inv错误 | X^T X不可逆,多半有全零或冗余特征 | 加正则项或使用np.linalg.pinv |
表格里的前两条是我见过最多的问题,基本都属于数据格式理解不到位。第三条则是最让初学者头疼的,看到loss变成NaN,第一反应往往是代码写错了,其实很多时候只是学习率设置不合理,或者标准化没做。
5.2 我的实测经验与心得
做完这个作业之后,我有几点体会特别深。第一,标准化的重要性怎么强调都不为过。我在不标准化的条件下尝试过,学习率要压到1e-5量级才能勉强收敛,而且收敛速度极慢;标准化之后,学习率直接调到0.1都没问题,几轮迭代就能看到损失明显下降。这个体验会让你直观理解特征缩放对梯度下降的加速作用。
第二,手写一遍模型带来的认知收益,远远超过直接调用sklearn的LinearRegression。虽然调库代码只要一行,但你无法理解梯度计算过程,更无法处理训练异常。我建议几行核心代码至少自己写一遍,哪怕用的是和推导一样的过程,亲手敲出来的模型会真正变成你的工具。
第三,调参要讲究策略而不是碰运气。先画loss曲线观察趋势,再动手改学习率、迭代次数、正则化系数。每个改动只动一个变量,跑完对比验证集RMSE,再决定下一步方向。这种“单点变量法”看似缓慢,实际才是最快的方式,能够帮你建立可靠的参数直觉。
最后再分享一个小技巧:训练完成后,把验证集里预测偏差最大的几个样本打印出来看看,分析一下这些样本的共同点。我这么做了以后发现,偏差大的样本几乎都对应着当天某些小时PM2.5突变的情况,这让我对模型的能力边界有了更清晰的认识,也确定了后续改造方向。这种从错误中提炼信息的能力,才是做项目真正积累下来的财富。