This experiment was carried out at the Key Laboratory of Agricultural Internet of Things, Ministry of Agriculture and Rural Affairs, Northwest A&F University, from December 2020 to January 2021. Cucumber (Cucumis sativus L. cv Bonai 14-3) was germinated and grown in 72-hole seed trays. The seedlings were transferred to plastic nutrition cups (10×10×9 cm³) with the Danish substrates (Pindstrup, Denmark) when they had two leaves and one heart. An artificial climatic chamber (RGL-P500D, Hefei) with three vertical full spectrum LEDs and an external ultrasonic humidifier provided a growing environment for the seedlings. The environmental settings in the chamber were as follows: the light intensity was 140 μmol·m^-2·s^-1 and the period was 14 h/10 h (light/dark). The air temperature during the light and dark were 25°C and 16°C, respectively. The relative humidity during the light and dark were 60% and 50%, respectively. The vapor pressure deficit during the light and dark were 1.1~1.5 kPa and 0.4~0.5 kPa, respectively. After seven days of acclimatization, the cucumber seedlings were separated into four groups with 24 cups each. Different groups of cucumber seedlings were placed in four different climate chambers. The air temperature of the four different climate chambers were 8°C, 10°C, 12°C, and 14°C, respectively. Other environmental parameters were consistent with the acclimatization period.

All plants were exposed to a low-temperature environment for five days. Subsequently, to evaluate their recovery, the two groups of cucumber seedlings at 8°C and 14°C were transferred to their previous temperature treatment (light/dark: 25°C/16°C) for 5 days.

Eight plants were randomly chosen as samples from each group to measure the ChlF induction curve. Four plants were randomly chosen as samples to measure the ChlF image and the activity and energy conversion of PSI and PSII. Twelve plants were selected for malondialdehyde (MDA) content measurement after each day of low-temperature light exposure.

The measurements were performed between 3 h and 6 h after the beginning of the photoperiod (10:00 h) every day. Eight plants were randomly selected as the samples from each group to measure the ChlF induction curve, and the second fully unfolded leaf from the growth point down was taken as the leaf to be tested. Cucumber seedlings were dark-adapted for over 20 mins prior to measurement. MINI-PAM-II measuring systems (Heinz Walz, Effeltrich, Germany) were utilized to measure the ChlF induction curve. F _o and F _m were determined with a very weak light (470nm,<0.1 μmol·m^-2·s^-1) and a transient saturated light pulse (470nm, >3000 μmol·m^-2·s^-1), respectively (Chen et al., 2021). The samples were then illuminated for 300 s with a continuous actinic light (900 μmol·m⁻²·s⁻¹) to obtain a steady-state fluorescence yield (F). Subsequently, saturating light pulses were applied to achieve F m , . Actinic light was then turned off and a far-red illumination was turned on to measure F o , (Lima et al., 2002). Other ChlF parameters in the induction curve are displayed in Table 1 . Four plants were randomly chosen as samples to measure the ChlF false-color image and the activity and energy conversion of PSI and PSII, and the second fully unfolded leaf from the growth point down was taken as the leaf to be tested. A blue version of IMAGING-PAM-MAX/B measuring systems (Heinz Walz, Effeltrich, Germany) was utilized to collect the images. The false-color images (F _v/F _o, F _o, and F _m) were measured using the standard measurement procedure of the ImagingWin software (Dong et al., 2019) and created by mapping the one-dimensional value to the three-dimensional RGB value using a false-color system. DUAL-PAM-100 measuring systems (Heinz Walz, Effeltrich, Germany) were utilized to measure the rapid light response curves (RLCs) using the standard measurement procedure of the Dual-PAM-100 software (Yang et al., 2020). A saturation pulse after illumination with far-red light for 10 s was used to measure the Pm after dark-adaption. The gradients of light intensity were: 0, 8, 34, 92, 170, 270, 419, 609, 921, 1453, and 2256 μmol·m^-2·s^-1. The parameters used in this paper are displayed in Table 1 .

