DNA子序列动态规划问题

Question

我正在尝试解决 DNA 问题，这更像是 LCS 问题的改进（？）版本。 在这个问题中，有一个字符串是字符串和半子字符串，它允许部分字符串跳过一个或不跳过一个字母。例如，对于字符串“desktop”，它有半子字符串{"destop", "dek", "stop", "skop","desk","top"}，所有这些都跳过一个或不跳过一个字母。

现在，我得到了两个由{a,t,g,c} 组成的 DNA 字符串。我正在尝试查找最长的半子字符串 LSS。如果有多个 LSS，则以最快的顺序打印出一个。

例如，两个 dnas {attgcgtagcaatg, tctcaggtcgatagtgac} 打印出 "tctagcaatg"

aaaattttcccc, cccgggggaatatca 打印出"aattc"

我正在尝试使用常见的 LCS 算法，但无法用表格解决它，尽管我确实解决了没有跳过字母的问题。有什么建议吗？

Answer 1

这是 LCS 动态编程解决方案的一种变体，用 Python 编写。

首先，我正在为可以使用跳过规则从每个字符串生成的所有子字符串构建一个Suffix Tree。然后我与后缀树相交。然后我正在寻找可以从该交叉树中获得的最长字符串。

请注意，这在技术上是 O(n^2)。最坏的情况是两个字符串都是同一个字符，一遍又一遍地重复。因为你最终会得到很多逻辑上类似的东西，“一个字符串中第 42 位的 'l' 可能与另一个字符串中第 54 位的位置 l 匹配”。但实际上它将是O(n)。

def find_subtree (text, max_skip=1):
    tree = {}
    tree_at_position = {}

    def subtree_from_position (position):
        if position not in tree_at_position:
            this_tree = {}
            if position < len(text):
                char = text[position]
                # Make sure that we've populated the further tree.
                subtree_from_position(position + 1)

                # If this char appeared later, include those possible matches.
                if char in tree:
                    for char2, subtree in tree[char].iteritems():
                        this_tree[char2] = subtree

                # And now update the new choices.
                for skip in range(max_skip + 1, 0, -1):
                    if position + skip < len(text):
                        this_tree[text[position + skip]] = subtree_from_position(position + skip)

                tree[char] = this_tree

            tree_at_position[position] = this_tree

        return tree_at_position[position]

    subtree_from_position(0)

    return tree


def find_longest_common_semistring (text1, text2):
    tree1 = find_subtree(text1)
    tree2 = find_subtree(text2)

    answered = {}
    def find_intersection (subtree1, subtree2):
        unique = (id(subtree1), id(subtree2))
        if unique not in answered:
            answer = {}
            for k, v in subtree1.iteritems():
                if k in subtree2:
                    answer[k] = find_intersection(v, subtree2[k])
            answered[unique] = answer
        return answered[unique]


    found_longest = {}
    def find_longest (tree):
        if id(tree) not in found_longest:
            best_candidate = ''
            for char, subtree in tree.iteritems():
                candidate = char + find_longest(subtree)
                if len(best_candidate) < len(candidate):
                    best_candidate = candidate
            found_longest[id(tree)] = best_candidate
        return found_longest[id(tree)]

    intersection_tree = find_intersection(tree1, tree2)
    return find_longest(intersection_tree)


print(find_longest_common_semistring("attgcgtagcaatg", "tctcaggtcgatagtgac"))

Answer 2

令g(c, rs, rt)代表字符串中最长的公共半子串S和T，以rs和rt结尾，其中rs和rt是该字符的排序出现次数，c，分别在S 和T 中，K 是允许的跳过次数。然后我们可以形成一个递归，我们必须对 S 和 T 中的所有 c 对执行。

JavaScript 代码：

function f(S, T, K){
  // mapS maps a char to indexes of its occurrences in S
  // rsS maps the index in S to that char's rank (index) in mapS
  const [mapS, rsS] = mapString(S)
  const [mapT, rsT] = mapString(T)
  // h is used to memoize g
  const h = {}

  function g(c, rs, rt){
    if (rs < 0 || rt < 0)
      return 0
    
    if (h.hasOwnProperty([c, rs, rt]))
      return h[[c, rs, rt]]
      
    // (We are guaranteed to be on
    // a match in this state.)
    let best = [1, c]
    
    let idxS = mapS[c][rs]
    let idxT = mapT[c][rt]

    if (idxS == 0 || idxT == 0)
      return best

    for (let i=idxS-1; i>=Math.max(0, idxS - 1 - K); i--){
      for (let j=idxT-1; j>=Math.max(0, idxT - 1 - K); j--){
        if (S[i] == T[j]){
          const [len, str] = g(S[i], rsS[i], rsT[j])

          if (len + 1 >= best[0])
            best = [len + 1, str + c]
        }
      }
    }
  
    return h[[c, rs, rt]] = best
  }

  let best = [0, '']
  
  for (let c of Object.keys(mapS)){
    for (let i=0; i<(mapS[c]||[]).length; i++){
      for (let j=0; j<(mapT[c]||[]).length; j++){
        let [len, str] = g(c, i, j)
        
        if (len > best[0])
          best = [len, str]
      }
    }
  }
  
  return best
}

function mapString(s){
  let map = {}
  let rs = []
  
  for (let i=0; i<s.length; i++){
    if (!map[s[i]]){
      map[s[i]] = [i]
      rs.push(0)
    
    } else {
      map[s[i]].push(i)
      rs.push(map[s[i]].length - 1)
    }
  }
  
  return [map, rs]
}

console.log(f('attgcgtagcaatg', 'tctcaggtcgatagtgac', 1))

console.log(f('aaaattttcccc', 'cccgggggaatatca', 1))

console.log(f('abcade', 'axe', 1))