高晨曦,盧秋萍,歐年青,胡清萍,林 雪,鮑玲鑫

(福建農(nóng)林大學計算機與信息學院/生態(tài)與資源統(tǒng)計福建省高校重點實驗室 福州 350002)

CO2等溫室氣體所產(chǎn)生的溫室效應加速了全球氣候變化的進程,世界各國已通過全球協(xié)議達成減少溫室氣體排放的共識,中國也積極參與在全球碳減排行動之中,并于2020年宣布“CO2排放力爭于2030年前達到峰值,努力爭取2060年前實現(xiàn)碳中和”。農(nóng)業(yè)生產(chǎn)作為碳排放的重要來源[1],約占全球溫室氣體排放總量的14%[2],僅次于電熱生產(chǎn)。中國作為農(nóng)業(yè)大國,農(nóng)業(yè)碳排放占總碳排放量的17%左右[3]。因此,農(nóng)業(yè)碳減排是中國實現(xiàn)碳達峰、碳中和目標的重要方面。

目前,農(nóng)業(yè)碳排放的研究主要集中在碳排放測算[4-7]、時空演變[8-9]、脫鉤效應[10-12]和碳減排策略[13-14]等方面,農(nóng)業(yè)碳排放影響因素及預測分析有待深究[15]。在農(nóng)業(yè)碳排放影響因素分析方法中,可拓展的隨機性的環(huán)境影響評估模型(Stochastic Impacts by Regression on Population,Affluence,and Technology,STIRPAT)可以充分考慮社會、經(jīng)濟和技術方面的影響因素[16],但忽略了區(qū)域城鎮(zhèn)化進程對農(nóng)業(yè)發(fā)展的影響。在農(nóng)業(yè)碳排放預測分析方法中,目前普遍采用傳統(tǒng)預測方法,如:環(huán)境庫茲涅茨曲線、IPAT 方程、STIRPAT 模型、投入產(chǎn)出模型、灰色預測模型和LEAP 模型等[17],這些傳統(tǒng)預測方法往往存在模型計算復雜或過擬合等缺陷,且針對復雜的非線性農(nóng)業(yè)碳排放系統(tǒng),由于缺乏學習和測試過程,易導致非線性關系模糊、泛化能力不足等諸多弊端。支持向量回歸(Support Vector Regression,SVR)作為一種基于統(tǒng)計學習理論的機器學習方法,計算的復雜度只取決于支持向量的數(shù)目,不敏感損失函數(shù)的引入與核技巧的應用也使得其能夠避免過擬合和較好地解決非線性問題,適用于農(nóng)業(yè)碳排放預測研究。

河南省作為農(nóng)業(yè)大省,2020年農(nóng)業(yè)生產(chǎn)總值占全國的7.23%,位居全國第二。本文以河南省為例,采用排放因子法測算2001—2020年河南省農(nóng)業(yè)碳放量。針對目前農(nóng)業(yè)碳排放影響因素及預測分析模型存在的不足之處,構建擴展STIRPAT 模型,對農(nóng)業(yè)碳排放的各種影響因素進行定性與定量分析。建立基于徑向基函數(shù)(Radial Basis Function,RBF)核ε-SVR 的農(nóng)業(yè)碳排放預測模型,預測不同情景下2021—2025年河南省農(nóng)業(yè)碳排放量及趨勢,以期為河南省農(nóng)業(yè)碳減排及相關政策的制定提供量化參考和科學依據(jù),同時也為區(qū)域農(nóng)業(yè)碳排放研究提供新的方法借鑒。

1 研究方法與數(shù)據(jù)來源

1.1 農(nóng)業(yè)碳排放測算

采用聯(lián)合國政府間氣候變化專門委員會(Intergovernmental Panel on Climate Change,IPCC)指南的排放因子法測算農(nóng)業(yè)碳排放量。將農(nóng)業(yè)溫室氣體排放源劃分為4 類:水稻種植、土地利用、動物腸道發(fā)酵和糞便管理,并以此建立農(nóng)業(yè)碳排放測算框架。確定各溫室氣體排放系數(shù)及對應的百年尺度CO2變暖潛力的倍數(shù)。最終測算出農(nóng)業(yè)CO2排放總量,測算公式如下:

式中:E為農(nóng)業(yè)CO2排放總量;g i為第i類農(nóng)業(yè)碳排放源的量;e ij為第i類農(nóng)業(yè)碳排放源第j類溫室氣體的排放系數(shù);GWPj為第j類溫室氣體對應的百年尺度CO2變暖潛力的倍數(shù)。根據(jù)IPCC 第5 次評估報告,溫室氣體CO2、CH4和N2O 的百年尺度CO2變暖潛力分別為1 倍、28 倍和265 倍。

農(nóng)業(yè)溫室氣體排放源與排放系數(shù)如表1所示。

表1 農(nóng)業(yè)溫室氣體排放源與排放系數(shù)Table 1 Agricultural greenhouse gas emission sources and emission coefficients

1.2 農(nóng)業(yè)碳排放影響因素分析模型

STIRPAT 模型是York 等[19]在Ehrlich 等[20]提出的IPAT 模型的基礎上改進而來的用于分析人文因素對環(huán)境影響的隨機模型。它將環(huán)境影響因素分解為人口規(guī)模、富裕程度和技術水平3 類,不僅克服了IPAT 模型在假設檢驗方面的局限性[21],還允許通過增減或因素分解進行擴展改進[22],且便于定性和定量地分析人文指標對碳排放的影響程度。STIRPAT 模型的基本形式為:

式中:I 、 P、 A、 T分別代表環(huán)境影響、人口規(guī)模、富裕程度和技術水平;a為模型系數(shù);b、 c、 d分別為P、 A、 T的彈性指數(shù);e為隨機誤差。

考慮到區(qū)域城鎮(zhèn)化進程對農(nóng)業(yè)發(fā)展的影響,本文在STIRPAT 模型的基礎上引入城鎮(zhèn)化率指標對其進行擴展,建立以農(nóng)村人口、農(nóng)作物播種面積、大型畜牧數(shù)、人均農(nóng)業(yè)GDP、農(nóng)村人均可支配收入、農(nóng)業(yè)機械化水平和城鎮(zhèn)化率為農(nóng)業(yè)碳排放影響因素的擴展STIRPAT 模型。模型表達為:

式(3)兩邊同時取對數(shù),將非線性模型轉化成線性模型,以便進行回歸分析確定模型參數(shù)。對數(shù)變換后模型表達為:

式中:E為農(nóng)業(yè)CO2排放總量(萬t);a為常數(shù)項;R為農(nóng)村人口(萬人);S為農(nóng)作物播種面積(萬hm2);L為大型畜牧數(shù)(萬頭),代替原模型中的人口規(guī)模項;A為人均農(nóng)業(yè)GDP(元);I為農(nóng)村人均可支配收入(元),代替原模型中的富裕程度項;M為農(nóng)業(yè)機械化水平,以農(nóng)業(yè)機械總動力(萬kWh)表示,代替原模型中的技術水平項;U為城鎮(zhèn)化率(%),反映城鎮(zhèn)化進程;β1、β2、β3、β4、β5、β6、β7為模型的彈性系數(shù),根據(jù)彈性系數(shù)概念,當R、S、L、A、I、M、U每發(fā)生1%的變動時,將分別引起E發(fā)生β1%、β2%、β3%、β4%、β5%、β6%、β7%的變動。

1.3 農(nóng)業(yè)碳排放預測模型

SVR 是Drucker 等[23]在支持向量機的基礎上擴展而來的一種能夠進行回歸計算的機器學習方法。由于SVR 適合小樣本的學習,理論上可獲得全局最優(yōu)解,且能夠通過核技巧較好地解決非線性回歸問題。不敏感損失函數(shù)ε的引入允許其輸出值f(x)與真實值y之間存在至多ε的偏差,在一定程度上能夠避免過擬合,提高泛化能力。計算的復雜度也只取決于支持向量的數(shù)目,而與樣本空間維數(shù)無關,在高維空間中具有一定優(yōu)勢。本文所研究的系小樣本非線性問題,為避免過擬合和提高泛化能力,采用如下ε-SVR 模型:

輸入訓練集數(shù)據(jù):

式中:x i為輸入向量;y i為輸出值,i=1,2,···,n。

設R m上的一個線性表達式為:

式中:ω ∈R m為權值向量;b∈R為偏置。

