西藏茶葉穩定同位素特征與區域差異及環境因子關系分析

李 梁 李春霖 王貞紅 劉振東 聶 晶 張永志 袁玉偉,*

(1 西藏農牧學院食品科學學院/西藏特色農牧資源研發省部共建協同創新中心，西藏 林芝 860000；2 浙江省農業科學院農產品質量安全與營養研究所/農業農村部農產品信息溯源重點實驗室，浙江 杭州 310021；3 西藏農牧學院資源與環境學院，西藏 林芝 860000)

我國西藏人民飲茶歷史悠久，茶葉對于藏區人民的日常飲食起著重要作用[1]。過去西藏茶葉產量較低，主要依賴從四川等地運輸。近年來，西藏高度重視茶產業發展，把茶產業作為促進經濟發展和助推脫貧攻堅的重要產業之一[2]。西藏茶葉主產于林芝市，2019年林芝市茶園種植面積為4.5萬畝，2020年春季干毛茶產量約1.55萬公斤，帶動1 400余人人均增收3 445元，茶產業已成為西藏人民致富的“金葉子”[3]。西藏主要生產紅茶、綠茶和黑茶類產品，由于高海拔和獨特的地理環境，西藏茶葉品質較優，香氣豐富，滋味醇厚，具有高原特色。然而，西藏茶葉的產地特征尚不清楚，市場中可能出現的假冒標識問題，對西藏茶葉的貿易公平和品牌保護產生不利影響。因此，挖掘西藏茶葉的產地特征，開展溯源機理研究對西藏茶葉質量監管和產業發展具有重要意義。

植物穩定同位素值(δ13C、δ15N、δ2H和δ18O)與栽培地的氣候環境和栽培方式等因素密切相關[4-5]。已有研究表明，茶葉中穩定同位素存在一定的地域差異，具有特征指示性。Liu等[6]研究發現，浙江和山東茶葉中穩定同位素比率有明顯差異，結合多元素和化學計量學方法，能夠對產地進行有效判別。Jin等[7]分析了福建安溪和武夷山茶葉中穩定同位素比率，發現兩地樣品差異顯著。植物穩定同位素在不同地域表現的差異主要是由環境和氣候因子的變化導致的，這些因子通過調控植物生理過程從而引起同位素特征發生改變[8]。目前已有相關研究探索農產品中穩定同位素與環境因子的關聯性。如Camin等[9]測定了意大利葡萄酒樣本的穩定同位素比率和環境氣候因子，發現δ18O值與氣候和地理位置信息相關性最強；Akamatsu等[10]統計分析了日本水稻中δ13C和δ18O值與環境因子之間的關系，結果表明灌漿期的日最低氣溫是影響大米同位素的主要因素；Wang等[11]研究了氣候、地形和土壤等生長環境對香蕉同位素組成的影響，發現降水和氣溫對δ18O值影響較大。已有研究發現，茶葉中穩定同位素受到環境溫度和光照強度的影響[12]，δ2H值與氣溫、光照和降水量呈現顯著相關性[13]。上述研究對探索茶葉同位素特征與環境因子關系具有一定的借鑒意義。

西藏茶葉作為我國特色農產品，其產地同位素特征和環境影響因素尚不明確。本研究通過分析西藏不同海拔地區茶葉的同位素特征，研究西藏茶葉與我國四大主產區(江北茶區、江南茶區、華南茶區、西南茶區)茶葉的同位素差異，挖掘環境因子對茶葉中同位素的影響，旨在為西藏茶葉品牌保護提供技術支撐和科學依據。

1 材料與方法

1.1 材料與試劑

本試驗材料為93份茶葉樣品，收集于2017年4月，來自全國12個省，包括27個地市及自治州，覆蓋我國江北(山東、陜西)、江南(安徽、浙江)、華南(福建、廣東、海南)和西南(貴州、四川、云南、重慶)四大茶區，其中江北茶區樣品9個、江南茶區樣品21個、華南茶區樣品15個、西南茶區樣品30個。西藏自治區樣品18個，取自林芝市不同海拔的3個取樣點，其中低海拔(750 m)樣品5個、中海拔(1 250 m)樣品8個、高海拔(2 120 m)樣品5個。所有樣品均于2017年春季按照綠茶工藝加工制成。樣品分布情況見圖1。

圖1 樣品分布及氣象信息收集圖Fig.1 Collection map of tea samples and climate information

