基于衛(wèi)星遙感的臨沂市冬小麥種植時(shí)空變化特征分析

李莉，劉志紅，韓晨琛，申高航

(1.成都信息工程大學(xué)，四川 成都 610225；2.臨沂市氣象局，山東 臨沂 276000)

引言

臨沂地區(qū)是山東省的糧食主產(chǎn)區(qū)之一，一年兩熟，主要以冬小麥、夏玉米、棉花、花生為主，其中冬小麥產(chǎn)量占糧食總產(chǎn)量的40%～50%[1-2]。近年來(lái)，受經(jīng)濟(jì)作物和政策影響，臨沂地區(qū)冬小麥種植區(qū)域出現(xiàn)波動(dòng)，因此，快速、準(zhǔn)確地獲取臨沂冬小麥種植區(qū)域信息，把握其種植結(jié)構(gòu)的空間變化，對(duì)臨沂種植結(jié)構(gòu)調(diào)整、糧食安全和農(nóng)業(yè)發(fā)展具有重要的數(shù)據(jù)支撐作用[3-4]。

遙感監(jiān)測(cè)作物種植區(qū)域，了解農(nóng)作物的類型、結(jié)構(gòu)和分布特征是宏觀把握農(nóng)作物布局，開(kāi)展農(nóng)業(yè)減災(zāi)等的重要依據(jù)，也是進(jìn)行農(nóng)業(yè)結(jié)構(gòu)調(diào)整的基礎(chǔ)[5-8]。準(zhǔn)確、連續(xù)、多年遙感監(jiān)測(cè)作物空間分布對(duì)面積統(tǒng)計(jì)[9]起著重要的作用。目前關(guān)于農(nóng)作物時(shí)空格局變化的研究主要是利用農(nóng)作物生育期內(nèi)增強(qiáng)型植被指數(shù)EVI和歸一化植被指數(shù)NDVI的特點(diǎn)，依據(jù)空間分析原理和數(shù)理統(tǒng)計(jì)方法，結(jié)合統(tǒng)計(jì)年鑒和衛(wèi)星遙感數(shù)據(jù)，研究水稻、玉米、棉花、小麥等主要農(nóng)作物的種植分布信息，為研究區(qū)的農(nóng)業(yè)結(jié)構(gòu)以及生產(chǎn)管理提供有效依據(jù)[10-14]。多集中于較大尺度和空間面積，對(duì)于以地級(jí)市為分析單元的研究還較少，但農(nóng)業(yè)種植不同地區(qū)政策導(dǎo)向不同，除了需要在國(guó)家層面、大區(qū)層面了解農(nóng)業(yè)種植，更需要對(duì)地市縣進(jìn)行剖析，及時(shí)調(diào)整基層單元的種植結(jié)構(gòu)。

本文以臨沂地區(qū)冬小麥為研究對(duì)象，利用Terra/MOD13Q1遙感植被數(shù)據(jù)，充分考慮遙感的光譜特征和冬小麥的物候特征，快速提取冬小麥種植區(qū)域信息，在此基礎(chǔ)上分析2009—2019年臨沂地區(qū)冬小麥種植的時(shí)空分布特征，穩(wěn)定度種植特征，以及變化波動(dòng)的原因，對(duì)臨沂地區(qū)的農(nóng)業(yè)產(chǎn)業(yè)政策制定及結(jié)構(gòu)調(diào)整起到積極的參考作用。

1 數(shù)據(jù)與方法

1.1 研究區(qū)域

臨沂地區(qū)位于山東省東南部，地跨34°22′～36°13′N，117°24′～119°11′E，總面積17 191.2 km2，轄蘭山等3區(qū)和郯城等9縣，區(qū)域內(nèi)北部主要以山地、丘陵為主，南部以平原為主。此地區(qū)屬暖溫帶大陸性季風(fēng)氣候，陽(yáng)光充足，雨量豐沛，有利于農(nóng)作物、植被的生長(zhǎng)發(fā)育，是糧食作物和經(jīng)濟(jì)作物的重要原產(chǎn)地[15]。

1.2 數(shù)據(jù)收集及預(yù)處理

