## Abstract

It is shown that Boltzmann's methods from statistical physics can be applied to a much wider range of systems, and in a variety of disciplines, than has been commonly recognized. A similar argument can be applied to the ecological models of Lotka and Volterra. Furthermore, it is shown that the two methodologies can be applied in combination to generate the Boltzmann, Lotka and Volterra (BLV) models. These techniques enable both spatial interaction and spatial structural evolution to be modelled, and it is argued that they potentially provide a much richer modelling methodology than that currently used in the analysis of ‘scale-free’ networks.

## 1. Introduction

It is well known that the methods of statistical mechanics introduced by Boltzmann provide powerful methods for model building. His methods have been generalized and applied widely, in particular to predict flows between zones in spatial systems. It is also the case that Lotka's and Volterra's prey–predator (PP) model has been applied in a variety of contexts beyond ecology. In this paper, we draw attention to the way in which the two approaches can be combined in the building of models to represent urban structure. It can then be seen that this combined approach is a fruitful one for a much wider class of dynamical systems than has hitherto been realized. The resulting models can be characterized as the Boltzmann, Lotka and Volterra (BLV) models. The purposes of this paper, therefore, are threefold: (i) to elucidate this methodology in general terms using illustrative simple assumptions, (ii) to explore the range of application, and (iii) to show that the BLV methodology potentially provides a much more powerful basis for generating and analysing networks.

We first outline the BLV method in general terms (§2). In §3, we outline the range of application. In §4, we show how this methodology overtakes much of the scale-free modelling enterprise. In §5, we outline the challenges for future development.

## 2. The BLV methodology

The deployment of the BLV methods in urban and regional geography is well established (Wilson 1970*a*,*b*, 2000), and it has recently been shown that they can be extended to ecology (Wilson 2006). Here, we present the methods in a general format, though with one example, an urban retail system, to help fix ideas.

Consider a spatial system characterized by a set of origin locations {*i*}, such as zones (and zone centroids) or points on a network and a set of destination locations {*j*}. Let {*X*_{i}} be a set of measures of activity at each origin and {*Z*_{j}} a set of measures of activity at each destination. In terms of the urban retail example, the {*X*_{i}} might be demands for a consumer good and the {*Z*_{j}} the availability of those goods in retail centres. Then, let {*Y*_{ij}} be a matrix whose elements represent the flows (or ‘interaction’) between zones *i* and *j*. In the specific example, it could be the flow of expenditure from the consumers in *i* to shops in *j*. In the first instance, we neglect to represent the underlying network explicitly but can do so implicitly through a set of measures of the impedance of the interaction between each *i* and each *j*. This can be a matrix {*c*_{ij}}. Each *c*_{ij} will usually be the sum of different kinds of costs associated with carrying flows along the routes of the real network—the so-called ‘generalized’ cost. We emphasize that, notwithstanding the example, this is a general formalism.

There are many possible models that can be generated within the framework of this methodology. Let us first assume that the {*X*_{i}} are fixed and given, and that the modelling task is (i) to estimate the flows {*Y*_{ij}}, while also given the {*Z*_{j}}, and (ii) to model the evolution of the {*Z*_{j}}. It is the {*Z*_{j}} in this case that we are taking to represent spatial structure. We assume that the modelling of the interactions {*Y*_{ij}} represents a fast dynamic, with a rapid return to equilibrium after a change, and the modelling of the structural variables {*Zj*} is a slow dynamic.

By a suitable choice of units, we can assume that the following constraints hold:(2.1)and also that(2.2)where *C* is a known total—analogous to an ‘energy’ in the Boltzmann terms—to support the interactions. We also add constraints(2.3)where *b*_{j} represents the benefit of going to *j*. To link with the structural variables, we assume that(2.4)That is, if *Z*_{j} is a measure of ‘size’, then the benefit of interaction with *j* is proportional to the log of the size. Hence, (2.3) can be written(2.5)Note that if we can estimate {*Y*_{ij}}, then the set of total flows into a destination *j*, say *D*_{j}, can be calculated as(2.6)

