基于轉(zhuǎn)錄組數(shù)據(jù)分析藥用大黃的密碼子使用偏好性

高夢(mèng)琦，鄒建珍，霍小位，李依民，張明英，高 靜，張 崗*，劉蒙蒙*

1.河北大學(xué)中醫(yī)學(xué)院，河北 保定 071002

2.河北大學(xué)藥學(xué)院，河北 保定 071002

3.陜西中醫(yī)藥大學(xué)藥學(xué)院 陜西省秦嶺中草藥應(yīng)用開發(fā)工程技術(shù)研究中心，陜西 西安 712046

藥用大黃Rheum officinaleBaill.屬于蓼科大黃屬多年生草本植物，為《中國藥典》2020年版收錄的大黃3種基原植物之一，具有瀉下攻積、清熱瀉火、涼血解毒、逐瘀通經(jīng)、利濕退黃等功效[1]。現(xiàn)代中藥化學(xué)研究表明大黃中主要含有蒽醌類、酚類等物質(zhì)，其中蒽醌類化合物是主要的活性物質(zhì)，具有抗炎、抗腫瘤作用[2-3]。目前藥用大黃的資源主要來源于人工種植，隨著過度采挖，野生資源已經(jīng)瀕臨滅絕。因此，開展蒽醌類化合物生物合成機(jī)制研究，探索蒽醌類化合物生物合成路徑，是解決藥用大黃野生資源瀕危的思路之一。

遺傳密碼子是生物體DNA與蛋白之間信息傳遞的信息載體。大多數(shù)氨基酸不僅由一種密碼子編碼，而是可以同時(shí)由2～3種密碼子編碼，生物對(duì)于不同的密碼子的使用具有一定偏好性[4]。研究發(fā)現(xiàn)，密碼子使用偏好性受到GC含量（鳥嘌呤和胞嘧啶含量）[5]、基因表達(dá)水平[6]、環(huán)境壓力[7]、突變壓力[8]、基因長度[9]、蛋白結(jié)構(gòu)[10]及tRNA豐度[11]的影響。隨著高通量測(cè)序技術(shù)的發(fā)展，藥用植物的密碼子使用偏好性研究有了顯著的成果[12-13]，但是目前還沒有關(guān)于藥用大黃密碼子偏好性的研究。生物體選擇最優(yōu)密碼子可以體高基因的翻譯效率及蛋白的表達(dá)量。因此，研究密碼子的使用偏好性，有利于基因工程研究中密碼子優(yōu)化、最適外源宿主選擇，并對(duì)生物的系統(tǒng)進(jìn)化研究具有重要意義。

本研究利用藥用大黃的轉(zhuǎn)錄組數(shù)據(jù)對(duì)其編碼序列進(jìn)行密碼子偏好性分析，揭示藥用大黃基因密碼子使用偏好性的特征，為藥用大黃基因基因工程研究中進(jìn)行密碼子優(yōu)化及選擇最佳的外源載體提供理論基礎(chǔ)。

1 數(shù)據(jù)

藥用大黃轉(zhuǎn)錄組（SRR10838870）來源于本課題組前期選取藥用大黃種苗根、根莖、樣品，利用Illumina HiSeq 2500 測(cè)序平臺(tái)進(jìn)行轉(zhuǎn)錄組測(cè)序并拼接獲得Uingenes序列，利用BLASTx及ESTScan軟件對(duì)上述Uingenes序列進(jìn)行編碼序列（coding sequence，CDS）分析。使用perl程序篩選編碼序列中長度大于300 bp并且同時(shí)具有起始密碼子及終止密碼的序列，用于藥用大黃密碼子使用偏好性分析，最后共得到4 733條CDS進(jìn)行分析。

2 方法

2.1 編碼序列中堿基含量及中性繪圖分析

使用Codon W軟件分析藥用大黃編碼序列中的有效密碼子數(shù)（effective number of codons，ENc）、密碼子適應(yīng)指數(shù)（codon adaptation index，CAI）及同義密碼子相對(duì)使用度（relative synonymous codon usage，RSCU），計(jì)算序列中的GC含量，利用編寫的perl程序計(jì)算GC1、GC2及GC3的含量，即密碼子第1、2、3位的GC含量。GC12為GC1、GC2的平均值。

