劉春霖, 夏建新
(中央民族大學生命與環境科學學院,北京 100081)
隨著城市化進程的快速發展,人類對土地資源大力開發利用,導致區域范圍內地表覆被時空格局發生巨大變化,不但改變了城市地表景觀結構,而且影響陸地生態系統的物質循環、能量流動及其穩定性。因此,為了實現區域生態環境的科學保護,必須對城市土地變化規律有充分的認識和研究。
近年來元胞自動機模型(cellular automata, CA)在時空動態模擬領域受到廣泛關注,CA提供了時空動態模擬運算框架,利用自下而上的模擬方式,由局部規則演化出全局的變化模式,實現了受宏觀地理條件和土地利用局部變化雙重作用下的土地利用變化模擬[1-7]。Von Neumann[8]首次將CA定義為一個離散的動力學模型,由個體之間局部的行為演化出時間與空間上全局的變化模式。CA模型可與統計學方法和人工智能方法結合,實現更加精確、智能的復合模型。土地利用覆被變化過程是一個復雜的非線性變化過程,有學者指出需要建立轉換規則自動獲取方法應對復雜的土地利用變化過程。他們對利用統計學理論自動提取元胞轉換規則的研究取得了突破性進展,例如多層感知機(multi-layer perceptron, MLP)、支持向量機 (support vector machine, SVM)和隨機森林(random forest, RF)模型[9-10]。
隨著深度學習的不斷發展,越來越多的深度學習模型被應用于土地利用分類及變化模擬中,卷積神經網絡是深度學習中常用的一種網絡模型,它通常使用卷積濾波器來提取局部區域的隱含特征。He等[11]考慮到空間效應,利用卷積神經網絡提取驅動因素的空間特征來模擬土地利用變化。然而這些方法都是在一定假設理論基礎上進行的,即某個單元的狀態只與上一個時間步驟的狀態有關,與其他時刻狀態無關[9,12-14],但是土地利用變化是一個長期過程,因此這些模型對于歷史土地利用變化規則并未得到充分挖掘。
在深度學習模型中遞歸神經網絡用于建立時間序列數據之間的依賴性,它可以在時間步長之間傳遞歷史信息,但普通的循環神經網絡(recurrent neural network, RNN)在實際應用中很難處理長距離的依賴。長短時記憶網絡(long short term memory network, LSTM)成功解決了原始RNN的問題,能夠對具有長時間依賴性的時間序列進行建模,能夠充分考慮歷史時序土地分類的時間變化特征,被廣泛應用于土地覆蓋分類中[15-16]。因此,本文以張家口市中心城區(包括橋東區、橋西區及萬全區)為研究對象,將1995年、2000年、2005年、2010年、2015年歷史土地利用數據作為基礎數據,提出將LSTM和CA模型結合來模擬2020年土地利用的動態變化,最后將模擬結果與實際2020年土地利用分布結果比較,證明本文提出的LSTM-CA模型能夠有效提高土地利用變化模擬準確性。
近年來隨著城市化的不斷發展,張家口市中心城區(橋東區、橋西區及萬全區)面積迅速擴張,該區域土地類型變化明顯,因此本文以張家口市中心城區為主要研究區,如圖1所示,該區域位于N40.645°~41.059°,E114.440°~114.996°之間,其中橋東區及橋西區以清水河為界,研究區域主要是山地與河谷地,其海拔范圍為800~1 500 m。

