1 Introduction

The global population is projected to be 8.5 billion people in the 2030s (Forman and Wu 2016). Large amounts of phosphorus (P), a nonrenewable resource, are derived from phosphate rock reserves to meet the growing demand for food. Global reserves of phosphate rock will be depleted within several hundred years with the increase in the use of P, especially P fertilizers (Fixen and Johnston 2012; Johnston et al. 2014). According to the statistics of FAO in 2021, the consumption rate of P fertilizer has been higher than the growth rate of global grain production, resulting in low P-utilization efficiency (FAOSTAT 2021). Most crops only use 15-20% of the applied P fertilizer annually (Zhang et al. 2008). Excessive P fertilizers will account for lower P-use efficiency, leading to eutrophication of freshwater and marine ecosystems through leaching and erosion (Vorosmarty et al. 2010; Elser and Bennett 2011). In recent years, an increasing number of studies have explored methods to reduce P-fertilizer consumption and improve P-use efficiency through developing novel fertilizers, microbiome inoculations, and breeding desirable genotypes and cropping systems (Cong et al. 2020). Intercropping is an economical and environmentally-friendly approach to achieve the goal by increasing crop diversity (Chien et al. 2009; Jing et al. 2022).

Intercropping refers to cultivation practices where two or more crop species are planted in the same field during at least part of the same growing period (Vandermeer 1992), and has been repeatedly reported to enhance productivity, soil fertility, nutrient-use efficiency, and temporal stability (Gong et al. 2020; Li et al. 2021b; Wang et al. 2021b; Yu et al. 2022). Based on a meta-analysis (Tang et al. 2021) and case studies (Eichler-Loebermann et al. 2021; Song et al. 2021), intercropping also increased P uptake compared with monocultures. The productivity and P-uptake advantages through crop diversification can be explained by two processes: the selection effect and the complementarity effect (Loreau and Hector 2001). The selection effect (SE) in P uptake is defined as greater P uptake in intercropping by the dominant crop species. The complementarity effect (CE) in P uptake indicates greater P uptake of both components of a mixture relative to that of monocultures (Li et al. 2018; Clark et al. 2019). The positive CE of P uptake in crop diversification can be divided into resource partitioning (Hinsinger et al. 2011) and facilitation (Li et al. 2007; Li et al. 2014; Dissanayaka et al. 2015). Facilitation occurs between two or more crop species in mixture with a beneficial outcome (e.g., increase biomass) for one or more crop species, which plays an important role in sustainable intercropping systems and has received increasing attention as an important mechanism driving P-uptake advantages (Wright et al. 2017; Brooker et al. 2021; Yu et al. 2021a). However, several studies indicate that most intercropping systems increases P uptake, but not all intercropping combinations do (Li et al. 2019; Gitari et al. 2020; Tang et al. 2021). Therefore, we surmise the greater performance of intercropping requires desirable crop combinations mediated by functional traits of intercropped species (Yu et al. 2021a).

Plant species are diverse in their ability to mobilize P (Lyu et al. 2016; Lambers 2022). Efficient P-mobilizing species are defined as those with a large capacity to mobilize P, releasing exudates (phosphatases, carboxylates, and/or protons), as so-called physiological root traits to convert unavailable P into plant-available P in the rhizosphere. In contrast, inefficient species exhibit weak physiological root traits to mobilize sorbed and organic soil P (Lambers et al. 2006). As an efficient P-mobilizing species, chickpea is more efficient at mobilizing organic P into available forms by releasing phosphatases; and faba bean releases carboxylates and protons to mobilize insoluble P into dissolved phosphate (Li et al. 2004; Li et al. 2007). As inefficient P-mobilizing species, maize and wheat usually exhibit weaker physiological root traits to mobilize sorbed P; instead, these plant species exhibit greater specific root length and root branching intensity than efficient P-mobilizing species, i.e. root morphological traits, to acquire soil P (Wen et al. 2019).

When plant species with different P-acquisition traits are grown together, the advantage of P uptake may occur through niche differentiation and facilitation (Yu et al. 2020a; Navarro-Cano et al. 2021). Intercropping efficient P-mobilizing species with inefficient species enhances the biological potential to mobilize sorbed soil P, which may also save P-fertilizer inputs (Li et al. 2014). However, most studies focused on how efficient P-mobilizing species facilitate species that are inefficient at P mobilization, whereas few explored how inefficient P-mobilizing species respond to their neighbors. Furthermore, most experiments were conducted in the greenhouse to explore belowground processes and P uptake in plant species grown together (e.g., Li et al. 2004, 2007). The links among root-root interactions, P uptake, and P-fertilizer saving have largely been underexplored under field conditions.

To explore the effects of cropping system (i.e. monoculture and intercropping) and P-application rate on P uptake, a field experiment with nine cropping systems (faba bean/maize, soybean/maize, chickpea/maize, oilseed rape/maize, and corresponding five monocultures) and three P-fertilizer application rates (0 kg P ha−1, 40 kg P ha−1, and 80 kg P ha−1) was established in 2009 (see Baiyun-2 Experiment in Li et al. 2021b; Fig. 1a, b). Our previous study calculated the apparent recovery efficiency of applied P to reflect the utilization efficiency of P fertilizer. In the early stage of the field experiment (the second to the third year after establishment), we found that intercropping had a greater apparent recovery efficiency of applied P and shoot P content in the four intercropping systems than those of monocultures (Xia et al. 2013). Agronomic efficiency of applied P and internal utilization efficiency of P are also used to reflect P-use efficiency (Adiele et al. 2020; Bayuelo-Jimenez and Ochoa-Cadavid 2014); while how root-root interactions mediate P-use efficiency and sustainable P management in intercropping systems remains largely unknown under field conditions.

Fig. 1
figure 1

An aerial view of experimental plots at Baiyun Experiment Station in Wuwei city, Gansu province (a). The four intercropping combinations in the growing season (June), including faba bean/maize intercropping system, chickpea/maize intercropping system, oilseed rape/maize intercropping system, and soybean/maize intercropping system (b). Conceptual diagram of the main topic of the article. We explored the way root traits affect productivity, shoot phosphorus (P) content and P-use efficiency indicators in intercropping with different P-application rates. We also link changes in productivity, shoot P content, and P-use efficiency indicators with sustainable P management in agroecosystems (c).

Full size image

