水稻种质资源稻瘟病抗性全基因组关联分析

周海平 张 帆 陈 凯 申聪聪 朱双兵 邱先进 徐建龙,,5,*

水稻种质资源稻瘟病抗性全基因组关联分析

周海平1张 帆2陈 凯3申聪聪3朱双兵3邱先进4,*徐建龙2,3,5,*

1温州市农业科学研究院 / 浙南作物育种重点实验室, 浙江温州 325000;2中国农业科学院作物科学研究所, 北京 100081;3岭南现代农业科学与技术广东省实验室深圳分中心 / 中国农业科学院农业基因组研究所, 广东深圳 518210;4长江大学农学院, 湖北荆州 434025;5海南省崖州湾种子实验室, 海南三亚 572024

稻瘟病是一种对全球水稻生产威胁极大的真菌性病害, 鉴定抗稻瘟病基因并将其导入到现有感病品种改良品种的抗性是控制这种病害的有效途径。本研究利用5个稻瘟病菌株鉴定了212份籼稻和235份粳稻种质资源的苗瘟抗性, 分别筛选到8个和12个抗全部5个菌株的籼稻和粳稻种质材料。采用全基因组关联分析在籼粳混合群体、籼稻和粳稻种质资源中共定位到43个影响水稻苗瘟抗性的QTL, 抗GD00-193、GD08-T19、GD17-CQ16、HB1708和HLJ13-856菌株的QTL分别为9、4、14、14和2个。其中, 12个抗病QTL仅在籼稻亚群中检测到, 7个仅粳稻亚群中检测到, 1个为籼粳2个亚群共同检测到, 说明籼稻抗稻瘟病总体好于粳稻, 而且稻瘟病抗性存在明显的籼粳分化。同时影响水稻对2个及2个以上菌株的抗性或在2个及2个以上群体中同时被定位到的QTL共计11个, 利用候选区间关联分析和单倍型分析鉴定到23个抗病候选基因, 不同抗病候选基因在籼、粳群体中的分布频率不同。研究结果为水稻品种稻瘟病抗性分子改良提供种质资源和有利基因信息及不同抗病基因的育种利用策略。

稻瘟病抗性; 种质资源; 全基因组关联分析; 数量性状位点; 有利等位基因

水稻是世界上最重要的粮食作物之一, 是全世界一半以上人口的主食。由稻瘟病菌()引起的稻瘟病每年都给水稻生产造成至少10%左右的产量损失, 严重时减产幅度高达40%～ 50%[1]。目前防治稻瘟病主要采用化学防治, 但该方法不仅成本高, 而且会造成环境污染。研究水稻抗稻瘟病的遗传基础, 挖掘影响水稻抗稻瘟病的重要有利基因, 利用分子育种技术选育抗稻瘟病的水稻品种是解决这一问题最经济有效的方法。

