Satellite-Based Determination of the Water Footprint of Carrots and Onions Grown in the Arid Climate of Saudi Arabia

Abstract

Increasing demand for food, climate change, and other human interventions are leading to significant increases in water consumption by the agricultural sector. This requires rationalizing the water used for the production of agricultural crops through improved irrigation management practices. Therefore, this study aimed to estimate the water footprint (WF) of onion (Allium cepa L.) and carrot (Daucus carota) crops using the CROPWAT model and the SSEB (Simplified Surface Energy Balance) algorithm. Experiments were carried out at two center-pivot irrigated fields belonging to Tawdeehiya Commercial Farms in the southeastern region of the Riyadh governorate, Saudi Arabia. Individual bands and vegetation indices (VIs) were retrieved from Sentinel-2 satellite data, including the normalized difference vegetation index (NDVI), soil adjusted vegetation index (SAVI), optimized soil adjusted vegetation index (OSAVI), renormalized difference vegetation index (RDVI), and enhanced vegetation index (EVI), and the land surface temperatures (LST) extracted from Landsat-8 data were used to estimate crop productivity (CP), crop water use (CWU) (i.e., evapotranspiration—ET_a), and crop WF. Crop growth/phenology stages and georeferenced biophysical parameters were recorded during the growth period, and crop yield samples were collected randomly from predetermined sampling locations. It was found that the NIR band was appropriate for predicting onion yield (R² = 0.68; p > F = 0.02) and carrot yield (R² = 0.77; p > F = 0.02). The results also showed the feasibility of using the RDVI and EVI to estimate the yields of onion and carrot crops, with bias values of 15% and –17%, respectively. The CWU has also been successfully estimated using the SSEB algorithm, with an overall accuracy of 89%. The SSEB-estimated CWU was relatively high compared to the applied amounts by 10.6% (onions) and 12.6% (carrots). Finally, the crop WF was successfully estimated at 312 m³ t⁻¹ and 230 m³ t⁻¹ for carrots and onions, respectively, with an overall accuracy of 71.11%. The outcomes of this study can serve as a reference for crop irrigation management practices in the study region and areas with similar environmental conditions.

1. Introduction

The demand for water has increased dramatically around the world, with particular attention being paid to the water used for the irrigation of agricultural crops and related industries. Therefore, the sustainable use of water in agricultural production, particularly in arid regions such as Saudi Arabia, has become one of the key priorities of modern agricultural strategies. Due to the drastic shortage of freshwater resources and harsh environmental conditions, Saudi Arabia is making great efforts to develop sustainable agriculture. In this regard, improved agricultural monitoring and water footprint (WF) estimation techniques provide important information on the areas under different crops, crop conditions, yield forecasts, and the amounts of water used. These data are very important for optimizing decision-making strategies related to agricultural production and food security measures [1].

The Water Footprint Analysis (WFA), described by Hoekstra and Hung [2], is a common tool for estimating the freshwater consumed by different products, and has gained increasing applicability in estimating freshwater consumption by agricultural crops. The WF, however, refers to the amount of water applied to produce a unit of the produce [1,3] and is commonly estimated as m³ kg⁻¹, m³ t⁻¹, etc. The crop WF has three main components, namely, the green WF, the blue WF, and the grey WF. The green (irrigation) and blue (rain) WFs represent the amount of water used by plants (evapotranspiration, ET_a), while the grey WF represents the amount of water used to absorb pollutants from agrochemicals or water quality levels. What distinguishes the WF method is that it focuses primarily on the actual water consumed by plants, instead of the total amount of applied water [4,5]. Therefore, application of the WFA tool to estimate the amount of freshwater used by crops has become essential for the effective management of irrigation water [6,7,8,9].

Evapotranspiration can be either measured in situ by specialized devices or indirectly by mathematical models using climatic, soil, and crop data. Advances in remote sensing and GIS technologies have contributed to the improvement of agricultural management systems by integrating weather, soil, and crop data [10,11,12]. Satellite remote sensing allows for estimating crop water use (CWU, i.e., ET_a) by means of various algorithms developed based on the surface energy balance components, such as the simplified surface energy balance (SSEB), surface energy balance algorithm over land (SEBAL), surface energy balance system (SEBS), mapping evapotranspiration with internal calibration (METRIC), two-source model (TSM), two-source time-integrated model (TSTIM), and simple algorithm for evapotranspiration retrieval (SAFER). However, the SSEB algorithm is characterized by its simplicity, minimal data requirements, ease of implementation, and acceptable accuracy, and it is therefore widely used to estimate ET_a for agricultural fields [13,14,15].

