UAV Video Collection

This study used a UAV with a hyperspectral sensor to collect remote sensing images of soybean fields. In order to ensure the quality and accuracy of the data, a detailed design was carried out for different flight altitudes, speeds, and image resolutions to obtain wide coverage and clear and detailed image data.

The choice of flight altitude was set according to the overall situation of the soybean field. A higher flight altitude of 120-150 m is selected in relatively open and uniform fields. This altitude can ensure a larger range of image coverage and effectively reduce local image distortion caused by individual differences in soybean plants during flight. In areas with more complex terrain or uneven plant distribution, the flight altitude is appropriately lowered to 60-80 m to improve the spatial resolution of the image. By flexibly adjusting the flight altitude, the image coverage is ensured to be extensive, and the details are clear, thereby improving the accuracy of subsequent data analysis. In terms of the choice of flight speed, a moderate flight speed is set according to the requirements of image resolution and acquisition accuracy. A moderate flight speed ensures the continuity of image acquisition and data stability and avoids image blur caused by flying too fast. Moderate flight speed control also provides a sufficient time window for multi-band data collection of hyperspectral sensors, which helps fully record data in different bands and ensure data diversity. The UAV flight parameter design is shown in Table 1.

Table 1 is the design of the UAV flight parameters, including the basic parameters such as the model and specifications of the UAV, and clearly states that the maximum take-off weight is 25 kg and the flight time is up to 40 min. The selection of the hyperspectral sensor is crucial. This study uses a hyperspectral sensor that can cover visible light, near-infrared, and short-wave infrared bands, with a band range of 400-2500 nm, ensuring the multidimensionality of the image data. The hyperspectral sensor has a high spectral resolution and can effectively capture the reflectance characteristics of soybean plants in different spectral bands. Hyperspectral data can more accurately analyze the physiological status of soybeans at different growth stages and distinguish the effects of CD stress. In order to ensure the timeliness and representativeness of the data, the image acquisition plan is strictly carried out according to the growth cycle of soybeans. The physiological changes of soybeans at different growth stages have different effects, and the impact of CD stress on soybean growth is especially obvious. In this regard, multiple flight missions were arranged to collect high-frequency images at key growth stages of soybeans, such as seedling, bud, flowering, and maturity. This collected data covers different ecological states of soybean growth and can better reflect the impact of CD stress on soybean crops at various stages. Table 2 shows the UAV payload equipment.

Table 2 lists the UAV-mounted equipment used in studying the effects of CD stress on soybean physiology and yield, including hyperspectral sensors, multispectral sensors, thermal imaging sensors [17], RGB sensors [18], and positioning and navigation equipment. The frequency and timing of data collection were adjusted according to climate and soybean growth conditions. Before and after the physiological changes that may be caused by CD stress occur, the acquisition frequency can be increased to obtain as much data as possible that can reflect the stress effect. In good weather conditions, image errors caused by strong sunlight and cloud changes can be avoided, and flight missions can be carried out as much as possible to ensure the clarity and accuracy of the images. During the flight mission and data collection process, in order to ensure the stability of the flight path and the continuity of data collection, an advanced flight control system [19] is used to ensure that the UAV flies according to the plan through preset flight routes and real-time monitoring, avoiding deviations from the route and repeated collection. The execution of the flight plan is combined with GPS (Global Positioning System) and inertial measurement unit data [20] to perform real-time correction of the image’s geographic location to ensure that each image’s data has accurate geographic coordinate information.

Through precise flight strategies and coordinating hyperspectral sensors, high-quality remote sensing images covering soybean fields were successfully obtained, providing a reliable basis for subsequent hyperspectral data processing, feature extraction, and physiological and ecological analysis. The diversity and high resolution of the data provided sufficient information to support the analysis of the effects of CD stress on soybean growth and laid a solid data foundation for the subsequent establishment of accurate physiological and ecological models and yield estimation models.

Data Acquisition and Calibration of CD Stress

