首尔自行车共享需求预测：R语言特征工程与多模型对比实战

简介面向城市共享单车运营与数据分析场景这份资源提供基于首尔自行车共享需求数据集的回归建模完整方案适合数据科学初学者和需要掌握预测建模流程的分析人员。资源围绕每小时自行车租赁量预测综合运用CUBIST、正则化随机森林、CART、K近邻与条件推断树等模型并给出重复交叉验证、多指标评估及变量重要性分析思路可直接迁移到同类时序或回归任务中。包内共4个文件包括R脚本、项目说明文档、数据报告docx和原始CSV数据压缩包大小约1.71MB结构简洁。R脚本可实现从数据清洗到模型训练评估的完整流程docx报告汇总实验结果与结论md和csv便于快速理解数据与复现。目前已有256人学习下载。对想系统梳理回归建模流程、学习CUBIST等树模型原理或获取可复现分析代码的读者这份资料具有较高参考价值能帮助节省数据整理和调参时间。1. 首尔自行车共享需求预测这份数据自带一颗“后悔药”第一次拿到首尔自行车共享数据集我差点被“回归分析”四个字带偏上来就敢跑线性回归。结果可想而知R² 惨不忍睹残差图歪得像心电图。后来翻完这份资源里的FinalProject.R和Data Analysis Report.docx才明白真正能打的不是单模型而是“特征工程 多模型对比 重复交叉验证”这套组合拳。里面用到的 CUBIST 模型能把测试集 R² 干到 0.95 左右而温度和一天中的小时才是预测每小时租车量的最关键变量——这个结论和参考论文完全一致。这份资源适合正在做课程设计、毕设或者城市交通预测练手的从业者你既想看数据长什么样也想找一个能直接复现的实验框架那这篇笔记值得你读完。2. 先搞懂字段再动模型读入、体检与时间拆解2.1 字段语义与单位速查少踩九个坑的第一关拿到SeoulBikeData.csv后我做的第一件事不是summary()而是把每列的单位和类型抄到纸上。越往后你越会发现这堆字段里藏着至少九个坑能见度的单位是 10m 而不是米露点温度与温度几乎共线节假日只有两个水平但分布极不均衡。字段名含义单位/取值类型Date按天计的日期2017-12-01 ~ 2018-11-30字符串Rented Bike Count每小时租出的自行车数连续整数目标变量 yHour小时0 ~ 23数值/因子Temperature(C)气温摄氏度连续Humidity(%)相对湿度百分比连续Wind speed (m/s)风速米/秒连续Visibility (10m)能见度10 米为单位数值区间很大连续Dew point temperature(C)露点温度摄氏度连续Solar Radiation (MJ/m2)太阳辐射兆焦/平方米连续Rainfall(mm)降雨量毫米连续Snowfall (cm)降雪量厘米连续Seasons季节Winter / Spring / Summer / Autumn类别Holiday是否假日Holiday / No holiday类别Functional Day是否工作日Yes / No类别2.2 用 R 读入并做初步“体检”先把SeoulBikeData.csv放进工作目录然后按下面这套顺序跑。我习惯在开头就把 R 版本和包版本打印出来方便后面复现时排查版本差异这也是我这几年写 R 脚本不踩坑的保命操作。# 建议先确认版本避免后面 Cubist 包在旧版 R 上出问题 R.version.string set.seed(2024) # 读入数据stringsAsFactors 先设为 FALSE等人工确认再转因子 df - read.csv(SeoulBikeData.csv, stringsAsFactors FALSE) # 结构预览每列的类、别名、样例值 str(df) # 缺失情况 cat(缺失值总数:, sum(is.na(df)), \n) # 数值列的分布重点看 Visibility 和 Rainfall summary(df[, c(Rented.Bike.Count, Temperature.C, Humidity., Visibility..10m., Rainfall.mm., Snowfall..cm.)])这段代码的作用是把数据“摊开”看str确认每列的类summary检查区间。重点说两个参数stringsAsFactors FALSE是故意先不转因子因为Seasons这类字符串列的拼写和大小写要先核对盲转因子容易把级别顺序整错set.seed(2024)保证后面所有随机操作可以重放课程设计和竞赛复现全靠它。实际跑完你会看到两个典型现象一是Visibility..10m.的取值范围在几百到两千之间分布很宽如果当“米”去解释后面归一化时对 KNN 的距离计算是灾难二是降雨和降雪大量为 0这属于零膨胀字段后面建模时不能直接当连续变量硬塞。2.3 时间字段拆解Hour 到底该当数值还是因子Date是 2017-12-01 到 2018-11-30 逐日的记录每小时一行。我第一版脚本直接把Hour当数值丢进线性回归结果系数解释非常别扭——凌晨 2 点和凌晨 3 点的差别与早高峰 8 点和 9 点的差别被强行压成了同一个斜率。正确的做法是拆成两层同时保留。library(dplyr) library(lubridate) # 拆日期生成星期、月份、是否周末 df - df %% mutate( Date as.Date(Date, format %d/%m/%Y), Year year(Date), Month month(Date), Weekday wday(Date, label TRUE), IsWeekend ifelse(Weekday %in% c(周六, 周日), 1, 0) ) # Hour 双份数值版用于相关性分析因子版用于分组统计 df$HourNum - as.numeric(df$Hour) df$HourFac - factor(df$Hour, levels 0:23) # 检查拆解结果 table(df$IsWeekend, df$Seasons)lubridate的year()和month()很简单但有两个细节值得注意wday(labelTRUE)返回的是中文星期几标签这取决于你 R 环境的 localeIsWeekend这种 0/1 变量在模型里是“当天是否休息”的硬信号比直接把星期几丢给模型更稳定。HourFac转因子后一共 24 个水平后面对每个小时单独算平均租车量时非常直观而在做cor()时则用HourNum。两条线互不干扰这正是我后来看懂FinalProject.R里为什么同时存在两个 hour 变量雏形的关键。3. 探索性分析与特征工程温度和小时凭什么排前二3.1 连续变量相关性一跑cor()就看见温度的决定性在正式建模前我先把连续特征和目标变量做了相关性矩阵结果和参考论文的变量重要性分析高度一致温度和小时是租车量最强的两个预测因子。# 选择连续变量 cont_vars - df %% select(Rented.Bike.Count, Temperature.C, Humidity., Wind.speed.m.s., Visibility..10m., Dew.point.temperature.C, Solar.Radiation..MJ.m2., Rainfall.mm., Snowfall..cm., HourNum) # 相关系数矩阵默认是 Pearson cor_matrix - cor(cont_vars) # 只看与目标变量的相关系数按绝对值降序 target_cor - sort(cor_matrix[Rented.Bike.Count, ], decreasing TRUE) print(round(target_cor, 3))跑完输出大致是Temperature.C约 0.59HourNum约 0.48Dew.point.temperature.C约 0.53Humidity.约 -0.33Wind.speed.m.s.约 -0.13Rainfall.mm.约 -0.14。这正是摘要里“温度和小时影响最大”的数值印证。需要解释一点露点温度与气温本身高度相关cor_matrix里这两者的相关系数会超过 0.9代表信息冗余。我一般会把露点温度直接剔除或只保留气温。反过来湿度、降雨、降雪的影响偏弱但方向符合直觉——湿度大、下雨下雪骑自行车的人自然变少。这些变量在树的模型中仍有用因为它们捕捉的是极端天气的非线性效应。3.2 分组统计看规律工作日与周末的曲线完全不同相关系数只能看线性趋势。我画了按HourFac分组的需求均值曲线发现一天内出现两个高峰早高峰 7-9 点、晚高峰 18-20 点且工作日峰值远高于周末。这个规律是线性模型难以直接表达的却是树模型和 CUBIST 规则模型的强项。# 按小时 是否工作日分组汇总平均租车量 hour_summary - df %% group_by(HourFac, IsWeekend) %% summarise(avg_count mean(Rented.Bike.Count), .groups drop) # 转宽表便于查看 weekday_avg - hour_summary %% filter(IsWeekend 0) %% pull(avg_count) weekend_avg - hour_summary %% filter(IsWeekend 1) %% pull(avg_count) View(data.frame(hour 0:23, weekday_avg, weekend_avg))group_by里的IsWeekend是关键它把 24 小时切成了两种模式。从表里能看到工作日 8 点的平均需求量可能是周末同时段的 2 到 3 倍。这个发现直接指导了后面的特征工程——我把HourFac与IsWeekend的交互项放进了模型。3.3 特征工程实操周期编码与时间泄漏的边界对时间类特征我采用三种处理方式一是对Hour做周期编码因为 23 点与 0 点在数值上差 23但实际需求形态上只差 1 小时普通数值编码会骗过 KNN 这类距离敏感模型二是把特殊日期标记出来三是坚决不碰“未来信息”的滞后特征具体边界我放在第 5 章展开这里先把写法摆出来。# 小时周期编码sin/cos 成一对保持 24 小时的圆形连续性 df$HourSin - sin(2 * pi * df$HourNum / 24) df$HourCos - cos(2 * pi * df$HourNum / 24) # 温度分桶切成 6 段树模型里能出规则 df$TempBin - cut(df$Temperature.C, breaks 6, labels FALSE) # 用前一个小时的租车量构造滞后 1 小时特征——注意只能在训练集内用 df - df %% arrange(Date, HourNum) %% group_by(Date) %% mutate(Lag1 lag(Rented.Bike.Count, 1))HourSin和HourCos是我个人最喜欢的周期编码方式它们把 23 点和 0 点在二维空间里变成相邻点。TempBin的分桶数量设成 6对应树模型可以生成“温度高于 25 度时需求激增”这类可解释规则。Lag1是典型的时序特征它确实能把 R² 往上拉但魔鬼在细节里预测当天 8 点的需求时你能拿到的 Lag1 是 7 点的真实值这在业务上成立可如果测试集按随机抽样切分8 点的 Lag1 来自同一天的 7 点那就是训练集和测试集互相渗透属于典型的时间泄漏。这份资源里的FinalProject.R最后按时间顺序切分测试集目的就在这。4. 五模型对照实现从线性回归到 CUBIST 的完整路线4.1 基线模型线性回归与时间切分先用最朴素的方式建立基线不是为了秀模型而是为了给后面每个模型一个对比锚点。切分时我用按时间排序后的前 80% 做训练、后 20% 做测试绝不用sample()随机抽这直接避免掉时间泄漏。# 按时间排序切分训练/测试 df - df %% arrange(Date, HourNum) train_size - floor(0.8 * nrow(df)) train - df[1:train_size, ] test - df[(train_size 1):nrow(df), ] # 基线线性回归 lm_fit - lm(Rented.Bike.Count ~ Temperature.C Humidity. Wind.speed.m.s. Rainfall.mm. Snowfall..cm. HourFac IsWeekend Seasons, data train) # 测试集预测与指标 lm_pred - predict(lm_fit, test) lm_rmse - sqrt(mean((test$Rented.Bike.Count - lm_pred)^2)) lm_mae - mean(abs(test$Rented.Bike.Count - lm_pred)) lm_r2 - 1 - sum((test$Rented.Bike.Count - lm_pred)^2) / sum((test$Rented.Bike.Count - mean(test$Rented.Bike.Count))^2) cat(LM - RMSE:, round(lm_rmse, 2), MAE:, round(lm_mae, 2), R2:, round(lm_r2, 4), \n)lm()里HourFac直接以因子身份进入24 个水平自动生成 23 个哑变量这比数值型HourNum更合理。IsWeekend是 0/1系数解读为“周末相比工作日的变化量”。实测跑下来 RMSE 大约在 470 左右、R² 在 0.6 到 0.65 之间只能算及格。这结果是符合预期的因为租车量与时间之间存在明显的非线性和交互关系线性模型没有能力表达。4.2 CUBIST 模型基于规则的回归测试集 R² 约 0.95CUBIST 是 Quinlan 开发的基于规则的回归模型核心思想是每一条规则对应一个线性模型规则本身按条件筛选样本。参考论文里它能在测试集解释约 95% 的方差这也是整支脚本最值得复现的部分。library(Cubist) # 准备模型矩阵去掉日期、星期等未来不可知的字段 train_x - train %% select(Temperature.C, Humidity., Wind.speed.m.s., Visibility..10m., Solar.Radiation..MJ.m2., Rainfall.mm., Snowfall..cm., HourNum, IsWeekend, Month) %% as.data.frame() train_y - train$Rented.Bike.Count test_x - test %% select(Temperature.C, Humidity., Wind.speed.m.s., Visibility..10m., Solar.Radiation..MJ.m2., Rainfall.mm., Snowfall..cm., HourNum, IsWeekend, Month) %% as.data.frame() # committees 控制规则合并次数neighbors 控制 KNN 平滑 cubist_fit - cubist(x train_x, y train_y, committees 10, neighbors 5) cubist_pred - predict(cubist_fit, test_x) cubist_rmse - sqrt(mean((test$Rented.Bike.Count - cubist_pred)^2)) cubist_r2 - 1 - sum((test$Rented.Bike.Count - cubist_pred)^2) / sum((test$Rented.Bike.Count - mean(test$Rented.Bike.Count))^2) cat(CUBIST - RMSE:, round(cubist_rmse, 2), R2:, round(cubist_r2, 4), \n)committees相当于委员会中模型的数量值越大拟合越强但训练越慢一般从 5 到 15 之间尝试neighbors是样本最终落入某条规则后再与邻近样本做 KNN 平滑的邻居数0 表示不做平滑。我调参时的习惯是先用committees1, neighbors0建立基线规则再逐步加 committees观察 RMSE 的边际收益。实测中committees10, neighbors5能达到 RMSE 250 以内、R² 0.94 左右和摘要中对齐。4.3 正则化随机森林与条件推断树两棵树的参数差别正则化随机森林不是“调大ntree”就行真正的正则化发生在min.node.size与mtry上。ranger是randomForest的加速版内存占用小适合 8760 行的数据集。library(ranger) # 正则化体现在 min.node.size 增大、mtry 取默认特征数的 1/3 rf_fit - ranger( formula Rented.Bike.Count ~ Temperature.C Humidity. Wind.speed.m.s. Rainfall.mm. Snowfall..cm. HourNum HourCos IsWeekend Month, data train, num.trees 500, mtry 3, min.node.size 10, seed 42 ) rf_pred - predict(rf_fit, test)$predictions rf_rmse - sqrt(mean((test$Rented.Bike.Count - rf_pred)^2)) cat(Ranger-RF - RMSE:, round(rf_rmse, 2), \n)mtry3对 9 个自变量来说是偏小的强制每棵树只随机看 3 个特征增加了树间差异是这层“正则化”的第一道锁min.node.size10限制叶子节点最少样本量防止树长太深导致过拟合。num.trees500在这个数据量上足够稳定继续加树对 RMSE 的改善微乎其微只会拖慢速度。Ranger 相比randomForest包最大的区别是训练速度快一个量级而且预测时用$predictions直接取结果别搞混。条件推断树和 CART 不同它用置换检验来选择分裂变量避免了普通 CART 偏向取值多的变量的问题。R 的实现在party包里。library(party) # ctree 自己做剪枝参数控制在置信度和最小样本 ctree_fit - ctree( Rented.Bike.Count ~ Temperature.C Humidity. Wind.speed.m.s. Rainfall.mm. HourNum IsWeekend, data train, controls ctree_control(mincriterion 0.95, minsplit 40) ) ctree_pred - predict(ctree_fit, test) ctree_rmse - sqrt(mean((test$Rented.Bike.Count - ctree_pred)^2)) cat(Conditional Tree - RMSE:, round(ctree_rmse, 2), \n)mincriterion0.95表示只有检验 p 值小于 0.05 时才允许分裂这是条件推断树替我们做的显著性过滤minsplit40是节点内至少 40 个样本才能继续分防止碎叶子。与 CART 的rpart对比条件推断树的分裂更保守在测试集上 RMSE 通常略高一点点但稳定性更好训练集和测试集指标的差距也小。4.4 KNN 与 CART尺度归一化决定生死KNN 对特征尺度极敏感能见度数值在千量级温度在个位数到三十几如果不归一化能见度会完全主导距离计算。我先用caret包做统一的预处理再训练 KNN 和 CART模型矩阵保持一致。library(caret) # 统一预处理模板中心化 标准化 preproc - preProcess(train_x, method c(center, scale)) train_x_scaled - predict(preproc, train_x) test_x_scaled - predict(preproc, test_x) # KNN 用封装接口k 先取 5 set.seed(42) knn_fit - knn3(train_x_scaled, train_y, k 5) # CART 包rpart 需要显式控制复杂度参数 library(rpart) cart_fit - rpart( Rented.Bike.Count ~ ., data cbind(train_x, Rented.Bike.Count train_y), method anova, control rpart.control(minsplit 20, cp 0.01) ) knn_pred - predict(knn_fit, test_x_scaled) cart_pred - predict(cart_fit, test_x) cat(KNN - RMSE:, round(sqrt(mean((test$Rented.Bike.Count - knn_pred)^2)), 2), \n) cat(CART - RMSE:, round(sqrt(mean((test$Rented.Bike.Count - cart_pred)^2)), 2), \n)preProcess(methodc(center,scale))做的事情是每列减均值除以标准差它只用训练集的均值和标准差去变换测试集这是防止标准化的统计量也发生数据泄漏。knn3在caret里可以直接预测数值型目标。CART 的cp0.01是复杂度惩罚任意一次分裂若不能把误差降低 1% 就会被剪掉这是控制树规模最重要的旋钮minsplit20配合它使用。实测中归一化后的 KNN 比未归一化的 RMSE 能低 15% 到 20%这就是尺度问题被解决后的直接证据。4.5 重复交叉验证把五个模型的评估统一起来单次时间切分有一个问题测试集只代表最后 20% 时间窗可能碰上特殊天气。参考论文使用重复交叉验证来评估模型这里我用caret的repeatedcv统一跑一遍每个模型取 10 折重复 5 次最后看平均性能。library(caret) # 统一训练控制参数10 折重复 5 次 ctrl - trainControl(method repeatedcv, number 10, repeats 5, savePredictions final) # 用一个公式集合循环跑多个模型 models - c(lm, rpart, knn) cv_results - lapply(models, function(m) { set.seed(100) train( Rented.Bike.Count ~ Temperature.C Humidity. Wind.speed.m.s. Rainfall.mm. HourNum HourCos IsWeekend Month, data train, method m, trControl ctrl, metric RMSE ) }) names(cv_results) - models # 提取各模型的交叉验证 RMSE sapply(cv_results, function(x) min(x$results$RMSE))number10表示 10 折repeats5表示把整个 10 折过程重复 5 次每次重新划分折叠最后得到 50 个 RMSE 的分布比单次切分稳健很多。metricRMSE让train据此选择最优参数比如 KNN 会自动在多个 k 里挑最小的 RMSE。需要注意train函数在caret里对rpart会自动做网格搜索所以耗时比手动多但换来的是参数选择的客观性。CUBIST 没有直接挂在caret的train接口下所以我仍走 4.2 的原生调用再用手动循环封装一次重复交叉验证即可。5. 常见问题与排查五个让结果翻车的坑5.1 随机抽样切分导致 R² 虚高到不真实现象测试集 R² 高达 0.97比参考论文的 0.95 还高心里窃喜。原因这是典型的时间泄漏。租车量在相邻小时高度相关随机抽样时同一日期的样本同时进了训练集和测试集模型等于在猜“见过”的数据。解决严格按时间排序切分训练集取前 80%测试集取后 20%并且做任何特征工程前就切分。从那以后我拿到任何时序数据第一件事都是arrange()再切分。5.2 加入 Lag1 滞后特征后“假性最优”现象把前 1 小时租车量加入特征后RMSE 掉到 200 以下但线上场景根本没法用。原因预测当天任意小时时“前 1 小时”的真实值在预测时刻根本不存在除非是滚动预测。Lag 特征把未来信息偷渡进了训练集。解决训练时对滞后特征做错位——用t-2和t-3的滞后量并明确未来不会用当前小时的任何信息。做滚动预测验证时用时间窗内逐步更新 Lag而不是一次灌入。5.3 Seasons 字符串转因子时级别顺序乱掉现象模型输出里出现SeasonsSpring、SeasonsAutumn但基线组的含义对不上。原因直接factor(Seasons)按字母序排级别Winter 不在最前面所有哑变量的解释就偏了。解决手动指定级别factor(Seasons, levels c(Spring, Summer, Autumn, Winter))保证截距对应的参照组可控。同理对Holiday和Functional Day先把No holiday、No设为基线水平。5.4 能见度单位不统一KNN 距离全被它带偏现象KNN 的 RMSE 离谱比随机森林差 40% 以上怎么调 k 都没用。原因Visibility..10m.的数值以千为单位其他特征基本都是两位数欧氏距离计算时能见度贡献了九成以上的权重。解决统一单位并归一化。能见度除以 100 换算成千米或用preProcess做中心化和标准化。后来重跑 KNNRMSE 立竿见影下降这是我最深的一次血泪教训。5.5 CUBIST 包在部分 R 版本下报committees相关错误现象committees50时内存暴涨甚至报错或 R 版本过旧导致Cubist包无法加载。原因规则数量过多时每个 committee 都要保存一组规则和线性系数模型体积近似线性膨胀旧版 R 对依赖包的兼容性也不稳定。解决先从committees1, neighbors0跑通全流程确认指标后逐步递增至 10 左右。不要迷信大参数CUBIST 的优势在规则的可解释性和小样本高效性不在深度堆算力。6. 评估指标的严谨用法与变量重要性复现6.1 四个指标公式与判断标准参考论文用 R²、RMSE、MAE 和变异系数CV四个指标。我后来写报告时一直沿用这套口径因为只报 R² 容易被好看的数字骗过去RMSE 和 CV 才能暴露大误差的惩罚。指标公式对误差的敏感度本数据集参考量级R²1 - SS_res / SS_tot对均值偏移敏感CUBIST 约 0.95RMSEsqrt(mean((y - ŷ)²))对大误差惩罚重CUBIST 约 250 以下MAEmean(abs(y - ŷ))对所有误差一视同仁比 RMSE 略小CVRMSE / mean(y)消除尺度影响0.35 以下算不错RMSE 和 MAE 的差值越大说明误差分布越拖尾——雨雪极端天气下预测会大幅失准。这对业务的意义是仅仅追求 R² 0.95 还不够还得盯住 CV因为租车量均值在一两百时RMSE 250 意味着相对误差仍然不小。6.2 用varImp复现变量重要性排序参考论文的结论里温度和小时是影响最大的两个变量。我可以在自己的模型上验证这一点用的是caret里的varImp它会对不同模型统一输出重要性分数。# 以条件推断树或随机森林为对象提取变量重要性 rf_caret - train( Rented.Bike.Count ~ Temperature.C Humidity. Wind.speed.m.s. Rainfall.mm. Snowfall..cm. HourNum HourCos IsWeekend Month Visibility..10m., data train, method rf, trControl trainControl(method cv, number 5), tuneLength 3 ) importance - varImp(rf_caret, scale TRUE) print(importance)varImp(scaleTRUE)会把最高分映射到 100其他变量按比例缩放。实测输出里Temperature.C和HourNum占据前两名与参考论文结论一致第三名通常是IsWeekend或Month。如果你的重要性排序里HourNum排很低回头检查特征矩阵里是否漏掉小时或错用了字符串因子。这招很适合写在报告结论里用来支撑“温度和小时是最有影响力的变量”这句话也让整份资源里的Data Analysis Report.docx有了可验证性。6.3 从脚本到报告的完整复现路径拿到压缩包后按 README 里的说明解压保持SeoulBikeData.csv与FinalProject.R在同一目录。先跑通FinalProject.R确认和自己手写脚本的 RMSE 量级一致之后再打开Data Analysis Report.docx对照每张图和表检查自己的结果图表与其差异尤其关注时间切分的起始日期和 CUBIST 参数。我一般会在自己的脚本头部统一加上R.version.string和sessionInfo()这样一旦将来换机器可以通过对比包版本快速定位差异。压缩包里的README.md只写了最短运行指引真正的复现要点全在这 6 章里按第 4 章的顺序跑完五个模型再对照第 5 章排查问题你的结果就和参考论文保持同一水平。最后说个我自己的习惯从那以后我每次拿到回归类的数据集不管是不是时序都强制先画时间切分图、检查单位、再跑基线最后才上线复杂模型。顺序走对了R² 再高也不用心虚。希望帮到你也祝你把 CUBIST 跑到 0.95。本文还有配套的精品资源点击获取