A remote sensing approach to estimate variable crop coefficient and evapotranspiration for improved water productivity in the Ethiopian highlands

Proper and reasonable estimation of actual evapotranspiration is critical for the design, operation, and management of irrigation systems. However, the density of climatic stations is low in many parts of Ethiopia to estimate spatially reasonable reference evapotranspiration (ETo), and the lack of spatial variability of crop coefficient (Kc) is clear barrier to the proper management of irrigation water. Therefore, the main objective of this study was to estimate reasonable crop evapotranspiration (ETc) by deriving spatially and temporally varying crop coefficients using remote sensing products in 10 wheat plots at the Koga irrigation scheme. The moderate-resolution imaging spectroradiometer (MODIS) potential evapotranspiration was calibrated at two class-I climate stations (Bahir Dar and Dangila) based on Penman–Monteith estimates and Normalized Difference Vegetation Index (NDVI) was derived from Sentinel 2B, which was later used to derive Kc. The correlation between reference evapotranspiration from MODIS Penman Monteith was acceptable for both Dangila (R2 = 0.64) and Bahir Dar (R2 = 0.74) stations. The MODIS calibration constant for Koga irrigation schemes was 0.27 and 1.99 (regression slope and constant) on an 8-day basis. Similarly, a strong correlation (R2 = 0.95) was found between Sentinel-based NDVI and FAO crop coefficient, which indicated an alternative pathway for estimating crop coefficient. The value of Kc varies in space (across the 10 plots) from 0.16 to 0.42 at the initial stage and from 1.20 to 1.32 at the mid-stage. Similarly, the mean value of Kc varies in time from 0.29 at the initial stage to 1.26 at the mid-stage. On the other hand, evaluation of MODIS and WaPOR ETc found a significant difference (p < 0.05) with the calibrated MODIS-Sentinel 2B derived ETc. This indicated the need to calibrate both MODIS and WaPOR for the proper estimation of crop water needs. Underestimation of ET was observed from MODIS, and the reverse is true for WaPOR. Our findings showed that calibrating and integrating MODIS with Sentinel 2B would be a feasible approach to estimate Kc and hence ETc that varies in time and space. This would assist water managers in estimating crop water needs for better productivity in the region.


