A density-functional-theory based approach to efficiently compute numerically exact vibrational free energies - including anharmonicity - for chemically complex multicomponent alloys is developed. It is based on a combination of thermodynamic integration and a machine-learning potential. We demonstrate the performance of the approach by computing the anharmonic free energy of the prototypical five-component VNbMoTaW refractory high entropy alloy.