歐美楊人工林生態系統凈碳交換對環境因子響應的時滯*

馮鑫煒 張志強 許 行 律 江 張海泉 孟祥雪

(1.北京林業大學水土保持學院 北京 100083;2.北京市共青林場 北京 101300)

大氣中CO2濃度持續升高引起的全球變暖和極端氣候頻率增加已成為事實(IPCC,2013;Kaushaletal.,2014)。森林生態系統作為陸地生態系統中最主要的植物碳庫(戴巍等,2017)，其植被對減緩氣候變化進程發揮著不可替代的作用。隨大范圍造林，我國人工林面積不斷增加，已占森林總面積的33.17%(國家林業局,2014)，人工林生態系統對碳收支的潛在影響不可忽略。分析人工林凈生態系統碳交換(NEE)的變化特征及其與環境因子間的相互關系，對準確評估我國人工林碳匯功能具有重要意義(唐祥等,2013;吳亞叢等,2013)，可為全球范圍內森林生態系統碳循環研究提供重要支撐(鄒文濤等,2017)。

NEE受到土壤、植被、水、氣、光、熱等多種因素的共同作用，存在復雜的響應機制(Lasslopetal.,2010)。大量研究證明，輻射是日尺度下驅動NEE變化的主要因子(吳志祥等,2014;譚麗萍等,2015;馬小紅等,2017)。也有研究表明，土壤含水量是NEE日變化的主導因子(徐勇峰等,2018)。然而NEE對環境因子的響應存在時滯現象(Pinginthaetal.,2010;Jiaetal.,2018)，會影響判斷環境因子對NEE的作用。定量研究時滯現象有助于分析氣候變化對陸地生態系統碳循環過程的影響機制、對提高陸地生態系統碳循環模型的參數化和驗證精度(Dietzeetal.,2011;Jiaetal.,2018)具有十分重要的意義，但目前相關研究只停留在定性或定量描述時滯現象(Pinginthaetal.,2010;高翔,2013;Ouyangetal.,2014;Jiaetal.,2018)，很少定量分析忽略時滯現象對判斷NEE與環境影響因子間的相關性以及NEE數據擬合模型精度的影響。

本研究以北京市順義區歐美楊(Populus×euramericana)人工林為研究對象，揭示在日尺度下NEE與環境因子間的時滯現象，并分析消除時滯現象前后NEE與環境因子間的定量關系以及NEE擬合精度的變化情況，以期為提高估算生態系統碳源/匯量的準確性提供依據。

1 研究區概況

研究區位于北京市順義區共青林場(116°42′41″E，40°06′27″N)，海拔29 m，地下水位約2 m。氣候類型屬于暖溫帶半濕潤大陸性季風氣候，年均氣溫11.5 ℃。1月平均氣溫4.9 ℃，最低氣溫-19.1 ℃；7月平均氣溫25.7 ℃，最高氣溫40.5 ℃。年日照2 750 h，全年無霜期195天左右。年均空氣相對濕度50%，年均降水量576 mm，為華北地區降水量較均衡的地區之一。全年75%的降水集中在夏季，每年生長季為4—10月。該地區屬于海河水系潮白河沖積扇下段，地勢平坦，土壤以砂土、亞砂土為主，具有高滲透性和低持水能力。研究區為歐美楊人工林，種植于1996年，2014年胸徑為(25.7 ± 1.6) cm，樹高為(17.5 ± 1.6) m，栽植密度為4 m × 3 m，葉面積指數為2.15 ～ 3.41，林下灌木草本較少，主要有紅瑞木(Swidaalba)和珍珠梅(Sorbariasorbifolia)等。

2 研究方法

2.1 環境因子測定與計算方法

