在给定结构中寻找RNA碱基配对

Question

我有一个RNA序列，看起来像这样。每个字符代表一个核苷酸(程序员可以忽略这一点，您可以将它们视为元素)：

         (((((((..((((.....(..)))).((((.........)))).....(((((..)....))))))))))))....

我会用这个名称来回答这个问题：

          ( = lhb
          ) = rhb 
          . = dot

因此，本质上，具有lhb的元素与具有rhb的元素相连接，而所有的点都是自由的，它们的连接方式是复杂的。很难用文字表达，所以为了方便起见，我会把数字放在下面一些元素的下面，这些元素是相连的：

   (((((((..((((.....(..)))).((((.........)))).....(((((..)....))))))))))))....
   1                 2  2       3         3           45  5    4          1   

我想这会给你一个关于它们是如何联系的概念。我有兴趣找出连接的元素和免费的区域的位置。

(例如，要素1与要素72连接，要素8至9是自由的)。

我选择了C来编写代码，但我并不接近逻辑。

而且，用C编写程序变得越来越困难，我觉得在python中使用正则表达式可以很容易地做到这一点，或者可能是perl，但我对这些语言没有多少经验。所以，如果有人能提供一个更简单的方法，那将是一个巨大的帮助。改进C代码的想法也受到欢迎。这是我的C代码:-

#include <stdio.h>

int main() {

char dot[500];
int i = 0, j = 0;
int count = 0, count1 = 0, count2 = 0;
int lhb[100];
int rhb[100];
int dots[100];
int pair_1[100];
int pair_2[100];
int pair_3[100];
FILE * fp;

fp = fopen("structure.txt", "r");

while (fscanf(fp, "%c", & dot[i]) != EOF) {

  i++;
}

fclose(fp);

for (i = 0; dot[i] != '\0'; i++) {

  if (dot[i] == '(') {
    lhb[count] = dot[i];
    pair_1[count] = i;
    count++;
    } 
  else if (dot[i] == '.') {
    rhb[count1] = dot[i];
    pair_2[count1] = i;
    count1++;
    }
  else {
    dots[count2] = dot[i];
    pair_3[count2] = i;
    count2++;
    }
}


printf("Base-pair details :\n");

for (j = 0; j < count; j++)

  printf("%d--%d\n", pair_1[j] + 1, pair_3[count - j - 1] + 1);

printf("Loop details :\n");

// for(j=0;j<=count;j++)

// printf("--%d-",pair_2[j]+1);

return 0;

}

Answer 1

这是一个潜在的解决方案。Free包含一个对每个点进行索引的列表，而对包含一个包含每对索引的元组的列表。这是假设您的数据是完美的(即相同数量的左对和右对，在左对之前没有右对)。然而，这可以修改，以检查边缘的情况。根据你的问题，我从1开始索引，而不是0。

data = '(((((((..((((.....(..)))).((((.........)))).....(((((..)....))))))))))))....'

left = []
group = []
prev = ''

free = []
pairs = []

for index, elem in enumerate(data, 1):
    if elem == '.' and prev == '.':
        group.append(index)
    elif elem == '.':
        group = [index]
    else:
        if len(group) >= 1:
            free.append(group)
            group = []
        if elem == '(':
            left.append(index)
        elif elem == ')':
            pairs.append([left.pop(), index])
    prev = elem
if len(group) > 0:
    free.append(group)
pairs.sort()

输出：

[[1, 72], [2, 71], [3, 70], [4, 69], [5, 68], [6, 67], [7, 66], [10, 65], [11, 25], [12, 24], [13, 23], [19, 22], [27, 43], [28, 42], [29, 41], [30, 40], [49, 64], [50, 63], [51, 62], [52, 61], [53, 56]]
[[8, 9], [14, 15, 16, 17, 18], [20, 21], [26], [31, 32, 33, 34, 35, 36, 37, 38, 39], [44, 45, 46, 47, 48], [54, 55], [57, 58, 59, 60]]

代码的工作方式是在每个元素上迭代一次。如果遇到左paren，它会将索引追加到列表中。当它遇到一个右paren时，它会弹出最后一次看到的左paren的索引，创建一个对。对于自由点，跟踪前一个元素将告诉您如何处理每个点。如果已经看到一个列表，则继续添加到当前列表中，否则启动一个新列表。

Answer 2

