蔡燕秋 唐翔宇 張建強 程 軍 唐家良 章熙峰

(1.西南交通大學(xué)地球科學(xué)與環(huán)境工程學(xué)院，四川 成都 611756；2.中國科學(xué)院、水利部成都山地災(zāi)害與環(huán)境研究所，四川 成都 610041)

河流是流域內(nèi)污染物遷移進入湖泊、水庫等水體的主要途徑。河流的水質(zhì)指標和徑流量是污染物流失通量計算的重要依據(jù)。由于高頻次水質(zhì)監(jiān)測的成本較高，我國各地區(qū)水質(zhì)監(jiān)測頻率一般為1～4次/月，在污染物流失通量計算時需要通過線性插值法或非線性插值法進行估算。在人類活動強度高的區(qū)域內(nèi)，水質(zhì)與水文情況復(fù)雜，線性插值法嚴重偏離實際情況[1]。已有多種模型如SWAT、雨水管理模型(SWMM)等[2-6]被應(yīng)用于較大流域污染物流失通量的估算，但此類模型并不適用于小流域。

河流污染物流失通量與徑流量間一般存在著較好的相關(guān)關(guān)系，美國地質(zhì)調(diào)查局(USGS)據(jù)此開發(fā)了LOADEST模型。該模型能利用連續(xù)的徑流量數(shù)據(jù)和離散的水質(zhì)數(shù)據(jù)，建立污染物流失通量與徑流量之間的非線性回歸方程，因其符合污染物波動式流失的常態(tài)，并且操作簡便、數(shù)據(jù)需求量小，已被廣泛應(yīng)用[7-10]。宋方方[11]通過LOADEST模型對趙家溪污染物流失通量進行估算，為趙家溪污染控制方案的制定提供數(shù)據(jù)支撐；DWIVEDI等[12]研究發(fā)現(xiàn)，LOADEST模型用于估算低水平范圍的大腸桿菌負荷有較好的效果；JHA等[13]用LOADEST模型對美國紐斯河的氮磷流失通量進行估算，發(fā)現(xiàn)估算值與實測值相比總體偏高；PARK等[14]指出回歸方程不僅與水質(zhì)參數(shù)有關(guān)，而且與暴雨事件樣本所占比例有關(guān)；高信娟[15]在九龍江的研究結(jié)果表明，水質(zhì)數(shù)據(jù)缺失越少，污染物流失通量的估算值越接近實測值。

本研究以四川省鹽亭縣萬安丘陵小流域截流、大興及萬安3個河流斷面為對象，比較不同監(jiān)測時間間隔條件下總氮(TN)、總磷(TP)日流失通量估算值相比實測值的相對誤差，以確定既經(jīng)濟可行又結(jié)果可靠的最低監(jiān)測頻次，并揭示小流域?qū)傩缘目臻g尺度差異所產(chǎn)生的影響。

1 材料與方法

1.1 研究區(qū)

研究區(qū)位于四川省鹽亭縣萬安小流域(105°27′E，31°16′N)，海拔為389～426 m，屬中亞熱帶季風(fēng)氣候。2015年、2016年和2017年的年降雨量分別為1 029.3、886.7、628.7 mm，其中5—10月降雨量占全年總降雨量的比例分別為87%、83%和76%。該小流域內(nèi)截流、大興、萬安集水面積分別為0.35、4.80、12.36 km2。研究區(qū)概況見圖1。

圖1 萬安小流域概貌與監(jiān)測點位置Fig.1 Map of Wan’an catchment and monitoring locations

1.2 LOADEST模型

LOADEST模型參數(shù)的估計方法有極大似然估計(MLE)法、校正極大似然估計(AMLE)法和最小絕對偏差(LAD)法3種。一般采用AMLE法，但當(dāng)殘差不符合正態(tài)假設(shè)時，則考慮采用LAD法。當(dāng)源文件中將模型選擇設(shè)置為0時，LOADEST模型運用Akaike信息準則(AIC)[16-18]和Schwarz后驗概率準則(SPPC)[19]，自動從9個回歸方程(見式(1)至式(9))中取得最小AIC值和SPPC值的方程為最優(yōu)方程，并進行自變量中心化消除多重共線性[20]。

