Abstract. Parameter Sensitivity analysis plays a critical role in efficiently determining main parameters, enhancing the effectiveness of estimation of parameters, and uncertainty quantification in hydrologic modeling. In this paper, we demonstrate uncertainty and sensitivity analysis technique for the holistic SWAT+ model, coupled with new gwflow module, spatially distributed, physically based groundwater flow modeling. Main calculated groundwater inflows and outflows include boundary exchange, pumping, saturation excess flow, groundwater–surface water exchange, recharge, groundwater–lake exchange, and tile drainage outflow. We present the method for four watersheds located in different areas of the United States for 16 years (2000–2015), emphasizing regions of extensive tile drainage (Winnebago River, Minnesota, Iowa), intensive surface–groundwater interaction (Nanticoke River, Delaware, Maryland), groundwater pumping for irrigation (Cache River, Missouri, Arkansas), and mountain snowmelt (Arkansas Headwaters, Colorado). The main parameters of coupled SWAT+gwflow model are estimated utilizing the parameter estimation software (PEST). The monthly streamflow of holistic SWAT+gwflow is evaluated based Nash–Sutcliffe efficiency index (NSE), percentage bias (PBIAS), determination coefficient (R2), and Kling–Gupta efficiency coefficient (KGE), whereas groundwater head is evaluated using mean absolute error (MAE). The Morris method is employed to identify the key parameters influencing hydrological fluxes. Furthermore, the iterative ensemble smoother (iES) is utilized as a technique for Uncertainty Quantification (UQ) and Parameter Estimation (PE) and to decrease the computational cost owing to the large number of parameters. Depending on the watershed, key identified selected parameters include aquifer specific yield, aquifer hydraulic conductivity, recharge delay, streambed thickness, streambed hydraulic conductivity, area of groundwater inflow to tile, depth of tiles below ground surface, hydraulic conductivity of the drain perimeter, river depth (for groundwater flow processes); runoff curve number (for surface runoff processes); plant uptake compensation factor, soil evaporation compensation factor (for Potential and actual evapotranspiration processes); soil available water capacity, percolation coefficient (for Soil water processes). The presence of gwflow parameters permits for the recognition of all key parameters in the surface/subsurface flow processes, with results substantially differing if the base SWAT+ models are utilized.