In the present study, we used data in the 11th (2020) and 12th (2021) year in the long-term intercropping experiment (Fig. 1) to explore the effect of intercropping on productivity, aboveground P content, and P-use efficiency (i.e. apparent recovery efficiency of applied P, agronomic efficiency of applied P, and internal utilization efficiency of P). We also measured physiological and morphological traits, namely: leaf manganese concentrations (leaf [Mn]), as a proxy for rhizosheath carboxylate concentrations (Lambers et al. 2021; Lambers et al. 2015); acid phosphatase activity; alkaline phosphatase activity; pH; soil microbial biomass P, as an indicator for microbial activity; root surface area, specific root length, root tissue density, root branching intensity, average root diameter and proportion of fine roots. We address two key questions: (1) How do long-term intercropping and P-application rates affect productivity, P-use efficiency, and shoot P content of intercropping systems compared with monocultures? (2) How do root traits affect P-use efficiency through belowground interspecific interactions under field conditions (Fig. 1c)?

2 Materials and methods

2.1 Site, climate, and description

The sampling was carried out in 2020 and 2021 in a long-term field experiment, which was established in 2009, at Baiyun Experiment Station in Wuwei city, Gansu province (38°37′N, 102°40′E, at an altitude of 1504 m above sea level). The site has a typical arid climate with an average annual temperature of 7.7 °C, mean annual precipitation of 150 mm and mean potential annual evaporation of 2021 mm. Total solar radiation is 5988 M J m−2 yr−1, and the frost-free period is 170-180 days. The soil type is classified as calcareous Aridisol. The initial soil (0-20 cm) had a pH of 8.0, soil organic matter content of 19.1 g kg−1, total nitrogen (N) concentration of 1.08 g kg−1, Olsen-P concentration of 20.3 mg kg−1, and plant-available potassium (K) concentration of 233 mg kg−1 before the establishment of the experiment (Xia et al. 2013). After 12 years, the soil Olsen-P concentration in the soils (0-20 cm) was 4.81 mg kg−1 without P application, 14.8 mg kg−1 with P40, and 29 mg kg−1 with P80, respectively.

2.2 Experimental design and treatments

The experiment was a split-plot design with three replicates. The main-plot factor refers to three P-fertilizer application rates: no P (P0), 40 kg P ha−1 (P40), and 80 kg P ha−1 (P80). The sub-plot factors are the identity of nine cropping systems, which included (a) sole maize (Zea mays L. cv. Xianyu No. 335), (b) sole faba bean (Vicia faba L. cv. Lincan No. 5), (c) sole chickpea (Cicer arietinum L. cv. Longying No. 1), (d) sole soybean (Glycine max L. cv. Zhonghuang No. 33), (e) sole oilseed rape (Brassica campestris L. cv. Longyou No. 5), (f) faba bean/maize intercropping, (g) chickpea/maize intercropping, (h) soybean/maize intercropping, (i) oilseed rape/maize intercropping. There were 81 experimental plots in total, and each plot experienced the same treatment for 12 years, e.g., in faba bean/maize intercropping plot, the two plant species were continuously intercropped on the same place every year.

Faba bean/maize, chickpea/maize, soybean/maize, and oilseed rape/maize intercropping systems were planted in alternating 1.4-m strips, including a 0.8-m wide maize strips and a 0.6-m wide companion crops strips, and each intercropping plot comprised four strips. Each intercropping strip included two maize rows (0.4 m inter-row and 0.25 m inter-plant distance) and three companion crop rows (0.2 m inter-row and 0.2 m inter-plant distance). Thus, maize occupied 57% and the companion crops occupied 43% of the intercropped area. Maize and companion crops monoculture plots were planted with the same inter-row and inter-plant distance in monoculture and intercropping systems. Faba bean and chickpea were sown in late March and harvested in early July. Oilseed rape was planted in late April and harvested in early July. Soybean and maize were sown in late April and harvested in late September and early October, respectively.

All P (as superphosphate) and K (60 kg K ha−1 as potassium sulfate) fertilizers were evenly broadcast before sowing. For maize and intercropping systems, one-third of N fertilizer (75 kg N ha−1 as urea) was applied as basal fertilizer, and the rest was evenly divided into two parts applied at the stem-elongation stage (75 kg N ha−1 as urea) and the pre-tasseling stage (75 kg N ha−1 as urea) of maize. In total, maize and intercropping systems received 225 kg N ha−1. The amounts of fertilization in intercropping were recommended by agronomists to meet the nutrient requirement of the maize. For monoculture legumes and oilseed rape, half of N fertilizer (75 kg N ha−1 as urea) was applied as basal fertilizer, and the other half was applied at the stem-elongation stage (75 kg N ha−1 as urea), which received 150 kg N ha−1 during the growth stage. Six irrigations were applied on mid-May, early June, late June, early July, late July, and mid-August. Each irrigation was 100 mm, both in monocrops and in intercrops.

2.3 Sample collection

Soil, mature leaf, and root samples were collected at the flowering stage of crops. The root samples were only collected in 2021, and the other samples were collected in 2020 and 2021. The samples of faba bean, chickpea, and oilseed rape were collected in mid-June, and those of soybean and maize in mid-July. Rhizosheath soil was collected from 5 to 10 individuals to ensure sufficient soil. Roots were dug out of the 0-30 cm soil layer with the area of 40×25 cm for maize and the area of 20×20 cm for companion crops which comprised most roots and ensured a consistent sampling range for the root layer. The soil attached to the roots was defined as rhizosheath soil and brushed into a collector (Pang et al. 2017). In the monoculture plots, five soil cores were collected by soil auger (diameter 5 cm) at 0-20 cm in a given plot. In the intercropping plots, five soil cores were separately collected from corresponding strips. Then soil samples were mixed as bulk soil samples. All the soil samples were passed through the 2-mm sieve and divided into two parts, then stored at 4 °C. Root samples were washed carefully to remove the remaining soil and stored at −20 °C to measure the root morphology. Mature intact leaf samples were oven-dried and ground for leaf [Mn] measurements. At maturity, three continuous rows of companion crops and two rows of maize were harvested and separated into stems and grains. The samples were air-dried to measure shoot P concentration.

Total root length, total root surface area, root volume, root length of the fine roots (diameter less than 0.2 mm; Bergmann et al. 2017), and average root diameter were measured with a scanner (Epson Perfection V750 Pro, Epson, Suwa, Japan) and analyzed with a WinRHIZO scanner-based system (WinRHIZO system; Regent Instruments Inc., Quebec City, QC, Canada). Soil pH was measured in a 1:2.5 soil: CaCl2 ratio (Seven Compact pH meter S210; METTLER TOLEDO, Switzerland, Zurich). Soil microbial biomass P was fumigated with CHCl3, and extracted by 0.5 M NaHCO3 (Brookes et al. 1982). Acid phosphatase and alkaline phosphatase activity were measured using a spectrophotometric method with p-nitrophenol (PNP) as substrate at pH 5.2 and 8.5, respectively (Tabatabai and Bremner 1969). Shoot P concentration was measured by the molybdo-vanadophosphate method (Bao 2005). Leaf [Mn] was determined using inductively coupled plasma optical-emission spectroscopy (ICP-OES; OPTIMA 3300 DV, Perkin-Elmer, Waltham, MA, USA).

