Subsurface temperature data is usually only accessible as point information with a very limited number of observations. To spatialize these isolated insights underground, we usually rely on interpolation methods. Unfortunately, these conventional tools are in many cases not suitable to be applied to areas with high local variability, like densely populated areas, and in addition are very vulnerable to uneven distributions of wells. Since thermal conditions of the surface and shallow subsurface are coupled, we can utilize this relationship to estimate shallow groundwater temperatures from satellite-derived land surface temperatures. Here, we propose an estimation approach that provides spatial groundwater temperature data and can be applied to natural, urban, and mixed environments. To achieve this, we combine land surface temperatures with anthropogenic and natural processes, such as downward heat transfer from buildings, insulation through snow coverage, and latent heat flux in the form of evapotranspiration. This is demonstrated for the city of Paris, where measurements from as early as 1977 reveal the existence of a substantial subsurface urban heat island (SUHI) with a maximum groundwater temperature anomaly of around 7 K. It is demonstrated that groundwater temperatures in Paris can be well predicted with a root mean squared error of below 1 K by means of satellite-derived land surface images. This combined approach is shown to improve existing estimation procedures that are focused either on rural or on urban conditions. While they do not detect local hotspots caused by small-scaled heat sources located underground (e.g., sewage systems and tunnels), the findings for the city of Paris for the estimation of large-scale thermal anomalies in the subsurface are promising. Thus, the new estimation procedure may also be suitable for other cities to obtain a more reliable insight into the spatial distribution of urban ground and groundwater temperatures.

1. Introduction

Natural in situ temperatures usually do not substantially vary at a depth of more than 10-20 m. Here, ground and groundwater temperatures are close to the annual mean values in the atmosphere, and commonly only a marginal attenuated influence of the coupled seasonal temperature variation in the atmosphere can be detected [1]. During recent years, attention has been growing to the thermal conditions beneath cities, which were revealed to be completely different when compared to the undisturbed rural surrounding [29]. However, the knowledge of the spatial and temporal variability of temperature in urban ground is crucial for proper planning of geothermal energy use. Elevated temperatures offer enhanced opportunities for geothermal heating applications, whereas warmer groundwater is less useful for cooling applications [1017]. Aside from its role for the geothermal potential, thermal anomalies can influence chemical transport in shallow urban groundwater that often serves as a freshwater resource [18, 19]. Finally, anthropogenic accumulation of heat threatens the stability of groundwater ecosystems [20, 21].

In most studied cases, temperature in urban ground and groundwater is elevated, which manifests in a large-scale heat carpet underneath a city. This so-called subsurface urban heat island (SUHI) is highly case specific, often with the highest temperatures beneath city centers and strong local spatial variations [7, 2228]. Ferguson and Woodbury [29] demonstrated that elevated groundwater temperatures below cities are mainly caused by heat losses from buildings. Accordingly, the investigation of SUHIs enables a determination of the initial period of urbanization in a certain region [30, 31]. In addition to building basements, groundwater flow, surface cover, and man-made climate change influence groundwater temperatures [3236]. Moreover, anthropogenic heat sources like subsurface infrastructure, power plants, landfills, or geothermal installations can have an impact on the spatial distribution of subsurface temperatures [14].

Various studies have reported SUHIs below large cities worldwide, for example, in Moscow [24], London [37], Cologne [7, 28, 38], Osaka [27, 39, 40], and Ankara [41]. Within these cities, temperatures are measured in a limited number of boreholes or groundwater wells and are commonly several degrees higher than in surrounding rural soil and pristine groundwater bodies. In the city of Paris, measurements from the 18th century onward are reported in a 28 m deep cellar at the city center [42, 43]. The long-term temperature increase of 0.1-0.7 K per decade is attributed to global warming and the developing atmospheric urban heat island (UHI) in Paris. Related studies such as by Pal et al. [44] primarily focus on the atmospheric and surface layer of Paris’ UHI. An analysis of surface temperatures shows that the intensity of Paris’ UHI, which is the difference between temperatures in the city and undisturbed conditions, can be up to 6 K [45]. Escourrou [46] confirms the presence of a strong atmospheric heat island in Paris and its effect on the local climate (e.g., influence on the rainfall regime and the occurrence of local winds). Cantat [47] underlines that the size and intensity of Paris’ UHI depend on various factors such as seasonal variability, climate, or weather conditions.

