植物耐盐性的全基因组关联及差异表达分析陆地棉L在萌发阶段

背景

盐胁迫是严重影响作物产量的主要非生物胁迫。开发和利用耐盐品种是解决土壤盐渍化最经济的途径。陆地棉是盐碱地上的主要纤维作物和先锋植物，其耐盐遗传结构有待深入研究。

结果

本研究采用全基因组关联分析和RNA测序技术，对196个陆地棉基因型萌发期耐盐品质性状位点(qtl)和候选基因进行了检测。利用四种环境下的耐盐性综合评价值，鉴定出33个显著单核苷酸多态性(snp)，其中至少在两种和四种环境下分别有17个和7个显著snp。17个稳定snp位于13个qtl的98个候选基因内或附近，其中35个基因被功能注释为参与盐胁迫反应。RNA-seq分析表明，98个候选基因中有13个稳定差异表达。此外，13个候选基因中有12个通过qRT-PCR进行了验证。RNA-seq分析在三个采样时间点检测到6640、3878和6462个差异表达基因，其中869个是共享的。

结论

这些结果，包括具有准确耐盐性评价的棉花优质材料、显著SNP标记、候选基因和耐盐途径，有助于我们进一步了解棉花耐盐胁迫的分子调控机制和遗传改良。

盐度是一种重要的非生物胁迫，在世界范围内降低作物产量和质量[1]。全球8亿农业用地中有6%以上受到盐碱化的影响[2]。全面了解盐响应分子机制和探索耐盐基因将有助于提高作物的耐盐性[1，3.，4]。陆地棉(Gossypium是天然纤维、植物油和蛋白质的重要来源，也是一种中等耐盐的先锋作物，可在盐碱地种植。然而，随着土壤盐碱化程度的增加，其产量将急剧下降[5]。棉花种质间的耐盐性差异显著。因此，筛选优质高耐盐种质是选育耐盐棉花的关键，也是鉴定耐盐分子机制和关键基因的关键。

目前已有几种评价种质资源的方法。主成分分析的因子分析(Factor analysis of principal component analysis, PCA)是一种常用的通过分析大量样本和主要相关指标来评价一组中每种材料状态的方法[6，7，8，9]。隶属函数分析(SFA)常用于评价应力耐受性[8，10]。然而，当仅使用PCA或SFA时，评估可能是片面的[10，11]。结合PCA和SFA的综合评价值可以将各指标转化为独立的因子，在保持原始信息的前提下相互比较，从而对植物耐受性进行更全面的评价[10]。该综合评价值比直观估算干枯面积比例的等级评价具有更高的准确性和效率。这种综合方法已被用于评估甘蔗的抗逆性[10，12]，黄瓜[13]，番茄[11]，苜蓿[8]，小麦[14]和棉花[15]。

耐盐性是一种多基因控制性状，受环境因素影响较大。基于连锁不平衡(LD)的关联图谱是识别与表型特征特定变异相关的基因组区域的一种常用且强大的技术[16基于它比连锁图谱分析更多的等位基因的能力。然而，利用关联定位鉴定棉花耐盐qtl /基因的研究[7，16，17]甚至联动映射[18，19是有限的。此外，利用高密度单核苷酸多态性(SNP)标记对棉花耐盐性进行关联定位的研究结果也有报道。共有23个显著snp和280个可能的候选基因被发现与两种耐盐相关性状相关，其中大多数与转录因子、转运体和酶有关[20.]。随着二代测序技术的成熟和普及，RNA测序(RNA-seq)已成为挖掘候选基因、构建分子调控通路和潜在调控网络的主要手段。一些对盐有反应的mrna [21，22]， miRNAs [3.]，备选剪接[23]，或长链非编码rna (lncRNAs) [24利用RNA-seq技术在棉花中检测到]，初步探索部分基因潜在的分子调控途径或调控网络。在棉花中，通过关联或连锁定位和RNA-seq检测到几个盐胁迫诱导基因，包括GhNHX1［25]、金属硫蛋白(GhMT3a) [26]，GhERF2-GhERF6［27，28]，GhDREB1［29]、ccch型锌指(GhZFP1) [30.]，GhNAC1-GhNAC631［31]，GhMPK2［32]，GhMKK1［33]，GhSOD1和GhCAT1［34]，GhWRKY17［35),而GhAnn1［36]。

在本研究中，我们结合综合评价、关联定位和RNA-seq来探索棉花萌发期耐盐候选染色体区域/基因。本研究为剖析棉花耐盐遗传机制和品种选育提供了候选QTL区域和基因。

耐盐性能及评价

在萌发期，对196个陆地棉基因型的10个性状进行了测定，结果见表1。结果表明，盐胁迫对10个性状均有显著影响。盐胁迫下各性状的均值和极值均显著低于正常条件。相比之下，GP(发芽势)是受影响最大的性状，而SFM(芽鲜质量)受影响最小。0.3% (200 mmol/L) NaCl胁迫不仅降低了总发芽率，而且降低了萌发速度和随后的生长。在正常条件下，根系干质量(RDM)的cv(变异系数)为7.1% ~活力指数(VI)的20.6%。在0.3% (200 mmol/L) NaCl胁迫下，除VI和RFM外，其余10个性状的cv均显著升高。茎干质量(SDM)、GP和萌发率(GR)分别排在前三位。这些结果表明，这组陆地棉对盐胁迫的大部分性状响应是高度多样化的。

