肉苁蓉属植物叶绿体基因组密码子偏好性分析

冯 展，江 媛，郑 燕，缪雨静，黄林芳*，罗光明*

1.江西中医药大学药学院，江西 南昌 330004

2.中国医学科学院北京协和医学院 药用植物研究所，北京 100193

3.大理大学药学院，云南 大理 671000

肉苁蓉属CistancheHoffmanns.&Link 是列当科（Orobanchaceae）的多年生寄生草本植物，主要分布于欧、亚洲，中国主要有4 种肉苁蓉属植物，分布在内蒙古、宁夏、甘肃、青海以及新疆等地[1]。本属植物在中国药用植物中占有及其重要的地位，其中肉苁蓉Cistanche deserticolaMa 的药用价值最高，因其主要生长在沙漠地区，素有“沙漠人参”的美誉。然而，肉苁蓉属药用植物面临着植物分类混乱[2]、市售品种混用[3]的困境。

叶绿体是大多数绿色植物的光合作用场所，参与发育和次级代谢活动[4]，并协调细胞器和核基因组之间的基因表达[5]。叶绿体具有自主遗传的基因组，被广泛用于植物系统发育分析，物种鉴定和遗传多样性表达等研究。近年来，随着叶绿体基因组高通量测序技术逐渐成熟，肉苁蓉属植物肉苁蓉、盐生肉苁蓉C.salsa(C.A.Mey.)G.Beck、沙苁蓉C.sinensisG.Beck、管花肉苁蓉C.tubulosaWight等多种植物均已开展叶绿体测序研究，并对其系统发育和遗传多样性进行了分析[2,6-7]，但尚未发现有关肉苁蓉属植物绿体基因组密码子偏好性方面的研究报道。

密码子（codon）又称遗传密码，是核酸与蛋白质联系的桥梁，是识别和传递生物遗传信息的重要载体，在生物遗传和变异中至关重要的一部分[8]。由于不同的物种其蛋白质翻译过程存在差异，在翻译过程中倾向于使用一种或几种特定的同义密码子，这种现象称为密码子偏好性（codon usage bias，CUB）[9]，密码子偏好性在mRNA 翻译、DNA 转录、蛋白的结构、表达、功能和共翻译折叠等细胞代谢过程中发挥着重要作用[10]。石岩硕等[11]通过分析4 种人参属PanaxLinn.植物叶绿体基因组密码子偏好性，推断出同属植物间亲缘关系更近；宋芸等[12]证明可基于密码子偏好性对ICE1 基因进行优化，使得其在低温胁迫下表达；李显煌等[13]发现密码子偏好性可反映出物种间的进化关系；张俊焱等[14]证明突变和自然选择共同影响獐牙菜

Swertia bimaculate(Sieb.et Zucc.)Hook.f.et Thoms.ex C.B.Clark 叶绿体基因组密码子偏好性，为獐牙菜引种驯化提供依据。因此，研究植物叶绿体基因组的密码子使用模式，可为提高基因表达载体构建效率、探讨物种进化关系、理解生物与环境适应的分子机制、改良植物品种等方面提供数据支持[15]。

本研究基于高通量测序对肉苁蓉、盐生肉苁蓉、沙苁蓉、管花肉苁蓉四种肉苁蓉属植物的叶绿体基因组密码子的使用模式进行了分析，确定4 种肉苁蓉属的叶绿体基因组密码子偏好性及其影响因素，为后续肉苁蓉属叶绿体基因组的基因表达检测、适应性及进化关系研究提供参考依据。

1 材料

肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉植株由不同地区采集（表1），经北京协和医学院药用植物研究所黄林芳研究员鉴定为肉苁蓉属植物肉苁蓉C.deserticolaMa、盐生肉苁蓉C.salsa(C.A.Mey.)G.Beck、沙苁蓉C.sinensisG.Beck 和管花肉苁蓉C.mongolicaBeck 的全株植物，储存于北京协和医学院药用植物研究所植物标本室，凭证标本为CMPB16201、CMPB16202、CMPB16203和CMPB16204。叶绿体全基因组序列已上传至GenBank，收录号分别为MN614127（肉苁蓉）、MN614128（盐生肉苁蓉）、MN614129（沙苁蓉）和MN614130（管花肉苁蓉）。肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉的叶绿体基因组大小分别为109 495、111 710、111 500、75 375 bp，分别包括60、61、60、53 个蛋白质编码基因序列（coding DNA sequence，CDS）。由于短序列无法准确地计算有效密码子数，剔除了长度小于 300 bp 的CDS[16]，又因为终止密码子UAA、UGA 和UAG不编码任何氨基酸，并且UGG 和AUG 分别是色氨酸和蛋氨酸的唯一密码子，这些密码子不存在偏好性[17]，故选择的CDS 以ATG 为起始[18]，以TAA、TAG、TGA 为结尾，最后为了降低结果的误差，剔除了重复CDS。处理后每个物种留存18 条CDS，并用于后续分析。

表1 肉苁蓉样品信息Table 1 Information of samples

2 方法

2.1 密码子偏好性分析

将挑选出的72 条CDS 整合到一个fasta 文件中，利用线上软件CUSP（https://www.bioinformatics.nl/ cgi-bin/emboss/cusp）计算每个基因的编码区总GC 含量（GCall）以及密码子第1 位（GC1）、第2位（GC2）、第3 位（GC3）的GC 含量。利用CodonW1.4.2 软件[19]计算有效密码子数（effective number of codon，ENC）和同义密码子的相对使用频率（relative synonymous codon usage，RSCU），并进行对应性分析（correspondence analysis，COA）。然后通过Microsoft Excel 和IBM SPSS19.0 软件对分析数据进行整理和绘图。

2.2 中性绘图分析

中性绘图分析是判断密码子偏好性影响因素的方法之一[20]。同义密码子的突变位点通常为第3 碱基，而非同义密码子的突变位点通常为第1 或第2碱基。以每个基因的GC1 和GC2 的平均值（GC12）为纵坐标、以GC3 为横坐标在Microsoft Excel 软件中绘制散点图，剖析GC3 与GC12 的相关性，其中每个散点则代表1 个基因。如果所有的点都沿对角线分布，说明密码子3 个位置的碱基没有明显差异，只受突变压力的影响；如果GC12 和GC3 变异的相关性很低，则表明自然选择是密码子使用模式的主要影响因素[21]。

2.3 ENC-plot 分析

ENC 代表了密码子偏离随机选择的程度[22]，以ENC 为纵坐标，GC3 为横坐标在Microsoft Excel软件中进行二维散点图的绘制，以公式ENC=2＋GC3＋29/[GC32＋(1-GC3)2]绘制标准曲线[23]，每个散点代表1 个基因。散点图可反映出密码子使用偏好性和基因碱基组成之间的关系，如果密码子使用偏好性主要受到突变压力的影响，散点将位于标准曲线上或略低于标准曲线，相反，如果基因低于标准曲线，则主要受到选择和其他因素的影响[24]。

2.4 PR2-plot 分析

PR2 偏倚分析（PR2-bias plot analysis）对密码子第3 位碱基上的A、T、C、G 含量进行分析，避免了密码子第3 位碱基A/T 和C/G 之间的突变不平衡。分别以G3/(G3＋C3)为横坐标、以A3/(A3＋T3)为纵坐标在Microsoft Excel 软件中制作散点图。通过判断散点与中心点（A＝T，C＝G）的矢量距离，来确定碱基偏移的程度和方向[25]。

2.5 COA

COA 从多维空间中提取最具影响的轴和方向，被广泛用于研究基因间密码子使用变化的主要趋势[26]，通过分析COA 轴分离出的基因能识别影响密码子使用偏性的主要原因[27]。利用CodonW1.4.2 软件进行COA，根据肉苁蓉属四种植物的RSCU 值，将所有基因密码子绘制在59 维的向量空间中，以反映密码子使用变化的趋势。轴1（Axis 1）的占比代表对密码子使用频率变化的最大影响，其余58 个轴代表了逐渐减小的因素。根据基因在多维空间中的同义密码子使用情况，可以分析其主要变异来源[28]。