Several previous crop WF studies have focused on either global [16,17] or local [18,19] scales. However, most crop WF estimations performed at regional, local, or global scales, based on crop, soil, and climate data, return crop WF results of varying accuracy [20,21]. Furthermore, the use of remote sensing data for the estimation of crop water requirement, irrigation scheduling, and water footprint (WF) mapping is still limited in the Middle East. Hence, a knowledge gap in the area of WF still exists in this region, especially in Saudi Arabia. Therefore, this research was undertaken to address the spatial assessment of the WF of carrot and onion crops using free satellite images, such as Landsat-8 and Sentinel-2 images, along with site-specific data, such as crop phenology, yield, climate, and irrigation data. However, the objectives covered in this research include: (1) quantification of the WF of carrot and onion crops according to local climates, and (2) a comparison between the crop WF determined using the satellite-based SSEB and that obtained from field data combined with the CROPWAT algorithm.

The WF of carrot and onion crops was estimated using the CROPWAT and remote sensing approaches. The crop water requirement was calculated using the CROPWAT approach as defined by Allen et al. [22]. The remote sensing estimation of the crop water use was accomplished by assessing the fraction of ET and then the actual ET (ET_a) by applying surface energy balance principles to the Landsat-8 TIRS bands and Sentinel-2 visible bands. However, the lack of rainfall and its irregularity in the study area led to an absence of Landsat-8 images coinciding with rainy days. Hence, this study mainly coped with the blue component of the WF.

2. Materials and Methods

2.1. Experimental Area

The research study was carried out on two center-pivot-irrigated fields belonging to Tawdeehiya Commercial Farms in the southeastern region of the Riyadh governorate, Saudi Arabia, between latitudes 24.172992° and 24.210347°N and longitudes 47.937389° and 48.085711°E. As shown in Figure 1, a 27-ha field (ID: N6) and a 2-ha field (ID: T/PVT) were grown with carrot and onion crops, respectively. The study area is characterized by an arid climate, where temperatures range from 11 °C in the winter to 49 °C in the summer, and there is an average rainfall of 98 mm in the period from November to early March. The texture of experimental pivots is characterized as sandy loam soil. The main crops cultivated in this establishment include onions, carrots, and other vegetable crops.

2.2. Experimental Details

The experimental work was carried out in the period between February and July 2020. The onion crop (cultivar: red onion) was planted on 26 February 2020, on soil beds of 1.4 m width, 30 cm row spacing, and 10 cm between plants in a row. The carrot crop (cultivar: scarlet red) was sown on 22 March 2020, on soil beds of 2.5 m width, 30 cm row spacing, and a 15 cm between plants in a row. Irrigation water was applied to the experimental fields using a center-pivot irrigation system. Harvesting of onion and carrot crops took place on 26 April 2020 and 12 July 2020, respectively.

2.3. Field Sampling and Data Collection

A GeoXH6000 (Trimble, Westminster, CO, USA) handheld GPS receiver was used to georeference 120 randomly selected sampling locations from the onion field (30 sampling locations) and carrot field (90 sampling locations). Yield samples were collected from the same geo-referenced locations by weighing the harvest collected from an area of 2.5 m² at each sampling site. After removing damaged/irregularly shaped carrots or onions, the yield samples were weighed and presented in a common harvest unit (t ha⁻¹).

A manual soil auger with a bit diameter of 62 mm and a steel rod length of 1.0 m was used for the collection of soil samples at a depth of 0 to 10 cm from the soil surface. The collected samples were characterized for bulk density (g m⁻³), field capacity (%), and wilting point (%), adopting the standard laboratory methods. The Kc values for the studied crops were determined as described by Allen et al. [22].

2.4. Water Footprint Assessment

The overall methodological flow is given in Figure 2 and described in subsequent subsections. The blue and green water use, in terms of depth (mm), was computed based on the crop water requirement (CWR) with the use of CROPWAT (Land and Water, Food and Agriculture Organization, Hot Springs, VA, USA) (Ver. 8.0, FAO), as described by Allen et al. [22].

2.4.1. CROPWAT Input Data

