河相關(guān)系的隨機(jī)微分方程建模與研究

宋曉龍，鐘德鈺，王光謙

（1.清華大學(xué) 水沙科學(xué)與水利水電工程國(guó)家重點(diǎn)實(shí)驗(yàn)室，北京 100084；2.天津大學(xué) 水利工程仿真與安全國(guó)家重點(diǎn)實(shí)驗(yàn)室，天津 300072）

1 研究背景

河流生態(tài)系統(tǒng)，是陸地與海洋聯(lián)系的紐帶，在生物圈的物質(zhì)循環(huán)中起著主要作用。隨著近代以來全球化的推進(jìn)，借助基礎(chǔ)河流理論和數(shù)字分析手段的發(fā)展，河流資源開發(fā)力度逐步加大，同時(shí)與河流保護(hù)的矛盾問題也日益凸顯。在制定河流保護(hù)和水資源開發(fā)計(jì)劃之前，應(yīng)該最大限度地掌握河流相關(guān)信息，以做到合理科學(xué)地規(guī)劃河流。在中國(guó)等發(fā)展中國(guó)家，河流動(dòng)力學(xué)和河流管理始終是相關(guān)基礎(chǔ)研究領(lǐng)域的重要課題，對(duì)實(shí)現(xiàn)資源的可持續(xù)利用意義重大。不過研究目標(biāo)需要分階段逐步實(shí)現(xiàn)，河流系統(tǒng)的非線性和突變性等特征始終是擺在專家學(xué)者面前的一道難題。

在很長(zhǎng)一段時(shí)間內(nèi)，線性模型在研究河灣演化規(guī)律中扮演重要的角色，揭示了一些河灣現(xiàn)象的主要特征［1-4］，例如:河道幾何特征如河寬和彎曲率等的局部波動(dòng)可以引起河道中縱向沙壩和橫向上河道平面的大尺度變化；小尺度自由沙壩在小寬深比條件下傾向于向下游移動(dòng)，而在大寬深比條件下有可能向上游移動(dòng)；在小擾動(dòng)下，河流可在長(zhǎng)期的發(fā)展過程中達(dá)到準(zhǔn)穩(wěn)定狀態(tài)，等等。雖然各種研究層出不窮，但遠(yuǎn)遠(yuǎn)稱不上完美。事實(shí)上，在洪水影響下，河流通常呈現(xiàn)出季節(jié)性和不規(guī)則波動(dòng)。在極端天氣增多增強(qiáng)的今天，突發(fā)強(qiáng)降雨等事件無疑增加了我們的疑慮:是否河流系統(tǒng)可以與時(shí)時(shí)變化的水沙輸入保持平衡，達(dá)到所謂的“動(dòng)態(tài)平衡狀態(tài)”；用一個(gè)單獨(dú)的造床流量是否可以代表影響河道形態(tài)的各種非恒定因素［5］；是否可以另外認(rèn)為河道形態(tài)是一種歷史擾動(dòng)的綜合反映，而不是同外界時(shí)時(shí)保持平衡的一種準(zhǔn)穩(wěn)定狀態(tài)［6］；文獻(xiàn)［7-9］中有一些對(duì)氣候變化下的水文響應(yīng)特征研究，通過自然實(shí)驗(yàn)或衛(wèi)星遙感分析了極端降雨對(duì)河道的潛在影響，但至今缺少對(duì)隨機(jī)擾動(dòng)下的河流特征關(guān)系的機(jī)制性研究。

近幾十年來，眾多物理學(xué)家轉(zhuǎn)向現(xiàn)代系統(tǒng)理論拓展本領(lǐng)域的研究。其方法與以往經(jīng)典科學(xué)的顯著不同之處是，把自己研究的對(duì)象不再看作某個(gè)單獨(dú)事物，而是看作一個(gè)整體系統(tǒng)，要求做整體綜合性的研究和把握。特別是，表示線性對(duì)象的輸入與輸出間關(guān)系的傳遞函數(shù)法獲得了廣泛應(yīng)用。我們關(guān)注的自然河流系統(tǒng)，具有高度的規(guī)律性和交互性，其非平衡態(tài)特征正獲得越來越多的關(guān)注，而簡(jiǎn)單的確定性方程并不能滿足相關(guān)研究的要求。1980年代以來，越來越多的研究基于擴(kuò)展的傳遞函數(shù)進(jìn)行河流演化過程的追蹤和識(shí)別，典型案例包括:用傳遞函數(shù)噪聲模型預(yù)測(cè)產(chǎn)匯流過程［10］，用隨機(jī)微分方程代替確定性方程來數(shù)值模擬水質(zhì)量因子的動(dòng)態(tài)概率演化過程［11］，研究河岸帶寬度和植被密度隨流量變化的隨機(jī)響應(yīng)［12］，研究土壤水平衡動(dòng)態(tài)模型在隨機(jī)降雨下的隨機(jī)變化［13］，研究地表徑流中懸浮顆粒的隨機(jī)運(yùn)動(dòng)［14］等等。這些方法為河流特征研究和流域管理提供了新的思路:將水流輸入項(xiàng)視為隨機(jī)源，通過概率化的研究來有效處理氣候突變等外界影響下的河流響應(yīng)問題。

