Land Use / Land Cover Change Detection and Urban Sprawl Analysis in the Mount Makiling Forest Reserve Watersheds and Buffer Zone , Philippines

Ecologically Valuable Areas play an important role in providing ecosystem services, however, human activities such as land conversion and urban sprawl pose pressures and threats to these areas. The study assessed the land use/land cover and urban sprawl in the Mount Makiling Forest Reserve (MMFR) Watersheds and Buffer Zone from 1992 to 2015 using remote sensing and Geographic Information System (GIS). Results showed that the land use/cover within the MMFR buffer zone has changed from 1992 to 2015 with built-up areas increasing by 117% despite Proclamation 1257, s. 1998 which regulates human activities in the zone. Based on the Shannon entropy analysis the land development in the MMFR buffer zone tends to be dispersed and sprawling. However, when the magnitude of change of urban sprawl in the buffer zone from 2002 to 2015 was calculated, a decrease in the entropy value was observed which implies a compacting pattern as the human settlement in the buffer zone increases over time. Proclamation 1257, s. 1998 needs to be strengthened to protect MMFR and its buffer zone from further encroachment and pressure. Moreover, remote sensing and GIS proved to be useful tools for assessing urban sprawl in ecologically valuable areas such as MMFR.


Introduction
Ecologically valuable areas (e.g., forest reserve, national parks, and protected areas) [1] within the borders of urban areas provide many ecosystem services for urban [2] and neighboring population.These areas are especially valuable for the protection of biodiversity ( [3][4][5]; as cited in [6]).Despite their importance, ecologically valuable areas are threatened with the continuing urbanization [1].
The impacts of urbanization along the fringes of these areas include fragmentation of habitats, loss of rare species, edge effects, the introduction of alien species, degradation of the environment, reduced tourism potential and the loss of ecosystems services vital for the wellbeing and health of urban residents ( [1,7,8]; as cited in [9]).According to Mutuga [9], the harmful effects of urbanization is particularly caused by urban sprawl.Urban sprawl is the expansion of human settlement or development outside of an urban city center [10].As urban areas expand they take up more land and extends further to ecologically valuable areas [9].A disproportionate increase in the population growth near the boundaries of ecologically valuable areas may threaten their ability to conserve biodiversity [11].Urbanization near the borders of ecologically valuable areas creates a conflict between biodiversity protection and demands for housing and economic development.These conflicts are more frequent and intensive in areas with a greater number of rare and endangered species and in areas with severe population and development pressures ([12]; as cited in [9]).
In response to these problems, buffer areas have been designed to protect and enhance the conservation value of ecologically valuable areas.However, Rathore [13] notes that most forest buffers have already been reduced to a degraded status due to the massive pressures from adjacent settlements.McDonald, Kareiva, and Forman [14] observed the shrinking distance between protected areas and cities in some parts of the world.Considering these threats and pressures, periodic assessment [15] and monitoring of the urban development are important for conservation [16] and protection of ecologically valuable areas from further degradation, to evaluate natural resources and the effectiveness of management [15] and to make effective and appropriate measures using remote sensing and Geographic Information System (GIS) [16].Remote sensing is a cost-and-time efficient tool [17] that can be used for monitoring urban expansion [16,18] and the impact of human activities on the landscape [17] of forest reserves, national parks or protected areas.Remote sensing is a science and art of acquiring information about the Earth's surface and atmosphere using sensors designed to detect electromagnetic radiation either ground-based, aerial or spaceborne platforms [19,20].It is a reliable source of data, provides spatially consistent coverage of large areas with temporal frequency, and high spatial detail, which is useful for analyzing time-dependent phenomenon such as urban expansion [21,22].
This study assessed the land use/land cover change and urban sprawl in the Mount Makiling Forest Reserve (MMFR) Watersheds and its buffer zone from 1992 to 2015.Shannon's entropy was calculated to characterize and measure the relative trends in the dispersion or compaction of urban sprawl in the MMFR watersheds and its buffer zone.The magnitude of change of urban sprawl in the watersheds and its buffer zone over time was also analyzed.

