Insufficient nitrogen supply from symbiotic fixation reduces seasonal crop growth and nitrogen mobilization to seed in highly productive soybean crops.

Nitrogen (N) supply can limit the yields of soybean [Glycine max (L.) Merr.] in highly productive environments. To explore the physiological mechanisms underlying this limitation, seasonal changes in N dynamics, aboveground dry matter (ADM) accumulation, leaf area index (LAI) and fraction of absorbed radiation (fAPAR) were compared in crops relying only on biological N2 fixation and available soil N (zero-N treatment) versus crops receiving N fertilizer (full-N treatment). Experiments were conducted in seven high-yield environments without water limitation, where crops received optimal management. In the zero-N treatment, biological N2 fixation was not sufficient to meet the N demand of the growing crop from early in the season up to beginning of seed filling. As a result, crop LAI, growth, N accumulation, radiation-use efficiency and fAPAR were consistently higher in the full-N than in the zero-N treatment, leading to improved seed set and yield. Similarly, plants in the full-N treatment had heavier seeds with higher N concentration because of greater N mobilization from vegetative organs to seeds. Future yield gains in high-yield soybean production systems will require an increase in biological N2 fixation, greater supply of N from soil or fertilizer, or alleviation of the trade-off between these two sources of N in order to meet the plant demand.

INTRODUCTION

n class="Species">Soybean [pan> class="Species">Glycine max (L.) Merr.] is the world's most important legume crop and the primary vegetable protein source for food and feed (Wilson, 2008). Its nitrogen (N) requirement per unit of photosynthate produced is the highest among all food crops (Sinclair & de Wit, 1975). On average, soybean requires 80 kg N in aboveground dry matter (ADM) per metric ton of seed produced (Salvagiotti et al., 2008; Tamagno et al., 2017). This value is about three times greater than the N requirement per unit of grain in cereal crops such as rice, wheat, or maize (Barraclough et al., 2010; Setiyono, Walters, Cassman, Witt, & Dobermann, 2010a; Yin et al., 2019). Unlike cereal crops, fertilizer N is rarely used on soybean, except for an occasional small amount applied as ‘starter’ at sowing. Hence, soybean relies on two major N sources to meet its large N requirement: (a) biological N2 fixation (BNF) and (b) N from soil organic matter mineralization, inorganic soil N left by previous crop, atmospheric dry and wet deposition, water table and irrigation water (hereafter collectively called ‘indigenous soil N supply’). However, it is unclear if N from BNF and indigenous soil N supply are sufficient to meet soybean N requirement as average yield continues to increase due to improved cultivars and agronomic practices (Grassini, Specht, Tollenaar, Ciampitti, & Cassman, et al., 2014; Specht et al., 2014).

There is a well‐known trade‐off between BNF and indigenous sn class="Chemical">oil N supply, that is, BNF decreases as the contribution from inpan>digenous span> class="Chemical">oil N supply increases and vice versa (Santachiara, Borrás, Salvagiotti, Gerde, & Rotundo, 2017; Streeter & Wong, 1988). As a result, application of fertilizer N reduces BNF, resulting in a small marginal increase in crop N uptake and yield (Salvagiotti et al., 2008). While this trade‐off has often been reported, it is not clear whether the combined N supply from BNF and soil is sufficient to meet the N requirement in soybean in highly productive environments (Ciampitti & Salvagiotti, 2018). To address this question, Cafaro La Menza, Monzon, Specht, and Grassini (2017); Cafaro La Menza et al. (2019) developed a protocol to assess N limitation across a wide range of environments, where seed yields ranged from 2.5 to 6.7 Mg ha−1. The protocol consisted of a side‐by‐side comparison between a ‘control’ treatment (hereafter called ‘zero‐N’) that forced the crop to rely on BNF and indigenous soil N supply, and a ‘full‐N’ treatment designed to provide the crop with fertilizer N supply to optimally match the expected seasonal plant N demand. The full‐N treatment aimed to eliminate N limitation at any time of the growing season, independent of the contribution from BNF. Average seed yield was higher in the full‐N than in the zero‐N treatment, but the magnitude of the yield difference was larger in highly productive environments, confirming the existence of a N limitation in high‐yield soybean.

The studies by Cafaro La Menza et al. (2017, 2019) also documented that higher seed yieln class="Chemical">ds inpan> the full‐N versus zero‐N treatment were associated with greater end‐of‐season shoot pan> class="Disease">dry matter, total N, seed number and seed mass. Despite the other well‐known trade‐off between seed yield and seed N concentration (Chung et al., 2003), it was remarkable that the full‐N treatment also exhibited higher seed N concentration, which is of interest for soybean processors who seek high protein concentration when producing meal for animal feed (Brumm & Hurburgh Jr, 2006). In contrast, there were no differences in harvest index, N‐use efficiency, N‐harvest index or seed oil concentration between the full‐N and zero‐N treatments. While these prior studies advanced knowledge of N limitation in high‐yield soybean crops, they did not provide insight into the physiological mechanisms over a seasonal timeframe that eventually led to the differences in seed yield and seed N concentration. Such an assessment would be useful for a more mechanistic understanding of the N limitation in high‐yield soybean production systems.

Considering that n class="Species">soybean yield has to inpan>crease to meet future demanpan>d while avoidinpan>g further conversion of natural ecosystems inpan>to croplanpan>d (Cassmanpan> & Grassinpan>i, 2020), anpan>d that future yield gainpan>s should not have negative effects on seed proteinpan> concentration, it is importanpan>t to designpan> strategies to overcome the N limitation inpan> high‐yield pan> class="Species">soybean systems. With that goal in mind, the objective of this study was to better understand seasonal N dynamics by which key physiological mechanisms (e.g., BNF, ADM and N accumulation, leaf area index (LAI), photosynthesis, and N mobilization) account, separately or in concert, for the observed differences in seed yield and seed N concentration between soybean crops growing under contrasting N supply. For that comprehensive assessment of N limitation in soybean, we used original data on accumulated ADM and N, BNF and other physiological processes collected during the entire crop season from field experiments conducted over 2 years in a high‐yield environment.

MATERIALS AND METHODS

Experimental sites and design

