基于ArcGIS ModelBuilder的LUCC轉移矩陣自動計算模型研究

王育紅

(江蘇師范大學 地理測繪與城鄉規劃學院，江蘇 徐州 221116)

自1995年，“國際地圈生物圈計劃(IGBP)”和“全球環境變化人文計劃(HDP)”兩大國際科學計劃組織共同擬定發表《土地利用/覆被變化科學研究計劃》以來，土地利用覆被/變化(Land Use/Cover Change，LUCC)一直是世界各國學術界持續關注和廣泛研究的熱點。人們從不同目的出發構造多種LUCC分析與評價模型[1]，如空間統計模型、系統動力學模型、元胞自動機模型、智能體模型以及綜合模型等[2]。

轉移矩陣作為其中最為核心的基礎性模型，更是被眾多專家學者所使用。例如，史培軍等利用轉移矩陣分析了深圳市1980—1994年間的土地利用變化空間過程[3]；段增強等基于轉移矩陣分析挖掘了北京市海淀區1991—2001年間的土地利用變化信息[4]；謝金開等通過轉移矩陣研究了烏魯木齊河下游地區1989— 2007年間的土地利用/覆蓋時空動態變化規律[5]；廖諶婳等憑借轉移矩陣對中老緬泰交界地區1990—2010年間的土地利用變化特征進行了分析[6]；任東風等運用轉移矩陣對于彰武縣2008—2018年間的土地沙漠化問題進行了研究[7]。

從文獻資料上看，現有的大多數研究都是直接給出研究區的轉移矩陣結果，對其計算方法與過程往往缺乏深入詳細的闡述。針對這一情況，文中在界定總結轉移矩陣概念及其計算原理的基礎上，描述ArcGIS平臺下人工手動計算轉移矩陣的過程及工具，并采用ModelBuilder技術通過集成現有工具設計實現一個轉移矩陣自動批量計算工具模型。

1 轉移矩陣及其計算原理

LUCC轉移矩陣來源于系統分析中對系統狀態與狀態之間轉移的定量描述，可以從數量和結構上，綜合描述某一地區在期初和期末限定的研究時期內，各種土地利用/覆被類型之間的總體變化轉移情況。在實際應用中，轉移矩陣有兩種形式：即以變化轉移面積絕對數量為統計指標的轉移面積矩陣(通常稱為轉移矩陣)和以變化轉移面積相對數量為統計指標的轉移概率矩陣。

1.1 定義與性質

假設某一研究區的土地利用/覆蓋類型(簡稱地類)在研究期初(T0)為m種，其編碼值分別為1，2，3，…，m；在研究期末(T1)為n種，其編碼值分別為1，2，…，n。如果n=m，表示研究時期內地類種類沒有變化；如果n>m，表示研究期末地類種類增加；如果n

表1 地類轉移面積矩陣

表2 地類轉移概率矩陣

1.1.1 地類轉移面積矩陣及其基本性質

1)aii(i≤min(m,n))表示研究期間地類i沒有發生變化轉移的土地面積。

2)aij(i≠j,1≤i≤m,1≤j≤n)表示研究期初地類i在研究期末向地類j變化轉移的土地面積。

3)每一行的面積值之和(記作Ai)表示研究期初地類i對應的土地總面積，即：

(1)

4)每一列的面積值之和(記作Ai)表示研究期末地類j對應的土地總面積，即：

(2)

5)由于在期初和期末研究區的土地總面積保持不變，即式(3)恒成立。

(3)

1.1.2 地類轉移概率矩陣及其基本性質

2)每一行的轉移面積比例值之和為1，即式(4)恒成立。

(4)

為全面反映地類變化轉移數量和程度，建議將兩種轉移矩陣結合起來使用，避免誤判或低估某些地類的變化轉移情況。例如，地類變化轉移絕對面積值相對較小而相對變化轉移比例(概率)值卻很大的情況。但在土地利用/覆蓋空間數據集采用經緯度地理坐標系，不便于計算真實面積值的情況下，應使用轉移概率矩陣。通過轉移矩陣可以進一步挖掘獲取地類轉移變化的凈變化量、交換變化量、總變化量及其持續性、隨機性、系統性等特征[8-9]，轉移概率矩陣還可用于地類未來構成及數量的估算與預測[10]。有關這方面的具體模型與方法，這里不再詳述，請參考相應文獻資料。

1.2 基本計算原理

