转录组学分析强调了合子胚胎发育过程中的表观遗传和转录调控松果体松树

背景

在胚胎发生期间，植物体计划被确定，负责所有胚胎后生长的分生组织被指定。控制针叶树胚胎发生的分子机制在很大程度上仍然是未知的。它们的阐明可能为阐明被子植物和裸子植物胚胎发育的不同特征是否由差异基因调控引起提供有价值的信息。为了解决这个问题，我们首次对针叶树物种(松果体松树)，重点研究在植物胚胎发育过程中起重要作用的调控基因，即表观遗传调控因子和转录因子。

结果

基因芯片分析p .松树合子胚胎从早期发育到成熟胚胎分五个阶段进行。结果表明，转录本水平的变化主要发生在胚胎的第一个阶段和最后一个阶段，即早期到前子叶胚胎和子叶到成熟胚胎。通过胚胎发生差异表达的基因的功能分类分析强调了几种表观遗传调控机制。虽然与转座因子和重复序列机制相关的转录本的假定同源物在胚胎发生早期强烈表达，但prc2介导的基因抑制似乎在胚胎发生晚期更为相关。另一方面，在胚胎发育的各个阶段，与sRNA通路相关的功能似乎受到差异调节，miRNA功能在胚胎发育中后期普遍存在。鉴定出在连续胚胎阶段之间可能存在差异调控的转录因子基因，强烈暗示了在胚胎发生早期生长素反应和生长素载体调控的相关性。这种反应可能与建立胚胎模式有关。在发育后期，与调节生长素流动和决定正反极性的基因同源的转录本被上调，分生组织形成和功能以及器官边界建立所需的基因的同源基因也被上调。比较分析答:芥胚胎发生也突出了参与生长素介导反应的基因，以及表观遗传调控，表明两个物种之间的转录谱高度相关。

结论

本文首次报道了针叶树合子胚胎发生的时间过程转录组学分析。综上所述，我们的研究结果表明，与生长素运输和应答相关的表观遗传调控和转录调控在松木胚胎发生的早期到中期是至关重要的，胚胎发生过程中的重要事件似乎是由被子植物中主要发育调控因子的同源物协调的。

胚胎发生是大多数植物生命周期的关键时期。在模式被子植物中，生殖生物学和胚胎发育的分子方面已经得到了广泛的研究，模式被子植物在3亿多年前从裸子植物中分化出来，在精子植物中遵循一条独特的进化途径[1]。因此，在生殖和胚胎发生过程中存在着显著的差异，其中被子植物的双受精现象最为明显与裸子植物的单受精就是一个主要的例子。在针叶树胚胎发育过程中观察到的重要差异还包括:(1)前胚胎发育过程中没有细胞分裂的核复制，而不是被子植物合子中常见的初始不对称细胞分裂;(2)多胚现象的频繁发生;(3)胚胎早期发育过程中管细胞的分化;(4)胚发育后期多子叶的形成[2]。这些差异暗示两类植物在胚胎发育的分子调控上一定存在差异。松果体代表针叶树科中最大的属，在裸子植物中也是北半球分布最广的树木属。在各种数据库中都有关于几种松树的转录基因库的大量信息[3.，4]。然而，迄今为止使用这些信息的大多数转录组学研究都集中在抗逆性/耐受性和木材发育上[5- - - - - -8]。最全面的成绩单分析研究松果体胚胎已经在p . taeda，其中约有68,700项无害环境技术是由体细胞和受精卵产生的[9]。作者认为被子植物和裸子植物胚胎发育途径的差异主要是由于空间和时间基因表达的差异控制以及独特蛋白的表达。然而，在胚胎发育的不同阶段，转录基因的表达动态尚未研究。在答:芥斯宾塞等。[10]得出结论，就几个胚胎阶段的总体转录谱而言，时间表达差异比空间差异更显著。最近，挪威云杉体细胞胚胎发生过程中的差异基因表达(p .冷杉属)已经用微阵列进行了探测p . taeda序列，它允许识别调节与模式形成和分化相关的假定过程的分子事件[11]。迄今为止，对针叶树的研究大多依赖于体细胞胚胎，虽然体细胞胚胎发生是研究针叶树胚胎学的一个有用的实验模型系统，但人们也认识到，体细胞胚胎发生过程提供的条件在体外培养，如合成生长素的提供，可以对转录谱产生影响[12]。大多数对生长素的发育反应似乎是通过基因表达的改变介导的，外用生长素对植物的生长发育产生深远的影响[j]。13]。此外，体细胞胚胎形态异常也有报道p .松树，通常在含有生长素的培养基上诱导[14]。虽然合子胚胎发生是通常比较体细胞胚胎发生的模式，但很少研究合子胚胎的发育，因为从未成熟的针叶树种子中分离合子胚胎在技术上具有挑战性，特别是在胚胎发育的早期阶段[15]。

在本研究中，我们鉴定了p .松树从早期胚胎发生到胚胎成熟的不同发育阶段分离的合子胚胎，使用含有约25,000个独特cDNA的火炬松cDNA微阵列[7]，目的是识别与特定发育阶段相关的转录本和生物过程，并强调早期胚胎发育。据我们所知，这是第一次对松树合子胚胎发生的转录谱进行全基因组研究。我们的方法揭示了在关键发育过程的表观遗传和转录控制中可能起作用的主要调控基因。对照转录组学分析答:芥胚胎发生模型[16]进一步凸显了这些监管功能。

微阵列分析松果体松树发育中的胚胎转录组

转录水平动态p .松树利用PtGen2 cDNA芯片(GPL11184)对代表胚胎发育5个阶段的样本进行杂交，分析合子胚胎的发生情况。根据先前监测海松胚胎发育的研究[15，17]，我们在五个时间点分离了优势合子胚胎，代表胚胎发育的连续阶段，分为早期、前子叶、早期子叶、子叶和成熟胚胎(图2)1）.PtGen2微阵列包含25848个(除去缓冲空白和重复点共26496个)cDNA克隆扩增子，这些扩增子来自于火炬松(p . taeda)根和针组织;在其构建中未使用胚胎组织[7，18]。利用该阵列与来自不同的目标样品进行跨种杂交松果体物种,包括p .松树已被证明[7，18]。事实上，火炬松cDNA阵列也已成功用于其他针叶树属的基因表达分析[11，19]。我们的研究共使用了30张微阵列载玻片(图2)1A).归一化前后表达式值的箱形图证实数据归一化成功(附加文件1）.通过使用Pearson相关和平均连锁的分层聚类分析验证重复之间的可重复性，证明了微阵列数据集的质量1B).除第15天外，同一时间点采集的样本聚集在一起，与其他时间段的样本分开。技术重复总是聚集在一起，其相关值在0.7 ~ 0.85之间。总体而言，观察到技术重复与同一时间点采集的样本(生物重复)之间的密切关系，而来自遥远时间点的样本则表现出更大的变异性。

