Penman–Monteith Reference Evapotranspiration Estimation Models, Using Latitude–Temperature Data, in the State of Sinaloa, Mexico

1. Introduction

Historically, to guarantee the feeding of the world population, agriculture has been the activity that has consumed the greatest amount of water [1]. Approximately 40% of the world’s food depends on activities inherent to agricultural irrigation [2]. This constantly increasing water demand [3] can trigger significant meteorological droughts [4,5,6], which are accentuated in arid regions [6, 7], where the incident solar radiation (SR) is more intense [8, 9]. Intense SR causes approximately 60% of precipitation to return to the atmosphere in the form of evapotranspiration [10, 11], causing these regions to be classified as vulnerable to desertification [12]. For example, in semi-arid regions, agricultural irrigation is a parameter that should trend towards intelligent irrigation [13,14,15,16]. To develop intelligent irrigation, valuable information must be available that establishes the relationship between crop growth and water balance [17], in which reference evapotranspiration (ETo) is essential [18]. According to [19] and [20], ETo is the potential evapotranspiration of a hypothetical grass surface, with uniform height, well-watered and with active growth, and which depends entirely on climatological variables [9, 21]. According to [9, 22–24], it is always advisable to use empirical equations to estimate ETo by the Penman–Monteith (PM_C) method, even when data is lacking, mainly because it remains the most precise method. Of alternative methods, Hargreaves (HA_C) continues to be the most used, mainly due to its high accuracy/number of variables used ratio [25,26,27]. However, [26] and [27], state that another possible way to estimate ETo by Penman–Monteith is through simple linear regressions (SLR) and simple nonlinear regressions (SNR); PM_R (dependent variable) vs HA_C (independent variable), which more accurately calculates the hydric requirements of crops.

In Mexico, approximately 77% of the volume of the total water resource is allocated to the agricultural sector, and two thirds of the national territory is characterized by an aridity index ranging from arid to semi-arid [13]. In particular, the state of Sinaloa has a predominantly semi-arid climate [7], and according to [13] and [23], this condition predisposes it to focus efforts on the characterization of PM_R, as well as the subsequent design and administration of intelligent irrigation systems [13, 23]. Intelligent irrigation could improve the volumes of yields of Sinaloan crops, as well as encourage the conservation of water resources [23, 28].

In this study, daily series (1979–2017) of minimum (Tmin) and maximum (Tmax) temperatures were obtained from seven weather stations in Sinaloa from National Water Commission (CONAGUA) [29]. PM_C, HA_C and PM_R were calculated. At another weather station, Acaponeta, observed daily series (2005–2008, 2011–2013 and 2015–2017) were obtained of mean temperature (Tmea_O), incident solar radiation (SRg), average relative humidity (${\mathrm{R}\mathrm{H}}_{\mathrm{O}}$) and average wind speed at a height of 10 m (${\mathrm{U}}_{\mathrm{z}\mathrm{O}}$). The data for the eight stations were provided by the CONAGUA [29] and CONAGUA–National Meteorological Service (SMN) (CONAGUA–SMN) [30]. At Acaponeta, daily observed Penman–Monteith ETo (PM_O) was calculated. For validation, three SLR were obtained: SR vs SRg, PM_C vs PM_O and PM_R vs PM_O. A hypothesis test was applied: Pearson correlation (Pr) vs Pearson critical correlation (Pcr). In the three SLR, the condition Pr > Pcr was met; that is, all Pr were significantly different from zero [31].

The goal was to create PMR estimation models using latitude–temperature data for the state of Sinaloa, Mexico.

Although most of the weather stations for public use in Sinaloa lack the full set of climate variables necessary for the calculation of PM_O [7], in this study, predictive models of PM_R are provided using the variables latitude–temperature. These models can help ensure the feeding of “the breadbasket of Mexico” through intelligent irrigation [13, 23].

2. Materials and Methods

2.1. Study Area

Sinaloa is in the northwest of Mexico (Figure 1), and because it is the most important agricultural state in Mexico, it is called “the breadbasket of Mexico” [32]. Furthermore, this state is the main producer of export-oriented crops [33] cited by [32]. Due to the planted area and sensitivity to extremes of RH_O–Tmax–Tmin, two of the most important crops in Sinaloa are corn and beans [28].

