Land use/land cover changes and their impact on land surface temperature using remote sensing technique in district Khanewal, Punjab Pakistan

ABSTRACT

The aim of this research was to assess the land use/land cover (LULC) changes and its impact on land surface temperature (LST) using remote-sensing (RS) technique in the district Khanewal, Punjab, Pakistan. Data were pre-processed using ERDAS imagine 15 and Arc GIS 10.4 software for layer stacking, mosaicking, and sub-setting of Landsat images. After pre-processing, the supervised classification scheme was applied for the years 1980, 2000, and 2020, which explains the maximum likelihood algorithm to identify LULC changes observed in the study area. “Built-up area” in 1980 occupied 1.75% but in 2020, the build-up area increased (5.27%) compared to 2020. Vegetation area was decreased by 4.12% from 1980 to 2020 in district Khanewal. It was observed that there has been a rapid change in vegetation area to build-up area. The LST values were increased by 0.50°C from 1980 to 2020 due to the increasing build-up area from East to West direction in district Khanewal. Maximum and minimum normalized difference vegetation index (NDVI) values were observed 0.72 and −0.2 for 1980 and 2020. The regression line produced a definitive explanation, showing a strong negative correlation with NDVI and LST. The outcomes of this study indicated that a dramatic transformation took place in district Khanewal regarding a decrease in greenness due to a rapid increase in population density, urban growth, and other infrastructural developments. Thus, these results will be used in regional and urban planning and will be used for managing agriculture in the coming years of rapid environmental changes.

KEYWORDS:

• Land use/land cover
• climate change
• remote sensing
• land surface temperature
• normalized difference vegetation index
• GIS

1. Introduction

Impact of climate change directly and indirectly on agriculture worldwide (Akram et al., 2018; Hashim et al., 2020; Kawo et al., 2021; Nasim et al., 2018; Sabagh et al., 2020). Climatic and environmental changes negatively influenced the health of the crops with the passage of time in the following ways: variations in the rainfall (Kimuku & Ngigi, 2017), extreme heat stress, change in the land surface temperature (LST), and migration of the peoples from rural to urban areas (Hussain et al., 2020a; Mukherjee & Singh, 2020). The LST patterns provide helpful information about the environment and, therefore, help understand climate change (Cai et al., 2018; Ding & Shi, 2013; Sadiq Khan et al., 2020). Several factors influence it, such as soil content, vegetation changes, building material utilized inroads, and impervious and pervious surfaces (Athick et al., 2019; Chaudhuri & Mishra, 2016; Yamamoto & Ishikawa, 2020). The LST has inversely impacted green cover, and vegetation can expressively decrease LST in urban climate and make urban regions less vulnerable to climate change (Ahmed et al., 2009; Choudhury et al., 2018; Das & Angadi, 2020; Mubeen et al., 2021; Zahoor et al., 2019).

Land use/land cover (LULC) changes play a significant role in ecosystem services, impacting ecosystem functioning, local toward global scales (Sahana et al., 2016), and human variables like environment and government planning (Athick & Shankar, 2019; Sinha et al., 2012). The LULC change studies provide useful information for a better understanding of previous practices, current LULC patterns, and future LULC trajectory (Almazroui et al., 2017; Bhagyanagar et al., 2012; Zhang et al., 2016). The change detection technique is categorized in three levels based on pixel (M. J. Khan et al., 2020), features, and objects level image processing and presents a LULC classification using different techniques for change detection (Akhtar et al., 2017; Abdullah et al., 2019; Hussain et al., 2020b). It is noted by previous research, LULC changes have the prospective to have a significant environmental impact changing the urban climate (M. Ali et al., 2019). Different research has also confirmed the environmental effects on LULC changes (Cheema & Bastiaanssen, 2010; Das & Sarkar, 2019; Tan et al., 2020; B. Zoungrana et al., 2015).

