国兰叶色突变体根状茎差异表达基因分析

蒋 彧 何俊蓉

(四川省农业科学院园艺研究所，四川 成都 610000)

兰科是单子叶植物第一大科，包含了5个亚科801个属，品种数量近3万种，约占所有种子植物种类的10%[1]。我国十大名花中的兰花也称国兰，属于兰科(Orchidaceae)兰属(Cymbidium)小花型地生种，在我国有悠久的栽培历史，是世界上最早栽培的兰花。兰花是高产值的花卉，随着我国经济及出口贸易的发展，全球兰花市场逐渐繁荣，因此对兰花的育种工作也提出了更高的要求。

兰花育种的常用手段包括选择育种、杂交育种、倍性育种、体细胞无性系变异以及分子育种等。我国兰花品质相关的育种以选择育种及杂交育种为主，但这两种育种方式的育种周期长、工作量大，近年来也开始采用如倍性育种、体细胞无性系变异等方法。随着测序技术及基因编辑技术的发展，分子育种技术在水稻、玉米、棉花等农作物中被广泛使用，并取得了很好的育种效果[2]。兰花分子育种起步较晚，兰科植物如小兰屿蝴蝶兰(Phalaenopsisequestris)[3]、铁皮石斛(Dendrobiumofficinale)[4]、深圳拟兰(Apostasiashenzhenica)[5]、白花蝴蝶兰(Phalaenopsisaphrodite)[6]的基因组发布将有效促进兰花分子育种工作的开展。近年来，基于多组学数据的分析，兰花重要性状相关基因的挖掘取得了重要进展。兰花中受低温抑制的CgSVP[分离自春兰(cymbidumgoeringii)]参与了兰花低温开花调控，该调控机制与模式植物中的低温响应机制存在显著差异[7]。兰花中的MADS-box基因家族参与了国兰特异性唇瓣分化和重瓣化的调控，通过表达区域扩增产生花型的变异[8-9]。通过CRISPR/Cas9基因编辑技术，研究者也成功获得了兰科植物铁皮石解C3H、C4H[10]以及蝴蝶兰MADS8/36/44的突变株系[11]。上述研究表明分子育种技术在兰花育种中能够发挥重要作用。

叶艺是国兰特有的一种体征，表现为叶片呈现金银色的条纹或斑点，具有较高的观赏价值，叶艺兰的价格也在逐步提高，其消费市场具有巨大潜力。目前市场上销售的叶艺兰品种大部分属于野生的变异种(也称为下山兰，即直接从山上挖回的兰花)，另有少部分来自于栽培品种芽变，但前者会对兰花野生资源造成破坏，后者变异率较低。人工诱变技术较栽培品种芽变具有更高的变异率，对兰花的育种工作起到了重要的推动作用。叶艺突变体成因复杂，目前的研究认为叶绿体发育缺陷和叶绿素过量降解是叶艺性状产生的直接生理原因[12]。墨兰叶片颜色在由红变黄，最后转化为绿色的过程中，有28种代谢物在叶片变绿后含量降低到近零水平，还有15种花色素苷类代谢物含量呈下降趋势[13]。CsERF2则通过乙烯信号途径参与调控墨兰(Cymbidumsinense)叶色的变化[14]。但总的来说，导致国兰叶艺表型背后的分子机制还有待进一步研究。

研究认为，60Co-γ射线可有效地对植物进行诱变产生突变体[15]，通过60Co-γ射线对隆昌素根状茎进行诱变后可获得能产生叶艺苗的叶艺隆昌素。隆昌素与叶艺隆昌素根状茎分别表现为绿色和黄色，具有非常明显的差异，这两种材料中的叶绿素与类胡萝卜素含量有较大差别[2]。本研究采用RNA-Seq技术对这两种根状茎进行分析，鉴定在隆昌素与叶艺隆昌素根状茎中差异表达的基因，旨在为解析叶艺隆昌素叶艺成因的分子机制奠定研究基础。

1 材料与方法

1.1 试验材料与试剂