河相關(guān)系是一種描述河流形態(tài)要素特征關(guān)系的經(jīng)典工具，文獻(xiàn)［15-18］中有大量的關(guān)于其非線性特征的研究。然而問題是，大多研究只是將焦點(diǎn)聚焦在利用實(shí)測(cè)數(shù)據(jù)來求解其冪函數(shù)的系數(shù)和指數(shù)本身，而很少考慮過該函數(shù)關(guān)系在隨機(jī)事件的影響下可能發(fā)生怎樣的變化。因此，本文將利用隨機(jī)微分方程的方法，對(duì)河相關(guān)系的隨機(jī)擾動(dòng)提出新的見解，推動(dòng)相關(guān)基礎(chǔ)理論的發(fā)展，有效提高河流管理水平。

2 河相關(guān)系的隨機(jī)微分方程模型

2.1 理論模型在河床演變學(xué)中，平灘水位是指河流灘地（唇）高程齊平，河道發(fā)生漫灘的臨界水位。水位上升到平灘水位時(shí)，由于相應(yīng)水流的流速大，輸沙能力高，造床能力強(qiáng)，因此對(duì)應(yīng)的所謂“平灘流量”成為反映河道排洪輸沙能力和設(shè)計(jì)穩(wěn)定河道幾何形態(tài)的重要指標(biāo)。經(jīng)典研究普遍認(rèn)為平灘流量可以由一個(gè)具體數(shù)值代表，無論來水來沙是否發(fā)生趨勢(shì)性變化，都會(huì)隨著河床的沖淤演變而隨時(shí)間變化。如吳保生［19］基于河流自動(dòng)調(diào)整的基本原理，根據(jù)河床受到擾動(dòng)后調(diào)整速率與其當(dāng)前狀態(tài)和平衡狀態(tài)之間差值成正比的基本規(guī)律，建立了一種平灘流量的滯后響應(yīng)模型:

式中:Qt為t時(shí)刻的平灘流量；β為平灘流量的變化速率，是一個(gè)與水流能量大小和河床邊界可動(dòng)性有關(guān)的參數(shù)；Qe為相對(duì)某一給定水沙條件下，河床調(diào)整達(dá)到相對(duì)平衡時(shí)的平灘流量，可由式（2）表示:

式中:Qf為汛期平均流量；ξf=CtfQf，稱作來沙系數(shù)；Ctf為汛期平均懸移質(zhì)含沙量；K、b、c分別為待定系數(shù)和指數(shù)。

在一定的來沙水沙條件下，河床將不斷調(diào)整其比降和斷面形態(tài)，力求挾沙能力與上游來水來沙相適應(yīng)。在此調(diào)整過程中河床所能達(dá)到的最適宜狀態(tài)，或者說平衡狀態(tài)，我們稱為河相關(guān)系。它是通過某一特征流量，如平灘流量式（1），將河床不同斷面的比降、河寬、水深、流速等特征聯(lián)系在一起，用以描述沖積河床橫斷面形態(tài)以及河床的總輪廓［20］。如馬卡維耶夫公式和Leopold［21］等都是從穩(wěn)定斷面形態(tài)開始，得到一組描述天然河流河床形態(tài)的經(jīng)驗(yàn)關(guān)系式:

式中:St、Bt、Dt、Ut分別為平灘流量下t時(shí)刻某河段平均比降及入口斷面的河寬、平均水深、平均流速；下標(biāo)不同的4種α和m為待定系數(shù)和指數(shù)。

根據(jù)式（3）的冪函數(shù)形式，我們可以用一概化變量X來代表河相關(guān)系特征參量:

冪函數(shù)形式的河相關(guān)系滿足非線性系統(tǒng)進(jìn)行線性化處理的條件，可對(duì)式（4）求導(dǎo)得:

將式（1）、（2）帶入式（5），并轉(zhuǎn)換為微分形式，有:

如果取平灘流量滯后響應(yīng)模型的多步遞推模式［19］來表示Qt的通解:

式中:Q0為t=0時(shí)Qt的值，即初始平灘流量；Δt為時(shí)段長(zhǎng)度；n為迭代時(shí)段數(shù)；i為時(shí)段編號(hào)。

