中国北方砂栖纤毛虫的空间分布及其驱动因素

陈欣怡，粱月琴，郑伊珊，吴成倩，许媛*

( 1. 华东师范大学 河口海岸学国家重点实验室，上海 200241；2. 华东师范大学 生命科学学院，上海 200241)

1 引言

了解生物在时空上的分布模式和驱动机制是生态学研究的主要目标之一，也是制定保护规划的基础[1]。对于大型动植物，现有研究已明确环境、历史和空间距离在其生物群落组成中的重要作用[2-4]。反观微型生物，由于自身体积微小、种群密度高，其分布模式和驱动机制可能与大型生物有着根本不同[5-8]。众多研究表明，局地环境变量（如pH、食物资源和栖息地的可利用性等）对微型生物多样性和群落结构的驱动比空间距离更为重要[9-11]。而在沿海生态系统中情况则更为复杂，海水的流动性和连通性促使生活在其中的微型生物群落在邻近区域高度分散。这种“强扩散”效应协同环境过滤和空间距离，在沿岸海域共同影响着微型生物群落结构[12]。

砂质海岸遍布全球，是大多数温带和热带海岸的主要类型[13]。在环境污染和全球变化的压力下，沙滩对人类生产生活的重要性与日俱增[14]。明确其中重要生物组分的地理分布模式，对评估沙滩生态系统功能和健康至关重要。由于临海界面不断发生着风、浪、潮汐和沉积物的交互作用[15]，海滩是一个动态系统，其不稳定的生境直接形成了一个天然“环境过滤器”。同时，砂质潮间带在沿海区域的连续分布以及海水流动交换促进微型生物在邻近海域的扩散，为进一步区分微型生物空间分布的影响因子提供了理想研究条件。砂栖纤毛虫作为砂质潮间带优势物种[16-18]，是微型食物网能量传递的关键环节[19]，并且其细胞体积相对较大，分类学和生物学特征的文献记载相对丰富[20-21]，有利于物种鉴定和生物多样性的计算。以往研究虽已总结出纤毛虫群落与沉积物粒度[22]、水环境盐度[23]和富营养化[24-25]间的相关性，但一贯认为物理环境（海滩形态）是砂质海滩群落的主要决定因子[26]。除沉积物粒度外，表征海滩形态动力学特征的潮汐范围和海滩坡度等因子均较少被考虑，仅在亚北极潮间带的微型生物群落研究中有所涉及[27]。然而，在较大地理尺度上对海洋砂栖纤毛虫空间分布模式的研究则更加匮乏[28-31]。

我国北方砂质海滩分布广泛，跨越黄海、渤海两个海区，其间大量河流纵横交汇，驱动着悬沙、营养物质及有机质的分布，不同海滩表层沉积物粒度和水动力条件存在明显差异[32-34]，为砂栖纤毛虫群落研究创造了天然条件。在以往的黄海、渤海纤毛虫群落生态学调查中，研究对象往往是浮游类纤毛虫，底栖类群得到的关注相对较少[35-37]。对砂栖纤毛虫的研究尺度也多局限于单个海滩，未从海域尺度对其群落的空间分布和驱动机制进行过深入探究。因此，本研究对黄海、渤海沿岸14个砂质潮间带进行采样调查，比较了不同海域砂栖纤毛虫群落的多样性和物种组成，明确了纤毛虫群落的空间分布模式。通过进一步分析其与环境、空间距离的相关性，揭示区域尺度上环境过滤和空间距离对微型生物群落组成的相对重要性，为砂质潮间带环境保护与修复提供必要的基础数据。

2 材料方法

2.1 研究区域及样品采集

本 课 题 组 于2020年10月13-23日，在34°46′～40°36′N，119°07′～122°36′E范围 的14个砂质潮间带进行底栖生物调查（图1）。采样方法主要遵循文献[23]。

关于海滩形态动力学特征，最大潮差（maxTR）通过潮汐记录获取，为每个海滩全年的月最大潮差均值。海滩坡度（Slope）定义为一个无量纲比值[27]，公式为

式中，H（单位：m）是从潮汐记录中获得的中潮位与低潮位潮差（图1）；中潮线与低潮线之间的距离通过卷尺现场测得（D1，图1），而水平方向的长度为D2=（D12-H2）0.5（单位：m）（图1）。

