简介：GEFCOM2014-EPFL 能源负荷预测资料包面向电力行业数据科学学习者、R 语言用户及时间序列预测课程实践者，用于借助 GEFCOM2014 历史小时级负荷数据练习消耗模式分析与未来负荷预测。压缩包整体约 111.53MB，页面文件总数标注为 0、类型明细暂未提供，实际目录与文件构成需下载后核对。页面显示已有 1426 人学习或下载，具备一定参考热度。资料内容围绕多地区负荷数据展开，涉及时间序列对象构建、季节性趋势分解、ARIMA 与状态空间模型、随机森林/支持向量机/神经网络等预测方法，也覆盖温度与节假日等特征工程、MSE/RMSE/MAE/R² 评估、交叉验证调参及集成学习思路，并提及 ggplot2、forecast、caret、dplyr、caretEnsemble 等 R 工具链，适合需要复现能源负荷预测流程、对照多种模型或完成课程设计的数据分析人员。

1. 从一次备用容量争论说起：GEFCOM2014 能源负荷预测到底在预测什么

做负荷预测的人迟早会遇到这个问题：模型输出的那条线，到底该往上留多少余量。前几年帮一个园区做日前负荷，点预测的 MAPE 已经压到 3% 以内，调度那边还是不敢用。他们真正想要的是“明天下午三点有九成把握不超过多少兆瓦”，一条能覆盖风险的曲线，而不是一个数。GEFCOM2014 能源负荷预测赛题把这个需求直接写进了评分规则：每个待预测时刻要交 99 个分位数，从 0.01 到 0.99，用 pinball loss 打分，偏保守和偏激进都要受到处罚。

gefcom2014-epfl 这份仓库就是围绕这个赛题写的一整套 R 实现，从原始 CSV 读入、气象站温度映射、特征构造、分位数模型训练，一直到提交文件生成，整个流程能跑通。它不负责教统计，它负责展示一条能跑起来的路径。想搞清楚概率预测整套流程，或者手里已经有一份自己的负荷数据、想补上区间输出，从这套代码入手比查阅文档快得多。

2. 项目骨架与技术选型：R 里搭一条从原始 CSV 到分位数提交的完整流程

2.1 概率预测和点预测的差别不只是一个区间

先把评分规则说明白，后面所有设计决策都来自这里。赛题用的 pinball loss 对单个分位数 q 是这样定义的：

# y 是真实值，yhat 是 q 分位数的预测值 # d > 0 说明真实值高过预测，惩罚 q * d（低估要罚） # d < 0 说明真实值低于预测，惩罚 (q - 1) * d = (1 - q) * |d|（高估也要罚） pinball <- function(y, yhat, q) { d <- y - yhat mean(ifelse(d >= 0, q * d, (q - 1) * d)) }

这段逻辑看着简单，但它决定了一件事：分位数回归不能只拟合均值。q=0.9 的模型会主动往高预测，因为低估的代价是高估的 9 倍；q=0.1 的模型反过来。所以同一套特征要训 99 次，每次目标函数不一样。这也是为什么不能拿一个点预测模型加个正态假设去凑区间——负荷的残差分布明显是偏的，冬天早晨的尾部特别厚。

常见的手法是先训中心分位数（q=0.5）拿到中位预测，再用它去初始化其他分位数的模型，收敛快一些。不过这里为了代码结构清晰，99 个分位数是各自独立训的，方便并行化，出问题也好定位是哪个 q 出问题。

2.2 为什么是 R、quantreg 和 gbm，而不是深度学习

数据规模是个容易忽略的细节。7 个完整年份的逐小时数据，20 个分区，展开成“分区 × 时刻”的长表大概一百多万行；特征加上小时哑变量、星期、温度 hinge、滞后项之后，也就几十维。这个体量下，深度学习的优势根本发挥不出来，反而会带来一堆调参负担和不可解释性。

