Change dynamics of land-use and land-cover for tropical wetland management

Sustainable use and management of tropical coastal wetlands requires assessment of direct and indirect human impact on these ecosystems. Humans exert pressure on vegetation and agricultural land in wetland areas, turning to agriculture as a sustainable resource for the population. The ecosystem's resilience and recovery capacity must be determined after such impacts. They are often found in sensitive areas that are hard to reach due to the existence of habitats that provide shelter, flooded conditions and thick vegetation. Thus, remote sensing and geographical information systems (GIS) are useful tools for monitoring the spatial and temporal tendencies of land-use and land-cover (LULC) dynamics in tropical basin wetlands (Dahdouh-Guebas 2002). An assessment of the displacement of the island's estuary coast and sand barrier along the Virginia, Maryland, and New Jersey coasts, using GIS techniques, documented the impacts of Hurricane Sandy such as long-term changes, which include the impacts of past storms (Plant et al. 2018).

Different change detection techniques based on pixels – the basic image analysis unit – are classified as: (i) pre-classification for measuring change/no change: image differencing, image ratioing, regression analysis, vegetation index differencing, normalized differencing vegetation index (NDVI), change vector analysis (CVA), and principal component analysis (PCA); (ii) comparison based on classification for measuring detailed change (post-classification and multi-temporal direct or compound classification) (Hussain et al. 2013).

To assess wetland change in Kafue Flats, Zambia, Munyati (2000) used Landsat MSS images and Thematic Mapper™, as well as the change detection technique based on post-classification comparison, carrying out a supervised classification using a maximum likelihood algorithm. This methodology seems to be applicable for monitoring continental wetland systems in southern Africa.

According to Hayes & Sader (2001), the Guatemalan Mayan Biosphere Reservoir has experienced high rates of deforestation corresponding to human migration and agricultural frontier expansion. Techniques for detecting change with multi-temporal Landsat TM satellite images were used. Three pre-classification methods were assessed: (i) NDVI image differencing, (ii) principal component image differencing and (iii) RGB-NDVI change detection; NDVI values NDV(j) = (TM4 – TM3)/(TM4 + TM3) from three dates were assigned the colors red, green and blue, and categorized high, medium or low. The highest accuracy was accomplished by the RGB-NDVI method (85%).

In the Bharthapuzha River Basin, India, LULC changes were analyzed using Landsat™ multispectral images. The change detection technique was post-classification comparison using a supervised classification algorithm. An 8.7% depletion of the wetland's agricultural area was observed, and urban expansion in the basin increased by 32%. The study highlights the need for a scientific management plan for river basin sustainability (Nikhil & Azeez 2010).

For seasonal variations of the Malinda wetland, Tanzania, Kuria et al. (2014) used orthophotos from unmanned aerial vehicle (UAV) images fused with Spotlight TerraSAR-X dual-polarized radar data, combining two techniques for change detection: image differencing based on unsupervised classification by applying maximum likelihood diffuse inference and change vector analysis (CVA). The accuracy of determination of land covers was satisfactory.

Mwita (2016) assessed restoration of the Usangu wetland, Tanzania, using Landsat time-series images, and the NDVI change detection method of image differencing. The land cover and swamp size results show a need to implement measures for wetland protection and sustainability.

For the Surjpur wetland, Uttar Pradesh, India, Saha et al. (2017) studied LULC changes as a component of human interactions with the environment, using Landsat 7 and Landsat 8 OLI images, and comparing them with Google Earth images. LULC changes were detected by post-classification comparison. Saha et al. recommend wetland management and restoration.

Spruce et al. (2020) conducted a LULC mapping study for the Lower Mekong Basin, Asia, using GIS along with LULC change maps from 1997 to 2010 to be used as inputs to the SWAT hydrological model. The main purpose was the evaluation of the results and support for the management and planning of water and terrestrial resources.

The purpose of this study was to analyze the dynamics of bi-temporal and spatial changes in land use (LU) and land cover (LC) in a tropical basin with a natural wetland – the Urama River, Venezuela. The wetland has been impacted because of land demand for agricultural and livestock production (Petrochemistry of Venezuela (PEQUIVEN) 2014). The study was carried out using time periods 1986–1991, 2000–2008, and 2015–2017. The methodology included digital pre-processing of Landsat MSS, TM, ETM+ and OLI images to generate classification maps LU and LC type, by applying pre- and post- classification techniques. This novel method enables greater accuracy in determining the dynamics of LU and LC bi-temporal and spatial changes in tropical wetlands, thus contributing to their management. The change/no change of cover use is determined by the qualitative pre-classification technique, and bi-temporal change variation by classes is determined by the quantitative post-classification techniques, and thus the magnitude and direction. This enables definition of the proportionality of dynamic factors for management models for wetlands in tropical basins.

