Laurentide Ice Sheet basal temperatures during the last glacial cycle as inferred from borehole data

Thirteen temperature–depth profiles (≥ 1500 m) measured in boreholes in eastern and central Canada were inverted to determine the ground surface temperature histories during and after the last glacial cycle. The sites are located in the southern part of the region that was covered by the Laurentide Ice Sheet. The inversions yield ground surface temperatures ranging from−1.4 to 3.0 C throughout the last glacial cycle. These temperatures, near the pressure melting point of ice, allowed basal flow and fast flowing ice streams at the base of the Laurentide Ice Sheet. Despite such conditions, which have been inferred from geomorphological data, the ice sheet persisted throughout the last glacial cycle. Our results suggest some regional trends in basal temperatures with possible control by internal heat flow.


Introduction
The impact of future climate change on the stability of the present-day ice sheets in Greenland and Antarctica is a major concern of the scientific community (e.g., Gomez et al., 2010;Mitrovica et al., 2009).Satellite gravity measurements performed during the Gravity Recovery and Climate Experiment (GRACE) mission suggest that the mass loss of the Greenland and Antarctic glaciers has accelerated during the decade 2002-2012 (Velicogna and Wahr, 2013).Over the past 2 decades, mass loss from the Greenland Ice Sheet has quadrupled and contributed to a fourth of global sea level rise from 1992 to 2011 (Church et al., 2011;Straneo and Heim-bach, 2013).In 2014, two teams of researchers noted that the collapse of the Thwaites Glacier basin, an important component holding together the West Antarctic Ice Sheet, was potentially underway (Joughin et al., 2014;Rignot et al., 2014).The collapse of the entire West Antarctic Ice Sheet would lead to rise in sea level by at least 3 m.The present observations of the ice sheet mass balance are important.However, to predict the effects of future climate change on the ice sheets, it is necessary to fully understand the mechanisms of ice sheet growth, decay, and collapse throughout the past glacial cycles.The models of ice sheet dynamics during past glacial cycles show that the thickness and elevation of the ice sheets and the thermal conditions at their base are key parameters controlling the basal flow regime and the evolution of the ice volume (Marshall and Clark, 2002;Hughes, 2009).
During the last glacial cycle (LGC), ∼ 120 000-12 000 yr BP, large ice sheets formed in the Northern Hemisphere, covering Scandinavia and almost all of Canada (Denton and Hughes, 1981;Peltier, 2002Peltier, , 2004;;Zweck and Huybrechts, 2005).The growth and decay of ice sheets are governed mainly by ice dynamics and ice-sheet-climate interactions (Oerlemans and van der Veen, 1984;Clark, 1994;Clark and Pollard, 1998).Ice dynamics are also strongly controlled by the underlying geological substrate and associated processes (Clark et al., 1999;Marshall, 2005).(Clark and Pollard, 1998) studied how the ice thickness and the nature of the geological substrate control flow at the base of the glacier.They showed that soft beds, beds of unconsolidated sediment at relatively low relief, result Published by Copernicus Publications on behalf of the European Geosciences Union.
in thin ice sheets that are predisposed to fast ice flow when possible.Hard beds, beds of high-relief crystalline bedrock, on the other hand, provide the ideal conditions for the formation of larger, thick ice sheets that experience stronger bed-ice-sheet coupling and slow ice flow.Basal flow rate is also affected by the basal temperatures, which play a key role in determining the velocity.However, the ice sheet evolution models use basal temperatures that are poorly constrained because of the lack of direct data pertaining to thermal conditions at the base of ice sheets.The objective of the present study is to use borehole temperature-depth data to estimate how the temperatures varied at the base of the Laurentide Ice Sheet during the last glacial cycle.
Temporal variations in ground surface temperatures (GST) are recorded by Earth's subsurface as perturbations to the "steady-state" temperature profile (see, e.g., Hotchkiss and Ingersoll, 1934;Birch, 1948;Beck, 1977;Lachenbruch and Marshall, 1986).With no changes in GST, the thermal regime of the subsurface is governed by the outflow of heat from Earth's interior resulting in a profile where temperature increases with depth.In homogeneous rocks without heat sources, the "equilibrium" temperature increases linearly with depth.When changes in ground surface temperature occur and persist, they are diffused downward and recorded as perturbations to the semi-equilibrium thermal regime.These perturbations are superimposed on the temperature-depth profile associated with the flow of heat from Earth's interior.For periodic oscillations of the surface temperature, the temperature fluctuations decrease exponentially with depth such that where z is the depth, t is time, κ is the thermal diffusivity of the rock, ω is the frequency, and T is the amplitude of the temperature oscillation.Long time-persistent transients affect subsurface temperature to great depths.(Hotchkiss and Ingersoll, 1934) were the first to attempt to infer past climate from such temperature-depth profiles; specifically, they estimated the timing of the last glacial retreat from temperature measurements in the Calumet copper mine in northern Michigan.(Birch, 1948) estimated the perturbations to the temperature gradient caused by the last glaciation and suggested a correction to heat flux determinations in regions that had been covered by ice during the LGC.A correction including the glacial-interglacial cycles of the past 400 000 years was proposed by (Jessop, 1971) to adjust the heat flow measurements made in Canada.It was only in the 1970s that systematic studies were undertaken to infer past climate from borehole temperature profiles (Cermak, 1971;Sass et al., 1971;Beck, 1977).The use of borehole temperature data for estimating recent (< 300 years) climate changes became widespread in the 1980s because of concerns about increasing global temperatures (Lachenbruch and Marshall, 1986;Lachenbruch, 1988).High-precision borehole temper-ature measurements have been mostly made for estimating heat flux in relatively shallow (a few hundred meters) holes drilled for mining exploration.Such shallow boreholes are suitable for studying recent (< 500 years) climate variations and the available data have been interpreted in many regional and global studies (e.g., Bodri and Cermak, 2007;Jaupart and Mareschal, 2011, and references therein).However, very few deep (≥ 1500 m) borehole temperature data are available to study climate variations on the timescale of 10 to 100 kyr.Oil exploration wells are usually a few kilometers deep but are not suitable because the temperature measurements are not made in thermal equilibrium and lack the required precision.Nonetheless, a few deep mining exploration holes have been drilled, mostly in Precambrian shields, where temperature measurements can be used for climate studies on the timescale of the LGC.In Canada, (Sass et al., 1971) measured temperatures in a deep (3000 m) borehole near Flin Flon, Manitoba, and used direct models to show that the surface temperature during the Last Glacial Maximum (LGM), ∼ 20 000 yr BP, could not have been more than 5 K colder than present.Other Canadian deep-borehole temperature profiles have since been measured revealing regional differences in temperatures during the LGM.From a deep borehole in Sept-Îles, Québec, (Mareschal et al., 1999b) found surface temperatures to have been approximately 10 K colder than in the present day.This was confirmed by (Rolandone et al., 2003b) who studied four deep holes and suggested that LGM surface temperatures were colder in eastern Canada than in central Canada.(Chouinard and Mareschal, 2009) examined eight deep boreholes located in central to eastern Canada and observed significant regional differences in heat flux, temperature anomalies and ground surface temperature histories.(Majorowicz et al., 2014) have studied a 2400 m deep well near Fort McMurray, Alberta, penetrating the basement.They interpreted the variations in heat flux as being due to a 9.6 K increase in surface temperature at 13 ka.(Majorowicz and Šafanda, 2015) have revised this estimate to include the effect of the variations in heat productions with depth which is non-negligible in granitic rocks with high heat production.They concluded that temperatures at the base of the ice sheet in Northern Alberta were about −3 • C during the LGM.On the other hand, studies of deep boreholes in Europe lead to different conclusions for the Fennoscandian Ice Sheet which covered parts of Eurasia during the LGC (see, e.g., Demezhko and Shchapov, 2001;Kukkonen and Jõeleht, 2003;Šafanda et al., 2004;Majorowicz et al., 2008;Demezhko et al., 2013;Demezhko and Gornostaeva, 2015) .(Kukkonen and Jõeleht, 2003) analyzed heat flow variations with depth in several boreholes from the Baltic Shield and the Russian Platform and found a 8 ± 4.5 K temperature increase following the LGM.(Demezhko and Shchapov, 2001) studied a ∼ 5 km deep borehole in the Urals, Russia, and found a postglacial warming of 12-13 K, with basal temperatures below the melting point of ice during the LGM.This was confirmed by recent work indicating that temperatures in the Urals were ∼ −8 • C at the LGM (Demezhko and Gornostaeva, 2015).
In this study, we shall examine all the deep-borehole temperature profiles measured in central and eastern Canada in order to determine the temperature at the base of the Laurentide Ice Sheet, which covered the area during the LGC.The geographical extent of the study is confined to the southern portion of the Laurentide Ice Sheet because deep mining exploration boreholes have only been drilled in the southernmost part of the Canadian Shield.

