Examining the relationship between land surface temperature and landscape features using spectral indices with Google Earth Engine.

Land surface temperature (LST) is strongly influenced by landscape features as they change the thermal characteristics of the surface greatly. Normalized Difference Vegetation Index (NDVI), Normalized Difference Water Index (NDWI), Normalized Difference Built-up Index (NDBI), and Normalized Difference Bareness Index (NDBAI) correspond to vegetation cover, water bodies, impervious build-ups, and bare lands, respectively. These indices were utilized to demonstrate the relationship between multiple landscape features and LST using the spectral indices derived from images of Landsat 5 Thematic Mapper (TM), and Landsat 8 Operational Land Imager (OLI) of Sylhet Sadar Upazila (2000-2018). Google Earth Engine (GEE) cloud computing platform was used to filter, process, and analyze trends with logistic regression. LST and other spectral indices were calculated. Changes in LST (2000-2018) range from -6 °C to +4 °C in the study area. Because of higher vegetation cover and reserve forest, the north-eastern part of the study region had the greatest variations in LST. The spectral indices corresponding to landscape features have a considerable explanatory capacity for describing LST scenarios. The correlation of these indices with LST ranges from -0.52 (NDBI) to +0.57 (NDVI).

Introduction

Monitoring land use and land cover (LULC) dynamics of an area at various temporal scales allow for the assessment of landscape dynamics and environmental health [1]. The dynamic change of items or continuous events in a specific region across time is described by spatiotemporal data [2, 3]. Urbanization is rapidly changing the landscape all across the world by covering up larger surface areas [4]. Dynamic LULC patterns have a variety of effects on the thermal environment, specifically changes in land surface temperature (LST) due to the build-up of impervious surfaces [5, 6]. Land Surface Temperature (LST) is one of the most important aspects of surface energy and water balance on local to global scales and is also considered a key metric when analyzing the exchange of comprised material, energy balance, and biophysical and chemical processes of the land surface [7, 8, 9, 10, 11, 12, 13, 14].

Satellite imagery and data obtained by Remote Sensing (RS) and Geographic Information Systems (GIS) have enhanced our ability to notice intricate changes on the Earth's surface [15]. With this revolution, it has been considerably easier to identify changes across a large geographical area and over time [16]. The thermal conditions of cities and their surroundings are strongly influenced by human-oriented development and related land-use conversion of the natural landscape into impermeable surfaces [15, 17, 18, 19].

Due to the surface reflectance and roughness of different LULC types, the LST of different surface areas varies greatly [20]. Researchers have discovered that remote sensing plays a vital role in determining ecological conditions and tracking change at both geographical and temporal scales [21, 22, 23, 24]. Due to increasing urbanization, land surface types have been changing in recent years [25]. Using remotely sensed satellite imagery, a change in surface temperature may be observed, which is an indicator of impervious surface build-up [26]. The distribution of LST is greatly influenced by the presence of natural vegetation [15, 27, 28, 29]. To study variations in LST, the Normalized Difference Vegetation Index (NDVI) is often used [12, 30, 31]. Different scientists from different historical periods have studied the linear correlation between LST and other landscape features using spectral indices for different research regions, including Addis Ababa [32], Brisbane [33], Florence and Naples [34], Mexico [35], Philadelphia [36], Raipur [37], and Shanghai [38].

Some researchers have discovered that both natural and socioeconomic variables have an impact on the LST pattern [39, 40, 41]. Furthermore, recent research has demonstrated that remotely sensed Landsat imagery from satellites was successfully used to construct LULC and LST maps to analyze land cover changes that replaced natural vegetated surfaces with manufacturing facilities [42, 43, 44]. Landsat's multi-spectral images were used to perform our research. LST was estimated using thermal imagery from Landsat 5 TM and Landsat 8 TIRS data. In addition, the relationship between LST and various land covers has been the subject of discussion among researchers [45, 46, 47]. Differential urban-rural LULC composition, heat conductivities of urban surfaces, vegetation coverage, anthropogenic discharge, and built-up density are all contributors to LST intensity [48, 49]. Bala et al. [46] discovered that water and vegetation contribute negatively to the urban heat island, while barren ground and built-up areas contribute positively.