R 在这个场景的优势很具体：quantreg::rq直接给你的就是线性分位数回归，一行公式能跑；gbm的distribution = list(name = "quantile", alpha = q)支持分位数目标，配合树模型能自动捕捉温度和小时的交互；data.table处理百万行级别的按组聚合不吃力；lubridate处理夏令时和节假日也算顺手。整套流程都在这几个包里完成，不需要跨语言传递数据。

选型上我一般会采用两条方案并行的策略：线性分位数回归当基准，因为它稳定、可解释、跑得快，出问题能当参照物；分位数 GBDT 当主力，因为它能自己学出“周六下午 + 高温”这种组合效应。两个模型的预测结果可以做加权平均，权重在验证集上调——这比单押一个模型安全得多。

需要说明的是，树模型对分位数目标的训练成本比回归高不少：99 个分位数意味着 99 次训练。实践中会做两个妥协，一是把分位数抽稀（比如只训 0.05 到 0.95 的 19 个，中间线性插值），二是降低树的复杂度（interaction.depth控制在 4 以内，n.trees800 左右，shrinkage0.03）。

2.3 目录结构和一条能跑通的执行链路

仓库的组织方式大致是这样，按职责切分而不是按模型切分：

gefcom2014-epfl/ ├── R/ │ ├── load_data.R # 原始 CSV 读入、时区统一、缺失值标记 │ ├── mapping.R # 气象站到分区的温度映射 │ ├── features.R # 日历特征、温度 hinge、滞后项 │ ├── models_linear.R # quantreg 分位数线性回归 │ ├── models_gbm.R # gbm 分位数提升树 │ ├── evaluate.R # pinball loss、分位数交叉检查 │ └── run_all.R # 主入口，按分区循环调度 └── data/ ├── raw/ # 赛题原始负荷与温度文件 └── processed/ # 中间宽表，方便反复调试

run_all.R是主入口，逻辑很直白：读数据 → 建映射 → 造特征 → 按分区循环训模型 → 拼提交表。按分区循环这点值得说，因为赛题里 20 个分区的负荷形态差别很大，有的以居民为主，双休日负荷塌得厉害；有的以工业为主，周末几乎不掉。混在一起训一个全局模型，树模型会被强行平均掉，尾部全糊。分开训虽然要多跑 20 遍，但每个模型都能吃满自己的数据。

代价是内存。20 个分区 × 99 个分位数 × 若干模型对象，如果不在循环里及时清理，R 的 GC 会明显拖慢速度。我一般的做法是每个分区处理完就存 RDS 到data/processed/，主进程只保留提交用的预测矩阵。

2.4 环境依赖和第一次运行

依赖装齐之后，第一次运行先从单个分区、单个分位数开始，确认链路通，再扩到全量。这一步能省掉大量排查时间。

# 依赖安装，版本用近几年的 R 都行 install.packages(c( "data.table", # 长表处理 "lubridate", # 时区与节假日 "quantreg", # 线性分位数回归 "gbm", # 分位数提升树 "doParallel" # 分区并行 )) # 先跑通一个分区 source("R/load_data.R") source("R/features.R") dt <- load_zone("data/raw/load_zone01.csv", "data/raw/temp.csv") feat <- build_features(dt, zone_id = 1) nrow(feat) # 确认行数对得上，先别急着训模型

load_zone的第二个参数是温度文件，第一个是单分区负荷。先把行数和时间跨度打出来核对，比直接扔进模型里强。特征表一旦有重复时戳或者缺行，后面所有分位数的结果都是错的，而且错得很隐蔽——pinball loss 会整体偏高一点点，但不至于明显到让你立刻发现。

3. 数据层拆解：分区负荷、气象站温度和日历特征怎么对齐成一张长表

3.1 原始数据的形状和第一道坎