In order to accurately construct a quantitative model of the effects of CD stress on the physiological and ecological characteristics and yield of soybeans, the soil cadmium concentration and soybean plant physiological indicators were systematically collected and calibrated. The sensors carried by drones are primarily designed for remote sensing monitoring of the ground surface and vegetation, but their measurement range does not extend directly below the soil surface. Soil sample collection and analysis were performed independently to ensure an accurate assessment of soil cadmium concentrations. The sampling areas were divided according to the degree of cadmium pollution in the fields, and sampling points were set in each area. The sampling points were evenly distributed, which could reflect the spatial heterogeneity of soil cadmium distribution. The boundary points of the cadmium pollution gradient were appropriately increased, and soil samples were collected at two depths, 0-20 cm and 20-40 cm, at each sampling point. These soil samples will be used for subsequent analysis of cadmium concentrations in the soil. This data is the basis for understanding the impact of cadmium stress on soybeans and will be combined with drone remote sensing image data to build and verify the model. It was also necessary to ensure that the outside world did not contaminate the soil during the sampling process. The soil samples were air-dried, ground, and sieved through a 100-mesh sieve, and the microwave digestion method [21] was used for sample pretreatment. The cadmium concentration of each treated sample was measured by inductively coupled plasma mass spectrometry [22]. This method has high sensitivity and can effectively eliminate interference. Standard reference materials were used for quality control. A blank sample and a standard sample were inserted every 20 samples to ensure data accuracy. Fig. 1 shows the sampling and sample pretreatment diagram.

Fig. 1 describes the collection, processing, and analysis process of soil and plant samples. Representative soybean plants were selected from different cadmium pollution gradient areas, and the samples covered the seedling stage, bud stage, flowering stage, and maturity stage. No less than 30 plants were randomly selected from each area, and their growth status and obvious stress symptoms were recorded. Leaf and stem samples were collected from the plants, and the functional leaves of the stems were measured leaf by leaf using a portable chlorophyll meter. The average value was taken to represent the chlorophyll level of the individual plant. The measured value can directly reflect the effect of CD stress on photosynthesis. The activities of superoxide dismutase [23], catalase, and peroxidase were determined by spectrophotometry [24]. Specifically, it includes making leaf homogenate, centrifuging and separating the supernatant, and reacting with the corresponding substrate system. The spectrophotometer can monitor changes in enzyme activity at a specific wavelength. The thiobarbituric acid reaction method determines the malondialdehyde content, which reflects the degree of lipid peroxidation in the plant and evaluates the effect of CD stress on cell membrane stability. The soluble sugar content is determined by the sulfuric acid-anthrone colorimetric method, and the free amino acid content is determined by the ninhydrin method [25] to reveal changes in the plant’s osmotic regulation ability.

Fig. 1.

Sampling and sample pretreatment.

To ensure the data’s accuracy and relevance, special attention was paid to the positional correspondence between plant and soil samples during sampling. When selecting the 30 plants to obtain chlorophyll and other physiological data, we ensured that their sampling locations were the same as the locations of the soil samples used for ICP analysis. The physiological data of each plant can be directly correlated with the soil cadmium content of its growing environment, providing a solid foundation for subsequent data analysis.

After preliminary inspection of the collected soil cadmium concentration and plant physiological index data, extreme values and outliers were removed. The z-score standardization method was used to normalize the data to eliminate the dimension differences of different indicators. With soil cadmium concentration as the independent variable and plant physiological index as the dependent variable, the Pearson correlation coefficient was calculated to screen the significantly correlated indicators. The focus was on the relationship between cadmium concentration and chlorophyll content, antioxidant enzyme activity, and malondialdehyde content. The nonlinear relationship between soil cadmium concentration and plant physiological indicators was explored by combining multi-distance linear regression and a generalized additive model. The rationality of the model was verified by model fitting residual analysis. R² evaluated the model’s performance and root mean square error (MSE). A segmented regression model was further applied for significant nonlinear relationships to capture the threshold effect of soil cadmium concentration.

