基于Crystal Ball的泥沙來(lái)源指紋判別及不確定性評(píng)價(jià)

史忠林，陳佳村,2，嚴(yán)冬春，文安邦，龍 翼，吳 毅

(1.中國(guó)科學(xué)院/水利部 成都山地災(zāi)害與環(huán)境研究所，四川 成都 610041；2.中國(guó)科學(xué)院大學(xué)，北京 100049)

流域泥沙來(lái)源的研究由來(lái)已久，方法眾多。傳統(tǒng)的研究方法主要有徑流小區(qū)法、水沙資料分析法和大面積人工或遙感調(diào)查法等[1]。20世紀(jì)70年代以來(lái)，基于泥沙性質(zhì)的指紋識(shí)別技術(shù)興起并快速發(fā)展，目前已成為流域泥沙來(lái)源研究的熱點(diǎn)技術(shù)之一。指紋技術(shù)的基本原理為利用某一種或一組在不同泥沙源地間具有顯著差異的理化性質(zhì)作為指示因子，建立源地物質(zhì)與目標(biāo)泥沙之間的對(duì)應(yīng)關(guān)系，定量研究各源地對(duì)目標(biāo)泥沙的相對(duì)貢獻(xiàn)。指紋技術(shù)的發(fā)展大體經(jīng)歷了單因子識(shí)別與復(fù)合因子識(shí)別兩個(gè)階段。早期研究使用礦物性質(zhì)、地質(zhì)特征、磁性、顆粒顏色及植物孢粉等作為示蹤因子。20世紀(jì)90年代以來(lái)，放射性核素137Cs、210Pb、7Be等相繼用于侵蝕泥沙示蹤研究，指紋因子的種類(lèi)不斷豐富[2]?？紤]到單因子在診斷多物源時(shí)固有的局限性，研究者開(kāi)始嘗試采用多因子組合代替單一指紋因子，復(fù)合指紋技術(shù)應(yīng)運(yùn)而生[3]。復(fù)合指紋法目前在國(guó)外研究較多且較為深入[4]，在國(guó)內(nèi)則相對(duì)較少且以應(yīng)用性探索為主，對(duì)技術(shù)本身的研究相對(duì)欠缺[5-7]。

盡管復(fù)合指紋法歷經(jīng)了20多年發(fā)展已取得重要進(jìn)展，但人們?cè)谘芯恐兄饾u認(rèn)識(shí)到該技術(shù)的應(yīng)用仍存在諸多不確定性，例如對(duì)潛在泥沙源地的預(yù)判、樣品采集的代表性、指紋因子在源匯轉(zhuǎn)化過(guò)程中的穩(wěn)定性、物源與泥沙指紋因子間的可比性及混合模型的選擇與優(yōu)化等。早前的研究多利用經(jīng)統(tǒng)計(jì)檢驗(yàn)篩選得到的指紋因子的平均值或中位數(shù)來(lái)表征不同物源與目標(biāo)泥沙的差異，然而土壤本身固有的空間變異性及對(duì)侵蝕強(qiáng)度響應(yīng)程度的不同，造成利用采集的少量樣品來(lái)代表某一物源或泥沙時(shí)存在較大的不確定性，這就極大地影響到物源貢獻(xiàn)判別結(jié)果的可靠性。為此，近年來(lái)研究者嘗試采用不同的分布函數(shù)(如正態(tài)分布[8]、t分布[9])來(lái)描述物源或泥沙中的指紋因子，然后通過(guò)蒙特卡洛隨機(jī)抽樣方法實(shí)現(xiàn)混合模型的優(yōu)化求解和物源貢獻(xiàn)結(jié)果的不確定性評(píng)價(jià)。

已有的泥沙來(lái)源研究中，蒙特卡洛模擬運(yùn)算主要通過(guò)Matlab[10]、R[11]等軟件平臺(tái)進(jìn)行，要求具有計(jì)算機(jī)基礎(chǔ)，這對(duì)一般研究者來(lái)說(shuō)具有一定挑戰(zhàn)性。為此，本研究介紹一種簡(jiǎn)便快捷的分析軟件Crystal Ball (V.11.1)，并結(jié)合實(shí)例分析其在流域泥沙來(lái)源判別中的應(yīng)用。

1 蒙特卡洛方法及Crystal Ball概述

