| 17 Geoplanning Vol 5, No. 1, 2018, 17-34 Journal of Geomatics and Planning E-ISSN: 2355-6544 http://ejournal.undip.ac.id/index.php/geoplanning doi: 10.14710/geoplanning.5.1.17-34 MONITORING AND PREDICTING LAND USE-LAND COVER (LULC) CHANGES WITHIN AND AROUND KRAU WILDLIFE RESERVE (KWR) PROTECTED AREA IN MALAYSIA USING MULTI-TEMPORAL LANDSAT DATA J. Gambo a,b, H. Z. M. Shafri a , N. S. N. Shaharum a, F. A. Z. Abidin c, M. T. A. Rahman c a Department of Civil Engineering and Geospatial Information Science Research Centre (GISRC), Universiti Putra Malaysia (UPM), Malaysia b School of General Studies, Binyaminu Usman Polytechnic, Nigeria c Department of Wildlife and National Parks (DWNP), Malaysia Abstract: Natural and anthropogenic activities surrounding a Protected Area (PA) may cause its natural area to change in terms of Land Use-Land Cover (LULC). Thus, there is need of environmental change monitoring within and around PA because of its significant values to ecosystem at conservation scales. Effects and influences of local community within and around PA turn into the major problems for natural resource and conservations management as well as environmental impact assessment. Ascertaining the complex interface in relations to changes and its driving factors over period of time within and around PA is significant in order to predict future LULC changes, build alternative scenarios and serve as tools for decision making. The main objective of this work was to evaluate temporal change detection and prediction of LULC as well as the trends of changes from 1989 to 2016 within and around Krau Wildlife Reserve (KWR). The cloud issues were mitigated by producing cloud free image and object-based image analysis (OBIA) was adopted after a comparison with pixel-based analysis for overall accuracy and kappa statistics. The comparison of classified maps had produced a satisfactory results of overall accuracies of 91%, 86% and 90% for 1989, 2004 and 2016 respectively. The natural/dense forest between periods of 1989-2016 was decreased whereas built-up and agricultural/sparse forest were increased. The simulation model of Land Change Modeler (LCM) was utilized with digital elevation model (DEM) and past LULC maps to project future LULC pattern using Markov chain. The predicted map trend showed an increase of dense forest converted to agricultural/sparse forest in the north-western, and urban/built-up in east-southern part of KWR. The study is important for the conservation of habitat species and monitoring the current status of the KWR. Copyright © 2018 GJGP-UNDIP This open access article is distributed under a Creative Commons Attribution (CC-BY-NC-SA) 4.0 International license. Gambo, J., et al. (2018). Monitoring and Predicting Land Use-Land Cover (LULC) Changes within and around Krau Wildlife Reserve (KWR) Protected Area in Malaysia using Multi-Temporal Landsat Data. Geoplanning: Journal of Geomatics and Planning, 5(1), 17-34. doi: 10.14710/geoplanning.5.1.17-34 1. INTRODUCTION Protected areas (PAs) represent a massive investment around the world both at national and international levels in looking after our environment. Awareness and dialogue about protection and conservation of these environmentally sensitive areas be well informed and shared understanding among all the beneficiaries involved within and around. Natural areas have been effected greatly in many world locations by natural and human activities. Thus, monitoring, detecting and forecasting of land features modifications are significant for sustainable management, biodiversity, conservation, and development of PA (Bozkaya et al., 2015). Geographic spaces which, because of their particular environmental values for conservation purposes, deserved a special forms of safety ranging from total closure, except for protection purposes, to various forms of intervention required to maintain or restore habitats, to direct human use, remove OPEN ACCESS Article Info: Received: 20 August 2017 in revised form: 06 March 2018 Accepted: 30 April 2018 Available Online: 30 April 2018 Keywords: LULC, OBIA, Protected Area, Krau Wildlife Reserve, Land Change Modeler Corresponding Author: Helmi Zulhaidi Mohd Shafri Coordinator of Remote Sensing and GIS programme Department of Civil Engineering Faculty of Engineering Universiti Putra Malaysia (UPM) 43400 Serdang, Selangor, Malaysia Email: hzms04@gmail.com https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 mailto:hzms04@gmail.com Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 18 | invasive species, re-introduce extirpated spaces, or facilitate visitation by scientists or the public for purposes of research, monitoring and mapping, recreation and education. Tourism was also considered to consist of the facilitation of recreational visitation: availability of food services, guide to access roads, accommodation, services, and water supply and regular sanitation is referred to PA (Dudley & Stolton, 2008). PA is also clearly defined by International Union for Conservation of Nature (IUCN) 2008, as a geographical space, dedicated, acknowledged and managed, through legal or other effective means, to achieve the long-term future conservation of nature with associated ecosystem services and cultural values or an area of land and/or sea especially dedicated to the protection and maintenance of biological diversity, and of natural and associated cultural resources, and managed through legal or other effective means (Dudley & Stolton, 2008). UNEP (2004) described PA as an area of land or sea especially dedicated to the protection and maintenance of biological diversity and conservation of natural and related cultural resources, managed through legal or other positive means. Natural forests are the unique land cover types of ecosystems especially in most of the wildlife reserve that provide important ecological services for habitat species within the reserve area. Rural community along the PA relies on forest resources for their livelihood (Despot Belmonte & Bieberstein, 2016). Frequent changes in PA has led to degradation and fragmentation of wildlife habitat. In Malaysia reserved areas have been traditionally established with creation of Chior Wildlife Reserve in 1903. However, based on the Master List there are 490 protected areas listed for Malaysia: 271 PAs for Peninsular Malaysia, 173 for Sabah (including 3 in FT Labuan) and 46 for Sarawak. These encompass land area (3,510,239 ha) with a total size of 4,586,273 ha (Table 1) (Interim Master List of Protected Area in Malaysia). Some PAs are administered by the department within the Ministry of Natural Resources and Environment (NRE) and others are administered and managed at the states and NGOs level. According to the NRE, terrestrial PAs currently cover over 1.8 million ha in Peninsular Malaysia. They can be divided into four legal categories: Areas reserved for a public purpose under the land laws; Permanent reserved forests (PRFs) under the forestry laws; National parks and state parks under the parks’ laws; and Sanctuaries or reserves under the wildlife laws (UNDP, 2012) KWR is one of the threatened PA in Malaysia due to frequent illegal logging (Lin, 2016) of important forest tree which alter the condition of natural forest over a decade (Ahmad, Abdullah, & Jaafar, 2012). Geospatial technologies have potentials for mapping changes in Pas (Willis, 2015), environmentally sensitive areas, and biosphere reserve area. Mapping and monitoring PAs and their surrounding areas at both local and regional scales are crucial given that the vulnerability to anthropogenic activities, including climatic change, and important for conservation and biodiversity management. Monitoring using geospatial technologies and field information, can play a significant role in developing baselines for understanding condition of habitats and related species diversity (Bush et al., 2017) as well as measures the gain and losses, associated with specific activities around the protected area (Nagendra et al., 2013). The challenges for detecting and monitoring of LULC using optical remote sensing data especially in tropical region were the cloud cover and haze (Nagendra et al., 2013), but this study utilized image patching using multi-date Landsat data. KWR, for example has been particularly input for reason such as recent report about unstop illegal logging and mining around KWR, the use of geospatial data has been used before (Ahmad et al., 2012), but limited and outdated information about the loss future prediction and less number of LULC classes by (Ahmad et al., 2012), also understanding the current status of natural forest and conditions of habitat species. This could then be of use as updated study of and also utilized the products for proper land use planning, biodiversity and conservation management of KWR based on recent report of illegal logging in KWR (Norawi, 2017) and also illegal gold mining Lakum Forest Reserved within the river streams Sg. Teris around KWR. Thus, this study aims to provide the most up-to-date study on the status of KWR and evaluate changes that have occurred over the period and prediction (Mishra, Rai, & Mohan, 2014; Kumar et al., 2015; Reveshty, 2011). Technique of OBIA will be investigated and compared with traditional pixel-based method to ascertain on the proper method to be used in generating the required information. https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 | 19 Table 1. Terrestrial and Marine PA Coverage in Malaysia State Terrestrial (ha) Marine (ha) Johor 235,407.80 61,869.80 Kedah 2.00 10,720.10 Kelantan 127,946.70 0.00 Melaka 106.80 2,401.80 Negeri Sembilan 57,323.90 1,893.80 Pahang 855,160.90 55,800.80 Perak 286,673.20 0.00 Perlis 4,441.20 0.00 Penang 1,414.10 1,339.30 Selangor 106,673.10 0.00 Terengganu 139,844.10 110,931.40 Federal Territory 0.00 0.00 Kuala Lumpur 156.30 0.00 Putrajaya 0.00 0.00 Labuan 0.00 9,288.30 Sabah 1,555,022.70 118,768.50 Sarawak 616,288.20 226,914.00 Malaysia 3,986,345.60 599,927.80 2. DATA AND METHODS Study area base map was prepared by using subset of landsat satellite image of 2016 (Islam et al., 2018). The systematic workflow indicate entirely steps conducted and adopted throughout the research period, starting with selection of study area, data collection (landsat 1989, 2004 and 2016). Data analysis were also applied to both imageries, initially from image correction (atmospheric correction and geometric correction). Due to presence of cloud cover in 2016 image, image patching was also applied using SmartGEO Fill Tool. The detailed description of the entire methodology workflow was done in sub headings of data and methods as included; data and preprocessing, dealing with cloud and image patching, image classification and change detection approaches. The general systematic flowchart of the steps conducted in this study is shown in Figure 1. 2.1. Study Area The research area of interest is known as Krau Wildlife Reserve (KWR), located nearby Mountain Benom with a streams/tributaries drained to Lompat, Teris and Krau River in the district of Temerloh and Jerantut of Pahang, Malaysia. It is geographically bounded to the south-east of Taman Negara forest which is the largest natural forest in. KWR covers approximately 62,395 ha as the largest wildlife reserve which makes it the third largest PA in Peninsular Malaysia (Figure 2) with elevation of 2,107 metres ranging from the top of Benom Mountain in Kuala Lompat, to 43 metres to the reserve area. The Krau Wildlife Reserve Office, Institute of Biodiversity and National Elephant Conservation Centre located in the southern part with an entrance through Lanchang while Jenderak Wildlife Conservation Centre is located in the east of the reserve bordering Felda Jenderak Selatan (DANCED & Jabatan Perlindungan Hidupan Liar dan Taman Negara, 2001). KWR was managed and controlled by the Department of Wildlife and National Parks. Because of its diversity of the landscape and biodiversity fullness within such a compact area makes KWR a unique centre of habitat species, natural forest with different flora and fauna (DANCED & Jabatan Perlindungan Hidupan Liar dan Taman Negara, 2001). In International Union for Conservation of Nature (IUCN) Management https://doi.org/10.14710/geoplanning.5.1.17-34 https://en.wikipedia.org/wiki/Department_of_Wildlife_and_National_Parks Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 20 | Categories for protected area, KWR in list Category 1a: Strict Nature Reserve/Wilderness Area its meant protected area mainly for science or wilderness protection. Nevertheless, activities surrounding the area lead KWR to faces many problems related to conservation purposes including encroachment and conversion of natural forest, illegal harvesting of non-timber products, degazettment, and over-hunting (Ahmad et al., 2012). 1989 Imageries 2004 Imageries 2016 Imageries Mosaic Image Patching Geo Smart Fill Cloud Free Images Subset 1989 Subset 2016Subset 2004 Image Segmentations Segmented Image for 1989 Segmented Image for 2004 Segmented Image for 2016 Training Object Selection for Segmented Images (1989,2004 and 2016) Nearest Neighbor Image Classification for Segmented Images (1989,2004 and 2016) Time Series Classified Maps (1989,2004 and 2016) Reference points for Subset (1989.2004,2016) Accuracy Assessment (1989,2004 and 2016) RESULTS (MAPS, TABLES & FIGURES Change Analysis (1989-2004 and 2004-2016) predictions (2028 and 2040) Landsat TM Imagery (1989) Landsat TM Imagery (2004) Landsat OLI Imagery (2016) Atmospheric and Geometric Corrections STUDY AREA (Krau Willdlife Reserve) Figure 1. Methodology Workflow Adopted in this Study https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 | 21 Figure 2. The KWR Location in Peninsular Malaysia 2.2. Datasets The downloaded Landsat 8 OLI and Landsat TM Level-I imagery of the area from U.S Geological Survey (USGS) website (http://earthexplorer.usgs.gov/) are shown in Table 2. To analyze LULC changes in on a yearly basis, the Landsat imageries of the year 1989, 2004 and 2016 were obtained with an interval of 15 and 12 years. Due to the heavy rainfall and cloud cover in Malaysia as tropical region, obtaining optical sensor imageries of 10% cloud free is quite difficult, for this reason in this study four different scene of 2016 Landsat 8 OLI and two different scene of 2004 Landsat TM were utilized and produced 0-1% cloud free image of study area within the same dry seasons between June, July and August, the preprocessing operation started from atmospheric correction and image patching of Landsat TM 2004 and Landsat OLI 2016 due to presence of cloud in area of interest of this study using PCI Geomatics Smart Geofill tools. Using region of interest (ROI) created in Google Earth Pro, the area of interest was subset following the image patching. The completion of pre-processing stages, the study progressed to the image analysis process, first stage image segmentation was conducted by object-based image classification. The segmentation algorithm of multiresolution in eCognition Developer 9.0 applied to each individual image for generating image objects. Table 2. Satellites Data Used in this Study SENSOR ID DATE ACQUIRED Time Path/Row Resolution TM 4 LANDSAT 6/15/1989 3:01:32 127/057 30 TM 5 LANDSAT 6/16/1989 2:49:51 126/057 30 TM 5 LANDSAT 8/2/2004 3:07:27 126/057 30 TM 5 LANDSAT 7/18/2004 3:10:35 127/057 30 OLI 8 LANDSAT 7/3/2016 3:28:03 127/057 30 OLI 8 LANDSAT 6/26/2016 3:21:49 126/057 30 OLI 8 LANDSAT 6/26/2016 3:22:12 126/058 30 OLI 8 LANDSAT 6/1/2016 3:27:51 127/057 30 The fusing of neighboring segments together to ensure a heterogeneity threshold is stretched to one-pixel image segment within segmentation algorithm. Nearest Neighbor image classification was applied by following three stages that included creating class hierarchy, training data sets, and training samples were https://doi.org/10.14710/geoplanning.5.1.17-34 http://earthexplorer.usgs.gov/ Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 22 | taken for each LULC type to be classified in the image. Testing samples were randomly selected from three different Landsat image of the study (Waiyasusri, Yumuang, & Chotpantarat, 2016). With reference to our previous training samples (Hackman, Gong, & Wang, 2017).Second stage was the classified images of different dates were used for change detection and landcover projection in Land Change Modeler (LCM) in Idrisi Selva software. Five different LULC types comprising of water, dense forests, urban/built-up, agriculture/sparse forest and bare soil were recognized as the final classes in Table 3. Table 3. Level LULC within and around KWR, as Derived from Landsat Satellitesin 1989, 2004 and 2016 Land cover classes Descriptions Water River, Dam, Pond, Stream, Reservoir, Pool etc. Dense Forest Natural and dense Forest Built-up Building, Road etc. Bare soil Open land, Harvested land etc. Agriculture/Sparse Forest Rubber plantation, oil palm plantation, Banana plantation, shrubs/Grass, Cropland, Orchard and low density forest etc. 2.3. Dealing with Cloud and Image Patching It was possible to utilize uncorrupted Landsat scenes within the same season with corrupted one. In each scene, some pixels in the areas with high altitude regions of KWR had been covered by clouds, haze, predominantly in the northern part of the study area. It was easy to identify scenes with cloudy areas because of the available cloud percentage information. However, we had to visually look for areas where pixels had been corrupted by haze or live fires and the resulting thick smoke using false composites (Hackman et al., 2017). We used Smart Geofill tool which allowed user to copy a specified area of an image layer with non-cloud effects, make changes to it and then paste the selection to another layer of the image with cloud effects using color balancing method either overlap or histogram trim (Smart GeoFill Geomatica 2015 Tutorial), and can also adjust settings for color balance, blend width, contrast, and brightness of the selected area to enhance or adjust its appearance in the destination layer (Hruby et al., 2016). Meanwhile (Hackman et al., 2017) has masked out areas with corrupted pixels before the classification which had an extremely thick layer of cloud. In this study image patching was applied for all available 2004 and 2016 multi-date Landsat data images with heavily cloudy (Figure 3) showed the 2004 multi-date landsat before image patching and after applied using Smart GeoFill). Nonetheless, except for the landsat scenes taken on June 15/16, 1989. Figure 3. (A) Before Image Patching, (B) After Image Patching A B https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 | 23 2.4. Image Classification and Change Detection Approaches The main aim of image classification was to automatically categorize all pixels in an image into land cover classes (Figures 4). The classification legend was made based on spectral characteristics. Image segmentation was conducted by object-based image classification. Traditionally, pixel-based image analysis was utilized for image classification on both low and moderate resolution satellites data, OBIA incorporates not only the spectral information, but also included the shape, size, spatial, texture and contextual of the data. OBIA merged both spatial and spectral information about the features to extract land use for specific objects (Kindu et al., 2013), which grouping many pixels in to one image objects during segmentation stage to avoid the salt-pepper effects (Desclée, Bogaert, & Defourny, 2006). While pixel-based directly focus to one single image objects. In recent research, many studies utilized these technologies such as OBIA in change detection at different scales, with both low and high resolution satellites sensor such as landsat data (Dutta, Reddy, Sharma, & Jha, 2016; Waiyasusri et al., 2016; Balaji, Geetha, & Soman, 2016; Ranjan et al., 2016) and data like ALOS (AVNIR2) was used for land use changes (Munthali & Murayama, 2011) and the OBIA approach the same of work by Zhang et al (2017) maintained OBIA approach for change analysis in florida everlades water conservation area using landsat data. The object-based image analysis showed the expansion/reduction of land use types when applied the classified image in change detection algorithms as done by Son et al. (2015). The multiresolution segmentation algorithm in eCognition Developer 9.0 was applied to generate image objects for each individual image. The segmentation algorithm starts with making homogeneous object clusters with one-pixel image segment, and considerably merges neighboring segments together until a heterogeneity threshold is reached. The heterogeneity threshold identification depended on user-defined scale parameter, as well as the shape and compactness weights. The scale of the segmentation determined the quality of segmentation, and classification. The image segmentation is scale-dependent. For this study different scales were applied before the selection of appropriate classification. Figure 4 shows two segmented results with different scales assigned. Image (A) has the scale parameter of 40, and maintained the default values of shape 0.1 and compactness equal to 0.5. For image (B) we tried to change different segmentation scale from default of shape, compactness and scale parameter to 70, but the segment between the features have an overlap. In this study we adopted the image (A) segmentation scale because each land use and land cover was segmented properly. Figure 4. Segmentation Testing Results using Different Scale Parameters A B 0.1, 05 and 70 0.1, 05 and 70 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 24 | Nearest Neighbor Image classification was done by following three stages that included creating class hierarchy, training data sets, and accuracy assessment. Training samples were taken for each LULC type to be classified in the image. Following the adaptation of Nearest Neighbor object based image analysis in this study, a pixel-based image analysis was applied to 2016 OLI landsat image using support vector machine (SVM) as to compare with OBIA approach. The comparison in terms of the overall and kappa statistics showed good results for both OBIA and pixel-based techniques (for OBIA 90%, 0.87 and SVM 98%, 0.87). However, the visualization of two classified image indicated that there were misclassifications between dense forest and agriculture/sparse forest as well as bare soil and built-up area. Figure 5. (A) OBIA Classified Map, (B) Pixel-based SVM Classified Map The pixel-based SVM result produced the highest value of overall accuracy and OBIA classified image obtained lower accuracy. However, the OBIA classification results showed better accuracy and realistic representations when compared to the topographic map of KWR collected from Department of Survey and Mapping Malaysia (JUPEM) and also Google Earth Map. Because of these comparison and validation (figure 5). This study adopted the OBIA image classification for all three temporal Landsat satellite data classifications. 2.5. Markov model The stochastic model that the model output is depending on the probabilities of a transition of current change scenarios of Pi – j, between states i and j. The multiple land covers land uses categories in a landscape with transition probability Pi j , would be the land-cover type (pixels) i probability in time t0 A B https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 | 25 changes to land-cover type j probability in time t1 is called Markov chain model (Bozkaya et al., 2015). The Markov transitions probabilities, expressed as; (1) The derived transition probabilities is from a sample of transitions that occurred between two certain time intervals data and these probabilities was shown through the matrix P of following transition. The equation below proportion probability of land cover of the second date, and calculated using the equation. (2) Where vi j x Pi j is the proportion of land cover of the later date, Pi j is the matrix of the probability of landcover transition, Vi is the proportion of land cover of the current date (vector), i is the type of land cover of the first date, j is the type of land cover of the second date, P11 is the probability that land cover 1 at the first date will change into land cover 1 by the second date, P12 is the probability that land cover 1 at the first date will change into land cover 2 by the second date and so on, and m is the number of land-cover types in the study area (Bozkaya et al., 2015). 3. RESULTS AND DISCUSSION 3.1. Accuracy Assessment of Classified Images The Nearest Neighbor OBIA image classification of the temporal images, generated land cover maps and accuracy assessment report of confusion matrix were performed on (1989, 2004 and 2016) classified images indicated a satisfactory overall accuracy and a kappa statistics as work of (Kindu et al., 2013; Son et al., 2015; Yu, et al, 2016) in (Tables 4, 5 and 6). An overall Kappa Statistics of 0.88, 0.82 and 0.87 was achieved for 1989, 2004 and 2016 LULC with Classification accuracy of 91%, 86% and 90% respectively. In Table 4 built-up has a lower producer accuracy, followed by dense forest with 71.43% user accuracy. Bare soil has producer’s accuracy of 73.33% while dense forest has the highest producer’s accuracy of 100% showed in Table 5. All the remaining LULC classes were having their accuracies above 60%. The user’s accuracies of all the LULC types were above 60% with water and bare soil having the highest accuracy of 100%. Figure 6. Classified LULC Maps (1989, 2004 and 2016) 1989 2004 2016 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 26 | Table 4. Accuracy Results for the Landsat-5 1989 Image Derived from OBIA Classification Methods Note: Descriptions of LULC classes; U= Unclassified, W= Water, BS= Bare soil, DS=Dense forest, BU=Built- up, AFS= Agriculture/Sparse forest. Table 5. Accuracy Results for the Landsat-8 Image 2004 Derived from OBIA Classification Methods classified Reference Producer's accuracy % User's accuracy % U W F BU BU ASF Total Unclassified 0 0 0 0 0 0 0 0.00 0.00 Water 0 23 0 0 0 0 23 76.67 100.00 Dense forest 0 2 40 1 0 2 45 100.00 88.89 Built-up 0 1 0 22 1 1 25 73.33 88.00 Bare soil 0 2 0 7 27 3 39 90.00 69.23 Agriculture/ Sparse forest 0 2 0 0 2 34 38 85.00 89.47 Total 0 30 40 30 30 40 170 Overall Accuracy 86% Kappa coefficient 0.82 Note: Descriptions of LULC classes; U= Unclassified, W= Water, BS= Bare soil, DS=Dense forest, BU=Built- up, AFS= Agriculture/Sparse forest. Table 6. Accuracy Results for the Landsat-8 Image 2016 Derived from OBIA Classification Methods Classified Reference Producer's Accuracy % User's Accuracy % U W BS DF BU ASF Total Unclassified 0 0 1 0 1 0 2 0.00 0.00 Water 0 27 0 0 0 0 27 90.00 100.00 Bare land 0 0 23 0 1 0 24 76.66 95.83 Dense forest 0 1 1 39 0 3 44 97.50 88.64 Built-up 0 1 4 0 27 0 32 90.00 81.82 Agriculture/ Sparse forest 0 1 1 1 1 37 41 92.50 90.24 Total 0 30 30 40 30 40 170 Overall Accuracy 90% Kappa coefficient 0.87 Note: Descriptions of LULC classes; U= Unclassified, W= Water, BS= Bare soil, DS=Dense forest, BU=Built- up, AFS= Agriculture/Sparse forest. Classified Reference Producer's Accuracy % User's Accuracy % U W BS DF BU ASF Total Unclassified 0 0 0 0 0 1 1 0.00 0.00 Water 0 27 0 0 0 0 27 90.00 100.00 Bare soil 0 0 40 0 0 0 40 97.56 100.00 Dense forest 0 3 0 25 7 0 35 86.21 71.43 Built-up 0 0 0 2 23 0 25 76.67 92.00 Agriculture/ Sparse forest 0 0 1 2 0 39 42 97.50 92.85 Total 0 30 41 29 30 40 170 Overall Accuracy 91% Kappa coefficient 0.88 https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 | 27 3.2. Change Detection Trend analysis of the KWR (within and around) reveals changes in area of the five LULC of 27 year period of the study in Table 7 and Figure 7 and 8. The land cover changes being taken by both natural activities and community activities surrounding KWR between the study periods were measured by using modifications from the late classified map to early classified map as normally applied by all researches including (Islam et al., 2018). Table 7 below showed the results changes in hectares and percentages that has been revealed in the past three distinct years of the study generated and computed in terms of maps and tables through LCM. Areas covered with dense forest showed an immense changes of -7.05% between 2004 -2016 sooner than 1989-2004 with interval of 15 years period with an only 5.63%. The significant proportions of changes in agriculture/sparse forest from 4137.93 ha in between 1988-2004 to 18772.92 ha between 2004-2016 with relation to earlier trend changes of dense forest to bare soil between the period of 1989-2004. Furthermore, there was an increase of built up areas and water. The land cover change detection maps in figure 3 and 4, indicated the changes amongst five classes recognized in this study. 1989 t0 2004 changes result showed that the highest changes was between dense forest to bare soil/open land around the KWR area despite all environmental and biodiversity management measures. However, between the periods of 2004 to 2016 the results indicated that the bare soil area and dense forest were converted more to agriculture/sparse forest with little encroachment around the PA boundary (Ahmad et al., 2012; de Oliveira et al., 2017; Zhang et al., 2017). Table 7. LULC Changes Between the Study Periods in Hectares and Percentages Land cover classes 1989-2004 2004-2016 1989-2004 2004-2016 Area (ha) Area (ha) Area (%) Area (%) Water -983.43 2140.92 -0.31 0.67 Forest -17996.4 -22629.78 -5.63 -7.05 Built-up -2334.78 14436.99 -0.73 4.50 Bare soil 17289.18 -12718.80 5.39 -3.96 Agriculture 4137.93 18772.92 1.28 5.85 Figure 7. LULC Changes Between 2004-2016 of KWR (within and around) https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 28 | Figure 8. LULC Changes Between 2004-2016 of KWR (within and around) 3.3. Gain and Losses between 1989 to 2004 and 2004 to 2016 It is clearly indicated from Figure 9 that there were significant negative changes and transitions within and around the boundary of KWR for various LULC classes during the period from 1989 to 2004 and 2004 to 2016. The main gain and losses occurred basically between the dense forest to bare soil from 1989 to 2004 while between 2004 to 2016 showed that the losses in bare soil are positively changed to gain in agriculture/sparse forest. This make the analysis valid since changes started from removal of forest at the beginning then planting agriculture products as started by previous study that agricultural activities and illegal logging is the major land use activities taking place around the Krau Wildlife Reserve. Figure 8 illustrated the increase and decline that occurred between LULC adopted in this study in hectares from 1989 to 2004 and 2004 to 2016. The green bars represent the gain per class measured in hectares, and the left side brown bars describe the loss (decline) of each class in the same unit. In terms of net changes between the periods of study Figure 9 also indicated. Between 1989 and 2004 there is increase in the amount of bare soil (20,803 ha), more of the dense forest were proportionally lost to about 39,706 ha. The agriculture/sparse forest had the maximum extent of gains (49,000 ha) between 2004 and 2016 while about 38,324 ha was lost for dense forest within the same period (Figure 8). Figure 10 and figure 11 indicate the net changes contributions of each land use categories between 1989 to 2004 and 2004 to 2016. 0 20000 40000-20000-40000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C las se s Gains and losses between 1989 and 2004Gains and losses between 1989 and 2004 -1355358 -39706 21678 -83796036 -3524 20803 -33841 37930 0 20000 40000-20000-40000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C las se s Gains and losses between 2004 and 2016Gains and losses between 2004 and 2016 -4042538 -38324 15752 -4593 19007 -20560 7856 -30281 49008 Figure 9. Land Use/Land Cover Gain and Losses in (Ha) from 1989 to 2004 and 2004 to 2016 https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 | 29 From figure 8 it is clear that there are significant changes and transitions among various LULC categories during the period from 1989 to 2016. The main changes and transitions are mostly occurred among dense forest, bare soil and agriculture/sparse forest (Reddy et al., 2017). 0 6000 12000 18000-6000-12000-18000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Net Change between 1989 and 2004Net Change between 1989 and 2004 -997 -18029 -2343 17280 4090 0 10000 20000-10000-20000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Net Change between 2004 and 2016Net Change between 2004 and 2016 2135 -22572 14414 -12704 18727 Figure 10. Land Use/Land Cover Net Changes in (Ha) at KWR from 1989 to 2004 and 2004 to 2016 The contributions of other categories to their net change is presented in figure 11 and figure 12 below. It has been clearly shown that dense forest contribute about 6,925 ha to bare soil between 1989 to 2004 and dense forest for both period of study explained the majority of the total increase in agricultural/sparse forest areas (11,868 ha and 16,269 ha). Nevertheless, other LULC contribute to the changes throughout the study period as shown in figure 11 and figure 12 respectively. 0-40-80-120-160-200-240-280-320-360-400 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Contributions to Net Change in WaterContributions to Net Change in Water 0 -381 -80 -155 -381 0-2000-4000-6000-8000-10000-12000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Contributions to Net Change in Dense forestContributions to Net Change in Dense forest 381 0 384 -6925 -11868 0-200-400-600-800-1000-1200 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Contributions to Net Change in Built-upContributions to Net Change in Built-up 80 -384 0 -1286 -753 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 30 | 0 1000 2000 3000 4000 5000 6000 7000 8000 9000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Contributions to Net Change in Bare soilContributions to Net Change in Bare soil 155 6925 1286 0 8913 0 3000 6000 9000 12000-3000-6000-9000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Contributions to Net Change in Agriculture/Sparse forestContributions to Net Change in Agriculture/Sparse forest 381 11868 753 -8913 0 Figure 11. Contribution to Net Changes for All LULC Classes from 1989 to 2004 0 200 400 600 800 1000 1200 1400 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C las se s Contributions to Net Change in WaterContributions to Net Change in Water 0 1446 -36 94 631 0-3000-6000-9000-12000-15000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Contributions to Net Change in Dense forestContributions to Net Change in Dense forest -1446 0 -4380 -477 -16269 0 1000 2000 3000 4000 5000 6000 7000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Contributions to Net Change in Built-upContributions to Net Change in Built-up 36 4380 0 3139 6859 0-2000-4000-6000-8000-10000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C la ss es Contributions to Net Change in Bare soilContributions to Net Change in Bare soil -94 477 -3139 0 -9948 0 4000 8000 12000 16000-4000-8000 Water Dense forest Built-up Bare soil Agriculture/Sparse forest LU LC C las se s Contributions to Net Change in Agriculture/Sparse forestContributions to Net Change in Agriculture/Sparse forest -631 16269 -6859 9948 0 Figure 12. Contribution to Net Changes in Agriculture/Sparse Forest from 2004 to 2016 https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 | 31 3.4. Prediction The simulation model of LCM was used to simulate LULC modifications pattern (Mishra et al., 2014) within and around KWR for 2028 and 2040 based on Markov transition probabilities showed in appendix A and B in appendix respectively. The probability of diagonal cells (Table 8 and Table 9) represent an area which remain under the same class (Areendran et al., 2017). However, the predicted LULC (Reddy et al., 2017) showed the forest area keep reduced in the next projected time period (2028 and 2040) and has the highest probability of changes from dense forest converted to agriculture/sparse forest and built-up in Figure 13. Table 8. Markov Probability of changes for 2028 Figure 13. Predicted Land Cover Maps of KWR Area over 22 Years with Interval of 12 Years (2028 and 2040) Markov Probability of changes for 2028 LULC Class Water Forest Built-up Bare soil Agriculture Water 0.52 0.14 0.08 0.00 0.25 Forest 0.01 0.78 0.03 0.01 0.17 Built-up 0.00 0.07 0.42 0.05 0.45 Bare soil 0.00 0.07 0.17 0.04 0.72 Agriculture 0.01 0.12 0.09 0.05 0.73 https://doi.org/10.14710/geoplanning.5.1.17-34 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 32 | The simulated projected LULC maps in Figure 13 can be utilized as tools of decision making toward the protection, conservation and implementation of law enforcement in KWR, because the results indicate that if proper actions not taking the rate of degradation of dense forest will keep increasing and can lead to decline of both plants and animals habitat species (Conservation and Environmental Management Division, 2006) Malaysia. Table 9. Markov Probability of Changes for 2040 Markov Probability of changes for 2040 LULC Class Water Forest Built-up Bare soil Agriculture Water 0.2742 0.2195 0.1056 0.0202 0.3805 Forest 0.0129 0.6395 0.0516 0.0189 0.2771 Built-up 0.0081 0.1419 0.2319 0.0464 0.5717 Bare soil 0.0092 0.1575 0.1458 0.0455 0.642 Agriculture 0.011 0.1931 0.1191 0.0434 0.6334 4. CONCLUSION Mapping and predicting the LULC changes in a PA is very important for monitoring the activities within and around it. This can minimize the negative impact and help to plan for future managements to safeguard the KWR. The mapping analysis of KWR using multi-temporal satellites data showed and predicted the gradual loss in natural forest area within and around from 1989-2016 and verified by field visits and interviews with KWR officials. Oil palm and rubber plantations are one of the main factors leading to the encroachment around the KWR boundary. Moreover, the size of low dense forest/agriculture land, of the analyzed years (1989, 2004 and 2016) was increased and was found in predicted results of 2028 and 2040. The changes within and around KWR showed a massive degradation and if left unattended through current situation based on projected land use and land cover, it will be detrimental for biological conservation of wildlife in the PA. The changes and encroachment around KWR boundary have a link with the dynamics of political and social issues of local communities surrounding the wildlife reserve. This study also suggested the implementation of buffer zones which may be one of the key solutions for a better conservation of PAs to protect from the negative effects of illegal activities within and around the KWR. It is also recommended that it increase the number of forest rangers in the KWR to monitor the encroachments by local and indigenous communities. Future work can improve the work further by utilizing more data from various systems such as radar, lidar and very high resolution data. Furthermore the results for the PA analysis in this study can be imported to a GIS for further analysis and model development. 5. ACKNOWLEDGMENTS The authors would like to thank UPM for the facilities and funding for research and travel in completing this task. In addition, we thank the Department of Wildlife and National Parks, Peninsular Malaysia (Perhilitan) and the Institute of Biodiversity (IBD) for providing useful data and information. 6. REFERENCES Ahmad, C. B., Abdullah, J., & Jaafar, J. (2012). Community Activities Around Protected Areas and the Impacts on the Environment at Krau Wildlife Reserve, Malaysia. Procedia-Social and Behavioral Sciences, 68, 383–394. [Crossref] Areendran, G., Raj, K., Mazumdar, S., & Sharma, A. (2017). Land Use and Land Cover Change Analysis for Kosi River Wildlife Corridor in Terai Arc Landscape of Northern India: Implications for Future Management. Tropical Ecology, 58(1). https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.1016/j.sbspro.2012.12.235 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 | 33 Balaji, S. A., Geetha, P., & Soman, K. P. (2016). Change Detection of Forest Vegetation using Remote Sensing and GIS Techniques in Kalakkad Mundanthurai Tiger Reserve - (A Case Study). Indian Journal of Science and Technology, 9(30), 1–6. [Crossref] Bozkaya, A. G., Balcik, F. B., Goksel, C., & Esbah, H. (2015). Forecasting land-cover growth using remotely sensed data: a case study of the Igneada protection area in Turkey. Environmental Monitoring and Assessment, 187(3). [Crossref] Bush, A., Sollmann, R., Wilting, A., Bohmann, K., Cole, B., Balzter, H., … Yu, D. W. (2017). Connecting Earth observation to high-throughput biodiversity data. Nature Ecology & Evolution, 1(7), 176. [Crossref] Conservation and Environmental Management Division. (2006). Biodiversity in Malaysia. DANCED, & Jabatan Perlindungan Hidupan Liar dan Taman Negara. (2001). Krau Wildlife Reserve Management Plan. Perhilitan. de Oliveira, S. N., de Carvalho Júnior, O. A., Gomes, R. A. T., Guimarães, R. F., & McManus, C. M. (2017). Deforestation analysis in protected areas and scenario simulation for structural corridors in the agricultural frontier of Western Bahia, Brazil. Land Use Policy, 61, 40–52. [Crossref] Desclée, B., Bogaert, P., & Defourny, P. (2006). Forest change detection by statistical object-based method. Remote Sensing of Environment, 102(1–2), 1–11. [Crossref] Despot Belmonte, K., & Bieberstein, K. (2016). Protected Planet Report 2016. How Protected Areas Contribute to Achieving Global Targets for Biodiversity. Dudley, N., & Stolton, S. (2008). Defining protected areas: An international conference in Almeria, Spain Mayo 2007. IUCN Protected Areas Categories Summit. Dutta, K., Reddy, C. S., Sharma, S., & Jha, C. S. (2016). Quantification and Monitoring of Forest Cover Changes in Agasthyamalai Biosphere Reserve, Western Ghats, India (1920-2012). Current Science, 110(4), 508. [Crossref] Hackman, K. O., Gong, P., & Wang, J. (2017). New land-cover maps of Ghana for 2015 using Landsat 8 and three popular classifiers for biodiversity assessment. International Journal of Remote Sensing, 38(14), 4008–4021. [Crossref] Hruby, F., Melamed, S., Ressl, R., Stanley, D., Balancing, C., Imagery, S., & Data, B. (2016). Mosaicking Mexico - The Big Picture of Big Data, XLI (July), 407–412. [Crossref] Islam, K., Jashimuddin, M., Nath, B., & Nath, T. K. (2018). Land use classification and change detection by using multi-temporal remotely sensed imagery: The case of Chunati wildlife sanctuary, Bangladesh. The Egyptian Journal of Remote Sensing and Space Science, 21(1), 37–47. [Crossref] Kindu, M., Schneider, T., Teketay, D., & Knoke, T. (2013). Land Use/Land Cover Change Analysis Using Object-Based Classification Approach in Munessa-Shashemene Landscape of the Ethiopian Highlands. Remote Sensing, 5(5), 2411–2435. [Crossref] Kumar, K. S., Valasala, N. V. A. S. S., Subrahmanyam, J. V. V, Mallampati, M., Shaik, K., & Ekkirala, P. (2015). Prediction of Future Land Use Land Cover Changes of Vijayawada City Using Remote Sensing and Gis. International Journal of Innovative Research in Advanced Engineering (IJIRAE), 2(3), 91–97. Lin, T. S. (2016). Sumber Air Orang Asli Tercemar Akibat Pembalakan. Pahang: Astro Awani. Retrieved May Tuesday, 2017, from http://www.astroawani.com/berita-malaysia/sumber-air-orang-asli- tercemar- akibat-pembalakan-111904. (2017). Mishra, V., Rai, P., & Mohan, K. (2014). Prediction of land use changes based on land change modeler (LCM) using remote sensing: A case study of Muzaffarpur (Bihar), India. Journal of the Geographical Institute Jovan Cvijic, SASA, 64(1), 111–127. [Crossref] Munthali, K. G., & Murayama, Y. (2011). Land use/cover change detection and analysis for Dzalanyama forest reserve, Lilongwe, Malawi. Procedia - Social and Behavioral Sciences, 21, 203–211. [Crossref] Nagendra, H., Lucas, R., Honrado, J. P., Jongman, R. H. G., Tarantino, C., Adamo, M., & Mairota, P. (2013). Remote sensing for conservation monitoring: Assessing protected areas, habitat extent, habitat condition, species diversity, and threats. Ecological Indicators, 33, 45–59. [Crossref] Norawi, M. F. (2017). Alam dakwa 2,000 balak dicuri di Kuala Krau. Sinar Online, pp. 4–7. Ranjan, A. K., and Akash Anand, S, V., & Singh, R. K. (2016). LU/LC Change Detection and Forest Degradation Analysis in Dalma Wildlife Sanctuary Using 3S Technology: A Case Study in Jamshedpur-India. AIMS Geosciences, 2(4), 273–285. [Crossref] https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.17485/ijst/2016/v9i30/99022 https://doi.org/10.1007/s10661-015-4322-z https://doi.org/10.1038/s41559-017-0176 https://doi.org/10.1016/j.landusepol.2016.10.046 https://doi.org/10.1016/j.rse.2006.01.013 https://doi.org/10.18520/cs/v110/i4/508-520 https://doi.org/10.1080/01431161.2017.1312619 https://doi.org/10.5194/isprsarchives-XLI-B2-407-2016 https://doi.org/10.1016/j.ejrs.2016.12.005 https://doi.org/10.3390/rs5052411 https://doi.org/10.2298/ijgi1401111m https://doi.org/10.1016/j.sbspro.2011.07.035 https://doi.org/10.1016/j.ecolind.2012.09.014 https://doi.org/10.3934/geosci.2016.4.273 Gambo et al. / Geoplanning: Journal of Geomatics and Planning, Vol 5, No 1, 2018, 17-34 doi: 10.14710/geoplanning.5.1.17-34 34 | Reddy, C. S., Singh, S., Dadhwal, V. K., Jha, C. S., Rao, N. R., & Diwakar, P. G. (2017). Predictive modelling of the spatial pattern of past and future forest cover changes in India. Journal of Earth System Science, 126(1). [Crossref] Reveshty, M. A. (2011). The Assessment and Predicting of Land Use Changes to Urban Area Using Multi- Temporal Satellite Imagery and GIS: A Case Study on Zanjan, IRAN (1984-2011). Journal of Geographic Information System, 03(04), 298–305. [Crossref] Son, N.-T., Chen, C.-F., Chang, N.-B., Chen, C.-R., Chang, L.-Y., & Thanh, B.-X. (2015). Mangrove Mapping and Change Detection in Ca Mau Peninsula, Vietnam, Using Landsat Data and Object-Based Image Analysis. IEEE Journal of Selected Topics in Applied Earth Observations and Remote Sensing, 8(2), 503–510. [Crossref] Waiyasusri, K., Yumuang, S., & Chotpantarat, S. (2016). Monitoring and predicting land use changes in the Huai Thap Salao Watershed area, Uthaithani Province, Thailand, using the CLUE-s model. Environmental Earth Sciences, 75(6). [Crossref] Willis, K. S. (2015). Remote sensing change detection for ecological monitoring in United States protected areas. Biological Conservation, 182, 233–242. [Crossref] Yu, W., Zhou, W., Qian, Y., & Yan, J. (2016). A new approach for land cover classification and change analysis: Integrating backdating and an object-based method. Remote Sensing of Environment, 177, 37–47. [Crossref] Zhang, C., Smith, M., Lv, J., & Fang, C. (2017). Applying time series Landsat data for vegetation change analysis in the Florida Everglades Water Conservation Area 2A during 1996{\textendash}2016. International Journal of Applied Earth Observation and Geoinformation, 57, 214–223. [Crossref] https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.14710/geoplanning.5.1.17-34 https://doi.org/10.1007/s12040-016-0786-7 https://doi.org/10.4236/jgis.2011.34026 https://doi.org/10.1109/jstars.2014.2360691 https://doi.org/10.1007/s12665-016-5322-1 https://doi.org/10.1016/j.biocon.2014.12.006 https://doi.org/10.1016/j.rse.2016.02.030 https://doi.org/10.1016/j.jag.2017.01.007