陳瑞卿，陶建斌，徐 猛

（1.地理過程分析與模擬湖北省重點實驗室/華中師范大學城市與環(huán)境科學學院，武漢 430079；2.中國農(nóng)業(yè)科學院農(nóng)業(yè)資源與農(nóng)業(yè)區(qū)劃研究所/農(nóng)業(yè)部農(nóng)業(yè)信息技術(shù)重點實驗室，北京 100081）

0 引言

耕地是進行農(nóng)業(yè)生產(chǎn)的基本資源和條件，是人類物質(zhì)文明和精神文明的來源。近年來，由于我國城市化水平和工業(yè)化水平不斷提高，人口總數(shù)持續(xù)增長，導致耕地被不斷蠶食，耕地資源總量和質(zhì)量都面臨著嚴峻的挑戰(zhàn)［1］。同時，整體上耕地資源缺乏科學管理和合理配置，導致耕地資源利用效率長期處于較低水平。因此，提高耕地的集約利用水平對于農(nóng)業(yè)的可持續(xù)發(fā)展和國家糧食安全保障工作十分重要［2］。

在耕地集約利用研究中，國內(nèi)學者通常使用復種指數(shù)來表示耕地資源的利用程度。復種指數(shù)是指耕地在年內(nèi)種植農(nóng)作物的次數(shù)，復種指數(shù)越大，表明人類農(nóng)業(yè)活動對耕地資源的利用程度越高［3-6］，在單產(chǎn)不明顯變化時，對應耕地的糧食產(chǎn)量也越高［7］。同時，國外學者通常使用種植頻率(Cropping Frequency)和種植強度(Cropping Intensity)作為耕地資源集約利用的評價指標，其中種植頻率表示一年內(nèi)耕地種植農(nóng)作物的次數(shù)，種植強度表示年內(nèi)耕地種植的農(nóng)作物收獲的次數(shù)［8-12］。這兩個指標和國內(nèi)使用較多的復種指數(shù)的含義很相近，都可以作為耕地資源利用程度的評價指標，揭示耕地資源的利用狀況并為區(qū)域糧食安全提供重要保障［13］。遙感具有大面積同步觀測的能力，更新速度快，時效性強，在區(qū)域尺度耕地資源監(jiān)測中得到廣泛應用。江漢平原所處的長江中下游地區(qū)是我國經(jīng)濟和人口快速發(fā)展的區(qū)域，也是耕地資源被各種建設(shè)項目大量侵占的區(qū)域，更是我國耕地資源長期保護的重點區(qū)域［14］。該文基于時間序列MODIS數(shù)據(jù)中提取的種植頻率數(shù)據(jù)，構(gòu)建了3種指標（雙季種植強度、平均種植強度、總體種植強度）對耕地種植強度進行表達。江漢平原耕地資源的集約利用研究，對于及時準確地掌握農(nóng)作物分布格局以及該區(qū)域的農(nóng)業(yè)發(fā)展和糧食安全具有重要的現(xiàn)實意義［15］。

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

1.1 研究區(qū)概況

江漢平原位于湖北省中南部（圖1），屬于長江中下游平原的一部分，是由長江及其最大干流漢江沖積形成。江漢平原海拔高度小于500 m、地形起伏平緩，水熱條件較好，耕地資源集中［16］，是我國重要的糧食基地之一。主要的種植作物包括水稻、小麥、油菜、棉花等，農(nóng)作物種植頻率主要分為一年一季種植和一年兩季種植。

圖1 研究區(qū)行政區(qū)劃Fig.1 Administrative map of research area

1.2 數(shù)據(jù)源

本文采用的數(shù)據(jù)包括EOS/TERRA衛(wèi)星的MOD13Q1植被指數(shù)產(chǎn)品數(shù)據(jù)和Landsat8數(shù)據(jù)。MOD13Q1植被指數(shù)產(chǎn)品數(shù)據(jù)來源于美國國家航空航天局（NASA）地球資源觀測系統(tǒng)數(shù)據(jù)中心網(wǎng)站（https：//ladsweb.modaps.eosdis.nasa.gov/）。MOD13Q1數(shù)據(jù)為陸地專題植被指數(shù)產(chǎn)品數(shù)據(jù)，其空間分辨率為250 m，時間分辨率為16 d，包括12種數(shù)據(jù)集圖層，主要有歸一化植被指數(shù)、增強型植被指數(shù)、質(zhì)量圖層等。本文選取的MOD13Q1產(chǎn)品歸一化植被指數(shù)（Normalized Difference Vegetation Index，NDVI）數(shù)據(jù)的時間為 2000—2017年，18年間共計411期，Landsat8數(shù)據(jù)時間為2015年3月6日，條帶號和行編號分別為124、39，用于制作驗證區(qū)土地覆蓋分類圖。

