1
$\begingroup$

I have a fasta file containing the sequence of a gene across different species. In total there are around 900 samples and 12 species. (Each sequences is over multiple lines and longer than 100bp.)

My fasta file looks like:

>Species-1-samplenameA
CTATCCTTAAACGCATATCTCGCACAGTAACTCCCCAATATGTGAGCATCTGATGTTGCCCGGGCCGAGTTAGTCTTGTGCTCACGGAACTTATTGTATG
>Species-2-samplenameB
AGTAGTGATTTGAAAGAGTTGTCAGTTAGCTCGTTCAGGTAATGGTTCCTCACACTACGTCAAAATAAGAGAGCGGTCGTGACATTATCCGTGATTTTCT
>Species-3-samplenameC
CACTACTATCAGTACTCACGACTCGATTCTGCCGCAGCCACGTATCGCCAGAAAGCCAGTCAGCATTAAGGAGTGCTCTGGGCAGGACAACTCGCATAGT
>Species-3-samplenameD
GAGAGTTACATGTTCGTTGGGCTCTTCCGACACGAACCTCAGTTGGCCTACATCCTACCTGAGGTCTGTGCCCCGGTGGTGAGAAGTGCGCATTTCGTTC

I want to split this file in one fasta file per species. I think it’s possible to use the awk function for this but I’m stuck. Does anyone have a script/code that might help me?

Thanks a lot.

Improve this question
$\endgroup$

1 Answer 1

Reset to default
3
$\begingroup$

A simple Biopython solution- iterate over the sequence records, identify the species, open the file handle using append mode to ensure no data is overwritten, and write the record:

from Bio import SeqIO

for record in SeqIO.parse("myfile.fa", "fasta"):
    species = record.id.split('_')[0]
    with open(f"{species}.fa", "a") as f:
        SeqIO.write(record, f, "fasta")
Improve this answer
$\endgroup$
5
  • $\begingroup$ Thanks. I tried: from Bio import SeqIO for record in SeqIO.parse("test.fasta", "fasta"): species = record.id.split('Species-')[0] with open(f"{species}.fasta", "a") as f: SeqIO.write(record, f, "fasta") But I got as output every sample in a separate file and not combined all species one in one file. $\endgroup$
    – Julie
    Jul 26, 2020 at 14:38
  • $\begingroup$ @Julie your code on the line species = does not match my code, try what I suggested (if your input format is as described) $\endgroup$
    – Chris_Rands
    Jul 26, 2020 at 14:41
  • $\begingroup$ Sorry. I made a mistake when posting the question, my fasta headers actually look like: >Species-1-samplenameA ... I edited my post. $\endgroup$
    – Julie
    Jul 26, 2020 at 14:44
  • $\begingroup$ @Julie species = '-'.join(record.id.split('-', 2)[:2]) $\endgroup$
    – Chris_Rands
    Jul 26, 2020 at 14:46
  • $\begingroup$ Works perfectly! Thanks a lot. $\endgroup$
    – Julie
    Jul 26, 2020 at 14:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.