Based on the comprehensive analysis of soil cadmium concentration and plant physiological response, the CD stress levels were divided into no stress, mild, moderate, and severe. Each level corresponds to a specific range of soil cadmium concentration and plant physiological index thresholds to ensure that the CD stress level has clear physiological significance. The data was divided into a training set and a validation set. The 10-fold cross-validation method [26] was used to evaluate the stability of the CD stress model and compare the prediction errors of the training set and the validation set to avoid overfitting. To verify the repeatability of the measured data, the sampling was repeated three times in different pollution gradient areas, and the soil cadmium concentration and plant physiological indicators were measured, respectively. The error of the results was controlled within 5% to ensure the reliability of the experimental data.

Hyperspectral Data Preprocessing

Remote sensing image data has errors caused by sensor noise and environmental interference. In order to improve the data quality, the wavelet transform and principal component analysis method were used to remove the noise. The Daubechies wavelet basis is selected to decompose the image signal [27], and the image is divided into sub-band coefficients of different frequencies. In the high-frequency sub-band, the threshold can be estimated according to the noise variance, and the high-frequency coefficients can be clipped using the soft threshold function to retain the main signal information and suppress the noise signal. The corrected high-frequency and low-frequency coefficients are subjected to inverse wavelet transform to reconstruct the image and generate denoised image data. It can reduce the impact of random noise and avoid excessive loss of important details.

The signal-to-noise ratio is further improved through principal component analysis; the covariance matrix of the image data is calculated, the eigenvalues and eigenvectors are extracted, and the principal components of each band are obtained. The principal components are sorted according to their contribution rate, and the principal components with a contribution rate of more than 99% are retained. The rest are considered noise components and removed. Combining the main components to reconstruct the image, the multi-band data is generated after noise reduction. This can remove the secondary components to reduce the redundancy of image data, optimize the correlation of image data, and ensure that the hyperspectral data reflects the vegetation characteristics more accurately.

The geometric deformation of hyperspectral images can directly affect the subsequent geographic information analysis and feature extraction. In order to ensure the accuracy of the image’s geographic location, ground control points and image registration technology are used for geometric correction. Differential GPS is used to obtain high-precision ground control point coordinates, and easily identifiable features in the image are selected as control points, and the corresponding pixel positions are manually marked. Based on the ground control point data, a quadratic polynomial model is used to establish the transformation relationship from the image coordinate system to the geographic coordinate system. The multivariate regression model can consider the image between multiple independent variables (pixel coordinates of the image) and one or more dependent variables (geographic coordinates) to more accurately describe and predict the relationship between them. By fitting this relationship, the conversion coefficient can be obtained to accurately locate each pixel in the image into the geographic coordinate system. The model form is:

In the formula, (A, B) are geographic coordinates, (M, N) are image pixel coordinates, and α₁ and β₁ are fitting coefficients. Compared with the reference image, the nearest neighbor interpolation method is used to resample the pixel values to ensure that the image is not distorted after geographic correction and that the geographic coordinate system can be calibrated at the same time. This process can effectively eliminate the image’s tilt, rotation, and local stretching errors, ensure the image’s spatial accuracy, and provide a basis for comparative analysis of multi-temporal images.

In order to eliminate the interference of atmospheric scattering and absorption on hyperspectral data and make the image data closer to the true reflectance of the ground, the ATCOR (Atmospheric Correction) model [28] is used for atmospheric correction. According to the collection area’s climate conditions and terrain characteristics, the mid-latitude summer atmospheric model was selected, and the parameters such as aerosol optical thickness, water vapor content, and initial ground reflectivity of the study area were input. The radiative transfer model (RTM) [29] was used to calculate the effects of atmospheric molecular scattering and aerosol scattering. The path radiation correction formula is:

Fig. 2.

Hyperspectral image processing flow chart.