Twelve plants were selected for malondialdehyde (MDA) content measurement after each day of low-temperature light exposure. The thiobarbituric acid (TBA) colorimetry method was employed (Xu et al., 1993; Ruan et al., 2002). For each plant, 0.5 g fresh leaf material was homogenized with 5 ml of 10% trichloroacetic acid (TCA). 2 ml of extract solution was combined with 2 ml of 0.6% TBA solution and heated in a boiling water bath for 15 mins and then cooled and centrifuged for 10 mins. A UV spectrophotometer (LAMBDA 365, America) was used to measure the absorbance of the reaction mixture at 450 nm, 510 nm, 532 nm, and 560 nm, respectively. The MDA content was calculated according to the following formula.

Where C represents the MDA content, in nmol·g^-1. D 450 , D 510 , D 532 , D 560 represent the absorbance values at wavelengths of 450 nm, 510 nm, 532 nm, and 560 nm, respectively. N represents the volume of extraction liquid. W represents the fresh plant tissue weight.

To facilitate subsequent data analysis, the data set was normalized in the range [0,1] using a linear normalization function defined by:

Where X, Y represents the data to be normalized and normalized, and X _min, X _max represent the minimum and maximum in the data sequence. Meanwhile, we apply the 3σ rule to filter the data for outliers (Lehmann, 2013).

Pearson correlation coefficient is an index that accurately measures the correlation between two variables. In this study, it could be used to examine the relationship between different ChlF parameters. For two different ChlF variables A = [ a 1 , a 2 , … , a n ] and B = [ b 1 , b 2 , … , b n ] , the formula was given by (Wiedermann and Hagmann, 2016):

The value range of P A , B is [-1,1], that is | P A , B | = [ 0 , 1 ] . The correlation between A and B is judged by the following range of values: | P A , B | is extremely strong at 0.8 and 1.0, strong at 0.6 and 0.8, moderate at 0.4 and 0.6, weak at 0.2 and 0.4, and uncorrelated at 0.0 and 0.2.

Statistical analyses of the ChlF parameters were conducted using Levene’s test and the one-way ANOVA test in SPSS (SPSS Inc., Chicago, IL). The Levene’s test can be employed for variance homogeneity. If its p-value is less than 0.05, then we can reject the null hypothesis of variance homogeneity and conclude that there is significant heterogeneity among the groups. The one-way ANOVA test was performed to examine the differences in means among the groups. The F-value obtained from the ANOVA test represents the ratio of between-group variability to within-group variability. A larger F-value suggests a greater difference in means among the groups, indicating a significant effect. Conversely, a smaller F-value indicates a smaller difference in means, suggesting a non-significant effect.

The PCA method is primarily used to decrease the number of variables and discover the relationship structure between similarly classified variables. In doing so, the main input is transformed into principal components (PCs) that are non-correlated. Most of the variations in the original data may be explained by fewer indicators (Dittrich et al., 2021). The percentage of variance in the data that the PCs could explain was set at 95% (Sun et al., 2021).

Uniform manifold approximation and projection (UMAP) was proposed in 2018 to reduce dimensionality and improve visualization (Mcinnes and Healy, 2018). Based on Riemannian geometry and algebraic topology, UMAP not only has the speed advantages of PCA but also retains the local and global structure. The hyper-parameters of UMAP were set as follows: The dimensionality of the target embedding, n = 4; the number of neighbors, k = 5; the minimum allowed distance between points in the embedding space, d = 0. In this study, PCA and UMAP were applied to process ChlF parameters and select the most important variables to assess chilling injury in cucumber seedlings.

The fuzzy C-means clustering algorithm (FCM) is an improvement of the conventional K-means algorithm. When used in conjunction with fuzzy theory, FCM calculates the similarity between the samples and the cluster centers and estimates the probability that a sample fits into a particular category according to the membership (Kamolov and Park, 2021). The FCM algorithm iteratively optimizes the objective function J by updating the membership u_ij and the cluster center V. The objective function J could be expressed as:

Where m represents the fuzzy coefficient, m>1. N represents the sample size. d_ij represents the Euclidean distance from the sample to cluster center.

The traditional FCM randomly selects a sample as the initial centroid, however, this may cause the model to fall into a local optimum. Genetic algorithm (GA) is a popular optimization algorithm inspired by the idea of biological evolution that generates high-quality optimal solutions based on selection, crossover, and mutation processes (Reddy et al., 2020). Thus, the initial centroid of FCM was optimized by GA in the proposed chilling injury classification model. In the optimization process, genetic parameters were set as follows: population 100, crossover probability 0.8, mutation probability 0.2. Taking selection into account, individuals were chosen based on the Euclidean distance between the sample set and the chromosome of the population.

Analyzing the correlations among ChlF parameters provides a preliminary understanding of their mathematical relationships. Figure 4 depicts the correlation coefficient between ChlF parameters under experimental conditions. Only data from the lower triangle was displayed, as the matrix is symmetric. F _v/F _m is strongly positively correlated with F _v/F _o, F _v, Y(II), F m , , Y(NPQ), qN, NPQ, F _m, F, but not with F o , and qL. F _o was only strongly correlated with F o , , but weakly or not correlated with other parameters. Y(NO) was extremely strongly correlated with many parameters, particularly Y(NPQ) and NPQ. There were large or small correlations between all the ChlF parameters, and each parameter had at least a significant or extremely significant correlation with another parameter. It may not be reasonable to manually select one or multiple parameters to characterize the chilling injury. If all ChlF parameters were employed for subsequent data analysis, data redundancy would occur, resulting in intensive calculation and strong data fluctuation. It is necessary to extract the effective components of ChlF parameters.

All the clustering algorithms have limitations of dimensionality because high dimensional data requires more observed samples to produce much density (Allaoui et al., 2020). Methods of dimension reduction (PCA and UMAP) make it easier for these algorithms to cluster the data. Public information on ChlF parameters was extracted using the PCA method, and several comprehensive and independent indicators containing most of the original information were constructed ( Figure 5 ).

The first four PCs accounted for 95.639% of the original information. The first axis (PC1) explained 72.413% of the variance and was dominated by the damage of PSII in cucumber leaves at low temperature. The parameters with a high positive contribution in PC1 are F _v/F _m, F _v, F _v/F _o, and F _m. These parameters were obtained before continuous actinic light. The damage of low temperature to photosynthetic organs can lead to significant changes in these parameters. The reduced photosynthesis in low temperature conditions can lead to the accumulation of excessive energy. Excessive excitation energy reached the reaction center and dissipated through the antenna chlorophyll, which caused a sharp decline in PSII activity, and caused extensive damage to the PSII complex (Zhang et al., 2020). Then, the second axis (PC2), which was related to the PSII activity in cucumber leaves as shown by the loadings of qP, qL, F, and F o , , explained 12.297% of the variance. F and F o , were negatively distributed on the PC2, while qP and qL were positively distributed. Low temperature treatment reduces the PSII activity and increases the risk of photoinhibition, resulting in a decrease in qP (the approximate fraction of open PSII reaction center) and an increase in NPQ. Thermal energy dissipation might not be sufficient to prevent a decline in qP at low temperatures, hence the negative correlation between qP and NPQ was manifested in PC3 (Kitao et al., 2004). The third axis (PC3) accounted for 6.389% of the variance and included the photoprotection ability of PSII in cucumber leaves. PC3 was dominated by F _o and Y(NPQ). Y(NPQ) indicates the degree of PSII photoprotection (Dias et al., 2018). Plants reduced the impact of low temperature stress by increasing the ratio of energy dissipation [Y(NPQ) and Y(NO)]. The fourth axis (PC4) accounted for 4.54% of the variance and was dominated by the efficiency of PSII in cucumber leaves, as it mainly integrated the information of Y(NO) and F _o. In PC3 and PC4, F _o had the largest positive contribution. In PC1 and PC4, Y(NO) had the largest negative contribution. F _o is related to the chlorophyll concentration and the transfer efficiency of excitation energy (Georgieva and Lichtenthaler, 1999). Low temperature treatment reduces the energy transfer efficiency of antenna chlorophyll a to the PSII reaction center (inactive), which could explain the rise of F _o. Considering the result of correlation analysis, F _v/F _m, qP, Y(NO), and F _o were selected as key ChlF parameters sensitive to chilling injury of cucumber seedlings. This finding is supported by the research conducted by Dong et al. (2019) and Aazami et al. (2021). The aforementioned analysis showed the distribution of different ChlF parameters in the PCA-extracted general information and indicated several fluorescence parameters that could better characterize chilling injury.