The study area is in the central North Coastal Region, Venezuela, whose boundaries are between the Carabobo and Yaracuy states. There are two rivers, the Temerla and Canoabo, in the Urama Basin. On the alluvial plain, a natural wetland is created by discharge to the Caribbean (Figure 1).

Seven multispectral images acquired by Landsat sensors were used as data sources in this study. They are available on the Earth Explorer website (USGS 2017). The criterion used to analyze LULC change detection was the 31-year data collection period from 1986 to 2017, including 1991, 2000, 2008, 2015, 2016 to 2017 (Table 1).

Preliminary Landsat image processing comprised making absolute and relative geometric, radiometric, topographic and atmospheric corrections corresponding to each image (Jensen 2009).

The selected methods are: (1) pre-classification: (A) algebraic methods (reflectance image differencing) and (B) transformation method (principal components image differencing); (2) post-classification method. The procedures are as follows.

Two multi-temporal images co-registered accurately are used corresponding to times t₁ and t₂. A pixel by pixel difference operation is applied to produce a residual image representing the change between t₁ and t₂ (Singh 1989; Hussain et al. 2013).

This involves a linear transformation of variables corresponding to rotation and translation of the original coordinate system (Singh 1989; Marquez et al. 2019), to reduce the dimensionality of the reflectance variable in the spectral bands to one or two components before applying the LULC change detection techniques.

Change dynamics is estimated as the attribute difference between the initial and final years (Jensen 2009), and the supervised classification technique is used (Marquez et al. 2019). Once classified, the map is analyzed to elaborate the recoded version by comparing the reference image from websites like Google Earth corresponding to each time period studied, in conjunction with field reconnaissance, which enables identification of spectrally homogeneous classes (PEQUIVEN 2014; Marquez et al. 2018a, 2018b).

The results are expressed by interval limits for the pixel distribution of change/no change areas in the river basin for the period 1986–2017. The limit criterion for change/no change areas is based on an interval defined by the mean (μ) +/− n standard deviation (σ), assuming that the distribution of the reflectance differencing image pixels in the near-infrared spectral region (NIR) approaches a normal distribution function and that the change/no change area limits are equidistant from the mean. The first limit is no change area equal to the mean (μ) +/−1 standard deviation (σ). The second is equal to the change area obtained as <μ−1σ,> μ+1σ; respectively; giving the lower limit (LL) and upper limit (UL).

Bi-temporal image differencing shows a slight tendency towards symmetry with respect to the mean (μ) for the time periods 1986–2017, 2000–2017, 2008–2017, 2015–2017 and 2016–2017. A bias or asymmetry was found for the period 1991–2017 (Table 2). The standard deviation is in the order of 10¹–10².

Three categories are defined for the results of the bi-temporal reflectance image differencing as expressed by the error matrix (Table 3) in the change detection accuracy for the images observed: decrease (DEC), no change (NC), and increase (INC). The reflectance differencing image pixels associated with reflectance DEC are negative because the NIR reflectance observed in the image year for time t₂ is less than the NIR reflectance observed for t₁. The NC pixels are taken from the limit mean region +/− 1 standard deviation centered at zero, whose value shows that the image reflectances at t₁ are equal to those at t₂. The INC image pixels are those where the reflectance difference is greater for the reflectance image for t₂ than the observed reflectance differencing image (ORDI) for t₁. Thus, for the period 2016–2017: (1) DEC in the predicted reflectance differencing image (PRDI) corresponding to DEC in ORDI: 109, NC in ORDI: 32, INC in ORDI: 0. The accuracy and Kappa indices in the error matrix are 79% and 0.77, respectively.

Comparison of PRDI with ORDI shows that 236 out of 300 pixels match, justifying the 79% accuracy obtained. The differences between predicted and observed pixels are based on the change/no change interval selection criterion, which depends on the number of times assigned to the standard deviation, defined by the change/no change threshold. For example, the NC pixels selected in the PRDI showed a user accuracy of 100% in the ORDI, for a sample of 31 randomly distributed points and K = 1, while, in the pixels where change occurred – increase or decrease – user accuracy varied between 75 and 77%, and K between 0.6 and 0.71 (Table 3). The Kappa statistic takes into account the non-diagonal elements by comparing them with the total accuracy index. Sinha & Kumar (2013) put forward a two-step thresholding approach for threshold value determination for the spectrally increased and decreased area using NDVI differencing image.

The results show that there is no significant difference in accuracy when distribution is nearly symmetric, confirming the importance of considering normal data distribution in LC change/no change analysis. Accuracies exceeding 90% and K 0.89 were achieved from all the change/no change determination methods used. Comparison with the accuracies in this study has shown that there are no significant differences between the accuracy and Kappa indices for the reflectance differencing image.

