A Machine-Learning Model Based on the Fusion of Spectral and Textural Features from UAV Multi-Sensors to Analyse the Total Nitrogen Content in Winter Wheat

Abstract

Timely and accurate monitoring of the nitrogen levels in winter wheat can reveal its nutritional status and facilitate informed field management decisions. Machine learning methods can improve total nitrogen content (TNC) prediction accuracy by fusing spectral and texture features from UAV-based image data. This study used four machine learning models, namely Gaussian Process Regression (GPR), Random Forest Regression (RFR), Ridge Regression (RR), and Elastic Network Regression (ENR), to fuse data and the stacking ensemble learning method to predict TNC during the winter wheat heading period. Thirty wheat varieties were grown under three nitrogen treatments to evaluate the predictive ability of multi-sensor (RGB and multispectral) spectral and texture features. Results showed that adding texture features improved the accuracy of TNC prediction models constructed based on spectral features, with higher accuracy observed with more features input into the model. The GPR, RFR, RR, and ENR models yielded coefficient of determination (R²) values ranging from 0.382 to 0.697 for TNC prediction accuracy. Among these models, the ensemble learning approach produced the best TNC prediction performance (R² = 0.726, RMSE = 3.203 mg·g⁻¹, MSE = 10.259 mg·g⁻¹, RPD = 1.867, RPIQ = 2.827). Our findings suggest that accurate TNC prediction based on UAV multi-sensor spectral and texture features can be achieved through data fusion and ensemble learning, offering a high-throughput phenotyping approach valuable for future precision agriculture research.

1. Introduction

Nitrogen is a critical nutrient for crop growth, influencing development, yield, and quality [1]. TNC is a primary indicator of a crop’s nitrogen status [2]. Monitoring TNC can provide valuable insights into a crop’s nutritional status and facilitate effective field management practices. Traditional nitrogen content determination methods involve destructive sampling, which is a time-consuming and labour-intensive process requiring significant resources [3]. While non-destructive methods such as chlorophyll meters have emerged to estimate nitrogen content, they do not fully reflect the plant’s overall condition.

In several areas, quantitative remote sensing plays an important role [4,5,6]. With the development of UAV remote sensing technology, several studies have been conducted in the agricultural field using UAV-mounted sensors for soil and crop TNC monitoring in a high throughput and non-destructive manner [7,8,9]. For example, Lopez-Calderon [8] demonstrated the effectiveness of estimating the whole nitrogen content of forage maize based on UAV multispectral imagery data, and Liu [10] demonstrated the effectiveness of the inversion of the leaf nitrogen content of winter wheat based on UAV RGB imagery. The commonly used RGB sensor provides information in three high-resolution bands, while the multispectral sensor has five bands with more sensitive band information [11]. These two sensors are widely used in agriculture because of their small size, low cost, simple data processing, and easy disassembly and installation for portability [12,13]. Red, near-red, and thermal infrared bands have been found to perform well for crop monitoring using UAV remote sensing techniques [14]. However, the obtained spectral information and vegetation indices may perform poorly due to soil information and large canopy biomass [15]. Therefore, combinations of spectral features with varying sensitivities are chosen to achieve highly accurate prediction data. Furthermore, most studies have used only a single sensor to demonstrate the effectiveness of UAV remote sensing technology in predicting nitrogen content. Nonetheless, there are a limited number of studies that have explored the fusion of multi-source sensor data from UAV remote sensing to determine total crop nitrogen content.

Texture information is an essential complement to remotely sensed imagery that helps identify important features of objects or regions in an image. It is commonly used in image classification [16,17]. Different nitrogen treatments can affect crop growth, resulting in plant height differences, structural differences, and changes in leaf size and colour, ultimately leading to changes in texture features in spectral imagery [18,19]. Texture features have been used in vegetation identification and classification [20], nitrogen inversion [21] and condition detection [22]. However, most of the previous studies have only analysed RGB texture features. There are fewer studies that comprehensively assess nitrogen content using both RGB and multispectral texture features.

