Python+Stata:如何获取中国气象历史数据

发布时间:2021-11-27 阅读 1764

Stata连享会   主页 || 视频 || 推文 || 知乎 || Bilibili 站

温馨提示: 定期 清理浏览器缓存,可以获得最佳浏览体验。

New! lianxh 命令发布了:
随时搜索推文、Stata 资源。安装:
. ssc install lianxh
详情参见帮助文件 (有惊喜):
. help lianxh
连享会新命令:cnssc, ihelp, rdbalance, gitee, installpkg

课程详情 https://gitee.com/lianxh/Course

课程主页 https://gitee.com/lianxh/Course

⛳ Stata 系列推文:

PDF下载 - 推文合集

作者:秦利宾 | 许梦洁

编者按:本文主要参考自「2000-2020年中国地面气象数据:从NOAA到分省面板」「空气质量与气象网站 (网友整理)」,特此致谢。


目录


1. 数据介绍

中国气象历史数据来源是美国国家气候数据中心 (NCDC),其数据存放网址为 https://www.ncei.noaa.gov/data/。具体来看,下图左半部分美国国家气候数据中心的数据清单,右半部分为我们感兴趣的全球气象小时数据。

可以看到,global-hourly 目录下有三个文件夹,其中 access 存放的是全球每年每个站点的所有监测数据。例如,我们打开 access 目录下的 2020 年文件夹,则可以看到 2020 年所有站点的数据。文件 00841599999.csv 命名规则是,前 6 位为站点 ID。

然后,我们打开厦门高崎站点 2020 年的监测数据表 (关于如何知道厦门高崎站点 ID,详见下文),可以看出该表提供了站点 ID、监测时间、站点名称等基础信息,以及相关监测信息 (下图仅列示部分信息)。

archive 文件夹下分别存放的是 csv 格式和 isd 格式的压缩文件。例如,我们打开 archive 目录下的 csv 文件夹,则可以看到按年打包的压缩文件。

doc 文件夹下存放的是关于数据的帮助文档。

当然国内更为方便的途径是空气质量与气象网站 (网友整理,不定期更新)空气质量与气象网站主要提供以下几类数据:

  • 中国空气质量历史数据 (2014-05-13 至今);
  • 北京市空气质量历史数据 (2013-12-06 至今);
  • 中国气象历史数据 (1942-07 至今)。

其中,我们关心的中国气象历史数据选取了中国区域 (含港澳台) 的观测站点数据,并重新按年打包。其主要特征如下:

  • 时间范围:1942 年至今;
  • 时间精度:近年的数据大多为 3 小时数据,少量站点有 1 小时数据;
  • 站点数量:近年为 400 多个;
  • 气象要素:气温、气压、露点、风向风速、云量、降水量。

本次案例演示,我们使用的数据是空气质量与气象网站的中国气象历史数据。该原始数据可通过以下几种方式获取:

下图为中国气象数据百度网盘 (空气质量与气象网站开发者提供) 中的数据目录,该目录下主要有两类文件:一是存放气象数据的压缩文件 china_isd_lite_****.zip,解压后是几百个站点数据文件,每个文件是单个站点全年的数据。二是数据说明文件,主要是关于数据和指标的介绍,其中比较重要的是 中国地面气象站基本要素观测资料表.pdf

接着,我们以 china_isd_lite_2020.zip 为例来查看文档结构。在下图中,左半部分为压缩文件中打包的文件。其数据格式为 ISD-Lite,是简化的 ISD (Integrated Surface Data) 数据。文件 450070-99999-2020 的命名规则是,第 1 段数字代表站点 ID,第 3 段数字代表年份。

右半部分为文件内容,主要有 9 列,依次为观测年份、观测月份、观测日期、观测小时、空气温度、露点温度、海平面气压、风向、风速、云层厚度、液体渗透深度 (1小时)、液体渗透深度 (6小时)。关于指标的详细定义、单位、和换算规则,请参考原始数据中的 isd-lite-format.txt 文档。

