Groundwater storage in aquifers has become a vital water source due to the water scarcity inthe last years. However, aquifer systems are full of uncertainties, which inevitably propagate throughout the calculations, mainly reducing the reliability of the model output. This study develops a novel two-dimensional stochastic confined groundwater flow model. The proposed model is developed by linking the stochastic governing partial di↵erential equations through their one-to-one correspondence to the nonlocal Lagrangian-Eulerian form of the Fokker- Planck equation (LEFPE). In the form of the LEFPE, the resulting deterministic governing equation describes the spatio-temporal evolution of the probability density function of the state variables in the confined groundwater flow process. This probability evolution is performed by one single numerical realization rather than requiring thousands of simulations in the Monte Carlo simulation. Consequently, the ensemble groundwater flow process’s mean and standard deviation behavior can be modeled under uncertainty in the transmissivity field and recharge and/or pumping conditions. In addition, an appropriate numerical method for its solution is subsequently devised. Then, LEFPE’s solution is presented, discussed, and illustrated through numerical examples. Furthermore, they are compared against the results obtained through the Monte Carlo simulations to evaluate the performance of the LEFPE framework. Results suggest that the proposed model appropriately characterizes the ensemble behavior in confined groundwater systems under uncertainty in the transmissivity field.