基于ＭＷＰ 檢驗的黃河內蒙古河段水沙過程時變點分析

張羽 丁夢霞 董明家 劉明瀟

摘 要：運用統計學方法確定水沙過程發生變異的原因及標志性節點，是理解河勢演變制約機制與發展趨勢的一種重要研究手段。近年提出的ＭＷＰ 非參數檢驗法對時序數據總體分布的假定要求較寬，具有避免嚴重分析誤差、適用性強的特點。采用這一方法對黃河內蒙古河段３ 個水文站１９７３—２００４ 年實測徑流泥沙資料進行時間序列分析，獲得了以含沙量、流量為指標的多級水沙時變點信息。研究表明，河流特殊水文事件會引起水沙時間序列時變點的響應，Ｕｔ曲線的波峰與波谷（ＫＴ ）體現了時間序列中對水沙過程影響最大的節點即非參數檢驗統計突變點。一級時變點是長時段水沙過程出現趨勢性變化的標志性節點，上游龍羊峽水庫運用就在一級時變點中得到映射。二級時變點多為短期偶發事件引起水沙過程局部調整的標志性節點。借助ＭＷＰ 檢驗法確定河流水沙時間序列的多級時變點，可確定影響河流水沙過程長期變化趨勢與局部時段調整的關鍵因子，追溯水文情勢變化誘因，有助于認識近年黃河水沙變異的態勢。

關鍵詞：水沙過程；ＭＷＰ 非參數檢驗法；時變點；響應；黃河內蒙古河段

中圖分類號： Ｐ３３３．４；ＴＶ８８２．１ 文獻標志碼：Ａ ｄｏｉ：１０．３９６９／ ｊ．ｉｓｓｎ．１０００－１３７９．２０２３．０７．００５

引用格式：張羽，丁夢霞，董明家，等．基于ＭＷＰ 檢驗的黃河內蒙古河段水沙過程時變點分析［Ｊ］．人民黃河，２０２３，４５（７）：２４－２９．

近年來，黃河內蒙古河段受上游大型水庫運用、氣候及水文異常等復雜環境因素影響，來水來沙條件發生變化，進而導致該河段河勢不斷調整。內蒙古河段水沙關系的變異及其作用過程是該河段河勢變化的主要動力因素，研究長序列水沙過程的變化特性，可以探尋河床演變的發生機理，尋找兩者在時空上的響應關系［１－２］ 。沖積性河流上大型水利工程的修建及其多種方式的運用會在不同時間段引起下游水沙過程發生變化。通過統計學方法確定水沙過程中發生顯著變化的時間臨界點（時變點），對于理解河勢演變的成因及規律是至關重要的。

為了探尋水文序列中徑流泥沙的變化特征，目前常用的變異點分析方法有Ｒ／ Ｓ 分析法［３－４］ 、Ｔ 檢驗法［５－６］ 、Ｆ 檢驗法［７－８］ 、秩和檢驗法［９－１０］ 、Ｍａｎｎ－Ｋｅｎｄａｌｌ檢驗法［１１］ 、Ｍａｎｎ－Ｗｈｉｔｎｅｙ－Ｐｅｔｔｉｔｔ（簡稱ＭＷＰ） 檢驗法［１２－１５］ 、有序聚類法等［１６］ 。相比其他方法，ＭＷＰ 檢驗法不要求樣本所屬的總體符合某種理論分布，檢驗不受總體分布形態的限制，適用范圍較為廣泛，對于原始數據不要求有很精確的計量值，可以使用分類數據和順序數據，幾乎可以處理任何類型的數據，達到要求的檢驗精度所需要的樣本容量相比其他方法小。

