multiFASTA Dateiverarbeitung

https://stackoverflow.com/questions/1789307

22-09-2019
|

Frage

Ich war neugierig zu wissen, ob es ein Bioinformatik-Tool ist da draußen die Lage, eine multiFASTA Datei zu verarbeiten mir Infos wie die Anzahl der Sequenzen geben, Länge, Nukleotid / Aminosäuregehalt usw. und vielleicht automatisch beschreibendes Diagramme zu zeichnen. Auch eine R Bioconductor Lösung oder ein BioPerl Modul tun würde, aber ich habe es nicht schaffen, etwas zu finden.

Können Sie mir helfen? Vielen Dank: -)

Lösung

Einige der Prägungs Tools sind eine Sammlung von kleinen Tools, die Ihnen helfen können.

seqstats kehrt Sequenzlänge
pepstats sollten Sie Inhalt usw. Aminosäure Einige der Werkzeuge bieten auch Funktionen Plotten. Sehr praktisch. http://emboss.sourceforge.net/apps/release/5.0 /emboss/apps/groups.html

Zur Anzahl der fasta Einträge zu zählen, die ich benutze: grep -c '^>' mySequences.fasta.

Um sicherzustellen, dass keiner der Einträge werden doppelte, überprüfe ich, dass ich die gleiche Nummer erhalten dabei ist: grep '^>' mySequences.fasta | sort | uniq | wc -l

Andere Tipps

Sie können auch in faSize interessiert sein, die ein Werkzeug aus der

Es sollte beachtet werden (für jedermann auf diese Stolpern, wie ich gerade), dass es eine robuste Python-Bibliothek ist speziell für diese Aufgaben genannt zu behandeln entworfen biopython . In wenigen Zeilen Code, können Sie schnell Zugriff Antworten für alle oben genannten Fragen. Hier sind einige sehr einfache Beispiele, vor allem aus der Verbindung angepasst. Es gibt Kessel-Platte GC% Diagramme und Sequenzlänge Graphen im Tutorial auch.

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(), '*'))

Lizenziert unter: CC-BY-SA mit Zuschreibung

Nicht verbunden mit StackOverflow