基于Voronoi圖的林分空間模型及分布格局研究

劉 帥，吳舒辭，* ，王 紅，2，張 江，李建軍，王傳立

(1.中南林業科技大學，長沙 410004;2.湖南省科學技術廳，長沙 410013)

與林木空間位置有關的林分結構統稱為林分空間結構，主要包含樹種混交、大小分化和空間格局3個方面的內容。林分空間結構決定了林木間的競爭優勢及其空間生態位，在很大程度上影響著林分的健康和穩定。采用何種手段從而更加準確、有效地表征林分空間結構特征一直是現代森林經理學研究熱點，對于近自然森林經營意義重大［1-2］。

林分空間格局源于種群空間格局，指林木在水平空間上的配置狀況和分布特征，是林木生長特性、種內種間關系、環境因子等綜合作用的結果［3］。研究林木空間格局有助于理解森林的建立、更新和發展過程，是實際營林所必須測定的空間結構因子。目前，林分空間格局研究多采用與距離有關的格局分析方法，主要包括聚集指數、精確最近鄰體分析［4］、Ripley's K(d)函數分析［5］。其中，精確最近鄰體分析和Ripley's K(d)函數分析以林木坐標點圖為基礎進行格局分析，亦稱為點格局分析［4，6］，已被廣泛應用于描述林木的空間分布，可分析不同距離尺度下的種群格局及不同種群的種間關聯關系，具有顯著的優越性。但該類方法提供的是種群的統計格局，很難從中導出與森林空間經營直接相關的具體信息，且方法復雜，計算量較大。惠剛盈［7-8］提出基于相鄰木角度測量的角尺度，該方法簡單可靠，既能分析林分整體空間格局，又可表征林分精細結構即林木點結構的分布格局，是林分經營中最常用的格局分析方法［2-3］。但角尺度與距離無關，無法反映林分疏密程度及相鄰木距離等空間要素。

Voronoi圖是計算幾何的基本數據結構，按照對象集合元素的最近鄰原則將空間劃分為許多不相重疊的凸多邊形(Voronoi多邊形)。每個凸多邊形都有唯一的集合元素與之相關聯，其內任一點到該元素的距離總小于到其他元素的距離。Voronoi圖具有鄰接性、唯一性、空間動態等特性，適用于表征空間目標的復雜鄰域關系［9］。由于可根據需要對Voronoi圖定義進行擴展，Voronoi圖在計算機圖形學、分子生物學、空間規劃等眾多領域具有廣闊的應用背景［10-11］。Brown［12］最早將 Voronoi圖引入林學領域，提出林木競爭分析的 APA(area potentially available)指數;湯孟平等將Voronoi圖用于群落優勢樹種的種內種間分析［13］和混交度研究［14］;趙春燕［15］等人通過 Voronoi圖研究紅樹林空間結構指數。但應用Voronoi圖深入分析林分空間結構特征及林木分布格局的研究尚不多見。本研究選擇南洞庭湖龍虎山林場典型常綠闊葉林為對象，建立基于Voronoi圖的林分空間模型，量化林木空間結構，包括林木水平位置、空間近鄰關系、空間密度等信息。在該模型的基礎上，分析樣地林木整體空間格局及主要喬木種群的格局動態。

1 研究區概況和樣地調查

研究地點位于南洞庭湖區龍虎山林場，地理位置為東經 112°17'2″—112°18'30″，北緯 28°54'4″—28°54'56″。該地屬中亞熱帶季風濕潤氣候，年均氣溫16.7℃。林場內有木本植物100余種，其代表性林分是青椆混交林，主要喬木樹種有青椆(Cyclobalanopsis myrsinaefolia)、香樟(Cinnamomum camphora)、青岡櫟(Cyclobalanopsis glauca)、毛豹皮樟(Litsea coreana var.lanuginosa)、苦櫧(Castanopsis sclerophylla)、鉤栗 (Castanopsis tibetana)、馬褂木(Liriodendron chinense)、木荷(Schima superba)、楓香(Liquidambar formosana)等。當地資料和林分調查表明，研究地森林群落是受過輕微干擾而形成的天然次生林，青椆、香樟、青岡櫟等為當地優勢樹種［16-17］。