试验材料来自四川省农业科学院园艺研究所组培室。培养温度为23℃(白天)/18℃(晚上)，光照周期为14 h(光照)/10 h(黑暗)，光照强度为2 000 lx。无菌条件下随机选取叶艺隆昌素(YR)和隆昌素(GR)根状茎各4瓶，分别从每个培养瓶中随机挑选1～2 cm根状茎各4条，置于液氮中冷冻后，-80℃超低温冰箱贮存备用。

多糖多酚植物总RNA提取试剂盒(离心柱型)，北京天根生化科技有限公司；mRNA富集用的磁珠mRNA Capture BeadsRNA，南京诺唯赞生物科技有限公司；逆转录试剂盒 TUREscript 1st Strand cDNA Synthesis Kit，北京艾德莱生物科技有限公司；2×SYBR®Green Supermix实时荧光定量PCR试剂盒，德国DBI®Bioscience。

1.2 仪器与设备

Eppendorf高速冷冻离心机(Eppendorf，德国)；Nanodrop2000超微量分光光度计(Thermo Scientific，美国)；qTOWER2.2实时PCR系统(DBI®Bioscience，德国)

1.3 试验方法

1.3.1 总RNA提取 根据植物总RNA提取试剂盒说明书，从根状茎组织样品中提取总RNA，利用Nanodrop2000超微量分光光度计检测RNA的浓度和纯度，1.5%琼脂糖凝胶电泳检测RNA完整性。

1.3.2 cDNA文库的构建 使用Oligo dT磁珠法从总RNA中分离出mRNA，用超声波破碎法将mRNA随机断裂成200 bp左右的小片段。在逆转录酶的作用下合成第一链cDNA，在DNA聚合酶作用下进一步合成第二链cDNA。

1.3.3 转录组测序及差异表达基因筛选 将构建好的cDNA送至上海美吉生物医药科技有限公司使用Illumina Hiseq2500进行转录组测序，测序方式为双末端法，测序长度为150 bp。测序数据经质控，截除测序接头、引物序列，并对低质量值数据进行过滤，最终获得高质量数据(reads)。使用Trinity软件[16]对高质量Reads进行denovo拼接获得单基因(unigene)。采用RSEM软件(V 1.2.4)中的FPKM(fragments per kilobase per million reads)方法对unigene表达量水平进行计算，使用Bioconductor软件包中的DESeq(V 1.14.0)进行基因的差异表达分析。差异表达基因(different expression genes, DEGs)筛选阈值为|log2 value|>1且FDR≤0.05。

1.3.4 unigene功能注释及基因结构分析 使用BlastX将unigene与NR、GO、KOG数据库进行比对，根据基因的相似性进行功能注释，得到与给定unigene具有最高序列相似性的蛋白，从而得到该unigene的蛋白功能注释信息。使用TransDecoder软件(transdecoder.github.io/)对unigene的编码区及其对应的氨基酸序列进行预测。

1.3.5 差异表达基因富集分析 分别使用GOseq与KOBAS对差异表达基因进行GO功能显著性及KEGG Pathway富集分析。

1.3.6 实时荧光定量PCR(quantitative real-time PCR, qRT-PCR) 提取根状茎RNA，用Aidlab公司反转录试剂盒(TUREscript 1st Strand cDNA Synthesis Kit)进行cDNA的合成。从转录组数据中筛选出准备验证的unigene,分别设计这些unigene的特异引物，采用qTOWER2.2实时PCR系统进行qRT-PCR分析，所有反应均设置3个生物学重复和技术重复。qRT-PCR反应体系为10 μL：5 μL 2×SYBR®Green Supermix，0.5 μL 200 nmol·L-1正向引物，0.5 μL 200 nmol·L-1反向引物，1 μL cDNA模板，3 μL ddH2O。qRT-PCR反应程序为：95℃预变性3 min；95℃变性10 s，60℃扩增30 s，循环39次；溶解曲线分析，即随后将温度从60℃逐渐升至95℃(每次提高1℃，每提高1℃停留4 s)。ACTIN-7基因作为内参，通过2-△△Ct法计算目的基因相对表达量。

