# What is the full script to extract sequences from fasta file by using Ids in text file in python 3 and pycharm?

I have a fasta file (gene.fasta)with sequences with names patterns:

>Hsc_gene_9627.t1
ATGGCACGCATTTTCCTCATTCTTTTATTGCTGCACAACATCTGCTGTGCAGCCGCTTCA
TTGCTCATTTTGAATGCCGTTACATTGGAGAAGGATGCTAATGATTATGCCGTTGGCGAT


and I have the Ids in text file (seq.txt) which are not the exact as in the fasta file:

HSC_gene_996
HSC_gene_9734
and some of the names came as
HSC_gene_996|HSC_gene_9734


How can extract the sequences? I am not experienced in python so please use python for dummies language :) Thank you in advance

• Does it have to be python? There are far easier ways using only bash, etc. Jun 3, 2020 at 11:47
• I just have python in my computer, if there is please help me, I tried to find some online tools but it didn't work Jun 3, 2020 at 12:37
• You may want to try Joe Healey's bioinfo-tools repository. This script looks like the one to use. Jun 3, 2020 at 14:07
• To second @RamRS, that script should work indeed - provided your IDs are exact full or substring matches (excluding the >). It's been a while since I used it, so if it could be tweaked to be more useful, I'm happy to take suggestions. Whilst I love bash etc., I'd probably disagree with NatWH that it is easier than using something like biopython which is purpose built for this. Jun 3, 2020 at 15:08
• OP's problem seems to be their platform. Python excels in being a cross-platform solution, hence OP looking for Python based solutions makes sense. Jun 3, 2020 at 15:24

## 1 Answer

For example, the fasta are like these:

gene.fasta:

>Hsc_gene_9627.t1
ATGGCACGCATTTTCCTCATTCTTTTATTGCTGCACA
ACATCTGCTGTGCAGCCGCTTCA
>Hsc_gene_9627.t2
GGGGTTTTCCCC
>Hsc_gene_962.t1
AGTCGTCAGTCAGTAGTCGC

seq.txt:

HSC_gene_9627
HSC_gene_9734
HSC_gene_996|HSC_gene_9734


We use a few modules, read in the fasta:

import re
from Bio import SeqIO
records = list(SeqIO.parse("gene.fasta","fasta"))


Read in the list:

seq = [i.rstrip().lower() for i in open('seq.txt').readlines()]


Here I use two functions, one to write fasta files, and the other to format the header for comparison:

def fasta_out(rec):
return(">"+str(rec.description)+"\n"+str(rec.seq))

def format(rec):
return re.sub("[.][^ ]*\$","",rec.description).lower()


Define the records to keep:

keep = [fasta_out(i) for i in records if format(i) in seq]


Write them out:

f=open('matched.fasta','w')
f.writelines("\n".join(keep))
f.close()


We can check:

open('matched.fasta').readlines()
Out[50]:
['>Hsc_gene_9627.t1\n',
'ATGGCACGCATTTTCCTCATTCTTTTATTGCTGCACAACATCTGCTGTGCAGCCGCTTCA\n',
'>Hsc_gene_9627.t2\n',
'GGGGTTTTCCCC']


Most likely there are some simple while, for loops that does not require modules..

• Thank you very much @StupidWolf Jun 10, 2020 at 12:06