Hydrological Response to Meteorological Droughts in the Guadalquivir River Basin, Southern Iberian Peninsula

: Drought is an extreme phenomenon that will likely increase in frequency and severity in the current context of climate change. As such, it must be studied to improve the decision-making process in affected areas. As a semi-arid zone, the Guadalquivir River basin, located in the southern Iberian Peninsula, is an interesting area to perform this study. The relationship between meteorological and hydrological droughts is studied using drought indices with data from 1980 to 2012. The chosen indices are the Standardized Streamﬂow Index (SSI) and the Standardized Precipitation Evapotranspiration Index (SPEI). Their correlations are calculated, based on SPEI accumulation periods, and these values are analyzed with a principal component analysis to ﬁnd spatial patterns in drought behavior inside the basin. This analysis was performed for the continuous series and also for monthly series, to account for seasonal changes. It has been found that the relationship of drought types occurs at different time scales depending mainly on orography and catchment area. Two main patterns were found. Generally, for low altitudes and small catchment areas, accumulation periods are shorter indicating that hydrological system in this area respond rapidly to meteorological conditions. In mountainous parts of the basin, longer accumulation periods have a stronger inﬂuence due to effects such as snowmelt.


Introduction
Several regions around the world, their populations, and their economies, are currently affected by the extreme phenomenon of drought [1]. Climate change projections using regional models show that some arid and semiarid areas, such as the Iberian Peninsula, will likely suffer an increase in the number and severity of droughts [2,3]. Therefore, such extreme events must be studied in order to understand their mechanisms and how different parameters affect their development.
One of the mechanisms that can be studied is the relationship between different drought types. The definition of drought type depends on the variable that determines the occurrence of drought. For example, meteorological drought happens when rainfall is lower than a certain threshold; in contrast, hydrological drought appears when streamflow and runoff are lower than normal. Depending on the purpose of the definition, the timeframe may span from weeks or months used in policies to years used in megadrought studies. The relationship between hydrological and meteorological droughts has been studied in different areas of the world [4][5][6][7][8], and it depends on factors such as climate characteristics and orography [9,10]. Therefore, it is not an immediate conclusion that hydrological drought will happen shortly after meteorological drought. On the contrary, it requires a detailed quantitative analysis on each basin.
Drought indices are one of the tools that can be used to perform the aforementioned quantitative analyses. There are a plethora of indices, and each of them has a different Regarding the orography, the Guadalquivir basin consists of mountainous are the north and southeast, with a high presence of limestone and conglomerate. The no ern mountain range, Sierra Morena, is composed of low permeability soils, while so the southeastern mountains of Sierra Nevada and Sierra de Segura show moderate meability [18]. A long valley is located in the center of the basin, crossing from ea west, where the Guadalquivir flows. The river ends in a marsh area in the southwe the basin. As can be inferred from this description, the area cannot be defined as hom neous, increasing its interest as the place to study the relationship between drought ty Concerning the relationship between precipitation and streamflow, during spring mo (from March to May), streamflow will increase due to snowmelt, that is, due to preci Regarding the orography, the Guadalquivir basin consists of mountainous areas to the north and southeast, with a high presence of limestone and conglomerate. The northern mountain range, Sierra Morena, is composed of low permeability soils, while soils in the southeastern mountains of Sierra Nevada and Sierra de Segura show moderate permeability [18]. A long valley is located in the center of the basin, crossing from east to west, where the Guadalquivir flows. The river ends in a marsh area in the southwest of the basin. As can be inferred from this description, the area cannot be defined as homogeneous, increasing its interest as the place to study the relationship between drought types. Concerning the relationship between precipitation and streamflow, during spring months (from March to May), streamflow will increase due to snowmelt, that is, due to precipitation from past months. Therefore, different precipitation regimes can be found in this basin: rainfall regime, snowmelt and rainfall regime, and human-induced regime due to regulation [27].
This study used high quality observational data for the variables required in the calculation of drought indices: precipitation, potential evapotranspiration, and streamflow. Precipitation was obtained from the Spanish Precipitation at Daily Scale (SPREAD) dataset [28]. This dataset comprises gridded data based on previous observations. One of the reasons to choose this dataset is that the dimension of grid squares is 5 km, which is enough to describe precipitation in the entirety of the Guadalquivir basin at a fine and appropriate resolution. Potential evapotranspiration was extracted from the Spanish Reference Evapotranspiration (SPETo) gridded database [29]. This dataset contains weekly data at a resolution of 1.1 km, and it has been calculated using temperature, wind speed, dewpoint temperature, and sunshine duration data by using the Penman-Monteith formulation [30]. Regarding streamflow data, the main source of observational data is the Guadalquivir hydrographic authority (Confederación Hidrográfica del Guadalquivir) website [31], where data are freely available. Particularly, monthly reservoir data have a high quality with respect to series duration and low presence of gaps and inconsistencies. Nevertheless, some stations still showed gaps or had a short data timespan. These stations were not included in the analysis. Furthermore, time series had to be trimmed so that their data were simultaneous. SPREAD and SPETo end in December 2012, while most streamflow data begin or improve their quality in January 1980. Therefore, the shared dates of the three databases range from January 1980 to December 2012, and this is the period of study.
When managing streamflow data from reservoirs, there are two main values to consider: inflow and outflow. Inflow depends on the physical characteristics of the catchment basin, while outflow is totally dependent on dam operation. Certainly, the outflow is higher when extreme precipitation events occur, since dam sluices are opened to avoid floods. It is, however, more interesting for this study to consider only inflow data. The reason is the important human interaction component present in outflow data, which distorts the correlation between precipitation and streamflow. Therefore, for the purpose of this work, only inflow has been used to calculate the aforementioned relationships. Note, however, that inflow data are also affected by regulation from upstream reservoirs. Nevertheless, the impact of such reservoirs is lowered the further they are from the data source. Considering this explanation, the reservoirs that provided no useful information because of this effect, i.e., reservoirs that had another reservoir few kilometers upstream, were excluded from the analysis.
After discarding 19 stations with insufficient data quality or quantity, data from 43 reservoir stations for the period 1980-2012 were used to perform this study. The designation of these stations, as shown in Figure 1 and Table 1, is taken from the aforementioned database of observational reservoir data. For these stations, drought indices were calculated. Geographically, the most restrictive variable was the streamflow, since it stands as a series of non-grid points. Consequently, the 43 closest grid points from SPREAD and SPETo datasets were chosen for the calculation of SPEI.
The drought indices calculation is slightly different for SPEI and SSI. When computing SPEI, the hydrological balance (i.e., the difference between precipitation and evapotranspiration) is calculated as a monthly series, adding values from previous months depending on the accumulation period, which in this study ranges from 1 to 24. These data are then adjusted to a log-logistic distribution [12] in order to compute the probability series. These values must then be normalized [32], giving as a result the values of the drought index. Regarding SSI, streamflow data cannot be accumulated when calculating index values [13], so the accumulation period concept is not relevant for SSI. Additionally, the distribution to which the data must be adjusted has to be chosen between the following six [13]: Generalized Extreme Values, Pearson type III, Log-logistic, Log-normal, Pareto, and Weibull.
The best-fitting distribution needs to be chosen for each station, and after obtaining its probability values, the final index values are calculated following the same methodology of SPEI. By correlating the 1-24 months accumulated SPEI indices with SSI, which accounts only for the current month, the propagation from meteorological to hydrological drought can be analyzed for different areas of the river basin [5,33]. This response has been analyzed in two separate ways: on the one hand, the continuous SSI and SPEI series from 1980 to 2012 have been compared; on the other hand, both series were divided into months, showing how the relationship varies in different moments of the hydrological year. Additionally, the second monthly analysis could be helpful to identify the moment when the hydrological year starts in the Guadalquivir basin, since its definition is related to the correlation between precipitation and streamflow [34].
The proposed methodology, that correlates monthly SSI with multi-scale SPEI, intends to easily evaluate the response time of hydrological drought [35] and has been successfully applied in other basins in the Iberian Peninsula [5].
After computing correlations, a principal component analysis (PCA) was performed in order to analyze the different spatial patterns of correlation between the SSI and the SPEI at different time scales. PCA is widely used in climate research to find patterns in the data from a complex dataset [36,37]. The purpose of this PCA was to check if the basin could be divided in sub-basins with similar hydrological responses. It was performed using both the continuous series and the monthly series, and the divisions were made using the PCA scores and loadings.
PCA was performed in S-mode, so that the PCA explores relationships between the time series of the grid points, as opposed to T-mode, which explores relationships for a given event. Since the aim is to study spatial patterns of stations with a similar behavior of correlations between SPEI and SSI, S-mode is used [38]. PCA is able to reduce a large number of interrelated variables to a few independent principal components (PCs) that capture much of the variance of the original dataset. Moreover, it produces a few major spatial variability patterns (or empirical orthogonal functions, EOFs). Additionally, a varimax rotation was applied to the components, thus maximizing the variance of the squared correlation coefficients, and therefore better capturing the physically meaningful and simplified spatial patterns. This procedure allows the extraction of representative stations for each of the areas obtained in the regionalization process. Consequently, by performing this transformation, the found patterns were physically more stable [39], and it was easier to understand the information given by loadings [40]. Thus, scores of a given principal component show what accumulation periods dominate the correlation of stations with high loadings for the same principal component. By describing the basin with more than one principal component, different accumulation periods can be considered in the analysis.
Finally, the last part of the analysis includes a multilinear regression that tries to find if the patterns found during the PCA can be explained by the geographical characteristics of the basin. To this end, the loading obtained for each PC and for each station is introduced together with its values for the geographical variables in a linear regression model. Then, p-values of each variable are compared to find which ones are more relevant for each PC. The geographical variables were station altitude (Figure 1), catchment area (Figure 2a), and percentage of permeable soil in that catchment area ( Figure 2b). This way, it is possible to elucidate if the similarities in drought behavior (obtained by means of the PCA) can be explained by the physical characteristics of the stations. The variables corresponding to each station are detailed in Table 1.
Water 2022, 14, x FOR PEER REVIEW 7 of 20 PCA) can be explained by the physical characteristics of the stations. The variables corresponding to each station are detailed in Table 1. For clarity purposes, the methodology has been represented as a flowchart in Figure  3. For clarity purposes, the methodology has been represented as a flowchart in Figure 3.

