基于ARIMA-小波分析的加拿大氣溫時(shí)空變化趨勢(shì)研究*

肖 旋 楊新凱

（上海師范大學(xué)信息與機(jī)電工程學(xué)院 上海 201418）

1 引言

IPCC 第五次氣候評(píng)估得出結(jié)論為1880 年～2012年地表均溫呈現(xiàn)上升趨勢(shì)，漲幅大約0.85℃［1］。很多科學(xué)家持有觀點(diǎn)是地球變暖可能致使更多的極端天氣事件的發(fā)生。季飛、吳召華使用最新開發(fā)的非線性趨勢(shì)檢測(cè)方法-多維集合經(jīng)驗(yàn)?zāi)Ｊ椒纸鈦碓u(píng)估了全球陸地表面氣溫［2］。樸世龍與團(tuán)隊(duì)使用北半球高緯度區(qū)域大氣CO2含量觀測(cè)數(shù)據(jù)、遙感數(shù)據(jù)、大陸生態(tài)系統(tǒng)和海域碳往復(fù)運(yùn)動(dòng)變化模型，并聯(lián)合大氣輸導(dǎo)，系統(tǒng)分析了1980年～2010年北半球高緯度地域生態(tài)系統(tǒng)春季碳匯因氣溫變化相應(yīng)的動(dòng)態(tài)改變及其機(jī)制［3］。探索溫度時(shí)空演變規(guī)律，對(duì)未來地區(qū)氣候預(yù)測(cè)、保護(hù)和改善自然環(huán)境及地區(qū)的可持續(xù)發(fā)展提供科學(xué)支撐。

針對(duì)氣候變化的數(shù)據(jù)：1）根據(jù)年代際的距平變化和ARIMA 時(shí)間序列分析的綜合結(jié)果可以得出，從1940 年～2010 年加拿大溫度的時(shí)間變化趨勢(shì)為升高、降低、升高的趨勢(shì)。2）選取ArcGIS 繪圖軟件，通過反距離權(quán)重插值法在地圖上繪制溫度的空間分布得出東西海岸平均溫度上漲了大約2℃左右，內(nèi)陸地區(qū)大部分上升幅度都大于了4℃，北極三省的上漲幅度較大，基本都大于了6℃。特別北極圈內(nèi)，漲幅超過了10℃。3）通過小波變換分析得出1980 年～2010 年地球海洋表面氣溫的改變趨向?yàn)樽兣厔?shì)。大約10 年為一個(gè)全球海表溫度的變化周期，一個(gè)較強(qiáng)的擾動(dòng)出現(xiàn)在約5 年～6 年的尺度上。

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

加拿大劃分為十個(gè)省和三個(gè)特別行政區(qū)，本文所用的加拿大各地天氣的歷史數(shù)據(jù)均來自加拿大氣象網(wǎng)站，預(yù)處理的數(shù)據(jù)是溫度、降雨量、降雪量等403 個(gè)氣象站點(diǎn)1940 年～2010 年月值數(shù)據(jù)，具體氣象站點(diǎn)如圖1 所示，落在研究區(qū)域的這些氣候數(shù)據(jù)符合均勻分布、完整性、典型性的可用性數(shù)據(jù)特征。

圖1 氣象站點(diǎn)分布

對(duì)加拿大403 個(gè)站點(diǎn)1940 年～2010 年各氣象資料進(jìn)行處理加工，得到月平均值及年平均值。

3 數(shù)據(jù)準(zhǔn)備

首先將下載的數(shù)據(jù)采用Python 的pandas 進(jìn)行數(shù)據(jù)清洗，主要包括以下流程。

3.1 缺失值處理