我国稻瘟病菌主要分为8个亚群和50余个生理小种, 不同生理小种的遗传差异很大[2-3], 而且稻瘟病菌的致病性变异很快, 从而导致抗病品种的稻瘟病抗性丧失频繁[4]。稻瘟病抗性分为水平抗性和垂直抗性, 垂直抗性由少数主效基因控制, 存在生理小种专化性, 通常表现为抗性强但很容易丧失。水平抗性由微效多基因控制, 无生理小种专化性[5], 表现为有一定程度的发病但对产量影响不大。因此, 聚合垂直抗性和水平抗性基因对提高水稻品种的稻瘟病抗性, 延缓抗性丧失更有优势[6]。随着遗传学和分子标记技术的发展, 遗传学家利用F2、重组自交系和导入系等遗传群体已定位了超过500个影响水稻稻瘟病抗性的QTL (https://www.gramene.org/), 其中37个基因已被成功克隆[7]。

虽然基于双亲遗传群体的QTL定位方法是挖掘水稻稻瘟病抗性基因的重要途径, 但这种方法也存在一定的局限性。例如, 构建遗传群体需要花费大量的时间, 且由于群体样本有限, 导致QTL定位的精度较低[8]。此外, 由于群体只携带双亲2种基因型, 因而只能比较每一QTL亲本2种等位基因的优劣, 无法挖掘该位点种质资源中存在的优异等位基因。近年来, 基于连锁不平衡的全基因组关联分析已被成功应用于各种复杂性状的遗传解析[9-11]。利用全基因组关联分析策略, Wang等[12]采用336份来源于中国的籼稻种质资源和80万的SNP基因型定位到了抗16个稻瘟病菌株的30个QTL; Kang等[13]利用420份水稻种质资源和44,000个SNP基因型定位到影响水稻对5个稻瘟病菌株抗性的97个QTL; Mgonja等[8]利用70份籼稻和92份温带粳稻种质资源以及44,000个SNP基因型定位到影响水稻对8个抗非洲稻瘟病菌株抗性的31个QTL; Raboin等[14]利用150份热带粳稻品种和190份籼稻品种分别定位到2个和1个影响水稻田间抗性的主效QTL; Mgonja等[15]利用190份非洲品种定位到影响6个稻瘟病病菌抗性的25个QTL。

本研究鉴定了来源于全球的212份籼稻和235份粳稻种质资源对我国5个稻瘟病菌株的苗期稻瘟病抗性, 并利用全基因组关联分析方法分别定位了籼粳混合群体(212份籼稻和235份粳稻合在一起的群体)、籼稻和粳稻群体中抗稻瘟病的QTL, 并对同时影响抗2个及2个以上菌株或在2个及2个以上群体中同时定位到的重要QTL进行单倍型分析, 鉴定优异等位基因。研究结果将为分子育种改良水稻稻瘟病抗性提供重要的抗病种质资源和有利等位基因。

1 材料与方法

1.1 试验材料

试验材料为从菲律宾国际水稻研究所种质库引进的212份籼稻和235份粳稻种质资源, 其来源于全球56个国家和地区。其中籼稻种质资源主要来源于中国(36份)、菲律宾(34份)、印度(20份)、斯里兰卡(17份)、孟加拉国(14份)、马达加斯加(9份)、越南(7份)和印度尼西亚(5份)等国家; 粳稻种质资源主要来源于印度尼西亚(22份)、科特迪瓦(20份)、巴西(19份)、马达加斯加(18份)、泰国(12份)、哥伦比亚(9份)、中国(8份)、美国(8份)、印度(7份)和不丹(7份)等国家。

1.2 稻瘟病抗性鉴定

选用5个强致病力稻瘟病小种代表菌株, 采用人工接种的方法鉴定种质资源的苗期稻瘟病抗性, 其中GD00-193 (中B13)和GD08-T19菌株(中B03)分别源于中国广东省海丰县和雷州市, GD17-CQ16菌株(中B15)来源于重庆涪陵市, HB1708菌株(中C13)来源于湖北恩施市, HLJ13-856菌株(中B05)来源于黑龙江856农场, 前4个稻瘟病菌株收集于籼稻品种, 最后1个菌株收集于粳稻品种。参照Park等[16]描述的方法培养水稻幼苗和稻瘟病菌株。水稻种子催芽后穴播于30 cm × 20 cm × 5 cm规格的搪瓷盆里, 每盆播24个材料, 每个材料播种量为8～10 粒种子。稻苗采取旱育栽培, 长至一叶一心期每盆施用0.5 g硫酸铵, 接种前共施3次。待稻苗长至3.5～4.0叶龄, 用配置好的孢子液进行人工喷雾接种, 接种菌液量为20 mL盆–1。接种后置于遮光密闭的培养箱中, 在25℃下保湿24 h, 之后置于温室, 在25～28℃下保湿培植至稻苗发病, 接种7 d后进行调查。菌株的小种鉴定及致病性测定均设2个重复。病级调查按国际水稻研究所稻瘟病圃苗瘟分级标准进行[17], 按0～9级记载(0级为高抗, 9级为高感), 其中0～3级为抗病类型、4～9级为感病类型。试验重复2次, 取2次的平均值作为该品种的稻瘟病病级。

1.3 表型数据分析

采用Microsoft Excel 2010进行种质资源表型数据进行整理, 利用Stataistca 5.5[18]进行统计描述、显著性分析和不同稻瘟病菌株抗性的相关性分析。

1.4 种质资源基因型鉴定与全基因组关联分析

每个品种随机挑选20个单株的叶片混合, 采用MATAB法提取DNA[19], 并将DNA浓度调整至100 ng μL–1。利用HDRA 700K SNP芯片[20]对籼粳2套种质资源进行基因型分析。获得全基因组700K的基因型数据后, 利用PLINK1.9软件[21]去掉最低等位基因频率均小于5%和缺失率大于20%的SNP, 最终籼粳混合群体、籼稻和粳稻种质资源分别鉴定出209,493、175,812和115,157个高质量SNP。利用SVS8.4.0软件包中EMMAX模块中的单位点混合线性模型进行SNP与表型之间的关联分析[11], 阈值设置为1.0×10–5。水稻的LD平均约为200 kb[20], 本研究中最高峰值SNP两侧各100 kb内的QTL即认定为相同的QTL。QTL命名遵循+菌株前2个字母代号+菌株前2个数字编号+所在染色体号。

1.5 重要QTL的候选基因鉴定

针对同时影响2个及2个以上菌株抗性或2个及2个以上群体同时定位到的重要QTL进行候选区间关联分析。以日本晴作为参考基因组, 在700K SNP基因型数据中挑选出每个QTL最显著的峰值SNP上下游各100 kb区间内所有基因上的SNP, 利用SVS8.4.0软件包中EMMAX模块中的单位点混合线性模型进行基于候选基因的关联分析[11], 阈值设置为1.0×10–5。随后筛选出位于基因启动子区、5¢UTR、3¢UTR、外显子中导错义突变和内含子中导致剪切模式发生改变的显著SNP进行相关性状的单倍型分析, 采用单因素方差分析和多重比较检验不同单倍型之间的差异显著性, 阈值设置为<0.01。

2 结果与分析

2.1 种质资源的稻瘟病抗性分析

447份种质资源对5个稻瘟病菌株的抗性表现见表1。混合群体、籼稻和粳稻种质资源对GD08-T19菌株的抗性最强, 抗级最低, 分别为2.87、2.87和2.83; 对GD17-CQ16菌株的抗性相对较差, 抗级最高, 分别为4.58、4.72和4.35。其中, 大部分籼稻品种对GD08-T19和HLJ13-856表现为抗病, 对其他3个稻瘟病菌株表现为感病; 大部分粳稻品种对GD08-T19表现为抗病, 而对其他4个稻瘟病菌株表现为感病。籼粳种质资源对GD08-T19、GD17-CQ16和HB1708的抗感反应相似, 抗级无显著差异, 而籼稻对GD00-193的抗性显著差于粳稻, 对HLJ13-856的抗性显著优于粳稻。显然, 籼稻或粳稻对不同稻瘟病菌株的抗性反应不同, 即使对同一稻瘟病小种, 籼、粳之间的抗性反应也不完全相同, 稻瘟病抗性可能存在一定的籼粳分化。此外, 籼粳种质资源对5个稻瘟病菌株的抗性均表现为连续变异, 表明种质资源对5个稻瘟病菌株的抗性受多基因控制。其中8份籼稻(BIRAIN 360、CHUA DAU、IR 32453-20-3- 2-2、KARAYAL、KHAO PON、MA WAINE OHN、ZI GAN NAN GU和SAHEL 108)和12份粳稻(CHAHORA 144、IRAT 177、LLANERO 501、OS 6、GBUAPU 1、SENG、TOS 5790、CURINCA、IRAT 212、CIRAD 403、CT 13582-15-5-M和VIETNAM 1)对全部5个稻瘟病菌株均表现为抗病。

3个群体对不同稻瘟病菌株的抗性之间的相关性见表2。在3套种质资源群体中, 品种对3个菌株GD00-193、GD17-CQ16和HB1708的抗性之间均表现为显著正相关, 说明品种抗某一个稻瘟病菌株的同时很可能抗另外2个菌株。此外, 籼稻种质资源对粳型小种HLJ13-856的抗性与对其余4个籼型小种的抗性均不相关, 粳稻群体对粳型小种HLJ13-856的抗性与对2个籼型小种(GD08-T19和GD17-CQ6)的抗性不相关, 对另2个籼型小种(GD00-193和HB1708)呈弱的显著正相关。

2.2 抗稻瘟病QTL定位

利用3套群体共定位到43个影响5个稻瘟病菌株的抗病QTL, 它们分布于水稻全部12条染色体上。

2.2.1 抗GD00-193和GD17-CQ16的QTL定位

定位到9个影响水稻对GD00-193抗性的QTL, 分布于1号、3号、4号、6号、8号、9号和11号染色体上(表3和附图1)。其中7个QTL (、、、、、和)只在混合群体中被定位到, 对表型的贡献率为1.72～17.94%;只在籼稻核心种质中被定位到, 对表型的贡献率为13.49%;在混合群体和粳稻中同时被定位到, 分别解释表型变异的15.86%和15.90%。

共鉴定到14个影响水稻对GD17-CQ16抗性的QTL (表3和附图1)。其中9个QTL (、、、、、、、和)只在混合群体中被定位到, 最高能解释9.70%的表型变异; 2个QTL (和)只在籼稻核心种质中被定位到, 平均解释12.70%的表型变异;在混合群体和粳稻种质资源中同时被定位到, 分别解释4.93%和4.23%的表型变异;在混合群体和籼稻种质资源中稳定表达, 对表型的贡献率最高为10.64%;在3个群体中同时被定位到, 对表型的贡献率最低为3.69%。

2.2.2 抗GD08-T19、HB1708和HLJ13-856的QTL定位 在6号、11号和12号染色体上共定位到4个影响水稻对GD08-T19抗性的QTL (表4和附图1)。其中3个QTL (、和)只在籼稻核心种质中被定位到, 对表型的平均贡献率为9.96%;只在混合群体中被定位到, 解释5.89%的表型变异。

表1 种质资源对5个稻瘟病菌株的苗期抗性表现

a籼粳种质资源间稻瘟病抗级差异。

aDifference in blast resistance betweenandsubpopulations.

表2 种质资源对不同稻瘟病菌株之间的相关性

第1、2和3行分别代表混合群体、籼稻和粳稻种质资源对5个稻瘟病菌株抗性之间的相关系数;*和**分别代表显著水平为< 0.05和< 0.01。

Data in the first, second, and third rows represent correlation coefficients inmix population,andsubpopulations, respectively.*and**indicate significant difference at< 0.05 and< 0.01, respectively.

共鉴定到14个影响水稻对HB1708抗性的QTL, 位于除2号、9号和12号染色体外的9条染色体上(表4和附图1)。其中7个QTL (、、、、、和)只在混合群体中被定位到, 对表型的贡献率为1.81%～10.87%;只在籼稻种质资源中被定位到, 能解释表型变异的17.50%;、和只在粳稻稻种质资源中被定位到, 对表型的贡献率最大为10.94%;和在混合群体和籼稻种质资源中同时被定位到, 对表型的平均贡献率为14.95%;在混合群体和粳稻种质资源中同时被定位到, 解释1.08%的表型变异。

共有2个QTL与水稻对HLJ13-856的抗性关联(表4和附图1)。其中籼稻核心种质只定位到1个主效QTL (), 对表型变异的贡献率为18.93%;只在混合群体中被定位到, 对表型的贡献率为2.02%。

稻瘟病抗性在不同群体间存在较大的差异, 对5个菌株的抗性在混合群体中单独定位到25个抗病QTL, 其在籼或粳群体中均未检测到, 说明这25个抗病QTL在籼或粳亚群体中的变异较小; 籼稻亚群中检测到12个抗病QTL, 明显多于粳稻亚群中检测到的7个, 说明籼稻抗稻瘟病总体上好于粳稻, 其中籼粳共同定位到的仅针对GD17-CQ16菌株的1个QTL (), 表明稻瘟病抗性存在明显的籼粳分化。

2.3 重要位点的有利等位基因鉴定

在43个影响水稻对5个稻瘟病菌株抗性的QTL中, 有11个QTL同时影响水稻对2个及2个以上菌株的抗性, 或在2个或2个以上群体中同时被定位到。利用候选区间关联分析和单倍型分析在这11个QTL共鉴定到23个候选基因。

2.3.1 单个候选基因的QTL的有利等位基因鉴定

1号染色体28.92～29.12 Mb区间()、6号染色体10.32～10.52 Mb区间()和22.60～ 22.80 Mb区间()、7号染色体4.73～4.93 Mb区间()、8号染色体17.40～17.60 Mb区间()和20.66～20.86 Mb区间()和11号染色体6.92～7.12 Mb区间()均只鉴定到1个候选基因(图1)。其中编码细胞色素P450, 分为C和T两种单倍型, 其中单倍型T对GD00-193和GD17-CQ16的抗性均显著优于单倍型C, 为有利等位基因(图1-A)。编码一个含有NBS-LRR结构域的蛋白, 分成A和C两种单倍型, 单倍型A在混合群体和籼稻种质资源中对HB1708的抗性均显著优于单倍型C, 为有利等位基因(图1-B)。编码一个表达蛋白, 分成A和T两种单倍型, 单倍型A对GD17- CQ16和HB1708的抗性均显著优于单倍型T, 为有利等位基因(图1-C)。编码IPP转移酶, 分成C和T两种单倍型, 单倍型T在混合群体中对GD17-CQ16和粳稻种质资源中对GD17- CQ16和HB1708的抗性均显著优于单倍型C, 为有利等位基因(图1-D)。编码一个反转录转座子, 分成单倍型A和单倍型G, 其中单倍型A在混合群体和粳稻种质资源中对GD00-193的抗性显著优于单倍型G, 为有利等位基因(图1-E)。编码14-3-3蛋白, 分成单倍型A和单倍型G, 其中单倍型A在混合群体中对GD17- CQ16和HB1708的抗性显著优于单倍型G, 为有利等位基因(图1-F)。编码受体激酶5的前体, 有A和G两种单倍型, 其中单倍型A在混合群体中对GD17-CQ16和HLJ13-856和在籼稻种质资源中对GD17-CQ16的抗性显著优于单倍型G, 为有利等位基因(图1-G)。

2.3.2 多个候选基因的QTL的有利等位基因鉴定

有4个位点鉴定到多个候选基因。在1号染色体1.47～1.67 Mb区间内()共有105个SNP位于24个候选基因上, 其中有3个基因上有SNP超过阈值(图2-A)。编码甘露糖-6-磷酸异构酶, 有A和G两种单倍型, 单倍型G对GD00-193、GD17-CQ16和HB1708的抗性均显著优于单倍型A, 为有利等位基因。编码核酸酶PA3, 共分为3种单倍型, 其中单倍型TG对GD17-CQ16的抗性显著优于GA, 为有利等位基因。编码一个表达蛋白, 分成了A和C两种单倍型, 单倍型A对GD17-CQ16的抗性显著优于C, 为有利等位基因(图2-A)。

在1号染色体33.05～33.25 Mb区间内()有143个SNP位于30个候选基因上。有4个基因(、、和)内有SNP超过阈值(图2-B),其中没有启动子区、5¢UTR、3¢UTR、外显子中导错义突变和内含子中导致剪切模式发生改变的显著SNP。编码一个SAND结构域蛋白质, 分成了3种单倍型, 在粳稻种质资源中单倍型GT对HB1708的抗性均显著优于单倍型AC和GC, 为有利等位基因。编码一个表达蛋白, 在混合群体中将其分为3种单倍型, 在粳稻种质资源中分为2种单倍型。单倍型AGA在混合群体对GD00-193和HB1708和粳稻种质资源中对HB1708的抗性显著优于单倍型CAC, 为有利等位基因。编码rp1, 分为A和G两种单倍型。其中单倍型A在混合群体对GD00-193和HB1708和粳稻种质资源中对HB1708的抗性显著优于单倍型G, 为有利等位基因(图2-B)。

图1 7个单候选基因位点的候选区间关联分析及单倍型分析

单倍型分析中不同小写字母代表性状在< 0.01水平差异显著。

Different lowercase letter represents significantly difference at< 0.01 in haplotype analysis.

图2 4个多候选基因位点的候选区间关联分析及单倍型分析

单倍型分析中不同小写字母代表性状在< 0.01水平差异显著。

Different lowercase letter represents significantly difference at< 0.01 in haplotype analysis.

在8号染色体6.36～6.56 Mb区间内()共有187个SNP位于27个候选基因上, 9个基因(、、、、、、、和)内有SNP超过阈值(图2-C)。其中和没有启动子区、5¢UTR、3¢UTR、外显子中导错义突变和内含子中导致剪切模式发生改变的显著SNP, 其余基因(、、、、、和)均编码反转录转座子。、、、、和均有2种单倍型, 其中单倍型A、A、T、A、T和T在混合群体和籼稻种质资源中对HB1708的抗性显著优于单倍型G、T、C、G、C和G, 为有利等位基因。分为AA和TC两种单倍型, 单倍型AA在混合群体和籼稻种质资源中对HB1708的抗性显著优于单倍型TC, 为有利等位基因(图2-C)。

利用104个位于21个候选基因上的SNP对9号染色体10.66～10.86 Mb区间()进行关联分析, 5个基因(、、、和)内有SNP超过阈值(图2-D), 其中和没有启动子区、5¢UTR、3¢UTR、外显子中导错义突变和内含子中导致剪切模式发生改变的显著SNP。编码一个表达蛋白,编码捕光叶绿素/结合蛋白基因,编码一个含短链脱氢还原酶家族结构域的氧化还原酶。这3个基因均有2种单倍型, 其中单倍型A、G和C对GD17-CQ16的抗性显著优于单倍型C、A和T, 为有利等位基因。此外,和的单倍型C和T对GD00-193的抗性显著优于单倍型C和T (图2-D)。

3 讨论

3.1 水稻对稻瘟病菌抗性分化与籼粳亚种间稻瘟病抗性差异

稻瘟病是稻瘟病菌()感染引起的一种真菌性病害, 我国稻瘟病菌主要分为8个亚群和50余个生理小种, 且不同生理小种的遗传差异较大[2]。李旭升等[17]利用3K种质资源中的1217份种质资源筛选到144份抗稻瘟病的种质资源。对其中30份种质资源进行苗期30个稻瘟病菌接种发现仅17份种质资源的抗性频率高于70%。杨楠等[22]对云南的46份地方品种进行了8个生理小种抗性鉴定, 结果表明同一品种对不同菌株的抗性差异很大。此外, 目前已经克隆的37个稻瘟病抗性基因中, 仅、和等少数基因具有广谱抗性, 其余基因均只能抗少数几种生理小种, 表明同一抗性基因对不同生理小种的抗性不同[7,23-24]。本研究利用212份籼稻和235份粳稻种质资源对5个稻瘟病菌株进行了抗性鉴定, 结果表明大部分籼稻品种对GD08-T19和HLJ13-856表现为抗病, 对其他3个稻瘟病菌株表现为感病; 大部分粳稻品种对GD08-T19表现为抗病, 而对其他4个稻瘟病菌株表现为感病。此外, 品种对GD00-193、GD17-CQ16和HB1708这3个菌株的抗性之间均表现为显著的正相关, 而对其他2个菌株的抗性均无显著相关性, 表明品种对不同稻瘟病菌株的抗性存在分化。关联分析定位到的43个稻瘟病抗性QTL中, 有9个QTL同时影响对2个或2个以上菌株的抗性, 其中7个位点均同时影响水稻对GD00-193、GD17-CQ16和HB1708中2个菌株的抗性, 1个位点同时影响水稻对这3个菌株的抗性, 1个位点同时影响水稻对GD17-CQ16和HLJ13-856的抗性, 表明水稻对不同稻瘟病菌株的抗性也存在遗传分化。因此, 在抗稻瘟病抗性育种中, 要尽量选择同时抗不同菌株的抗病QTL, 以提高品种对稻瘟病的抗谱。

另一方面, 籼粳亚种对同一菌株的抗性也表现出显著差异。李旭升等[17]在湖北恩施进行了1217份种质资源的抗稻瘟病鉴定, 并筛选到的144份不同时期抗稻瘟病的抗病材料中, 一半以上为籼稻, 温带粳稻仅占7.6%, 说明在籼稻对恩施环境下的籼型生理小种的抗性显著优于粳稻。杨楠等[22]对云南地方品种进行8个籼型菌株的稻瘟病抗性评价, 结果也表明籼稻的稻瘟病抗性整体高于粳稻。本研究所用的5个菌株中, 除HLJ13-856为粳型菌株外, 其余4中均为籼型菌株。籼稻对GD00-193的抗性显著差于粳稻, 对HLJ13-856的抗性显著优于粳稻。结果表明, 对同一稻瘟病菌株, 籼粳之间的抗性存在差异。籼稻群体对粳型小种HLJ13-856的抗性与对其余4个籼型菌株的抗性均不相关, 粳稻群体对粳型小种HLJ13-856的抗性与对2个籼型小种(GD00-193和HB1708)呈弱的显著正相关, 也说明对稻瘟病菌的抗性存在着一定的籼粳分化。此外, 在本研究定位的43个抗5个稻瘟病菌株中, 有3个QTL在混合群体与籼稻种质资源中同时被定位到, 3个QTL在混合群体与粳稻种质资源中同时被定位到,仅1个QTL在3个群体中均被定位到, 表明籼粳群体对稻瘟病抗性的遗传基础不同。因此, 在稻瘟病抗性育种中, 尽量选用在不同遗传背景下定位到的QTL进行稻瘟病抗性改良, 这些抗病QTL的抗性表达受遗传背景的影响可能较小。

3.2 本研究定位的抗稻瘟QTL与前人研究的比较

本研究利用全基因组关联分析共定位到43个影响水稻对5个稻瘟病菌株抗性的QTL, 部分QTL与前人克隆的稻瘟病抗性基因/QTL位于相同或相近位置。例如, 1号染色体的、、和分别与()和位于相同区间[25-28], 2号染色体的与位置相同[29], 5号染色体的与()处于相近区间[30], 6号染色体的、和分别与()、//位置相同[23,30], 位于7号染色体的和与qDI-7位置相同[31], 位于11号染色体的、和与和物理位置一致[32-34], 位于12号染色体的与()位于相同位置[30]。将来需要对这些相同或相近的QTL进行精细定位和克隆验证它们是否为同一QTL。

3.3 重要抗稻瘟位点的候选基因分析

本研究利用全基因组关联分析的方法共定位到43个影响5个菌株抗性的QTL, 其中11个QTL同时影响2个以上菌株的抗性或在2个以上群体中同时被定位到, 为重要的QTL。在这些QTL中, 有7个位点挖掘到了1个候选基因, 有3个位点都挖掘到了3个候选基因, 1个位点挖掘到了7个候选基因。其中, 1号染色体1.47～1.67 Mb区间内编码核酸酶PA3的是抗GD17- CQ16菌株的最有可能候选基因, 由于其在混合群体中的3种单倍型对GD00-193和HB1708的抗性均无显著差异, 推测不同于抗病QTL1和的候选基因, 也即、和虽然定位在同一个区间, 但不是同一抗病基因, 可能存在紧密连锁的关系。水稻中已发现核糖核酸酶III[35]与水稻稻瘟病抗性相关, 表明核酸酶参与到水稻对稻瘟病的抗性, 因此最可能是该位点的候选基因。1号染色体28.92～29.12 Mb区间的候选基因编码细胞色素P450, 其同源基因[36]编码细胞色素P450单加氧酶, 可催化色胺转变成血清素, 增强水稻对稻瘟病均侵染的抗性。因此,可能是该位点的候选基因。1号染色体33.05～33.25 Mb区间内,编码一个SAND结构域蛋白质, 提高其表达量将导致颖花发育异常。水稻另一个基因与稻瘟病抗性有关, 该基因突变后也导致颖花发育异常[37-38]。因此,是该位点最有可能的候选基因。6号染色体10.32～10.52 Mb区间的编码一个含有NBS-LRR结构域蛋白, 前人克隆的稻瘟病基因中大多数均为含有NBS-LRR的蛋白质[7], 因此可能是该位点的候选基因。8号染色体20.66～20.86 Mb区间的编码14-3-3蛋白, 为已克隆的稻瘟病抗性基因[39]。在9号染色体10.66～10.86 Mb区间内,编码一个含短链脱氢还原酶家族结构域的氧化还原酶。前人研究发现[40]、[41]和[42]等多个编码氧化还原酶类的基因具有稻瘟病抗性, 可以推测植物体内的氧化还原反应可以影响植物对稻瘟病的抗性, 因此是该位点最可能的候选基因。上述这些可能的候选基因还有待通过基因编辑敲除和转基因过表达进行功能验证, 这些候选基因的有利等位基因也需要通过分子标记辅助选择的方法导入到感病品种中, 明确其在水稻抗病育种中的利用价值。

3.4 对水稻抗稻瘟病分子育种的启示

稻瘟病对全球水稻生产造成了巨大威胁, 尤其在亚洲和非洲, 稻瘟病每年都导致巨大的产量损失[40]。选育抗稻瘟病的水稻品种能极大地帮助这些地区降低稻瘟病造成的损失, 保障粮食安全。然而, 抗稻瘟病新品种的选育也极大地依赖于抗稻瘟病种质资源的发掘和利用。本研究鉴定了447份籼稻和粳稻种质资源对5个稻瘟病菌株的抗性, 从中筛选出了8份和12份抗全部5个菌株的籼稻和粳稻种质资源。这些种质资源能作为抗稻瘟病的抗源直接应用于育种。

在本研究定位到的43个影响水稻抗稻瘟病QTL中, 其中籼、粳亚群体分别定位到12个和7个主效抗病QTL, 籼稻中抗瘟基因多于粳稻, 与水稻育种实践中观察到的籼稻对稻瘟病抗性总体强于粳稻是一致的。本研究从籼、粳种质资源中鉴定到一批效应较大(解释表型贡献率籼稻背景>15%、粳稻背景>10%)且仅抗1个菌株的抗病QTL, 如仅在籼稻资源中鉴定到的抗GD17-CQ16菌株的, 抗HLJ13-856菌株的, 抗HB1708菌株的和。这4个QTL中,的有利等位基因为A, 其在籼稻种质资源中的频率为78.33%, 但其在粳稻种质资源中的频率为2.34%;的有利等位基因为C, 其在籼稻种质资源中的频率为20.57%, 但其在粳稻种质资源中的频率为0.45%;的有利等位基因为T, 其在籼稻种质资源中的频率为31.75%, 但其在粳稻种质资源中的频率为0;的7个候选基因在籼稻种质资源中的频率分别为36.77%、38.24%、19.59%、76.19%、33.33%、48.21%和35.23%, 而它们在粳稻种质资源中的频率分别为1.82%、3.20%、2.26%、1.94%、2.25%、0.46%和1.81%。结果表明, 这些QTL有利等位基因在籼稻中的比例显著高于粳稻, 因此在水稻抗稻瘟病育种实践中, 可以利用籼稻资源作为抗源将这4个抗瘟基因导入粳稻来提高粳稻品种的抗稻瘟病水平。同样, 仅在粳稻资源中鉴定到抗HB1708菌株的和抗GD00-193菌株的, 其中的候选基因在粳稻种质资源中的有利等位基因频率为23.41%, 而其在籼稻种质资源中的频率仅为0.49%;的候选基因在粳稻种质资源中的有利等位基因频率为75.69%, 而其在籼稻种质资源中的频率仅为0.96%。这2个QTL的候选基因有利等位基因在粳稻中的频率显著高于籼稻, 也可以利用粳稻资源作为抗源将这2个抗瘟基因导入籼稻来提高籼稻品种的抗稻瘟病水平。此外, 在1号染色体1.47～1.67 Mb区间定位到抗HB1708、GD00-193和GD7-CQ16的3个抗病QTL (、和),是最可能的候选基因, 该候选基因在籼稻种质资源的有利等位基因频率分别为59.89%, 显著高于粳稻中的40.20%。在后续的籼稻、粳稻和籼粳交品种选育过程中, 可以针对这些大效应位点和同时影响水稻对多个稻瘟病菌株抗性的重要位点, 根据等位基因的差异开发KASP标记, 通过分子标记辅助选择提高水稻品种的稻瘟病抗性; 此外, 还可以将抗不同稻瘟病菌株的QTL进行聚合, 以达到拓宽水稻抗谱的目的。

4 结论

212份籼稻和235份粳稻种质资源苗期对5个稻瘟病菌株的抗性变异十分丰富, 且对不同菌株的抗性存在分化, 对同一菌株的抗性也表现出籼粳差异。利用全基因组关联分析共定位到43个影响水稻对5个稻瘟病菌株抗性的QTL, 其中11个位点同时影响水稻对2个以上菌株抗性或在2个以上群体中同时被定位到, 为广谱较广的抗稻瘟病QTL。挖掘到23个候选基因和23个优异等位基因。根据抗病QTL在籼粳群体中的不同分布频率, 可以针对性地选择抗病QTL用于籼粳品种的稻瘟病抗性改良, 研究结果为水稻抗稻瘟病育种提供了抗源和抗病基因信息及不同抗病基因利用的策略。

致谢 广东省农业科学院植物保护研究所朱小源团队帮助完成稻瘟病接种鉴定, 在此深表感谢！

附图 请见网络版: 1) 本刊网站http://zwxb. chinacrops.org/; 2) 中国知网http://www.cnki.net/; 3) 万方数据http://c.wanfangdata.com.cn/Periodical- zuowxb.aspx。