The reflectance differencing images during the six periods – 1986–2017, 1991–2017, 2000–2017, 2008–2017, 2015–2017, 2016–2017 – show the occurrence of positive, negative and no change values (Figure 2). The changes occur mainly during the longer periods 1986–2017 and 1991–2017. Most of the Urama Basin is covered by vegetation. Reflectance decrease (white dots in Figure 1) could be due to change from vegetation to degraded soil. Reflectance increase (black) could arise from a change to agricultural use.

The pre-classification method results based on principal component (PC) image differencing are derived by transforming the reflectance variables (ρ) in the principal component for each Landsat image. The Landsat 8 OLI image transformation for 2017 is shown in Table 4, for example, and the subsequent PC1 differences in Tables 4 and 5. The correlated reflectance dataset contained in the seven 2017 Landsat image spectral bands is transformed into a dataset consisting of variables identified as PCs, which are uncorrelated linear combinations of the original variables represented by the reflectances in the seven spectral bands.

As observed in Equation (1), the reflectance variable coefficients (covariances) are similar in the first four bands, where reflectances have been detected by sensors in the optical region. This is confirmed as the reflectance dataset correlation coefficients are not statistically different between the PCs PC1 to PC4 (Table 4), but differ slightly in the others (PC5 to PC7).

The similarity in the reflectance covariance coefficients in the optical region bands, and the difference with respect to band 5 (NIR) covariance might arise because the predominant landscape cover is vegetation and relates to the spectral response from this type of cover. The spectral reflectance profile obtained from the type of vegetation found in the Urama Basin is associated with slight and mild reflectance variations in the optical region, ranging from 5 to 7%. In the NIR, reflectance increases from 5 (optical region) to 30%. According to Jensen (2009), healthy green leaves absorb radiant energy in the blue and red portions of the spectrum where incident light is required for photosynthesis. In the NIR, healthy green vegetation is generally characterized by high reflectance (40–60%), high transmittance (40–60%) through the leaf to its underside, and relatively low absorption (5–10%).

The PC transformation results are also expressed by reflectance image eigenvalues from 1986 to 2017. Thus, the variance vector of the PCs or eigenvalues in the Landsat 8 OLI image from 10/03/2017 is: Var [ρ(PC_1–7)] = (544.938, 59.492, 4.384, 0.257, 0.216, 0.087, 0.059); while the population variance is represented by the sum of eigenvalues = 604.43 (Table 5). The highest proportion of the total population variance is represented by PC1, varying between 81.75 and 93.45% for images from 1986 to 2017 (Table 5). The PC1 selection criterion as the image to be included in the change detection technique based on principal component image differencing is that the first PC expresses the maximum possible proportion of variance in the original dataset (Ingebritsen & Lyon 1985).

The PC1 differencing images covering 1986–2017 and 2016–2017 are shown in Figure 3. The changes are shown as in Figure 2. Most changes are associated with reflectance increases through time. Like the reflectance differencing image method, the principal change could be replacement of the vegetation by agricultural use (Hayes & Sader 2001).This method highlights the reflectance increase in the reservoir in the middle of the river basin from 2000 to 2017, in particular. This could be due to a change from clear water to sediment and algae caused by eutrophication. As the concentration of suspended sediments and algae-laden increase in clear water, reflectance increases in all wavelengths for both clayey and silty soils (Han & Rundquist 1997; Lodhi et al. 1997).

The post-classification method results include procuring the statistics associated with two stages: (1) generation of the LULC classification map based on reflectance images and (2) recoding the classified map based on the acquired reference image from Google Earth.

The post-classification comparison method results using the image classification error matrix for 1986 are shown in Table 6. Some 12,637 randomly selected pixels on the classified map (CM), (left column), were compared with the uses observed in the 2017 Google Earth reference image (RI) (top row). The classes are: 1. Vegetation (C1), 2. Waterbody (C2), 3. Agriculture (C3), 4. Rural area (C4), 5. Bare soil (C5), 6. Clouds (C6), and 7. Shadows (C7).

As for the 1986 Urama Basin classified map (Table 6), vegetation (C1) was associated with 8,419 pixels, of which 8,411 coincided with C1 in RI and 8 corresponded to agriculture (C3) in RI, resulting in 99.9% user accuracy. In the waterbody (C2): 623 pixels corresponded exactly with C2 in RI –user accuracy 100%. For agriculture (C3): of 910 pixels, 133 corresponded to vegetation (C1) in RI, 775 with agriculture class (C3) and 2 with bare soil (C5). Accuracy was thus 97% and K = 0.942.