Correlations between SSI and SPEI
When correlating the indices, it was considered that soil memory could introduce a delay in the correlation values between SPEI and SSI [41]. This situation would be particularly interesting for short accumulation periods. Figure 4 shows the correlation between SPEI and SSI for station 5018 (located to the northeastern area) accounting with different SPEI delays and accumulation periods. As can be seen, correlation decreases as the delay increases in this station for all accumulation periods. The same effect appeared in all stations of the basin. This is in accordance with the results of Koster and Suarez [10], where it is shown that the value of this lag is between 1 and 2 months in arid and semi-arid regions. The fact that the data used in this study are monthly accumulated instead of daily explains that correlations are higher for simultaneous SPEI and SSI series.

Correlations between SSI and SPEI
When correlating the indices, it was considered that soil memory could introduce a delay in the correlation values between SPEI and SSI [41]. This situation would be particularly interesting for short accumulation periods. Figure 4 shows the correlation between SPEI and SSI for station 5018 (located to the northeastern area) accounting with different SPEI delays and accumulation periods. As can be seen, correlation decreases as the delay increases in this station for all accumulation periods. The same effect appeared in all stations of the basin. This is in accordance with the results of Koster and Suarez [10], where it is shown that the value of this lag is between 1 and 2 months in arid and semi-arid regions. The fact that the data used in this study are monthly accumulated instead of daily explains that correlations are higher for simultaneous SPEI and SSI series. Once the effect of soil memory was checked, propagation from meteorological to hydrological drought was analyzed. It was discovered that the studied stations have  Once the effect of soil memory was checked, propagation from meteorological to hydrological drought was analyzed. It was discovered that the studied stations have different behaviors in their correlations. Figure 5 shows the SPEI accumulation period for which correlations between SSI and SPEI reach the maximum value. As can be seen, maximum correlations occur at short accumulation periods (around 2 or 3 months) for most stations. There are, however, several exceptions. Most notably, stations 5019 and 5025, located in the south, and 5032, in the center, show that maximum correlations appear with large accumulation periods of 24 months. This is the first sign of the special treatment that these three stations require, which will be discussed later in this study. This differentiation is required because of their distinct physical characteristics, mainly their catchment area, larger than 4000 km 2 , as shown in An interesting part of the basin appears in Sierra Nevada, to the southeast. There are two stations with maximum correlations at 21 and 24 months of accumulation, while the surrounding stations are between 10 and 13 months. Note that the highest altitudes of the whole basin are reached in this area, so water received in these stations comes mainly from meltdown. Thus, in this area, streamflows are not created rapidly after precipitation, they are formed several months thereafter. Figure 2a.
Three detailed examples of the evolution of correlation with the accumulation period are shown in Figure 6. Station 5016 is the westernmost station in the study area (Figure 1), and it is located in the header of its sub-basin, so it has no upstream regulation. Its Water 2022, 14, 2849 9 of 19 correlation peak appears in the 2-month accumulation period, and it decreases thereafter. The explanation for this behavior is that the catchment area is small and relatively irregular, causing water from precipitation events to rapidly find the stream. Station 5026 is located in the center and southern part of the Guadalquivir basin, and it has one of the largest catchment areas (Figures 1 and 2a). It shows a peak for the 13-month accumulation period. Note that, for this station, a significant part of its inflow comes from Sierra Nevada, which means that its streamflow is delayed and influenced by precipitations-snow-of several previous months. Lastly, station 5035, located in the northwestern part of the basin, is shown here since it presents a particular feature, with a local maximum in the 2-month accumulation period. However, its correlation later grows between the 3 and the 12-month accumulation periods.
Water 2022, 14, x FOR PEER REVIEW 9 of 20 different behaviors in their correlations. Figure 5 shows the SPEI accumulation period for which correlations between SSI and SPEI reach the maximum value. As can be seen, maximum correlations occur at short accumulation periods (around 2 or 3 months) for most stations. There are, however, several exceptions. Most notably, stations 5019 and 5025, located in the south, and 5032, in the center, show that maximum correlations appear with large accumulation periods of 24 months. This is the first sign of the special treatment that these three stations require, which will be discussed later in this study. This differentiation is required because of their distinct physical characteristics, mainly their catchment area, larger than 4000 km 2 , as shown in Figure 2a. An interesting part of the basin appears in Sierra Nevada, to the southeast. There are two stations with maximum correlations at 21 and 24 months of accumulation, while the surrounding stations are between 10 and 13 months. Note that the highest altitudes of the whole basin are reached in this area, so water received in these stations comes mainly from meltdown. Thus, in this area, streamflows are not created rapidly after precipitation, they are formed several months thereafter.
Three detailed examples of the evolution of correlation with the accumulation period are shown in Figure 6. Station 5016 is the westernmost station in the study area (Figure 1), and it is located in the header of its sub-basin, so it has no upstream regulation. Its correlation peak appears in the 2-month accumulation period, and it decreases thereafter. The explanation for this behavior is that the catchment area is small and relatively irregular, causing water from precipitation events to rapidly find the stream. Station 5026 is located in the center and southern part of the Guadalquivir basin, and it has one of the largest catchment areas (Figures 1 and 2a). It shows a peak for the 13-month accumulation period. Note that, for this station, a significant part of its inflow comes from Sierra Nevada, which means that its streamflow is delayed and influenced by precipitations-snow-of several previous months. Lastly, station 5035, located in the northwestern part of the basin, is shown here since it presents a particular feature, with a local maximum in the 2-month accumulation period. However, its correlation later grows between the 3 and the 12month accumulation periods.

