1961-2015 年中巴經(jīng)濟(jì)走廊逐日氣象數(shù)據(jù)集

陳金雨，陶輝，劉金平

1. 中國科學(xué)院新疆生態(tài)與地理研究所荒漠與綠洲生態(tài)國家重點(diǎn)實(shí)驗(yàn)室，烏魯木齊830011

2. 中國科學(xué)院大學(xué)，北京 100049

3. 華北水利水電大學(xué)測(cè)繪與地理信息學(xué)院，鄭州 450046

引 言

高時(shí)空分辨率的氣象數(shù)據(jù)對(duì)氣候變化研究至關(guān)重要[1]。在發(fā)展中國家，由于基礎(chǔ)設(shè)施的建設(shè)和維護(hù)成本高，地面氣象觀測(cè)站稀少，難以獲得能夠滿足研究需求的高空間分辨率氣象數(shù)據(jù)集[2]。中巴經(jīng)濟(jì)走廊（China-Pakistan Economic Corridor，CPEC）地處南亞次大陸西北部，介于24°N-40°N 和60°E-80°E，氣候類型多樣；北起中國喀什地區(qū)，南至巴基斯坦的瓜達(dá)爾港，是“一帶一路”的重要組成部分[3]。高時(shí)空分辨率氣象數(shù)據(jù)對(duì)中巴經(jīng)濟(jì)走廊地區(qū)氣候變化研究和氣象災(zāi)害風(fēng)險(xiǎn)評(píng)估具有重要意義。

國際上已相繼建立了多氣象要素和多時(shí)空分辨率的全球格點(diǎn)數(shù)據(jù)集，如CHIRPS[4]、MSWEP[5]、PGFMD[6]、CPC[7]。但大多數(shù)據(jù)集空間尺度較大、時(shí)間序列不一，在研究中小區(qū)域尺度氣候變化時(shí)存在偏差[8-9]。尤其在地形復(fù)雜的地區(qū)，表達(dá)區(qū)域氣候特征能力有限[10]。因此，為了得到適合中小區(qū)域尺度長時(shí)間序列、高精度的格點(diǎn)化氣象數(shù)據(jù)，通常會(huì)使用空間插值技術(shù)[11]。由澳大利亞國立大學(xué)開發(fā)的氣象數(shù)據(jù)空間插值軟件ANUSPLIN 可以有效地模擬地形對(duì)降水的影響，其基于薄盤光滑樣條插值技術(shù)在許多研究中得到了應(yīng)用，并被證明是可靠的空間插值方法[12-15]。ANUSPLIN 已為許多氣象數(shù)據(jù)集的構(gòu)建做出了貢獻(xiàn)，且在世界各地得到了廣泛的應(yīng)用[16-18]。然而，中巴經(jīng)濟(jì)走廊地區(qū)目前還沒有一套完整的高時(shí)空分辨率的氣象數(shù)據(jù)集，這使得在該地區(qū)開展氣候變化相關(guān)研究具有一定困難。

本數(shù)據(jù)集以中巴經(jīng)濟(jì)走廊及其周邊地區(qū)日降水、日最高氣溫和日最低氣溫站點(diǎn)數(shù)據(jù)為基礎(chǔ)，結(jié)合研究區(qū)DEM 數(shù)據(jù)，利用ANUSPLIN 軟件進(jìn)行數(shù)據(jù)處理，經(jīng)過數(shù)據(jù)重采樣和空間插值，生成中巴經(jīng)濟(jì)走廊地區(qū)0.25°×0.25°空間分辨率氣象數(shù)據(jù)集，并利用廣義交叉驗(yàn)證和統(tǒng)計(jì)學(xué)方法對(duì)數(shù)據(jù)集進(jìn)行質(zhì)量評(píng)估，得到結(jié)果可為中巴經(jīng)濟(jì)走廊地區(qū)氣候變化研究提供參考。

1 數(shù)據(jù)采集和處理方法

