Как эффективно перебрать большую строку и записать индекс каждой трехбуквенной подстроки?

Question

У меня есть очень большая последовательность, которую я читаю в виде строки (250 000 000 букв). Буквы G, A, C, T.

Например:

'GACTCGTAGCTAGCTG'

Я хотел бы создать способ хранения индекса каждой трехбуквенной подстроки в моей последовательности для использования в другой функции позже.

Например:

{'GAC': 1, 'ACT': 2 'CTC':3, 'TCG':4 ...}

В моем текущем подходе моя проблема заключается в том, что я не нашел эффективного способа сохранить индексы каждой трехбуквенной подстроки в моей последовательности. Как только я знаю индексы для каждой подстроки, я случайным образом выберу некоторые из них на основе имеющихся у меня вероятностей и заменю их на другую известную подстроку.

Я попытался выполнить итерацию, используя для l oop Трехбуквенная подстрока за раз в скользящем окне и сохранение позиций в виде значений словаря с трехбуквенной подстрокой в ​​качестве ключа. Также, когда я сохраняю этот файл словаря, я использую pd.to_csv, но он кажется очень неэффективным и создает файл размером 8 ГБ. Однако это займет слишком много времени для очень большой строки (250 000 000 букв).

Answer 1

Ваш пример выходных данных предполагает, что вам нужен только один индекс для каждой подпоследовательности, но ваше описание подразумевает, что вам, вероятно, понадобятся все индексы.

Вот функция, которая создаст словарь со всеми индексами для подпоследовательностей ( заданной длины), которые присутствуют:

from collections import defaultdict
def sequenceDict(seq,count):
    result = defaultdict(list)
    for i,subseq in enumerate(zip(*(seq[i:] for i in range(count)))):
        result["".join(subseq)].append(i)
    return result

r = sequenceDict('GACTCGTAGCTAGCTG',3)
print(r)

# {'GAC': [0], 'ACT': [1], 'CTC': [2], 'TCG': [3], 'CGT': [4], 'GTA': [5], 'TAG': [6, 10], 'AGC': [7, 11], 'GCT': [8, 12], 'CTA': [9], 'CTG': [13]})

Если вы действительно хотите только первый индекс каждой трехбуквенной подпоследовательности, словарь можно получить намного быстрее, используя одно словарное понимание:

from itertools import product
{ ss:sequence.index(ss) for p in product(*["ACGT"]*3)for ss in ["".join(p)] if ss in sequence}

Я провел тест производительности для случайной последовательности из 250 миллионов букв, и словарь с одним индексом можно получить за несколько микросекунд. Получение всех индексов занимает чуть больше минуты (с помощью функции выше):

import time

size = 250_000_000
print("loading sequence...",size)
start = time.time()
import random
sequence = "".join(random.choice("ACGT") for _ in range(size))
print("sequence ready",time.time()-start)


start = time.time()
from itertools import product
seqDict = { ss:sequence.index(ss) for p in product(*["ACGT"]*3)for ss in ["".join(p)] if ss in sequence}
print("\n1st index",time.time()-start)

start = time.time()
r = sequenceDict(sequence,3)
print("\nall indexes",time.time()-start)

output:

loading sequence... 250000000
sequence ready 193.82172107696533

1st index 0.000141143798828125

all indexes 71.74848103523254

Учитывая, что время для загрузки последовательностей составляет намного больше, чем время на создание индекса, вы можете просто go сохранить этот словарь индекса и просто перестраивать его каждый раз из исходных данных (которые вам все равно нужно загружать для вашего процесса)

Вы также можете просто сохранить словарь счетчиков и извлекать индексы по требованию по мере необходимости:

Эта функция получает количество вхождений для каждой подпоследовательности:

from collections import Counter
def countSubSeqs(seq,size):
    return Counter("".join(s) for s in zip(*(seq[i:] for i in range(size))))

Он выполняется примерно в то же время, что и функция sequenceDict, но создает гораздо меньший словарь.

Чтобы получить индексы указанной подпоследовательности c (включая перекрывающиеся позиции), вы можете использовать это:

subSeq  = "ACT"
indexes = [ i for i in range(len(sequence)) if sequence[i:i+3]==subSeq ] 

Если вам не нужны сразу все индексы для всех подпоследовательностей, вы, вероятно, могли бы соответствующим образом структурировать свой код и получить индексы, которые en, необходимый (и возможно сохраните их в словаре для повторного использования)

Answer 2

Вы можете разбить длинную строку на 3-буквенные подстроки:

string='GACTCGTAGCTAGCT'
substrings=[string[3*x:3*x+3] for x in range(int(len(string)/3))]

substrings будет:

['GAC', 'TCG', 'TAG', 'CTA', 'GCT']

Получить указатели этих в другом списке:

indices=[x for x in range(int(len(string)/3))]

Это просто даст:

[0, 1, 2, 3, 4]

Элемент n в списке подстрок будет соответствовать элементу n в списке индексов.

Как поместить файл в строковую переменную, на которую вы можете посмотреть: Как прочитать текстовый файл в строковую переменную и удалить символы новой строки?