Normalized difference vegetation index (NDVI) has a great ability to signal across the annual and seasonal variation in the vegetation cover activity and vegetation’s response to climate change (Hussain, 2018; Zaidi et al., 2017). The NDVI values were interlinked with the plant’s biological activities, and change in the NDVI represents the different biological activities in plants (Aredehey et al., 2018). The daily variations in LST can be well characterized easily through NDVI value, which indicates that vegetation cover state (Xu et al., 2016). The NDVI also helps to study the different vegetation phenology cycles at global and regional scales (Olmanson et al., 2016; Ullah et al., 2019). The vegetation response to climate variations, plant processes of seasonal cycles, and health of vegetation cover all of these can be monitored by NDVI (Forkel et al., 2013; Mohammed et al., 2019). The green cover growth and development stages can be monitored with the help of NDVI (Kidane et al., 2019; Omran, 2012; Yuan et al., 2015). The NDVI values depend on maximum reflectance in the near-infrared; reflectance in the red region (Naz & Rasheed, 2017). The productivity of vegetation covers relates to evapotranspiration and precipitation, and NDVI is a valuable tool to study vegetation cover performance at the landscape scale (Gilani et al., 2015; I. Khan et al., 2019; Rizvi et al., 2020; B. J. Zoungrana et al., 2018).

Geographic information system (GIS) and RS have proved to be very significant in analyzing and assessing LULC changes as well as LST (Ayele et al., 2018; Hc et al., 2020). Satellite-based RS has the ability to provide synoptic LULC data at a specific time and location (Chew et al., 2016; Saleem et al., 2020). The RS and GIS have a remarkable ability to map and identify LULC and LST changes (Aboelnour & Engel, 2018; Malik et al., 2020). From the last few decades, research studies highly relied on the satellite remote sensed spatial information to map individual vegetation species (Huyen et al., 2016) and describe the variation in the vegetation types (Kharazmi et al., 2018; Pal & Ziaul, 2017; Romaguera et al., 2018). The important objective of this research is to analyze the variations in LULC and LST patterns in district Khanewal. The main objectives of this research are:

1. Assess urban change by calculating LULC for the years 1980, 2000, and 2020 in the study area.

2. Calculate Land Surface Temperature (LST) in the study area from 1980, 2000, and 2020.

3. Examine vegetation change by calculating NDVI from 1980 to 2020.

4. Examine the relationship between LST and NDVI from 1980 to 2020 in district Khanewal.

2. Materials and methods

2.1. Study area

The study area was bounded by district Khanewal, located in the Southern Punjab Pakistan, selected for LULC and LST changes using RS and GIS techniques. The area of Khanewal lies between longitude 71° 56ʹ 41ʹ^’ E to 72° 50ʹ 17ʹ^’ E and latitude 29° 89ʹ 30ʹ^’ N to 30° 72 ^’85ʹ^’ N. District Khanewal boundary lies on the east by Sahiwal, on the west by Multan and Muzaffargarh, on the south by Vehari and Lodhran and the north by districts of Jhang and Toba tek singh. District Khanewal has the second largest railway station in the country, the Khanewal Junction. The area of district Khanewal is an alluvial plain and flat, which is an ideal area for agricultural activities. The district Khanewal has a total area of 4,290.41 square kilometers divided into four tehsils, namely Jahanian, Mian Channu, Kabirwala, and Khanewal (Census Report 2017). Many canals are present here that irrigate this district. This irrigation network of canals makes the land very fertile. In the monsoon season, the land area close to the Chenab River is usually flooded. All areas with sweltering summers and mild winter’s features are an arid climate. The average rainfall is roughly 186 millimeters (7.3 inches). The area witnessed some of the most extreme weathers in the country with the highest recorded temperature approximately 52°C (126°F), and the lowest recorded temperature approximately −1°C (30°F). The study area map was shown in Figure 1.

Figure 1. Study area map of district Khanewal.

2.2. Landsat data

