A time-averaged 30-m-depth shear wave velocity ( V S30) map is developed for New Zealand as a weighted combination of a geology-based and a terrain-based model. A Bayesian updating process allows local V S30 measurements to control model estimates where data exist and uses model estimates developed for other parts of the world where local data are sparse or nonexistent. Geostatistical interpolation is performed on the geology-based and terrain-based models using local V S30 measurements to constrain the model in the vicinity of data. Conventional regression kriging is compared with a flexible multivariate normal (MVN) approach that allows for arbitrary assumptions regarding measurement uncertainty at each data location. A modification to the covariance structure in the MVN application allows for more realistic estimates by reducing undesirable extrapolation across geologic boundaries. The results of kriging and MVN approaches are compared. The geology-based and terrain-based MVN models are combined to produce a final model suitable for engineering applications. The 100-m resolution map outputs are publicly available.