We explore the range of applicability of the nuclear ensemble method (NEM) for quantitative simulations of absorption spectra and their temperature variations. We formulate a “good practice” for the NEM based on statistical theory. Special attention is paid to proper treatment of uncertainty estimation including the convergence with the number of samples, which is often neglected in the field. As a testbed, we have selected a well-known chromophore, (E)-azobenzene. We measured its temperature difference UV–vis absorption spectra in methanol, which displayed two dominant features: a moderate increase in the intensity of the nπ* band and a pronounced decrease in intensity of the low-energy part of the ππ* band. We attributed both features to increasing non-Condon effects with temperature. We show that the NEM based on the path integral molecular dynamics combined with range-separated hybrid functionals provides quantitatively accurate spectra and their differences. Experimentally, the depletion of the absorption in the ππ* band showed a characteristic vibrational progression that cannot be reproduced with the NEM. We show that hundreds of thousands of samples are necessary to achieve an accuracy sufficient for the unambiguous explanation of the observed temperature effects. We provide a detailed analysis of the temperature effects on the spectrum based on the harmonic model of the system combined with the NEM. We also rationalize the vibrational structure of the spectrum using the Franck–Condon principle.