翟曉燕，張永勇

(1.中國水利水電科學(xué)研究院水利部防洪抗旱減災(zāi)工程技術(shù)研究中心，北京 100038；2.中國科學(xué)院地理科學(xué)與資源研究所陸地水循環(huán)及地表過程重點(diǎn)實(shí)驗(yàn)室，北京 100101)

流域水質(zhì)惡化已成為流域管理中的突出水問題之一，可導(dǎo)致水體功能、生態(tài)和環(huán)境系統(tǒng)遭受災(zāi)害性影響，直接損害居民公共健康和經(jīng)濟(jì)財(cái)產(chǎn)，嚴(yán)重威脅了水安全[1-4]。2018年我國約34%的水功能區(qū)水質(zhì)不達(dá)標(biāo)[5]，其中淮河流域約39%的河段水質(zhì)劣于Ⅲ類水，省界斷面達(dá)標(biāo)測次比例約為52%[6]。據(jù)統(tǒng)計(jì)，每年因水污染導(dǎo)致的經(jīng)濟(jì)損失高達(dá)2 400億元[7-8]。2011年我國政府提出了水資源管理的“三條紅線”和“三項(xiàng)制度”，2015年頒布實(shí)施了《水污染防治行動條例》，主要江河湖泊水功能區(qū)的水質(zhì)狀況得到明顯改善。然而，水污染仍然是淮河流域重要的水問題之一[9-11]。開展流域水質(zhì)指標(biāo)時空變化診斷與影響因素分析，可為探索流域水質(zhì)演變規(guī)律和成因提供依據(jù)，也對流域水污染治理具有重要的現(xiàn)實(shí)意義。

受工業(yè)市政廢水排放、城鎮(zhèn)生活污水、營養(yǎng)物流失、土地利用變化、閘壩修筑及調(diào)控等多重因子影響，流域水質(zhì)分布呈現(xiàn)顯著的時空異質(zhì)性[12-16]。統(tǒng)計(jì)分析是一種簡便有效的信息挖掘技術(shù)，可根據(jù)大量監(jiān)測站長序列實(shí)測數(shù)據(jù)挖掘潛在的流域水質(zhì)問題及其成因。其中，趨勢檢測技術(shù)，如Mann-Kendall趨勢檢驗(yàn), Sen趨勢度分析法, Spearman’s Rho檢驗(yàn)等，已被廣泛用于水情變化和水污染等問題，尤其是非參數(shù)檢驗(yàn)方法，因其對于數(shù)據(jù)結(jié)構(gòu)的要求較少，具有更好的適用性[17-20]。已有研究多未考慮水質(zhì)指標(biāo)的空間相關(guān)性，而是去除空間臨近、可能存在自相關(guān)性的站點(diǎn)，人為減少了流域樣本數(shù)量，影響了流域水質(zhì)指標(biāo)時空特征分析的準(zhǔn)確性[21]。此外，聚類分析方法可將大量斷面/站點(diǎn)尺度多個水質(zhì)指標(biāo)的變化識別為具有相似特征的典型水質(zhì)類型，從流域尺度綜合刻畫水污染狀況，已逐漸用于辨識水污染關(guān)鍵要素及其主要影響因子[22-23]。Cao等[24]結(jié)合歐式距離和動態(tài)時間規(guī)整距離，通過動態(tài)k均值聚類將養(yǎng)殖水塘中溶解氧(DO)序列分為2類，DO預(yù)測結(jié)果的均方根誤差相比于未聚類序列減小了7.6%；任婷玉等[25]采用自組織映射神經(jīng)網(wǎng)絡(luò)，將63個湖泊2006—2016年的9種水質(zhì)指標(biāo)分為3類，分別有6、27和30個湖泊污染程度較嚴(yán)重、中等和較輕。

