李守娟,楊 磊,陳利頂,趙方凱,孫 龍
1 中國科學院生態環境研究中心城市與區域生態國家重點實驗室,北京 100085 2 中國科學院大學, 北京 100049
土壤中的碳和氮共同調節和維持著生態系統的生產力和穩定性,碳氮循環是生態系統的重要過程[1]。城市化通過改變城市周邊的土地利用,改變了土壤中碳氮循環在內的一系列生物地球化學過程[2]。作為城市生態系統和農業/自然生態系統之間的過渡區域,城郊是受城市化過程影響最為直接和強烈的區域,城郊地區土地利用方式的改變,也帶來了土壤中碳氮儲量的變化。例如,城市化造成森林等生態用地減少,森林轉變為農地可使土壤碳庫減少25%,此外,長期農業耕作的影響使得城郊土壤中碳、氮儲量逐漸減少[3- 5]。系統闡明城市化過程中的土地利用變化及其對土壤碳氮儲量的影響,可有效揭示土壤肥力提供等關鍵生態系統服務對城市化過程的響應,也是快速城市化背景下土壤安全研究的重要內容[2]。
土壤中碳、氮含量的變化過程非常緩慢,傳統的定位實驗等研究方法在時間和空間上具有很大的局限性,難以揭示城市化對土壤碳、氮儲量的影響過程[6]。隨著研究手段的發展,土壤碳氮模型已經成功應用于長期定位實驗觀測數據的整合預測[7]。其中,DNDC(DeNitrification—DeComposition)模型的模擬精度較高,并且輸入參數簡單,能準確模擬溫室氣體排放[8]和土壤中碳、氮變化過程[9-10],適用于多種氣候條件下的多尺度農業生態系統[11],在全球許多地區得以驗證和應用。在DNDC模型基礎上開發出來的Forest-DNDC模型可應用于多種森林生態系統的生長和土壤碳、氮過程模擬[12-13]。雖然利用DNDC模型已經開展了若干景觀和流域尺度土壤碳、氮儲量的研究[8,14],但多限于單一土地利用下的土壤碳氮儲量模擬,而土地利用變化對土壤碳、氮儲量影響的研究亟需開展?;贒NDC模型模擬城郊地區土地利用變化驅動下的土壤碳、氮儲量響應,可有效揭示快速城市化過程對土壤生態系統關鍵過程的影響,對認識生態系統服務的演變機制具有重要意義。
長三角地區是我國城市化發展最快的區域之一,本研究以長三角典型城郊小流域浙江寧波樟溪流域為研究區,通過20世紀70年代以來土地利用變化過程的解析并結合DNDC模型的模擬,系統研究了城郊土地利用變化及其土壤碳、氮儲量的響應過程,以期為快速城市化過程中的土壤安全維護和生態系統服務提升提供科學依據。
研究區位于浙江省寧波市城郊的樟溪流域(29°44′—29°51′N,121°13′—121°21′E),流域面積91.6 km2(圖1)。研究區屬亞熱帶季風氣候,年平均氣溫16.3 °C,年日照時數1800 h,無霜期228 d。流域位于城鎮向自然生態系統的過渡地帶,土地利用類型多樣且結構復雜,主要土地利用類型為農地、林地、園地、水域、城鎮建設用地,主要土壤類型為水稻土、紅壤和黃紅壤。其中,農地主要種植貝母、花生、蔬菜等經濟作物,園地主要種植果樹及觀賞性苗木,林地是亞熱帶長綠闊葉林,以香樟、麻櫟等為主。