Landsat-based images of 30 m × 30 m spatial resolution covering the study area were collected from for the years 1980, 2000, and 2020 downloaded from the USGS obtained by NASA’s website (https://earthexplorer.usgs.gov/.) with mixed vegetation, an agricultural area, and bare surface as mapped LULC types (Rani et al., 2018). For identification of LULC, LST and vegetation changes, Landsat image of Operational Land Imager (OLI), Landsat 8, Enhanced Thematic Mapper Plus (ETM+), Landsat 7, and Landsat 4, 5 Thematic Mapper (TM) was used (Table 1). The LULC pattern is characterized into four classes: vegetation area, build-up area, bare land, and waterbodies of the study area.

Table 1. Specification of Landsat satellite data.

No.	Satellite/Sensor	Pixel size	Spectral resolution	Band used	Path/Row	Date
1	LANDSAT 5	30 m	Multispectral (8 bands)	1,2,3,4,5,7	149/039	March, 1980
					149/040
2	LANDSAT 5	30 m	Multispectral (8 bands)	1,2,3,4,5,7	149/039	March, 2000
					149/040
3	LANDSAT 8	30 m	Multispectral (11 bands)	1,2,3,4,5,6,7,9	149/039	March, 2020
					149/040

2.3. Image processing and classification

For the analysis of LULC and LST, multispectral Landsat images from satellite observatory were used. Landsat images consist of eight discrete bands; for this study, only 1–5 and 7 bands were used as band 6 is a thermal band used for LST and band 8 so it was not used for further analysis to confirm the bands and then write more. Landsat images were processed in software ERDAS imagine 15 for layer stacking (stacking the method used to produce a multiband image from discrete bands.), mosaicking (to combine the two stacked images) and sub setting (After stacking study area was extracted) of the image according to the study area (Rahman et al., 2017). After pre-processing, the supervised classification scheme was applied on temporal satellite datasets for the years 1980, 2000, and 2020 which explains the maximum likelihood algorithm to identify LULC changes observed in the study area. For each LULC types, training samples were selected by delimiting polygons around different demonstrative sites using software Arc GIS 10.4 (Nayak & Fulekar, 2017). Spectral signatures of different LULC class’s resultant from the Landsat images were noted using the pixels enclosed by these polygons. Detail methodology was represented in Figure 2.

Figure 2. Flow chart for methodology.

2.4. Accuracy assessment

Accuracy assessment was considered an important phase in the valuation of several image processing processes in classification. Error matrix was the greatest mutual and general means to present accuracy assessment. Producer accuracy, user accuracy, and overall accuracy were calculated by using Equation (1)(1) $\mathrm{O}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{a}\mathrm{l}\mathrm{l}\mathrm{a}\mathrm{c}\mathrm{c}\mathrm{u}\mathrm{r}\mathrm{a}\mathrm{c}\mathrm{y}=\frac{\begin{array}{c}\mathrm{n}\mathrm{u}\mathrm{m}\mathrm{b}\mathrm{e}\mathrm{r}\mathrm{o}\mathrm{f}\mathrm{s}\mathrm{a}\mathrm{m}\mathrm{p}\mathrm{l}\mathrm{i}\mathrm{n}\mathrm{g}\mathrm{c}\mathrm{l}\mathrm{a}\mathrm{s}\mathrm{s}\mathrm{e}\mathrm{s}\\ \mathrm{c}\mathrm{l}\mathrm{a}\mathrm{s}\mathrm{s}\mathrm{i}\mathrm{f}\mathrm{i}\mathrm{e}\mathrm{d}\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{l}\mathrm{y}\end{array}}{\begin{array}{c}\mathrm{n}\mathrm{u}\mathrm{m}\mathrm{b}\mathrm{e}\mathrm{r}\mathrm{o}\mathrm{f}\mathrm{r}\mathrm{e}\mathrm{f}\mathrm{e}\mathrm{r}\mathrm{e}\mathrm{n}\mathrm{c}\mathrm{e}\\ \mathrm{s}\mathrm{a}\mathrm{m}\mathrm{p}\mathrm{l}\mathrm{i}\mathrm{n}\mathrm{g}\mathrm{c}\mathrm{l}\mathrm{a}\mathrm{s}\mathrm{s}\mathrm{e}\mathrm{s}\end{array}}$(1) .

(1) $\mathrm{O}\mathrm{v}\mathrm{e}\mathrm{r}\mathrm{a}\mathrm{l}\mathrm{l}\mathrm{a}\mathrm{c}\mathrm{c}\mathrm{u}\mathrm{r}\mathrm{a}\mathrm{c}\mathrm{y}=\frac{\begin{array}{c}\mathrm{n}\mathrm{u}\mathrm{m}\mathrm{b}\mathrm{e}\mathrm{r}\mathrm{o}\mathrm{f}\mathrm{s}\mathrm{a}\mathrm{m}\mathrm{p}\mathrm{l}\mathrm{i}\mathrm{n}\mathrm{g}\mathrm{c}\mathrm{l}\mathrm{a}\mathrm{s}\mathrm{s}\mathrm{e}\mathrm{s}\\ \mathrm{c}\mathrm{l}\mathrm{a}\mathrm{s}\mathrm{s}\mathrm{i}\mathrm{f}\mathrm{i}\mathrm{e}\mathrm{d}\mathrm{c}\mathrm{o}\mathrm{r}\mathrm{e}\mathrm{c}\mathrm{t}\mathrm{l}\mathrm{y}\end{array}}{\begin{array}{c}\mathrm{n}\mathrm{u}\mathrm{m}\mathrm{b}\mathrm{e}\mathrm{r}\mathrm{o}\mathrm{f}\mathrm{r}\mathrm{e}\mathrm{f}\mathrm{e}\mathrm{r}\mathrm{e}\mathrm{n}\mathrm{c}\mathrm{e}\\ \mathrm{s}\mathrm{a}\mathrm{m}\mathrm{p}\mathrm{l}\mathrm{i}\mathrm{n}\mathrm{g}\mathrm{c}\mathrm{l}\mathrm{a}\mathrm{s}\mathrm{s}\mathrm{e}\mathrm{s}\end{array}}$(1)

The kappa (K) values measure how well RS classification agrees or is accurate with the reference data. A statistical representation of kappa (K) is as follows: A Statistical representation of kappa (K) is as follows (Safder, 2019):

(2) $\mathrm{K}=\frac{\mathrm{O}\mathrm{b}\mathrm{s}\mathrm{e}\mathrm{r}\mathrm{v}\mathrm{e}\mathrm{d}\mathrm{A}\mathrm{c}\mathrm{c}\mathrm{u}\mathrm{r}\mathrm{a}\mathrm{c}\mathrm{y}-\mathrm{C}\mathrm{h}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}\mathrm{A}\mathrm{s}\mathrm{s}\mathrm{e}\mathrm{s}\mathrm{s}\mathrm{m}\mathrm{e}\mathrm{n}\mathrm{t}}{1-\mathrm{C}\mathrm{h}\mathrm{a}\mathrm{n}\mathrm{c}\mathrm{e}\mathrm{A}\mathrm{g}\mathrm{r}\mathrm{e}\mathrm{e}\mathrm{m}\mathrm{e}\mathrm{n}\mathrm{t}}$(2)