站點(diǎn)觀測(cè)數(shù)據(jù)主要來源于巴基斯坦氣象局、中國氣象局和美國國家環(huán)境預(yù)報(bào)中心的逐日降水、最高和最低氣溫站點(diǎn)數(shù)據(jù)。其中巴基斯坦地區(qū)有74 個(gè)站點(diǎn)，中國段有8 個(gè)站點(diǎn)。剔除缺測(cè)率超過50%的17 個(gè)站點(diǎn)，使用剩下的65 個(gè)氣象站點(diǎn)進(jìn)行數(shù)據(jù)制作。DEM 數(shù)據(jù)來自美國航空航天局的SRTM GRID 數(shù)據(jù)處得到的成品數(shù)據(jù)，使用ArcMap 軟件將DEM 數(shù)據(jù)進(jìn)行重采樣為0.25°×0.25°，并轉(zhuǎn)換成ANUSPLIN 軟件能夠識(shí)別的ASCII 格式數(shù)據(jù)。

1.1 ANUSPLIN 原理

薄盤樣條函數(shù)插值方法最早是Wahba 于1979 年提出，Hutchinson 等于1984 年對(duì)其改進(jìn)能夠適用于更大的數(shù)據(jù)集，Bates 等于1987 年將其進(jìn)一步拓展為局部薄盤光滑樣條法[19-20]。為了方便薄盤樣條函數(shù)法的使用，Hutchinson 等基于普通薄盤和局部薄盤樣條函數(shù)的插值理論，開發(fā)了專業(yè)氣候數(shù)據(jù)空間插值軟件ANUSPLIN，它除了可引入自變量外，還允許引入?yún)f(xié)變量（海拔、海岸線等）[21]。ANUSPLIN 軟件的核心是局部薄盤光滑樣條算法，其理論統(tǒng)計(jì)模型公式為：

其中，Zi為位于空間點(diǎn)i的因變量；f(xi)為關(guān)于xi的未知光滑函數(shù)；xi為獨(dú)立變量；bT為yi的p維系數(shù)；yi為p維獨(dú)立協(xié)變量；ei為隨機(jī)誤差；N為觀測(cè)值數(shù)量。當(dāng)式（1）缺少第一項(xiàng)f(xi)時(shí)，該統(tǒng)計(jì)模型簡化為簡單多元線性回歸模型，但是在ANUSPLIN 軟件的實(shí)際使用中不允許出現(xiàn)這種情況；當(dāng)式（1）缺少第二項(xiàng)bTyi時(shí)，即不存在協(xié)變量（p=0），該統(tǒng)計(jì)模型就簡化為普通的薄盤光滑樣條模型。式（1）中，函數(shù)f和系數(shù)bT通過最小二乘估計(jì)來確定：

其中，Jm(f)是函數(shù)f(xi)的粗糙度測(cè)度函數(shù)，為函數(shù)f的m階偏導(dǎo)（也稱為樣條次數(shù)）；ρ為正的光滑參數(shù)，主要用來平衡插值數(shù)據(jù)的保真度以及擬合曲面的粗糙度。當(dāng)ρ →0時(shí)，函數(shù)f為精確內(nèi)插式；當(dāng)ρ →+∞時(shí)，函數(shù)f為最小二乘多項(xiàng)式。在ANUSPLIN 軟件中通常以廣義交叉驗(yàn)證GCV 和最大似然法GML 的最小化來確定。GCV 的計(jì)算原理主要為逐個(gè)移除數(shù)據(jù)點(diǎn)，在同樣的ρ下利用其他數(shù)據(jù)點(diǎn)來估算該點(diǎn)的殘差，并且在ANUSPLIN 軟件中的log 日志文件中有記錄。