Study Site
The Mount Makiling Forest Reserve (MMFR) is an important landscape and "training laboratory for the advancement of science and technical knowledge for the conservation and development of forests and associated ecosystems" [23].Mount Makiling is also an important catchment area for Laguna de Bay [24,25] providing a water supply for industrial, agricultural and residential sectors of the CALABARZON (Cavite, Laguna, Batangas, Rizal, and Quezon) region [26].MMFR is also classified as a national park ( [27]; as cited in [28]) and was declared as the 33rd ASEAN (Association of Southeast Asian Nations) Heritage Park on 3 October 2013 [29].By virtue of Republic Act 6967 (Section 3), the exclusive control, jurisdiction, and management of MMFR was handed to the University of the Philippines Los Baños (UPLB) [23].MMFR covers a total land area of 4224 hectares and lies between 14 • 06 to 14 • 15 north latitude and 121 • 09 to 121 • 16 east longitude [24].The forest reserve covers parts of the municipalities of Bay, Calamba City and Los Baños in the province of Laguna, and Santo Tomas in the province of Batangas [30].There are four watersheds in MMFR and these are the Molawin-Dampalit Watershed, Cambantoc Watershed, Tigbi Watershed and the Greater Sipit Watershed.In this study, they are collectively referred to as MMFR Watersheds.
Proclamation 1257 or "Designating Certain Areas Surrounding the Makiling Forest Reserve as Buffer Zone and Prescribing Guidelines Thereof" institutionalize the coverage of MFFR buffer zone.Accordingly, areas adjacent to the MMFR boundary with 18% slope and above (Figure 1. Map of Mount Makiling Forest Reserve Buffer Zone) forms part of the park's buffer zone.Proclamation 1257 was passed in 1998 "to guide the local government units (LGUs), entrepreneurs, developers, businessmen, landowners and other interested individuals on the use of the buffer zone surrounding the MMFR.It was also intended to avoid further encroachment on the forest reserve that would create disturbances or damage to the ecology of the area and reconcile socio-economic undertakings and environmental protection on the use of the buffer area to minimize negative off-site effects of land developments therein" [31].However, despite the Proclamation, MMFR is continually being threatened by the expansion of human settlements, agricultural areas, and urban and industrial zones [24].The urbanization along the borders of MMFR and the development of residential subdivisions on its lower slopes also adds pressure on the forest reserve [32], which may eventually have tremendous impacts on the communities depending on it [33].
The buffer zone covers portions of the alienable and disposable land of the municipalities of Bay, Los Baños, and Santo Tomas, and the city of Calamba.The population density in these areas averages as more than 1000 individuals per km 2 (Table 1.Population, Land Area and Population Density of Areas portions of which are covered by the Mount Makiling Forest Reserve).Among the municipalities, Bay is the smallest in terms of land area with only 42.66 km 2 .It is a second-class municipality with a population of 62,143 [34].The municipality of Los Baños is a first-class municipality with a land area of 54.22 km 2 and a population of 112,008 [34].Los Baños is also considered as a Special Science and Nature City of the Philippines [35].Calamba City is a first-class city and the most populated city in the province of Laguna.Calamba City is the center of commerce and industry in the province and has attracted many job seekers from nearby provinces which have led to a significant increase in its population [36].The population in Calamba City has increased from 389,377 in 2010 [37] to 454,486 in 2015 [34].Santo Tomas is also a first-class municipality and one of the oldest towns of Batangas [38].The municipality has a land area of 95.41 km 2 and a population of 179,884 [34].

Urban Sprawl Detection and Assessment
Various procedures were utilized to assess the urban sprawl in the buffer zone of MMFR, Philippines (Figure 2).These include data acquisition, image preprocessing, classification and analysis, post-processing, accuracy assessment, and Shannon's entropy calculation.

Spatial Data Acquisition
The study utilized information in vector and raster format including those acquired through remote sensing platforms.The different data used in the study are summarized in Table 2. Vector data include the physical boundary of MMFR watersheds, MMFR boundary, and MMFR buffer zone boundary which were requested from Makiling Center for Mountain Ecosystems (MCME) in UPLB.The Landsat satellite images were acquired from 26 January 1992 (Landsat Thematic Mapper (TM)) [40], 4 April 2002 (Landsat 7 Enhanced Thematic Mapper Plus (ETM+)) [41] and 29 July 2015 (Landsat 8 Operational Land Imager/Thermal Infrared Sensor (OLI/TIRS)) [39].These were downloaded from the USGS Earth Explorer website [39][40][41] in October 2017.The selected images were either cloud-free or present cloud cover of less than 10% to achieve more accurate results of the study.All vector and raster files were projected into Universal Transverse Mercator (UTM) System, Zone 51 North and World Geodetic System (WGS) 84 datum.

