Analytical study of land surface temperature with NDVI and NDBI using Landsat 8 OLI and TIRS data in Florence and Naples city, Italy

ABSTRACT The present study focuses on determining the relationship of estimated land surface temperature (LST) with normalized difference vegetation index (NDVI) and normalized difference built-up index (NDBI) for Florence and Naples cities in Italy using Landsat 8 data. The study also classifies different land use/land cover LU–LC) types using NDVI and NDBI threshold values, iterative self-organizing data analysis technique and maximum likelihood classifier, and analyses the relationship built by LST with the built-up area and bare land. Urban thermal field variance index was applied to determine the thermal and ecological comfort level of the city. Several urban heat islands (UHIs) were extracted as the most heated zones within the city boundaries due to increasing anthropogenic activities. The difference between the mean LST of UHI and non-UHI is 3.15°C and 3.31°C, respectively, for Florence and Naples. LST build a strong correlation with NDVI (negative) and NDBI (positive) for both the cities as a whole, especially for the non-UHIs. But, the strength of correlation becomes much weaker within the UHIs. Moreover, most of the UHIs (85.21% in Naples and 76.62% in Florence) are developed within the built-up area or bare land and are demarcated as an ecologically stressed zone.

Italian urban areas have some special characteristics in the evolution of built-up surfaces. Italy is considered as one of the most developed European countries with artificial land cover or built-up surface. High LST is normally found within the UHI and specifically, within the UHS. There are a number of analytical research works available regarding the relationships between LST and LU-LC indices, especially in the Italian urban environment. A new linearregression-based downscaling method of LST has been developed for Landsat 5 TM from 120 to 30 m resolution in Florence city and it has been compared with the other downscaling techniques (Bonafoni, Anniballe, Gioli, & Toscano, 2016). An interesting and innovative attempt was made in Milan, Rome, Bologna and Florence cities by using MODIS data products (1 km) to observe annual and seasonal pattern of spatial relationships between LST and built-up surfaces (Morabito et al., 2016). Another unique research work was performed to develop urban heat-related elderly risk index in 11 major Italian cities using MODIS daytime and night-time LST (Morabito et al., 2015). In Italy, inland and coastal environment can play different roles in the LST generation and in the relationship of generated LST with green areas and built-up surfaces within the urban areas.
A number of research works have been carried out to introduce some thermal comfort indices for measuring the effect of UHI intensity, i.e. temperature humidity index, physiological equivalent temperature, wet-bulb globe temperature and urban thermal field variance index (UTFVI). (Kakon, Nobuo, Kojima, & Yoko, 2010;Matzarakis, Mayer, & Iziomon, 1999;Willett & Sherwood, 2012;Zhang, 2006). Among these indices, UTFVI is the most widely used index for the ecological evaluation of urban environment due to its direct relation to LST (Alfraihat, Mulugeta, & Gala, 2016;Li et al., 2013;Liu & Zhang, 2011;Mackey, Lee, & Smith, 2012;Nichol, 2005). Determining the ecological comfort level of each and every important city is a very important task, especially for developed Italian cities.
The specific objectives of this study are: (1) to retrieve LST from Landsat 8 OLI and TIRS data and prepare spatial LST distribution maps for two Italian cities: (a) Florenceconsidered as an inland city and (b) Naplesconsidered as a coastal city; (2) to identify the UHI and non-UHI, based on the retrieved LST over the whole Florence and Naples cities; (3) to correlate LST with normalized difference vegetation index (NDVI) and normalized difference built-up index (NDBI) for the whole Florence and Naples and for the UHI and non-UHI existed within the city; (4) to build a relationship between UHI with built-up area and bare land portions of both the cities; and (5) to evaluate the thermal comfort level for both the cities using UTFVI.