2014年生長季(5—10月)環境因子由微氣象梯度觀測系統直接測定或間接計算所得。其中，空氣溫濕度傳感器(HC253,Campbell Scientific Inc,USA)測量空氣溫度(Ta)和空氣濕度(RH)，觀測高度分別為0.5、1.5、5、15和25 m；土壤溫度傳感器(TCAV107,Campbell Scientific Inc,USA)測量土壤溫度(Ts)，采樣深度分別為5、10和25 cm；光量子傳感器(LI-190SB,LI-COR Inc)測量光合有效輻射(PAR)，觀測高度為30 m；土壤水分觀測儀TDR(CS616,Campbell Scientific Inc,USA)測量土壤體積含水量(VWC)，采樣深度分別為5、25、50、100、150和200 cm。所有數據使用CR1000(Campbell Scientific Inc,USA)采集。飽和水汽壓差(VPD)定義為：

(1)

式中：Ta和RH分別為15 m高處的空氣溫度和空氣濕度，e為自然常數。

本研究土壤溫度采用5 cm深處傳感器數據(T5)，土壤體積含水量采用25 cm深處傳感器數據(VWC25)。

2.2 凈生態系統碳交換計算方法與數據預處理

2014年生長季(4—10月)凈生態系統碳交換(NEE)由渦度相關觀測系統測得的湍流原始數據計算所得。渦度相關觀測系統由開路式H2O/CO2紅外氣體分析儀(EC150,Campbell Scientific Inc.,Logan,Utah,USA)和三維超聲風速儀(CSAT3,Campbell Scientific Inc.,USA)組成，觀測高度25 m，所有10 Hz湍流原始數據采用CR5000(Campbell Scientific Inc,USA)采集。利用Eddypro 6.2.0軟件對所采集的10 Hz湍流原始數據進行統計測試與再處理，處理過程包括去除異常值(約占6%)、坐標轉換(平面擬合法)、WPL校正和通量數據質量分析。將存在10%以上野點或超出臨界值的數據、湍流穩態測試質量等級為8和9的數據以及超出通量貢獻區的數據剔除，最終獲得NEE的30 min原始數據，其中缺失數據占全年數據總量的22%(Xuetal.,2018)。不同因子缺失數據的時間段有可能不同步，故在生長季選取的201天中至少有165天數據有效。根據PAR將NEE和各環境因子數據分為白天(PAR >4 μmol·m-2s-1)和夜間(PAR≤ 4 μmol·m-2s-1)兩部分(Kangetal.,2015)，在以下分析中除描述NEE與環境因子動態變化及簡單線性回歸分析使用全天30 min數據，其余分析分別使用白天和夜間30 min數據，且所有分析均使用未插補的數據。

2.3 回歸分析

通過建立因變量與自變量間相關關系的回歸曲線可反映2個變量間的相關性，以回歸方程是否通過F檢驗和決定系數R2來判斷是否存在相關性及相關程度。利用SPSS 22.0軟件分別對NEE和環境因子進行簡單線性回歸分析，從而初步判斷所選擇的環境因子是否對NEE有影響，并對消除時滯前后NEE與環境因子進行曲線估計得到最優回歸模型。

2.4 錯位平移

按照觀測時間序列順序，以NEE數據為基準，將對應時刻的環境因子數據提前kh(k>0)或滯后-kh(k<0)平移，k=0時不做平移。其中，k為時滯時間，取0，±0.5，±1，±1.5，…，正值表示NEE先于環境因子變化，負值則相反。

2.5 時滯分析方法

小波互相關分析法可以對非平穩的時間序列進行時滯分析，通過分別計算NEE與各環境因子在不同時滯下小波互相關系數的大小來判斷時滯關系(桑燕芳等,2010)。對于時間序列xn(t)和yn(t)(n=1，2，3，…，N)，兩者之間的小波互相關系數定義為：

(2)

(3)

(4)

(5)

(6)

(7)

(8)

式中：c為無量綱時間；ω0為無綱量時間序列頻率，一般取ω0=6(Grinstedetal.，2004)；i為虛數單位。

本研究使用RStudio 1.1.447軟件，對公式(2) ～ (7)編程，將NEE和環境因子的半小時原始數據分別作為時間序列xn(t)和yn(t)輸入程序中，計算結果即為NEE與各環境因子的小波互相關系數數據。由于時滯分析使用白天和夜間2組數據，任意相鄰2天的數據并不連續，在時滯時間k下需要分別計算每天的小波互相關系數值，即最終對于任意時滯k將獲得至少165個小波互相關系數數據。