The four PCs and sensitive ChlF parameters were extracted by PCA, which may be important for early diagnosis of cucumber chilling injury. However, the loose distribution of data transformed by PCA may compromise the precision of subsequent clustering models ( Figure 5 ). Similarly, we analyzed the distribution of ChlF parameters [F _v/F _m, Y(NO), and qP] that contributed significantly to the PCs, but it has the same shortcomings as PCA ( Figure 6A ). UMAP seemed more appropriate for dimension reduction. Compared with PCA, the visualization result was superior when we applied the UMAP ( Figure 6B ). All samples were reduced to four dimensions using UMAP, and the new features of the samples were Feature 1, Feature 2, Feature 3, and Feature 4, respectively. It should be highlighted that UMAP successfully reflects much of the large-scale global structure that is well represented by PCA while also preserving the local fine structure. However, UMAP lacks the robust interpretability of PCA. As a result, we applied PCA to explain the significance of various ChlF parameters in chilling injury and used UMAP to reduce the sample dimensionality.

Different models, including FCM (15 ChlF parameters as input), ChlF-FCM (F _v/F _m, F _o, Y(NO), and qP as input), PCA-FCM (4 PCs as input), and UMAP-GA-FCM (Feature 1, Feature 2, Feature 3, and Feature 4 as input), were constructed for comparison and selection ( Figure 7 ). Using the GA algorithm to optimize the first centroid can reduce the initial movement error and the iterations to a certain extent. Significant differences were observed in the initial movement error of ChlF-FCM and FCM, but they did not require a lot of iterations. The PCA, as a linear method, has the tendency to neglect some important local information. If all ChlF parameters for the datasets in this paper were used, it would result in a dimensional disaster. Correspondingly, the outcomes of models that used them as inputs were not the best. All models can ultimately reach the set error threshold, but the UMAP-GA-FCM method has the smallest initial error and iteration steps. Therefore, the UMAP-GA-FCM method was employed to evaluate the chilling injury in this study.

Taking the four new features obtained by the UMAP method as the input, an evaluation model of cucumber seedlings chilling injury was constructed based on the UMAP-GA-FCM algorithm. The corresponding relationship between the new features and centroids is illustrated in Figure 8 . The UMAP-GA-FCM algorithm was successful at pulling together clusters corresponding to similar samples. Samples from different classes exhibit distinct feature values and centroids, while the centroids and the average values of the features at different classes were generally consistent.

All the samples were categorized into four chilling injury classes according to the classification model. Figure 9 illustrates the severity of cucumber seedlings chilling injury under different low temperatures and time periods. The chilling injury levels of cucumber seedlings subjected to the same treatment may be different. After 1 d of low-temperature treatment, all samples were observed to have slight chilling injury. 25% of the samples treated at 14°C showed signs of slight chilling injury for 1 d or more. The samples treated at 8°C for 2 d and 37.5% of those treated at 10°C for 4 d were categorized as suffering from severe chilling injury. Under the other environmental treatments, the samples were categorized as having moderate chilling injury. Many studies have broken down the extent of chilling injury by temperature and duration, which may be applicable to the lengthening of time at a temperature or the lowering of temperature at a time, but this evaluation standard is constrained when both exist because of individual differences and subjective judgment. The classification model for evaluating cucumber chilling injury that uses a clustering algorithm is more objective and effectively reduces the influence of individual differences.

