IPCC第五次评估报告指出, 由于海平面上升, 全球沿海地区将遭受越来越多的不利影响, 发生洪泛与海岸侵蚀的风险递增^[1]。世界上1/3以上的人口生活在滨海地区, 一旦这些地区遭到破坏, 将会给当地的居民和经济发展造成重大损失^[2]。因此, 制定可以适应长期气候变化的滨海防护策略是一项十分紧迫的任务^[3-4]。

盐沼、红树林、海草床等湿地是滨海生态系统的重要组成部分, 可以提供固碳、净化水质、生物多样性保育、消浪等多种生态系统服务^[5]。研究表明, 即使在风暴潮情景下滨海湿地也可以发挥消浪、缓流作用^[6-7]。传统结构型工程的建造和维护成本较高, 并且对自然生态系统有负面影响^[2]。2004年以后, 许多学者主张将基于生态系统的防护策略作为对结构工程方法的重要补充, 融入到未来海岸带保护蓝图中^{[2, 8-9]}。然而, 生态系统的消浪效度在不同时空环境下表现不尽相同, 因此, 需要深入理解滨海湿地消浪服务空间分布, 才能对海岸带空间发展规划、生态修复资源分配等重要管理决策形成有效指导^{[3, 10-12]}。

目前人类对生态系统服务的产生缺乏足够认识, 加之海陆交界处基础资料匮乏, 一些大尺度的滨海防护功能评价以湿地面积或单位面积波高衰减度、经济价值代替指示^[13-15], 这种核算方法在进行尺度转移时会产生误差。尽管海岸带工程领域学者较早开展了植物消波研究^[16], 但目前生态学家对于海洋物理环境的认知是不够的^[17]。滨海湿地的消波过程具有明显的非线性特征^[18-21], 这种非线性体现在植被带边缘通常可以消散绝大部分的入射波能, 例如在长江口地区只需要80 m宽的互花米草植被带即可消散完入射波高, 这个数值在世界其他地区甚至更低^[21]。这意味着先锋植被在其所属景观中对入射波浪的消减贡献最大, 位于先锋植被后方的其他植被对消波的贡献则较小。在评估滨海湿地消浪服务时融合其非线性特征十分必要^[22], 因为基于线性思维的评估结果会诱导管理者做出全部保留湿地或全部转化为其他土地类型的决策, 这对土地资源紧张的滨海地区来说并不是最佳的资源利用方式^[23]。

生态系统服务制图是制定符合区域生态环境特征与可持续发展政策的有利工具^[24]。已开展的湿地消浪研究多是小尺度的野外观测实验^{[10, 18]}。同时, 也有学者根据线性波理论提出植物消波模型^[16]；如在刚性植物和规则波前提下, 提出Kobayashi指数与Darlymple幂函数波高衰减模型^[24-25], 并在此基础上发展出不规则波与柔性植物相互作用的改进模型^{[19, 26]}。这些模型主要应用于断面尺度的相关研究^[23], 仅有少数学者开展了区域尺度以上的消浪功能评估^{[12, 27]}, 如Marsooli等^[12]利用MDO模型探究了美国牙买加湾盐沼湿地消浪功能, 该模型计算过程复杂, 需要充足的数据支持, 不便大范围推广^[28]。因此, 需要发展一套相对简易的消浪服务制图手段来提高公众对自然生态系统的重视。

以位于长江口的崇明岛自然滨海湿地为例, 基于线性波理论植物消波模型, 提出一个相对简易的GIS空间评估模型。从断面向整个区域扩展, 分别模拟了大、小潮升情景下崇明岛环岛滨海湿地消波服务的空间分布。研究尝试回答以下科学问题：(1)崇明环岛滨海湿地消波服务空间分布特征及其影响因素；(2)实际参与到消波服务的湿地景观规模在不同潮升条件下的变化与空间分布特征；(3)不同湿地景观类型的消波服务有何差异。研究结果可为沿海地区自然资产核算、滨海湿地保护与修复、土地利用政策制定提供科技支撑。

1 研究区概况

