Influence of Climate Change on the Thermal Condition of Yakutia’s Permafrost Landscapes (Chabyda Station)

s: This paper presents the results of 39 years of observations conducted at the Chabyda station to monitor the thermal state of permafrost landscapes under current climatic warming. The analysis of long-term records from weather stations in the region has revealed one of the highest increasing trends in mean annual air temperature in northern Russia. The partitioning of the energy balance in different landscape units within the study area has been analyzed. Quantitative relationships in the long-term variability of ground thermal parameters, such as the ground temperature at the bottom of the active layer and seasonal thaw depth, have been established. The ground temperature dynamics within the depth of zero annual amplitude indicates that both warm and cold permafrost are thermally stable. The short-term variability of the snow accumulation regime is the main factor controlling the thermal state of the ground in permafrost landscapes. The depth of seasonal thaw is characterized by low interannual variability and exhibits little response to climate warming, with no statistically significant increasing or decreasing trend. The results of the ground thermal monitoring can be extended to similar landscapes in the region, providing a reliable basis for predicting heat transfer in natural, undisturbed landscapes.


Introduction
In the context of research into global climatic change, interest in the problems of the response of permafrost to these changes has increased. Research on the impact of climate warming on permafrost involves a wide variety of issues, including the thermal evolution of upper permafrost in natural and disturbed landscapes under anthropogenic effects. This research has become a priority in geocryological science, with theoretical and practical significance worldwide.
In the face of the major climate changes in recent decades in Central Yakutia [1]-a very densely populated and promising region in terms of agricultural and industrial development-obtaining information on the response of permafrost to modern warming is critical.
Pioneering systematic observations of temperatures in 10-15 m depth bore holes were conducted in 1935. The study revealed peculiarities of temperature variations and ground thaw depth under the effects of vegetation and snow cover, the "cultural layer", and geological and geomorphological parameters [2,3]. In the 1940s, year-round experimental observations were conducted to study ground temperature regime variations under thermal insulation, snow and vegetation cover removal, and in natural landscapes of experimental sites of the Yakut permafrost research station [4][5][6]. In the 1950s-1960s, more detailed integrated thermal physics observations of t h e shallow soil formation thermal regime were performed, applying updated methods at experimental sites of the North-Eastern department of the Permafrost Institute [7].

