I need to convert contigs into their respective protein sequences given a reference genome (i.e. I need to take a substring, whose position is already known on the string, and I need to locate the nearest start and stop codons).
This is tricky because sometimes the first position of contig will not be a multiple of three (i.e. the first three nucleotides in the contig may not perfectly match a codon). Also, sometimes the contigs could be located in the intergenic region (i.e. between genes). The goal is to separate both coding and non-coding DNA.
This is my code so far:
from Bio.Seq import Seq from Bio.Alphabet import generic_dna, generic_protein start_codons = ['ATG'] stop_codons = ['TAG', 'TAA', 'TGA'] string = 'GG*TAG*CCAATT*ATG*AACGAA*TAG*GAC' #remove '*', just for visual contigs = ['CCAA', 'TGAAC', 'GAA', 'GGAC'] positions = [5, 12, 17, 22] #position indices for each contig on string extended_contigs =  extended_position_contigs =  intergenic_contigs =  intergenic_position_contigs =  for i in contigs: extended_contigs.append(#some code) extended_positions_contigs.append(#some code) intergenic_contigs.append(#some code) intergenic_positions_contigs.append(#some code)
I should get
extended_contigs = ['ATGAAC', 'ATGAACGAA'] and
extended_positions_contigs = [12, 17]. These are the contigs that are located within genes. In order to code them into peptides, I need to reach back in the string until I find the start codon and extend the initial contig (e.g.
I should also get
intergenic_contigs = ['CCAA', 'GGAC'] and
intergenic_positions_contigs = [5, 22]. When the missing code is run, the computer searches to the left of the string and finds a stop codon (e.g.
TAG) before the start codon. Thus, the contig is located in between two genes and nothing needs to be added. These intergenic contigs are just stored in a new list.
My code continues:
prot_contigs =  for i in extended_contigs: my_dna = Seq(i, generic_dna) my_prot = my_dna.translate() prot_contigs.append(str(my_prot))
Here, no new code needs to be added. After the above is run,
prot_contigs = ['MN', 'MNE'].
The last step of code (which I need help with) will transform
new_prot_contigs = ['MN', 'E'].
How? If for any contig (e.g.
'TGAAC'), the beginning or end is part of another codon (not perfect multiple of 3), the extra codon on either end will be maintained (e.g.
'MN'). Otherwise, if the contig (e.g.
'GAA') perfectly matches a codon, anything that was added onto it will be removed (e.g.
I would try to solve the 2 portions of code myself, but I am unsure on how to take a position on a string (i.e. beginning points of contig) and look along the string to find the nearest start/stop codons, so I could determine the DNA's function and accurately sequence the protein coding contigs into peptides.
Any help would be appreciated!