劉春霖， 夏建新

(中央民族大學(xué)生命與環(huán)境科學(xué)學(xué)院，北京 100081)

0 引言

隨著城市化進(jìn)程的快速發(fā)展，人類對(duì)土地資源大力開發(fā)利用，導(dǎo)致區(qū)域范圍內(nèi)地表覆被時(shí)空格局發(fā)生巨大變化，不但改變了城市地表景觀結(jié)構(gòu)，而且影響陸地生態(tài)系統(tǒng)的物質(zhì)循環(huán)、能量流動(dòng)及其穩(wěn)定性。因此，為了實(shí)現(xiàn)區(qū)域生態(tài)環(huán)境的科學(xué)保護(hù)，必須對(duì)城市土地變化規(guī)律有充分的認(rèn)識(shí)和研究。

近年來元胞自動(dòng)機(jī)模型(cellular automata, CA)在時(shí)空動(dòng)態(tài)模擬領(lǐng)域受到廣泛關(guān)注，CA提供了時(shí)空動(dòng)態(tài)模擬運(yùn)算框架，利用自下而上的模擬方式，由局部規(guī)則演化出全局的變化模式，實(shí)現(xiàn)了受宏觀地理?xiàng)l件和土地利用局部變化雙重作用下的土地利用變化模擬[1-7]。Von Neumann[8]首次將CA定義為一個(gè)離散的動(dòng)力學(xué)模型，由個(gè)體之間局部的行為演化出時(shí)間與空間上全局的變化模式。CA模型可與統(tǒng)計(jì)學(xué)方法和人工智能方法結(jié)合，實(shí)現(xiàn)更加精確、智能的復(fù)合模型。土地利用覆被變化過程是一個(gè)復(fù)雜的非線性變化過程，有學(xué)者指出需要建立轉(zhuǎn)換規(guī)則自動(dòng)獲取方法應(yīng)對(duì)復(fù)雜的土地利用變化過程。他們對(duì)利用統(tǒng)計(jì)學(xué)理論自動(dòng)提取元胞轉(zhuǎn)換規(guī)則的研究取得了突破性進(jìn)展，例如多層感知機(jī)(multi-layer perceptron, MLP)、支持向量機(jī) (support vector machine, SVM)和隨機(jī)森林(random forest, RF)模型[9-10]。

隨著深度學(xué)習(xí)的不斷發(fā)展，越來越多的深度學(xué)習(xí)模型被應(yīng)用于土地利用分類及變化模擬中，卷積神經(jīng)網(wǎng)絡(luò)是深度學(xué)習(xí)中常用的一種網(wǎng)絡(luò)模型，它通常使用卷積濾波器來提取局部區(qū)域的隱含特征。He等[11]考慮到空間效應(yīng)，利用卷積神經(jīng)網(wǎng)絡(luò)提取驅(qū)動(dòng)因素的空間特征來模擬土地利用變化。然而這些方法都是在一定假設(shè)理論基礎(chǔ)上進(jìn)行的，即某個(gè)單元的狀態(tài)只與上一個(gè)時(shí)間步驟的狀態(tài)有關(guān)，與其他時(shí)刻狀態(tài)無關(guān)[9,12-14]，但是土地利用變化是一個(gè)長(zhǎng)期過程，因此這些模型對(duì)于歷史土地利用變化規(guī)則并未得到充分挖掘。