中性繪圖分析主要用于分析GC12與GC3之間的相關(guān)性，以GC12為縱坐標(biāo)，以GC3為橫坐標(biāo)制作散點(diǎn)圖，并進(jìn)行線性擬合回歸分析?？捎糜谠u(píng)價(jià)自然選擇壓力及突變對(duì)密碼子使用偏好性造成的影響。

2.2 GC3關(guān)聯(lián)及奇偶偏好性分析

ENc-GC3關(guān)聯(lián)分析主要用于分析密碼子使用偏好性與堿基組成之間的關(guān)聯(lián)。以ENc為縱坐標(biāo)、GC3s為橫坐標(biāo)繪制二維散點(diǎn)圖，并進(jìn)行非線性擬合回歸分析。如果基因的點(diǎn)分布在擬合曲線附近，則表明密碼子的使用偏好性主要是受到突變作用的影響；如果基因點(diǎn)分布在擬合曲線下方，則表明密碼子的使用偏好性主要是受到選擇壓力的影響。奇偶偏好性分析以A3/(A3＋T3)為縱坐標(biāo)，以G3/(G3＋C3)為橫坐標(biāo)制作二維散點(diǎn)圖，以（0.5，0.5）坐標(biāo)為中心點(diǎn)，從這個(gè)中心點(diǎn)制作的失量表示了奇偶偏倚的程度和方向[14]。

2.3 最優(yōu)密碼子分析

對(duì)于最優(yōu)密碼子的篩選，本研究利用藥用大黃編碼書寫的CAI值進(jìn)行由高到低排列，然后取前后各5%的編碼序列，分別代表高表達(dá)與低表達(dá)的基因序列。之后分析各個(gè)密碼子的平均RSCU值。當(dāng)密碼子在高低表達(dá)基因中RSCU值的差異＞0.3，且在高表達(dá)基因中RSCU≥1、在低表達(dá)基因中的RSCU≤1，即可被確定為最優(yōu)密碼子[15-16]。

2.4 對(duì)應(yīng)分析

對(duì)應(yīng)分析是指利用多元統(tǒng)計(jì)分析的方法研究不同編碼序列密碼子的變異趨勢(shì)。利用Codon W軟件中對(duì)應(yīng)分析功能，將所有編碼序列分布到59維向量空間，進(jìn)而獲得密碼子使用偏好性的特點(diǎn)。通過分析基因位置研究基因向量及基因間值的變異程度，預(yù)測(cè)使密碼子使用發(fā)生偏好性的原因[17]。

3 結(jié)果

3.1 密碼子GC含量組成分析

利用Codon W軟件計(jì)算藥用大黃編碼序列中GC含量，將藥用大黃編碼序列作為一個(gè)整體進(jìn)行計(jì)算。結(jié)果如圖1所示，GC含量范圍分布在28.6%～71.9%，主要分布在40%～50%，GC的平均含量為45.6%。GC1、GC2、GC3的平均含量分別為48.63%、40.89%、44.37%。對(duì)比發(fā)現(xiàn)，GC2的含量最小，GC1與GC3的值相對(duì)比較接近，說明藥用大黃編碼序列更傾向于使用鳥嘌呤（G）和胞嘧啶（C）這2種堿基。

圖1 藥用大黃編碼序列中GC含量組成Fig.1 Distribution of GC content ofcoding sequences(CDS) in R.officinale

3.2 中性繪圖分析

藥用大黃編碼序列的中性圖分析如圖2所示，其中GC12的數(shù)值范圍為0.15～0.87，GC3的數(shù)值范圍為0.09～0.87。GC12與GC3的相關(guān)性系數(shù)為0.155，表明這兩者之間的相關(guān)性比較弱，說明選擇壓力在藥用大黃密碼子使用偏好性的選擇具有一定影響。

圖2 中性圖分析Fig.2 Neutrality plot

3.3 ENc-plot繪圖分析