In this work, the focus is set on the shallow subsurface thermal regime in the metropolitan area of Paris. The objective is to examine the SUHI of Paris spatially based on repeated temperature measurements in wells located in and around the city. Due to the low well density within the studied urban area, especially within the city center, we investigate the auxiliary use of satellite-derived land surface temperatures to estimate spatial groundwater temperature distribution. In Method and Data, the measured groundwater temperature data is first described in detail. This information is contrasted to the satellite-based land surface temperatures for the city. Based on this, a new spatial estimation procedure for groundwater temperatures in urban, rural, and mixed environments is explained and applied to the case study of Paris.

2. Method and Data

2.1. Geographical Setting, Geology, Climate, and Hydrology

Paris is the capital of France and located in the northern part of the country in the region Île-de-France (Figure 1(a)). The focus of investigation here is the entire metropolitan area of Paris, which includes the city center district, the suburbs, and the surrounding rural areas. Paris has a circumference of nearly 55 km and also includes two huge woodlands: Vincennes in the east and Boulogne in the west. The average elevation of the city is 30 m a.s.l. The Seine River with a total length of 13 km inside of Paris crosses the city from northwest to southeast, dividing it into two sections [48].

In general, Paris’ city climate is moderate with relatively warm summers and mild winters. The mean annual air temperature is 10.8°C. Temperatures within Paris are slightly higher compared to the suburbs and the surrounding rural areas; the temperature difference typically amounts to 1-2 K. Elevated temperatures within the city center are caused by the increasing urbanization (e.g., higher building density) of the capital [47, 49].

Paris is eponymous for the Paris Basin (French: Bassin parisien). The Paris Basin is an elliptical geological basin filled with continental, epicontinental, and marine Mesozoic and Tertiary sediments. Postvariscan crustal extension caused subsidence processes during the Permian that formed the Paris Basin alongside other central European basins. It is present in the northeast of France, in the west of Belgium, and in the southeast of England. The center of the basin is located 80 km east of Paris. Aquiferous formations in the Paris Basin are mainly built up by sandstones, limestones, chalk, and carbonates. The aquifer units are intermitted by less permeable claystone and marl layers which results in a complex layered aquifer system [50].

2.2. Measured Groundwater Temperatures

There have been several groundwater measurement campaigns in the metropolitan area of Paris. Data from 1990 to 2015 were measured by local authorities in wells of the Paris metropolitan area and are available in the French national database on groundwater resources [51]. The well coordinates and measuring depths are accessible at InfoTerre [52]. Even though a rich dataset of 3018 measurements exists, the records are heterogeneous, often incomplete over time and information on the specific sampling methods is sparse. The majority of data were measured on site while pumping. Therefore, the given depths represent the mean depth of the screened well interval. During preprocessing, we excluded outliers exceeding a 6 K offset to the median for a single well location for further analysis. Moreover, measured values with an offset above 2.5 times the interquartile range (IQR) were expelled if the IQR was above 1 K. This was necessary to exclude outliers most probably caused by poor handling during sampling and typos. In addition, four outliers that are directly influenced by a waste disposal site north of Paris were eliminated. A total of 74 outliers (2.3% of the available data) were removed from the data pool. The remaining data covers 377 sampling wells in the outer and inner suburbs of Paris. The derived spatial coverage is nonuniform and there are no data points in the city center (Figure 1).

Diffre et al. [53] provide data that include measurements in the city center (Figure 1(b)) from February 1977. As the corresponding coordinates of the observation wells were not given in that report, they had to be extracted from the online database InfoTerre [52]. Unfortunately, none of the wells measured in 1977 are listed in the ADES [51] database (Figure 1). The 1977 GWTs focus on the inner suburbs and the city center, while the more recent ones reveal the thermal conditions in the surrounding.