Study area and data
In the present study, the spatial distribution of LST and the correlation of LST with NDVI and NDBI have been analysed in cities having Mediterranean and humid subtropical climate. Florence and Naples of Italy were selected for the research work. False colour composite (FCC) imageries of Florence and Naples city are shown in Figure 1 (Source: Unites States Geological Survey). Red colour shows the green vegetation, while blue or cyan colour presents the built-up areas. Florence is an inland city of northern Italy far away from the coastal environment, while Naples is situated on the western coast of southern Italy.
Florence is the capital city of Tuscany and metropolitan city of Italy. It is located in the Arno river basin and interlocked by a number of hills, like Careggi, Fiesole, Settignano, Arcetri, Poggio Imperiale and Bellosguardo. The latitudinal and longitudinal extension of Florence is from 43°42ʹ 36ʺ N to 43°50ʹ 24ʺ N and from 11°08ʹ 24ʺ E to 11°20ʹ 24ʺ E, respectively. Geologically, Florence is characterized by Quaternary alluvial and lacustrine deposits. The city has an elevation of 22-350 m. The soil is mainly alluvial, rich in iron oxides, clay minerals and organic matters with shallow water content. Florence is characterized by humid subtropical (Cfa) climate blended with Mediterranean (Csa) climate. Summer season is hot, while winter season is cool and moist. July and August are the hottest months, while December and January are the coldest months. The temperature of Florence often crosses the coastal temperature due to the absence of prevailing wind. Rainfall is convectional (in summer) and orographic (in winter). The mean annual temperature varies from 6°C to 25°C and average annual rainfall is about 850-900 mm [Source: World Meteorological Organisation (United Nations)]. The total area of Florence city is about 63 km 2 . Florence experiences a very high density of population (3700 persons/km 2 ) [Source: Italian National Institute of Statistics (Italian: Istituto Nazionale di Statistica; Istat)]. Most of the residents are Italian while the majority of the immigrants come from Romania, Albania, China and the Philippines. Tourism, manufacturing, food and wine production are the dominant industries of Florence. Florence is well known for Italian fashion hub and architectural heritage.
Naples is the capital of Campania and it is the third-largest municipality in Italy. It is located on the Gulfs of Naples and attains an elevation of 467 m from the mean sea level. It extends from 40°47ʹ 24ʺ N to 40°55ʹ 12ʺ N latitude and 14°08ʹ 24ʺ E to 14°21ʹ 02ʺ E longitude. Mount Vesuvius and Campi Fiegrei are the two nearest volcanoes. Geologically, Naples is characterized by extrusive igneous rocks, escarpment and included valleys. Volcanic soils with clay materials, iron oxides and alluvial deposits are the major soil types in Naples region. It has Mediterranean (Csa) and humid subtropical climate (Cfa) characterized by hot, dry summer and cool, wet winter. Mean annual temperature is 8-23°C and average annual rainfall is about 900-1000 mm [Source: World Meteorological Organisation (United Nations)]. Hottest months are July and August, while minimum temperature is found in December and January. Naples has an area of 107 km 2 . Naples is considered as the most densely populated city in Italy with >8000 persons/km 2 [Source: Italian National Institute of Statistics (Italian: Istituto Nazionale di Statistica; Istat)]. Percentage of foreign immigrants is low compared to the cities of northern Italy. Naples is considered as the fourth-largest economy in Italy and also is listed among the UNESCO world heritage sites.
Landsat 8 OLI and TIRS data of 10 July 2016 for Florence (Figure 2(a) shows the FCC image of Florence) and of 21 July 2016 for Naples (Figure 2 (b) shows the FCC image of Naples) were used to determine the LST and to detect the UHI and UHS during the summer season. ASTER GLOBAL DEM data of 17 October 2011 for Florence (ASTGDEMV2_0N43E011) and for Naples (ASTGDEMV2_0N40E014) were used to get the information regarding elevation (Figure 2(c) and (d)). The ASTER DEM data sets (Source: Unites States Geological Survey) were georeferenced and resampled at 30 m resolution to match the Landsat data set. These elevation data were very helpful to understand the relationships between LST-NDVI and LST-NDBI.
The Landsat 8 OLI and TIRS data sets were created by the US Geological Survey and were obtained in geographically tagged image file format (GeoTIFF). The data specification of Landsat 8 OLI and TIRS sensor is shown above (Table 1). Landsat 8 TIRS data set has two bands (band 10 and band 11) with thermal characteristics. Due to the larger calibration uncertainty associated with band 11, it is recommended that users should work with band 10 data as a single spectral band and should not attempt a split-window algorithm using both bands 10 and 11. The optical bands have been used in developing NDVI and NDBI. High-resolution Google Earth image was used as the reference data for LU-LC identification. ERDAS Imagine 9.1 and ArcGIS 9.3 were used as remote sensing and GIS software to analyse the data and obtain the final output throughout the entire study.