本研究的所有样品采集统一在潮位退至中潮线后两小时开展，并在半小时内完成。在每个海滩上，平行于海岸线设置3个间距为100 m，面积为1 m2的样方。每一样方内，使用内径约2.2 cm的改良注射器采集表层2 cm芯样5管（共计20 mL），与等体积甲醛溶液（终浓度为1%）混合，用于砂栖纤毛虫的鉴定与计数[38]。同时，每个样方收集表层0～2 cm砂质沉积物约80 g，-20°C保存，用于测量沉积物特性，即中值粒径（Sediment Grain Size, SED）及总有机碳（Total Organic Carbon, TOC）、总氮（Total Nitrogen, TN）和叶绿素a（chlorophylla, Chla）含量。使用孔径为0.45 μm的Millipore注射器滤膜过滤每个平行样方内现场挖穴收集的渗出海水（间隙水），滴加1.5‰（v/v）饱和HgCl2，以进一步分析间隙水环境的营养盐含量（氨氮、硝酸盐、亚硝酸盐、磷酸盐和硅酸盐）。此外，间隙水盐度（S）、pH和氧化还原电位（Oxidation Reduction Potential, ORP）通过折光盐度计（MASTER-S/Millα，ATAGO）、便携式手持pH计（pH400型）和氧化还原电位仪（FJA-6）现场测量获取，表层沉积物温度（T）由便携式手持pH计所记录。本次研究共采集生物样品、沉积物样品及水样各42份。

图1 研究站位分布及海滩坡度计算示意图Fig. 1 Schematic illustration of study station distribution and beach slope computation

2.2 样品处理及环境参数测定

使用 Ludox-QPS 法对砂栖纤毛虫进行制片观察[38]，依据文献资料完成鉴定与计数工作[20,39]。

在超声波分散仪分散后使用激光衍射粒度分析仪（BECKMAN COULTER LS13320）测量沉积物SED。利用元素分析仪（Vario EL CHNOS）测定沉积物TOC和TN含量。Chla含量通过叶绿素荧光仪（Trilogy Laboratory Fluorometer 7200-000）进行测定。连续流动分析仪（SKALAR Sanplus System, 荷兰）用于测定孔隙水样品中的5种营养盐（）含量，并以亚硝酸盐、硝酸盐及氨氮含量之和作为孔隙水的溶解无机氮（Dissolved Inorgani c Nitrogen，DIN）含量。

2.3 数理统计与分析

数据分析前，将各沙滩3个采样点的环境参数和纤毛虫数据进行平均，环境数据除ORP和S左偏，使用反向 lg（xmax-x+1） 转换，其余均为正向 lg（x+1） 转换，同时进行标准化；物种数据进行二次方根转换。

2.3.1 多样性参数

本研究使用PRIMER v7.0计算得到各样点的物种数（s）、个体数（N）和Pielou均匀度指数（J′）对群落多样性状况进行描述。

2.3.2 环境因子在黄海、渤海差异分析

运用单因素方差分析（One-way ANOVA）比较环境因子在两海域间的差异；采用非度量多维尺度分析（NMDS）和相似性分析（ANOSIM）对不同海域沙滩纤毛虫群落结构进行比较。

2.3.3 驱动因素分析

将环境因子与生物多样性指数进行皮尔逊相关性检验，通过生物-环境（BIO-ENV）分析，探寻影响纤毛虫群落组成的最佳环境因子组合。Mantel和偏Mental检验用于确定群落组成的非相似性与空间距离矩阵的关系以及群落组成的非相似性与环境变量差异矩阵间的关系。与Mental分析不同的是，偏Mental分析是在控制环境距离矩阵（或空间距离矩阵）的z轴后，度量纯空间距离（或纯环境效应）对物种组成非相似性距离矩阵的影响。本研究中局域环境变量由 BIO-ENV 分析筛选出的最优组合构成，而空间变量矩阵则采用各沙滩采样点的经纬度均值。在进行Mental与偏Mental分析时，bioenvdist函数与earth.dist函数分别用于计算环境距离矩阵和空间距离矩阵。上述分析方法除皮尔逊相关性检验和单因素方差分析于SPSS v24.0中进行，Mental及偏Mental检验利用R v4.0.0中的“vegan”软件包完成，其余分析均在PRIMER v7.0中进行。