本文利用淮河流域2008—2018年22個水質(zhì)斷面周監(jiān)測數(shù)據(jù)，綜合采用seasonal Mann-Kendall(SMK)趨勢檢驗(yàn)、莫蘭指數(shù)、動態(tài)k均值聚類、回歸分析等多種統(tǒng)計(jì)分析技術(shù)，檢測關(guān)鍵水質(zhì)指標(biāo)的時空變化特征，辨識典型水質(zhì)類型及其與不同時空尺度土地利用類型的關(guān)系，明晰淮河流域水污染關(guān)鍵要素，以期為淮河流域水質(zhì)演變特征分析和水污染防治提供科學(xué)依據(jù)。

1 研究區(qū)概況與數(shù)據(jù)來源

淮河流域(111°55′E～121°25′E，30°55′N～36°36′N)是我國水污染最嚴(yán)重的地區(qū)之一?；春痈闪靼l(fā)源于河南省桐柏縣主峰太白頂西北側(cè)，流經(jīng)湖北、河南、安徽、山東和江蘇5省，自西向東匯入長江。研究區(qū)地處我國南北氣候過渡帶，年平均氣溫為11～16℃，多年平均降水量為883 mm，多年平均徑流深為230 mm。流域總面積約27萬km2，以廢黃河為界，分為淮河水系(71%)和沂沭泗水系(29%)，多年平均降水量分別為910 mm和836 mm，多年平均徑流深分別為237 mm和215 mm。受點(diǎn)源排污和農(nóng)業(yè)非點(diǎn)源等影響，流域內(nèi)水環(huán)境狀況不容樂觀，主要受污染影響的水質(zhì)指標(biāo)包括高錳酸鹽指數(shù)(CODMn)和氨氮(NH3-N)、DO質(zhì)量濃度等。