蒙特卡洛方法是一種以概率統(tǒng)計(jì)為理論基礎(chǔ)的隨機(jī)抽樣模擬方法，通過(guò)估計(jì)不確定變量的概率分布，多次抽樣模擬試驗(yàn)得到結(jié)果變量的概率分布，從而進(jìn)行風(fēng)險(xiǎn)評(píng)價(jià)及敏感度分析。

假定函數(shù)Y=f(x1,x2,…,xn)，其中x1,x2,… ,xn為隨機(jī)變量，其概率分布已知。通過(guò)重復(fù)抽樣取出每一組隨機(jī)變量的值(x1,x2,…,xn)，然后按函數(shù)關(guān)系式確定函數(shù)值Y=f(x1,x2,…,xn)。反復(fù)獨(dú)立抽樣若干次，可得到函數(shù)Y的一批抽樣數(shù)據(jù)Y1,Y2,…,Yk，當(dāng)達(dá)到合理的模擬次數(shù)時(shí)，樣本特征可以近似反映總體的特征。若誤差在可接受范圍內(nèi)，這時(shí)可以得出與實(shí)際情況相近的函數(shù)Y的概率分布及其數(shù)字特征[12]。

Crystal Ball中文名為水晶球，是美國(guó)Oracle公司基于Windows平臺(tái)開(kāi)發(fā)的一款商業(yè)分析和評(píng)估軟件，廣泛應(yīng)用于工程投資和公共安全等領(lǐng)域。Crystal Ball是加載在Excel電子表格上的第三方程序，其界面友好，操作簡(jiǎn)便，分析結(jié)果以圖表和報(bào)告的形式呈現(xiàn)，直觀易懂。用戶可以通過(guò)定義假設(shè)變量的概率分布，用計(jì)算機(jī)產(chǎn)生的隨機(jī)數(shù)進(jìn)行多次模擬計(jì)算，獲得滿足約束條件下決策變量的一組解。

2 泥沙來(lái)源指紋判別分析案例

采集云南省元謀縣涼山鄉(xiāng)某小流域農(nóng)耕地、林地、溝岸等三種潛在泥沙源地的土壤樣品和流域出口處泥沙樣品，測(cè)定樣品理化性質(zhì)，運(yùn)用復(fù)合指紋技術(shù)分析各物源對(duì)流域出口處泥沙的相對(duì)貢獻(xiàn)率。測(cè)試指標(biāo)經(jīng)范圍檢驗(yàn)、Kruskal-Wallis H檢驗(yàn)和多元判別分析后，確定S、Rb、U、Y四種元素構(gòu)成最佳指紋組合，統(tǒng)計(jì)結(jié)果見(jiàn)表1。

2.1 混合模型

利用多元混合模型定量計(jì)算各源地泥沙貢獻(xiàn)百分比，模型結(jié)構(gòu)[3]為

式中：Res為殘差平方和；Cssi為流域出口處泥沙中指紋因子i的濃度；Csi為泥沙源地s中指紋因子i的濃度；Ps為泥沙源地s的泥沙貢獻(xiàn)百分比；m為泥沙源地?cái)?shù)量，m=3；n為指紋因子的數(shù)量，n=4。

模型應(yīng)用需滿足兩個(gè)基本前提，即各源地泥沙貢獻(xiàn)率Ps非負(fù)且總和為1。

2.2 變量定義

Crystal Ball軟件中需要用戶定義的變量有三個(gè)，即假設(shè)變量、決策變量和預(yù)測(cè)變量。假設(shè)變量是Crystal Ball中基本的輸入量，可以是單一的數(shù)值，也可以是服從某一分布的不確定變量。本例中假設(shè)變量為目標(biāo)泥沙中指紋因子的濃度Cssi和各源地指紋因子的濃度Csi。Crystal Ball 11.1分布庫(kù)中提供了正態(tài)分布、三角形分布、均勻分布等21種預(yù)先設(shè)定的分布類(lèi)型和一種自定義分布，用戶可根據(jù)需要自行選擇，也可根據(jù)實(shí)測(cè)值進(jìn)行擬合。但用實(shí)測(cè)值擬合時(shí)，要求實(shí)測(cè)值的數(shù)量至少為15個(gè)。本例假設(shè)目標(biāo)泥沙和源地土壤中各指紋因子濃度均服從正態(tài)分布，用均值和標(biāo)準(zhǔn)差依次定義各指紋因子的分布函數(shù)，并將其下限設(shè)置為0[13]。

