Abstract

We propose a discrete agent-based model to investigate the migration dynamics of heterogeneous individuals. Compatibility among agents of different types is expressed in terms of homophily parameters capturing the extent to which similar individuals are attracted to, or dissimilar individuals are repelled by, each other. Based on agent-based simulations, we establish the connection between emerging spatiotemporal patterns and the homophily parameters. Key results are presented in a novel phase diagram, which reveals a wide range of spatial patterns including the cell, worm, herd, amoeba, and swarm modes under the dynamic regime and the separation, ghetto, and integration modes under the stationary one. Our model thus provides a generalized framework encompassing both static equilibrium and nonstationary systems to investigate the impact of agent heterogeneity on population dynamics. We demonstrate potential applications of our model to social systems using sexual segregation of ungulate habitats as a case study.

1. Introduction

Spatiotemporal patterns are ubiquitous in natural and social systems. Many of these patterns are the results of bottom-up, self-organizing behavior without the presence of a top-down, central control. Numerous studies have explored whether such complex patterns, either spatial or temporal, can be explained with a simple set of rules defined between interacting individuals, or agents, in the system. Discrete multiagent models are widely used to pursue such a quest in physics, biology, and social sciences. In such a model, the rules of interaction between agents can be defined in unlimited ways, depending on the system of interests.

One approach is to employ a set of rules borrowed from the laws of nature. For example, when biologists try to identify a fission-fusion process, the chemical and physical forces acting on each cell are identified and the equations of motion are applied to describe cell motility in the most precise way [1, 2]. Another example is the application of conventional physics laws, such as Newton’s law of motion, which is relatively straightforward when we deal with particle-like individuals. But when the constituting individuals behave autonomously in a way that cannot be predicted with certainty, the laws of physics alone are inadequate. For example, when simulating the migratory patterns of a school of fish [3, 4], we cannot directly identify the forces that lead to the spontaneous emergence of collective behavior. Each individual has an internal degree of autonomy [5], in contrast to particle-like agents that follows simple, mechanistic rules. Pedestrians [6], bird flocks [7], and schools of fish are examples of autonomous agents that react collectively to environmental changes. Often the external forces driving the behavior of autonomous agents are cultural, economical, psychological, or sociodemographical. In this case physical laws alone are insufficient, and new models need to be developed in order to explain the emerging collective behavior with a simple set of rules. In [8] for example, the flight pattern of a flock of birds is successfully reduced to just three simple rules. Most of these examples dealing with the spatial patterns of natural systems share at least one thing in common, namely, that the rules, be it physical or not, generate a continuous path of motion for every agent in the system. Thus, even though the rules are not directly derived from physical laws, the resulting continuous trajectory is nevertheless consistent with the laws of physics.

Another approach is to assume that laws of physics are of secondary importance or even ignore them completely altogether. The well-known Conway’s game of life [9] is a good example of this class of models. Depending on the initial setup and rules, distinct spatial patterns can arise from a two- or three-dimensional cellular automaton. There are those who have even argued that many automaton creatures generated using simple rules are indistinguishable from real-life organisms [7]. However, cellular automaton rules are typically highly abstract, not paying much attention to how they can be realized in the real world. Between the two extremes of true physical and purely conceptual models, there exist hybrid ones. For example, consider the migration of households [10, 11] across neighborhoods in a city. In this case, the migration’s intermediate trajectory is usually considered of little importance. Instead, only the origin and the final destination matter. The migration of households could be explained in terms of simple rules as in the Schelling model [10]. Sometimes such a migration can be explained in terms of the physics-like rules. For example, we can define the potential field as a function of the spatial distribution of all agents in the system [11]. As agents relocate and situation unfolds, the potential field also changes accordingly, which in turn affects agents’ migratory behavior in the subsequent period. This coevolving system is known to be associated with adaptive fitness landscapes or adaptive potentials [12].

We can similarly define the fitness landscape for households in a city. A household migrates to another neighborhood in order to improve its welfare, that is, its fitness level. For example, a household relocates because of its preference for same-ethnic neighbors [11]. The migration dynamics is comparable to the dynamics of Brownian particles [5]. A Brownian particle also moves on an adaptive landscape to improve its potential, but its trajectory is a continuous one. The typical trajectory of Brownian particles may not be appropriate to describe the movement of migrating households. A particle cannot pass through a zone that has higher potential than its present one. By contrast, any migrating households can pass through Beverly Hills even though housing there is unaffordable to most households. Similarly, but in opposite sense, buffalos frequently move between groups during mating periods at the risk of greater exposure to predation [13]. Clearly, conventional laws of physics cannot describe the rules of discrete migration since the individual trajectory is not continuous. But those rules are still not purely conceptual, in the sense that rules are based on the physics-like ansatz. In describing discrete migration dynamics, the rules are not defined in terms of acceleration or velocity but in terms of position. In addition to migration of households, the fission-fusion of animal groups is another example [13, 14]. And though not transpiring on physical space, the dynamics of human networks is an example in which social distance between groups replaces the role of geometrical distance [15].

