The integration of distributed energy resources (DERs) and unpredictable loads has increased uncertainty in power systems. Traditional methods struggle to assess performance under these uncertainties, and existing probabilistic methods face challenges with complexity and accuracy. This paper introduces a new combined analytical–numerical probabilistic method to assess the impact of DERs on voltage stability. Using Bayesian Parameter Estimation (BPE), the method derives the analytical properties of random variables (RVs) associated with DERs and loads, obtaining posterior distributions. The Metropolis–Hastings sampling technique then estimates these posteriors numerically, enabling accurate predictions of DERs and load outputs. Voltage stability analysis was performed using the continuation power flow method and validated on the IEEE 59-bus test system in MATLAB/Simulink. The results show that integrating DERs significantly improves voltage stability. The proposed method outperforms the Monte Carlo simulation (MCS)-based method in accuracy and computational speed, increasing DERs penetration and voltage stability limits by 3%. It closely matches MCS voltage estimates but requires fewer iterations (500 per loading increment) compared to MCS’s 1000, leading to faster computation times (a few hours to one day versus up to three days for MCS). This method provides an efficient solution for managing uncertainties in power systems.