周 建，張鳳榮，徐 艷，邱孟龍，謝 臻

（1. 陜西師范大學(xué)西北國(guó)土資源研究中心，西安 710119；2. 中國(guó)農(nóng)業(yè)大學(xué)土地科學(xué)與技術(shù)學(xué)院，北京 100193）

0 引 言

半干旱沙區(qū)是指在干濕氣候區(qū)劃上為半干旱區(qū)且土壤質(zhì)地以沙質(zhì)為主的區(qū)域[1]，荒漠化是該區(qū)域最突出的生態(tài)環(huán)境問題[2]。相關(guān)研究表明，自然因素和人類活動(dòng)因素共同導(dǎo)致了土地荒漠化，但是，各區(qū)域引起土地荒漠化的主導(dǎo)因素存在差異[3-4]。就半干旱沙區(qū)而言，土地資源廣袤，自然降水為半干旱，種植業(yè)與養(yǎng)殖業(yè)在經(jīng)濟(jì)上存在互補(bǔ)[5-6]，以上因素使得廣種薄收的土地利用方式是土地荒漠化的主要原因。因此，研究該區(qū)域耕地利用特征對(duì)實(shí)現(xiàn)土地荒漠化的進(jìn)一步防治具有重要意義。

當(dāng)前，相關(guān)研究基于區(qū)域視角，對(duì)半干旱沙區(qū)耕地面積的變化及其空間差異、驅(qū)動(dòng)機(jī)制等進(jìn)行了大量研究[7-11]，指出在人口增長(zhǎng)背景下該區(qū)域耕地面積呈現(xiàn)增長(zhǎng)的趨勢(shì)。在此基礎(chǔ)上，相關(guān)文獻(xiàn)對(duì)半干旱沙區(qū)土地耕作的適宜性及其空間差異[12]、耕地質(zhì)量的變化[13]、耕地資源的利用特征[14]及如何實(shí)現(xiàn)高效利用進(jìn)行了研究[15]，指出不適宜開發(fā)利用為耕地的土地被開發(fā)利用引起了土地荒漠化[16]。同時(shí)，相關(guān)研究對(duì)半干旱沙區(qū)耕地利用的穩(wěn)定性及不穩(wěn)定耕地的類型進(jìn)行了分析[1,17]。盡管相關(guān)文獻(xiàn)對(duì)半干旱沙區(qū)耕地變化及其利用狀況進(jìn)行了較多的分析，但是，很少有研究基于耕地地塊的尺度，研究半干旱沙區(qū)的耕地利用特征。尤其是在中國(guó)土地荒漠化防治呈“整體逆轉(zhuǎn)、局部擴(kuò)張”的背景下[18-19]，研究耕地地塊的利用特征，對(duì)進(jìn)一步防治土地荒漠化具有重要參考。同時(shí)，當(dāng)前中國(guó)農(nóng)用地分等均以耕地圖斑作為評(píng)價(jià)單元[20-22]，研究半干旱沙區(qū)的耕地利用特征，也可為該區(qū)域耕地分等評(píng)價(jià)單元的確定提供參考。

基于此，本研究基于地塊尺度分析半干旱沙區(qū)耕地地表溫度特征，在此基礎(chǔ)上探討耕地地表溫度特征形成的原因，據(jù)此提出進(jìn)一步防治土地荒漠化、確定該區(qū)域耕地分等評(píng)價(jià)單元的建議。

1 研究區(qū)概況

