MAPPING PERMAFROST DISTRIBUTION IN THE PARVATI VALLEY, KULLU USING LANDSAT 8 DERIVED LAND SURFACE TEMPERATURE

Permafrost covers much of the Parvati valley, but mapping its distribution is difficult due to the scarcity of ground observed data sets and the region's high spatial heterogeneity. We show permafrost distribution maps for the Parvati valley from 2013 to 2021 based on the Mean Annual Air Temperature. We created maps of Mean Annual Air Temperatures (MAAT) using Landsat 8 Land Surface Temperature (LST) products in Google Earth Engine. Google Earth Engine (GEE) is a cloud-based platform that was utilised to quickly and efficiently obtain the spatial and temporal variations of permafrost distribution. Permafrost is defined as ground that maintains a temperature below zero degrees Celsius for at least two consecutive years. To justify the definition, we examined Mean Biennial Air Temperature (MBAT), Mean Triennial Air Temperature (MTAT), Mean Quadrennial Air Temperature (MQAT), Mean Quinquennial Air Temperature (MQnAT), Mean Septennial Air Temperature (MSpAT), Mean Octennial Air Temperature (MOAT), and Mean Novennial Air Temperature (MNOAT) (MNAT). Our findings demonstrate that the percentage of Permafrost distribution in the Parvati valley is about the same in all situations, accounting for 22 percent of the overall study area excluding glaciers. Our maps were divided into four categories: Continuous Permafrost Zone, Discontinuous Permafrost Zone, Sporadic Permafrost Zone, and No Permafrost Zone.


INTRODUCTION
Permafrost is a vital part of the cryosphere, and it is expected to rapidly degrade as a result of present climate warming (Batbaatar et al. 2020). It covers approximately 22% of the land area in the Northern Hemisphere (Brown et al. 1997). Permafrost is extremely vulnerable to climate change. It has temperatures below freezing and contains ice-rich material near the surface, which could thaw as a result of climate warming, causing enormous landscape rearrangement via thermokarst (Panda, et al., 2014).
Temperature measurements are necessary for assessing the state of permafrost, which is at risk of thawing as the temperature warms (Osterkamp et al. 2009). Annual maps of ground surface Temperature (GST) could be utilised as inputs for regional models of the permafrost thermal regime, allowing for monitoring of climate change at the ground surface level (Hachem et al., 2009). Permafrost's thermal condition has traditionally been represented by the air temperature (Ta). Because of the thermal effect of surface characteristics and substrates, different MAAT (Mean Annual Air Temperature) isotherms are used to determine the continuation of permafrost in different areas (Luo et al. 2018). Monitoring near-surface air and GSTs using in-situ approaches across broad areas is difficult and expensive logistically. Remote sensing, which collects spatially dense data over broad areas, is used to solve this challenge. Remote sensing can be used to quantify permafrost extent and identify changes by remote identification and mapping of surface landmasses that validate the existence of permafrost and remote sensing of physical parameters that correspond to thermal subsurface conditions (Westermann et al. 2015).
Over the last decade, remote sensing data processing has shifted from traditional workstations outfitted with cutting-edge (and often very expensive) hardware and satellite image processing software to cloud-based platforms that enable users to instantly access and analyse massive pre-processed geospatial data via user-friendly, web-based interfaces and powerful scripting languages. Google Earth Engine (GEE) is one of such platforms. Google Earth Engine (GEE) is a web-based platform that allows remote sensing users to easily do massive data analysis without the need for computational resources (Gorelick et al. 2017). It enables users to tackle the fundamental difficulties associated with the administration of very huge volumes of data, such as storage, integration, processing, and analysis, in a very effective manner (Tassi and Vizzari 2020;Gorelick et al. 2017).
There is a paucity of data on the extent and consequences of permafrost thawing across the Himalayan region as a whole (Gruber et al. 2017). Permafrost mapping in the Indian Himalayan region has largely relied on spatial interpolation between near-surface air temperatures measured 1-3 metres above ground at widely dispersed meteorological stations (Hachem et al., 2009) or downscaling medium resolution satellite data, such as MODIS satellite data, to 30m resolution (Haq and Baral 2019;Khan et al. 2021;Allen et al. 2016). To date, no studies have used LSTs from the high-resolution Landsat 8 platforms to map surface temperatures and derive regional thermal conditions for permafrost distribution in the Indian Himalaya.
In this study, we used the Google Earth Engine (GEE) platform to map the Permafrost distribution (Mean air temperature below 0°C for at least 2 years) of Parvati Valley, Kullu, in Northwest Himalaya utilising Landsat 8 derived Land Surface Temperature (LST).