近年來也有學者采用ＭＷＰ 非參數檢驗法開展氣候與水文變化規律的研究。呂季蓉［１７］ 建立臺灣南部地區各雨量站的標準化降雨指標，對照歷史干旱記錄，使用ＭＷＰ 檢驗法，檢驗了近百年長期雨量與氣象資料的變化趨勢，并研究了干旱發生的頻率。樊輝等［１８］利用１９５０—２００７ 年黃河利津站水沙數據，采用ＭＷＰ檢驗法，對黃河入海水沙通量變化情況進行分析，發現徑流泥沙過程出現多個時序轉折點，與時變點對應的黃河入海徑流通量與泥沙通量顯著減少，汛期變化幅度大于非汛期的。丁明軍等［１９］ 根據鄱陽湖周邊１０ 個氣象觀測站的日均氣溫、日最高氣溫以及日最低氣溫資料，利用Ｍａｎｎ－Ｋｅｎｄａｌｌ 趨勢分析法和ＭＷＰ 檢驗法，分析了１９６１—２００７ 年鄱陽湖周邊地區氣溫時序變化特征，發現日均氣溫和日最低氣溫的冬季均溫突變時間開始得早。

雷紅富等［２０］ 對上述檢測方法進行了評價，通過設計的水文變異序列生成器生成滿足Ｐ 型分布的統計參數，用統計試驗的方法比較分析了１０ 種常用水文序列變異點檢驗方法對３ 類變異序列的檢驗性能，認為秩和檢驗、ＭＷＰ 檢驗法效果較好，特別是ＭＷＰ 檢驗法，具有不受總體分布形態限制的優越性，可以更為清晰、準確地捕捉水文序列中的特征時變點。

根據上述徑流泥沙變化特征研究方法的討論，筆者采用ＭＷＰ 非參數檢驗法對黃河內蒙古河段實測的歷史水沙數據進行分析，找出水沙序列的顯著變化時刻，將整個研究時間段分為若干時段，再與引起河勢變化的各種因素發生變化的時間進行關聯，探尋上述因素與黃河內蒙古河段河勢變化的響應關聯，以期為深入了解河道演變趨勢與影響因素提供基礎性數據依據。

１ 基本原理與研究方法

１．１ ＭＷＰ 非參數檢驗法的基本原理

ＭＷＰ 檢驗法是統計學中對基于時間序列的連續數據變化特征進行檢測的一種方法。該方法認為在連續數據出現顯著差異的前后，其特征值會有明顯的變化，可以通過對多級時變點的尋找，發現連續數據的變化規律。ＭＷＰ 檢驗的基本原理如下。

１．２ ＭＷＰ 非參數檢驗法的應用特點

ＭＷＰ 非參數統計與傳統的參數統計相比，要求的假定條件比較少，適用范圍更廣，對時序數據總體分布的假定要求條件很寬，因此不會因總體分布假定而導致較大誤差。當研究問題推論多于３ 個時，ＭＷＰ 非參數檢驗方法具有優越性。

２ 研究對象與計算結果

２．１ 研究對象

２．１．１ 黃河內蒙古河段形態及水沙特征

黃河內蒙古河段是黃河上游的沖積性河段，全長６７３ ｋｍ。黃河內蒙古河段出石嘴山峽谷段后進入內蒙古平原，比降逐漸減小，平均坡降０．０１％。自巴彥高勒至頭道拐河段河流寬淺、散亂，彎道達６９ 處，河槽擺動幅度３～７ ｋｍ／ ａ。汛期水量占全年的７０％，最大與最小年徑流量之比為２．６８，年徑流變差系數ＣＶ 為０．１６，天然徑流具有豐枯交替的特征，每個周期約１７ ａ。以內蒙古河段末端的頭道拐水文站為例，多年平均徑流量與輸沙量見表１。內蒙古河段出口頭道拐水文站的水沙變化過程如圖１ 所示。

內蒙古河段有３ 個重要的水文站，分別為巴彥高勒站、三湖河口站和昭君墳站，分別位于該河段的上游、中游和下游，３ 個水文站所處的河段具有不同的沖淤規律，有時呈現上沖下淤的沿程沖淤變化特點。３個水文站均具有１９７３—２００４ 年的月均徑流泥沙數據，可以用來分析內蒙古河段水沙過程變化特點。

２．１．２ 黃河內蒙古河段上游的水電工程及影響