1.3 數(shù)據(jù)預處理

本文利用Modis Reprojection Tool批處理工具對MOD13Q1數(shù)據(jù)進行批量拼接、重采樣、重投影等處理，重采樣方法為Nearest Neighbor，重投影方式為UTM（Universal Transverse Mercator Projection）投影。由于 Landsat土地覆蓋數(shù)據(jù)的空間分辨率為 30 m，遠高于MODIS數(shù)據(jù)的250 m分辨率，因此需要進行MODIS數(shù)據(jù)重采樣，使之與Landsat數(shù)據(jù)的空間分辨率相同，以便進行兩種數(shù)據(jù)的空間對比和精度評價。

在遙感數(shù)據(jù)獲取過程中，傳感器易受到云、水汽、氣溶膠、太陽高度角等因素的影響，造成了植被指數(shù)時間序列曲線不穩(wěn)定，例如時間序列呈現(xiàn)鋸齒狀波動狀態(tài)、存在異常值等狀況，這一問題給物候信息的提取和趨勢分析造成了一定的困難。因此需要對時間序列數(shù)據(jù)進行平滑處理，以消除噪聲的影響。本文采用Savizky-Golay濾波法對時序數(shù)據(jù)進行擬合重建。

2 研究方法

時間序列NDVI曲線圖可以用來描述植被生長期內(nèi)的生長與衰落的周期性變化過程。對農(nóng)作物而言，單季農(nóng)作物的NDVI曲線在年內(nèi)顯示為一個波峰，雙季農(nóng)作物的NDVI曲線在年內(nèi)形成兩個波峰。因此，通過S-G濾波法對時間序列NDVI曲線進行擬合處理，得到重建的農(nóng)作物生長曲線，以此作為區(qū)分單/雙季農(nóng)作物的依據(jù)。

江漢平原種植強度時空格局分析的過程分為4個步驟：首先，從時間序列NDVI中提取了生長季起始期、生長季結(jié)束期、生長季長度3個物候參數(shù)；然后運用非監(jiān)督分類法對提取的物候參數(shù)進行聚類處理，得到種植頻率數(shù)據(jù)；接下來基于種植頻率數(shù)據(jù)定義了3種指標：雙季種植強度、平均種植強度、總體種植強度，來表達耕地的種植強度；最后，利用自組織神經(jīng)網(wǎng)絡(luò)（Self-Organizing Map，SOM）進行種植強度分析，并得到種植強度分區(qū)圖。

2.1 物候參數(shù)提取

植被物候是指植被在完整生長周期中的生長特點，主要包括植被的發(fā)芽生根、出苗、開花、結(jié)果、落葉等過程［17-18］，反映的是植被對各種自然環(huán)境、社會人文等因子動態(tài)變化的響應過程［19］。隨著遙感技術(shù)的發(fā)展，利用遙感數(shù)據(jù)反演得到的植被指數(shù)數(shù)據(jù)可以模擬植被隨時間變化的生長全過程，因此通過提取植被生長模擬過程中的關(guān)鍵信息，如生長季起始期、生長季結(jié)束期、生長季長度等，即可進行植被物候信息提取。

TIMESAT（Time-series Satellite Data Analysis Tool，TIMESAT）［20］是一款強大的時間序列遙感影像數(shù)據(jù)處理分析軟件。TIMESAT集成了Savitzky-Golay濾波法、雙線性邏輯函數(shù)擬合法、非對稱高斯擬合法等3種不同的數(shù)據(jù)擬合算法，能夠?qū)Σ煌瑫r間分辨率的時間序列數(shù)據(jù)進行擬合處理，并提取豐富的物候信息。

利用TIMESAT軟件對時間序列NDVI曲線進行擬合重建，并添加了質(zhì)量控制數(shù)據(jù)集來提高擬合效果。同時，通過重建的時間序列NDVI曲線提取得到植被物候信息，主要包括生長季起始期、生長季結(jié)束期、生長季長度、基準值、最大值、振幅等物候參數(shù)。本文選取比較有代表性的生長季起始期、生長季結(jié)束期和生長季長度等3種物候參數(shù)（如圖2所示）。圖3為2017年的3種物候參數(shù)分布數(shù)據(jù)?？梢钥吹剑p季農(nóng)作物的第一個生長季的起始期大約在一年中的第33d至第49d之間，結(jié)束期大約在第113d至第129d之間，第一個生長季長度大約為80—96d；單季作物的生長季起始期大約在第145d至第161d之間，結(jié)束期大約在第257d至273d之間，生長季長度大約為96d至112d。

