The probability distribution functions (PDFs) of momentum flux and zonal flow formation in ion-temperature-gradient (ITG) turbulence are investigated in two different models. The first is a general five-field model (ni, ϕ, Ti, Te, vi∥) where a reductive perturbation method is used to derive dynamical equations for drift waves and a zonal flow. The second is a reduced two-field model (ϕ, Ti) that has an exact non-linear solution (bipolar vortex soliton). In both models the exponential tails of the zonal flow PDFs are found with the same scaling (), but with different coefficients cZF. The PDFs of momentum flux is, however, found to be qualitatively different with the scaling (PDF ∼ exp{−cMRs}), where s = 2 and s = 3/2 in the five and two-field models, respectively.