基于二項分布檢驗法的水質(zhì)達標評價(jià)方法研究

　1 引言

　　水質(zhì)達標評價(jià)(Water Quality Compliance Assessment)是根據一定時(shí)段內的水質(zhì)監測數據，評價(jià)水體狀態(tài)是否符合其功能對應的水質(zhì)標準的決策過(guò)程，是流域水污染防治和水質(zhì)管理的重要科學(xué)依據.區別于水質(zhì)評價(jià)研究中對水質(zhì)綜合指標、時(shí)間序列趨勢、季節性特征等的關(guān)注，水質(zhì)達標評價(jià)的關(guān)注點(diǎn)是水質(zhì)是否達標，并可以由此評估已有的流域治理措施效果、制定新的治理規劃(如TMDL).在流域水污染防治資金緊缺的條件下，準確識別水體是否達標至關(guān)重要：如果評價(jià)方法過(guò)于保守，會(huì )導致水生態(tài)系統健康受損;若評價(jià)方法過(guò)于嚴格，會(huì )導致水體過(guò)于保護，水體自?xún)裟芰�無(wú)法得到有效利用及高的治理資金投入.科學(xué)合理的水質(zhì)達標評價(jià)方法是進(jìn)行水質(zhì)達標評價(jià)的基礎.我國地表水環(huán)境質(zhì)量評價(jià)采用單因子評價(jià)方法，對斷面時(shí)空監測數據求算術(shù)平均值判斷其水質(zhì)類(lèi)別，根據斷面水質(zhì)類(lèi)別比例對流域水質(zhì)進(jìn)行分級.美國清潔水法303(d)條款規定對受損水體實(shí)施TMDL計劃，EPA導則推薦使用監測數據的90%分位數評價(jià)水體是否受損，即監測數據的超標個(gè)數不能超過(guò)10%.在歐盟，采用每年12個(gè)月監測數據的算術(shù)平均值作為評價(jià)依據，少數成員國將監測數據的變異性考慮在內，給出統計置信區間.

　　水質(zhì)指標監測值存在不確定性已經(jīng)被廣泛地認識和研究.我國水質(zhì)達標評價(jià)方法未考慮水質(zhì)指標的不確定性，美國和歐盟已經(jīng)將不確定性納入到水質(zhì)達標評價(jià)中.水質(zhì)指標的不確定性包括固有的變異性(如監測數據的季節性差異)和可減小的不確定性(如測量誤差)2部分.通過(guò)合理選擇監測斷面、監測指標以及質(zhì)量保證和質(zhì)量控制(QA/QC)可以減小不確定性，卻不能完全消除.忽視不確定性可能會(huì )帶來(lái)巨大的決策風(fēng)險，科學(xué)合理的水質(zhì)達標評價(jià)方法必須能夠表征監測數據的不確定性，將由水質(zhì)指標不確定性帶來(lái)的決策風(fēng)險定量化，為決策提供有用信息.

　　在統計學(xué)視角下，將水質(zhì)指標視為隨機變量，可有效地表征不確定性.水質(zhì)指標具有時(shí)空分異性，在同一時(shí)間獲得空間所有點(diǎn)位的監測數據或者在同一點(diǎn)位進(jìn)行連續的監測是不可能的，只能通過(guò)有限的監測數據對水質(zhì)進(jìn)行評價(jià)，即根據樣本數據估計總體特征.統計學(xué)假設檢驗方法已經(jīng)廣泛應用于環(huán)境系統中，假設檢驗的2類(lèi)錯誤可以定量表征水質(zhì)指標的不確定性帶來(lái)的決策風(fēng)險.針對美國EPA提出的最大超標率為10%，Smith等和McBride等認為10%應該指10%的時(shí)間，即總體分布的90%分位數，并采用二項分布檢驗法(Binomial Test)進(jìn)行了水質(zhì)達標評價(jià).該方法是一種針對二元數據的非參數方法，通過(guò)假設檢驗給出在一定允許超標率和樣本容量時(shí)的最大超標個(gè)數，并與實(shí)際監測得到的超標個(gè)數對比，評價(jià)水質(zhì)是否超標.盡管研究者指出二項分布檢驗法存在信息損失的問(wèn)題，并提出了參數檢驗方法，然而McBride的研究表明，在大部分情況下，二項分布檢驗法由于不受限于原始數據的分布而更具有魯棒性.本文以基于二項分布檢驗法的水質(zhì)達標評價(jià)方法為基礎，分析其決策過(guò)程及決策風(fēng)險的定量表征，并以海河流域為例，對比不同評價(jià)方法的差異性，分析其適用性和靈活性.

