Estimation of Maize Yield and Protein Content under Different Density and N Rate Conditions Based on UAV Multi-Spectral Images

: In the ﬁeld of precision agriculture research, it is very important to monitor crop growth in time so as to effectively conduct ﬁeld diagnosis and management and accurately predict yield and quality. In this experiment, the relationship between the vegetation index of Zhengdan 958 and Suyu 41 and their yield and quality when reducing N application (25 and 50% N reduction compared to local conventional N application rate) under low, medium and high planting densities (60,000, 75,000 and 90,000 plants · ha − 1 ) during 2018–2020 was investigated using multispectral images obtained from UAV monitoring. The results showed that under different density treatments, the normalized vegetation index (NDVI) and ratio vegetation index (RVI) decreased with the decrease in nitrogen application, while the plant senescence reﬂectance index (PSRI) increased. Through principal component analysis (PCA) and subordinate function analysis, the comprehensive score of each treatment can reﬂect the maize yield and total protein content under each treatment. Based on the vegetation index, predictive models of maize yield and protein content were established. The best prediction period for grain yield and protein content were physiological maturity and 35 days after silking (R4), respectively. The R 2 of the predictive models are greater than 0.734 and 0.769, respectively. Multi-period and multi-vegetation indexes can better monitor crop growth and help agricultural ﬁeld management.


Introduction
Maize is a multifunctional crop integrating grain and forage crops.In China, approximately 70% of the maize is used as feed [1].However, with the continuous improvement in maize yield, the protein content of hybrids decreased gradually, which greatly reduced the feed efficiency [2].It is essential to improve the grain yield and protein quality simultaneously to balance the demand for maize production.The main way to improve yield is to select density tolerant varieties and adopt high-density cultivation [3][4][5].In addition, reasonable nitrogen (N) application is also important to the yield and quality of maize.Nevertheless, farmers used to apply excessive N fertilizer to increase the yield, which was not conducive to yield and quality improvement and led to resource waste [6].Therefore, it is imperative to determine the suitable N and density combination to improve yield and maintain protein quality.
As a short-term crop, maize has a tall plant type, difficult sampling and poor representation in the field.When the micro-environment of crops alters, the vegetation coverage, biomass accumulation and translocation changes and the plant canopy reflection spectrum is also adjusted [7].Existing research shows that the canopy reflectance spectra of plant leaves with different biochemical compositions is also different [8].Different physical and chemical properties of crops will show different absorption of frequency doubling and Agronomy 2023, 13, 421 2 of 16 frequency combination at the spectral level; thus we can use visible light and near-infrared spectroscopy to analyze the substance qualitatively or quantitatively [9].Based on Unmanned Aerial Vehicle (UAV) multispectral images, the N status of crops can be monitored quickly and non-destructively, which is conducive to the fine management of crops [10].High-throughput phenotypic analysis of crops can also be carried out according to the obtained plant height, canopy coverage and vegetation index, so as to indirectly estimate the physiological and environmental stress status of crops [11], which is of great significance for realizing agricultural precision, informatization and intelligence.Multispectral imaging technology is a new, rapid, nondestructive testing technology developed in recent years.It combines spectral information with imaging data, uses chemometrics technology to screen and model the data, and then constructs a fast, nondestructive, efficient and accurate evaluation method of the grain yield [12] and quality [13].It can timely and effectively guide maize production and promote the development of maize production to precision.However, the accuracy of the predictive model established using single-band reflectance to monitor crop growth is not as good as that of the vegetation index [14] because the vegetation index not only weakens the influence of external environmental factors on spectral signals, but also contains information about canopy structure and material characteristics.The vegetation index has been used to predict nitrogen accumulation [15], crown-layer chlorophyll content [16] and leaf area index [17].Previous studies have shown that NDVI, RVI, DVI and other vegetation indexes have a positive correlation with crop yield [18].There is significant correlation between crop protein and some bands of the canopy spectrum [19].In theory, yield and protein content can be estimated by spectral characteristics of these sensitive bands.Previous researchers have also established yieldprediction models using NDVI [20], RVI [21], DVI [22], and protein-content-prediction models using NDVI [23], RVI [24] and SIPI [25].
The subordinate function analysis is often used to evaluate crop resistance [26] and variety selection [27].However, there are many indicators involved in the evaluation and there is a certain correlation between the indicators, so there are some limitations in the comprehensive evaluation of crops using only the subordinate function method [28].Principal component analysis (PCA) is an analysis method that clusters multivariate data and reduces dimensions.The advantage of this analysis method is that many original indicators can be integrated into several relatively independent comprehensive indicators by taking advantage of the deep relationship between indicators.Therefore, following PCA, subordinate function analysis is more suitable for comprehensive evaluation.However, this assessment method is rarely used to assess crop growth among treatments.Therefore, this study tried to use this evaluation method to evaluate the difference in maize growth among treatments.
At present, most studies on monitoring maize yield and quality using UAV hyperspectral imaging technology are reflected by a single or several vegetation indexes [29].Few studies have combined and synthesized these vegetation indexes as new indicators for monitoring maize yield and quality.ew studies have investigated the monitoring effect of UAV on maize yield and quality, taking into account the impact of cultivation measures.Therefore, we conducted a field experiment in 2019-2020.In this experiment, Zhengdan 958 and Suyu 41 were used as materials, and three planting densities of low, medium and high (60,000, 75,000 and 90,000 plants•ha −1 ) with conventional N application, N reduction by 25% and N reduction by 50% (300, 225 and 150 kg•ha −1 ) were set.The purpose of our study was to (1) determine the response effect of the vegetation index on maize growth, (2) analyze maize growth under different treatments by comprehensive scoring of the vegetation index, (3) explore the predictive effect of multispectral imaging on maize yield and protein content prediction, (4) expand the application scope and improve the precision of UAV monitoring of crop growth and ( 5) create a practical model for predicting the grain yield and protein content of summer maize.

