Land Use/Land Cover (LULC), change Detection, and Simulation Analysis of Manila Bay’s Dolomite Mining Site in Cebu, Philippines Using Sentinel–2 Satellite

Ordoñez, C. A.^1* and Silva, D.²

¹School of Graduate Studies, Mapua University, 658 Muralla St., Intramuros, Manila 1002, Metro Manila Philippines

E-mail: cayeamariezordonez@gmail.com^*

²School of Civil, Environmental, and Geological Engineering, Mapua University, 658 Muralla St.,

Intramuros, Manila 1002, Metro Manila, Philippines

E-mail: dantesilva2000@yahoo.com

*Corresponding Author

Abstract

Dolomite mining in Barangay Pugalo – Pasol, Alcoy, Cebu, and the Manila Bay beach project along Roxas Boulevard, Metro Manila, Philippines has long been controversial due to the ecological, environmental, and economic implications. From the threat of rising sea levels to questionable management solutions, this sand project cannot be fully measured because of the absence of an environmental impact study. The primary objective of this research is to develop a geospatial study of the Cebu mining site using Sentinel–2 satellite data utilizing remote sensing techniques and Quantum GIS where 2018, 2020, 2030, and 2050 classification and simulation analysis shows an increase in Soil and a decrease in Vegetation and Built-up classes. Dolomite sand sustainability and longevity in Manila Bay's Dolomite Beach, are also assessed using related literature. With the help of this analysis, it’s possible to identify specific changes and predict which land will be impacted by upcoming years when current practices in the area do not change. It can also be utilized as resource management, environmental policy, and regulation support tools in identifying ecological problems while serving as a source of historical information for future land management and support for the long-term use of natural resources and expanding populations.

Keywords: Geographic Information System (GIS), Geospatial, Land Use/Land Cover (LULC), Remote Sensing, and Satellite Data

1. Introduction

The purpose of the paper is to increase understanding of the use of remote sensing techniques to measure local levels of growth in environmental sustainability and longevity. In relation, the 2030 Agenda for Sustainable Development Goals (SDG) particularly the Geo4SDGs or GEO-for-SDGs, which outlines the implementation, monitoring, and tracking of various sectors and indices, is effective in achieving these goals is a valuable concept of SDG in transforming the world into a sustainable future. Earth observation data, remote sensing, and geospatial information are proven to be relevant aspects in this regard. The geospatial study is a powerful communication tool that cuts across language and cultural barriers. It enables countries to design visualization tools that support accurate assessment and evaluation of the development impact across the seventeen (17) SDGs consistently when combined with other statistical data. With low utilization, the potential of the geospatial study to be a vital tool in accomplishing development goals is yet to be maximized.

The Department of Public Works and Highways (DPWH) and the Department of Environment and Natural Resources (DENR) aim to complete the beach nourishment, coastal restoration, and enhancement of the Manila Bay walk area as part of the Manila Bay Rehabilitation initiative, which Secretary Roy Cimatu lead, in January 2019, as seen in [1]. Work on these projects began in August 2020 when the government authorized Philippine Mining Service Corporation to transport crushed dolomite from Alcoy, Cebu to Manila as stated at [2]. The beach is 60 meters wide, spanning 5.4 hectares. According to [3] the project progress currently covers 500 meters of its proposed 900 meters length.

The paper's primary goal is to create a geospatial analysis of Barangay Pugalo – Pasol, Alcoy, Cebu, Philippines mining site using Sentinel-2 satellite data and understanding the sustainability and longevity of the dolomite sand in Manila Bay's Dolomite Beach Project along Roxas Boulevard, Malate, Metro Manila, Philippines adopting Remote Sensing (RS) techniques and Geographic Information System (GIS) software like Quantum Geographic Information System (QGIS) and related literature, respectively. The present study will pull from related resources that indicate how the sea level changes and how it affects the Dolomite Beach Project. According to Oceana [4], an international organization that focuses on ocean conservation, the sea had already reclaimed the initial dolomite load dumped in Manila Bay’s project between September and December 2020 with the dolomite beach degrading by at least 300 square meters between December 2020 and February 2021.

