Regional drought assessment using improved precipitation records under auxiliary information

Abstract

Changes in the climate and weather conditions, as well as rising earth’s average temperature are likely to escalate deterioration of global drought occurrence. Drought is considered an interwoven natural disaster composed by a number of different factors, as for example agricultural, meteorological or hydrological. Hydrological drought estimation with regional accuracy is the most problematic and challenging issue. In order to monitor and characterize drought conditions, using Standardized Drought Indices (SDI) is recently the most frequently used practice. In this research article, we suggest an improved hydrological drought index that incorporates upgraded monthly rainfall estimation records, which play an important role in defining regional drought conditions, with regard to the global temperature rise. Rainfall is highly changeable even at a low distance and therefore should be also considered in precipitation estimation records because temporal rainfall records play a significant role in determining long-term rainfall shortages. Thus, the integration of regional aspect to the amount of rainfall is essential for accurate regional drought assessment. This research article proposes adding auxiliary data such as regional weights in order to make monthly rainfall records more accurate in relation to the dependency characteristics of temperature and rainfall records under regression and product estimation settings. Subsequently, we propose an innovative method of hydrological drought evaluation, a so-called Regionally Improved Weighted Standardized Drought Index (RIWSDI). We evaluated hydrological drought with the usage of RIWSDI at seven various meteorological regions situated in climatologically different areas in Pakistan. We assessed and compared the results using RIWSDI, Standardized Precipitation Index (SPI) on 3 and 12-month interval period on the basis of Pearson correlation. Under both parametric and non-parametric standardization, we discovered that there is a high positive correlation between RIWSDI and current methodology (SPI). To sum up, we proved that the upgraded estimations of rainfall are able to improve systems for monitoring droughts.

Keywords:

• drought
• standardized precipitation index
• drought indices
• regional drought
• auxiliary information

1. Introduction

Drought is a complicated natural phenomenon that occurs basically due to the lack of precipitation over a time period. And its occurrence usually results in great costs on various parts of nature and society. Drought is usually believed as a prolonged lack of precipitation. A worldwide definition based on precipitation quantity and period is: A drought implies that a spot has less precipitation (rain or snow) than ordinary over a couple of months or significantly more or drought is usually a situation of adversity because of the insufficient water due to surprising meteorological conditions. It could be categorized into disparate types such as; Meteorological, Hydrological, Agricultural and Socio-economic. Meteorological drought is related to precipitation for assessing dryness and dry period for a specific region where average precipitation may diverge spatially. Hydrological drought occurs when stream flow, reservoir, soil moisture, groundwater recharge, and lake levels are affected due to a decrease in precipitation. Agricultural drought occurs when there is not sufficient moisture to maintain average crop production for specific land area. Socioeconomic drought relates to the supply and demand of some economic goods with fundamentals of hydrological, meteorological, and agricultural drought. Climate and weather conditions in Pakistan differ across its large territory primarily due to the variety of altitude, longitude and latitude, wind flows and distance from the sea. A number of climate classification schemes have been introduced to account for spatial variability. One of the distinguished is the Köppen-Geiger scheme (Kottek and Rubel, 2007; Peel et al., 2007), which on the basis of precipitation and temperature separates climate into five classes, such as polar and alpine, continental, mild, dry and tropical climate. Since the late 20th century, a trend of so-called global warming can be observed around the world (Hansen et al., 2010; Rohde et al., 2013). Continuance of the trend is expected by numerous climate models under most emission scenarios and concentration pathways, which calculate the rise of the surface temperature between 0.7 and 2.4 °C by the year 2050 in comparison with the reference period of 1986–2005 (Collins et al., 2013). Field et al. (2012) assumes that the rise of overall temperature might cause unprecedented climate disasters (such as droughts and floods). Few severe droughts occurring in the early 21st century are supporting this presumption, as for example the 2010 drought in Australia (Cai et al., 2014), the 2011 drought in southern China Sun and Yang (2012) and the 2011–2012 drought in the USA (Seneviratne et al., 2017).

Several further environmental factors are engaged in the occurrence of drought, such as temperature, wind flows, relative humidity and intensity, duration and severity of rain (Wilhite et al., 1994). However, the long term precipitation and temperature factors play the leading role in methods for calculating drought indices (Coffel and Horton, 2015). In recent decades, the frequency, intensity and influencing area of the drought have significantly increased, which have drawn the attention of many researchers. It was mainly caused by climate changes and human activities.

Not only precipitation and temperature on a global level but also their regional distribution is of high importance for accurate monitoring of climate changes and other natural disasters. Furthermore, the documentation of rainfall and its regional distribution is one of the main responsibilities of the governmental meteorological department. Therefore, it is essential to track drought conditions with the use of regional statistics and data. Moreover, accurate regional precipitation estimates are indispensable for a wide range of research fields, such as hydrology, meteorology and others. Understanding of drought characteristics on a regional level is necessary for mitigating drought risks, moderating potential impacts on various socio-economic sectors and adopting proper measures and strategies (Hirabayashi et al., 2008; Svoboda et al., 2016).

Spatial distribution of precipitation and temperature belong to the major watershed factors playing a significant role in advanced hydrological researches. Consequently, taking into account the intricacy of the temporal structure of the regional climate, several authors Coles and Tawn (1991), Guler et al. (2007) and Mahdian et al. (2009) jointly constructed two methods based on geospatial tools and advanced statistical models. Yet, those techniques are built on temporal data collected from single stations, which means that it only covers a single realization at a continuous spatial domain. That deprives the findings of those methods from the effect of spatial prevalence in climate. Additionally, increased the uncertainty of the prediction might have negative consequences on climate shifts policies and reliability of environmental, climate and weather conditions prognoses. Therefore, the incorporation of regional rainfall records combined with the temperature already at the entry stage might provide a substantial contribution to the accuracy, efficiency, and reliability of drought moderation policies. A regional rainfall may be defined as ‘the average of all stations’ rainfall (monthly/annually) in a region’.

However, contemporary advances in estimation methods which add auxiliary variables are accessible in sampling theory and environmental statistics. Cochran (2007) proposes comprehensive theory and methods connected to the usage of auxiliary data in order to improve the assessment of unknown features of unsystematic variables. Several other researchers also incorporate additional supplementary data to enhance the estimation of an examined variable under auxiliary information based statistical methods such as regression and kriging (Zhu and Lin, 2010). Paloscia et al. (2013) works with Australian Landsat images having Landsat (path = 92/row = 84) that are similar to ENVISAT data in order to observe the effects of auxiliary information on vegetation. Apaydin et al. (2011) deals with altitude as a basis of auxiliary data connected to climatic inputs (precipitation and temperature) for interpolation under co-kriging settings.

Keeping in mind the significance of auxiliary information in an estimation procedure, the aim of this research article is to incorporate and apply regional temperature as auxiliary information in the process. By regional temperature we mean the average of all stations’ monthly mean temperature of a region. We propose another drought index, its framework is constructed on the basis of regional temperature used as an auxiliary variable to improve the accuracy and reliability of average precipitation estimates. Calculated precipitation values are further utilized to acquire Standardized Drought Index (SDI) values (Ali et al., 2019). In contrast, temperature and rainfall can be perceived as global representatives of a specific natural homogenous optimized regional catchment area. Thus, using regional temperature as the auxiliary information might improve the validity of precipitation data.

Although hydrological characteristics of droughts have been examined broader than ecological and socio-economic aspects, deeper understanding of hydrological aspects of drought is still necessary to assess potential effects of drought more precisely and to design and implement adequate and effective moderation measures. In this regard, further research on drought conditions taking into account hydrological and regional disparity factors seems crucial. Indeed, regional spatial variability of drought affects the definition of the drought itself. As Yevjevich (1967) indicates, the term ‘drought’ can be used for a large covered area together with a long term severe deficit of water, whereas expressions as ‘water shortage’ or ‘deficit’ describe occasions with less harmful impacts. Ali et al. (2019) proposed Locally Weighted Standardized Drought Index (LWSDI). They utilized auxiliary information as a local weight for the improvement of monthly precipitation record. Ali et al. (2019) used regression estimator as a tool to incorporate auxiliary information in the study variable. One disadvantage is that they didn’t cover the scenario for the negative correlation between actual rainfall and average temperature (auxiliary variable). And the other disadvantage is that they improved only the specified stations but not all the regions of Pakistan. Our study is designed to cover all the disadvantages in a precise manner.

