展开

Kraken序列分类

最后发布时间 : 2024-05-25 08:48:08 浏览量 :

学习资料

淘宝商品: 转录组自动化分析流程

淘宝商品: 转录组自动化分析流程

淘宝商品: 转录组自动化分析流程

背景

宏基因组学是研究直接从环境中获得的基因组序列。在研究了海水、酸性矿井排水和人体等各种环境的项目中,宏基因组学使研究人员能够创建环境微生物生命的图景,而无需分离和培养单个微生物。结合快速对DNA进行测序的能力,宏基因组学项目可以生成大量的序列数据来描述这些以前不可见的世界。对于许多宏基因组样品,在测序时,样品中存在的物种、属甚至门在很大程度上是未知的,测序的目的是尽可能精确地确定这种微生物组成。当然,如果一个生物体与以前见过的任何东西完全不同,那么它的DNA序列就不能被表征,只能将其标记为新奇的。然而,许多物种与已知物种具有一些可检测的相似性,并且这种相似性可以通过灵敏的比对算法来检测。最著名的此类算法,也是将分类标签分配给未知序列的最佳方法之一,是BLAST程序。它可以通过找到与大型基因组序列数据库的最佳比对来对序列进行分类。虽然BLAST不是为宏基因组序列设计的,但它很容易适应这个问题,它仍然是最好的方法之一。

其他序列分类方法,利用序列比对和机器学习技术,试图提高BLAST的准确性

  • 在 MEGAN程序中,针对多个数据库搜索序列(使用 BLAST),并将针对每个数据库的最佳匹配的最低共同祖先 (LCA) 分配给该序列。
  • PhymmBL将 BLAST 的结果与插值马尔可夫模型产生的分数相结合,实现了比单独使用 BLAST 更高的精度。
  • 朴素贝叶斯分类器(NBC)将贝叶斯规则应用于基因组中k-mers的分布。

然而,所有这些程序的运行速度都比BLAST慢,而BLAST本身需要大量的CPU时间来对齐典型的Illumina测序运行产生的数百万个序列。这种处理负担要求如此之高,以至于它提出了另一种更快的宏基因组序列分析方法:丰度估计(MetaPhlAn )。

丰度估计程序的工作原理是创建一个比所有基因组集合小得多的数据库,这些数据库经过精心设计,包含“标记”基因(几乎存在于所有微生物中的单拷贝基因),或已发现对某些分支具有特异性的基因。由于数据库仅包含每个基因组的非常小的样本,因此这些程序只能从典型的宏基因组学样本中对一小部分序列进行分类。它们旨在用于表征给定样品中存在的生物体的分布,而不是标记每个读数。尽管丰度估计程序提供了宏基因组的摘要级表征,但它们无法帮助需要有关样本的更多细节的分析。例如,它们不能用于估计样品中的基因含量,因为这需要将每次读数与已知基因进行比较。

Kraken,一种新的序列分类工具,其准确性可与最佳序列分类技术相媲美,其速度远远超过分类器和丰度估计程序。这种速度优势很大程度上源于对 k-mer 的精确匹配数据库查询的使用,而不是序列的不精确对齐。它的准确性之所以成为可能,是因为测序的微生物基因组数量非常大且仍在增长,目前数量超过8,500个。Kraken 能够实现与最快的 BLAST 程序 Megablast 获得的非常相似的属级灵敏度和精度。

原理

为了对序列进行分类,序列中的每个 k-mer 都映射到数据库中包含该 k-mer 的基因组的最低共同祖先 (LCA)。也就是说数据库中的每个 k-mer都有一个分类标签。
与序列的 k-mers 相关的分类群以及分类群的祖先构成了一般分类树的修剪子树,用于分类。将序列剪切为k-mer后,用这些短的k-mer与数据库中的k-mer进行比对,reads的k-mer会得到一个分类标签,并且reads的k-mer都有权重(分类下reads k-mer的数量)。这些reads的k-mer的分类标签与分类群的祖先构成了修建的子树
在分类树中,每个节点的权重等于与节点分类单元关联的序列中的 k-mer 数量。(分类下reads k-mer的数量)
分类树中的每个根到叶 (RTL) 路径通过添加路径中的所有权重来评分,分类树中的最大 RTL 路径是分类路径(以黄色突出显示的节点)。
此分类路径的叶子(分类树中最左边的橙色叶子)是用于查询序列的分类。

