We present an upscaled Lagrangian approach to predict the plume evolution in highly heterogeneous aquifers. The model is parameterized by transport independent characteristics such as the statistics of hydraulic conductivity and the Eulerian flow speed. It can be conditioned on the tracer properties, the conductivity data at the injection region, and is able to account for mobile-immobile mass transfer. Thus, the model is transferable to different solutes and hydraulic conditions. It captures the large scale non-Gaussian features for the evolution of the longitudinal mass distribution observed for the bromide and tritium tracer plumes at the macrodispersion experiment (MADE) site (Columbus, Ohio, USA), which are characterized by a slow moving peak and pronounced forward tailing. These large scale features are explained by advective tracer propagation due to a broad distribution of spatially persistent Eulerian flow speeds as a result of spatial variability in hydraulic conductivity.