Methodology
The overall methodology of the present study has been presented in a flowchart (Figure 3).

Image preprocessing
The satellite data acquired from the Landsat 8 sensor were subset to limit the data size. First of all, the satellite imageries were geometrically and radiometrically corrected to enhance the image quality. The thermal infrared band (band 10) for Landsat 8 TIRS image has a spatial resolution of 100 m. This thermal band was resampled using the nearest neighbour algorithm with a pixel size of 30 m to match the optical bands. To analyse the changes in temperature and LU-LC types in the study region, the NDVI and NDBI were applied for determining the correlation with derived LST.

Extraction of different LU-LC types using NDVI and NDBI
The NDVI (Purevdorj, Tateishi, Ishiyama, & Honda, 1998;Tucker, 1979) is the most common and widely used index for vegetation extraction which was applied in this study. NDBI (Zha, Gao, & Ni, 2003) was also applied in this study to detect the built-up area. These two indices can be applied to categorize different types of LU-LC (Chen et al., 2006) by the suitable threshold values. To get more accurate classification, boolean operators may be used on the spectral bands of the indices. For example, NDVI > 0.2 and NDBI < 0 may be used together to extract vegetation. Similarly, NDVI < 0 and NDBI < 0 may be used together to extract water bodies, whereas 0 < NDVI < 0.2 and NDBI > 0.1 may be used together to extract built-up area and bare land. But these threshold values may differ due to atmospheric conditions. These values may also be integrated for LU-LC classification (Chen et al., 2006). An iterative self-organizing data analysis technique classifier with the unsupervised clustering method and the maximum likelihood classification were used later to get a better result in extracting the built-up area and bare land. The values of kappa

Retrieving LST from Landsat thermal band
LST was retrieved from band 10 of the Landsat 8 OLI and TIRS image of Florence and Naples cities using the following algorithm (Artis & Carnahan, 1982): where L λ is the spectral radiance in Wm −2 sr −1 mm −1 .
where T B is the brightness temperature in Kelvin (K), L λ is the spectral radiance in Wm −2 sr −1 mm −1 ; K 2 and K 1 are calibration constants. For Landsat 8 OLI, K 1 is 774.89 and K 2 is 1321.08. The surface emissivity,ε, was estimated using the NDVI thresholds method (Sobrino, Munoz, & Paolini, 2004;Sobrino, Raissouni, & Li, 2001). The fractional vegetation,F v , of each pixel was determined from the NDVI using the following equation (Carlson & Ripley, 1997): where NDVI m in is the minimum NDVI value (0.2) where pixels are considered as bare soil and NDVI m ax is the maximum NDVI value (0.5) where pixels are considered as healthy vegetation. dε is the effect of the geometrical distribution of the natural surfaces and internal reflections. For heterogeneous and undulating surfaces, the value of d ε may be 2%.
where ε v is vegetation emissivity, ε s is soil emissivity, F v is fractional vegetation and F is a shape factor whose mean is 0.55 (Sobrino et al., 2004).

Mapping UHI
UHI and non-UHI were identified by the range of LST determined by the following equations (Ma, Kuang, & Huang, 2010;Guha, Govil, & Mukherjee, 2017): where μ and δ are the mean and standard deviation of LST in the study area, respectively.

Delineating the UHS
In this study, LST maps were used in delineating the UHS over Florence and Naples to provide special emphasis for continuous monitoring. These small pockets are too hot and unfavourable for human settlement and mostly developed within the UHI. These UHS were delineated by the following equation (Guha, Govil, & Mukherjee, 2017): The urban thermal field variance index (UTFVI) A number of thermal comfort indices are available for evaluating the UHI impacts on the quality of urban life. In this study, the UTFVI has been used for the ecological evaluation of UHI zones of Florence and Naples cities in a dry summer month. UTFVI has been estimated using the following equation (Zhang, 2006): where UTFVI= Urban Thermal Field Variance Index T s = LST (°C) T mean = Mean LST (°C)