This study links discrete migration dynamics to an agent-based model in order to explain the emergence of distinct spatial patterns in complex systems. To that end, a dynamical system composed of two types of migrating individuals is defined and investigated. In that system, compatibility between two individuals is described in terms of differential attraction coefficients or homophily parameters [16]. The original concept of “homophily” developed in the sociology literature is based on the principle that similarity breeds compatibility [17]. Similarities in race, ethnicity, age, religion, education, and gender all create strong homophilous relationships. The result is agents that are attracted to same-type individuals sharing similar characteristics, while at the same time repelled by dissimilar individuals, or perhaps attracted but only to a lesser degree. Under this pretext, positive homophily arises from the interaction of two same-type individuals that benefit from the presence of each other; while negative homophily is generated when two incompatible individuals of different types attempt to evade each other. A more general notion may also allow positive homophily to result from the interaction of dissimilar individuals while negative one could stem from the interaction of same-type individuals. Similar homophily principles have been applied in biology to examine the preferential mixing behavior [18] of animals and in sociology to understand in-group attraction and out-group avoidance [19]. The general setup of our homophily-augmented agent-based model is briefly described next.

Starting from an initial random distribution of agents on a two-dimensional square lattice, individuals relocate under a single move condition defined in terms of adaptive potential and frictional cost. Specifically, an agent migrates to a new location if its increase in potential net of travel costs is positive. We show below that the emerging spatial patterns of the system population depend on specific values of the homophily parameters, which capture the extent to which individuals are attracted or repelled by each other. A key result of this study is the phase-diagram representation of the emerging patterns on the homophily plane. We then demonstrate how our model can be used as a tool to study the grouping behavior of ungulates.

2. The Model

Agents in our model interact on a two-dimensional lattice system composed of sites. A total, fixed number of heterogeneous agents of multiple types are distributed on the lattice system in such a way that, at any given time, a location is either empty or occupied by a maximum of one agent. Thus, the average population density is = . Although potentially there is no limit on the number of agent types that our model could accommodate, in this study we focus mainly on a two-group system. Specifically, the agent population can be broken down into two subgroups denoted by set and set , respectively. Set comprises of agents of type while set comprises of agents of type . Hence, = + . For completeness, let denote the set of empty locations, comprising of = sites.

The potential or fitness level of an agent depends on the spatial distribution of all the agents in the system. Specifically, every agent in the system contributes to an agent’s potential. The spatial distribution evolves due to migration, which means every time an agent moves; it affects the adaptive potential of all other agents. Next, define as the homophily parameter that comes into effect when agent i of type and agent of type interact. We thus have four homophily constants corresponding to interactions between and within the two types of agents, namely, , , , and . We assume that empty sites also contribute to an agent’s potential, which can be viewed as agents’ preference for empty space. Denote the contribution of an empty site to the potential of a type- and a type- agent as and , respectively. Note that empty sites do not have preferences, which mean that and are undefined in our model. The potential of an agent i is then defined as in the following equation: In (2.1), denotes agent i’s type and denotes either agent j’s type or an empty cell , while measures the Euclidian distance between agents i and j. The summation is carried over all sites in the system except itself. It is important to note that the power exponent is strictly positive implying that the transmission of influence across individuals decays with distance. The exponent can take any real value, depending on the nature of the interactions. It is shown in [11] that, in general, agent-agent interactions are global when , while local when is greater than 2. In the present study we assume global interactions and fix = 1, which implies that every individual is affected by everybody else in the system. We recognize that in complex-systems studies local interactions are often assumed. For example, the classic Schelling’s model of residential segregation assumes purely local interactions, so that only the nearest neighbors have an impact on the focal agent. Conway’s game of life [9] and the voids model [8] are other well-known examples of purely local interactions.

