Abstract: Soil salinization has threatened the ecological stability and sustainable agriculture under the global environmental and climate change. Among them, the Yellow River Delta region can exhibit the seriously soil salinization and ecological heterogeneity with the complex features, compared with the arid areas. However, conventional two-dimensional feature spaces have limited to only two environmental variables under such complex environments. Furthermore, only a few specific indices cannot systematically construct and evaluate a index pool to identify the optimal parameters of the soil salinization. It is often required for the full extraction of the salinization information. In contrast, three-dimensional feature spaces can be expected to effectively utilize the multi-band remote sensing data and multiple environmental variables, indicating the high-accuracy monitoring. In this study, the three-dimensional feature spaces were constructed with the spectral indices using feature selection and Bayesian-optimized XGBoost models, in order to estimate the tidal soil salinization in the Yellow River Delta. Landsat 9 satellite imagery was employed to build a spectral index pool. A total of 35 spectral indices were extracted, including the vegetation, salinity and water indices. A XGBoost model with Bayesian optimization was utilized to evaluate the features for the modelling efficiency and parameter screening, according to the built-in Gain metrics. The top two most important indices were retained from each category. Multiple three-dimensional feature space models were constructed to combine the representative indices. The three coordinate axes were represented the different index types in the three-dimensional spaces. Any point (x, y, z) in the feature space was corresponded to the values of the three indices for a specific pixel in the remote sensing image. Simultaneously, the multiple two-dimensional feature spaces were built using the single most important index from each category. The accuracy evaluation metrics were compared with the field-measured data. The optimal three- and two-dimensional feature space models were determined for the soil salinization inversion in the Yellow River Delta. Regional salinization spatial analysis was then conducted after optimization. The results show that: 1) The XGBoost model with Bayesian optimization was effectively screened the most relevant indices. Salinity indices were achieved the highest modelling accuracy (R² = 0.921, RMSE = 0.964 g/kg, and RPIQ = 8.422), with the salinity index 7 (SI7), indicating the highest feature importance (0.341). Ultimately, eight of the most informative feature indices were selected after comparison. 2) Three-dimensional feature spaces were more fully exploited the spectral information, compared with the two-dimensional feature spaces. The optimal three-dimensional model was improved 0.059 in R² and 1.191 in RPIQ, with a reduction of 0.069 g/kg in RMSE, indicating that the three-dimensional approach was improved the high-precision prediction of the soil salinization. 3) Among the three-dimensional feature space models, the SI8-Albedo-WI was achieved in the highest accuracy (R² = 0.922, RMSE = 0.863 g/kg, RPIQ = 7.645, and Kappa coefficient = 86%), whereas the ERVI-WI-Albedo model performed the worst (R² = 0.519, RMSE = 3.464 g/kg, and RPIQ = 1.087). 4) The moderately salinized areas were accounted for the largest proportion (29.7%) in the Yellow River Delta region, primarily distributed in the central-western part of the Kenli District and Lijin County; The severely salinized areas were constituted the smallest proportion (9.8%), which were located mainly in the eastern part of Kenli District. The findings can also provide crucial references and decision-making support to prevent and remediate the soil salinization in the Yellow River Delta.