PCA Analysis
When performing the PCA, the first step is to determine the number of principal components retained. In the case of continuous series, two components explain 96% of the total variance. After applying the varimax rotation, the first principal component explains

PCA Analysis
When performing the PCA, the first step is to determine the number of principal components retained. In the case of continuous series, two components explain 96% of the total variance. After applying the varimax rotation, the first principal component explains 67.5% of the total variance, while the second component explains an additional 28.5%.
As shown in Figure 7, PC1 mainly represents stations where short-time responses are the most important, while stations where PC2 has the highest loadings are more affected by accumulation periods longer than 9 months. The loadings of every station are shown in Figure 8

PCA Analysis
When performing the PCA, the first step is to determine the number of principal components retained. In the case of continuous series, two components explain 96% of the total variance. After applying the varimax rotation, the first principal component explains 67.5% of the total variance, while the second component explains an additional 28.5%.
As shown in Figure 7, PC1 mainly represents stations where short-time responses are the most important, while stations where PC2 has the highest loadings are more affected by accumulation periods longer than 9 months. The loadings of every station are shown in Figure 8    Principal component loadings show that there is a clear distinction between north and south in the Basin, with Guadalquivir River as the dividing line. Northern stations are generally represented by PC1 (Figure 8a), while southern stations, with some exceptions, present higher PC2 loadings (Figure 8b).
A possible explanation is that tributaries on the northern bank of the Guadalquivir are generally shorter than those on the southern bank. Therefore, shorter accumulation periods dominate in the northern area. On the contrary, the course of rivers in the southern part such as Guadiana Menor, Guadajoz, or Genil is longer, and large accumulation periods gain influence in this area.
The area belonging to Sierra Nevada, south of the basin, is also interesting. Even though there are no reservoirs directly on the mountains, the stations in this area are located at high altitudes, in a mountainous region which includes the highest peaks in the Iberian Peninsula. For Sierra Nevada and Sierra de Segura stations, high PC2 loadings are related to longer accumulation periods. This result is reasonable, since streamflow happens here mostly when snow melts. It is, therefore, realistic to assume that it is highly correlated to precipitation (or snow) from previous months. This result is confirmed with those presented in [42], even though accumulation periods in the latter are larger than in the Guadalquivir basin. However, it is explainable because of the higher altitudes of the studied basin in [42]. Principal component loadings show that there is a clear distinction between north and south in the Basin, with Guadalquivir River as the dividing line. Northern stations are generally represented by PC1 (Figure 8a), while southern stations, with some exceptions, present higher PC2 loadings (Figure 8b).
A possible explanation is that tributaries on the northern bank of the Guadalquivir are generally shorter than those on the southern bank. Therefore, shorter accumulation periods dominate in the northern area. On the contrary, the course of rivers in the southern part such as Guadiana Menor, Guadajoz, or Genil is longer, and large accumulation periods gain influence in this area.
The area belonging to Sierra Nevada, south of the basin, is also interesting. Even though there are no reservoirs directly on the mountains, the stations in this area are located at high altitudes, in a mountainous region which includes the highest peaks in the Iberian Peninsula. For Sierra Nevada and Sierra de Segura stations, high PC2 loadings are related to longer accumulation periods. This result is reasonable, since streamflow happens here mostly when snow melts. It is, therefore, realistic to assume that it is highly correlated to precipitation (or snow) from previous months. This result is confirmed with those presented in [42], even though accumulation periods in the latter are larger than in the Guadalquivir basin. However, it is explainable because of the higher altitudes of the studied basin in [42].
In addition to Sierra Nevada, there are other interesting results in the basin. Stations 5035 and 5042, both located close to each other and to the river estuary (Figure 1), do not follow the same trend as other stations in their vicinity. Station 5042 shows medium to In addition to Sierra Nevada, there are other interesting results in the basin. Stations 5035 and 5042, both located close to each other and to the river estuary (Figure 1), do not follow the same trend as other stations in their vicinity. Station 5042 shows medium to high loadings for both principal components since maximum correlations occur here at 11 months ( Figure 5), with similar although lower correlation values for 9 and 10 months. Following Figure 7, this is the accumulation period range in which PC1 score is surpassed by PC2 score, so it is logical that this station is represented by both components. However, station 5035 is exclusively represented by PC2, which implies longer accumulation periods despite the fact that the maximum correlation value in this station is reached for an accumulation period of 2 months ( Figure 5). This seems to contradict the score results, but it can be explained by looking at the correlation values of each accumulation period for this station in Figure 6. While 2 months is the maximum and the following short accumulation periods present lower correlation values, they rise again between 9 and 12 months, which are the accumulation periods for which PC2 has the maximum scores. Therefore, correlation values of station 5035 are better represented by PC2.
Stations 5019, 5025 (central-southern part of the basin), and 5032 (north-central part of the basin) are not represented by either of the two principal components. Moreover, according to Figure 5, SPEI accumulation periods for maximum correlation are the longest of the basin on these points. As stated before, these stations are particular because they are located in the largest catchment areas in the basin (Figure 2a), showing a different behavior that is not represented by the PC1 or PC2. For this reason, they will be analyzed separately.
In addition to the analysis of the continuous series, another PCA was performed on the twelve monthly series. These twelve PCAs were performed with two principal components each. The percentage of total variance represented by the first two rotated PCs ranges from 66% (September) to 98% (January). The specific value for each month appears in Figures 9 and 10 for PC1 and PC2, respectively. There are some stations where the streamflow dataset had monthly gaps, mostly during the summer months. These gaps are provoked most probably due to errors in the stream gauging stations, due to the low or inexistent streamflow during those months. For such stations, the results of the monthly PCA were not reliable enough, so they were marked with crosses in Figures 9 and 10, and excluded from this analysis. Figure 11 shows the scores of both principal components regarding the accumulation period. At monthly scale, it is clear that, except for May and December, PC1 is representative of stations with low accumulation periods, while PC2 is related to longer accumulation periods. This is similar to the result that appeared in the continuous series analysis. Figures 9 and 10 show the corresponding loadings for every month.
Generally, monthly principal components follow the same general behavior that appeared in the continuous series analysis. There are, however, two important exceptions: May and December. In May, PC1 represents longer accumulation periods, and PC2 denotes shorter accumulation periods. Moreover, spatial loading distribution is not as clear as in other months. Longer periods appear mostly in Sierra Nevada and the northwestern part of the basin, with shorter periods gaining importance only in some stations in the northeastern part of the basin. This result could be explained by considering that longer periods in May include snowmelt that happened during the spring, together with the normal precipitation increase that occurs in March and April [43]. For December, PC1 explains the south of the basin, which means that longer periods are dominant in this part, although with low scores. Shorter periods, represented by PC2, appear in the northern part of the basin. Since December is the month with the highest precipitation values [18], this result means that precipitation is converted into streamflow shortly after it happens, while mountainous areas accumulate snow that will turn into streamflow in future months. A similar result is discussed in [44], where nival regimes require longer accumulation periods, driven by snow-melt processes.
In order to analyze the differentiated behavior found between meteorological and hydrological droughts in detail, an exhaustive examination of the physical characteristics of the sub-basin has been carried out. To this end, it was studied how PCA results are impacted by the physical characteristics of each sub-basin where stations are located. For PC1 and PC2 separately, all stations with positive loadings are introduced in a multiple linear regression model together with the values of altitude, percentage of permeable soil, and catchment area. Figure 12 represents how well the stations represented by each PC are adjusted by the linear model. The line represents the perfect adjustment, so the horizontal distance between the line and each point characterizes the residuals. For PC1 (Figure 12a), the regression has R 2 = 0.71, with p-value = 5.05 × 10 −9 . Based on p-values for each variable, shown in Table 2, the only variable that seems to impact the loadings is the catchment area, while permeability and altitude are not statistically significant. This result confirms the hypothesis that stations 5019, 5025, and 5032, not shown in Figure 12, had to be studied separately; due to their exceptionally large catchment area, they can be considered as outliers in the regression analysis. These results are in accordance with the hypotheses studied by Heudorfer and Stahl [45] and van Langen et al. [46]. Water 2022, 14, x FOR PEER REVIEW 13 of 20 Figure 9. PC-1 loadings for the monthly series. Crossed stations indicate insufficient data for that monthly analysis. Crossed stations indicate insufficient data for that monthly analysis. Figure 10. PC-2 loadings for the monthly series. Crossed stations indicate insufficient data for that monthly analysis. Figure 11 shows the scores of both principal components regarding the accumulation period. At monthly scale, it is clear that, except for May and December, PC1 is representative of stations with low accumulation periods, while PC2 is related to longer accumulation periods. This is similar to the result that appeared in the continuous series analysis. Figures 9 and 10 show the corresponding loadings for every month. Generally, monthly principal components follow the same general behavior that appeared in the continuous series analysis. There are, however, two important exceptions: May and December. In May, PC1 represents longer accumulation periods, and PC2 denotes shorter accumulation periods. Moreover, spatial loading distribution is not as clear as in other months. Longer periods appear mostly in Sierra Nevada and the northwestern part of the basin, with shorter periods gaining importance only in some stations in the northeastern part of the basin. This result could be explained by considering that longer periods in May include snowmelt that happened during the spring, together with the normal precipitation increase that occurs in March and April [43]. For December, PC1 explains the south of the basin, which means that longer periods are dominant in this part, although with low scores. Shorter periods, represented by PC2, appear in the northern part of the basin. Since December is the month with the highest precipitation values [18], this result means that precipitation is converted into streamflow shortly after it happens, while mountainous areas accumulate snow that will turn into streamflow in future months. A similar result is discussed in [44], where nival regimes require longer accumulation periods, driven by snow-melt processes.
In order to analyze the differentiated behavior found between meteorological and hydrological droughts in detail, an exhaustive examination of the physical characteristics of the sub-basin has been carried out. To this end, it was studied how PCA results are impacted by the physical characteristics of each sub-basin where stations are located. For PC1 and PC2 separately, all stations with positive loadings are introduced in a multiple linear regression model together with the values of altitude, percentage of permeable soil, and catchment area. Figure 12 represents how well the stations represented by each PC are adjusted by the linear model. The line represents the perfect adjustment, so the horizontal distance between the line and each point characterizes the residuals. For PC1 (Figure 12a), the regression has R 2 = 0.71, with p-value = 5.05 × 10 −9 . Based on p-values for each variable, shown in Table 2, the only variable that seems to impact the loadings is the catchment area, while permeability and altitude are not statistically significant. This result confirms the hypothesis that stations 5019, 5025, and 5032, not shown in Figure 12, had to be studied separately; due to their exceptionally large catchment area, they can be considered as Figure 11. Scores of the two principal components obtained in the monthly series analysis for each accumulation period .  Concerning PC2, the regression results show R 2 = 0.73 and p-value = 1.9 × 10 −2 ( Figure  12b). According to the p-values of each variable (Table 2), this regression excludes permeability, while altitude is the most significant factor, followed by catchment area." Table 2. p-values of each variable and PC for the multiple linear regression analysis.  Concerning PC2, the regression results show R 2 = 0.73 and p-value = 1.9 × 10 −2 (Figure 12b). According to the p-values of each variable (Table 2), this regression excludes permeability, while altitude is the most significant factor, followed by catchment area."