2 结果与分析

2.1 叶艺隆昌素和隆昌素根状茎形态

叶艺隆昌素根状茎呈淡黄色，形状较为粗短(图1-A)；隆昌素根状茎则呈深绿色，形状为不规则长条形(图1-B)。经过培养后，淡黄色的根状茎分化出叶艺苗，绿色的则只分化成正常的绿色试管苗。

图1 叶艺隆昌素根状茎(A)和隆昌素根状茎(B)观察Fig.1 Observation of Yeyi Longchangsu (A) and Longchangsu (B) rhizome

表1 测序结果Table 1 Sequencing statistics

2.2 测序数据组装与分析

分别对叶艺隆昌素和隆昌素的根状茎进行转录组测序，每个样本分别进行4次重复，8个文库得到的原始数据(raw reads)条数介于13 193 936～20 145 078，质控后的高质量数据(clean reads)条数介于10 686 612～16 676 096，clean reads占raw reads的比例为70.15%～84.12%，clean reads数据量大小为3.96～6.08 Gb(表1)。使用Trinity和TGICL软件对clean reads进行从头组装，共得到162 407条unigene(表2)。unigene的N50长度和平均长度分别为1 049和678 bp。unigene长度分布结果显示大于1 000 bp的unigene有26 107个(图2)。

2.3 基因功能注释及分类

基因功能注释与分类为阐明细胞内特定基因的分子功能和生物代谢途径提供了有价值的信息。将162 407 条unigene与Nr、GO和KOG数据库进行比对，获得unigene的注释信息。结果显示有26 473条unigene注释到Nr数据库，有20 119条unigenes注释到KOG数据库，有38 926条unigene注释到GO数据库(图3)，其中有14 380条unigene在3个数据库中均有注释。在GO部分，注释到细胞组分、分子功能、生物学过程的unigene数目分别为12 114、13 170、13 642 条。

表2 转录组数据的组装及unigene预测Table 2 Statistics of transcriptome assembly and predicted unigenes

图2 unigene长度分布图Fig.2 Length distribution of unigene

图3 GO、Nr和KOG数据库注释结果统计Fig.3 GO, Nr and KOG database annotation results statistics

图4 相关性热图(A)和差异表达基因火山图(B)Fig.4 Correlation heat map (A) and volcano map of DEGs (B)

2.4 差异表达基因的鉴定及功能分类

在阈值为|log2 value|>1和FDR≤0.05情况下，鉴定出644个unigene在叶艺隆昌素和隆昌素根状茎中差异表达，其中497个unigene的转录本在叶艺隆昌素中表达丰度较高，147个unigene的转录本在隆昌素中表达丰度较高(图4)。

对差异表达基因进行GO聚类，结果表明有42个类别显著富集，其中属于生物学过程的有20个，属于细胞组分的有11个，属于分子功能的有11个(图5)。在生物学过程类别中，代谢过程富集度最高，节律过程富集度最低；在细胞组分类别中，细胞富集度最高，膜封闭腔富集度最低；在分子功能类别中，催化活性富集度最高，金属伴侣活性富集度最低。

对差异表达基因进行KEGG聚类，结果显示有8个代谢过程被显著富集(P<0.05)(表3)。这些过程涉及植物-病原物互作(P=0.000 26)、矿物质吸收(P=0.001 34)、光形态建成(P=0.001 99)、倍半萜和三