为了改进后续的标注，p . taeda与PtGen2阵列上发现的cdna的3 '和5 '端对应的ESTs被用于BLASTX对SustainPineDB的搜索，其中包含一组不重复的转录本p .松树(Unigenes)。使用3′或5′端序列将10922个点与同一基因对齐。另外6911个位点只有一个单端序列(大多数情况下是3 '端)与单个单基因对齐。根据使用的是3′端序列还是5′端序列，有3105个位点与不同的unigenes对应，因此这些位点与6210个unigenes相关。在所有情况下，重复的点都与相同的单基因对齐，在少数情况下，不同的克隆与相同的单基因对齐。总共有20938个位点与来自SustainPineDB (Additional file)的14996个不同的单基因比对2）.相比之下，在5294个微阵列位点上，在SustainPineDB中没有发现明显的比对。在这种情况下，p . taeda如果有3 ' EST序列，则在注释期间使用。或者，使用5 ' EST序列。14996个基因加上5294个基因p . taedaEST末端序列用BLASTX与NCBI蛋白数据库(nr)比对进行注释(e值< 1e-10)。共找到13280个同源基因和3482个cDNA克隆(附加文件)3.）.同一性(匹配/比对长度)分布峰值在75%。在最热门的点击中，47.4%对应于云杉sitchensis,10.3%,至葡萄3.4%萝藦。对于任何一个物种，剩余的命中率都不超过2.5%(附加文件)4）.基因本体术语(http://www.geneontology.org/) [20.]随后使用Blast2GO与12659个单基因和EST序列关联3.）.标注的GO术语等级从2级到11级，集中在7级左右。大部分无法标注的序列都小于1kb(附加文件)4）.我们还使用BLASTX (e值< 1e-10)来查看答:芥对应于每一个p .松树unigenes。共有13265个unigenes与7732个对齐答:芥从而提供基于TAIR映射的植物本体、途径和基因家族注释(http://www.arabidopsis.org)(附加文件3.）.

胚胎发育过程中基因功能类别的时程分析

为了全面了解与合子胚胎发育不同阶段显著相关的过程和功能，我们首先使用maSigFun中的回归模型对表达转录本进行了功能评估，以分析时间过程微阵列实验。在这项分析中，随着时间的推移，基因以相同的模式显著改变转录物水平的功能类别表明，它们内部存在高水平的共表达，并与胚胎发生过程有关[21]。差异表达(FDR < 0.01)在103个功能类别中被注意到，这些功能类别聚集在9个谱中(图2)2）.在相同的发育阶段，一些基因的表达趋势相似。例如，基因型1、6和7在胚胎发生早期被上调，而基因型3、4和5在同一时期被下调。然而，这些剖面在这些阶段都表现出不同的褶皱变化，使它们彼此区别开来。在附加文件中描述了根据它们的本体论注释分配到每个类别的基因3.。每个概要文件中的类别数量从2到29不等。除了氧化石墨烯的生物过程、功能和细胞成分外，鉴定的功能类别还包括18个EC编号、6个途径、1个基因家族(细胞骨架图1)和7个植物本体术语。然而，应该指出的是，在所有情况下，EC编号对于GO术语都是多余的。

氨基酸转运和代谢以及核苷酸代谢是胚胎早期最普遍的功能类别，但在随后的发育阶段呈逐渐下降的趋势(图1和图6)。位于碳和氮代谢之间重要分支点的酶，如谷氨酸脱氢酶，或参与碳和能量代谢中心代谢产物的生物合成的酶，如乙酰辅酶a，或在谷胱甘肽代谢中，如谷胱甘肽硫酯酶，也普遍存在于这些转录谱中。与GDP-甘露糖代谢相关的功能，特别是与GDP- l -聚焦生物合成相关的GDP-甘露糖4,6-脱水酶活性，在胚胎发生早期也被上调。细胞骨架基因家族成员在胚胎发生早期被上调，细胞壁成分的生物合成途径，如dTDP/ udp - l -鼠李糖也被上调(图2)2）.

图谱9，其中功能类别在胚胎早期到中期高度表达，但在成熟胚胎(第25天)时急剧下调，特别令人感兴趣。在该图谱中发现的所有共表达类别似乎都与表观遗传调控机制有关，包括染色质沉默的维持、组蛋白乙酰化或甲基化的调控以及DNA甲基化的调控。相反，在图3中，小rna介导的过程作为共表达类别的一个子集普遍存在，其表达趋势在胚胎发生后期增加。包括miRNA代谢过程、siRNA与miRNA结合、miRNA沉默基因等。大量共表达类别遵循类似的趋势，包括与嘌呤分解代谢循环活动相关的功能，包括黄嘌呤脱氢酶和黄嘌呤氧化酶，以及鸟嘌呤核苷酸降解途径。DNA复制和修复过程，这表明高DNA复制率，也被确定在这个剖面与生物合成的结构成分，如二磷酸植基。在图谱2中，与图谱3不同的是，它在胚胎成熟阶段的表达持续增加，与脂肪酸代谢相关的过程明显被过度代表，但其他机制，如与极化生长、染色质组织和细胞分裂素的精细调节相关的机制也存在。

胚胎发育过程中的转录谱分析

对于单个基因的分析，使用maSigPro提取差异转录序列[22]。该方法首先调整全局回归模型来识别与时间相关的转录差异序列，然后应用变量选择策略来研究组间差异并找到显著不同的谱。在胚胎发育过程中，共有3081个位点被差异转录p .松树(fdr < 0.001;额外的文件5）.其中，与2个单基因相关的位点有384个，与1个单基因相关的位点有2210个，与任何单基因不相关的位点有487个。然而，由于部分位点映射到相同的单基因，并且唯一的单基因数量为2633个，因此我们分析的差异转录序列总数为3120个(2633个单基因和487个ESTs)。其中2814个序列在NCBI nr数据库中发现同源物，2161个序列在NCBI nr数据库中发现同源物答:芥TAIR10数据库。我们能够将GO术语与这些序列中的2485个相关联(附加文件)5）.

根据5个发育时间点的转录水平，3081个位点可分为6个簇，分别由796个、631个、812个、5555个、125个和162个位点组成。这些位点分别映射到837、646、846、594、134和170个unigenes和ESTs，总共3227个差异转录序列(图2)3.）.映射到多个位点的一些独特基因分在多个簇中。6个转录谱可进一步按模式分为向上(簇1和3)、向下(簇2和4)、向上-向下-向上(簇5)和向下-向上-向下(簇6)，描述了从早期胚胎发生到成熟胚胎的变化。在整个胚胎发育过程中，有1683个转录本(52.2%)被上调，1240个转录本(38.4%)被下调。其余170个(5.3%)和134个(4.2%)转录本分别在中间时间点(簇5和6)上调或下调，对应于前子叶和早期子叶发育阶段。

在发育过渡期间，转录谱的变化在最极端的时间点，即第0天→第5天和第15天→第25天，比中间时间点，特别是第11天→第15天的变化更为明显，涉及的基因数量也更多。考虑到第11天和第15天之间的时间窗口较短，以及早期子叶胚胎和子叶胚胎之间的发育接近性，这并不特别令人惊讶，它们的区别主要在于预先形成的器官的扩大。

我们评估了生物过程GO项在各聚类中的分布(图1)3.）.当转录本数小于20时，词条的连接水平较优。氧化石墨烯项分析表明，氧化石墨烯的代谢过程氧化还原在大多数聚类中是过度代表的，其次是对压力的反应。例外是簇5和6，它们包含相对较少的注释基因。然而，某些氧化石墨烯术语仅在特定的集群中被确定，这表明它们与胚胎发育的特定时期有关。例如,解剖结构发育(包括器官和拍摄形态发生,光形态发生),多细胞生物发育(包括胚胎发育)都与早期胚胎发生有关(簇4)。Post-embryonic发展，与繁殖有关的发育过程;和生殖结构发育(包括种子、果实、胚和花的发育)被限制在第1和第3簇上，代表转录本主要从第11-15天积累到成熟胚。Organ发展成绩单(包括post-embryonic器官茎、根发育、叶片衰老)也存在于第1和第3簇，第6簇也存在。虽然更一般的过程细胞对刺激的反应，对化学物质的反应,压力或非生物刺激出现在所有集群中，对激素刺激的反应仅与中期和晚期胚胎发生阶段(集群1和6)相关，而对无机物质的反应83个差异表达转录本仅在裸子植物中被发现(注释在Additional file ?5）.这些转录本与53个单基因和26个ESTs相对应，因为一些转录本与多个单基因相关联。平均而言，每簇中约有3%的序列对应于裸子植物特异性序列，它们在胚胎发生早期和晚期之间分布均匀。这些序列很少与NCBI序列相似。

与表观遗传调控相关的差异转录谱

我们将差异转录基因的分析重点放在鉴定可能驱动胚胎转录组表达的主调控因子上，重点放在可能影响发育的表观遗传调控因子和转录因子上。表中列出了24个被选中的基因，它们的注释表明它们参与了表观遗传调控1。一些与染色质重塑相关的基因，包括组蛋白翻译后修饰、DNA甲基化和小RNA生物发生和加工的调控，在胚胎发育的不同阶段都被上调。然而，表中71%的转录本1在簇1和簇4中发现，它们在胚胎发育的早期和晚期表现出相反的转录模式。

四组蛋白去乙酰化酶(HDA)基因在差异表达的转录本中被鉴定出来。其中两个，HDA8和HDA9，从胚胎早期到子叶胚期，转录增加(簇1);HD2C，属于只存在于植物中的一类[23]显示相反的转录谱(聚类4)。HDA2在胚胎发生中期，特别是在子叶前期，表达上调(聚类6)。

与组蛋白甲基化相关的基因也被鉴定出不同的转录谱。的假定同源词答:芥polycomb-group (Pc-G)基因卷曲的叶子（CLF)是多梳抑制复合体2 (PRC2)的一部分，在胚胎发育中后期被上调(簇1)Ash1同源物2基因(ASHH2)，该基因编码一种具有组蛋白赖氨酸n -甲基转移酶活性的蛋白质，该蛋白质通过H3K36三甲基化调节基因表达[24]，表明成熟过程中转录增加(簇3)SU(VAR)3-9同源物（SUVH1)在早期子叶胚胎中被特异性上调(簇6)。几个被认为与SWI2/SNF2染色质重塑atp酶同源的基因调节基因组区域对转录机制的可及性[25]是差异表达的。其中两个在集群1中被发现，即CHC1和RAD5，而第三个，浓密的增长（BSH)是酵母SNF5 (SWI/SNF染色质重塑复合体的亚基)的一个假定的同源基因，在胚胎发生晚期被上调(集群3)INO80它是SNF2 atp超家族的一员。

