2017, Vol. 37
文章信息
- 刘金福, 薛凡, 黄嘉航, 林志玮, 兰思仁, 洪伟
- LIU Jinfu, XUE Fan, HUANG Jiahang, LIN Zhiwei, LAN Siren, HONG Wei
- 基于生态值的戴云山保护区地块选择研究
- Study on site selection of Daiyun Mountain Reserve based on ecological value
- 森林与环境学报,2017, 37(1): 1-7.
- Journal of Forest and Environment,2017, 37(1): 1-7.
- http://dx.doi.org/10.13324/j.cnki.jfcf.2017.01.001
-
文章历史
- 收稿日期: 2016-08-23
- 修回日期: 2016-11-26
2. 福建农林大学计算机与信息学院, 福建 福州 350002
2. College of Computer and Information Science, Fujian Agriculture and Forestry University, Fuzhou, Fujian 350002, China
建立自然保护区是保护生物多样性最有效的方法，其规划设计的科学性和实用性直接影响着保护区管理和生物多样性保护功能[1]。 我国自然保护区规划设计主要基于体系规划上，依条例和规划原则，重点结合行政区域或地理热点分区，定性分析边界划定，存在空缺区域。
随人口数量增多，人类活动范围扩大，资源减少，导致不得不用最合理面积保护最大生物多样性，合理的保护区设计成为保护区管理有效性的重要保障。 定量最优化设计能在有限地域面积中最大化地保护生物多样性。 如何科学规划和合理设计自然保护区亟待解决[2]。
运筹学为自然保护区空间最优化选择设计提供了新思路[3]，如地块选择优化设计方法，能将研究区划分为规则或不规则地块，依设定指标或规则建立数学模型，划定选择保护区域[4]。 国内早在十几年前就着手保护区最优设计研究，徐海根等[5]就内蒙古东北部地区丹顶鹤(Grus japonensis)保护区网络设计中应用最优化方法，仅以湿地面积为约束，未包含保护区空间特征。 空间特性的忽略，尤其是空间连续性，将加剧空间破碎化并造成生物多样性损失[6-7]。 保护区地块间的连续性，能保证物种在区内的正常迁移、 扩散和繁衍，提高物种的存活率[8]。 目前国内森林生态系统类型保护区地块优化设计尚无讨论，而植物相对于动物不存在因迁移导致分布不确定的问题，使得以植物分布为主的森林生态系统可能更适合利用确定的分布数据进行多物种地块选择设计，因而有必要探讨基于确定数据的地块选择法在我国森林类型保护区设计的适用性并进行实证研究。
戴云山国家级自然保护区是闽中重要的水源涵养地与发源地，现有保护区面积134.724 km2。 区内有较高的植物多样性保护水平，生态系统结构和功能较为稳定，但缺乏完善的规划建设方案[9-10]，无法满足许多珍稀物种的活动、 繁殖、 发展的生态环境要求。 为此，以福建戴云山为对象，首次将地块选择算法运用于我国森林生态系统类型自然保护区设计，旨在探讨戴云山自然保护区规划的优化设计，保证在有限资源下最大程度地保护生物多样性，为我国森林生态系统类型自然保护区设计提出新思路。
1 研究区概况戴云山国家级自然保护区位于福建省泉州市德化县境内(东经118°05′22″～118°20′15″，北纬25°38′07″～25°43′40″)，所辖赤水、 雷峰、 南埕、 桂阳、 上涌和大铭6个乡镇22个行政村，总面积134.724 km2，森林植被覆盖率93.4%，主峰海拔1 856 m，有“闽中屋脊”之称。 该区属海洋季风气候，年均气温15.6～19．5 ℃，极端最低气温-16.8 ℃，极端最高气温36.6 ℃，年均降水量1 700～2 000 mm； 该区成土母质以花岗岩为主，土壤类型主要有赤红壤、 红壤、 黄壤、 黄红壤和泥炭沼泽土等； 为我国东南沿海典型的山地森林生态系统，分布着全球面积最大的黄山松(Pinus taiwanensis)林、 丰富的兰科(Orchidaceae)植物资源及濒危动植物资源。 植被类型有温带针叶林、 针阔叶混交林、 暖性针叶林、 落叶阔叶林、 常绿阔叶林、 竹林、 灌丛、 沼泽和水生植被9种。 本区域具有南、 中亚热带过渡特征，地位特殊且重要，具有较高的保护和研究价值[11-13]。
2 研究方法 2.1 资料收集与处理图形资料主要是戴云山自然保护区所在6个乡镇林业基本图，其经营小班数据主要来自2012年度德化县二类调查数据。 采用ArcGIS10.2软件处理林分地形图及经营小班数据。 通过乡镇行政图数据转换，采用高斯克吕格投影坐标系，考虑数据和模型设计对计算效率影响，以绘制4 km2网格划分研究区，以小班地理坐标为准对网格进行生态值赋值编辑，利用连接功能将6个乡镇18 762个小班及其属性与网格匹配，获得261个网格内所有土地属性。 即把每一个包含土地属性的网格视为地块，其也是若干小班的集合。
2.2 地块选择算法保护区地块选择算法设计源于运筹学领域最优化选址的2个经典问题，即集合覆盖问题(SCP)和最大覆盖问题(MCP)，常用于城市应急设施选址优化研究[14]，后被学者引入并应用于自然保护区规划设计[15]。 SCP和MCP是以地块数量最小化和物种丰富度最大化为目标，不包含连续、 聚集和边界等空间特性[16-17]，忽略空间特征导致选择地块彼此孤立且相距很远，不利于物种迁移、 繁衍，最终影响物种的长期存活。 为克服原始覆盖模型的缺陷，OLDFIELD et al[18]和ÖNAL et al[19]基于图论提出设计连续自然保护区的“尾长法”； 王宜成[20]和DISSANAYAKE et al[21]提出Önal-Briers优化新模型，运算效率高于原模型，一定程度上解决了保护区设计需考虑空间特征问题，但未顾及地块的生态属性差异，忽视了除物种分布外的其他生态因素。 为此，基于王宜成改进模型，新增反映戴云山植物多样性的生态指标，构建线性整数规划模型：
| $\text{Mini}\sum\limits_{i}{{{U}_{i}}{{e}_{i}}}$ | (1) |
| $\begin{matrix} s.t.\text{ }\sum\limits_{i}{{{\delta }_{si}}{{U}_{i}}\ge {{k}_{s}}} & \text{对所有}s\in R \\ \end{matrix}$ | (2) |
| $\begin{matrix} \begin{align} & {{X}_{ij}}+{{X}_{ji}}\le {{U}_{i}} \\ & {{X}_{ij}}+{{X}_{ji}}\le {{U}_{j}} \\ \end{align} & \text{对所有弧线}\left( i,j \right),i＜j \\ \end{matrix}$ | (3) |
| $\begin{matrix} \sum\limits_{j}{{{X}_{ij}}}\le {{U}_{i}} & \text{对所有}i \\ \end{matrix}$ | (4) |
| $\sum\limits_{i,j}{{{X}_{ij}}}=\sum\limits_{i}{{{U}_{i}}-1}$ | (5) |
| $\begin{matrix} {{Z}_{ij}}\ge {{W}_{i}}+1-m\left( 1-{{X}_{ij}} \right) & \text{对所有}i,j \\ \end{matrix}$ | (6) |
| $\begin{matrix} {{W}_{i}}=\sum\limits_{j}{{{Z}_{ji}}} & \text{对所有}i \\ \end{matrix}$ | (7) |
| $\begin{matrix} {{U}_{i}},{{X}_{ij}}=0,1 & \text{并且}{{Z}_{ij}}, \\ \end{matrix}{{W}_{i}}\ge 0$ | (8) |
式中： i，j为独立的备选地块(网格)编号； s为物种编号； R为物种集合； Ui为二元变量，若地块i被选为自然保护区地块，Ui=1，否则Ui=0； ei为地块i生态值ei′的变换量，为满足目标函数最小化，取ei=100-ei′； δsi为指定参数，当物种s在地块i中出现时，δsi=1，否则δsi =0； ks为指定参数，表示选择保护区包含物种s最少地块数； Xij为二元变量，若有弧线从i指向j，Xij=1，否则Xij =0； Zij表示i加到j上的尾长； Wi表示i的尾长； m为任意远大于选择自然保护区地块数量的常数值。
目标函数(1) 使所选保护区面积最小且生态价值最大。 约束条件(2) 满足树种分布要求，对于物种s，选择的保护区必须至少有ks个地块包含该种。 约束条件(1) ，(2) 构成传统SCP模型。 约束条件(3) ，(4) 约束弧线形成，其中约束条件(3) 表示如果地块i没有被选择为保护区，即Ui =0，则无弧线连接地块i与其他地块； 反之，如果地块i被选择为保护区，即Ui =1，则最多有一条弧线连接i与相邻地块； 对于地块j亦如此。 约束条件(4) 说明若地块i没有被选中，即Ui =0，则没有弧线以i为起点； 反之，如果地块被选为保护区一部分，即Ui =1时，则最多有一条弧线从它发出指向相邻的地块。
约束条件(5) ，(6) 和(7) 考虑保护区空间特性，使选择的保护区形成连续图。 由图论原理可知，弧线数量等于点数量减1是图连续的条件之一，约束条件(5) 为该条件显式表达。 约束条件(6) 计算地块i尾长，等于所有与地块i相邻地块的尾长之和。 相邻地块i，j尾长Zij与Xij关系由约束条件(7) 控制。 约束条件(8) 定义Ui，Xij为二元变量，Zij，Wi为非负值。
每个网格内所代表地块的生态值e′是由每一个小班的生态值vi和小班面积比重求得，即(9) 式。 生态值是反映生态功能特征，考虑植被为森林生态系统重要组成部分，以植被作为生态值赋值对象是可行的。 依现有数据构建生态值赋分标准体系[22]，并对地块进行赋分评价，见表 1。
| ${{e}^{\prime }}=\frac{{{v}_{1}}{{s}_{1}}+{{v}_{2}}{{s}_{2}}+{{v}_{3}}{{s}_{3}}+\ldots +{{v}_{n}}{{s}_{n}}}{{{s}_{1}}+{{s}_{2}}+{{s}_{3}}+\ldots +{{s}_{n}}}$ | (9) |
| 分值Score | 地类Land | 树种组成Varieties of trees | 林分郁闭度Stand closure | 平均树高Mean height/m | 平均胸径Mean DBH/cm | 龄组Age composition | 立地等级Site quality |
| 10 | 乔木林地Arbor land | 优势种明显，多个伴生种Obvious dominating species with multiple accompanying species | ＞0.8 | ＞16 | ＞20 | 中龄林Middle-aged stand | 肥沃Fertile |
| 8 | 经济林地、 竹林地Commercial and bamboo forests | 有优势种多个伴生种Dominating species with multiple accompanying species | 0.7～0.8 | 10～16 | 14～20 | 近熟林Near-mature forest | 较肥沃级Relative fertile |
| 6 | 疏林、 灌木林地Open woodland and shrub land | 一个优势种一个伴生种Single dominating species with single accompanying species | 0.6～0.7 | 7～10 | 10～14 | 成熟林Mature forest | 中等肥沃级Mid-fertile |
| 4 | 宜林地、 非规划林地Suitable and unplanned lands | 无优势种No dominating species | 0.4～0.6 | 4～7 | 8～10 | 过熟林Over-mature forest | 贫瘠级Infertile |
| 2 | 未利用荒地、 非林地Unutilized wasteland and non-forest land | 单种Single species | 0～0.4 | 0～4 | 0～8 | 幼龄林Young-aged stand | — |
| 0 | 非林地Non-forest land | 无树种No species | 0 | 0 | 0 | — | — |
式中，e′为网格生态值，vi(i=1，…，n)为网格中所包含小班生态值，si(i=1，…，n)为小班面积。 保护区保护对象为地带性及代表性的生态系统、 珍稀濒危物种，受限于调查资料，依戴云山地区6个乡镇森林资源调查数据，选择10类树种作为规划对象： 地带性植被类型以黄山松、 罗浮栲(Castanopsis fabri)、 马尾松(Pinus massoniana)、 杉木(Cunninghamia lanceolata)、 其他硬阔叶树种为代表，濒危珍稀植物以福建柏(Fokienia hodginsii)为代表，考虑物种多样性，又选择甜槠(Castanopsis eyrei)、 米槠(Castanopsis carlesii)、 青冈(Cyclobalanopsis glauca)等为保护对象。 鉴于保护区现有面积为134.724 km2，占用35个网格，建立保护区覆盖物种地块数应至少不比原保护区少，取ks=35，即选择戴云山自然保护区至少包含35个地块。 同时，对于树种分布网格数小于35的树种及其他如闽楠(Phoebe bournei)、 秃杉(Taiwania cryptomerioides)、 南方红豆杉(Taxus wallichiana var. mairei)等珍稀树种所在的地块，可由手动添加到模型运算结果中。 设定m值为500，采用交互式通用优化求解器软件LINGO12.0，在内存6G计算机运行程序。
3 结果分析 3.1 生态值计算依公式(9) 计算戴云山地区生态值分布情况，结果见图 1。 若不考虑空间特征，地块选择法只会选中生态值高的地块，其分布极为离散，难以完整保护生物多样性。 因此设计自然保护区需考虑空间连续等特征。 基于数学规划软件LINGO编辑程序并运行，全局最优解为4 872.210，模型迭代8 262次，运行时间为84 s。 模型求解共包含135 720个变量，203 332个约束，占用内存34 092 k，输出结果见图 2。 考虑m的取值对模型求解效率有较大影响。 m取50，略大于保护物种数ks=35，在10 min内迭代超过2万次仍无结果； 而m取100，500，1 000，5 000时计算时间均在80 s左右，表明m可取远稍大于保护区地块数量值，否则可能会造成计算效率的反复。 图 2为依戴云山地区地块选择算法结果，包括了地块数小于35个树种所分布的地块，共选择88个地块，主要集中在北纬25°37′30″～25°43′30″，与戴云山自然保护区范围北纬25°38′07″～25°43′40″极为接近。 该结果是满足最小破碎化和生态价值最大化条件下的最优地块组合，地块之间并不完全连续的原因在于模型运行结果与手动添加的地块不完全重叠且存在间隔。 将生态值分布与地块选择结果叠加得到图 3。 最终选择的保护区地块集合与生态值较大的地块集合并不完全重合，即模型选择最小破碎化结果，在生态值差异不大情况下，模型优先考虑空间连续性。
|
图 1 戴云山地区生态现状 Fig. 1 Ecological status of Daiyun Mountain district |
|
图 2 戴云山地区地块选择算法结果 Fig. 2 Result of network selection method in Daiyun Mountain district |
|
图 3 生态值分布与模型计算结果 Fig. 3 Distribution of ecological value and modeling result |
自然保护区规划设计应考虑保护区区域面积和形状以减弱边缘效应，保护生态系统边缘内的多样性。 将原戴云山保护区与图 3叠加可得图 4，其中深色部分即为原保护区行政区划，边界形态很不规则，与矩形或圆形相比，将极大增加边界周长，形成明显的边缘效应，易受到破坏。
|
图 4 戴云山自然保护区模型计算结果 Fig. 4 Modeling result of Daiyun Mountain Nature Reserve |
依CorelDRAW软件绘制图 5，其中标注的5个区域为戴云山保护区备选区域，其范围确定从各局部生态值最高点出发，选取连续区域并考虑整体区域的均匀性(即边缘效应的影响)。 保护区网络优化方案依文献[23]保护区设计原则并结合树种保护状况进行设计。 从空间分布来看，区域2、 3和4间的间隔较小，其整体布局较为均匀，与区域1、 5的间隔都较远。 同时，2、 3、 4区域覆盖了10类树种分布的地块，能起到基本的保护作用。 考虑保护区网络整体布局的完整性和紧凑性以减少边缘效应作用，建议将区域2、 3和4中间环绕的区域也纳入保护区网络，同时舍掉1、 5区域及最左边下方的1个地块，不影响树种保护水平，因此共包含50个备选地块，其中44个地块为选择结果，比原保护区区域新增了13个地块，见图 6虚线内区域。
|
图 5 戴云山地区地块选择算法5个备选区域 Fig. 5 Five candidate regions of Daiyun Mountain district using block selection method |
|
图 6 戴云山地区最优地块组合 Fig. 6 Optimal combination blocks of Daiyun Mountain district |
优化模型设计中不考虑土地成本约束，主要基于我国土地归国家或集体所有，因此不存在土地购买问题的境况。 国外在规划自然保护区时通常会考虑经济预算及备选地块是否均质(单位地价相同)的问题[24-25]。 依我国实际情况，可考虑在模型中结合生态补偿替代土地购买成本规划设计保护区网络。
自然保护区体系规划是生物多样性保护的统筹设计，目的是在时间和空间大尺度上保护遗传多样性、 物种多样性和生态系统多样性。 地块选择算法根据设定的指标或规则建立数学模型，通过软件计算模型选择保护区域，探讨最大程度保护戴云山生物多样性下土地资源的最优化配置。 顺利引用国外自然保护区优化设计模型对戴云山地区进行整体规划，其设计方案只作为保护区优化设计的初始方案，最终方案的确定仍需结合实际境况与保护要求作出选择，如在保护区网络中手动添加其他珍稀树种所在的地块，考虑选择结果排除居民区地域等。 同时，在数据选择与模型计算方面尚需进一步斟酌与探讨： 1) 实地调查和3S技术相结合是自然保护区优化设计需关注的问题。 2) 后续研究可结合更详尽的科考调查资料，丰富物种类型，拓展与完善生态值赋分标准体系。 3) 考虑自然地貌、 人造设施对保护区规划的影响，可依保护区优化设计范围精确划分核心区、 缓冲区、 实验区。 4) 若将更多空间特性，即紧凑、 形状等公式化加入地块选择模型，优化设计可能更加合理，有待进一步探讨。 5) 随着计算机软件和运筹学的发展，在大规模自然保护区地块选择问题上可考虑运用并行计算技术提高模型计算效率。 6) 可考虑利用除规则网格外的其他不规则形状拟合小班，探讨精细设计自然保护区的可能性。
5 结论首次运用结合空间特性的地块选择技术设计基于森林生态系统的戴云山国家级自然保护区，选择共88个最优地块组成5个自然保护区备选区域，确定编号2，3和4为建设戴云山自然保护区最优区域组合。 88个地块并没有完全连续且有间隔，为满足约束条件下目标函数最小值，即最小破碎化自然保护区。 选择编号2，3和4备选区域拓展保护区，共包含50个备选地块，其中44块地为选择结果，总面积220 km2。 通过戴云山自然保护区优化设计，新增13个地块，北部和东部均为选择备选区域，因空间间距问题，优化保护区可向东部编号4备选区域发展。
| [1] | 薛达元, 武建勇, 赵富伟. 中国履行《生物多样性公约》二十年:行动、进展与展望[J]. 生物多样性, 2012, 20(5): 623–632. |
| [2] | 崔国发. 自然保护区学当前应该解决的几个科学问题[J]. 北京林业大学学报, 2004, 26(6): 102–105. |
| [3] | 胡晓东, 袁亚湘, 章祥荪. 运筹学发展的回顾与展望[J]. 中国科学院院刊, 2012, 27(2): 145–160. |
| [4] | 王宜成. 自然保护区空间特征和地块最优化选择方法[J]. 生态学报, 2011, 31(14): 4094–4106. |
| [5] | 徐海根, 王连龙, 包浩生. 我国丹顶鹤自然保护区网络设计[J]. 农村生态环境, 2003, 19(4): 5–9. |
| [6] | 武正军, 李义明. 生境破碎化对动物种群存活的影响[J]. 生态学报, 2003, 23(11): 2424–2435. |
| [7] | 徐高福, 胡广, 斯幸峰, 等. 千岛湖生境片段化与生物多样性研究综述[J]. 林业调查规划, 2015, 40(5): 42–48. |
| [8] | ÖNALH, BRIERSR A. Optimal selection of a connected reserve network[J]. Operations Research, 2006, 54(2): 379–388. |
| [9] | 薛凡, 刘金福, 兰思仁, 等. 戴云山国家级自然保护区植物多样性评价与可持续发展策略研究[J]. 武夷科学, 2013, 29: 6–15. |
| [10] | 梁文慧. 福建戴云山自然保护区建设初探[J]. 福建林业科技, 2004, 34(3): 132–135. |
| [11] | 刘金福, 朱德煌, 兰思仁, 等. 戴云山黄山松群落与环境的关联[J]. 生态学报, 2013, 33(18): 5731–5736. |
| [12] | 任国学, 刘金福, 徐道炜, 等. 戴云山国家级自然保护区黄山松群落类型与物种多样性分析[J]. 植物资源与环境学报, 2011, 20(3): 82–88. |
| [13] | 朱德煌, 刘金福, 洪伟, 等. 戴云山黄山松群落主要树种更新生态位研究[J]. 热带亚热带植物学报, 2012, 20(6): 561–565. |
| [14] | CHURCH R, REVELLE C. The maximal covering location problem[J]. Papers of Regional Science Association, 1974, 32(1): 101–118. |
| [15] | UNDERHILL L G. Optimal and suboptimal reserve selection algorithms[J]. Biological Conservation, 1994, 70(1): 85–87. |
| [16] | MARIANOV V, REVELLE C. The queueing maximal availability location problem:a model for the siting of emergency vehicles[J]. European Journal of Operational Research, 1996, 93(1): 110–120. |
| [17] | SHIRABE T. Districting modeling with exact contiguity constraints[J]. Environment and Planning B:Planning and Design, 2009, 36(6): 1053–1066. |
| [18] | OLDFIELD T E E, SMITH R J, HARROP S R, et al. A gap analysis of terrestrial protected areas in England and its implications for conservation policy[J]. Biological Conservation, 2004, 120(3): 303–309. |
| [19] | ÖNALH, BRIERSR A. Designing a conservation reserve network with minimal fragmentation:a linear integer programming approach[J]. Environmental Modeling & Assessment, 2005, 10(3): 193–202. |
| [20] | 王宜成. 最优化设计连续的自然保护区[J]. 生态学报, 2011, 31(17): 5033–5041. |
| [21] | DISSANAYAKE S T M, ÖNAL H, WESTERVELT J D, et al. Optimal selection of clustered conservation lands using integer programming:the case of Fort Stewart in Georgia, USA[M]//NUNES P A L D, KUMAR P, DEDEURWAERDERE T. Handbook on the economics of ecosystem services and biodiversity. Cheltenham, UK:Enward Elgar Publishing, 2014:193-207. |
| [22] | 赵惠勋, 周晓峰, 王义弘, 等. 森林质量评价标准和评价指标[J]. 东北林业大学学报, 2000, 28(5): 58–61. |
| [23] | DIAMOND J M. The island dilemma:lessons of modern biogeographic studies for the design of natural reserves[J]. Biological Conservation, 1975, 7(2): 129–146. |
| [24] | TÓTH S F, HAIGHT R G, SNYDER S A, et al. Reserve selection with minimum contiguous area restrictions:an application to open space protection planning in suburban Chicago[J]. Biological Conservation, 2009, 142(8): 1617–1627. |
| [25] | ÖNALH, WANGY C. A graph theory approach for designing conservation reserve networks with minimal fragmentation[J]. Networks, 2008, 51(2): 142–152. |