[1] 杨义强, 朱林峰, 李晓芳, 付杰, 黄道强, 邱先进, 周少川, 王重荣. 抗稻瘟病基因的基因特异性KASP标记开发与应用. 植物遗传资源学报, 2021, 22: 1314–1321. Yang Y Q, Zhu L F, Li X F, Fu J, Huang D Q, Qiu X J, Zhou S C, Wang C R. Development and application of KASP marker specific for rice blast resistancegene., 2021, 22: 1314–1321 (in Chinese with English abstract).

[2] 杨小林, 施仕胜, 张舒, 吕亮, 喻大召. 湖北省稻瘟病重发区病菌群体致病性分化的研究. 湖北农业科学, 2016, 55: 4169–4171. Yang X L, Shi S S, Zhang S, Lyu L, Yu D Z. Population pathotype ofin rice blast epidemic areas of Hubei province., 2016, 55: 4169–4171 (in Chinese with English abstract).

[3] 钟宝玉, 黄德超, 朱小源, 陈玉托, 邹寿发, 杨伟新, 赖新红. 近十年广东稻瘟病菌生理小种变化分析. 仲恺农业工程学院学报, 2018, 31: 24–29. Zhong B Y, Huang D C, Zhu X Y, Chen Y T, Zou S F, Yang W X, Lai X H. Analysis of physiological races ofin Fuangdong during recent decade., 2018, 31: 24–29 (in Chinese with English abstract).

