The substantial spatial and temporal variability of pesticides has led to large uncertainties when determining their peak aqueous concentrations. There is however a lack of large-scale studies dealing with accurate determination of annual maximum daily concentration (AMDC) across the landscape and over time based on the publicly available monitoring data. We developed a novel data-driven approach that firstly used time series modeling to generate AMDCs for qualified water monitoring sites in the conterminous U.S. With feature variables such as pesticide use and land cover compiled into the dataset, machine learning models using eXtreme Gradient Boosting (XGBoost) and Random Forest Regressor (RF) were then developed to estimate AMDCs in surface waters across the U.S. Both models exhibited significant predictability, while a hybrid model consisting of the average predictions by XGBoost and RF model had the highest prediction accuracy (mean absolute error (MAE): 1.23; R2: 0.61). The analysis of permutation variable importance indicated that pesticide use and drainage area were the two most important drivers. Partial dependence analysis revealed that pesticide use, precipitation, cultivated crop land cover and solubility exhibited concentration-promoting effects, whereas drainage area and molecular weight had concentration-demoting effects. Soil adsorption coefficient (Koc) showed nonmonotonic effects. The hybrid model was used to predict and map AMDCs of four example pesticides, including 2,4-dichlorophenoxyacetic acid (2,4-D), atrazine, glyphosate and imidacloprid during 2016–2019 at national scale. The predictive capability was validated using independent monitoring datasets. The fully evaluated approach significantly reduced the uncertainties in modeling annual peak concentrations and served as a valuable solution for conducting geographically oriented, highly refined exposure assessments for pesticides.