There are, however, complex systems in which purely local interactions alone are not sufficient to understand the observed communal behavior. For example, a select group of living cells behave as if they are “aware” of the state of their entire system [20], fish adjust its swimming velocity to stay united with all other fish in its school [21], and traveling vertebrates appear to utilize global information in order to move collectively together in certain directions [22]. There have been also hybrid models of self-organizing systems that combine local with global interactions, for example, to study the motility-proliferation cycle of malignant brain tumor cells [23]. The distance-decay function that we employ here reflects the combination of short-range, local interactions with long-range influences of agents located farther away that weaken with distance but never become completely negligible.

The potential of an empty location can also be calculated using (2.1), which can be thought of as the hypothetical potential that an evaluating agent would realize if it chooses to move to an empty site. An empty location’s potential therefore depends on the type of the evaluating agents. When an agent migrates to enhance its potential, the same empty location acquires a different value (potential) depending on the type of the relocating agent. Relocation therefore alters not only the potential of the migrating agent but also that of all other agents and empty sites in the system as well. In this way agents in our model climb up the hills on an adaptive fitness landscape.

We specify next the condition under which migration occurs. The key behavioral assumption here is agents’ desire to find locations that yield higher potentials. However, since a site’s maximum carrying capacity is one agent at any given time, and since we do not allow swapping between agents, the only possible relocation of an agent is into an empty site. Consider the case in which agent a evaluates the attractiveness of empty site e. That agent migrates to empty site e if the following condition is satisfied: where f is the friction factor and the distance between agent a’s current location and the empty site e. The right-hand term corresponds to moving cost that is proportional to the distance between the agent’s origin location and its destination. Thus, according to condition (2.2) an agent relocates if its increase in potential obtained by doing so is greater than the moving cost.

The dynamics of the system in our model is then completely defined by the sequence of matching-evaluating-moving. At every single, discrete time step, only one agent is allowed to move, if ever. Specifically at any given time, an agent and an empty cell are randomly selected and matched. The selected agent then evaluates the attractiveness of the empty cell. Migration occurs if condition (2.2) is satisfied. If it is not, no action is taken, and the simulation advances to the next time step. Starting from a random initial distribution of agents, we repeat the simple process of matching-evaluating-moving until no agents have any incentives to change location anymore, or until a long-run spatial pattern can be identified.

3. Simulations and Results

We demonstrate next the dependence of the emerging spatiotemporal pattern on the model’s parameters. Since our goal is to establish the connection between the extent of agent heterogeneity and population dynamics, we focus on the impact of changing the homophily parameter values while holding fixed the other system parameters. The reported results here are based on agent-based simulations of a system in which the lattice size is 100 100, the subgroup populations are = = 500, and the friction factor is = 0.1, unless stated otherwise. With two types of agent, there are six agent-specific parameters whose numerical values remain to be determined, namely, and . The first three parameters influence directly the migratory behavior of type- agents: represents the in-group homophily among the type- agents and the extent to which a type- agent affects a type- one. More specifically, a positive implies that type- agents benefit from the presence of type- ones, while a negative value indicates that the former suffers a reduction in potential from the presence of the latter. The parameters are not constrained to be symmetric across types so that, in general, . Of course, symmetric homophily, , is automatically included as a special case.

We can easily show that, without loss of generality, the set of parameters is shift invariant (see the appendix), which implies that the set of parameters is equivalent to the set for arbitrary constant . Using this fact, we can standardize the parameter set such that . Thus we are left with four nonzero homophily parameters , and . This reduction in the number of nonzero parameters substantially reduces the computation time required to compute agent potentials as defined by (2.1). The reward of standardizing the homophily parameters is that we no longer have to compute empty sites’ contribution to agent potentials since , and hence at any given time we need to keep track of only 10 percent of total cells that are occupied by agents. This results in 90 percent savings of the original computation time without standardization.

Although the number of homophily parameters to be investigated has been reduced to four, an exhaustive search is still impractical. We can further constrain the values of these four parameters to those satisfying 1.0, 1.0. These constraints can be thought of as a second standardization that normalizes the sum of the homophily parameter values associated with a particular agent type to unity. With these two constraints, we have only two independent parameters left. We select (,) to be the independent parameters.

We show next the distinct spatiotemporal patterns that emerge from the system corresponding to different combinations of the two independent homophily parameters. In general, starting from a random initial distribution, the long-run behavior of the system converges toward two broad classes of temporal outcomes. Under some parameter setting, the system settles down to a stationary equilibrium pattern. Every agent in the system is satisfied with its current location or, equivalently, there is no possible relocation that would make an agent profited from a move. This case therefore corresponds to the static regime. Under some other parameter setting, however, the system never settles down. Instead, the system remains out of equilibrium even in the long run, with agents constantly moving in a state of perpetual motion. This case corresponds to the dynamic nonequilibrium regime in which the prevailing spatial pattern is ever-changing. Figures 1 and 2 show examples of spatial patterns that emerge under both stationary and dynamic regimes. A select number of distinct patterns were already reported in a previous study [24]. Nevertheless, a comprehensive study establishing the link between patterns and the model parameters has not been done until now.

