https://peerj.com/articles/13902/#

摘要

作物雄性不育具有重要的理论研究和育种应用价值。HTS-1 的雄蕊转化为雌蕊或雌蕊状结构，是春季三雌蕊 （CSTP） 小麦中重要的雄性不育材料。然而，HTS-1 中雌蕊发育的分子机制仍然是一个谜。11 个小麦组织的 RNA-seq 数据来自美国国家生物技术信息中心 （NCBI），包括 CSTP 的雄蕊和 HTS-1 的雌蕊和雌蕊。Salmon程序用于量化 11 种小麦组织的基因表达水平;基因定量结果通过每百万转录本 （TPM） 进行标准化。共表达网络分析 （WGCNA） R 包总共使用了 58,576 个基因来构建块级网络。我们获得了与 11 种小麦组织显著相关的所有模块。使用 AgriGO V2.0 进行基因本体论 （GO） 富集分析;从 GO 富集结果中鉴定出这些关于小麦雌蒲发育的重要模块中的基因和转录因子 （TFs）。使用基本的局部比对搜索工具 （BLAST） 将 HTS-1 蛋白与已发表的雌蕊相关蛋白和 TFs 进行比对。通过 qRT-PCR 分析和验证有关小麦雌蕊发育的基因。MEturquoise、MEsaddlebrown、MEplum、MEcoral1、MElightsteelblue1 和 MEdarkslateblue 模块与雌蕊显著相关 （相关性 p < 0.05）。此外，通过基因本体论 （GO） 分析，仅在 MEturquoise 模块中鉴定出 206 个与心皮发育 （GO 0048467：0048440） 或股蓝发育相关的基因，206 个基因中有 42 个是 MEturquoise 模块中的枢纽基因。qRT-PCR 结果显示，42 个 hub 基因中有 38 个在雌蕊和雌蕊中高表达于雄蕊。共有 15 种 pistillody 发育相关蛋白通过 BLAST 验证。还在 MEturquoise 模块中分析了转录因子 （TFs），鉴定出 618 个 TFs。总共有 11 个科的 56 个 TFs 被认为可调节雌蕊的发育。共表达网络显示，在 42 个 hub 基因中鉴定出 6 个 HB 基因和 3 个 BES1 基因。这表明 TFs 在小麦雌蕊发育中起重要作用。此外，还有 11 个乙烯相关基因与 TFs 或 hub 基因相连，表明乙烯相关基因在雌蕊发育中的重要作用。这些结果为雌蕊发育的分子相互作用提供了重要的见解。

介绍

