《Frontiers in Plant Science》：A comparative study of phyllostachys edulis and its dwarf variant phyllostachys edulis ‘Tubaeformis’ at the anatomical, transcriptomic, and DNA methylation levels
编辑推荐：
本研究通过整合解剖学、转录组学和全基因组甲基化测序（WGBS）技术，系统比较了毛竹（Phyllostachys edulis）及其矮化变种圣音竹（P. edulis ‘Tubaeformis’）的节间发育差异。研究发现圣音竹节间缩短源于细胞分裂与伸长异常，并首次揭示DNA甲基化（CG/CHG/CHH）通过动态调控生长调节因子（GRF）、细胞周期蛋白（Cyclin）和扩张蛋白（EXPA）等关键基因的表达时序，进而影响竹子节间伸长的表观遗传新机制。该研究为竹类植物矮化育种提供了重要的理论依据。
引言
毛竹（Phyllostachys edulis）作为禾本科（Poaceae）竹亚科（Bambusoideae）植物，以其惊人的生长速度著称，竹秆从笋芽到高达15–20米仅需45至60天。作为中国最重要的竹资源，其分布面积达473万公顷，占全国竹林总面积的70%。圣音竹（P. edulis f. Tubaeformis）是毛竹的一个矮化变种，其节间显著缩短，竹秆下部粗壮并向基部逐渐扩张呈喇叭形，因其独特的形态和用于制作乐器时发出的空灵音色而成为珍贵的观赏竹种。
竹笋的生长主要依赖于节间细胞的伸长和分裂。在节间内部，不同位置的组织处于不同的生长状态：基部组织以细胞分裂为主，中部组织以细胞伸长为主，而顶部组织则停止生长。先前研究表明，外源施加油菜素内酯合成抑制剂PPZ（propiconazole）或赤霉素合成抑制剂多效唑（PAC）会显著降低毛竹实生苗的株高和节间长度，而外源施加赤霉素（GA3）或生长素（NAA）则能促进其节间伸长。比较转录组学研究曾揭示，圣音竹中与细胞壁松弛相关酶、纤维素和木质素生物合成途径以及激素生物合成和信号转导相关基因的异常表达是其矮化的关键因素，但调控这些基因表达的表观遗传机制尚不清楚，且筛选出的关键功能基因尚未进行功能验证。
DNA甲基化作为一种重要的表观遗传修饰，在植物的非生物胁迫响应、笋芽生长发育、竹鞭发育和开花等方面已有研究。例如，在箣竹（Bonia amplexicaulis）茎秆生长过程中，快速生长阶段（ST2-ST4）的CG和CHG甲基化水平显著低于潜伏期（ST1）和平台期（ST5）。对不同高度毛竹笋单节间（基部和中部）的甲基化分析发现，2米高的竹笋在RNA甲基化水平上N6-甲基腺嘌呤（m6A）修饰比例较高，而4米高的竹笋在DNA甲基化水平上5-甲基胞嘧啶（5mC）水平较高。用RNA甲基化抑制剂（DZnepA）和DNA甲基化抑制剂（5-azaC）处理毛竹实生苗，其根长变短，但侧根数量增加。然而，圣音竹的节间缩短是否由DNA甲基化调控引起仍属未知。
本研究首先在解剖学水平比较了圣音竹与毛竹的节间，发现圣音竹节间缩短由细胞分裂和伸长异常引起。随后，在转录组水平鉴定出圣音竹中大量差异表达基因（DEGs），包括与DNA甲基化、细胞周期和细胞壁合成相关的基因，其异常表达可能导致圣音竹节间缩短。接着，在DNA甲基化水平比较了两种竹子，证明DNA甲基化在调控基因表达中扮演重要角色。最后，对在DNA甲基化调控下于毛竹中高表达的关键基因GRF10（Growth-regulating factor 10）进行了功能验证。该研究为理解圣音竹矮化的表观遗传调控机制提供了新的见解和证据。
材料与方法
样品采集自浙江省安吉县中国竹子博览园，毛竹与圣音竹的种植和管理条件一致。研究选取基部第15节间作为实验材料。具体而言，毛竹样品包括发育早期（M1，代表细胞分裂期，长度4厘米）和发育晚期（M2，代表细胞伸长期，长度15厘米）的节间基部。圣音竹则采集对应的S1和S2样品。部分样品用于石蜡切片（保存于FAA固定液），部分用于转录组和甲基化测序（液氮速冻后于-80°C保存）。
扫描电子显微镜（SEM）观察显示，成熟毛竹和圣音竹秆经2.5%戊二醛固定、PBS清洗、梯度乙醇脱水、醋酸异戊酯置换、临界点干燥和喷金后，在2 kV加速电压下利用JSM-6360LV扫描电镜成像。
石蜡切片制备与观察过程包括FAA固定超过48小时、梯度乙醇脱水（30%至100%）、二甲苯透明、石蜡包埋，使用Leica RM 2235旋转切片机以7-8微米厚度切片，番红-固绿染色，中性树胶封片后于Leica DM2500光学显微镜下观察拍照。
转录组测序使用RNAprep pure Plant Kit提取总RNA，经DNase I处理去除基因组DNA残留。通过琼脂糖凝胶电泳、Nanodrop 2000和Agilent 2200评估RNA质量。合格的总RNA用于构建链特异性RNA-seq文库（dUTP法），在Illumina Novaseq 6000平台上测序。原始数据经fastp进行质控（QC）和接头序列去除，使用Hisat2将质控后的reads比对到毛竹参考基因组，仅使用唯一比对到基因组的reads进行下游分析。基因表达量以FPKM（Fragments Per Kilobase Million）计算，采用DESeq2方法进行差异表达分析（≥1倍变化，FDR <0.01）。差异表达基因通过BLASTX比对Swiss-Prot、Pfam、Gene Ontology（GO）和KEGG数据库进行注释和分类。
全基因组重亚硫酸氢盐测序（WGBS）使用DNeasy Plant Mini Kit提取总基因组DNA（gDNA），经琼脂糖凝胶电泳和紫外分光光度法评估质量。gDNA通过BioRuptor超声破碎成100–500 bp片段，末端修复、腺苷化、接头连接后，使用Zymo Methylation Kit进行重亚硫酸氢盐处理。纯化后的DNA用于构建测序文库，在Illumina Novaseq 6000平台上测序。原始测序数据去除含接头、未知碱基或低质量碱基的reads后，使用Bowtie2将clean reads唯一比对到毛竹参考基因组。利用BS-Seq数据检测单个胞嘧啶的甲基化状态、甲基化胞嘧啶位点数量及每种基因组背景下的甲基化比例。甲基化水平计算为每个甲基化胞嘧啶（mC）的reads数除以该胞嘧啶的总reads数。差异甲基化区域（DMRs）定义为至少包含5个CG（或CHG/CHH）位点且甲基化水平有2倍变化的区域（Fisher精确检验，P ≤ 0.05）。所有测序数据已上传至NCBI，登录号为PRJNA1371697。
定量实时PCR（qRT-PCR）分析使用RNA Easy Fast Kit从冷冻叶片中提取总RNA，PrimeScript? RT Reagent Kit合成第一链cDNA。qRT-PCR反应使用SYBR? Green Premix Pro Taq HS qPCR Kit在ABI StepOne系统上进行。基因特异性引物使用Oligo 7软件设计，以液泡膜内在蛋白41（TIP41）作为内参基因，采用2?ΔΔCt法计算各基因的相对表达水平。
PheGRF10（Ped09CXg23240）的过表达通过从cDNA中扩增PheGRF6a的全长开放阅读框，克隆至pMD18-T载体，再通过同源重组方法插入到由CaMV 35S启动子驱动的Cambia1301::PhePheGRF6a载体中。重组载体通过农杆菌（Agrobacterium tumefaciens，EHA105）介导的转化法导入水稻（Oryza sativa L. subsp. japonica cv. Zhonghua 11）。经过广泛的表型鉴定、PCR确认和卡那霉素筛选，成功建立了15个独立的T3代转基因株系。收集转基因株系和野生型植株茎秆从基部到顶部的第一、第二、第三和第四节间，选取各节间中部进行石蜡切片，方法同前。
结果
解剖结构比较显示，在节间发育早期，毛竹节间基部的细胞核出现率显著高于圣音竹。发育后期，毛竹进入细胞伸长期，其细胞长度显著长于圣音竹。茎秆成熟时，毛竹节间细胞长度是圣音竹的两倍以上。这表明圣音竹节间缩短是细胞分裂和伸长异常共同作用的结果。
转录组测序鉴定出导致圣音竹矮化的关键异常表达基因。基于基因表达水平的主成分分析（PCA）表明，代表毛竹节间发育早期的M1与发育晚期的M2以及圣音竹的两个节间样品（S1和S2）之间存在显著差异。差异表达基因（DEGs）统计显示，M1与S1之间共有13,754个基因差异表达，其中7,518个在毛竹中高表达，6,236个在圣音竹中高表达。在节间发育晚期，M2与S2之间共有7,037个基因差异表达，包括3,238个在毛竹中高表达和3,799个在圣音竹中高表达的基因。DEGs的功能注释（GO）表明，与核糖体生物发生、响应生长素、rRNA代谢过程和DNA甲基化相关的基因在M1与S1之间差异表达，而与响应刺激、次级代谢过程和响应非生物刺激相关的基因在M2与S2之间差异表达。此外，响应生长素、细胞壁组织或生物发生以及响应激素的基因也表现出差异表达。GO注释结果表明，细胞分裂相关基因（如与核糖体生物发生相关的基因）和DNA甲基化相关基因在毛竹笋生长早期发挥关键作用，而与内源激素和细胞壁合成相关的基因在晚期扮演重要角色。
WGBS分析表明，转录组结果证明大量DNA甲基化相关功能基因在毛竹和圣音竹之间差异表达，提示两者笋芽DNA甲基化水平的差异可能是导致圣音竹矮化的关键因素。因此，比较了毛竹和圣音竹的全局DNA甲基化水平。配对比较（M1对S1，M2对S2）的DMRs以毛竹中的甲基化增益为主。然而，毛竹在早期的甲基化增益显著高于晚期。DMRs数量统计显示，在M1和S1中，最常见的差异甲基化类型是CHG，其次是CG。相反，在M2和S2之间，最常见的差异甲基化类型是CG，其次是CHG。
CG甲基化在蛋白质编码基因及其±2 kb侧翼区域呈钟形分布。甲基化水平在转录起始位点（TSSs）和转录终止位点（TTSs）最低，TTSs周围的水平略高于TSSs。CHG甲基化的分布模式与CG甲基化相似，但在转录区域内波动较小。相比之下，CHH甲基化在转录区域的分布几乎均匀。然而，蛋白质编码基因的上游启动子区域具有比下游区域更高的CHH甲基化水平。此外，在转座子（TEs）主体内检测到高水平的DNA甲基化，而其侧翼区域的甲基化水平较低。另外，笋芽生长晚期（M2和S2）基因主体、TE主体及其上下游区域的DNA甲基化水平高于早期（M1和S1）。
将基因表达水平从低到高分为四个类别（Flank1-4）。总体而言，基因启动子和基因主体的甲基化水平与基因表达呈负相关，即基因表达水平越高，DNA甲基化水平越低。其中，启动子中的CHG和CHH甲基化水平对基因表达的影响大于CG甲基化，而基因编码区CG甲基化水平对基因表达的影响更为显著。DNA甲基化类型（CG, CHG, CHH）与不同基因组区域基因表达的相关性分析显示，大多数基因的基因间区、启动子和内含子区的甲基化水平与基因表达呈显著正相关或负相关，而外显子区的甲基化水平与基因表达无明显相关性。这些结果表明基因及其周围区域的甲基化水平显著促进或抑制基因表达。
基于差异表达基因的功能注释，选取了几个关键基因家族进行进一步的基因表达和启动子甲基化水平分析。这些家族包括参与DNA甲基化调控的（Domains Rearranged Methyltransferase, DRM；DNA Methyltransferase, DNMT）、细胞周期控制的（Cyclin；E2F transcription factor/Dimerization Partner, E2F/DP）、细胞壁合成与生长的（Expansin-like A, EXPA；Cellulose Synthase-Like D, CSLD；4-Coumarate: CoA Ligase, 4CL）以及一般细胞分裂与生长的（GRF）。DRM、DNMT、cyclin和GRF家族的大多数基因在节间生长早期于毛竹中的表达水平显著高于圣音竹，而在后期未见显著差异。这些家族中大多数基因的启动子甲基化水平与其相应的基因表达呈负相关。在GRF家族中，Ped09CXg23240（PheGRF10）在两个发育阶段于毛竹中的表达水平均显著高于圣音竹，且差异显著。此外，这种表达模式与其启动子甲基化水平呈负相关。E2F/DP在节间发育的早期和晚期均于毛竹中表现出更高的表达模式，且早期毛竹与圣音竹之间的表达差异更为显著。启动子甲基化水平分析显示，绝大多数基因的DNA甲基化水平在早期于圣音竹中更高，与基因表达呈负相关。相比之下，晚期基因启动子甲基化水平与基因表达无显著相关性。EXPA家族的大多数成员在节间发育晚期于毛竹中高表达。绝大多数成员在节间发育早期表现出启动子DNA甲基化水平与基因表达之间的显著负相关，而晚期甲基化与基因表达无显著相关性。例如，Ped01Dg09700和Ped04Cg23400在节间伸长晚期于毛竹中的表达均高于圣音竹。然而，Ped01Dg09700的表达与其启动子甲基化水平呈正相关，而Ped04Cg23400的表达则与其启动子甲基化水平呈负相关。此外，还分析了两个纤维素合成相关基因家族（CSLD和XTH）。这两个家族的基因在节间发育晚期于毛竹中的表达高于圣音竹，并且超过一半的CSLD和XTH成员的基因表达与其启动子的DNA甲基化水平呈负相关。
PheGRF10（Ped09CXg23240）的功能验证基于转录组数据，筛选出可能参与毛竹节间细胞伸长的关键基因PheGRF6a（Ped09CXg23240）。经过潮霉素筛选和半定量PCR分析，共获得15个转基因水稻株系。在水稻品种‘中花11’中过表达PheGRF10显示，与野生型（WT）相比，过表达株系主茎的株高显著增加，具有显著的统计学差异。此外，过表达株系的穗、种子和叶片均显著长于WT。qRT-PCR结果显示，PheGRF10在转基因水稻中高表达，而在WT中未检测到表达。而且，PheGRF10的高表达抑制了转基因水稻中水稻KNOX家族基因（KN1和KN2）的表达水平。转基因茎秆解剖结构观察显示，过表达株系的节间细胞长度长于WT。
讨论
DNA甲基化是调控植物生长发育的关键表观遗传修饰，通过调节基因表达发挥作用。WGBS揭示了毛竹和圣音竹之间不同的DNA甲基化模式，为竹子矮化的表观遗传调控提供了直接证据。本研究中观察到的CG、CHG和CHH甲基化在蛋白质编码基因和转座子（TEs）中的分布模式与被子植物中保守的甲基化谱一致。CG和CHG甲基化在基因主体呈钟形分布，在转录起始位点（TSSs）和转录终止位点（TTSs）处最低——这种模式与转录调控相关。相比之下，CHH甲基化在启动子区域富集，这可能与抑制转座子活性和维持基因组稳定性有关。在毛竹发育晚期，基因主体、转座子（TE）主体及其侧翼区域的DNA甲基化水平较高。然而，圣音竹早期和晚期之间的差异不如毛竹显著。这表明DNA甲基化在节间发育过程中受到动态调控，并可能参与从细胞分裂到细胞伸长的转变。
DNA甲基化模式的改变可导致染色质结构和基因可及性的变化，最终影响基因表达。基因的CG甲基化水平与其表达水平呈正相关，而许多高表达基因表现出低CG甲基化水平。在本研究中，我们得到了类似的结果：与CHG和CHH相比，基因主体内CG甲基化水平对基因表达的抑制作用更为显著。这表明CG水平在毛竹笋节间伸长过程中的基因表达调控中扮演重要角色。
研究表明，竹笋的节间伸长是细胞分裂和细胞伸长共同作用的结果。我们分别比较了毛竹和圣音竹节间发育的早期和晚期。我们发现毛竹在节间发育中表现出更明显的阶段特异性特征：其生长在早期以细胞分裂为主，在晚期以细胞伸长为主。不同功能的基因家族在毛竹节间伸长的早期和晚期履行不同的职责。在水稻中，当OsGRF1的表达被RNA干扰（RNAi）敲低时，会导致叶片变小、植株矮化和抽穗期延迟。在二倍体杨树中过表达GRF5会导致转基因株系叶片更大、生长更快；此外，转基因株系的茎粗和株高均增加。在本研究中，发现大多数GRF基因在毛竹节间伸长早期高表达，与甲基化水平呈负相关；此外，在水稻中异源过表达PheGRF10显著促进了茎秆高度生长，共同表明这些基因在DNA甲基化调控下，增强了毛竹的节间伸长。另外，与细胞增殖和细胞周期相关的E2F/DP基因家族成员在节间发育的早期和晚期均于毛竹中表现出比圣音竹更高的表达模式，且在早期观察到的差异更为显著。这表明诸如E2F/DP、GRF和cyclin等与细胞分裂和增殖相关的基因在毛竹发育早期的节间细胞分裂中发挥重要作用。相反，这些基因在圣音竹中的异常表达导致了节间细胞分裂异常。
在节间伸长晚期，与细胞伸长相关的EXPA基因家族在毛竹中的表达水平高于圣音竹。此外，总体来看，圣音竹早期和晚期之间的EXPA表达水平无显著差异，且其早期表达水平甚至高于毛竹。这与解剖观察结果一致：在节间伸长早期，毛竹的节间细胞长度实际上短于圣音竹。类似地，参与细胞壁合成和纤维素合成的XTH和CSLD家族也显示出相同的表达模式。这表明在节间伸长晚期，在DNA甲基化调控下，与细胞壁合成和细胞生长相关的基因在毛竹中高表达，从而促进了节间生长从以细胞分裂（早期）为主转变为以细胞伸长（晚期）为主。
细胞分裂相关基因（Cyclin, GRF）和细胞伸长相关基因（EXPA）在毛竹中的较高表达伴随着较低的启动子甲基化，而在圣音竹中情况则相反。所有这些表明DNA甲基化在调控竹笋节间伸长中起着重要作用。在毛竹中，参与节间细胞分裂的基因和参与细胞伸长的基因分别在节间发育的早期和晚期履行其功能。然而，在圣音竹中，基因表达显得更为紊乱。许多与细胞伸长和细胞壁合成相关的基因在早期高表达。这种不规则的表达间接影响了圣音竹的节间伸长生长——在细胞分裂生长完成之前就启动了细胞伸长生长。
结论
本研究通过整合解剖学、转录组学和全基因组重亚硫酸氢盐测序（WGBS）分析，探讨了圣音竹（毛竹变种）矮化的分子机制。解剖学观察揭示，与毛竹相比，圣音竹节间缩短源于细胞分裂和伸长异常，毛竹表现出明显的阶段特异性生长（早期以细胞分裂为主，晚期以细胞伸长为主）。WGBS分析鉴定出两种竹子之间不同的DNA甲基化模式：CG、CHG和CHH甲基化在基因主体和转座子（TEs）中呈现保守分布，启动子和基因主体甲基化与基因表达呈负相关。关键生长相关基因家族（GRF、E2F/DP、Cyclin、EXPA、CSLD、XTH）在毛竹中有序地时序表达，并受DNA甲基化调控，而在圣音竹中其表达紊乱（例如，细胞伸长相关基因早期高表达）。功能验证证实过表达PheGRF
