基于RNA-seq技术的茶树CsbHLH62生物信息学分析

2019-03-01金孝芳马林龙刘艳丽郭桂义

西南农业学报 2019年1期

曹丹，金孝芳*，马林龙，刘艳丽，郭桂义

(1.湖北省农业科学院果树茶叶研究所，湖北武汉 430064；2.信阳农林学院/河南省豫南茶树资源综合开发重点实验室，河南信阳 464000)

【研究意义】硒是人和动物必需的微量元素之一，过量或不足对机体都有不可忽视的影响[1]。通过调查发现我国有近66.7 %的人口不同程度缺硒[2]，因其有益与有害的阈值较低，故开发安全有效的硒资源非常有必要。茶是一种非常流行的饮品，茶树的富硒能力较强，但是其硒含量因品种、部位和生长环境的不同而表现出较大的差异[3]，茶树硒累积的分子机理至今尚并不明确。bHLH转录因子参与植物多种生物学过程，对茶树bHLH基因进行生物信息学分析，有助于理解其在茶树硒吸收转运过程中的生物学功能，进一步对利用生物技术手段调控茶叶硒含量以及茶树新品种的培育都具有非常重要的作用。【前人研究进展】转录因子bHLH(basic helic-loop-helic，bHLH)广泛存在于真核生物中[4]，是植物体内第二大家族转录因子。bHLH结构域约含60个氨基酸，由2个功能不同的区域组成，一个是能与DNA结合的碱性区域[5]，另一个是α螺旋-环-α螺旋结构，其中的α螺旋亲水又亲脂，被不同长度的连接区分开[6]。bHLH转录因子内部或者转录因子之间的2个α螺旋可以相互作用，形成同源或者异源二聚体，并与靶基因启动子区域的顺式作用元件结合，对下游基因发挥调控作用[7]。目前，通过对一些重要物种的研究发现，拟南芥(ArabidopsisthalianaL.)、烟草(NicotianatabacumL.)、水稻(OryzasativaL.)、葡萄(VitisvinferaL.)、玉米(Zeamays)、杨树(Populustrichocarpa)等中分别含有164、190、180、191、276、183个bHLH转录因子[8-10]。在植物中该类转录因子主要参与了与种子萌发[11-12]、雄性不育[13]、分枝发育[14]、离子转运[15]、激素应答[16-17]、逆境胁迫[18-19]等过程。茶树bHLH转录因子的研究也有一些进展，如CsbHLH2从抗逆性较强的品种‘陕茶一号’中被成功克隆出来[20]，赵磊等[21]发现CsbHLH2和CsbHLH24与茶树类黄酮合成相关，然而通过对bHLH的进化树分析发现，在功能已知的第IIIf亚组中，没有茶树基因组中的bHLH转录因子，认为更多的茶树bHLH需要进一步挖掘；孙彬妹等[22]克隆了2个与拟南芥高度同源bHLH转录因子GL3和EGL3，分析发现其参与调控了花青素的合成。迄今为止，bHLH转录因子在茶树硒富集过程中的作用少有报道。【本研究切入点】本研究团队对茶树品种进行了转录组测序，通过分析筛选出一个表达丰度较高、差异倍数较大的bHLH转录因子(Unigene0020663)，结合基因功能注释，将其命名为CsbHLH62，通过对其进行生物信息学分析，探讨其生物学功能。【拟解决的问题】旨在为进一步验证该基因功能奠定基础。

1 材料与方法

1.1 实验材料与测序分析

1年生茶苗(鄂茶1号)于2016年夏季在湖北省农业科学院果树茶叶研究所温室大棚进行水培培养(处理条件为内部资料，未发表)，分别对茶树的根部及叶片进行转录组测序，由广州基迪奥生物科技有限公司完成。

1.2 生物信息学分析方法