为更好地分析砂栖纤毛虫功能群组在环境梯度上的变化趋势，根据Xu和Soininen[9]在2019年对潮间带纤毛虫功能群组所做的分类，本研究将纤毛虫划分为藻食性、菌食性和肉食性，并使用模糊编码程序来说明单个分类群出现杂食性的情况[21]。有关食性的本底数据主要来源于描述物种的原始资料、专家意见和文献记录[20,40]。根据物种丰度数据，使用R v4.0.0的“FD 1.0.12”软件包计算每个海滩各食性类别纤毛虫的群落加权平均数指数（Community Weighted Means, CWM），而后分别对每个沙滩的DIN含量与藻食性和菌食性纤毛虫功能组的CWM进行线性拟合。

3 结果

3.1 环境特征

本次调查的14个砂质潮间带绝大多数为典型海洋栖息地（盐度约为30，pH呈碱性，在7.1～8.6之间），个别采样站位（1号、6号和14号）由于周围河流淡水注入的影响，盐度降至20～25。此次调查样点均处温带地区，且为一次性连贯取样，表层沉积物的采样瞬时温度取值在13～21°C。除12号站位的沉积环境为弱氧化环境（ORP＜400 mV）外，其余站位均为氧化环境（ORP＞400 mV）。各站位沉积物的有机质和叶绿素a含量普遍较低，但2号和8号站位的Chla、TOC和TN含量较高。在Pearson相关性分析中，发现TN、TOC及Chla含量存在极显著的正相关关系（p＜0.01）。间隙水中的5种营养盐含量走势基本一致，其中溶解无机氮含量主要受氨氮含量的控制，4号和14号站位的DIN、氨氮、硅酸盐和磷酸盐含量明显高于其他沙滩。在各形态动力学因子中，SED和坡度的变化趋势趋于一致，而各海滩的最大潮差从1.1 m至5.3 m不等（图2）。

单因素方差分析结果表明，亚硝酸盐、磷酸盐和最大潮差在两个海域间具有显著差异（p＜0.05），其中亚硝酸盐和最大潮差均为黄海区显著高于渤海区，而磷酸盐则相反（图2）。

3.2 黄海、渤海群落组成比较

本次调查共检获纤毛虫105种，隶属26目65属，各站位种数从7种至50种不等。其中帆口目（34.0%）、小胸目（22.8%）、环毛目（10.9%）、斜头目（7.4%）和侧口目（7.1%）纤毛虫的种数最多，由类群组成图也可看出这5类纤毛虫是此次调查中构成砂栖纤毛虫群落最主要的类群，达到总量的82.2%（图3）。

从生物多样性的空间分布格局可获悉，2号和8号站位的砂栖纤毛虫物种数较高，分别为50种和35种，其个体数也远高于多数样点，达到了158 个/mL和200 个/mL。而纤毛虫个体数最高值（239 个/mL）则出现在4号站位。在群落均匀度方面，则是物种数和个体数均较低的3号（15种，10 个/mL）、7号（14种，7 个/mL）和9号（11种，4 个/mL）站位物种分布较为均匀（J′＞0.9）（图4）。

图2 环境因子趋势图Fig. 2 Environmental factors trend graph

图3 各样点砂栖纤毛虫类群组成Fig. 3 The community composition of the benthic ciliates at various sampling sites

从群落结构的研究看，砂栖纤毛虫的丰度优势种依次为 简单小胸虫（Microthorax simplex，10.8%）、蚤状中缢虫（Mesodinium pulex，10.7%）、尖梭刺膜袋虫（Acucyclidium atractodes，10.3%）、柠檬镰膜袋虫（Falcicyclidium citriforme，9.5%）和乳突盘毛虫（Discotricha papillifera，6.5%）。这5个物种的丰度累计值达47.8%，将近所观测纤毛虫总丰度值的一半。且这5种纤毛虫并非在各站位均匀分布，而是在某特定站位聚集性出现。例如简单小胸虫在8号站位能够达到95个/mL；蚤状中缢虫在10号及12号站位分别能够达到53 个/mL和38 个/mL；尖梭刺膜袋虫、柠檬镰膜袋虫和乳突盘毛虫更是在4号站位大量出现，是导致4号站位纤毛虫个体数在所有站位中位列第一的主要贡献种。此外，瞬闪膜袋虫（Cyclidium glaucoma）、细纤刺叶虫（Kentrophyllum fibrillatum）和斯坦楯纤虫（Aspidisca steini）虽然丰度值较低，但在此次所调查的站位中广泛分布，在超过70%的站位中均有发现。