2.2. Data

2.2.1. Daily Maximum ($\mathrm{T}\mathrm{m}\mathrm{a}\mathrm{x}$) and Minimum Temperature ($\mathrm{T}\mathrm{m}\mathrm{i}\mathrm{n}$)

Using data from CONAGUA (https://smn.conagua.gob.mx/es/climatologia/informacion-climatologica/informacion-estadistica-climatologica)[29], daily series of Tmax and Tmin were obtained from 70 weather stations in Sinaloa for the period 1942–2019. These same series were previously obtained by [34]. Through a review of the availability of recent information (< 5% missing data), in this study it was decided to work with seven weather stations (Culiacán, El Playón, Las Tortugas, Rosario, La Concha, Ixpalino and Sanalona II (Figure 1), for the period 1979–2017.

2.2.2. Imputation of Missing Data, Homogenization of Series and Determination of the Mean Daily Temperature (Tmn)

Using RStudio software, with the Climatol library [35] and the orthogonal regression method, missing daily data of Tmax and Tmin were estimated by imputation. Using the standard normal homogeneity test (SNHT) [36] method, with Climatol, the series were also homogenized. By means of the semi-sum of the complete and homogeneous series of Tmax and Tmin, the daily series of Tave were determined.

In general, the greatest thermal extremes were registered in Ixpalino (Tmax = 46.50 °C day⁻¹), El Playón and Las Tortugas (Tmin = −6.00 C day⁻¹) and El Playón (Tmn = 38.00 C day⁻¹, Table 1).

2.2.3. Wind Speed at 10 m Height (U_z)

Through the National Oceanic and Atmospheric Administration (NOAA)

[37] (https://downloads.psl.noaa.gov/Datasets/ncep.reanalysis2/Monthlies/gaussian_grid/), the monthly series (Jan–Dec) of wind speed at a height of 10 m (U_z) were obtained for the period 1979–2017. Due to the availability of satellite information, U_z was obtained for only two coordinates in the state of Sinaloa: 1) 25° 43′ 14″ N by 108° 45′ 00″ W and 2) 23° 55′ 42″ N by 106° 46′ 48″ W.

2.3. Empirical Equations to Estimate Penman–Monteith Reference Evapotranspiration, Calculated with Missing Data (PM_C) and Observed Data (PM_O)

2.3.1. Wind Speed at 2 m Height (U₂)

Although [9, 38], state that wind speed is not very relevant for estimating PM_C in semi-arid regions, in this study, using Equation 1, wind speed was obtained at a height of 2 m (U₂) [19, 22].

$${\mathrm{U}}_2={\mathrm{U}}_{\mathrm{z}}·\left[\frac{4.87}{\mathrm{l}\mathrm{n}·\left(67.8·\mathrm{z}-5.42,\right)}\right],$$

(1)

where U₂ = monthly average wind speed at a height of 2 m (m s^–1), U_z = average wind speed measured at a height of 10 m (m s^–1) and z = measurement height of U_z (m).

Since PM_C is the international standard because of its greater measurement accuracy [9], in this study, PM_C was estimated daily using Equations 2–10. These equations, which are recommended by [19, 22] when there are missing data, are given as follows:

$${\mathrm{e}}_{\mathrm{a}\mathrm{C}}=0.6108·\mathrm{e}\mathrm{x}\mathrm{p}\left(\frac{17.27·\mathrm{T}\mathrm{m}\mathrm{i}\mathrm{n}}{\mathrm{T}\mathrm{m}\mathrm{i}\mathrm{n}+237.3}\right),$$

(2)

where e_aC = calculated actual vapor pressure and Tmin = daily minimum air temperature (°C).

$${\mathrm{e}}_{\mathrm{s}}=\left(0.6108·\mathrm{e}\mathrm{x}\mathrm{p}\frac{17.27·\mathrm{T}\mathrm{m}\mathrm{n}}{\mathrm{T}\mathrm{m}\mathrm{e}\mathrm{a}+237.3}\right),$$

(3)

where e_s = saturation vapor pressure (kPa) and Tmn = daily mean air temperature (°C).

$$∆=\frac{4098·\left[0.6108·\mathrm{e}\mathrm{x}\mathrm{p}\frac{17.27·\mathrm{T}\mathrm{m}\mathrm{n}}{\mathrm{T}\mathrm{m}\mathrm{n}+237.3}\right]}{{\left(\mathrm{T}\mathrm{m}\mathrm{n}+237.3\right)}^2},$$

(4)

where ∆ = slope of the saturated vapor pressure curve (kPa °C^–1).

$$\mathrm{S}\mathrm{R}={\mathrm{K}}_{\mathrm{R}\mathrm{S}}·{\left(\mathrm{T}\mathrm{m}\mathrm{a}\mathrm{x}-\mathrm{T}\mathrm{m}\mathrm{i}\mathrm{n}\right)}^{0.5}·\mathrm{R}\mathrm{a},$$

(5)

where SR = calculated incident solar radiation (MJ m^–2 day^–1), KRS = solar radiation adjustment coefficient (dimensionless), with a value of 0.16 for continental conditions and 0.19 for coastal conditions (in this study, KRS = 0.16 was used), Tmax = maximum daily air temperature (°C) and Ra = extraterrestrial solar radiation (obtained through tabulated values, with respect to latitude, MJ m^–2 day^–1).

$$\mathrm{S}\mathrm{R}\mathrm{o}=0.75·\mathrm{R}\mathrm{a},$$

(6)

where SRo = incident solar radiation with clear sky (MJ m^–2 day^–1).

$$\mathrm{R}\mathrm{n}\mathrm{l}=\left[(\mathsf{\sigma }·{\mathrm{T}\mathrm{m}\mathrm{n}\mathrm{K}}^4)·\left(0.34-0.14·{\mathrm{e}}_{\mathrm{a}}^{0.5}\right)·\left(1.35·\frac{\mathrm{S}\mathrm{R}}{\mathrm{S}\mathrm{R}\mathrm{o}}-0.35\right)\right]$$

(7)

where Rnl = net longwave radiation (MJ m^–2 day^–1), σ = Stefan–Boltzmann constant (0.4903 × 10–8 MJ K^–4 m^–2 day^–1) and TmnK = mean daily air temperature (°K⁴).

$$\mathbf{R}\mathbf{n}\mathbf{s}=0.77·\mathbf{S}\mathbf{R}$$

(8)

where Rns = net shortwave radiation (MJ m^–2 day^–1).

$$\mathrm{R}\mathrm{n}=\mathrm{R}\mathrm{n}\mathrm{s}-\mathrm{R}\mathrm{n}\mathrm{l},$$

(9)

where Rn = net radiation (MJ m^–2 day^–1).

$${\mathrm{P}\mathrm{M}}_{\mathrm{C}}\mathrm{a}\mathrm{n}\mathrm{d}{\mathrm{P}\mathrm{M}}_{\mathrm{O}}=\frac{0.408·∆·\left(\mathrm{R}\mathrm{n}-\mathrm{G}\right)+\mathsf{\gamma }·\frac{900}{\mathrm{T}\mathrm{m}\mathrm{n}+273}·{\mathrm{U}}_2·({\mathrm{e}}_{\mathrm{s}}-{\mathrm{e}}_{\mathrm{a}\mathrm{C}})}{∆+\mathsf{\gamma }·(1+0.34·{\mathrm{U}}_2)},$$

(10)

where PM_C = Penman–Monteith reference grass evapotranspiration (mm day^–1, calculated with missing data), PM_O = Penman–Monteith observed grass reference evapotranspiration (mm day^–1, at the Acaponeta station), G = soil heat flux density (MJ m^–2 day^–1, null for daily estimates) and γ = psychrometric constant (0.0677 kPa °C^–1).

2.4. Calculated Hargreaves Reference Evapotranspiration (HAC, Alternative Method Used)

When the absence of data does not allow Equation 10 to be used, [25] recommend the use of Expression 11 to estimate HAC, which is widely recommended worldwide, due to the high ratio accuracy/number of variables used.

$${\mathrm{H}\mathrm{A}}_{\mathrm{C}}=0.0023·\mathrm{R}\mathrm{a}·\left(\mathrm{T}\mathrm{m}\mathrm{n}+17.8\right)·{\left(\mathrm{T}\mathrm{m}\mathrm{a}\mathrm{x}-\mathrm{T}\mathrm{m}\mathrm{i}\mathrm{n}\right)}^{0.5},$$

(11)

where HA_C = Hargreaves reference evapotranspiration (mm day^–1).

PM_C and HA_C were also calculated as monthly (Jan–Dec), seasonal (Mar–Aug) and annual (Jan–Dec) averages.

2.5. Pre-Validation

2.5.1. Normality Test and Correlation Coefficients

A Shapiro–Wilk normality test was applied to all PM_C and HA_C series [39]. To find out whether PM_C and HA_C were significantly correlated, a Pr was applied to the series that presented normality and a Spearman correlation (Sr) was applied to the series that did not present normality.

2.5.2. Simple Linear Regressions (SLR) and Simple Nonlinear Regressions (SNR)

To generate sensitive models [7] to predict PM_R (dependent variable) based on HA_C (independent variable), SLR were initially fitted (Equation 12). A Shapiro–Wilk normality test was applied to the SLR residuals. When the residuals were not normal, a SNR (10 different functions) was applied, fitting a curvilinear estimate. Of the 10 functions, the following were chosen: a) exponential function (monthly series, Equation 13) and potential function (seasonal series, Equation 14), which were selected due to the highest R² recorded.

$${\mathrm{P}\mathrm{M}}_{\mathrm{R}}=\mathrm{a}+\mathrm{b}·{\mathrm{H}\mathrm{A}}_{\mathrm{C}},$$

(12)

$${\mathrm{P}\mathrm{M}}_R=\mathrm{a}·{\mathrm{e}}^{\mathrm{b}{·\mathrm{H}\mathrm{A}}_{\mathrm{C}}},$$

(13)

$${\mathrm{P}\mathrm{M}}_R=\mathrm{a}·{{\mathrm{H}\mathrm{A}}_{\mathrm{C}}}^{\mathrm{b}},$$

(14)

where e = Euler number (2.7182) and a, b = regression coefficients which describe the relationship between PM_R and HA_C.

2.5.3. Hypothesis Test

For each SLR and SNR, the Pr and Sr were obtained by the square root of R² (Sections 2.5.2). To find out if each Pr and Sr were significantly different from zero, hypothesis tests were applied [31, 40]. Each Pr and Sr (Section 2.5.1) were compared with a Pcr = |0.316| (Equation 15) and a critical Spearman correlation coefficient (Scr = 0.318, Equation 16).

$$\mathrm{P}\mathrm{c}\mathrm{r}=\sqrt{\frac{\frac{{\mathrm{t}}_{\mathrm{c}}^2}{\mathrm{d}\mathrm{f}}}{\frac{{\mathrm{t}}_{\mathrm{c}}^2}{\mathrm{d}\mathrm{f}}+1}},$$

(15)

where t_c = critical value of the student t statistic and df = degrees of freedom (n–2).

$$\mathrm{S}\mathrm{c}\mathrm{r}=\pm \mathrm{z}\sqrt{\mathrm{n}-1},$$

(16)

where z = 1.96, n = 39 (for the period 1979–2017).

The design of the hypotheses is shown in Equations 17–18:

Also, the root mean square error (RMSE) between PM_C and PM_R was calculated.

2.6. Validation

Using the CONAGUA–SMN database (https://smn.conagua.gob.mx/tools/GUI/sivea_v3/sivea.php)[30], the following observed data were obtained from the Acaponeta station: Tmn_O, U_ZO, SRg and RH_O, for the periods 2005–2008, 2011–2013 and 2015–2017. PM_O was calculated, reapplying Equations 1, 3–4, 6–10 and 19.

$${\mathrm{e}}_{\mathrm{a}\mathrm{O}}=\frac{{\mathrm{R}\mathrm{H}}_{\mathrm{O}}}{100}·{\mathrm{e}}_{\mathrm{s}},$$

(19)

where e_aO = observed actual vapor pressure and RH_O = observed mean daily relative humidity (%).

Three SLR were applied: 1) SR (La Concha) vs SRg (Acaponeta), 2) PM_C (La Concha) vs PM_O (Acaponeta) and 3) PM_R (La Concha) and PM_O (Acaponeta). A Shapiro–Wilk normality test was applied to the residuals of the three SLR. From each SLR, Pr = (R²)^0.5 was obtained. To find out if Pr were significantly different from zero, another hypothesis test was carried out between Pr vs Pcr = |0.576| (for n = 12). Finally, the RMSE values were calculated between the calculated and observed values of the three SLR. The pre-validation and validation were an adaptation of the development by [7].

