Predicting Rice Grain Yield Based on Dynamic Changes in Vegetation Indexes during Early to Mid-Growth Stages

Predicting the grain yield during early to mid-growth stages is important for initial diagnosis of rice and quantitative regulation of topdressing. In this study, we conducted four experiments using different nitrogen (N) application rates (0–400 kg N·ha−1) in three Japonica rice cultivars (Wuyunjing24, Ningjing4, and Lianjing7) grown in Jiangsu province, Eastern China, from 2015–2016. Spectral reflectance data were collected multiple times during early to mid-growth stages using an active mounted sensor (RapidScan CS-45, Holland Scientific Inc., Lincoln, NE, USA). Data were then used to calculate optimal vegetation indexes (normalized difference red edge, NDRE; normalized difference vegetation index, NDVI; ratio vegetation index, RVI; red-edge ratio vegetation index, RERVI), which were used to develop a dynamic change model and in-season grain yield prediction model. The NDRE index was more stable than other indexes (NDVI, RVI, RERVI), showing less standard deviation at the same N fertilizer rate. The R2 of the relationships between leaf area index (LAI), plant nitrogen accumulation (PNA), and NDRE also increased compared to other indexes. These findings suggest that NDRE is suitable for analysis of paddy rice N nutrition. According to real-time series changes in NDRE, the resulting dynamic model followed a sigmoid curve, with a coefficient of determination (R2) >0.9 and relative root-mean-square error <5%. Moreover, the feature platform value (saturation value, SV) of the NDRE-based model accurately predicted the differences between treatments and the final grain yield levels. R2 values of the relationship between SV and yield were >0.7. For every 0.1 increase in SV, grain yield increased by 3608.1 kg·ha−1. Overall, our new dynamic model effectively predicted grain yield at stem elongation and booting stages, providing real-time crop N nutrition data for management of N fertilizer topdressing in rice production.


Introduction
Rice (Oryza sativa L.) is a staple food for billions of people worldwide.China is the largest rice producer in the world, with approximately 31 million hectares of paddy rice cultivation [1].To ensure sustainable food production, agricultural producers invest in a significant amount of nitrogen (N) fertilizer [1].Zhang et al. reported that the N use efficiency (NUE) of China's major food crops was less than 30% [2], representing 33% of the global NUE, and 40% of that in developed countries in Europe and America [3].Scientific management of N fertilizer use is, therefore, crucial in realizing the goals of high yield, quality, and efficiency in rice production in eastern China.confirmed a relatively stable relationship between NDVI and grain yield in double-cropping rice during major growth periods, with R 2 values ranging from 0.6 to 0.65 [25].
Xue et al. also described the relationship between NDVI and N accumulation from tillering to early grain filling stages, and used the obtained sufficiency index (SI) to predict potential yield and guide topdressing nitrogen management [26].The R 2 value was 0.8, resulting in effective improvements in NUE.Furthermore, Evert et al. studied the relationship between the weighted difference vegetation index (WDVI) and the N status of potato using MSR87 (CropScan Inc., Rochester, MN, USA) [27], revealing a significant correlation between WDVI and N uptake.Accordingly, they were able to establish a potential yield model for potato.Similarly, Morier et al. used FieldSpec HandHeld (Analytical Spectral Devices (ASD), Inc., Boulder, CO) to obtain hyperspectral data (325-1075 nm and a spectral resolution of 3.5 nm at 700 nm) for monitoring of potato populations, revealing that the best relationship existed between the index incorporating the red-edge wave (CIred-edge) and yield.As a result, they were able to explain 76% of the variability in total tuber yield at 55 days after planting [28].
Lukina et al. and Raun et al. also proposed an N fertilization optimization algorithm based on spectral index and the corresponding time of the growing degree days (GDD), which allowed calculations of the estimated yield coefficient in-season estimated yield (INSEY) [29,30].This method eliminated the effect of different growth stages, compared with approaches based on vegetation indexes, and improved accuracy by 30% in terms of grain yield predictions in winter wheat [31].However, studies also showed that NDVI and other spectral indexes have a significant saturation effect under high vegetation coverage [32].In rice, accuracy decreases when the plant dry matter is greater than 3736 kg•ha −1 [33].As a result, the prediction accuracy of the algorithm also decreased, with an R 2 of only 0.5 [34].Moreover, in rice crops in particular, R 2 was only about 0.4-0.5 at stem elongation and booting stages [8,35].However, Thompson et al. revealed that red-edge bands based on normalized difference red edge (NDRE) could effectively identify crop populations, with a higher correlation with indicators such as plant nitrogen accumulation, helping solve the problem of saturation [36].Cao et al. used multispectral data of Crop Circle-ACS470 (Holland Scientific, Lincoln, NE, USA) to fit a variety of vegetation indexes including the normalized NIR index (NNIR), green soil-adjusted vegetation index (GSAVI), and modified enhanced vegetation index (MEVI); then, based on these indexes, they determined the relationship between INSEY and rice yield at stem elongation and booting [37].These results revealed an R 2 value of 0.8 at stem elongation, 21-26% higher than that based on INSEY-NDVI and INSEY-RVI [38].However, few studies examined the use of NDRE-related models in accurate management of crop N fertilizer topdressing.
Previous studies suggested problems of saturation with models using vegetation indexes such as the single spectral index, resulting in low prediction accuracy when used to predict nutrition or grain yield in paddy rice [38].Screening of effective vegetation indexes and development of a stable accurate method for predicting rice grain yield is, therefore, required.This paper, therefore, aimed to (1) evaluate spectral index data obtained using the active multi-spectrometer RapidScan CS-45 (Holland Scientific Inc., Lincoln, NE, USA) during rice growth; (2) construct a dynamic model using the vegetation index data; and (3) develop a prediction algorithm of rice grain yield based on the dynamic model.