图4 砂栖纤毛虫多样性的空间分布Fig. 4 Spatial distribution of diversity of benthic ciliates

黄海、渤海砂栖纤毛虫群落相对丰度累计值达90%的优势种名录由表1列出。从表中可以看出，渤海的优势种为尖梭刺膜袋虫、 乳突盘毛虫和柠檬镰膜袋虫，黄海的优势种为蚤状中缢虫和简单小胸虫。单因素方差分析发现，两海域砂栖纤毛虫群落多样性指数并无显著差异（p＞0.05）。且基于Bray-Curtis距离的NMDS和ANOSIM结果显示，两海域的砂栖纤毛虫群落相似性较高，并不根据所处海域形成鲜明分组，其群落组成也无显著差异（p＞0.05）（图5）。

3.3 砂栖纤毛虫多样性及群落结构的驱动因素解析

皮尔逊相关分析显示，纤毛虫群落的物种数和个体数都与Chla含量呈显著正相关（p＜0.05），此外，物种数还与TOC含量呈极显著正相关（p＜0.01），与SED呈极显著负相关（p＜0.01）（表2）。

BIO-EVN分析结果显示，纤毛虫群落结构与环境因子相关性较高，能够解释其物种组成的最优环境因子组合是S、SED和Slope，相关系数为0.613（p＜0.05）。次要的影响因子还有DIN含量、含量和maxTR （表3）。

Mental 检验表明，纤毛虫群落物种组成的非相似性距离与环境距离矩阵显著相关（R=0.626，p＜0.01），并且控制空间距离的效应对二者的相关系数影响不大。但空间距离矩阵对纤毛虫群落的变异影响并不显著（R=0.167，p＞0.05），在控制环境异质性后也无明显关联（表4）。

4 讨论

4.1 局域环境影响

4.1.1 沉积物特性及沙滩形态学特征

在所有局域环境因子中，沉积物特性（SED、TOC含量、Chla含量）在预测纤毛虫多样性方面发挥着重要作用，而海滩形态动力学因子（Slope、maxTR）则更多地影响了砂栖纤毛虫的群落组成。

沉积物颗粒大小一直是预测沙滩物种丰富度的主要因素，因为粒径的大小和形状直接决定着砂栖纤毛虫的可生存空间。如若空间过小，大量沉积物堵塞会阻碍溶解氧扩散，进一步减少生存空间；而沙隙过大，淤泥、黏土内容物含量降低，不足以维持大量纤毛虫的代谢消耗[41]。已有研究表明，SED在120～250 μm的海洋细砂沉积物孕育着最多样化和最丰富的纤毛虫群落[22]，而本研究结果也支持这一结论（图6a，图6b）。

SED成为影响砂栖纤毛虫多样性关键因素的另一个原因，可能是与其他环境因子的相关关系。例如，沉积物颗粒大小总是与TOC和TN含量呈显著负相关[42-43]。在本研究中，SED与TOC含量的显著负相关关系也符合这一结论（R=-0.56，p＜0.05）。这是由于粗砂在退潮后的保水能力和对有机质的吸附力较差，附着其上的有机碎屑残余较少；而细砂表面覆有水膜，在退潮后能够保持浸润状态，因此有机质不容易

被冲刷殆尽，能够留存较长时间，为硅藻和其他光合作用生物的生长提供了条件，以保证纤毛虫食物来源的供应[41,44]。这也在一定程度上解释了TOC和Chla含量（表征底栖藻类生物量）与纤毛虫群落多样性指数的显著正相关关系。

表1 砂栖纤毛虫群落丰度优势物种在黄海、渤海的组成及分布Table 1 Composition and distribution of benthic ciliates abundance-dominant species in the Yellow Sea and Bohai Sea

