Based on a stochastic approach, drainage in a layered soil profile with macropores was predicted using a one-dimensional numerical flow model. Soil hydraulic functions, theta(psi) and K(psi), were obtained from laboratory measurements on undisturbed 100 cm(3) soil cores. The cores were collected in 60 different locations for each of three soil horizons along a 31 m transect on a sandy loam soil. The effect of a stochastic description of theta(psi) and K(psi) on predicted fluxes and water content profiles was evaluated with the aid of Monte-Carlo simulations. Initially simulations were carried out with the actual field data. Alternatively, 500 sets of the multi-normal parameter vector P = {theta(s), alpha, n, K-s} were generated and used in a Monte-Carlo simulation approach. In addition, deterministic calculations were carried out using (1) arithmetic mean values of untransformed hydraulic parameters, (2) arithmetic mean values obtained after backtransformation of Gaussian distributed parameters, and (3) alternatively using a set of mean bimodal retention and conductivity curves. Predictions obtained with the two stochastic methods and the deterministic cases are compared with measured cumulative outflow and moisture content profiles from fifteen 1 m long undisturbed soil columns collected along the same transect. All simulations based on unimodal hydraulic functions underestimated the mean observed drainage. If actual field data are used, the mean outflow is 74% of the mean measured cumulative outflow, compared to 70% for the Monte-Carlo simulation. The range mean +/- two standard deviations of the simulated outflow for both methods is considerably smaller compared to the measured range of variation. Deterministic predictions of the total drainage amounted to 83, 83, and 127% of the observed drainage for the cases 1, 2, and 3, respectively. Observed outflow curves displaying macropore flow could not be predicted with the unimodal Richards' flow equation. The results indicate that the field-scale unsaturated flow behaviour in soils with macropores cannot be accurately predicted assuming a unimodal retention characteristic. This is true for both the mean drainage but even more so for the extremes. The use of bimodal hydraulic functions considerably increased the total drainage, resulting in an overprediction of the mean measured drainage. Predicted soil water contents in the 1 m deep soil profile are mostly higher than the measurements. Contrary to the drainage, the range of variation of the predicted water contents is always larger than the measured variability, except for the surface values.