Investigating the inter-annual precipitation changes of Iran

Precipitation is the most important climate element and is usually investigated from the two perspectives of time and place. Emphasizing the perspective of time determines temporal change and emphasizing the perspective of place determines spatial change. In climatological studies, these two aspects are always interdependent (Alijani et al. 2012), because climate is an element with both spatial and temporal aspects. Therefore, understanding temporal and spatial variability is very important in environmental planning. Knowledge of the spatial and temporal distribution of precipitation is necessary for making precipitation forecasts and monitoring weather conditions (Chappell et al. 2013). Precipitation is one of the basic meteorological elements which affects human life directly and indirectly as an essential component of the global hydrological cycle and is a key parameter in ecology, hydrology and meteorology (Goovaerts 2000; Langella et al. 2010; Li & Shao 2010; Yatagai et al. 2012). The existence of observational data and its spatial and temporal distribution over vast areas is very important in many studies of water balance, flood forecasting, and precipitation (Fotovatikhah et al. 2018; Ghorbani et al. 2018; Moazenzadeh et al. 2018; Mosavi et al. 2018; Qasem et al. 2019; Homsi et al. 2020). Precipitation measurement stations which constantly monitor rainfall are often not available. These stations have a sporadic distribution, especially in mountains, and the close interaction of these components with the climate adds to the complexity of precipitation estimation, due to the complexity of topography of such areas (Oki et al. 1991; Barry 1992; Hofinger et al. 2000; Sturman & Wanner 2001; Sotillo et al. 2003). Several methods might be used to estimate precipitation in areas with no stations; among these methods geostatistical techniques should be noted.

The geostatistic method has been confirmed as an appropriate way to assess precipitation data (Delhomme 1978; Goovaerts 2000; Asakreh 2008). Geostatistics can be used specifically for precipitation variability (Barancourt et al. 1992; Berne et al. 2004). In this study, the spatial structure of Iran's spatio-temporal precipitation variability has been investigated using geostatistical techniques. An accurate estimation of the spatial precipitation distribution requires a dense and regular cellular network (Goovaerts 2000). A static spatial model must have a fixed mean, variance and direction across the study area (Fortin & Dale 2005). Considering the climate change debate and the temporal and spatial variations of precipitation, several researchers have conducted studies using parametric and non-parametric methods, including (Kahya & Kalaycı 2004; Zhang et al. 2006; Liang et al. 2010; Tang et al. 2011; ShiftehSome et al. 2012; Duhan & Pandey 2013). Jia et al. (2011) used spatial statistical methods to perform a fine scale view of an exponential algorithm of Tropical Rainfall Measuring Mission (TRMM) Precipitation Network, based on the Normalized Difference Vegetation Index (NDVI) and digital elevation model in China's cadmium basin. In this study, the Global Moran index was used to analyse the pattern of NDVI index. The results indicate that precipitation's fine scale viewing of data using the mentioned algorithms were different at six studied stations based on the data. Allard & Soubeyrand (2012) used practical spatial statistics and identified the areas vulnerable to climate change in research in the Colmar area in eastern France and epidemiology dispersion models of plant species. Robeson et al. (2014) analysed meteorological stations of the world in research using the K Ray ply function and have introduced the best pattern of this function. As temporal-spatial changes are among the most important issues of applied climatology, the main purpose of this study was to monitor spatio-temporal changes in Iran's precipitation at a monthly scale.

Iran is a rugged land, which emphasizes the variability of precipitation in the country. The need for studies on precipitation and its behaviour is strongly felt, due to the climatic characteristics of Iran and the importance of precipitation in water resource management. Therefore, a precipitation data investigation study would be useful, in which long-term precipitation data are available to the researchers. Accessing long-term precipitation data is not possible in Iran. There are currently more than 380 synoptic stations but only 50 or 60 of them have long-term (30 years or more) statistics of meteorological variables which is not sufficient to investigate the spatial autocorrelation precipitation with high spatial and temporal fluctuations given Iran's large area. In addition to the long-term data accessing problems, desert areas and internal deserts and areas with an altitude of over 2,600 m (country's pond centres) lack precipitation measurement stations. Further adjacent areas with high climate variability, such as the north-south slopes of Alborz and eastern and western slopes of Zagros, lack stations (Asakreh 2008). The data obtained from all available precipitation databases such as Global Precipitation Climatology Project (GPCP), Global Precipitation Climatology Centre (GPCC) and Asian Precipitation – Highly-Resolved Observational Data Integration Towards Evaluation of Water Resources (APHRODITE) were reviewed. On GPCP and GPCC databases, the rainfall data were available on a monthly scale, thus not showing extremes of rainfall events. In addition, the high spatial resolution of the APHRODITE database (0.25 by 0.25 degree spatial resolution) was compared to the GPCP database (2.5 by 2.5 degree spatial resolution). Finally its allocation was compared to that of Asia. Therefore, the outputs obtained from the APHRODITE database evaluated the temporal-spatial changes of Iran's seasonal precipitation. Figure 1 shows synoptic stations in the country and APHRODITE precipitation database with a spatial resolution of 0.25 × 0.25°.