The four post-classification method accuracy indices can be interpreted in the seven classified maps using statistics for a sample size equal to 12,637 pixels (Table 6, 1986 image): interval (I), mean (M) and standard deviation (SD): (1) accuracy of change/no change detection: I: 72 and 99%; M: 94.11% and SD: 9.81%. (2) Kappa Coefficient: I: 0.51 and 0.99; M: 0.89 and SD: 0.17. Trodd (1995) reviewed the methods used to assess image classifications in a literature survey of 84 classifications reported in 25 papers published in major journals from 1994 to 1995; finding that the mean accuracy was 67.66%. Foody (2002) discussed accuracy classifications below 85%, contemplating a wide range in the accuracy with which individual classes were classified.

The proportions of areas associated with LULC classes of classified maps in the Urama Basin from 1986 to 2017 show that the order of class occurrence prevalence is: a) vegetation (50–85%), b) agriculture (15–20%), c) rural areas (5–10%), d) bare soil (5–10%), e) waterbodies (5%), f) clouds (1–10%) and g) shadows (<1%). Comparing this method's results with those of the PCs method demonstrates that the pattern found in the spectral band co-variances in the optical infrared region (OIR) matches the spectral profile corresponding to vegetation or agricultural areas.

The post-classification method results using the change detection matrix corresponding to the recoded LULC classification maps, based on reflectance images in the basin from 1986 to 2017, are shown in Table 7. The total pixel area is 818,707 (736.84 km²). The no change areas were: C1: 294,983 pixels (265 km², 36%), C2: 8 pixels (0.0072 km², 0.00097%), C3: 44,141 pixels (39.72 km², 5.39%), C4: 4,328 pixels (4 km², 0.54%), C5: 6,196 pixels (5.57 km², 0.75%). The LULC changes from C1 (vegetation) from 1986 to 2017 were: 93,940 pixels to agriculture (84.85 km²), 36,030 to rural areas (32.42 km²), and 52,559 to bare soil (47.3 km²). In general, the bi-temporal change area varies between 36.89 and 57.28%; the no change area between 42.71 and 63.1%; and the change area average is C: 46.43% and NC: 53.55%.

For example, the post-classification results using LULC recoded maps in 2000 and 2016 are shown in Figure 4. The year-by-year issues are as follows.

1986: most of the area is covered by vegetation, the reservoir near the middle is relatively new, there are agricultural areas in some high and low zones, with bare soil in an upper sub-basin and around the drainage network, connecting the reservoir with the Caribbean in the north of the basin.

1991: most of the area is covered by vegetation and, to a lesser extent, reservoirs and bare soil.

2000: the reservoir's supply zone is now protected, a rural use, as the agricultural area has expanded and there is bare soil in the lower part of the basin in wetland 1, close to the Caribbean.

2008: the water supply rural area and the agricultural area within wetland 1 have expanded.

2015: the area of bare soil in the reservoir's water supply protection zone has increased, with a decrease in the extent of the occupied area by the water reservoir, bare soil with agricultural use in the lower part of the basin (wetland 1).

2016: the extent of bare soil in the reservoir's resource protection zone, and in wetland 2 (near the middle of the basin) and wetland 1 (close to the Caribbean) has increased.

2017: the extent of the agricultural areas, and of bare soil throughout the basin, including areas around the two wetlands has increased.

The difference results from the LULC are presented in order of importance. Vegetation cover tends to increase by between about 5 and 25%, it is the most sensitive of the classes to change. Agricultural use decreased to a range of 5–15%. Rural use showed an annual increase between 2 and 5%. The only other significant variation was a decrease of 10% in degraded soil between 2000 and 2016, followed by an increase of less than 5% between 2016 and 2017.

For the tropical wetland studied, the NIR reflectance difference can be used to estimate LULC vegetation, agricultural use, and bare soil changes where reflectance is reasonably high. PC1, which includes the greatest reflectance variance in the satellite image n spectral bands, was the best option in the principal component difference method. It enabled the detection of changes in five LULC classes: waterbodies, rural areas, vegetation, agricultural use and bare soil. The post-classification method enables validation of classes where change/no change were detected using pre-classification methods.

The combination of methods described enables higher levels of accuracy in the LULC change predictions of the PC method in scenarios where one class is predominant, one of the characteristics of tropical wetlands. The Urama Basin scene has extensive vegetation cover. PC1 encompasses between 80 and 90% of the reflectance variance detected by Landsat sensors and was recorded in seven spectral bands. The eigenvalue vector results from the reflectance image dataset for the time period confirmed this. This is a significant difference from the reflectance differencing method, where reflectance corresponds to a specific spectral band, which is restrictive for detecting LULC change. The post-classification method has limitations because the supervised classification results depend on producer/user accuracy in classifying specific LULC classes.

Land-use and land-cover change dynamics will enable proposal of a management model for the Urama Basin wetland, identifying vegetation, waterbodies, agricultural use and bare soil as dynamic factors, to guide land management unit planning.