生信小木屋

介绍

Kraken 是一个分类序列分类器,分配短 DNA 读取分类标签。 它通过检查 read 中的 k-mers 并使用这些 k-mers 查询数据库来实现这一点。这个数据库包含了Kraken基因组文库中的每一个 k-mer 到包含 k-mer 的所有基因组的分类树中的最近公共祖先(LCA)的映射。然后分析与读取中的 k-mers 对应的 LCA 分类群集合,为读取创建单个分类标签;这个标签可以是分类树中的任何节点。Kraken被设计成快速、敏感和高度精确。我们对各种实际和模拟数据的测试表明,Kraken 的灵敏度略低于 Megablast,精度略高。在一组模拟的100bp 读数上,Kraken 在正常操作下每分钟在一个核上处理超过130万个读数,在快速操作下每分钟处理超过410万个读数。

系统要求

注意: 关心磁盘或内存需求的用户应该阅读下面关于 MiniKraken 的段落。

磁盘空间: 到2017年10月,建造 Kraken 的标准数据库至少需要500GB 的磁盘空间。自定义数据库可能需要或多或少的空间。构建完成后,所需的最小数据库文件需要大约200GB 的磁盘空间。使用的磁盘空间与不同 k-mers 的数量成线性关系; 截至2017年10月,Kraken 的默认数据库包含大约140亿(1.4 e9)不同 k-mers。

此外,用于存储数据库的磁盘应该是本地附加存储器。将数据库存储在网络文件系统(NFS)分区上会导致 Kraken 的操作非常慢,或者完全停止。由于 NFS 访问比本地磁盘访问慢得多,因此使用 NFS 会降低预加载和数据库构建的速度。

内存: 为了有效地运行,Kraken 需要足够的空闲内存来将数据库保存在 RAM 中。虽然这可以通过使用内存磁盘完成,Kraken 提供了一个实用程序,通过操作系统缓存将数据库加载到内存中。默认数据库大小为174 GB (截至2017年10月) ,因此如果您想使用默认数据库构建或运行,那么至少需要这么大的 RAM。

依赖关系: Kraken 目前广泛使用 Linux 实用程序,如 sed、 find 和 wget。许多脚本是使用 Bashshell 编写的,而主要脚本是使用 Perl 编写的。构建数据库和运行分类器所需的核心程序是用 C + + 编写的,需要使用 g + + 进行编译。使用 OpenMP 处理多线程。NCBI 数据的下载由 wget 执行,在某些情况下由 rsync 执行。大多数安装了任何类型的开发包的 Linux 系统都可以使用上面列出的所有程序和库。

最后,如果您想构建自己的数据库,则需要安装 Jellyfish k-mer 计数器。请注意,Kraken 只支持使用水母版本1。水母版本2与Kraken不兼容。

网络连接: Kraken 的标准数据库构建和下载命令希望能够不受限制地访问 NCBI FTP 服务器的 FTP 和 rsync。如果您在代理后面工作,可能需要设置某些环境变量(如 ftp _ xy 或 RSYNC _ PROXY) ,以使这些命令正常工作。

MiniKraken: 为了让使用低内存计算环境的用户使用 Kraken,我们提供了一个简化的标准数据库,可以从 Kraken 网站下载。当Kraken 使用简化的数据库运行时,我们称之为“MiniKraken”。
我们提供的数据库只有4GB 和8GB 的大小,在只有8GB 和16GB RAM (分别)的计算机上应该能够很好地运行。每个 MiniKraken 数据库所需的磁盘空间也只有4GB 或8GB。

安装