圖1 研究區位置及土地利用(2015年)Fig.1 Location of study area and land use in 2015
系統收集了研究區1974至2015年MSS、TM、ETM遙感影像,根據野外實地調查并參考Google Earth在ArcGIS 10.2中目視解譯,識別不同時期研究區土地利用類型及其空間分布特征,獲取了1974、1977、1981、1985、1990、1995、2000、2005、2010、2015年土地利用數據。
在研究區內根據土地利用、海拔、土壤類型、管理措施等布設樣點82個(圖1),采樣點包括農地樣點41個、園地樣點11個和林地樣點30個,城鎮建設用地和水域未做采樣。在各樣點采用“S”形布點采集0—10 cm土壤樣品,測定質地、有機質、銨態氮、硝態氮、pH和容重。其中,土壤質地采用激光粒度儀測定,有機質采用重鉻酸鉀法測定,銨態氮和硝態氮采用連續流動分析儀測定,pH采用電位法(水土比為2.5∶1)測定,容重采用環刀法測定[15]。土壤數據用于DNDC模型的校正和參數輸入[16]。
DNDC模型包括兩個部分,由6個子模型組成。第一部分包括土壤-大氣模型、作物生長發育模型、有機質分解模型,其作用是利用輸入的環境因素和人類活動數據來反映各土壤環境因子的動態變化;第二部分包括硝化模型、反硝化模型、發酵模型,其作用是估算生態系統中N2O、CH4、CO2氣體的排放量。應用DNDC模型模擬生物地球化學過程時,需要根據當地種植耕作情況輸入氣象、土壤以及作物等數據,便可進行一年至多年的模擬[17]。DNDC模型的模擬以日為單位,根據每個樣點的數據計算點位尺度的土壤碳氮含量,并可從點位尺度擴展到區域尺度。由點位尺度擴展到區域尺度是將模擬區域劃分為多個單元,每一單元內部各種條件都是均勻的,DNDC對所有單元進行逐一模擬,最后將各單元疊加在一起[18]。
樟溪流域是浙貝母主要產地,貝母是研究區農地種植的主要作物。在2017年5月—7月間,本研究通過入戶式問卷調查和查閱文獻資料的方式獲取農戶貝母種植模式信息。問卷包含貝母種植模式中各項農地管理措施的基礎信息,由此建立貝母種植模式數據庫,作為分析計算的數據基礎。問卷調查對象包括當地政府農業部門、種植大戶和以家庭為單位的種植戶等。其中農業部門獲取的資料用以代表區域貝母種植的平均狀況,具有專業農地管理技術的種植大戶提供了相對先進的貝母種植技術數據,而以家庭為單位的種植戶則提供了目前普遍使用的貝母管理措施信息。文獻資料則包括貝母的基本生理參數和科學的種植方法,用以對調查問卷獲取的數據進行輔助和修正。
表1為DNDC模型所需參數和本研究中的數據來源,本研究采用DNDC模擬輸出的作物產量與實測值進行比較,來驗證DNDC模型在研究區的模擬精度。因DNDC模擬的作物產量單位為kgC/hm2,為便于分析實測結果與模擬結果之間的差異,采用每公頃作物產量中碳含量進行校正。各土地利用類型的土壤基本理化性質及主要管理措施如表2、表3所示。
采用決定系數R2、相對誤差RE和標準化的均方根誤差NRMSE來評價模型模擬的準確性[23],其中,相對誤差和標準化的均方根誤差的計算公式如下:

表1 DNDC模型所需輸入參數
注: DNDC模型為生物地球化學模型,DeNitrification- DeComposition model.

表2 不同土地利用類型土壤基本理化性質
*:同一列中字母不同表示不同土地利用類型間存在顯著差異,相同則表示無顯著差異,顯著性P<0.05

表3 農地和園地管理措施設置
(1)相對誤差法,公式為:
式中,E為相對誤差,Oi為觀測值,Pi為預測值,n為觀測樣本數。
(2)標準化的均方根誤差法

決定系數R2越接近于1,表明實測值與模擬值線性相關性越好。RE絕對值小于5%,表明模擬的結果準確,RE絕對值小于10%,表明模擬的結果可行。NRMSE的值越小,表明模擬值與實測值的擬合度越高。一般情況下,NRMSE值小于25%,即模擬值與實測值一致性非常好,在25%—30%之間表明模擬效果一般,大于30%則表明模擬值與實際值偏差大,模擬效果不理想[24]。
土壤有機碳儲量的計算公式為:

式中,SOCi為第i個樣點單位面積有機碳含量,S為土地利用類型面積,n為觀測樣本數。
土壤氮儲量的計算公式為:

