7

我很想知道是否有任何生物信息学工具能够处理 multiFASTA 文件,为我提供序列数量、长度、核苷酸/氨基酸含量等信息,并可能自动绘制描述图。也可以使用 R BIOconductor 解决方案或 BioPerl 模块,但我没有找到任何东西。

你能帮助我吗?非常感谢 :-)

4

4 回答 4

7

一些浮雕工具是可以帮助您的小工具的集合。

要计算 fasta 条目的数量,我使用: grep -c '^>' mySequences.fasta.

为确保所有条目都不重复,我在执行此操作时检查是否得到相同的数字:grep '^>' mySequences.fasta | sort | uniq | wc -l

于 2009-11-24T14:53:14.240 回答
2

您可能还对faSize感兴趣,它是Kent Source Tree中的一个工具,尽管这比仅使用 grep 需要更多的努力(您必须 dload 和编译)...这里是一些示例输出:

me@my-lab ~/data $ time faSize myfile.fna
215400419 bases (104761 N's 215295658 real 215295658 upper 0 lower) in 731620 sequences in 1 files
Total size: mean 294.4 sd 138.5 min 30 (F5854LK02GG895) max 1623 (F5854LK01AHBEH) median 307
N count: mean 0.1 sd 0.4
U count: mean 294.3 sd 138.5
L count: mean 0.0 sd 0.0
%0.00 masked total, %0.00 masked real

real    0m3.710s
user    0m3.541s
sys     0m0.164s
于 2010-01-16T17:42:55.497 回答
0

应该注意(对于任何偶然发现这一点的人,就像我刚才所做的那样),有一个强大的 python 库专门设计用于处理这些称为Biopython的任务。只需几行代码,您就可以快速获得上述所有问题的答案。以下是一些非常基本的示例,主要改编自链接。教程中还有样板 GC% 图和序列长度图。

In [1]: from Bio import SeqIO

In [2]: allSeqs = [seq_record for seq_record in SeqIO.parse('/home/kevin/stack/ls_orchid.fasta', """fasta""")]

In [3]: allSeqs[0]
Out[3]: SeqRecord(seq=Seq('CGTAACAAGGTTTCCGTAGGTGAACCTGCGGAAGGATCATTGATGAGACCGTGG...CGC', SingleLetterAlphabet()), id='gi|2765658|emb|Z78533.1|CIZ78533', name='gi|2765658|emb|Z78533.1|CIZ78533', description='gi|2765658|emb|Z78533.1|CIZ78533 C.irapeanum 5.8S rRNA gene and ITS1 and ITS2 DNA', dbxrefs=[])

In [4]: len(allSeqs) #number of unique sequences in the file
Out[4]: 94

In [5]: len(allSeqs[0].seq) # call len() on each SeqRecord.seq object
Out[5]: 740

In [6]: A_count = allSeqs[0].seq.count('A')
        C_count = allSeqs[0].seq.count('C')
        G_count = allSeqs[0].seq.count('G')
        T_count = allSeqs[0].seq.count('T')

        ​print A_count # number of A's

        144

In [7]: allSeqs[0].seq.count("AUG") # or count how many start codons
Out[7]: 0

In [8]: allSeqs[0].seq.translate() # translate DNA -> Amino Acid
Out[8]: Seq('RNKVSVGEPAEGSLMRPWNKRSSESGGPVYSAHRGHCSRGDPDLLLGRLGSVHG...*VY', HasStopCodon(ExtendedIUPACProtein(), '*'))
于 2016-01-27T19:30:12.273 回答
0

python 中的熨平板很棒。

import screed
for record in screed.open(fastafilename):
    print(len(record.sequence))

https://screed.readthedocs.io/en/v1.0/screed.html

于 2021-10-04T07:52:46.480 回答