文献解读|Dev Cell-植物单细胞多组学揭示玉米穗发育过程

文章题目：Single-cell RNA sequencing of developing maize ears facilitates functional analysis and trait candidate gene discovery

研究单位：上海交通大学医学院

发表时间：2021年发表杂志：Developmental Cell

影响因子：13.417

关键词：scRNA-seq，ATAC-seq，Chip-seq, GWAS

研究背景

玉米的玉米穗的发育是与玉米产量紧密相关的生物学过程，在玉米穗的发育过程中，分生组织的细胞分裂、分化与发育构建花序结构，包括由花序分生组织（IM）形成的小穗对分生组织（SPM）和由SPM分生组织的分枝形成的小穗分生组织（SMS）。这些过程涉及重要的转录调控因子的激活与调控。

对玉米突变体的研究揭示了一些关键基因在玉米穗发育与进化中的生物学功能与作用机制。但是，这些研究受限于遗传冗余和多效性。因此，了解玉米穗发育过程中特定细胞类型和区域的表达图谱，从而进一步解析控制发育的基因网络至关重要。

本篇文章运用单细胞转录组测序技术（scRNA-seq）研究了玉米穗，构建玉米穗发育过程中scRNA-seq图谱，并通过与其他组学技术联合分析，展示scRNA-seq如何预测遗传冗余、整合转录调控网络和识别与作物产量性状相关的候选基因来促进遗传研究。

技术路线

本研究采用sc-RNA技术研究发育中的玉米穗，通过已知或预测的marker基因进行细胞注释，并通过mRNA原位杂交和FACS RNA-seq数据集的比较验证了scRNA-seq结构，从而构建了一个发育中的玉米穗的scRNA-seq图谱。此外，作者通过与其他组学联合分析分别从三个方面去展示scRNA-seq的应用范围，首先，scRNA-seq数据的高细胞分辨率可以预测冗余；其次，scRNA-seq可以与ChIP-seq结合来预测直接调控的TF；最后，scRNA-seq的marker基因可以提供对作物产量很重要的性状的调节因子的信息。

研究结果

1. 玉米穗发育过程的scRNA-seq图谱构建

本研究选取的是5-10 mm 阶段的玉米穗进行scRNA-seq（图1A），获得来自3个样本的12,525个单细胞，并使用玉米V3基因组检测到表达28,899个基因，与同一组织的bulk RNA-seq检测到的数量相当。接着，作者通过MetaNeighbor来测试细胞类型在3个scRNA-seq数据集中是否一致，最终聚类得到12个meta-cluster（图1b）。

图1 分离玉米穗原生质体构建单细胞转录图谱

2. mRNA原位杂交和FACS RNA-seq验证scRNA-seq结果

为了确定每个meta-cluster的类型，作者整理了一个已知或未知的花序发育marker基因列表，已在玉米或拟南芥中研究过其表达模式。其中，74个marker基因在玉米中具有突变表型的功能特征。在本实验中，检测到其中73个基因的表达，而每个基因在12个meta-cluster中的一个或多个中表达量高（图1D-N）。

为了验证鉴定的meta-cluster的细胞类型，作者首先使用差异表达分析来识别至少一个重复中AUROC ≥ 0.7的marke基因，我们确定了813个候选marke基因。根据表达标记的细胞百分比进一步选择每个meta-cluster的顶部标记，并显示出高表达（图2A）。然后，根据这些marke基因对发育的作用进行了优先排序后，使用mRNA原位杂交对其进行验证（图2B-2M）。验证结果证实meta-cluster的细胞注释可靠。

然后，作者在ATAC-seq实验中，通过FACS分选的pZmYAB14-TagRFPt表达细胞来研究染色质可及性在耳朵中确定的侧向器官分化过程中如何变化。

图2 通过mRNA原位杂交和FACS RNA-seq验证scRNA-seq

3. scRNA-seq可以预测遗传冗余并有助于预测转录调控网络

为了解决基因冗余常常掩盖单基因敲除表型从而导致区分冗余和非冗余旁系同源物的挑战，作者采用如下方法进行研究：

（1）在对玉米突变体ramosa3 (ra3) 的研究中，作者确定了ra3增强子作为其旁系同源物ZmTREHALOSE PHOSPHATE PHOSPHATASE 4 (ZmTPP4)；

    （2）在12个玉米ZmTPP基因中，其中ZmTPP4和ZmTPP12在ra3突变体中上调，这两个基因是补偿冗余旁系同源物的常见预测因子；

（3）CRISPR敲除ZmTPP12后发现并不会影响穗的发育，也没有增强ra3突变体。因此，推测ZmTPP4充当冗余补偿的角色；

