The experiment was conducted in an experimental base within Wuming County (22°59′58″N, 107°49′26″E), Guangxi, China, in 2019 and 2020. The climate there is a subtropical monsoon climate, with an average annual precipitation of 1,233 mm and an average annual temperature of 21.7 °C (Figure 1). The field was left fallow for 2 years before the experiment was conducted, and its soil physical and chemical characteristics were well described in Li et al. (2023).

A split plot design was adopted for the field experiment (Supplementary Figure S2). Tillage treatment was designated as the main factor, and fertilization regimes were set as the secondary factor. Two tillage methods, CT and FL, were applied, and two blocks were first split along the long dimension of the experimental field for CT and FL, respectively. To ensure stable and high-yielding cassava production with improved fertilizer use efficiency, a balanced fertilization experiment was also included. The recommended NPK application rate for cassava, as proposed by local research institutions, was N 358.80 kg·ha⁻¹, P₂O₅ 89.10 kg·ha⁻¹, and K₂O 187.50 kg·ha⁻¹ (Gao et al., 2011), representing a nutrient ratio of approximately N: P₂O₅: K₂O ≈ 4: 1: 2. This recommended rate was defined as the 100% NPK treatment (100 NPK). Based on this benchmark, three reduced fertilization levels were established: 33% NPK (33NPK), 66% NPK (66NPK), and a zero-fertilizer control (0 N0P0K, CK). To address the negative consequences of nitrogen over-application on the environment and soil, reduced-nitrogen treatments were incorporated into the experimental design. In total, eight fertilization rates were then set, replicated three times, and randomly arranged within the tillage treatment blocks, resulting in a total of 48 plots (Table 1). The individual plot size was 14 m × 5 m, ending up with 70 plants in each plot with a row spacing of 1 m × 1 m.

Cassava variety South China 205 was used as plant material, which was one of the most popular cultivated varieties in China. Tillage was implemented for field preparation on April 15, 2019. Given the high soil disturbance associated with FL tillage, it was commonly practiced on an alternate-year basis. Consequently, on April 9, 2020, all plots underwent only two passes of raking to a depth of 18 cm. Cassava was sown on the date after field preparation and harvested on January 13th 2020 and January 17th 2021, respectively. No base fertilizer was used, and the designed fertilizers were applied as topdressing, 60% were applied on the 60 days of planting (DAP) and the rest were applied on the 120 DAP. No artificial watering was added. The planting method, fertilization amount and field management in 2020 were all consistent with those in 2019.

Five plants were randomly selected from each plot. The plant height and the stem diameter (the diameter of the cassava stem at about 15 cm above the ground) of these plants were tracked and measured. The plant material sampling was conducted during the cassava root formation period (July 14), expansion period (September 18), maturity period (December 15), and harvest time (January 13, 2020) in 2019, while the samples for 2020 were collected on July 11, September 17, December 18, 2020, and January 17, 2021, respectively. The fresh weight (FW) of cassava tuberous roots and stems were measured and transported to the laboratory, cleaned, chopped into small fractions (<2cm³), and thoroughly homogenized. Quartering method was then applied, and about 500 g of fresh sample was taken for each plot, oven-dried at 105 °C for 15 min, dried to constant weight at 60 °C, ground to pass through an 80-mesh sieve, and then preserved in plastic bags before wet chemical analysis. The root and stem yield of cassava were calculated using the same formula: (1) yield (kg.hm⁻²) = 10,000 (plant population/hm²) × FW (kg). To evaluate the leaf yield, 10 cassava plants were randomly selected for each plot. The average leaf number of each plant (ALN) was estimated by calculating the remain green leaves and the leaf scars. Ten mature green leaves were taken and oven-dried to calculate the average weight of a single leaf (ASL). Leaf dry matter yield (LDY) was calculated using formula (2) LDY (kg.hm⁻²) = 10,000 (plant population/hm²) × ALN× ASL (kg).

