程強強,王 哲,3,易發成,孫 卿,張怡萍,蔣 虎,王成霞
(1西南科技大學環境與資源學院,四川綿陽621010;2固體廢物處理與資源化教育部重點實驗室,四川綿陽621010;3中國科學技術大學地球和空間科學學院,合肥230026)
關鍵字:茶園;潛在毒性元素;內梅羅指數;地積累指數;潛在生態風險指數;多元分析;空間分析
中國作為茶的故鄉,人工種植茶葉可以追溯到6000多年前的浙江余姚田螺山一帶,而制茶和飲茶也有上千年歷史。近年來,中國茶葉發展迅速,截止2019年底,全國茶園面積為306.52萬hm2,茶葉年產量約260萬t,分別占世界的60%和45%,穩居世界第一位。茶葉中富含許多對人體健康有益的礦物元素和化學組分,具有降血壓、降血糖、抗氧化、抗衰老和緩解疲勞等功效[1]。近年來,由于茶園中普遍使用化肥、除草劑、農藥等農業生產資料以及茶園周邊工業發展,茶園土壤、茶葉等不同程度地受到重金屬的污染,茶園環境中重金屬污染日益加劇,導致茶葉中也含有一些對人體有害的潛在毒性元素(Potentially toxic elements,PTEs),PTEs通過土壤-植物生態系統在茶葉中富集并當累積達到安全閾值時,可通過飲茶途徑進入人體,從而構成健康風險。相關研究表明,中國大多茶園土壤和茶葉中PTEs含量處在限量范圍內,部分茶園受到不同PTEs不同程度的污染,茶園土壤以Cr、Cd、As和Hg污染為主,商品茶中以Pb和Cd污染為主[2]。隨著人們對茶飲料功效的認識不斷提高,關于茶葉PTEs與飲茶健康的問題也備受國內外學者關注[3]。學者們已圍繞茶葉中PTEs含量與茶樹器官、大氣沉降、土壤環境質量和施肥[4]等開展了深入的研究工作。以往相關研究主要側重于茶園土壤污染現狀的分析,其主要采用內梅羅綜合指數法[5-6]、地積累指數法、富集因子、潛在生態危險指數、土壤和農產品綜合質量指數法等評價方法[7]。此外,隨著地理信息科學技術在土壤污染領域的應用,利用地統計分析和地理信息分析系統ArcGIS成為PTEs空間分布與污染程度評價的主流分析手段[8]。
以往的工作大多局限于污染程度方面的研究,對土壤中PTEs污染的研究大多數以元素的分布、賦存形態和生態環境危害與修復等為主要內容,而對土壤中PTEs的污染來源解析方面的研究較少。而調查土壤PTEs的來源,是從源頭上控制并解決土壤PTEs污染問題最有效的途徑之一。現有的土壤PTEs源解析方法主要包括同位素的比值分析、多元統計分析及基于GIS空間分析。在土壤PTEs來源分析方法中,以同位素示蹤法最為精準,但該方法成本較高,且不適于大范圍分析[9],其他常用分析方法為主成分分析(Principal Component Analysis,PCA),聚類分析方法(Cluster Analysis,CA),地統計分析[10]等方法解釋PTEs污染來源[11]。葉宏萌等[12]通過主成分分析對武夷山土壤重金屬來源進行分析,得到潛在生態風險較高的元素Cd和Hg來源于農業污染源;張曉茹等[13]對亞青會期間南京的大氣降塵中重金屬進行聚類分析,得到Cd、Pb和Sn來源于工業排放,Fe、As、Mn和Sb來源于燃煤排放,Co、Cu、Cr、V、Ni和Zn來源于道路揚塵和機車尾氣排放;戴彬等[14]通過地統計分析對山東萊蕪市鋼城進行重金屬來源分析,得到As、Cd、Hg、Cu、Pb和Zn的污染范圍與工業企業的位置相一致。可見,多元統計分析是一種綜合性的統計分析方法,能夠分析多指標和對象相互關聯的規律,為分析樣點中各類PTEs間相互關系提供了有力工具,可以較為可靠的分析土壤中PTEs的來源[15]。而空間分析是基于地理對象位置和PTEs總量數據的分析技術,分析土壤中PTEs的污染情況,可以揭示出多元分析等研究方法難以發現的規律,應用空間分析技術分析異常空間分布與污染源的關系可以直觀地判斷出可能的異常成因[16]。
因此,為更好的揭示宜賓市屏山縣某茶園土壤PTEs污染現狀與來源,通過對研究區樣品釆集和檢測分析,以多元統計和基于ArcGIS空間分析方法相結合的手段對研究區茶園土壤中7種PTEs:Cd、As、Cr、Cu、Pb、Ni和Zn進行污染評價和來源分析,可避免利用單一方法來進行土壤PTEs污染來源分析造成較大的局限和誤差,為研究區茶園土壤PTEs污染綜合防治提供可靠的科學依據。
研究區位于四川省宜賓市屏山縣某茶園,屬中亞熱帶濕潤型季風氣候區,垂直氣候差異大。年均氣溫為14.9℃,年均日照為950.7 h,無霜期為300天,降水量為1066 mm。在現場踏勘基礎上,按照《土壤環境監測技術規范》(HJ/T 166—2004)[17],采用系統隨機布點與分塊隨機布點相結合的方法,對選定研究區內2處茶園開展土壤樣品采集。為了保證采集的土壤樣品具有代表性,2個茶園(A區和B區)布置了31+43個采樣點,共采集74個土壤樣品,土壤采樣點分布圖如圖1所示。土壤樣品采集0~20 cm深表層土,每個采樣點采用梅花布點法分別采取5個1 kg的子樣本混合成一個土壤樣品。