f(x)的求解問題可轉化為式(7)、式(8)的最優(yōu)化問題:

1.4 數(shù)據(jù)來源及處理

本研究中河南省農(nóng)業(yè)碳排放原始數(shù)據(jù)根據(jù)2002—2021年《中國統(tǒng)計年鑒》和《河南省統(tǒng)計年鑒》整理得出。其中:化肥施用量采用折純量表示,農(nóng)藥、農(nóng)膜、農(nóng)用柴油均以實際使用量表示,翻耕數(shù)據(jù)以當年實際農(nóng)作物播種面積為準,大型畜牧數(shù)均以年末存欄數(shù)量為準,人均農(nóng)業(yè)GDP 為種植業(yè)和畜牧業(yè)生產(chǎn)總值與鄉(xiāng)村人口的比值,農(nóng)業(yè)CO2排放強度為農(nóng)業(yè)CO2排放總量與種植業(yè)和畜牧業(yè)生產(chǎn)總值的比值,農(nóng)業(yè)機械化水平以農(nóng)業(yè)機械總動力表示。

2 結果與分析

2.1 農(nóng)業(yè)碳排放結構及時序特征分析

2001—2020年河南省農(nóng)業(yè)碳排放結構表明(圖1):土地利用對農(nóng)業(yè)CO2排放總量的平均貢獻率最大,為48.31%;其次是動物的腸道發(fā)酵、糞便管理和水稻種植,分別貢獻了31.28%、16.05%和4.37%。水稻種植和土地利用的農(nóng)業(yè)CO2排放量(貢獻率)總體呈上升趨勢,分別從2001年的275.66 萬t(3.02%)、3443.34 萬t(37.75%)上升至2020年的408.97 萬t(5.62%)、4357.12 萬t(59.88%)。 動物的腸道發(fā)酵、糞便管理的農(nóng)業(yè)CO2排放量(貢獻率)總體呈下降趨勢,分別從2001年的3857.85 萬t(42.30%)、1543.79萬t(16.93%)下降至2020年的1444.25 萬t(19.85%)、1066.23 萬t(14.65%)。 總體來看,河南省農(nóng)業(yè)CO2主要排放源已經(jīng)從畜牧業(yè)的腸道發(fā)酵和糞便管理轉移至種植業(yè)的土地利用和水稻種植;農(nóng)業(yè)CO2排放總量主要呈現(xiàn)“上升-下降”的趨勢,年均增長率為—1.18%,當前全省農(nóng)業(yè)碳排放總體已呈降低趨勢,特別是在2006年和2017年有大幅下降,分別較上年下降9.40%和16.67%,且農(nóng)業(yè)CO2排放量于2005年達峰值,為10 256.69 萬t,而2020年已降至 7276.56 萬t,較峰值降低29.06%。

圖1 2001—2020年河南省農(nóng)業(yè)CO2 排放結構Fig.1 Structure of agricultural CO2 emissions in Henan Province from 2001 to 2020

農(nóng)業(yè)碳排放強度是指每單位農(nóng)業(yè)生產(chǎn)總值的增長所帶來的碳排放量[24],能有效反映碳排放與農(nóng)業(yè)經(jīng)濟增長之間的聯(lián)系。2001—2020年河南省農(nóng)業(yè)CO2排放強度主要呈現(xiàn)“下降-上升-下降”的趨勢(圖2),年均增長率為—8.69%,2001—2020年河南省農(nóng)業(yè)CO2排放強度呈現(xiàn)先降低后升高的趨勢,且于2003年達到峰值4.86 t·萬元—1,此后逐年下降,2004年的下降幅度最大,相比上年排放強度降低24.27%,2020年已下降至0.80 t·萬元—1,較峰值降幅為83.54%。這些結果表明近年來河南省農(nóng)業(yè)產(chǎn)業(yè)轉型優(yōu)化對農(nóng)業(yè)碳減排效果顯著,農(nóng)業(yè)碳排放與農(nóng)業(yè)經(jīng)濟增長已逐漸脫鉤。

圖2 2001—2020年河南省農(nóng)業(yè)CO2 排放強度變化Fig.2 Changes of agricultural CO2 emission intensity in Henan Province from 2001 to 2020

2.2 農(nóng)業(yè)碳排放影響因素分析

