《Molecular Ecology Resources》：Estimating Organism Abundance Using Within-Sample Haplotype Frequencies of eDNA Data
编辑推荐：
本篇综述提出了一种创新方法，通过分析环境DNA（eDNA）宏条形码数据中的样本内单倍型频率与群体水平单倍型频率之间的偏差，来估算贡献给单个eDNA样本的个体数量（即“贡献者”数）。该方法首次建立了直接从eDNA数据近似估计群体单倍型频率的理论框架，并引入了最大似然估计量进行贡献者数量推断。研究表明，在单倍型变异足够丰富、群体频率特征明显且样本量充足的情况下，能够获得准确的估算结果。这为从eDNA宏条形码数据中定量评估物种丰度提供了新的路径，弥补了现有分子方法的不足。
引言：环境DNA（eDNA）研究的定量挑战
环境DNA研究主要聚焦于单物种PCR（用于特定类群的精确检测与定量）和宏条形码（用于获取群落组成的广泛快照，但通常仅能提供比例数据）。尽管这两种方法极大地增进了我们对生物多样性的理解，但从eDNA数据中获取的定量信息有限，阻碍了我们在更大尺度上准确监测物种丰度和动态的能力。宏条形码中的读取计数虽然偶尔被定量解读，并且测序读取数与物种丰度之间可能存在表面相关性，但即使在相同的测序深度下，分配给特定生物的读取数也受到多种观测误差源的影响，例如扩增偏倚（由于引物序列不匹配导致）以及数据的组成性本质。即使DNA数量被精确估计，连接总eDNA与生物丰度之间的生物学变异性仍然存在，这涉及生物量、代谢、系统特异性生物和非生物因素以及eDNA衰变和运输的时空动态。
然而，宏条形码数据提供了超越物种组成和读取计数的更多信息。它们还包含种内遗传变异，因为单个环境样本可以包含来自同一物种不同个体的多个单倍型（或扩增子序列变异，ASV）。在随机交配的无结构群体（即泛交群体）中，eDNA样本中检测到的单倍型代表了从群体单倍型频率分布中的多项分布抽取。随着贡献给样本的个体数量（以下简称“贡献者”）增加，检测到的等位基因可能更多，并且样本内观察到的等位基因频率预计会收敛于这些等位基因的群体频率。因此，我们假设eDNA样本中测量的单倍型频率与群体预期的频率之间的差异反映了贡献者的数量，差异越大表明贡献者越少，差异越小表明贡献者越多。
方法：理论基础与估算框架
该方法利用eDNA观测中单倍型频率相对于群体频率的变异性来提取关于每个样本贡献者数量（N）的信息，分两步进行：首先，使用完整的、空间重复的eDNA数据集推导群体水平单倍型频率向量π；其次，提供一种正态近似最大似然估计量（MLE）来获得N的点估计和置信区间。
设π = (p1, …, pK) 表示同一基因座上K个单倍型的群体频率。假设收集了独立的eDNA样本j，用Nj表示样本j的贡献者总数。在泛交条件下，贡献给每个样本的单倍型可被视为多项分布抽取。观察到的频率fij的方差反映了从个体DNA释放后所有下游过程的累积效应。将所有增加的方差合并为一个单一的乘法膨胀因子ξ，作用于pi(1 - pi)/Nj项，使得随着贡献者数量减少，等位基因频率中观察到的方差成比例增加。
关键假设与要求包括：需要已知的等位基因频率；观测噪声不过度 overwhelming；群体处于泛交状态（即研究区域内无群体结构），且所用基因座处于哈迪-温伯格平衡；假设eDNA数据中观察到的ASV对应于群体中的真实单倍型，而非错误等位基因。
从eDNA样本推导群体单倍型频率
在所述多项分布模型下，由于每个观察到的fij都是以pi为中心的无偏抽取，因此跨多个独立样本取平均值提供了真实群体频率的自然估计量。为了避免循环性，当估计给定样本j的贡献者数量时，我们采用留一法（LOO）计算π(-j)，即从除j之外的所有样本推导出的群体频率估计。
估计单个样本的贡献者数量
一旦群体单倍型频率已知，就可以根据观察到的频率的变异性推导个体贡献者的数量。一种统计上严格的方法将是使用层次模型，明确考虑图2中突出的每个采样阶段增加的方差。然而，在实践中，大多数数据集很少有所需的高水平重复，这限制了模型的可识别性和适用性。因此，我们通过忽略这些中间未观察阶段及其相关方差来简化模型，直接连接观察到的单倍型频率与现实世界中贡献者数量的生物学方差。
正态近似对数似然法
根据中心极限定理（CLT），每个观察到的频率fij可以近似为正态分布。在该假设下，对数似然函数可以导出，并且通过微分可以得到贡献者数量Nj的封闭形式解（点估计）。为了评估N的不确定性，我们采用轮廓似然法。根据标准正则条件，似然比统计量渐近服从自由度为1的卡方分布。置信区间定义为使得似然比统计量小于或等于卡方分布特定分位数的N值集合。
基于倍性的最小贡献者数量
除了似然法，我们还根据所分析基因座的倍性施加样本的最小贡献者数量。对于单倍体基因座，要求eDNA样本的估计贡献者数量至少等于观察到的单倍型数量。对于二倍体基因座，最小贡献者数量定义为观察到的单倍型数量的一半（向上取整）。在无法推断倍性的情况下（例如受异质性影响的基因座），可以禁用此约束。
模拟验证
为了验证统计方法，我们开发了一系列R函数来模拟eDNA数据，并将其作为概念验证应用。这些函数允许我们生成复制eDNA采样关键过程的合成数据集，包括个体贡献、衰变、等位基因频率和采样偏倚。我们进行了两种类型的模拟：一种在理想化的“无误差”条件下（ξ = 1），另一种使用更具生物学真实性的参数（ξ > 1）。
结果
eDNA数据模拟
eDNA模拟导致样本内单倍型数量可变，随着贡献者数量增加，观察到的单倍型也更多。此外，当我们将观察到的单倍型频率fij的方差作为贡献者数量的函数进行观察时，即使在个体贡献异质的情况下（ξ > 1），图2中描述的模式也会出现。观察到的频率与所有单倍型的群体频率之间的这种差异本质上是所提出的定量方法用来推导最可能贡献者数量的依据。
从eDNA样本快速估计群体频率
我们评估了仅使用LOO平均法从模拟eDNA样本中恢复群体水平单倍型频率（π）的准确性。模拟在两种场景下进行：一种是假设的超可变标记（包含14个单倍型），另一种是保守标记（一个单倍型占主导，共6个单倍型）。正如预期，跨空间重复的eDNA样本聚合单倍型频率减少了随机方差，导致向真实群体频率收敛。增加模拟样本数量提高了准确性，且具有更多个体贡献者的样本收敛更快。对于两种标记类型，频率估计都迅速稳定：无论贡献者数量或单倍型频率如何，20个样本足以在平均平方误差（MSE）阈值内恢复群体频率。当eDNA样本平均包含超过10个贡献者时，使用少至100个有目标检测的样本即可将群体频率估计误差控制在1%以内。在所有情况下，eDNA推导的频率比组织样本更快地近似真实群体频率。
模拟数据中准确预测贡献者数量
一旦群体单倍型频率已知，就可以应用正态近似最大似然法来估计每个样本的贡献者数量。为了评估我们方法的最佳性能，我们首先在理想条件下模拟eDNA数据（所有个体向当地eDNA库贡献等量DNA，且未引入观测误差，ξ = 1）。使用具有17个等位基因的假设超可变标记，该方法提供了与真实贡献者数量密切匹配的准确估计。与真实值有很强的相关性，并且90%的估计在其95%置信区间内捕获了真实值。然后，我们使用具有保守单倍型特征的标记（一个主要等位基因，仅6个总单倍型）重复模拟。在这些条件下，该方法产生了严重偏差且准确性较低的估计，对于较低贡献者数量往往被高估。与真实贡献者数量的相关性下降，平均相对误差（MRE）增加，覆盖率下降。这些结果强调，即使在理想条件下，贡献者估计的准确性也 fundamentally 与基因座上的单倍型多样性相关。随着单倍型分布变得越来越偏斜，特别是当一个等位基因占主导时，估计误差显著增加。
额外方差对估计的影响
然后，我们使用两组等位基因频率重复模拟，这次引入了个体eDNA贡献的高度可变性和一些观测偏倚，以评估这些方法在更现实的误差源下的表现。这导致贡献者单倍型频率与观察到的eDNA单倍型频率之间存在显著差异。我们随后重复估计量，使用一个或三个生物学重复（从同一当地eDNA库模拟的独立瓶装样品）。正如预期，在这种场景下，该方法低估了贡献者数量，尤其是在贡献者最多的样本中。这是由于观察到的单倍型频率中的误差被误解为较少的贡献者。尽管如此，在可变标记下估计仍然合理。当分析单个生物学重复时，与真实贡献者数量的相关性下降，MRE增加，但覆盖率下降更为显著。当存在观测噪声时，理论上生物学重复可以帮助缓解其影响，因为下游噪声（除eDNA生态学外的所有部分）可以平均掉。因此，当考虑三个生物学重复时，MRE和相关系数仅略有恢复。这是因为模拟中贡献者单倍型频率与观察到的频率之间的大部分差异是由于eDNA生态学（每个个体的可变贡献）造成的，而不是重复共享的实验室和观测过程。当使用噪声和保守标记的组合时，估计要差得多。与真实贡献者数量的相关性下降到很低，覆盖率也显著下降。与保守标记类似，添加生物学重复对改善估计几乎没有影响。
讨论
本研究引入了一个新框架，利用样本内单倍型频率和已知或推断的群体水平单倍型分布来估计环境DNA样本的个体贡献者数量。通过利用观察到的和群体单倍型频率之间的自然统计收敛性，我们表明可以从eDNA宏条形码数据中推导出有意义的丰度估计，而无需直接依赖eDNA分子的数量。模拟证实，该方法在理想化条件下表现良好，并且在纳入生物学上现实的噪声源（如贡献者大小、距离和释放率的变化）时，只要标记足够可变，仍然能提供信息。总之，这些结果表明，单倍型频率模式为eDNA数据提供了一个很大程度上尚未开发的定量信息来源，弥合了宏条形码的传统存在/缺失解释与更具物种特异性的定量方法之间的差距。
一个关键的方法学进步是使用样本内单倍型频率来估计贡献者丰度。早期的模型完全依赖于单倍型的身份或数量，而忽略了观察到的频率。因此，我们的方法充分利用了数据中存在的信息，从而得到更准确和细致的估计。值得注意的是，我们没有显式地对单倍型存在或组合排列进行建模；相反，这些效应隐式地捕获在频率模式中。这需要假设观察到的频率以贡献者的真实频率为中心。然而，我们的模拟（包含了大量且现实水平的噪声）表明，即使这个假设被适度违反，该方法也是稳健的。
我们还表明，该方法广泛适用于不同单倍型多样性的基因座，但在应用于高变标记时表现明显更好。这种模式反映了先前提出的仅基于单倍型计数或身份的方法。对标记变异性的依赖给将该方法应用于许多现有的eDNA数据集带来了挑战，因为常用标记（如12S和COI）相对保守。这种权衡反映了eDNA宏条形码内部的权衡：为了检测广泛类群，引物必须靶向保守区域，这限制了分类分辨率和定量信号的丰富性。然而，随着测序技术的进步和从eDNA中可靠地获得更长的读取，将可以恢复更多的可变位点，从而增加基于单倍型推断的能力。
除了估计单个eDNA样本的贡献者数量外，我们还提供了一个直接从eDNA数据估计群体单倍型频率的框架。虽然可以使用现有的组织衍生参考数据集，但从eDNA估计π提供了几个方法学优势。首先，它确保了内部一致性：用于估计贡献者的相同eDNA数据也定义了参考频率。这种对齐自然捕获了任何数据集特定的扩增或测序偏倚——这些偏倚是组织衍生的单倍型频率可能遗漏的，特别是在涉及引物错配或多个测序运行的情况下。其次，eDNA经常恢复组织数据集中罕见或缺失的单倍型。当此类单倍型出现在后续样本中时，除非群体频率是从eDNA数据本身推断出来的，否则它们缺乏参考频率将无法使用。虽然没有单个样本能够捕获完整的单倍型多样性，但跨空间和时间重复样本的汇总频率可以作为一个生物学和方法学上一致的代理。这种方法类似于传统的群体调查，其中可靠性来自于重复，而不是任何单一的观察。
实际考虑与方法局限性
在实际应用中，我们建议在以下条件下应用此框架：首先，群体单倍型频率应通过跨空间或时间尺度的重复，通过多次eDNA检测（> 20次）合理地近似。其次，测序深度和模板DNA应足够，以使单倍型频率保留有关贡献者的信息，而不是被噪声主导。第三，应选择具有足够种内变异性的标记，对于提供有限分辨率的保守基因座（如12S）应谨慎。这些技术考虑与关于群体结构和独立性的更深层生物学假设同时运作。
重要的是，我们的方法从根本上依赖于这样的假设：eDNA贡献者代表来自非结构化群体的独立随机抽取（泛交）。然而，当自然群体具有广泛的地理范围并表现出群体结构时，这一假设很容易被违反。例如，在具有显著地理结构的群体中，本地样本可能反映与整体群体频率系统性偏差的单倍型频率，如果与全局频率比较，将导致有偏的贡献者估计。因此，如果可用，本地群体频率更可取。同样常见的是，入侵或扩张物种会在其分布边缘显示出空间排序，这将打破泛交假设。此外，在具有强烈社会联系或家庭群体的物种中，eDNA样本可能会捕获来自相关个体的遗传物质，这些个体的单倍型不是独立的。这在分析线粒体标记时尤其成问题，因为线粒体标记遵循严格的母系遗传模式。在这种情况下，我们的模型将无法区分具有相同遗传组成的不同场景。这些生物学现实强调了仔细选择标记和深思熟虑解释结果的重要性。尽管如此，需要指出的是，这种对泛交的依赖并非我们框架独有，而是反映了许多经典群体遗传模型（如F统计和Kingman溯祖模型）的简化基线。这样的假设在自然界中很少完美成立，但这些模型仍然有用，因为它们提供了一个易于处理的参考点：偏差会引入偏倚但不会消除所有的推断价值。该框架的未来扩展同样可以放宽这一假设以适应结构化群体。
正如引言中所讨论的，估计贡献者数量的方法的常见限制，包括本文提出的方法以及Andres等人开发的方法，是它们对群体水平单倍型频率的依赖。在本研究中，我们建议使用跨多个eDNA样本的单倍型分布来近似这些频率，将汇总频率作为群体水平值的代理。虽然我们表明这种方法能快速近似真实频率，但这必须谨慎进行。首先，必须意识到可能包含由PCR或测序错误引起的虚假或人工单倍型。为了在将本文方法应用于真实数据集时降低这种风险，重要的是使用过滤方法以确保观察到的序列变异确实与单倍型相关。即使采用严格的生物信息学流程，eDNA数据中也可能存在错误等位基因。幸运的是，像lulu这样的工具可以通过基于相似性和共现对序列进行聚类来折叠错误变异，特别是在使用模拟群落校准时。最近，专用于eDNA数据的R包，如TombRaider和gmmDenoise，也引入了统计过滤功能，可以移除虚假ASV。总之，这些方法使得从环境样本中移除人工产物并恢复可靠的单倍型变得切实可行。
迄今为止，大多数关于eDNA群体推断的实证研究都是在水生环境中进行的。我们提出的框架是否适用于陆地或空中eDNA仍不确定，因为这些系统中DNA的生态学特征远未得到充分了解。特别是空气中的样本通常产生异质性且相对较低的DNA浓度，与水相比，这些条件可能会限制跨多个贡献者恢复单倍型频率信息。在这种情况下，像本文提出的这种基于频率的方法可能信息量较少，而依赖单倍型身份或计数的方法可能更合适。阐明这些替代框架在不同样本类型中的适用性将需要有针对性的实证研究。
最后，在解释结果时，必须强调推导出的贡献者数量不应被解释为一个地区个体生物的确切计数。相反，它们代表了其DNA通过采样和测序过程被捕获和保留的贡献者的保守估计。这个值仅间接地，并且可能通过一个复杂的非线性函数，与真实的生物丰度相关。例如，如果两个湖泊被采样进行eDNA分析，一个湖泊估计有10个贡献者，另一个有100个，可以合理地推断后者总体上拥有更多个体，但不一定是正好100个个体，甚至是前者的10倍。该方法最好被解释为提供了贡献者数量的保守下限，排除了关键的现实世界方差来源，如差异释放率、环境降解和测序偏倚。
未来步骤与验证
我们在此仅将该方法应用于模拟的eDNA样本，但需要并鼓励在广泛方法应用之前，使用野外和真实的eDNA样本进行进一步的方法验证，因为真实eDNA样本中添加的“噪声”水平可能过高，以至于该方法在实践中不适用。验证本文方法的理想数据集应包括：已知单倍型频率的群体、已知或独立估计的贡献个体数量、空间分布的采样以及使用高变标记进行的eDNA测序。据我们所知，没有公开可用的数据集满足所有这些标准，生成这样一个数据集作为应用该方法的 proof of concept 超出了本研究的范围。大多数现有的eDNA宏条形码工作依赖于相对保守的线粒体基因座（如12S、COI或CytB），在这些基因座上，所提出的方法可以工作但准确性非常低。因此，虽然我们选择使用模拟数据来呈现和应用该方法以展示概念，但未来的工作应使用经验数据来验证模型的有效性。
该方法的未来改进或迭代还应探索更可靠地估计ξ的方法，以便校正噪声场景下的低估偏倚，并使该方法对于大量贡献者更可靠。此外，未来的研究可以通过将基于单倍型的推断与其他独立数据流相结合来改进丰度估计。一个很有前景的途径是将这种方法与总eDNA定量（例如来自qPCR或dPCR）相结合，后者反映了总生物量但受时空偏倚的影响。因为这两种方法捕获了不同的生物学信号（个体数量 vs. 生物量），它们可以提供互补且独立的丰度估计。原则上，一个具有高eDNA数量但低单倍型多样性的样本表明来自少数个体的许多分子，而两者都具有高值则表明来自许多个体的贡献。此外，跨多个不连锁的基因座运行该方法可以进一步提高稳健性，特别是在二倍体或重组的生物中。这可以通过跨基因座估计值相乘似然分布来实现。虽然我们在此没有探索这种集成框架，但它们代表了从eDNA更准确且更有生物学意义地估计群体丰度的明确下一步。