筆者于2010年8月在研究區內設置面積為100m×100m的標準樣地1塊。采用全站儀等設備對樣地內胸徑大于5cm的活立木進行定位、編號及每木檢尺，檢測內容為林木位置、樹種信息、胸徑、樹高、冠幅等。以樣地西南角為坐標原點，記樣地內林木的相對空間坐標為(x，y)，其中x表示東西方向坐標，y表示南北方向坐標。該次研究共涉及林木982株。

2 研究方法

2.1 Voronoi圖定義

設二維平面空間中有一離散點集P={p1，p2，…，pn}，其中任意兩點不相等，任意四點不共圓。點集中任一元素pi(1≤i≤n)根據空間鄰接關系確定Voronoi多邊形 tpi(1 ≤ i≤ n)［10］:

式中，符號‖·‖表示兩點間的歐氏距離，點pi稱為Voronoi多邊形tpi的生成元，該定義確保Voronoi多邊形tpi內任意點s到生成元pi的距離小于到其他生成元pj的距離，則由tpi(1≤i≤n)組合成的平面圖形稱為離散點集P的Voronoi圖。至此，離散點集完成對二維平面的空間剖分。

2.2 基于Voronoi圖的格局分析方法

離散點集在平面空間通常有3種格局分布形式:均勻分布、隨機分布和聚集分布，如圖1所示。點集的Voronoi圖可以直觀地反映點狀目標在平面空間的分布情況［18］。當平面點均勻分布時，相鄰Voronoi多邊形大小相當;聚集分布時，Voronoi多邊形大小各異，聚集區內的較小，聚集區外的較大;隨機分布時的情況則介于均勻分布和聚集分布之間。可見，通過量化Voronoi多邊形面積的變化關系則可確定平面點的分布格局。因此，本研究引入變異系數(Coefficient of Variation，CV)來衡量Voronoi多邊形面積的相對變化，其表達式為:

式中，σ是n個彼此鄰接的Voronoi多邊形面積的標準差，ˉs是其面積均值，xi表示第i個Voronoi多邊形的面積，ˉxi表示第i個Voronoi多邊形及其相鄰Voronoi多邊形面積的均值。按照公式(2)計算圖1的Voronoi圖的變異系數，結果依次為:0.125、0.273、0.781。

圖1 點集的空間分布Fig.1 Spatial distribution of point set

變異系數反映了單位均值上點狀目標的離散程度和分布狀態，當取值從低到高時，空間點將依次經歷均勻、隨機、聚集3種分布格局。為明確3種不同分布格局下變異系數的取值范圍，采用Monter Carlo模擬方法對平面點的分布進行顯著性檢測，具體步驟如下:模擬本研究林分樣地在100m×100m的二維平面隨機布點982個，根據點位置對二維平面進行Voronoi分割，計算每個Voronoi多邊形的面積;為研究格局分布與空間尺度的關系，采用種-面積關系研究中的巢式取樣［23］，始終以坐標原點(西南角)為起點，取樣尺度按每次4m遞增，一共獲取25個方形小樣地:4m×4m，8m×8m，12m×12m，…，100m×100m，計算每個小樣地內Voronoi多邊形的變異系數，以上過程重復99次，按99%的置信水平來確定置信區間。再根據樣地林木實際分布情況計算變異系數，計算結果如果在置信區間內，則為隨機分布;如果在置信區間外，則顯著偏離隨機分布。以上計算過程在MATLAB 7.1平臺上完成。

2.3 邊緣校正

受人為劃分樣地的影響，樣地邊緣林木的相鄰木可能位于樣地以外，因而樣地邊緣木的Voronoi多邊形是不完整的。計算這部分林木的變異系數，會使整體結果產生較大誤差，進而影響林木空間格局判斷的精度。

