Predicting rice (Oryza sativa L.) canopy temperature difference and estimating its environmental response in two rice cultivars, ‘Koshihikari’ and ‘Takanari’, based on a neural network

ABSTRACT Canopy photosynthesis is an important component of biomass production in field-grown rice (Oryza sativa L.). Although canopy temperature differences (CTD) provide important information for evaluating canopy photosynthesis, the measurement of CTD is still a labor-intensive task. Therefore, we designed this study to establish a model for predicting CTD under different field conditions using meteorological data and evaluated the environmental response of CTD using the established model. Our study collected 2,056,264 CTD data points from two rice cultivars having different photosynthetic capacities, ‘Koshihikari’ and ‘Takanari’, and then used these data to create a novel model using a neural network (NN). The input variables were limited to meteorological data, and the output variable was set to CTD. The established NN model produced a prediction accuracy of R2 = 0.792 and RMSE = 0.605°C. We then used this NN model to simulate the CTD response of the Koshihikari and Takanari cultivars in response to various environmental changes. These predictions revealed that Takanari had a lower CTD than Koshihikari when exposed to high relative humidity (RH) or low to moderate solar radiation (Rs ). In contrast, the CTD of Koshihikari tended to be lower than that of Takanari under lower RH or higher Rs . This result implies that the advantages of the single-leaf gas exchange system in Takanari can be mitigated under extremely high-VPD conditions. Thus, our new method may provide a powerful tool to gain a better understanding of gas exchange, growth processes, and varietal differences in rice cultivated under field conditions. Graphical Abstract