2.5. Estimation of NDVI

The vegetation cover, for example, NDVI was calculated for the monitoring of vegetation area. The NDVI values were estimated as follows (Mohammed et al., 2019):

(3) $NDVI=\frac{NIR-RED}{NIR+RED}$(3)

where NIR is near-infrared band (TM and ETM band 4, OLI band 5) as well as R is the red band (TM and ETM band 3, OLI band 4).

The NDVI was estimated from Landsat images as well as have characteristic values between −1 and 1.

2.6. Estimation of LST

The LST is a general term describing the joint temperature of complete objects existing on the terrestrial. Different researchers used defined measurement on Landsat data to compute the LST. The LST values were estimated from the TM, TM+, and OLI thermal band like (band 6 and band 10). All steps for calculation of LST are given below (Bento et al., 2018).

Firstly, Lλ the values were estimated to spectral radiance by Equation (4)(4) $L\lambda =gain\times QCAL+offset$(4) .

(4) $L\lambda =gain\times QCAL+offset$(4)

where Lλ represents spectral radiance and QCAL is quantized calibrated pixel value in DN.

Second step, the spectral radiance value was changed to temperature by using Equation (5)(5) $T=\frac{K_2}{In\left(\frac{K_1}{L\lambda }+1\right)}$(5) (Orimoloye et al., 2018).

(5) $T=\frac{K_2}{In\left(\frac{K_1}{L\lambda }+1\right)}$(5)

where K1 and K2 values are 607.76 and 1260.56 for Landsat 4/5 (TM), and 772.88 and 1321.07 for Landsat 8 (OLI), respectively.

