利用 Tax4Fun2 对 16S 微生物组数据进行功能预测
16S rRNA 基因的高通量测序已被广泛用于研究各种海洋,地表和宿主相关环境中微生物群落的组成和结构。但许多生物学问题更需要我们研究其功能变化,而不仅仅是微生物分类组成。近年来,有不少研究团队开发了几款预测工具,例如 PICRUSt,Tax4Fun,Piphillin,Faprotax 和 paprica。尽管这些工具并不能替代宏基因组测序,但它们仍在一定程度上为我们提供了独特见解。
这些工具的预测能力往往取决于公共数据库中基因组的功能信息。但这些公开的基因组并不一定代表所研究的生态系统中存在的总功能多样性。鉴于公开基因组的数量迅速增加,尤其是一些通过宏基因组分箱得到的数据,未发表的数据,亦或是某些特定环境的基因组信息,所以纳入这些数据可以大大提高功能预测的准确性。
基于这些想法,Tax4Fun 迎来了新的版本升级 Tax4Fun2。Tax4Fun2 是一个快速且易用的 R 包,默认的参考数据集包含 275 个古细菌和 12,002 个细菌基因组。Tax4Fun2 最大的一个升级之处在于可以合并我们的自定义数据集,以增强功能预测的鲁棒性和特异性。
开发团队还将 Tax4Fun2 与 Tax4Fun、PICRUSt 进行比较。结果表明在所有数据集中,Tax4Fun2 的预测性能均优于 PICRUSt 和 Tax4Fun。
Tax4Fun2 预测的功能谱与从宏基因组数据分析得到的功能谱高度相关。尽管 Kelp 相关群落的预测特征与功能特征显著相关,但中位 Spearman 相关系数仅为 0.72,这表明缺少合适的参考数据可能限制了 Tax4Fun2 的性能。
为了解决这个问题,开发团队使用了来自 90 个海带宏基因组的 68 个 MAGs ,构建特定于海带的宏基因组数据集。结果表明这大大提高了预测的准确性(中位 Spearman 相关系数为 0.86),也减少了预测中未使用的序列比例。
此外,添加了自定义的海带数据集后,甚至还可以成功预测那些用参考数据集预测失败的样本。这也证实了添加特定环境参考数据库的优势之处,这种方法也使 Tax4Fun2 更加灵活,并使之在一众类似软件中脱颖而出。
微生物生态学中另一个主要的研究问题是,微生物群落是否包含功能冗余的成分,及其在何种程度上包含功能冗余,从而使之在面对多变的环境时保持生态系统的相对稳定。在 Tax4Fun2 中,针对单个功能通路引入了功能冗余指数(FRI)的概念。我们可以基于包含特定功能的物种比例及其相互之间的系统发育关系来计算 FRI。高 FRI 表示该功能几乎在所有成员中普遍存在,而低 FRI 则表示该功能仅存在于一些密切相关的物种中或仅在一个物种中被检测到,FRI 为 0 表示该功能不存在。Tax4Fun2 会计算相对 FRI(rFRI)和绝对 FRI(aFRI),前者会用环境中的平均系统发生距离进行归一化,后者则用 Tax4Fun2 参考数据中所有原核生物的平均系统发生距离进行归一化。简而言之,rFRI 可用于比较某一项研究中的样本,而 aFRI 允许比较不同生态系统中的功能冗余指数。
为了检测 FRI 的准确性,开发团队模拟了 1,000 个环境,每个环境都包含 100 个原核基因组。接着从每个模拟环境中提取了 16S rRNA 基因序列,将它们以 97% 的相似度聚类并计算了 FRI 值。随后根据模拟环境的实际基因组信息将这些值与 FRI 值进行了比较。结果表明, Tax4Fun2 对微生物群落中功能冗余的进行了较好的估计(Spearman 等级相关性> 90%)。此外,用海水样本进行进一步验证。其中六个样本是为赤潮期间采集,而另外三个为正常海水样本,作为参照数据。在参照样本中,将近 7,000 个功能显示出更高的功能冗余,而只有 1,468 个功能在赤潮样本中有更高的冗余度。
这表明在赤潮期间功能冗余发生了很大变化。浮游植物水华通常以基质控制的演替为特征,即在水华期间和之后的不同阶段,不同的细菌进化枝在浮游植物群落中占主导地位。因此,在特定阶段参与某些底物更新的环境成员占主导地位,因此它们的基因组和相关功能将更加冗余,而对于其他所有环境成员,情况则相反。
综上所述,Tax4Fun2 不但具有较高的预测能力,还可通过结合用户自定义数据集来进一步提高其准确性。除此之外,Tax4Fun2 还可以计算特定功能的冗余度,这对于预测在多变环境中定功能丢失的可能性非常重要。
Tax4Fun2 安装及使用
GitHub 地址:https://github.com/bwemheu/Tax4Fun2
Step1. 安装 Tax4Fun2 和下载参考数据库
安装 Tax4Fun2
下载:
wget https://github.com/bwemheu/Tax4Fun2/releases/download/1.1.5/Tax4Fun2_1.1.5.tar.gz
用 install.packages()
安装:
install.packages(pkgs = "Tax4Fun2_1.1.5.tar.gz", repos = NULL, source = TRUE)
加载 Tax4Fun2:
library(Tax4Fun2)
下载参考数据库
在 Tax4Fun2 中我们可用自带的 buildReferenceData()
命令来安装配置参考数据集:
buildReferenceData(path_to_working_directory = ".", use_force = FALSE, install_suggested_packages = TRUE)
参数
•path_to_working_directory
表示数据库安装位置,默认为当前目录;•use_force
表示是否覆盖,默认为 FALSE
;•install_suggested_packages
表示是否安装依赖包 ape
和 seqinr
,默认为 TURE
。
在 Windows 下可能出现数据库已下载但无法解压的情况,我们也可以手动下载然后解压,数据库下载地址:https://owncloud.gwdg.de/index.php/s/pii5Hoitr5SgDqK/download
安装依赖包
此命令将下载 blast(v2.9.0),并将文件夹放在 Tax4Fun2 参考数据库文件夹中。此外,它将测试 R 包 ape 和 seqinr 是否已正确安装。
buildDependencies(path_to_reference_data = "Tax4Fun2_ReferenceData_v2", install_suggested_packages = TRUE)
也可以直接手动下载 blast,然后将
blast_bin
目录放到Tax4Fun2_ReferenceData_v2
文件夹下,blast 下载地址:ftp://ftp.ncbi.nlm.nih.gov/blast/executables/blast+/LATEST/
下载测试数据
getExampleData(path_to_working_directory = ".")
Step 2. 生成自定义参考数据集
1. 提取 SSU 序列(16S rRNA 和 18S rRNA)
# Option A) 从单个基因组中提取 SSU 序列
extractSSU(genome_file = "OneProkaryoticGenome.fasta", file_extension = "fasta", path_to_reference_data = "Tax4Fun2_ReferenceData_v2")
# Option B) 从多个基因组中提取 SSU 序列
extractSSU(genome_folder = "MoreProkaryoticGenomes", file_extension = "fasta", path_to_reference_data = "Tax4Fun2_ReferenceData_v2")
参数
•genome_file
或 genome_folder
:指定单个基因组文件或包含多个基因组的文件夹。必须将多个基因组放在一个文件夹中(一个基因组一个文件),且所有文件为相同的文件扩展名;•file_extension = "fasta"
:指定 fasta
格式的基因组文件;•path_to_reference_data = "Tax4Fun2_ReferenceData_v2"
:指定参考数据集路径。
2. 为原核基因组进行功能注释
# Option A) 对单个基因组进行功能注释
assignFunction(genome_file = "OneProkaryoticGenome.fasta", file_extension = "fasta", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", num_of_threads = 1, fast = TRUE)
# Option B) 对多个基因组进行功能注释
assignFunction(genome_folder = "MoreProkaryoticGenomes/", file_extension = "fasta", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", num_of_threads = 1, fast = TRUE)
参数
•num_of_threads = 1
:设置 diamond 运行的线程数;•fast = TRUE
:用默认参数运行 diamond;若设置为 FALSE
则表示运行 diamond --sensitive
,可提高灵敏度,但速度要慢得多。
Mac 可能需要自定义 diamond 路径:
# Option A) 对单个基因组进行功能注释
assignFunction(genome_file = "OneProkaryoticGenome.fasta", file_extension = "fasta", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", num_of_threads = 1, fast = T, path_to_diamond_binary_mac = "diamond")
# Option B) 对多个基因组进行功能注释
assignFunction(genome_folder = "MoreProkaryoticGenomes/", file_extension = "fasta", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", num_of_threads = 1, fast = T, path_to_diamond_binary_mac = "diamond")
3. 生成自定义参考数据集
提取 SSU 序列并注释功能后,每个基因组都至少包括两个文件:SSU序列文件和功能谱文件。下一步就是用 generateUserData()
命令生成 自定义参考数据集:
# 1) 根据单个基因组生成参考数据集
generateUserData(path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_user_data = ".", name_of_user_data = "User_Ref0", SSU_file_extension = "_16SrRNA.ffn", KEGG_file_extension = "_funPro.txt")
# 2) 根据多个基因组生成参考数据集
generateUserData(path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_user_data = "MoreProkaryoticGenomes", name_of_user_data = "User_Ref1", SSU_file_extension = "_16SrRNA.ffn", KEGG_file_extension = "_funPro.txt")
# 3) 根据多个基因组生成参考数据集并用 uclust 去冗余
generateUserDataByClustering(path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_user_data = "MoreProkaryoticGenomes", name_of_user_data = "User_Ref2", SSU_file_extension = "_16SrRNA.ffn", KEGG_file_extension = "_funPro.txt", use_force = T)
参数
•path_to_refernce_data
:指定参考数据路径;•path_to_user_data
:指定自定义数据集路径;•name_of_user_data
:自定义数据集名称。
Step 3. 进行功能预测
这里需要把你的 OTU 表格式化成和示例数据中一样的格式。
PlanA 只用默认的参考数据集进行功能预测
# 1. 先进行 blast
runRefBlast(path_to_otus = "KELP_otus.fasta", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Kelp_Ref99NR", database_mode = "Ref99NR", use_force = T, num_threads = 6)
# 2. 再进行功能预测,normalize_pathways = FALSE
makeFunctionalPrediction(path_to_otu_table = "KELP_otu_table.txt", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Kelp_Ref99NR", database_mode = "Ref99NR", normalize_by_copy_number = TRUE, min_identity_to_reference = 0.97, normalize_pathways = FALSE)
# 或 normalize_pathways = TRUE
makeFunctionalPrediction(path_to_otu_table = "KELP_otu_table.txt", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Kelp_Ref99NR", database_mode = "Ref99NR", normalize_by_copy_number = TRUE, min_identity_to_reference = 0.97, normalize_pathways = TRUE)
参数
•path_to_otus
:指定 OTU 的 fasta 文件;•path_to_otu_table
:指定 OTU 矩阵路径;•path_to_refernce_data
:指定参考数据库路径;•path_to_temp_folder
:指定临时目录,结果将会存放在该目录下;•database_mode = "Ref99NR"
:选择参考数据库,可选用 ' Ref99NR 或 Ref100NR,数值表示用 uclust 聚类的阈值(分别为 99% 和 100%);•num_threads
:设置 diamond 运行的线程数;•use_force
:表示是否覆盖;•normalize_by_copy_number
:是否进行归一化;•normalize_pathways = FALSE
:将每个 KO 的相对丰度关联到它所属的每个通路。将其设置为 true,相对丰度将平均分配到所属的每个通路。
PlanB 结合默认数据库和自定义数据库(未聚类)进行功能预测
# 1. 生成自定义数据集,命名为 KELP1
generateUserData(path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_user_data = "KELP_UserData", name_of_user_data = "KELP1", SSU_file_extension = ".ffn", KEGG_file_extension = ".txt")
# 2. 用 include_user_data = TRUE,path_to_user_data,name_of_user_data 三个参数添加自定义数据集进行 blast
runRefBlast(path_to_otus = "KELP_otus.fasta", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Kelp_Ref99NR_withUser1", database_mode = "Ref99NR", use_force = T, num_threads = 6, include_user_data = T, path_to_user_data = "KELP_UserData", name_of_user_data = "KELP1")
# 3. 添加自定义数据集进行功能预测
makeFunctionalPrediction(path_to_otu_table = "KELP_otu_table.txt", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Kelp_Ref99NR_withUser1", database_mode = "Ref99NR", normalize_by_copy_number = T, min_identity_to_reference = 0.97, normalize_pathways = F, include_user_data = T, path_to_user_data = "KELP_UserData", name_of_user_data = "KELP1")
参数
•include_user_data = T
:包括自定义数据集;•name_of_user_data
:自定义数据集的名称;•path_to_user_data
:自定义数据集的路径。
PlanC 结合默认数据库和自定义数据库(聚类)进行功能预测
# 1. 生成自定义参考数据集并用 uclust 去冗余,命名为 KELP2
generateUserDataByClustering(path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_user_data = "KELP_UserData", name_of_user_data = "KELP2", SSU_file_extension = ".ffn", KEGG_file_extension = ".txt", similarity_threshold = 0.99)
# 2. 用 include_user_data = TRUE,path_to_user_data,name_of_user_data 三个参数添加自定义数据集进行 blast
runRefBlast(path_to_otus = "KELP_otus.fasta", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Kelp_Ref99NR_withUser2", database_mode = "Ref99NR", use_force = T, num_threads = 6, include_user_data = T, path_to_user_data = "KELP_UserData", name_of_user_data = "KELP2")
# 3. 添加自定义数据集进行功能预测
makeFunctionalPrediction(path_to_otu_table = "KELP_otu_table.txt", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Kelp_Ref99NR_withUser2", database_mode = "Ref99NR", normalize_by_copy_number = T, min_identity_to_reference = 0.97, normalize_pathways = F, include_user_data = T, path_to_user_data = "KELP_UserData", name_of_user_data = "KELP2")
Step 4. 计算(多)功能冗余指数 FRI
计算 KEGG 功能通路的系统发育分布(高 FRI -> 高冗余,低 FRI -> 低冗余度):
# 1. 进行 blast
runRefBlast(path_to_otus = "Water_otus.fna", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Water_Ref99NR", database_mode = "Ref99NR", use_force = T, num_threads = 6)
# 2. 计算 FRI
calculateFunctionalRedundancy(path_to_otu_table = "Water_otu_table.txt", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Water_Ref99NR", database_mode = "Ref99NR", min_identity_to_reference = 0.97)
# New in the latest pre-release (v1.1.6): prevalence_cutoff (see comment on pre-release)
calculateFunctionalRedundancy(path_to_otu_table = "Water_otu_table.txt", path_to_reference_data = "Tax4Fun2_ReferenceData_v2", path_to_temp_folder = "Water_Ref99NR", database_mode = "Ref99NR", min_identity_to_reference = 0.97, prevalence_cutoff = 1.0)
Have Fun!
参考链接
[1]
文章地址:https://www.biorxiv.org/content/10.1101/490037v1
[2]
GitHub 地址:https://github.com/bwemheu/Tax4Fun2
推荐阅读
-
利用 Tax4Fun2 对 16S 微生物组数据进行功能预测
-
windows下进程间通信的(13种方法)-摘 要 本文讨论了进程间通信与应用程序间通信的含义及相应的实现技术,并对这些技术的原理、特性等进行了深入的分析和比较。 ---- 关键词 信号 管道 消息队列 共享存储段 信号灯 远程过程调用 Socket套接字 MQSeries 1 引言 ---- 进程间通信的主要目的是实现同一计算机系统内部的相互协作的进程之间的数据共享与信息交换,由于这些进程处于同一软件和硬件环境下,利用操作系统提供的的编程接口,用户可以方便地在程序中实现这种通信;应用程序间通信的主要目的是实现不同计算机系统中的相互协作的应用程序之间的数据共享与信息交换,由于应用程序分别运行在不同计算机系统中,它们之间要通过网络之间的协议才能实现数据共享与信息交换。进程间通信和应用程序间通信及相应的实现技术有许多相同之处,也各有自己的特色。即使是同一类型的通信也有多种的实现方法,以适应不同情况的需要。 ---- 为了充分认识和掌握这两种通信及相应的实现技术,本文将就以下几个方面对这两种通信进行深入的讨论:问题的由来、解决问题的策略和方法、每种方法的工作原理和实现、每种实现方法的特点和适用的范围等。 2 进程间的通信及其实现技术 ---- 用户提交给计算机的任务最终都是通过一个个的进程来完成的。在一组并发进程中的任何两个进程之间,如果都不存在公共变量,则称该组进程为不相交的。在不相交的进程组中,每个进程都独立于其它进程,它的运行环境与顺序程序一样,而且它的运行环境也不为别的进程所改变。运行的结果是确定的,不会发生与时间相关的错误。 ---- 但是,在实际中,并发进程的各个进程之间并不是完全互相独立的,它们之间往往存在着相互制约的关系。进程之间的相互制约关系表现为两种方式: ---- (1) 间接相互制约:共享CPU ---- (2) 直接相互制约:竞争和协作 ---- 竞争——进程对共享资源的竞争。为保证进程互斥地访问共享资源,各进程必须互斥地进入各自的临界段。 ---- 协作——进程之间交换数据。为完成一个共同任务而同时运行的一组进程称为同组进程,它们之间必须交换数据,以达到协作完成任务的目的,交换数据可以通知对方可以做某事或者委托对方做某事。 ---- 共享CPU问题由操作系统的进程调度来实现,进程间的竞争和协作由进程间的通信来完成。进程间的通信一般由操作系统提供编程接口,由程序员在程序中实现。UNIX在这个方面可以说最具特色,它提供了一整套进程间的数据共享与信息交换的处理方法——进程通信机制(IPC)。因此,我们就以UNIX为例来分析进程间通信的各种实现技术。 ---- 在UNIX中,文件(File)、信号(Signal)、无名管道(Unnamed Pipes)、有名管道(FIFOs)是传统IPC功能;新的IPC功能包括消息队列(Message queues)、共享存储段(Shared memory segment)和信号灯(Semapores)。 ---- (1) 信号 ---- 信号机制是UNIX为进程中断处理而设置的。它只是一组预定义的值,因此不能用于信息交换,仅用于进程中断控制。例如在发生浮点错、非法内存访问、执行无效指令、某些按键(如ctrl-c、del等)等都会产生一个信号,操作系统就会调用有关的系统调用或用户定义的处理过程来处理。 ---- 信号处理的系统调用是signal,调用形式是: ---- signal(signalno,action) ---- 其中,signalno是规定信号编号的值,action指明当特定的信号发生时所执行的动作。 ---- (2) 无名管道和有名管道 ---- 无名管道实际上是内存中的一个临时存储区,它由系统安全控制,并且独立于创建它的进程的内存区。管道对数据采用先进先出方式管理,并严格按顺序操作,例如不能对管道进行搜索,管道中的信息只能读一次。 ---- 无名管道只能用于两个相互协作的进程之间的通信,并且访问无名管道的进程必须有共同的祖先。 ---- 系统提供了许多标准管道库函数,如: pipe——打开一个可以读写的管道; close——关闭相应的管道; read——从管道中读取字符; write——向管道中写入字符; ---- 有名管道的操作和无名管道类似,不同的地方在于使用有名管道的进程不需要具有共同的祖先,其它进程,只要知道该管道的名字,就可以访问它。管道非常适合进程之间快速交换信息。 ---- (3) 消息队列(MQ) ---- 消息队列是内存中独立于生成它的进程的一段存储区,一旦创建消息队列,任何进程,只要具有正确的的访问权限,都可以访问消息队列,消息队列非常适合于在进程间交换短信息。 ---- 消息队列的每条消息由类型编号来分类,这样接收进程可以选择读取特定的消息类型——这一点与管道不同。消息队列在创建后将一直存在,直到使用msgctl系统调用或iqcrm -q命令删除它为止。 ---- 系统提供了许多有关创建、使用和管理消息队列的系统调用,如: ---- int msgget(key,flag)——创建一个具有flag权限的MQ及其相应的结构,并返回一个唯一的正整数msqid(MQ的标识符); ---- int msgsnd(msqid,msgp,msgsz,msgtyp,flag)——向队列中发送信息; ---- int msgrcv(msqid,cmd,buf)——从队列中接收信息; ---- int msgctl(msqid,cmd,buf)——对MQ的控制操作; ---- (4) 共享存储段(SM) ---- 共享存储段是主存的一部分,它由一个或多个独立的进程共享。各进程的数据段与共享存储段相关联,对每个进程来说,共享存储段有不同的虚拟地址。系统提供的有关SM的系统调用有: ---- int shmget(key,size,flag)——创建大小为size的SM段,其相应的数据结构名为key,并返回共享内存区的标识符shmid; ---- char shmat(shmid,address,flag)——将当前进程数据段的地址赋给shmget所返回的名为shmid的SM段; ---- int shmdr(address)——从进程地址空间删除SM段; ---- int shmctl (shmid,cmd,buf)——对SM的控制操作; ---- SM的大小只受主存限制,SM段的访问及进程间的信息交换可以通过同步读写来完成。同步通常由信号灯来实现。SM非常适合进程之间大量数据的共享。 ---- (5) 信号灯 ---- 在UNIX中,信号灯是一组进程共享的数据结构,当几个进程竞争同一资源时(文件、共享内存或消息队列等),它们的操作便由信号灯来同步,以防止互相干扰。 ---- 信号灯保证了某一时刻只有一个进程访问某一临界资源,所有请求该资源的其它进程都将被挂起,一旦该资源得到释放,系统才允许其它进程访问该资源。信号灯通常配对使用,以便实现资源的加锁和解锁。 ---- 进程间通信的实现技术的特点是:操作系统提供实现机制和编程接口,由用户在程序中实现,保证进程间可以进行快速的信息交换和大量数据的共享。但是,上述方式主要适合在同一台计算机系统内部的进程之间的通信。 3 应用程序间的通信及其实现技术 ---- 同进程之间的相互制约一样,不同的应用程序之间也存在竞争和协作的关系。UNIX操作系统也提供一些可用于应用程序之间实现数据共享与信息交换的编程接口,程序员可以通过自己编程来实现。如远程过程调用和基于TCP/IP协议的套接字(Socket)编程。但是,相对普通程序员来说,它们涉及的技术比较深,编程也比较复杂,实现起来困难较大。 ---- 于是,一种新的技术应运而生——通过将有关通信的细节完全掩盖在某个独立软件内部,即底层的通讯工作和相应的维护管理工作由该软件内部来实现,用户只需要将通信任务提交给该软件去完成,而不必理会它的具体工作过程——这就是所谓的中间件技术。 ---- 我们在这里分别讨论这三种常用的应用程序间通信的实现技术——远程过程调用、会话编程技术和MQSeries消息队列技术。其中远程过程调用和会话编程属于比较低级的方式,程序员参与的程度较深,而MQSeries消息队列则属于比较高级的方式,即中间件方式,程序员参与的程度较浅。 ---- 4.1 远程过程调用(RPC)