Conclusions
The relationship between hydrological and meteorological droughts in the Guadalquivir basin has been analyzed using two drought indices, SSI and SPEI, through the correlation analysis.
Two main correlation types have been found in the basin. These types are short-time and long-time correlations. Short time scale correlations, with accumulation periods of maximum correlation ranging generally from 1 to 6 months, appear in the northern part of the basin, comprised mainly by flat areas, or short tributaries basins. In these areas, hydrological drought occurs shortly after the appearance of meteorological drought. In mountain areas, however, longer time scale correlations are present. Here, the snowmelt has a significant influence on the development of droughts. Accumulation periods above 9 months cause the highest correlations between SSI and SPEI. Additionally, there is a differenced behavior that happens in stations located downstream of another close reservoir. For these cases, since their streamflow is regulated, the correlation between streamflow and precipitation is nullified.
Through a PCA with varimax rotation, the basin has been divided into different areas with similar correlation patterns between the SSI and the 1-24 months accumulated period SPEI. The PCA was performed both on the continuous correlation series and separately on each of the monthly correlation series.
The results from the continuous series analysis showed that a first division can be made between North and South of the Guadalquivir River, with few stations showing a differentiating behavior from their neighboring stations. Short accumulation periods dominate the northern part, where tributaries are shorter and streamflow rates are lower, while longer accumulation periods have more influence in the south. The orography is also important in this basin, since accumulation periods presenting maximum correlations in mountainous areas on the southeastern part of the basin are different from those in flat areas. Generally, longer accumulation periods show higher correlation in the former, with the contrary happening in the latter.
The monthly PCA analysis shows that the behavior observed in the continuous series analysis is maintained during each separate month, except for May and December. In May, snowmelt and increased precipitations provoke that the spatial distribution that has been explained is not maintained. Meanwhile, in December, precipitation increases, so streamflows are formed in the flatter areas of the north with short accumulation periods, while the mountainous areas are represented by longer accumulation periods since snow will increase streamflows in later months.
For stations presenting low accumulation periods, catchment area is the most important parameter in order to describe the correlations between PCs and the values of altitude, percentage of permeable soil, and catchment area. Nevertheless, altitude is the variable with highest influence for long accumulation periods. This is confirmed by the spatial distribution of stations in the PCA, with PC2 mostly describing points in mountainous areas.
Results from this study can be useful to understand the mechanism behind hydrological droughts, and the importance of river regulation in its development. However, in this study, there are no data available in the area close to the mouth of the Guadalquivir, to the