Thermostability enhancement is a fundamental aspect of protein engineering as a biocatalyst's half-life is key for its industrial and biotechnological application, particularly at high temperatures and under harsh conditions. Thermostability changes upon mutation originate from modifications of the free energy of unfolding (ΔG u ), making thermostabilization extremely challenging to predict with computational methods. In this contribution, we combine global conformational sampling with energy prediction using AlphaFold and Rosetta to develop a new computational protocol for the quantitative prediction of thermostability changes upon laboratory evolution of acyltransferase LovD and lipase LipA. We highlight how using an ensemble of protein conformations rather than a single three-dimensional model is mandatory for accurate thermostability predictions. By comparing our approaches with existing ones, we show that ensembles based on AlphaFold models provide more accurate and robust calculated thermostability trends than ensembles based solely on crystallographic structures as the latter introduce a strong distortion (scaffold bias) in computed thermostabilities. Eliminating this bias is critical for computer-guided enzyme design and evaluating the effect of multiple mutations on protein stability.