ANUSPLIN 軟件的log 日志文件中提供了一系列用于判別誤差來源和插值質(zhì)量的參數(shù)：觀測(cè)數(shù)據(jù)統(tǒng)計(jì)量（均值、方差、標(biāo)準(zhǔn)差等）、廣義交叉驗(yàn)證（GCV）、最大似然法誤差（GML）、擬合曲面參數(shù)的信號(hào)自由度（Signal）和剩余自由度（Error）、均方殘差（MSR）、光滑參數(shù)（RHO）、期望真實(shí)均方誤差（MSE）等。log 日志文件中的統(tǒng)計(jì)結(jié)果還給出了均方根殘差（RMSR，Root mean square residual）的數(shù)據(jù)點(diǎn)序列，可以用來控制數(shù)據(jù)質(zhì)量，檢驗(yàn)并消除原始數(shù)據(jù)在位置和數(shù)值上的錯(cuò)誤。

對(duì)于log 日志文件中數(shù)據(jù)擬合表面的結(jié)果，RHO 過小和Signal 大于觀測(cè)站點(diǎn)的一半或RHO 過大都表明在擬合過程中找不到最優(yōu)的光滑參數(shù)，說明數(shù)據(jù)點(diǎn)過于稀疏、存在短相關(guān)或擬合函數(shù)過于復(fù)雜，所選模型不適合用于插值，這些情況在ANUSPLIN 軟件的log 日志文件中以符號(hào)（*）標(biāo)出。ANUSPLIN 軟件插值過程中最佳模型的選擇標(biāo)準(zhǔn)：log 日志文件中GCV 或GML 最小，模型殘差比（MRR）或信噪比（SNR）最小，Signal 小于站點(diǎn)的一半，文件中無*號(hào)指示[22]。

1.2 數(shù)據(jù)處理流程

數(shù)據(jù)處理流程主要包括4 個(gè)部分（圖1）：原始數(shù)據(jù)輸入、數(shù)據(jù)處理、數(shù)據(jù)輸出（符合要求的數(shù)據(jù)格式）和空間插值（編寫批處理代碼）。輸入數(shù)據(jù)主要包括1961-2015 年中巴經(jīng)濟(jì)走廊地區(qū)氣象要素（日降水、日最高氣溫和日最低氣溫）站點(diǎn)數(shù)據(jù)、氣象臺(tái)站信息資料和DEM 數(shù)據(jù)。數(shù)據(jù)處理部分分別把氣象要素站點(diǎn)數(shù)據(jù)和DEM 數(shù)據(jù)處理成ANUSPLIN 軟件需要的數(shù)據(jù)格式。其中，將氣象要素站點(diǎn)數(shù)據(jù)樣本量小于50%的站點(diǎn)作為無效站點(diǎn)進(jìn)行剔除，用反距離加權(quán)法（IDW）對(duì)剩下站點(diǎn)的缺測(cè)值進(jìn)行插補(bǔ)，以保證插值過程和結(jié)果的可信度，然后輸出為ANUSPLIN 軟件需要的數(shù)據(jù)格式；另外，將中巴經(jīng)濟(jì)走廊地區(qū)DEM 數(shù)據(jù)進(jìn)行重采樣，根據(jù)插值目標(biāo)把空間分辨率重采樣為0.25°×0.25°，然后以ASCII 碼數(shù)據(jù)格式類型輸出。空間插值部分主要在ANUSPLIN 軟件中完成，通過編寫批處理腳本文件，進(jìn)行空間插值。

為保證每個(gè)擬合表面的插值精度和模型的穩(wěn)定性，并使之在連續(xù)的時(shí)間序列上具有可比性，在對(duì)3 個(gè)氣象要素（降水、日最高和最低氣溫）連續(xù)55 年逐日站點(diǎn)數(shù)據(jù)進(jìn)行曲面插值過程中，首先選取1979 年進(jìn)行實(shí)驗(yàn)（該年為平水年）。實(shí)驗(yàn)?zāi)Ｐ蜑楸”P樣條和局部薄盤樣條函數(shù)的6 個(gè)spline 模型（獨(dú)立變量、協(xié)變量和樣條次數(shù)多種組合，表1）。根據(jù)最佳模型的選擇標(biāo)準(zhǔn)，初步選出每個(gè)氣象要素的最優(yōu)待用模型，再用這些待用模型對(duì)不同氣象要素進(jìn)行連續(xù)55 年插值。對(duì)于個(gè)別模型不符的月份，利用殘差分析，剔除個(gè)別殘差較大的站點(diǎn)以使模型能夠使用。

