旧版本｜1980年1月1日～2023年7月31日中国各区县逆温数据&每年各区县逆温天数数据（空气污染影响相关研究的工具变量）

1558 字 约 5 分钟 数据资料 空气质量旧版本

新版本在这里：https://rstata.duanshu.com/#/brief/course/6716f6a463b641a88b7308f86b9753d0

一直以来都有小伙伴跟我要逆温数据，我之前觉得这个数据处理起来很麻烦就一直没弄，最近处理了下发现还挺简单的，所以就处理了份分享给大家～

首先我们先了解下逆温。正常情况下，海拔越高，气温越低，但是有时候会出现海拔高的地方气温反而更高，这个时候冷空气在下方，空气不易流通，污染物扩散慢。这种现象称为逆温，由于逆温对污染物扩散的外生性和相关性，所以逆温是研究空气污染影响很好的工具变量。

逆温的衡量计算方法有很多，这里我们采用 Shuai Chen（2018）中的方法进行计算：

不过文中对这种逆温的衡量计算方法并不完全准确，NASA 的 MERRA-2 数据集提供了 0.5 度×0.625度分辨率、42 个气压层（不同的海拔有着不同的气压）的气温等数据，但是并不是每个格点在 42 个气压层都有数据，例如青藏高原的海拔在 5000 m 以上，底部很多气压层是没数据的：

library(terra)
rast("data/MERRA2_100.inst6_3d_ana_Np.19900718.SUB.nc") -> rst
cairo_pdf("levtemp.pdf", width = 10, height = 10)
plot(rst)
dev.off()

所以我们计算逆温的应该采用每个格点近地的第一层和第二层（也就是最开始有数据的两层）。

因此我们可以通过下面的流程计算中国各区县的逆温：

对逆温数据计算感兴趣的小伙伴可以学习这个课程：

×

  1. 从 NASA 下载 1980 年至今的 MERRA-2 数据：https://disc.gsfc.nasa.gov/datasets/M2I6NPANA_5.12.4/summary （需要登陆注册）；
  2. 从下载到的 netCDF 数据提取近地的第一、二、三层数据（第三层数据可以用来进行稳健性分析），例如 1990 年 7 月 18 日的：

  1. 把栅格数据准换成区县面板数据（使用的是 2021 年中国县级行政区划）：

  1. 比较第一层和第二层的平均气温就可以判断逆温情况了，如果第二层的气温高于第一层就是逆温，差值就是逆温强度，反之就不是逆温，逆温强度为 0（每天都有四个时点的数据，只有有一个时点的结果是逆温当天就视为逆温天）；
  2. 另外可以使用一三层的温度进行比较得到逆温2 和逆温强度2 用于稳健性检验。

经过上述的计算就得到了 1980年1月1日～2023年7月31日的区县逆温数据：

1980年1月～2023年7月区县逆温数据.rds（13 GB）

由于这个文件高达 13 GB （使用 R 语言的 readr::read_rds() 函数读取，使用该文件请确保你的电脑有 32G 以上的内存），为了方便大家使用，我把数据拆分成了每年的 xlsx 文件：

该数据包含如下变量，各个变量的含义如下：

  • date：日期
  • year：年
  • month：月
  • day：日
  • 省：省份名称
  • 省代码：省份行政区划代码
  • 市：城市名称
  • 市代码：城市行政区划代码
  • 县：区县名称
  • 县代码：区县行政区划代码
  • T1_layer1：0 点近地第一层的平均气温
  • T1_layer2：0 点近地第二层的平均气温
  • T1_layer3：0 点近地第三层的平均气温
  • T2_layer1：6 点近地第一层的平均气温
  • T2_layer2：6 点近地第二层的平均气温
  • T2_layer3：6 点近地第三层的平均气温
  • T3_layer1：12 点近地第一层的平均气温
  • T3_layer2：12 点近地第二层的平均气温
  • T3_layer3：12 点近地第三层的平均气温
  • T4_layer1：18 点近地第一层的平均气温
  • T4_layer2：18 点近地第二层的平均气温
  • T4_layer3：18 点近地第三层的平均气温
  • T1_thermal_inv12：0 点根据一二层计算的逆温强度
  • T2_thermal_inv12：6 点根据一二层计算的逆温强度
  • T3_thermal_inv12：12 点根据一二层计算的逆温强度
  • T4_thermal_inv12：18 点根据一二层计算的逆温强度
  • T1_thermal_inv13：0 点根据一三层计算的逆温强度
  • T2_thermal_inv13：6 点根据一三层计算的逆温强度
  • T3_thermal_inv13：12 点根据一三层计算的逆温强度
  • T4_thermal_inv13：18 点根据一三层计算的逆温强度
  • thermal_inv_yesno12：根据一二层判断当天是否是逆温天
  • thermal_inv_yesno13：根据一三层判断当天是否是逆温天

进一步我们可以根据这份详细数据计算各个区县每年的逆温天数：

read_rds("1980年1月～2023年7月区县逆温数据.rds") -> df

df %>%
  group_by(year, 省, 省代码, 市, 市代码, 县, 县代码) %>%
  summarise(days1 = sum(thermal_inv_yesno12),
            days2 = sum(thermal_inv_yesno13)) %>%
  ungroup() -> yeardf
yeardf

#> # A tibble: 123,711 × 9
#>     year 省     省代码 市     市代码 县     县代码 days1 days2
#>    <dbl> <chr>   <dbl> <chr>   <dbl> <chr>  <chr>  <dbl> <dbl>
#>  1  1980 安徽省 340000 安庆市 340800 大观区 340803   203   144
#>  2  1980 安徽省 340000 安庆市 340800 怀宁县 340822   213   120
#>  3  1980 安徽省 340000 安庆市 340800 潜山市 340882   212   108
#>  4  1980 安徽省 340000 安庆市 340800 宿松县 340826   160   101
#>  5  1980 安徽省 340000 安庆市 340800 太湖县 340825   191   113
#>  6  1980 安徽省 340000 安庆市 340800 桐城市 340881   223   151
#>  7  1980 安徽省 340000 安庆市 340800 望江县 340827   189   112
#>  8  1980 安徽省 340000 安庆市 340800 宜秀区 340811   206   133
#>  9  1980 安徽省 340000 安庆市 340800 迎江区 340802   204   137
#> 10  1980 安徽省 340000 安庆市 340800 岳西县 340828   167   109
#> # … with 123,701 more rows

其中 days1 表示根据一二层判断的逆温天数，days3 表示根据一三层判断的逆温天数。

下图展示了 2021 年各个区县的逆温天数分布：

还挺有意思，原来华北地区的空气状况不好也跟逆温天太多有关系！

点击这里跳转到 RStata 短书平台获取附件：旧版本｜1980年1月1日～2023年7月31日中国各区县逆温数据&每年各区县逆温天数数据（空气污染影响相关研究的工具变量）

相关文章

评论