CROPWAT is considered an important decision-support model in the field of agricultural production, created by the Food and Agriculture Organization (FAO) through the Land and Water Development Division, used to calculate the water required to irrigate different crops according to crop parameters, soil characteristics, and climatic conditions. CROPWAT requirements, such as mean values of precipitation, wind speed, minimum and maximum temperatures, sunshine hours, and relative humidity were collected from the records of the weather station facility of the experimental farm (Figure 2). The crop coefficients, used to determine the actual evapotranspiration (ET_c), vary based on local conditions, the crop planting date, crop development stages, and the length of the growing season [23]. This approach is used to calculate the crop evapotranspiration (ET_c) by multiplying the reference evapotranspiration (ET_o), which reflects the atmospheric water demand based on the meteorological conditions, by the crop coefficient (K_c), which integrates crop water requirements and farmland characteristics. The crop coefficients (K_c) have been estimated for the majority of commercial crops, and are available in tables with global values covering the different growth stages of the crops [22]. Despite the fact that the global Kc values are valid only for unstressed crops grown under optimum agricultural conditions, these are the values most often used when running the CROPWAT model (Table 1). The crop water footprint (WF) is then estimated according to Equations (1) and (2) [3].

$${\mathrm{ET}}_{\mathrm{c}}=\mathrm{CWU}$$

(1)

$$\mathrm{WF}=\frac{\mathrm{CWU}}{\mathrm{CY}}$$

(2)

where:

WF = crop water footprint (m³ t⁻¹),

CWU = crop water use (m³ ha⁻¹),

CY = crop yield (t ha⁻¹).

Table 1. Soil and crop data used for the execution of the CROPWAT model.

Parameter		Carrots	Onions
Root depth (cm)	Min.	18.00	15.00
	Max.	35.00	25.00
Crop coefficients (length of the growth period in days)	Kc—initial	0.35 (30)	0.40 (15)
	Kc—mid	1.15 (50)	1.05 (30)
	Kc—end	0.65 (60)	0.60 (15)
Soil characteristics	Bulk density (g cm−3)	1.27	1.29
	Field capacity	0.34	0.34
	Wilting point	0.25	0.25

Table 1. Soil and crop data used for the execution of the CROPWAT model.

Parameter		Carrots	Onions
Root depth (cm)	Min.	18.00	15.00
	Max.	35.00	25.00
Crop coefficients (length of the growth period in days)	Kc—initial	0.35 (30)	0.40 (15)
	Kc—mid	1.15 (50)	1.05 (30)
	Kc—end	0.65 (60)	0.60 (15)
Soil characteristics	Bulk density (g cm−3)	1.27	1.29
	Field capacity	0.34	0.34
	Wilting point	0.25	0.25

2.4.2. Satellite Data and Image Analysis

Sentinel-2 (S2) and Landsat-8 satellite images for the study period (January–July 2020) were downloaded (

Table 2. Summary of satellite and source data.

Satellite	Particulars	Source
Sentinel-2	Level 2A (MSI) BOA	https://scihub.copernicus.eu (accessed on 30 May 2022)
Landsat-8	Level 2A (OLI, TIRS)	https://earthexplorer.usgs.gov (accessed on 28 August 2022)
Weather data	Historical monthly weather data (2000–2020 *)	https://Worldclim.org (accessed on 12 September 2022)

* 2019 and 2020 are the projected and downscaled datasets.

), and the image analysis and interpretation were performed with the use of the SentiNel Application Platform (SNAP) developed by the European Space Agency (ESA) and the ESRI GIS programs. Sentinel-2 multispectral pictures were utilized to generate land use/land cover (LU/LC) maps, vegetation indices (VI), albedo, and other data. Landsat-8 images were utilized to estimate the land surface temperature (LST) of the experimental fields. The extracted region of interest (ROI), i.e., the area corresponding to the experimental fields, along with all the datasets, was managed in the Universal Transverse Mercator (UTM) map projection, WGS84 datum, North zone-38. Subsequently, all the generated datasets were utilized as input for developing crop productivity (CP), evapotranspiration factor (ET_f), crop water consumption (CWU), and crop water footprint (WF).

2.5. Estimation of Crop Productivity