Moreover, growing vulnerability and dependence of a high-quality agricultural production on the accessibility of water and the development of large-scale multi-purpose water-supply systems indicate the fact, that the analysis of drought only at a single station is insufficient and a regional scale would be more beneficial. Lastly, the necessity for working with the whole array of existing regional data in measuring stations has arisen and it can be similarly used in different kinds of hydrological researches (e.g. flood analysis) as well. Historically it has been recorded that almost all Pakistani provinces faced repeated droughts: Khyber Pakhtunkhwa (KPK) in 1902 and 1951; Punjab province in 1899, 1920 and 1935 and Sindh province in 1871, 1881, 1899, 1931, 1947 and 1999 (Drought Bulletin of Pakistan, 2015). Therefore, it is necessary to assess the regional impact of drought on a particular area in order to better understand the situation and subsequently adopt appropriate measures.

There are numerous number of estimation procedures such as product, regression and ratio estimators, which integrate complementary information on single/multiple auxiliary variables. Detailed overview, as well as the mathematical structure of regression, product and ratio estimators can be found in (Cochran, 2007). In each respective method, validation of the proposed estimator depends on the incorporation of the auxiliary variables. For instance, when it is perceived that there is a negative correlation between study variable and auxiliary variable, then product estimator provides precise estimates for the population characteristics. On the contrary, when there is positive and perfect correlation then ratio estimator is useful. Preliminary applications based on the auxiliary information employed in 1973 Tarima and Pavlov (2006), when Pugachev (1973) incorporated the auxiliary information using correlation effects. As, in many surveys and records keeping modules, a collection of some extra information related to study variables are common practice. Ali et al. (2019) proposed a new method for the assessment of drought-the Locally Weighted Standardized Precipitation Index (LWSDI). They used LWSDI on ten different stations having different regions of Pakistan. Different probability distributions have been used to calculate LWSDI. For each stations, the CDF of those distributions having smallest value of BIC, are then selected for standardization.

This research article suggests a new way to characterize annual drought conditions, in which precipitation data specified by the auxiliary variable are incorporated in Standardized Drought Index (SDI) process. Hereby we recommend adding a regional average of monthly temperature data as weights in order to improve dependent precipitation estimates under regression and product estimation settings. Thus, we are able to diminish the sampling defects in the projected rainfall quantity records and to take account of the global warming effect while monitoring the drought conditions.

This paper is organized as follows: Section 2 consists of materials and methods which includes the overview of the data, study area as well as methodology which provides a comprehensive exposition of the theoretical details behind the product and regression estimation, SPI, and RIWSDI. We also present the mathematical formula employing average temperature as auxiliary information under regional settings. Here we proposed a drought index employing improved precipitation estimates under parametric and non-parametric approaches. In Section 3, we will provide a detailed discussion and results on how we improved the precipitation estimates and used in SDI procedure. Here, some temporal plots showing the difference between simple and improved records will also be presented. Finally, a conclusion on the proposed method will be provided in Section 4.

2. Materials and methods

2.1. Data and study area

Pakistan is located in Southern Asia which is the junction of Middle East and Central Asia having 23°–37° N latitude and 61°–78°E longitude (Ahmed et al., 2018). It hosts the triple point (junction) of three world’s famous mountain ranges Karakoram, Hindukush and Himalayas in its north. Pakistan lies between 0 to 8611 m altitude. Mean temperature as well as precipitation data are considered as the climatic inputs in the current study. The highest precipitation (1038.6 mm) was recorded in Islamabad in July, 2001, and minimum precipitation record of (0 mm) in several regions were also observed during the study. There are both anthropogenic and natural reasons of change in the climate but earlier is the most leading with a constantly increasing trend since 1940 (the industrial revolution) (Anwar, 2011). Pakistan shares her border with four countries; in the west with Afghanistan and Iran, in the east with India, in the north with China, and the Arabian sea is located in the south of Pakistan. According to 2016 census, Pakistan has a total population of more than 200 million. Where typically most of the people concerned with agriculture sectors either directly or indirectly.

This research consists of seven meteorological regions, and these seven regions contains overall 50 stations, situated in different climatic regions of Pakistan. Figure 1 shows the chosen meteorological regions. Detailed clustering of the regions of Pakistan can be seen (Hussain et al., 2011). The map is generated through Geographic Information System (GIS). The proposed method, required a time series data of long-term monthly precipitations, as well as maximum and minimum temperatures. Consequently, a data ranges from January 1967 to December 2016 is obtained from – Karachi Data Processing Centre (KDPC) – via Pakistan Meteorological Department (PMD). The regions have significant high variability in rainfall and temperature during different seasons. Table 1 reveals the statistics of the selected stations for 50 years (1967–2016).

Fig. 1. Map of the selected regions.

Table 1. Statistics of the selected stations.

			Precipitation (mm)			Temperature (°C)			Coordinates
Regions		Stations	Mean	Min	Max	Mean	Min	Max	Lat (N)	Lon (E)	Altitudes (m)
		Astor	39.4	21.6	72.8	9.9	8.3	11.3	35.3570°	74.8624°	2546
		Bunji	13.3	5.1	28.3	17.6	15.9	19.4	35.6431°	74.6342°	1532
C1		Chilas	15.8	3.0	47.4	20.3	19.0	21.4	35.4222°	74.0946°	1265
		Gilgit	11.7	3.4	22.3	15.8	14.6	16.9	35.8819°	74.4643°	1500
		Gupis	15.9	0.4	56.3	12.6	10.9	14.1	36.2274°	73.4421°	2176
		Skardu	19.1	5.7	41.3	11.8	10.1	13.2	35.3247°	75.5510°	2228
		Cherat	52.1	16.0	109.4	17.2	14.8	18.9	33.8215°	71.8883°	892
		Chitral	37.9	17.8	64.9	16.0	14.7	17.4	35.7699°	71.7741°	1494
		Dir	115.6	73.6	179.1	15.5	14.7	18.9	35.1977°	71.8749°	1420
		Drosh	47.8	25.1	79.3	17.6	15.3	19.6	35.5684°	71.8038°	1359
C2		Kohat	47.8	19.5	83.3	23.3	22.0	24.4	33.5889°	71.4429°	489
		Parachinar	72.2	44.2	184.8	14.4	10.8	16.3	33.9011°	70.0860°	1705
		Peshawar	40.2	15.8	75.4	22.9	21.8	24.0	34.0151°	71.5249°	331
		Risalpur	55.9	20.6	89.8	22.2	21.2	23.1	34.0751°	71.9876°	309
		Zhob	23.5	9.1	41.3	19.3	17.2	20.8	31.3497°	69.4665°	1426
		Balakot	131.5	89.1	207.4	18.6	17.4	19.8	34.5482°	73.3532°	3212
		Garhidupatta	124.2	77.7	175.9	19.2	17.4	20.6	34.2264°	73.6157°	819
		Islamabad	120.6	56.3	192.1	21.7	20.6	23.4	33.6844°	73.0479°	540
C3		Kakul	111.2	78.3	146.6	16.9	15.2	18.7	34.1875°	73.2618°	1300
		Kotli	103.9	58.9	160.8	21.9	20.1	22.9	33.5008°	73.9007°	707
		Murree	145.1	103.4	202.8	13.0	11.2	14.8	33.9070°	73.3943°	2291
		Muzaffarabad	124.1	76.5	179.4	20.6	19.2	21.9	34.3551°	73.4769°	737
		Bahawalnagar	19.9	4.8	43.6	25.1	23.8	26.6	30.0025°	73.2412°	163
		Bahawalpur	15.2	0.9	55.9	25.7	24.6	26.8	29.3544°	71.6911°	214
		DI Khan	26.1	11.6	63.0	23.7	22.4	24.8	31.8626°	70.9019°	165
		Faisalabad	32.5	14.4	67.2	23.5	22.3	24.7	31.4504°	73.1350°	184