Theory
For calculating the temperature-depth profile, we assume that heat is transported only by vertical conduction and that the temperature perturbation is the result of a time-varying, horizontally uniform, surface temperature boundary condition.The temperature at depth in the Earth, for a homogeneous half-space with horizontally uniform variations in the surface temperature, can be written as where T o is the reference ground surface temperature (steady-state or long-term surface temperature), Q o is the reference heat flux (steady-state heat flux from depth), λ(z) is the thermal conductivity, z is depth, H is the rate of heat generation, and T t (z) is the temperature perturbation at depth z due to time-varying changes to the surface boundary condition.R(z) is the thermal resistance to depth z, which is defined as Thermal conductivity is measured on core samples, usually by the method of divided bars (Misener and Beck, 1960).Radiogenic heat production is also measured on core samples.
The temperature perturbation at depth z resulting from surface temperature variations can be written as (Carslaw and Jaeger, 1959) where t is time before present, κ is thermal diffusivity and T o (t) is the surface temperature at time t (Carslaw and Jaeger, 1959).For a stepwise change T in surface temperature at time t before present, the temperature perturbation at depth z is given by (Carslaw and Jaeger, 1959) where erfc is the complementary error function.
If the ground surface temperature variation is approximated by a series of constant values T k during K time intervals (t k−1 , t k ), the temperature perturbation can be written as follows: The T k represent the departure of the average GST during each time interval from T o .

