Summary: We explore the evolution of the probability density function (PDF) for an initially deterministic passive scalar diffusing in the presence of a uni-directional, white-noise Gaussian velocity field. For a spatially Gaussian initial profile we derive an exact spatio-temporal PDF for the scalar field renormalized by its spatial maximum. We use this problem as a test-bed for validating a numerical reconstruction procedure for the PDF via an inverse Laplace transform and orthogonal polynomial expansion. With the full PDF for a single Gaussian initial profile available, the orthogonal polynomial reconstruction procedure is carefully benchmarked, with special attentions to the singularities and the convergence criteria developed from the asymptotic study of the expansion coefficients, to motivate the use of different expansion schemes. Lastly, Monte-Carlo simulations stringently tested by the exact formulas for PDFs and moments offer complete pictures of the spatio-temporal evolution of the scalar PDFs for different initial data. Through these analyses, we identify how the random advection smooths the scalar PDF from an initial Dirac mass, to a measure with algebraic singularities at the extrema. Furthermore, the Péclet number is shown to be decisive in establishing the transition in the singularity structure of the PDF, from only one algebraic singularity at unit scalar values (small Péclet), to two algebraic singularities at both unit and zero scalar values (large Péclet).