A Deterministic Cellular Automata Model for Simulating Rural Land Use Dynamics : A Case Study of Lake Chad Basin

ISPRS Archive Vol. XXXVIII, Part 4-8-2-W9, "Core Spatial Databases - Updating, Maintenance and Services – from Theory to Practice", Haifa, Israel, 2010 A Deterministic Cellular Automata Model for Simulating Rural Land Use Dynamics: A Case Study of Lake Chad Basin a A. P. Ozah a, *, F. A. Adesina b, A. Dami b, Regional Centre for Training in Aerospace Surveys (RECTAS), PMB 5545, Off Road 1, OAU Campus, Ile-Ife, NIGERIA - azukaozah@yahoo.com b Department of Geography, Obafemi Awolowo University, Ile-Ife, NIGERIA – (tonidamy@yahoo.com, faadesin@yahoo.com) Commission VIII, WG VIII/1 and WG VIII/8 KEY WORDS: cellular automata, grid cell, transition model, land-use, potential model, markov chains, transition rule ABSTRACT Cellular automata (CA) models (deterministic, stochastic or hybrid) have recently garnered tremendous popularity as spatial simulation techniques in a wide range of urban modelling domains and, as such, the vistas of research in this direction are rapidly expanding. Over the past few years, CA models have found application in spatial simulation involving a plethora of themes, including population dynamics, polycentricity, urban land-use evolution, gentrification, urban sprawl and a host of others. Compared to conventional mathematical tools of spatial simulation such as differential equations, partial differential equations and empirical equations, CA models are relatively simple yet produce results that are stunningly meaningful and useful to support decision making in a planning context. Operating in synergy with other planning models and such other cutting-edge technologies as Geographic Information Systems (GIS) and digital image processing, CA can help to portray the dynamics and patterns of growth in a given spatial context. We present a synergistic approach featuring an integration of a deterministic CA model with Markov chains transition model to determine the patterns and dynamics of land-use change in a rural setting. The site chosen for this study is the Lake Chad Basin, an endorheic basin supporting a population of over twenty million people living in four African countries (Chad, Cameroon, Niger and Nigeria) where the people are among Africa's most chronically vulnerable to food insecurity due to the drastic impact of natural and anthropogenic agents of ecological transformation in the basin. Using remotely sensed images obtained at three different dates with at least ten years separating any two dates (1975, 1987 and 1999) together with other supporting attribute data, simulation runs were executed to construct two different future scenarios (2011 and 2023) of rural land-use dynamics. In the simulation analysis, Water, Wetland, Openland and Forest land-use classes were predicted to register net losses while the Road, Settlement and Farm land-use classes were predicted to register net gains in both 2011 and 2023 predicted land-use scenarios. Based on the simulation results, the Lake Chad ecosystem was noted to undergo extensive land-use/land-cover transformation in the future. This study demonstrates that the proposed methodological approach of integrated rural land-use scenario building and analysis relying on the CA-based land-use simulation model possesses an encouraging and exciting prognosis as a technique to support rural land-use planning and policy for sustainable development. The spatio-temporal dynamics of rural ecosystems is an invariably complex phenomenon that involves a complex nexus of interacting forces of causal agents. It has been demonstrated that, in order to disentangle the complex suite of socio-economic and bio-physical forces that influence the rate, spatial pattern and distribution of land use change and to estimate the impacts of such changes, spatially-explicit land use models are indispensable (Costanza, R. et al, 1998). As reproducible tools, such spatial simulation models have the potential to expand the planner’s knowledge domain and to support the exploration of future land use changes under different scenario conditions, thus supplementing his existing mental capabilities to analyze land use change and to make more informed decisions. Such tools can help to predict ecological responses to changing landscape heterogeneity and to gain insights into the variability of landscape patterns and processes over time. 1. INTRODUCTION The issues of land management and land-use change (rural or urban) dominate the development agenda of all countries of the world and, over the years, have remained highly political and contentious. This is expected in the context in which optimal land-use planning is perceived as an indispensable factor for ensuring food security, environmental sustainability and economic development. One of the recent studies on the impact of increased pressure on land and the effects of land-use and land management practices on the dynamic character of rural ecosystems show that a strong correlation exists between balanced, sustainable land development and, human, food and environmental security (Houet, T,. et al, 2006). Other recent studies on land management indicate that, in many parts of the world, the rural landscape is experiencing rapid land-use/landcover (LULC) changes (Kamusoko, C. et al, 2008). It is not surprising therefore that the issue of rural land use change occupies the front-burner of the development initiatives of responsible governments the world over. The desire of planning authorities and municipal governments is to articulate policies and programmes capable of maintaining a balanced ecosystem or to mitigate or prevent, on a sustainable basis, the devastating consequences of severe land-use change when they occur. Achieving this goal requires that such authorities must adopt responsible, holistic and sustainable development strategies which must axiomatically be carried out using spatial decision support tools and methodologies. The research on dynamic landscape modelling using spatial simulation models is quite extensive and varied (Batty and Xie, 1994). These models have been successfully applied in such domains as land-use allocation and land-use planning. The most recent of these spatial simulation models are based on simulated annealing, genetic algorithms, cellular automata (CA), or agent based models. Until this recent development, spatial simulation models were exclusively based on techniques such as differential equations, partial differential equations and empirical equations. CA models (deterministic, stochastic or hybrid) have recently garnered tremendous popularity as a spatial simulation technique in a wide range of urban and rural * Corresponding author. 75 ISPRS Archive Vol. XXXVIII, Part 4-8-2-W9, "Core Spatial Databases - Updating, Maintenance and Services – from Theory to Practice", Haifa, Israel, 2010 ƒ The objects move from one state to the next according to the transition probabilities which depend only on the current state (previous history is not taken into account). ƒ The transition probabilities do not change over time. A homogeneous Markov chain can be represented mathematically as follows (Equations 1 and 2). land use simulation and modelling domains and, as such, the vistas of research in this direction are rapidly expanding. Over the past few years, CA models have found application in spatial simulation involving a plethora of themes, including population dynamics, polycentricity, urban land-use evolution, gentrification, urban sprawl and a host of others. Compared to conventional mathematical tools of spatial simulation such as differential equations, partial differential equations and empirical equations, CA models are relatively simple yet produce results that are quite meaningful and useful to support decision making in a planning context. Although Geographical Information Systems (GIS) are powerful tools to collect, store, manage and analyze spatial data, current GISs have shown considerable weakness and limitation in spatial decision making which are due, in great part, to their lack of sophisticated analytical and spatial modeling tools (Malczewski, 1999; Park et al, 1997). Many studies have shown that the integration of Geographic Information Systems (GIS), cellular automata (CA) models, land use allocation models (such as multi-criteria evaluation) and statistical simulation models (such as Markov chains) provides a powerful environment to simulate and predict dynamic phenomena such as urban and rural spatial growth (Clarke et al, 1997; Fengyun, M., et al, 2005). Q t +1 = Q tT P n T  q1   p11 p11 ... p11   q1         q 2   p 21 p 22 ... p 2m   q2   ...  =  ...   ... ... ... ...        q  p q  p mm   m  t +1  m  t  m1 p m 2 where n m (1) n (2) = number of time steps; = number of states; = vector of initial states at an initial time, t; Qt Q t +1 = vector of states at the next time, t+1; P = transition probabilities matrix. A common practical application of the Markovian simulation model is the prediction of a future land-use scenario given a past history of land-use situations. In this analytical situation, the states are represented by the land-use classes (eg., the classes in a classified satellite image). An example of a Markov transition probability is the probability of one land-use class at an initial time, t changing to another (or remaining in the same) class at a later time, t+1. The crossing of two raster-structured maps representing the terrain situation at two different times is a contingency table (class occupation statistics) of the classes from the two maps. The table is known as the transition cells matrix and can be used to compute the transition probabilities matrix by dividing every entry in a row (class at an initial time) by the total number of cells in that row. In this study, we present a methodological framework that leverages the suitability-based cellular automata land-use simulation model. The proposed approach involves loosecoupling Geographical Information Systems with Markov chains, Multi-criteria Evaluation and Cellular Automata models to model and simulate the rural land-use dynamics using spatial and thematic data sets covering a section of the Lake Chad Basin, an endorheic basin supporting a population of over twenty million people living in four African countries (Chad, Cameroon, Niger and Nigeria) where the people are among Africa's most chronically vulnerable to food insecurity due to the devastating impact of natural and anthropogenic drivers of ecological transformation in the basin. Our proposed methodology features a workflow-centric, three-step process starting from change quantification (transition model using Markov Chains) through change location (potential model using MCE) to change differentiation (cellular automata). A timeseries of multi-scale, multi-resolution and multi-temporal historical satellite imagery and other existing 1:50,000 topographic maps of the study site served as input spatial data sets to determine and delineate the landscape features and the land-use/land-cover trends over the period from 1975 through 1987 to 1999. Bio-physical data (digital elevation model and precipitation data) and accessibility data (distance to roads, distance to water and distance to settlements) served as input suitability spatial factors that were compiled, analyzed and assessed to quantify spatial dependencies using raster-based map algebra and spatial statistical techniques. 1.2 Multi-criteria Evaluation Land-use change is often modelled as a function of a selection of socio-economic and bio-physical variables that act as the “driving forces” or factors of land-use change. Driving forces are generally categorized into three main groups: socioeconomic drivers (e.g. population pressure, income levels and agricultural production) bio-physical drivers (e.g. climatic factors such as rainfall, temperature, humidity and topographic variables such as altitude, slope and aspect) and proximate or accessibility factors (e.g. distance to road, distance to settlement and distance to water). Multi-criteria decision-making (MCDM) or Multi-criteria Evaluation (MCE) problems involve a set of alternatives that are evaluated on the basis of a set of evaluation criteria (Jacob, N et al., 2008). Using MCE, different factors or change drivers can be combined using appropriate weights assigned to the factors. The result of such a combination is a numerical value map representing the land “transition suitability” or “transition potential”, values (scores) that describe the potential of cells to undergo transition from a current state (e.g. “forest”) to a new state (e.g. “built-up”). A number of methods have been proposed for weighting the factors. Examples include direct weighting, rank-order weighting and Analytical Hierarchy Process (AHP). Of these approaches, AHP has been identified as a weighting strategy that can overcome the problem of weighting bias which are obvious short-comings of the direct and rank-order methods (Deekshatulu et. al. 1999). AHP can be used to determine the relative importance of a set of activities or criteria through a pair-wise comparison of the various factors. The first step of the AHP is to form a hierarchy of objectives, criteria and all other 1.1 Markov Chain Model A Markov chain model is defined by a set of states and a set of transitions with associated probabilities where the transitions emanating from a given state define a distribution over the possible next states. The controlling factor in a Markov chain is the transition probability, a conditional probability for the system to undergo transition to a new state, given the current state of the system. One special class of Markov chains that finds wide application in practical simulation problems is the homogeneous Markov chain. A problem can be considered a homogeneous Markov chain if it has the following properties (Anthony, T. L., 2000): ƒ For each time period, every object in the system is in exactly one of the defined states. 76 ISPRS Archive Vol. XXXVIII, Part 4-8-2-W9, "Core Spatial Databases - Updating, Maintenance and Services – from Theory to Practice", Haifa, Israel, 2010 elements involved in the problem. Once the hierarchical structure has been formed, comparison matrices are developed as results of evaluations made by the decision-makers on the intensity of difference in importance, expressed as a rank number on a given numerical scale for each level in the hierarchy. This forms the basis of the final computation of the factor weights. An elaborate description of the MCE process by AHP can be found in Jacob N., et al.(2008). In a typical practical application, the factors (input) into the process are spatial factors or factor maps (eg. aspect map, slope map, distance map, population surface map, etc), spatial constraints (e.g. a binary map of the exclusion zones such as maps of underground treatment plants, water canal, landfill sites, etc) and non-spatial constraints. Some MCE modelling methods require that the factors be standardized prior to the computation of the final transition suitability (potential) map (ILWIS, 2007). t N i , j = neighbourhood of a cell, Ci , j at time t; t R i , j = transition rule applied to cell, Ci , j at time t. Transition rules: Transition rules guide and control the dynamic evolution of CA. In classical CA, transition rules are deterministic and unchanged during evolution. In several recent studies, however, they are modified into stochastic and fuzzy logic controlled methods (Wu, 1998). Neighbourhood: This is defined by the local neighbours of a cell. In a two-dimensional cellular automata model there are two common types of neighbourhood: the Von Neumann neighbourhood with four neighbouring cells and the Moore neighbourhood with eight neighbours (see Figure 1 below). 1.3 Cellular Automata Cellular Automata (CA) are dynamic models that can be employed to simulate the evolution or dynamics of a wide variety of natural and human systems. They are processing algorithms that were originally conceived by Ulam and Von Neumann in the 1940s to study the behavior of complex systems (von Neumann, 1966). CA models present a powerful simulation environment represented by a grid of space (raster), in which a set of transition rules determine the site attribute of each given cell taking into account the attributes of cells in its vicinities. These models have been very successful in view of their operationality, simplicity and ability to embody both logic and mathematics-based transition rules, thus enabling complex global patterns to emerge directly from the application of simple local rules. A cellular automaton system consists of a regular grid of cells, each of which can be in one of a finite number of k possible states, updated synchronously in discrete time steps according to a local, identical interaction (transition) rule. The state of a cell is determined by the previous states of a surrounding neighborhood of cells. The types of spatial problems that can be approached using CA models include spatially complex systems (e.g., landscape processes), discrete entity modeling in space and time (e.g., ecological systems, population dynamics) and emergent phenomena (e.g., evolution, earthquakes). From the application perspective, CA are dynamic models that inherently integrates spatial and temporal dimension. (a) (b) Figure 1. 3 x 3 Neighbourhood Kernels showing: (a) Moore neighbourhood and (b) Von Neumann neighbourhood The future state of a cell in a CA is dependent on its current state, neighborhood states, and transition rules which are setup and fine-tuned using transition suitability or potential scores of individual cells. Iterative local interaction between cells within the neighborhood finally produce the global pattern. 2. STUDY SITE The site for this study is a section of the Lake Chad Basin, an endorheic basin located approximately between latitudes 12oN and 14o 30’ N and Longitudes 13o E and 15o 30’ E. The basin is shared by four West African countries namely, Chad, Cameron, Nigeria and Niger and supports a population of over twenty million people. The Lake basin comprises five bio-climatic zones, namely, Saharan, sahelo-saharan, sahelo-sudanian, dudano-sahelian and sudano-guinea ecological zones. The southwest humid Atlantic (monsoon) and the north-east Egyptian hot and dry (harmattan) currents influence the climate and consequently the ecological zonation of the basin. The sudanoguinean climate in the south for example has annual rainfall of over 950 mm, a rainy season of six to seven months (May November) with an average annual temperature at Sarh of 28oC (absolute minimum 10oC, absolute maximum 450) and annual Piche-recorded evaporation of 2027mm in 1961 (Thambyahphillay, G. G. R., 1987). The topography of the lake basin can be described as generally flat with a few shallow depressions and a few widely scattered elevated spots. In terms of hydrogeology, the basin lies in a tectonic zone with an extensive sedimentary basin where depositional events, resulting in the formation of four aquifers, had taken place in tertiary and quaternary times (LCBC, 2007). The soil characteristics of the Chad basin region are Ferruginous Tropical and undifferentiated semi-arid brown soils. These cover about two-fifths of the basin while the remaining 60 percent is covered by a zonal vertisols, regosols and mixtures of alluvial and vertisols characterized by a high shrink-swell potential. The Nigerian section of the basin has been rated as soil having about 90 percent potentiality of medium to high fertility. The predominant vegetation of the basin is comprised of the woodland and pseudo-steppe types populated with trees and shrubs. Figure 2 shows a map of the study site depicting a section of Lake Chad Basin chosen for this study. CA is composed of a quadruple of elements as defined in the following equation (White and Engelen, 2000). CA = {X, S, N, R} (3) where CA = cellular automaton; X = CA cell space; S = CA states; N = CA cell neighbourhood; R = CA transition rule; Cell space: The cell space is composed of individual cells. Although these cells may be in any geometric shape, most CA adopts regular grids to represent such space, which makes CA very similar to the cellular structure of raster GIS. Cell states: The states of each cell may represent any spatial variable, e.g., the various land-use types. The state transition of a CA is defined by the following relation: t +1 S i , j = f (( t S i , j ),( t N i , j ),( t R i , j )) where t +1 (4) S i , j = new (next) state of a cell, Ci , j at time t+1; t S i , j = initial state of a cell, Ci , j at time t; 77 ISPRS Archive Vol. XXXVIII, Part 4-8-2-W9, "Core Spatial Databases - Updating, Maintenance and Services – from Theory to Practice", Haifa, Israel, 2010 Data Landsat TM Landsat TM Landsat ETM+ Topographical map ASTER GDEM CGIAR-CRU Rainfall Source GLCF1 GLCF1 GLCF1 OSGOF2 ASTER3 CGIAR-CRU4 Scale/Resolution 30m 30m 28.5m 1:50,000 30m 0.5° Date 1975 1987 1999 1965 2009 1999 Table 3. Spatial data sets employed in the study CHAD 1975 Landsat image 1987 Landsat image 1999 Landsat image Figure 3. Source Landsat images of 1975, 1987 and 1999. CHAD 3.2 Materials CAMEROUN NIGERIA Although several off-the-shelf GIS and digital image processing software packages exist and have been successfully used by several researchers to implement CA-based spatial simulations (see Houet, T,. et al, 2006), such tools were not available to us during the course of this study. However, an open-source digital image processing and GIS package, Integrated Land and Water Information System (ILWIS, 2007) was found adequate for the execution of most of the data processing, conversion, integration and presentation tasks. Special tasks and functionalities required for the simulation but are not supported in the ILWIS software were implemented using in-house programs developed in Visual Basic 6.0 environment. Figure 2. Map of the study site showing a section of Lake Chad Basin (Source: http://www.answers.com) 3. DATA, MATERIALS AND METHODS 3.1 Data Sources 3.3 Methodology The simulation of rural land-use dynamics undertaken in this study for the chosen site was intended to statistically and spatially associate future land-use scenarios with historical growth patterns in the study area by employing the site attributes (bio-physical data) covering the area. Consequently, the following spatial data types were acquired: a set of satellite image data acquired at three different dates; topographic maps; digital elevation data; and climatic data. A set of three Landsat image scenes acquired at three different dates (1975, 1987 and 1999) were downloaded from the Global Land Cover Facility (GLCF) digital image archive at http://www.glcf.com. The 1975 Landsat TM image was obtained in the mosaicked format with a spatial resolution of 30m. The 1987 Landsat TM had a spatial resolution of 30m while the 1999 Landsat ETM+ image had a spatial resolution of 28.5m. Scanned copies of the 1:50,000 topographic maps covering the study site were obtained from the Office of the Surveyor-General of the Federation (OSGOF) of Nigeria. An appropriate image tile of the Advanced Spaceborne Thermal Emission and Reflection Radiometer (ASTER) GDEM was downloaded from global data server at http://www.gdem.aster.ersdac.or.jp. Similarly, a tile of CGIARCRU rainfall data covering the study site was downloaded from the CGIAR-CRU Global Climate Database at http://cru.csi.cgiar.org. Table 3 lists the source data sets with their scale/resolution properties and dates of acquisition while Figure 3 shows the Landsat images as obtained. Data preparation, conversion and processing: Our proposed framework for simulating rural land-use dynamics requires five input maps as summarised in Table 4 below. Input map Initial configuration Transition suitability map Suitability threshold map Predominant count map Predominant class map Description Map showing the initial (seed) scenario for the simulation. Map showing transition potential values for each land-use class. Attribute map depicting suitability cut-off values for each land-use class. Map depicting the count of the predominant class within the neighbourhood of each cell of the original (initial configuration) map. Map showing the class with the predominant count within the neighbourhood of each cell of the original (initial configuration) map. Table 4. Input maps for CA-based simulation The first step in our methodological approach was the preprocessing of input maps and tables for the simulation. All the required data sets were first imported into the ILWIS environment and converted into ILWIS format. The CGIARCRU Rainfall data (point raster format) was then interpolated using the Inverse Distance Weighted (IDW) method and resampled to obtain a 30-m resolution surface map. All the three satellite images, including the scanned topo maps, the ASTERGDEM and the CGIAR-CRU rainfall surface map were then geo-referenced to the UTM (Zone 33) projection on WGS84 Ellipsoid using ground control (tie-points) read off from the topo maps. The required slope and aspect maps were then derived from the ASTER-GDEM (altitude) map using raster map calculation. The 1975 (mosaicked) source image was observed to show marked differences in spectral characteristics due 1 1 Office of the Surveyor-General of the Federation, Nigeria. Global Land Cover Facility (www.glcf.com). Advanced Spaceborne Thermal Emission and Reflection Radiometer (http://www.gdem.aster.ersdac.or.jp) 4 Consultative Group for International Agriculture Research - Consortium for Spatial Information: Global Climate Database (http://cru.csi.cgiar.org). 2 3 78 ISPRS Archive Vol. XXXVIII, Part 4-8-2-W9, "Core Spatial Databases - Updating, Maintenance and Services – from Theory to Practice", Haifa, Israel, 2010 this operation are as summarized in Tables 7 (a), (b) and (c) below. From Tables 7(b) and (c), the suitability threshold map was generated by interpolating the values in these tables (column F) from the cumulative histograms of the individual classes for the two scenarios (Table 6). This operation was performed using an in-house program developed in Visual Basic 6.0. A composite suitability threshold map was finally generated by assigning the computed threshold values to each of the land-use classes. probably to disparities in image acquisition conditions. Since this could lead to erroneous classification results, this image was first unglued to obtain its component parts. In the study area, seven land-use classes (Road, Settlement, Water, Wetland, Openland, Farm and Forest) were identified and considered as candidate “states” for the simulation. However, all the source images were classified into only five land-use classes (Water, Wetland, Openland, Farm and Forest). The Road and Settlement classes were omitted from the classification since they showed very similar spectral characteristics with the Openland class. The separate components of the 1975 image were first classified in the ILWIS environment using the maximum likelihood classifier and finally glued back together to obtain a seamless image. The same digital classification parameters were then applied to classify the 1987 and the 1999 images. To include the Road and Settlement classes, the two feature types were extracted from the topo maps by on-screen digitizing and overlaid on each of the three original images. The corresponding polygons (Settlement class) and segments (Road class) were then updated to reflect the situation at each period. Each of the Road and Settlement layers so delineated were then rasterized and merged into the corresponding classified image using raster map calculation function available in ILWIS. Finally, a spatial constraint map layer was created by digitizing the Water Canal features from the topo maps. Similarly, three accessibility maps (distance to road, distance to settlement and distance to water) were generated from the road, settlement and water layers using the distance calculation function in ILWIS. An extract of the defined study area was thereafter made from each of the resulting images and stored for further analysis. The final 1975, 1987 and 1999 classified images are shown in Figure 4 below: (a) (b) The next step in the simulation process was the computation of the potential model (using MCE) to generate the transition suitability (potential) maps. In this study, seven spatial factors (altitude, slope, rainfall, distance to road, distance to settlement and distance to water) and one spatial constraint (water canal) were considered to affect the suitability of cells for conversion into other land-use classes. Socio-economic factors such as population density and agricultural productivity were not available during the course of this study. For each of the available factors, a factor map was prepared. This process was effected in ILWIS by first creating a factor group for each landuse class and then populating the group with the maps corresponding to each of the factors. This was followed by a process of standardization of each map considering the overall goal of the evaluation in relation to the contribution of each factor towards the goal. The factors per group were then weighted using the Analytical Hierarchy Process (AHP) by pairwise comparison of factors (see Saaty, 1980 for details). This resulted in specific weights assigned to each factor (see Table 5). The AHP was assessed using the inconsistency ratios and the computed weights were accepted and used in the computation of a composite index map (suitability map) for each factor. Thus seven different suitability maps corresponding to the seven different land-use classes were obtained from the MCE-AHP process. To obtain a single suitability map, all the seven suitability maps were combined using the map calculation functionality in ILWIS. The final transition suitability map was then generated using cell-by-cell multiplication of the composite suitability map and the binary constraint map. This operation assigned a suitability value of zero (unsuitable) to all cells corresponding to the water canal feature. (c) Figure 4: Classified images of : (a) 1975; (b) 1987 and (c) 1999 In this study, we adopted the classified 1999 image (latest image) as the initial (seed) configuration for the CA simulation. Based on this classified image, the neighbourhood predominant count map and the neighbourhood predominant class maps were generated by applying the neighbourhood map calculation functionalities supported in the ILWIS software based on the Moore neighbourhood (Figure 1(a)). To obtain the transition cells matrices, two successive classified images were crossed, giving the class occupation statistics of the two image pairs. The 1975-1987 and 1987-1999 transition cells matrices so obtained were then employed to compute the corresponding transition probabilities matrices as described in Section 1.1 of the paper. Using the method described in Berchtold A. (1998), a homogeneous transition matrix representing the general trend of land-use evolution over the period between 1975 and 1999 was computed based on the two separate transition matrices. The homogeneous transition matrix and the class populations (initial states vector) from the seed image constituted input variables for the computation of the transition model. Factor F1 F2 F3 F4 F5 F6. F7 TW IR C1 0.140 0.071 0.071 0.045 0.316 0.231 0.127 1.000 0.099 C2 0.036 0.073 0.085 0.025 0.278 0.158 0.346 1.000 0.089 Land-use Classes C3 C4 C5 0.045 0.116 0.053 0.070 0.089 0.053 0.029 0.054 0.053 0.315 0.325 0.333 0.096 0.046 0.109 0.131 0.044 0.109 0.315 0.325 0.291 1.000 1.000 1.000 0.090 0.049 0.059 C6 0.032 0.049 0.028 0.335 0.154 0.116 0.287 1.000 0.099 C7 0.047 0.053 0.026 0.431 0.115 0.192 0.134 1.000 0.090 Table 5. Multi-Criteria Evaluation factor weights for different Landuse classes computed based on Analytical Hierarchy Process (AHP) (Classes: C1=Road, C2=Settlement, C3=Water, C4=Wetland, C5=Openland, C6=Farm, C7=Forest; Factors: F1= Altitude, F2=Slope, F3=Aspect, F4=Rainfall, F5=Distance to road, F6=Distance to settlement, F7=Distance to water; TW=Total Weight; IR=Inconsistency Ratio) Suitability Threshold (2023) LU Class Suitability Threshold (2001) 0.916626547 0.912322272 Road 0.962713488 0.956168315 Settlement 0.757331349 0.746932254 Water 0.73871287 0.67372059 Wetland 0.815529404 0.788445542 Openland 0.787391162 0.777337315 Farm 0.737500963 0.730205614 Forest Table 6. Predicted transition suitability threshold values for 19992011 and 1999-2023 simulations The transition (quantification of land-use changes) model for the simulation was then executed to obtain simulated rates of change corresponding to proposed future scenarios (2011 and 2023). Adopting the latest period (1999) as the initial (seed) period, the homogeneous Markov chain model (Equation 2, Section 1.1) was applied (substituting the homogeneous transition matrix for P and the 1999 class population values for Q) to compute simulated transition rates (cells) for 2011 (n=1 or first time step) and 2023 (n=2 or second time step). The results obtained from 79 ISPRS Archive Vol. XXXVIII, Part 4-8-2-W9, "Core Spatial Databases - Updating, Maintenance and Services – from Theory to Practice", Haifa, Israel, 2010 Class populations 2011 45216 2186508 900545 1900899 6584558 6196237 4582088 LU Classes Road Settlement Water Wetland Openland Farm Forest 1999 40596 1661531 984230 2977357 7191009 4202185 5339142 Class C1 C2 C3 C4 C5 C6 C7 A 40596 1661538 984234 2977370 7191041 4202204 5339166 B 5401 599371 189062 1551940 1056422 4606912 2424867 C 781 74394 272747 2628398 1662872 2612860 3181922 Class C1 C2 C3 C4 C5 C6 C7 A B 4059 6 1661538 984234 2977370 7191041 4202204 5339166 11415 1187798 307434 1587111 2040933 4783369 2445938 Algorithm 1: Basic-like CA-based transition rule for the simulation --------------------------------------------------------------------------------------If (Suitability_Map = 0) Then Simulated_Map = Initial_Map; ‘No change in cell state Else If ((Initial_Map = "Settlement") Or (Initial_Map = "Road")) Then Simulated_Map = Initial_Map; ‘No change in cell state Else If ((Suitability_Map >= Suitability_Threshold_Map) And (Predominant_Count >= 3)) Then Simulated_Map = Predominant_Class; ‘ Update cell state Else Simulated_Map = Initial_Map; ‘No change in cell state End if End if End if 2023 50465 2706862 825165 2018867 6541046 6138675 4114970 (a) D 4620.1 524976.4 -83685.2 -1076458.1 -606450.8 1994052.0 -757054.4 E 13.3 36.1 19.2 52.1 14.7 109.6 45.4 F 1.9 4.5 27.7 88.3 23.1 62.2 59.6 G 11.4 31.6 -8.5 -36.2 -8.4 47.5 -14.2 C D E F G 1546 142467 466499 2545601 2690895 2846879 3670110 (c) 9868.9 1045331.2 -159065.3 -958490.2 -649962.5 1936490.2 -1224172.3 28.1 71.5 31.2 53.3 28.4 113.8 45.8 3.8 8.6 47.4 85.5 37.4 67.7 68.7 24.3 62.9 -16.2 -32.2 -9.0 46.1 -22.9 ----------------------------------------------------------------------------Simulation Evaluation: In order to determine the degree of reliability of the simulation results, a numerical evaluation was conducted using the historical image data sets for the three periods (1975, 1987 and 1999). This task was done by simulating transition rates for 1987 and 1999 using the 1975 and 1987 class populations respectively based on Markov chain prediction. For each simulation scenario, the accuracy in the estimation of each land-use class was computed using the following formula:   Simulated Value - Real Value   Accuracy(%) = 1001 - abs   Real Value    The results obtained for the two scenarios are as presented in Table 8 below. (b) Table 7. Simulation statistics: (a) Class populations for seed and simulated configurations (b) 1999-2011 and (c) 1999-2023 (A=number of cells in the seed map, B=gain in cells between initial and simulated maps, C= loss in cells between initial and simulated maps, D= net gain in cells between initial and simulated maps, E= % gain in cells between initial and simulated maps, F=% loss in cells between initial and simulated maps, G= % net gain in cells between initial and simulated maps; C1=Road, C2=Settlement, C3=Water, C4=Wetland, C5=Openland, C6=Farm, C7=Forest) LU Class C1 C2 C3 C4 C5 C6 C7 Cellular Automata Simulation Runs: In keeping with the goal of simulating the rural land-use dynamics for our study area, CA-based transition rules (algorithms) were implemented in an ILWIS-based script. In designing the CA transition rule, the Moore neighbourhood kernel (see Section 1.3) was adopted with a kernel threshold of 3. A cell would therefore undergo state transition to the state of the predominant cell in its 8-cell neighbourhood if its transition suitability value is greater than zero and if it is not “Settlement” or “Road” and if its transition suitability value is greater than the suitability threshold value and if the count of the predominant cell is greater than or equal to the neighbourhood kernel threshold value of 3. It is to be noted that several other constraints can be integrated into the transition rule to achieve some set configuration. Algorithm 1 shows a simple Basic-like CA-based transition algorithm developed for the simulation. The script was designed to reference the five raster-structured maps generated in the previous sub-section as input maps to compute final (simulated) raster maps corresponding to the 1999-2011 and 1999-2023 scenarios based on neighbourhood map calculation. Figure 5 shows the resulting simulated maps for the two proposed scenarios. Simulation for 1987 Real Predicted Acc(%) 38378 35177 91.7 1135212 941912 83.0 882703 1074995 78.2 1148271 2051885 21.3 8016597 7109876 88.7 8561286 6642613 77.6 2613603 4539592 26.3 Simulation for 1999 Real Predicted Acc (%) 40596 39781 98.0 1661531 1564405 94.2 984230 957563 97.3 2977357 2109228 70.8 7191009 7016688 97.6 4202185 6429666 47.0 5339142 4278720 80.1 Table 8: Simulation evaluation results (C1=Road, C2=Settlement, C3=Water, C4=Wetland, C5=Openland, C6=Farm, C7=Forest) 4. RESULTS AND DISCUSSION The quantitative results of the integrated spatial simulation of rural land-use dynamics undertaken in this study are presented in Tables 6 and 7 while its graphical (image) outputs are presented in Figure 5. As shown in Tables 7 (b) and (c), Water, Wetland, Openland and Forest classes were predicted to register net losses of 8.5%, 36.2%, 8.4% and 14.2% respectively for the 2011 scenario and 16.2%, 32.2%, 9.0% and 22.9% respectively for the 2023 scenario, while Road, Settlement and Farm classes were predicted to register net gains of 11.4%, 31.6%, and 47.5% respectively for the 2011 scenario and 24.3%, 62.9% and 46.1% respectively for the 2023 scenario. As can be discerned from Table 8, the accuracies of the simulation of the various land-use change rates are better for the 1999 scenario (long-tern projection) than that for the 1987 scenario (short-tern projection). This situation needs to be further investigated. From the initial and simulated images presented in Figure 5, the morphological changes in the land-use classes between the original image and its simulated versions are clearly visible. 5. CONCLUSION AND RECOMMENDATIONS Seed image (1999) Simulated (2011) An integrated methodological approach featuring the coupling of GIS with Markovian, MCE and CA models for modelling and simulating the spatio-temporal dynamics in a rapidly changing Simulated (2023) Figure 5. Seed image (1999) and simulated images (2011 and 2023). 80 ISPRS Archive Vol. XXXVIII, Part 4-8-2-W9, "Core Spatial Databases - Updating, Maintenance and Services – from Theory to Practice", Haifa, Israel, 2010 ecosystem (the Lake Chad Basin) was presented in this study. The proposed integrated model was successfully used to model, analyze and construct future land-use scenarios based on empirical, ground-truth spatial data sets acquired over a period of twenty-four years (1975-1999). The results of the simulation analysis indicate that the Lake Chad ecosystem is steadily undergoing land-use/land-cover changes. In particular, the study reveals that the water stock of the lake is rapidly shrinking. LCBC, 2007. Regional Roundtable on sustainable development of the Lake Chad Basin, University of Maiduguri, Nigeria. Malczewski, J., 1999. Spatial multicriteria decision analysis. In: J.-C. Thill (ed.), Multicriteria Decision Making and Analysis: A Geographic Information Sciences Approach. Brookfield, VT, Ashgate Publishing, pp. 11-48. Mu Fengyun and Zhang Zengxiang, 2005. Cellular automata model based on GIS and urban sprawl dynamics simulation, Proceedings of SPIE, the International Society for Optical Engineering. This study was conducted on only a section of the basin. To perform a more comprehensive analysis of the causes and consequences of the land-use dynamics in the basin, the study area needs to be extended to cover the entire basin. Socioeconomic factors that were not available for the determination of the transition potential values also need to be integrated in future studies to enhance the realism of the simulation. Park, S., and Wagner, D.F., 1997, Incorporating cellular automata simulators as analytical engines in GIS, Transactions in GIS 2, pp. 213-231. The results obtained from this study demonstrate that integrated rural land-use scenario building and analysis relying on the CAbased land-use simulation model can support land-use planning and policy for sustainable land development. However, issues concerning simulation evaluation, calibration and validation need to be further considered and investigated. Saaty, T. , 1980. The Analytical Hierarchy Process. New York, McGraw Hill. Tamara Lynn Anthony, 2000. Markov Chains, http://ceee.rice.edu/Books/LA/markov/markov2.html URL: Thambyahphillay, G. G. R., 1987. Meteorological and Climatological Perspective of Drought and desertification in the Lake Chad Basin of Sahelo-Soudan Nigeria. Paper presented to the Chad Basin Commission’s International Seminar on “Water Resources in the Lake Chad Basin: Management and Conservation” N’Djamena (Republic of Chad): 3rd-5th June, 1987. REFERENCES Batty M, Couclelis H, Eichen M, 1997. Urban systems as cellular automata, Environment and Planning B: Planning and Design 24, pp. 159-164. Berchtold A., 1998. Chaînes de Markov et modèles de transition: Application, aux sciences sociales, Hermes, 284 pages. von Neumann, John, and A.W. Burks, eds. 1966. Theory of selfreproducing automata, Urbana-Champaign: University of Illinois Press. Clarke K C, Hoppen S, Gaydos L, 1997. A self-modifying cellular automaton model of historical urbanization in the San Francisco Bay area, Environment and Planning B: Planning and Design 24, pp.247-261. Wu F,Webster C J, 1998. Simulation of land development through the integration of cellular automata and multi-criteria evaluation, Environment and Planning B: Planning and Design 25, pp.103-126 Costanza R. and Ruth M., 1998. Using dynamic modelling to scope environmental problems and build consensus, Environmental Management 22, pp. 183–195. Deekshatulu, B.L., Krishnan, R., Novaline, J., 1999. Spatial Analysis and Modelling Techniques – A Review, Proceedings of Geoinformatics - Beyond 2000 conference, Published by Indian Institute of Remote Sensing (1999), DehraDun, India, pp.268275. Houet Thomas, Hubert-Moy Laurence, 2006. Modelling and Projecting Land-Use and Land-Cover Changes with a Cellular Automaton in Considering Landscape Trajectories: An Improvement for Simulation of Plausible Future States, EARSeL eProceedings 5, pp 63 – 76. ILWIS, 2007. ILWIS Users Guide, ITC, The Netherlands. Jacoba Novaline, Krishnan R, Prasada Raju PVSP, Saibaba J, 2008. Spatial and Dynamic Modeling Techniques for Land Use Change Dynamics Study, The International Archives of the Photogrammetry, Remote Sensing and Spatial Information Sciences. Vol. XXXVII. Part B2. Beijing. Kamusoko Courage, Aniya Masamu, Adi Bongo and Manjoro Munyaradzi, 2008. Rural sustainability under threat in Zimbabwe – Simulation of future land use/cover changes in the Bindura district based on the Markov-cellular automata model, URL: www.sciencedirect.com 81