与DNA甲基化相关的基因转录本，DNA甲基化减少1（DDM1),DNA甲基转移酶1相关蛋白（DMAP1)，在胚胎发生早期通过前子叶胚胎阶段被上调(集群4)甲基化变异1（ORTH2 / VIM1)在群集5中被发现。

最后，一些与已知的小RNA生物发生、加工和功能调节因子同源的转录本在胚胎发育过程中也受到差异调节。在集群1中鉴定出了一个假定的DAWDLE (DDL)编码转录本。的假定同源物偏下性的LEAVES1（HYL1),DICER-LIKE1（DCL1)分别出现在聚类1和聚类2中。另外，两个假定的ARGONAUTE成绩单、AGO1和AGO9，在集群1中被发现。根据相反的转录记录，我们还发现了一个开花地点ca（葬礼)簇4中假定的同源物。

参与连续胚胎阶段到阶段转换的转录因子

研究了在连续胚胎发育阶段之间显示≥2倍差异的转录本，以确定可能与特定发育时期相关的基因(附加文件)6）.与转录谱和功能分类的时间过程分析一致，从第0天到第5天的表达变化最显著，有173个转录本特异性下调，78个转录本特异性上调，从第15天到第25天的表达变化最显著，分别有280个转录本特异性上调，139个转录本特异性下调(图2)4）.在第11天→第15天期间，只有4个转录本被特异性差异调节(上调)4)，没有基因在第0天→第5天、第5天→第11天、第11天→第15天同时上调，说明早期和晚期胚胎的转录组存在较大差异。