2.3. Satellite-Based Estimation of Groundwater Temperatures

Recently, satellite data was suggested to investigate the thermal ground conditions in cities and SUHIs. Zhan et al. [54] estimated ground temperatures in Beijing, China, from MODIS land surface temperatures. The results were compared with recorded temperatures (measuring depths of 0.05 m, 0.40 m, and 3.20 m). A time delay between the maximum temperatures of the atmosphere and the ground was detected depending on the measuring depth. Hafner and Kidder [55] used Advanced Very High-Resolution Radiometer (AVHRR) satellite data to determine the albedo of the land surface as well as soil thermal and moisture properties around Atlanta, USA. These parameters were utilized for a numerical simulation of Atlanta’s atmospheric temperature dynamics. The spatial extension and intensity of atmospheric urban heat islands in 18 Asian megacities were analysed by Tran et al. [56] using satellite-derived land surface temperature () data. This study demonstrated that the spatial heterogeneity of the investigated heat islands mainly depends on vegetation cover and building density, whereas its total extent and intensity are related to the population size.

Benz et al. [38] analysed whether the spatial distribution of interpolated is linked to the spatial distribution of satellite-derived . In this context, the four German cities Berlin, Munich, Cologne, and Karlsruhe were analysed. In order to quantify the relation between and , the Spearman correlation coefficient was determined and correlations of up to 80% were found. Benz et al. [38] also propose a method for calculating groundwater temperatures with the help of mean annual , building density (), and basement temperature (). This reflects that accelerated ground heat flux in cities is caused not only by land surface changes and modified above-ground temperatures as detected from satellites but also by other heat sources such as heat release from basements, buildings, and underground networks.

Even for undisturbed rural areas, LSTs are generally colder than adjacent s. This offset is caused by a variety of processes in the soil transition zone and varies with local climate, biological activity, and the hydrogeological setting [57, 58]. Benz et al. [59] introduce an approach to calculate this offset on a global scale. This approach empirically relates evapotranspiration () and the number of days with snow cover () to the offset between and . For a global dataset, they could reduce the root mean squared error () of prediction by 0.5 K to 1.4 K and increase the squared Spearman correlation coefficient by 0.5 to 0.95 compared to a prediction relying solely on .

Due to the fact that we have a strong spatial variety of our measured in Paris (Figure 1), in the following method, we combine the individual approaches of Benz et al. [38] and Benz et al. [59] to estimate in the region of Paris, both for rural and urban areas at the same time (Figure 2). This combined approach has been suggested by Benz [60] for German cities, and it has the specific advantage of covering natural and urban processes within one estimate.

Estimated groundwater temperatures () are calculated in three steps:

The first step consists of estimating in rural areas. Here, the offset between and is dominated by the latent heat flux cooling down the surface in the form of evapotranspiration and by snow cover insulating the subsurface during the winter period. An empirical best-fit approximation for estimating rural groundwater temperatures () from satellite-derived data on a global scale is given by equation (1) [2]. The second step is estimating groundwater temperatures in urban environments () by using equation (2), which takes into account building density and basement temperatures [38]. The third step combines rural and urban estimations by replacing by in equation (2). The resulting equation (3) is intended to estimate groundwater temperatures in rural, urban, and mixed settings ().

Building densities were calculated at a resolution of from OpenStreetMap [61] building polygons. Basement temperatures were estimated to be following guidelines by the International Organization for Standardization (ISO) [62] and assuming similar conditions for Central Europe. This value range is also consistent with average basement temperatures obtained from borehole temperature profile inversion in the area of Zurich, Switzerland [30].