Δ is the total radiation received by the sensor, Δ_p is the path radiation, and Δ_s is the surface reflected radiation. After correcting for atmospheric scattering, the surface reflectivity is inverted according to the radiation transfer formula:

P_d is solar radiation, T is atmospheric transmittance, and T_atm is atmospheric transmission coefficient. Interpolation is used to replace bands with strong water vapor absorption in the image to eliminate the influence of abnormal absorption in the hyperspectral data. After inspection, the preprocessed hyperspectral image has a significantly improved signal-to-noise ratio, the noise level is reduced to less than 5% before acquisition, and the geometric positioning error is controlled within 1 pixel. After atmospheric correction, the root MSE between the reflectance and the measured ground spectral data is less than 0.02, which verifies the accuracy and reliability of the data processing process. It provides a high-quality data foundation for subsequent feature extraction and model analysis. Fig. 2 shows the hyperspectral image preprocessing process. Fig. 2 shows the logical relationship of the hyperspectral data processing process, including each step’s input, output, and processing methods.

Feature Extraction and Physiological and Ecological Analysis

Multiple parameters that can reflect the growth status of vegetation and the impact of stress are selected, including the structural sensitive vegetation index, the soil-regulated vegetation index, and a new index based on the blue edge and green edge reflectance characteristics. Through time series retrospective analysis, the changing trend of the index is monitored at different growth stages to capture the specific response under CD stress. In the index construction, the weight of the normalization factor in the NDVI calculation is optimized to reduce background soil reflection and other non-target interference, given the CD stress effect. By calculating the temporal change rate before and after stress by segmenting the time window function, the stress response dynamics at different stages were evaluated, and it was found that the exponential changes in the seedling and mature stages were most sensitive to the cadmium concentration gradient, laying the foundation for subsequent feature analysis.

In order to break through the limitations of single-band reflectance, the cross-fusion method between bands was used to extract CD stress sensitivity features. Combined with the spectral difference method, the gradient calculation of the continuous reflectance in the range of 400-2500 nm was performed to enhance the differences in specific bands. The band combination screening method based on mutual information gain optimizes the sensitivity of different band combinations to CD stress characteristics. The screening results show that the index after the fusion of 680 nm red light, 740 nm red edge, and 810 nm near-infrared bands is the most significant in distinguishing CD stress.

The dynamic spectral index is used to capture the difference in reflectance changes of soybean plants at different times. The characteristic bands of multi-temporal images were differentiated to construct a dynamic difference matrix. The segmented threshold optimization method based on curve fitting was used to capture the inflection points of stress characteristic changes, and the dynamic characteristics of spectral response were accurately extracted. The vegetation reflectance data was transformed from the time domain to the frequency domain for analysis, and the low-frequency characteristic components caused by CD stress were extracted with the help of the Fourier transform model. The original spectral curve is decomposed by sliding the window to remove high-frequency noise during the processing. Through spectrum reconstruction, the key low-frequency components reflecting the stress effect are retained, and the changes in their energy concentration are analyzed. The effects of CD stress on the reflectance spectrum of soybeans are mainly reflected in the enhancement of low-frequency components and the contraction of spectrum width.

Using the physiological index data sampled in the field, a multivariate nonlinear model is constructed and coupled with the hyperspectral reflectance characteristics for analysis. Using the kernel function optimized support vector regression algorithm, the relationship between reflectance and physiological data was nonlinearly fitted based on the radial basis function. After model training and verification, the R² value between the red edge reflectance feature and the antioxidant enzyme activity was obtained, proving its significant correlation. In order to explore the interaction between multiple indicators, the partial least squares regression method was applied to map the reflectance features to the low-dimensional latent variable space and calculate the feature weight contribution rate. To verify the significance of the extracted features, the non-parametric test method Mann-Whitney U [30] was used to conduct a significance analysis of the feature distribution under different stress levels. Based on the significance test results, the critical threshold of the effect of cadmium concentration on vegetation index was defined in combination with segmented regression analysis, and the distribution range of the characteristic parameters under cadmium concentration stress of 0, 0.5, 1, 2, and 5 mg/kg was clarified.