2.6 偏相關分析

分析2個指定變量在去除其他變量影響后的線性相關性，并通過雙尾t檢驗和偏相關系數大小來可反映是否具有相關性及相關程度。利用SPSS 22.0軟件進行偏相關分析，分別計算NEE與不同環境因子的t檢驗P值和偏相關系數，以判斷影響NEE的主導因子以及NEE與各因子的相關程度。

2.7 主成分分析

通過降維把多個具有一定相關性的指標轉化為少數幾個相互獨立的綜合指標，從而簡化數據集。NEE受多種環境因子的影響，且這些環境因子間具有一定的相關性(周媛媛等,2017)。使用RStudio 1.1.447軟件的PCA函數對環境因子進行主成分分析得出主成分表達式，根據貢獻率加權平均得到綜合環境因子，以進一步分析NEE與綜合環境因子間的關系并進行模型模擬。

2.8 模型精度評價指標

本研究選取決定系數(R2)、均方根誤差(RMSE)和林氏調和系數(LCCC)作為評價模型模擬精度的指標(Lin,1989;陳亮等,2018)：

(9)

(10)

(11)

3 結果與分析

3.1 主要環境因子月均日動態變化特征

環境因子月均日變化具有一定的規律性(圖1)。Ta在各月均存在1個峰值(14:00—15:30)和1個谷值(5:00—6:30)，整個生長季7月出現最大值，10月出現最小值。T5的變化趨勢與Ta相近，最低點出現在5:30—6:30，最高點出現在12:00—13:30，7月出現最大值，10月出現最小值。兩者最小值出現時間相近，但Ta的峰值明顯落后于T5，這是由于研究區的土壤類型主要為砂土，其土壤熱容量比較小。PAR呈倒“U”型曲線，5:00前后開始上升，中午12:00前后達到峰值，隨后開始下降，直到19:00前后趨于平緩，最大值出現在8月，最小值出現在10月。VPD各月在凌晨4:30 — 6:00處于最低值，下午14:30—16:00達到峰值，6月出現最大值，10月出現最小值。VWC25各月日變化平緩，6:00前后達到最大值，18:00前后下降到最低值，整個生長季最大值出現在9月，最小值出現在8月。

圖1 Ta、T5、PAR、VPD和VWC25的月平均日動態變化Fig.1 Diurnal variations of monthly mean Ta,T5,PAR,VPD,and VWC25

3.2 NEE月均日動態變化特征

生長季NEE各月24 h日均值呈現出相似的變化趨勢(圖2)。各月份白天NEE均為負值，植物吸收大氣中的CO2，表現為碳匯；夜間均為正值，生態系統向大氣排放CO2，表現為碳源。NEE各月的日變化均呈單峰“U”型分布，從5:30前后開始，隨著光合作用加強，NEE由正值開始向負值逐漸變化，于中午12:00前后出現吸收峰，午后NEE吸收減緩，漸趨于0，在18:00前后由負值轉變為正值，隨后一直為正值，處于碳排放狀態。整個生長季，4—10月NEE的平均值分別為-4.96、-5.67、-7.26、-5.40、-5.27、-4.33和-2.10 μmol·m-2s-1，波動范圍分別為 -19.59～4.28，-23.02～7.18，-28.47～9.66，-23.72～7.68，-21.86～5.79，-21.29 ～ 5.23和-14.47～3.47 μmol·m-2s-1，6月NEE變化幅度最大，10月最小。

圖2 2014年生長季凈生態系統碳交換(NEE)月平均日動態變化Fig.2 Diurnal variation of the monthly mean net ecosystem exchange of CO2 (NEE) during the growing season in 2014

3.3 NEE與環境因子的線性回歸分析

NEE與Ta、T5、PAR和VPD均極顯著負相關(圖3)，P值均小于0.000 1，且與PAR的相關性最強(R2=0.465 9)，其次是T5(R2=0.227 5)、Ta(R2=0.224 2)和VPD(R2=0.173 0)。由于NEE與VWC25的回歸方程P值(0.151 5)大于0.05而未通過F檢驗，不具有顯著的相關性，但為了避免因忽略時滯現象而誤判VWC25對NEE沒有影響，還需進一步結合時滯分析結果來判斷。

