A multi-agent model system for land-use change simulation

Environmental Modelling & Software 42 (2013) 30e46

Contents lists available at SciVerse ScienceDirect

Environmental Modelling & Software journal homepage: www.elsevier.com/locate/envsoft

A multi-agent model system for land-use change simulation Célia G. Ralha a, *, Carolina G. Abreu a, b, Cássio G.C. Coelho a, Alexandre Zaghetto a, Bruno Macchiavello a, Ricardo B. Machado c a

Computer Science Department, Institute of Exact Sciences, University of Brasília, P.O. Box 4466, Zip Code 70.904-970 Brasília, DF, Brazil The Brazilian Institute of Environment and Renewable Natural Resources (IBAMA), SCEN Trecho 2, Ed. Sede, P.O. Box 9566, Zip Code 70.818-900 Brasília, Brazil c Zoology Department, Institute of Biological Sciences, University of Brasília, Campus Darcy Ribeiro, Zip Code 70.910-900 Brasília, Brazil b

a r t i c l e i n f o

a b s t r a c t

Article history: Received 14 June 2012 Received in revised form 1 November 2012 Accepted 4 December 2012 Available online 12 February 2013

This paper presents a multi-agent model system to characterize land-use change dynamics. The replicable parameterization process should be useful to the development of simulation frameworks, important to environmental policy makers to analyze different scenarios during decision making process. The methodological two-fold approach intends to form a solid backbone based on: (i) the systematic and structured empirical characterization of the model; and (ii) the conceptual structure definition according to the agent-based model documentation protocol e Overview, Design concepts and Details. A multiagent system for land-use change simulation was developed to validate the model, which is illustrated with a case study of the Brazilian Cerrado using LANDSAT ETM images. The simulation results prove the model importance with a figure of merit greater than 50%, what means the amount of correctly predicted change is larger than the sum of any type of error. The results are very good compared with nine popular peer-reviewed land change models. Ó 2012 Elsevier Ltd. All rights reserved.

Keywords: Agent-based model Brazilian Cerrado Environmental simulation Model description Multi-agent system MASE

Software availability Name of software: MASE e Multiagent System for Environmental Simulation Developer: Cássio Giorgio Couto Coelho (InfoKnow-Computer Systems for Information and Knowledge Treatment Group, ComNet e Computer Networks Lab) Contact address: Computer Science Department, University of Brasília, Campus Universitário Darcy Ribeiro - Asa Norte e Edifício CIC/EST, Caixa postal 4466 e CEP 70.910-900 Brasília e DF e Brazil Email:[email protected], Ralph. [email protected], [email protected] Availability online for free download at: https://sourceforge.net/p/ mase-unb Year first available: 2013 Hardware required: Dual-core Processor 2 GB RAM e 1 GB HDD (minimum) Software required: JRE 7.0 or superior (tested on Windows 7 Home Premium 64bit)

* Corresponding author. Tel.: þ55 61 31073675; fax: þ55 61 31070487. E-mail address: [email protected] (C.G. Ralha). 1364-8152/$ e see front matter Ó 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.envsoft.2012.12.003

Programming language: JAVA 7.0 Agent Development Framework: JADE e Java Agent Development Framework, version 4.0, release date 04/20/10 Executable and Libraries: 32.2 MB Resources (images): 181 MB 1. Introduction Brazil is a nature privileged country. It houses various biomes including the Amazon, the largest tropical rainforest in the world and the largest biome in Brazil, occupying almost half (49.29%) of the country (4,196,943 km2) (IBGE, 2004). Despite this privilege, historically environmental information in Brazil has not been taken into account during planning and decision making processes. This common critique is usually addressed to traditional solutions in the way decision makers treat Land Use/Cover Change (LUCC), among other environmental issues. Recent changes in the Brazilian Forest Code, demonstrates the effects of policies on the environment are barely identified, assessed or even discussed by the competent authorities in the decision making process. A clear example of this gap between environmental information and political action is the Brazilian Cerrado (woodland savanna). Cerrado is the second largest biome in South America and covers approximately 22% of the Brazilian territory, or 204.7 million hectares

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

(Sano et al., 2010; Oliveira, 2002). According to Cavalcanti and Joly (2002) and Myers et al. (2000), the Cerrado harbors outstanding biodiversity, being included in the world’s 25 principal hotspot areas, with great endemism and less than 30% remaining natural vegetation. Unfortunately, the biologically diverse, more than 12,356 catalogued plant species as presented in Rezende et al. (2008), is poorly protected. Only 2.2% of Cerrado’s original area is protected by formal reserves (Klink and Machado, 2005). The Cerrado savanna of central Brazil has been undergoing rapid transformation to cattle ranching and more recently to soy production, requiring urgent action for the conservation of biological richness in the face of threats of destruction (McAlpine et al., 2009). For the cited reasons, the Federal District was selected for our case study since it is the single state of Brazil fully covered by the Cerrado biome. Changes in LUCC are relevant not only for the Cerrado Biome, since they are amongst the most pervasive and important sources of recent alterations of the Earth’s land surface (Houet et al., 2010; Smith et al., 1998). According to Vitousek (1997), LUCC is one of the most profound human-induced alterations on the environment. Its research aims to support insightful management of land resources in order to avoid irreversible damage (Le et al., 2008). There is now more than a decade of international research on land-use processes, and as a result, a much clearer understanding of these processes and a much better appreciation of their complexity and comprehensiveness (Schaldach et al., 2011; Lambin and Geist, 2006). LUCC can be summarized as a complex process caused by the interaction between natural and social systems at different temporal and spatial scales (Valbuena et al., 2008; Rindfuss et al., 2004; Lambin and Geist, 2001). Environmental policy and management demand the integration of cross-disciplinary knowledge of socio-ecological processes, where agent-based model (ABM) has become a popular investigatory technique to simulate LUCC dynamics (Matthews et al., 2007; Verburg et al., 2004; Parker et al., 2003). ABM use is driven by increasing demand from decision makers to provide support for understanding the potential implications of decisions in complex situations (Smajgl et al., 2011). ABM explicitly represent human decision making processes by means of agents, presented as autonomous computer entities interacting directly with themselves and the environment, in order to achieve goals (Naivinit et al., 2010). Agents in the model can represent individuals or group of people and also can be designed with heterogeneous, autonomous and dynamic characteristics. For this inherent possibilities, ABM explicitly deal with the diversity of land-use decisions (Valbuena et al., 2008). Simulation models of LUCC in space and time can inform policy setting and decision making processes on the use and management of land resources (Le et al., 2008). By mimicking the causal mechanisms and feedback loops of LUCC, these models can be learning tools for understanding the dynamics and driving forces of the land-use system and show how landowners’ choices might affect the direction the future may take (Verburg, 2006). To successfully achieve these ambitions, environmental models and software need to be scientifically and technically sound, reliable, usable and cost effective (Mcintosh et al., 2011). However, LUCC modeling involves the complexity of human drivers and natural constraints, since land-use change emerges from the interactions of these components (Le et al., 2008). To overcome limitations of classical mathematical models, in the Computer Science field, the sub-area of Artificial Intelligence (AI), considering AI modern approach of Russell and Norvig (2009) where agents are the focus, specifically the Multi-Agent Systems (MAS) approach, as studied by Wooldridge (2009), emerges as a candidate technique for solving this class of problems. MAS can be defined as a collection of autonomous entities, called agents, interacting with each other and with their environment (Ferber, 1999). Contrary to traditional modeling techniques, MAS are not expressed in terms of variables,

31

functions and equations, but in terms of agents, objects and environments. They provide the software resources necessary to implement in software ABM descriptions, capturing emergent phenomena resulting from the interactions of individual entities. This paper presents the use of an ABM and a MAS methodology to characterize LUCC dynamics that can be used in the development of different simulation frameworks. The methodology validation was done through the development of MASE, a Multi-Agent System for Environmental LUCC simulation. MASE aims to assist analyzing LUCC dynamics using technical information to aid the decision making process. MASE was used to simulate a LUCC model of the Brazilian Cerrado. The rest of the paper presents: in Section 2, the empirical characterization, the conceptual structure and MASE overview; in Section 3, the case study with simulations and results; in Section 4, a short discussion of the adopted method and MASE results; and finally, in Section 5, the conclusions and future works. 2. Methods In this section, we describe the multi-agent model system used to characterize LUCC dynamics, focusing the parameterization process to help in the development of simulation frameworks. We present the validation of the model defined with the development of MASE, including the architecture and implementation aspects of the application. 2.1. Empirical characterization The purpose of a LUCC simulation framework based in ABM is to represent and simulate the human decision-making processes through the use of agents behavior, the physical environment and the understanding of how the policy decisions affect the dynamics of LUCC. By observing how these variables interact, it is possible to explore causes, consequences and formulate useful scenarios during the decisionmaking process and planning. The results should not be used predictively, but as a tool to develop sustainable strategies for land use. The implementation of human decision-making processes is the main strength of ABM, but the behavioral response functions that represent these processes require knowledge support from qualitative and quantitative empirical sources. Unfortunately, there are no standard approaches to documenting empirical support that underlies modeling and design decisions in agent-based models (Robinson et al., 2007; Berger and Schreinemachers, 2006). Thus, in this work we have used a set of specific methods, commonly used in practice, based on the methodological parameterization process for human behavior in ABM proposed by Smajgl et al. (2011). According to Smajgl et al. (2011) there are two fundamental steps in which empirical data are required: the development of behavioral categories and the scaling to the whole population of agents. The challenge involves the representation of larger populations, identifying suitable empirical methods for eliciting behavioral data, where the sample-based data need to be translated into behavioral representations for the whole population (up-scaling). For this challenge, we have used a proportional up-scaling approach assuming that the sample of agents represent the whole population. Smajgl et al. (2011) build a comprehensive framework with the set of the most commonly used methods for the challenge of parameterize behavioral responses of humans empirically: (i) expert knowledge, (ii) participant observation, (iii) social Table 1 IBGE agricultural census of 2006. Seq.