In recent years, the use of machine-learning methods to automatically detect patterns from data and make predictions about unknown data has become increasingly common in data-intensive fields [23]. These algorithms can effectively solve multivariate non-linear agricultural problems with good results [24]. For example, Li [25] demonstrated that a hyper-spectral inversion model based on the random forest algorithm was interpretable, generalisable, required few samples, did not overfit, and had high accuracy (validation area test accuracy R² = 0.73) for estimating rice canopy nitrogen content. Berger [26] combined machine learning regression to estimate crop nitrogen content and found that the Gaussian process regression model accurately simulated aboveground nitrogen. Zhang [27] demonstrated the applicability of ridge regression in the field by introducing ridge regression analysis to spectral detection methods for crop nitrogen nutrient monitoring. Mahajan [28] demonstrated the feasibility of ENR in estimating the nutrient status of mango leaves through machine learning modelling. These studies provide ample evidence that RFR, GPR, RR, and ENR have good accuracy in agricultural monitoring. Compared to individual machine learning models, ensemble learning models have better accuracy. As one of these models, stacking regression combines multiple weak learning models to obtain a more comprehensive model that performs well in datasets with different sample sizes [29]. Stacking regression is an ensemble learning model that improves accuracy by combining multiple individual learners and capturing their best features [30]. The variety and adequacy of individual learner selection ensure that the information between learners complements each other, which is crucial to obtaining correct results [31]. Stacking regression methods are widely used in agriculture, for example, to estimate chlorophyll content in potatoes [32], to evaluate nitrogen content in citrus leaves [33], and to estimate alfalfa yield [30], where accuracies higher than those of individual machine learning models were obtained. In particular, for understanding the spatial distribution of soil organic carbon (SOC) content in different climatic regions, the stacking method outperforms the single machine learning model for estimating nitrogen content [34]. This method can compensate for the deficiencies of the basic learning model effectively. Although the interpretability of the stacking method may decrease when using multiple basic stacking tools, the increase in model prediction accuracy is significant, so further research on the performance of the stacking model is needed. To date, there are no studies on predicting winter wheat nitrogen content using ensemble learning methods with stacked spectral and texture features from multiple sources.

In summary, this study aims to achieve two main objectives: (1) to assess the efficacy of spectral and texture features obtained from multi-source sensors mounted on UAVs, through data fusion, for predicting TNC during the heading stage of winter wheat, and (2) to develop ensemble learning models that can enhance the accuracy of TNC prediction compared to individual machine learning models.

2. Materials and Methods

2.1. Experimental Area and Design

This study was conducted at Qiliying Comprehensive Experimental Base of the Chinese Academy of Agricultural Sciences, located in Xinxiang City, Henan Province (113°45’38”E, 35°8’10”N) (Figure 1), from 2020 to 2021. The study site has a temperate continental monsoon climate. Figure 2 illustrates the variations in average daily temperature, rainfall, and radiation during the winter wheat growing season. As depicted in the figure, the highest temperature and radiation intensity are observed in May, whereas the highest rainfall occurs in March. The lowest temperature is recorded in January, and the radiation intensity is lowest in November but increases as the winter wheat grows. Rainfall is mainly concentrated in November and March, followed by December, April, and May.

The trial area consisted of 180 plots with three treatments, N1 (300 kg·hm⁻²), N2 (180 kg·hm⁻²), and N3 (60 kg·hm⁻²), applied at two fertility stages: the jointing and the heading stage (

Table 1. Information about the fertiliser.

Treatments	Jointing Stage (kg·hm−2)	Heading Stage (kg·hm−2)	Fertiliser Types
N1	200	100	Urea
N2	120	60	Urea
N3	40	20	Urea

). The total fertiliser application for each plot was proportionally divided into three parts: two at the jointing stage and one at the heading stage. Each nitrogen treatment comprised 60 plots arranged in rows 20 cm apart, with a size of 3 m long and 1.5 m wide, totalling an area of 4.5 m². To ensure objectivity, thirty wheat varieties were selected for this trial, with two replications in each treatment. Pesticide, fertiliser, and irrigation amounts were based on local management practice standards in the field experiment. The TNC data were obtained from collected wheat samples at the heading stage on 23 April 2021. Six representative wheat plants were taken from each plot of uniform growth as wheat samples, leaving only the aboveground parts with a scissor treatment, resulting in 180 wheat samples. The wheat samples were dried and weighed at 80 °C. The sample plants were ground and sieved, and the TNC of the wheat samples was obtained using concentrated sulphuric acid with hydrogen peroxide digestion and a Kjeldahl nitrogen analyser.

2.2. Acquisition and Processing of Spectral Data

This experiment used an M210 (DJI Technology Co., Ltd., Shenzhen, China) UAV with a Red-Edge MX multispectral sensor and a Phantom 4 Pro (DJI Technology Co., Ltd., Shenzhen, China) UAV with an RGB sensor to obtain UAV multi-sensor image data collection (Figure 3).