[4] Khan M A I, Latif M A, Khalequzaman M, Tomita A, Ali M A, Fukuta Y. Genetic variation in resistance to blast (Cavara) in rice (L.) germplasms of Bangladesh., 2017, 67: 493–499.

[5] Ashikani S, Rafii M Y, Rahim H A, Latif M A. Genetic dissection of rice blast resistance by QTL mapping approach using an F3population., 2013, 40: 2503–2515.

[6] Zenbayashi K, Ashizawa T, Tani T, Koizumi S. Mapping of the QTL (quantitative trait locus) conferring prtial resistance to leaf blast in rice cultivar Chubu 32., 2002, 104: 547–552.

[7] 杨德卫, 王莫, 韩立波, 唐定中, 李生平. 水稻稻瘟病抗性基因的克隆、育种利用及稻瘟菌无毒基因研究进展, 植物学报, 2019, 54: 265–276. Yang D W, Wang M, Han L B, Tang D Z, Li S P. Progress of cloning and breeding application of blast resistance genes in rice and avirulence genes in blast fungi., 2019, 54: 265–276 (in Chinese with English abstract).

[8] Mgonja E M, Blimponya E G, Kang H X, Bellizzi M, Park C H, Li Y, Mabagala R, Sneller C, Correll J, Opiyo S, Talbot N J, Mitchell T, Wang G L. Genome-wide association mapping of rice resistance genes againstisolates from four African countries., 2016, 106: 1359–1365