A project entitled APHRODITE was started by the Human and Nature Research Center and Meteorological Research Center of Japan in order to create daily networked precipitation databases with high spatial resolution and based on rain gauges across Asia in 2006. This database covers a period of more than 57 years and has been created with the ADW (angular-distance weighting) interpolation method (Yatagai et al. 2012). The data of the APHRODITE stations with the spatial interval of 0.25 × 0.25 and 0.5 × 0.5 arc degree for the Middle East and South East Asia and a temporal interval from January 1957 until December 2007 exist at the following website (http://www.chikyu.ac.jp/precip/index.html) (Yatagai et al. 2012).

In the present research, the data available in the Middle East (APHRO_ME) were taken from the latest product of the APHRODITE database entitled v1101 with spatial resolution 0.25 × 0.25 and degree of arc 0.5 × 0.5 with the NETCDF obtained from the above website. Then, the programming capabilities of Grads and Matlab software packages were used to elicit the data from all existing data in the database (APHRO ME). The cells pertaining to the area of Iran were extracted using the capabilities of ArcGIS software. As the mentioned data exist at a spatial resolution 0.25 × 0.25, for the next stage, using geostatistical techniques and Kriging interpolation methods, the cells to the dimensions of 15 × 15 km² were spanned on the study region. The Kriging interpolation method refers to an optimal technique that provides regional observations in the regions without any value (Asakreh & Razmi 2012). Also, Phillips et al. (1992) showed that the Kriging geostatistical model is the best method to predict rainfall for areas without data. Instead of fulfilment of the analysis process on a matrix with data representing the spatial resolution 0.25 × 0.25 and acquisition of the results from the study area, the micro scale data were used to the dimensions 15 × 15 km. From this, boundaries of climate regions and spatial patterns were revealed.

To examine the prevailing pattern in precipitation in Iran during the past half-century, the modern spatial methods such as spatial auto-correlation (Local Moran and Global Moran) and hotspot were used. Recognition of patterns and the discovery of existing processes in spatial data are of importance (Diggle 2003; Waagepetersen & Schweder 2006), because, before any analysis of the maps of spatial statistics, we must engage in pre-judgement over how the data have been distributed in space and what patterns and rules are used for distribution of the data in the space (Illian et al. 2008). In the present research analysis methods of Moran, the local pattern has been used for a better understanding of data and accurate decision-making about the confidence level. In recent decades, a variety of scenarios concerning analysis of spatial data patterns have been expanded. For this the Global Moran index, also known as Global Moran I, is used. This index obtains numerical statistics (z-score) through which it can measure the degree of scattering in spatial data in space (Getis & Ord 1992; Levine 1996; Mitchell 2005; Wheeler 2007; Illian et al. 2008).

The Global Moran's I statistic provides a measure of the degree of spatial auto-correlation based on both the locations of events and the values associated with the events at the same time. It indicates the degree of spatial concentration or dispersion for a given point pattern (Scott & Janikas 2010). In general, if the value of the Moran index is close to +1, data will experience spatial auto-correlation, and if the amount of the index is close to −1, data will experience a discrete pattern. Cluster maps, a function mentioned in the analysis of pattern, imply that the general statistics which seek to give a response to the question ‘whether a significant spatial dispersion exists between the data or not’ help to identify where the clusters have been formed at the region under study and where are the borders in this region. In this study, cluster and outlier analyses (Anselin Local Moran's I) and hotspot analysis have been used to study patterns and spatio-temporal variations in precipitation. The cluster and outlier analysis which is well known as Anselin Local Moran's I is regarded as an optimal pattern to display statistical distribution of phenomena in space (Getis & Ord 1992; Anselin 1995; Getis & Aldstadt 2004; Anselin et al. 2009; Wheeler & Paéz 2009).

