The location of the interface of geological formations is an important piece of information for tunneling construction. As site investigation data are usually limited, the uncertainties in locating geological interfaces for the sections between boreholes can be large and challenging to estimate. A suitable geostatistical method is thus needed for spatial prediction of the geological interfaces. In this paper, the performance of three commonly used spatial prediction methods, namely the multivariate adaptive spline regression (MARS), conditional random field (CRF) method, and thin-plate spline interpolation (TPSI) methods, are evaluated for two-dimensional cases using the boreholes data from three sites in Singapore. The prediction accuracies, patterns of the predicted surfaces, and prediction uncertainties obtained from the three methods are compared. A zonation is also proposed to improve the prediction accuracy of the MARS method. The results indicate that the MARS method can show the spatial trend of the geological interface more clearly than the other two methods. The TPSI method produces undesirable oscillations of the surface of geological interfaces and the CRF method may underestimate the extreme values of the geological interface elevations. In general, the prediction accuracy of the MARS method is similar to that of the CRF method, but higher than that of the TPSI method. For cases with very limited data in geologically complex areas, the MARS may have larger errors than the CRF method. However, the accuracy of the former can be significantly improved if a reasonable zonation is performed.