Previous and recent studies as seen in [5] are focused on the controversial budget of this ‘Rehabilitation Program’ with an estimated three hundred eighty-nine hundred million Philippine pesos (₱ 389, 000, 000) in its 1^st phase of construction and a two hundred sixty-five million Philippine pesos (₱ 265, 000, 000) budget for the 2^nd phase with a proposed total of six hundred fifty-four million Philippine pesos (₱ 654, 000, 000). This financial budget raises a high concern for the citizens on its financial management considering this is not the final figure for the project as it continues to expand for the continuation of phases. Dolomite mining in Cebu and the Manila Bay beach has long been controversial due to its ecological, environmental, and economic implications. From the threat of rising sea levels to the questionable management solutions, this sand project cannot be fully measured because of the absence of an environmental impact study particularly, spatial imagery data to assess the environmental change. As stated, [6] the beach was originally scheduled to be finished by December 2020 but was delayed to 2021 as a result of the suspension of dolomite mining activities in September 2020 for several reasons, one of which is the opposition from activists, environmentalists, and historical conservationist groups according to [7]. For instance, the University of the Philippines Marine Science Institute (MSI) [8] opposed using crushed dolomite sand to build the beach, indicating that this would not improve Manila Bay's water quality and that ongoing sand replenishment would be costly. The UP Institute of Biology [9] also stated that using mangroves rather than crushed dolomite rock would be better for rehabilitation, adding that the International Union for Conservation of Nature (IUCN) has similarly recommended the alternative.

Further criticism for the program mounted following the deluge of trash left behind by Typhoon Ulysses in November. The DENR [10] refuted claims that the dolomite had been washed away, claiming that black sand had instead washed onto the beach. Following a string of typhoons in late 2020, the DENR [11] replenished the beach with fresh supplies of crushed dolomite rock during the construction of phase one (1) which stretches the budget furthermore. Satellite data, specifically the Sentinel–2 data from the Copernicus Program, was used due to its wide range of services and applications, especially in the LULC study. However, when selecting the timestamp for each piece of data, it was important to keep in mind that the availability of the requested imagery for the subject region was limited due to issues with cloud coverage, date, year, undated spatial data, etc. In this study, the effects of geospatial analysis on environmental management and sustainability will be shown. Remote sensing and GIS will be used to modify large-scale geographic data quickly and accurately. Management of LULC and research that aids in regulating the allotment of land for uses are included in this paper to provide auxiliary data for eco-centered government policies as well as a source of historical information for ongoing and future land development research.

2. Methods and Methodology

2.1 Methodology and Data Acquisition

This chapter presents the methodology used in understanding the two study areas. Quantitative research, divided into three (3) phases with seven (7) steps is expected to be completed for the research results and study. QGIS software will be used as the main research material, produced by a GNU General Public License (GPL)-compliant project from the Open-Source Geospatial Foundation (OSGeo). Other research materials are available at data processing satellite imagery websites such as USGS Earth Explore, a geospatial information system analytical cloud-based platform that may be used for visualization and satellite image analysis in RS research and administration. Using these technologies, the researcher calculated and analyzed the statistical data. Maximum Likelihood Classification Spectral Signature - (Jeffries-Matusita Distance) is used in classification analysis. Pearson's Correlation Coefficient is used for the correlation and area changes analysis section, and Cohen's Kappa Statistic is used in both the research classification and Cellular Automata simulation process.

The study has gathered primary and secondary data. Primary data are actual spatial data analysis from the satellite data that are divided into three phases, shown in Figure 1.

2.2 Phase 1