在深度學(xué)習(xí)模型中遞歸神經(jīng)網(wǎng)絡(luò)用于建立時(shí)間序列數(shù)據(jù)之間的依賴性，它可以在時(shí)間步長(zhǎng)之間傳遞歷史信息，但普通的循環(huán)神經(jīng)網(wǎng)絡(luò)(recurrent neural network, RNN)在實(shí)際應(yīng)用中很難處理長(zhǎng)距離的依賴。長(zhǎng)短時(shí)記憶網(wǎng)絡(luò)(long short term memory network, LSTM)成功解決了原始RNN的問題，能夠?qū)哂虚L(zhǎng)時(shí)間依賴性的時(shí)間序列進(jìn)行建模，能夠充分考慮歷史時(shí)序土地分類的時(shí)間變化特征，被廣泛應(yīng)用于土地覆蓋分類中[15-16]。因此，本文以張家口市中心城區(qū)(包括橋東區(qū)、橋西區(qū)及萬全區(qū))為研究對(duì)象，將1995年、2000年、2005年、2010年、2015年歷史土地利用數(shù)據(jù)作為基礎(chǔ)數(shù)據(jù)，提出將LSTM和CA模型結(jié)合來模擬2020年土地利用的動(dòng)態(tài)變化，最后將模擬結(jié)果與實(shí)際2020年土地利用分布結(jié)果比較，證明本文提出的LSTM-CA模型能夠有效提高土地利用變化模擬準(zhǔn)確性。

1 研究區(qū)概況及數(shù)據(jù)處理

近年來隨著城市化的不斷發(fā)展，張家口市中心城區(qū)(橋東區(qū)、橋西區(qū)及萬全區(qū))面積迅速擴(kuò)張，該區(qū)域土地類型變化明顯，因此本文以張家口市中心城區(qū)為主要研究區(qū)，如圖1所示，該區(qū)域位于N40.645°～41.059°，E114.440°～114.996°之間，其中橋東區(qū)及橋西區(qū)以清水河為界，研究區(qū)域主要是山地與河谷地，其海拔范圍為800～1 500 m。

圖1 研究區(qū)概況圖(Landsat8 B7(R),B5(G),B2(B)合成遙感影像)Fig.1 Overview of the study area

本文以1995年、2000年、2005年、2010年、2015年、2020年Landsat5/8遙感數(shù)據(jù)為基礎(chǔ)，參考《城市用地分類與規(guī)劃建設(shè)用地標(biāo)準(zhǔn)》，通過卷積神經(jīng)網(wǎng)絡(luò)將研究區(qū)劃分為耕地、建設(shè)用地、水體、林地、草地5類土地利用類型，根據(jù)研究實(shí)測(cè)調(diào)研數(shù)據(jù)進(jìn)行精度檢驗(yàn)，總體分類精度達(dá)90%以上。

土地利用變化的概率往往取決于一系列的距離變量和單元的自然屬性等。例如，某一模擬單元越接近市中心及交通要道，其轉(zhuǎn)變?yōu)槌鞘杏玫氐母怕试礁?。由于研究區(qū)所在區(qū)域主要是山地與河谷地，土地利用分布受地形影響顯著，因此河流水系是影響該城市土地分布的主要因素之一。與此同時(shí)，隨著城市化的不斷發(fā)展，其土地利用變化與城市交通緊密相連，土地利用變化受人類活動(dòng)影響顯著，對(duì)于地勢(shì)相對(duì)平緩的區(qū)域，人類活動(dòng)相對(duì)頻繁，城市交通分布更加復(fù)雜，因此，本文在選擇驅(qū)動(dòng)因子時(shí)主要考慮城市水系、地形及交通等指標(biāo)，其中地形因素包括高程、坡度、坡向等信息，對(duì)于交通因素主要選取市中心點(diǎn)、高速公路、火車站、城市主干道、城鎮(zhèn)點(diǎn)、村莊點(diǎn)等要素，越是靠近城市中心，交通分布就越復(fù)雜，因此，利用ArcGIS軟件距離分析工具，計(jì)算研究區(qū)內(nèi)每一點(diǎn)到城市交通要素的最近歐式距離，以此作為土地利用變化另一空間指標(biāo)，并對(duì)這些空間變量進(jìn)行歸一化處理，結(jié)果如圖2所示。除此以外，鄰近現(xiàn)有土地利用類型的數(shù)量對(duì)于土地利用變化也有一定的影響，例如當(dāng)鄰近范圍內(nèi)存在同一種土地類型越多，該單元在下一時(shí)刻不變的概率就較高，在本文中將與當(dāng)前土地利用單元狀態(tài)相同的鄰域單元數(shù)量加入空間變量中。

