Teleconnection analysis of runoff and soil moisture over the Pearl River basin in southern China

This study explores the teleconnection of two climatic patterns, namely the El Niño–Southern Oscillation (ENSO) and the Indian Ocean Dipole (IOD), with hydrological processes over the Pearl River basin in southern China, particularly on a sub-basin-scale basis. The Variable Infiltration Capacity (VIC) model is used to simulate the daily hydrological processes over the basin for the study period 1952–2000, and then, using the simulation results, the time series of the monthly runoff and soil moisture anomalies for its ten sub-basins are aggregated. Wavelet analysis is performed to explore the variability properties of these time series at 49 timescales ranging from 2 months to 9 yr. Use of the wavelet coherence and rank correlation method reveals that the dominant variabilities of the time series of runoff and soil moisture are basically correlated with IOD. The influences of ENSO on the terrestrial hydrological processes are mainly found in the eastern sub-basins. The teleconnections between climatic patterns and hydrological variability also serve as a reference for inferences on the occurrence of extreme hydrological events (e.g., floods and droughts).


Introduction
Regional responses of hydrological processes to climate variability are often different from one river basin to another (Eltahir, 1996).Accordingly, the different sub-basins of a large-scale river basin may play their own roles in modulating the influences of climatic patterns (Chen andKumar, 2002, 2004;Marković and Koch, 2005) and in manifesting the tendency of flood or drought occurrence due to extreme precipitation.Therefore, an important step in the assessment and prediction of floods or droughts is to understand the relationships between the key hydrological variables (e.g., precipitation, runoff, soil moisture) and climate patterns.
Generally speaking, observations of hydrological variables over large-scale basins are scarce, except for precipitation.For instance, even though streamflow gauging stations currently exist along most large rivers around the world, long-term flow observations are still not available, especially for sparsely populated regions and complex terrains.In addition, while soil moisture is a key component of the hydrological cycle, it is also one of the least observed hydrological variables in terms of spatial and temporal coverage.Because of the lack of long-term and large-scale observations of runoff and soil moisture, utilization of hydrological variables forced by observed precipitation and other near-surface meteorological data is considered as a sensible and viable alternative in hydrological and water resources studies.For instance, using a macroscale hydrological model, the Variable Infiltration Capacity (VIC) model, Niu andChen (2009, 2010) simulated the daily terrestrial hydrological processes over the Pearl River basin in southern China, and routed the simulated runoff to the downstream control stations of three major tributaries (i.e., the West River, the North River, and the East River).
In the specific context of the Pearl River basin, many other recent studies have also helped advance our understanding of the hydrological variability over the basin and the surrounding regions, including those on precipitation (e.g., Cui et al., 2007;Yang et al., 2010;Fischer et al., 2012;Chen et al., 2011;Zhang et al., 2012) and runoff (S.R. Zhang et al., 2008;Niu and Chen, 2010), among others.For instance, several studies have reported that the significant fluctuations in precipitation over the basin are associated with atmospheric circulation patterns (Li and Zhou, 2012;Niu, 2013) influencing the monsoon system over the region (Zhang et al., 2010;Zhou et al., 2012).Zhang et al. (2013) examined the meteorological drought tendencies based on the standardized precipitation index by employing copula functions and the Mann-Kendall trend test.The Pearl River basin was categorized into six regions based on latitude, longitude, elevation, and the mean annual precipitation using the L-moments approach (Yang et al., 2010), and four homogenous regions based on the changing properties of precipitation using the k mean clustering algorithm (Zhang et al., 2012).Zhang et al. (2009) analyzed the scaling and persistence features of the long, daily streamflow series for four hydrological stations on the mainstream of the East River.The changes of observed water discharge and sediment load (1950s-2004) were analyzed by S. R. Zhang et al. (2008) for nine stations in the Pearl River basin, although performed only on annual timescale.Despite these developments, it must be recognized that these studies were mainly performed based on observed data, the exceptions being those of Niu andChen (2009, 2010).
Following up on the studies by Niu andChen (2009, 2010), the present study evaluates the sub-basin hydrological processes over the Pearl River basin in response to two large-scale climatic patterns, the El Niño-Southern Oscillation (ENSO) and the Indian Ocean Dipole (IOD), in order to improve our understanding of basin-featured flood and drought occurrences.Compared with the previous studies over the Pearl River basin, our focus in this study is on runoff and soil moisture processes, particularly on a sub-basin scale basis.Although runoff and soil moisture processes are generally driven by precipitation, details of the analysis of precipitation on sub-basin scale over the basin are not presented here.For such details, the interested reader is directed to the studies by Niu (2013), which present extensive investigations on the precipitation process over the basin as well as its teleconnections with climatic patterns.
In this study, we attempt to address the following specific question: are there some characteristic periods and variability processes of different hydrological variables within the Pearl River basin that may be related to climatic patterns, which may be a possible trigger for major floods or droughts (i.e., extreme hydrological events) in the basin?We address this question through (1) investigation of the possible variability linkages between the runoff/soil moisture and the above two climatic patterns, and (2) discussion of the underlying causes for extreme hydrological events from this perspective.For analysis, we consider data over the period 1952-2000, similar to the period studied by Niu andChen (2009, 2010).To identify the patterns and capture the differences in continuous timescales of fluctuation embedded in the hydrological time series (i.e., streamflow and soil moisture), we employ wavelets, in particular continuous wavelet transform.
Torrence and Compo (1998) developed a significance test for the wavelet power spectrum to provide quantitative results and reliable insights, and make the applications of continuous wavelet transform extensive.Subsequently, Grinsted et al. (2004) discussed the cross-wavelet transform and wavelet coherence for the applications to geophysical time series in examining relationships between time series in a time-frequency domain.Meanwhile, Maraun and Kurths (2004) suggested that the interpretation of nonnormalized cross spectra for the interrelation between two processes is not straightforward if not including appropriate confidence intervals, as it may show misleading peaks when just one of the time series has a strong peak.Further, Schaefli et al. (2007) reviewed recent developments in wavelet spectral analysis, recommended the calculation of wavelet coherence in the applications, and provided a valuable example to illustrate the usefulness of wavelet spectral analysis for detecting the flood-triggering situations in the Alpine catchment.In addition, the importance of the pervious knowledge about the physical processes and a methodologically rigorous significance testing are presented for the proper application and interpretation of the wavelet spectral analysis (Maraun and Kurths, 2004;Schaefli et al., 2007;Maraun et al., 2007).
The rest of this paper is organized as follows.Section 2 gives details of the Pearl River basin study area and data considered for analysis.Section 3 describes the methods, including Morlet wavelet transform, the wavelet coherence and Spearman's rank correlation.Section 4 presents the results and discussions from the wavelet application to runoff and soil moisture data in the study basin.Section 5 discusses the underlying causes of major flood/drought events in the basin.Conclusions are drawn in Sect.6.