Considering that it is itself a kind of a Z score, Z's re-calculation has been avoided.

Table 1 shows some of the statistical characteristics of monthly precipitation. According to this table, the maximum precipitation mean of 57.88 mm occurred in January; the minimum of 3.30 mm occurred in August. However, the highest precipitation variation coefficient of 410.98 occurred in September. The distribution coefficient of skewness (G1) of all months' precipitation has been positive. A positive skewness indicates that the frequency of precipitation events is higher than the mean, and is higher than the frequency of precipitation events lower than the mean. In other words, the frequency of the data greater than the mean value is higher than the ones with a value lower than the mean. Rain-bearing systems, due to their dynamic and thermodynamic conditions and depending on their geographic locations, can cause different areas of rain when dealing with local conditions (Ghayour et al. 2011). Therefore, the amount of precipitation will have different indices. For every month of the year, the difference between the median, mean and standard deviation of the mean indicates that the data do not follow a normal distribution. As also specified in Table 1, the highest value of G2 (peaking coefficient) is related to September (98.32) and its lowest value is related to April (6.62). The peaking factor coefficient represents the spatial variation of precipitation in the country. The reason behind the higher peaking factor coefficient in September is due to the high spatial variation of precipitation in Iran at this time of year. The maximum rainfall, 258.59 mm, occurred in Anzali on the Caspian Sea coastline in northern Iran and the minimum rainfall, 0.00002 mm, occurred in Sarbaz city in south-eastern Iran. In winter, the difference in precipitation is less, due to rainfall in most parts of Iran. In some places, a small amount of rainfall will cause a relative regularity of precipitation rate in Iran and, to some extent, the low value of the peak coefficient.

Table 2 represents features of precipitation spatial statistics for different months. According to the ID index, September and October have shown the highest numerical values. This index is a good indicator of the country's precipitation spatial anomaly in a way that when Babolsar, Ramsarand Astara stations experience precipitation of 200 mm in the mentioned months, there is no record of precipitation in the eastern, southeastern and central parts of Iran. As shown in Table 2, the ICS index has been positive and above 10 for all 12 studied months which determines the cluster pattern of precipitation in the country. The GI index determines the same feature. The ICF index shows the size of the cluster. According to this index, the largest precipitation clusters are created during winter which indicates the relative order of the country's precipitation. The IP and GI indexes indicate the spatial difference of rain clusters, and summer and fall have the highest levels of cluster displacement.

Geostatistical methods will be optimal in cases where data are normally distributed, thus before performing these methods, examination of the histograms of the data is required to check the normality and identification of outlier data. In this research, data were obtained from the APHRODITE database. These data were quality controlled before being displayed on the website. The results of the Kolmogorov–Smirnov test showed that the H₀ theory is rejected with the minimum certainty of 0.86 compared to the H₁ theory. Therefore the normality of the data has been confirmed as well.

To reduce adverse effects on variogram estimation, the observations need to be removed. An overview of the process in the data must be examined through analysis of the process in geostatistical analysis. Then, by reducing the process at each cell with the observed value, the residuals which lack the process are used instead of actual data for fitness of variogram. To eliminate the trend, GS+ software, the de-trending tool was used, and as a result values of Global (0.46) and Local (0.53) were obtained. These values were selected utilizing the least root mean square deviation (RMSE) and testing various values. The scatter plot of cells before and after removal of the trend is represented in Figure 2. With the fitness of the process to scatter plot, anisotropy data were examined. The results indicate that the scatter plot demonstrates kurtosis in one direction, and with regard to one-directed kurtosis, the scatter plot will be conveyed with anisotropic geometry (Hassani Pak 2007). Anisotropy is interpreted in this way, that largeness of the area in which the data express spatial structure and associate together are different in various directions, regarding equal variation in various directions.