科爾沁左翼后旗（科左后旗，121°30′ E－123°42′ E，42°40′ N－43°42′ N）位于中國(guó)面積最大的科爾沁沙地的東南部，屬于半干旱沙區(qū)。面積為1.1×104km2，海拔在88.5～208.4 m之間。氣候?yàn)闇貛Т箨懶约撅L(fēng)氣候，多年平均降水量為415 mm，主要集中在6－8月份。年平均風(fēng)速3～4 m/s之間，且大風(fēng)天氣主要集中在冬、春季節(jié)。土壤類型以風(fēng)沙土為主，占土壤總面積的68.9%，沙性土壤使土地易風(fēng)蝕沙化?？谱蠛笃斓匦纹骄彛潭ㄉ城?、半固定沙丘、流動(dòng)沙丘相結(jié)合，坨甸相間分布是其主要的地貌特點(diǎn)，坨甸相間的微地形導(dǎo)致了自然降水的重新分布，使得土壤水分地形分異顯著?？谱蠛笃鞂儆趪?guó)家級(jí)貧困縣，農(nóng)牧業(yè)仍為當(dāng)?shù)剞r(nóng)民的主要收入來源。2013年人均土地面積和耕地面積分別為2.81、0.74 hm2，均遠(yuǎn)高于中國(guó)0.70、0.09 hm2的平均水平。耕地作物以旱地玉米為主，2016年全旗糧食作物播種面積1.89×105hm2，其中玉米播種面積1.62×105hm2，占85.7%。2016牧業(yè)年度牲畜頭數(shù)304.29 萬頭（只、口）。50年代末期至80年代末期，荒漠化土地面積占土地總面積的比例由 22%增至48%[2]。圖1為2013年科左后旗遙感影像。

圖1 2013年科左后旗遙感影像及其擦除范圍Fig.1 Remote sensing images of Horqin Left Back Banner in 2013 and deleted parts

2 數(shù)據(jù)與研究方法

2.1 數(shù)據(jù)與處理

本文所用數(shù)據(jù)包括遙感影像、土地利用/覆被、大氣水汽含量、地形高程和土壤化驗(yàn)數(shù)據(jù)。

利用2013年9月12日（行編號(hào)條帶號(hào)為030/120）和2013年9月5日（行編號(hào)條帶號(hào)為030/119）兩幅Landsat－8遙感影像反演科左后旗的地表溫度，兩幅遙感影像覆蓋科左后旗的范圍如圖1。由于2013年9月5日的030/119幅遙感影像存在一定的云量，為避免云對(duì)地表溫度反演產(chǎn)生影響，對(duì)有云的區(qū)域進(jìn)行了刪除[23]。

2013年土地利用/覆被數(shù)據(jù)來源于對(duì)2013年9月12日（行編號(hào)條帶號(hào)030/120）、2013年9月5日（行編號(hào)條帶號(hào)030/119）、2013年7月3日（行編號(hào)條帶號(hào)030/119）共3幅遙感影像的解譯[5]，土地利用/覆被類型包括耕地、林地、草地、建設(shè)用地、沙地、水域。此外，土地利用數(shù)據(jù)還有2012年科左后旗土地利用變更調(diào)查數(shù)據(jù)。

大氣水汽含量為 NASA網(wǎng)站（https://ladsweb.nascom.nasa.gov/data/search.html）的Modis05_L2數(shù)據(jù)，利用 HEGTool（Hdf-Eos to GIS Conversion Tool）將Modis05_L2近紅外大氣水汽含量數(shù)據(jù)的坐標(biāo)轉(zhuǎn)換為UTM投影，并將其1 km空間分辨率利用最近鄰法重采樣為30 m。

地形高程為GDEMV2 30 m（http://www.gscloud.cn/）數(shù)據(jù)，經(jīng)過鑲嵌、裁剪等操作，獲得科左后旗的DEM。

2015年 10月采集科左后旗邊布拉村耕地地塊不同位置的表層土壤樣品，經(jīng)實(shí)驗(yàn)室分析獲得土壤有機(jī)質(zhì)含量數(shù)據(jù)。2014年10月8日至2014年10月16日、2015年7月20日至2015年7月31日對(duì)科左后旗進(jìn)行了有關(guān)土地利用、農(nóng)業(yè)生產(chǎn)、牲畜養(yǎng)殖等方面情況的實(shí)地調(diào)研。

2.2 研究方法

2.2.1 地表溫度遙感反演

對(duì)沙區(qū)而言，沙的比熱容小，受熱后溫度上升較快。因此，沙區(qū)地表覆蓋狀況對(duì)地表溫度影響明顯，地表溫度對(duì)植被長(zhǎng)勢(shì)敏感[24-26]。科左后旗玉米生長(zhǎng)期為 5月 5日至 9月 17日[27]，本文選用玉米生長(zhǎng)末期（2013年 9月12日、2013年9月5日）的兩幅Landsat－8遙感影像反演地表溫度，以更好地體現(xiàn)玉米生長(zhǎng)差異。本文利用Rozenstein 等提出的劈裂窗算法（式1）反演地表溫度[28]。