The Pearl River basin
The Pearl River basin in southern China has a total area of about 450 000 km 2 , and consists of four main components: the East River, the North River, the West River and the Pearl River delta.The basin is located in a tropical and subtropical monsoon climate region (Pearl River Water Resources Commission, 2005), with a long-term average annual precipitation of about 1480 mm yr −1 and mean annual temperature ranging from 14 (in the western part) to 22 • C (in the eastern part).As the largest river in terms of the rate of annual runoff generation (about 740 000 m 3 km −2 yr −1 ) among the seven major rivers in China (i.e., the Yangtze River, the Yellow River, the Songhua River, the Pearl River, the Huai River, the Hai River, and the Liao River), the Pearl River experiences an abundant amount of variability in the terrestrial hydrological processes.For instance, the runoff ratio of the  wettest year to the driest year is as high as 6-7, and about 80 % of the total annual runoff occurs in the wet season of April-September (Niu and Chen, 2009).As a consequence, extreme hydrological events frequently occur over the basin (Pearl River Hydraulic Research Institute, 2007;Cui et al., 2007).For example, a basin-scale flood occurred in 1994 and a severe drought event in 1963.
With due consideration to hydrological factors, water resources management, topographic features, and climate zones, this study divides the basin into ten sub-basins for analysis, similar to the approach by Niu (2013).These subbasins are the East River, the North River, and the eight sub-basins in the West River basin (i.e., the Guihe, Lijiang, Hongshui, Beipan, Nanpan, Youjiang, Zuoyu, and the West River's lower reach) (see Fig. 1).The sub-basin areas range from 29 000 to 58 600 km 2 , and the mean annual precipitation varies from 1072 to 1772 mm yr −1 .

Hydrological model
To obtain runoff and soil moisture data over the Pearl River basin, we adopt the VIC hydrological model to simulate such data.The VIC model (Liang et al., 1994) is a macroscale hydrological model, usually applied with spatial resolutions ranging from 0.125 to 2 • latitude/longitude and from hourly to daily temporal resolutions (Nijssen et al., 2001a, b;Maurer et al., 2002;Sheffield et al., 2004;Troy et al., 2008).The model uses the variable infiltration curve (Zhao et al., 1980), which allows for the representation of subgrid variability in soil storage capacity and infiltration capacity, and the parameterization of baseflow as a nonlinear recession curve from the bottom soil layer (Dümenil and Todini, 1992;Todini, 1996).The land surface of the study area is represented by a specified number of tiled land cover classes (e.g., cropland, grassland) and the soil column is divided into three soil layers.The total runoff of each grid cell is the sum of the overland runoff and baseflow obtained for each grid cell (Nijssen et al., 1997).

Forcing and output data
The terrestrial hydrological processes over the Pearl River basin are simulated by the VIC model for the period 1952-2000, with 1 • × 1 • spatial resolution at daily time steps (Niu and Chen, 2010).The GTOPO30 DEM (digital elevation model) data set with 1 km spatial resolution (Verdin and Verdin, 1999) is used to delineate the Pearl River basin.Meteorological forcings, including daily precipitation, maximum and minimum temperature, wind speed, and relative humidity, reported in Feng et al. (2004), are used.The vegetation and soil data are extracted from two global data sets provided by Nijssen et al. (2001a).The vegetation data are derived from the 1 km vegetation database of the University of Maryland, on the basis of Advanced Very High Resolution Radiometer (AVHRR) data (Hansen et al., 2000).The soil texture and soil bulk density are derived from the 5 min FAO-UNESCO (Food and Agriculture-United Nations Educational, Scientific, and Cultural Organization) digital soil map of the world and the WISE (World Inventory of Soil Emission Potentials) Pedon database; the remaining soil parameters are obtained based on the study of Cosby et al. (1984).The thicknesses of the three soil layers over the Pearl River basin are 0.3 m (the upper layer), 1.0-2.5 m (the second layer), and 0.2 m (the bottom layer), respectively.The VIC-model's daily simulation results are compared with the observations at a monthly scale after aggregation of daily data.The available streamflow observations at 10 gauging stations (i.e., Gaoyao, Shijiao, Boluo, Zhedong, Baise, Nanning, Liuzhou, Wuzhou, Xinfengjiang, and Longchuan) in the basin (see Fig. 1) are used (Niu and Chen, 2010).The monthly scale is chosen for (1) the more "acceptable" accuracy of the model simulation results, and (2) the compatibility with the large-scale climate indices data.Table 1 presents a comparison between the observations and the model simulations of streamflow data for these 10 stations, by employing three objective functions (i.e., relative bias -RB, relative root mean square error -RRMSE, and Nash-Sutcliffe efficiency -NSE -coefficient).The comparison suggests that the VIC model can simulate reasonably well the streamflow from the North River and the West River with its several sub-basins.However, there is a considerable discrepancy between the streamflow simulations and the observations for the Boluo station in the East River.The studies by Niu and Chen (2010) and Wu and Chen (2012) have revealed that the Boluo streamflow is heavily regulated by the reservoirs in the upstream of the station.Therefore, the observations without the reservoir influences from two upstream stations (Xinfengjiang and Longchuan) are also used herein to evaluate the streamflow simulations in the East River basin.The new results for the East River basin are comparable to the results for the other stations in the North or West River basins (see Table 1).Based on the performance ratings suggested by Moriasi et al. (2007), the streamflow simulations at a monthly time step in the Pearl River basin in this study are evaluated as "satisfactory" (NSE > 0.5 and RB < ±0.25).With the validation of the VIC model performance, the monthly time series of runoff and soil moisture are further aggregated for the 10 sub-basins (see Fig. 1 for their locations), and the waveletbased analysis is carried out for the monthly runoff and soil moisture anomalies for the sub-basins for the study period 1952-2000.It is also noted that the routing effects are effective within 2-3 weeks for the Pearl River basin (Niu and Chen, 2010), and thus the monthly scale we choose, as the basis of the teleconnection study, is larger than the effective timescale of the routing effects.Table 2 lists the basic statistical properties for the observed precipitation and simulated runoff and soil moisture over the 10 sub-basins in the basin.

Large-scale climatic patterns
Two prominent large-scale climatic patterns influencing southern China (Saji et al., 1999;Guan and Yamagata, 2003;Zhou and Chen, 2007;Zhang et al., 2007), ENSO and IOD, are selected for examining the influence of largescale climatic patterns on the hydrological processes in the Pearl River basin.The ENSO indices used are obtained from the monthly data of the National Oceanic and Atmospheric Administration (NOAA) Niño 3.4 (5 • N-5 • S, 170-120 • W) SST (sea surface temperature; Trenberth and Stepaniak, 2001;Labat, 2010a) for the base climatological period 1952-2000.The intensity of the IOD is characterized by an anomalous SST gradient between the western equatorial Indian Ocean (10 • S-10 • N, 50-70 • E) and the southeastern equatorial Indian Ocean (10 • S-0 • N, 90-110 • E).The IOD indices used are the combination of monthly IOD indices  from the Hadley Centre Global Sea Ice and Sea Surface Temperature (HadISST) data set and monthly IOD indices (1997)(1998)(1999)(2000) from the NOAA Optimum Interpolation (OI) SST data set (Saji et al., 1999).

Morlet wavelet transform
As the continuous wavelet transform is suitable for analyzing nonstationary signals (such as streamflow and soil moisture) and identifying the temporal localization of dominant events, such wavelets are considered here.In particular, the Hydrol.Earth Syst.Sci., 18, 1475-1492, 2014 www.hydrol-earth-syst-sci.net/18/1475/2014/ ) is obtained from the monthly raw time series for its corresponding hydrological variable.The values of σ 2 (variance) (unit: (mm month −1 ) 2 ) and α (lag-1 coefficient of the assumed red-noise process) are computed from the monthly anomaly time series obtained by removing the respective 49 yr monthly mean.
continuous Morlet wavelet transform is performed, as (1) it offers a good balance between time and frequency on the interpretation of variability properties; and (2) its complex wavelet function consists of more oscillation waves, thus, making it capable of capturing the oscillatory characteristics.
The Morlet wavelet consists of a complex exponential wave modulated by a Gaussian function.
To realize the continuous wavelet transform (for a given wavelet scale s) at all time index n simultaneously and efficiently (Torrence and Compo, 1998), the wavelet transform can be calculated as where "hat" denotes the Fourier transform, * indicates the complex conjugate, k is the frequency index and ω k is the angular frequency equal to ± 2 π k N δt .For convenience, the scales are written as fractional powers of two: where s 0 is the smallest resolvable scale and J determines the largest scale.In this study, for the period 1952-2000, the values of δt = 1/12 yr, s 0 = 2 δt, δj = 0.12, and J = 48 give a total of 49 scales, ranging from 1/6 yr (2 months) up to 9 yr.The wavelet power spectrum is defined as W n (s j ) 2 and gives a local measure of the time series variance at time and frequency domains.The total energy is conserved under the wavelet transform, and the variance of the analyzed time series is given as where C δ is a reconstruction factor specific to different wavelet functions (C δ = 0.776 for the Morlet) (Torrence and Compo, 1998).
The time-averaged wavelet power spectrum over all the local wavelet spectra is defined as a global wavelet spectrum (GWS): The chi-square (χ 2 ) test is performed to evaluate the significance of the timescale observed by the GWS.The scaleaveraged wavelet power spectrum is the weighted sum of the wavelet power spectrum over scales j 1 to j 2 : In addition, to display the variance evolution on the dominant variability band, the band-pass reconstruction is performed.
The reconstructed time series is the sum of the real part of the wavelet transform (Torrence and Compo, 1998) over the above band: where the factor ψ 0 (0) removes the energy scaling, Re W n (s j ) is the real part of the transform, and s 1/2 j converts the wavelet transform to an energy density.