Final step, temperature which was in Kelvin was converted to Celsius (C°) through Equation (6)(6) $T\left(C^O\right)=T\left(K\right)-273.15$(6) (Cristóbal et al., 2018).

(6) $T\left(C^O\right)=T\left(K\right)-273.15$(6)

3. Results and discussion

3.1. Land use/and cover (LULC) changes

The analysis of LULC changes using supervised classification technique from 1980 to 2020 presented that district Khanewal was covered with different land features (built-up area, vegetation area, waterbodies, and bare soil) shown in Figure 3. In the year 1980, the vegetation area was 94.15% while waterbodies covered by (1.67%); the built-up area was covered by 1.75% while bare soil covered about 2.43%. In the year 2000, vegetation area was covered by 92.8%, followed by waterbodies 1.51%; the built-up area was covered by 3.5% while the area was covered by bare soil 2.19% (Table 2). On the other hand, in the year 2020, vegetation area was covered by 90.03%, followed by water (1.2%); the built-up area was covered by 7.02% while the area covered by bare soil 1.75%. The area under the built-up indicated the increasing trend from 2000 to 2020 in the study area.

Table 2. The LULC areal distribution from 1980 to 2020 in district Khanewal.

LULC classes	1980		2000		2020
	Ha	%	Ha	%	Ha	%
Vegetation	410622.21	94.15	404730.4	92.8	392654.65	90.03
Build-up area	7640.81	1.75	15266.18	3.5	30623.24	7.02
Bare soil	10578.01	2.43	9572.67	2.19	7625.3	1.75
Water bodies	7293.57	1.67	6565.35	1.51	5231.41	1.2
Total	436134.6	100	436134.6	100	436134.6	100

Figure 3. The LULC maps for the years 1980, 2000 and 2020 of district Khanewal.

“Built-up area” in 1980 was covered 1.75% in district Khanewal. On the other hand, in 2020, the build-up area was increased (5.27%) as compared to 2020. Although there is a large increase in “building area” has practiced an important increase during 1980 to 2020 shown in Table 3. “Vegetation area” in 1980 was covered with 94.15%, but in 2000 vegetation area decreased (1.35%) as compared to 1980 and similarly “vegetation” in 2000 occupies the class with 92.8%, but in 2020 vegetation area decreased (4.12%) as compared to 2000. The total area in district Khanewal was 10578.01 ha (2.43%) under bare soil in 1980 that decreased by 2952.71 ha (10.68%) in 2020. The waterbodies decreased by 2062.16 ha (0.47%) from 1980 to 2020 in district Khanewal, as shown in Figure 3. These changes occurred due to an increase in urban areas near the city in the study area. On the other word, a slight increase was perceived in the eastern and western sides. The built-up area was increased in the study area due to the construction of several new colonies in the east and western sides.

Table 3. Change analysis from 1980 to 2020 in district Khanewal.

LULC classes	Change 1980–2000		Change 2000–2020		Change 1980–220
	Ha	%	Ha	%	Ha	%
Vegetation	−5891.81	−1.35	−12075.8	−2.77	−17967.6	−4.12
Build-up area	7625.37	1.75	15357.06	3.52	22982.43	5.27
Bare soil	−1005.34	−0.23	−1947.37	−0.45	−2952.71	−0.68
Water bodies	−728.22	−0.17	−1333.94	−0.30	−2062.16	−0.47

Our results indicated that vegetation area and bare soil were converted to roads, commercial, and residential areas. In 1981, the total population was 1 Million, which increased to 2 Million in 2017, raising more than 1 Million populations from 1981 to 2017 in Khanewal district (Table 4). Due to the increasing population, the build-up area was increased in different regions of district Khanewal. Our results showed that almost 5% of the vegetation and bare soil have converted to built-up areas for 40 years. Different vegetation areas have been converted into residential colonies in the study area. The commercial area is also growing with the densification of a residential area. The objective of change detection is to know which LULC was decreased or increased from the previous few years. Change detection gives precise information about the rate at which the area is changing in LULC. By A. Ali et al. (2018), the built-up area has increased, while the populated area’s increase rate was a bit less during the last few years. According to Hussain et al. (2020c) found that a rapid increase in urban areas is expected; significantly impacting vegetation areas in the approaching years in district Multan. According to Usman et al. (2015), LULC classes with minimum and maximum changes under a LULC were nominated as upper and lower baselines to detect the vegetation changes during the last few years.

