蛹虫草转录组分析及类胡萝卜素生物合成相关基因的挖掘

娄海伟，赵 玉，赵 逸，林俊芳,3,*，赵仁勇，叶志伟,3，郭丽琼,3,*

（1.河南工业大学粮油食品学院，河南 郑州 450001；2.华南农业大学食品学院，广东 广州 510642；3.广东省微生态制剂工程技术研究中心，广东 广州 510642）

蛹虫草（Cordyceps militaris）是一种珍贵的食药两用真菌，含有多种生物活性成分，如虫草素[1]、免疫调节蛋白[2]、虫草多糖[3]、喷司他丁[4]等，使蛹虫草具有抗炎[5]、增强免疫[6]、抗肿瘤[7]、抗癌[8]等功效。2014年，蛹虫草被国家卫生计生委批准为新食品原料。这使得蛹虫草的应用领域和市场前景更加广阔。

类胡萝卜素作为一类重要的萜类色素，广泛分布于光合细菌、真菌、植物和藻类中[9]，具有抗氧化[10]、抗炎[11]、抗癌[12]、预防肥胖[13]等生理功能，是维持人体健康不可或缺的生物活性物质。近年来，在蛹虫草中发现了新型的类胡萝卜素，该新型色素具有良好的亲水性和抗氧化活性[14-16]，这亦使得蛹虫草的市场需求量逐年增多。但是，蛹虫草中新型类胡萝卜素的含量较低，不能满足市场的需求，主要是因为其生物合成途径未知，生物合成相关基因的研究报道较少[17]。前人的研究结果表明，培养基显著影响蛹虫草类胡萝卜素的生物合成[18-19]。基于此，本研究旨在对不同培养基培养的蛹虫草进行转录组测序和生物信息学分析，挖掘蛹虫草类胡萝卜素的生物合成相关基因，这将有利于揭示蛹虫草类胡萝卜素的生物合成途径，为通过生物技术方法提高蛹虫草类胡萝卜素的产量奠定基础。

1 材料与方法

1.1 材料与试剂

1.1.1 菌株

蛹虫草菌株CM10，4 ℃保存于PDA斜面培养基上。

1.1.2 培养基

马铃薯葡萄糖水（potato dextrose broth，PDB）培养基：马铃薯200.0 g、葡萄糖20.0 g、MgSO41.5 g、KH2PO43.0 g，自来水定容至1 000 mL。

马铃薯葡萄糖琼脂（potato dextrose agar，PDA）培养基：PDB培养基中加入20.0 g琼脂，自来水定容至1 000 mL。

蛹虫草生长所需培养基的最适碳氮比（C/N）为20∶1[20]，根据简锦辉[18]的研究结果，蛹虫草在（C/N为20∶1）培养基W上产类胡萝卜素的量较高，而在（C/N为41.78∶1）对照培养基R上产类胡萝卜素的量较低，且蛹虫草在上述2 种培养基上产类胡萝卜素的量差异显著。因此，本研究选择了培养基W和培养基R。

培养基W：小麦14.6 g、麸皮4.0 g、玉米1.4 g，粒度均为10 目，合计20.0 g放入一个直径为9 cm的培养皿，加入25.0 mL营养液W。营养液W配方为：MgSO41.0 g、KH2PO41.0 g、VB140 mg，自来水定容至1 000 mL。

培养基R：大米20.0 g（粒度10 目）放入一个直径9 cm的培养皿，加入25.0 mL营养液R。营养液R配方为：MgSO41.0 g、KH2PO41.0 g、VB140 mg、柠檬酸1.0 g、葡萄糖12.0 g、胰蛋白胨8.0 g，自来水定容至1 000 mL。

1.1.3 试剂

葡萄糖、MgSO4、KH2PO4、柠檬酸 广州化学试剂厂；胰蛋白胨 北京奥博星生物技术有限公司；VB1上海源叶生物科技有限公司；TRIzol试剂 美国Thermo Fisher Scientific公司；Fungal RNA Kit 北京Omega Bio-tek公司；TransScript All-in-One First-Strand cDNA Synthesis SuperMix for qPCR、TransStart Tip Green qPCR SuperMix 北京全式金生物技术有限公司。