2.4 Calculations

The aboveground biomass produced by intercropping can be compared with the weighted means of corresponding monoculture crops based on their proportions of the area in the intercropping systems (Li et al. 2021b). The weighted means of aboveground biomass were calculated using the following equation:

$$\mathrm{weighted}\ \mathrm{means}\ \mathrm{of}\ \mathrm{aboveground}\ \mathrm{biomass}={\mathrm{B}}_{\mathrm{m}}\times {\mathrm{Area}}_{\mathrm{m}}+{\mathrm{B}}_{\mathrm{B}}\times {\mathrm{Area}}_{\mathrm{B}}$$
(1)

where m represents maize, while B represents the companion crops (oilseed rape, chickpea, faba bean, and soybean); Bm and BB indicate the biomass of maize and companion crops, respectively. Aream and AreaB refer to the proportions of the area in the intercropping systems, where Aream = 57% and AreaB = 43%.

The above formulas were also used to calculate the weighted means of grain yield.

The weighted means of shoot P content were calculated using the following equation:

$$\mathrm{weighted}\ \mathrm{means}\ \mathrm{of}\ \mathrm{shoot}\ \mathrm{P}\ \mathrm{content}={\mathrm{P}}_{\mathrm{m}}\times {\mathrm{B}}_{\mathrm{m}}\times {\mathrm{Area}}_{\mathrm{m}}+{\mathrm{P}}_{\mathrm{B}}\times {\mathrm{B}}_{\mathrm{B}}\times {\mathrm{Area}}_{\mathrm{B}}$$
(2)

where Pm and PB are the shoot (stem + grain) P concentration of maize and companion crops.

The apparent recovery efficiency of applied P (AREP) was calculated using the following equation (Mei et al. 2012; Xia et al. 2013):

$$\mathrm{the}\ \mathrm{weighted}\ \mathrm{means}\ {\mathrm{of}\ \mathrm{AREP}}_{\mathrm{monoculture}}=\frac{\left({\mathrm{PC}}_{\mathrm{monoculture}\_\mathrm{f}}-{\mathrm{PC}}_{\mathrm{monoculture}\_0}\right)\times 100}{\mathrm{f}}$$
(3)
$${\mathrm{AREP}}_{\mathrm{intercropping}}=\frac{\left({\mathrm{PC}}_{\mathrm{intercropping}\_\mathrm{f}}-{\mathrm{PC}}_{\mathrm{monoculture}\_0}\right)\times 100}{\mathrm{f}}$$
(4)

where PCmonoculture _ f and PCmonoculture _ 0 represent the weighted means of shoot P content in the monoculture systems with P application or without P application, respectively; f refers to the amount of P fertilizer applied. PCintercropping _ f represents the P content of intercropping systems when P fertilizer was applied. Interspecific facilitation plays an important role in increased P uptake in intercropping system (Li et al. 2014); the use of PCintercropping _ 0 as a baseline value to calculate AREPintercropping probably eliminated the facilitation effect and did not reflect the effects of intercropping. Therefore, we used PCmonoculture _ 0 rather than PCintercropping _ 0 as a baseline value when calculating AREPintercropping.

The agronomic efficiency of applied P (AEP) was calculated using the following equation (Dobermann 2007):

$$\mathrm{the}\ \mathrm{weighted}\ \mathrm{means}\ \mathrm{of}\ {\mathrm{AEP}}_{\mathrm{monoculture}}=\frac{\left({\mathrm{yield}}_{\mathrm{monoculture}\_\mathrm{f}}-{\mathrm{yield}}_{\mathrm{monoculture}\_0}\right)\times 100}{\mathrm{f}}$$
(5)
$${\mathrm{AEP}}_{\mathrm{intercropping}}=\frac{\left({\mathrm{yield}}_{\mathrm{intercropping}\_\mathrm{f}}-{\mathrm{yield}}_{\mathrm{monoculture}\_0}\right)\times 100}{\mathrm{f}}$$
(6)

where yieldmonoculture _ f and yieldmonoculture _ 0 represent the weighted means of grain yield in the monoculture systems with P application or without P application, respectively; f refers to the amount of P fertilizer applied. yieldintercropping _ f represents grain yield of intercropping systems when P fertilizer was applied. We used yieldmonoculture _ 0 as a baseline value to calculate AEPintercropping which reflects the effects of interspecific interactions on agronomic efficiency of applied P of the intercropping system.

The internal utilization efficiency of P (IEP) was calculated using the following equation (Dobermann 2007):

$$\mathrm{the}\ \mathrm{weighted}\ \mathrm{means}\ \mathrm{of}\ \mathrm{IEP}=\frac{\mathrm{yield}}{{\mathrm{PC}}_{\mathrm{grain}}}$$
(7)

where PCgrain represent the weighted means of grain P content.

The total land-equivalent ratio for aboveground biomass was calculated using the following equation:

$${\mathrm{LER}}_{\mathrm{B}}=\frac{{\mathrm{B}}_{\mathrm{intercropping}\_\mathrm{m}}}{{\mathrm{B}}_{\mathrm{monoculture}\_\mathrm{m}}}\times {\mathrm{A}\mathrm{rea}}_{\mathrm{m}}+\frac{{\mathrm{B}}_{\mathrm{intercropping}\_\mathrm{B}}}{{\mathrm{B}}_{\mathrm{monoculture}\_\mathrm{B}}}\times {\mathrm{A}\mathrm{rea}}_{\mathrm{B}}$$
(8)

where Bintercropping _ m and Bintercropping _ B, and Bmonoculture _ m and Bmonoculture _ B refer to the aboveground biomass of maize and companion crops in the intercropping and monoculture systems, respectively. The above formula was also used to calculate the total land-equivalent ratio based on grain yield.

The biodiversity effect for P content was calculated using the following equation (Loreau and Hector 2001):

$$\mathrm{Complementarity}\ \mathrm{effect}\ \left(\mathrm{CE}\right)=\mathrm{n}\times \mathrm{mean}\ \left(\Delta \mathrm{RY}\right)\times \mathrm{mean}\ (M)$$
(9)
$$\mathrm{Selection}\ \mathrm{effect}\ \left(\mathrm{SE}\right)=\mathrm{n}\times \mathrm{covariance}\ \left(\Delta \mathrm{RP},M\right)$$
(10)
$$\mathrm{Net}\ \mathrm{effect}=\mathrm{Observed}\ \mathrm{P}\ \mathrm{content}-\mathrm{Expected}\ \mathrm{P}\ \mathrm{content}=\mathrm{CE}+\mathrm{SE}$$
(11)
$$\Delta \mathrm{RY}={\mathrm{RP}}_{\mathrm{O}}-{\mathrm{RP}}_{\mathrm{E}}$$
(12)