要开始使用 Kraken,首先需要安装它,然后下载或创建一个数据库。

Kraken 包括两个主要的脚本(“ Kraken”和“ Kraken-build”) ,以及一些程序和较小的脚本。作为安装过程的一部分,所有脚本和程序都安装在同一个目录中。安装之后,您可以将主脚本移动到其他地方,但是移动其他脚本和程序需要编辑脚本并更改“ $KRAKEN _ DIR”变量。

一旦选择了一个目录,您需要在解压 Kraken 源代码的目录中运行以下命令:

./install_kraken.sh $KRAKEN_DIR

(将上面的“ $KRAKEN _ DIR”替换为您想要安装 KRAKEN 的程序/目录的目录。)
Install _ kraken.Sh 脚本应该编译所有的 Kraken 的代码和设置您的 Kraken 数据目录。如果您看到消息“ Kraken 安装完成”,则安装成功

一旦安装完成,您可能希望将两个主 Kraken 脚本复制到 PATH 变量中的目录(例如,“ $HOME/bin”) :

cp $KRAKEN_DIR/bin/kraken $HOME/bin
cp $KRAKEN_DIR/bin/kraken-build $HOME/bin

安装完成后,就可以创建或下载数据库了。

Kraken Databases

Kraken 数据库是一个至少包含4个文件的目录:

database.kdb: Contains the k-mer to taxon mappings
database.idx: Contains minimizer offset locations in database.kdb
taxonomy/nodes.dmp: Taxonomy tree structure + ranks
taxonomy/names.dmp: Taxonomy names
其他文件可能作为数据库生成过程的一部分出现。

在与 Kraken 交互时,您不必直接引用这些文件中的任何一个,只需提供存储它们的目录的名称即可。Kraken 允许使用标准数据库和自定义数据库; 这些分别在下面的Standard Kraken Database和Custom Databases 库章节中描述。

Standard Kraken Database

注意: 建立标准的 Kraken 数据库下载和使用所有完整的细菌,古细菌和病毒基因组在 Refseq 建立时。截至2017年10月,这包括约25,000个基因组,需要33 GB 的磁盘空间。然后,构建过程将需要大约450GB 的额外磁盘空间。在构建标准数据库之后,使用数据库将只需要用户保留 database. idx、 database. kdb 和 taxonomy/files,这需要大约200GB 的磁盘空间。在此数据库上运行示例时,用户将需要175GB 的 RAM。如果您没有这些计算资源,或者需要针对约25,000个基因组的 Refseq 数据库进行测试,我们建议只使用应用程序所需的基因组构建自定义数据库。

要创建标准的 Kraken 数据库,可以使用以下命令:

kraken-build --standard --db $DBNAME

(将上面的“ $DBNAME”替换为您首选的数据库名称/位置。)

这将下载 NCBI 分类信息,以及 RefSeq 中细菌、古细菌和病毒域的完整基因组。下载完所有这些数据后,开始构建过程; 这是最耗时的步骤。如果您有多个处理核心,您可以使用多个线程运行此进程,例如:

kraken-build --standard --threads 24 --db $DBNAME

2017年10月,在一台拥有244GB 内存的计算机上使用24个线程,构建过程花费了大约5个小时(带星号的步骤启用了一些多线程)。请注意,建立数据库所需的时间取决于基因组序列的数量:

  24m50s  *Step 1 (create kmer set)
     n/a  Step 2 (reduce database, optional and skipped)
2h34m53s  *Step 3 (sort set)
     n/a  Step 4 (GI number to sequence ID map)
   0.17s  Step 5 (Sequence ID to taxon map)
 2h7m28s  *Step 6 (set LCA values)
--------
 5h7m11s  Total build time

注意,如果任何步骤(包括初始下载)失败,构建过程将中止。但是,Kraken-build 将在整个安装过程中产生检查点,如果您试图在部分构建的数据库上再次运行相同的命令,则将在最后一个不完整的步骤中重新启动构建。

