1

大家好，我有一个生物信息学问题，我可以帮忙解决。它很长，但我会尝试将其分解为更小的部分，任何帮助都很棒。

我有一个由 4 个字母 A、U、C、G 组成的 RNA 长度“n”序列，它作为字符串导入 Python，可以折叠形成一个循环。一个循环是通过匹配序列中的一对字母来形成的，这样 A 与 U，C 与 G 和 G 与 U 一起，因此字符串自身折叠起来。

问题是必须有三个或更多的字母彼此相邻形成一对，超过或等于 3 个字母连续形成一对，并且至少 3 个字母的部分之间也必须有间隙.

我试图发布一张图片，但我没有足够的声望点:(

在我引用的期刊中，作者谈到了一种嵌套循环方法，可以在可能的情况下找到所有可能的组合，然后将它们包含在一个组中，以便稍后调用。

我的问题是编写嵌套循环，因为我是编程和 python 的新手。以及以可以识别对并可能将它们添加在一起的方式存储序列。

再次任何帮助都会很棒，如果有任何不清楚的地方，请告诉我

编辑：

一个例子是 seq='aggcuugaguuu' ，其中一个输出显示 seq[0:2] 与 seq[9:11] 的配对，这意味着代码形式像 U 形。

如果你把绳子想象成一根物理的绳子，在 3 个点上握住它，在三个不同的点上握住它，然后将这些点接触在一起，它会导致绳子形成一个循环。我正在寻找确定使用的 6 点。

我不是在寻找为我编写的代码，我只是想知道一种编写代码的方法。

我尝试了一种方法，其中 seq1=input code 和 seq2=reverse input code 并沿着 seq1 移动 seq2 寻找三个相邻对，但这并没有给我正确的输出。

4

2 回答 2

4

您是否考虑过使用itertools的产品。然后您可以迭代结果并仅选择您喜欢的这些结果。

于 2012-04-24T23:00:44.543 回答
3

如果你的 RNA 不是很长（一千个碱基可能还可以；几十万个肯定不行），你可以使用简单的 O(n^3) 算法。O(n^3) 意味着执行时间在最坏的情况下与碱基数的立方成正比。作者提到嵌套循环在很大程度上暗示了这种简单但相当缓慢的方法。

def find_loops(rna, min_pairs=3, min_loop=3):
    n = len(rna)
    result = []
    for loop_start in xrange(min_pairs, n - min_pairs - min_loop + 1):
        for loop_end in xrange(loop_start + min_loop, n - min_pairs):
            if (loop_end - loop_start < min_loop + 2 or 
                    not base_pair(rna[loop_start], rna[loop_end - 1])):
                max_pairs = min(loop_start, n - loop_end)
                for k in xrange(max_pairs):
                    if not base_pair(rna[loop_start - k - 1], rna[loop_end + k]):
                        break
                else:
                    k = max_pairs
                if k >= min_pairs:
                    result.append((loop_start - k, k, loop_end - loop_start))
    return result

def base_pair(x, y):
    return (x == 'A' and y == 'U' or
            x == 'C' and y == 'G' or
            x == 'G' and y == 'C' or
            x == 'U' and y == 'A')

这会遍历所有可能的 RNA 环的开始和结束，然后在两个方向上从潜在环的末端走开，只要碱基仍然配对。当它到达一对不匹配的碱基时，它会停下来检查它是否至少有最少的碱基对。如果有，它将循环添加到结果列表中。

第一个if是为了避免列出可能被“压缩”得更紧的循环。如条件所示，如果循环太短（少于五个碱基）或末端不匹配，则无法将循环拉得更紧。

结果是一个元组列表，每个可能的循环一个，形式为(start_pos, pair_count, loop_length)。这意味着pair_count从碱基编号开始的碱基序列start_pos后面是一个loop_length碱基环，然后是反向的互补序列。该序列的反义拷贝从碱基开始start_pos + pair_count + loop_length。第一个基数是数字 0，而不是 1（我们是这里的程序员）。

一个例子可能会更清楚：print find_loops('GGGGAUUACAGCGUGUAAUCAAUA')returns [(4, 3, 13), (3, 7, 3)]，也就是说，它找到了两个循环：

  • 在位置 4 处，三个碱基 ,AUU围成一个 13 个碱基的环，并与AAU位置 20 结合；
  • 在位置 3，七个碱基，GAUUACA包围三个碱基的环，并与UGUAAUC位置 13 结合。

如果没有第一个if，该函数还将返回类似 (3, 6, 5) 的循环（即GAUUAC在位置 3 包含五个碱基的循环并绑定到GUAAUC位置 14），这与 (3, 7, 3) 相同) 上面，只是没有拉得那么紧。

希望这可以帮助。如果您需要更快的算法，我相信有一个动态编程解决方案可以处理更长的字符串。让我知道，我会考虑的。不过，这不会那么容易理解......

于 2012-04-25T00:17:42.523 回答