Satellite-derived data represent the decadal mean from 2005 to 2014 and include land surface temperatures, evapotranspiration, and the percentage of snow days. Data was retrieved and processed with Google Earth Engine according to the procedure that is explained in detail in Benz et al. [59]. The MODIS Aqua and Terra Land Surface Temperature and Snow Cover Daily Global 500 m products were retrieved from Google Earth Engine, courtesy of the NASA EOSDIS Land Processes Distributed Active Archive Center (LP DAAC), USGS/Earth Resources Observation and Science (EROS) Center, Sioux Falls, South Dakota [63, 64]. In contrast to Benz et al. [59], the updated MODIS data version 6 was used for determining values. Snow days were computed with MODIS version 5 due to different snow cover algorithms and product availability in the version 6 data, which could not be used with the empirical variables in equation (1) without biasing . For a detailed discussion of the given raster data, we refer to Benz et al. [59]. Estimated groundwater temperatures were computed by equations (1), (2), and (3) using ArcPy and ArcMap 10.5.1 [65]. Note that extracted raster data on point locations were bilinearly interpolated between the middle of adjacent raster cells. All raster data was exported at a resolution of approximately (Figure 3).

Land surface temperatures () range from 10.4°C in the outer suburbs to 14.9°C within the city center (Figure 3(a)). This dataset exhibits two central cold spots (12.7 and 12.5°C), caused by the two large forests in the eastern and western parts of the city center. The annual percentage of snow days () ranges from 0 to 7.3% and as expected s are more frequent in the outer suburbs than in the city center. Again, the two forests manifest as anomalies in the spatial distribution (Figure 3(b)). In contrast to and , the data on has a lower resolution of (approx. ). The city center is almost fully covered by one cell of low with a value of 10 mg m-2 s-1. Towards the rural areas, increases up to 17 mg m-2 s-1 (Figure 3(c)). Finally, building density () is calculated from building polygons of OpenStreetMap at a resolution of . is highest in the city center (up to 60%) and decreases towards the suburbs to nearly zero in areas of low population.

2.4. Misfit and Correlation

The estimate of is assessed by calculating the root mean squared error (), mean absolute errors (), and the mean error () between measured and estimated . and are used to evaluate the misfit of the predicted . Calculating is probably the most common way to evaluate the predictive capabilities of a model. However, it is more sensitive than to large prediction errors or outliers, because errors are squared before being averaged. is used to describe the error orientation, as positive and negative errors cancel each other out when calculating . It is thus useful to test if observed values are over- or underestimated by the prediction. In addition, a linear least square regression is performed to correlate measured and estimated s by using the Spearman method. Hereby, the obtained Spearman correlation coefficients () quantify the link between estimated and measured s, but are less useful as (mis)fit measures.

3. Results and Discussion

3.1. Spatial and Vertical Distribution of

For the first step, we characterize the spatial and vertical distribution of s based on the measured data (Figure 4). For this, the available early measurements from 1977, which focused on the city center, and averages from the years 1990 to 2015 for the greater Paris region were spatially combined.

In February 1977, the reported beneath the city center ranged between 14.1°C and 18.3°C. Towards the suburbs, the measured temperatures decreased yielding a temperature range from 12.2°C to 16.6°C. The measurement depths range from 17 m to 100 m below ground level. For the shallower measurements, a slight seasonal bias can be expected. For the 87 observation wells measured in February 1977, no repeated or more recent groundwater temperature measurements are available. From the period between 1990 and 2015, a total of 297 sampling points are located above a depth of 80 m below ground level showing a temperature range between 9.6°C and 22.8°C. Above 15 m below ground level, groundwater temperatures are strongly affected by seasonal influences. Here, they show a standard deviation of 2.82 K compared to a standard deviation of 1.23 K for measurements between 15 m and 80 m. For the latter interval, temperatures in the inner suburbs range between 11.8°C and 18.1°C and between 11.7°C and 17.2°C in the outer suburbs.

While no clear spatial trend can be deduced for the shown data, the lower limit of clearly declines towards the outskirts. Nevertheless, the lower limit of 14.1°C for the groundwater temperatures in the city center in 1977 is remarkable, being more than 3 K higher than the average air temperature of 10.8°C during that time. On average, they are also warmer than temperatures in the inner and outer suburbs. If we ignore any temporal and or vertical effects, the maximum temperature in the city center (measured in 1977) exceeds rural groundwater temperatures (measured after 1990) by 5 to 6 K. This offset is typical for SUHIs in Western Europe and is on the upper end of observed intensities in previously studied German cities [7].

