Bi-Temporal Analysis of Vegetation Index on Land Surface Temperature in Kottayam, Kerala

The impact of NDVI (Normalized Difference Vegetation Index) on the LST (Land Surface Temperature) as well as on the genesis of surface heat islands in urban areas during two different time periods was assessed in Kottayam district, Kerala, India. Landsat TM, Landsat OLI and TIRS imagery from the years 1988 and 2020 were employed to scrutinize the relationship between NDVI and LST. The area covered under different NDVI classes were quantified. The finding indicated that NDVI of the research region decreased from 0.77 in 1988 to 0.59 in 2020, resulting in an increase in LST max from 34.46 °C in 1988 to 40.63 °C in 2020. Decrease in NDVI resulted in an increase in the high UHI class from 20.83 km2 in 1988 to 660.59 km 2 and from 7.26 km2 to 181.35km 2 in the very high UHI class. An inverse relationship was observed between NDVI and LST, with Pearson coefficients of 0.5737 and 0.5199 for 1988 and 2020, respectively, which indicates that NDVI could serve as a crucial metric for evaluating LST and UHI effects. Future research will explore the effect of seasonal variability in LULCC on LST, day and night time UHI and their impacts on human health and energy consumption.


Introduction
Accelerated urbanization of a particular area is always associated with the cost of the transformation of its original Land Use and Land Cover (LULC) pattern.This results in the natural green spaces to undergo a considerable transition into artificial impermeable surfaces like concrete, asphalt and metals etc. 1 Escalation of impervious surfaces increases absorption of solar radiation, storage, and reduction of long-wave radiation loss, lower albedo, thereby affecting air humidity and atmospheric temperature. 2All these factors form major contributors for increased Land Surface Temperature (LST) and for the creation of Urban Heat Island (UHI). 3,4reen cover provides shade, thus lowering surface temperature thereby reducing surface temperature through evapotranspiration process.But, LULC modification mediated through urbanization reduces green cover, there by altering the energy balance, thus rendering certain areas warmer than their surroundings, triggering the formation of UHI. 5 UHI studies can be performed either by measuring air temperature or by measuring the LST. 6LST is a response of different types of land covers.Therefore, investigating variations in LST is useful for comprehending the sustainability and ecological health of urban areas.Remote sensing (RS) sensors have the excellent potential to measure LST, and therefore, their successful adoption in UHI studies.Understanding the pattern and intensity of UHI requires an understanding of the differences in surface temperatures of land areas between well-established towns and areas with sparse development or rural settings.As a result, LST turns out to be a critical indicator to detect the level and strength of UHI. 7,8Various vegetation indices are also generated employing remote sensing technology to assess vegetation cover.NDVI (Normalized Difference Vegetation Index) is one such index, which is an indicator of LULC change (LULCC) and is employed widely for vegetation extraction.NDVI values are suggestive of greenness intensity of vegetation. 9rusal of literature reveals several reports on LULC changes and associated genesis of UHI and corresponding changes in vegetation index by using satellite images.A temporal analysis of NDVI and LST has been conducted to investigate LULC changes in Iberia.1 0 Spatial variability in the LST across Delhi, with the highest LST in urban areas has been reported. 11,12An increase in LST with land cover (LC) dynamics was observed in Bangalore. 13eports showed contrast between diurnal and various LC indices on LST over Ahmedabad. 14ndependent studies to assess the impact of LULC change for Surat showed that LULCCs had a significant impact on surface temperature. 15,168][19] The association between land use change, urban expansion and LST for Jaipur and Chennai have been investigated. 20,21Significant research is being conducted across the globe on study of Urban Heat Island (UHI) dynamics employing remote sensing indices, particularly utilizing Landsat imagery.3][24][25][26] Landsatderived indices such as the Index Based Built-up Index (IBI) has been demonstrated as valuable tools for precise UHI estimation as has been reported in Istanbul. 27A comparative study was conducted on the NDVI and UHI of Delhi and Mumbai, employing thermal satellite data for investigating the accelerated increase in urban heat and its correlation with NDVI. 28The significance of NDVI and NDBI (Normalized Difference Built-up Index) in unravelling effects of SUHI (Surface Urban Heat Island) in lasi, employing Landsat 8 (L 8 ) imagery was demonstrated. 29The study unveiled significant spatial heterogeneity, disparities across seasons, and the influence of correlations between LST and NDVI on SUHI dynamics. 30Seasonal discrepancies in RS indices as well as heat islands on urban surfaces has been witnessed along the boundary lines between the rural and urban areas in Dalian, Northeast China. 31NDVI and NDBI analysis on SUHI for Bangalore and New Delhi, indicated that NDBI was a better indicator of SUHI intensity than NDVI for both study areas. 32The urbanization in Jaipur had led to changes in LULC with significant increase in surface temperature resulting in significant negative impact on vegetation. 33The study demonstrated the importance of considering elevation in UHI studies and indicated the potential of using NDVI and Enhanced Vegetation Index (EVI) as reliable UHI indicators.Changes in LULC based on NDVI and Normalized Difference Moisture Index (NDMI) have been reported in Charaideu district of Assam. 34However, there is a paucity of reports on the LULC index and LST and UHI aspects in Kerala, India.][37][38] Reports on NDVI and UHI in Kerala's Kottayam district is inadequate.Therefore, main aspect of the present work was to analyze the potential influence of NDVI on LST and associated UHI phenomena in Kottayam district from 1988 to 2020, using geospatial techniques.

