Recent Climate Change in Western Black Sea Region of Turkey by Paleoclimatic Reconstruction of Borehole Temperatures

The impact of recent climate change varies around the world, and understanding the regional variations is important for making accurate deductions and generating future climate predictions and mitigation plans. Borehole temperature reconstruction is one of the common methods to determine the ground surface temperature (GST) history of a region. In this study, high-precision (< 0.01 K) borehole temperature-depth data from five different locations in Western Black Sea Region of Turkey were used for reconstruction of the GST changes for the last century. Four of the measurement sites are located in rural areas in the inland section and free from the urban heat island interference. The reconstructions reveal that GSTs have dropped by an average of 0.1 °C between 1900 and 1970 and have risen by an average of 0.99 °C between 1970 and 2010 in the study area. Results from inland sites show average cooling of 0.24 °C between 1900 and 1975 and average warming of 0.94 °C until 2010. The results in the coastal area show no cooling period. The rapid warming trend in 1990–2010 revealed by GST reconstructions indicates high sensitivity of the region to present-day global warming.


Introduction
Climate change is a global phenomenon that impacts many aspects of the environment.Climate change increases drought risk (Diffenbaugh et al., 2015), affects lakes (Wang et al., 2018), influences the extinction rates of species (Roma ´n-Palacios & Wiens, 2020), can create conditions for disease outbreaks (Bryson et al., 2020) and increases wildfire risks (Jolly et al., 2015).The effects of climate change vary in different parts of the world.Therefore, understanding of local conditions is crucial for mitigation plans.Turkey has diverse geographical features and is influenced by different teleconnection (connection of meteorological events that occur a long distance apart) patterns.These characteristics lead to different types of extreme weather events occurring around the country, and climate change is increasing the frequency of those extreme weather events (Baltaci & Arslan, 2022).
Surface air temperature (SAT) around the globe has been rising for the last century as a result of climate change (Alexander et al., 2006;Brunet et al., 2007;Hansen et al., 2006;Ji et al., 2014;Ren et al., 2017).Rock formations at the subsurface level are affected by temperature changes occurring at the surface of the earth and hold information on surface temperature history (Pollack, 1993).This information is particularly important for understanding the climatic conditions of the Earth prior to instrumental temperature recordings.
Temperature variation at the Earth's surface can occur in different time periods.Night and day cycle creates the shortest temperature oscillations whereas orbital movement of the earth dictates seasonal and annual temperature trends.Surface temperature variations penetrate to subsurface layers as thermal wave and strength of this effect decrease exponentially with depth as dictated by the theory of heat conduction (Birch, 1948;Pollack, 1993).Thus, ground surface temperature (GST) variations influenced by long-term climate change create temperature anomaly zones that GST history information can be extracted from by using temperature-depth profile analysis (Eppelbaum et al., 2006).
The physical conditions of the boreholes are crucial for data quality (Eppelbaum et al., 2014).Regional groundwater activity is one of the most important subjects that limits the availability of boreholes for climate reconstructions.Boreholes must have a uniform conductive heat flow originating from the depths.Vertical groundwater flows within boreholes (Erkan, 2015) also disturb this heat flow locally and prevent obtaining useful temperature-depth data from the boreholes for climate reconstructions.Canopy and humanmade structures around the borehole are also important.Open space around the borehole with minimal vegetation is the ideal condition for climate studies (Pollack, 1993).
Paleoclimatic reconstruction of borehole temperatures has been applied in many different parts of the world to determine regional scale climate change.For example, Akkiriaju and Roy (2011) collected temperature-depth data from a 210-m borehole and showed 0.5 ± 0.1 °C increase in GST over the past 92 ± 7 years in Hyderabad District of South India.Pickler et al. (2018) found a cooling period between 1850 and 1980 in their study in Northern Central Chile.According to the study, the cooling period was followed by 1.9 °C increase in the last 20-40 years.In their study in northern Quebec, Canada, Chouinard et al. (2007) showed 1.4 °C increase between the mid-eighteenth century and 1940s, followed by 0.4 °C cooling period that lasted for 40-50 years.Following that, results showed 1.7 °C increase in GST in 15 years.Correia and Safanda (2001) found in their study that about 1 °C increase in GST occurred between 1950 and the mid-1990s, with an increase over the last 15 years in Evora Town, located in Southern Portugal.In Turkey, Erkan et al. (2019) utilized borehole paleoclimate method to determine GST values for Western and Central Anatolia.Results showed 0.33 and 0.42 °C cooling between 1900Results showed 0.33 and 0.42 °C cooling between and 1990Results showed 0.33 and 0.42 °C cooling between and warming of 1.47 and 1.66 °C between 1990Results showed 0.33 and 0.42 °C cooling between and 2010 for the western and central Anatolian sites, respectively.
In this study, reconstruction of recent GST history was done in six provinces of Western Black Sea Region of Turkey using borehole temperature data.First, boreholes, which are available for measurements and have adequate depth, were identified.Then, temperature-depth (T-D) data were collected from these boreholes in the field studies.Collected T-D profiles were further analyzed to check whether they were suitable for GST history reconstruction.Finally, GST histories were created by utilizing an inverse modeling algorithm.