where n is the number of crop species in the intercropping, i.e. two in the present study. ΔRY indicates the deviation from expected relative shoot P content of crop species in the intercropping system. The observed relative shoot P content of crops in the intercropping system (RPO) is the shoot P content in the intercropping system divided by the shoot P content in monoculture (M). The expected relative shoot P content of crop species in the intercropping system (RPE) is the area proportion of a crop species in intercropping.

To maintain consistency with the sample size of the apparent recovery efficiency of applied P and agronomic efficiency of applied P, which only had data at the P40 and P80 levels, we processed all rest indices relative to P0 in correlation analyses and stepwise regression. The formula was as follows:

$$\mathrm{Relative}\ \mathrm{P}0\ \mathrm{data}=\frac{{\mathrm{data}}_{\mathrm{f}}-{\mathrm{data}}_0}{{\mathrm{data}}_0}$$
(13)

The indices included acid and alkaline phosphatase activity in rhizosheath soil, pH in rhizosheath soil, leaf [Mn], microbial biomass P in bulk soil, root surface area, specific root length, root tissue density, root branching intensity, average root diameter, and proportion of fine roots of maize and companion crops; we also calculated the relative P0 data of shoot P content and internal utilization efficiency of P (i.e. intercropping systems and the weighted means of corresponding monocultures) using the above formula.

Root morphological traits were calculated as follows:

$$\mathrm{specific}\ \mathrm{root}\ \mathrm{length}=\frac{\mathrm{total}\ \mathrm{root}\ \mathrm{length}\ }{\mathrm{root}\ \mathrm{dry}\ \mathrm{weight}}$$
(14)
$$\mathrm{specific}\ \mathrm{root}\ \mathrm{surface}\ \mathrm{area}=\frac{\mathrm{total}\ \mathrm{root}\ \mathrm{surface}\ }{\mathrm{root}\ \mathrm{dry}\ \mathrm{weight}}$$
(15)
$$\mathrm{the}\ \mathrm{proportion}\ \mathrm{of}\ \mathrm{fine}\ \mathrm{root}\mathrm{s}=\frac{\ \mathrm{root}\ \mathrm{length}\ \mathrm{of}\ \mathrm{fine}\ \mathrm{root}\mathrm{s}\times 100}{\mathrm{total}\ \mathrm{root}\ \mathrm{length}}$$
(16)
$$\mathrm{root}\ \mathrm{tissue}\ \mathrm{density}=\frac{\mathrm{total}\ \mathrm{dry}\ \mathrm{weight}\ }{\mathrm{root}\ \mathrm{volume}}$$
(17)
$$\mathrm{root}\ \mathrm{branching}\ \mathrm{intensity}=\frac{\ \mathrm{number}\ \mathrm{of}\ \mathrm{root}\ \mathrm{tips}\ }{\mathrm{total}\ \mathrm{root}\ \mathrm{length}}$$
(18)

2.5 Statistical analyses

Linear-mixed effect models were used in this study, using the ‘nlme’ package (Pinheiro et al. 2022). To examine the effects of P-application rate and cropping system (i.e. monoculture and intercropping) on aboveground biomass, grain yield, shoot P content, internal utilization efficiency of P, agronomic efficiency of applied P, the apparent recovery efficiency of applied P, and P-related parameters during sampling years (i.e. 2020 and 2021), we treated the factor ‘block’ as a random effect in statistical models as in previous studies (Yang et al. 2022; Yu et al. 2020b). First, aboveground biomass, grain yield, shoot P content, internal utilization efficiency of P, agronomic efficiency of applied P, and the apparent recovery efficiency of applied P of the monoculture and intercropping systems were tested using year, cropping system, P-application rate, and crop combination as fixed effects, and block as a random effect. Second, we examined the effects of year, cropping system, P-application rate, species, and their interaction effect on aboveground biomass, grain yield, shoot P content, internal utilization efficiency of P, agronomic efficiency of applied P, and the apparent recovery efficiency of applied P of each component species using year, cropping systems, P-application rate, and species as fixed factors; block was treated as a random factor. Third, the same method was applied to detect the effects of year, cropping system, P-application rate, and their interaction effect on rhizosheath pH, rhizosheath acid and alkaline phosphatase activity, leaf [Mn], MBP in bulk soil of maize, with block and neighbor identity as random factors. The root morphological traits of maize were tested with cropping systems and P-application rate as fixed effects, and block and neighbor identity as random effects. The effect of P-application rate and species on physiological traits in monoculture species of each year were tested with block as a random effect. We also tested the effects of P-application rate and cropping system on morphological root traits of companion crops, with block as a random effect. Complementarity effect calculated on the basis of P content, selection effect calculated on the basis of P content, LER based on aboveground biomass, and LER based on grain yield among crop combinations of each year were accessed with P-application rate as a fixed effect, and block as a random effect. A t-test was used to assess the difference of complementarity effect and selection effect from 0, and the difference of LER from 1. Tukey’s post-hoc HSD test was conducted at the 5% probability level in ANOVAs.

Principal component analyses (PCA) were conducted using physiological traits (rhizosheath acid phosphatase activity, rhizosheath alkaline phosphatase activity, rhizosheath pH, MBP in bulk soil, and leaf [Mn]) of the species (maize, faba bean, chickpea, soybean, and oilseed rape), using the ‘ggbiplot’ package (Vincent 2011). A PERMANOVA test was performed to calculate the P value between crop species using the ‘pairwiseAdonis’ package (Pedro 2017). Before correlation analyses and stepwise regressions, to be consistent with the sample size of the apparent recovery efficiency of applied P and agronomic efficiency of applied P of intercropping systems and the weighted means of corresponding monocultures, we turned the 11 P-related parameters associated with maize, shoot P content and internal utilization efficiency of P of intercropping systems and the weighted means of corresponding monocultures into relative P0 data (see Eq. 13). Correlation analyses were conducted to test potential correlations among P-related agronomic indicators (shoot P content, internal utilization efficiency of P, agronomic efficiency of applied P, and the apparent recovery efficiency of applied P of intercropping systems and the weighted means of corresponding monocultures) and 11 P-related parameters associated with maize, using the ‘Matrix’ package (Bates and Maechler 2010). We divided 11 P-related parameters that may affect P-efficiency into three groups, including root morphological traits, physiological traits, and microorganisms. Stepwise regressions were performed to calculate the relative contribution of maize root morphological traits (i.e. root surface area, specific root length, root tissue density, root branching intensity, average root diameter, and proportion of fine roots), physiological traits (i.e. rhizosheath acid phosphatase activity, rhizosheath pH, and leaf [Mn]), and microorganism (i.e. rhizosheath alkaline phosphatase activity and microbial biomass P in bulk soil) on the apparent recovery efficiency of applied P, agronomic efficiency of applied P, internal utilization efficiency of P, and shoot P content of intercropping systems and the weighted means of the corresponding monocultures in 2021, using the ‘hier.part’ package (Nally and Walsh 2004). Finally, we also examined the correlation between productivity (aboveground biomass and yield), P-related agronomic indicators (i.e. shoot P content, internal utilization efficiency of P, agronomic efficiency of applied P, and the apparent recovery efficiency of applied P of intercropping systems and the weighted means of the corresponding monocultures) and 11 P-related parameters associated with maize by correlation analyses at each P fertilizer level. All statistical analyses were performed with R version 4.1.3 (R Development Core Team 2022).