2.6 最优密码子分析

将4 种植物密码子的ENC 值作为参考指标，取ENC 值最高（10%）和最低（10%）的基因序列分别建库，作为高表达组和低表达组。利用CodonW 1.4.2 软件计算2 个表达库的RSCU 值和ΔRSCU 值（ΔRSCU＝RSCU 高表达－RSCU 低表达）。其中，ΔRSCU≥0.08 的密码子作为高表达密码子，RSCU＞1 的为高频密码子，同时满足2 个条件的密码子为最优密码子[29]。

3 结果与分析

3.1 叶绿体基因的密码子偏好性分析

肉苁蓉属4 个物种叶绿体基因组的GC 含量见图1。肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉的18 个叶绿体基因整个编码区的GC 含量（GCall）范围在36.60%～37.22%，叶绿体基因组第1、第2和第3 位密码子的GC 含量分别为 44.03%～44.75%，38.97%～39.82%和25.80%～27.08%。结果显示，4 个物种的平均GC 含量均小于50%，这说明4 个物种的叶绿体基因组倾向于使用A/T 碱基和A/T 结尾的密码子。ENC 值从20（极偏好，每个氨基酸只使用1 个密码子）到61（无偏好，每个氨基酸均匀使用所有同义密码子），当ENC 值≤35时，认为密码子使用有非常显著的偏好[22,30]，4 个物种的叶绿体基因组编码区的ENC 值均大于35，说明肉苁蓉属物种绿体基因密码子偏好性较弱。

图1 肉苁蓉属4 个物种的叶绿体基因组密码子不同位置的GC 含量Fig.1 GC content of different positions of codon in chloroplast genome of Cistanche

4 个肉苁蓉属物种各叶绿体基因的GC1、GC2、GC3、GCall、ENC 和密码子数（codon counts，CC）相关性分析见图2。4个物种的GCall均与GC1和GC2呈极显著相关（P＜0.01），其中盐生肉苁蓉与管花肉苁蓉的GC1 和GC2 呈显著相关（P＜0.05），表明这2 个物种密码子第1 位碱基和第2 位碱基组成相似，但与第3 位碱基组成存在显著差别；肉苁蓉的ENC值与GC3 呈显著相关，盐生肉苁蓉的ENC 值与GC1和GC2 呈显著相关，说明在肉苁蓉中，密码子第3位碱基的组成对密码子使用偏性有重要影响，在盐生肉苁蓉中密码子第1、3 位碱基的组成对密码子使用偏性有重要影响，而在沙苁蓉和管花肉苁蓉中，碱基的组成对密码子使用偏性并无显著影响；4 个物种的基因密码子数CC与ENC值之间的相关系数介于0.18（盐生肉苁蓉）～0.27（沙苁蓉），相关性均未达到显著水平，这说明CC 对ENC 的影响较弱，叶绿体基因的长度对密码子使用偏性无影响。

图2 肉苁蓉属叶绿体基因组相关参数的相关性分析Fig.2 Correlation analysis of chloroplast genome-related parameters of Cistanche

为了进一步探究密码子使用模式，分别计算四个物种的RSCU 值（表2），结果显示4 个物种的RSCU 值相差无异。RSCU＞1 表明该密码子使用偏性较高，反之，则为非偏好密码子。肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉中，RSCU＞1 的氨基酸占比分别为48.4%、50%、51.6%和46.9%，其中以A/U 结尾的氨基酸的个数分别为27，27，28 和26，这说明肉苁蓉属叶绿体基因组中A 或U 出现频率较高，是该属叶绿基因组偏好密码子。RSCU＜1的氨基酸中，多数以G/C 结尾，表明这些是肉苁蓉属叶绿体基因的非偏好密码子。