关注参与转录调控的基因，我们在表中确定了被注释为可能转录因子(TFs)的转录本2。在鉴定的23个tf中，约有2/3在特定的发育过渡时期上调，其中第15天→第25天差异表达的tf数量最多。bHLH，其次是NAC和MYB转录因子家族，在我们的分析中最具代表性。在胚胎发生早期(第0天→第5天)发现的最具特征的TF是GARP转录因子家族成员KANADI 2 (KAN2)。在胚胎发生的早期和中期，也观察到一个假定的同源基因的上调AINTEGUMENTA（蚂蚁)(表2)，而假定的YABBY2转录本在早期子叶期到子叶期的过渡过程中似乎很重要。假定的三个同源词南汽记录,即不结盟运动/NARS2，拟南芥nac结构域含有蛋白75（ANAC075),ATAF1，在第15天→D25过渡期间表达上调(表2)2)，而假定的bHLH成绩单,抗亮氨酸3（ILR3)，在第15天→第25天强烈上调。值得注意的是，从第0天到第5天的强烈下调(20.9倍)生长素反应因子，ARF16这可能有助于强调生长素反应机制在早期胚胎中的相关性。在这个早期阶段，一个假定的bZIP TF的同源编码，分生组织3的花过渡（FTM3），假设WRKY28被发现是下调的。在第15天→第25天的过渡中特异性下调的tf包括bHLH超家族的成员，例如推定的法玛同系物。此外，植物特异性TFs的PLATZ家族同源物和HAP5B作用不明确或尚未明确的基因在这一时期也出现了下调。

胚胎发生过程中差异表达转录物的比较时程分析紫菀与拟南芥

我们比较了参与胚胎发育的基因的转录谱答:芥和p .松树为了寻找推测的同源转录本水平之间的相关性。对于我们认为物种之间相等的每个发育阶段，基因被绘制在散点图中答:芥和p .松树作为坐标的表达式值(附加文件7）.线性回归分析表明，相关性最高的是第1期(Day0和球形胚)和最后期(Day25和成熟胚)，相关性最低的是第3期(Day11和鱼雷胚)。

为了比较两个物种之间的转录谱，在每个阶段的转录水平都发生了变化与每个物种在胚胎发生过程中的平均值被量化为每个数据系列。224个基因的转录谱p .松树和答:芥Pearson相关性大于0.90(附加文件8）.Agrigo富集分析[j]26在206个基因中有GO注释答:芥TAIR10数据库显示，7种生物学功能被显著过度代表(FDR < 0.001;额外的文件9),具体分解代谢的过程，细胞成分组织，发展的过程，多细胞有机体过程和几个代谢过程类别。

在两个物种之间具有高度相关转录谱的基因中，我们发现:(1)6个胚胎缺陷基因,即GNOM/EMB30(AT1G13980),S-腺苷-同型半胱氨酸水解酶/EMB1395(AT4G13940),秘书同源物2/EMB2289(AT2G31530)， EMB2296 (AT2G18020)，核糖体蛋白1 / EMB2207 (AT2G18020)和细胞分裂缺陷1 / emb101(AT2G39770);(ii)几种差异表达的基因发展的过程一类先前被证明会影响胚胎发育的答:芥通过转录和调控过程，比如SKP1-Cullin/CDC53-F-box(也叶片卷曲反应，电感电容电阻测量),ETHYLENE-INSENSITIVE 2（EIN2;AT5G03280);(三)两种相关表观遗传调控因子;IRA 1的多拷贝抑制因子（MSI1;AT5G58230)和泛素载体蛋白（UBC2;AT2G02760)。

微阵列数据验证

为了验证微阵列表达数据，一组10个差异转录基因被推测参与表观遗传和转录调控(表2)1和2)进行RT-qPCR分析。所选基因在合子胚胎发育过程中表现出不同的转录谱。聚类1有5个基因，聚类6有2个基因，聚类2、3、4各有1个基因(图2)3.）.其中4个基因(Unigenes 1030、7237、27680和17529)为转录因子。获得每个基因的微阵列绿/红强度比(M值)和RT-qPCR转录水平(附加文件)10）.对两个数据序列的大小进行均衡化或归一化，计算每个时间点与平均值之间的折射率变化，然后进行对数₂转换为如图所示5。微阵列与RT-qPCR数据之间存在较高的Pearson相关系数，其范围为0.69 (Unigene 25311)至0.99 (Unigene 17529)。

合子胚胎分离是大多数植物物种面临的主要挑战，特别是在胚胎发育的早期阶段[16]。我们的方法允许按阶段快速分离胚胎，并在冷冻前验证结构完整性，从而能够收集足够数量的组织用于RNA分离。通过这种方式，我们能够进行第一个全基因组转录谱分析研究，该研究涵盖了针叶树物种合子胚胎发育的广泛时间窗口。

