編寫Python代碼解析病毒基因組

該畫作描繪了新冠狀病毒進入肺組織的過程。此時病毒被呼吸細胞所分泌的粘液所包繞，並且周圍還有分泌的抗體、小免疫系統的多種蛋白質。病毒由包膜包裝，包膜上含有S蛋白、M蛋白、E蛋白。病毒內部為病毒RNA。回答：病毒（virus）是介於生命和非生命之間的一種物質形式。病毒游離於細胞之外時，沒有複製能力，不表現任何生命特徵，只是一種有機物質的存在的形式；但進入細胞後就表現出生命的特徵，能利用細胞體系進行複製，感染鄰近細胞，是一種非常小的感染性細胞內（分子）專性寄生物。言歸正傳，今天就讓我們來編寫Python程序解析病毒基因組。

首先，下載病毒基因組信息。

登錄 GenBank，網址：

www.ncbi.nlm.nih.gov/genbank/

查詢輸入：MN908947

MN908947 為我國復旦大學上海公共衛生臨床中心和公共衛生學院於今年1月5日上傳的新冠狀病毒全基因組數據。下載fasta文件，另存為：MN908947.fa。

構建讀入基因組函數：

def loadCoronaVirusDNA(filename):
    f = open(filename)
    dna=''
    for line in f:
        if line[0] != '>':
            dna += line.rstrip()
        else:
            header = line.split()
            name = header[0][1:]
    return name, dna

coronavirusname, coronavirusdna = loadCoronaVirusDNA('MN908947.fa')

其次，我們要病毒的部分蛋白質的序列。這就需要我們將DNA/RNA的信息「翻譯」成蛋白質，即三聯密碼子的識別。我們定義一個字典（dictionary），key是不同的三聯密碼子，value是20中胺基酸，且是一個三項元組（tuple），第一項是胺基酸的全稱，第二項是三字符縮寫，第三項是單字符縮寫：

codon = {'ATT': ('Isoleucine', 'Ile', 'I'),
         'ATC': ('Isoleucine', 'Ile', 'I'),
         'ATA': ('Isoleucine', 'Ile', 'I'),

         'CTT': ('Leucine', 'Leu', 'L'),
         'CTC': ('Leucine', 'Leu', 'L'),
         'CTA': ('Leucine', 'Leu', 'L'),
         'CTG': ('Leucine', 'Leu', 'L'),
         'TTA': ('Leucine', 'Leu', 'L'),
         'TTG': ('Leucine', 'Leu', 'L'),

         'GTT': ('Valine', 'Val', 'V'),
         'GTC': ('Valine', 'Val', 'V'),
         'GTA': ('Valine', 'Val', 'V'),
         'GTG': ('Valine', 'Val', 'V'),

         'TTT': ('Phenylalanine', 'Phe', 'F'),
         'TTC': ('Phenylalanine', 'Phe', 'F'),

         'ATG': ('Methionine', 'Met', 'M'),

         'TGT': ('Cysteine ', 'Cys', 'C'),
         'TGC': ('Cysteine ', 'Cys', 'C'),         

         'GCT': ('Alanine', 'Ala', 'A'),
         'GCC': ('Alanine', 'Ala', 'A'),
         'GCA': ('Alanine', 'Ala', 'A'),
         'GCG': ('Alanine', 'Ala', 'A'),

         'GGT': ('Glycine ', 'Gly', 'G'),
         'GGC': ('Glycine ', 'Gly', 'G'),
         'GGA': ('Glycine ', 'Gly', 'G'),
         'GGG': ('Glycine ', 'Gly', 'G'),

         'CCT': ('Proline', 'Pro', 'P'),
         'CCC': ('Proline', 'Pro', 'P'),
         'CCA': ('Proline', 'Pro', 'P'),
         'CCG': ('Proline', 'Pro', 'P'),

         'ACT': ('Threonine', 'Thr', 'T'),
         'ACC': ('Threonine', 'Thr', 'T'),
         'ACA': ('Threonine', 'Thr', 'T'),
         'ACG': ('Threonine', 'Thr', 'T'),

         'TCT': ('Serine', 'Ser', 'S'),
         'TCC': ('Serine', 'Ser', 'S'),
         'TCA': ('Serine', 'Ser', 'S'),
         'TCG': ('Serine', 'Ser', 'S'),
         'AGT': ('Serine', 'Ser', 'S'),
         'AGC': ('Serine', 'Ser', 'S'),

         'TAT': ('Tyrosine', 'Tyr', 'Y'),
         'TAC': ('Tyrosine', 'Tyr', 'Y'),

         'TGG': ('Tryptophan', 'Trp', 'W'),

         'CAA': ('Glutamine', 'Gln', 'Q'),
         'CAG': ('Glutamine', 'Gln', 'Q'),

         'AAT': ('Asparagine', 'Asn', 'N'),
         'AAC': ('Asparagine', 'Asn', 'N'),

         'CAT': ('Histidine ', 'His', 'H'),
         'CAC': ('Histidine ', 'His', 'H'),

         'GAA': ('Glutamic acid', 'Glu', 'E'),
         'GAG': ('Glutamic acid', 'Glu', 'E'),

         'GAT': ('Aspartic acid', 'Asp', 'D'),
         'GAC': ('Aspartic acid', 'Asp', 'D'),

         'AAA': ('Lysine', 'Lys', 'K'),
         'AAG': ('Lysine', 'Lys', 'K'),

         'CGT': ('Arginine', 'Arg', 'R'),
         'CGC': ('Arginine', 'Arg', 'R'),
         'CGA': ('Arginine', 'Arg', 'R'),
         'CGG': ('Arginine', 'Arg', 'R'),
         'AGA': ('Arginine', 'Arg', 'R'),
         'AGG': ('Arginine', 'Arg', 'R'),

         'TAA': ('Stop', 'Stop', '___'),
         'TAG': ('Stop', 'Stop', '___'),
         'TGA': ('Stop', 'Stop', '___')
         }

構建翻譯函數：

def translate(dna):
    aminoacid = ''
    amino = ''
    a = ''
    i = 0
    dna_len = len(dna)
    DNA = dna.upper()

    while(1):
        if i+3 <= dna_len:
            c = DNA[i:i+3]
            residue = codon[c][0]
            res = codon[c][1]
            r = codon[c][2]

            if i == 0:
                aminoacid = residue
                amino = res
                a = a + r

            else:
                aminoacid = aminoacid + '-' + residue
                amino = amino + '-' + res
                a = a + r

            i = i + 3
        else:
            break
    return aminoacid, amino, a

這時，讓我們輸出S蛋白的序列。在Genbank中，S蛋白的基因是：

說明定位是 21563 到 25384,。注意，在Python中，由於 slice 是從 0 開始計數，所以引用 index 則要減 1，則為 21562 到 25383：

spike_dna = coronavirusdna[21562:25383]
spike_protein_tuple = translate(spike_dna)
spike_protein = spike_protein_tuple[2]

最後，新冠病毒 S 蛋白的空間結構也於2月26日在Science發表了（DOI: 10.1126/science.abb2507），PDB 號 6vsb。卡通圖如下：

應用BioPython解析和可視化蛋白質的結構

傳染病模型的Python實現

如何安裝和使用免費的蛋白質空間結構可視化軟體 PyMol

關注Chuan Workstation

獲得更多臨床醫學和

生物信息學新技術和新進展