3.4 NEE與環境因子時滯現象分析

小波互相關分析表明，白天NEE與各環境因子的小波互相關系數變化幅度較大(圖4)。NEE與Ta、T5和VPD的小波互相關系數分別在2.5、2和2.5 h達到最大值，表明NEE比Ta和VPD提前2.5 h達到峰值，比T5提前2 h達到峰值。而NEE與PAR的小波互相關系數在時滯為0 h時處于最大值，說明NEE與PAR同步變化。NEE與VWC25的小波互相關系數在任何時滯下沒有明顯大小關系，無法確定準確的時滯時間。夜間NEE與各環境因子在任何時滯下的小波互相關系數均沒有十分明顯的大小關系。因此，夜間NEE與環境因子間的時滯關系不明確，同樣無法找到確切的時滯時間(圖5)。

NEE與VWC25在白天和夜間均不存在時滯現象，結合回歸分析(圖3)證明VWC25不是影響NEE的環境因子，其余環境因子在夜間同樣不存在時滯關系，故在后面的分析中將VWC25剔除且不考慮其他因子夜間的時滯情況。

3.5 忽略時滯現象產生的影響

3.5.1 忽略時滯現象對判斷NEE驅動因子的影響 消除時滯后，NEE與Ta和T5的偏相關系數絕對值小幅提高，而與PAR和VPD的偏相關系數絕對值有所下降，但與PAR的偏相關系數絕對值始終最大，說明PAR是影響NEE的主導因子，時滯現象并不會影響對NEE主導因子的判斷。此外，消除時滯前T5的偏相關系數P值(0.224 8)大于0.05而未通過顯著性檢驗，即NEE與T5的線性關系不顯著。但消除時滯后的偏相關系數P值(0.005 1)比0.05小而通過檢驗，即NEE與T5有顯著的線性關系。這種變化表明不考慮時滯現象會忽略白天T5對NEE的真實影響(表1)。

3.5.2 忽略時滯現象對NEE與環境因子擬合模型的影響 對消除時滯前后的環境因子數據進行主成分分析，計算得到特征值、貢獻率、累計貢獻率(表2)和主成分載荷(表3)。

圖3 凈生態系統碳交換(NEE)與Ta、T5、PAR、VPD和VWC25日均值的回歸曲線Fig.3 Regression curves of diurnal mean NEE with Ta,T5,PAR,VPD,and VWC25

圖4 白天凈生態系統碳交換(NEE)與Ta、T5、PAR、VPD和VWC25間的時滯分析Fig.4 Time-lag anal ysis of NEE with Ta,T5,PAR,VPD,and VWC25 during the daytime,respectively矩形盒中間黑線表示小波互相關系數的中位數，上下邊緣分別表示第3分位數Q3和第1分位數Q1，兩者的差值為四分位距IQR，由矩形盒向上延伸的豎線端點為最大值(Q3+1.5IQR)，向下延伸的豎線端點為最小值(Q1-1.5IQR)，最大值與最小值以外的黑色點為異常值。下同。The middle black line of the rectangular box represents the median value of the cross-correlation coefficient.The upper and lower edges represent the third quartile (Q3) and the first quartile (Q1),respectively.The difference value between the two quartiles is the interquartile range (IQR).The vertical endpoint extending upward from the rectangular box indicates the maximum value (Q3+1.5IQR),while the vertical endpoint extending downward shows the minimum value (Q1-1.5IQR).And the black points beside the two endpoints are the outliers.The same below.

圖5 夜間凈生態系統碳交換(NEE)與Ta、T5、VPD和VWC25間的時滯分析Fig.5 Time-lag analysis of NEE with Ta,T5,VPD,and VWC25 during the night time,respectively

表1 凈生態系統碳交換(NEE)與Ta、T5、PAR和VPD的偏相關系數Tab.1 Partial correlation coefficient of NEE with Ta,T5,PAR and VPD,respectively

