How Will Rwandan Land Use/Land Cover Change under High Population Pressure and Changing Climate?

In recent decades, population growth and economic development have greatly influenced the pattern of land use/land cover (LULC) in Rwanda. Nevertheless, LULC patterns and their underlying change mechanisms under future climate conditions are not well known. Therefore, it is particularly important to explore the direction of LULC transfer in the study area, identify the factors driving the transfer of different types of LULC and their changes, and simulate future LULC patterns under future climate conditions. Based on LULC analyses of Rwanda in 1990, 2000, 2010, and 2015, the LULC pattern of Rwanda in the next 30 years was simulated using an LULC transition matrix, random forest sampling, the Markov chain model, and the PLUS model. The results showed that LULC change in the study area primarily comprised a decrease in forest area and expansion of cropland area, accompanied by a small increase in grassland area and an annual increase in urban land area. Prior to 2000, the LULC in Rwanda was mainly converted from forest and grassland to cropland, with the ratio being 0.72:0.28. After 2010, the LULC was mainly converted from forest to grassland and cropland, with the ratio being 0.83:0.17. Changes in forests, grasslands, and cropland are driven by multiple factors, whereas changes in wetlands, water, urban land, and unused land are more likely to be driven by a single factor. The existing trend of LULC change will continue for the next 30 years, and the future LULC pattern will exhibit a trend in which cropland area will increase in the west and grassland area will decrease, whereas grassland area will increase in the east and cropland area will decrease.


Introduction
Over the years, sustainable development in Africa has been the focus of global attention. Although the population density of Sub-Saharan Africa is still relatively low (43 persons/km 2 ) [1], the population distribution in East Africa is extremely uneven because of the differences in topography and climate brought about by the Great Rift Valley. The population density in some areas exceeds 300-500 persons/km 2 [2,3]. The land resources in this area have been under long-term pressure from the ever-increasing population. The rapid population growth has given rise to numerous challenges, such as land degradation, loss of biodiversity, excessive deforestation, water shortages, and difficulty in ensuring food security. The emergence of these problems undoubtedly has a negative effect on sustainable development [4,5]. The land resources, which are critical for human survival, are the carriers of these problems. Anthropogenic activities and the natural environment can both induce and restrict changes in LULC [6,7]. The expansion of cities caused by rapid economic development and the expansion of arable land caused by the rapidly expanding population occur at the expense of forests, wetlands, and grasslands [8].
Rwanda has witnessed rapid population growth and economic development since the 1990s, which have brought about dramatic changes in the LULC pattern of the country [9]. Rwanda relies heavily on rain-fed agriculture, which makes agricultural production and The population of Rwanda in 2015 was estimated to be about 11.37 million, of which 28.81% resided in urban areas [22]. By 2030, the population is expected to reach about 17 million, with 41.55% living in cities. By 2050, the population is projected to exceed 25 million, with 52% living in urban areas. The topography of Rwanda is characterized by strong landscape heterogeneity; its 30 administrative districts exhibit significant differences in topography, as shown in Figure 1. The terrain of Rwanda can be described as "high west and low east". High mountains and hills (e.g., the Congo-Nile Ridge, Volcanic Highland, and Buberuka Highland) are present in the northwest, and plains (e.g., the Eastern Savanna Region, Eastern Plateau Region, Central Plateau Region, and Bugesera-Mayaga Region) are present in the southeast.

Data Sources
The LULC data used in this study was obtained from the Rwanda LULC product generated by the Regional Centre for Mapping of Resources for Development (RCMRD) in 2015 [23]. The data spans four periods (1990, 2000, 2010, and 2015), and has a spatial resolution of 30 × 30 m 2 , and an average inspection accuracy of 77.50% for all LULC types. The original data contains 14 types of LULC. These 14 types were reclassified to obtain the following land use types: forest, grassland, cropland, wetland, water, urban land, and unused land. The reclassification rules are shown in Table 1. The digital elevation model The population of Rwanda in 2015 was estimated to be about 11.37 million, of which 28.81% resided in urban areas [22]. By 2030, the population is expected to reach about 17 million, with 41.55% living in cities. By 2050, the population is projected to exceed 25 million, with 52% living in urban areas. The topography of Rwanda is characterized by strong landscape heterogeneity; its 30 administrative districts exhibit significant differences in topography, as shown in Figure 1. The terrain of Rwanda can be described as "high west and low east". High mountains and hills (e.g., the Congo-Nile Ridge, Volcanic Highland, and Buberuka Highland) are present in the northwest, and plains (e.g., the Eastern Savanna Region, Eastern Plateau Region, Central Plateau Region, and Bugesera-Mayaga Region) are present in the southeast.