This section will focus on installing QGIS software and acquiring spatial satellite data in a raster format for the LULC classification analysis using Semi-Automatic Classification Plug-in (SCP). Figure 2 shows a flowchart of the process of the micro-integration of SCP using the maximum likelihood algorithm. It is an effective spectral signatures calculator that compares the random closest neighborhood pixels and assumes that the statistics for each class in each band are normally distributed with the highest probability and calculates the chance that a given pixel belongs to a specific class.

Figure 1: The conceptual framework highlights the workflow of the classification, forecasting, and modeling process/ stages of the LULC data

Figure 2: Flowchart of the classification process using SCP

If the highest probability is less than the threshold's permitted value, the pixel is left unclassified. Whereas spectral signatures are the distinctive way a given type of land cover absorbs and reflects light. To get the combination of wavelengths that best identifies the spectral signature of the desired land cover, researchers manipulate the colors acquired by the satellite by using various RGB values. This allows the opportunity to discern between different types of land cover and their characteristics. Jeffries-Matusita distance is effective at determining the spectral distance between signatures when utilizing the maximum likelihood technique; it is also chosen as the spectral distance algorithm for this research. Where if the results of the signature are fully distinct, this procedure asymptotically approaches the value of two (2) and tends to show a value of zero (0) when the signatures are identical as stated in [12]. The researcher wants to know if the Maximum Likelihood Classification – Jeffries-Matusita method shows a distinct signature towards each class. The study also used Cohen’s kappa statistic to gauge the inter-rater dependability which is referred to as accuracy that occurs when data raters give the same score to the same data item. By comparing the classified map layer to the reference map layer, Cohen’s kappa statistic tabulates the error matrix and validates the classification result. As in Equation 1.

Equation 1

where:

Kis Cohen’s kappa statistic

P₀ is relative observed reliability,

P_e is hypothetical probability of chance reliability

2.3 Phase 2

This section is the application of the Modules for Land Use Change Evaluation (MOLUSCE) plug-in shown in the flowchart in Figure 3. The first step is the input of the output of phase one (1) data and the extraction of the Digital Elevation Model (DEM) and slope. The following steps are the correlations, detection of area changes, Transition potential modeling, and the Cellular Automata (CA) simulation. Pearson’s Correlation coefficient will be used to summarize the characteristics of the dataset and identify the strength and direction of the linear relationship between the quantitative variables of the research as stated in [13]. In this paper, the correlation is going to be conducted between the data of LULC 2018, LULC 2020, DEM, and Slope as they are the main variables. Pearson correlation coefficient results that are between 0 and 1 are considered a positive correlation, results equal to 0 are considered as no correlation, and results between 0 to -1 are a negative correlation. To determine how effectively LULC was evaluated, the kappa coefficient was frequently used. The validation between the predicted and reference map is carried out by calculating the entire kappa value. Cellular automata simulation is used to carry out the LULC change prediction method after the kappa value satisfies the assessment requirement. This approach is frequently used in geographic information systems to forecast changes in land use. This approach is based on the idea that each cell or pixel adapts to the circumstances of its nearest neighbors as well as several other spatial factors.

Figure 3: Flowchart of the detection changes and simulation process using the MOLUSCE plug-in

2.4 Phase 3

The last phase of the study will extract important LULC raster-classified data maps from the previous phases of the methodology to develop a layout map and a three-dimensional (3D) model of the study area using the Qgis2threejs tool in QGIS.

3. Results and Discussion

3.1 Study Area

The researcher’s first objective is to understand the paper's study areas which are done by backtracking each historical background and topographic detail. Figure 4 shows the map of the Philippines wherein the study area is highlighted in different colors. Metro Manila, where Manila Bay is located is highlighted by colored lines particularly pink and Cebu province in red color lines where Barangay Pugalo – Pasol, Alcoy is located – mining site of the dolomite sand.

Figure 4: a) Philippine Map (b)Metro Manila in color pink lines where Manila Beach project is located (c) Cebu Province in red color lines where Barangay Pugalo – Pasol is located

