生物信息學(xué)方法預(yù)測病原體抗原蛋白序列中多肽疫苗候選表位

趙靜靜,澹小秀,王廣志,歐陽儉,簡星星,3*,謝 鷺*

（1.上海海洋大學(xué)食品學(xué)院,中國上海 201306;2.上海市生物醫(yī)藥技術(shù)研究院上海生物信息技術(shù)研究中心,中國上海 201203;3.中南大學(xué)湘雅醫(yī)院生物信息中心,中國湖南 長沙 410008）

疫苗是當(dāng)前人類對抗病原生物最經(jīng)濟(jì)有效的方式。現(xiàn)階段,針對病原體蛋白中能夠引起免疫系統(tǒng)效應(yīng)功能的表位而設(shè)計的多肽疫苗是新型疫苗研發(fā)的熱點(diǎn)。多肽疫苗一般是指利用化學(xué)合成技術(shù),按照病原體抗原基因中具有免疫原性的一段或多段氨基酸序列構(gòu)成而制備的疫苗[1]。其因成分簡單、特異性強(qiáng)、安全性高、易于保藏等優(yōu)點(diǎn)被廣泛應(yīng)用于傳染病和腫瘤的預(yù)防和治療中,例如:預(yù)防人乳頭瘤狀病毒的疫苗[2]、預(yù)防乙型肝炎病毒的HBsAg/preS疫苗[3]、新型艾滋病疫苗RV144[4]及個性化腫瘤新抗原疫苗[5]等。

近年來,得益于生物信息學(xué)技術(shù)的迅猛發(fā)展,反向疫苗學(xué)在疫苗設(shè)計領(lǐng)域顯示出廣闊的應(yīng)用空間。以“序列-結(jié)構(gòu)-功能”思想為依據(jù),借助生物信息學(xué)工具對病原體高通量的組學(xué)信息（基因組學(xué)、轉(zhuǎn)錄組學(xué)、蛋白質(zhì)組學(xué)等）預(yù)先進(jìn)行嚴(yán)格的預(yù)測和篩選,能極大地提高研究人員發(fā)現(xiàn)病原體蛋白質(zhì)序列中候選表位的效率。隨著候選表位的持續(xù)產(chǎn)出和收集,相關(guān)表位數(shù)據(jù)庫應(yīng)運(yùn)而生,如IEDB（Immune Epitope Database）[6]及本課題組自行構(gòu)建的dbPepNeo[7]等。研究人員利用這些數(shù)據(jù)集開發(fā)了諸多表位預(yù)測和篩選的工具,但其中的多數(shù)均只針對特定的功能進(jìn)行開發(fā),如:NetMHC-pan EL 4.0[8]主要用于預(yù)測與主要組織相容性復(fù)合體（major histocompatibility complex,MHC）結(jié)合的抗原表位;本課題組開發(fā)的INeo-Epp[9]則用于預(yù)測具有免疫原性的抗原表位。因此,表位預(yù)測工具功能的多樣化,對于候選表位的鑒定以及多肽疫苗的設(shè)計具有重要指導(dǎo)意義[10～12]。

本研究基于反向疫苗學(xué)策略,綜合利用多種生物信息學(xué)工具,構(gòu)建了一個完整的、易于操作的表位預(yù)測和篩選流程,并以當(dāng)下正在全球蔓延的嚴(yán)重急性呼吸綜合征冠狀病毒2（severe acute respiratory syndrome coronavirus 2,SARS-CoV-2）為實例,成功地檢驗了該流程的實用性。文中共篩選到34條關(guān)于SARS-CoV-2的T細(xì)胞候選表位,其中20條候選表位與數(shù)據(jù)庫中經(jīng)過驗證的表位高度同源,且能夠被T細(xì)胞受體（T cell receptor,TCR）特異性識別。此外,為了驗證流程的適用性,我們還對鼠類肉瘤病毒癌基因KRAS與大腸桿菌基因OmpC所編碼蛋白質(zhì)的抗原表位進(jìn)行了示例性分析。本研究不僅為新冠疫苗設(shè)計提供了候選表位,而且還建立了一個廣泛適用于多肽疫苗候選表位預(yù)測和篩選的操作流程。

