基于盒維數(shù)的雅魯藏布江流域水系分形維數(shù)及影響因素研究

費俊源，吳鵬飛，劉金濤，2

（1.河海大學水文水資源學院，南京210098；2.河海大學水文水資源與水利工程科學國家重點實驗室，南京210098）

0 引 言

地貌水文過程研究的一大難點是定量化地描述水系結(jié)構(gòu)特征［1］，而水系結(jié)構(gòu)特征會對水資源空間分布［2］、聚居點規(guī)劃［3］、地貌演化［4］產(chǎn)生影響，這又使得此特征受到了廣泛關(guān)注。由于水系的局部與整體間具有自相似性，目前主要使用針對自相似結(jié)構(gòu)的分形理論描述水系結(jié)構(gòu)特征［5，6］。

在分形理論中，分形維數(shù)被用于量化描述非規(guī)則客體［7］。目前的水系分形特征研究主要依托數(shù)字高程模型（Digital Elevation Model，DEM），使用Horton 定律或計盒法獲取水系分維數(shù)［7，8］。除了針對中小流域的研究［5，9］，以較大流域為對象的水系分形研究大多將流域劃分為少量子區(qū)域或子流域，單獨分析各子區(qū)域的分形特征。例如，竇明等［2］依照水資源分區(qū)將淮河流域劃分為13個區(qū)域，發(fā)現(xiàn)各區(qū)域的水系盒維數(shù)能夠反映區(qū)域的人類活動、土地利用等情況；鄭楠炯等［8］將韓江流域劃分為9個子流域，發(fā)現(xiàn)子流域與整個流域間的分維數(shù)接近。但是，上述研究使用的分區(qū)較大，且分維數(shù)樣本數(shù)量較少，由于分維數(shù)與氣象［7］、地形地質(zhì)［10］間聯(lián)系緊密，少量的分區(qū)不足以充分反映大型流域內(nèi)氣象和地形條件強異質(zhì)性帶來的水系分維數(shù)的差異。

本研究以雅魯藏布江（以下簡稱雅江）流域為研究區(qū)，采用自動化提取方法得到流域內(nèi)大量不同面積尺度的子流域并通過計盒法獲取子流域的盒維數(shù)?；诖撕芯S數(shù)數(shù)據(jù)集分析雅江流域內(nèi)不同尺度的子流域水系分維數(shù)與降水、坡度、起伏度等氣象、地形特征的關(guān)系。

1 研究區(qū)及數(shù)據(jù)源

雅江流域位于西藏自治區(qū)南部，涉及拉薩市、山南市、日喀則市、林芝市4 個地級市部分地區(qū)，流域面積24.8 萬km2，絕大多數(shù)區(qū)域海拔＞3 000 m，是世界上平均海拔最高的流域［11］。流域隨雅江自西向東延伸，南北方向相對狹窄。流域內(nèi)降水的空間分布差異大，同時存在干旱、半干旱半濕潤、濕潤地區(qū)［11］。

圖1 雅魯藏布江流域概況Fig.1 Overview of Yarlung Zangbo river basin

本研究使用1″（約30 m）分辨率的AW3D30 DEM，數(shù)據(jù)下載自日本宇宙航空研究開發(fā)機構(gòu)網(wǎng)站（https：//www.eorc.jaxa.jp/ALOS/en/aw3d30）。雅江流域范圍由更大范圍的DEM，經(jīng)填洼和流向識別［12，13］，依照雅江出口點提取。研究使用1 km 分辨率青藏高原地區(qū)（2000-2015）降水數(shù)據(jù)，下載自國家青藏高原科學數(shù)據(jù)中心（http：//data.tpdc.ac.cn），并以此計算出流域各位置的平均年降水量數(shù)據(jù)。

2 研究方法

2.1 水系及子流域提取

本研究首先對DEM 進行填洼預處理，移除DEM 內(nèi)因采樣誤差導致的洼地和平地，再計算每個柵格單元的流向。然后，依據(jù)流向累積得到每個單元的上游匯水面積。通過設(shè)置臨界源面積（Critical Source Area，CSA），提取匯水面積大于CSA的單元作為雅江的主要水系。最后，以水系的全部源點和部分主要干支流交匯節(jié)點作為子流域窗口，提取其完整的上游區(qū)域作為不同尺度的子流域以供分析。具體方案如下：