圖2 基于時間序列NDVI數(shù)據(jù)的物候參數(shù)提取a.生長季起始期；b.生長季結(jié)束期；c.生長季長度Fig.2 Phenological parameter extraction based on time series NDVI data a.start of growing season；b.end of growing season；c.length of growing season

圖3 2017年江漢平原植被物候信息分布Fig.3 Distribution of vegetation phenological information in Jianghan Plain in 2017

2.2 種植頻率提取

在沒有樣本先驗知識的情況下，利用3種物候參數(shù)采用ISODATA算法進行聚類，目標分類數(shù)為5類：單季農(nóng)作物、雙季農(nóng)作物、城鎮(zhèn)建設(shè)用地、水體、天然植被。將天然植被、水體、城鎮(zhèn)建設(shè)用地等合并為其他一類。最終得到3個類別：單季作物、雙季作物、其他。2000—2017年江漢平原部分年份種植頻率分布見圖4。

2.3 種植強度提取

基于上述種植頻率數(shù)據(jù)，本文定義了3種指標：雙季種植強度、平均種植強度、總體種植強度，來表征耕地的種植強度。

2.3.1 雙季種植強度

將雙季種植強度定義為研究期內(nèi)某像元為雙季農(nóng)作物種植的累積年數(shù)。綜合歷年江漢平原種植頻率分布圖進行重分類，得到雙季農(nóng)作物種植強度分布圖。其中，像元值為“0”表示該地物類型為非農(nóng)用地，像元值域為“1～7”表示雙季種植強度為“低”，值域“8～14”表示雙季種植強度為“中”，值域“15～18”表示雙季種植強度為“高”。

2.3.2 平均種植強度

將歷年江漢平原農(nóng)作物種植頻率分布圖中的單季農(nóng)作物區(qū)域像元值賦值為“1”，雙季農(nóng)作物區(qū)域像元值賦值為“2”，其他區(qū)域像元值賦為“0”，然后進行多年平均種植強度計算。

2.3.3 總體種植強度

將農(nóng)作物總體種植強度定義為2000—2017年江漢平原耕地中進行過農(nóng)作物種植活動的像元年數(shù)。利用SOM自組織神經(jīng)網(wǎng)絡(luò)對雙季種植強度、平均種植強度、總體種植強度3個指標進行聚類分析，得到江漢平原耕地種植強度分區(qū)圖。該方法將總體種植強度數(shù)據(jù)進行重分類，其中像元值為“0”的區(qū)域代表非農(nóng)用地，值域“1～7”的區(qū)域表示該地區(qū)總體種植強度為“低”，值域“8～14”的區(qū)域表示該地區(qū)總體種植強度為“中”，值域“15～18”的區(qū)域表示該地區(qū)總體種植強度為“高”。

2.4 種植強度分區(qū)

SOM自組織神經(jīng)網(wǎng)絡(luò)又稱為Kohonen網(wǎng)絡(luò)，是芬蘭學者Kohonen在1981年提出的一種自組織競爭人工神經(jīng)網(wǎng)絡(luò)［21］。自組織神經(jīng)網(wǎng)絡(luò)模仿人腦的特點，能夠在沒有先驗知識的情況下自動尋找神經(jīng)網(wǎng)絡(luò)中存在的內(nèi)在規(guī)律和本質(zhì)特征，其本質(zhì)上是一種非監(jiān)督分類［22］。本文利用SOM自組織神經(jīng)網(wǎng)絡(luò)方法對雙季農(nóng)作物種植強度、平均種植強度和總體種植強度進行聚類分析，得到江漢平原耕地種植強度分區(qū)。

3 種植強度空間格局

3.1 雙季種植強度空間格局

2000—2017年江漢平原雙季種植強度分布見圖5、圖6。分析可知，江漢平原雙季種植強度的高值區(qū)域主要位于江陵、潛江、公安北部、松滋東北部、當陽中南部、枝江西部、天門西南部、仙桃西部、云夢西部等地，表明以上區(qū)域為雙季農(nóng)作物穩(wěn)定種植的主要分布區(qū)；中值區(qū)域集中分布于松滋中部以及枝江、荊州、天門等部分區(qū)域；低值區(qū)域主要分布于應城、監(jiān)利中部、仙桃中東部、天門北部、云夢中東部、漢川北部等地。同時，低、中、高種植強度區(qū)域面積占農(nóng)用地總面積的比例分別為42.84%、33.09%、24.07%。其中研究期間均種植雙季農(nóng)作物的區(qū)域僅占5.08%，主要分布于江陵、當陽中南部、天門西南部等地。