First-order estimate of the GST history
We have used variations in the long-term surface temperature with depth directly to obtain a first-order estimate of time variations in surface temperature.To interpret subsurface anomalies as records of GST history variations, we must separate the climate signal from the quasi steady-state temperature profile.The quasi steady-state thermal regime is estimated by a least-squares linear fit to the lowermost 100 m section of the temperature-depth profile.The slope and the surface intercept of the fitted line are interpreted as the steady-state temperature gradient o and the long-term surface temperature T o , respectively.The bottom part of the profile is the section least affected by the recent changes in surface temperature and dominated by the steady-state heat flow from Earth's interior.When we are using a shallower section of the profile, the estimated long-term temperature and gradient are more affected by recent perturbations in the surface temperature.The shallower the section, the more recent the perturbation.Thus, we determine how the estimated long-term surface temperature varies with time by iterating through different sections of the temperature-depth profile.
Using the scaling of Eq. ( 5), the time it takes for the surface temperature variation to propagate into the ground and the depth of the perturbations are related by where κ is thermal diffusivity, z is the depth, and t is time.This approach yields only a first-order estimate because the temperature perturbations are attenuated as they are diffused downward.However, the temporal variation in the long-term surface temperature is obtained by extrapolating the semi-equilibrium profile to the surface.The effect of a small change in gradient on the long-term surface temperature is proportional to the average depth where the gradient is estimated.Variations in the calculated long-term surface temperature are thus less attenuated than the perturbations in the profile and may well provide a gross approximation of the GST history but be independent of the model assumptions in an inversion procedure.
C. Pickler et al.: Laurentide Ice Sheet basal temperatures during the last glacial cycle

Inversion
In order to obtain a more robust estimate of the GST histories, we have inverted the temperature-depth profiles for each site, and whenever profiles from nearby sites are available, we have inverted them jointly.The inversion consists of determining T o , Q o , and the T k from the temperaturedepth profile (Eq.6).For each depth where temperature is measured, Eq. ( 6) yields a linear equation in T k .For N temperature measurements, we obtain a system of N linear equations with K + 2 unknowns, T o , Q o , and K values of T k .Even when N ≥ K + 2, the system seldom yields a meaningful solution because it is ill-conditioned.This means that the solution is unstable and a small error in the data results in a very large error in the solution (Lanczos, 1961).Different authors have proposed different inversion techniques to obtain the GST history from the data (see, e.g., Vasseur et al., 1983;Shen andBeck, 1983, 1991;Nielsen and Beck, 1989;Wang, 1992;Mareschal and Beltrami, 1992;Clauser and Mareschal, 1995;Mareschal et al., 1999b).Here, to obtain a solution regardless of the number of equations and unknowns and to reduce the impact of noise and errors on this solution, singular-value decomposition is used (Lanczos, 1961;Jackson, 1972;Menke, 1989).This technique is well documented, and further details can be found in (Mareschal and Beltrami, 1992), (Clauser and Mareschal, 1995), (Beltrami and Mareschal, 1995), and (Beltrami et al., 1997).

Simultaneous inversion
As meteorological trends remain correlated over a distance on the order of 500 km (Hansen and Lebedeff, 1987), boreholes within the same region are assumed to have been affected by the same surface temperature variations and their subsurface temperature anomalies are expected to be consistent.This holds only if the surface conditions are identical for all boreholes.If these conditions are met, jointly inverting different temperature-depth profiles from the same region will increase the signal-to-noise ratio.Singular-value decomposition was used to jointly invert the sites with multiple boreholes.A detailed description and discussion of the methodology can be found in papers by (Beltrami and Mareschal, 1995) and (Beltrami et al., 1997).