2.7. Software Used

To carry out this research, the following programs were used: RStudio version 4.3.0, Past version 4.08, XLstat version 2023, Panoply version 5.2.6 and CorelDRAW version 2019.

3. Results

3.1. Calculated Monthly Average Reference Evapotranspiration: Penman–Monteith (PM_C) and Hargreaves (HA_C) Methods

The average ETo ranged from PMC = 1.483 mm day^–1 in 1992 (Culiacán–Jan, Figure 2a) to PMC = 6.656 mm day^–1 in 1982 (El Playón–May, Figure 2b); and from HAC = 2.256 mm day^–1 in 1991 (Culiacán–Dec, Figure 2a) to HAC = 8.133 mm day^–1 in 2002 (Sanalona II–May, Figure 2g).

3.2. Normality Test for the Calculated Average Reference Evapotranspiration: Penman–Monteith (PM_C) and Hargreaves (HA_C) Methods

3.2.1. Monthly (Jan–Dec), Seasonal (Mar–Aug) and Annual (Jan–Dec) Series

For PM_C–Ixpalino in all months (Jan–Dec), p(normal) and W ranged from 0.090 to 0.623 and from 0.951 to 0.978, respectively (Figure 3a). In total, 37 monthly series did not present normality; PM_C = 15 series and HA_C = 22 series (Figure 3a). The seasonal series (Mar–Aug) did not present normality [p(normal) < 0.05] were PM_C–El Playón [p(normal) = 7.8 × 10–5], HA_C–El Playón [p(normal) = 1.7 × 10–6], HA_C–Las Tortugas [p(normal) = 0.003] and HA_C–Rosario [p(normal) = 0.008, Figure 3b]. The annual series (Jan–Dec) without normality were PM_C–El Playón [p(normal) = 2.3 × 10–4], HA_C–El Playón [p(normal) = 9.4 × 10–7], HA_C–Las Tortugas [p(normal) = 0.006] and HA_C–Ixpalino [p(normal) = 0.018, Figure 3b].

