寫一個(gè)腳本寓盗，完成以下任務(wù)：
有一個(gè)文本文件A如下：barcode.fa灌砖，其中包含了若干長度8bp至11bp的DNA序列璧函。

>bc1_0
GTTTGTTT
>bc1_1
ACCGTGTTT
>bc1_2
GATAGTGTTT
>bc1_3
TGAGGCGGTTT
>bc1_4
GATCGTTT
>bc1_5
ATCACGTTT
>bc1_6
GATGTAGTTT
>bc1_7
TGACACAGTTT
>bc1_8
CTTTCTTT
>bc1_9
AGCCTCTTT
>bc1_10
GACGGGCTTT

另有一個(gè)fastq文件B，fastq文件不做介紹了基显。
要求對此fastq文件進(jìn)行處理蘸吓，輸出滿足以下條件的序列：
1）序列的前8bp-11bp與前述文本文件A中的DNA序列hamming distance不大于2；
2）能夠唯一匹配到文本文件A中的某一條DNA序列
（例如撩幽，如果fastq中某條序列的前8bp-11bp在文本文件中沒有完全匹配的DNA序列库继，而在A文件中有兩條或以上的DNA序列的hamming distance為1，則拋棄該序列）窜醉。
注意：
主要是在沒有0的情況下宪萄，多于1個(gè)barcode的hamming distance等于1，或者在沒有0和1的情況下榨惰，多于一個(gè)barcode的hamming distance等于2拜英，都是不應(yīng)該輸出的。
有多個(gè)barcode跟同一個(gè)read的hamming distance都在2以內(nèi)琅催，這個(gè)也分很多種情況居凶，比如，沒有barcode的distance是0藤抡，但是又1個(gè)barcode的distance 是1侠碧，n(n>1)個(gè)barcode的distance是2，這個(gè)時(shí)候最小的distance是1缠黍，且只跟1個(gè)barcode有這個(gè)最小值弄兜，那么就應(yīng)該輸出。另一個(gè)例子瓷式，如果沒有barcode的distance是0挨队，有2個(gè)barcode的distance是1，那就不該輸出蒿往。

計(jì)算hamming 距離
漢明距離是使用在數(shù)據(jù)傳輸差錯(cuò)控制編碼里面的，漢明距離是一個(gè)概念湿弦，它表示兩個(gè)（相同長度）字對應(yīng)位不同的數(shù)量瓤漏，我們以d（x,y）表示兩個(gè)字x,y之間的漢明距離。對兩個(gè)字符串進(jìn)行異或運(yùn)算颊埃，并統(tǒng)計(jì)結(jié)果為1的個(gè)數(shù)蔬充，那么這個(gè)數(shù)就是漢明距離。from 維基百科https://zh.wikipedia.org/wiki/%E6%B1%89%E6%98%8E%E8%B7%9D%E7%A6%BB
直接上代碼：

import gzip
from Bio import SeqIO
import itertools
# 定義函數(shù)計(jì)算hanming distance班利。
def hamming(str1, str2):
  return sum(itertools.imap(str.__ne__, str1, str2))
# 處理fasta文件饥漫，將id與seq存儲(chǔ)為dict
def deal_dna_file(a,):
    dna_dict = {}
    for record in SeqIO.parse(a,"fasta"):
        dna_dict[record.id] = record.seq
    return dna_dict

A_DNA_file = sys.argv[1]
B_fastq_file = gzip.open(sys.argv[2],"r")

# B_fastq_file = open(sys.argv[2],"r")
dna_dict = deal_dna_file(A_DNA_file)
# print dna_dict
# 遍歷fastq文件去處理每行序列
for record in SeqIO.parse(B_fastq_file, "fastq"):
    a = 0
    b = 0
    # 遍歷dict，去判斷hamming距離罗标，分0,1,2三種情況庸队，記錄距離為1积蜻，和2的次數(shù)，根據(jù)次數(shù)去判斷彻消。
    for k,v in dna_dict.items():
        if a >1:
            break
        if hamming(v,record.seq[:len(v)]) == 0:
            print record.seq
        elif hamming(v,record.seq[:len(v)]) == 1:
            a +=1
        elif hamming(v,record.seq[:len(v)]) ==2:
            b +=1
    print a,b
    if a == 1 and b > 1: # one barcode distance is 1 and more than one barcode distance are 2
        print record.seq
    if a == 1 and b == 0: # only barcode distance is 1. 
        print record.seq
    if a == 0 and b == 1: # only barcode distance are 2.
        print record.seq

這樣就符合要求了竿拆，有bug請反饋。
生信學(xué)習(xí)者練習(xí)題宾尚；