Research Sites and Methods
The Chabyda station is located in the middle taiga zone and is characterized by both low and high-temperature permafrost zones. A detailed description of the observation sites and their location in relation to the station is given in a number of publications [12,13].
The relief of the study area is characterized by a combination of lowered and elevated sections of a hilly-steep plain. Lower areas are the bottoms of streams and gullies, while dry elevated areas are slopes of various steepness and watershed spaces. On the slopes, soils are represented by fine and medium-grained sands, soils on the watersheds in the upper horizons are sandy loam, and those below are sand. The objects of research are soils within the layer of annual heat rotation; in other words, the layer of annual temperature fluctuations (the upper 10-20 m).
Depending on the terrain conditions, soil texture, moisture content, vegetation, and surface cover, the active layer thickness varies over a wide range, from 0.4 to 4.0 m, and the mean annual ground temperature at a depth of 10 m varies from -0.2 to -5.0 °C [13].
Landscape and monitoring methods were utilized to study the ground thermal regime. Landscape studies involved remote sensing and ground surveys of landscapes and their classification and mapping. The information obtained was used for site selection and for the design of a thermal monitoring program.
The experimental sites (S) were organized into two landscape types of terrain. In the stream valley: on the mari (S-1, S-3a); on the grassy lowland (S-8a); in the larch forest (S-8). On slopes: in pine forest on a gentle slope (S-4, S-5, S-6b); in pine (S-7b) and larch (S-9) forests on a watershed slope; in pine forest on moderately steep slopes of the Northern (S-10) and southern (S-11) expositions.
The ground thermal regime within the depth of annual temperature fluctuations is determined both by external (solar radiation, air temperature, and precipitation) and internal (moisture conditions and lithology) factors. Snow cover, vegetation, and organic mat are also important controlling factors [10][11][12][13].
The radiation balance is the result of the incoming and outgoing solar energy at the Earth's surface and is described by the following equation: where Q is the total incoming radiation; S is the reflected radiation; Jef is the effective radiation (net long-wave radiation); and R is the net radiation. The components of the heat and moisture exchange of the Earth's surface with the atmosphere are expressed by the heat balance equation or the energy balance equation on the Earth's surface [10]: where P is the sensible heat flux; LE is the latent heat flux (E is the evaporation rate and L is the latent heat of evaporation of water); and B is the ground heat flux. Long-term studies at a thermal monitoring station in Yakutsk suggest that all components of the energy budget are half or one order of magnitude lower during the cold season compared to the warm season [10]. Differences in the energy balance characteristics among landscape types are therefore more distinctive in summer. Summertime net radiation can be determined from instrumental measurements or by calculation. Pavlov [10] proposed a simple model to estimate the energy balance using incident solar radiation and surface albedo. We measured albedo at the Chabyda research sites during summer months ( Table 1). The measured А and ho were used to calculate the monthly sums of surface net radiation over the warm season, LE and B were derived from instrumental measurements, and Р was estimated as a residual of the energy balance equation. Meteorological and actinometric observations were made using standard instruments of Roshydromet. The permafrost thermal monitoring program included measurements of ground temperature and thaw depth, as well as the major controlling variables, such as ground material, cryostructure, moisture content, soil unit weight, thermal conductivity of soils and surface covers, and snow depth and density.
Ground temperatures were measured using MMT-4 semiconductor thermistors with an accuracy of 0.1 °C. Until 1990, the measuring devices used were the potentiometer PP-63 and the resistance bridge MO-62; since 1990, the multimeters EDM-169S and EDM-89S have been used. Thermistor cables placed in backfilled boreholes recorded temperatures at depths of 1, 2, 3, 4, 6, 8 and 10 m. Active layer depths were determined by probing with a steel rod and by hand boring. The analysis and summary of the observation results at various monitoring stages were presented by Skryabin et al. [12] and Varlamov et al. [13].
Observations of the thermal regime of the active layer at the Chabyda station were carried out in the warm seasons of 1981 and 1982 daily, in seven shifts, and in the cold season of 1982-1983 in four shifts. In the annual cycles of 1983-1986, in the summer, observations were carried out once in a pentad, in four shifts. In 1981-1987, soil temperature measurements in the layer of annual temperature variations were carried out once a decade. Since 1987, thermal observations have been carried out according to a more abbreviated program on the 15th of each month. Snow cover is also observed monthly, and at the end of the warm season, the depth of seasonal thawing is observed. The methodology of geothermal monitoring has shown its reliability and can be successfully used in various climatic conditions.