The DJI M210 is a quadcopter drone with a maximum take-off weight of 6.14 kg and an average flight endurance of about 30 min. It has a full horizontal flight speed of 18 m/s and is equipped with a Red-Edge MX multispectral sensor that has five channels: red, green, blue, NIR, and red edge. These channels have centre wavelengths of 668 nm, 560 nm, 475 nm, 840 nm, and 717 nm, respectively, with bandwidths of 10 nm, 20 nm, 20 nm, 40 nm, and 10 nm. Each channel has a resolution of 1280 × 960 and a field of view of 47.2°. To convert the DN value of the multispectral sensor into reflectance, a calibration plate is required to calibrate the sensor before and after the mission during post-image processing. The DJI Phantom 4 Pro is a quadrotor drone with a maximum take-off weight of 1.38 kg, a full horizontal flight speed of 20 m/s, a full ascent speed of 6 m/s, and a flight endurance of approximately 30 min. Mounted on this drone is an RGB sensor with a resolution of 3000 × 4000 and a 94°FOV of lens. Both UAV missions took place on 23 April 2021 between 11:00 and 14:00, a period of clear and cloudless weather that allowed the avoidance of shadows as much as possible. The missions were flown at an altitude of 30 m, with a heading overlap of 85% and a collateral overlap of 80%. Each sensor uses GNSS (Global Navigation Satellite System) with millimetre accuracy to accurately record the set ground control points (GCPs) location for later geo-correction. The sensor uses a photo mode with equal intervals of vertical ground photography.

2.3. Pre-Processing of UAV Images

In this study, separately acquired UAV multispectral and RGB images during the same heading period were imported into Pix4DMapper Pro software (Pix4D SA, Switzerland) and aligned using a feature point matching algorithm. Firstly, a sparse point cloud of the flight area was generated based on UAV image and position data. Then, a spatial grid was created based on the sparse point cloud, and the spatial coordinate information of the ground control points (GCPs) was added. Thirdly, a sparse point cloud with precise positions was generated, and the surface geometry and spatial texture information of the flight area were created. Finally, high-definition digital orthophotos (DOM) and digital surface models (DSM) of the flight area were generated, and the processed images were exported as TIFF images. For RGB images, histogram equalisation was used to enhance the contrast and brightness, and noise was removed and sharpening was used to enhance detail and clarity. For multispectral images, radiometric calibration images with known reflectance were used to correct the images for radiometric calibration, and DN values were converted to reflectance. Illumination differences between different bands were corrected to improve data quality. ArcMap 10.5 (Environmental Systems Research Institute, Inc., RedLands, CA, USA) was used to divide the high-definition digital orthophotos into plots, create shapefile files, and divide them into 180 areas with IDs to separately obtain spectral information for the corresponding ID areas. The raster calculator in ArcMap 10.5 was used to calculate the corresponding spectral index according to the band information of RGB and multispectral, respectively. The vegetation index of the corresponding cell was extracted according to the ID areas. To minimise edge effects on the image, the shapefile was created to omit the image edge areas and cropped to obtain the required image for the experiment, which was then imported into ENVI 5.3 (Exelis Visual Information Solutions, Inc., Boulder, CO, USA) for texture feature extraction. The mean of all feature pixel values extracted according to the ID was used as the corresponding feature.

2.4. Spectral and Textural Features

Texture information of the multispectral and RGB images was extracted using the widely used greyscale co-occurrence matrix (GLCM), based on the wavelength information of the spectral images. ENVI 5.3 software was utilised to extract texture information for both types of images, which included mean (ME), variance (VA), homogeneity (HO), contrast (CO), dissimilarity (DI), entropy (EN), second moment (SE), and correlation (COR). Additionally, 21 TNC-sensitive vegetation indices were calculated from the spectral reflectance of the multispectral images. The average DN values of the digital images’ three channels were normalised, and the three channels of red, green, and blue were labelled as R, G, and B. Subsequently, the variables r, g, and b were obtained by normalising the DN values of the digital images’ three channels, and six more TNC-sensitive vegetation indices were calculated based on these three digital image variables.

Table 2. Spectral and textural features of the RGB sensor.

Data Type	Feature	Formula	Source
RGB	r	$\mathrm{r}=R/(R+G+B)$	/
	g	$g=G/(R+G+B)$	/
	b	$b=B/(R+G+B)$	/
	Visible atmospherically resistant index	$VARI=(g-r)/(g+r+b)$	[35]
	Ground-level image index	$GLA=(2\times g-r+b)/(2\times g+r+b)$	[36]
	Green, red vegetation index	$GRVI=(g-r)/(g+r)$	[37]
	Excess red index	$EXR=2\times g-r-b$	[38]
	Normalised difference index	$NDI=(r+g)/(r+g+0.01)$	[39]
	g/r	$\mathrm{g}/r=g/r$	[40]
	r/b	$r/b=r/b$	[40]
	Grey-level co-occurrence matrix	ME, HO, DI, EN, SE, VA, CO, COR	[41]

/ empirical visible vegetation index, ME mean, HO homogeneity, DI dissimilarity, EN entropy, SE second moment, VA variance, CO contrast, COR correlation.

shows the spectral and textural characteristics of the RGB images, while

