This study proposes a surrogate-assisted evolutionary optimization framework for small dataset that integrates machine learning–based surrogate models with evolutionary algorithms for the aerodynamic optimization of a spiked blunt body in supersonic flow. A database of simulated cases covering a range of Mach numbers, spike length ratios (L/D), and diameter ratios (d/D) was used to train regression models and identify optimal geometries. Among the tested algorithms, the Gradient Boosting Regressor (GBR) achieved the best predictive performance (R² = 0.8909, RMSE = 0.00775), accurately capturing the nonlinear dependencies of the drag coefficient (Cd). Evolutionary optimization methods, including Differential Evolution (DE), Covariance Matrix Adaptation Evolution Strategy (CMA-ES), and Genetic Algorithm (GA), consistently converged to near-optimal configurations, with DE exhibiting the most stable behavior across Mach regimes. SHapley Additive exPlanations (SHAP) analysis revealed that (L/D) is the most influential parameter on Cd, followed by Mach number and (d/D), highlighting the dominant effect of geometric slenderness in drag reduction. The integration of data-driven modeling with evolutionary computation provides a robust framework for aerodynamic optimization, offering both predictive accuracy and physical interpretability. These results demonstrate the potential of hybrid Machine Learning-Evolutionary algorithms and CFD approaches to accelerate the design of high-performance configurations in supersonic applications.