穩定同位素標準物質：IAEA-CH-6(蔗糖，δ13CV-PDB= -10.45‰±0.03‰)、IAEA-N-2(硫酸銨，δ15Nair=20.30‰±0.20‰)，購于國際原子能機構(international atomic energy agency, IAEA，奧地利)；B2155(酪蛋白，δ13CV-PDB=-26.98‰±0.13‰，δ15Nair=5.94‰±0.08‰)，購于英國Elemental Microanalysis公司；USGS64(甘氨酸，δ13CV-PDB=-40.81‰±0.04‰)、USGS40(δ15Nair=-4.52‰±0.06‰)、USGS54(加拿大松，δ2HV-SMOW=-150.40‰±1.10‰，δ18OV-SMOW=17.79‰±0.15‰)、USGS56(南非紅象牙木，δ2HV-SMOW=-44.0‰±1.8‰，δ18OV-SMOW=27.23‰±0.03‰)，購于美國地質勘探局Reston同位素實驗室。

1.2 儀器與設備

XP6型天平，瑞士Mettler-Toledo公司；Vario PYRO cube、Vario Isotope cube型元素分析儀，Isoprime 100型、Biovision穩定同位素比率質譜儀，德國Elementar公司；HR2864粉碎機，中國飛利浦電子公司。

1.3 試驗方法

1.3.1 樣品預處理 茶葉樣品烘干至含水量不變，經粉碎機粉碎處理，過80目網篩，避光常溫保存，待測。

1.3.2 碳、氮穩定同位素檢測 稱取茶葉樣品3.0～4.0 mg，放入錫箔舟中包樣，按順序放入元素分析儀固體樣品自動進樣盤。樣品中的碳、氮元素經燃燒還原后分別轉化為純凈的CO2和N2氣體，進入同位素比率質譜儀進行檢測。元素分析儀的燃燒管主要填料為WO3,溫度保持1 150℃，還原管填料為Cu，溫度保持850℃，載氣為高純氦氣(He, 99.999%)，流量為250 mL·min-1，參考氣體為CO2和N2。元素分析儀和熱導檢測器(thermal conductivity detector, TCD)檢測流速分別為230 mL·min-1和40～50 mL·min-1。氧氣流量為40 mL·min-1，注氧時間70 s。C模式質譜參數：加速電壓為3 964 V，捕集電流為200 μA，磁場強度為4 000 mA；N模式質譜參數：加速電壓為4 162 V，捕集電流為400 μA，磁場強度為3 000 mA。樣品檢測前進行氣密性檢測、離子源真空度和穩定性測試，CO2參考氣壓力為12 psi，N2參考氣壓力為5 psi。檢測方法參照文獻[6]和[13]。

1.3.3 氫、氧穩定同位素檢測 稱取茶葉樣品0.5 mg，放入銀舟中包樣，按順序放入元素分析儀固體樣品自動進樣盤。樣品中的氫、氧元素經高溫裂解后分別轉化為H2和CO氣體，進入同位素比率質譜儀進行檢測。元素分析儀裂解爐溫度保持1 450℃，氦氣吹掃流量為150 mL·min-1，參考氣體為CO和H2。H模式質譜參數：加速電壓為4 497 V，捕集電流為600 μA，磁場強度為890 mA；O模式質譜參數：加速電壓為 4 191 V， 捕集電流為200 μA，磁場強度為3 000 mA。樣品檢測前進行氣密性檢測、離子源真空度和穩定性測試，H2參考氣壓力為20 psi，CO參考氣壓力為25 psi。檢測方法參照文獻[6]和[13]。

1.3.4 穩定同位素比率計算

(1)

式中，R樣品為檢測樣品中重同位素與輕同位素豐度之比，即13C/12C、15N/14N、2H/1H、18O/16O；R標準為國際標準品中重同位素與輕同位素豐度比[14]。測試數據均采用兩點法進行校正。

1.4 氣象信息采集

在試驗樣品采集點中，選取16個采樣點(圖1)進行氣象信息采集，涵蓋了西藏和我國茶葉主產區，收集當地環境氣象信息，包括平均相對濕度、最小相對濕度、日降雨量、日照時間、平均氣溫、最高氣溫、最低氣溫，每個環境因子均包括采樣當月(2017年4月)和采樣前3個月(2017年1月至4月)的2組平均數據，共計14項環境數據。

1.5 數據處理

采用正態檢驗和方差齊性檢驗判斷數據是否符合正態分布和方差齊性，采用單因素方差分析中的Duncan法對不同海拔樣品的同位素比率進行差異顯著性分析，采用Tamhane’s T2方法對不同產區樣品的同位素比率進行差異顯著性分析，在IBM SPSS Statistics 21軟件中進行。

通過偏最小二乘法(partial least squares，PLS)分別建立4種同位素比率與環境因子之間的關聯性模型，以樣品14項環境因子數據為自變量，穩定同位素比率為因變量，建立定量預測模型。模型共包含16個樣本，隨機選取其中14個樣本為訓練集，2個樣本為預測集。模型結果根據所有樣品的決定系數(R2)、均方根誤差(root mean square error，RMSE)和相對預測偏差(relative prediction deviation，RPD)表示，計算方法見公式(2)～(4)[15-16]。結合連續變量投影法進行變量重要性分析(variable importance for projection，VIP)，并計算自變量的標準化相關系數。PLS和VIP在XLSTAT 2019軟件中進行。

