Damage Assessment of Rice Crop after Toluene Exposure Based on the Vegetation Index (VI) and UAV Multispectral Imagery

: Chemical spill accidents lead to environmental problems, especially for plants. Plant vegetation assessment is necessary after a chemical accident; however, conventional methods can be inaccurate and time-consuming. This study used the vegetation index (VI) extracted from unmanned aerial vehicle (UAV) multispectral imagery for crop damage assessment after chemical exposure. The chemical accident simulations were conducted by exposure of rice at ﬁve growth stages to four levels of toluene. The VI was measured at ﬁve days after damage and 67 days after planting. Physiological characteristics (chlorophyll content and grain yield) were also measured. As a result, the mean normalized difference VI (NDVI) of toluene-exposed rice was signiﬁcantly decreased with respect to toluene exposure concentration increases at most growth stages. Recovery after toluene exposure was lower in rice exposed to higher concentrations at the earlier growth stages. The chlorophyll content and grain yield were also decreased after toluene exposure with respect to increasing toluene concentrations and showed positive correlations with the NDVI. It indicates that the NDVI is capable of reﬂecting the plant response to chemical exposure. Thus, the results demonstrated that the VI based on UAV multispectral imagery is feasible as an alternative for crop monitoring, damage assessment after chemical exposure, and yield prediction.


Introduction
Environmental concern regarding the use of toxic chemicals has increased as these chemicals are increasingly used industrially. There are major hazards with the industrial use of chemicals that are potentially dangerous, most of which are toxic, flammable, or explosive [1]. The chemical industry market share of Korea is the sixth largest share in the world. Approximately 2500 kinds of hazardous chemicals are used in South Korea, and approximately 400 kinds of chemicals are newly introduced annually [2]. The amount of toxic gases used in South Korea is increasing by approximately 39% annually [3]. Due to the increase in the amount of domestic chemical usage, the annual number of chemical accident occurrences has greatly increased. Although the annual number of chemical accident occurrences has decreased since 2016, these numbers are still high, and the damages are increased due to the complexity and diversity of the chemical industries [2].
The release of toxic chemicals leads to environmental problems, especially damage to crops and vegetation. Approximately 8 tons of hydrofluoric acid was spilled from a chemical leak accident in Gumi city of South Korea in 2012, leading to severe damage to plants. More than 190 ha of crops, including rice, 18 ha of vegetation, and 57,000 landscaping trees, were damaged due to leakage of hydrofluoric gas [4]. The economic damage was estimated to be greater than USD 20 million [3]. Environmental risk assessment after chemical accidents is necessary, but it is challenging due to the complexity of chemical accidents and The paddy field was segmented by pots. Twenty rice seeds were planted in one pot 30 cm in diameter. Each pot was transplanted to a paddy field 39 days after seeding (DAS). The rice planted in each pot formed one hill, a unit of exposure in this study. The rice paddy and field plots are described in Figure 1.

Toluene Exposure Simulation
To simulate the chemical accident, rice was exposed to toluene vapor. Toluene exposure simulations were conducted in five different growth stages, i.e., early tillering (ET), late tillering (LT), stem elongation (SE), booting (B), and flowering (F). The decision of the rice growth stage followed the BBCH (Biologische Bundesanstalt, Bundessortenamt und Chemical industry) scale. The BBCH scale is a coding system that identifies the plant growth stages based on phenological development [16]. The details of the treatment period are described in Table 1. For each simulation, four concentrations of toluene exposure (45.1 ppm (T1), 98.9 ppm (T2), 164.9 ppm (T3), and 227.4 ppm (T4)) were simulated and compared to the control (0 ppm (Ctrl)). The control group was not exposed to toluene but moved to the laboratory with toluene-exposed rice to be placed under the same conditions as the toluene-exposed rice. Table 1. Toluene vapor exposure simulation and UAV multispectral imagery acquisition periods presented as days after planting (DAP) and days after damage (DAD). The first multispectral capture was for damage assessment, and the second was for recovery assessment. The toluene exposure concentrations were selected based on the toxicity of toluene to rice. Since there are no data on the toxicity of toluene in the air to rice, the toxicity of

