張學(xué)良，王華東，肖鵬峰，鄭子賢，馮學(xué)智

1.南京大學(xué)地理與海洋科學(xué)學(xué)院，自然資源部國土衛(wèi)星遙感應(yīng)用重點(diǎn)實(shí)驗(yàn)室，江蘇省地理信息技術(shù)重點(diǎn)實(shí)驗(yàn)室，南京 210023

引 言

積雪是冰凍圈最活躍的存在形式，也是全球氣候系統(tǒng)的重要組成部分，在全球地表能量平衡、水文循環(huán)和人類活動中占據(jù)重要地位[1-4]。全球約有98%的季節(jié)性積雪覆蓋在北半球，覆蓋面積約占冬季北半球陸地面積的 50%[5]，積雪融水匯集河流，不僅可為全球六分之一的人口提供淡水資源，同時(shí)影響河流下游地區(qū)的農(nóng)業(yè)和生態(tài)系統(tǒng)[6]。在中國，多年平均積雪面積約占陸地面積的 94%，尤其是在我國的三大季節(jié)性積雪穩(wěn)定區(qū)[7]。在全球變暖背景下，積雪對氣候變化具有高度敏感性和重要反饋?zhàn)饔肹8-9]，積雪及其屬性發(fā)生變化，將造成冰凍圈和其他圈層發(fā)生改變，繼而影響生態(tài)環(huán)境及人類社會生活。

最初的積雪研究主要依賴有限的人工觀測和地面臺站數(shù)據(jù)，數(shù)據(jù)準(zhǔn)確性高，但難以實(shí)現(xiàn)大規(guī)模積雪要素動態(tài)監(jiān)測，尤其是在高海拔的山區(qū)。遙感以其大范圍、多尺度、多時(shí)相等特點(diǎn)，逐漸成為積雪調(diào)查研究的重要手段。遙感與地面調(diào)查手段的聯(lián)合為了解和掌握中國積雪特性的時(shí)空分布特征提供了有力保障和支撐。而積雪物候、積雪深度/雪水當(dāng)量、積雪反照率、積雪密度是積雪的重要參數(shù)，與氣候變化、水資源管理和自然災(zāi)害等密切相關(guān)[10-12]。

積雪要素專題地圖可直觀、定量表達(dá)積雪的時(shí)空分布，已有的專題地圖包括積雪分布[13]、積雪穩(wěn)定性[7]、積雪災(zāi)害[13-14]、積雪覆蓋日數(shù)[15]、最大積雪深度[15]等。在國家科技基礎(chǔ)資源調(diào)查專項(xiàng)的支持下，基于“中國積雪特性及分布調(diào)查”項(xiàng)目生產(chǎn)的中國區(qū)域不同積雪要素逐日遙感產(chǎn)品和地面調(diào)查數(shù)據(jù)，編制中國積雪系列圖件，包括4個(gè)子集，分別為1980-2020年中國積雪日數(shù)（積雪物候）圖集、1980-2020年中國積雪深度/雪水當(dāng)量圖集、2000-2020年中國積雪反照率圖集和2013-2020年中國積雪密度圖集，共計(jì)1244張電子地圖及相應(yīng)制圖數(shù)據(jù)。從不同時(shí)間頻度（逐月、逐年、多年逐月平均等）定量表達(dá)積雪要素的空間分布及其年內(nèi)、年際變化，對逐日積雪要素?cái)?shù)據(jù)的進(jìn)一步整合歸納與再生產(chǎn)，為認(rèn)識和掌握中國積雪特性的空間分布差異和時(shí)間變化特征提供數(shù)據(jù)支撐。

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

1.1 數(shù)據(jù)來源

從國家冰川凍土沙漠科學(xué)數(shù)據(jù)中心網(wǎng)站（http://www.ncdc.ac.cn/）獲取的中國1980-2020年積雪面積5 km逐日無云產(chǎn)品[16]、中國1980-2020年雪水當(dāng)量25 km逐日產(chǎn)品[17]和中國2000-2020年積雪反照率逐日產(chǎn)品[18]，分別編制1980-2020年中國積雪日數(shù)圖集、1980-2020年中國積雪深度/雪水當(dāng)量圖集和2000-2020年中國積雪反照率圖集。