(2)

(3)

(4)

式中，N為樣品數量,ypred為同位素預測值,yref為同位素檢測值,yref為同位素檢測值的平均數，SD為同位素檢測值的標準差。

2 結果與分析

2.1 西藏不同海拔地區茶葉穩定同位素特征

本研究共采集西藏茶葉樣品18份，穩定同位素比率分析結果見表1。由于西藏海拔差異較大，樣品取自3個不同海拔的取樣點，不同海拔的樣品穩定同位素比率差異整體顯著。δ13C隨海拔變化最為明顯，茶葉樣品中δ13C值隨海拔高度的升高而減小。低海拔地區(750 m)樣品δ13C值介于-28.8‰ ～ -27.8‰ 之間，平均值為-28.3‰；中海拔地區(1 250 m) 樣品δ13C值介于-30.1‰ ～ -28.9‰之間，平均值為-29.6‰；高海拔地區(2 120 m)樣品δ13C值介于-31.9‰ ～ -30.9‰之間，平均值為-31.2‰。δ15N 值在低海拔地區樣品中最高，平均值為3.0‰，顯著高于中海拔(0.9‰)和高海拔地區(0.0‰)的樣品。δ2H和δ18O隨海拔變化趨勢一致，均為低海拔樣品最高，分別為-48.8‰和27.3‰，中海拔樣品最低，分別為-100.9‰和23.4‰。

2.2 西藏茶葉與四大主產區茶葉穩定同位素分布差異

西藏由于其獨特的地理位置，與我國茶葉主產區在環境氣候條件上差異明顯，穩定同位素比率也呈現一定的特征差異。將西藏樣品分別與江北、江南、華南、西南(除西藏)四大主產區代表性樣品進行差異顯著性分析，結果見圖2。西藏樣品中δ13C值特征最為明顯，介于-31.9‰ ～ -27.8‰之間，平均值為-29.7‰，顯著低于其他主產區樣品。δ15N值介于-2.2‰ ～ 5.0‰之間，平均值為1.2‰，與其他產區樣品無顯著差異。δ2H值介于-111.5‰ ～ -40.5‰ 之間，平均值為-83.3‰，顯著低于西南、江南和江北產區樣品，與華南產區樣品無顯著差異。δ18O值介于22.2‰ ～ 27.6‰之間，平均值為25.3‰，顯著高于華南產區樣品，與其他產區樣品均無顯著差異。

注：**代表P<0.01，*代表P<0.05，ns代表P>0.05。Note: **represent P<0.01, * represent P<0.05, ns represent P>0.05.圖2 西藏與四大主產區茶葉中穩定同位素箱線圖Fig.2 Boxplot of stable isotopes in tea samples from Tibet and four main producing area

2.3 茶葉穩定同位素與環境因子關系模型

為探索西藏茶葉同位素特征形成的氣候環境影響，本研究選取西藏和我國茶葉主產區范圍內16個取樣點進行環境因子數據采集和茶葉同位素比率檢測，分別將δ13C、δ15N、δ2H、δ18O與14項環境因子數據進行PLS分析，定量預測模型結果見表2。其中R2和RPD越高，RMSE越低，模型效果越優[15-16]。4個模型的R2均高于0.55，RMSE均低于5.09，RPD介于1.51～3.08之間。δ2H模型預測效果最優，R2和RPD均最高，分別為0.90和3.08，表明δ2H與14項環境因子關聯性更強。由于δ2H絕對值較其他同位素比率偏大，而RMSE受自變量絕對值影響，導致δ2H模型RMSE最高，為5.09。其次為δ13C模型，R2為0.84，RPD為2.63，RMSE最低，為0.48。δ18O模型R2為0.75，RMSE為0.64，RPD為2.13。δ15N模型預測效果最差，R2和RPD均最低，分別為0.55和1.51，RMSE為1.20，表明δ15N與14項氣候環境因子關聯性較弱。

表1 西藏不同海拔茶葉樣品穩定同位素分布表Table 1 Distribution of stable isotopes in tea samples from Tibet with different altitudes /‰

