The experiment was performed in the test site of Shanghai Academy of Agricultural Sciences in Fengxian District, Shanghai City, China (30.891°N, 121.359°E), as shown in Figure 1 . Fengxian District is located in the alluvial plain of the Yangtze River Delta. It has a subtropical marine monsoon climate and the average annual temperature is about 15.8°C.

A total of 209 distinct cultivars of lettuce, consisting of 160 green and 49 red varieties, were randomly assigned to 209 plots. Each plot, measuring approximately 4 m² (4m x 1m), contained approximately 24 plants of each cultivar, as illustrated in Figure 1 . All the lettuces were sown on September 24, 2020, and were transplanted to the field on October 22, 2020. The harvest period commenced around December 29, 2020.

Figure 2 shows the minimum and maximum temperatures from Dec 15, 2020 to Jan 8, 2021. The temperatures from Dec 15, 2020 to Dec 28, 2020 and from Jan 2, 2021 to Jan 5, 2021 range from 0 to 15°C. Lettuce leaves remain undamaged at temperatures near freezing, but are susceptible to damage at temperatures below freezing (Turini et al., 2011). The first day when temperature went below 0°C after field-planting was Dec 29, 2020, on which day the temperature dropped to -5°C. The lowest temperature, which was -9°C, appeared on Jan 7, 2021. During these two periods, from Dec 29, 2020 to Jan 1, 2021 and from Jan 6, 2021 to Jan 8, 2021, different cultivars of lettuce suffered from frost damage of different degrees, as shown in Figure 3 .

The field investigation of FDI of each lettuce cultivar plot was carried out on Jan 8, 2021. The FDI was defined by referring to the statistical method of other damage indexes, such as chilling injury index(Fernandez-Trujillo et al., 1998; Porat et al., 2000; González-Aguilar et al., 2004; Zhao et al., 2006; Yang et al., 2011; Pan et al., 2016) and leaf disease index(Wang et al., 2015), which involved in a variety of crop stress investigations. In this research, the severity of frost damage was graded according to the characteristics of frost damage, which was shown in Supplementary Table 1 . Then the FDI was calculated using the following equation:

where i is the damage level value, n is the number of damage levels counted in this plot, P N i is the number of plants at damage level of i in this plot, M A X _ D L is the maximum damage level, T P N is the total number of plants in this plot.

Two aerial surveys were performed on Dec 15, 2020, when the temperature had not dropped below 0°C and there was no frost damage, and Jan 8, 2021, when the lettuce had suffered from low temperature and frost damage happened. In each survey, the UAV-based RGB and multispectral images of lettuce field were acquired respectively.

The RGB images were captured by a quadrotor named DJI Phantom 4 RTK (SZ DJI Technology Co., Shenzhen, China), which is a compact and lightweight UAV with a 20-megapixel RGB camera. The flight height of flight route was set as 30 m and the corresponding ground resolution of images was 0.012 m, and the forward and side overlaps were 80% and 70%, respectively.

The multispectral images were acquired using a five-band multispectral camera with a resolution of 1280 × 980 pixels, RedEdge-M (Micasense Inc., Seattle, WA, USA), mounted on a DJI M600 Pro UAV (SZ DJI Technology Co., Shenzhen, China). The central wavelength of each band with the corresponding bandwidth was 475nm (20nm), 560nm (20nm), 668nm (10nm), 717nm (10nm), and 840nm (40nm). The flight height of flight route was set as 30 m and the corresponding ground resolution of images was 0.018 m, and the forward and side overlaps were 85% and 80%, respectively. The images of reference panel which is a gray board with 50% reflectance by the size of 15.5 cm × 15.5 cm were captured before flight missions for radiometric calibration.

The RGB and multispectral images of the whole lettuce field were generated by mosaicking the originally acquired images within the aerial survey area using Pix4Dmapper Pro (PIX4D, Lausanne, Switzerland). The multispectral images were radiometrically corrected before mosaicking according to the images of reference panel. All the mosaiced images were geometrically corrected based on the RGB image of Dec 15, 2020 using ArcGIS (ESRI, Redlands, CA, USA).

The growth of lettuces almost stopped during Dec 15, 2020 to Jan 8, 2021 and their changes were mainly caused by frost damage. On Jan 8, 2021, many lettuces had reduced their coverage and lost biometric features because of frost damage, making it difficult to separate them from the background. Therefore, we chose the image on Dec 15, 2020 to extract pure lettuce regions for image features calculation from images on Jan 8, 2021. As NDVI is one of the most sensitive indexes to vegetation cover(Liu et al., 2020; Murphy et al., 2020), decision tree classification was performed on the multispectral image (MSI) of Dec 15, 2020, when the lettuces were in health status, to separate the lettuce plants from the background by setting NDVI greater than 0.5 as the rule using ENVI (Harris Geospatial Solutions, Inc., Broomfield, CO, America). Majority/Minority analysis was then applied to the classified results to reduce small plaques. Regions of interest of lettuce plants of each plot were generated based on the classification results and converted to shapefiles which were used to extract image and spectral features from the RGB and multispectral images of Jan 8, 2021 through zonal statistics tool.