1 材料和方法

1.1 數(shù)據(jù)下載

本研究涉及的SARS-CoV-2的結(jié)構(gòu)蛋白包括刺突糖蛋白（spike glycoprotein,S）、膜糖蛋白（membrane glycoprotein,M）、包膜蛋白（envelope protein,E）與核衣殼蛋白（nucleocapsid phosphoprotein,N）,它們的蛋白質(zhì)序列與鼠類肉瘤病毒癌基因KRAS及大腸桿菌基因OmpC相關(guān)的蛋白質(zhì)序列文件均從NCBI（https://www.ncbi.nlm.nih.gov/）數(shù)據(jù)庫獲取,檢索ID分別為QHR63290.2、QHD43419.1、QHD-43418.1、QHD43423.2、CAC5395073.1 和 ADU34-074.2。

1.2 表位預(yù)測識別

本研究使用生物信息學(xué)工具NetMHCpan EL 4.0[8]（http://tools.immuneepitope.org/mhci/）,對SARSCoV-2的S蛋白、M蛋白、E蛋白與N蛋白及KRAS和OmpC蛋白序列中潛在的T細(xì)胞表位分別進(jìn)行預(yù)測。然后,依據(jù)“抗原肽-MHC”復(fù)合物間的親和力（通常,%rank<0.5被定義為強(qiáng)結(jié)合,%rank<2為弱結(jié)合,其他為不結(jié)合）進(jìn)行篩選,保留%rank在0～2的短肽作為結(jié)合表位,進(jìn)一步評估其作為多肽疫苗候選表位的潛力。

1.3 表位免疫原性與抗原性預(yù)測

為了篩選出具有免疫原性的候選表位,我們使用IEDB中提供的免疫原性預(yù)測工具[6,13]（http://tools.iedb.org/immunogenicity/）及本課題組自行開發(fā)的INeo-Epp[9]（http://www.biostatistics.online/INeo-Epp/）兩個生物信息學(xué)工具,對S蛋白、M蛋白、E蛋白與N蛋白及KRAS和OmpC蛋白結(jié)合表位的免疫原性分別進(jìn)行預(yù)測。IEDB與INeo-Epp工具是通過不同的算法對“抗原肽-MHC”復(fù)合物錨定位置處氨基酸的理化特性進(jìn)行建模,預(yù)測效果表現(xiàn)良好[前者的曲線下面積（area under the curve,AUC）=0.69,后者AUC=0.81]。對于免疫原性評分（score）,一般認(rèn)為評分越高表示該表位引發(fā)免疫反應(yīng)的可能性越大[14]。為了降低假陽性率,通常以0.2為標(biāo)準(zhǔn)進(jìn)行篩選。因此,本研究僅保留IEDB中score大于0.2和INeo-Epp中輸出結(jié)果為強(qiáng)陽性（positive-high,PH）的表位。

隨后,我們使用Vaxijen v2.0[15]（http://www.ddgpharmfac.net/vaxijen/）工具預(yù)測結(jié)合表位的抗原性。Vaxijen是第一個利用非序列比對策略預(yù)測表位抗原性的生物信息學(xué)工具,也是反向疫苗學(xué)中常用的工具。該模型設(shè)置了0～1的閾值,以預(yù)測病毒、細(xì)菌和腫瘤蛋白質(zhì)序列中的保護(hù)性抗原。當(dāng)閾值為0.4時,該模型在病毒測試數(shù)據(jù)集的效果較好（AUC=0.74）。因此,本研究在0.4的閾值下獲取S蛋白、M蛋白、E蛋白、N蛋白、KRAS蛋白及OmpC蛋白結(jié)合表位的抗原性分值。

1.4 表位理化特性分析