Individual bands (blue, green, red, and NIR) and vegetation indices (Vis) (e.g., normalized difference vegetation index—NDVI, soil adjusted vegetation index—SAVI, optimized soil adjusted vegetation index—OSAVI, reorganized difference vegetation index—RDVI, enhanced vegetation index—EVI, and the simple ratio vegetation index—SRVI) were used to predict the crop productivity (CP). Linear regression analysis was performed for the CP prediction models. The collected yield samples were divided into training (60%) and validation (40%) datasets. Dependent variables (i.e., individual bands and the generated VIs) were regressed against crop yield as an independent variable as a training set for prediction models. The accuracy of the obtained linear models was evaluated, and then the most accurate models were used for mapping the CP. The remaining set of independent variables (40% of yield samples) was used to validate the developed models. Then, the model accuracy was assessed using standard statistical indicators: the coefficient of determination (R²), mean bias error (MBE), and root mean square error (RMSE). On the other hand, potential yield distribution maps and the corresponding error maps were also developed and examined.

2.6. Estimation of the Crop Water Use

The green (CWU_G) and blue (CWU_B) portions of the crop water use (CWU) were calculated by accumulating the daily recorded evapotranspiration (ET, mm d⁻¹) for the entire growth period, following the procedure described by Allen et al. [23]. The evaporative demand of the atmosphere (ET_o) was determined using the standard Penman–Monteith method revised and recommended by the FAO. The CROPWAT model, which was used to calculate evapotranspiration, offers two different options for calculating evapotranspiration, i.e., the option of crop water requirements considering optimum constraints, and the option of irrigation practices, such as determining in real time the actual irrigation supply. For a remote sensing quantification of ET, the SSEB technique reported by Senay et al. [14] was used to estimate the CWU as actual evapotranspiration (ET_a). The ET_a was obtained in two steps, namely, the estimation of reference ET fraction (ET_f) and the reference ET (ET_o) as shown in Equations (3)–(5).

$$\mathrm{CWU}=10\times {\sum }_{\mathrm{d}=1}^{\mathrm{lgp}}{\mathrm{ET}}_{\mathrm{a}}$$

(3)

$${\mathrm{ET}}_{\mathrm{a}}={\mathrm{ET}}_{\mathrm{f}}\times {\mathsf{\alpha }\mathrm{ET}}_0,$$

(4)

where α is the scale element, usually 1.2. The reference evapotranspiration (ET_o) was determined after the FAO–Penman–Monteith method (Equation (5)), according to Allen et al. [23]:

$${\mathrm{ET}}_{\mathrm{o}}=\frac{0.408\Delta \left({\mathrm{R}}_{\mathrm{n}}-\mathrm{G}\right)+\mathsf{\gamma }\frac{900}{\mathrm{T}+273}{\mathrm{u}}_2\left({\mathrm{e}}_{\mathrm{s}}-{\mathrm{e}}_{\mathrm{a}}\right)}{\Delta +\mathsf{\gamma }\left(1+0.34{\mathrm{u}}_2\right)},$$

(5)

where:

ET_o = reference evapotranspiration (mm day⁻¹),

R_n = net radiation at the crop surface (MJ m⁻² day⁻¹),

G = soil heat flux density (MJ m⁻² day⁻¹),

T = air temperature at 2 m height (°C),

u₂ = wind speed at 2 m height (m s⁻¹),

e_s = saturation vapor pressure (kPa),

e_a = actual vapor pressure (kPa),

e_s = e_a saturation vapor pressure deficit (kPa),

Δ = slope vapor pressure curve (kPa °C⁻¹),

γ = psychrometric constant (kPa °C⁻¹)

ET_f is the key variable in the SSEB approach because it takes into account the effect of soil moisture on ET_a. ET_o determines a potential ET under unconstrained watering conditions. ET_f was computed using temperature datasets (LST and air), considering that hot pixels (T_h) experience little or no ET [24,25]. Cold pixels (Tc) represent the highest ET. Both S2 and L8 images were utilized to map the ET_f. In this study, however, vegetation coverage was calculated using the S2 data, and the LST was calculated using the TIRS bands of the L8 data.

Cloud-free Landsat-8 OLI/TIRS images were subjected to radiometric correction (linear contrast stretching) to reduce interference errors. The processed images were then used to calculate the crop evapotranspiration (ET_c) through several computational stages including the brightness temperature (T_b), surface temperature (T_s), net radiation (R_n), geothermal flux (G), air heat flux (H), latent heat flux (LE), and finally the evapotranspiration using a digital image processing module (i.e., the SSEB).