C4		Lahore	55.3	27.8	103.9	24.0	22.7	25.0	31.5204°	74.3587°	217
		Mianwali	45.6	11.4	90.5	23.7	22.8	24.7	32.6645°	71.4774°	210
		Multan	17.9	6.9	42.8	24.9	23.8	26.9	30.1575°	71.5249°	122
		Sargodha	39.9	21.1	63.9	23.8	22.8	25.0	32.0740°	72.6861°	190
		Sialkot	83.1	45.4	157.2	22.4	21.1	23.5	32.4945°	74.5229°	256
		Badin	18.6	0.0	76.2	26.8	25.8	27.8	24.6459°	68.8467°	10
		Chhor	19.1	0.4	45.8	26.5	25.4	27.6	25.5114°	69.7823°	4
		Hyderabad	13.2	0.4	43.7	27.7	26.8	28.9	25.3960°	68.3578°	13
		Jacobabad	11.2	1.2	48.6	27.2	26.2	28.1	28.2823°	68.4472°	56
C5		Khanpur	10.5	0.7	32.9	25.4	23.4	26.7	28.6332°	70.6574°	80
		Nawabshah	12.3	0.0	55.3	26.9	25.9	28.2	26.2447°	68.3935°	29
		Padidan	10.1	0.2	45.6	26.5	25.4	27.8	26.7724°	68.2922°	46
		Rohri	9.1	0.2	37.7	29.6	27.1	30.7	27.6687°	68.8943°	62
		Sibbi	14.4	0.0	31.7	30.3	26.9	31.9	29.5532°	67.8808°	130
		Dalbandin	6.7	0.3	17.0	22.6	20.4	23.9	28.8854°	64.3964°	843
		Kalat	15.3	0.0	81.4	14.6	12.9	16.6	29.0523°	66.5879°	2007
		Khuzdar	21.3	4.4	49.6	21.8	19.4	23.3	27.8165°	66.6057°	1237
C6		Lasbela	15.6	0.7	61.6	26.6	25.3	28.0	25.8700°	66.7129°	241
		Nokkundi	2.9	0.0	15.6	24.5	22.6	26.1	28.8259°	62.7500°	679
		Panjgur	8.3	1.4	25.4	22.4	20.4	23.7	26.7303°	64.1478°	980
		Quetta	21.4	5.2	79.2	16.8	14.4	18.9	30.1798°	66.9750°	1679
		Jiwani	8.6	0.3	32.2	25.7	24.1	26.9	25.0538°	61.7707°	57
C7		Karachi	16.3	0.0	59.4	26.6	25.6	27.7	24.8607°	67.0011°	8
		Pasni	8.3	0.0	29.7	25.3	22.6	26.5	25.2510°	63.4154°	10

2.2. Methodology

In order to assess and compute the comparative ascertainment of drought indices based on RIWSDI, secondary time series data of average temperature and monthly total precipitation (1967–2016) is used.

In the current study, we evaluate 32 probability distributions using an R package i.e. Propagate (Tellinghuisen and Spiess, 2014). Here, various types of fitness criterion e.g. Anderson-Darling, Kolmogorov-Smirnov and Chi-Square tests were used in order to figure out the candidate distribution among others, for each individual indicator. For standardization of each study region and all the stations included in the regions, a CDF of the distributions, having a minimum value of BIC, are subsequently selected. To portrait different results in tables some regions has been chosen randomly i.e. Region-1,3,6,7. In this research, we are handling the problem of updating and improving regional precipitation estimates. A step-by-step procedure can be observed from the flow chart by analyzing Fig. 2. Here, the average temperature is suggested as an auxiliary information. In previous work, various surveys indicate that there is a positive correlation between rain and temperature. Zhao and Khalil (1993) examine the relationship between precipitation and temperature for eight regions including the USA. Their survey indicates that there is a positive correlation between these variables in all the seasons. At the Guliya ice core, detailed analyses of the precipitation index (glacier accumulation) and the temperature proxy recorded in since 300 years BP show that precipitation correlates with temperature in this region (Yang et al., 2006). Rajeevan et al. (1998) found that temperature and rainfall were positively correlated during January and May but negatively correlated during July. Sneva (1977) found a positive month-wise correlation between temperature and rainfall in southeastern Oregon.

Fig. 2. Flow chart. Precipitation (P), Average Temperature (AT), Correlation coefficient (r), Product Estimator (PE), Regression Estimator (RE), Improved Precipitation (IP), Probability Distributions (PDs).

As the temperature is a globally representative environmental variable and homogeneous in nature has a strong association with precipitation. Thus, the use of temperature as auxiliary data is logically valid. Here the auxiliary variable is the mean monthly temperature. The following equations were used to incorporate auxiliary information, (1) $${\overline{y}}_r={\overline{y}}_i+b_1({\overline{X}}_j-{\overline{x}}_i)$$(1) and, (2) $${\overline{y}}_p={\overline{y}}_i\left(\frac{{\overline{x}}_i}{{\overline{X}}_i}\right)$$(2)

In Eqs. (1)(1) $${\overline{y}}_r={\overline{y}}_i+b_1({\overline{X}}_j-{\overline{x}}_i)$$(1) and (2)(2) $${\overline{y}}_p={\overline{y}}_i\left(\frac{{\overline{x}}_i}{{\overline{X}}_i}\right)$$(2) ${\overline{y}}_r$ and ${\overline{y}}_p$ are the updated regression and product means of the study variable respectively, ${\overline{X}}_j$ is the overall mean auxiliary variable of $$ station, ${\overline{x}}_i$ is the sample mean of the $i\text{th}$ month of the auxiliary variable. In Eq. (1)(1) $${\overline{y}}_r={\overline{y}}_i+b_1({\overline{X}}_j-{\overline{x}}_i)$$(1) $b_1$ is the regression slope between the study variable and auxiliary variable. After assessing theoretical support about the positive and negative correlations between precipitation and average temperature, this study suggests average temperature as a piece of auxiliary information to improve annual meteorological records of precipitation. So, before defining drought characteristics and precipitation deficient, we utilized the concepts of regression and product estimator to improve the annual estimates of precipitation using average temperature as auxiliary information. The mathematical structures of regression and product estimator employing monthly mean temperature as auxiliary information for the estimation of the total monthly amount of precipitation are as follows, (3) $$IW{P}_r={\overline{p}}_i+b_1({\overline{T}}_j-{\overline{t}}_i)$$(3) and, (4) $$IW{P}_p={\overline{p}}_i\left(\frac{{\overline{t}}_i}{{\overline{T}}_j}\right)$$(4) where, $IW{P}_r$ and $IW{P}_p$ are the improved weighted precipitation estimates of mean using regression and product estimator respectively, for the improvement of original series of precipitation, where, simple precipitation records are weighted by the dependence characteristics of average temperature with precipitation. From Eqs. (3)(3) $$IW{P}_r={\overline{p}}_i+b_1({\overline{T}}_j-{\overline{t}}_i)$$(3) and (4)(4) $$IW{P}_p={\overline{p}}_i\left(\frac{{\overline{t}}_i}{{\overline{T}}_j}\right)$$(4) ${\overline{p}}_i$ and ${\overline{t}}_i$ is the mean precipitation and mean of average temperature of the $i\text{th}$ region, respectively, whereas ${\overline{T}}_j$ is the overall mean of the average temperature of a region for the $j\text{th}$ month. In Eq. (3)(3) $$IW{P}_r={\overline{p}}_i+b_1({\overline{T}}_j-{\overline{t}}_i)$$(3) b₁ is the regression slope between the study variable (precipitation) and auxiliary variable (average temperature).

2.2.1. Comparative statistics and quality measures

In this study, we encompass a well-known correlation statistic ‘Pearson Product-Moment Correlation coefficient r’, commonly called the correlation coefficient, for the comparison of the outcomes of our proposed index with that of existing indices. Correlation coefficient r is most widely used test statistic and measures the collinearity between two series. Formula for r is in Eq. (5)(5) $$r=\frac{\sum \left(x-\overline{x}\right)\left(y-\overline{y}\right)}{\sum {\left(x-\overline{x}\right)}^2{\left(y-\overline{y}\right)}^2}$$(5) (5) $$r=\frac{\sum \left(x-\overline{x}\right)\left(y-\overline{y}\right)}{\sum {\left(x-\overline{x}\right)}^2{\left(y-\overline{y}\right)}^2}$$(5) where, in Eq. (5)(5) $$r=\frac{\sum \left(x-\overline{x}\right)\left(y-\overline{y}\right)}{\sum {\left(x-\overline{x}\right)}^2{\left(y-\overline{y}\right)}^2}$$(5) x and y representing the two series with n number of elements, $\overline{x}$ and $\overline{y}$ representing mean values of the two series. The range –1 to 1 contains the values of r. Positive values close to 1 indicates a strong positive correlation between two series, on the other hand negative values indicate an inverse correlation.