Producer legal status

Amount

Units or hectares

1

Individual owner e No. establishments Individual owner e Area establishments Condominium, trust or company e No. Condominium, trust or company e Area Cooperative e Number establishments Cooperative e Area establishments Corporation or limited liability e No. Corporation or limited liability e Area Public utility institution e No. Public utility institution e Area Government (federal, state, municipal) e No. Government (federal, state, municipal) e Area Other conditions e No. Other conditions e Area Total No. establishments Total Area establishments

3407 182,031 426 33,564 3 1092 93 21,432 5 1087 7 11,536 14 578 3955 251,320

Units Hectares Units Hectares Units Hectares Units Hectares Units Hectares Units Hectares Units Hectares Units Hectares

2 3 4 5 6 7

32

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46 Sample Population

Agent attributes M2

Theory

Model Formulation M1

M4a Agent attributes

Agent Classes & Structure of agent behaviour class

Agent types M5 M4b M3

Behavioural parameters

Behavioural components

Fig. 1. MASE agents parameterization defined according to Smajgl et al. (2011). surveys, (iv) interviews, (v) census data, (vi) field or lab experiments, (vii) role-playing games, (viii) cluster analysis, (ix) dasymetric mapping, and (x) Monte Carlo method. Relating to this challenge, in our work we have used methods (i), (v) and (vi):  The expert knowledge (EK) covering a variety of formal and informal methods for capturing the understanding of experts in the environmental decisionmaking processes from the Brazilian Institute of Environment and Renewable Natural Resources (IBAMA) and the Institute of Biology of the University of Brasília. According to Cooke and Goossens (2008), the applications of structured expert judgment knowledge can also be employed to quantify uncertainties associated with experts judgment using the TU Delft model (The Delft University of Technology);  Aggregated census data (CD), including GIS data, to summarize responses for groups of agents within enumeration districts (Rees et al., 2002). For this we have used the agriculture census of 2006 provided by the Brazilian Institute of Geography and Statistics (IBGE), as presented in Table 1 (IBGE, 2006);  Lab experiments (LE) designed to observe how the change in independent variables affects a dependent variable, for example, changing the number of agents in the simulation grid, improving the degree of land exploration by different types of agents, and other experiments to establish a more realistic behavior to be observed in the model (Patzer, 1996). Smajgl et al. (2011) cites that ABM requires the systematic representation of three main phenomena: agents, their social networks and the agent environment, but authors discuss only the parameterization of human agents, leaving the development of empirical frameworks for social networks and the agent environment to other researchers. In order to achieve agent parameterization authors identified six different methodological steps, each with multiple empirical methods available. Fig. 1 presents the agent parameterization used in our work according to Smajgl et al. (2011) framework for ABM parameterization. The six steps and the empirical methods used for MASE parameterization is not all inclusive, but we posit that it is complete enough to build an adequate framework according to Smajgl et al. (2011) recommendation for behavioral responses of humans empirically, as presented in Fig. 1 and summarized in Table 2:  Method 1 (M1) e identify different classes of agents, assuming that agents that share the same sequence of actions, or represent the same modeled behavior, are defined in an agent class. The agent classes were defined based on expert knowledge, laboratory experiments and census data. Within each agent class, agents are grouped into agent types (i.e., types are the equivalent of subclasses), according to the classification of Russell and Norvig (2009): simple reflex agents, model-based reflex agents, goal-based agents, and utility-based agents. We have defined two different groups of simple reflex agents (subclasses): (i) service agents, to configure the simulation grid, the MATLAB software1 and window agents; and (ii) transformation agent (TA), that act over the simulation space, composed by farmer, rancher, urbanization and conservative agents. There are also goal-based agents: GRID manager, spatial manager and transformer manager agents. Associated to the spatial manager agent there is the cell agent, a simple reflex agent to deal with the state of the grid cells (colors). For more details see Section 2.3.1, Table 3 and Fig. 12.  Method 2 (M2) e specify values for agents attributes, considering that the parameters governing the magnitude of the actions and the degree to which various sources of information influence those actions can vary. Agent types can be developed from agent attributes or from data on agent behaviors. We considered that data was elicited to measure real-world attributes and

1

http://www.mathworks.com/products/matlab/.









behavioral responses. Behavioral responses were observed relationships between qualitative contextual variables based on expert knowledge, and quantitative contextual variables based on census data from IBGE. For example, the sub-class of TA has different behavioral responses distinguished by attributes like the default tax of land use, which is different for the farmers and urban developers. Method 3 (M3) e specify parameters for the agents behavioral rules. Agents with the same behavior and similar or identical parameters of behavior belong to the same agent type or sub-classes. The sub-class of transformation management agents coordinates the TA behavior (farmers, ranchers, urban and conservative) using a finite state machine (FSM) to define the exploration behavior and trade between the types of land use (i.e., planting or breeding). The FSM used in our work (Fig. 2) was defined based in the work of Smith et al. (1998) as presented in Fig. 11, including land use classes and process transformation in a flow chart, which was based on expert knowledge. Method 4a (M4a) e as noticed agent types can be defined from agent attributes with different methods. In our work we have combined a mix of correlation and expert knowledge with expert knowledge individually to develop agent types from behavioral data. To define agent types of TA we have used the Federal District Spatial Plan (Plano de Ordenamento Territorial do Distrito Federal e PDOT), that sets out the specific areas that can accommodate agriculture or urbanization for the agents behavior. Method 4b (M4b) e an alternative method to define agent types is participant observation from observed behavioral responses, which was done also by expert knowledge in the Federal District of Brazil. Method 5 (M5) e after defining agent types, their attributes and characterize behavioral response functions, agents in the whole population need to be assigned to the appropriate agent type. For this method we have used a upscaling approach as the sub-population and the population are identical as presented in Fig. 1. Since we have used IBGE agriculture census data of the Federal District (Table 1), there is some assurance that the sample used in our model is representative of the population, so the scaling is proportionally.

Smajgl et al. (2011) framework for ABM parameterization is cited to provide an initial structure that can be used within rigorous discussions of parameterization methods for human behavior in ABM. We consider that with the use of methods M1eM5 to ABM parameterization, it was possible to have a more robust empirical model development, where the sequence characterization facilitates the description and documentation of a ABM applications. Also we expect that the model description developed would be important to allow scientific peers to make good judgments about the value of the MASE model results.

2.2. Conceptual structure As criticized by Lorek and Sonnenschein (1999), the documentation of individual or agent-based models (I/ABM) lacks of consistency in the text describing such models. In particular, such descriptions were often incomplete not providing enough information to replicate the results, or were difficult to read as it was not presented in a logical order. To address this problem, Grimm et al. (2006) originally developed the Overview, Design concepts and Details (ODD) protocol for documenting I/ABM. Since then, the ODD protocol has been widely adopted in ecological models. According to Grimm et al. (2010), about 70% of articles citing Grimm et al. (2006) are in the ecological modeling area. Similar proposals for documentation standards to agent-based social simulations are done by Triebig and Klégl (2009), Janssen et al. (2008) and Richiardi et al. (2006). MASE is going to be presented using the ODD sequence. According to the ODD protocol, the I/ABM information follows a structure sequence of seven elements grouped into three blocks (i.e., the ODD sequence) as presented in Sections 2.2.1e2.2.3.

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

33

Table 2 Overview of the parameterization of behavioral traits of human agents. Methods M1

M2

M3

M4a

M4b

M5

Expert knowledge Aggregated census data Lab experiments

Expert knowledge Aggregated census data

Expert knowledge

Expert knowledge

Expert knowledge

Aggregated census data

Table 3 MASE agent classes. Agent class

Type

No. of instances

Functions

GRID Manager (GRIDM)

Goal-based

1

Spatial Manager (SM)

Goal-based

1

Transformer Manager (TM)

Goal-based

1

Cell Agent (CA) Transformer Agent (TA)

Reflexive state record Reflexive state record

User defined User defined

Promote interface parameterizations defined by users Manage start, pause and end of agents Receive agents state for the visualization Promote agents state visualization for the user Instantiate the amount of cells to simulate Get orders from GRIDM and replicates to cells Receive the states of cells and replicates to GRIDM Instantiates number of TA for simulation Get orders from GRIDM and replicates for TA Receive TA states and replicates to GRIDM Perceptions and actions defined in Table 4 Perceptions and actions defined in Table 4

Native Vegetation

Conservative (natural ecosystems protection)

Farmers (agriculture)

Ranchers (cattle raising)

Urban (urbanization)

Fig. 2. FSM deterministic diagram.

