如何爲生物信息學查詢改進一個python腳本

我對python相當陌生，如果可能，我將不勝感激。我正在比較兩個密切相關的生物體[E_C & E_F]的基因組並嘗試識別一些基本的插入和缺失。我使用兩種生物體的序列進行了FASTA配對比對（glsearch36）。如何爲生物信息學查詢改進一個python腳本

以下是我的python腳本的一部分，我已經能夠在一個序列（數據庫）中識別與另一個序列（查詢）中的空位相對應的7個核苷酸（七聚體）。這是什麼，我有一個例子：

ATGCACAA-ACCTGTATG # query 
ATGCAGAGGAAGAGCAAG # database 
9 
GAGGAAG

假設缺口是在9號位置。我試圖改進腳本以選擇20個核苷酸或更多分開放在兩個序列且僅當週圍的核苷酸差距也匹配

ATGCACAAGTAAGGTTACCG-ACCTGTATGTGAACTCAACA 
       ||| ||| 
GTGCTCGGGTCACCTTACCGGACCGCCCAGGGCGGCCCAAG 
21 
CCGGACC

這是我的腳本部分，上半部分涉及打開不同的文件。它還打印一個字典，並在每個序列的末尾計數。

list_of_positions = [] 

for match in re.finditer(r'(?=(%s))' % re.escape("-"), dict_seqs[E_C]): 
    list_of_positions.append(match.start()) 

set_of_positions = set(list_of_positions) 
for position in list_of_positions: 
    list_no_indels = [] 
    for number in range(position-20, position) : 
     list_no_indels.append(number) 
    for number in range(position+1, position+21) : 
     list_no_indels.append(number) 
    set_no_indels = set(list_no_indels) 
    if len(set_no_indels.intersection(set_of_positions))> 0 : continue 

    if len(set_no_indels.intersection(set_of_positions_EF))> 0 : continue 


    print position 
    #print match.start() 

    print dict_seqs[E_F][position -3:position +3] 

    key = dict_seqs[E_F][position -3: position +3] 

    if nt_dict.has_key(key): 
     nt_dict[key] += 1 
    else: 
     nt_dict[key] = 1 


print nt_dict

從本質上講，我試圖編輯兩兩比對的結果，試圖找出核苷酸同時在查詢和數據庫序列的空白相反，以進行一些基本的插入/缺失分析。

爲了減少噪音，我能夠通過將間隙「 - 」之間的距離增加到20 nt來解決我之前的一個問題，這改善了我的結果。上面編輯的腳本。

這是我的結果的一個例子，最後我有一個字典來計算每個序列的出現次數。

ATGCACAA-ACCTGTATG # query 
ATGCAGAGGAAGAGCAAG # database 
9 (position on the sequence) 
GAGGAA (hexamer) 


ATGCACAAGACCTGTATG # query 
ATGCAGAG-AAGAGCAAG # database 
9 (position) 
CAAGAC (hexamer)

但是，我仍然試圖修復腳本，在那裏我得到的間隙核苷酸完全匹配如此，在哪裏|只是爲了顯示每個序列的匹配核苷酸：

GGTTACCG-ACCTGTATGTGAACTCAACA # query 
    ||| || 
CCTTACCGGACCGCCCAGGGCGGCCCAAG # database 

9 
ACCGAC

任何幫助，將不勝感激！

來源

2013-10-31 sheaph

對於生物信息學來說，這是+2還是至少有2個人能夠理解你的問題？看來我不能。你能否澄清你的問題，即你正在尋找什麼答案？一些測試用例？ – alko

@alko我提供了一些更多的細節。任何建議都會很棒。謝謝！ – sheaph

稍微解釋一下您的術語可能會有幫助。例如，「查詢」和「數據庫」是什麼意思？你的「查詢」實際上是否被髮送到數據庫中？在你的第一個例子中，你有'ATGCACAA-ACCTGTATG'（查詢），'ATGCAGAGGAAGAGCAAG'（數據庫），然後拉出'GAGGAAG'（hexamer），它出現在你的數據庫字符串中，但不在你的查詢字符串中。如果您將該模式從數據庫中提取出來，則會留下「ATGCA-AGCAAG」，這似乎與您的查詢無關。你能否澄清數據庫，查詢和haxamer子集之間的關係？ – jeremiahbuddha

我想我明白你要做什麼，但正如@alko所說 - 在你的代碼中的評論肯定會有很大的幫助。

至於發現周圍間隙的精確匹配，你可以運行一個字符串比較：

線沿線的東西：

if query[position -3: position] == database[position -3: position] and query[position +1: position +3] == database[position +1: position +3]: 
    # Do something

您需要替換「查詢」和「數據庫」與你稱之爲你想要比較的字符串。

來源

2013-11-06 15:43:13 Qui

如何爲生物信息學查詢改進一個python腳本

回答

相關問題