崇明岛位于长江入海口(121°09′30″—121°54′00″E, 31°27′00″—31°51′15″ N, 图 1), 总面积1267 km²。全岛90%以上的土地高程在3.21 m到4.20 m之间(吴淞零点基准, 下同)。环岛主海塘总长220.9 km(含江苏省兴隆沙和永隆沙海塘)^[29], 受海洋和河流动力交互影响, 南北大潮期间潮高差异接近30 cm, 大潮时长江口南、北支表层水体盐度差异超过26‰^[30]。岛屿南岸边滩多属于侵蚀或微侵蚀型淤泥滩与草滩, 而崇明岛北岸与东滩中北部的绝大部分边滩为稳定和堆积型的粉砂-淤泥滩涂^[31]。岛上已调查到的植物种类达到31科74属89种, 典型代表为芦苇(Phragmites australis)、海三棱藨草(Scirpus mariqueter)与互花米草(Spartina alterniflora), 崇明东滩是全球候鸟迁徙路上的重要驿站, 每年有290种约100万只迁徙水鸟在此过境^[32]。研究的评估范围包括环岛0 m等深线与海塘之间的所有自然潮滩湿地(图 1)。

2 研究方法与数据 2.1 光滩与植被带内波浪传播模拟

波浪在光滩上的衰减率根据Yang等^[18]和史本伟等^[33]的野外观测实验, 将大、小潮升期间光滩上的波浪衰减率分别设为0.06%/m与0.12%/m。对于盐沼部分, 采用Kobayashi指数形式波高衰减模型来模拟波高在植被带中的变化^[25]：

(1)

式中, H_o与H_x为植被带边缘的入射有效波高与距离植被带边缘x m处的波高(m)；k_v为衰减系数(m^-1)。研究假设同种生态系统所处的高程相对一致, 不同生态系统所处高程可能存在差异, 暂不考虑波浪破碎以及陡坎、潮沟等复杂因素对消浪的影响^{[11, 34]}。k_v值受多种因素影响^[34], 综合数据的可获取性, 本研究仅考虑水深与株高比值对k_v的影响。通过收集全球各地不同盐沼的野外观测结果(植株高度、水深、以及波浪经过植被带前后的有效波高), 经公式(1)推算出不同k_v值, 将k_v值与对应的水深/株高比做非线性拟合, 获取k_v值与水深/株高比的经验公式(R²=0.34, n=64)：

(2)

式中, water_depth与plant_height分别表示水深与植物株高(m)。不同植物所处水深为潮高与植物生长高程之差, 水深数据来源于近期长江口附近海图中发布的多个水文站点资料, 通过克里金插值可获取沿岸任意位置的大潮升(spring rise)与小潮升(neap rise)数据。依据前人的调查^[35], 芦苇、互花米草与莎草科植物生长最低高程分别为3.0、2.73、2.58 m, 3种植物生长季的平均株高分别为2.8、2.03、0.41 m。

2.2 不同植被带与光滩宽度提取

选取2018年长江口附近云量小于10%且处于较低潮位的Landsat 8 OLI/TIRS(2018-7-26)遥感影像, 基于ENVI 5.3中的监督分类识别不同湿地类型(莎草科植物、互花米草、芦苇与光滩), 最终分类精度大于80%。运用ArcGIS 10.4中Feature To Line工具进行面转线, 再经人工矫正得到不同植被带和光滩向海最远分界线(图 1b)。采用数字岸线分析系统(DSAS 4.0)每隔200 m生成垂直于海塘的断面(Transects), 经过筛选后保留918条, 通过断面统计不同类型植被带的宽度^[36]。

2.3 断面尺度的消浪服务核算 2.3.1 断面初始入射有效波高估算

通过2.1、2.2中获取的断面和沿岸潮升数据, 提取每个断面末端端点处的潮升数值。浅水环境中, 每条断面的入射有效波高通过公式(3)估算^{[20, 34]}:

(3)

式中, H_offshore和D_{water_offshore}为初始入射有效波高与水深(m)。

2.3.2 断面消波高度与实际参与消波断面宽度估算