[9] Agrama H A, Eizenga G C, Yan W. Association mapping of yield and its components in rice cultivars., 2007, 19: 341–356.

[10] Huang X H, Wei Z H, Sang T, Zhao Q, Feng Q, Zhao Y, Li C Y, Zhu C R, Lu T T, Zhang Z W, Li M, Fan D L, Guo Y L, Wang A H, Wang L, Deng L W, Li W J, L u Y Q, Weng Q J, Liu K Y, Huang T, Zhou T Y, Jing Y F, Li W, Lin Z, Buckler E, Qian Q, Zhang Q F, Li J Y, Han B. Genome-wide association studies of 14 agronomic traits in rice landraces., 2010, 42: 961–967.

[11] Qiu X J, Yang J, Zhang F, Niu Y N, Zhao X Q, Shen C C, Chen K, Teng S, Xu J L. Genetic dissection of rice appearance quality and cooked rice elongation by genome-wide association study., 2021, 9: 1470–1480.

[12] Wang C H, Yang Y L, Yuan X P, Xu Q, Feng Y, Yu H Y, Wang Y P, Wei X H. Genome-wide association study of blast resistance inrice., 2014, 14: 311.

[13] Kang H X, Wang Y, Peng S S, Zhang Y L, Xiao Y H, Wang D, Qu S H, Li Z Q, Yan S Y, Wang Z L, Liu W D, Ning Y S, Korniliev P, Leung H, Mezey J, McCouch S R, Wang G L. Dissection of genetic architecture of rice resistance to blast fungus., 2016, 17: 959–972.