　　2 研究方法

　　2.1 基于二項分布檢驗法的水質(zhì)達標評價(jià)方法

　　在二項分布檢驗法之前，需將監測數據轉化為0~1變量，方法為：對于某一特定監測數據x，將其與水質(zhì)標準θ對比，對比后的結果記為y，如果x>θ，即超過(guò)水質(zhì)標準，記y=1;如果x≤θ，記y=0.令p(0<p<1)為成功概率，表示x>θ的概率，則y服從成功概率為p的0~1分布，其分布律為：

(1)

　　假設收集到的N個(gè)監測數據用于評價(jià)水質(zhì)是否達標，記隨機變量X=x1，x2，...，xN，將X中的每個(gè)元素轉化為0~1變量，便得到隨機變量Y=y1，y2，...，yN，則Y服從參數為(N，p)的二項分布，記為Y~B(N，p).令n為Y成功的總次數，即為N個(gè)樣品中超標總個(gè)數，其分布律為：

(2)

　　二項分布的分布律表征了水質(zhì)變量轉化為二元數據后超標總個(gè)數不確定性.欲檢驗總體分布的90%分位數是否達標，即超標率是否小于10%，可令原假設為水質(zhì)達標，即H0:p≤0.1，備擇假設為H1:p>0.1.原假設為真時(shí)，可得p=0.1時(shí)超標個(gè)數n的分布函數F;實(shí)際監測數據中超標個(gè)數nc不應太大，在給定的顯著(zhù)性水平α下，應有F(nc)=P(n≤nc)≤1-α;否則拒絕原假設，判定為水質(zhì)超標.

　　2.2 棄真錯誤和取偽錯誤

　　在經(jīng)典統計學(xué)的假設檢驗過(guò)程中，不可避免地要討論2類(lèi)錯誤：棄真錯誤和取偽錯誤.在水質(zhì)達標評價(jià)中，這2類(lèi)錯誤是由于水質(zhì)變量的不確定性導致的.棄真錯誤為當原假設為真時(shí)，放棄原假設而犯的錯誤(falsely inferring a breach of standard)，棄真錯誤的概率記為α;取偽錯誤為當原假設為假時(shí)，相信原假設而犯的錯誤(falsely inferring compliance)，棄偽錯誤的概率記為β.對于二項分布，$\alpha =P\left( n>{{n}_{s}} \right)=\sum\limits_{j={{n}_{s}}+1}^{N}{C_{N}^{j}{{P}^{j}}}{{\left( 1-P \right)}^{N-j}}$.計算β時(shí)，需要給定效應值η(備擇假設和原假設之間的差距)，二項分布檢驗法一般取η=0.15，即p′=p+η=0.25來(lái)計算β，計算公式為β=G-1(ns)，其中G(n)=$\sum\limits_{k=1}^{N}{C_{N}^{k}{{P}^{'k}}}{{\left( 1-P' \right)}^{N-k}}$，為成功概率為p′的二項分布的分布函數，G-1(n)為G(n)的反函數.

　　按照上述α和β的計算公式，可以得到N一定時(shí)，以不同的ns(0≤ns≤N)作為決策準則(即nc>ns判定為超標，nc≤ns判定為達標)時(shí)，α和β的變化情況.由圖 1可知N一定時(shí)，對于特定的ns存在唯一的(α，β)組合;α隨ns減小，β隨ns增加.圖中豎直虛線(xiàn)假設變量服從對稱(chēng)分布，算術(shù)平均值等于50%分位數，則平均值法等價(jià)于最大允許超標率為50%.圖中豎直實(shí)線(xiàn)代表采用10%超標樣品數目(以下稱(chēng)為“超標比例法”)得到的ns，當N=12時(shí)，ns=1，α=0.341;當N=50時(shí)，ns=5，α=0.384.可見(jiàn)對于相同樣本容量，平均值法棄偽錯誤的概率很大，超標比例法棄真錯誤的概率較大，而二項分布檢驗法則可根據決策者需要選擇不同的(α，β)組合.相對于選擇特定的顯著(zhù)性水平值做決策，(α，β)組合能量化決策準則的風(fēng)險，為決策者提供更大的決策空間.α和β隨N增加均有減小的總體趨勢，參見(jiàn)Smith等.

　　圖 1N=12(a)和50(b)時(shí)2類(lèi)錯誤概率隨最大允許超標個(gè)數的變化

　　2.3 決策過(guò)程