Experimental Site
The field experiment was conducted in the Experimental Farm of Lianyungang Academy of Agricultural Sciences, Jiangsu Province, China (34.57• N, 119.32 • E) from 2019 to 2020.The altitude of this area is 1 m.The average content of organic matter in the tillage layer was 10.2 g•kg −1 and the total nitrogen (N), rapid phosphorus (P), and rapid potassium (K) were 0.84 g•kg −1 , 94.2 mg•kg −1 and 141.5 mg•kg −1 , respectively.The available phosphorus and potassium were determined with molybdenum blue colorimetry and flame photometry, respectively.The type of soil is clay with a field capacity of 28%, PH of 7.3% and bulk density of 1.3 g/cm3.The total rainfall during the crop growth period was 257.45 mm and 437.2 mm in 2019 and 2020, respectively.The average temperature during the growth period in 2019 and 2020 was 25.5 • C and 22.7 • C, respectively.The location diagram of the test area is shown in Figure 1.
yield and protein content prediction, (4) expand the application scope and improve the precision of UAV monitoring of crop growth and ( 5) create a practical model for predicting the grain yield and protein content of summer maize.

Experimental Site
The field experiment was conducted in the Experimental Farm of Lianyungang Academy of Agricultural Sciences, Jiangsu Province, China (34.57°N, 119.32 °E) from 2019 to 2020.The altitude of this area is 1 m.The average content of organic matter in the tillage layer was 10.2 g•kg −1 and the total nitrogen (N), rapid phosphorus (P), and rapid potassium (K) were 0.84 g•kg −1 , 94.2 mg•kg −1 and 141.5 mg•kg −1 , respectively.The available phosphorus and potassium were determined with molybdenum blue colorimetry and flame photometry, respectively.The type of soil is clay with a field capacity of 28%, PH of 7.3% and bulk density of 1.3 g/cm3.The total rainfall during the crop growth period was 257.45 mm and 437.2 mm in 2019 and 2020, respectively.The average temperature during the growth period in 2019 and 2020 was 25.5 °C and 22.7 °C, respectively.The location diagram of the test area is shown in Figure 1.