赛题给的负荷数据是分区 × 逐小时的结构，温度是气象站 × 逐小时，两边的时区标注不一定完全一致。这是第一个容易出问题的地方：R 里as.POSIXct默认按系统时区解析，如果服务器是 UTC 而数据带夏令时偏移，每年的切换日就会多出一小时或者少一小时，整个分区错位 24 行。这种错位在验证集上可能看不出来，但在提交时会导致某一天的 24 个预测值整体偏了一格，pinball loss 直接变差很多。

稳妥的写法是显式指定时区，并且读完立刻做一次完整性检查：

library(data.table) library(lubridate) load_zone <- function(load_path, tz = "UTC") { dt <- fread(load_path) # 显式给时区，别依赖系统默认 dt[, timestamp := as.POSIXct(timestamp, tz = tz)] dt[, load := as.numeric(load)] # 完整性检查：逐小时序列不应该有缺口或重复 full <- seq(min(dt$timestamp), max(dt$timestamp), by = "hour") stopifnot(!anyDuplicated(dt$timestamp)) stopifnot(length(setdiff(full, dt$timestamp)) == 0) dt[order(timestamp)] }

stopifnot这两行是关键。anyDuplicated抓重复时戳，setdiff抓缺口。真实数据里缺口几乎必然存在，多半是传感器掉线或者采集系统重启，缺个几小时。直接删掉会让序列不连续，滞后特征跟着错位；填 0 更糟，负荷曲线会出现一排凹坑，模型学进去之后每逢类似季节就往下压。

我一般会按缺口长度分开处理：连续缺 1 到 2 小时的用线性插值，zoo::na.approx一行搞定；缺一整天以上的就标记成 NA 并在训练时整段跳过，不硬填。赛题数据里跨年的长缺口不多，但真碰上就别勉强。

3.2 气象站到分区的温度映射

温度是负荷预测里最重要的外生变量，但赛题只给气象站数据，不给站点和分区的对应关系。这是整个项目里最费心思的一步，也是很多人直接跳过、导致后面怎么调都上不去的地方。

思路是用训练期的相关性给每个分区选站。负荷对温度的响应在冬天是负相关（越冷负荷越高，取暖），夏天是正相关（越热负荷越高，制冷），所以要用绝对值排序，不能直接按相关系数从大到小排。

# 给每个分区挑相关性最强的 top_k 个站点 pick_stations <- function(load_dt, temp_dt, top_k = 3) { zones <- unique(load_dt$zone) rbindlist(lapply(zones, function(z) { l <- load_dt[zone == z, .(timestamp, load)] m <- merge(l, temp_dt, by = "timestamp") # 用绝对值，因为制冷和取暖的符号相反 cors <- m[, .(rho = cor(load, temp, use = "complete.obs")), by = station] setorder(cors, -abs(rho)) cors[seq_len(min(top_k, .N)), .(zone = z, station, rho)] })) }

top_k = 3是个经验值。取 1 个站容易受单点传感器故障影响，取太多又会把远处不相关站点的噪声带进来。选完之后对这几个站的温度做加权平均，权重可以用相关系数的绝对值归一化。

要注意的是，选站必须只用训练期数据。用全量数据算相关性，等于把未来信息泄进了特征构造，验证集分数会虚高。这个坑在时间序列比赛里很常见，稍不注意就会影响最终结果。

还有一层：温度本身也需要处理。原始温度有离群值，比如传感器在阳光直射下读数偏高，某几个小时突然飙到 45 度。这些点会严重干扰 hinge 特征。我一般会按分区做一次分位数截断，把 0.1% 和 99.9% 以外的值夹到边界上，不动中间的分布。

3.3 日历特征与温度衍生特征

特征构造这块，做了和做对差别很大。基础的是小时、星期几、月份这三个类别变量，R 里用factor或者直接展开成哑变量都行，树模型对 factor 处理得还可以，线性分位数回归就必须展开。

温度这块要做的是 hinge 特征，把非线性的温度响应掰成几条直线。pmax(0, temp - k)是标准写法，每个 knot 一个特征：

add_temp_hinge <- function(dt, knots = c(0, 8, 16, 22, 28)) { # 平滑温度，压掉逐小时的抖动 dt[, temp_smooth := frollmean(temp, n = 24, align = "right", fill = NA)] for (k in knots) { nm <- paste0("t_hinge_", k) dt[, (nm) := pmax(0, temp_smooth - k)] } dt[] }

knots 的选择要贴合当地气候。0/8/16/22/28 这组分段是按“严寒—冷—舒适—暖—酷热”切，覆盖了负荷曲线拐点最明显的位置。如果数据在热带，可以改成 18/22/26/30/34，或者干脆用分位数来定 knot，保证每段都有足够样本。

滞后特征要格外小心。shift(load, 24)这种 24 小时前负荷，在训练期是有的，在预测期第一天可能就没有。如果训练时全用了滞后项，模型会高度依赖它们，到了预测期一旦缺失，预测精度立刻下降。我一般采用一种折中策略：只保留“预测发起时刻已知”的特征，比如预测明天全天时，今天同一时刻的负荷是已知的，可以做 lag 24 的基准；但明天 12 点的 lag 24（即今天 12 点）也是已知的，所以逐小时递推其实可行。只是复杂度上去了，容易写成循环依赖。简单项目的做法是干脆不用负荷滞后，纯靠温度和日历外生变量，虽然损失一点精度，但工程上干净。

3.4 时间序列切分：验证集绝不能随机抽

这条单独拎出来说，因为它带来的影响最隐蔽。负荷数据有强自相关，随机切分会让验证集和训练集在时间上交错，模型见过“验证样本的邻居”，分数虚高十几个点是常事，一提交就暴露问题。

正确做法是最后一年做验证，前面若干年做训练；如果要调参，再从训练集里往前切一段做开发集。三段时间严格不重叠。

# 按时间切，绝不随机 split_by_time <- function(feat, valid_from = as.POSIXct("2011-01-01", tz = "UTC")) { list( train = feat[timestamp < valid_from], valid = feat[timestamp >= valid_from] ) }

顺带提一句，滚动预测场景下还要做“多步递推验证”：用训练集训模型，预测下周，再把预测值当输入往后再推一周。这会暴露很多单步验证看不出来的问题，比如误差累积、分位数越推越宽。赛题的待预测时段是分批下发的，每批之间模型需要重新适配，所以这一步不能省。

4. 模型层拆解：温度分段响应加 99 个分位数的提升树

4.1 用 hinge 特征把温度-负荷曲线掰弯

负荷和温度的关系是典型的 V 型（或者叫 U 型）：低温段随温度升高负荷下降，高温段随温度升高负荷上升，中间有个谷底。直接用原始温度做线性项，拟合出来是一条直线，两头全错。用多项式又容易在极端温度处发散，冬天零下十几度的时候预测出离谱的值。

hinge 特征的好处是把这条曲线拆成几段直线，每段一个斜率，连续但不光滑——这在负荷预测里反而合适，因为温度响应本身在拐点处就有突变（比如集中供暖开始的那条线）。线性分位数回归用 hinge 特征可以直接解释每个分段的斜率；树模型用 hinge 特征相当于帮它预先切好了分裂点候选，收敛更快。

这里有个取舍：knots 固定还是让树自己找。固定 knots 更稳，特别是在数据量小的分区；让树自己找更灵活，但样本少的时候容易过拟合到某几年的极端天气上。我的习惯是混合——固定一组粗 knots，再加一维原始温度让树自由发挥。

4.2 分位数回归与分位数 GBDT 的两条腿

两条腿的代码结构类似，都是把公式和数据喂进去，输出一批预测值。线性分位数回归更短：

library(quantreg) fit_quantile_lm <- function(train, formula, q) { rq(formula, data = train, tau = q, method = "fn") }

method = "fn"是 Frisch-Newton 算法，比默认的br在中等规模数据上快不少，尤其是分位数靠边（0.01 或 0.99）的时候，单纯形法容易慢。

GBDT 那边要注意参数对分位数目标的影响：

library(gbm) fit_quantile_gbm <- function(train, formula, q, n_trees = 800) { gbm( formula, data = train, distribution = list(name = "quantile", alpha = q), n.trees = n_trees, interaction.depth = 4, shrinkage = 0.03, bag.fraction = 0.7, n.minobsinnode = 20, verbose = FALSE ) }

shrinkage = 0.03配n.trees = 800是个偏保守的组合，适合分位数目标——分位数损失函数在尾部不太平滑，学习率大了容易在极端分位上报废。n.minobsinnode = 20是防止小分区在叶子上只有几个样本，尤其是 q=0.99 这种只有 1% 数据落在上面的情况。

bag.fraction = 0.7比默认的 0.5 高一点，是为了让每次迭代看到更多样本，减小分位数估计的方差。分位数回归本身比均值回归方差大，靠 bagging 压一压有用。

预测的时候记得传n.trees，否则gbm用的是默认值，你会得到一个完全不对的结果而且不一定报错：

predict(fit, newdata = test, n.trees = 800, type = "response")

4.3 分位数交叉和单调重排

99 个分位数独立训，几乎必然出现交叉：q=0.7 的预测值低于 q=0.6。理论上不可能，实际上因为每个模型是分开优化的，小事。

交叉不修，pinball loss 不会崩，但会浪费掉一部分信息，而且下游如果按分位数去算区间宽度会很别扭。简单的处理是逐样本做一次排序：

# 每行（一个时刻）的 99 个分位数强制单调不减 enforce_monotone <- function(pred_mat) { t(apply(pred_mat, 1, sort)) }

sort是最粗暴的重排方式，会把交叉的那两个分位数的值互换，等于承认“模型认为 0.7 分位应该更低”，然后强行掰回来。更讲究的做法是保序回归isoreg，最小化重排带来的损失，但在这个场景下差别不大，sort够用而且快。

要检查交叉有多严重，直接在评估脚本里统计：

crossing_rate <- function(pred_mat) { # pred_mat: 行是样本，列是按 q 升序排列的分位数 mean(rowSums(diff(t(pred_mat)) < 0) > 0) }

crossing_rate超过 5% 就该回头看看模型了，通常是某个分位数的训练样本太少，或者特征里有共线性太强的东西。

4.4 评估函数与提交格式必须同源

评估脚本和提交脚本用的分位数顺序必须是同一套，否则分数对不上还查不出原因。我的做法是只维护一个全局的QUANTILES常量，评估和提交都从它取：

QUANTILES <- seq(0.01, 0.99, by = 0.01) # 平均 pinball loss，跨所有分位数 mean_pinball <- function(y, pred_mat, qs = QUANTILES) { stopifnot(ncol(pred_mat) == length(qs)) mean(vapply(seq_along(qs), function(i) { pinball(y, pred_mat[, i], qs[i]) }, numeric(1))) } # 生成提交表：一行一个时刻，列是各分位数 to_submission <- function(ts, pred_mat, qs = QUANTILES) { out <- as.data.table(pred_mat) setnames(out, paste0("q", sprintf("%02d", round(qs * 100)))) out[, timestamp := ts] setcolorder(out, c("timestamp", setdiff(names(out), "timestamp"))) out[] }

stopifnot(ncol(pred_mat) == length(qs))这一行拦住过不止一次事故。列数对不上时报错，比生成一个列名错位的文件、传上去发现分数只有随机水平强。

5. 避坑与常见问题排查：五个让 pinball loss 变差的细节

5.1 现象一：本地验证很好，提交分数下降明显

原因几乎都是未来信息泄漏。最常见的三种：用全量数据算温度映射的相关系数；在特征里放了全局标准化时用到的均值和标准差；负荷滞后项跨过了预测起点。

解决方式是所有“需要拟合”的步骤都严格放在训练集内完成——选站、标准化参数、knots 的分位数、模型权重，一个都不能用验证期数据。我一般的做法是在run_all.R里把数据切成 train/valid 之后，所有拟合函数只接受 train 参数，验证集只进predict。

5.2 现象二：pinball loss 出现负数或者大得离谱

负数基本只可能来自一个地方：分位数顺序和列顺序对不上。QUANTILES是升序，但预测矩阵可能被某个apply或者t()转置过，导致 q=0.01 那一列拿到了 q=0.99 的预测值。这时候低估的惩罚变成 99 倍，loss 直接冲到天上。

解决办法是在mean_pinball里加一行断言，检查每个分位数列的均值是不是随 q 单调递增——正常情况下 q 越高预测值越大，这个性质很容易验证。

5.3 现象三：冬季负荷曲线整体塌陷

温度缺失值填 0 的典型后果。0 度在冬天算“不算太冷”，模型据此判断这是个温和的日子，负荷预测偏低。而实际那几小时可能是零下十几度，只是传感器掉线了。

正确做法是先对温度序列做缺失标记，短缺口插值，长缺口用同期其他站点的中位数补。补完之后再看一眼冬季月份的温度分布，如果补出来的值明显偏离历史同期，说明补法有问题。

5.4 现象四：夏令时切换那天预测整体错位

时间索引按 UTC 解析但原始数据带本地时区偏移，切换日会多出或少掉一小时。整条序列往后错一格，那一天的 24 个预测全错位。

排查办法很简单：对每个分区算一遍diff(timestamp)的分布，正常应该全是 1 小时，出现 0 或者 2 就说明有问题。这个检查放在load_zone里，比出错之后再回头找省事得多。

5.5 现象五：某个分区的分数就是上不去，其他分区正常

多半是这个分区的负荷形态特殊。比如以工业负荷为主的分区，周末和节假日负荷几乎不掉，用统一的日历特征会让模型学偏；或者这个分区只有一两个气象站覆盖，温度的代表性差。

处理方式是单独给这个分区加特征或者调参。节假日特征按分区的重要性加权，温度换一组站点。这种情况只能一个分区一个分区看残差图，没有通用解法。

6. 二次开发：换数据源、换模型、把评估脚本当回归测试用

这套代码的改造价值在于它的结构是“数据 → 特征 → 模型 → 评估”四段式，每一段都能单独替换。换成自家园区的负荷数据，主要动load_data.R和mapping.R两个文件，把列名和时区对上，其他不用碰。如果只有一个总表没有分区概念，把zone列固定成常量，run_all.R里的循环自然退化成一次执行。

换模型的时候有个接口约束要注意：新的模型函数必须接受(train, formula, q)三个参数并返回一个带predict方法的对象。守住这个约定，上层代码一行不用改。想试 LightGBM 的分位数目标，写个适配函数包一层就行。

改造之后怎么验证没改坏？我的习惯是把evaluate.R当成回归测试用：固定一份小样本，跑一遍，把mean_pinball的数值记下来，改完代码再跑一遍，差得超过一个阈值就回头查。这个方法比肉眼看预测曲线靠谱，因为曲线的变化往往很细微，而 loss 的数值变化是客观的。

还有一个更省事的技巧：在run_all.R里加个开关，只跑一个分区、三个分位数（0.1、0.5、0.9），几十秒就能出一版结果。改代码的过程中全程用这个快速模式，确认没问题再跑全量。

血泪经验是别在最后一天才跑全量。99 个分位数 × 20 个分区，光训练就能跑掉大半天，中间任何一个分区报错都要重来。从那以后我每次动模型代码，都强制走一遍“单分区三分位数”的快速验证，再走全量。希望帮到你。

本文还有配套的精品资源，点击获取