2.2.2. Methodology for standardized precipitation index (SPI)

The SPI is a common indicator of drought that does not require information about land surface conditions and needs only precipitation data to compute drought properties. According to McKee et al. (1993) and Wu et al. (2007), the SPI can be calculated in a given year ‘o’, a calendar month ‘p’ and for time scale q, by the following steps;

• The long-term record of precipitation is fitted to a probability distribution, which is then transformed into a normal distribution.

• The 1st step is the computation of the cumulative precipitation data, $X_{op}^q(o=1,2,3,\dots ,n)$ for a period of interest p.

• The 2nd step is to fit a cumulative probability distribution (commonly gamma distribution function), but in the current study different distributions have been applied. The PDF of gamma is defined as; (6) $$f\left(x\right)=\frac{1}{{\beta }^{\alpha }\Gamma \left(\alpha \right)}{x}^{\alpha -1}{e}^{-\frac{x}{\beta }},\mathit{for}\mathrm{ }x\ge 0$$(6)

• A cumulative probability distribution for the particular time scale and given month of the observed precipitations event, (7) $$F\left(x\right)={\int }_0^{x}f\left(x\right)dx={\int }_0^x\frac{1}{{\beta }^{\alpha }\Gamma \left(\alpha \right)}{x}^{\alpha -1}{e}^{-\frac{x}{\beta }}dx$$(7)

• The cumulative probability of each observed precipitation event x_i can be derived by Eq. (7)(7) $$F\left(x\right)={\int }_0^{x}f\left(x\right)dx={\int }_0^x\frac{1}{{\beta }^{\alpha }\Gamma \left(\alpha \right)}{x}^{\alpha -1}{e}^{-\frac{x}{\beta }}dx$$(7) . An equiprobability transformation is then made from the cumulative probability to the standard normal random variable Z with zero mean and unit variance, where the SPI takes on the value of Z. (8) $$\mathit{SPI}={\phi }^{-1}\left[F\left(X\right)\right]$$(8)

• The distribution of precipitation may contain zeros. For instance, taking a special case of Gamma distribution, suppose the probability of all zero values in a time series of IWP_(r,p) is denoted by q. Let m be the number of all zero values and n is the total number of observations contained in IWP_(r,p) time series, therefore, q can be estimated by the ratio of m and n i.e. m/n, (9) $$H\left(x\right)=q+\left(1-q\right)F\left(x\right)$$(9)

• As precipitation is not distributed normally, so an equiprobability transformation is done from the CDF of the mixed distribution to CDF of a standard normal distribution, with zero mean and variance 1, which is given by; (10) $$\mathit{SPI}={\phi }^{-1}\left[H\left(X\right)\right]$$(10)

• Because the SPI is normalized, wetter and drier climates can be represented in the same way, and wet periods can also be monitored using the SPI.

• The 3rd step is to show the adequacy of the selected distribution, using some numerical or graphical methods.

• The 4th step is to verify the normality of SPI using numerical or graphical techniques.

2.3. Proposing a hydrological drought index: regionally improved weighted standardized drought index (RIWSDI)

In this study, we used IWP(r,p) weights in the placement of generally simple precipitation series in the SDI technique. As the estimates IWP(r,p) are representative more regionally which encounter direct effect of maximum and minimum temperatures in the estimation stage. Consequently, the inclination of using IWP(r,p) estimates in place of simple precipitation records is logically adequate. Therefore, to achieve the SDI, follow the instruction of Stagge et al. (2015) on the parametric standardization approach, the current study integrates the probability distributions, that fits well on the specified temporal series of IWP(r,p) estimates. Furthermore, this research incorporates different probability plotting (PP) methods for the assessment of the coherence and validity of IWP(r,p) estimates in non-parametric approach. The following are the brief descriptions of both the methods.

2.3.1. Parametric mechanism: selection of optimal probability distributions

In parametric approach, we arrange a list of different, 32 candidate distributions among various commonly available probability distributions. Although, the study suggests the probability distributions having multiple-parameters, for instance, in place of a Weibull and Gamma distributions having two-parameters, a goodness of fit must be applied to a four-parameter Weibull and three-parameter Gamma, respectively. Where the lowest value of Bayesian Information Criteria (BIC) Schwarz (1978) indicates the optimal probability distribution.

In the experimental and computational analysis, we incorporate Tellinghuisen and Spiess (2014) an r package to achieve the optimum probability distribution BIC for all SPI-3 and 12 – and RIWSDI-3 and 12. The study consist of 32 highly parametric distributions, for example; gumbel distribution, generalized extreme value, and generalized normal distribution, etc. Chi-square, Anderson Darling, and Shapiro Wilk tests are used for the nomination of the optimal distributions. Whereas, a Levenberg-Marquardt algorithm using minpack.lm Elzhov et al. (2010) an r package, is used to estimate the parameters of each chosen distribution. Furthermore, the Cumulative Distribution Function (CDF) of the selected optimal distribution is then converted using the method in Eq. (9)(9) $$H\left(x\right)=q+\left(1-q\right)F\left(x\right)$$(9) .

In researches related to hydrology, particularly in proposing a new drought index as well as their comparative assessment, r is the most frequently used statistical technique (Tsakiris and Vangelis, 2005; Ali et al., 2017). Although, the existing various drought indices bounds to choose optimum and most pertinent drought indices. In the previous researches, numerous authors proposed new drought indices related to hydrology including Cumbie-Ward and Boyles (2016); Jain et al. (2015); Naumann et al. (2014) and Ye et al. (2016) and compared those indices with the well-known index: Standardized Precipitation Index (SPI) (McKee et al., 1993). RIWSDI is mainly based on a long-term rainfall (precipitation) series as well as on average temperature as auxiliary variable, and is used to highlight the insufficiency in the amount of precipitation for various time scales (1, 3, 6, 9, 12, 24, 48) at a particular chosen station. Initially, the RIWSDI method is based on the standardization of CDF of the chosen distribution from 32 different candidate distributions. Where negative and positive values of RIWSDI respectively show less than or greater than the median precipitation.

Here, the study finds out that for defining drought related to hydrology, RIWSDI at 3 and 12-month time scale procedure uses, respectively, three- and twelve-months’ average data of monthly precipitation records as well as the auxiliary variable. Therefore, a comparative analysis of RIWSDI with SPI at 3 and 12-month time scale, is considered. Habibi et al. (2018) stated that most of the researches related to drought monitoring and hydrology, SPI-3 and SPI-12 are most extensively, efficient and effective used indicators for highlighting drought related to hydrology. Furthermore, characterization of drought based on 3 and 12-month timescale show an overview of the seasonal drought and an overall behaviour of hydrological conditions related to regions (Gumus and Algin, 2017).

In computational and experimental analysis, it has been shown that the utilization of Gamma distribution is appropriate for modelling precipitation records at different time scales across all accumulation periods and regions within Europe (Stagge et al., 2015). That study is for Europe, as a result, a PDF is required among the enlisted probability functions for the regions under study. In the current study, we are following guidelines of Stagge et al. (2015) in the process of estimation. Hence, a general expression of RIWSDI index can be written as follows, (11) $$IW{P}_{\left(r,p\right)}\sim \mathit{PDF}\left(\mathit{Parameter}\mathrm{ }1,\mathit{Parameter}\mathrm{ }2,\dots ,\mathit{Parameter}\mathrm{ }\mathrm{n}\right)$$(11) where IWP_(r,p) is the monthly cumulative total of either using regression or product type estimator for the improvement of rainfall data, and PDF shows the optimum probability function with n parameters. Estimation of the values of RIWSDI can be done by normalizing CDF of the pertinent selected PDF, which are fitted to the time series records of improved monthly cumulative precipitations.