The size of grid is required to perform geostatistical methods (Masoodian 2008). The size of grid is under the influence of a variety of factors including scale, ability for computer processing, spatial accuracy, density of samples and spatial autocorrelation structure (Hengl 2006). In the present research, size and number of steps regarding the used method (Kriging) have been considered. The size of the steps has been found to have a major effect on the semivariogram. If the size of the steps is large, autocorrelation at short ranges is ignored.

If samples are taken as grids in this case, the distance between grids will be a suitable measure for the size of steps. Yet if the data are taken at random, the determination of step size is not simple; in this case the size of data must be considered equal to half of the maximum distance between all the samples (Johnston et al. 2003). To create an information layer at the station distance, a spatial analysis function has been used and maximum distance between stations has been calculated equal to 0.25°, thus the step size was calculated at 0.125°.

Outputs of Global Moran's precipitation spatial auto-correlation are listed in Table 3 and Figure 3. In general, if the Moran index is close to +1, data will have spatial auto-correlation and a clustering pattern, and if the Moran index is close to −1, data will be discrete.

Graphical output represents scattering or clustering data. According to the Global Moran index, the null hypothesis is based on the assumption that there is no spatial clustering between values of the element associated with geographical features. Yet when the P-value is so small and the calculated Z-value is so large, the null hypothesis can be rejected.

If the Moran index is greater than 0, the data will represent a type of spatial clustering. If the Moran index is less than 0, the geographical features under study will have a scattered pattern. As shown in Table 2, the Global Moran index for all 12 months is above 0.70. This point indicates that precipitation of the country during all the 12 months is at the 95 and 99% levels; this indicates a high clustering pattern. Nonetheless, the highest Global Moran index with a value of 0.975363 has been associated with January, which indicates the relative pattern in rainfall over Iran. Z-value for the 12 months is high, ranging from 193 to 263. Hence, in general, based on the Global Moran index, it can be deduced that inter-annual changes of precipitation over the country follow a high clustering pattern. Hence, with regard to high Z-value and low P-value, the null hypothesis can be rejected based on the lack of spatial auto-correlation between data for all 12 months. If spreading the precipitation in a normal way during various months over Iran was assumed, the Global Moran index will be equal to −0.000139.

As observed, the spatial auto-correlation of the Global Moran index just specifies the type of pattern. For this reason, a Local Moran index for the months under study has been used to display the spatial distribution pattern of precipitation in Iran. The results indicate whether the features have been distributed randomly, scattered or clustered in the space. If the I-value is positive, this will imply that the considered condition has been dominated by similar features. Hence, the considered feature is a part of that cluster. If the I-value is negative, this will imply that the feature has been enclosed by dissimilar features. These are called outliers. The value of these statistics has been calculated in the framework with the standard score, and P-value can be analysed. In these statistics, HH represents clusters with positive spatial auto-correlation at 99% confidence level, LL represents clusters with negative spatial auto-correlation at 99% confidence level, HL represents outliers in which a high value has been enclosed via low values, and LH represents single cells in which the condition enjoys a low value, enclosed via high values. In Figure 3, areas with dark black colour indicate high clusters (HH) and areas with a lighter colour indicate low clusters (LL) which are statistically significant at the 99% level and are distinct from each other. The grey areas are not statistically significant. In other words, in this region there is no dominant pattern and therefore it lacks any spatial autocorrelation pattern. This condition in January as representative of winter is seen in more than 50% of the country. Only about 20% of the total area of the country, mostly on the southern shores of the Caspian Sea and the Zagros belt, has a high cluster model (positive spatial autocorrelation) (Table 4). Low cluster models (negative spatial autocorrelation) exist in almost 23% of the country's area as models with low precipitation in central Iran, north of Sistan and the south-east of the country. In April, as the representative month of spring, low cluster model (negative spatial autocorrelation), high cluster model (positive spatial autocorrelation) and areas with no spatial autocorrelation are dominant over 40.62 and 84.62% of the country's areas, respectively. In November, as the representative month of autumn, the high cluster and low cluster models are dominant over 21.15% and 29.24% of the country, respectively, and areas with no spatial autocorrelation have been allocated 49.61% of the country's total area (see Table 4). As can be seen, top clusters (HH) cover the north coast (especially Babolsar to Anzali), with the centre and northwest of the Zagros as islands. In just one month (July) in the southern part of the country, a top cluster model can be seen which is due to the impact of the Indian Ocean monsoon precipitation in this region of the country. Due to the shift in Tibet's high-pressure location to the west in July, these monsoons create favourable synoptic weather conditions that allow the wet Indian monsoon air system to enter Iran from the southeast and intensify rainfall in these areas. The results of this section are consistent with those of Shukla & Misra (1977) which showed that with the increase in sea surface temperature (SST) in the Arabian Sea during July, rainfall will increase significantly in the western India Ocean (Indian monsoons), including south-eastern Iran. The lower cluster areas are mostly in the central, south-eastern, southern coastal regions of Iran and the coasts of the Oman Sea to the Strait of Hormuz, due to the location of the Zagros Mountains in the western strip (northwest to southeast direction) and the Alborz Mountains in the northern strip (west-eastern direction), which prevents the entry of moist air masses into central and south-eastern Iran, thus creating low-rainfall areas.