Specific band combinations from the satellite imagery were used to calculate spectral indices, such as the Normalized Difference Vegetation Index (NDVI), Normalized Difference Water Index (NDWI), and Normalized Difference Built-up Index (NDBI), and the Normalized Difference Bareness Index (NDBAI). NDVI corresponds to vegetation cover, NDWI corresponds to water bodies, NDBI corresponds to impervious build-ups, and NDBAI corresponds to bare lands. These indices are then utilized to demonstrate the relationship between multiple landscape features and LST.

The main objective of this research is to investigate the relationship between LST and landscape features by their corresponding spectral indices. We used Google Earth Engine (GEE) and the raster datasets available to calculate LST and four spectral indices and then perform logistics regression to visualize correlations between the aforementioned indices and LST. Jupyter notebook was used to calculate the Pearson's correlation coefficient and visualize the time series.

Materials and methods

Google Earth Engine was used to filter out the cloud-free images of Landsat 5 TM and Landsat 8 TIRS, which were used to calculate LST for the first epoch (2001–2012) and second epoch (2013–2018) using thermal bands of Landsat 5 TM and Landsat 8 TIRS, respectively. We selected this period (2000–2018) for our study to visualize the changes that occurred since the turn of the 21st century. The temporal extent was set to 18 years as it provides six-year equal intervals (2000, 2006, 2012, and 2018).

Study area

Sylhet Sadar Upazila was the study area for this research. It is located in the south-eastern part of Bangladesh at 24.8917° North and 91.8833° East. Around 829,103 people live in this area, 70% of whom are in urban areas [50]. It is the eighth-largest city corporation by population as well as the third most important city after Dhaka and Chattogram. Sylhet Sadar Upazila is situated on the banks of the Surma River. It has a moderate subtropical climate with highlands. The digital elevation derived from USGS (United States Geological Survey) Earth Explorer is shown in Figure 1.

Figure 1

The Digital elevation map of the study area.

In Figure 1, we can see the study area of Sylhet Sadar Upazila, which is located in north-eastern Bangladesh on the north bank of the Surma River and serves as the administrative center for the Sylhet division [51]. Sylhet features a lush highland topography with a subtropical climate [52]. The elevation map of the study area is depicted in Figure 1. The majority of the hilly region is located in the eastern half of the research area.

Landsat image filtering

Google Earth Engine was used to filter cloudless images from 2000 to 2018 from Landsat 5 and Landsat 8 OLI/TIRS. A total of ten images (at biennium intervals) were used for this study. Images from 2004 and 2012 have minor cloud covers; however, as these clouds are outside of the study area, no corrections were performed. Characteristics of data are shown in Table 1.

Table 1

Characteristics of collected data.

Year	Satellite	Sensor	Resolution	Cloud Cover
2000	Landsat 5	Multispectral Scanner (MSS) & Thematic Mapper (TM)	Band 1–5 and 7: 30 m; Band 6: resampled to 30 m from 60 m	0
2002				0
2004				0.6
2006				0
2008				0
2010				0
2012				0.6
2014	Landsat 8	Operational Land Imager & Thermal Infrared Sensor	Band 1–9: 30 m; Band 10–11: resampled to 30 m from 100 m	0
2016				0
2018				0

Land surface temperature estimation

Landsat thermal infrared measurements were utilized to estimate LST using the single-channel (SC) method [53, 54, 55, 56]. Landsat 5 and 8 both have thermal bands; however, Landsat 8 is the only one with two [57]. Using GEE, we utilized the following equations (Eqs. (1), (2), (3), (4), (5), (6), (7), and (8)) to calculate the LST.

From landsat 5 TM

Spectral radiance was calculated using Eq. (1) [58].where = Spectral radiance, = Quantized calibrated pixel value in DN, and = Spectral radiance scaled to QCALMIN and QCALMAX, and = Minimum and Maximum quantized calibrated pixel value (corresponding to LMINA) in DN.

After calculation spectral radiance, temperature or LST was calculated using Eq. (2) [58].where = Calibration constant one of band 6, = Calibration constant two of band 6, = Spectral radiance.

From landsat 8 OLI/TIRS

Landsat 8 collections, given by the USGS and included in the GEE data catalog, are the primary source of data for this study.

Digital numbers (DN) were converted to Top of Atmospheric radiance (TOA) using Eq. (3) [58].where = Multiplicative rescaling factor of band 10, = Band 10, = Additive rescaling factor of band 10.

To calculate brightness temperature, TOA and two thermal conversion constants were used in Eq. (4) [58].where = Thermal conversion constant one of band 10, = Thermal conversion constant two of band 10.