lnL=a0+a1lnQ

(1)

lnL=a0+a1lnQ+a2lnQ2

(2)

lnL=a0+a1lnQ+a2td

(3)

lnL=a0+a1lnQ+a2sin(2πtd)+a3cos(2πtd)

(4)

lnL=a0+a1lnQ+a2lnQ2+a3td

(5)

lnL=a0+a1lnQ+a2lnQ2+a3sin(2πtd)+
a4cos(2πtd)

(6)

lnL=a0+a1lnQ+a2sin(2πtd)+a3cos(2πtd)+a4td

(7)

lnL=a0+a1lnQ+a2lnQ2+a3sin(2πtd)+
a4cos(2πtd)+a5td

(8)

(9)

式中：L為污染物流失通量；a0～a6為回歸系數(shù)；Q為徑流量；td為分數(shù)化的時間與中心化時間之差；參數(shù)單位均視具體情況而定。

LOADEST模型主要通過以下3種方式檢驗回歸方程的有效性：(1)判定系數(shù)(R2)越接近1，證明方程的擬合程度越好；(2)殘差序列相關(guān)系數(shù)(SCR)越接近0，證明方程的各殘差之間越獨立；(3)概率曲線相關(guān)系數(shù)(PPCC)越趨近于1，表明殘差越接近正態(tài)分布。

1.3 方法比較與優(yōu)化

于2015—2017年，在上述3個河流斷面進行徑流量的連續(xù)監(jiān)測，并對河水TN和TP兩種污染物指標開展采樣監(jiān)測，頻次為雨季(5—10月)每周1次、非雨季(11月至次年4月)每10天1次。

采用LOADEST模型以不同時間間隔條件下的水質(zhì)數(shù)據(jù)和連續(xù)的徑流量數(shù)據(jù)(日數(shù)據(jù))分別對3個河流斷面的TN和TP日流失通量進行估算。為避免監(jiān)測日期選擇的偶然性，將用于流失通量估算的污染物濃度數(shù)據(jù)起算日期錯開得到時間間隔相同的不同數(shù)據(jù)分組。采用式(10)計算2015—2016年較大與較小時間間隔條件下污染物總流失量估算值的相對誤差：

(10)

式中：RE為相對誤差，%；Li和L10為用LOADEST模型估算得到的時間間隔為i及10 d時某斷面的污染物兩年總流失量，kg；i為時間間隔，d。

2 結(jié)果與討論

2.1 污染物流失通量的時空變化特征

由圖2可見，在時間分布特征上，雨季各斷面的TN、TP日流失通量的平均值明顯高于非雨季，且在雨季的波動范圍比非雨季大。從空間尺度來看，截流斷面的TN、TP日流失通量的平均值明顯低于大興和萬安斷面；大興斷面的TN日流失通量略高于萬安斷面，其TP流失通量平均值與萬安斷面在非雨季基本持平，而在雨季略高于萬安斷面；截流斷面在雨季與非雨季的TN、TP日流失通量波動最大，動態(tài)變化特征最為突出。

2.2 回歸方程的有效性

表1給出了不同監(jiān)測頻次下3個河流斷面的TN和TP流失通量回歸方程的R2。對于水質(zhì)監(jiān)測時間間隔相同的不同數(shù)據(jù)分組應(yīng)選擇不同的最優(yōu)回歸方程?？梢钥闯?，大部分方程的擬合程度均比較好，除大興斷面外，R2在0.75以上。SCR為0.01～0.57，表明部分模型擬合受到共線性影響，且受影響的模型為水質(zhì)監(jiān)測時間間隔為30～40 d的數(shù)據(jù)分組。PPCC均在0.9以上，表明各污染物流失通量回歸方程的參數(shù)估計均可采用AMLE法。

