李娜，徐慧芳，黄国勤

江西农业大学生态科学研究中心，南昌 330045

水稻是世界三大粮食作物之一。2019 年，我国水稻收获面积为6 092.19 万hm2，居世界首位，约占全球水稻收获面积的1/3（http：//www.fao.org/home/en/）。稻田土壤长期处于淹水厌氧状态，反硝化作用强烈，导致稻田氮素气态损失严重，氮肥利用率低下，约为30%~40%[1］。反硝化产物氧化亚氮（N2O）是一类重要的温室气体[2］。因此，稻田反硝化作用过程成为人们关注的热点问题。

水稻田属于典型的人工湿地生态系统[3］，土壤生境长期处于淹水厌氧条件下，其中N2O 的产生和稻田的氮素转化受到反硝化作用影响。稻田氮素的积累、迁移和流失等循环转化过程会导致施入稻田中氮肥的流失、氮肥利用率下降，直接或间接地影响稻田生态系统健康和土壤可持续发展等问题。近年来，国内外针对稻田反硝化过程的研究主要在稻田的不同施肥配比[4］、稻田水分调控[5］、稻田根际土壤[6］、反硝化抑制剂[7］、稻虾共作[8］、秸秆还田[9］及种植绿肥作物[10］等方面，但对不同轮作休耕模式下稻田的特殊生境研究较少。稻田长期大量化肥的投入及土壤高强度连作种植，在短期内会增加粮食的产量，但从长期来看，会产生很多对生态环境的不利影响[11］。连作和单一作物种植导致根系分泌的有机酸和酚酸积累、土壤养分枯竭、病原微生物滋生和生物多样性降低等[12］，最终影响作物生产及土壤环境的可持续性发展。与连作或单一种植模式相比，轮作休耕是改善土壤环境质量和提高作物生产力的替代模式，在增加作物多样性的同时储蓄了土壤养分[13］。因此，轮作休耕农业管理制度在世界许多国家越来越受欢迎，轮作休耕对稻田土壤反硝化过程的影响值得探索。

本试验通过研究不同轮作休耕模式下土壤理化性质、土壤反硝化潜势及土壤反硝化功能基因nirS和nirK 群落结构的变化，分析土壤反硝化功能微生物种群的丰度和群落组成，并明确其变化的重要影响因子及关键菌属，旨在为提高农田生态系统氮素利用效率提供科学参考。

1 材料与方法

1.1 试验地概况

于2020 年10 月至2021 年11 月在 江西省鹰潭市余江区农业科学研究所试验田（28°14′8″N，116°51′22″E）进行，该试验地属亚热带气候，周年总太阳辐射量为4 542.7 kJ/m2，周年平均日照时长1 852.4 h；周年平均气温17.7 ℃，极端最高温40.5 ℃，极端最低温5.1 ℃；周年平均降雨量1 796.8 mm；试验地为泥沙淤土，土质肥沃，偏微酸性。种植前的土壤理化性质为pH 4.28；碱解氮159.6 mg/kg，有效磷31.98 mg/kg，速效钾26.67 mg/kg，全氮1.86%，有机质36.19 g/kg，NH4+-N 3.73 mg/kg，NO3−-N 8.79 mg/kg，容重2.23 mg/kg，含水率22.56%。

1.2 试验设计

本研究5 种种植模式，A（CK）为传统轮作模式，B、C、D为3个多种作物轮作模式，E处理为休耕模式（表1）。每个处理3 次重复，每个小区面积为66.7 m2，随机区组排列，小区周边设置同步种植同种品种的水稻保护区。小区间的田埂、保护行和环沟的宽度分别为0.5、1.0 和0.5 m。试验作物品种：紫云英（余江大叶籽），油菜（华赣油1号），早稻（中早33），晚稻（黄华占），春大豆（沪鲜豆6 号），秋大豆（奎鲜二号），玉米（赣1号），甘蔗（赣紫皮果蔗），甘薯（赣南瓜红薯）。紫云英播种量约为22.5 kg/hm2，播种时间2020 年10 月，2021 年4 月紫云英在盛花期直接翻压还田；油菜种植密度约为11.1 万株/hm2，移栽时间在2020年12月，收获时间2021年5月；早稻播种量约为45 kg/hm2，移栽时间在2021 年4 月，收获时间在2021 年7 月；晚稻播种量约为42 kg/hm2，移栽时间在2021年7月，收获时间在2021年11月；大豆种植密度约为6.7万株/hm2，春大豆移栽时间在2021 年4 月，收获期在6 月，秋大豆移栽时间在2021 年6 月，收获时间在2021 年11 月；玉米播种量为6.7 万株/hm2，播种时期在2021 年7 月，收获时间2021 年10 月；甘蔗种植密度约为8 230 株/hm2，移栽时间在2021 年5月，收获时间在2021 年12 月；甘薯种植密度约为5.6万株/hm2，移栽时间在2021 年7月，收获时间在2021年10 月。作物施肥量采用当地习惯施肥量，土壤取样时间为试验收获后进行取样。

表1 试验设计Table 1 Test design

1.3 土壤理化性质分析

用5点取样法分别取每小区0~20 cm耕层土，混合土样，带回实验室自然风干后分析土壤理化性质。用pH 计测定土壤pH，采用重铬酸钾法-浓硫酸外加热法测定有机质，NH4OAc 浸提-火焰光度法测定土壤速效钾含量，NaHCO3浸提-钼锑抗比色法测定土壤有效磷含量，半微量凯氏法测定土壤全氮，KCl 浸提-靛酚蓝比色法测定铵态氮，紫外分光光度法测定硝态氮，碱解扩散法测定土壤碱解氮，环刀法测定土壤容重和土壤含水率。

1.4 土壤反硝化潜势测定

采用乙炔抑制法测定土壤反硝化潜势，参照Pell等[14］的方法进行。

1.5 土壤反硝化功能基因测定

从每个小区随机挑选5株作物，用无菌勺子刮取其根系表面土壤，装入50 mL 无菌离心管中，加入液氮，用干冰将土壤样品运到实验室于−80 ℃冰箱储存。再将取好的土样委托上海美吉生物医药科技有限公司扩增测序。根据biowest agArose（biowest，E.S.）说明书进行土壤微生物群落总DNA 抽提。nirS基因采用的引物为：cd3aF_R3cdR（5' -GTSAAC⁃GTSAAGGARACSGG - 3'、5'- AGTTCTGS⁃GTRGGCTTSAG-3'），nirK基因引物为：F1aCu_R3 Cu（5' - TTAGCATGGAATAATRRAATAGGA -3'、5' - TCTGGACCTGGTGAGTTTCC - 3'）[15］。扩增条件（PCR 仪：ABI GeneAmp® 9700 型）为：95 ℃预变性3 min；95 ℃30 s，55 ℃30 s，72 ℃30 s，共27个循环；然后在72 ℃延伸10 min；最后于4 ℃保存。PCR 反应体系为：5×TransStart FastPfu缓冲液4 μL，2.5 mol/L dNTPs 2 μL，上游引物（5 μmol/L）0.8 μL，下游引物（5 μmol/L）0.8 μL，TransStart Fast⁃PfuDNA 聚合酶0.4 μL，模板DNA 10 ng，ddH2O 补足至20 μL，每个样本3 个重复。利用Illumina 公司的Illumina Miseq PE300 平台（上海美吉生物医药科技有限公司）进行测序。

1.6 数据分析

根据条形码的精确匹配，对不同样品进行nirS、nirK基因配对端读码解复，与条形码序列出现1个或多个错配，或至少出现2个与引物序列错配的读数被丢弃[16］。然后利用FLASH（fast length adjusting of short reads）将得到的读数进行合并，合并后的读数至少有10 bp 的重叠区域[17］。删除长度超出预期200~500 bp 的序列或含有任何模糊碱基的序列，获得的干净读取数分别为：421 362和232 177。经过筛选和嵌合体去除后，使用UCLUST 以97%的序列同源性进行de novo操作分类单元（OTU）的筛选[18］。最后，分别获得nirS和nirK基因OTUs 1 642和1 653。

采用Microsoft Excel 2019 处理数据；采用SPSS软件（version18.0，Chicago，Illinois，USA）对土壤理化性质进行方差分析（ANOVA），Duncan’s 法进行多重比较（α=0.05）。采用QIIME（Version 1.9.0）软件估算nirS、nirK基因的α-多样性，采用未加权的UniFrac 距离来计算多样性[19］。采用非度量多维度尺度分析（NMDS）可视化样本间群落结构的差异[20］。相似性分析（ANOSIM）被用来评估群落结构的相似性和差异性[21］。用群落barplot分析图来表示不同样品间nirS、nirK基因的组成[22］。利用Canoco 5.0（Microcoputer Power，Ithaca，NY，USA），采用基于Bray-Curtis 距离的冗余分析（RDA）来解释环境因子对反硝化群落的影响。

2 结果与分析

2.1 不同轮作休耕模式下土壤理化性质

如表2 所示，土壤pH、碱解氮、速效磷、速效钾、NO3−-N 均为A（CK）处理最低（P<0.05），且均为B、C 处理的养分指数较高。而全氮、有机质、NH4+-N 含量中，A（CK）、B 和E 处理显著高于其他处理（P<0.05）。土壤化学性质均为B、C 处理更有利于土壤肥力的增加，利于土壤氮含量的积累。在土壤物理性质中，A（CK）处理的土壤容重显著高于其他处理（P<0.05），分别高出0.81%、0.76%、0.86%、0.90%。E 处理的土壤含水率显著高于其他处理，其次是A 处理，E 处理显著高出了32.14%、101.86%、83.58%、780.87%。E 处理能有效降低土壤容重，改善土壤团粒结构，且保水性较好。结果表明，不同轮作休耕模式均能改善土壤肥力和土壤结构，并抑制了土壤氮含量的流失，其中B（紫云英-早稻-玉米||甘薯）、C（油菜-甘蔗||春大豆）和E（休耕）处理效果较好。

表2 不同轮作休耕模式下土壤理化性质Table 2 Soil physicochemical properties

2.2 不同轮作休耕模式下土壤反硝化潜势

由图1 可见，各处理的土壤反硝化潜势中A（CK）处理最高，其次是E 处理，C 处理最低。各处理 分 别 低 于A（CK）处 理82.49%、93.23%、89.43%、70.08%。通过Pearson 相关性分析（表3），可知土壤反硝化潜势与pH 呈极显负相关（r2=−0.935，P<0.01），而与AP、AK 和NO3−-N呈现显著负相关，相关系数分别为r2=−0.992、r2=−0.938 和r2=−0.944（P<0.05）。上述结果表明，除传统轮作模式A（CK，紫云英-早稻-晚稻）以外，其他轮作休耕模式均减少了氮素的气态损失，pH、AP、AK、NO3−-N与土壤反硝化潜势呈显著性相关。

表3 土壤反硝化潜势与理化性质相关性分析Table 3 Correlation analysis between soil denitrification potential and physical chemical properties

图1 不同轮作休耕模式对土壤反硝化潜势的影响Fig.1 Effects of different crop rotation-fallow patterns on denitrification potential

2.3 不同轮作休耕模式下土壤nirS 和nirK 基因型微生物的多样性

由表4可知，土壤反硝化功能基因的α-多样性指数受轮作休耕模式影响存在显著差异（P<0.05）。nirS基因除Simpson 指数是B 处理最高（0.110 2）以外，其他指数均为E 处理最高，分别高出A（CK）处理11.73%、64.50%、43.52%。nirK基因α-多样性指数各处理间无显著性差异。

表4 土壤反硝化功能基因α-多样性Table 4 Study on soil denitrification function gene α-diversity

由图2 可以看出，在β-多样性指数方面，nirS基因除B、C 和D 处理有部分重叠以外，其他处理间无重叠，通过ANOSIM 相似性检验β-多样性指数，P值为0.013，表明各处理间有显著差异（P<0.05）。nirK基因中除A（CK）和E处理没有与其他处理间有重叠以外，各处理均有重叠，进一步进行ANOSIM相似性检验，P值为0.015，表明各处理间有显著性差异（P<0.05）。

图2 土壤反硝化功能基因nirS（A）和nirK（B）的β-多样性指数Fig.2 Study on β-diversity of soil denitrification function gene nirS（A）and nirK（B）

结果表明，不同轮作休耕模式的反硝化功能基因的多样性指数有显著性差异（P<0.05），且nirS基因的多样性指数及群落结构受影响程度大于nirK基因，E处理（休耕）更有利于土壤nirS基因α-多样性指数的提高，B、C、D 处理的nirS和nirK基因群落结构相似。

2.4 稻田不同轮作休耕模式下土壤nirS 和nirK 基因菌属相对丰度

属水平下，各处理的反硝化功能基因的相对丰度具有显著性差异（P<0.05）。由图3 可见，不同轮作休耕模式的土壤反硝化基因菌属具有一致性，但不同处理下菌属相对丰度不同。nirS基因的优势菌属 有unclassified_k_norank_d、unclassified_p_Proteo⁃ba、红杆菌属Rhodanobacter。其中B 处理（72.27%）在unclassified_K_norank_d 菌属相对丰度最高，C（20.21%）和E 处理（21.76%）在unlasified_p_Proteo⁃ba相对丰度高于其他处理，而D 处理（12.70%）在红杆菌属Rhodanobacter丰度最高。nirK基因优势菌属有unclassified_k_norank_d_Bacteria、norank_d_Bacte⁃ria、大豆根瘤菌Bradythizobium、unclassified_p_Pro⁃teobacteria、亚硝化螺菌属Nitrosospira。D（24.26%，14.70%）和E 处理（18.48%，14.97%）中norank_d__Bacteria和大豆根瘤菌Bradythizobium的相对丰度高于其他处理，B 处理（15.46%，12.12%）的unclassi⁃fied_p_Proteobacteria和亚硝化螺菌属Nitrosospira相对丰度高于其他处理。A 处理（55.97%）仅有un⁃classified_k_norank_d_Bacteria的相对丰度高于其他处理。综上，不同轮作方式下nirS、nirK 基因在属水平的相对丰度存在差异（P<0.05），其中B、D、E处理的菌属相对丰度高于其他处理。

图3 nirS基因（A）和nirK基因（B）群落组成相对丰度Fig.3 nirS gene（A）and nirK gene（B）community composition abundance

2.5 nirS 和nirK 反硝化功能基因的优势菌属与土壤反硝化潜势的相关性

由表5 可知，通过Pearson 相关性分析发现，nirS基因的优势菌属与土壤反硝化潜势无显著相关性（P>0.05）。nirK基因的大豆根瘤菌Bradythizobi⁃um和unclassified_p_Proteobacteria是影响土壤反硝化潜势的关键菌属（P<0.05），其中大豆根瘤菌Bra⁃dythizobium与土壤反硝化潜势呈现负相关性（r2=−0.85，P<0.05），unclassified_p_Proteobacteria与土壤反硝化潜势呈现正相关性，（r2=0.88，P<0.05）。

表5 土壤反硝化潜势与优势菌属相关性分析Table 5 Correlation analysis of soil denitrification poten⁃tial and dominant bacterial genus

2.6 nirS 和nirK 反硝化功能基因细菌群落结构与土壤理化性质的相关性

基于RDA 分析，采用Bray-Curtis 距离冗余度分析对群落组成排序的环境变量进行分析（图4），结果显示，反硝化功能基因的群落结构与土壤理化性质之间有显著相关性（P<0.05），SWC 是影响nirS、nirK基因群落结构的重要影响因子，相较于其他土壤理化性质，SWC 的贡献度较高，显著性影响较大。不同轮作休耕模式下，SWC 是土壤反硝化功能基因群落结构的重要关键因子，水田转旱作的种植模式改变了稻田土壤的含水率，从而影响了土壤反硝化功能基因群落结构。

图4 土壤反硝化基因nirS（A）、nirK（B）群落与土壤环境因子的相关性Fig.4 Correlation analysis of soil denitrification gene nirS（A），nirK（B）communities with soil environmental factors

3 讨 论

本研究探讨了不同轮作休耕模式下土壤反硝化潜势和反硝化功能微生物（nirS/nirK）群落的多样性和组成，结果表明，不同轮作休耕处理改变了土壤肥力及土壤结构，土壤反硝化潜势及反硝化功能基因群落结构均受到显著性影响。传统轮作模式A（CK）（紫云英-早稻-晚稻）对土壤综合肥力及土壤结构的改善较为缓慢，土壤反硝化潜势也较大，使得土壤氮素气态损失严重。稻田适当采用水旱轮作和休耕模式可以通过不同作物相互作用及土壤生境的改变，增加nirK和nirS功能基因的物种多样性及菌属的相对丰度，不仅有利于土壤环境质量的提高，还促进了稻田土壤氮素循环的可持续发展。

3.1 稻田不同种植模式下土壤反硝化潜势及土壤理化性质相关性

5 种稻田轮作模式下土壤反硝化潜势存在显著差异，pH、AP、AK 和NO3−-N 是影响土壤反硝化潜势的关键因子。含水量较高的A（CK）（紫云英-早稻-晚稻）处理反硝化潜势较强，可能由于长期淹水条件导致土壤容易酸化，土壤结构不良，加快了土壤氮的流失。研究表明，稻田土壤淹水上层（氧化层）中的氮在硝化细菌的作用下转化为硝态氮，而NO3−-N不能被土壤胶粒所吸附，随水渗漏于下边的土壤层，逐渐在反硝化细菌作用下还原成水稻难以吸收利用的气体逸失于大气中，土壤反硝化潜势增强[23］。其他轮作休耕模式的反硝化潜势较低，原因可能是土壤环境由厌氧变成好氧，改善了土壤结构，有利于土壤养分的吸收，有效阻止土壤次生潜育化和酸化。Shcherbak 等[24］研究发现，土壤pH 增加有效抑制了N2O 的排放，降低了反硝化潜势。陈盟[25］研究表明，当土壤含水量稍高时，形成低氧或厌氧条件，土壤中N2O 的排放主要是受土壤的反硝化过程所主导。C处理（油菜-甘蔗||春大豆）的N2O 释放量和反硝化潜势最低，原因可能是大豆作为养地作物有效减少了甘蔗施肥量，从而减少蔗田土壤N2O 的释放量和反硝化潜势。这与龙鹏宇等[26］的研究结果相似，他们认为合理的滴灌施肥措施能有效提高甘蔗产量且减少土壤N2O 的排放。邢肖毅等[27］也发现pH 是影响土壤反硝化作用的最重要因素，认为较低的pH 会限制土壤反硝化微生物的生长。综上，稻田长期淹水下NO3−-N 不容易被土壤胶粒所吸附，并还原成水稻难以吸收的气体，导致N2O 大量排放、土壤结构遭到破坏、土壤养分降低、土壤酸化，加快了土壤反硝化的作用。

3.2 不同稻田种植模式下土壤nirS 和nirK 基因群落结构

本试验中稻田不同轮作休耕模式对nirS 基因群落结构的影响大于nirK 基因，可能是试验年限较短，nirK 基因不敏感，不容易受短时间环境的改变所影响。朱杰等[8］发现由于稻虾种植年限较短，nirK 基因微生物多样性指数无显著性差异。尽管nirK 和nirS 行使的功能相同，但nirK 的丰度与磷或钾的供应量有关，而nirS 的丰度由氮的供应量决定[28］。谢婉玉等[29］研究表明秸秆还田处理提供了提供氮底物，对nirS 基因丰度有显著影响，对nirK 基因丰度没有显著影响。本研究除E 处理（休耕）外均种植豆科作物，豆科作物可固氮，所以nirS 基因的物种多样性出现显著性差异。试验中E 处理更有利于土壤nirS基因α-多样性指数的提高，其原因可能是休耕处理的土壤容重较低，保水性较好，土壤养分不易流失。在本研究中，种植豆科绿肥作物及施肥处理nirS、nirK 功能基因丰度均有显著性差异。本研究中unla⁃sified_p_Proteoba 是nirS 和nirK 型反硝化细菌的共同优势菌属，这与Snaidr 等[30］对土壤反硝化微生物的研究结果吻合。nirS、nirK 功能基因中的各菌属丰度均为B 处理（紫云英-早稻-玉米||甘薯）较高，主要是unclassified_k_norank_d、unclassified_p_Proteobac⁃teria、亚硝化螺菌属Nitrosospira 菌属的丰度较高，可能由于该处理作物类型多且采用了水旱轮作模式，使得根系分泌物及作物的凋落物增加，相应增加了细菌相对丰度。D处理（紫云英-春大豆-秋大豆）和E处理（休耕）红杆菌属Rhodanobacter 和大豆根瘤菌Bradythizobium 相对丰度较高，且只有D处理中有博斯氏菌属（Bosea），可能由于该模式有养地作用（种植豆科作物和休耕），所以增加了有益菌属的丰度。衷炜华等[31］研究也表明，红杆菌喜光合好氧细菌，所以，红杆菌属Rhodanobacter 在旱地种植模式中相对丰度较高。

通过Pearson 相关性分析发现，大豆根瘤菌Bra⁃dythizobium 和unclassified_p_Proteobacteria 是 影 响土壤反硝化潜势的关键菌属（P<0.05），其中大豆根瘤菌Bradythizobium 与土壤反硝化潜势呈现负相关，unclassified_p_Proteobacteria 与土壤反硝化潜势呈正相关，原因可能是大豆根瘤菌Bradythizobium能够固氮，可以减少土壤中的氮素损失。程玉柱等[32］研究表明禾本科作物可以利用豆科作物所固定的氮素，促进对氮素的吸收。unclassified_p_Proteo⁃bacteria 属于变形菌纲，有研究表明反硝化菌主要为原核生物，大量存在于在α-，β-和y-变形菌纲中[23］。冗余度分析表明，影响土壤反硝化功能基因群落结构的重要关键因子是土壤含水率，水田转旱作或者休耕的种植模式改变了稻田土壤生境，土壤水分减少，对土壤反硝化功能基因群落结产生了影响，周维[33］研究也证实土壤水分对反硝化微生物群落组成的影响最大且最敏感。