1. Introduction

According to ^{Flores et al. (2012)}, the area in the state of Sinaloa is predominantly classified as subtropical-intertropical, with zones of low topographic relief near the Pacific Ocean and extreme elevations in the mountain system known as the Sierra Madre Occidental. The climatic variability-randomness of this area can be numerically modeled using criteria based on significant trends (ST), adjusted return periods (ARP) and fitted probability distribution functions (^{Cohen et al., 2012}; ^{Rustom 2012}; ^{Franzke 2015}; ^{Kienzle 2018}; ^{Chapman et al., 2019}; ^{Garry et al., 2019}; ^{Latif & Mustafa 2020}; ^{Alves et al., 2021}). These mathematical tools were used to characterize the climate variability (^{Llanes et al., 2022b}; Llanes, 2023) and randomness of Sinaloa; i.e., indicators of hot extremes (HE, Llanes et al., 2022a), since determination of the STs can be helpful for identifying zones where extreme climate changes occur, and knowledge about such changes has been demonstrated to be very useful for designing adaptation/prevention strategies to face natural disasters (^{Blanco et al., 2014}). Globally, the HE maximum maximorum (highest maximum) temperatures (MmT) are some of the most important indicators of change climate and can provide valuable information about the climate conditions to which the crops in a region are subjected (^{Li et al., 2018}; ^{Ruiz et al., 2005}). The major crops produced in the area of Sinaloa state, nicknamed "the breadbasket of Mexico" (^{Rincón 2012}), are maize, sorghum, wheat, and beans, which are particularly vulnerable to temperature extremes and droughts (SAGARPA 2015; ^{Gudko et al., 2021}; ^{Gurung et al., 2021}; ^{Habibi & Meddi 2021}). The negative effects of a high frequency of maximum extreme temperature (MeT) events can include human deaths as a consequence of gastrointestinal diseases and heat stroke (^{Murage et al., 2017}; OPS 2019).

With ARPs, the indicators of climate extremes can be estimated with very low uncertainty and the highest efficiency and reliability (^{Aguilera 2007}; ^{Campos 2014}; ^{Lustenberger et al., 2014}; ^{Deng & Zhu 2021}). Worldwide, there are multiple studies that address ST using non-parametric Mann-Kendall tests (^{Kumar et al., 2023}), as well as the use of fitted probability distribution functions, due principally to the importance, reliability, and efficiency of modeling (^{González et al., 2011}; ^{Fontanelli et al., 2019}).

In this study, based on HE indicators (^{Li et al., 2018}) of the MeT and MmT (^{Ruiz et al., 2005}), the negative climate effects in Sinaloa were calculated according to the criteria established by Li et al. (2018). According to ^{Wang et al. (2016)} and ^{Hochman et al. (2022)}, many studies indicate that extreme climate events may intensify or become more frequent and long-lasting with climate change, especially in semi-arid regions (^{Demisse et al., 2022}) such as the state of Sinaloa (^{Llanes et al., 2022b}). Due to the high interconnection with drought (^{Horton et al., 2016}; ^{Yuan et al., 2016}), the HE indicators calculated were frequency (FR), annual average duration (AAD), annual daily duration (ADD) and intensity (IN) (^{Seneviratne et al., 2012}).

To learn about previous research related to this study, a review was done of the literature (documents available locally and worldwide). ^{Jaiswal et al. (2023)} calculated extreme climate indices in India (hot days and daily maximum temperature) and ^{Llanes et al. (2022a)} calculated indicators of MeT (HE and MmT) in Sinaloa, but did not calculate their trends or return periods through fitted probability distribution functions.

The modeling tools of this study (ST, fitted probability distribution functions and ARP) provide a new method to estimate the frequency of occurrence of MeT indicators HE (FR, AAD, ADD and IN) and MmT. This method can help to prevent impacts on agricultural yields through climate-smart practices such as drip irrigation, conservation agriculture, precision fertilization, protected agriculture, and crop substitution, which will withstand future hot temperature extremes (World Bank et al., 2014), not only in Sinaloa, but throughout Mexico. For this reason, the goal of this investigation was to calculate MeT indicators based on STs and ARPs, to design climate-smart agricultural adaptation/prevention strategies.

2. Material and methods

2.1 Study area

