Development of Algorithm for Determining N Fertiliser Requirements of Winter Wheat Based on N Status Using APSIM Modelling

Abstract

The determination of optimum nitrogen (N) fertilisation rates, which maximise yields and minimise N losses, remains problematic due to unknown upcoming crop requirements and near-future supply by the soil. Remote sensing can be used for determining the crop N status and to assess the spatial variability within a field or between fields. This can be used to improve N fertilisation, provided that the optimal fertilisation rate at the time of fertiliser application for an expected yield is known. Using the APSIM-wheat model, we developed an algorithm that relates the N status of the plants at early development stages to the yield response to N. Simulations were performed for winter wheat under growth conditions in Denmark. To obtain a range of different N status in the biomass at early growth stages, the soil N in autumn was varied from 20 to 180 kg N ha⁻¹, and at BBCH23, fertiliser was applied at a rate of 50 kg N ha⁻¹. In a full factorial setup, additional N fertiliser was applied ranging from 0 to 150 kg N ha⁻¹ during three different development stages (BBCH30, 32, and 37). The algorithm was evaluated by comparing model outputs with a standard N application of 50 kg N ha⁻¹ at BBCH23 and 150 kg N ha⁻¹ at BBCH30. The evaluation showed that, depending on the N status of the soil, the algorithm either provided higher or lower optimal N fertilisation rates when targeting 95% of the maximum yield, and these affected the grain yield and the grain N, as well as the amount of N leaching. Split application of fertiliser into three applications was generally beneficial, with decreased product-related N leaching of up to nearly 30%. Further testing of the model under different environmental conditions is needed before such an algorithm can be used to guide N fertilisation.

1. Introduction

Nitrogen fertilisation rates for grain crops are typically based on expected yields and respective N demands, as anticipated from average environmental conditions on a given site. Recommendations either provide a single value for a crop or vary depending on climatic conditions, soil type and N credits from preceding crops (e.g., grain legumes or service crops), and manure applications [1]. However, such blanket fertiliser rates over areas with similar soils and climates do not account for differences in soil N supply, which generally varies across years and different fields, and even within fields. This variability is due to variations in soil properties and differences in mineralisation potential from soil organic matter, manure and crop residues, leftover N from previous crops, and uneven application of fertilisers [2,3,4,5,6]. The recommendations of fertiliser applications without considering the potential supply of N from the soil can pose an environmental risk, as well as economic [7] and societal costs [8].

Nitrate leaching from agricultural soils is of considerable global environmental concern [9,10]. In Denmark, strict fertiliser regulations and various mitigation measures have been implemented to improve N use efficiency and reduce N losses to aquatic environments. Despite this, N loads still exceed the Water Framework Directive thresholds, especially in coastal areas and vulnerable groundwater bodies [11]. Estimated N leaching across Denmark varies widely, depending on crop sequence, soil characteristics, annual percolation, the residual N left after harvest, and the soil mineralisation rate. Zhao et al. [12] reported variations in N leaching from 3 to 92 kg N ha⁻¹ based on field experiments from 44 site x years under optimal N fertilisation rates. Similarly, using the Daisy model with 20 different crop rotations, Rashid et al. [11] found substantial differences in simulated N leaching under Danish conditions, varying from 16 to 85 kg N ha⁻¹. Winter wheat in Denmark is mostly grown on the more fertile soils [13], with common N fertilisation rates around 200 kg N ha⁻¹; a applied in two split applications in the middle of March and April. Average grain yields and harvested N in grain of 8.74 t ha⁻¹ and 154 kg N ha⁻¹ have been reported based on measurements on farmers’ fields across Denmark between 2010 and 2015 [13]. In a two-year fertiliser study with winter wheat, Rasmussen et al. [14] measured around 40 kg residual N before sowing (with pre-crops of oats and winter barley), and similar amounts after the harvest of the winter wheat fertilised with 150 kg N ha⁻¹. Increasing the N rate to 250 kg N ha⁻¹ substantially increased the amount of soil N at harvest. Thus, accounting for the previous fertilisation practice and cropping sequence is important for reducing N losses to the environment, especially after the breakup of grasslands substantial N is released during mineralization of crop residues [9,15].