为了减少环境因素的影响，达到较高的评价精度，根据每个性状的9个(3年× 3个重复)耐盐指数计算出BLUPed (BLUP: the best linear unbiased prediction) STI(耐盐指数)，结果见附加文件1:表S1 [37，38]。3年间性传播感染及各性状的BLUPed性传播感染无显著差异。然而，BLUPed性传播感染及其cv在各性状之间存在显著差异。GR、RL(根长)和PH(株高)表现出较大的BLUPed sti和较小的cv，而VI表现出最小的BLUPed sti和较大的cv。GP和SFM的BLUPed sti较小，cv较大。高耐盐性性状(高BLUPed STI)的CV值通常较低，反之亦然，但SDM的BLUPed STI值次之，CV值中等。为了进一步揭示盐胁迫下所测性状之间的关系，我们估计了Pearson相关系数，并将其列在附加文件中2表S2。除RFM-SDM和RFM-SFM呈极弱的不显著相关外，其余性状对均呈显著正相关。

耐盐性综合评价

为了综合评价该小组的耐盐性，对sti进行了PCA。在主成分分析前，采用KMO (Kaiser-Mayer-Olkin)测度和Bartlett检验检验sti是否适合进行因子分析。4个数据集的KMO值(2014:0.825;2015: 0.856;2016: 0.852;BLUP: 0.857)均> 0.5,Bartlett’s检验值均表现良好(χ²>χ²_{0.01, 45}， Sig = 0.000 < 0.05) [6]，表明STI数据适合进行PCA分析。四个STI数据集在PCA中得到一致的结果。对于每组STI数据，根据特征值-大于1.0规则(附加文件)，获得了至少占总方差86.85%的两个主成分3.:图S1和附加文件4表S3)。如PCA图所示(附加文件)5:图S2)，因子1代表除STI_RFM(根鲜质量耐盐指数)为因子2外的所有耐盐性状。

进而根据两个主因子(附加文件)计算的从属函数值(U)，得到各基因型的综合评价值(D)6:表S4)。通过K-mean聚类分析，将196个基因型根据代表基因型2014、2015、2016年和BLUPed数据的4组D值分为4组。196个基因型中，高耐盐27个，D值在0.669 ~ 0.934之间;中耐盐35个，D值在0.486 ~ 0.633之间;盐敏感67个，D值在0.320 ~ 0.474之间;高耐盐67个，D值在0.053 ~ 0.312之间(另附文件)6:表S4)。

SNP标记与耐盐性的关系

在GWAS(全基因组关联研究)中，综合评价值(D)作为表型数据来检测显著的耐盐qtl /基因。为此，在混合线性模型(MLM)中引入协变Q矩阵，减少群结构的负面影响，以减少误报。为了检验结构Q矩阵和亲缘关系矩阵(MLM [Q + K])的混合线性模型和关联结果的质量，我们采用了QQ (Quantile Quantile)图。数字1结果表明，MLM [Q + K]模型在检测显著SNP标记时略过于严格，说明假阳性的概率要低得多。因此，共有33个snp显著(p< 0.001)与耐盐性相关，使用MLM [Q + K]检测(图2)。2和表2）.其中，至少在两种环境(3年和BLUPed数据)中检测到17个snp，在所有四种环境中均检测到7个snp，解释了5.98 ~ 10.76%的表型变异，平均为8.21%。因此，这17个snp被用来鉴定候选qtl。当第一个SNP与相邻SNP之间的距离小于LD在成对r处的衰减距离时²=每条染色体上0.1级(附加文件7表5)或r²第一个SNP与相邻SNP之间> 0.1时，这些SNP标记作为候选QTL的置信区间。最终在10条染色体上获得13个qtl，其中4个qtl位于同源染色体A07-D07上。

GWAS检测候选基因

根据候选qtl的物理位置参考g .分子(TM-1)基因组[39]，共检测到98个候选基因(表2)2)，所有这些基因都被注释了进去拟南芥(附加文件8表6)。在UniProtKB中检索候选基因以注释其生物学过程(附加文件)8表6)。功能注释显示35个候选基因与耐盐性相关(表2)3.）.其中，11个基因涉及8个转录因子家族，包括CO-like、MYB、bZIP、ERF、TALE、SBP、HD-ZIP和ARR-B;17个基因与酸性pH、高温、过氧化氢、缺水和高渗盐胁迫等“应激反应”或“防御反应”相关;9个基因与水杨酸(SA)、赤霉素酸(GA)、茉莉酸(JA)、钙介导的信号传导、细胞分裂素和脱落酸(ABA)等“信号传导”或“对信号因子的反应”相关;6个基因参与“离子稳态”或“离子运输”，1个基因与“脂肪酸生物合成过程”相关，5个基因与氨基酸合成或运输相关。此外，这35个候选基因中有12个涉及两到三个功能类别(表2)3.）.