在获取到站点 ID 和站点监测数据后,我们如何将其对应到具体地址,例如省份。这里就需要使用到上文提及的 中国地面气象站基本要素观测资料表.pdf。该文件来自中国国家气象科学数据中心的中国气象数据网,里面的区站号末尾加 0 即为站点 ID。

在此基础上,我们就可以整理得到中国某地的气象数据。为了帮助大家更好的掌握数据处理流程,我们将以「省份年度气象数据」的整理为案例进行演示。

2. 数据处理

对于 ISD-Lite 格式的文件,Stata 并不能直接处理,因此我们需要借助 Python 来转换合并。Python 的主要工作包括以下几项:

  • 读取各站点监测数据,并按年度合并为 txt 文件 (在这一过程中,需要新生成站点 ID 变量);
  • 由于 中国地面气象站基本要素观测资料表.pdf 是 pdf 格式,我们需要借助 Python 的 pdfplumber 提取表格。
# 将各站点监测数据按年合并为txt文件
# 导入包
import os
import pandas as pd
import zipfile
import csv
import shutil

os.chdir(r"D:\中国气象数据")
os.getcwd()

# 解压文件
files = os.listdir()
files = [file for file in files if ".zip" in file]
for file in files:
    zip_file = zipfile.ZipFile(file) 
    zip_list = zip_file.namelist() # 得到压缩包里所有文件
    for subfile in zip_list:       # 循环解压文件
        zip_file.extract(subfile)  

# 将年度所有文件写入一个txt中
files = os.listdir() 
files = [file for file in files if os.path.isdir(file) and "china" in file]

for file in files:
    with open(file+".txt", "w", encoding="utf-8") as f:
    	f.write("")
    subfiles = os.listdir(file)
    for subfile in subfiles:
        with open(r"./{}/{}".format(file,subfile),"r") as f:
            content = f.read()
        station = subfile.split("-")[0]
        content = content.replace("\n", " {}\n".format(station)) # 纵向追加气象站id
        with open(file+".txt", "a", encoding="utf-8") as f:
            f.write(content)
    print("整理完文件夹{}!".format(file))  

# 删除解压的文件夹
for file in files:
    if os.path.exists(file):  
        shutil.rmtree(file)
# 将中国地面气象站基本气象要素观测资料台站表pdf表转换为csv
# !pip3 install pdfplumber #安装并导入pdfplumber
import pdfplumber
pdf =  pdfplumber.open("./_数据说明/中国地面气象站基本气象要素观测资料台站表.pdf") # 读取pdf文件

with open("./_数据说明/台站表.csv", "w", newline="", encoding="gb18030") as f:
    writer = csv.writer(f)
    writer.writerow("")
    
for i in range(0,len(pdf.pages)):
    page = pdf.pages[i]
    table = page.extract_table()
    for row in table:
        with open("./_数据说明/台站表.csv", "a", newline="", encoding="gb18030") as f:
            writer = csv.writer(f)
            writer.writerow(row)
    print("整理完第{}页!".format(i))

需要注意的是,由于 中国地面气象站基本要素观测资料表.pdf 中有几个站名长度问题,pdfplumber 识别并不是太准确,例如黄泛区农试站、科尔沁右翼中旗、乌鲁木齐牧试站等。为此,大家可以手工修改一下。

在完成上述工作后,我们就可以用 Stata 来合并站点 ID 数据和气象数据,并根据年度和省份分组计算「省份年度气象数据」。

**# 1 中国地面气象站基本气象要素观测资料台站表
    import delimited using "./_数据说明/台站表.csv", clear
    tostring 区站号, replace
    replace 区站号 = 区站号 + "0" //里面的区站号加上个 0 即站点 ID
    rename 省份 province
    rename 区站号 station
    rename 站名 location
    rename 纬度度分 altitude
    rename 经度度分 latitude
    rename 气压传感器拔海高度米 sensor_sea_level
    rename 观测场拔海高度米 station_sea_level
    save 台站表.dta, replace