The vertical resolution of the measured data is depicted in Figure 4. Most sampling points were obtained in shallow depths of a few tens of meters, with some measurements also taken in depths of more than 100 m. However, no clear evidence of the geothermal gradient can be observed. Still, for the comparison between satellite-derived and the measured , the lower depth limit for this analysis is set to 80 m in order to zoom in on the regime influenced by land surface effects. Additionally, the shallowest 15 m of the subsurface are disregarded to avoid any influence of seasonal atmospheric temperature variation. The usual seasonal temperature variation in a humid climate is below 1 K for depths of more than 10 m [1].

In the next step, equation (3) is applied to estimate groundwater temperatures using the combined method, . This resulting map in Figure 5(b) represents at shallow depths with a temperature range of 11.3°C to 16.6°C. Consistent with the available groundwater temperature measurements from different points in time, the city center is characterized by warm temperatures, which are decreasing towards the suburbs.

3.2. Comparison of Measured and Estimated

To assess the quality of predicted groundwater temperatures from satellite imagery, , they are compared with well measurements in the field. It is important to note the time difference: data in the city center was only observed in 1977, about three decades earlier than the satellite data considered for calculating . However, as s are expected to have increased during that time [30, 66, 67], these measurements are included in the following comparison viewing them as minimum values. In order to inspect the role of temporal variability and trends, we separately compare measured s from (i) 1977, (ii) 1990-2015, and (iii) 2005-2015 with satellite-derived estimates. The latter subset is extracted as it covers roughly the same time period as the satellite-based product. Since satellite data covers the time period from 2005/01 till 2014/12 and groundwater temperatures are recorded till 2015/11, the data contains a negligibly small amount of 85 unsynchronized single measurements (4.4% of the single measurements between 2005 and 2015).

Figure 6 shows the measured s versus and along with misfit and Spearman correlation indices. The city center wells measured in 1977 (Figure 6(a)) yield misfits (, ) between and measured s around 0.5 K higher than the error observed between and . Additionally, is higher (Figures 6(a) and 6(d)) when comparing measured groundwater data with (1.27 K) than with just (0.09 K). This means that subsurface temperatures are overestimated by when using recent values of and old well temperatures in equation (3). This is because is expected to have increased by around 1 K from the 1970s to the period between 2005 and 2015 due to climate change and atmospheric urban heating, which is not surprising. At the same time, s in 1977 were as high as 18.3°C, which is warmer than the assumed basement temperature of 17.5°C in equations (2) and (3). Therefore, other heat sources such as heat release from subsurface infrastructures may play a significant role but cannot be reproduced by this method. This is supported by the low sensitivity of the measured to as revealed by Figure 6(a).

The picture is different when analysing the recent datasets from 1990 to 2015 and 2005 to 2015: when using instead of , s decrease by 0.53 K and 0.49 K and decreases by 0.51 K and 0.46 K for the time periods of 1990 to 2015 and 2005 to 2015, respectively (Figures 6(b), 6(c), 6(e), and 6(f)). and relative to are below 1 K for both time periods, with both misfit values of the reference time period (2005 to 2014) being slightly lower compared to the entire dataset from 1990 to 2015. We can therefore conclude that the estimation approach is most applicable for overlapping time periods between subsurface and satellite measurements.

Overall, the combined approach as given in equation (3) delivered groundwater temperature estimates with a below 1 K. However, even for this time period, extreme temperature anomalies cannot be predicted. It is very likely that these wells are impacted by local heat sources other than buildings. When applying this method to the German cities of Berlin, Cologne, and Karlsruhe, Benz [60] found similar results with s between and of 0.96, 0.86, and 0.82 K, respectively.

3.3. Impact of Urban and Rural Estimation Approach on s

To evaluate the suitability of the different estimation procedures (equations (1), (2), and (3) and only), we use the correlation coefficient (), , , and . Figure 7 shows the correlation to the measured s of the different equations for calculating , , and and the sole use of .