In Eq. (9)(9) $$H\left(x\right)=q+\left(1-q\right)F\left(x\right)$$(9) , a small amendment is made in the CDF to fully adjust the adverse effect of non-precise values of the series. Furthermore, following McKee et al. (1993) and Ali et al. (2017), quantitative records of RIWSDI are then classified in accordance to the severity of drought. The drought characterization for SPI and RIWSDI can be seen in Table 2.

Table 2. Categories of SPI and RIWSDI.

S. No.	SPI and RIWSDI values	Categories
1.	2.00 and above	Extremely wet
2.	1.50 to 1.99	Very wet
3.	1.00 to 1.49	Moderate wet
4.	–0.99 to 0.99	Near normal
5.	–1.00 to –1.49	Moderate drought
6.	–1.50 to –1.99	Severe drought
7.	–2.00 and less	Extremely drought

2.3.2. Non-parametric mechanism: incorporation of graphical methods

In each of the probabilistic models, as a result ambiguity repeatedly exists in the precise and accurate estimation procedures (Parker, 2014). Furthermore, the selection of the precise probability distributions considering each indicator is exclusively subjective in its nature. To avoid these kind of problems, for non-parametric drought monitoring, Hao and AghaKouchak (2014) gave an idea by using Probability Position formulas (PP-formula) of Gringorten (1963) as an alternative technique of Gamma distribution for obtaining the SDI. The fact behind the use of graphical technique is to fascinate the extreme events, and as a result, to reduce the errors, in precise as well as accurate estimation of the specified drought index. Farahmand and AghaKouchak (2015); Ghamghami et al. (2017) and Zhang et al. (2018) also used a non-parametric approach for drought monitoring. Although, the behaviour of data is varying, as it is varied by place to place, therefore, it is not enough to the incorporate only one probability plotting position. Stagge et al. (2015) studied that it is not enough to use only the Gamma distribution to capture and observe the behaviours of different climatic regions to acquire the drought indices. Hence, the deployment and use of different PP-formulas are mandatory for evaluating different behaviours of varies specified probability distributions (Cunnane, 1978; Vogel, 1986; Shukri Yah et al., 2012). In the current study rather than utilizing Gingorten PP-formula, we use other six well known non-parametric PP-formulas for the computation of RIWSDI.

3. Results and discussion

3.1. Basic statistics, deviations, and temporal behaviour

Statistics of each region’s stations is presented in Table 1. It contains the information of the precipitation and temperature of all the stations included in different regions. Coordinates and Altitude(m) can also be seen from Table 1. Mean, Minimum and Maximum precipitation as well as temperature can be obtained from this table. For example, in C1 (Region-1), Astor station has a mean, minimum and maximum annual precipitation of 39.5(mm), 21.6(mm) and 72.8(mm), respectively. Similarly, it has a mean, minimum and maximum annual temperature of 9.9(°C), 8.3(°C) and 11.3(°C), respectively. Astor has a latitude (N) and longitude (E) of 35.3570° and 74.8624°, respectively. Astor is 2546(m) above the sea level. Similarly, we can see the statistics of each and every station from Table 1.

SPI-12 is the most commonly used drought index for annual monitoring and characterization of hydrological drought (Habibi et al., 2018). In previous research, several studies proposed hydrological drought indices and compared it with SPI-12 (Naumann et al., 2014; Jain et al., 2015; Cumbie-Ward and Boyles, 2016). In this study, SPI-3 and SPI-12 with the comparison of RIWSDI-3 and RIWSDI-12 have been computed, before the standardization of IWP_(r,p), a little graphical analysis is done by assessing the temporal behaviour and deviations in improving precipitation records with those which are used in the SPI-3 and SPI-12. SPI-3 and SPI-12 uses three- and twelve-month average of monthly precipitation records, respectively. Analogous to SPI, the proposed structure of precipitation records used in the IWP_(r,p) model has the same mathematical structure and rationale. Therefore, it is necessary to show how auxiliary information plays a role in the temporal estimation of precipitation records. Figure 3 shows the graphical representation of temporal precipitation records of regionally weighted and the usual records for all 7 regions (clusters) along with their stations as well. From Fig. 3 we observed that there are significant changes, particularly in upper precipitation record of Astor, Dir, Parachinar, Badin and Dalbandin stations. These contrasts reveal how drought characterization and analysis of regional meteorology of climates and their phenomena can be modified for the betterness by the introduction of advanced estimation techniques. It can be seen that the improved records in all the stations shows a better trend. Furthermore, before the standardization of the precipitation records, we may obtain the assorted proposition of the appropriate probability distributions for improved and simple precipitation records. Table 3 summarizes the statistics of the selected regions (clusters). It is found that the mean precipitation for simple and improved precipitation are almost the same. But, there is a significant difference, between the standard deviation (SD) and coefficient of variation (CV), of the saying records.

Fig. 3. Comparison of precipitation records.

Table 3. Statistics of the selected clusters (Regions).

Regions	Mean	P(MA) (mm) S.D	C.V	Mean	IP(MA) (mm) S.D	C.V	Mean	AT(MA)(°C) S.D	C.V
1	19.2	5.8	30.5	18.3	5.0	27.5	14.7	0.5	3.7
2	55.7	10.0	17.9	51.9	9.1	17.5	18.7	0.5	2.7
3	120.4	17.4	14.5	117.4	15.1	12.9	18.8	0.5	2.8
4	37.3	9.5	25.5	34.9	7.4	21.3	24.7	0.5	2.0
5	13.2	8.1	61.9	12.4	7.3	58.6	26.8	0.5	1.7
6	13.1	5.9	45.0	12.9	5.6	43.6	24.7	0.5	2.0
7	11.0	6.7	61.3	11.3	7.3	64.7	25.9	0.4	1.6

Precipitation Monthly Annual (P(MA)), Improved Precipitation Monthly Annual (IP(MA)), Average Temperature Monthly Annual (AT(MA)), Standard Deviation (SD), Coefficient of Variation (CV).

3.2. Parametric computation

Analogous to RIWSDI this study also calculates SPI. Subsequently, for each of the indicator, some sort of list of probability distributions is equipped in order to check their optimum fitness. Goodness-of-fit (GoF) test statistics are usually used for confirming validity and for choosing the best-fit among different distributions for a specific data set. In the current study BIC is used for choosing the best-fit. Numerical (BIC) as well as graphical performance (pdf) were used in the current study to assess and select the best fitted probability distribution among other 32 candidate distribution. A pdf with lowest value of BIC is considered as a best fit. The Bayesian Information Criterion (BIC) Konishi and Kitagawa (2008); Sakamoto et al. (1986) is applied here. A distribution, for which the value of BIC calculated by means of the following Eq. (12)(12) $$\mathit{BIC}=-\frac{2\hspace{0.17em}\text{ln}\hspace{0.17em}L}{N}+\frac{k\hspace{0.17em}\text{ln}\hspace{0.17em}N}{N}$$(12) is the lowest, is considered as the best fit, (12) $$\mathit{BIC}=-\frac{2\hspace{0.17em}\text{ln}\hspace{0.17em}L}{N}+\frac{k\hspace{0.17em}\text{ln}\hspace{0.17em}N}{N}$$(12) where L, k and N is the function of likelihood, number of the estimated parameters and total number of observations respectively, for the analyzed data (Kotowski and Kaźmierczak, 2013).

Tables 4–7 summarize the CDF of the optimum probability functions for all the indicators at all study regions and their associated stations. It also shows the chosen probability function and their estimated parameters as well as BIC values. Selection of each probability distribution is based on (weighted) residual sum-of-squares as the minimization criterion based on the Levenberg-Marquardt algorithm. Where, the estimation phase consists of method of moments, method of maximum likelihood estimation and method of L-moments. All these methods were implemented using lmom R packages. Especially in SPI and RIWSDI, fitness and selection of different probability distribution validate the finding of Stagge et al. (2015). A number of modifications to SPEI and SPI methodology, and also for assessing SPEI and SPI an updated procedure based on Shapiro-Wilk test were proposed by Stagge et al. (2015). They found gamma distribution and generalized extreme value distribution for SPI and SPEI respectively as the best fitted distribution in their study.