**# 2 获取省级年度气象数据
    cap mkdir tmp //创建临时文件夹
    cap fs *.txt
    foreach file in `r(files)'{
        import delimited using "`file'", clear
        replace v1 = stritrim(v1)
        split v1, p(" ")
        drop v1
		
        gen year = real(v11)
        gen month = real(v12)
        gen day = real(v13)
        gen hour = real((v14))
        gen air_temperature = real(v15)
        gen dew_point_temperature = real(v16)
        gen sea_level_pressure = real(v17)
        gen wind_direction = real(v18)
        gen wind_speed = real(v19)
        gen cloud_cover = real(v110)
        gen precipitation_one_hour = real(v111)
        gen precipitation_six_hour = real(v112)
        gen station = v113
		
        keep year-station
        merge m:1 station using 台站表
        keep if _merge == 3
        drop _merge
		
       *每个观测站每日平均
        foreach var of varlist air_temperature-precipitation_six_hour{
            replace `var' = . if `var' == -9999 //将 -9999 替换成缺失值.
            bys year month day station: egen `var'_day = mean(`var')
        }
        bys year month day station: keep if _n == _N
		
        *每个省份年度平均
        foreach var of varlist air_temperature_day-precipitation_six_hour_day{
            local varname = subinstr("`var'", "_day", "",.)
            bys province: egen `varname'_year = mean(`var')
        }
        bys province: keep if _n == _N
		
        keep year province air_temperature_year-precipitation_six_hour_year
        order year province air_temperature_year-precipitation_six_hour_year
        local dtaname = subinstr("`file'", ".txt", "",.)
        save "./tmp/`dtaname'", replace    
    }
	
    *合并数据
    clear 
    fs "./tmp/*.dta"
    foreach file in `r(files)'{
        append using "./tmp/`file'"
    }
    drop precipitation_one_hour_year //缺失数值太多
    sort province year
    save 省级气象数据.dta, replace
    ! rmdir /s/q "./tmp" //删除临时文件夹

3. 结束语

本文仅是根据中国气象数据,整理出年度省份数据。推而广之,我们还可以整理出每个站点和每个城市的日度、月度、以及年度的气象数据。如果需要更为丰富的中国气象数据,大家可从美国国家气候数据中心 (NCDC) 获得。最后,欢迎大家分享获取中国气象数据经验,文中有不足之处还请谅解。

4. 相关推文

Note:产生如下推文列表的 Stata 命令为:
lianxh 数据, m
安装最新版 lianxh 命令:
ssc install lianxh, replace

相关课程

免费公开课

最新课程-直播课

专题 嘉宾 直播/回看视频
最新专题 文本分析、机器学习、效率专题、生存分析等
研究设计 连玉君 我的特斯拉-实证研究设计-幻灯片-
面板模型 连玉君 动态面板模型-幻灯片-
面板模型 连玉君 直击面板数据模型 [免费公开课,2小时]
  • Note: 部分课程的资料,PPT 等可以前往 连享会-直播课 主页查看,下载。

课程主页

课程主页

关于我们

  • Stata连享会 由中山大学连玉君老师团队创办,定期分享实证分析经验。
  • 连享会-主页知乎专栏,700+ 推文,实证分析不再抓狂。直播间 有很多视频课程,可以随时观看。
  • 公众号关键词搜索/回复 功能已经上线。大家可以在公众号左下角点击键盘图标,输入简要关键词,以便快速呈现历史推文,获取工具软件和数据下载。常见关键词:课程, 直播, 视频, 客服, 模型设定, 研究设计, stata, plus, 绘图, 编程, 面板, 论文重现, 可视化, RDD, DID, PSM, 合成控制法

连享会小程序:扫一扫,看推文,看视频……

扫码加入连享会微信群,提问交流更方便

✏ 连享会-常见问题解答:
https://gitee.com/lianxh/Course/wikis

New! lianxhsongbl 命令发布了:
随时搜索连享会推文、Stata 资源,安装命令如下:
. ssc install lianxh
使用详情参见帮助文件 (有惊喜):
. help lianxh