The spatial distribution of geochemical elements is complex, and traditional regional geochemical surveys and analyses consider only spatial information elements, which do not reflect the heterogeneity of other elements. Therefore, this study used hierarchical regionalization labels, such as geological, geographical, spatial, and ecological elements along with anthropogenic activities, as covariates of the machine learning model to predict the spatial distribution of selenium elements. By comparing Random Forest (RF), extreme gradient lifting tree (XGboost), Feature Contribution-Random Forest (FC-RF), Deep learning Neural Network (DNN) models, and traditional spatial interpolation techniques; this study showed that the FC-RF model has higher precision in selenium prediction, which indicates that multilevel interpolation technology based on the FC-RF model is more feasible for predicting the spatial distribution of soil selenium. Based on the predicted results, a Bayesian Laser Interferometer Space Antenna (LISA) spatial correlation analysis was conducted on the spatial distribution of selenium on Naozhou Island, indicating that the key factors affecting the spatial distribution of selenium are basic geological conditions and anthropogenic activities. This study further expands the service application dimensions of geochemical data, thus making the use of soil geochemical data more accurate and scientific.