Pythonで文字列検索を最適化する

Question

大規模な50 MBのDNA配列と15文字程度の小さなものを与えたPythonプログラムを書く必要があります。これは15文字すべての配列のリストを返します。 1つは与えられただけでなく、大きなもののどこにあるのか。

私の現在のアプローチは、最初にすべてのサブシーケンスを取得することです。その後、

def get_subsequences_of_size(size, data): 
    sequences = {} 
    i = 0 
    while(i+size <= len(data)): 
     sequence = data[i:i+size] 
     if sequence not in sequences: 
      sequences[sequence] = data.count(sequence) 
     i += 1 
    return sequences 

と（私はポジションを得るために忘れてしまった）問題が尋ねた内容に従って辞書のリストにそれらをパック：

def find_similar_sequences(seq, data): 
    similar_sequences = {} 
    sequences = get_subsequences_of_size(len(seq), data) 
    for sequence in sequences.keys(): 
     diffs, muts = calculate_similarity(seq,sequence) 
     if diffs not in similar_sequences: 
      similar_sequences[diffs] = [{"Sequence": sequence, "Mutations": muts}] 
     else: 
      similar_sequences[diffs].append({"Sequence": sequence, "Mutations": muts}) 
     #similar_sequences[sequence] = {"Similarity": (len(sequence)-diffs), "Differences": diffs, "Mutatations": muts} 
    return similar_sequences 

私の問題は、この実行方法が遅すぎることです。 50MBの入力では、処理を終了するのに30分以上かかります。

Answer 1

次のアプローチについてはどう：あなたの長いシーケンスを超えると、すべてのサブシーケンスのための長さ15のスライディングウィンドウと

ゴー：

• 店舗長いシーケンスの開始位置
• 計算し、類似

import re 
from itertools import islice 
from collections import defaultdict 

short_seq = 'TGGCGACGGACTTCA' 
long_seq = 'AGAACGTTTCGCGTCAGCCCGGAAGTGGTCAGTCGCCTGAGTCCGAACAAAAATGACAACAACGTTTATGACAGAACATT' +\ 
      'CCTTGCTGGCAACTACCTGAAAATCGGCTGGCCGTCAGTCAATATCATGTCCTCATCAGATTATAAATGCGTGGCGCTGA' +\ 
      'CGGATTATGACCGTTTTCCGGAAGATATTGATGGCGAGGGGGATGCCTTCTCTCTTGCCTCAAAACGTACCACCACATTT' +\ 
      'ATGTCCAGTGGTATGACGCTGGTGGAGAGTTCCCCCGGCAGGGATGTGAAGGATGTGAAATGGCGACGGACTTCACCGCA' +\ 
      'TGAGGCTCCACCAACCACGGGGATACTGTCGCTCTATAACCGTGGCGATCGCCGTCGCTGGTACTGGCCCTGTCCACACT' +\ 
      'GTGGTGAGTATTTTCAGCCCTGCGGCGATGTGGTTGCTGGTTTCCGTGATATTGCCGATCCCGTGCTGGCAAGTGAGGCG' +\ 
      'GCTTATATTCAGTGTCCTTCTGGCGACGGACTTCACGCGTCAGCCCGGAAGTGGTCAGTCGCCTGAGTCCGAACAAAAAT' 


def window(seq, n=2): 
    "Returns a sliding window (of width n) over data from the iterable" 
    " s -> (s0,s1,...s[n-1]), (s1,s2,...,sn), ...     " 
    # from https://docs.python.org/release/2.3.5/lib/itertools-example.html 
    it = iter(seq) 
    result = tuple(islice(it, n)) 
    if len(result) == n: 
     yield ''.join(result) 
    for elem in it: 
     result = result[1:] + (elem,) 
     yield ''.join(result) 

def hamming_distance(s1, s2): 
    if len(s1) != len(s2): 
     raise ValueError("Undefined for sequences of unequal length") 
    return sum(ch1 != ch2 for ch1, ch2 in zip(s1, s2)) 

k = len(short_seq) 
locations = defaultdict(list) 
similarities = defaultdict(set) 

for start, subseq in enumerate(window(long_seq, k)): 
    locations[subseq].append(start) 
    similarity = hamming_distance(subseq, short_seq) # substitute with your own similarity function 
    similarities[similarity].add(subseq) 

with open(r'stack46268997.txt', 'w') as f: 
    for similarity in sorted(similarities.keys()): 
     f.write("Sequence(s) which differ in {} base(s) from the short sequence:\n".format(similarity)) 
     for subseq in similarities[similarity]: 
      f.write("{} at location(s) {}\n".format(subseq, ', '.join(map(str, locations[subseq])))) 
     f.write('\n') 

を保存

これは、指定されたシーケンスにどれくらい近いかによって順序付けられたサブシーケンスのリストを出力します。

Sequence(s) which differ in 0 base(s) from the short sequence: 
TGGCGACGGACTTCA at location(s) 300, 500 

Sequence(s) which differ in 5 base(s) from the short sequence: 
TGGCGATCGCCGTCG at location(s) 362 

Sequence(s) which differ in 6 base(s) from the short sequence: 
TGGCAACTACCTGAA at location(s) 86 
TGGTGAGTATTTTCA at location(s) 401 
TGGCGAGGGGGATGC at location(s) 191 

Sequence(s) which differ in 7 base(s) from the short sequence: 
ATGTGAAGGATGTGA at location(s) 283 
AGGGGGATGCCTTCT at location(s) 196 
TGACAACAACGTTTA at location(s) 53 
CGCTGACGGATTATG at location(s) 154 
TTATGACCGTTTTCC at location(s) 164 
TGGTTGCTGGTTTCC at location(s) 430 
TCGCGTCAGCCCGGA at location(s) 8 
AGTCGCCTGAGTCCG at location(s) 30, 536 
CGGCGATGTGGTTGC at location(s) 422 

[... and so on...] 

また、50 MB FASTAファイルでスクリプトを実行しました。私のマシンでは、結果を計算するには42秒かかり、結果をファイルに書き出すにはもう30秒かかりました。（印刷すると時間がかかりました）