The state of Sinaloa is located between 105° 23' 32" and 109° 26' 52" N, and 27° 02' 32" and 22° 28' 02" W (Figure 1). It comprises approximately 2.9% of the total territory of the country (INEGI 2016).

Figure 1 Location of the study area. State of Sinaloa Source: ^{Llanes et al. (2021)}  

2.2 Data collection and processing

Daily maximum temperature data for April-October were obtained for 19 weather stations in the study area from the CLImate COMputing database (CLICOM, 2021) at the website http://clicom-mex.cicese.mx/mapa.html for the period 1982-2014 (Table 1). Missing data were estimated though the multiple imputation method (^{Little & Rubin, 1987}) with the Markov chain Monte Carlo algorithm (MCMC; ^{Austin & Van Buuren, 2022}) since, according to ^{Puerta (2002)} and ^{Salgado & Largo (2018)}, with this method it is possible to obtain unbiased estimators that can be used efficiently. For this calculation, the XLstat version 2014 program was used. This program works in the Windows Excel environment and is widely used worldwide in climatology. These data were obtained from ^{Llanes et al. (2021)}

Table 1 Location of the weather stations. 

HE indicators and MmT

To obtain some of the most important variables used worldwide to characterize extreme weather, the following HE indicators were calculated following ^{Li's method (Li et al., 2018}): FR, AAD, ADD and IN. These authors define a HE event as occurring when the value of daily maximum temperature exceeds the 85^th and 90^th percentiles (P85 and P90, respectively) for three consecutive days (Stefanon et al., 2012, Zschenderlein et al., 2019 and ^{Llanes et al., 2022a}). The values of MmT were obtained according to ^{Ruiz et al. (2005)}, who defined these variables as the maximum monthly value in a series of maximum daily temperatures. To calculate the HE indicators and MmT, the XLstat version 2014 program was used. The authors used the data of the HE indicators FR, AAD, ADD and IN, previously calculated by Llanes et al. (2022a).

ST of MeT indicators

In order to determine whether local climate trends are cooling or warming, the STs of the annual average values of FR, AAD, ADD, IN and MmT were calculated.

The ^{Mann-Kendall test for trend analysis (Mann, 1945}; Kendall, 1975, Table 2) and ^{Sen's slope test (Sen 1968)} for magnitude of trend (Figure 3) were applied. To identify the presence of STs of MeT indicators, an absolute statistical threshold based on the standard normal distribution (Zstd) = |1.96| was set with significance level a = 0.05 (^{Ghasemi, 2015}). For a better interpretation of the ST variation data, the results were multiplied by ten; i.e., the time scales used and interpreted in this study were decadal. Since single-year changes ST may be small, decadal trends were chosen to make the trends clearer (^{Rojas et al., 2010}). To calculate the STs of the MeT indicators (Mann-Kendall and Sen's slope tests), Past version 3.25 and XLstat version 2014 were used.

Table 2 P-values for the fitted probability distribution functions of MeT. 

FR-frequency, AAD- annual average duration, ADD-annual daily duration, I- intensity, MmT-maximum maximorum temperature and (boxed and bold)-fitted probability distribution functions.

Source: ^{Llanes et al. (2021)}

Statistical analysis of the fitted probability distribution functions

To obtain the fitted probability distribution functions of HE indicators and MmT, in this study the maximum likelihood estimation method was used. This mathematical method was applied due to its adaptability to the different HE indicators and its precision in the estimation of large samples (data series for the period 1982-2014; ^{Gómez, 2020}). The p-values of 16 probability distribution functions were calculated for the annual average values of each MeT indicator. This statistical method was used because is considered one of the most consistent, sufficient and asymptotically unbiased ways to identify the weather parameters of the fitted probability distribution functions (^{Strupczewski et al., 2001}). In addition to being more complex than the method of moments, it allows for greater rigor for all distributions (XLstat, 2022), and a greater reliability of fit for each indicator of HE and MmT. In this study, the fitted probability distribution functions used were logarithmic (Log), log-normal (Lognormal), and generalized extreme values (GEV) (^{Chow et al., 1994}; ^{Campos, 2016}). Log, lognormal and GEV functions are shown in Equations 1, 2 and 3, respectively. These calculations were carried out using XLstat version 2014.

Adjusted return periods

To obtain the cumulative distribution function F(x) of each fitted probability distribution function, the variable (x) and the following parameters were substituted into Equations 1, 2 and 3: the location (μ) and scale (σ) parameters in the log and lognormal functions; and the shape (k), scale (β) and location (μ) parameters into the GEV distribution. Equation 4 was used to calculate the adjusted return periods (ARPs, ^{Chow et al., 1994}; ^{Woong et al., 2003}). To calculate the frequency of occurrence (ARPs, Equation 4) of extreme events for each indicator of MeT, with Equations 1, 2 and 3, implicit in the function F(x), the values of the parameters μ, σ, k, and β , (obtained with the XLstat version 2014 software) were aggregated, respectively.

The modeling tools of this study (ST, fitted probability distribution functions and ARP), associated with the indicators of MeT; HE (Fr, AAD, ADD and IN) and MmT, serve to estimate the frequency of MeT occurrence. This methodology has also been applied by ^{Chow et al. (1994)}, ^{Ariza (2013)}, ^{Gomáriz (2015)} and ^{Campos (2016)}. To calculate the ARPs, XLstat version 2014, Surfer version 10.0 and corelDRAW version X7 were used.

3. Results and Discussion

3.1 Mann-Kendall test for STs

Figure 2 shows the weather stations with warming or cooling STs (ST positive or ST negative, respectively) of the indicators of MeT; FR, AAD, ADD, IN and MmT. The number of negative and positive STs for indicators of MeT were FR 6 negative and 3 positive; AAD 5 negative and 3 positive; ADD 6 negative and 3 positive; IN 5 negative and 4 positive; and MmT 6 negative and 3 positive. The weather stations that showed negative STs in the five MeT indicators were JAI, OTA and ROS, and the highest value of ADD was recorded at the ROS station with ST=-3.60. Only one station showed a highest positive ST in the five MeT indicators; namely, GUT with ST=3.90. These results are consistent with those of ^{Llanes et al. (2015)} in that although these authors calculated trends for the annual average maximum temperature, they also found a negative ST at the JAI station and a positive ST at the GUA station. These results provide essential information to characterize the evolution of regional warming-cooling, which should be attributed, according to ^{Alghamdi and Harrington (2019)}, not only to local and regional changes (elevation, latitude, and proximity to water bodies), but also to synoptic scale changes (atmospheric circulation) such as the climate indices Pacific decadal oscillation, Atlantic multidecadal oscillation, oceanic El Niño index and North Atlantic oscillation (Llanes, 2023).

Figure 2 Interpolation of Mann-Kendall test (STs = black line, dimensionless). 

Sen's slope test for STs

Figure 3 shows the magnitudes of the ST indicators. The highest magnitudes of negative ST were recorded at stations EPL (FR = -2.81 HE dec^-1, AAD = -9.38 day HE^-1 dec^-1 and MmT = -1.88 °C day^-1 dec^-1), ACA (ADD = -20.94 day dec^-1) and BAD (IN = -12.92 °C dec^-1). The highest positive magnitudes of ST were recorded at stations CUL (FR = 3.44 HE dec^-1), GUT (AAD = 6.15 day HE^-1 dec^-1 and IN = 13.62 °C dec^-1), IXP (ADD = 35.00 day dec^-1) and POT (MmT = 2.50 °C day^-1 dec^-1), as shown in Figure 3. If Li's criteria are considered, which are applied in the presence of extreme droughts, the positive trends of STs om FR, AAD, ADD, IN, and MmT suggest that extreme drought events may occur near stations CUL, GUT, IXP and POT. In this situation, the strength of these events could affect the yield of priority rainfed crops grown in the study area; beans, sorghum (Sagarpa, 2015) and corn (^{Llanes et al., 2022b}). This is verified with the results obtained by Llanes et al. (2022b), who established that in Sinaloa, the indices of meteorological drought that require evapotranspiration and maximum temperature as input data for their calculation (i.e., reconnaissance drought index and effective reconnaissance drought index) registered the highest significant correlations with rainfed maize yield; i.e., in Las Tortugas (station with the highest average maximum temperature for the period 1982-2013). This could put the socioeconomic and agrifood sustainability of Sinaloa and the whole country at risk (^{Robeson, 2004}; ^{Ruiz et al., 2005}; ^{Sheridan & Lee, 2018}), especially considering that this is one of the most important agricultural states in Mexico (^{Rincón, 2012}). In summer, Sinaloa has extremely high temperatures (up to 48 °C) and is often affected by diseases caused by high MeTs (^{Riojas et al., 2006}). It is important to note that the positive values of ST for MmT suggests that the people near stations CUL, GUT, IXP and POT could be at risk of sporadic events of heatstroke and, as a consequence, potential loss of human life (OPS, 2019).

Figure 3 Interpolation of Sen's slope test (STs = black line). Measurement units: a (HE dec^-1), b (day HE^-1 dec^-1), c (day dec^-1), d (°C dec^-1) and e (°C day^-1 dec^-1). 

The results of the analysis of FR, AAD, ADD, IN and MmT presented in Table 3 show a prevailing cooling trend, due to the number of indicators of MeT with negative STs (28 stations) being greater in absolute value than those with positive values (15 stations).

Table 3 Location, scale, and shape parameters of the fitted probability distribution 

MeT-maximum extreme temperature, FR-frequency, AAD- annual average duration, ADD-annual daily duration, I-intensity, MmT-maximum maximorum temperature, Log-logarithmic, Lognormal-Log-normal, GEV-generalized extreme values, u-location, o-scale, k-shape and P-scale.

Source: ^{Llanes et al. (2021)}

These results are consistent with the historical evidence reported by ^{Llanes et al. (2015)}, where they found a higher number of weather stations with negative STs for the minimum temperature in the states adjacent to Sinaloa (Baja California Sur, Durango, Chihuahua, and Sonora). Of 38 weather stations studied, 20 showed negative STs and 18 positive STs. Llanes et al. (2015), found that in some stations near the state of Sinaloa there are negative trends for the minimum and average temperatures, but not for the maximum temperature, which presents a positive trend (i.e., in the stations Puerto San Carlos, El Rosario and Ojito de Camellones), which can be attributed to the high randomness between the spatial and temporal influence of the Pacific decadal oscillation, Atlantic multidecadal oscillation, oceanic El Niño index, North Atlantic index, and geopotential height at the level of 700 hPa (Llanes et al., 2018; Llanes et al., 2020; Llanes, 2023), and the magnitudes of the maximum, minimum and average temperatures, and the meteorological droughts of the northern states of Mexico (Llanes et al., 2018).

The difference between Figures 2 and 3 is that Figure 2 only identifies which indicators are statistically significant and Figure 3 shows the amount of significant trend.

Statistical analysis of the fitted probability distribution functions

Table 2 shows the p-values of the fitted probability distribution functions, where the maximum values correspond to the Log, GEV and Lognormal functions. The highest p-values of the MeT indicators were FR (Log, p=0.9944), AAD (GEV, p=0.9795), ADD (Log, p=0.9728), IN (Lognormal, p=0.9517 and MmT (Log, p=0.9088). The lowest p-values (p < 0.00001) corresponded to the fitted probability distribution functions Erlang, Fisher-Tippett (1), Gumbel, Normal-standard, Student and Weibull (1). If the log function were chosen for all the indicators (FR, AAD, ADD, IN and Mmt), then even the indicators AAD and IN (which do not present large p-values) can show good fit (p = 0.9573 and p = 0.9070 respectively).

These results are in agreement with ^{Aguilera (2007)} and with ^{Papacharalampous et al. (2018)}, who states that hydroclimatological data are rarely normally distributed; MeT indicators theoretically fit distributions that represent seasonal extreme values (^{Linsley et al., 1977}; ^{Pliego & Ruiz 2004}; ^{Campos 2014}; ^{Wang et al., 2016}; ^{Alghamdi & Harrington, 2019}).

Adjusted return periods

Table 3 shows the magnitudes of the shape, scale and location parameters of the fitted probability distribution functions used to calculate the adjusted return periods (ARPs). These parameters, based on the average values of MeT, are shown in Figure 4 for the ranges FR, 1.07 to 377.84 years; AAD, 1.03 to 24.88 years, ADD 1.06 to 816.94 years, IN, 1.02 to 16.59 years; and MmT, 1.03 to 118.15 years.

The P85 of MeTs as a function of the ARPs are FR (P85=3.65 HE year^-1; Figure 4a), AAD (P85= 5.09 day HE^-1 year^-1; Figure 4b),ADD (P85=19.34 day^-1 year^-1; Figure 4c), IN (P85=36.99 ^oC year^-1; Figure 4d), and MmT (P85=40.85 °C day^-1; Figure 4e).

Figure 4 Adjusted return periods (ARPs) of the maximum extreme temperature (MeT) indicators: hot extremes (FR, frequency; AAD, annual average duration; ADD, annual daily duration; IN, intensity) and maximum maximorum temperature (MmT). Fitted probability distribution functions: logarithmic (Log), Lognormal (Log-normal) and generalized extreme values (GEV).  

According to Figure 4, the average occurrence frequency of extreme events (values with the highest probability of exceedance) per 100 years are FR = 1.16 HE dec^-1 (ARP = 1.07 years; prob. = 0.93 = 93 times ), AAD = 2.24 day HE^-1 dec^-1 (ARP = 1.03 years; prob = 0.97 = 97 times), ADD = 4.16 day dec^-1 (ARP = 1.06 years; prob = 0.94 = 94 times), IN = 20.42 °C dec^-1 (ARP = 1.02 years; prob = 0.98 = 98 times) and MmT = 39.18 °C day^-1 dec^-1 (ARP = 1.03 years; prob = 0.98 = 98 times).

The average frequency of occurrence of extreme events (values with less probability of exceedance) per 100 years are FR = 6.11 HE dec^-1 (ARP = 377.84 years; prob. = 0.01 = 1 time), AAD = 6.64 day HE^-1 dec^-1 (ARP = 24.88 years; prob = 0.04 = 4 times), ADD = 38.68 day dec^-1 (ARP = 816.94 years; prob = 0.01 = 1 time), IN = 39.09 °C dec^-1 (ARP = 16.59 years; prob. = 0.06 = 6 times) and MmT = 41.95 °C day^-1 dec^-1 (ARP = 118.15 years; prob. = 0.01 = 1 time).

These results are especially useful for carrying out adaptation measures in sectors such as the agrifood industry, which is the main economic activity in the state of Sinaloa (^{Fiscal et al., 2017}).

If these results are not considered in the decision making for seeding and harvesting dates, significant reductions in crop yields could occur in the region, resulting in significant negative effects on the economy of the region, which can generate high costs (^{Ruiz et al., 2013}; ^{Llanes et al., 2022b}). These unusual temperature fluctuations (indicators of MeT) that will increase by the year 2030 (World Bank et al., 2014), could trigger significant decreases in agricultural crop yields in Sinaloa (Llanes et al., 2022b), which will aggravate the food sustainability of Sinaloa and of Mexico.

Through the results of this study, the association among the indicators of MeT can be observed; HE (FR, AAD, ADD and IN) and MmT and the modeling tools (STs, fitted probability distribution functions and ARPs), which satisfy the requirements for calculating estimates of the frequency of occurrence of MeT indicators.

Source: ^{Llanes et al. (2021)}

4. Conclusions

Based on the ST results of these parameters, a greater number of negative values were observed than positive, which shows a trend toward cooling. The fitted probability distribution functions used to calculate the indicators of MeT were for FR, Log; for AAD, GEV; for ADD, Log; for IN, Lognormal; and for MmT, Log. This contribution provides a useful methodological tool to design adaptation strategies for climate change and natural disasters that may affect the agrifood sector. With the strategies that can be devised, non-reduction of the agricultural yield should be guaranteed for the most important crops in Sinaloa; maize, sorghum, wheat, and beans. These products are a mainstay in feeding Mexico and a great part of the harvest comes from the study area, which is still called "the breadbasket of Mexico." Limitations of this study are the small number of weather stations studied and the small number of probability distribution functions used for the fit. It is recommended that the MeT indicators be recalculated with data from a greater number of meteorological stations, as well as using a greater number of probability distribution functions to try to increase the p-values and obtain more accurate ARP. Some proposals for climate-smart agricultural strategies (adaptation/prevention) are drip irrigation, protected agriculture, and crop substitution, which can serve to minimize the effects of high temperatures in the agricultural activities of Sinaloa state.