Results and discussion
Spatial distribution of LST, NDVI and NDBI

Validation of derived LST with respect to MODIS data
Before performing any kind of application process, a validation of derived LST is very necessary with in situ measurement or with another type of satellite sensor. In the present research, MODIS data were used to validate the LST values as a reference image.
MODIS and any Landsat sensors do not provide quality images for a same study area in a specific acquisition date. Thus, MOD11A1 data (1 km resolution) of 11 July 2016 have been taken instead of 10 July 2016 for Florence city, and MOD11A1 data (1 km resolution) of 22 July 2016 have been taken instead of 21 July 2016 for Naples city. No precipitation or atmospheric disturbances occurred in between the acquisition date of both Landsat 8 and MOD11A1 imageries. Spatial resolution of thermal infrared band is 1 km for MOD11A1 data and 100 m for Landsat 8 data. Without performing any upscaling of Landsat 8 data or downscaling of MOD11A1 data, a strong positive correlation (0.61 for Florence city and 0.56 for Naples city) was determined between the LST values estimated from Landsat 8 data and MOD11A1 data.

Spatial distribution of UHI and non-UHI
The intensity of UHI may be defined as the difference between the average temperature of UHI and non-UHI (Table 3). In Florence, the UHI zones consistently extend from west to east along the central part of the city, while in Naples UHI zones are developed in the entire eastern portion along with the western parts ( Figure 5). The threshold value of UHI is estimated at 33.81°C for Florence and 34.24°C for Naples. For UHI in both the cities, the standard deviation values of LST show more variability than non-UHI. On an average, the calculated mean LST values of UHI are 3.15°C and 3.31°C more than the calculated mean LST values of non-UHI for Florence and Naples, respectively.

Identification of UHS
UHS were more abundant in the built-up areas along the western and central parts of Florence and eastern parts of Naples due to lack of both vegetation and shadows despite the higher albedo of exposed surface ( Figure 6). The UHS were identified by a threshold value of 36.72°C and 37.61°C for Florence and Naples, respectively. Parking area, roadways, power plants, metal roofs and industrial factories are most suitable places for the development of UHS. Almost all such hot spots have a very little or negligible amount of vegetation and water bodies.    Figure 5. Spatial extents of UHI and non-UHI in Florence (a) and Naples (b).

Relationship of LST-NDVI and LST-NDBI within the whole city, within the UHI and within the non-UHI
In general, LST presents a positive relationship with NDBI and an inverse relationship with NDVI. NDVI shows a strong negative correlation with LST for the whole city (−0.71 for Florence and −0.57 for Naples) and for the non-UHI (−0.53 and −0.42, respectively, for Florence and Naples). This strength becomes weaker for UHI (−0.38 for Florence and −0.19 for Naples). NDBI presents a highly positive correlation with LST for the whole city (0.71 for Florence and 0.61 for Naples) and for the non-UHI (0.61 and 0.50, respectively, for Florence and Naples). But in UHI, the value of correlation coefficient becomes much lower (0.21 for both the cities). This phenomenon is observed due to the presence of complexity in landscape composition. Thus, LST-NDVI and LST-NDBI both build stronger correlation in large natural landscapes, while it tends to be weaker in small built-up areas.