胚胎发育所需的许多基因不是胚胎特异性的，因为它们的基本功能在植物的整个生命周期中都需要[27]。等。[28报道称，在整个过程中，16000个基因中有活性答:芥种子发育，只有289种是种子特异性的。然而，基因活性的显著定量变化发生在特定的发育阶段。事实上，种子发育的每个阶段都有一组特有的基因，相对于其他阶段活跃的基因，这些基因要么是特异性的，要么是上调的[28]。我们发现氧化-还原代谢过程在大多数簇中都有过度的代表，这表明在发育中的胚胎中发生了高代谢活动。通过氧化还原代谢物(如谷胱甘肽)维持细胞氧化还原稳态，可能在胚胎发育中起重要作用。事实上，谷胱甘肽代谢在我们对差异表达功能类别(谷胱甘肽硫酯酶活性)和大量的转录本编码谷胱甘肽转移酶。这些观察结果似乎与先前的报告一致，即增殖细胞中谷胱甘肽的丰度通过在生长素运输和信号传导中的作用对茎和根分生组织的发育至关重要[29]。此外，它已经表明操纵谷胱甘肽代谢在在体外胚胎发生会影响胚胎的产量和质量，即在体细胞胚胎发生诱导过程中，还原型谷胱甘肽水平升高的环境会导致未成熟胚胎数量增加，而氧化程度更高的环境则会促进胚胎发育[30.]。

针叶树是裸子植物中的一个主要类群，与被子植物相比，针叶树的独特特征可能揭示出独特的基因和基因网络，从而进一步阐明植物胚胎发育及其进化意义，是研究胚胎发生的有趣对象。在我们的分析中，我们发现大约3%的差异调节转录本似乎是裸子植物所特有的，其中大多数编码未知或未表征的蛋白质。虽然这些未知的蛋白质可能在针叶树胚胎发育中起重要作用，但似乎胚胎的发生主要是由被子植物和裸子植物中一组相似的转录本的协调活动完成的，正如先前提出的那样[9]。这两个植物类群胚胎发育的差异特征可能主要是由于基因调控的差异。因此，我们的研究重点是分析在植物发育过程中起重要作用的基因表达调控机制，即表观遗传调控和转录因子的转录调控。

胚胎发育的表观遗传调控途径

组蛋白的共价修饰、DNA甲基化、染色质重塑酶和小rna，以及其他因素，通过调节DNA的获取和定义不同的染色质状态，最终决定基因组序列的选择性读出，在基因表达中发挥核心作用[31]。

在我们对松树胚胎发生功能类别的时间过程分析中，染色质沉默的维持、特定的组蛋白翻译后修饰以及DNA甲基化和转位的调节在早期胚胎发生过程中被认为是共同调节的功能。在胚胎发生早期，co-调控的功能类别以及鉴定的差异调控转录本，指出了与转座因子(TE)控制相关的基因沉默机制的重要性。事实上，DNA甲基化和异染色质维持在co .的分析中都得到了强调-规范的功能类别和向上-对假定的规则DDM1这是一个关键的调节异色形成答:芥TE所需-特异性DNA甲基化[32]。在假设的FCA中观察到相同的转录谱，最近涉及与转座子、反转录转座子和通常被RNA沉默的分散重复相关的RNA序列的调控-定向DNA甲基化途径[33]。也是一个假定的ORTH2/VIM1在细胞分裂过程中，DNA甲基化状态与染色质结构的建立/维持有关[34)。-尽管在整个胚胎发育过程中表现出明显的转录谱下降，但在早期胚胎发生过程中受到调控-在胚胎发生中期受到特殊调控。

DCL1［35]，它的同系物是up-在Day0样品中受到调节，最近被认为在TAS中发挥作用-衍生小干扰RNA-通过直接加工触发DNA甲基化助教基因转录物[36]。DCL1是早在8岁时细胞分化事件所必需的-细胞阶段答:芥胚胎(37]，在那里它参与了早期的胚胎模式。DCL1通过其对miRNA生物发生的作用，阻止了在胚胎发生后期促进分化的miRNA靶标，即转录因子的积累[37]。

同时发现组蛋白H4乙酰化负调控和H3K9甲基化正调控的证据-功能分类的过程分析表明，在胚胎发生早期有转录抑制的趋势。具体来说，在高表达基因中发现H3K9三甲基化，而不是H3K9甲基化或二甲基化[38]。相反，高乙酰化组蛋白通常与基因激活有关，而低乙酰化组蛋白则与基因抑制有关[23]。我们的研究结果表明，在胚胎发生过程中，不同的基因群通过不同类型的去乙酰化酶的作用成为表观遗传调控的目标，例如在早期胚胎发生中推测的HD2C与推测HDA8和HDA9在胚胎发生后期的作用。组蛋白去乙酰化酶(hdac)的抑制答:芥是否已知会影响胚胎发育以及种子的表达-相关基因，包括转录因子[39，40]。在另一种针叶树，挪威云杉，Uddenberg等。[41报道了类似的行为，这表明乙酰化模式的变化与胚胎发生水平之间存在联系-相关基因表达。

的同系物SUVH1它编码H3赖氨酸-9特异性组蛋白-赖氨酸N-甲基转移酶上升了-在早期子叶胚(簇6)中受到特异性调控，表明可能在转录抑制染色质的组织中起作用[42]。总之，这些结果表明，不同的井-已知的基因沉默机制在松树的早期胚胎发育过程中是活跃的。

从胚胎发生的中期到晚期，大的染色质重塑事件似乎也发生在针叶树胚胎中，这表明几个表观染色质的转录显著增加-重塑atp酶，如推测的CHC1, RAD5和BSH。另一个基因出现了-中期调控-到胚胎发生晚期(簇1)是同源的CLFpolycomb的用法和样例-集团(电脑-G)蛋白质，通过抑制细胞中的基因表达来调节植物和动物的许多发育过程-通过三甲基化组蛋白H3赖氨酸27 (H3K27me3)的特异性方式[43]。CLF最初被认为是花同源基因的抑制因子[44]，但其他基因如杯-形子叶2（CUC2),PIN1被描述为Pc-G靶基因[45，46]。中联科利例如，基因是由PIN1，对胚芽顶端分生组织的建立和两个子叶的形成至关重要，可能是通过阻止子叶间区域的细胞增殖和子叶的生长[47]。的Picea abies PaNAC01的同源词中联科利，也受极性生长素运输的调节，并与茎尖分生组织的分化和分离子叶的形成有关[48]。与这些发现一致，H3K27me3的缺失拟南芥被证明与生长素反应相关的基因调控不当有关[49]。似乎在我们的结果中，观察到的差异调节CLF和其他染色质-基因重塑与其中一些过程的调控有关。

在关于小RNA途径的研究中，一些与已知的小RNA生物发生、加工和功能调节因子同源的转录本在松树胚胎发育过程中存在差异表达，包括一个假定的DDL，大部分是向上的-在胚胎发生中后期受到调控，可以参与mirna和内源性sirna的生物发生。DDL不直接影响mirna的转录，而是通过其他蛋白质，如DCL，通过促进它们对pri的获取或识别而起作用-microrna [50]。来自较长初级转录本(pri)的mirna的加工-miRNAs)需要几种蛋白质的活性，包括DCL1和double-链RNA-HYPONASTIC LEAVES1结合蛋白[35，51，52]，发现了一个假定的同源词-在与假定的DDL相同的发育时期受到调节。

