论文部分内容阅读
摘 要： 研究辣椒花药中sRNA分布，筛选育性相关miRNA，通过miRNA及其靶基因的表达分析探讨miRNA对雄性育性的调控。利用高通量测序技术对辣椒细胞核雄性不育两用系处于小孢子单核靠边期的花药进行sRNA测序，同时分析两材料间miRNA的表达差异。通过qRT-PCR技术验证差异分析结果，并对miRNA及其靶基因在小孢子发育不同时期的表达模式进行分析。在构建的不育株和可育株两个文库中共发掘出857个可信度高的miRNA；筛选得到42个表达差异显著的miRNA；通过靶基因预测与注释，预测出的差异表达miRNA的靶基因中与繁殖有关的基因有152 个，与生殖过程有关的基因有151个。通过qRT-PCR验证了选出的9条miRNA的存在，其表达差异与测序分析结果基本一致；单核靠边期多数miRNA负调控其靶基因，不同发育时期其调控模式会发生变化。本研究为揭示辣椒miRNA与雄性育性的关系提供了重要信息。
关键词： 辣椒； 细胞核雄性不育； miRNA； 高通量测序； qRT-PCR
Identification and expression analysis of fertility-related miRNAs for the genic male sterile line in pepper （Capsicum annuum L.）
LIU Siyue， LIU Chen， GUO Jinju， DA Fengjiao， SHEN Huolin
（Beijing Key Laboratory of Growth and Developmental Regulation for Protected Vegetable Crops·College of Agriculture and Biotechnology，China Agricultural University，Beijing 100193，China）
Abstract： The objective of this study is to investigate the distribution of sRNAs from anthers of pepper，and to screen miRNAs related to genic male sterility. By expression analysis of miRNAs and their target genes，we discussed the regulation mechanisam of male sterility by miRNAs. Using high-throughput sequencing technology，sRNA populations from anthers in late-uninucleate stage for GMS line in pepper were sequenced. qRT-PCR were performed to validate the results of sRNA sequencing and to analyze the expression profile in different developmental stages of microspore. Of the millions of high-quality sRNAs from 2 libraries，857 miRNAs of high credibility were identified，42 significant differently expressed miRNAs were screened. Prediction and annotation of target genes indicated that there are 152 targets involved in reproduction and 151 involved in reproductive process. qRT-PCR results confirmed the existence of the 9 selected miRNAs，and the expression files of most of the miRNAs is consistent with sequencing results. In late-uninucleate stage，most of miRNAs negatively regulate their targets，but the regulation file differs in different development phase. This study provides important information for uncovering the relation between miRNAs and male fertility in pepper.
Key words： Pepper； Genic male sterility； miRNA； High-throughput sequencing； qRT-PCR
辣椒是世界各地普遍栽培的蔬菜和调味品。目前，人工授粉是生产辣椒杂交种的主要方式，但这种方式成本高并且种子纯度难以保证，而利用雄性不育系则可以克服这些缺点。本研究通过初步分析细胞核雄性不育辣椒两用系miRNA和其靶基因的表达模式，为辣椒miRNA与育性关系的研究提供参考。辣椒细胞核雄性不育类型（genetic male sterile，GMS）因败育彻底，恢复源广，在杂交制种过程中发挥重要作用。迄今为止，世界各国学者已发现的辣椒细胞核雄性不育基因近20个，其中由国内学者发现的不育基因有2个[1]。克隆分析得到的育性相关基因有CaMF2、CaPME1等[2-3]。 sRNA（small RNA）是调节植物生长发育、抵御逆境的多种小分子RNA的总称，miRNA（microRNA）是sRNA的一种，是广泛存在于生物体内的一类长度为20～24 nt的内源非编码小RNA分子，它们能通过沉默靶基因或者抑制翻译来调控基因的表达。Lee等[4]1993年在秀丽隐杆线虫（Caenorhabditis elegans）体内发现了第一个miRNA lin-4，它参与线虫发育进程的调控。最早发现的植物miRNA是2002年由Reinhart等[5]从拟南芥（Arabidopsis thaliana）sRNA文库中获得。植物中许多物种都发现了大量的miRNA，并且大多数植物miRNA在物种间是高度保守的。据估计，编码miRNA的基因可能占基因总数的1%～5%[6]。miRNA主要通过与靶基因的特异性碱基互补配对，从而降解靶基因或抑制其翻译，在转录后水平对动植物的基因表达进行调控[7-9]。
目前有关miRNA的研究主要包括调控植物的生长发育、逆境生理以及响应激素信号等方面[10-15]，其中对miRNA与花粉发育关系的研究主要集中在模式植物拟南芥中[16-17]，而对于miRNA与育性的关系研究目前还较少。Schwab等[18]研究发现，过量表达miR156会导致植株在短日照条件下开花延迟和育性降低，其作用机制可能是通过靶向作用于转录因子SPL；Mallory等[19]在含有5个沉默突变位点的ARF17 mRNA转基因植株中观察到，植株表现出多种异常表型，包括畸形叶片产生、开花期提前、育性降低等。这是由于miR160作用的靶标——ARF17 mRNA含有5个沉默突变位点，导致miR160不能正常介导ARF17 mRNA的切割。虽然这些研究初步证明了miRNA与植物育性的相关性，但是要揭示具体的调控机制，目前的研究成果还相当匮乏，辣椒中更是未见报道。
研究表明，辣椒的小孢子单核靠边期是花粉败育的重要时期[20-22]，并且是许多已报道育性相关基因表达差异大的时期[23-24]。因此，本试验通过对辣椒GMS小孢子单核靠边期建立miRNA差异表达文库，研究其与育性的关系。通过建立辣椒小孢子单核靠边期差异表达的miRNA文库，筛选差异表达miRNA并预测其靶基因，利用qRT-PCR技术验证差异分析结果并鉴定miRNA与靶基因相关关系，为辣椒miRNA育性调控功能研究提供信息。
1 材料与方法
1.1 材料
本试验所用GMS两用系辣椒材料‘629AB’由中国农业大学辣椒课题组提供。2014年春在中国农业大学上庄实验站大棚内种植，常规管理。在开花期，分别收集不育株和可育株的花蕾，根据小孢子发育时期与花器官形态的相关性，将收集到的花蕾分为四分体时期、单核早中期、单核靠边期和双核期[25]，剥取花药并立即用液氮冷冻，-80 ℃保存待用。
1.2 small RNA文库的测序及生物信息学分析
1.2.1 总RNA的提取，sRNA文库的建立并测序
Trizol（Invitrogen公司，美国）法分别提取不育株和可育株单核靠边期花药RNA。用NanoDrop紫外分光光度计检测浓度和纯度（RNA浓度≥500 mg·L-1，OD260/280为1.8～2.2，260 nm处有正常峰值），1%琼脂糖凝胶电泳检测总RNA完整性（RNA 28S/18S≥1.5）。
样品检测合格后，利用sRNA的3’及5’端特殊结构（5’端有完整的磷酸基团，3’端有羟基）的特征，以总RNA为起始样品，在sRNA的3’端和5’端添加接头，反转录合成cDNA。经过PCR扩增，PAGE凝胶电泳分离目标DNA片段，切胶回收得到cDNA文库。经检测合格后，利用HiSeqTM 2500基因分析仪测序。
1.2.2 miRNA鉴定 原始reads经过过滤低质量并除去无接头后长度小于18 nt和大于30 nt的reads，最终得到clean reads。将clean reads分别与GenBank库（http：//www.ncbi.nih.gov./genebank/）、Rfam库（http：//rfam.sanger.ac.uk/）比对，去除其他ncRNA（非编码RNA，包括rRNA、tRNA、snRNA、snoRNA）序列，剩下的与参考基因组比对获得比对到基因组的序列。对于比对到参考基因组的序列，利用miRDeep2软件进行miRNA前体预测，软件对star序列、序列数目、MFE（最小折叠自由能））、cons.seed等分别进行打分，最终分数大于0即为前体序列。
预测到了前体序列则可以确定前体序列中的成熟序列、环状结构、star序列的位置，将比对到成熟序列位置上的序列提取出来即为成熟的miRNA。
将预测到的前体序列和成熟序列与miRBase数据库中的所有植物做比对，来定义保守miRNA和非保守miRNA。
1.2.3 microRNA差异表达分析 整合2个样品中比对到各miRNA上的reads数，利用IDEG6，选取FDR<0.01且Fold Change>=2作为标准，进行差异表达分析。FDR（False Discovery Rate）表示错误发现率，它是通过对p-value进行校正得到的。转录组中的差异miRNA分析是对大量的miRNA进行独立的统计假设检验，它会导致总体假阳性偏高的问题，因此在差异软件进行差异分析过程中，采用Benjamini Hochberg校正方法[26]对原有假设检验得到的p-value进行校正，最终采用FDR作为差异基因筛选取得的关键指标。
1.2.4 靶基因预测及注释 利用预测的miRNA序列信息和辣椒基因组中的序列信息文件，通过TargetFinder软件进行靶基因预测。靶基因预测的罚分系统如下：1）错配、单碱基缺失和凸起分别罚1分；2）出现G：U配对罚0.5分；3）microRNA 5’端的第2～第13位碱基出现G：U以外的错配罚1分。若配对出现以下情况则认定其不是该miRNA的靶基因：1）多于1个单碱基缺失或凸起；2）错配、G：U配对、缺失或凸起的总数大于7；3）多于4个错配数或多于4个G：U配对。最终罚分小于或等于3则预测其为该miRNA的靶基因。 使用BLAST软件将预测得到的靶基因序列与GO数据库[27]比对，获得靶基因序列的注释信息。
1.3 qRT-PCR分析
1.3.1 miRNA的qRT-PCR检测 Trizol法分别提取小孢子四分体时期、单核早中期、单核靠边期、双核期4个时期的总RNA。检测浓度和完整性后，置于-80 ℃ 贮存备用。利用miRNA cDNA Kit（康为世纪，中国）进行miRNA对应的cDNA第一链合成：首先在miRNA3’末端加多聚A尾的方法来使miRNA具有Poly（A）尾，之后再使用Oligo（dT）-Universal tag通用逆转录引物进行逆转录反应，最终合成miRNA对应的第一链cDNA。最后用荧光定量 PCR 检测 miRNA 的表达情况。实时荧光定量利用miRNA Real-Time PCR Assay Kit（康为世纪，中国），按照说明冰上配制下列反应体系：2×miRNA qPCR premix（With SYBR and ROX） 10 μL，Forward primer（10 μmol·L-1） 0.4 μL，Reverse primer（10 μmol·L-1） 0.4 μL，miRNA第一链cDNA模板2.5 μL，RNase-Free Water补足20 μL反应体系。利用ABI 7500荧光定量PCR仪执行以下反应条件：预变性95 ℃ 10 min，变性95 ℃ 15 s，退火延伸60 ℃ 1 min，45个循环。反应结束后确认Real Time PCR的扩增曲线和溶解曲线，记录Ct值并计算2-ΔΔCt。Reverse primer使用试剂盒自带引物，Forward primer序列列于表1。
1.3.2 靶基因qRT-PCR检测 小孢子发育4个时期总RNA的提取方法和检测方法同上。利用PrimeScript? II 1st strand cDNA Synthesis Kit（Takara公司，日本）按照说明书进行cDNA第一链的合成，然后利用Promega GoTaq qPCR Master Mix（Promega公司，美国）冰上配制反应体系：2×GoTaq? qPCR Master Mix 10 μL，Forward primer（10 μmol·L-1） 0.4 μL，Reverse primer（10 μmol·L-1） 0.4 μL，cDNA模板2 μL，RNase-Free Water 7.2 μL共20 μL的反应体系。利用ABI 7500荧光定量PCR仪执行以下反应条件：预变性95 ℃ 10 min，变性95 ℃ 15 s，退火延伸60 ℃ 1 min，40个循环。反应结束后确认Real Time PCR的扩增曲线和溶解曲线，记录Ct值并计算2-ΔΔCt。引物序列列于表2。
2 结果与分析
2.1 测序结果统计
分别对小孢子发育单核靠边期的不育株和可育株的花药进行Illumina测序。2个材料通过测序，过滤低质量reads、接头及过长过短序列，最终得到clean reads分别为19 852 077、16 503 277条。通常，sRNA库的测序结果在成分上比较复杂，除了包括大量miRNA和siRNA，还有来自其他编码和非编码序列的降解片段，因此将这些clean reads分别与GenBank库、Rfam库比对，提取注释信息 （表3）。
不育株中匹配到参考基因组（http：//peppersequence.genomics.cn/page/species/index.jsp）的序列有9 489 233条，占总序列数的47.80%，可育株中可以匹配到的序列有7 464 111条，占总序列数的45.23%（表3）。
表3 clean reads分类注释结果
2个文库的sRNA结果显示，不育株和可育株中的sRNA都是以24 nt最多，分别占总数的51.28%、51.63%。miRDeep2预测得到不育株和可育株的miRNA长度主要在21、22、24 nt。分析miRNA位点对碱基的偏好发现，miRNA成熟体序列首位碱基对碱基U具有很强的偏向性。2个库共发掘出857个可信度高的新的miRNA，与miRBase数据库比对后得到保守的miRNA有262个，非保守的有595个。
2.2 差异表达的miRNA
2个库共筛选出42个差异显著表达的miRNA，其中可育株相对于不育株表达上调的有16个，下调的有26个。上调的16个miRNA中，有3个保守的miRNA和13个非保守miRNA。下调的miRNA中，有11个保守miRNA，15个非保守miRNA。
2.3 靶基因预测及注释
找到miRNA所调控的mRNA靶标对理解miRNA的生物功能至关重要。在857个预测得到的miRNA中，有733个可以预测到靶基因，共预测到12 374个靶基因，其中223个保守的miRNA共预测到4 850个靶基因，510个非保守的miRNA共预测到7 524个靶基因（表4）。
对预测得到的miRNA的靶基因进行GO分类注释，共注释到6 879个靶基因，其中，差异表达的miRNA靶基因共注释得到2 176个，与总的靶基因GO分类相比，差异miRNA的靶基因在细胞过程和免疫系统过程方面有较多富集。
所有靶基因中与繁殖有关的注释到527个，与生殖过程有关的基因有523个。差异miRNA的靶基因在2种分类中分别为152和151个。经过进一步的分级注释，筛选得到6个注释与生殖有关的靶基因（表5）。其他还包括一系列调控分生组织营养生长向生殖生长转变的基因。
表5 GO注释与生殖有关的部分
2.4 qRT-PCR
为了验证新预测的miRNA在辣椒中的真实性、在小孢子发育的不同时期表达量变化情况及其与靶基因的关系，选取9条miRNA与其相应的9条靶基因，根据序列设计引物，用qRT-PCR对辣椒小孢子发育的4个时期的miRNA和靶基因进行检测（表6）。为了验证测序结果的准确性，选取了Capana08g000340、Capana08g000534、Capana00g000886为测序数据显示有显著表达差异的miRNA对应的靶基因。Capana04g001683、Capana06g002954、Capana12g002894为表达差异大并在GO分类中与生殖相关的基因，CaSEP1、CaPROF、CaNAC为本实验室筛选得到的育性相关候选基因。 [注] CaSEP1、CaPROF、CaNAC在基因组中的编号分别为Capana02g003249、Capana06g001250、Capana03g003315。
荧光定量融解曲线显示各miRNA及靶基因扩增产物为单峰，表明9条miRNA在辣椒中真实存在。unconservative_Chr00_6109657、conservative _Chr00_6253770、unconservative_Chr05_2660285在两样本中表达差异大，在单核靠边期，3条miRNA在不育株中的表达量分别是可育株的2.74倍、3.83倍、2.74倍，与测序结果一致。
9个miRNA在辣椒小孢子发育的不同时期中均有不同程度的表达，CaSEP1、CaPROF、CaNAC这3个实验室筛选到的基因对应的miRNA呈现一致的表达趋势，即前期表达量较低，单核靠边期表达量显著增加，达到最高，双核期表达量又迅速下降，不育株和可育株中趋势变化一致。GO分类中与育性相关的靶基因对应的miRNA表达趋势各不相同，但从总量来看，单核靠边期都是这3条miRNA大量表达的时期，同时又是unconservative_Chr06_2966220、unconservative_Chr00_6146461这2条miRNA的表达量显著差异的时期，在不育株中的表达量分别是可育株的13.21倍和2.31倍（图1）。
分析这些miRNA在4个时期表达趋势变化，发现最多的一类就是四分体时期、单核早中期表达量相对低，单核靠边期表达量显著增高，双核期又明显下降。
靶基因方面，在小孢子发育的不同时期，选取的9条靶基因呈现各异的表达趋势，并且在不同阶段，与其相应的miRNA存在的相关关系各不相同。从本试验选出做qRT-PCR的miRNA相对表达丰度最高且表达差异比较大的单核靠边期来看，有 7个miRNA与其相应的靶基因呈现负相关关系。但这种负相关关系从各个时期来看并不严格，会因为时期的不同和表达量变化而呈现其他相关关系。例如，unconservative_Chr00_6146461与其靶基因Capana12g002894在单核靠边期、双核期表现负相关，但在四分体时期和单核早中期却呈现正相关（图1）。
3 讨 论
近年来对于植物miRNA的研究越来越深入，但主要集中在拟南芥[28-29]、水稻[30]等模式植物上，对辣椒的miRNA研究少有报道，关于辣椒miRNA的测序研究主要关注miRNA在不同组织中的表达[31]。本研究第一次从细胞核雄性不育两用系的角度利用高通量测序技术对辣椒miRNA进行了分析，通过qRT-PCR技术探讨了miRNA及其靶基因在不同发育时期的表达模式以及二者的关系，为阐明miRNA对辣椒雄蕊育性的影响提供了依据。
本试验预测得到的sRNA的长度分布以24 nt最多，这与许多已报道的物种拟南芥[32]、小麦[33]和棉花[34]中24 nt的sRNA分布最多情况一致。已报道的辣椒miRNA长度多在21 nt[35]，本研究的结果表明，miRNA的序列长度主要分布在21、22、24 nt，这是对前人研究结果的验证和补充。miRNA成熟序列首位碱基对U具有很强的偏向性，这一结果符合了前体miRNA（pre-miRNA）在Exportin-5帮助下转运到细胞核外，由细胞质Dicer酶对特异位点进行酶切的特性[36]。2个材料共预测出857个miRNA，其中差异表达的有42个，差异表达miRNA对应的靶基因多富集在辣椒的新陈代谢过程中，说明在辣椒小孢子单核靠边期，miRNA可能是通过调节新陈代谢来调控育性。
通过对挑选出的9个miRNA与其靶基因进行qRT-PCR分析发现，多数miRNA在单核靠边期表达差异最大，并负调控其靶基因。之前的研究表明unconservative_Chr11_4964464的靶基因与雄蕊发育有关，unconservative_Chr12_5488192的靶基因与花粉的萌发和花粉管的生长有关[23]，这2个miRNA在单核靠边期均大量表达并负调控其靶基因。这一时期多数miRNA与其靶基因的负相关关系一定程度上反映了miRNA与其靶基因的负调控模式。miRNA负调控其靶基因此前已多有研究[7-9]，其中育性调控方面也是如此。例如植物过量表达miR156后，其靶标的一类转录因子基因SPL表达量减少，出现短日照条件下开花延迟、育性降低的现象[18]。miR160合成减少时，其靶标ARF17会积累，植株出现育性降低的现象[19]。另外，本试验的结果显示，有少数miRNA与其靶基因并不存在严格的负相关关系，同一miRNA在不同时期可能呈现不同的相关关系，这些结果可能与之前报道的miRNA与其靶基因呈现“一对多”（一个miRNA对应多个靶基因）、“多对一”（多个miRNA对应一个靶基因）的对应关系有关[19]。因此我们进一步推测不同时期miRNA侧重调控的靶基因不同，同一时期发挥主要作用的miRNA不同。
本研究证明了前期发现的育性相关基因和本试验中预测的部分育性相关基因与其对应miRNA间的负相关关系，说明miRNA可能是通过负调控其靶基因进而影响辣椒育性。但是，由于目前辣椒育性相关基因的报道较少，育性决定基因也没有被克隆，辣椒细胞核雄性不育的分子机理尚不明确，miRNA与育性的确切关系也还需要进一步的验证。
4 结 论
本研究分析了辣椒细胞核雄性不育两用系花药中sRNA的分布情况，筛选得到了一系列在两材料中表达差异显著的miRNA，这些miRNA的靶基因多富集在新陈代谢过程中，并且与生殖直接相关的靶基因功能涉及生殖结构发育、多细胞生物生殖相关、生殖中涉及到的细胞过程、有性生殖等方面，这些miRNA与其靶基因为辣椒育性研究提供了重要信息。本研究初步分析了细胞核雄性不育辣椒中的miRNA和其相应靶基因的表达模式，为辣椒育性相关miRNA的功能研究提供了参考。 参考文献
[1] Wang D，Bosland P W. The genes of Capsicum[J]. HortScience，2006，41（5）： 1169-1187.
[2] Chen C，Chen G，Hao X，et al. CaMF2，an anther-specific lipid transfer protein （LTP） gene，affects pollen development in Capsicum annuum L.[J]. Plant Science，2011，181（4）： 439-448.
[3] Chen C M，Liu S Q，Hao X F，et al. Characterization of a pectin methylesterase gene homolog，CaPME1，expressed in anther tissues of Capsicum annuum L.[J]. Plant Molecular Biology Reporter，2011，30（2）： 402-413.
[4] Lee R C，Feinbaum R L，Ambros V. The C.elegans heterochronic gene lin-4 small RNAs with antisense complementarity to Lin-14[J]. Cell，1993，75（5）： 843-854.
[5] Reinhart B J，Weinstein E G，Rhoades M W，et al. MicroRNAs in plants[J]. Genes
关键词： 辣椒； 细胞核雄性不育； miRNA； 高通量测序； qRT-PCR
Identification and expression analysis of fertility-related miRNAs for the genic male sterile line in pepper （Capsicum annuum L.）
LIU Siyue， LIU Chen， GUO Jinju， DA Fengjiao， SHEN Huolin
（Beijing Key Laboratory of Growth and Developmental Regulation for Protected Vegetable Crops·College of Agriculture and Biotechnology，China Agricultural University，Beijing 100193，China）
Abstract： The objective of this study is to investigate the distribution of sRNAs from anthers of pepper，and to screen miRNAs related to genic male sterility. By expression analysis of miRNAs and their target genes，we discussed the regulation mechanisam of male sterility by miRNAs. Using high-throughput sequencing technology，sRNA populations from anthers in late-uninucleate stage for GMS line in pepper were sequenced. qRT-PCR were performed to validate the results of sRNA sequencing and to analyze the expression profile in different developmental stages of microspore. Of the millions of high-quality sRNAs from 2 libraries，857 miRNAs of high credibility were identified，42 significant differently expressed miRNAs were screened. Prediction and annotation of target genes indicated that there are 152 targets involved in reproduction and 151 involved in reproductive process. qRT-PCR results confirmed the existence of the 9 selected miRNAs，and the expression files of most of the miRNAs is consistent with sequencing results. In late-uninucleate stage，most of miRNAs negatively regulate their targets，but the regulation file differs in different development phase. This study provides important information for uncovering the relation between miRNAs and male fertility in pepper.
Key words： Pepper； Genic male sterility； miRNA； High-throughput sequencing； qRT-PCR
辣椒是世界各地普遍栽培的蔬菜和调味品。目前，人工授粉是生产辣椒杂交种的主要方式，但这种方式成本高并且种子纯度难以保证，而利用雄性不育系则可以克服这些缺点。本研究通过初步分析细胞核雄性不育辣椒两用系miRNA和其靶基因的表达模式，为辣椒miRNA与育性关系的研究提供参考。辣椒细胞核雄性不育类型（genetic male sterile，GMS）因败育彻底，恢复源广，在杂交制种过程中发挥重要作用。迄今为止，世界各国学者已发现的辣椒细胞核雄性不育基因近20个，其中由国内学者发现的不育基因有2个[1]。克隆分析得到的育性相关基因有CaMF2、CaPME1等[2-3]。 sRNA（small RNA）是调节植物生长发育、抵御逆境的多种小分子RNA的总称，miRNA（microRNA）是sRNA的一种，是广泛存在于生物体内的一类长度为20～24 nt的内源非编码小RNA分子，它们能通过沉默靶基因或者抑制翻译来调控基因的表达。Lee等[4]1993年在秀丽隐杆线虫（Caenorhabditis elegans）体内发现了第一个miRNA lin-4，它参与线虫发育进程的调控。最早发现的植物miRNA是2002年由Reinhart等[5]从拟南芥（Arabidopsis thaliana）sRNA文库中获得。植物中许多物种都发现了大量的miRNA，并且大多数植物miRNA在物种间是高度保守的。据估计，编码miRNA的基因可能占基因总数的1%～5%[6]。miRNA主要通过与靶基因的特异性碱基互补配对，从而降解靶基因或抑制其翻译，在转录后水平对动植物的基因表达进行调控[7-9]。
目前有关miRNA的研究主要包括调控植物的生长发育、逆境生理以及响应激素信号等方面[10-15]，其中对miRNA与花粉发育关系的研究主要集中在模式植物拟南芥中[16-17]，而对于miRNA与育性的关系研究目前还较少。Schwab等[18]研究发现，过量表达miR156会导致植株在短日照条件下开花延迟和育性降低，其作用机制可能是通过靶向作用于转录因子SPL；Mallory等[19]在含有5个沉默突变位点的ARF17 mRNA转基因植株中观察到，植株表现出多种异常表型，包括畸形叶片产生、开花期提前、育性降低等。这是由于miR160作用的靶标——ARF17 mRNA含有5个沉默突变位点，导致miR160不能正常介导ARF17 mRNA的切割。虽然这些研究初步证明了miRNA与植物育性的相关性，但是要揭示具体的调控机制，目前的研究成果还相当匮乏，辣椒中更是未见报道。
研究表明，辣椒的小孢子单核靠边期是花粉败育的重要时期[20-22]，并且是许多已报道育性相关基因表达差异大的时期[23-24]。因此，本试验通过对辣椒GMS小孢子单核靠边期建立miRNA差异表达文库，研究其与育性的关系。通过建立辣椒小孢子单核靠边期差异表达的miRNA文库，筛选差异表达miRNA并预测其靶基因，利用qRT-PCR技术验证差异分析结果并鉴定miRNA与靶基因相关关系，为辣椒miRNA育性调控功能研究提供信息。
1 材料与方法
1.1 材料
本试验所用GMS两用系辣椒材料‘629AB’由中国农业大学辣椒课题组提供。2014年春在中国农业大学上庄实验站大棚内种植，常规管理。在开花期，分别收集不育株和可育株的花蕾，根据小孢子发育时期与花器官形态的相关性，将收集到的花蕾分为四分体时期、单核早中期、单核靠边期和双核期[25]，剥取花药并立即用液氮冷冻，-80 ℃保存待用。
1.2 small RNA文库的测序及生物信息学分析
1.2.1 总RNA的提取，sRNA文库的建立并测序
Trizol（Invitrogen公司，美国）法分别提取不育株和可育株单核靠边期花药RNA。用NanoDrop紫外分光光度计检测浓度和纯度（RNA浓度≥500 mg·L-1，OD260/280为1.8～2.2，260 nm处有正常峰值），1%琼脂糖凝胶电泳检测总RNA完整性（RNA 28S/18S≥1.5）。
样品检测合格后，利用sRNA的3’及5’端特殊结构（5’端有完整的磷酸基团，3’端有羟基）的特征，以总RNA为起始样品，在sRNA的3’端和5’端添加接头，反转录合成cDNA。经过PCR扩增，PAGE凝胶电泳分离目标DNA片段，切胶回收得到cDNA文库。经检测合格后，利用HiSeqTM 2500基因分析仪测序。
1.2.2 miRNA鉴定 原始reads经过过滤低质量并除去无接头后长度小于18 nt和大于30 nt的reads，最终得到clean reads。将clean reads分别与GenBank库（http：//www.ncbi.nih.gov./genebank/）、Rfam库（http：//rfam.sanger.ac.uk/）比对，去除其他ncRNA（非编码RNA，包括rRNA、tRNA、snRNA、snoRNA）序列，剩下的与参考基因组比对获得比对到基因组的序列。对于比对到参考基因组的序列，利用miRDeep2软件进行miRNA前体预测，软件对star序列、序列数目、MFE（最小折叠自由能））、cons.seed等分别进行打分，最终分数大于0即为前体序列。
预测到了前体序列则可以确定前体序列中的成熟序列、环状结构、star序列的位置，将比对到成熟序列位置上的序列提取出来即为成熟的miRNA。
将预测到的前体序列和成熟序列与miRBase数据库中的所有植物做比对，来定义保守miRNA和非保守miRNA。
1.2.3 microRNA差异表达分析 整合2个样品中比对到各miRNA上的reads数，利用IDEG6，选取FDR<0.01且Fold Change>=2作为标准，进行差异表达分析。FDR（False Discovery Rate）表示错误发现率，它是通过对p-value进行校正得到的。转录组中的差异miRNA分析是对大量的miRNA进行独立的统计假设检验，它会导致总体假阳性偏高的问题，因此在差异软件进行差异分析过程中，采用Benjamini Hochberg校正方法[26]对原有假设检验得到的p-value进行校正，最终采用FDR作为差异基因筛选取得的关键指标。
1.2.4 靶基因预测及注释 利用预测的miRNA序列信息和辣椒基因组中的序列信息文件，通过TargetFinder软件进行靶基因预测。靶基因预测的罚分系统如下：1）错配、单碱基缺失和凸起分别罚1分；2）出现G：U配对罚0.5分；3）microRNA 5’端的第2～第13位碱基出现G：U以外的错配罚1分。若配对出现以下情况则认定其不是该miRNA的靶基因：1）多于1个单碱基缺失或凸起；2）错配、G：U配对、缺失或凸起的总数大于7；3）多于4个错配数或多于4个G：U配对。最终罚分小于或等于3则预测其为该miRNA的靶基因。 使用BLAST软件将预测得到的靶基因序列与GO数据库[27]比对，获得靶基因序列的注释信息。
1.3 qRT-PCR分析
1.3.1 miRNA的qRT-PCR检测 Trizol法分别提取小孢子四分体时期、单核早中期、单核靠边期、双核期4个时期的总RNA。检测浓度和完整性后，置于-80 ℃ 贮存备用。利用miRNA cDNA Kit（康为世纪，中国）进行miRNA对应的cDNA第一链合成：首先在miRNA3’末端加多聚A尾的方法来使miRNA具有Poly（A）尾，之后再使用Oligo（dT）-Universal tag通用逆转录引物进行逆转录反应，最终合成miRNA对应的第一链cDNA。最后用荧光定量 PCR 检测 miRNA 的表达情况。实时荧光定量利用miRNA Real-Time PCR Assay Kit（康为世纪，中国），按照说明冰上配制下列反应体系：2×miRNA qPCR premix（With SYBR and ROX） 10 μL，Forward primer（10 μmol·L-1） 0.4 μL，Reverse primer（10 μmol·L-1） 0.4 μL，miRNA第一链cDNA模板2.5 μL，RNase-Free Water补足20 μL反应体系。利用ABI 7500荧光定量PCR仪执行以下反应条件：预变性95 ℃ 10 min，变性95 ℃ 15 s，退火延伸60 ℃ 1 min，45个循环。反应结束后确认Real Time PCR的扩增曲线和溶解曲线，记录Ct值并计算2-ΔΔCt。Reverse primer使用试剂盒自带引物，Forward primer序列列于表1。
1.3.2 靶基因qRT-PCR检测 小孢子发育4个时期总RNA的提取方法和检测方法同上。利用PrimeScript? II 1st strand cDNA Synthesis Kit（Takara公司，日本）按照说明书进行cDNA第一链的合成，然后利用Promega GoTaq qPCR Master Mix（Promega公司，美国）冰上配制反应体系：2×GoTaq? qPCR Master Mix 10 μL，Forward primer（10 μmol·L-1） 0.4 μL，Reverse primer（10 μmol·L-1） 0.4 μL，cDNA模板2 μL，RNase-Free Water 7.2 μL共20 μL的反应体系。利用ABI 7500荧光定量PCR仪执行以下反应条件：预变性95 ℃ 10 min，变性95 ℃ 15 s，退火延伸60 ℃ 1 min，40个循环。反应结束后确认Real Time PCR的扩增曲线和溶解曲线，记录Ct值并计算2-ΔΔCt。引物序列列于表2。
2 结果与分析
2.1 测序结果统计
分别对小孢子发育单核靠边期的不育株和可育株的花药进行Illumina测序。2个材料通过测序，过滤低质量reads、接头及过长过短序列，最终得到clean reads分别为19 852 077、16 503 277条。通常，sRNA库的测序结果在成分上比较复杂，除了包括大量miRNA和siRNA，还有来自其他编码和非编码序列的降解片段，因此将这些clean reads分别与GenBank库、Rfam库比对，提取注释信息 （表3）。
不育株中匹配到参考基因组（http：//peppersequence.genomics.cn/page/species/index.jsp）的序列有9 489 233条，占总序列数的47.80%，可育株中可以匹配到的序列有7 464 111条，占总序列数的45.23%（表3）。
表3 clean reads分类注释结果
2个文库的sRNA结果显示，不育株和可育株中的sRNA都是以24 nt最多，分别占总数的51.28%、51.63%。miRDeep2预测得到不育株和可育株的miRNA长度主要在21、22、24 nt。分析miRNA位点对碱基的偏好发现，miRNA成熟体序列首位碱基对碱基U具有很强的偏向性。2个库共发掘出857个可信度高的新的miRNA，与miRBase数据库比对后得到保守的miRNA有262个，非保守的有595个。
2.2 差异表达的miRNA
2个库共筛选出42个差异显著表达的miRNA，其中可育株相对于不育株表达上调的有16个，下调的有26个。上调的16个miRNA中，有3个保守的miRNA和13个非保守miRNA。下调的miRNA中，有11个保守miRNA，15个非保守miRNA。
2.3 靶基因预测及注释
找到miRNA所调控的mRNA靶标对理解miRNA的生物功能至关重要。在857个预测得到的miRNA中，有733个可以预测到靶基因，共预测到12 374个靶基因，其中223个保守的miRNA共预测到4 850个靶基因，510个非保守的miRNA共预测到7 524个靶基因（表4）。
对预测得到的miRNA的靶基因进行GO分类注释，共注释到6 879个靶基因，其中，差异表达的miRNA靶基因共注释得到2 176个，与总的靶基因GO分类相比，差异miRNA的靶基因在细胞过程和免疫系统过程方面有较多富集。
所有靶基因中与繁殖有关的注释到527个，与生殖过程有关的基因有523个。差异miRNA的靶基因在2种分类中分别为152和151个。经过进一步的分级注释，筛选得到6个注释与生殖有关的靶基因（表5）。其他还包括一系列调控分生组织营养生长向生殖生长转变的基因。
表5 GO注释与生殖有关的部分
2.4 qRT-PCR
为了验证新预测的miRNA在辣椒中的真实性、在小孢子发育的不同时期表达量变化情况及其与靶基因的关系，选取9条miRNA与其相应的9条靶基因，根据序列设计引物，用qRT-PCR对辣椒小孢子发育的4个时期的miRNA和靶基因进行检测（表6）。为了验证测序结果的准确性，选取了Capana08g000340、Capana08g000534、Capana00g000886为测序数据显示有显著表达差异的miRNA对应的靶基因。Capana04g001683、Capana06g002954、Capana12g002894为表达差异大并在GO分类中与生殖相关的基因，CaSEP1、CaPROF、CaNAC为本实验室筛选得到的育性相关候选基因。 [注] CaSEP1、CaPROF、CaNAC在基因组中的编号分别为Capana02g003249、Capana06g001250、Capana03g003315。
荧光定量融解曲线显示各miRNA及靶基因扩增产物为单峰，表明9条miRNA在辣椒中真实存在。unconservative_Chr00_6109657、conservative _Chr00_6253770、unconservative_Chr05_2660285在两样本中表达差异大，在单核靠边期，3条miRNA在不育株中的表达量分别是可育株的2.74倍、3.83倍、2.74倍，与测序结果一致。
9个miRNA在辣椒小孢子发育的不同时期中均有不同程度的表达，CaSEP1、CaPROF、CaNAC这3个实验室筛选到的基因对应的miRNA呈现一致的表达趋势，即前期表达量较低，单核靠边期表达量显著增加，达到最高，双核期表达量又迅速下降，不育株和可育株中趋势变化一致。GO分类中与育性相关的靶基因对应的miRNA表达趋势各不相同，但从总量来看，单核靠边期都是这3条miRNA大量表达的时期，同时又是unconservative_Chr06_2966220、unconservative_Chr00_6146461这2条miRNA的表达量显著差异的时期，在不育株中的表达量分别是可育株的13.21倍和2.31倍（图1）。
分析这些miRNA在4个时期表达趋势变化，发现最多的一类就是四分体时期、单核早中期表达量相对低，单核靠边期表达量显著增高，双核期又明显下降。
靶基因方面，在小孢子发育的不同时期，选取的9条靶基因呈现各异的表达趋势，并且在不同阶段，与其相应的miRNA存在的相关关系各不相同。从本试验选出做qRT-PCR的miRNA相对表达丰度最高且表达差异比较大的单核靠边期来看，有 7个miRNA与其相应的靶基因呈现负相关关系。但这种负相关关系从各个时期来看并不严格，会因为时期的不同和表达量变化而呈现其他相关关系。例如，unconservative_Chr00_6146461与其靶基因Capana12g002894在单核靠边期、双核期表现负相关，但在四分体时期和单核早中期却呈现正相关（图1）。
3 讨 论
近年来对于植物miRNA的研究越来越深入，但主要集中在拟南芥[28-29]、水稻[30]等模式植物上，对辣椒的miRNA研究少有报道，关于辣椒miRNA的测序研究主要关注miRNA在不同组织中的表达[31]。本研究第一次从细胞核雄性不育两用系的角度利用高通量测序技术对辣椒miRNA进行了分析，通过qRT-PCR技术探讨了miRNA及其靶基因在不同发育时期的表达模式以及二者的关系，为阐明miRNA对辣椒雄蕊育性的影响提供了依据。
本试验预测得到的sRNA的长度分布以24 nt最多，这与许多已报道的物种拟南芥[32]、小麦[33]和棉花[34]中24 nt的sRNA分布最多情况一致。已报道的辣椒miRNA长度多在21 nt[35]，本研究的结果表明，miRNA的序列长度主要分布在21、22、24 nt，这是对前人研究结果的验证和补充。miRNA成熟序列首位碱基对U具有很强的偏向性，这一结果符合了前体miRNA（pre-miRNA）在Exportin-5帮助下转运到细胞核外，由细胞质Dicer酶对特异位点进行酶切的特性[36]。2个材料共预测出857个miRNA，其中差异表达的有42个，差异表达miRNA对应的靶基因多富集在辣椒的新陈代谢过程中，说明在辣椒小孢子单核靠边期，miRNA可能是通过调节新陈代谢来调控育性。
通过对挑选出的9个miRNA与其靶基因进行qRT-PCR分析发现，多数miRNA在单核靠边期表达差异最大，并负调控其靶基因。之前的研究表明unconservative_Chr11_4964464的靶基因与雄蕊发育有关，unconservative_Chr12_5488192的靶基因与花粉的萌发和花粉管的生长有关[23]，这2个miRNA在单核靠边期均大量表达并负调控其靶基因。这一时期多数miRNA与其靶基因的负相关关系一定程度上反映了miRNA与其靶基因的负调控模式。miRNA负调控其靶基因此前已多有研究[7-9]，其中育性调控方面也是如此。例如植物过量表达miR156后，其靶标的一类转录因子基因SPL表达量减少，出现短日照条件下开花延迟、育性降低的现象[18]。miR160合成减少时，其靶标ARF17会积累，植株出现育性降低的现象[19]。另外，本试验的结果显示，有少数miRNA与其靶基因并不存在严格的负相关关系，同一miRNA在不同时期可能呈现不同的相关关系，这些结果可能与之前报道的miRNA与其靶基因呈现“一对多”（一个miRNA对应多个靶基因）、“多对一”（多个miRNA对应一个靶基因）的对应关系有关[19]。因此我们进一步推测不同时期miRNA侧重调控的靶基因不同，同一时期发挥主要作用的miRNA不同。
本研究证明了前期发现的育性相关基因和本试验中预测的部分育性相关基因与其对应miRNA间的负相关关系，说明miRNA可能是通过负调控其靶基因进而影响辣椒育性。但是，由于目前辣椒育性相关基因的报道较少，育性决定基因也没有被克隆，辣椒细胞核雄性不育的分子机理尚不明确，miRNA与育性的确切关系也还需要进一步的验证。
4 结 论
本研究分析了辣椒细胞核雄性不育两用系花药中sRNA的分布情况，筛选得到了一系列在两材料中表达差异显著的miRNA，这些miRNA的靶基因多富集在新陈代谢过程中，并且与生殖直接相关的靶基因功能涉及生殖结构发育、多细胞生物生殖相关、生殖中涉及到的细胞过程、有性生殖等方面，这些miRNA与其靶基因为辣椒育性研究提供了重要信息。本研究初步分析了细胞核雄性不育辣椒中的miRNA和其相应靶基因的表达模式，为辣椒育性相关miRNA的功能研究提供了参考。 参考文献
[1] Wang D，Bosland P W. The genes of Capsicum[J]. HortScience，2006，41（5）： 1169-1187.
[2] Chen C，Chen G，Hao X，et al. CaMF2，an anther-specific lipid transfer protein （LTP） gene，affects pollen development in Capsicum annuum L.[J]. Plant Science，2011，181（4）： 439-448.
[3] Chen C M，Liu S Q，Hao X F，et al. Characterization of a pectin methylesterase gene homolog，CaPME1，expressed in anther tissues of Capsicum annuum L.[J]. Plant Molecular Biology Reporter，2011，30（2）： 402-413.
[4] Lee R C，Feinbaum R L，Ambros V. The C.elegans heterochronic gene lin-4 small RNAs with antisense complementarity to Lin-14[J]. Cell，1993，75（5）： 843-854.
[5] Reinhart B J，Weinstein E G，Rhoades M W，et al. MicroRNAs in plants[J]. Genes