Pearson correlation analysis method was used to select the RGB and multispectral image features. The correlation analysis was employed between these features and FDI in each group. The features that reached a highly significant level (p<0.01) were selected. Three regression algorithms, namely MLR, SVM and neural network (NN), were used to establish the estimation models of lettuce FDI by taking selected color, texture and spectral features as independent variables.

MLR is a basic method in multiple regression analysis and widely used in remote sensing monitoring because of their good theoretical basis(Zhang et al., 2021). In this study, the MLR models were constructed based on the following calculation formula:

where Y i ,   p and ϵ i are the FDI of the ith sample, the number of independent variables and the ith independent identically distributed normal error; X i j and β j are respectively the jth independent variable and its coefficient of the ith sample ( j =1, 2, …, p ).

SVM is one of the commonly used regression methods to predict physiological parameters of crops(Shah et al., 2018; Zhang et al., 2022). The essence of SVM is to construct a set of planes or hyperplanes in a high or infinite dimensional space(Cortes and Vapnik, 1995). In this study, linear kernel function and sequential minimal optimization were chosen to construct the SVM models. The value of kernel scale was 1 and the approximations of box constraint and epsilon ranged from 0.18 to 0.2 and 0.018 to 0.02, respectively.

NN is a mathematical model that simulates the brain for information processing(Agatonovic-Kustrin and Beresford, 2000). The NN models are powerful predictive tools for crop growth status(Romero et al., 2018). The network structure of NN includes the number of hidden layers, the number of nodes in each layer, the initialization of weights, the training algorithm, and the learning rate. The NN models in this study used 10 hidden layers, and the training algorithm was the Levenberg-Marquardt algorithm.

The experimental cultivars contained lettuce in both green and red colors. There may be some influence of lettuce color on the predicted results. To understand this influence, the data was divided into three groups according to the color of lettuce: a group of all the lettuce (ALL), a group of green lettuce (GREEN) and a group of red lettuce (RED). Each group was divided into a training set and a test set by stratified random sampling according to the sample ratio of 7:3. Table 3 describes the statistical characteristics of the FDI of the samples.

The procedure of statistical analysis is summarized by the flowchart in Figure 4 . Coefficient of determination (R-squared, R²), Root Mean Square Error (RMSE) and Mean Absolute Error (MAE) were used to evaluate the accuracy of the estimation models established with FDI as the dependent variable. The value of R² (0≤R² ≤ 1) determines the degree of closeness of correlation. The larger R² is, the closer the relationship between dependent and independent variables is. RMSE and MAE are used to measure the deviation between the predicted and actual values. The smaller RMSE and MAE are, the closer the predicted values are to the actual values. Therefore, the closer R² is to 1, RMSE is to 0 and MAE is to 0, the higher the accuracy of model will be. The formulas are as follows:

where y i , y ^ i   and y ¯ are the actual FDI, the predicted FDI and the average of actual FDI, respectively; n is the number of samples.

Supplementary Figure 1 illustrates the correlation between the RGB and multispectral image features and FDI for three distinct lettuce groups. The results reveal that certain CIs demonstrated a significant correlation with FDI, while all texture features exhibited no noticeable correlation with FDI. The majority of MSI features displayed a significant correlation with FDI, and their correlation trends with lettuce FDI remained largely consistent within each group. Notably, the correlations between certain features, primarily CIs, and FDI exhibited opposing trends for GREEN and RED. Subsequent to the correlation analysis (p<0.01), features with absolute Pearson correlation coefficients (absolute r) greater than 0.181, 0.202, and 0.358 were selected respectively for ALL, GREEN and RED for further analysis, as depicted in Figure 5 .

For the ALL group, a total of 32 features were selected, including 15 features from the RGB images and 17 features from the multispectral images. The CIs derived from the RGB images exhibited low correlations with FDI (absolute r< 0.55), with ExR and ExB demonstrating absolute r values higher than 0.5, and the highest absolute r recorded at 0.549. The correlations between various MSI features and FDI displayed considerable disparity (absolute r arranged from 0.183 to 0.717), with R 840 , the four FD_VIs and EVI surpassing an absolute r of 0.6, and the highest absolute r observed for EVI.