Photosynthesis can be evaluated using thermal imaging techniques as leaf/canopy transpiration is reflected in leaf/canopy surface temperature (Gates, 1968;Jones, 2014). Given this, many studies have evaluated changes in the canopy surface temperature to evaluate the crop canopy status of open fields. Takai et al. (2010) revealed that canopy surface temperature is strongly related to the photosynthetic rate and stomatal conductance of the outer leaves, while Horie et al. (2006) reported that canopy diffusive conductance under field conditions could be estimated based on the surface temperature of both shaded and sunlit canopies. The canopy transpiration rate was estimated in short-time intervals based on the canopy surface temperature and the heat balance model (Monteith & Unsworth, 2013) for rice (Kondo et al., 2021), soybean (Hou et al., 2019), and cotton (Jones et al., 2018). All these techniques assume that the leaf surface temperature is strongly affected by the latent heat flux of transpiration. Therefore, canopy temperature difference (CTD), represented by the difference between air and canopy surface temperature, is fundamental to the use of thermal imaging techniques. However, the measurement of CTD remains a challenging task because of the labor-intensive and technically sophisticated evaluations required to complete these evaluations. Furthermore, the use of thermal imaging techniques has spatiotemporal limitations. The canopy surface temperature measured using thermal imaging devices are typically evaluated as a single image, and the capturing range of thermal imaging devices is limited. Even though this spatial limitation can partly be overcome by combining with unmanned aerial vehicle, it is difficult to take continuously in shorttime intervals. To apply thermal imaging techniques for evaluations of plants cultivated in a field with greater practicality, these limitations need to be overcome.
Micrometeorological conditions strongly affect canopy gas exchange (Jones, 2014). For instance, favorable temperatures and strong radiation generally accelerate photosynthesis (Choudhury, 1987;Yamori et al., 2014), which leads to stomatal opening and greater transpiration. Furthermore, air humidity affects stomatal aperture and transpiration rates (Monteith, 1995;Morison & Gifford, 1983). Thus, the activity of canopy gas exchange affects the canopy surface temperature and CTD, and it can be assumed that CTD can be estimated from micrometeorological data. Meteorological data are universally utilized in the field of agriculture, easy to access, can be collected in short-time intervals, and cover wide area of lands. Given this, the evaluation of meteorological data may be advantageous to realize easier estimation of CTD and transpiration, thus reducing costs and overcoming the spatiotemporal limitations in the use of thermal imaging techniques.
Neural networks (NNs) are intensively being applied in the field of crop science. For example, NNs have been used to help detect crop diseases (Mujahidin et al., 2021), predict biomass production (Jin et al., 2020;Ma et al., 2019), and grain yield (Das et al., 2018;Haghverdi et al., 2018;Nevavuori et al., 2019;Tanaka et al., 2021). In fact, some reports have even described using NN to estimate crop transpiration rates (Nam et al., 2019;J. Fan et al., 2021). However, these studies were conducted at a limited scale and had difficulties translating these results to large-scale field conditions. The use of thermal imaging techniques may solve these problems to some extent. In addition, to the best of our knowledge, ours is the first study designed to evaluate the application of NN to predict CTD in rice (Oryza sativa L.). The major reason for this is that the accumulation of data describing changes in the canopy surface temperature is a labor-intensive task. We collected the canopy surface temperatures of rice and the corresponding micrometeorological data at 2,056,264 points by establishing a system for estimating rice canopy transpiration rate over short time intervals (Kondo et al., 2021). We then used these data to develop an NN-mediated model to estimate CTD in an effort to increase accuracy and simplify these kinds of evaluations. Besides, the CTD prediction model based on big data may cancel noises/ errors caused by fluctuations of weather conditions (e.g. wind velocity, air temperature). This advantage may contribute to extracting physiological responses to the variations of weather conditions more clearly in rice.
The gas exchange of and genotypic differences in rice under either steady state or simple environmental conditions have been reported numerous times (Ikawa et al., 2017;Qu et al., 2016;Yamori et al., 2020). Many previous studies have reported that the saturated photosynthetic rate of a single indica leaf from the 'Takanari' cultivar was significantly higher than that of various japonica cultivars, including 'Koshihikari' (Adachi et al., 2019b;Hirasawa et al., 2010;Takai et al., 2010;Taylaran et al., 2011). However, little information has been published regarding its canopy level gas exchange and genotypic differences in response to various meteorological conditions. This study focused on two rice cultivars, Koshihikari, the most popular cultivar in Japan, and Takanari, a cultivar having high photosynthetic capacities. With the model established to predict CTD from micrometeorological data, the environmental responses of CTD in these two cultivars under virtual meteorological conditions can be evaluated or simulated. Thus, our model may be useful in revealing differences in the environmental response of the canopy gas exchange in both Koshihikari and Takanari.
The aim of this study was to establish a less laborious method to predict CTD and then to evaluate the difference of environmental responses in two rice cultivars, Koshihikari and Takanari. We first developed an NN model to predict CTD using micrometeorological data for two rice cultivars, Koshihikari and Takanari, and the accuracy of the model was evaluated by comparison with two popular regression methods: multiple linear regression (MLR) and partial least squares regression (PLSR). Different CTD values produced under different meteorological conditions were then simulated for both cultivars and the results of these simulations were then used to evaluate genotypic differences in the response of these cultivars to changing environmental conditions. These experiments allowed us to produce a novel methodology for the easy estimation of CTD and evaluation of genotypic differences in environmental responses using micrometeorological data.

Plant materials
Koshihikari and Takanari were cultivated in a paddy field at the Graduate School of Agriculture, Kyoto University, Kyoto, Japan (35° 2ʹ N, 135° 47ʹ E, 65 m altitude) in 2017, 2018, and 2019 with seeds being sown on April 20, 23 April 2017, 2018, and 19 April 2019, and the seedlings transplanted on May 16, 18 May 2017, 2018, and 17 May 2019. The planting density was 22.2 plants·m −2 and we used 60 kg N ha −1 , 47 kg P ha −1 , and 56 kg K·ha −1 as the basal dressing across all 3 years. In 2018 and 2019, 40 kg N·ha −1 , 31 kg P·ha −1 , and 37 kg K·ha −1 were additionally applied as top dressings at the beginning of July. Weeds, diseases, and insects were strictly controlled, and the field was fully irrigated throughout the growing season.