1)時(shí)間序列的Terra/MOD13Q1遙感植被數(shù)據(jù)產(chǎn)品，時(shí)間覆蓋范圍為2009—2019年冬小麥生長(zhǎng)發(fā)育關(guān)鍵期(3—5月)，共計(jì)有55景數(shù)據(jù)，空間分辨率為250 m，從MODIS官網(wǎng)(https://ladsweb.nascom.nasa.gov/search/)獲得。2)DEM數(shù)據(jù)，30 m分辨率，來(lái)源于地理空間數(shù)據(jù)云網(wǎng)站(http://www.gscloud.cn/)。3)Google Earth目視解譯冬小麥數(shù)據(jù)，主要用來(lái)驗(yàn)證冬小麥面積提取結(jié)果。4)冬小麥物候歷數(shù)據(jù)，來(lái)源于中華人民共和國(guó)農(nóng)業(yè)部種植業(yè)管理司(http://www.zzys.moa.gov.cn /)。5)冬小麥種植面積統(tǒng)計(jì)數(shù)據(jù)，以市域?yàn)閱挝?，時(shí)間范圍為2009—2019年，來(lái)源于臨沂地區(qū)統(tǒng)計(jì)年鑒。

對(duì)遙感數(shù)據(jù)進(jìn)行一些預(yù)處理，包括拼接、裁剪、投影轉(zhuǎn)換(選擇正軸等積雙標(biāo)準(zhǔn)緯線圓錐投影，即Albers)、去除異常值、重分類等流程，最終生成便于研究的柵格格式的NDVI時(shí)間序列數(shù)據(jù)。

1.3 冬小麥種植信息遙感提取

不同的農(nóng)作物有不同的生長(zhǎng)規(guī)律，比如：播種、出苗的日期不同，什么時(shí)候生長(zhǎng)最旺盛，什么時(shí)候開(kāi)始成熟、收獲，以及同一種作物在同一地區(qū)具有相對(duì)穩(wěn)定的生長(zhǎng)發(fā)育規(guī)律，各個(gè)發(fā)育期的時(shí)間段相對(duì)固定等[16]。由于不同作物的物候期不同，同一時(shí)期不同農(nóng)作物的植被指數(shù)是不同的，導(dǎo)致反映作物生長(zhǎng)的植被指數(shù)(NDVI)的時(shí)間曲線形態(tài)不同，并按一定的規(guī)律變化[17]。因此，根據(jù)冬小麥生長(zhǎng)過(guò)程的典型特征，結(jié)合臨沂地區(qū)冬小麥物候歷和平均NDVI時(shí)間序列曲線，可以得到臨沂地區(qū)冬小麥空間分布信息提取模型。

1.4 精度驗(yàn)證

從以下兩個(gè)方面進(jìn)行精度驗(yàn)證：

1)空間位置上的精度驗(yàn)證。為有效地對(duì)提取結(jié)果進(jìn)行驗(yàn)證，采用地面樣點(diǎn)和Google Earth影像目視采樣點(diǎn)兩種方式對(duì)提取結(jié)果進(jìn)行空間分布上的驗(yàn)證。

2)數(shù)量上，主要用統(tǒng)計(jì)數(shù)據(jù)來(lái)進(jìn)行驗(yàn)證。以臨沂地區(qū)統(tǒng)計(jì)年鑒數(shù)據(jù)作為評(píng)價(jià)標(biāo)準(zhǔn)，對(duì)提取出來(lái)的2009—2019年臨沂地區(qū)冬小麥種植面積進(jìn)行精度評(píng)價(jià)。

1.5 時(shí)空演變?cè)u(píng)價(jià)指標(biāo)

采用變化幅度和動(dòng)態(tài)度2個(gè)指標(biāo)[18]評(píng)價(jià)臨沂地區(qū)冬小麥3個(gè)時(shí)期、2個(gè)時(shí)間段時(shí)空演變規(guī)律。變化幅度，即變化率，主要用于評(píng)價(jià)在研究時(shí)段內(nèi)不同土地類型的面積變化，反映其變化的幅度和速度[19]。本文的研究對(duì)象為冬小麥，所以在一定時(shí)段內(nèi)其面積的變化速度又可以用動(dòng)態(tài)度來(lái)表示，也就是單一土地利用類型的動(dòng)態(tài)度[20]。動(dòng)態(tài)度和變化幅度的計(jì)算公式如下：