The synchronisation of crop N demand and supply by soil in both time and space is the most effective way to increase N fertilisation responses and reduce N losses. Various approaches have been tested for refining N fertilisation, including measurements of soil mineral N or potential mineralisation rates and N budgets [16,17]. The critical N concentration curve approach has also been used for refining in-season N fertilisation rates, which is based on critical crop N content as a function of crop biomass, which results in optimal growth [18].

Remote sensing is increasingly being used as a timely and nondestructive tool for mapping the N nutrition status of plants and to rapidly assess the spatial variability within a field based on the canopy reflectance. In-season measurements of the crop N status, linked with local production information, offer promise for fine-tuning N fertilisation rates and are increasingly being used for maize and wheat [19]. A drawback of these approaches is the need for an expected yield. When considering spatial, in-field variability, expected yield differences can be obtained from historical crop maps [20]. To include year-to-year variability, Raun et al. [21] developed an approach for an in-season estimation of yield based on the normalised difference vegetation index NDVI and growing degree days. By using this approach, the N use efficiency of winter wheat could be increased by 15%. Various sensor-based algorithms have been developed to estimate the N requirements based on the in-season crop N status. A detailed review of these has been performed by Franzen et al. [22], who showed that these algorithms depend on user-specified optimum N rates, and either require a high-N reference strip or a virtual reference, as in the Holland–Schepers algorithm. The virtual reference is based on a statistical approach, with a frequency distribution to identify crops in more fertile parts of a field with adequate N.

The use of dynamic simulation models, which simulate soil and crop processes, for adjusting the N fertiliser rates based on simulated soil or plant N has also been shown to be promising [23,24]. Such dynamic simulation models have also been used in conjunction with sensor technology [25], but they are not designed as decision-support tools for supporting in-season fertilisation management. Specific simulation models for in-season N management have been developed [26], and Cichota et al. [27] used the APSIM model to develop an algorithm for fine-tuning the N requirements for grasslands based on plant N status and potential growth. A similar approach would also be relevant for winter wheat and other cereals, and such an algorithm could then be coupled with remote sensing and potentially be developed into a decision support tool to optimise fertiliser management and to reduce production cost and environmental risk.

The objective of this simulation study is to develop an algorithm that can guide N fertilisation in winter wheat depending on the plant N status in spring and evaluate the algorithm regarding grain yield and area- and product-based N leaching. The algorithm was developed with the APSIM model, which has previously been calibrated for the environmental conditions in Denmark regarding phenology, grain yield, N uptake and biomass development under 13 different fertiliser management strategies.

2. Materials and Methods