圖1 宜賓某茶園土壤采樣點分布
將采集的土壤風干后磨細通過100目的尼龍篩網,并測試土壤中的PTEs元素的總含量和pH值。本文采用電熱恒溫石墨消解儀進行濕法消解:使用分析天平稱取100目土樣樣品0.1000 g,放入50 mL石墨消解儀配套的四聚氟乙烯消解罐中;先加入5 mL密度為1.42 g/mL的HNO3浸潤0.5 h去除土壤中的有機物,再加入2 mL密度為1.49 g/mL的HF,最后加入2 mL HClO4,將石墨消解儀的升溫階段設置為100℃加熱2 h,再以180℃恒溫消解24 h以上,直至有淡黃色粘稠狀液體產生;將其冷卻至室溫后,加入2%的HNO3后用25 mL容量瓶定容;用0.45 μm孔徑的微孔過濾膜對定容的樣品進行過濾并低溫保存,樣品中所含總砷,鎘,鉻,鎳,鋅,銅和鉛總量采用西南科技大學分析檢測中心的ICP-OES進行測定。
在多元統計中對樣品數據的正態分布研究是必不可少的,為避免指標之間的水平差異太大,需將其轉換為無量綱指標,并通過變異系數反映樣本數據集的分散情況。同時,為反映出正態分布樣品的偏態和尾重,可通過箱線圖來觀測多組數據中異常值。本文采用Z-Score標準化來確保數據結果的穩定性,并應用SPSS 25.0軟件開展茶園土壤PTEs含量的相關性分析、PCA和CA,尋找區域內潛在污染源與土壤PTEs分布特征之間的聯系。其中,相關性分析主要用于分類和分析土壤中各個變量之間的Pearson相關性;PCA用于解釋土壤中PTEs的來源,并采用Varimax旋轉減少每個因素的高負荷下變量數量,可便于結果的解釋;CA主要通過聚類樹形圖對土壤中七個PTEs進行分類,以進一步分析PTEs的來源。
空間分析是以空間數據為基礎、以地學原理為依托的分析技術。地理信息系統(GIS)是空間分析技術的主要平臺[16]。本次研究利用ArcGIS10.0軟件進行反距離加權插值法(IDW)可得到研究區茶園土壤中PTEs濃度空間分布及污染程度特征圖,為開展污染評價和來源分析提供基礎數據。
參照《土壤環境監測技術規范》(HJ/T 166—2004),將土壤中PTEs環境質量評估模型分為單因素污染指數法和綜合污染指數法,通過單因子指數法可以確定某種重金屬污染物污染程度。計算見公式(1):

在本文研究中,單因子指數只能反映單個元素的污染程度,并不能全面的反映土壤的整體污染情況。因此采用內梅羅綜合污染指數法進行分析,也稱為Nemerow指數,它可以充分反映各種污染物對土壤的共同作用影響,并突出顯示高濃度污染物對環境質量的影響。Nemerow綜合指數由公式(2)進行廣泛計算:

地積累指數是用來評價PTEs污染的定量指標,該指標彌補了其他方法的不足,既考慮環境自然背景值與自然成巖作用對污染評價結果的影響,同時也考慮人為污染因素。地積累指數為公式(3):

由Hakanson提出的潛在生態風險指數是一種根據重金屬性質及環境行為特點,從沉積學角度對土壤重金屬污染進行評價的方法。該方法綜合考慮了多元素協同作用、毒性水平、污染濃度以及環境對重金屬污染敏感性等因素,以反映特定研究區域內多種污染物對環境的綜合影響。其計算見公式(4)~(6):

研究區茶園土壤pH與PTEs(As、Cd、Cr、Cu、Ni、Pb、Zn)進行分析,其結果如表1所示。首先,繪制了研究區茶園土壤PTEs與土壤pH的箱線圖,為避免因實驗環節導致的誤差,使結果更加接近研究區實際水平,將箱線圖中的2個異常值點去掉(如圖2);其次,將研究區茶園土壤PTEs含量通過Z-Score標準化處理,并將標準化結果通過K-S檢驗,結果表明標準化后的pH和PTEs濃度值符合正態分布。
由表1可以看出,研究區茶園土壤中As、Cd、Cr、Pb、Cu、Ni、Zn 的平均值分別為 22.75、3.21、89.45、37.78、34.24、41.23、101.88 mg/kg,分別是四川省土壤背景值(A層)的2.15、15.28、0.95、1.06、0.38、1.28、3.33倍。由此可見As、Cd、Ni、Zn在茶園土壤中有明顯的富集情況,而Cr、Pb、Cu未見明顯富集。

表1 各個茶園研究區土壤PTEs富集含量

圖2 砷、鎘、鉻、鎳、銅、鋅、鉛和pH的箱線圖
通過對研究區茶園土壤中PTEs和土壤pH的Pearson相關性分析(如圖3所示),結果表明,Ni分別與Cd、Cr元素存在強相關性(γ2=0.82、0.83,P<0.01);Zn與Cd,Cr,Ni之間存在顯著性相關(γ2=0.62、0.71、0.64,P<0.01);Cd與Cr也同樣存在顯著性相關(γ2=0.76,P<0.01);As,Cu,Pb,Zn之間未見顯著性相關關系。這說明Cd,Cr,Ni和Zn的富集可能存在相同的來源。

圖3 Pearson相關性分析
利用公式(2)和(3)對研究區的采樣點Nemerow指數和地積累指數進行計算,如圖4所示。其中,Nemerow指數最大值9.05,最小值5.5,平均值7.65,污染指數值均大于3.0,根據表2的污染程度分級標準可知研究區茶園土壤PTEs污染屬于重度污染;地積累指數最大值3.87,最小值3.03,平均值3.38,地累積指數(Igeo)介于3~4之間,屬于重度污染。Nemerow指數與Igeo的評價結果相同,表現出較高的一致性。

表2 污染程度、地積累指數及潛在風險指數分級表

圖4 內梅羅指數與地積累指數散點分布圖
利用公式(4)~(6)對研究區茶園土壤PTEs進行潛在生態風險評價,其中As、Cd、Cr、Cu、Ni、Pb、Zn的單元素的污染指數(Cf)分別為Cd(40.63)?As(2.19)>Pb(1.27)>Cu(1.22)>Zn(1.18)>Cr(1.13)>Ni(1.10),單個元素環境風險指數(Er)分別為Cd(1218.80)?As(21.876)>Pb(6.323)>Cu(6.11)>Ni(5.51)>Cr(2.27)>Zn(1.178),而研究區整體的潛在生態風險指數為RI=1262.1,如圖5所示。

圖5 污染指數與潛在生態風險評價指數圖
為了能夠直觀看出研究區茶園土壤中PTEs空間分布情況和污染程度,輔助判定PTEs來源的類型與方位。本次研究利用ArcGIS10.0軟件進行反距離加權插值法(IDW)可得到了研究區茶園土壤中As、Cd、Cr、Cu、Ni、Pb、Zn濃度空間分布及污染程度特征圖,如圖6,圖7,圖8所示。

圖6 研究區(A區)PTEs含量空間分布圖

圖7 研究區(B區)PTEs含量空間分布圖