圖1 數(shù)據(jù)處理流程圖

表1 6 個(gè)spline 模型詳細(xì)列表

對(duì)于降水?dāng)?shù)據(jù)、日最高/低氣溫?cái)?shù)據(jù)，用初定的6 個(gè)模型對(duì)1979 年的數(shù)據(jù)進(jìn)行實(shí)驗(yàn)，結(jié)果顯示選擇以高程作為協(xié)變量的三變量局部薄盤光滑樣條函數(shù)、樣條次數(shù)為2 的LLD2 模型能保證大部分插值結(jié)果最為精確。

2 數(shù)據(jù)樣本描述

1961-2015 年中巴經(jīng)濟(jì)走廊逐日氣象數(shù)據(jù)集共包含1961-2015 年60 264 個(gè)數(shù)據(jù)文件，命名方式為CPEC_XXX_YYYYMMDD。其中CPEC 為中巴經(jīng)濟(jì)走廊；XXX 為氣象要素，包括日降水量PRE、日最高氣溫TMAX 和日最低氣溫TMIN；YYYY 為數(shù)據(jù)年份；MM 為月份；DD 表示天。圖2 為中巴經(jīng)濟(jì)走廊地區(qū)1979 年8 月23 日降水?dāng)?shù)據(jù)，圖3 為中巴經(jīng)濟(jì)走廊地區(qū)1979 年8 月23 日最高（a）和最低（b）氣溫?cái)?shù)據(jù)。

圖2 中巴經(jīng)濟(jì)走廊地區(qū)1979 年8 月23 日降水

圖3 中巴經(jīng)濟(jì)走廊地區(qū)1979 年8 月23 日最高和最低氣溫

3 數(shù)據(jù)質(zhì)量控制和評(píng)估

為了驗(yàn)證本數(shù)據(jù)集的精度和可靠性，采用了研究區(qū)內(nèi)3 個(gè)未進(jìn)行插值的臺(tái)站的氣象要素作為驗(yàn)證數(shù)據(jù)（表2）。本數(shù)據(jù)的質(zhì)量控制并未涉及觀測(cè)站點(diǎn)的搬遷、觀測(cè)儀器變更和觀測(cè)規(guī)范變更等信息。

表2 驗(yàn)證站點(diǎn)

同時(shí)，本研究制作的格點(diǎn)化氣象數(shù)據(jù)集（以下簡稱CPEC-P、CPEC-T）與國際上較為常用的逐日降水?dāng)?shù)據(jù)集（表3）與逐日最高、最低氣溫?cái)?shù)據(jù)集（表4）進(jìn)行了對(duì)比。

表3 常用逐日降水?dāng)?shù)據(jù)集

表4 常用逐日最高、最低氣溫?cái)?shù)據(jù)集