根據(jù)微分方程式（6）得到一組描述河相關(guān)系的確定性的迭代解。

河流系統(tǒng)是一個(gè)自然結(jié)構(gòu)、生態(tài)環(huán)境和經(jīng)濟(jì)社會(huì)相互耦合的開放系統(tǒng)，由于水體的流動(dòng)性，系統(tǒng)與外界不斷進(jìn)行物質(zhì)和能量的交換以及信息的傳遞，從而表現(xiàn)出明顯的耗散性和協(xié)同性［22］。河流系統(tǒng)不能完全消除外界影響而獨(dú)立存在，因此普通微分方程在真實(shí)反映河流演化過程的效果要大打折扣。考慮到河床形態(tài)因子受各種因素的影響，比如氣候或人為環(huán)境擾動(dòng)，引起水沙條件、河床邊界條件等產(chǎn)生不確定性，因此可參照金融領(lǐng)域的Black-Shloes模型等經(jīng)典的隨機(jī)擴(kuò)散模型［23-25］，假設(shè)河相關(guān)系的動(dòng)態(tài)行為具有隨機(jī)性，在穩(wěn)態(tài)模型式（6）的基礎(chǔ)上增加噪聲模型來對(duì)其進(jìn)行分析:

經(jīng)典隨機(jī)微分方程的噪聲項(xiàng)通常由高斯白噪聲表示。然而其缺點(diǎn)是高斯白噪聲只能表現(xiàn)連續(xù)性強(qiáng)、較穩(wěn)定的隨機(jī)擾動(dòng)，而不能控制突發(fā)情況的發(fā)生。實(shí)際的白噪聲就像白色光包含所有可見光譜一樣，也包含各種類型的噪聲，依據(jù)概率分布形式不同，可主要分為高斯白噪聲和非高斯白噪聲。其中，常用泊松噪聲（也稱泊松跳）來模擬的非高斯噪聲項(xiàng)，起著重要的傳遞突發(fā)信息的作用。其他領(lǐng)域的隨機(jī)研究中已經(jīng)發(fā)現(xiàn)了泊松跳的重要價(jià)值，并進(jìn)行了大量的推廣應(yīng)用［26-28］。此外，傳統(tǒng)的線性思維模式可能并不符合現(xiàn)實(shí)情況:天然河流系統(tǒng)就像分子生物系統(tǒng)，并不是純粹的機(jī)械系統(tǒng)，而是一個(gè)生命系統(tǒng)，其可表現(xiàn)出某些自發(fā)性和特殊的復(fù)雜性。研究表明，真實(shí)的系統(tǒng)輸入具有長(zhǎng)相關(guān)和自相似性特征，而分子布朗運(yùn)動(dòng)是描述這種有機(jī)隨機(jī)過程的典型工具。因此，我們有必要將分形理論也應(yīng)用進(jìn)隨機(jī)微分方程以更好地刻畫河相關(guān)系的演變特征。考慮到環(huán)境變化可能給河相關(guān)系帶來較穩(wěn)定的擾動(dòng)，或類似脈沖的隨機(jī)沖擊，而這些隨機(jī)過程不管是發(fā)生頻率、發(fā)生時(shí)間或發(fā)生強(qiáng)度都是隨機(jī)的，因此筆者使用三種流行的噪聲模型，包括高斯白噪聲、泊松噪聲和分?jǐn)?shù)白噪聲，進(jìn)行三種組合，評(píng)估隨機(jī)微分方程（8）描述河相關(guān)系隨機(jī)演變特征的適用性。詳情如下:

2.1.1 簡(jiǎn)單的高斯白噪聲模型

式中:{Wt:0≤t≤T}為標(biāo)準(zhǔn)維納過程；m和σ分別為漂移項(xiàng)和擴(kuò)散項(xiàng)參數(shù)；X0為初值。

方程式（9）的解可表示為:

2.1.2 高斯白噪聲+泊松噪聲組合模型

式中:Vj為跳躍強(qiáng)度；Nj(λj)為兩種相互獨(dú)立的泊松過程；λj為跳躍事件發(fā)生頻率；j=u*，d*分別表示上跳和下跳；跳躍強(qiáng)度的對(duì)數(shù)分布規(guī)律(Y=ln(V))可由下面的雙指數(shù)分布函數(shù)表示:

式中:ηu＞1，