3.3. Pearson (rP) and Spearman (rS) Correlations of Calculated Average Reference Evapotranspiration: Penman–Monteith (PMC) and Hargreaves (HAC) Methods

3.3.1. Monthly Correlations (Jan–Dec)

As shown in Table 2, the correlations ranged from rP = 0.443 (El Playón–Jul) to rP = 0.929 (Las Tortugas–Jan). All rP and rS were significantly different from zero (rP > rcP = |0.316| and rS > rcS = |0.318|).

3.3.2. Seasonal (Mar–Aug) and Annual (Jan–Dec) Correlations

As shown in Table 3, all seasonal (Mar–Aug) and annual (Jan–Dec) rP and rS were significantly different from zero (rP > rcP = 0.316 and rS > rcS = 0.318;). Seasonal correlations (Mar–Aug) ranged from rP = 0.693 (Sanalona II) to rP = 0.907 (La Concha). The annual correlations (Jan–Dec) ranged from rS = 0.831 (El Playón) to rP = 0.921 (La Concha).

3.4. Linear (SLR) and Simple Nonlinear Regressions (SNR) of Calculated Average Reference Evapotranspiration: Penman–Monteith (PM_R, Dependent Variable) and Hargreaves (HA_C, Independent Variable) Methods

3.4.1. Normality Test Of Monthly (Jan–Dec), Seasonal (Mar–Aug) and Annual (Jan–Dec) Residuals