為了降低多肽疫苗輸注后患者發(fā)生不良反應(yīng)的概率,本研究使用AllerTOP v2.0[16]（https://www.ddg-pharmfac.net/AllerTOP/）、ClanTox[17]（http://www.clantox.cs.huji.ac.il/）工具對結(jié)合表位的理化特性進(jìn)行分析。AllerTOP和ClanTox擁有非常簡潔的操作界面,不需要額外設(shè)置閾值即可得到預(yù)測結(jié)果。此外,AllerTOP還是當(dāng)前較被認(rèn)可的預(yù)測蛋白質(zhì)致敏性的工具,它通過自協(xié)方差和交叉協(xié)方差變換的方法對氨基酸進(jìn)行編碼,綜合利用邏輯回歸、決策樹、隨機(jī)森林、樸素貝葉斯、多層感知器、K-近鄰等算法進(jìn)行建模,最后選擇結(jié)果最好的K-近鄰算法（AUC=0.85）完成工具的構(gòu)建。而ClanTox是一種蛋白質(zhì)毒素分類器,可以根據(jù)蛋白質(zhì)的一級結(jié)構(gòu)計算出該蛋白質(zhì)是否為毒性蛋白質(zhì)。在本研究中,我們僅保留AllerTOP和ClanTox輸出結(jié)果為非致敏（non-allergen）與非毒性（nontoxin）的表位作為SARS-CoV-2多肽疫苗研發(fā)的候選表位。

1.5 人群覆蓋率評估

根據(jù)先前的研究[18～19],我們構(gòu)建了12種在中國人群中常見的人類白細(xì)胞抗原（human leukocyte antigen,HLA）的等位基因（allele）分型。它們分別為HLA-A*01:01、HLA-A*02:01、HLA-A*03:01、HLA-A*11:01、HLA-A*23:01、HLA-A*24:02、HLA-B*07:02、HLA-B*08:01、HLA-B*35:01、HLA-B*40:01、HLA-B*44:02、HLA-B*44:03（表1）。此外,我們使用IEDB提供的人口覆蓋率計算工具Population Coverage[20]（http://tools.immuneepitope.org/population/）評估了每條候選表位在中國及世界人群中的覆蓋率。需要注意的是,使用Population Coverage工具評估表位在特定人群中的覆蓋率時,需要提交一份包含該表位所對應(yīng)HLA等位基因分布頻率的信息文件。

表1 12種HLA等位基因的分布頻率Table 1 Distribution frequency of 12 HLA alleles

1.6 同源比對及特異性TCR檢索

IEDB[21]（http://www.iedb.org/）、VDJdb[22]（https://vdjdb.cdr3.net/）、McAPS-TCR[23]（http://friedmanlab.weizmann.ac.il/McPAS-TCR/）是當(dāng)前表位及相關(guān)TCR數(shù)據(jù)儲存最多的數(shù)據(jù)庫。首先,本研究收集和整理了3個數(shù)據(jù)庫中有關(guān)human、CD8+T、TCRβ的表位信息及其相關(guān)的TCR序列信息。然后,應(yīng)用BLAST（Basic Local Alignment Search Tool）[24]（https://blast.ncbi.nlm.nih.gov/Blast.cgi/）工具,將本研究篩選的候選表位與數(shù)據(jù)庫中的表位進(jìn)行同源比對,記錄比對結(jié)果中所有同源表位的序列信息及E value和identity值。最后,在已整理的TCR數(shù)據(jù)中手動檢索與該同源表位匹配的TCR序列信息。

2 結(jié)果

2.1 SARS-CoV-2結(jié)構(gòu)蛋白中T細(xì)胞識別的候選表位

本研究綜合多個生物信息學(xué)工具（圖1A）,系統(tǒng)地篩選了SARS-CoV-2結(jié)構(gòu)蛋白序列中可激發(fā)T細(xì)胞免疫應(yīng)答的表位,即T細(xì)胞表位。首先,我們構(gòu)建了12種在中國人群普遍存在的HLA等位基因分型:HLA-A*01:01、HLA-A*02:01、HLAA*03:01、HLA-A*11:01、HLA-A*23:01、HLA-A*24:02、HLA-B*07:02、HLA-B*08:01、HLAB*35:01、HLA-B*40:01、HLA-B*44:02、HLAB*44:03;然后,預(yù)測了SARS-CoV-2結(jié)構(gòu)蛋白中長度為8～12個氨基酸的T細(xì)胞表位。通過“抗原肽-MHC”結(jié)合預(yù)測工具NetMHCpan EL,共得到117 180條預(yù)測表位。之后,根據(jù)%rank分值篩選出0～2的結(jié)合表位,僅有2 593條表位顯示能與MHC分子結(jié)合。