真正解决这个问题的最好方法是维护一个堆栈。对于遇到的每个(，您都可以将索引值推到堆栈中，对于每个)，您需要用插入上一个(索引的数字来弹出堆栈。这意味着(的索引和)的索引形成一对。

这可以通过做这样的事情来实现

seq = '(((((((..((((.....(..)))).((((.........)))).....(((((..)....))))))))))))....'
xStack = []
for i, x in enumerate(seq):
    if x == '(':
        xStack.append(i)
    if x == ')':
        o = xStack.pop()

现在已经完成了基本步骤，除了维护括号的索引之外，还需要一些其他的东西。在pop操作之后，您需要存储匹配的对，为此，让我们引入另一个变量，在遇到.时基本上什么也不做

seq = '(((((((..((((.....(..)))).((((.........)))).....(((((..)....))))))))))))....'

xStack = []

resultingPairs = []

for i, x in enumerate(seq):
    if x == '(':
        xStack.append(i)
    if x == ')':
        o = xStack.pop()
        tempPair = [o, i]
        resultingPairs.append(tempPair)
    if x == '.':
        pass

现在，我们有了如下所示的结果对

[[18, 21], [12, 22], [11, 23], [10, 24], [29, 39], [28, 40], [27, 41], [26, 42], [52, 55], [51, 60], [50, 61], [49, 62], [48, 63], [9, 64], [6, 65], [5, 66], [4, 67], [3, 68], [2, 69], [1, 70], [0, 71]]

我们需要找出所有的空闲空间，这可以很容易地做一些事情，如下所示

spacesList = [i for i in range(len(seq)) if seq.startswith('.', i)]

这会导致

[7, 8, 13, 14, 15, 16, 17, 19, 20, 25, 30, 31, 32, 33, 34, 35, 36, 37, 38, 43, 44, 45, 46, 47, 53, 54, 56, 57, 58, 59, 72, 73, 74, 75]

现在，您可以轻松地编写一个函数，在其中传递spacesList和resultingPairs，并获得每个可能的对之间的空闲空格数。有很多可能的优化，但这应该能够使您开始朝着正确的方向。

def getSpacesCount(spacesList, resultingPairs):
    for pair in resultingPairs:
        a = pair[0]
        b = pair[1]
        spacesCount = 0
        for val in spacesList:
            if a < val < b:
                spacesCount+=1
        print a,b,spacesCount

seq = '(((((((..((((.....(..)))).((((.........)))).....(((((..)....))))))))))))....'
xStack = []
resultingPairs = []
for i, x in enumerate(seq):
    if x == '(':
        xStack.append(i)
    if x == ')':
        o = xStack.pop()
        tempPair = [o, i]
        resultingPairs.append(tempPair)
    if x == '.':
        pass

spacesList = [i for i in range(len(seq)) if seq.startswith('.', i)]

getSpacesCount(spacesList, resultingPairs)

你有开始括号的位置，结束括号的位置，以及它们之间的自由空格数。

>>> getSpacesCount(spacesList, resultingPairs)
18 21 2
12 22 7
11 23 7
10 24 7
29 39 9
28 40 9
27 41 9
26 42 9
52 55 2
51 60 6
50 61 6
49 62 6
48 63 6
9 64 28
6 65 30
5 66 30
4 67 30
3 68 30
2 69 30
1 70 30
0 71 30

编辑似乎无法理解如何编辑堆栈溢出的注释，将函数更新为

def getSpacesCount(spacesList, resultingPairs):
    for pair in resultingPairs:
        a = pair[0]
        b = pair[1]
        spacesCount = 0
        spaces = []
        for val in spacesList:
            if a < val < b:
                spaces.append(val)
                spacesCount+=1
        print a,b,spacesCount,spaces

这会给你数量和位置。你可以留着你真正喜欢的任何一个。

>>> getSpacesCount(spacesList, resultingPairs)
18 21 2 [19, 20]
12 22 7 [13, 14, 15, 16, 17, 19, 20]
11 23 7 [13, 14, 15, 16, 17, 19, 20]
10 24 7 [13, 14, 15, 16, 17, 19, 20]
29 39 9 [30, 31, 32, 33, 34, 35, 36, 37, 38]
28 40 9 [30, 31, 32, 33, 34, 35, 36, 37, 38]
27 41 9 [30, 31, 32, 33, 34, 35, 36, 37, 38]
26 42 9 [30, 31, 32, 33, 34, 35, 36, 37, 38]
52 55 2 [53, 54]
51 60 6 [53, 54, 56, 57, 58, 59]
50 61 6 [53, 54, 56, 57, 58, 59]
49 62 6 [53, 54, 56, 57, 58, 59]
48 63 6 [53, 54, 56, 57, 58, 59]
9 64 28 [13, 14, 15, 16, 17, 19, 20, 25, 30, 31, 32, 33, 34, 35, 36, 37, 38, 43, 44, 45, 46, 47, 53, 54, 56, 57, 58, 59]
6 65 30 [7, 8, 13, 14, 15, 16, 17, 19, 20, 25, 30, 31, 32, 33, 34, 35, 36, 37, 38, 43, 44, 45, 46, 47, 53, 54, 56, 57, 58, 59]
5 66 30 [7, 8, 13, 14, 15, 16, 17, 19, 20, 25, 30, 31, 32, 33, 34, 35, 36, 37, 38, 43, 44, 45, 46, 47, 53, 54, 56, 57, 58, 59]
4 67 30 [7, 8, 13, 14, 15, 16, 17, 19, 20, 25, 30, 31, 32, 33, 34, 35, 36, 37, 38, 43, 44, 45, 46, 47, 53, 54, 56, 57, 58, 59]
3 68 30 [7, 8, 13, 14, 15, 16, 17, 19, 20, 25, 30, 31, 32, 33, 34, 35, 36, 37, 38, 43, 44, 45, 46, 47, 53, 54, 56, 57, 58, 59]
2 69 30 [7, 8, 13, 14, 15, 16, 17, 19, 20, 25, 30, 31, 32, 33, 34, 35, 36, 37, 38, 43, 44, 45, 46, 47, 53, 54, 56, 57, 58, 59]
1 70 30 [7, 8, 13, 14, 15, 16, 17, 19, 20, 25, 30, 31, 32, 33, 34, 35, 36, 37, 38, 43, 44, 45, 46, 47, 53, 54, 56, 57, 58, 59]
0 71 30 [7, 8, 13, 14, 15, 16, 17, 19, 20, 25, 30, 31, 32, 33, 34, 35, 36, 37, 38, 43, 44, 45, 46, 47, 53, 54, 56, 57, 58, 59]