Table 4. District Khanewal population during 1981 to 2017 (GOP, 2017).

Census	Urban	Rural	Total
1981	152,426	915,567	1,067,993
1998	229,928	1,139,838	1,369,766
2017	364,261	1,704,229	2,068,490
Change 1981–2017	211,835	788,667	1,000,497

3.2. Accuracy assessment

Detail of user’s and producer’s accuracy as well as kappa (K) values for different LULC classes for the years 1980, 2000, and 2020 is shown in Table 5. Average producer and user accuracies were 88.1% and 91% for 1980, 89.5% and 89.6% for 2000 and 89.7% and 87.3% for 2020, respectively. The highest and lowest producer’s accuracy as well as user’s accuracy values observed for “building area” was 95.5% (2020), 86.7% (1980) and 100% (2000), 90.8% (1980), respectively. Similarly, the highest and lowest producer’s accuracy as well as user’s accuracy values observed for “bare soil” were 88.2% (2020), 84.6% (1980), and 86% (1980) and 83.7% (2000), respectively. Overall accuracy for classification for the studied years was also reasonable 0.83%, 0.88%, and 0.86% for the years 1980, 2000, and 2020, respectively. Kappa (K) values were noted 0.80%, 0.85% and 0.82% for the years 1980, 2000 and 2020, respectively, in the study area. The producer’s accuracy is the map accuracy (the producer of the map), but the user’s accuracy is the accuracy from the map user’s point of view, not the mapmaker. In our results, the two accuracies were in a fairly well range, which showed fewer errors in classification.

Table 5. Summary of producer’s and user’s accuracy and kappa (K) coefficients in district Khanewal.

LULC Classes	Season & Class										Overall accuracy	K
	Producer’s accuracy (%)					User’s accuracy (%)
	1	2	3	4	Avg.	1	2	3	4	Avg.
1980	89	86.7	84.6	80.2	88.1	92.5	90.8	86	94.7	91	0.83	0.80
2000	90.4	91.3	87	89.1	89.5	90.7	100	83.7	84.2	89.6	0.88	0.85
2020	87.7	95.5	88.2	87.4	89.7	82.5	96.7	84.7	85.5	87.3	0.86	0.82

Where: 1 = Vegetation area 2 = Build-up area 3 = Bare soil 4 = Water bodies

3.3. Land surface temperature (LST) changes

Figure 4 indicates a spatial pattern and areal distribution of LST in threeyears, for example, 1980, 2000 and 2020 in district Khanewal. The LST spatial patterns temporal and concentration shift LST highlight the changes in LULC classes rapidly. Generally, LST was estimated the range of 10.3–44.6°C during 1980, 10.8–45°C during 2000, and 11.2–44.7°C during 2020. This growth is wholly arithmetic; however, more accurate growth of the temperature was measured using the spatial average, and it designates that since 1980–2020, approximately 0.5°C LST has increased. The LST was increased due to increasing the build-up area in the study area from 1980 to 2020. The eastern part of the study area exhibits a lower temperature due to higher vegetation and agricultural land. In contrast, the central part (like district Khanewal) exhibits the rise in LST due to rapid urban expansion and declining waterbodies and vegetation cover. The highest deviation was noticed in maximum and minimum temperatures for the year 2020 (44.7°C) and (11.2°C), respectively. In 1980, maximum and minimum deviation in temperature was (44.6°C) and (10.3°C), respectively, shown in (Table 6). Considering all the limitations of RS derived LST estimation, the deviation between estimated and recorded LST is acceptable and can be used further analysis in the future like LST simulation and temperature condition index in the study area. Furthermore, the maximum LST values in Kanewala were estimated in the city’s central part (Khanewal), which is severely urbanized and quickly increasing, representing a positive relationship between LST and built-up areas.

Table 6. Maximum and minimum values of LST and NDVI.

Years	LST			NDVI
	Min	Max	Average	Min	Max	Average