普通小麦 （Triticum aestivum L.） 是人类蛋白质、维生素和矿物质的主要来源（Ramírez-González 等人，2018 年）。尽管近年来全球小麦产量稳步增长，但一些研究表明，小麦产量几乎无法满足人类日益增长的需求（国际小麦基因组测序联盟，2014 年;Marcos-Barbero et al.， 2021）。因此，许多研究人员努力提高谷物的生产力和品质。小麦花的发育会影响谷物的产量和品质（Yang et al.， 2015;Ullah et al.， 2022）。在过去的十年中，已经报道了各种小麦花突变材料（Murai & Tsunewaki， 1993;Murai 等人，2002 年;Murai，2013 年）。突变体（cr）-N26和（cr）-CSdt7BS是雄性不育材料，其雄蕊转化为雌蕊或雌蕊状结构，这种现象称为雌蕊状（Murai & Tsunewaki，1993;Murai 等人，2002 年;Murai，2013 年）。（cr）-N26 和 （cr）-CSdt7BS 突变体包含山羊杂草的外源基因;因此，（cr）-N26 和 （cr）-CSdt7BS 的雌蕊是由异源细胞质雄性不育引起的 （Murai等人，2002）。2013 年，Peng 等人从中国春三雌蕊 （CSTP） 小麦中发现了一种名为 HTS-1 的新雌蕊突变体（Peng，2003 年;Peng et al.， 2013）。CSTP 是常见小麦品种中国春的近等基因系，其 Pis1 基因来源于三雌蕊突变体（Peng et al.， 2013）。HTS-1 的遗传基础与 （cr）-N26 和 （cr）-CSdt7BS 的遗传基础不同 （Murai， 2013;Peng et al.， 2013）。HTS-1 突变体不含外源细胞质;因此，HTS-1 的雌蕊受核基因控制 （Peng et al.， 2013）。以前的报告表明，小麦雄蕊转化为雌蕊的分子机制非常复杂（Murai等人，2002 年;Zhu et al.， 2008;Yamamoto et al.， 2013;Yang et al.， 2018）。Murai 等人 （2002） 发现小麦APETALA3同源物 （WAP3） 表达降低可能与 （cr）-CSdt7BS 中雌蕊的诱导有关。Aegilops crassa 细胞质的线粒体基因 orf25 被认为是与异质小麦中雌腈相关的细胞质因子的候选者（Hama等人，2004 年）。Zhu et al. （2008） 报道了一个称为 orf260 的线粒体基因CRA，这可能导致 （cr）-CSdt7BS 的雄蕊中产生雌蕊。原位杂交分析表明，钙调蛋白结合蛋白基因小麦钙调蛋白结合蛋白 1 （WCBP1） 在发育早期至晚期在雌蕊样雄蕊中高度表达（Yamamoto et al.， 2013）。这一发现表明 WCBP1 基因在雌蕊样雄蕊发育中起重要作用（Yamamoto et al.， 2013）。HTS-1 小花中雄蕊向雌蕊的同源转化是由两个隐性核基因 hts1 和 hts2 的相互作用决定的（Peng等人，2013 年）。Yang et al. （2015） 使用 RNA-seq 分析确定了 206 个与 HTS-1 雄蕊和雌蕊发育相关的差异表达基因 （DEG）。其中，123 个基因下调，83 个基因上调。但这项研究从未发现过关于雌蕊发育的基因。Sun et al. （2019） 报道了表皮模式因子样 （TaEPFL1） 基因在雌蕊和雄蕊中的表达高于雌蕊和雄蕊。TaEPFL1 基因在拟南芥中的异源表达导致细丝和花梗缩短，这可能诱导雄蕊转变为雌蕊或雌蕊状结构。此外，MADS-box 的一些转录因子 （TFs）（例如，酵母中的 MCM1、拟南芥中的 AG、抗裂芥中的 DEFICIENCYCIENS 和哺乳动物中的 SRF）和 YABBY 家族的下垂叶 （DL） 基因调节小麦雌蕊的发育（Murai，2013;Hama et al.， 2004;Meguro et al.， 2003;Ishikawa 等人，2009 年;Zhang et al.， 2015;Wang et al.， 2015;Zhu et al.， 2015）。表达模式分析表明，TaDL1、TaDL2 和 TaDL3 可能导致 HTS-1 的雄蕊完全或部分转化为雌蕊（Zhang等人，2015）。这些研究报告了小麦雌蒲的发展。然而，雌蕊的机制仍不清楚。