本研究关注每条断面最终消波高度以及实际参与消浪的景观规模, 在评估中遵循两个原则：一, 潮升条件是否能淹到断面中每种植被带。二, 是否有波浪进入特定植被带中。以光滩与芦苇构成的断面为例介绍消波服务各要素核算方法：当光滩宽度足以消散所有波高时, 无论潮高是否可以淹到芦苇群落, 芦苇占据的断面宽度不纳入消波核算。当光滩宽度不满足消波需求时, 则要考虑潮升是否可以淹到芦苇带, 即当潮高可以淹到芦苇带时, 则通过公式(1)、(2)计算消散完剩余波高所需最低植被带宽度, 若该宽度大于芦苇实际占据宽度, 则认为该断面不足以消散完所有波高, 且该断面所有景观都参与消波。当潮高淹不到芦苇带时, 则认为波浪只能在光滩上衰减, 最终消散波高依据光滩波高衰减率计算, 而剩余波浪可能直接与植被带和光滩之间的陡坎碰撞, 并且剩余植被带不参与消浪。对于由3个以上不同生态系统构成的断面按上述算法进行类推。

2.3.3 生态系统服务综合制图

使用断面初始端点将全岛岸线分割成为多个相互连接的岸段, 采用ArcGIS 10.4中Near工具获取每个独立岸段邻近断面的消浪信息和景观结构组成。采用公式(4)、(5)计算每个断面的实际消波效率、实际参与消浪断面的宽度占比。采用基本描述性统计分析所有断面的潮高、入射有效波高、实际消波高度与实际消波断面宽度、消波效度与实际消波断面宽度比的数据分布。采用自然断裂点法将所有的岸段消波指标分成5个等级, 并统计处于该等级的岸段长度占环岛岸线总长的比例。

(4)

(5)

式中, H_i^attenuation与H_i^offshore分别表示断面i的实际消散波高与入射有效波高(m), k_i为断面i的实际消波效率。WD_i^attenuation与TW_i分别表示断面i中实际参与消波的景观总宽度与该断面总长(m), WD_i为断面i中的实际参与消波景观宽度占比。

2.4 不同滨海湿地景观的消波服务特征

主要受潮滩高程和盐度的限制, 本地湿地植被具有由光滩-海三棱藨草/互花米草群落—芦苇群落/互花米草群落演替的普遍特征^[32]。本研究由海向陆将环岛湿地景观归为7类, 每一类景观由单一或多个植被群落呈带状组合(表 1)。使用箱线图揭示7种滨海湿地景观在大、小潮升条件下的实际消波高度、消波效度、实际参与消波断面宽度分布及其差异。

3 研究结果 3.1 崇明环岛湿地消浪服务的空间分布特征

大潮升下环岛滨海湿地沿岸潮高与入射波高明显大于小潮升, 实际消波高度、实际消波效度比小潮升时分别提升0.4 m和23.2%(表 2)。大、小潮升下, 不同岸段实际消弱波高的空间分布相似, 表现为崇明岛北岸的实际消波高度要大于南岸, 小潮升时仅有38.3%的岸段消波高度处于中等等级(0.82—1.07 m), 而大潮升期间绝大部分(80.3%)岸段的消波高度均大于0.82 m(图 2)。小潮升时, 消波效率与消波高度的空间分布十分相似。在大潮升下, 南部除少数无植被覆盖岸段外, 其余岸段湿地的消波效率有明显提升, 但北岸提升不明显, 环岛消波效率在80%以上的岸段长度比小潮升时多出约1倍。

3.2 崇明环岛湿地实际参与消浪的宽度空间分布特征

大潮升时崇明环岛实际参与消浪的湿地宽度比小潮升时增加219 m, 断面宽度占比也有所增加(表 3)。崇明岛南岸实际参与消波的景观宽度总体低于崇明东滩北部, 而北部岸段实际参与消浪的景观宽度空间分布表现为两端宽中间窄(图 3)。尽管大潮升时整体入射波高增加, 但崇明北滩附近实际参与消波景观宽度增加并不明显, 大潮升期间实际参与消浪断面宽度在840 m以上的岸段长度比例比小潮升增加30%(图 3)。大、小潮升下环岛实际参与消波的景观宽度占比空间分布大体相似, 大潮升时崇明岛南部岸段的实际参与消波的景观宽度占比比小潮升有明显提高(图 3), 大潮升时断面湿地整体利用率在90%以上的岸线总长占比相对小潮升时增加15%(图 3)。

