寻找优雅的类通配符DNA字符串扩展方法

2 投票
5 回答
672 浏览
提问于 2025-04-15 12:46

我正在尝试对一组DNA字符串进行类似通配符的扩展,这些字符串有多个可能的碱基。

我的DNA字符串的碱基包含字母A、C、G和T。不过,我可以有一些特殊字符,比如M,它可以表示A或C。

举个例子,假设我有这个字符串:

ATMM

我想把这个字符串作为输入,然后输出四个可能的匹配字符串:

ATAA ATAC ATCA ATCC

我不想用暴力破解的方法来解决这个问题,我觉得应该有一些优雅的Python、Perl或正则表达式的技巧可以做到这一点。

谢谢你的建议。

编辑,感谢cortex提供的乘积运算符。这是我的解决方案:

我还是Python的新手,所以我觉得处理每个字典键的方式可能还有比再用一个循环更好的方法。任何建议都非常欢迎。

import sys
from itertools import product

baseDict = dict(M=['A','C'],R=['A','G'],W=['A','T'],S=['C','G'],
                  Y=['C','T'],K=['G','T'],V=['A','C','G'],
                  H=['A','C','T'],D=['A','G','T'],B=['C','G','T'])
def glob(str):
    strings = [str]

    ## this loop visits very possible base in the dictionary
    ## probably a cleaner way to do it
    for base in baseDict:
        oldstrings = strings
        strings = []
        for string in oldstrings:
            strings += map("".join,product(*[baseDict[base] if x == base 
                                 else [x] for x in string]))
    return strings

for line in sys.stdin.readlines():
    line = line.rstrip('\n')
    permutations = glob(line)
    for x in permutations:
        print x

5 个回答

0

这其实不是一个“扩展”的问题,而且几乎肯定用常规的正则表达式是无法解决的。

我觉得你想要的其实是“如何生成排列”。

1

你可以在Python中使用yield这个操作符来做类似的事情。

def glob(str):
      if str=='':           
          yield ''
          return      

      if str[0]!='M':
          for tail in glob(str[1:]): 
              yield str[0] + tail                  
      else:
         for c in ['A','G','C','T']:
             for tail in glob(str[1:]):
                 yield c + tail                 
      return

编辑:正如大家指出的,我之前犯了一些错误。这里是我尝试过并且有效的版本。

2

我同意其他人的看法,这听起来确实是个奇怪的想法。当然,如果你真的想这么做,Python(2.6及以上版本)总是有一种优雅的方法来实现:

from itertools import product
map("".join, product(*[['A', 'C'] if x == "M" else [x] for x in "GMTTMCA"]))

这是一个完整的解决方案,包括输入处理:

import sys
from itertools import product

base_globs = {"M":['A','C'], "R":['A','G'], "W":['A','T'],
              "S":['C','G'], "Y":['C','T'], "K":['G','T'],

              "V":['A','C','G'], "H":['A','C','T'],
              "D":['A','G','T'], "B":['C','G','T'],
              }

def base_glob(glob_sequence):
    production_sequence = [base_globs.get(base, [base]) for base in glob_sequence]
    return map("".join, product(*production_sequence))

for line in sys.stdin.readlines():
    productions = base_glob(line.strip())
    print "\n".join(productions)

撰写回答