Plant samples were digested using H₂SO₄-H₂O₂ method (Bao, 2000). Contents of nitrogen and phosphorus were then determined by using the automatic intermittent chemical analyser (SmartChem 200, Zeal Quest Equipments), while potassium content was measured using a flame photometer (WGH6400, Shanghai Changxi Instrument & meters Co., Ltd.). The nutrient accumulation (NA) of N, P, and K of each plot was calculated by using formula (3) NA (kg.hm⁻²) = root and stem dry-matter yield (kg.hm⁻²) × nutrient content (%). The FUE of each element was then calculated accord to formula (4) (Fageria and Baligar, 2003):

where N_f and N_u are the nutrient accumulation of fertilized plots (kg) and unfertilized plots (kg), respectively, while N_a is the quantity of nutrient applied (kg).

A comprehensive evaluation of each treatment in terms of cassava yield and FUE was conducted using principal component analysis (PCA). The model included only principal components with eigenvalues greater than 1. The composite score for all fertilizer treatments under FL tillage was calculated using the equation as follows (Zhang et al., 2023):

J = 1,2,3,4,5,6.

In the equation, T1 _j and T2 _j represent the scores of each variable on principal component 1 and 2, respectively, derived from the PCA component score matrix. X _j denotes the original measured value of the yield of root and stem, and FUE in 2019 to2020. The coefficients 0.46196 and 0.23947 are the variance percentages explained by component 1 and component 2, respectively, which correspond to their initial eigenvalues of 2.772 and 1.437.

To evaluate how tillage treatments influence the soil bulk density and porosity, five sampling points of each tillage treatment were randomly chosen according to the S-shaped route (Xiao et al., 2020) and the samples of soil profile with three layers of 0 ~ 10 cm, 11 ~ 20 cm, and 21 ~ 30 cm of each point were collected before planting and harvesting in 2019 and 2020, respectively. Soil samples for all layers were collected by ring knives, marked, brought back to the laboratory, and measured (Li et al., 2022).

To measure the chemical properties, soil samples were collected on July 14 (the cassava root formation stage), September 18 (root expansion stage), and December 15 (root maturity stage) in 2019 and 2020, respectively. Five sampling points were randomly selected with an S-shaped pattern within each plot. The impurities on the soil surface were removed before sampling. Soils from the 0 ~ 20 cm depth at these five points were collected using a soil sampler and thoroughly mixed. The quartering method was subsequently applied to obtain a representative soil sample of approximately 500 g per plot, and then air-dried. To measure the soil microorganisms, soil samples collected during the root maturity stage of 0 N, 25 N, 50 N, and 100 N treatments were immediately divided into two portions after mixing and quartering. One portion (100 g) was immediately frozen with liquid nitrogen and stored at −80 °C for the determination of soil microorganisms, while the other part was air-dried for soil chemical properties detection.

The pH value was detected using a PHS-2F meter (Shanghai Inesa Scientific Instrument Co., Ltd.) with a soil-water ratio of 1: 2.5 (W/V). The content of organic matter (OM) was measured by the potassium dichromate volumetric method. Alkali-hydrolysable nitrogen (AN) was measured by alkaline hydrolysis diffusion method, while available phosphorus (AP) was measured by molybdenum acid colorimetry, and available potassium (AK) was measured by NH₄OAc extraction and a flame photometer (WGH6400, Shanghai Changxi Instrument & Meters Co., Ltd.). Details of these methods were well described by Bao (2000).