　　基于二項分布檢驗法的水質(zhì)達標評價(jià)的決策過(guò)程包括2個(gè)環(huán)節：首先確定ns，然后根據監測樣品的超標個(gè)數nc判斷水質(zhì)是否達標.一般而言，決策過(guò)程應該使得α和β最小.在N一定時(shí)，由于α、β隨ns單調性相反(圖 1)，對于特定的ns不能同時(shí)使α和β最小.作為決策者，在選擇ns時(shí)：一方面，需要根據風(fēng)險偏好和科學(xué)判斷，選擇最大允許超標率p;另一方面，需要權衡各相關(guān)方的利益，如：作為受體人群和生態(tài)系統的及作為治理投資者的政府或企業(yè).在同等條件下，ns隨著(zhù)p的增大而不嚴格單調遞增;在p和N一定時(shí)，根據決策者的利益權衡，可以得到不同的ns.

　　下面介紹3種典型的(α，β)組合選擇策略：①保證α在一定的低水平，例如α≤0.05，使得β最小;②保證β在一定的低水平，例如β≤0.05，使得α最小;③控制α≈β，即minα-β.圖 2給出了當p=0.1，N為4~100時(shí)，這3種策略得到的最大允許超標個(gè)數ns和最大允許超標率η=nsN(散點(diǎn)圖，連續線(xiàn)為點(diǎn)過(guò)于密集的結果).可見(jiàn)當N較小時(shí)，控制α≤0.05得到的ns最大，控制β≤0.05得到的ns最小，且當N≤10時(shí)，無(wú)論如何選擇ns均不能使β≤0.05;隨著(zhù)N變大，控制α≤0.05得到的ns最小，控制β≤0.05得到的ns最大;策略③的結果介于①和②之間.總體上，控制α≤0.05得到的η呈遞減趨勢，由于假設檢驗是p=0.1時(shí)二項分布的95%分位數，因此η恒大于0.1(圖 2d虛線(xiàn));控制β≤0.05得到的η呈遞增趨勢，由于假設檢驗是p=0.25時(shí)二項分布的5%分位數，因此η恒小于0.25(圖 2e虛線(xiàn)).綜上，基于二項分布檢驗法的水質(zhì)達標評價(jià)決策過(guò)程可表示為圖 3.決策者可以根據風(fēng)險偏好和利益權衡，選擇p值和(α，β)組合，進(jìn)而確定N一定時(shí)的ns，結合水質(zhì)監測數據中的實(shí)際超標個(gè)數nc，對水體的達標狀況進(jìn)行判斷;通過(guò)p和(α，β)可將決策的各種風(fēng)險定量化.

　　圖 2(Fig. 2)

　　圖 2最大允許超標個(gè)數和最大允許超標率隨樣本容量的變化

　　圖 3基于二項分布檢驗法的水質(zhì)達標評價(jià)決策過(guò)程

　　3 案例研究(Case study)3.1 研究對象和方法

　　選擇海河流域的CODMn作為達標評價(jià)的對象，數據來(lái)源為環(huán)保部網(wǎng)站，收集2014年7月22日-2015年7月14日的周監測數據，共有北京密云古北口(51)、北京門(mén)頭溝沿河城(51)、天津三岔口(51)、天津果河橋(51)、河北張家口八號橋(45)、河北石家莊崗南水庫(51)、山東聊城秤鉤灣(42)等7個(gè)監測斷面，括號中為不同監測斷面的有效監測樣品數.

　　當N=42、45、51時(shí)，采用平均值法、二項分布檢驗法的3種策略、超標比例法確定α、β、ns(對于平均值法，求各個(gè)樣本中大于平均值的樣品個(gè)數，作為平均值法的ns)，對比監測數據與地表水環(huán)境質(zhì)量標準中各水質(zhì)類(lèi)別的濃度限值θi(i=1，2，...，5)，得到超過(guò)θi的樣品個(gè)數nci，水質(zhì)類(lèi)別判定遵循如下過(guò)程：當i=1時(shí)，若nc1≤ns則水質(zhì)類(lèi)別判定為I類(lèi);當i=5時(shí)，若nc5>ns則判定為劣V類(lèi);若nci≤ns<nc(i-1)則判定為i類(lèi).最后，根據地表水環(huán)境質(zhì)量評價(jià)辦法(試行)對流域水質(zhì)進(jìn)行分級.本文的所有統計方法和運算均在r軟件(version:r-3.2.3，https: p="" )中實(shí)現.<="" www.r-project.org="">

　　3.2 水質(zhì)評價(jià)結果

