In this study, a first-order, nonstationary stochastic model of transient flow is developed which is applicable to the entire domain of a bounded vadose zone in the presence of sink/source. We derive general equations governing the statistical moments of the flow quantities by perturbation expansions. Owing to the mathematical complexity of the equations, in general we need to solve them numerically. The numerical moment equation approach, however, has the flexibility in handling different boundary conditions, flow configurations, input covariance structures, and soil constitutive relationships. The moment equations presented in this study are simpler and easier to solve than those in the literature for transient unsaturated flow. We solve these moment equations by the method of finite differences and demonstrate the developed model through some one- and two-dimensional examples under various transient conditions. ¿ 1999 American Geophysical Union