Data Sources
The LULC data used in this study was obtained from the Rwanda LULC product generated by the Regional Centre for Mapping of Resources for Development (RCMRD) in 2015 [23]. The data spans four periods (1990, 2000, 2010, and 2015), and has a spatial resolution of 30 × 30 m 2 , and an average inspection accuracy of 77.50% for all LULC types. The original data contains 14 types of LULC. These 14 types were reclassified to obtain the following land use types: forest, grassland, cropland, wetland, water, urban land, and unused land. The reclassification rules are shown in Table 1. The digital elevation model (DEM) downloaded from the United States Geological Survey (USGS) EarthExplorer [24] database has a 30 × 30 m 2 spatial resolution and was selected to derive the topographic factor. Precipitation data was taken from the Climate Hazards Group InfraRed Precipitation with Station data (CHIRPS), which has been successfully used for scientific research in East Africa [25,26]. CHIRPS data for the period 2000-2015 was used; the spatial resolution of the data was 0.05 × 0.05 • . The soil data, which had a spatial resolution of 250 m, was sourced from the Africa Soil Information Service (AfSIS) [27]. The soil data used in this study included soil texture (percentage of sand, silt, and clay) and soil organic carbon content. Population data were obtained from WorldPop [28]. The resolution of the population data was 1 arc-second (about 30 m at the equator). The vector data for Rwanda Waterways was provided by the National Institute of Statistics of Rwanda (NISR); it was published in 2015 and covers the main waterways in Rwanda [29]. The vector point data for major cities was based on the coordinates (longitude and latitude) of 15 major cities in Rwanda, published by the RCMRD in 2008 [30]. The road data for Rwanda were downloaded from the website of the African Development Bank Group [31]. The data includes national roads in Rwanda, which are further classified as paved and unpaved according to surface conditions. The LULC transfer matrix helps reflect the internal transfer of regional LULC change within a period of interest. This method can comprehensively describe the structural characteristics of regional LULC change and the direction of regional LULC change. This method is derived from the quantitative description of the system state and transition state in system analysis. The transfer matrix is not an index; rather, it only expresses the area of a certain LULC type that has changed in the form of a matrix, which is the basis of LULC structure analysis and changes the square analysis. The mathematical form is as follows [32]: · · · · · · · · · · · · · · · A n1 A n2 A n3 · · · A nn (1) where A is the area, n is the number of LULC types, i (i = 1, 2, . . . n) is the LULC type before transfer, j (j = 1, 2, . . . n) is the LULC type after transfer, and A ij is the area of j when it is transferred from i. of occurrence of events based on Markov process theory. It has the feature of no aftereffect, that is, the current state is only related to the state of the previous moment and has nothing to do with other factors. The land use status can be predicted using the following formula [34]:

Generating Land
where S(T), S(T 0 ) is the land use status at time T and T 0 , respectively and P ij is the transfer probability matrix of land use change.

Land Use Pattern Simulation Using Patch-Generating Land Use Simulation Model (PLUS)
The patch-generating land use simulation (PLUS) model is a multi-type LULC simulation model developed by Liang et al. [35]; it integrates the influences of anthropogenic and natural factors. This model proposes a framework for rule mining based on Land Expansion Analysis Strategy (LEAS). On this basis, a Cellular Automata (CA) model based on multi-type Random Patch Seeds (CARS) was proposed to better simulate the fine-scale patch growth of various LULC types. The working principle of the model is shown in Figure 2. LEAS is used to obtain the land expansion map in this period by superimposing the LULC maps for two periods. Then, the land expansion map and selected driving factors ( Figure 3) are sampled, and the sampled data are used to develop dual-state Random Forest Models. The relationship between the growth of each LULC type and a variety of driving factors is discussed, and a map of the growth probabilities of each land use type and the contribution rate of each driving factor to land expansion is generated. CARS is a preset scenario-driven LULC simulation model that simulates LULC from the macroscale (target LULC amount) to the microscale (each grid LULC type) using a series of competitive mechanisms.      The model simulates LULC changes based on the driving factors behind the dynamics of land use transition, taking into account the background conditions of climate during land use simulation, the interaction and competition between different land use types during simulation training, and other related factors (nature reserves, cities, transportation routes, etc.) on land use simulation. This model can better simulate the LULC pattern with scattered patches. Moreover, research has shown that the PLUS model can simulate LULC and landscape patterns well, and that the simulation accuracy is better than that of models such as the Future Land Use Simulation (FLUS) model, the Conversion of Land Use and its Effects (CLUE) model, and Artificial Neural Networks-Cellular Automata (ANN-CA) [35,36].

