Plant species diversity (PSD) has always been an essential component of biodiversity and plays an important role in ecosystem functions and services. However, it is still a huge challenge to simulate the spatial distribution of PSD due to the difficulties of data acquisition and unsatisfactory performance of predicting algorithms over large areas. A surge in the number of remote sensing imagery, along with the great success of machine learning, opens new opportunities for the mapping of PSD. Therefore, different machine learning algorithms combined with high-accuracy surface modeling (HASM) were firstly proposed to predict the PSD in the Xinghai, northeastern Qinghai-Tibetan Plateau, China. Spectral reflectance and vegetation indices, generated from Landsat 8 images, and environmental variables were taken as the potential explanatory factors of machine learning models including least absolute shrinkage and selection operator (Lasso), ridge regression (Ridge), eXtreme Gradient Boosting (XGBoost), and Random Forest (RF). The prediction generated from these machine learning methods and in situ observation data were integrated by using HASM for the high-accuracy mapping of PSD including three species diversity indices. The results showed that PSD was closely associated with vegetation indices, followed by spectral reflectance and environmental factors. XGBoost combined with HASM (HASM-XGBoost) showed the best performance with the lowest MAE and RMSE. Our results suggested that the fusion of heterogeneous data and the ensemble of heterogeneous models may revolutionize our ability to predict the PSD over large areas, especially in some places limited by sparse field samples.