Image Preprocessing
Image preprocessing prepares and restores band data by recalculating Digital Number (DN) values and eliminating external and space-borne platforms errors.A standard radiometric calibration was applied to enhance DN values and extract the Top-of-Atmosphere (TOA) reflectance [42] of the satellite images covering the MMFR that were taken at different periods.The preprocessing of the Landsat images was conducted using the ENVI 5.2 software (Exelis Visual Information Solutions, Boulder, CO, USA).Only a subset containing the MMFR boundary was extracted and processed for a more focused analysis and better data handling.
The surface reflectance free from clouds and other atmospheric components was also determined which is vital in many remote sensing applications [19,43].Pan-sharpening of the Landsat 7 and 8 satellite images was completed using the Gram Schmidt method in the ENVI 5.2 software to maximize the available high-resolution panchromatic band (15-m) by combining it with a low-resolution multispectral image (30-m).

Classification and Analysis
The supervised classification method employs the Maximum Likelihood Classification (MLC) to classify and segment the pre-processed Landsat satellite images.MLC quantitatively analyze remote sensing data [19,[44][45][46] by classifying each cell based on the highest probability of being a member and spectral pattern represented in class mean vectors and covariance matrices [46].For Supervised classification, representative pixels or training data [46] for each land cover classes to be segmented were selected and utilized in the classification algorithm.Ten (10) training points per land cover class were identified and randomly selected from Google Earth [47] online [48,49] for Landsat 8 (Table 3).The land cover classes of these training points were determined using the Food and Agriculture Organization (FAO) land cover classes and description [50].Likewise, the land cover classes of these training points were identified from Landsat 5 (81 random points), Landsat 7 (104 random points) and Landsat 8 (72 random points) satellite images for change analysis.

Accuracy Assessment
Classification accuracy refers to the extent of correspondence between remotely sensed data and reference data [51].The generated land cover data from Landsat 8 were assessed through a selection of 30 sampling points [45] per class using the Google Earth [47] online [44,52].It was assumed that the 2015 classified map accuracy was enough to justify the overall accuracy of the study area in the different periods [53].The accuracy assessment of the classified images was performed using the ArcGIS 10.5 software (Environmental Systems Research Institute, Redlands, California, CA, USA) and the results were shown in a confusion matrix.The overall accuracy was generated which shows the overall quality of the data by dividing the total correct pixels by the total number of pixels.Kappa coefficient was also used to measure the classification accuracy.Kappa coefficient (K) is a discrete multivariate technique commonly used in accuracy assessment that provides a more unbiased estimate of overall agreement [48,54].The K interpretation values range from poor to excellent agreement ranging from 0 to 1.The closer the value of K to one (1.0) the more acceptable the classification [48].The error of commission and omission were also calculated.The error of commission represents pixels or the reliability of the map since it indicates how accurate the maps are to represent what is actually seen on the ground while the error of omission represents pixels that belong to the truth class but fail to be classified into the proper class.User and producer accuracy was also generated using the correct classified pixel per category (row and column) divided by the total pixel per class (row and column).

Classification Post Processing
Final data editing, consolidation, and integration are important procedures to enhance the quality of the final map [44].The post-classification change detection technique, which detects the nature, rate, and location of changes was used in this study.Vectorization using unsimplified polygons [44,55] of the classified images was executed to generate the total area and land cover change per class in hectare per satellite scene.Final data editing and integration were implemented for doubtful areas (e.g., clouds using secondary data such as land cover maps [56], Google Earth [47] and personal knowledge in the study site).The MMFR landscape patterns were determined based on the post edited final results of the image classification.The overlay analysis e.g., intersection of the final MMFR watersheds' land cover maps and the MMFR buffer zone was conducted to extract the landcover patterns inside the MMFR buffer zone for 1992, 2002 and 2015.
Using an overlay procedure in GIS, spatial changes in the Land Use/Land Cover (LULC) of the Mount Makiling Forest Reserve watersheds and buffer zone in two phases (i.e., [1992][1993][1994][1995][1996][1997][1998][1999][2000][2001][2002][2002][2003][2004][2005][2006][2007][2008][2009][2010][2011][2012][2013][2014][2015] were determined.This technique resulted in a cross-tabulation matrix describing the main types of change in the study area.The cross-tabulation obtained was analyzed on a pixel-by-pixel basis and this procedure consequently facilitated the determination of the number of conversions from a particular class to other land use classes and their corresponding area and weight of change (in percent).As a result, a new thematic map stating different combinations of "from to" change classes were produced.

Calculation of Shannon's Entropy
Quantifying urban sprawl is crucial to monitor urban development especially in areas adjacent to Ecologically Valuable Areas like MMFR.Spatial metrics are numeric measurements that quantify spatial structure and pattern of a landscape ( [57,58]; as cited in [59]).Shannon's entropy is a spatial metric commonly used to efficiently identify and characterize the degree of spatial concentration or dispersion in a geographic area [60].The Shannon's entropy (H n ) was computed to analyze the degree of dispersion or compactness of the spatial expansion in the MMFR watersheds and buffer zone.Shannon's entropy (H n ) can be used to measure the degree of spatial concentration or dispersion of a geographical variable (X i ), e.g., built-up area [60][61][62].The value of entropy ranges from 0 to log(n) [60,63].The value of 0 [63] or closer to 0 indicates that the distribution of the built-up areas is compact while values closer to log(n) indicate that the distribution of built-up areas is dispersed [60,63].The high values of entropy indicate the occurrence of sprawl.The halfway mark of the log(n) value is generally considered as the threshold and if the value crosses the threshold the area is considered sprawling [63].
The land cover classes from the three periods (1992, 2002, and 2015) were re-classified into two categories: built-up and non-built-up areas.In the analysis, the MMFR watersheds' boundary was divided into 500 m by 500 m grids (Figure 3) [64] in ArcGIS 10.5 software using the Fishnet Tool.To overlay the fishnet with the MMFR watersheds and the land cover data per period, the overlay toolset (i.e., intersect) was used.Afterward, the land cover area for each class in the grids was calculated using the field calculator.The output table was converted into an Excel matrix to complete the computation of the entropy values [65].The same procedure was used for the calculation of the entropy values for the MMFR buffer zone.The computation of the entropy was done using this formula [60]: where p i is the probability or proportion of occurrence of a phenomenon in the ith zone, and thus, is given by: where X i is the area of built-up at the ith zone and n is the total number of zones.Moreover, the relative entropy was measured to scale the entropy value into a value that ranges from 0 to 1.The formula for relative entropy (H' n ) [62] is: To determine the magnitude of change of urban sprawl in the study area over time, the difference of entropy values between time t and t + 1 was measured [60,66].A decrease in entropy over time shows that the urban development is compact while an increase in entropy shows a dispersed urban development [65].The formula [60] for measuring the difference is: where ∆H is the change in the entropy between time (t + 1) and (t).

Land Classification Change Analysis
The study assessed urban sprawl in the buffer zone for MMFR from 1992 to 2015.The land cover data during the years 1992, 2002 and 2015 were initially examined to describe the historical land cover change in the MMFR watersheds and its buffer zone.The land cover classifications of MMFR watersheds and buffer zone include forest, agricultural areas, built-up, and water (Figures 4  and 5).As shown in Table 4, the forest area in the MMFR watersheds had significantly decreased over time.The forest area of the MMFR watersheds decreased by as much 35% from 7884.90 has in 1992 to 5149.76 has in 2015.Moreover, there were only 1235.88 hectares of built-up areas in 1992 which more than doubled to 3368.63 hectares by 2015.During these periods, the agricultural area also increased from 5526.36 in 1992 to 6130.60 hectares in 2015.The decrease in the forest area and increase in agricultural areas and built-up areas in the MMFR watersheds can be attributed to the growing population in Bay, Los Baños, Calamba, and Sto.Tomas (Table 1).Population growth may consequently increase the need for a higher income to meet the household basic needs and in some cases, landowners may decide to give-off and trade their lands to support the growing needs of the family.Gradinaru et al. [67] stated that the landowner's anticipation of higher profits from the sale of land with urbanization may lead to the abandonment of agricultural land to give way to built-up development.

Land Cover Change Detection
Table 5 shows the land cover change in the buffer zone for MMFR.Results show that the forest area in the buffer zone decreased by 28% from 1193 hectares in 1992 to 860.15 hectares in 2015 (Figures 6 and 7).The significant increase (96%) in the agricultural area in the buffer zone occurred between 1992 and 2002 when the agricultural land area increased from 225.06 ha to 441.98 ha.During these periods, 95% (from 22.67 has to 44.18 has) increased in the built-up area was also observed in the buffer zone.Panaguiton, Eustaquio, and Campo [68] stated that the increasing settlements near the forest reserve constantly threaten its stability and sustainability.They further explained that the slow but steady encroachment of communities near the borders of the forest reserve could have resulted from unclear government policies dealing with these communities and occupants.In 1998, Proclamation 1257 was issued as a response to the increasing settlements and to limit the human activities in the buffer zone of MMFR.Despite the Proclamation 1257, s. 1998, the built-up area still increased by 11% and the forest area decreased by 10% from 2002 to 2015.The authors believe that the government needs to "strengthen" Proclamation 1257 and not just a "law enforcement" considering that a legislative act is required for this.

Classification Accuracy Assessment
The indispensable procedure after various image processing techniques is the accuracy assessment [19,44,46].Accuracy refers to the correctness of data [69].It measures the agreement between the reference data and the output of the image classification of unknown quality [44,46,69,70].It is ought to be done since it gives the degree of confidence and reliability of the data output [19,46] In essence, it aims to determine the percentage of pixel values that have correctly classified against the reference data.Many researchers incorporated the accuracy assessment in their studies [44,52,69,[71][72][73][74].
It is important and ought to be done in the conduct of image classification for the following reasons [70,75]: 1.
The legal standing of maps and reports derived from remotely-sensed data.

2.
The operational use of such data for decision making e.g., watershed management.

3.
The validity as input for scientific research.
In this study, a confusion matrix was utilized based on 93 sampling points with the help of ArcGIS 10.4.1 [44][45][46]74].After image classification, at least 30 samples per class [45] were generated using the random sampling to target at least 90% overall accuracy [46].The sampling points were the input to the accuracy assessment of the Landsat 8 classified image.Jensen [45] emphasized that the sampling points and signature data should be collected during the same period.Due to the time and resources constraints, the accuracy of the classified image was based on the available higher resolution online maps (e.g., Google Earth) [44,52,55] and other relevant secondary information available such as land cover maps from National Mapping and Resource Information Authority (NAMRIA) [56].The matrix (Table 6) tallied the correctness of the classification generated through the application of sampling points [44][45][46]48,70,74].The matrix was constructed based on the reference data as represented in the columns and compared against to the classified data in rows (Table 6).The major diagonal of the table shows the agreement between these two datasets.
Land cover classes such as forest, agricultural areas, built-up and water were extracted from Landsat TM, Landsat ETM+, and Landsat OLI/TIRS.The kappa and overall accuracy values of 0.91 and 93.55% were generated for the 2015 land cover map.Table 6 shows the confusion matrix with the user, producer, overall accuracies and kappa coefficient.The overall assessment result suggests that the classified land cover of MMFR Watersheds and Buffer Zone has very high accuracy and acceptability.The accuracy of at least 85% is required for satisfactory use of land cover data for resource management [52,76].Due to the absence of higher resolution images for 1992 and 2002, only data in 2015 were assessed accordingly [52,53].The assumption was 2015 accuracy was sufficient to shed light on the overall accuracy of the land cover change mapping procedures [53] conducted in the study.Thus, the statistical values of the following parameters were created to evaluate the quality of the 2015 classified image [44][45][46]70,74]: Overall accuracy-shows the overall quality of the data and computed as: Overall Accuracy = Total Correct (sum of the major diagonal)/Total number of pixels in the matrix 2.
Producer's accuracy-measures the error of omission or samples that are omitted from the correct classification since it indicates the probability that a reference sample will be correctly classified.
It is called the producer's accuracy since the analyst is concerned in mapping the Earth surface correctly.Computed as: Producer's Accuracy = Total number of correct pixels in a category/Total number of pixels of that category as derived from the reference data (column total)