While and are lowest for the combined method, indicating an improvement of the estimation accuracy, correlation is highest for . According to these results, this approach which estimates naturally occurring offsets also has an impact on the predicted groundwater temperatures which is higher than that of the urban approach (). This is a consequence of the location of the observation points that are solely in the inner and outer suburbs of Paris and have values below 0.3 (30%). In the city center, values are more than doubled, which also doubles the amount of positive temperature correction (equation (2)). At the same time, has a lower impact on the city center because both and are significantly lower (Figure 3).

All predictions underestimate measured s with a minimum of -0.23 K for . However, we could improve the negative offset by 0.9 K compared to the misfit between measured and . In contrast, Benz [60] found that the combined method overestimates temperatures in the city centers of Karlsruhe, Cologne, and Berlin. This overestimation can in part be attributed to the insulation of basements, which is not considered in this estimation. Overall, we expect to work best in the city centers of larger urban areas, due to higher and lower and values. For the given groundwater data in the outskirts of Paris, is only slightly favorable over the use of .

Although the comparison between the results of the individual equations correlates well with measured groundwater, the used approach also has some limitations. Particularly, local thermal disturbances caused by underground anthropogenic heat sources (e.g., geothermal devices, sewer networks, and underground constructions) cannot be detected from above-ground data or their thermal footprints are too small to be resolved by the resolution of of . Moreover, the resolution of the dataset is very coarse (about 28 km) and cannot detect small-scale or even city-scale features. Applying the combined method to Karlsruhe and Cologne in Germany, the dataset did not display the reduction of commonly associated with urban areas, which added to an overestimation of s in those areas [60]. Unfortunately, no recent s within the Paris city center are available here, and these observations can thus not be tested for such a large-scale city. However, future measurement campaigns could help to eliminate this limitation.

Another limitation of the method used is that as primary input data already indicates a heat anomaly in the city center, and the data used for and are also affected by urbanization. Accordingly, it cannot be differentiated between natural and anthropogenic thermal contributions by the use of this method, whereas the processes itself can fairly easily be attributed to one or the other.

4. Conclusions

For the first time, the large-scale subsurface urban heat island (SUHI) in the groundwater of Paris is spatially investigated. The minimum and maximum groundwater temperatures (s) are 11.4°C and 17.2°C measured between 2005 and 2015 for the shallow subsurface at depths between 15 m and 80 m. These groundwater temperatures represent the inner and outer suburbs of Paris. They should be interpreted as regional trends and are not meant to represent local temperature disturbances caused by anthropogenic heat sources which may be of a much greater magnitude. The minimum and maximum s observed in 1977 below the city center ranges between 14.1°C and 18.3°C. Despite the temporal gap of 28 years between the two measurement campaigns, this indicates a characteristic difference between the center and the suburbs. Based on the examined data, the maximum anomaly reaches 6.9 K in Paris.

In the present study, an approach is applied that estimates groundwater temperatures () from satellite-derived data and building footprints combining existing methods for estimating rural () and urban () groundwater temperatures. Estimates are validated with temperatures measured in the wells of Paris. Predicted s have a very good fit () and only slightly underestimate measured with a mean error of -0.23 K. In addition to that, a coherent spatial estimation of the SUHI of Paris is obtained. Even so, measured extreme s that cannot be reproduced as subsurface heat sources other than in basements are not considered. In an ideal case, the estimation techniques displayed here would supplement existing temperature measurements in wells, which would allow the use of this method as a tool to spatially connect point information obtained from wells.

In future studies, the presented procedure needs to be further validated by estimating s in other rural and urban environments. This is ideally accomplished in regions, where data with a very dense spatiotemporal resolution is available. When applying this approach to other cities, it should be kept in mind that the anthropogenic fingerprint varies among different regions, countries, and climates. The role of this variability will be scrutinized by a comparison of different case studies in the future.

Data Availability

Groundwater temperature data is available in the French national database on groundwater resources [51].

Conflicts of Interest

The authors declare that there are no conflicts of interest.


We thank Enzio Schoffer and Jakob Michael for language editing. This work is funded by the German Research Foundation (grant number BA2850/3-1).