Table 3. Spectral and textural features of multispectral sensor.

Data Type	Feature	Formula	Source
MS	Chlorophyll vegetation index	$CVI=(NIR\times R)/G^2$	[42]
	Colouration index	$CI=(R-B)/R$	[43]
	Canopy chlorophyll content index	$CCCI=(NIR-RE)/(NIR+RE)/(NIR-R)/(NIR+R)$	[44]
	Chlorophyll index Red-edge	$CIRE=(NIR/RE)-1$	[45]
	Green difference vegetation index	$GDVI=NIR-G$	[46]
	Normalised difference vegetation index	$NDVI=(NIR-R)/(NIR+R)$	[47]
	Green NDVI	$GNDVI=(NIR-G)/(NIR+G)$	[48]
	Normalised difference red-edge	$NDRE=(NIR-RE)/(NIR+RE)$	[49]
	Green soil adjusted vegetation index	$GSAVI=1.5\times (NIR-G)/(NIR+G+0.5)$	[50]
	Green optimised soil adjusted vegetation index	$GOSAVI=(NIR-G)/(NIR+G+0.16)$	[51]
	Nitrogen reflectance index	$NRI=(G-R)/(G+R)$	[52]
	Green ratio vegetation index	$GRVI=NIR/G$	[53]
	Normalised red-edge index	$NREI=RE/(NIR+RE+G)$	[54]
	Normalised NIR index	$NNI=NIR/(NIR+RE+G)$	[55]
	Modified normalised difference index	$MNDI=(NIR-RE)/(NIR-G)$	[54]
	Difference vegetation index	$DVI=NIR-R$	[37]
	Renormalised difference vegetation index	$RDVI=(NIR-R)/SQRT(NIR+R)$	[56]
	Soil-adjusted vegetation index	$SAVI=1.5\times (NIR-R)/(NIR+R+0.16)$	[57]
	Optimised SAVI	$OSAVI=1.16\times (NIR-R)/(NIR+R+0.16)$	[58]
	MERIS terrestrial chlorophyll index	$MTCI=(NIR-RE)/(RE-R)$	[59]
	ModifiedNon-linear index	$NLI=1.5\times (NI{R}^2-R)/(NI{R}^2+R+0.5)$	[60]
	Grey-level co-occurrence matrix	ME, HO, DI, EN, SE, VA, CO, COR	[41]

MS Multi-spectral, ME mean, HO homogeneity, DI dissimilarity, EN entropy, SE second moment, VA variance, CO contrast, COR correlation.

presents the corresponding characteristics for the multispectral images.

2.5. Model Framework

To improve the TNC prediction accuracy of ensemble models based on multiple source sensors, this study proposes a stacking-based approach which involves two steps. First, four individual machine learning TNC prediction models were constructed based on multiple source sensor data trained separately: Gaussian Process Regression (GPR), Random Forest Regression (RFR), Ridge Regression (RR), and Elastic Network Regression (ENR). Second, multiple predictions were stacked by an RR learner. These four machine learning models have been evaluated for their applicability in many studies, can be used for TNC prediction, and can be supplemented with more useful information by multiple machine learning models for outcome prediction, which is essential for the construction of ensemble machine learning models. The four individual machine-learning models are briefly described below. GPR is a supervised learning process that estimates the parameters of a regression model by sample learning. It can theoretically approximate any continuous function in a tight space and can be used to solve a variety of engineering problems [61]. RFR is a machine-learning model that contains multiple decision trees and can simulate the relationship between dependent and independent variables based on decision rules. It can handle large numbers of input variables, assess the importance of variables when deciding on categories, produce higher accuracy, balance errors, and mine data quickly [62]. RR is a biased estimation regression method that yields more realistic and better fitting results at the expense of losing some information and reducing accuracy [63]. ENR is a combination of Ridge and Lasso regression, an iterative method that maintains the canonical nature of Ridge, produces reasonable solutions, and does not produce cross-paths [64]. In this study, the grid search method was used to optimise the hyperparameters. The model’s best performance was achieved by exhaustively searching a set of predefined hyperparameter value combinations.

The stacking regression model is an ensemble learning approach that utilises ensemble methods to learn different data features for improved prediction results [31]. Figure 4 depicts the construction of the stacking ensemble learning model. This study employed a 5-fold cross-validation method to partition the datasets into five random and equal parts 80 times. The same division method was applied to different input features, with each one serving as the validation set, while the remaining four parts were the training set. This process was repeated five times, and all the data obtained were utilised as training and validation samples. After constructing predictions for the base machine learning models based on the initial dataset, 5 sets of validation data corresponding to the 5 training sets were generated. These 5 sets of validation data were stacked vertically to obtain the test set prediction matrix, which was further utilised as the test set for the secondary machine learning models. The results of the validation set predictions were averaged to obtain the prediction accuracy of each base machine learning model. The RR model served as an ensemble machine learning model to blend the predictive power of each machine learning model, and a 5-fold cross-validation method was employed to train the RR model. Five validation results were obtained based on the test set prediction matrix, and their mean value was utilised to obtain the final prediction accuracy. Dividing the dataset multiple times in line with the 5-fold cross-validation method facilitates the interpretation of the prediction accuracy of different models and improves the reliability of the predictions.