Modern Climate Changes
Since the second half of the 1960s, according to Skachkov [1], Central Yakutia has experienced the highest increasing trend in mean annual air temperature in northern Russia (up to 0.07 °С/year). In 2001-2019, the highest ever temperatures in the history of weather observations were registered in Yakutsk, at −7.7 °C as compared with an average of −10.0 °C. The evaluation of predicted trends of temperature variations in the 21st century by various researchers is ambiguous. For example, predictions of climatic change performed by the Voeikov General Geophysical Observatory (GGO) based on climatic feature extrapolation results suggest that the observed increasing trend in the Russian mean annual air temperature will continue to 2030 by 1.5 °C as compared with values from the year 2000. A further increase of the mean amount of precipitation, predominately in the cold season, is predicted. In the north of eastern Siberia, an increase of precipitation amounts by 7%-9% and of snow accumulation by 2%-4% in winter is expected [27]. In the late 21st century in eastern Siberia, the increase of the near surface mean annual air temperature is expected to be 6.1±1.6 °C compared to the values of the last two decades of the 20th century [ 2 8 ] . The expected near surface mean annual air temperatures in Yakutsk up to 2015, as estimated at the Permafrost Institute of the Siberian Branch, Russian Academy of Sciences, applying the frequency analysis method, almost coincide with the predictions made by GGO. Furthermore, a temperature decrease of 2-3 °C is expected up to the first half of the 21st century [29,30]. According to predictions made by Neradovky and Skachkov [31], by 2050, an increase of the mean annual air temperature will pass the achieved climatic level of no more than 0.7-1.0 °C. Researchers are increasingly concerned about an expected increase of the amount precipitation in winter and snow accumulation, both of which play a significant role in the increase of soil temperatures.
The variability of the main climate elements (air temperature and precipitation) over the study period at Chabyda station can be determined from the data from Yakutsk weather station-the nearest one to the study area (Figures 2 and 3). We should mention the high correlation between Yakutsk and other weather stations in Central Yakutia [1]. The mean annual air temperature (1961-1990 climate norm) in Yakutsk is -10.0 °C and the mean annual precipitation is 235 mm.  As shown in Figure 2, the trend of the mean annual air temperature increase is significant and exhibits stable growth. This increase was mainly caused by the growing temperatures in winter (October to April). In summer (May to September), warming was less evident.
The data presented in Figure 3 indicate that the recent decades have seen major year-to-year variations of annual precipitation totals (October to September), but these have not increased in Yakutsk on the whole. The precipitation total by months and for entire years observed during the last 40 years is close to normal. Nevertheless, it should be noted that, in some years, abnormally large (1989,1993,(2005)(2006)(2007)2013) and abnormally small (1986,2001,2004) levels of precipitation occurred.
Thus, the variability of the main climate parameters (air temperature, precipitation, and snow cover) over the past decades is manifested in different ways. The positive trend of the average annual air temperature is the most significant. Precipitation and snow cover height experienced short-period fluctuations without a clear trend. Annual precipitation (mm)