　　各種方法對應的α、β、ns結果見(jiàn)表 1，各種方法的嚴格程度為：超標比例法>B(2)>B(3)>B(1)>平均值法(“a>b”表示a嚴格于b，B(i)表示二項分布檢驗法的第i種策略);采用平均值法β接近于1，超標比例法得到的α很大(>0.4)，二項分布檢驗法的3種策略得到的2類(lèi)錯誤之和小于0.3，平均在0.2左右，且均處于相對較低水平.水質(zhì)評價(jià)結果見(jiàn)表 2，在7個(gè)斷面中，石家莊崗南水庫水質(zhì)最好，聊城秤鉤灣水質(zhì)最差;石家莊崗南水庫和門(mén)頭溝沿河城2個(gè)斷面采用各種方法的評價(jià)結果相同，另外5個(gè)斷面有差異.根據斷面的水質(zhì)類(lèi)別得到的流域分級結果也有差異.斷面水質(zhì)類(lèi)別的評價(jià)結果，除了與評價(jià)方法的嚴格程度有關(guān)，也與水質(zhì)監測數據有關(guān).采用相對嚴格方法評價(jià)達標的斷面，在不嚴格方法下自然達標;采用相對不嚴格方法評價(jià)達標的斷面，在嚴格的方法下不一定達標，例如密云古北口采用平均值法評價(jià)為I類(lèi)水質(zhì)，而采用其它方法評價(jià)時(shí)卻為II類(lèi)水質(zhì).總體上看，評價(jià)結果符合嚴格方法得到的水質(zhì)類(lèi)別較低的規律.

　　表 1 3種樣本容量對應的各種評價(jià)方法的2類(lèi)錯誤概率和最大允許超標個(gè)數

　　表 2 海河流域CODMn達標評價(jià)結果

　　綜合表 1和表 2，對于二項分布檢驗法的3種方法，控制β≤0.05有利于保護受體的利益，對水體水質(zhì)要求嚴格，流域分級為“輕度污染”;控制α≤0.05有利于減少治理投資，對流域分級為“良好”;二項檢驗的3種策略的評價(jià)結果對比，說(shuō)明在采用二項分布檢驗法進(jìn)行水質(zhì)達標評價(jià)時(shí)，即使對于完全相同的水質(zhì)監測數據，決策時(shí)傾向于不同利益群體，可能得出不同的水質(zhì)達標評價(jià)結果.二項分布檢驗法為水質(zhì)達標評價(jià)的決策提供了較為廣闊的決策空間.

　　案例研究中，取最大超標率p=0.1，采用二項分布檢驗法的3種典型策略進(jìn)行水質(zhì)達標評價(jià);在實(shí)踐中，決策過(guò)程有更大的靈活性：對于p，決策者可以根據風(fēng)險偏好選擇不同值.對于策略①和②，決策可設定不同的錯誤概率限值，例如令α≤0.1使得β最小，再比如，可同時(shí)限定α，β的最大值，例如令α≤0.2同時(shí)β≤0.1，即可得到多組(α，β)組合;對于策略③，可根據需要設定限制條件，例如出于保護受體的目的令α<β，或者給定最大可接受的差值δ，令α-β≤δ，亦可得到不同的(α，β)組合.對得到的多組(α，β)組合，決策者可結合其它約束條件進(jìn)行優(yōu)選確定最終的ns.

　　4 討論

　　4.1 先驗知識與假設檢驗

　　研究中采用的原假設是水質(zhì)達標，α代表水質(zhì)達標而判定為不達標的概率，β代表水質(zhì)未達標而判定為達標的概率，二項分布檢驗法的策略①實(shí)質(zhì)上是求成功概率為p的二項分布的95%分位數;Gibbons指出如果有很強的先驗知識認為水體超標，則假設檢驗的原假設應為水質(zhì)未達標，此時(shí)α和β的含義互換，而二項分布檢驗法的策略①變?yōu)榍蟪晒Ω怕蕿閜的二項分布的5%分位數.由計算可得，當先驗知識認為水質(zhì)超標時(shí)，得到的ns較小.如當效應值為0.15，樣本容量N=50，最大超標率p=0.25時(shí)，按策略①，原假設為水質(zhì)達標時(shí)ns=18，原假設為水質(zhì)超標時(shí)ns=7，則對于超標個(gè)數在8~18的樣本按照2種原假設得到的評價(jià)結果不同，例如當nc=10時(shí)，若原假設為水質(zhì)達標，則判定為達標，若原假設為水質(zhì)超標，則判定為超標.可見(jiàn)，水質(zhì)現狀的先驗知識對評價(jià)結果具有顯著(zhù)影響，在先驗知識充足時(shí)，水質(zhì)達標評價(jià)具有“寬進(jìn)嚴出”的規律，即將先驗認為達標的水體列為不達標狀態(tài)時(shí)對超標個(gè)數的要求較為保守，而將先驗認為未達標水體列為達標狀態(tài)時(shí)對超標個(gè)數的要求較為嚴格.β值與效應值的選取有關(guān)，本文參照文獻做法取效應值為0.15，即備擇假設的p=0.25.在同等條件下，盡管根據不同效應值可以得到不同的β，但是其隨N和ns的變化規律一致.假設檢驗的效果隨樣本容量的增加而增大，對海河流域的研究顯示，采用二項分布檢驗法可以控制各個(gè)斷面的總錯誤概率平均值在0.2左右，評價(jià)結果具有很高的可信度.