建立数据库后,若要删除任何不必要的文件(包括不再需要的库文件) ,请运行以下命令:

kraken-build --db $DBNAME --clean

若要创建自定义数据库或使用其他来源的数据库,请参见 Custom Databases。
内存数量较少的用户须知:

如果您遇到 Jellyfish 无法在您的系统上分配足够的内存来运行构建过程的问题,您可以使用 Kraken-build 的—— Jellyfish-hash-size 开关为 Jellyfish 提供一个更小的 hash 大小。哈希表中的每个空间大约使用6.9个字节,因此使用“—— jellyfish-hash-size 6400M”将使用64亿个空间的哈希表,并且需要44.3 GB 的 RAM。

Kraken 的构建过程通常试图通过分配大块 RAM 并在其中操作来最小化磁盘写入,直到需要将数据写入磁盘。但是,这种额外的 RAM 使用可能会超出您的容量。在这种情况下,您可能需要使用 Kraken-build 的—— work-on-disk 交换机。这将最大限度地减少 RAM 使用量,并导致 Kraken 的构建程序执行大多数磁盘文件外的操作。这个开关对于构建在内存磁盘或固态驱动器上的人也很有用。请注意,在某些计算机上,处理磁盘文件的速度可能相当慢,导致构建需要几天甚至几周的时间。

Custom Databases

我们意识到标准的数据库可能不适合每个人的需要。 Kraken 也允许创建定制的数据库。

安装一个分类法。通常,你只需要使用 NCBI 分类法,你可以很容易地下载使用:

kraken-build --download-taxonomy --db $DBNAME

这将从 NCBI 下载序列 ID 到分类群地图,以及分类名称和树信息。这些文件可以在 $DBNAME/taxonomy/中找到。如果需要修改分类法,可以编辑该目录中的 names.dmp 和 nodes.dmp 文件; gi _ taxid _ nucl.Dmp 文件也需要适当地更新。

安装一个基因组文库。四套标准基因组可以很容易地通过kraken-build获得:
bacteria: RefSeq complete bacterial/archaeal genomes
plasmids: RefSeq plasmid sequences
viruses: RefSeq complete viral genomes
human: GRCh38 human genome

要下载并安装其中任何一个,请使用—— download-library 开关,例如:

kraken-build --download-library bacteria --db $DBNAME

FASTA/multi-FASTA 文件中的其他基因组也可以添加:
如果从 NCBI 下载,基因组可以直接使用—— add-to-library 开关添加,例如:

kraken-build --add-to-library chr1.fa --db $DBNAME
kraken-build --add-to-library chr2.fa --db $DBNAME

请注意,如果您有要添加的文件列表,您可以在 bash 中执行类似的操作:

for file in chr*.fa
do
    kraken-build --add-to-library $file --db $DBNAME
done

或者甚至添加在基因组目录中找到的所有 * . fa 文件:

find genomes/ -name '*.fa' -print0 | \
    xargs -0 -I{} -n1 kraken-build --add-to-library {} --db $DBNAME

(如果有多个处理器,xargs 的-P 选项对于并行添加许多文件也很有用。)
没有从 NCBI 下载的复制品可能需要显式地分配它们的分类信息。这可以使用字符串 kraken: taxid | XXX 在序列 ID 中完成,用所需的分类单元 ID 替换 XXX。例如,要将已知的适配器序列放入分类单元32630(“合成构造”) ,可以使用以下方法:

>sequence16|kraken:taxid|32630  Adapter sequence
CAAGCAGAAGACGGCATACGAGATCTTCGAGTGACTGGAGTTCCTTGGCACCCGAGAATTCCA

Taxid 字符串必须以序列 ID 开始,或者前面立即加上一个管道字符(|)。以这种方式显式分配分类 ID 将覆盖 NCBI 提供的序列 ID 映射。

