Refractory multi-principal element alloys (RMPEAs) draw great interest with their superior mechanical properties and extremely high melting points, yet the strengthening mechanism remains unclear. Here, we calculate the critical resolved shear stress (CRSS) for a single dislocation to move in RMPEAs consisting of 4 or 5 elements with or without short-range order (SRO) to represent strength by a machine learning-based interatomic potential. The increased CRSS is then attributed to high lattice distortion, elastic mismatch, and SRO strengthening, all of which originate from the solid solution strengthening theory. After detailed research of the CRSS across many RMPEAs systems with different composition ratios, we construct an XGBoost model to predict the CRSS from a few parameters and rank their importance. We find that lattice distortion strongly influences both dislocation types and reduces the screw-to-edge ratio in CRSS, while the elastic mismatch has a more significant impact on the screw dislocation than the edge one.