Crop Management and Experimental Design
The maize varieties Zhengdan 958 and Suyu 41 were selected as experimental materials.They are erectophile maize hybrids that are suitable for dense planting.Their growth period was about 95-110 d in the test area.Their plant height is 230-240 cm, and the ear length is 18-20 cm.Both the two varieties are widely grown in Jiangsu Province.Three planting densities of low (local conventional density), medium and high (60,000, 75,000 and 90,000 plants•ha −1 ) and three nitrogen application levels of local conventional nitrogen (N) application, 25 and 50% N reduction compared with conventional N application (300, 225 and 150 kg•ha −1 ) were set.The experiment was sown at the end of June and harvested at the beginning of October every year.Due to the different weather conditions between years, the specific dates varied.The average moisture content of the grain at harvest was 25-30%.The forecrop was wheat.Before sowing, a compound fertilizer with N-P-K of 18%-18%-18% was applied in each plot.At the V6 and V12, urea was topdressing.The

Crop Management and Experimental Design
The maize varieties Zhengdan 958 and Suyu 41 were selected as experimental materials.They are erectophile maize hybrids that are suitable for dense planting.Their growth period was about 95-110 d in the test area.Their plant height is 230-240 cm, and the ear length is 18-20 cm.Both the two varieties are widely grown in Jiangsu Province.Three planting densities of low (local conventional density), medium and high (60,000, 75,000 and 90,000 plants•ha −1 ) and three nitrogen application levels of local conventional nitrogen (N) application, 25 and 50% N reduction compared with conventional N application (300, 225 and 150 kg•ha −1 ) were set.The experiment was sown at the end of June and harvested at the beginning of October every year.Due to the different weather conditions between years, the specific dates varied.The average moisture content of the grain at harvest was 25-30%.The forecrop was wheat.Before sowing, a compound fertilizer with N-P-K of 18%-18%-18% was applied in each plot.At the V6 and V12, urea was topdressing.The experiments were established in a split-plot design consisting of three replicates, with planting densities as the main plot and N rates as the sub-plot.The size of each subplot was 12.0 m × 5.0 m (rows spaced 60 cm apart) and the interval between each plot was 1 m.During the growth period, we ensured a good water supply, timely management of diseases, pests and weeds and a good growth environment for the plants.

Grain Yield
Three rows of 5 m ears were harvested at maturity.The moisture was measured with a moisture meter and dried naturally.The GY was expressed at a moisture content of 14.0%.

Grain Protein Content
Grain samples were ground and then digested by adding 10 mL H 2 SO 4 -H 2 O 2 to each sample.The N content of grain was measured using the semi-micro Kjedahl method (KjeltecTM 8200 Auto Distillation Unit, Foss, Hillerød, Denmark, Yuen and Pollard 1953; Bremner 1960).The grain protein content is determined by multiplying the nitrogen content by 6.25.

Multi-Spectral Data Acquisition
In this study, the Ecodroneuas-8 multi-spectral UAV (EcoDroneUAS-8) was used with a multispectral camera (rededge MX, Micasense, Raleigh, NCS, USA) for data acquisition (Table 1).The field test data were obtained at 14 days (R2) after silking, 35 days after silking (R4) and physiological maturity (R6).The flight time of the UAV was about 12:00 am, with good lighting conditions, no wind or breeze and no clouds.The UAV flew at 20 m, the heading overlap rate was 75% and the sideway overlap rate was 75%.Pix4dMapper software was used to Mosaic the UAV image, and the UAV digital orthophoto was obtained.The ENVI software was used to extract the data.The reflectivity values of the maize canopy in 5 bands were extracted through the processing and analysis of UAV multi-spectral images and the data were statistically analyzed.According to the band characteristics of multispectral sensors and the results of previous studies, several spectral indexes were calculated (Table 2).

Statistical Analysis
Microsoft Excel 2013 and SigmaPlot 12.5 were used to process the test data and plot, and SPSS 20.0 was used for data statistics and analysis and significance analysis (LSD method).
The vegetation index value was used for principal component analysis (PCA) to evaluate the differences among the 9 treatments, and the final total score was calculated through the subordinate function analysis to represent the advantages and disadvantages of the treatment.Pearson correlation was calculated to determine the relationship between the score and its yield and protein content.
According to the subordinate function formula of fuzzy mathematics, the parameters and indicators in the test process are quantitatively transformed to obtain the subordinate function value, which is used to compare the yield and quality under different treatments.At present, subordinate function analysis is a common method to evaluate abiotic resistance, and it is calculated using the following formula [30,31]: u(Xj) = (Xj-Xmin)/(Xmax-Xmin) j = 1, 2, 3, ..., n where u(Xj) is the subordinate function value of the j-th comprehensive indicator, Xj is the j-th comprehensive indicator value, Xmax is the maximum value of the j-th comprehensive indicator, and Xmin is the minimum value of the j-th comprehensive indicator.
where Wj is the index weight of the j-th comprehensive indicator among all comprehensive indicators and Pj is the contribution rate of the j-th comprehensive indicator.
where D is the comprehensive evaluation value of the different treatments.The higher the D value, the better the treatment; the lower the D value, the worse the treatment.