2.2.1. Overview According to the ODD sequence, the context and general information is provided first e overview, followed by more strategic considerations e design concepts, and finally more technical details. The purpose of the overview block is to give a very quickly idea of the model’s focus, resolution and complexity. The overview consists of three elements: purpose, state variables and scales, process overview and scheduling. 2.2.1.1. Purpose. MASE purpose is to understand how the behavior of human agents, the physical environment and policy decisions affect the dynamics of use and land cover. The observation of the interaction of these variables will lead to the exploration of different scenarios, which are useful for the planning and development of sustainable strategies for land use. 2.2.1.2. State variables and scales. The state variables and scales is concerned with the following questions: What is the structure of the model system? What kind of low-level entities (e.g., individuals, habitat units) are described in the model? How are they described? What hierarchical levels exist? How are the abiotic and biotic environments described? What is the temporal and spatial resolution and extent of the model system? Considering MASE model, it addresses the issues of LUCC using agents of transformation (representing the socioeconomic changes) and agents of space (representing the changes in the physical environment). MASE model allow the use of policies to land use. In the case study we have adopted the cited PDOT to illustrate where agents interactions to land use will be influenced by rules associated to their specific perceptions and actions over time and space. The physical environment is spatially represented by agents, named agents of space. The agents of space represent the whole simulation grid, divided into different cells, each one represented by an agent. The physical state of the agents corresponds to the set of real spatial data including: (i) water courses (rivers); (ii) water bodies (lakes); (iii) buildings; (iv) highways; (v) railways; (vi) streets; and (vii) Protected Areas (PA). A sub model of land use transitions allows each cell to change its attributes over time, resulting from the action of agents of transformation in every cell and in its neighborhood. Thus, the pattern of landscape change is a result of the composition of the ecological changes that occur in the range of pixels. The human factor is considered in terms of TA: farmers, ranchers, urban developers and conservative agents. A TA is the minimum unit for the representation of

human variables and they are divided in classes and types (sub-classes), according to the description in Section 2.1. The decision model comes from the perceptions of the state of the physical environment and the behavior defined for each type of agent. The decision model is influenced by policies that are being tested in the model. The decision model is probabilistic and universal for every TA. However, decisions are taken in a particular manner for each class or sub-classes of agents, their behavior and state of the environment. The interactions between the landscape and humans come from the cycles of perception and action between the agents of space and TA. The perception corresponds to the physical state of the cell and its vicinity, as well as the benefits and allowances provided by policies related to land use. The action is the effect of the decision of the transforming agent in changing the environment for land use. Fig. 3 presents the different interactions according to Haggith et al. (2003) and Le et al. (2008). The temporal and spatial scale of the model are configurable. The model structure is generic to allow the exploration of space locally, regionally or globally, within the appropriate time to the problem. In the case study, the model operates on a local scale of regions, with areas of the Federal District of Brazil, chosen to illustrate the simulations since it is the single unit in Brazil covered by the Cerrado (Section 3.1). We choose to adapt the land use model defined by Smith et al. (1998), presented in Fig. 11, from the more usual land use Cerrado classes as shown in Fig. 4. A cell in the grid represents one hectare (ha). Each time-step equals one week of chronological time. The simulations were carried out involving the period 2002e 2008, requiring 365 time-steps for the seven years. We have used two time instants (2002 and 2008) with satellite images of LANDSAT ETM, classified by the PROBIO software (IBAMA/Environment Ministry (MMA) or Remote Sensing Center (CSR)-IBAMA/MMA),2 as presented in Fig. 13. 2.2.1.3. Process overview and scheduling. The ODD sequence element is concerned with the process overview: a verbal, conceptual description of each process and its effects is sufficient, since the idea is to give a concise overview. In the implementing of MASE model, the interactions succeed every time-step. The landscape of the initial grid is given by images obtained from GIS raster files. The

2

http://siscom.ibama.gov.br/monitorabiomas/metodologia.htm.

34

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

Land Use Policies Performace indicator

Policy Intervention

Performace indicator

Changes in the physical environment

Changes in Socioeconomic Patterns and processes Interactions

Transformer Agents

Land use Land cover

Spatial Agents

Typology (attributes and behaviours) Decision making

Physical state perception

(GIS data)

action

LUCC submodels

Human System

Landscape System

Fig. 3. Diagram of the conceptual structure.

initial population of agents is generated by the SMA, described in Section 2.3. Users can set the initial population of each type of agent, as well as the policies and parameters that will be explored in the model. Initially there is the preprocessing of the sub model of LUCC and the GIS information that define the physical environment (proximal variables). In each step, it is updated the state of the environment, followed by the TA decision-making. This ODD sequence is also concerned with the scheduling, that deals with the order of the processes and the order in which the state variables are updated. As mentioned in Section 2.2.1 e State variables and scales, the time in MASE is modeled by discrete intervals of time (time-steps). All activities take place in an atomic fashion: only after all actions have been deliberate and executed by agents occurs the transition from time-step and execution of group activities. The variables are updated synchronously in any order of execution within a single time-step. 2.2.2. Design concepts The design concepts provide a common framework for designing and communicating I/ABM. The purpose is to link model design to general concepts identified in the field of Complex Adaptive Systems (Grimm and Railsback, 2005; Railsback, 2001). These concepts include questions about emergence, the type interactions among individuals, whether individuals consider predictions about future conditions, or why and how stochasticity is considered. Here we provide a short checklist of five items that apply to MASE model description as suggested by the ODD protocol Grimm et al. (2006):  emergence e the human agents population dynamics emerge from the behavior of each agent at the land use decisions characteristics that characterize each agent class. However, the life cycle of each agent, as well as their behavior is completely described by empirical rules. The landscape changes in the system emerge solely from the conversion process and land-use change caused by the TA.  sensing e each TA belongs to a specific class type and have a corresponding capacity conversion/exploitation of the environment, so that they apply, for example, specific probabilities depending on the environment condition in which they are inserted. The cells also have internal state that controls how much each piece of land can be exploited or preserved.  fitness (goal-based) e the achievement of objectives is explicitly modeled in the TA decision making process of land use. Each agent checks what area will give a greater amount of return by land exploitation using probability maps, and guides his actions in those directions, so that agents will get the best return from the land use.  interaction e agents interact directly and indirectly. Indirect interactions involve the fact that changes in land use caused by an agent may lead to changes in the decisions of other agents in the next time-step. Direct interaction happens in the race for space, when two agents have the same exploration probability about the same cell space. In this case, rules of succession between different agent classes (i.e., farmers, ranchers) are applied or an agent is chosen randomly among agents of the same class (i.e., farmers).  observation e includes maps of successive coverage and land use for each time-step.

2.2.3. Details According to the ODD sequence, the more technical details include three elements: initialization, input, sub models, to present the details that were omitted in the overview. The sub models implementing the model’s processes are described in detail and the information required to completely re-implement the model and run the baseline simulations is provided in this block. 2.2.3.1. Initialization. The initial landscape model used in MASE simulation uses images (from the classified images of raster-GIS programs obtained by satellite monitoring LANDSAT ETM, as cited in Section 2.2.1), state variables and scales. The physical environment state is also imported into a set of images that correspond to the physical environment interest variables, or proximal variables. In our study case, the time-step t ¼ 0 starts with the map of 2002. TA are initially created with a different number of farmers and ranchers. These agents are allocated spatially in places of greater likelihood of exploitation. Thus, the initial state is always the same in all simulations. The number of TA was parameterized according to the Federal District population provided by the IBGE rural census. 2.2.3.2. Input. The dynamics of many I/ABM models are driven by some environmental conditions which change over space and time. The environmental conditions are input to MASE model, that are imposed dynamics of certain state variables. The model output gives the response of the model to the input. Here we present the input data used in MASE, detailing how they were generated and also how they can be obtained in order to achieve reproducibility.3 The input variables to the simulation process are:  Images of the study area classified for time-step t ¼ 0;  Set of images for each variable proximal set;  Set of policy information to be exploited (the parameters of spatialized images PDOT). In the study case, we used as input two images of raster-GIS obtained by satellite monitoring LANDSAT ETM (IBAMA/MMA) of the Federal District region of Brazil (2002 and 2008), as presented in Fig. 13. Note that the green color is the remaining area, the yellow the anthropic activity area and the blue the water courses (Fig. 13 presents colors in the web version). 2.2.3.3. Proximal variables. We consider that the exploitation of land use depends on several social, economic, geographic and even politic aspects. Many of these factors are difficult to be modeled and predicted; however, certain information from the environment can be taken into account to create land use simulation, that we refer as proximal variables. The proximal variables considered in this study are: railways, highways, rivers, lakes, streets, buildings and PA. For each type of land use stimulus relationships were

3 Online archives of the used input maps from the Federal District of Brazil can be downloaded from the site https://svn.infoknow.comnet.unb.br/mase/model/inputs.

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

T ¼

5 X

Mi  pi

35

(2)

i¼0

Fig. 10 shows an example of how T and Taux can be used to determine the movement of a transformer agent. Take Fig. 10(a) as a reference, where the circle indicates a transformer agent (Fig. 10 presents colors in the web version), consider the yellow areas (in web version) as being or had been explored areas and the green areas, the ones that could potentially be explored. Once there is the need for the agent to move, the Taux map, presented in Fig. 10(b) is consulted and the neighboring cells marked with the logical value 1 are identified as candidates, while the cells marked with logic value 0 (zero) are discarded for exploration. Based on the T and Taux, the new destination for the transformer agent is determined in a probabilistic manner. Fig. 10(c) shows the superimposition of T and Taux inside the agent’s moving range. The first step to determine the agent’s next location is to build a vector V using only the values of T that fall inside its neighboring cells. In Fig. 10(c) example, these values are sequentially scanned from left to right, top to bottom, as presented in vector 3.

