Impact of Urbanization on Urban Heat Island Dynamics in Shillong City, India Using Google Earth Engine and CA-Markov Modeling

Abstract

Growth in urban areas contributes to environmental degradation through increased land surface temperature (LST), exacerbating the urban heat island (UHI) effect. This study examined how land use and land cover (LULC) characteristics of Shillong City are linked to the UHI phenomenon. The LULC was classified into five broad categories: agricultural land, barren land, settlement, vegetation, and water bodies. The results show that the study area experienced notable changes in the LULC pattern from 1993 to 2023, with settlement areas increasing by 10.96%, transforming previously barren lands. The emergence and growth of settlements (and/or built-up areas) and impervious surfaces have led to a steady increase in LST. The settlement land use class had an average LST of 17.45 °C in 1993, 21.56 °C in 2003, 21.37 °C in 2013, and 21.75 °C in 2023. From 1993 to 2023, surface temperatures in settlement areas rose by a maximum of 4.3 °C, while barren land and vegetated areas also saw an increase of 4.9 °C and 4.0 °C, respectively. The relationship between LULC and the LST has been evaluated to identify hotspot areas. The highest temperatures are found in crowded and dense built-up areas, while the lowest temperatures are found in vegetated areas and water bodies. The findings also reveal a clear warming trend over the 30-year period, marked by a substantial decrease in areas with LST below 12 °C and between 12–17 °C, highlighting a shift towards warmer temperatures. Projected LULC changes indicate that urban areas will experience significant growth, increasing from 17.36% of the total area in 2023 to 21.39% in 2030, and further to 28.56% by 2050. The results suggest that the settlement land use class will increase by 11.2%, accompanied by a decrease in agricultural lands, vegetation, and water bodies.

1. Introduction

Over the past few decades, urban areas worldwide have experienced rapid growth [1] The key factors contributing to urbanization include economic development and population growth. It was estimated that global urban areas will expand by over a million square kilometers by 2030 [2,3,4,5,6,7]. According to the Census of India 2011, India’s urban population has significantly increased over the last two decades, rising from 217 million in 1991 to 377 million in 2011. The United Nations reports that India’s urban population is projected to grow by over 400 million between 2018 and 2050 [8]. Urbanization can positively impact people’s lives by improving their living standards and reducing energy consumption. However, the expansion of cities due to human activities has resulted in both positive and negative effects [6,9].

The rapid growth in the urban population and the expansion of urban areas have led to an increase in settlement areas and a decrease in green cover within cities. Urbanization results in the conversion of forests and agricultural land into settlement areas. Studies indicate that changes in land use and land cover (LULC) lead to the loss of agricultural and forest lands while increasing barren lands and impermeable surfaces [10,11]. The rate of LULC changes is more significant in cities than in other areas, mainly due to population growth and industrialization. Such trends could lead to irreversible changes in landscapes and pose threats to human well-being, natural resources, and ecosystem services [11,12].

The physical changes in the landscape, such as a decrease in green cover and an increase in settlement areas, likely increase the land surface temperature (LST) [13,14]. LST is defined as the temperature felt when the land surface is touched or as the skin temperature of the ground, and it is the measured in Kelvin [15]. The spatial variation in LST in urban areas is directly linked to both large- and small-scale urban features [16]. For example, forest cover, concrete buildings, roads, and parking lots exhibit specific thermal, radiative, and aerodynamic properties that are influenced by their surroundings. The reduction in agricultural lands and water bodies, coupled with an increase in LST due to the absorption of solar radiation, results in adverse effects on people’s daily lives due to climate change. LST is an indicator of the urban heat island (UHI) effect, which impacts sustainable development across various sectors [10,11].

Shillong City has undergone multiple transformations in its LULC characteristics, leading to significant impacts on LST. The emergence of concrete structures and increased urban development have resulted in warmer temperatures, particularly in densely populated areas. Therefore, it is essential to conduct a thorough investigation. While several studies have been conducted in the region, some have focused on LULC changes [17,18,19], while a few have emphasized urban heat island effects [20]. However, there remains a scarcity of research specifically exploring the prediction of LULC dynamics and their effects on surface temperature within the city.