收集研究區(qū)22個水質(zhì)斷面(圖1)2008—2018年周監(jiān)測數(shù)據(jù)，其中淮河水系18個斷面，沂沭泗水系4個斷面(徐州李集橋、棗莊臺兒莊大橋、邳州邳蒼艾山西大橋、臨沂澇溝橋)，選取pH值和DO、CODMn、NH3-N質(zhì)量濃度作為主要水質(zhì)指標(biāo)進(jìn)行分析，數(shù)據(jù)來源于中國環(huán)境監(jiān)測總站(http://www.cnemc.cn/)。收集研究區(qū)2010年和2015年土地利用類型數(shù)據(jù)，主要的土地利用類型均為旱地(55.85%和55.24%)、水田(17.89%和17.58%)和農(nóng)村居民點(diǎn)(8.66%和8.72%)，數(shù)據(jù)來源于中國科學(xué)院地理科學(xué)與資源研究所資源環(huán)境科學(xué)與數(shù)據(jù)中心(https://www.resdc.cn/)。基于ArcGIS 9.3提取各水質(zhì)監(jiān)測斷面5種尺度圓形緩沖帶(半徑分別為0.5 km、1 km、2 km、5 km、10 km)內(nèi)主要的土地利用類型面積占比，用于分析水質(zhì)類型與2008—2012年和2013—2018年土地利用類型的關(guān)系。

2 研究方法

2.1 趨勢分析

采用SMK檢測各水質(zhì)指標(biāo)的周尺度和年尺度變化趨勢，適于水質(zhì)指標(biāo)序列存在潛在的季節(jié)性、相關(guān)性和缺失值等問題[9,26]。設(shè)X={Xi|i=1,2,…,n}為n年實(shí)測水質(zhì)指標(biāo)序列，Xi={xij|j=1,2,…,m}為子樣本m周水質(zhì)指標(biāo)序列(本文中m取值為52)。假設(shè)實(shí)測水質(zhì)指標(biāo)為服從相同分布的隨機(jī)變量，且不隨時間單調(diào)變化。檢驗(yàn)的統(tǒng)計(jì)值S由式(1)計(jì)算，若序列中存在缺失值，sgn(xig-xih)=0。

(1)

標(biāo)準(zhǔn)化后的統(tǒng)計(jì)值Z服從漸進(jìn)的標(biāo)準(zhǔn)正態(tài)分布，顯著性水平p=0.05，相應(yīng)的臨界值為1.96，若|Z|>1.96，拒絕原假設(shè)，認(rèn)為在該顯著性水平下變化趨勢是顯著的。進(jìn)一步采用seasonal Kendall坡度值B估計(jì)水質(zhì)序列變化趨勢的幅度，該指標(biāo)不易受極值的影響[27]。

(2)

2.2 空間自相關(guān)性分析

采用全局莫蘭指數(shù)I診斷各水質(zhì)指標(biāo)序列的總體空間分布結(jié)構(gòu)，采用局部莫蘭指數(shù)Ii識別斷面i各水質(zhì)指標(biāo)序列的空間分布模式[28]，I和Ii取值范圍均為[-1,1]，計(jì)算公式如下：

(3)

(4)

水質(zhì)空間分布模式包括“高-高”“低-低”“高-低”“低-高”，分別表示斷面為顯著的高聚集中心且周邊指標(biāo)值較高、斷面為顯著的低聚集中心且周邊水質(zhì)指標(biāo)值較低、斷面為顯著的高聚集中心但周邊指標(biāo)值較低、斷面為顯著的低聚集中心但周邊指標(biāo)值較高。

2.3 聚類及回歸分析

采用主成分分析法和動態(tài)k均值聚類法辨識典型水質(zhì)類型。各水質(zhì)斷面的不同水質(zhì)指標(biāo)間可能存在相關(guān)性，采用主成分分析法將各水質(zhì)指標(biāo)降維轉(zhuǎn)換為獨(dú)立的綜合指標(biāo)，按主成分累積方差貢獻(xiàn)率不低于85%的準(zhǔn)則確定主成分和樣本矩陣[22,29]。以聚類組內(nèi)水質(zhì)指標(biāo)的離差平方和最小為目標(biāo)函數(shù)(式(5))，采用動態(tài)k均值聚類法將樣本矩陣分為k個互斥的水質(zhì)類型，通過迭代計(jì)算使得最終聚類分組不再變化，即各聚類組內(nèi)的水質(zhì)指標(biāo)特征較為接近，聚類組間的水質(zhì)指標(biāo)特征差異較大。

(5)

式中：J為目標(biāo)函數(shù)；x′i為第i個綜合指標(biāo)；cj為第j個聚類中心；Cj為第j個聚類組的綜合指標(biāo)數(shù)；k為聚類組數(shù)；‖·‖為歐式距離。

(6)

其中

式中：si為斷面i水質(zhì)指標(biāo)的輪廓系數(shù)；bi為斷面i的水質(zhì)指標(biāo)與其最鄰近聚類組斷面水質(zhì)指標(biāo)的平均歐式距離，表征聚類組間分離度；ai為斷面i的水質(zhì)指標(biāo)與其所在聚類組內(nèi)其他斷面水質(zhì)指標(biāo)的平均歐式距離，表征聚類組內(nèi)凝聚度。

采用線性回歸模型探索各水質(zhì)類型與不同尺度土地利用類型的相關(guān)性，其中不同尺度土地利用類型用各水質(zhì)斷面5種尺度圓形緩沖帶內(nèi)主要的土地利用類型面積占比表示，并對Pearson相關(guān)系數(shù)r進(jìn)行顯著性檢驗(yàn)，檢驗(yàn)選取的顯著性水平p=0.05。

3 結(jié)果與討論

3.1 水質(zhì)指標(biāo)時空變化趨勢

3.1.1時間變化趨勢

淮河流域2008—2018年22個水質(zhì)斷面的年變化趨勢如圖2所示。5個斷面的pH值呈顯著減少趨勢(-0.03～-0.10 a-1)；2個斷面的DO質(zhì)量濃度呈顯著減少趨勢(-0.14～-0.12 mg/(L·a))，2個斷面呈顯著增加趨勢(0.35～0.62 mg/(L·a))；9個斷面的CODMn質(zhì)量濃度呈顯著減少趨勢(-1.37～-0.17 mg/(L·a))，2個斷面呈顯著增加趨勢(0.20～0.23 mg/(L·a))；12個斷面的NH3-N質(zhì)量濃度呈顯著減少趨勢(-0.61～-0.01 mg/(L·a))，駐馬店班臺呈顯著增加趨勢(0.04 mg/(L·a))。

(a)pH值

各斷面周尺度顯著性變化趨勢主要集中在非汛期(10月至次年5月，即第41周至次年第22周)，部分?jǐn)嗝嫠|(zhì)指標(biāo)周尺度變化坡度如圖3所示(圖中“周序數(shù)”表示一年52周中的序數(shù))。淮河流域pH值基本在6～9之間，且呈減少趨勢(水體酸化)，其余3個水質(zhì)指標(biāo)整體呈輕微改善趨勢；淮河水系各水質(zhì)指標(biāo)周尺度變化坡度略大于沂沭泗水系，CODMn質(zhì)量濃度的變化坡度絕對值最大，NH3-N、DO質(zhì)量濃度和pH值次之?；春铀祊H值和DO質(zhì)量濃度的最大變幅分別為-0.25 a-1和-1.51 mg/(L·a)，均出現(xiàn)在永城黃口(第39周和第41周)，CODMn和NH3-N質(zhì)量濃度的最大變幅分別為-2.81 mg/(L·a)和-2.43 mg/(L·a)，均出現(xiàn)在亳州顏集(第12周和第21周)；沂沭泗水系pH值和DO、CODMn、NH3-N質(zhì)量濃度的最大變幅分別為-0.26 a-1、1.03 mg/(L·a)、-0.50 mg/(L·a)、-0.31 mg/(L·a)，分別出現(xiàn)在棗莊臺兒莊大橋(第1周)、邳州邳蒼艾山西大橋(第10周)、臨沂澇溝橋(第21周)、徐州李集橋(第28周)。