图5 黄海、渤海砂栖纤毛虫群落NMDS排序图Fig. 5 A ranking map of benthic ciliate communities based on the NMDS in the Yellow Sea and Bohai Sea

表2 砂栖纤毛虫多样性指数与环境因子的皮尔逊相关分析Table 2 Correlation investigation of benthic ciliate diversity indices with environmental conditions using Pearson’s correlation analysis

在砂栖纤毛虫的群落构建中，环境因子SED、Slope和maxTR起到了重要作用。这验证了有关沙滩群落生态学的两个假说。一是“生境恶劣性假说”，主要与SED和Slope这两个环境因子有关。海滩较陡的坡度和较大的沙粒都会增强潮间带的水动力强度，生物体需花费更多的时间和能量成本来躲避急流，维持生存，从而摄食时间减少[26,45]。因此，陡峭的粗砂海滩是一个严格的“环境过滤器”，只有具备足够吸附力、身体外壳坚韧（如 楯纤虫属（Aspidisca））或形状相对圆整、扁平（如帆口虫属（Pleuronema）、膜袋虫属（Cyclidium）和急游虫属（Strombidium））的纤毛虫才能够定植[27]。在本研究中7号及9号站位坡度较大，其群落的物种数和个体数也较低。二是“生境可利用性假说”，与潮汐范围（即潮差）有关，本研究中以maxTR指代。较大的潮差拓宽了海滩，改变了潮间带轮廓[26]，增加了生境的异质性和可用性，以支持更多物种的定居和共存[46-47]。此外，更大的潮间带面积也能够支持更大的种群规模，可以减少物种灭绝的概率[26]。因此，海滩的形态动力学特征与微型生物群落结构密切相关。

表3 与群落结构相关性最密切的环境因子组合Table 3 The most closely connected environmental parameters with community structure

表4 纤毛虫群落的非相似性与环境和空间距离间的Mantel检验和偏Mantel检验Table 4 Mantel test and partial Mantel test for non-similarity between ciliate communities and environmental and spatial variables

4.1.2 间隙水环境状况

间隙水环境状况主要指间隙水的盐度和营养状态。在本研究中，间隙水盐度虽然与群落的多样性并没有显著的相关性，但却是决定纤毛虫群落物种组成的一个关键因子。这与Mazei等[48]在白河河口的实验结果一致。盐度差异会影响各种生物体的能量消耗和代谢[49]，室内实验也表明纤毛虫可能具有遗传上固定的耐盐范围[50]。顺盐度梯度变化，不同耐盐范围的纤毛虫进行物种更替[51]，因此在多样性上并无明显变化。有关盐度对群落结构的主导作用多见于河口区域的研究中[52]，因为河口地带通常涵盖了强烈的盐度梯度。本研究中海滩的盐度均值为29.6±7.5，由于部分海滩周围有淡水河流注入，盐度变化范围为22～33，跨度较大，因此群落结构对盐度变化的响应较为明显。

图6 本研究中沉积物中值粒径、孔隙水的溶解无机氮含量与砂栖纤毛虫群落物种数、个体数、藻食性和菌食性功能群关系的散点图Fig. 6 Scatter plots of SED, pore water DIN in relation to the number of species, individuals, algivores, and bactivores functional groups in this study’s benthic ciliate community

另一方面，间隙水DIN含量是构建砂栖纤毛虫群落的次要因素。以往研究中也发现了氮对沿海地区纤毛虫群落组成的重要作用，这主要是通过改变不同食性纤毛虫所占比例来实现的[23,53-54]。间隙水DIN含量与藻食性和菌食性纤毛虫功能组的线性拟合结果也表明，藻食性纤毛虫随着富营养化的发展呈现负向变化，而菌食性纤毛虫则相反（图6c，图6d）。在本研究中除个别站位DIN含量极高，其余站位DIN含量的变化范围不大，取值区间不够连续可能是导致拟合结果不够显著的原因之一，但该图仍能体现出明显的变化趋势。在DIN含量较高的4号和14号站位可观测到明显的优势类群：4号站位的优势种为隶属于小胸目的兼顾藻食性和菌食性的乳突盘毛虫；而14号站位的优势种为隶属于帆口目的以菌食性为主的纤毛虫如维氏帆口虫（Pleuronema wiackowskii） 和艾氏裂纱虫 （Schizocalyptra aeschtae）。一般认为，氮是沿海地区藻类生长和富营养化的关键限制因子[55]。在富营养化初期，藻类暴发性增殖，形成了有利于藻食性纤毛虫生长的环境；而当富营养化进行到一定程度时，大量藻类死亡分解，则促进了菌食性纤毛虫的增长。在本研究中，整体上随着DIN含量的增加，纤毛虫群落的食性组成呈现出菌食性逐渐替代藻食性功能群的趋势。因此我们认为，沿海富营养化程度会通过影响藻食性和菌食性纤毛虫的比例而改变潮间带砂栖纤毛虫群落的组成，并且随着富营养化的发展，纤毛虫组成逐渐从藻食性转向菌食性。