式中 Ts地表溫度，K；A0、A1、A2是由大氣透射率（atmospheric transmittance）和地表比輻射率（land surface emissivity）決定的系數(shù)，T10、T11是TIRS傳感器10波段和11波段的亮度溫度。TIRS傳感器10波段和11波段的大氣透射率利用中緯度夏季大氣剖面的透射率公式進(jìn)行計(jì)算。對(duì)于TIRS傳感器10波段和11波段的地表比輻射率，首先根據(jù)土地利用/覆被類型將地表劃分為水體、自然表面和城鎮(zhèn)3類，水體和沙地（裸土）的地表比輻射率按照宋挺等給定的參數(shù)進(jìn)行取值[29]，城鎮(zhèn)和自然表面的比輻射率按覃志豪等提出的方法進(jìn)行計(jì)算[30]。TIRS傳感器10波段和11波段的亮度溫度是在計(jì)算其輻射亮溫的基礎(chǔ)上進(jìn)行計(jì)算。地表溫度的具體計(jì)算過程詳見文獻(xiàn)[23]。

2.2.2 耕地地塊地表溫度及其邊緣地表溫度

從2012年土地利用變更調(diào)查數(shù)據(jù)庫中提取耕地，并與 2013年遙感解譯的耕地進(jìn)行空間疊加，以獲得 2013年耕地地塊邊界。利用ArcGIS軟件的區(qū)域統(tǒng)計(jì)功能計(jì)算各地塊溫度范圍。

在分析耕地地塊邊緣溫度時(shí)，為避免耕地地塊面積過小對(duì)結(jié)果造成影響。首先，選擇已獲得2012年耕地邊界的面積大于8100 m2的2013年遙感解譯的耕地地塊。其次，利用緩沖功能向地塊內(nèi)部緩沖30 m，以獲得地塊邊緣部分。最后，利用ArcGIS的區(qū)域統(tǒng)計(jì)功能分析地塊邊緣部分的地表溫度。

2.2.3 高程指數(shù)

利用高程指數(shù)展示坨甸相間的微地形狀態(tài)（式 2）。首先，將科左后旗劃分為2 km×2 km的格網(wǎng)，然后統(tǒng)計(jì)每個(gè)格網(wǎng)的高程平均值，最后計(jì)算高程指數(shù)。某柵格高程指數(shù)越小，表示該柵格在某一格網(wǎng)中地形越低。

式中EI為高程指數(shù)；E為高程，m；ˉE為格網(wǎng)的高程平均值，m。

3 結(jié)果與分析

3.1 耕地地表溫度總體狀況分析

行編號(hào)條帶號(hào)為 030/120部分的耕地地表溫度最小值為294.03 K，最大值為308.61 K，相差14.58 K。就耕地地表溫度的分布而言，地表溫度分布在[299 K，300 K)范圍內(nèi)的面積最大，占耕地總面積的 26.278%（表 1）；占耕地總面積的比例在10%以上的溫度范圍包括[297 K，298 K)、[298 K，299 K)、[299 K，300 K)、[300 K，301 K)、[301 K，302 K)共 5個(gè)溫度范圍，面積占耕地總面積的92.492%，跨度5 K。030/119部分的耕地地表溫度最小值為291.43 K，最大值為310.62 K，相差19.19 K。就耕地地表溫度的分布而言，地表溫度分布在[302 K，303 K)范圍內(nèi)的耕地最多，面積占耕地總面積的30.582%（表1）；占耕地總面積的比例在10%以上的溫度范圍包括[300 K，301 K)、[301 K，302 K)、[302 K，303 K)、[303 K，304 K)共 4個(gè)溫度范圍，面積占耕地總面積的 82.519%，跨度4 K。以上結(jié)果表明，耕地地表溫度分布范圍較大。

表1 耕地地表溫度Table 1 Land surface temperature of cultivated land

3.2 耕地地塊溫差分析