注: 1: 生物调控; 2: 细胞成分组织或起源; 3: 细胞过程; 4: 发育过程; 5: 定位系统的建立; 6: 生长; 7: 免疫系统过程; 8: 定位; 9: 代谢过程; 10: 多机体过程; 11: 多细胞组织过程; 12: 生物过程负调控; 13: 生物过程正调控; 14: 生物过程调控; 15: 再生; 16: 再生过程; 17: 刺激应答; 18: 节律过程; 19: 信号; 20: 单一的生物过程; 21: 细胞; 22: 胞间连丝; 23: 细胞要素; 24: 细胞间区域; 25: 大分子复合物; 26: 膜; 27: 膜要素; 28: 膜结合腔体; 29: 细胞器; 30: 细胞器要素; 31: 共质体; 32: 抗氧化剂活性; 33: 结合剂活性; 34: 催化剂活性; 35: 电荷载体活性; 36: 酶调控因子活性; 37: 金属伴侣活性; 38: 分子感应器活性; 39: 核苷酸结合转录因子活性; 40: 蛋白质结合转录因子活性; 41: 结构分子活性; 42: 转运因子活性。Note: 1: Biological regulation. 2: Cellular component organization or biogenesis. 3: Cellular process. 4: Developmental process. 5: Establishment of localization. 6: Growth. 7: Immune system process. 8: Localization. 9: Metabolic process. 10: Multi-organism process. 11: Multicellular organismal process. 12: Negative regulation of biological process. 13: Positive regulation of biological process. 14:Regulation of biological process. 15: Reproduction. 16: Reproductive process. 17: Response to stimulus. 18: Rhythmic process. 19: Signaling. 20: Single-organism process. 21: Cell. 22: Cell junction. 23: Cell part. 24: Extracellular region. 25: Macromolecular complex. 26: Membrane. 27: Membrane part. 28: Membrane-enclosed lumen. 29: Organelle. 30: Organelle part. 31: Symplast. 32: Antioxidant activity. 33: Binding activity. 34: Catalytic activity. 35: Electron carrier activity. 36: Enzyme regulator activity. 37: Metallochaperone activity. 38: Molecular transducer activity. 39: Nucleic acid binding transcription factor activity. 40: Protein binding transcription factor activity. 41: Structural molecule activity. 42: Transporter activity.图5 不同表达基因的GO注释分析Fig.5 The GO annotation statistical of DGEs

表3 DEG基因KEGG分析Table 3 The KEGG path way of DEG

萜生物合成(P=0.047 25)。

2.5 叶绿素与类胡萝卜素相关基因的鉴定及表达分析

叶绿素、类黄酮、类胡萝卜素是决定植物色泽的三大主要物质。本试验通过unigene注释结果进行筛选，分别获得了15编码参与叶绿体与12个编码参与类胡萝卜素生物合成的关键酶的基因。对15个叶绿素相关基因的表达分析发现，其中5个基因在两种实验材料中表达水平基本一致，10个基因表达水平存在一定差异，但未达到显著水平(表4)。与类胡萝卜素相关的12个unigene中(表5)，编码类胡萝卜素解离双加氧酶7(CDD7)的基因c1920_g1在叶艺隆昌素中表达下调。

表4 叶绿素相关unigene鉴定Table 4 Identification of chlorophyll related unigene

表5 类胡萝卜素相关unigene鉴定Table 5 Identification of carotenoid related unigene

2.6 叶绿体发育相关基因的鉴定及表达分析

从转录组测序结果中筛选出一些编码光系统I/II、 细胞色素b6f复合物、ATP合酶的unigene，其中光系统Ⅱ中的psbC、psbH以及光系统I中的Ycf3基因在叶艺隆昌素根状茎中均表现为表达水平下调(表6)。叶绿体自身编码的第二类基因是和叶绿体DNA复制与叶绿体基因表达相关的基因。DEG分析结果显示，编码RNA聚合酶亚基的glgC、rbcL以及编码核糖体蛋白亚基的rpoB、rpsR基因在叶艺隆昌素根状茎中均表现为表达下调(表7)。

2.7 qRT-PCR验证

为了验证转录组数据的可靠性，本试验选择了9个涉及叶绿素生物合成、光合系统和叶绿体发育的基因设计引物(表8)进行qRT-PCR反应，对转录组分析

表6 光合作用相关的叶绿体编码基因筛选Table 6 Identification of chloroplast structure related unigenes