(a)pH值

3.1.2空間分布模式診斷

2008—2018年淮河流域pH值和DO、CODMn、NH3-N質(zhì)量濃度的全局莫蘭指數(shù)分別為0.38、0.45、0.34和0.29，均達(dá)到0.05的顯著性水平。4個指標(biāo)在淮河流域各斷面間均呈現(xiàn)顯著的空間正相關(guān)性，即鄰近斷面間的水質(zhì)指標(biāo)總體呈現(xiàn)相同的變化趨勢。各斷面水質(zhì)指標(biāo)的局部空間自相關(guān)性診斷結(jié)果如圖4所示。信陽淮濱水文站、阜南王家壩和駐馬店班臺為3個低pH值聚集中心(p分別為0.00、0.02、0.00)，主要分布在淮河干流上游和洪汝河下游；阜陽張大橋、周口鹿邑付橋閘和亳州顏集為3個低DO質(zhì)量濃度聚集中心(p分別為0.00、0.00、0.01)，主要分布在沙潁河和渦河；阜陽張大橋和亳州顏集為兩個高CODMn(p分別為0.01和0.00)和NH3-N(p均為0.00)質(zhì)量濃度聚集中心，主要分布在沙潁河和渦河。

圖4 水質(zhì)指標(biāo)空間分布模式

受近年來點(diǎn)源污染治理、水量水質(zhì)聯(lián)合調(diào)度等影響，全流域水質(zhì)指標(biāo)空間分布受外部干擾的程度有所減弱，空間異質(zhì)性降低，尤其是DO質(zhì)量濃度由隨機(jī)分布(1994—2005年)逐漸變?yōu)榭臻g正相關(guān)模式，淮河流域水污染問題的局部性緩解主要受區(qū)域人類活動和自然因素的影響[9,21]。沙潁河和渦河是受人類活動影響(如點(diǎn)源排污、閘壩調(diào)控等)最為劇烈的兩條支流，易于形成水污染聚集中心[9]，與本文研究結(jié)果一致。淮河干流上游和洪汝河下游地區(qū)土地利用以水田為主，受氮、磷肥等營養(yǎng)物流失影響，水體逐漸酸化，逐漸形成低pH值聚集中心[31]。

3.2 典型水質(zhì)類型特征辨識