随着高通量测序和计算机分析技术的快速发展，可以使用各种高通量测序数据系统地研究复杂的生物学问题。各种网络分析方法，如神经网络分析（Greenside等人，2018年）、蛋白质-蛋白质相互作用网络分析（Stelzl等人，2005年）和WGCNA（Zhang & Horvath，2005年），已经构建并广泛应用于使用高通量测序数据揭示分子机制。WGCNA是研究高度相关基因的相关网络的有效方法（Langfelder & Horvath，2008）。根据不同组织中的基因表达水平，具有相似表达模式的基因可以通过WGCNA聚集到同一个模块中（Zhang & Horvath，2005）。因此，WGCNA 被广泛用于探索各种植物中花发育的枢纽基因、模块-性状关系和基因相互作用网络（Hollender等人，2014 年;Roberts & Roalson， 2017;马、丁和李，2017 年）。特别是，Ramírez-González 等人（2018 年）使用来自 22 种组织类型的 209 个 RNA-seq 样本，分析了小麦 Azhurnaya 品种的基因表达模式;并建立了组织和胁迫特异性共表达网络，揭示了整个小麦发育过程中同源物表达的广泛协调。这些网络，包括小穗发育网络和详细的基因表达图谱，被存放在小麦电子荧光象形文字 （eFP） 浏览器 （http://bar.utoronto.ca） 中（Ramírez-González等人，2018 年）。但是，wheat eFP 浏览器中的网络不会揭示和提供有关雌螨发育基因机制的信息。Azhurnaya 栽培品种不是雌螨小麦。此外，WGCNA 对与雌蕊发育相关的基因的系统鉴定尚未报道。

本研究使用 WGCNA 分析了来自小麦组织的 11 个 RNA-seq 数据，包括 CSTP 的正常雄蕊和 HTS-1 的雌蕊和雌蕊，以探索与小麦雌蕊发育相关的基因及其相互作用机制 （Yang et al.， 2015）。探讨了模块与小麦 RNA-seq 数据中特定花组织之间的关系。然后，通过 GO 分析从与 pistillody 显着相关的模块中鉴定出与 pistillody 发育相关的基因。采用实时定量 PCR （qRT-PCR） 验证 HTS-1 中雌蕊发育相关的候选基因的表达模式，以验证 WGCNA 结果。这项工作为小麦雌蕊发育的分子机制提供了进一步的见解。

材料和方法

植物材料和 RNA 提取

CSTP 及其雌蕊突变体 HTS-1 于 2020 年 10 月的最后一周种植在中国南充市西师范大学的实验田中。HTS-1 的雌蕊和雌蕊于 2021 年 2 月中旬在开花初期采集，同时收获 CSTP 的正常雄蕊。从 30 株 HTS-1 和 CSTP 植株中收集 3 种花材料。每种材料都有两个生物学重复。特别是，仅收集了中央穗状花序的第一朵小花。立即用液氮将花材冷冻并保存在超低温冰箱中（Forma 900 系列;Thermo Scientific，马萨诸塞州沃尔瑟姆）在 −80 °C 下。 使用植物 Easy Spin RNA 小量制备试剂盒 （Biomiga， San Diego， CA） 提取 HTS-1 的雌蕊和雌蕊以及 CSTP 雄蕊的总 RNA。在进行互补 DNA （cDNA） 合成之前进行 DNase 处理。通过凝胶电泳测定 cDNA 的质量。计算 cDNA 浓度，并使用分光光度计 （Thermo， NanoDrop 2000） 调节至 500 ng/μL。

数据源

以下普通小麦的 RNA-seq 数据来自美国国家生物技术信息中心 （NCBI） 序列读取档案数据库 （SRA， http://www.ncbi.nlm.nih.gov/sra）：根 （ERR424732）、茎 （ERR424762）、旗叶 （SRR3068439）、孕穗 （SRR6802610）、穗穗 （SRR6802611）、轴 （SRR6802608）、早期卵巢 （SRR6802613）、晚期卵巢 （SRR6802612）、HTS-1 雌蕊 （SRR1175868）、HTS-1 雌蕊 （SRR1177760） 和 CSTP 雄蕊 （SRR1177761）。CSTP 是中国春小麦和三雌蕊小麦的回交后代，每个 CSTP 小花都有三个雌蕊和三个雄蕊（Peng et al.， 2013）。HTS-1 是来自 CSTP 的雌蕊突变体 （Peng et al.， 2013）。一些 HTS-1 雄蕊转化为雌蕊或雌蕊状结构（Peng等人，2013 年）。因此，HTS-1 雌蕊 （SRR1175868） 和 CSTP 雄蕊 （SRR1177761） 的 RNA-seq 数据为野生型，HTS-1 雌蕊 （SRR1177760） 为小麦花突变体。所有 RNA-seq 数据均由 Illumina HiSeq 系列平台生成，布局为成对末端读长。小麦编码序列 （CDS） 和蛋白质序列从 Ensembl Plants 数据库 （版本 41， http://ftp.ensemblgenomes.org/pub/plants/release-41/） 下载。共获得 133,346 个 CDS 和蛋白质序列。