[14] Raboin L M, Ballini E, Tharreau D, Ramanantsoanirina A, Frouin J, Courtois B, Ahmadi N. Association mapping of resistance to rice blast in upland field conditions., 2016, 9: 59.

[15] Mgonja E M, Park C H, Kang H X, Blimponya E G, Opiyo S, Bellizzi M, Mutiga S K, Ganeshan V D, Mabagala R, Sneller C, Correll J, Zhou B, Talbot N J, Mitchell T K, Wang G L. Genotyping-by-sequencing-based genetic analysis of African rice cultivars and association mapping of blast resistance genes againstpopulations in Africa., 2017, 107: 1039–1046.

[16] Park C H, Chen S B, Shirsekar G, Zhou B, Khang C H, Songkumarn P, Afzal S J, Ning Y S, Wang R Y, Bellizzi M, Valent B, Wang G L. Theeffector AirPiz-t targets the RING E3 Ubiquitin ligase APIP6 to suppress pathogen-associated molecular pattern-triggered immunity in rice., 2012, 24: 4748–4762

[17] 李旭升, 向小娇, 申聪聪, 杨隆维, 陈凯, 王小文, 邱先进, 朱小源, 邢丹英, 徐建龙. 水稻重测序核心种质资源的稻瘟病抗性鉴定与评价. 作物学报, 2017, 43: 795–810. Li X S, Xiang X J, Shen C C, Yang L W, Chen K, Wang X W, Qiu X J, Zhu X Y, Xing D Y, Xu J L. Identification and evaluation of blast resistance for resequenced rice core collections., 2017, 43: 795–810 (in Chinese with English abstract).