(1)

(2)

式(1)、(2)中：KT為冬小麥面積在某一時(shí)間段內(nèi)的變化幅度或變化率；Ua和Ub分別表示研究時(shí)段初期、末期的冬小麥面積(單位：km2)；K為冬小麥面積的動(dòng)態(tài)度，即年均變化率；T為研究時(shí)段的時(shí)間間隔(單位：a)。

2 結(jié)果與分析

2.1 模型參數(shù)確定與結(jié)果驗(yàn)證分析

2.1.1 參數(shù)確定

為了準(zhǔn)確提取冬小麥NDVI日序曲線，從研究區(qū)選取10個(gè)冬小麥地面純像元點(diǎn)，并保證樣點(diǎn)250 m范圍內(nèi)均為冬小麥，利用研究時(shí)段內(nèi)每年11月至次年5月的Terra/MOD13Q1植被數(shù)據(jù)產(chǎn)品計(jì)算NDVI，提取各個(gè)樣點(diǎn)的NDVI曲線，應(yīng)用所有樣點(diǎn)的NDVI平均值來(lái)構(gòu)造平均NDVI時(shí)序曲線(圖1)。

圖1 臨沂地區(qū)冬小麥平均NDVI日序曲線

利用MODIS數(shù)據(jù)獲取冬小麥種植面積、空間分布信息需要幾個(gè)關(guān)鍵發(fā)育期(返青、拔節(jié)、孕穗、抽穗期)的遙感影像，這時(shí)候同一時(shí)期的林地、草地等植被雖然已經(jīng)返綠，但NDVI值都低于冬小麥，同時(shí)期的非植被的NDVI值更是遠(yuǎn)低于冬小麥，冬小麥處于拔節(jié)、抽穗期的時(shí)候，棉花、水稻等其他作物也才剛剛播種，冬小麥信息在遙感影像上最為突出[21]。所以，可以選擇冬小麥的關(guān)鍵生育期作為研究時(shí)段[22]。為了模型便于推廣，用NDVIxy表示不同時(shí)段的NDVI值，其中x取值范圍為1～12，代表月份，y取值范圍為1～3，分別表示上、中、下旬，如NDVI33表示3月下旬的NDVI值[23]。冬小麥時(shí)間序列遙感監(jiān)測(cè)提取的一般模型(一般提取耕地中同時(shí)滿足以下條件的晴空天氣的像元為冬小麥)為：NDVI31>T1，NDVI32>NDVI31，NDVI41>NDVI32，NDVI43>NDVI41，NDVI43>T2。

因?yàn)镸OD13Q1的植被數(shù)據(jù)產(chǎn)品，時(shí)間分辨率是16 d，所以一個(gè)月只有兩景影像，不是每旬都有數(shù)據(jù)。模型中T1、T2為不同生育期的閾值(表1)。由于受播種時(shí)間、物候差異、氣候變化、作物長(zhǎng)勢(shì)差異等影響，NDVI值具有區(qū)域特性，不同監(jiān)測(cè)年份、不同區(qū)域閾值不同，需要根據(jù)每年的物候歷或農(nóng)情野外監(jiān)測(cè)數(shù)據(jù)進(jìn)行修正。

表1 2009—2019年臨沂地區(qū)冬小麥不同生育期閾值

2.1.2 結(jié)果驗(yàn)證分析

2.1.2.1 空間分布位置的驗(yàn)證

地面樣點(diǎn)數(shù)據(jù)包括2009—2019年臨沂市氣象局農(nóng)業(yè)氣象觀測(cè)數(shù)據(jù)、野外調(diào)查觀測(cè)數(shù)據(jù)以及Google Earth目視采樣點(diǎn)數(shù)據(jù)，其中農(nóng)業(yè)氣象站的觀測(cè)數(shù)據(jù)主要用于冬小麥發(fā)育期的確定。野外觀測(cè)中共采集了5個(gè)GPS樣點(diǎn)，鑒于人力、物力有限，選取的野外調(diào)查點(diǎn)數(shù)量較少，故又在Google Earth 高分辨率影像中隨機(jī)選取了35個(gè)樣點(diǎn)，基本覆蓋了冬小麥主要種植區(qū)(圖2)。將這些數(shù)據(jù)進(jìn)行空間化和數(shù)字化處理，用來(lái)輔助驗(yàn)證臨沂地區(qū)冬小麥種植信息空間分布的提取結(jié)果。結(jié)果表明：40個(gè)冬小麥樣點(diǎn)中，38個(gè)在影像中被正確分類為冬小麥，空間匹配精度為95%。