1.2 仪器与设备

LRH-70生化培养箱 上海一恒科学仪器有限公司；GXZ智能光照培养箱 宁波江南仪器厂；YX-280手提式压力蒸汽灭菌器 合肥华泰集团股份有限公司；SW-CJ-1FD洁净工作台 苏州安泰空气技术有限公司；MQD-S2R振荡培养箱 上海旻泉仪器有限公司；U410超低温冰箱 德国Eppendorf公司；QuantStudio 3荧光定量聚合酶链式反应（polymerase chain reaction，PCR）仪美国Thermo Fisher Scientific公司。

1.3 方法

1.3.1 菌丝制备和类胡萝卜素含量的测定

取4 ℃保存的斜面菌种CM10接种至PDA培养基上，25 ℃避光培养21 d，然后光照培养7 d，即得活化的蛹虫草菌种。在洁净工作台中把活化的菌种CM10用手术刀切成1 cm×1 cm大小的菌块，每100 mL PDB培养基接种3 块CM10菌块，25 ℃避光振荡培养（150 r/min）4 d，即得蛹虫草液体菌种。接种培养基W和培养基R，每个培养皿接种5 mL液体菌种，在相同培养条件（25 ℃避光培养10 d，然后25 ℃光照培养8 d）下培养，分别得到CM10_WL和CM10_RL菌丝样品。每个样品进行6 个生物学重复，采用液氮收集菌丝样品，同一样品的6 个重复样混合，存放于-80 ℃超低温冰箱中备用。取出一部分菌丝进行冷冻干燥，根据Zhang Jiaojiao等[21]报道的方法测定冻干菌丝中类胡萝卜素的含量。

1.3.2 cDNA文库构建和转录组测序

-80 ℃保存的样品送至华大基因有限公司提取样品的总RNA，根据NEB Next Ultra RNA Library Prep Kit试剂盒说明书构建CM10_WL和CM10_RL样品的cDNA文库，然后采用Illumina HiSeq 4000测序平台进行测序。原始测序数据经质控以除去低质量的reads（接头reads、未知碱基N含量大于5%的reads），最后得到clean data，已存储于美国国家生物技术信息中心（National Center of Biotechnology Information，NCBI）的生物项目数据库（BioProject）检索号PRJNA541819。

1.3.3 基因组比对和差异表达基因（differentially expressed gene，DEG）分析

使用HISAT软件将每个样品的clean reads与蛹虫草基因组序列（NCBI reference sequence: NZ_AEVU00000000.1）进行比对[22-23]。以样品CM10_RL作为对照样品，与样品CM10_WL进行基因表达量的比较。将错误发现率（false discovery rate，FDR）≤0.001且差异倍数（fold change）≥2作为DEG的筛选标准，满足此筛选标准的基因即为DEG。

1.3.4 DEG的功能注释

将筛选到的所有DEG分别与NCBI非冗余蛋白数据库（Non-redundant Protein Sequence Database，NR）、京都基因与基因组百科全书（Kyoto Encyclopedia of Genes and Genomes，KEGG）数据库、基因本体论（Gene Ontology，GO）数据库中的蛋白序列进行比对，从而获得与DEG对应的蛋白功能分类统计和功能注释。

1.3.5 实时PCR（real-time PCR）分析

为了验证转录组数据的可靠性，随机挑选4 个基因（CCM_05119、CCM_06728、CCM_07824、CCM_09155），采用real-time PCR检测上述4 个基因在样品CM10_RL和CM10_WL中的表达情况。采用Fungal RNA Kit提取总RNA，采用TransScript All-in-One First-Strand cDNA Synthesis SuperMix for qPCR（One-Step gDNA Removal）试剂盒合成第1链cDNA。以tef1基因（GenBank: DQ070019）作为内参基因[24-25]，real-time PCR反应体系和反应条件参见TransStart Tip Green qPCR SuperMix试剂盒说明书，每个样品进行3 次生物学重复。以CM10_RL作为对照样品，采用2−ΔΔCt法计算基因的相对表达量[26]。所用引物见表1。