4.2 纤毛虫地理分布驱动机制

微型生物群落空间分布驱动机制的争论焦点在于局域环境和空间距离的相对重要性[56]。某一物种在特定生境的缺失，可能是由于当地环境不利于其定居繁衍，也可能是由于自身扩散限制未能够到达该生境。物种在生境中的生存能力主要取决于其自身特性与环境因子的适配程度，而扩散能力则与观测尺度和生态系统类型有关[9]。例如，淡水和陆地生态系统存在明显的地理屏障隔离，因此其中的生物往往表现出强烈的限制性分布模式和局域性的空间结构[57-58]。而海洋是相互连通的，这使得海洋生物有可能在全球范围内扩散，因此预测海洋物种比陆地和淡水物种更不容易受到扩散限制，即地理障碍和历史因素对海洋群落的影响较弱[56]。一般而言，随着研究尺度的逐渐增大，局域环境因子的解释量会逐渐降低，空间距离的解释量将有所提高。而在本研究中，站位虽然涵盖了较大的地理范围，但局域环境因子的变化在解释砂栖纤毛虫群落结构的变异中依然占据主导地位，且解释率高达62.6%，远高于空间距离的解释率（16.7%）。另外，尽管ANOVA显示渤海和黄海的部分站位环境因子存在显著差异，但两海域纤毛虫的生物多样性和群落组成却表现出相似的模式。例如瞬闪膜袋虫、细纤刺叶虫以及斯坦楯纤虫等物种，在黄海、渤海区域相距较远、环境条件并不一致的站位中均有分布。我们猜想这可能与两点原因有关。其一，由于相同空间尺度下，单细胞微型生物受到的扩散限制相比多细胞生物较弱[59]，海域可能并非解释砂栖纤毛虫群落的最佳空间变量。在更大尺度的世界海洋生物地理分区中，黄海、渤海隶属于同一个海洋生态区，因此我们可以推测二者共享着同一个区域物种库[60]。许多物种在黄海、渤海区域广泛分布，未受扩散限制，具体站位的砂栖纤毛虫群落更多地是当地环境条件筛选后的结果。其二，在海洋沿岸系统中存在“强扩散”现象，即海水的连通性和流动性有利于活动力较强的物种扩散，高扩散率使得物种可以出现在环境条件不够理想的地方，从而导致邻近斑块之间部分同质化[61-63]。扩散限制通常随着样点之间空间距离的增加而增加，与之相反，“强扩散”效应的重要性可能随着样点之间距离的减小而增加[12]。根据本研究尺度观测结果，“强扩散”效应和局域环境条件在中国北方砂栖纤毛虫群落的组成中都相当重要，且前者在一定程度上掩盖了环境异质性影响，能够解释局域环境和空间距离所不能解释的部分。

5 结论

本研究分析了中国北方砂栖纤毛虫的空间分布格局及其驱动因素，揭示了：（1）沙滩沉积物特性和形态动力学特征是预测纤毛虫物种多样性和群落结构的良好因子；（2）沿海富营养化程度通过改变潮间带砂栖纤毛虫群落中藻食性和菌食性纤毛虫的比例，进而影响其群落结构；（3）黄海、渤海区域间，环境过滤作用和“强扩散”效应的影响强于空间距离。本研究结果为未来全球变化背景下，潮间带生物多样性的保护和生态系统管理提供了更微观精细的视角。

致谢：感谢张嘉伟、张淑坤和王苗勋同志在野外采样中的帮助，感谢范鑫鹏教授在物种鉴定工作中的指导。