Coupling Mechanism of Markov Chain and the PLUS Model
Although the Markov chain can effectively reflect the temporal changes in LULC, it cannot effectively characterize the spatial attributes of the data, nor can it display the spatial change distribution of these elements. Therefore, the method based on the coupling of the Markov chain and the PLUS model can be used to simulate and predict the temporal and spatial distribution of future LULC in the study area. The coupling mechanism is shown in Figure 4. The total amounts of different LULC types in Rwanda in different years (2015,2020,2030,2040, and 2050) were forecasted. Then, the total LULC data for different years were input into the PLUS model, and the spatial distribution of different LULC types in Rwanda in each year was simulated using 2010 as the starting point of the simulation. The simulation results for 2015 were used to verify the model. of the Markov chain and the PLUS model can be used to simulate and predict the temporal and spatial distribution of future LULC in the study area. The coupling mechanism is shown in Figure 4. The total amounts of different LULC types in Rwanda in different years (2015,2020,2030,2040, and 2050) were forecasted. Then, the total LULC data for different years were input into the PLUS model, and the spatial distribution of different LULC types in Rwanda in each year was simulated using 2010 as the starting point of the simulation. The simulation results for 2015 were used to verify the model.

Model Validation
Cohen's kappa coefficient is commonly used to verify data consistency [37]. Theoretically, the calculated results of Cohen's kappa coefficient fall between −1 and 1, but values typically range is between 0 and 1. Values of Cohen's kappa coefficient that are close to 1 indicate a high consistency between the two maps. Generally, values of Cohen's kappa coefficient that are greater than 0.6 indicate substantial consistency between two images. However, because Cohen's kappa coefficient has certain limitations related to the verification of consistency, the figure of merit (FoM) as well as Cohen's kappa coefficient were used to characterize the accuracy of the LULC model in this work.
Proposed by Pontius et al. [38], the FoM is an indicator that reflects unit level consistency and pattern level similarity. It is the ratio of the intersection of the observed and predicted changes to the union of the observed and predicted changes. The calculation formula is as follows: where A is the area of error due to observed change predicted as persistence, B is the area of correct due to observed change predicted as change, C is the area of error due to observed change predicted as wrong gaining category, and D is the area of error due to observed persistence predicted as change. The FoM index is superior to the commonly used kappa coefficient for evaluating the accuracy of simulated changes [38]. Through verification, it has been found that the range of FoM index values generally falls between