2.6. Parameters for Model Accuracy Evaluation

This study divided the initial data set into training and validation sets 80 times, using a 5-fold cross-validation approach to train the model. Four hundred test results were obtained after 80 divisions, and the mean values of these test results were used as model accuracy evaluation parameters, including coefficient of determination (R²), root mean square error (RMSE), mean square error (MSE), ratio of performance to deviation (RPD) and ratio of performance to quartile distance (RPIQ). The larger the R², RPD and RPIQ and the smaller the RMSE of a predictive model, the better the predictive ability of the model. The equations for the above four model accuracy assessment parameters are as follows:

$$R^2=1-\frac{{\displaystyle {\sum }_{i=1}^n{({\stackrel{⌢}{y}}_i-y_i)}^2}}{{\displaystyle {\sum }_{i=1}^n{(y_i-\overline{y})}^2}}\begin{array}{ccc}& & \end{array}$$

(1)

$$RMSE=\sqrt{\frac{{\displaystyle {\sum }_{i=1}^n{(y_i-{\widehat{y}}_i)}^2}}{N}}\begin{array}{ccc}& & \end{array}$$

(2)

$$MSE=\frac{{\displaystyle {\sum }_{i=1}^n{(y_i-{\widehat{y}}_i)}^2}}{N}\begin{array}{ccc}& & \end{array}$$

(3)

$$RPD=\frac{SD}{RMSE}\begin{array}{ccc}& & \end{array}$$

(4)

$$RPIQ=\frac{Q_3-Q_1}{RMSE}\begin{array}{ccc}& & \end{array}$$

(5)

where $y_i$ is the observed value, $\widehat{y_i}$ is the predicted value, $\stackrel{-}{y}$ is the mean of the measured values, N is the sample size, SD is the standard deviation of the measured values of the prediction set, is the lower limit of the third quartile, and is the upper limit of the first quartile [65].

In this study, the importance of each basic learning model in the ensemble model was calculated using the importance function. Each basic learning model was assigned a percentage representing its contribution to the explanatory power of the ensemble learning model. The higher the importance percentage, the greater the contribution of the basic learning model.

3. Results

3.1. Sampling Statistics

Table 4. Statistics on the characteristics of TNC in samples from each plot at the heading stage (mg·g⁻¹).

Category	Observations	Min	Max	Mean	SD	Q25	Q50	Q75	CV
All datasets	180	8.26	31.63	20.07	5.70	16.12	19.28	24.74	0.28
N1 dataset	60	15.33	31.63	23.66	4.38	20.44	24.50	26.97	0.19
N2 dataset	60	12.20	30.72	21.28	5.00	17.40	20.68	25.47	0.23
N3 dataset	60	8.26	26.34	15.28	4.07	11.66	15.56	18.01	0.27

SD standard deviation, Q25 lower quartile, Q50 median quartile, Q75 upper quartile, CV coefficient of variation.

shows the TNC values for all test plots and plots under the three nitrogen treatments in the experiment. The mean TNC value for all test plots sampled in this experiment was 20.07 mg·g⁻¹. The mean TNC value of the three nitrogen treatments differed, with the N1 treatment having a significantly higher TNC value than the N2 and N3 treatments, at 23.66 mg·g⁻¹, while the N3 treatment had the lowest TNC value at 15.28 mg·g⁻¹. The range, standard deviation (SD), quantile statistics, and coefficient of variation (CV) for all plots and plots under each N treatment showed significant differences in TNC among the N treatments, and good data separation.

3.2. Analysis of TNC Prediction Accuracy

This study employed four individual machine-learning methods and one ensemble machine-learning method to predict TNC based on RGB and multispectral images of the winter wheat heading stage. The prediction results are presented in

Table 5. TNC prediction accuracy based on different machine learning methods.

Sensor Type	Feature Type	Metrics	GPR	RFR	RR	ENR	Stacking (RR)
RGB	Spectral	R2	0.493	0.382	0.481	0.479	0.511
		RMSE (mg·g−1)	4.273	4.591	4.303	4.401	4.216
		MSE (mg·g−1)	18.259	21.077	18.516	19.369	17.775
		RPD	1.386	1.279	1.374	1.342	1.384
		RPIQ	2.083	1.962	2.069	2.026	2.125