對河南省農(nóng)業(yè)碳排放量及影響因素數(shù)據(jù)做自然對數(shù)處理后進行普通最小二乘回歸分析,發(fā)現(xiàn)模型存在嚴重的多重共線性。為確保模型的穩(wěn)定性和精確度,選擇主成分分析對變量進行降維處理,將線性相關的農(nóng)業(yè)碳排放各影響因素通過正交變換合并成少數(shù)幾個線性無關的主成分,達到消除多重共線性的目的。

將自然對數(shù)處理后的數(shù)據(jù)進行標準化處理以消除量綱影響,分別以ZlnR、ZlnS、ZlnL、ZlnA、ZlnI、ZlnM、ZlnU表示標準化后的變量。運用SPSS 23.0 對處理后的變量進行主成分分析,提取兩個主成分F1、F2,可解釋原始變量的98.347%(表2)。主成分方程系數(shù)其中:r ij為成分矩陣各分量,λi為主成分特征值(表3)。F1、F2為:

表2 農(nóng)業(yè)碳排放影響因素主成分分析的總方差解釋Table 2 Total variance interpretation of principal component analysis for influencing factors of agricultural carbon emissions

表3 農(nóng)業(yè)碳排放影響因素的成分矩陣和主成分方程系數(shù)Table 3 Component matrix and principal component equation coefficients of influencing factors of agricultural carbon emissions

利用樣本數(shù)據(jù)計算F1、F2值, 并與ZlnE做最小二乘回歸, 可得:

2.3.1 支持向量回歸結果分析

運用MATLAB,按照基于RBF 核ε-SVR 的農(nóng)業(yè)碳排放預測模型的步驟,對樣本數(shù)據(jù)進行歸一化處理,從歸一化后的樣本數(shù)據(jù)中隨機選取10年數(shù)據(jù)組成訓練集,剩余10年數(shù)據(jù)組成測試集,設定ε的初始值為0.01,采用網(wǎng)格搜索和5 折交叉驗證對懲罰參數(shù)C和核函數(shù)參數(shù) γ進行尋優(yōu),C和 γ的取值范圍均為[2—10,210],網(wǎng)格寬度為0.2。系統(tǒng)將在訓練集中尋找出訓練誤差最小的一組超參數(shù)組合。參數(shù)選擇效果見圖3,最佳模型超參數(shù)(C,γ)的值為(1024,0.000 976 562 5),二乘回歸,可得:

圖3 參數(shù)選擇效果圖(3D)Fig.3 The effect drawing of the parameter selection(3D)

可決系數(shù)R2為0.940,調整后的R2為0.933,表明模型具有較高的擬合度。Sig.F=0.000＜0.05,表明在P＜0.05 的顯著性水平下模型線性關系顯著。Sig.t=0.000＜0.05,表明在P＜0.05 的顯著性水平下模型回歸系數(shù)顯著(表4)。將式(17)、(18)代入式(19)中可得:

表4 農(nóng)業(yè)碳排放影響因素的主成分回歸系數(shù)和顯著性檢驗Table 4 Coefficient and significance test of principal component regression of influencing factors of agricultural carbon emissions

將式(20)進行反標準化處理,可得:

由式(21)得到最終的河南省農(nóng)業(yè)碳排放擴展STIRPAT 模型:

模型表明:鄉(xiāng)村人口(R)、農(nóng)作物播種面積(S)、大型畜牧數(shù)(L)、人均農(nóng)業(yè)GDP(A)、農(nóng)村人均可支配收入(I)、農(nóng)業(yè)機械化水平(M)、城鎮(zhèn)化率(U)每發(fā)生1%的變動,將分別引起河南省農(nóng)業(yè)碳排放量0.162%、 0.175%、 0.130%、 —0.018%、 —0.029%、0.120、—0.071%的變動;影響河南省農(nóng)業(yè)碳排放量的7 個因素中,農(nóng)作物播種面積的促進作用最大,其次是鄉(xiāng)村人口和大型畜牧數(shù),農(nóng)業(yè)機械化水平的促進作用最小,城鎮(zhèn)化率的抑制作用最強,其次是農(nóng)村人均可支配收入,人均農(nóng)業(yè)GDP 的抑制作用最小。