This article presents new data of seasonal dynamics of n class="Chemical">ADM, N, LAI, anpan>d other physiological variables from the set of experiments conducted inpan> Nebraska (NE), USA by Cafaro La Menza et al. (2017, 2019) (Table 1). These previous studies addressed the occurrence of N limitation inpan> pan> class="Species">soybean at different seed yield levels and the influence of indigenous soil N supply on the magnitude of N limitation by looking at end‐of season productivity and total N. Here, we investigated the underlying mechanisms that explain the nature of the N limitation in high‐yield production systems by looking at original data collected during the growing season from the same set of experiments. These NE irrigated production systems provide a suitable background to evaluate N limitation in high‐yield soybean, considering that NE has the largest irrigated soybean area in the world (ca. 2.6 million ha) where producers achieve average irrigated yields of ca. 4.4 Mg ha−1 (USDA‐NASS, 2014–2018; https://www.nass.usda.gov/).

TABLE 1

Description of the seven field experiments conducted in Nebraska, USA.

Site	Crop season	Cultivar (and MG)	Sowing date	Seeding rate (m−2)	Yield potential (Mg ha−1)	Tillage method	pH a	SOC a (g kg−1)	Soil P a (mg kg−1)	Soil K a (mg kg−1)	Applied fertilizer (kg ha−1) b
Mead (41.244 N; 96.502 W)	2016	AG2723 (MG 2.7)	May 8	25	6.4	Strip till	6.3	22	12	335	N (20); P (21); S (17); Zn (2)
	2017	AG2723 (MG 2.7)	May 6	35	6.4	Strip till	6.7	24	25	268	N (12); P (25)
Saronville (40.601 N; 97.966 W)	2016	AG2431 (MG 2.4)	April 26	25	5.5	Disk	6.2	21	17	395	N (17); P (35); S (10); Zn (1)
	2017	P33T72R (MG 3.3)	April 23	35	5.6	Strip till	7.1	20	33	438	N (6); P (20)
Smithfield (40.538 N; 99.683 W)	2016	P24T19 (MG 2.4)	May 13	32	5.7	No till	7.5	16	12	487	P (38); S (1)
	2017	P27T59 (MG 2.7)	May 25	44	4.8	No till	7.2	18	25	373	N (24); P (48); S (27)
Atkinson (42.658 N; 99.029 W)	2017	AG2723 (MG 2.7)	April 24	40	6.3	Disk	7.0	14	7	184	N (20); P (22); K (50); Ca (62); Mg (8); S (33); Zn (2)

Note: K: 1 M NH4‐acetate extractable potassium; P: Bray‐1 phosphorus; MG: maturity group; SOC: soil organic carbon; S: sulfur; Ca: calcium; Mg: magnesium; Zn: zinc.

Average concentration in the 0.60 m depth measured at crop emergence (VE stage).

Parenthetic values indicate rates of elemental nutrient applied to both treatments (full‐N and zero‐N) at or before sowing in each environment.

n class="Chemical">Note: K: 1 M n class="Chemical">NH4‐acetate extractable potassium; P: Bray‐1 phosphorus; MG: maturity group; SOC: soil organic carbon; S: sulfur; Ca: calcium; Mg: magnesium; Zn: zinc.

Experiments were conducted in pivot‐irrigated producer fieln class="Chemical">ds located at four sites across NE (Mead, Saronville, Smithfield anpan>d Atkinpan>son) durinpan>g 2016 anpan>d 2017. All fielpan> class="Chemical">ds at these sites have consistently produced high soybean yields (>5 Mg ha−1) in previous years. This combination of site‐years portrayed well the range of weather, soil, and management (e.g., tillage methods, cultivar maturity groups [MGs]) used in NE and the rest of the US North‐Central region (Table 1). We excluded the experiment in Atkinson in 2016 due to a severe infestation of powdery mildew (Microsphaera diffusa). A portable weather station was installed within 50 m of each experiment in both years to monitor hourly solar radiation, soil and air temperature, rainfall, relative humidity and wind speed. A well‐validated soybean crop model (SoySim; Setiyono et al., 2010b) was used to simulate the yield potential for each site‐year using local weather and recorded sowing date and cultivar MG at each site‐year (Table 1). In all cases, simulated yields exceeded the yield threshold of 4.5 Mg ha−1 proposed in the literature to define high‐yield soybean production environments where N limitation is likely to occur (Salvagiotti et al., 2008).

Each of the seven site‐year combinations is referred to hereafter as an ‘environment’ that included two N treatments (full‐N and zero‐N) in a complete randomized design with four replicates per treatment (size: 176 m2 each). Experimental plots were purposely placed in areas within each field where the highest (n class="Species">maize anpan>d/or pan> class="Species">soybean) yields had been achieved in previous years. Experiments were conducted in fields following a 2‐year rotation with maize, which is the dominant practice in the US North Central region (Grassini et al., 2014). The practice of coating seeds with rhizobia inoculum (commonly referred to as ‘inoculation’) is not needed for proper nodulation in fields with prior soybean history as documented by a number of field studies (e.g., Carciochi et al., 2019; de Bruin et al., 2010; Leggett et al., 2017 and references cited therein). Hence, in our experiments, seeds were treated with fungicide and insecticide, but received no inoculant. Row spacing was 0.76 m in all environments, with seeding rates adjusted to maximize soybean yields (De Bruin & Pedersen, 2009). Soil water content in the upper meter was monitored using Watermark® sensors (Irmak, Payero, VanDeWalle, Rees, & Zoubek, 2014), and maintained above 50% of plant available water throughout the entire growing season, except for a short period of time at Atkinson in 2017 (Supplementary Figure S1). Several prophylactic foliar applications of herbicide, fungicide and insecticide and pre‐sowing nutrient applications (based on soil test results) kept the crops free from biotic and nutrient stresses (aside from N in the zero‐N treatment). Inorganic soil N in the upper 0.60 m at sowing ranged from 25 to 58 kg N ha−1 across environments, which is within the range expected for soybean grown after a maize crop (Farmaha et al., 2016).

Nitrogen treatments

Two n class="Chemical">N treatments were compared inpan> each environment: (a) a ‘zero‐pan> class="Chemical">N’ treatment in which the crop relied on indigenous soil N supply and BNF and (b) a ‘full‐N’ treatment designed to provide the crop with sufficient fertilizer N to optimally match seasonal crop N demand. The N fertilizer was applied as urea and broadcast between plant rows in the full‐N treatment. A total seasonal amount of 870 kg N ha−1 was applied at all sites based on (a) site‐specific yield potential simulated using the SoySim model (Setiyono, et al., 2010b), (b) N uptake requirement of 80 kg N per Mg−1 seed yield (Salvagiotti et al., 2008; Tamagno et al., 2017), and (c) an extra 40% of fertilizer N to compensate for potential N losses through volatilization and leaching resulting from mismatches between irrigation or rainfall events and fertilizer N application. Total N fertilizer amount was split into five applications to approximate the expected increase in crop N requirement during the crop season (Bender, Haegele, & Below, 2015; Thies, Singleton, & Bohlool, 1995). Of the total N fertilizer amount, 10%, 10%, 20%, 30%,and 30% were applied at the V2, V4, R1, R3 and R5 stages, respectively (Figure 1).

FIGURE 1

Scheme showing time of N applications in the full‐N treatment and respective amounts expressed as percentage of total amount of applied N fertilizer. Phenological time is shown in the bottom x‐axis using a dimensionless scale adapted from Lindquist, Arkebauer, Walters, Cassman, and Dobermann (2005) which allows comparisons to be made among environments with different temperature, sowing date, and cultivar maturity group. The major development stages defined by the Fehr & Caviness (1977) scale are shown on the upper x‐axis; VE: emergence, V3: three fully developed leaves at main stem nodes 1, 2, & 3, R1: beginning of flowering, R3: beginning of pod setting, R5: beginning of seed filling, R6: full seed, and R7: physiological maturity. Photos illustrating plants at each of these crop stages are shown at the top

Field measurements

Detailed measurements of phenology, leaf area, photosynthesis and accumulated n class="Chemical">ADM anpan>d N were made at regular inpan>tervals over the entire crop season. Phenological stages definpan>ed by Fehr & Cavinpan>ess (1977) (Figure 1) were recorded from 10 consecutive planpan>ts withinpan> one row inpan> each replicate for both treatments every seven days, startinpan>g at VE anpan>d concludinpan>g at R7. Followinpan>g Linpan>dquist et al. (2005), phenological stages were made comparable across environments (with different temperature regimes, sowinpan>g date anpan>d cultivar pan> class="Chemical">MG) by defining each development stage (DS) using a dimensionless scale where 0, 1 and 2 correspond to the VE, R3 and R7 stages, respectively (Figure 1). In our case, each DS was calculated based on daily mean air temperature using a beta function as defined by Wang & Engel (1998), with appropriate cardinal temperatures for each phase as reported by Setiyono et al. (2007). We did not account for photoperiod in our DS calculation as the four locations were located within a narrow latitudinal band (from 40.5°N to 42.6°N).

Seasonal dynamics of n class="Chemical">ADM anpan>d accumulated N were assessed inpan> each treatment‐replicate by collectinpan>g all planpan>ts withinpan> a 1‐m section of row (inpan>clusive of the 10 planpan>ts from which phenological development was assessed) surrounded by two rows receivinpan>g the same N treatment inpan> order to avoid edge effects. We did not attempt to collect root biomass, which accounts for only ca. 10% of total planpan>t biomass at R7 (Setiyono et al., 2010b). Samples were collected weekly from VE until R7 anpan>d separated inpan>to green leaves, stems, seed, pod walls anpan>d senesced leaves. Green leaves (definpan>ed as anpan>y leaf with >50% green area) were scanpan>ned to determinpan>e the LAI inpan> each pan> class="Chemical">ADM sample (LAI‐3100 area meter LI‐COR, Lincoln, NE). Abscised leaves were also collected every week within a 1‐m long net placed between rows in each treatment‐replicate. All samples were oven‐dried at 70°C until reaching constant weight. Total ADM was calculated as the dry matter sum of all plant organs (included collected abscised leaves). Each plant organ sample was separately ground in a Wiley mill (1‐mm screen mesh), and N concentration was determined with a dry combustion‐based analyzer (LECO Corporation, St Joseph, MI). Larger plant samples (4.6 m2) were collected shortly after R7 from the two central rows in each plot, surrounded by two rows receiving the same N treatment, to obtain an end‐of‐season post‐R7 estimate of seed yield. Two sub‐samples of 200 seeds from the threshed seed of each plot were weighed to estimate the mean individual seed dry mass and to derive the number of seeds per harvested area. Seed yield and seed mass were adjusted to 130 g H2O kg−1 seed, which is the market‐based standard for soybean seed moisture content.

Seasonal dynamics of BNF were determined using the natural n class="Chemical">15N abundanpan>ce method (Shearer & Kohl, 1986) inpan> the zero‐N treatment. Due to logistic constrainpan>ts, BNF was measured inpan> five of the seven environments. We present here data on BNF inpan> the zero‐N treatment only; determinpan>ation of BNF inpan> the full‐N treatment was not possible due to the high isotopic fractionation as a result of inpan>creased denitrification with fertilizer addition (Mathieu et al., 2007). The natural pan> class="Chemical">15N abundance method requires a ‘reference crop’ that does not fix N2 and that is grown with the same N management of the legume crop (Peoples, Unkovich, & Herridge, 2009; Unkovich et al., 2008). Following previous studies (Collino et al., 2015; Pate, Unkovich, Armstrong, & Sanford, 1994), we used maize grown in ‘N‐omission’ plots as the reference crop to measure BNF. A N‐omission plot (9.1 × 11 m) located adjacent (<30 m) to each soybean experiment was sown with maize on the same date as soybean and did not receive any N fertilizer (but received P and K fertilizer). Calculation of BNF following the natural 15N abundance method requires that timing of N uptake between soybean and the reference crop is similar, regardless of the total amount of N absorbed (Witty, 1983). In our case, the maize plot was sown (and emerged) on the same date as soybean; maize relative maturity was purposely selected so that both maize and soybean crops reached physiological maturity around the same date. Also, rate of root depth over time and distribution of root length with depth are similar between the two crops, as has been documented in field studies conducted in the US Corn Belt (Nichols et al., 2019; Ordóñez et al., 2018). To summarize, we believe that our experimental setup was appropriate to ensure that the timing of N uptake coincided between the reference maize crop and soybean as it was required to calculate BNF following the natural 15N abundance method.

Every week, we collected aboveground plant samples consisting of (a) 0.5‐m row of consecutive n class="Species">soybean planpan>ts inpan> the zero‐N treatment, anpan>d (b) two pan> class="Species">maize plants in the N‐omission maize plot. These samples were collected on the same dates as other plant samples for ADM and accumulated N determination. The natural 15N abundance was measured using an elemental analyzer interfaced to a continuous flow isotope ratio mass spectrometer (UC Davis, California, USA). The fraction of the total N in ADM derived from BNF was calculated as:where δ15Nref and δ15Nsoy are the natural 15N abundance of the reference crop (maize grown in the N‐omission plot) and soybean, respectively, and B is the natural 15N abundance in soybean that relies only on BNF. The B value used in this work was −1.75, which corresponds to the average of the B values reported in the literature (Balboa & Ciampitti, 2020; Unkovich et al., 2008). Sensitivity analysis revealed that ±20% variation in the B value would have resulted in a range of BNF that fall within the experimental error, which justifies the use of an average B value from the literature. In the case of δ15Nref, we adjusted a quadratic model to smooth the observed variation among sampling times, with r of fitted models ranging from 0.67 to 0.93 across experiments (Supplementary Figure S2). The decline in δ15Nref over time was likely associated with the root system exploring soil layers with different 15N enrichment as the crop season advances (Högberg, 1997; Shearer & Kohl, 1986). For a given sampling time, BNF (in kg N ha−1) was calculated based on the fraction of BNF and accumulated N; indigenous soil N supply was calculated as the difference between accumulated N and BNF. We do not expect results from our analysis to be influenced by the method selected for BNF determination as average BNF at R7 stage measured using the natural 15N abundance method was identical to that estimated independently using the ‘difference method’ (i.e., absorbed N in soybean minus absorbed N in the reference maize crop). However, we reported here only the results on BNF derived from the natural 15N abundance technique as it has been reported in the literature to be less prone to errors compared with the difference method (Unkovich et al., 2008).

Incident and absorbed photosynthetically active radiation (PAR and APAR, respectively) were measured in one or two replicates in each treatment and environment. Measurements were taken every second and recorded as a 30‐minute average, starting soon after VE and ending at R7. All sensors were calibrated by the manufacturer and a cross‐calibration among sensors was performed every year before placing them in the field and again after removal at harvest. Sensors were leveled (if needed) and cleaned every 3 to 5 days during the season. The PAR was measured above the canopy using a point quantum sensor facing up (LI‐190SA, LI‐COR, Lincoln, NE). Transmitted PAR was measured with a singleline quantum sensor (LI‐191SA, LI‐COR, Lincoln, NE) placed at the sn class="Chemical">oil surface diagonally across rows. Total (canpan>opy plus span> class="Chemical">oil) reflected PAR was measured using an inverted point quantum sensor (ibid) placed 2 m above crop canopy. Reflected PAR from soil was measured using an inverted line quantum sensor (ibid) placed 5 cm above soil surface and diagonally across rows. Daily canopy reflected PAR was calculated as total reflected PAR minus soil reflected PAR. Daily APAR was calculated as incident PAR minus transmitted PAR and canopy reflected PAR; and expressed as the fraction of daily incident PAR (fAPAR). Finally, radiation‐use efficiency (RUE) was estimated as the slope of the relationship between accumulated ADM sampling points and their corresponding accumulated APAR from VE to R7 stages.

Leaf photosynthesis was measured at four stages during the growing season (R1, R3, R5 and R6) in all treatments and environments in 2017. All photosynthesis measurements were performed between 10 a.m. and 2 p.m., using only plants that were purposely selected based on any given above‐stated crop stage so that all measured plants were at the same average crop stage recorded on that day. Photosynthesis was measured on the central leaflet of the third most recently developed leaf of one plant in each treatment‐replicate. Light response curves were generated by varying the photosynthetic photon flux density (PPFD), from 1800 to 0 μmol photons m−2 s−1 of a red/blue LED light source of an open‐flow gas exchange system (LI‐6400, Lincoln, NE). n class="Chemical">Carbon dioxide (pan> class="Chemical">CO2) concentration inside the chamber of the open‐flow gas exchange system was maintained at 400 ppm, leaf temperature was set at 25°C, and the minimum measurement waiting time was 60 s or until reaching a coefficient of variation (CV) ≤ 3% of the CO2 assimilation rate.

Data analysis

A logistic model (France & Thornley, 1984) was fitted to the weekly data of n class="Chemical">ADM and accumulated N collected from VE to R7:where W0 is the estimated n class="Chemical">ADM (g m−2) or N (kg ha−1) at emergence, Wf is the maximum ADM or accumulated N during the growing season, t is time in days after emergence, and u is a constant of proportionality of plant growth and its decaying with time. The derivative of Equation 2, with respect to time, represents the crop growth rate (CGR; g m−2 d−1) or N accumulation rate (NArate; kg N ha−1 d−1). Estimating CGR and NArate using this approach helps remove the measurement error associated with specific sampling times (e.g., Bange, Hammer, & Rickert, 1997; Hall, Connor, & Sadras, 1995; Lindquist et al., 2005). Similarly, a Gaussian distribution model was fit to the seasonal dynamics of LAI:where LAImax is estimated maximum LAI, t is time at which at LAImax is reached, and SD is a parameter of the exponential equation.

Seasonal dynamics in fAPAR were compared between N treatments. Because we did not have light sensors in every plot, we used the following approach to estimate daily fAPAR for each experimental unit. First, we generated a relationship between measured fAPAR and LAI using all available dates of LAI sampling (Supplementary Figure S3). The extinction coefficient (k) was estimated to be 0.54 across environments and treatments, without differences between N treatments. Second, we used that relationship to derive the daily fAPAR for each treatment in each experiment based on daily LAI obtained from the fitted models (Equation 3).

Sources of n class="Chemical">carbon for seed pan> class="Disease">dry matter accumulation during the reproductive phase include new photo‐assimilates created on a daily basis (Yamagata, Kouchi, & Yoneyama, 1987), and dry matter mobilization from vegetative organs (Egli, Guffy, & Leggett, 1985; Stephenson & Wilson, 1977). Apparent dry matter mobilization to seeds was estimated as the difference between non‐seed ADM (including stems, green, senesced and abscised leaves, and pod walls) at R5 and R7. We used the term ‘apparent’ because the amount of non‐seed ADM mobilized to the seeds was calculated as the difference in non‐seed ADM between R5 and R7, rather than being directly measured, and also to recognize that it does not account for the conversion efficiency (Borrás, Slafer, & Otegui, 2004). In the case of seed N accumulation, major sources include crop soil N uptake and BNF during seed filling as well as mobilization from non‐seed ADM to the seeds. This mobilized N was estimated following the same approach as for dry matter mobilization. We did not attempt to account for N mobilized from belowground biomass, which we estimated to represent a very low fraction (ca. 3%) of the total N mobilized to seed, based on root‐to‐shoot ratio, root N concentration, and N mobilization fraction reported in the literature (Amthor et al., 1994; Connor, Loomis, & Cassman, 2011; Setiyono et al., 2010b).

Following Connor et al. (2011), the following function was used to quantify the response of net photosynthesis to incident light:where A is the photosynthetic rate (μmol n class="Chemical">CO2 m−2 s−1), Amax is the maximum photosynthetic rate at light saturation (μmol pan> class="Chemical">CO2 m−2 s−1), PPFD is photosynthetic photon flux density (μmol photons m−2 s−1), Ic is the light compensation point (μmol photons m−2 s−1) and α is the initial slope of the response curve. Dark respiration (Rd; μmol CO2 m−2 s−1) can be estimated from the fitted parameter values as –α*Ic.

We used a combined analysis of variance (ANOVA) for the effect of N treatments across environments where environments and N treatments were treated as fixed effects (SAS® PROC MIXED v.9.3; Moore & Dixon, 2015). The combined ANOVA assumes homogeneity of the error variances within a given environment and when pooled over all seven environments. Our analysis indicated that the error variances did not deviate much from homogeneity (F max <6), which led us to consider our analysis to be robust (Milliken & n class="Disease">Johnson, 2009). The followinpan>g directly observed or calculated parameters were evaluated: seed yield, seed number, seed mass, pan> class="Chemical">ADM at R7, accumulated N at R7, seed N concentration and mobilized ADM and N. A similar analysis was performed to identify the crop stage interval with the largest difference between the two N treatments relative to CGR and NArate. To do this, the CGR and NArate obtained per replicate were averaged by treatment over four crop phases (VE‐R1, R1‐R3, R3‐R5, R5‐R7). Variation in seed yield is largely driven by differences in seed set, which, in turn, is modulated by the crop growth rate during key reproductive phases; the latter has been referred to as ‘critical period’ (Egli & Bruening, 2006; Fischer, 1975). Following Kantolic, Peralta, & Slafer (2013), accumulated ADM during the critical period for soybean (between R3 and R6 stages) was used to understand seed number differences between N treatments. Least square means were computed for each treatment, and treatment differences evaluated using Dunnett's test.

Differences in the parameters of the light response curve (Amax, Ic, α, Rd) between N treatments were evaluated using mixed models (InfoStat; Di Rienzo et al., 2011). Linear regression analysis was used to investigate relationships and trade‐offs between the measured and calculated variables. Slope and intercept were calculated and their differences between treatments were tested with F tests. To remove the confounding effect of differences in phenology across environments, seasonal patterns in n class="Chemical">ADM, accumulated N, CGR, pan> class="Chemical">NArate and LAI are shown as a function of the DS calculated for each environment (Figure 1), but with data across environments pooled for each N treatment to facilitate the comparison.

RESULTS

Impact of soybean nitrogen limitation on seed yield components

n class="Species">Soybean seed yield ranpan>ged from 5.3 to 5.8 pan> class="Chemical">Mg ha−1 under zero N and from 5.5 to 6.7 Mg ha−1 with non‐limiting N across the seven environments (Table 2). In all cases, measured yields in the full‐N treatment were within 15% of the simulated yield potential (Table 1). Average seed yield was 0.6 Mg ha−1 greater in the full‐N than in the zero‐N treatment (Table 3). Total accumulated N was 51 kg N ha−1 greater in the full‐N than in the zero‐N treatment (13%; 446 versus 395 kg N ha−1, respectively), which translated to 10% greater ADM at R7 in the full‐N compared with the zero‐N treatment (13.2 versus 12.0 Mg ha−1, respectively). The full‐N treatment exhibited 7% greater seed number (3,301 versus 3,096 seeds m−2), and 4% greater seed mass (185 versus 179 mg) than the zero‐N treatment. Ordinarily, in absence of strong environmental variation, seed number and seed mass are negatively correlated, and that was the case across our irrigated soybean experiments in NE (Figure 2). However, the slope of the linear regression between seed mass and seed number was significantly different between the zero‐N and full‐N treatments (−0.02 versus −0.05 mg per additional seed, respectively; P = 0.005), revealing that the seed mass versus seed number trade‐off was substantively alleviated in the case of the full‐N treatment (Figure 2).

TABLE 2

Mean soybean seed yield, number and mass, aboveground dry matter (ADM), accumulated nitrogen (N) in ADM, seed N concentration at physiological maturity (R7 stage), and mobilized ADM and N from non‐seed ADM to seeds in the two N treatments across seven environments

Location	Year	Treatment	Seed yield (Mg ha−1)	Seed number (m−2)	Seed mass (mg)	ADM (Mg ha−1)	Accumulated N (kg N ha−1)	Seed N concentration (g N kg −1)	Mobilized ADM (Mg ha−1)	Mobilized N (kg N ha−1)	BNF (kg N ha−1) a
Mead	2016	Zero	5.3	3,320	161	12.8	419	56.8	1.0	154	264
		Full	6.2	3,606	172	13.7	460	57.1	1.5	174	n.c.
	2017	Zero	5.4	2,992	179	11.3	370	56.7	1.9	147	242
		Full	6.7	3,557	187	13.1	459	60.0	2.5	179	n.c.
Saronville	2016	Zero	5.8	3,405	171	12.1	408	57.7	0.9	134	270
		Full	6.4	3,476	184	13.4	436	58.0	1.0	140	n.c.
	2017	Zero	5.7	3,317	173	13.1	401	58.0	2.1	167	n.c.
		Full	6.4	3,584	179	14.7	485	62.7	1.6	196	n.c.
Smithfield	2016	Zero	5.3	3,031	174	10.4	364	54.7	1.7	130	293
		Full	5.7	3,065	186	11.5	384	57.5	1.1	157	n.c.
	2017	Zero	5.5	2,676	204	11.2	367	62.3	1.5	151	218
		Full	5.5	2,813	197	12.4	432	64.5	1.0	172	n.c.
Atkinson	2017	Zero	5.5	2,933	188	13.1	436	61.2	1.5	165	n.c.
		Full	5.7	3,009	191	13.3	466	61.6	1.8	214	n.c.

Note: Also shown is biological N fixation (BNF) measured for the zero N treatment.

BNF measured for a subset of five environments. BNF was not calculated (n.c.) for the full‐N treatment given the high isotopic fractionation due to the large amount of N applied (see Section 2.3).

TABLE 3

Analysis of variance for effect of nitrogen (N) treatment on soybean seed yield, number and mass, aboveground dry matter (ADM), accumulated nitrogen (N) in ADM, and seed N concentration at physiological maturity (R7 stage), and mobilized ADM and N from non‐seed ADM.

		Seed yield	Seed number	Seed mass	ADM	Accumulated N	Seed N concentration	Mobilized ADM	Mobilized N
Fixed effects	d.f.	F	F	F	F	F	F	F	F
E	6	7.3***	20.4***	59.2***	16.1***	5.9**	4.6**	5.3**	5.0**
E	6	7.3***	20.4***	59.2***	16.1***	5.9**	4.6**	5.3**	5.0**
N	1	42.6***	18.7***	24.8***	16.8***	27.7***	5.6*	0.2	16.8***
N x E	6	3.0*	2.2	4.1**	0.5*	1.2	0.3	2.9*	0.6
Random effects		MS	MS	MS	MS	MS	MS	MS	MS
Rep (E)	21	0.08	31,207	15.7	0.5	1,011	11	0.25	566
N x rep (E)	21	0.12	31,437	24.1	1.3	1,326	10	0.17	578
Estimated means		(Mg ha−1)	(m−2)	(mg)	(Mg ha−1)	(kg N ha−1)	(g N kg ADM−1)	(Mg ha−1)	(kg N ha−1)
Full‐N		6.1	3,301	185	13.2	446	60	1.5	176
Zero‐N		5.5	3,096	179	12.0	395	58	1.5	150
Difference (full‐N minus zero‐N)		0.6***	205***	7***	1.2***	51***	2*	Nil	26***

Note: Each experiment was considered to be a separate environment (E) for the N × E interaction in this analysis. Also shown are the F‐test values and probabilities for the fixed effects, mean squares (MS) for the random effects, estimated means for each treatment, and contrast between N treatments means. Asterisks indicate statistical significance at *P < 0.05, **P < 0.01, and ***P < 0.001.

FIGURE 2

Observed seed number and seed mass at physiological maturity in the full‐N (red squares) and zero‐N (blue circles) treatments across the seven environments. Connecting grey lines relate paired treatments from the same environment. Dotted lines indicate different seed yields across the seed number and mass ranges. Parameters of the fitted linear regressions (solid blue and red lines for zero‐N and full‐N treatment, respectively) and coefficients of determination (r) are also shown. Note that regression lines do not imply causality (with respect to one trait versus the other); instead, they are shown to illustrate the trade‐off between seed number and mass for each N treatment found across the seven irrigated experiments in Nebraska

Mean n class="Species">soybean seed yield, number and mass, aboveground n class="Disease">dry matter (ADM), accumulated nitrogen (N) in ADM, seed N concentration at physiological maturity (R7 stage), and mobilized ADM and N from non‐seed ADM to seeds in the two N treatments across seven environments

n class="Chemical">Note: Also shown is biological n class="Chemical">N fixation (BNF) measured for the zero N treatment.

Bn class="Chemical">NF measured for a subset of five environments. Bn class="Chemical">NF was not calculated (n.c.) for the full‐N treatment given the high isotopic fractionation due to the large amount of N applied (see Section 2.3).

Analysis of variance for effect of n class="Chemical">nitrogen (N) treatment on n class="Species">soybean seed yield, number and mass, aboveground dry matter (ADM), accumulated nitrogen (N) in ADM, and seed N concentration at physiological maturity (R7 stage), and mobilized ADM and N from non‐seed ADM.

n class="Chemical">Note: Each experiment was considered to be a separate environment (E) for the n class="Chemical">N × E interaction in this analysis. Also shown are the F‐test values and probabilities for the fixed effects, mean squares (MS) for the random effects, estimated means for each treatment, and contrast between N treatments means. Asterisks indicate statistical significance at *P < 0.05, **P < 0.01, and ***P < 0.001.

Observed seed number and seed mass at physiological maturity in the full‐N (red squares) and zero‐N (blue circles) treatments across the seven environments. Connecting grey lines relate paired treatments from the same environment. Dotted lines indicate different seed yieln class="Chemical">ds across the seed number and mass ranges. Parameters of the fitted linear regressions (solid blue and red lines for zero‐N and full‐N treatment, respectively) and coefficients of determination (r) are also shown. Note that regression lines do not imply causality (with respect to one trait versus the other); instead, they are shown to illustrate the trade‐off between seed number and mass for each N treatment found across the seven irrigated experiments in Nebraska

Differences in seasonal accumulated ADM and N between treatments

Accumulated n class="Chemical">ADM anpan>d N followed sigmoidal patternpan>s (Figure 3). The sigmoidal inpan>flection poinpan>ts corresponded to ca. R5 anpan>d to a pan> class="Chemical">DS value of ca. 1.25, which are more precisely documented in the inset plots that show the peaks of the derivative variables CGR and NArate. It was also evident in the sigmoid patterns that greater ADM and accumulated N became visibly different just before R1 in the full‐N compared with the zero‐N treatment. However, the post‐peak downward slopes after R5 for CGR and NArate were coincident between N treatments. Maximum CGR was slightly higher in the full‐N than in the zero‐N treatment (22.4 versus 21.7 g m−2 d−1; P = 0.046), but the growth stage at which the peak of CGR occurred did not differ between N treatments (P = 0.212; Figure 3a, inset). In contrast, maximum NArate was similar between N treatments (7.7 versus 7.4 kg N ha−1 d−1 in the full‐N and zero‐N treatments, respectively; P = 0.180); however, the peak of NArate occurred earlier in the full‐N compared with the zero‐N treatment (P = 0.001; Figure 3b, inset). On average, CGR and NArate between VE‐R5 were 11% and 22% higher in the full‐N compared with the zero‐N treatment, respectively (Supplementary Table S1 and S2). These differences led to an additional 676 kg ADM ha−1 and 45 kg N ha−1 accumulated at R5 in the full‐N versus zero‐N treatment.

FIGURE 3

(a) Soybean aboveground dry matter (ADM), and (b) nitrogen (N) accumulation in the full‐N (red squares) and zero‐N (blue circles) treatments as a function of development stage (DS) on the bottom x‐axis. (c) Same as (b) but also showing accumulated N derived from indigenous soil N (ISN; brown downward triangles) and biological N fixation (BNF; green upward triangles) in the zero‐N treatment. Vn and Rn stages based on Fehr & Caviness (1977) are shown in the top x‐axis. Solid red and blue lines represent the fitted Equation 2 for the full‐N and zero‐N treatment means computed from pooled data across environments. Insets show sigmoid derivative daily crop growth rates (CGR), N accumulation rates (NArate), and rates of ISN and BNF. Coefficient of determination of fitted models was >0.90 in all cases [Colour figure can be viewed at wileyonlinelibrary.com]

(a) n class="Species">Soybean aboveground pan> class="Disease">dry matter (ADM), and (b) nitrogen (N) accumulation in the full‐N (red squares) and zero‐N (blue circles) treatments as a function of development stage (DS) on the bottom x‐axis. (c) Same as (b) but also showing accumulated N derived from indigenous soil N (ISN; brown downward triangles) and biological N fixation (BNF; green upward triangles) in the zero‐N treatment. Vn and Rn stages based on Fehr & Caviness (1977) are shown in the top x‐axis. Solid red and blue lines represent the fitted Equation 2 for the full‐N and zero‐N treatment means computed from pooled data across environments. Insets show sigmoid derivative daily crop growth rates (CGR), N accumulation rates (NArate), and rates of ISN and BNF. Coefficient of determination of fitted models was >0.90 in all cases [Colour figure can be viewed at wileyonlinelibrary.com]

On average, indigenous sn class="Chemical">oil N anpan>d BNF accounted for 33% anpan>d 67% of the accumulated N at R7 stage (Figure 3c, Table 2). However, relative contribution of inpan>digenous span> class="Chemical">oil N and BNF to plant N accumulation varied during the season. For example, indigenous soil N accounted for the largest portion of accumulated N between VE and R5 stages (65%), while BNF supplied most of the N during the seed filling (90%). Maximum rates of N accumulation from indigenous soil N and BNF occurred around R3 (DS = 1.07) and R5.5 stages (DS = 1.56), respectively (Figure 3c, inset).

Accumulated n class="Chemical">N inpan> the full‐pan> class="Chemical">N treatment can be taken as a measure of plant N demand when N supply is not limiting. Comparison of rates of accumulated N in full‐N treatment versus BNF in the zero‐N treatment is of interest as to discern when (and the extent to which) the combined N supply from BNF and indigenous soil N was not sufficient to meet plant N demand. Daily rates of plant N demand and BNF increased gradually during the season until reaching a peak, declining subsequently during the seed filling (Figure 3c, inset). However, the peak of BNF rate occurred later compared with plant N demand (DS = 1.54 versus 1.28) and daily BNF was consistently lower than the plant N demand during the entire crop cycle until the middle of the seed filling phase. During the same period, NArate was consistently lower in the zero‐N than in the full‐N treatment, suggesting that BNF was not sufficient to fill in the ‘N gap’ between plant N demand and indigenous soil N supply.

Drivers for differences in seed number between N treatments

Seed number was positively associated with accumulated n class="Chemical">ADM between the R3 anpan>d R6 stages (Figure 4a). No differences inpan> slope or inpan>tercept of the relationship between seed number anpan>d accumulated pan> class="Chemical">ADM during R3‐R6 period were detected between the two N treatments (P > 0.60), so N treatment data were pooled. Inferentially, this means that number of seeds set per unit of accumulated ADM between R3 and R6 remained unchanged between N treatments. Accumulated ADM between R3 and R6 was greater in the full‐N than in the zero‐N treatment (816 versus 785 g m−2, respectively; paired T‐test, P = 0.045) (Figure 4b), leading to higher seed number in the full‐N compared with the zero‐N treatment (Tables 2 and 3). Differences in R3‐R6 duration between N treatments was not significant (P = 0.356), indicating that differences in accumulated ADM during the R3‐R6 phase between N treatments was associated with differences in CGR (Supplementary Table S1, and S2).

FIGURE 4

(a) Relationship between seed number and accumulated aboveground dry matter (ADM) between the beginning of pod setting (R3) to full seed (R6) stages in the full‐N (red squares) and zero‐N (blue circles) treatments. Connecting lines relate paired treatments from the same environment. Parameters of the fitted linear regression (solid black line) and coefficient of determination (r) are shown. (b) Comparison of accumulated ADM between the R3–R6 phase between the full‐N and zero‐N treatments for the seven environments. Solid black line labeled 1:1 indicates y = x null hypothesis of no difference [Colour figure can be viewed at wileyonlinelibrary.com]

(a) Relationship between seed number and accumulated aboveground n class="Disease">dry matter (pan> class="Chemical">ADM) between the beginning of pod setting (R3) to full seed (R6) stages in the full‐N (red squares) and zero‐N (blue circles) treatments. Connecting lines relate paired treatments from the same environment. Parameters of the fitted linear regression (solid black line) and coefficient of determination (r) are shown. (b) Comparison of accumulated ADM between the R3–R6 phase between the full‐N and zero‐N treatments for the seven environments. Solid black line labeled 1:1 indicates y = x null hypothesis of no difference [Colour figure can be viewed at wileyonlinelibrary.com]

Mechanisms explaining differences in accumulated ADM between N treatments

The two n class="Chemical">N treatments differed in terms of their impact on the seasonal dynamics of LAI (Figure 5a). First, LAI early in the season was greater in the full‐n class="Chemical">N than in the zero‐N treatment, which led to greater leaf area duration before R5 (i.e., the integral of the LAI curve before R5; P = 0.010). Second, while the maximum LAI value was nearly identical in the full‐N and zero‐N treatments (5.8 and 5.7, respectively; P = 0.804), the LAImax value was reached sooner in the full‐N compared with the zero‐N treatment (1.22 versus 1.37 DS units, respectively; P = 0.007). Third, the rate of (post‐max) decline in LAI values was greater in the full‐N than in the zero‐N treatment, though LAI trends of the two N treatments converged towards the end of the season.

FIGURE 5

(a) Soybean leaf area index and (b) estimated fraction of absorbed photosynthetically active radiation (fAPAR) in the full‐N (red squares) and zero‐N (blue circles) treatments as a function of developmental stage (DS). Stages based on Fehr & Caviness (1977) are shown in the top x‐axis. Solid lines represent the fitted models for the full‐N (red) and zero‐N treatment (blue) based on the pooled data across environments. Data for each N treatment were pooled across environments. See Section 2.4 for detailed explanation on fAPAR estimation [Colour figure can be viewed at wileyonlinelibrary.com]

(a) n class="Species">Soybean leaf area inpan>dex anpan>d (b) estimated fraction of absorbed photosynthetically active radiation (fAPAR) inpan> the full‐N (red squares) anpan>d zero‐N (blue circles) treatments as a function of developmental stage (pan> class="Chemical">DS). Stages based on Fehr & Caviness (1977) are shown in the top x‐axis. Solid lines represent the fitted models for the full‐N (red) and zero‐N treatment (blue) based on the pooled data across environments. Data for each N treatment were pooled across environments. See Section 2.4 for detailed explanation on fAPAR estimation [Colour figure can be viewed at wileyonlinelibrary.com]

Seasonal patterns in estimated fAPAR followed the observed dynamics in LAI (Figure 5b). The full‐N treatment exhibited consistently greater fAPAR between VE and R5, reaching ca. 95% of full interception 3 days earlier (P < 0.001) compared with the zero‐N treatment. These differences resulted in larger total estimated cumulative APAR from VE to R7 stages in the full‐N than in the zero‐N treatment (700 versus 679 MJ m−2; P = 0.014).

Canopy‐level radiation‐use efficiency was significantly greater in the full‐N compared with the zero‐N treatment in six of the seven environments (Figure 6). The RUE was 8% greater in the full‐N than in the zero‐N treatment (2.00 versus 1.86 g MJ−1, respectively). In contrast, we could not detect significant differences in leaf‐level photosynthesis (i.e., net n class="Chemical">CO2 assimilation rate) between N treatments (Figure 7). Crop stage (i.e., pan> class="Chemical">DS) significantly influenced Amax and Ic (P < 0.04), with Amax increasing until reaching a maximum (40.8 μmol CO2 m−2 s−1) at ca. R5 (DS = 1.3) stage, but with Ic and Rd gradually decreasing from the first sampling date just prior to stage R1 onward to the last sampling date at R6 stage. Changes in α were small, without a clear pattern. Overall, the above four photosynthetic parameters, and also leaf N concentration and specific leaf weight (SLW), were not statistically significant between N treatments (P > 0.1), but the full‐N treatment tended to have greater leaf N concentration and SLW compared with the zero‐N treatment (Supplementary Figure S4).

FIGURE 6

Aboveground biomass plotted against cumulative absorbed photosynthetically active radiation (APAR) in the full‐N (red squares) and zero‐N (blue circles) treatments across the seven environments. The upper left panel shows the pooled data from the seven environments. Slope of the fitted linear regression represents the radiation‐use efficiency (RUE; g MJ−1). Significance of the statistical test for the null hypothesis of no difference between slopes between N treatments is shown. Asterisks indicate significance at *P < 0.10, **P < 0.05, or ***P < 0.001 [Colour figure can be viewed at wileyonlinelibrary.com]

FIGURE 7

Leaf net photosynthesis (A) as a function of photosynthetic photon flux density (PPFD) in the full‐N (red squares) and zero‐N (blue circles) treatments at four development stages (DS). Parameters of the fitted models are shown: maximum photosynthesis in μmol CO2 m−2 s−1 (Amax), light compensation point in μmol m−2 s−1 (Ic), initial slope of light response curve (α), and dark respiration in μmol CO2 m−2 s−1(Rd). Data for each N treatment were pooled across environments [Colour figure can be viewed at wileyonlinelibrary.com]

Leaf net photosynthesis (A) as a function of photosynthetic photon flux density (PPFD) in the full‐N (red squares) and zero‐N (blue circles) treatments at four development stages (n class="Chemical">DS). Parameters of the fitted models are shown: maximum photosynthesis inpan> μmol pan> class="Chemical">CO2 m−2 s−1 (Amax), light compensation point in μmol m−2 s−1 (Ic), initial slope of light response curve (α), and dark respiration in μmol CO2 m−2 s−1(Rd). Data for each N treatment were pooled across environments [Colour figure can be viewed at wileyonlinelibrary.com]

Drivers for differences in seed mass and seed N between N treatments

Rate of seed n class="Disease">dry matter accumulation during the seed filling period (R5‐R7 phase) was greater in full‐N than in the zero‐N treatment (P = 0.032). There was a strong relationship between mobilized N from non‐seed n class="Chemical">ADM to seed and the amount of N in non‐seed ADM at R5 stage (Figure 8a). Differences in mobilized N between the treatments were associated with greater accumulated N at R5 in the full‐N compared with the zero‐N treatment (281 versus 242 kg N ha−1, respectively; P < 0.001), without changes in the fraction of non‐seed N that was mobilized between N treatments (62%; P = 0.525). Indeed, mobilized N was 17% greater in the full‐N than in the zero‐N treatment (176 versus 150 kg N ha−1, respectively, P = 0.001; Table 3, Figure 8b). Apparent dry matter mobilization from non‐seed ADM to seed was smaller than N mobilization (ca. 14% of non‐seed ADM at R5) and not affected by the N treatments (P = 0.283; Table 3).

FIGURE 8

(a) Relationship between N mobilized from non‐seed dry matter to seed during seed filling and the amount of accumulated N in non‐seed dry matter at R5 in the full‐N (red squares) and zero‐N (blue circles) treatments. Parameters of the fitted linear regression (solid black line) and coefficient of determination (r) are shown. (b) Comparison of mobilized N between the full‐N and zero‐N treatments. Solid black line labeled 1:1 indicates y = x null hypothesis of no difference. Parameters of the fitted linear regression (dashed line) and coefficient of determination (r) are also shown [Colour figure can be viewed at wileyonlinelibrary.com]

(a) Relationship between N mobilized from non‐seed n class="Disease">dry matter to seed durinpan>g seed fillinpan>g anpan>d the amount of accumulated N inpan> non‐seed pan> class="Disease">dry matter at R5 in the full‐N (red squares) and zero‐N (blue circles) treatments. Parameters of the fitted linear regression (solid black line) and coefficient of determination (r) are shown. (b) Comparison of mobilized N between the full‐N and zero‐N treatments. Solid black line labeled 1:1 indicates y = x null hypothesis of no difference. Parameters of the fitted linear regression (dashed line) and coefficient of determination (r) are also shown [Colour figure can be viewed at wileyonlinelibrary.com]

DISCUSSION

In highly productive n class="Species">soybean systems, the combinpan>ed N supply from BNF anpan>d inpan>digenous span> class="Chemical">oil N supply cannot meet the crop N requirement as indicated by the differences in NArate between the N treatments before R5 stage (Figure 3b). There was an asynchrony between BNF and plant N demand: BNF lagged behind plant N demand (i.e., started and peaked later) and it was not sufficient to meet crop N demand after accounting for the indigenous soil N supply (Figure 3c). Sequential N fertilizer application induced a greater NArate, as indicated by the difference between the full‐N and zero‐N treatments, that in turn supported faster leaf area expansion which, combined with greater fAPAR and CGR, allowed a larger seed set during the R3‐R6 phase (Figures 3, 4, and 5; Supplementary Table S1 and S2). Similarly, the extra N accumulated in the ADM before R5 in the full‐N compared with the zero‐N treatment (+45 kg N ha−1) resulted in greater N mobilization from non‐seed ADM to seed (+26 kg N ha−1) (Figure 8, Tables 2 and 3). This is consistent with Sinclair, Farias, Neumaier, and Nepomuceno (2003) and Sinclair & Rufty (2012), who postulated that the degree to which the plant can supply N to meet seed requirements will determine its capacity to meet the potential seed filling rate. These results also suggest that strategies to increase N supply in soybean should aim to increase NArate before seed filling, with the goal of increasing the CGR during the critical period for seed number determination, and increase the amount of N in non‐seed ADM at R5 to support a greater seed mass and seed N concentration. Finally, leaf‐level photosynthesis did not differ among N treatments. In contrast, canopy‐level RUE was greater in the full‐N compared with the zero‐N treatment, which might be associated with changes in root‐to‐shoot ratio, reduced costs due to lower BNF as a result of N fertilizer application, and/or variation in leaf N distribution within the canopy (Bonelli & Andrade, 2020; Cassman, Whitney, & Stockinger, 1980; Pate & Layzell, 1990).

The zero‐n class="Chemical">N treatment seed yield varied withinpan> a narrow ranpan>ge (from 5.3 to 5.8 pan> class="Chemical">Mg ha−1) suggesting that high‐yield soybean relying exclusively on N supply from soil and BNF has an upper yield limit near 5.5 Mg ha−1 in Nebraska (Figure 2, Table 2). In contrast, in four of the seven environments, seed yield in the full‐N treatment reached ca. 6.5 Mg ha−1 because of a simultaneous increase in both seed number and mass. In the other three environments, the full‐N treatment yields were lower (5.5–5.7 Mg ha−1), and thus closer to the zero‐N treatment yields, which was attributable to the high indigenous soil N supply (ca. 100–150 kg N ha−1) at these three site‐years as documented previously in Cafaro La Menza et al. (2019) together with the relatively lower yield potential at these three site‐years (Table 1). Similarly, it was remarkable also to find a simultaneous increase in seed yield and seed N concentration given the well‐documented trade‐off between these two variables (Chung et al., 2003). The degree to which seed N demand is met in soybean depends on (a) NArate during seed filling (R5‐R7), which includes N from soil and BNF, and (b) N mobilized from non‐seed ADM to the growing seeds (Egli et al., 1985; Stephenson & Wilson, 1977; Yamagata et al., 1987). In our study, the NArate differed between N treatments, except for the R5‐R7 seed filling phase. In turn, the amount of mobilized N was greater in the full‐N compared with the zero‐N treatment as a result of higher accumulated N in the non‐seed ADM at R5, with no change in the fraction of N mobilized compared with the zero‐N treatment. Clearly, the extra accumulated N in non‐seed ADM at R5 in the full‐N versus zero‐N treatment should not be seen as ‘luxury N consumption’ as it helped increase seed mass and maintain seed N concentration (Staswick, 1994). Even when soil N was non‐limiting in the full‐N treatment, greater remobilized N (rather than NArate) was the underlying mechanism explaining greater seed mass and seed N concentration compared with the zero‐N treatment. These results are consistent with the CROPGRO model (Boote, Jones, Hoogenboom, & Pickering, 1998) in which mobilized N from non‐seed organs to seed is calculated as a function of thermal time and unaffected by soil N availability during the seed filling. This preference for remobilizing N from non‐seed ADM may be associated with the lower cost of protein breakdown and re‐synthesis compared with de novo protein synthesis from indigenous soil N supply and BNF (De Vries, Van Laar, & Chardon, 1983). It may also indicate a loss of root functionality during seed filling as has been reported for other crop species (Lisanti, Hall, & Chimenti, 2013; Thibodeau & Jaworski, 1975). Understanding the underlying drivers behind the apparent limited capacity of the plant to make use of available soil N during seed filling deserves further research.

In a global context, crop yieln class="Chemical">ds need to inpan>crease ca. 50% by year 2050 inpan> order to meet food demanpan>d anpan>d avoid a massive conversion of natural ecosystems inpan>to croplanpan>d (Cassmanpan> & Grassinpan>i, 2020). For pan> class="Species">soybean, the important issue is the degree to which N limitation may or not allow such a yield achievement. Given the current average yield of ca. 3 Mg ha−1 in major producing areas, the goal would have to be an average field yield of 4.5 Mg ha−1 by year 2050, which entails a N requirement of 360 kg N ha−1. Indigenous soil N supply may be able to cover half of that N requirement considering that fertile agricultural soils in USA and Argentina can provide ca. 100–150 kg N ha−1 (Cafaro La Menza et al., 2019). The remaining N requirement (ca. 200–250 kg N ha−1) can be readily supplied by BNF (Salvagiotti et al., 2008). Thus, adoption of improved agronomic practices and cultivars can be cost‐effective options to increase yields in these intermediate‐yield production systems (e.g., Rattalino Edreira et al., 2017; Di Mauro et al., 2018). However, much higher BNF would be required to sustainably advance yield gains in high‐yield irrigated soybean production areas such as Nebraska and other areas in the Central US Great Plains where producers are already achieving average yields of ca. 4.5 Mg ha−1. In these environments, a 50%‐yield increase would imply an average yield of 6.8 Mg ha−1 by year 2050, which has an associated N requirement of 540 kg N ha−1. Assuming the same level of indigenous soil N (150 kg N ha−1), BNF would then have to increase to 400 kg N ha−1. None of the recent reviews on BNF on soybean provide evidence that reaching such a high level of BNF and indigenous soil N supply is possible (Ciampitti & Salvagiotti, 2018; Salvagiotti et al., 2008). Our N fertilizer treatment was successful at increasing both yield and seed N concentration. However, it was far from being cost‐effective and obviously not an environmentally sound practice to be adopted in commercial farms. Increasing BNF and/or indigenous soil N supply, improving the synchrony between N fixation and plant N demand, and alleviating the trade‐off between the two sources of N are avenues worth exploring, even though the associated probability of success and timeline for impact are unknown (Denison & Kiers, 2005; Giller & Cadisch, 1995; Van Kessel & Hartley, 2000).

CONFLICT OF INTEREST

The authors declare that there is no conflict of interest regarding the publication of this article.

Appendix S1: Supporting Information

Click here for additional data file.