Data description
We have used 13 deep boreholes (≥ 1500 m) across eastern and central Canada to determine the temperature throughout and after the LGC.All the sites are located in the southern portion of the Canadian Shield, which was covered by the Laurentide Ice Sheet, which extended over most of Canada during the LGC (Fig. 1).Borehole locations and depths are summarized in Table 1.All the holes that we logged were drilled for mining exploration.The only exception is Flin Flon that was drilled to be instrumented with low-noise seis- mometers for monitoring nuclear tests.Detailed descriptions of the measurement techniques as well as the relevant geological information can be obtained from heat flow publications (Sass et al., 1971;Mareschal et al., 1999b, a;Rolandone et al., 2003a, b;Perry et al., 2006Perry et al., , 2009;;Jaupart et al., 2014).We shall briefly summarize the main steps of the measurement technique, which has been described in several papers (e.g., Mareschal et al., 1989;Perry et al., 2006;Lévy et al., 2010).Temperature is measured at 10 m intervals along the hole by lowering a calibrated probe with a thermistor.The precision of the measurements is better than 0.005 K with an overall accuracy estimated to be of the order 0.02 K Thermal conductivity is measured on core samples.Samples are collected in every lithology with an average of one sample per 80 m.The thermal conductivity is measured with a divided bar apparatus on five cylinders of the core with thickness varying between 0.2 and 1.0 cm.This method based on five measurements on relatively large core samples provides the best estimate of the thermal conductivity of the bulk rock and is unaffected by small heterogeneities.Samples of the core have also been analyzed for heat production following the method described by (Mareschal et al., 1989).
A description of eight sites -Flin Flon, Pipe, Manitouwadge 0610, Manitouwadge 0611, Balmertown, Falconbridge, Lockerby, and Sept-Îles -can be found in (Rolandone et al., 2003b) and (Chouinard and Mareschal, 2009)  We found systematic variations in thermal conductivity at Flin Flon, Thompson (Owl), and Matagami, and we corrected them accordingly (Bullard, 1939).We calculated the thermal resistance and obtained a temperature vs. thermal resistance profile that is almost linear.We calculated the heat flux as the slope of the temperature resistance and found no discontinuity in heat flux along the profile (Table 1).For all the other sites that show no systematic variations in conductivity, we have used the mean thermal conductivity to calculate heat flux.Heat production was measured and found to only be significant at two sites -Lockerby (3 µW m −3 ) and Victor Mine (0.9 µW m −3 ) -and was therefore only taken into account at these sites.Variations in heat production with depth may affect the temperature profiles and the GST (Majorowicz and Šafanda, 2015).Because systematic variations in heat production with depth were not present at these sites, they were not included in the inversion.In the absence of heat production and in steady state, the heat flux does not vary with thermal resistance.Variations in heat flux with thermal resistance (or depth) are thus a diagnostic of departure from a 1-D steady-state thermal regime.A decrease in heat flux toward the surface is associated with surface warming, and enhanced heat flux is due to cooling.The heat flux profiles that we have calculated for all the sites (Fig. 2) exhibit clear departures from the 1-D steady-state condition.The heat flux was calculated as the product of the temperature gradient and the thermal conductivity within each interval.No smoothing was applied.As expected, the gradient profile contains highfrequency variations as the gradient always amplifies noise and errors in the temperature measurements.Some holes appear to be noisier at depth near the exploration targets because of small-scale conductivity variations due to the presence of mineralization.Furthermore, the inclination of some of the holes decreases markedly at depth resulting in larger errors in the gradient.For example, the inclination at Val d'Or, 85 • at the collar, was only 20 • near the end of the hole.Most of the profiles show a very pronounced increase in heat flux with depth at shallow depth (< 200 m) and a clear trend of increasing heat flux between 500 and 1500 m.The increase at shallow depth is related to very recent (< 300 years) climate warming.The trend between 500 and 1500 m is the result of the surface warming that followed the glacial retreat at ca. 10 ka.Heat flux is calculated as the product of thermal conductivity by the temperature gradient calculated over three points.The Flin Flon, Owl and Matagami profiles have been corrected to account for thermal conductivity variations with depth as shown in Table 1.

