In this Communication, we study the effect of spurious oscillations in the profiles of energy derivatives with respect to nuclear coordinates calculated with density functional approximations (DFAs) for formaldehyde, pyridine, and furan in their ground and electronic excited states. These spurious oscillations, which can only be removed using extensive integration grids that increase enormously the CPU cost of DFA calculations, are significant in the case of third- and fourth-order energy derivatives of the ground and excited states computed by M06-2X and ωB97X functionals. The errors in question propagate to anharmonic vibronic spectra computed under the Franck-Condon approximation, i.e., positions and intensities of vibronic transitions are affected to a large extent (shifts as significant as hundreds of cm-1 were observed). On the other hand, the LC-BLYP and CAM-B3LYP functionals show a much less pronounced effect due to spurious oscillations. Based on the results presented herein, we recommend either LC-BLYP or CAM-B3LYP with integration grids (250, 974) (or larger) for numerically stable simulations of vibronic spectra including anharmonic effects.