土壤环境DNA（eDNA）破译景观尺度植物覆盖与生物多样性：与传统调查方法的比较与整合

本研究通过对比789个土壤样本的高通量扩增子测序（eDNA）与传统实地调查，评估了土壤eDNA在描述植物群落、多样性指标、指示类群及植被分类预测中的潜力。结果表明，eDNA方法在属级分类分辨率下能提供与实地调查相似的物种丰富度和组成信息，但在指示类群和植被群落分类预测方面敏感性较低。研究揭示了eDNA方法在景观尺度生物多样性监测中的优势与局限，为未来整合分子与传统方法提供了重要依据。

　　

1 引言

面对日益加速的气候变化，景观尺度的生物多样性监测对于检测物种组成、生态系统功能和栖息地健康的变化至关重要。传统的国家尺度生物多样性调查需要大量野外生态学家对数百个地点进行访问，以进行一致且准确的植物物种测量。然而，由于资源可用性（人力和财力）的限制，调查的规模和采样强度往往受到制约。环境DNA（eDNA）监测技术具有简化野外工作、实现大批量样品处理和高分类分辨率的潜力，为大规模生物多样性调查提供了新的解决方案。

eDNA分析已在许多难以监测的环境中提供了深入的见解，并协助评估了微生物和宏观生物的多样性与分布。然而，基于土壤eDNA进行植物群落分析的准确性在较大尺度上尚未得到充分验证。本研究基于Fahner等人的工作，将性能最佳的分类标记应用于国家景观尺度上，评估土壤eDNA表征当地植物群落的能力。研究分析了789个土壤样本，使用高通量扩增子测序技术，并与同地进行的实地调查数据进行了比较，包括DNA为基础的多样性指标、指示类群、预测植被类别和植物覆盖度。

2 材料与方法

2.1 植被调查

调查人员按照已发布的指南，作为2007年英国乡村调查的一部分进行了植被调查。该调查涵盖了英国各土地类别的多个地点，使用了1?m²植物物种记录（巢级0）。调查在距离土壤样本采集位置至少1米、最多2.5米的范围内进行。

2.2 土壤采集

土壤采样与植被调查同时进行。使用直径5厘米的干净、未使用过的塑料管从每个采样点采集顶部15厘米的土壤芯样本，并将芯样密封在预标记的塑料袋中以防止样品间土壤残留物的转移。芯样在采集当天转移到实验室，并进行多项分析。用于分子工作的芯样在-20°C下冷冻保存，以供后续处理。芯样经轻微解冻后，从有机层下方（排除细根）收集子样本，均质化后在-20°C下存档，用于后续DNA提取和植物ITS2扩增子测序。

2.3 植被分类

基于调查员对植物物种覆盖度的估计，对每个采样点应用了聚合植被分类（AVC）。AVC的确定遵循Bunce等人的方法；简要来说，植物调查数据经过整合并进行多变量分析，使用DECORANA进行排序。样本在八个组中的聚类形成了分类的基础，这八个AVC类别定义为：（1）“作物和杂草”，（2）“高草和草本”，（3）“肥沃草地”，（4）“贫瘠草地”，（5）“低地林地”，（6）“高地林地”，（7“沼泽草地镶嵌”，和（8）“荒原和沼泽”。

2.4 植物ITS2的分子分析

从0.25克存档土壤中提取DNA。土壤使用预先灭菌的装置称重，并放入Powersoil DNA 384分离试剂板中，按照制造商的说明提取DNA。样本随机分布在提取板中，每个提取板包含阴性提取对照。使用两步扩增方法生成扩增子，使用Illumina Nextera标记的ITS2区域引物，正向引物ITS2-S2 ATGCGATACTTGGTGTGAAT和反向引物ITS4 TCCTCCGCTTATTGATATGC，遵循Fahner等人的方法。每个引物在5'端修饰了Illumina预适配器和Nextera测序引物序列。使用高保真DNA聚合酶生成扩增子。PCR条件如下：初始变性95°C 2分钟，之后95°C变性15秒，55°C退火30秒，72°C延伸30秒，共35个循环。最后在72°C进行10分钟的最终延伸步骤。使用Zymo ZR-96 DNA纯化试剂盒纯化PCR产物。在第二次PCR扩增步骤中添加了MiSeq适配器和8核苷酸双索引条形码序列。PCR条件如下：初始变性95°C 2分钟，之后95°C变性15秒，55°C退火30秒，72°C延伸30秒，共8个循环，最后在72°C进行10分钟的最终延伸。通过凝胶电泳验证提取和PCR阴性对照样本为阴性。使用Agilent 2200 TapeStation系统确定扩增子大小。使用SequalPrep标准化板试剂盒对文库进行标准化，使用Qubit dsDNA HS试剂盒进行定量，并合并文库。合并后的文库在变性和中和后稀释至400 pM。变性通过0.2N NaOH处理5分钟实现，随后用0.2N HCl中和。然后用HT1缓冲液和10%变性的PhiX对照文库将文库稀释至12 pM的加载浓度。最终变性通过96°C加热2分钟实现，随后在碎冰中冷却。使用Illumina MiSeq和V3 600循环试剂进行测序。789个样本随机分成三个测序文库，每个文库在独立的流动池上测序，每个流动池产生超过1700万条原始读数。