首先，自行編程實現(xiàn)由Garbrecht 和Martz［12］提出的算法進行填洼，然后使用Wu 等［13］提出的iFAD8 單流向算法獲取單元流向并進一步得到上游匯水面積。

為了綜合考慮不同尺度流域間的分形維數(shù)差異，以250 km2作為CSA 閾值，提取出雅江流域的主要水系。將所有水系源點對應的上游區(qū)域作為最小子流域。再選取2個上游入流匯水面積差異較?。ㄝ^小入流的匯水面積不得少于較大入流的1/3）的單元作為主要干支流交匯點，提取對應的完整上游流域作為較大面積尺度子流域。值得注意的是，不同尺度子流域間存在嵌套現(xiàn)象，即較小子流域可能被某一較大子流域覆蓋。這一方面是為了保證用于分析的是完整流域，另一方面則是為了保證充足的數(shù)據(jù)量。對此，本研究將所有子流域劃分至不同面積區(qū)間，主要針對尺度近似的子流域進行分析。盡管依然存在少量處于同面積區(qū)間但相互嵌套的子流域，但由于后續(xù)使用的都是面平均參數(shù)，可以認為這些子流域依然具有代表性。

對于每個子流域，通過設(shè)置新的CSA 閾值得到其內(nèi)部的水系，而本研究使用拐點法確定該閾值［7，14］。拐點法通過檢查不同匯水面積閾值下盒維數(shù)變化趨勢，選擇盒維數(shù)變化趨勢的拐點作為最終閾值。由于本研究使用的子流域數(shù)目較多，難以逐個確定每個子流域的最佳閾值，因此使用不同面積閾值下所有子流域的平均盒維數(shù)的變化趨勢確定拐點，將確定的拐點閾值作為所有子流域的最佳匯水面積閾值。

2.2 盒維數(shù)計算方法

首先將柵格DEM 區(qū)域劃分成大量邊長為r的正方形窗口，然后檢索出存在水系的窗口數(shù)量N（r）。r、N（r）與盒維數(shù)D間存在如下關(guān)系［15］：

根據(jù)上述關(guān)系，使用一組r和N（r）數(shù)據(jù)點繪ln（r）-ln［N（r）］關(guān)系圖，由最小二乘法擬合的直線斜率的絕對值即為所求盒維數(shù)值。由于當r的擴大倍數(shù)超過2 時分維數(shù)會出現(xiàn)波動［16］，因此研究將r的擴大倍數(shù)固定為2，以5 倍柵格單元邊長（約150 m）作為初始的正方形窗口邊長r。

3 結(jié)果與分析

3.1 子流域劃分及匯水面積閾值確定

本研究從雅江流域中提取出了356 個子流域，其中依據(jù)259 個水系源點提取得到259 個最小子流域，并依據(jù)3.1 節(jié)的方法選取97 個主要水系節(jié)點提取了97 個較大尺度的子流域。然后，以250 個柵格為間距，共計算了9 個CSA 閾值（500～2 500 之間）對應的平均盒維數(shù)，結(jié)果如圖2（a）所示?？梢园l(fā)現(xiàn)，平均盒維數(shù)隨著CSA 的增大而減小。當以1 500 個柵格（灰色點）為分界點時，匯流柵格數(shù)目＜1 500 個（黑色點）和＞1 500 個（空心點）對應的平均盒維數(shù)呈現(xiàn)出2 種不同的線性減小趨勢，且二者的擬合效果極佳。此外，如圖2（b）所示，各CSA 對應的所有子流域盒維數(shù)間平均絕對誤差和均方根誤差都較小，整體與平均值接近，可見采用平均盒維數(shù)確定閾值具備代表性。因此，將拐點對應的1 500 個柵格面積（約1.35 km2）作為匯水面積閾值提取水系。

圖2 平均盒維數(shù)隨上游匯流柵格數(shù)目的變化關(guān)系和所有子流域盒維數(shù)的平均絕對誤差和均方根誤差隨匯流柵格數(shù)目的變化關(guān)系Fig.2 The relationship between the average box dimension and the number of upstream grids and The relationship between the mean absolute error and root mean square error of all sub-basins’box dimension and the number of upstream grids