表7 叶绿体结构相关的基因筛选Table 7 Identification of chloroplast structure related unigenes

结果进行验证，如图6所示。qRF-PCR结果显示验证基因的表达趋势均与转录组数据相符合，表明转录组数据具有准确性与可信性。

2.8 其他调控相关的差异基因

除了上述与叶绿素或类胡萝卜素相关的差异基因，转录组测序结果还获得了大量与植物激素信号传导、蛋白转运、编码细胞色素P450蛋白、次生物质代谢、内源激素相关的差异表达基因(表9)。这些基因也可能影响叶艺隆昌素的叶色表型。

3 讨论

目前已有研究对兰科植物进行了基因组测序，其花色变异及分子标记得到了部分解析[8,17-18]，但大多数兰花物种的遗传信息尚不清楚。转录组代表组织特定发育时期或阶段的所有转录本的总和，对于从宏观层面揭示分子调控网络有着独特的优势[19-21]。隆昌素和叶艺隆昌素试管苗根状茎的颜色呈现明显不同，前者呈墨绿色，后者呈黄色。通过对叶艺兰分化的研究发现，黄色根状茎可以分化出叶艺苗，绿色根状茎则分化正常绿色试管苗，是研究国兰叶艺形成机制的重要材料[2]。本研究通过转录组分析鉴定到644个在隆昌素和叶艺隆昌素根状茎中差异表达的unigene，通过对这些差异表达基因的筛选，初步分析了影响叶艺隆昌素根状茎变黄的因素,为解析国兰叶艺形成分子机制提供了一定的科学依据。

植物的色泽与类胡萝卜素、类黄酮、叶绿素含量有关[22]，植物叶绿素的合成过程要经历大约16个步骤，其中以谷氨酰t-RNA为起点，过程中涉及16种酶和20多个基因[23]，这些基因会影响叶绿素合成，进而使叶绿体中各种色素的比例与含量发生变化，最终改变叶色。水稻中已有35个叶色相关基因得到克隆，如水稻ygl1黄绿叶突变体幼叶中叶绿体因为叶片内四吡咯中间体含量的变多而发育迟缓，导致其叶绿素含量远低于野生型[24]。相关研究结果显示，国兰叶片颜色变异由叶绿素代谢途径异常而非合成途径异常导致[12]。本研究中筛选出的编码叶绿素合成相关酶的基因在叶艺隆昌素与隆昌素根状茎中也未检测出显著差异，与前人研究结果相似。

注:*表示具有显著差异(P<0.05), ** 表示具有极显著差异(P<0.01)。Note: * and ** indicate significant at 0.05 level and very significant difference at 0.01 level, respectively.图6 unigene qRT-PCR验证结果Fig.6 qRT-PCR identification of unigenes

表8 qRT-PCR引物设计Table 8 qRT-PCR primers

CCD7编码一个ABA生物合成途径中的关键限速酶。研究表明，植物体内源ABA的变化与CCD7表达变化存在显著相关性[25]。另一方面，花瓣颜色的形成与CCDs的表达水平存在显著相关性[25]，菊花白色花瓣中CCD4表达水平远高于黄色花瓣，CCD4沉默后的白色花瓣会转变为黄色，原因是白色花瓣合成类胡萝卜素之后进行了降解。在叶艺隆昌素根状茎内，CCD7的表达量下调，可使胡萝卜素降解速度变快，从而使得根状茎颜色发黄，进一步说明CCD的表达参与了颜色调控。

