文章阅读目录大纲
宏基因组测序所处理的对象是直接对环境样本中的所有DNA进行测序。达到无需培养即可揭示微生物群落的组成和功能潜力的目的。在数据处理中，一个核心任务是从海量短读序列中估算物种丰度（即每个物种在样本中的相对含量）和基因丰度（即每个基因或功能单元的相对含量）。传统的基于序列比对的方法计算成本高昂，而基于k-mer的方法通过利用固定长度的子序列（k-mer）信息，能够在不依赖完整比对的情况下快速估算丰度。
k-mer是指长度为k的连续子序列，例如在k=2的时候，DNA序列“ATCG”包含的2-mers有“AT”、“TC”、“CG”。通过统计读序列中k-mer的出现频率，并将其与参考数据库中的k-mer频率进行比较，我们可以推断出样本中各物种或基因的丰度。这种方法具有计算速度快、内存效率高的优势，并且无需对每个读进行精确比对，因此在处理大规模宏基因组数据时非常实用。
k-mer方法在宏基因组学中的应用
k-mer分析是现代基因组学和宏基因组学中的一项基础技术，其应用范围广泛，包括序列组装、错误校正、物种分类和丰度估算等。在宏基因组学中，k-mer方法主要用于快速分类和定量微生物群落，以及功能注释和丰度估算。与传统的基于比对的方法相比，k-mer方法具有以下优势：
- 速度与效率：k-mer统计可以通过哈希表或布隆过滤器等数据结构高效实现，避免了耗时的序列比对。例如，Kraken等工具通过将读序列的k-mer与预构建的k-mer数据库进行快速匹配，实现对数百万条读的快速分类。这种查找表的方式比动态规划比对快几个数量级。
- 内存优化：通过概率数据结构（如布隆过滤器）和紧凑的编码（如2-bit编码核苷酸），k-mer计数可以在有限内存下处理海量数据。例如，khmer工具利用Count-Min Sketch数据结构，在内存中仅存储k-mer的计数而非k-mer本身，从而大幅降低内存需求。
- 无需参考基因组：某些k-mer方法（如基于k-mer频率的组成分类）不依赖完整的参考基因组，仅利用k-mer组成模式即可对读进行分类或聚类，适用于未知微生物的发现。
- 内在生物学意义：k-mer的分布和频率本身携带了序列的生物学信息。例如，某些k-mer可能在不同物种中具有特异性频率，可用于物种鉴定；又如，k-mer的缺失（nullomers）或重现（neomers）在疾病检测和疫苗设计中具有潜在应用。
在宏基因组学中，k-mer方法通常用于物种丰度估算和基因丰度估算这两个层面。前者关注哪些物种存在及其相对数量，后者关注哪些基因或功能存在及其相对丰度。虽然两者都基于k-mer统计，但在数据库选择和算法策略上存在差异。下面我们将分别讨论这两类应用，然后重点介绍基因丰度估算所需的方法学。
物种丰度估算与基因丰度估算的差异
物种丰度估算和基因丰度估算是宏基因组分析的两种常见需求，它们在目标和应用上有所不同，因此在k-mer数据库的构建和使用上也存在差异：
- 目标不同：物种丰度估算旨在确定样本中各微生物物种（或菌株）的相对含量，通常以分类学单位（如种、属）为单位报告结果。基因丰度估算则关注样本中各基因或功能单元（如蛋白家族、代谢通路）的相对丰度，以功能单位为单位报告结果。
- 参考数据库不同：物种丰度估算通常使用基因组或标记基因数据库。例如，Kraken使用完整的微生物基因组构建k-mer数据库，每个k-mer标记其所属物种；MetaPhlAn等工具则使用物种特异性标记基因（如单拷贝的管家基因）数据库，通过比对读的k-mer到标记基因来推断物种存在和丰度。相比之下，基因丰度估算使用基因或蛋白家族数据库。例如，HUMAnN等工具将读比对到整合的泛基因组（pangenome）或功能基因数据库（如UniRef、KEGG Orthology）来量化基因家族和通路的丰度。因此，物种估算数据库侧重于分类标识，而基因估算数据库侧重于功能标识。
- k-mer数据库规模和复杂性：物种估算数据库往往包含成千上万个微生物基因组，k-mer空间巨大，但每个k-mer通常只对应一个物种（理想情况下），因此可以通过精确匹配或最小公共祖先（LCA）策略进行分类。基因估算数据库则可能包含数百万个基因序列，许多k-mer在不同基因中重复出现，因此需要处理多映射（ambiguity）问题。例如，同一个k-mer可能存在于多个同源基因中，需要算法来分配其贡献。
- 丰度计算策略：物种丰度估算通常通过读分配或k-mer计数直接得到。例如，Kraken统计每个物种的k-mer命中数，Bracken则利用贝叶斯模型校正这些计数以获得更准确的物种丰度。基因丰度估算则更复杂，因为一个基因可能跨越多个k-mer，且读可能只覆盖基因的一部分。因此，基因丰度往往通过映射读的覆盖度或k-mer频率的加权来估算。例如，HUMAnN首先将读比对到物种的泛基因组，然后对未比中的读进行翻译搜索，最后根据比对结果量化基因家族和通路的丰度。
简而言之，物种丰度估算关注“谁在那里”，其k-mer数据库侧重于分类标记；基因丰度估算关注“它们在做什么”，其k-mer数据库侧重于功能基因。开发一个基因丰度估算工具，需要针对功能数据库的特点设计算法，例如处理同源基因的k-mer共享、考虑基因长度差异、以及将k-mer计数转换为基因相对丰度等。下面我们将重点讨论基因丰度估算工具所需的k-mer数据库和算法原理。
k-mer数据库的选择与构建
要开发一个基因丰度估算工具，首先需要构建或选择一个合适的k-mer参考数据库。该数据库应包含目标功能基因的k-mer信息，用于与样本读序列的k-mer进行比较。以下是构建此类数据库时需要考虑的关键因素：
- 数据库内容：根据分析目标选择基因集合。例如，如果关注代谢功能，可以使用KEGG Orthology (KO)或MetaCyc数据库中的基因；如果关注蛋白家族，可以使用Pfam或eggNOG等。每个基因序列将作为k-mer提取的来源。
- 基因序列来源：可以来自公共基因数据库（最常见的核酸数据来源就是NCBI上的nt库。顺便提一句，假若做宏蛋白组分析，常用的数据库则是UniProt数据库）或整合的泛基因组。泛基因组包含一个物种所有已知基因，有助于捕获物种内部的基因多样性。HUMAnN2就采用了物种泛基因组数据库，将读比对到已知物种的泛基因组以获得更全面的功能谱。
- k-mer长度选择：k值的选择需要权衡特异性和敏感性。较长的k-mer（如k=31）具有更高的特异性，能够区分相似基因，但可能因测序错误而丢失匹配；较短的k-mer（如k=15）更鲁棒于错误，但特异性降低，容易产生假阳性匹配。通常，对于Illumina短读（如150 bp），常用k=31或k=35，以在特异性和覆盖度之间取得平衡。
- 数据库构建方法：对于每个基因序列，提取其所有可能的k-mer（通常考虑正反链），并记录这些k-mer与基因的对应关系。构建数据库时，可以采用哈希表（将k-mer映射到基因列表）或布隆过滤器（快速判断k-mer是否存在于数据库中）。为了节省内存，可以使用压缩或计数 sketches（如Count-Min Sketch）来存储k-mer频率信息。此外，需要处理冗余k-mer：如果一个k-mer出现在多个基因中，数据库应记录其关联的所有基因，以便后续分配权重。
- 更新与定制：理想情况下，数据库应允许用户根据需要添加自定义基因序列。例如，用户可能希望将新测序的微生物基因纳入分析。因此，工具应提供数据库更新接口，允许用户从FASTA文件构建新的k-mer数据库。
构建完成后，该k-mer数据库将作为基因丰度估算的参考。在分析时，工具会将样本读序列的k-mer与数据库进行比对或查询，以推断哪些基因存在及其丰度。那，在接下来这里，我们将介绍如何利用k-mer数据库进行基因丰度估算的算法原理。
基于k-mer的基因丰度估算算法原理
基于k-mer的基因丰度估算通常包括以下步骤：k-mer提取与计数、k-mer与数据库匹配、基因丰度推断和归一化。下面将详细讨论每个步骤涉及的算法原理和数学模型。
1. k-mer提取与计数
首先，需要对宏基因组样本的读序列进行k-mer提取和计数。这一步的目的是获得样本中每个k-mer的出现频率。算法流程如下：
- 读序列预处理：去除低质量序列和接头污染（可选），以减少错误k-mer的干扰。
- 滑动窗口提取k-mer：对每条读序列，以滑动窗口方式提取所有长度为k的子序列。例如，读序列“ATCGAT”在k=3时包含“ATC”、“TCG”、“CGA”、“GAT”四个3-mers。
- 正反链处理：由于DNA是双链的，通常考虑正反链的k-mer。可以采用最小字典序策略：对每个k-mer，取其本身和反向互补中字典序较小者作为代表，以合并正反链的计数。
- 计数数据结构：使用高效的数据结构统计k-mer频率。传统方法是使用哈希表（键为k-mer，值为出现次数），但内存消耗高。为节省内存，可采用布隆过滤器（Bloom Filter）或Count-Min Sketch等概率数据结构。布隆过滤器可以快速判断一个k-mer是否在数据库中，但无法精确计数；Count-Min Sketch则可以给出k-mer计数的近似值，且内存占用远低于哈希表。例如，khmer工具使用Count-Min Sketch实现了对海量k-mer的在线计数。
- 错误k-mer过滤：由于测序错误，样本中会存在大量仅出现一次的k-mer（singleton k-mers），这些通常不提供有用信息且占用内存。可以在计数阶段设置阈值，忽略或过滤掉出现次数过低的k-mer，从而降低内存使用并提高后续分析的准确性。
经过上述步骤，我们得到样本中每个k-mer的计数（或存在性信息）。这些k-mer计数将用于与参考数据库进行匹配，以推断基因丰度。
2. k-mer与数据库匹配
接下来，将样本k-mer与参考基因k-mer数据库进行匹配，以确定样本中可能包含哪些基因。匹配策略可以分为精确匹配和模糊匹配两类：
- 精确匹配：对于每个样本k-mer，查找其在数据库中是否存在。如果存在，则将其与数据库中记录的基因关联起来。精确匹配可以通过哈希表快速实现，但要求k值足够长以避免随机匹配。例如，Kraken使用k=35的精确k-mer匹配来分类物种。在基因丰度估算中，精确匹配可以确保我们只考虑高度可信的k-mer，但可能因测序错误而漏掉一些真实存在的基因。
- 模糊匹配：允许一定程度的错配或使用更短的k-mer，以增加对错误的鲁棒性。例如，可以使用 spaced seeds（间隔种子）或最小哈希（MinHash）技术来捕捉序列相似性。MinHash通过比较k-mer集合的Jaccard相似度，可以快速估计序列相似度，被用于Mash等工具进行基因组距离估算。在基因丰度估算中，模糊匹配可以检测与参考基因高度相似但不完全相同的序列，但需要权衡特异性。
- 多映射处理：一个样本k-mer可能匹配数据库中的多个基因（同源基因共享k-mer）。此时需要决定如何分配该k-mer的计数贡献。简单策略是平均分配（将计数均分给所有匹配基因），或加权分配（根据基因长度、k-mer特异性等因素加权）。更高级的方法使用统计模型来推断最可能的来源基因。例如，HUMAnN2在将读比对到泛基因组后，对未比中的读进行翻译搜索，以捕获与参考基因有差异的同源基因。
- 数据库查询优化：为了加速匹配过程，可以预先对数据库建立索引。例如，使用布隆过滤器快速过滤掉不在数据库中的k-mer。对于大型数据库，可采用分层或分块策略：先将k-mer分配到不同的子数据库（如按物种或功能分类），再分别查询，以减少单次查询的数据量。
经过匹配，我们得到每个基因在样本中被匹配到的k-mer计数（或读数）。这些原始计数需要进一步处理才能转换为基因丰度。
3. 基因丰度推断
从k-mer匹配结果推断基因丰度是算法的核心。这里涉及将k-mer计数转换为基因丰度，并考虑各种偏差和不确定性。主要考虑因素包括：
- 基因长度偏差：较长的基因会包含更多的k-mer，因此在相同表达水平下，其k-mer计数会更高。为了公平比较不同基因的丰度，需要对基因长度进行归一化。常见做法是计算每千碱基每百万读（RPKM）或每百万映射读的转录本每百万（TPM）等指标。TPM是更推荐的方法，它先对基因长度进行归一化，再对测序深度进行归一化，使得各样本的TPM总和相同，便于跨样本比较。
- 测序深度偏差：不同样本的测序数据量不同，需要进行归一化。TPM已经包含对测序深度的归一化，也可以使用每百万读（RPM）等简单归一化方法。
- k-mer多映射偏差：如果采用平均分配策略，多映射k-mer可能稀释高丰度基因的信号。可以采用期望最大化（EM）算法来迭代估计各基因的丰度，使得多映射k-mer的分配更符合实际丰度分布。EM算法假设每个k-mer的观测计数是各基因丰度的混合，通过迭代更新基因丰度估计，使似然最大化。这种方法在RNA-seq定量中也有应用（如RSEM类似思路）。
- 零膨胀和过离散：宏基因组数据中，许多基因在样本中可能完全未被观测到（计数为零），而观测到的计数往往具有过离散（方差大于均值）。直接使用泊松分布假设可能低估方差，导致统计推断偏差。为此，可以采用零膨胀模型或负二项分布模型来拟合计数数据。例如，零膨胀泊松（ZIP）模型将数据视为“零”和“非零”两部分，分别建模，能够更好地处理大量零值。负二项分布则引入离散参数，允许方差大于均值，适用于过离散计数数据。这些模型可以用于差异丰度分析，也可以用于改进丰度估计的准确性。
- 统计推断：在获得基因丰度估计后，如果需要比较不同条件下的基因丰度差异，应使用考虑计数分布特性的统计检验。例如，使用零膨胀负二项回归来识别差异丰度基因，可以提高统计功效并控制假阳性。
一个完整的工具并非孤立使用所有这些模型，而是将它们有机组合。综合以上因素，一个典型的基因丰度推断流程如下：
- 原始k-mer计数汇总：将匹配到每个基因的所有k-mer计数相加，得到每个基因的原始k-mer计数（或读数）。
- 长度和深度归一化：计算每个基因的RPKM或TPM值，以消除基因长度和测序深度的影响。
- 可选：EM迭代：如果存在大量多映射k-mer，可运行EM算法调整基因丰度估计，使多映射k-mer的分配更合理。
- 可选：统计模型校正：使用零膨胀或负二项模型对丰度数据进行建模，以获得更准确的丰度估计和置信区间。
- 输出：得到每个基因的相对丰度（如TPM）以及可能的统计不确定性（如标准误或置信区间）。
通过上述算法，我们可以从k-mer匹配结果推断出样本中各基因的丰度。接下来，我们将讨论如何将这些算法整合为一个完整的数据处理流程。
基因丰度估算算法详细介绍
假若我们需要自己开发这样的一个用于从宏基因组测序数据中估算基因丰度的计算工具，工具的输入理论上通常是宏基因组测序的FASTQ文件（单端或双端读序列）以及一个参考基因k-mer数据库。这个工具应该向用户输出是每个基因的丰度估计，可以以表格形式输出（如CSV或TSV），包含基因标识符、丰度值（TPM）以及可能的统计度量。在基因丰度估算的核心算法模块，会涉及一些关键的数学算法和模型，下面我们来进行详细说明：
1. Count-Min Sketch
Count-Min Sketch方法是一种概率数据结构，用于高效统计k-mer频率。它由多个哈希函数和二维数组组成：对于每个k-mer，使用多个哈希函数映射到数组的不同位置，取这些位置中的最小值作为该k-mer的计数估计。Count-Min Sketch能够以高概率保证计数估计不超过真实计数，且内存占用远低于哈希表。在实现时，需要选择合适的哈希函数数量和数组宽度，以平衡准确性和内存。这种概率数据结构，其核心优势在于用可控的内存开销来近似统计海量元素的频率，并且能够提供一个计数上限的估计。
但是需要我们进行注意的是，这种概率方法其估计值永远是有偏的（例如，总是高估），但能让我们以可控的内存误差换取高效率。我们会需要在使用这个概率模型的时候通过增加哈希函数和数组宽度的调整，以可以在准确性和内存之间进行权衡。
2. TPM计算
与转录组表达数据处理方法一样，在宏基因组数据处理中同样也可以按照相同的计算方法来计算出可比较的基因丰度数据。在这里面，TPM（Transcripts Per Million）是一种常用的基因丰度归一化指标。计算步骤为：
- (1) 对每个基因，计算其k-mer计数（或读数）除以基因长度（以kb为单位），得到RPK值；
- (2) 将所有基因的RPK值求和，得到总RPK；
- (3) 对每个基因，计算其TPM = (RPK / 总RPK) * 1,000,000。
TPM确保每个样本的基因丰度总和相同，便于跨样本比较。实现时需要注意基因长度的计算（通常以有效碱基对数计，排除低复杂度区域）。
3. 期望最大化（EM）算法
我们可以通过期望最大化（EM）算法来处理多映射k-mer的分配问题。假设有 G 个基因， K 个k-mer，nk 为k-mer k 的观测计数，ag为基因 g 的真实丰度（待估计）。每个k-mer k 可能属于多个基因，记其所属基因集合为Sk。EM算法的步骤如下：
E步：对于每个k-mer k，计算其由基因 g 产生的期望计数：
# e_kg = n_k * (a_g * w_kg) / sum(a_h * w_kh for h in S_k)
# 这是一个概念性的 R 表达式，展示计算逻辑
# 假设我们有以下数据结构:
# n: 一个向量，存储每个 k-mer 的观测计数, n[k] 对应 n_k
# a: 一个向量，存储每个基因的丰度估计, a[g] 对应 a_g
# w: 一个矩阵，w[k,g] 对应 w_kg (k-mer k 对基因 g 的权重)
# S_k_list: 一个列表，S_k_list[[k]] 包含 k-mer k 所属的所有基因的索引
# E步: 计算期望计数 e_kg
e <- matrix(0, nrow = length(n), ncol = length(a)) # 初始化 e 矩阵
for (k in 1:length(n)) {
S_k <- S_k_list[[k]] # 获取 k-mer k 的所属基因集合
if (length(S_k) > 0) {
# 计算分母: sum_{h in S_k} a_h * w_kh
denominator <- sum(a[S_k] * w[k, S_k])
# 避免除以零
if (denominator > 0) {
# 计算分子: n_k * a_g * w_kg
numerator <- n[k] * a[S_k] * w[k, S_k]
# 更新 e 矩阵中对应的行
e[k, S_k] <- numerator / denominator
}
}
}
# e[k,g] 现在存储了 e_kg 的值其中wkg表示k-mer k 对基因 g 的权重（例如，若k-mer在基因 g 中出现多次，则权重高）。
M步：更新基因丰度估计：
# a_g = sum(e_kg for all k) / L_g
# L_g: 基因 g 的长度（以k-mer数量计）
# 假设我们有:
# L: 一个向量，存储每个基因的长度, L[g] 对应 L_g
# M步: 更新基因丰度 a
# 对每个基因 g，将 e 矩阵中对应列的所有 e_kg 求和
sum_e_kg <- colSums(e)
# 更新丰度 a
a <- sum_e_kg / L
# a[g] 现在存储了更新后的 a_g 值其中Lg为基因 g 的长度（以k-mer数量计），相当于对每个基因的期望k-mer计数求和并归一化。
迭代E步和M步，直到收敛（例如，丰度估计变化小于阈值）。EM算法能够将多映射k-mer按照当前丰度估计进行分配，从而逐步逼近真实的基因丰度分布。需要注意的是，EM算法对初始值敏感，容易陷入局部最优。实践中常使用更简单的方法（如将唯一映射的k-mer计数分配给基因）的结果作为EM算法的初始丰度估计，以改善收敛结果。
目前，我们主要可以是通过上述的数学算法和模型的实现，来构建出一个可以准确地从k-mer数据推断基因丰度，并提供合理的统计不确定性度量的计算工具。
在这里，EM算法是核心丰度估算算法。EM算法是本篇文章中的方法主角。它的任务是利用k-mer的观测计数，推断出每个基因的丰度（TPM）。它输出的是一个点估计。在基于EM方法完成了基因丰度的估算之后，我们可以得到包含有出多个样本的基因丰度矩阵。假若我们需要进行下游的差异丰度统计分析，例如我们想知道“基因A在疾病组和健康组的丰度是否有显著差异？”。由于基因丰度计数数据具有过离散和大量零值的特点，负二项分布（如DESeq2、edgeR工具）和零膨胀负二项分布（如metagenomeSeq工具）成为了DA分析的标准模型，在这里我们一同介绍下下游相关分析的统计算法原理：
4. 零膨胀模型
这个算法可以用于处理大量零计数的统计模型。零膨胀泊松（ZIP）模型假设观测计数 Y 由两部分混合产生：以概率 π 来自零点分布（恒为0），以概率 1-π 来自泊松分布（均值 λ）。其概率质量函数为：
# P(Y=y) = { pi + (1-pi)*exp(-lambda) if y = 0
# { (1-pi) * (exp(-lambda) * lambda^y) / y! if y > 0
# 定义一个函数来计算 ZIP 的概率质量函数
zip_pmf <- function(y, pi, lambda) {
ifelse(y == 0,
pi + (1 - pi) * exp(-lambda),
(1 - pi) * (exp(-lambda) * lambda^y) / factorial(y))
}
# 参数:
# y: 观测计数值 (或值的向量)
# pi: 零膨胀概率
# lambda: 泊松分布的均值其中 π 为零膨胀概率，λ 为泊松均值。类似地，零膨胀负二项（ZINB）模型将泊松替换为负二项分布，以引入过离散参数。在实现时，可以使用最大似然估计或贝叶斯方法估计模型参数。例如，可以使用期望-条件最大化（ECME）算法迭代估计 π 和 λ，或使用MCMC采样进行贝叶斯推断。零膨胀模型可以用于改进丰度估计（通过区分“真正的零”和“未观测到的非零”）以及差异丰度分析（提供更准确的p值）。
5. 负二项分布
负二项分布主要用于建模过离散计数数据。负二项分布的概率质量函数为：
# P(Y=y) = C(y+r-1, y) * (1-p)^r * p^y
# 其中 C(n, k) 是组合数
# 定义一个函数来计算负二项分布的概率质量函数
# 在 R 中，有内置函数 dnbinom，但这里展示如何根据公式实现
nb_pmf_from_formula <- function(y, r, p) {
# choose(n, k) 计算 C(n, k)
choose(y + r - 1, y) * (1 - p)^r * p^y
}
# 参数:
# y: 观测计数值
# r: 成功次数参数
# p: 单次试验的成功概率其中 r 为成功次数参数，p 为成功概率，
# E[Y] = (p * r) / (1 - p)
mean_nb <- (p * r) / (1 - p)
# Var(Y) = (p * r) / (1 - p)^2
variance_nb <- (p * r) / (1 - p)^2当 r 固定时，方差随均值增大而增大，能够捕捉过离散特性。在实现时，可以使用矩估计或最大似然估计估计参数(r,p)，或使用贝叶斯方法（如Gamma-Negative Binomial共轭先验）进行推断。
参考文献
- Moeckel, C., Mareboina, M., Konnaris, M. A., et al. (2024). A survey of k-mer methods and applications in bioinformatics. Computational and Structural Biotechnology Journal, 23, 2289–2303.
- Rozday, T. J. (2025). Metagenomics distilled: new k-mer-based methods. Nature Reviews Microbiology, 23, 473.
- Shaw, J., & Yu, Y. W. (2024). Rapid species-level metagenome profiling and containment estimation with sylph. Nature Biotechnology, 43, 1348–1359.
- Shen, C., Wedell, E., Pop, M., & Warnow, T. (2025). TIPP3 and TIPP3-fast: Improved abundance profiling in metagenomics. PLoS Computational Biology, 21(4), e1012593.
- Lu, J., Breitwieser, F. P., Thielen, P., & Salzberg, S. L. (2017). Bracken: estimating species abundance in metagenomics data. PeerJ Computer Science, 3, e104.
- Lindner, M. S., & Renard, B. Y. (2012). Metagenomic abundance estimation and diagnostic testing on species level. Nucleic Acids Research, 41(1), e10.
- Nayfach, S., Bradley, P. H., Wyman, S. K., et al. (2015). Automated and Accurate Estimation of Gene Family Abundance from Shotgun Metagenomes. PLoS Computational Biology, 11(11), e1004573.
- Franzosa, E. A., McIver, L. J., Rahnavard, G., et al. (2018). Species-level functional profiling of metagenomes and metatranscriptomes. Nature Methods, 15(11), 962–968.
- Jonsson, V., Österlund, T., Nerman, O., & Kristiansson, E. (2019). Modelling of zero-inflation improves inference of metagenomic gene count data. Statistical Methods in Medical Research, 28(12), 3712–3728.
- Melsted, P., & Pritchard, J. K. (2011). Efficient counting of k-mers in DNA sequences using a bloom filter. BMC Bioinformatics, 12, 333.
- Zhang, Q., Pell, J., Canino-Koning, R., et al. (2014). These Are Not the K-mers You Are Looking For: Efficient Online K-mer Counting Using a Probabilistic Data Structure. PLoS ONE, 9(7), e101271.
- Ai, D., Pan, H., Huang, R., & Xia, L. C. (2018). CoreProbe: A Novel Algorithm for Estimating Relative Abundance Based on Metagenomic Reads. Genes, 9(6), 313.
- 宏基因组测序数据基因丰度估算方法理论 - 2025年12月8日
- 最低共同祖先（Lowest Common Ancestor, LCA）算法讲解 - 2025年12月2日
- 宏基因组去除宿主序列的主流算法原理与基于kmer的方法详解 - 2025年11月29日
No responses yet