It is increasingly recognized that the assumption of stationarity is not always appropriate for estimating return periods from meteorological time series under the effects of climate change and climate variability. Here, we assessed the capability of a non-stationary framework for modeling seasonal precipitation series using Generalized Additive Models for Location, Scale and Shape (GAMLSS), a well-established approach for modeling hydro-meteorological variables under non-stationary conditions. Seasonal precipitation series were considered from 95 stations covering the period of 1970–2021 across Turkey. Four widely used homogeneity tests showed that the seasonal precipitation series were generally reliable, with no obvious anthropogenic influences or errors. The parameters of the fitted distributions were modeled as a function of large-scale oscillation indices known to control precipitation across Turkey, namely the North Atlantic Oscillation (NAOI), North Sea-Caspian Pattern (NCPI), Mediterranean Oscillation (MOI), and Southern Oscillation (SOI). The model with oscillation indices performed better at 85%, 79%, 76%, and 54% of sites for autumn, winter, spring, and summer, respectively, than the model with no covariates. In particular, the NCPI was seen as a significant predictor during the winter, while the MOI captured precipitation variability well across the country during autumn and summer. The NAOI appeared as another important predictor during all seasons except summer, while the SOI appeared as a significant explanatory variable in certain regions. Using the non-stationary models, we then computed seasonal precipitation estimates for different return periods (i.e., 20, 50, and 100 years) and considered the minimum and maximum possible extreme precipitation scenarios at each site. We show how the use of simple minimum/maximum values derived from the non-stationary models can help provide water resource managers and policy makers with a plausible range of extreme values, rather than the single deterministic value obtained from the traditional stationary approach.