這里用一個(gè)混合方法來模擬跳躍發(fā)生的時(shí)間點(diǎn)(τ1，τ2，…) 。其中，時(shí)間差Ri+1=τi+1-τi可根據(jù)參數(shù)為1λ的指數(shù)分布律計(jì)算，而λ可根據(jù)歷史數(shù)據(jù)進(jìn)行估計(jì)［29］。Xt?從一個(gè)跳躍點(diǎn)到下一個(gè)跳躍點(diǎn)之間按幾何布朗運(yùn)動(dòng)方式演化。假設(shè)一時(shí)間點(diǎn)t位于兩跳躍點(diǎn)τ1，τi+1之間(τ1＜t＜τi+1) ，帶跳隨機(jī)微分方程（11）的離散解可表示為:

2.1.3 分?jǐn)?shù)白噪聲+泊松噪聲組合模型

該模型與模型2的不同之處在于，含有Hurst參數(shù)H的{WtH，t≥0}是一個(gè)分?jǐn)?shù)布朗運(yùn)動(dòng)（fBm），它是一個(gè)連續(xù)中心化高斯過程，有協(xié)方差函數(shù):

當(dāng)H=1/2時(shí)，fBm是標(biāo)準(zhǔn)布朗運(yùn)動(dòng)，模型3退化為模型2。

假設(shè)跳躍時(shí)間點(diǎn)為τ1，τ2，…，方程式（14）可以用如下改進(jìn)的歐拉方法進(jìn)行離散求解:

其中，為滿足“準(zhǔn)確與效率”的雙重標(biāo)準(zhǔn)，使用一種快速的一維fBm生成器［30］模擬分?jǐn)?shù)白噪聲WtH，而用Yin［31］的方法處理fBm的差值:

2.2 模型參數(shù)的非參數(shù)估計(jì)法基于短期序列離散測(cè)量值進(jìn)行參數(shù)估計(jì)是隨機(jī)微分方程相關(guān)研究中的一項(xiàng)重要工作。國(guó)內(nèi)外文獻(xiàn)中已有大量方法，例如最大似然估計(jì)法（MLE）、廣義矩估計(jì)法、模擬矩估計(jì)法和MCMC方法等，每一種都有自己的獨(dú)特特征。S?rensen［32］經(jīng)過詳細(xì)調(diào)查和比較，證明MLE估計(jì)法由于其連續(xù)性、漸進(jìn)正態(tài)性、漸進(jìn)有效性等特征，在弱正則條件下表現(xiàn)最好。并且，MLE在從擴(kuò)散項(xiàng)中分離跳躍項(xiàng)方面也具有優(yōu)越性。因此，針對(duì)本研究提出的隨機(jī)微分方程，下面將使用一種非參數(shù)估計(jì)方法，充分發(fā)揮MLE法的優(yōu)越性，并排除傳統(tǒng)參數(shù)估計(jì)法帶來的限制，特別是重點(diǎn)將放在轉(zhuǎn)移密度函數(shù)的有效構(gòu)建上。

假設(shè){xi，i=0，1，…，N} 為實(shí)測(cè)值序列；θ為待定參數(shù)集向量（包括確定項(xiàng)河道參數(shù)及噪聲項(xiàng)參數(shù)）；Pθ(ti，xi； (ti-1，xi-1)，θ)是從xi-1到xi的轉(zhuǎn)移密度函數(shù)；那么θ的極大似然函數(shù)就為下面根據(jù)如下步驟采用蒙特卡洛模擬方法對(duì)L()θ進(jìn)行估計(jì):

（1）假設(shè)時(shí)間范圍{Ft:t∈[0 ，T]}。實(shí)測(cè)值時(shí)間間隔:考慮時(shí)間區(qū)間為[ ]ti-1，ti，把該區(qū)間插入M個(gè)長(zhǎng)度為的子區(qū)間，那么以Xti-1為初值，利用Euler-Maruyama算法，就可以通過隨機(jī)微分方程的計(jì)算得到Xti的近似值。重復(fù)R次，可以得到R維向量

其中，核函數(shù)取標(biāo)準(zhǔn)正態(tài)形式，帶寬hpi取流行的Silverman［33］形式:

其中si表示的標(biāo)準(zhǔn)差。

對(duì)于帶跳隨機(jī)微分方程的情況，假設(shè)lr=ln(Xr)-ln(x)(r=1，2，…，R)，用m，n分別表iii-1ud示時(shí)間間隔Δ內(nèi)發(fā)生上跳和下跳情況的數(shù)目，可以用如下4種條件概率密度函數(shù)的泊松加權(quán)和來表示區(qū)間[ti-1，ti]內(nèi)的非條件概率密度函數(shù)，進(jìn)而得到轉(zhuǎn)移密度函數(shù):