傳統的邊緣校正方法，如設置8鄰域樣地［19］、建立緩沖區［20］等不適用本研究。本研究的處理辦法是標記所有以樣地邊界為邊的Voronoi多邊形，在林木空間格局分析時不納入變異系數的計算中。

3 結果與分析

3.1 林分空間Voronoi圖模型

點格局分析將林木的空間位置視為平面點，從點的數量特征入手分析其空間格局［6］。然而，實際情況是森林中每株林木都要通過與周邊林木競爭來獲取必要的資源和生長空間，僅僅將林木抽象成平面點會引起林木空間信息的缺失。因此，本研究建立樣地林木水平空間的Voronoi圖模型，借此反映林木及其“影響圈”的空間分布，作為進一步分析林木空間格局的基礎。

圖2 樣地林木空間分布Voronoi圖Fig.2 Voronoi diagram of spatial distribution in sample plot

為了建模需要，忽略樣地的地形、地貌特征，視其為二維平面，樣地內所有起徑5cm以上的林木為該二維平面的點狀目標，在MATLAB平臺上編程繪制樣地林木Voronoi圖并計算其相關參數，結果如圖2所示。該Voronoi圖模型不但顯示了林木詳細的空間位置，可用于定性分析林木空間格局，還能通過該圖提取林分空間結構信息，具體如下:

(1)平面點即生成元代表林木在樣地中的空間相對位置 (xi，yi)，1≤ i≤ 982。

(2)Voronoi多邊形代表林木的“影響圈”，即林木所占據的空間生態位，反映了林木間相互作用及對環境適應的結果，其面積或周長度量了該空間的大小，其面積的倒數則表示該林木的點密度。一定范圍內Voronoi多邊形平均面積ˉs的倒數則表示該區域內林木的空間密度ρ，即:

(3)林木的生長狀況、分布格局等情況與其直接鄰體的數量、規模、距離有密切關系。兩株林木相鄰的依據是兩者對應的Voronoi多邊形共享一條邊，該邊稱為公共邊。對于非樣地邊緣林木，其Voronoi多邊形公共邊的條數代表了該林木的相鄰木數目。

(4)Voronoi圖具有唯一的對偶結構Delaunay三角網，如圖3所示，實線代表Delaunay三角網，虛線表示相對應的Voronoi圖。單個Delaunay三角形包含相鄰木間的距離信息和角度信息，其邊長長度等于參照樹與其相鄰木的距離，其內角可直接計算參照樹的角尺度［7］。

經樣地邊緣校正，本研究樣地林木Voronoi圖的有關參數及取值如表1所示。

圖3 Delaunay三角網示例Fig.3 An example of Delaunay Triangulation

3.2 林木空間格局分析

3.2.1 林分整體分布格局

本文采用基于Voronoi圖變異系數的格局分析方法研究不同空間尺度下的林木分布格局，結果如圖4所示。當小尺度取樣時，空間格局呈顯著聚集分布，20—40m為強度聚集，觀察圖2可看到該尺度下有明顯的林木聚生區域。此后，隨著尺度的增大，變異系數逐漸減小并趨于穩定。當取樣尺度大于76m時，空間格局呈接近隨機分布的輕度聚集，說明在大尺度范圍林木的空間格局有從聚集分布向隨機分布過渡的趨勢。由此可知，基于Voronoi多邊形變異系數的格局分析與取樣尺度有密切關系，尺度越大所包含樣本越多，則分析結果越接近林分真實的空間格局。

表1 樣地Voronoi圖信息Table 1 Information of Voronoi diagram in the plot

為對比基于Voronoi圖的格局分析方法，本研究選取結構化森林經營中有代表性的角尺度［7-8］方法分析同一塊樣地的空間格局。圖4為樣地中角尺度各個取值的統計情況，其均值為0.532，表明樣地林木整體上呈輕度聚集分布。可見，基于Voronoi圖變異系數方法與角尺度方法對于林分整體空間格局的判斷是一致的。圖4顯示了25種不同空間尺度下變異系數和角尺度的相關性分析結果，其相關系數R2等于0.7714(顯著正相關)，則進一步說明兩種格局分析方法可以相互替代。