目前關(guān)于積雪密度的產(chǎn)品有積雪密度再分析數(shù)據(jù)集，它是將各種觀測數(shù)據(jù)同化到預(yù)測模型中，提供積雪密度網(wǎng)格數(shù)據(jù)集[19]，例如歐洲的ERA-5產(chǎn)品，但與中國站點(diǎn)數(shù)據(jù)對比后發(fā)現(xiàn)驗(yàn)證精度較低。本積雪密度圖集利用2013-2020年中國氣象站觀測數(shù)據(jù)（http://data.cma.cn），2017-2019年中國典型積雪區(qū)線路積雪觀測數(shù)據(jù)集[20]（http://www.ncdc.ac.cn/），和 ERA-5再分析數(shù)據(jù)集（https://cds.climate.copernicus.eu），包括風(fēng)速、氣溫、融雪量、降雪量、積雪表層溫度、積雪蒸發(fā)量、降水量、地表氣壓、葉面積指數(shù)，以及數(shù)字高程模型（http://www.gscloud.cn/），編制了2013-2020年中國積雪密度圖集。

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

1.2.1 積雪日數(shù)

利用中國1980-2020年積雪面積5 km逐日無云遙感產(chǎn)品，提取逐年積雪日數(shù)（積雪持續(xù)日數(shù)、積雪初日、積雪終日），其中積雪持續(xù)日數(shù)在本文中由積雪覆蓋日數(shù)指代，指一個(gè)積雪年內(nèi)（9月1日-翌年8月31日）觀測到的累計(jì)積雪總天數(shù)。計(jì)算公式如下：

SCD(snow cover days)表示積雪持續(xù)日數(shù)，N表示一個(gè)完整的積雪年(9月1日-翌年8月31日)，Si表示第i天的逐日無云積雪面積產(chǎn)品，Si=1表示存在積雪，Si=0表示無雪。

為了避免短暫積雪的影響，積雪初日定義為在一個(gè)積雪年內(nèi)，第一個(gè)連續(xù)5天積雪的第一天；積雪終日定義為在一個(gè)積雪年內(nèi)，最后一個(gè)連續(xù)5天積雪的最后一天。如果在一個(gè)積雪年內(nèi)不存在連續(xù)5天積雪，則表示沒有積雪初日和積雪終日，計(jì)算公式如下：

其中SCOD(snow cover onset date)表示積雪初日，SCED(snow cover end date)表示積雪終日，t為儒略日，N表示一個(gè)完整積雪年，St表示第t天的逐日無云積雪面積產(chǎn)品，St=1表示存在積雪，St=0表示無雪，k表示連續(xù)5天。

利用逐年積雪持續(xù)日數(shù)、積雪初日、積雪終日，計(jì)算多年平均積雪持續(xù)日數(shù)、積雪初日和積雪終日，以及利用一元線性回歸進(jìn)行趨勢傾向估計(jì)積雪日數(shù)的多年變化率，并進(jìn)行α=0.05的顯著性檢驗(yàn)。

對得到的所有制圖數(shù)據(jù)進(jìn)行分級和制圖符號設(shè)計(jì)，采用分層設(shè)色法表示積雪日數(shù)，獲得積雪日數(shù)專題地圖。積雪日數(shù)制圖流程如圖1所示。

圖1 積雪日數(shù)制圖流程圖Figure 1 Flow chart of snow phenology mapping

1.2.2 積雪深度/雪水當(dāng)量

利用中國1980-2020年雪水當(dāng)量25 km逐日遙感產(chǎn)品，計(jì)算逐年平均、逐月平均和多年逐月平均積雪深度/雪水當(dāng)量，再根據(jù)逐年平均和逐月平均積雪深度/雪水當(dāng)量，分別計(jì)算得到積雪深度/雪水當(dāng)量多年變化率和多年逐月變化率，并進(jìn)行α=0.05的顯著性檢驗(yàn)。利用逐月平均和多年逐月平均積雪深度/雪水當(dāng)量，計(jì)算得到積雪深度/雪水當(dāng)量多年逐月標(biāo)準(zhǔn)偏差。

對得到的所有制圖數(shù)據(jù)進(jìn)行分級和制圖符號設(shè)計(jì)，采用分層設(shè)色法表示積雪深度和雪水當(dāng)量，獲得積雪深度/雪水當(dāng)量專題地圖。積雪深度/雪水當(dāng)量制圖流程如圖2所示。