ENc的取值范圍為20～61，平均值為47.58，其中ENc≤35的基因有386個(gè)（8.15%），ENC＝61的基因有251個(gè)（5.30%），說明藥用大黃的編碼序列除了少數(shù)密碼子具有較強(qiáng)的偏好性，其整體的密碼子使用偏好性較弱。以ENc為縱坐標(biāo)、GC3s為橫坐標(biāo)繪制散點(diǎn)圖（圖3），圖中各點(diǎn)代表基因分布情況。圖3顯示大部分基因遠(yuǎn)離期望曲線，說明密碼子的偏好性主要來源于選擇壓力的影響，少部分基因均勻的分布在期望曲線兩側(cè)，說明部分基因收到突變壓力的影響。

圖3 藥用大黃轉(zhuǎn)錄組編碼序列ENC-GC3s相關(guān)性分析Fig.3 Neutrality plot (ENC vs GC3)

利用(ENCexp－ENCobs)/ENCexp統(tǒng)計(jì)分析期望ENc與觀測(cè)ENc之間的區(qū)別[18]。如圖4所示，(ENCexp－ENCobs)/ENCexp的取值主要分布在0～0.2，說明期望ENc與觀測(cè)ENc兩者之間的存在很小差異。進(jìn)而說明在藥用大黃編碼序列密碼子使用偏好性形成過程中，突變壓力不是唯一的影響因素，但卻是一個(gè)重要的影響因素。

圖4 ENc比值頻率分布Fig.4 Frequency distribution of ENc ratio

3.4 對(duì)應(yīng)性分析

對(duì)應(yīng)性分析結(jié)果見圖5，GC含量不同的基因被標(biāo)記上不同的顏色，編碼序列中GC含量大于60%的被標(biāo)記為紅色，GC含量在45%～60%的被標(biāo)記為藍(lán)色，GC含量小于45%的被標(biāo)記為綠色。其中，GC含量大于60%的編碼序列分散的分布在坐標(biāo)系中，而GC含量小于60%的編碼序列則在坐標(biāo)系中分布比較集中。

圖5 對(duì)應(yīng)性分析Fig.5 Correspondence analysis

此外，選取6個(gè)編碼序列相關(guān)的參數(shù)（GC12、GC3、GCall、ENC、CAI以及Axis 1）進(jìn)行相關(guān)性分析，結(jié)果如表1所示，Axis 1與GC12、GC3這2個(gè)重要參數(shù)存在顯著的負(fù)相關(guān)性（r＝-0.418，r＝-0.728，P＜0.01）。結(jié)合圖5與表1，可以說明GC含量與GC3含量高的編碼序列取向于分布在Axis 1左側(cè)。因此，推測(cè)以G/C結(jié)尾的密碼子應(yīng)集中分布在Axis 1負(fù)軸而以A/U結(jié)尾的密碼子集中分布在Axis 1正軸。

表1 藥用大黃編碼序列對(duì)應(yīng)參數(shù)的相關(guān)性分析Table 1 Correlation analysis of R.officinale gene-related parameters

3.5 偏倚性分析

偏倚性分析結(jié)果如圖6所示，大部分藥用大黃轉(zhuǎn)錄組編碼序列樣本點(diǎn)分布于中心點(diǎn)附近，說明其堿基中嘌呤與嘧啶出現(xiàn)的頻率相對(duì)一致，少部分編碼序列偏離中心點(diǎn)較遠(yuǎn)，進(jìn)一步說明藥用大黃密碼子使用偏好性形成過程中主要受到突變壓力的影響。

圖6 編碼序列偏倚性分析Fig.6 Parity rule 2 (PR2 ) -bias plot

3.6 基因表達(dá)水平對(duì)密碼子使用偏好性的影響

藥用大黃編碼序列的CAI的范圍為0.261～0.835，表明藥用大黃的基因表達(dá)水平比較高。同時(shí)，對(duì)CAI與其他幾個(gè)重要的參數(shù)進(jìn)行了相關(guān)性分析，包括ENc、GC3、GC。結(jié)果如表1及圖7所示，CAI與ENc及GC之間都存在顯著的正相關(guān)性（r＝0.07，P＜0.01；r＝0.695，P＜0.01）。因此，基因表達(dá)水平在藥用大黃密碼子使用偏好性形成過程中起到關(guān)鍵作用，基因表達(dá)水平高并且基因中GC含量高，密碼子使用偏好性程度越高。

圖7 中性圖分析Fig.7 Neutrality plot

3.7 最優(yōu)密碼子分析