表 1 real-time PCR引物Table 1 Sequences of primers used for real-time PCR

2 结果与分析

2.1 蛹虫草菌株的活化和液体菌种的制备

蛹虫草CM10菌株在PDA培养基上活化结果见图1a，其菌丝为橙黄色，表明蛹虫草产橙黄色色素[27]。制备的液体菌种见图1b，蛹虫草呈米粒状的菌丝球，PDB培养基澄清，液体菌种具有蛹虫草特有的清香。

图 1 蛹虫草CM10的活化菌种（a）和液体菌种（b）Fig. 1 C. militaris CM10 activated on PDA plate (a) and its cell suspension (b)

2.2 蛹虫草菌丝的制备

在相同培养条件下、不同培养基上培养的蛹虫草菌丝样品CM10_RL和CM10_WL见图2，样品间色泽表型差异显著，样品CM10_RL呈较淡的橙黄色，而样品CM10_WL呈较深的橙红色。

如图3所示，样品CM10_WL的类胡萝卜素含量显著高于样品CM10_RL的类胡萝卜素含量，表明培养基显著影响蛹虫草类胡萝卜素的产生[19,28]，这可能是由于培养基中C/N的差异引起的，C/N为20∶1的培养基W更适合蛹虫草产类胡萝卜素，而C/N偏离20∶1的培养基R（C/N为41.78∶1）不利于蛹虫草产类胡萝卜素，这与前人的研究结果一致[18,29]。这亦可能是导致样品CM10_WL较CM10_RL的色泽更深的主要原因。

图 3 蛹虫草菌丝中类胡萝卜素的含量Fig. 3 Contents of carotenoids in C. militaris mycelia

2.3 蛹虫草转录组测序数据统计分析和参考基因组比对

表 2 蛹虫草菌丝转录组测序统计Table 2 Transcriptomic sequencing statistics of C. militaris mycelia

由表2可知，样品CM10_RL和CM10_WL分别获得32 388 522 个和32 387 910 个raw reads，质控后，样品CM10_RL和CM10_WL分别获得30 595 116 个和30 261 506 个clean reads，2 个样品质控后的碱基数均大于4.5 Gb，测序数据量足够大，能准确反映碱基序列，Q20均大于97%，Q30均大于94%，测序数据的质量较高，可用于后续分析。

将质控后的clean reads与参考基因组比对，样品CM10_RL和CM10_WL的比对率分别为79.83%和79.94%，样品间均匀的比对率表明样品之间的数据具有可比性。同时，从样品CM10_RL中检测到8 923 个表达基因（包括8 786 个已知基因和137 个新基因），从样品CM10_WL中检测到8 807 个表达基因（包括8 670 个已知基因和137 个新基因），与已知报道的蛹虫草蛋白编码基因数量9 684 个很接近[22]。

2.4 蛹虫草DEG检测和功能注释分析

以样品CM10_RL作为对照样品，检测到样品CM10_WL有936 个DEG，其中318 个基因上调表达和618 个基因下调表达。

图 4 DEG的GO富集柱状图Fig. 4 GO enrichment analysis of DEGs

由图4可知，共有995 个DEG归类到生物过程（453 个DEG）、细胞成分（176 个DEG）、分子功能（366 个DEG）；在生物过程类中，富集最多DEG的条目是代谢过程；在细胞组分类中，富集最多DEG的条目是细胞膜；在分子功能类中，富集最多DEG的条目是催化活性。因此，富集到“代谢过程”、“细胞膜”、“催化活性”条目中的DEG可能负责蛹虫草色素的生物合成。

