{'fraction_non_alphanumeric': 0.03949441420796333, 'fraction_numerical': 0.04615439702091091, 'mean_word_length': 4.763309946347503, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Assessment and Prediction of Land Use/Land Cover Change in the National Capital of Burundi Using Multi-temporary Landsat Data and Cellular Automata-Markov Chain Model

Gitega District has experienced significant land use and land cover changes due to human activity. This has increased land degradation and environmental issues. However, there is no data on LULC change to guide land-use planning. This study assessed the rate and magnitude of LULC change over the last 35 years and also simulated future scenarios using Geoinformatics. In the first step, five LULC classes were extracted from satellite images from 1984, 2002, and 2019 using the supervised classification method. Overall accuracy and Kappa statistics of more than 85% and 82% respectively were achieved with 30 reference samples. Change analysis highlighted by Land Change Modeler (1984-2019) indicated a significant increase in Agriculture of 94 km 2 , a slight increase in Shrub Land and Built-up Area of 5.5 km 2 and 2 km 2 , respectively; and a steep decrease in Trees Cover and Grass Land of 62.5 km 2 and 39 km 2 , respectively. Markov Chain and CA-Markov models were further calibrated to simulate LULC changes in 2038 and 2057 using the 2019 base map. Evaluation and analysis of 2019-2057 simulation results showed a moderate agreement of 75% for Kappa and the same trends of LULC change: Trees Cover, Grass Land, and Shrub Land will decrease by 11.5 km 2 , 13 km 2 , 11.5 km 2 respectively, whereas Agriculture and Built-up Area will increase by 30 km 2 and 6 km 2 respectively in 2057. These study outcomes can support decision-making towards restoration measures of land degradation and long-term environmental conservation in the region. Change Modeler/


INTRODUCTION
Land cover change denotes alteration in physical land types such as forests, vegetation, and so on, whereas land-use change refers to changes in a specific area of land used or managed by humans for their provisions (Liu and Shi, 2019;Patel et al., 2019). Land use and land cover (LULC) change is a result of the human activities and subsequent development taking place across the world (Bai et al., 2008). Land use pattern is triggered by population dynamics which is directly responsible for the changes in land cover (IPBES, 2018;Serneels and Lambin, 2001;Verma and Raghubanshi, 2019).
The impact of LULC change on the Earth's environment has increased as a result of rapid growth in human population which leads in turn to deforestation, intensive agricultural activities, urbanization, and other unplanned human settlement (Boissière et al., 2009;Lambin et al., 2000;Meyer et al., 1994). According to a global report published by FAO, 24% of the world's land is degraded due to anthropogenic influences, primarily in Africa, Asia, Australia, and North America (Bai et al., 2008). In recent years, assessing and modelling LULC change has become an essential aspect since the LULC change has influences on socioeconomic, and environmental structure (Halmy et al., 2015;Twisa and Buchroithner, 2019).
The study of assessing and modelling LULC change is an expensive and time-consuming process. However, through Remote Sensing (RS) and Geographic Information System (GIS) techniques, it is possible to examine LULC change on large spatial and temporal scales (Dewan and Yamaguchi, 2009;Nijimbere et al., 2019). R.S advancements allow us to obtain synoptic information about LULC at a specific time and location (Anderson, 1976;Patil et al., 2012). Multi-temporary satellite images have been widely used for LULC change detection (Kumar et al., 2010;Vila and Barbosa, 2010). The GIS approach, namely Geoinformatics, through its ability to analyze, update and backdate LULC information with object-oriented or pixel-based image classification, has been shown to produce accurate results for LULC monitoring (Mishra et al., 2014;Shen, 2019). Such classification accuracy is performed by both the spatial and spectral resolution of satellite image data (Poursanidis et al., 2015).
Moreover, the future LULC change prediction requires the integration of past trends and current landscape, and plausible assumptions. Therefore, the literature describes various forecasting models based on research purposes, i.e., landscape dynamics, ecological modelling, disaster and deforestation assessment, agriculture projects and urban sprawl planning, etc., (Baker, 1989). Verburg (2006) argues that no single land-use model can be superior to another. But, some other researchers have expressed their confidence in the significant results obtained when integration of Cellular Automata and Markov Chain models is used in the LULC modelling process (Mondal et al., 2016;Nadoushan et al., 2015). Several national land cover mapping programs have been developed to continuously monitor LULC changes and produce national and regional land cover maps (Gómez et al., 2016). This has contributed to a better understanding of natural-human interactions and advanced monitoring and modelling of landscape dynamics (Ghosh et al., 2017;Turner et al., 2007).
Burundi is a small landlocked country with an area of 27,834 km 2 , and the most densely populated country in Africa with 480 people/km 2 (UNdata, 2020). Due to an increasing number of householders who rely on agriculture to meet their needs (Moore, 2007), the country is vulnerable to major landscape dynamics and environmental problems (Kamungi et al., 2005;Nzabakenga et al., 2013). It is currently agreed that access to arable land is a top priority for any household. With 90% of the 12 million as workforces in farming land, land resources are being much more depleted and scarcer than in the past (Nzabakenga et al., 2013).
Gitega District, one of the most densely populated regions in Burundi with 476 people/km 2 (Niyuhire, 2018), is prone to land degradation and environmental vulnerabilities as a result of long-term major landscape dynamics. This region has attracted a large number of immigrants due to its geographic location in the center of the country, which is an important aspect of being a National Capital (Guichaoua, 1982). Significant LULC change started occurring in 2007, when the government of Burundi began the process of restoring that district's former status as National Capital (Patel et al., 2019). Thus, the rapid population increase in that area results in agricultural expansion along with built-up area and reduction of forest and grassland in the Gitega's landscape (Guichaoua, 1982;Niyuhire, 2018). As per El-Hassanin et al. (1993) and Marathianou et al. (2000), LULC changes reduce normalized difference vegetation index (NDVI) of land. This provokes other extreme impacts on the environment such as climate change, extreme radiative forcing, pollution and quality reduction of natural ecosystems, changes in hydrological regimes, runoff, soil loss and depletion of soil fertility (IPBES, 2018;Verma and Raghubanshi, 2019). Certain recent studies have only addressed agricultural issues confronting farmers at the national and district levels. Kamungi et al. (2005), land access, Nzabakenga et al. (2013), agriculture and farmers' livelihoods, Nijimbere et al. (2019), soil erosion, Niyuhire (2018), soil depletion and integrated soil fertility. Though, they did not study land degradation with LULC change patterns.
However, Burundi has neither national nor regional mapping platforms for generating basic information about the LULC change and environmental dynamics that could aid in land use planning (Burrell et al., 2018;Ndzabandzaba, 2015). Therefore, the main objective of this study was to assess LULC changes that occurred in the Gitega District over the last 35 year by calculating the rate and magnitude of change for the LULC category. This study also aimed to simulate the future scenario of LULC changes under an assumption of the continuation of past and current trends using the combination of RS, GIS, and Cellular Automata and Markov Chain approach.

METHODOLOGY
The Figure 1 shows the overall methodology used for this study in three main parts: data source, data processing, and results analysis.

Study area
Gitega District is located in the center of Burundi at a specific geographic grid reference of "03° 25' 35" South Latitude and 29° 50' 37" East Longitude ( Figure 2). It has an area of 315 km 2 , 150,151 people, and a density of 476 inhabitants/km 2 (Niyuhire, 2018). Its topography dominated by plateau dispersed by hills, valleys, and moderate plains rises between 1,600-2,000 m. Climate is subtropical highland and tropical savanna climate depicting summer and winter (Vassolo et al., 2019). The average annual temperature is 20°C with lower and higher temperatures of 10°C and 26°C respectively, while the average annual rainfall is 1,130 mm with lower and higher rainfall in July and April respectively (Vassolo et al., 2019). Due to those mild climate and plentiful rainfall, the area is ideal for forests, vegetation, and plantations. However, high population pressure on land and poor farming practices are unsuitable for environmental preservation (El-Hassanin et al., 1993;Kamungi et al., 2005). With two harvests per year, cultivation of subsistence foods such as bananas, corn, manioc, sweet potatoes, Irish potatoes, beans, peas, wheat, peanuts, vegetables, etc., is intensified even on steep slopes and hilltops (Kamungi et al., 2005;Nijimbere et al., 2019). This study area, which includes both urban and rural land use characteristics, is likely to experience the most change in its landscape. Since the project to restore its former status as the national capital was launched in 2007, several buildings and road patterns have been constructed to meet the county's second city functions.

Data acquisition and processing
Three Landsat Images with 30 m resolution were downloaded from the United States Geological Survey (http://earthexpolorer.usgs.gov/) and important information is shown in Table 1. All these images were taken on satellite track path/row:173/062 and projected in UTM with WGS-84 datum 36 N.

Data Source
Landsat imagery (1984, 2002, and  Prior to classifying the above satellite images from different dates, an image-based atmospheric correction was performed by subtracting a common minimum digital number (DN) value from the image using Idrisi software. Aside from noise removal, another image preprocessing based on image optical properties, such as sun elevation, angle, and spectral was applied using image metadata file. False-color composite and contrast stretching techniques were further applied to improve brightness and tonal differentiation between various features (Song et al., 2001).

Image classification and mapping
Several classification methods were developed to extract useful information from imagery. Two of the most popular land cover classification machines were also used in this research (Pulighe et al., 2016;Tilahun and Islam, 2015): Supervised classification, which involves manually defining land cover classes and forcing signatures based on field identification and local knowledge of land cover types (Eastman, 2009;Patil et al., 2012). In this process, the common classification scheme in the LULC system (Anderson, 1976) and field survey data were used to form five LULC categories, for what training sits" samples were created.
Maximum Likelihood classifier was then used, whose function is to assign the pixel to the class with the highest probability. The training site samples generated during the supervised process were fully distributed on the entire image in terms of spectral information for each LULC class. The signature file containing spectral information was thereby generated (Eastman, 2009;Patil et al., 2012). Finally, ArcGIS software was used to map the five LULC categories: Agriculture, Built Area, Grass Land, Shrub Land, and Tree Cover. Figure 3 shows LULC maps for 1984, 2002, and 2019, respectively, and results of image classification and the area distribution is given in Table 2.

Classification evaluation process
The high-resolution contemporary satellite imagery available on Google Earth pro was used to collect the ground truth data for 1984, 2002, and 2019 maps. Thirty pixels were selected using a stratified sample random method with ArcGIS to examine image classification accuracy (Pulighe et al., 2016;Tilahun and Islam, 2015). The confusion matrix Boundary of Gitega Class 1. Agriculture 2. Built-up area 3. Grass land 4. Shrub land 5. Trees cover method reporting overall accuracy (OA) user's accuracy (UA), producer's accuracy (PA) and Kappa statistics (K) were generated (Dewan and Yamaguchi, 2009;Rwanga and Ndambuki, 2017). OA measures the percentage of total classified pixels that were truly labeled into the specific land cover and is computed by dividing the total correctly classified pixels (TCS or the sum of the diagonals) by the number of reference pixels (TS) in the error matrix using equations (1), (2), and (3): Where; TCSij is the total number of the correctly classified pixels in row i and column j, TS is the total reference sample, TScj is the total number of pixels in column j and TSri is the total number of pixels in the row i. The quantitative measure of the level of agreement was done with utilization of kappa statistic (K) assumed that a K of 1 indicates ideal agreement, whereas a kappa of 0 indicates agreement equivalent to chance to truly classify the pixels as computed in the equation (4): Where; TCSij is the total number of the correctly classified pixels in row i and column j, TS is the total reference sample, TScj is the total number of pixels in column j and TSri is the total number of pixels in the row i. Moderate agreement of more than 85% and 82% for overall accuracy and Kappa statistics respectively was obtained as presented in Table 3.

Cellular Automata and Markov Chain for simulation the future LULC Change
The combination of Cellular Automata and Markov Chain models have been the most popular integrated model used for modelling spatial and temporal changes (Borana and Yadav, 2017;Sang et al., 2011). Markov Chain is a dynamic process that calculates the probability of changes from a particular object (e.g., vegetation) to another object (e.g., agriculture) in the next state depending on the former state. It has been used equally in LULC studies due to its ability to quantify conversion rates between two categories (Ghosh et al., 2017). Cellular Automata (CA) is a discrete dynamic system consisting of a normal finite cell network that can change in various directions in neighboring (Singh, 2003;Torrens, 2000). CA is also widely used in simulating LULC change because it is fairly close to GIS raster data (Pijanowski et al., 2002;Surabuddin et al., 2013). Therefore, the combination of the two above models was also applied in this study to simulate the LULC in 2019, 2038, and 2057 using the following equations (5) and (6) (Mondal et al., 2016;Nadoushan et al., 2015): Where; S (t) is the system status at time t; S (t+1) is the system status at a time of t+1, and Pij is the transition probability obtained in equation (6) below: P = |P ij | = | P 1,1 P 1,2 … P 1,N P 2,1 P 2,2 … P 2,N ⋯ ⋯ ⋯ P N,1 P N,2 … P N,N | 0 ≤ P ij ≤ 1 Where; P=Transition probability; Pij=stands for probability of changing from particular state i to another state j; P_N=state probability of any time, and N=Land cover type (Eastman, 2009;Mishra et al., 2014).
In this study, 1984, 2019LULC maps were used to run Markov Chain. Firstly, 1984 LULC maps were used as earlier and later images to project the 2019 LULC map (Figure 7), along with transition probability matrix shown in Table 4. This method produced a suitability map matching the current rate and quantity of 2019 LULC which will be further used to simulate future scenarios (Sang et al., 2011). The 1984 LULC map was used in CA-Markov as base map and Cellular automata iterations were set to 35 as the time interval between 1984 and 2019 to predict the 2019 LUC map. This simulation outcome was compared with the existing LULC to evaluate the model performance (Figure 4) (Gupta and Sharma, 2020). The 2019 LULC simulation was achieved using 1984-2002 transition probability area matrix, transition suitability map, and a set window technique 5*5 kernels called standard contiguous filter for influencing each cell center. After successful validation of 2019 simulation result, predictions of 2038 and 2057 LULC maps were done (Figure 7). At this stage, we used the 2019 real LULC map as the simulation's input, along with the relevant transition probability area matrices, and transitional suitability maps generated from Markov Chain (Nadoushan et al., 2015;Singh, 2003). Based on the modelling theory, we kept a default set of 5*5 kernel size contiguity filter, while the cellular automata iteration was set to 19 to simulate the 2038 LULC map, and the same method was used to simulate 2057 LULC map.

Simulation accuracy and model evaluation
Simulation accuracy and model validation were determined by comparing the predicted LULC maps for 2019, 2038, and 2057 to the actual LULC map for 2019. The Kappa statistics is the most widespread method agreed to quantify the power and suitability of the simulation model (Maingi et al., 2002). This study utilized the submodule in Idrisi namely, GIS Analysis, to generate the Kappa index of Agreement shown in Table 5. The validation was done based on the assumption that the higher the kappa values, the better the model (Borana and Yadav, 2017). All kappa index values surpass the minimum acceptable standard and they range from 65% to 89%. This indicates a high degree of agreement between projected and actual LULC maps, thus a successful validation of the CA-Markov Chain model (Kundel and Polansky, 2003;van Vliet et al., 2011).

LULC Change Analysis
The LULC changes analysis enables an understanding of physical modification or loss of features in the natural landscape such as vegetation and forests clearing, agricultural land, waterbody as worthwhile information useful for planning land use and environmental conservation (Wang et al., 2021). We applied the comparative method to analyze five LULC classes from image classifications based on remote sensing techniques (Dewan and Yamaguchi, 2009). The land use and land cover change analysis were accomplished by using land cover data with an interval of 35 years, from 1984 to 2019 which depicts Based on Figure 5, it is revealed that the landscape of Gitega in 1984 was extensively covered by Grass Land and Trees Cover: 112 km 2 (36.97%) and 81 km 2 (26.73%), respectively. Shrub Land is the third most dominant class, accounting for 61 km 2 (20.13%), while Agriculture and Built-up Area account for 48 km 2 (15.84%) and 1 km 2 (0.33%), respectively. In the year 2002, the study area was largely covered with Grass Land and Agriculture:121 km 2 (39.84%) and 111.5 km 2 (36.92%), respectively. However, Trees Cover and Shrub Land decreased to 38 km 2 (12.24%) and 31 km 2 (10.21%), respectively, and Built-up Area had slightly increased. In 2019, the largest land cover was Agriculture with 142 km 2 , an increase of 10% from 2002 to 2019. Grass Land took a second areal position with 73 km 2 due to a decrease of nearly 5%. Shrub Land rose to the third position due to an increase by 11.74% with an area of 66,5 km 2 , while Trees Cover steeply decreased by 6.44% and occupied 18.5 km 2 . Built-up Area has a steady increase of 1.5 km 2 from 2002-2019. Over past 35 years, the dynamic in LULC change patterns is intensely high. Trees Covers and Grass Land were continuously reduced by 62.5 km 2 at an annual rate of 1.79 km 2 and 39 km 2 at the rate of 1.12 km 2 per year, respectively. In contrast, a net change indicates a large extension in Agriculture of 94 km 2 at the rate of 2.68 km 2 per year. A slight increase is also detected in both Shrub Land and Built Area, which is 5.5 km 2 with a rate of 0.15, and 2 km 2 with a rate of 0.8 km 2 per year, respectively. For getting a better understanding, LULC change detection was investigated using Land Change Modeler (LCM), and the model outputs shown in Figure 6 are consistent with certain results of land cover change detection obtained using the same model. For instance, Gupta and Sharma (2020) have acknowledged interesting results of LULC change monitoring by LCM model in LULC dynamics study (over past 50 years) conducted in India. These decreasing areas in Trees Cover and Grass Land (by about 101.5 km 2 in total with an annual rate of 2.9 km 2 ) have been converted mostly into Agriculture (94 km 2 ), and little to Shrub Land and Built-up Area (7.5 km 2 together) from 1984 to 2019. Similar trends of LULC change, especially deforestation and agricultural expansion, were globally addressed (FAO, 2016), and in many other developing countries (Alawamy et al., 2020;Berakhi, 2013;Henry et al., 2011;Islam et al., 2018). Interesting LULC change issues in terms of urban sprawl were recently dealt with by Patel et al. (2019) in India. Thereby, a general perception is that forest and vegetation reductions lead to dramatic degradation of productive ecosystems and loss of biodiversity (Wang et al., 2021). Primarily, the land degradation will radiate outwards within the global environment (IPCC, 2019;Marathianou et al., 2000;Niyogi et al., 2009). Deforestation alters hydrological processes and affects water conductivity like surface runoff (El-Hassanin et al., 1993) and, therefore, the occurrence of soil erosions and soil losses as observed in the study area (Henry et al., 2011;Nijimbere et al., 2019;Niyuhire, 2018;Nzabakenga et al., 2013). A review of the trends in population growth and density (476 people/km 2 ) in the Gitega District can reveal that the most LULC changes are likely the shortcomings of human encroachment on natural land resources and ecosystems (Kamungi et al., 2005). In this study, the results of the LULC change analysis for the last thirty-five years were used to make predictions for the next thirty-eight years.

Simulation result analysis
CA-Markov Chain simulation model requires earlier and later LULC maps to predict the future scenario (Sang et al., 2011). Future LULC simulations were projected for 2057 with a time interval of 38 years, from 2019 to 2057, which also depict two subtime intervals: 19 years (from 2019 to 2038) and 19 years (from 2019 to 2057). According to the Markov Chain concept, changes that occurred in the previous period can be used to predict changes in the next period (Ghosh et al., 2017;Wang et al., 2021). In this study, the simulation inputs were the real and projected LULC maps for 2019 and the probability matrix for the transition to land-use change in 2019-2038 and 2019-2057.
Based on the theory of transition probability, the Markov model performs the likelihood of a change from one land cover category to another within a dimensioned matrix of LULC classes (Nadoushan et al., 2015). For example, the analysis of the results from Markov Chain application in this study given in the cross-tabulation showed revealed that from 2019 to 2038, the probability of change for Agriculture to Agriculture is 68.24%, while the probability of future change of Agriculture to Grass Land is 20.76%. Trees Cover has a probability as low as 5.9% to remain as they are but has a probability of 37.06% to change to Agriculture and the same process for other LULC classes. In the second prediction scenario, from 2019 to 2057, Agriculture has the highest probability of 71.09% to remain the same in 2057, whereas Trees Cover indicates the most declining probability of 3% to remain same in 2057. Built-up Area, Grass Land, and Shrub Land have the probability of 49.01%, 30.84%, and 26.24%, respectively, to remain as they are in 2019. Table 4   Cellular Automata and Markov Chain models were spatially integrated to predict change location (Verburg, 2006). This simulation process for 2019-2057, which combines Cellular Automata and the Markov model, is similar to the simulation process for 1984-2019. The difference was in the time interval of input data used as the basis for simulation, which was due to the need for high-quality image data to successfully detect and model Gitega District LULC changes. The transition probability and changing area prediction from the Markov Chain process served as the input basis for Cellular Automata to determine the change location represented in pixel (Nadoushan et al., 2015;Singh, 2003). The simulation results are presented as predicted LULC maps. Table 5 displays the index values that indicate reasonable agreement for the model validation and the prediction results accuracy. The Figure 7 depicts maps of LULC prediction results in the Gitega District for 2038 and 2057 respectively simulated using estimated LULC map for 2019.
Based on Table 6, it is projected that Agriculture will continuously expand with the increase of 14 km 2 and 30 km 2 in 2038 and 2057, respectively. Built-up Area will increase by 6 km 2 in 2057. However, trees cover, grass land, and shrub land will be decreasing by 11.5 km 2 , 13 km 2 and 11.5 km 2 , respectively, in 2057. These predicted trends of LULC change showing such continuous agricultural expansion at expenses of forest and vegetation are also similar to those reported by other LULC modelling studies (Henry et al., 2011;Mondal et al., 2016;Nadoushan et al., 2015;Wang et al., 2021). Eventually, a better decision is needed to avoid these projected major LULC changes in Gitega District. Table 6 displays the results of change analysis of predicted LULC maps.  Shrub land 54 55 -12.5 -11.5 Trees cover 10 7 -8.5 -11.5

CONCLUSION
This paper assessed the rate and trends of change in LULC patterns quantified and predicted using Geoinformatics in Gitega District. The Change analysis was performed using 5 LULC classes obtained from images classification of years 1984, 2002, and 2019 and was highlighted with Land Change Modeler. Markov Chain and CA-Markov models were applied to project future trends of LULC changes in 2038 and 2057. Satisfactory accuracy with good agreements of more than 85% and 82%, respectively, for overall accuracy and Kappa statistics was achieved. RS, GIS, and Cellular Automata and Markov Chain models were therefore proven to be Boundary of Gitega 1. Agriculture 2. Built-up area 3. Grass land 4. Shrub land 5. Trees cover Class effective tools for LULC change monitoring and generating multi-temporary land cover dataset to assist decision-makers towards land use planning and environmental protection.
Overall findings reveal dramatic decreasing area in Trees Cover and Grass Land of 101.5 km 2 together, and which was likely converted mostly to Agriculture (which increased by 94 km 2 at a rate of 2.68 km 2 per year), and to Shrub Land and Built-up Area (7.5 km 2 together) over the past 35 years. These results are reflective of the current rate, trends, and magnitude of LULC changes and local policies in the study area. If Gitega attempts to avoid further irreversible land degradation and associated environmental problems, the government and policymakers should urgently apply agroecological approaches and reforestation.
Similar trends in LULC change patterns are revealed by predictive result analysis. Agriculture has 71.09% probability of remaining as Agriculture in 2057, with an expected area of 172 km 2 (56.76%). On the other hand, Trees Cover have the lowest probability of 3% to remain as Trees in 2057, with an area of 7 km 2 . These study outcomes are the potential to support decision-making to undertake restoration measures of land degradation and future sustainable land use management and environmental presservation .However, future studies should consider the correlation of Gross Domestic Products (GDP) and Population Growth with these LULC changes analysis results.