Introduction
Global consumption of freshwater resources is increasing nowadays, especially in the agricultural sector (Leguen 2019). The increasing population in the world leads to the intensification of irrigation for increased food production resulting in additional pressure on freshwater resources (Kumar et al. 2017). The current trend in irrigation intensification is unsustainable (Sappington and Rifai 2018). Sustainable intensification of irrigated agriculture and proper estimation of crop water need is essential to efficiently utilize water resources (Ewaid et al. 2019). The estimation of crop water demand is dependent on the estimation of evapotranspiration (ET) for irrigation planning, design, and efficient 1 3 168 Page 2 of 15 irrigation water management of irrigation schemes (Maina et al. 2014). Evapotranspiration estimation is essential for various applications including the assessment of irrigation feasibility in arid and semi-arid regions (Allen et al. 2007), drought monitoring (Djaman et al. 2018), and conservation of irrigation water through proper irrigation scheduling (Subedi and Chávez 2015;Koech and Langat 2018). The spatial and temporal values of evapotranspiration estimates would help to determine crop water use which is vital for improving the performance of irrigation schemes (Kharrou et al. 2021).
Several studies were conducted on the performance of irrigation schemes in Ethiopian highlands and found insufficient to address the designed command area (Asres 2016;Abera et al. 2019;Haile et al. 2020;Tesfaye 2020). Koga, the major irrigation scheme in the Ethiopian highlands, was designed to irrigate a maximum command area of 7000 ha (MacDonald 2004). However, the scheme was so far able to irrigate a maximum size of about 5000 ha in 2014 (Abiyu and Alamirew 2015). Tewabe and Dessie (2020) reported that the irrigation ratio of the scheme was around 76.6% in the 2016 and 2017 irrigation seasons. The main reasons were poor irrigation water management at the field level, insufficient reservoir release to the command area, and traditional water application techniques resulting in high water loss throughout the scheme. This demands a proper estimation of crop evapotranspiration in the scheme and hence irrigation scheduling (Kassie and Alemu 2021), which would help to improve the performance of the Koga irrigation scheme and productivity in general.
Several approaches can be used to measure and/or estimate actual evapotranspiration (ETc) including eddy covariance, soil water balance approach, Bowen Ration Energy Balance System, pan evaporation, and weighting lysimeters. However, these methods are expensive, timeconsuming, cannot be easily automated, and provides point-based measurements (Genanu 2017;Zhao et al. 2013Zhao et al. , 2006Enku et al. 2011). To overcome such challenges, remote sensing approaches are practiced worldwide as efficient alternatives (Reyes-Gonzalez et al. 2015). This includes Surface Energy Balance Algorithm for Land, SEBAL, , Two Source Energy Balance, TSEB, (Norman et al. 1995), Simplified Surface Energy Balance Index, S-SEBI, (Roerink et al. 2000), Surface Energy Balance System, SEBS, (Basit et al. 2018), Mapping Evapotranspiration at High Resolution with Internalized Calibration, METRIC, (Allen et al. 2007), and the Simplified Surface Energy Balance, SSEB, (Chen et al. 2016). These models require detailed standard metrological data (Kitsara et al. 2013;Najmaddin et al. 2017;Allam et al. 2021), which can be tricky in region where climate stations are scarce (Najmaddin et al. 2017). One of the most commonly used methods for ETc estimation is the crop coefficient approach which is the product of reference evapotranspiration (ETo) and crop coefficient (Kc) (Yimam et al. 2020(Yimam et al. , 2021(Yimam et al. , 2021. Several have developed equations to quantify reference evapotranspiration from the basic fundamental physical principles of energy balance and mass transfer (Maeda et al. 2011). The Penman-Monteith method is among the standard approach for determining ETo, though it requires detailed meteorological data inputs (Zhang et al. 2018). Remote sensing approaches are nowadays becoming popular in estimating evapotranspiration (ET) at various temporal and spatial scales (Kharrou et al. 2021;Maeda et al. 2011) which is a feasible solution from various perspectives for in-situ estimation (Wanniarachchi and Sarukkalige 2022). Remote sensing-based estimation of ET helps water managers in the irrigation sector to manage both water and crop stress conditions . Of the various remote sensing indices used to estimate ETc, The Normalized Difference Vegetation Index (NDVI) is the most commonly used indices used to estimate ETc (Carpintero et al. 2020;Pôças et al. 2020). In the crop coefficient approach, Kc is typically taken from the Food and Agricultural Organization, FAO (Singh and Irmak 2009), which does not capture the spatial and temporal variability in estimating the crop evapotranspiration. The spatial and temporal variability of Kc can be easily derived from NDVI considering actual crop growth conditions and stages (Reyes-González et al. 2018). Similarly, open-access remote sensing datasets would help in accessing reference evapotranspiration products at varying spatial and temporal resolutions during the estimation of crop evapotranspiration. The moderate-resolution imaging spectroradiometer (MODIS) evapotranspiration product (Mu et al. 2011) and the FAO portal product to monitor water productivity through open-access portal (WaPOR) developed for Africa and the Near East (Blatchford et al. 2020) are among the open sources to estimate ET. However, the applications of those open-source remote-sensing products to estimate ET and help in irrigation scheduling are an unexploited approach in the region that needs refined studies.
Recently, Sentinel and Landsat satellites are among the open-source remote sensing products used for site-specific management of irrigation. Even though the sentinel -2 satellite was launched recently 7 years ago, it has an advantage over the Landsat satellite in terms of better spatial and temporal resolution. According to Morave et al., (2021), atmospheric correction is essential to increase the NDVI values.
Therefore, the main objective of this research was to estimate crop evapotranspiration using the combined application of open-access remote sensing ETo and NDVI-derived Kc which was generated from atmospherically calibrated Sentinel-2 images for the major irrigated crop in the Koga irrigation scheme (wheat).
We also did a comparison of remote sensing-based ET for their applications in irrigated water management in the same irrigation scheme. In Ethiopia, wheat grain is the flour for the preparation of traditional as well as modern processed food products and industrially processed products (Aklilu et al. 2015). Wheat in Ethiopia is a significant stable and cash crop in increasing the income of the households, food security, and national Gross Domestic Product (GDP) increment (Anteneh and Asrat 2020). At present, the Ethiopian government has started planning and working to produce better wheat production using modern irrigation technology to reduce the import of wheat from foreign countries and ensure food security. Ukraine and Russia are the most important producer and exporters of wheat in the world (Fellmann et al. 2014). The start of the Russian-Ukraine war has significantly undermined grain exports due to raising the global price and the lack of access to Black Sea ports. As a direct importer from those countries, Ethiopia is the most affected country by the wheat shortage (WFP 2022). The results from this study would help water managers in proper irrigation scheduling to monitor irrigation application trends at the farm level in a Koga irrigation scheme. In general, the output of this research will be used, to improve water use efficiency by the farmers and for decision-makers to plan and design irrigation infrastructures in the region. It is also vital to support policy decisions on water allocation.. Besides, the study can provide a feasible alternative way to obtain temporally and spatially variable evapotranspiration estimates for hydrological analysis and other water-related studies in the Ethiopian highlands.

Study area description
This study was conducted in Koga's large-scale irrigation scheme located in the west Gojam Zone of Amhara Regional State in the Ethiopian highland (11.35 °N latitude and 37.14 °E longitude, at an altitude of 1900 m above sea level) (Fig. 1). Its catchment area is about 220 km 2 including a reservoir area of 17 km 2 with an impounding capacity of 83.1 million cubic meters (MCM). The scheme is planned to irrigate a maximum area of 7000 ha for 10,000 smallholder households (Asres 2016). The irrigation scheme networks include 19.7 km of primary canals, 42.3 km of secondary canals, 117 km of tertiary canals, 783 km of quaternary canals, and a 12-night storage structure. Agriculture in Koga includes rain-fed in rainy periods and irrigation in dry monsoon phases where maize and wheat are, respectively, the most common cereal crops. Irrigated wheat covers about 70% of irrigated command area annually while potato and onion are the second and third major irrigated crops in the Koga irrigation scheme.

Data sources and description
This research was conducted using open-source remote sensing products and meteorological data (Table 1). Sentinel images with 10 m spatial resolution were obtained every 10 days to develop the NDVI map. Similarly, MODIS potential evapotranspiration (PET) products were obtained to estimate reference evapotranspiration. Also, actual evapotranspiration was taken from the MODIS sensor and WaPOR database to extract wheat crop evapotranspiration. The MODIS and WaPOR actual evapotranspiration (ETa) was obtained from late December 2020 to mid of April 2021. A total of 12 ETa image products were downloaded from each website considering the same temporal resolution as Sentinel-2B images. The level-3 dataset of WaPOR actual evapotranspiration for the Koga irrigation scheme was obtained at 30 m spatial resolution every 10 days. In WaPOR, the pixel value represents the average daily ETa in a given 10 days (FAO 2019). Meteorological data (precipitation, mm; relative humidity, %; sunshine hour, hr; and wind speed, m/sec) were collected from two meteorological

Overall framework
Wheat crop evapotranspiration (ETc) was estimated for the Koga irrigation scheme by combining open-source remote sensing products and meteorological data (Fig. 2). The daily meteorological data were collected from nearby stations of Dangila and Bahir Dar from the year 2016 to 2021. The Penman-Monteith method (Eq. 1) was used to estimate ETo at the stations and the result was compared with the MODIS PET product. Then, crop evapotranspiration for the Koga irrigation scheme was computed by multiplying the reference evapotranspiration from MODIS with the crop coefficient in the Geographic Information System, GIS, and spatial analysis. The crop coefficient was derived by developing a regression equation between the NDVI-derived Kc from sentinel 2B and FAO Kc for the wheat crop. Reconnaissance surveys were conducted to identify wheat field plots along the command area during the first stage of sawing. A total of ten wheat field plots were selected from the entire Koga irrigation scheme. Depending on the size and shape of each field, point coordinates were taken at representative points per plot. Atmospherically corrected sentinel-2B Level-1C images were used to generate the NDVI value over the Koga irrigation scheme. In addition, a correlation analysis was performed between MODIS-based crop evapotranspiration with WaPOR and MODIS crop evapotranspiration products to validate uses for the metrological data-scarce region.
where ET o is reference evapotranspiration (mm/day), R n is net radiation at the crop surface (MJ m −2 /day), G is soil heat flux density (MJ m −2 /day), T is mean daily air temperature at 2 m height (°C), u 2 is the wind speed at 2 m height (ms −1 ), e s is saturation vapor pressure (kpa), e a is actual vapor pressure (kpa), e s − e a is saturation vapor pressure deficit (kpa), ∆ is Slope vapor pressure curve (kpa °C −1 ), and is the psychometric constant (kpa °C −1 ).

Comparison of MODIS PET product and Penman-Monteith ETo
The moderate-resolution imaging spectroradiometer (MODIS)-based evapotranspiration (MOD16A2)  The MOD16 global evapotranspiration datasets have a spatial resolution of 500 m and temporal resolution of 8-day, monthly and annual intervals. The 8-day potential evapotranspiration (PET) is the sum of potential evapotranspiration during these 8-day periods. The MOD16 algorithm is based on the logic of the Penman-Monteith approach. The MODIS PET was obtained from January first, 2016 to April 2021 at 500 m spatial resolution every 8 days, and a total of 244 PET MODIS image products were obtained. (1) The true color of Sentinel 2B image: A Top-of-atmosphere reflectance and B bottomof-atmosphere reflectance The remote sensing potential evapotranspiration product was calibrated using the Penman-Monteith reference evapotranspiration estimates. The MODIS potential evapotranspiration image products were extracted by points to the exact meteorological stations (i.e., Dangila and Bair Dar). A regression equation was developed between the Penman-Monteith evapotranspiration and MODIS potential evapotranspiration for each meteorological station. Based on the goodness of fit, a recommendation was made to use a remote sensing tool for evapotranspiration estimation, particularly in data-scarce regions. where ETo c and ET m are calibrated reference evapotranspiration and 8-day MODIS potential evapotranspiration, respectively. The calibration constants a and b are average regression constants (Table 1) computed from the regression equation between the remote sensing evapotranspiration and evapotranspiration from the meteorological station at Dangila and Bahir Dar stations.

Computation of satellite-derived crop coefficient
The information obtained from the satellite images needs to undergo various corrections (i.e., atmospheric, geometric, and radiometric) to minimize distortion and enhance quality data (Young et al. 2017). The Sentinel-2 Level-1C products had geometrically registered and radio-metrically calibrated (Phiri et al. 2020), and thus only atmospheric correction was done. The atmospheric correction for the bottom of atmosphere reflectance for Sentinel-2 images was obtained from the European space agency (ESA) Sentinel Application Platform (SNAP). The surface reflectance values were obtained for red and infrared bands to calculate NDVI.
The Sentinel images were obtained at the different growth stages of the wheat crop. A total of 12 images were obtained (2) ETo C = a * ET m + b where NIR and RED are near-infrared (Band 8) and red band (Band 4) of the sentinel image, respectively.

Computation of crop evapotranspiration and comparisons with remotely sensed ETc
Crop evapotranspiration (ETc) was computed and mapped as a product of the NDVI-derived crop coefficient and calibrated ETo from 8-day MODIS PET. The 8-day calibrated MODIS reference evapotranspiration map was resampled to 10 by 10 spatial resolution to fit with the Kc map. The crop evapotranspiration map was extracted by field plot polygon to determine wheat water requirement for each growing season. A correlation analysis was done between crop evapotranspiration obtained from the Penman-Monteith method and MODIS for validation. Similarly, the WaPOR crop evapotranspiration and MODIS ETc product were obtained at 30 m and 500 m resolution, respectively, and resampled to 10 m spatial resolution. The MODIS and WaPOR crop evapotranspiration product was compared with the ones derived from the MODIS PET. A two-tailed student t-test was used to examine the mean difference between calibrated ETc derived from MODIS and Sentinel 2B, and the direct

Atmospheric correction
The reflectivity difference (Fig. 3) between the true color of the top of the atmosphere and the bottom of the atmosphere radiance presents the big difference between atmospherically calibrated and uncorrected images. To determine the effects of atmospheric correction on the raw data, changes in Sentinel 2B red and NIR bands for NDVI calculation were examined. The visual comparison of Fig. 1A and B showed that atmospheric correction is essential for Sentinel-2B image quality enhancement. David Morave et al., (Moravec et al. 2021) revealed that atmospheric correction is used to increase the NDVI values and reported an increased NDVI value of 0.497 from an atmospherically corrected image than an uncorrected image (0.421).

Evaluation of MODIS reference evapotranspiration
The MODIS potential evapotranspiration values were evaluated with the Penman-Monteith reference evapotranspiration estimated for Bahir Dar and Dangila climatic stations. The root means square error and means absolute error between remote sensing evapotranspiration product and the Penman-Monteith method for both climatic stations were within the acceptable limit, 0.2 to 0.5, based on Yu et al., (2016) data on an 8-day basis (Table 2 and Fig. 4). Similarly, the coefficient of determination (R 2 ) showed a good correlation, more than 50%, based on Taylor (1990) for both stations. The result would help to estimate the temporal and spatial variation of evapotranspiration at a regional and local scale in areas where climate stations are scarce or unavailable. The calibration parameters in Table 2 were obtained from the regression equations of each of the trend lines in Fig. 4. The average of the 8-day regression constants (i.e., 0.27 for a and 1.99 for b) was taken to map the calibrated MODIS reference evapotranspiration for further estimation of crop evapotranspiration in the Koga irrigation scheme. The monthly time series of the calibrated reference evapotranspiration (Fig. 5) showed an increasing trend in the dry period and reached a maximum in April for both stations with 4.7 mm and 4.5 mm for Bahir Dar and Dangila, respectively. Whereas, the minimum value was around 2.6 mm for both stations and occurs in July. This is because the region experiences an annual maximum temperature in April resulting in the high value of ET. The reverse is true in July, the air temperature is low resulting in a reduced ET rate. The ranges of evapotranspiration for Bahir Dar and Dangila were 2.6 mm to 4.7 mm, and 2.6 mm to 4.5 mm, respectively. The mean daily reference evapotranspiration for Bahir Dar and Dangila, respectively, were found 3.9 and 3.6 mm. Mengistu and Amente (2020) found 3.67 mm and 3.59 mm of mean daily reference evapotranspiration for Bahir Dar and Dangila stations, respectively, which supports our findings.

Deriving NDVI for the various growth stages of wheat
The NDVI values for ten wheat field plots along the crop growth stage are presented in Fig. 6. At the initial stage, the value of NDVI is low because the ground cover is small at this stage. More precisely, at the initial stage, the reflectance values of the near-infrared (NIR) wave band are low, whereas the red wave band has high reflectance at the initial stage (Panda et al. 2010). In addition, the red wave band is strongly absorbed by chlorophyll whereas the NIR wave band is highly reflected from the mesophyll structure of the crop leaves. According to Yin et al. (2012) the high values of NDVI are related to dense vegetation with high  reflectivity values of NIR band. Maximum values of NDVI were obtained at the mid-stage (48-88 DAP) in each wheat field plot which ranges from 0.61 to 0.83. Whereas, the minimum values rages from 0.13 to 0.45 for the initial and final stages of the crop cycle. The lowest average value of the NDVI for the ten plots was found at 0.21 for the initial stage, the 8th day after plantation in late December (Table 3). The average NDVI values gradually increased, proportionally to their growth stage. On the 78th day after the plantation, the maximum average NDVI was found 0.8 (in early March). A relatively high range of average NDVI was observed at the mid-stage, early February to mid-March, (day 50 to day 88 after plantation). At the end of the growing stage, the average NDVI values gradually decreased to 0.27 in the middle of April (118th days after planting). The value of NDVI across all plots ranged from 0.13 to 0.83 (Table 3). Our estimate of NDVI was found consistent with similar research findings with a range of maximum NDVI values of 0.67 to 0.9 (Hunsaker et al. 2005;Farg et al. 2012;Lei and Yang 2014;. The spatial and temporal variability of the NDVI in the study area is presented in Fig. 10.

Deriving crop coefficient for the various growth stages of wheat
The correlation between the average NDVI and FAO crop coefficient of wheat along the crop growth stage was found strong, R 2 = 0.95, (Fig. 7). This depicted crop coefficient could be determined from remote sensing indices (NDVI). The regression equation was used to develop a crop coefficient, and the resulting crop coefficient for the wheat crop is presented in Table 4. The average crop coefficient was found about 0.83 for the crop base period. The minimum crop coefficient was found 0.16 during the initial period (on the 8th day after plantation) and the maximum value was 1.32 at mid-stage (on the 78th day after plantation). In the early stage, the plant height and surface area of crop leaves are small resulting in a lower Kc value and the reverse is true as the cop grows. As the crop maturity increases, the greenness of the crop reduced which resulted in a reduction in crop coefficient. Several research findings reported a closer estimate of the crop coefficient for the wheat crop from NDVI (Gontia and Tiwari 2010;Er-Raki et al. 2013;Adamala et al. 2013;). The spatial and temporal variability of the crop coefficient is presented in Fig. 11. The FAO-Kc values were taken from Food Agricultural Organization Irrigation and Drainage Paper-56 for the wheat crop, while calibrated Kc values were derived from NDVI values of Sentinel 2B at the various growth stages (Fig. 8). The result illustrated that the Kc values estimated from NDVI were higher than FAO Kc at an early stage (DAP 1-18) and mid-stage (DAP 58-88). From mid-development

Crop evapotranspiration
The 8-day MODIS PET and derived Kc values from Sentinel 2B were used to generate the crop evapotranspiration (ETc) in the study area with 10 m spatial resolution (Fig. 12). The minimum crop evapotranspiration was 0.59 mm/day during the initial period, while the maximum crop evapotranspiration was 5.66 mm/day during the midstage (Table 5). The average crop evapotranspiration of wheat for the crop base period was found to be 3.2 mm/ day. This indicated that a small amount of irrigation water was needed at the initial and late stages, whereas, the highest crop water need was observed at the mid-stage. Crop evapotranspiration at middle stage of the crop was larger than other growing stages, which agreed with (Reyes- Gonzalez et al. 2015;Gontia and Tiwari 2010;Er-Raki et al. 2013;Farg et al. 2012;. This study argued that the estimation of crop evapotranspiration using remotely sensed data is very important for designers to better understand the amount of irrigation for various crop growth stages and for better irrigation scheduling. In Fig. 12

Comparison between MODIS-sentinel 2B derived ETc, WaPOR ETc, and MODIS ETc
The comparison between the MODIS-Sentinel 2B derived wheat crop evapotranspiration, direct MODIS ETc observation, and WaPOR is presented in Table 7 and Fig. 9. The disparity in variance between data sets was checked and the result showed equal variance. A two-tailed, t-test assuming equal variances with a 5% significance level was used to examine the significant difference between the MODISbased calibrated crop evapotranspiration and the direct remotely sensed crop evapotranspiration (WaPOR ETc and MODIS ETc) of the wheat crop (  (Fig. 9). On the other hand, the advantage of combing MODIS and Sentinel 2B is its higher spatial resolution (10 m) capable of capturing even small plots of fields. The point scale daily reference evapotranspiration estimates of the remote sensing products are presented in Table 7. The mean daily point reference evapotranspiration for Bahir Dar and Dangila, respectively, was 3.9 and 3.6 mm ( Table 7). The range of evapotranspiration for Bahir Dar and Dangila was from 2.6 to 4.7 mm, and 2.6 to 4.5 mm, respectively. The result showed a closer estimate of reference evapotranspiration for both stations with a slight increment in Bahir Dar station. The table caption presented here is as presented in the revised manuscript.

Conclusion
This study aimed to develop an NDVI-based crop coefficient and examine crop evapotranspiration of wheat that varies in time and space using MODIS and Sentinel 2B products for the Koga irrigation scheme. We found a good correlation (R 2 > 0.6) between MODIS potential evapotranspiration and Penman Monteith's estimate of reference evapotranspiration and developed a regression coefficient that serves for the calibration of MODIS. Similarly, a strong correlation (R 2 = 0.95) between NDVI derived from Sentinel 2B and FAO Kc values, and regression coefficients were used to generate time and space variable Kc from NDVI. The mean value of Kc from the 10 plots significantly varied in time; from 0.29 at the initial stage to 1.26 at the mid-stage. Similarly, a significant Kc variation was found across the 10 plots; 0.16 to 0.42 at the initial stage and 1.20 to 1.32 at the mid-stage. On the other hand, the evaluation of MODIS and WaPOR ETc indicated that both need to be calibrated before being directly used for the estimation of crop water needs. The combination of MODIS and Sentinel 2B to derive time and space variable Kc and hence ETc was found a feasible alternative for proper water management and improvement of crop productivity.
Author contributions DWM contributed to the experimental design, data acquisition, analysis, and interpretation, and drafted the manuscript; TTA contributed to the experimental design, data acquisition, and analysis, and revised the manuscript for the scientific content; AYY contributed to data analysis and interpretation, and revised the manuscript. All authors have read and approved the final version of the manuscript; SAB contributed to data analysis and interpretation and revised the manuscript.
Funding This research was funded by Kansas State University -Agricultural Dean's Office, Contract No: 41339. The opinions expressed herein are those of the author(s) and do not necessarily reflect the views of Kansas State University or any other institution.

Conflict of interest The authors declare no conflict of interest.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.