Study Location
Kottayam, a central district in Kerala, India is located between 9° 23' 20" and 9° 50' 45" N latitude and 76 o 21' 40" and 77 o 58' 20" E longitude with an area of 2204 km 2 and 3 m altitude.It is Kerala's only district that does not border the Arabian Sea or any other state.The district is surrounded by the Idukki district in the east, the Vembanad Lake and Alappuzha district in the west, the Pathanamthitta district in the south and Ernakulam in the North.The Köppen-Geiger classification system categorizes the climate in this region as "Am." (https://en.climate-data.org/asia/philippines/zambales-1832/).Kottayam enjoys a tropical climate with intense hot season in the plains and ample rains throughout the year.The district experiences relatively constant temperatures, averaging between 25°C (77°F) and 32°C (90°F), owing to its proximity to the equator.Following the hot season, which spans from March to May, the southwest monsoon persists from June to September.Gradual increase in daytime temperatures occur throughout the post-monsoon months of October and November, when the heat is almost as intense as in summer.Northeast monsoon is from December to February, and the rainfall ceases in early January.Typically, the region experiences an annual precipitation of 3130.33 millimeters.The district comprises 5 taluks namely, Vaikom, Kottayam, Meenachil, Changanassery and Kanjirapally (Fig. 1).

Methodology Remote Sensing (RS) Data
Landsat images (path/row:144/53) procured from the website of United States of Geological Survey (USGS) (https://earthexplorer.usgs.gov/)over two different periods was used for analysis using ArcGIS 10.5 software.To minimise atmospheric and seasonal effects, Landsat 5 TM (L5) (1988)  and Landsat 8 (L8) OLI/TIRS (2020) with a 30 meter spatial resolution for the same month (January) were used.The two imageries were registered to a common UTM (Universal Transverse Mercator) coordinate system, Zone 43N, with UTM WGS (World geocoded system) 1984 Projection parameters.ArcGis software was employed for standard image pre-processing techniques like extraction, georeferencing, atmospheric correction, layer stacking, and subsetting.A detailed description of the two Landsat data is presented in Table 1.

Retrieval of LST
Data from a satellite image is shown as Digital Number (DN), indirectly representing the brightness value of ground objects.Thermal bands 6 and 10 respectively for L 5 and L 8 are employed for LST computation.Even though Landsat 8 has two thermal spectral bands, 10 and 11, only 10 is employed in this study's LST computation. 39,40Band 11 is not used according to USGS recommendation due to the high calibration uncertainty.However, the extraction process of LST from Landsat TM and Landsat OLI differs slightly in terms of calculating spectral radiance (L λ ).

Conversion of the DN to L λ
Converting the DN of terrestrial objects to spectral radiance adopting equation 1 in the TIRS sensor is the basic step in LST retrieval. 39,41The L λ from thermal bands of L 5 and L 8 is calculated using the following equation.The values for L MIN , L MAX , QCal MIN and Qcal MAX were acquired from metadata file attached with Landsat images.
For L 8 thermal band, top of atmospheric radiance (L λ ) was computed using the approach described below. 42r Landsat 8, Here, M L represents band-specific multiplicative rescaling factor (value = 0.0003342), Q cal signifies digital numbers of band 10, A L stands for the bandspecific additive rescaling factor (0.1) Transformation of L λ into At-satellite brightness temperature (TB): The equation ( 3) was employed to deduce At-satellite brightness temperature from spectral radiance.

LST calculation
The obtained temperature values referred to as black body temperatures were corrected for spectral emissivity (Ɛ) to determine LST.Emissivity correction is done according to the nature of land cover by using NDVI values for every individual pixel. 43The emissivity corrected LST have been computed following, 44 ... For L 5 , Float (4-3)/ Float (4+3) For L 8 , Float (5-4)/ Float (5+4).

Estimation of NDVI
NDVI, a widely used vegetation index, is a descriptor of vegetation phenology that quantifies the disparity between near-infrared and red reflectance by summing these two components and dividing them. 45DVI extraction was accomplished with the following procedure. 46,47Values of NDVI ranges from -1 to 1.
The year 2020 also observed inclusion of 0.40 km 2 area under very high LST class in Meenachil taluk, and 0.07 km 2 in Kanjirapally taluk, which was not prevalent in 1988 (Table 3).The spatial pattern shows migration of Vaikom, Meenachil and Kottayam districts from low temperature zones in 1988 to medium and high temperature zones in 2020.These variations in LST may be attributed to changes in LULC as the developmental activities progresses in those areas.This ultimately reduces vegetative cover, thereby resulting in increased LST.Meenachil, is known for its rapid development and prosperity due to its well established education system, being well-connected to major cities and towns and with excellent access to basic amenities.It is a hub for small and medium sized enterprises, and also a popular tourist destination known for its scenic beauty, green landscapes, temples and social cohesion.Kanjirapally is also a rapidly developing taluk in Kottayam district.People from Kanjirapally have easy access to various parts of Kerala, facilitated by its well-developed infrastructure which includes a railway station and bus stand.Easy access to basic amenities like electricity, water supply etc. are also significant in this taluk, all of these which attracts human settlement in the area, resulting in increased housing and commercial structures.Table 3 summarises the spatial distribution of LST in taluks of Kottayam district in the years 1988 and 2020.Fig. 4 clearly illustrates that LST was high in urban areas whereas areas with vegetative covers and water bodies recorded lower temperatures.
Quantitative research indicates that rapid expansion of urban and/or built-up areas along with impervious materials such as materials used for construction and other infrastructure developments are key factors for increase in urban temperatures and for the resultant changes in the urban and regional climate. 50Maximum area under high LST class was observed in Kottayam taluk.Kottayam, being the district headquarters is known for its developmental activities when compared to other taluks in the district.
During the years, it has established a well-developed transport infrastructure, bloom in educational institutions and associated infrastructure like hostels etc., healthcare facilities, commercial activities, and well-established tourist destinations, and all those activities which might have accelerated the constructional activities for infrastructural facilities, thereby increasing the concrete and impervious surfaces, eventually resulting in high LST.Enhanced infrastructure, ease in access to raw materials, concentration of good administrative and entrepreneurial and technical expertise, better transport facilities creates a tendency for industries to cluster in Kottayam, thus accelerating urbanisation. 51Similar changes in spatial distribution of LST due to urbanisation has been reported for several places like Delhi, with high LST in urbanised areas, South Brazil and for Tianjin, China. 28,11,52,48ur results corroborates with the findings in Kerala's Ernakulam district where increase in LST owing to LULC changes was demonstrated.Their results showed an increase in built-up area by 2% and decline in waterbodies and wetlands by 1%.

Analysis of NDVI
The NDVI scale is used to evaluate the vegetation status of a particular area.The values of NDVI within the research region ranged from -0.33 to 0.77 in 1988 and -0.14 to 0.59 in 2020 (Fig. 5).Very high vegetation cover was observed at all the taluks in the year 1988, which disappeared in the year 2020.The NDVI values were classified under five groups (Table 4 and Fig. 7): (i) < 0 (no-vegetation); (ii) 0-0.3).This may be due to the conversion of native forests into agroecological  zones, specifically for rubber plantations, which attract farmers due to its economic importance.Kottayam comes under the core rubber producing zone, with the highest production reported from Kanjirapally taluk. 53Conversion of native forests into monoculture plantations such as rubber imposes detrimental effects on the temperature and microclimate of that area. 53,54Kerala has encountered a massive increase in urban land over the years 1973-2016 55 with expansion of built-up land by replacing vegetation and agricultural land.

UHI effects
Developmental activities in any area will be associated with changes in LULC, resulting in loss of vegetative cover and increased surface temperatures.This subsequently triggers the genesis of UHIs. Figure 8 depicts the spatial distribution of UHI in Kottayam district for the year 1988 and 2020.
The year 2020 was characterized by intense growth in UHI areas, where moderate, high and very high level UHIs increased (Table 6).The occurrence of high UHI areas was more concentrated in the Kottayam taluk, followed by Meenachil and Changanaserry, which showed good rate of development in the kottayam district indicating urban expansion.Intensity of UHI was associated with high LST areas.Rapid growth of industries and urbanisation coupled with loss of vegetation has proved to accelerate UHI formation. 57In this context, Kottayam district, being the prime rubber producer of the country provides large scale employment to people thus increasing settlement in the area.
Reduction in natural vegetation, and increase in impervious surfaces like pavements, buildings, roofs etc., provides less shade and moisture, contributing to higher temperatures.Moreover, heat generated from human activities like movement of vehicles and industrial activities emit heat into urban environment, which can contribute to heat island effects. 58,59All these form a major reason for increase in high and very high UHI classes from 1.34 % in 1988 to 40% in 2020, which is a tremendous increase.Conversely, areas under low and very low level UHIs decreased (Table 7) indicating urbanization process.Genesis of UHI associated with land cover change has been reported for Skopje, Macedonia. 25Correlation analysis among LST and NDVI and NDBI suggested weaker UHI effects in green areas and strong UHI effects in the built-up areas which are similar to the results of the present study.

Analysis of UTFVI
UTFVI was utilized to evaluate UHI impacts on ecological stability.UTFVI was classified as (Fig. 9) (i) excellent (ii) good (iii) normal (iv) bad (v) worse and (vi) worst based on ecological indexing.High UTFVI values corresponds to high UHI intensity. 60oth the study periods, 1988 and 2020 witnessed six levels of UTFVI (Table 8).However, when areas under excellent, normal and bad thermal comfort levels recorded a minimal decrease in 2020, areas under worst thermal comfort zone increased by 250% (Table 8), which was observed in the rapidly developing Meenachil and Kanjirapally taluks (Table 9).This was in good connection with LST, NDVI and UHI analysis of both the taluks.Such a remarkable increase in worst UTFVI class was also observed for Tianjin city, China, where increased UHI attributed to industrialisation and urbanisation . 48Analysis of UTFVI is considered as a very important aspect since alterations in the UTFVI impose detrimental effects on native climatic conditions, indirect economic losses, decreased comfort, and increased mortality rates. 61[64][65]

For
equation, L λ represents Spectral Radiance, L MAXλ corresponds to maximum spectral radiance (15.303 for L 5 and 17.04 for L 8 ), L MIN signifies for the minimum spectral radiances (1.238 for TM), Qcal MIN ) denotes the minimum DN value (1), Q calMAX designates maximum DN Value (255), QCal represents the DN of band 6.
Kelvin (K) to Celsius ( 0 C) degrees ...(4) (5)   Where, LST in 0 C (Celsius degrees), λ denotes the wavelength of emitted radiance in meters (11.45 for L 5 (Band 6)) and 10.8 (L8 (Band 10)).p = h*c/s = 1.4388*10 -2 m K ...(6) h stands for Planck's constant (6.626*10 -34 J s), s signifies the Boltzmann constant, equal to 1.38 * 10 -23 J/K, c represents the velocity of light, which is approximately 2.998*10 8 m/s.Land surface emissivity (Ɛ) = 0.004*Pv+0.986...(7) Applying the equation in the raster calculator, a correction value of 0.986 corresponds to the equation's adjustment factor.Where the PV (Proportion of vegetation) can be determined as;...(8)Generally, the NDVI minimum and NDVI maximum values can be revealed directly in the image.NDVI is calculated as:NDVI= Float (NIR-RED) /Float (NIR+RED) ...(9) negative value denotes water and positive values denotes vegetation.High values indicates dense greenery.NDVI= NIR-Red / (NIR+Red) ) ...(10) Mapping of Urban Heat Island (UHI): LST range value was used to identify UHI by following Equation 11. 48UHI= (T-T min )/T min ...(11) T denotes LST raster value, and T min denotes minimum LST value of the study region.Mapping of Urban Thermal Field Variance Index (UTFVI): The urban thermal ecology and the thermal comfort of the Kottayam region was defined in terms of UTFVI.The value of the LST is proportional to the amount of heat generated.The environmental impact of Kottayam's UHI zones was assessed using UTFVI, which was calculated as follows (equation 12). 49UTFVI= (T s -T mean )/T mean ...(12) Where, T s = LST ( o C); and T mean = Mean of the LST ( o C).Association between LST and NDVI: The influence of the NDVI on LST was assessed using Pearson's Product Moment correlation.The correlation was calculated in ArcGIS and Microsoft Excel software using scatter plot tool.Results and Discussion Analysis of LST Mean LST of Kottayam area was 24.82 0 C in 1988, which was elevated by 2.55 0 C reaching 27.37 0 C in 2020 (Fig. 2).Highest temperature was observed at Meenachil taluk which recorded a maximum of 40.63 o C during 2020, whereas low LST of 33.36 o C was recorded in Changanaserry (Fig. 3).LST was categorised into five classes: (a) very low, (b) low, (c) medium, (d) high, and (e) very high.The spatial distribution of LST for the years 1988 and 2020 is illustrated in Fig. 4. Between the years 1988 and 2020, a major shift was observed in areas categorized under medium LST which surged from 740.91 km 2 in 1988 to 1968.25 km 2 in 2020 (Table

Fig. 7 :
Fig. 7: Spatial distribution of NDVI of Kottayam District for the years 1988 and 2020

Fig. 8 :
Fig. 8: Spatial distribution of UHI for the years 1988 and 2020

Fig. 9 :
Fig. 9: Spatial distribution of UTFVI for the years 1988 and 2020 Fig.10and 11 confirms that NDVI shares a negative correlation with LST for the years 1988 and 2020 respectively [(R 2 = 0.5737 (1988) and 0.5199 (2020)], indicating that a dense vegetative cover lowered the surface temperature.High surface temperatures were observed in areas of inactive or less vegetation, whereas low surface temperatures were associated with dense green vegetative areas, which prevents the earth's surface from absorbing more radiation.Such negative correlation between LST and NDVI has been reported by several authors for various study regions like Wuhan, Lucknow, Vila Velha in Brazil, Skopje, Macedonia.24,49,66,25

Table 9 : Talukwise spatial distribution of EEI of Kottayam District for the years 1988 and 2020.
Correlation between LST and NDVI NDVI and LST are very sensitive to changes and variations in NDVI may initiate alterations in LST and vice versa.