Construction of Hyperspectral Recognition Model

To accurately classify soybeans under CD stress, a recognition model based on hyperspectral data was constructed, using machine learning algorithms such as Transformer, support vector machine, and random forest suitable for hyperspectral data processing. This model was combined with feature screening and optimization technology to improve classification effect and model stability.

Hyperspectral data contains hundreds of continuous bands with a lot of redundant information. In order to reduce the computational complexity and highlight the characteristics of CD stress, the maximum correlation and minimum redundancy method is used for feature screening. The mutual information value between each band and the CD stress classification label is calculated to evaluate the correlation between the band and the classification target. Then, the mutual information value between each band is used to evaluate the redundancy between bands. The most representative feature bands are selected by maximizing the correlation between the band and the target and minimizing the redundancy between bands.

On this basis, further dimensionality reduction is performed, and the high-dimensional feature space is mapped to the low-dimensional principal component space by combining the principal component analysis method, retaining the principal components with more than 99% cumulative variance contribution rate. The principal component analysis method can not only effectively reduce data redundancy but also retain the main information structure of hyperspectral features, providing concise and effective input features for subsequent model training. When constructing a random forest model, the performance is optimized by adjusting key parameters. The number of decision trees is set at intervals of 50-300 steps to determine the balance between classification accuracy and computational efficiency. The maximum tree depth can be set to a dynamic adjustment mode, and the error curve of the validation set can be used to avoid overfitting. The number of small sample splits can be adjusted to balance the generalization ability and accuracy of the tree. The random forest model constructs multiple decision trees by randomly selecting feature subsets and sample subsets to enhance the robustness of the model to hyperspectral data noise.

The model is fully evaluated using a confusion matrix, classification accuracy, recall rate, and F1 value indicators. The real and predicted distribution of each type of sample in the classification results can be calculated to clarify the classification effect of the model on different CD stress levels. The diagonal elements in the confusion matrix can be observed to evaluate the model’s accuracy for the main classification target. It can be defined as the ratio of the number of correctly classified samples to the total number of samples, which can measure the overall classification ability of the model. The off-diagonal elements analyze the misclassification pattern of the model. For small sample categories, such as the severity level of CD stress, the recall rate is calculated to reflect the model’s ability to identify small categories. Combining precision and recall rate, the F1 value comprehensively evaluates the balanced performance of the model in each classification task. Based on cross-validation, the model’s generalization ability is further verified using an independent test set. Comparing the performance with the training set ensures that the model does not overfit or underfit.

Yield Estimation Model

In order to accurately establish the relationship between CD stress and soybean yield, it is necessary to construct a comprehensive input feature set, combining multi-source data, including key spectral features of hyperspectral images, soil cadmium concentrations, and plant physiological and ecological indicators. The feature extraction process focuses on representativeness and simplicity to ensure the effectiveness of model training and generalization capabilities. Band combinations and derived indices sensitive to CD stress, such as the structure-sensitive vegetation index and red edge spectral features, were screened from hyperspectral images. These features can reflect the physiological changes and stress levels of vegetation. Soil cadmium concentration was used as the core variable to characterize stress intensity, and physiological and ecological characteristics such as soybean chlorophyll content, antioxidant enzyme activity, and soluble sugar content were integrated to enhance the characterization of crop stress responses. Multicollinearity analysis was used to calculate the variance inflation factor between features and eliminate redundant variables with significant collinearity to optimize the feature combination. A recursive feature elimination algorithm further screened the most important variables to ensure that the input feature set had minimal redundancy and maximum representativeness.

