With a population estimated to reach 222 million by 2050 (Streatfield et al., 2008; Farid et al., 2011), Bangladesh is a developing country small in size with a large population that has increased the threat to the nations vegetation cover (Mukul, 2016; Hasnat et al., 2018; Singh et al., 2020). The nations current land use patterns have long been brought to light by the countrys accelerated population boom (Biswas and Choudhury, 2007; Reddy et al., 2016). Furthermore, in the twenty-first century, there has been a rapid shift in land uses and land cover change (LULC) due to economic and industrial development, high population rates, and other factors. This has resulted in the conversion of landscapes into residential areas and resulting in increase of impervious surface area globally. Buildings, roads, and industrial zones are examples of impermeable surfaces that absorb solar radiation in the short wavelength range and decrease solar radiation in the long wavelength range that the earth emits (Adiguzel et al., 2022a and b; Adiguzel et al., 2020). One well-known technique to evaluate ecological and environmental deterioration is to identify patterns of changes in land cover and land use (Beevi et al., 2015; Hadeel et al., 2011; Giri et al., 2005). For the responsible and lasting use of natural resources, environmental protection, and food security, it is crucial to comprehend the complexity of LULC changes and to analyze and monitor them (Drummond et al., 2012; Foley et al., 2005). Research on LULC changes can also be useful in forecasting future trends and informing planning decisions for natural resource management (Prenzel, 2004; Bekere et al., 2023).

There are two standard methods for identifying the UHI effect. According to Hejazizadeh et al. (2019) and Monteiro et al. (2016), the first method uses ground-based air temperature measurements in microstudies based on the modeling of meteorological data. When conducting macrostudies using land surface temperature (LST) data, the second approach is employed. LST is determined using remote sensing techniques that track the thermal energy that the earth releases into the atmosphere using data from satellite images. Since it allows us to measure the energy emitted from the Earth into the atmosphere, LST measurements using thermal satellite photos are simple, quick, continuous, and highly accurate (Matzarakis, 2002; Cetin, 2020a; Cetin, 2019; Matzarakis, 2007).  The primary factors influencing LST values are changes in LULC and greenhouse gas emissions. Because industrialization mostly relied on fossil fuels, it increased emissions of greenhouse gases, Leading to a rise in greenhouse gas concentrations in the troposphere, the lowest layer of the atmosphere. The earths surface temperature rises due to absorption and reemission of solar radiation in the troposphere (Zeren Cetin and Sevik, 2020; Zeren Cetin et al., 2020; Zhong and Chen, 2019; Cetin, 2020a; b; Cetin, 2019). The LULC variations affect the radiation emissions from the planet surface. Radiation is reflected by impermeable surfaces like highways and buildings, but is absorbed primarily by vegetation and surfaces coated in permeable materials. As a result, LULC changes alter LST values by altering the land surfaces rate of radiation absorption.  Furthermore, a number of studies have highlighted how evaporation from dense, healthy vegetation influences the LST and UHI. Thus, vegetation indicators and vegetation evaluation are done using remote sensing techniques. To characterize vegetation patterns, one of those indices the normalized difference vegetation index, or NDVI is frequently employed in research (De Freitas, 2003; Lise and Tol, 2002; Cetin, 2020a and b; Cetin, 2019; Lin and Matzarakis, 2008). Since the 1870s, Bangladeshs forest area has decreased, covering under 16% of the countrys total land, equivalent to 2.33 million hectares (Mukul et al., 2016; Nesha et al., 2021).

On the contrary, a lot of researchers looked at how urban impervious surfaces, such building roofs and roadways, affected LST using the normalized difference built-up index (NDBI). Furthermore, research on LST and UHI demonstrates that identifying the variables influencing LST and the ways in which these variables interact are critical to developing strategies that mitigate the urban heat island effect (Matzarakis, 2006; Zhong and Chen, 2019).

