Infrared (IR) spectroscopy has been widely utilized for the study of protein folding, unfolding, and misfolding processes. We have previously developed a theoretical method for calculating IR spectra of proteins in the amide I region. In this work, we apply this method, in combination with replica-exchange molecular dynamics simulations, to study the equilibrium thermal unfolding transition of the villin headpiece subdomain (HP36). Temperature-dependent IR spectra and spectral densities are calculated. The spectral densities correctly reflect the unfolding conformational changes in the simulation. With the help of isotope labeling, we are able to capture the feature that helix 2 of HP36 loses its secondary structure before global unfolding occurs, in agreement with experiment.