For the group of GREEN, 16 RGB image features and 21 MSI features were selected. The absolute r of RGB image features were all above 0.21 (mostly surpassing 0.55) with ExG recording the highest absolute r of 0.738. Specifically, GNI, ExG, GLA, and RGBVI exhibited absolute r values exceeding 0.72. As for the MSI features, FD_VI3, FD_VI4, EVI, FD_VI1, FD_VI2, NDVI, SIPI, ARVI, and R_840 demonstrated absolute r values greater than 0.65. Among them, FD_VI3 attained the maximum absolute r of 0.773.

For the group of RED, a total of 19 features were selected, including 3 RGB image features and 16 MSI features. The selected RGB image features were B, G and R with absolute r of 0.558, 0.523, and 0.492, respectively. They all showed positive correlations with FDI. The MSI features with absolute r above 0.6 were FD_VI1, FD_VI2, ARVI and NDVI, and FD_VI1 had the highest absolute r of 0.626. Additionally, R 668 was the only feature that displayed a positive correlation with FDI.

Notably, the four newly proposed FD_VIs were explicitly correlated with FDI and displayed exceptional performance across the three groups of lettuce data, highlighting their broad utility in accurately evaluating the impact of frost damage.

The estimation models were constructed with selected RGB image features as independent variables and FDI as a dependent variable by using MLR, SVM and NN algorithms. The accuracy of models was demonstrated in Table 4 and Figure 6 . For the models of the group of ALL, the R² and RMSE of training set ranged from 0.465 to 0.56 and 0.014 to 0.132 with MAE between 0.094 and 0.104, and the R² and RMSE of test set ranged from 0.437 to 0.545 and 0.019 to 0.142 with MAE between 0.111 and 0.115. The FDI estimation accuracy of models for GREEN were highest among the three groups with the R² of training and test sets between 0.517 and 0.637, RMSE between 0.012 and 0.131 and MAE between 0.090 and 0.094. The performance of FDI estimation models for RED was neither sufficient nor stable with R² ranging from 0.223 to 0.386 for the training set, but from 0.357 to 0.737 for the test set. For all the three groups, the models using NN algorithm always achieved better accuracy and stability than others with higher R², lower RMSE and lower MAE, and the difference of R², RMSE and MAE between the training and test set were smaller.

The FDI estimation models based on selected MSI features were constructed using MLR, SVM and NN, respectively. Table 5 and Figure 7 illustrate the accuracy of FDI prediction. All the models based on MSI features had better performance than those based on RGB image features. The R² of models by different algorithms for the group of ALL in the training set were relatively close, ranging from 0.6 to 0.651; but the R² of test set varied from 0.484 to 0.639. The models for GREEN had better accuracy than other groups with R² up to 0.718 and 0.665, RSME down to 0.01 and 0.014 and MAE down to 0.080 and 0.090 for the training and test set, respectively. The models for RED were improved by using MSI features, but still not good enough for FDI estimation with R² lower than 0.5. As for the algorithms, NN was still superior to MLR and SVM.

In order to make full use of the information from different data sources, the estimation models were constructed using all selected features from both RGB and multispectral images, as shown in Table 6 and Figure 8 . Again, NN remained the best modeling algorithm for each group. In the training sets, The R² of the NN models for ALL, GREEN, and RED raised to 0.653, 0.722, and 0.592; while the RMSE and MAE decreased to 0.011, 0.009, and 0.009, and to 0.086, 0.074 and 0.082, respectively. Correspondingly, The R² of NN models in the test sets for each group increased to 0.694, 0.715, and 0.575, respectively; the RMSE dropped to 0.014, 0.014, and 0.018; and the MAE reduced to 0.097, 0.093 and 0.113.

Table 7 illustrates the average accuracy of FDI estimation models with different sensors and groups. It can be seen that the predictive performance of models with multisource features were better than those with the single-source features. Therefore, it can be concluded that the multisource features are helpful in improving the accuracy of FDI prediction.

The damage caused by frost on lettuce will change the way how solar radiation interacts with lettuce leaf cells. When frost damage happened, chlorophyll broke down and Photosynthesis was weakened(Choudhury et al., 2019), which would cause the decrease of absorption of visible light, especially the blue and red light which were mainly used for photosynthesis(Gates et al., 1965). In consequence, the reflectance of red and blue bands increased with the FDI. As can be seen in Supplementary Figure 1B and R of RGB image as well as R 475 and R 668 of MSI were positively correlated with FDI. For the same reason, most of the VIs (including NDVI, SR, EVI, ARVI, and SIPI) that had been proven to be related to chlorophyll content on(Susantoro et al., 2017; Chen et al., 2021; Shu et al., 2022) had significant negative correlation with FDI.

