Hi everyone, I recently started using python with biopython. I'm trying to practice to get the translate ORF using this gene taken from Genbank as input: NM_100684.3
However, my output does not show me the correct ORF and I get a different amino acid sequence both in composition and length.
What am I doing wrong?
These are the scripts used by me
>>>from Bio import SeqIO
>>>record = SeqIO.read("sequence.fasta", "fasta")
>>> table = 1
>>> min_pro_len = 100
>>>for strand, nuc in [(+1, record.seq), (-1, record.seq.reverse_complement())]:
for frame in range(3):
length = 3 * ((len(record)-frame) // 3) #Multiple of three
for pro in nuc[frame:frame+length].translate(table).split("*"):
if len(pro) >= min_pro_len:
print("%s...%s - length %i, strand %i, frame %i" \
% (pro[:30], pro[-3:], len(pro), strand, frame))
YSDIDQINLNQISNLQRNLKYFITMGDSTG...NNV - length 554, strand 1, frame 2
SSPGDKGHNCKGGSASSLCPHREEHHSHNG...ILT - length 162, strand -1, frame 1
IEHQDSHDDVQPTGYKEGDPPGREGCGTAA...HNW - length 216, strand -1, frame 1
TKVTGNVQATIITPIHVSPCSVVKCEVEKK...SDA - length 122, strand -1, frame 2
This above is my output but isn't corrected and do not start with methionine, in Genbank the correct protein has 530 a.a. and start with "MGDSTGEPGSSMHGVTGREQ ..."
2 answers
As far as I can see your code does not have a check for the starting methionine, you are retrieving ORFs that generate at least 100 amino acid long peptides.
And I am pretty sure the 530 a.a. long "Genbank peptide" is within your 554 a.a. long peptide prediction. Check from 25th amino acid on in your first result, that goes like MGDSTG...
I was tryng something similar, try this:
from Bio import SeqIO
record = SeqIO.read("NC_005816.1.fna", "fasta")
table = 11
min_pro_len = 50
x= 0
for strand, nuc in [(+1, record.seq), (-1, record.seq.reverse_complement())]:
for frame in range(3):
length = 3 * ((len(record)-frame) // 3) #Multiple of three
for pro in nuc[frame:frame+length].translate(table).split("*"):
splitlocal = pro.find('M')
seq_final = pro[splitlocal:]
if len(seq_final) >= min_pro_len:
print("%s...%s - length %i, strand %i, frame %i" \
% (seq_final[:10], pro[-3:], len(seq_final), strand, frame))
x = x+1
print("Numero de ORFs:",x)
MVAHRFTCSL...MTI - length 54, strand 1, frame 0
MLKQPTATVC...KIW - length 58, strand 1, frame 0
MLWMPTSRRP...SGF - length 66, strand 1, frame 0
MPHRCVRRTC...KNM - length 56, strand 1, frame 0
MKKSSIVATI...YRF - length 312, strand 1, frame 0
MMELQHQRLM...NPE - length 259, strand 1, frame 1
Log in to answer this question.
From the docs you're following:
As it stands, all you've assessed is regions uninterrupted by stop codons, you haven't gone to the next step of identifying starts etc.