Long-term surface temperatures
Estimated long-term surface temperatures as a function of time and depth (i.e., ∝ depth squared) were determined for each borehole (Fig. 3).The time in these plots represents the time it took the signal to propagate and not the time that the surface temperature perturbation occurred.The range of long-term surface temperatures for each borehole was estimated over its sampled depth and revealed the persistent long-term climate trends.These trends are the mirror image of the heat flow trends.A decreasing temperature trend with time and depth is apparent in all the boreholes in Manitoba except Pipe, which does not show a clear trend.Variations in surface temperature are not consistent, with almost no trend at Pipe and a weak signal at Owl.The Flin Flon borehole is the deepest available for this study and provides a history 4 times longer than that of the two shallower boreholes at Owl and Pipe.It is consistent with a colder period coinciding with the LGM.In western Ontario, Manitouwadge 0611 exhibits very strong oscillations at depth.The source of noise is difficult to ascertain because of the complicated geological structure and the absence of thermal conductivity data.Borehole 0610 at Manitouwadge is consistent with colder temperatures during the LGM, but the Balmertown hole does not show any variation in long-term surface temperature.Four profiles at Sudbury, are consistent with colder temperature during the LGM, but the amplitude of the trend varies between sites.Similar trends are observed for the three easternmost holes in Québec with long-term temperatures on average 5 K lower near the bottom of the holes than near the surface.
These first-order estimates suggest colder surface temperatures during the LGM at most of the sites.In order to better quantify the surface temperature changes, we must turn to inversion and obtain the GST histories.

Individual inversions
The GST histories at all the studied sites for the time period of 100 to 100 000 yr BP were inverted from the temperaturedepth profiles.The time span of the GST history model consists of 16 intervals whose distribution varies logarithmically because the resolution decreases with time.A singular-value cutoff of 0.08 was used for all the individual profile inversions (Figs.4-7).The singular-value cutoff eliminates the part of the solution that is affected by noise and effectively introduces a smoothing constraint on the solution (Mareschal and Beltrami, 1992).A summary of the inversion results can be found in Table 2. Two main episodes can be recognized in the GST histories: one is associated with a minimum temperature that occurred around the LGM at ca. 20 ka.The second is a warming observed at ca. 2-6 ka coinciding with the Holocene climatic optimum (HCO), a warm period that followed the deglaciation (Lamb, 1995).
For the purpose of the discussion, we have grouped the sites that are from the same geographical region.We shall thus distinguish between Manitoba, western Ontario, the Sudbury area, and Québec.
The sites from Manitoba -Flin Flon, Owl and Pipehave not recorded a very strong signal and do not exhibit common regional trends as observed in (Fig. 4).The lack of regional trends is expected as Thompson (Owl and Pipe) and Flin Flon sites are ∼ 300 km apart and the present-day ground surface temperatures differ by 3 K.The weak signal recorded may in part be due to the present ground surface temperature being very low; it is close to 0 • C in Thompson where intermittent permafrost is found.At Flin Flon, where the present ground temperature is near 3 • C, we found that ground surface temperature variations were small, confirming previous studies by (Sass et al., 1971).The surface temperature was at a minimum around the LGM and was near the melting point of ice (−0.3 • C).For Pipe, little to no signal was recorded.We found that there was minimal change in the ground surface temperature (∼ 2.5 K in amplitude) over the past 100 000 years.For Owl, the amplitude of the temperature changes is double (∼ 5 K), with the minimum temperature (−2.4 • C) around the LGM and a warming around the HCO.Although this result is plausible, some uncertainty remains as (Guillou-Frottier et al., 1996) noted that high heat flux is correlated with high thermal conductivity in the Thompson Belt.The elevated thermal conductivity is due to the presence of vertical slices of quartzites, which increase thermal conductivity by a factor of 1.7.Furthermore, the site is affected by a poorly resolved conductivity structure as thermal conductivity measurements vary in the deepest part of the borehole, between 2.21 and 5.14 W m −1 K −1 .The lateral heat refraction effects due to the thermal conductivity contrast affect the temperature profiles and alter the GST history.2.
For these reasons, we have little confidence in the robustness of the Owl GST history reconstruction.In western Ontario (Balmertown and Manitouwadge 0610 and 0611), results at all sites show minimum temperatures around the LGM, with a very weak minimum at Balmertown (Fig. 5).However, the amplitude of the temperature change is much larger at Manitouwadge 0611 (∼ 10 K) than at Balmertown (∼ 2 K) and Manitouwadge 0610 (∼ 3 K).As the two Manitouwadge sites are only ∼ 40 km apart, this difference is surprising.While Manitouwadge 0611 yields a plausible GST history, it appears to have an amplified signal.The site is located in a complex geological structure and lacks thermal conductivity data.There is also a change in the temperature gradient at 500 m, which cannot be explained.In the absence of thermal conductivity data, we cannot consider the GST history for Manitouwadge 0611 as reliable.
All four sites in the Sudbury region -Craig Mine, Falconbridge, Lockerby, and Victor Mine -have recorded minimal temperatures around the LGM (Fig 6).However, the minimum past temperatures for the region do vary between sites.The coldest minimum temperature was found at Craig Mine.The amplitude of the temperature changes for the site (∼ 12 K) is much larger than those of the other sites, which vary between ∼ 5 and 7 K.This difference is unexpected as these sites are all within the Sudbury structure and should have recorded similar histories.The Craig Mine signal appears amplified, which could be the result of water flows induced by pumping at levels below 2000 m in the mine.We thus believe that the GST history for Craig Mine is not reliable.The minimum temperatures at Lockerby and Victor Mine, 2.8 and 3.0  The temperature profile at this site may be distorted by horizontal contrasts in thermal conductivity. 2The temperature profile in the lowermost part of the hole may be affected by subvertical layering and thermal conductivity contrasts. 3The temperature profile may be affected by water flow caused by pumping in the nearby mine.
Figure 5. Ground surface temperature history for the western Ontario boreholes: Balmertown and Manitouwadge 0610 and 0611.
The temperatures have been shifted with respect to the reference surface temperature of the site, T o , as shown in Table 2.
over, these are the only two sites with non-negligible heat production: 3 and 0.9 µW m −3 , respectively.The corrections for heat production produce an increase in the temperature gradient proportional to depth and result in an amplification of the warming signal in the profile.Consequently, the minimum and maximum temperatures would be higher at these sites -Lockerby and Victor Mine -than at those with negligible heat production.
The GST histories from the three boreholes in Québec, Matagami, Val d'Or and Sept-Îles, display regional differences (Fig. 7).However, for all three sites, the minimum temperatures are synchronous and occurred around the LGM.For all the sites used in this study, the lowest minimum temperature occurred in this region (−1.4 • C at Sept-Îles).
At all sites, excluding Pipe, we found that the ice retreat was followed by a warm episode that can be associated with the HCO, a warm period whose maximal temperatures have been dated to 4.4-6.8ka with palynological reconstructions from northern Ontario and northern Michigan (Boudreau et al., 2005;Davis et al., 2000).
We have compared the ranges of temperature in the inverted GST histories with those of the long-term surface temperature variations for all the sites (Table 4).Although the total range varies between sites from ∼ 2 to ∼ 8 K, the two methods yield consistent values that differ by less than 1 K at most of the sites.The warming trend in the long-term surface temperature variations is consistent with the inverted GST histories.This correlation between the long-term surface temperature ranges and the persistent long-term GST history trends suggests that our results are robust.The inversion of borehole profiles has a very limited resolving power, and short-period oscillations can seldom be recovered.In practice, the duration of an episode must be proportional (1/3 to 1/2) to the time when it occurred.The last glacial period is easily identified, but shorter-period events such as the HCO that lasted for 1-2 kyr around 5 ka are just beyond the threshold of resolution.It is thus possible that constraining the GST   to include the HCO would result in colder GST during the LGM, as suggested by Will Gosnold in his comment on the discussion paper and in Gosnold et al. (2011).This would require obtaining proxy data from sites close to Flin Flon and Thompson and requires documenting the temperature conditions at the bottom of Lake Agassiz, which covered most of Manitoba after the glacial retreat.

