ARTICLE DETAIL

资讯详情

深耕编程入门与网站建设的一线实战洞察。

线性回归预测PM2.5:特征工程、数据清洗与模型评估全流程解析

线性回归预测PM2.5:特征工程、数据清洗与模型评估全流程解析 简介这是一份基于线性回归的PM2.5预测系统Python源码适合机器学习初学者、数据分析课程设计者参考。资源共19个文件压缩包大小2.14MB其中12个csv文件包含训练集、测试集以及经特征工程处理后的中间数据如特征拼接表、标准答案5个png图片直观展示数据分布或预测效果另有model.npy保存训练得到的模型参数主程序.py实现从数据读取、清洗、特征构造到模型训练与预测输出的完整流程。源码以繁体中文编码的空气质量观测数据为输入逐步骤演示了big5文件读取、numpy矩阵运算、最小二乘求解及csv结果导出等关键操作并附带预测结果样例便于逐行对照理解线性回归在空气污染预测中的实际应用。已有594人学习下载可据此快速复现实验也可结合自身数据替换修改用于课程报告或毕业设计。1. 用线性回归做PM2.5预测这个老模型为什么仍值得跑一遍去年给某地级市做空气质量数据分析时对方上来就问能不能上深度学习我反问一句你们现在有能用的基线吗对方摇头。于是我先把历史监测数据喂给了线性回归结果这个十几分钟跑完的模型在预报次日均值上的表现直接超过了他们之前人工估算的方案。线性回归预测PM2.5的意义在于它是所有预测系统的基线也是检验特征工程是否有效的试金石。模型本身没有秘密秘密全在特征和数据处理上。这套源码适合两类人——想入门机器学习、需要一个能跑通全流程的完整案例的新手以及需要快速搭建业务基线、给领导一个可解释结论的工程师。2. 把监测数据整理成特征表时间戳、滞后项与缺失值处理2.1 先看数据和字段这套源码最常见的输入格式PM2.5预测系统无论代码怎么写第一步永远是数据整理。绝大多数监测平台导出的原始数据是CSV或Excel表格每行是一个整点记录列包含时间戳、PM2.5浓度、PM10、SO2、NO2、CO、O3以及气象要素如温度、湿度、风速、气压。这套源码里我一般会把数据文件命名为data.csv字段名做成英文小写避免中文列名在pandas里反复出幺蛾子。import pandas as pd df pd.read_csv(data.csv, parse_dates[time]) df df.sort_values(time).reset_index(dropTrue) print(df.head()) print(df.dtypes)逻辑说明parse_dates把时间列解析成 datetime 类型sort_values按时间排序。很多原始数据是乱序的如果跳过排序后面做滞后特征时shift会把不同日期的数据混在一起模型立刻就废了。这一步是整个数据预处理的第一道关卡。值得提醒的是现场对接的数据经常混着小数点和千分位逗号或者PM2.5浓度列里有—这种占位符。读取时应当加na_values[—, , NULL]把这些统一转成缺失值而不是让pandas把它们读成字符串。源码包里如果附带了数据说明文档第一件事永远是确认字段含义和单位而不是急着跑模型。2.2 构造时间与滞后特征让线性回归「看见」昨天的风线性回归没有记忆能力它看到的每一行都是独立的。如果只给模型当前时刻的PM10和气象数据模型完全不知道昨天发生了什么预测效果会非常差。因此要做滞后特征lag features把过去几小时的浓度、过去24小时的平均值作为新列加入。for lag in [1, 3, 6, 24]: df[fpm25_lag_{lag}] df[pm25].shift(lag) df[pm25_lag24_mean] df[pm25].rolling(24).mean().shift(1) df[hour] df[time].dt.hour df[weekday] df[time].dt.weekday df df.dropna().reset_index(dropTrue)逻辑说明shift(lag)表示取前 lag 行的数值比如pm25_lag_1就是上一小时的PM2.5浓度。rolling(24).mean().shift(1)是过去24小时滑动平均值注意shift(1)必不可少否则当前时刻的浓度会被包含进去造成数据泄漏。hour和weekday让模型能捕捉早晚高峰和周末效应。时刻、滞后浓度、滑动均值这几类特征组合起来线性回归的效果会有明显跃升。参数怎么调滞后阶数的选择取决于业务场景。做小时级预测时lag 取 1、3、6、24 基本够用做次日日均值预测时lag 应该取前一天相同时刻的值并用过去7天的均值做平滑。这个参数没有固定最优解可以先用少量组合跑一轮看RMSE变化再决定。2.3 缺失值与异常值处理坏数据才决定RMSE下限监测数据最常见的缺失发生在设备校准时段和通信中断时表现为连续几小时为NaN。而异常值更典型冬季某天PM2.5实测值忽然飙到800可能是秸秆焚烧也可能是传感器故障。线性回归对异常值极其敏感——因为它最小化的是平方误差一个极端值就能把回归系数拉偏不少。df[pm25] df[pm25].replace(0, np.nan) df[pm25] df[pm25].interpolate(methodlinear, limit6) df df[df[pm25] 500]逻辑说明第一行把0值替换成缺失因为PM2.5浓度出现0基本是仪器归零错误线性插值适合短时间缺失超过6小时连续缺失建议用前后一天的同一时刻均值填充而不是插值df[pm25] 500是过滤异常大值阈值按地区情况调整华北地区冬季400多并不罕见但500以上基本都是数据质量问题。处理缺失值时我的习惯是保留一个is_missing标记列给模型让它自己学习缺失时段与其他特征的关系这对线性回归也是有效特征。3. 写最小可跑的线性回归代码sklearn 版与手写梯度下降版3.1 sklearn 版本十分钟跑通基线拿到处理好的特征表之后建模本身反而是最快的一步。这套源码里train.py 核心代码不超过二十行。使用 scikit-learn 的LinearRegression是最靠谱的起点它默认使用最小二乘法求解数据量在几万行内几乎瞬时完成。import pandas as pd from sklearn.model_selection import train_test_split from sklearn.linear_model import LinearRegression from sklearn.metrics import mean_squared_error, r2_score df pd.read_csv(features.csv) feature_cols [c for c in df.columns if c not in [time, pm25]] X df[feature_cols].values y df[pm25].values X_train, X_test, y_train, y_test train_test_split( X, y, test_size0.2, shuffleFalse ) model LinearRegression() model.fit(X_train, y_train) y_pred model.predict(X_test) rmse mean_squared_error(y_test, y_pred, squaredFalse) print(RMSE:, round(rmse, 2)) print(R2:, round(r2_score(y_test, y_pred), 4))逻辑说明这段代码的关键在shuffleFalse。时间序列数据不能随机打乱后切分否则测试集里混着训练集前后时间段的数据模型等于开卷考试RMSE会虚低。squaredFalse让mean_squared_error直接返回RMSE均方根误差单位与PM2.5浓度一致比较好解释。如果你刚装好Python环境这一步需要先安装依赖pip install pandas numpy scikit-learn。版本上 scikit-learn 1.0 以上即可不需要追新。遇到ModuleNotFoundError时先确认是不是装进了别的虚拟环境这个坑我见得太多了。3.2 手写梯度下降版理解学习率和收敛sklearn 把细节都封装了但作为工程落地理解参数怎么调仍然绕不开底层原理。这套源码如果附带手写版本一般会包含一个用梯度下降训练线性回归的文件它能让读者直观看到学习率如何影响收敛。import numpy as np class LinearRegressionGD: def __init__(self, learning_rate0.01, epochs1000): self.lr learning_rate self.epochs epochs self.weights None self.bias 0 def fit(self, X, y): n_samples, n_features X.shape self.weights np.zeros(n_features) for _ in range(self.epochs): y_pred X self.weights self.bias dw (2 / n_samples) * X.T (y_pred - y) db (2 / n_samples) * np.sum(y_pred - y) self.weights - self.lr * dw self.bias - self.lr * db def predict(self, X): return X self.weights self.bias逻辑说明核心是每次迭代计算损失函数对所有参数的偏导数然后沿负梯度方向更新参数。dw是权重梯度db是偏置梯度learning_rate决定每一步走多远。数据量较大时矩阵乘法X.T (y_pred - y)仍然是全量计算这是批梯度下降。参数怎么设学习率 0.01 对标准化后的数据通常能稳定收敛学习率 1.0 时损失函数会震荡甚至爆炸。特征没有做标准化时不同特征的量纲差异会导致梯度方向偏斜收敛慢。因此使用手写版本前务必对特征做标准化或者把学习率调到 0.0001 级别否则loss曲线会像心电图一样乱跳。3.3 源码工程的目录组织train.py / predict.py / features.py一个能交付的预测系统源码文件组织至少要让人一眼看出「数据从哪来、特征在哪造、训练怎么跑、预测怎么调」。这套源码包里的典型结构是文件职责load_data.py读取原始CSV字段清洗返回干净的DataFramefeatures.py滞后特征、滑动平均、时间特征构造train.py训练模型输出RMSE、R²保存模型到本地predict.py加载模型对最新一条数据做下一时刻预测model.pkl训练产出的模型文件predict 阶段反序列化使用我一般会在predict.py里也放一份特征构造代码保证预测时的输入格式与训练时完全一致。很多翻车现场就是训练时用了pm25_lag_1预测时忘了构造这一列导致维度对不上直接报错。如果源码包的README里没有强调这一点你接手后要主动把它对齐。4. 模型评估与参数调优R²、RMSE 与三个必调参数4.1 评价指标怎么选RMSE 惩罚大偏差MAE 更稳预测系统的评价指标决定了调参方向。PM2.5预测里最常用的是RMSE和MAE。RMSE把误差平方后求平均再开方因此它对那些偏差特别大的预测点非常敏感。比如某天实际浓度250模型预测150这个100的误差在平方之后会主导整个RMSE数值。MAE则只看绝对误差的平均值更平稳但对大偏差不敏感。from sklearn.metrics import mean_absolute_error, mean_squared_error, r2_score rmse mean_squared_error(y_test, y_pred, squaredFalse) mae mean_absolute_error(y_test, y_pred) r2 r2_score(y_test, y_pred) print(fRMSE: {rmse:.2f} MAE: {mae:.2f} R2: {r2:.4f})逻辑说明三个指标一起打印能看出模型误差分布的特征。如果RMSE明显大于MAE说明存在少数预测得很离谱的样本如果两者接近说明误差分布比较均匀。对于PM2.5业务而言污染峰值时段的预测偏差才是考核重点因为重污染天的应急响应才是刚需。所以我在项目里会额外关注「实际浓度200时模型的RMSE」这个指标比全局RMSE更能说明问题。评估代码里还要关注r2_score。R²接近1说明模型解释了大部分方差但R²低不一定意味着模型没用——如果测试集正好赶上一次极端沙尘天气R²会很难看可模型在常规天气下其实表现不差。所以判断模型质量时不要只看单个指标要结合时间序列图观察预测曲线与实际曲线的跟随程度。4.2 三个必调参数滞后窗口、训练集比例、是否加正则线性回归的参数不多但这三个地方值得花时间调。首先是滞后窗口的长度。滞后项太少模型缺少历史信息滞后项太多特征之间相关性变强数值稳定性变差。我的经验是小时级预测用[1, 3, 6, 24]四档日均值预测用[24, 48, 72]加7日滑动平均在RMSE上会有明显差异。第二个是训练集与测试集的比例和切分方式。时间序列模型不能用随机切分而是按时间顺序切。常见做法是前80%训练、后20%测试这样测试集模拟的是「用历史预测未来」的真实场景。如果你想严谨一点可以再用最后两周数据做验证集用于调参避免在测试集上反复调参造成信息泄漏。第三个是是否使用正则化。当特征数量多、共线性严重时普通线性回归的系数方差会变大。此时把LinearRegression换成Ridge加上 L2 正则往往能让测试集上的RMSE下降几个点。正则化系数alpha的典型值从 0.01 到 10 之间网格搜索选验证集上RMSE最小的值这是性价比最高的调参操作。from sklearn.linear_model import Ridge from sklearn.model_selection import GridSearchCV param_grid {alpha: [0.01, 0.1, 1.0, 10.0]} ridge Ridge() grid GridSearchCV(ridge, param_grid, cv5, scoringneg_mean_squared_error) grid.fit(X_train, y_train) print(best alpha:, grid.best_params_)逻辑说明GridSearchCV在训练集内部做交叉验证来选超参数而不是直接用测试集评估。因为时间序列的连续性这里cv5用的是默认的KFold划分严格来说会在交叉验证中造成时间错位。更严谨的做法是用TimeSeriesSplit代替默认KFold这点后面避坑章节会细说。neg_mean_squared_error是sklearn的常规写法sklearn的评分规则是越大越好所以MSE要取负号。4.3 用交叉验证验证稳定性别让一次划分决定结果单次划分训练集和测试集具有偶然性——测试集恰好落在重污染过程或恰好全是好天结果天差地别。我一般会用TimeSeriesSplit做滚动验证把数据切成多段每段都用前面的数据训练、后面的数据测试综合多次结果看稳定性。from sklearn.model_selection import TimeSeriesSplit tscv TimeSeriesSplit(n_splits5) rmse_list [] for train_idx, test_idx in tscv.split(X): X_tr, X_te X[train_idx], X[test_idx] y_tr, y_te y[train_idx], y[test_idx] m LinearRegression() m.fit(X_tr, y_tr) pred m.predict(X_te) rmse_list.append(mean_squared_error(y_te, pred, squaredFalse)) print(RMSE mean:, round(np.mean(rmse_list), 2), std:, round(np.std(rmse_list), 2))逻辑说明TimeSeriesSplit保证每次划分中训练集永远在测试集之前符合时间序列预测的真实约束。输出5次RMSE的均值和标准差均值是模型的真实水平标准差是稳定性。如果标准差很大说明模型在不同时间段表现波动大需要回头检查特征是否遗漏了季节性信息。5. PM2.5预测避坑指南这些坑我都替你踩过5.1 数据泄漏把未来数据混进了特征现象训练集上R²高达0.98测试集也接近0.97模型表现好得不真实部署到线上后预测值却一直滞后。原因构造滞后特征时用了当前时刻的数据算滑动平均比如用当天的PM2.5参与计算当天的预测特征。测试阶段没问题因为特征是从同一时刻数据里算出来的但线上预测时你还没有当前时刻的PM2.5真实值特征构造不完整模型只能使用残缺输入。解决所有特征必须来自「预测时刻之前」的数据。滑动平均构造完要shift(1)确保特征矩阵里第 i 行的数据只包含第 i 行以前的信息。写代码时养成一个习惯构造完特征后随机抽一行手动验证一下该行特征是否都能用历史数据独立算出。5.2 用 R² 判断模型好坏会误判现象模型A的R²为0.90模型B的R²为0.82团队选择上线A。实际运行后发现A在高浓度时段的误差巨大。原因R² 衡量的是模型解释的方差比例对均值附近的拟合程度非常敏感。PM2.5数据大部分时间在100以下模型只要把低浓度时段拟合好R²就会很高高浓度时段样本少对R²贡献小模型即使在这些时段完全失效也不影响R²。解决用RMSE作为调参主指标同时单独统计重污染时段实测值150的RMSE。只有这两个指标同时优化模型才有业务价值。在汇报时也习惯把两个数字一起放出来避免其他人被单一指标带偏。5.3 时间序列数据不能用随机切分现象train_test_split默认test_size0.2且shuffleTrue跑出来的RMSE远好于预期模型看起来完美知道未来趋势。原因随机打乱后训练集和测试集的数据来自同一时间段测试集里包含了与训练集相邻时刻的信息。比如测试集取了下午3点的数据训练集里可能就有下午2点的数据。时间序列存在强自相关性模型相当于在「默写」测试集。解决一律使用shuffleFalse按时间先后切成前80%和后20%。更严格的方案是用TimeSeriesSplit滚动验证并保证训练集和测试集之间留出至少一天的间隔gap避免滑动平均特征跨越切分边界造成泄漏。5.4 缺失值填充错误让滞后特征「空转」现象数据缺失率约15%直接用dropna()删掉所有含缺失的行结果剩下的数据集中在秋冬季模型在夏季表现一塌糊涂。原因缺失不是均匀分布的——设备检修常常集中在固定的时间段直接删行相当于重塑了数据分布。另外如果先做滞后特征再删行一行缺失会导致后面多个滞后特征全变NaN删除量比预期大得多。解决先填充缺失值再做滞后特征。短时间缺失用线性插值长时间缺失用前后同时间点的均值填充实在不行就保留NaN让树模型类算法处理。线性回归不能接受NaN所以必须在特征构造前把缺失值解决干净。5.5 极端浓度把回归系数拉偏现象模型各项指标看起来可以但某次沙尘天气后重新训练RMSE整体变大且预测值系统性偏高。原因一次极端事件产生的样本点在损失函数里被平方放大为了拟合这少数几天模型牺牲了大量常规天气样本的精度。我在实际项目中遇到过PM2.5实测值飚到900多的记录事后确认是传感器故障但模型已经为此付出代价。解决建模前画分布图把明显超过物理上限比如500的样本挑出来人工核对。确属真实污染事件时可以保留但要用截断或加权方式降低其影响比如把极端值截断到某个百分位数。对业务预测来说保住常规天气的精度比死磕极端值更有价值。6. 进阶从单步预测到多步预测的滞后窗口策略前面的模型预测的是下一小时的PM2.5浓度属于单步预测也是这套源码的核心功能。但实际业务中常常要的是未来24小时或未来3天的浓度变化曲线这时候单纯调参解决不了需要换预测策略。两种常见做法一种是递推预测即把上一时刻的预测值当作历史值放回特征里继续预测下一步——优点是代码简单、依赖少缺点是误差会逐级累积预测周期越长越不靠谱另一种是直接多步预测直接训练24个独立模型分别预测未来1小时、2小时到24小时的浓度——每个模型各自构建对应的滞后特征缺点是训练成本增大但避免了误差累积问题我在项目里更倾向这种方式。def train_direct_multistep(df, horizons[1, 6, 24]): models {} for h in horizons: df[ftarget_{h}] df[pm25].shift(-h) tmp df.dropna() X tmp[[c for c in tmp.columns if c.startswith(feat_)]].values y tmp[ftarget_{h}].values m LinearRegression() m.fit(X, y) models[h] m return models逻辑说明对每个预测步长h把目标列向上移动h行即第 t 行的目标是 th 时刻的浓度。注意shift(-h)后尾部会多出 h 行NaN需要删掉。每个模型都用自己对应的目标列训练预测时分别调用不同模型得到未来曲线上的离散点。我个人的经验是做小时级预测系统时递推预测能用就用省事做决策支撑比如是否启动重污染应急预案时直接多步预测的可靠性更值得信赖。这套源码拿来在实际业务中跑最重要的事不是追求R²再高零点几而是把特征构造、数据切分、评估口径这三件基本功做扎实。很多团队翻车不是模型不够强是数据管线和评估方法出了问题再花哨的算法也救不回来。希望帮到你。本文还有配套的精品资源点击获取
返回列表