圖4 林分空間分布格局Fig.4 Stand spatial pattern

3.2.2 主要樹種的空間格局動態

林木分布格局不僅可以定量描述林分的水平空間結構，而且能夠說明林分種群的動態變化。目前多以林木徑級結構代替年齡結構分析種群格局動態［21］。本研究按照胡小兵等人［22］對木本植物種群的劃分標準，以胸徑10cm為界將樣地林木分為小徑木(小樹和中樹，胸徑小于等于10cm)和大徑木(大樹，胸徑大于10cm)，分析林分內主要樹種的空間格局動態(圖5)，圖中橫坐標表示取樣尺度，縱坐標對應具體樹種的Voronoi圖變異系數值。

圖5 主要樹種的空間分布格局Fig.5 Spatial pattern of main tree species

青椆屬地帶性常綠闊葉優勢樹種，具有水源涵養和水土保持作用，更新能力強，兼有實生和萌生兩種繁殖方式，其小徑木呈強度聚集，且隨尺度增大聚集程度進一步加深。由于在群落中分布最廣，數量最多，該樹種所承受的競爭壓力也最大，因而其大徑木雖呈顯著聚集，但聚集程度相比小徑木明顯降低(圖5)。香樟是群落中另一種優勢樹種，根系發達，具有較強的適應性，其種子輕，結實量大，散播范圍廣，在群落中與青椆全方位競爭，小樹在林分內呈輕度聚集分布，而成年大樹則因自然稀疏其格局過渡到隨機分布(圖5)。青岡櫟與青椆同屬山毛櫸科，林下更新良好，種子粒大而重，萌發率高，小樹多聚生在一起。但青岡櫟在群落中的優勢地位不及青椆和香樟，由于自稀疏及環境影響，使得其成年大樹的空間格局在隨機與聚集分布之間震蕩(圖5)。毛豹皮樟是群落中的非優勢種(限于篇幅，非優勢種僅列出毛豹皮樟)，其空間格局波動較大，分布特征不明顯，這可能與該種群數量較少有關(圖5)。總的來看，樣地內小徑木的空間格局多為聚集分布，大徑木主要呈輕度聚集或隨機分布;優勢種空間分布較明顯，非優勢種分布特征不明顯或呈隨機分布;林分空間格局主要由優勢種的格局決定，各樹種分布格局的動態變化與林分整體呈聚集態，且有向隨機分布演變的趨勢保持一致。

以上分析表明，受隨機自稀疏、個體生長差異等多重因素影響，同一樹種在不同的發育階段，其空間格局會發生動態變化。同時，受空間尺度改變的影響，林木在小尺度上形成的分布格局可能會在大尺度上復合進而引起分布狀態的改變，使得種群分布動態表現出強烈的空間及時間特征。

4 結論與討論

本研究以南洞庭湖龍虎山林場典型常綠闊葉林群落為研究對象，根據林木點位置剖分樣地平面空間，建立基于Voronoi圖的林分空間模型。該模型不但能顯示出林木的空間位置和疏密程度，而且量化了林木的空間鄰近關系，具有明確的林學及生態學意義。對于以參照樹為中心的空間結構單元，其相鄰木可以由參照樹的Voronoi多邊形公共邊唯一確定，解決了長期以來相鄰木選取尺度模糊的問題。

本研究提出基于Voronoi圖的格局分析方法，應用該方法分析林分的整體空間格局。與角尺度方法的對比研究顯示，Voronoi圖變異系數和角尺度存在顯著的正相關關系，表明基于Voronoi圖的格局分析方法同樣能有效刻畫林木空間分布特征。進一步使用該方法從不同的時空尺度上描述了林分主要樹種的空間格局動態，并結合林木生長特性闡釋了格局動態變化的原因:在時間尺度上，林木的自然稀疏是引起種群格局動態的主要原因;在空間尺度上，自然稀疏的作用范圍會因尺度的改變而相應改變，使得小尺度上形成的格局會因尺度增大而復合成更大規模的格局分布，從而改變種群格局面貌。