Argonaute (AGO)蛋白是rna诱导沉默复合体(RISC)的一部分，它结合小rna并导致基因沉默。两个公认的前类似于AGO1和AGO9，在晚松胚转录组中有高度的代表性。AGO9属于一个所有成员识别24个核苷酸的小干扰rna并在转录水平上沉默te和其他重复序列的系统发育分支[53]。相比之下，AGO1的同源物主要是miRNA活性的介质[54，55]，尽管AGO1也参与了rdr6依赖性sirna的产生[56，57]。综上所述，不同的基因沉默机制似乎在松木胚胎发育的不同阶段更为活跃。在胚胎发生早期，以TE和重复序列为目标的机制似乎占主导地位，而prc2介导的参与特定发育过程(如子叶形成)的基因抑制似乎在胚胎发生晚期更为相关。相反，与sRNA通路相关的基因被发现在胚胎发育的各个阶段受到不同的调控。

转录调控因子和生长素介导的事件

当染色质结构允许表达的tf进入其结合位点时，这些蛋白质在基因表达的调节中起主要作用。的假定同源词ARF16从早期胚胎到前子叶发育阶段，转录量急剧下降，表明它可能在早期松树胚胎发生中起主要作用。ARFs是生长素调控基因表达的关键调控因子[58通过结合早期生长素反应基因的启动子来激活或抑制目标基因。ARF16和ARF10表现出抑制作用WOX5转录并将其限制在答:芥根静止中心[59]需要维持多能小柱干细胞的地方[60]。有趣的是，除了在早期胚胎发育中表达外，假定ARF16在松木(簇5)中，从早期子叶到成熟胚的转录谱也呈现出增加的趋势，这与上文描述的相似ARF16在答:芥胚胎(61]。另一个在早期松树胚胎中上调的基因被认为是一种同源基因生长素RESISTANT1（AUX1)，它编码生长素内流载体。与ATP结合盒亚家族B (ABCB)转运蛋白和PIN蛋白一起，AUX1载体是生长素极性转运的主要协调者。的同系物N-myc下调样1（NDL1)，它可能通过调节生长素运输载体蛋白如PIN2和AUX1来调节生长素运输[62]，在早期松树胚胎中也上调表达。

如果这些基因产物在裸子植物胚胎发生中具有保守功能，那么生长素与具有明确时空表达模式的转录因子之间的相互作用对于松木胚胎模式的建立至关重要。例如，假定的ARF16在松树胚胎中的表达模式似乎与已描述的一致ARF16在答:芥在胚胎发生的早期阶段，它通过参与根尖分生组织形成的启动，参与根尖-基部模式的建立。

KANADI蛋白(KAN2)在被子植物早期胚胎发生中起重要作用，可能是通过调节PIN蛋白的极性表达来调节生长素的流动。63]。我们的研究表明，KAN2同源物可能在早期松树胚胎发生到前子叶胚胎阶段的过渡中起重要作用。有趣的是，PIN3同源物的差异调控(附加文件)6)，一种生长素外排载体，是GNOM的同源物(簇4，附加文件)5)，对内体腔室和质膜之间PIN蛋白的再循环很重要[64]，在同一发育时期观察到的结果与假定的KAN2的这种作用是一致的。在心脏阶段答:芥胚胎发生,KAN2以及KAN1和KAN3,在子叶初生叶的背面基部也有类似的表达模式[63]。同样的作者提出，沿中央-外周轴的模式形成是生长素与KANADI和III类HD-Zip转录因子相互作用的结果。