The surface temperature (T_s) of each pixel was examined, and the determined hot and cold pixels were used for the calculation of the fraction of ET. Hot pixels were selected, using the NDVI map as a guide, by locating dry bare land (or sparse vegetation) with very low NDVI values. Similarly, cold pixels were selected from areas of high moisture, healthy, completely covered with vegetation, and with maximum values of NDVI. The ET fraction (ET_f,x) of an individual pixel (x) was also calculated using Equation (6):

$${\mathrm{ET}}_{\mathrm{f},\mathrm{x}}=\frac{{\mathrm{dT}}_{\mathrm{h}}-{\mathrm{dT}}_{\mathrm{s}}}{{\mathrm{dT}}_{\mathrm{h}}-{\mathrm{dT}}_{\mathrm{c}}},$$

(6)

where:

dT_h = the temperature difference between the L8 estimated surface temperature (T_s) and the hot pixel air temperature (T_a).

dT_c = the temperature difference between L8 estimated T_s and cold pixel T_a.

dT_x = the temperature difference between T_s and T_a.

Ground temperatures of the six selected pixels (three hot and three cold) were estimated using the ArcGIS software, ver. 10.7.1 (ESRI, Redlands, CA, USA). The database file of the 9.0 software results was exported to an Excel (Microsoft, Redmond, WA, USA) spreadsheet, and then, the hot and cold pixels were averaged. Images containing the ET_f of individual pixels were then utilized to calculate the ET_a across the study period.

The green CWU (CWU_G) and blue CWU (CWU_B) for each crop were calculated following the method described by Hoekstra et al. [3], as in Equations (7) and (8). ET_G and ET_B represent the green and blue water evapotranspiration, respectively, from the first day (d = 1) to the end of the growing season (lgp). The factor of 10 is to convert the water depth to millimeters, and then to water volume per land area (m³ ha⁻¹).

$${\mathrm{CWU}}_{\mathrm{G}}=10\times {\sum }_{\mathrm{d}=1}^{\mathrm{lgp}}{\mathrm{ET}}_{\mathrm{G}}$$

(7)

$${\mathrm{CWU}}_{\mathrm{B}}=10\times {\sum }_{\mathrm{d}=1}^{\mathrm{lgp}}{\mathrm{ET}}_{\mathrm{B}}$$

(8)

ET_G and ET_B are the evapotranspiration (mm) by crops with respect to the green and blue components of the crop water footprint, respectively.

2.7. Assessment of Crop Water Footprint

The total water consumption (WF, m³ t⁻¹) of a crop is the sum of the blue and green components of the WF, as shown in Equations (9)–(11) (Hoekstra et al., 2011). Both the blue and green WFs for a given crop was calculated by dividing the crop water consumption (CWU, m³ ha⁻¹) by the crop productivity (CP, t ha⁻¹):

$$\mathrm{WF}={\mathrm{WF}}_{\mathrm{B}}+{\mathrm{WF}}_{\mathrm{G}}$$

(9)

$${\mathrm{WF}}_{\mathrm{B}}=\frac{{\mathrm{CWU}}_{\mathrm{B}}}{\mathrm{CP}}$$

(10)

$${\mathrm{WF}}_{\mathrm{G}}=\frac{{\mathrm{CWU}}_{\mathrm{G}}}{\mathrm{CP}}$$

(11)

where WF_B and WF_G are the blue and green water footprint, CWU_G and CWU_B are the crop consumption of green (precipitation) and blue (surface and groundwater) water, and CP is the crop productivity (yield) based on the crop’s water demand and actual transpiration output from the CROPWAT/SSEB model.

WF_B is calculated based on field data as the amount of irrigation water applied is considered blue water, as given in Equation (12) [26]:

$${\mathrm{WF}}_{\mathrm{MB}}=\frac{10\times \left(\mathrm{Ir}-\left(\mathrm{DP}+\mathrm{RO}\right)\right)}{{\mathrm{Y}}_{\mathrm{f}}}$$

(12)

