An iterative geostatistical inverse approach is developed to estimate conditional effective unsaturated hydraulic conductivity parameters, soil-water pressure head, and degree of saturation in heterogeneous vadose zones. This approach is similar to the classical cokriging technique, and it uses a linear estimator that depends on covariances and cross covariances of unsaturated hydraulic parameters, soil-water pressure head, and degree of saturation. The linear estimator is, however, improved successively by solving the governing flow equation and by updating the residual covariance and cross-covariance functions in an iterative manner. As a result, the nonlinear relationship between unsaturated hydraulic conductivity parameters and head is incorporated in the estimation and the estimated fields are approximate conditional means. The ability of the iterative approach is demonstrated through some numerical examples.¿ 1997 American Geophysical Union