表2 肉苁蓉属叶绿体同义密码子的使用频率Table 2 RSCU analysis of protein coding region in chloroplast of Cistanche

3.2 中性绘图分析

肉苁蓉属叶绿体基因中性绘图见图3，肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉的GC12 与GC3的相关系数分别0.083、0.084、0.033 和0.245，且相关性均未达到显著水平（P＝0.05），说明突变压力不会对对密码子使用偏性产生明显的影响；中性绘图斜率显示，4 个叶绿体基因组密码子使用模式的突变压力占4.44%～32%，这说明4 个叶绿体基因组密码子使用模式的突变压力占比小，自然选择是影响密码子使用偏性的主要因素。

图3 肉苁蓉属中性绘图分析Fig.3 Analysis of neutrality plot

3.3 ENC-plot 绘图分析

肉苁蓉属4 个物种的ENC-plot 绘图见图4，大部分基因与标准曲线有一定的距离，只有小部分基因分布在标准曲线上，这说明ENC 的实际值与理论值存在一定的偏差且偏性较弱。大部分基因分布在标准曲线下方，说明基因偏性主要受到自然选择的影响。

图4 ENC-plot 分析Fig.4 ENC-plot analysis

3.4 PR2-plot 绘图分析

肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉中各基因不均匀分布在PR2-plot 绘图的4 个平面区域中（图5），结果显示，大多数基因远离中心，分布在右下角，值得注意的是几乎没有基因分布在左下角，这表明GC3 的使用主要受到自然选择影响，且4 种植物的G3 存在T/G使用偏性。

图5 PR2-plot 分析Fig.5 PR2-plot analysis

3.5 对应性分析

为了探究肉苁蓉属植物叶绿体密码子变异源头，对肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉叶绿体基因组密码子进行了基于RSCU 值的对应性分析（图6），前4 轴的累计可解释变异分别占总变异的52.97%、53.96%、53.29%和52.18%，第1 轴的可解释变异分别占总变异的 17.89%、18.70%、17.39%和17.57%，第2～4 轴的可解释变异分别为肉苁蓉（12.87%、11.35%和10.87%）、盐生肉苁蓉（12.91%、12.00%和10.34%）、沙苁蓉（13.10%、12.07%和 10.73%）和管花肉苁蓉（12.12%、11.65%和10.83%）。4 个物种的第2、3、4 轴的可解释变异均小于第1 轴，因此第1 轴对密码子偏好性的贡献率最大，但第1 轴仅代表肉苁蓉属叶绿体基因组中密码子使用偏性的部分变异，这说明肉苁蓉属叶绿体的密码子使用偏性的影响因素并非一个，还可能与突变、自然选择、基因长度和功能等因素有关。

图6 对应性分析Fig.6 Corresponding analysis

3.6 最优密码子分析

肉苁蓉属4 个物种的最优密码子分析见表3。结果显示，肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉各有10、11、13 和12 个。其中肉苁蓉与盐生肉苁蓉最为相似，共有9 个同样的密码子（图7），4 个物种共有的密码子有4 个（GUA、UCA、CCU和ACA）。在共同的最优密码子中，ΔRSCU 均大于0.5 的仅有1 个，为GUA。

表3 最优密码子分析Table 3 Optimal codons in chloroplast genome of Cistanche

续表3

图7 最优密码子个数韦恩图Fig.7 Venn diagram of optimal codons

4 讨论

密码子偏好性指的是在蛋白质的翻译过程中不同同义密码子使用频率的高低现象，是生物基因组进化的一个重要特征，密码子偏好性与GC 含量、tRNA 丰度、基因表达水平和蛋白质结构等因素相关[15]，研究植物密码子的偏性对深入探究分子进化和外源蛋白表达有着及其重要的意义[31-32]。在影响密码子偏好性的诸多因素中，突变和自然选择为主导因素，自然选择使得植物在翻译中优先编码最优密码子，突变会导致植物中部分非偏好性密码子的存在[33-34]，所以在此次研究中，密切关注了基因突变压力和自然选择压力。