Crop 24%

Native Vegetation 38%

Pasture 20%

Urban Area 17%

Reforestation

V½i ¼ ½ 96 21 98 98 110 98 164 74 53 110 175 173 148 148 164 :

1% Fig. 4. Cerrado LUCC model.

 Vnorm i ¼ created, discouraging, or neutral for each of the proximal variables. It is considered that rail, road, river, lake, street and building stimulate anthropic use of areas in their vicinity. In other words we understand that areas close to a highway, for example, are more likely to be occupied than areas farther from roads. To model the ability to stimulate the occupation of each of these variables, the images that represent them were filtered using an isotropic Gaussian filter (circularly symmetric) bi dimensional h(x, y), defined by Equation (1), where hsize and s represent the filter size and the standard deviation, respectively, considering they are parameters that can be adjusted by the system user.

hðx; yÞ ¼

1 x2 þy2 2 e 2s 2ps2

(1)

Fig. 5 presents the Gaussian filter used in generating the images in Fig. 6 (hsize ¼ 11 and s ¼ 2). Fig. 7 shows a magnified portion of a proximal variable. Note that the lighter regions in the images represent the places where there is a greater probability of land occupation for anthropic use. The darker areas, which are surrounded by lighter pixels indicate regions that cannot be used during the simulation process. The darker regions represent areas that are already occupied (i.e., roads, buildings, or areas of the water courses). The variable PA, shown in Fig. 8, has the effect of inhibiting the use of human activity. In that example, the blank region represents the PA that cannot be occupied. After the filtered images are generated, they are used to compute a weighted average, according to Equation (2); where T is the resulting map, Mi represents the filtered image related to the proximal variable i and pi indicates how the variable i will contribute in the generation of the map T. All pi weights were considered equal to 1/6. In addition, there has been generated an assistant map Taux showing the disabled regions for any use, since they are PA or areas that are already occupied. Fig. 9 shows the T and Taux maps for the Federal District region of Brazil.

0.05

h(x,y)

0.04 0.03 0.02 0.01 0 5 5 y

0

0 −5

−5

(3)

Then, the values are normalized according to Equation (4), where n is the number of elements in V[i].

x

Fig. 5. Gaussian filter used to generate images of Fig. 6 (for the color figure please refer the online version).

V½i ; n1 X V½j

i ¼ 0.n  1;

(4)

j¼0

The next step is to calculate a C[i] function according to Equation (5). i X  C i ¼ Vnorm ½j;

i ¼ 0.n  1:

(5)

j¼0

Finally, a pseudorandom number N is drawn from the standard uniform distribution on the open interval (0, 1) and our algorithm searches for the first value of i where C[i] is greater then N. The retrieved index represents the location to where the agent will move. Fig. 10(e) and (f) show Vnorm[i] and C[i] for Equation (3). If, for instance, N ¼ 0.6, the value of i would be 10 (C[10] ¼ 0.63). The final agent’s movement is shown in Fig. 10(d). 2.2.3.4. Sub models. Some model features presented are essential to understand the architectural decisions implemented in the prototype. The FSM, presented in Fig. 2, is used in the management of the transition states of the grid cells. It is necessary to emphasize that the focus of the state machine used in this work is the transition between the native vegetation for agriculture and other exploration. The approach adopted for the SMA prototype, using the proximal variables to determine the movement patterns of the agents and the patterns of exploitation of the environment are all represented by cells in the grid. 2.3. MASE software Considering MAS methodologies, we adopted TROPOS (Dam and Winiko, 2003), which allows to incorporate system requirements to the ODD protocol for documenting ABM (see Section 2.2). The TROPOS diagrams are not included in this paper, since the development of MASE is not the focus, but the methodological aspects involved in the definition of the model that can be replicated in other LUCC simulation packages. One of the important definitions during a MAS project is the environment characterization. In this project we have used the environment characteristics proposed by Russell and Norvig (2009) and Wooldridge (2009). Although our approach deals with real environment (i.e., the geographical area over which the changes occur in coverage and use of land), the simulation grid by the agents perceptions is partially observable, stochastic, sequential, dynamic, continuous, multiagent and competitive, where several restrictions are necessary to be made in order to deal with computable model simulations. To the case study developed in this research work the environment has the following characteristics:  partially observable e each agent has a restricted field of perception related to the grid cells neighborhood (i.e., adjacent cells);  deterministic e the next state of the environment is determined by current state and the actions taken by the agent;  episodic e no time is treated continuously, there are atomic time-steps that are considered for the simulation execution;  static e the environment does not change while an agent is acting;  discrete e the possible transitions are defined by the FSM (Fig. 2);  multi-agent e a set of agents with different roles and behaviors are used in the system;  competitive e agents have interests that are competing, while the grid of resources is limited and has to be shared by agents.

36

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

Fig. 6. Proximal variables: (a) railways; (b) highways; (c) rivers; (d) lakes; (e) streets; and (f) buildings.

2.3.1. Architectural aspects For the definition of MASE architecture different classes and sub-classes of agents where defined as presented in Section 2.1. Each agent class is related to the definition of entities responsible for specific decision making, execution of actions, perceiving the environment and the execution of the time-steps in the simulation. For the first version of MASE prototype five classes of agents were created as shown in Table 3. Note that these classes can be expanded, or new instances of these classes can be created, as the modeling of different behaviors of agents is found necessary. The good aspect of a SMA project is the fact that we adopted an incremental approach, like other software methodologies, what makes easy the development of new features in the system. Linked to each class of agents there are images with the spatial domain to the native class, for example, the map of urban areas are associated with urban agents, thus forming a specific level of simulation. Fig. 12 presents the defined architecture for the first version of MASE prototype system. According to the classification of agent types defined by Russell and Norvig (2009), and cited in Section 2.1, we are going to detail two reflexive agents: cell and TA. These agents play central roles in the simulation, according to the details of perceptions and actions presented in Table 4. To facilitate the coordination of agents, cooperation and conflict resolution, we opted for a hierarchical organization, using different layers. In these layers, agents of higher levels have greater control over agents at lower levels, as shown in Fig. 12. Note that a layer has been designed for the user interface, in order to allow the configuration of the simulation model and also to present the simulation results. All user-supplied configuration settings are translated into the Extensible Markup Language (XML)4 files and loaded into MASE, in order to start the action of agents, linking them to the other simulation layers. Through the user interface layer, user can configure which agents will be part of the simulation, adding behaviors to agents from a library available. User can also associate specific FSM to agents in order to define their behaviors. Moreover,

4

http://www.w3.org/XML/.

additional rules as guidelines for global simulation are available at this stage. The flexibility of configuration extends to the creation of the number of agents that compose each class, which interact in the space defined by the user simulation layer of each respective class. The control layer uses the previously defined rules and settings to produce the simulation, where the logical simulation mechanisms are defined. At this layer, there is some layers overlapping, since the storage of the intermediate results is shared. The physical layer is responsible for loading the actual images and provides image attributes to agents in synchronized approach in order to guarantee threads safety and avoid image inconsistence. 2.3.2. Implementation aspects The development of MASE prototype was performed using the Java Agent Development Framework (JADE), version 4.0 (release date 04/20/2010), which implements the architecture proposed by the Foundation for Intelligent Physical Agents (FIPA) for a definition of SMA. FIPA is an organization connected to the Institute of Electrical and Electronics Engineers (IEEE), responsible for specifying standards for the development of technologies based on intelligent agents. The JADE framework, developed by Telecom Italia LAB (TILAB), is a middleware for developing and executing applications based on intelligent agents, developed in Java. JADE provides many resources such as white and yellow pages, communication and interaction protocols for agents, and tools for agent mobility in distributed platforms (Bellifemine et al., 2007). For the calculation of the influence maps we used MATLAB, version 7.10.0.499 (release R2010a), exported by the matlabcontrol library, which provides an interface and allows JADE agents to interact with systems developed with MATLAB, as presented in Fig. 12 (image processing module). Since we understand that the use of MATLAB could impose a restriction in the use of MASE, we developed the same routines using the open source ImageJ Lybrary (Abramoff et al., 2004). Thus, the references to MATLAB can be shifted by the ImageJ Lybrary without lose of functionalities. Although not part of the prototype first version scope, we plan to improve the actual implementation allowing it to perform in a distributed infrastructure using a computational grid or a cloud infrastructure. To allow MASE

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

37

3. Case study: the Cerrado

Fig. 7. Magnified portion of a proximal variable.