3 Results and discussion

3.1 Intercropping enhanced productivity, shoot P content, and P-use efficiency under an appropriate P-application rate

The average aboveground biomass and grain yield of intercropping were significantly greater than those of the weighted means of corresponding monocultures by 31% and 36% in 2020, and by 22% and 24% in 2021 (Fig. 2a, b; Table S1; Fig. S1). The aboveground biomass of chickpea/maize, faba bean/maize, oilseed rape/maize, and soybean/maize intercropping across two years was increased by 35%, 30%, 21%, and 15%, and grain yields of those intercropping were enhanced by 39%, 32%, 28%, and 14%, compared with the weighted means of corresponding monocultures, respectively (Fig. 2c, d; Fig. S1). Phosphorus application significantly increased aboveground biomass and grain yields, regardless of monoculture or intercropping system (Fig. 2e, f; Fig. S1). The LERs of chickpea/maize, faba bean/maize, and oilseed rape/maize intercropping systems were greater than 1 under the three P-application rates in 2020 and 2021 (Fig. S2). Several experiments and a meta-analysis also showed the positive effects of intercropping on productivity (Daryanto et al. 2020; Li et al. 2021b; Tang et al. 2021; Ndayisaba et al. 2021). Our recent study also shows that the enhanced productivity might be related to nitrogen (N) acquisition, as intercropping significantly increases the amount of N derived from atmospheric N2 fixation of legumes which indicates that the N-fertilization rate was appropriate in this experiment (Xing et al. 2023). Here, we mainly focus on the acquisition and utilization of P in intercropping.

Fig. 2
figure 2

Effects of year (Y) and cropping system (i.e. monoculture and intercropping, Cs) on aboveground biomass (a) and grain yield (b) in 2020 and 2021. Effects of crop combination (Cc) and cropping system (Cs) on aboveground biomass (c) and grain yield (d) across two years. Effects of phosphorus (P)-application rate and cropping system (Cs) on aboveground biomass (e) and grain yield (f) across two years. Aboveground biomass and grain yield of the monocultures were calculated as the weighted means of corresponding monoculture crops based on their land proportions in intercropping. P0: 0 kg P ha−1; P40: 40 kg P ha−1; P80: 80 kg P ha−1. FM: faba bean/maize intercropping system; SM: soybean/maize intercropping system; CM: chickpea/maize intercropping system; RM: oilseed rape/maize intercropping system. Uppercase letters refer to differences among P-application rates. Lowercase letters indicate differences among treatments if the interaction effect was significant. The same letter means there was no significant difference (Tukey HSD). *, P < 0.05; **, P < 0.01; ***, P < 0.001; ns, not significant.

Full size image

Compared with the weighted means of corresponding monocultures, intercropping significantly increased shoot P content by 12.1 kg ha−1 and 8.8 kg ha−1, increased agronomic efficiency of applied P (reflecting P-fertilizer-use efficiency on productivity) by 58 kg kg−1 and 42 kg kg−1, and increased the apparent recovery efficiency of applied P (an indicator of P-fertilizer-use efficiency on shoot P content) by 28% and 21% in 2020 and 2021, respectively (Fig. 3a, c, d). This is consistent with previous studies and indicates substantial increases in P-use efficiency with intercropping (Darch et al. 2018; Li et al. 2018; Eichler-Loebermann et al. 2021). Oilseed rape/maize combination exhibited greater shoot P content and the apparent recovery efficiency of applied P than other combinations; chickpea/maize exhibited greater agronomic efficiency of applied P than the other three crop combinations (Fig. 3e, g, h).

Fig. 3
figure 3

Effects of year (Y) and cropping system (i.e. monoculture and intercropping, Cs) on shoot phosphorus (P) content (a), internal utilization efficiency of P (b), agronomic efficiency of applied P (c), and apparent recovery efficiency of applied P (d) in 2020 and 2021. Effects of crop combination (Cc) and cropping system (Cs) on shoot phosphorus (P) content (e), internal utilization efficiency of P (f), agronomic efficiency of applied P (g), and apparent recovery efficiency of applied P (h) across two years. Effects of P-application rate and cropping system (Cs) on shoot P content (i), internal utilization efficiency of P (g), agronomic efficiency of applied P (k), and apparent recovery efficiency of applied P (l) across two years. Shoot P content, internal utilization efficiency of P, agronomic efficiency of applied P, and apparent recovery efficiency of applied P of the monocultures were calculated as the weighted means of corresponding monoculture crops based on their land proportions in intercropping. P0: 0 kg P ha−1; P40: 40 kg P ha−1; P80: 80 kg P ha−1. FM: faba bean/maize intercropping system; SM: soybean/maize intercropping system; CM: chickpea/maize intercropping system; RM: oilseed rape/maize intercropping system. Uppercase letters refer to differences among crop combinations. Lowercase letters indicate differences among treatments if the interaction effect was significant. The same letter means there was no significant difference (Tukey HSD). *, P < 0.05; **, P < 0.01; ***, P < 0.001; ns, not significant.

Full size image

Shoot P content was positively affected by P-application rates (P< 0.001), and the increase of P content in intercropping with P40 and P80 was greater than in P0 (Fig. 3i; Table S1). Shoot P content of intercropping systems under 40 kg P ha−1 supply was even greater than that of monoculture systems under 80 kg P ha−1 supply, which indicates that less P-fertilizer input can obtain the same or a greater P content through intercropping.

The apparent recovery efficiency of applied P and agronomic efficiency of applied P are widely used approaches to estimate P-use efficiency (Dobermann 2007; Adiele et al. 2020); the apparent recovery efficiency of applied P and agronomic efficiency of applied P under P40 were greater than that in P80, regardless of intercropping systems or the weighted means of corresponding monocultures, and the increase of those in intercropping with P40 was greater than that in P80 (Fig. 3k, l), indicating that intercropping produced more P content or yield under per unit of P application under 40 kg P ha−1 than 80 kg P ha−1, i.e. exhibited the potential to save P fertilizer and contribute to the sustainable management.