A variety of machine learning and statistical modeling methods were used to explore the quantitative relationship between CD stress and soybean yield and optimize the model’s performance. When constructing the benchmark regression model, the input feature set was used as the independent variable, and soybean yield was used as the dependent variable. The model fitted the regression coefficients using the least squares method to explain the direct linear contribution of the explanatory variables to yield. In order to enhance the explanatory power, the interaction term regression was further introduced to capture the interaction effect between cadmium concentration and plant physiological characteristics. Given the nonlinear relationship in the data, the radial basis kernel function [31] was used to support vector regression modeling. The method maps high-dimensional space and transforms complex nonlinear relationships into linear fitting problems, which can optimize hyperparameters and improve model fitting accuracy. Random forest constructs multiple regression trees in an integrated learning manner and builds a robust model through random sampling of features and sample subsets. Its advantage is that it can automatically capture the nonlinear relationship between features and yield and evaluate the importance of features. The model complexity and generalization ability are balanced by adjusting the number of regression trees and the maximum depth.

A nonlinear regression model based on a multi-layer perceptron [32] with input, hidden, and output layers is constructed. Nonlinear activation functions are used to improve the model’s ability to fit complex relationships. Stochastic gradient descent is used to optimize weight parameters, combined with an early stopping strategy to avoid overfitting. The grid search method [33] and the Bayesian optimization method [34] are used to jointly adjust key parameters. The optimal parameter combination is determined through the validation set results. The results of multiple models are integrated through weighted averaging or stacked regression methods, combining the advantages of linear and nonlinear models to improve prediction accuracy and stability. Stacked regression uses the secondary model to learn the prediction results of the primary model and further optimizes the yield estimation.

Evaluation of the Fit of the Hyperspectral Recognition Model

The coefficient of determination measures the model’s ability to fit the data. The closer the R² value is to 1, the stronger the model’s explanatory ability is. A lower value indicates that the relationship between hyperspectral data and soybean characteristics has not been fully explored. The R² value can well reflect the degree of explanation of hyperspectral characteristics on soybean physiological and ecological indicators and stress responses.

In the formula, α_i is the actual value, α^i is the predicted value, and α¯i is the mean. The root MSE measures the deviation between the model’s predicted value and the true value. The smaller the value, the more accurately the model fits the data. It is suitable for the quantitative relationship between different physiological indicators and spectral characteristics.

Performance Evaluation of CD Stress Yield Estimation Model

The mean absolute error (MAE) measures the average absolute deviation between the predicted value and the true value. The smaller the MAE, the higher the model’s accuracy in predicting yield. It can evaluate the overall error margin of the model in predicting soybean yield and reflect the practical applicability of the method.

The relative error represents the ratio of the prediction error to the true value, which can explain the model’s adaptability to different production ranges and provide a relative measure of the prediction error, making it easier to compare model performance under different conditions.

The cross-validation method evaluates the model’s robustness to assess its generalization ability on unseen data, avoiding overfitting and underfitting problems. The 10-fold cross-validation can comprehensively test the model’s stability under different data distributions to ensure high reliability.

Data Appropriateness Evaluation

The signal-to-noise ratio [35] evaluates the noise level to assess the quality of the denoised hyperspectral data. The noise level should be less than 5% of the noise level before acquisition. Too high a noise level can interfere with extracting CD stress features. The purpose is to ensure that the raw data is fully processed and provides high-quality input. The geometric error is the spatial registration error of the remote sensing image. The geometric correction error of the remote sensing image is controlled within 1 pixel to avoid the deviation between the spectral characteristics and the actual distribution of objects and ensure the accuracy of the spatial positioning of the data. The spectral reflectance after atmospheric correction must be consistent with the ground-measured data to reflect the accuracy of the correction algorithm. An error that is too large can lead to failure of trench feature extraction.

Comprehensive Performance Evaluation

The model interpretability is to analyze the feature weight distribution of the model to clarify that the hyperspectral band and remote sensing features contribute the most to the detection and yield estimation of CD stress, which can help understand the data rules and is beneficial to further optimize the model. Time efficiency evaluates the time cost of data collection, preprocessing, and model training to ensure the operability of the method in practical applications, especially its practicality in large-scale agricultural management.