Such LULC pattern analysis may be quickly completed and used to display the affected areas and their effects on the surrounding environment by integrating Geographic Information Science (GIS) and Remote Sensing (RS) techniques. Numerous studies have been conducted on land use and land cover (LULC) changes using different satellite data, such as Landsat, MODIS, and SPOT (Mondal et al., 2021, 2022; Thakur et al., 2020; Thakur et al., 2020, 2020). However, Landsat satellite imagery is particularly valuable for LULC change detection due to its freely accessible, moderate-resolution data available in multitemporal time series dating back to 1972 (Lu et al., 2019). According to various sources (Mondal & Bandyopadhyay, 2016, 2022; Chamling and Bera, 2020; Lai, 2020; Mondal et al., 2019; Mondal et al., 2016), analyzing LULC change dynamics is logical for evaluating the environmental transformations occurring in a given area and aids in the development of an efficient management plan which can assist in achieving sustainable development and appease both local and global environmental changes. Given all of these advancements, it is crucial to look into the implications of variations in the LST, NDVI, and NDBI values in Coxs Bazars central district. Thus, in order to assess the change in the UHI effect, this study first examined the effects of changes in LULC, LST, NDVI, and NDBI values in the Coxs Bazar district between 1990 and 2020.

Objectives

1. To calculate the Land Surface Temperature (LST), Normalized Difference Vegetation Index (NDVI), and Normalized Difference Built-up Index (NDBI) for the study area in the years 1990, 2000, and 2020.
2. To evaluate the changes in land use and land cover over the past three decades.

Study area

The Coxs Bazar district, covering an area of 2,121 square kilometers (158,000 hectares), had a population of 2,906,281 as of the 2022 census. The total district area is 2,492 square kilometers, with a population density of 1,166 individuals per square kilometer in 2022. The annual population growth rate from 2011 to 2022 is 1.8%. The rural population stands at 1,231,639, while the urban population is 1,591,629, as reported by the Bangladesh Bureau of Statistics (BBS, 2021).

Fig. 1: Location of study area.

Geographically, Coxs Bazar is bordered by the Chattogram district to the north, the Bay of Bengal to the south and west, and the Bandarban district, Arakan (Myanmar), and the Naf River to the east. It is renowned for having the longest sea beach in the world. The climate, typical of the tropical monsoon region, features high temperatures, substantial rainfall, excessive humidity, and clear seasonal variations. Coxs Bazars coastal location influences its climate, with average annual temperatures ranging from a high of approximately 34.8°C to a low of 16.1°C, and an average annual rainfall of 4,285 mm .

The transformation of land use in Coxs Bazar over recent decades have been driven by several key factors. Rapid population growth and urbanization are the main reasons of them. They have led to the expansion of residential, commercial, and infrastructural developments (Hossain et al., 2019). The influx of Rohingya refugees since 2017 has further intensified land use pressures, with large areas of forest and agricultural land cleared for settlements and supporting infrastructure (UNHCR, 2018). Also, the tourism industry centered around the districts famous sea beach, has significantly increased, prompting the construction of hotels, resorts, and other tourist facilities. Additionally, agricultural practices have shifted, with some areas being converted to accommodate these urban developments. These changes have resulted in a decrease in vegetation cover and agricultural land, contributing to environmental degradation and altering the natural landscape (Rahman & Islam, 2016).

Data collection 

This study uses secondary data for analysis where multi-spectral satellite imagery from three different years (1990, 2000, and 2020) was downloaded from the United States Geological Survey (USGS) as the secondary source to detect the changes in LULC, LST, NDVI and NDBI. Landsat images with less than 5% cloud coverage and 30 m resolution were acquired from USGS for the year 1990 (Landsat 5 Thematic Mapper), 2000 (Landsat 5 Thematic Mapper), 2020 (Landsat 8 Operational Land Imaget_TIRS). To detect land use changes, supervised classification method (by ArcGIS 10.2) was used. 

Table 1: The characteristics of the images.

LULC

A false color composition was used to classify the land uses into four classes e.g. Agriculture, Vegetation, Water, and Aquaculture. Only pre-monsoon (March to May) sunmmer seasons images were included in this study due to the because of the occurrence of cloud-free sky and least rainfall in this period that aided the actual detection of classified areas.  The data were inbuilt georeferenced to UTM zone 46 North projections with WGS-84 datum. Before the analysis, radiometric corrections and image enhancement procedures were executed with Arc-GIS 10.8 software. LULC maps were classified using the supervised classification method.  

Fig. 2: The schematic framework of methodology.

The supervised classification was processed by ArcGIS 10.8. Bands 1–5 and 7 were used for preparing LULC maps with Landsat 5–TM images, but the 6 was excluded because it was a thermal band. The bands 1–7 were considered for Landsat 8-OLI images. All the bands were compounded in ArcGIS 10.8 software using the image analyst tool. More than 200 training samples were collected randomly using the training sample manager tool for each year.