免疫原性是指抗原誘導(dǎo)機(jī)體發(fā)生特異性免疫應(yīng)答的性能,抗原性則是指抗原與其所誘導(dǎo)產(chǎn)生的抗體或致敏淋巴細(xì)胞特異性結(jié)合的能力。簡而言之,表位的免疫原性與抗原性越強(qiáng),表明其被T細(xì)胞捕獲并引發(fā)免疫應(yīng)答的潛力越高?？紤]到表位的免疫原性和抗原性均是決定該抗原能否被T淋巴細(xì)胞識別的關(guān)鍵特征,我們進(jìn)一步對結(jié)合表位的免疫原性及抗原性進(jìn)行了分析。首先,使用IEDB和INeo-Epp兩個工具分別對2 593條SARS-CoV-2結(jié)合表位的免疫原性進(jìn)行預(yù)測,保留在兩個工具中結(jié)果都較好的表位（IEDB,score>0.2;INeo-Epp,score>0.5即為PH）,共獲得131條;隨后,使用Vaxijen工具分析這131條表位的抗原性,最終僅獲得71條既有免疫原性又有抗原性的表位。

疫苗設(shè)計的前提是對該病毒候選表位有充分的理解。為此,我們使用蛋白質(zhì)理化特性分析工具AllerTOP和ClanTox評估了這71條表位的致敏性與毒性,發(fā)現(xiàn)其中非致敏和非毒性的候選表位有49條。此外,由于HLA等位基因的多態(tài)性,我們發(fā)現(xiàn)候選表位中存在同一表位與多個HLA等位基因結(jié)合的情況,例如:表位M135～144和表位S1110～1121以不同的%rank 分值與 HLA-B*40:01、HLA-B*44:02及HLA-B44:03三種分型結(jié)合。候選表位與HLA等位基因分型的結(jié)合率如圖1B所示。

最后,在忽略同一表位與不同HLA等位基因結(jié)合的情況下,我們?nèi)栽?9條表位中找到了34條獨(dú)特的表位（表2）,其中15條來自S蛋白序列、9條來自M蛋白序列、5條來自E蛋白序列、5條來自N蛋白序列（圖1C）?？傮w來看,這34條候選表位經(jīng)過了系統(tǒng)和嚴(yán)格的篩選,可以與MHC分子結(jié)合,并由抗原提呈細(xì)胞加工后呈遞到細(xì)胞表面,供TCR識別,繼而誘導(dǎo)T細(xì)胞發(fā)揮效應(yīng)功能。

表2 SARS-CoV-2的候選表位信息Table 2 Information about candidate epitopes of SARS-CoV-2

圖1 多肽疫苗候選表位預(yù)測流程及其在SARS-CoV-2中的應(yīng)用（A）SARS-CoV-2 T細(xì)胞候選表位的預(yù)測流程圖;（B）候選表位與HLA等位基因的結(jié)合率;（C）不同預(yù)測過程中得到的SARSCoV-2表位數(shù)量。Fig.1 The workflow for predicting candidate epitopes for peptide vaccines and its application in SARS-CoV-2（A）Screening process of SARS-CoV-2 T cell candidate epitopes;（B）Binding rates of candidate epitopes and HLA alleles;（C）The number of SARS-CoV-2 epitopes retained during the screening process.

2.2 候選表位的人群覆蓋率