黃河內蒙古河段的上游修建有多座水庫，其中較大的有劉家峽水庫（１９６８ 年） 和龍羊峽水庫（１９８６年），它們在調節徑流與洪水過程方面起了很大作用，雖然離內蒙古河段很遠，但也給水庫下游河道的水沙運動與河床演變帶來很大影響，例如河床活動性增強，主槽排洪能力降低［２１］ 。上游水庫運用、本河段支流（孔兌）高含沙洪水入匯等都會對該河段的水沙情勢造成不同的影響。２．２ 內蒙古河段水沙過程時變點計算結果本研究選擇實測含沙量、流量兩個參數作為主要的隨機因子，運用ＭＷＰ 檢驗法對３ 個代表水文站月平均流量序列及含沙量序列進行分析，探尋水沙時變點及其反映的河流信息。

選取有年內水沙過程詳細資料的典型時間序列（１９７３—２００４ 年）進行內蒙古河段水沙過程時變點計算。這一時段接近兩個水文循環周期，期間有大型水利工程投入運用。首先統計整理了內蒙古河段３ 個典型水文站在所選研究時段的月平均含沙量和流量，同時按年、月排序得到３ 站的月平均含沙量、流量的３２ ａ時間序列（共計３８４ 個月）。對內蒙古河段這一時間序列的水沙代表數據進行了非參數檢驗計算，以下給出對各站月均含沙量、流量數據進行ＭＷＰ 檢驗的計算結果。

從１９７３ 年１ 月至２００４ 年１２ 月，有３８４ 個月數據的時間序列。以３ 站對應的月均含沙量、流量為數據樣本，以時刻ｔ 為分割點將該時間序列分為前后兩部分，計算特征指標Ｕｔ值。Ｕｔ計算采用基于Ｍａｔｌａｂ 運算程序，獲得分析時段一系列的Ｕｔ ，做出Ｕｔ 在分析時段的變化曲線，如圖２ 和圖３ 所示，體現為具有不同震蕩特征的波峰與波谷過程。波峰（波谷）過程曲線的最大（小）值對應的Ｕｔ 極值點就是對水沙過程具有重要影響作用的時變點（統計突變值）。

以圖２（ｂ）為例，可以看到Ｕｔ 曲線總體呈現波谷過程，但有局部波動震蕩，｜Ｕｔ ｜最大值位置大約在時間序列（月序）２００ 與２５０ 之間；于是極值點２２０ 即為此時間序列的一級時變點，對應時間為１９９１ 年４ 月，這是３２ ａ 水沙序列中體現三湖河口站月均含沙量變化特征的時變點。利用式（３）、式（４），可以獲得最顯著改變時間點對應的｜Ｕｔ ｜最大值即ＫＴ ，通過Ｍａｔｌａｂ 計算程序獲得ＫＴ 為１６ ３８２，對應發生概率ｐ 為０．９６５。同理也可以得到３ 個水文站其他水沙指標Ｕｔ 曲線的時變點和對應的峰值ＫＴ以及ＫＴ發生概率，見表２。

３ 結果分析與討論

３．１ 一級時變點分析

根據對３ 站水沙序列數據進行的ＭＷＰ 非參數檢驗，得到了每一個水沙指標的時變點。為了更清晰地觀察時變點在整個時間序列中的態勢和作用，將所有月均流量、含沙量與整個時間序列的特征流量、特征含沙量進行比對，并分析與時變點的關聯關系。這里將其比值設為縱坐標，整個時間序列為橫坐標，分析３ 個水文站３２ ａ 的相對水沙過程及其水沙指標時變點位置，見圖４～ 圖６。豐水年到枯水年流量Ｕｔ 為波峰曲線，枯水年到豐水年流量Ｕｔ 為波谷曲線，Ｕｔ 絕對值越大整個時間序列上水沙特征值的變化越劇烈，Ｕｔ 絕對值越小整個時間序列上水沙特征值的變化越平緩。

分析圖４～圖６ 及表２ 中含沙量、流量的時間序列特征值和一級時變點位置，結合３ 個站上游的大型水電工程龍羊峽水庫運用情況，可以解讀出如下信息：

１）１９８６ 年龍羊峽水庫運用后，河道徑流泥沙過程發生很大改變，流量過程調勻、漫灘洪水大大減少、小水期大大延長。而經過ＭＷＰ 非參數檢驗，發現３ 站中最上游的巴彥高勒站的流量一級時變點出現在１９８６ 年，含沙量一級時變點出現在１９８８ 年。時變點反映出上游水庫運用改變了水沙同步的一般徑流過程，使得流量和含沙量一級時變點出現了錯位現象，即大水未必大沙，大沙未必大水。