在進(jìn)行問題探究和收集數(shù)據(jù)中發(fā)現(xiàn)，數(shù)據(jù)集中數(shù)值缺失較多。經(jīng)過數(shù)據(jù)庫檢測(cè)發(fā)現(xiàn)多數(shù)的變量都或多或少的缺失，而這些數(shù)據(jù)的缺失對(duì)加拿大的氣象站點(diǎn)的數(shù)據(jù)分析有較大的影響，所以這一部分的處理是非常重要的，不可忽視。對(duì)于這種缺失數(shù)據(jù)值，我們決定采用數(shù)據(jù)補(bǔ)齊的方式，首先可以先根據(jù)其他的相關(guān)字段進(jìn)行預(yù)估，比如地區(qū)的變量等等。若上述方法無法實(shí)現(xiàn)填充，缺失的樣本值則通過相關(guān)性尋找最近距離的K個(gè)樣本值，再將其加權(quán)平均的結(jié)果值來補(bǔ)全數(shù)據(jù)［4］。

3.2 數(shù)據(jù)標(biāo)準(zhǔn)化

在處理數(shù)據(jù)前，通常需要規(guī)范化數(shù)據(jù)。因?yàn)橐话阍诙嘀笜?biāo)評(píng)價(jià)機(jī)制中，每一個(gè)指標(biāo)它的衡量性質(zhì)并不完全相同，數(shù)量級(jí)和量綱一般有差異。若是這些指標(biāo)之間的情況相差非常小，影響會(huì)較小；但如果這些指標(biāo)之間相差水平很大，一般來說會(huì)增大數(shù)值大的指標(biāo)在整體分析中的影響，而數(shù)值小的指標(biāo)相對(duì)應(yīng)的影響就會(huì)被減弱。本文在進(jìn)行數(shù)據(jù)標(biāo)準(zhǔn)化的時(shí)候，主要采用了min-max 標(biāo)準(zhǔn)化，通常也將之稱為規(guī)范化方法，下面對(duì)其進(jìn)行詳細(xì)介紹。

對(duì)序列x1,x2,x3,···,xn進(jìn)行變換：

3.3 數(shù)據(jù)初步展示

經(jīng)過數(shù)據(jù)預(yù)處理后，我們先用雷達(dá)圖簡(jiǎn)單的呈現(xiàn)全部數(shù)據(jù)的整體走向趨勢(shì)的概況，預(yù)先了解溫度的大體變化走向，我們將數(shù)據(jù)的csv 文件導(dǎo)入python 程序，經(jīng)數(shù)據(jù)統(tǒng)計(jì)，用雷達(dá)圖逐年代展示溫度變化趨向。

如圖2 所示越往里溫度值越低，越往外溫度值越高，每年溫度12 個(gè)月份的變化趨勢(shì)，6 月～9 月溫度全年最高，12 月～2 月溫度全年最低。1940 年～2010年每月溫度的變化趨勢(shì)基本一致。

圖2 溫度變化雷達(dá)圖

4 研究方法

4.1 線性趨勢(shì)估計(jì)

本文運(yùn)用的主要是根據(jù)魏鳳英老師的《現(xiàn)代氣候統(tǒng)計(jì)診斷與預(yù)測(cè)技術(shù)》一書中推薦的線性趨勢(shì)評(píng)估法，線性趨勢(shì)估計(jì)法是用一條合理的直線指示天氣變量與時(shí)間變量的聯(lián)系。［5］用xi表示樣本量為n的某一氣象變化值，用ti表示xi所對(duì)應(yīng)的時(shí)間，構(gòu)建xi與ti之間的一元線性回歸方程：

其中，a為回歸常數(shù)，b為回歸系數(shù)，a和b可用最小二乘法進(jìn)行估量［5］。

對(duì)于觀測(cè)數(shù)據(jù)xi及相應(yīng)的時(shí)間ti，回歸常數(shù)a和系數(shù)b的最小二乘法估值為

時(shí)間ti與變量xi之間的相關(guān)系數(shù)借助回歸系數(shù)b與相關(guān)系數(shù)的關(guān)聯(lián)性計(jì)算：

4.2 相關(guān)分析與回歸分析

4.3 移動(dòng)平均

