A Gaussian mixture model approach is proposed for accurate uncertainty propagation through a general nonlinear system. The transition probability density function, is approximated by a finite sum of Gaussian density functions whose parameters (mean and covariance) are propagated using linear propagation theory. Two different approaches are introduced to update the weights of different components of a Gaussian mixture model for uncertainty propagation through nonlinear system. The first method updates the weights such that they minimize the integral square difference between the true forecast probability density function and its Gaussian sum approximation. The second method uses Fokker-Planck equation error as a feedback to adapt for the amplitude of different Gaussian components while solving a quadratic programming problem. The proposed methods are applied to a variety of problems in the open literature, and argued to be an excellent candidate for higher dimensional uncertainty propagation problems.