Experimental Design 2.2.1. Toluene Exposure Simulation
To simulate the chemical accident, rice was exposed to toluene vapor. Toluene exposure simulations were conducted in five different growth stages, i.e., early tillering (ET), late tillering (LT), stem elongation (SE), booting (B), and flowering (F). The decision of the rice growth stage followed the BBCH (Biologische Bundesanstalt, Bundessortenamt und Chemical industry) scale. The BBCH scale is a coding system that identifies the plant growth stages based on phenological development [16]. The details of the treatment period are described in Table 1. For each simulation, four concentrations of toluene exposure (45.1 ppm (T1), 98.9 ppm (T2), 164.9 ppm (T3), and 227.4 ppm (T4)) were simulated and compared to the control (0 ppm (Ctrl)). The control group was not exposed to toluene but moved to the laboratory with toluene-exposed rice to be placed under the same conditions as the toluene-exposed rice. Table 1. Toluene vapor exposure simulation and UAV multispectral imagery acquisition periods presented as days after planting (DAP) and days after damage (DAD). The first multispectral capture was for damage assessment, and the second was for recovery assessment. The toluene exposure concentrations were selected based on the toxicity of toluene to rice. Since there are no data on the toxicity of toluene in the air to rice, the toxicity of toluene in the air was estimated using the toxicity of toluene in the soil to rice and the volatilization factor (VF) of toluene. The 168 h LC50 for toluene in the soil to rice was approximately  (1) and (2) [17]. Considering the toluene exposure duration, rice should be exposed to approximately 358.68 ppm of toluene for 2 h. However, the primary exposure pathways of toluene in the air to rice are dry deposition and wet deposition, which affect rice leaves more directly than exposure to toluene in soil. Therefore, this experiment determined that the highest toluene concentration is 227.4 ppm (T4), which is lower than 358.68 ppm. The toluene concentration was measured using a GC-MSD (Agilent 7890A GC system equipped with an Agilent 5975C with Triple-Axis-HED-EM Detector, Santa Clara, CA, USA) and DB-5MS capillary column (30 m length, 0.25 mm i.d., 0.25 µm film thickness, Agilent Technologies (Santa Clara, CA, USA). Toluene sampling was performed every 20 min during exposure by withdrawing 500 µL of gas at the top and bottom of the chamber using a 1 ml gas-tight syringe. The toluene samples were equilibrated for 5 min, and 50 µL of samples were manually injected into the GC instrument. The toluene concentration was expressed as the maximum of the average concentration of the top and bottom of the chamber.

Exposure
where D A , d b , and T are the apparent diffusivity (cm 2 /s ), dry soil bulk density (g/cm 2 ), and exposure interval (exposure duration) (s), respectively. Toluene vapor exposures were conducted in the chamber. The rice pots were moved to the laboratory at each growth stage and exposed to toluene vapor in the chamber for 2 h ( Figure A1). After toluene vapor exposure, the rice pots were placed outdoors until UAV multispectral imagery was captured.

Damage and Recovery Assessment after Toluene Exposure
To assess the damage of toluene vapor exposure to rice, UAV multispectral imagery was taken five days after damage (DAD) from toluene vapor exposure simulation at each growth stage. Table 1 summarizes the multispectral image-capturing period.
Control and toluene-exposed rice were retransplanted in the field to investigate the recovery of rice after toluene exposure. UAV multispectral imagery of retransplanted rice was taken 67 days after planting (DAP). The details of the imaging periods of recovery assessment are described in Table 1.

Instrument Setup and Flight Mission
The high-resolution aerial multispectral imageries were captured by a multispectral camera (Parrot Sequoia, Parrot SA, Paris, France) mounted on UAV (Inspire 2, DJI, Shenzhen, China). The Parrot Sequoia multispectral camera is composed of a 16 MPIX RGB camera and four 1.2 MPIX single-band spectral cameras that detect single spectral bands of green (530-570 nm), red (640-680 nm), red edge (730-740 nm), and NIR (770-810 nm). The sunshine sensor connected to the Parrot Sequoia multispectral camera provided ambient sunshine intensity, GPS, and IMU information during the flight missions. The Parrot Sequoia multispectral camera simultaneously captures 4608 × 3456 pixels of RGB images and 1280 × 960 pixels of single-band images during the flight missions. Each image contains sunshine intensity, GPS, and IMU information from the sunshine sensor to calibrate and construct reflectance maps.
The flight missions were operated automatically by a predefined flight plan defined on Pix4D capture (Pix4D S.A., Lausanne, Switzerland). The flight altitude was 40 m above ground level (AGL), and the capturing gap was 5.6 m. The calculated front overlap and side overlap achieved 85% and 80%, respectively. Multispectral imagery was taken on Remote Sens. 2021, 13, 25 5 of 19 a clear day and near noon (11:00 A.M.~12:00 P.M.) to minimize variation in solar zenith angle and other environmental conditions, such as ambient sunshine intensity. White reflectance calibration panel (MicaSense Inc., Seattle, WA, USA) images were captured to calibrate reflectance before UAV takeoff to correct variations in environmental conditions. The detailed parameters of flight missions are described in Table 2.

Image Processing and Analysis
The visible symptoms of toluene-damaged rice are chlorosis, necrosis, and tip burn on rice leaves. Crop damage can be measured using VIs, which have great potential for index-based crop insurance when damaged. VIs measure plant vigor and greenness and provide a description of overall crop health [15,18]. This study estimated the damage to rice leaves using VIs. This study calculated three VIs: the normalized difference vegetation index (NDVI) [19], the green normalized difference vegetation index (GNDVI) [20], and the soil-adjusted vegetation index (SAVI) [21]. The VIs were selected by considering their plant vigor description capacity and their limitations. Ray et al. (2006) demonstrated that the NDVI and SAVI had the highest correlation coefficient with the leaf area index (LAI) [22]. NDVI is the most widely used vegetation index and describes the vegetation health well. However, it can be insensitive to canopy with a high concentration of chlorophyll [6,23]. GNDVI is more sensitive to the wide range variation of chlorophyll [20]. SAVI is capable of compensating for the influence of soil [6]. Therefore, this study selected these three VIs. Each VI was calculated using Equations (3)- (5).
where NIR, R, and G are the reflectances of near-infrared, red, and green, respectively. L is the soil adjustment factor, where an L value of 0.5 is generally used to reduce soil noise [21]. This study used an L value of 0.5. Parrot Sequoia captures four multispectral image sets of green, red, red edge, and NIR bands during flight missions. These high-resolution aerial images were processed to calculate VIs. The image processing schematic is described in Figure 2. This study used Pix4Dmapper Pro (Pix4D S.A., Lausanne, Switzerland) and R studio version 1.1.456 (R Studio Inc., Boston, MA, USA). The key steps of the image processing algorithm are as follows: (1) radiometric calibration of multispectral images using calibration panel images, (2) generation of single-band reflectance maps by ortho-rectification of the images, (3) extraction of VIs using Equations (3)-(5), (4) background and soil parts of reflectance map removal by removal of pixels that have VI values lower than the threshold value, and (5) mean VI and standard deviation calculation. Among the image processing steps, (1)- (2) are processed in Pix4D mapper Pro, and (3)-(5) are processed in R studio.
follows: (1) radiometric calibration of multispectral images using calibration panel images, (2) generation of single-band reflectance maps by ortho-rectification of the images, (3) extraction of VIs using Equations (3)-(5), (4) background and soil parts of reflectance map removal by removal of pixels that have VI values lower than the threshold value, and (5) mean VI and standard deviation calculation. Among the image processing steps, (1)-(2) are processed in Pix4D mapper Pro, and (3)-(5) are processed in R studio. This study performed a sensitivity analysis to select VIs that can demonstrate crop damage with respect to toluene vapor exposure damage. The sensitivity of each VI is calculated using the slope of the linear relationship between the toluene exposure level and the corresponding VI. This study calculated the slope of the VI using Equation (6) [15,23], where M is the slope of the VI, is the VI value of control rice, and is the VI value of toluene-exposed rice.

Physiological Characteristic Data Collection
Physiological characteristics of rice were measured to observe physiological changes of rice with respect to toluene exposure levels and to evaluate the potential of VIs based on UAV multispectral imagery for crop damage assessment. This study measured the leaf chlorophyll contents, grain yield, and yield components as physiological characteristics of rice.

Leaf Chlorophyll Contents
This study measured leaf chlorophyll contents at 5 DAD from each toluene exposure simultaneously with NDVI measurements at each growth stage. Leaf chlorophyll content was measured in two ways: using a chlorophyll meter (SPAD 502 Plus; KONICA MI-NOLTA Inc., Osaka, Japan) and using a microplate reader (MQX 200; BIO-TEK Instruments INC., Winooski, VT., USA). This study randomly sampled three leaves from each rice pot five days after toluene exposure simulations for leaf chlorophyll content measurements. This study performed a sensitivity analysis to select VIs that can demonstrate crop damage with respect to toluene vapor exposure damage. The sensitivity of each VI is calculated using the slope of the linear relationship between the toluene exposure level and the corresponding VI. This study calculated the slope of the VI using Equation (6) [15,23], where M is the slope of the VI, V I control is the VI value of control rice, and V I t is the VI value of toluene-exposed rice.

Physiological Characteristic Data Collection
Physiological characteristics of rice were measured to observe physiological changes of rice with respect to toluene exposure levels and to evaluate the potential of VIs based on UAV multispectral imagery for crop damage assessment. This study measured the leaf chlorophyll contents, grain yield, and yield components as physiological characteristics of rice.

Leaf Chlorophyll Contents
This study measured leaf chlorophyll contents at 5 DAD from each toluene exposure simultaneously with NDVI measurements at each growth stage. Leaf chlorophyll content was measured in two ways: using a chlorophyll meter (SPAD 502 Plus; KONICA MINOLTA Inc., Osaka, Japan) and using a microplate reader (MQX 200; BIO-TEK Instruments INC., Winooski, VT, USA). This study randomly sampled three leaves from each rice pot five days after toluene exposure simulations for leaf chlorophyll content measurements.
The SPAD 502 Plus meter measures the leaf transmittance of red (650 nm) and infrared (940 nm) regions and translates the measurements to a relative SPAD value, which is proportional to actual chlorophyll contents [24]. This study measured ten SPAD values per leaf sample evenly distributed over the leaf except near the leaf edge. The microplate reader analysis used the same sample as the chlorophyll meter analysis. The leaf chlorophyll contents were calculated using the equations published in previous research [25]. The leaf samples were frozen, dried and homogenized by shaking for 1 min at 30 Hz (TissueLyser II; Qiagen, Hilden, Germany). Chlorophylls were extracted from samples with 1 mL of methanol and shaken for 1 min at 30 Hz, and then the supernatants were isolated by centrifuging the samples for 2 min at 10,000 RCF (Relative Centrifugal Force). The supernatants were transferred to another tube. The remaining chlorophylls in the pallets were re-extracted by adding 1 mL of methanol and repeating the prior sequence. The re-extracted chlorophylls were added to prior supernatants. Two hundred microliters of samples were transferred to 96-well plates. The absorbances of the samples were measured at 652 nm and 665 nm using a microplate reader. The chlorophyll content calculation follows Equations (7)-(10) published formulae of Warren (2008) [25].
where A 652, 1cm , A 665, 1cm , A 652, microplate , and A 665, microplate are absorbance at 652 nm with 1 cm pathlength, absorbance at 665 nm with 1 cm pathlength, absorbance of 200 µL sample at 652 nm, and absorbance of 200 µL sample at 665 nm, respectively. This study calculated the chlorophyll contents per unit area. The leaf area was measured using ImageJ software.

Grain Yield and Yield Components
This study measured grain yield and yield components at 130 DAP. The yield components include the number of panicles per hill, spikelet number per panicle (SNPP), filled grain percentage, and 1000-grain weight. The grain weight was measured after drying at room temperature until the moisture content was 14%. The adjusted grain weight calculation follows Equation (11), provided by Yoshida et al. (1976), where M is percent moisture content of grain and W is grain weight in grams [26].

Data Analysis
This study performed statistical analysis (ANOVA) to evaluate rice damage with respect to toluene vapor exposure. ANOVA and post-hoc analysis were performed on the NDVI, leaf chlorophyll contents, grain yield, and yield components at a 5% significance level using R studio and SigmaPlot version 12.5 (SSI, San Jose, CA, USA).
Correlation analysis between the NDVI and physiological characteristics was performed to evaluate the feasibility of the usage of the VI based on UAV multispectral imagery as a crop damage assessment tool after chemical exposure. The Pearson correlation coefficient (r) between the NDVI measured at 5 DAD and associated physiological characteristics was analyzed using R studio and SigmaPlot version 12.5. For the significance of the correlation analysis, correlation analysis was conducted at a 5% significance level; the null hypothesis was rejected at a p value < 0.05.

The Availability of UAV Multispectral Imagery to Crop Monitoring
The five flights for damage assessment after the toluene exposure using UAV showed the availability of UAV multispectral imagery for crop monitoring. Figure 3 illustrates the orthomosaic RGB images, NDVI map, GNDVI map, and SAVI map at five growth stages of toluene exposed rice. For all VI maps at all five stages, each rice pot had higher VIs compared to the background, so it was possible to be identified well. The spatial resolution of the VI maps was varied between 4.15 cm and 5.22 cm GSD (ground sampling distance) in spite of the actual flight condition variation. Considering the diameter of rice pots are 30 cm, the spatial resolution was sufficiently small to distinguish the rice and background and observe the damage on rice leaves. the availability of UAV multispectral imagery for crop monitoring. Figure 3 illustrates the orthomosaic RGB images, NDVI map, GNDVI map, and SAVI map at five growth stages of toluene exposed rice. For all VI maps at all five stages, each rice pot had higher VIs compared to the background, so it was possible to be identified well. The spatial resolution of the VI maps was varied between 4.15 cm and 5.22 cm GSD (ground sampling distance) in spite of the actual flight condition variation. Considering the diameter of rice pots are 30 cm, the spatial resolution was sufficiently small to distinguish the rice and background and observe the damage on rice leaves. Figure 3. The orthomosaic RGB images, normalized difference vegetation index (NDVI) map, green NDVI (GNDVI) map, and soil-adjusted vegetation index (SAVI) map at five growth stages of toluene exposed rice.

Sensitivity Analysis
Based on the sensitivity analysis results, this study selected the VI that best reflects the change in rice after toluene exposure simulations. The sensitivity of the VI was described as the slope of the VI before and after toluene exposure (Table A1). The higher slope value of the VI indicates the higher sensitivity of that VI to damage [15]. The negative slope value of the VI indicates that, compared to control rice, the corresponding VI of rice decreases after toluene vapor exposure. Figure 4 shows the negative slope values of VIs except for the slope of the GNDVI measured at the SE and B stages. In general, the NDVI slope values were negatively higher than those of the GNDVI and SAVI at most of the growth stages, which means the NDVI best reflects the toluene damage on rice among the three VIs. For this reason, this study used the NDVI for the toluene exposure damage assessment. Figure 3. The orthomosaic RGB images, normalized difference vegetation index (NDVI) map, green NDVI (GNDVI) map, and soil-adjusted vegetation index (SAVI) map at five growth stages of toluene exposed rice.

Sensitivity Analysis
Based on the sensitivity analysis results, this study selected the VI that best reflects the change in rice after toluene exposure simulations. The sensitivity of the VI was described as the slope of the VI before and after toluene exposure (Table A1). The higher slope value of the VI indicates the higher sensitivity of that VI to damage [15]. The negative slope value of the VI indicates that, compared to control rice, the corresponding VI of rice decreases after toluene vapor exposure. Figure 4 shows the negative slope values of VIs except for the slope of the GNDVI measured at the SE and B stages. In general, the NDVI slope values were negatively higher than those of the GNDVI and SAVI at most of the growth stages, which means the NDVI best reflects the toluene damage on rice among the three VIs. For this reason, this study used the NDVI for the toluene exposure damage assessment.

NDVI-Based Rice Damage and Recovery Assessment
The mean NDVI of the control and four toluene-exposed rice (T1~T4) in five growth stages were measured ( Figure 5). The mean NDVI of control rice increased until the B stage and then decreased in the F stage. The mean NDVI of the T1 level (45.1 ppm) of

NDVI-Based Rice Damage and Recovery Assessment
The mean NDVI of the control and four toluene-exposed rice (T1~T4) in five growth stages were measured ( Figure 5). The mean NDVI of control rice increased until the B stage and then decreased in the F stage. The mean NDVI of the T1 level (45.1 ppm) of toluene-exposed rice did not significantly decrease at all growth stages except the LT stage (Figure 5b). The mean NDVI of the T2 or higher levels of toluene-exposed rice significantly decreased at all five growth stages. Figure 6 shows the canopy changes of toluene-exposed rice at each growth stage. Chlorosis, necrosis, and tip burn occurred on toluene-exposed rice. Compared to the mean NDVI of the control of each growth stage, Figure 5f shows the relative mean NDVI. Compared to the control rice, the extent of decreases in the mean NDVI of toluene-exposed rice increases as the toluene exposure concentration increases until the B stage. The mean NDVI of toluene-exposed rice in the F stage was not significantly changed.

NDVI-Based Rice Damage and Recovery Assessment
The mean NDVI of the control and four toluene-exposed rice (T1~T4) in five growth stages were measured ( Figure 5). The mean NDVI of control rice increased until the B stage and then decreased in the F stage. The mean NDVI of the T1 level (45.1 ppm) of toluene-exposed rice did not significantly decrease at all growth stages except the LT stage (Figure 5b). The mean NDVI of the T2 or higher levels of toluene-exposed rice significantly decreased at all five growth stages. Figure 6 shows the canopy changes of toluene-exposed rice at each growth stage. Chlorosis, necrosis, and tip burn occurred on toluene-exposed rice. Compared to the mean NDVI of the control of each growth stage, Figure 5f shows the relative mean NDVI. Compared to the control rice, the extent of decreases in the mean NDVI of toluene-exposed rice increases as the toluene exposure concentration increases until the B stage. The mean NDVI of toluene-exposed rice in the F stage was not significantly changed. The same capital letters refer to the same growth stages. The small letters refer to the significant differences in the mean NDVI among the control and four toluene-exposed rice (α = 0.05). (f) shows the relative mean NDVI of tolueneexposed rice at each growth stage, which is expressed as the proportion to the control group (100% of the control group). The same capital letters refer to the same growth stages. The small letters refer to the significant differences in the mean NDVI among the control and four toluene-exposed rice (α = 0.05). (f) shows the relative mean NDVI of toluene-exposed rice at each growth stage, which is expressed as the proportion to the control group (100% of the control group).
The recovery assessment of toluene-exposed rice was conducted on rice exposed to toluene at the LT, SE, and B growth stages. The control and toluene-exposed rice at LT, SE, and B were retransplanted to paddy fields after toluene exposure ( Figure A2). Then, the NDVI of retransplanted rice was measured at 67 DAP. Figure 7 shows the mean NDVI of retransplanted rice measured at 67 DAP. The NDVI of rice exposed to the T4 level of toluene at the LT stage was not measured because it completely withered after toluene exposure. Except for the T4 level of toluene-exposed rice at the LT stage, all retransplanted rice recovered, and the mean NDVI increased. The toluene-exposed rice at the LB and SE stages showed a similar NDVI pattern: the NDVI significantly decreased with increasing toluene exposure concentration. The mean NDVI of toluene-exposed rice at the B stage measured at 67 DAP increased, but there was no significant difference in the mean NDVI among the control rice and toluene-exposed rice. Comparing each growth stage, the mean NDVI of rice exposed to toluene in the earlier growth stage was lower than the NDVI of rice exposed to toluene at later growth stages with the same toluene exposure concentration.
Remote Sens. 2021, 13, x FOR PEER REVIEW 11 of 20 Figure 6. RGB image of the control group and toluene-exposed rice canopy. All images were captured five days after damage (DAD) at each growth stage.
The recovery assessment of toluene-exposed rice was conducted on rice exposed to toluene at the LT, SE, and B growth stages. The control and toluene-exposed rice at LT, SE, and B were retransplanted to paddy fields after toluene exposure ( Figure A2). Then, the NDVI of retransplanted rice was measured at 67 DAP. Figure 7 shows the mean NDVI of retransplanted rice measured at 67 DAP. The NDVI of rice exposed to the T4 level of toluene at the LT stage was not measured because it completely withered after toluene exposure. Except for the T4 level of toluene-exposed rice at the LT stage, all retransplanted rice recovered, and the mean NDVI increased. The toluene-exposed rice at the LB and SE Figure 6. RGB image of the control group and toluene-exposed rice canopy. All images were captured five days after damage (DAD) at each growth stage.  Figure 7. Recovery of rice after toluene exposure. The left mean NDVI was measured at 5 days after damage at each growth stage, and the right mean NDVI (recovery) was measured at 67 days after planting. The small letters refer to the significant differences in the mean NDVI of 14 retransplanted rice (α = 0.05). The mean NDVI within the same small letter is not significantly different at Tukey's HSD (p > 0.05).

Physiological Characteristic-Based Damage Assessment
The leaf chlorophyll contents after toluene exposure were analyzed to demonstrate the leaf color changes after toluene exposure. The measurements were conducted at 5 DAD simultaneously with NDVI measurements at each growth stage. Figure 8 shows that the leaf chlorophyll decreases after toluene exposure at each growth stage. Similar leaf chlorophyll patterns were observed in chlorophyll contents measured using a microplate reader ( Figure 8a) and a chlorophyll meter (Figure 8b). The SPAD value was highly correlated with chlorophyll content (R = 0.8381; p < 0.0001; n = 23) ( Figure A3), as described in previous studies [22,27]. In the early three growth stages (ET, LT, and SE stages), the leaf chlorophyll contents and SPAD value significantly decreased with respect to the toluene exposure concentration increases. Compared with the two tillering stages (ET and LT stages) and SE stage, the decreases in chlorophyll contents regarding toluene exposure concentration increase were relatively greater at the two tillering stages. In the latter two growth stages (B and F stages), leaf chlorophyll contents did not change significantly after toluene exposure.
(a) (b) Figure 8. Leaf chlorophyll content-based damage assessment of rice after toluene exposure: (a) leaf chlorophyll content measured by photometric method; (b) SPAD value measured by chlorophyll meter. The same capital letters refer to the same growth stages. The small letter refers to the significant difference in leaf chlorophyll content among the control group and toluene exposed rice (α = 0.05).

Figure 7.
Recovery of rice after toluene exposure. The left mean NDVI was measured at 5 days after damage at each growth stage, and the right mean NDVI (recovery) was measured at 67 days after planting. The small letters refer to the significant differences in the mean NDVI of 14 retransplanted rice (α = 0.05). The mean NDVI within the same small letter is not significantly different at Tukey's HSD (p > 0.05).

Physiological Characteristic-Based Damage Assessment
The leaf chlorophyll contents after toluene exposure were analyzed to demonstrate the leaf color changes after toluene exposure. The measurements were conducted at 5 DAD simultaneously with NDVI measurements at each growth stage. Figure 8 shows that the leaf chlorophyll decreases after toluene exposure at each growth stage. Similar leaf chlorophyll patterns were observed in chlorophyll contents measured using a microplate reader ( Figure 8a) and a chlorophyll meter (Figure 8b). The SPAD value was highly correlated with chlorophyll content (R = 0.8381; p < 0.0001; n = 23) ( Figure A3), as described in previous studies [22,27]. In the early three growth stages (ET, LT, and SE stages), the leaf chlorophyll contents and SPAD value significantly decreased with respect to the toluene exposure concentration increases. Compared with the two tillering stages (ET and LT stages) and SE stage, the decreases in chlorophyll contents regarding toluene exposure concentration increase were relatively greater at the two tillering stages. In the latter two growth stages (B and F stages), leaf chlorophyll contents did not change significantly after toluene exposure.
Remote Sens. 2021, 13, x FOR PEER REVIEW 12 of 20 Figure 7. Recovery of rice after toluene exposure. The left mean NDVI was measured at 5 days after damage at each growth stage, and the right mean NDVI (recovery) was measured at 67 days after planting. The small letters refer to the significant differences in the mean NDVI of 14 retransplanted rice (α = 0.05). The mean NDVI within the same small letter is not significantly different at Tukey's HSD (p > 0.05).

Physiological Characteristic-Based Damage Assessment
The leaf chlorophyll contents after toluene exposure were analyzed to demonstrate the leaf color changes after toluene exposure. The measurements were conducted at 5 DAD simultaneously with NDVI measurements at each growth stage. Figure 8 shows that the leaf chlorophyll decreases after toluene exposure at each growth stage. Similar leaf chlorophyll patterns were observed in chlorophyll contents measured using a microplate reader ( Figure 8a) and a chlorophyll meter (Figure 8b). The SPAD value was highly correlated with chlorophyll content (R = 0.8381; p < 0.0001; n = 23) ( Figure A3), as described in previous studies [22,27]. In the early three growth stages (ET, LT, and SE stages), the leaf chlorophyll contents and SPAD value significantly decreased with respect to the toluene exposure concentration increases. Compared with the two tillering stages (ET and LT stages) and SE stage, the decreases in chlorophyll contents regarding toluene exposure concentration increase were relatively greater at the two tillering stages. In the latter two growth stages (B and F stages), leaf chlorophyll contents did not change significantly after toluene exposure.  The grain yield and yield components of retransplanted rice were measured at 130 DAP to investigate the long-term effect of toluene damage to rice. The grain yield of rice exposed to the T3 and T4 levels of toluene at the ET stage and the rice exposed to the T4 level of toluene at the LT stage were not measured because they were completely withered. In general, rice grain yield decreased with respect to the increasing toluene concentrations (Figure 9), except for the T1 level of toluene exposure at the LT and SE stages. The yield component changes with respect to toluene exposure damage are summarized in Table 3. The panicle number per hill and 1000-grain weight decreased after toluene exposure at all growth stages. Filled grain percentage decreased after toluene exposure at early growth stages (ET and LT stages), but it was not changed with respect to toluene exposure at later growth stages (SE, B, and F stages). The SNPP decreased after toluene exposure, but a decreasing tendency with increasing toluene exposure concentration was not observed.

Relationship between the NDVI and Physiological Characteristics at Different Growth Stages
The correlation between the NDVI and physiological characteristics was evaluated, and the Pearson correlation coefficients of each correlation are summarized in Table 4. The NDVI of each growth stage measured at 5 DAD showed a positive correlation with the leaf chlorophyll content and SPAD value at the 5% significance level (n = 25; R = 0.550; p = 0.0046 and n = 23; R = 0.652; p = 0.0007, respectively) (Figure 10a,c). The NDVI-chlorophyll content plot is divided into two clusters, early growth stages (ET and LT stages) and later growth stages (SE, B, and F stages). Significantly, the rice in the earlier growth stages is highly clustered in the left of the plot ( Figure 10). The linear relationship between the NDVI and chlorophyll contents showed different correlations concerning the growth stages (Figure 10b,d). Rice exposed to toluene at both earlier and later growth stages showed a positive correlation between the NDVI and chlorophyll content, but only rice in later growth stages showed a significant correlation (R = 0.941; p < 0.0001) (Figure 10b). Similar patterns were observed in the correlation between the NDVI and the SPAD value (R = 0.771; p = 0.0008) (Figure 10d). Figure 11 shows a positive correlation between the NDVI measured at 67 DAP and the grain yield measured at 130 DAP (n = 14; R = 0.645; p = 0.0127). Table 4. Correlation between the NDVI and physiological characteristics. The E growth stage refers to the earlier growth stages, ET and LT stages. The L growth stage refers to the later growth stages; SE, B, and F stages. n, R, and p refer to the number of plots, the Pearson correlation coefficient, and the p-value, respectively (* p < 0.05; ** p < 0.01; *** p < 0.001).   Table 3. Grain yield and yield components of the control group and toluene-exposed rice. For spikelet number per panicle (SNPP), the same capital letters refer to the same growth stage. The small letters refer to the significant differences in mean SNPP (α = 0.05).

Relationship between the NDVI and Physiological Characteristics at Different Growth Stages
The correlation between the NDVI and physiological characteristics was evaluated, and the Pearson correlation coefficients of each correlation are summarized in Table 4. The NDVI of each growth stage measured at 5 DAD showed a positive correlation with the leaf chlorophyll content and SPAD value at the 5% significance level (n = 25; R = 0.550; p = 0.0046 and n = 23; R = 0.652; p = 0.0007, respectively) (Figure 10a,c). The NDVIchlorophyll content plot is divided into two clusters, early growth stages (ET and LT stages) and later growth stages (SE, B, and F stages). Significantly, the rice in the earlier growth stages is highly clustered in the left of the plot (Figure 10). The linear relationship between the NDVI and chlorophyll contents showed different correlations concerning the growth stages (Figure 10b,d). Rice exposed to toluene at both earlier and later growth stages showed a positive correlation between the NDVI and chlorophyll content, but only rice in later growth stages showed a significant correlation (R = 0.941; p < 0.0001) (Figure 10b). Similar patterns were observed in the correlation between the NDVI and the SPAD value (R = 0.771; p = 0.0008) (Figure 10d). Figure 11 shows a positive correlation between the NDVI measured at 67 DAP and the grain yield measured at 130 DAP (n = 14; R = 0.645; p = 0.0127). Table 4. Correlation between the NDVI and physiological characteristics. The E growth stage refers to the earlier growth stages, ET and LT stages. The L growth stage refers to the later growth stages; SE, B, and F stages. n, R, and p refer to the number of plots, the Pearson correlation coefficient, and the p-value, respectively (* p < 0.05; ** p < 0.01; *** p < 0.001).  Figure 10. Correlation between the NDVI and (a) leaf chlorophyll content for all growth stages; (b) leaf chlorophyll contents for the earlier and later growth stages; (c) SPAD value for all growth stages; (d) SPAD value for the earlier and later growth stages. The blue line refers to the 95% confidential interval. R and p refer to the Pearson correlation coefficient and p-value, respectively (* p < 0.05; ** p < 0.01; *** p < 0.001).  Figure 11. Correlation between the NDVI and grain yield. The blue line refers to the 95% confidential interval. R and p refer to the Pearson correlation coefficient and p-value, respectively (* p < 0.05; ** p < 0.01; *** p < 0.001).

Discussion
In this study, we demonstrated the availability of UAV multispectral imagery to crop monitoring and damage assessment. The UAV-based multispectral imagery exhibited suf- R=0.645 p=0.0127* Figure 11. Correlation between the NDVI and grain yield. The blue line refers to the 95% confidential interval. R and p refer to the Pearson correlation coefficient and p-value, respectively (* p < 0.05; ** p < 0.01; *** p < 0.001).

Discussion
In this study, we demonstrated the availability of UAV multispectral imagery to crop monitoring and damage assessment. The UAV-based multispectral imagery exhibited sufficient spatial resolution (4.15 cm to 5.22 cm GSD) which was suitable for detailed crop monitoring.
To select the optimized VI that best reflects the change in rice with respect to the toluene damage, we conducted the sensitivity analysis. There are no strict rules for selecting the VI. Huete (1988) recommended the use of the SAVI, which can reduce soil noise [21]. Zhou et al. (2016) assessed hail damage in potatoes using the GNDVI [15]. Several studies have recommended using the most sensitive VI for vegetation monitoring or damage assessment in terms of the signal-to-noise ratio or change in the VI after the disaster [28,29]. In this study, we observed that the NDVI slope values were negatively higher than other VIs, so we assessed the toluene damage on rice based on the NDVI.
The damage on rice after toluene exposure can be assessed based on NDVI extracted from UAV imagery. The mean NDVI of toluene exposed rice was significantly decreased after severe levels of toluene exposure ( Figure 5). These results reflect the canopy changes or rice after toluene exposure ( Figure 6). The NDVI is closely related to canopy structures, such as LAI, fractional vegetation cover, and chlorophyll contents [27,30]. This study claims that the changes were caused by toluene exposure (i.e., chlorosis, necrosis, and tip burn). The mean NDVI of the rice exposed to T1 level (45.1 ppm) of toluene was not decreased which indicates that the rice was not damaged when exposed to the low level of toluene. However, the mean NDVI of toluene-exposed rice in the F stage was not significantly changed. In the flowering stage, anthers form at the top of the panicle, and leaves turn yellow [16]. These natural changes in rice are difficult to distinguish from toluene exposure damage through the NDVI. Further research is needed to distinguish between natural changes and chemical exposure damage. Based on the recovery assessment of tolueneexposed rice, we demonstrated that the rice exposed to the higher concentration of toluene at the earlier growth stage had a lower NDVI value (Figure 7). These results indicate that toluene exposure at the earlier growth stage is more critical to rice regrowth and recovery.
The damage on rice after toluene exposure was also assessed based on physiological characteristics. The leaf chlorophyll contents and SPAD values were significantly decreased after toluene exposure at early growth stages, however, it was not changed after toluene exposure at later growth stages ( Figure 8). These results indicate that the rice damage from toluene exposure decreases as the rice grows; the toluene exposure damage to rice cannot be quantified through leaf chlorophyll content after the B stage. The grain yield decreased with respect to the increasing toluene exposure levels. The decrease in grain yield might result from a decrease in chlorophyll contents. Since more than 90% of crop biomass is derived from photosynthesis, the decrease in crop yield is closely related to the decrease in photosynthesis due to a decrease in chlorophyll contents [31].
The correlation between NDVI and physiological characteristics was evaluated. The NDVI and leaf chlorophyll contents showed a positive relationship. The NDVI-leaf chlorophyll content plot is clustered in two, early and later growth stages ( Figure 10). The rice in the earlier growth stages has lower NDVI values than the rice in the later growth stages, even if it has similar chlorophyll content because LAI and biomass are smaller. In general, the NDVI has a strong linear or exponential relationship with LAI and biomass and leaf chlorophyll content [32,33]. The rice in the earlier growth stages is also more susceptible to toluene exposure, therefore, the leaf chlorophyll content decrease after toluene exposure was greater than the rice in the later growth stages. Due to the lower NDVI values and small NDVI decreases in the earlier growth stages, the decrease in leaf chlorophyll content after toluene exposure is much greater than the decrease in NDVI. Based on these characteristics of the rice in the earlier growth stages, the plots of two earlier growth stages are highly clustered in the left of the plot. In addition, only rice in later growth stages showed a significant correlation (Figure 10b,d). This is because the NDVI cannot perfectly distinguish vegetation and soil if the NDVI of vegetation is too low in the earlier growth stages. The grain yield also showed a positive correlation between NDVI ( Figure 11). These results were consistent with previous research, demonstrating the positive correlation between VIs and leaf chlorophyll content and crop yield [13,15,30,34]. These significant correlations between the NDVI and physiological characteristics indicate that the NDVI based on UAV multispectral imagery can be alternatives for crop monitoring, such as crop damage assessment and yield prediction for long-term monitoring after chemical exposure. However, considering the imperfect distinguishment of vegetation and soil and the non-significant correlation between NDVI and chlorophyll content in early growth stages, there is a limitation to assessing damage based on NDVI for crops in early growth stages. The damage assessment of rice in earlier growth stages might be better to use the physiological characteristics rather than VIs. The applicability of NDVI to damage assessment has to be further developed in terms of the removal of background noise.

Conclusions
This study assessed toluene damage to rice at four different exposure levels with five growth stages using the NDVI based on UAV multispectral imagery. The rice exposed to higher toluene concentrations and at the earlier growth stage was severely damaged and showed a lower mean NDVI and lower recovery after exposure. Physiological characteristics, such as chlorophyll content and grain yield, showed similar decreasing patterns with decreasing NDVI patterns after exposure of rice to toluene, especially a larger decrease at the earlier growth stage with exposure to increasing concentrations of toluene. Positive correlations were observed between the NDVI and physiological characteristics, indicating that the NDVI based on UAV multispectral imagery can reflect chemical exposure damage to rice. Although this study performed only for rice exposed to toluene using NDVI, which showed the highest sensitivity to toluene damage, there might be better VIs than NDVI for other crop species and other chemicals. Sensitivity analysis can be a useful tool for the selection of optimal VI. Overall, this study demonstrated that the VI based on UAV multispectral imagery is capable of being an alternative tool for crop monitoring of damage assessment by chemical exposure.
Author Contributions: H.K. contributed to the remote sensing data and physiological data acquisition, methodology development, data processing and analysis, and writing the original draft. W.K. contributed to rice cultivation and the toluene exposure experiment. S.D.K. contributed to the overall experiment design, review and editing the manuscript, and funding acquisition. All authors have read and agreed to the published version of the manuscript.  Figure A1. Images of (a) the chamber for toluene exposure simulation; (b) the rice five days after toluene exposure.  Figure A1. Images of (a) the chamber for toluene exposure simulation; (b) the rice five days after toluene exposure.  Figure A2. Field plot description for the retransplanted rice. LT, SE, and B indicate the growth stages late tillering, stem elongation, and booting, respectively. Figure A2. Field plot description for the retransplanted rice. LT, SE, and B indicate the growth stages late tillering, stem elongation, and booting, respectively. Figure A2. Field plot description for the retransplanted rice. LT, SE, and B indicate the growth stages late tillering, stem elongation, and booting, respectively. Figure A3. Correlation between chlorophyll contents measured using a microplate reader and the SPAD value. R and p refer to the Pearson correlation coefficient and p-value, respectively (* p < 0.05; ** p < 0.01; *** p < 0.001).