2.1.2.2 基于統(tǒng)計(jì)數(shù)據(jù)的數(shù)量驗(yàn)證

用以上方法得到2009—2019年臨沂地區(qū)冬小麥空間分布(圖3)信息，且提取面積與統(tǒng)計(jì)面積具有較高的吻合度，提取精度在97%以上(表2)。

圖2 驗(yàn)證點(diǎn)空間位置分布

圖3 2009—2019年臨沂地區(qū)冬小麥種植的空間特征信息(a. 2009年，b. 2010年，c. 2011年，d. 2012年，e. 2013年，f. 2014年，g. 2015年，h. 2016年，i. 2017年，j. 2018年，k. 2019年)

表2 基于統(tǒng)計(jì)數(shù)據(jù)的精度驗(yàn)證

2.2 冬小麥種植的空間特征信息

圖3為基于遙感植被數(shù)據(jù)產(chǎn)品提取出來(lái)的2009—2019年臨沂冬小麥種植的空間特征信息。

2009年，臨沂地區(qū)冬小麥種植面積2 838 km2，占行政面積的16.5%，南部平原的蘭陵、臨沭和郯城為主要分布區(qū)，合計(jì)占冬小麥總面積的59%。西北部費(fèi)縣、蒙陰、平邑、沂水4個(gè)縣丘陵山地較多，冬小麥相對(duì)較少，占比均小于10%。

2019年，臨沂地區(qū)冬小麥種植面積2 985 km2，占行政面積的17.4%?？臻g分布格局與2009年一致，蘭陵、臨沭和郯城冬小麥面積合計(jì)占冬小麥總面積的66%。西北部4縣的冬小麥依然相對(duì)較少。

2.3 冬小麥種植穩(wěn)定度分析

從不同區(qū)域冬小麥種植年數(shù)(圖4)來(lái)看，2009—2019年，持續(xù)種植冬小麥10 a及以上的區(qū)域總面積為1 638 km2，主要分布在蘭陵和郯城，占比均超過(guò)30%，穩(wěn)定度高；連續(xù)種植冬小麥7～9 a的區(qū)域總面積為976 km2，主要分布在郯城、蘭陵和臨沭，較為穩(wěn)定；持續(xù)種植4～6 a的區(qū)域主要分布在蘭陵和臨沭，面積955 km2，具有一定的穩(wěn)定度；持續(xù)種植小于4 a的區(qū)域主要分布在費(fèi)縣和莒南，面積為1 914 km2，種植不穩(wěn)定。近10年來(lái)，蘭陵、郯城冬小麥種植分布格局相對(duì)穩(wěn)定，連續(xù)種植10 a及以上的區(qū)域面積占各自冬小麥種植總面積的50%左右。由此可見(jiàn)，平原地區(qū)冬小麥種植區(qū)域相對(duì)穩(wěn)定，種植時(shí)間相對(duì)較長(zhǎng)；而丘陵山地、城中心地區(qū)，冬小麥的種植區(qū)域不穩(wěn)定，種植時(shí)間較短，變化較大。

圖4 臨沂地區(qū)2009—2019年不同地塊冬小麥種植年數(shù)分布

2.4 冬小麥種植面積的時(shí)空演變

2.4.1 時(shí)間變化

臨沂地區(qū)三區(qū)九縣冬小麥種植面積近10年期間均有不同程度的變化，但起伏都不大(圖5)。從圖上可以明顯看出，2009—2019年，郯城和蘭陵地區(qū)的冬小麥種植面積一直遠(yuǎn)高于其他地區(qū)，其次是臨沭，是臨沂冬小麥的主產(chǎn)區(qū)；雖然這三個(gè)縣的冬小麥面積變化起伏不一，但總體上都呈不同程度的增加趨勢(shì)。而蒙陰和沂水則相反，種植面積一直低于其他縣區(qū)，且與2009年相比，2019年的面積均略有減少。