Change in Vegetation Index
From 2019 to 2020, variety, density and N rate all had significant effects on NDVI, RVI and PSRI at different growth stages (Supplementary Tables S1 and S2).At the same time, NDVI, RVI and PSRI also showed a high correlation with yield and total protein content (Figures 2 and 3).Under different densities, the NDVI and RVI of the two varieties decreased with the decrease in N rate.NDVI and RVI increased with the increase in density under the same N rate.The trend of PSRI with nitrogen dose and density was opposite to that of NDVI and DVI (Figures 4-6).Take the example of the change in NDVI, RVI and PSRI of Zhengdan 958 in the R2 period in 2019.At a low density, the NDVI and RVI of Zhengdan 958 decreased by 1.3% and 1.7%, respectively, and the PSRI increased by 3.0%, when the N rate was decreased by 25%.The NDVI and RVI decreased by 6.8% and 8.7%, respectively, and the PSRI increased 4.3%, when the N rate was decreased by 50%.At a medium density, the NDVI and RVI decreased by 0.33% and 0.62%, respectively, and the PSRI increased by 1.3%, when the N application was decreased by 25%.The NDVI and RVI decreased by 3.1% and 1.5%, respectively, and the PSRI increased by 5.5%, when the N rate was reduced by 50%.At a high density, when the N rate was decreased by 25%, the NDVI and RVI decreased by 1.6% and 0.46%, respectively, and the PSRI increased by 0.58%.When the N rate was decreased by 50%, the NDVI and RVI decreased by 2.6% and 3.8%, respectively, and the PSRI increased by 3.8%.

PCA of Vegetation Index of Different Treatments of Reduced Nitrogen Application at Different Densities
In order to eliminate the factors with less influence and greater interference and improve the accuracy of measurement data analysis, the single index introduced above was converted into a smaller number of more effective indicators [32].Therefore, in the experiments in 2019 and 2020, a total of 87 vegetation indexes were used for PCA based on three periods.Taking the PCA of Suyu 41 in 2019 as an example, the contribution ratio and eigenvalue of the first principal component were the largest, at 35.8% and 31.1, respectively (Table 3).The vegetation indexes that made greater contributions to the first principal component were NDVI, RVI, DVI and other vegetation indexes closely related to coverage (Supplementary Table S3).Similarly, the contribution ratios of the second, third and fourth principal components with higher eigenvalues were 19.4%, 17.5% and 12.9%, respectively.The cumulative contribution of the total variation of the first to fourth principal components is 85.7%, which is enough to represent a large part of the original indicator information.Therefore, we can objectively analyze the growth of maize with reduced nitrogen application under different densities through four independent comprehensive indicators.

PCA of Vegetation Index of Different Treatments of Reduced Nitrogen Application at Different Densities
In order to eliminate the factors with less influence and greater interference and improve the accuracy of measurement data analysis, the single index introduced above was converted into a smaller number of more effective indicators [33].Therefore, in the experiments in 2019 and 2020, a total of 87 vegetation indexes were used for PCA based on three periods.Taking the PCA of Suyu 41 in 2019 as an example, the contribution ratio and eigenvalue of the first principal component were the largest, at 35.8% and 31.1, respectively (Table 3).The vegetation indexes that made greater contributions to the first principal component were NDVI, RVI, DVI and other vegetation indexes closely related to coverage (Supplementary Table S3).Similarly, the contribution ratios of the second, third and fourth principal components with higher eigenvalues were 19.4%, 17.5% and 12.9%, respectively.The cumulative contribution of the total variation of the first to fourth principal components is 85.7%, which is enough to represent a large part of the original indicator information.Therefore, we can objectively analyze the growth of maize with reduced nitrogen application under different densities through four independent comprehensive indicators.