The primary goal of this study is to examine the changes in LULC and their effects on LST in Shillong City, India. This study analyzed the effects of urbanization on Shillong City’s surface temperatures using a Markov chain model. Our main objectives are to: (1) analyze the LULC changes in the years 1993, 2003, 2013, and 2023; (2) analyze the changes in LST during the same years; (3) determine the relationship between LULC characteristics and LST; and (4) predict future LULC changes for the years 2030 and 2050.

2. Study Area

Shillong, the capital city of Meghalaya, is located in the East Khasi Hills at 25°35′N latitude and 91°53′E longitude, occupying the Shillong valley with an area of about 45 km² (Figure 1). It is situated at an altitude of 1520 m above mean sea level (msl). Shillong is an agglomeration of the six urban centers: Shillong Municipality, Shillong Cantonment, Madanriting, Mawlai, Nongthymmai, and Pynthorumkhrah. According to the 2001 Census, Shillong has a population of about 267,662 [21]. Shillong underwent a transformation in 1864 when it became the new civil station for the Khasi and Jaintia hills.

The climate of Shillong is generally very humid and is directly influenced by the southwest monsoon and the northeast winter winds [22]. Shillong experiences four main seasons: spring (March to April), monsoon (May to September), autumn (October to November), and winter (December to February). Temperatures begin to warm at the onset of spring, peaking during the summer months. Winters are quite severe, with temperatures dropping to as low as 4 °C. April and May are the warmest months, while January is the coldest month [23].

3. Material and Methods

3.1. Acquisition of Satellite Imageries

Landsat satellite data (USGS Collection 2, Level 2 archive) were downloaded from the United States Geological Survey (USGS) Earth Explorer (https://earthexplorer.usgs.gov/, accessed on 15 March 2024). The datasets include Landsat 5 Thematic Mapper (TM) for 1993, Landsat 7 (ETM+ SLC) Enhanced Thematic Mapper Plus (ETM+) and Scan Line Corrector (SLC) for 2003, Landsat 8 (OLI/TIRS) Operational Land Imager (OLI) for 2013, and Landsat 9 (OLI/TIRS) for 2023 (

Table 1. Details of the satellite imagery used in this study.

Satellite	Sensor Type	Spatial Resolution (m)	Year
Landsat 5	TM	30	1993
Landsat 7	ETM+	30	2003
Landsat 7	OLI/TIRS	30	2013
Landsat 9	OLI/TIRS	30	2023

). Additionally, field surveys and high-resolution images from Google Earth Pro 7.3.4 were used to determine the ground truth conditions.

3.2. Image Preprocessing and LULC Classification

An integrated flowchart depicts a series of preprocessing steps used to generate the results (Figure 2). First, we processed the datasets in Google Earth Engine (GEE) to create a false color composite (FCC). A georeferenced map of Shillong City was employed to extract and clip the study area from all the datasets. Landsat images for the years 1993, 2003, 2013, and 2023 were used to map the land use and land cover (LULC) characteristics. The supervised maximum likelihood classification (MLC) algorithm was used to classify different LULC classes. The MLC algorithm calculates the probability that each pixel belongs to a specific LULC class. The classification was performed in ERDAS Imagine 2014. All the satellite image pixels were categorized into five LULC classes: agricultural land, barren land, settlement, vegetation, and water bodies (

Table 2. Description of the various LULC classes considered in this study.

LULC Class	Description
Agricultural land	Agricultural lands include arable land (cropland and fallows), land under permanent crops, pastures, and hayfields.
Barren land	These include the area in which less than one third of the land has vegetation or other cover. Barren land has thin soil, sand, or rocks. Managed area with exposed rock and gravel pits.
Settlement	Settlement includes those areas covering urban areas, built-up areas, residential areas, roads, and commercial and educational buildings.
Vegetation	Vegetation includes dense as well as sparsely covered vegetated areas ranging from forests, dense bushes, urban green spaces, scrubs, and grasses.
Waterbody	Water bodies include wetlands, ponds, streams, rivers, and lakes. Lakes are the only feature detected in this study.

). Following classification, an accuracy assessment was conducted utilizing at least 200 stratified random sampling points for ground truth data. Upon achieving acceptable accuracy, the LULC maps of 1993, 2003, and 2013 were used to predict the LULC map for 2023 using a CA-Markov model. This predicted map was then validated against the classified LULC map for 2023, derived from the satellite imagery. After successful validation, the predicted maps for 2030 and 2050 were produced.