（3）對(duì)每一個(gè)xi重復(fù)以上步驟，就可以得到極大似然函數(shù)。經(jīng)過對(duì)數(shù)變換，我們可以用Matlab中的fmincon函數(shù)求解如下負(fù)對(duì)數(shù)極大似然函數(shù)的極小值，從而得到需要的參數(shù)估計(jì)值:

3 實(shí)例分析

3.1 研究區(qū)域黃河下游高村至孫口河段（如圖1）經(jīng)常出現(xiàn)過流能力明顯小于其上下游河段的駝峰現(xiàn)象，且在高村附近河底高程抬升幅度最大［34］。筆者也加入眾多對(duì)此駝峰河段重點(diǎn)關(guān)注的學(xué)者［34-38］行列，在此選取1952年到2013年間黃河下游高村至孫口段的深泓線縱比降，以及高村站河流橫斷面形態(tài)要素（河寬、水深、流速）資料［39］作為原始數(shù)據(jù)，對(duì)本研究提出的隨機(jī)微分方程模型進(jìn)行非參數(shù)估計(jì)以及演變趨勢(shì)研究。

根據(jù)吳保生方法［19］，利用高村站歷年平灘流量、汛期平均流量、來沙系數(shù)等資料（表1），通過最小二乘法得到流量公式（1）中的系數(shù)和指數(shù)K、b、c和β分別為30.26、-0.47、0.42和0.346。根據(jù)該公式計(jì)算高村站平灘流量，并與實(shí)測(cè)值比較，如圖2所示，得出結(jié)論符合較好。

圖1 黃河下游平面示意圖

圖2 高村站的平灘流量實(shí)測(cè)值與計(jì)算值對(duì)比

表1 高村站歷年汛期平均流量、來沙系數(shù)、實(shí)測(cè)及計(jì)算平灘流量

3.2 模型參數(shù)估計(jì)對(duì)于未知參數(shù)集θ，采用2.2節(jié)介紹方法進(jìn)行參數(shù)估計(jì)。本例中，時(shí)間區(qū)間為［1952，1982］，初值取1952年河相關(guān)系各特征數(shù)據(jù)，時(shí)間間隔Δ=1，步長(zhǎng)取在單次估值試算中，由于Matlab中的fmincon函數(shù)對(duì)初值依賴性較強(qiáng)，導(dǎo)致結(jié)果失真嚴(yán)重。經(jīng)過多次嘗試確定的解決方法是:先給模型一個(gè)經(jīng)驗(yàn)初值進(jìn)行計(jì)算，之后每次計(jì)算前先計(jì)算已得結(jié)果的均值作為新的初值，經(jīng)過多次重復(fù)計(jì)算取全體估值的平均值為最終結(jié)果。當(dāng)求解帶跳模型的跳躍發(fā)生頻率λj(j=u，d)時(shí)，可根據(jù)其歷史數(shù)據(jù)計(jì)算［29］，例如:在離散時(shí)間區(qū)間［0，T］內(nèi)，將極端事件（X≥xu或X≤xd）發(fā)生次數(shù)定義為nj，則λj=njT。對(duì)于本例中的河相關(guān)系特征量，使用1000個(gè)樣本重復(fù)計(jì)算5000次，從而得到估計(jì)結(jié)果。此外為便于對(duì)比分析，應(yīng)用最小二乘法估計(jì)計(jì)算了確定性方程（6）的參數(shù)m。

3.3 模擬計(jì)算基于以上的模型參數(shù)估計(jì)結(jié)果，可以通過多次重復(fù)計(jì)算（5000次），近似得到3種模型條件下全時(shí)段的河相關(guān)系概率分布演化結(jié)果。其中特別是，在一給定跳躍強(qiáng)度條件下，通過最大化隨機(jī)分布均值與實(shí)測(cè)值的相關(guān)系數(shù)，從而可近似估計(jì)出最優(yōu)跳躍發(fā)生時(shí)間點(diǎn)。如圖3—圖6為不同條件下的河相關(guān)系概率演化結(jié)果與實(shí)測(cè)數(shù)據(jù)的平面對(duì)比關(guān)系圖。

從圖3—圖6可以看出，模擬的概率穩(wěn)定分布圖與實(shí)測(cè)值趨勢(shì)總體相同，在一定程度上證明了隨機(jī)微分方程描述河相關(guān)系動(dòng)態(tài)演化趨勢(shì)的有效性。然而對(duì)于模型1來說，問題是，在1990年左右隨機(jī)平均解大幅度偏離真實(shí)值；似乎跟預(yù)想的一樣，高斯白噪聲驅(qū)動(dòng)的隨機(jī)過程連續(xù)性有余，突變性不足，局部模擬值過大或過小，這無疑凸顯了模型對(duì)突發(fā)情況模擬能力的不足。實(shí)際動(dòng)態(tài)河流系統(tǒng)可能需要更加復(fù)雜的噪聲模型來刻畫，現(xiàn)在的白噪聲還不足以滿足要求。就像有足夠的離心力才能使得汽車平衡保持在賽道，白噪聲隨機(jī)方程還需要進(jìn)一步加強(qiáng)。