對(duì)于降水?dāng)?shù)據(jù)（CPEC-P），從不同數(shù)據(jù)集評(píng)估結(jié)果（表5）和月降水?dāng)?shù)據(jù)驗(yàn)證散點(diǎn)圖（圖4）中可以看出，本研究制作的CPEC-P 能夠較好反映出真實(shí)的降水水平，其中在德拉·伊斯梅爾·汗站點(diǎn)評(píng)估結(jié)果最好，回歸系數(shù)為0.87，R2=0.73，均方根誤差RMSE=18.76 mm。對(duì)于德拉·伊斯梅爾·汗和曼迪·巴奧丁兩個(gè)氣象站，CPEC-P 與PGFMD 評(píng)估結(jié)果相一致，但分別低估了100 mm（德拉·伊斯梅爾·汗）與200 mm（曼迪·巴奧?。┮陨系慕邓?，而CHIRPS 與MSWEP 兩個(gè)數(shù)據(jù)集均整體低估了降水。對(duì)于吉德拉爾站點(diǎn)，CPEC-P 整體高估了該站點(diǎn)的降水，但其他數(shù)據(jù)集同樣不能夠很好反映出該站點(diǎn)的真實(shí)降水，這可能與該站點(diǎn)高程（1500 m）有關(guān)，高程較高的站點(diǎn)插值出來的結(jié)果誤差較大。

表5 不同降水?dāng)?shù)據(jù)集評(píng)估結(jié)果比較

圖4 月降水?dāng)?shù)據(jù)驗(yàn)證散點(diǎn)圖

對(duì)于最高、最低氣溫?cái)?shù)據(jù)（CPEC-T），從不同日最高、最低氣溫?cái)?shù)據(jù)集評(píng)估結(jié)果（表6）和月平均最高、最低氣溫驗(yàn)證散點(diǎn)圖（圖5）中可以看出，本研究制作的CPEC-T 比其他數(shù)據(jù)集能夠更好反映出站點(diǎn)的真實(shí)氣溫。對(duì)于德拉·伊斯梅爾·汗和曼迪·巴奧丁兩個(gè)氣象站，3 個(gè)數(shù)據(jù)集都能很好反映出站點(diǎn)的真實(shí)氣溫，CPEC 數(shù)據(jù)R2均在0.98 以上且RMSE 均在1℃以內(nèi)。對(duì)于吉德拉爾站點(diǎn)，3 個(gè)數(shù)據(jù)集都低估了該站點(diǎn)的氣溫，但3 個(gè)數(shù)據(jù)集評(píng)估結(jié)果表現(xiàn)為與觀測(cè)數(shù)據(jù)擬合程度較好，R2均在0.9以上，這可能因?yàn)檎军c(diǎn)的高程影響了插值的效果，造成數(shù)據(jù)結(jié)果的低估。

表6 不同日最高、最低氣溫?cái)?shù)據(jù)集評(píng)估結(jié)果比較

圖5 月平均最高(a-c)、最低氣溫(d-f)驗(yàn)證散點(diǎn)圖

4 數(shù)據(jù)使用方法和建議

本數(shù)據(jù)集為tif 文件格式，解壓后可使用Matlab 或ArcMap 軟件打開、顯示、查看、統(tǒng)計(jì)分析等。因?yàn)閿?shù)據(jù)量較大，建議使用Matlab 軟件進(jìn)行批處理，提取數(shù)據(jù)經(jīng)緯度代碼已上傳至網(wǎng)站。

致 謝

感謝巴基斯坦氣象局（PMD）、國家氣候中心、美國國家環(huán)境預(yù)報(bào)中心（GSOD）提供站點(diǎn)觀測(cè)數(shù)據(jù)。

數(shù)據(jù)作者分工職責(zé)

陳金雨（1998—），男，河南省信陽市人，碩士研究生，研究方向?yàn)闅庀笏臑?zāi)害風(fēng)險(xiǎn)評(píng)估。主要承擔(dān)工作：論文撰寫，數(shù)據(jù)質(zhì)量控制和評(píng)估。

陶輝（1981—），男，新疆昌吉市人，研究生學(xué)歷，副研究員，研究方向?yàn)闅夂蜃兓c風(fēng)險(xiǎn)評(píng)估。主要承擔(dān)工作：數(shù)據(jù)制作、評(píng)估整體思路的設(shè)計(jì)。

劉金平（1990—），男，河南省商丘市人，研究生學(xué)歷，講師，研究方向?yàn)槿蜃兓难h(huán)。主要承擔(dān)工作：數(shù)據(jù)整理、插值。