Study Area
The Parvati Valley (31° 58' 41'' N to 32° 05' 51" N Latitudes and 77° 14' 23" E to 77° 27' 08'' E Longitudes), located in the Kullu district of Himachal Pradesh along the Parvati River, encompasses an area of 1,766 sq km ( Figure 1). It is a narrow valley with steep slopes. The Parvati River is the main drainage of the catchment, and it is fed by tributaries such as Malana nallah, Tosh nallah, Grahan nallah (Sharma and Samant 2017). The elevation varies from 1100 metres at Bhuntar town, where the Parvati river meets Beas river, to 4100 metres in Mantalai Lake. The valley has warm and temperate (900-1800 m), cool and temperate (1100-2400 m), and cold alpine climatic distribution. It also has glacial climate (2400-4800 m) in the north and east. In the alpine sections of the valley, significant rainfall (800-900 mm) and substantial snowfall (8-9 m) are common (Gupta et al. 2017). The glaciers cover an area of 361.91 sq km of the total Parvati valley ( Figure 1).

Data
In this study we used Landsat-8 Surface Reflectance Tier-1 and Landsat-8 Top of the Atmosphere (TOA) brightness temperature Tier-1 images from 2013 through 2021. Landsat 8 has a spatial resolution of 30 m and a temporal repetivity of 16 days. The Landsat 8 SR (Surface Reflectance) T1 (Tier 1) dataset is an atmospherically corrected surface reflectance data product obtained from the Landsat 8 OLI/TIRS sensors using the LaSRC algorithm. These images contain five visible and near infrared (VNIR) bands, two short wave infrared (SWIR) bands, and two thermal infrared (TIR) bands. In addition, we used Randolph Glacier inventory (RGI) data from GLIMS (Global Land Ice Measurements from Space) to extract glaciers from our study area. The GLIMS project gathers data from the Landsat Enhanced Thematic Mapper Plus (ETM+), the Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER), and publicly available historical data. The data of this study is freely available on Google Earth Engine (GEE), which is given and maintained by the USGS.