Data collection
Micrometeorological data from the paddy field were measured using a meteorological data acquisition system from July 18th to 25th over 3 years at 1-s intervals, except for 24 July 2019, because of serious problems with the device. The air temperature (T a ) and relative humidity (RH) were recorded in 2017, using a temperature and relative humidity probe (CS215, Campbell Scientific, Inc., USA) with an aspirated radiation shield. Solar radiation (R s ) was recorded using an albedo meter (PCR-3, Kipp & Zonen, Netherlands). These instruments were set approximately 2 m above the ground (1 m above the canopy) and connected to a data logger (CR-1000; Campbell Scientific, Inc.). In 2018 and 2019, T a and RH were recorded using a composite meteorological sensor (WS500; EKO Instruments, Japan) while R s was recorded using a pyranometer (MS802; EKO Instruments, Japan). All sensors were set at the same height as in 2017 and were connected to a data logger (GL840, GRAPHTEC Co., Japan).
In all 3 years, the canopy surface temperature (T c ) was recorded using an infrared thermal imaging camera (InfReC S30, Nippon Avionics Co., Ltd., Japan) with a resolution of 160 × 120 pixels, with a spectral range of 8-13 μm, and thermal sensitivity of less than 0.2°C at 30°C. The camera was placed 16 m from the paddy field and 7 m above the ground with an elevation of 24°. Thermal images were then recorded every second and all plots were placed in a single image. T c was corrected using reference temperatures and a temperature reference board was placed adjacent to the experimental field. The white felt was attached to a -60 × 40 cm wooden board and the surface temperature of the reference board was simultaneously recorded using a thermometer (TR-52i, T&D Corporation, Japan) and thermal imaging camera. The difference in these two data points was then used to correct the temperature drift of the thermal imaging camera. Our modeling assumed that this corrected for changes in canopy surface emissivity. Once corrected, we extracted the mean T c value from each thermal image with four, two, and three replications in 2017, 2018, and 2019, respectively. CTD was then calculated using the following equation:

Data preprocessing, model establishment, and prediction
The distribution of all observed CTD data is shown in Figure 1(a) with each datasets classified as one of three categories: training, validation, and prediction, as shown in Figure 1(b). All parameters were normalized as follows: where V norm is the value after normalization, V is the observed value, α is the correction term, V max is the maximum observed value, and V min is the minimum observed value. The parameters for normalizing each variable on a 0-1 scale are shown in Table 1. After this preprocessing, the model was established using three computational methods: MLR, PLSR, and NN, using the training and validation datasets. The input variables were T a , RH, and R s , and the output variable was CTD. The MLR and PLSR models were established using the scikit-learn library in Python, version 3.8.5, while the NN model was established using the Neural Network Console software (version 2.0 (Sony Network Communications Inc., Japan). First, the model structure was automatically determined using a structure search function in the Neural Network Console software using the Koshihikari dataset. The model structure with the least validation error was then selected (Figure 2). We then used this structure to train the model for both the Koshihikari and Takanari datasets. In both the structure search and training, the loss function was expressed as a squared error, and the optimizer was Adam. The batch size and learning rate were set as 512 and 0.0001, respectively. The learning process was used to minimize the errors between the estimated and observed CTDs in the training and validation datasets. The epoch size was set to 50, and the validation error was calculated every two epochs. This training revealed that the lowest validation error was observed in the 12th epoch in Koshihikari and the 20th epoch in the Takanari datasets. The model with the least validation error was then saved for each cultivar and used in all the downstream analyses. We then confirmed the efficacy of these models by using them to predict the CTD values for each strain on 22 July 2017, and 21 July 2018, respectively.

Analysis
Prior to completing these predictions, we exposed our data to a correlation analysis focusing on the relationships between T a , RH, R s , CTD in Koshihikari, and CTD in Takanari. We then used the predictions from the MLR, PLSR, and NN models to determine the coefficient of determination (R 2 ) and root mean squared error (RMSE) for the observed versus predicted CTD values. We then evaluated CTD prediction under virtual meteorological conditions using our NN model. These simulations were completed using RH and R s values ranging from 30% to 100% and 0 to 1500 W m −2 , respectively. The T a was fixed at one of three positions: 25°C, 30°C, or 35°C, and then used to predict CTD for both Koshihikari and Takanari using these parameters. Once completed, the simulated CTD values were compared and the difference in CTD between Koshihikari and Takanari (T Taka -T Koshi ) was defined by the equation below:  where CTD Taka and CTD Koshi represent the CTDs of Takanari and Koshihikari, respectively. All analyses were conducted in Python (Van Rossum & Drake, 2009).

Meteorological conditions
The mean daytime T a , mean RH, and cumulative R s from 18 to 25 July in 2017, 2018, and 2019 are shown in Figure 3. The mean T a ranged from 25.43°C on 19 July 2019, to 33.48°C on 18 July 2018. The mean RH ranged from 53.21% on 24 July 2018, to 93.55% on 22 July 2019. The cumulative R s ranged from 6.56 MJ·m −2 ·d −1 on 22 July 2019, to 28.00 MJ·m −2 ·d −1 on 24 July 2018. Generally, 2018 was characterized as hot and dry, while 2019 was cool and humid. 2017 displayed a larger variance in weather and ran the full gambit of these conditions. The correlation matrix calculated for the meteorological data and CTD is shown in Figure 4. This data revealed a strong negative correlation between T a and RH (r = −0.96) and more moderate correlations between T a and R s (r = 0.60) and RH and R s (r = −0.62). Of the three meteorological factors evaluated only CTD, T a , and RH were moderately correlated with CTD in Koshihikari (r = −0.48, r = 0.47, respectively) and Takanari (r = −0.55, r = 0.53, respectively). In contrast, R s was not correlated with CTD in either cultivar (r = 0.014 and r = −0.028, respectively).
The differences in CTD between Takanari and Koshihikari (T Taka -T Koshi ) under various meteorological conditions are shown in Figure 7. The RH, which is set on the x-axis, was restrained to 30% to 100%, while R s was set on the y-axis and restrained to 0 W m −2 to 1500 W m −2 . The T a was set to 25°C, 30°C, and 35°C in Figure 7(a-c), respectively. The plots describing the observed weather conditions in response to a T a between 24.5°C and 25.5°C, 29.5°C and 30.5°C, and 34.5°C and 35.5°C are shown in Figure 8(a-c), respectively. Figure 8 used the same x-and y-axes as Figure 7 to allow for direct comparisons. In Figure 7(a), CTD in Takanari was generally lower than Koshihikari when of RH was over 75% and R s was under 900 W m −2 . A similar weather condition was wholly observed when T a was from 24.5°C to 25.5°C (Figure 8(a)). When the T a was set to 30°C, CTD in Takanari was generally lower than Koshihikari when the RH was over 60% and R s under 1050 W m −2 , while the CTD for Koshihikari was generally lower than Takanari when RH was over 60% and R s over 1050 W m −2 (Figure 7(b)). However, both conditions could be observed when T a was from 29.5°C to 30.5°C (Figure 8(b)). When T a was set to 35°C, CTD was generally decreased in the Koshihikari   Heatmaps showing the varietal difference in the simulated values of canopy surface temperature in Koshihikari and Takanari (T Taka -T Koshi ). In each panel, RH was retrained to between 30% and 100% and R s was sequenced from 0 to 1500 W m −2 . The T a was set to 25°C (a), 30°C (b), and 35°C (c). The color bar represents the T Taka -T Koshi values. Figure 8. Plots describing the weather conditions observed over the course of the entire study. In each panel, RH was restrained between 30% and 100% and R s was sequenced from 0 to 1500 W m −2 . The panels (a), (b), and (c) show the plots of each of these conditions in response to a T a of between 24.5 and 25.5°C, 29.5 and 30.5°C, and 34.5 and 35.5°C, respectively. variant compared to the Takanari cultivar when the RH was under 50% or over 75%, or R s was under 900 W m −2 . However, the CTD in Koshihikari was higher than Takanari when RH ranged from 50% to 75% and R s was over 900 W m −2 (Figure 7(c)). Both conditions were also observed when T a was from 34.5°C to 35.5°C (Figure 8(c)).

Discussion
CTD was positively correlated with T a and negatively correlated with RH (Figure 4), this is likely because higher temperatures are associated with lower RH and higher vapor pressure deficits (VPD), which induces higher latent heat flux, producing increased transpiration and reduced CTD. This result implies that the acceleration of transpiration by VPD exceeds the effect of stomatal closure, although a high VPD generally induces stomatal closure in most plants (Jalakas et al., 2021;Kimm et al., 2020;Ohsumi et al., 2008). R s moderately correlates with T a and higher R s values are believed to contribute to greater transpiration, leading to higher T a and photosynthetic rates. However, there was no direct correlation between R s and CTD in our evaluations. The response of the stomata to sudden fluctuating light is relatively slow (Ohkubo et al., 2020;Sakoda et al., 2021;Tanaka et al., 2019;Taniyoshi et al., 2020), and this may have resulted in a time lag between light fluctuation and observable CTD. The difference between the R s values recorded by the sensor and solar irradiation in the canopies may also have influenced.
The relationships between the observed and predicted CTD values based on MLR, PLSR, and NN modeling were much better than the direct correlation analyses estimating the relationships between various meteorological components and CTD (Figures 4-6). This is because T c is more obviously affected by complex combinations of multiple meteorological components than single meteorological conditions, thus the application of large datasets accumulated over long periods allows a better understanding of these complex relationships. This means that our data models were better equipped to evaluate the relationships between CTD and specific weather conditions.
The predicted CTD values fit the observed values for 2017 and 2018 (Figures 5 and 6); however, the prediction accuracy on 21 July 2018, was lower than that on 22 July 2017. This is likely because 21 July 2018, produced CTD values ranging from −8 to 4°C, reducing the range of these values when compared to other studies, where the minimum CTD value was approximately −4°C (Fukuda et al., 2018;Zheng et al., 2020). It is worth noting that 21 July 2018, was extremely hot and dry, which may have led to a remarkably high VPD (Figure 3). This may have caused a considerable increase in the transpiration rate and thus a lower CTD (Figure 6), which could not be accurately predicted using any of the three computational methods described in this study. This may go some way to explaining the reduced accuracy for extremely low CTD values in this study, as these values fall outside of 94.8% of the total data sets, making these data points very rare (Figure 1(a)).
On 22 July 2017, and 21 July 2018, the slopes of the linear regression were closer to 1 in the prediction model created using NN than those produced using MLR and PLSR (Figures 5 and 6). In addition, the R 2 and RMSE values of the NN were better than those of the other two methods. The advantage of the NN model is its prediction range, where NN accommodated a CTD range of −4 to +2°C, while the MLR and PLSR models only supported a CTD range of −4 to 0°C. Because the majority of the CTD values over 0°C were observed under field conditions, it is thought that the NN prediction model is likely to be the best method for predicting CTD under these conditions. This is further supported by the fact that this NN model can predict CTD, an important indicator of crop growth and physiological condition, using only open-field scale meteorological data. The established NN model had better prediction accuracy within the CTD range from −4 to +2°C (Table S1), which covered 94.8% of all the observed data. Although the CTD prediction model in the present study is thought to be practical under most of weather conditions, further accumulation of data following extreme conditions (i.e. CTD under −4°C or over +2°C) are necessary to improve the robustness of our model.The simulated CTD values for the Takanari cultivar were generally lower than those of Koshihikari when R s was under 900 W m −2 (Figure 7), which was largely consistent with the observed data ( Figure 8) and previous studies showing that the greater leaf photosynthetic capacity of Takanari is supported by high stomatal conductance and leaf nitrogen content (Hirasawa et al., 2010;Ohsumi et al., 2007a), stomatal density (Ohsumi et al., 2007b), hydraulic conductivity (Taylaran et al., 2011), and consequently, low leaf/canopy surface temperatures Ikawa et al., 2017;Takai et al., 2010). However, the simulated CTD values for the Takanari cultivar were significantly increased when compared to the CTD values of Koshihikari when T a was 30°C and R s was over 1200 W m −2 (Figure 7(b)) or when T a was 35°C, RH was between 50% and 75%, and R s was over 900 W m −2 (Figure 7(c)). Many of these behaviors were confirmed in the observed data ( Figure 8) and were largely consistent with the fact that reduced RH increases VPD. The simulated values on CTD in Koshihikari and Takanari under the same set of environments are shown in Figure S1.
Generally, the stomatal aperture in plants decrease under high VPD conditions (Inoue et al., 2021). A previous study reported that stomatal conductance decreased in Takanari when compared with Nipponbare, a typical japonica cultivar, under drought stress (Ohsumi et al., 2008). This is also supported by our previous observations that Takanari demonstrate a midday depression in leaf photosynthesis (Adachi et al., 2019a) and canopy transpiration (Kondo et al., 2021). Thus, the stomata of Takanari are thought to be more sensitive to extremely high VPD conditions than those of Koshihikari. Another possible reason for the reversal phenomenon may be the difference of aerodynamic characteristics in Koshihikari and Takanari. Takanari has larger leaf area, lower canopy height, and greater aerodynamic resistance under windless conditions (r a *) than Koshihikari (Kondo et al., 2021). In other words, canopies of Takanari have structures that tend to prevent heat transfer to and from the atmosphere. In some cases of dry conditions, stomata of plants are closed, where the genotypic differences of photosynthetic traits were almost canceled. Under these conditions, the genotypic difference of heat diffusion attributed to canopy structures is thought to be a dominant factor affecting CTD. The reversal phenomenon was simulated under the conditions of R s over 900 W m −2 and T a at 30°C, but not in T a at 35°C (Figure 8). Further physiological research on Koshihikari and Takanari may be needed to reveal reasons for this inconsistency.
Despite the fact that our novel NN model enables fairly accurate prediction of CTD using only micrometeorological data, it remains limited in terms of genotypic and environmental diversity. This is because our model is specific to only two cultivars, both cultivated in Kyoto, Japan: Koshihikari and Takanari. One option to overcome this limitation is to extend the dataset by collecting CTD data for other cultivars and environments. However, collecting a comparable amount of data for new cultivars and environments is a time-consuming task. Another option for solving this problem is fine-tuning. In fine-tuning, we apply all the layers from the pretrained NN model except for the last one and replace the final layer in the model with a new layer that has weights trained for its new targets. Previous studies have reported that the same level of prediction performance as full training can be achieved using a reduced dataset (Cetinic et al., 2018;Howard & Ruder, 2018;Tajbakhsh et al., 2016). Therefore, this technique may be useful in improving the versatility of this model moving forward.

Conclusion
Our NN model was able to predict CTD under field conditions using only meteorological data with practical accuracy under the conditions of CTD from −4°C to +2°C, covering 94.8% of all the observed data. This model accurately described the CTD values for Koshihikari and Takanari under moderate meteorological conditions and clearly modeled the varietal differences in response to specific environmental conditions. Thus, this newly established method may enable the evaluation of canopy gas exchange and its varietal characteristics in field-grown rice while reducing both labor and cost and expanding area of measurement with high timeresolution. Further collection of datasets under extreme conditions, in novel genotype, and at new location is needed to improve the versatility of this method.