With regard to what has been mentioned above, those areas of the country with positive and negative spatial auto-correlations were specified, and the hotspots index or GI* used from the areas with clusters with high and low values with the results displayed in Figures 4 and 5 and listed in Table 5. The GI* statistic, which is calculated for existing conditions in the data, is a Z-value. For a positive value, the higher the score, the greater the degree of clustering showing a hotspot. For a negative score, the smaller value l, there will be fewer clusters, indicating cold spots. For the group of Z-scores with statistically significant negative values, a smaller Z-score indicates a tighter clustering of low values which indicates areas with low precipitation. On this basis, high and mountainous parts of Iran include most parts of East and West Azerbaijan, Kurdistan, Zanjan and Qazvin, Bijar, the West Highlands, Khansar, Golpaygan, Tuyserkan, Isfahan, Arak, Hamedan, Kermanshah as well as the Alborz mountain areas from Taleghan to Alamut and Damavand, Northern Khorasan, Khorasan Razavi and northern parts of Southern Khorasan which are under the dominance of heavy precipitation patterns. According to Figures 4 and 5, the south-eastern and south-western areas of the country are dominated by low precipitation patterns. These areas include most parts of the Dasht Yari coastal plain of the Strait of Hormuz, central Iran, the hinterlands and the coastal plain of the Persian Gulf. Table 5 shows the percentage of areas under the influence of precipitation patterns in each month of the year at 90, 95 and 99% confidence levels. As can be seen in Figures 4 and 5, the effects of altitude can be seen quite well in the precipitation patterns. Adjacent areas that are of ecological divergence like the north-south and slopes of the Alborz and eastern-western slopes of Zagros (Asakreh 2008) in the Western System (cold periods) are displayed in strips across the northwest and southwest of Iran. At the same time as the sub-tropical high pressure (STHP) affects Iran during the warm period of the year, the spatial precipitation anomalies increase; a typical example of the effect of this system can be seen in September. Moreover, with the arrival of this system in Iran, dry conditions prevail throughout the country, but the precipitation of the Caspian Sea coast in northern Iran may continue because of the separate system of precipitation in this region compared to the central parts of Iran. The Alborz Mountains prevent STHP entering the Caspian region. Therefore, due to the presence of moist resources in the coastal areas of the Caspian Sea, in late summer (natural convection) and early autumn (horizontal convection) separate precipitation systems are formed in this region of Iran and create high precipitation patterns (positive spatial auto-correlation).

Limited movement of convective precipitation in the Caspian region is due to the existence of high Alborz Mountains in the south of this region, as they prevent the entry of rainfall into the central part of Iran and the only displacement in the mountain belt is possible during convective precipitation systems (Figure 5).

The findings of the present research support the work of Asakereh & Seifipour (2013), which indicated that the Zagros Mountains, the northwest mountains and the Caspian regions are considered as the country's high precipitation regions. Nazemosadat et al. (1995, 2006) and Tabari & Hosseinzadeh Talaee (2011) indicated also that the country's northern parts and the Caspian Sea coasts have precipitation consistency but the southern parts have precipitation dispersion.