The only series of monthly residuals (Jan–Dec) that did not present normality was Sanalona II–Oct [p(normal) = 0.046 and W = 0.942, Figure 4a]. In the normal monthly series, the p(normal) values ranged from 0.059 (El Playón–May) to 0.951 (Culiacán–Nov, Figure 4a). As seen in Figure 4b, Las Tortugas for the seasonal period (Mar–Aug) [p(normal) = 0.028 and W = 0.936] was the only series that did not register normality. In the normal seasonal series, the p(normal) values ranged from 0.237 (Ixpalino) to 0.445 (Rosario). In the normal annual series, the p(normal) values ranged from 0.221 (Ixpalino) to 0.964 (La Concha).

3.4.2. Monthly Coefficients and Goodness of Fit (Jan–Dec)

The fit ranged from R² = 0.196 (rP = 0.443, El Playón–Jul) with RMSE = 0.274 mm day^–1 to R² = 0.863 (rP = 0.929, Las Tortugas–Jan) with RMSE = 0.218 mm day^–1 (Table 4). For SNR–exponential function (Sanalona II–Oct, Table 4), R² = 0.706 (rS = 0.840 > rcS = |0.318|).

3.4.3. Coefficients and Seasonal (Mar–Aug) and Annual (Jan–Dec) Goodness of Fit