3.2 分形維數(shù)提取結(jié)果統(tǒng)計

根據(jù)流域的面積大小，研究以500、1 000、2 000、4 000 km2為面積尺度分區(qū)節(jié)點，分5 個面積區(qū)間對356 個子流域的氣候、地形特征以及分形維數(shù)分別進行了統(tǒng)計（表1）。其中起伏度的定義為流域最大高程和最小高程之差。可以看出，提取出的面積250～500 km2的子流域數(shù)目最多，為259 個，也就是有259 個由水系源點提取的最小子流域。值得注意的是，由于DEM單元匯水面積增長的躍進特征，由CSA 得到的水系源點的匯水面積不一定等于CSA 對應的250 km2。其他4 個面積區(qū)間內(nèi)的子流域較少，但也均超過10個。本研究得到的分區(qū)盒維數(shù)均值最小為1.04，最大為1.17。這一數(shù)值范圍與鄰近流域的研究成果接近［5，16］，屬于較小的盒維數(shù)范圍。水系發(fā)育程度越高的區(qū)域盒維數(shù)越大［16］，由此體現(xiàn)出雅江流域水系發(fā)育程度較低。雅江流域的盒維數(shù)D≤1.6，表明雅江流域的河流地貌處于侵蝕發(fā)育的幼年期［17］。根據(jù)表1 中的數(shù)據(jù)可以發(fā)現(xiàn)，盒維數(shù)均值存在隨面積區(qū)間增大而增大的情況，但增幅并不明顯。這是由于較大流域包含更多下游河谷區(qū)域，該區(qū)域水系發(fā)育程度優(yōu)于上游山丘區(qū)。此外，不同面積子流域的平均坡度、平均海拔相近，起伏度呈現(xiàn)出隨流域面積增大而增加的趨勢。

表1 不同面積范圍的流域降水、地形特征及分維數(shù)統(tǒng)計Tab.1 Precipitation，topographic features and fractal dimension statistics of catchments with different area scales

3.3 分形維數(shù)和地形及氣象特征的關(guān)系

本研究分別分析了5 個面積區(qū)間內(nèi)盒維數(shù)隨平均坡度、起伏度以及年平均降水的變化關(guān)系，并用線性關(guān)系、指數(shù)型關(guān)系、對數(shù)型關(guān)系擬合了關(guān)系函數(shù)。對于所有數(shù)據(jù)組，3 種擬合關(guān)系的確定系數(shù)差值均小于0.04，因此下文均以最簡單的線性關(guān)系進行描述分析。

圖3 展示了盒維數(shù)隨平均坡度的變化。在5 種面積區(qū)間內(nèi)，盒維數(shù)均呈現(xiàn)隨平均坡度增大而減小的趨勢。這體現(xiàn)出雅江流域中高地形起伏的上游區(qū)域河流處于發(fā)育更早期的階段，越平坦的下游區(qū)域河流發(fā)育情況越好。此外，隨著面積區(qū)間的增大，盒維數(shù)與平均坡度間擬合的確定系數(shù)同樣增加。這是由于小流域的水系發(fā)育更容易受各類局地因素的影響，大流域更能體現(xiàn)盒維數(shù)（或河網(wǎng)密度）與平均坡度間的相關(guān)性。

圖4展示了盒維數(shù)與起伏度之間的相關(guān)關(guān)系。大致可以判斷出雅江流域內(nèi)子流域存在盒維數(shù)隨著起伏度的增加而減小的趨勢。但整體看來，盒維數(shù)與起伏度間的相關(guān)性不佳。由于高起伏度地區(qū)的平均坡度傾向于更大，這里的變化趨勢可以視作圖3中盒維數(shù)與流域平均坡度相關(guān)性的體現(xiàn)。

圖3 不同面積（s）范圍下分形盒維數(shù)與流域平均坡度的關(guān)系Fig.3 Relationship between box dimension and average slope of catchments with different area scales（s）

圖4 不同面積范圍（s）下分形盒維數(shù)與起伏度的關(guān)系Fig.4 Relationship between box dimension and average topographic relief of catchments with different area scales（s）