不同的基因在不同的生物过程中常具有协同作用，KEGG通路富集分析有助于进一步了解基因的生物学功能[30]。DEG的KEGG通路富集结果见图5，图中的圆点越大，说明富集的DEG越多，Q值越小表明富集越显著[31]。由图5可以看出，“代谢途径”富集的DEG最多，这些DEG可能是造成样品CM10_WL产生更多类胡萝卜素的原因，此结果亦表明培养基对代谢途径影响最大，与前人的研究结果一致[32]。

图 5 DEG的KEGG pathway富集Fig. 5 KEGG pathway enrichment analysis of DEGs

2.5 蛹虫草类胡萝卜素生物合成相关基因的挖掘

前人研究报道蛹虫草黄色素为类胡萝卜素[14,33]，因此本研究所得的所有转录本与KEGG通路数据库中的类胡萝卜素生物合成途径（map00906）进行比对，结果发现蛹虫草CCM_06728、CCM_09155基因与类胡萝卜素生物合成途径相关，参与了由圆酵母素到链孢霉黄素的生物合成（图6），圆酵母素在圆酵母素双加氧酶（CCM_06728）的作用下生成β-阿朴-4ʹ-胡萝卜醛，β-阿朴-4ʹ-胡萝卜醛在胡萝卜醛加氧酶（CCM_09155）的作用下生成链孢霉黄素。这2 个类胡萝卜素生物合成基因（CCM_06728和CCM_09155）在样品CM10_RL和CM10_WL中表达差异不显著，这可能是因为它们的表达不受培养基差异的影响。前人的研究结果表明，基因CCM_06728和CCM_09155的表达水平也不受光照因素的影响，而且这2 个基因分布在类胡萝卜素生物合成途径的其中一个支路的末端，对整个类胡萝卜素生物途径的影响较小[17]。

图 6 蛹虫草CCM_06728和CCM_09155基因参与类胡萝卜素生物合成的途径Fig. 6 CCM_06728 and CCM_09155 genes involved in carotenoid biosynthesis in C. militaris

类胡萝卜素生物合成途径的第一步是异戊烯基焦磷酸和二甲基丙烯基焦磷酸在香叶基香叶基焦磷酸合成酶（geranylgeranyl pyrophosphate synthase，GGPPS）催化下生成香叶基香叶基焦磷酸（geranylgeranyl pyrophosphate，GGPP）[34]，而香叶基香叶基焦磷酸合成酶基因GGPPS被认为是蛹虫草类胡萝卜素生物合成途径中的第1个关键酶基因[35]。蛹虫草GGPPS基因（CCM_03697、CCM_03059、CCM_06355）已被成功克隆，且基因CCM_03059被认为参与蛹虫草类胡萝卜素的生物合成[36]。但是基因CCM_03059在样品CM10_RL和CM10_WL中的表达水平无显著差异，这表明基因CCM_03059的表达不受培养基差异的影响，且GGPPS不是蛹虫草类胡萝卜素生物合成过程中的限速酶。GGPP在八氢番茄红素合成酶的作用下生成八氢番茄红素，八氢番茄红素在八氢番茄红素脱氢酶的作用下生成番茄红素，番茄红素通过去饱和、环化、羟基化等反应进入类胡萝卜素生物合成的分支途径。但是，在样品CM10_RL和CM10_WL的转录组中均未发现八氢番茄红素合成酶基因和八氢番茄红素脱氢酶基因，前人亦采用转录组测序技术研究了蛹虫草的转录组，也未发现上述2 个关键酶基因[17,37]。有关蛹虫草类胡萝卜素生物合成途径的研究，目前仅有一条假设的蛹虫草类胡萝卜素生物合成途径，在这条假设的生物合成途径中，基因CCM_00310、CCM_03203、CCM_03697、CCM_04994、CCM_05367、CCM_05380、CCM_06728、CCM_07659、CCM_07937、CCM_08582、CCM_09155被认为参与蛹虫草类胡萝卜素的生物合成，但是这11 个基因在样品498L（光照条件下培养的菌丝）和498D（避光条件下培养的菌丝）中的表达水平无显著差异[17]，而且这11 个基因在样品CM10_RL和CM10_WL中的表达水平也无显著差异，这表明上述11 个基因的表达水平不受光照因素和培养基差异的影响，因此上述11 个基因的功能有待鉴定。