藥用大黃編碼序列高、低表達(dá)序列的RSCU值如表2所示。從中篩選到29個(gè)密碼子作為藥用大黃編碼序列的最優(yōu)密碼子，在這些最優(yōu)密碼子中，以U（15個(gè)）及A（13個(gè)）結(jié)尾的較多，表明藥用大黃轉(zhuǎn)錄組編碼序列的最優(yōu)密碼子偏好以U或A結(jié)尾。

表2 藥用大黃轉(zhuǎn)錄組高、低表達(dá)基因的RSCU比較分析Table 2 Optimal codons of R.officinale genes based on the RSCU values

4 討論

在生物長期的進(jìn)化過程中，形成了對(duì)密碼子的使用的固定模式。今年來，隨著二代測(cè)序技術(shù)的發(fā)展，探索天然產(chǎn)物合成路徑并進(jìn)行異源合成，成為了研究熱點(diǎn)。因此，探索藥用大黃基因密碼子使用偏好性，可為研究其密碼子使用特點(diǎn)，篩選最佳蛋白異源表達(dá)載體及密碼子優(yōu)化奠定理論基礎(chǔ)。

本研究對(duì)藥用大黃轉(zhuǎn)錄組數(shù)據(jù)中4733條完整編碼序列的密碼子使用偏好性進(jìn)行分析，結(jié)果顯示，藥用大黃轉(zhuǎn)錄組編碼序列的GC、GC3平均含量分別為45.6%、44.73%，表明藥用大黃轉(zhuǎn)錄組中的編碼序列偏好使用富含A、U且以A/U結(jié)尾的密碼子，這與小麥Triticum aestivumL.[19]、水稻Oryza sativaL.[20]、玉米Zea maysL.[21]密碼子使用模式不一致，這些植物的密碼子富含GC堿基，并且大多數(shù)密碼子以G/C結(jié)尾，藥用植物川貝母轉(zhuǎn)錄組中的編碼序列的GC和GC3平均含量都接近于50%，并沒有表現(xiàn)出密碼子使用偏好性[15]，表明堿基組成在藥用大黃響密碼子偏好性形成過程中起到重要作用。

此外，堿基突變和遺傳選擇也影響著密碼子使用偏好性形成的重要因素。本研究中，藥用大黃轉(zhuǎn)錄組GC12與GC3含量具有顯著相關(guān)性，ENc-plot繪圖及偏倚性分析分析結(jié)果表明，大部分基因遠(yuǎn)離于標(biāo)準(zhǔn)曲線及坐標(biāo)軸中心點(diǎn)，基因的期望ENc值與觀測(cè)ENc值具有較大差異，說明這部分基因密碼子偏好性的形成是以遺傳選擇為主導(dǎo)因素；少數(shù)的基因分布于標(biāo)準(zhǔn)曲線附近及坐標(biāo)軸中心點(diǎn)，基因的期望ENc值與觀測(cè)ENc值基本保持一致，說明這部分基因密碼子形成偏好性的形成是以突變?yōu)橹鲗?dǎo)因素。

在構(gòu)建蛋白表達(dá)載體時(shí)，不同的基因密碼子使用偏好性差異很大，因此分析藥用大黃的密碼子使用特點(diǎn)對(duì)于載體的密碼子優(yōu)化具有重要意義。本研究篩選得到藥用大黃29個(gè)最優(yōu)密碼子，大多數(shù)密碼子以A/U結(jié)尾，從藥用大黃編碼序列中的GC含量可以看出, 其偏好使用A或U堿基結(jié)尾的密碼子。由此推斷藥用大黃基因的密碼子偏好性的產(chǎn)生與其基因進(jìn)化過程中AU到GC突變壓力要比GC到AU突變壓力高有關(guān)。

本研究對(duì)藥用大黃轉(zhuǎn)錄組數(shù)據(jù)中編碼序列的密碼子使用特點(diǎn)進(jìn)行了研究，探討了影響其密碼子使用偏好性形成的影響因素，并篩選出藥用大黃基因的29個(gè)最優(yōu)密碼子，可為蒽醌類化合物的異源合成及其合成路徑研究提供理論依據(jù)。

利益沖突所有作者均聲明不存在利益沖突