(a) 高程 (b) 坡向 (c) 坡度 (d) 到城市中心距離 (e) 到高速公路入口距離

(f) 到火車站距離 (g) 到主干道距離 (h) 到城鎮(zhèn)距離 (i) 到村莊距離 (j) 到河流水系距離圖2 地形、距離變量空間分布Fig.2 Spatial distribution of terrain and distance variables

2 研究方法

將LSTM網(wǎng)絡(luò)、CA和地理信息系統(tǒng)(geographic information system，GIS)相結(jié)合進(jìn)行土地利用動(dòng)態(tài)模擬，流程如圖3所示，圖中x1，x2，…，xT表示時(shí)間序列的輸入樣本數(shù)據(jù)，h1，h2，…，hT表示時(shí)間序列的輸出特征。主要流程包括： ①基礎(chǔ)數(shù)據(jù)獲取，包括計(jì)算各歷史時(shí)期初始化的土地分類概率，統(tǒng)計(jì)相同土地利用狀態(tài)的鄰域單元數(shù)量，生成驅(qū)動(dòng)因子； ②根據(jù)基礎(chǔ)數(shù)據(jù)在研究區(qū)內(nèi)隨機(jī)生成訓(xùn)練樣本數(shù)據(jù)； ③基于長(zhǎng)時(shí)間序列土地分類數(shù)據(jù)，利用LSTM網(wǎng)絡(luò)開展土地利用發(fā)展概率計(jì)算； ④通過添加擾動(dòng)因子、限制因子和鄰域因子進(jìn)行土地利用全局轉(zhuǎn)化概率計(jì)算，開展土地利用動(dòng)態(tài)模擬研究。

圖3 LSTM-CA流程Fig.3 Process of LSTM-CA

2.1 基于LSTM網(wǎng)絡(luò)的土地利用發(fā)展概率計(jì)算

LSTM(圖4)作為一種特殊RNN，可以解決長(zhǎng)時(shí)間序列數(shù)據(jù)在訓(xùn)練中的梯度爆炸和消失問題。對(duì)于普通的RNN，權(quán)重在各個(gè)時(shí)間上是共享的，由于梯度被近時(shí)間梯度主導(dǎo)，導(dǎo)致模型難以學(xué)到遠(yuǎn)距離的依賴關(guān)系。對(duì)于LSTM網(wǎng)絡(luò)，由于總的遠(yuǎn)距離梯度等于各條路徑的遠(yuǎn)距離梯度之和，因此，LSTM能夠在更長(zhǎng)時(shí)間序列數(shù)據(jù)中具有更好的表現(xiàn)。RNN通常由一系列循環(huán)單元組成，相對(duì)于普通RNN，LSTM在循環(huán)單元的基礎(chǔ)上增加了門限限制，這樣可以使得信息進(jìn)行有效過濾，在LSTM中每一個(gè)單元都有輸入門、輸出門和遺忘門3部分組成，其結(jié)構(gòu)如圖4所示。圖中σ為sigmoid激活函數(shù)；tanh為tanh激活函數(shù)； 紅色線表示遺忘門; 綠色線表示輸入門; 藍(lán)色線表示輸出門;xt為當(dāng)前時(shí)刻輸入數(shù)據(jù);ht表示當(dāng)前時(shí)刻的特征輸出; ⊕和?分別表示加操作和乘操作。