圖5 2009—2019年臨沂地區(qū)各縣區(qū)冬小麥種植面積年際變化

近10年來(lái)，臨沂地區(qū)冬小麥種植面積變化率最大的是平邑(圖6)，但面積增加并不多，只增加了34 km2，由2009年占全市行政面積的1%增長(zhǎng)到2019年的3%。2009—2019年冬小麥種植面積增加最多的是臨沭，由2009年占全市行政面積的33%增長(zhǎng)到2019年的48%；其次是蘭陵，增加了111 km2。

圖6 2009—2019年臨沂地區(qū)冬小麥種植面積變化

圖7為臨沂市三區(qū)九縣冬小麥面積動(dòng)態(tài)度空間分布。2009—2019年，臨沂冬小麥種植總面積的變化幅度和年均變化率分別為5.5%和0.5%；各縣區(qū)變化趨勢(shì)不一致，年均變化率也相差較大，冬小麥面積增加主要在平邑和臨沭，年均變化率為13%、5%；冬小麥面積減少主要在蘭山和費(fèi)縣，年均變化率分別為-4%、-3%。2009—2012年，三區(qū)九縣冬小麥面積整體呈增加趨勢(shì)，增加主要在平邑、費(fèi)縣、蒙陰；2012—2019年，冬小麥面積各縣區(qū)整體呈下降趨勢(shì)，年均變化率在-9%～-1%之間。

圖7 2009—2019年臨沂地區(qū)冬小麥面積動(dòng)態(tài)度空間分布

2.4.2 空間變化

2009—2019年臨沂地區(qū)冬小麥種植空間位置上略有轉(zhuǎn)移，但總體分布情況幾乎不變，仍呈南多北少的格局。圖8直觀表示了臨沂地區(qū)2019年冬小麥種植分布較2009年的變化情況。其中，冬小麥種植面積新增了367 km2，占總增多面積的87%，新增區(qū)域主要分布在蘭陵、臨沭和莒南，呈明顯的擴(kuò)張趨勢(shì)。冬小麥面積減少達(dá)220 km2，主要分散分布在臨沂三區(qū)和費(fèi)縣。

圖8 2009—2019年臨沂地區(qū)冬小麥空間分布變化(a. 2009—2012年，b. 2012—2019年，c. 2009—2019年)

2009—2012年冬小麥種植面積新增最多的是蘭陵，占總增加面積的18%，其次是臨沭、費(fèi)縣、莒南；新增最少的是蒙陰，只占總增加面積的0.9%；種植區(qū)域消失最多的是臨沭，其次是郯城、河?xùn)|、莒南、蘭陵；消失的最少的是蒙陰。4年來(lái)，持續(xù)種植冬小麥的區(qū)域面積共計(jì)2 274 km2，其中蘭陵所占比例最大，為26%，蒙陰最小，僅占0.1%。2009—2012年，冬小麥種植總面積減少的只有河?xùn)|區(qū)，其他各區(qū)縣的冬小麥種植面積均呈不同程度的增加趨勢(shì)。

2012—2019年冬小麥種植面積新增最多的區(qū)縣是臨沭，占總增加面積的27%，其次是蘭陵、郯城、莒南，增加最少的是蒙陰，只占總增加面積的0.5%。2012—2019年冬小麥種植區(qū)域減少最多的是費(fèi)縣，最少的是蒙陰。7年來(lái)，持續(xù)種植冬小麥的區(qū)域面積共計(jì)2 207 km2，主要分布在蘭陵。

3 冬小麥種植區(qū)域時(shí)空變化特征的原因分析

3.1 空間分布原因分析