通过GO (Gene Ontology)富集分析进一步推断这98个候选基因的功能。在一个P-value < 0.05且基因数> 3时，将98个候选基因分为8个GO项(图2)。3.）.生物过程(biological process, BP)和分子功能(molecular function, MF)这两个主要类别都包含四个GO术语。本研究中BP类包括“几丁质分解代谢过程”、“碳水化合物代谢过程”、“代谢过程与转录调控”，MF类包括“几丁质酶活性”、“水解酶活性”、“转录因子活性与转运体活性”。

根据物理位置，GWAS中检测到的98个候选基因中有21个接近显著snp，其中6个(Gh_A07G0622-TM18816、Gh_A07G0729-TM19028、Gh_A12G0874-TM41811、Gh_D07G0251-TM63245、Gh_D10G0640-TM73567和Gh_D10G0647-TM73579)在编码序列(CDS)中具有显著snp(另附文件)9表S7)。此外，根据它们的功能注释或生物过程(表1)3.)， 6个基因中有3个(Gh_A07G0622、Gh_A07G0729和Gh_D07G0251)与“应激或防御反应”、“信号或对信号因子的反应”和“转录因子”的功能类别相关，暗示它们参与了耐盐反应。

转录组测序

在0.3% NaCl胁迫和正常条件下，分别于处理后3、24和72 h收获萌发期的Han682根组织进行RNA分离。从H3_1、H3_2、S3_1、S3_2、H24_1、H24_2、S24_1、S24_2、H72_1、H72_2、S72_1和S72_2的12个RNA测序文库中，共获得50.5、45.9、46.6、52.6、46.8、49.4、45.1、44.5、43.6、55.8、4530和4060万个原始reads(附加文件)11表8)。至少88%的clean reads(没有低质量reads和适配器序列的原始reads)被映射到g .分子(TM-1)基因组，其中唯一映射的reads约占干净reads的81%。共有6640、3878和6462个deg(差异表达基因)(P-值< 0.05)，其中3956(59.6%)、2238(57.7%)和3226(49.9%)表达上调(图2)。4）.在三个采样时间点中，我们确定了869个共有的deg，其中562个基因持续上调，307个基因持续下调。

为了验证我们的DEG结果，我们对随机选择的20个DEG采用了qRT-PCR，其中包括13个持续上调的基因和7个下调的基因(另附文件)10:图S3)。使用ΔΔCt方法计算相对表达量。转录组测序结果与我们的qRT-PCR结果一致。

为了监测耐盐基因的表达，利用氧化石墨烯富集分析(P-值经FDR < 0.05校正)对562个连续上调基因和307个连续下调基因进行了分析(图2)。5）.对于持续上调的基因(图2)。5a) BP类包括“单生物代谢过程”、“氧化还原过程”、“碳水化合物代谢过程”、“光合作用过程”、“有机氮化合物分解代谢过程”、“低聚糖代谢过程”和“芳香族化合物分解代谢过程”;细胞组分(CC)类别主要与“类囊体”和“光系统”有关;MF类别包括“氧化还原酶活性”和“o -甲基转移酶活性”。对于持续下调的基因(图2)。5b) BP类主要涉及“氧化应激反应”、“废GTP分解代谢过程”和“多元醇代谢过程”;CC类主要涉及“病毒衣壳”、“细胞器膜整体成分”和“细胞器膜内在成分”;MF类主要涉及“转运蛋白活性”、“过氧化物酶活性”和“无机二磷酸酶活性”。869个共有deg的所有富集GO项都与应激反应相关。

对869个连续上调/下调的deg进行KEGG(京都基因与基因组百科全书)途径富集(图2)。6）.持续上调的基因(P-value < 0.05)在9条KEGG通路中富集(图2)。6A)，包括“光合作用”、“类黄酮生物合成”、“淀粉和蔗糖代谢”、“光合作用-天线蛋白”、“丙酮酸代谢”、“植物激素信号转导”、“β -丙氨酸代谢”、“半乳糖代谢”、“花生四烯酸代谢”。持续上调基因的KEGG通路对应于与应激反应相关的GO富集结果和基因功能。持续下调的基因主要与9种KEGG通路相关(图2)。6B)，包括吞噬体、苯丙氨酸代谢、苯丙氨酸生物合成、淀粉和蔗糖代谢、戊糖和葡萄糖酸盐相互转化、抗坏血酸和醛酸盐代谢、硫代葡萄糖苷生物合成、氧化磷酸化和丙氨酸、天冬氨酸和谷氨酸代谢。

TFs在植物逆境响应中起关键作用[40]。在本研究中，共有62个(7.13%)编码tf的基因被鉴定为共享deg，其中41个基因持续上调，21个基因持续下调(表1)4）.所有62个基因均属于19个TF家族，包括ARF、bHLH、bZIP、C2H2、ERF、HD-ZIP、MYB、MYB_related、NAC和WRKY等几个响应非生物和生物胁迫的关键调控基因家族。

