肖 旋 楊新凱
(上海師范大學信息與機電工程學院 上海 201418)
IPCC 第五次氣候評估得出結論為1880 年~2012年地表均溫呈現上升趨勢,漲幅大約0.85℃[1]。很多科學家持有觀點是地球變暖可能致使更多的極端天氣事件的發生。季飛、吳召華使用最新開發的非線性趨勢檢測方法-多維集合經驗模式分解來評估了全球陸地表面氣溫[2]。樸世龍與團隊使用北半球高緯度區域大氣CO2含量觀測數據、遙感數據、大陸生態系統和海域碳往復運動變化模型,并聯合大氣輸導,系統分析了1980年~2010年北半球高緯度地域生態系統春季碳匯因氣溫變化相應的動態改變及其機制[3]。探索溫度時空演變規律,對未來地區氣候預測、保護和改善自然環境及地區的可持續發展提供科學支撐。
針對氣候變化的數據:1)根據年代際的距平變化和ARIMA 時間序列分析的綜合結果可以得出,從1940 年~2010 年加拿大溫度的時間變化趨勢為升高、降低、升高的趨勢。2)選取ArcGIS 繪圖軟件,通過反距離權重插值法在地圖上繪制溫度的空間分布得出東西海岸平均溫度上漲了大約2℃左右,內陸地區大部分上升幅度都大于了4℃,北極三省的上漲幅度較大,基本都大于了6℃。特別北極圈內,漲幅超過了10℃。3)通過小波變換分析得出1980 年~2010 年地球海洋表面氣溫的改變趨向為變暖趨勢。大約10 年為一個全球海表溫度的變化周期,一個較強的擾動出現在約5 年~6 年的尺度上。
加拿大劃分為十個省和三個特別行政區,本文所用的加拿大各地天氣的歷史數據均來自加拿大氣象網站,預處理的數據是溫度、降雨量、降雪量等403 個氣象站點1940 年~2010 年月值數據,具體氣象站點如圖1 所示,落在研究區域的這些氣候數據符合均勻分布、完整性、典型性的可用性數據特征。

圖1 氣象站點分布
對加拿大403 個站點1940 年~2010 年各氣象資料進行處理加工,得到月平均值及年平均值。
首先將下載的數據采用Python 的pandas 進行數據清洗,主要包括以下流程。
在進行問題探究和收集數據中發現,數據集中數值缺失較多。經過數據庫檢測發現多數的變量都或多或少的缺失,而這些數據的缺失對加拿大的氣象站點的數據分析有較大的影響,所以這一部分的處理是非常重要的,不可忽視。對于這種缺失數據值,我們決定采用數據補齊的方式,首先可以先根據其他的相關字段進行預估,比如地區的變量等等。若上述方法無法實現填充,缺失的樣本值則通過相關性尋找最近距離的K個樣本值,再將其加權平均的結果值來補全數據[4]。
在處理數據前,通常需要規范化數據。因為一般在多指標評價機制中,每一個指標它的衡量性質并不完全相同,數量級和量綱一般有差異。若是這些指標之間的情況相差非常小,影響會較小;但如果這些指標之間相差水平很大,一般來說會增大數值大的指標在整體分析中的影響,而數值小的指標相對應的影響就會被減弱。本文在進行數據標準化的時候,主要采用了min-max 標準化,通常也將之稱為規范化方法,下面對其進行詳細介紹。
對序列x1,x2,x3,···,xn進行變換:

經過數據預處理后,我們先用雷達圖簡單的呈現全部數據的整體走向趨勢的概況,預先了解溫度的大體變化走向,我們將數據的csv 文件導入python 程序,經數據統計,用雷達圖逐年代展示溫度變化趨向。
如圖2 所示越往里溫度值越低,越往外溫度值越高,每年溫度12 個月份的變化趨勢,6 月~9 月溫度全年最高,12 月~2 月溫度全年最低。1940 年~2010年每月溫度的變化趨勢基本一致。

圖2 溫度變化雷達圖
本文運用的主要是根據魏鳳英老師的《現代氣候統計診斷與預測技術》一書中推薦的線性趨勢評估法,線性趨勢估計法是用一條合理的直線指示天氣變量與時間變量的聯系。[5]用xi表示樣本量為n的某一氣象變化值,用ti表示xi所對應的時間,構建xi與ti之間的一元線性回歸方程:

其中,a為回歸常數,b為回歸系數,a和b可用最小二乘法進行估量[5]。
對于觀測數據xi及相應的時間ti,回歸常數a和系數b的最小二乘法估值為

時間ti與變量xi之間的相關系數借助回歸系數b與相關系數的關聯性計算:




移動(滑動)平均可用來模擬趨勢,在訓練量為n時,給定時間序列x,此處為氣候變化序列x,取其平均值描繪分析波動趨勢[5],表達其滑動平均序列是:

式(8)中,滑動長度為k,j=1,2,…,n-k+1。觀察七十年的時間序列以及分析周期,將滑動窗口k設置為10,即觀測氣候序列的變化趨勢是通過10年一個滑動分離出變化趨勢以探究規律。
ARIMA 模型(Autoregressive Integrated Moving Average model),也就是所謂的差分整合移動平均自回歸模型(移動也可稱作滑動),預測時間序列的趨勢[10]。ARIMA(p,d,q)里,p 是自回歸項數,用來解釋時間序列模型中數據自身的滯后數,即AR項,d 是要達到平穩序列需要做的差分次數,即Integrated 項,q 是滑動平均項數,表示該模型在預測誤差范圍內產生的滯后數,也就是MA 項[10]。“差分”是關鍵步驟。
此模型被選擇來預測溫度隨時間推移的變化序列。由ARMA(p,q)模型的延伸擴充獲得ARIMA(p,d,q)模型。ARIMA(p,d,q)模型表達式是:

1)小波概念
小波分析來自于Fourier,并且實現了Fourier所不能實現的一些艱難問題,適合于處理非穩定信號[12]。小波分析可以用于時間-尺度分析和多分辨率分析。其波形同時具有衰減性和波動性,定義公式如下:

2)小波變換
小波變換是將函數和數據等切分成不同的片段,然后分析相應的尺度分量,小波變換有兩種變換方式,一種是通過連續形式實現變換,另一種則是采用離散手段達到變換;其中,離散小波又有著多種不同形式,比如:正交、雙正交小波。[13]由于海表溫度的數據有一定的連續震蕩的特征,所以采用一維連續小波變換分析其變化特點。一維連續小波變化的公式如下:

3)小波方差
小波方差即小波變換系數的平方值在b域上的積分[13],其計算公式為

小波方差圖是小波方差值伴隨時間變化形成的曲線,該圖可以解釋在不同時間尺度上小波能量的分布狀況。
距平是顯示氣候要素偏離正常值的最常見的量。數據集中的某個數xi與平均數之間的差值就是距平,即:

其中,一組氣候變量數據x1,x2,xn與其平均數的差組成距平序列。

在氣象研究中,經常選取距平序列來替換氣候指標本身的觀測值。
將數據按照時間序列輸出距平值,反映了10年一個單位時間尺度上及整體的溫度變化特征[12],從圖3 能夠發現,加拿大各地區從1940 年起,70 年來每10 年單位時間尺度之間差異性較大,整體上溫度有上升趨勢。

圖3 溫度變化趨勢
通過統計各月的變化幅度,即圖4可知,1月份每單位時間尺度實際地表溫度與1940 年~2010 年地表溫度均值相差度數分別是-0.64℃、0.76℃、-0.54℃、-1.66℃、-0.29℃、0.99℃以 及1.31℃,表明20 世紀40 年代至80 年代地表溫度偏低,90 年代以及21 世紀00 年代地表溫度偏高,總體表現為“升高-降低-升高”的變化趨勢,與2、4、9月份各年代際地表溫度變化趨勢相似,都是70 年代氣溫偏低最嚴重。3 月的溫度變化趨勢“下降-上升”,類似于5 月份各年代際地表溫度變化趨勢,60 年代氣溫偏低最嚴重。6、10、11 月總體表現為“上升-下降-上升”的變化趨勢,80 年代溫度偏低最嚴重。7、8、12 月總體表現為“上升-下降-上升”的變化趨勢,60 年代溫度偏低最嚴重。除了12 月之外,溫度升高最大的都是00年代,12月溫度升高最大的是90年代。

圖4 各月實際地表氣溫和地表氣溫均值
經過觀察分析,我們可以看出月份之間存在相似性,于是我們做了相關性分析,采用斯皮爾曼相關性分析做了進一步的相關性矩陣,具體如圖5 所示。

圖5 相關矩陣
相關性系數分布在0~1區間,越靠近1,說明相關程度越高,越靠近0,相關程度越低,從圖中可以看出,5 月和10 月相關關系最低,相關系數低至0.021。
通過ARIMA 模型對1940 年~2010 年的數據進行STL 分解,得出1940 年~2010 年隨時間的溫度變化趨勢,如圖6所示。