２）１９８９ 年內蒙古河段的十大孔兌（支流）相繼出現高含沙洪水過程，在入黃口門附近發生嚴重淤堵，影響了附近三湖河口河段的河床邊界條件，再加上龍羊峽水庫運用影響，造成距巴彥高勒下游２００ 多ｋｍ 的三湖河口站水沙時變點均出現在１９９１ 年４ 月。時變點位置也體現了該站水沙響應的滯后效應。

３）昭君墳水文站位于十大孔兌之一的黑賴溝下游，由于最臨近十大孔兌，因此受孔兌影響最為直接。１９８６ 年前后十大孔兌經常發生高含沙洪水，直接導致附近河段水沙及河勢變化。采用ＭＷＰ 檢測方法得出昭君墳水文站在１９８６ 年８ 月和１９８５ 年１１ 月迅速響應分別形成水、沙的一級時變點。反映出上一年汛末河床強烈沖淤變化引起的含沙量一級時變點對次年汛期河槽行洪能力流量一級時變點出現的影響。

３．２ 二級時變點分析

以一級時變點為界將原時間序列分為兩個時間序列，分別以兩個時間序列為基礎再次進行ＭＷＰ 非參數檢驗，對水沙時間序列進行更精細的二級時變點分析。用ＭＷＰ 檢驗方法可計算出每個水文站水沙指標對應的兩個二級時變點，見表３。其中數據序列起始點與一級時變點之間的時變點稱為上變點，一級時變點與數據序列結束點之間的時變點稱為下變點。二級時變點可以從更精細時間尺度和更精準的層面研究水沙過程的變化情況。

分析表３ 中含沙量及流量的二級時變點位置，從水沙時序的二級時變點分析可以解讀出如下信息：

１）巴彥高勒站月均流量的二級時變上變點出現在１９８１ 年６ 月，這與１９８１ 年汛期持續５０ ｄ 的特大洪水有關。由于上游劉家峽水庫洪水調節作用有限，特大洪水造成內蒙古河段巴彥高勒至三湖河口段發生沖刷，突變的時間點與洪水時間基本吻合，說明水沙時序的二級時變點可準確刻畫上游洪水對該河段水沙情勢的影響。

２）昭君墳站月均流量的二級時變上變點雖然也同樣出現在１９８１ 年６ 月，但在這一時期昭君墳至頭道拐河段實際上發生了淤積，河床變形的性質與巴彥高勒至三湖河口段不同。昭君墳河段相應發生水沙條件的改變，出現二級時變點也是必然的，這也是在１９８６年６ 月受龍羊峽水庫運用影響發生大趨勢改變（一級時變點）下的局部調整。

３）１９７６ 年和１９７８ 年十大孔兌中的黑賴溝與西柳溝均有高含沙洪水災害的記錄，昭君墳站恰在這兩個孔兌之間，昭君墳站月均含沙量的二級時變上變點就在這一期間發生（１９７７ 年１０ 月）。昭君墳站含沙量二級時變點（上變點）是對應該站一級時變點１９８５ 年１１月（１９８５ 年底龍羊峽水庫影響）的局部調整。

４）１９８５ 年到１９８６ 年三湖河口上游出現了一次中水大沙過程，對三湖河口河段沖淤產生重要影響，造成了含沙量的二級時變點。１９９４ 年、１９９５ 年連續兩年中水沖刷了三湖河口河段河槽，１９９６ 年后連年枯水，形成了三湖河口流量的二級時變點。

４ 結論

基于ＭＷＰ 檢驗方法對黃河內蒙古河段水沙過程進行了分析研究，得到以下幾點主要結論：

１）ＭＷＰ 非參數檢驗法具有發現長時間序列隨機事件重要影響點、追溯河流水文情勢變化關鍵誘因的功能。通過確定水沙因子的多級時變點，可以了解長時間序列中的河流水沙變化趨勢與局部時段的調整特點。

２）通過對黃河內蒙古河段巴彥高勒、三湖河口及昭君墳３ 個水文站３２ ａ 水沙數據進行ＭＷＰ 非參數檢驗，建立了水沙因子的Ｕｔ 過程線，極值分析獲得了一級時變點和二級時變點。表明水沙時間序列時變點對特殊水文事件有所響應的特征。