一旦您的库完成,您需要构建数据库。根据您的尺寸要求,您可能需要从默认值调整 k-mer 和/或最小化长度。除了一些小的簿记字段,Kraken 数据库将使用 sD + 8(4M)字节,其中 s 是用于存储 k-mer/taxon 对的字节数(通常为12,但对于较小的 k-mers 更低) ,D 是库中不同 k-mers 的数量,M 是最小化子的长度(以 bp 为单位)。虽然 D 确实随着 k 的增加而增加,但是如果不实际执行计数,就不可能确切地知道对于给定的 k,库中将存在多少个不同的 k-mer。默认情况下,k = 31,M = 15。

最小化器用于使查询序列中相邻的 k-mers 在数据库中彼此靠近,这使得 Kraken 可以利用 CPU 缓存。改变 M 的值可以显著影响Kraken的速度,增加或减少 M 都不能保证速度更快或更慢。

要构建数据库,您将使用—— build 开关:

kraken-build --build --db $DBNAME

如上所述,您可能还需要使用任何—— thread、—— kmer-len 或—— miniizer-len 来调整数据库构建时间和/或最终大小。

收缩数据库: “收缩”任务允许您获取一个现有的 Kraken 数据库,并从中创建一个较小的 MiniKraken 数据库。使用这个选项可以删除除了指定数量的 k-mer/taxon 对之外的所有数据,从而创建一个新的、更小的数据库。例如:

kraken-build --shrink 10000 --db $DBNAME --new-db minikraken

这将创建一个名为 minikraken 的新数据库,其中包含从原始数据库($DBNAME)中选择的10000k-mers。

收缩任务只能在已完成的数据库上运行。但是,如果在创建数据库之前就知道只能使用一定数量的内存,那么可以使用—— build 任务的—— max-db-size 开关为数据库提供最大大小(以 GB 为单位)。这允许您创建一个 MinKraken 数据库,而不必首先创建一个完整的 Kraken 数据库。

使用 Kraken-build —— help 可以获得 Kraken-build 的完整选项列表。

构建数据库后,如果想减少数据库的磁盘使用量,可以使用 kraken-build 的—— clean 开关从数据库目录中删除所有中间文件。

Classification

要对一组序列进行分类(读取) ,请使用 Kraken 命令:

kraken --db $DBNAME seqs.fa

默认情况下,输出将被发送到标准输出。应该在命令行中指定包含要分类的序列的文件。序列也可以通过使用特殊文件名/dev/fd/0的标准输入提供。

注意,为了获得最佳速度,Kraken 的数据库应该首先加载到 RAM 中。如果您有超级用户权限,可以通过使用内存磁盘来实现这一点。如果没有,您可以使用—— preload 开关到 Kraken,例如:

kraken --preload --db $DBNAME seqs.fa

数据库文件将在使用此开关进行分类之前加载。有关详细信息,请参阅内存使用情况和效率。

kraken项目允许几种不同的选择:

多线程: 使用—— thread NUM 开关来使用多个线程。
快速操作: 不是按顺序搜索所有 k-mers,而是在第一个数据库命中后停止分类; 使用 --quick 启用此模式。请注意—— min-hit 允许您在声明序列分类之前需要多次命中,在测试序列是否属于特定的基因组时,这对定制数据库尤其有用。

序列过滤: 已分类或未分类的序列可以分别使 --classified-out and --unclassified-out交换机发送到文件以便稍后处理。

输出重定向: 可以使用标准 shell 重定向(| 或 >)或使用—— Output 开关定向输出。

FASTQ 输入: 输入通常以 FASTA 格式进行,但是您可以使用 --fastq-input 开关对 FASTQ 数据进行分类。

压缩输入: Kraken 可以将 gzip 和 bzip2压缩文件作为输入处理,方法是指定适当的—— gzip 压缩或—— bzip2压缩开关。

输入格式自动检测: 如果在命令行中指定常规文件作为输入,Kraken 将尝试在分类之前确定输入的格式。可以通过适当地显式指定—— fast-input、—— fastq-input、—— gzip 压缩和/或—— bzip2压缩来禁用它。注意,使用字符设备文件/dev/fd/0从标准输入(又名 stdin)读取将不允许自动检测。