式中,SN為第i個樣點單位面積總氮含量,S為土地利用類型面積,n為觀測樣本數。
研究區位于城鎮向自然生態系統的過渡地帶,土地利用以農地、城鎮建設用地和林地為主。由表4和圖2可以看出,1974至2015年研究區林地和農地面積均在不斷減小,城鎮建設用地面積在不斷增加,園地自2000年左右出現以后面積呈不斷增加的趨勢。研究區農地受地形條件限制,主要分布在山間平原,并且受城鎮化的影響,面積在不斷減小。城鎮建設用地主要沿河流分布在平原地區,地形的限制使得城鎮面積擴張受到嚴重制約,但2015年城鎮建設用地面積也是1974年的2.3倍。園地主要為茶園、果園、苗圃等,2000年左右開始出現,主要由村鎮附近低山丘陵的林地轉換而來。2000年以后,伴隨著園地面積的不斷擴大,林地面積在不斷減小。研究區丘陵山地面積較大且山體較陡,土壤層發育較薄,是林地的主要分布范圍。由于山地地形和土壤條件不適宜農業耕作,因而研究區內林地作為自然生態系統在幾十年中沒有顯著的變化。

表4 研究區1974—2015年土地利用變化/%

圖2 研究區1974—2015年土地利用圖Fig.2 Land use map of study area from 1974 to 2015

圖3 模擬值與實測值的相關性 Fig.3 Correlations between simulated values and observed values
選取研究區38個典型農地樣點用于DNDC模型校正,按DNDC點位模式要求輸入校正數據后對貝母產量(單位:kgC/hm2)進行模擬。數據表明,實測值和模擬值的決定系數為0.93,RE和NRMSE分別為-1.44%和7.69%,表明DNDC模型在本研究區有較好的模擬結果(圖3)。
DNDC模型在計算園地土壤的碳循環機制方面已經得到了較為廣泛的應用[25],并且已經成功應用在與研究區土壤和氣候條件相似的上海地區園地,用于模擬和定量分析農業生產過程中氮素的遷移轉化規律,相關研究表明DNDC模型可用于對本研究區園地的土壤碳、氮模擬。Forest-DNDC模型在我國貢嘎山、川中丘陵區和喀斯特石質山區等環境中成功進行了土壤碳、氮循環過程的模擬[26-28],研究區的氣候、土壤、植被等在相似環境中已有相關模擬研究,因而采用Forest-DNDC模型對本研究區林地土壤碳、氮循環進行模擬。
DNDC模擬表明,1974年至1980年農地單位面積有機碳含量略有增長,1980年至1995年有所波動,2000年以后則逐漸降低(圖4)。園地單位面積有機碳含量呈波動變化,而林地單位面積有機碳含量則持續增加。不同土地利用類型單位面積總氮變化表明,1974年至1985年農地單位面積總氮含量持續增加,1990年有明顯的減少,1995年以后農地單位面積總氮含量逐漸降低(圖4)。園地和林地單位面積總氮含量變化與單位面積有機碳含量變化較為相似。

圖4 不同土地利用類型單位面積有機碳和總氮的變化Fig.4 Changes of soil organic carbon and total nitrogen per unit area in different land use types
從圖5土壤單位面積有機碳含量同降雨和氣溫的相關分析可以看出,除農地外,園地和林地土壤有機碳含量均隨著溫度的升高而增加,不同土地利用類型有機碳含量均隨著降雨量的增加而減少。圖6單位面積土壤總氮含量與氣象因子的相關分析表明,除農地外,土壤總氮含量均隨著溫度的升高而增加,但隨著降雨量的增加并無明顯變化。

圖5 單位面積有機碳含量與年均溫和年降雨量之間的相關性Fig.5 Correlations between soil organic carbon per unit area and annual average temperature and annual precipitation