CsbHLH62基因的开放阅读框采用NCBI中的ORF Finder(https://www.ncbi.nlm.nih.gov/orffinder/)预测，功能结构利用ScanProsite程序(http://www.expasy.org/tools/scanprosite/)分析，编码产物理化性质采用ProtParam程序(https://expasy.org/tools/protparam.html)分析，信号肽分析、跨膜结构域与亚细胞定位分别由在线软件SignalP4.1(http://www.cbs.dtu.dk/services/SignalP/)、TMHMM(http://www.cbs.dtu.dk/services/TMHMM/)和Cell-PLoc 2.0(http://www.csbio.sjtu.edu.cn/bioinf/Cell-PLoc-2/)预测，磷酸位点预测利用在线软件NetPhos 3.1 Server(http://www.cbs.dtu.dk/services/NetPhos/)进行，二级、三级结构分别由在线程序SMOPA(http://npsa-pbil.ibcp.fr/cgi-bin/npsa_automat.pl?page=/NPSA/npsa_sopma.html)和Phyre2(http://www.sbg.bio.ic.ac.uk/phyre2/html/page.cgi?id=index)预测，采用软件DNAman和MEGA6进行多序列比对和同源性分析。

2 结果与分析

2.1 CsbHLH62基因的开放阅读框分析

该基因全长933 bp，共编码277个氨基酸，开放阅读框(ORF)长831 bp，5’非编码区(5’UTR)长102 bp。

2.2 CsbHLH62基因编码蛋白的理化性质

利用在线程序ProtParam分析蛋白质的理化性质，得到该基因编码蛋白的化学方程式为C1318H2122N390O438S9，预测分子量30728.24 kD，共包括4277个原子，理论等电点7.73，脂肪系数62.64，不稳定系数36.46(不稳定指数大于40时，表示不稳定)，表明该蛋白稳定。

该基因编码276个氨基酸，组成中最多的氨基酸是丝氨酸(Ser)，所占比例为11.6 %，其次是亮氨酸(Leu)和赖氨酸(Lys)，均占7.6 %(图1)。另外，负电荷残基总数(Asp+Glu)是36，正电荷残基总数(Arg + Lys)是37。

利用Protscale程序预测CsbHLH62编码产物的疏水性/亲水性(图2)，分析氨基酸残基的分值发现，多肽链的第191位具有最高的分值1.422，疏水性最强；第104和106位具有最低的分值-2.933，亲水性最强，总平均亲水性(Grand average of hydropathicity )为-0.851，整条多肽链表现为亲水性，故推测该蛋白是一种可溶性蛋白这有利于其发挥转录因子的功能。

2.3 CsbHLH62编码蛋白的信号肽分析、跨膜结构域与亚细胞定位

利用软件SignalP4.1进行信号肽(Signal peptide)分析，结果如图3所示。C-score是原始信号肽裂解位点的分值，数值越大表明该点出现裂解位点的概率就越大；S-score是信号肽分值，数值越大表明该氨基酸位于信号肽区域的可能性就越大，数值越低表示该氨基酸不含信号肽或者位于成熟蛋白部分；Y-score是基于C值和S值的斜率得出的几何平均数，是最大可能的信号肽裂解位点；D-score是S均值和Y最大值的加权平均值，用于区分信号肽和非信号肽[23]。C、S、Y的最大值值分别是0.109、0.116、0.105，S-mean和D值分别是0.097、0.101，另外Signal peptide结果显示为No，表明CsbHLH62编码的蛋白不存在信号肽，是一种非分泌蛋白。

图1 氨基酸含量分布Fig.1 Distribution of amino acid content

通过TMHMM分析CsbHLH62蛋白质的跨膜结构，结果如图4所示，CsHLH62编码的1～277位氨基酸都位于细胞膜表面，另外，蛋白N端位于膜内的可能性为0.01881，非跨膜螺旋区预测值是0，该值大于18表示预测蛋白是跨膜蛋白[23]，综合分析表明CsbHLH62基因编码产物不是跨膜蛋白质，在细胞中处于游离状态。采用Cell-PLoc 2.0[24]程序进行亚细胞定位分析，结果表明该基因编码的蛋白定位于细胞核。

图2 CsbHLH62编码蛋白疏水性/亲水性分析Fig.2 Hydrophilia/hydrophobicity analysis of CsbHLH62 encoding protein

图3 CsbHLH62编码产物信号肽预测Fig.3 Signal peptide prediction of CsbHLH62 encoding protein

2.4 CsbHLH62与其他植物间的多序列比对和进化分析

利用NCBI上Blast在线工具，对CsbHLH62编码蛋白进行同源性检索，发现其与胡桃(JuglansregiaXP_018844170.1)、绿豆(Vignaradiatavar.RadiataXP_014496615.1)、大豆(GlycinemaxXP_003532257.1)、桃(PrunuspersicaXP_007201182.1)、菜豆 (PhaseolusvulgarisXP_007140690.1)、短小蛇根草 (OphiorrhizapumilaBAU67840.1)、木豆(CajanuscajanXP_020234271.1)、欧洲越橘(VacciniumcorymbosumAOY34395.1)、枣(ZiziphusjujubaXP_015892293.1)的相似性在50 %～70 %，其中与欧洲越橘的序列同源性最高，相似性为68 %。通过对bHLH62编码氨基酸序列进行同源性分析(图5)，发现该基因在不同植物间结构上保守性较强。

图4 CsbHLH62编码产物跨膜结构域预测Fig.4 Transmembrance domain prediction of CsbHLH62 encoding protein

图5 bHLH62同源基因氨基酸序列对比Fig.5 Homologous alignments of the bHLH62 genes

采用MEGA6.0邻接法构建系统进化树(图6)，结果表明CsbHLH62基因与欧洲越橘bHLH024亲缘关系最近，推测与其具有相似的功能。

2.5 CsbHLH62编码产物的功能结构域和磷酸化位点的预测分析

将CsbHLH62基因编码序列输入ScanProsite在线程序，进行功能结构域搜索，结果表明该基因编码蛋白有1个Myc型的螺旋-环-螺旋(basic helix-loop-helix，bHLH)结构域，该结构域分布在第131～181位氨基酸之间。

数字表示支持率，分枝长短代表亲缘关系远近图6 不同物种的bHLH62编码产物的系统进化树Fig.6 Phylogenetic tree of bHLH62 genes in different plants

应用NetPhos 3.1 Server在线软件对CsbHLH62基因翻译后磷酸化情况进行分析，结果如图所示，该基因编码的蛋白可能发生磷酸化的位点有35个，其中丝氨酸(Ser)、苏氨酸(Thr)和酪氨酸(Tyr)可能发生磷酸化修饰的位点分别有24、10、1个(图7～8)。

2.6 CsbHLH62编码产物的二级和三级结构预测

利用SMOPA程序预测CsbHLH62编码产物的二级结构。结果如图9所示，134个氨基酸形成了无规则卷曲(Random coil)，占整体的48.38 %；86个氨基酸组成α-螺旋(Alpha helia)结构，占整个二级结构的31.05 %；延伸链(Extended strand)和β转角(Beta turn)分别由35、22个氨基酸残基构成，在所占比例分别为12.64 %和7.94 %。以上结果表明无规则卷曲和α-螺旋是CsbHLH62编码产物的主要结构。

图7 CsbHLH62编码产物结构域分析Fig.7 Domain analysis of CsbHLH62 encoding protein

图8 CsbHLH62编码产物磷酸位点分析Fig.8 Phosphorylation sites analysis of CsbHLH62 encoding protein

图9 CsbHLH62编码产物二级结构预测Fig.9 Secondary structure prediction of CsbHLH62 encoding protein

图10 CsbHLH62编码产物三级结构预测Fig.10 3D structure of prediction of CsbHLH62 encoding protein

将 CsbHLH62氨基酸序列提交至Phyre2中进行预测，得到了该基因编码产物的三维模型(图10)。

3 讨论

茶树是一种重要的经济型作物，培育优良茶树品种对于茶产业的发展具有重要的意义。前期通过转录组测序发现一个在根部表达丰度较高的bHLH转录因子，该结果与黄芪响应硒的研究相似[25]，但其在硒吸收及累积过程中的生物学功能还不明确，本研究通过对其进行生物信息学分析，为后续的研究做好准备。

CsbHLH62基因编码是一种可溶性蛋白，无信号肽，二级结构主要由无规则卷曲和α-螺旋组成，这与在苦荞、苹果和丹参等物种中的bHLH转录因子的分析结果一致[26-28]，该基因定位在细胞核中，这与在茶树中另一bHLH转录因子的研究类似[20]。其中二级结构中的无规则卷曲所占比例最大，前人认为其可用于分子识别，在蛋白质特性和功能中发挥关键作用[29]；而对于该基因的三维结构，虽然SWISS-MODEL同源建模应用较广，但可能会出现目标基因预测模型和模板不一致的情况，而Phyre2使用远程同源检测的方法，较好地弥补了这一缺陷[30]；此外，CsbHLH62因含有MYC型的bHLH结构域而归属于bHLH类转录因子家族。以往研究发现MYC类转录因子是植物茉莉酸信号转导途径中的核心转录因子，在植物生长发育、次级代谢和抵御反应中具有重要的作用[31]。比如拟南芥中的MYC2与JA诱导的许多抗虫基因正相关，对提高其抵抗虫害的能力中发挥重要作用[32]。在烟草中MYC类转录因子也参与了尼古丁的生物合成[33]。在本研究中，转录组分析发现茶树受到硒的刺激后，茉莉酸类相关基因的表达发生显著变化，根据bHLH转录因子的结构特点，猜测CsbHLH62是与茉莉酸类相关基因的DNA结合，进而调控其下游的响应基因，该结论尚需进一步的试验验证。