3.
User's Accuracy-measures the error of commission or the reliability of the map since it indicates how accurate the maps to represent what actually seen on the ground.Computed as: User's Accuracy = Total number of pixels in a category/Total number of pixels that were actually classified in that category (row total) 4. Kappa Coefficient (K)-provides a more unbiased estimate of the overall agreement [48,54].The K interpretation values range from poor to excellent agreement ranging from 0 to 1.The closer the value of K to one (1.0) the more acceptable the classification [48] and computed using the formula [45]: where: The confusion matrix of the classified image (Landsat 8) was exported into table format (Excel file) as the basis of the quality of the classification [44] and analyzed the needed parameters (e.g., Overall Accuracy) [77].

Urban Sprawl in the MMFR Watersheds
The Shannon's entropy (H n ) was computed to analyze the degree of dispersion or compactness of the spatial expansion of built-up areas in the MMFR watersheds and buffer zone from the three periods (1992, 2002, and 2015).The results of Shannon's Entropy analysis are summarized and illustrated in Table 7 and Figure 8, respectively.Table 7 shows that the highest value of Shannon's entropy [Log(670)] is 2.83.The entropy results obtained for the three study periods were 2.34, 2.50, and 2.49, respectively.Moreover, the entropy values for each period is greater than the threshold value [Log(n)/2] of 1.41.Several authors stated that when the entropy value is greater than the threshold value the land development is dispersed [60,63] and sprawling [63].
However, when the magnitude of change of urban sprawl over time was calculated, a different pattern of urban sprawl in the MMFR Watersheds was observed.The results showed a difference in entropy values of 0.16 from 1992 to 2002 and −0.01 from 2002 to 2015.The decrease of entropy values was also observed in the study of Chong [68], where she explained that a decrease in entropy values over time implies a reduction on the proportion of low-intensity urban land cover across the zones resulting from the increase of high or medium intensity urban land cover or conversion to a different land cover type.An urban development is dispersed when there is an increase in entropy values over time and compact when there is a decrease in entropy values [68].

