We report on the first fully differential calculation of the next-to-next-to-leading-order (NNLO) QCD radiative corrections to the production of bottom-quark pairs at hadron colliders. The calculation is performed by using the qT subtraction formalism to handle and cancel infrared singularities in real and virtual contributions. The computation is implemented in the Matrix framework, thereby allowing us to efficiently compute arbitrary infrared-safe observables in the four-flavour scheme. We present selected predictions for bottom-quark production at the Tevatron and at the LHC at different collider energies, and we perform some comparisons with available experimental results. We find that the NNLO corrections are sizeable, typically of the order of 25–35%, and they lead to a significant reduction of the perturbative uncertainties. Therefore, their inclusion is crucial for an accurate theoretical description of this process.