2.5 分子生物信息学

三个测序文库的Illumina解 multiplexed 序列使用HONEYPI（一种为处理英国国家蜂蜜监测计划的ITS2序列而开发的生物信息学管道）独立处理。生成扩增子序列变异（ASV）表和序列分类学。简要来说，原始序列使用HONEYPI中的DADA2进行质量过滤，参数如下：在质量分数下降低于Q=30时截断读取，并移除短于150 bp的读取。平均每个样本有约17,605条读取通过初始质量（Q≥30）过滤步骤，标准差为13,566条读取 per sample。使用HONEYPI中实现的默认DADA2参数进行去噪和ASV推断。随后，使用在更新的NCBI核苷酸数据库上定制训练的自然贝叶斯分类器进行分类分配，分配基于97%的序列相似性。输出文件包括分类学和ASV表。通过HONEYPI后，使用R按ASV序列合并三个文库的分类学和ASV表，生成所有样本的复合表。读取数少于1000的样本被视为超出正态分布范围，因此被移除。读取数大于1000的样本（n=798）被保留，所有样本均稀释至此读取深度以进行后续分析。

2.6 分析

进行测试以比较（i）属的多样性度量，（ii）记录的属的丰富度，（iii）AVC指示属，和（iv）每种调查方法对AVC预测的准确性。

检查的地点总数为789个，其中126个为AVC-1“作物和杂草”，35个为AVC-2“高草和草本”，173个为AVC-3“肥沃草地”，190个为AVC-4“贫瘠草地”，19个为AVC-5“低地林地”，28个为AVC-6“高地林地”，69个为AVC-7“沼泽草地镶嵌”，149个为AVC-8“荒原和沼泽”。

从ASV和调查表中移除了所有非开花植物数据；因此，苔藓植物、藻类、裸子植物、裸地、落叶和岩石被排除在外。为了尽量减少虚假物种级别分配的风险，我们在属级别进行了所有下游分析。例如，源自油菜（Brassica napus）的序列被错误地分配为其一个或两个亲本谱系（Brassica rapa和Brassica oleracea），因此需要将所有Brassica检测结果合并到属级别。这种保守的方法确保了一致性，并避免了因不确定的物种鉴定而引入人为误差。

然后使用decostand计算每个样本分子数据的比例丰度，并在后续比较和分析之前从每个数据集中移除每个样本的稀有属（<5%丰度）。

为了评估分子数据和调查数据在分类学观察上的相似性，使用每个地点的AVC分类来确定指示属、香农多样性指数以及每种调查方法记录的属的丰富度。通过Shapiro检验确定来自任一种调查方法的数据均非正态分布；因此，使用非参数Spearman相关来评估调查方法之间的关系。计算Spearman's Rho统计量以估计调查方法之间基于等级的关联度量，并使用Benjamini-Hochberg方法调整p值。具体来说，计算了每个AVC内每个地点共同记录的属的丰度之间的相关性，以及每个地点和每个AVC内属级别多样性度量之间的相关性。

我们使用机器学习比较了eDNA丰度和1?m²植物覆盖调查在分配样本地点AVC分类方面的预测能力。为此，我们使用了R包xgboost，采用4:1的训练测试分割，设置如下：method="xgbTree", tuneGrid=expand.grid(nrounds=c(50,100), max_depth=c(2,4,6), eta=c(0.1,0.3), gamma=c(0,1), colsample_bytree=c(0.7), min_child_weight=c(1), subsample=c(0.8))，并进行5折交叉验证检查，以评估每个数据集在预测敏感性（我们能够正确预测的实际阳性案例的数量）、特异性（正确预测的案例中实际为阳性的数量）和准确性（分类器正确预测的频率）方面的表现。

3 结果

3.1 测量的属级别多样性和丰富度以及方法间的相关性

