A surrogate model based groundwater optimization model was developed to solve the non-aqueous phase liquids (NAPLs) contaminated groundwater remediation optimization problem. To illustrate the impact of sampling method improvement to the surrogate model performance improvement, aiming at a nitrobenzene contaminated groundwater remediation problem, optimal Latin hypercube sampling (OLHS) method was introduced to sample data in the input variables feasible region, and a radial basis function artificial neural network was used to construct a surrogate model. Considering the surrogate model's uncertainty, a chance-constrained programming (CCP) model was constructed, and it was solved by genetic algorithm. The results showed the following, for the problem considered in this study. (1) Compared with the Latin hypercube sampling (LHS) method, the OLHS method improves the space-filling degree of sample points considerably. (2) The effects of the two sampling methods on surrogate model performance were analyzed through comparison of goodness of fit, residual and uncertainty. The results indicated that the OLHS-based surrogate model performed better than the LHS-based surrogate model. (3) The optimal remediation strategies at 99%, 95%, 90%, 85%, 80% and 50% confidence levels were obtained, which showed that the remediation cost increased with the confidence level. This work would be helpful for increasing surrogate model performance and lowering the risk of a groundwater remediation strategy.