NDVI is used to calculate emissivity. NDVI and proportion of vegetation are calculated using Eqs. (5) and (6) respectively [58].

Emissivity is the radiation capacity of a surface compared to that of a black body [58, 59]. It is calculated by using Eq. (7) [58].Finally, LST is calculated using Eq. (8) [26, 58]

Spectral indices calculation

Spectral indices (e.g., NDVI, NDWI, NDBI, and NDBAI) were calculated with Eq. (9).

The band combination for each index is mentioned in Table 2.

Table 2

Band combination of spectral indices.

Index	Landsat 5 TM		Landsat 8 OLI
	First Band	Second Band	First Band	Second Band
NDVI	4	3	5	4
NDWI	2	4	3	5
NDBI	5	4	6	5
NDBAI	5	6	6	10

Results

Changes in the overall LST range from −6 °C to +4 °C; This change in overall LST is correlated to NDVI, NDBAI, and NDWI, while mixed correlation is overserved in NDBI. The study's detailed findings are explained below.

Changes in LST from 2000 to 2018

LST is heavily influenced by air surface temperature, although it is a good indicator of heat-retaining or heat-reflecting surfaces. Built-up and urban regions have greater temperatures because they reflect more heat than the Earth's surface [60, 61]. According to Figure 2, the LST difference for the period of the study in the area ranges from −6 °C to +4 °C. The north-eastern part of the area shows relatively low temperature due to higher vegetation and reserve forest.

Figure 2

Changes in LST from 2000 to 2018.

Furthermore, temperatures in the western part of the area fluctuate from 0 to 4 °C. LST, on the other hand, increased in the south-eastern section of the region as a result of urban activities and built-up areas. Urbanization increases built-up regions, which retain more heat than surrounding areas due to their impervious nature. Mean LST and Average air temperature describe in Table 3.

Table 3

Comparison of mean LST and Mean air temperature.

Type	2000	2006	2012	2018
Air Temperature (°C)	21.97	23.33	23.24	23.62
Land Surface Temperature (°C)	21.12	24.95	23.45	23.67

Dynamics of spectral indices

In Figure 3, we can see the graphical relationship of all spectral indices with LST. These values were normalized for better visualization. The NDVI, NDWI, and NDBAI all exhibit a rapidly increasing trend versus the LST, as shown in Figure 3(a, c, d), respectively. However, the NDBI indicates a mixed decrease tendency with LST (see Figure 3(b)). NDVI was the most effective indicator for categorizing native forests, and it was also useful for classifying water bodies, fallow land, and artificially planted woods. Since it is less affected by the soil and the impacts of the atmosphere, NDVI is a useful metric for places with medium to high vegetation density [62]. The NDWI is best suited for bodies of water that emphasize certain water features while minimizing others.

Figure 3

Dynamics of all four spectral indices (NDVI, NDBI, NDWI, NDBAI) and LST from 2000 to 2018 (biennium interval).

Discussion

Figure 4 shows the Pearson's correlation coefficient between all four spectral indices and LST, while the color scheme indicates the level of correlation. In Figure 4(a), LST shows a moderate correlation with the spectral indices. NDBI shows a correlation of −0.52, indicating that a moderate negative correlation exists between LST and build-ups; LST is positively correlated to the three other indices.

Figure 4

Correlation between all four spectral indices and LST; (a) mean from 2000 to 2016, (b) for 2018.

NDWI and NDBAI have correlation coefficients of 0.98 and 0.81, respectively, with NDVI, showing that they are each strongly positively correlated with NDVI. NDVI and NDBI, on the other hand, correlate at −0.98, indicating that they are significantly negatively associated. Furthermore, the NDBI has a substantial negative connection with all of the indicators. The correlation between NDWI and NDBAI is 0.73, denoting a strong positive correlation between these two indicators.

Vegetation type and growth phase might affect NDVI [63]. As NDVI shows a significant increase with a clear upwards trend, vegetation and canopy cover also increased during this period. High NDVI values typically indicate the presence of green vegetation, while high NDBI values indicate the presence of built-up areas and impervious land features. LST increases with the more built-up area and barren ground but lowers with increased forest, farmland, wetland, and water bodies [34] due to the thermal characteristics of these different surfaces. Research demonstrates that NDVI varies seasonally and is heavily influenced by the landscape, with over half of all land regions worldwide exhibiting a strong seasonal trend [64]. Increased NDBI signifies increased build-ups of impervious surfaces or urbanization. As impervious surfaces show the highest temperatures in built-up areas [26], LST should be positively correlated with NDBI. Figure 4(b) shows a substantial positive association between NDBI and LST in 2018, suggesting built-up areas increased in 2018 compared to 2000–2016, altering the thermal characteristics. The amount of bare land diminishes as built-up areas grow or vegetation reclaims lands [65, 66].