Intercropping decreased the internal utilization efficiency of P (i.e. grain yield: grain yield P content ratio) in 2020 and 2021 (Figs. 3b, S3, S4). The internal utilization efficiency of P of intercropping was lower than that of the weighted means of corresponding monocultures in chickpea/maize and oilseed rape/maize, while it did not change that in the other two combinations across the three P-application rates (Fig. 3f). In addition, intercropping decreased internal utilization efficiency of P in P40 and P80, but not in the P0 treatment compared with the weighted means of corresponding monocultures across the four crop combinations (Fig. 3j). The decreased internal utilization efficiency of P in intercropping systems was mainly affected by maize; internal utilization efficiency of P of intercropped maize decreased by 8.6% compared with monoculture maize (Fig. S6) which probably reduced internal utilization efficiency of P in intercropping systems (Fig. 3b, f, j). In addition, companion crops could meet P requirement of intercropped maize via increased availability of soil P in the rhizosphere which may also decrease internal utilization efficiency of P (Li et al. 2014). The decreased internal utilization efficiency of P in intercropping under P application indicate that intercropping systems could meet crop P requirement under P application or exhibit a lower degree of P deficiency stress, while other abiotic stresses may become a limiting factor for grain yield (Dobermann 2007; Adiele et al. 2020).

Our results indicated that P40 was sufficient to maintain relatively high productivity and shoot P content, with the highest apparent recovery efficiency of applied P and agronomic efficiency of applied P of all intercropping systems. The average P-application rate on the North China Plain is 92 kg P ha−1 yr−1 (Vitousek et al. 2009); our findings indicate that intercropping with 40 kg P ha−1 yr−1 is a desirable rate to achieve high P-use efficiency, which may have the potential to save P fertilizer compared with conventional P-application rates in monocultures.

Linear-mixed effect models showed that the productivity, shoot P content, agronomic efficiency of applied P, and the apparent recovery efficiency of applied P varied among crop species depending on cropping system and P-application rate (Figs. S5, S6; Table S2). Intercropped maize exhibited greater aboveground biomass and grain yield than monoculture (Fig. S5). The results also showed that intercropped maize had greater shoot P content, the apparent recovery efficiency of applied P, and agronomic efficiency of applied P than monoculture maize (Fig. S6). The shoot P content of intercropped maize, intercropped oilseed rape, and intercropped chickpea were 16.7 kg ha−1, 5.3 kg ha−1, and 3.3 kg ha−1 greater than those of corresponding monocultures across two years (Fig. S6e). Therefore, we highlight that the advantage of P content and P-use efficiency in intercropping systems was mainly due to the greater shoot P content and P-use efficiency of maize.

3.2 The positive complementarity effect partly drives enhanced P content in intercropping

For P content, the complementarity effect in faba bean/maize, chickpea/maize, and oilseed rape/maize systems were significantly greater than zero under three P-application rates in 2021 and under P40 and P80 in 2020; the soybean/maize system exhibited a positive complementarity effect in P40 and P80 in two years which may be due to facilitation and niche differentiation. At all P-application rates, the selection effect calculated on the basis of P content was significantly greater than zero in the soybean/maize system (Table 1). The complementarity effect of P content in P40 and P80 was greater than that without P application in all four crop combinations.

Table 1 Complementarity effect and selection effect of shoot P content under different phosphorus (P) application rates in different crop combinations in 2020 and 2021.
Full size table

For P content, the mean value of selection effect, which was significantly greater than zero, was greater than complementarity effect in soybean/maize and oilseed rape/maize intercropping systems without P application (Table 1), reflecting the potential interspecific competition by crop species with greater biomass, i.e. maize (Li et al. 2018). Oilseed rape is highly dependent on soil P availability and grows poorly without P application and is affected by competition from maize, resulting in a positive selection effect calculated on the basis of P content without P application. Additionally, maize exhibited stronger interspecific competition in intercropping systems because of its greater root extension and light interception, resulting in shading of legumes and limited legume root distribution (Gao et al. 2010a; Gao et al. 2010b; Li et al. 2021a; Zhang et al. 2021). Another reason for the greater selection effect in soybean/maize intercropping was that the temporal niche differentiation (TND) in soybean/maize was lower than that in other crop combinations in this study, leading to greater competition for maize when intercropped with soybean. The longer co-growth period (lower TND) limits crops from acquiring light and nutrients, leading to competition (Huang et al. 2017; Li et al. 2020). For P content, the higher TND may partly explain the positive complementarity effect in other crop combinations.

For P content, the mean value of complementarity effect was higher than that of selection effect in faba bean/maize and chickpea/maize intercropping systems under all P levels (except for faba bean/maize under P0 in 2020) and oilseed rape/maize intercropping systems at P40 and P80 in both years (Table 1). For P content, the complementarity effect of soybean bean/maize was higher than selection effect at P40 and P80 in 2020 and P40 in 2021. Thus, the complementarity effect was more important in these combinations, suggesting that niche differentiation and facilitation played an important role in P uptake (Barry et al. 2019).

The sole crops showed significant differences in physiological traits. In the principal component analysis, the two principal axes explained 60.9% of the variance in physiological traits across crop species (Fig. 4). Soybean, oilseed rape, and faba bean exhibited higher leaf [Mn] than maize (Fig. 4; Table S5), which partly reflects the higher release of P-mobilizing carboxylates and protons in the rhizosheath (Yu et al. 2020a; b; Wen et al. 2021; Yu et al. 2023). Oilseed rape and faba bean showed greater acid phosphatase activity and alkaline phosphatase activity in rhizosheath (Tables S3, S4), which indicates a strong ability to mobilize soil organic P (Nannipieri et al. 2011). Chickpea and soybean showed greater soil microbial biomass P in the bulk soil than maize did, indicating a potentially greater competition of soil microorganisms to take up soil P (Fig. 4; Table S6). Soil microbial biomass P in the bulk soil is potentially available to plant species, which can maintain P in labile forms away from reactions from soil (Olander and Vitousek 2004; Richardson and Simpson 2011; Peng et al. 2021). Overall, the companion crops (oilseed rape, faba bean, chickpea, and soybean) exhibited stronger P-mobilizing capacity than maize did (Fig. 4; Tables S3-S7); they were efficient in mobilizing sorbed soil P and increased P availability. The advantage of P uptake may occur through niche differentiation and interspecific facilitation in intercropping systems resulting in greater complementarity effect in four intercropping systems.

Fig. 4
figure 4

Principal component analysis (PCA) of phosphorus (P)-related functional traits of five tested crop species in monoculture in 2020 and 2021. PC1 represents the first axis, PC2 represents the second axis, and the percentage number represents the proportion of variation the axis could explain. ACr, rhizosheath acid phosphatase activity; AKr, rhizosheath alkaline phosphatase activity; pHr, rhizosheath pH; MBP, soil microbial biomass P in bulk soil; leaf [Mn]: leaf manganese concentration. P value means the differences between crop species by PERMANOVA test.