圖1 研究區概況圖(Landsat8 B7(R),B5(G),B2(B)合成遙感影像)Fig.1 Overview of the study area
本文以1995年、2000年、2005年、2010年、2015年、2020年Landsat5/8遙感數據為基礎,參考《城市用地分類與規劃建設用地標準》,通過卷積神經網絡將研究區劃分為耕地、建設用地、水體、林地、草地5類土地利用類型,根據研究實測調研數據進行精度檢驗,總體分類精度達90%以上。
土地利用變化的概率往往取決于一系列的距離變量和單元的自然屬性等。例如,某一模擬單元越接近市中心及交通要道,其轉變為城市用地的概率越高。由于研究區所在區域主要是山地與河谷地,土地利用分布受地形影響顯著,因此河流水系是影響該城市土地分布的主要因素之一。與此同時,隨著城市化的不斷發展,其土地利用變化與城市交通緊密相連,土地利用變化受人類活動影響顯著,對于地勢相對平緩的區域,人類活動相對頻繁,城市交通分布更加復雜,因此,本文在選擇驅動因子時主要考慮城市水系、地形及交通等指標,其中地形因素包括高程、坡度、坡向等信息,對于交通因素主要選取市中心點、高速公路、火車站、城市主干道、城鎮點、村莊點等要素,越是靠近城市中心,交通分布就越復雜,因此,利用ArcGIS軟件距離分析工具,計算研究區內每一點到城市交通要素的最近歐式距離,以此作為土地利用變化另一空間指標,并對這些空間變量進行歸一化處理,結果如圖2所示。除此以外,鄰近現有土地利用類型的數量對于土地利用變化也有一定的影響,例如當鄰近范圍內存在同一種土地類型越多,該單元在下一時刻不變的概率就較高,在本文中將與當前土地利用單元狀態相同的鄰域單元數量加入空間變量中。

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

(f) 到火車站距離 (g) 到主干道距離 (h) 到城鎮距離 (i) 到村莊距離 (j) 到河流水系距離圖2 地形、距離變量空間分布Fig.2 Spatial distribution of terrain and distance variables
將LSTM網絡、CA和地理信息系統(geographic information system,GIS)相結合進行土地利用動態模擬,流程如圖3所示,圖中x1,x2,…,xT表示時間序列的輸入樣本數據,h1,h2,…,hT表示時間序列的輸出特征。主要流程包括: ①基礎數據獲取,包括計算各歷史時期初始化的土地分類概率,統計相同土地利用狀態的鄰域單元數量,生成驅動因子; ②根據基礎數據在研究區內隨機生成訓練樣本數據; ③基于長時間序列土地分類數據,利用LSTM網絡開展土地利用發展概率計算; ④通過添加擾動因子、限制因子和鄰域因子進行土地利用全局轉化概率計算,開展土地利用動態模擬研究。