Table 4. Top-notch distributions with their parameters for Region-1.

RIWSDI/SPI		SPI3	SPI12
RIW1	Function parameters BIC	Generalised Normal α = 33.40, ξ = 45.92, κ = −0.57 –1099	Gumbel µ = 195.04, β = 53.50 –736.81
Astor	Function parameters BIC	Gamma k = 1.77, θ = 0.014 –1289.66	Triangular a = 160.30, b = 469.96, c = 766.18 –913.50
Bunji	Function parameters BIC	Generalised-normal α = 27.88, ξ = 26.21, κ = −0.70 –827.59	Laplace µ = 151.19, b = 77.90 –692.36
Chilas	Function parameters BIC	4P Beta α1 = 0.98, α2 = 6.93, a = 2.50, b = 410.24 –815	Laplace µ = 187.10, b = 102.14 –1335.20
Gilgit	Function parameters BIC	Gamma k = 1.11, θ = 0.03 –836.32	Inverse Gamma α = 7.28, β = 934.80 –1261.47
Gupis	Function parameters BIC	Gamma k = 0.78, θ = 0.02 –1270.42	Johnson SB ξ = 98.74, λ = 51.22, γ = −0.60, δ = 0.79 –940.08
Skardu	Function parameters	Generalised extreme value µ = 27.35, ξ = 28.08, σ = 0.38	Gamma k = 5.30, θ = 0.02
	BIC	–656.22	–1032.46

Table 5. Top-notch distributions with their parameters for Region-3.

RIWSDI/SPI		SPI3	SPI12
RIW3	Function parameters	Gumbel µ = 254.36, β = 160.12	Inverse Gaussian µ = 1412.11, λ = 0.00
	BIC	–1187.93	–1119.06
Balakot	Function parameters	Gumbel µ = 283.18, β = 186.72	Generalised extreme value µ = 1428.91, ξ = 249.33, σ = 0.19
	BIC	–1780.05	–1224.08
Garhidupatta	Function parameters BIC	Logistic µ = 331.46, σ = 115.38 –1728.47	Cosine µ = −33873.65, σ = 803.06 –1030.12
Islamabad	Function parameters	Gamma k = 1.52, θ = 0.004	Normal µ = 1426.21, σ = 368.42
	BIC	–1046.01	–1435.98
Kakul	Function parameters BIC	Triangular a = −31.70, b = 201.11, c = 796.58 –1176.07	Inverse Gaussian µ = 1344.88, λ = 0.00 –1476.64
Kotli	Function parameters BIC	Generalised extreme value µ = 200.89, ξ = 150.87, σ = 0.28 –1487.17	Cosine µ = 1237.87, σ = 724.98 –1721.01
Murree	Function parameters BIC	Gumbel µ = 312.41, β = 222.11 –1609.62	Inverse Gaussian µ = 1756.46, λ = 0.00 –1188.48
Muzaffarabad	Function parameters	Gumbel µ = 267.78, β = 178.14	Gumbel µ = 1372.88, β = 232.52
	BIC	–1547.78	–1072.52

Table 6. Top-notch distributions with their parameters for Region-6.

RIWSDI/SPI	SPI3		SPI12
RIW6	Function	Gamma	Trapezoidal
	Parameters	k = 1.14, θ = 0.03	a = −5.40, b = 153.4, c = 102.5, d = 325.37
	BIC	−906.97	−655.33
Dalbandin	Function	3P Weibull	Gumbel
	Parameters	µ = 0.82, α = 0.48, β = 42.48	µ = 51.24, β = 45.21
	BIC	–820.47	–1031.25
Kalat	Function	4P Beta	Trapezoidal
	Parameters	α1 = 0.38, α2 = 2.56, a = 2.42, b = 302.62	a = 6151, b = −6170.14, c = 145, d = 457.65
	BIC	–918.54	–1365.03
Khuzdar	Function	4P Beta	Gamma
	Parameters	α1 = 0.97, α2 = 8.15, a = 2.50, b = 714.26	k = 4.95, θ = 0.02
	BIC	–933.31	–818.74
Lasbela	Function	4P Beta	Gamma
	Parameters	α1 = 0.45, α2 = 77.80, a = 1.49, b = 9769.92	k = 2.01, θ = 0.01
	BIC	–1204.54	–848.40
Nokkundi	Function	Generalised normal	Gamma
	Parameters	α = 226.14, ξ = 71.56, κ = −3.20	k = 0.86, θ = 0.02
	BIC	–873.31	–882.77
Panjgur	Function	3P Weibull	Logistic
	Parameters	µ = 0.99, α = 0.84, β = 50.57	µ = 84.07, σ = 32.69
	BIC	–767.37	–588.99
Quetta	Function	3P Weibull	Laplace
	Parameters	µ = 2.49, α = 0.79, β = 130.59	µ = 215.85, b = 134.21
	BIC	–1364.30	–1075.69

Table 7. Top-notch distributions with their parameters for Region-7.

RIWSDI/SPI		SPI3	SPI12
RIW7	Function	4P Beta	Rayleigh
	Parameters	α1 = 0.46, α2 = 2.96, a = 1.18, b = 247.03	µ = −40.46, σ = 131.59
	BIC	–752.57	–775.79
Jiwani	Function	Generalised Normal	Gamma
	Parameters	α = 268.56.ξ = 94.80, κ = −2.91	k = 1.13, θ = 0.01
	BIC	–734.53	–739.74
Karachi	Function	3P Weibull	3P Weibull
	Parameters	µ = 2.01, α = 0.44, β = 94.20	µ = 25, α = 0.98, β = 221.30
	BIC	–1285.18	804.69
Pasni	Function	3P Weibull	Exponentail
	Parameters	µ = 0.99, α = 0.57, β = 147.63	λ = 0.008
	BIC	–1045.91	–713.748

Figure 4 shows temporal behaviours of RIWSDI-3 and SPI-3 indices, Fig. 5 shows temporal behaviours of RIWSDI-12 and SPI-12 indices. Figures 6 and 7 shows the histograms, Q-Q plots and also provides the values of correlation coefficients between RIWSDI and SPI, respectively. In almost all regions the correlation coefficient is high, which shows the significance of the correlation between RIWSDI and SPI. Whereas some stations show low correlation i.e. Region-2-Chitral (0.32), Drosh (0.46) for SPI3, Region-2-Chitral (0.46), Drosh (0.48), Zhob (0.34) for SPI-12, and Region-6-Nokkundi (0.45) for SPI-12 as well. This shows that RIWSDI can be recommended as an alternate drought index which incorporates auxiliary information-based precipitation data for characterization of hydrological drought. Figures 8–11 shows the plots of the chosen distributions among 32 distributions on the basis of the lower value of BIC.

Fig. 4. Comparison of RIWSDI-3 and SPI-3.

Fig. 5. Comparison of RIWSDI-12 and SPI-12.

Fig. 6. Correlations of RIWSDI-3 and SPI-3.

Fig. 7. Correlations of RIWSDI-12 and SPI-12.

Fig. 8. Observed probability distribution functions for RIWSDI-3 and SPI-3 (Regions 1-4).

Fig. 9. Observed probability distribution functions for RIWSDI-3 and SPI-3 (Regions 5-7).

Fig. 10. Observed probability distribution functions for RIWSDI-12 and SPI-12 (Regions 1-4).

Fig. 11. Observed probability distribution functions for RIWSDI-12 and SPI-12 (Regions 5-7).

3.3. Non-parametric computation

Beside of cautious determination of probability distribution, our observed findings exhibit that outliers as well as extreme values cannot be fully capture by probability function. Figures 8–11 reveals that for Chilas station 4 P Beta distribution is the best fit with low BIC value, does not capture uncertainty in the significant parts of the data for SPI-3 model in Region-1. Similar results can be seen from Figs. 8–11 in some stations. All chosen distributions appeared to inefficient for coverage of extreme records. Hence, in order to check and validate furthermore, the current study incorporates six Probability Plotting Position (PPP) formula, see Table 8.

Table 8. Probability plotting position techniques.