表1 截流、大興及萬安斷面污染物流失通量回歸方程的檢驗結(jié)果1)

2.3 不同監(jiān)測頻次污染物總流失量估算值的相對誤差分析

為驗證水質(zhì)監(jiān)測時間間隔為10 d條件下小流域污染物流失通量LOADEST模型估算值的可靠性，將其與流失通量實測值(污染物濃度實測值與日徑流量實測值的乘積)的時間變化特征對比，結(jié)果表明：污染物日流失通量的估算值與實測值吻合較好(估算值相對誤差在-25%～25%的數(shù)據(jù)占49.59%)，動態(tài)變化特征基本一致(見圖3)。因此，可采用水質(zhì)監(jiān)測時間間隔為10 d的污染物流失通量估算值作為參比值，驗證較大的監(jiān)測時間間隔條件下污染物流失通量估算值的相對誤差。

圖3 監(jiān)測時間間隔為10 d條件下污染物流失通量LOADEST模型估算值與實測值的時間動態(tài)變化Fig.3 Dynamics of the estimated pollutant loss flux in the river using LOADEST with the measured pollutant loss flux at a time interval of 10 d

較大的監(jiān)測時間間隔條件下污染物兩年總流失量估算值(大興斷面部分數(shù)據(jù)由于無法收斂而缺乏參數(shù)估算結(jié)果)的相對誤差見圖4。結(jié)果表明，對于同一個河流斷面，在校準文件中輸入不同時間間隔條件下的污染物濃度及日徑流量實測數(shù)據(jù)時，所得到的污染物總流失量估算值存在差異。相對誤差絕對值最小為0.07%，最大為148.08%，說明不同時間間隔條件下獲得的兩年總流失量估算值之間存在極大差異。從相同時間間隔的不同數(shù)據(jù)分組來看，在40 d間隔條件下，截流斷面TN總流失量估算值的相對誤差變化范圍很大，但萬安斷面TN總流失量估算值卻較為穩(wěn)定。盡管截流斷面在不同監(jiān)測時間間隔條件下獲得的TN流失通量最優(yōu)回歸方程R2均在0.94以上，但其兩年總流失量的模型估算值之間卻存在較大差異，這表明影響截流斷面TN流失的相關(guān)因素及過程隨時間推移呈現(xiàn)隨機、非單調(diào)性變化。

圖4 不同監(jiān)測時間間隔條件下兩年總流失量的LOADEST模型估算結(jié)果及相對誤差Fig.4 The estimated 2-year total exports of pollutants based on monitoring data obtained at different time intervals and their relative errors

從不同的監(jiān)測頻次來看，時間間隔為20 d條件下，不同河流斷面的污染物總流失量估算值的平均相對誤差為10.70%；時間間隔為30 d條件下，平均相對誤差為15.02%；時間間隔為40 d條件下，平均相對誤差為31.76%。這說明水質(zhì)監(jiān)測時間間隔越長，污染物總流失量的估算值與實測值的偏差可能越大。

注：STN、STP分別為單位集水面積TN、TP日流失通量，g/(d·hm2)。圖2 2015—2016年河流斷面的非雨季與雨季單位集水面積實測污染物日流失通量Fig.2 Daily pollutant loss flux per unit drainage area in non-rainy and rainy season during 2015-2016

從不同污染物來看，TN的平均相對誤差為28.93%，TP的平均相對誤差為20.90%。這說明LOADEST模型對TP流失通量的估算更為可靠，但對TN流失通量估算的適用性相對較差，此結(jié)果與KIM等[21]得出的結(jié)果相似。

以25%作為污染物總流失量LOADEST模型估算的允許誤差，由表2可以看出，截流斷面TN和TP流失通量監(jiān)測最大時間間隔為20 d；大興斷面TN監(jiān)測最大時間間隔可設(shè)定為30 d，而TP流失通量則不能用LOADEST模型估算；萬安斷面TN和TP監(jiān)測最大時間間隔為40 d。