Surface Energy Balance
A brief review of the radiation regime in the region is presented below based on radiation measurement data from the Yakutsk weather station [12]. The incoming short-wave radiation is largely controlled by sun angle, cloud coverage, and day length. The incident solar radiation varies seasonally over a wide range from 17 (December) to 620-712 MJ/m 2 (June). The annual total of radiation is 3710 MJ/m 2 . Approximately 34% of the incoming radiation is reflected back into the atmosphere. The energy loss is greater during the winter period (75%-85%) than in the summer (18%-22%). About one-third of the annual incoming radiation is lost by net long-wave radiation. The heat loss by net long-wave radiation increases in proportion to the increase in surface temperature from 29-63 MJ/m 2 in the winter months to 167-188 MJ/m 2 in the summer. During the winter, the net radiation is negative due to the loss of heat by terrestrial radiation and the reduced solar radiation. Positive net radiation values occur between late March and mid-October. The average annual total of net radiation is 904 MJ/ 2 , or 33% of incoming solar radiation.
The landscape units displayed different energy balance characteristics depending on surficial factors. The partially inundating riparian area covered with mosses and grasses (site 2) showed little variation in albedo during the summer months (18%-20%). At the creek valley bottom (site 3), a 0.1m-thick icing melted out by late May; therefore, the heat loss by reflection was 1.4 times higher than at the riparian site. From June to September, the reflectivity of the surface at site 3 with a moss-shrub cover and hummocky topography was consistently low (15%-16%). At the lower and upper slopes (sites 4 and 5), the albedo of bearberry patches with dead plant spots under the open pine canopy varied only slightly over the summer. At the bottom of a drainage trough (site 8), the albedo of the moss-ledum-lingonberry cover under the closed birch and larch canopy was reduced to 12%-13 %.
The net radiation (R) of the moss-shrub cover on the drainage trough bottom was 1.04 times larger than the R of the riparian site due to the lower albedo.  (Table 2). The evapotranspiration at each site was determined from measurements with two to four evaporation pans to account for vegetation type ( Table 3). Rates of soil water evaporation depend on energy factors and moisture conditions. The evaporation rates on the creek valley bottom were higher between May and July, when the monthly totals of net radiation were largest. In August and September, the reduced net radiation resulted in a decrease in evapotranspiration by a factor of 1.5-2 despite the sufficiently high soil moisture contents. At the forest bog site with a moss-grass-shrub cover, the evapotranspiration averaged over two seasons was 230 mm, which was 1.5 times the mean evapotranspiration from forb meadow at the Yakutsk monitoring station. The seasonal average latent heat flux LE at the forest bog site comprised 42% of R.
At the open-canopy forest site on the lower slope (site 4), the evaporative capacity of the vegetation covers was not the same. For example, the evapotranspiration Eс from the lingonberry cover was 2.3 times higher than that from the lichen cover. The evapotranspiration under the pine forest was greatly reduced and averaged 110 mm due to both the reduced R and the lower moisture storage in the sandy soils. Compared to the forest bog (site 3), the evapotranspiration under the pine forest was approximately twice as low. The ratio of LE to R was 0.26. In Central Yakutia, the LE/R values reported by Pavlov [10] were 0.18 and 0.35 for open-canopy pine and larch forests, respectively.
An accumulation of heat in the ground at the riparian site occurred between May and September. The ground heat flux was largest in May and constituted 20% of R. In August and September, ground heat fluxes comprised only 3%-7% of R. Overall, 159 MJ/m 2 of heat, or 12% of R, entered the ground during the thaw season. In October and November, when the snow cover was shallow, the heat flux was strongly negative, comprising half of the annual heat loss (Figure 4). Between January and April, heat losses were low due to the thick snow cover. The amount of heat utilized in the snowmelt process was 20.9 MJ/m 2 . Because snow cover depths were greater than the long-term average, the input components of the ground heat flux through the surface exceeded the outputs. The difference between the inputs and outputs comprised 28% of the magnitude of the annual heat exchange.
At the forest bog site, the seasonal sum of sensible and latent heat fluxes, P + LE, was 60 MJ/m 2 greater than at the bank. Therefore, the ground heat flux was, on a seasonal average, 8% lower despite the higher R values. The B/R was 11%. During the cold season, the heat flux to the atmosphere was strongly reduced due to sustained groundwater flow down the gentle slope onto the bog surface. On an annual basis, the heat input was 1.56 times the heat loss.
The drier sands with a low thermal conductivity under the open canopy (site 4) accumulated 1.7 times less heat during the warm season compared to the forest bog. At site 5, heat accumulation was strongly reduced in the dry sands due to the small temperature gradients near the surface. The B/R ratio was 6%, which was also due to the lower moisture content of the sands. Comparing the absolute values of the energy balance components among the sites, it can be noted that they are of the same order for the gentle slope sites. The only difference is in the magnitude of the sensible heat fluxes, which is attributed to site differences in canopy shading. At the drainage trough bottom (site 8), the LE and B values are much higher due to the higher soil moisture content and the significant reduction in sensible heat flux under the larch canopy. As a result, the input component of heat exchange at this site is 1.7-1.9 times larger than at the slope sites.
In summary, investigations in the Chabyda area provided important information about the radiation and heat balance characteristics for different landscape units. During the warm season, the albedo is highest for the riparian site with a moss-grass cover