（4）接着，作者提出问题：为什么是ZmTPP4充当冗余补偿器而非ZmTPP12？所以，接下来就是查询这两个基因的共表达网络，在89个玉米bulk-seq数据集中，ra3及其两个旁系同源物具有相似的共表达分数。而在scRNA-seq数据中，ZmTPP4与ra3高度共表达，符合实际研究结果。

为了进一步验证上述研究结果的可靠性，作者接着进行了如下研究：

（1）鉴定了一个VASCULAR PLANT One-ZINC-FINGER（ZmVOZ）基因，它的同源基因调控拟南芥的开花时间。其中ZmVOZ4和ZmVOZ5在scRNA-seq数据中表现出高度的共表达相关性；

（2）作者确定了一组共表达的高置信度基因，它们与ZmVOZ4和ZmVOZ5显示出一致的共表达相关性；

（3）为了了解这些旁系同源物是否具有冗余作用，作者采用CRISPR-Cas9敲除所有ZmVOZ基因（包括未在scRNA-seq数据集中未检测到的ZmVOZ1和ZmVOZ2）。实验结果如预测一样，单个ZmVOZ突变体没有发生明显的表型变化，但是ZmVOZ1,2,4,5四重突变体的开花时间严重延迟（图3A）。

因此，基于上述两个研究结果，作者得出结论：玉米scRNA-seq数据可以预测遗传冗余。

由于转录因子（TF）的直接靶基因应该在相同的细胞类型中共表达，因此，作者想验证scRNA-seq能否有助于构建转录调控网络，作者使用scRNA-seq数据来计算KN1与其已知的直接靶基因的共表达关系。结果发现，与使用所有玉米基因的对照组相比，它显著高于预期（p < 0.01），支持作者的假设（图3B）。

因此，作者采用ChIP-seq技术得到两个数据集：ZmHOMEODOMAIN LEUCINE ZIPPER IV6 (ZmHDZIV6)，在表皮中唯一表达（图3C）。以及ZmMADS16 (ZmM16)，在特定的花序器官中表达（图3D）。每个TF在ChIP-seq的生物学重复中均有显著重叠。作者确定了ZMMADZIV6和ZmM16的907个高置信度峰，并对这两个基因的motif进行进一步的研究。在ZmHDZIV6结合峰的motif分析发现了类似同源结构域亮氨酸拉链IV家族的成员类似的序列（图3E），类似的，对ZmM16结合峰的motif分析发现类似的结果。基于上述发现，作者对ZmHDZIV6和ZmM16结合靶点预测的准确性有信心。

由于调控TF靶点通常需要比较ChIP-seq结合靶点和突变RNA-seq来预测。但考虑到ZmHDZIV6和ZmM16基因家族成员众多，因而RNA-seq并不适用。因此，作者采取基于与每个TF的scRNA-seq共表达来预测靶点，分别识别了ZmHDZIV6和ZmM16的79和55个候选调控靶点。在ZmHDZIV6的候选预测靶点之中，作者选定了ZmHDZIV家族的另外5个基因，其中一些在玉米SAM表皮中表达相似并可能形成转录级联反应来调节表皮分化。验证了候选预测靶点的表皮表达，包括ZmNOD26样膜内在蛋白1(ZmNIP1A、GRMZM2G041980)和ZmPRECURSOR ELCITOR PEPTIDE1(ZmPROPEP1、GRMZM5G899080)(图3G-3L)。

图3 scRNA-seq网络预测冗余并构建基因调控网络

4. scRNA-seq识别与玉米产量性状相关的基因

玉米穗形态与产量相关指标有关。为了验证scRNA-seq中鉴定的细胞类型和marker基因是否与玉米产量的候选调控因子一致。作者使用靶向GWAS方法，比较来自分生组织、侧向器官和维管组织的scRNA-seq的marker基因与281个玉米穗部与产量相关的形态性状的GWAS的结果集合（图4A-4D）。使用scRNA-seq的marker 基因2kb内的单核苷酸多态性（SNP），作者发现meta-cluster3的marker基因ZmYABBY9 (ZmYAB9)对于穗轴重（CW）有显著的作用（图4B），还发现两个与穗直径（ED）有显著的SNP（图4C 和4D），这与meta-cluster 9,10和11的marker 基因GRMZM2G361210有关。

图4  scRNA-seq的marker基因与玉米穗性状相关

研究结论

本研究通过scRNA-seq为理解玉米穗发育过程提供了有价值的见解，刻画的细胞图谱可以为发育遗传学研究和育种提供信息。开发的方法可以应用于其他物种的其他组织，并且随着植物相关scRNA-seq数据集的增加，在单细胞分辨率下进行跨物种或跨组织比较分析可以提供有关在进化过程中如何选择基因特征从而促进农业生产的信息。

参考文献：

Xu, X.S., Crow, M., Rice, B.R., et al. 2021. Single-cell RNA sequencing of developing maize ears facilitates functional analysis and trait candidate gene discovery. Dev. Cell 56, 557-568.