Iran has special precipitation circumstances due to its wide range of latitude and longitude, the existence of ripples' patterns and exposure to the prevailing air mass. The general structure of precipitation in Iran is influenced by latitude, longitude and air masses in a way that with the change of each of these factors precipitation will change as well. In other words, the general conditions of precipitation are a function of latitude and altitude, and other factors such as water and surface characteristics which are referred to as local factors, are involved in the formation of Iran's precipitation. Spatio-temporal analysis has been performed in the present study using modern methods of spatial statistics. In this context, Global and Local Moran, hotspots and cluster and non-cluster analysis methods have been used. The present study has concentrated on the assumption that Iran's precipitation is a function of the cluster model and the pattern of precipitation itself is also a function of the internal and external conditions. To achieve this goal, the APHRODITE precipitation database of monthly mean for the period 1977–2007 has been used. According to the statistical outputs, September has the highest precipitation coefficient of variation (410.98) in Iran. The peaking factor and Gates indexes also confirm this result. The highest precipitation anomalies based on spatial distribution belong to September (75.26) and October (45.02). Cluster size index outputs also indicate that Iran's largest precipitation clusters are created during winter, which indicate the relative order of Iranian precipitation. The results of the Global Moran method showed that inter-annual changes in rainfall over Iran follow the high cluster pattern. Accordingly, in September, Iran has the highest difference in precipitation patterns with the lowest numerical index value of 0.712121 (between +1 and −1) and January showed the highest Global Moran index value of 0.973363. The reason for this is the entry of western systems into Iran and its complete domination over the country. This system, which is controlled by the Mediterranean trough, provides most of Iran's instability and winter rainfall during the cold period of the year (January, February and March) by conveying cyclones and wet Mediterranean air masses.

The results of the Global Moran method indicated that yearly precipitation changes in Iran follow the above cluster model. On this basis, September with the lowest index number (0.712114) of the Moran statistics (between +1 and −1) has the highest difference in Iran's precipitation patterns and January, with the highest Global Moran index of 0.975363 shown as well. The reason for this is the westerly weather system entering Iran and having full control over the whole country. This system is controlled by the Mediterranean Basin in the cold period of the year (January, February and March) and causes winter precipitation and instability in most parts through cyclonic activity and moist air masses. Based on the analysis of spatial autocorrelation, areas with negative spatial autocorrelation are related to the south-east and coast of Oman to Abadan and parts of the north-east of the country throughout the year. Areas with positive spatial autocorrelation are mostly in the southern confines of the Caspian Sea and Zagros. In all 12 months of the year, a large area of the country lacks significant spatial autocorrelation. Areas with positive spatial autocorrelation, or in other words, areas of Iran with high precipitation, from June to August in south and south-eastern Iran, examples of which are Sarbaz, Sarvan, Iranshahr and Nikshahr in Sistan and Baluchestan, Haji Abad in Hormozgan and Lar in Fars province, indicate the effect of the warm humid tropical air in these regions which enters Iran through the following ways: (1) Oman and the Persian Gulf breeze (to a limited extent) impact the warm season precipitation that is to be expected in the Fars province; (2) seasonal weather which creates Iran's southern and south-eastern precipitation systems through low thermal stress which is created in Pakistan and India in summer. However, it should be emphasized that the impact of this system is discovered when the uplift factor is also available. To sum up, spatial analysis has shown that Iran's precipitation models are classified into two forms of rainfall patterns in the south region (LL, low precipitation pattern), northern region (the Caspian Sea shores) and west and northwest (HH, high precipitation pattern). The investigation indicated that in the study period, low precipitation models (negative spatial autocorrelation) have had a much higher frequency compared to high precipitation models. To sum up, given that Iran is a rugged land, the variability of rainfall is very high and the management of water resources is strongly felt. The main problem in managing water resources is changes in rainfall, which were difficult to detect due to the lack of sufficient rain gauges and synoptic stations over the country. Therefore, in this study, by using the APHRODITE database, high and low precipitation patterns and clusters were identified within the country. This can be used in water resources management such as supplying agricultural water resources, dam construction, etc. Therefore, in future research, different regions of Iran can be examined in terms of the cultivation of various agricultural products (such as the cultivation of drought-tolerant plants in low-rainfall areas and vice versa), location of dam construction and so on.

The authors would like to thank Peter Smithson, formerly Senior Lecturer in Climatology, University of Sheffield, UK, for the revision of the article from a linguistic and scientific point of view. The authors would also like to acknowledge the University of Isfahan for providing access to their data.

All relevant data are included in the paper or its Supplementary Information.