Correlation and Model Construction between Vegetation Index and Grain Yield and Protein Content
From 2019 to 2020, the vegetation indexes of Zhengdan 958 and Suyu 41 were significantly correlated with the final GY (Figures 2 and 3).The GY-prediction model was constructed using regression analysis.During modeling, two-thirds of all data were randomly selected as the training set and the remaining one-third comprised the test set.Comparing the determination coefficients of the regression equations, it was found that the higher the correlation between the vegetation index and GY, the better the GY-prediction model constructed with the vegetation index (Table 7).Comparing the predictive effect and accuracy of the models, it was found that the GY-prediction model constructed with the vegetation index during the R6 period of the two varieties had a higher R 2 and lower mean relative error (MRE), and therefore a good predictive effect during the two-year experiment.The precision of the yield-prediction model for Suyu 41 was better than that for Zhengdan 958.There was a high correlation between RVI, DVI and protein content per unit area in 2019 (Figure 3).The correlation between NDVI and PSRI and the protein content per unit area was high in 2020 (Figure 4).This difference is related to different weather conditions between years.A protein-content-prediction model was established using stepwise regression analysis.By analyzing the R 2 of the regression equations, it was found that the protein-content-prediction model oforf Zhengdan 958 established using the vegetation index was more accurate than that for Suyu 41 (Table 8).Through the regression analysis, it was found that the protein-content-prediction model constructed with the vegetation index during the R4 period of the two varieties had a good predictive effect during the two-year experiment.Compared with the predictive model established for 2020, the predictive model in 2019 had a higher accuracy.

Discussion
Maize yield and total protein content are closely related to its growth.In the past, destructive sampling was used with high investment and low efficiency.Previous studies have shown that due to the growth and development of crops, the spectral characteristics of images of the field will change, which indirectly reflects the growth and development process of crops [33].Crop-growth monitoring has had good progress, but most methods focus on the quantitative relationship between remote-sensing parameters and agronomic indicators in a single period [34,35].In this study, we analyzed the changes in the NDVI, RVI and PSRI of summer maize in the R2, R4 and R6 periods under different treatment conditions.The NDVI [36], RVI [37] and PSRI [38], which have been verified by predecessors and can be used to monitor plant growth, were selected for the detection of plant growth.Under different densities, NDVI and RVI decreased with the decrease in the nitrogen application rate, while PSRI increased with the decrease in the nitrogen application rate.Therefore, we can infer that the vegetation index represented by NDVI, RVI and PSRI will change regularly with the change in treatment conditions under different treatments.
The value of the vegetation index is closely related to crop growth, but it is easily affected by weather and other external environments.Therefore, when weather conditions change, the performance of the same vegetation index is different.As a result, it is difficult to accurately reflect the growth of crops with a single vegetation index.Multiple vegetation indexes are used to comprehensively evaluate the growth of crops with reduced nitrogen application under different densities.However, there is a certain correlation between some vegetation indexes, which makes the crop growth information reflected by them cross and overlap when nitrogen is reduced at different densities, and the importance of comprehensive evaluation of each indicator is also different.Previous studies have also used principal component analysis to transform multiple indicators into a few unrelated comprehensive indicators [39] when evaluating the deep-sowing tolerance [40], drought resistance [41] and cold tolerance [42] of crops.On this basis, the subordinate function method was used for comprehensive evaluation, and the results are more scientific and reliable [43].Hence, in this study, the original multiple vegetation indexes were converted into several independent comprehensive indexes using principal component analysis.On this basis, the corresponding subordinate function value is obtained by using the comprehensive index value to obtain the comprehensive score of reduced nitrogen application at different densities.According to the comprehensive score, it can be seen that at different densities, the growth of crops will be poor with the reduction in nitrogen application.The comprehensive score of each treatment was highly correlated with its yield and protein content.Therefore, by scoring, we can infer that summer maize can maintain a good yield and protein content while reducing N rate by 25% at a high density, which is also consistent with our final measured value.For these reasons, we infer that we can directly use the comprehensive scores of vegetation index under different treatments to compare the yield and quality of maize under different treatments, providing a basis for field management.
The canopy structure of summer maize varied among growth periods, and the different spectral information obtained with the UAV resulted in the precision of predictive models changing a lot.The formation of maize yield and total protein content is a long-term accumulation process.Although the yield and protein content can be estimated using the vegetation index of a single key growth period or a single vegetation index [44,45], multiple vegetation indexes of multiple growth periods in the maize growing season can more fully reflect the formation process of yield and total protein content [46].Therefore, based on the above research, we selected vegetation indicators highly related to yield and protein content for further analysis in the R2, R4 and R6 periods.It was found that the yield-prediction model established in the R6 period and the protein-content-prediction model established in the R4 period can achieve a better prediction.In addition, compared with the predictive model of yield and protein content established by predecessors, the determination coefficient of the model is significantly improved [47,48].The predictive model obtained in this study is greatly affected by weather changes, and the interval between periods is long, so it still has great limitations.In the future, the detection interval should be shortened, and many years of experiments should be carried out to reduce the impact of weather on the prediction results, so as to make the predictive effect more accurate.
The research methods proposed in this study can be realized under different treatments of different crops.The main advantage of this program is that it uses different vegetation indexes as comprehensive indicators for analysis, which improves the application scope and accuracy of UAVs in monitoring crop growth in the field.The predictive model of yield and protein content was established in different periods, which improves the accuracy of the model by defining the best modeling period.It provides an opportunity to monitor the land, and supports the management of agricultural systems and decision-making actions.