圖8 研究區土壤中PTEs污染內梅羅、地積累、潛在生態風險指數分布圖
2.3.1 主成分分析 通過Varimax和Kaiser進行歸一化旋轉,得到三個特征值大于1的主成分,其總方差解釋矩陣結果如表3所示。其中,第一組主成分(C1)為Cd,Cr,Ni和Zn,第二組主成分(C2)為Cu和pH,第三組主成分(C3)為As和Pb,如圖9所示。從表3可知,三組主成分分別占總因素所解釋方差的42.76%,16.95%,15.19%,三者累積占總因素所解釋方差的74.91%,為推測研究區土壤中PTEs的來源提供依據。

圖9 3種主要因子載荷圖
2.3.2 聚類分析 為了進一步分析PTEs的來源,尋找區域內潛在污染源與茶園土壤PTEs分布類型之間的聯系,本次研究將標準化后的各個因子采用離差平方和法(Ward Linkage Method)進行聯接,計算每個因子之間的歐氏距離(Euclidean Distance),并根據距離將因子進行聚類分析,其結果如圖10所示(其中:X軸為歐氏距離,Y軸為Z-Score標準化后的變量因子)。

圖10 PTEs和pH聚類分析結果
將研究區中各PTEs的平均值與《土壤環境質量農用地土壤污染風險管控標準(試行)》(GB15618—2018)中的篩選值相比,Cd的濃度遠超其風險篩選值,表明研究區土壤有被該元素污染的風險;而其他元素未超過風險篩選值,未構成污染風險。各元素變異系數均在0.1~1之間,表明數據分布為中等變異程度,其中As和Pb元素變異程度最大,說明可能受人類活動影響較大。通過對比我國部分地區茶園土壤中PTEs含量(見表1),可明顯看出研究區茶園土壤中Cd、Cr、Ni、Cu的含量顯著高于其他研究區水平,表明這4種元素存在過度富集,可能會造成茶園土壤污染。而研究區茶園土壤中As、Zn和Pb元素富集水平接近其他研究區水平(除四川地區),未見明顯富集。同時,從四川、浙江、江西、山東、貴州、廣東、湖北等地區茶園土壤中Cd的污染情況來看,其富集濃度均超過中國平均水平,表明Cd是大部分茶園研究區土壤中的主要污染元素。
研究區兩處茶園土壤中PTEs的濃度分布如圖6(A區)和7(B區)所示。研究區兩處茶園土壤中Cr、Ni、Cd和Zn的濃度空間分布圖中顏色分層分布情況大致相同,其中:A區茶園土壤中污染最嚴重區集中在該區中部和東南部區域,而西北方向污染較輕;B區茶園土壤中污染最嚴重的區域集中在該區南部區域,可能為來自同一污染源。研究區兩處茶園土壤中As和Pb濃度空間分布圖中相同顏色深度的區塊面積較大,其中:A區顏色較深的區域主要集中在該區中部和西北區域;B區除西北側外,其余區域顏色分層分布情況相同,可能為同一污染源導致。研究區A區茶園土壤中Cu污染最嚴重區主要集中在該區中部;B區茶園土壤中Cu污染最嚴重區主要集中該區北部。
研究區茶園土壤PTEs污染的內梅羅污染指數、地積累指數和潛在生態風險指數分布如圖8所示。3種指數分布圖中顏色分層分布情況與圖6和7中Cd的顏色分層分布情況相似,說明了重金屬Cd對整個研究區土壤中PTEs的污染情況、污染程度及潛在生態風險等級起決定性作用。同時,由表2及潛在生態風險等級劃分標準可知,除Cd為重度環境風險(Er≥320),其他元素均為輕度環境風險(Er<40);研究區綜合潛在生態污染風險(RI=1262.06),表明茶園土壤PTEs為極強污染風險等級,其中Cd的潛在生態污染風險的貢獻率為84.4%,說明了重金屬Cd對于整個研究區茶園土壤PTEs潛在生態風險評價起到控制性影響。
此外,相關研究表明茶樹中不同生長部位的Cd含量差異較大,向素雯等[32]認為茶樹各器官部位中Cd含量大小為莖>老葉>新葉,其中莖中Cd含量是老葉的2.4倍,是新葉中的8倍。雖然茶園土壤中Cd的潛在生態污染風險高,但新葉作為成品茶葉的原料,茶園土壤中重金屬Cd對成品茶的污染有限。因此,為了減少Cd元素在茶園土壤中富集及確保茶葉生產安全,應對樹莖和老葉進行及時清理,可避莖葉腐爛后重金屬Cd重新進入土壤-植物的循環系統。
3.3.1 主成分分析 PTEs的組分矩陣和旋轉組分矩陣分析結果表明(如表3所示)。第一組(C1)中Cd、Cr、Ni、Zn正載荷較高,分別為0.90、0.92、0.92和0.80,元素之間相關性較高且變異系數相近(0.13~0.16)。段續川等[33]認為同一組分中載荷較高的元素可判斷為同一來源,則研究區茶園土壤中Cd、Cr、Ni、Zn可判斷為同一污染源。PTEs的來源分為人類活動和自然背景積累,4種污染物濃度高于四川省土壤背景值,說明自然背景積累并非是主要來源。而研究區區域內主要的礦業生產活動為硫鐵礦的采選冶[34],其早期礦業生產過程中廢氣、廢水、廢渣等導致Cd、Cr、Ni、Zn在水體與大氣中富集,并經引水灌溉與大氣降塵形式等在茶園土壤中積累,這與李志濤等[35]研究認為Cr、Cd、Ni、Zn在農用地土壤中富集與硫鐵礦礦區人類活動有關的結論相一致。