Wavelet coherence
The wavelet coherence analysis is used to identify regions in a time-frequency domain where the two time series covary, but may not have a common high power (Grinsted et al., 2004).The software package on wavelet coherence by Grinsted et al. (2004) is employed here to provide both correlation analysis and rough information on the phase difference between two time series in a time-frequency domain.For two given time series, x n and y n , with corresponding wavelet transforms W X n (s) and and the wavelet coherence of two time series is defined as where S is a smoothing operator, defined as (Grinsted et al., 2004): where S scale is the smoothing along the wavelet-scale axis and S time is the smoothing in time.A suitable smoothing operator for the Morlet wavelet is given by Torrence and Webster (1999) as follows: where c 1 and c 2 are normalization constants, and is the rectangle function.The value 0.6 represents the decorrelation length of the empirically determined scale for the Morlet wavelet.The statistical significance level of the wavelet coherence is estimated based on a large ensemble of surrogate data set pairs with the same lag-1 autoregressive (AR1) coefficients given by Grinsted et al. (2004).The global wavelet coherence spectrum can be obtained by averaging the local wavelet coherence spectrum in time, and is helpful to examine the characteristic scales (Torrence and Compo, 1998;Labat, 2010a).

Rank correlation
The Spearman's rank correlation coefficient r S is applied to measure the association between the ranks of two time series of scale-averaged wavelet power (x n and y n ).The rank correlation analysis is used here because the transformed waveletpower time series, for either an individual timescale or a timescale band, normally does not follow Gaussian distribution (Torrence and Webster, 1999).The score d is calculated after ranking both x n and y n , which is the difference between the rank on x n and the rank on y n .
The significance test of the Spearman's correlation has t and z versions as where significance statistic t o has (n − 2) degrees of freedom, and where significance statistic z o has (n − 1) degrees of freedom.Following the now-familiar rules (Asquith, 2008), the t test is used when the cases or paired data points are fewer than 30, while z is used for 30 cases or more.

