新型冠狀病毒Omicron變異株密碼子偏愛性及其進化分析

李亞飛，羅春雨，石哲芳，劉 奇

2019年12月，由新型冠狀病毒(Severe Acute Respiratory Syndrome Coronavirus 2，SARS-CoV-2)引起的新型冠狀病毒肺炎(Corona Virus Disease 2019，COVID-19)暴發(fā)，截至2022年2月，SARS-CoV-2大流行已在全球造成超過4億人感染和超過600萬人死亡，對全球健康和經(jīng)濟發(fā)展構(gòu)成了嚴重威脅。SARS-CoV-2為單股正鏈RNA病毒，主要編碼核衣殼蛋白(nucleocapsid protein，N蛋白)、刺突糖蛋白(spike glycoprotein，S蛋白)、膜蛋白(membrane protein，M蛋白)和高度疏水的包膜蛋白(envelope protein，E蛋白)4種結(jié)構(gòu)蛋白和ORF1ab、ORF1a、ORF3a、ORF6、ORF7a、ORF7b、ORF8、ORF10共8種非結(jié)構(gòu)蛋白(nonstructural protein，NSP)[1]。

冠狀病毒的突變率高，復(fù)制容易出錯[2]。隨著感染人口基數(shù)的增多，SARS-CoV-2變體也隨之出現(xiàn)。世界衛(wèi)生組織(World Health Organization，WHO)將能引起大范圍流行的變異株歸類為關(guān)切變異株(Variants of Concern，VOC)[3]，目前主要包括Alpha 變體(B.1.1.7 首次發(fā)現(xiàn)于英國)，Beta變體(B.1.351 首次發(fā)現(xiàn)于南非)，Gamma變體(P.1 首次發(fā)現(xiàn)于巴西)、 Delta變體(B.1.617.2 首次發(fā)現(xiàn)于印度)及2021年11月發(fā)現(xiàn)的Omicron最新變異體(B.1.1.529/BA.1 首次發(fā)現(xiàn)于南非)[4]。Omicron自發(fā)現(xiàn)以來傳播迅速，目前在非洲已經(jīng)成為SARS-CoV-2感染的主導(dǎo)病毒株[5]。

不同種生物在蛋白質(zhì)編碼過程中同義密碼子的使用頻率是非均衡的，存在特定傾向使用的現(xiàn)象，即其具有密碼子使用的偏愛性[6]，與其基因表達或者其他選擇壓力如自然選擇或者突變壓力等有關(guān)。研究證明，病毒關(guān)鍵蛋白編碼基因中特異性的密碼子使用偏倚與其復(fù)制能力、適應(yīng)性免疫逃避能力和致病性等都有著密切關(guān)系[7]。目前Omicron的來源和現(xiàn)有疫苗抗性等尚不明確，進行Omicron與其他VOC毒株密碼子使用偏愛性分析可為探究其復(fù)制能力、適應(yīng)性免疫逃避能力和致病性等方面的差異研究提供分子機制基礎(chǔ)。因此，我們分析了Omicron不同蛋白編碼基因密碼子使用特性的改變及影響因素，并與標準株Wuhan-Hu-1及其他4種VOC(Alpha、Beta、Gamma、Delta)進行比較分析，研究其密碼子偏愛性改變及方向性變化并探索之間的進化關(guān)系。

1 材料與方法

1.1 材 料