Regional Setting
Western Black Sea Region of Turkey includes Kastamonu, Zonguldak, Karabuk, Bartin, Bolu and Duzce provinces (Fig. 1).In the region, mountains are oriented west to east, running parallel to the Black Sea.This situation creates an orographic effect that leads to microclimatic divisions in the area.The northern sides of the mountains have warm and humid conditions, whereas the southern sides have warm but step-like climatic conditions (Go ¨ktu ¨rk et al., 2011).Combination of cool northern air masses with the Black Sea acts as a moisture source and dictates the regional climate.The climatic difference between southern and northern sides also affects ecological structure in such a way that the moist summer conditions in the north allow dense forestation and prevent wildfires.Climatic and ecological diversity of the region makes it an important agriculture and forestry zone for Turkey.In our study, four of the boreholes (Akoren, Devrekani, Dizdarli and Hasirli) are located in the southern (inland) parts of region, and one borehole (Bartin) is located in the northern (coastal) part.Borehole locations in the southern part have average elevation of 1052 m (varying between 865 and 1185 m), and the location of the borehole in the northern part has an elevation of 42 m.

Data Collection
GST history reconstructions require high-quality data obtained from boreholes.Abandoned, low-yield boreholes without any mechanical attachments (pumps, etc.) are readily available for temperaturedepth measurements.However, only a few boreholes can provide high-quality data although they fit the description.The most common problem is the groundwater flow.Second, even the groundwater activity is negligible; intra-borehole vertical flow due to differences in pressure heads can disturb background thermal conditions along the borehole (Erkan, 2015).Many measurement sites had to be eliminated because of these effects.According to the class definitions given in Erkan (2015), boreholes were classified into four types, denoted as A, B, C and D, to qualify their data availability for GST modeling.As shown in Fig. 1, boreholes in Kastamonu Province are located in Dizdarli, Devrekani and Hasirli Districts.These districts are located in inland parts of Kastamonu Province.Each location has its distinct topographic characteristics.Dizdarli is located in the western region of Kastamonu.The studied borehole is located near a village on a hill.Devrekani and Hasirli Districts are located about 20 km north of Kastamonu city center, south of Kure Mountains.The borehole in Hasirli is located near a village, with a vast open space around it with very light vegetation.The borehole in Devrekani is located on a small hill with no vegetation, and it is 500 m away from the farmyards.The borehole in Bartin Province is located about 8 km south of the central district in an area that is mostly grassland.It is inside an orchard that was planted in 2017.Finally, the borehole in Bolu Province is located near Akoren village.The village is located on an open area.There are no lakes near any borehole.Most of the borehole locations were obtained from Provincial Special Administrations (local administrations that are in charge of development and support of rural settlements), and few of them were obtained by contacting private owners through drilling contractors.Our investigations show that areas around the boreholes were not part of a forest, and no dried lakes had existed nearby in the last century.No farming activity was recorded near the boreholes in the inland section (Akoren, Devrekani, Dizdarli, Hasirli).Area around the Bartin borehole was used for low density unirrigated (rainfed) farming until 2017 and converted to orchard plantation later.Our study covers the 1900-2010 period, which predates this land use change.The properties of the boreholes are shown in Table 1.
The temperature logging device includes a thermistor sensor that can obtain high-resolution temperature-depth data (Erkan et al., 2017).Detection limit of the temperature logging device is \ 0.01 °K; however, the practical temperature resolution in the field was determined to be 0.02 °K.This value is used as the error bound for the T-D data in the model input.Signal strength depth that can be measured is 300 m.Clow et al. (1996) showed that the precision of the temperature instrument has a significant effect on the resolution of GST reconstructions.Our sensitivity test for various incremental precisions indicated that our measurement precision (20 mK) prohibits GST reconstructions beyond 110 ± 10 years.During the measurement process, temperature logging probe was released into the boreholes, and temperature and depth data were recorded with 2.5-5 m logging intervals.When the probe is in water, about 20 s was employed to reach equilibrium at each logging depth.On the other hand, measurements above static water level needed to be approached differently as the sensitive temperature logging probe took longer to reach equilibrium in air.Therefore, a time-lapse measurement procedure was necessary for measurements above the water level.This type of measurement was done in Hasirli and Devrekani boreholes for 30-60 and 30-80 m, respectively.During these measurements, temperatures at each 30-s interval were recorded for 5 min at each depth (Fig. 2).
For extrapolation of the time-lapse temperatures to the final equilibrium temperature, the following equation was used (Costain, 1970): where t is the time of recording, T is the temperature at time t, T f is the final temperature, T 0 is the initial temperature, and s is the time constant of the probe (taken as 3.5 min for our probe).Equilibrium temperature value (T f ) was then obtained by applying generalized least squares (GLS) method (Chapra and Canale, 2011) between the data and the equation above.For the detection and elimination of the problematic measurements, we applied a two-step procedure.First, we accepted the results at a certain depth satisfying |T f -T 0 |B 0.02 as this upper bound is within the practical limit of temperature measurements in the field.For the measurements in the remaining depths, we enforced a standard error (SE) of the GLS fit to be B 5% in accordance with the measurement precision of our logging instrument.
Figure 2 shows example temperature measurements from Devrekani that were taken above the water table level.For 30 m (Fig. 2a), |T f -T 0 | value Figure 2 Time-lapse temperature measurements and data fitting at various depths above water level.T obs is the recorded temperature, and T fit is the temperature calculated by extrapolation for the given time is 0.04, which is above the threshold, but SE is found as 4%, which is below the 5% limit.Therefore, T f value of 10.79 °C calculated by our algorithm is acceptable.As another example, Fig. 2-b shows temperature measurements and curve fit data from 45-m temperature-depth log data.|T f -T 0 | value here is found to be 0.02, which is below the threshold so calculated value of 10.88 °C is acceptable for T f .Temperature-depth data used as the input for GST modeling are shown in Fig. 3. Except Deverakani and Hasırlı, all data were recorded below the water table level.In the subsurface, short-term annual variations in surface temperatures are generally effective down to depths of * 20 m (Erkan et al., 2019;Harris & Chapman, 1998).Our repeated T-D measurements in Bartin, which are 1 year apart, also showed that seasonal variations can be effective down to around 25 m.Therefore, temperature logs influenced by short-term annual temperature variations were omitted from GST model inputs.A unit temperature change that occurs during the last century can influence the temperature-depth profile of a borehole in up to 150 m (Pollack, 1993;Safanda et al., 1997).However, such strong climate signals have not been observed in the world except at high latitudes.Observations at the mid-latitudes show that the climatic component of the T-D data occur at shallower depths (Erkan et al., 2019;Gosnold et al., 1997).Some of the boreholes in our study are relatively shallow becuase of characteristics of boreholes in the region; however, linearity of T-D profiles at the bottom of these boreholes indicates steady state conditions at the lower portions.
Variations of the topography can have significant static effect on the background geothermal gradients and hence the background heat flow.The boreholes used for GST modeling are located in areas with insignificant topography except Dizdarlı.When generating the transient data (data after the static conditions have been removed), the effect of topography on the background conditions is automatically eliminated.Hence, no topographic correction was applied for Dizdarli site.