Urban Sprawl in the MMFR Buffer Zone
The results of Shannon's entropy (H n ) for MMFR buffer zone is presented in Table 8 and Figure 9.A similar pattern of urban development was observed in the buffer zone of MMFR.The entropy values for the three periods is 1.29, 1.47, and 1.46, respectively.All H n values were greater than the halfway mark of the log(n) value (1.085) which implies that the land development in the zone is dispersed and sprawling.The magnitude of change of urban sprawl over time in the buffer zone of MMFR was likewise analyzed.The magnitude of change from 1992 to 2002, that is prior to the enactment of Proclamation 1257, s. 1998 up to its early years of implementation, showed a dispersed pattern of urban sprawl with a difference value of 0.18.However, from 2002 to 2015, a difference value of −0.01 was observed, which implies a compacting pattern as the human settlement in the zone increases over time.Chong [65] stated that a decrease in entropy values over time indicates that the land development in an area is becoming compact or concentrated while an increase in entropy values indicate that it is developing in a dispersed manner.

Conclusions
The study analyzed the urban sprawl in the MMFR watersheds and its buffer zone using time series assessment based on Landsat data.The land cover was characterized through a supervised classification approach.The Shannon's entropy (H n ) was also calculated to describe the spatial patterns of urban sprawl in the area.Results showed that the land cover within the buffer zone of MMFR changes from 1992 to 2015 with a 117% increase in built-up areas despite the existence of policies concerning the management of the zone.In 2015, the total built-up area in the buffer zone was 49.22 ha representing 3.4% of the total area of the buffer zone.From 2002 to 2015, a few years after the enactment of Proclamation 1257, s. 1998, the magnitude of change of urban sprawl in the buffer