The response of green channel of RGB image (G) and reflectance at green band of MSI ( R 560 ) responded quite differently for green and red lettuce. The absolute r of G and R 560 of red lettuce were much higher than that of green one. The possible reason is that red lettuce contained much more anthocyanin which absorbed the green light(Yang et al., 2016); when anthocyanins of lettuce decomposed due to frost damage, the absorption of green light was decreased and the reflection was enhanced; therefore, the positive relations between FDI and G and R 560 of red lettuce were more significant. For the same reason, the anthocyanin reflectance indexes of MSI, ARI1 and ARI2, which reflected anthocyanin content, was more sensitive to FDI of red lettuce. For green lettuce, the change of color from green to yellow caused by frost damage was obvious. Although G and R 560 did not achieve a significant correlation with FDI of the green lettuce, this change could still make the CIs(including GNI, ExG, VARI, ExGR, NGRDI, MGRVI, GLA, RGBVI, and VEG) which enhanced the green component show much closer relation to FDI of the green lettuce than that of the red one.

The destruction of cell structure due to frost damage would cause the decrease of reflectance at NIR band(Gates et al., 1965). As a result, R 840 was obviously negative correlation with FDI of all lettuces. The result is similar to that obtained in an experiment evaluating maize frost damage, where it was found that the frost damage caused a sharp decline of reflectance between 720 and 1350 nm(Choudhury et al., 2019).

Although Pearson correlation analysis showed that R 717 and R 560 were weakly correlated with FDI, the four proposed FD_VIs also showed high correlations with FDI. These improvements were mainly due to the involvement of red edge band (717nm), green band (560nm) and blue band (475nm), on the basis of the near-infrared band (840nm) and red band (668nm). In particular, the correlations between FD_VI4 of all five bands and FDI of green and red lettuce reached the maximum. The possible reason is that canopy spectrum is the result of comprehensive influence of multiple factors such as internal components of leaves and canopy structure, and there are synergistic changes between FDI and the spectrum of these bands.

The correlations between the texture features and the FDI were not significant. One possible reason is that there were many cultivars with different structures in this experiment, and the difference in texture features among cultivars exceeds the difference in texture caused by frost damage.

The lettuce FDI estimation models based on MSI features had better performance than those based on RGB image features, which is consistent with previous studies on the vegetation coverage monitoring and nitrogen accumulation estimation of rice (Zheng et al., 2018; Furukawa et al., 2021). It was because that multispectral imager typically captured more information than RGB cameras(Shu et al., 2022). While RGB cameras capture only three color channels, multispectral imagers capture more pronounced changes in reflectance at multiple bands across the electromagnetic spectrum, including near-infrared bands. This allows for the calculation of vegetation indexes that are highly sensitive to changes in plant health and vigor. This additional information can be leveraged to better detect and quantify frost damage. Some relevant studies have proved that RGB images are generally used to monitor the early growth of crops, and the information about near-infrared band provided by multispectral images is more suitable for the later growth of crops(Marcial-Pablo et al., 2019). In addition, the use of multispectral images in some studies enhanced the monitoring of biodiversity(Tait et al., 2019; Wolff et al., 2023). In this research, the monitored field lettuces were already in the harvesting period, and the changes were no longer obvious as in the early growth stage. As a result, the changes captured by the RGB camera were limited. At the same time, many cultivars of lettuce were selected in this study, the diversity of which was suitable for monitoring with a multispectral camera. However, although multispectral imagers had the advantage of high accuracy, RGB cameras could also be an alternative selection for low-cost monitoring of frost damage, for green lettuce at least. The FDI predictive accuracy was further improved by taking the combination of both RGB and multispectral image features as model independent variables, in which way all the information related to frost damage was fully used and the problem of spectrum saturation can be solved(Zhou et al., 2021; Shu et al., 2022). The conclusion that multisource data fusion can improve model accuracy has also been demonstrated in previous studies on crop monitoring(Jiang et al., 2019; Liu et al., 2019; Zhu et al., 2021).

Compared to the models trained based on RED and ALL, the models trained based on GREEN worked better. The major reason for the low accuracy of RED models was that there were only 49 red lettuce plots, and the sample size of red lettuce was too small which led to severe overfitting in most models for RED. In the future studies, the number of red lettuce cultivars should be increased in the selection of experimental cultivars to improve the sample size of the red lettuce dataset. Since green and red lettuces had different levels of secondary metabolites, they presented different colors and produced different responses after frost damage, as manifested in section 3.1. Therefore, the accuracy of models for ALL were not as good as those for GREEN. But even in the best FDI estimation model for GREEN, the RGB and multispectral image features could explain only about 70% of the variation. The differences between lettuces were not only caused by the different frost damage degrees, but also influenced by the different plant morphologies and green levels among different lettuce cultivars. The accuracy of models built by different regression algorithms did not show significant distinction for ALL and GREEN. NN models tended to be more stable than MLR and SVM with less difference of R² between training and test set.