MDA contents and F _v/F _m, F _m, and F _o false-color images of the samples with different chilling injury classes were analyzed ( Figure 10 ). Low temperature stress causes changes in cell membrane permeability, transforming the membrane from a liquid crystalline state to a gel state, resulting in an increase in MDA content (Amin et al., 2021). Compared to unstressed samples, the average values of MDA contents with slight, moderate, and severe chilling injury samples increased 19.07%, 28.14%, and 67.24%, respectively ( Figure 10E ). The variations were minor at low chilling injury classes, indicating that the relationship between the MDA content and the degree of low temperature stress was nonlinear. Using MDA content to accurately evaluate the chilling injury classes requires the classification model proposed in this research. ChlF images are well suited to visualize spatio-temporal heterogeneity in plants’ responses under stress conditions (Fathi-Najafabadi et al., 2021). For F _m false-color images, they exhibited a color progression from cyan to green and eventually to yellow as the chilling injury intensified. With regard to the F _o false-color images, their color changes were not obvious under different chilling injury classes and different areas of cucumber leaves. The changes in F _v/F _m were attributed to a decrease in F _m and an increase in F _o. For F _v/F _m false-color image, it is uniformly blue for the unstressed sample ( Figure 10A ); it appeared cyan when the cucumber seedling was affected by low-temperature ( Figure 10B ); it gradually became yellow-green, the main leaf veins appeared orange-red chilling spots with the deepening of chilling injury ( Figure 10C ); it turned to red, and some areas were identified as background areas because the F _v/F _m values were less than 0.038 ( Figure 10D ). Spatial-temporal heterogeneity in the three parameters clearly illustrated the onset of cold sensitive symptoms in the vein and edge of leaves that progressed toward the center of leaves. Further research could be performed on the relationship between areas and classes of chilling injury using ChlF images. According to the above results, our model could be validated from MDA content and ChlF images.

The effect of low temperature on the photosynthetic reaction center of cucumber seedlings is not only reflected in PSII but also in PSI, which may be more cold-sensitive and not easily rectified. CEF is an important photoprotective mechanism that plays a crucial role in regulating the distribution of light energy between photosystems, balancing photoprotection, and facilitating the processes of photochemical reactions. Consequently, it is essential to investigate the changes in PSI and CEF under different chilling injury classes based on the classification model. We further analyzed the RLCs of Y(I), Y(ND), and Y(CEF). Figure 11 depicts their distributions on different days at four low-temperature treatments. The parameter values of cucumber leaves remained stable at 0 d but changed rapidly when exposed to low temperatures. Under low-temperature treatment, the photosynthetic assimilation rate of cucumber seedlings decreased significantly and the photosynthetic electron transport chain became over-reduced, with captured light energy exceeding the absorption capacity of the electron transport chain. Excess light energy obstructed electron transfer, resulting in severe photoinhibition in both PSII and PSI. Thus, Y(I) and Y(CEF) decreased significantly. The high thermal dissipation capacity on the donor side [Y(ND)] was maintained to protect PSI from photoinhibition induced by strong light. Y(I) and Y(ND) changed linearly at low light intensity and reached a transition point at 420 to 610 μmol·m^-2·s^-1, and tended to stabilize at 921 μmol·m^-2·s^-1, indicating that the photoinhibition was deepened and Y(ND) was rapidly induced under strong light. Y(ND) differed substantially under different light intensities due to PSI’s limitation on accepting additional electrons for photosynthesis electron transfer on the donor side. Y(CEF) rose steadily with increasing light intensity and reached the maximum at 2250 μmol·m^-2·s^-1. In light of this, four light intensities of the RLCs (420, 610, 921, and 2250 μmol·m^-2·s^-1) were chosen to analyze the changes in PSI and CEF of cucumber seedlings under different chilling injury classes.