圖6 1940年~2010年隨時間的溫度變化趨勢
根據分解結果Trend 圖的結論,溫度仍然是出現先略微下降再升高的趨勢。
1)反距離權重插值法
反距離權重插值法是常用的、易于理解的、方便實用的、比較好的空間插值方式之一,是一種精準插值,源于相近相似原理,具體表現為距離小,屬性越接近,類似程度越高,相反,間距遠,特性差別較大[14]。基于插值點與樣本點之間的距離賦予權重求平均。間隔越小給予該點的權重值就會越大。經過查閱文獻資料,我們選取ArcGIS 繪圖軟件,用反距離權重插值法在地圖上繪制溫度的空間分布[15]。
2)溫度空間分布分析
加拿大幅員遼闊,主要城市分布在5 個氣候區:太平洋地區(西海岸區)、平原地區、中部地區(大湖-圣勞倫斯蒂地區)、大西洋地區和北方地區,各個城市氣候有著不一樣的特點[16]。
我們可以從圖7 中分析溫度空間變化趨勢,在這段時間中,加拿大所有地區的平均氣溫都有著不同程度的上升。上升比較慢的是東西海岸,大約在2℃~3℃之間。結合相關地理知識易于理解,靠近海洋的氣候本來就是比較溫和。

圖7 溫度空間變化趨勢
內陸地區與沿海地區差別較大,大部分平均溫度上升幅度都大于了4℃,北極三省的上漲幅度最大,基本都大于了6℃。顯而易見,特別突出明顯的是北極圈內,漲幅甚至超過了10℃。通過新聞可知,溫室氣體增加,由此導致溫度上漲,成為了主要原因。
首先根據收集的數據對海表溫度進行初步展示,并且對1981 年~2010 年一共30 年的數據進行小波分析,全球海表溫度變化趨勢、距平值及平均溫度分布如圖8、圖9、圖10所示。

圖8 全球海表溫度變化趨勢

圖9 全球海表溫度距平值曲線圖

圖10 全球海表平均溫度分布圖
全球海表溫度在近30 年間總體呈現上升的趨勢,如圖8 所示,1981 年海表溫度13.80 左右,在2006年出現最高溫度,溫度大約是13.99,增幅約為0.19。1981 年之后全球海表的距平值開始上升,1994 年到達第一個高峰,距平值約為0.35,在2010年出現最大距平。從圖8 中可以看出,海表溫度出現緩慢增長趨勢,圖9 呈現海表溫度的距平在波動曲折中上升,圖10 直觀顯示全球海表平均溫度分布[17]。
生成函數選擇墨西哥小帽小波,用一維相關連續小波對30 年這個時間段的數據集做變換操作,得到的灰度圖直觀顯示小波變換系數絕對值的大小程度[12]。橫軸為時間變換,縱軸為溫度變化,圖中越高的值,表征此處有越大的小波變換系數絕對值,意味著該時間段下海表溫度變化越明顯;小波變換系數為0 的坐標點是該尺度下海表氣溫變化的拐點[18]。為了準確地剖析30 年間海表溫度的變化和突變特點,繪制了絕對值灰度圖所對應的等值線圖。通過等值線圖11 可以得出,全球海表溫度的變化周期約為10 年一個周期,在約5 年~6 年的尺度上有一個較強的擾動。

圖11 等值線圖
首先借助距平值揭示地表溫度隨時間變化的動態規律,得到12 個月的地表溫度變化。通過ARIMA對1940年~2010年的數據進行STL分解,得到整體的氣溫趨勢是略下降再上升的趨勢[19]。對空間的溫度趨勢分析時,采用了反距離權重插值法和基于小波分析的方式對海表溫度規律探究,能清楚的用數據和圖表展示出1980 年~2010 年全球海洋表面溫度的變化趨勢可以得出全球海表溫度的變化趨勢為變暖的趨勢。
由于收集到的數據中指標很多,初始樣本量更是高達十萬條,給數據量化處理工作帶來很大困難,因此在篩選過程中會遺漏部分可能存在重要作用的變量信息,并且我們在求解權重的時候,可能方法不盡合理,會對結果產生一定程度的影響。最后,由于指標提取后,丟棄一些指標會造成小信息的丟失,所以得到的結論可能會存在一些可接受的誤差。
收集的數據庫中數據缺失很多,而這些缺失的數據對加拿大的氣象站點的數據分析也是有影響的,雖然我們在進行建模和數據分析的時候,采用了K 最近距離鄰法,對數據進行了填充,但是對最后的結果可能造成一定的誤差。
在分析空間和時間的溫度變化趨勢時,影響地表溫度變化的因素中沒有包括地形因素。由于氣候的形成往往受到地形的影響,地表溫度的空間分布因此發生了變化[20]。因此,若獲取到更詳細的地形數據集加入到考慮因素中分析時空分布的變化趨勢,結果的準確度會顯著提高。