[18] Qiu X J, Chen K, Lyu W K, Ou X X, Zhu Y J, Xing D Y, Yang L W, Fan F J, Yang J, Xu J L, Zheng T Q, Li Z K. Examining two sets of introgression line reveals background-independent and stably expressed QTL that improve grain appearance quality in rice (L.)., 2017, 130: 951–967.

[19] Qiu X J, Pang Y L, Yuan Z Y, Xing D Y, Xu J L, Dingkuhn M, Li Z K, Ye G Y. Genome-wide association study of grain appearance and milling quality in a worldwide collection ofrice germplasm., 2015, 10: e0145577.

[20] McCouch S R, Wright M H, Tung C W, Maron L G, McNally K L, Fitzgerald M, Singh N, DeClerck G, Agosto-Perez F, Korniliev P, Greenberg A J, Naredo M B N, Mercado A M Q, Harrington S E, Shi Y, Branchini D A, Kuser-Falcao P R, Leung H, Ebana K, Yano M, Eizenga G, McClung A, Mezey J. Open access resources for genome-wide association mapping in rice., 2016, 7: 10532.

[21] Purcell S, Neale B, Todd-Brown K, Thomas L, Ferreira M A, Bender D, Maller J, Sklar P, de Bakker P I, Daly M J, Sham P C. PLINK: a tool set for whole-genome association and population-based linkage analyses., 2007,81: 559–575.

[22] 杨楠, 陈恺茜, 杨勤忠, 唐文强, 张文龙, 何平, 杨暮英, 李勇成, 韩光煜. 云南地方籼粳稻稻瘟病抗性和农艺性状差异分析及优异稻种筛选. 南方农业学报, 2021, 52: 2680–2689. Yang N, Chen K X, Yang Q Z, Tang W Q, Zhang W L, He P, Yang M Y, Li Y C, Han G Y. Blast resistance and agronomic traits of local Indica rice and Japonica rice in Yunnan and screening of elite rice varieties., 2021, 52: 2680–2689 (in Chinese with English abstract).

[23] Deng Y W, Zhai K R, Xie Z, Yang D Y, Zhu X D, Liu J Z, Wang X, Qin P, Yang Y Z, Zhang G M, Li Q, Zhang J F, Wu S Q, Milazzo J, Mao B Z, Wang E T, Xie H A, Thareau D, He Z H. Epigenetic regulation of antagonistic receptors confers rice blast resistance with yield balance., 2017, 355:962–965.

[24] Jiang J F, Mou T M, Yu H H, Zhou F S. Molecular breeding of thermos-sensitive genic male sterile (TGMS) lines of rice for blast resistance usinggene., 2015, 8: 11.

[25] Hayashi K, Yoshida H. Refunctionalization of the ancient rice blast disease resistance geneby the recruitment of a retrotransposon as a promoter., 2009, 57: 413–425.

[26] Barman S R, Gowda M, Venu R C, Chattoo B B. Identification of a major resistance gene in the rice cultivar “Tetep”., 2004, 123: 300–302.

[27] Takahashi A, Hayashi N, Miyao A, Hirochika H. Unique features of the rice blast resistancelocus revealed by large scale retrotransposon-tagging., 2010, 10: 175.

[28] Lin F, Chen S, Que Z Q, Wang L, Liu X Q, Pan Q H. The blast resistance geneencodes a nucleotide biding site-leucine-rich repeat protein and is a member of a resistance gene cluster on chromosome 1., 2007, 177: 1871–1880.

[29] Wang Z X, Yano M, Yamanoouchi U, Iwamoto M, Monna L, Hayasaka H, Katayose Y, Sasaki T. Thegene for rice blast resistance belongs to the nucleotide binding and leucine-rich repeat class of plant disease resistance genes., 1999, 19: 55–64.

[30] Sallaud C, Lorieux M, Roumen E, Tharreau D, Berruyer R, Svestasrani P, Garsmeur O, Ghesquiere A, Notteghem J L. Identification of five new blast resistance genes in the highly blast-resistant rice variety IR64 using a QTL mapping strategy.,2003, 106: 794–803.

[31] Guo L Y, Zhao H W, Wang J G, Liu H L, Zheng H L, Sun J, Yang L M, Sha H J, Zou D T. Dissection of QTLs alleles for blast resistance based on linkage and linkage disequilibrium mapping inrice seedlings., 2016, 45: 209–218.

[32] Cesari S, Thilliez G, Ribot C, Chalvon V, Michel C, Jauneau A, Rivas S, Alaux L, Kanzaki H, Okuyama Y, Morel J B, Fournier E, Tharreau D, Terauchi R, Kroj T. The rice resistance protein pair RGA4/RGA5 recognizes theeffectors AVR-Pia and AVR1-CO39 by direct binding., 2013, 25: 1463–1481.

