For spatial prediction of soil variables where the local mean can be expressed as a local trend, or a linear function of some auxiliary (external drift) variable, the state-of-the-art is to estimate the spatial covariance parameters of the residual variation by residual maximum likelihood (REML). Method of moments estimators, as used to compute point estimates of the variogram in most geostatistical studies, are biased when applied to residuals from separately fitted trend or external drift models. Both REML and method of moments estimates are susceptible to the effects of extreme values such as may appear in soil data due to a second process (e.g. pollution) superimposed on continuous background variation. While robust point estimators of the variogram exist, the problem of robustly estimating spatial covariance parameters by REML has not been tackled. A robust REML estimator is described that is less sensitive to outliers. In tests upon contaminated, second order stationary, data sets the robust REML estimator estimates the underlying variogram more accurately than the standard REML estimator, the method of moments and robust point estimators of the variogram. However the primary advantage of the robust REML algorithm is that it is applicable in the presence of external drift. The performance of the REML estimator depends upon a damping parameter, the optimal value of which varies according to the contaminating distribution. A diagnostic statistic for selecting this parameter is suggested. The robust REML estimator performs well with damping parameters suggested by this method although these values tend to be larger than the optimal value. The robust REML estimator is applied to a survey of the lead content of soil in the Swiss Jura.