移動(dòng)（滑動(dòng)）平均可用來模擬趨勢(shì)，在訓(xùn)練量為n時(shí)，給定時(shí)間序列x，此處為氣候變化序列x，取其平均值描繪分析波動(dòng)趨勢(shì)［5］，表達(dá)其滑動(dòng)平均序列是：

式（8）中，滑動(dòng)長(zhǎng)度為k，j=1,2,…,n-k+1。觀察七十年的時(shí)間序列以及分析周期，將滑動(dòng)窗口k設(shè)置為10，即觀測(cè)氣候序列的變化趨勢(shì)是通過10年一個(gè)滑動(dòng)分離出變化趨勢(shì)以探究規(guī)律。

4.4 ARIMA時(shí)間序列分析

ARIMA 模型（Autoregressive Integrated Moving Average model），也就是所謂的差分整合移動(dòng)平均自回歸模型（移動(dòng)也可稱作滑動(dòng)），預(yù)測(cè)時(shí)間序列的趨勢(shì)［10］。ARIMA（p，d，q）里，p 是自回歸項(xiàng)數(shù)，用來解釋時(shí)間序列模型中數(shù)據(jù)自身的滯后數(shù)，即AR項(xiàng)，d 是要達(dá)到平穩(wěn)序列需要做的差分次數(shù)，即Integrated 項(xiàng)，q 是滑動(dòng)平均項(xiàng)數(shù)，表示該模型在預(yù)測(cè)誤差范圍內(nèi)產(chǎn)生的滯后數(shù)，也就是MA 項(xiàng)［10］。“差分”是關(guān)鍵步驟。

此模型被選擇來預(yù)測(cè)溫度隨時(shí)間推移的變化序列。由ARMA（p，q）模型的延伸擴(kuò)充獲得ARIMA（p，d，q）模型。ARIMA（p，d，q）模型表達(dá)式是：

4.5 小波分析

1）小波概念

小波分析來自于Fourier，并且實(shí)現(xiàn)了Fourier所不能實(shí)現(xiàn)的一些艱難問題，適合于處理非穩(wěn)定信號(hào)［12］。小波分析可以用于時(shí)間-尺度分析和多分辨率分析。其波形同時(shí)具有衰減性和波動(dòng)性，定義公式如下：

2）小波變換

小波變換是將函數(shù)和數(shù)據(jù)等切分成不同的片段，然后分析相應(yīng)的尺度分量，小波變換有兩種變換方式，一種是通過連續(xù)形式實(shí)現(xiàn)變換，另一種則是采用離散手段達(dá)到變換；其中，離散小波又有著多種不同形式，比如：正交、雙正交小波。［13］由于海表溫度的數(shù)據(jù)有一定的連續(xù)震蕩的特征，所以采用一維連續(xù)小波變換分析其變化特點(diǎn)。一維連續(xù)小波變化的公式如下：

3）小波方差

小波方差即小波變換系數(shù)的平方值在b域上的積分［13］，其計(jì)算公式為

小波方差圖是小波方差值伴隨時(shí)間變化形成的曲線，該圖可以解釋在不同時(shí)間尺度上小波能量的分布狀況。

4.6 距平

距平是顯示氣候要素偏離正常值的最常見的量。數(shù)據(jù)集中的某個(gè)數(shù)xi與平均數(shù)之間的差值就是距平，即：

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

在氣象研究中，經(jīng)常選取距平序列來替換氣候指標(biāo)本身的觀測(cè)值。

5 結(jié)果分析

5.1 年代際溫度變化趨勢(shì)

將數(shù)據(jù)按照時(shí)間序列輸出距平值，反映了10年一個(gè)單位時(shí)間尺度上及整體的溫度變化特征［12］，從圖3 能夠發(fā)現(xiàn)，加拿大各地區(qū)從1940 年起，70 年來每10 年單位時(shí)間尺度之間差異性較大，整體上溫度有上升趨勢(shì)。

圖3 溫度變化趨勢(shì)

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