Analysis of Long-Term Variability of Soil Temperature
To date, there is a database for sites 1 and 5 for 39 years of observations, for site 8 for 37 years, for site 9 for 35 years, and for other sites there are 34 years of observations. Such representative material allows a qualitative analysis of the variability of the main thermal parameters of the layer of annual temperature variations and allows us to draw objective conclusions. The idea of the longterm variability of quantity ξ, Tξ, and Т10 at the experimental sites is given in Table 4.
The greatest fluctuations of Тξ and T10 occurred in landscapes of a shallow-valley type of terrain (see Table 4 and Figure 5a, b). These changes were mainly determined by the influence of two winter factors: snow accumulation conditions and the sum of air temperatures during the cold period. Moreover, the first factor is prevailing in the sharply continental climate of Central Yakutia [13,32,33].
Over the entire observation period, the lowest and highest average annual temperatures at the sites of the shallow-valley terrain type were noted in hydrological years 2003/04 and 2006/07, respectively (October-September) (see Figure 5). In the slope type of terrain, the peak of the lower temperature of the soil also occurred in 2003/04, and the peak of the temperature increase was noted in 2007/2008. The exception here is site 7b, where the lowest soil temperature was observed in 1987/88 ( Figure 6).  The winter in 2002/03 featured abnormally light snow and had an abnormally late period of formation of stable snow cover (Figure 7). This was the main reason for the strong cooling of the soil, despite a fairly warm year. Subsequent years were characterized by snowy winters and higher than normal rainfall, which led to a sharp increase in soil temperature. In the period from 2002/03 to 2006/07 (2007/08), the temperature of the soils at the bottom of the active layer in these types of terrain increased by 0.5-6.0 °С and at the depth of zero annual amplitude by 0.3-2.7 °С (see Figures 5  and 6).
The warming of the permafrost and deepening of the active layer have accelerated thermokarst, thermal erosion and other permafrost processes in the region. Water ponding and paludification have been observed in the low lying areas. These changes have resulted in increased water levels in the existing lakes, the filling of the previously dry lake basins, and in the rise of the water table in the active layer. In the water-logged topographic lows, residual thaw layers have developed due to the lowering of the permafrost table and significant warming of the near-surface ground [20,23]. Thus, 2007 and 2008 were critical years for the permafrost thermal regime in the 40-year period of observations or even in its century-long history.
The winters of 2008/09 and 2009/10 were characterized as abnormally snowy with late periods of formation of stable snow cover; thus, the soils were much cooled to almost the level of 2002/03, despite the abnormally warm years. Subsequent winters were characterized by close to long-term average values of average winter snow depths, and therefore, mean annual ground temperatures at the active-layer bottom and at the depth of zero annual amplitude were close to long-term average values. It should be noted that the correlation between the temperatures at the bottom of the active layer and the average winter snow depths is quite high (correlation coefficients are 0.55-0.74) [13].
According to long-term observations, trends in the increasing average annual temperatures of the layer of annual heat rotations are extremely weakly expressed. Noticeable trends in temperature increase were observed only at two sites: 7b (0.23°C/10 years) and 8a (0.33°C/10 years) ( Table 5). This was largely caused by the abnormally snowy winters mentioned above, and site 8a was flooded during several warm seasons with meltwater. In addition, the first years of the analyzed period (1986)(1987)(1988) were quite cold. It should also be noted that in these types of terrain at sites 8 and 9, there are decreasing trends of T10. This is due, first of all, to an increase in the shading of the surface under the canopy of the stand with its crowns and the growth of shrubs and shrubs. Table 5. Trends in mean annual ground temperature at the depth of zero annual amplitude (10 m). According to long-term observations, the thermal effect of abnormal winters has been experimentally established. Thus, one warm and less snowy winter can lower the temperature of the soil more than a cold and snowy one. An abnormally cold winter can produce a stronger effect than several abnormally warm ones following one after another [32].

Terrain Type Ground, Experimental Site (S) and Number Observation Periods
The most important conclusion from the results of long-term observations is that the temperature of the soils experiencing significant interannual and short-period fluctuations associated mainly with the great variability of the snow accumulation regime has a very weak tendency to increase. Against the background of significant changes in the average annual air temperature (in Yakutsk for the period 1981-2019, the growth trend was 0.064 °C/year), the thermal state of the layer of annual temperature variations at the slope-type sites of the Chabyda station area remains stable.