圖3 LSTM-CA流程Fig.3 Process of LSTM-CA
LSTM(圖4)作為一種特殊RNN,可以解決長時間序列數據在訓練中的梯度爆炸和消失問題。對于普通的RNN,權重在各個時間上是共享的,由于梯度被近時間梯度主導,導致模型難以學到遠距離的依賴關系。對于LSTM網絡,由于總的遠距離梯度等于各條路徑的遠距離梯度之和,因此,LSTM能夠在更長時間序列數據中具有更好的表現。RNN通常由一系列循環單元組成,相對于普通RNN,LSTM在循環單元的基礎上增加了門限限制,這樣可以使得信息進行有效過濾,在LSTM中每一個單元都有輸入門、輸出門和遺忘門3部分組成,其結構如圖4所示。圖中σ為sigmoid激活函數;tanh為tanh激活函數; 紅色線表示遺忘門; 綠色線表示輸入門; 藍色線表示輸出門;xt為當前時刻輸入數據;ht表示當前時刻的特征輸出; ⊕和?分別表示加操作和乘操作。
遺忘門限控制當前單元狀態中丟棄的信息,利用當前時間的輸入和前一個時間的輸出來通過sigmoid函數來使得單元狀態乘以這個sigmoid函數的輸出。若sigmoid函數輸出0則該部分信息需要被丟棄或遺忘,反之該部分信息繼續在單元狀態中繼續傳播。輸入門限層是控制更新舊的單元狀態。之前的遺忘門限中sigmoid 層決定哪些信息需要更新,通過 tanh 激活函數計算用來更新的內容,把這2部分聯合起來,對單元狀態進行更新。輸出門限是對單元狀態的限制,從而決定輸出的信息。假設X=(x1,x2,…,xt,…,xT)表示時間序列的輸入,H=(h1,h2,…,ht,…,hT)表示經過LSTM后對應時間序列的輸出特征,則X與H之間的映射關系可以表示為:
(1)
(2)
(3)
(4)
(5)
(6)
式中:it為t時刻輸入門限;σ為sigmoid激活函數;W為全連接網絡的權重;b為全連接網絡的偏置;Ct為t時刻單元狀態;ft為t時刻遺忘門限;ot為t時刻輸出門限。
在本文中,LSTM采用了歷史土地利用分類及驅動因子作為輸入向量,xt=[xt1,xt2,…,xt16],其中xt1~xt5表示t時刻土地利用分類概率,xt6表示t時刻與當前土地利用單元狀態相同的鄰域單元數量,xt7~xt16分別表示高程、坡向、坡度、到市中心距離、到城鎮距離、到村莊距離、到高速公路入口距離、到火車站距離及到主干道距離等驅動因子。在模型訓練階段,將1995—2010年4期土地利用分類及驅動因子作為LSTM模型的輸入,通過實際2015年土地利用構建交叉熵損失,采用Adam算法優化網絡模型。在預測階段,采用2000—2015年4期土地利用及驅動因子數據作為輸入,模擬2020年土地利用單元發展概率。
LSTM計算的發展概率只考慮各種空間變量對土地轉化的影響,在CA中鄰域元胞狀態對土地轉化同樣至關重要,例如鄰域元胞有較多轉變為建設用地,則該元胞轉化為建設用地的概率也越大。因此,本文選擇3×3鄰域窗口計算土地類型所占的比例作為鄰域因子Ω,公式為:
(7)
式中:si為第i個元胞3×3鄰域中單元狀態的土地利用類型;con(si=c)表示第i個元胞鄰域狀態為c時則返回1,否則返回0。
土地利用轉換的過程中受到自然因素、人類活動等影響,使土地利用變化過程更加復雜。為了使模型結果更加符合實際情況,反映土地變化的不確定性,因此CA中引入隨機擾動因子R[17],其表達公式為:
(8)
式中:γ為0~1范圍內的均勻隨機變量;α為控制隨機擾動大小的參數,在本文中α取值為5。
對于土地變化模擬指定相應的規則至關重要,有必要對其地表變化進行限制,轉換概率主要與研究區的城市發展水平相關,例如,城市建設用地轉化為草地的成本相對較高,而農用地轉為城市建設用地的成本為相對較低。Liu等[18]提出每個土地利用對的轉換成本(cost)是根據當地專家經驗和城市規劃者確定的(表1)。轉換成本的值在 [0,1] 的范圍內變化。較大的值表示較大的轉換難度,值為 1 表示幾乎不可能轉換。

表1 土地利用各類型轉換成本Tab.1 Conversion costs of various types of land use
因此,本文根據1-cost作為限制因子,根據歷史時期的土地利用狀態分別計算每個時期對應的限制因子,最終將多個時期限制因子進行平均計算。
考慮到近鄰范圍的影響、不確定因素影響及限制因素影響,通過隨機擾動因子R、單元發展概率Pr、鄰域因子Ω3×3(8鄰域)及限制性因子L計算最終的土地利用概率Pt,其計算公式為:
(9)
本文根據提出的LSTM-CA模型,以1995—2015年土地利用及驅動因子數據作為基礎,模擬2020年土地利用變化,并與MLP-CA模型結果進行對比,模擬對比如表2所示。從表2中可以看出,整體模擬結果與2020年真實土地利用分類結果相似,且2020年土地利用中建設用地呈現出增加的趨勢,在實際2020年實際土地利用分布中也可以看出,建設用地面積增加,主要表現在圍繞城市中心向四周擴張。為了進一步驗證方法的性能,通過局部細節進行對比分析,從第3行中可以看出,LSTM-CA模擬的道路信息更加精細,相對于MLP-CA建設用地的分布更加緊湊,這主要是由于LSTM-CA充分考慮了不同時間各土地類型之間的轉換關系; 另外也有效利用了路網等交通因素,同時通過觀察可以發現,LSTM-CA模擬的林地噪聲更少,主要是由于林地的分布相對穩定,通過多期土地利用分類數據能夠有效減少林地分布的不確定性。