圖2 積雪深度/雪水當(dāng)量制圖流程圖Figure 2 Flow chart of snow depth / snow water equivalent mapping

1.2.3 積雪反照率

利用中國2000-2020年積雪反照率逐日遙感產(chǎn)品，計(jì)算逐旬（上旬、中旬、下旬）平均、逐月平均和多年逐月（9月-翌年5月）平均積雪反照率，然后利用逐月平均和多年逐月平均積雪反照率，計(jì)算積雪反照率多年逐月標(biāo)準(zhǔn)偏差圖。

對得到的所有制圖數(shù)據(jù)進(jìn)行分級和制圖符號設(shè)計(jì)，采用分層設(shè)色法表示積雪反照率圖集，得到積雪反照率專題地圖。積雪反照率制圖流程如圖3所示。

圖3 積雪反照率制圖流程圖Figure 3 Flow chart of snow albedo mapping

1.2.4 積雪密度

首先建立時(shí)空加權(quán)神經(jīng)網(wǎng)絡(luò)模型得到2013-2020年逐日積雪密度，模型輸入變量包括積雪特性、氣象條件和地形特征，輸出變量為積雪密度。輸入變量的數(shù)據(jù)來自于遙感數(shù)據(jù)和再分析數(shù)據(jù)。使用的遙感數(shù)據(jù)包括逐日積雪深度、積雪反照率、積雪面積，以及地形特征數(shù)據(jù)（海拔、坡度、坡向），利用逐日積雪面積提取積雪累積日數(shù)，即從某個(gè)積雪日向前延伸到積雪開始出現(xiàn)的累積日數(shù)。使用的ERA-5再分析數(shù)據(jù)包括風(fēng)速、氣溫、融雪量、降雪量、積雪表層溫度、積雪蒸發(fā)量、降水量、地表氣壓、葉面積指數(shù)。地面觀測積雪密度數(shù)據(jù)包括氣象站觀測數(shù)據(jù)和地面調(diào)查數(shù)據(jù)。

時(shí)空加權(quán)神經(jīng)網(wǎng)絡(luò)模型主要包括兩部分[21]：一是時(shí)空加權(quán)，二是廣義回歸神經(jīng)網(wǎng)絡(luò)。時(shí)空加權(quán)考慮積雪密度的時(shí)空異質(zhì)性，可以得到預(yù)測點(diǎn)和其周圍樣本點(diǎn)的時(shí)空權(quán)重，而廣義回歸神經(jīng)網(wǎng)絡(luò)考慮積雪密度與影響因素之間的非線性關(guān)系，可以得到預(yù)測點(diǎn)和其周圍樣本點(diǎn)的影響因素權(quán)重，時(shí)空加權(quán)神經(jīng)網(wǎng)絡(luò)模型計(jì)算公式如下：

其中，dST表示時(shí)空距離，(j,k,t)表示樣本點(diǎn)s和預(yù)測點(diǎn)p的空間位置和時(shí)間，φ是調(diào)節(jié)空間和時(shí)間距離重要性的參數(shù)，hST表示自適應(yīng)帶寬，即不同預(yù)測點(diǎn)選擇相同數(shù)量的周圍樣本參與模型計(jì)算，而不是確定一個(gè)固定范圍（固定帶寬），這可以解決樣本分布不均勻問題，最后根據(jù)公式（4）和（5）得到不同樣本點(diǎn)和預(yù)測點(diǎn)的時(shí)空權(quán)重WGT。

其中，x表示樣本點(diǎn)s和預(yù)測點(diǎn)p的影響因素，D表示歐幾里得距離，spread是控制擬合函數(shù)平滑度的參數(shù)，N表示預(yù)測點(diǎn)周圍樣本點(diǎn)的數(shù)量，yi表示樣本點(diǎn)i的真實(shí)積雪密度值。根據(jù)公式（4）和（5）計(jì)算得到第i個(gè)樣本點(diǎn)和預(yù)測點(diǎn)的時(shí)空權(quán)重，公式（6）計(jì)算得到第i個(gè)樣本點(diǎn)和預(yù)測點(diǎn)的廣義回歸神經(jīng)網(wǎng)絡(luò)權(quán)重，最后根據(jù)公式（7）得到預(yù)測點(diǎn)的積雪密度(snow density)估計(jì)值。