就行編號(hào)條帶號(hào)為 030/120部分的耕地地塊溫差的分布而言，地塊溫差小于1 K的耕地面積占耕地總面積的比例僅為9.687%，地塊溫差小于2 K的耕地面積占總面積的35.656%（表2）。但是，地塊溫差大于2 K的耕地面積占總面積的64.344%，接近耕地總面積的2/3；地塊溫差大于3 K的耕地面積占總面積的36.895%，超過耕地總面積的1/3；地塊溫差大于4 K的耕地面積占總面積的16.841%。就行編號(hào)條帶號(hào)為030/119部分的耕地地塊溫差的分布而言，地塊溫差小于1 K的耕地面積占耕地總面積的比例僅為3.299%，地表溫差小于2 K的耕地面積占總面積的19.526%。相反，地表溫差大于2 K的耕地面積占總面積的80.474%；地表溫差大于3 K的耕地面積占總面積的57.971%；地表溫差大于4 K的耕地面積占總面積的41.571%。反映出科左后旗耕地地塊溫差較大，且地塊溫差較小的耕地面積占總面積的比例較小，地塊溫差較大的耕地面積所占比例較大。

表2 耕地地塊溫差Table 2 Temperature difference of cultivated land parcel

3.3 耕地地塊邊緣溫度分析

在行編號(hào)條帶號(hào)為030/120部分面積大于8100 m2的22107塊耕地中，有13268塊耕地的邊緣地表平均溫度較地塊平均溫度高，數(shù)量占耕地總塊數(shù)的60.02%，面積占耕地總面積的73.71%（表3）。在行編號(hào)條帶號(hào)為030/119部分的面積大于8100 m2的1644塊耕地中，有866塊耕地的邊緣地表平均溫度較地塊平均溫度高，數(shù)量占耕地總塊數(shù)的 52.68%，面積占耕地總面積的 53.80%（表 4）。反映出耕地地塊邊緣的地表平均溫度高于地塊平均溫度。

表3 行編號(hào)條帶號(hào)030/120部分耕地地塊邊緣平均溫度與地塊平均溫度Table 3 Average temperatures at cultivated land parcel margin and cultivated land parcel in row/path 030/120 part

表4 行編號(hào)條帶號(hào)030/119部分耕地地塊邊緣平均溫度與地塊平均溫度Table 4 Average temperatures at cultivated land parcel margin and cultivated land parcel in row/path 030/119 part

4 討 論

4.1 耕地地表溫度范圍較大原因分析

科左后旗位于半干旱沙區(qū)，耕地利用以旱地玉米為主，自然降水量決定了玉米的生長(zhǎng)狀況。根據(jù)相關(guān)文獻(xiàn)，科左后旗自然降水量由東部的 462.8mm下降到西部的358 mm，降水量空間分異明顯[23]，如此便導(dǎo)致了玉米生長(zhǎng)空間分異顯著，耕地地表溫度范圍較大。

4.2 耕地地表溫度非均質(zhì)性分析

坨甸相間是科左后旗主要的地形，坨地為海拔相對(duì)較高的高地，甸地為海拔較低的洼地。沙質(zhì)土壤持水性能差，坨甸微地形對(duì)自然降水進(jìn)行再分配，導(dǎo)致微地形不同部位的土壤水分狀況差異顯著[31]，甸地土壤水分條件相對(duì)最多，坨地土壤水分條件相對(duì)最少[32]。地形相對(duì)越高，作物生長(zhǎng)所需水分滿足程度越低，生長(zhǎng)越差，地表裸露，地表溫度相對(duì)越高。就行編號(hào)條帶號(hào)為030/120部分的耕地地表溫度與高程指數(shù)而言，在[294 K，302 K）范圍內(nèi)，隨著高程指數(shù)的增加，耕地地表溫度逐漸升高（圖2a），此部分耕地面積占030/120部分耕地總面積的93.729%（表1）。就030/119部分的耕地地表溫度與高程指數(shù)而言，在[291 K，307 K）范圍內(nèi)，隨著高程指數(shù)的增加，耕地地表溫度逐漸升高（圖 2b），此部分耕地占030/119部分耕地總面積的99.948%（表2）。表明，在坨甸微地形下，隨著相對(duì)海拔的升高，耕地地表溫度在升高。坨甸微地形的普遍存在，使得耕地地塊并非完全平坦而是存在一定的高差，造成耕地地塊水分條件和作物生長(zhǎng)存在差異，使得耕地地塊溫差較大。

在行編號(hào)條帶號(hào)為030/120部分當(dāng)?shù)乇頊囟却笥?02 K和行編號(hào)條帶號(hào)為030/119部分當(dāng)?shù)乇頊囟却笥?07 K時(shí)，隨著高程指數(shù)的降低，耕地地表溫度反而出現(xiàn)下降的趨勢(shì)。此部分耕地占耕地總面積的比例較小，分別為6.271%和0.052%（表1）。這主要是由于大陸性季風(fēng)氣候條件下坨甸微地形低洼處匯集水分較多，耕地作物生長(zhǎng)在一定程度上受瀝澇影響，作物生長(zhǎng)較差[1]，植被覆蓋度低，地表溫度較高。