Full size image

3.3 The higher P content of maize when intercropped was associated to its altered root physiological traits by the neighbor species

When growing maize and efficient P-mobilizing species together, compared with monoculture, rhizosheath alkaline phosphatase activity and leaf [Mn] of maize increased across two years (Fig. 5; Table S8). The physiological traits of maize were also affected by P-application rates. Rhizosheath acid phosphatase activity and rhizosheath alkaline phosphatase activity without P application were significantly greater than those at P80 (Fig. 5); soil microbial biomass P in bulk soil and leaf [Mn] without P application was significantly lower than that at P80 (Fig. 5d, e). The physiological traits of maize were weaker than those of companion legumes and oilseed rape; therefore, we surmise the enhanced performance of maize was related to the direct facilitation by neighboring species. Phosphorus facilitation by efficient P-mobilizing species has been shown in intercropping systems, e.g., in faba bean/maize and white lupin/wheat intercropping systems, efficient P-mobilizing legumes mobilize inorganic P (Cu et al. 2005; Li et al. 2007). In a chickpea/wheat intercropping system, chickpea mobilizes organic P and enhances the system’s P content (Li et al. 2003). A recent meta-analysis also shows that compared with monocultures, plant species mixtures enhance soil P availability and facilitate P uptake in natural ecosystems and agroecosystems (Chen et al. 2022).

Fig. 5
figure 5

Effects of phosphorus (P)-application rate and cropping system (i.e. monoculture and intercropping, Cs) on acid phosphatase activity in rhizosheath soil (a), alkaline phosphatase activity in rhizosheath soil (b), pH in rhizosheath soil (c), soil microbial biomass P in bulk soil (d) and leaf manganese concentration (leaf [Mn], e) of maize across two years. P0: 0 kg P ha−1; P40: 40 kg P ha−1; P80: 80 kg P ha−1. Lowercase letters indicate differences among treatments if the interaction effect was significant. Uppercase letters refer to differences among P-application rates. The same letter means there was no significant difference (Tukey HSD). #, 0.05 < P < 0.1; *, P < 0.05; **, P < 0.01; ***, P < 0.001; ns, not significant.

Full size image

The shoot P content of intercropping systems and the weighted means of corresponding monocultures were positively correlated with rhizosheath alkaline phosphatase activity, leaf [Mn], and soil microbial biomass P in bulk soil of maize at flowering stages, indicating the greater P-mobilizing ability of intercropped maize was stimulated by companion crops and increased P uptake (Fig. 7a). The apparent recovery efficiency of applied P and agronomic efficiency of applied P of intercropping systems and the weighted means of corresponding monocultures showed a significantly positive correlation with maize physiological traits, especially with acid phosphatase activity and pH in rhizosheath, indicating that the increase of physiological traits of maize may also enhance the apparent recovery efficiency of applied P and agronomic efficiency of applied P (Fig. 7a). By contrast, the physiological traits of companion crops did not correlate with P-related agronomic indicators (Fig. S11). Acid phosphatase may provide a pathway to explain the correlation between rhizosphere heterogeneity and P acquisition in diverse crops (Bargaz et al. 2017). In addition to physiological traits (i.e. exudates), root morphological traits should also be integrated into P-acquisition efficiency (Campos et al. 2018).

3.4 The higher P content was also associated with thinner roots of intercropped maize compared with monocultures

Compared with monoculture, root surface area, specific root length, root branching intensity of maize increased, while root tissue density decreased in intercropping (Fig. 6). This indicates that maize, a crop species that is inefficient at P mobilization, exhibited thinner roots and more acquisitive root morphological traits in intercropping systems (Cahill et al. 2010; Wen et al. 2020). The shifts of physiological and morphological traits of maize and companion species were also affected by P-application rates (Figs. 6, S7-S10).

Fig. 6
figure 6

Effects of phosphorus (P)-application rate and cropping system (i.e. monoculture and intercropping, Cs) on root surface area (a), root tissue density (b), average root diameter (c), specific root length (d), root branching intensity (e), and proportion of fine roots (f) of maize in July 2021. P0: 0 kg P ha−1; P40: 40 kg P ha−1; P80: 80 kg P ha−1. Lowercase letters indicate differences among treatments if the interaction effect was significant. Uppercase letters refer to differences among P-application rates. The same letter means there was no significant difference (Tukey HSD). #, 0.05 < P < 0.1; *, P < 0.05; **, P < 0.01; ***, P < 0.001; ns, not significant.

Full size image

Specific root length and root tissue density without P application were significantly higher than those at P40 and P80 (Fig. 6b, d); root branching intensity and average root diameter were greater in P40 and P80 than in P0 (Fig. 6c, e). For companion crops, intercropping decreased both root surface area and specific root length in faba bean and chickpea (Figs. S7-S8). Similarly, intercropping also decreased the proportion of fine roots in oilseed rape (Fig. S9). Intercropped soybean had a greater proportion of fine roots compared with the monoculture (Fig. S10). The variation in root morphological traits in intercropped companion crops, except for soybean were not conducive to soil nutrient acquisition. Intercropped maize showed greater root morphological plasticity than intercropped companion crops, indicating that maize may be more effective at resource-acquisition in intercropping. Similar responses were found in peanut/maize and soybean/maize intercropping systems under different N-applications, where intercropped maize showed greater root length density, root weight density, and total root surface than the monoculture (Yang et al. 2022). A global root trait database also showed that domestication of crops, which tend to have higher specific root length and lower root diameter than their wild relatives confers an advantage for successful resource acquisition in organic agroecosystems (Isaac et al. 2021). Our evidence shows that maize exhibited greater morphological plasticity in intercropping which may enhance maize benefits via interspecific facilitation and lead to intercropping advantages (Yu et al. 2020a).

Plant species with a greater specific root length and root branching intensity, but lower root tissue density, exhibit a fast resource-acquisitive strategy and facilitate P acquisition (Wen et al. 2022). We found that the apparent recovery efficiency of applied P and agronomic efficiency of applied P of intercropping systems and the weighted means of corresponding monocultures were tightly correlated with morphological root traits (e.g., specific root length, root surface area, the proportion of fine roots, and average root diameter) rather than physiological traits of maize (Fig. 7). By contrast, the morphological traits of companion crops had insignificant effects on P-related agronomic indicators (Fig. S11), indicating that maize had a greater impact on efficient P-uptake than companion crops by changes of root morphological traits. Root morphological traits (64%) showed a greater relative contribution to the apparent recovery efficiency of applied P of intercropping systems and the weighted means of corresponding monocultures than root physiological traits (27%) and microorganisms (9%, Fig. 7b). Similarly, the relative contribution of maize morphological traits to the agronomic efficiency of applied P, internal utilization efficiency of P, shoot P content of intercropping systems and the weighted means of corresponding monocultures was 60%, 46%, and 40%, respectively (Fig. 7c-e). We also found that root morphological traits of maize were positively correlated with greater productivity, the apparent recovery efficiency of applied P, agronomic efficiency of applied P, internal utilization efficiency of P, and shoot P content of intercropping systems and the weighted means of corresponding monocultures under each P-application rate (Fig. S12).