MS	Spectral	R2	0.541	0.465	0.515	0.505	0.551
		RMSE (mg·g−1)	4.013	4.205	4.149	4.174	3.978
		MSE (mg·g−1)	16.104	17.682	17.214	17.422	15.824
		RPD	1.468	1.373	1.420	1.405	1.468
		RPIQ	2.194	2.068	2.113	2.104	2.198
RGB + RGB	Spectral + textural	R2	0.494	0.531	0.509	0.507	0.562
		RMSE (mg·g−1)	4.179	3.955	4.138	4.158	3.947
		MSE (mg·g−1)	17.464	15.642	17.123	17.289	15.579
		RPD	1.401	1.466	1.413	1.395	1.469
		RPIQ	2.156	2.262	2.178	2.165	2.280
MS + MS	Spectral + textural	R2	0.625	0.650	0.630	0.625	0.672
		RMSE (mg·g−1)	3.610	3.536	3.584	3.608	3.415
		MSE (mg·g−1)	13.032	12.503	12.845	13.018	11.662
		RPD	1.641	1.686	1.657	1.645	1.738
		RPIQ	2.483	2.543	2.490	2.478	2.625
RGB + MS + RGB	Spectral + spectral + textural	R2	0.570	0.554	0.554	0.546	0.597
		RMSE (mg·g−1)	3.936	3.881	3.942	3.991	3.788
		MSE (mg·g−1)	15.492	15.062	15.539	15.928	14.349
		RPD	1.504	1.504	1.484	1.468	1.540
		RPIQ	2.256	2.278	2.239	2.224	2.337
RGB + RGB + MS	Spectral + textural + textural	R2	0.651	0.651	0.671	0.662	0.686
		RMSE (mg·g−1)	3.599	3.553	3.468	3.495	3.386
		MSE (mg·g−1)	12.953	12.624	12.027	12.215	11.465
		RPD	1.689	1.680	1.742	1.719	1.765
		RPIQ	2.508	2.494	2.576	2.544	2.628
RGB + MS + MS	Spectral + spectral + textural	R2	0.659	0.675	0.668	0.664	0.699
		RMSE (mg·g−1)	3.487	3.466	3.433	3.438	3.300
		MSE (mg·g−1)	12.159	12.013	11.785	11.820	10.890
		RPD	1.720	1.714	1.745	1.734	1.802
		RPIQ	2.568	2.534	2.579	2.562	2.668
MS + RGB + MS	Spectral + textural + textural	R2	0.666	0.675	0.671	0.675	0.710
		RMSE (mg·g−1)	3.504	3.404	3.435	3.416	3.257
		MSE (mg·g−1)	12.278	11.587	11.799	11.669	10.608
		RPD	1.719	1.713	1.738	1.734	1.802
		RPIQ	2.605	2.637	2.639	2.643	2.746
RGB + MS + RGB + MS	Spectral + spectral + textural + textural	R2	0.670	0.697	0.700	0.692	0.726
		RMSE (mg·g−1)	3.456	3.365	3.352	3.362	3.203
		MSE (mg·g−1)	11.944	11.323	11.236	11.303	10.259
		RPD	1.735	1.769	1.822	1.798	1.867
		RPIQ	2.647	2.731	2.724	2.708	2.827

. Among the individual machine learning models, GPR performed the best when spectral indices of RGB were used as input features (R² = 0.493, RMSE = 4.273 mg·g⁻¹, MSE = 18.259 mg·g⁻¹, RPD = 1.386, RPIQ = 2.083), and GPR also performed the best when spectral features of multiple spectra were used as input variables (R² = 0.541, RMSE = 4.013 mg·g⁻¹, MSE = 16.104 mg·g⁻¹, RPD = 1.468, RPIQ = 2.194). To explore the effect of adding texture information on improving the model’s prediction accuracy, we added texture information of RGB images and multispectral images, respectively. As seen from Table 4, adding texture features to the input variables improved the accuracy for all four individual machine learning models, with the most significant improvement seen for the RFR model, where R² improved from 0.382 to 0.531. Adding multispectral texture features also improved the accuracy for all four individual machine learning models, with RFR performing the best, and R² improved from 0.465 to 0.65. In this study, the spectral and texture features of RGB images and multispectral images were used again as input variables, and all four individual machine learning models performed best when based on the spectral and texture features of multispectral images and texture features of RGB images. The RFR and ENR models had the largest R² of 0.675, and the RFR model had the smallest RMSE and MSE at 3.404 mg·g⁻¹ and 11.587 mg·g⁻¹, respectively, indicating that the RFR model performed the best. In this study, the spectral and texture features based on RGB and multispectral images were used as input variables again, and the results showed that all four individual machine-learning models had the highest accuracy, and the RR model was the best TNC prediction model with an R² of 0.7, RMSE of 3.352 mg·g⁻¹, MSE of 11.236 mg·g⁻¹, RPD of 1.822, and RPIQ of 2.724. Furthermore, the prediction results of individual machine learning models were used to build an ensemble machine learning model using the Stacking (RR) method. As shown in Figure 5, the ensemble machine learning models were more accurate than the four individual machine learning models when constructed with the same input features. Among the ensemble machine learning models, the model constructed based on the spectral and texture features of RGB and multispectral images had the highest accuracy (R² = 0.726, RMSE = 3.203 mg·g⁻¹, MSE = 10.259 mg·g⁻¹, RPD = 1.867, RPIQ = 2.827) and was the best TNC prediction model. The distribution of the degree of importance of the results of the four individual machine learning models corresponding to an ensemble machine learning RR model when constructed based on the seven input features is shown in Figure 6. The results of the RFR model had the highest weight in all seven ensemble machine learning models, indicating the higher importance of individual machine learning models with high accuracy when building ensemble machine learning models and highlighting the close relationship between the performance of individual machine learning models and ensemble machine learning models.