The DNA of bacteria and fungi in soil samples was extracted using the FastDNA® Spin Kit for Soil (MP Biomedicals, USA) and the E.Z.N.A.® Soil DNA Kit (Omega Bio-Tek, Norcross, GA, USA) according to the manufacturer’s instructions, respectively. For bacteria, the V3-V4 hypervariable regions of the 16S rRNA gene were targeted using primers 338F (5’-ACTCCTACGGGAGGCAGCAG-3′) and 806R (5′-GGACTACHVGGGTWTCTAAT-3′), while primers SSU0817F (5′-TTAGCATGGAATAATRRAATAGGA-3′) and 1196R (5′-TCTGGACCTGGTGAGTTTCC-3′) were used to amplify the target fragment of the 18S rRNA sequence in the SSU0817F_1196R region for fungi (Borneman and Hartin, 2000). The DNA was amplified by PCR according to the method described Li et al. (2023). The PCR products were separated using a 2% agarose gel, then purified and quantified using the AxyPrep DNA Gel Extraction Kit (Axygen Biosciences, Union City, CA, USA) and the Quantus™ Fluorometer (Promega, USA), respectively. Library was constructed using the NEXTFLEX Rapid DNA-Seq Kit. Paired-end sequencing (2 × 300 bp) was conducted using the Illumina MiSeq platform (Illumina, San Diego, CA, USA) by Majorbio Bio-Pharm Technology Co., Ltd. (Shanghai, China). All raw reads of bacteria and fungi obtained in this study have been deposited in the National Center for Biotechnology Information (NCBI) with accession numbers PRJNA1188335 and PRJNA1188458, respectively.