配对阅读: Kraken 不查询含有不明确核苷酸(非 ACGT)的 k 聚体。如果你有成对的读数,你可以利用这个事实,通过在序列之间用一个单一的 N 连接成对来提高 Kraken 的准确性。在运行 Kraken 时使用—— double 选项将自动为您做到这一点; 只需在命令行上指定两个配对文件。我们发现这样可以提高灵敏度,比将序列分类为单端读数提高大约3个百分点。有关配对读取输入/输出的更多信息,请参见配对读取

要得到完整的选项列表,请使用Kraken --help。

Output Format

Kraken分类的每个序列都产生一行输出。输出行包含五个以制表符分隔的字段; 从左到右,它们是:

  • “ C”/“ U”: 一个字母代码,表示该序列是classified or unclassified。
  • 序列 ID,从 FASTA/FASTQ 头获得。
  • Kraken使用的记序列的分类 ID ; 如果序列未分类,则为0。
  • 序列的长度(以 bp 为单位)。
  • 一个空格分隔的列表,指示序列中每个 k-mer 的 LCA 映射。例如,“562:13 561:4 A:31 0:1 562:3”表示:
    • 映射到分类 ID # 562的前13k-mers
    • 接下来的4个 k-mers 映射到分类 ID # 561
    • 接下来的31个 k 聚合体含有一个模棱两可的核苷酸
    • 下一个 k-mer 不在数据库里
    • 最后3个 k-mers 映射到分类 ID # 562

对于希望获得与每个输入序列相关联的完整分类名称的用户,我们提供了一个名为 kraken-trans 的脚本,该脚本为分类序列生成两种不同的输出格式。该脚本操作的kraken的输出,如下所示:

kraken --db $DBNAME sequences.fa > sequences.kraken
kraken-translate --db $DBNAME sequences.kraken > sequences.labels

(用于运行 Kraken 的数据库应该用于翻译输出; 请参阅下面的 Kraken Environment Variables 以了解减少命令行冗余的方法。)

上面的示例生成的文件 Sequences.label 是一个文本文件,其中包含两个以制表符分隔的列,并且在 Sequences.fa 中每个分类序列都有一行; 不分类的序列不会由 kraken-trans 报告。Kraken 翻译输出的第一列是分类序列的序列 ID,第二列包含序列的分类。例如,Kraken 的输出一行:

C     SEQ1    562     36      562:6

将产生相应的 Kraken 翻译输出一行:

SEQ1  root;cellular organisms;Bacteria;Proteobacteria;Gammaproteobacteria;Enterobacteriales;Enterobacteriaceae;Escherichia;Escherichia coli

或者,Kraken-trans 接受这个选项—— mpa-format,它只报告分类层次和标准等级分配(超级王国,王国,门,类,目,科,属,物种) ,并使用管道来划分分类层次。例如,kraken-translate --mpa-format --db $DBNAME 与上面的 Kraken 示例输出将导致以下输出行:

SEQ1  d__Bacteria|p__Proteobacteria|c__Gammaproteobacteria|o__Enterobacteriales|f__Enterobacteriaceae|g__Escherichia|s__Escherichia_coli

superkingdom (d _ _)级别以上的分类分配在使用—— mpa-report 选项和 kraken-trans 时仅表示为“ root”。

Paired Reads

当用户指 --paired 项时,Kraken 将对成对读取进行分类,方法是首先使用 | 对读取进行连接,然后再根据 Kraken 数据库对组合读取进行分类。
Kraken v1.0包含了许多其他选项,简化了配对读数的分析。下面描述这些选项,并列出这些选项的可能组合及其应用时的行为。请注意,所有选项都要求指定——成对选项,并提供两个输入 FASTA/FASTQ 文件。