淮河流域DO質(zhì)量濃度與CODMn、NH3-N質(zhì)量濃度間存在顯著的負(fù)相關(guān)關(guān)系(Pearson相關(guān)系數(shù)r≤-0.50)，CODMn質(zhì)量濃度與NH3-N質(zhì)量濃度間存在顯著的正相關(guān)關(guān)系(r=0.82)，如表1所示。

表1 水質(zhì)指標(biāo)間Pearson相關(guān)系數(shù)

表2 水質(zhì)指標(biāo)聚類評估結(jié)果

各類型水質(zhì)指標(biāo)分布如圖5和圖6所示。類型1的pH值和DO、CODMn、NH3-N質(zhì)量濃度平均值分別為7.69、7.60 mg/L、4.43 mg/L和0.70 mg/L，主要特征為弱堿性(pH值略大于7)、CODMn和NH3-N質(zhì)量濃度均偏低、DO質(zhì)量濃度介于類型2和3之間；類型2的各水質(zhì)指標(biāo)平均值分別為8.07、4.84 mg/L、10.63 mg/L和5.00 mg/L，主要特征為偏堿性(相比于類型1，pH值平均偏高4.97%)、DO質(zhì)量濃度偏低、CODMn和NH3-N質(zhì)量濃度偏高；類型3的各水質(zhì)指標(biāo)平均值分別為8.10 mg/L、8.49 mg/L、5.77 mg/L和0.58 mg/L，主要特征為偏堿性(相比于類型1，pH值平均偏高5.33%)、DO質(zhì)量濃度偏高、CODMn和NH3-N質(zhì)量濃度均偏低。各水質(zhì)類型的水質(zhì)指標(biāo)分布存在顯著差異，反映不同水質(zhì)類型的變化特征。

圖5 不同水質(zhì)類型空間分布

(a)pH值

類型1各斷面的DO質(zhì)量濃度在Ⅰ、Ⅱ類水標(biāo)準(zhǔn)限值內(nèi)(占比均為50%)，CODMn和NH3-N質(zhì)量濃度主要在Ⅱ、Ⅲ類水標(biāo)準(zhǔn)限值內(nèi)(占比均為30%～50%)，類型1主要分布在淮河上游和淮河干流；類型2各斷面的DO質(zhì)量濃度在Ⅲ、Ⅳ類水標(biāo)準(zhǔn)限值內(nèi)(占比均為50%)，CODMn質(zhì)量濃度在Ⅳ、Ⅴ類水標(biāo)準(zhǔn)限值內(nèi)(占比均為50%)，NH3-N質(zhì)量濃度均在劣Ⅴ類水標(biāo)準(zhǔn)限值內(nèi)，類型2主要分布在沙潁河；類型3各斷面的DO質(zhì)量濃度主要在Ⅰ類水標(biāo)準(zhǔn)限值內(nèi)(占比為80%)，CODMn質(zhì)量濃度主要在Ⅲ、Ⅳ類水標(biāo)準(zhǔn)限值內(nèi)(占比均為40%～50%)，NH3-N質(zhì)量濃度主要在Ⅱ、Ⅲ類水標(biāo)準(zhǔn)限值內(nèi)(占比均為40%～50%)，類型3主要分布在淮河中游和沂沭泗水系。

3.3 土地利用對不同類型水質(zhì)時空變化影響分析

類型2僅有2個水質(zhì)斷面，未進(jìn)行土地利用影響分析。類型1各斷面水質(zhì)指標(biāo)(表3)中，pH值與2013—2018年不同尺度緩沖區(qū)內(nèi)水田存在顯著正相關(guān)關(guān)系(r≥0.66,p<0.05)，與部分緩沖區(qū)內(nèi)旱地存在顯著負(fù)相關(guān)關(guān)系(r≤-0.64,p<0.05)；DO質(zhì)量濃度與2008—2012年10 km緩沖區(qū)內(nèi)其他林地存在顯著正相關(guān)關(guān)系(r=0.72,p<0.05)；CODMn質(zhì)量濃度與5～10 km緩沖區(qū)內(nèi)旱地存在顯著正相關(guān)關(guān)系(r>0.65,p<0.05)；NH3-N質(zhì)量濃度與2013—2018年5～10 km緩沖區(qū)內(nèi)旱地存在顯著正相關(guān)關(guān)系(r>0.65,p<0.05)。