GST Reconstructions
Ground surface temperature reconstructions were performed by utilizing the functional space inversion (FSI) algorithm (Shen & Beck, 1991).FSI algorithm is based on functional space analysis instead of conventional discrete formulation; the solution is provided by using the iterative gradient method.As described by Shen et al. (1995) FSI algorithm selects the model that minimizes the total misfit function, which is: where d is the calculated borehole temperatures, d 0 is the a priori borehole temperatures, C d is the covariance matrix that defines uncertainties in the a priori temperatures, m is the calculated model parameters, m 0 is the a priori model parameters, and C m is the covariance matrix that defines uncertainties in the a priori parameters.FSI model requires temperaturedepth profile, thermal properties of subsurface layers, background heat flow and undisturbed surface temperature as input.In our study, temperature-depth profiles were determined by instrumental loggings as   1).Background heat flow (qb) was calculated by taking the deepest section of the temperature-depth log where temperature-depth profile shows a linear trend and can be defined as the undisturbed zone.According to Fourier's law: where d 1 (lowest) and d 2 (deepest) are two depth points in linear portions of the temperature-depth curve, T 1 and T 2 are temperatures measured at depths d 1 and d 2 and, finally, k is the thermal conductivity (W/m/K) of the medium.Thermal conductivity values for Hasirli and Devrekani were determined by instrumental measurements of collected rock samples (Table 1).Our thermal conductivity instrument (QTM-500) can only measure consolidated (bulk) rock samples.When the samples were in the form of drill cuts prohibiting direct measurements, we determined the thermal conductivity values based on literature values after specifying its lithological properties (Balkan et al., 2017;Gul & Maqsood, 2006).In Bartin, Akoren and Dizdarli boreholes, thermal conductivity values were taken as range of values because of unavailability of thermal conductivity measurements.Undisturbed surface temperatures (''pre-observational mean GST'', Harris & Chapman, 1997) were determined by extending background T-D profiles (see Fig. 3) to the surface.Thermal diffusivity, which dictates the distance of temperature signal, is calculated as: where j is thermal diffusivity (m 2 /s), k is the thermal conductivity (W/m/K) and c is the volumetric specific heat (J/m 3 /K).Volumetric specific heat was selected as 2.30 MJ/m 3 /K.FSI model processes these data as inputs and produces a posteriori temperature-depth and background heat flow values.These products are used to test the model accuracy by comparing them with the observed values.If a priori and a posteriori profiles are in agreement, reconstructed GST values are accepted as the confirmed model.Standard deviation (SD) values for a priori thermal conductivity and temperature data are also required for GST solutions.SD value for thermal conductivity was selected as 0.5 W/m/K, and SD value for borehole temperature data was selected as 20 mK.The resulting propagated error is provided by the FSI algorithm along with the GST reconstructions.