--out-fmt legacy: [ default ]使用 N 作为序列分隔符,如果分类/非分类读取是使用 --classified-out or --unclassified-out標签打印的。--out-fmt legacy遗留问题目前不支持 FASTQ 输出。
--out-fmt legacy --classified-out C_reads.fa: 打印分类成对读数 N 连接两个成对读数。
--out-fmt paired:当使用--classified-out or --unclassified-out 标记时,将配对序列分离成两个独立的 FASTA 文件。
--out-fmt paired --fastq-output:当使用--classified-out or --unclassified-out 标记时,将成对序列分离成两个独立的 FASTQ 文件。FASTQ 报头将包括原始 FASTQ 报头中第二个空白字符的所有内容。
--out-fmt paired --classified-out C_reads: 将分类的成对读取打印到 FASTA 文件 C_reads_R1.fa and C_reads_R2.fa
--out-fmt interleaved:将配对序列打印到一个 FASTA 文件中,而不连接配对读取; 相反,将一个接一个地打印配对读取。

Memory Usage and Efficiency

Kraken的执行需要很多随机访问一个非常大的文件。为了获得最大的速度,这些访问需要尽快完成。这意味着在执行期间数据库必须在物理内存中。尽管我们为不能使用 ramdisk 的用户提供了—— preload 选项,ramdisk 可能是最简单的选项,并且非常适合在 Kraken 大部分时间都在运行的计算机上安装。此外,使用内存磁盘可以更快地完成 Kraken 的初始启动。如果使用了内存磁盘,则不应使用—— preload 开关。

我们还注意到,在某些情况下,可能不需要—— preload (甚至是建议的)。如果您知道数据库已经在内存中(例如,如果它最近被读取或解压缩,那么它应该在驻留在物理内存中的操作系统缓存中) ,那么就不需要执行这个步骤。我们已经注意到,在低内存(大约8GB)的情况下,预加载 MiniKraken DB 实际上比简单地使用 cat MiniKraken/database 要慢得多。* >/dev/null.选择将数据库存入内存的最佳方式取决于几个因素,包括 RAM 总量、操作系统和当前可用内存。出于这个原因,您可能需要对自己的设置进行试验,以便为您找到一个好的解决方案。

要创建一个 ramdisk,您需要有超级用户(root)权限。作为 root 用户,您可以使用以下命令创建一个 ramdisk:

mkdir /ramdisk
mount -t ramfs none /ramdisk

或者,您可能有一个希望能够将数据库复制到此目录中的受信任用户。在这种情况下,您需要通过 chown 使该用户成为目录的所有者。

要将数据库放到 ramdisk 上,只需将数据库目录复制到 ramdisk 目录:

cp -a $DBNAME /ramdisk

然后你可以通过指定数据库拷贝在内存盘上使用它,例如:

注意,如果 Ramdisk 被卸载或计算机重新启动,那么复制到 ramdisk 中的任何内容都将被删除,因此请确保在硬盘(或其他非易失性存储器)上有数据库的副本。

请注意,当使用——成对选项时,Kraken 将不会(默认情况下)尝试确保您指定的两个文件确实匹配成对结束读取集。为了验证每个读取的名称是否确实匹配,可以将--check-names 选项与 --paired 选项组合使用。

Sample Reports

为了了解整个样本中Kraken的结果,我们提供了Kraken报告脚本。它是这样使用的:

kraken-report --db $DBNAME kraken.output

请注意,所使用的数据库必须与生成输出文件所使用的数据库相同,否则报表脚本可能会遇到问题。输出被发送到标准输出。

Kraken-report 的输出以制表符分隔,每个分类单元一行。从左到右的输出字段如下:

  • 植根于该分类群的进化枝所覆盖的读数百分比
  • 植根于此分类单元的进化枝所覆盖的读数
  • 直接分配给该分类群的读数
  • 一种等级代码,表示(U)未分类,(D) omain,(K) Kingdom,(P) hylum,(C) lass,(O) rder,(F) family,(G) enus,或(S) peces。所有其他等级都是简单的“-”
  • NCBI生物分类编号
  • 缩进的科学名