表3 類型1水質(zhì)指標(biāo)與土地利用類型Pearson相關(guān)系數(shù)

類型3各斷面水質(zhì)指標(biāo)(表4)中，除2013—2018年10 km緩沖區(qū)外，DO質(zhì)量濃度與各時期城鎮(zhèn)用地存在顯著負(fù)相關(guān)關(guān)系(r<-0.80,p<0.05)；CODMn質(zhì)量濃度與2013—2018年其他林地存在顯著負(fù)相關(guān)關(guān)系(r<-0.75,p<0.05)；除10 km緩沖區(qū)外，NH3-N質(zhì)量濃度與各時期城鎮(zhèn)用地存在顯著正相關(guān)關(guān)系(r≥0.79,p<0.05)；pH值與各時期土地利用均無顯著相關(guān)性。

表4 類型3水質(zhì)指標(biāo)與土地利用類型Pearson相關(guān)系數(shù)

類型1和類型3的水質(zhì)指標(biāo)與不同時空尺度土地利用類型關(guān)系密切。類型1水質(zhì)狀況主要受農(nóng)田影響，尤其是旱地施肥造成的農(nóng)業(yè)面源污染已成為淮河上游和淮河干流地區(qū)主要污染源之一，與已有研究結(jié)果[14]基本一致。此外，水田和旱地對水體酸堿度的影響與其施肥品種有關(guān)，氮肥等酸性肥料中的游離酸可造成土壤和水體酸化，而有機(jī)肥可增加土壤有機(jī)質(zhì)和養(yǎng)分，緩沖酸化[31-32]。類型3水質(zhì)狀況主要受城鎮(zhèn)用地影響，隨著城鎮(zhèn)發(fā)展水平逐漸提高，入河污染物負(fù)荷大幅增加，加劇了水污染的惡化態(tài)勢，成為淮河中游和沂沭泗水系地區(qū)DO和NH3-N污染主要來源之一。

4 結(jié) 論

a.淮河流域各斷面水質(zhì)指標(biāo)周尺度顯著性變化趨勢主要集中在非汛期(第41周至次年第22周)，分別有23%、9%、41%和55%斷面的pH值顯著減少、DO質(zhì)量濃度顯著增加、CODMn和NH3-N質(zhì)量濃度顯著減少，淮河水系水質(zhì)指標(biāo)變化坡度大于沂沭泗水系。

b.淮河流域各斷面間水質(zhì)指標(biāo)均呈現(xiàn)顯著的空間正相關(guān)性，3個低pH值聚集中心主要分布在淮河干流上游和洪汝河下游，3個低DO質(zhì)量濃度聚集中心和2個高CODMn和NH3-N質(zhì)量濃度聚集中心主要分布在沙潁河和渦河。

c.淮河流域共劃分3種典型水質(zhì)類型。類型1為弱堿性、低CODMn和NH3-N質(zhì)量濃度斷面，分布在淮河上游和淮河干流；類型2為偏堿性、低DO質(zhì)量濃度、高CODMn和NH3-N質(zhì)量濃度斷面，分布在沙潁河；類型3為偏堿性、高DO質(zhì)量濃度、低CODMn和NH3-N質(zhì)量濃度斷面，分布在淮河中游和沂沭泗水系。

d.類型1水質(zhì)指標(biāo)主要與2013—2018年水田和旱地等顯著相關(guān)，其中水田對pH值影響較大，旱地對其余指標(biāo)影響較大；類型3的DO和NH3-N質(zhì)量濃度與各時期不同緩沖區(qū)城鎮(zhèn)用地顯著相關(guān)，CODMn質(zhì)量濃度與2013—2018年其他林地顯著相關(guān)。