圖4 各月實(shí)際地表氣溫和地表氣溫均值

經(jīng)過觀察分析，我們可以看出月份之間存在相似性，于是我們做了相關(guān)性分析，采用斯皮爾曼相關(guān)性分析做了進(jìn)一步的相關(guān)性矩陣，具體如圖5 所示。

圖5 相關(guān)矩陣

相關(guān)性系數(shù)分布在0～1區(qū)間，越靠近1，說明相關(guān)程度越高，越靠近0，相關(guān)程度越低，從圖中可以看出，5 月和10 月相關(guān)關(guān)系最低，相關(guān)系數(shù)低至0.021。

5.2 基于ARIMA 分析溫度序列變化趨勢(shì)

通過ARIMA 模型對(duì)1940 年～2010 年的數(shù)據(jù)進(jìn)行STL 分解，得出1940 年～2010 年隨時(shí)間的溫度變化趨勢(shì)，如圖6所示。

圖6 1940年～2010年隨時(shí)間的溫度變化趨勢(shì)

根據(jù)分解結(jié)果Trend 圖的結(jié)論，溫度仍然是出現(xiàn)先略微下降再升高的趨勢(shì)。

5.3 溫度的空間分布變化

1）反距離權(quán)重插值法

反距離權(quán)重插值法是常用的、易于理解的、方便實(shí)用的、比較好的空間插值方式之一，是一種精準(zhǔn)插值，源于相近相似原理，具體表現(xiàn)為距離小，屬性越接近，類似程度越高，相反，間距遠(yuǎn)，特性差別較大［14］。基于插值點(diǎn)與樣本點(diǎn)之間的距離賦予權(quán)重求平均。間隔越小給予該點(diǎn)的權(quán)重值就會(huì)越大。經(jīng)過查閱文獻(xiàn)資料，我們選取ArcGIS 繪圖軟件，用反距離權(quán)重插值法在地圖上繪制溫度的空間分布［15］。

2）溫度空間分布分析

加拿大幅員遼闊，主要城市分布在5 個(gè)氣候區(qū)：太平洋地區(qū)（西海岸區(qū)）、平原地區(qū)、中部地區(qū)（大湖-圣勞倫斯蒂地區(qū)）、大西洋地區(qū)和北方地區(qū)，各個(gè)城市氣候有著不一樣的特點(diǎn)［16］。

我們可以從圖7 中分析溫度空間變化趨勢(shì)，在這段時(shí)間中，加拿大所有地區(qū)的平均氣溫都有著不同程度的上升。上升比較慢的是東西海岸，大約在2℃～3℃之間。結(jié)合相關(guān)地理知識(shí)易于理解，靠近海洋的氣候本來就是比較溫和。

圖7 溫度空間變化趨勢(shì)

內(nèi)陸地區(qū)與沿海地區(qū)差別較大，大部分平均溫度上升幅度都大于了4℃，北極三省的上漲幅度最大，基本都大于了6℃。顯而易見，特別突出明顯的是北極圈內(nèi)，漲幅甚至超過了10℃。通過新聞可知，溫室氣體增加，由此導(dǎo)致溫度上漲，成為了主要原因。

5.4 基于小波分析的海表溫度規(guī)律研究

首先根據(jù)收集的數(shù)據(jù)對(duì)海表溫度進(jìn)行初步展示，并且對(duì)1981 年～2010 年一共30 年的數(shù)據(jù)進(jìn)行小波分析，全球海表溫度變化趨勢(shì)、距平值及平均溫度分布如圖8、圖9、圖10所示。

圖8 全球海表溫度變化趨勢(shì)

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

圖10 全球海表平均溫度分布圖