研究林木空間關系的最終目的是采取合理的經營措施調控林分空間結構，而具體的營林效果可通過經營前后的林分Voronoi圖體現出來。例如，擇伐部分林木后，原本由被伐林木所占據的空間單元并入周邊鄰木，由此引發的空間格局變化可從林木Voronoi圖變異系數上反映出來。對于整體呈聚集分布的林分，其中必然存在某些林木聚集區域和稀疏區域，根據林分Voronoi圖的空間密度信息可以找出這些待調控區域并施加不同的經營措施(擇伐或補植)。因此，本研究結果對于天然林空間結構重構及調控具有一定指導意義。

通過相鄰木連線的垂直平分線來構造Voronoi圖，所隱含的林學假設是彼此相鄰的兩株林木處在同等競爭地位，然而實際情況中這些林木多為不對稱競爭狀態。胸徑粗大、根系發達的林木將處于優勢地位，占據更多的資源及生存空間。為使所建模型能反映林分真實情況，后續研究可將胸徑、冠幅等林木屬性作為權值建立林分加權Voronoi圖，給林木空間單元劃分提供更多科學依據。

此外，還可深入挖掘森林的空間結構信息，建立更為完善的空間模型。例如，考慮林木對光能的利用，可引入樹高因子，對林分三維空間進行Voronoi剖分;結合GIS等先進技術，建立更大研究尺度、多層次的森林空間模型。

［1］Shan Y，Bettinger P，Cieszewski C J，Li R T.Trends in spatial forest planning.Mathematical and Computational Forestry＆Natural Resource Sciences，2009，1(2):86-112.

［2］Tang M P.Advances in study of forest spatial structure.Scientia Silvae Sinicae，2010，46(1):117-122.

［3］Zhang H R，Tang S Z.Research of the Northeast Natural Forest Sustainable Management Technology.Beijing:China Forestry Publishing House，2011:81-82.

［4］Moeur M.Characterizing spatial patterns of trees using stemmapped data.Forest Science，1993，39(4):756-775.

［5］Ripley B D.Modelling spatial patterns.Journal of Royal Statistical Society.Series B，1977，39(2):172-212.

［6］Zhang J T.Analysis of spatial point pattern for plant species.Acta Phytoecologica Sinica，1998，22(4):344-349.

［7］Hui G Y，Von Gadow K，Alert M.The neighbourhood pattern-a new structure parameter for describing distribution of forest tree position.Scientia Silvae Sinicae，1999，35(1):37-42.

［8］Hui G Y，Von Gadow K，Hu Y B，Chen B W.Characterizing forest spatial distribution pattern with the mean value of uniform angle index.Acta Ecologica Sinica，2004，24(6):1225-1229.

［9］Charles D，Gilles G.Voronoi tessellation to study the numerical density and the spatial distribution of neurones.Journal of Chemical Neuroanatomy，2000，20(1):83-92.

［10］Skare ?，M?ller J，Jensen E B V.Bayesian analysis of spatial point processesin the neighbourhood ofVoronoinetworks.Statistics and Computing，2007，17(4):369-379.

［11］Kristensen L，Olsen J，Weiner J，Griepentrog H W，N?rremark M.Describing the spatial pattern of crop plants with special reference to crop-weed competition studies.Field Crops Research，2006，96(2/3):207-215.

［12］Brown G S.Point density in stems per acre.New Zealand Forestry Service Research Notes，1965，38:1-11.

［13］Tang M P，Chen Y G，Shi Y J，Zhou G M，Zhao M S.Intraspecific and interspecific competition analysis of community dominant plant populations based on Voronoi diagram.Acta Ecologica sinica，2007，27(11):4707-4716.

［14］Tang M P，Zhou G M，Chen Y G，et al.Mingling of evergreen broad-leaved forests in Tianmu mountain based on Voronoi diagram.Scientia Silvae Sinicae，2009，45(6):1-5.