基因定量分析

Salmon （Patro et al.， 2017） 用于量化上述 11 种小麦组织的基因表达水平。首先，采用Salmon指数程序使用默认参数创建小麦指数，并使用小麦的 133,346 个 CDS 作为输入数据。接下来，使用基于准映射的模式，通过 Salmon 定量程序针对小麦指数对 11 个 RNA-seq 数据进行定量。基因定量结果通过每百万转录本 （TPM） 进行标准化。过滤掉 TPM 值小于 0.5 的雌蕊雄蕊中的基因，以减少低表达基因的干扰。总共保留了 58,576 个基因。

WGCNA

WGCNA R包（Langfelder & Horvath，2008）被用来构建11个小麦组织（58,576个基因）的共表达网络。设置一个软阈值 （power） 以使网络中的基因符合无标度拓扑分布。第一个软阈值（功率）用于使用 0.9 的无标度拓扑拟合指数构建 WGCNA（图 1A）。然后，应用块级 Modules 函数构建块级网络，参数如下：power = 21， networkType = “unsigned”， maxBlockSize = 10000， TOMType = “unsigned”， minModuleSize = 30， mergeCutHeight = 0.15。图 1B 描述了分层集群树和合并的模块。如果一个模块与其他模块高度相关，则应将它们合并为一个模块。模块特征基因 （MEs） 被定义为模块中的第一个主成分，代表该模块中的基因表达谱。应用模块 eigengenes 函数计算 MEs。然后，使用模型矩阵函数为每个性状构建一个二进制矩阵，包括 11 个向量，每个向量有 11 个元素。在这个矩阵中，当位置等于 1 时，它表示所陈述的性状;否则，位置为 0。之后，通过模块-性状关系分析分析 MEs 与性状之间的相关性。模块-特征关系的 P 值是使用 corPvalueStudent 函数估计的。然后，绘制模块-性状关系热图（图 2）。当模块的 P 值小于 0.05 时，该模块被视为与特征显著相关。之后，选择与孕穗、早期子房、晚期子房、穗尖、雌蕊、雄蕊和雌蕊雄蕊显著相关的模块作为候选模块进行进一步分析。

候选模块基因的 GO 富集分析

AgriGO V2.0 （http://systemsbiology.cau.edu.cn/agriGOv2/index.php） （Tian et al.， 2017） 用于 GO 富集分析。选择奇异富集分析 （SEA） 工具和 Triticum aestivm 物种进行此分析。查询 ID 是来自候选模块的小麦的集成 ID，这些模块与 7 个小麦花组织显著相关。AgriGO 的其他参数设置如下：Fisher’s test作为统计检验方法，Yekutieli （依赖下的 faslse 发现率 （FDR）），显着性水平为 0.05，一个最小映射条目数，以及完整的 GO 类型。表 S1 显示了上述候选模块中小麦花发育的重要 GO 项，FDR 值小于 0.05。然后，应用 ggplot2 构建不同模块中花发育 GO 项的气泡图 （图 3）。

枢纽基因鉴定

Hub 基因是模块中具有高连通性 （前 10%） 的基因。使用基于特征基因的连接 （kME） 和基因显着性 （GS） 的绝对值确定枢纽基因。kME 测量从基因表达谱到模块特征基因的距离，WGCNA 中的 signedKME 函数用于计算 kMEs。GS 是每个基因与性状的关联。使用以下阈值确定来自与 pistillody 显着相关的模块的枢纽基因： GS >0.2 和 kME >0.8 的绝对值。然后，将枢纽基因的信息与花发育相关基因整合，并显示在表 S1 中。

