Abstract: Abstract: More effective methods are needed to evaluate the water demand to improve water resource management as the population and economy grow leading to water shortage. Time series models are useful tools in the estimation and forecasting of reference evapotranspiration series and their changes. However, due to the dynamic nature of reference evapotranspiration, accurate estimation of variance has been being a challenging task and requires new modelling approaches in application. The Generalized Autoregressive Conditional Heteroscedasticity (GARCH) family models, which provide an appropriate framework for focusing on the conditional variance remaining in the residuals of the time series models, were applied in this study to comply with the estimation. Daily climate data including average air pressure, air temperature, relative humidity, wind speed, and saturation deficit were provided by the Chinese National Climatic Data Center (CDC). The observed meteorology dataset from 1953-2007 in Yichang station in Hubei province in China was chosen as the selected station for its data available. The FAO Penman-Monteith method was applied to obtain the potential reference evaporation series. Then, for deseasonalization of the seasonal variation in the series, a Seasonal Autoregressive Moving Average (SARMA) model was set up to estimate the conditional mean of the potential reference evaporation. Besides, four types of GARCH models (GARCH, Threshold GARCH, Exponential GARCH, and Power GARCH models) were investigated to take the advantages of GARCH family models to simulate the conditional variance of the potential reference evaporation series. Also, six criteria were utilized to compare the model performance. It was shown that the data of potential reference evaporation series in Yichang Station had skewness, and right tail on an annual cycle, indicating the necessary of SARMA model used in this study. The final SARMA model was chosen by the minimum of AIC values. The results of SARMA model showed that it was efficient for modelling the monthly mean total daily reference evaporation series as a small redisudal mean square error (RMSE) between the observed and estimated values (0.089 mm). However, the heteroscedasticity was present in the residuals of SARMA model according to the Engle test, which suggested the necessary of GARCH models used for modeling of reference evaporation series. Results of GARCH models showed its ability to remove the heteroscedasticity from the reference evaporation residuals. The asymmetric effects of solar radiation series were also confirmed by the application of GARCH family models for estimating the residuals of SARMA model. Among the SARMA and GARCH models, the EGARCH model was best to predict the series because the prediction of the EGARCH model had a narrower confidence level. The better interval estimation would provide more useful information for analysis of hydrological processes and was undoubtedly more favorable to further uncertainty analysis in water resource management. The multi-criteria evaluation for model performance also proved that the EGARCH model was best among SARMA and GARCH models, and it was recommended to use the EGARCH model to estimate the potential reference evaporation series in practice. In conclusion, the results indicated that the application of GARCH family models would be a promising alternative over the traditional approaches in the estimation of potential reference evaporation series, and can be a useful tool to existing water resources management. Further studies about the comprehensive application of GARCH models should focus on the model performances with different observed dataset among stations. It also suggests to use more GARCH models to estimate the hydrological series and the results of accuracy estimation should be linked with water risk analysis.