configuration template we choose to use XML as a language to allow future integration to other Web semantic based applications. The prototype implementation and simulations were performed locally, but we intend to distribute MASE into different JADE containers to increase computational performance and possibly integrate in a Web platform. For the tests of the prototype, we used the Computer Networks Lab (COMNET) infrastructure of the Computer Science Department at the University of Brasília, using the following machine: Intel Core i5 2.27 GHz 4 GB RAM, 64-bit operating system. Since we are using JADE framework we adopted the FIPA specifications for the Agent Communication Language (ACL) to provide the exchange of messages among agents. Considering the agents behavior aspects, MASE implements simple approaches, such as OneShotBehaviour and CyclicBehaviour, and composed ones, such as SequentialBehaviour and FSM-Behaviour. Related to the behavior execution aspects, MASE implements a proper scheduling procedure through a distributed approach, where the activation of agents pass from waiting to running states of their life cycle not in a preemptive way. The procedure was defined to allow as an attribute the defined number of simultaneous running TA agents. In our experiments we have used five TA agents as default to be active. We also put in waiting state the other TA agents until one of the five agent finishes its execution, in order to rationalize resources. Note that the other control layer agents are active running in parallel way by the JADE preemptive control. Considering the interaction protocol MASE implements a hybrid approach using cooperative and competitive means of relations. Agents interact through ACL protocol of communication and act in a common grid environment. Since MASE uses a hierarchical organization with different layers, the GRIDM manager creates the TM and SM managers as soon as the pre-processing module finishes, while the TM and the SM managers creates the TA and the CA agents respectively by user-supplied configuration settings. The cooperative aspects are presented among the GRIDM, SM and TM managers relating to the CA and TA agents since they control their activation and specific behaviors, all registered in JADE Agent Management System (AMS). For example, when a simulation step begins the GRIDM inform the start event to the SM and TM managers that inform their related agents CA and TA. During the pre-processing phase MASE registers in the JADE Directory Facilitor (DF) the services, such as the matlabcontrol agent for image processing, the config agent for the user-supplied configuration settings and the grid service that provides the access for the physical layer. The competitive aspects are presented by the TA self-interested agents since they compete for grid resources during the execution of the simulation process.

Fig. 8. Protected areas (strict protection only): inhibit the anthropic use.

The Federal District Cerrado has 68.11% of 5789 km2 of native vegetation cleared and 90% covered by PA of strict protection or sustainable use, corresponding to the IeIII and IVeVI IUCN management categories (IUCN Commission On National Parks And Protected Areas and World Conservation Monitoring Centre, 1994). It is the major Brazilian savanna-like ecosystem and the second largest in Brazil. The Cerrado is characterized as a neotropical savanna with marked climatic seasonality and strong dynamics of human occupation. The region is located between the equatorial zone and 23 south latitude bordered by the Amazon forest to the north, by the Atlantic forest to the south and southeast, and by the Caatinga (deciduous xerophytic vegetation) of the semiarid region to the northeast (Klink and Moreira, 2002). The Cerrado is unique in the sense it serves as corridor for species inhabiting neighboring biomes. The attention of conservationists has focused on rainforest (Amazon, Atlantic), but the Cerrado is currently one the most threatened biomes of South America due to the rapid expansion of agriculture (Oliveira, 2002). According to Sano et al. (2008), it is considered one of the last agricultural production frontiers in the world. The main activities in the region focus on soybean production and large-scale grazing. The conversion of savannas for soy production is accelerating to meet the demand for feed rations for cattle and swine feedlots (McAlpine et al., 2009). Soybeans and soy products are amongst the largest of Brazils export commodities, and the Cerrado support the largest cattle herd in the country (Klink and Moreira, 2002). The expansion of these activities is driven by a series of interconnected socioeconomic factors, often encouraged by government policy. Data from IBAMA show a cumulative loss of 47.8% of the Cerrado natural vegetation cover (around three decades). The intensity and speed of the process of destruction require urgent action to slow down landscape fragmentation, loss of biodiversity, biological invasion, soil erosion, water pollution, land degradation, and the heavy use of chemicals (Cavalcanti and Joly, 2002). To reconcile land use and conservation of this biome is the main challenge for the next decade (Klink and Machado, 2005), where understanding the dynamics of these phenomena is essential for an adequate policy of sustainable use of resources. Although the important aspects of the Brazilian Cerrado biome, we may state that other region could be used to illustrate the multi-agent model system presented in this paper. 3.1. Simulations and results Like any application that predicts change between points in time, we used two maps of the Federal District Cerrado to compare from an initial time (2002 e t0) to a subsequent time (2008 e t6). Each time-step in our simulation corresponds to a week of real time, so 365 steps were necessary to cover the length of seven years. Fig. 13 presents both time instants of the region maps, where the yellow part represents the anthropic used land, the green represents the native vegetation and the blue the water courses (Fig. 13 presents colors in the web version). The total area of the Federal District, representing the physical environment with 5789 km2, is divided into cells, that together form the total space of the simulation grid. Every set of four cells represents one hectare, which is occupied by a different agent. The physical state of the cells corresponds to the set of real spatial data including six proximal variables: (i) water courses (rivers); (ii) water bodies (lakes); (iii) buildings; (iv) highways; (v) streets; and (vi) PA. As described in Section 2.2.3, during the initialization phase of MASE each of the proximal variables are filtered using an isotropic Gaussian filter with MATLAB software. After the filtering process

38

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

Fig. 9. T and Taux maps of the Federal District of Brazil.

the proximal variables are recomposed adding values to water courses, water bodies and highways, and removing buildings and PA to form the primary influence matrix. Also the PDOT polygon is considered to generate the final influence matrix, where each pixel of the simulation grid is compared to the influence matrix entries. Where the PDOT indicates prohibited area for use, the value of the influence matrix was altered to 99,999, in order to maintain the land unchanged. But, where the PDOT indicates encouraged area for use the value is increased in 10%. In the simulation process the human factor over the land is represented by two different types of TA: farmers and ranchers. The instantiation of the farmers and ranchers were random in the simulation grid using a dice to decide from zero to the total number of the matrix entries. In the case the drawn cell of the grid is of native vegetation, included in PDOT as prohibited area or it is already occupied by an agent, then the dice is thrown again until a cell can host a TA. Considering the agents movement, a fixed neighborhood routine was used based on the final influence matrix. For each movement the neighborhood values of the matrix were added without considering the negative values and the occupied cells. All the neighborhood values were processed being divided by the total sum of them and multiplied by 100 in order to have the absolute probability percentage. From the second neighborhood cell on the total value is being added to the previous cell total value turning the neighborhood into a cumulative probabilistic percentage space. Whenever the TA has to move a dice is throw from one to 100 to indicate the probabilistic space interval where it can move. In the case all cells are zero, what means the neighborhood area is totally explored or very close to a forbidden area for exploitation, then the transformer manager relocates the agent to a new anthropic area (i.e., a grid cell). Given the simulation assumptions, we are going to present one simulation to illustrate the use of MASE. As simulation parameters, we adopted that the native vegetation have 1500 units of exploitation potential and the anthropic cells 500 units. Whenever an exploitation happens in a cell, the adjacent cells decreases 150 exploitation units, representing the environmental impact caused by the transformation activity. The TA (farmers and ranchers) explore at the rate of 500 exploration units and the number of agents vary considering ten of each type (ten farmers and ten ranchers). Not only the number of farmer agents varies, but also their exploitation potential to demonstrate the difference of agent groups related to LUCC. The number of TA with different exploitation potential was calibrated from the 2006 agricultural census published by IBGE according to data in Table 1. There is also a fixed number of manager agents (grid, spatial and transformer) and cell agents (see Table 3). Each set of five steps were logged in a MATLAB batch file to be processed afterward. After 365 steps, the simulation grid was reset, the influence matrix recalculated and TA also reset, in order to add ten more agents of each type to

repeat the simulation process, which was done up to 150 agents of each type. In order to compare quantitative and qualitative aspects of the multi-agent model proposed we applied two scientifically rigorous statistical techniques of map comparison to land change models. The methods were developed by Pontius et al. (2008) to compare: (i) a reference map of initial time (2002), (ii) a reference map of subsequent time (2008), and (iii) a prediction map of the subsequent time (2008). According to the authors, the three possible two-map comparisons for each application characterize: (i) the dynamics of the landscape, (ii) the behavior of the model, and (iii) the accuracy of the prediction. The authors also propose the threemap comparison for each application to specify the amount of the predictions accuracy that is attributable to land persistence versus land change. According to the method the applications can be summarized and compared using two statistics: the null resolution and the figure of merit. Considering the figure of merit, the more accurate applications are the ones where the amount of observed net change in the reference maps is larger. 3.1.1. Simulation with three possible two-map comparisons The analytic procedure of the three possible two-map comparisons considering the Federal District Cerrado application has basically two categories: anthropic used land (Fig. 13, yellow area, in web version) and the native vegetation (Fig. 13, green area, in web version), since the water courses cannot be explored (blue area, in web version). According to Pontius et al. (2008) method, Figs. 14 and 15 show the three possible two-map overlay with 90 agents of each type of TA (farmers and ranchers), since in previous simulations good prediction results were reached with this amount of agents (see Fig. 17 for the figure of merit percentage). As cited before, the number of TA with different exploitation potential was calibrated from the 2006 agricultural census (Table 1). We defined two different exploitation potential groups: the individual (seq. 1) and group ones (seq. 2e7). The individual TA represents 3407 units (86%) of 182,031 ha (72%) and the group 548 units (14%) of 69,289 ha (28%). As cited, the group formed by the individual TA explores 86% units (3407 of 3955) of the total Cerrado area (251,320 ha), and we defined for them the exploitation rate of 500 units per step. For the group formed by condominium, cooperatives, corporation, public utility, government and other institutions we defined 1500 units of exploration per step of simulation. Fig. 14 illustrates the three possible two-map overlay for the 90 TA, with 50% of agents exploiting at a 500 unit rate and 50% of agents exploiting at a 1500 unit rate. Fig. 14(a) examines the difference between the reference map of initial time (2002) and the reference map of subsequent time (2008), which Fig. 15 quantifies by the length of the top bar labeled Observed. Fig. 14(b) examines the difference between the reference map of initial time (2002) and the prediction map of the subsequent time (2008), which Fig. 15

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