Methodology
We imported the Landsat 8 SR (Surface Reflectance) and Landsat 8 TOA (Top of the Atmosphere) data from 2013 to 2021 on the Google Earth Engine and clipped it to our study area. We applied cloud masking to remove cloudy pixels from the image collection. For the calculation of LST, values of brightness temperature and surface emissivity are required. Brightness temperature was calculated from the TIR band (Band 10) of the Landsat 8 TOA data and surface emissivity was calculated through NDVI. To calculate NDVI, we utilized the Red (Band 4) and NIR (Band 5) bands of the Landsat 8 Surface Reflectance product using equation 1. The fractional vegetation cover (FVC) was calculated with the help of NDVI using equation 2. Surface Emissivity values were calculated using the equation 3. We used the Statistical Mono Window (SMW) Algorithm proposed by (Duguay-Tetzlaff et al. 2015) and utilised by (Ermida et al. 2020) to calculate LST.
In the equation (4) Tb is the TOA brightness temperature of the TIR channel, and ε is the surface thermal emissivity in the same TIR channel and the values of the coefficients Ai, Bi, and Ci were provided by (Ermida et al. 2020).
A time series analysis was done on Land Surface Temperature values obtained for all scenes from 2013 to 2021. The LST value was missing for several days at certain locations due to cloud masking. The cloud-free LST values were taken and outliers were identified and removed using the box plot to account for any possible bias. We validated the derived LST values against the IMD (Indian Meteorological Department) observed air temperature for various IMD stations of Himachal Pradesh.
Using linear regression analysis, we formulated the empirical equations between Mean Annual Air Temperature (MAAT) and Land Surface Temperature (LST) for each location. As the Parvati valley is located in the Kullu district, so we applied the empirical equation derived for Kullu district to our study area and calculated MAAT values from the satellite derived LST values.
As the definition of Permafrost states that it is ground that has remained at or below 0° C for at least two years. This suggests that consecutive two or more than two years could be considered for identifying the permafrost areas. So, to better understand permafrost distribution and changes in its characteristics between 2013 and 2021, we considered all the possible combinationof consecutive years such as at least 2 years, 3 years, 4 years, up to for at least 9 years (as from 2013  We then removed the glaciated regions from the study area by clipping the GLIMS Randolph Glacier Inventory data so as to compute and carry out the long-term changes that have happened in the permafrost areas. The detailed methodology in the flowchart is given in figure 2. We applied a threshold of Air temperature < 0°C to all the 36 maps. We assumed MAAT < -5 °C to be Continuous Permafrost, -5°C < MAAT < -2°C to be Discontinuous Permafrost and -2°C < MAAT < 0°C to be Sporadic Permafrost. So, we reclassified the maps as Continuous Permafrost (CP), Discontinuous Permafrost (DP), Sporadic Permafrost (SP) and No Permafrost. We analysed the pattern of permafrost distribution and the changes occurred over the years.

RESULTS
The Our result indicates approx. 22% area (excluding glaciers) of the Parvati valley is under Permafrost throughout these time period.
In every combination, a cyclic pattern can be seen in the distribution of Permafrost. Whether we do our analysis for at least two years or for nine years, the % area of permafrost is nearly equal. This can be seen in Table (1-8).   (SP), it was observed that continuous and sporadic were following the inverse trend. That is when CP was increasing, SP was decreasing and when SP was increasing CP was decreasing while DP remained more or less stagnant. This indicates that during this time period due to cyclicity behaviour, CP was converting to SP thus decrease in CP lead to increase in SP and vice versa as shown in figure 4 and indicated by table 1. The explanation for this is the cyclic nature of climate warming and cooling. When the climate warms, the air temperature rises, melting the permafrost. As a result, sporadic permafrost develops. The amount of SP increases while the CP drops. When the environment cools, the water in the Permafrost freezes and produces CP. We carried out similar work for other consecutive time periods as explained earlier in the text and we observed similar trends and results. We prepared all the figures for these, but to remove the redundancy, we are not showing those maps here but just giving the data in the tabular form (Table 2 to Table 8). The cyclic pattern of permafrost distribution is very clear and the conversion of CP to SP and SP to CP is also clearly seen in the data and its graph (not shown here) in all the data sets.        The work carried out in this paper suggest that for Parvati valley even if we take minimum of 2 consecutive years or upto 9 consecutive years, there is no significant change in the area of permafrost distribution. The changes from CP to SP and vice versa is very prominently visible. This indicates that the degradation of permafrost is a slow process in this area. However similar work needs to be carried out in other parts of the global permafrost to assess the temporal variation of the permafrost areas.

DISCUSSION
Nearly one fourth of the overall study area is underlain by Permafrost between the years 2013 and 2021. Cyclicity effect plays a major role in the distribution of Permafrost. There was no major difference in the distribution of Permafrost in all the 8 cases. Rise or fall in the percentage area of permafrost is due to new formation or thawing of Permafrost respectively. Out of total Permafrost area around 36% is continuous Permafrost, 37% is discontinuous Permafrost and 27% is Sporadic Permafrost. We validated our result with the help of Google Earth Pro, by taking rock glaciers as a proxy for the presence of Permafrost as shown in the figure 6. The response time of rock glaciers to climate warming and cooling is longer than that of glaciers. Because the ice in a rock glacier is protected by a debris mantle. The insulated core allows an ice record to accumulate and be stored for an extended length of time (Burger, Degenhardt Jr, and Giardino 1999). Active rock glaciers have a substantially lower average yearly mean specific discharge than glaciers (Krainer and Mostler 2002). A debris-covered glacier can convert into a considerably shorter rock glacier as the environment warms (Anderson et al. 2018). Occurrence of many rock glaciers in the permafrost area shows that the changes happening in this terrain is a slow process and needs to be studied in detail in the Himalayan region. To acquire a better understanding of temporal variation of permafrost, the same study must be conducted for other global permafrost regions.

CONCLUSION
This study used better-resolution Landsat 8 images to determine mean LST for the Parvati Valley to construct better permafrost maps than any previously created maps in this region by interpolating data from climate stations or resampled satellite data. Permafrost distribution is modelled using MAAT calculated from spatially continuous LST. Our results demonstrate around 22% area of the Parvati valley is Permafrost region, out of which approx. 36% area comes under Continuous Permafrost zone, 37% under Discontinuous Permafrost zone and 27% as Sporadic Permafrost zone. During these years, Continuous Permafrost is converted to Sporadic Permafrost when there is warming, and when there is less warming or enhanced colder period then Sporadic Permafrost is converted to Continuous Permafrost. This probable permafrost map developed in this study will be useful for future research since it identifies areas in need of more indepth ground investigations.