Eshed等[65提出最初的不对称叶片发育，主要由KANADI和PHABULOSA (PHB)样基因之间的相互拮抗调节，被翻译成极性的YABBY表达，随后有助于背面细胞的命运和背面/正面并置介导的叶片扩张。一个假定的淡水螯虾基因家族成员，YABBY2，在松木中，从子叶早期开始显著上调，这可能与参与了背面细胞命运的决定和沿背面-正面边界的叶面生长相一致[66，67]。另一个文本，一个假定的AS2 / LOB，编码一种植物特异性蛋白，参与决定叶片原基的正面-背面极性[68- - - - - -70]和1类KNOX同源盒基因的抑制[71]，在第15天→第25天过渡期间下调。可能参与顶端分生组织形成和功能的松基因也有差异表达，包括一个推测的NARS2 /南，在松木胚胎发生晚期显著上调。NARS2 /南在器官边界和胚芽顶端分生组织的形成中都有涉及。例如，拟南芥［72),CUC2［73还有牵牛花无顶端分生组织(名称)［74]，表达于花器官原基之间的边界和子叶之间的边界。然而，基因突变CUC1，CUC2,不结盟运动也影响茎尖分生组织的形成[72- - - - - -75]。AINTEGUMENTA(蚂蚁)参与了茎尖分生组织功能的调控[76]。推定的松树蚂蚁同源基因可能在子叶前胚向早期子叶胚过渡的过程中起重要作用，从这一时期同源基因的上调可以看出。含有基因突变的植物蚂蚁基因对极性生长素运输中断的敏感性增加[77这表明它在胚胎发生过程中的作用可能与生长素有关。

假定同源物剖面之间的相关性答:芥和p .松树胚胎发生

通过对两个远亲种的胚胎阶段进行匹配，对模式被子植物进行比较分析答:芥和p .松树我们认为被子植物和裸子植物胚胎发育的几个重要过程是相同的。事实上，一些已知的基因对于答:芥胚胎发育，包括胚胎缺陷我们的分析强调了基因和各自假定的松树同源物。其中一个基因，GNOM（EMB30)，在早期起着关键作用答:芥胚胎形成需要胚胎发生[78]。它在PIN1(和PIN2)到基基质膜的内体再循环中的功能是胚胎发生过程中顶端-基部极性建立的基础。通过组织特异性表达GNOM在gnom突变胚胎，Wolters等。[78表明顶端和基部的胚胎组织都依赖于GNOM并提出gnom依赖的PIN重定位(细胞自主)和sink驱动的生长素运输(非细胞自主)触发胚胎模式形成过程中生长素运输路线的启动。本研究通过对比分析确定了生长素介导反应的其他相关基因。例如，LCR已被证明在叶片发育过程中负调控几种生长素响应基因[79]， EIN2与AUX1和gnm一起参与了根尖生长素浓度梯度的建立，这对根毛在表皮层内的定位很重要[80]。

最后，有两个相关的表观遗传调控因子，其转录水平随着时间的推移以高度相似的方式下降，表明在发育中具有相同的作用答:芥和p .松树胚胎:MSI1和UBC2。MSI1是植物正常种子发育所需的PRC2复合体的核心蛋白[81]。需要MSI1来维持同源基因的正确的时间和器官特异性表达，包括无性生殖的和APETALA2，以及在SOC1染色质中建立表观遗传标记，如H3K4二甲基化和H3K9乙酰化[82- - - - - -84]。UBC2参与调节组蛋白H2B的泛素化和控制开花时间答:芥通过调控开花位点C/ MADS影响开花染色质[85]。

综上所述，我们的结果表明，特征转录变化与每个发育时期有关，正如之前观察到的答:芥胚胎发生(10，28]，并表明胚胎发生过程中的主要事件是由被子植物中这些过程的调节剂的假定同源物精心安排的。我们的分析在一定程度上是有限的，因为本研究中使用的胚胎发育阶段没有涵盖区分被子植物和裸子植物胚胎发生的主要事件的某些时期，例如受精后的第一次合子分裂和多胚。因此，在确定优势胚胎之前，进一步的工作应该集中在松树合子胚胎发育的最早阶段。扩大针叶树基因组的知识对于进一步了解裸子植物胚胎发生特征的分子基础也至关重要。进一步深入了解这些物种胚胎发育的分子调控，将对针叶树体细胞胚胎发生的改良及其营养繁殖具有重要意义。

植物材料

未成熟雌球果松果体松树河中的小岛。1970 ~ 1975年，在葡萄牙Escaroupim国家森林(经度8°44′w，北纬39°4′n)嫁接建立的无性繁殖果园中，随机采集了开放授粉的树木。在2007年7月至8月的5个不同的时间点进行了锥体收获，每个日期从3-5棵树中选择。每次收获球果后，立即从球果中取出种子，打开，露出巨形植物。在体视显微镜下用细钳和手术刀切除优势胚胎。分离的胚胎被迅速评估发育阶段，并立即在液态氮中冷冻₂根据舞台的不同，进入不同的水池。五个时间点分别为第一次采收日期的第0天，第5、11、15和25天，根据第0天之后的采收日期。根据先前描述的胚胎分期系统[15]， Day0包括发育阶段T0、T1和T2的胚胎(早期胚胎期)，Day5包括T3和T4期(前子叶胚胎)，Day11包括T4B期(早期子叶胚胎)，Day15包括T5期(子叶胚胎)，Day25包括T7期胚胎(成熟胚胎)(图1)1A).为每个时间点准备几个单独的池，每个池中每个阶段含有30-60个受精卵，并保存在单独的试管中，作为后续分析的生物重复。

RNA的分离和扩增

根据制造商的说明，使用RNeasyPlant Mini Kit (Qiagen, Valencia CA, USA)从每个冷冻受精卵池中提取总RNA，在分离过程中使用柱上DNase I酶切和RNase-Free DNase Set (Qiagen)。采用NanoDrop 1000 (Thermo Scientific, Waltham, Massachusetts, USA)分光光度法测定RNA浓度，琼脂糖凝胶电泳测定RNA质量和完整性。

在检查RNA完整性并在进行RNA扩增之前排除任何DNA污染的可能性后，根据制造商的说明，在使用TURBO DNase (Ambion, Foster City CA, USA)进行扩增之前，立即对样品进行额外的DNA酶消化处理。使用Amino Allyl MessageAmp II aRNA扩增试剂盒(Ambion)进行RNA扩增。简而言之，该反应包括一个反转录步骤，在42°C下反应2小时生成含有T7启动子序列的第一链cDNA，然后在16°C下合成2小时的第二链cDNA。cDNA纯化后，在体外在37°C下转录12 - 14h，生成反义扩增RNA。根据制造商的说明，在两轮连续扩增中扩增了总共100 ng的总RNA，其中氨基烯丙基UTP核苷酸的掺入仅在第二轮扩增中进行。

微阵列实验与信号采集

本研究中使用的PtGen2 cDNA微阵列先前已被描述[7]，有关该阵列的更多信息可从NCBI的基因表达综合数据库(GEO)获得，平台登录号为GPL11184。实验样品由图中所示的五个时间点的合子胚胎扩增的RNA组成1A.每个时间点分析3个生物重复，每个样本进行2个技术重复。杂交参考样本由一个包含来自所有五个时间点的等量总RNA的池组成。荧光标记、杂交、杂交前后清洗、扫描和数据处理均按照Lorenz等人的方法进行。86]，除了使用碎片化试剂(Ambion, Cat.)进行碎片化后，用2 μg的cy标记的扩增RNA与阵列杂交。# AM8740)。微阵列扫描使用ProScanArrayTM共聚焦扫描仪(Perkin Elmer, Waltham, MA)，配备532 nm和635 nm激光器。使用ImaGene (Bio-Discovery Inc.， El Segundo, CA, USA)处理原始荧光数据。

低质量的原始数据和空白对照按上文所述进行过滤[7，87，88]。

微阵列数据规范化和注释

过滤后的原始数据使用Bioconductor Limma软件包进行分析[89]。根据Smyth和Speed对阵列内(Print-tip黄土)和阵列间(循环黄土)的数据进行归一化处理[90]，强度比转化为对数₂值。在归一化过程中，不应用背景校正，重复探针独立处理。

使用MeV程序对样本进行聚类[91]。纳入更新的成绩单信息为p .松树，将PtGen2 cDNA微阵列(GPL11184)重新注释到p .松树数据库，SustainPineDB (Version 2);http://www.scbi.uma.es/sustainpine)，其中包含92 478份来自海松的独特转录本。代表用于生成该阵列的cDNA克隆的3 '和5 '端EST序列[5]用BLASTX (e值1e-10)进行比对识别p .松树每个点对应的转录本。

转录档案和功能分类分析