Figure 2 .
Figure 2. Flowchart for urban sprawl detection and assessment.

Figure 4 .
Figure 4. Land Cover Maps of MMFR Watersheds in 1992, 2002 and 2015 generated using the Supervised Classification Approach.

Figure 5 .
Figure 5. Land Cover Maps of MMFR Buffer Zone in 1992, 2002 and 2015 generated using the Supervised Classification Approach.

Figure 6 .
Figure 6.Change Maps of the MMFR Watersheds from 1992-2002 and 2002-2015.The two insets show the same zoomed area for its land cover change.

Figure 7 .
Figure 7. Change Maps of MMFR Buffer Zone from 1992-2002 and 2002-2015.The two insets show the same zoomed area for its land cover change.

Table 1 .
[34]lation, Land Area and Population Density of Areas portions of which are covered by the Mount Makiling Forest Reserve[34].

Table 2 .
Different data used in the study.

Table 4 .
Land Cover and Percentage of Land Cover Change in Mount Makiling Forest Reserve Watersheds in three periods: 1992, 2002, and 2015.

Table 5 .
Land cover and Percentage of Land Cover Change in Buffer Zone for Mount Makiling Forest Reserve Watersheds in three periods: 1992, 2002, and 2015.

Table 6 .
Confusion matrix of the 2015 classification maps derived from the Landsat OLI data of MMFR Watersheds.

Table 7 .
Shannon's entropy values for three years in the study area (MMFR Watersheds).

Table 8 .
Shannon's entropy values for three years in the study area (MMFR buffer zone).