由于與表位結(jié)合的HLA等位基因在不同國家的分布頻率差異明顯,我們評估了SARS-CoV-2的34條候選表位在中國及世界范圍的覆蓋率。首先,我們獲取了HLA等位基因在中國及世界范圍內(nèi)的分布頻率[18～19]。然后,基于HLA等位基因的分布頻率,使用Population Coverage工具計算了每條候選表位在中國和世界人群的覆蓋率。結(jié)果顯示,本研究篩選的SARS-CoV-2候選表位在中國人群廣泛覆蓋（1%～91%）（表2）。其中,表位S1110～1121與表位 M135～144在中國人群的適用范圍高達(dá)91%,同時在世界人群的適用范圍也達(dá)到71%。

2.3 同源表位及其誘導(dǎo)的TCR序列

為了檢驗本研究所構(gòu)建預(yù)測流程的實用性,我們使用BLAST工具將本文篩選的34條SARSCoV-2表位與數(shù)據(jù)庫中經(jīng)過驗證的表位進(jìn)行同源比對。結(jié)果顯示,有20條SARS-CoV-2候選表位與 IEDB、VDJdb、McAPS-TCR 3個數(shù)據(jù)庫中經(jīng)過驗證的表位同源（表3）。另外,當(dāng)比對結(jié)果的參數(shù)E value調(diào)整為 0.05 時,E29～38和 E30～37這兩條表位與數(shù)據(jù)庫中LLAILTYYV（blast epitope）表位的identity高達(dá)100%。與此同時,我們還在數(shù)據(jù)庫中檢索到6條與該同源表位相關(guān)的特異性TCR序列（CASSLVRDRHTEAFF、CASSPTGTGGSDTQYF、CASSQAGEQYF、CASSWVGGADTQYF、CASTVRQGSNQPQHF 和 CSASFHNGFWGGTEAFF）。

表3 SARS-CoV-2的同源表位及相關(guān)TCR序列信息Table 3 Information about homologous epitopes and relevant TCR sequences of SARS-CoV-2

2.4 其他病原體蛋白質(zhì)序列中候選表位的預(yù)測與篩選

為了進(jìn)一步驗證該流程的適用性,我們在NCBI中檢索了鼠類肉瘤病毒癌基因KRAS及大腸桿菌基因OmpC相關(guān)的蛋白質(zhì)序列,并對蛋白質(zhì)序列中適合多肽疫苗研發(fā)的表位進(jìn)行了預(yù)測和篩選。從表4可知,本研究共預(yù)測到5條有關(guān)KRAS的表位,即 KLVVVGAGGVGK、VVVGAGGVGK、LVVVGAGGVGK、QYMRTGEGFL和YMRTGEGF。基于這5條表位的%rank分值、免疫原性、抗原性、致敏性及毒性,我們發(fā)現(xiàn)7～16區(qū)段的表位VVVGAGGVGK有潛力成為KRAS多肽疫苗研發(fā)的候選表位。另外,從表5的同源比對結(jié)果可以看出,VVVGAGGVGK與數(shù)據(jù)庫中YMDDVVLGA表位的identity達(dá)到80%,同時還有一條特異性TCR序列CASSYLTGEGDYGYTF與該同源表位相關(guān),這進(jìn)一步表明表位VVVGAGGVGK具有免疫原性,能夠誘導(dǎo)T細(xì)胞免疫應(yīng)答。

表4 KRAS的候選表位信息Table 4 Information about candidate epitopes of KRAS

表5 KRAS的同源表位及相關(guān)TCR序列信息Table 5 Information about homologous epitopes and relevant TCR sequences of KRAS

我們在基因OmpC相關(guān)的蛋白質(zhì)序列中共預(yù)測到12條表位（表6）。從預(yù)測表位的各預(yù)測分值及同源比對結(jié)果可知,表位YEGFGIGGAI（205～214）、VLPEFGGDTY（127～136）、TDVLPEFGGDTY（125～136）、DVLPEFGGDTY（126～136）可以作為OmpC多肽疫苗研發(fā)的候選表位（表6～7）。而且,對比分析 125～136、126～136、127～136 區(qū)段可以發(fā)現(xiàn),125～136區(qū)段是OmpC蛋白中適用于多肽疫苗研發(fā)的優(yōu)勢區(qū)段。由此,我們認(rèn)為本研究所提出的預(yù)測流程可適用于病原體抗原蛋白序列中多肽疫苗候選表位的預(yù)測和篩選。

