Laurentide Ice Sheet basal temperatures at the Last Glacial Cycle as inferred from borehole data

Laurentide Ice Sheet basal temperatures at the Last Glacial Cycle as inferred from borehole data Carolyne Pickler1, Hugo Beltrami2,3, and Jean-Claude Mareschal1 1GEOTOP, Centre de Recherche en Géochimie et en Géodynamique, Université du Québec à Montréal, Canada 2Climate & Atmospheric Sciences Institute and Department of Earth Sciences, St. Francis Xavier University, Antigonish, Nova Scotia, Canada 3Centre ESCER pour l’étude et la simulation du climat à l’échelle régionale, Université du Québec à Montréal, Canada Correspondence to: Hugo Beltrami(hugo@stfx.ca)


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 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 two 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 Heimbach, 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), ∼ 120000-12000 yrs 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 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 (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 and 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 temperature 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 or 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 time scale of 10 to 100 kyr.Oil exploration wells are usually a few km 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 time scale 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 yrs BP, could not have been more than 5 K colder than present.Other Canadian deep boreholes temperature profiles have since been measured revealing regional differences in temperatures during the LGM.From a deep borehole in Sept-Iles, Québec, Mareschal et al. (1999b) found surface temperatures to be approximately 10 K colder than present.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 from 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 2,400m deep well near Fort McMurray, Alberta, penetrating in the basement.They interpreted the variations in heat flux as due to a 9.6 K increase in surface temparature at 13ka.Majorowicz and Šafanda (2015) have revised this estimate to include the effect of the variations of 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 o 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 (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 o 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"/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 directly used variations in the long-term surface temperature with depth 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 equation 5, the time it takes for the surface temperature variation to propagate into the ground and the depth of the perturbations are related by: with κ 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 of 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 independent of the model assumptions in an inversion procedure.

Inversion
In order to obtain a more robust estimate of the GST histories, we have inverted the temperaturedepth 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 temperature-depth profile (equation 6).For each depth where temperature is measured, equation 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 (e.g., Vasseur et al., 1983;Shen and Beck, 1983;Nielsen and Beck, 1989;Wang, 1992;Shen and Beck, 1991;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 (Beltrami et al., 1997), 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 thirteen 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 that extended over most of Canada during the LGC (Figure 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 seismometers for monitoring nuclear tests.Detailed description of the measurement techniques as well as the relevant geological information can be obtained from the 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 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.005K with an overall accuracy estimated to be of the order 0.02K Thermal conductivity is measured on core samples.Samples are collected in every lithology with an average of one sample per 80m.The thermal conductivity is measured with a divided bar apparatus on five cylinders of the core with thickness varying between 0.2 and 1.0cm.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 Victor Mine operated in 1959 and 1960 but exploration and engineering work is presently underway to prepare for reopening the mine at greater depth.
We found systematic variations of thermal conductivity at Flin Flon, Thompson (Owl), and Matagami, and we corrected 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 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 of heat production with depth were not present at these sites, they were not included in the inversion.In 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) is thus a diagnostic of departure from 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 (Figure 2) exhibit clear departures from 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 high frequency 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 Val d'Or, 85 degrees at the collar, was only 20 degrees near the end of the hole.Most of the profiles show a very pronounced increase of 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.
4 Analysis and Results

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 (Figure 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 reveal 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 that does not show a clear trend.Variations of 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 four 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 years BP was inverted from the temperature-depth 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 (Figures 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 20ka.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 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 Pipe, have not recorded a very strong signal do not exhibit common regional trends are observed (Figure 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 3K.The weak signal recorded may be in part because the present ground surface temperature is very low; it is close to 0 o C in Thompson where intermittent permafrost is found.At Flin Flon, where the present ground temperature is near 3 o C, we found that ground surface temperature variations were small, confirming previous studies by Sass et al. (1971).
The surface temperature was minimum around the LGM and was near the melting point of ice (-0.3 o 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 has doubled (∼5 K), with the minimum temperature (-2.4 o 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 high heat flux 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.For these reasons, we have little confidence in the robustness of the Owl GST history reconstruction.In western Ontario (Balmertown, Manitouwadge 0610, and 0611), results at all sites show minimum temperatures around the LGM, with a very weak minimum at Balmertown (Figure 5).However, the amplitude of the temperature change is much larger at Manitouwadge 0611 (∼10 K) than 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 accounted for.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-7 K.This difference is unexpected as these sites are all within the Sudbury craton 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 o C and 3.0 o C, are also the highest of the study.These are also the only two sites with non-negligible heat production, 3µW m −3 and 0.9µW m −3 respectively.
The corrections for heat production produce an increase of 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 those with negligible heat production.
The GST histories from the three boreholes in Québec, Matagami, Val d'Or and Sept Iles, display regional differences (Figure 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 o C at Sept Iles.
At all sites, excluding Pipe mine, 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 at 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 K 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 1-2ky around 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.This would require obtaining proxy data from sites close to Flin Flon and Thompson.This also 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 simultaneously the boreholes of Thompson (Owl and Pipe), Manitouwadge (0610 and 0611) and Sudbury (Craig Mine, Falconbridge, Lockerby, and Victor Mine) to observe regional trends and to improve the signal to noise ratio (Figure 8).The simultaneous inversion results are summarized in Table 3.For simultaneous inversion, the temperature-depth profiles were truncated to ensure a common depth for all the boreholes.This insures 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.Altough Owl and Manitouwadge 0611 have questionable individual inversions, we have still performed simultaneaous inversions for the regions of Thompson and Manitouwadge in an attempt to decrease the signal to 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 slightly the questionable signal recorded at Owl.The amplitude difference for the Manitouwadge GST history is ∼8 K, 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 11K, 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 o 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 (e.g., Dahl-Jensen et al., 1998;Pattyn, 2010;Schneider et al., 2006).The 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, ∼54mW m −2 , because crustal heat production is higher than average (Perry et al., 2009).The lowest basal temperature is recorded at Sept Iles, where the heat flux is the lowest of the studied regions, ∼34mW m −2 .These correlations suggest a link between heat flux and basal temperatures.This is further supported by modelling 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 Iles basal temperature has been also been linked to its proximity to the edge of the ice sheet and area of thinner ice (Rolandone et al., 2003b).Our study is consistent with a possible link between basal temperature and heat flux along with ice dynamics in the Laurentide Ice Sheet during the LGC but further modelling 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 increase 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 associated melt persisted prior to and throughout the LGM over 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 collapse, it demonstrates that they are a key parameter in ice sheet evolution, 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 alone this 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 ca 10ka.The corresponding temperatures at the base of of the Laurentide Ice Sheet range from −1.4 to 3.0 o 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 parameters as basal temperature.Despite the suggestion that melting took place at its base, the Laurentide ice sheet persisted throughout the LGM 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.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.

Figure 1 .
Figure 1.Map of central and eastern Canada and adjoining US 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 .Figure 3 .
Figure 2. Heat flux variation as a function of depth.Heat flux is calculated as the product of thermal conductivity by the temperature gradient calculated over 3 points.The Flin Flon, Owl and Matagami profiles have been corrected to account for thermal conductivity variations with depth as shown in Table1.

Figure 4 .Figure 5 .Figure 6 .
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, To, as shown in Table2.

Figure 7 .Figure 8 .
Figure 7. Ground Surface Temperature History for the boreholes in Quebec, Matagami, Val d'Or and Sept-Iles.The temperatures have been shifted with respect to the reference surface temperature of the site, To, as shown in Table2.

Table 1 .
Technical information concerning the boreholes used in this study.

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).The temperature profile at this site may be disturbed by horizontal contrats in thermal conductivity, ‡ The temperature profile in the lowermost part of the hole may be affected by subvertical layering and thermal conductivity contrasts, § The temperature profile may be affected by waterflow caused by pumping in the nearby mine.