最近整理欧洲地面气象数据时看到一个很有代表性的案例法国在 7 月同时刷新了高温与干旱纪录。这类新闻每年都会出现但如果只停留在“今年又热了”的层面就很难理解纪录究竟是怎样被打破的。是单站偶然偏高还是区域整体异常与气候平均状态相比偏了多少降水到底缺了多少这些问题都需要用数据来回答。本文将围绕“法国 7 月高温与干旱破纪录”这个场景用 Python 走一遍完整的气象数据分析流程从原始日值数据读取、清洗、月尺度聚合到计算 1991-2020 气候基准期距平、识别历史纪录、计算简化 SPI 干旱指数最后用 Matplotlib 输出可视化结果。文章适合有 Python 基础、想入门气象数据处理或者想积累时间序列分析经验的读者。1. 背景与核心概念1.1 什么是“破纪录”的气象数据“破纪录”听起来很直观但在气象统计中需要明确口径。媒体所说的“法国 7 月打破高温纪录”可以有几种理解方式是全国 7 月平均气温的历史最高值还是某一站点最高气温突破历史极值或者是整个 7 月的平均气温比过去所有 7 月都高不同指标对应的结论完全不同。在本文的数据分析流程中我们主要关注“月平均气温打破历史同期纪录”。即把每年 7 月的月平均气温排成一条时间序列若某一年的数值超过此前所有年份则该年 7 月就算“破纪录”。这种判断方法用累积最大值cummax就能实现逻辑简单也能体现气象纪录最核心的统计含义。同样干旱也需要量化。降水异常、土壤湿度、河流径流都是干旱的表现形式。本文使用最基础的降水序列先计算 7 月降水距平当月降水与气候平均值的差值再计算简化版的标准化降水指数SPI用来衡量某年 7 月到底有多干。1.2 高温与干旱如何量化距平与指数距平Anomaly是气象分析里最常用的量。它表示某个时刻的数值相对于气候平均状态的偏差气温距平 当月平均气温 - 气候基准期同月平均气温降水距平百分比当月降水 - 基准期同月平均降水 / 基准期同月平均降水 × 100%其中“气候基准期”通常使用世界气象组织推荐的 30 年标准气候平均值。过去常用 1961-1990后来更新为 1981-2010目前官方推荐使用 1991-2020。基准期不同距平结果也会不同这是分析时必须注意的地方。SPIStandardized Precipitation Index标准化降水指数则是把历史降水序列拟合成概率分布再转换成标准正态分布的分位数。SPI 为负表示偏干小于 -1 属于轻度干旱小于 -2 属于极端干旱。正式研究中 SPI 还需要处理“零降水”概率、按不同时间尺度1 个月、3 个月、6 个月分别计算本文只做教学演示。1.3 为什么用 Python 做气象数据分析气象数据本质上就是带时间戳的观测表格非常适合用 pandas 处理空间栅格数据可以用 xarray绘图用 Matplotlib 就能覆盖大部分需求。Python 的优势在于生态完整、脚本可复现、方便批量处理多个站点和多年数据。再加上 Jupyter Notebook 可以边写边看非常适合做这类探索性分析。2. 环境准备与数据来源2.1 运行环境与依赖库本文示例在以下环境验证操作系统Windows / Linux / macOS 均可Python 版本3.9 或 3.10 以上核心依赖pandas、numpy、scipy、matplotlib可选依赖netCDF4、xarray用于读取栅格数据版本不需要完全一致但建议使用较新的稳定版本。创建独立环境可以避免依赖冲突conda create -n climate python3.10 -y conda activate climate pip install pandas numpy scipy matplotlib netCDF4 xarray如果使用 pip 安装 netCDF4 失败可以改用 condaconda install -c conda-forge netcdf4 xarray2.2 数据来源说明真实的法国气象观测数据可以从以下渠道获取Météo-France 开放数据平台提供法国本土多个气象站的日值数据ECAD欧洲气候评估与数据集项目提供欧洲范围的站点数据和 E-OBS 栅格数据NOAA GHCN-D 全球历史气候网络数据E-OBS 是覆盖欧洲的观测插值网格数据分辨率约为 0.1° 和 0.25°适合做区域平均分析但下载和读取需要 xarray 与 netCDF4。本文为了让读者能直接复现先生成一个模拟站点日值数据再完成后续分析。真实数据下载后只需要把列名对齐成date、tmean、precip即可套用同一套代码。2.3 示例数据字段设计我们模拟法国南部一个气象站 1950 年至 2023 年的日平均气温与日降水量。字段如下字段含义单位date日期年-月-日tmean日平均气温°Cprecip日降水量mm这样设计的好处是贴近真实站点数据后面所有分析都围绕这 3 列展开。2.4 项目结构建议按下面结构组织文件climate-analysis/ ├── data/ │ └── france_south_daily.csv ├── scripts/ │ ├── 01_generate_sample_data.py │ ├── 02_monthly_aggregate.py │ ├── 03_record_and_spi.py │ └── 04_plot.py └── output/代码先按脚本拆分最终也可以合并成一个脚本运行。3. 数据读取与预处理3.1 生成模拟数据为了让大家零依赖复现这里先用固定随机种子生成示例数据。随机种子固定后每次运行得到的结果一致方便对照。# 文件路径scripts/01_generate_sample_data.py import numpy as np import pandas as pd np.random.seed(42) # 生成 1950-01-01 到 2023-12-31 的日期序列 dates pd.date_range(1950-01-01, 2023-12-31, freqD) # 用正弦函数模拟季节温度变化再叠加随机波动 tmean 15 10 * np.sin(2 * np.pi * dates.dayofyear / 365.25) np.random.normal(0, 2, len(dates)) # 使用 Gamma 分布模拟日降水夏季略少 precip np.random.gamma(2, 3, len(dates)) precip[dates.month.isin([6, 7, 8])] * 0.6 df pd.DataFrame({ date: dates, tmean: np.round(tmean, 2), precip: np.round(precip, 1) }) df.to_csv(data/france_south_daily.csv, indexFalse) print(df.head()) print(df.info())运行后可以看到前 5 行数据和字段信息。这个模拟序列中温度存在明显年周期夏季平均在 25°C 左右冬季 5°C 左右符合法国南部地中海气候的大致特征。3.2 读取 CSV 并解析日期读取时最核心的是日期解析。如果日期列不是标准的YYYY-MM-DD一定要单独处理。import pandas as pd df pd.read_csv(data/france_south_daily.csv, parse_dates[date]) df[year] df[date].dt.year df[month] df[date].dt.month df[day] df[date].dt.day print(df.head()) print(df[date].min(), df[date].max())parse_dates[date]会让 pandas 自动将这一列转换为datetime64类型。之后用dt.year、dt.month、dt.day提取年、月、日能极大方便后续分组聚合。注意数据若包含 2 月 29 日在闰年处理上是正常的但如果原始文件里的日期格式是29/02/2023这种不符合实际的日期pandas 会报错此时需要检查数据源而不是盲目跳过。3.3 缺失值与异常值检查真实气象数据几乎不可能没有缺失值。处理原则是先统计缺失比例再决定填充还是删除。对于日值数据最稳妥的方式是缺多少就标记多少聚合到月尺度时设置阈值。print(缺失值统计) print(df.isna().sum()) # 去除气温明显不合理的记录例如超过 -50°C 或低于 -100°C 这类物理异常 df df[(df[tmean] -50) (df[tmean] 50)] # 降水量非负 df df[df[precip] 0]这里只是最基本的物理范围检查。正式的科研流程还会检查站点迁移、仪器更换、相邻站点一致性等问题本文不展开。3.4 月尺度聚合日值数据噪音大判断“某月破纪录”通常使用月平均值。聚合逻辑是气温取月平均降水取月总和。monthly df.groupby([year, month]).agg( tmean_month(tmean, mean), precip_month(precip, sum), days_count(tmean, count) ).reset_index() # 如果某月有效观测天数太少直接舍弃避免月平均值失真 monthly monthly[monthly[days_count] 25] print(monthly.head()) print(monthly.tail())days_count 25是一种简单的质量控制某月如果缺测超过 5 天该月的平均气温就不具备代表性。实际项目中阈值可以设为 28 天甚至 29 天具体取决于数据质量要求。4. 核心分析气温距平与破纪录判断4.1 计算气候基准期平均值先取出 1991-2020 年的数据按月份计算 30 年平均气温。这个值就是“气候平均值”也是距平对比的基准。clim ( monthly[monthly[year].between(1991, 2020)] .groupby(month)[tmean_month] .mean() .rename(tmean_clim) ) print(clim)注意between(1991, 2020)是闭区间包含 1991 年和 2020 年。如果用官方标准气候值应以 1991 年 1 月 1 日到 2020 年 12 月 31 日完整 30 年为准。4.2 合并基准值并计算距平把气候平均值合并回原表然后用“当月气温 - 当月气候平均值”得到距平。monthly monthly.merge(clim, onmonth, howleft) monthly[tmean_anom] monthly[tmean_month] - monthly[tmean_clim] # 只看 1991 年以后避免基准期内数据自身参与对比造成误导 recent monthly[monthly[year] 1991] print(recent.tail())距平为正说明比气候平均热为负说明偏凉。这条时间序列可以直接用来观察“7 月高温破纪录”到底偏高了多久。4.3 识别 7 月气温破纪录年份破纪录的判断逻辑很朴素当前年份的月气温是否超过此前所有年份的月气温。用 pandas 的cummax累积最大值一条语句就能实现。july monthly[monthly[month] 7].sort_values(year).copy() july[cummax] july[tmean_month].cummax() july[is_record] july[tmean_month] july[cummax] record_years july[july[is_record]] print(record_years[[year, tmean_month, cummax, tmean_anom]])输出结果会显示哪些年份的 7 月刷新了历史纪录。由于示例数据是模拟的实际年份不一定对应真实新闻但整套判断逻辑与真实研究一致。真正的法国观测数据中近年如 2018、2019、2022 年7 月多次刷新纪录这正是标题所述“破纪录”在数据层面的体现。4.4 破纪录判断的注意事项这里要提醒几点cummax从序列开始位置计算。如果数据从 1950 年开始那么 1950 年本身就是纪录这是历史序列的自然结果不代表异常。站点资料长度不同会影响“历史纪录”的判定。一个 1950 年开始观测的站点和一个 2000 年开始观测的站点可比性完全不同。破纪录与“异常程度”不是一回事。某年打破了纪录但可能只比上一年高 0.01°C而有些年份没有破纪录距平却非常显著。因此破纪录判断通常要和距平、百分位放在一起看。5. 干旱分析降水距平与简化 SPI5.1 降水距平计算干旱分析的第一步同样是看降水序列。这里以 7 月为例先算降水距平百分比。clim_precip ( monthly[monthly[year].between(1991, 2020)] .groupby(month)[precip_month] .mean() .rename(precip_clim) ) monthly monthly.merge(clim_precip, onmonth, howleft) monthly[precip_anom_pct] (monthly[precip_month] - monthly[precip_clim]) / monthly[precip_clim] * 100 july_precip monthly[monthly[month] 7] print(july_precip.tail())降水距平为负表示偏旱。百分比形式比绝对差值更容易理解-50% 意味着当月降水只有气候平均的一半。5.2 用 Gamma 分布计算简化 SPISPI 的基本思路把某月的历史降水序列拟合为 Gamma 分布然后求累积概率再转换成标准正态分布的分位数。SPI 为 0 表示接近中位水平负值越小越干旱。from scipy import stats import numpy as np # 取完整历史 7 月降水序列 july_rain july_precip[precip_month].values # 过滤掉全零月份避免分布拟合失败 positive july_rain[july_rain 0] # 用极大似然估计拟合 Gamma 分布floc0 表示位置参数固定为 0 alpha, loc, beta stats.gamma.fit(positive, floc0) # 计算每个 7 月降水的累积概率 p stats.gamma.cdf(july_rain, aalpha, scalebeta) # 映射到标准正态分位数 spi_july stats.norm.ppf(np.clip(p, 1e-6, 1 - 1e-6)) july_spi july_precip[[year, precip_month]].copy() july_spi[spi] np.round(spi_july, 2) print(july_spi.tail(10))这段代码是教学简化版。正规的 SPI 计算还需要处理“零降水概率”问题即降水为 0 的月份不能直接参与 Gamma 拟合而要用混合分布处理。实际科研可以直接使用SPEIPython 包或climate_indices库。5.3 SPI 等级划分与解读SPI 常用等级如下SPI 值干旱等级0 到 -0.99轻度干旱-1.00 到 -1.49中度干旱-1.50 到 -1.99严重干旱小于 -2.00极端干旱如果某年 7 月 SPI 达到 -1.5 以下同时气温距平显著为正就会出现“高温叠加干旱”的复合极端事件。这正是标题里“heat and drought records”同时出现的数据表现。需要强调SPI 只反映降水偏差不反映温度对干旱的加剧作用。要更全面地评估农业干旱和水文干旱应进一步使用标准化降水蒸散指数SPEI它需要气温、潜在蒸散发等数据。这篇文章先走到 SPI 这一层已经足够说明分析方法。6. 可视化与结果解读6.1 绘制 7 月气温距平时间序列可视化能直观展示“破纪录”的过程。先画一条从 1950 年至今的 7 月气温距平曲线并标注最近年份。# 文件路径scripts/04_plot.py import matplotlib.pyplot as plt import pandas as pd july pd.read_csv(output/july_monthly.csv, parse_dates[date]) fig, ax plt.subplots(figsize(12, 5)) ax.plot(july[year], july[tmean_anom], color#d62728, lw1.2, label7月气温距平) ax.axhline(0, colorblack, lw0.8, ls--) ax.set_xlabel(年份) ax.set_ylabel(气温距平 (°C)) ax.set_title(法国南部某站 7 月气温距平基准期 1991-2020) ax.legend() fig.tight_layout() plt.savefig(output/july_temp_anomaly.png, dpi150) plt.show()从图中可以快速判断近十年往往出现较长一段正距平正距平峰值年份大概率就是“破纪录”年份。6.2 用双轴图叠加降水与气温高温和干旱经常同时出现把 7 月降水距平百分比和气温距平画在同一个图里能直观看出两者是否同步。fig, ax1 plt.subplots(figsize(12, 6)) ax1.bar(july[year], july[precip_anom_pct], color#9ecae1, alpha0.7, label7月降水距平百分比) ax1.set_xlabel(年份) ax1.set_ylabel(降水距平 (%), color#3182bd) ax1.tick_params(axisy, labelcolor#3182bd) ax2 ax1.twinx() ax2.plot(july[year], july[tmean_anom], color#d62728, markero, ms3, label7月气温距平) ax2.axhline(0, colorgray, lw0.8, ls--) ax2.set_ylabel(气温距平 (°C), color#d62728) ax2.tick_params(axisy, labelcolor#d62728) plt.title(法国南部某站 7 月降水距平与气温距平) fig.tight_layout() plt.savefig(output/july_precip_temp.png, dpi150) plt.show()图中如果出现“降水柱状图明显偏低 气温曲线明显偏高”的年份基本就是复合极端事件的候选年份。6.3 柱状图突出破纪录年份也可以只关注 7 月平均气温本身用柱状图标出每个破纪录年份效果更接近新闻报道里的图表。colors [#d62728 if r else #bbbbbb for r in july[is_record]] fig, ax plt.subplots(figsize(12, 5)) ax.bar(july[year], july[tmean_month], colorcolors) ax.set_xlabel(年份) ax.set_ylabel(7月平均气温 (°C)) ax.set_title(法国南部某站 7 月平均气温红色为破纪录年份) fig.tight_layout() plt.savefig(output/july_record_bar.png, dpi150) plt.show()红色柱子出现的频率能直观反映“纪录被不断刷新”的节奏如果破纪录年份主要出现在近期说明近期 7 月确实频繁偏热。6.4 结果解读的边界可视化只是分析工具解读时要注意模拟数据的结果不能直接等同于真实法国气候结论真实分析必须使用官方观测数据。单个站点的结果不能代表整个法国。法国南北气候差异很大巴黎、马赛、布列塔尼的 7 月气温和降水特征完全不同。“破纪录”反映的是统计事实不自动等于“气候变化导致的”归因需要单独的事件归因研究。本文的范围是描述性统计和可视化。7. 常见问题与排查思路问题现象常见原因解决思路读取 CSV 报ParserError: day is out of range for month日期格式不规范或数据包含非法日期打印原始日期列排查统一为YYYY-MM-DD必要时用errorscoercegroupby().mean()结果全为 NaN列类型还不是数值类型或缺失值未处理用pd.to_numeric(..., errorscoerce)转换先处理缺失值ImportError: No module named netCDF4依赖未安装或安装到了错误环境确认先conda activate climate再conda install -c conda-forge netcdf4stats.gamma.fit返回 nan数据全为 0 或样本量不足过滤正降水加入最小样本数判断样本过少时改用经验分布距平符号和预期相反单位或基准期错误确认气温是摄氏度而不是开尔文确认基准期是 1991-2020matplotlib 中文乱码系统缺少中文字体设置plt.rcParams[font.sans-serif] [SimHei]或改用英文标签排查时建议从数据源头开始先看原始数据前 100 行确认日期格式、单位、缺失值表现再逐步向后排查。不要一上来就怀疑算法。8. 最佳实践与工程建议8.1 数据管理规范气象分析项目建议统一命名规则原始观测数据放在data/raw/不做任何修改清洗后的中间数据放在data/processed/图表输出到output/每个处理步骤单独存成脚本运行顺序用文件名编号控制这样做的最大好处是几个月后回看项目时能清楚知道每个文件是怎么来的也能快速定位某个结果对应哪一步处理。8.2 合理选择气候基准期世界气象组织建议使用连续的 30 年作为气候标准期当前官方推荐期是 1991-2020。选择基准期时要注意基准期不能和研究时段重叠过多否则“距平”和“破纪录”会失真不同文献使用的基准期不同比较数值前先确认口径如果是长期趋势研究还可以同时计算多个基准期结果观察稳定性8.3 单位、缺失值与质量控制气象数据最容易出错的地方就是单位温度可能是摄氏度也可能是开尔文降水可能是毫米也可能是英寸。项目一开始就要在元数据中写明单位并在代码里增加单位校验。月尺度聚合前务必设置有效观测天数阈值。一个月的观测如果缺了 15 天算出来的月平均气温几乎没有参考价值。质量控制规则应当写在文档里并随数据一起保存。8.4 破纪录与归因结论要谨慎分析师很容易把“破纪录”直接写成“气候变化导致”这在科学上是不严谨的。破纪录是统计事实归因需要回答“如果没有人类活动影响这个事件还会不会发生”这属于事件归因研究通常需要气候模式模拟。在业务报告或博文中建议的表述方式是可以写“该年 7 月平均气温为 26.4°C比 1991-2020 平均值高 1.8°C为有记录以来最高。”不要写“气候变化导致今年 7 月高温破纪录。”保持统计结论和因果推断的边界是数据分析的基本职业素养。8.5 可复现性可复现性体现在三个层面代码层面固定随机种子如np.random.seed(42)保证模拟数据一致环境层面用requirements.txt或 conda env 记录依赖版本数据层面真实数据保存下载日期、版本号、来源链接这样即使几个月后重新运行也能得到相同结果或者能快速排查出结果变化的原因。8.6 大数据量下的性能优化站点日值数据通常只有几十 MBpandas 足够处理。但如果换成 E-OBS 这类覆盖全欧洲、时间跨度几十年的栅格数据建议使用 xarray 读取 netCDF利用维度标签自动对齐数据对长时间序列做区域平均时先用sel裁剪研究区域内存不足时使用chunk和 dask 延迟计算避免在循环里逐日读写文件9. 总结与拓展通过这篇文章我们完成了一条完整的气象数据分析路径从模拟日值数据出发经过日期解析、缺失值检查、月尺度聚合计算了 1991-2020 基准期的气温距平和降水距平用累积最大值识别出 7 月破纪录年份并用简化 Gamma 分布计算了 SPI 干旱指数最后用 Matplotlib 画出了三种不同类型的可视化图。这套流程稍加修改就能推广到其他场景把月份从 7 月改成任意月份分析“史上最热 X 月”是否成立把单站改成多站用groupby([station, year, month])做区域对比把日值换成小时值分析“最热一天”或“最长热浪持续天数”用 xarray 直接读取 E-OBS 栅格数据做空间分布分析如果继续深入可以学习 SPEI 干旱指数、热浪识别算法HWMId、气候趋势显著性检验Mann-Kendall 检验等内容。这些方向都需要本文的基础功底数据清洗要扎实、统计口径要清楚、结论边界要把握住。建议你下载一份真实的法国或中国气象站数据把本文代码跑一遍再试着改成自己感兴趣的城市。数据分析和写代码一样只有亲手处理过缺失值、亲手画过破纪录年份的柱状图才能真正理解那些新闻标题背后的统计细节。如果本文对你有帮助可以收藏备用后续遇到气象数据处理问题也能随时回来查阅。 
长河君 · 2026/8/28 14:33:01 · 阅读 1.2k 用Python分析气象数据：从距平计算到破纪录判断与SPI干旱指数
最近整理欧洲地面气象数据时看到一个很有代表性的案例法国在 7 月同时刷新了高温与干旱纪录。这类新闻每年都会出现但如果只停留在“今年又热了”的层面就很难理解纪录究竟是怎样被打破的。是单站偶然偏高还是区域整体异常与气候平均状态相比偏了多少降水到底缺了多少这些问题都需要用数据来回答。本文将围绕“法国 7 月高温与干旱破纪录”这个场景用 Python 走一遍完整的气象数据分析流程从原始日值数据读取、清洗、月尺度聚合到计算 1991-2020 气候基准期距平、识别历史纪录、计算简化 SPI 干旱指数最后用 Matplotlib 输出可视化结果。文章适合有 Python 基础、想入门气象数据处理或者想积累时间序列分析经验的读者。1. 背景与核心概念1.1 什么是“破纪录”的气象数据“破纪录”听起来很直观但在气象统计中需要明确口径。媒体所说的“法国 7 月打破高温纪录”可以有几种理解方式是全国 7 月平均气温的历史最高值还是某一站点最高气温突破历史极值或者是整个 7 月的平均气温比过去所有 7 月都高不同指标对应的结论完全不同。在本文的数据分析流程中我们主要关注“月平均气温打破历史同期纪录”。即把每年 7 月的月平均气温排成一条时间序列若某一年的数值超过此前所有年份则该年 7 月就算“破纪录”。这种判断方法用累积最大值cummax就能实现逻辑简单也能体现气象纪录最核心的统计含义。同样干旱也需要量化。降水异常、土壤湿度、河流径流都是干旱的表现形式。本文使用最基础的降水序列先计算 7 月降水距平当月降水与气候平均值的差值再计算简化版的标准化降水指数SPI用来衡量某年 7 月到底有多干。1.2 高温与干旱如何量化距平与指数距平Anomaly是气象分析里最常用的量。它表示某个时刻的数值相对于气候平均状态的偏差气温距平 当月平均气温 - 气候基准期同月平均气温降水距平百分比当月降水 - 基准期同月平均降水 / 基准期同月平均降水 × 100%其中“气候基准期”通常使用世界气象组织推荐的 30 年标准气候平均值。过去常用 1961-1990后来更新为 1981-2010目前官方推荐使用 1991-2020。基准期不同距平结果也会不同这是分析时必须注意的地方。SPIStandardized Precipitation Index标准化降水指数则是把历史降水序列拟合成概率分布再转换成标准正态分布的分位数。SPI 为负表示偏干小于 -1 属于轻度干旱小于 -2 属于极端干旱。正式研究中 SPI 还需要处理“零降水”概率、按不同时间尺度1 个月、3 个月、6 个月分别计算本文只做教学演示。1.3 为什么用 Python 做气象数据分析气象数据本质上就是带时间戳的观测表格非常适合用 pandas 处理空间栅格数据可以用 xarray绘图用 Matplotlib 就能覆盖大部分需求。Python 的优势在于生态完整、脚本可复现、方便批量处理多个站点和多年数据。再加上 Jupyter Notebook 可以边写边看非常适合做这类探索性分析。2. 环境准备与数据来源2.1 运行环境与依赖库本文示例在以下环境验证操作系统Windows / Linux / macOS 均可Python 版本3.9 或 3.10 以上核心依赖pandas、numpy、scipy、matplotlib可选依赖netCDF4、xarray用于读取栅格数据版本不需要完全一致但建议使用较新的稳定版本。创建独立环境可以避免依赖冲突conda create -n climate python3.10 -y conda activate climate pip install pandas numpy scipy matplotlib netCDF4 xarray如果使用 pip 安装 netCDF4 失败可以改用 condaconda install -c conda-forge netcdf4 xarray2.2 数据来源说明真实的法国气象观测数据可以从以下渠道获取Météo-France 开放数据平台提供法国本土多个气象站的日值数据ECAD欧洲气候评估与数据集项目提供欧洲范围的站点数据和 E-OBS 栅格数据NOAA GHCN-D 全球历史气候网络数据E-OBS 是覆盖欧洲的观测插值网格数据分辨率约为 0.1° 和 0.25°适合做区域平均分析但下载和读取需要 xarray 与 netCDF4。本文为了让读者能直接复现先生成一个模拟站点日值数据再完成后续分析。真实数据下载后只需要把列名对齐成date、tmean、precip即可套用同一套代码。2.3 示例数据字段设计我们模拟法国南部一个气象站 1950 年至 2023 年的日平均气温与日降水量。字段如下字段含义单位date日期年-月-日tmean日平均气温°Cprecip日降水量mm这样设计的好处是贴近真实站点数据后面所有分析都围绕这 3 列展开。2.4 项目结构建议按下面结构组织文件climate-analysis/ ├── data/ │ └── france_south_daily.csv ├── scripts/ │ ├── 01_generate_sample_data.py │ ├── 02_monthly_aggregate.py │ ├── 03_record_and_spi.py │ └── 04_plot.py └── output/代码先按脚本拆分最终也可以合并成一个脚本运行。3. 数据读取与预处理3.1 生成模拟数据为了让大家零依赖复现这里先用固定随机种子生成示例数据。随机种子固定后每次运行得到的结果一致方便对照。# 文件路径scripts/01_generate_sample_data.py import numpy as np import pandas as pd np.random.seed(42) # 生成 1950-01-01 到 2023-12-31 的日期序列 dates pd.date_range(1950-01-01, 2023-12-31, freqD) # 用正弦函数模拟季节温度变化再叠加随机波动 tmean 15 10 * np.sin(2 * np.pi * dates.dayofyear / 365.25) np.random.normal(0, 2, len(dates)) # 使用 Gamma 分布模拟日降水夏季略少 precip np.random.gamma(2, 3, len(dates)) precip[dates.month.isin([6, 7, 8])] * 0.6 df pd.DataFrame({ date: dates, tmean: np.round(tmean, 2), precip: np.round(precip, 1) }) df.to_csv(data/france_south_daily.csv, indexFalse) print(df.head()) print(df.info())运行后可以看到前 5 行数据和字段信息。这个模拟序列中温度存在明显年周期夏季平均在 25°C 左右冬季 5°C 左右符合法国南部地中海气候的大致特征。3.2 读取 CSV 并解析日期读取时最核心的是日期解析。如果日期列不是标准的YYYY-MM-DD一定要单独处理。import pandas as pd df pd.read_csv(data/france_south_daily.csv, parse_dates[date]) df[year] df[date].dt.year df[month] df[date].dt.month df[day] df[date].dt.day print(df.head()) print(df[date].min(), df[date].max())parse_dates[date]会让 pandas 自动将这一列转换为datetime64类型。之后用dt.year、dt.month、dt.day提取年、月、日能极大方便后续分组聚合。注意数据若包含 2 月 29 日在闰年处理上是正常的但如果原始文件里的日期格式是29/02/2023这种不符合实际的日期pandas 会报错此时需要检查数据源而不是盲目跳过。3.3 缺失值与异常值检查真实气象数据几乎不可能没有缺失值。处理原则是先统计缺失比例再决定填充还是删除。对于日值数据最稳妥的方式是缺多少就标记多少聚合到月尺度时设置阈值。print(缺失值统计) print(df.isna().sum()) # 去除气温明显不合理的记录例如超过 -50°C 或低于 -100°C 这类物理异常 df df[(df[tmean] -50) (df[tmean] 50)] # 降水量非负 df df[df[precip] 0]这里只是最基本的物理范围检查。正式的科研流程还会检查站点迁移、仪器更换、相邻站点一致性等问题本文不展开。3.4 月尺度聚合日值数据噪音大判断“某月破纪录”通常使用月平均值。聚合逻辑是气温取月平均降水取月总和。monthly df.groupby([year, month]).agg( tmean_month(tmean, mean), precip_month(precip, sum), days_count(tmean, count) ).reset_index() # 如果某月有效观测天数太少直接舍弃避免月平均值失真 monthly monthly[monthly[days_count] 25] print(monthly.head()) print(monthly.tail())days_count 25是一种简单的质量控制某月如果缺测超过 5 天该月的平均气温就不具备代表性。实际项目中阈值可以设为 28 天甚至 29 天具体取决于数据质量要求。4. 核心分析气温距平与破纪录判断4.1 计算气候基准期平均值先取出 1991-2020 年的数据按月份计算 30 年平均气温。这个值就是“气候平均值”也是距平对比的基准。clim ( monthly[monthly[year].between(1991, 2020)] .groupby(month)[tmean_month] .mean() .rename(tmean_clim) ) print(clim)注意between(1991, 2020)是闭区间包含 1991 年和 2020 年。如果用官方标准气候值应以 1991 年 1 月 1 日到 2020 年 12 月 31 日完整 30 年为准。4.2 合并基准值并计算距平把气候平均值合并回原表然后用“当月气温 - 当月气候平均值”得到距平。monthly monthly.merge(clim, onmonth, howleft) monthly[tmean_anom] monthly[tmean_month] - monthly[tmean_clim] # 只看 1991 年以后避免基准期内数据自身参与对比造成误导 recent monthly[monthly[year] 1991] print(recent.tail())距平为正说明比气候平均热为负说明偏凉。这条时间序列可以直接用来观察“7 月高温破纪录”到底偏高了多久。4.3 识别 7 月气温破纪录年份破纪录的判断逻辑很朴素当前年份的月气温是否超过此前所有年份的月气温。用 pandas 的cummax累积最大值一条语句就能实现。july monthly[monthly[month] 7].sort_values(year).copy() july[cummax] july[tmean_month].cummax() july[is_record] july[tmean_month] july[cummax] record_years july[july[is_record]] print(record_years[[year, tmean_month, cummax, tmean_anom]])输出结果会显示哪些年份的 7 月刷新了历史纪录。由于示例数据是模拟的实际年份不一定对应真实新闻但整套判断逻辑与真实研究一致。真正的法国观测数据中近年如 2018、2019、2022 年7 月多次刷新纪录这正是标题所述“破纪录”在数据层面的体现。4.4 破纪录判断的注意事项这里要提醒几点cummax从序列开始位置计算。如果数据从 1950 年开始那么 1950 年本身就是纪录这是历史序列的自然结果不代表异常。站点资料长度不同会影响“历史纪录”的判定。一个 1950 年开始观测的站点和一个 2000 年开始观测的站点可比性完全不同。破纪录与“异常程度”不是一回事。某年打破了纪录但可能只比上一年高 0.01°C而有些年份没有破纪录距平却非常显著。因此破纪录判断通常要和距平、百分位放在一起看。5. 干旱分析降水距平与简化 SPI5.1 降水距平计算干旱分析的第一步同样是看降水序列。这里以 7 月为例先算降水距平百分比。clim_precip ( monthly[monthly[year].between(1991, 2020)] .groupby(month)[precip_month] .mean() .rename(precip_clim) ) monthly monthly.merge(clim_precip, onmonth, howleft) monthly[precip_anom_pct] (monthly[precip_month] - monthly[precip_clim]) / monthly[precip_clim] * 100 july_precip monthly[monthly[month] 7] print(july_precip.tail())降水距平为负表示偏旱。百分比形式比绝对差值更容易理解-50% 意味着当月降水只有气候平均的一半。5.2 用 Gamma 分布计算简化 SPISPI 的基本思路把某月的历史降水序列拟合为 Gamma 分布然后求累积概率再转换成标准正态分布的分位数。SPI 为 0 表示接近中位水平负值越小越干旱。from scipy import stats import numpy as np # 取完整历史 7 月降水序列 july_rain july_precip[precip_month].values # 过滤掉全零月份避免分布拟合失败 positive july_rain[july_rain 0] # 用极大似然估计拟合 Gamma 分布floc0 表示位置参数固定为 0 alpha, loc, beta stats.gamma.fit(positive, floc0) # 计算每个 7 月降水的累积概率 p stats.gamma.cdf(july_rain, aalpha, scalebeta) # 映射到标准正态分位数 spi_july stats.norm.ppf(np.clip(p, 1e-6, 1 - 1e-6)) july_spi july_precip[[year, precip_month]].copy() july_spi[spi] np.round(spi_july, 2) print(july_spi.tail(10))这段代码是教学简化版。正规的 SPI 计算还需要处理“零降水概率”问题即降水为 0 的月份不能直接参与 Gamma 拟合而要用混合分布处理。实际科研可以直接使用SPEIPython 包或climate_indices库。5.3 SPI 等级划分与解读SPI 常用等级如下SPI 值干旱等级0 到 -0.99轻度干旱-1.00 到 -1.49中度干旱-1.50 到 -1.99严重干旱小于 -2.00极端干旱如果某年 7 月 SPI 达到 -1.5 以下同时气温距平显著为正就会出现“高温叠加干旱”的复合极端事件。这正是标题里“heat and drought records”同时出现的数据表现。需要强调SPI 只反映降水偏差不反映温度对干旱的加剧作用。要更全面地评估农业干旱和水文干旱应进一步使用标准化降水蒸散指数SPEI它需要气温、潜在蒸散发等数据。这篇文章先走到 SPI 这一层已经足够说明分析方法。6. 可视化与结果解读6.1 绘制 7 月气温距平时间序列可视化能直观展示“破纪录”的过程。先画一条从 1950 年至今的 7 月气温距平曲线并标注最近年份。# 文件路径scripts/04_plot.py import matplotlib.pyplot as plt import pandas as pd july pd.read_csv(output/july_monthly.csv, parse_dates[date]) fig, ax plt.subplots(figsize(12, 5)) ax.plot(july[year], july[tmean_anom], color#d62728, lw1.2, label7月气温距平) ax.axhline(0, colorblack, lw0.8, ls--) ax.set_xlabel(年份) ax.set_ylabel(气温距平 (°C)) ax.set_title(法国南部某站 7 月气温距平基准期 1991-2020) ax.legend() fig.tight_layout() plt.savefig(output/july_temp_anomaly.png, dpi150) plt.show()从图中可以快速判断近十年往往出现较长一段正距平正距平峰值年份大概率就是“破纪录”年份。6.2 用双轴图叠加降水与气温高温和干旱经常同时出现把 7 月降水距平百分比和气温距平画在同一个图里能直观看出两者是否同步。fig, ax1 plt.subplots(figsize(12, 6)) ax1.bar(july[year], july[precip_anom_pct], color#9ecae1, alpha0.7, label7月降水距平百分比) ax1.set_xlabel(年份) ax1.set_ylabel(降水距平 (%), color#3182bd) ax1.tick_params(axisy, labelcolor#3182bd) ax2 ax1.twinx() ax2.plot(july[year], july[tmean_anom], color#d62728, markero, ms3, label7月气温距平) ax2.axhline(0, colorgray, lw0.8, ls--) ax2.set_ylabel(气温距平 (°C), color#d62728) ax2.tick_params(axisy, labelcolor#d62728) plt.title(法国南部某站 7 月降水距平与气温距平) fig.tight_layout() plt.savefig(output/july_precip_temp.png, dpi150) plt.show()图中如果出现“降水柱状图明显偏低 气温曲线明显偏高”的年份基本就是复合极端事件的候选年份。6.3 柱状图突出破纪录年份也可以只关注 7 月平均气温本身用柱状图标出每个破纪录年份效果更接近新闻报道里的图表。colors [#d62728 if r else #bbbbbb for r in july[is_record]] fig, ax plt.subplots(figsize(12, 5)) ax.bar(july[year], july[tmean_month], colorcolors) ax.set_xlabel(年份) ax.set_ylabel(7月平均气温 (°C)) ax.set_title(法国南部某站 7 月平均气温红色为破纪录年份) fig.tight_layout() plt.savefig(output/july_record_bar.png, dpi150) plt.show()红色柱子出现的频率能直观反映“纪录被不断刷新”的节奏如果破纪录年份主要出现在近期说明近期 7 月确实频繁偏热。6.4 结果解读的边界可视化只是分析工具解读时要注意模拟数据的结果不能直接等同于真实法国气候结论真实分析必须使用官方观测数据。单个站点的结果不能代表整个法国。法国南北气候差异很大巴黎、马赛、布列塔尼的 7 月气温和降水特征完全不同。“破纪录”反映的是统计事实不自动等于“气候变化导致的”归因需要单独的事件归因研究。本文的范围是描述性统计和可视化。7. 常见问题与排查思路问题现象常见原因解决思路读取 CSV 报ParserError: day is out of range for month日期格式不规范或数据包含非法日期打印原始日期列排查统一为YYYY-MM-DD必要时用errorscoercegroupby().mean()结果全为 NaN列类型还不是数值类型或缺失值未处理用pd.to_numeric(..., errorscoerce)转换先处理缺失值ImportError: No module named netCDF4依赖未安装或安装到了错误环境确认先conda activate climate再conda install -c conda-forge netcdf4stats.gamma.fit返回 nan数据全为 0 或样本量不足过滤正降水加入最小样本数判断样本过少时改用经验分布距平符号和预期相反单位或基准期错误确认气温是摄氏度而不是开尔文确认基准期是 1991-2020matplotlib 中文乱码系统缺少中文字体设置plt.rcParams[font.sans-serif] [SimHei]或改用英文标签排查时建议从数据源头开始先看原始数据前 100 行确认日期格式、单位、缺失值表现再逐步向后排查。不要一上来就怀疑算法。8. 最佳实践与工程建议8.1 数据管理规范气象分析项目建议统一命名规则原始观测数据放在data/raw/不做任何修改清洗后的中间数据放在data/processed/图表输出到output/每个处理步骤单独存成脚本运行顺序用文件名编号控制这样做的最大好处是几个月后回看项目时能清楚知道每个文件是怎么来的也能快速定位某个结果对应哪一步处理。8.2 合理选择气候基准期世界气象组织建议使用连续的 30 年作为气候标准期当前官方推荐期是 1991-2020。选择基准期时要注意基准期不能和研究时段重叠过多否则“距平”和“破纪录”会失真不同文献使用的基准期不同比较数值前先确认口径如果是长期趋势研究还可以同时计算多个基准期结果观察稳定性8.3 单位、缺失值与质量控制气象数据最容易出错的地方就是单位温度可能是摄氏度也可能是开尔文降水可能是毫米也可能是英寸。项目一开始就要在元数据中写明单位并在代码里增加单位校验。月尺度聚合前务必设置有效观测天数阈值。一个月的观测如果缺了 15 天算出来的月平均气温几乎没有参考价值。质量控制规则应当写在文档里并随数据一起保存。8.4 破纪录与归因结论要谨慎分析师很容易把“破纪录”直接写成“气候变化导致”这在科学上是不严谨的。破纪录是统计事实归因需要回答“如果没有人类活动影响这个事件还会不会发生”这属于事件归因研究通常需要气候模式模拟。在业务报告或博文中建议的表述方式是可以写“该年 7 月平均气温为 26.4°C比 1991-2020 平均值高 1.8°C为有记录以来最高。”不要写“气候变化导致今年 7 月高温破纪录。”保持统计结论和因果推断的边界是数据分析的基本职业素养。8.5 可复现性可复现性体现在三个层面代码层面固定随机种子如np.random.seed(42)保证模拟数据一致环境层面用requirements.txt或 conda env 记录依赖版本数据层面真实数据保存下载日期、版本号、来源链接这样即使几个月后重新运行也能得到相同结果或者能快速排查出结果变化的原因。8.6 大数据量下的性能优化站点日值数据通常只有几十 MBpandas 足够处理。但如果换成 E-OBS 这类覆盖全欧洲、时间跨度几十年的栅格数据建议使用 xarray 读取 netCDF利用维度标签自动对齐数据对长时间序列做区域平均时先用sel裁剪研究区域内存不足时使用chunk和 dask 延迟计算避免在循环里逐日读写文件9. 总结与拓展通过这篇文章我们完成了一条完整的气象数据分析路径从模拟日值数据出发经过日期解析、缺失值检查、月尺度聚合计算了 1991-2020 基准期的气温距平和降水距平用累积最大值识别出 7 月破纪录年份并用简化 Gamma 分布计算了 SPI 干旱指数最后用 Matplotlib 画出了三种不同类型的可视化图。这套流程稍加修改就能推广到其他场景把月份从 7 月改成任意月份分析“史上最热 X 月”是否成立把单站改成多站用groupby([station, year, month])做区域对比把日值换成小时值分析“最热一天”或“最长热浪持续天数”用 xarray 直接读取 E-OBS 栅格数据做空间分布分析如果继续深入可以学习 SPEI 干旱指数、热浪识别算法HWMId、气候趋势显著性检验Mann-Kendall 检验等内容。这些方向都需要本文的基础功底数据清洗要扎实、统计口径要清楚、结论边界要把握住。建议你下载一份真实的法国或中国气象站数据把本文代码跑一遍再试着改成自己感兴趣的城市。数据分析和写代码一样只有亲手处理过缺失值、亲手画过破纪录年份的柱状图才能真正理解那些新闻标题背后的统计细节。如果本文对你有帮助可以收藏备用后续遇到气象数据处理问题也能随时回来查阅。