3.3 不同滨海湿地景观的消浪服务特征

结果显示, 在小潮升时, “光滩-芦苇”类消散的波高偏低, 而除“光滩”和“光滩-芦苇”类型外其他类景观消波高度均较高(>0.75 m)且消波效度接近100%, 同时“光滩”与“光滩-芦苇”实际参与消波的断面宽度普遍较低, 但“光滩”参与消波宽度占比是所有类型最高的, 接近75%(图 4)。

大潮升时, 不同景观类型的实际消波高度有明显增加, 其中“光滩-芦苇”类岸段平均消波高度比小潮升时高出0.7 m以上。除“光滩”外的其他滨海湿地景观消波效度接近100%。大潮升下, 部分景观类型实际参与消波宽度的增幅与实际消弱波高的增幅不成正比；相比小潮升, 大潮升时“光滩-芦苇”消弱波高增幅接近200%, 远高于其参与消浪景观宽度比例的增加(约30%), 这一现象在“光滩-互花米草”类中也有明显体现。除“光滩”和“光滩-芦苇”两种类型外, 其他景观断面的实际消浪利用率普遍低于50%(图 4)。

4 讨论

目前, 只有少数学者开展了区域乃至全球尺度的湿地消浪功能评估工作^{[14, 37]}, 这些研究更关注湿地在极端事件下的作用。滨海湿地对常规情况下波浪的响应似乎比极端事件更加敏感^[38], 而针对常规状态下大尺度滨海湿地消浪功能空间分布评估并不多见, 本研究通过集合GIS空间分析, 进行了初步尝试。

滨海地区物理环境因子与生态系统组成的空间分布及其相互作用异质性导致湿地消浪服务的空间分布也具有异质性^[5]。潮汐作为滨海湿地中重要的环境因子, 是影响不同岸段生态系统消波高度和效度的重要因素^{[20, 39]}。潮汐条件直接影响潮滩被淹到时的水深、范围和持续时间^[40], 水深与植株高度比对单位距离的植被消波效率有很大影响；在植物未被淹没时, 入射波高越大、植被宽度越大, 对波浪的衰减作用越明显^[41]。长江口潮升表现为口外大于口内, 崇明岛南岸潮升小于北岸^[42], 岛屿南岸多为光滩与芦苇组成的较窄潮滩, 小潮升条件下这里的潮高达不到芦苇的生长高程, 对应的岸段可能只有光滩发挥消浪作用, 但光滩对波能消耗十分有限, 因此, 在小潮升时崇明南岸的消波效度偏低。但当水深足够淹到芦苇带时, 芦苇群落将会发挥高效的消波作用。崇明环岛湿地在常规状态下平均消波效度与Narayan等^[43]的全球性荟萃分析结果接近, 表明在常规情况下, 自然湿地的消波效度与人工潜堤相当, 但对应的生态修复成本最多是修筑潜堤的1/2, 说明基于生态系统的海岸带防护是一种具有较高效益的方案^[43]。

研究表明, 一个断面上的所有生态系统并非都参与消波, 大潮升时全岛断面参与消浪的平均宽度比例为71%, 具体到不同景观类型中这个比例差异性更加明显, 如“光滩-互花米草”景观在大潮升时的实际参与比例不到37%。邻近地点的野外观测研究表明, 即使是海三棱藨草这种株高较低的植物, 消弱80%的入射波高也只需要185 m^[44], 而现实中崇明东滩海三棱藨草群落宽度远大于185 m, 这说明消浪服务在空间上具有非线性特征^[22]。植被带前缘的边界效应(boundary effect)是形成这种非线性特征的原因, 即当波浪与植被带前缘碰撞时, 波高会急剧衰减, 但随后植物对波高的衰减强度会逐渐变弱^{[16, 45]}。因此, 就消浪服务而言并不是所有的滨海湿地生态系统都会发挥作用^[5], 在进行自然资本核算时需要考虑生态系统空间位置属性对其实际价值评估的影响。

