如何使用 Python 基于较长序列提取较短序列的位置?

How to extract the location of shorter sequence based on longer sequence using Python?

我有一个带有序列ID的文件A,还有结合位点位置的信息。我想只提取位置信息而没有 A、T、C、G 信息。较长序列上方的较短序列显示其位置和左侧的每个数字,例如在文件 A 中,值 451 是左侧的位置值。我想在以 453(起始站点)开始的较长序列上获取短序列的位置,并获取大小为 21 的较短序列的长度,然后添加到 453 以获得 474 的结束站点。任何人都可以帮助我?

文件A.txt

chr1:152806601-152807450

        TTCAGCACCATGGACAGCGCC          
451  GGCTTCAGCACCACGGACAGCGCCCCACCCGCGGCCCTCCCCCCGGCGGCGCGCTCCAGCCGGTGTAGGCGAGGC

            TTCAGCACCATGGACAGCGCC            
751  AGAGCCCCCCGGGACTGCAGAGAGCACCTGGGAGGCTGGACTGGGAACGAGACATACTCGAAGGAGTAAGTGAAG


chr10:125364276-125364825

                       TTCAGCACCATGGACAGCGCC
301  CAGTAATGTGGGGTTGTGGTCAGCACCATGGACAGCTCCCCTGTTGCTTCATATTGAGGAATAGGAAAGCGCCGC

        TTCAGCACCATGGACAGCGCC
376  TATCTCCGGATCCTGGCTAGCTCCAGCCACTGCAGGTAACTGTCTTGAATGGGCTTAGAAACATGGTGATGTCTG

期望输出

chr1:152806601-152807450 453 474
chr1:152806601-152807450 757 778
chr10:125364276-125364825 318 339
chr10:125364276-125364825 378 399

示例代码

import re
with open("A.txt", "r") as f:
    lines = f.readlines()


label_ptrn = re.compile("")   # insert regular expression sequence ID
line_ptrn = re.compile("")    # insert regular expression start site
inner_ptrn = re.compile("")   # insert regular expression end site

all_matches = []
for line in lines:
    m = label_ptrn.match(line)
    if m:
        label = m.groupdict().get("label")
        continue
    m = line_ptrn.match(line)
    if m:
        start = m.groupdict().get("start_value")
        sequence = m.groupdict().get("sequence")
        mi = inner_ptrn.search(sequence)
        if not mi:
            continue
        span = mi.span()
        all_matches.append((label, int(start)+span[0], int(start)+span[1]))

with open("A_ouput.bed", "w+b") as f:
    for m in all_matches:
        f.write('%s\t%i\t%i\n' % m)

所以在我看来,你想要的输出的起始位置偏移了一位。

        TTCAGCACCATGGACAGCGCC          
451  GGCTTCAGCACCACGGACAGCGCCCCACCCGCGGCCCTCCCCCCGGCGGCGCGCTCCAGCCGGTGTAGGCGAGGC

较短序列中的第一个 T 看起来像是在较长序列的第四个字符上方。如果较长序列的第一个字符位于位置 451,则较短序列的第一个字符位于位置 454。

如果文件结构是常量这里是一个非正则表达式的解决方案。

result = []
with open('file.txt') as f:
    for line in f:
        if line.startswith('chr'):
            label = line.strip()
        elif line[0] == ' ':
            # short sequence
            length = len(line.strip())
            # find the index of the beginning of the short sequence
            for i, c in enumerate(line):
                if c.isalpha():
                    short_index = i
                    break
        elif line[0].isdigit():
            # long sequence
            n = line.split(' ')[0]
            # or
            # n = line[:line.index(' ')]
            # find the index of the beginning of the long sequence
            for i, c in enumerate(line):
                if c.isalpha():
                    long_index = i
                    break
            start = int(n) + short_index - long_index
            # start -= 1
            end = start + length
            result.append('{} {} {}'.format(label, start, end))
            offset, n, start, length = 0, 0, 0, 0

结果

['chr1:152806601-152807450 454 475',
 'chr1:152806601-152807450 758 779',
 'chr10:125364276-125364825 319 340',
 'chr10:125364276-125364825 379 400']

如果我误解了您的示例数据,请取消注释 start -= 1