0.12

39

1 0.9

0.1

0.8 0.7

0.08 C[i]

0.06

V

norm

[i]

0.6 0.5 0.4 0.04

0.3 0.2

0.02

0.1 0 0

2

4

6

i

8

10

12

14

0 0

2

4

6

i

8

10

12

14

Fig. 10. TA movement: (a) original cells; (b) Taux map; (c) values of T and Taux inside the agent’s moving range; (d) T map showing the agent’s movement to the next location; (e) Vnorm[i]; and (f) C[i] (for the color figure please refer the online version).

quantifies by the length of the bar labeled 10TA e Predicted to 150TA e Predicted, representing the predicted using 10e150 TA. Fig. 14(c) examines the difference between the reference map of subsequent time (2008) and the prediction map of the subsequent

time (2008), which Fig. 15 quantifies by the length of the bar labeled 90 TA e 5% e Error to 90 TA e 50% e Error, representing the accuracy of the prediction using 90 TA with 5%e50% grouped agents exploiting at a 1500 unit rate.

40

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

PL PLANTATION

EF AG

CROP

EXT

CSV

AP

NATIVE VEGETATION

EXT

CHANGED AREA

DM DM

PLANTED PASTURE

UR

PC

RE

LEGEND UR

PC

EF

UR

URBAN AREA

ABANDONMENT AREA UR UR EF

Transformation Processes Land Use Classes AG agriculture CSV conservation DM deforestation PC Cattle raising UR Urbanization

Fig. 11. Flow chart of land use classes and process transformation.

Fig. 14 allows the reader to assess visually the nature of the observed, predicted and error according to the different gray shades of the legend. For the two-map comparisons maps (a) e Observed and (b) e Predicted, the black pixels show the proximal variables in different layers for rivers, lakes, buildings, highways, streets and PA, since they are not to be considered for the LUCC. The dark gray pixels show the anthropic gain, the medium gray pixels show the anthropic persistence and the light gray pixels show the native vegetation persistence. For the two-map comparisons map (c) e Error, the black pixels show the layers for the proximal variables. The dark gray pixels show the error due to observed native

vegetation which was predicted as anthropic. The medium gray pixels show where the LUCC model predicts anthropic use correctly. The light gray pixels show where the LUCC model predicts native vegetation correctly. The white pixels show the error due to observed anthropic predicted as native vegetation. A series of previous papers describe how to budget the total disagreement between any two maps that share a categorical variable in terms of separable components (Pontius, 2000, 2002; Pontius et al., 2004). According to the authors, the two most important components are quantity disagreement (i.e., net change) and location disagreement (i.e., swap change), which sum to the

Fig. 12. MASE architecture.

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46 Table 4 Perception and action roles of Cell Agent (CA) and Transformer Agent (TA). Agent

Perception

Action

CA

Receive tasks from SM

TA

Receive tasks from TM

Inform state to SM agent Begin land/vegetation recovering or stop it Signal whether or not land is occupied by TA Inform state to TM agent Request position change to TM agent Moving from one cell to another Explore the cell Identify if cell has exhausted its ability to be exploited

total disagreement. While the quantity disagreement derives from differences between the maps in terms of the number of pixels for each category, the location disagreement is the disagreement that could be resolved by rearranging the pixels spatially within one map. If the location disagreement can be resolved by swapping the pixels over small distances, then it budgets the error as “near” location disagreement, otherwise it budgets the error as “far” location disagreement. We adopted in this paper what Pontius et al. (2008); Pontius (2002) suggested for “near” location disagreement: the location disagreement that can be resolved by swapping within 64-row  64-column clusters of pixels of the raw data (raster data). In this paper we did not adopt the coarsening procedure used by Pontius et al. (2008), where the side of each coarse pixel is 64 times larger than the spatial resolution, since all of our maps are in the same resolution. The bars of Fig. 15 are helpful to compare the application results. The graph presents a single technique to show important characteristics of the observed and predicted maps using bars in order to interpret the error bar properly, where the observed bar is the error of a null model that predicts pure persistence (i.e., no change between 2002 and 2008). So if the observed bar is smaller than the error bar, then the null model is more accurate than the LUCC model. Fig. 15 quantifies the differences in each of the two-map overlays for each run of the simulation, with 90 TA with variable exploration rates. According to Pontius et al. (2008), if the model were to predict the observed change perfectly, then Fig. 14(a) and (b) would be identical, and the Observed bar would be identical to the Predicted bar in Fig. 15. The error bars presented in Fig. 15 are much smaller than the observed change bar, meaning that the MASE LUCC model is more accurate than its null model. 3.1.2. Simulation with one possible three-map comparison According to Pontius et al. (2008), there is an additional validation technique that considers the overlay of all three maps: the reference map of initial time (2002), the reference map of subsequent time (2008), and the prediction map of the subsequent time (2008). We consider this three-map comparison very good, since it

41

allows one to distinguish visually the pixels that are correct due to persistence versus the pixels that are correct due to change. Fig. 16 allows the reader to assess visually the nature of the prediction errors using different gray shades. The black pixels show the proximal variables in different layers for rivers, lakes, buildings, highways, streets, railways and PA, since they are not to be considered in the land transformation. The dark gray pixels show the LUCC model predicts change correctly. Medium gray pixels show error where change is observed at locations where the model predicts persistence. Light gray pixels show error where persistence is observed at locations where the model predicts change. White pixels show locations where the LUCC model predicts persistence correctly. If there are more light gray pixels than medium gray pixels, means that the LUCC model predicts more change than is observed. Also a more trained user can distinguish between “near” and “far” location disagreement, for example if there is a medium gray pixel next to a light gray pixel means the error is of near location. Figs. 17 and 18 presents a summary of the Federal District application according to the logic of the legend for Fig. 16 as suggested by Pontius et al. (2008). Each bar is a rectangular Venn diagram where the solid cross-hatched central segments represent the intersection of the observed change and the predicted change; the central solid black segment is change that the model predicts correctly. The union of the segments on the left and center portions of each bar represents the area of change according to the reference maps, and the union of the segments on the center and right portions of each bar represents the area of change according to the prediction map. If a prediction were perfect, then its entire bar would have exactly one solid segment, which would have a length equal to both the observed change and the predicted change. The “figure of merit” is a statistical measurement that derives from the information of the bars in Fig. 18. The figure of merit is the ratio of the intersection of the observed change and predicted change to the union of the observed change and predicted change, which can range from 0% (no overlap between observed and predicted change) to 100% (perfect overlap between observed and predicted change, a perfect accurate prediction). Equation (6) defines mathematically the figure of merit, where A is the area of error due to observed change predicted as persistence, B area of correct due to observed change predicted as change, C area of error due to observed change predicted as wrong gaining category, and D area of error due to observed persistence predicted as change.

Figure of merit ¼ B=ðA þ B þ C þ DÞ;

(6)

Figs. 17 and 18 can also be used to show two types of conditional accuracy: producers accuracy and users accuracy, which are expressed as percents in the figures. Equation (7) gives the producers accuracy, which is the proportion of pixels that the model predicts accurately as change, given that the reference maps

Fig. 13. LANDSAT ETM satellite images from two different years of the Federal District of Brazil, classified by PROBIO software (IBAMA/MMA) (for the color figure please refer the online version).

42

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

Fig. 14. The Federal District maps of: a) observed change 2002e2008, b) predicted change 2002e2008, and c) prediction error 2008 for the simulation with behavior diversity among the same transformer class of agents.

indicate observed change. Equation (8) gives the users accuracy, which is the proportion of pixels that the model predicts accurately as change, given that the model predicts change.

Producers accuracy ¼ B=ðA þ B þ CÞ;

(7)

Users accuracy ¼ B=ðB þ C þ DÞ;

(8)

In summary, Figs. 17 and 18 summarize the most important results in a manner that facilitates cross application comparison at the resolution of the raw data. Considering the figure of merit, if it is greater than 50% the amount of correctly predicted change is larger than the sum of the various types of error, what means that the producers and users accuracy is greater than 50% and the application is more accurate than the null model at the resolution of the raw data. The more accurate applications are the ones where the amount of observed net change in the reference maps is larger. Considering the simulations discussed, Fig. 17 shows that with more than 100 TA the figure of merit is greater than 50%. Considering the producers accuracy it is greater than 50% with more than 50 TA. The users accuracy is greater than 50% for all number of TA, what means that the proportion of pixels that the model predicts accurately as change is very good. Fig. 18 shows the figure of merit

greater than 50% with all number of TA, presenting the best result with 90 TA and 50% exploiting at 1500 units rate, achieving 53%, 75% and 65% for the figure of merit, the producers accuracy and the users accuracy, respectively. Simulations with 100 TA and 5% exploiting at 1500 units rate were conducted and presented even better results, achieving 55%, 88% and 60% for the figure of merit, the producers accuracy and the users accuracy, respectively. The reader can obtain color version of the Federal District simulation maps by contacting the first author or visiting https://svn. infoknow.comnet.unb.br/mase/model/results with “anonymous” for user and password. Pontius et al. (2004) also describe an additional statistical method of validation that considers all three maps simultaneously, the technique compares the accuracy of the LUCC model to the accuracy of its null model at multiple resolutions, but we did not use this method since our case study maps are not with multiple resolutions. 4. Discussion In the real world, change in time and space patterns are produced by the interaction of biophysical and socioeconomic processes. Although ABM and MAS have been cited as natural

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