圖3 3種模型條件下比降的概率分布演化與實(shí)測(cè)值對(duì)比

圖4 3種模型條件下河寬的概率分布演化與實(shí)測(cè)值對(duì)比

包含高斯白噪聲和對(duì)數(shù)正態(tài)分布泊松噪聲的隨機(jī)微分方程最早由Merton［40］用于研究股票期權(quán)波動(dòng)?？紤]到氣候條件突變或人為擾動(dòng)可能給河相關(guān)系帶來類似脈沖的隨機(jī)影響，而這些隨機(jī)過程不管是發(fā)生頻率、發(fā)生時(shí)間或發(fā)生強(qiáng)度都是隨機(jī)的，因此筆者考慮在高斯白噪聲模型基礎(chǔ)上，增加泊松跳來對(duì)河相關(guān)系進(jìn)行研究。在該部分，河相關(guān)系的演化過程由幾何布朗運(yùn)動(dòng)表示的連續(xù)性過程與泊松跳表示的非連續(xù)性過程組成?？梢钥闯?，跟圖3—圖6中各子圖（a）相比，各子圖（b）的計(jì)算結(jié)果有了大幅度提升，特別是，平面概率穩(wěn)定帶明顯從平緩變得更加精致化，隨機(jī)平均解也更加貼近真實(shí)解的變化。這證明帶有泊松跳的隨機(jī)微分方程能顯著提高模擬水平，是基礎(chǔ)隨機(jī)理論的一個(gè)重要擴(kuò)展。唯一的遺憾是圖3—圖6的（b）中局部突出的變量大大超出了自定義的有效區(qū)域解范圍。某些時(shí)刻，較大（或較?。┑暮酉嚓P(guān)系特征量會(huì)發(fā)生過大的突變。因此帶跳隨機(jī)方程還不完美。除了突發(fā)性，可能還需要增加考慮非線性。

分?jǐn)?shù)布朗運(yùn)動(dòng)的研究最早可追溯到1940年代。Mandelbrot and Aizenman［41］命名并大大推動(dòng)了分形理論的發(fā)展。如今其已廣泛應(yīng)用于自然科學(xué)和工程科學(xué)的廣泛領(lǐng)域，例如水力學(xué)中的湍流研究。因此，筆者建立了分?jǐn)?shù)-帶跳隨機(jī)微分方程。整體來看，圖3—圖6中各子圖（c）比相應(yīng)子圖（b）表現(xiàn)更好。對(duì)于河寬和流速算例，是由于其Hurst參數(shù)（0.3490和0.3014）距離0.5較遠(yuǎn)，因此促成了顯著地時(shí)增擴(kuò)散現(xiàn)象。換句話說，河寬和流速過程時(shí)間序列是顯著負(fù)相關(guān)的，在相鄰時(shí)間點(diǎn)之間長(zhǎng)期發(fā)生高低值轉(zhuǎn)換現(xiàn)象。同時(shí)，對(duì)于比降和水深，該模型（3）計(jì)算結(jié)果與帶跳模型（2）計(jì)算結(jié)果十分接近，均成功預(yù)測(cè)出了局部跳躍點(diǎn)，但也漏掉了一些。還有令人遺憾的是，仍然有大量散點(diǎn)遠(yuǎn)離概率圖中的核心有效區(qū)間。筆者懷疑可能是由于真實(shí)的泊松項(xiàng)要遠(yuǎn)比我們這里的簡(jiǎn)單概化方法更復(fù)雜，其跳躍過大或過小都不合適?？赡苡行У母倪M(jìn)手段包括在跳躍強(qiáng)度或跳躍頻率中增加縮放因子，或者將本文隨機(jī)模型中的互不相關(guān)噪聲改為有相關(guān)關(guān)系的噪聲。

4 模型討論與應(yīng)用