表2 2020年土地利用模擬結果Tab.2 Land use simulation results in 2020
為了定量檢測方法的模擬結果,將模擬的各土地類型與2020年實際土地類型進行疊加分析,統計得到各類像元一致性對比情況,如表3所示。從總體預測像元數量上來看,LSTM-CA模擬的結果與實際土地利用更為接近,建設用地、耕地、草地、林地及水體的準確率分別94.03%,97.38%,66.61%,97.78%和73.83%,而MLP-CA模擬結果的準確率分別為79.99%,96.27%,67.46%,95.42%和71.32%,除草地結果基本持平外,LSTM-CA其他類別模擬精度均高于MLP-CA結果,可見LSTM-CA方法兼顧了多時期的土地利用類型變化,對于模擬的各土地利用類型像元數量更符合實際情況。

表3 2020年土地利用實際柵格與模擬柵格對比Tab.3 Comparison between actual grid and simulated grid of land use in 2020 (個)
為了更好地驗證土地類型在空間上的合理性,本文將采用全局點對點的Kappa系數和變化差異對比的FoM(figure of merit)指標作為評價標準。Kappa落在(0,1)間,通常情況下認為Kappa達到0.8以上表現出較高的精度。FoM指標的計算公式為:
(10)
式中:A為實際變化但模擬結果未發生變化的單元;B為實際變化同時模擬結果也發生變化的單元;C為實際變化同時模擬結果也變化,但是變化方向不一致的單元;D為實際未變化但模擬結果變化的單元。當FoM指標大于或接近于0.21時,說明模擬結果具有一定的可信度,認為此模型具有較強的適用性。
通過對2種模擬方法進行評價指標計算,結果如表4所示。從表4中可以看出,2種土地利用模擬方法Kappa系數均達到了0.85以上,說明2種方法模擬的結果全局精度較高,而本文提出的LSTM-CA模型相比于MLP-CA方法Kappa精度提高了0.03,說明本方法在全局模擬結果上更具有優勢。而FoM指標均大于0.21,并且LSTM-CA相比于MLP-CA提升了0.17,除此以外, LSTM-CA模擬結果在A,B,C,D這4個指標上均優于MLP-CA,說明本文提出的LSTM-CA方法具有較高的模擬精度,說明LSTM-CA在充分挖掘歷史土地利用變化之間的關系,從而有效提升模擬精度。

表4 土地利用模擬結果精度評價Tab.4 Accuracy evaluation of land use simulation results
本文在CA模型基礎上提出了基于LSTM-CA模型用于土地利用變化模擬,以張家口市為研究區,基于1995年、2000年、2005年、2010年、2015年時序土地利用數據,結合距離變量、鄰域變量及單元自然屬性預測了2020年的土地利用分布情況,與MLP-CA進行精度對比分析,結果表明本文提出的LSTM-CA模型Kappa系數可達0.90,FoM指標為0.39,達到了更高精度。說明該模型能夠較充分地挖掘時間序列土地利用之間的時空變換特征,通過歷史時序土地利用分類數據及驅動因子構建神經網絡模型計算單元發展概率,模擬未來土地利用變化。
然而LSTM-CA也存在一定的局限性,即限制規則根據前人經驗獲取,因此,未來將繼續研究如何更加靈活地獲取限制因子,同時考慮更多的驅動因素(如土地政策和社會經濟等)是下一步研究的重點。