A good estimate of {*Y*_{ij}} can then be obtained by Boltzmann's methods, essentially a statistical averaging procedure (Wilson 1970*a*). We can work directly by analogy with statistical mechanics. The system can be treated as a microcanonical ensemble and the most probable {*Y*_{ij}}, and hence the model estimate is obtained by maximizing an entropy term. The entropy measure based on the microcanonical ensemble used here can be shown to be formally equivalent to the information theory measure (see Wilson 1970*b*). Thus, maximize(2.7)subject to the constraints (2.1), (2.2) and (2.5). It can easily be shown that this gives(2.8)where(2.9)It can be seen that, from a statistical physics perspective, 1/*A*_{i} plays the role of a partition function. *α* and *β* are parameters, derived from the Lagrangian multipliers associated with (2.5) and (2.2), respectively. In the physics model, they represent inverses of different kinds of temperatures. In the applications envisaged here, the parameters *α* and *β* can be expected to be positive because the size *Z*_{j} is assumed to be a benefit and travel impedance, and *c*_{ij} is assumed to be a cost.

Note that can be written as exp(*α* log *Z*_{j}) and hence(2.10)showing explicitly how log *Z*_{j} as a benefit measure is set against the impedance or ‘cost’ *c*_{ij}—weighted by *α* and *β*, respectively.

Equations (2.8) and (2.9) represent the most probable estimates of {*Y*_{ij}} given {*X*_{i}}, {*Z*_{j}} and {*c*_{ij}}. It is reasonable in this case to assume that, because this is overwhelmingly the most probable state, this represents an equilibrium state, and that, for any change in the *X*, *Z* or *c* variables, there will be a rapid return to a new equilibrium, and hence this represents the ‘fast dynamics’.

To model the {*Y*_{ij}}, we assumed the spatial structural variables {*Z*_{j}} to be fixed. We now consider the possibility of these changing on the basis of a ‘slow dynamics’ hypothesis. Note that, from equation (2.6), we can calculate the total flow into each *j*, *D*_{j}. These flows have been attracted by a pulling power at *j* of the structural variable *Z*_{j}, raised to a power . For illustrative purposes, we can now hypothesize (again assuming that everything is measured in commensurate units) that if *D*_{j}>*Z*_{j}—that is if the inflow at *j*, *D*_{j}, exceeds the notional capacity *Z*_{j}—then *Z*_{j} should grow, and vice versa—if *D*_{j}<*Z*_{j}, then *Z*_{j} should decline (Harris & Wilson 1978). Of course, there are many other possible hypotheses but this assumption serves to illustrate the method: how Boltzmann can be combined with Lotka and Volterra. In reality, of course, it would be important to incorporate stochastic fluctuations in the dynamics, but we retain the simple assumption again for clarity of illustration. This hypothesis can be represented by the following set of difference equations:(2.11)for a suitable parameter *ϵ*. (An extra factor *Z*_{j} has been added, which determines the form of growth near *Z*_{j}=0.) These equations are recognizably related to the Lotka and Volterra equations, albeit, in this case, with the ‘populations’ being spatially distributed but a single species. These equations, in various forms, turn up in many fields, particularly ecology, and have been extensively studied.

In this general case, we can note first that—setting Δ*Z*_{j}=0 in (2.11)—the equilibrium condition for {*Zj*} is(2.12)and it is instructive to make the substitutions for *D*_{j} from (2.6), (2.8) and (2.9) to give an explicit set of equations in {*Z*_{j}}—combining Boltzmann, Lotka and Volterra—the total inflows predicted from the Boltzmann component and the evolution of the structural variables from the Lotka–Volterra hypothesis(2.13)

Equations (2.13) are complicated nonlinear equations in the {*Z*_{j}}. They can be solved numerically, and it is possible to obtain some explicit insights into the way in which spatial structures such as {*Z*_{j}} evolve. The l.h.s. of (2.13) can be interpreted in a straightforward way. Each *i* term on the l.h.s. of (2.13) is the share of the origin total at *i* that flows to *j*. The factor(2.14)in this expression is the proportion of *X*_{i} going to *j* and is determined by the ‘attractiveness’ of *j*, measured as and the impedance of connecting *i* to *j*, measured as exp(−*βc*_{ij}). It is inversely proportional to(2.15)which represents the ‘competition’ of other destination zones in an obvious way.