全球海表溫度在近30 年間總體呈現(xiàn)上升的趨勢(shì)，如圖8 所示，1981 年海表溫度13.80 左右，在2006年出現(xiàn)最高溫度，溫度大約是13.99，增幅約為0.19。1981 年之后全球海表的距平值開始上升，1994 年到達(dá)第一個(gè)高峰，距平值約為0.35，在2010年出現(xiàn)最大距平。從圖8 中可以看出，海表溫度出現(xiàn)緩慢增長(zhǎng)趨勢(shì)，圖9 呈現(xiàn)海表溫度的距平在波動(dòng)曲折中上升，圖10 直觀顯示全球海表平均溫度分布［17］。

生成函數(shù)選擇墨西哥小帽小波，用一維相關(guān)連續(xù)小波對(duì)30 年這個(gè)時(shí)間段的數(shù)據(jù)集做變換操作，得到的灰度圖直觀顯示小波變換系數(shù)絕對(duì)值的大小程度［12］。橫軸為時(shí)間變換，縱軸為溫度變化，圖中越高的值，表征此處有越大的小波變換系數(shù)絕對(duì)值，意味著該時(shí)間段下海表溫度變化越明顯；小波變換系數(shù)為0 的坐標(biāo)點(diǎn)是該尺度下海表氣溫變化的拐點(diǎn)［18］。為了準(zhǔn)確地剖析30 年間海表溫度的變化和突變特點(diǎn)，繪制了絕對(duì)值灰度圖所對(duì)應(yīng)的等值線圖。通過等值線圖11 可以得出，全球海表溫度的變化周期約為10 年一個(gè)周期，在約5 年～6 年的尺度上有一個(gè)較強(qiáng)的擾動(dòng)。

圖11 等值線圖

6 結(jié)語

6.1 模型評(píng)價(jià)

首先借助距平值揭示地表溫度隨時(shí)間變化的動(dòng)態(tài)規(guī)律，得到12 個(gè)月的地表溫度變化。通過ARIMA對(duì)1940年～2010年的數(shù)據(jù)進(jìn)行STL分解，得到整體的氣溫趨勢(shì)是略下降再上升的趨勢(shì)［19］。對(duì)空間的溫度趨勢(shì)分析時(shí)，采用了反距離權(quán)重插值法和基于小波分析的方式對(duì)海表溫度規(guī)律探究，能清楚的用數(shù)據(jù)和圖表展示出1980 年～2010 年全球海洋表面溫度的變化趨勢(shì)可以得出全球海表溫度的變化趨勢(shì)為變暖的趨勢(shì)。

由于收集到的數(shù)據(jù)中指標(biāo)很多，初始樣本量更是高達(dá)十萬條，給數(shù)據(jù)量化處理工作帶來很大困難，因此在篩選過程中會(huì)遺漏部分可能存在重要作用的變量信息，并且我們?cè)谇蠼鈾?quán)重的時(shí)候，可能方法不盡合理，會(huì)對(duì)結(jié)果產(chǎn)生一定程度的影響。最后，由于指標(biāo)提取后，丟棄一些指標(biāo)會(huì)造成小信息的丟失，所以得到的結(jié)論可能會(huì)存在一些可接受的誤差。

收集的數(shù)據(jù)庫中數(shù)據(jù)缺失很多，而這些缺失的數(shù)據(jù)對(duì)加拿大的氣象站點(diǎn)的數(shù)據(jù)分析也是有影響的，雖然我們?cè)谶M(jìn)行建模和數(shù)據(jù)分析的時(shí)候，采用了K 最近距離鄰法，對(duì)數(shù)據(jù)進(jìn)行了填充，但是對(duì)最后的結(jié)果可能造成一定的誤差。

6.2 模型改進(jìn)

在分析空間和時(shí)間的溫度變化趨勢(shì)時(shí)，影響地表溫度變化的因素中沒有包括地形因素。由于氣候的形成往往受到地形的影響，地表溫度的空間分布因此發(fā)生了變化［20］。因此，若獲取到更詳細(xì)的地形數(shù)據(jù)集加入到考慮因素中分析時(shí)空分布的變化趨勢(shì)，結(jié)果的準(zhǔn)確度會(huì)顯著提高。