使用每个地点开花植物的群落水平相对丰度数据（合并至属级别）计算每次调查的香农多样性指数。按AVC计算得分平均值，分子和1?m²调查的结果分别为：作物和杂草：0.57和0.16，高草和草本：0.48和0.69，肥沃草地：0.65和0.8，贫瘠草地：0.85和1.32，低地林地：0.56和0.78，高地林地：0.35和0.73，沼泽草地镶嵌：0.64和1.22，荒原和沼泽：0.38和1.93。在所有样本中，分子和1?m²调查的平均香农多样性指数分别为0.61和0.89。

除了AVC-1“作物和杂草”外，分子调查记录的平均香农多样性指数较低。调查方法间香农得分的Spearman等级相关系数分别为：作物和杂草：0.00，高草和草本：0.23，肥沃草地：0.31，贫瘠草地：0.12，低地林地：0.13，高地林地：0.39，沼泽草地镶嵌：0.28，荒原和沼泽：0.38。其中两个AVC（肥沃草地和荒原和沼泽）显示出两种调查方法之间存在显著关联，但Rho系数未超过0.38，表明关联为“弱”或“中等”。基于群落水平相对丰度构建的多样性度量在两种调查方法之间似乎相关性不佳。

类似地，评估了每个地点分子和1?m²调查检测到的属的丰富度。除“作物和杂草”生境外，分子调查的平均丰富度得分较低。计算了每种方法在每个AVC内记录的总属数，并统计了两种调查类型共同记录的属数。各AVC的独特属数值，分子和1?m²调查分别为：作物和杂草：69和38，高草和草本：42和42，肥沃草地：70和38，贫瘠草地：88和68，低地林地：28和28，高地林地：29和39，沼泽草地镶嵌：41和42，荒原和沼泽：23和29。在AVC“作物和杂草”中，分子方法记录了更高的丰富度；这可能是由于先前作物的残留和短暂杂草的检测所致。分子方法在“肥沃”和“贫瘠草地”中记录的总属数也高于传统方法，因为在调查茂密草甸时可能记录到较低的丰富度。在所有样本和所有AVC中，分子和1?m²调查的总属丰富度分别为158和151，其中110个属被共同记录。

3.2 调查方法间共同记录的属丰度相关性

为了评估从分子丰度数据预测开花植物覆盖量的能力，与1?m²植物调查数据进行了比较。在此评估中，仅使用每种AVC中两种调查类型共同记录的属，通过Spearman等级相关分析每个属的总百分比丰度。丰度相关性系数分别为：作物和杂草：0.51，高草和草本：0.37，肥沃草地：0.66，贫瘠草地：0.75，低地林地：0.39，高地林地：0.75，沼泽草地镶嵌：0.58，荒原和沼泽：0.19。五个AVC（作物和杂草、肥沃草地、贫瘠草地、高地林地、沼泽草地镶嵌）显示出调查方法之间的显著关联，Rho系数>0.50。丰度关联为“中等”至“良好”，表明在 taxa 被共同记录的情况下，丰度通常具有良好的相关性。

3.3 聚合植被分类的属级别指示剂

我们评估了调查方法间指示属的共性，以确定分子数据是否可以提供相关的AVC栖息地分类判别器。为每个AVC从每种调查类型中确定了具有显著性的指示属。列出了按指示值排序的每个AVC的显著指示剂。粗体显示的属名是两种调查方法均记录为显著指示剂的属。与先前的多样性和丰富度度量一致，分子数据比1?m²植物调查数据记录了更多的“作物和杂草”指示剂，其中农业或受干扰土地的常见杂草也被记录为指示剂。分子数据还记录了相同的“肥沃草地”指示剂分类单元和更多的“贫瘠草地”指示剂分类单元。在其余的AVC类别中，1?m²植物调查数据产生了更多的植被类型指示剂。

3.4 按调查方法记录的属清单

属级别丰度数据摘要显示了每种调查方法记录的开花植物属的完整清单，以及每个AVC的摘要统计信息，包括样本点数量、记录的总属丰富度及其同记录数量、丰度的Spearman相关得分以及平均AVC香农多样性指数。每种调查方法的丰度得分提供了热图的基础。

3.5 通过机器学习进行聚合植被分类预测

使用属级别丰度数据来训练和验证XGBoost模型，以从eDNA或1?m²植物覆盖调查衍生数据中预测每个样本的AVC分类。在所有AVC中，分子eDNA输出的交叉验证混淆矩阵准确度为0.61，而1?m²植物覆盖调查的准确度为0.74。通过敏感性、特异性和准确性测量的预测能力在调查类型之间大致相当，分子调查数据的预测能力较低，尤其是在低地林地和