Nonlinear equations exhibit characteristic properties and this is true of the BLV equations. They often have *multiple stable solutions*. Typically, particular solutions are *path dependent*, i.e. the solutions at a point in time are dependent on the *initial conditions* at earlier points in time and this has implications for our understanding of the evolution of structures. And there are *phase transitions*—sudden jumps from one kind of structure to another. The solutions to equations (2.13) can be shown to exhibit these properties (Clarke & Wilson 1985).

The equilibrium solution gives the spatial structure at a point in time—if equilibrium has been achieved. It is more likely in practice of course that stochastic fluctuations in exogenous variables mean that we should think in terms of ‘progress towards equilibrium’. As parameters (or indeed the {*X*_{i}} or the {*c*_{ij}}) change over time, the {*Z*_{j}}, as represented by the equilibrium state, will obviously change (through equation (2.11)) and this is the sense in which it is possible to model the evolution of spatial structure. We have mentioned phase transitions: there will be critical parameter values at which sudden changes in structure can take place. One such transition has been demonstrated empirically for the ‘corner shop’ to supermarket transition in retailing (Wilson & Oulton 1983).

We can take the interpretation of the {*Z*_{j}} equations (2.13) a stage further by focusing on the structures that can arise from different combinations of *α* and *β* values (Clarke & Wilson 1985). Model simulations and analytical interpretations concur with intuition: larger *α* implies greater attractiveness as measured by the size, and smaller *β* implies higher probabilities of interaction over longer ‘distances’ to achieve the benefits of size. Hence (high *α*, low *β*) combinations generate structures with a small number of large *Z*_{j}'s, ‘centres’, and vice versa.

## 3. The range of application

We can now review the range of actual and potential applications. Recall the characterization of our generalized network: the core is a matrix of interactions between two sets of points {*Y*_{ij}}, and the structures are represented by the scale of activity attracted, {*D*_{j}}, and the associated infrastructure {*Z*_{j}}. The flows are carried on an underlying network. If the elements that make up the flows are weakly interacting, like the individuals going on shopping trips, then a Boltzmann statistical mechanics methodology can be used to model the interaction matrix. If the attracting nodes are competing with each other to achieve their levels of activity, then a Lotka–Volterra kind of equation is likely to be a good starting point for modelling structural evolution. In each particular case, of course, the bare bones of the models outlined above have to be developed further. We will show in the rest of this section that many varieties of models can be developed within the formalism.

If it happens that the conditions are not satisfied for a Boltzmann-like interaction model, say for protein–protein interactions in biology, then the generalized framework may still be valuable. The hypotheses that have been used in our illustrative models can be amended. Some of the examples cited below use the techniques of network analysis, which we believe could be fruitfully replaced by the BLV methods. We explicitly address the scale-free networks (SFN) enterprise in §4.

Retailing, our earlier example, is the archetypal BLV system from urban *geography*. Consider equations (2.8), (2.9), (2.11) and (2.13). These can be taken as representing an idealized retail system if the variables are defined as follows. The {*X*_{i}} represent the (spatially distributed) spending power within the system, the {*Z*_{j}} represent the attractiveness of retail centres (measured by size) and, hence, the {*Y*_{ij}} are the flows of expenditure on retail goods from the residents of zones *i* to retail centres *j*. The *c*_{ij} represent (generalized) travel costs. Equations (2.8) and (2.9) predict the interaction pattern, (2.11) the dynamics of retail centre spatial structures and (2.13) the equilibrium retail structure. These models, particularly the interaction submodels, have been widely applied in both academic and commercial environments (e.g. Birkin *et al*. 1996).

These methods can be extended to modelling all aspects of the urban and regional systems. The classical models that represent different aspects of geographical structure are based on underlying interactions with structural variables that can be modelled with greater generality with the BLV methods (Wilson 2000). The methods can also be deployed in multi-regional demographic modelling (Rees & Wilson 1977). For multi-regional economic modelling, see Wilson (1970*a*,*b*); building on Leontief & Strout (1963), for the theory; and Roy & Hewings (2005) for an example.