遺忘門限控制當(dāng)前單元狀態(tài)中丟棄的信息，利用當(dāng)前時(shí)間的輸入和前一個(gè)時(shí)間的輸出來通過sigmoid函數(shù)來使得單元狀態(tài)乘以這個(gè)sigmoid函數(shù)的輸出。若sigmoid函數(shù)輸出0則該部分信息需要被丟棄或遺忘，反之該部分信息繼續(xù)在單元狀態(tài)中繼續(xù)傳播。輸入門限層是控制更新舊的單元狀態(tài)。之前的遺忘門限中sigmoid 層決定哪些信息需要更新，通過 tanh 激活函數(shù)計(jì)算用來更新的內(nèi)容，把這2部分聯(lián)合起來，對(duì)單元狀態(tài)進(jìn)行更新。輸出門限是對(duì)單元狀態(tài)的限制，從而決定輸出的信息。假設(shè)X=(x1,x2,…,xt,…,xT)表示時(shí)間序列的輸入，H=(h1,h2,…,ht,…,hT)表示經(jīng)過LSTM后對(duì)應(yīng)時(shí)間序列的輸出特征，則X與H之間的映射關(guān)系可以表示為：

(1)

(2)

(3)

(4)

(5)

(6)

式中：it為t時(shí)刻輸入門限；σ為sigmoid激活函數(shù)；W為全連接網(wǎng)絡(luò)的權(quán)重；b為全連接網(wǎng)絡(luò)的偏置；Ct為t時(shí)刻單元狀態(tài)；ft為t時(shí)刻遺忘門限；ot為t時(shí)刻輸出門限。

在本文中，LSTM采用了歷史土地利用分類及驅(qū)動(dòng)因子作為輸入向量，xt=[xt1,xt2,…,xt16]，其中xt1～xt5表示t時(shí)刻土地利用分類概率，xt6表示t時(shí)刻與當(dāng)前土地利用單元狀態(tài)相同的鄰域單元數(shù)量，xt7～xt16分別表示高程、坡向、坡度、到市中心距離、到城鎮(zhèn)距離、到村莊距離、到高速公路入口距離、到火車站距離及到主干道距離等驅(qū)動(dòng)因子。在模型訓(xùn)練階段，將1995—2010年4期土地利用分類及驅(qū)動(dòng)因子作為L(zhǎng)STM模型的輸入，通過實(shí)際2015年土地利用構(gòu)建交叉熵?fù)p失，采用Adam算法優(yōu)化網(wǎng)絡(luò)模型。在預(yù)測(cè)階段，采用2000—2015年4期土地利用及驅(qū)動(dòng)因子數(shù)據(jù)作為輸入，模擬2020年土地利用單元發(fā)展概率。

2.2 基于LSTM-CA模型的土地利用模擬

LSTM計(jì)算的發(fā)展概率只考慮各種空間變量對(duì)土地轉(zhuǎn)化的影響，在CA中鄰域元胞狀態(tài)對(duì)土地轉(zhuǎn)化同樣至關(guān)重要，例如鄰域元胞有較多轉(zhuǎn)變?yōu)榻ㄔO(shè)用地，則該元胞轉(zhuǎn)化為建設(shè)用地的概率也越大。因此，本文選擇3×3鄰域窗口計(jì)算土地類型所占的比例作為鄰域因子Ω，公式為：

(7)

式中：si為第i個(gè)元胞3×3鄰域中單元狀態(tài)的土地利用類型；con(si=c)表示第i個(gè)元胞鄰域狀態(tài)為c時(shí)則返回1，否則返回0。

土地利用轉(zhuǎn)換的過程中受到自然因素、人類活動(dòng)等影響，使土地利用變化過程更加復(fù)雜。為了使模型結(jié)果更加符合實(shí)際情況，反映土地變化的不確定性，因此CA中引入隨機(jī)擾動(dòng)因子R[17]，其表達(dá)公式為：

(8)

式中：γ為0～1范圍內(nèi)的均勻隨機(jī)變量；α為控制隨機(jī)擾動(dòng)大小的參數(shù)，在本文中α取值為5。