Results
Reconstructed GST histories are shown in Fig. 4. Due to elimination of near-surface seasonal variations, very recent changes in the GST histories cannot be modeled.Therefore, GST models cover the 1900-2010 period.Boreholes investigated in Kastamonu (Devrekani, Dizdarli and Hasirli) Province show a cooling period starting from 1900, followed by a warming period.Durations of cooling and warming periods differ in each location.Devrekani district (Fig. 4a) shows 0.4 °C cooling between 1900 and 1960, followed by 0.96 °C warming until early 2000s.GST model shows a slight 0.15 °C cooling in the late 2000s but net change from 1900 to 2010 still adds up to 0.56 °C warming.In Dizdarli (Fig. 4b), a longer cooling period is found as 0.27 °C until 1990 for k = 2 W/m/K and 0.4 °C until the late 1990s for k = 3 W/m/K.After cooling period, 0.64 and 0.9 °C warming is found until 2010 for k = 2 W/m/K and k = 3 W/m/K conditions, respectively.Hasirli district (Fig. 4c) shows 0.83 °C cooling period between 1900 and early 1990s.Following the cooling period, 1.38 °C warming is found until 2010.Borehole investigated in Bartin Province (Fig. 4d) shows a different warming trend.Here, GST starts rising with the beginning of the twentieth century, and 1.4 and 1.67 °C warming is found for k = 3 and k = 2 W/m/ K conditions, respectively.The result in Akoren (Fig. 4e) shows similarities with Hasirli and Dizdarli in terms of cooling prior to warming.In Akoren, 0.1 °C cooling until the 1940s and 0.3 °C cooling until the 1970s were found for k = 2 W/m/K and k = 3 W/m/K conditions, respectively.Following this To compare SAT and GST results, SAT archive data were obtained from Turkish State Meteorological Service (TSMS).Comparison of GST model results and SAT history from Kastamonu, Bartin and Bolu TSMS stations is shown in Fig. 6.SAT values are expressed as anomalies regarding the mean of the 1961-1990 period, and GST values are given as  from the 1990s.Dizdarli and Hasirli GST histories show similar trends; however, late warming starts earlier in Devrekani, and a slight cooling period is observed near 2010 (Fig. 6a).In Bartin, SAT recordings show a sharp increase starting from the late 1980s, and accordingly, GST histories show that 50% of the total warming occurred in that period (Fig. 6b).Similarly, SAT recordings show an increase in the late 1980s in Bolu and warming in Akoren GST, especially at k = 3 W/m/K conditions, which also increases in that period (Fig. 6c).

