Temperature distribution at depth is of key importance for characterizing the crust, defining its mechanical behavior and deformation. Temperature can be retrieved by heat flow measurements in boreholes that are sparse, shallow, and have limited reliability, especially in active and recently active areas. Laboratory data and thermodynamic modeling demonstrate that temperature exerts a strong control on the seismic properties of rocks, supporting the hypothesis that seismic data can be used to constrain the crustal thermal structure. We use Rayleigh wave dispersion curves and receiver functions, jointly inverted with a transdimensional Monte Carlo Markov Chain algorithm, to retrieve the VS and VP/VS within the crust in the Italian peninsula. The high values (〉1.9) of VP/VS suggest the presence of filled-fluid cracks in the middle and lower crust. Intracrustal discontinuities associated with large values of VP/VS are interpreted as the α−β quartz transition and used to estimate geothermal gradients. These are in agreement with the temperatures inferred from shear wave velocities and exhibit a behavior consistent with the known tectonic and geodynamic setting of the Italian peninsula. We argue that such methods, based on seismological observables, provide a viable alternative to heat flow measurements for inferring crustal thermal structure.