[33] Okuyama Y, Kanzaki H, Abe A, Yoshida K, Tamiru M, Saitoh H, Fujibe T, Matsura H, Shenton M, Galam D C, Undan J, Ito A, Sone T, Terauchi R. A multifaceted genomics approach allows the isolation of the rice Pia-blast resistance gene consisting of two adjacent NBS-LRR protein genes., 2011, 66: 467–479.

[34] Hua L X, Wu J Z, Chen C X, Wu W H, He X Y, Lin F, Wang L, Ashikawa I, Matsumoto T, Wang L, Pan Q H. The isolation of, an allele at thelocus which confers broad spectrum resistance to rice blast., 2012, 25: 1047–1055.

[35] Zhang D D, Liu M X, Tang M Z, Dong B, Wu D X, Zhang Z G. Repression of microRNA by silencing ofactivities the basal resistance toin rice., 2015, 237: 24–32.

[36] Fujiwara T, Maisonneuve S, Isshiki M, Mizutani M, Chen L T, Wong H L, Kawasaki T, Shimamoto K. Sekiguchi lesion gene encodes a cytochrome P450 monooxygenase that catalyzes conversion of tryptamine to serotonin in rice., 2010, 85: 11308–11313.

[37] Lee J R, Park S C, Kim M H, Jung J H, Shin M R, Lee D H, Cheon M G, Park Y, Hahm K S, Lee S Y. Antifungal activity of rice Rex5p, a receptor for peroxisomal matrix proteins., 2007, 359: 941–946.

[38] You X M, Zhu S S, Zhang W W, Zhang J, Wang C M, Jing R N, Chen W W, Wu H M, Cai Y, Feng Z M, Hu J L, Yan H G, Kong F, Zhang H, Zheng M, Ren Y L, Lin Q B, Cheng Z J, Zhang X, Lei C L, Jiang L, Wang H Y, Wan J M.regulates rice spikelet development through modulating jasmonic acid biosynthesis., 2019, 224: 712–724.

[39] Lu T, Diao Z J, Yang D W, Wang X, Zheng X X, Xiang X Q, Xiao Y P, Chen Z W, Wang W, Wu Y K, Tang D Z, Li S P. The 14-3-3 proteinpositively regulates immunity by modulating the protein homoeostasis of the GRAS protein OsSCL7 in rice., 2022, 45: 1065–1081.

[40] Mei C S, Qi M, Sheng G Y, Yang Y N. Inducible overexpression of a rice allene oxide synthase gene increases the endogenous jasmonic acid level,gene expression, and host resistance to fungal infection., 2006, 19: 1127–1137.

[41] Ke Y G, Yuan M, Liu H B, Hui S G, Qin X F, Chen J, Zhang Q L, Li X H, Xiao J H, Zhang Q F, Wang S P. The versatile functions of OsALDH2B1 provide a genic basis for growth-defense trade-offs in rice., 2020, 11: 3867–3873.

[42] Jiang C J, Shimono M, Maeda S, Inoue H, Mori M, Hasegawa M, Sunano S, Takatsuji H. Suppression of the rice fatty-acid desaturase geneenhances resistance to blast and leaf blight diseases in rice., 2009, 22: 820–829.

Identification of rice blast resistance inandgermplasms by genome-wide association study

ZHOU Hai-Ping1, ZHANG Fan2, CHEN Kai3, SHEN Cong-Cong3, ZHU Shuang-Bing3, QIU Xian-Jin4,*, and XU Jian-Long2,3,5,*

1Wenzhou Academy of Agricultural Sciences / South-Zhejiang Crop Breeding Key Laboratory, Wenzhou 325000, Zhejiang;2Institute of Crop Sciences, Chinese Academy of Agricultural Sciences, Beijing 100081, China;3Shenzhen Branch, Guangdong Laboratory for Lingnan Modern Agriculture, Agricultural Genomics Institute at Shenzhen, Chinese Academy of Agricultural Sciences, Shenzhen 518120, Guangdong, China;4College of Agriculture, Yangtze University, Jingzhou 434025, Hubei, China;5Hainan Yazhou Bay Seed Laboratory, Sanya 572024, Hainan, China

Rice blast is one of the major fungal diseases that threaten rice production worldwide. To improve rice blast resistance, identifying blast resistant genes and introgressing them into elite rice varieties is an effective way. In this study, a panel of 212accessions and 235accessions collected worldwide were evaluated for resistance against five blast isolates at seedling stage. All of them showed large variations in resistance against five isolates, and 8and 12accessions were detected to present resistance to all five blast isolates. Using genome-wide association strategy, a total of 43 QTLs were identified for resistance to five isolates in mix population (subpopulation andsubpopulation), including 9, 4, 14, 14, and 2 QTLs for GD00-193, GD08-T19, GD17-CQ16, HB1708, and HLJ13-856, respectively. Among them, 12 resistant QTLs were detected only inricesub-population, 7 only detected inrice sub-population, and 1 simultaneously detected in bothsub-populations, indicating that blast resistance was generally better inthan inrice, and there was obvious differentiation in blast resistance betweenandrice. A total of 11 QTLs affected resistance to two or more trains or were simultaneously identified in two or more populations, and 23 candidate genes were identified by candidate interval association analysis and haplotype analysis. Different resistance candidate genes had different frequencies inandpopulations. The results provide germplasm resources and favorable genes information for molecular improvement of blast resistance in rice varieties and the breeding and utilization strategies of different resistance genes.

blast resistance; germplasm; genome-wide association study; quantitative trait locus (QTLs); favorable allele

10.3724/SP.J.1006.2023.22024

本研究由海南崖州湾种子实验室揭榜挂帅项目(B21HJ0216)和温州市农业新品种选育协作组专项(2019ZX006)资助。

This study was supported by the Hainan Yazhou Bay Seed Laboratory (B21HJ0216) and the Wenzhou Agricultural New Variety Breeding Collaborative Project (2019ZX006).

徐建龙, E-mail: xujianlong@caas.cn; 邱先进, E-mail: xjqiu216@yangtzeu.edu.cn

E-mail: zhouhaiping7583@126.com

2022-04-24;

2022-09-05;

2022-09-15.

URL: https://kns.cnki.net/kcms/detail/11.1809.S.20220915.0835.002.html

This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).