1980	10.3	44.6	27.45	−0.32	0.72	0.2
2000	10.8	45	27.9	−0.3	0.67	0.186
2020	11.2	44.7	27.95	−0.2	0.55	0.175

Figure 4. The LST maps for the years 1980, 2000 and 2020 of district Khanewal.

3.4. Normalized difference vegetation index (NDVI) changes

Assessment of NDVI using Landsat image is based on the spectral properties of vegetation, absorbs visible light on vegetation, energy used in photosynthesis, and reflects near-infrared (NIR) radiation. The NDVI model derived through Landsat images of 1980, 2000, and 2020 as shown in Figure 5. ERDAS Imagine 2015 software was used to identify vegetation cover. The NDVI values ranges in 1980 were −0.32 to + 0.72 in the study area, while in 2000, NDVI values changed (minimum −0.3 and maximum +0.67). In 2020, NDVI values changed (minimum −0.2 and maximum +0.55) in the study area (Table 6). Higher NDVI values showed that there are most productive and productive areas, such as vegetation and forest. On the contrary, lower values of NDVI showed that there are least and less productive areas, for example, bare soil, waterbodies, and build-up area. Therefore, the NDVI map showed a major decrease in the most productive and productive area in the study area. Overall, the NDVI value is higher in forest lands than that for bare soil, and thus, the enlarged vegetation area may have its influence on satellite-observed greenness of vegetation throughout the study area. It was noted that there was a large (NDVI values) change in 1980 as compared to 2020. The finding by Ahmad (2012) presented that there was an increase in NDVI values in the district Multan. According to him, one of the greatest successful methods is established on tracking the historical change of a VI like NDVI values.

Figure 5. The NDVI maps for the years 1980, 2000, and 2020 of district Khanewal.

3.5. Relationship between LST and NDVI

Figure 6 shows the correlation between NDVI and LST. The regression line produced dramatic clarification, showing a strong negative relationship between NDVI and LST. The linear regression analysis (R²) for 3 years of the analysis indicated that LST and NDVI are negatively correlated to each other like 0.81, 0.74, and 0.78 for the years 1980, 2000, and 2020, respectively (Figure 6). These results show that vegetation-covered areas can decrease due to the effects of LST. The negative relationship between NDVI and LST identifies that the greater biomass of vegetation cover has lower LST. The LST and NDVI have a direct impact on LULC changes. It also implies that the territories showing most reduced NDVI values have less vegetal spread due to the urban extension, while the higher NDVI values have thick vegetal spread, and in this way, LST increase with the abatement in vegetal density. A solid relationship is seen among LST and NDVI, which guarantees the possibility of utilizing direct relapse in foreseeing LST if the estimations of NDVI are known in the study area. From this time forward, the exact LST can be calculated by utilizing NDVI. It is perceived that the NDVI values have decreased from 1980 to 2020 due to increasing urban growth and decreasing vegetation area. The maximum LST values show the littler NDVI and vice versa as well as it additionally distinguishes the territories with high LST being bare soil and regions with low temperature secured with more vegetation.

Figure 6. Relationship between LST and NDVI for the years 1980, 2000, and 2020 of district Khanewal.

4. Conclusion

The important objective of this study was to assessment LULC changes and its impact on LST in district Khanewal using RS and GIS software. “Built-up area” in 1980 covered 1.75%, but in 2020 build-up area increased (5.27%) as compared to 2020. Vegetation area was decreased by 4.12% from 1980 to 2020 in district Khanewal. It was observed that there has been a rapid change in vegetation area to build-up area. The LST values were increased 0.5°C from 1980 to 2020 due to the increasing build-up area in district Khanewal. The regression analysis is also showing a strong negative relationship between NDVI and LST. Our results show that decrease in vegetation area and water resources has also been one of the main reasons for reducing vegetation areas. Due to these changes, we lost natural ecology and biodiversity, due to increasing build-up areas may lead to a lot of environmental problems. This study will help in improving the capacity of government in the study area to device complete strategies at the local and national level for land management. Increase in agricultural production would improve economy of the people as well as our country will progress quite rapidly.

Disclosure statement

No potential conflict of interest was reported by the author(s).