This paper proposes a new method to propagate uncertainties undergoing nonlinear dynamics using the Koopman Operator (KO). Probability density functions are propagated directly using the Koopman approximation of the solution flow of the system, where the dynamics have been projected on a well-defined set of basis functions. The prediction technique is derived following both the analytical (Galerkin) and numerical (EDMD) derivation of the KO, and a least square reduction algorithm assures the recursivity of the proposed methodology.