Conclusions
According to the results of this study, the vegetation index will change regularly with the change in treatment.NDVI, RVI and PSRI are important vegetation indexes for monitoring crop growth.Moreover, the comprehensive score of each treatment can directly reflect the yield and protein content of maize under different treatments.Furthermore, the NDVI, RVI and DVI obtained at R6 can be used to predict production (R 2 > 0.734 **).The NDVI and RVI acquired at R4 were valuable for predicting kernel protein content (R 2 > 0.769 **).However, the predictive model is easily affected by weather conditions and the cycle interval is long, resulting in significant differences in the accuracy of the predictive model during the experimental years.

Figure 1 .
Figure 1.Location of the study area, the Experimental Farm of Lianyungang (China).

Figure 1 .
Figure 1.Location of the study area, the Experimental Farm of Lianyungang (China).

Figure 2 .
Figure 2. The correlation matrix between all parameters at different growth stages in 2019.

Figure 2 .
Figure 2. The correlation matrix between all parameters at different growth stages in 2019.

Figure 3 .
Figure 3.The correlation matrix between all parameters at different growth stages in 2020.

Figure 3 .
Figure 3.The correlation matrix between all parameters at different growth stages in 2020.

Figure 5 .
Figure 5. Effects of decreasing N rates under different planting densities on RVI in 2019 and 2020.

Figure 5 .
Figure 5. Effects of decreasing N rates under different planting densities on RVI in 2019 and 2020.

Figure 5 .
Figure 5. Effects of decreasing N rates under different planting densities on RVI in 2019 and 2020.

Figure 6 .
Figure 6.Effects of decreasing N rates under different planting densities on PSRI in 2019 and 2020.

Figure 6 .
Figure 6.Effects of decreasing N rates under different planting densities on PSRI in 2019 and 2020.

Table 2 .
Vegetation indexes, formulas, and references used in this study.

Table 3 .
Eigenvalue and contribution of each comprehensive index.

Table 3 .
Eigenvalue and contribution of each comprehensive index.

Table 4 .
Comprehensive index value (Cl), subordinate function value u(X), comprehensive evaluation value (D) and correlation analysis in 2019.

Table 5 .
Comprehensive index value (Cl), subordinate function value u(X), comprehensive evaluation value (D) and correlation analysis in 2020.

Table 6 .
Effects of densities and decreasing nitrogen rates on yield and protein content of summer maize in 2019 and 2020.
Note: The different alphabet in each column of the same cultivar means the significant difference between treatments at the 0.05 level.ns, not significance, ** Significance at 0.01 probability level, the same as below.

Table 7 .
Construction of GY model for different growth periods of summer maize in 2019 and 2020.

Table 8 .
Construction of protein-content-prediction model for different growth stages of Zhengdan 958 in 2019 and 2020.