根据分类法指定的树结构,科学名称使用空格缩进。
默认情况下,分类群如果没有分配给(或在)它们的读数,将不会产生任何输出。但是,如果希望显示所有分类群,可以使用—— show-zeros 开关进行显示。如果您希望对报告进行进一步的下游分析,并且希望对样本进行比较,那么这将非常有用。按分类 ID (使用 sort-nf5)排序可以在报表之间提供一致的行排序。

此外,我们还提供了 Kraken-mpa-report 程序; 该程序提供的输出格式类似于 MetaPhlAn 的 tab 分隔输出。对于 Kraken-mpa-report,可以在命令行上指定多个 Kraken 输出文件,每个文件将被视为一个单独的示例。对于标准等级的每个分类群(从域到种) ,显示分配给该分类群的进化枝中任何节点的每个样本中的读取计数。Kraken-mpa-report 的运行方式与 Kraken-report 相同,其输出也被发送到标准输出。

Confidence Scoring

目前,我们还没有建立一个可靠的海怪概率解释的信心分数。然而,我们已经开发了一个简单的计分方案,已经为我们产生了良好的结果,我们已经在 Kraken-filter 脚本中提供了这个方案。我们使用的方法允许用户在[0,1]间隔内指定一个阈值评分; 然后 Kraken-filter 脚本将向上调整标签,直到标签的评分(如下所述)满足或超过该阈值。如果分类树根部的标签的得分没有超过阈值,则序列被 Kraken-filter 调用为未分类。

序列标签的得分是分数 C/Q,其中 C 是映射到根于标签的进化枝中 LCA 值的 k 聚体的数目,Q 是序列中缺乏模糊核苷酸的 k 聚体的数目(即,它们针对数据库进行了查询)。考虑前面给出的 Kraken 输出中的 LCA 映射的例子:

“562:13561:4 A: 310:1562:3”表示:

  • the first 13 k-mers mapped to taxonomy ID #562
  • the next 4 k-mers mapped to taxonomy ID #561
  • the next 31 k-mers contained an ambiguous nucleotide
  • the next k-mer was not in the database
  • the last 3 k-mers mapped to taxonomy ID #562

在本例中,ID # 561是 # 562的父节点。在这里,这个序列的 # 562标签的得分是 C/Q = (13 + 3)/(13 + 4 + 1 + 3) = 16/21。一个 # 561的标签的 C/Q 值为(13 + 4 + 3)/(13 + 4 + 1 + 3) = 20/21。如果用户指定的阈值超过16/21,Kraken-filter 会将原始标签从 # 562调整到 # 561; 如果阈值大于20/21,序列将变得不分类。

Kraken 过滤器是这样使用的:

kraken-filter --db $DBNAME [--threshold NUM] kraken.output

如果未指定,则阈值为0。Kraken-filter 的输出类似于 Kraken 的输出,但是在长度和 LCA 映射列表之间出现了一个新的字段,表示新标签的得分(或根标签的得分,如果序列未分类)。

为了给选择一个合适的阈值提供一些指导,我们在这里展示了 Kraken 论文中 MiSeq 宏基因组上不同阈值的结果(更多细节见本文; 注意这里使用的数据库比本文中使用的数据库更新)。精确度,灵敏度和 F 评分是按属级测量的:

Thres	Prec	Sens	F-score
0	95.43	77.32	85.43
0.05	97.28	76.31	85.53
0.10	98.25	75.13	85.15
0.15	98.81	73.87	84.54
0.20	99.13	72.82	83.96
0.25	99.38	71.74	83.33
0.30	99.55	70.75	82.71
0.35	99.61	69.53	81.90
0.40	99.66	68.35	81.09
0.45	99.70	66.93	80.09
0.50	99.71	65.49	79.06

可以看到,没有阈值(比如,Kraken 的原始标签) ,Kraken 的精确度相当高,但是它确实随着阈值的增加而增加。然而,在决定自己的项目所使用的阈值时,必须考虑到敏感性报酬递减的损失。