由于蛹虫草含有新型的色素[14-15]，且被认为与已知的粗糙脉孢菌中的类胡萝卜素生物合成途径完全不同[37]，又基于类胡萝卜素属于萜类化合物，结合GO功能富集和KEGG通路富集分析，从DEG中挖掘到5 个参与萜类生物合成的基因CCM_03738、CCM_06092、CCM_05423、CCM_08567、CCM_09041。相对于样品CM10_RL，基因CCM_03738、CCM_06092、CCM_05423在样品CM10_WL中上调表达，表明这3 个基因促进蛹虫草色素的合成；而基因CCM_08567和CCM_09041在样品CM10_WL中下调表达，表明这2 个基因抑制蛹虫草色素的合成。基因CCM_03738的注释结果为4-羟基苯甲酸聚戊烯基转移酶，该酶的功能是催化聚戊烯基二磷酸和4-羟基苯甲酸反应生成4-羟基-3-聚戊烯基苯甲酸。基因CCM_06092的注释结果为长链脂酰辅酶A合成酶，该酶的功能是催化长链脂肪酸和辅酶A反应生成长链脂酰辅酶A。基因CCM_05423的注释结果为甲基丙二酰辅酶A合成酶，该酶的功能是催化生物素和阿朴-甲基丙二酰辅酶A反应生成甲基丙二酰辅酶A。基因CCM_08567的注释结果是雀稗灵合成酶，但该酶的催化反应式未被注释，因此其催化功能有待被鉴定。基因CCM_09041的注释结果是丙酰辅酶A合成酶，该酶的功能是催化丙酸和辅酶A反应生成丙酰辅酶A。下一步，将鉴定上述5 个基因在蛹虫草中的功能。

2.6 real-time PCR验证转录组数据

由表3可知，基因CCM_05119和CCM_07824的FDR＜0.001，log2foldChange（CM10_WL/CM10_RL）≥1，即fold change（CM10_WL/CM10_RL）≥2，表明基因CCM_05119和CCM_07824显著上调表达；而基因CCM_06728和CCM_09155的FDR＞0.001，log2foldChange（CM10_WL/CM10_RL）＜1，即fold change（CM10_WL/CM10_RL）＜2，表明基因CCM_06728和CCM_09155的表达水平在样品CM10_RL和CM10_WL中无显著差异。由表3中的real-time PCR实验结果可知，相对于样品CM10_RL，基因CCM_05119和CCM_07824在样品CM10_WL中的表达水平均大于2 倍，显著上调表达，而基因CCM_06728和CCM_09155在2 种样品中的表达水平差异不显著。real-time PCR实验结果与转录组实验结果一致，说明转录组数据可靠。

表 3 real-time PCR检测4 个基因的相对表达水平Table 3 Relative expression levels of four genes measured by real-time PCR

3 讨论与结论

不同培养基培养的蛹虫草菌丝样品（CM10_RL和CM10_WL）具有不同的颜色表型和不同的类胡萝卜素含量，样品CM10_WL较样品CM10_RL具有更深的色泽和更高的类胡萝卜素含量。在此基础上，采用高通量测序技术对蛹虫草菌丝样品CM10_RL和CM10_WL进行转录组测序和生物信息学分析，结果发现样品CM10_RL和CM10_WL之间有936 个DEG（包括318 个上调表达基因和618 个下调表达基因）；通过GO功能富集分析，发现在生物过程、细胞组分、分子功能类中富集最多DEG的条目分别是代谢过程、细胞膜、催化活性；通过KEGG通路富集分析，发现富集DEG最多的通路为代谢途径；把本研究所有蛹虫草转录本与KEGG通路数据库中的类胡萝卜素生物合成途径进行比对，结果发现基因CCM_06728和CCM_09155参与蛹虫草类胡萝卜素的生物合成；基于类胡萝卜素属于萜类化合物，从DEG中筛选到5 个参与萜类生物合成的基因（CCM_03738、CCM_06092、CCM_05423、CCM_08567、CCM_09041）。此外，real-time PCR实验结果与转录组结果一致，转录组数据可靠。