NDVI

NDVI quantifies vegetation cover by calculating the difference between the near-infrared and red portions of the electromagnetic spectrum. Healthy and lush vegetation absorbs visible light and reflects most of the nearinfrared light. In contrast, unhealthy vegetation reflects more visible light and less near-infrared light. NDVI has values ranging from + 1 to − 1 (+ 1 for the healthy vegetation and − 1 for areas with no vegetation). NDVI was calculated in ArcGIS 10.8 by using the following equation –

NDBI

The Normalized Difference Built-up Index (NDBI) ranges from −1 to +1 and is used to distinguish built-up areas from other land uses or land surfaces. Built-up regions typically exhibit higher reflectance in the mid-infrared (MIR) band compared to the near-infrared (NIR) band. Consequently, high NDBI values indicate built-up zones. NDBI was calculated using ArcGIS 10.8 software, employing MIR band 6 and NIR band 5 from Landsat-8 data, and MIR band 5 and NIR band 4 from Landsat-5 data, as represented in the equation below.

LST

LST refers to the temperature measured by the remote sensor. As the LST provides essential data about a climate system, the procedure has been addressed and clarified in many studies. In this study, Landsat 5 images from 1990 and 2000, along with Landsat 8 images from 2020, were used to conduct the Land Surface Temperature (LST) analysis as outlined below:

Step 1: Conversion to Top of Atmosphere (TOA) Radiance - The formula in Equation 3 was applied to convert thermal infrared digital numbers (DN) into TOA radiance values for Landsat 8 data.

 ML = the band-specific multiplicative rescaling factor = 0.0003342; F = re-scaling factor = 0.1; Oj = correction value = 0.29. 

The formula in the equation below was applied to the Landsat-5 TM data.

 LMAXλ = is the radiance that is scaled to QCALMAX; LMINλ is the radiance that is scaled to QCALMIN, QCALMIN is the lowest calibrated value, and LMIN is the highest calibrated value LMAX λ

Step 2: Brightness temperature: radiance values were converted to brightness temperature using the equation below

 T = brightness temperature (0 C); K1 = calibration constant 1 (Table 2); K2 = calibration constant 2 (Table 2)

Step 3: Calculation of LST applying the equation

  

Table 2: Thermal constant of Landsat TM and OLI thermal imageries.

T = brightness temperature (°C), λ = wavelength of emitted radiance = 10.8 for band, E = emissivity of land surface, C = 14,388mK

Obtaining the Land Surface Emissivity (E) is essential for calculating LST. Emissivity represents a surfaces ability to emit and absorb radiation in the long-wave spectrum and typically varies depending on the type of land cover.

E = land surface emissivity; Pv = proportion of vegetation; 0.986 is the correction value

To determine the land surface emissivity (E) value, the vegetation proportion (Pv) parameter equation proposed is shown below

Pv is vegetation proportion; NDVImin is the lowest value of NDVI; NDVImax is the highest value of NDVI

Accuracy assessment

The classification accuracy assessment was conducted with the reference of the raw satellite images. The entire process was executed by comparing the reference images with the classified images with some random points following stratified random sampling procedure. Another widely used metric for assessing image classification accuracy is the Kappa coefficient (K). Generally, its values range between 0 and 1, with higher values indicating greater agreement and accuracy (Billah et al., 2021; Islam et al., 2021; Hasnat, 2021). In this study, the Kappa statistic was also calculated using the formula below:

Where TP is total pixels, and TCP is total corrected pixels. The result reveals that the value of overall classification accuracy varies from 79.22% to 90.79% and Kappa statistics varies from 0.76 to 0.87 that indicates high accuracy.

LULC

The land use of Coxs Bazar district was categorized into four: waterbody, vegetation, bare lands, aand built-up area. Fig. 2 shows the spatial distribution of land use in Bartin in 1990 and 2020. Table 3 below presents the LULC and change statistics. A minus (−) sign indicates a decrease compared to the previous period, while a plus (+) sign denotes an increase. The analysis results showed that between 1990 and 2020 , built-up areas increased significantly from 467.1 to 945.7 km2 by 22.6%, but bare lands have decreased from 842.64 (38.80%) to 148.44 (6.86%) by 694.2 (-31.90%). Here we can also see that the vegetation cover has increased a little from 842.16 km2 (38.78%) to 1067.70 km2 (49.17%) by 225.54 km2 (+10.39%). Land use maps provide essential information about the spatial change in urban regions.