Simultaneous inversion
We have inverted the boreholes of Thompson (Owl and Pipe), Manitouwadge (0610 and 0611) and Sudbury (Craig Mine, Falconbridge, Lockerby, and Victor Mine) simultaneously to observe regional trends and to improve the signalto-noise ratio (Fig. 8).The simultaneous inversion results are summarized in Table 3.For simultaneous inversion, the temperature-depth profiles were truncated to ensure a com-C.Pickler et al.: Laurentide Ice Sheet basal temperatures during the last glacial cycle mon depth for all the boreholes.This ensures consistency and facilitates comparison because we are examining the subsurface temperature anomalies for the same time period.For all inversions, the minimum temperature occurs around the LGM and the deglaciation is followed by a warming associated with the HCO.Although Owl and Manitouwadge 0611 have questionable individual inversions, we have still performed simultaneous inversions for the regions of Thompson and Manitouwadge in an attempt to decrease the signalto-noise ratio.For Thompson, the width of the GST history temperature range is ∼ 4 K.As Pipe did not appear to record a signal, the simultaneous inversion appears to have damped the questionable signal recorded at Owl slightly.The amplitude difference for the Manitouwadge GST history is ∼ 8 K, which lies between that of 0610 (∼ 3 K) and 0611 (∼ 10 K).The Sudbury simultaneous inversion was performed with and without Craig Mine.This was done to check whether the inclusion of Craig Mine, a site which we suspect to have been affected by water flow, affects the reconstructions.Both inversions display similar trends; however, there is a difference in the range of temperature variations.The inversion excluding Craig Mine yielded a T of 7 K, similar to those of the individual inversions of Falconbridge (∼ 7 K), Lockerby (∼ 7 K) and Victor Mine (∼ 5 K).Upon inclusion of Craig Mine, T was increased to 11 K, demonstrating that the presence of a signal due to water flow in the temperature profile at this site has a strong effect on the results of the simultaneous inversion.