3.3. Analysis of TNC Observations and Predictions

Figure 7 displays the observed and predicted values of the best TNC prediction model constructed based on each of the seven input features. The R² value was 0.511 when RGB spectral features were used as input features, which improved to 0.562 with the addition of RGB texture features. Similarly, the R² value was 0.551 when multispectral spectral features were used as input features, and it improved to 0.672 with the addition of multispectral texture features. When three of the spectral and texture features of RGB and multispectral were combined as model input features, it was discovered that the TNC prediction model constructed based on multispectral spectral and texture features and RGB texture features had the highest R² value of 0.71. In contrast, the TNC prediction model built based on RGB spectral and texture features and multispectral spectral features had the smallest R² value of 0.597. Notably, the TNC prediction model based on the spectral and textural features of RGB and multispectral had the largest R² value. Finally, the R² value of the TNC yield prediction model based on all three feature combinations was higher than that of the TNC yield prediction model based on a single feature or a combination of two features.

This paper presents a comparative analysis of the accuracy of TNC prediction models constructed based on multiple features, and it is found that the ensemble machine learning model (RR) constructed based on spectral and texture features of RGB and multispectral achieves the best TNC prediction accuracy. Therefore, the model was used to generate the predicted TNC distribution, as shown in Figure 8. The t-test analysis of TNC between different N treatments is presented in

Table 6. t-test under different nitrogen treatments.

Feature	t	p-Value
N1 vs. N2	3.847	0.000
N1 vs. N3	9.416	0.000
N2 vs. N3	5.654	0.000

, and all p-values were less than 0.001, indicating significant differences in TNC between the three N treatments in the order of N1 > N2 > N3. The predicted yields’ distribution showed that the TNC of the N1 treatment ranged from 15 to 31 mg·g⁻¹. Based on the measured TNC results, it was observed that the N1 treatment had the highest TNC, ranging between 15 and 32 mg·g⁻¹, followed by the N2 and N3 treatments. These results are consistent with the predicted TNC distribution of the ensemble machine learning model (RR), and demonstrate that the model can be used for winter wheat TNC estimation.

4. Discussion

4.1. Analysis Based on Multi-Source Spectral Features and Texture Features

In this study, four individual machine learning models (GPR, RFR, RR, and ENR) and a stacking (RR) ensemble learning model were constructed for predicting TNC based on UAV RGB and multispectral image data. The results showed that prediction models constructed based on multispectral spectral features had higher accuracy than those constructed based on RGB spectral features. Similarly, models constructed based on multispectral spectral and texture features had higher accuracy than those constructed based on RGB spectral and texture features. This is because multispectral images have five bands and more near-red and red-edge bands than RGB images [66], which provide richer data information and higher accuracy, consistent with Furukawa’s findings [67]. Texture features were also added to the spectral features in this study, and the accuracy of each model improved with the addition of texture features, consistent with Wang’s findings [68]. Both RGB spectral features and multispectral features, when combined with corresponding texture features, outperformed TNC prediction models constructed using spectral features alone. The improvement in accuracy was greater after adding multispectral texture features to multispectral spectral features than after adding RGB texture features to RGB spectral features. This may be because multispectral texture features extracted in this study based on the grey-scale co-occurrence matrix consisted of texture information in five bands, which contain more information sensitive to TNC. Additionally, Zhang and Liu [18,69] found that texture features were unsuitable as independent remote sensing variables for constructing prediction models, so texture features were combined with spectral features in this study. This study also used a combination of three spectral and texture features of RGB and multispectral as input variables for the prediction model, and found that the accuracy of the model constructed based on the combination of all three features was higher than the accuracy of the model constructed based on a single feature or a combination of two of the three features, consistent with Fei’s findings [31]. The prediction model based on the combination of multispectral spectral features and texture features and RGB texture features had the highest accuracy. In contrast, the prediction model based on the combination of RGB spectral features and texture features and multispectral texture features had higher accuracy than the model based on the combination of RGB spectral features and texture features and multispectral spectral features, indicating that texture features can significantly improve model accuracy compared with spectral features [70]. Constructing these five models by combining all four of the RGB and multispectral spectral features and texture features as input features showed the highest accuracy, demonstrating that data fusion from multiple sensors can produce higher prediction accuracy due to the fact that data acquired by different sensors all contribute to TNC prediction in unique and complementary ways [71].