圖5 2000—2017年江漢平原雙季種植強度Fig.5 Double-season cropping intensity in Jianghan Plain from 2000 to 2017

圖6 2000—2017年江漢平原雙季種植強度重分類 Fig.6 Reclassification of double-season cropping intensity in Jianghan Plain from 2000 to 2017

3.2 平均種植強度空間格局

2000—2017年江漢平原平均種植強度分布圖如圖7所示。分析發(fā)現(xiàn)，高值區(qū)域和江漢平原雙季種植高強度區(qū)域范圍比較一致，低值區(qū)域和江漢平原水體、天然植被、人工建筑等類別分布范圍比較一致，中間值區(qū)域和江漢平原單季農(nóng)作物分布范圍比較一致。其中平均種植強度高值區(qū)集中在江陵、公安、松滋、石首等地的沿長江區(qū)域，以及當陽東部、潛江、天門西南部、仙桃西部等區(qū)域；中間值區(qū)域集中在枝江西部、松滋西部、監(jiān)利、應城、云夢、仙桃中東部、洪湖中部、荊州區(qū)、沙市區(qū)、公安南部等區(qū)域。

圖7 2000—2017年江漢平原農(nóng)作物平均種植強度空間分布Fig.7 Distribution of average cropping intensity of crops in Jianghan Plain from 2000 to 2017

3.3 總體種植強度空間格局

2000—2017年江漢平原總體種植強度分布見圖8、圖9。分析發(fā)現(xiàn)，江漢平原耕地中絕大部分區(qū)域均位于“高”值區(qū)，其面積占全部耕地面積的比例為83.12%；“中”值區(qū)域零星，其面積占全部耕地面積的比例為10.13%；“低”值區(qū)域主要分布于水域、城鎮(zhèn)、天然植被等非耕地周邊，其面積占全部耕地面積的比例為6.75%?？傮w上，2000—2017年江漢平原農(nóng)用地中有57.11%的區(qū)域每年均有農(nóng)作物種植。

3.4 種植強度分區(qū)

圖10為2000—2017年江漢平原耕地種植強度分區(qū)圖。其中類別C1表示該區(qū)域耕地種植強度高，類別C2表示該區(qū)域耕地種植強度中等，類別C3表示該區(qū)域耕地種植強度較低。分析可知，江陵、潛江、公安北部、松滋東北部、當陽東南部、枝江東部、天門西南部、仙桃西部等地耕地種植強度較高；松滋中部、枝江中東部、荊州中部等地耕地種植強度中等；監(jiān)利中部、公安南部、天門北部、應城、云夢中東部、漢川北部、仙桃東部等地耕地種植強度較低。其中，耕地種植強度較高的類別C1面積占農(nóng)用地總面積的比例為35.64%，耕地種植強度的中等類別C2面積占農(nóng)用地總面積的比例為33.85%，耕地種植強度較低的類別C3面積占農(nóng)用地總面積的比例為30.51%。

圖8 2000—2017年江漢平原農(nóng)作物總體種植強度Fig.8 Distribution of total cropping intensity of crops in Jianghan Plain from 2000 to 2017

圖9 2000—2017年江漢平原農(nóng)作物 總體種植強度重分類 Fig.9 Reclassification of total cropping intensity in Jianghan Plain from 2000 to 2017

圖10 2000—2017年江漢平原耕地種植強度空間分區(qū)Fig.10 Spatial division of cropping intensity in Jianghan Plain from 2000 to 2017

4 結(jié)論

本文利用近20年MODIS數(shù)據(jù)對江漢平原耕地種植強度的空間格局進行了分析。首先對MODIS時序數(shù)據(jù)進行擬合，提取植被物候參數(shù)，并運用非監(jiān)督分類方法得到2000—2017年江漢平原農(nóng)作物種植頻率分布圖；在此基礎(chǔ)上，通過定義雙季種植強度、平均種植強度、總體種植強度3個指標，運用SOM自組織神經(jīng)網(wǎng)絡(luò)分析法得到江漢平原耕地種植強度分區(qū)圖。結(jié)果表明，江漢平原絕大部分地區(qū)均為高種植強度區(qū)域，具備良好的農(nóng)業(yè)生產(chǎn)條件，農(nóng)業(yè)活動強度高。同時傳統(tǒng)的雙季農(nóng)作物種植分布區(qū)的耕地種植強度較高，傳統(tǒng)的單季農(nóng)作物種植分布區(qū)的耕地種植強度較低，表明農(nóng)作物種植頻率與耕地種植強度之間存在明顯的正相關(guān)性。本文構(gòu)建的耕地種植強度的指標，對于及時準確了解耕地種植強度的時空格局，提高耕地資源的集約利用水平具有重要意義。