本研究中肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉叶绿体基因组中GCall 为36.60%～37.22%，GC1、GC2 和GC3 的含量分别为44.03%～44.75%，38.97%～39.82%和25.80%～27.08%，含量均在50%以下，表明4 种肉苁蓉属植物的叶绿体基因第3 位碱基多以A/T 结尾或更偏向于A/T 碱基，这与前人报道的大戟科[32]、茄科[35]、人参属[11]、豆蔻属[36]等植物叶绿体基因组研究一致，说明不同物种间叶绿体基因密码子使用偏好性存在着一定的相似性。

密码子第3 位碱基的变化通常不会导致编码氨基酸的改变，因此研究密码子第3 位碱基组成对研究密码子偏好性具有重要意义[37]。对肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉叶绿体基因组的第3位碱基上的A、T、C、G 含量进行PR2 分析,结果显示：T＞A、G＞C，这与小麦Triticum aestvumLinn.[38]、剑麻Agave sisalanaPerr.ex Engelm.[24]、黄芩Scutellaria baicalensisGeorgi[39]、陆地棉Gossypium hirsutumLinn.[40]等植物表现出的T＞A、C＞G 和胡杨Populus euphraticaOliv.[41]、禾本科（Poaceae Barnhart）植物[42]表现出的A＝T、G＝C的结果不一致。由此可见，不同植物的叶绿体基因组中碱基组成情况存在明显差异，这导致了密码子偏好性有所不同。ENC 值是判断判断密码子偏好性强弱的标准，当ENC 低于35 时，密码子使用偏性较大，当ENC 高于35 时，偏性较小。在肉苁蓉属的四种植物种，叶绿体基因组编码区的ENC 值均大于35，说明肉苁蓉属物种绿体基因密码子偏好性较弱。在相关性分析中，4 个物种的GCall 与GC1 和GC2 均呈极显著相关，除此之外，其余数值的相关性各有不同，由于这4 种肉苁蓉属植物的生境不尽相同，可以推测生长环境的不同可能会导致不同物种间的GC 差异。

中性绘图分析、ENC-plot 分析、PR2-plot 的结果均说明肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉叶绿体基因组密码子更偏向于受自然选择的影响，这与之前报道的豆蔻属AmomumRoxb.[36]、含笑属MicheliaLinn.[43]和人参属[11]植物结果一致，但蒺藜苜蓿Medicago truncatulaGaertn.[17]、拟南芥Arabidopsis thaliana(L.)Heynh.和杨树Populus przewalskiiMaxim.[41]等植物的叶绿体基因组密码子偏好性则主要受突变的影响，这表明不同植物的密码子偏好性是多因素综合的结果。

最优密码子分析表明，肉苁蓉、盐生肉苁蓉、沙苁蓉和管花肉苁蓉的最优密码子数为10、11、13和12 个，共有的最优密码子数个数为4。值得注意的是，肉苁蓉与盐生肉苁蓉共享9 个相同的最优密码子，是肉苁蓉属4 个物种中最为相近的2 个物种，可以推断出在肉苁蓉属内，肉苁蓉与盐生肉苁蓉具有相近的密码子使用偏好性，亲缘关系更为接近，这与缪雨静等[44]和杨俏俏等[2]通过肉苁蓉属的系统发育分析所得结论一致。

本研究首次基于叶绿体基因组密码子偏好性对国产肉苁蓉属进行了分析，系统性的阐明了影响密码子偏好性的因素，揭示了影响其密码子偏好性的因素与自然选择压力相关，筛选出了最优密码子，并且发现了属内肉苁蓉与盐生肉苁蓉亲缘关系更为接近。为后续的肉苁蓉属物种基因工程研究、资源鉴定、外源表达基因构建等提供了科学依据和理论支撑。

利益冲突所有作者均声明不存在利益冲突