Table 3: LULC and Change Statistics of Coxs Bazar.

LST

Fig. 3 shows the LST maps of the study area in 1990, 2000 and 2020. As seen in the figure, the average temperature in most of the Coxs Bazar district was 20–25 °C in the 1990s, while some part were less than 20°C; in 2000 we see whole of the area having avarage LST of 20-25 °C. But it rose to 25-30°C in this region in 2020, while some part having LST of more than 30°C. Similarly, there is a significant increase in urban settlements LST values. Table 4 shows a statistical overview of LST values. 

Fig. 3: Spatial distribution of LULC in 1990 and 2020.

In general, there has been an increase in the maximum and average temperature values of the study area over the 30 years. The average temperature rose from 20.28 to 25.19 °C from 1990 to 2020.

NDVI

The normalized difference vegetation index (NDVI) identifies and monitors vegetation. Several studies have also examined the effects of NDVI on the urban heat island (UHI) phenomenon (Zhong & Chen, 2019; Fletcher & Morakabati, 2008). NDVI values between 0.1 and 0.75 typically indicate vegetation cover. NDVI maps are presented in Fig. 4, where the lowest values correspond to bare land and built-up areas.

Fig. 4: LST maps for the years 1990, 2000 and 2020.

 However, as highlighted by Lin and Matzarakis, (2008) and Scott and Lemieux, (2009) it is difficult to distinguish between bare lands and built-up areas solely based on NDVI values. NDVI values above 0.1 generally indicate vegetation cover. In this study, NDVI values were categorized into low (0.1–0.2), medium (0.2–0.3), and high (0.4 and above) classes. 

Table 4: The statistical summary of LST values for Coxs Bazar.

As shown in Fig. 4, the most significant changes occurred in three main regions. The first was Bartin city center, where NDVI values decreased substantially. Changes were also noticeable in agricultural lands to the north and upland settlements to the south. In particular, the upland settlements in the mountainous areas, which typically have rich vegetation, contributed to the overall decline in vegetation cover across the study area. 

Table 5 shows the NDVI density classes for 1990 and 2021, while Table 6 provides a statistical summary of NDVI values. As indicated in Table 5, decreases occurred in the low, medium, and high-density vegetation classes, while there was a marked increase in built-up areas.

 Fig. 5: NDVI maps for the years 1990, 2000 and 2020.

NDBI

NDBI has been applied in many studies as it helps distinguish the built-up areas from other land uses and land cover. The NDBI density maps are given in Fig. 5 and statistical data are shown in Table 7. 

Table 6: NDVI and Change Statistics for Coxs Bazar.

Negative NDBI values correspond to vegetation, small positive values indicate bare land, and larger positive values represent built-up areas. However, the NDBI values associated with residential areas have varied across different studies. However, the NDBI values indicating residential differed in some studies (Zhong and Chen, 2019; De Freitas, 2003). NDBI value representing the residential areas was 0.07 for the year 1990 and, 0.16 for the year 2021. 

In this study, LST values of 1990, 2000 and 2020 were analyzed in order to determine the change in the urban heat effect in Coxs Bazar over 30 years. The study results showed that the average temperature of Coxs Bazar was 20.28 °C in 1990, but it rose to 25.19 °C in 30 years. In order to determine the spatial changes in the LST values, LULCC maps for 1990 and 2021 were formed. As a result of the LULCC maps, it was found that the land use of Coxs Bazar district had changed rapidly from 1990 to 2020 due to urbanization and agricultural and industrial practices. A negative correlation also found between NDBI and LST. The built areas has increased in the last 30 years, but the vegetation has decreased resulting in the rise of LST.

The data used in this study, including Geographic Information Science (GIS) and Remote Sensing (RS) datasets, are available upon request from the corres-ponding author. Researchers interested in accessing the data are encouraged to contact the corresponding author for further information and assistance.

I would like to express my sincere gratitude to Honorable Maam Meher Neegar Neema for her invaluable guidance and continuous support throughout this work. Her insightful advice and encouragement have been instrumental in the successful completion of this study. I am deeply thankful for her time and dedication.

The authors declare that they have no known financial conflicts of interest or personal relationships that could have influenced the work reported in this paper.