分別提取觀測站點(diǎn)的積雪密度及其對應(yīng)位置上的影響因素作為樣本，共得到16933個(gè)樣本。將所有樣本隨機(jī)平均分成十份，采用十折交叉驗(yàn)證的方法確定最優(yōu)模型，估算2013-2020年逐日積雪密度。進(jìn)而計(jì)算逐月（9月-翌年5月）平均和多年逐月平均積雪密度，以及積雪密度多年逐月標(biāo)準(zhǔn)偏差。

對得到的所有制圖數(shù)據(jù)進(jìn)行分級和制圖符號設(shè)計(jì)，采用分層設(shè)色法表示積雪密度，獲得積雪密度專題地圖。積雪密度制圖流程如圖4所示。

圖4 積雪密度制圖流程圖Figure 4 Flow chart of snow density mapping

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

中國積雪特性時(shí)空分布電子地圖集共有四種積雪要素圖集，分別為積雪日數(shù)圖集、積雪深度/雪水當(dāng)量圖集、積雪反照率圖集和積雪密度圖集，共生產(chǎn)1244張jpg格式的電子地圖和1244張geotiff格式的制圖數(shù)據(jù)，每個(gè)圖集的詳細(xì)信息如表1所示。

表1 圖集介紹Table 1 Introduction of electronic atlas

2.1 積雪日數(shù)圖集

積雪日數(shù)圖集共包含5個(gè)子圖集：（1）1980-2020年中國逐年積雪持續(xù)日數(shù)圖；（2）1980-2020年中國逐年積雪初日圖；（3）1980-2020年中國逐年積雪終日圖；（4）1980-2020年中國平均積雪日數(shù)圖；（5）1980-2020年中國積雪日數(shù)變化率圖。

命名規(guī)則為：（1）yyyy-yyyy年中國積雪持續(xù)日數(shù)圖；（2）yyyy-yyyy年中國積雪初日圖；（3）yyyy-yyyy年中國積雪終日圖；（4）1980-2020年中國平均積雪日數(shù)要素圖；（5）1980-2020年中國積雪日數(shù)要素變化率圖。其中圖5和圖6展示了中國積雪持續(xù)日數(shù)、積雪初日、積雪終日的多年平均和多年變化率圖。

圖5 多年平均積雪日數(shù)圖（審圖號：GS(2022)922號）Figure 5 Map of average snow phenology

圖6 積雪日數(shù)多年變化率圖（審圖號：GS(2022)922號）Figure 6 Map of snow phenology change rate

2.2 積雪深度/雪水當(dāng)量圖集

積雪深度/雪水當(dāng)量圖集共包括8個(gè)子圖集：（1）1980-2020年中國逐月平均積雪深度圖；（2）1980-2020年中國逐月平均雪水當(dāng)量圖；（3）1980-2020年中國積雪深度逐月標(biāo)準(zhǔn)偏差圖；（4）1980-2020年中國雪水當(dāng)量逐月標(biāo)準(zhǔn)偏差圖；（5）1980-2020年中國積雪深度逐月變化率圖；（6）1980-2020年中國雪水當(dāng)量逐月變化率圖；（7）1980-2020年中國積雪深度變化率圖；（8）1980-2020年中國雪水當(dāng)量變化率圖。

命名規(guī)則為：（1）1980-2020年中國逐月平均積雪深度圖_mm月；（2）1980-2020年中國逐月平均雪水當(dāng)量圖_mm月；（3）1980-2020年中國積雪深度逐月標(biāo)準(zhǔn)偏差圖_mm月；（4）1980-2020年中國雪水當(dāng)量逐月標(biāo)準(zhǔn)偏差圖_mm月；（5）1980-2020年中國積雪深度逐月變化率圖_mm月；（6）1980-2020年中國雪水當(dāng)量逐月變化率圖_mm月；（7）1980-2020年中國積雪深度變化率圖；（8）1980-2020年中國雪水當(dāng)量變化率圖。其中圖7和圖8分別展示中國積雪深度和雪水當(dāng)量不同統(tǒng)計(jì)特征的一張示例圖。

2.3 積雪反照率圖集