圖5 展示了不同條件下水系盒維數(shù)與年平均降水量的關(guān)系。雖然不同面積區(qū)間內(nèi)盒維數(shù)與平均年降水量的確定性系數(shù)存在波動，但可以發(fā)現(xiàn)整體上不同面積區(qū)間的流域水系盒維數(shù)均呈現(xiàn)出隨年平均降水量增大而減小的趨勢。研究顯示的雅江流域內(nèi)部的盒維數(shù)隨年平均降水量的變化趨勢與王博等［16］統(tǒng)計全國多個流域得到的濕潤半濕潤地區(qū)分形維數(shù)普遍高于干旱半干旱地區(qū)的結(jié)果不同。結(jié)合分形維數(shù)與河網(wǎng)密度的對應關(guān)系，此差異也許可以類比河網(wǎng)密度與降水量的U 型關(guān)系。Abrahams 等［18］證明了河網(wǎng)密度隨干旱程度呈U 型關(guān)系，即以某一干旱情境為節(jié)點，在相對此情境越干旱或越濕潤的情況下，河網(wǎng)密度都會增加。因此，可以推斷分形維數(shù)與降水量可能也呈類似U 型關(guān)系，本研究和已有研究發(fā)現(xiàn)的特征規(guī)律大致分屬U 型的兩端。在本研究使用的子流域中，當年平均降水量較大（＞600 mm）時，部分子流域并未遵從隨年平均降水量增大而遞減的趨勢，因此可以觀察到大量盒維數(shù)增大的子流域［如圖5（a）、圖5（b）、圖5（c）］，正是這些子流域?qū)е铝讼嚓P(guān)區(qū)間內(nèi)線性擬合效果較差，也間接佐證盒維數(shù)與平均年降水可能存在U型關(guān)系。隨著面積尺度的增加，受限于流域數(shù)目及降水量分布，難以觀測到水系分形維數(shù)與降水的U型關(guān)系。

圖5 不同面積范圍（s）下盒維數(shù)與年平均降水的關(guān)系Fig.5 Relationship between box dimension and average precipitation of catchments with different area scales（s）

4 結(jié)果與討論

針對已有水系分維研究中樣本數(shù)量過少的問題，本研究使用自動化處理方法從雅江流域選取了合適的節(jié)點，提取出了356個不同面積的子流域，經(jīng)拐點法確定1.35 km2為最適合雅江流域分形研究的匯流面積閾值。然后使用計盒法計算了子流域的水系分形維數(shù)，發(fā)現(xiàn)雅江流域的河流處于地貌發(fā)育的幼年期。研究還劃分了5 個流域面積區(qū)間，分別分析了不同面積尺度下盒維數(shù)與流域平均坡度、起伏度以及平均年降水量間的關(guān)系。結(jié)果顯示，使用盒維數(shù)與平均坡度和平均年降水量間相關(guān)性較好，大致呈線性下降的關(guān)系，不同面積區(qū)間下擬合函數(shù)的確定系數(shù)則呈現(xiàn)出隨流域面積增大而增大的趨勢。而流域的分形維數(shù)與起伏度間的關(guān)系則較差。此外，根據(jù)本研究發(fā)現(xiàn)流域的盒維數(shù)與平均年降水量間關(guān)系呈現(xiàn)出非單調(diào)趨勢，推斷盒維數(shù)應與河網(wǎng)密度一致，隨流域的干旱程度呈U 型關(guān)系。其中，較小面積尺度流域的盒維數(shù)與平均降水一定程度上呈U 型關(guān)系，而此關(guān)系在較大面積尺度的流域則不明顯，未來仍需要增加不同降水量的大面積尺度流域樣本數(shù)量做進一步分析驗證。

本研究主要關(guān)注雅江流域水系分形結(jié)構(gòu)與地形、氣象要素的關(guān)系。但是部分已有研究表明，水系的分形結(jié)構(gòu)與植被、土地利用等多種其他要素同樣存在關(guān)聯(lián)［2，19］，且這些要素與社會、經(jīng)濟關(guān)系密切。因此，后續(xù)研究需要關(guān)注相關(guān)領(lǐng)域，將水系分形特征在生產(chǎn)生活中的推廣應用。 □