3.1.1 Manila Bay’s area (MBA) Dolomite beach project along Roxas Boulevard, Malate Manila, Metro Manila, Philippines

This study area is also known as the nation’s biggest Bay and the sewer systems discharge destination of some regions, such as the majority of Central Luzon along with the Laguna Bay Area, CALABARZON, and Metro Manila, see Figure 4. According to the Manila Bay Sustainable Development Master Plan (MBSDMP) [14] of the National Economic and Development Authority (NEDA), MBA's environmental status has been

deteriorating over time for a variety of reasons, including environmental and biodiversity loss, and pollution. In addition, the natural environment has been altered by both natural and man-made situations over several decades, notably on the shore of Metro Manila. However, the project’s long-term viability became an issue given the questionable number of uncontested funds committed for the acquisition of dolomite and the lack of established data on how long the sand expenditure will stay on the coast given that the Bay’s substrate spans from mud to sand as stated in. Since then, the DENR has stated that two-meter-diameter geotubes were placed to contain the dolomite sand and prevent continuous erosion. However, the Mines and Geosciences Bureau (MGB) [15] explicates in its news brief that sand manually deposited on beaches would shift due to wind patterns, tidal variations, waves, currents, and other human and natural occurrences. Furthermore, any materials used for the Roxas Boulevard overlay, including the dolomite sand, will likely be spread either north to south near Paranaque or the Navotas, Obando, Meycauayan areas, or perhaps the Cavite region. It was also mentioned in [16] that the purpose of beach nourishment is too short-term prevent erosion rather than to halt it altogether. Given the probability of coastal erosion, DENR has not yet clarified how regularly the sand would be supplied.

The Ecosystems Research and Development Bureau (ERDB) of the DENR gave coastal communities in the National Capital Region (NCR), a high Coastal Vulnerability Index (CVI) rating in 2019, ranging from 62.01% to 74.70%. The CVI is a gauge of how fragile a site is to coastal erosion; hence, places with greater susceptibility are assigned a higher CVI. The ongoing alterations to the coastline, the lack of natural ecosystems like coral reefs, limited vegetation, and seagrass all contributed to Manila’s high CVI, which is 71.67% - 74.70%. The city’s coast is exposed to stronger waves and uneven sediment transport, which is especially susceptible to storm surges and coastal erosion, see Figure 5 for a reference of the area. The first mention of storm surges in Manila Bay dates to 1589 when they were said to have reached heights of four (4) to six (6) meters. The barrier along Roxas Boulevard was destroyed by typhoon “Pedring” in 2011, which led to the flooding of nearby regions. These damaging surge heights have been observed in previous years by tropical cyclones such as “Yoling,” which dates from 1970, and “Bebeng,” which dates from 1983 [17]. The MBA Sea level rises 13.2 millimeters each year, which is four times faster than the global average, according to the National Mapping and Resource Information Authority (NAMRIA) [18], a subsidiary of DENR. The effects of these have been seen in places close to Bulacan, Obando, Meycauayan, Marilao, Navotas, and Malabon. All these climatic and coastal threats call for a thorough and open evaluation of the program’s viability.

In response, The DENR [19] responds by saying that the dolomite sand won't just be "washed out". However, maintenance of the beach replenishment project is anticipated to be non-strategic and financially unsustainable given the site's susceptibility to coastal and climate challenges.

3.1.2 Barangay Pugalo – Pasol, Alcoy, Cebu, Philippines

Dolomite mining in Cebu began in 1981, and owing to demand, a continuation in the form of a Dolomite Mining Expansion Project was granted by the DENR last November 25, 2002, with Environmental Compliance Certificate (ECC) No. 0107-556-302. According to [20], this project is under the management of two big companies: Philippine Mining Service Corporation (PMSC), which oversees processing and production operations, and Dolomite Mining Corporation (DMC), which handles quarrying operations. These corporations are also liable for supplying the MBA Beach Project. Located in Central Visayas Region, with a total area of 66.53 square kilometres, this will serve as the main study area of the paper. The primary location of the existing mining site of Dolomite material is in Barangay Pugalo, Alcoy, Cebu, Philippines as shown in Figure 4. However, as seen in Figure 6, it is set to expand to a portion of its neighboring barangay named Pasol. The researcher decided to use both barangay as the main study area.