積雪反照率圖集共包含4個(gè)子圖集：（1）中國逐旬平均積雪反照率圖；（2）中國逐月平均積雪反照率圖；（3）2000-2020年中國逐月平均積雪反照率圖；（4）2000-2020年中國積雪反照率逐月標(biāo)準(zhǔn)偏差圖。

命名規(guī)則為：（1）中國逐旬平均積雪反照率圖_yyyy年mm月t旬；（2）中國逐月平均積雪反照率圖_yyyy年mm月；（3）2000-2020年中國逐月平均積雪反照率圖_mm月；（4）2000-2020年中國積雪反照率逐月標(biāo)準(zhǔn)偏差圖_mm月。其中圖9展示了中國積雪反照率圖集中每個(gè)子圖集的一張示例圖。

圖7 積雪深度圖集示例圖（審圖號：GS(2022)922號）Figure 7 Sample map of snow depth atlas

2.4 積雪密度圖集

積雪密度圖集共包括2個(gè)子圖集：（1）2013-2020年中國逐月平均積雪密度圖；（2）2013-2020年中國積雪密度逐月標(biāo)準(zhǔn)偏差圖。

命名規(guī)則為：（1）2013-2020年中國逐月平均積雪密度圖_mm月；（2）2013-2020年中國積雪密度逐月標(biāo)準(zhǔn)偏差圖_mm月。其中圖10為中國積雪密度圖集子圖集的示例圖。

圖8 雪水當(dāng)量圖集示例圖（審圖號：GS(2022)922號）Figure 8 Sample map of snow water equivalent atlas

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

積雪反照率圖集使用中國2000-2020年積雪反照率逐日產(chǎn)品編制，逐日積雪反照率產(chǎn)品的均方根誤差（root mean squared error，RMSE）達(dá)到0.11[18]，積雪深度/雪水當(dāng)量圖集使用中國1980-2020年雪水當(dāng)量25 km逐日產(chǎn)品編制，逐日雪水當(dāng)量產(chǎn)品的RMSE在15 mm以內(nèi)[17]。

利用中國 1980-2020年積雪面積逐日無云產(chǎn)品[16]提取得到的積雪日數(shù)遙感數(shù)據(jù)，經(jīng)與 7282個(gè)氣象站觀測積雪日數(shù)（雪深≥2 cm視為積雪日）數(shù)據(jù)進(jìn)行對比驗(yàn)證，結(jié)果顯示：積雪持續(xù)日數(shù)平均日數(shù)誤差為3.5天，積雪初日和積雪終日的平均日數(shù)誤差分別為3.8天和5.4天。

圖9 積雪反照率圖集示例圖（審圖號：GS(2022)922號）Figure 9 Sample map of snow albedo atlas

圖10 積雪密度圖集示例圖（審圖號：GS(2022)922號）Figure 10 Sample map of snow density atlas

經(jīng)十折交叉驗(yàn)證，積雪密度的決定系數(shù)（R2）、平均絕對誤差（mean absolute error，MAE）、RMSE分別為0.517、0.028 g/cm3和0.043 g/cm3。在三大積雪區(qū)中，北疆的估計(jì)精度最高，R2達(dá)到0.626，其次為東北-內(nèi)蒙古（R2=0.566）和西藏-青海（R2=0.502），其MAE和RMSE都在0.06 g/cm3以下。

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

本數(shù)據(jù)集共享了長時(shí)間序列積雪要素的時(shí)空分布特征及其變化規(guī)律，是對逐日積雪要素?cái)?shù)據(jù)的進(jìn)一步整合歸納與再生產(chǎn)。數(shù)據(jù)文件按積雪要素屬性文件存儲，分別采用geotiff格式存儲制圖數(shù)據(jù)，jpg格式存儲電子地圖。用戶可按照實(shí)際情況選擇不同時(shí)間頻度、不同積雪屬性的文件下載。

致 謝

感謝中國氣象局（http://data.cma.cn）提供的 2013-2020年中國氣象站觀測數(shù)據(jù)，以及 ERA-5（https://cds.climate.copernicus.eu）和地理空間數(shù)據(jù)云（http://www.gscloud.cn/）提供的再分析數(shù)據(jù)集和數(shù)字高程模型數(shù)據(jù)。