3.3. CA-Markov Prediction Model Analysis

The CA-Markov model uses a stochastic probability matrix to predict the transition from one state to another [6]. The model employs a conditional probability formula to estimate trend lines, as shown in Equations (1)–(3):

$$S \left(t+1\right)=P_{ij}\times S\left(t\right)$$

(1)

$$P_{ij}=\left[\begin{array}{ccc}{P}_{11}& P_{12}& P_{1n}\\ P_{21}& P_{22}& P_{2n}\\ P_{n1}& P_{n2}& P_{n3}\end{array}\right]$$

(2)

However,

$$\left(0\le P_{ij}<1 and \sum _{j=1}^N{P}_{ij}=1,(i,j=\mathrm{1,2},3\dots ..n)\right)$$

(3)

A Markov chain and cellular automata modelling techniques are used to predict LULC changes for the years 2030 and 2050. These predictions are calculated using Terrset’s (Clark Labs TerrSet 18.31, Worcester, MA 01610, USA) land use change modeler (LCM). The accuracy of predicting LULC can be enhanced using the CA-Markov model, provided previous LULC characteristics remain consistent. However, the CA-Markov model may not yield precise spatial predictions for raster datasets [24]. Influential factors between the CA-Markov model and other variables can be directly identified, as the model is based on a probability matrix [25]. Different variables play important roles in identifying the most crucial factors. The CA-Markov model creates training patterns and initiates training automatically when it receives inputs from the strata. The input parameters are not given specific weights based on established standards. Changes such as urbanization, deforestation, and rising temperatures are mainly driven by human actions and decisions at regional and metropolitan levels, making accurate prediction difficult [6]. Dynamic models have limitations, but they remain useful in creating assumptions and making decisions about LULC changes in any given area, despite lacking strict rules.

3.4. Retrieving Land Surface Temperature