The simulations were conducted using the APSIM modelling framework (version 7.10). APSIM is a process-based deterministic crop model that simulates crop phenology, crop growth and development, and carbon (C) and N dynamics in the soil and plants at a daily time step as a function of climatic, agronomic, and soil characteristics inputs. The key APSIM modules used in this study were SoilWat for simulating water movement and SurfaceOM and SoilN, which simulate the dynamics of N and C, with manager scripts accounting for management such as sowing, harvesting and fertiliser application. For simulating winter wheat, the cultivar ‘Dan_winter’ was used, which has been calibrated and evaluated by Kumar et al. [28] based on field data from seven locations across Denmark, five years, two sowing dates, and with 7 to 13 fertiliser treatments. Field data included the phenology, N status of the biomass during early phenological stages, and grain yield and grain N at harvest maturity. A detailed description of the plant process and parameters involved in the simulation of C and N dynamics in the APSIM wheat model is available online (https://www.apsim.info/wp-content/uploads/2019/09/WheatDocumentation.pdf; accessed on 24 July 2023).

APSIM simulations were set up for climate and growth conditions at the Flakkebjerg location, Denmark, for 2018/2019. The soil at the site is a sandy loam, with the setup of the soil profile characteristics as provided by Kumar et al. [28]. The simulations were initialised in March 2018, and a generic crop (Canola) was sown to adjust to environmental conditions, thus reducing the effects of the initial conditions. No N fertiliser was applied to the canola crop, which was harvested at the end of August. The winter wheat was sown according to common practice at the end of September (20 September 2018). Under common agricultural practices, the N fertilisation rate to winter wheat is around 200 kg ha⁻¹, with mineral N fertiliser surface-applied in two split applications in the middle of March (50 kg N ha⁻¹) and April (150 kg N ha⁻¹). Recently, a split application with three application timings has also been promoted to reduce environmental impacts due to overfertilization with N. Thus, different fertilisation schemes were set up within APSIM.

For the development of the algorithm, fertilizer rates ranged from 0 to 250 kg N ha⁻¹ (interval of 50 kg ha⁻¹), and these were applied at the BBCH (a phenological stage, see [22] scales of 30, 32, and 37. Five different algorithms (Algorithm 1 to 5) were developed, which differed in the timing (BBCH stage) of fertiliser application (

Table 1. Fertiliser scenarios for developing the fertiliser algorithm, with N_r = N fertilisation rate (kg N ha⁻¹) and BBCH = phenological stage of the wheat according to Lancashire et al. [30]. Nmin are the initial values of soil mineral N.

Fertilisation Scheme	Nr BBCH23	Nr BBCH30	Nr BBCH32	Nr BBCH37
Algorithm 1	50	0–250	-	-
Algorithm 2	50	-	0–250	-
Algorithm 3	50	-	-	0–250
Algorithm 4	50	50	0–250	-
Algorithm 5	50	100	-	0–250

). Additionally, each simulation received a single dose of 50 kg N ha⁻¹ at BBCH23. To obtain a range of plant N status during the early growth stages (at the time when decisions regarding N fertilisation rates are made), the mineral soil N at the time of sowing in autumn was varied between 20 and 180 kg N ha⁻¹ (with an interval of 20 kg N ha⁻¹). This range in mineral soil N has been chosen to obtain realistic N uptake rates, which align with plant N uptake amounts measured spring across various locations in Denmark [29]. The soil organic carbon and nitrogen were constant (at 1.5% organic carbon) across all simulations and did not change with the initial soil mineral N. This implies that the difference in mineral soil is only due to leftover fertiliser from the previous year and not due to differences in organic matter mineralisation. Within a full factorial setup, these initial soil N were simulated with the various fertiliser schemes at the three BBCH stages, resulting in a combination of 3087 simulations.

The algorithm, a three-dimensional surface response function, is based on the Mitscherlich yield response function:

$$Y=Y_{max}-\left(Y_{max}-Y_0\right)\exp -\beta N_r$$

(1)

where Y is the grain yield (kg DM ha⁻¹), Y_max is the maximum yield under the climatic and edaphic conditions (kg DM ha⁻¹), N_r is the rate of N applied (kg N ha⁻¹), Y₀ is the yield when no external N is applied (N_r = 0), and β is an ‘activity’ coefficient, which is a measure of the availability of the applied nutrient to the crop. Following Vogeler et al. [31], it is assumed that both Y₀ and β are dependent on the N uptake during early development:

$$Y_0=a+b N_{uptake}$$

(2)

$$\beta =c+d N_{uptake}$$

(3)

Simulated grain yield responses (to soil mineral N and fertilisation rates) were fitted to the developed 3D model (Equations (1)–(3)) for the different BBCH stages using Table Curve 3D (v. 4.0; SYSTAT Software Inc., Richmond, CA, USA) and using a Y_max range between 9300 and 9400 kg DM ha⁻¹, as the maximum obtained in the simulations for the algorithm development.

APSIM simulations with different fertilisation management were set up to evaluate the performance of the fertiliser algorithm. The simulations were derived from the simulation described above and comprised five different soil mineral N contents (25, 50, 75, 100, and 140 kg N ha⁻¹) and either used the standard site and crop fertilisation rates or were based on one of the algorithms (

Table 2. Fertilisation rates (N_r; kg ha⁻¹) scenarios for testing the performance of the fertiliser algorithm. For each of the fertilisation schemes, the soil mineral N in autumn at sowing of the winter wheat was either 25, 50, 75, 100 or 140 kg N ha⁻¹.

Fertilisation Scheme	Nr BBCH23	Nr BBCH30	Nr BBCH32	Nr BBCH37
Standard	50	150
Algorithm 1	50	Algorithm 1	-	-
Algorithm 2	50	-	Algorithm 2	-
Algorithm 3	50	-	-	Algorithm 3
Algorithm 4	50	50	Algorithm 4	-
Algorithm 5	50	100	-	Algorithm 5

). The algorithm was used at different BBCH stages to apply the required fertilisation rate while targeting a maximum yield of 95% of Y_max. To avoid excessive N application, the maximum N fertilisation rate at any BBCH stage was limited to 200 kg ha⁻¹.

3. Results

3.1. Yield Response Curves

To obtain yield estimates depending on the N uptake of the wheat during the early development and the fertilisation rate, Equation (1) (with 2 and 3) was fitted to the three-dimensional yield data (yield, N uptake, and fertiliser rate), providing a response surface that can be used to guide N fertilisation, as shown below. The fitting was performed separately for the three different BBCH stages and the five different algorithms, which included single and split applications of N. For each of the BBCH stages, only simulations were used in which the fertilisation rates at the other BBCH stages were the same (e.g., for BBCH32, only simulations with N fertilisation of 0 kg N ha⁻¹ at BBCH30 and BBCH 37 were used). The developed equation describes the response curve reasonably well for both the algorithms with a single N application (Figure 1) and the algorithms with split N fertilisation (Figure 2). Fitted standard errors for the yield range between 167 and 451 kg DM ha⁻¹. The values for the response surface parameters for the different Algorithms and BBCH stages are given in

Table 3. Parameter values for the response surface function algorithms for different BBCH stages and fertilisation rates (N_r; kg ha⁻¹) used in the simulations for the respective fertiliser schemes, Y_max is the maximum yield under the climatic and edaphic conditions (kg DM ha⁻¹).

Fertilisation Scheme	BBCH	Nr BBCH30	Nr BBCH32	Nr BBCH37	Ymax	a	b	c	d
Algorithm 1	30	0–250	-	-	9400	3849	49.6	−1.0 × 10−3	−1.9 × 10−4
Algorithm 2	32	-	0–250	-	9400	3367	48.1	−1.6 × 10−3	−1.9 × 10−4
Algorithm 3	37	-	-	0–250	9397	3000	44.1	−2.2 × 10−3	−1.8 × 10−4
Algorithm 4	32	50	0–250	-	9402	3000	55.7	−1.5 × 10−3	−2.2 × 10−4
Algorithm 5	37	100	-	0–250	9330	3031	42.4	−1.2 × 10−2	−2.7 × 10−4

.

By substitution of Equations (2) and (3) into Equation (1) and rearranging, the amount of fertiliser required for a targeted yield (Y_T) at or below Y_max can then be calculated, and for any N uptake and for the different BBCH stages:

$$N_{r,NuptakeBBCH}=ln \left(\frac{\left(Y_{max}-Y_T\right)}{\left(Y_{max}-\left(a+b N_{uptake}\right)\right)}\right)/\left(c-d N_{uptake}\right)$$

(4)

The fertiliser requirement for 95% of the maximum yield based on the various algorithms for the different BBCH stages is shown in Figure 3, which shows an exponential decrease in N fertiliser requirement with increased N uptake at any BBCH stage. For example, at BBCH30, an N uptake of 60 kg ha⁻¹ requires a N fertilisation of 103 kg N ha⁻¹, while at an uptake of 80 kg N/ha⁻¹, only 38 N ha⁻¹ needs to be applied for achieving 95% of the maximum yield. At any N uptake, the required N fertilisation rate increases with increasing phenological development. For example, at an N uptake of 60 kg ha⁻¹, the required N fertilisation rate at BBCH30, BBCH32, and BBCH37 are 103, 123, and 161 kg ha⁻¹. Split application, as in Algorithm 4 and Algorithm 5, reduces the required N fertilisation rate, as the future supply would be increased, as not all the applied N fertiliser would have been taken up by the crop.

3.2. Algorithm Evaluation

The comparison between simulation results, when N fertilisation at BBCH30 was either based on the standard application (50 kg N ha⁻¹ at BBCH23 and 150 kg N ha⁻¹ at BBCH30) or the N status and using the algorithm (Equation (4), with the values for Algorithm 1 provided in Table 3), shows that at most N status in spring, the algorithm would apply more N (

Table 4. Comparison of using standard fertilisation schedule with fertilisation based on Alg1, in which the fertilisation rate is based on the plant N status at BBCH30. The colour scaling indicates the ranking of each variable, with green and red being the lowest or highest depending on the parameter. N_min is the soil mineral nitrogen (kg N ha⁻¹) in autumn at the sowing of the winter wheat. N uptake (grain and straw), N_r (fertilisation rate), grain N, and N leaching in kg N ha⁻¹, grain yield in kg ha⁻¹.

Fertilisation Scheme	Nmin	Nuptake BBCH30	Nr	N Uptake	Yield	Grain N	N Leaching	Nleach/ kg DM
Standard	25	47	200	165	8389	135	28	0.0033
Algorithm 1	25	47	250	202	9345	163	37	0.0039
Standard	50	52	200	173	8645	141	32	0.0037
Algorithm 1	50	52	250	210	9311	162	42	0.0045
Standard	75	57	200	181	8852	149	36	0.0040
Algorithm 1	75	57	249	216	9265	161	47	0.0051
Standard	100	62	200	189	9072	156	39	0.0043
Algorithm 1	100	62	225	205	9241	161	46	0.0050
Standard	140	70	200	202	9237	161	47	0.0051
Algorithm 1	140	70	189	195	9199	160	43	0.0047

). This higher application results in a higher total N uptake at maturity and an increase in both grain yield and grain N. Leaching is, however, also substantially increased, by 17 to 32% when area-based and by 15 to 27% when product-based (N leaching/kg grain DM). This shows that increased N fertilisation above the standard of 200 kg N ha⁻¹ at stage BBCH30 can increase grain yields, at least in the year simulated. This would, however, come at a high cost regarding N leaching. As in these simulation setups, only the soil mineral N was altered; additional N from increased N mineralisation in areas with higher soil organic C is not accounted for. This would likely increase the N uptake in spring and thus reduce the required N fertilisation rates in such areas further.

Delaying fertilisation until BBCH32 or BBCH37 and using Alg2 or Alg3 similarly resulted in higher N fertilization rates compared to the standard application with higher grain yield, grain N and N leaching (

Table 5. Model outputs for using fertiliser algorithms, in which the fertilisation rate (N_r) is based on the plant N status at different BBCH stages. N_min is the soil mineral N in autumn (kg N ha⁻¹). N uptake (BBCH), N_r, N uptake (grain and straw), grain N, and N leaching are in kg N ha⁻¹; yield is in kg DM ha⁻¹. The number in brackets indicates the BBCH stage for the N uptake. The N fertilisation rate includes the total amount applied at the various BBCH stages.

Nmin	Fert Scheme	Nuptake (BBCH)	Nr	N Uptake	Yield	Grain N	N Leach	N Leach/kg DM	Diff N Leach/kg DM
25	standard	47 (30)	200	165	8389	135	28	0.0033
25	Alg 1	47 (30)	250	202	9345	163	37	0.0039	17%
25	Alg 2	49 (32)	250	203	9177	166	40	0.0044	32%
25	Alg 3	52 (37)	250	215	8264	141	40	0.0048	45%
25	Alg 4	66 (32)	234	190	9214	159	34	0.0037	11%
25	Alg 5	113 (37)	218	178	8831	148	30	0.0034	3%
50	standard	52 (30)	200	173	8645	141	32	0.0037
50	Alg 1	52 (30)	250	210	9312	162	42	0.0045	22%
50	Alg 2	55 (32)	250	209	9243	164	46	0.0050	36%
50	Alg 3	60 (37)	250	219	8907	155	46	0.0051	39%
50	Alg 4	74 (32)	207	179	8825	147	33	0.0037	1%
50	Alg 5	122 (37)	198	172	8615	141	31	0.0036	−2%
75	standard	57 (30)	200	181	8852	149	36	0.0040
75	Alg 1	57 (30)	249	216	9265	161	47	0.0051	27%
75	Alg 2	62 (32)	235	206	9289	162	47	0.0050	25%
75	Alg 3	68 (37)	247	219	9255	161	51	0.0055	37%
75	Alg 4	80 (32)	189	176	8675	144	32	0.0037	−8%
75	Alg 5	128 (37)	183	174	8615	142	30	0.0034	−14%
100	standard	62 (30)	200	189	9072	156	39	0.0043
100	Alg 1	62 (30)	225	205	9241	161	46	0.0050	15%
100	Alg 2	68 (32)	206	192	9244	159	43	0.0047	9%
100	Alg 3	77 (37)	212	202	9309	161	44	0.0048	11%
100	Alg 4	85 (32)	174	176	8661	144	30	0.0035	−18%
100	Alg 5	133 (37)	170	174	8587	141	29	0.0034	−21%
140	standard	70 (30)	200	202	9237	161	47	0.0051
140	Alg 1	70 (30)	189	195	9199	160	43	0.0047	−7%
140	Alg 2	80 (32)	163	180	8758	147	35	0.0041	−20%
140	Alg 3	92 (37)	165	183	8913	150	36	0.0041	−20%
140	Alg 4	93 (32)	149	174	8574	141	31	0.0036	−29%
140	Alg 5	136 (37)	160	179	8736	146	34	0.0039	−23%

). The product-related N leaching was also increased with the use of Algorithms 1, 2 and 3. Targeting a lower maximum yield would obviously reduce the amount of fertiliser applied via these algorithms and consequently yield N leaching. Only at high soil, N_min was the N rate reduced, and this was at the cost of reduced yield and grain N but also reduced N leaching. The use of a more split application of N fertiliser at BBCH30 and BBCH32 based on the N status, as performed using Algorithms 4 and 5, can be beneficial with similar yields with reduced N leaching, especially at high soil Nmin, where product-related reductions in N leaching were up to 29%.

4. Discussion

This simulation study was set up to develop algorithms for guiding N fertilisation rates for winter wheat based on the N uptake at different early development stages in spring. The algorithms were then tested and compared to a standard blanket fertilisation rate. Results indicated that when targeting 95% of the maximum yield, the use of the algorithms resulted in higher N application rates, with higher yields, but also increased N leaching. Using a lower yield target would obviously decrease the N fertilisation rate and also N leaching but decrease the yield. The use of a split application with three applications resulted in a substantial decrease in product-related N leaching at medium to high soil mineral N. This shows that the developed algorithms is promising for improving N fertilisation and controlling N losses in cereal production systems. So far, the algorithms have only been developed for one year and one location, and it needs to be tested and parameterised for other environments as well as tested in the field. It should also be noted that the variation in crop growth here was limited to N supply, while in-field variation can also be due to other constraints, such as water availability and other macro- and micro-nutrients—although various studies have shown the latter to be of minor importance for the often observed large in-field yield variation [32].

The developed algorithms can be combined with remote sensing of the N status and developed into a decision support system. A dataset, including various winter wheat trials conducted in different locations in Denmark [29] and with a range of N fertilisation, shows a good correlation between the drone-estimated NDRE (Normalised Difference Red Edge) index at growth stages 28–57 and N uptake (Figure 4), with an R² of 0.62. As the NDRE saturates at about 0.6, this relationship should only be used up to a N uptake of about 75 kg N ha⁻¹. Such information obtained from remote sensing can then be used to calculate the amount of fertiliser required for a targeted yield (Y_T) using the developed algorithm (Equation (4)), with the parameter values for the corresponding BBCH stage (Table 3). The targeted yield can be the maximum, but also the one corresponding to the economic optimum, the fertilisation rate at which the marginal cost of the fertilisation corresponds to the marginal revenue [33]. This would also comply with the maximum allowable N fertilisation rates in Denmark [34]. So, for example, an NDRE value of 0.4 equates to a N uptake of 34.5 kg ha⁻¹ (Figure 4).

If such an NDRE value of 0.4 would be measured at BBCH30, the N fertilisation rate required to achieve 95% of the maximum yield would be 285 kg ha⁻¹ (Figure 3). Similarly, an NDRE value of 0.6 would equate to a N uptake of 75.2 kg ha⁻¹. If this were to be measured at BBCH30, BBCH 32 or BBCH37, the required N rates would be 52.5 kg ha⁻¹, 74.5 and 104.9, respectively (Figure 3). The algorithms could also be used to identify the lowest product-related N leaching. Especially when the residual soil N in the autumn, prior to sowing of the winter wheat, is high, a split application of N fertiliser at BBCH30 and BBCH32 based on the N status, as performed using Algorithms 4 or 5, can be beneficial with similar yields and low N leaching, resulting in reductions of product related N leaching of up to 29% compared with the standard application rate.