The classical Lotka–Volterra models in *ecology* are the PP and competition-for-resources (CFR) models. We demonstrate how an *N* species version with populations {*Z*^{m}} can be built as a BLV model. In this case, we are adding the ‘B’ to the ‘LV’. *m* (and later, *n*) is a species label and we introduce this as a superscript to facilitate the introduction of spatial zone labels, *i* and *j*, as subscripts later. An obvious notation can now be developed. *Note, however, in this case, that the {X} and {Z} populations are one and the same*. They are both the same *species* populations, either on a PP basis or on a CFR basis, interacting with each other. In the CFR case, each ‘resource’ also counts as a ‘species’. The dynamics, following the argument in Wilson (2006), can be presented in terms of the carrying capacities {*D*^{n}} for each population as(3.1)The {*D*^{n}} are shown as formally dependent on all the population variables. We can be more specific by setting(3.2)to give(3.3)which resembles equation (2.11). Different assumptions for the signs of *α*'s will generate appropriate models—pure PP, CFR or mixed. The details are given by Wilson (2006). The dynamical properties of these equations—the general and the special cases—are well known and so will not be articulated in detail here.

We can now absorb the −*Z*^{n} term in (3.3) into the *α* coefficients by redefining *α*^{mn} as *α*^{mn}+1 without the loss of generality and then introduce space into the model through our zone systems {*i*} and {*j*}. An obvious generalization of (3.3) is(3.4)We have now introduced the spatial zone labels, *i* and *j* (in ecological terms, ‘patches’), as subscripts to bring the notation nearer to that of the core BLV model.

We can then make the spatial interaction terms explicit using the Boltzmann methods and this gives a new ecological model (Wilson 2006). There is only one empirical example in the literature known to the author (Smith 1983). There are two elements of knowledge transfer from urban to ecological modelling here: first, an explicit and potentially more effective way of handling spatial interaction and, second, the ways of analysing the dynamics of spatial patterns.

*Economics* is something of a special case. The retail model, for example, is clearly a model of consumer behaviour and this is an important element of microeconomics. However, economists prefer to view the rationale for the model within their own theoretical framework, maximizing a utility function. They can achieve this by defining a utility function that has the properties of the entropy function used in the model development here and, of course, deriving equivalent results. How the model is perceived is perhaps a matter of preference. However, there is one key issue in comparing the entropy-based models with more traditionally formulated economic models: they represent, from a market perspective, a suboptimal solution. Consumers do not travel the shortest distance to achieve their maximum utility. The shortest distance is achieved by letting *β*→∞ (Evans 1973). A finite value of *β* indicates longer trips on average. This is likely to be closer to reality, of course, because there will be dispersion from the ‘optimum’, partly from the lack of perfect information and partly because there will be preferences not captured by the model.

The {*Z*_{j}} mechanism has an attractiveness factor and, when *α*>1, there are *positive* returns to scale. This turns out to be very important in creating hierarchical structures and networks but is, of course, contrary to the usual assumption in the economics of *diminishing* returns to scale. We begin to have a model that explains the SFN distributions and takes us beyond empiricism. In the case of all the urban models, translated into an economic framework, it is possible to add evolutionary equations for prices and rents also using the Lotka–Volterra principles. These principles have been widely used in economics in a number of other ways, most directly and obviously in resource management, especially fisheries management (Rosser 1991). They have also been applied in a very interesting way in organizational settings (e.g. Allen *et al*. 2006).

An example in *epidemiology* relates to the spatial diffusion of infection. There are network models of this process examining the links from a host to a recipient (who may or not become infected) and then an infected person can either be taken out of the population (in isolation) or in turn become a host. There are many examples, for instance, exploring the spread of viruses across the web, as well as the more traditional epidemiology of infectious diseases (Moreno & Vazquez 2003). There should be a way of formulating the interaction elements of this modelling task in the Boltzmann form. We then need to model the sequence of events through time: the {*Z*} ‘structures’ would represent the concentrations of infection.