No	Formula (Authors)	Year	Description
PPP-1	Hazen	1914	P(xi) = i−0.5/n
PPP-2	Weibull	1939	P(xi) = i/n + 1
PPP-3	Beard	1943	P(xi) =i−0.31/n + 0.38
PPP-4	Blom	1954	P(xi) =i−0.375/n + 0.25
PPP-5	Chegodayov	1955	P(xi) = i−0.3/n + 0.4
PPP-6	Tukey	1962	P(xi) =3i−1/3n + 1

After standardization the vector of time series based on the PPP, we see that RIWSDI remains aligned with SPI. Table 9 summarizes the correlation between RIWSDI and SPI of the non-parametric approach. There is a highly positive correlation between the values of SPI and RIWSDI for all the stations and their respective regions, under different methods of non-parametric approach, 3 out of 50 stations have weak but positive correlation with their respective regions.

Table 9. Non-parametric determination: correlation of RIWSDI-12 and SPI-12.

Regions	Stations	Hazen	Weibull	Beard	Blom	Chegodayov	Tukey
	Astor	0.81	0.81	0.81	0.81	0.81	0.81
	Bunji	0.74	0.74	0.74	0.74	0.74	0.74
C1	Chilas	0.74	0.74	0.74	0.74	0.74	0.74
	Gilgit	0.75	0.75	0.75	0.75	0.75	0.75
	Gupis	0.64	0.64	0.64	0.64	0.64	0.64
	Skardu	0.68	0.68	0.68	0.68	0.68	0.68
	Cherat	0.75	0.75	0.75	0.75	0.75	0.75
	Chitral	0.45	0.45	0.45	0.45	0.45	0.45
	Dir	0.54	0.54	0.54	0.54	0.54	0.54
	Drosh	0.47	0.47	0.47	0.47	0.47	0.47
C2	Kohat	0.56	0.56	0.56	0.56	0.56	0.56
	Parachinar	0.67	0.67	0.67	0.67	0.67	0.67
	Peshawar	0.77	0.77	0.77	0.77	0.77	0.77
	Risalpur	0.72	0.72	0.72	0.72	0.72	0.72
	Zhob	0.34	0.34	0.34	0.34	0.34	0.34
	Balakot	0.71	0.71	0.71	0.71	0.71	0.71
	Garhidupatta	0.60	0.60	0.60	0.60	0.60	0.60
	Islamabad	0.65	0.65	0.65	0.65	0.65	0.65
C3	Kakul	0.82	0.82	0.82	0.82	0.82	0.82
	Kotli	0.72	0.72	0.72	0.72	0.72	0.72
	Murree	0.77	0.77	0.77	0.77	0.77	0.77
	Muzaffarabad	0.78	0.78	0.78	0.78	0.78	0.78
	Bahawalnagar	0.66	0.66	0.66	0.66	0.66	0.66
	Bahawalpur	0.66	0.66	0.66	0.66	0.66	0.66
	DI Khan	0.57	0.57	0.57	0.57	0.57	0.57
	Faisalabad	0.56	0.56	0.56	0.56	0.56	0.56
C4	Lahore	0.73	0.73	0.73	0.73	0.73	0.73
	Mianwali	0.67	0.67	0.67	0.67	0.67	0.67
	Multan	0.59	0.59	0.59	0.59	0.59	0.59
	Sargodha	0.63	0.63	0.63	0.63	0.63	0.63
	Sialkot	0.66	0.66	0.66	0.66	0.66	0.66
	Badin	0.85	0.85	0.85	0.85	0.85	0.85
	Chhor	0.82	0.82	0.82	0.82	0.82	0.82
	Hyderabad	0.78	0.78	0.78	0.78	0.78	0.78
	Jacobabad	0.77	0.77	0.77	0.77	0.77	0.77
C5	Khanpur	0.65	0.65	0.65	0.65	0.65	0.65
	Nawabshah	0.90	0.90	0.90	0.90	0.90	0.90
	Padidan	0.83	0.83	0.83	0.83	0.83	0.83
	Rohri	0.71	0.71	0.71	0.71	0.71	0.71
	Sibbi	0.68	0.68	0.68	0.68	0.68	0.68
	Dalbandin	0.72	0.72	0.72	0.72	0.722	0.72
	Kalat	0.69	0.69	0.69	0.69	0.69	0.69
	Khuzdar	0.82	0.82	0.82	0.82	0.82	0.82
C6	Lasbela	0.62	0.62	0.62	0.62	0.62	0.62
	Nokkundi	0.44	0.44	0.44	0.44	0.44	0.44
	Panjgur	0.73	0.73	0.73	0.73	0.73	0.73
	Jiwani	0.78	0.78	0.78	0.78	0.78	0.78
	Quetta	0.52	0.52	0.52	0.52	0.52	0.52
C7	Karachi	0.80	0.80	0.80	0.80	0.80	0.80
	Pasni	0.71	0.71	0.71	0.71	0.71	0.71

4. Conclusion

The new proposed index shows significant advantages over the former index (SPI) by incorporating apart from the record of precipitation, an additional meteorological parameter as an auxiliary variable, the average temperature. It is figured that although the RIWSDI usually responds in an identical fashion towards the SPI, it is more delicate and appropriate in regional drought analysis. Drought is generally acknowledged as a regional phenomenon. Even so, facts are accumulated by determined meteorological stations, which can also be regarded as representing the regions related to them. This study advantages drought-observing component by the incorporation of improved precipitation series in standard strategies of SDI technique. Whereas, improvement of the precipitation series is derived from the use of auxiliary information in the estimation phase of mean rainfall amount. In this article, the simple regression and product estimator approaches are utilized to weight the rainfall amount of each region. Therefore, the current study proposed a new regional hydrological drought index: The Regionally Improved Weighted Standardized Drought Index (RIWSDI). Performance of RIWSDI is evaluated by verifying the direction of the relationship, the form (shape) of the relationship, and the degree (strength) of the relationship among the different values of SPI through Pearson correlation statistics. To check out the consistency and efficiency of the proposed index RIWSDI, this study contains seven meteorological regions dispersed in several climatic settings of Pakistan, see Fig. 1. To compare and compute the values of RIWSDI and SPI for time scale 3 and 12, estimation methodology comprises of both the parametric Stagge et al. (2015) and non-parametric Hao and AghaKouchak (2014) approaches.

Comparative analysis of improved precipitation records and simple precipitation records can be observed from Fig. 3. A significant difference between the probability distributions of simple and improved precipitation series can be observed from Tables 4–7, and Figs. 8–11. As we described previously, parametric distributions which are used to compute SPI and RIWSDI should be capable of providing normally distributed SPI and RIWSDI series. An SPI and RIWSDI series are consider non-normal if the criteria given below is satisfy simultaneously: (a) p-values ≤ 0.10 and Shapiro–Wilk (W) statistic lower than 0.96; (b) an absolute value of median > 0.05. Further information on W test may be found in numerous studies including Razali and Wah (2011). W test result can be found in Tables 10 and 11. By analyzing Tables 10 and 11, it can be concluded and verified that all RIWSDI and SPI for time scale-3 and 12 are normally distributed. For example, taking the values of RIWSDI-3 and RIWSDI-12 of Region-1, all the criterion of distribution to be normal are satisfied. The confirmation of normality of RIWSDI and SPI for all other regions and the stations included in the regions can be confirmed by analyzing the results of W, p-value and absolute median of the respective regions and stations. Moreover, on the contrary with the standardization based on various distribution, the pattern of the records of RIWSDI are very closed to that of SPI (time scale 3 and 12), see Figs. 4 and 5. It can be observed that the simple precipitation records within the regions are relatively different, whereas it can be seen from Figs. 6–7 that the correlations among the stations are significantly low, whereas the correlation between RIWSDI and SPI (time-scale 3 and 12) for different regions and the stations within the regions are significantly high. Generally, comparative evaluation indicates that RIWSDI is significantly high correlated with SPI (time-scale 3 and 12) in both standardization i.e. parametric and nonparametric. Overall, comparative evaluation indicates that RIWSDI is strongly correlated with SPI in each parametric and non-parametric standardization. On the other hand, some of the discrepancies can be observed in parametric standardization and show low but positive correlation.

Table 10. Shapiro–Wilk test for Region (1–5).

		Time scale-3			Time scale-12