決策變量為每個(gè)潛在物源的泥沙貢獻(xiàn)百分比(Ps)，依次定義P1、P2、P3，將其下限和上限分別設(shè)為0和1，變量類(lèi)型設(shè)置為連續(xù)型。

在完成假設(shè)變量和決策變量定義后，需要定義預(yù)測(cè)變量。預(yù)測(cè)變量單元格為包含一個(gè)或多個(gè)假設(shè)變量和決策變量單元格的公式。本例中預(yù)測(cè)變量為各指紋因子濃度的殘差平方和，根據(jù)混合模型的公式編輯單元格。在運(yùn)行模擬時(shí)，預(yù)測(cè)單元格會(huì)自動(dòng)關(guān)聯(lián)模型中的各個(gè)單元格并計(jì)算生成最終結(jié)果。

2.3 模型運(yùn)行

點(diǎn)擊運(yùn)行Crystal Ball軟件中的OptQuest優(yōu)化求解器。由于混合模型的求解要求參與運(yùn)算的所有指紋因子濃度殘差平方和最小，因此將“目標(biāo)”設(shè)置為最小化預(yù)測(cè)單元格的最終值。同時(shí)，在“約束”選項(xiàng)中將決策變量之和指定為1。

優(yōu)化求解時(shí)，需要對(duì)模擬參數(shù)進(jìn)行設(shè)置，主要包括試驗(yàn)次數(shù)和抽樣方法等。軟件中隨機(jī)數(shù)抽樣方法有蒙特卡洛抽樣和拉丁超立方抽樣兩種，其中前者所產(chǎn)生的隨機(jī)數(shù)之間相互獨(dú)立，而后者把假設(shè)變量的概率分布分成幾個(gè)等概率區(qū)間，區(qū)間數(shù)量可根據(jù)需要設(shè)定，模擬試驗(yàn)時(shí)將為每個(gè)區(qū)間產(chǎn)生不同系列的隨機(jī)數(shù)[14]。本例中將模擬次數(shù)設(shè)為1 000次，抽樣方法采用蒙特卡洛隨機(jī)抽樣。

2.4 優(yōu)化求解結(jié)果與不確定性分析

將混合模型模擬運(yùn)算1 000次后，統(tǒng)計(jì)各物源泥沙貢獻(xiàn)率的頻率分布，見(jiàn)圖1。由圖1可知，三種物源的泥沙貢獻(xiàn)頻率呈尖峰分布，集中趨勢(shì)顯著，表明模型計(jì)算結(jié)果的不確定性較小。

指紋因子濃度的殘差平方和(預(yù)測(cè)變量)與各物源對(duì)目標(biāo)泥沙的相對(duì)貢獻(xiàn)百分比(決策變量)統(tǒng)計(jì)結(jié)果列于表2。三種源地中，溝岸對(duì)流域出口處泥沙的貢獻(xiàn)率最大，達(dá)到78.3%±17.1%；其次為農(nóng)耕地，為17.1%±11.6%；林地的貢獻(xiàn)率最小，僅為4.6%±16.1%。據(jù)調(diào)查，泥沙來(lái)源判別結(jié)果符合該流域?qū)嶋H情況。

3 結(jié) 論

運(yùn)用Crystal Ball軟件的優(yōu)化求解功能，結(jié)合蒙特卡洛模擬，成功實(shí)現(xiàn)了某流域泥沙來(lái)源及其相對(duì)貢獻(xiàn)率的定量判別，模型運(yùn)算結(jié)果與流域?qū)嶋H情況吻合。Crystal Ball通過(guò)對(duì)指紋因子的分布擬合，極大地降低了指紋因子單一取值帶來(lái)的模型求解結(jié)果的不確定性。同時(shí)，模型輸出結(jié)果給出了不同物源對(duì)流域泥沙貢獻(xiàn)率的平均值和標(biāo)準(zhǔn)差，相比傳統(tǒng)的Excel規(guī)劃求解只給出唯一解而言，提高了判別結(jié)果的科學(xué)性和合理性。此外，Crystal Ball操作簡(jiǎn)便，不需要復(fù)雜的計(jì)算機(jī)語(yǔ)言，具有較強(qiáng)的實(shí)用性，未來(lái)將在流域泥沙來(lái)源研究中得到更廣泛的應(yīng)用[15]。