*Chemistry* offers a different kind of example. The physical chemistry of mixtures provides examples of the application of equations that are structurally similar to the Lotka–Volterra equations, and perhaps the best known example is the Belousov–Zhabotinsky reaction that involves two liquid chemicals oscillating in proportions in the mixture (Gray & Scott 1990). It is also possible to model chemical reactions using a more general form of the Lotka–Volterra equations with particular reference to the role of catalysts (Jain & Krishna 2001). Interestingly, in their case, they see the modelling of the (chemical) populations as the fast dynamics and the graph of interactions between them as the slow dynamics.

There are many examples in *biology*, for example using the methods of network analysis in cell biology and proteomics (Albert 2005). In these cases, the networks are very different and, because there will be high levels of interdependence between elements, the Boltzmann conditions are unlikely to apply.

A very different kind of application of network analysis is offered by *geomorphology* (Rinaldo *et al*. 2004). There are explorations of the structure of networks in the landscape, particularly river networks. This in turn connects to much earlier work on such networks (Haggett & Chorley 1969).

*Physics* provides examples that give us an alternative approach for formulating the evolutionary model. This involves deploying a Hamiltonian framework and is both interesting in its own right and leads to the possibility of incorporating a entropy term on the structural variables. The techniques of interest arise out of modelling populations at the nodes of a lattice, whether spins or numbers of molecules.

We begin with the Ising model (following Thurner 2005). The model is concerned with spin systems and the alignment of spins at certain temperatures that produce magnetic fields. The model is usually based on a lattice. In the case of the one-dimensional model with only neighbour-to-neighbour interaction, and spins *σ*_{i} at each lattice point *i*, then the Hamiltonian is(3.5)where *J* measures the strength of the interaction, and the probability of a particular configuration *r* occurring is(3.6)We can translate this into our problem as follows. Consider our *j* nodes to be the nodes of a lattice. Let *r* label a configuration of the structural variables {*Z*_{j}}. The task, then, is to find a Hamiltonian *H*_{r}, as a function of the structural vector. We can then write (3.6) as(3.7)and we have to find the {*Z*_{j}} that maximizes *p*_{r}. Suppose we take the measure of profit from our BLV methodology(3.8)where *K* is a unit cost of provision and *D*_{j} can be obtained from equations (2.6), (2.8) and (2.9). It can then be shown that this produces a result that is equivalent to equation (2.13).

## 4. Scale-free networks

At this point, it is convenient to confront the SFN literature and to contrast it with the BLV methodology. We will indicate that the latter has the potential to overtake the former. The key variable that is used in the SFN analysis is the number of edges at each node, and *p*(*k*) is taken as the probability that a node has *k* edges connected to it. It is found empirically that this distribution often takes the form of a power law(4.1)for some parameter *γ* and a normalizing constant *a*.

A ‘network’ is defined to be a *single* set of nodes, and this is a much more restricted concept of a network than that which underpins the BLV models. To develop the argument in the terminology of this paper, we can consider a single set of nodes either as our {*i*} or {*j*} sets. For definiteness, consider the {*j*} set. Then, if we measured *Z*_{j} by the number of edges—which might be flows above some threshold—at *j*, then the size distribution of the *Z*_{j} would be equivalent to the *p*(*k*) distribution.

We can immediately see how the BLV models are much more general than the usual SFN models in a number of ways.

They are based on a more detailed characterization of the system of interest, in particular of the network and associated flows.

The flows are being estimated explicitly, and this is richer than simply counting edges. A method for articulating network structure based on flows was given by Nystuen & Dacey (1961) and used by Rihll & Wilson (1987) in their work on settlement patterns in Ancient Greece: edges can be counted if that is thought appropriate within the BLV formulation.

Because there is an explicit hypothesis for the evolution of {

*Z*_{j}} in the BLV formulation, and for the equilibrium structures that emerge, (*a*) an explanation is offered for the spatial structure and the size distribution (and this is typically not the case in the SFN analysis) and (*b*) if a power law is found for the distribution of activity levels at nodes, whether measured by edges or otherwise, then the BLV model will offer an insight into how this is generated.The SFNs are not really networks except in a topological sense. The ‘edges’ are more like ‘significant interactions’. In many cases, as we have seen, there are real networks that underpin the flow models—road networks being an obvious example. In this case, it is necessary to model both the interactions and network loads (which are not the same thing) and the relationships between them.