３）一級時變點是長時段水沙過程出現趨勢性變化時的標志性節點，例如內蒙古河段水沙過程的一級時變點反映了上游龍羊峽水庫運用引起重大水文情勢變化影響；二級時變點多為偶發水文事件引起河段水沙局部調整的標志性節點，例如內蒙古河段支流高含沙洪水入匯對干流輸沙干擾的影響就在二級時變點得到明確體現。

４）分析時變點的特征及對應的水文事件，有助于了解河流系統的水文情勢變化調整原因，有助于認識近期黃河水沙變異及發展趨勢。

參考文獻：

［１］ 張邀丹，李志軍，李春江，等．黃河內蒙古段河冰微結構特性及影響因素分析［Ｊ］．水利學報，２０２１，５２（１２）：１４１８－１４２９．

［２］ 孫東坡，楊真真，張立，等．基于能量耗散關系的黃河內蒙古段河床形態調整分析［Ｊ］．水科學進展，２０１１，２２（５）：６５３－６６１．

［３］ 李保，周豐年，高健．基于Ｒ／ Ｓ 分析法的長江口總磷變化趨勢分析［Ｊ］．水利水電快報，２０１９，４０（１０）：５４－５８．

［４］ 譚曉慧．邊坡穩定的非線性有限元可靠度分析方法研究［Ｊ］．巖石力學與工程學報，２００８，２７（８）：１７２８．

［５］ 李春永，楊中蘭．前處理方式對化學需氧量測定的影響［Ｊ］．人民長江，２０２１，５２（增刊２）：４１－４４．

［６］ 王文周．改良ｔ 檢驗法［Ｊ］．西華大學學報（自然科學版），２００８，２７（６）：９１－９４，１２５．

［７］ 郭鴻飛．Ｆ 檢驗法和Ｔ 檢驗法在方法驗證過程中的應用探究［Ｊ］．山西冶金，２０１９，４２（４）：１１４－１１６．

［８］ 胡偉，張蕾．基于改進Ｆ 檢驗法對實際化工過程的穩態判斷［Ｊ］．廣東化工，２０２０，４７（６）：５３－５５．

［９］ 吳吉玲，黃一珍，姜鵬飛，等．排序法在大西洋鯖魚脫腥工藝篩選中的應用［Ｊ］．食品研究與開發，２０２０，４１（２３）：７３－７９．

［１０］ 拜存有，張升堂．渭河關中段年徑流過程變異點的診斷［Ｊ］．西北農林科技大學學報（自然科學版），２００９，３７

（１０）：２１５－２２０．

［１１］ ＣＵＲＩＡＣ Ｃ Ｄ，ＭＩＣＥＡ Ｍ Ｖ．Ｉｄｅｎｔｉｆｙｉｎｇ Ｈｏｔ Ｉｎｆｏｒｍａｔｉｏｎ Ｓｅ?ｃｕｒｉｔｙ Ｔｏｐｉｃｓ Ｕｓｉｎｇ ＬＤＡ ａｎｄ Ｍｕｌｔｉｖａｒｉａｔｅ Ｍａｎｎ?ＫｅｎｄａｌｌＴｅｓｔ［Ｊ］．ＩＥＥＥ Ａｃｃｅｓｓ，２０２３，１１：１８３７４－１８３８４．

［１２］ ＷＡＮＧ Ｘ Ｊ，ＹＵＡＮ Ｐ，ＹＵＡＮ Ｘ Ｍ，ｅｔ ａｌ．Ｐｏｗｅｒ Ｌａｗ Ｄｉｓｔｒｉ?ｂｕｔｉｏｎ Ｃｈａｒａｃｔｅｒｉｓｔｉｃｓ ｏｆ Ｄａｉｌｙ Ｓｔｒｅａｍｆｌｏｗ ｉｎ ｔｈｅ ＹｅｌｌｏｗＲｉｖｅｒ Ｂａｓｉｎ Ｕｎｄｅｒ ａ Ｃｈａｎｇｉｎｇ Ｅｎｖｉｒｏｎｍｅｎｔ［Ｊ］．Ｊｏｕｒｎａｌ ｏｆＷａｔｅｒ ａｎｄ Ｃｌｉｍａｔｅ Ｃｈａｎｇｅ，２０２０，１１（４）：１６０３－１６１８．