Figure 3 unifies the various pattern-parameters pairs using compact phase diagrams. In that figure, patterns are classified based on visual inspection. Figure 3(a) shows the case of and . Only the upper triangular part is shown since the lower part is symmetric with respect to the diagonal line where . Figures 3(b) and 3() corresponds to the cases in which at least one summation over the homophily parameters, or , is negative.

In the phase diagrams, dynamic nonequilibrium regimes are represented as shaded areas. An interesting question is whether the phase transition lines between the dynamic and static regimes can be estimated. We will not attempt to fully solve the question here. Instead, we will present a simple approximation to help understand the phase transition that separates the two regimes. The starting point is the recognition that dynamic regimes can be obtained when the two types of agents interact as if they were predators and preys. This predator-prey relationship can be approximated by the following condition: For example, if and , then type- agents are the predators that benefit from closer proximity to type- agents, which are the preys since they suffer from such proximity. If we apply (3.1) to cases shown in Figure 3(a) where and , then we obtain the dynamic regime almost surely when and . This approximation seems to hold, at least for those dynamic patterns plotted in the phase diagram. It lacks precision, however, because the actual zone covered by the dynamic regime is smaller than that predicted by the approximation. This is due to the friction factor that is not reflected in approximation (3.1). In particular, the dynamic zone gets smaller as the friction factor gets larger, and vice versa. If the friction factor is increased to a level greater than some critical value, the system is “frozen” in its initial state [16], which means the static regime trivially rules the entire phase plane.

Table 1 lists the long-run patterns of the population spatial distribution found in the present study. It is important to keep in mind that the list is not meant to be exhaustive. A number of patterns were already identified and reported in a previous study [24], including the worm and the cell modes. However, the list also includes newly identified patterns such as the herd, the amoeba, and the stripes modes. One contribution of this paper is the pinning down of the patterns’ relative positions vis-à-vis others on the phase diagram. The cell pattern in particular is fascinating since it can be thought of as representing self-replicating cells in a live organism. In contrast to cellular automatons which may be considered the most abstract form of life,the cell mode here can be viewed as a life-form at an intermediate level of abstraction. It is also interesting that such complex life-form lies just on the edge of the static regime in the phase diagram.

4. Applications

One potential application of our model is to study the population dynamics of social systems. Many social groups experience continual change in terms of both size and memberships. Changes occur when groups split into smaller ones, or when small clusters merge to form a larger aggregation. For example, select ungulates form groups of sexually segregated habitats [13, 14] except in the mating season. During the intermittent periods, the two groups intermingle only when they are grazing. When they rest, the two groups again self-segregate themselves sexually. Our model can be applied to explain such a dynamically changing communal behavior. Denote male and female ungulates by type- and agents, respectively. Our grid system in this case represents their habitats in which the two types of agents are spatially distributed. We can then simulate the emerging spatiotemporal patterns of the heterogeneous agent population using the present model.

During rest time, an ungulate agent prefers the company of another same-sex agent to a different one. This can be represented by the conditions and . If during rest time all agents prefer to live next to a neighbor than to be isolated, then a location next to an empty grid would be the least desirable. In this case, the homophily parameters could be set equal to, for example, = (0.8, 0.2, 0.0) for males and = (0.2, 0.8, 0.0) for females. The parameter set in this case is (0.8, 0.8) in Figure 3(a). The corresponding pattern is a stationary separation, like shown in Figure 1(a). During grazing however, agents consider a location next to an empty grid—interpreted as the ungulates’ common grassland—to be the most desirable. Thus the homophily parameters during grazing may be set equal to = (0.8, 0.2, 1.0) for males. Because of the shift invariance property, this configuration is equivalent to = (−0.2, −0.8, 0.0). The parameter configuration associated with grazing period therefore corresponds to the distributed pattern, such as that shown in Figure 1(d), which again is shown in Figure 3(b).