where WF_MB (m³ t⁻¹) is the WF_B calculated from the recorded irrigation data, Y_f is the crop yield (carrots and onions), I_r is the amount of water used throughout the irrigation season (mm), DP is the deep percolation of water out of the root zone (mm), and RO is the surface runoff. Since measuring losses in such large fields is challenging, an average irrigation efficiency of 70% was used for center pivot irrigation systems to account for total losses due to runoff, surface flow, and deep percolation [27]. The total crop water requirement for the studied crops was estimated at 1940 mm and 824 mm for carrots and onions, respectively. However, the total volume of water used to irrigate the crops was 2658 mm and 1199 mm for carrots and onions, respectively.

2.8. Statistical Analysis

The crop WF estimates based on the CROPWAT/SSEB models were compared with the actual field data using different accuracy indicators, namely, root mean square error (RMSE), mean bias error (MBE), and relative error (RE). The use of these accuracy indicators allows for a more precise assessment of the generated models’ performance [19]. Statistical analysis was performed using the analysis of variance (ANOVA) statistical tool in the Statistical Analysis System (SAS) for Windows (v. 9.4, SAS Institute, Inc. (Cary, NC, USA)).

3. Results

The average daily minimum and maximum air temperatures were 9 °C (January) and 47 °C (August), respectively, with an annual mean temperature of 34 °C. The annual mean relative humidity was 22.6%, with peak values (36–47%) from November to February and minimum values (10–17%) from May to September. The highest (58 km h⁻¹) and lowest (6 km h⁻¹) wind speeds in the study area were recorded during July and October, respectively. The sunshine hours in the study area ranged between 9.5 and 13.5 h, with an annual mean value of 11 h.

The reference evapotranspiration and rainfall data (Figure 3 and Figure 4) showed that a mean annual rainfall of 68 mm occurred mainly between February and April. The mean peak reference evapotranspiration (ET_o) of about 20 mm d⁻¹ was recorded from May to August, with an average daily ET_o of 13.2 mm d⁻¹ and an average monthly ET_o of 375 mm month⁻¹, ranging between 138 mm (January) and 588 mm (August).

3.1. Crop Yields

The field-collected (i.e., sampled) harvest, total yield (Y), commercial harvest (CH) and aboveground biomass (AGB) of the tested crops (onions, carrots) are given in

Table 3. Commercial harvest, aboveground biomass, and total yields of carrot and onion crops.

Crop	Yield (t ha−1)	NDVI			Mean
		0.18–0.35	0.35–0.50	0.50–0.62
Carrots	Total yield	61.85	74.90	85.61	74.12
	Commercial yield	47.22	54.64	67.21	56.36
	Aboveground biomass	13.30	15.47	17.34	15.37
Onions	Total yield	36.54	44.75	51.63	44.31
	Commercial yield	30.93	36.05	44.16	37.05
	Aboveground biomass	3.84	5.22	5.81	4.96

. They were subsequently analyzed with respect to NDVI classes (Table 3). The mean carrot yield ranged from 47.2 to 74.1 t ha⁻¹, with an average commercial harvest (CH) of 56 t ha⁻¹, while the onion bulb yield was estimated at 44 t ha⁻¹; with a commercial yield of 37 t ha⁻¹.

3.2. Yield Prediction Models

A total of nine satellite images over the study area, which showed less than 1% cloud cover (i.e., cloud-free), were selected for the study period of 1 February to 18 July 2020. The results indicated that the reflectance at 36 days after sowing (DAS) increased in the NIR band and decreased in the R, G, and B bands (Figure 5). The amplitude decreased across the R, G, and B bands at 40 DAS in combination with the highest growth stage of carrot roots [28,29,30].

A linear regression analysis was performed with the Sentinel-2 bands and vegetation indices (RDVI, NDVI, SAVI, OSAVI, EVI, and RVI) layers against the actual yield to generate carrot and onion yield prediction models. Subsequently, the NDVI layer was used along with L8 TIRS bands processed for evapotranspiration (ET) mapping (Figure 6 and Figure 7). The most accurate crop productivity models are presented in

Table 4. Onion and carrot yield prediction models (* significant at 0.05; ** significant at 0.01).

Crop	Model No.	Prediction Model	Model	Cross-Validation
			R2	R2	MBE (%)	RMSE (%)
Carrots	M1	1143.6 × NIR 2212 459.51	0.77 **	0.64 **	7.82	13.41
	M2	973.1 × EVI − 226.51	0.69 **	0.62 **	–17.46	9.21
	M3	962.86 × RDVI − 219.74	0.58 **	0.59 **	5.98	10.43