1.1.1 目的序列來源 收集截至2022年2月底NCBI SARS-CoV-2數(shù)據(jù)庫(https://www.ncbi.nlm.nih.gov/sars-cov-2/)中VOC所有高質(zhì)量完整序列：Alpha(B.1.1.7)65 143條、Beta(B.1.351)295條、Gamma(P.1)4 010條、 Delta(B.1.617.2)6 018條、Omicron 34 396條(B.1.1.529/BA.1)及SARS-CoV-2標準株Wuhan-Hu-1(NC_45512.2)序列。

1.1.2 計算分析軟件 EMBOSS子程序CUSP(https://www.bioinformatics.nl/cgi-bin/emboss/cusp)用于計算目的基因密碼子GC3值和各蛋白密碼子使用頻率(Frequency)等數(shù)據(jù)；Codon W軟件計算分析有效密碼子數(shù)和同義密碼子相對使用頻率等；Sigmaplot 14.0軟件繪圖分析；MEGA 11.0軟件構(gòu)建系統(tǒng)進化樹；SPSS 22.0軟件進行蛋白聚類分析等。

1.2 方 法

1.2.1 有效密碼子數(shù)分析 利用Codon W軟件計算目的病毒株各蛋白的有效密碼子數(shù)(Effective Number of Codon，ENC)。ENC取值范圍為20～61，用于描述密碼子使用偏離隨機選擇的程度，能反應(yīng)同義密碼子非均衡使用的偏好程度[8]。ENC值越小使用偏性越強(每個氨基酸只使用一個密碼子時值為20；各密碼子使用概率相等時值為61)。ENC<35說明具有較強烈的密碼子使用偏性，ENC>55表明使用偏性弱。

1.2.2 相對同義密碼子使用度 相對同義密碼子使用度(Relative Synonymous Codon Usage，RSCU)表示某一特定的密碼子在編碼對應(yīng)氨基酸的同義密碼子間的相對概率[9]。若某一密碼子的RSCU=1，表示沒有使用偏性；RSCU>1，表示該密碼子為使用相對較多的密碼子；RSCU>1.5，說明該密碼子為使用高頻密碼子。

1.2.3 ENC-plot分析 ENC-plot是以ENC值為縱坐標，GC3(密碼子第3位為G或C的密碼子含量)為橫坐標繪制圖形，用于反應(yīng)突變壓力對密碼子使用偏倚的影響程度[10]。標準曲線顯示的是密碼子使用偏性僅受突變壓力影響的情況。當(dāng)對應(yīng)點落在標準曲線附近時，提示突變是影響觀察到的密碼子偏性的主要因素；遠離曲線時，提示密碼子偏倚還受自然選擇等因素的影響。

1.2.4 中性分析 以GC3為橫坐標，GC12為縱坐標作中性繪圖分析(Neutrality Plot)和回歸曲線分析，二者的相關(guān)性分析結(jié)果可反映突變壓力和自然選擇對密碼子使用模式的相對影響[11]。當(dāng)回歸系數(shù)為1或接近1時，則該蛋白編碼基因GC3和GC12相關(guān)性顯著，其密碼子3個位置上的堿基組成無明顯差異，突變壓力為影響密碼子偏性的主要因素。若回歸系數(shù)為0或者遠離1，則說明自然選擇為影響密碼子使用的主要因素。

1.2.5 奇偶偏好分析 奇偶規(guī)則(Parity Rule 2，PR2)用于評估自然環(huán)境選擇或突變壓力對基因密碼子使用模式的影響。若2個互補鏈DNA之間不存在任何突變或選擇效應(yīng)上的偏倚，則預(yù)測其密碼子A=T、C=G[12]。此規(guī)則以目的基因的GC偏向值[G3/(G3+C3)]為橫坐標，AT偏向值[A3/(A3+T3)]為縱坐標進行PR2-plot繪圖。中心坐標點(0.5，0.5)代表A=T、C=G(PR2)的位置，表示此處目的基因內(nèi)鏈密碼子使用沒有偏差，不受突變和翻譯選擇影響。若點在四個象限分布均勻，則表示密碼子使用偏好由突變引起，分布不均則表示還受自然選擇等其他因素影響。

1.2.6 基于密碼子偏性的聚類分析和發(fā)育樹分析

利用MEGA 11.0軟件對比分析6種SARS-CoV-2病毒株全基因組序列，采用最大似然法(Maximum Likelihood Method，MLM)構(gòu)建其系統(tǒng)進化樹，描述其之間的相關(guān)關(guān)系和可能的進化過程[13]。采用EMBOSS 子程序CUSP計算不同病毒株各蛋白CDS基因序列的密碼子使用頻數(shù)，再利用SPSS 22.0軟件系統(tǒng)進行聚類分析。

2 結(jié) 果

2.1 有效密碼子數(shù)(ENC)分析 Wuhan-Hu-1及5種VOC(Alpha、Beta、Gamma、Delta、Omicron)全基因組結(jié)果(圖1)顯示，各蛋白ENC分布范圍為26.38～57.89。Omicron蛋白ENC均值為46.05±7.80，較Wuhan-Hu-1(46.02±8.03)略高，其總體密碼子偏性較弱。進一步對12種蛋白編碼基因的ENC值分別進行分析發(fā)現(xiàn)，Omicron的S(t=6 327.577,P<0.01)、E(t=17 358.633,P<0.01)、ORF1ab(t=3 859.225,P<0.01)、ORF1a(t=37 046.417,P<0.01) 4種蛋白的ENC較Wuhan-Hu-1出現(xiàn)了增高，表明其密碼子偏性減弱； M(t=-10 313.811,P<0.01)、N(t=-2 209.230,P<0.01)、ORF3a(t=-4487.525,P<0.01)、ORF7b(t=-24 736.158,P<0.01) 4種蛋白的ENC較Wuhan-Hu-1密碼子偏性增強。對比其他4種VOC，Omicron變異株 S、ORF1ab、ORF1a蛋白ENC值最大，M、ORF3a、ORF7b蛋白ENC值最小。

圖1 SARS-CoV-2 及其VOC各蛋白ENC分布Fig.1 ENC values of SARS-CoV-2 and VOC

2.2 RSCU分析 為進一步探索Omicron各蛋白基因組密碼子使用情況的具體變化情況，參考Gun Li[14]等人計算方法，分析比較了Omicron和Wuhan-Hu-1各蛋白RSCU與人基因組RSCU 差異的變化(結(jié)果見表1)。與Wuhan-Hu-1相比，Omicron RSCU更接近人基因組密碼子的，以“+”表示；更遠離人基因組的以“-”表示。結(jié)果顯示Omicron病毒株的S(靠近密碼子28個；遠離18個)、ORF1ab(靠近密碼子14個；遠離9個)、ORF1a(靠近密碼子15個；遠離13個)等主要蛋白均有更多數(shù)量的密碼子與人基因組密碼子偏愛性趨于一致，尤其以S蛋白最為明顯。S和ORF1ab編碼基因覆蓋了整個基因組的近4/5，說明Omicron密碼子使用模式相較Wuhan-Hu-1更靠近人類基因組。

表1 Omicron突變株RSCU與人基因組RSCU的關(guān)系Tab.1 Relationship of RSCU between Omicron variant and human genome

2.3 ENC-plot分析 為了探究Omicron與Wuhan-Hu-1及其他VOC病毒株密碼子使用偏性改變的影響因素，對其12種蛋白編碼基因CDS進行ENC-plot繪圖分析?？傮w結(jié)果(圖2A)顯示Omicron總體相較其他VOC和Wuhan-Hu-1更接近標準曲線，說明Omicron的總體密碼子使用偏倚改變受突變壓力的影響更多。各蛋白分析結(jié)果可以看出，Omicron的 S蛋白和E蛋白相較于其他VOC更接近于標準曲線，表明這兩種蛋白較其他VOC受突變壓力的影響更大。而M蛋白和ORF3a相較其他VOC距離標準曲線更遠，表明相較其他VOC，這些蛋白受突變壓力的影響小。

注：A為SARS-CoV-2 Wuhan-Hu-1及其5種VOC 12種蛋白CDS總圖；B、C、D、E分別S、E、M、ORF3a蛋白ENC-plot圖。圖2 Omicron 和其他VOC病毒株ENC-Plot分析Fig.2 ENC-Plot analysis of Omicron variant and other VOC

2.4 中性繪圖分析 利用中性分析探究5種VOC之間密碼子使用偏倚差異的影響因素，各蛋白中性繪圖分析結(jié)果(圖3)顯示S、E、ORF1ab、ORF1a、ORF3a、ORF6、ORF7b等大多數(shù)蛋白斜率b均遠離1，回歸模型均具有統(tǒng)計學(xué)意義。GC3和GC12之間的弱相關(guān)性表明，這幾種蛋白之間密碼子使用偏倚之間的差異則更多受自然選擇壓力的影響，突變壓力起次要作用。該結(jié)果提示，盡管冠狀病毒基因具有高突變性，但在SARS-CoV-2不斷變異進化過程中，各蛋白經(jīng)受的選擇壓力是不一致的，自然選擇作用導(dǎo)致了不同VOC之間的S、E、ORF1ab等蛋白密碼子使用偏愛性差異的變化，這可能與世界范圍的疫苗接種有關(guān)。

圖3 Omicron 和其他VOC病毒株部分蛋白中性分析Fig.3 Neutral plot analysis of partial protein of Omicron and other VOC

2.5 奇偶偏好分析 PR2分析結(jié)果顯示(圖4)，Omicron與SARS-CoV-2 Wuhan-Hu-1及其他VOC各CDS主要集中在左下象限，密碼子使用不平等。多數(shù)點位于y<0.5以下，提示各VOC基因組第3位密碼子C/T的使用頻率仍舊較高，與前期對SARS-CoV-2的研究[16]一致。該結(jié)果還進一步證明了除突變壓力外，自然選擇等其他因素也在影響VOC密碼子使用偏好性。

圖4 SARS-CoV-2及VOC病毒株的所有CDS的PR2分析Fig.4 PR2 analysis of all CDS of SARS-CoV-2 and VOC

2.6 系統(tǒng)發(fā)育分析和聚類分析 Mahmoud Kandeel等[13]采用鄰位歸并法(Neighbor-joining Method，NJ)方法分析出Omicron與Alpha關(guān)系密切，同為一個分支。我們利用MEGA 11.0分析軟件，基于最大似然法(Maximum Likelihood Method，MLM)分析系統(tǒng)發(fā)育，結(jié)果顯示Omimon已形成了一個新的分支群(圖5A)。而基于密碼子偏愛性的聚類分析結(jié)果顯示，Omicron的S、ORF1a蛋白(圖5B、C)與其他VOC差異最大，單獨聚為一類；E、N、ORF3a等(圖5D、E、F)大多數(shù)蛋白與其他VOC關(guān)系密切。

注：A為VOC及Wuhan-Hu-1 系統(tǒng)發(fā)育分析；B、C、D、E、F分別為基于Wuhan-Hu-1及VOC病毒株的S、ORF1a、E、N、ORF3a蛋白編碼序列密碼子偏愛性的聚類分析。圖5 SARS-CoV-2 VOC病毒株的系統(tǒng)發(fā)育分析和部分蛋白聚類分析Fig.5 Phylogenetic analysisand of partial protein's cluster analysisof SARS-CoV-2 VOC virus strain

3 討 論

2021年11月9日，南非首次報道了SARS-CoV-2 Omicron的變異株(B.1.1.529/BA.1)，11月25日感染Omicron病例數(shù)的占比就高達新增病例的75%，其傳播速度比Delta毒株更快[17]，傳染性比原始毒株高500%以上[18]。全球共享流感數(shù)據(jù)倡議組織(Global influenza data sharing initiative，GISAID)顯示Omicron 目前已在163個國家傳播，成為全球優(yōu)勢流行毒株。2022年1月8日，中國內(nèi)陸也出現(xiàn)了Omicron首例病例[18]。前人研究發(fā)現(xiàn)，SARS-CoV-2的突變主要發(fā)生在N、S以及ORF1a、ORF3a和ORF8等蛋白中[19]。Omicron 50多處突變位點同樣集中在S、ORF1ab蛋白[18]，其中S蛋白高達 32 處及以上，是 Delta 的兩倍[20]。Omicron氨基酸突變位點的增多，被認為是其傳播性增加的主要原因[21]。本研究則主要以O(shè)micron各蛋白的密碼子偏愛性改變作為研究對象，探究其與流行病學(xué)之間的聯(lián)系，為Omicron提供基因水平上的相關(guān)研究基礎(chǔ)。

研究發(fā)現(xiàn)，Omicron 關(guān)鍵蛋白如S、E、ORF1ab蛋白等密碼子使用偏愛性較其他VOC進一步減弱，這種使用偏性的減弱能使其在不競爭有限的tRNA資源的情況下保持有效復(fù)制，從而更易于在宿主細胞中表達[22]。這可能與Omicron具有更強的傳播速度快[5]有關(guān)。N、M、ORF3a等蛋白使用偏性增強則可能會降低該蛋白在宿主體內(nèi)的表達量，降低病毒對宿主的損傷，逃避宿主免疫防御，這是與宿主協(xié)同進化的結(jié)果。該結(jié)果從基因水平上解釋了Omicron出現(xiàn)傳播速率快，癥狀較輕，多數(shù)為輕型患者或者無癥狀感染者[5]這一流行病學(xué)特點的原因。

此外，為適應(yīng)宿主，病毒密碼子使用偏性往往會朝宿主方向進化[23-24]。RSCU結(jié)果顯示，Omicron毒株S、ORF1ab等關(guān)鍵蛋白的密碼子使用模式較早期的Wuhan-Hu-1更靠近人類基因密碼子，兩者使用模式更趨于一致。這意味著較Wuhan-Hu-1相比，Omicron在宿主體內(nèi)基因利用和翻譯效率增加，導(dǎo)致S、ORF1ab編碼表達的蛋白量的提高。鑒于S蛋白主要負責(zé)介導(dǎo)病毒進入，ORF1ab編碼產(chǎn)生的一系列非結(jié)構(gòu)蛋白主要參與病毒的復(fù)制和釋放[24]，因此，這種與人類基因組密碼子使用模式一致性的提升，可能是Omicron感染性增強的另外一個原因。

我們進一步研究了Omicron 密碼子偏愛性的影響因素。其中，中性分析結(jié)果提示自然選擇是導(dǎo)致SARS-CoV-2密碼子使用偏性改變的主要動力，這可能與世界范圍的疫苗接種有關(guān)。而ENC-plot 分析進一步顯示Omicron毒株S、E蛋白受到了比其他VOC更大的突變壓力，推測這是Omicron為逃避免疫而在其他VOC毒株突變基礎(chǔ)上進一步突變積累的結(jié)果。此外，系統(tǒng)發(fā)育分析顯示，Omicron已形成了一個新的分支群。而基于密碼子偏愛性的聚類分析發(fā)現(xiàn)，Omicron毒株N、E、ORF3a蛋白分別與其他VOC聚為一類，而S、ORF1a蛋白的密碼子偏愛性已經(jīng)發(fā)生較大改變并單獨聚為一類。

綜上可見，關(guān)鍵蛋白密碼子偏愛性的改變可能是SARS-CoV-2突變株，尤其是Omicron感染性增強的重要原因。這提示我們在進行SARS-CoV-2流行監(jiān)測及研究中，除了關(guān)注氨基酸的突變外，也應(yīng)進一步加強對其密碼子偏愛性變化的關(guān)注。

利益沖突：無

引用本文格式：李亞飛，羅春雨，石哲芳，等. 新型冠狀病毒Omicron變異株密碼子偏愛性及其進化分析[J].中國人獸共患病學(xué)報，2022，38(7)：559-565，576. DOI:10.3969/j.issn.1002-2694.2022.00.085