Figure 5:Philippine Map from Google Satellite:(MBA) (a) September 6, 2020 (b) December 23, 2020

Figure 6:Topographic map of Dolomite mining site

Though both barangays have ample land space, focusing on this area alone can cause a LULC error during classification due to its size. To ensure a higher accuracy assessment for the classification process, a vast space, or a higher number of pixels to classify is a must - to receive a positive or a reasonable classification. For these matters, refer to the classification section.

PMSC Dolomite Processing Facility [20] conducted a Socio-Economic and Perception Survey of the residents. The survey garnered five hundred eighty-eight hundred (588) respondents which unraveled the changes and experiences they have encountered since the mining activities started on their land and its effect on their community. The results show that in PARAMETER B natural disasters and calamities are more or less likely to occur 56% of the time. In contrast, man-made disasters and calamities are more or less likely to occur 49% of the time. When natural and man-made disasters and calamities are combined, they are more or less likely to occur 52% of the time.

The effects of the operation on the environment over the past ten years in the barangay revealed that the level of the quality of the water had changed little to not at all, noise has increased, the quality of the marine or seawater has not changed, the catch of fish has changed little to not at all, the quarry area has not changed, developed agricultural fields have shrunk, wildlife populations have declined, public green spaces have shrunk or declined, tree cover has improved, the quality of air has dropped, the corals and seaweeds have significantly reduced, settlement and dwelling areas have grown, and the scope of mangroves has expanded. Only 12% of the respondents reported that businesses or contractors such as PMSC were their primary source of income. In comparison, 4% of residents work for the government, 5% are employed by other private enterprises, and 21% are involved in farming and fishing. Notably, almost 25% reported not having a job, being jobless, or lacking a stable source of income . Regarding PMSC Dolomite Processing Operations' benefits, 98% of respondents benefit from the organization, as opposed to the two percent of respondents who disagree. 98% of these are attributable to the employment of locals and the provision of help for local livelihood, local tax income, and local employment. Following are the detrimental consequences of PMSC activities on the environment and community that respondents listed: With an eight percent chance of pollution and contamination of waterways, an 11% chance of accidents that could cause damage to life and property, a 16% chance of noise, a 22% of ground vibration due to blasting operations, an 18% threat to community health and safety, and a 24% flooding and landslides.

3.2 SCP: Land Use/Land Cover Analysis

The researcher acquired satellite data using SCP tools of the mining site area. PARAMETER A, undergone a classification process where land classes were identified automatically after creating training input with ROI and signature with assigned macro classes using the SCP using Maximum likelihood Classification and runs into a post-classification of accuracy assessment. This assessment calculates classification kappa value results that validate the LULC of the data produced by the SCP, in which the result shows a value greater than 80%, a reliable classification.

3.2.1 Satellite data and Region of Interest (RIO) & signature

The satellite imagery selection process paramount consideration on the acquisition date and its cloud cover. In this case, the satellite data chosen are from 2018 and 2020, see Figure 11, an estimated thirty-seven (37) months apart of spatial imagery, with the lowest cloud cover percentage of 0.85 and 1.70, respectively. With a minimum latitude and longitude of 8.96 and 123, respectively. A maximum latitude and longitude of 9.95 and 124, respectively. A zone/ path of 51PWL. A cloud covers higher than ten percent can cause a lower accuracy assessment. The selection of the color composite RGB values helps create ROI and spectral signature, which are the major components of the classification process. These ROIs are defined by their Macro class name and Class ID, which are each assigned to a specific land cover area; see Table 1. The display color composites of the selected RGB value are related to the band numbers in the band set, see Figure 7. Some RGB values used in the study are Band 3-2-1, 8-4-3, and 4-3-2. The researcher created several ROIs for each macro class.