Onions	M1	915.78 × NIR − 316.2	0.68 *	0.61 *	–17.19	12.65
	M2	1756.4 × RDVI − 360.81	0.72 **	0.69 **	15.21	17.67
	M3	1314.3 × EVI − 253.29	0.52 **	0.49 **	–6.19	11.24

. The NIR band was found to be appropriate for predicting the productivity of onions with an R² value of 0.68 (p > F = 0.02) and carrots with an R² value of 0.77 (p > F = 0.02). On the other hand, the findings of this study have proven that RDVI and EVI can be used for the prediction of onion and carrot productivity at bias values of 15% (R² = 0.72) and −17% (R² = 0.69), respectively.

3.3. Evapotranspiration Mapping

The temporal changes in the reference ET (ET_o) during the experiments was found to be 272 mm in February and more than 600 mm in July.

Table 5. Crop water use (ET_a, mm) values were estimated based on remote sensing models. MBE = mean bias error (%); RMSE = Root mean square error (%).

Month	Carrots			Onions
	ETa	RMSE (%)	MBE (%)	ETa	RMSE (%)	MBE (%)
February				272	−3.63	−13.2
March	234.0	−1.97	−3.9	351	−2.76	−7.6
April	391.6	7.82	61.1	378	−3.33	−11.1
May	560.2	−0.77	−0.6
June	522.0	−6.06	−36.7
July	612.1	−4.67	−21.8
Overall	2319.9	−1.13	−12.6	1001	−3.24	−9.0

summarizes the results of the L8 estimated carrot and onion evapotranspiration (ET_c) and the ET values determined from remote sensing data are depicted in Figure 6 and Figure 7 for onion and carrot crops, respectively. The SSEB-based ET_a was 2320 mm and 1001 mm for carrots and onions, respectively (Table 5).

4. Discussion

With the use of Sentinel-2-generated empirical models, a yield map of tested crops was created, as illustrated in Figure 6a (onions) and Figure 8a (carrots). The crop water use (CWU, i.e., ET_c) maps generated with the use of Landsat-8 data are depicted in Figure 6b (onions) and Figure 8b (carrots). In order to examine the accuracy of predicted maps, an error assessment of model performance was performed. As illustrated in Figure 9 and

Table 6. Accuracy assessment of the generated CWU (i.e., ET_a) map against the amount of water actually applied (AWA).

		Reference Data
		ETo	ETa (Onions)	AWA (Onions)	ETa (Carrots)	AWA (Carrots)	SUM	User’s Accuracy
Map data	ETo	18	1	2	0	0	21	85.71%
	ETa (onions)	3	19	1	2	1	26	73.07%
	AWA (onions)	0	5	17	4	1	27	62.96%
	ETc (carrots)	4	2	1	19	6	32	59.37%
	AWA (carrots)	1	1	2	2	23	29	79.31%
	SUM	26	28	23	27	31	135	72.08%
	Producer’s accuracy	69.23%	67.85%	73.91%	70.37%	74.19%	71.11%
	Overall accuracy	(18 + 19 + 17 + 19 + 23)/135 = 71.11%

, the performance of prediction models showed a good to moderate ability to predict yield. Furthermore, the highest concentration of samples was 136 kg/pixel (carrots) and 96 kg/pixel (onions), corresponding to an average yield of 68 t ha⁻¹ and 48 t ha⁻¹ for carrots and onions, respectively. The obtained yield range is consistent with the actual average yield per hectare of the studied crops. The coefficient of determination (R²) during model cross-validation was found to be 0.69 and 0.62 for carrot and onion crops, and this is a good indicator for the validity of the predicted yield.

Empirical equations for crop productivity forecasting were generated by regressing the harvested yield against the vegetation indices (VIs) calculated with the datasets provided by the Sentinel-2A satellites. Landsat-8 datasets were used for the estimation of crop ET by the SSEB model. For the carrot crop, the correlation between the actual ET_c and the predicted ET_c resulted in an R² value of 0.97, and a slope of 0.971. The obtained results are in agreement with the findings of Bezerra et al. [31], who estimated the ET_c employing the SSEB approach; their results showed R² values of 0.94 for the calibration and 0.90 for the validation. Additionally, these results are in agreement with Saggi and Jain [32] and Granata [33], who estimated ET_c using machine learning models, with R² values ranging between 0.95 and 0.99. The producer accuracy, overall accuracy, and user accuracy were used to determine the effectiveness of the comparison between the CWU (i.e., ET_a) map and the actual applied irrigated water (Table 6). The overall accuracy of the generated CWU (i.e., ET_a) map employing the SSEB algorithm was 71.11%.