结合关联分析和转录组测序

我们结合GWAS和RNA-Seq结果进一步筛选耐盐候选基因。在GWAS的98个候选基因中，有13个基因的表达水平显著不同(p< 0.005)，在多个采样时间点进行RNA-seq分析5）.在13个推测的基因中，8个基因仅在一个时间点表现出显著差异，4个基因(Gh_A03G1740、Gh_A12G0877、Gh_D07G0263和Gh_D07G0500)在两个时间点表现出显著差异，1个基因(Gh_A07G0622)在三个时间点都表达上调。此外，13个推测的DEGs中有6个(Gh_A01G1563、Gh_A02G1100、Gh_A07G0622、Gh_A12G0877、Gh_D07G0251和Gh_D07G0500)在GWAS中接近显著snp, 6个基因中有2个(Gh_A07G0622和Gh_D07G0251)在编码序列区域内存在显著snp(附加文件)9表S7)。此外，如上表所示，Gh_A07G0622和Gh_D07G0251参与“应激或防御反应”、“信号或对信号因子的反应”和“转录因子”(表1)3.）.

为了进一步验证推测的基因，除Gh_A12G0877无法与同源基因区分外，所有13个推测的deg在盐或水处理后3、12、24、48和72 h，通过qRT-PCR对每个样品进行三组重复分析。使用ΔΔCt方法计算相对表达量。数字7结果表明，盐处理后Gh_A01G1563、Gh_A07G0622、Gh_D07G0243和Gh_D07G0251等推定基因持续上调。Gh_D07G0623、Gh_D07G0250、Gh_D07G0258和Gh_D07G0500基因在盐胁迫初期(处理后3 h或12 h)对盐胁迫有响应。Gh_A07G0623基因在处理后72 h对盐有反应。候选基因Gh_A03G1740、Gh_A02G1100和Gh_D07G0249在盐胁迫下持续下调。qRT-PCR结果显示，这12个基因可能都参与了耐盐反应。

耐盐性综合评价

陆地棉耐盐性状复杂，因品种、发育阶段和组织而异[j]。41]。因此，采用精确的方法对种质耐盐性进行准确评价具有重要意义。棉花耐盐性评价指标包括种子萌发指标(GR、GP、GI和VI)，植物形态指标(PH、SFM、SDM和RL)，生理生化指标(Na⁺K⁺和甜菜碱)，以及产量性状(铃数、铃重和皮棉产量)[7，20.，42]。目前普遍认为多指标、多方法相结合的耐盐性综合评价更为真实可靠，如PCA、SFA、综合评价D值[6，7，8，9，10，11]。综合评价D值具有较高的准确性，已广泛应用于不同作物或蔬菜种质的抗逆性评价[j]。8，10，11，12，13，14，15]。本研究选取了10个与耐盐性相关的性状，并将其综合评价D值应用于GWAS中，有助于提高该品系的耐盐性评价，进而提高GWAS的耐盐性评价。

GWAS检测候选qtl /基因

耐盐性是棉花重要而复杂的性状。以前的一些研究利用不同的标记类型和定位群体对耐盐分子标记进行了研究[7，16，17，18，19，43，44]。此外，还报道了棉花耐盐性qtl的meta分析[45]。与以往报道相比，未发现与本文检测到的QTL共定位的QTL。这种低一致性可能与以下因素有关，包括棉花耐盐机制的复杂性，以及不同研究群体、耐盐评估性状、标记类型和标记密度的差异。本研究采用高密度SNP阵列和4组D值(3年和BLUPed D值)进行GWAS分析，在不少于2种环境中检测到13个控制萌发期耐盐性的候选qtl。此外，这些qtl区域包含候选基因，其功能被认为与盐耐受性有关，并受盐胁迫调节。

这13个qtl中的许多候选基因很可能与耐盐性有关拟南芥和功能注释(表1)3.和附加文件8表6)。7个TF家族共涉及11个基因，分别为MYB、bZIP、ERF、TALE、SBP、HD-ZIP和ARR-B。TFs在植物生物和非生物胁迫响应中发挥重要作用[46]。MYB [46，47，48]， bZIP [49，50，51]及ERF [28，52，53]也参与对生物和非生物应激的反应。co样基因(COL4)正调节非生物(盐和渗透)胁迫耐受性拟南芥通过依赖aba的方法[54]。SBP是参与植物发育的植物特异性TF，可能通过调节植物生长来增强植物的抗逆性[j]。55]。ARR-B基因ARR12(本研究Gh_D07G0251)负向调控拟南芥干旱[56]。HD-ZIP tf参与植物的发育和应激(干旱或盐)反应答:芥［57，58]，大米[59]，马尼奥特·埃斯库伦塔·克兰茨［60),Craterostigma plantagineum［61]。此外，HD-ZIP基因(GhHB1)g .分子对根系发育、脱落酸和盐的反应[62]，表明HD-ZIP基因在植物发育和逆境反应中发挥重要作用。如表所示3.在GWAS中检测到的这17个候选基因对酸性pH、高温、过氧化氢、缺水或高渗盐胁迫作出反应或参与“防御反应”，极有可能与耐盐性有关。此外，Gh_A01G1564响应GA和JA, Gh_A02G1099和Gh_A02G1100响应SA, Gh_D07G0251和Gh_D07G0256响应细胞分裂素(CK)， Gh_D07G0500响应ABA, Gh_D10G0643响应生长素(IAA)。GA、JA、SA、CK、ABA和IAA等植物激素在植物对外界刺激和环境变化的适应中发挥着重要作用[j]。3.，22，63]。这些参与“激素刺激反应”或“信号传导”的基因与以往盐胁迫下植物表达谱的研究结果一致[5，64，65]。Gh_A07G0729也参与钙介导的信号传导，在适应盐胁迫中发挥重要作用[3.，66，67]。此外，有6个基因参与“离子稳态”或“离子运输”，从而增强植物的耐盐性[68，69，70]。