圖6 單位面積總氮與年均溫和年降雨量之間的相關性Fig.6 Correlations between total nitrogen per unit area and annual average temperature and annual precipitation
對不同土地利用土壤碳、氮儲量的分析表明,流域農地土壤有機碳儲量呈逐漸降低的趨勢,而園地和林地土壤有機碳儲量則不斷增加。圖4表明農地單位面積有機碳含量在模擬期間呈小幅度波動變化,而流域農地土壤總有機碳儲量則呈不斷降低的趨勢(圖7),這主要是由于流域農地面積的不斷減少,使得農地土壤有機碳總儲量逐漸降低,另一方面,氣象條件的變化也在一定程度上影響了農地土壤總有機碳的積累。園地土壤有機碳總儲量的變化則源于園地面積的迅速增加。相比而言,林地單位面積有機碳含量處于不斷增加的過程,再加上林地面積在這一時間內變化并不十分明顯,從而使得林地土壤有機碳總儲量也不斷增加。流域不同土地利用類型土壤總氮儲量的變化趨勢和有機碳較為一致,其中農地土壤總氮的減少和園地土壤總氮儲量的增加均與面積的改變顯著相關,而林地則同樣受單位面積土壤氮含量的影響,呈現不斷增加的趨勢(圖7)。

圖7 不同土地利用有機碳和總氮儲量變化Fig.7 The soil organic carbon and the soil total nitrogen changes in different land use types
DNDC模擬表明,41年間非城鎮用地(農地+園地+林地)面積逐漸減少(圖8),但流域土壤有機碳和總氮儲量整體逐年增加,由1974年的1.38×108t增至2015年的5.49×108t,總氮由1974年的0.25×108t增至2015年的0.36×108t。對于土壤有機碳和總氮來說,由于林地在面積和單位面積儲量上均遠高于農地和園地,使得流域土壤有機碳和總氮儲量的變化趨勢與林地的變化趨勢基本一致。