盡管《農(nóng)用地分等規(guī)程（GB/T 28407—2012）》中提出耕地質(zhì)量的評(píng)價(jià)單元可以采用疊置法、地塊法、網(wǎng)格法、多邊形法等方法進(jìn)行確定。但是，目前在農(nóng)用地分等的實(shí)際工作中，耕地質(zhì)量評(píng)價(jià)單元的確定基本上采用地塊法，即以土地利用現(xiàn)狀數(shù)據(jù)庫中的耕地圖斑作為耕地質(zhì)量的評(píng)價(jià)單元[20-22]。本研究發(fā)現(xiàn)在半干旱沙區(qū)，耕地地塊溫差較大，耕地地塊作物生長(zhǎng)差異顯著，耕地質(zhì)量存在明顯差別。如果采用地塊法（耕地圖斑）勢(shì)必會(huì)消除耕地地塊內(nèi)部的質(zhì)量差異，建議在半干旱沙區(qū)進(jìn)行耕地質(zhì)量評(píng)價(jià)時(shí)采用網(wǎng)格法確定評(píng)價(jià)單元。

圖2 耕地地表溫度與高程指數(shù)Fig.2 Land surface temperature of cultivated land and elevation index

4.3 耕地地塊邊緣效應(yīng)分析

通常情況下地塊邊緣種內(nèi)競(jìng)爭(zhēng)較地塊內(nèi)部弱，作物生長(zhǎng)獲得的養(yǎng)分等資源多，生長(zhǎng)狀況一般較好[33]，地塊邊緣溫度較地塊平均溫度低。但是，對(duì)科左后旗的實(shí)地調(diào)研發(fā)現(xiàn)，科左后旗耕地承包分配時(shí)均是土壤條件較好的耕地，此后農(nóng)戶普遍以承包地地塊為基礎(chǔ)，開墾周邊的土地以擴(kuò)大耕種面積[6]。根據(jù)耕層土壤有機(jī)質(zhì)分析結(jié)果，1號(hào)點(diǎn)土壤有機(jī)質(zhì)含量為 44.33 g/kg，2號(hào)點(diǎn)為25.92 g/kg。3、4號(hào)點(diǎn)分別為24.29 g/kg、23.41 g/kg，5、6號(hào)點(diǎn)分別為1.88 g/kg、2.56 g/kg（圖3），反映出耕地邊緣的土壤有機(jī)質(zhì)含量較耕地內(nèi)部的土壤有機(jī)質(zhì)含量低，即以承包地地塊為基礎(chǔ)向外擴(kuò)展的耕地土壤質(zhì)量較差。其次，質(zhì)量較高的承包地一般位于地形較低的區(qū)域[34]，向外開墾的土地地勢(shì)較高[35]，耕地的土壤水分條件較差。第三，由于承包地周邊開墾的耕地質(zhì)量差，農(nóng)戶對(duì)開墾耕地的利用管理更加粗放[6]。以上因素導(dǎo)致了半干旱沙區(qū)耕地地塊邊緣作物生長(zhǎng)較內(nèi)部差，地塊邊緣溫度較地塊內(nèi)部溫度高。由于耕地地塊周邊作物生長(zhǎng)差，地表裸露，土地易于風(fēng)蝕沙化，阻礙了生態(tài)環(huán)境的進(jìn)一步恢復(fù)。為此，承包地地塊周邊開墾的新增耕地應(yīng)當(dāng)進(jìn)行退耕，以進(jìn)一步防治土地荒漠化和恢復(fù)生態(tài)。

