Abstract

Water is a primary element for human life on Earth. Fresh water, which includes rivers, lakes, streams, and ponds, contributes less than one thousandth of a percent of the total water on Earth, but it is critical for the environment and human life. Change in land use and land cover (LULC) is a foremost concern in global environment change. Rapid changes in LULC lead to the degradation of ecosystems and have adverse effects on the environment. There is an urgent need to monitor changes in LULC and evaluate the effects of these changes in order to inform decision makers on how to support sustainable development. This study used Moderate Resolution Imaging Spectroradiometry images to detect and investigate changes in LULC patterns in Gilgit-Baltistan, Pakistan, between 2008 and 2017. Six types of LULC were used to explain the major changes of LULC in the study area. The results showed that there was a reduction of barren lands and an increase of urban areas. It also showed an inconsistent behavior of water bodies during the study. Snow area, which also increased, needs further investigation.

HIGHLIGHTS

  • Comprehensive overview of six different LULC (barren, snow, urban, forest, water bodies, and wetlands) in Gilgit Baltistan were studied.

  • Representation of the variation in spatial, temporal, compositional, and configurational patterns of different LULC were assessed using MODIS data.

  • Trend analysis of LULC changes over the last decade was measured using the Mann–Kendall test and Sen's Slope Estimator.

Graphical Abstract

Graphical Abstract
Graphical Abstract

INTRODUCTION

Freshwater resources are a primary necessity for human life on Earth; however, comprehensive evidence about variations in freshwater bodies and the capacity to store the water from these water bodies is surprisingly limited (Khandelwal et al. 2017). Despite the significance of research into long-term variations and trends of freshwater sources, the current global awareness is surprisingly inadequate (Döll et al. 2016). Water flows on the Earth's surface, water storage on the land surface, and water storage under the Earth's surface have key roles in sustaining human life and maintaining the ecological system. They have been gradually transformed by human activities including urbanization, deforestation, greenhouse gase production changes in land use, and the construction of dams and barrages, many of which cause severe and rapid flooding and can affect atmospheric processes, sea level increment, and global biogeochemical cycles. These transformations affect whole ecosystems across most regions of the world, and rapid development by humans is leading to water scarcity. Surface water observation is a functional necessity to study the processes of ecology and hydrology (Huang et al. 2018). It is important to characterize the natural and human component of freshwater systems for environmentally friendly land and water management systems and for a better ecosystem.

One of the major issues of global environmental change is land use and land cover (LULC) change. The scientific research community first met for the purpose of productive study of changes in land use at the Stockholm Conference on the Human Environment in 1972, and the second meeting was held 20 years later in 1992 at the United Nations Conference on Environment and Development. In the meantime, the International Human Dimension Programme and International Geosphere and Biosphere Programme (IGBP) have established a group to organize a research program and encourage research activities into LULC changes (Fan et al. 2007).

In recent years, the continuous monitoring of water resources with respect to natural and man-made changes has been a concern all over the world. Study of spatio-temporal dynamic of surface water, configurational and compositional features of water areas and their link with different LULC types are vital for managing and protecting ecosystems of specific types because these type a are often highly vulnerable to human activities like tourism and urbanization. Zhichao Li and Yujie Feng recently presented a study about the monitoring of water body patterns in a Mediterranean lagoon using remote sensing techniques. The spatio-temporal dynamics of water surface by water frequency index was studied on a yearly basis (Li et al. 2019).

Remotely sensed data using satellites has been widely adopted to determine the amount of LULC change across the world (Ehlers et al. 1990; Treitz et al. 1992; Harris & Ventura 1995; Xia 1996; GAR-ON YEH & LI 1999; Weng 2000). Remotely sensed images from a Moderate Resolution Imaging Spectroradiometer (MODIS) have been extensively used to study large-scale areas over a long time period (Hu & Zhang 2013; Mertes et al. 2015). MODIS have great potential for consistent large-area land surface monitoring because of its frequent revisit times and wide swaths. Nearly daily observations have enabled reliable spatial analyses of land surface phenology more reliable (Ganguly et al. 2010). Multi-temporal MODIS normalized difference vegetation index (NDVI) composites and a standard normal distribution statistical approach was applied to detect land cover changes on a yearly basis. Successful results were achieved for detection of changes in non-agricultural areas using a 4-year MODIS time series from 2002 to 2005 (Lunetta et al. 2006). Annual mapping of LULC changes using MODIS time series for Inner Mongolia for 2000–2011 was carried out successfully (Yin et al. 2014). LULC maps for the Lower Mekong Basin were developed using MODIS imagery, which showed an 87% accuracy for this area (Spruce et al. n.d.).

However, the spatial resolution of MODIS is not great for detecting small patches of water. Satellite images with high spatial resolution have been useful in such cases to obtain detailed geographical information about LULC. However, their resolution does not have great temporal resolution and they cannot perform a continuous observation at a comparatively large spatial scale. Fast and simple extraction of water bodies from remote sensing images are often done by multi-spectral thresholding methods. It is especially useful for large-scale and long time series studies. In general, to extract the water body by MODIS imagery, single band thresholds or index methods cannot efficiently handle the effect of cloud shadows and mountains, and the spectral relationship method will cause the reduction of spatial resolution of data (Huang et al. 2012).

Water patches are important for wildlife because they offer a food source, breeding areas, a safe refuge and they harbor many animal and plant species that would not survive in the adjacent landscape. Obviously, monitoring the surface water dynamics for water management is very important. However, the monitoring of water bodies dynamics is also important for ecosystem valuation (Ahamad et al. 2020) and biodiversity preservation over a long period of time (Li et al. 2019).

However, the spatial and temporal dynamics of surface water over a long period in this area of research has not given much attention to the researchers. In addition, studies have not yet been conducted on the linkage of surface hydrodynamics with the types of LULC. In this regard, this study focuses on a comprehensive overview of surface water dynamics in Gilgit-Baltistan, Pakistan, as well as the diagnosis of the effects of natural and human activities on frequency of water changes during the period 2008–2017, including the following:

  • Monitoring all frequency maps on a yearly basis using all MODIS (MCD12Q1) images available from 2008 to 2017.

  • Studying the variation in spatial, temporal, compositional, and configurational patterns of water areas using MODIS data.

  • Analyzing the connection between different LULC types and yearly surface water dynamics based on maps created by MCD12Q1.