3.2.2 Digital Elevation Modelling (DEM) and slope

The researcher also acquired the DEM and slope data masked to the shapefile of the study area from two separate plug-ins: the SRTM downloader plug-in and the Slope plug-in. The DEM imagery extracted from the plug-in is located at exactly 09°N, 123°E with a pixel size of ± 10, see Figure 9. These DEM data are then clipped to the shapefile of the study area. Moreover, using the same DEM clipped file, the researcher extracted a slope analysis of the same area using the Raster Terrain Analysis plug-in.

Table 1:MOLUSCE – Class Statistics

Figure 7: Study Area in its selected color composite (RGB)

Figure 8: LULC Classification – PARAMETER A (a) 2018 (b) 2020

Figure 9: STRM Downloader: DEM Data Located at 09°N, 123°E with (a) PARAMETER B (b) Shape File PARAMETER B DEM (c) Slope

3.2.3 Classification and accuracy assessment

The Maximum Likelihood approach is used throughout the classification phase to assume that the statistics compute the likelihood that a particular pixel belongs to a certain class and is normally distributed for each class in each band. After the classification of “PARAMETER A” in the year 2018 and 2020, See Figure 8. A second (2^nd) classification is executed where the researcher masks the area with the Barangay Pugal – Pasol, Alcoy, Cebu Philippines shapefile, see Figure 12 and Figure 14. This classification displays an output result named “PARAMETER B”, which only shows the 65.54 square kilometers area and a macro class of three (3) in the classification process. This classified data will be used for the paper’s change detection analysis. Due to the reclassification of the LULC “PARAMETER A,” some of the changes experienced in “PARAMETER B” are Macro class named Water bodies are removed, which causes the area of the Macro class Unclassified to increase. Thus, the Spectral Signature Plot and Details, and Spectral Distances show minimal changes. Nonetheless, both calculations of PARAMETER A and PARAMETER B show good classification results as the Maximum Likelihood output is ranging from 1.90 -2.0 which is described as each spectral distance between classes being fully distinct. Under the post-processing window of the same plug-in, the researcher used the generated LULC named “PARAMETER A” for the accuracy assessment of both years to assess the information's quality based on remotely sensed data. Comparing the ROIs and the classification is used to evaluate accuracy. With each value representing a class of comparison between the classification and reference file. With a 95% confidence interval, results show a 98.76% and 97.76% overall accuracy for 2018 and 2020, respectively, indicating that the reference data were proportioned and mapped correctly. The table also calculates the Kappa Classification value in which the result shows a value of 0.96, which is a perfect agreement result and a measure of how the categorization results correlate to values assigned by chance.

3.3 MOLUSCE: Change Detection and Simulation Analysis

The researcher input LULC PARAMETER B, DEM, and slope raster as spatial variables in the MOLUSCE tool. These data pass through Pearson's Correlation which shows the correlation from the data and tabulates each class's LULC area, its area changes from 2018 to 2020, and kappa validation with the 2030 and 2050 simulations.

3.3.1 Evaluation correlation and area changes

The researcher used Pearson’s Correlation as the main method for evaluating the relationship of geographic or spatial variables between the two LULC raster images. The output shows a strong correlation between LULC 2018 and 2020 with values of 0.95. An output value of 1 and 0.95 to the correlation of DEM to LULC 2018 and 2020, respectively. And a weak correlation with slope to LULC 2018 and 2020 with values of 0.25. The class statistics show the change in each category from the year 2018 to 2020, where there has been a relative increase in soil categories in the area and a definite decrease in vegetation and built-up areas for 2018 and 2020, see Table 1.