Model Validation
Cohen's kappa coefficient is commonly used to verify data consistency [37]. Theoretically, the calculated results of Cohen's kappa coefficient fall between −1 and 1, but values typically range is between 0 and 1. Values of Cohen's kappa coefficient that are close to 1 indicate a high consistency between the two maps. Generally, values of Cohen's kappa coefficient that are greater than 0.6 indicate substantial consistency between two images. However, because Cohen's kappa coefficient has certain limitations related to the verification of consistency, the figure of merit (FoM) as well as Cohen's kappa coefficient were used to characterize the accuracy of the LULC model in this work.
Proposed by Pontius et al. [38], the FoM is an indicator that reflects unit level consistency and pattern level similarity. It is the ratio of the intersection of the observed and predicted changes to the union of the observed and predicted changes. The calculation formula is as follows: where A is the area of error due to observed change predicted as persistence, B is the area of correct due to observed change predicted as change, C is the area of error due to observed change predicted as wrong gaining category, and D is the area of error due to observed persistence predicted as change. The FoM index is superior to the commonly used kappa coefficient for evaluating the accuracy of simulated changes [38]. Through verification, it has been found that the range of FoM index values generally falls between 0 and 0.59 [36], and it is most often below 0.3. When the FoM exceeds 0.3, the consistency between the two images is very high and is almost perfect. Figure 5 indicates that the area of forest, grassland, cropland, wetland, water, urban land, and unused land in Rwanda accounted for 29.64%, 18.13%, 41.51%, 3.77%, 6.06%, 0.82%, and 0.08% of the total area, respectively. Among them, forest, grassland, and cropland are the main LULC types; together, they account for 89.28% of the study area (Appendix A Table A1). As shown in Figure 5, the areas of forest, grassland, and cropland changed significantly from 1990 to 2015 and of the urban land area showed an increasing trend, whereas the area of other land use types showed no obvious change trend. From 1990 to 2015, the forest area gradually decreased with time, the relationship between cropland and forest was a declining relationship, and the grassland area first decreased and then increased. In 1990, the forest area in Rwanda was 10,894.89 km 2 , which was much higher than the cropland area in the same period (6241.79 km 2 ). However, in 2000, the cropland area (11,262.8 km 2 ) was higher than the forest area (7471.05 km 2 ) ( 1990 to 2015, the forest area gradually decreased with time, the relationship between cropland and forest was a declining relationship, and the grassland area first decreased and then increased. In 1990, the forest area in Rwanda was 10,894.89 km 2 , which was much higher than the cropland area in the same period (6241.79 km 2 ). However, in 2000, the cropland area (11,262.8 km 2 ) was higher than the forest area (7471.05 km 2 ) ( Table 2). From 1990 to 2000, the cropland area increased by 5021.01 km 2 , accounting for 19.80% of the national area; conversion from forest accounted for 3525.38 km 2 of this area, while conversion from grassland accounted for 1387.17 km 2 of this area. The ratio of area converted from forest to that converted from grassland was 0.72:0.28. Between 2000 and 2010, 2415.23 km 2 of forest area was converted to cropland, whereas 2632.37 km 2 of cropland area was converted to forest during the same period. The area of cropland converted to grassland was 744.76 km 2 , and the area of grassland converted to cropland was 921.74 km 2 . The area of grassland converted into forest was 335.90 km 2 . From 2010 to 2015, the forest area was reduced to 4322.12 km 2 , resulting in the cropland area reaching 13,446.27 km 2 in 2015 and the grassland area increasing to 4816.66 km 2 . The area converted from grassland to forest was 335.90 km 2 .    Compared with 2010, the cropland area increased by 2280.07 km 2 , the forest area reduced by 3066.17 km 2 , and the grassland area increased by 819.45 km 2 . The net area converted from forest to cropland was 2497.92 km 2 , and the net area converted from forest to grassland was 512.67 km 2 , at a ratio of 0.83:017.

Driving Factors and Changes in Land Use/Land Cover
In this study, 11 indicators (Figure 3) of five driving factors, namely location, soil, population, precipitation, and topography, were selected. Changes in the driving factors of LULC during three phases, namely Phase A (1990-2000), Phase B (2000-2010), and Phase C (2010-2015), were analyzed. As shown in Figure 6, changes in forest, grassland, and cropland are driven by multiple factors such as location, soil, population, precipitation, and topography; By contrast, changes in wetland, water, urban land, and unused land are mainly driven by a single factor, and supplemented by other factors.
In this study, 11 indicators (Figure 3) of five driving factors, namely location, soil, population, precipitation, and topography, were selected. Changes in the driving factors of LULC during three phases, namely Phase A (1990)(1991)(1992)(1993)(1994)(1995)(1996)(1997)(1998)(1999)(2000), Phase B (2000-2010), and Phase C (2010-2015), were analyzed. As shown in Figure 6, changes in forest, grassland, and cropland are driven by multiple factors such as location, soil, population, precipitation, and topography; By contrast, changes in wetland, water, urban land, and unused land are mainly driven by a single factor, and supplemented by other factors. During Phases A and B, precipitation had the highest contribution rate to forest change, followed by population, topography, location, and soil. During Phase C, topography had the greatest contribution rate (25.5%) to forest change, followed by population (23.00%), precipitation (21.84%), location (16.89%), and soil (12.73%). During Phases A and B, population had the highest contribution rate to grassland change, followed in decreasing order by topography, precipitation, soil, and location. During Phase C, precipitation became the largest force driving grassland change, with a contribution rate of 26.45%; it was followed by population (23.97%), topography (21.82%), location (16.54%), and soil (11.22%). The driving factors contributing to cropland change were population, precipitation, topography, location, and soil, in decreasing order of their contribution rate. Population had the highest contribution rates in Phases A (26.89%), B (27.90%), and C (29.19%); its contribution exhibited an increasing trend. Precipitation had the second-highest contribution rates to cropland change in Phases A (20.98%), B (21.65%), and C (23.51%); its contribution also exhibited an increasing trend with time. Topography had the third-highest contribution rates to cropland change in Phases A (23.72%), B (22.91%), and C (18.72%); its contribution rates exhibited a decreasing trend with time. The contribution rates of topography to wetland change during Phases A, B, and C were as high as 41.56%, 43.32%, and 45.54%, respectively; they demonstrated an increasing During Phases A and B, precipitation had the highest contribution rate to forest change, followed by population, topography, location, and soil. During Phase C, topography had the greatest contribution rate (25.5%) to forest change, followed by population (23.00%), precipitation (21.84%), location (16.89%), and soil (12.73%). During Phases A and B, population had the highest contribution rate to grassland change, followed in decreasing order by topography, precipitation, soil, and location. During Phase C, precipitation became the largest force driving grassland change, with a contribution rate of 26.45%; it was followed by population (23.97%), topography (21.82%), location (16.54%), and soil (11.22%). The driving factors contributing to cropland change were population, precipitation, topography, location, and soil, in decreasing order of their contribution rate. Population had the highest contribution rates in Phases A (26.89%), B (27.90%), and C (29.19%); its contribution exhibited an increasing trend. Precipitation had the second-highest contribution rates to cropland change in Phases A (20.98%), B (21.65%), and C (23.51%); its contribution also exhibited an increasing trend with time. Topography had the third-highest contribution rates to cropland change in Phases A (23.72%), B (22.91%), and C (18.72%); its contribution rates exhibited a decreasing trend with time. The contribution rates of topography to wetland change during Phases A, B, and C were as high as 41.56%, 43.32%, and 45.54%, respectively; they demonstrated an increasing trend, indicating that topographic factors were the most important factors limiting wetland change. Topography, precipitation, and location were among the factors driving the change in water. Of these, the contribution rates of topography during Phases A (43.68%), B (55.54%), and C (25.45%) were the highest. Precipitation had the second-highest contribution rates during Phases A (16.32%), B (10.17%), and C (21.17%). The contribution rate of location during Phase C increased sharply to 29.21%, making location the main factor driving the change in water during Phase C. Urban land change was mainly influenced by population. During Phases A, B, and C, the contribution rates of population to urban change were 72.54%, 79.38%, and 55.24%, respectively. Moreover, the contribution rates of precipitation and soil to urban change increased gradually with time. From Phase B to Phase C, there was a change in the main factors driving changes in unused land. The main factor limiting land use change during Phases A and B was precipitation, which had contribution rates of 50.22% and 55.75%, respectively. During Phase C, soil and topography became the main factors limiting the change in unused land. The contribution rates of soil to the change in unused land were 6.35%, 13.23%, and 35.65%, respectively, during Phases A, B, and C, showing a gradually increasing trend.

Simulation of Future Land Use/Land Cover Pattern
Driven by the demand for land use in the future, the PLUS model was used to simulate LULC changes in different years. The simulation results are shown in Figure 7. The results indicate that the actual land use in Rwanda in 2015 (Figure 7a) and the simulated land use (Figure 7b) were very similar. The kappa coefficient of the two LULC images was 0.7325, indicating that the two images were highly consistent. The FoM value was 0.3584, which was much higher than the FoM value obtained by Liu et al. [36] (0.1962), who used the PLUS model to simulate LULC in mainland China. The FoM value obtained herein was also higher than that obtained by Liang et al. [35] (0.2514), who simulated LULC in Wuhan. This indicated that the simulated LULC pattern was almost the same as the actual pattern, and the simulation results were credible. In terms of temporal evolution, the forest area in Rwanda was projected to decrease gradually from 4326.64 km 2 in 2015 to 2332.72 km 2 in 2050. Both grassland and cropland showed an increasing trend, increasing, respectively by 1180.01 km 2 and 999.09 km 2 . Wetlands and water decreased by 343.36 km 2 and 117.66 km 2 , respectively, and urban land expanded from 363.29 km 2 in 2015 to 637.97 km 2 in 2050.
In this study, the agro-ecological zones proposed by Iiyama et al. [1] were used as the standard to calculate the future LULC distribution in Rwanda (Appendix B Figure A1). Iiyama et al. [1] divided Rwanda into six eco-agricultural zones, namely the Eastern Savanna (ES), Eastern Plateau (EP), Bubureka Highland (BH), Central Plateau and Collines (CP), Volcanic Highland (VH), and Congo-Nile Crest (CNC). The main LULC types in ES and EP were grassland and cropland, and both regions demonstrated a trend of increasing grassland and decreasing cropland ( Figure A2a,b). BH and CP primarily consisted of cropland, which demonstrated an increasing trend in both regions. It is noteworthy that the forest land in BH and CP, although not large, exhibited a decreasing trend in both regions ( Figure A2c,d). VH mainly consisted of forest and cropland, followed by urban land; the woodland and cropland area in this region decreased, but the urban land area increased significantly ( Figure A2e). The main LULC types in CNC were cropland and forest. The forest area showed a decreasing trend, whereas the cropland area In terms of temporal evolution, the forest area in Rwanda was projected to decrease gradually from 4326.64 km 2 in 2015 to 2332.72 km 2 in 2050. Both grassland and cropland showed an increasing trend, increasing, respectively by 1180.01 km 2 and 999.09 km 2 . Wetlands and water decreased by 343.36 km 2 and 117.66 km 2 , respectively, and urban land expanded from 363.29 km 2 in 2015 to 637.97 km 2 in 2050.
In this study, the agro-ecological zones proposed by Iiyama et al. [1] were used as the standard to calculate the future LULC distribution in Rwanda (Appendix B Figure A1). Iiyama et al. [1] divided Rwanda into six eco-agricultural zones, namely the Eastern Savanna (ES), Eastern Plateau (EP), Bubureka Highland (BH), Central Plateau and Collines (CP), Volcanic Highland (VH), and Congo-Nile Crest (CNC). The main LULC types in ES and EP were grassland and cropland, and both regions demonstrated a trend of increasing grassland and decreasing cropland ( Figure A2a,b). BH and CP primarily consisted of cropland, which demonstrated an increasing trend in both regions. It is noteworthy that the forest land in BH and CP, although not large, exhibited a decreasing trend in both regions ( Figure A2c,d). VH mainly consisted of forest and cropland, followed by urban land; the woodland and cropland area in this region decreased, but the urban land area increased significantly ( Figure A2e). The main LULC types in CNC were cropland and forest. The forest area showed a decreasing trend, whereas the cropland area exhibited an increasing trend ( Figure A2f). To further identify the spatial change characteristics of the main LULC types in the study area, the hotspots of the changes in area for forest, grassland, cropland, and urban land were identified (Figure 8)  Overall, over the next 30 years, most of the areas of forest loss in Rwanda were projected to be located on either side of the Congo/Nile watershed in the western part of the country and in the VH region in the north. The increase in grassland was projected to occur primarily in the eastern ES and EP regions. The increase in cropland was mainly projected to occur in the BH, PC, and CNC regions in the northwest.

Analysis of Historical Land Use/Land Cover Change
Forest land, grassland and cropland are the main LULC types in the study area, and the total amount between 1990 and 2015 remained essentially constant. However, analysis Overall, over the next 30 years, most of the areas of forest loss in Rwanda were projected to be located on either side of the Congo/Nile watershed in the western part of the country and in the VH region in the north. The increase in grassland was projected to occur primarily in the eastern ES and EP regions. The increase in cropland was mainly projected to occur in the BH, PC, and CNC regions in the northwest.

Analysis of Historical Land Use/Land Cover Change
Forest land, grassland and cropland are the main LULC types in the study area, and the total amount between 1990 and 2015 remained essentially constant. However, analysis of the conversion between different LULC types indicated that the conversion between the three LULC types was quite pronounced ( Table 2).
Between 1990 and 2000, the LULC change consisted mainly of the conversion of grassland and a large amount of forest land to cropland. Of the newly transferred cropland, 72% was converted from forest land and 28% from grassland. This was caused by the surge in population during that period. Wars and turbulence during the early 1990s in Rwanda caused the loss of a significant amount of population by death and displacement. After the wars, the displaced population gradually moved back, and the population began to recover. However, the economic impacts of the war meant that the survival of the returning population could not be effectively guaranteed [39]. Therefore, large portions of forest were cut down for firewood and charcoal to cope with the energy shortage caused by the insufficient production capacity [40]. At present, the vast majority of the rural population (93%) uses firewood for cooking. In cities, charcoal is used by 51% of households, followed by firewood (45%) [41]. The precipitation and temperature conditions unique to Africa cause the rapid deterioration of the surface humus in woodland areas after logging [42], producing fertile soil, which is then converted to cropland. Even in barren, underdeveloped shrublands or woodlands, fire is used to clear the ground in order to produce farmland with high organic matter content [41]. Between 2000 and 2010, there were mutual conversions between forest and cropland, and between grassland and cropland. However, the circulation among the three land types achieved an overall balance, and the total amount of change was not significant. The population of Rwanda continued to increase between 2000 and 2010, but there was no significant reduction in forest area. This was considerably influenced by afforestation during the early part of the 21st century. A large reforestation effort increased the overall forest cover by an average of 8% annually between 2000 and 2005, although virtually none of the resulting forest is classified as primary forest. The demographic structure during this period also contributed to this phenomenon. Although the population growth rate began to increase gradually (from 1.2% to 2.6%) after the post-war population migration [3], the migrating population primarily consisted of women and children; compare with young adult men, such populations do not have a high production capacity [43]. Therefore, there was no major change in the amount of forest land during this period. After 2010, there was a massive conversion of forest to cropland and grassland; this was accompanied by a conversion of cropland to grassland. Of the forest area that was transferred, 83% became cropland and 17% became grassland. This was attributed to the dual effects of a rapidly growing population and a rapidly developing economy.
Owing to the high population pressure, Rwanda is facing an energy crisis, food security problems, and other similar social challenges. The development undertaken to overcome these challenges causes deforestation. The abundant surface vegetation can be used directly as energy material. The newly transformed cropland (from forest) exhibits high nutrient content, good soil structure, and high yield. Moreover, as part of the 'Vision 2020 of Rwanda, its national economy has also witnessed rapid growth. Between 2004 and 2014, the gross agricultural product accounted for 34% of the total annual GDP [44]. This means that the crops produced on the limited arable land not only need to meet the survival needs of the large population, but also need to support the GDP growth as part of the rapid development vision.
The massive conversion of forest into cropland and grassland can also be attributed to the lack of sustainable and effective land management measures. Although the first land policy (i.e., the National Land Policy) of the new century was promulgated in 2004, it focused on land law reform, land use rights, and land fees and taxation management; it did not provide clear provisions for sustainable land development [45]. The lack of management policies for sustainable land development has led to several land problems, such as the excessive fragmentation of agricultural land by families and land degradation [12]. In 2004, the average area of cropland possessed by a Rwandan rural household was about 0.6 ha, and this value could be as low as 0.25 ha. The very high number of families has caused severe fragmentation, which has accelerated the degradation of the land and intensified soil erosion [46]. This kind of unsustainable development of cropland means that forest land is continuously being converted to other LULC types.

Analysis of Factors Driving Land Use Change
The area of forest, grassland, and cropland changed over time, and the factors driving this change in area have also changed gradually. In Rwanda, changes in forest, grassland, and cropland were driven by multiple factors, while changes in wetlands, waters, urban land, and unused land were driven to a greater degree by a single factor than the first three. In 2010 (from Phase B to C), the driving factor of forest changed from 'precipitationpopulation' to 'topography-population'. This change indicated that precipitation was no longer the main factor limiting the change in forest land. Instead, topography and the high population pressure together drove the change in forest. In East Africa, population was the main factor driving LULC change, but the main factor driving forest change did not change to population [6]. This is because low-altitude, scattered forest patches were gradually converted into cropland. Most of the newly added forest land was located in high-altitude areas that were not suitable for human habitation, or in areas such as nature reserves (Nyungwe Forest National Park); population pressure could not be expressed effectively in areas with few people. Moreover, Rwanda has a temperate tropical plateau climate with abundant rainfall and rapid vegetation growth, so precipitation became the main driving factor. The driving factor of grassland changed from 'population-topography' to 'precipitation-population'. Under high population pressure, grassland suitable for cultivation was transferred to cropland (Table A2). Therefore, during Phases A and B, population was the main factor driving grassland change, and topography was the factor limiting the expansion of grassland. Owing to the lack of sustainable land management measures, part of the cropland underwent severe deterioration and was abandoned because it was unsuitable for cultivation (Table A3). After a short period of natural succession, the vegetation in these abandoned cropland areas was restored to grassland vegetation. The natural succession of grasslands was mostly observed in areas with abundant rainfall and dense populations. Therefore, precipitation and population became the new driving forces for grassland change. The driving factor of cropland changed from 'populationtopography' to 'population-precipitation'. Owing to the increase in population pressure, croplands gradually expanded. The expansion of cropland caused a reduction in the forest distributed in high-altitude and high-latitude areas. The topography could no longer limit the expansion of cropland, and thus, topography was no longer the main factor restricting the expansion of cropland. Furthermore, owing to the increase in the elevation of cropland, which was more abundant in the high-altitude areas of Rwanda, precipitation replaced population as the main factor driving cropland change. Wetlands were primarily distributed along river courses, and were therefore primarily found in the low-altitude and relatively flat eastern regions with developed water systems; the fluctuations to their area were small. Thus, topography was the main factor limiting wetland changes.
The driving factor of water area changed from 'topography-population' to 'locationtopography'. Rwanda's water bodies are primarily natural rivers and lakes, which serve as the sources of water used by Rwandans for various activities. Thus, the change in water area is driven by the population. [47].
The cropland in the river valleys is also an important part of the arable land in Rwanda. Rwanda receives abundant rainfall and frequent heavy rains during the rainy seasons, which raises level of water in its rivers. The existence of drainage measures maintains the stability of the arable land area in the river valleys. However, with the growth of population and the intensification of urbanization, the demand for water among the population has increased, and a portion of the agricultural land in the river valleys has been occupied by river channels, or even used for water storage. This is the reason why 14.81 km 2 of cultivated land was converted into water area during 2010-2015 (Table A4).
The driving factor of urban land changed from 'population' to 'population-soil'. Population was undoubtedly the most important factor driving urbanization, but soil became a secondary driving factor because the locations of new cities were not determined solely by population gathering, as before; rather, these locations were the results of more scientific planning. Soil properties are an important factor in the development of cropland. Therefore, the growth of urban land was driven by soil factors.
The driving factor of unused land changed from 'precipitation-location' to 'soiltopography'. This was related to the degradation of cropland owing to the lack of effective cropland management measures under the high population pressure. As a result of early unreasonable reclamation, part of the low-fertility land was developed as cropland. However, with the passage of time, such land was quickly abandoned because it could not meet the demand for food crops. Furthermore, no natural grass or other vegetation was derived before it was transformed into grassland, and this land therefore became unused land, which demonstrated that this change was driven by soil factors.

Analysis of Future Land Use Change
In the next 30 years, LULC changes in Rwanda will mainly consist of reductions in forest area and increments in grassland and cropland area, followed by urban expansion. This is closely related to the gradual increase in the population of Rwanda. Given a certain unit agricultural output, reclamation of forest becomes the best approach to ensure food security [48]. Rwanda's post-war recovery was achieved in this manner in the past, and this approach is expected to continue for the foreseeable future.
Spatially, the future LULC of Rwanda exhibits 'grassland moving eastward' and 'cropland moving westward' trends. The 'westward movement' of cropland is the result of a combination of several factors. The coupling of precipitation, topography, and farming patterns is an important factor. The precipitation and topography of Rwanda, which is located in the interior of the East African Plateau, are high in the west and low in the east (Figure 3i). Moreover, the agricultural system consists primarily of rain-fed agriculture. Most food crops need to be planted and harvested during the two rainy seasons. However, heavy rainfall and frequent floods occur during the rainy seasons, necessitating the implementation of agricultural drainage systems in Rwanda [49,50]. However, the agricultural drainage system in Rwanda is not well developed, and therefore, the natural advantages of the topography are employed to solve this problem. In addition, the western region consists mainly of forest, and cropland that has been converted from forest has more abundant soil nutrients and produces higher yields.
The 'eastward movement' of grassland is caused by natural and socio-economic factors. In terms of natural conditions, the road network in the ES and EP regions in the east is not well developed (Figure 3e,f), and the vegetation consists mainly of sparse grassland, which is not suitable for reclaiming high-quality cropland. With the rapid development of society and economy, the diet of the population has changed, and the demand for meat, dairy products, and animal husbandry products has increased [51,52]. The eastern region has flat topography and numerous pastures. The water system in this region is well developed and the water supply is adequate, which favors the development of animal husbandry.

Conclusions
Based on the analysis of LULC transfer in Rwanda during the historical period , this study found that the LULC in Rwanda was dominated by the decrease in forest land and the expansion of cropland. Moreover, the wetland, water and unused land areas fluctuated with the annual increase in grassland and urban land. Before 2000, LULC transfer in Rwanda primarily consisted of the conversion of forest and grassland to cropland, with 72% of the newly generated cropland coming from forest land, and 28% coming from grassland. From 2000 to 2010, the land transfer was balanced. After 2010, 83% of the transferred forest land was converted to cropland, and 17% was converted to grassland. Changes in forest, grassland, and cropland areas were driven by multiple factors, whereas changes in wetland, water, urban land and unused land areas were dominantly influenced by a single factor. Moreover, the factors driving different LULC types changed with the passage of time.
The PLUS model was applied to simulate the LULC pattern in Rwanda during the next 30 years. The simulation results showed that the PLUS model has high simulation accuracy and generates a reliable LULC pattern. In Rwanda, the LULC will continue to exhibit a trend of decreasing forest and increasing grassland and cropland. Moreover, the model indicated that the cropland in the western part has increased and the forest has decreased. The grassland in the eastern part has increased, but the cropland has decreased. Thus, the trend observed is that of cropland 'moving westward' and grassland 'moving eastward'.  Institutional Review Board Statement: Not applicable.

Informed Consent Statement: Not applicable.
Data Availability Statement: All data used in this article are public record. All data sources have been listed in Section 2.2.

Acknowledgments:
The authors would like to thank Bo Ma and all the reviewers who participated in the review.

Conflicts of Interest:
The authors declare no conflict of interest.  Appendix B Figure A1. Map of the six agro-ecological zones [53]. Figure A1. Map of the six agro-ecological zones [53].