臨沂地區(qū)冬小麥種植區(qū)域相對(duì)集中，呈南多北少的分布格局，空間分布變化主要有以下2個(gè)特點(diǎn)：1)從城區(qū)到郊外，冬小麥種植面積減小幅度趨于緩和，離城區(qū)越近的地方冬小麥種植面積減小越迅速，而且更加明顯；2)海拔越高的地方冬小麥種植面積減小的速度越慢，反之，海拔越低的地方冬小麥種植面積減少得越快。該結(jié)論與前人的研究結(jié)果[24]一致。這主要是因?yàn)榕R沂地區(qū)北部的蒙陰、沂水以山地、丘陵為主，不利于種植冬小麥，而處于南部的郯城、蘭陵和臨沭，以平原為主，適宜農(nóng)作物種植。

3.2 種植穩(wěn)定度原因分析

城市周邊地區(qū)冬小麥種植穩(wěn)定度相對(duì)較差，這主要是因?yàn)榻陙?lái)我國(guó)經(jīng)濟(jì)持續(xù)快速發(fā)展，城鎮(zhèn)化水平不斷提高，城鎮(zhèn)人口比例增加，城區(qū)不斷向四周擴(kuò)展，耕地逐漸被占用，大量耕地被轉(zhuǎn)化為建設(shè)用地，導(dǎo)致冬小麥種植區(qū)域大幅減少。而處于城市邊緣的郯城、蘭陵種植穩(wěn)定度則相對(duì)較好，50%左右的種植區(qū)域連續(xù)種植冬小麥達(dá)到10 a及以上，這主要是由于該區(qū)域陽(yáng)光充足、雨量充沛、地勢(shì)平坦，有利于農(nóng)作物、植被的生長(zhǎng)發(fā)育。

3.3 種植面積變化趨勢(shì)分析

臨沂地區(qū)三區(qū)九縣冬小麥種植面積近10年期間均有不同程度的變化，作為臨沂主城區(qū)的蘭山、河?xùn)|和羅莊，近10年來(lái)的種植總面積呈減少趨勢(shì)，這主要是因?yàn)槌鞘谢M(jìn)程加快，主城區(qū)人口不斷增多，城市擴(kuò)張。費(fèi)縣和沂南種植面積減少也比較明顯，這可能是由于該區(qū)鄉(xiāng)鎮(zhèn)企業(yè)發(fā)達(dá)，且冬小麥種植費(fèi)時(shí)費(fèi)力，農(nóng)作物灌溉費(fèi)用較高，經(jīng)濟(jì)效益下降，導(dǎo)致農(nóng)民種植的積極性下降，出現(xiàn)“棄耕”的現(xiàn)象，原來(lái)冬小麥種植用地轉(zhuǎn)為種植蔬菜瓜果等經(jīng)濟(jì)效益較高的作物或轉(zhuǎn)化為建設(shè)用地[25-26]。而蘭陵和臨沭總體上都呈不同程度的增加，呈明顯的擴(kuò)張趨勢(shì)，這可能與土地整理復(fù)墾、中低產(chǎn)田改造有關(guān)。其余地形地貌相近的縣區(qū)冬小麥種植面積近10年情況相似，均變化不大。

4 結(jié)論

1)基于Terra/MOD13Q1遙感植被數(shù)據(jù)產(chǎn)品，通過(guò)分析冬小麥的遙感光譜特征、物候歷特征及其生物學(xué)特性，建立臨沂地區(qū)冬小麥時(shí)間序列遙感監(jiān)測(cè)提取的一般模型，空間匹配精度為95%，基于統(tǒng)計(jì)數(shù)據(jù)的提取精度達(dá)到97%，精度較高，滿足研究需要。

2)近10年來(lái)，郯城和蘭陵冬小麥種植穩(wěn)定度較高，冬小麥種植區(qū)域中有近一半的面積連續(xù)種植達(dá)到10 a及以上；而費(fèi)縣和莒南穩(wěn)定度較低，冬小麥種植區(qū)域連續(xù)種植年數(shù)大部分在4 a以下。

3)2009—2019年，臨沂地區(qū)冬小麥種植在空間分布上略有變化，但總體變化較小，仍呈南多北少的分布格局，種植區(qū)域相對(duì)集中，主要分布在南部蘭陵、臨沭、郯城等平原地區(qū)。

4)2009—2019年臨沂地區(qū)冬小麥種植面積略有增加，由2009年占全市行政面積的16.5%增長(zhǎng)到2019年的17.4%，變化幅度和年均變化率分別為5.2%、0.5%。