此外，在该GWAS中富集候选基因的氧化石墨烯术语包括碳水化合物代谢过程、代谢过程、水解酶活性、转录调控、转录因子活性和转运蛋白活性等，这些术语也主要富集在前期棉花根系对盐胁迫的响应研究中[j]。3.，21，22，65，71]。转基因烟草中几丁质酶的过度表达增强了其对生物(真菌和细菌病原体)和非生物(盐和重金属)胁迫的耐受性[j]。72]，重金属胁迫诱导植物几丁质酶亚型的积累[73]。在该GWAS中，与盐胁迫相关的富集氧化石墨烯术语包括“几丁质分解代谢过程”和“几丁质酶活性”。

转录组测序中与耐盐性相关的候选基因和途径

迄今为止，对棉花耐盐机制的研究有限，大多数研究集中在盐反应基因和途径的转录组测序上。徐等人(2013)[65研究了将植物暴露在200 mM NaCl中3、12、72或144小时后根系基因表达的变化，发现富集的氧化石墨烯与细胞成分有关，包括“膜固有”、“细胞质囊泡”和“膜部分”。Peng et al. (2014) [3.在盐胁迫(200 mM NaCl)处理后4和24 h，棉花叶片中的氧化石墨烯含量得到了鉴定，并揭示了“转录因子活性”、“胁迫响应”和“生物过程调控”等氧化石墨烯富集的术语。郭等(2015)[5报道称，在150 mM NaCl处理下，耐盐和盐敏感品种的转录本均上调，富集了与“刺激响应”、“转录因子活性”、“过氧化物酶体”和“脯氨酸代谢过程”相关的氧化石墨烯术语。张等(2016)[22鉴定出了DEGsg . davidsonii在盐胁迫(200 mM NaCl)后12、24、48、96和144 h的根和叶中，发现了丰富盐响应型氧化石墨烯类别的deg，包括“对氧化应激的响应”、“对渗透胁迫的响应”、“离子运输”、“对各种激素刺激的响应”、“对蔗糖刺激的响应”和“代谢过程”。此外，与代谢过程相关的基因还涉及“碳水化合物”、“激素”、“蛋白质”、“脂质”、“氨基”、“氧化还原”和“有机物”。Shu et al. (2017) [21]报道了NaCl/CK的deg与氧化石墨烯的“氧化还原”、“低聚糖代谢过程”、“光合作用”、“类囊体”和“氧化还原酶活性”相关。

本研究在高耐盐旱地品种汉682中，分别在盐胁迫后3、24和72 h筛选到6640、3878和6462个基因。其中562个基因持续上调，307个基因持续下调。持续上调的DEGs富集的氧化石墨烯术语与“代谢过程”、“氧化还原过程”、“碳水化合物代谢过程”、“光合作用”、“低聚糖代谢过程”、“类囊体”和“氧化还原酶活性”有关，这在以往的研究中也有发现[5，21，22，74]。被持续下调的deg富集的氧化石墨烯术语包括“膜固有”、“细胞质囊泡”和“膜部分”、“氧化应激反应”和“转运蛋白活性”，与以往的研究结果一致[3.，5，22，65]。与“氧化还原”相关的氧化石墨烯术语充满了连续上调和下调的deg，这表明氧化还原系统引发了更复杂的应激反应。此外，与“膜”、“转运体活性”和“类囊体”相关的氧化石墨烯术语在根和叶中都富含显著的deg，这表明与耐盐性相关的一些机制可能在不同的植物组织和器官中共享[3.，5，21，22，65]。

在前人对棉花幼苗根系耐盐机理的研究中，也发现了与耐盐相关的途径，包括类黄酮生物合成、淀粉和蔗糖代谢、植物激素信号转导、淀粉和蔗糖代谢、苯丙素生物合成、苯丙氨酸代谢和吞噬体[qh]5，21，22，74，75]。这些途径可能在植物适应逆境中起重要作用。然而，丙酮酸代谢、半乳糖代谢和花生四烯酸代谢途径是苗期特有的。在幼苗期盐胁迫下的叶片中也发现了丙酮酸代谢途径[75]。盐条件下棉花种子的萌发需要适宜的盐离子稳态、充足的能量供应和去除有害物质的能力等生理状态。这些特定的途径，如丙酮酸代谢，可能能够提高盐胁迫下幼苗期种子的存活率，尽管其潜在的机制尚不清楚。黄酮类化合物生物合成途径中因子含量最高，具有极显著性(P-值< 0.0001)在持续上调的deg富集途径中。黄酮类化合物是一种多酚类次生代谢化合物，在植物生长发育、繁殖和应激防御等方面具有重要作用[j]。76，77]。Petrussa et al. (2013) [77研究表明，在严重或长期的胁迫条件下，黄酮类化合物在植物中构成了二级ros消除系统。Petrussa et al. (2013) [77]和Fini等人(2011)[76]提示液泡在ROS稳态中的关键作用可能是由类黄酮介导的。Ma等人(2014)[78研究了小麦中与类黄酮途径相关的基因表达和与耐旱性相关的类黄酮积累。黄酮类化合物还在生物应激的防御反应中发挥关键作用[79]。

类似的MYB- bhlh - wdr (MBW)复合物和MYB小蛋白家族(R3-MYB)已被证明在类黄酮生物合成的调控中发挥关键作用[47，80]。这证实了TFs在植物的生物和非生物胁迫响应中起着重要作用。在目前的研究中，有19个TF家族始终与盐胁迫反应有关[3.，5，22]，连续差异表达(表2)4）.其中，ARF、ERF、C2H2下调3个，C3H、CO-like、HSF、LBD、M-type_MADS、NF-YA、NF-YB上调7个。排名前六的tf分别是bZIP、MYB (MYB相关)、HD-ZIP、NAC、bHLH和HSF。

假定的盐反应基因

根据UniProtKB的生物学过程和功能，通过关联图谱和转录组测序检测到的13个推定deg中，大多数可能对盐胁迫有反应(表1)6和附加文件8表6)。其中，8个(Gh_A01G1563、Gh_A02G1100、Gh_A07G0622、Gh_A07G0623、Gh_A12G0877、Gh_D07G0249、Gh_D07G0251和Gh_D07G0500)参与了“应激或防御反应”、“信号或对信号因子的反应”或“转录因子”，很可能与盐胁迫反应有关(表2)3.和5）.这8个盐调控deg与与耐盐性相关的基因同源，如参与脂肪酸生物合成过程的Gh_A01G1563(CUT1)，编码壁相关受体激酶2的Gh_A02G1100 (WAK2) [81，82]，并响应ABA，编码热休克蛋白的Gh_A07G0622 (CLPB1)和Gh_A07G0623 (HSP18.2) [3.]， Gh_A12G0877 (ERF2)编码乙烯反应性转录因子[53，56， Gh_D07G0249 (CHIT1)编码壳三酸苷酶1 [5，73]， Gh_D07G0251 (RR23)参与细胞分裂素激活信号通路[3.，22，64]和Gh_D07G0500 (HVA22E)对ABA和应激的响应[83，84]。此外，Gh_A03G1740 (BGAL3)参与碳水化合物代谢过程，而碳水化合物代谢过程一直是盐胁迫下植物主要富集的氧化石墨烯项[qh]3.，5，22];Gh_D07G0250与甲基化有关，甲基化在棉花耐盐性中起重要作用[46，85];Gh_D07G0263 (GAPN)参与与盐胁迫相关的甘油醛-3-磷酸脱氢酶(NADP+)(非磷酸化)活性[86，87]。

qRT-PCR结果还显示，Gh_A01G1563、Gh_A07G0622、Gh_D07G0243、Gh_D07G0251、Gh_D07G0623、Gh_D07G0250、Gh_D07G0258、Gh_D07G0500、Gh_A07G0623、Gh_A03G1740、Gh_A02G1100和Gh_D07G0249基因的表达在萌发期受到盐胁迫的调控。这些结果为进一步研究棉花耐盐机理提供了候选基因。这些基因在棉花耐盐性中的具体功能和分子调控有待进一步研究[88]。

本研究利用10个盐相关性状的综合D值对196份材料的耐盐性进行了综合评价。在此基础上，进行了耐盐GWAS试验。在GWAS中，从至少两种环境中的17个显著snp的13个候选qtl中获得了98个候选基因。功能注释显示，98个候选基因中有35个参与耐盐反应。此外，在盐处理或对照处理后的三个时间点对高耐盐品种Han682进行转录组测序，以验证GWAS的结果。结合GWAS和RNA-seq结果，鉴定出13个推测基因，并利用qRT-PCR对其中12个基因的表达进行了验证。这些结果将加深我们对棉花耐盐性的分子遗传调控的认识，并有助于棉花耐盐性相关性状的改良。

植物材料和SNP标记

利用来自中国5个棉花产区的169个基因型和来自其他国家27个基因型的196个不同陆地棉花材料，进行了萌发期耐盐性关联图谱的研究。所有196份材料在本研究使用前均进行了至少3年的自交系杂交。关于这个面板的详细信息在附加文件中描述12:表S9 [89]。

从77774个SNP (CottonSNP80K， [90])进行种群结构分析和GWAS。将196个基因型划分为2个亚群，采用4种方法(UPGMA (unweighted pair group method with arithmetic means))进行系统发育、主成分分析(PCA)、结构分析(STRUCTURE)和亲缘关系矩阵(kinship matrix)进行鉴定13:图S4) [89]。snp、多态性信息含量(PIC)、基因多样性和LD衰减的汇总也按照Yuan等人(2018)的描述进行了计算7:表S5) [89]。

发芽期耐盐性评价

对盐胁迫最敏感的萌发期耐盐性评价[j]91]，于2017年进行了三重随机区组实验，种子来自三个日历年(2014-2016)。用0.1% HgCl对硫酸脱色后的棉籽表面消毒15 min，选择大小均匀的棉籽进行萌发试验，在28±1℃、80%±2%相对湿度的培养箱中进行萌发试验。种子(100粒)播种于含有750 g干砂的萌发箱(13 × 19 × 12 cm)中，种子上方均匀覆盖250 g干砂，然后加入250 mL 200 mmol/L NaCl溶液或蒸馏水作为对照。每个基因型的每个重复(处理或对照)由两个萌发箱中的200颗种子组成。