表6 OmpC的候選表位信息Table 6 Information about candidate epitopes of OmpC

表7 OmpC的同源表位及相關(guān)TCR序列信息Table 7 Information about homologous epitopes and relevant TCR sequences of OmpC

3 討論

快速、安全、靈活、高效且低成本的新型疫苗研發(fā)路線是未來人類對抗病原體最有效的手段[25]。近年來,隨著生物信息學(xué)技術(shù)的發(fā)展和應(yīng)用,研究人員已經(jīng)可以從病原體高通量組學(xué)信息中挖掘出大量的免疫原性表位序列,同時實現(xiàn)對T/B細(xì)胞靶點(diǎn)處表位的精準(zhǔn)識別,從而為新型疫苗的研發(fā)提供理論指導(dǎo)[26]。基于此,借助生物信息學(xué)工具初篩病原體蛋白質(zhì)序列中的候選表位已經(jīng)成為疫苗研究中關(guān)鍵的一步。在此次新型冠狀病毒肺炎（COVID-19）發(fā)生之初,Grifoni團(tuán)隊[19]就使用Net-MHCpan EL web工具對SARS-CoV-2蛋白質(zhì)序列中潛在的T、B細(xì)胞表位進(jìn)行了預(yù)測,經(jīng)同源比對發(fā)現(xiàn),這些表位與IEDB中記錄的SARS表位有較高的同源性（B細(xì)胞表位為69%～100%;T細(xì)胞表位為17%～100%）。由此,我們綜合利用當(dāng)前的生物信息學(xué)工具,對SARS-CoV-2結(jié)構(gòu)蛋白中潛在的T細(xì)胞候選表位進(jìn)行了預(yù)測和篩選,共獲得34條符合MHCⅠ類分子的T細(xì)胞候選表位信息,其中候選表位 S1110～1121和 M135～144覆蓋高達(dá)91%的中國人群和71%的世界人群,可為COVID-19的疫苗設(shè)計提供參考。

目前,利用機(jī)器學(xué)習(xí)算法,相關(guān)人員已經(jīng)開發(fā)了系列針對抗原提呈過程的生物信息學(xué)工具。比如,利用人工神經(jīng)網(wǎng)絡(luò)算法開發(fā)的NetMHCpan系列工具,既擁有操作簡潔的web頁面,也有適用于Linux系統(tǒng)的版本。另外,本課題組利用隨機(jī)森林算法開發(fā)的INeo-Epp工具可以預(yù)測抗原肽的免疫原性,其AUC達(dá)到0.81,而IEDB中提供的同類型工具的AUC僅為0.69。除本研究中涉及的工具外,還有諸多類似功能的工具可用于疫苗候選表位的預(yù)測和篩選,例如:NetCTL 1.2（http://www.cbs.dtu.dk/services/NetCTL/）可用于預(yù)測蛋白質(zhì)序列中潛在的細(xì)胞毒性T淋巴細(xì)胞（cytotoxic T lymphocyte,CTL）表位[27];MixMHC2pred[28]可用于預(yù)測HLAⅡ類分子呈遞的抗原肽。這些生物信息學(xué)預(yù)測工具在抗原特定功能的預(yù)測中均顯示出良好的性能,但在實現(xiàn)疫苗候選表位的完整預(yù)測和分析方面仍存在一定難度。

綜上所述,本研究綜合利用多種生物信息學(xué)工具,結(jié)合自主研發(fā)的算法工具,整理并提出了一個可廣泛應(yīng)用于病原體多肽疫苗候選表位預(yù)測和分析的流程。通過該工作流程可以預(yù)測出不同抗原蛋白質(zhì)序列中高可信度的疫苗候選表位,從而為預(yù)防和控制感染性疾病與惡性腫瘤的新型疫苗的研制提供參考。