The results of the SSEB for the ET and green and blue components of WF of carrot and onion fields are given in Table 7. The crop water use of the green portion (i.e., rain water evapotranspiration—ET_G) for the carrot field was 68 (±11) mm. The blue component (i.e., irrigation water—ET_B) ranged between 2187 mm and 2394 mm. However, the mean green, blue, and WF_G+B of carrots were found to be 15.6, 283, and 312 m³ t⁻¹, respectively. In the case of onions, however, the WF_G, WF_B, and WF_G+B were 8.6, 227, and 230 m³ t⁻¹, respectively.

The WF result of the carrot crop (312 m³ t⁻¹) was 31% lower than what was reported (450 m³ t⁻¹) by Multch et al. [34]. The WF of onions (230 m³ t⁻¹) is comparable with the values (280 m³ t⁻¹) reported by Penaloza-Sanchez et al. [35]. Moreover, it was found that the green and blue WF components calculated using the SSEB model in this study exceeded the global averages reported by Mekonnen and Hoekstra [36] by 18% and 133% for onion and carrot crops, respectively. These differences may be because the global statistics use both irrigated and rain-fed crop data to monitor the WF of global crops.

Variations in the monthly ET due to changes in weather conditions significantly affected the crop WF, as higher temperatures in the summer season increased crop water use due to higher evaporative demand, while lower temperatures in the autumn and the winter seasons resulted in a longer growing season. Despite the scarcity of water resources in Saudi Arabia, vegetable production depends mainly on irrigation water, and to ensure food security, optimal management to ensure high efficiency of crop water use and yield may contribute significantly to reducing water use in food production.

Table 7. Water footprint (WF) of studied crops calculated with SSEB model and published data.

Crop	Predicted Yield-YP (kg ha−1)	CWU (i.e., ETa, mm)	Water Footprint—WF (m3 t−1)			Reference
			WFG	WFB	WFG+B
Carrots	74	2320	16	283	312	This study
			23	427	450	[34]
			106	28	134	[36]
			104	12	116	[37]
			0	114	114	[38]
Onions	44	1001	9	227	230	This study
			18	243	261	[34]
			176	53	229	[35]
			192	88	280	[36]

Table 7. Water footprint (WF) of studied crops calculated with SSEB model and published data.

Crop	Predicted Yield-YP (kg ha−1)	CWU (i.e., ETa, mm)	Water Footprint—WF (m3 t−1)			Reference
			WFG	WFB	WFG+B
Carrots	74	2320	16	283	312	This study
			23	427	450	[34]
			106	28	134	[36]
			104	12	116	[37]
			0	114	114	[38]
Onions	44	1001	9	227	230	This study
			18	243	261	[34]
			176	53	229	[35]
			192	88	280	[36]

5. Conclusions

This study was carried out to estimate the productivity and water footprint of carrot and onion crops grown in the central regions of Saudi Arabia under a center-pivot irrigation system. Empirical equations for yield forecasting were generated by regressing the harvested yields against the vegetation indices (VIs) extracted from datasets provided by Sentinel-2A satellites. Landsat-8 datasets were used for the estimation of crop ET by the SSEB model. The results revealed that the NIR band was appropriate for yield prediction of both onions (R² = 0.68, p > F = 0.02) and carrots (R² = 0.77, p > F = 0.02). RDVI and EVI showed the best results for predicting the crop yields, with bias values of 15% for onions and –17% for carrots. The SSEB-estimated CWU amounts are higher than the actually applied quantities by 10.6% (onions) and 12.6% (carrots). The crop water footprint values obtained in this study are somewhat higher than earlier published results. However, the outcomes of this study can be a reference for crop irrigation management practices in the study region and areas with similar environmental conditions. Under the arid climatic conditions of Saudi Arabia, the return flow of irrigation water and recharge of ground water is almost nil due to sparse rainfall. Thus, water pollution through agricultural inputs (fertilizers, pesticides, etc.) is essential and should be considered grey WF for the maintenance of soil quality when quantifying the WF.