Discussion
In our study, averaged GST history shows a slight (B 0.1 °C) cooling until 1970 and 0.99 warming between 1970 and 2010, resulting in a net 0.91 warming between 1900 and 2010 (Fig. 7a).According to our study, averaged GST histories of Dizdarli, Devrekani, Hasirli and Akoren (boreholes located in inland sections of study area) show two phases as 0.24 °C cooling between 1900-1975 and 0.94 °C warming afterwards (Fig. 7b).In a similar study in Turkey, Erkan et al. (2019) focused on Western and Central Anatolia.According to this study, GST history shows three phases.First, a small warming occurs between 1900 and the 1950s, which is \ 0.25 °C.Following that, cooling of an average * 0.5 °C occurs until the 1990s.Finally, * 1.2 °C warming occurs from 1990 to the 2010s.Although GST histories of Western Black Sea Region show a different warming/cooling trend than Western and Central Anatolian Regions, both similarly show significant warming in the late twentieth century and early twenty-first century.
Results include two different GST history profiles for Akoren, Bartin and Dizdarli where thermal conductivity values were selected as k = 2 W/m/K and k = 3 W/m/K.This difference is caused by thermal diffusivity values as it determines distance of temperature signal and is dictated by thermal conductivity value.Although the results have nuances in cooling and warming magnitudes and in the time points where transition between cooling and warming occur, results still clearly indicate the rapid warming since the late twentieth century.
All boreholes located in the inland section of the study area show net warming.There are also variations in cooling/warming periods and magnitudes.Notably, Devrekani borehole shows different GST trends and a late cooling period.Regional variations in climatic trends also occur in SAT studies.Tu ¨rkes ȩt al. (2002) showed that annual minimum SAT values showed cooling trends between 1940 and 1970 in Western Black Sea region.This cooling trend is associated with increased seasonal northerly circulations in that period.Study also showed that SAT records in Black Sea region show the highest amount of diversity in terms of cooling and warming periods among all regions.Tu ¨rkes ¸et al. (2016) showed that three different SAT trend clusters (humid Black Sea cluster, semi-humid Western Anatolia cluster, semiarid Central Anatolia cluster) exist in our study site.Areas of these clusters change between 1950-1980 and 1981-2010 periods.Late cooling period observed in Devrekani in our study can be attributed to these micro-climate variations.These climate events are most likely due to topographic properties of Western Black Sea Region and local effects of the Black Sea (Go ¨ktu ¨rk et al., 2011).Further studies in the inland section of Black Sea Region may bring more information.
Mean GST model history of Bartin shows a different warming trend than other sites (Fig. 7c).It shows slow and steady 1.53 °C warming starting from the early twentieth century, increasing greatly in the last quarter of the century.Fifty percent of this warming occurred between 2000 and 2010.Main geographical difference between Bartin and other locations in this study is that Bartin is located on the Black Sea side of the orographic effect.The results can also be explained by urban heat island effect due to its proximity to the center of the district as observed in other provinces of Turkey (Tayanc ¸& Toros, 1997).
Our results show a rapid increase in GST starting from the late 1990s.Projections on regional effects of climate change on the Eastern Mediterranean show that Anatolian Peninsula is among regions that will have a dramatic increase in surface temperatures (Gao and Giorgi, 2008).Studies predict that surface air temperatures in Western Black Sea Region will increase 4.5-5 °C during 2070-2100 compared to 1961-1990 average temperature values and are in agreement with our results (O ¨nol et al., 2009).