表3 總方差解釋矩陣和旋轉解釋矩陣
第二組(C2)中Cu和pH所占方差貢獻率較小,因子載荷為0.75和0.79。而據相關報道[36],研究區依托當地竹產業發展起來的1118家土法造紙廠,生產過程中產生的廢水排放到小龍溪內,導致水體及周邊土壤污染,2016年當地環保部門對土法造紙廠進行了關停。龔鑫梅等[37]研究認為造紙廠廢水等污染物的排放在一定程度上是土壤PTEs的主要來源,且土壤pH與Cu存在明顯的正相關性。本次研究表明同一組分Cu和pH具有較高的載荷,可以判斷Cu的富集和土壤酸化為同一污染來源,與上述分析相一致。
第三組(C3)中As和Pb所占方差貢獻率較小,因子載荷分別為0.57和0.90。煤炭資源也是研究區的主要礦產資源之一[34],Pb和As是煤礦中標志性元素[38]。劉海彪等[39]研究表明四川地區的煤中重金屬含量Pb較高,As其次。研究區茶園土壤中Pb載荷較高,表明工礦企業燃煤引起產生的大氣降塵是導致Pb在土壤中富集的主要原因;而As在C3中載荷較低且在C1組分中載荷為0.46,表明茶園土壤中As的富集推斷由大氣降塵與其他污染源疊加導致,而茶園施用有機肥與化肥也是導致茶園土壤As富集的原因之一。這與麻萬諸等[40]認為土壤中As的來源主要為大氣沉降,其次為有機肥與化肥的使用的研究結論相一致。
3.3.2 聚類結果討論 以歐氏距離15為標準,8種因子被分為三類且與上述PCA分析結果相同;以歐氏距離2為標準,Cr與Ni為一組,Cd,Zn,As,Pb,Cu和pH分別各自為一組。按查瓦里茨基分類[41]分為親硫元素和鐵族元素,Cr和Ni為鐵族元素,離子半徑相近,所以容易通過同晶置換反應存在于硅酸鹽礦物中;根據戈爾德施密特的分類[41],Cd、Zn、As、Pb、Cu為親硫元素,與硫的親和力強,從側面證明這幾種元素為早期硫鐵礦山開發和燃煤使用所導致。
茶園土壤中7種PTEs中僅有Cd超過農用地土壤篩選值,且As、Ni和Zn的濃度值均超過土壤背景值,表明這4種元素在研究區土壤中均存在富集現象。從污染程度特征來看,茶園土壤內梅羅污染指數和地積累污染指數平均值分別為7.65和3.38,為重度污染;而綜合潛在生態污染風險指數平均值為1262.06,表明茶園土壤PTEs潛在生態風險高。多元分析結果表明,Cd、Cr、Ni和Zn的富集可能為硫鐵礦開采導致;Pb主要為鐵礦的采選冶和大量燃煤使用導致;As除了燃煤使用外,還與農田中肥料的使用有關;造成茶園土壤中Cu富集及酸化與區內小龍溪沿岸土法造紙廠廢水排放密切相關。