资讯详情

资讯详情

气象数据分析全流程:数据清洗、插值与趋势检验实战

简介面向数据分析入门者与Python爬虫学习者这份气象分析实战资料以中国天气网某一城市的真实天气数据为对象完整演示从爬虫采集、数据整理到可视化分析的闭环流程覆盖数据获取、清洗、图表呈现与结论输出等多个关键环节。源代码可用Jupyter Notebook直接打开运行配合整理好的实验报告读者能对照学习如何生成雷达图、条形图等常用图表并针对气温、风力、降水等维度写出有价值的分析说明快速建立数据分析项目的基本框架。压缩包共39个文件大小仅1.69MB其中包含20张PNG格式的分析图表以及15个XML和4个REL类型的文档结构文件——这些是实验报告Word文档的内部组成用于保存正文内容与图片关联关系便于理解报告排版机制。该资料已有1853人学习下载代码注释清晰、报告结构完整特别适合用作数据分析课程设计、爬虫可视化练手或实验报告撰写的参考模板。1. 气象分析数据分析到底在解决什么不是画曲线是处理时空错位气象分析数据分析有一个反直觉的结论这个方向最难的从来不是算法而是把数据“摆正”。A同学第一次接手某市一年的逐时气象观测资料画出温度曲线后以为拿到了完整时序可在每年1月1日零点都出现一个-80℃的怪点——那是仪器结冰后回传的乱码。更麻烦的是城市里明明有几十个站点能连续三年不缺测且质量达标的往往只剩三四个。气象数据分析真正要解决的是站点稀疏、时间错位、统计口径不一致这一串问题然后才能谈得上气候特征和预测模型。适合读这篇文章的是做气象服务、环境评估、农业保险、新能源功率预测的工程师和数据从业者以及想找一个有挑战性的数据分析练习对象的人。2. 数据选型与解析站点观测和再分析资料怎么影响你的分析结果气象分析数据分析的第一步不是跑模型而是选数据。我见过不少项目栽在数据源上有人拿着区域再分析资料去分析某个山谷的逐时低温有人拿站点观测去描述整个流域的降水分布结论自然站不住。常见做法是先按分析尺度选数据源再按数据源定清洗方案。这一章把两类数据源的差异讲透并给出一条最小可复现的解析路径。2.1 站点观测与再分析资料的五个差异维度站点观测来自真实仪器时间连续、数值可信但空间上只是“点”。同一城市的两个站热岛效应强的时候温差可达3-5℃。再分析资料是把历史观测同化进数值模式后生成的全球格点场空间连续、系统误差小但本质是模式输出局地的极端值会被网格平均抹平。选型时要看时间分辨率、空间覆盖、获取成本和分析尺度。下面这个表格是我常用的选型参照。差异维度站点观测再分析资料空间形态离散点依赖站点密度全球格点常用0.25°或0.5°时间连续性受仪器故障、迁站影响逐小时连续极少缺测极端值表现保留真实极值极端值被网格平均削弱数据体积小适合整段分析大需按区域裁剪下载适用场景单点时序、站点评估区域背景、缺测回填、驱动模型选型的核心判断是你要回答“这个点发生了什么”还是“这个区域发生了什么”。前者优先站点观测后者优先再分析资料。如果做单站缺测回填通常用再分析资料和邻近站做回归而不是直接填格点值。还有一个常见误用拿站点观测画区域等值线站点密度不够时插值结果会很好看但从气象上说不通。2.2 用Python解析逐时观测CSV从原始报文到整洁表格无论数据源是什么落地第一步都是把原始文件变成“一列时间戳、一列要素、一列质量标记”的长表或宽表。我一般先写一个最小解析脚本把原始CSV的列名对齐、单位统一、时间统一到UTC或本地时区。下面这段代码可以处理大多数逐时站点数据。import pandas as pd # 读取原始逐时观测列名按常见中文表头写实际以数据说明为准 df pd.read_csv( obs_hourly.csv, encodingutf-8, parse_dates[观测时间], # 直接解析成DatetimeIndex na_values[, 9999, -9999], # 把缺测标记统一读成NaN ) # 统一列名只留分析需要的字段 df df.rename(columns{ 观测时间: time, 气温(℃): temp, 相对湿度(%): rh, 风速(m/s): wind, 气压(hPa): pres, 降水量(mm): precip, })[[time, temp, rh, wind, pres, precip]] # 时间标准化排序、去重、检查间隔 df df.sort_values(time).drop_duplicates(subsettime).reset_index(dropTrue) df[temp] df[temp].astype(float) print(df.head())这段代码首先要做的是把“观测时间”转成pandas的datetime类型后面所有重采样和插值都依赖这一步。na_values是关键参数不同数据源的缺测编码不一样常见的有9999、-9999、32744最好先把原始文件里的缺测标记全部看一遍统一写进这个列表。排序和去重是为了防止原始文件里出现重复时次这在报文拼接后很常见。单位转换也需要在进入正式分析前完成比如气压从mmHg转成hPa、风速从km/h转成m/s一旦进入建模阶段再改单位特征和系数全要重解释。2.3 质量控制三个基本动作极值、时次、要素一致性拿到整洁表格不代表数据可信。气象数据的质量问题不像用户日志那样只表现为空值和异常大数更多是“合理但错误”。我常用的质量控制只有三件事。第一是极值校验。不要用固定的-60℃这类硬边界而要用该站历史分位数把逐时气温按月份分组超过99.9分位数或低于0.1分位数且有物理常识冲突的记录标出来人工确认后再决定删还是留。第二是时次校验。检查时间戳是否连续日均应1440条逐时记录缺时次要记录缺失位置而不是直接填充。第三是要素一致性校验。比如降水量大于0时相对湿度应接近饱和且气温不会太高如果某条降水记录对应的湿度只有20%那一定是仪器故障或报文错位。这三种校验做完后建议把每条被标记记录的原始值、校验规则、处理结果单独留一列而不是直接覆盖。气象分析里“可追踪”比“好看”重要得多后期做趋势和预测时任何一版清洗规则变更都能追溯结果差异来自哪里。3. 缺测与插值把气象数据洗到能上模型的程度站点观测最磨人的不是格式而是缺测。气象分析数据分析里有个说法缺测是常态连续完整才是意外。迁站、仪器维护、通信中断都会造成缺失而且缺测往往不是随机发生的——夏季强对流天气常常伴随通信故障于是“下雨时缺测”就成了系统性偏差直接分析会偏向晴天。这一章的处理原则只有一条先判断缺测能不能插再选择插法最后做质量验证。3.1 缺测率与缺失形态先判断能不能插再决定怎么插拿到缺测结果先算缺失率和缺失形态。常规判断流程缺失率小于5%且分布随机直接对邻近时间点做线性插值即可缺失率在5%到20%之间按要素选插法温度用时间插值或站点间插值都行降水不适合时间插值只能用邻站回归缺失率超过20%且呈季节性集中比如每年7月缺半个月就不要硬插了用再分析资料在该站位置做回归回填然后注明回填范围。缺失情形建议处理风险点5%随机零散线性插值或邻近时次复制插值会压低逐时方差5%-20%连续数小时温度用三次样条降水用邻站比例法样条在强对流时段可能过冲20%季节性集中再分析资料回归回填回填值不能参与极端事件统计站点完全断档数月不插按缺失去评估补出的序列不可信提示连续缺测超过72小时哪怕缺失率只有2%也不要直接线性插值。这里有一个容易翻车的判断连续缺测超过72小时哪怕缺失率只有2%也不要直接线性插值。72小时足够一次完整的天气过程过境线性插值会把一次降温过程生硬地拉直成缓慢过渡。遇到这种情况应该用邻近站点的同步观测做回归或者直接用再分析资料回填。原则是时间插补应对“短的随机缺失”空间回填应对“长的系统缺失”。3.2 时间重采样到日值三个聚合函数和一个闰年坑逐时数据整理干净后下一步往往是聚合成日值。温度取日均温、降水取日总量、风速取日平均看起来简单但有三个细节需要提前定清楚。第一个细节日均温到底怎么算。业务上常用的日平均有“四时次平均”和“24小时平均”两种。如果你手里的资料是逐时的直接用24小时平均更合理但做多年对比时要注意基线资料可能用的是四时次平均两者可能差0.5-1℃后面做趋势分析时要保持一致。第二个细节风向不能平均只能说“日主导风向”需要做16方位或8方位的频次统计这又是一个容易被人忽略的点。第三个细节是闰年下面是聚合代码。# 按日重采样先处理跨年与闰年 daily df.resample(D).agg( temp(temp, mean), # 24小时平均温度 temp_max(temp, max), # 日最高温 temp_min(temp, min), # 日最低温 precip(precip, sum), # 日降水量 wind(wind, mean), # 日平均风速 ).reset_index() # 降水量为0和缺测要区分开先把NaN保留下来 daily[precip] daily[precip].where(daily[precip].notna(), pd.NA)resample(D)会严格按照自然日聚合闰年多出来的一天会自动作为独立的一天参与统计不需要手工处理。需要注意这里precip直接sum的前提是原始数据中缺测已经被前面的步骤处理干净如果某一天缺了12个小时sum出来的日降水量是偏小的而不是缺测。所以要先把逐时缺测解决掉再聚合否则会出现日降水总量系统性偏低。聚合完成后最好再检查一遍日序的连续性确认没有生成空日期。处理闰年的正确姿势是信任DatetimeIndex不要自己写按年份循环的脚本那样反而容易漏掉2月29日。3.3 空间插值选型反距离加权和克里金的边界在哪里当分析要从点扩展到面时空间插值就躲不开了。站点少的时候反距离加权IDW通常比克里金稳因为克里金的变差函数拟合需要足够多的点对站点少于10个时拟合出来的变差函数基本靠猜。反过来站点密度高且分布均匀时克里金能把空间自相关结构利用起来插值更圆滑。下面给一个可以直接跑的IDW实现适合手头只有三五个站的情况。import numpy as np def idw_interp(obs_lon, obs_lat, obs_val, grid_lon, grid_lat, p2): 反距离加权插值obs为站点坐标和值grid为目标格点 out np.full(grid_lon.shape, np.nan) for i in range(grid_lon.shape[0]): for j in range(grid_lon.shape[1]): dist np.hypot(obs_lon - grid_lon[i, j], obs_lat - grid_lat[i, j]) if np.any(dist 1e-6): # 目标点与站点重合直接取站值 out[i, j] obs_val[np.argmin(dist)] continue w 1.0 / np.power(dist 1e-9, p) # 距离越近权重越大 out[i, j] np.sum(w * obs_val) / np.sum(w) return outIDW唯一要调的参数是幂指数pp2是默认选择适合大多数气象要素。p越大插值结果越接近最近站点的值极端值保留得好但表面不光滑p越小结果越接近区域均值空间平滑但会把降水极值磨平。对于降水这类偏态分布要素我一般先把降水量做开方或对数变换插值后再变换回去否则IDW会在强降水站四周画出很大的中雨区——这是最常见的插值翻车现场。克里金比IDW多了一步变差函数拟合参数里range和sill对结果影响很大站点少时调参基本是玄学不建议一上来就上克里金。4. 气候态、距平与趋势统计指标怎么选才不会主观气象分析数据分析做到这个阶段才开始真正回答业务问题今年冬天算不算冷、这个站的年雨量在减少吗。这两个问题看似简单统计口径选择不当会得出完全相反的结论。这一章讲两个核心工具气候态与距平、非参数趋势检验并说明各自避开统计陷阱的原理。4.1 气候态与距平定义一个冷冬需要哪几步“今年冬天偏冷”必须有参考基准。气象上常用最近三个年代的平均作为气候态比如当前常用的1991-2020年。计算距平就是用今年冬季均值减去气候态均值这里有一个必须处理的跨年细节冬季一般指12月到次年2月12月的年份与1月、2月不同必须先把冬季的年份标签对齐否则会把12月的值归错年份。下面代码是一个可直接套用的流程。winter df.set_index(time)[df[time].dt.month.isin([12, 1, 2])].copy() # 给冬季做“气候年”标签12月归到下一年 winter[season_year] winter.index.year winter.loc[winter.index.month 12, season_year] 1 # 计算每年冬季均值按season_year聚合 seasonal winter.groupby(season_year)[temp].mean() # 气候态取1991-2020年冬季均值的平均 clim seasonal.loc[1991:2020].mean() anomaly seasonal - clim # 距平序列 std_anomaly anomaly / seasonal.loc[1991:2020].std() # 标准化距平这段代码里最关键的逻辑是把12月归到下一个气候年12月的气温属于那一年的冬天但不是那一年的日历年season_year在12月加1其余月份不动。标准化距平是用气候态期间的标准差去除距平得到无量纲的数±1以内算正常超过±1.5才容易被认为是明显异常。温度这种近似正态的要素可以直接这样算降水是偏态分布标准化距平前要先做开方变换。这样定义出来的“冷冬”才经得起复查。4.2 用Mann-Kendall和Sen斜率检验趋势避开线性回归的过度自信判断年降水在减少很多人会画一条线性回归线然后看p值。这个做法对气象序列来说过于自信气象数据有很强的自相关今天的温度会影响明天的温度样本之间不独立普通线性回归的显著性检验会低估方差也就是更容易给出显著的结论。Mann-Kendall是气象分析里最常用的非参数趋势检验它只比较序列的相对次序不假设分布对异常值也更稳健。配合Sen斜率可以给出趋势的稳健估计不受个别极端年份影响。from scipy.stats import theilslopes import itertools import numpy as np # 计算Sen斜率每上升1年序列的变化量 res theilslopes(annual_series.values, annual_series.index.values) trend_slope res[0] # 简化的Mann-Kendall统计量S n len(annual_series) diff [annual_series.iloc[j] - annual_series.iloc[i] for i, j in itertools.combinations(range(n), 2)] S sum(np.sign(d) for d in diff)theilslopes返回的斜率是每一年序列的变化量比如斜率是-3.2mm/年说明年降水大约以每年3.2mm的速度在减少这个幅度比线性回归系数稳健。MK统计量S是所有两两比较中的“上升对”减去“下降对”S正值代表上升趋势。实际应用时不要只用S的符号判断趋势显著性还要计算方差并用标准正态近似得到p值序列有自相关时先把序列做预白化处理再检验否则p值依然会偏小。业务报告里我一般同时给MK的p值和Sen斜率前者回答“趋势是否存在”后者回答“趋势有多大”。统计误用方式后果正确做法用线性回归直接判趋势显著性虚高Mann-Kendall 预白化年际序列不做自相关检查p值不可信先算lag-1相关系数用日历年拼冬季序列冬季均值错位生成season_year标签直接平均降水年总量偏态分布下均值失真先变换或报告中位数5. 气象分析的五个高频排查点现象、成因、解决前面的流程走完大多数项目会进入“看起来能算”的状态。实战里我复盘过不少项目最花时间的不是建模和画图反而集中在五个反复出现的问题上。这些问题单独看都像数据细节但每一个都能让整条结论反转值得单独开一章讲。排查顺序建议按数据处理的先后走先查数据源和拼接再查清洗和插值最后查统计口径和绘图。5.1 站点迁移导致的温度断崖现象某站2015年前后温度序列出现约3℃的整体跳变同期邻站没有同向变化直接做距平和趋势分析都不对。原因气象站在2014年底迁址新址比旧址海拔高约200米两地气候背景不同但原始文件里站号没变。解决先查站址变更记录按记录把序列按时段拆开用两站重叠期做线性调整归一到新址再拼接。没有重叠期的用再分析资料在新旧两个位置分别取值做差值校正。5.2 插值把区域性暴雨磨平现象用IDW把某次暴雨过程的降水量插到格点降水中心从大于250mm变成了不足100mm的斑块。原因IDW对离群大值天然不友好它把周围小值站的权重分摊进来中心极值被压低了。解决对降水量做Box-Cox变换后插值再变换回来或者改用站点实测值画散点加阈值圈选。业务上如果目标是了解哪些地方超警戒宁可用实测站插值加阈值圈选也别强求平滑好看。5.3 日平均算法不一致制造虚假增温趋势现象同一个站一种资料的日均温比另一种高1℃左右拼在一起画长期趋势后出现明显增暖。原因早期资料用的是四时次平均后续资料用24小时平均夜里偏冷时段被剔除导致早期序列系统性偏高。解决向前看齐把早期时次资料重采样成24小时平均或者把整段序列统一成四时次平均再分析。报告里注明算法同一条序列内永远只用一种口径。5.4 自相关让Mann-Kendall检验误报显著现象降水序列MK检验p0.02看似显著但做了预白化后p0.13趋势不显著。原因年降水序列本身有年代际波动正自相关导致MK检验的方差被低估p值失真。解决先检验lag-1自相关系数大于0.1就先做预白化把序列减去lag-1回归后的残差再做MK。进一步可以用block bootstrap给出置信区间。这个问题的教训是不要只看一个p值趋势检验的结论必须能经得起自相关修正。5.5 忘了格点掩膜把无资料区画成降水大值区现象区域降水图上本该是干涸的地带出现漂亮的等值线和大范围强降水中心。原因再分析资料或插值结果在无观测区域仍会给出数值绘图时没加掩膜把无资料区当成真实低值区画出来了。解决下载数据时连同陆地掩膜一起取或读取通用的掩膜文件生成布尔数组绘图前把所有无资料点置为NaN再上色。这个小问题最常见的出现位置是区域边界外扩的一两圈检查图像时优先看四角。6. 从小时序列到温度预测一个最小建模流程的验证方法到这里数据和统计基础都具备了。气象分析数据分析的常见终点是预测明天温度多少、未来一小时有没有强降水。最后一个章节用一个温度预测小流程演示建模与验证的完整习惯不需要复杂模型但要把“验证”这道工序做扎实。6.1 特征工程构造滞后特征与滑动平均模型能学到什么取决于喂给它什么特征。温度预测常用的是滞后特征、滑动平均以及同期气压湿度。下面是把逐时序列构造成监督学习样本的常用方法。import numpy as np df[temp_lag1] df[temp].shift(1) # 上一小时温度 df[temp_lag24] df[temp].shift(24) # 昨天同一时刻温度 df[temp_sma3] df[temp].rolling(3).mean() # 近3小时滑动平均 df[hour] df[time].dt.hour # 周期性时间特征 df[sin_hour] np.sin(2 * np.pi * df[hour] / 24) df[cos_hour] np.cos(2 * np.pi * df[hour] / 24)滞后特征是时间序列预测里最重要的一组特征。temp_lag1抓状态的连续性temp_lag24抓住日循环规律滑动平均则平滑短时波动hour本身是数值型但0和23距离很远直接喂模型会误导树模型的分裂因此同时构造sin_hour和cos_hour把时间编码成圆环坐标。构造完后要dropna因为序列头部会因滞后产生空值。6.2 时间序列交叉验证随机切分会把验证变成自欺欺人训练和预测必须按时间顺序进行。把样本随机打散再切训练集和测试集模型会偷看未来数据验证指标虚高。常见做法是使用时间序列切分每次用前一段数据训练后一段数据验证滚动多次取平均误差。模型我用的是梯度提升树它不需要完备的特征缩放对缺失值和异常值也相对耐受。from sklearn.model_selection import TimeSeriesSplit from sklearn.ensemble import HistGradientBoostingRegressor tscv TimeSeriesSplit(n_splits5) for train_idx, val_idx in tscv.split(X): X_train, X_val X.iloc[train_idx], X.iloc[val_idx] y_train, y_val y.iloc[train_idx], y.iloc[val_idx] model HistGradientBoostingRegressor(max_iter200, learning_rate0.05) model.fit(X_train, y_train) pred model.predict(X_val) # 记录RMSE和MAE最后取5折平均TimeSeriesSplit默认按时间顺序取前段训练、后段验证不做随机洗牌。这里的关键参数是n_splits5代表划分成5组训练集越来越大验证集是固定的时间块模拟的是“用过去预测未来”的真实业务场景。如果换成train_test_split(random_state42)验证结果就会虚高这是时间序列建模里最常见的坑。特征重要性要在验证之后再解释且要小心滞后特征带来的虚假重要。6.3 把这个流程交给下一个人用配置文件管理数据源和参数分析能跑通只是第一步交付时真正值钱的是可复现性。我自己的习惯是数据读取、清洗、特征构造、训练、验证各自独立成函数所有数据源路径、聚合方式、阈值参数全部放进一个配置文件。下一个人接手时只改配置和路径不用改算法代码。这样做的底气来自一个血泪教训A同学曾因为临时改动插值参数没有记录重跑后结果和提交报告不一致最后花了三天排查才发现是清洗脚本里一处写死的缺测编码改掉了。希望这段经验能帮到你少走一次这个弯路。本文还有配套的精品资源点击获取
觉得有用,分享给同行:

为您的企业打造数字门面

稳重轻奢商务风格,端正雅致视觉,长效耐看不易过时。

立即咨询 →