Analysis of the Long-Term Variability of the Depth of Seasonal Thawing
As is known, the depth of seasonal thawing (ξ) depends on the mechanical composition of soils, their moisture content, the nature of the vegetation and soil cover. The greatest long-term variability ∆ξ is typical for soils of a shallow-valley type of terrain and makes up 30%-55 % of their maximum extremum. In the slope type of terrain, ∆ξ varies within the range of 14%-24 % (see Table 4 and Figure 8). It was previously believed that the main factors determining the long-term variability of the thickness of the seasonally thawed layer are the sum of positive air temperatures and summer precipitation. However, recent developments on this issue show that long-term changes in the depth of seasonal thawing and the sums of summer air temperatures do not correlate sufficiently well [34].   Most likely, the ground conditions-namely the moisture content-can majorly affect the thawing depth. An analysis of the data obtained at the station shows that the largest long-term variabilities of ξ are characteristic of landscapes of a shallow-valley type of terrain (up to 80 cm). On sites in the slope type of terrain, these changes were 30-60 cm (see Table 4).
In the shallow-valley type of terrain at sites 3a and 8, a significant increase in the depth of seasonal thawing was noted. Significant decreasing tendencies are observed at sites of the slope type of terrain (sites 5, 6b, 9,11). Additionally, if in the watershed area (Site 9) this decrease can be attributed to the intensive growth of shrubs, then for sites 5 and 6b, fluctuations in the level of permafrost waters of the seasonally thawed layer most likely play a significant role. At other sites (7b, 10), with significant interannual changes in the depth of seasonal thawing, insignificant growth tendencies are observed ( Table 6). Table 6. Trends in the depth of seasonal thawing of soils of the Chabyda station. The maximum depths of seasonal thawing at individual sites were observed in different years. In the shallow-valley type of terrain (sites 1, 3a, 8, 8a), they were noted in 2007-2008. In the slope type of terrain, the maxima ξ were revealed in the following years: site 5-1984, site 6b-1995, site 7b-1989, site 9-1991, and site 10 and 11-1988. Thus, during climate warming, even within local undisturbed areas, an opposite trend can be detected in modern changes at the depth of the seasonal thawing of soils. This serves as strong evidence that the thickness of the active layer is not a sensitive indicator of modern climate change.
The analysis of a large amount of experimental data confirmed the conclusion about the weak response of the depth of seasonal thawing to modern climate changes [35]. The observations convincingly show that under the conditions of modern climate warming, an increase in the depth of seasonal thawing of soils does not always occur. Maxima and minima of this value at various sites of the station are observed, most often in different years. This indicates that the depth of seasonal thawing depends not only on long-term changes in air temperature in the warm period but also on other meteorological factors.

Conclusions
At the end of the last century, we concluded that the thermal regime of permafrost soils was stable and their reaction to climate variability was weak [36]. Subsequent years have confirmed this conclusion.
Based on the analysis of the 39-year monitoring, the main conclusions are as follows: 1. The radiation and heat balance characteristics have been determined for different landscape conditions in the Chabyda area.
2. The long-term dynamics of the thermal state of the layer of annual temperature variations during climate warming indicates the thermal stability of both high-temperature and lowtemperature permafrost. The main regulatory factor in the dynamics of the thermal state of permafrost soils is the snow accumulation mode.
3. The depth of seasonal thawing, despite significant interannual fluctuations, weakly responds to climate warming and has no significant trends. 4. The results of studies of the thermal regime of soils can be extended to similar landscapes of Central Yakutia. 5. The practical value of the materials obtained at the Chabyda station also lies in the fact that the observation results can be used to model heat transfer processes in natural landscapes.
Author Contributions: All authors-Varlamov S.P., Skachkov Y.B., and Skryabin P.N.-made a significant contribution to the published work. All authors have read and agreed to the published version of the manuscript.
Funding: This research received no external funding.