對模型進一步進行變量重要性分析和標準化相關系數分析，VIP值大于1的環境因子對同位素比率影響較大[17]，結果見圖3。δ13C主要受日照和氣溫影響，重要性最高的3個因子為4月(采樣當月)日照時數、平均氣溫和最低氣溫，3個變量因子的標準化相關系數均為正值，表明這3個因子與δ13C成正相關。在δ15N預測模型中，4月日照時數、1-4月(采樣前3個月)日照時數和1-4月平均相對濕度影響較大。δ2H模型顯示，日照、氣溫和相對濕度均會對δ2H影響較大，VIP大于1的因子以相對濕度為主，VIP最高的3個因子分別為4月日照時數、1-4月最高氣溫和1-4月最小相對濕度。δ18O則主要受氣溫影響，VIP大于1的因子主要包括1-4月和4月的最高氣溫、平均氣溫和最低氣溫，6個因子的標準化相關系數均為負值，表明這6個因子與δ18O成負相關。

表2 偏最小二乘模型結果Table 2 Results of partial least squares models

圖3 基于偏最小二乘的變量重要性分析與標準化相關系數Fig.3 Variable importance for projection and standardized coefficients based on partial least squares

3 討論

3.1 環境因子對碳、氮穩定同位素的影響

本研究結果發現，相較于其他同位素比率，西藏茶葉中δ13C值隨海拔變化規律最明顯，茶葉中δ13C值隨海拔升高而降低。經PLS-VIP分析表明，δ13C值主要受采樣當月日照和氣溫影響，且δ13C值與日照時數和平均氣溫成正相關。而氣溫隨海拔的升高而降低，因此茶葉中δ13C值隨海拔的升高而降低。在西藏與四大主產區茶葉同位素比率比較研究中發現，西藏茶葉δ13C值顯著低于四大主產區茶葉。而西藏林芝市4月平均氣溫僅10℃，較其他主產區偏低。該結果與PLS分析中δ13C值與氣溫成正相關的結果一致。已有研究分析了廬山不同海拔茶樹光合響應的差異，結果表明茶樹在低海拔處的光合作用和光系統Ⅱ活性顯著低于高海拔處[18]。此外，氣溫和CO2濃度升高能夠明顯提高茶樹葉片的光合效率[19]，而光合作用與植物13C穩定同位素分餾密切相關[20-21]，進一步解釋了茶葉δ13C值隨海拔升高而降低的原因。

經單因素方差分析，西藏不同海拔樣品中δ15N值差異不如δ13C值明顯，而西藏與不同產區茶葉樣品中δ15N值差異不顯著，并且在PLS分析中δ15N值與環境氣象因子關聯性較弱。這說明相較其他同位素比率，δ15N值受氣候影響較小。已有研究表明，植物中δ15N值與栽培和土壤關系更為密切，主要受到栽培方式[22]、土壤類型[23]、根系情況[24]等多方面影響。本研究有待進一步收集西藏茶葉種植環境樣品，探究土壤背景、栽培方式以及肥料種類對茶葉同位素的影響。

3.2 環境因子對氫、氧穩定同位素的影響

不同海拔茶葉樣品中δ2H和δ18O值分布規律較為一致，但西藏茶葉與主產區茶葉同位素比較中，δ2H和δ18O值分布規律有較大差異，推測與環境影響因素有關。已有研究表明，δ2H值與降雨、氣孔開合等密切相關，具有顯著的緯度效應、海拔效應和內陸效應[25]。本研究PLS結果表明，δ2H與環境氣象因子關聯性最強，利用14項環境因子能對δ2H進行定量預測。VIP結果顯示，光照、氣溫、濕度多種環境因素均會影響茶葉中δ2H值，VIP大于1的變量中以相對濕度為主。該結果與已有研究結果較為一致，如Deng等[13]利用同位素和多元素技術對我國綠茶進行產地溯源，探究了不同產地綠茶的同位素特征及其與氣候環境的相關性，發現茶葉中δ2H值與氣溫、光照和降水量相關性顯著；Liu等[26]研究了不同產地綠茶的穩定同位素比率，發現湖北的綠茶樣品中δ2H和δ18O值顯著高于貴州和廣東省樣品，推測與環境相對濕度和蒸騰速率等因素有關。本研究PLS-VIP分析結果表明，氣溫是茶葉中δ18O值的主要影響因素，影響δ18O值的氣溫包括了采樣當月和采樣前3個月的平均氣溫、最低氣溫和最高氣溫。此外，除了氣溫對δ18O值的影響外，也已有研究表明植物纖維素中δ18O值與羰基和水之間的平衡同位素效應等因素相關[27]。

4 結論

本研究結果表明，西藏茶葉中穩定同位素比率具有一定的海拔規律和區域特征，其中δ13C值特征最為明顯，隨海拔升高而降低，顯著低于其他主產區樣品。茶葉穩定同位素比率均受到環境因子影響，其中δ2H值與氣象環境關聯性最強，主要與相對濕度等因素相關。在同位素特征分析的基礎上，后續研究可結合代謝組學、光譜學、質譜學等技術，進一步挖掘西藏茶葉品質特征，助力西藏茶葉品牌建設與推廣，推動西藏茶產業快速發展。