［１３］ ＣＨＥＮ Ｓ Ｔ，ＬＩ Ｙ Ｈ．Ｌｏｎｇ?Ｔｅｒｍ Ｐａｎ Ｅｖａｐｏｒａｔｉｏｎ Ｔｒｅｎｄ ａｎｄＰｒｅｌｉｍｉｎａｒｙ Ｉｎｖｅｓｔｉｇａｔｉｏｎ ｏｆ Ｅｖａｐｏｒａｔｉｏｎ Ｃｏｍｐｌｅｍｅｎｔａｒｉｔｙ ｉｎＴａｉｗａｎ［ Ｊ］． Ｊｏｕｒｎａｌ ｏｆ Ｔａｉｗａｎ Ａｇｒｉｃｕｌｔｕｒａｌ Ｅｎｇｉｎｅｅｒｉｎｇ，２０１６，６２（１）：１２－２８．

［１４］ ＫＵＮＤＵ Ｓ，ＫＨＡＲＥ Ｄ，ＭＯＮＤＡＬ Ａ，ｅｔ ａｌ．Ａｎａｌｙｓｉｓ ｏｆ Ｓｐａｔｉａｌａｎｄ Ｔｅｍｐｏｒａｌ Ｖａｒｉａｔｉｏｎ ｉｎ Ｒａｉｎｆａｌｌ Ｔｒｅｎｄ ｏｆ ＭａｄｈｙａＰｒａｄｅｓｈ，Ｉｎｄｉａ（１９０１－２０１１） ［Ｊ］．Ｅｎｖｉｒｏｎｍｅｎｔａｌ Ｅａｒｔｈ Ｓｃｉ?ｅｎｃｅｓ，２０１５，７３（１２）：８１９７－８２１６．

［１５］ ＷＡＮＧ Ｙ，ＷＡＮＧ Ｒ Ｙ，ＭＩＮＧ Ｊ，ｅｔ ａｌ．Ｅｆｆｅｃｔｓ ｏｆ Ｄｕｓｔ ＳｔｏｒｍＥｖｅｎｔｓ ｏｎ Ｗｅｅｋｌｙ Ｃｌｉｎｉｃ Ｖｉｓｉｔｓ Ｒｅｌａｔｅｄ ｔｏ Ｐｕｌｍｏｎａｒｙ Ｔｕ?ｂｅｒｃｕｌｏｓｉｓ Ｄｉｓｅａｓｅ ｉｎ Ｍｉｎｑｉｎ，Ｃｈｉｎａ［Ｊ］．Ａｔｍｏｓｐｈｅｒｉｃ Ｅｎｖｉ?ｒｏｎｍｅｎｔ，２０１６，１２７：２０５－２１２．

［１６］ 拜存有，蘇瑩，郭旭新．流域水文過程變點分析研究綜述［Ｊ］．水資源與水工程學報，２０１０，２１（１）：８３－８６．

［１７］ 呂季蓉．臺灣南部地區長期干旱趨勢分析之研究［Ｄ］．臺南：臺灣成功大學，２００６：５０－８１．

［１８］ 樊輝，劉艷霞，黃海軍．１９５０—２００７ 年黃河入海水沙通量變化趨勢及突變特征［Ｊ］．泥沙研究，２００９，３４（５）：９－１６．

［１９］ 丁明軍，鄭林，楊續超．１９６１—２００７ 年鄱陽湖周邊地區氣溫變化趨勢分析［Ｊ］．中國農業氣象，２０１０，３１（４）：５１７－５２１．

［２０］ 雷紅富，謝平，陳廣才，等．水文序列變異點檢驗方法的性能比較分析［Ｊ］．水電能源科學，２００７，２５（４）：３６－４０．

［２１］ 孫東坡，陳丹，張羽，等．基于系列水沙條件的黃河內蒙古段河床演變試驗研究［Ｊ］．泥沙研究，２０１２，３７（３）：７３－７９．

【責任編輯 張 帥】