Fig. 7
figure 7

Correlations between P-related parameters of maize and apparent recovery efficiency of applied phosphorus (P) (system-AREP), agronomic efficiency of applied P (system-AEP), internal utilization efficiency of P (system-IEP), and shoot P content (system-shoot P) of intercropping systems and the weighted means of corresponding monoculture systems (a) in 2021. The 11 P-related parameters of maize, system-IEP, and system-shoot P were turned into relative P0 data to maintain consistency with the sample size of system-AREP and system-IEP. The red and blue boxes indicate positive and negative correlations, respectively. Relative contribution of maize morphological traits (including SRL, SRA, RBI, PFR, ARD, RTD), physiological traits (including ACr, pHr, Mn), and microorganisms (including AKr, MBP) to system-AREP (b), system-AEP (c), system-IEP (d), and system-shoot P (e) in 2021. Mn: leaf manganese concentration; pHr: rhizosheath pH; ACr: rhizosheath acid phosphatase activity; AKr: rhizosheath alkaline phosphatase activity; MBP: soil microbial biomass P in bulk soil; SRL: specific root length; SRA: root surface area; RBI: root branching intensity; PFR: proportion of fine roots; ARD: average root diameter; RTD: root tissue density. *, P < 0.05; **, P < 0.01; ***, P < 0.001.

Full size image

The shifts in morphological root traits of maize may be more effective to take up inorganic P from soil (Lyu et al. 2016). According to the root economic spectrum, changes in root morphological traits may incur more return than shifts of root physiological traits of plant species that are inefficient at P-mobilization (Wen et al. 2022). Therefore, we surmise that thinner roots of maize in intercropping are more effective at exploring P resources and enhancing the apparent recovery efficiency of applied P and agronomic efficiency of applied P than the changes in root physiological traits, which are partly mediated by interspecific facilitation.

3.5 Implications for sustainable P management in crop diversification

The enhanced apparent recovery efficiency of applied P and shoot P content of the present study (2020-2021) in the four intercropping systems compared with monocultures were greater than those in the early stage of the experiment (2010-2011, Xia et al. 2013), which suggests that long-term intercropping increases the magnitude and stability of P uptake advantage. Several long-term experiments showed that grain yields increase significantly over time with greater year-to-year stability in intercropping systems than in monocultures (Li et al. 2021b). A meta-analysis also showed that diversity effects on productivity increased over time in other terrestrial ecosystems (Cardinale et al. 2007). Therefore, long-term field crop diversification experiments are needed to explore the fertilizer-saving potential and to establish sustainable agroecosystems.

Management practices such as crop diversification (e.g., intercropping or rotational cropping) and selection of crop genotypes can reduce negative plant-soil feedback (e.g., accumulation of plant pathogens) or enhance the positive effect of legacies (e.g., greater soil biodiversity, better soil structure, and lower pest density) to benefit the performance of the crops (Jing et al. 2022). Most studies focused on short-term legacies, while long-term effects on soil legacies have received little attention. A greenhouse study was conducted to test microbial effects on growth of crop species by collecting soil from a 10-year field experiment. The results showed that positive microbial legacies drove overyielding by repressing pathogens and increasing beneficial microorganisms with long-term crop diversification (Wang et al. 2021a). The legacy P assessment model also exhibits that long-term legacy soil P decreased P demand without crop production loss for 69 years (Yu et al. 2021b). Therefore, long-term positive soil legacies may improve P uptake and reduce the demand for P fertilizer which is urgently needed to understand how long-term intercropping contributes to sustainable P management (Fig. 1c).

Plant species show different strategies in P acquisition; combining P-efficient crop species or genotypes (e.g., P-acquisition efficiency and P-utilization efficiency) is a promising way to tighten the P cycle in agroecosystems (Cong et al. 2020). The new diversified cropping systems can be designed by combing species and genotype diversity to enhance P-use efficiency. In addition, the phenotypic plasticity of root traits varies among crop species; selecting crop species or genotypes that exhibit a greater trait plasticity to respond to soil nutrient availability that is enhanced by neighboring species may be better facilitated and enhance P-use efficiency in intercropping (Yu et al. 2021a).

Overall, the morphological and physiological root traits of maize are facilitated by companion crops, enhancing productivity, P uptake, and P-use efficiency of intercropping systems. The thinner roots of maize in intercropping mainly drive the increase of the P-use efficiency in intercropping systems and save P-fertilizer inputs. Crops tend to recover only 15%-20% of applied P fertilizer in China during the current growing season, and most of the P fertilizer is sorbed in the soil (Zhang et al. 2008). Compared with intensive monocultures, intercropping, with desirable crop combinations, may unleash species’ potential to acquire P efficiently via interspecific facilitation, consequently increasing productivity, P uptake, agronomic efficiency of applied P, and the apparent recovery efficiency of applied P. This may enhance P-use efficiency and improve the sustainability of agriculture. The greater performance of maize in intercropping suggests an adaptive trait suited to neighboring crop species as well as a strategy that can be applied to impact root growth to improve crop nutrient-use efficiency in agroecosystems.

4 Conclusions

We suggest intercropping with an appropriate P-application rate maintains relatively high productivity (i.e. aboveground biomass and grain yield), shoot P content, agronomic efficiency of applied P, and the apparent recovery efficiency of applied P. Our findings concern not only the role of physiological traits but also the significant role of morphological traits of maize on enhanced P-use efficiency under intercropping systems. Enhanced shoot P content, agronomic efficiency of applied P, and the apparent recovery efficiency of applied P by 9.9 kg ha−1, 47.6 kg kg−1, and 23.5% were observed in four intercropping systems than the weighted means of corresponding monocultures, partly accounted for by facilitation underlying a positive complementarity effect. Root morphology and physiology of maize rather than companion crops were related to the increase of the apparent recovery efficiency of applied P and agronomic efficiency of applied P in intercropping systems, in which morphological traits of maize mainly contributed to a greater the apparent recovery efficiency of applied P and agronomic efficiency of applied P. Our findings highlight that root traits are important in understanding P facilitation and designing diversified crop systems for sustainable P management in agroecosystems. Desirable crop combinations increase the P-use efficiency through belowground facilitation, which may decrease environmental impact and improve the sustainability of agriculture in terms of P management globally.