4.1 模擬比較從圖3—圖6可以看出，應(yīng)用隨機(jī)微分方程，傳統(tǒng)的確定性的河相關(guān)系解進(jìn)化成為一種概率穩(wěn)定帶，從而為河流控制和河流管理提供了一種重要的新思路。顯然，如果外輪廓令人滿意，模型好壞評(píng)判準(zhǔn)則主要是與有效區(qū)域解，或者，隨機(jī)平均解與測(cè)量值相關(guān)度（如，相對(duì)誤差-δ；Nash-Sutcliff系數(shù)-NSE）等有關(guān)。也就是說，假如我們定義離散值在概率分布平面演化圖中的發(fā)生概率為有效概率，內(nèi)限值為有效穩(wěn)定寬度；那么可以說，有效概率越大，有效穩(wěn)定寬度越小，隨機(jī)平均解與測(cè)量值的相對(duì)誤差越小、Nash-Sutcliff系數(shù)越大，模型整體評(píng)價(jià)越高。所以，考慮到以上因子對(duì)模型調(diào)查和比較的重要意義，在圖7—圖8中基于上文的概率分布演化結(jié)果，給出相對(duì)有效穩(wěn)定寬度、相對(duì)誤差、Nash-Sutcliff系數(shù)的演化和比較結(jié)果圖。結(jié)合圖3—圖6，接下來我們對(duì)各河相關(guān)系特征量給出針對(duì)性的模型比較和最優(yōu)建議。

圖5 3種模型條件下水深的概率分布演化與實(shí)測(cè)值對(duì)比

圖6 3種模型條件下流速的概率分布演化與實(shí)測(cè)值對(duì)比

圖7 基于隨機(jī)微分方程計(jì)算的河相關(guān)系有效概率穩(wěn)定寬度時(shí)變演化

圖8 基于隨機(jī)微分方程計(jì)算的河相關(guān)系隨機(jī)平均解與實(shí)測(cè)值比較

（1）比降:圖3—圖6表明泊松跳能夠大幅度提高系統(tǒng)動(dòng)態(tài)演化的模擬效果。簡(jiǎn)單高斯白噪聲模型理應(yīng)被拋棄。同時(shí)根據(jù)圖7—圖8可以判斷，分?jǐn)?shù)-泊松模型由于較小的有效穩(wěn)定寬度、較大的有效概率和NSE系數(shù)，要?jiǎng)龠^簡(jiǎn)單帶跳隨機(jī)模型。

（2）河寬:與比降算例類似，根據(jù)圖3—圖6，我們可以直觀發(fā)現(xiàn)后兩種泊松模型勝過簡(jiǎn)單白噪聲模型的地方。然后，再分析圖7—圖8可以發(fā)現(xiàn)，簡(jiǎn)單帶跳隨機(jī)模型并沒有大幅度落后分?jǐn)?shù)-泊松模型，其有效概率穩(wěn)定寬度在后半段大部分時(shí)間，和NSE系數(shù)在全程的表現(xiàn)，均領(lǐng)先對(duì)手。分?jǐn)?shù)-泊松模型的優(yōu)勢(shì)是相對(duì)誤差-δ較小、隨機(jī)均解的波動(dòng)情況模擬較好。綜合來看，兩種泊松模型應(yīng)該結(jié)合在一起作為河流管理的參考資源。

（3）水深:與以上兩算例不同，圖3—圖6和圖7—圖8均顯示出高斯白噪聲+泊松噪聲模型在模擬水深動(dòng)態(tài)演進(jìn)過程中的優(yōu)勢(shì)，換句話說，水深主要受到突發(fā)隨機(jī)事件的線性影響。

（4）流速:從圖3—圖6和圖7—圖8可以看出，分?jǐn)?shù)-泊松模型要優(yōu)越于高斯-泊松模型，特別是在隨機(jī)均解波動(dòng)情況的模擬、相對(duì)誤差和NSE系數(shù)上的表現(xiàn)。顯然分?jǐn)?shù)-泊松模型更令人信服。

總結(jié)來看，分?jǐn)?shù)-泊松模型較適合模擬比降和流速的演進(jìn)，兩種泊松模型最好結(jié)合起來共同管控河寬，而簡(jiǎn)單高斯-泊松模型對(duì)于探索水深特征已經(jīng)足夠。

4.2 系統(tǒng)應(yīng)用系統(tǒng)穩(wěn)定性的評(píng)估始終都是理解河流的必要工作。有一些很好的方法已經(jīng)用于從整體上描述河流系統(tǒng)，包括能夠區(qū)分河型的河床穩(wěn)定指數(shù)（Zw）［42］（游蕩型:Zw＜5； 分汊型:5＜Zw＜15；彎曲型:Zw＞15），代表河岸侵蝕強(qiáng)度的寬深比系數(shù)（ζ），以及表示水流強(qiáng)度的河流功率（Ω）（式（24））。接下來，我們將使用分?jǐn)?shù)-泊松隨機(jī)微分方程對(duì)沿高村站向下游的Zw、ζ、Ω進(jìn)行概率分布演化的統(tǒng)計(jì)計(jì)算，結(jié)果如圖9所示。