表2 Ta、T5、PAR和VPD 4個環境因子各主成分的特征值和貢獻率Tab.2 Eigenvalue and contribution rate of four environmental factors (Ta,T5,PAR,and VPD)

由表2可知，消除時滯前后，前2個主成分的累計貢獻率均超過90%，因此可以選取前2個主成分來反映原評價對象。根據表3可以得到消除時滯前環境因子第一主成分C1和第二主成分C2的表達式：

C1=0.523Ta+0.533T5+0.419PAR+0.517VPD，

(12)

C2=0.466Ta+0.390T5-0.772PAR-0.228VPD。

(13)

(14)

(15)

C=(74.92C1+17.65C2)/92.57,

(16)

(17)

表3 Ta、T5、PAR和VPD 4個環境因子的主成分載荷Tab.3 Principal component loads of four environmental factors (Ta,T5,PAR,and VPD)

利用NEE與消除時滯前后的綜合環境影響因子分別進行回歸分析，經過多種曲線估計后結果見表4。由表4可知，所有回歸模型的決定系數R2在消除時滯后均得到了提升。消除時滯前后，NEE與綜合環境因子的最優回歸模型均為對數回歸模型，R2分別為0.450(P<0.000 1)和0.531(P<0.000 1)。消除時滯前后NEE與綜合氣象因子的相關性均為極顯著，消除時滯后R2較消除時滯前提高了0.081。進一步利用消除時滯前后的原始環境因子數據通過最優對數回歸模型模擬NEE值，結果如圖6所示，可以看出NEE觀測值與模擬值均分布在1∶1直線兩側，消除時滯后RMSE減小了0.417 μmol·m-2s-1，LCCC提升了0.02。綜合R2、RMSE和LCCC 3個指標，忽略時滯現象影響了白天NEE數據模型的擬合精度，消除時滯后模型的擬合效果更好，精度更高。

表4 凈生態系統碳交換(NEE)與綜合環境因子(C和C’)的回歸分析Tab.4 Regression analysis between net ecosystem exchange of CO2 (NEE) and comprehensive environmental factors(C and C’)

圖6 凈生態系統碳交換(NEE)實際觀測值與模型模擬值的精度評價Fig.6 Accuracy evaluation of the observation values and the simulation values of NEE

4 討論

4.1 驅動NEE變化的主導因子

結合NEE與環境因子間的回歸分析(圖3)與偏相關分析(表1)可知，PAR對NEE日變化的解釋度(46.59%)和偏相關系數絕對值均大于其余因子，在不考慮時滯和考慮時滯2種情況下這種關系不會發生改變，所以Ta、T5和VPD對NEE沒有足夠的影響力，PAR是驅動NEE變化的最主要因子，這一結論與吳志祥等(2014)和Kwon等(2010)的結論一致。對于農田(Lietal.,2006;Guoetal.,2013)、濕地(Polsenaereetal.,2013;王莉莉等,2017)和草甸(王海波等,2014;王婧等,2015)等生態系統，輻射是驅動NEE變化的主導因子，說明陸地生態系統中主要由輻射影響植物的光合作用，進而決定生態系統的碳源/匯功能。此外，本站點NEE與VWC25沒有相關性(P=0.151 5)和明確的時滯關系，不是影響NEE變化的環境因子。這是由于研究區靠近潮白河，水分充足，地下水位在2 m左右，水分已不再是制約NEE變化的環境因子(Xuetal.,2018)。

4.2 NEE對環境因子響應產生時滯的原因

森林生態系統NEE受不同環境因子的控制而產生不同的響應機制(Borchardetal.,2015)，在這些復雜的環境因子作用下，NEE與環境因子產生了非同步變化(徐同慶等,2017)，這種非同步變化是由環境因子對植物生理過程發生作用而引起的。