在播种后第3天至第10天记录种子发芽数。然后计算GR、GP、GI、VI，公式如下:GR =\(\frac{G7}{TS}\times 100\% \)和GP =\(\frac{G3}{TS}\times 100\% \),在那里TS每个重复的种子总数(200个)和G_3.或G₇为播种后第1天至第3天或第7天两个萌发箱中总发芽种子数;和GI =\(\sum \frac{Gt}{Dt} \)， vi = gi x年代,在那里t是种植后的天数，Gt发芽种子的数量在t播种后的第二天，Dt种植后的天数是否对应Gt,年代是一株幼苗的新鲜重量。

除了发芽特性外，还研究了其他几种耐盐性相关性状。第10天，每个萌发箱随机选取10株，测定PH、SFM、SDM、RL、RFM和RDM。每个重复两个萌发箱的平均值用于进一步的数据分析。

各性状STI的计算公式如下:\(\mathrm{STI}=\left(Y\mathrm{non}\ \mathrm{salt}\ \mathrm{stress}\times Y\mathrm{salt}\ \mathrm{stress}\right)/{\left(\overline{Y}\mathrm{non}\ \mathrm{salt}\ \mathrm{stress}\right)}^2 \),在那里Y_{non-salt压力}为无盐胁迫时的表型值，Y_盐胁迫盐胁迫下的表型值是多少\(\overline{Y} \)_{non-salt压力}无盐胁迫时的平均表型值为[92，93]。