圖3 耕地表層土壤采樣點(diǎn)Fig.3 Soil sampling point at cultivated land surface

4.4 遙感反演地表溫度的準(zhǔn)確性及其反映作物生長(zhǎng)狀況的有效性分析

本文分析了耕地的地表溫度特征及其形成原因。在遙感反演地表溫度的準(zhǔn)確性方面，本文利用遙感反演的方法計(jì)算的是科左后旗2013年的地表溫度，由于無法對(duì)2013年的地表溫度進(jìn)行實(shí)測(cè)，使得無法利用實(shí)測(cè)地表溫度對(duì)遙感反演的地表溫度進(jìn)行驗(yàn)證。但是，相關(guān)研究表明本文所采用的劈裂窗算法反演地表溫度的精度較高，Rozenstein等采用本方法反演的地表溫度的均方根誤差為0.93 ℃[28]。宋挺等對(duì)利用Landsat－8遙感影像反演地表溫度的精度進(jìn)行分析，表明Jiménez-Mun?oz劈裂窗算法和Rozenstein劈裂窗算法的精度較高，誤差僅為0.7 K左右[29]。遙感反演地表溫度反映作物生長(zhǎng)差異的有效性方面，首先沙區(qū)土壤質(zhì)地為沙質(zhì)，受熱后地表溫度上升快。其次植被生長(zhǎng)狀況對(duì)沙區(qū)地表溫度上升具有顯著影響[24-26,36]。植被覆蓋度高的土地較植被覆蓋度低的土地地表溫度上升小。相關(guān)文獻(xiàn)對(duì)科爾沁沙地坨甸地形耕地作物生長(zhǎng)狀況與耕地地表溫度進(jìn)行研究，表明隨著坨甸地形的升高土壤水分條件變差，耕地作物生長(zhǎng)也在變差，耕地地表溫度隨之升高[23]。因此遙感反演地表溫度可以表征半干旱沙區(qū)耕地作物的生長(zhǎng)差異。但是，一方面本文所采用的Landsat－8遙感影像的獲取時(shí)間為2013年，另一方面由于未在遙感影像獲取的時(shí)間進(jìn)行作物生長(zhǎng)狀況的采樣，導(dǎo)致無法對(duì)遙感反演地表溫度與作物生長(zhǎng)狀況的關(guān)系進(jìn)行分析，未來將對(duì)沙區(qū)遙感反演的地表溫度與作物生長(zhǎng)狀況之間的關(guān)系進(jìn)行深入研究。

5 結(jié) 論

1）就耕地地表溫度的總體狀況而言，耕地地表溫度范圍較大。在行編號(hào)條帶號(hào)為030/120部分，占耕地總面積比例在 10%以上的溫度范圍包括[297 K，298 K)、[298 K，299 K)、[299 K，300 K)、[300 K，301 K)、[301 K，302 K)，跨5個(gè)溫度范圍，面積占耕地總面積的92.492%；在行編號(hào)條帶號(hào)為030/119部分，占耕地總面積的比例在10%以上的溫度范圍包括[300 K，301 K)、[301 K，302 K)、[302 K，303 K)、[303 K，304 K)，跨4個(gè)溫度范圍，面積占耕地總面積的82.519%。

就耕地地塊溫差而言，耕地地塊溫差較大，且地塊溫差較小的耕地面積占總面積的比例較小，地塊溫差較大的耕地面積所占比例較大。

2）耕地地塊邊緣的地表平均溫度高于地塊平均溫度。

3）基于以上研究，認(rèn)為在半干旱沙區(qū)進(jìn)行耕地質(zhì)量評(píng)價(jià)時(shí)，評(píng)價(jià)單元的確定宜采用網(wǎng)格法；建議對(duì)耕地地塊邊緣的耕地進(jìn)行退耕，以進(jìn)一步防治土地荒漠化和恢復(fù)生態(tài)。