Archives

Categories

合并和拆分双端序列:seqtk_interleave和seqtk_deinterleave

seqtk_interleave和seqtk_deinterleave 为 seqtk_utils 系列工具中的一组对应工具,即合并PE双端序列与拆分Interleaved格式的fasta/q文件。

1. seqtk_interleave

很多应用程序需要 Interleaved 格式序列作为输入文件,比如 SortMeRNA 很多工具都可以解决这个问题,比如: seqtk mergepe 可以很方便的将两个fq文件的R1和R2转换成下列文件格式:

@HISEQ:573:HLWMHBCXY:1:1101:11739:2308 1:N:0:TGGACTAT CCTATGGGGTGCAGCAGTGGGGAATATTGGACAATGGGGGGAACCCTGAT + IIIIIIIIIIIIIIIIIIIIIIHIIIIIIIIIIIIIIIIIIIIIIIIIII @HISEQ:573:HLWMHBCXY:1:1101:11739:2308 2:N:0:TGGACTAT GCTGTGAGGACTACTTGGGTATCTAATCCTGTTTGCTACCCACGCTTTCG + BDDDDHHIIIIIIIIIIIIIIIIIIIIIIIIIIIIHHIIIIIIIIIIIII

但是有时我们需要下面格式。

@HISEQ:573:HLWMHBCXY:1:1101:11739:2308/1 CCTATGGGGTGCAGCAGTGGGGAATATTGGACAATGGGGGGAACCCTGAT + IIIIIIIIIIIIIIIIIIIIIIHIIIIIIIIIIIIIIIIIIIIIIIIIII @HISEQ:573:HLWMHBCXY:1:1101:11739:2308/2 GCTGTGAGGACTACTTGGGTATCTAATCCTGTTTGCTACCCACGCTTTCG + BDDDDHHIIIIIIIIIIIIIIIIIIIIIIIIIIIIHHIIIIIIIIIIIII

seqtk_interleave 可以输出这样old fashion格式。

命令行接口:

$ seqtk_interleave Usage: seqtk_interleave [options] <in1.fq> <in2.fq> Options: -o old fashion interleave fastq format -v print version number […]

OTU表子集操作:otutable_subset 和 otutable_subsamples

otutable_subset 和 otutable_subsamples 程序设计的目的是解决获取OTU表子集的操作方式, 比如:

1. 提取给定OTU,并按照指定顺序生成新的OTU; 2. 提取特定样本,并按照特定顺序生成新的OTU; otutable_subset

Usearch 的 otutab_otu_subset 可以完成类似的事情, 但是不能很好解决一些问题, 比如不支持带有注释的OTU表,支持streaming 模式。

命令行接口:

$ otutable_subset Usage: otutab_subset [options] <otu-table> Options: -r reorder the OTUs by OTUs in list -v print version number

默认按照 OTU表中出现的顺序输出, 如指定 -r 可选项可以按照输入文件顺序。

实例:

文件 otu_ids.txt 为我们需要提取的OTU列表:

OTU_2 OTU_1 OTU_3

命令

otutable_subset -r otutable_counts_ann.txt otu_ids.txt otutable_subset […]

生成不同水平系统分类丰度表:otutable_abundance 和 otutable_level

这篇文章介绍 otutable_abundance 和 otutable_level, 对注释的OTU表统计不同系统分类水平丰度信息表, 支持OTU频数表和相对丰度表,注释信息按照单字母缩写模式:

比如:d:Bacteria,p:Cyanobacteria/Chloroplast,c:Chloroplast,f:Chloroplast,g:Streptophyta

1. otutable_level

命令行接口

$ otutable_level Usage: otutable_level [options] <otu-table:annotated> Options: -l taxonomy level in [dkcofgs] Default:[‘g’] -v print version number

通过 -l 指定分类水平, 默认是 genus 水平,映射表如下:

k Kingdom o Order d Domain f Family p Phylum g Genus c Class s Species

程序接口简单接受一个参数: OTU注释表。

实例

otutable_level -l g […]

OTU表生成 Krona 输入文件:otutable_krona

这篇文章介绍如何根据注释的 OTU Table 进行 Krona 可视化

1. Krona介绍

Krona是一种基于HTML5的交互式数据可视化方式,兼容大部分现代浏览器,比如 Chrome,Firefox,IE等。

参考文件:Interactive metagenomic visualization in a Web browser Github仓库地址: Krona: Interactively explore metagenomes and more from a web browser.

KronaTools 可本地安装,并包含各种程序进行导入不同工具的输出结果,比如:RDP、BLAST、MEGAN等, Krona的各种使用我们通过另外一篇文章进行描述,现在主要使用Krona的基本功能。 ktImportText

主要命令行接口:

$ ktImportText ktImportText \ [options] \ text_1[,name_1] \ [text_2[,name_2]] \ … [-o <string>] Output file name. [Default: ‘text.krona.html’] [-n <string>] Name […]

对OTU表进行注释:otutable_annotation

这篇文章介绍小工具 otutable_annotation, 这个工具设计比较通用,可以对 usearch sintax 或者 usearch utax 、RDP、 DADA2的分类结果进行OTU, 不过需要需要一些适配器才能完成工作,下面介绍如何结合其它工具完成各种OTU表注释。

1. 输入文件格式:

OTU表格式为标准QIIMME OTU 表格式, 包含以 # 开头的注释行;

#OTU ID F_5 F_6 F_7 F_8 OTU_1 1763 8408 1420 2674 OTU_2 1582 7132 1234 2317 OTU_3 0 10 1 15

注释文件需要这样格式:

OTU_1 d:Bacteria,p:Cyanobacteria/Chloroplast,c:Chloroplast,f:Chloroplast,g:Streptophyta

我们可以通过特定适配器灵活完成注释。

otutable_counts.txt OTU表; otu.sintax.txt sintax分类结果; otu.rdp.txt 分类结果; dada2.taxa.txt 分类结果;

2. 命令行接口:

$otutable_annotation […]

使用SINTAX对16S序列进行分类:SINTAX

一、SINTAX介绍

先前介绍过RDP分类器的使用:使用RDP对16S序列进行分类:rdp_classifier, 现在介绍另外一个基于非朴素贝叶斯分类算法: SINTAX(SImple Non-Bayesian TAXonomy) , 文章见:

SINTAX: a simple non-Bayesian taxonomy classifier for 16S and ITS sequences.

SINTAX使用 Kmer (默认k=8, 放回式采样32 kmers)去计算和参考库共享Kmer, 确定最佳Hits, 迭代100次, 确定每个Level的分类(出现次数最多的分类)和可信度(出现的频率)。由此可见并不需要训练参考序列库。

后面会继续介绍 DADA2 assignTaxonomy 和 assignSpecies QIIME2 的 q2-feature-classifier,新文章: Optimizing taxonomic classification of marker gene sequences

二、SINTAX 实践:

1 数据库下载

SINTAX提供了 RDP training set 16 (13k seqs, with species […]

多tsv文件组合操作: tabtk_join 和 tabtk_agg

今天介绍两个合并 tsv文件 的工具: tabtk_join 和 tabtk_agg

1. tabtk_join 组合多个文件, key设置为第一列 2. tabtk_agg 指定Key列(可多列)和目标列,将多个文件的目标列合并到同一个文件

下面以下图为示例介绍如何使用tabtk_join和tabtk_agg完成类似操作

一、tabtk_join 表合并操作

我们经典使用场景:

组合不同16S多样性实验的结果(比如不同文章的不同16S分区的物种分类信息),不推荐使用Closed-OTU方式构建新的OTU表, 最合适的方式是单独构建OTU表,然后合并不同水平的分类结果。

命令行接口:

$tabtk_join Usage: tabtk_join [options] Options: -p CHAR placehold for missing value, default: [‘-‘]; -v print version number

可选项:

-p 占位符,数据缺失时使用占位符代替, 数值型可以使用’0’, 文本型使用’-‘

示例:

$cat A-phylum.txt #level F_1 F_2 Fusobacteria 0.001872 0.000248 Firmicutes 0.85635 0.827091 […]

计算OTU表相对丰度:otutable_freq

本文介绍otutable-utils系列的一个小工具:otutable_freq (Github 链接 ), 就是计算OTU表的相对丰度(频数转频率).

1. 命令行接口 $ otutable_freq Usage: otutable_freq <matrix> version 0.0.2 2. 命令行接口 otutable_freq otutable_counts.txt >otutable_freqs.txt 3. 其它实现

3.1 otutab_counts2freqs

otutab_counts2freqs 是 Usearch 9.2 版本才加入的命令行工具,接口也比较简单

usearch -otutab_counts2freqs otutable_counts.txt -output otutable_freqs.txt

本文材料为 BASE (Biostack Applied bioinformatic SEies ) 课程 Linux Command Line Tools for Life Scientists 材料, 版权归 上海逻捷信息科技有限公司 所有.

Last Upate: 2017-10-18 […]

通过 minhash 指纹快速进行菌种鉴定:MASH 篇

前面一个文章提到了如何下载数据,如何并行下载NCBI拼装数据库的Genbank和Refseq数据, 现在开始要讲讲如何使用这些数据。

一、 mash介绍

发表在 Genome Biology 杂志 Mash: fast genome and metagenome distance estimation using MinHash。 可以通过Github Mash 仓库 下载可执行文件;

wget https://github.com/marbl/Mash/releases/download/v2.0/mash-Linux64-v2.0.tar tar xvf mash-Linux64-v2.0.tar mv mash-Linux64-v2.0 mash-2.0

MinHash介绍可以参考:Wikipedia MinHash

Mash原理图

通过构建序列(基因组,原始reads)Kmer集合,并计算出kmer 的Hash集合, 然后通过杰卡德指数Jaccard Index 估计相似度, 即 A和B共享的Hash集合(Hash值可以通过Hash函数计算)与A和B所有的不同Hash集合的比值。

另外: 传统的DDH杂交算法也可以很好的进行菌株鉴定,但是速度没有使用MinHash快,下次介绍几种DDH相关工具。

二、命令行接口

主程序为单一可执行文件程序,分发比较方便,注册主程序命令接口:

$ mash Mash version 2.0 Type ‘mash –license’ for license and […]

行计数 和 (组合)列 uniq 并统计计数工具: tabtk_nlines 和 tabtk_uniq

本文介绍两个小工具:tabtk_nlines 和 tabtk_uniq

1. tabtk_nlines

其实使用 wc -l 或者 grep -P ‘^#’ | wc -l 可以很方便的解决行计数问题,但是处理数据时候需要做两个事情:

添加样本信息 过滤掉注释行

使用 tabtk_nlines 可以干净的解决这个问题:

命令行接口:

$ tabtk_nlines Usage: tabtk_nlines [options] <text> <samples> Options: -i flag to ignore lines start with ‘#’; -v print version number

需要两个参数, 一个文本文件(支持压缩包)和样本信息

示例:

统计F1 样本 F1-ko.tsv.gz 注释到ko的基因个数,有表头信息,以#开头为注释信息

$ tabtk_nlines -i F1-ko.tsv.gz F1 F1 […]