2.3 農(nóng)業(yè)碳排放預測分析

交叉驗證的MSE 為0.0395。將訓練集和測試集的擬合值與實測值進行線性回歸,其中訓練集線性回歸的R2為0.96,MSE 為0.01,表明模型具有較高的擬合度。測試集線性回歸的R2為0.97,MSE 為0.01,表明模型具有較強的泛化能力。反歸一化后訓練集和測試集的擬合值和實測值對比見圖4。

圖4 訓練集(a)和測試集(b)農(nóng)業(yè)CO2 排放擬合值與實測值結果對比Fig.4 Comparing the fitting values of training set(a)and test set(b)with the measured values of agricultural CO2 emissions

同時,與2.2 節(jié)中擴展STIRPAT 模型的R2為0.94,MSE 為0.067 相比,基于RBF 核ε-SVR 的農(nóng)業(yè)碳排放預測模型的測試集和訓練集均具有更高的擬合度和更小的均方誤差。可見其在區(qū)域農(nóng)業(yè)碳排放預測研究方面有更好的學習和推廣能力。

2.3.2 影響因素預測分析

根據(jù)《河南省國民經(jīng)濟和社會發(fā)展第十四個五年規(guī)劃和二〇三五年遠景目標綱要》(以下簡稱規(guī)劃目標),并結合河南省農(nóng)業(yè)發(fā)展情況,設置基準情景、低碳情景Ⅰ和低碳情景Ⅱ 3 種情景模式,以期能夠全面分析未來河南省農(nóng)業(yè)碳排放趨勢。其中:農(nóng)村人均可支配收入基于規(guī)劃目標中人均可支配收入年均增長率6.0%設定;人均農(nóng)業(yè)GDP 基于規(guī)劃目標地區(qū)生產(chǎn)總值年均增長率設定,同為6.0%;城鎮(zhèn)化率根據(jù)規(guī)劃目標要在2025年達60%以上,2020年河南省城鎮(zhèn)化率為55.43%,年均最低增長率約為2.5%;鄉(xiāng)村人口同城鎮(zhèn)化率具有相反的增長率;農(nóng)作物播種面積、大型畜牧數(shù)和農(nóng)業(yè)機械化水平3 項影響因素沒有具體的發(fā)展目標,建立自回歸移動平均模型

(Autoregressive Integrated Moving Average Model,ARIMA)進行預測。對3 項影響因素數(shù)據(jù)進行差分平穩(wěn),確定其為非白噪聲序列后,根據(jù)自相關系數(shù)和偏自相關系數(shù)并結合赤池信息量和貝葉斯信息量選定各模型階數(shù)。各影響因素預測結果歸納為表5。

表5 農(nóng)業(yè)碳排放影響因素預測分析Table 5 Prediction and analysis of influencing factors of agricultural carbon emission

2.3.3 不同情景下農(nóng)業(yè)碳排放預測分析

分別對基準情景和平衡低碳情景下的河南省農(nóng)業(yè)CO2排放量進行預測,得到各情景下農(nóng)業(yè)CO2排放量預測值(圖5)。

圖5 不同情景模式下河南省農(nóng)業(yè)CO2 排放量預測Fig.5 Prediction of agricultural CO2 emissions in Henan Province under different scenarios

結果表明:在基準情景下,2021—2025年河南省農(nóng)業(yè)CO2排放量仍呈持續(xù)下降的變化趨勢,2025年預測值為6483.80 萬t,與2020年7276.56 萬t 相比減少10.89%;在低碳情景Ⅰ下,呈現(xiàn)出較基準情景更快的下降速率,2025年預測值為6369.19 萬t,較2020年減少12.47%;在低碳情景Ⅱ下,河南省農(nóng)業(yè)CO2排放量的下降速率最高,2025年預測值為6338.32 萬t,較2020年減少12.89%。綜上所述,2021—2025年河南省農(nóng)業(yè)CO2排放量仍呈持續(xù)下降的變化趨勢,且相比基準情景,低碳情景下河南省具有更大的碳減排潛力,在實現(xiàn)農(nóng)業(yè)碳減排的同時,也能兼顧農(nóng)村經(jīng)濟發(fā)展和城鎮(zhèn)化進程。

3 討論和結論

3.1 討論