Seasonal fit (Mar–Aug, Table 5) ranged from R² = 0.480 (rP = 0.693, Sanalona II) with RMSE = 0.156 mm day^–1 to R² = 0.823 (rP = 0.907, La Concha) with RMSE = 0.117 mm day^–1. Annual fit (Jan–Dec) ranged from R² = 0.719 (rP = 0.848, Sanalona II) with RMSE = 0.112 mm day^–1 to R² = 0.848 (rP = 0.921, La Concha) with RMSE = 0.082 mm day^–1. For SNR–potential function (Las Tortugas–seasonal, Table 5), the fit was R² = 0.699 (rS = 0.836 > rcS = |0.318|).

3.5. Validation

3.5.1. Simple Linear Regressions (SLR) between Calculated and Observed Values from: 1) Incident Radiation (SR vs SRg), 2) Penman–Monteith Reference Evapotranspiration, Calculated with Equations (PM_C vs PM_O) and 3) Calculated with Regressions (PM_R vs PM_O)

All the monthly average SLR (Figures 5a–5c) recorded rP significantly different from zero (rP > rcP = |0.576|, for n = 12). Specifically for SR vs SRg, the measures of fit were: R² = 0.905, rP = 0.951 and RMSE = 0.684 mm day^–1 (Figure 5a). In PM_C vs PM_O, the measures of fit were R² = 0.350, rP = 0.592 and RMSE = 0.590 mm day^–1 (Figure 5b). For PM_R vs PM_O, the measures of fit were R² = 0.391, rP = 0.625 and RMSE = 0.578 mm day^–1 (Figure 5c). The residuals of the three SLR presented normality: p(normal) = 0.193 and W = 0.907 (Figure 5a), p(normal) = 0.344 and W = 0.927 (Figure 5b) and p(normal) = 0.464 and W = 0.937 (Figure 5c).

4. Discussion

The results of Figure 2a are similar to those reported by [41], who found a range from PM_C = 3.0 mm day^–1 to PM_C = 5.8 mm day^–1 for the Culiacán valley, in the period 2013–2014. The variation between PM_C vs HA_C ranged from RMSE = 1.861 mm day^–1 (El Playón, Figure 2b) to RMSE = 1.972 mm day^–1 (Culiacán, Figure 2a), that is, ETo presents a tendency towards underestimation of PM_C and overestimation of HA_C (RMSE > 0.3 mm day^–1) [9].