圖8 流域有機碳、總氮年儲量變化Fig.8 The organic carbon and total nitrogen changes at watershed scale
研究區處于多山丘陵區,由于受地形限制,林地面積占流域面積約80%,農地主要分布山間平原地區,園地零星分布于農地中,城鎮建設用地基本沿河分布。1974—2015年土地利用變化主要體現在城鎮建設用地侵占農地而持續增加,園地主要靠開墾地勢較低的林地。這與已有的長江三角洲地區城鎮化的時空動態格局變化相一致[29]。
對于DNDC模型和Forest-DNDC模型,在點位模擬中,默認單元內的土壤屬性是均一的,而實際上土壤屬性在空間上具有一定的差異性,在進行儲量計算時會導致模擬結果與實際情況有一定的偏差[30],這也是模型模擬普遍存在的問題。另一方面,為了維持土壤肥力和作物產量,施肥和農地管理措施在不斷發展和更新,土壤有機碳的模擬結果受土壤有機碳初始含量的影響,而土壤氮平衡過程除了受有機碳的影響,也受到氮肥投入的影響[31]。本研究重點探討了城郊流域土地利用變化對土壤碳氮儲量的影響,未考慮多年施肥和管理措施的演變對土壤碳氮儲量的影響,在一定程度上影響了模擬結果,在今后的研究中需加以考慮,以達到更為精準的模擬。
影響土壤有機碳、總氮的因素較為復雜,如溫度、濕度、水土流失和生物量等。已有研究表明,土壤有機碳含量隨溫度和降雨量的變化而變化,其影響因素因區域環境不同而存在差異[28]。本文中農地單位面積有機碳含量、總氮含量與溫度和降雨均呈負相關,這與周文強等[32]和曹宏杰等[33]的結論相同。氣溫升高和降雨增加能促使土壤微生物的呼吸作用增加,加速土壤有機碳的消耗,從而降低土壤的固碳能力[34],并且長期翻耕能夠加速土壤微生物分解活動,降低土壤固碳能力[35],使得農地土壤有機碳含量隨著溫度和降雨的增加而降低。而且降雨影響氮淋溶和流失,與溫度均會影響NOx和NH3的揮發速率以及揮發量[36],并影響土壤硝化-反硝化的反應速率[37],從而影響土壤氮含量的變化。另一方面,施肥對農地NOx和NH3的排放具有明顯作用,且肥料種類不同,氣體排放量不同[8-39]。研究區農地每年三耕兩耱,按照植物生長周期定期翻耕和施肥,土質疏松孔隙大,降雨和溫度會增加農地土壤氮的流失以及提高NOx和NH3的揮發量,從而降低了土壤的固氮能力。
林地單位面積有機碳含量、總氮含量與溫度呈顯著正相關,與降雨量呈負相關,這與畢珍等[40]在四川盆地森林的研究結果相似,但與其溫度的相關性結果有差異。森林喬木層郁閉度高,生長季結束后在地面形成較厚的枯枝落葉層,以及更茂密的地下根系,積累了大量的有機質。并且溫度升高除了促使COx和CH4的排放量增加,同時也會加速枯枝落葉和地下根系的分解并補給土壤碳庫,并隨著樹齡的增長,也會促進土壤有機碳的累積[41],從而使得園地和林地土壤有機碳隨溫度變化與農地有所不同。同溫度對林地有機碳含量的影響一樣,溫度升高也會加速枯枝落葉和地下根系的分解并補給土壤氮庫,增加土壤氮含量。周濤等[42]的研究也表明,土壤有機碳含量與溫度和降雨的關系在不同的氣候帶有顯著的差異,當溫度為10—20℃時,土壤有機碳與溫度呈顯著正相關(P<0.01),與降水呈正相關,并且降水與溫度之間存在顯著的正相關,這種強的相關性可能導致土壤有機碳含量與某一因子之間的關系受到第三個因子的影響。本研究中降雨與溫度之間存在正相關,且溫度對有機碳的影響也高于降雨的影響,土壤有機碳與降水量呈負相關的關系,有可能是溫度效應較強的結果。研究區林地單位面積總氮含量變化與溫度呈顯著正相關,與降雨呈負相關,但不顯著。這是由于研究區屬亞熱帶季風氣候,光照充足,雨量豐沛,溫度升高和降雨增加促進NOx和NH3排放量和氮淋溶量,造成土壤氮的損失[43]。同時,這一氣候條件也能促使枯枝落葉加速分解,并補給土壤氮庫[44]。園地單位面積有機碳含量、總氮含量與溫度和降雨的相關性與林地的規律類似,但由于園地出現時間較短,數據量較少,并不具代表性。
土地利用變化會改變土壤中的碳、氮含量,其中,林地轉變為農地、農地轉變為城鎮用地等均會大幅度降低土壤有機碳和總氮含量[45]。農地因總面積和單位面積有機碳含量逐年降低,成為流域中土壤有機碳儲量流失的主要區域,說明農地面積的變動和人為擾動都將導致土壤有機碳儲量發生變化[46-47]。而研究區林地面積約占流域總面積的80%,并且其單位面積有機碳含量處于逐漸增加的狀態,從而影響了流域土壤有機碳總儲量及其動態,這說明固碳區域在土地利用結構上的優勢能決定流域整體有機碳儲量的增減[48]??偟獌α康淖兓f明總氮的輸入與輸出是逐年動態變化的,并且變化幅度較大,與已有研究結果相似[49],這些均表明土地利用的變化驅動了土壤有機碳、總氮儲量的變化。
研究表明,在1974—2015年間,典型城郊流域的林地和農地面積均在不斷縮減,而城鎮建設用地面積不斷增加,園地自2000年左右出現后面積不斷擴大。應用DNDC模型和Forest-DNDC模型對農地、園地和林地土壤碳氮含量的模擬發現,農地單位面積有機碳、總氮含量逐年降低,而林地則不斷增加,園地呈波動變化,不同土地利用類型的單位面積有機碳、總氮含量對溫度和降雨的變化有不同程度的響應。隨著城鎮化對土地利用的改變,使得流域農地土壤有機碳總儲量降低,園地土壤有機碳總儲量迅速增長,林地土壤有機碳總儲量不斷增加,土壤總氮含量與有機碳含量變化趨勢相似,而流域土壤總有機碳和總氮儲量則受優勢土地利用類型的影響較大。