4.2. Potential for Ensemble Learning Models

The study used four individual machine learning models to predict TNC [72]. However, the calibration of these models required large sample sizes, and in cases of small sample sizes the results obtained by machine learning models are often inconsistent [73]. To address this issue, this study used the Bootstrap method [74] to calculate the mean of 400 times the model prediction results as the accuracy of the machine learning model. This not only solved the problem of the small sample size, but also improved the generalisation ability of the model. As a single machine learning model cannot effectively perform with multiple sources of data [31], this study constructed a Stacking (RR) ensemble learning model by combining four individual learner models. The results showed that when prediction models were constructed based on the seven input features selected for this study, the ensemble learning models outperformed the individual models, with the highest R² of the models improving to 0.726 and the lowest RMSE reaching 3.203 mg·g⁻¹. To construct ensemble models, individual models need to be selected based on their adequacy and diversity [75]. The four individual learning models chosen for this study had good predictive power and low similarity between models, which effectively complemented information during the ensemble [76]. The study found that the ensemble learning model could collect the advantages of several individual learning models, compensate for the limitations and shortcomings of individual models, and produce results with better robustness and generalisation ability in regression prediction [77]. An ensemble learning model, which combines different base models, can reduce the risk of overfitting of a single model, avoid the impact of the curse of dimensionality, and improve the overall generalisation ability of the ensemble model [33]. Overall, each basic learning model can be learned as an independent model or a model of the same type with different parameters or subsets of data for training [65]. After evaluating the contribution of four models (GPR, RFR, RR, and ENR) in the ensemble model, it was found that RFR had the highest contribution to the ensemble learning model among all the models constructed based on different input features. RFR is one of the most commonly used and important machine learning methods in current research. The model can utilise various input features to obtain optimal and balanced results, and has demonstrated good predictive performance in several studies [78,79], further demonstrating the effectiveness of base model selection. While linear regression models are often used as secondary learners in stacking models and have shown good accuracy in most studies, Bayesian model averaging, decision layer fusion, and other secondary learners have also been shown to be effective in various domains [80,81], improving model accuracy. In this study, RR as a secondary learning model also showed good performance.

4.3. Implications and Reflections

In this study, we constructed a TNC prediction model based on spectral image data acquired by UAV-carried RGB and multispectral sensors at the heading stage of the winter wheat canopy. The model achieved reliable prediction accuracy, but there is still room for improvement. To enhance the accuracy of TNC prediction, future studies could consider fusing data from different sources such as thermal infrared and hyperspectral [31,82]. Additionally, this study only used TNC data for winter wheat at the heading stage. Therefore, including trial data for other fertility stages can test the stability and practicality of the proposed method. Furthermore, we constructed a stacking (RR) ensemble learning model based on four individual learning models in this study, which significantly improved the prediction accuracy. The next step of the study should consider two things: (1) adopting multiple ensemble learning models to predict TNC and comparing the accuracy and applicability of each model; (2) adding other well-performing individual learning models as the base model of the ensemble learning model. Finally, this study was limited to 30 varieties of winter wheat from the Yellow and Huaihe River wheat regions, which had certain regional variety limitations. Including winter wheat varieties from other areas as research subjects could further demonstrate the applicability of the proposed method to different regions and materials.

5. Conclusions

In this study, the spectral and textural features of RGB and multispectral data were fused as input features using UAV remote sensing technology. The study investigated the prediction accuracy of GPR, RFR, RR, and ENR models with different input features, as well as the potential for an ensemble of the four base models in predicting winter wheat TNC. The results indicated that the fusion of multi-source spectral and texture features, along with the ensemble learning method, improved the prediction accuracy of winter wheat TNC. The method successfully estimated the TNC of winter wheat at the heading stage under different nitrogen treatments, providing a basis for future evaluations of winter wheat TNC. Future studies should consider conducting trials at different fertility stages and in various growth environments to improve the stability and practicality of the method.