［15］Zhao C Y，Li J P，Li J J.Quantitative analysis of forest stand spatialstructure based on Voronoidiagram and Delaunay triangulated network.Scientia Silvae Sinicae，2010，46(6):78-84.

［16］Lü Y，Zhang H，Wang X Y，Deng Z A，Li J J.Storey structure study of Cyclobalanopsis myrsinaefolia mixed stand based on storey index.Forest Resources Management，2012，(3):81-84.

［17］Li J J，Zhang H R，Wang C L，Zang H，Cao X P.Discussion on spatial optimization modelling of water resource conservation forests and management practice of forest functions.Journal of Central South University of Forestry and Technology，2012，32(3):23-28.

［18］Qin X W，Zhang S Q，Li X F，Na X D，Liu C Y.Spatial pattern of red-crowned crane nest-sites based on Voronoi diagram.Chinese Journal of Ecology，2008，27(12):2118-2122.

［19］Haase P.Spatial pattern analysis in ecology based on Ripley’s K-function:Introduction and methods of edge correction.Journal of Vegetation Science，1995，6(4):575-582.

［20］HuiG Y， Gadow K， HuY B.Structure-BasedForest Management.Beijing:China Forestry Publishing House，2007:153-154.

［21］Zhang H R，Tang M P，Shu Q T.Theory and Practice of Forest Ecological Harvesting.Beijing:Publishing House ofChina Forestry，2006:225-226.

［22］Hu X B，Yu M J.Size structure and distribution pattern of Cyclobalanopsisglauca population in evergreen broad-leaved forests.Journal of Zhejiang University:Science Edition，2003，30(5):574-579.

［23］Tang Z Y，Qiao X J，Fang J Y.Species-area relationship in biological communities.Biodiversity Science，2009，17(6):549-559.

參考文獻

［2］湯孟平.森林空間結構研究現狀與發展趨勢.林業科學，2010，46(1):117-122.

［3］張會儒，唐守正.東北天然林可持續經營技術研究.北京:中國林業出版社，2011:81-82.

［6］張金屯.植物種群空間分布的點格局分析.植物生態學報，1998，22(4):344-349.

［7］惠剛盈，Von Gadow K，Alert M.角尺度——一個描述林木個體分布格局的結構參數.林業科學，1999，35(1):37-42.

［8］惠剛盈，Gadow K V，胡艷波，陳伯望.林木分布格局類型的角尺度均值分析方法.生態學報，2004，24(6):1225-1229.

［13］湯孟平，陳永剛，施擁軍，周國模，趙明水.基于Voronoi圖的群落優勢樹種種內種間競爭.生態學報，2007，27(11):4707-4716.

［14］湯孟平，周國模，陳永剛，，趙明水，何一波.基于Voronoi圖的天目山常綠闊葉林混交度.林業科學，2009，45(6):1-5.

［15］趙春燕，李際平，李建軍.基于Voronoi圖和Delaunay三角網的林分空間結構量化分析.林業科學，2010，46(6):78-84.

［16］呂勇，臧顥，萬獻軍，鄧志昂，李建軍.基于林層指數的青椆混交林林層結構研究.林業資源管理，2012，(3):81-84.

［17］李建軍，張會儒，王傳立，臧顥，曹旭鵬.水源涵養林多功能經營結構優化模型初探.中南林業科技大學學報，2012，32(3):23-28.

［18］秦喜文，張樹清，李曉峰，那曉東，劉春悅.基于Voronoi圖的丹頂鶴巢址空間格局分析.生態學雜志，2008，27(12):2118-2122.

［20］惠剛盈，克勞斯·馮佳多，胡艷波.結構化森林經營.北京:中國林業出版社，2007:153-154.

［21］張會儒，湯孟平，舒清態.森林生態采伐的理論與實踐.北京:中國林業出版社，2006:225-226.

［22］胡小兵，于明堅.青岡常綠闊葉林中青岡種群結構與分布格局.浙江大學學報:理學版，2003，30(5):574-579.

［23］唐志堯，喬秀娟，方精云.生物群落的種-面積關系.生物多樣性，2009，17(6):549-559.