Reconstructing historical climate change from deep ground temperature measurements in cold regions is often complicated by the presence of permafrost. Existing methods are typically unable to account for latent heat effects due to the freezing and thawing of the active layer. In this work, we propose a novel method for reconstructing historical ground surface temperatures (GST) from borehole temperature measurements that accounts for seasonal thawing and refreezing of the active layer. Our method couples a recently developed fast numerical modeling scheme for two-phase heat transport in permafrost soils with an ensemble-based method for approximate Bayesian inference. We evaluate our method on two synthetic test cases covering both cold and warm permafrost conditions as well as using real data from a 100m deep borehole on Sardakh Island in northeastern Siberia. Our analysis of the Sardakh Island borehole data confirms previous findings that ground surface temperatures in the region have likely risen by 5 to 9°C between the pre-industrial period of 1750–1855 and 2012. We also show that latent heat effects due to seasonal freeze-thaw have a substantial impact on the resulting reconstructed surface temperatures. We find that neglecting the thermal dynamics of the active layer can result in biases of roughly -1 to -1.5°C in cold conditions (i.e. mean annual ground temperature below -5°C) and as much as -2 to -3°C in warmer conditions where substantial active layer thickening (>200cm) has occurred. Our results highlight the importance of considering seasonal freeze-thaw in GST reconstructions from permafrost boreholes.