Steel circular hollow section (CHS) members are one of the commonest axial force-resisting structural members of steel structures. Due to the working conditions of the steel CHS members, axial buckling strength is the critical parameter to determine steel CHS members’ performance. The traditional prediction methods of axial buckling strength can’t gain a satisfactory result due to the variabilities in geometry, material, and initial imperfections, which have a significant impact. To address this research gap, a prediction model is created to predict the nominal axial buckling strength of steel CHS members. This model is developed after the comparison among ten different machine learning algorithms. The database is grounded in extensive numerical simulations where the accuracy of the numerical simulation method is confirmed by experimental results. Additionally, the prediction model is interpreted by Shapley Additive exPlanations (SHAP) method to investigate the relationship between the input parameters and the probabilistic axial buckling strength. Then, the probabilistic axial buckling strength prediction model is further established based on the developed prediction models, where the Latin hypercube sampling method is applied to address the variability of geometries, material, and initial imperfections. The distribution of the prediction model’s results exhibits a high degree of concordance with the numerical experimental results. With the developed prediction model, the probabilistic axial buckling strength of CHS members can be calculated quickly and efficiently.