Python/Biopython: 如何在有缺口的序列中查找参考序列(字符串)?
Python/Biopython: How to search for a reference sequence (string) in a sequence with gaps?
我遇到了以下问题,还没有找到解决办法:
我正在开发一种序列分析工具,它使用带有参考序列的文件并尝试在测试序列中找到这些参考序列之一。
问题是测试序列可能包含间隙(例如:ATG---TCA
)。
我希望我的工具找到一个特定的参考序列作为测试序列的子字符串,即使参考序列在测试序列中被间隙 (-
) 打断。
例如:
我的参考序列之一:
a = TGTAACGAACGG
我的测试顺序:
b = ACCT**TGT--CGAA-GG**AGT
(来自参考序列的相应部分以粗体给出)
我考虑过正则表达式并尝试将自己投入其中,但如果我没记错的话,正则表达式只能以相反的方式工作。所以我需要将间隙位置作为正则表达式包含到参考序列中,而不是将其映射到测试序列。
但是,我不知道测试序列中的位置、长度和间隙数。
我的想法是将测试序列字符串中的间隙位置(因此所有 -
)交换为某种正则表达式或代表参考序列中任何其他字符的特殊字符。比我将未修改的参考序列与我修改的测试序列进行比较......
不幸的是,我没有在 python 中找到用于字符串搜索的函数或一种可以做到这一点的正则表达式。
非常感谢!
你可以这样做:
import re
a = 'TGTAACGAACGG'
b = 'ACCTTGT--CGAA-GGAGT'
temp_b = re.sub(r'[\W_]+', '', b) #removes everything that isn't a number or letter
if a in temp_b:
#do something
有好消息也有坏消息...
首先是坏消息:您尝试做的事情并不容易,正则表达式确实不是实现它的方法。在一个简单的情况下,正则表达式可以工作(也许),但它效率低下并且无法扩展。
然而,好消息是这是生物信息学中很好理解的问题(例如参见 [=11=]
编辑
从下面的讨论来看,您似乎是在说 'b' 可能很长,但假设 'a' 仍然很短(上面示例中的 12 个碱基)我认为您可以通过遍历每个'b' 中的 12 聚体。 IE。将 'b' 分成 12 个碱基长的序列(显然你最终会得到很多!)。然后您可以轻松地比较这两个序列。如果你真的想使用正则表达式(我仍然建议你不要)那么你可以用'.'替换'-'。并做一个简单的匹配。例如。
import re
''' a is the reference '''
a = 'TGTAACGAACGG'
''' b is 12-mer taken from the seqence of interest, in reality you'll be doing this test for every possible 12-mer in the sequence'''
b = 'TGT--CGAA-GG'
b = b.replace('-', '.')
r = re.compile(b);
m = r.match(a)
print(m)
我遇到了以下问题,还没有找到解决办法:
我正在开发一种序列分析工具,它使用带有参考序列的文件并尝试在测试序列中找到这些参考序列之一。
问题是测试序列可能包含间隙(例如:ATG---TCA
)。
我希望我的工具找到一个特定的参考序列作为测试序列的子字符串,即使参考序列在测试序列中被间隙 (-
) 打断。
例如:
我的参考序列之一:
a = TGTAACGAACGG
我的测试顺序:
b = ACCT**TGT--CGAA-GG**AGT
(来自参考序列的相应部分以粗体给出)
我考虑过正则表达式并尝试将自己投入其中,但如果我没记错的话,正则表达式只能以相反的方式工作。所以我需要将间隙位置作为正则表达式包含到参考序列中,而不是将其映射到测试序列。
但是,我不知道测试序列中的位置、长度和间隙数。
我的想法是将测试序列字符串中的间隙位置(因此所有 -
)交换为某种正则表达式或代表参考序列中任何其他字符的特殊字符。比我将未修改的参考序列与我修改的测试序列进行比较......
不幸的是,我没有在 python 中找到用于字符串搜索的函数或一种可以做到这一点的正则表达式。
非常感谢!
你可以这样做:
import re
a = 'TGTAACGAACGG'
b = 'ACCTTGT--CGAA-GGAGT'
temp_b = re.sub(r'[\W_]+', '', b) #removes everything that isn't a number or letter
if a in temp_b:
#do something
有好消息也有坏消息...
首先是坏消息:您尝试做的事情并不容易,正则表达式确实不是实现它的方法。在一个简单的情况下,正则表达式可以工作(也许),但它效率低下并且无法扩展。
然而,好消息是这是生物信息学中很好理解的问题(例如参见 [=11=]
编辑 从下面的讨论来看,您似乎是在说 'b' 可能很长,但假设 'a' 仍然很短(上面示例中的 12 个碱基)我认为您可以通过遍历每个'b' 中的 12 聚体。 IE。将 'b' 分成 12 个碱基长的序列(显然你最终会得到很多!)。然后您可以轻松地比较这两个序列。如果你真的想使用正则表达式(我仍然建议你不要)那么你可以用'.'替换'-'。并做一个简单的匹配。例如。
import re
''' a is the reference '''
a = 'TGTAACGAACGG'
''' b is 12-mer taken from the seqence of interest, in reality you'll be doing this test for every possible 12-mer in the sequence'''
b = 'TGT--CGAA-GG'
b = b.replace('-', '.')
r = re.compile(b);
m = r.match(a)
print(m)