用Python代码查找带有起始和终止密码子的DNA

Question

所以最近我一直在尝试编写一个程序，根据起始密码子和终止密码子来检测和切断DNA序列的编码部分。在

最终的目标是比较两个240个核苷酸的序列，但是其中一个会导致镰状细胞疾病，所以你想看看这两个结果的编码部分之间的差异。在

这是我到目前为止编写的代码，它实际上可以处理随附的序列。在

sequence = "CCATGCTTGATCA"
sequence_list = list(sequence)
codon_list = ["ATG", "TAA", "TAG", "TGA"]
position_list = []

length_sequence = len(sequence)
length_codon = len(codon_list)
length_position = len(position_list)

n = length_sequence-1
while n >= 0:

    for i in range(0, len(codon_list)):
        codon_sub_list = list(codon_list[i])

        if sequence_list[n] == codon_sub_list[2] and sequence_list[n-1] == codon_sub_list[1] and sequence_list[n-2] == codon_sub_list[0]:
            position_list.append(n-2)
            print(sequence_list[n], "@", n)
            print(sequence_list[n-1], "@", n-1)
            print(sequence_list[n-2], "@", n-2)  

    n-=1

print(len(position_list))
print(sequence[position_list[length_position-1]:(position_list[0]+3)])

现在，我两天前做这个实验的结果很有希望。正如预期的那样，结果如下：

A at location 9

G at location 8

T at location 7

G at location 4

T at location 3

A at location 2

[7,2]

ATGCTTGA

然而，今天我尝试了一个不同的序列来继续这项工作，这次是240个核苷酸长序列中的一个。 下面是这两个序列以及它是哪一个。在

镰状细胞病序列号：gaccatcattgcttacattgtctgagaactgttcactgcaacctcagcaagacaggtgctctcagcaagagagagagagagtgtgttggttgctctgtgtgaggaagtgaggatgaaggtggtgagtgtgtgtgtgtgtgttagttggttagttagagtagagagagagagagttagttagagagagagagagttggttggttagagagagagagagagttggtttttcagagagagagagagagagttggtttttct在

正常序列：gaccatcattgcttacattgtctgagagagagcttgttggttgctctcagcaagagagtctgttagtctgttagagagagtctgttagtcttaggcagtgtgagcaagagtgagatgagtgaaggtgggtggtgagtgtgtgtgttagttagtagagagagagagagttagttagagagagagagagttggttagagagagagagagttggtttcttagagagagagagagagtgtttttct

不过，这是我从执行它得到的结果，我会很快列出核苷酸和它们的位置，因为它们大部分都是无关的，而最重要的是最后一个。在

[G,211] [T,210] [A,209] [G,199] [A,198] [T,197] [A,187] [A,186] [T,185] [A,145] [G,144] [T,143] [A,133] [G,132] [T,131] [G,132] [T,131] [A,130] [A,123] [G,122] [T,121] [A,78] [G,77] [T,76] [G,68] [T,67] [A,66] [G,47] [A,46] [T,45] [A,29] [G,28] [T,27] [A,1] [G,0] [T,-1]

[209, 197, 185, 143, 131, 130, 121, 76, 66, 45, 27, -1]

No sequence, just an empty line

现在很明显，当最后一个密码子，TGA，它注意到T的位置为-1时，问题就出现了，但是我不知道是什么原因导致了这个问题，并尝试调整了几个值来使它工作，但在任何情况下都没有做到。在

我想知道是什么引起的，该怎么办？另外，我两天前做了这个，主要是作为一个初稿开始的，而且可能还有其他更好的东西，所以如果有什么东西看起来有点草率，那就找借口吧，在我看来整个while循环可以做得更好，但在那一刻，我选择了它，因为一个不同的循环方法不起作用，已经记不清到底是什么了。在

注：我制作了一个空闲输出的屏幕截图，让您了解：

Answer 1

我发现你的原始代码有点不清楚（从右边开始搜索对我来说是违反直觉的），所以我试图编写一个替代方案。最重要的变化是，我现在从左到右遍历序列，并通过一次比对子序列来搜索密码子，而不是逐个核苷酸。这是我的代码，希望有一些有用的注释。这能满足你的需要吗？如果没有，请告诉我。在

sequence = "GAGCCATCTATTGCTTACATTTGCTTCTGACACAACTGTGTTCACTAGCAACCTCAAACAGACACCATGGTGCACCTGACTCCTGTGGAGAAGTCTGCCGTTACTGCCCTGTGGGGCAAGGTGAACGTGGATGAAGTTGGTGGTGAGGCCCTGGGCAGGTTGGTATCAAGGTTACAAGACAGGTTTAAGGAGACCAATAGAAACTGGGCATGTGGAGACAGAGAAGACTCTTGGGTTTCT"
codon_list = ["ATG", "TAA", "TAG", "TGA"]

# store the starting positions of the codons
found_codon_positions = []

# note that we can use len() and [] with strings as well, no need to
# convert to list first
n = len(sequence)
k = 0
while k < n-2:
    # extract a three-nucleotide subsequence
    possible_codon = sequence[k:k+3]
    if possible_codon in codon_list:
        found_codon_positions.append(k)
    k += 1

print('found codons at indices {}'.format(found_codon_positions))

print('extracted sequence:')
print(sequence[found_codon_positions[0]:found_codon_positions[-1]+3])

输出：

^{pr2}$

Answer 2

我不确定我是否遵循您的所有逻辑，但如果您想存储所有索引并找到开始的第一个和最后一个匹配的子序列：

def find_seq(s, cd):
    od = dict((s, []) for s in codon_list)
    mn, mx = None, None
    for n in range(len(s) - 1):
        seq = sequence[n:n + 3]
        if seq in od:
            od[seq].append((seq, n))
            if  mn is None:
                mn = n
            mx = n + 3
    return mn, mx, od

dict将包含找到的所有子序列以及每个子序列开始位置的索引：

^{pr2}$