溫度影響著植物的新陳代謝過程，主要通過對酶活性、微生物活性和植物根系生長等產生作用來直接或間接影響植物的光合作用與呼吸作用(Lasslopetal.,2012;溫旭丁等,2014;馬濤等,2017)。本站點NEE比Ta提前2.5 h達到峰值，比T5提前2 h達到峰值，這可能是因為中午輻射最大，光合作用也最大，NEE達到峰值，但溫度還因受光照繼續加熱，并未達到最大值(李小梅等,2015;韋志剛等,2016)，因此中午NEE依然領先于Ta和T5達到峰值。

在溫度等環境條件適宜的情況下，PAR對植物光合作用起主要作用(Krauseetal.,2013)。本站點NEE與PAR同步變化，Ouyang等(2014)也得出了同樣的結論，這是由于光直接作用于植物冠層，葉子接受光的同時開始光合作用，所以通常情況下NEE與PAR同步變化(彭鎮華等,2009;藥靜宇等,2016)。但謝瀟(2011)的研究發現，由于下墊面植株間距小，且植被十分繁茂，冠層導度小，使原本較大的湍流被分成若干較小的湍流，不能及時進行氣體交換，CO2在冠層內發生滯留，儀器需要更長的時間才能觀測到NEE數據，最終導致崇明東灘濕地CM1站點出現NEE滯后PAR 0.5 h的現象。而Jia等(2018)研究發現，NEE比PAR提前1.25 h達到峰值，這是因為其研究區域水資源有限，上午比中午具有更高的葉水勢和更強的光合能力，隨著光照增強，水分蒸發增多，阻礙光合作用，NEE先于PAR達到峰值。

VPD主要影響植物氣孔開閉行為，繼而影響植物的光合作用(譚麗萍等,2015)。本站點NEE比VPD提前2.5 h達到峰值，這是由于當VPD升高時氣孔關閉，導致氣孔導度減小，植物冠層的光合作用也隨之降低，從而影響了NEE的變化(Sulmanetal.,2016)。此外，由于VPD受到Ta的影響產生變化(Wagleetal.,2014)，所以NEE比VPD提前達到峰值的時間與NEE和Ta的時滯一致，均為2.5 h。

4.3 忽略時滯現象對分析結果的影響

目前大多數研究認為土壤溫度主要影響土壤呼吸過程(付雨龍等,2013;魏書精等,2013)，而土壤呼吸是NEE的重要組成部分，NEE的變化也應該受到土壤溫度的影響，但本研究未考慮時滯現象時NEE與T5的線性關系被忽略。因此，消除時滯現象對明確NEE的環境影響因子具有重要意義。

此外，不考慮時滯現象還會影響NEE與環境因子回歸模型的擬合精度，不利于插補缺失數據以及模擬NEE動態變化。本研究中考慮時滯后的NEE擬合模型解釋度可以提高8.1%，RMSE降低5.6%，LCCC增加2.0%，使用此模型對白天NEE數據進行擬合，模擬值與實際觀測值較不考慮時滯時更精確。土壤呼吸和枝干呼吸作為NEE的組成部分，也有研究表明忽略時滯對分析土壤呼吸和枝干呼吸與環境因子間關系有影響。Vargas等(2008)發現，在模擬土壤呼吸動態變化時，如果不考慮時滯現象的影響，會導致高估或低估土壤呼吸速率。韓風森等(2015)的研究表明，枝干呼吸與溫度存在明顯的時滯，這種時滯現象會影響計算敏感系數Q10的準確性，而且這種影響隨著測量枝干的高度降低而有所增加。因此，定量NEE及其相關變量的數值時，忽略時滯因素會影響結果的準確性，需要先排除時滯的影響。

5 結論

PAR直接作用于植物冠層而與NEE同步變化，是NEE的主導因子。由于研究區水分充足，VWC25不是制約NEE變化的因子。NEE還受到Ta、T5和VPD調控，分別比三者提前2.5、2和2.5 h變化。不考慮時滯現象會忽略白天NEE與T5的線性關系，還會使NEE的擬合精度降低，影響分析結果的準確性。考慮時滯現象并剔除其影響可以明確NEE的環境影響因子并提高NEE數據插補和模擬準確性，對估算生態系統碳源/匯量具有重要意義。