圖9 基于分?jǐn)?shù)-泊松隨機(jī)微分方程計(jì)算的河床穩(wěn)定指數(shù)、寬深比系數(shù)和河流功率隨機(jī)演化圖

從圖9（a）中可以看出，河床穩(wěn)定指數(shù)從1986年左右開始大幅下滑，在此期間，其隨機(jī)均值幾乎跌破Zw=5的游蕩閾值線，直到接近2000年才又逐漸增大并回到分汊區(qū)間（5＜Zw＜15）。與此同時(shí)，高村站的水流強(qiáng)度經(jīng)歷緩慢增加過程（圖9（b）），且在2000年之后經(jīng)受的沖深作用大于展寬作用（圖9（c））。其原因，經(jīng)過查找文獻(xiàn)，我們可以了解到，是由于上游小浪底水庫(kù)的運(yùn)行化解了下游長(zhǎng)期少水多沙的負(fù)擔(dān)。這也與相關(guān)研究［35］做出的結(jié)論，認(rèn)為駝峰現(xiàn)象與水流流量大小成反比相吻合。整個(gè)圖9顯示出，河流穩(wěn)定性的低谷已經(jīng)在1999年被挽救，然而2008年左右出現(xiàn)的危機(jī)依然值得關(guān)注。并且，概率分布穩(wěn)定寬度與隨機(jī)均值間的正相關(guān)關(guān)系，也提醒我們未來需要提高對(duì)黃河下游河道演變的警戒水平。此外，考慮到測(cè)量解并不能完全反應(yīng)真實(shí)變化過程，例如2000年左右的寬深比系數(shù)（ζ），我們也因此可以看出隨機(jī)微分方程的優(yōu)越性。換句話說，對(duì)河相關(guān)系特征量個(gè)別單獨(dú)的測(cè)量值并不能幫助我們很好地了解河流系統(tǒng)性的變化，而在僅有少量汛期流量和懸移質(zhì)沙量的有限條件下，我們可以通過本研究的隨機(jī)方法表示的，例如寬深比系數(shù)（ζ），來進(jìn)行河流的系統(tǒng)性評(píng)估。通過Zw、ζ、Ω的隨機(jī)均值與平灘流量相關(guān)系數(shù)的比較（表2），我們可以得知流量主導(dǎo)的黃河下游適合用隨機(jī)狀態(tài)空間上的河流功率（Ω）來系統(tǒng)表征。近年我們收集的最新資料顯示，受水利與水土保持工程攔截、上游水庫(kù)調(diào)節(jié)、沿河用水增加及氣候變化等多重因素影響，黃河干流水沙形勢(shì)發(fā)生了明顯變化，表現(xiàn)為來水來沙量的持續(xù)減小。在新的水沙形勢(shì)和河道狀況下，河道自身過流能力顯著增強(qiáng)，同時(shí)中小流量出現(xiàn)的機(jī)率增加。這些都降低了洪水漫灘的機(jī)會(huì)。而長(zhǎng)期相關(guān)的影響還需要進(jìn)行系統(tǒng)性地測(cè)量、識(shí)別和評(píng)估。本研究具有重大的理論和現(xiàn)實(shí)意義。

表2 Zw、ζ、Ω的隨機(jī)平均解與平灘流量的相關(guān)系數(shù)

5 結(jié)論

河流是一個(gè)由流量驅(qū)動(dòng)，以河相關(guān)系為重要表現(xiàn)形式的動(dòng)態(tài)系統(tǒng)。本文通過在確定性方程基礎(chǔ)上加入隨機(jī)噪聲模型（包括簡(jiǎn)單高斯白噪聲模型、高斯-泊松模型和分?jǐn)?shù)-泊松模型），成功地建立并應(yīng)用隨機(jī)微分方程研究了黃河下游高村-孫口段的河相關(guān)系隨時(shí)間演化的概率分布動(dòng)態(tài)演化過程，為河流基礎(chǔ)理論做出了重要貢獻(xiàn)，同時(shí)為相關(guān)河流決策提供了一個(gè)新思路。本方法可以幫助我們?cè)谟邢薜臈l件下系統(tǒng)性地掌握河流變化特征。通過模型比較，本文認(rèn)定同時(shí)考慮非線性和突發(fā)性擾動(dòng)特征的分?jǐn)?shù)-泊松擴(kuò)散模型是適宜研究河相關(guān)系隨機(jī)演化特征的較優(yōu)模型。未來將把該方法進(jìn)行更多的推廣應(yīng)用和進(jìn)一步的完善，包括避免使用簡(jiǎn)化的噪聲模型，而使用互相相關(guān)的更真實(shí)的復(fù)雜噪聲模型等。