BLUPed STI采用rlmer4包计算每个性状的9个STI(3个重复× 3年)。然后，在随后的计算中使用了四组sti(三年sti和BLUPed sti)。各主成分因子的权重计算为W我=国际扶轮/∑国际扶轮（我= 1, 2，…，n)，和R_我表示的贡献率我主成分因子。

从属函数分析如下:U(西) = (西−X最小值)/ (X马克斯−X分钟)(我= 1,2，…，n)，其中X_我的值是我主成分因子，和X_最小值和X_马克斯的最小值和最大值我主成分因子[10，11]。

耐盐性综合评价值(D)计算为D =∑(U)我×W我) (i = 1, 2，…，n) [10，11]。使用Excel (Office 365)和SPSS(23.0版本，RRID:SCR_002865)进行数据处理。

关联映射

采用TASSEL版本5.2.40 (RRID:SCR_012837)软件进行标记-性状关联，阈值设为p< 0.001表示显著的标记-性状关联。对于每条染色体，LD衰减距离作为其上检测到的候选QTL的置信区间。与g .分子AD1基因组NAU-NBI组装体v1.1 [39]，每个QTL获得候选基因。

RNA-seq和DEG

在对196个基因型的耐盐性进行综合评价的基础上，选择高耐盐品种Han682(即H1)进行盐胁迫下的RNA-seq分析。用0.1%盐酸表面消毒后，种子在萌发箱中湿沙中播种。当胚根长到2 ~ 3 mm时，选择胚根均匀的种子，在其他萌发箱中以0.3%(重量百分比，约200 mmol/L)或0% NaCl(对照)的沙床中种植。在种植后3、24和72 h收获根组织，并立即保存在液氮中。盐处理植株的根样编号为S3、S24和S72，对照植株的根样编号为H3、H24和H72。每次治疗重复2次。按照制造商的建议，使用NEBNext®Ultra™RNA Library Prep Kit for Illumina®(New England Biolabs, Ipswich, MA, USA)生成RNA-seq库。文库制备在中国北京Novogene生物信息学研究所Illumina HiSeq平台上测序。