43

Fig. 15. Observed change, predicted change, and prediction error for the 15 runs of the simulation (i.e. 90 TA 50%: Run with 90 TA exploiting at a rate of 1500 units).

metaphors for the analysis of models and theories of interactivity of human societies, the representation of such interactions in a limited computer model is not a trivial task, translating into a dilemma for many environmental research works. We agree that ecological models are adequate to be described by ABM, since they involve the interaction of different agents in the environment, relating socioeconomic issues to biophysical, temporal and spatial aspects, but MAS can effectively provide the software resources to implement those models thought powerful simulation frameworks. In this direction, we consider this papers approach a good place to begin the definition of a simulation framework for a variety of reasons. The presented multi-agent model system has a replicable

Fig. 16. Validation map for the simulation obtained by overlaying the three maps. Simulation with 90 TA, 50% exploiting 500 and 50% exploiting 1500 units rate.

parameterization process that can be useful to the development of simulation frameworks. The methodological approach adopted, which is based in the empirical characterization model of Smajgl et al. (2011) done in a systematic and structured way, what forms a solid backbone together with the conceptual model definition according the ODD agent-based documentation protocol of Grimm et al. (2006). Thus, we consider that this paper presents a generally applicable model useful to the modelers community. During the development of this research several simulations were performed, but in Section 3.1 we present only one for space constraints. The simulation was performed in the machine cited in Section 2.3.2, a Intel Core i5 2.27 GHz 4 GB RAM, 64-bit operating system. The simulation considered variable exploration rates and percentage of group agents varying from 5% to 50% presenting the execution time of 60 min for the 365 steps, corresponding to seven years, from 2002 to 2008 (Figs. 14 and 15). According to the scientifically rigorous method of Pontius et al. (2008), our simulation results indicate the potential of the presented multi-agent model system. Considering the accuracy of the simulations using MASE, the application results were better than the null model, what examines both the behavior of the model and the dynamics of the landscape. In Pontius et al. (2008) we find the method to multiple resolution map comparison to quantify characteristics applied for 13 LUCC modeling applications of 9 different popular peer-reviewed LUCC models (Geomod, SLEUTH, Land Transformation ModeleLTM, CLUE, CLUE-S, Land Use Scanner, Environment Explorer, Logistic Regression, SAMBA). The results show that 12 of the 13 applications contain more erroneous pixels than pixels of correctly predicted land change at the fine resolution of the raw data. The models compared by Pontius et al. (2008) used a variety of techniques including statistical regression, cellular automata, machine learning and ABM, but the figure of merit statistical analysis shows that only for one of the 13 applications, the figure of merit was greater than 50%. The multi-agent model system presented in this paper could simulate the observed scenario correctly, being able to overcome the null resolution model using a multi-agent approach with the right calibration (90 TA with 50% exploiting at 1500 units rate).

44

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

Fig. 17. Sources of percent correct and percent error in the validation with variable number of TA exploiting at the same rate (500).

Fig. 18. Sources of percent correct and percent error in the validation with 90 TA, 50% exploiting 500 and 50% exploiting 1500 units rate.

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

MASE was able to use the correct or nearly correct net quantities for the categories in the prediction map. The simulation presented demonstrates that the variety among the same class of agents presents a more accurate model than static agents with fixed behavior. Furthermore, it is reasonable to think that models to be able to simulate accurately various possible dynamics according to numerous alternative scenarios would need the models underlying mechanisms to be validated under a wide range of circumstances. This paper does not compare those underlying mechanisms according to numerous alternative scenarios, since the focus of the paper was not this. But, we consider the results of the defined model system promising to be considered by change modelers since they can benefit from the reported experience. Agarwal et al. (2002) developed a detailed study of environmental tools by analyzing tendencies and methodologies. As a result they present a review and assessment of land-use change models considering three dimensions: dynamics of space, time and human choice. According to this study, the ultimate goal for modeling the dynamics between man and the environment involves high complexity in the three dimensions. As a result they define six levels of complexity for human decision-making: (i) without human decision-making, only the biophysical variables in the model; (ii) human decision-making determined by population size or population density; (iii) human decision-making seen as a probability function depending on socioeconomic or biophysical variables beyond population variables, without feedback from the environment to the function that determines the decision; (iv) similar to (iii), but with feedback; (v) a type of agent, whose decisions are modeled in relation to choices made and based on variables that affect other processes and outcomes; and (vi) similar to (v), but with multiple agents, where the model may also be able to make decisions on domains, as the time steps are processed, or about interactions between agents decision makers. Considering Agarwal et al. (2002) study, we believe MASE achieved the sixth level of complexity for modeling the dynamics between man and the environment when considering the three dimensions: dynamics of space, time and human choice. MASE multiple agents are able to make decisions of LUCC exploration considering different domains, including areas that are already occupied (i.e., roads, buildings), areas of water courses or PA. Also agents in MASE maintain their autonomy to move and explore using different exploitation rates in the simulation space keeping interactions between different agent types (individuals) and groups of agents (e.g., condominium). 5. Conclusions Different environmental simulation tools use multiple techniques, such as Markov models, cellular automata and models based in individuals or agents; consider the space in a realistic way representing it implicitly or explicitly; and the adjacent cells are computed with static or variable radius. The LUCC models are normally top-down or bottom-up as cited in Section 1. The approach presented in this paper is a hybrid one with a multi-agent model system to simulate LUCC dynamics, using multiple agents to represent the interaction between different types of agents with autonomy. The space is represented in a realistic way using real LANDSAT ETM satellite images and the simulation grid can be computed with static or variable exploration rates, according to the behavior of individual or multiple agents, achieving the sixth level of complexity for modeling LUCC dynamics (Agarwal et al., 2002). Considering the experimental results presented in Section 3.1, and the figure of merit statistical analysis of simulations, we note that the application results overcome the null model obtaining figure of merit greater than 50% and also good percentage of

45

producers and users accuracy. These results examine both the behavior of the LUCC model proposed and the dynamics of the landscape. Considering Pontius et al. (2008) results, where only one of the 13 LUCC modeling applications with nine different popular peer-reviewed LUCC models presented the figure of merit greater than 50%, we consider the multi-agent model system presented represents an interesting model for LUCC simulation developers. The model is done in a spatially-explicit, integrated and multi-scale manner, being important for the projection of alternative pathways for conducting experiments that test human understanding of key processes in land-use changes. As future work we may cite the increase of agents rationality since this involves the most challenging aspects of the research focusing AI, MAS and ABM. But, we consider that the definition of MASE using MAS methodology would help with of a more comprehensive rational model to agents to represent behavior and social interactions, e.g. beliefedesireeintention model. Also other techniques to improve expressivity of LUCC complex process caused by the interaction between natural and social systems can be used allowing different temporal and spatial scales, as cited in Section 1. This can be achieve for example using fuzzy logics, foundational ontology, qualitative spatial/spatiotemporal reasoning, what would increase complexity of the proposed model demanding to prove the validity for policy makers. Acknowledgments The first author acknowledges the scholarship received from the Brazilian Coordenação de Aperfeiçoamento de Pessoal de Nível Superior (CAPES) e Proc. no BEX 6113/12  5. We are grateful to the Brazilian Institute of Environment and Renewable Natural Resources (IBAMA) for the Federal District LANDSAT ETM satellite images, being specially grateful to Carolina G. Abreu. Also, we would like to thank Alexey Voinov and the anonymous referees for their valuable comments on preliminary version of this article. References Abramoff, M.D., Magelhaes, P.J., Ram, S.J., 2004. Image processing with ImageJ. Biophotonics International 11 (7), 36e42. Agarwal, C., Green, G.M., Grove, J.M., Evans, T.P., Schweik, C.M., 2002. A Review and Assessment of Land-use Change Models: Dynamics of Space, Time, and Human Choice. U.S. Dept. of Agriculture, Forest Service, Northeastern Research Station. Bellifemine, F.L., Caire, G., Greenwood, D., 2007. Developing Multi-Agent Systems with JADE. John Wiley & Sons. Berger, T., Schreinemachers, P., 2006. Creating agents and landscapes for multiagent systems from random samples. Ecology And Society 11 (2), 19. Cavalcanti, R.B., Joly, C.A., 2002. Biodiversity and conservation priorities in the Cerrado region. In: Marquis, R.J., Oliveira, P.S. (Eds.), The Cerrados of Brazil: Ecology and Natural History of a Neotropical Savanna, first ed. Columbia University Press, New York, pp. 351e367 (Chapter 18). Cooke, R.M., Goossens, L.L.H.J., 2008. TU Delft expert judgment data base. Reliability Engineering & System Safety 5 (93), 657e674. Dam, K.H., Winiko, M., 2003. Comparing Agent-oriented Methodologies. Ferber, J., 1999. Multi-Agent System: an Introduction to Distributed Artificial Intelligence. Addison-Wesley Professional. Grimm, V., Berger, U., Bastiansen, F., Eliassen, S., Ginot, V., Giske, J., Goss-Custard, J., Grand, T., Heinz, S.K., Huse, G., 2006. A standard protocol for describing individualbased and agent-based models. Ecological Modelling 198 (1e2), 115e126. Grimm, V., Berger, U., DeAngelis, D.L., Polhill, J.G., Giske, J., Railsback, S.F., 2010. The ODD protocol: a review and first update. Ecological Modelling 221 (23), 2760e2768. Grimm, V., Railsback, S.F., 2005. Individual-based Modeling and Ecology: (Princeton Series in Theoretical and Computational Biology). Princeton University Press. Haggith, M., Muetzelfeldt, R.I., Taylor, J., 2003. Modelling decision-making in rural communities at the forest margin. Small-scale Forest Economics. Management and Policy 2 (2), 241e258. Houet, T., Verburg, P.H., Loveland, T.R., 2010. Monitoring and modelling landscape dynamics. Landscape Ecology 25, 163e167. IBGE, 2004. Mapa de biomas do brasil (map of biomes of brazil). Instituto Brasileiro de Geografia e Estatística e IBGE (brazilian institute of geography and statistics). Escala 1:5.000.000. Available at: http://mapas.ibge.gov.br/ (accessed 25.05.12.).