捕光色素蛋白复合体由类囊体蛋白与分布在叶绿体类囊体膜上的光合色素组成，这些色素蛋白复合体均参与了光反应，起到转换能量与传递电子的作用[26]。本研究中光系统Ⅱ中的psbC、psbH以及光系统Ⅰ中Ycf3基因在叶艺隆昌素根状茎中表达量下调。这些基因的表达下调可能影响了叶绿体基粒片层结构，造成光合电子传递速率下降，光合能力下降，进而影响根状茎增殖及分化效率。光合作用碳固定和光呼吸碳氧化过程中的关键酶为1,5-二磷酸核酮糖羧化酶/加氧酶(Rubisco)，其大亚基编码通过叶绿体基因实现，是植物叶片内的关键储存蛋白。有研究显示，水稻(OryzasativaL.)幼苗白化大多由Rubisco小亚基因的突变导致[27-28]。编码Rubisco的unigene在叶艺隆昌素内的表达量会减少也可能是导致叶艺隆昌素根状茎变黄的原因之一。前人研究发现，大亚基与小亚基是核糖体的主要组成要素，烟草(NicotianatabacumL.)叶片因为缺失叶绿体核糖体小亚基编码基因rps18而出现了畸形并变白的现象，存活率极低；Fleischmann等[29]通过研究发现，叶片形态会因为缺失核糖体大亚基编码基因rpl36而发生变化，光合速率和植株生长速

表9 其他相关基因筛选Table 9 Identification of others related unigenes

表9 (续)

度明显变慢，本研究中编码核糖体蛋白大小亚基的rpoB、rpsR基因叶艺隆昌素根状茎均表达下调，而这些基因是否突变或缺失引起叶绿体结构破坏，有待进一步研究确定。

ABC转运蛋白参与了细胞器(如液泡、过氧化物体、线粒体、叶绿体等)间次生代谢产物、植物激素、蛋白质、氨基酸、糖类等物质的转运，是一种跨膜运输蛋白，存在于大部分生物膜上[30-31]。本研究中，ABC转运蛋白家族基因(ABCB1与ABCG2)表达下调，由此推测由于叶艺隆昌素根状茎颜色变黄、叶绿体结构破坏，导致光合作用减弱、产生的次生代谢物质减少，从而通过信号转导途径反馈抑制ABC转运蛋白基因的表达。线粒体与叶绿体在植物细胞中都具有重要作用，线粒体内膜中会出现氧化磷酸化的现象，在ATP合酶催化作用下，电子在线粒体内膜上一系列电子传导作用下传递到氧，促使ADP和磷酸合成ATP，线粒体内膜上的ATP合成酶、细胞色素氧化酶、辅酶Q-细胞色素bc1还原酶、琥珀酸脱氢酶、NADH脱氢酶五种蛋白复合体都属于电子传递体，其中线粒体电子传递链内以NADH脱氢酶和NADH-泛醌氧化还原酶最为重要，电子传递速率受到这些酶活性的影响。本试验中，NADH脱氢酶、NADH-泛醌氧化还原酶和细胞色素氧化酶均表达下调，可以推测由于叶艺隆昌素根状茎叶绿体结构不完整，导致体内新陈代谢减弱，通过信号传导途径抑制NADH脱氢酶、NADH-泛醌氧化还原酶和细胞色素氧化酶基因的表达，减少ATP的合成。作为折叠水解酶超家族成员的BDG，其主要在皮细胞中表达，拟南芥纯合突变植株矮小、叶片异常、表皮毛细胞数不断变小且细胞萎陷，存在许多畸型细胞，植株发育缓慢，推测其表皮细胞组织是在某种未知调节机制的作用下实现分化和增殖的[32]。本研究中，叶艺隆昌素根状茎类联苯水解酶蛋白显著上调，但该蛋白的作用机制尚不清楚，可知其功能可能与角质的生物合成或与细胞壁本身交联酶相关，同时该蛋白是否与根状茎颜色变黄相关有待进一步证实。

4 结论

本研究采用转录组分析发现，隆昌素与叶艺隆昌素根状茎差异基因在植物-微生物互作、矿物质吸收、光形态建成与倍半萜和三萜生物合成通路中显著富集。此外，参与ABA合成的CDD7、光系统Ⅰ中的Ycf3、光系统Ⅱ中的psbC与psbH、编码RNA聚合酶亚基的glgC与rbcL、编码核糖体蛋白亚基的rpoB与rpsR等基因表达量的改变可能通过影响植株叶绿体发育从而影响植株的叶艺表型。