Note that if *p*(*Z*) is the probability of a destination having size *Z*, then it is this distribution that is to be compared with the *p*(*k*) of SFN. It is clearly a more sensitive measure and it has been modelled in such a way that it is predicted by a theory—that which is being modelled—about the system of interest. If this distribution turns out to be a power distribution, as in the case with SFN, this will be a by-product of the analysis rather than a prime focus.

The conventional literature on network evolution relies on network-generating algorithms. It is possible to classify these into three groups: random, scale-free and small world. For the first category, edges are added on a random basis (Erdos & Renyi 1959). There is less of a hierarchical structure than in the second (Barabasi & Albert 1999) where the algorithm is modified to make it more probable that edges are added to nodes that already have more than the average number of connections—the ‘rich get richer’ algorithms. It is not surprising in this case that hierarchies emerge. What is more surprising is their ubiquitous scale-free power-law properties. The third class, the small-world networks, is based on the algorithms that generate a higher propensity of clustering (Watts & Strogatz 1998).

It is not the purpose of this paper to offer an exposition of these algorithms: there are many such. Rather it is to review the extent to which there are systems being modelled using these network algorithms for which greater depth of explanation and understanding could be achieved by employing the BLV methodology or a variant of it. First, any system founded on matrices of flows or interactions—whether people, goods, money, or whatever is likely to be more effectively modelled using the Boltzmann-based spatial interaction models. Second, care should be taken in understanding the real network structures in the system of interest. For most real networks, there will be a set of edges, with intermediate network nodes, which constitute a route, connecting any origin–destination pair. The modelling of interactions can then be separated from the task of loading flows onto the real network—what transport engineers call the assignment problem.

Third, great care should obviously be taken in formulating hypotheses for the modelling of structural dynamics. We have shown that successful hypotheses can be developed for retail and other urban systems and ecosystems, and the range of other examples that have been considered suggests a much wider potential. These models are the basis for creating hierarchical structures. This can be argued for flows of goods, at different scales, from the local up to international trade (Baskaran & Bruck 2005), and can also be applied to utilities flows such as water, gas and electricity. These are carried on real networks, supplied from origins, used by the consumers at destinations. It can be applied to communications technologies, including the World Wide Web (Pastor-Satorras & Vespigniani 2004). It has been estimated that this has 800 000 nodes and it seems to be a clear case where it is potentially important to model origin–destination flows and to recognize that these flows are carried on an underlying physical (or wireless) network. The fact that significant hubs develop is not surprising: it connects to the ‘economics’ argument about scale economies. The challenge is to assemble enough data to enable a sensible model to be constructed. This area has been explored using entropy-maximizing models (Tomlin 2003), but focusing on the Boltzmann component, not, as yet, the dynamics.

## 5. Challenges for the future

This brief analysis provides the basis for indicating some of the challenges for the future.

The spatial interaction models based on entropy-maximizing principles could be applied more widely in different kinds of systems.

The matrices of flows could be analysed to provide richer portraits of network structures.

It would be possible as a matter of routine to note the size distributions of the destination flow totals thus produced to examine the extent to which these matrices represented, in the narrow sense, SFN.

The slow dynamics models could be explored for a wider range of systems, and it has been argued that this provides a new generating mechanism for networks. To what extent, in different circumstances, are these scale free?

There are two aspects of network analysis that need to be borne in mind for future priorities. First, most real networks have hierarchical structures—motorway links relative to country lanes provide an obvious everyday example—and these are worthy of analysis. Second, while we have used LV principles to model the dynamics of nodal evolution, we have not considered edge evolution on the same basis. This can be done in principle (Wilson 1983).

There is a mathematical challenge: to find a way of explicitly connecting the {

*Z*_{j}} size distributions that arise in the BLV models to the statistical distributions used as measures of network structure in the scale-free literature. We need to find analytical expressions for the size distributions of {*Z*_{j}}, which are themselves the solutions of equations such as (2.13), and which cannot be solved analytically. The beginning of a way forward can be found (Rosser 1991; Batty 2006).

## Footnotes

- Received October 31, 2007.
- Accepted November 19, 2007.

- © 2007 The Royal Society