Soil particle size fractions (PSFs) are vital for understanding soil functions and are widely used in environmental and land surface modeling. Accurate spatial mapping of PSFs is crucial and relies on robust prediction methods. This study, conducted in the typical loess plateau, integrates multidimensional auxiliary data-including remote sensing bands, topography, vegetation indices, socioeconomic indicators, and soil properties-to enhance prediction accuracy. We employ a two-point machine learning (TPML) model that incorporates spatial autocorrelation and attribute similarity into a unified framework for predicting PSF distribution. TPML enhances local predictions and overcomes the curse of dimensionality, ensuring robust performance with limited training samples. The TPML model is evaluated against random forest (RF), random forest regression kriging (RFRK), inverse distance weighting, and ordinary kriging under various sample sizes. Accuracy is assessed using mean absolute error (MAE), root mean square error (RMSE), and coefficient of determination (R-2). TPML consistently outperforms other models across sample sizes (50-200), with optimal results at 150 samples (e.g., silt: MAE 4.21, RMSE 7.61, R-2 0.62). Its performance improves with larger training data, particularly for silt and sand, where R-2 increased from 0.59 to 0.69 and 0.63 to 0.73, respectively. Key environmental variables-latitude, precipitation, valley depth, vegetation type, and SAVI-strongly influenced spatial PSF patterns, showing a gradient from coarse textures in the north to fine textures in the south. The results demonstrate the effectiveness of TPML in enhancing PSF mapping by integrating spatial and attribute-based similarities within a high-dimensional prediction framework.