p .松树使用unigenes或ESTs作为BLASTX查询，鉴定NCBI蛋白(nr)或同源序列(e值< 1e-10)答:芥TAIR10数据库。使用Blast2GO收集相应的基因本体、酶编码、植物本体、途径和基因家族术语[92]，也用于富集分析。

使用Bioconductor maSigPro包从不同时间点的数据和复制中鉴定差异转录基因[22]。我们定义了一个二次回归模型来识别五个时间点上的差异转录基因，前提是该模型应该有足够的能力来分析减少的时间点[22]。错误发现率(FDR)为0.01，用于鉴定显著差异转录基因。maSigPro逐步回归调整为使用“两种方法”。α α 0.01, R²截断值拟合为0.7。高R的使用²已经提出了基因模型来捕捉生物学上有意义的表达趋势和谱差异[22]。使用基于相关距离的分层方法对显著差异转录基因进行聚类。在开花植物(或其它更远的植物)中不存在的独特基因分类单元)，但存在于裸子植物中，被注释为裸子植物特有。为了进行这种鉴定，在NCBI蛋白(nr)数据库中找到了属于裸子植物的完整蛋白列表分类单元生成(税号58021、3312、58020和58022)(GI list)，并从NCBI Entrez数据库(104709个蛋白)下载。我们使用了BLASTX选项negative_gilist或gilist，分别排除或限制与裸子植物的比较。使用maSigFun进行功能评估，maSigFun适合用于注释到类似功能类别的基因组的回归模型[21]。在这种情况下，分组是基于基因本体论术语、酶代码号、植物本体论术语、基因家族和途径，以及我们有重要模型的其他类别的选择。

时间过程微阵列分析的强度值答:芥从Xiang等人的工作中获得了胚胎发生。http://www2.bri.nrc.ca/plantembryo/) [16]，并进行归一化，通过减去一个恒定的背景值得到表达式值，如作者所述。对于在我们的研究中进行的比较分析答:芥球形、心形、鱼雷、弯曲和成熟的胚胎阶段被认为与第0天、第5天、第11天、第15天和第25天的胚胎样本相当p .松树,分别。对于每个基因，在每个阶段的转录水平值之间的折叠变化与计算各物种各阶段的平均值。自p .松树单一的假定同源物答:芥使用Pearson相关性比较转录谱。

使用实时RT-qPCR验证微阵列数据

RT-qPCR用于独立验证从微阵列数据中鉴定出的一组基因的转录谱。引物采用Primer3程序设计(http://primer3.sourceforge.net/)，使Tm为58°C，结合基因保守区域，扩增约100 bp(附加文件)10）.使用美国罗氏诊断公司(Roche Diagnostics, USA)的转录高保真cDNA合成试剂盒(Transcriptor High Fidelity cDNA Synthesis Kit)和Random Hexamer引物，用1.25 μg的总RNA合成cDNA进行反转录。随后在Roche LightCycler 480系统中按照标准循环条件进行定量PCR (qPCR)。采用美国罗氏诊断公司SYBR Green I Master试剂盒制备qPCR反应混合物，每个反应混合物含有1.8 μl cDNA，每个引物400 nM，总体积为20 μl。RNA样品经qPCR扩增后无信号，即排除DNA污染泛素引物。阴性对照(无模板)包含在所有运行中，并且在熔化曲线(解离曲线)中只检查一个特定峰的存在。从所有五个时间点提取的含有等量总RNA的池作为校准器和阳性对照。每个样品分析2个生物重复，每个生物重复分析3个技术重复。利用LinReg 11.3程序对PCR反应动力学进行回归分析，计算每对引物的效率[93]，通过Roche Lightcycler软件(附加文件)测定Ct值10）.每个转录本的相对数量通过delta-delta-Ct方法计算[94)使用泛素作为内参基因，该基因作为内源调控已被证实可靠p .松树体细胞胚胎发生[95]。

支持数据的可用性

支持本文结果的数据集可在NCBI GEO数据库中获得，登录号为GSE32551 (http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE32551）.

FCT通过SFRH/BD/32037/2006(给M. Simões)和SFRH/BD/79779/2011(给A.S. Rodrigues)、PTDC/AGR- GLP/102877/2008、P-KBBE/AGR- gpl /0001/2009和PEst-OE/EQB/LA0004/2011项目获得了资金支持。这项工作的一部分也得到了葡美合作与发展基金会的支持。来自国家研究机构 Agrária e Veterinária (INIAV)的Alexandre Aguiar和Isabel Carrasquinho因提供植物材料而获得认可。

作者指出

从属关系

1. iBET -生物实验研究所Tecnológica，厄伊拉斯，2780-901，葡萄牙

   josise J de Vega-Bartol, Marta Simões, Andreia S Rodrigues和camelia M Miguel

2. 技术研究所Química e Biológica，葡萄牙新里斯本大学República，葡萄牙奥伊拉斯，2780-157

   josise J de Vega-Bartol, Marta Simões, Andreia S Rodrigues和camelia M Miguel

3. 美国佐治亚大学沃内尔林业与自然资源学院，佐治亚州雅典，30602

   沃尔特·洛伦兹和杰弗里·迪安

4. 孟山都公司，邮编:700 Chesterfield Parkway West, Chesterfield, MO, 63017, USA

   罗伯·阿尔巴

相应的作者

对应到米格尔先生。

相互竞争的利益

作者宣称他们没有竞争利益。

作者的贡献

MS参与了实验设计、样品采集和制备、微阵列杂交、结果解释和论文撰写。JVB进行了生物信息学分析，参与了结果的解释和手稿的准备。WWL监督样品标记，微阵列杂交，并参加了手稿的关键审查。ASR有助于微阵列数据的验证。RA对实验设计做出了贡献，并为结果的解释提供了输入。JFDD提供了材料和设备，并参与了手稿的批判性审查。CM构思并监督项目，参与结果解读，撰写并编辑稿件。所有作者都阅读并批准了最终稿件。

josise J de Vega-Bartol, Marta Simões对这项工作也做出了同样的贡献。

开放获取本文由BioMed Central Ltd.授权发表。这是一篇开放获取的文章，在知识共享署名许可(https://creativecommons.org/licenses/by/2.0)，允许在任何媒体上不受限制地使用、分发和复制，前提是正确引用原创作品。

转载及权限