Conclusion
In this study temperature-depth data from five boreholes were used to reconstruct ground surface temperature history in the Western Black Sea Region.The results of this study can be used in regional climate projections for this region.
Figure 1 Map of the study area.Labels marked with red squares show borehole locations; italic labels with blue markings show province centers, and yellow labels with yellow triangle markings show locations of Turkish State Meteorology Service (TSMS) stations used in this study

Figure 3 T
Figure 3 T-D curves of boreholes used in climate reconstruction.Water tabel levels are given in Table 1.Dashed lines show background conditions.Both 2021 and 2022 measurements for Bartin are shown.Measurements from depths [ 25 m are omitted

Figure 4
Figure 4 Reconstructed GST histories for Western Black Sea Region study sites.Akoren, Bartin and Dizdarli have an upper and lower range for results due to unavailability of thermal conductivity measurement samples.In their GST plots, orange lines represent k = 3 W/m/K and blue lines represent k = 2 W/m/K conditions.Vertical bars show SD values provided by the FSI model

Figure 5
Figure 5 Measured transient temperatures (diamonds) and modeled transient temperatures (lines).For Akoren, Bartin and Dizdarli blue lines show k = 2 W/m/K and orange lines show k = 3 W/m/K conditions.Pearson correlation coefficient (r) of modeled and measured temperature profiles for each borehole is shown at the top right corner.Vertical lines at 0 °C correspond to the background thermal conditions as a reference point for each borehole

Figure 6
Figure 6 Comparison of reconstructed GST histories and SAT history records from TSMS stations

Four
of the boreholes are located in the inland section, and one of them is located in the coastal section.Reconstructed GST histories of all boreholes starting from the year 1900 show B 0.1 °C cooling until 1970 and 0.99 warming between 1970 and 2010, adding up to 0.91 °C net increase on average until 2010.The averaged GST models of inland section show a 0.24 °C cooling period between 1900 and 1975 prior to 0.94 °C average warming in Akoren,

Figure 7
Figure 7 Mean GST histories.Horizontal lines below the temperature values show correlation length

Table 1
Properties of the studied boreholes and sites H elevation, LD logging date as month.year,BD borehole depth, WT depth of water table, D borehole diameter, T 0 undisturbed GST at time = 0, k thermal conductivity, qb background heat flow

Table 1 .
Dashed lines show background conditions.Both 2021 and 2022 measurements for Bartin are shown.Measurements from depths [ 25 m are omitted