Hydrological extremes in the water cycle can significantly affect surface water engineering design, and represents the high-impact response of surface water and groundwater systems to climate change. Statistical analysis of these extreme events provides a convenient way to interpret the nature of, and interaction between, components of the water cycle. This study applies three probability density functions (PDFs), Gumbel, stable, and stretched Gaussian distributions, to capture the distribution of extremes and the full-time series of storm properties (storm duration, intensity, total precipitation, and inter-storm period), stream discharge, lake stage, and groundwater head values observed in the Lake Tuscaloosa watershed, Alabama, USA. To quantify the potentially non-stationary statistics of hydrological extremes, the time-scale local Hurst exponent (TSLHE) was also calculated for the time series data recording both the surface and subsurface hydrological processes. First, results showed that storm duration was most closely related to groundwater recharge compared to the other storm properties, while intensity also had a close relationship with recharge. These relationships were likely due to the effects of oversaturation and overland flow in extreme total precipitation storms. Second, the surface water and groundwater series were persistent according to the TSLHE values, because they were relatively slow evolving systems, while storm properties were anti-persistent since they were rapidly evolving in time. Third, the stretched Gaussian distribution was the most effective PDF to capture the distribution of surface and subsurface hydrological extremes, since this distribution can capture the broad transition from a Gaussian distribution to a power-law one.process-based models requiring intensive data that are typically not available for many study sites, statistical analysis of hydrological extremes is an attractive tool to interpret these extreme events within and across systems from simple measurements [7-10].There are two major challenges when applying basic statistics to analyze hydrologic extremes. First, hydrologic processes in the water cycle are interconnected, while basic statistical analysis of extremes often does not consider multiple systems and tends to oversimplify the complexity of the correlated processes like precipitation and groundwater table fluctuations [11,12]. Understanding how the subtle properties of one system's extreme events can affect the other interconnected systems requires more in-depth analysis and use of advanced statistical techniques. Second, these basic statistical techniques often rely on assumptions or major simplifications of properties for water systems, such as stationarity in both space and time, which may not always be valid for real world dynamics [13][14][15]. These issues motivated this study.One example of the assumption/simplification used by basic statistical studies in hydrological extremes is the well-known Gumbel distribution. Statistics of extremes is one of th...