Regions	Stations	W	p-value	Med	W	p-value	Med
	Region-1	0.99	0.52	0.03	0.99	0.07	0.04
	Astor	0.99	0.04	0.04	0.99	0.00	0.06
	Bunji	0.99	0.051	0.02	0.99	0.00	0.03
C1	Chilas	0.99	0.00	0.001	0.97	0.00	0.12
	Gilgit	0.99	0.29	0.02	0.99	0.002	0.05
	Gupis	0.99	0.00	0.007	0.97	0.00	0.06
	Skardu	0.99	0.002	0.04	0.99	0.03	0.03
	Region-2	0.99	0.003	0.04	0.99	0.09	0.03
	Cherat	0.99	0.000	0.02	0.99	0.001	0.06
	Chitral	0.99	0.06	0.05	0.99	0.05	0.04
	Dir	0.99	0.08	0.02	0.98	0.00	0.04
C2	Drosh	0.99	0.02	0.03	0.98	0.00	0.08
	Kohat	0.99	0.005	0.01	0.97	0.00	0.07
	Parachinar	0.97	0.00	0.08	0.87	0.00	0.02
	Peshawar	0.99	0.13	0.04	0.99	0.004	0.06
	Risalpur	0.99	0.73	0.016	0.99	0.000	0.04
	Zhob	0.99	0.26	0.008	0.99	0.000	0.12
	Region-3	0.99	0.78	0.03	0.99	0.001	0.001
	Balakot	0.99	0.61	0.001	0.97	0.000	0.15
	Garhidupatta	0.99	0.79	0.03	0.99	0.02	0.03
C3	Islamabad	0.99	0.11	0.005	0.99	0.01	0.08
	Kakul	0.99	0.09	0.05	0.99	0.04	0.001
	Kotli	0.99	0.49	0.05	0.99	0.001	0.03
	Murree	0.99	0.06	0.004	0.99	0.005	0.02
	Muzaffarabad	0.98	0.00	0.003	0.98	0.00	0.09
	Region-4	0.99	0.08	0.06	0.98	0.00	0.03
	Bahawalnagar	0.99	0.004	0.001	0.99	0.001	0.03
	Bahawalpur	0.99	0.01	0.04	0.98	0.00	0.04
	DI Khan	0.99	0.57	0.006	0.98	0.00	0.04
C4	Faisalabad	0.97	0.00	0.03	0.99	0.008	0.03
	Lahore	0.99	0.06	0.04	0.98	0.00	0.08
	Mianwali	0.99	0.88	0.03	0.98	0.00	0.04
	Multan	0.99	0.00	0.03	0.99	0.002	0.03
	Sargodha	0.99	0.22	0.34	0.98	0.00	0.01
	Sialkot	0.99	0.54	0.01	0.99	0.00	0.07
	Region-5	0.99	0.33	0.066	0.99	0.00	0.11
	Badin	0.98	0.00	0.07	0.97	0.00	0.019
	Chhor	0.98	0.00	0.15	0.95	0.00	0.06
	Hyderabad	0.98	0.00	0.03	0.99	0.00	0.06
C5	Jacobabad	0.99	0.001	0.02	0.98	0.00	0.08
	Khanpur	0.99	0.004	0.02	0.98	0.00	0.02
	Nawabshah	0.99	0.001	0.03	0.99	0.001	0.07
	Padidan	0.99	0.01	0.08	0.99	0.01	0.01
	Rohri	0.98	0.00	0.20	0.99	0.00	0.98
	Sibbi	0.98	0.00	0.03	0.99	0.00	0.02

Table 11. Shapiro–Wilk test for Region (6–7).

		Time scale-3			Time scale-12
Regions	Stations	W	p-value	Med	W	p-value	Med
	Region-6	0.99	0.83	0.02	0.99	0.006	0.02
	Dalbandin	0.99	0.00	0.23	0.97	0.00	0.07
	Kalat	0.98	0.00	0.068	0.97	0.00	0.06
	Khuzdar	0.99	0.12	0.04	0.98	0.00	0.04
C6	Lasbela	0.98	0.00	0.01	0.99	0.03	0.01
	Nokkundi	0.97	0.00	0.64	0.98	0.00	0.008
	Panjgur	0.98	0.00	0.15	0.98	0.00	0.04
	Quetta	0.99	0.02	0.02	0.97	0.00	0.04
	Region-7	0.99	0.33	0.06	0.98	0.00	0.12
	Jiwani	0.98	0.00	0.35	0.98	0.00	0.01
C7	Karachi	0.99	0.01	0.007	0.98	0.00	0.06
	Pasni	0.99	0.00	0.30	0.97	0.00	0.07

Some extreme values are addressed by incorporating non-parametric methods for the analysis and comparison of RIWSDI and SPI. Here, six different probability plotting position formulas are used to handle the extreme values and outliers, see Table 8. Table 9 shows different results of the correlation of RIWSDI and SPI under non-parametric approach. It can be observed that SPI (3–12) values of all the stations are significantly high correlated with RIWSDI (3–12) values of the respective regions.

The foremost advantage of RIWSDI is to signify hydrological drought primarily based on regionally improved series of precipitation. Though, the limitation of the current study is that RIWSDI cannot be generalized in the settings of multiscales (Edwards, 1997).

Acknowledgements

The authors are quite grateful to the editor and reviewers for their helpful and insightful comments.

Disclosure statement

No potential conflict of interest was reported by the authors.

Additional information

Funding

This research is supported by the National Natural Science Foundation of China (Grant No. 71901109, No. 11561026, No. 71861012), Natural Science Foundation of Jiangxi, China (No. 20181BAB211020). Jiangxi Double Thousand Plan, Postdoctoral Foundation of Jiangxi Province (No. 2018KY08), Scientific Research Fund of Jiangxi Provincial Education Department (Grant No. G180267, No. GJJ190264) and Human and Social Science Foundation of Jiangxi Province (No. T J19202).

Appendix

Appendix A: Standardization

This step considers a suitable transformation technique for the standardization of the chosen CDF, as well as all numeric vectors contains the time series records primarily based on the PP-formula and CDF. Hence, following the studies of McKee et al. (1993) and Ali et al. (2017), the current research employed an approximate transformation given in Abramowitz and Stegun (1965) for the transformation of the CDF into a standard normal distribution, which is as follows: (A1) $$\mathit{RIWSDI}=-\left(z+\frac{c_0+c_{1}z+c_2{z}^2}{1+d_{1}z+d_2{z}^2+d_3{z}^3}\right)$$(A1) for z, (A2) $$z=\sqrt{\hspace{0.17em}\text{ln}\hspace{0.17em}\left(\frac{1}{{\left\{T\left(x\right)\right\}}^2}\right)},0\mathrm{<}\mathrm{T}\left(x\right)\le 0.5$$(A2) (A3) $$\mathit{RIWSDI}=+\left(z-\frac{c_0+c_{1}z+c_2{z}^2}{1+d_{1}z+d_2{z}^2+d_3{z}^3}\right)$$(A3) and also for, (A4) $$z=\sqrt{\hspace{0.17em}\text{ln}\hspace{0.17em}\left(\frac{1}{{\left\{1-T\left(X\right)\right\}}^2}\right)},0.5\mathrm{<}\mathrm{T}\left(x\right)\le 1$$(A4) where c₀ = 2.515517, c₁ = 0.802853, c₂ = 0.010328, d₁ = 1.432788, d₂ = 0.985269, and d₃ = 0.001308 are constants.

In the above scheme, $$T(x)=\{{}_{Q\left(x\right)\mathrm{ }\text{if}\mathrm{ }\text{standardization}\mathrm{ }\text{is}\mathrm{ }\text{done}\mathrm{ }\text{under}\mathrm{ }\text{non}-\text{parametric}\mathrm{ }\text{approach}}^{P\left(x\right)\mathrm{ }\text{if}\mathrm{ }\text{standardization}\mathrm{ }\text{is}\mathrm{ }\text{done}\mathrm{ }\text{under}\mathrm{ }\text{parametric}\mathrm{ }\text{approach}}$$ here, P(x) and Q(x) are probability vectors, which are defined respectively in Eq. (9) and in Table 8. After standardization the records of RIWSDI will contain zero (0) mean as well as unit (1) variance.