由于根據轉移面積矩陣可計算派生出轉移概率矩陣，因此本節及以后內容將圍繞轉移面積矩陣(簡稱轉移矩陣)討論相關問題。目前，轉移矩陣的計算主要利用遙感影像解譯分類后所獲得的土地利用/覆被柵格數據集，其核心思想是地類代碼融合。為便于敘述與理解，文中采用具體的示意性數據對其進行解釋。

假設研究區由8×8個像元組成，反映研究期初與期末土地利用/覆被狀況的柵格數據集分別為R0和R1，兩期都有5種地類，分別用代碼1～5表示，前后兩期的具體像元地類代碼值分別如圖1所示，那么根據R0和R1計算轉移矩陣將由如下兩大步驟組成：

1)采用地圖代數按照式(5)對R0地類代碼進行放大，并將放大處理結果與R1相加，形成反映地類變化轉移情況的柵格數據，即:

R01=R0×s+R1.

(5)

其中，s為初地類代碼放大倍數，其值取決于期末地類代碼最大值的位數(b)，兩者之間的關系為s=10b。圖1中R1數據集相應的b為1，故放大倍數s為10，因此可得如圖1所示的R01。

如果土地分類柵格數據集的像元地類代碼信息以文本名稱、標準編碼等形式表示，計算前應按相應的映射規則將其轉換成連續的整數型代碼。

圖1 基于柵格數據的轉移面積矩陣計算核心思想示意圖

2)對R01中不同編碼值的像元個數進行統計，根據空間分辨率(柵格像元所代表的實地面積)計算不同編碼對應的土地總面積，并對統計表格進行相應轉換處理以生成轉移矩陣。假設圖1中像元的空間分辨率為10 m×10 m，則對R01進行統計,所得表格和最終生成的轉移矩陣可分別為表3和表4的形式。

表3 柵格數據R01統計表

表4 轉移矩陣計算結果

2 基于ArcGIS的計算方法與模型

ArcGIS是全球著名GIS技術與軟件提供商—美國ESRI公司發布的全系列GIS軟件平臺的統稱，主要包括桌面GIS、嵌入GIS、移動GIS以及服務端GIS等多種子系統。據ARC咨詢公司2019年調查研究報告的數據顯示，ArcGIS產品占全球GIS市場的45%以上[11]。初步估計，目前ArcGIS產品也至少占國內GIS市場的25%以上。基于ArcGIS在國內的應用廣泛性，本節進一步討論基于ArcGIS桌面系統的轉移矩陣具體計算方法、步驟及工具。

2.1 手動計算方法

根據上述原理，使用ArcGIS桌面系統中的ArcMap模塊即可通過人工手動方式實現轉移矩陣的計算，具體計算過程主要包括如下步驟。

2.1.1 生成柵格數據集R01

打開位于“空間分析工具箱地圖代數工具集”下的“柵格計算器”工具，在該工具中的計算表達式窗口中輸入形如“R0”*s+“R1”的表達式，并指定輸出柵格數據集(即R01)存儲格式和位置，確認無誤后單擊“確定”生成R01。其中，s為代碼放大系數。

2.1.2 添加和計算字段值

打開R01的屬性表，在“ObjectID”、“Value”與“Count”現有字段的基礎上，添加“期初代碼”(整型)、“期末代碼”(整型)和“面積”(浮點型)，并在相應字段依次選擇右鍵菜單中的“字段計算器”，分別使用形如:Left([Value], Len([Value]) -b)、Right([Value],b)、[Count]*r的VB腳本語言計算表達式計算表中記錄3個字段的相應值。其中，b為期末地類代碼最大值位數；r為R01的空間分辨率(面積值)。

2.1.3 復制和透視屬性表

打開位于“數據管理工具箱表工具集”下的“復制行”工具，通過該工具將處理后的R01屬性表，復制備份到指定存儲位置。打開位于“數據管理工具箱表工具集”下的“數據集透視表”工具，在該工具中選擇剛復制的表作為輸入表，并分別將輸入字段、透視表字段和值字段三參數選擇設置為“期初代碼”“期末代碼”和“面積”字段，最后指定輸出表位置單擊“確定”即完成轉移矩陣的計算。

由于ArcGIS表格不支持數值型字段名及單元格的合并，所生成的轉移矩陣還難以做到與上述的定義形式完全一致。如果需要可通過修改字段別名，做到除“期初(T0)地類”“期初(T1)地類”字樣之外的非常一致。