本文利用機器學習方法對河南省農(nóng)業(yè)碳排放問題進行了系統(tǒng)研究,建立區(qū)域農(nóng)業(yè)CO2排放測算指標體系,定性與定量分析農(nóng)業(yè)碳排放的影響因素與未來排放趨勢。在研究結果上,河南省農(nóng)業(yè)已實現(xiàn)碳達峰,農(nóng)業(yè)碳排放的主要來源已從畜牧業(yè)轉移至種植業(yè),其中農(nóng)作物播種面積為促進農(nóng)業(yè)碳排放的首要因素,這為未來全省實現(xiàn)農(nóng)業(yè)碳中和指明了政策方向。伴隨著城鎮(zhèn)化率和農(nóng)村經(jīng)濟的進一步提升,未來全省農(nóng)業(yè)碳排放將持續(xù)下降,助力國家“雙碳”目標的實現(xiàn)。在研究方法上,與現(xiàn)有的灰色預測模型、LEAP 模型等農(nóng)業(yè)碳排放研究方法相比,本文采用機器學習方法對農(nóng)業(yè)碳排放進行影響因素及預測分析,很大程度上避免了傳統(tǒng)模型的過擬合、非線性關系模糊和泛化能力不足等缺陷,且通過訓練和測試所得到的農(nóng)業(yè)碳排放預測模型也具有良好的學習和推廣能力,亦可應用于其他區(qū)域的農(nóng)業(yè)碳排放預測研究。本文與李波等[4]、田成詩等[7]、王劼等[12]對農(nóng)業(yè)碳排放的測算與影響因素研究相比,對區(qū)域農(nóng)業(yè)CO2排放的測算指標選取更加細致,從農(nóng)業(yè)產(chǎn)業(yè)、農(nóng)村經(jīng)濟、農(nóng)業(yè)技術和城鎮(zhèn)化進程4 方面更加全面地探究了農(nóng)業(yè)碳排放的影響因素。與黎孔清等[15]、邱子健等[21]對農(nóng)業(yè)碳排放的預測研究相比,所選用的預測模型更加符合復雜的非線性農(nóng)業(yè)碳排放系統(tǒng),且模型的泛化能力更強。本文為河南省農(nóng)業(yè)低碳化發(fā)展提供了量化參考和科學依據(jù),同時也為區(qū)域農(nóng)業(yè)碳排放研究提供了新的方法借鑒。

在今后的農(nóng)業(yè)碳排放研究中將考慮更大區(qū)域尺度的農(nóng)業(yè)碳排放測算及空間差異分析。同時,結合深度學習對農(nóng)業(yè)碳排放進行預測分析也是值得深入研究的方向。

3.2 結論

1)2001—2020年河南省農(nóng)業(yè)碳排放結構及時序特征分析表明:河南省農(nóng)業(yè)CO2主要排放源已經(jīng)從畜牧業(yè)的腸道發(fā)酵和糞便管理轉移至種植業(yè)的土地利用和水稻種植;排放量及排放強度總體均呈下降趨勢,且分別于2005年與2003年達峰,農(nóng)業(yè)碳減排效果顯著,農(nóng)業(yè)碳排放與農(nóng)業(yè)經(jīng)濟增長已逐漸脫鉤。

2)通過擴展的STIRPAT 模型對河南省農(nóng)業(yè)碳排放影響因素分析表明:河南省農(nóng)業(yè)碳排放的治理重點是農(nóng)作物的土地利用及大型畜牧的糞便管理,促進農(nóng)業(yè)碳減排的著重點是穩(wěn)妥推進城鎮(zhèn)化及農(nóng)村經(jīng)濟發(fā)展。

3)基于RBF 核ε-SVR 的農(nóng)業(yè)碳排放預測模型對不同情景下2021—2025年河南省農(nóng)業(yè)碳排放量預測分析表明:2021—2025年河南省農(nóng)業(yè)CO2排放量仍將呈現(xiàn)持續(xù)下降的變化趨勢,且較基準情景,低碳情景具有更大的碳減排潛力,并能實現(xiàn)農(nóng)村經(jīng)濟、城鎮(zhèn)化進程與低碳農(nóng)業(yè)兼顧的高效發(fā)展,有助于加速全省“雙碳”目標的實現(xiàn)。