對(duì)于土地變化模擬指定相應(yīng)的規(guī)則至關(guān)重要，有必要對(duì)其地表變化進(jìn)行限制，轉(zhuǎn)換概率主要與研究區(qū)的城市發(fā)展水平相關(guān)，例如，城市建設(shè)用地轉(zhuǎn)化為草地的成本相對(duì)較高，而農(nóng)用地轉(zhuǎn)為城市建設(shè)用地的成本為相對(duì)較低。Liu等[18]提出每個(gè)土地利用對(duì)的轉(zhuǎn)換成本(cost)是根據(jù)當(dāng)?shù)貙＜医?jīng)驗(yàn)和城市規(guī)劃者確定的(表1)。轉(zhuǎn)換成本的值在 [0,1] 的范圍內(nèi)變化。較大的值表示較大的轉(zhuǎn)換難度，值為 1 表示幾乎不可能轉(zhuǎn)換。

表1 土地利用各類型轉(zhuǎn)換成本Tab.1 Conversion costs of various types of land use

因此，本文根據(jù)1-cost作為限制因子，根據(jù)歷史時(shí)期的土地利用狀態(tài)分別計(jì)算每個(gè)時(shí)期對(duì)應(yīng)的限制因子，最終將多個(gè)時(shí)期限制因子進(jìn)行平均計(jì)算。

考慮到近鄰范圍的影響、不確定因素影響及限制因素影響，通過隨機(jī)擾動(dòng)因子R、單元發(fā)展概率Pr、鄰域因子Ω3×3(8鄰域)及限制性因子L計(jì)算最終的土地利用概率Pt，其計(jì)算公式為：

(9)

3 結(jié)果與分析

本文根據(jù)提出的LSTM-CA模型，以1995—2015年土地利用及驅(qū)動(dòng)因子數(shù)據(jù)作為基礎(chǔ)，模擬2020年土地利用變化，并與MLP-CA模型結(jié)果進(jìn)行對(duì)比，模擬對(duì)比如表2所示。從表2中可以看出，整體模擬結(jié)果與2020年真實(shí)土地利用分類結(jié)果相似，且2020年土地利用中建設(shè)用地呈現(xiàn)出增加的趨勢(shì)，在實(shí)際2020年實(shí)際土地利用分布中也可以看出，建設(shè)用地面積增加，主要表現(xiàn)在圍繞城市中心向四周擴(kuò)張。為了進(jìn)一步驗(yàn)證方法的性能，通過局部細(xì)節(jié)進(jìn)行對(duì)比分析，從第3行中可以看出，LSTM-CA模擬的道路信息更加精細(xì)，相對(duì)于MLP-CA建設(shè)用地的分布更加緊湊，這主要是由于LSTM-CA充分考慮了不同時(shí)間各土地類型之間的轉(zhuǎn)換關(guān)系； 另外也有效利用了路網(wǎng)等交通因素，同時(shí)通過觀察可以發(fā)現(xiàn)，LSTM-CA模擬的林地噪聲更少，主要是由于林地的分布相對(duì)穩(wěn)定，通過多期土地利用分類數(shù)據(jù)能夠有效減少林地分布的不確定性。

表2 2020年土地利用模擬結(jié)果Tab.2 Land use simulation results in 2020

為了定量檢測(cè)方法的模擬結(jié)果，將模擬的各土地類型與2020年實(shí)際土地類型進(jìn)行疊加分析，統(tǒng)計(jì)得到各類像元一致性對(duì)比情況，如表3所示。從總體預(yù)測(cè)像元數(shù)量上來看，LSTM-CA模擬的結(jié)果與實(shí)際土地利用更為接近，建設(shè)用地、耕地、草地、林地及水體的準(zhǔn)確率分別94.03%，97.38%，66.61%，97.78%和73.83%，而MLP-CA模擬結(jié)果的準(zhǔn)確率分別為79.99%，96.27%，67.46%，95.42%和71.32%，除草地結(jié)果基本持平外，LSTM-CA其他類別模擬精度均高于MLP-CA結(jié)果，可見LSTM-CA方法兼顧了多時(shí)期的土地利用類型變化，對(duì)于模擬的各土地利用類型像元數(shù)量更符合實(shí)際情況。