Discussion
The minimum temperatures of the GST histories occur around the LGM, representing the basal temperatures of the Laurentide Ice Sheet.These temperatures vary spatially and range from −1.4 to 3.0 • C, near or above the pressure melting point of ice.Such spatial variation is expected as numerous studies have demonstrated present-day spatial basal temperature variability beneath the Antarctic and Greenland Ice Sheets (see, e.g., Dahl-Jensen et al., 1998;Pattyn, 2010;Schneider et al., 2006).highest basal temperatures occur within the Sudbury basin at Lockerby and Victor Mine.The Sudbury region has the highest average heat flux of the Canadian Shield, ∼ 54 mW m −2 , because crustal heat production is higher than average (Perry et al., 2009).The lowest basal temperature is recorded at Sept-Îles, where the heat flux is the lowest of the studied regions, ∼ 34 mW m −2 .These correlations suggest a link between heat flux and basal temperatures.This is further supported by modeling work showing that heat flux influences the thermal structure and properties of ice sheets, including basal temperatures and ice flow (Greve, 2005;Pollard et al., 2005;Llubes et al., 2006).However, the Sept-Îles basal temperature has been also been linked to its proximity to the edge of the ice sheet and to an area of thinner ice (Rolandone et al., 2003b).Our study is consistent with a Table 4. Ranges in surface temperature variations estimated from the iteration of the long-term surface temperature as a function of depth (column 1) and from inversion of the GST history (column 2), along with the difference between the two (column 3).