表2 各河流斷面模型估算結(jié)果在25%允許誤差范圍內(nèi)的占比

丘陵小流域水文過程一般存在顯著的空間尺度變化規(guī)律，溪溝河流的徑流系數(shù)隨集水面積增加先增大而后趨于恒定，降雨徑流量和污染物流失通量隨時間的變化也逐漸趨于平緩。本研究的數(shù)據(jù)表明，當(dāng)集水面積過小，污染物流失通量的估算值相對于實測值的誤差不穩(wěn)定，因此需要加大監(jiān)測頻率才能增加模型輸出結(jié)果的可靠性。

2.4 模型外推可行性分析

利用2015—2016年污染物濃度數(shù)據(jù)與2015—2017年徑流量數(shù)據(jù)外推得出2017年氮磷流失通量，與用2017年水質(zhì)指標和徑流量實測數(shù)據(jù)計算出的結(jié)果相比較，驗證外推的可行性。由表3可以看出，外推估算值與實測值差異較大，相對誤差最高可達182.70%。因此，外推估算法不可行，對于該丘陵小流域，當(dāng)年徑流量與水質(zhì)實測數(shù)據(jù)都是準確估算當(dāng)年TN、TP總流失量的必要數(shù)據(jù)，氣象水文情況的年際變化會對估算結(jié)果產(chǎn)生不可忽視的影響。

表3 污染物年總流失量的LOADEST模型外推估算值和實測值

2.5 模型估算誤差的來源分析

2015—2016年各河流斷面污染物濃度的變異系數(shù)分析表明，大興斷面TN、TP濃度的變異系數(shù)均為最高，截流斷面次之，萬安斷面的變異系數(shù)最小(見表4)。濃度分布越離散，就需要越多的數(shù)據(jù)才能實現(xiàn)模型的良好擬合。因此，在相同監(jiān)測頻次下LOADEST模型在大興斷面表現(xiàn)較差，而在萬安斷面則表現(xiàn)良好。

表4 河流各斷面污染物濃度的變異程度及集水區(qū)基本屬性

結(jié)合集水區(qū)的基本屬性，可辨識污染物濃度的變異原因。一般而言：集水區(qū)面積越大，緩沖性能越好，濃度變異系數(shù)就越??；坡降越大，越有利于污染物遷移，其濃度變異系數(shù)就越大；農(nóng)地和居民地占比越大，施肥、人為排放的污染物越多，污染物濃度越易受到人類活動的影響。由表4可見，對于面積小、坡降大的集水區(qū)(截流斷面)污染物徑流流失過程(特別是河道徑流)漲退變化快、峰值多；而對于面積大、坡降小的集水區(qū)(萬安斷面)，河道徑流的漲退變化較慢、峰值少。因此，前者需要較高的監(jiān)測頻次才能準確刻畫污染物流失通量隨時間的變化規(guī)律。大興斷面則有所不同，處于河谷平壩區(qū)，受居民地和農(nóng)田兩方面的影響均較大，故而氮磷流失通量的變異系數(shù)較高。

3 結(jié) 論

為探究LOADEST模型在丘陵小流域的適用性，優(yōu)化污染物流失通量估算方法，本研究以四川省鹽亭萬安丘陵小流域為對象，估算基流條件下的氮磷流失通量及其年際變化，比較不同監(jiān)測時間間隔的TN、TP流失通量估算值與實測值的差異，并分析其成因，結(jié)果表明：對于集水面積大、坡降小的河流斷面，允許的最大監(jiān)測時間間隔為20 d；而對于集水面積大、坡降小的河流斷面，允許的最大監(jiān)測時間間隔為40 d；LOADEST模型在TP流失通量估算中的表現(xiàn)優(yōu)于TN；外推估算法不適用于該小流域的污染物流失通量估算。