Changes in Soil Bacterial and Nematode Communities during Long-term Continuous Cotton Cropping in an Arid Region

© 2024 Academia Sinica, Taiwan Open Access Changes in Soil Bacterial and Nematode Communities during Long-term Continuous Cotton Cropping in an Arid Region Hong Chen1,2, Fenghua Zhang1,2,*, Lei Yang1,2, Zaixin Li1,2, and Jiajie Liu1,2 1Agricultural College, Shihezi University, Shihezi, Xinjiang 832003, P.R.China. *Correspondence: E-mail: [email protected] (Zhang) E-mail: [email protected] (Chen); 1[email protected] (Yang); [email protected] (Li); [email protected] (Liu) 2Key Laboratory of Oasis Eco-agriculture, Xinjiang Production and Construction Corps, Shihezi University, Shihezi, Xinjiang 832003, P.R.China Received 4 October 2022 / Accepted 27 December 2023 / Published 14 May 2024 Communicated by Benny K.K. Chan Clarifying the effects of continuous cotton cropping (CC) on soil biological communities is essential for maintaining agricultural productivity. In this study, high-throughput sequencing was used to study the effects of different CC durations (0-yr, 5-yr, 10-yr, 15-yr, 20-yr, and 25-yr CC treatments) on soil microbial and nematode communities. The results showed that the dominant bacterial phyla were Actinobacteria and Proteobacteria, and the dominant nematode genus was Helicotylenchus in all CC treatments. The richness indexes (ACE and Chao1 index) and diversity index (Shannon index) of bacterial and nematode communities were the highest in the 15-yr and 10-yr CC treatments, respectively. Bacterial community was significantly correlated with soil pH and available potassium (AK), and nematode abundance was significantly correlated with microbial biomass carbon (MBC). Soil bacterial PICRUSt analysis results showed that carbon metabolism and amino acid metabolism were the main metabolic functions of bacteria in the CC treatments. The composition and diversity of soil nematode communities were significantly related to the structure of soil bacterial communities, and the niche breadth of soil bacteria was negatively correlated with that of nematodes. Panagrolaimus and Acrobeles were the main genera of bacterialfeeding nematodes affecting bacterial communities, and their relative abundances were significantly positively correlated with the relative abundance of bacterial communities. Overall, long-term (10–15 years) continuous cotton cropping negatively impacts soil biota and the microecological environment of cotton fields in arid regions. Key words: Cotton, Continuous cropping, Soil bacteria, Soil nematode, High-throughput sequencing Citation: Chen H, Zhang F, Yang L, Li Z, Liu J. 2024. Changes in soil bacterial and nematode communities during long-term continuous cotton cropping in an arid region. Zool Stud 63:05. doi:10.6620/ZS.2024.63-05. BACKGROUND The Xinjiang Uygur Autonomous Region is the largest cotton-growing area in China, and its cotton acreage and output account for 82.8% and 89.5% of the country’s total. respectively (National Bureau of Statistics 2021). Continuous cotton cropping (CC) is very common in Xinjiang. However, long-term continuous cropping brings a series of negative effects, causing declines in crop yield and quality, soil nutrient imbalance, soil salinization, autotoxicity of secondary metabolites, and increases in soil pathogenic bacteria (Tian et al. 2010 2011 2013; Fu et al. 2017; Qin et al. 2017). This greatly restricts sustainable agricultural development in Xinjiang (Qin et al. 2017). Microbes and nematodes are the most active components of the soil ecosystem and indicate soil environmental changes (Bongers and Bongers 1998; Pascual et al. 2000; Romaniuk et al. 2011). The abundances of microbes and nematodes are greatly impacted by soil management (Coleman 2008). Some researchers found that CC could greatly affect soil Zoological Studies 63: 5 (2024) doi:10.6620/ZS.2024.63-05 1 © 2024 Academia Sinica, Taiwan microbial diversity by impacting the quality and quantity of soil resources (Li et al. 2010a; Yang et al. 2018; Tian et al. 2020). Microbes and nematodes occupy key positions in the soil food web (Bardgett and Van der Putten 2014; Manuel et al. 2016; Yurkevich et al. 2020), and play important roles in soil biogeochemical cycling and ecosystem services (Bardgett and Van der Putten 2014). Therefore, microbes and nematodes greatly influence crop growth and development (Ingham et al. 1985). Interestingly, selective predation by nematodes could stimulate microbial activity (Hu et al. 2019; Liu et al. 2019; Poromarto et al. 2019; Pawlowski and Hartman 2020), change the abundance and structure of the microbial community (Djigal et al. 2004), and inhibit diseases caused by bacteria, fungi, and pathogenic nematodes (Yeates et al. 1993a). In recent years, niche breadth has been widely concerned and applied in farmland ecosystem research (Zhang et al. 2018; Jiao et al. 2020). Some studies have found that microbiota with broad niche breadth are metabolically flexible at the community level (Pandit et al. 2009; Wu et al. 2018). However, it is still unclear how soil bacterial and nematode communities respond to longterm continuous cropping. In recent years, the application of high-throughput sequencing (HTS) technology has provided valuable semiquantitative information for the study of nematode community ecology (Porazinska et al. 2012; Griffiths et al. 2018; Gao et al. 2021). On the one hand, the physiological structure of most soil nematodes meets the DNA extraction requirement for PCR reactions (Powers 2004; Li et al. 2017). On the other hand, the identification of soil nematodes by morphological characteristics requires special taxonomic techniques and a significant amount of time (Mao et al. 2004), which limits the large-scale ecological study of soil nematodes (Lawton et al. 1998; Geisen et al. 2018). Treonis et al. (2018) showed that HTS annotated more nematode taxa and provided higher taxonomic resolution than morphological identification. Geisen et al. (2018) further confirmed this view. Therefore, highthroughput sequencing has great potential in nematode research. Kerfahi et al. (2016) and Moroenyan et al. (2016) studied the diversity of nematodes in tropical rain forests, Arctic regions, and the mediterranean region using HTS. Du et al. (2020) also used HTS technology to quantitatively investigate the nematode community composition in the soil of a steppe-forest ecotone. However, HTS technology also has some problems. For example, different primers and databases may produce different results. This study used Illumina Miseq sequencing to study how soil bacterial and nematode communities respond to long-term CC in arid areas in Xinjiang, China. The objectives were to clarify: (1) the responses of bacterial and nematode communities to long-term CC, (2) the changes in soil properties leading to the changes in soil nematode and bacterial communities, and (3) the relationship between soil nematode and bacterial community composition during longterm CC. This study will provide a reference for soil biological regulation of agro-ecosystems, and is of great significance to maintain farmland ecosystem health. MATERIALS AND METHODS Study site This experiment was conducted at Shihutan Town, Xinjiang, China (44°37'N, 86°8'E). This area has a temperate continental climate, with an average annual temperature of 6.6°C, an average annual precipitation of 110–200 mm, and an average annual evaporation of 1500–2000 mm. Continuous cotton cropping is a common practice in this area. Cotton yield was about 5250 kg ha-1. Experimental design This experiment adopted randomized complete block design. Due to soil salinization, the experimental site was abandoned for 29 years. From 1996 to 2016, the abandoned farmlands were separately reclaimed to plant cotton, which were marked as 25-year CC treatment, 20-year CC treatment, 15-year CC treatment, 10-year CC treatment, and 5-year CC treatment, respectively. An area with no cotton planting or anthropogenic influences since 1996 was selected as the control plot (0-yr). Each treatment had three replicates/plots. The area of each plot was 5.2 × 7 m, and 2 m buffer rows were set around each treatment. Cotton seeds (cultivar Zhuangjiahan 902) were sown (2.4 × 105 plants/hm2) in mid-April. The row spacing was 66 cm + 10 cm, and the plant spacing was 9 cm. Plastic film mulch and drip irrigation were used. Irrigation were conducted 10–12 times during the whole growth period, with a total amount of 4500 m3/hm2. Urea (150 kg/hm2), calcium superphosphate (450 kg/ hm2), and K2O (75 kg/hm2) were basally applied, and topdressing was conducted six times by a drip irrigation system from bud stage to boll opening stage (urea: 300 kg/hm2; potassium dihydrogen phosphate: 200 kg/ hm2). Other field managements were the same as local practice. To reduce the adverse effects of pests and diseases on cotton growth and development and ensure the consistency of control variables in all treatments, each treatment adopted the same insecticide application page 2 of 15Zoological Studies 63: 5 (2024) © 2024 Academia Sinica, Taiwan time and dosage. Except for the duration of continuous cropping, all management measures were consistent across all treatments. Soil sampling On July 14, 2020 (full flowering stage of cotton), twenty sampling points were selected in each plot along a “S”-shaped route. Before sampling, stones and fallen leaves on the soil surface were removed, and then soil samples (0–20 cm in depth) with a distance of 0–3 cm from cotton roots were collected. The twenty samples were mixed and used as the soil sample for the plot (500 g). A total of 18 soil samples were obtained. All soil samples were stored at 4℃, transported to the laboratory, passed through a 2-mm sieve, and divided into two portions. One portion was air-dried and stored at room temperature for soil property measurement, and the other was stored at -80℃ for molecular analysis. Determination of soil physicochemical properties Soil pH was determined using a pH meter (PHS3C, INESA Scientific Instrument CO., Ltd, Shanghai, China) (water: soil = 5:1). Soil electrical conductivity (EC) was determined using a conductivity meter (water: soil = 2.5:1) (Ciavatta et al. 1991). Soil organic matter (SOM) content was determined using the Walkley-Black dichromate oxidation method (Díaz-Zorita 1999). Soil nitrogen content was determined using the Semimicro Kjeldahl method (Bremner and Mulvaney 1982). Soil available phosphorus (AP) content was determined using the 0.5 moL L-1 NaHCO3 method (Page 1982). Soil available potassium (AK) content was determined using flame photometry (Shen et al. 2013). Soil nitrate nitrogen (NO3 --N) and ammonium nitrogen (NH4 +-N) content were determined using a flow analyzer (Bran Luebbe, Germany) (Mcleod 1992). Microbial biomass carbon (MBC) was determined using the chloroform fumigation-extraction method (Vance et al. 1987). Table 1 shows soil properties under different CC durations. DNA extraction, PCR amplification, and Illumina MiSeq sequencing The total DNA of nematodes and bacteria (from 0.5 g soil sample) were extracted using an E.Z.N.A.® soil DNA kit (Omega Bio-tek, Norcross, GA, U.S.) according to the manufacturer’s instructions. Agarose gel (1%) was used to check the DNA extract, and the concentration and purity of DNA were measured using a UV-vis spectrophotometer (NanoDrop 2000, Thermo Scientific, Wilmington, USA). The forward primer NF1-F (5'-GGTGTGCATGG CCGTTCTTAGTT-3') and the reverse primer 18Sr2b-R (5'-TACAAAGGGCAGGGACGTAAT-3') were used to amplify the variable region of soil nematode18S rDNA gene V4 using PCR (Du et al. 2020). PCR amplification of the soil bacteria 16S rRNA gene was performed using the PCR primers specific for the 338–806 (V3–V4) regions, and the following primers were used to attach a barcode (an 8-base sequence unique to each sample): 338F (5'-barcod-e-ACTCCTACGGGAGGCAGCAG-3') and 806R (5'-GGACTACHVGGGTWTCTAAT-3'). The PCR mix (20 μL) contained 4 μL of 5x TransStart™ FastPfu reaction buffer, 2 μL of 2.5 mM dNTP, 0.8 μL of forward primer (5 μM), 0.8 μL of reverse primer (5 μM), 0.4 μL of TransStart™FastPfu DNA polymerase, and 10 ng of ddH2O template DNA. The visualization and purification of PCR products were conducted using 2% agarose gel electrophoresis and the AxyPrep DNA Gel Extraction Kit (Axygen Biosciences, Union City, CA, USA) following the manufacturer’s protocol, followed by the quantification using the QuantiFluor™ -ST (Promega, USA). Purified amplicons were pooled in equimolar and paired-end sequenced (2 × 300) on an Illumina MiSeq platform (Illumina, San Diego, USA) by Majorbio BioPharm Technology Co., Ltd. (Shanghai, China). Raw fastq files were demultiplexed and qualityTable 1. Effects of different continuous cropping spans on soil properties Group pH EC (µS cm-1) SOC (g kg-1) TN (g kg-1) AP (mg kg-1) AK (mg kg-1) NH4 +-N (mg kg-1) NO3 --N (mg kg-1) MBC (mg kg-1) 0-yr 9.15 ± 0.04a6855.56 ± 218.98a5.12 ± 0.31e0.30 ± 0.03e18.82 ± 0.48f509.81 ± 4.01a75.03 ± 1.94b51.67 ± 0.92b56.48 ± 3.29e 5-yr 8.39 ± 0.09f534.00 ± 56.07b7.43 ± 0.61cd 0.61 ± 0.00c25.41 ± 2.39d268.46 ± 2.03c73.06 ± 3.67b22.66 ± 1.75d143.24 ± 7.59b 10-yr 8.92 ± 0.02b141.78 ± 6.53d11.11 ± 0.49a0.66 ± 0.02b33.10 ± 0.90c213.19 ± 0.45d97.19 ± 19.60a52.28 ± 1.78b222.82 ± 10.27a 15-yr 8.55 ± 0.04e468.00 ± 10.16b9.44 ± 0.24b0.83 ± 0.01a83.69 ± 0.35a463.66 ± 2.74b85.13 ± 8.76ab 93.94 ± 1.83a112.49 ± 14.35c 20-yr 8.80 ± 0.02d256.11 ± 8.89c6.96 ± 0.09d0.57 ± 0.00d22.30 ± 1.23e207.19 ± 17.53d83.16 ± 10.10ab 45.91 ± 2.21c84.22 ± 15.69d 25-yr 8.85 ± 0.04c276.56 ± 11.33c8.02 ± 0.12c0.56 ± 0.00d38.68 ± 1.17b273.55 ± 1.36c81.76 ± 9.68ab 8.78 ± 0.47e59.68 ± 6.24e Notes: Values are means ± standard deviation (n = 3). Different lowercase letters in the same column indicate significant difference at p < 0.05. EC, Electrical conductivity; SOC, Soil organic carbon; TN, Total nitrogen; AP, Available phosphorus; AK, Available potassium; NH4 +-N, Ammoniacal nitrogen; NO3 --N, Nitrate nitrogen; MBC, Microbial biomass carbon. page 3 of 15Zoological Studies 63: 5 (2024) © 2024 Academia Sinica, Taiwan filtered using QIIME (Quantitative Insights into Microbial Ecology, version1.9.1). Operational taxonomic units (OTUs) with a 97% similarity cutoff were clustered using USEARCH (version 7.1), and chimeric sequences were identified and removed using UCHIME (Edgar 2013). The bacterial taxonomy of each 16S rRNA gene sequence was analyzed using RDP classifier (version 2.2) based on the SILVA ribosomal RNA gene database (Release 123) (Wang et al. 2007). The representative sequences of nematode 18S rDNA were classified by RDP classifier (version 2.2) based on NCBI NT database. The confidence threshold was 70% (Du et al. 2020). Statistical analysis The data were statistically analyzed using SPSS software (version 20.0) (SPSS Inc., Chicago, IL, USA). The Chao1 richness index, ACE richness index, Shannon diversity index, and Simpson diversity index were determined using Mothur software (version 1.30.1) (Schloss et al. 2009). The similarities and differences among the communities were described using an Upset-plot diagram with common and unique OTUs. Principal component analysis (PCA) based on z-score normalization was performed using R software (version 3.6.1) (Mardia et al. 1979). Redundancy analysis (RDA) was performed using the Vegan package in R software to determine the relationship between soil bacterial and nematode communities, soil bacterial functional genes, and soil physicochemical properties (Dixon 2003). Spearman correlation coefficients between the relative abundance of major genus and soil properties were also calculated to estimate the difference between CC treatments. Bacterial functional gene prediction was performed using Phylogenetic Investigation of Communities by Reconstruction of Unobserved States (PICRUSt) (Langille et al. 2013). Functional genes were identified based on the Kyoto Encyclopedia of Genes and Genomes (KEGG) database (Ogata et al. 1999). Cluster heatmap was used to analyze the differences of pathway-3 functional genes among different treatments. To reveal habitatand region-based species classification, dispersal limitations, and their impact on soil microbial communities, the “niche breadth” function in the “spaa” package in the R software was used to estimate the niche breadth (B) of microbes (Levins 1968; Zhang and Zhang 2013). Communitylevel B-value (Bcom) was calculated as the average of the Bcoms for all taxa in a community (Wu et al. 2018). The Mantel test was used to determine the correlation between nematode and bacterial communities (Sunagawa et al. 2015). SPSS 20.0 (SPSS Inc., Chicago, USA) was used for one-way ANOVA and Tukey’s HSD test for multiple comparisons. RESULTS Composition and diversity of soil bacterial and nematode communities A total of 867,354 and 384,208 optimized sequences were obtained from the 18 soil samples through HTS analysis of soil bacteria and nematodes, respectively, and 6290 and 51 OTUs were detected, respectively (Table 2). Actinobacteriota and Proteobacteria were the dominant bacteria phyla, and the relative abundance were 21.32%–31.94% and 18.95%–27.16%, respectively (Fig. 1a, Table S1, p < 0.05). CC significantly decreased the relative abundance of Actinobacteria but increased the relative abundance of Proteobacteria, especially in the 25-yr CC treatment (p < 0.05). Other phyla with a higher relative abundance were Acidobacteria (2.61%–24.07%), Chloroflexi (8.81%–14.05%), Firmicutes (2.93%–5.89%), and Bacteroidetes (1.98%– 11.7%). The dominant genera were Gp6, Arthrobacter, Gemmatimonas, and Gp4 (Fig. 1b, Table S2, p < 0.05). There were differences in the Chao1 and Shannon indexes among the six CC treatments, and the Chao1 and Shannon indexes began to decrease after the 15-yr CC (Table 2, p < 0.05). The number of unique bacteria (genus level) in the 0-yr, 5-yr, 10-yr, 15-yr, 20-yr, and 25-yr CC treatments were 57, 2, 1, 11, 4 and 1, respectively, and the total quantity of bacterial genera decreased after 15-yr CC (Fig. S1a). Helicotylenchus (relative abundance: 19.87%– 55.13%) was the dominant nematode genus in all CC treatments except for the 20-yr CC treatment. The relative abundance of unclassified nematode genera in the 20-yr treatment was 70.59%. In particular, the relative abundances of OTU693 (28.1%), OTU282 (18.677%), and OTU281 (15.12%) were higher that those of other unclassified nematode genera (Table S4), and also higher than that of Helicotylenchus (11.18%) in the 20-yr treatment (Fig. 1c, Table S3 and S4, p > 0.05). Leptonchus was only detected in the 5-yr, 10-yr, and 15-yr CC treatments, Panagrolaimus was only detected in the 0-yr and 5-yr CC treatments, and Acrobeles was detected in all CC treatments except for the 0-yr CC treatment. The relative abundance of Panagrolaimus in the 5-yr CC treatment was reduced by 85.81% compared with that of the control. Soil nematode Chao1 and Shannon indexes of the CC treatments increased compared with those of the control (Table 2, p < 0.05). The 10-yr CC treatment had the largest number of soil nematode genera, and the number of unique genera page 4 of 15Zoological Studies 63: 5 (2024) © 2024 Academia Sinica, Taiwan in the 0-yr, 5-yr, 10-yr, 15-yr, 20-yr, and 25-yr CC treatments were 0, 4, 3, 0, 0, and 1, respectively (Fig. S1b). The analysis of similarities (ADONIS) of bacterial and nematode communities at the genus level showed that there were significant differences among the CC treatments (bacteria: R2 = 0.85, p = 0.001; nematode: R2 = 0.50, p = 0.002), which was also confirmed by the PCA analysis results (Fig. 2). The PCA axis divided bacterial communities into three parts: (1) 0-yr CC Fig. 1. Composition of the soil bacterial community at the phylum (a) and genus (b) level of the nematode community at the genus level (c) in differen treatment groups. Table 2. Effects of different continuous cropping spans on the richness and diversity of soil bacteria and nematode communities (similarity level: 97%) Microbial and nematode community Group Sequencing results Chao1 index Shannon index (H) Total no. of sequences Total no. of OUTs Bacteria 0-yr 48696 3154 2928 ± 152e 6.14 ± 0.08c 5-yr 41751 3606 3395 ± 148cd 6.22 ± 0.05bc 10-yr 44746 3804 3533 ± 32bc 6.37 ± 0.04ab 15-yr 40982 3872 3780 ± 83a 6.42 ± 0.04a 20-yr 46884 3919 3662 ± 22ab 6.38 ± 0.07ab 25-yr 40117 3531 3343 ± 33d 6.12 ± 0.19c Nematode 0-yr 216 6 2.33 ± 1.53b 0.35 ± 0.35b 5-yr 500 24 13.00 ± 6.08a 1.49 ± 0.06a 10-yr 1952 28 14.00 ± 2.65a 1.32 ± 0.33a 15-yr 224 17 8.83 ± 7.52ab 1.25 ± 0.31a 20-yr 334 17 8.67 ± 3.51ab 1.31 ± 0.56a 25-yr 797 14 7.33 ± 2.89ab 1.16 ± 0.60a Notes: Values are means ± standard deviation (n = 3). Different lowercase letters in the same column indicate a significant difference at p < 0.05. OTUs, operational taxonomic units. page 5 of 15Zoological Studies 63: 5 (2024) © 2024 Academia Sinica, Taiwan treatment, (2) 15-yr CC treatment, and (3) other CC treatments (Fig. 2a). However, the nematode community in the 10-yr CC treatment was obviously separated from those in the other CC treatments (Fig. 2b). Redundancy analysis results showed that bacteria (p = 0.001, Fig. 3a) responded strongly to soil physicochemical properties. RDA1 and RDA2 explained 64.01% and 18.56% of the total variation, respectively. It is worth noting that both pH (p = 0.007) and AK (p = 0.001) had significant effects on bacterial communities, and soil MBC (p = 0.004) was the vital soil factor affecting nematode community composition (Table 3). Potential functions of soil bacteria The PICRUSt analysis results of the metabolic functions of bacteria in CC soils showed that bacterial metabolism was the main pathway (level 1), and the relative abundance of genes in the bacterial metabolism pathway in each CC treatment was more than 79% (Fig. S2). The global and overview maps, carbohydrate metabolism, amino acid metabolism, metabolism of Fig. 2. PCA analysis of soil bacteria (a) and nematodes (b) at the genus level based on the Bray-Curtis algorithm. Fig. 3. Redundancy analysis of the top 15 bacteria (a) and nematodes (b) (at the genus level) with soil physicochemical properties. The colored dots represent relative abundance. The bars represent the total proportion of explained variations in the physicochemical properties by key bacteria and nematodes. AP, Available phosphorus; AK, Available potassium; NH4 +-N, Ammoniacal nitrogen; NO3 --N, Nitrate nitrogen; MBC, Microbial biomass carbon. (a) (a) (b) (b) page 6 of 15Zoological Studies 63: 5 (2024) © 2024 Academia Sinica, Taiwan cofactors and vitamins, and energy metabolism were the major pathways enriched by functional genes (Level 2) (Fig. S2). Cluster analysis results showed that the metabolic functions of bacteria in the 0-yr and 25-yr CC treatments were similar (level 3) (Fig. 4a). There was no difference in the abundance of carbon metabolism genes among CC treatments (p > 0.05). The abundance of bacterial oxidative phosphorylation genes and carbon fixation genes in the 20-yr CC treatment was significantly higher than those in the 0-yr and 25-yr CC treatments. However, the relative abundances of methane metabolism, sulfur metabolism, and nitrogen metabolism genes showed an opposite trend, and those in the 20-yr CC treatment were the lowest (Fig. 4b). The two components in RDA analysis explained 94.05% of the total variation in the relationship between soil bacterial functional genes and soil physicochemical properties, and soil pH, AK, and NO3 --N were the main influencing factors (Fig. 4C, Table 4). Relationship between bacterial and nematode communities The niche breadths of soil bacteria and nematodes were different in different CC treatments (p < 0.05). The bacterial community in the 0-yr CC treatment had the highest Bcom value, while that in the 20-yr CC treatment had the lowest Bcom value (Fig. 5a). This was opposite to that of nematode community (Fig. 5b). The Mantel test showed that nematode community diversity, composition, and niche breadth were the main factors affecting bacterial community composition. However, the diversity and niche breadth of the nematode community also affected the niche breadth of the bacterial community (Fig. 5c, Table S5, p < 0.05). Additionally, RDA1 and RDA2 explained 79.80% of the total variances in the relationship between bacterialfeeding nematodes and bacterial communities, and Panagrolaimus and Acrobeles were the two main bacterial-feeding nematode genera that affected soil bacterial communities (Fig. 5d, Table 5, p < 0.05). The correlation analysis results showed that the abundance of Panagrolaimus was significantly positively correlated with those of Pontibacter, Euzebya, Nitriliruptor, Truepera, Blastococcus, Aliifodinibius, Rubellimicrobium, Sphaerobacter, and Gillisia, and the abundance of Acrobeles was significantly positively correlated with that of Saccharibacteria genera incertae sedis (Table S6, p < 0.05). DISCUSSION Many factors affect the structure and diversity of soil bacterial community, such as climate, vegetation, and soil type (Xu et al. 2013; Gavrichkova et al. 2018; Sun et al. 2020). In this study, 25-yr continuous cropping significantly reduced the relative abundance of Actinobacteriota, but increased the relative abundance of Proteobacteria. Actinobacteriota and Proteobacteria are the dominant bacteria in most arid soils (Cheng et al. 2019). Doumbou et al. (2001) reported that Actinobacteria was an important microbe in rhizosphere soil, which could regulate plant pathogen resistance and soil microecological balance. The increase in the relative abundance of Proteobacteria may be because the increase in soil nitrogen content increases the uptake of other resources by soil bacteria (Fierer et al. 2012; Zhang and Du 2018; Lin et al. 2019). Besides, continuous cropping significantly increased the diversity and richness of soil bacteria, especially in the 15-yr CC treatment. This may be due to the changes in the soil environment, especially the changes in the nutrients required by soil bacteria (Cheng et al. 2018). Plant parasitic nematodes can induce inhibition to continuous cropping (Mweke et al. 2008). This study showed that Table 3. Importance of soil physicochemical properties in explaining the bacterial and nematode community structure changes obtained by redundancy analysis Soil properties Bacteria Nematode RDA1 RDA2 R2Pr (> r) RDA1 RDA2 R2Pr (> r ) pH -0.9737 -0.2280 0.4786 0.007** -0.6547 -0.7559 0.0027 0.982 AP 0.6945 0.7195 0.0548 0.687 0.9735 -0.2285 0.0494 0.685 AK -0.9769 0.2138 0.6861 0.001** -0.9484 -0.3169 0.2300 0.128 NH4 +-N 0.9995 0.0315 0.0590 0.645 0.9647 -0.2632 0.0897 0.506 NO3 --N -0.5363 -0.8440 0.2075 0.175 -0.3785 -0.9256 0.2572 0.108 MBC 0.9087 -0.4174 0.2777 0.093 0.3761 -0.9266 0.5598 0.004** Notes: AP, Available phosphorus; AK, Available potassium; NH4 +-N, Ammoniacal nitrogen; NO3 --N, Nitrate nitrogen; MBC, Microbial biomass carbon. *, p ≤ 0.05, **, p ≤ 0.01. page 7 of 15Zoological Studies 63: 5 (2024) © 2024 Academia Sinica, Taiwan Fig. 4. Relative abundance (a) and cluster heatmap (b) of bacterial functional genes in different treatments based on pathway 3 of the KEGG database, and RDA analysis of soil physicochemical properties and bacterial functional genes (c). Different lowercase letters indicate significant differences (Tukey’s test, p < 0.05). Table 4. Importance of soil physicochemical properties in explaining the relative abundance variances of bacterial KEGG functional genes obtained by redundancy analysis Soil physicochemical properties RDA1 RDA2 R2Pr (> r) pH -0.4590 0.8885 0.3329 0.049* AP 0.5852 -0.8109 0.0396 0.770 AK -0.4711 0.8821 0.6076 0.001** NH4 +-N 0.1297 -0.9916 0.0345 0.793 NO3 --N 0.4992 0.8665 0.4026 0.020** MBC 0.9475 -0.3198 0.1908 0.203 Notes: *, p ≤ 0.05, **, p ≤ 0.01. (a) (b) (c) page 8 of 15Zoological Studies 63: 5 (2024) © 2024 Academia Sinica, Taiwan Fig. 5. Analysis of the niche width (BCom) of bacteria (a) and nematodes (b) (at the genus level) in each treatment based on correlation and optimal multiple regression model, Mantel test of community diversity and niche breadth of soil bacteria and nematodes (c), and RDA analysis of soil bacterial-feeding nematode and bacterial genera (d). Table 5. Importance of soil bacterivores in explaining the bacterial genus variances obtained by redundancy analysis Bacterivores RDA1 RDA2 R2Pr (> r) Panagrolaimus 0.9552 -0.2961 0.3596 0.018* Acrobeles -0.5157 0.8568 0.4198 0.025* Acrobeloides -0.3809 -0.9246 0.0095 0.954 Pristionchus 0.0432 -0.9991 0.0740 0.498 Prismatolaimus -0.4506 0.8927 0.2308 0.161 Eumonhystera -0.2083 -0.9781 0.1283 0.342 Notes: *, p ≤ 0.05, **, p ≤ 0.01. (a) (b) (c) (d) page 9 of 15Zoological Studies 63: 5 (2024)