表3 2020年土地利用實(shí)際柵格與模擬柵格對(duì)比Tab.3 Comparison between actual grid and simulated grid of land use in 2020 (個(gè))

為了更好地驗(yàn)證土地類型在空間上的合理性，本文將采用全局點(diǎn)對(duì)點(diǎn)的Kappa系數(shù)和變化差異對(duì)比的FoM(figure of merit)指標(biāo)作為評(píng)價(jià)標(biāo)準(zhǔn)。Kappa落在(0,1)間，通常情況下認(rèn)為Kappa達(dá)到0.8以上表現(xiàn)出較高的精度。FoM指標(biāo)的計(jì)算公式為：

(10)

式中:A為實(shí)際變化但模擬結(jié)果未發(fā)生變化的單元；B為實(shí)際變化同時(shí)模擬結(jié)果也發(fā)生變化的單元；C為實(shí)際變化同時(shí)模擬結(jié)果也變化，但是變化方向不一致的單元；D為實(shí)際未變化但模擬結(jié)果變化的單元。當(dāng)FoM指標(biāo)大于或接近于0.21時(shí)，說明模擬結(jié)果具有一定的可信度，認(rèn)為此模型具有較強(qiáng)的適用性。

通過對(duì)2種模擬方法進(jìn)行評(píng)價(jià)指標(biāo)計(jì)算，結(jié)果如表4所示。從表4中可以看出，2種土地利用模擬方法Kappa系數(shù)均達(dá)到了0.85以上，說明2種方法模擬的結(jié)果全局精度較高，而本文提出的LSTM-CA模型相比于MLP-CA方法Kappa精度提高了0.03，說明本方法在全局模擬結(jié)果上更具有優(yōu)勢(shì)。而FoM指標(biāo)均大于0.21，并且LSTM-CA相比于MLP-CA提升了0.17，除此以外， LSTM-CA模擬結(jié)果在A，B，C，D這4個(gè)指標(biāo)上均優(yōu)于MLP-CA，說明本文提出的LSTM-CA方法具有較高的模擬精度，說明LSTM-CA在充分挖掘歷史土地利用變化之間的關(guān)系，從而有效提升模擬精度。

表4 土地利用模擬結(jié)果精度評(píng)價(jià)Tab.4 Accuracy evaluation of land use simulation results

4 結(jié)論

本文在CA模型基礎(chǔ)上提出了基于LSTM-CA模型用于土地利用變化模擬，以張家口市為研究區(qū)，基于1995年、2000年、2005年、2010年、2015年時(shí)序土地利用數(shù)據(jù)，結(jié)合距離變量、鄰域變量及單元自然屬性預(yù)測(cè)了2020年的土地利用分布情況，與MLP-CA進(jìn)行精度對(duì)比分析，結(jié)果表明本文提出的LSTM-CA模型Kappa系數(shù)可達(dá)0.90，F(xiàn)oM指標(biāo)為0.39，達(dá)到了更高精度。說明該模型能夠較充分地挖掘時(shí)間序列土地利用之間的時(shí)空變換特征，通過歷史時(shí)序土地利用分類數(shù)據(jù)及驅(qū)動(dòng)因子構(gòu)建神經(jīng)網(wǎng)絡(luò)模型計(jì)算單元發(fā)展概率，模擬未來土地利用變化。

然而LSTM-CA也存在一定的局限性，即限制規(guī)則根據(jù)前人經(jīng)驗(yàn)獲取，因此，未來將繼續(xù)研究如何更加靈活地獲取限制因子，同時(shí)考慮更多的驅(qū)動(dòng)因素(如土地政策和社會(huì)經(jīng)濟(jì)等)是下一步研究的重點(diǎn)。