METHODS AND DATA

The MODIS Land Cover Type MCD12Q1 uses six land cover legends to deliver an annual collection of science data sets (SDS) for global land cover with a spatial resolution of 500 m. Data were provided based on classifications of spectro-temporal features data obtained from MODIS (Sulla-Menashe & Friedl 2018).

In this study, MODIS Land Cover Data MCD12Q1 v.6 and the IGBP classification layer (which is a yearly composite at spatial resolution of 500 m) is used. A total of 20 tiles from 2008 to 2017 (the latest data is available at the end of manuscript) was used to extract land cover information for Gilgit-Baltistan (downloaded from https://earthexplorer.usgs.gov; Yao et al. (2019).

Data-preprocessing

Covers a large area of ground, from which a subset of Gilgit-Baltistan images was extracted. MODIS land products are provided by the U.S. Geological Survey in hierarchical data format, and the projection of this data is a sinusoidal (SIN) projection. Conventional data-processing software cannot handle storing this format and projection, and therefore, for further use, each scene is re-projected as a more commonly used projection known as Universal Transverse Mercator (UTM, WGS84) (Huang et al. 2012). The study area consisted of two tiles from the MCD12Q1 images and so preprocessing of the data involved the seamless combination of the mosaic of tiles and the preparation of the subsets from the full scenes.

Study area

This paper covers the high elevation region of Gilgit-Baltistan in northeast Pakistan and includes the mountain ranges of Karakorum, Hindukush, and Himalayas (Bilal 2019). Gilgit-Baltistan has a small geographic range which covers an area of 68,601.28 km2 (75° 08′ 48.12′′E & 37° 00′ 47.33′′N to 77° 41′ 11.82′′E & 35° 27′ 26.06′′N). Elevation range of Gilgit-Baltistan is 950–8,538 m above sea level (ASL)–approximately 90% is situated above 3,000 m above ASL and 12% is above 5,500 m ASL (Bin & Yinsheng 2015). The Indus River Basin water catchment relies on Gilgit-Baltistan, and the, majority of Pakistan is dependent on the area for hydroelectrically and for irrigation. The high mountain region is dominated by winter rain and is the backbone for water in Pakistan (Raza et al. 2015). This region hosts the world's three longest glaciers outside the polar regions: the Biafo Glacier, the Bartoro Glacier, and the Batura Glacier. It also has the world's second highest mountain, K-2, and several high-altitude lakes (Figure 1).

Figure 1

Study area with shuttle radar topography mission, digital elevation model values.

Figure 1

Study area with shuttle radar topography mission, digital elevation model values.

METHODOLOGY

Several change detection techniques have been used by researchers to identify differences in phenomenon or in the state of object at different times (Fan et al. 2007; Bengal 2018), such as principal component analysis, rationing of images, image differencing, change vector analysis, semi-automated classification, vegetation index differencing, unilabiate image differencing, direct multi-date classification, and so on (Fan et al. 2007; Hu & Zhang 2013). Generally, change detection methods are categorized into two different groups. In the first one, multiple images obtained at different times are classified and the results compared and evaluated (i.e. classify and compare scheme) to detect increment or decrement in water areas. In the second, there is a direct comparison of images obtained at different times, and the results are classified into different change levels (i.e. compare and classify scheme). Classify first and then compare is the most widely used methodology in the study of spatio-temporal monitoring surface water (Huang et al. 2018).

Supervised classification with the maximum likelihood algorithm is also a widely used methodology to classify MODIS images. In addition, comparison of different classes after classification in change detection technique was used to evaluate the changed types. Maximum likelihood classification works on the assumption that the facts for each class in each band are normally distributed and computes the probability that a specific pixel belongs to a specific class. Each pixel is assigned to the class that has the maximum probability (Hu & Zhang 2013).

The MODIS Land Cover Type Product (MCD12Q1) delivers land cover at global scale on an annual basis. This data is available from 2001 to present date. Supervized classification of MODIS reflectance data was used to create the MCD12Q1 product (Friedl et al. 2002). A C4.5 decision tree algorithm was used to create the IGBP scheme in Collection 5 MCD12Q1 that incorporated a full year of 8-day MODIS Nadir bidirectional reflectance distribution function-adjusted reflectance data (MCD43A2 and MCD43A4) (Schaaf et al. 2002). Substantial changes to the algorithm to classify and preprocess the data, and to the input features used in the classifications in Collection 6 were made to the MCD12Q1 SDS. Primary among these changes was the advancement of a new legend based on a nested set of classifications (Sulla-Menashe & Friedl 2018). Another advancement was the gap-fill technique to enhance the quality of data. The IGBP full name for the SDS is Land Cover Type 1, and its short name is LC_Type1. It is 8-bit unsigned data with valid range 1–17 and fill value is 255 (Table 1).

Table 1

Information of the IGBP band value used in this study

Class nameBand valueClass description
Forest with needle leaf trees The canopy is greater than 2 m of evergreen conifer trees, the tree cover is greater than 60% area. 
Forest with broadleaf trees The canopy is greater than 2 m of evergreen with broad leaves and palmate trees, the tree cover greater than 60% area. 
Forest with deciduous needle leaf trees The canopy is over 2 m and the trees cover more than 60% of the deciduous needle leaf (larch) trees. 
Deciduous broadleaf forests The canopy is over 2 m and the trees cover more than 60% of the deciduous broad leaf trees. 
Forest with mixed trees Neither deciduous nor evergreen (40–60%) tree type, canopy more than 2 m and tree coverage of more than 60%. 
Wetlands with permanent water 11 Always swamped lands with 30–60% water cover and vegetal cover greater than 10%. 
Built-up lands and urban areas 13 Minimum 30% impervious surface area should have asphalt, other building materials and vehicles as well. 
Area with permanent snow and ice 15 60% area should be covered with snow and ice for at least 10 months of the year. 
Barren lands 16 At least 60% of barren areas made up of rocks, soil, and sand. The area of vegetation less than 10%. 
Water bodies 17 Stable water bodies covered the minimum 60% of area. 
Unclassified data 255 These are missing inputs and do not have map labels. 
Class nameBand valueClass description
Forest with needle leaf trees The canopy is greater than 2 m of evergreen conifer trees, the tree cover is greater than 60% area. 
Forest with broadleaf trees The canopy is greater than 2 m of evergreen with broad leaves and palmate trees, the tree cover greater than 60% area. 
Forest with deciduous needle leaf trees The canopy is over 2 m and the trees cover more than 60% of the deciduous needle leaf (larch) trees. 
Deciduous broadleaf forests The canopy is over 2 m and the trees cover more than 60% of the deciduous broad leaf trees. 
Forest with mixed trees Neither deciduous nor evergreen (40–60%) tree type, canopy more than 2 m and tree coverage of more than 60%. 
Wetlands with permanent water 11 Always swamped lands with 30–60% water cover and vegetal cover greater than 10%. 
Built-up lands and urban areas 13 Minimum 30% impervious surface area should have asphalt, other building materials and vehicles as well. 
Area with permanent snow and ice 15 60% area should be covered with snow and ice for at least 10 months of the year. 
Barren lands 16 At least 60% of barren areas made up of rocks, soil, and sand. The area of vegetation less than 10%. 
Water bodies 17 Stable water bodies covered the minimum 60% of area. 
Unclassified data 255 These are missing inputs and do not have map labels. 

Water bodies and LULC detection

This annual composite data shows effectively and accurately the spatial and annual variations of surface water and LULC change during 2008–2017. The mosaicing and reprojection of downloaded tiles were performed using an Modis reprojection tool. Tiles were reprojected to WGS 1984 UTM zone 43N projection system. A shape file from a survey of Pakistan was used to extract the study area as a mask in ArcMap 10.3. For further analysis, each image's pixel values were transferred to Microsoft Excel.

In MS Excel, the following formula was used to calculate the percentage of different LULC types:
formula
(1)
where T = type of LULC, S = pixel count value, and ΣS = total number of pixels. Total six types of LULC (forest, permanent wetlands, built-up and urban area, permanent snow and ice, water bodies, and barren lands) was extracted for this paper.

Forest was extracted with an IGBP band value of 1–5 as a combination of all types of forest in this study area. Permanent wetlands with a band value of 11, Urban with a value of 13, permanent snow with a value of 15, barren lands with a value of 16, and water bodies with a value of 17 were extracted and converted to the shape file for further processing. These types of LULC were extracted every year from 2008 to 2017. In this study, continuous changes in LULC types and their effect of water bodies have been monitored using percentage and pattern changes of each type for each year.

Trends analysis

Mann–Kendall test

Generally, due to abruptness in data, the time series is unevenly distributed, making non-parametric tests a better choice than parameter tests (Porter et al. 2002). The non-parametric Mann–Kendall test has been widely used in time series to test for randomness against trends (Farhan et al. 2015; Hasson et al. 2015; Tahir et al. 2015, 2016; Azmat et al. 2017; Atif et al. 2018). The Mann–Kendall test is a rank-based non-parametric test for identifying trends in time series data. To detect significant statitical trends in this study, the Mann–Kendall test was used. The mathematical equations to calculate the Mann–Kendall Statistics U and test statistics X are:
formula
(2)
formula
(3)
formula
(4)
formula
(5)
where YI and Ym are the observations of time series in the form of occurrence order; i is the time series length i, vb represents the number of given ties for bth value, and b is the representation of number of tied values. Statistically significant trend existance is estimated through the value of X. A positive value of X denotes an upward trend in the time series data; negative value of X denote a downward trend in time series data.

Sen's slope estimator

To measure the magnitude of trend the non-parametric Sen's slope method was used (Wang et al. 2018b; Ye et al. 2019). Equations (6) and (7) were used to calculate slope estimate Pt:
formula
(6)
where Yu represents the data values at time U and Yv represents data values at time V. If Yu holds total number of n values in the time series, we have as many as N = n(n − 1)/2 slope estimates Pt. Median of these N values of Pt is Sen's estimator of slope, where N values of Pt were ranked in increasing order, along with the Sen's estimator.
formula
(7)

An upward trend can be detected by the positive value of P, while a negative value of P represents a downward trend in time series.

RESULTS

As discussed in Methods, a total of six types of LULC (barren land, forest, snow, permanent wetland, urban, and water bodies) were extracted over 10 years using the band threshold method from MODIS (MCD12Q1) data. This study generated two types of results for the interpretation of the LULC changes and their connection with water bodies from 2008 to 2017. Geographic information system-based maps were drawn to present a qualitative comparison between the change in specific type of LULC for all 10 years. Statistical tables were created to give us a quantitative comparison between the changes of any type of LULC for the last 10 years.

The results showed that barren land is the largest landform in the study area, and over the last 10 years (2008–2017) there was no significant change observed in barren land, with an average decrease in kilometers of 37,326–35,999 km2 and percentage of 54.41–52.47%. The decrease in barren land was very smooth, and no major unexpected behavior change was observed. However, some unexpected results were found in the total number of land patches. During 2016, there was a surprising increase of 1,451 patches, while in years, the change was steady and smooth (1,163–1,367 patches). This shows that, over time, barren land is decreasing, but its patches values showed that the barren land is converting to different types of LULC. In relation to patch size, the smallest patch pattern had no unexpected results, but the largest patch size did have some unexpected values in 2009, 2010, and 2011, and in 2010 it was very unusual. Trend analysis suggests that the barren land in Gilgit-Baltistan has a decreasing tendency. The Mann–Kendall coefficient value suggests a downward trend in the annual observation period of 2008–2017. Figure 2 shows the trends and variation of barren land during 2008–2017.

Figure 2

Barren land 2008–2017, X = Mann–Kendall value and P = sen's slope.

Figure 2

Barren land 2008–2017, X = Mann–Kendall value and P = sen's slope.

The second largest landform in this area is snow and there was a significant increase observed during (2008–2017. The average increase in kilometers was 9,108–9,984 km2 while the percentage increase was 12.48–13.68%. During 2008–2015, the increase in this landform was smooth and constant; however, in 2016 it abruptly changed before returning to normal in 2017. In 2016 the total number of snow patches and total area of snow cover land was surprisingly high. The Patch sizes are also increasing continuously, which shows the glacier and permanent snow cover of the area is increasing. Yearly trend analysis of snow cover area showed an upward trend during 2009–2018 with the value of X = 3.557 and P = 0.163% per year. The snow cover trend changes are shown in Figure 3.

Figure 3

Yearly snow cover during 2008–2017.

Figure 3

Yearly snow cover during 2008–2017.

The forests are the third largest landform in the area. During 2008–2017, there was a constant substantial increase in total forest area observed (51.34–55.83 km2). During 2008–2011, the percentage increase remained consistent (0.07%); in 2011–2014 there was a slight increase (0.08%); and in 2016–2017 it returned to its normal rate (0.07%). Overall, the forest area showed a slight and significant upward trend (X = 1.61, P = 0.0015) in Gilgit-Baltistan during the study period (Figure 4). A constant increase in number of land patches has also been observed patches from 35 to 40. The increase in total patches shows that the forest is being cut down and that large forest patches are at risk, particularly in 2017 when the largest patch decreased by 50% – which is alarming. Despite this, the total forest area has increased.

Figure 4

Forest area during 2008–2017.

Figure 4

Forest area during 2008–2017.

The urban land was the fourth largest landform, and the results show a steady increase in that land type. The increase in total urban land area was 17.83–18.69 km2, the percentage increase was 0.025–0.027% (as shown in Figure 5), and the number of land patches increased by 17–19. Trend analysis also showed an increasing trend in urban area with the Mann–Kendall trend test output values of X = 1.9677, P = 0.000078 showing an upward trend. During 2008–2017, the increase in all these parameters was very slow. Urbanization was slow and consistent, which showed that change is LULC types is largely because of natural activities and climate change. Human activities are limited in the area and did not significantly affect the LULC. Although urbanization is increasing, it is in a controlled manner and not alarming.

Figure 5

Urban area variations during 2008–2017.

Figure 5

Urban area variations during 2008–2017.

The water bodies land form was the fifth largest LULC, with a significant continuous increase in total number of land patches of 50–61, and a fluctuating percentage increase with a final increase of 0.022–0.032%. The increase in water bodies land area (in km) showed abrupt fluctuations. As we can see in 2011, 2014, and 2015, the total area of water increased abruptly (Figure 6). No significat trend was determine in the change of water bodies areas. The Mann–Kendall trend test showed values of X = 1.61 and P = 0.0014 with an upward trend. In 2011, the area of largest patches was more than three times the normal trend in the study area. For the smallest patches, the values were very consistent. The increase in water bodies patches showed that the rain in this area was high–we have already seen the increase in snow cover, which suggests that the increase in water bodies is not due to the melting of ice. Other LULC changes also do not support the area increase of water bodies. Climate change and change in precipitation can be the major reason for an increase in area of water bodies.

Figure 6

Water bodies areas during 2008–2017.

Figure 6

Water bodies areas during 2008–2017.

The study area mostly consists of mountain and snow cover, and so the wetlands was the last contributor to the LULC. The wetland was smallest land form during 2008–2017 and showed a gradual increase in total wetland area of 4.01–20.3 km2, a percentage increase of 0.005–0.02% (Figure 7), and an increase in the number of land patches of 14–50. During 2014–2017, the rate of change increased more than in previous years, but there was no abrupt change in the smallest or largest patch. Only the total area and number of patches increased and this change ratio jumped to almost double than in previous years. Trend analysis of wetlands showed an upward trend in Gilgit-Baltistan during 2008–2017 with values X = 3.75 and P = 0.0020% per year. This trend in the wetlands refers to the total precipitation change during the study years and is related to climate change.

Figure 7

Wetlands during 2008–2017.

Figure 7

Wetlands during 2008–2017.

For a better understanding and a clearer view of the trends and results, the study generated tables and used GIS-based maps. Standard deviation (SD) and mean for all LULC area changes for every year are shown in Table 2. For water bodies in 2011, the mean and SD was abrupt with a there was a significant increase, barren land had a consistent decreasing pattern, and all other LULC had a consistent increasing pattern.

Table 2

Statistics of different LULC

YearWater bodies
Barren land
Forest
Urban
Snow
Wetlands
MeanSDMeanSDMeanSDMeanSDMeanSDMeanSD
2008 0.304 0.237 32.095 980.362 1.467 3.128 1.049 1.820 9.518 162.619 0.287 0.196 
2009 0.266 0.174 30.587 679.273 1.431 3.060 1.049 1.820 9.984 187.434 0.296 0.192 
2010 0.275 0.180 30.200 673.104 1.442 3.077 1.015 1.775 10.223 191.919 0.302 0.180 
2011 0.532 1.064 29.830 668.633 1.470 3.081 1.015 1.775 9.902 188.442 0.299 0.168 
2012 0.343 0.290 28.568 913.466 1.461 3.178 1.015 1.775 9.938 187.619 0.297 0.160 
2013 0.358 0.279 28.054 904.258 1.514 3.236 1.015 1.775 10.059 190.266 0.309 0.170 
2014 0.355 0.319 26.928 882.537 1.541 3.280 1.015 1.775 9.780 188.718 0.375 0.260 
2015 0.351 0.359 25.779 861.275 1.583 3.298 1.015 1.775 9.463 186.056 0.404 0.295 
2016 0.289 0.172 24.950 845.103 1.446 3.083 1.015 1.775 8.939 181.690 0.366 0.301 
2017 0.364 0.315 26.335 867.427 1.362 1.730 0.984 1.732 9.628 189.137 0.406 0.342 
YearWater bodies
Barren land
Forest
Urban
Snow
Wetlands
MeanSDMeanSDMeanSDMeanSDMeanSDMeanSD
2008 0.304 0.237 32.095 980.362 1.467 3.128 1.049 1.820 9.518 162.619 0.287 0.196 
2009 0.266 0.174 30.587 679.273 1.431 3.060 1.049 1.820 9.984 187.434 0.296 0.192 
2010 0.275 0.180 30.200 673.104 1.442 3.077 1.015 1.775 10.223 191.919 0.302 0.180 
2011 0.532 1.064 29.830 668.633 1.470 3.081 1.015 1.775 9.902 188.442 0.299 0.168 
2012 0.343 0.290 28.568 913.466 1.461 3.178 1.015 1.775 9.938 187.619 0.297 0.160 
2013 0.358 0.279 28.054 904.258 1.514 3.236 1.015 1.775 10.059 190.266 0.309 0.170 
2014 0.355 0.319 26.928 882.537 1.541 3.280 1.015 1.775 9.780 188.718 0.375 0.260 
2015 0.351 0.359 25.779 861.275 1.583 3.298 1.015 1.775 9.463 186.056 0.404 0.295 
2016 0.289 0.172 24.950 845.103 1.446 3.083 1.015 1.775 8.939 181.690 0.366 0.301 
2017 0.364 0.315 26.335 867.427 1.362 1.730 0.984 1.732 9.628 189.137 0.406 0.342 

Yearly graphical representation of trends of LULC changes are shown in the Figures 27. The abrupt behavior of water bodies change in the area is apparent, while all other LULC have an overall consistent pattern. The decrease in barren land can be seen and it is very smooth. Urban areas are consistent during 2010–2016, or there are very small changes which cannot be detected from MODIS resolution, but in the last year (2017) an increase can be seen in urban area. Snow cover showed a constant increase, except in 2017 when it showed a very small decrease. Forests gradually increased up to 2014, but then it gradually deceased. Wetlands increased in a continuous pattern, but the growth rate increased after 2014.

Figure 8 shows GIS-based maps for the LULC changes for 2008–2017. Maps are combined for all LULC and small LULC contributors are hardly visible in the combined map due to the vastness of the study area.

Figure 8

LULC for 2008–2017.

Figure 8

LULC for 2008–2017.

Reference data and validation

The different LULC can be visually distinguished using high-resolution imagery provided by Google Earth. Google Earth imagery was therefore used to verify the sample points extracted from Landsat7 imagery with a 30-m resolution. Google Earth imagery with high resolution is the optimum platform for verifying sample points (El-zeiny & Hala 2017) for any LULC extracted from Landsat7. MODIS Tera MOD09A1 v. 6, which provided surface reflectance data set with bands 1–7 with a temporal resolution of 8 days, was used in this study to verify surface water bodies. MODIS Terra Snow Cover 8-Day with spatial resolution of 500 m (MOD10A2) was used to verify the snow cover area as a reference. Water bodies extraction was done using modified normalized difference water index from MOD09A1 and snow cover area was directly extracted from MOD10A2. Both snow and water were extracted over 8 days and converted to a yearly average. We acquired the MOD09A1 and MOD10A2 data for 2010 and 2016 for whole years and there was a total of 184 images. Results from the reference data were similar to the actual data, although we only verified the results of two LULC. Overall accuracy for snow and water was 87% and 84%, respectively.

Accuracy assessment

Accuracy assessment was done with the ready-to-use Landsat7 imagery with 30-m resolution. Annual ready-to-use Landsat7 data with 30-m resolution of enhanced vegetation index, normalized difference vegetation index (NDVI), normalized difference snow index (NDSI), and normalized difference water index (NDWI) is available on Google Earth. Accuracy for forest was determined using Landsat7 NDVI data. The simplest technique to evaluate and present accuracy is the error matrix. The main goal in the study was overall accuracy from the error matrix. The random selection technique was used to collect sample locations for this study.

In total, 30 sample points for every year were examined. Overall accuracy of forest detection with MODIS was 92.3%. NDSI is a well-known method to detect snow, and the accuracy of MODIS was tested against ready-to-use NDSI images provided by Google Earth. The same method was applied to snow and 30 samples for each year were examined. For snow, the overall accuracy of MODIS was 93%. Water detection accuracy was determined by ready-to-use NDWI Landsat imagery. MODIS did not show a very good accuracy for detecting water in mountainous areas and the overall accuracy was 86.7%. During 2008–2011, accuracy was above 90%, but after 2013 MODIS accuracy decreased. These images were also used to check the accuracy for wetlands, and we used visual interpretation to differentiate between water bodies and wetlands. In case of wetlands, MODIS accuracy was inefficient and overall accuracy was 89.3%. Wetlands in these areas only small patches and so accuracy was not high. For barren lands and urban areas, we use Landsat7 land cover with a resolution of 30 m. Same method as forest and snow was adopted and 30 sample points were chosen for visual interpretation. For the urban areas, MODIS showed a good accuracy of 94.1%, and for barren lands accuracy was 93.9%. Clear Landsat7 images for July or August were selected to allow the maximum choice of sample points.

DISCUSSION

Spatial and temporal changes in surface water patches in Gilgit-Baltistan were monitored using all available MODIS (MCD12Q1) images during 2008–2017. A quantitative connection between different LULC types and surface water changing patterns was derived to understand the possible factors affecting the changes in surface water. This method can be used to map the spatial and temporal variations of water patches in other similar types of areas and could be helpful for evaluating the environmental vulnerability caused by human and natural activities.

To improve the extraction of surface water dynamics over a yearly basis, MODIS imagery was used. The method of extraction for the different LULC is quite simple and effective for MCD12Q1. The water area results were unexpected in the study area, and there was no continuous increase or decrease of water area. Results for water area and water patches for 2011, 2014, and 2015 were unexpected. In 2011, the water area was more than twice that in 2010, but the water patches only increased compared to ten. In 2014, the water area also increased compared to 2013, and the water patches count unexpectedly increased. In 2015, water patches and water area were at their highest point of the study duration. Behavior of water patches and water area was not normal, which showed that climate change affected the area. Climate change studies should be conducted because this area is important to the whole Indus Basin in Pakistan.

It should be noted that a data comparison of different satellite (i.e. Landsat and MODIS) was not performed in this study, but different satellite and sensors have different level of performances for extraction of surface water using spectral indices (Zhai et al. 2015). The relative accuracy of each index also depends on the compositional and configurational features of surface water patches, water typologies, and thresholding methods (Fisher et al. 2016).

This study was carried out in the Gilgit-Baltistan region of Pakistan where the geography and location is suitable for a variety of human and economic activities. To protect the environment in Gilgit Baltistan from human activities, it is important to understand the connection between LULC types and surface water bodies, particularly the connection between land-use types of different social and economic activities and surface water changes. We tried to establish a connection between surface water patches and different LULC types. The complete maps representing the different LULC were important to help explain the role of LULC type in the surface water changes (Liu et al. 2015a), especially for an area with such complex geographical structure. Moreover, although there were clear changes in surface water bodies for 2011, 2014, and 2015 in the study area, the impact of LULC changes on surface water dynamic was not qualitatively assessed, and this could be considered in further studies. Many other studies have highlighted climatic factors for LULC change factors (e.g. temperature, precipitation, and evaporation) that could potentially affect the changes in surface water (Wang et al. 2018a) and the link between surface water patches in wetlands (Barbaree et al. 2018), and this should be taken into account in further studies.

Moderate spatial resolution data (MCD12Q1) was used to analyze different types of LULC in Gilgit Baltistan. MODIS data is USED extensively to study LULC and related environmental impacts OVER a large area BECAUSE rapid and convenient processing is possible (Liu et al. 2015b). Using high-resolution imagery like Landsat is very time consuming for analyzing the impact of the LULC changes on surface water patches for Gilgit Baltistan. Higher spatial resolution satellite data also usually has lower temporal resolution. However, using moderate resolution data may compromise the accuracy compared to using higher resolution data like Landsat. Future studies should use a fusion of multi-sensor data (e.g. Landsat series and MODIS) to investigate the LULC changes in a specific city or a region (Shen et al. 2016).

Although all available 500-m MODIS images from 2008 to 2017 were used in this study for detecting surface water and different type of LULC, extraction of small water patches and other small patches of different LULC was not possible, and, so these were omitted from the analysis. In future, 30-m Landsat images and time series images from Sentinel-2 with a resolution of 10 m could be used for comparison with MODIS images. Moreover, in such types of study area, satellite imagery has drawbacks due to vegetation and cloud cover. Synthetic Aperture Radar images are known to overcome such limitations. They allow observation of the Earth under any condition and can identify the water in areas with vegetation (Catry et al. 2018). Data fusion can be used for further studies (e.g. 500-m MODIS and 30-m Landsat) for more appropriate results. Fusion of MODIS and Landsat data can overcome the temporal and spatial resolution, which was not used in this study.

This study revealed an increase of snow cover area and wetlands, which indicates an increase in precipitation and possibly of temperature in the study area. This uncertainty needs further research to better understand the ecological and hydrological system of the study area.

CONCLUSION

This study observed the spatial and temporal variations of surface water bodies. Both configurational and compositional features of the water patches were studied in the Gilgit Baltistan region of Pakistan during 2008–2017 using all available MODIS (MCD12Q1) images. The study investigated variations of LULC and surface water bodies, and a quantitative link was established between different LULC types and water dynamic. This study can play a vital role in managing surface water bodies. It could also be used in other similar areas, offer baseline knowledge for the protection of such areas, and help in environmental management.

ACKNOWLEDGEMENTS

We thank our teammates for their valuable suggestion. Especially, we are grateful to the editor and anonymous reviewers for providing numerous comments and suggestions, which helped to improve this manuscript.

CONFLICTS OF INTEREST

The authors declare no conflict of interest.

DATA AVAILABILITY STATEMENT

All relevant data are included in the paper or its Supplementary Information.

REFERENCES

Ahamad
M. I.
Song
J.
Sun
H.
Wang
X.
Mehmood
M. S.
Sajid
M.
Su
P.
Khan
A. J.
2020
Contamination level, ecological risk, and source identification of heavy metals in the hyporheic zone of the Weihe River, China
.
International Journal of Environmental Research and Public Health
17
(
3
),
1070
.
Atif
I.
Iqbal
J.
Mahboob
M. A.
2018
Investigating snow cover and hydrometeorological trends in contrasting hydrological regimes of the Upper Indus Basin
.
Atmosphere 9 (5), 162 https://doi.org/10.3390/atmos9050162
Azmat
M.
Liaqat
U. W.
Qamar
M. U.
Awan
U. K.
2017
Impacts of changing climate and snow cover on the flow regime of Jhelum River, Western Himalayas
.
Regional Environmental Change
17
,
813
825
.
https://doi.org/10.1007/s10113-016-1072-6
.
Barbaree
B. A.
Reiter
M. E.
Hickey
C. M.
Elliott
N. K.
Schaffer-Smith
D.
Reynolds
M. D.
Page
G. W.
2018
Dynamic surface water distributions influence wetland connectivity within a highly modified interior landscape
.
Landscape Ecol.
33
,
829
844
.
https://doi.org/10.1007/s10980-018-0638-8
.
Bengal
W.
2018
Urbanisation and changing waterscapes: a case study of New Town, Kolkata
.
Applied Geography
97
,
109
118
.
https://doi.org/10.1016/j.apgeog.2018.04.012
.
Bilal
H.
2019
Recent snow cover variation in the Upper Indus Basin of Gilgit Baltistan, Hindukush Karakoram Himalaya
.
Journal of Mountain Science
16
,
296
308
.
Bin
S.
Yinsheng
F.
2015
Hydrological regimes under the conjunction of westerly and monsoon climates: a case investigation in the Astore Basin, Northwestern Himalaya
.
Climate Dynamics
44
,
3015
3032
.
https://doi.org/10.1007/s00382-014-2409-9
.
Catry
T.
Li
Z.
Roux
E.
Herbreteau
V.
Gurgel
H.
Mangeas
M.
Seyler
F.
Dessay
N.
2018
Wetlands and malaria in the Amazon: guidelines for the use of synthetic aperture radar remote-sensing
.
International Journal of Environmental Research and Public Health
15
,
468
.
https://doi.org/10.3390/ijerph15030468
Döll
P.
Douville
H.
Güntner
A.
Müller Schmied
H.
Wada
Y.
2016
Modelling freshwater resources at the global scale: challenges and prospects
.
Surveys in Geophysics
37
,
195
221
.
https://doi.org/10.1007/s10712-015-9343-1
.
Ehlers
M.
Jadkowski
M. A.
Howard
R. R.
Brostuen
D. E.
Company
J. W. S.
Town
O.
1990
Application of SPOT data for regional growth analysis and local planning
.
Photogrammetric Engineering and Remote Sensing
56
,
175
180
.
El-zeiny
A. M.
Hala
A. E.
2017
Environmental monitoring of spatiotemporal change in land use/land cover and its impact on land surface temperature in El-Fayoum governate, Egypt
.
Remote Sensing Applications: Society and Environment
8
,
266
277
.
https://doi.org/10.1016/j.rsase.2017.10.003
.
Fan
F.
Weng
Q.
Wang
Y.
2007
Land use and land cover change in Guangzhou, China, from 1998 to 2003, based on Landsat TM /ETM+ imagery
.
Sensors
7
,
1323
1342
.
Farhan
S. B.
Zhang
Y.
Ma
Y.
Guo
Y.
Ma
N.
2015
Hydrological regimes under the conjunction of westerly and monsoon climates: a case investigation in the Astore Basin, Northwestern Himalaya
.
Climate Dynamics
44
,
3015
3032
.
https://doi.org/10.1007/s00382-014-2409-9
.
Fisher
A.
Flood
N.
Danaher
T.
2016
Comparing Landsat water index methods for automated water classification in eastern Australia
.
Remote Sensing of Environment
175
,
167
182
.
https://doi.org/10.1016/j.rse.2015.12.055
.
Friedl
M. A.
Gopal
S.
Muchoney
D.
Strahler
A. H.
2002
Global land cover mapping from MODIS: algorithm design and preliminary results
.
IGARSS 2000. IEEE 2000 International Geoscience and Remote Sensing Symposium. Taking the Pulse of the Planet: The Role of Remote Sensing in Managing the Environment. Proceedings (Cat. No.00CH37120)
Honolulu, HI, USA, 2000, vol. 2, 527–529. doi: 10.1109/IGARSS.2000.861618
.
Ganguly
S.
Friedl
M. A.
Tan
B.
Zhang
X.
Verma
M.
2010
Land surface phenology from MODIS: characterization of the collection 5 global land cover dynamics product
.
Remote Sensing of Environment
114
,
1805
1816
.
https://doi.org/https://doi.org/10.1016/j.rse.2010.04.005
.
Gar-On Yeh
A.
Li
X. I. A.
1999
Economic development and agricultural land loss in the Pearl River Delta, China
.
Habitat International
23
,
373
390
.
https://doi.org/https://doi.org/10.1016/S0197-3975(99)00013-2
.
Harris
P. M.
Ventura
S. J.
1995
The integtation to lmprove lmagery sensed remotely atea in an urban classification
.
Photogrammetric Engineering & Remote Sensing
61
,
993
998
.
Hasson
S.
Böhner
J.
Lucarini
V.
2015
Prevailing climatic trends and runoff response from Hindukush–Karakoram–Himalaya, Upper Indus Basin
.
Earth System Dynamics
6
,
579
653
.
https://doi.org/10.5194/esdd-6-579-2015
.
Hu
J.
Zhang
Y.
2013
Seasonal change of land-use/land-cover (LULC) detection using modis data in rapid urbanization regions: a case study of the Pearl River Delta region (China)
.
IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing
6
,
1913
1920
.
https://doi.org/10.1109/JSTARS.2012.2228469
.
Huang
S.
Li
J.
Xu
M.
2012
Water surface variations monitoring and flood hazard analysis in Dongting Lake area using long-term Terra/MODIS data time series
.
Natural Hazards
62
,
93
100
.
https://doi.org/10.1007/s11069-011-9921-6
.
Huang
C.
Chen
Y.
Zhang
S.
Wu
J.
2018
Detecting, extracting, and monitoring surface water from space using optical sensors: a review
.
Review of Geophysics
56
,
333
360
.
https://doi.org/10.1029/2018RG000598
.
Khandelwal
A.
Karpatne
A.
Marlier
M. E.
Kim
J.
Lettenmaier
D. P.
Kumar
V.
2017
An approach for global monitoring of surface water extent variations in reservoirs using MODIS data
.
Remote Sensing of Environment
202
,
113
128
.
https://doi.org/10.1016/j.rse.2017.05.039
.
Li
Z.
Dessay
N.
Delaitre
E.
Gurgel
H.
2019
Continuous monitoring of the spatio-temporal patterns of surface water in response to land use and land cover types in a Mediterranean lagoon complex
.
Remote Sensing
11
(
12
),
1425
.
https://doi.org/10.20944/preprints201905.0119.v1.
Liu
Y.
Liu
X.
Gao
S.
Gong
L.
Kang
C.
Zhi
Y.
Chi
G.
Shi
L.
2015a
Social sensing: a new approach to understanding our socioeconomic environments
.
Annals of the Association of American Geographers
105
,
512
530
.
https://doi.org/10.1080/00045608.2015.1018773
.
Liu
Y.
Wang
Y.
Peng
J.
Du
Y.
Liu
X.
Li
S.
Zhang
D.
2015b
Correlations between urbanization and vegetation degradation across the world's metropolises using DMSP/OLS nighttime light data
.
Remote Sensing
7
,
2067
2088
.
https://doi.org/10.3390/rs70202067
.
Lunetta
R. S.
Knight
J. F.
Ediriwickrema
J.
Lyon
J. G.
Worthy
L. D.
2006
Land-cover change detection using multi-temporal MODIS NDVI data
.
Remote Sensing of Environment
105
,
142
154
.
https://doi.org/https://doi.org/10.1016/j.rse.2006.06.018
.
Mertes
C. M.
Schneider
A.
Sulla-menashe
D.
Tatem
A. J.
Tan
B.
2015
Detecting change in urban areas at continental scales with MODIS data
.
Remote Sensing of Environment
158
,
331
347
.
https://doi.org/10.1016/j.rse.2014.09.023
.
Porter
P. S.
Rao
S. T.
Hogrefe
C.
2002
Linear trend analysis: a comparison of methods
.
Atmospheric Environment
36
(
27
),
4420
4421
.
https://doi.org/10.1016/S1352-2310(02)00546-0
Raza
M.
Hussain
D.
Rasul
G.
Akbar
M.
Raza
G.
2015
Variations of surface temperature and precipitation in Gilgit-Baltistan (GB), Pakistan from 1955 to 2010
.
Journal of Biodiversity and Environmental Sciences
6
,
67
73
.
Schaaf
C. B.
Gao
F.
Strahler
A. H.
Lucht
W.
Li
X.
Tsang
T.
Strugnell
N. C.
Zhang
X.
Jin
Y.
Muller
J.-P.
Lewis
P.
Barnsley
M.
Hobson
P.
Disney
M.
Roberts
G.
Dunderdale
M.
Doll
C.
d'Entremont
R. P.
Hu
B.
Liang
S.
Privette
J. L.
Roy
D.
2002
First operational BRDF, albedo nadir reflectance products from MODIS
.
Remote Sensing of Environment
83
,
135
148
.
https://doi.org/10.1016/s0034-4257(02)00091-3
.
Shen
H.
Meng
X.
Zhang
L.
2016
An integrated framework for the spatio-temporal-spectral fusion of remote sensing images
.
IEEE Transactions on Geoscience and Remote Sensing
54
,
7135
7148
.
https://doi.org/10.1109/TGRS.2016.2596290
.
Spruce
J.
Bolten
J.
Srinivasan
R.
Lakshmi
V.
n.d
.
Developing land use land cover maps for the Lower Mekong Basin to aid hydrologic modeling and basin planning
.
Remote Sensing
10
(
12
),
1910
.
https://doi.org/10.3390/rs10121910
.
Sulla-Menashe
D.
Friedl
M. A.
2018
User guide to collection 6 MODIS land cover (MCD12Q1 and MCD12C1) product
.
1
18
.
https://doi.org/10.5067/MODIS/MCD12Q1.
Tahir
A. A.
Chevallier
P.
Arnaud
Y.
Ashraf
M.
Bhatti
M. T.
2015
Snow cover trend and hydrological characteristics of the Astore River basin (Western Himalayas) and its comparison to the Hunza basin (Karakoram region)
.
Science of The Total Environment
505
,
748
761
.
https://doi.org/10.1016/j.scitotenv.2014.10.065
.
Tahir
A. A.
Franklin
J.
Pierre
A.
Ul
A.
Silvia
H.
2016
Comparative assessment of spatiotemporal snow cover changes and hydrological behavior of the Gilgit, Astore and Hunza River basins (Hindukush – Karakoram – Himalaya region Pakistan)
.
Meteorology and Atmospheric Physics
128
,
793
811
.
https://doi.org/10.1007/s00703-016-0440-6
.
Treitz
P. M.
Howarth
P.
Gong
P.
1992
Application of satellite and GIS technologies for land-cover and land-use mapping at the rural-urban fringe: a case study
.
Photogrammetric Engineering & Remote Sensing
58
,
439
448
.
Wang
X.
Wang
W.
Jiang
W.
Jia
K.
Rao
P.
Lv
J.
2018a
Analysis of the dynamic changes of the Baiyangdian Lake surface based on a complex water extraction method
.
Water
10
(
11
),
1616
.
https://doi.org/10.3390/w10111616
.
Wang
Y.
Zhang
T.
Chen
X.
Li
J.
Feng
P.
2018b
Spatial and temporal characteristics of droughts in Luanhe River basin, China
.
Theoretical and Applied Climatology
131
,
1369
1385
.
https://doi.org/10.1007/s00704-017-2059-z
.
Xia
L. I.
1996
A method to improve classification with shape information
.
International Journal of Remote Sensing
17
,
1473
1481
.
https://doi.org/10.1080/01431169608948718
.
Yao
R.
Cao
J.
Wang
L.
Zhang
W.
Wu
X.
2019
Urbanization effects on vegetation cover in major African cities during 2001–2017
.
International Journal of Applied Earth Observation and Geoinformation
75
,
44
53
.
https://doi.org/10.1016/j.jag.2018.10.011
.
Ye
L.
Shi
K.
Zhang
H.
Xin
Z.
Hu
J.
Zhang
C.
2019
Spatio-temporal analysis of drought indicated by SPEI over northeastern China
.
Water
11
(
5
),
908
925
.
Yin
H.
Pflugmacher
D.
Kennedy
R. E.
Sulla-menashe
D.
Hostert
P.
Mapping
A.
2014
Mapping annual land use and land cover changes using MODIS time series
.
IEEE Journal of Selected Topics In Applied Earth Observations and Remote Sensing
7
,
3421
3427
.
Zhai
K.
Wu
X.
Qin
Y.
Du
P.
2015
Comparison of surface water extraction performances of different classic water indices using OLI and TM imageries in different situations
.
Geo-Spatial Information Science
18
,
32
42
.
https://doi.org/10.1080/10095020.2015.1017911
.
This is an Open Access article distributed under the terms of the Creative Commons Attribution Licence (CC BY 4.0), which permits copying, adaptation and redistribution, provided the original work is properly cited (http://creativecommons.org/licenses/by/4.0/).