　　4.2 最大允許超標概率的選擇

　　對于研究中最大超標率設定的目的，研究者給出了不同的解釋?zhuān)篠mith等、McBride等認為水質(zhì)變量不確定性造成的某次監測數據超標是難以避免的，因此對全部水質(zhì)監測數據給出一個(gè)可以接受的超標率(10%)是合理的;此時(shí)10%可理解為“容忍度”.Gibbons則認為對90%分位數進(jìn)行假設檢驗是為了保證總體的均值小于標準值;此時(shí)，樣本均值與90%分位數之間的差異可理解為“安全因子”.Smith等針對的是一類(lèi)分位數標準(percentile-based standards)，而Gibbons針對的則是平均值標準，平均值代表流域水質(zhì)變量波動(dòng)性的平均水平，與室內實(shí)驗或野外試驗的穩定濃度相對應.本文建議，在使用二項分布檢驗法進(jìn)行水質(zhì)達標評價(jià)時(shí)，需要首先區分評價(jià)變量的標準類(lèi)型，如果是分位數標準，則根據要求選擇較小的p值(例如0.1)是合理的;如果是平均值標準，過(guò)大的安全因子會(huì )導致水質(zhì)目標難以實(shí)現，需要選擇較大的p值(例如0.25).

　　結合4.1和4.2小節，最大允許超標概率和假設檢驗的2類(lèi)錯誤率均是由于水質(zhì)監測數據的不確定性導致的.最大允許超標概率體現了由于監測數據波動(dòng)性，而對其采用的“容忍度”或者“安全因子”，決策者對最大允許超標概率的選擇體現了風(fēng)險偏好;假設檢驗的2類(lèi)錯誤體現了由于水質(zhì)變量的不確定性，使用樣本估計總體特征的可信程度，決策者對2類(lèi)錯誤的選擇體現了對不同群體的利益權衡.通過(guò)這3個(gè)參數，基于二項分布檢驗法的水質(zhì)達標評價(jià)方法將水質(zhì)監測數據的不確定性對決策過(guò)程的影響定量化，對決策的各種風(fēng)險進(jìn)行了定量表征，使得決策更具科學(xué)性.

　　5 結論

　　1) 水質(zhì)變量存在不確定性，平均值法忽視了水質(zhì)變量的不確定性，會(huì )造成很大的取偽錯誤概率，而超標比例法會(huì )導致很大的棄真錯誤概率;基于二項分布檢驗法的水質(zhì)達標評價(jià)方法可以將由水質(zhì)變量的不確定性導致的決策風(fēng)險定量化，為決策者提供科學(xué)依據.

　　2) 本文在基于二項檢驗的水質(zhì)達標評價(jià)方法的基礎上，納入決策者的風(fēng)險偏好和對利益的權衡，提出了水質(zhì)達標評價(jià)的決策框架.二項分布檢驗法對原始數據的分布沒(méi)有要求，是一種方便快捷的非參數方法.本文提出的決策框架，可為決策者提供較大的決策空間，決策者可以根據風(fēng)險偏好和利益權衡選擇或限制最大超標率、效應值、棄真錯誤概率和取偽錯誤概率組合，確定最大允許超標個(gè)數，進(jìn)行水質(zhì)達標評價(jià).

　　3) 對海河流域的水質(zhì)評價(jià)結果表明，采用基于二項分布檢驗法的水質(zhì)達標評價(jià)方法，相對于平均值法和超標比例法，其評價(jià)結果具有很高的可信度，不同的利益權衡會(huì )得到不同的評價(jià)結果，驗證了二項分布檢驗法的適用性和靈活性.