The distribution of Y(ND), Y(I), and Y(CEF) at each chilling injury class under the four light intensities are shown in Figure 12 . The three parameters did vary between cucumber seedlings before and after low-temperature treatment, albeit some of them were not significant under different chilling injury classes. When the light intensity of the RLCs was 420 μmol·m^-2·s^-1, compared to unstressed samples, the Y(I) average values of slight, moderate, and severe chilling injury samples decreased by 34.350%, 52.413%, 73.403%, respectively; however, Y(ND) increased by 68.071%, 103.351%, 114.266%, respectively; and Y(CEF) decreased by 30.396%, 48.467%, 56.246%, respectively. When the light intensity of the RLCs was 610 μmol·m^-2·s^-1, compared to unstressed samples, the Y(I) average values of slight, moderate, and severe chilling injury samples decreased by 43.890%, 56.033%, 74.053%, respectively; Y(ND) increased by 61.601%, 84.931%, 91.302%, respectively; Y(CEF) decreased by 37.504%, 54.693%, 60.541%, respectively. when the light intensity of the RLCs was 921 μmol·m^-2·s^-1, compared to unstressed samples, the Y(I) average values of slight, moderate, and severe chilling injury samples decreased by 49.223%, 59.058%, 73.935%, respectively; Y(ND) increased by 54.021%, 71.605%, 76.924%, respectively; Y(CEF) decreased by 32.881%, 56.478%, 60.590%, respectively. when the light intensity of the RLCs was 2250 μmol·m^-2·s^-1, compared to unstressed samples, the Y(I) average values of slight, moderate, and severe chilling injury samples decreased by 45.708%, 61.155%, 75.435%, respectively; Y(ND) increased by 34.094%, 46.387%, 51.869%, respectively; Y(CEF) decreased by 31.599%, 59.590%, 67.052%, respectively. PSI regulates low temperature effects through high heat dissipation, thus, compared to unstressed samples, there is an obvious increase in Y(ND) in samples with slight chilling injury. The distribution of three parameters under different chilling injury classes indicated the validity of our classification model.

Since Y(CEF) is calculated by Y(I), Y(II), and light intensity, a positive correlation between Y(I) and Y(CEF) was consistently observed, and the relatively smaller decrease in Y(CEF) suggested that PSI is subjected to greater photoinhibition and damage compared to PSII. Furthermore, with the increase in chilling injury classes, the average value of Y(I) decreased by 43.293%, 57.165%, and 74.206%, respectively; the average value of Y(ND) increased by 54.447%, 76.569%, and 83.590%, respectively; the average value of Y(CEF) decreased by 33.094%, 54.807%, and 61.107%, respectively.

Low-temperature treatment could decrease the activity of PSI and induce PSI non-photochemical energy dissipation because of the donor-side limitation. When PSI is photo-inhibited, PSII is also more susceptible to photoinhibition due to the reduction in CEF, which aggravates the chilling injury and makes it difficult for plants to recover. The PSI photoinhibition of cucumber seedlings would be aggravated under strong light intensity or low-temperature conditions (Huang et al., 2016). As a result, the Y(I) values of cucumber seedlings with slight chilling injury and above were generally small. The Y(ND) values remained at a high level after low-temperature treatment, and there was no significant difference between moderate and severe chilling injury, demonstrating the donor side of PSI had a high thermal dissipation capacity (Lu et al., 2020). The rates of CEF in cucumber seedlings with moderate and severe chilling injury were largely inhibited, therefore, the Y(CEF) values were not significant. Because of this, PSI is more fragile than PSII and is more likely to be severely damaged when exposed to low temperature. In addition, PSII has a more efficient and dynamic recovery mechanism compared to PSI. Similarly, the PSI related parameters and the MDA contents could distinguish whether cucumber seedlings had been subjected to chilling injury to a certain extent, but they were not reliable in categorizing the chilling injury classes because their mechanisms of damage and recovery involve complex physiological reactions. We recommend that in future investigations, the PSII related parameters should be used to evaluate the chilling injury, and the PSI related parameters and biochemical parameters should be used to analyze the chilling injury.