Experiment Design
Japonica rice is the main variety grown in Jiangsu Province, which has the highest average rice production in all China [2].This study was carried out in two typical Japonica rice-producing areas in Jiangsu Province, the lower reaches of the Yangtze River (Rugao) and the Huaihe River Basin (Huai'an and Sihong) in northern Jiangsu.In Rugao, the predominant soil type is loam, with a total N content of 1.71 g•kg −1 , olsen phosphorus (P) content of 13.3 mg•kg −1 , and available potassium (K) content of 95.7 mg•kg −1 .The soil types in Huai'an and Sihong are yellow-brown and lime concretion black soil, respectively, with a total N content of 1.35 and 1.28 g•kg −1 , olsen P content of 32 and 27.6 mg•kg −1 , and available K content of 85.3 and 75.2 mg•kg −1 , respectively.The climate data of each experimental site (Exp.1, Exp. 2, and Exp. 3) is shown in Figure 1.Four experiments were carried out in total using different N application rates (0-400 kg•N•ha −1 ) and three Japonica rice varieties, Wuyunjing24 (WYJ24), Ningjing4 (NJ4), and Lianjing7 (LJ7), as detailed in Table 2.

Data Acquisition and Determination
Spectral reflectance data were collected using the active portable sensor RapidScan CS-45 (Holland Scientific Inc., Lincoln, NE, USA), and contained red (R, 670 nm), red-edge (Re, 730 nm), and near-infrared (NIR, 780 nm) wavebands.It has its own light source; thus, dependency on sunlight was avoided, preventing sensitivity to weather conditions.All spectral measurements were performed during cloud-free periods between 9:30 and 10:30 a.m. with the sensor placed 0.7 m above the canopy.Three representative lines were selected in each plot.The sensor records one spectral reflectance value every second along the line at uniform velocity, giving approximately 70-80 values per line.The average spectral index and reflectance in each plot was then calculated.
Exp. 1 and Exp. 3 were carried out from mid-tillering, with measurements every seven days until the end of tillering to flowering, when testing was conducted every three days.In Exp. 2 and Exp. 4, measurements were carried out simultaneously with sampling at each growth stage.
In each experiment, three hills from each plot were sampled for growth analysis at different stages during vegetative growth.Plants were manually uprooted and then cut at ground level for determination of N concentration.Fresh plants were separated into green leaf blade (leaf) and culm plus sheath (stem) samples, heated for 30 min at 105 • C to halt metabolic processes, and then dried at 80 • C in a forced-draft oven until reaching a constant weight.Plant dry matter (PDM) and leaf dry matter (LDM) were then determined before grinding the samples and passing them through a 1-mm sieve in a Wiley mill.The green fresh leaf samples were immediately scanned using an LI-3000A (Li-Cor, Lnr, Lincoln, NE) and the green leaf area of each layer was obtained to calculate the green LAI for each plot (the sum of different layers).The samples were then stored in plastic bags at room temperature until further chemical analysis.Whole-plant and leaf N concentrations (PNC and LNC) were determined using an elemental analyzer (vario MACRO cube; Elementar, Hanau, Germany); then, plant N accumulation (PNA) and leaf N accumulation (LNA) were calculated using the PDM/LDM and PNC/LNC values.Grain yield was determined in a 2-m 2 area in each plot and adjusted to a moisture content of 14.5%.

Calculation of Relative Accumulated Growing Degree Days (RAGDD)
Meteorological data including temperature were collected using the automated weather station Dynameta-1K (Dynamax Inc., Houston, TX, USA), which was installed at each test site, and recorded using the EM50 data acquisition system (Decagon Devices Inc., Washington, USA) every 5 min.Data were then used to calculate the accumulation growing degree days (AGDD) and RAGDD.In this paper, RAGDD was used as a time variable and was calculated from the AGDD, while AGDD represents the sum of growing degree days (GDD) throughout each experiment [39].AGDD can be predicted using GDD or time [40], and in this study, it was calculated as follows: where T Max and T Min represent the maximum and minimum temperatures on a specific day, respectively, and TBASE is the base temperature, which is usually set at 12.5 • C for Japonica rice [41].RAGDD was then calculated as follows: where AGDDi represents the AGDD on a specific day, and AGDD harvest is the AGDD at harvest, representing the AGDD of the entire growth period.

Statistical Analyses
Mapping, profiling, curve fitting, and model building were carried out using OriginPro 9.0 (OriginLab Corporation, Northampton, MA, USA).Four commonly used spectral indexes were determined in this study: NDVI, NDRE, the ratio vegetation index (RVI), and red-edge ratio vegetation index (RERVI).To do so, three wavebands tests were carried out using RapidSCAN during N nutrition monitoring.The specific equations are shown in Table 3. Data-fitting processes were performed using Origin 9.0, choosing different equations based on convergence.Linear, quadratic, logarithmic, exponential, and rational models were evaluated; then, the model with the highest coefficient of determination (R 2 ) was adopted.The accuracy of the model was evaluated using the R 2 and relative root-mean-square error (RRMSE), and the 1:1 relationship between the experimental observations and model predictions were plotted as follows: where n is the number of test samples, P i is the model estimate, O i is the observed value, and O i is the average observed value.

N IR Re
LAI, Biomass, N status [37] Note: NIR, Re, and R refer to the reflectance of near-infrared, red-edge, and red wavelengths, respectively.LAI: leaf area index, HI: harvest index (unit-less), RUE: radiation use efficiency.

Dynamic Changes in Agronomic Parameters and Spectral Indexes before Flowering
The dynamic changes in rice agronomic parameters (LAI, PDM, and PNA) in relation to N nutrition were examined in Exp.1-4 (Figure 2).Similar dynamic characteristics were observed between varieties under different N rates.The results showed that LAI first increased then remained stable (Figure 2A-D).During booting, LAI continued to increase, with the biggest growth rate appearing at the panicle initiation stage.After entering the booting stage, LAI stabilized and remained essentially unchanged thereafter.PDM increased at an invariable rate from tillering to flowering (Figure 2E-H

Relationship between the Spectral Indexes and Agronomic Parameters
Using the results of Exp.1-4 (Table 4), the relationships between each vegetation index and agronomic parameter were also determined, revealing differing results.The relationships between the spectral indexes and LAI and leaf N accumulation were strongest, while those with N concentration were lowest.Moreover, all correlations were low at the tillering stage, increasing with growth until reaching a maximum at stage.During the entire pre-flowering period, the R 2 of the rededge band-based indexes were higher than bidirectional reflectance simulation of the red-band-based vegetation indexes.Using the NDRE and NDVI findings in Exp.1-4, NDVI (Figure 4A) saturation and an intermittent phenomenon were obvious, while the NDRE (Figure 4B) data were more clearly distributed and had a better relationship with LAI.This is also consistent with results in other crops As shown in Figure 3A-H, NDVI and NDRE increased slowly before stem elongation followed by a rapid rise till booting then a slow increase or static period.This was similar to the changes in LAI.Changes in the RVI and RERVI were similar to those of PDM, PNA, and LNA (Figure 3I-P); however, growth rates differed.The growth rate of PDM and PNA accelerated at the end of the stem elongation stage, while RVI remained relatively stable.PNC showed an opposite trend to all spectral indexes and consistency was also poor.All vegetation indexes showed different values under different N gradients.The red-edge-based vegetation indexes NDRE and RERVI showed higher distinction under different N fertilizer rates, while NDVI showed slight overlap between fertilizer rates.Compared with NDRE, RERVI showed a smaller range before flowering, making it hard to distinguish between growth stages and, therefore, less useful.In each single test, NDRE errors and fluctuations were smaller than those of NDVI, making it more stable.Moreover, ranges of all agronomic parameters were also small, suggesting that the environmental effect on NDVI was much bigger than the effect on NDRE.We also compared the relationships between agronomic parameters and final yield during the main growth period (Table 5).As pointed out previously, the relationship between agronomic parameters and yield differs at different growth stages [44,45].This study yielded similar results.The relationship between LAI or N accumulation and yield was best at certain growth stages.Moreover, LAI had the highest coefficient and was more stable than N accumulation, highlighting its potential use in yield predictions.
As mentioned earlier, there were differences in the ranges of NDRE and RERVI before flowering.Stability and distinction between the N gradients of NDRE were better than those of RERVI (Figure 3).Moreover, the dynamic changes in LAI and NDRE were the most consistent (Figures 2 and 3).These results suggest that NDRE can be used to effectively reflect the growth status of rice under different N treatments, and in different eco-sites and different varieties, highlighting its potential use in growth diagnosis and dynamic model construction.

Construction of a Dynamic Model of NDRE
In this study, NDRE values were stable at first before showing a period of rapid growth throughout the entire pre-flowering stage, followed by a second stable, consistent with sigmoidcurve variation (Figure 5).Therefore, we chose a sigmoid curve to construct the dynamic change model of NDRE (Equation 4). −

Relationship between the Spectral Indexes and Agronomic Parameters
Using the results of Exp.1-4 (Table 4), the relationships between each vegetation index and agronomic parameter were also determined, revealing differing results.The relationships between the spectral indexes and LAI and leaf N accumulation were strongest, while those with N concentration were lowest.Moreover, all correlations were low at the tillering stage, increasing with growth until reaching a maximum at stage.During the entire pre-flowering period, the R 2 of the red-edge band-based indexes were higher than bidirectional reflectance simulation of the red-band-based vegetation indexes.Using the NDRE and NDVI findings in Exp.1-4, NDVI (Figure 4A) saturation and an intermittent phenomenon were obvious, while the NDRE (Figure 4B) data were more clearly distributed and had a better relationship with LAI.This is also consistent with results in other crops [43].The relationships between red-edge band-based vegetation indexes and agronomic parameters were, therefore, better than those indexes based on the red band.platform value corresponding to processing of the saturation value (SV), and  is the midpoint of the sigmoid curve, that is, the point at which the model reaches inflection, indicating the different stages of dynamic change (Figure 5).Before  , rice is considered to be in the rapid growth stage, while, after  , the growth rate decreases and the population reaches a peak.Moreover, when  =  ,  = , the growth rate of the model is at its highest.Finally,  represents the time constant of the model.At early and mid-tillering stages, resource requirements including fertilizer and water were relatively small, resulting in a small growing base.Moreover, vegetation coverage was also small, causing the canopy spectra to contain a lot of water and soil background information.As a result, the NDRE values were low and fluctuated around BV.In the mid-stage of growth, particularly early stem elongation, the growth rate increased rapidly, leading to a rapid growth period.After this, in the late stem elongation and booting stages, growth gradually reached a peak, and spectral index growth slowed until flowering.
Results of low (N2T1, 60 + 40 kg•ha −1 ), medium (N5T3, 120 + 80 kg•ha −1 ), and high (N8T7, 240 + 160 kg•ha −1 ) N fertilizer averages are shown in Figure 6.At early growth stages, the differences in NDRE were not significant; however, SV values differed during mid-growth stages.The SV of N2T1 was low at approximately 0.23, while the values of N5T3 and N8T7 were 0.29 and 0.38.Meanwhile, NDVI arrived at its threshold earlier than NDRE as shown in Figure 6, and it was hard to distinguish between N2T1, N5T3, and N8T7 treatments.
The model parameters were mainly obtained by fitting the measured values obtained at a high frequency testing, that is, under Exp. 1 and 3.In addition, BV and SV parameters were also estimated from NDRE values obtained at corresponding growth stages.For example, in the stem elongation to booting stages, NDRE entered the mid-stage of growth in the SV period.At this time, regardless of measurement errors, the measured NDRE values were approximately equal to the SV value.The SVs of Exp. 2 and Exp. 4 were subsequently obtained using this method.At the latest, SV values are determined at the early date.We also compared the relationships between agronomic parameters and final yield during the main growth period (Table 5).As pointed out previously, the relationship between agronomic parameters and yield differs at different growth stages [44,45].This study yielded similar results.The relationship between LAI or N accumulation and yield was best at certain growth stages.Moreover, LAI had the highest coefficient and was more stable than N accumulation, highlighting its potential use in yield predictions.As mentioned earlier, there were differences in the ranges of NDRE and RERVI before flowering.Stability and distinction between the N gradients of NDRE were better than those of RERVI (Figure 3).Moreover, the dynamic changes in LAI and NDRE were the most consistent (Figures 2 and 3).These results suggest that NDRE can be used to effectively reflect the growth status of rice under different N treatments, and in different eco-sites and different varieties, highlighting its potential use in growth diagnosis and dynamic model construction.

Construction of a Dynamic Model of NDRE
In this study, NDRE values were stable at first before showing a period of rapid growth throughout the entire pre-flowering stage, followed by a second stable, consistent with sigmoid-curve variation (Figure 5).Therefore, we chose a sigmoid curve to construct the dynamic change model of NDRE (Equation ( 4)).
where x represents the time-variable RAGDD, y is the spectral index NDRE value, A 1 is the previous platform value corresponding to the baseline value (BV) in each treatment, A 2 is the medium platform value corresponding to processing of the saturation value (SV), and x 0 is the midpoint of the sigmoid curve, that is, the point at which the model reaches inflection, indicating the different stages of dynamic change (Figure 5).Before x 0 , rice is considered to be in the rapid growth stage, while, after x 0 , the growth rate decreases and the population reaches a peak.Moreover, when x = x 0 , x = A 1 +A 2 2 , the growth rate of the model is at its highest.Finally, dx represents the time constant of the model.At early and mid-tillering stages, resource requirements including fertilizer and water were relatively small, resulting in a small growing base.Moreover, vegetation coverage was also small, causing the canopy spectra to contain a lot of water and soil background information.As a result, the NDRE values were low and fluctuated around BV.In the mid-stage of growth, particularly early stem elongation, the growth rate increased rapidly, leading to a rapid growth period.After this, in the late stem elongation and booting stages, growth gradually reached a peak, and spectral index growth slowed until flowering.
Results of low (N2T1, 60 + 40 kg•ha −1 ), medium (N5T3, 120 + 80 kg•ha −1 ), and high (N8T7, 240 + 160 kg•ha −1 ) N fertilizer averages are shown in Figure 6.At early growth stages, the differences in NDRE were not significant; however, SV values differed during mid-growth stages.The SV of N2T1 was low at approximately 0.23, while the values of N5T3 and N8T7 were 0.29 and 0.38.Meanwhile, NDVI arrived at its threshold earlier than NDRE as shown in Figure 6, and it was hard to distinguish between N2T1, N5T3, and N8T7 treatments.The model parameters were mainly obtained by fitting the measured values obtained at a high frequency testing, that is, under Exp. 1 and 3.In addition, BV and SV parameters were also estimated from NDRE values obtained at corresponding growth stages.For example, in the stem elongation to booting stages, NDRE entered the mid-stage of growth in the SV period.At this time, regardless of measurement errors, the measured NDRE values were approximately equal to the SV value.The SVs of Exp. 2 and Exp. 4 were subsequently obtained using this method.At the latest, SV values are determined at the early date.To confirm the accuracy of the dynamic models, they were fitted using each experiment and N treatment (Table 6).The results of model fitting based on NDRE showed that the R 2 values of different N gradients were greater than 0.9, while the values of RRMSE were less than 5%, suggesting that model values showed very small differences compared to the measured values, also indicating that the modeling effect of the different gradients was good.The normal distribution test of standard residuals under each treatment model showed that the probability of normal distribution was greater than 0.05, revealing that the standard residuals of each treatment were within normal distribution and that the model errors were random and generated by the model itself.Model fitting based on NDVI resulted in a high R 2 (0.87-0.96); however, the range of RRMSE was 7.5-21% higher than that of NDRE.Overall, these findings suggest that the model was applicable under different N levels, with different varieties, in different eco-sites, and in different years, confirming the accuracy of the model.

Changes in BV and SV under Different Treatments and Their Relationship with Agronomic Parameters
The NDRE values changed little at early stages of growth, fluctuating around the BV value.The BV under each treatment was less than 0.25 and there were no differences between treatments (Figure 7, NDRE), while BV and SV estimated using NDVI showed no obvious regularity.SV calculated using NDRE increased with increasing N fertilizer rates, and the trends between treatments were consistent.The differences with topdressing under the same basal fertilizer also differed, and the SV of T0 decreased by 15-25%.Moreover, a 10% increase in SV was observed between N0T0 and N5T0, N5T0 and N7T0, and N7T0 and N8T0.The differences between N application and topdressing levels also differed, suggesting that SV has obvious potential in growth status analysis in rice.
The relationship between SV and agronomic parameters was subsequently examined using the data from Exp. 1 and Exp. 3 (Figure 8).The relationships between SV and LAI, PDM, PNA, and LNA were good at both stem elongation (fourth leaf from the top stage) and booting stage (second leaf from the top stage).For NDRE, the R 2 of each of these relationships was more than 0.6.At stem elongation, there were slight differences but essentially similar correlations between SV and LAI, PDM, PNA, and LNA.Correlations at booting stage were better than those at stem elongation.Overall, LAI had the best relationship with SV.However, with NDVI, poor correlation was observed between all parameters.

Prediction Algorithm for Rice Grain Yield
The relationship between BV, SV, and x 0 with grain yield was also studied (Figure 9).The results showed a poor relationship between BV and yield (Figure 9A,E), while x 0 value was negatively correlated with yield, and, although the correlation was poor, differences were observed among cultivars (Figure 9B,F).There were no differences in SV between experiments, suggesting that SV showed the similar correlations with yield between different years, eco-sites, and varieties (Figure 9C).Meanwhile, the NDVI-based SV was poorly correlated with grain yield.The relationship between BV, SV, and  with grain yield was also studied (Figure 9).The results showed a poor relationship between BV and yield (Figures 9A,E), while  value was negatively correlated with yield, and, although the correlation was poor, differences were observed among cultivars (Figures 9B,F).There were no differences in SV between experiments, suggesting that SV showed the similar correlations with yield between different years, eco-sites, and varieties (Figure 9C).Meanwhile, the NDVI-based SV was poorly correlated with grain yield.
The quantitative relationship between SV and grain yield was developed (Figure 10) using the data from Exp. 1 and Exp. 3. As shown in Figure 10, grain yield was forecast using linear regression analysis based on the SV values.Grain yield was evenly distributed from 6000 to 13000 kg•ha −1 and the model was good.The slope of the linear model was 36081 kg•ha −1 , suggesting that with every 0.1 increase in SV, yield increased by 3608.1 kg•ha −1 .The quantitative relationship between SV and grain yield was developed (Figure 10) using the data from Exp. 1 and Exp. 3. As shown in Figure 10, grain yield was forecast using linear regression analysis based on the SV values.Grain yield was evenly distributed from 6000 to 13000 kg•ha −1 and the model was good.The slope of the linear model was 36081 kg•ha −1 , suggesting that with every 0.1 increase in SV, yield increased by 3608.1 kg•ha −1 .

Prediction Algorithm for Rice Grain Yield
The relationship between BV, SV, and  with grain yield was also studied (Figure 9).The results showed a poor relationship between BV and yield (Figures 9A,E), while  value was negatively correlated with yield, and, although the correlation was poor, differences were observed among cultivars (Figures 9B,F).There were no differences in SV between experiments, suggesting that SV showed the similar correlations with yield between different years, eco-sites, and varieties (Figure 9C).Meanwhile, the NDVI-based SV was poorly correlated with grain yield.
The quantitative relationship between SV and grain yield was developed (Figure 10) using the data from Exp. 1 and Exp. 3. As shown in Figure 10, grain yield was forecast using linear regression analysis based on the SV values.Grain yield was evenly distributed from 6000 to 13000 kg•ha −1 and the model was good.The slope of the linear model was 36081 kg•ha −1 , suggesting that with every 0.1 increase in SV, yield increased by 3608.1 kg•ha −1 .The data from Exp. 2 and Exp. 4 were subsequently used to test the accuracy of the SV-based yield prediction model (Figure 11).As shown, in Figure 11A, SV was obtained from the fitted sigmoid curve model, while, in Figure 11B, values were estimated using measured NDRE values from the stem elongation to booting stages.Both methods gave good prediction results.The predicted yields in both study years showed no obvious abnormalities compared with the observed yield, and the the error ranges of different gradients in spectral indexes containing red bands was observed after stem elongation, preventing clarification of the differences between N gradients.Studies showed that the absorption of chlorophyll on the red-edge waveband is weaker than that of red light, with the red-edge region having stronger transmission ability with the crop canopy and leaf [25].Use of red-edge rather than red-light bands can, therefore, reduce the saturation phenomenon, enhancing monitoring and diagnosis of crop N nutrition.
In this study, NDVI did not distinguish between different N treatments during mid-growth and late growth stages, similar to the findings of Liu et al. [41].It was also previously revealed that, at an early growth stage (before stem elongation), NDVI resulted in a better relationship with rice grain yield; however, on approaching maturity, NDRE performed better [50].Red-edge-based NDRE reflects this phenomenon via differences in saturation values.Our results suggest that the shortcoming of NDVI is that saturation at an early stage leads to no differences in SV with varying N rates.Yumiko et al. also reported that NDVI reached saturation earlier than NDRE under different N [50].NDVI can, therefore, be used to determine low-land rice at early growth stages, but is less useful at mid-growth to late growth stages.Moreover, using vegetation indexes with the red-light band, correlations with agronomic parameters were poor, while the red-edge band-based indexes (NDRE and RERVI) showed a linear relationship with agronomic parameters, allowing differences between N fertilizer levels to be distinguished.However, changes in RERVI throughout the entire growth period were small, thereby preventing effective analysis between different stages.In the same plot, NDRE showed lower variability than RERVI.These findings suggest that NDRE can be used to reflect the growth status of different treatments, making it a suitable index for diagnosis of rice growth.

The Dynamic Model and Parameters of Rice Spectral Indexes
The results of this study suggest that NDRE conforms to a sigmoid curve.In early stages, NDRE fluctuated around background BV values at a low level.However, values began to rise rapidly after tillering.The critical point (x 0 ) of the dynamic model was then reached on entering stem elongation.Meanwhile, at flowering stage, the dynamic model became stable again and approached the feature platform SV.
Previous monitoring methods such as the saturation index (SI) method and NNI are based mainly on data collected at each growth stage [8].These methods mostly lack biophysical mechanisms, often fail to obtain indicators, and consist of a complex calculation method [25].Moreover, deviations in identification of growth stages, environmental factors, time selection, and so on have a significant impact on these single tests.Furthermore, single measured values are treated as a reflection of an entire stage in the calculations, resulting in greater deviation from the actual values [51].
In this paper, the parameters in the model have practical biological meaning.The range of x 0 reflects growth and development of rice at early stages.The findings also suggest that higher x 0 values occur under abnormally high temperatures in early stages compared with normal growing conditions.The x 0 value can be used as an indicator of the growth environment during early stages of growth; however, the poor correlation between x 0 and grain yield prevented quantification of yield predictions.
BV was also able to reflect growth status at early stages.However, when rice coverage was less than 50%, including the depth of the water layer, turbidity, soil background color, and so on, there was an impact on the spectral reflectance test [52,53].In this study, BV values were low under all treatments, and no obvious differences between gradients were observed.Meanwhile, BV did not distinguish between differences in basal fertilizers levels in early stages of growth, suggesting that water and soil background information under low coverage conditions during this period obscure growth information contained in the spectral data.BV is, therefore, not applicable for direct monitoring and diagnosis of rice growth, although differences in LAI and other agronomic parameters were obvious.In the future, we aim to increase accuracy by adding variables such as coverage, soil, and water to improve the signal-to-noise ratio.
In the middle of the growing season, SV values overcame the above problems in x 0 and BV.SV values differed between different rice populations under different treatments and in different groups, and there was a good linear relationship with agronomic parameters such as LAI and aboveground N accumulation.SV values could be obtained not only from the fitted model, but also from measured values.Overall, the accuracy and availability of the model were high, confirming applicability for monitoring and diagnosis.
The dynamic model was successfully used to determine growth status using different rice varieties, in different eco-sites, and under different N rates before flowering with 90% variation.Bonfil et al. also reported similar dynamic changes in wheat in 2016 [52].Moreover, our dynamic model also determined the growth status of paddy rice using RapidScan CS-45 at different grain yield levels.Of course, it can also be used to predict yield before flowering without any specific calibration.It should also be used on non-rainy days, since water droplets on the blade affect monitoring.In the future, the growth status after flowering stage will be examined further, and additional vegetation indexes will be used in developing a dynamic model.

SV-Based Yield Prediction Algorithm
In previous studies, spectral indexes were used to predict yield in crops such as wheat [54] and sugarcane [55] with good results.Bonfil used RapidScan CS-45 to monitor the dynamic changes in wheat, revealing a good linear relationship between NDRE and wheat yield in different varieties [54].In this study, the dynamic pattern of NDRE at pre-flowering stages was consistent, and predicted yield was achieved under optimal growing conditions other than limited N fertilizer.Kanke et al. also revealed that vegetation indexes containing the red-edge band had a better linear relationship with agronomic parameters and yield at main growth stages in rice [55].A linear model of yield prediction using different vegetation indexes as input variables was also established, within which NDRE and RERVI gave best correlations, with R 2 values greater than 0.8.However, the prediction equation was predicted using data from a single growth stage, and the parameters differed between stages, making it of no use between stages or deviations [56].Similarly, Cowley et al. established an NDVI-based model for forecasting of rapeseed yield [57], while Liu et al. proved the dynamic model could be continuously used to monitor N nutrition status in a single growing season [41].However, it could only be used when GDD was between 210 and 320 • C.
Thus, the problem of low prediction accuracy resulting from NDVI saturation when using traditional algorithms in rice under high yield conditions remains unsolved.In this study, SV showed a good linear relationship with agronomic parameters and yield.Moreover, SV values were the feature platform of each treatment, eliminating deviations between tests and different environments, thereby increasing applicability of the model, improving stability of the algorithm, and widening the available time range.In addition, SV values could be obtained not only from the fitted model, but also from measured values.Moreover, measurement and processing are relatively simple and do not require normalization or other transformations to maintain accuracy, simplifying our method compared with other algorithms [58].
Satellite remote sensing (RS) data were extensively applied to a wide range of research problems and practical applications, including yield predictions [59].Zhu et al. used Moderate Resolution Imaging Spectroradiometer Enhanced Vegetation Index (MODIS EVI) time-series data to predict yield in winter wheat, revealing a high R 2 (0.70) and low RMSE (343.34 kg•hm −2 ) [60].Furthermore, Fletcher reported continuity of multispectral high-resolution optical observations over global terrestrial surfaces using Sentinel 2 [61], providing detailed information within field variability [62].Our newly developed SV-based yield prediction algorithm must now be tested and applied under different RS platforms, such as unmanned aerial vehicles and Sentinel 2.

Interval System Errors between Experiments and Uncertainties of the Model
In different experiments, trends in indicators were basically the same, although the ranges of change were quite different, especially in Exp. 1 (2015) and Exp. 3 (2016).This was possibly due to the differences in meteorological conditions in 2015 and 2016 (Figure 13).From June to early July 2016, the temperature was slightly higher than in 2015; thus, growth was relatively fast in early and mid-tillering.However, in late July to August, unusually high temperatures and continuous rain occurred, which is not conducive to rice growth and development, thereby inhibiting growth.As a result, the agronomic parameters, spectral indexes, and other indicators decreased sharply in 2015.

Interval System Errors between Experiments and Uncertainties of the Model
In different experiments, trends in indicators were basically the same, although the ranges of change were quite different, especially in Exp. 1 (2015) and Exp. 3 (2016).This was possibly due to the differences in meteorological conditions in 2015 and 2016 (Figure 13).From June to early July 2016, the temperature was slightly higher than in 2015; thus, growth was relatively fast in early and mid-tillering.However, in late July to August, unusually high temperatures and continuous rain occurred, which is not conducive to rice growth and development, thereby inhibiting growth.As a result, the agronomic parameters, spectral indexes, and other indicators decreased sharply in 2015.
This phenomenon was also reflected in another way.In both years,  values remained relatively concentrated at the end of the tillering stage; however, there were differences between N rates.Moreover, the critical point shifted in 2016 compared to 2015.The average daily temperature was higher than 33 °C at the end of tillering and early stem elongation stage in Exp. 3, 2016.The optimum temperature of paddy rice at tillering is 25~30 °C, while the maximum temperature is 32 °C [42] (Figure 13).As a result, these high temperatures caused an abnormal increase in the RAGDD time variable, even though growth was inhibited and  increased.This phenomenon was also reflected in another way.In both years, x 0 values remained relatively concentrated at the end of the tillering stage; however, there were differences between N rates.Moreover, the critical point shifted in 2016 compared to 2015.The average daily temperature was higher than 33 • C at the end of tillering and early stem elongation stage in Exp. 3, 2016.The optimum temperature of paddy rice at tillering is 25~30 • C, while the maximum temperature is 32 • C [42] (Figure 13).As a result, these high temperatures caused an abnormal increase in the RAGDD time variable, even though growth was inhibited and x 0 increased.
Model uncertainty received considerable attention in recent studies [63].There is growing pressure on crop field production to conform to quality standards, which require evaluation and expression of the uncertainty of measurement results [64].It is, therefore, important to determine the degree of uncertainty associated with a model and its application in predicting crop N status and yield.In this study, the dynamic model was developed using a large database, diverse eco-sites, different rice varieties, and differing N rates; however, errors were revealed when the model output was compared with observations.Moreover, it is difficult to evaluate uncertainty in practice.For example, the uncertainty of certain components is often unknown in crop field experiments, such as those associated with environmental effects or different sensors.During the 2015-2017 growing season, the sensors were calibrated each year to avoid errors; however, additional field experiments are needed to expand the database and therefore improve the robustness of our dynamic model.

Conclusions
We used the active spectrometer RapidScan CS-45 to determine the dynamic changes in different spectral indexes and agronomic parameters in rice.The results suggest that NDRE is more stable than other indexes, with the distinctions between different treatment gradients most obvious.NDRE also had the best relationship with agronomic parameters.The dynamic changes in NDRE and LAI were similar, and followed a "slow-fast-slow" trend during the pre-flowering period.
A sigmoid curve was used in NDRE dynamic model construction, with high accuracy and an R 2 of 0.9 and RRMSE of 5%.The model contained critical point x 0 , background value BV, and feature platform value SV parameters.The relationships with yield showed that each of these parameters was correlated with yield, although SV was better than the other two parameters.SV showed a good linear relationship with in-season potential grain.Based on these findings, we constructed a yield prediction model.With every 0.1 increase in SV, grain yield increased by 3608.1 kg•ha −1 .The overall prediction accuracy of the model was greater than 80%.It can be used at stem elongation and booting, improving timeliness compared with existing forecasting algorithms.Our model is, therefore, useful for rice production management, purchasing and storage preparation, rice spot price trend forecasts, and food policy decision-making.

1 .
Experiment 4 (Exp.4): Experiment 4 was a repeat of Exp. 2, carried out in 2016.Detailed information on each experiment is provided in Table 2. Data from Exp. 1 and Exp. 3 were mainly used to develop the model, while data from Exp. 2 and Exp. 4 were used for model testing.Exp. 2 and 4 were performed to maintain uniformity across different study years.Remote Sens. 2019, 11, x FOR PEER REVIEW 5 of 25

1 .
Experiment 4 (Exp.4): Experiment 4 was a repeat of Exp. 2, carried out in 2016.Detailed information on each experiment is provided in Table 2. Data from Exp. 1 and Exp. 3 were mainly used to develop the model, while data from Exp. 2 and Exp. 4 were used for model testing.Exp. 2 and 4 were performed to maintain uniformity across different study years.

Figure 1 .
Figure 1.Daily temperatures during the rice growing seasons in 2015-2016.

Figure 1 .
Figure 1.Daily temperatures during the rice growing seasons in 2015-2016.
), while PNA was consistent with LNA, continually increasing until flowering (Figure 2M-T).PNC differed slightly from the other N indicators, decreasing up until flowering at an increasing rate (Figure 2I-L).Remote Sens. 2019, 11, x FOR PEER REVIEW 8 of 25

Figure 4 .
Figure 4. Relationship between (A) the normalized difference vegetation index (NDVI) and (B) normalized difference red edge (NDRE) and leaf area index (LAI) before flowering.

Figure 4 .
Figure 4. Relationship between (A) the normalized difference vegetation index (NDVI) and (B) normalized difference red edge (NDRE) and leaf area index (LAI) before flowering.

25 Figure 5 .
Figure 5. Dynamic characteristics of the sigmoid curve.NDRE: normalized difference red edge, RAGDD: relative accumulated growing degree day, A1: the platform value corresponding to the baseline value (BV) in each treatment, A2: the medium platform value corresponding to processing of the saturation value (SV),  : the midpoint of dynamic change, that is, the point at which the model reaches its inflection point, indicating the different stages of dynamic change.

Figure 5 .
Figure 5. Dynamic characteristics of the sigmoid curve.NDRE: normalized difference red edge, RAGDD: relative accumulated growing degree day, A1: the platform value corresponding to the baseline value (BV) in each treatment, A2: the medium platform value corresponding to processing of the saturation value (SV), x 0 : the midpoint of dynamic change, that is, the point at which the model reaches its inflection point, indicating the different stages of dynamic change.

25 Figure 5 .
Figure 5. Dynamic characteristics of the sigmoid curve.NDRE: normalized difference red edge, RAGDD: relative accumulated growing degree day, A1: the platform value corresponding to the baseline value (BV) in each treatment, A2: the medium platform value corresponding to processing of the saturation value (SV),  : the midpoint of dynamic change, that is, the point at which the model reaches its inflection point, indicating the different stages of dynamic change.

Figure 8 .
Figure 8. Relationships between saturation values (SV, calculated using NDRE or NDVI) and leaf area index (LAI), plant dry matter, plant nitrogen accumulation, and leaf nitrogen accumulation at stem elongation and booting.

Figure 8 .
Figure 8. Relationships between saturation values (SV, calculated using NDRE or NDVI) and leaf area index (LAI), plant dry matter, plant nitrogen accumulation, and leaf nitrogen accumulation at stem elongation and booting.

Figure 9 .
Figure 9. Relationships between each vegetation index (NDRE, NDVI) based on different model parameters (SV, BV, and x0) and grain yield.BV: baseline value, SV: saturation value, and x0: the midpoint of the sigmoid curve.(A,E) show the BV, (B,F) show the x0, and (C,G) show the SV.

Figure 10 .
Figure 10.Yield prediction linear model based on saturation values (SV).SV in (A) was calculated using NDRE, and in (B) was calculated by NDVI.

Figure 9 .
Figure 9. Relationships between each vegetation index (NDRE, NDVI) based on different model parameters (SV, BV, and x 0 ) and grain yield.BV: baseline value, SV: saturation value, and x 0 : the midpoint of the sigmoid curve.(A,E) show the BV, (B,F) show the x 0 , and (C,G) show the SV.

Figure 9 .
Figure 9. Relationships between each vegetation index (NDRE, NDVI) based on different model parameters (SV, BV, and x0) and grain yield.BV: baseline value, SV: saturation value, and x0: the midpoint of the sigmoid curve.(A,E) show the BV, (B,F) show the x0, and (C,G) show the SV.

Figure 10 .
Figure 10.Yield prediction linear model based on saturation values (SV).SV in (A) was calculated using NDRE, and in (B) was calculated by NDVI.

Figure 10 .
Figure 10.Yield prediction linear model based on saturation values (SV).SV in (A) was calculated using NDRE, and in (B) was calculated by NDVI.

Figure 13 .
Figure 13.Daily average temperatures in the study site from late tillering to early stem elongation in 2015 and 2016.

Table 3 .
Equations of each vegetation index.

Table 4 .
Coefficients of determination (R 2 ) of the correlations between vegetation indexes at different stages and the agronomic parameters.
Notes: LAI: leaf area index, PNC: plant nitrogen concentration, PNA: plant nitrogen accumulation, LNA: leaf nitrogen accumulation, NDVI: normalized difference vegetation index, NDRE: normalized difference red edge, RERVI: red-edge ratio vegetation index, RVI: ratio vegetation index.The bold word means higher R 2 value among four experiments.

Table 5 .
Coefficients of determination (R 2 ) of the correlations between the agronomic parameters at different stages and yield.

Table 6 .
Summary of fitting analysis of the NDRE dynamic model under each treatment.

Table 7 .
Distribution of the critical point of dynamic change under a single time series.