46

C.G. Ralha et al. / Environmental Modelling & Software 42 (2013) 30e46

IBGE, 2006. Brazilian Census of Agriculture e Brazilian Institute of Geography and Statistics. Instituto Brasileiro de Geografia e Estatística e IBGE. Available at: http://www.ibge.gov.br/home/estatistica/economia/agropecuaria/censoagro/ default.shtm (accessed 31.05.12). IUCN Commission On National Parks And Protected Areas and World Conservation Monitoring Centre, 1994. Guidelines for Protected Area Management Categories. IUCN/UICN. Available at: http://books.google.com.br/books? id¼4gi3vtvRJ9EC (accessed 10.06.12.). Janssen, M.A., Alessa, L.N., Barton, M., Bergin, S., Lee, A., 2008. Towards a community framework for agent-based modelling. Journal of Artificial Societies and Social Simulation 11 (2). Klink, C.A., Machado, R.B., 2005. Conservation of the Brazilian Cerrado. Conservation Biology 19 (3), 707e713. Klink, C.A., Moreira, A.G., 2002. Past and current human occupation, and land use. In: Marquis, R.J., Oliveira, P.S. (Eds.), The Cerrados of Brazil: Ecology and Natural History of a Neotropical Savanna, first ed. Columbia University Press, New York, pp. 69e88 (Chapter 5). Lambin, E.F., Geist, H.J., 2001. Global land-use and land-cover change: what have we learned so far? Global Change Newsletter 46, 27e30. Lambin, E.F. (Ed.), 2006. Land-use and Land-cover Change: Local Processes and Global Impacts, first ed.. In: Geist, H.J. (Ed.), Global Change e The IGBP Series (Closed) Springer. Le, Q.B., Park, S.J., Vlek, P.L.G., Cremers, A.B., 2008. Land-Use Dynamic Simulator (LUDAS): a multi-agent system model for simulating spatio-temporal dynamics of coupled human-landscape system. I. Structure and theoretical specification. Ecological Informatics 3 (2), 135e153. Lorek, H., Sonnenschein, M., 1999. Modelling and simulation software to support individual-based ecological modelling. Ecological Modelling 115 (2e3), 199e216. Matthews, R.B., Gilbert, N.G., Roach, A., Polhill, J.G., Gotts, N.M., 2007. Agent-based land-use models: a review of applications. Landscape Ecology, 1447e1459. McAlpine, C.A., Etter, A., Fearnside, P.M., Seabrook, L., Laurance, W.F., 2009. Increasing world consumption of beef as a driver of regional and global change: a call for policy action based on evidence from Queensland (Australia), Colombia and Brazil. Global Environmental Change 19 (1), 21e33. Mcintosh, B.S., Alexandrov, G., Matthews, K., Mysiak, J., van Ittersum, M., 2011. Environmental Modelling & Software Preface: thematic issue on the assessment and evaluation of environmental models and software. Environmental Modelling and Software 26 (3), 245e246. Myers, N., Mittermeier, R.A., Mittermeier, C.G., da Fonseca, G.A.B., Kent, J., 2000. Biodiversity hotspots for conservation priorities. Nature 403, 853e858. Naivinit, W., Le Page, C., Trébuil, G., Gajaseni, N., 2010. Participatory agent-based modeling and simulation of rice production and labor migrations in Northeast Thailand. Environmental Modelling and Software 25 (11), 1345e1358. Oliveira, R.J.M.P.S., 2002. The Cerrados of Brazil: Ecology and Natural History of a Neotropical Savanna, first ed. Columbia University Press, New York. Parker, D.C., Manson, S.M., Janssen, M.A., Hoffmann, M.J., Deadman, P., 2003. Multiagent systems for the simulation of land-use and land-cover change: a review. Annals of the Association of American Geographers 93 (2), 314e337. Patzer, G.L., 1996. Experiment-research Methodology in Marketing: Types and Applications. Greenwood Publishing Group. Pontius, R.G., 2000. Quantification error versus location error in comparison of categorical maps. Photogrammetric Engineering and Remote Sensing 66 (8), 1011e1016. Pontius, R.G., 2002. Statistical methods to partition effects of quantity and location during comparison of categorical maps at multiple resolutions. Photogrammetric Engineering and Remote Sensing 68 (10), 1041e1049.

Pontius, R.G., Boersma, W., Castella, J.-C., Clarke, K., de Nijs, T., Dietzel, C., Duan, Z., Fotsing, E., Goldstein, N., Kok, K., Koomen, E., Lippitt, C.D., McConnell, W., Sood, A.M., Pijanowski, B., Pithadia, S., Sweeney, S., Trung, T.N., Veldkamp, A.T., Verburg, P.H., 2008. Comparing the input, output, and validation maps for several models of land change. Annals of Regional Science 42, 11e37. Pontius, R.G., Shusas, E., McEachern, M., 2004. Detecting important categorical land changes while accounting for persistence. Agriculture, Ecosystems and Environment 101 (2e3), 251e268. Railsback, S.F., 2001. Concepts from complex adaptive systems as a framework for individual-based modelling. Ecological Modelling 139 (1), 47e62. Rees, P., Martin, D., Williamson, P., 2002. The Census Data System. John Wiley & Sons, Chichester, UK. Rezende, A.V., Walter, B.M.T., Fagg, C.W., Felfili, J.M., da Silva Júnior, M.C., Noqueira, P.E., de Mendonça, R.C., de Filgueira, T.S., 2008. Flora vascular do bioma cerrado: checklist com 12.356 espécies. In: Sano, S.M. (Ed.), Cerrado: ecologia e floraEmbrapa Cerrados, first ed., vol. 2, ISBN 978-85-7383-421-5, pp. 423e1279 (Chapter 15). Richiardi, M., Leombruni, R., Saam, N.J., Sonnessa, M., 2006. A common protocol for agent-based social simulation. Journal of Artificial Societies and Social Simulation 9 (1). Rindfuss, R.R., Walsh, S.J., Turner II, B.L., Fox, J., Mishra, V., 2004. Developing a science of land change: challenges and methodological issues. Proceedings of the National Academy of Sciences of the United States of America 101 (39), 13976e 13981. Robinson, D.T., Brown, D.G., Parker, D.C., Schreinemachers, P., Janssen, M.A., Huigen, M., Wittmer, H., Gotts, N., Promburom, P., Irwin, E., Berger, T., Gatzweiler, F., Barnaud, C., 2007. Comparison of empirical methods for building agent-based models in land use science. Journal of Land Use Science 2 (1), 31e55. Russell, S., Norvig, P., 2009. Artificial Intelligence: a Modern Approach, third ed. Prentice Hall. Sano, E.E., Rosa, R., Brito, J.L.S., Ferreira, L.G., 2008. Mapeamento semidetalhado do uso da terra do bioma cerrado. Pesquisa Agropecuária Brasileira 43 (1), 153e156. Sano, E.E., Rosa, R., Brito, J.L.S., Ferreira, L.G., 2010. Mapeamento do uso do solo e cobertura vegetal e bioma cerrado: ano base 2002. MMA/SBF, Brasília. Schaldach, R., Alcamo, J., Koch, J., Kölking, C., Lapola, D.M., Schüngel, J., Priess, J.A., 2011. An integrated approach to modelling land-use change on continental and global scales. Environmental Modelling and Software 26 (8), 1041e1051. Smajgl, A., Brown, D.G., Valbuena, D., Huigen, M.G.A., 2011. Empirical characterisation of agent behaviours in socio-ecological systems. Environmental Modelling and Software 26, 837e844. Smith, J., Winograd, M., Gallopin, G., Pachico, D., 1998. Dynamics of the agricultural frontier in the Amazon and savannas of Brazil: analyzing the impact of policy and technology. Environmental Modeling and Assessment 3, 31e46. Triebig, C., Klégl, F., 2009. Elements of a documentation framework for agent-based simulation models. Cybernetics and Systems 40, 441e474. Valbuena, D., Verburg, P.H., Bregt, A.K., 2008. A method to define a typology for agent-based analysis in regional land-use research. Agriculture, Ecosystems and Environment 128 (1e2), 27e36. Verburg, P.H., 2006. Simulating feedbacks in land use and land cover change models. Landscape Ecology 21 (8), 1171e1183. Verburg, P.H., Schot, P.P., Dijst, M.J., Veldkamp, A., 2004. Land use change modelling: current practice and research priorities. GeoJournal 61 (4), 309e324. Vitousek, P.M., 1997. Human domination of Earth’s ecosystems. Science 277 (5325), 494e499. Wooldridge, M., 2009. An Introduction to MultiAgent Systems, second ed. John Wiley & Sons.