Relationship of LST with LU-LC types
It is very interesting fact that LST distribution is very closely related to the distribution of NDVI and NDBI. Generally, LST is negatively related to NDVI and positively related to NDBI. But, this relationship may be varied due to spatial resolution, latitudinal extension or seasonal variation. Florence and Naples are the two major cities of Italy; one is located in the inland northern part while other lies in the coastal southern part. Both the cities support the general relationship of LST with NDVI and NDBI. NDVI and NDBI are the two LU-LC indices those are highly dependent on the LU-LC types in any region. Normally, high NDVI values indicate the presence of green vegetation and high NDBI values indicate the presence of built-up area and bare land. Basically, LST increases with the increase in built-up area and bare land whereas it decreases with the increase in forest, cropland, wetland and water bodies. Figure 7 describes the built-up area and bare land of Florence and Naples cities with the help of the separate bands of Landsat data, LU-LC indices and Google Earth image. Built-up area and bare land together cover 25.18 km 2 area of Florence, i.e. 39.91% of the total city area of Florence, and cover 61.88 km 2 area of Naples, i.e. 57.95% of the total city area of Naples. Thus, Naples is considered as more urbanized city than Florence. Vegetation, river, water bodies, wetlands and croplands have been described as the other lands. Figure 8 presents a very important relationship of UHI with the built-up area and bare land. In Florence, 35.93% area is under UHI while in Naples it is almost equal (35.78%). But, the percentage of UHI in the total built-up area and bare land is significantly high in Naples (85.21%) compared to Florence (76.62%). In Florence, 23.38% of UHI zones are developed outside the built-up area and bare land. Thus, the relationship of UHI with the built-up area and bare land is stronger in Naples. Table 4 reflects the correlation between elevation with LST and LU-LC indices. There is no significant difference found in the relationships between elevation and NDVI (positive correlation) or elevation and NDBI (negative correlation) for Florence and Naples. But, Florence reflects a stronger negative correlation between elevation and LST.

Ecological evaluation of Florence and Naples using UTFVI
The UTFVI values were divided into six categories (Zhang, 2006) according to six different ecological evaluation indices (Table 5). It is clear from the Figure 9 that both the cities have shared almost equal percentage of land for two extreme categories of ecological evaluation: the excellent category (UTFVI < 0) and the worst category (UTFVI > 0.020). The areas enjoy excellent thermal conditions (i.e. UTFVI < 0), where vegetation, water bodies and wetlands are abundant. Mainly, the northeastern and the south-western portions of Florence and the central and the south-western parts of Naples (grey-coloured portions) experience such thermal condition.
But, the worst category (i.e. UTFVI > 0.020) of the ecological evaluation index also exists in a large portion of both the cities (Redcoloured portions). The strip extends from northwest to south-east section of Florence and the entire eastern part of Naples, which fall under the worst category. Here, most of the lands are impervious in nature (either bare land with exposed rock surface or built-up areas). The good and the normal thermal conditions (i.e. 0 < UTFVI < 0.010) are found in some small patches surrounding the areas under excellent condition, while the bad and the worse conditions (i.e. 0.010 < UTFVI < 0.020) exist around the areas of built-up class.

Conclusion
In this article, Landsat 8 OLI and TIRS data were used to investigate the UHI intensity effect in Florence and Naples cities in Italy and to interpret the dynamic relationship between LST with NDVI and NDBI. UHI zones were identified through LST which were distributed along the central region extends from north-west to south-east portions of Florence while it develops throughout the eastern part of Naples. Bare land and built-up area are mostly responsible for LST generation. The    presence of vegetation and water bodies reduces the LST level. Some UHS were also delineated within the UHI zones which are characterized by high concentrated LST. Furthermore, the relationships between LST-NDVI and LST-NDBI were interpreted quantitatively by linear regression analysis at the pixel level. For whole Florence and Naples cities, LST shows strong negative correlation (−0.71and −0.57) with NDVI; and strong positive correlation (0.71 and 0.61) with NDBI. The relationships become weaker for UHI. It may be due to the presence of more heterogeneous landscapes within the built-up area.
In addition, spatio-temporal dynamics of the ecological evaluation index of Florence and Naples were computed by UTFVI. Almost equal percentage of land is in excellent condition and worst condition simultaneously. It indicates that non-UHI zones (green areas and water bodies) remain almost unchanged or very little changed. Only the UHI zones are under severe heat stress. With consistent urban development, the UHI zones may worsen the eco-environmental quality and fall under worst ecological index.
In future, many additional research works may be included. First, LST may be retrieved using another method or different spatial resolution. Second, the in situ LST data may be measured with the same overpass of satellites for the calibration and validation of LST estimation. Third, apart from linear regression, several new statistical methods can be applied to estimate the correlation between LST and different LU-LC indices. Finally, ecological evaluation of UHI zones may be analysed with the inclusion of more biophysical parameters.