Wavelet power spectrum
In this study, ten sub-basins over the Pearl River are analyzed to examine the relationships between large-scale climatic patterns and local hydrological processes.Among these ten sub-basins, we consider the Boluo station in the East River for the purpose of illustration.
For the Boluo station, the monthly runoff anomaly time series is shown in Fig. 2a-i, and the normalized waveletpower spectrum, |W n (s)| 2 /σ 2 , is shown in Fig. 2a-ii; the horizontal axis is the time position over the period 1952-2000.In Fig. 2a-i, the anomaly data is obtained by removing the monthly mean for the period 1952-2000.In Fig. 2a-ii, the vertical axis has 49 timescales of the wavelet transform ranging from 2 months to 9 yr, and the strong variability power with greater than 95 % confidence is highlighted using white contours.We can observe that many abnormal runoff events (i.e., monthly runoff being too high or too low) appear at short timescales (e.g., less than 1 yr).Most of the wavelet power is distributed at the 1-4 yr band, and the strong variability power, shown as white-contour enclosed regions, is found during the periods 1955-1970, 1980-1985 and 1990-1995.ii. iii. i. ii.
iii.  Figure 2a-iii shows the global wavelet power spectrum (GWS) of the monthly runoff anomaly data (see Eq. 5).As many geophysical time series can be modeled as red noise (Torrence and Compo, 1998), the figure also includes the mean and 95 % confidence level red-noise spectrum (left and right dashed lines) by using the lag-1 autocorrelation.It is observed that, over the 49 timescales, the highest runoff variability is identified at the 1.92 yr scale, which is above the 95 % confidence level.This is consistent with the result obtained for precipitation data in the East River basin (Niu, 2013), where the highest precipitation variability is detected around the 2 yr band.The detected significant variability indicates less randomness of the underlying processes due to local influences, and it may well be related to external climatic patterns (see Sect. 4.4 for details).
Figure 2b shows the corresponding results for the monthly soil moisture in the East River basin.A comparison between Fig. 2a-ii and b-ii reveals that the high variability power at short timescales for soil moisture almost disappears, while for larger timescales, the related power for soil moisture is enhanced.In the terrestrial hydrological processes, the river basin, like a buffer, reduces the high-frequency precipitation variability by runoff, and the long-term precipitation variability is memorized by soil moisture.Figure 2b-iii shows that the variability power at short timescales is reduced much more, but improved at larger timescales (note that the magnitude of the GWS is doubled compared to that in Fig. 2a-iii).Further, the detected highest variability is transferred from 1.92 (see Fig. 2a-iii) to 3.16 yr (see Fig. 2b-iii), with the large spikes (about 6 spikes) in Fig. 2a-ii appearing around the 2 yr scale, and the spikes with the largest area (in a timefrequency domain) exist around the 2 yr scale in Fig. 2b-ii.