通過上述描述可以看出，利用ArcMap手動計算轉移矩陣仍是一項比較繁瑣的工作，人機交互頻繁，所需工具及參數設置多，稍有不慎將出現錯誤，進而影響計算效率及結果的準確性。為解決這一問題，作者進一步采用Modelbuilder技術對上述計算步驟與工具進行組合與封裝，構建形成一個通用型自動批量計算工具模型。

2.2 自動計算模型

ModelBuilder是ArcGIS桌面系統中一個用來創建、編輯和管理模型的子程序，也可將其看作一種可視化編程語言。ModelBuilder可直接利用ArcGIS系統已有的各類工具，按一定的規則將其有機組合串聯起來，最終形成一個流程化的模型程序，從而實現高效的分析計算[12-14]。

根據上述手動計算流程，可通過如下主要步驟在ModelBuilder模型編輯窗口中創建實現轉換矩陣自動計算模型工具。

2.2.1 添加輸入變量

為提高所建模型工具的通用性，添加兩個柵格圖層類型和一個整數類型的變量，其名稱分別為“期初柵格數據集”“期末柵格數據集”和“期末地類代碼最大值位數”。

2.2.2 添加、設置、連接現有處理工具

通過拖拽方式依次從ArcMap系統工具箱里將“計算值”(位于數據管理工具箱常規工具集下，用于計算放大倍數)、“柵格計算器”“添加字段”“計算字段”(與前一個工具一樣都位于數據管理工具箱字段工具集下，都需要添加3個)，“獲取柵格屬性”(位于數據管理工具箱柵格工具集下，需要添加兩個分別用于獲取單個像元的寬度值和高度值)、“復制行”“數據透視表”等工具添加到模型編輯窗口中。

在添加工具之后，依次打開相應處理工具通過選擇已有變量及相關工具輸出結果或手動輸入的方式設置處理工具所需參數，并將工具連接起來形成工作流。為提高模型易于閱讀和理解的程度，還需進一步修改工具及其輸出結果的名稱。

2.2.3 設置輸出變量及模型屬性

為方便處理結果的管理、查看與對比，除“數據透視表”工具輸出的最終處理結果(即轉移矩陣)外，還將“柵格計算器”工具輸出的中間處理結果(即研究期內地類變化轉移柵格數據集)設置為輸出變量。

在驗證模型能無誤運行后，在模型屬性對話框的“參數”選項卡中調整變量參數顯示順序，將“期末地類代碼最大值位數”的過濾器參數設置為“值列表”，并添加1、2、3等3個可選值。另外，在“常規”選項卡中，為名稱、標簽、描述等參數設置準確的信息，以方便工具的查找和使用。

最后所得到的模型工具內部工作流程及外部運行界面效果分別如圖2與圖3所示。

圖2 模型內部工作流程

3 初步實驗與效率對比

為了對比分析上述計算方法與工具的效率性能，筆者以前期相關研究積累的南京、鄭州、徐州3地在不同時期的6幅土地利用分類柵格數據集為例進行初步實驗。其中，南京市包括城市用地、非城市用地以及未定義用地3種地類，鄭州、徐州兩地包括水體、林地、農田、菜地、旱地、建設用地6種地類。筆者邀請3位本科生使用這些數據在各自的筆記本電腦上，分別采用手動和自動方式對相應時期內的轉移矩陣進行10次計算。

圖3 模型外部運行界面效果圖

表5給出實驗數據的基本特征及每人的平均計算用時。從表中可以看出，通過模型自動計算用時平均可比人工手動計算縮短10倍以上。另外，該工具只需要選擇輸入5個參數，避免常規手動計算方法用戶全程參與、工具及參數設置多、人機交互頻繁、易出錯的不足。

表5 實驗數據集基本特征及計算用時統計表

4 結束語

文中在界定總結轉移矩陣概念、性質及其計算原理的基礎上，詳細闡述ArcGIS桌面系統下的人工手動計算方法與自動計算模型創建過程，為土地資源管理、規劃、科研以及相關開發人員提供有益的理論參考與技術支持。隨著全國土地利用調查及地理國情監測等國家重大戰略工程的常態化實施，我國土地利用/覆被矢量正以前所未有的速度快速增長，亟需開展深層次的挖掘分析。針對這一需要，研究創建基于矢量數據的轉移矩陣自動計算方法和模型工具將是下一步研究工作的重點。