RClimDex 1.0实战:27个极端气候指数计算与QC排错
简介RClimDex 极端气候指数计算软件配套用户手册面向气候研究、气象监测、生态与水文分析人员用于在 R 环境中完成极端气候指数计算、数据质量控制及结果解读。手册覆盖软件加载、缺失值与异常值检查、高温日数、低温日数、极端降雨量等指数计算涵盖 ETCCDMI 推荐的 27 项核心指数及用户自定义阈值指数并说明阈值与基准期设置、已知问题反馈方式附录汇总指数清单、输入数据格式、指数定义及 Windows 环境常见问题便于按主题查阅。压缩包仅含 1 个 doc 文件大小约 600KB属于轻量级工具文档适合作为实际分析时的速查资料。已有 3839 人学习下载可帮助读者理解极端气候事件监测流程并为气候变化检测分析提供可复用的软件说明与指标参考。1. 一份 2004 年的 R 脚本为什么今天还在算极端气候指数手里这份 RClimDex 1.0 的用户手册只有二十多页讲的是一个跑在 R 上的极端气候指数计算包。它不训练模型、不做可视化大屏只干一件事把逐日气象站观测降水量、最高温、最低温按统一口径换算成 27 个极端气候指数再顺手画一条趋势线。做气候变化检测、区域极端事件统计、台站资料整编的人绕不开这套口径——ETCCDMI 定义的 27 个核心指数很多论文和业务报告里的 FD0、SU25、R95p 就是按它算的。新手拿它当“算指数的计算器”能直接上手熟手更在意它的 QC 规则和百分位指数的 bootstrap 处理这两处决定了结果能不能跨站、跨期比较。下面按“跑起来 → 喂数据 → 设参数 → 排错 → 接自己的分析”走一遍。2. 装环境、加载脚本让 RClimDex 的图形菜单弹出来RClimDex 1.0 不是一个能install.packages()装下来的包它是一份rclimdex.r源码脚本source()进 R 之后弹出一个 TclTk 写的图形菜单。这个设计在 2004 年是合理的——当年 R 的包管理还不成熟作者要的是一个“拷过去就能跑”的单文件工具。理解这一点很关键后面所有“装不上”“找不到包”“菜单不弹”的问题根子都在这里。2.1 R 与 TclTkGUI 到底依赖什么RClimDex 的界面全部由 TclTk 绘制菜单、弹窗、文件选择框都是 TclTk 控件。R 从 1.9.0 版本开始把 TclTk 打进默认安装包所以只要装的是近十几年的 RWindows 用户基本不会缺这一块。但 Linux 发行版经常把 R 拆成多个包r-base装了tk没装结果就是source(rclimdex.r)执行完毫无反应控制台也不报错——这是最典型的翻车现场之一。判断方法很简单在 R 控制台里跑一行# TRUE 表示当前 R 编译时链接了 TclTkFALSE 则图形菜单一定弹不出来 capabilities(tcltk)返回FALSE就别往下折腾数据了先把 TclTk 补上。返回TRUE但菜单仍然不出现那问题在别处通常是脚本没找到、或者source()的路径不对。需要提醒的是RClimDex 1.0 是在 R 1.84 时代写的脚本里用到的一些老式写法在新版 R 上偶尔会有兼容性摩擦但核心的读数据、算指数、画图逻辑都能正常跑。我自己的习惯是把它放在一个独立的、路径全英文的工作目录里不跟别的项目混在一起出问题好排查。2.2 把 R 装到位Windows 与 Linux 两条路Windows 下最省事从 R 项目官网下载对应系统的安装程序一路默认下一步桌面会出现 R 图标双击进控制台。安装时不要改动默认组件TclTk 就在默认组件里。如果你用的是很老的 R 版本1.9.0 以前需要单独补 TclTk这种环境现在基本见不到了。Linux 下用包管理器装注意把图形依赖一起带上# Debian/Ubuntu 系r-base 是解释器tk/tcl 是 RClimDex 界面依赖的图形库 sudo apt-get update sudo apt-get install -y r-base r-base-dev tk tcl # 装完立刻验证输出 TRUE 才说明图形界面这条路通了 Rscript -e print(capabilities(tcltk))这里的r-base-dev不是必需的但编译一些辅助包时会用到一起装上省事。tk和tcl才是关键缺了它们 R 能跑、脚本能source就是没有窗口。Windows 上还有一个历史坑部分机器启动 R 时提示找不到HOME环境变量脚本读写临时文件会失败。遇到这种情况在系统环境变量里加一个HOME指向一个存在的目录即可比如C:\Users\你的用户名。这个问题在手册的附录里也提过属于老版本 Windows 的遗留问题新系统上很少见了。2.3 加载 rclimdex.r 的两种姿势与首次自检脚本加载有两条路。第一条是命令行方式最直接# 先把工作目录切到 rclimdex.r 所在的位置 # 路径里不要有中文和空格R 对这两样东西的容忍度很低 setwd(D:/climate/rclimdex) # 确认文件确实在当前目录下返回 TRUE 再继续 file.exists(rclimdex.r) # 加载脚本成功后会自动弹出 RClimDex 主菜单窗口 source(rclimdex.r)第二条路是图形操作R 的菜单里选“文件 → 载入 R 代码”在弹窗里挑中rclimdex.r。两种方式等价但排查问题时我更信命令行——source()报的错会直接打在控制台里图形方式有时候把错误吞掉只剩一个“没反应”。source()成功之后屏幕上会出现一个标题为 RClimDex 的小窗口里面大致是“载入数据并运行质量控制”“指数计算”“退出”这几项。看到这个窗口环境这一关就算过了。这里有个容易被忽略的细节source()只是把函数和菜单定义加载进当前会话并没有把数据读进来。如果你关掉菜单窗口函数还在重新敲一次source()或者重新调用入口函数又能弹出来。反过来如果重启了 R之前source的东西全没了必须重新加载。我见过有人换了台机器直接点“指数计算”然后困惑为什么没有任何反应——因为新会话里根本没加载过脚本。再补一个自检习惯加载完成后先别急着读数据在控制台里敲getwd()看一眼当前目录确认它就是rclimdex.r和你的站点文件所在的地方。RClimDex 的日志目录log、结果目录indices、图件目录plots都是相对当前工作目录创建的工作目录错了结果会散落到你意想不到的地方找起来非常费劲。3. 喂数据与质量控制-99.9、NA 和 n 倍标准差数据质量控制在 RClimDex 里不是可选项而是计算指数的前置步骤手册里写得很直白必须先做 QC再算指数。原因不复杂——极端指数大量依赖阈值判断比如“日最高温大于 25 度”一个荒唐的极值就能把整年的统计带偏而气象站原始数据里恰恰不缺这种东西传感器故障、传输丢包、人工录入错位都会留下痕迹。3.1 六列定宽文本输入文件的字段与格式错RClimDex 吃的输入是纯文本一行一天六列空格或制表符分隔没有表头按时间升序排列。列的顺序是固定的年、月、日、日降水量、日最高温、日最低温。缺测统一用-99.9占位。1961 1 1 0.0 -3.2 -11.5 1961 1 2 2.4 -1.8 -9.6 1961 1 3 0.0 -2.5 -10.1 1961 1 4 7.8 0.6 -6.4四行的含义分别是1961 年 1 月 1 日无降水最高温 -3.2℃最低温 -11.5℃后面依次类推。降水单位是毫米温度单位是摄氏度不要混入开尔文或者华氏度程序不会帮你换算。文件名建议按站号.txt的规则来比如demo01.txt。这不是强制要求但 RClimDex 会用文件名前缀给输出文件命名站号规范一点后面几十个站的结果摆在一起时你会感谢自己。格式错是读取失败的头号原因常见的几类文件里带了表头行、列数不是六列、日期顺序写成“日月年”、用逗号而不是空格分隔、文件编码带 BOM。R 控制台里一旦出现读取相关的报错第一件事就是把文件头和文件尾各拉几行出来跟上面的样例逐列对一遍九成的错在这一步就能看出来。还有一类隐蔽的错数据不是严格按日期升序。RClimDex 内部有按年、按日序取值的逻辑顺序乱了不会立刻报错但算出来的百分位阈值可能是错的。这种事最坑因为结果看起来“有数”不核对根本发现不了。3.2 质量控制做了什么四类清洗与一类标记点开“载入数据并运行 QC”选中站点文件程序会做四件事第一把所有-99.9替换成 R 内部的缺失值NA。这一步是格式转换不是判断-99.9在 RClimDex 里就是“缺测”的约定记号。第二把日降水量小于 0 的值置为NA。负降水在物理上不成立只可能是仪器或录入问题。第三把日最高温小于日最低温的记录置为NA。这个判断看起来很笨但极其有效能抓出一大批列错位、符号错、单位错的数据。第四标记离群值。规则是按“同日序的气候态”算均值和标准差落在[mean - n*std, mean n*std]之外的值标为可疑。这里的n由用户在弹窗里填默认是 3。手册里明确建议3 会标出非常多的点可以先从 4 开始试。这里必须说清楚一个关键区别前三类是“改”第四类是“标”。被标为离群的值不会从数据里删掉indcal.csv里这些值原样保留只是另外写进tepstdQC.csv供你人工核对。这个设计是合理的——统计上偏离不等于数据错误一次真实的极端高温事件本来就应该离群程序没有资格替你删掉它。n怎么定是个经验活。台站资料年头长、质量控制好n3问题不大资料年头短、站点海拔高、日较差大n3会标出一大片人工核对的工作量爆炸。我一般的做法是先用 4 跑一遍看看标记比例如果某个站被标了超过 1% 的日值就把n放到 4.5 或 5同时把tepstdQC.csv里标记最集中的月份挑出来回去查原始报表。3.3 log 目录里的四个文件怎么看QC 跑完当前工作目录下会多出一个log子目录里面固定生成四个 CSV 和四张 PDF 时序图。四个 CSV 各有各的用处别混文件名以 demo01 为例内容怎么用demo01tempQC.csv温度类不合理值清单查 Tmax Tmin、明显越界的记录demo01prcpQC.csv降水类不合理值清单查负降水记录demo01tepstdQC.csv温度离群值清单带具体日期逐条人工判断是真是假demo01indcal.csv清洗后的数据指数计算的真正输入想手工修补就在这里改四张 PDF 分别是降水、最高温、最低温、日较差的时序图缺测点用红点标出。这份图价值很高一眼就能看出资料有没有整段缺失、有没有台阶式跳变。整段缺失往往意味着台站迁址或者观测中断台阶跳变则可能是仪器更换这两种情况都会让指数序列出现虚假的“突变”做趋势分析前必须先心里有数。3.4 异常值只标不改改不改由你定tepstdQC.csv拿到手工作才刚开始。这份清单要逐条看判断标准是常识加原始记录盛夏午后 38℃ 在多数地区正常隆冬凌晨 38℃ 就一定是错的日最低温比前后几天高出十几度而最高温没动多半是录入时小数点错了位。确认有错的在indcal.csv里改。Windows 下用 Excel 打开改是最快的注意保存时选 CSV 格式别存成 xlsx。改完必须重新走一遍“载入数据并运行 QC”——手册特意强调过这件事指数是从 QC 后的数据算的原始输入文件全程不被修改你改了原始文件却不重跑 QC指数结果不会跟着变。这一点我吃过亏改完数据直接点“指数计算”出来的还是旧结果对着屏幕愣了半天。4. 算 27 个核心指数基准期、纬度与五组自定义阈值数据洗干净进入正题。RClimDex 能算 27 个 ETCCDMI 核心指数外加一批带用户自定义阈值的温度、降水指数。菜单里点“指数计算”先弹参数窗再弹指数勾选窗两步都不能跳。4.1 27 个指数分成四组先想清楚要哪几个不勾选全部指数是个务实的选择。百分位类指数要跑 bootstrap27 个全勾一个长资料站点跑十几分钟很正常而你写报告可能只需要五六个数。按计算逻辑分四组看组别代表指数单位计算依赖绝对阈值日数FD0、SU25、ID0、TR20天日最高温、日最低温百分位类TN10p、TX10p、TN90p、TX90p天百分比天数基准期 bootstrap持续期类GSL、CDD、CWD、WSDI、CSDI天连续日判定极值与降水总量TXx、TNn、Rx1day、Rx5day、SDII、R95p、R99p、PRCPTOT℃ / mm日降水、日温度FD0 数的是日最低温小于 0℃ 的天数SU25 数的是日最高温大于 25℃ 的天数这两个最直观也最适合拿来当“结果对不对”的第一道验证。GSL 生长季长度的定义稍微绕北半球按 1 月 1 日到 12 月 31 日先找到第一个连续 6 天日均温大于 5℃ 的时段再找 7 月 1 日之后第一个连续 6 天日均温小于 5℃ 的时段两者之间的天数就是 GSL南半球整体平移半年从 7 月 1 日算起。日均温由最高温和最低温取平均得到。百分位类指数是这套工具的核心竞争力所在。TN10p 是日最低温低于基准期第 10 百分位的天数占比问题在于如果整段资料共用一个固定的百分位阈值序列前后期的比较会带上“不均匀性”趋势可能是算法造出来的而不是气候本身的。正确做法是 bootstrap——对每个日历年用基准期做重采样逐年重算阈值再判定。这个流程在 Excel 里几乎没法实现也正是 RClimDex 存在的理由。4.2 基准期、纬度与五组阈值一个都不能填错参数窗里有五类输入填错任何一个都会让结果悄悄跑偏基准期起止年。百分位阈值和 bootstrap 都基于这段时间常用 1961—1990也有用 1981—2010 的选哪个取决于你的资料长度和研究口径但必须写清楚并在全文统一。基准期太短少于 20 年会让百分位阈值不稳定太长又可能跨越资料断档。台站纬度。北半球填正南半球填负。这个参数唯一的作用是决定 GSL 和部分季节相关指数走北半球分支还是南半球分支。填错符号不会报错只会让南半球站点的生长季长度算成完全错误的数——属于典型的静默错误。日降水阈值 P。用于 Rnn 指数即日降水量超过该阈值的日数。填 1 就是 R1mm填 10 就是 R10mm填 20 就是 R20mm。想同时要几个就分几次跑或者一次填一个主力阈值。四组自定义温度阈值。对应 SUmm、FDmm、TRmm、IDmm 四个指数日最高温上限超过即计入 SUmm、日最高温下限低于即计入 FDmm、日最低温上限超过即计入 TRmm、日最低温下限低于即计入 IDmm。这里的mm就是你填的数值比如最高温上限填 35生成的指数就叫 SU35。站名或站号。QC 窗口里也出现过这个输入框那里可以不填但指数计算这一步建议填上它会进到输出文件里方便后面合并多站结果时溯源。填完点确定程序会先做一轮数据预处理慢的机器上要等几十秒别以为卡死了。预处理完成后弹出指数勾选窗默认全选按需取消再点确定。4.3 indices 与 plots结果文件怎么读计算完成后工作目录下多出两个子目录。indices里是结果命名规则是站号_指数名.csv列之间用逗号分隔一行一年。plots里是每个指数的年序列图JPEG 格式文件名规则相同。图上有两条线实线是普通最小二乘线性趋势虚线是局部加权回归LOESS平滑图角上还标了线性拟合的统计量。实线看整体斜率虚线看年代际波动两者不一致的时候要留神——比如实线平缓但虚线有明显的“先降后升”说明趋势被前后两段抵消了只报一个斜率会掩盖真实结构。结果文件是一年一行的年序列。这里有个使用习惯值得养成拿到indices目录先别急着做趋势检验挑一个站、挑一个指数回到indcal.csv里手工数一遍。比如 SU25直接在清洗后的数据里筛出日最高温大于 25 的记录按年计数跟demo01_SU25.csv对一遍。对上了说明数据读取、QC、阈值判定这条链路是通的对不上问题多半出在缺测处理或者年份边界上。这个笨办法我每次都做花十分钟能省掉后面几天的返工。5. 避坑与排查五类翻车现场和它们的解法这套工具用了二十年坑基本被踩遍了。下面五条是我自己和周围人真金白银换来的按“现象 → 原因 → 解决”写遇到问题可以直接对号入座。5.1 首年等于基准期首年程序直接崩现象数据从 1961 年开始基准期也填 1961—1990点完确定R 控制台刷出一串报错计算中断indices目录里什么都没有。原因百分位指数在计算年初和年末的阈值时需要取基准期边界之外的数据做窗口平滑。手册里把这条列为已知 bug——当资料首年与基准期首年重合程序取不到边界外的记录直接抛错。这是 RClimDex 1.0 的固有问题不是你的数据错了。解决在输入文件最前面手工插一行基准期开始前一年的 12 月 31 日三个要素全填缺测1960 12 31 -99.9 -99.9 -99.9 1961 1 1 0.0 -3.2 -11.5 1961 1 2 2.4 -1.8 -9.6这一行不参与任何统计只是给程序一个“边界外还有记录”的信号。插完重新走一遍 QC 和指数计算即可。注意日期必须是基准期首年的前一年年末基准期换成 1981—2010 就插 1980 年 12 月 31 日。5.2 脚本加载了菜单不弹控制台也没报错现象source(rclimdex.r)回车光标回到提示符一切安静就是没有窗口。原因三种可能。一是 TclTk 没编进当前 Rcapabilities(tcltk)返回FALSE二是在无图形界面的环境里跑比如纯命令行的服务器会话根本没有显示服务三是脚本source成功但入口函数没被自动调用。解决先跑capabilities(tcltk)确认依赖返回FALSE就补装tk、tcl后重启 R。在远程服务器上跑的话这套图形菜单本来就不适用只能退回到读取indcal.csv自己写脚本算指数。第三种情况最少见看一眼source()有没有返回值、控制台有没有被警告信息刷屏通常能定位。5.3 读数据报错或者读进去了但年份不对现象选完文件控制台报错或者不报错但算出来的年序列缺了好几年。原因格式问题。表头没删、列数不对、分隔符用了逗号、日月顺序颠倒、文件里有空行、编码带 BOM都会导致读取失败或者错位。年份缺失则多是因为某些年份的天数不足 365 天被程序判为不完整年跳过。解决用文本编辑器打开文件确认第一行就是数据、每行六列、只有空格或制表符、严格按年月日升序。年份缺失要回去查原始资料看那几年是不是真有整段缺测。缺测比例过高的年份与其硬算不如在结果里标注出来。5.4 改了原始数据结果却没变现象发现某天最低温录错了改了原始文件重新点“指数计算”出来的值和改之前一模一样。原因指数是从 QC 后的indcal.csv算的不是从原始文件算的。原始文件只在“载入数据并运行 QC”这一步被读取你不重跑 QCindcal.csv还是旧的。解决改完数据必须完整重走“载入数据并运行 QC → 指数计算”。更稳妥的做法是直接改log/indcal.csv改完从指数计算那一步接着走。但要记住indcal.csv会在下次 QC 时被覆盖所以改动要同步回原始文件否则下次一跑 QC 全白改。5.5 输出文件扩展名和手册写的不一致现象手册里说结果文件叫21946_XXX.cvs实际生成的是.csv用 Excel 双击打开没问题但写脚本按.cvs去找就找不到文件。原因手册里的cvs是笔误实际输出是标准的.csv。解决批处理脚本一律按.csv匹配。这类文档笔误还有几处比如手册里对某些指数的单位描述前后不完全一致真正算的时候以输出文件里的数值和定义为准别死抠文档字面。6. 把结果接进自己的分析批量合并与一次手工验证单站跑通只是起点真正的工作量在几十上百个站。RClimDex 的图形界面不适合批量但它的输出结构非常规整——indices目录下全是站号_指数名.csv一行一年列名固定。这意味着批量处理可以完全脱离界面用脚本把结果拼起来。我一般这么干先把所有站跑完把各自的indices目录汇总到一个大目录下然后用一段 R 代码按指数名分组读取、纵向拼接顺便算 Sen 斜率做趋势# 汇总目录下所有站点、所有指数的结果文件 files - list.files(all_indices, pattern _SU25\\.csv$, full.names TRUE) # 逐个读取并加上站号列站号从文件名里截出来 res - do.call(rbind, lapply(files, function(f) { d - read.csv(f, header TRUE) d$station - sub(_.*$, , basename(f)) # 文件名前缀即站号 d })) # 按站号统计 SU25 的多年均值快速看空间分布是否合理 agg - aggregate(SU25 ~ station, data res, FUN mean, na.rm TRUE) # 趋势检验建议用 trend 包里的 Sen 斜率 Mann-Kendall别只看最小二乘 # library(trend); sens.slope(res$SU25[res$station demo01])这段代码的关键点在pattern和站号截取两处pattern用指数名做后缀匹配避免把温度指数和降水指数混在一起站号靠文件名下划线前的部分截取所以前面强调过文件名规范要统一。真要做得细还得处理站点资料长度不一、年份不齐的问题merge按年对齐比rbind更稳妥。验证方法我固定用一个挑一个资料完整、气候背景熟悉的站把demo01_SU25.csv里某一年的数值跟indcal.csv里该年日最高温大于 25 的记录条数手工对一遍再挑一个降水指数比如 R95p用该年日降水排序后第 95 百分位以上的降水量求和核对。两个指数都对上才认为这条流水线是可信的。血泪经验是不要相信“跑完没报错”就等于结果正确。这套工具最危险的错误全是静默的——纬度符号填反、基准期填错、数据没按年升序程序一句提示都没有出来的图还挺好看。从那以后我每次跑完新站点都强制走一遍手工核对再顺手看一眼plots里那四张 QC 时序图有没有整段红点。希望帮到你。本文还有配套的精品资源点击获取