Coherent modes of the GWSs
In addition to the analysis of the East River presented above, we further inspect their global variability distributions (i.e., GWSs) for the other nine sub-basins.Performing principal component analysis (PCA) on the GWSs over these 10 sub-basins, we find that there are three coherent modes for hydrological variabilities at the 49 timescales (see Fig. 3).The main steps in the PCA application in this study are as follows: (1) a covariance matrix is obtained by using the wavelet power spectra for s (equal to 49) timescales of runoff/soil moisture in 10 sub-basins; (2) the covariance matrix is decomposed to get eigenvectors of the matrix, which accounts for the maximum amount of the joint variability of anomalies of initial wavelet power spectrum at each sub-basin; (3) the new variables are produced after projecting the data on the obtained eigenvectors; (4) the first three principal components are selected for runoff and soil moisture, respectively; and (5) each sub-basin is then classified as belonging to the mode where its coefficient of the eigenvectors (i.e., the length of the vector projected on corresponding direction) is the largest absolute value.
Over the Pearl River basin, there are mainly three types of precipitation events: frontal type storms (during April-June), tropical cyclonic rains (during July-October), and convective rainstorms (during June-August) (Niu, 2013).Therefore, the rainy season essentially lasts for 7 months in a year, from April to October (Niu, 2010).This high variability of precipitation is partly transferred to the high runoff variability at the timescale of 0.59 yr (see Fig. 3a).Apart from the high variability for less-than 1 yr periods, the variability in the dominant low-frequency variability band (DLVB), 0.83-5.21yr, demonstrates the distinct features of the three multiscale variability modes in the Pearl River basin (see Fig. 3a).The DLVB is determined as follows: the local dominant variability band is defined as the longest projecting portion with its peak spectrum at least above the mean red-noise spectrum and most likely to reach a higher confidence level (e.g., 95 % level); and the DLVB is the band that covers most of the identified dominant variability bands of the 10 sub-basins.The three coherent modes of soil moisture in Fig. 3b show that the variability at longer timescales (longer than 0.59 yr) has the striking features within the DLVB range of 1.92-5.21yr.
Figure 4a and b shows the sub-basin clusters of multiscale variability of runoff and soil moisture, respectively; the three clusters correspond to the three coherent modes.The cluster of runoff is basically consistent with that of precipitation (see Niu, 2010), while the cluster of soil moisture displays dissimilar results.One rainstorm center in the Pearl River basin is located in the lower reach of the Hongshui River basin (Pearl River Water Resources Commission, 2005), and this results in the highest power at the 0.5-1 yr period for the third mode of precipitation and runoff in the central part.The Nanpan River basin is in the northwest region of the Pearl River basin.This region has the highest altitudes, and is influenced  by the topographic rain shadow with respect to the prevailing storm tracks.The east region of the Pearl River basin is close to the South China Sea, and is subjected to convective movement of water by semitropical typhoons (Niu, 2013).The higher variability in the 1-4 yr band is more obvious, compared to the third mode of runoff in the central part of the Pearl River basin.
The Beipan, Nanpan, and Youjiang river basins, located in the northwestern region with the highest altitudes, belong to the first mode of soil moisture that displays a single peak of high variability at longer timescales.In these basins, the high-frequency precipitation variability is kept less by runoff and reflected more by soil moisture at low frequencies.The second mode of soil moisture, including the North, Guihe, and Zuoyu river basins and the West River's lower reach, shows relatively flat variability power distribution (see Fig. 3b), which either has lower altitudes or is located at the downstream of the Pearl River.Meanwhile, the third mode of soil moisture, to which the Liujiang, Hongshui and East river basins belong, sits between the above two distinctive modes, which is partly consistent with elevation variations between the other two modes.