Site
T o T (GSTH) Difference ( The temperature profile at this site may be disturbed by horizontal contrasts in thermal conductivity. 2The temperature profile in the lowermost part of the hole may be affected by subvertical layering and thermal conductivity contrasts. 3The temperature profile may be affected by water flow caused by pumping in the nearby mine. possible link between basal temperature and heat flux along with ice dynamics in the Laurentide Ice Sheet during the LGC, but further modeling work is necessary to confirm such a relationship.Variations in these parameters (heat flux and ice dynamics) could account for the differences observed in the basal temperatures beneath the Fennoscandian and Laurentide ice sheets.The basal temperatures recorded, near the pressure melting point of ice, indicate the possibility of basal flow and ice streams, two important factors affecting ice sheet evolution.These processes have been suggested by geomorphological evidence presented by (Dyke et al., 2002) and predicted by the ICE-5G model (Peltier, 2004).Basal flow has the ability to transport large amounts of water from the interior of the ice sheet, leading to thinning of the ice sheet and climatically vulnerable ice.It is postulated to be a key factor in glacial terminations (Marshall and Clark, 2002).Furthermore, these temperatures demonstrate that the southern portion of the Laurentide Ice Sheet was not frozen to the bed, suggesting basal sliding.These conditions can lead to instability.Widespread basal sliding and increased surface meltwater could have been a factor resulting in the rapid collapse of the Laurentide Ice Sheet (Zwally et al., 2002).However, these basal temperatures and the associated melt persisted prior to and throughout the LGM for more than 30 000 years with deglaciation only occurring rapidly during the early Holocene (Carlson et al., 2008).While this indicates that basal temperatures near the pressure melting point of ice cannot be solely responsible for ice sheet instability and col-lapse, it demonstrates that they are a key parameter in ice sheet evolution and one that ice sheet evolution models must take into account.Elevated basal temperatures, near or above the pressure melting point of ice, have been recorded in the present-day ice sheets (Fahnestock et al., 2001;Pritchard et al., 2012).Our results indicate that this alone cannot be considered as an indication of ice sheet collapse.However, combined with other processes it could lead to instability and collapse.

Conclusions
Thirteen deep boreholes from eastern and central Canada were analyzed to determine the GST histories for the last 100 kyr.The long-term trends are consistent between sites.A warm period following the retreat of the ice sheet is inferred at ∼ 2-6 ka in the inverted GST histories, correlated to the Holocene climatic optimum.
The surface temperatures reached their minima during the LGM and post glacial warming started at ca. 10 ka.The corresponding temperatures at the base of the Laurentide Ice Sheet range from −1.4 to 3.0 • C and are all near or above the pressure melting point of ice.Such temperatures allow for basal flow and fast-flowing ice streams, two important factors affecting ice sheet evolution, illustrating the need for models of ice sheet evolution to account for such a key parameter as basal temperature.Despite the suggestion that melting took place at its base, the Laurentide ice sheet persisted throughout the LGM for over more than 30 000 years.This demonstrates that basal temperatures near the melting point of ice do not indicate that an ice sheet is on the verge of collapse.However, combined with other processes it could lead to instability and collapse.
The differences between GST histories at different sites raise other questions concerning the controls on temperatures at the base of ice sheets.The equilibrium between heat flow from the Earth's interior and heat advection by glacial flow determines the temperature at the boundary between the ice and the bedrock.The correlation between higher heat flux and higher basal temperatures in the Sudbury region suggests that variations in crustal heat flux may account for some of the regional differences in basal temperatures along with the dynamics of ice thickness controlled by the accumulation rate and the distance to the edge of the ice sheet.
It is also noteworthy that similar deep-borehole studies in Europe suggest that basal temperatures beneath the Fennoscandian Ice Sheet and in the Urals during the LGC were much colder (Kukkonen and Jõeleht, 2003;Demezhko et al., 2013) than those observed in Canada (Chouinard and Mareschal, 2009).Because of the geological similarities between the two regions, this contrast is likely to be due to differences in climate and ice dynamics between Europe and North America during the LGC.

Figure 1 .
Figure 1.Map of central and eastern Canada and adjoining USA showing the location of sampled boreholes.Thompson (Owl and Pipe), Manitouwadge (0610 and 0611) and Sudbury (Falconbridge, Lockerby, Craig Mine, and Victor Mine) have several boreholes present within a small region.The number of profiles available at locations with multiple holes is enclosed in parenthesis.

Figure 2 .
Figure2.Heat flux variation as a function of depth.Heat flux is calculated as the product of thermal conductivity by the temperature gradient calculated over three points.The Flin Flon, Owl and Matagami profiles have been corrected to account for thermal conductivity variations with depth as shown in Table1.

Figure 3 .
Figure 3.Long-term surface temperature variations over time (left y axis) and depth (right y axis) for all the boreholes.Time is determined from depth by Eq. (7).

Figure 4 .
Figure 4. Ground surface temperature history from the Manitoba boreholes at Flin Flon and Thompson (Pipe and Owl).The temperatures have been shifted with respect to the reference surface temperature of the site, T o , as shown in Table2.

Figure 6 .Figure 7 .
Figure 6.Ground surface temperature history for all the boreholes around Sudbury, Ontario (Victor Mine, Falconbridge, Lockerby, and Craig Mine).The temperatures have been shifted with respect to the reference surface temperature of the site, T o , as shown in Table2.

Figure 8 .
Figure 8. GST changes from simultaneous inversion with respect to the long-term temperature at 100 ka.The Sudbury GST changes include (black) and exclude (red) Craig Mine.
. Five additional profiles (Owl, near Thompson, Manitoba; Victor and Craig Mines, both near Sudbury, Ontario; Matagami and Val d'Or, Québec) were analyzed.The Val d'Or borehole was logged in 2010 to a depth of ∼ 1750 m.It is situated 15 km east of the mining camp of Val d'Or, Québec, in a flat forested area.The Matagami borehole is located near the mining camp of Matagami, some 300 km north of Val d'Or.

Table 1 .
Technical information concerning the boreholes used in this study.The Owl borehole, which was logged in 1999 and 2001, is located ∼5 km from the Birchtree Mine and ∼ 8 km south of the city of Thompson, Manitoba.The two other new boreholes, Craig Mine and Victor Mine, are located within the Sudbury structure, northeast of Lake Huron, in Ontario.The Craig Mine borehole, near the town of Levack, northwest of Sudbury, was logged in 2004.The deep mine was in operation when measurements were made and pumping activity was continuous to keep the deep mine galleries from flooding.The Victor Mine site was sampled in 2013, close to the community of Skead, northeast of Sudbury, Ontario.Victor Mine operated in 1959 and 1960, but exploration and engineering work is presently underway to prepare for reopening the mine at greater depth.

Table 2 .
Summary of GST history results where T o is the long-term surface temperature, Q o is the quasi-equilibrium heat flow, T min is the minimal temperature, T pgw is the maximum temperature attained during the postglacial warming, and t min and t pgw are the occurrence of the minimal temperature and maximum postglacial warming temperature.Parentheses indicate sites where the GST history is not reliable.

Table 3 .
Summary of GST history results for simultaneous inversions where t min and t pgw are the occurrence of the minimal temperature and maximal temperature associated with postglacial warming and T is the temperature range.