The plug-in created a transition matrix that shows the percentage of pixels shifting through one land use classification to another. Each non-diagonal value shows the magnitude of the transition from one class to another, whereas the diagonal elements display the amount of class stability. Values near 1 in diagonal entries indicate a category's stability. Table 2 shows vegetation is the most stable category probability following Soil and Built-up. Soil added 0.33 and 0.025 to vegetation and built-up, respectively. Vegetation contributed 0.008 to built-up and decreased transition probabilities with 0.168. Built-up decreased transition probabilities of 0.238 and 0.361 from Soil and Vegetation, respectively. This is also used to evaluate how the study area's temporal changes and shows how each LULC class was projected to change, the generated area change map for 2018 to 2020.

Table 2: MOLUSCE – transition matrix

Figure 10: “PARAMETER B” Layout Map: Brgy. Pugalo – Pasol, Alcoy, Cebu, Philippines

Figure 11: LULC 2018 Layout Map – Brgy. Pugalo – Pasol, Alcoy, Cebu, Philippines

Figure 12: 2018 3D Layout – PARAMETER B: (a) Satellite Image (b) LULC

Figure 13: LULC 2020 Layout Map – Brgy. Pugalo – Pasol, Alcoy, Cebu, Philippines

Figure 14: 2020 3D Layout – PARAMETER B: (a) Satellite Image (b) LULC

Figure 15: LULC Simulation (a) 2030 (b) 2050 Layout Map – PARAMETER B

Figure 16: 3D Layout – PARAMETER B LULC Simulation (a) 2030 (b) 2050

3.3.2 Transition potential model and cellular automata

Using the same inputs and generated data, the MOLUSCE tool generates a transition potential modelling for future simulation of LULC. Artificial Neural Network (ANN) – (Multilayer perception), one of the algorithms available in the plug-in, is chosen by the researcher to lead the prediction analysis [21] together with Pearson’s results from the spatial factors chosen for modelling calibration that have a strong correlation with the input of LULC [22]. The researcher used LULC data from 2018, 2020, and geographic factors to anticipate LULC for 2030 and 2050 using the CA-ANN approach for transitional potential modeling and prediction.

Following the standard settings of the ANN method, the sample mode is set at random with a number of 1000 samples, a neighborhood value of 1 pixel, a learning rate of 0.100, maximum iterations of 1000, 10 hidden layers, and a momentum value of 0.050. The minimum validation overall error was 0.01717, which was the minimum reached error on the validation set of samples. After training the neural network, the researcher obtained a validation kappa value of 0.95 and a delta overall accuracy value of ±0.00612 which contains the difference between the error of the minimum reached value and current error, the neural network learning curve from the ANN-MLP algorithm integration, where the iterations are shown on the x-axis as cited in [23].

The Monte Carlo method, commonly referred to as multiple probability simulation, governs the cellular automata simulation process. The created model is used to anticipate the likelihood of a range of outcomes when there is a possibility for random variables. This algorithm is a computerized mathematical tool that enables users to account for risk in forecasting and decision-making accurately. In this section, the researcher set the generation of the map to the year 2030 with an iteration number of one (1) for the first simulation and then proceeded to simulate a much more significant year gap for the second simulation dated at the year 2050 with two (2) iteration value, see Figure 16. After predicting the LULC for 2030 and 2050, the researcher validates the generated simulation to the LULC 2018 reference map with a value of five (5) validation iterations each. The LULC validation shows kappa statistics results for overall kappa with 0.83% of correctness with 90.65, kappa histogram of 0.91, and kappa location of 0.91, following simulated LULC simulation for both 2030 and 2050. Table 3 shows the result of the simulation process, which indicates an increase in Soil and a decrease in Vegetation and Built-up for the years 2018, 2020, 2030, and 2050. The results show minimal changes that might be due to the plug-in version used in the study and the spatial variable to simulate results and other software glitches that are further explained in the conclusion section.

Figure 17: (a) DEM (b) SLOPE Layout Map – PARAMETER B

Figure 18: 3D Layout Map – PARAMETER B: (a) DEM (b) Slope