将DNA链转换为正则表达式

Question

以下是我目前掌握的代码：

import re

with open("link_allenz.txt") as f:
    s = f.read()


baseDict={'A':'A','C': 'C','T':'T','G':'G','R':'[GA]','Y':'[CT]',
          'M':'[AC]','K':'[GT]','S':'[GC]','W':'[AT]',
          'B':'[CGT]','D':'[ACT]','V':'[ACG]','N':'[ACGT]',',':'|'}

pat1=r'<1>(.*)\n(?:.*\n){3}<5>(.*)\n'
re1=re.compile(pat1)
e2r_dict=dict(re1.findall(s))

def make_seq_2_re(strand):
    while True:
        try:
            x=e2r_dict[strand]
            x_split=list(x)
            print(x_split)
            for base in x_split:
                if baseDict.get(base):
                    x_split[x_split.index(base)]=baseDict[base]
                else:
                    x_split.remove(base)
            y=''.join(x_split)
            if len(y)>0:
                return(y)
            else:
                print("No known strand for "+strand+"!")
            break
        except KeyError:
            print("Was there a spelling error?")
            pass

对于我的类项目，这应该是读取一个充满DNA链的文件，并将我告诉它的转换成正则表达式。出于某种原因，它并不总是转换它需要转换的每个字母。例如，如果它应该转换GACNNNN^NNGTC，它将输出GAC[ACGT][ACGT][ACGT][ACGT][ACGT][ACGT][ACGT]NGTC，忽略最后的N。我如何更改它

Answer 1

您在对x_split进行迭代时正在修改它。这根本不起作用。相反，解决方案是在迭代输入时生成一个新序列

此代码还存在多个其他问题

1. while True:做什么？这似乎是不必要的，而且您总是在一次迭代后中断。您可以完全删除它

2. 你的try块意味着如果用户用错误的键调用函数，他们会得到一条消息而不是一个错误。创建一条可读的消息是可以的，但接受错误是而不是。改为引发自定义错误

3. try块太大：它应该只包含e2_dict子集，而不是整个函数。将其缩小可以使错误发生的位置更清晰，并使代码更具可读性

4. x_split=list(x)是不必要的：字符串是序列，您可以像遍历列表一样遍历它们

5. x_split.index(base)做了错误的事情：它返回了x_split中base的第一个索引。考虑一下当你有了输入{{CD10}}，并且你在最后一个基点时发生了什么：^ {< CD11}}将返回1，而不是2（正确的索引）。您不能在此处使用index。如果在序列上迭代时需要索引，请使用for i, base in enumerate(x_split)。但是，在更正的代码中，我们实际上并不需要它，请参见下文

6. 组合baseDict.get(base)…baseDict[base]执行一次冗余查找。理想情况下，我们希望避免这种不必要的工作

7. if len(y)>0:有点奇怪：只有当e2r_dict[strand]存在但为空或仅包含无效符号（例如^^）时，此条件才会是False。这是一个实际发生的案例吗？如果是这样，则应改进错误消息。同样，不要打印错误消息引发错误。如果我们反转条件，我们也不需要else情况

8. 名字strand很奇怪：它似乎指的是DNA以外的东西；我不知道你的确切域，所以也许这是正确的术语，但我从来没有听说过它在序列分析中指的是DNA链以外的任何东西

9. 类似地，像pat1、re1、baseDict、x和y这样的名字也不是好名字，它们要么根本不提供任何关于它们的用途的信息（x、y），要么提供不相关的信息（在re1中是什么？{对于来说是什么？）

   选择好的名字很重要，因为它大大提高了可读性make_seq_2_re更好，但仍然可以改进：make通常过于泛化，在名称中不太有用（尽管也有例外），而且此类名称通常被表述为“从B生成A”，而不是“从B生成A”，以符合英语语法

10. 熟悉the official Python style guide, PEP 8。特别是，更喜欢snake_case而不是mixedCase，并在操作符周围放置一致的空格

考虑到所有这些因素，我将编写如下函数，例如：

def seq_re(name):
    try:
        seq = e2r_dict[name]
    except KeyError:
        raise KeyError(f'No sequence found for {name}')

    re = ''
    for base in seq:
        pattern = base_pattern_table.get(base)
        if pattern:
            re += pattern

    if not re:
        raise RuntimeError(f'sequence for {name} was invalid')

    return re

然而，Pythonists的爱列表理解和从空列表迭代建立序列被视为一种反模式。所以我很想重写这一部分：

def seq_re(name):
    try:
        seq = e2r_dict[name]
    except KeyError:
        raise KeyError(f'No sequence found for {name}')

    re = ''.join([pattern for base in seq
        if (pattern := base_pattern_table.get(base))])

    if not re:
        raise RuntimeError(f'sequence for {name} was invalid')

    return re

请注意，由于这使用了walrus operator，因此需要最新版本的Python(≥ 3.8). 我通常敦促大家注意不要过度使用这个操作符。但我认为在这种情况下是合适的