将 DNA 数据转换为 phylip 格式的氨基酸

Question

我有 Phylip 格式的 DNA 数据，我想将其转化为氨基酸。我尝试搜索可以执行此操作的库（或模块），但所有这些似乎都以 FastA 格式翻译/生成文件。

这是输入数据的样子：

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()

关于问题所在有什么建议吗？或者也许是完成这项任务的更好方法？

谢谢！

Answer 1

这些文件路径是文件夹还是文件？如果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 文件句柄将始终关闭。

但这可能是因为文件路径指向一个文件夹。如果不是这样，可能是因为挥之不去的文件句柄打开了，但没有关闭。

Answer 2

您的 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]))