mothur summary.seqs 统计fasta文件中每条序列的长度

Posted 2020-09-30 庐州月光

在介绍summary.seqs的用法之前，我们首先需要搞清楚两个概念：

1）ambiguous bases

中文叫做模糊碱基，对于DNA序列来说，只有ATCG 4种碱基，在IUPAC定义的碱基标准中，出了上述4种碱基之外，还包括其他的碱基，可以代表不同类型的碱基

代码		英文含义	中文含义
G		Guanine	鸟嘌啉
A		Adenine	腺嘌啉
T	(U)	Thymine (Uracil)	胸腺嘧啶 (尿嘧啶)
C		Cytosine	胞嘧啶
R	(A or G)	PuRine	嘌啉
Y	(C or T or U)	Pyrimidine	嘧啶
M	(A or C)	Amino	腺嘌啉或胞嘧啶(氨基)
K	(G or T)	Ketone	鸟嘌啉或胸腺嘧啶(酮基)
S	(C or G)	Strong interaction	强相互作用碱基
W	(A or T)	Weak interaction	弱相互作用碱基
H	(A or C or T)	Not-G (H after G)	非鸟嘌啉
B	(C or G or T)	Not-A (B after A)	非腺嘌啉
V	(A or C or G)	Not-T/U (V after U)	非胸腺嘧啶
D	(A or G or T)	Not-C (D after C)	非胞嘧啶
N	(A or C or G or T)	Any	不确定

模糊碱基实际上就是除了A T C G 这4种碱基之外的其他碱基

2）homopolymer base

由1个碱基重复多次的序列，比如GCAGAAAAAAA 序列中，末端的一串A就是 homopolymer base

summary.seqs的基本用法：

mothur "#summary.seqs(fasta = "input.fasta")"

 

运行成功之后，会生成input.summary 文件，内容如下：

seqname	start	end	nbases	ambigs	polymer	numSeqs
1	1	24	24	0	2	1
2	1	25	25	10	10	1
3	1	25	25	2	1	1
4	1	24	24	0	18	1
5	1	24	24	0	2	1
6	1	24	24	0	1	1
7	1	24	24	0	1	1
8	1	25	25	0	2	1

 

共7列，每列表头含义如下：

seqname : 序列标识符

start      :  起始位置，从1开始

end       :   终止位置，

nbases  :  总碱基数，可以看做序列长度

ambigs  : ambiguous bases 模糊碱基的数目

polymer :  homopolymer 碱基的最大长度

numSeqs : 序列数，对于每条序列来说，其值总是为1

除了上述的基本用法外，summary.seqs 还有很多的参数；

processors  : CPU个数，mothur 是支持并行的，通过设置processors 参数可以并行执行程序，用法如下：

mothur "#summary.seqs(fasta = "input.fasta"， processors = 10)"