参考基因组的索引(陆地棉AD1基因组NAU-NBI组装v1.1) [39]是使用Bowtie v2.2.3 [94]，使用TopHat v2.0.12将clean reads与参考基因组比对[95]。使用HTSeq v0.6.1以每百万定位片段每千碱基外显子片段来估计基因表达水平[96]。采用DESeq R软件包(1.18.0)对两种条件下的deg进行鉴定P-value < 0.05 [97]。

GO和KEGG富集

使用GOseq R软件包对DEGs进行GO富集分析[98]。与…妥协P经FDR校正的-值< 0.05视为显著富集。利用KOBAS软件检测KEGG通路中deg的统计学富集程度[99]。

中存在的验证

采用qRT-PCR验证deg的表达。使用OmniPlant RNA试剂盒(DNase I)提取总RNA。首先使用HiFiScript cDNA合成试剂盒对分离的RNA (1500 ng)进行反转录生成cDNA。然后使用UltraSYBR混合物在Applied Biosystems®7500 Real-Time PCR系统上作为qRT-PCR模板。在此过程中，所有试剂盒和混合物均购自北京科文生物科技有限公司。目的基因的特异性引物使用NCBI中的Primer-BLAST和多倍体基因组特异性引物设计网络平台GSP进行设计[One hundred.]。的GhUBQ7作为内参基因[101]。所有特定的引物(附加文件14表S10)由上海三工生物科技有限公司(中国上海)合成。

与本研究有关的大部分数据已包含在手稿的表格/图表中。作者很乐意应要求分享剩余的原始数据。

阿坝:

脱落酸

BLUP:

最佳线性无偏预测

英国石油公司:

生物过程

答:

蜂窝组件

cd:

编码序列

简历:

变差系数

度:

差异表达基因

遗传算法:

赤霉酸

GI:

发芽指数

走:

基因本体论

医生:

发芽势

格:

发芽率

GWAS:

全基因组关联研究

是:

茉莉酸

KEGG:

京都基因与基因组百科全书

KMO:

Kaiser-Mayer-Olkin

LD:

连锁不平衡

lncRNAs:

长链非编码rna

MF:

分子功能

传销:

混合线性模型

主成分分析:

主成分分析中的因子分析

PH值:

株高

图片:

多态信息含量

QQ:

分位数分位数

QTL:

质量性状位点

经济价值:

质量性状基因座

RDM:

根干质量

RFM:

根鲜质量

RL:

根的长度

RNA-seq:

RNA序列

山:

水杨酸

长效磺胺:

射干质量

国家林业局:

从属功能分析

SFM:

释放新鲜物质

单核苷酸多态性:

单核苷酸多态性

STI:

耐盐指数

TF:

转录因子

UPGMA:

带算术平均数的非加权对群法

第六:

vVigor指数

我们要感谢中国农业科学院棉花研究所种质资源库(中国河南省安阳市)提供的27份外来棉花材料。我们感谢LetPub (www.letpub.com)，在编写本手稿期间提供语言协助。

国家重点研发计划项目(2018YFD0100303)、中国转基因育种重大专项项目(2017ZX08005-004-006)、现代农业产业技术体系项目(SDAIT-03-03/05)、自然科学基金项目(ZR2017MC057)、山东省农业种子项目(棉花品种培育，2014-2017)资助。山东农业大学盐碱地改良创新项目(2015)。我们感谢所有基础的经济支持。资助机构为研究项目提供资金支持，但不参与研究的设计、数据收集、数据分析或手稿的撰写。

作者指出

1. 袁彦超和邢惠贤对这一工作的贡献相同。

从属关系

1. 山东农业大学作物生物学国家重点实验室/农学院，山东泰安

   袁彦超，邢慧娴，曾文冠，徐家玲，毛丽丽，王丽媛，冯伟，陶金才，王浩然，张海军，王庆康，宋贤良，孙学珍

2. 青岛农业大学生命科学学院，山东省高校植物生物技术重点实验室，青岛市长城路700号

   验钞元

3. 菏泽市农业科学院，菏泽

   这张

贡献

XLS和XZS设计了实验。手稿由YY和XLS撰写。HX、WZ、JX、LM、LW、WF、JT、HW、HZ、QW和GZ有助于收集表型数据。YY和HX分析结果。YY和HX完成了大部分的实验，对这项工作贡献相同。所有作者都阅读并批准了最终稿件。

相应的作者

对应到Xianliang歌或Xue-Zhen太阳。

伦理批准并同意参与

在我们的研究中使用的植物标本的收集符合机构和国家的指导方针。根据当地法律进行了实地研究。伦理批准不适用于本研究。

发表同意书

不适用。

相互竞争的利益

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

出版商的注意

施普林格·自然对已出版的地图和机构关系中的管辖权要求保持中立。

开放获取本文遵循知识共享署名4.0国际许可协议(http://creativecommons.org/licenses/by/4.0/)，它允许在任何媒体上不受限制地使用、分发和复制，前提是您要适当地注明原作者和来源，提供到知识共享许可协议的链接，并注明是否进行了修改。创作共用公共领域免责声明(http://creativecommons.org/publicdomain/zero/1.0/)适用于本文中提供的数据，除非另有说明。

转载及权限