All the patterns that have been discussed so far assume a small friction factor of = 0.1. We do not investigate cases with larger friction factors here. The effect of the friction factor has been well studied for a homogeneous, single-type agent population [16] in which the system dynamics are driven by only one homophily parameter in addition to the friction factor. If we standardize the homophily parameter to 1.0, then the population dynamics are fully determined by the friction factor for a given population density. In this case, depending on the values of the friction factor, four distinct spatial patterns can be identified in terms of the population size distribution, namely, the random uniform, the exponential law, the power law, and the black-hole patterns [16]. The friction factor of 0.1 then corresponds to the black-hole pattern in which the entire population converges into a giant cluster encompassing all agents in the system. In particular, if we set representing essentially a homogeneous population, then all agents regardless of type collapse to a single group, that is, the black hole. The reason why in the present study we observe separated groups as shown in the phase diagram is due to negative, homophily parameters.

When the homophily parameters are all positive and the friction factor = 0.1, we obtain at most two groups in the long run. But at higher friction factors we obtain a distributed system of groups—even with all-positive homophily parameters—where each group can contain a mix of different types. An example is a multiracial city system. Cities in such a system form spatially distributed groups, each composed of different types, that is, races. In a multiracial system the degree of segregation between the different races is often of interest. Another example is associated with mammalian groups. During the grazing period, ungulates self-organize themselves into smaller groups, each of which is sexually segregated. In order to quantify the overall degree of gender partitioning among ungulates, we compute the segregation coefficient (SC) as follows [14] where and denote the total number of males and females in the system, respectively, and = + . There are -groups in the system. In the th group, the male and female subpopulations are denoted by and , respectively. Groups with only one member are excluded from the calculation of the segregation coefficient. Our model can be easily applied to see how segregation behavior responds to changes in the homophily parameters. Again consider a population of ungulates where the group behavior of males is represented by the parameters = (0.8, 0.2, 0.0) while for females = (0.2, 0.8, 0.0). The friction factor is set at = 0.45, much higher than before in order to focus on cases where ungulates migrate not solely because of homophily forces, in particular during grazing periods when they strongly prefer to live next to an empty grassland.

With higher friction factor, the groups remain spatially distributed as shown in Figure 4. But the homophily forces remain sufficiently strong that each group exhibits uneven distribution of the two sexes. Note that we defined a “group” here as the cluster of locations connected by at least one non-empty site. From simulations of a 200 200 lattice system with = 1600 and = 2400, we obtain the average SC equal to (mean, standard deviation) = (0.35, 0.045) and the average group size for is 3.2. The segregation coefficient clearly suggests the tendency for males and females to live in gender-partitioned groups. Changing the male homophily parameter values to (0.7, 0.3, and 0.0) reduces the segregation coefficient to (0.16, 0.051). Thus, we could calibrate the parameters in order to fit the available experimental data and then perform simulations to test whether sexual homophily causes gender segregation, but not the other way around. Finally, we can also apply the model to a heterogeneous population with more than two types of agents. In Figure 5, we present an interesting spatial pattern generated from a system with three agent types. With greater number of parameters, more complex patterns could be generated.

5. Concluding Remarks

We have developed an agent-based model to investigate how diverse spatial patterns can be generated from a minimal set of interaction rules. Numerical simulations of the model allow us to identify distinct stationary and dynamic patterns of the population spatial distribution, and we present phase diagrams on the planes of homophily parameters. The complex spatial patterns are generated primarily by agents’ simple adaptive behavior, and not by the internal complexity of the agents themselves. The phase transition between the dynamic and stationary regimes is also identified and explained. Our model can be applied to study the group forming behavior of social communities. As an example, we apply the model to explain the sexual segregation of ungulates in terms of the homophily parameters and the friction factor. Future extensions of this study could include the development of the generalized Schelling model of segregation in a heterogeneous population with n types of agents. Also it will be an interesting topic of research to classify social systems in terms of the various spatial patterns of the population associated with the dynamic regime presented in this study.

Appendix

The Shift Invariance of Homophily Parameters

This appendix shows that the set of parameters is indeed shift invariant, which implies that the set is equivalent to the set for an arbitrary constant .

By definition, the potential of an agent a of type is given by the following equation:

In (A.1), the sets , , and correspond to agents of type and , as well as empty cells, respectively. Now, subtracting a constant c from every homophily parameter yields

Similarly, the potential of an agent b of type is obtained as follows:

Since we are dealing with systems exhibiting periodic boundaries, it follows that where is a constant, which does not depend on the specific agent under consideration. Thus, it follows that where is just another constant. This proves that the parameter sets and are equivalent. Therefore, and are shift invariant. This completes the proof.

Acknowledgment

This research was supported by the Yeungnam University research grants in 2009.