Accurate prediction of thermodynamic properties requires an extremely accurate representation of the free-energy surface. Requirements are twofold—first, the inclusion of the relevant finite-temperature mechanisms, and second, a dense volume–temperature grid on which the calculations are performed. A systematic workflow for such calculations requires computational efficiency and reliability, and has not been available within an ab initio framework so far. Here, we elucidate such a framework involving direct upsampling, thermodynamic integration and machine-learning potentials, allowing us to incorporate, in particular, the full effect of anharmonic vibrations. The improved methodology has a five-times speed-up compared to state-of-the-art methods. We calculate equilibrium thermodynamic properties up to the melting point for bcc Nb, magnetic fcc Ni, fcc Al, and hcp Mg, and find remarkable agreement with experimental data. A strong impact of anharmonicity is observed specifically for Nb. The introduced procedure paves the way for the development of ab initio thermodynamic databases.