Google Earth Engine offers various Landsat-specific processing methods, including computation of at-sensor radiance, top-of-atmosphere (TOA) reflectance, surface reflectance (SR), cloud score, and cloud-free composites [6]. In this study, we used the SR product, which measures the fraction of incoming solar radiation reflected by the Earth’s surface. The Landsat SR product, available in GEE, is a replica of the USGS Collection 2, Level 2 archive. Level 2 scientific products are atmospherically corrected images (https://www.usgs.gov/landsat-missions/landsat-collection-2-level-2-science-products, accessed on 15 March 2024). The SR product of Landsat 5 and 7 are generated using the Landsat ecosystem disturbance adaptive processing system (LEDAPS) version 3.4.0; the algorithm was originally developed by NASA [26]. For Landsat 8 and 9, the SR product is generated using the Land surface reflectance code (LaSRC) version 1.5.0 algorithm [27].

3.5. Calculation of Land Surface Temperature

Thermal satellite sensors can provide surface temperature and emissivity information. The land surface temperature is derived from measured thermal infrared (TIR) radiation. For Landsat 5 TM images and Landsat 7 (ETM+) images, band 6 is used for surface temperature (ST), and for Landsat 9 (OLI) images, band 10 is used to convert the raw image into spectral radiance. The land surface temperature (LST) is calculated using the following formula:

$$LST={\displaystyle \frac{TB}{1+\left(\frac{\mathsf{\sigma }\times \mathsf{\lambda }\times \mathrm{T}\mathrm{B}}{h\times c}\right)\mathrm{l}\mathrm{n}(\mathsf{ϵ})}}$$

(4)

where TB is the brightness temperature in Kelvin (K), ∈ is the land surface emissivity, λ is the effective wavelength, c is the speed of light, σ is the Boltzmann constant, h is the Plank constant.

4. Results

4.1. Accuracy Assessment

The LULC classification for urban growth is sensitive to accuracy assessments. A method was defined to assess the accuracy of the LULC classification, utilizing 200 stratified random samples for each LULC map. The results show that the overall accuracy was 85% for 1993, 86% for 2003, 90% for 2013, and 89% for 2023. Additionally, the Kappa coefficients for the LULC maps were 0.76, 0.77, 0.85, and 0.81, respectively.

To ensure compatibility between the classification and reference data, Kappa coefficients should generally be greater than 0.75. The USGS recommends the use of Landsat satellite images for LULC mapping if the accuracy level is 85% [28]. Our accuracy assessment results align with these recommended values.

4.2. Transition of LULC Characteristics over the Years

The LULC classification results for 1993, 2003, 2013, and 2023 showed notable changes in urban development and other land use patterns (Figure 3 and

Table 3. The area of different LULC classes (in square kilometers, km²) was measured after classifying satellite images across different years.

LULC Class	1993		2003		2013		2023
	Area (km2)	Area (%)	Area (km2)	Area (%)	Area (km2)	Area (%)	Area (km2)	Area (%)
Agricultural land	12.37	7.24	14.13	8.27	24.36	14.26	32.08	18.79
Barren land	69.87	40.87	49.63	29.03	36.49	21.36	11.65	6.82
Settlement	10.94	6.40	15.34	8.97	17.92	10.49	29.63	17.36
Vegetation	77.67	45.43	91.69	53.64	91.78	53.73	97.22	56.95
Waterbody	0.10	0.06	0.15	0.09	0.27	0.16	0.13	0.08

). The settlement land use class has experienced a steady increase, rising from 6.40% in 1993 to 8.97% in 2003, further to 10.49% in 2013, and finally reaching 17.36% in 2023. Conversely, barren land decreased from 40.87% in 1993 to 29.03% in 2003, continued to decrease to 21.36% in 2013, and further to 6.82% in 2023. Interestingly, agricultural land experienced a slight growth, increasing from 7.24% in 1993 to 8.27% in 2003, then to 14.26% in 2013, and finally to 18.79% in 2023. An analysis of land use change from 1993 to 2023 revealed a transformation of barren land into urban areas. During the same period, a notable increase was observed in built-up areas (settlements) by 10.96%, agricultural land by 11.56%, and vegetation by 11.52%.

To understand the changing patterns of LULC classes, we analyzed the net changes in various LULC classes that transitioned into settlement areas from 1993 to 2023 (Figure 4). The settlement area is the most stable land use class, suggesting that the likelihood of settlement areas being transformed into different land use classes is minimal. Additionally, forested areas are maintaining some level of stability, likely due to recent plantation initiatives undertaken by the government. The analysis revealed that the settlement area exhibited positive transformations over the studied periods. The results from all LULC classified maps show that settlement areas have transitioned from agricultural land, barren land, and vegetation (Figure 4). However, urban activities can lead to land degradation, making it difficult for vegetation to grow and turning previously built-up areas into barren land. The swift urban growth over the last three decades can be attributed to the continuous influx of a sizable population to the city, driven by factors such as the availability of fertile lands, business and employment opportunities, better urban infrastructure, healthcare facilities, and education. However, the rapid development of built-up areas will lead to sharp declines in forest cover and agricultural lands, potentially jeopardizing ecosystem health, human well-being, and food security. Additionally, the conversion of forests could expose bare ground, making it vulnerable to gully erosion in such a fragile and rugged terrain [29]. Population growth could lead to increased demand for groundwater and other natural resources, potentially resulting in groundwater pollution, which threatens the availability of groundwater resources and causes health problems [30].

4.3. Changes in LST between 1993 and 2023

The LST values have been divided into five categories (Figure 5). These categories were used to estimate the area covered by each temperature range. The results indicate a notable increase in surface temperatures between 1993 and 2023. Variations in temperature were observed throughout the studied period, for example, 5.46 °C to 26.28 °C (average of 13.85 °C) in 1993, 9.03 °C to 33.52 °C (average of 18.04 °C) in 2003, 11.09 °C to 31.62 °C (average of 19.03 °C) in 2013, and 11.76 °C to 29.68 °C (average of 19.24 °C) in 2023 (

Table 4. Land surface temperature (LST) was determined by classifying satellite images across different years.

Year	Minimum (°C)	Maximum (°C)	Mean (°C)
1993	5.46	26.28	13.85
2003	9.03	33.52	18.04
2013	11.09	31.62	19.03
2023	11.76	29.68	19.24

).

During the study period, a noticeable shift from lower temperature areas to higher temperatures areas was observed. In 1993, temperatures below 12 °C accounted for 27.95% of the total area. This area significantly decreased to 1.15% by 2003, further reduced to 0.12% by 2013, and became negligible by 2023. The temperature between 12 and 17 °C accounted for 56.25% of the total area in 1993. This area decreased to 36.25% by 2003, further decreased to 25.77% by 2013, and became 18.22% by 2023. For the temperature range of 17–23 °C, the total area covered was 15.68% in 1993, which increased to 56.94% in 2003, 66.22% in 2013, and finally increased to 76.83% in 2023. For the temperature range of 23–27 °C, the total area was 0.12% in 1993, which increased to 5.29% in 2003, 7.45% in 2013, before slightly decreasing to 4.83% in 2023. For the temperature range >27 °C, no area was recorded in 1993; however, the area increased to 0.37% in 2003, 0.44% in 2013, and then decreased to 0.11% in 2023. The findings reveal a clear and consistent warming trend over the 30-year period. The substantial decrease in areas with LST below 12 °C and between 12–17 °C, coupled with the significant increase in areas with temperatures between 17 and 23 °C, underscores a shift toward warmer surface temperatures. The initial rise and subsequent decrease in the 23–27 °C and >27 °C temperature ranges indicate variability in extreme temperatures. However, further research is necessary to understand the underlying causes and potential long-term consequences of these temperature shifts.

4.4. Variations in LST across Different LULC Classes over the Years

Urbanization is one of the most dynamic forces driving changes in LULC characteristics and has a significant impact on the environment. Built-up (settlement) areas have the highest average temperatures, along with barren land, and agricultural land (Figure 6). Since LST was measured during the winter, most of the agricultural lands were barren. Water bodies and vegetation exhibit lower LST, ranging between 12.20 °C and 16.67 °C during the study years. In 1993, built-up areas had an average LST of 17.45 °C, which rose to 21.56 °C in 2003, slightly decreased to 21.37 °C in 2013, and then increased to 21.75 °C in 2023. Overall, urban areas experienced a rise in mean LST by 4.3 °C. This increase is accompanied by a higher standard deviation and broader minimum and maximum ranges, indicating dynamic variation in land use characteristics. This variation is particularly due to the mixture of different LULC classes resulting from the development of built-up areas, industrial zones, and road networks. Similarly, agricultural land and barren land show an increase in mean LST of 4.5 °C and 4.9 °C, respectively. The data indicate that the factors influencing LST for each LULC class are distinct.

A significant rise in LST is associated with the development of buildings and urban infrastructure. Conversely, the smallest rise in LST is observed in areas with vegetative cover. This difference is primarily attributed to the high evapotranspiration from vegetation, which has a cooling effect [16]. Existing research has demonstrated that green cover tends to lower temperatures by absorbing and reflecting solar radiation and regulating latent and sensible heat exchange [31]. On the other hand, impervious surfaces and urban structures, such as concrete, asphalt, and metal, play a significant role in raising temperatures [16]. The spatiotemporal variation of LST can influence the urban heat island effect by increasing surface and air temperatures within the city [31].

4.5. Modelling and Prediction of Future Land Uses

We utilized a CA-Markov model within the LCM module of TerrSet 2020 geospatial monitoring and modeling software to simulate LULC patterns for the years 2030 and 2050. This simulation is based on historical LULC data and transition matrices specific to the study area. To validate the simulation model, classified LULC maps of the study area for the years 2013 and 2023 were used to generate a probability matrix, transition area matrix, and a series of conditional probability images. These datasets facilitated the generation of predicted land use scenarios for 2030 and 2050 [32].

To validate the accuracy of the prediction, the CA-Markov model was used to estimate the areas of different LULC classes for 2023. The predicted and classified maps were then compared and the prediction accuracy was calculated (

Table 5. Comparison of predicted and classified areas of various LULC classes for the year 2023.

LULC Class	Area Classified (km2)	Area Predicted (km2)	Accuracy (%)
Agricultural land	32.08	36.18	87.22
Barren land	11.65	23.95	5.58
Settlement	29.63	21.95	74.08
Vegetation	97.22	88.16	90.68
Waterbody	0.13	0.25	7.69

). The data reveal varying degrees of accuracy in the predicted and classified maps across different LULC classes (Table 5). For agricultural land, the CA-Markov model demonstrated high accuracy at 87.22%, with the classified area of 32.08 km² closely aligning with the predicted area of 36.18 km². The prediction of vegetation was also highly accurate, with an accuracy of 90.68%, indicating that the model effectively captured the changes. Settlement areas were predicted with an accuracy of 74.08%, though the model underestimated the total area by 7.68 km². In contrast, the model’s performance for barren land and water bodies was largely inaccurate, with accuracies of 5.58% and 7.69%, respectively, highlighting the challenges in accurately predicting small and dynamic land uses.

Since the transition sub-model focuses on settlements, we prioritized the accuracy of settlement predictions. The years 2003 and 2013 were used to predict changes in 2023. The results showed that the rate of change in settlements from 2003 to 2013 gradually increased, whereas from 2013 to 2023, there was an accelerated increase in settlements. Other parameters, such as roads, elevation, and the trend map of transitions from all classes converted into settlements, were added to the transition model to confirm its reliability. The model was then used to predict LULC for 2030 and 2050.

4.6. Predicted LULC for 2030 and 2050

The predicted LULC maps indicate that urban areas will experience significant growth, with the total area increasing to 21.39% in 2030 and further to 28.56% in 2050 (

Table 6. Predicted areas of various LULC classes for the years 2030 and 2050.

LULC Class	2030		2050
	$\mathbf{Area} ({\mathbf{k}\mathbf{m}}^2$)	Area (%)	$\mathbf{Area} ({\mathbf{k}\mathbf{m}}^2$)	Area (%)
Agricultural land	31.57	18.48	27.01	15.81
Barren land	8.21	4.81	6.59	3.86
Settlement	36.54	21.39	48.79	28.56
Vegetation	94.42	55.26	88.38	51.72
Waterbody	0.12	0.07	0.10	0.06

). This increase is likely if the current development trend continues without intervention. According to the LULC scenario, settlements will expand by 4.03% by 2030 and 11.2% by 2050, accompanied by a decrease in vegetation, agricultural land, and barren land. These results highlight the potential impacts on ecosystem services, urban health, and thermal characteristics due to deforestation and increased urbanization. If uncontrolled urban growth persists, it will likely lead to considerable environmental, economic, and health issues.

To promote environmental sustainability, it is crucial to implement a well-designed land use plan, safeguard water bodies, and actively engage in reforestation efforts. These measures will contribute to a more sustainable future for the city. The projected spatial changes in LULC for 2030 and 2050 indicate that forest cover will undergo fragmentation (Figure 7). Smaller and narrower patches of forests will be scattered throughout the city, particularly along streams, rivers, and marshlands. Much of the barren and forested land from 2023 will transform into built-up areas. The presence of bare ground suggests possible soil excavations, likely resulting from construction activities and landfills. These findings align with the ongoing trend of urbanization and current government policies.

5. Discussion

Shillong has experienced rapid population growth due to its administrative and educational facilities. The city’s expansion had reached its peak by 2010, leading to the development of new urban areas known as New Shillong. Despite this expansion, other LULC categories have remained relatively stable over the last thirty years. The forest area, in particular, has not undergone significant changes during the study period. Shillong is a tribal town, and strict regulations prohibit indiscriminate tree-cutting, which has helped maintain a relatively constant forest area. In 1970, the forest covered 31.88 km², which increased to 48.38 km² by 1991, and further grew to 54.68 km², remaining constant until 2010, with only a marginal decrease of 0.05 km² [17]. However, the predicted changes in LULC characteristics could exacerbate UHI dynamics, as the development of more urban areas is likely to increase LST in the future. Adnan et al. (2024) [33] studied heatwave vulnerability in five major cities in Bangladesh, using MODIS data to capture spatial and temporal variations in heat exposure. They identified populations highly sensitive to heatwaves. Therefore, it is important to consider social dynamics when assessing the impact LULC and LST variations on urban populations.

The rapid growth of built-up areas negatively affects the natural cooling process, influenced by factors such as decreased vegetation shading and transpiration [6,16]. Studies using the normalized difference vegetation index (NDVI) and LST have demonstrated that vegetation cover acts as a cooling sink within urban settings [34,35]. Urbanization modifies construction practices, leading to reduced soil infiltration and increased surface runoff. As a result, groundwater levels decline [30]. Climate change further exacerbates the deterioration of the water balance. Changes in land use significantly impact climate variables, such as daily maximum and minimum temperatures [36]. Additionally, alterations in surface albedo caused by land use changes disrupt the Earth’s radiation balance [37]. Various physical factors, including topography, land use, and vegetation, influence LST in urban areas [12]. The arrangement of vegetation, built-up areas, open land, and water bodies plays a crucial role in determining surface temperature. Vegetation and surface water bodies act as temperature regulators, making regions with high vegetation density and water features relatively cooler compared to areas with lower vegetation and fewer water features.

Shillong, a prominent hill station in the eastern Himalayas, has experienced rapid population growth in recent years. This study offers valuable insights into the environmental consequences of urbanization in Shillong, providing a foundation for future urban planning and climate adaptation strategies in the region and other areas with similar geographical characteristics.

6. Study Limitations

Kamusoko et al. (2009) [38] highlighted that a key limitation of CA-Markov analysis is its inability to account for human activities and government policies that impact land use behavior, as well as the lack of high-resolution imagery, which reduces the model’s accuracy. Zhou et al. (2012) [39] noted that the absence of socio-economic data has led to inaccuracies in their study on land salinization. Yagoub & Al Bizreh (2014) [40] pointed out that CA-Markov models often assume that factors driving land use change will remain constant over time, which is unrealistic and leads to errors in simulations. Yang et al. (2014) [41] also noted limitations in their study related to the methodology used for subdividing land into smaller regions. During the subdivision process, inconsistencies can arise when the same land use classes transition into different ones, potentially leading to errors. These issues represent an area for future research. It is recommended to use modern classification algorithms, such as recurrent neural networks (RNN) and long short-term memory (LSTM) in deep learning. These methods can effectively capture complex patterns and relationships in data, leading to improved predictions [42].

7. Conclusions

This study examines how LULC change affects LST in Shillong City, India. The results showed that regions with natural cover, such as water bodies and vegetation, had the lowest surface temperatures, with mean values of 16.57 °C and 16.21 °C, respectively. In contrast, barren lands had average surface temperatures of 20.01 °C, agricultural land 21.05 °C, and built-up areas 21.75 °C during 2023. Future LULC projections indicate significant growth in built-up areas, expected to reach approximately 36.54 km² by 2030 and 48.79 km² by 2050, primarily concentrating around the city center and expanding towards the northern and southwestern regions. The integration of satellite data and GIS technology is crucial for obtaining essential spatial information and LULC characteristics necessary for effective land resource planning and management. This study demonstrates the effectiveness of remote sensing data and GIS technology in analyzing LULC changes and predicting future scenarios. However, further research is needed to incorporate various environmental and socio-economic factors into the simulation, including the use of modern classification algorithms such as RNN and LSTM. Urban development and planning should focus on mitigating climate change impacts by creating green spaces, parks, and preserving water bodies. Promoting ecological resources and environmental education is crucial for effective urban planning and implementing green policies to tackle rising thermal stress. A detailed quantitative analysis of these aspects is required to highlight the direct influence of urbanization on LST and any lagging effects, which may provide valuable insights for urban planning.