In this paper we use a double-porosity/double-permeability approach to study transient, saturated flow in heterogeneous, fractured porous media. The spatial variabilities in both the fracture and matrix continua motivate us to treat some of the fracture and matrix medium properties as stochastic processes. Hence flow through heterogeneous, fractured porous media is also amenable to stochastic analysis. In particular, we develop a moment equation-based stochastic theory for flow in fractured porous media. We derive general equations governing the statistical moments of the pressure heads in such media by perturbation expansions. The solutions of these moment equations are the first two moments of the pressure heads. These may be used to construct the confidence intervals of the fracture and matrix pressure heads, which are measures of uncertainties caused by incomplete knowledge of the fracture and matrix medium heterogeneities. The resulting moment equations of pressure heads are coupled and need to be solved sequentially in an iterative manner. The two-dimensional version of the equations is implemented numerically and illustrated with some examples. ¿ 2000 American Geophysical Union