2016-04-22 73 views
0

我有我想要翻譯成氨基酸的Phylip格式的DNA數據。我試過搜索可以做到這一點的庫(或模塊),但所有這些似乎都以FastA格式轉換/生成文件。將Phylip格式的DNA數據翻譯成氨基酸

這是輸入數據的外觀:

3 1500 

seq1 TTTGCTA... 

seq2 TTCGCAA... 

seq3 TTTGCCA... 

,其中1500是序列

這是代碼我有,但我得到的輸出文件的長度爲空:

#!/usr/bin/python 

import sys 

filename = '/path/to/phylip/data/' 
finalrst = open('/path/to/translated/phylip/data/','w') 


def translate_dna(sequence): 

    codontable = { 
    'ACC':'T', 'ACG':'T', 'ACT':'T', 
    'AAC':'N', 'AAT':'N', 'AAA':'K', 'AAG':'K', 
    'AGC':'S', 'AGT':'S', 'AGA':'R', 'AGG':'R', 
    'CTA':'L', 'CTC':'L', 'CTG':'L', 'CTT':'L', 
    'CCA':'P', 'CCC':'P', 'CCG':'P', 'CCT':'P', 
    'CAC':'H', 'CAT':'H', 'CAA':'Q', 'CAG':'Q', 
    'CGA':'R', 'CGC':'R', 'CGG':'R', 'CGT':'R', 
    'GTA':'V', 'GTC':'V', 'GTG':'V', 'GTT':'V', 
    'GCA':'A', 'GCC':'A', 'GCG':'A', 'GCT':'A', 
    'GAC':'D', 'GAT':'D', 'GAA':'E', 'GAG':'E', 
    'GGA':'G', 'GGC':'G', 'GGG':'G', 'GGT':'G', 
    'TCA':'S', 'TCC':'S', 'TCG':'S', 'TCT':'S', 
    'TTC':'F', 'TTT':'F', 'TTA':'L', 'TTG':'L', 
    'TAC':'Y', 'TAT':'Y', 'TAA':'*', 'TAG':'*', 
    'TGC':'C', 'TGT':'C', 'TGA':'*', 'TGG':'W', 
    'ATG':'M' 
    } 
    proteinsequence = '' 
    for n in range (0,len(sequence),3): 
      if sequence[n:n+3] in codontable: 
        proteinsequence += codontable[cds[n:n+3]] 
      sequence = '' 
    print proteinsequence 

for line in open(filename): 
    if line[0] == "3 1500": 
     finalrst.write(line) 
    elif line == '': 
     finalrst.write(line) 
    elif line.startswith('sequence'): 
      finalrst.write(line + translate_dna(line.replace('sequence', ''))) 

finalrst.close() 

有什麼問題的建議?或者更好的方式來完成這項任務?

謝謝!

回答

0

這些文件路徑是文件夾還是文件?如果data是文件,則將/關閉。如果它是一個文件夾路徑,那麼你需要指定你想打開的文件。

然後,只是理智的緣故改變這樣的:

for line in open(filename): 
    if line[0] == "3 1500": 
     finalrst.write(line) 
    elif line == '': 
     finalrst.write(line) 
    elif line.startswith('sequence'): 
     finalrst.write(line + translate_dna(line.replace('sequence', ''))) 

爲了這樣的事情:

with open(filename, 'r') as readfile: 
    for line in readfile: 
     line = line.strip() 
     # Check the full line, stripped instead of the first character. 
     if line == "3 1500": 
      finalrst.write(line + '\n') 
     elif line == '': 
      finalrst.write(line + '\n') 
     elif line.startswith('sequence'): 
      finalrst.write(line + translate_dna(line.replace('sequence', '')) + '\n') 

這樣,readfile文件句柄總是會得到關閉。

但這可能是由於文件路徑指向一個文件夾。如果不是這樣,那可能是因爲打開了延遲的文件句柄,但沒有關閉。

+0

謝謝你幫助!還有一個問題:我正在訴諸使用elif line.startswith('seq1'): \t \t \t \t finalrst.write('seq1 \ t'+ translate_dna(line.split()[ - 1])+'\ n ')..對於數據集中的所有序列,關於如何讓代碼識別序列名稱以及(1)不翻譯它並(2)在翻譯後的數據之前寫入它的任何建議? – Hia3

+0

@ Hia3,我會把它作爲一個單獨的問題來輸入。 – DuckPuncher

0

您的translate_dna似乎無法正常工作。這裏有一個工作方法,儘管不是最優化的

def translate_dna(sequence): 
    sequence = sequence.upper() 

    codontable = { 
    'ACC':'T', 'ACG':'T', 'ACT':'T', 
    'AAC':'N', 'AAT':'N', 'AAA':'K', 'AAG':'K', 
    'AGC':'S', 'AGT':'S', 'AGA':'R', 'AGG':'R', 
    'CTA':'L', 'CTC':'L', 'CTG':'L', 'CTT':'L', 
    'CCA':'P', 'CCC':'P', 'CCG':'P', 'CCT':'P', 
    'CAC':'H', 'CAT':'H', 'CAA':'Q', 'CAG':'Q', 
    'CGA':'R', 'CGC':'R', 'CGG':'R', 'CGT':'R', 
    'GTA':'V', 'GTC':'V', 'GTG':'V', 'GTT':'V', 
    'GCA':'A', 'GCC':'A', 'GCG':'A', 'GCT':'A', 
    'GAC':'D', 'GAT':'D', 'GAA':'E', 'GAG':'E', 
    'GGA':'G', 'GGC':'G', 'GGG':'G', 'GGT':'G', 
    'TCA':'S', 'TCC':'S', 'TCG':'S', 'TCT':'S', 
    'TTC':'F', 'TTT':'F', 'TTA':'L', 'TTG':'L', 
    'TAC':'Y', 'TAT':'Y', 'TAA':'*', 'TAG':'*', 
    'TGC':'C', 'TGT':'C', 'TGA':'*', 'TGG':'W', 
    'ATG':'M' 
    } 

    translated = '' 
    while len(sequence) >=3: 
     substring_3 = sequence[:3] 
     if substring_3 in codontable: 
      translated+= codontable[substring_3] 
      sequence = sequence[1:] 
     else: 
      sequence = sequence[3:] 

    return translated 

此外,還有其他問題。例如:

elif line.startswith('sequence'): 
    finalrst.write(line + translate_dna(line.replace('sequence', ''))) 

。在你輸入的字符串沒有「序」。使它:

elif line.startswith('seq'): 
    finalrst.write(line + '\t' + translate_dna(line.split()[-1])) 
+0

感謝這工作(與其他答覆一些補充),除了我必須改變序列=序列[1:]到序列=序列[3:],因爲它被讀取的方式,例如ATTGTC作爲ATT,TTG,TGT和如此,而不是ATT,GTC。 – Hia3

+0

當然。很高興工作 – Spade