培养基对蛹虫草的生长和代谢有显著影响[19-20]。杨莹等[38]研究了不同氮源培养基对蛹虫草生长及类胡萝卜素含量的影响，发现不同氮源培养基显著影响蛹虫草的菌落形态、菌丝生长速度、分生孢子产量和类胡萝卜素含量。Guo Mingmin等[19]以小麦和蛹作为培养基栽培蛹虫草，发现培养基显著影响蛹虫草中的虫草素、腺苷、多糖、甘露醇、类胡萝卜素的含量及超氧化物歧化酶的活性。宋仙妹等[20]发现蛹虫草生长所需的最适C/N为20∶1，简锦辉[18]在此基础上进行了培养基的优化，获得了适用于蛹虫草产类胡萝卜素的培养基，因此本实验采用其所优化的培养基。本实验的研究结果表明，蛹虫草在培养基W上产类胡萝卜素的量显著高于在培养基R上产类胡萝卜素的量，证明了培养基显著影响蛹虫草类胡萝卜素的生物合成。

蛹虫草中新型类胡萝卜素的发现[14]，使其成为研究的热点。蛹虫草类胡萝卜素的含量被推荐作为衡量其商品价值的一个重要质量标准[33]，因此提高蛹虫草类胡萝卜素的含量将有利于提高其商品价值。通过优化培养基和培养条件可提高蛹虫草类胡萝卜素的含量[28-29]，但提高量有限，不能从根本上解决蛹虫草类胡萝卜素含量低的问题，因此挖掘蛹虫草类胡萝卜素的生物合成相关基因、阐明蛹虫草类胡萝卜素的生物合成途径是亟待解决的问题。本实验把蛹虫草样品CM10_RL和CM10_WL的所有转录本与KEGG通路数据库中的类胡萝卜素生物合成途径（map00906）进行比对，仅比对出2 个基因（CCM_06728和CCM_09155）参与蛹虫草类胡萝卜素的生物合成，这表明蛹虫草类胡萝卜素的生物合成途径比较独特，与KEGG数据库中已知的类胡萝卜素生物合成途径（map00906）有较大差异，这可能是因为蛹虫草类胡萝卜素是新型的类胡萝卜素，其生物合成途径可能是新的生物合成途径。由于蛹虫草类胡萝卜素是新型的类胡萝卜素，又基于类胡萝卜素属于萜类化合物，所以从DEG中挖掘到5 个参与萜类生物合成的基因（CCM_03738、CCM_06092、CCM_05423、CCM_08567、CCM_09041），这5 个基因可能参与蛹虫草类胡萝卜素的生物合成。因此，为了揭示蛹虫草类胡萝卜素的生物合成途径，将通过基因敲除、基因超表达、代谢组检测等技术鉴定上述7 个基因（CCM_06728、CCM_09155、CCM_03738、CCM_06092、CCM_05423、CCM_08567、CCM_09041）的功能。

本实验研究不同培养基培养的蛹虫草菌丝（CM10_RL和CM10_WL）的转录组差异，证实了培养基能够显著影响蛹虫草类胡萝卜素的生物合成，获得了蛹虫草菌丝样品CM10_RL和CM10_WL之间的DEG，挖掘到7 个（CCM_06728、CCM_09155、CCM_03738、CCM_06092、CCM_05423、CCM_08567、CCM_09041）可能参与蛹虫草类胡萝卜素生物合成的基因，这些研究结果为揭示蛹虫草类胡萝卜素的生物合成途径提供一定的参考依据。