本研究可为海岸带地区生态系统管理与应对长期气候变化提供科技支持。崇明岛南岸大部分岸段长久处于侵蚀状态^[31], 迫切需要对特定岸段实施主动的生态修复措施。近10年来长江入海径流量的增加为提高植被修复成功率增添了有力保障^[46]。研究表明, “光滩-莎草科植物—芦苇”的消浪效度较高, 建议在特殊岸段的芦苇群落前扩种海三棱藨草或藨草等植物, 减轻海浪对自然岸线的冲击。在大潮升时, 断面利用率在90%以上的岸段总长度增加了15%, 这意味着一旦遭受风暴潮等极端天气, 这类岸段的滨海湿地可能无法满足防灾需求。对于此类岸段, 应该加强监控和景观布局优化, 可以在无植被岸段先通过辅助的工程措施抬升滩面后, 再逐步扩增海三棱藨草或藨草, 或在无植被覆盖的低盐度高潮滩引种一定宽度的高大乔木, 如池杉、江南桤木、旱柳等^[32]。整体而言, 环岛滨海湿地不适宜进行大规模的围垦或向外建造更高标准的海塘, 这是为确保在海平面上升时, 有足够的适宜空间供不同高程的湿地向陆迁移^[2]。建议重新评估现有堤外植被与旧海塘的结合是否可以应对极端事件^[9], 尽量减少破坏自然环境和浪费社会资源。

近海波浪多为不规则波, 在不规则波情景下, 模拟波浪在植物群落中的传播过程更加复杂^[16]。除水深/株高比外, 风速、风向、基底性质、生物量、植株密度、植株茎秆韧性等因素都会影响消浪的效果^[8]。其中, 茎秆柔韧性对消浪的影响是目前研究的热点与难点^[47], 相关报道多见于基于材料模型的室内水槽实验研究^[47-48], 受研究技术的限制, 目前难以将波浪要素的比例和植物模型的几何比例以及现实环境进行有效统一^[16], 当研究尺度扩展到区域乃至全球时, 许多关键数据获取和参数量化更为困难^{[14, 28]}。本研究在向区域尺度应用时, 改进了前人在使用模型时未考虑水深和景观组成对消波影响的缺憾。研究所采用的k_v经验公式综合了全球各地具有不同生物量、植株密度的盐沼野外观测数据, 以弥补单一物种观测样本不足对参数拟合结果的影响。此外, 研究区内3种植物的高度和生长高程差异十分明显^[35], 使用该模型可以在一定程度上体现3种植物群落在消浪效度上的相对差异。因此, 该方法值得在更大空间尺度的研究中进行改进、验证和应用。河口海岸地区水动力条件与湿地植物的生理生态特征在不同季节表现并不一致。研究中使用的植被参数均在夏季采集, 理论上是一年中植物消浪最有效的时期。目前仅有少数学者开展了植物在其他季节的消浪观测^{[4, 44]}；Schoutens等^[4]在德国易北河河口的研究表明, 当地冬季海岸具有更强的防浪需求。长江口地区冬季平均波高较其他季节要高^[42], 这对应用基于自然的滨海防护方案提出新的挑战, 需要后续的研究来填补这项空缺。

5 结论

本研究结合GIS与Kobayashi指数形式波浪衰减模型, 分别模拟了大、小潮升情景下崇明环岛湿地消浪服务的空间分布, 探究不同滨海湿地景观的消波贡献差异, 填补了前人忽略对常规状态下消波服务空间分布研究的空缺。结果表明, 湿地消浪服务的各项指标随潮升的增加而增强, 大、小潮升条件下, 环岛湿地的平均消波高度分别为0.94, 0.54 m与效度分别为83.6%、60.4%, 崇明环岛湿地消浪服务的空间分布是物理环境与生物因子的空分布及其相互作用异质性的结果, 湿地消波高度大体表现为北岸高于南岸, 消波效度与消波高度的空间分布仅在小潮升下相似, 大潮升时环岛绝大部分岸段消波效度达到80%以上。滨海湿地消浪服务在空间上具有明显的非线性特征；在大、小潮升条件下, 环岛各断面实际参与消波的平均断面宽度比分别为70.8%和55.1%, 该指标在不同景观类型之间的差异十分明显。研究结果为相关生态系统服务评估和自然资本的精细化核算提供了依据, 也为盐沼湿地生态修复规划提供了参考。