Evolution of bandpass variance
To highlight the variability processes with respect to time for the coherent mode of GWS, the variance evolution in the dominant band is inspected by using the scale-averaged wavelet power and bandpass-filtered time series.Figure 5 illustrates the flowchart of unfolding regional variance evolution of hydrological variables.We use Eq. ( 6) to compute the scale-averaged wavelet power, and Eq. ( 7) to generate the bandpass-filtered series.The representative of each mode is the time series average of the independent sub-basins that belong to the same mode.
Figure 6 shows the runoff series of both scale-averaged wavelet power and corresponding bandpass-filtered time series for three coherent modes between the timescales of 0.83 and 5.21 yr.The second mode (i.e., the East River, North River, Guihe River, and Liujiang river basins and the WR's (West River's)lower reach) shows larger variance processes in terms of the magnitude of average variance (2508 mm 2 ), with high variance occurring sometime around 1963, and also two other periods from 1983 to 1984 and from 1991 to 1996.The variability processes for the third mode (i.e., the Hongshui, Beipan, Youjiang and Zuoyu river basins) show an irregular cycle of an increasing amplitude of runoff events followed by a decreasing amplitude.The approximately 5 yr cycle is also clearly apparent from 1967 to 1977 and from 1987 to 1997.
Figure 7 shows the variance evolution of soil moisture within 1.92-5.21yr.Compared to the results obtained for runoff, the soil moisture variability processes are relatively smooth along the time and relatively balanced in terms of the magnitude of average variance for the three modes (i.e., 792, 1151 and 911 mm 2 ).A large variance can be identified for the first mode (i.e., the Beipan, Nanpan and Youjiang river basins) for the period 1955-1970.The variance processes for the second mode (i.e., the North, Guihe and Zuoyu river basins and the WR's lower reach) are relatively gentle, while the third mode (i.e., the East, Liujiang and Hongshui river basins) shows noticeable changes with two high-variance periods: 1960-1970 and 1995-2000.

Teleconnection to climatic patterns
To examine the possible linkages of climatic patterns with variations in runoff and soil moisture in the Pearl River basin, their wavelet coherences are inspected.The results are presented below.

Runoff
Figure 8 shows the local and global wavelet coherence of the three modes of runoff (modes one to three, from top to bottom) in the Pearl River basin with ENSO (left panels) and IOD (right panels).It is observed that, for the first-mode runoff and ENSO (Fig. 8a, left panel), there are significant regions around the 1 yr band within the periods 1970-1975 and 1990-1995, and also for the 2-4 yr band around 1995.The vectors in the regions outlined by thick contours indicate the close in-phase relationship around 1995 (for 2-4 yr band), which is also consistent with the rank correlation results in Sect.4.5.These observations seem to indicate a link (i.e., close in-phase variance change) between the firstmode runoff and ENSO in the 2-4 yr band around 1995.The right panel in Fig. 8a shows the correlations of the first-mode runoff and IOD.The noticeable in-phase relations are found within the 2-4 yr band around 1995.
The runoff in the eastern part of the Pearl River basin belongs to the second coherent mode.Its correlation with ENSO is shown in Fig. 8b (left panel).The identifiable correlations (i.e., the time-frequency region is with significant wavelet coherence) are within the 1-4 yr band during the 1970-1975 period and around 1983, and within the 2-4 yr band around 1995.The phase relations indicate ENSO leading the second-mode runoff during the 1970-1975 period.The in-phase relation is identified within the 2-4 yr band around 1995.For the correlations between the secondmode runoff and IOD in Fig. 8b (right panel), the identifiable correlations are mainly within the 0.5-2 yr band around 1983 (anti-phase) and within the 1-4 yr band around 1995 (in-phase).
The left panel in Fig. 8c shows that the identifiable correlation between the third-mode runoff and ENSO exists around the 4 yr band during the period 1985-1995.The vectors in the regions show close in-phase relations of variance processes.Figure 8c (right panel) shows the main identifiable correlations between the third-mode runoff and IOD, with anti-phase relations within the 1-2 yr band during the period 1970-1975 and in-phase relations within the 2-4 yr band during the period 1990-1995.The global wavelet coherence results in Fig. 8 also show the oscillation correlation around the 1-6 yr band, especially between the second and third modes of runoff and IOD (the right panels in Fig. 8b and c).

Soil moisture
The wavelet coherence between the first-mode soil moisture and ENSO in Fig. 9a (right panel) discloses high coherence at both short and relatively longer periods.The  for the first, second and third modes are 263.26mm 2 , 2,507.64 mm 2 and 1,096.32mm 2 , respectively.Fig. 6.The scale-averaged (0.83-5.21 yr) wavelet power for runoff (thick black curve) and the bandpass filtered time series of square form (runoff squared) (thin curve) from 0.83 to 5.21 yr using a wavelet filter.The upper, middle and lower panels are the first, second and third modes of runoff, respectively.The black shading indicates positive peaks in the filtered time series, while the white are negative peaks.The wavelet power and bandpass filtered time series (before squaring) are the averaged time series over the sub-basins belonging to the same mode.The averaged variances of initial monthly runoff time series σ 2 for the first, second and third modes are 263.26,2507.64 and 1096.32 mm 2 , respectively. for the first, second and third modes are 792.45mm 2 , 1,151.50 mm 2 and 911.17 mm 2 , respectively.Fig. 7. Same as Fig. 5, but for soil moisture in the 1.92-5.21yr band.The averaged variances of initial monthly soil moisture time series σ 2 for the first, second and third modes are 792.45,1151.50 and 911.17 mm 2 , respectively.major identifiable correlations are distributed within the 1-2 yr band (around 1965 and 1973, and during the period 1990-1995), the 2-4 yr band (around 1995), and the 4 yr band (during the period [1980][1981][1982][1983][1984][1985][1986][1987][1988][1989][1990].The results also show close in-phase relationships for the regions during the period 1980-1990 and around 1995.The major identifiable correlations between the first-mode soil moisture and IOD in Fig. 9a (right panel) are found within the 1-2 yr band around 1975 (anti-phase), and within the 2-4 yr band around the period 1990-1995 (close in-phase).
The left panel in Fig. 9b shows that the identifiable correlations between the second mode of soil moisture and ENSO,  The identifiable correlations between the third mode of soil moisture and ENSO in Fig. 9c (left panel) can be found within the 1-2 yr band during the period 1960-1975 and the 1-4 yr band during the period 1990-2000, in which the vectors indicate that ENSO leads soil moisture in the former region and close in-phase relations in the latter region.The right panel in Fig. 9c shows that the major identifiable correlations are within the 0.5-2 yr band during the period 1980-1985 (anti-phase) and within the 1-4 yr band during the period 1990-2000 (close in-phase).The global wavelet coherence in Fig. 9b and c shows the teleconnection peaks between regional (the eastern and central parts) soil moisture and climatic patterns at longer timescales.

Correlation at multiple timescales
As indicated above, the significant correlations of runoff and soil moisture with the two climatic patterns (ENSO and IOD) mainly exist within the DLVB.To highlight the teleconnection influences further, the correlations of wavelet power time series between a hydrological variable and a climatic index, within the corresponding DLVB, are computed by using Spearman's rank correlation.A positive correlation value indicates that the amplitudes of variance over the desired time period have consistent evolution processes, while a negative correlation value corresponds to opposite evolution processes (Niu, 2010).
Figure 10 shows the scale-averaged wavelet power for each runoff mode within the 0.83-5.21yr band and the two climatic indices (i.e., ENSO and IOD), and Table 3 lists the corresponding correlation results.It is observed that the first and third modes of runoff correlate with ENSO and IOD occasionally for the divided 10 yr periods.The second mode of runoff is basically depicted by IOD, which is reflected by the significant correlations for all the divided periods and  also by the high correlation value of 0.643 for the complete time period .The opposite evolution processes are found for ENSO and IOD for the period 1960-1969 (see Fig. 10), which, in turn, produces positive correlations (0.295 and 0.744) and negative correlations (−0.417 and −0.501) with the first and second modes of runoff, respectively.This also confirms the different temporal processes of variance for the two climatic patterns (Niu, 2013).Figure 11 indicates a good relationship between each mode of soil moisture and IOD within the 1.92-5.21yr band for the entire period, with the rank correlation values of 0.669 (for the first mode), 0.319 (for the second mode) and 0.832 (for the third mode), as presented in Table 3. Significant correlations are also found between soil moisture and ENSO for most of the divided 10 yr periods, but with lower correlation degree for the first mode of soil moisture (i.e., in the western part of the basin) and higher for the second mode (i.e., in the eastern part of the basin) after the 1970s, compared to the corresponding correlation level between soil moisture and IOD.

Inferences on droughts and floods
To illustrate the potential utility of the above analysis, possible inferences of the results as to the occurrence of basin-featured extreme hydrological events (i.e., the flood event in 1994 and the drought event in 1963) are explored.
Two Asian monsoon systems, namely the Indian summer monsoon and the East Asian summer monsoon, strongly influence eastern China (Ding and Chan, 2005).Torrence and Webster (1999) discussed the relationship between the Indian monsoon and ENSO, through analysis of a 122 yr-long data set .They reported that the Indian monsoon tends to be weaker during the El Niño event periods and that the weak monsoon statistically leads by about 4 months to the peak of equatorial Pacific SST.A strong monsoon is associated with cold SST, and its peak lags the monsoon by about 4 months (Torrence and Webster, 1999).Even though investigating the mechanisms and relationships between climatic indices and monsoon systems is out of the scope of the current study, the results reported by Torrence and Webster (1999)

Flood event
In 1994, a major flood occurred over almost the entire Pearl River basin.This flood event offers a reasonable basis to interpret the results from the wavelet analysis in the context of relationships that exist between large-scale climatic patterns and basin-featured hydrological processes.The significant common wavelet power and coherence in Fig. 8a and c show strong in-phase relations between runoff and IOD within the 2-4 yr band in 1994, but the correlation of runoff with ENSO is not that visible; it is relevant to note that there was a strong positive IOD event in 1994.Further, in 1994, the moisture transport resulting from the southeastern tropical Indian Ocean gave rise to an enhanced convection over the Indian subcontinent and East Asia (Behera et al., 1999).
Therefore, heavy rainfall and flood occurred over southern China, including in the Pearl River basin.
The impactful occurrence of small-scale activity (variability) is indeed often supported by larger-scale activity (see Kumar, 1996 for details).Flood events normally occur in response to abundant precipitation over a short period of time, especially when the land has already been in a wet state (runoff or soil moisture) over a long time period.The abundant precipitation in 1994 occurred during an already wet state that existed due to high runoff in the basin before 1994, including due to high runoffs in the central part in 1992 and 1993 (shown in Fig. 6b) and in the western part in 1993 (shown in Fig. 6c), thus resulting in an extreme flood over the whole basin in 1994.Although abundant precipitation also occurred around 1983, it did not cause any major flood since the land surface was not as wet during the period before 1983 as it was before 1994, as the second mode of runoff indicates.Further, the soil moisture condition before 1983 was also in a long-term dry state (see Fig. 7), especially for the first and the third modes of soil moisture.Such relatively dry runoff and soil moisture states play a significant role in lowering and slowing down the runoff over the basin and, consequently, the resulting flow would not be sufficient to cause major floods; the fact that there was no flood in 1983 despite abundant precipitation explains this.

Drought event
In 1963, about 5 % of the total area of the Pearl River basin was subjected to a severe drought; this drought also provides a reasonable basis for interpreting the results from the wavelet analysis above.As seen in Fig. 9, the soil moisture (especially for the second and third modes) has significant wavelet coherence at > 95 % confidence levels with the ENSO, but not IOD, within the 1-2 yr band around 1963; ENSO was in a weak warm phase (El Niño) in that year.Statistically, a weak Indian monsoon may be related to an El Niño period (Torrence and Webster, 1999).However, the lower SST during summer in the tropical Pacific, the weaker the convective activities may be around southern China (Huang and Wu, 1989).Both these conditions might have been responsible for the limited precipitation over the Pearl River basin in 1963 and, consequently, for the severe drought in the basin.
It is observed that soil moisture was in a dry condition around 1963 (see Fig. 7).The general geographic feature of the Pearl River basin is that elevation decreases from the northwest to the southeast.Moreover, the insufficient precipitation condition in 1963 is also noticeable downstream of the Pearl River (Niu, 2013).The precipitation deficit in conjunction with the basin geography features are responsible for the serious dry state in both the northwestern and the southeastern parts, resulting in the severe drought in 1963 over the basin.Drought events normally occur when precipitation is insufficient over a long period (or even a short one) in an already dry state of soil moisture over a long time period.The precipitation deficit in 1989 is comparable to that in 1963, especially for the first and third modes (Niu, 2010).However, the long-term dry condition around 1989 is not comparable to that around 1963 (see Fig. 7).Therefore, the probability of a severe drought event occurrence in 1989 was less than in 1963, as the long-term dry condition around 1963 is more noticeable.
In this study, VIC-model-simulated runoff and soil moisture data over the Pearl River basin, southern China, were analyzed using wavelets to examine the teleconnections between such processes and two large-scale climatic patterns, ENSO and IOD.The results indicate that the hydrological variabilities in the basin have undergone considerable changes in variance and coherence with these two climatic patterns over the study period 1952-2000.
The wavelet power spectrum shows that the distinctive mode of runoff variability is dominant at the low-frequency band (DLVB) of 0.83-5.21yr.The bandpass-filtered time series displays the intervals of high or low variance for different coherent regions in the basin.The variance evolution of runoff is essentially influenced by IOD, especially for the second and third modes (i.e., the eastern and the central parts of the basin), according to the rank correlation analysis for the scale-averaged wavelet power series.In addition, the second mode of runoff in the eastern part of the basin displays more significant coherence with ENSO within the DLVB.
The soil moisture at the headwater basin area (i.e., regions with relatively high altitudes) shows an isolated high variability at longer timescales (2-4 yr), which is well correlated with the low-frequency variability of IOD for both the divided 10 yr periods and the total studied period (49 yr).The soil moisture has the capability to memorize the variability of climatic signals, as had previously been demonstrated (Chen and Kumar, 2002).This capability of soil moisture results in more significant wavelet-based correlations at relatively longer timescales when compared to those for runoff.
The identifiable relations with significant coherence at > 95 % confidence level in a time-frequency domain suggest that the variability in local runoff and soil moisture variability is triggered by large-scale atmospheric circulation and influenced by climatic patterns.These underlying correlations, coupled with the bandpass-filtered time series of hydrological variables, provide valuable information to draw inferences regarding the occurrence of basin-featured extreme hydrological events (i.e., flood or drought).Such an understanding, as evidenced by the current study, is paramount for providing skillful long-term prediction of hydrological processes for large-scale river basins (such as the Pearl River basin) that are of interest at national, and even international, levels.With increasing concerns on the impacts of global climate change on our water resources in the future (see, for instance, Sivakumar, 2011, for a more recent account), this understanding becomes all the more important.

Fig. 2 .
Fig. 2. (a) i. Monthly runoff anomaly time series at Boluo station in the East River and ii.Its local wavelet power spectrum using the Morlet wavelet, normalized by 1/σ 2 (σ 2 = 3,668.94mm 2 ).The dashed curve depicts the cone of influence beyond which the edge effects become important.The shaded contours are at normalized variances of 1, 2 and 4. The white contour closes regions of greater than 95% confidence for a red-noise process with a lag-1 coefficient α of 0.26.iii.Global wavelet power (solid line with circles), mean red-noise spectrum (left dashed line), and 95% confidence spectrum (right dashed line).(b) same as (a) but for soil moisture, in which (σ 2 = 1,342.57mm 2 , α = 0.64).

Fig. 2 .
Fig. 2. (a) (i)Monthly runoff anomaly time series at Boluo station in the East River and (ii) its local wavelet power spectrum using the Morlet wavelet, normalized by 1/σ 2 (σ 2 = 3668.94mm 2 ).The dashed curve depicts the cone of influence beyond which the edge effects become important.The shaded contours are at normalized variances of 1, 2 and 4. The white contour closes regions of greater than 95 % confidence for a red-noise process with a lag-1 coefficient α of 0.26.(iii) Global wavelet power (solid line with circles), mean red-noise spectrum (left dashed line), and 95 % confidence spectrum (right dashed line).(b) Same as (a) but for soil moisture, in which σ 2 = 1342.57mm 2 , and α = 0.64.

Fig. 3 .
Fig. 3. (a) Coherent modes of temporal variability for the wavelet power spectrum of monthly runoff anomaly over the Pearl River basin.The bold line is the first mode, while the line with circles and the dashed line are the second mode and the third mode respectively.DLVB stands for the dominant low-frequency variability band (0.83-5.21 yr).(b) Same with (a) but for soil moisture with DVLB (1.92-5.21yr).

Fig. 3 .
Fig. 3. (a) Coherent modes of temporal variability for the waveletpower spectrum of monthly runoff anomaly over the Pearl River basin.The bold line is the first mode, while the line with circles and the dashed line are the second mode and the third mode, respectively.DLVB stands for the dominant low-frequency variability band (0.83-5.21 yr).(b) Same with (a) but for soil moisture with DLVB (1.92-5.21yr).

Fig. 4 .
Fig. 4. (a) Spatial distribution of coherent modes (in Fig. 3) on global wavelet spectrum of runoff.(b) Same as (a) but for soil moisture.

Fig. 4 .
Fig. 4. (a) Spatial distribution of coherent modes (in Fig. 3) on global wavelet spectrum of runoff.(b) Same as (a) but for soil moisture.

Fig. 5 .
Fig. 5.The flowchart of regional variance evolutions for different hydrological variables.

Fig. 5 .
Fig. 5.The flowchart of regional variance evolutions for different hydrological variables.

Fig. 6 .
Fig. 6.The scale-averaged (0.83-5.21 yr) wavelet power for runoff (thick black curve) and the bandpass filtered time series of square form (runoff 2 ) (thin curve) from 0.83 to 5.21 yr using a wavelet filter.The upper, middle and lower panels are the first, second and third modes of runoff, respectively.The black shading indicates positive peaks in the filtered time series, while the white are negative peaks.The time series of the wavelet power and the bandpass filtered (before squaring) are the averaged time series over the sub-basins belonging to the same mode.The averaged variances of initial monthly runoff time series 2

Fig. 7 .
Fig. 7. Same to Fig. 5, but for soil moisture in 1.92-5.21yr band.The averaged variances of initial monthly soil moisture time series 2

Fig. 8 .
Fig. 8. (a) Squared wavelet coherence (local and global) between the first mode runoff and ENSO time series (left panel), and between the first mode runoff and IOD time series (right panel).The 5% significance level against red noise is shown as a thick contour and the cone of influence (COI) where edge effects might become important is shown as a thin black line.The relative phase relationship is shown as arrows (with in-phase pointing right, anti-phase pointing left, runoff leading ENSO by 90° pointing straight down, and ENSO leading runoff pointing straight up).(b) Same as (a) but for the second model of runoff.(c) Same as (a) but for the third model of runoff.

Fig. 8 .
Fig. 8. (a) Squared wavelet coherence (local and global) of the first mode runoff and ENSO time series (left panel), and between the first mode runoff and IOD time series (right panel).The 5 % significance level against red noise is shown as a thick contour and the cone of influence (COI) where edge effects might become important is shown as a thin black line.The relative phase relationship is shown as arrows (with in-phase pointing right, anti-phase pointing left, runoff leading ENSO by 90 • pointing straight down, and ENSO leading runoff pointing straight up).(b) Same as (a) but for the second model of runoff.(c) Same as (a) but for the third model of runoff.
are helpful for understanding the coherence results obtained in the present study.The complex correlations between monsoon and climatic indices indicate a phase difference of the variability of runoff or soil moisture with climatic signals.

Fig. 10 .
Fig. 10.Comparison of the scale-averaged wavelet power for runoff and selected climate indices for the dominant low-frequency variability band 0.83-5.21yr.

Fig. 10 .Fig. 11 .
Fig. 10.Comparison of the scale-averaged wavelet power for runoff and selected climate indices for the dominant low-frequency variability band 0.83-5.21yr.

Fig. 11 .
Fig. 11.Comparison of the scale-averaged wavelet power for soil moisture and selected climate indices for the dominant low-frequency variability band 1.92-5.21yr.

Table 1 .
Runoff simulation results at a monthly time step for 10 gauging stations over the Pearl River basin (in Fig.1).

Table 2 .
Statistical properties of the observed precipitation and the simulated monthly time seriesof 10 sub-basins over the Pearl River basin.

Table 3 .
Spearman's rank correlations between the wavelet powers of runoff/soil moisture and selected climate indices in the low-frequency dominant variability band.DOF denotes the degrees of freedom.Bold indicates greater than 95 % confidence level.