Quality filtering and adapter trimming of the demultiplexed reads were performed using Trimmomatic software. Flash software was used to filter and splice the data with the following criteria: (1) reads were trimmed when the average quality score fell below 20 within a 50-bp sliding window; (2) sequences were merged if they overlapped by more than 10 bp with no more than 2 mismatches; (3) sequences of each sample were separated according to barcodes (exactly matching) and primers (allowing two nucleotide mismatching), and reads containing ambiguous bases were removed. Subsequently, the Upaste algorithm (Upaste v7.1) was used to cluster the validated tags across all samples, with sequences being grouped into operational taxonomic units (OTUs) under the standard 97% similarity criterion. The RDP classifier 2.11 was used to annotate the screened non-chimeric sequences for species classification, using a threshold of 0.7 against the Silva database (SSU138 for fungi, Release_138 for bacteria, http://www.arb-silva.de/). Alpha diversity analysis (including the Shannon, Simpson, Ace, and Chao indexes) was carried out using MOTHUR software. In addition, the R 3.3.1 was used for data mining and plotting.

Statistical analysis, such as analysis of variance (ANOVA), paired t-test, multiple comparison with the Tukey HSD, and PCA modeling was performed using IBM SPSS Statistics 20.0, while plotting was done using OriginPro 2022. The structural equation model was generated by the SPSSPRO online platform (https://www.spsspro.com/analysis/index, accessed on 15 August 2025).

Multiple comparison of soil physical properties revealed that there were differences in soil bulk density and porosity in all layers between the two tillage treatments in 2019 after the first cassava planting season (Table 2). Soil bulk density of the 11 ~ 30 cm layers (1.37 ~ 1.53 g‧cm⁻³) in 2019 and the 21 ~ 30 cm (1.67 g‧cm⁻³) layer in 2020 under the FL treatment were significantly lower than that of the CT group (1.50 ~ 1.74 g‧cm⁻³ of 11 ~ 30 cm and 1.81 g‧cm⁻³ of 21 ~ 30 cm), while the difference gap of the 21 ~ 30 cm layer between two tillage treatments became narrower after the secon cassava planting season (2019, 1.53 g·cm⁻³ of FL, 1.74 g·cm⁻³ of CT; 2020, 1.67 g·cm⁻³ of FL, 1.81 g·cm⁻³ of CT). Soil porosity of the FL treatment (2019, 45.25% ~ 47.35%; 2020, 32.12% ~ 43.94%) was generally higher than that of the CT group (2019, 34.22% ~ 43.53%; 2020, 31.47% ~ 39.68%). In addition, soil chemical properties were influenced by both tillage and fertilizer treatments. Notably, tillage showed stronger influences on the nutrient parameters than the fertilizer treatment according to the mean square of ANOVA results, especially in 2019 (Table 3, Supplementary Figure S3).

ANOVA and paired t-tests of cassava agronomic traits revealed that FL had a significant promoting effect on cassava plant height and stem diameter in 2019, and this effect persisted until 2020 (Supplementary Figures S5A,B). In the CT group, higher plant height and larger stem diameter were found in the 75N, 100N, 33NPK, and 66NPK treatments than other treatments in 2019, indicating that N did promote cassava growth, especially when equilibrium regimens of N, P, and K applied (Supplementary Figures S5A,B). Among all fertilized treatments in the FL group, the growth-promoting effects of FL on cassava agronomic traits in 2019 was stronger than fertilizer and thus fertilizer showed little effect on cassava. In 2020, equilibrium fertilization regimens 66NPK and 100N in CT group showed stronger positive effects on cassava plant height and stem diameter than other fertilization rates (Supplementary Figures S5C,D). In the meanwhile, plant height and stem diameter of CK and 25N treatments of the FL group were significantly lower than other higher N input treatments, indicating that N from the first year accumulated in 2020 and showed stronger effects on cassava than the first year (Supplementary Figure S5).

Tillage treatment also dominantly influenced cassava yields in 2019 and 2020 (Figure 5). In 2019, the stem yields of all the fertilized treatments within the CT group were all higher than CK, while the root and total yields of all fertilized treatments within the FL group were higher than CK. In 2020, the root yield of cassava was significantly influenced by the interaction of tillage treatment and fertilizer regime. The root and yields were increased with the inputs of N within the CT group. During the whole study period, cassava yields under FL were significantly higher than under CT by 9.04 to 135.81%, according to the paired t-test results. However, there were no significant difference of yields among all fertilization rates within the FL treatment, suggesting that it was possible to maintain the cassava yield with low fertilizer input when the FL tillage was applied.

Tillage dominantly influenced the nutrients accumulations of cassava in 2019 and 2020, while fertilization showed significant effects on the N accumulation in 2020 (Figure 6). Compared to the CT group, FL significantly increased the accumulation of N, P, and K in cassava by 11.51% ~ 62.07, 9.86% ~ 33.51, and 21.75% ~ 40.76%, respectively, in 2019 (Figures 6A–C), while the nitrogen accumulation was found increased with N application rates among all fertilized treatments within the CT or FL groups in 2020 (Figure 6D). The interaction between tillage methods and fertilization treatments significantly affected the K accumulation in cassava in 2020, which increased with N application rates within the CT group (Figure 6F). Variations of N and K accumulations in 2020 indicated that FL showed positive effects on cassava yield, leading to a higher nutrient accumulation with lower N inputs than that of CT (Figures 5F, 6D,F).

According to ANOVA results, FUE of N, P, and K of cassava under different treatments were dominantly influenced by tillage treatment, while fertilizer treatment showed significant effects on FUE of N and K in 2019 and FUE of P in 2020 (Figure 7). Paired t-test indicated that FUE of N (2019, t = −3.98, p < 0.01; 2020, t = −3.25, p < 0.01) and K (2019, t = −5.52, p < 0.001; 2020, t = −2.97, p < 0.01) of cassava of the FL group were generally higher than that of CT group within the whole study period. The FUE of N and P decreased with the N inputs. However, when the application amounts of nitrogen, phosphorus, and potassium fertilizers were balanced (33NPK and 66NPK treatment areas), the FUE of nitrogen increased by 19.49 to 117.32%, indicating that a reasonable fertilizer ratio can enhance the nitrogen utilization rate by cassava. Furthermore, in FL group in 2019, FUE of P in the balanced fertilizer treatment plots (33NPK and 66NPK) were 19.42 and 34.53%, respectively, which were higher than other fertilizer treatments (Figures 7B,C).

The results of structural equation modeling indicated a clear mechanistic pathway for the year 2019 (Figure 8A): FL improved soil porosity (β = 1.00, p < 0.001), which subsequently increased the relative abundance of aerobic bacteria (β = 0.94, p < 0.001), thus facilitating the mineralization of soil organic matter. These processes ultimately resulted in a significant increase in cassava yield (β = 2.37, p < 0.01) and FUE (β = 2.50, p < 0.05). Compared with the 2019 model, the 2020 structural equation model showed that the pathways from FL to yield and FUE were no longer significant (Figure 8B). This mechanistic shift can be attributed to the concurrent decline in the effect of FL and the increase in the effect of fertilizer application. Soil physical properties under FL in 2020 further explained this change: bulk density increased (2019, 1.37 ~ 1.53 g·cm⁻³; 2020, 1.55 ~ 1.67 g·cm⁻³) while porosity decreased (2019, 45.25% ~ 47.35%; 2020, 32.12% ~ 43.94%), narrowing the difference in soil aeration between FL and CT (Table 2). Moreover, ANOVA and paired t-test results on soil chemical properties in 2020 indicated that tillage methods had no significant influence on soil OM (2019, mean square, MS = 74.15**, t = −3.720***; 2020, MS = 3.312^ns, t = −1.041^ns), AN (2019, MS = 187.5**, t = −3.075**; 2020, MS = 454.8^ns, t = −1.206^ns), AP (2019, MS = 4973***, t = −9.053***; 2020, MS = 0.001^ns, t = 0.002^ns), P accumulation (2019, MS = 842*, t = −4.43***; 2020, MS = 12.36^ns, t = 0.678^ns), or FUE-P (2019, MS = 0.241*, t = −3.84**; 2020, MS = 0.0002^ns, t = 0.17^ns) (Table 3, Figures 6, 7). Together, the diminished influence of FL on soil properties and the enhanced effect of fertilizer application in 2020 ultimately led to the non-significance of the pathways from FL to yield and FUE. Furthermore, a PCA model was developed to identify the optimal fertilization regime under FL tillage. Composite scores were calculated for all fertilizer treatments, with the 33NPK treatment achieving the highest score (Table 4), thereby establishing it as the optimal fertilization regime under the FL system.

Tillage can directly change the physical properties of soil, thus affecting chemical properties and microbial communities, creating favorable environments for crop growth (Ahmadi et al., 2024; Adeleke et al., 2023; Chichongue et al., 2020). Applying optimum tillage techniques combined with fertilizer application does not only improve soil quality but also results in higher crop production (Vilakazi et al., 2022). In the current study, the FL treatment showed lower soil bulk density and higher porosity than the CT plots, especially at the 21 ~ 30 cm layer (Table 2). This finding was consistent with previous researches (Zhai et al., 2017, 2019). Our previous study showed that the clay and the silt of cassava rhizosphere soil were also influenced by FL, fertilization, and their interaction (Li et al., 2023). For bulk soil chemical properties, tillage showed stronger effects on nutrient parameters than the fertilizer treatment, and the contents of OM, AN, AP, and AK of soil samples under the FL treatment were significantly higher than that of the CT treatment (Table 3). The same phenomenon was also found in cassava rhizosphere soils under FL and CT treatments, and tillage treatment also showed stronger influence on the chemical properties of cassava rhizosphere soils than N applications (Li et al., 2023). According to the climate data (Figure 1), waterlogging occurred during the summer time, especially in June of 2020. A lower bulk density and higher porosity of the bulk soil under FL treatment enhanced the infiltration of the rainfall, thus increasing the resistance to waterlogging stress of cassava. The aerobic conditions in soil created by FL can also accelerate the oxidation of reduced N compounds to nitrate, which is susceptible to leaching out of the soil due to its high solubility. Loss of nitrate as coupled ions promoted cation (Ca²⁺ and Mg²⁺) leaching, which can accelerate the soil acidification in FL treatment, ending up with a lower pH value than CT in 2020 (Table 3, Supplementary Figure S3) (Zhou et al., 2024). However, the high solubility of nitrate in FL treatment would be also benefit the N absorbency for cassava under drought stress, leading to a higher N accumulation in cassava than that of CT (Figure 6).

Soil microbiomes encompass a great diversity of organisms that perform vital ecological functions, including nitrogen and carbon cycling (García-Serquén et al., 2024), and their diversity positively correlates with multiple ecosystem functions such as plant productivity, nutrient cycling, and decomposition (Wu et al., 2024). In the current study, soil physical and chemical factors varied with both tillage and the fertilization regimes (Tables 2, 3), leading to a shifting of the soil microbe structure (Figure 2), which was consistent with previous reports (Chi et al., 2024; Özbolat et al., 2023; Ouverson et al., 2021). Since the field had been abandoned for 2 years before the experiment, the initial soil microbial species were rich and diverse. The application of nitrogen fertilizer reduced the Shannon, Ace, and Chao indices while increasing the Simpson index of soil microbes (Supplementary Table S2). This suggests that nitrogen input altered soil properties such as pH, OM, and AN content, ultimately leading to decreased richness and diversity in both bacterial and fungal communities (Table 3, Supplementary Figure S3).

Bacteria and fungi are key decomposers in the soil ecosystem, directly driving the carbon cycle through the breakdown of organic residues (Kramer and Gleixner, 2008; Strickland et al., 2009). Their necromass subsequently constitutes a primary source of stable soil organic matter (Bardgett and van der Putten, 2014). The dominant phyla of the bulk soil bacteria among all treatments in 2019 and 2020 were Actinobacteria, Proteobacteria, Chloroflexi, Acidobacteria, and Firmicutes (Figures 2A,C), which was similar to previous results (Jia et al., 2024; Nimnoi et al., 2024). Actinobacteria and Proteobacteria are considered as key roles in the decomposition of organic matter in soil, while organic C provides energy resource for microbial growth in soil (Zhao et al., 2019; Xiao et al., 2016). Under both FL and CT treatments in 2019, the relative abundance of Actinobacteria and Proteobacteria increased with the amount of nitrogen applied. In 2020, however, the relative abundance of Actinobacteria in all fertilized treatments was consistently higher than that in 2019 (Figure 2). This pattern aligns with previous findings indicating that the richness of Actinobacteria is negatively correlated with soil pH and can serve as an indicator of soil acidification induced by long-term nitrogen fertilization (Ren et al., 2020). Proteobacteria, a key component of soil bacteria that includes nitrogen-fixing and pathogenic species, showed a decrease in relative abundance in 2020. This decline may be attributed to the concurrent increase in Actinobacteria, which likely outcompeted Proteobacteria for carbon resources to support their own population growth.

Furthermore, FL treatment consistently maintained higher soil porosity than CT (Table 2), creating a more favorable aeration environment for aerobic bacteria. Correspondingly, it promoted the enrichment of several predominantly aerobic taxa, including OTU2000 (Gemmatimonadaceae, Family) and OTU1279 (Tumebacillus, Genus) in 2019, and OTU5900 (Micromonosporaceae, Family), OTU4505 (Tumebacillus, Genus), OTU5623 (Conexibacter, Genus), and OTU5940 (Acidobacteriales, Family) in 2020 (Figures 3A,C). Previous studies have confirmed the aerobic nature of these taxa (Ivanitskaia et al., 1982; Lei et al., 2023; Mujakic et al., 2023; Zhang et al., 2023). Aerobic bacteria contribute significantly to soil carbon cycling (Angélica et al., 2024), which likely accelerated nutrient release and led to higher contents of AN, AP, and AK under FL treatment compared to CT (Table 3). This suggests that FL tillage enhances the soil C/N cycling pathway—a key mechanism for maintaining soil ecosystem stability (Gao et al., 2025). Notably, OTU5623, identified as a strictly aerobic Conexibacter, can reduce nitrate to nitrite and may play an important role in nitrification (Pukall et al., 2010; Monciardini et al., 2003). Its enrichment under FL treatment could further support the regulation of soil nitrogen cycling.

Fungi regulate the balance of carbon and nutrients in the soil ecosystem by metabolizing and decomposing complex organic matter (Žifčáková et al., 2016). Ascomycota and Mucoromycota were dominant phyla in the bulk soil fungal communities among all treatments in current study (Figures 2B,D), which can decompose refractory organic matter in soil and play important roles in nutrient cycling (Zhao et al., 2023; Das et al., 2019). Ascomycota growth is primarily regulated by soil nitrogen content, with high N levels generally suppressing its abundance (Zhang et al., 2024; Fontaine et al., 2011). Consistent with this, our study observed that high N inputs significantly inhibited Ascomycota growth in the CT group throughout the experimental period (Figures 2B,D). However, the relative abundance of Ascomycota remained relatively stable across fertilized treatments within the FL group. This stability may be attributed to the ability of FL to improve bulk soil physicochemical properties, particularly by maintaining higher AN content, thereby helping to preserve the structure of the fungal community. Furthermore, FL tillage simultaneously created more stable soil conditions for Mucoromycota. These fungi play a key role in soil carbon cycling due to their strong capacity for decomposing OM—a function essential for maintaining soil fertility and microbial diversity. As important regulators of nutrient and carbon cycles, Mucoromycota significantly influence terrestrial ecosystems by shaping soil structure and function (Van der Heijden et al., 2015). In summary, compared to CT, FL treatment promoted greater stability in the soil microbial structure. This enhanced stability improved the ecosystem’s resistance to external disturbances, such as nitrogen stress, supporting FL as a beneficial practice for sustaining soil ecological functions.

As we mentioned before, tillage directly changes the soil physical properties, while soil microbes directly and indirectly influence crop productivity and its nutrient acquisition ability by changing the nutrients or stimulating the plant growth (Wang et al., 2021). An increase of nitrogen fertilizer application will accelerate soil acidification, which could enhance the activity of ammonifying bacteria leading to an increase of the volatilization of ammonia nitrogen (Dancer et al., 1973). Additionally, the population of nitrifying bacteria would be influenced under the low pH environment, leading to a decrease of nitrate formation rate, thereby restricting plant’s nitrogen acquisition ability (Saratchandra, 1978). Thus, crop yield and N accumulation increased with the amount of applied nitrogen, while the FUE decreases. FL tillage significantly improved soil conditions for both cassava root development (Figure 5) and microbial activity, particularly that of aerobic microorganisms (Figure 3), by thoroughly loosening the soil and increasing its porosity. This effect was especially pronounced in the 11–30 cm soil layer (Table 2). The resulting soil structure facilitated deeper root penetration compared to CT, enabling more efficient uptake of nutrients and water from subsurface soil layers (Xiao et al., 2023). The well-developed cassava root system under FL treatment was conducive to recruiting rhizosphere microorganisms, which in turn affected soil enzyme activities (Huang et al., 2025; Li et al., 2023). It also reduced the emission of the greenhouse gases, such as CO₂, N₂O and CH₄, and thus improved the fixation abilities of nitrogen and carbon in soil (Xiao et al., 2023; Zhu B. et al., 2023, Zhu S. et al., 2023; Duan et al., 2022).

In addition, the contents of AN, AP, and AK in soil under FL treatment were higher than that of CT (Table 3), offering more nutrients for cassava leading to a higher yield and nutrient accumulation (Figures 5, 6), which was consistent with studies for rice (Zhang et al., 2021), sugarcane (Xiao et al., 2023), and tobacco (Nicotiana tabacum L.) (Zhu B. et al., 2023; Zhu S. et al., 2023). The FUE-N and FUE-K of cassava within the FL group were also generally higher than that of CT group (Figure 7), indicating that FL enhances the utilization of nitrogen and potassium fertilizers by cassava.

Furthermore, identifying the optimal fertilization rate under the FL tillage system is critical for sustainable cassava production. Under FL treatment, the highest yield and FUE-K were all found in the 33NPK plot in 2019, while its FUE-N was lower than the 25NPK treatment but slightly higher than others treatments (Figure 7, Table 4). These superior performance of the 33NPK treatment of FL group lasted to 2020, which indicate that FL with 33NPK treatment strikes an ideal balance between cassava yield and FUE among all treatments. Although we did not identify specific functional genes through metagenomics, a significant influence of tillage practices on soil enzyme activities has been demonstrated in our previous work, particularly the urease activity was found notably promoted under FL tillage (Huang et al., 2025). Considering the effects on soil physicochemical properties and microbial communities caused by FL (Figure 8), structural equation models confermed that FL can accelerates the mineralization process of OM by increasing the abundance of aerobic bacteria, thus releases the nutrients from the soil for cassava growth, which ultimately leads to concurrent improvements in both yield and FUE of cassava.