Due to the presence of complexity in landscape composition, LST-NDVI and LST-NDBI both build stronger correlations in large natural landscapes but tend to be weaker in small built-up areas [34, 45, 47, 67, 68, 69, 70].

Figure 4(b) shows the correlation between all indices in 2018. Based on Figure 4(b) only, NDBI (0.6) and NDBAI (0.51) show a moderately positive correlation with LST, while NDWI (−0.6) shows a moderate negative correlation. NDVI has a positive correlation with NDWI (0.13) and NDBAI (0.24) but shows a negative correlation with NDBI (−0.13). Furthermore, NDBI and NDBAI have a 0.87 correlation, showing a significantly positive correlation, while NDWI shows a strong negative correlation. However, NDWI and NDBAI have a correlation coefficient of -0.87, showing that they are substantially negatively correlated.

Figures 5 and 6 show the scatter plots of LST against four spectral indices for 2000, 2006, 2012, and 2018. Logistic regression was conducted with 500 samples each year for each index to determine the trends of correlation. The upper and lower error level was set as ± r2 score. According to Figures 5 and 6, NDVI, NDBI, NDWI, and NDBAI, respectively, plot with the LST in 2000, 2006, 2012, and 2018. LST-NDVI, LST-NDBI, and LST-NDBAI show a significant positive correlation, and LST had a significant negative correlation with NDWI. Figure 5(e, f, h) indicates a high positive connection, whereas Figure 5(g) shows a strong negative correlation. Furthermore, Figure 6(i, n, p) demonstrates a positive trend, whereas Figure 6(m, o) demonstrates a decreasing trend and Figure 6(j, k, l) demonstrates no trend. Some studies also found that NDBI and NDBAI had a positive correlation with LST, while NDWI had a negative correlation with LST [71].

Figure 5

Scatter plot for NDVI and LST (a–d); NDBI and LST (e–h).

Figure 6

Scatter plot for NDWI and LST (i–l); NDBAI and LST (m–p).

Conclusion

The interaction between land surface temperature and landscape features was investigated using spectral indices for eighteen years (2000–2018) using Landsat time-series data from Sylhet Sadar Upazila. For the estimation of LST and other spectral indices, multiple Landsat imageries were used, which were obtained via GEE. We observed changes in LST from −6 °C to +4 °C over eighteen years. The built-up and barren areas had the greatest increase in LST. Compared to other land coverings, bare soil and built-up areas showed higher LST also shown by [71, 72]. The NDVI, NDWI, and NDBAI all show a rapid increase in comparison to the LST, while NDBI indicates a mixed decreased tendency with LST. Guha et al. [34] also found a strong correlation between LST and other indices. Most spectral indices have changed considerably in the previous several years.

In this study, the data reveal that as the NDVI grew from 2000 to 2018, so did the canopy cover. Several studies also found moderate-high vegetation change in this region [73, 74, 75]. The NDVI and NDWI have been used successfully to distinguish between vegetation, land cover, and surface water characteristics. We also discovered a greater correlation (0.98) between NDVI and NDWI and a negative (−0.98) correlation with NDBI. The results also showed that most of the urban built-up regions were located in the south-eastern part of the study area, where we discovered rather high LST values, rather than in the central parts.

This study is limited by the low resolution of historical datasets. Further investigation should be undertaken utilizing more credible and robust datasets to acquire more precise results. The outcomes and findings of this study may assist future sustainable planning.

Declarations

Author contribution statement

Bishal Roy, B.Sc; Ehsanul Bari, B.Sc: Conceived and designed the experiments; Performed the experiments; Analyzed and interpreted the data; Contributed reagents, materials, analysis tools or data; Wrote the paper.

Funding statement

This research did not receive any specific grant from funding agencies in the public, commercial, or not-for-profit sectors.

Data availability statement

All data used for this study are publicly available. Exact sources can be made available upon request.

Declaration of interests statement

The authors declare no conflict of interest.

Additional information

No additional information is available for this paper.