According to [42], the results of PM_C–Ixpalino series (Figure 3a) presents normality. According to [43], in the results of Figure 3b, the seasonal series (Mar–Aug) that did not present normality were: PM_C–El Playón, HA_C–El Playón, HA_C–Las Tortugas and HA_C–Rosario, because did not present the condition of p(normal) > 0.05. According to [43] and [44], the annual series (Jan–Dec) that did not present normality were PM_C–El Playón, HA_C–El Playón, HA_C–Las Tortugas and HA_C–Ixpalino, this due that p(normal) < 0.05.

According to [31] and [40], the results of Table 2 establish the significant monthly relationship (Jan–Dec) of PM_C vs HA_C, so monthly modeling of PM_R is appropriate (Equations 12–14), applying SLR and SNR, as recommended by [19, 22] and applied by [27, 45].

Because all the correlations of Table 3 were significant [31, 40], SLR and SNR can be applied to estimate PM_R with seasonal (Mar–Aug) and annual (Jan–Dec) scale, with HA_C as the independent variable [27,28,29,30,31,32,33,34,35,36,37,38,39,40,41,42,43,44,45].

According to [43] in the results of Figure 4b, and for the seasonal period (Mar–Aug), the only series that did not present normality was Las Tortugas, because p < 0.05. All series that did present the condition of p > 0.05 were considered as normal series [42].

All RLS (Table 4) exceeded rcP = |0.316| [31, 40] (significant correlation) and did not register a trend towards underestimation or overestimation (RMSE < 0.300 mm day^–1)[9]. In Table 4 and for Sanalona II–Oct (SNR–exponential function), the results were rS = 0.840 > rcS = |0.318| [31, 40] (significant correlation) [31, 40]. The methodology of Table 4 was applied to obtain more accurate estimates [19, 22], and was previously applied by [26, 27].

All SLR (seasonal and annual, Table 5) exceeded rcP = |0.316| [31, 40] (significant correlation) and showed no trend towards underestimation or overestimation (RMSE < 0.300 mm day^–1) [9]. These results are in agreement with [19, 22], who state that PM_R models are more accurate than when only Equations 1–10 (PM_C) are used. [46] also state that HA_C estimation is the most recommended method when data is not available to estimate PM_C.

In validation, the three SLR (Figures 5a–5c) performed well (RMSE < 1.0 mm día^–1)[47]. In this study, SR was highly influential (approximately 90.5%) for estimating PM_C, which agrees with [8], who point out that SRg in Sinaloa is decisive for the estimation of PM_C. According to results of PM_C vs PM_O (Figure 5b), Equations 1–10 were reliable and sensitive for estimating PM_C, even when the series presented missing data [9, 48, 49]. The results of Figure 5c and according to [19, 22, 26, 27], the models of this study are also reliable and sensitive for predicting PM_R. Finally, because the residuals of the three SLR (Figures 5a–5c) presented normality, the SLR are an appropriate statistical tool to use for comparison of calculated and observed data [7].

5. Conclusions

Due to the lack of data variables from weather stations of Sinaloa, PM_C and HA_C were estimated with the use of equations. PM_C presented trends towards underestimation and HA_C presented trends towards overestimation. For the first time in Sinaloa, monthly (Jan–Dec), seasonal (Mar–Aug) and annual (Jan–Dec) SLR and SNR were generated to estimate PM_R (dependent variable), using HA_C (independent variable). Although the equations are a good tool to estimate PM_C, the use of PM_R estimation models is more precise (without trends of underestimation or overestimation). To try to improve the fit of PM_R vs PM_O, in future studies it is recommended to estimate PM_R using any other alternative method for ETo, for example, Thornwaite, Priestley–Taylor, Valiantzas, Makkink, Schendel, Jensen, or Turc, among other methods. Knowledge of PM_R in Sinaloa can contribute to facilitating the calculation of crop evapotranspiration, which can enable the design of intelligent irrigation plans that are efficient, sustainable, and affordable. The PM_R models of this study are also a valuable tool when complete climate series are lacking, which are necessary for the calculation of PM_O, since in this study to obtain PM_R only latitude–temperature is required. These predictive models can also help ensure, in the near future, the feeding of the population of “the breadbasket of Mexico,” specifically through the relationship of less irrigation water/greater sustainability of food production.