六倍体小麦蛋白与已发表的雌雄蕊蛋白之间的相似性分析

表 S2 列出了与雌蕊发育相关的已发表蛋白质，包括蛋白质种质和来源参考信息。采用 BLAST 的 BLASTP 程序研究已发表的雌蕊雄蕊蛋白与获得的小麦蛋白 （-b 1 -m 8 -F F） 之间的相似性。共 172 个已发表的雌蕊蛋白被用作查询序列，六倍体小麦的蛋白质被用作主题序列（或数据库）。对齐结果如表 S2 所示。表 S2 中的标识定义为比对序列数的匹配碱基数。标识值越接近 100%，查询和主题序列之间的匹配碱基对就越多。使用以下方法计算查询序列的覆盖率，以评估查询序列和主题序列之间的相似性。首先，从起始位置 （query start） 中减去对齐结束位置 （query end）。然后，差值相加 1。此外，查询序列的覆盖率值是通过将查询序列的长度除以上一步的总和来获得的。在表 S2 中，使用 AgriGO V2.0 对主题序列进行 GO 富集，但仅显示与小麦花发育相关的 GO 术语。主题序列用于进一步分析，其相似性和覆盖率值与查询序列相比高达 95%。最后，如果主题序列满足以下标准，则被视为雌蕊发育的候选基因：（i） GO 分析结果与小麦心皮发育有关 （GO：0048440），或 （ii） 它们是重要模块中的枢纽基因。

TF 分析

使用与 pistillody 显着相关的模块的基因 ID 提取小麦 Ensembl 蛋白的序列，以鉴定 pistillody 发育的 TFs。接下来，采用 HMMER V3.0 软件包的 hmmsearch 程序使用来自重要模块的蛋白质对 HMM 谱预测 TF。总共使用了属于 69 个 TF 家族的 75 个 HMM 图谱来识别 TF （Chen等人，2015 年）。使用了 hmmsearch 的默认参数。当蛋白质的 e 值小于 0.01 时，这些蛋白质被认为是 TFs。HTS-1 的预测 TF 列于表 S3 中。将 MEturquoise 模块的枢纽基因 ID 与预测的 TFs 整合。然后，使用上述相同参数应用 AgriGO V2.0 进行 GO 富集分析。鉴定出的 TFs 的蛋白质 ID 用作输入数据。GO 分析结果如表 S3 所示。调查了已发表的与雌蕊发育相关的 TFs（表 S4）。此外，利用 BLASTP 程序将已发布的 TFs 与 HTS-1 （-b 1 -m 8 -F F） 预测的 TFs 进行比较，以找出它们的相似性。BLAST 结果显示在表 S4 中。如果 TFs 满足以下标准，则认为它们是调节雌蕊发育的候选基因：（i） GO 分析结果与心皮发育 （GO：0048440） 或股蓝发育 （GO：0048467） 有关;（ii） BLAST 比对的相似度和覆盖度值小于 95%。

小麦雌蕊发育共表达网络的可视化

根据 GO 富集结果，从表 S1 中获得了心皮发育 （GO：0048440） 和股蓝发育 （GO：0048467） 的中心基因。然后，从雌蕊相关的显著相关模块中鉴定出与这些枢纽基因相连的花发育目标节点。Cytoscape V3.3 用于可视化小麦心皮发育的网络（Shannon et al.， 2003）。网络仅显示相应拓扑重叠大于 0.08 的连接。导入数据后，我们选择了 Tools → NetworkAnalyzer → Network Analysis → Analyze Network。接下来，从统计数据菜单中选择常规样式，以设置将节点大小映射到度数，并将边缘大小映射到权重。图 4 显示了与心皮发育 （GO：0048440） 和股蓝发育 （GO：0048467） 相关的枢纽基因网络。