Diffusive Resettlement: Irreversible Urban Transitions in Closed Systems

We propose a non-equilibrium framework for modelling the evolution of cities, which describes intra-urban migration as an irreversible diffusive process. We validate this framework using the actual migration data for the Australian capital cities. With respect to the residential relocation, the population is shown to be composed of two distinct groups, exhibiting different relocation frequencies. In the context of the developed framework, these groups can be interpreted as two components of a binary fluid mixture, each with its own diffusive relaxation time. Using this approach, we obtain long-term predictions of the cities’ spatial structures, which define their equilibrium population distribution.


Introduction
Modern cities are diverse in their spatial structures: Some cities are monocentric, with a single centre of business, retail, and other types of activity, while some exhibit polycentric patterns in which multiple activity clusters are distributed across space [1][2][3][4][5]. It is well known that a city's structure affects economic productivity, environmental conditions, and other aspects of human life [6][7][8][9][10]. Importantly, spatial structures of cities change over time in intricate ways, with the process of intra-urban migration being one of the main drivers of a city's evolution [11][12][13][14][15][16][17][18][19][20]. However, there is a lack of models that can quantitatively explain and accurately predict a city's evolution in terms of intra-urban migration and the resultant patterns of an urban spatial structure.
In existing models of cities, intra-urban migration is usually considered as a fast process in which an equilibrium is reached very quickly [19][20][21][22][23][24]. Such an equilibrium is typically defined by the spatial distribution of infrastructure and employment [19,21,23,[25][26][27], and the intra-urban dynamics are typically not considered as an example of a non-equilibrium process in the thermodynamic sense. A non-equilibrium approach is also hard to explicitly validate due to the difficulty of collecting the data.
In this work, we consider intra-urban migration as an irreversible process, and explicitly derive the dynamics of the corresponding non-equilibrium evolution. In doing so, we shall draw on an analogy with diffusive relaxation. This opens a way towards a systematic and coherent framework describing human migration within cities (both temporally and spatially), as opposed to an unconstrained evolution of cities through expansion.
Human relocation has been widely considered as diffusion in open systems [28,29]. This approach employs an apparent and direct analogy between human migration and molecular diffusion, and has shown good predictions at different scales [15][16][17][30][31][32], from city growth [23,33,34] and epidemic spread [31,35,36] to inter-continental migration [28]. These migration processes can be characterised as expansive, as they increase the area of human habitat, viewing the human society as an open system. In the modern world, however, most of the migration processes result in a redistribution of the population across already occupied locations, rather than expanding to non-occupied areas. This constraint essentially confines the migration to happening in a closed system with a fixed area and population. In this paper, we propose a diffusion model describing migration processes in a closed system from the perspective of non-equilibrium thermodynamics.
In general, migration from rural to urban areas, as well as expansion of metropolitan boundaries, is a relatively slow and long-term process. In contrast, intra-urban migration happens at a much faster rate [19,20,26,27,[37][38][39]. Nevertheless, we argue that intra-urban migration is a fundamental driving force shaping the long-term evolution of the urban spatial structure, with external migration and spatial expansion playing only a secondary role. Thus, we aim to model urban transitions as dynamics developing in a closed system, at least in the first approximation.
Typically, intra-urban human mobility has been considered as a process driven by certain attractiveness of various locations within a city, perceived in terms of proximity to schools, business centres, recreational facilities, etc. [19,24]. This notion of attractiveness was modelled both explicitly, using specific socio-economic indicators [40,41], and implicitly, being reconstructed directly from migration data [26]. The latter approaches can be classified as "microscopic" [42,43], as they focus on a specific mechanism of human relocation.
In this paper, we propose a concise "phenomenological" approach that considers the migration flows from the perspective of diffusion. Importantly, we do not make any assumptions about particular choices that may motivate individuals to relocate. Instead, we analyse their collective movement and show that this process is similar to diffusion. As a result, we reveal the trends of population movement analogously to the macroscopic movement of a fluid. In doing so, we introduce a rigorous definition of an equilibrium state as the spatial configuration to which the apparent evolution relaxes, and a decomposition of the population into two distinct groups with different relocation frequencies.
This paper is organised as follows. In Section 2, we present a general framework of the irreversible evolution of a city driven by residential relocation. In Section 3, we apply this framework to the Australian capital cities. In particular, we analyse the dynamic relocation patterns in Section 3.1, develop the analogy between intra-urban resettlement and diffusion in Section 3.2, and predict the equilibrium population distribution in Section 3.3. In Section 4, we analyse the robustness of the model. Finally, in Section 5, we summarise the findings of this work.

Dynamics of Intra-Urban Migration
We consider an urban area as a set of N suburbs i with a certain residential population x i (t) at time t. The total population at any time is fixed: ∑ N i=1 x i (t) = x. We assume that time is discrete, and the choice of the time step is dictated by the resolution of the available data. A migration flow T ij (t) is defined as a change in residential location from suburb i to suburb j. This flow is uni-directional, so that, in general, T ij (t) = T ji (t), and the net flow J ij (t) ≡ T ij (t) − T ji (t) = 0. Non-zero net flow indicates that the system is out of equilibrium. In a diffusive system, the net flow gradually decays to zero with time as the system evolves towards an equilibrium. Such an equilibrium state is stationary on the "macroscopic" level, showing no change in the population of each suburb. However, on the "microscopic" level, there still exists some movement of people, resulting in non-zero unidirectional flows T ij (t). In an equilibrium, these uni-directional flows satisfy a microscopic detailed balance, so that T ij (t) = T ji (t), resulting in zero net flow J ij (t) between each pair of suburbs.
The uni-directional flow matrix allows one to predict the future population of any suburb. In particular, the population at the next time step t + 1 can be expressed through the migration flow T ij (t) at the current time step t as x j (t + 1) where the sum includes the term T jj (t) accounting for immobile population. Introducing the fraction of relocated people as p ij (t) ≡ T ij (t) /x i (t), we can write the population evolution equation as where X is the (row) vector of the suburbs' population and P is the relocation matrix denoting the fractions of relocating people between each pair of suburbs, with the diagonal elements p jj denoting the fraction of non-relocating residents. The column sum for each row of the relocation matrix is equal to ∑ N j=1 p ij = 1, so P is a row-stochastic matrix [44]. The population evolution Equation (1) represents a simple Markov process, converging to a distinct stationary state X eq , which we identify as the equilibrium state. In equilibrium, the population of each suburb x i,eq does not change in time, such that X eq · P eq = X eq , where P eq ≡ lim t→∞ P(t). We assume that the urban area is a closed system with no external shocks and, thus, the relocation matrix does not change in time, i.e., P eq ≈ P(0) ≡ P.
Since matrix P is row-stochastic, it has a unit eigenvalue [44], and the vector X eq can be found as a left eigenvector of matrix P that corresponds to the unit eigenvalue. This eigenvector is unique (up to a constant multiplier) if some power of matrix P has strictly positive elements [44]. Our next step is to explicitly represent the relocation dynamics in both spatial and temporal terms. We decompose the relocation matrix according to the following structure: where is the share of people who relocate to a different suburb within a period of time. Such a decomposition is known as the mover-stayer model [45], which has been used to describe relocation phenomena in biology, economics, and social sciences [46][47][48][49].
Here, matrix H shows the relocation structure of those residents who moved to a different suburb (i.e., have not stayed in the same suburb). The matrix H shows the spatial structure of the system and characterises the variation in microscopic "attractiveness" between different suburbs. Without loss of generality, we assume that h ii = 0. Furthermore, the coefficient can be alternatively written as ≡ 1/τ, where τ is the characteristic relocation time. We will refer to it as the relocation frequency or the population mobility.
We next extend the simple decomposition (3), so that the population mobility has a more complex structure than a single frequency . Without loss of generality, we assume that the relocation dynamics are governed by a discrete set of relocation frequencies k , where k = 1, 2, . . . , K. In the context of the population structure, this would suggest that the urban population comprises several distinct groups, which differ in their mobility. There may exist a number of classifications that differentiate population groups by their mobility, based on their ownership status (renters and home-owners), family status (singles and families), and employment status (students, professionals, retirees). In this work, we abstract away from the specific nature of these groups, assuming only their existence.
Expanding the structure of the population mobility does not affect the spatial structure of the relocation, i.e., the matrix H. We therefore assume that the spatial structure of the relocation dynamics is the same for each population group. Although this is, in principle, a strong condition, we show in Section 4.3 that it does not affect the results in practice.
The equilibrium population of each component is obtained similarly to Equation (2), as X k,eq · P k = X k,eq .
We show in Appendix A.1 that the population of each component X k (t) converges to the equilibrium population structure where α k is the total fraction of the city population belonging to the component k, so that ∑ C k=1 α k = 1 and X eq is the total equilibrium population structure that is independent of k and α k . This is also illustrated in Figure A1. Thus, the total equilibrium X eq can be obtained using the full spatial matrix H, without the component-specific relocation matrices P k or even the component's fractions α k . This is very convenient, as the component structure of the population is not known a priori, and in practice, it is only matrix H that can be obtained from the data directly.

Results
In this section, we present the results of our framework for the Australian capital cities. First, in Section 3.1, we demonstrate that the model with homogeneous population fails to consistently describe the intra-urban migration dynamics, while a heterogeneous model resolves this issue. Second, in Section 3.2, we develop an analogy between intra-urban migration and diffusion. This allows us to interpret the heterogeneous dynamics of intraurban migration as diffusion in a multi-component fluid mixture. Third, in Section 3.3, we predict the equilibrium configuration of the considered cities.

Revealing the Two-Component Structure of Intra-Urban Evolution
We first analyse the human relocation flows in eight Australian Greater Capital Areas, which represent populated metropolitan areas with diverse cultural and economic activities. The model is calibrated using the data from the Australian Census [50], which are reported as the migration flows T ij between each pair of suburbs within one year and five years, denoted as T ij;1Y and T ij;5Y , respectively. This suggests the natural choice for the time step as one year. The data are available for two census years, 2011 and 2016, with the migration counted backwards. Intra-suburb migration is not considered in this analysis. The data resolution we use, Statistical Area 2, is the finest in the Australian Census for which the migration data are available.
A naive approach suggests calculating the one-year migration matrix directly as This, however, produces results that are inconsistent with the five-year migration data. Indeed, the five-year migration matrix is, by definition, p ij; 1Y = P 5 ij , where P is the matrix of migration rates p ij; 1Y , and P 5 ij stands for the element in row i and column j of matrix P 5 . The five-year migration flow extrapolated from the one-year migration flow isT ij; 5Y = p ij; 5Y x i (t). Comparing the five-year population obtained from actual migration flow ∑ j =i T ij; 5Y (2016) with the five-year population obtained from the predicted migration flow ∑ j =iTij; 5Y (2016), as shown in Figure 1, we observe a systematic disagreement: The predicted numbers of movers are consistently higher than the actual numbers. In particular, in all Greater Capital Areas, the average share of people who do not change their place of residence within one year is about 0.87-0.92. The analogous share within five years is about 0.7-0.78, while the predicted one is approximately 0.9 5 ≈ 0.6, as shown in Table 1.
We note that there exist, in principle, several alternative ways to calibrate the singlegroup model. They, however, give the same result: The single-group model is not capable of explaining the five-year migration patterns from the one-year migration patterns. We refer to Appendix C.2 for the details of these calibrations.
In order to resolve this problem, we extend the model, assuming that the population comprises two groups instead of one, while staying with the general framework (11). Each group is characterised by its own relocation frequency, 1 and 2 ; in general, these differ from each other. Furthermore, we restrict ourselves to the case where the population share of each group, α 1 ≡ α and α 2 = 1 − α, is the same across all suburbs in the short term and is equal to the total population share. If α is different for each suburb, the model will have an excessive number of parameters, which may improve the goodness of fit, but will reduce the calibration robustness. We point out that the number of population groups with distinct relocation frequencies does not have to be equal to two: It simply has to differ from one. This is a crucial departure from a homogeneous population model that is not capable of explaining the actual relocation dynamics.  Table A1. Within this framework, we deduce two parameters, 1 and 2 , from the data sets described above. There exist multiple estimation algorithms for similar models with a parametric structure of the matrices (see, e.g., [47,51,52] for more details). Here, we use a simple calibration technique by selecting parameters 1 , 2 , and α without specifying a parametric functional form for the elements of relocation matrix H.
We calculate migration flows T ij as the sum of two components: where P k = (1 − k )I + k H, P 5 k ij stands for the element in row i and column j of the matrix P 5 k . Matrix H is estimated as follows: where T ij; 1Y (2016) is the number of people that migrated from suburb i to suburb j within the one-year period of 2015-2016. Relaxation rates k and α can be found from the conditions where s 1Y is the average share of stayers within one year, and s 5Y is the average share of stayers within a five-years period, which are calculated from the Census data. The actual values of s 1Y and s 5Y for the Australian Capital Areas are such that the solution to Equation (8) exists and is unique. It should also be noted that the available data do not allow us to calibrate the value of α, and thus, it has to be fixed beforehand. The general question of existence and uniqueness of the solution to Equation (8) with respect to α and k is discussed in Section 4.2 in more detail. The magnitudes of the migration outflow for each suburb predicted by this model for α = 0.9 are plotted against their actual magnitudes in Figure 1 (numbers of movers, ∑ j =i T ij; 5Y (2016)). It is evident that the values predicted by the two-component model are in a stronger agreement with the actual data than the predictions obtained by the one-component model.
The same analysis can be performed with the 2011 migration data, and the corresponding predictions are shown in Figure A2. Here, we again observe that the naive model produces a systematic bias in its predictions, while the two-component model provides a good fit to the data. From this comparison, we can conclude that the described methodology predicts the migration flows with a high precision once the systematic bias produced by the one-component model is eliminated.

Intra-Urban Migration as Diffusion
In this section, we describe intra-urban migration as an irreversible process of diffusion. We first follow the general description in Section 2, building the analogy for a general multicomponent fluid. Next, we illustrate the analogy for a specific case of intra-urban migration in Sydney using the results from Section 3.1.
The difference U(t) ≡ X(t) − X eq between the actual population at time t and the equilibrium population X eq shows how far the system is away from equilibrium. Introducing the rate of population change Q(t) ≡ X(t + 1) − X(t) as the difference between two subsequent time steps and using Equation (2), we can rewrite the population evolution in Equation (1) as where L ≡ P − I and I are the identity matrix. We view Equation (9) as the central expression underlying the analogy between intra-urban migration and diffusion. Indeed, if we consider two reservoirs with different fluid concentrations connected by a thin channel, there will exist a flow of fluid through that channel until these concentrations equilibrate.
The rate of the concentration change in each of the reservoirs is proportional to the difference between the current concentration and the equilibrium concentration, with the proportionality coefficient being related to the diffusion coefficient, and following the same dependency as Equation (9). In general, Equation (9) has the form of a typical transport equation in non-equilibrium thermodynamics [53], which describes the irreversible evolution of a thermodynamic system. For a closed system, this corresponds to a relaxation phenomenon, with U(t) being the driving force, which drives the system towards equilibrium, Q(t) being the rate of material change. Furthermore, L is the matrix of transport coefficients, or simply the transport matrix, comprising the transport coefficients between each pair of suburbs. The transport matrix determines how fast the system relaxes towards equilibrium. In the context of urban dynamics, the irreversibility is ensured by constancy of the relocation matrix P: The fractions of residents migrating between two suburbs, p ij , remain constant during relaxation (while the flows T ij and populations x i keep changing).
In other words, once the equilibrium is reached, there is no driving force to reverse the relocation dynamics (1). The transport coefficients are central in describing the irreversible evolution of a thermodynamic system. Similarly, the knowledge of the transport matrix is central in predicting the relocation dynamics in an urban system. An essential property of a transport coefficient in a physical system is that, in a closed system, it does not depend explicitly on time. This reflects the microscopic reversibility of molecular motion. While this principle does not exist a priori for an urban system, demanding that the transport matrix L (and, therefore, the relocation matrix P ≡ L + I) does not change with time, it indicates the microscopic reversibility of intra-urban relocation. We will see below that this assumption is supported by actual data, helping us to derive the transport matrix from the Census data on relocation.
Substituting decomposition (3), we obtain so that the transport matrix factorises into two terms. The temporal term, the coefficient , shows the speed of relaxation towards equilibrium and characterises the rate of the system irreversibility. The spatial term, the matrix H, shows the spatial distribution of a migration potential. We next point out that the population groups introduced above correspond to distinct components in a multi-component fluid mixture. Indeed, extending Equation (10) to multiple population groups, we can write Here, the temporal term k is different for each population group, corresponding to a fluid component with a distinct relaxation rate. In contrast, the spatial term H is the same for all components, corresponding to an external potential field. Writing the transport equation for each component separately, we obtain: where L k is the component-specific transport matrix defined by Equation (11), while With this analogy, the overall dynamics of intra-urban evolution follow a profile of diffusive relaxation. Specifically, at large t, the asymptotic decay of the driving force should be exponential, with exponent λ k < 1 being proportional to the second largest eigenvalue of matrix H. This implies that near equilibrium (when the values of U k (t) are small), the rate Q k (t) is asymptotically proportional to the driving force: We illustrate the long-term relaxation dynamics for a specific case of intra-urban migration in Sydney. As revealed in the previous subsection, there exist two population groups in Sydney, which correspond to two components in a fluid mixture. Figure 2 shows the relaxation profiles for these components. Particularly, the left panel shows that that the driving force decays exponentially with time, as expected. Furthermore, the right panel shows that the near-equilibrium rate of relaxation is linearly proportional to the driving force, according to Equation (13). This is, indeed, in agreement with the framework of linear irreversible thermodynamics [53], where the near-equilibrium relaxation rate is linearly proportional to the driving force, and the coefficient of proportionality is characterised by the second eigenvalue of the relocation matrix. Figure 2. Exponential convergence of U k (t) for each of the population groups, k = 1, 2: (A) log U k (t) is plotted against time step t; (B) Q k (t) is plotted against U k (t) (thick dotted curves) and the tangential lines with the slope 1 − λ k (solid straight lines), where λ k is the second eigenvalue of the group relocation matrix. For illustration purposes, both U k (t) and Q k (t) are normalised by the total number of residents, α k x, in the corresponding group.

Predicting Equilibrium Population Distribution
We next build a long-term forecast for the spatial structure of the Australian cities. In doing so, we assume that the current migration flows remain stable in the following years. As has been mentioned above, the equilibrium structure X eq is independent of α, These results reveal that the equilibrium states of three out of eight capital cities (Sydney, Melbourne, and Perth) are more spread out compared with their current structure (shown in Figure A3). The equilibrium structures of Brisbane, Adelaide, and Darwin are similar to the current ones, while the structures of Hobart and Canberra are more compact than the current one. These qualitative observations can be quantified by the spreading index method [26,54,55], which determines the degree of polycentricity and dispersalas opposed to monocentricity and compactness-of the city. The values of the spreading index (calculated for both actual configurations of the Australian capital cities and for the predicted ones) are shown in Table 2. It is remarkable that our long-term prediction is independent of α, , and even the number of heterogeneous components (which can be larger than two in reality) as long as the spatial migration pattern described by matrix H is shared by the entire population.

Robustness Evaluation
In the previous section, we showed that our non-equilibrium framework of diffusive intra-urban relaxation explains the short-term migration data and is able to provide longterm predictions. The important element of the model is the assumption that the population comprises multiple components with different relocation frequencies, which, in the context of our framework, correspond to different relaxation rates. In this section, we investigate the robustness of this claim, analysing the extent of its applicability. In particular, in Section 4.2, we study the two-component model, arguing that this case is sufficient to consistently describe the short-term migration. In Sections 4.3 and 4.4, we explore the sensitivity of the equilibrium configuration to the spatial migration patterns (captured by matrix H), varying between different population components. In Section 4.3, we do this for an abstract city with extreme migration patterns, while in Section 4.4, we extend this analysis to a specific case (Sydney).

Heterogeneity of the Population
In Section 3.1, we reported that the baseline one-component model systematically predicts higher rates of five-year migration than those observed in reality. Here, we show that having a homogeneous population mobility can not produce consistent migration predictions; hence, the population has to be heterogeneous with respect to its mobility. The population, which is comprised of two groups with two distinct relocation frequencies, is a minimal realisation of such heterogeneity.
One may expect that the five-year relocation rate observed in reality is lower than the one predicted from the one-year relocation rate due to the low mobility of recently relocated people. Indeed, if an individual has recently moved into a new home, there might be not much incentive for them to move further relatively quickly. To represent this constraint, we assume that people do not relocate for τ years after their last relocation (immobility assumption), and calculate the share of those who have not relocated for different values of this parameter, τ. Figure 4 compares how this share changes in time for homogeneous populations (full solution of this model is provided in Appendix A.3) and a two-component population without an immobility assumption. It is evident from the figure that the immobility assumption does not improve the model. Indeed, no value of the immobility period τ can match the five-year relocation flow. In such a setting, the presence of the residents that do not relocate within a certain period implies a higher relocation frequency for those who do. As a result, the five-year relocation rate predicted by the homogeneous model with the immobility assumption is higher than that produced by the simple homogeneous model (the baseline model), leading to an outcome opposite to the one that was expected.  (6)). All models are calibrated to the Sydney relocation data. All solid curves pass through the actual one-year relocation rate s 1Y = 0.91, but go well below the corresponding five-year value, s 5Y = 0.73. This analysis shows that a lower mobility of the recently relocated people cannot be a valid explanation for the lower actual five-year relocation rate. Therefore, we conclude that having a homogeneous population comprised of only one component is not sufficient to make the model consistent with the data.

Solution Space of the Two-Component Model
In Section 3.1, we calculated 1 and 2 using the average shares of stayers within a one-year period (s 1Y ), within a five-year period (s 5Y ), and with conditions (8). The values of 1 and 2 depend on α, which is not known without specifying the nature of the groups. This, however, does not affect the possibility of splitting the population into two groups and obtaining consistent predictions of the five-year migration patterns from the one-year migration patterns.
The non-linear system of algebraic Equation (8) allows one to calculate the relocation frequencies 1 and 2 for a given composition α. It consists of two equations while containing three unknown variables ( 1 , 2 , and α), and for a given α, it can have up to five real roots. Some of the roots may not belong to the range from 0 to 1 and, therefore, cannot be valid solutions for 1 and 2 . Figure 5 demonstrates that, depending on s 1Y , s 5Y , and α, there exist two, one, or zero solutions. Furthermore, it is evident from Figure 5 that, for any s 5Y ranging from s 5 1Y to s 1Y , there exists at least one solution for the pair 1 and 2 . This means that it is always possible to calibrate the model (8) as long as s 5Y ≥ s 5 1Y (the other inequality, s 5Y ≤ s 1Y , holds automatically), which is the case for the actual migration data.
Although it is not feasible to estimate α directly from the current data set, it would be possible to do so with a longer record of internal migration. In particular, if we also knew people's places of residence 10 years ago, 15 years ago, etc., we would be able to determine the value of α, which predicts the share of movers and stayers with a higher precision. This idea is illustrated in Figure 6, which demonstrates the stayer share values predicted by the two-component model. Parameters 1 and 2 are calibrated to s 1Y = 0.91 and s 5Y = 0.73 (Sydney values are used as an example). The values of α vary from 0.7 to 0.95 (the solution of (8) exists and is unique if 0.05 ≤ α ≤ 0.33 and 0.67 ≤ α ≤ 0.95, but the solutions for α = a and α = 1 − a are equivalent due to symmetry). For the 30-year horizon, the predictions for the share vary from 0.2 (if α = 0.95) to 0.6 (if α = 0.7). This means that the two-component model can consistently calibrate a larger variety of data than the one-component model. If, however, the actual structure of the population is more complex, e.g., is made of a larger number of components, the two-component model would not be able to adequately account for the corresponding data.

Component-Specific Relocation Matrix
As has been shown previously, heterogeneity in mobility rates k does not affect the equilibrium structure of the city as long as all groups have the same relocation matrix H. This assumption, however, may not be valid if we do not observe H for each group directly. This might be important, as there exists empirical evidence that different social groups (such as renters and mortgagors) have different migration patterns [27]. Thus, it is important to assess the possibility for the matrices H k to be component-specific, or in other words, heterogeneous.
If matrices H k are heterogeneous, the equilibrium population structure X eq is no longer independent of the compositions α k and relocation rates k . In particular, matrices H k may have different stationary vectors X k,eq , in which case it is impossible to estimate the stationary population structure unless matrices H k are observed directly. The equilibrium structure X eq = ∑ C k=1 X k,eq depends on the individual matrices H k and cannot be expressed via the aggregated relocation matrixĤ = ∑ C k=1 α k H k . In this section, we show that the structure calculated using the aggregated matrix H can still give a reasonable approximation for the stationary population structure X eq , even if the individual vectors X k,eq differ drastically. To demonstrate this, we consider two artificial examples. In the first example, the population components 1 and 2 generate the migration flows with opposite directions. In the second example, there are two groups of suburbs (A and B), and the members of component 1 always relocate to the suburbs within A, while the members of component 2 always relocate to suburbs within B. These two examples are considered for a linear toy city comprising 99 suburbs. All suburbs are located along a line, such that suburb 50 is the "central" one.
In the first example, matrix H 1 consists of elements h 1; ij given by: where d i = |i − 50| is the distance from i to the "central" suburb (suburb 50), β = 0.1.
Elements of H 2 are defined as follows: This form of H 1 and H 2 means that the members of component 1 prefer to relocate to more central suburbs (that are close to suburb 50), while the members of component 2 relocate to the peripheral suburbs (which are far from suburb 50) more frequently. The equilibrium population distributions X 1,eq and X 2,eq are displayed in Figure 7 (left column). As one might have anticipated, the population of component 1 forms a monocentric structure around the "central" suburb, while the population of component 2 predominantly inhabits the peripheral suburbs. The corresponding total population structure X eq and its approximationX eq obtained from matrixĤ = ∑ C k=1 α k H k are shown in the right column. In all three cases, (A) α = 0.1, (B) α = 0.5, and (C) α = 0.9, the actual population structures X eq (green bars) lie very close to the corresponding approximationsX eq (red solid line).
In the second example, we fill columns 26-74 of matrix H 1 with positive random numbers, and the other columns are filled with zeros. In contrast, to fill matrix H 2 , we assign zero values to columns 26-74 and positive random numbers to columns 1-25 and 75-99. Each row in both matrices is normalised so that its elements sum to one.
It is natural to anticipate that, in the equilibrium, all members of component 2 will live in suburbs 26-74, while the members of component 1 will live in suburbs 1-25 and 75-99 (left column in Figure 8; cases A, B, and C correspond to α = 0.1, 0.5, 0.9, respectively). In the right column of Figure 8, we again observe that, regardless of α, the actual equilibrium structure X eq (green bars) does not deviate significantly from its approximationX eq (red solid line). For all values of α, the approximationsX eq obtained from the observable matrixĤ (red solid line) are almost indistinguishable from the ground-truth equilibria X eq (green bars).

Sydney Case Study
To demonstrate the robustness of the results presented in Section 3.3 with respect to the heterogeneity of relocation patterns, we extend this analysis to the Greater Sydney Capital Area. In a real city, the migration matrices H 1 and H 2 are not normally observed separately. Moreover, these matrices cannot be assigned arbitrarily, as they need to be consistent with the actual migration data. In particular, following the procedure suggested in Section 4.3, the component-specific matrices have to be defined such that αH 1 + (1 − α)H 2 = H, with H 1 accounting for the relocations flowing primarily into central districts, and H 2 corresponding to the relocations flowing primarily into the peripheral areas.
To accomplish this task, we choose a distance threshold d, and select suburb groups A(d) and B(d) so that A(d) contains only suburbs with the distance to the central business district being less than d, while B(d) contains the rest of the suburbs. Next, we assume that when relocating, the members of component 1 almost always choose suburbs from set A(d), while group 2 members choose suburbs from set B(d). Finally, we calibrate H 1 and H 2 to actual relocation data denoting a i ≡ ∑ j:d j ≤d h ij in each row i, and define elements h 1; ij of H 1 as follows: if a i ≤ α, and if a i > α. The elements of H 2 are then given by: In other words, all members of the first component move to areas inside A(d) and all members of the second component move to areas inside B(d), but the total share a i of people from i who move to suburbs inside A(d) might differ from α. If a i ≤ α, we assume that all the people who relocate from i to A(d) belong to the first component, and so does the proportion (α − a i )/(1 − a i ) of the people migrating to the other suburbs, while the rest of i's residents belong to the second component. Conversely, if a i > α, we assume that only the proportion α/a i of those who relocate to A(d) belong to the first component, while others belong to the second component.
It is easy to see that, in that case, αH 1 + (1 − α)H 2 = H, and that all elements h 1; ij and h 2; ij are positive; in each row i, we have ∑ N j=1 h 1; ij = 1 and ∑ N j=1 h 2; ij = 1, which is required by construction. The resulting equilibrium structure of the population density is shown in Figure 9 for α = 0.9, d = 22 km (median distance to the central business district). Similarly to the previous examples, the approximatedX eq ( Figure 9D) does not differ significantly from the actual value X eq (Figure 9C), although X 1,eq and X 2,eq do differ drastically ( Figure 9A,B) .
From these examples, we can conclude that the heterogeneity in matrices H k has a limited effect on the long-term population structure, and it is possible to obtain an accurate prediction by using only the aggregate relocation matrix H = ∑ C i=1 α k H k . Figure 9. Equilibrium population density in Sydney for the case of heterogeneous relocation matrices H k . (A) First group's equilibrium structure X 1,eq ; (B) second group's equilibrium structure X 2,eq ; (C) total equilibrium structure X eq ; (D) approximationX eq obtained from the overall relocation matrix H. The scale bar in the lower-left corner indicates a distance equivalent to 20 km.

Conclusions
We have introduced a diffusive migration framework that describes intra-urban migration as an irreversible evolution of the urban population. The results have been tested with residential relocation data available from the Australian Census for eight Greater Capital areas over 10 years.
Using this framework, we were able to explain the medium-term (five years) migration patterns using the short-term (one year) migration patterns. We have shown that this is possible to achieve only if the population is not homogeneous and has an internal structure.
In particular, such a population should be comprised of at least two components, with each component having a distinct relocation frequency. Such a relocation frequency corresponds to a particular relaxation time of the component.
This heterogeneity of migration frequencies has an intuitive interpretation. For example, the group of residents that migrate more often can be interpreted as renters (who are less attached to their places of residence, and are relatively free to change them as soon as they identify a better option), and the other group could be interpreted as homeowners (for whom it may be more problematic to change the place of residence due to the transaction costs, peculiarities of the housing market, and individual circumstances).
Using this diffusive migration framework, we produced a long-term prediction for the Australian capital cities' structures based on the short-term migration data, with only an assumption about the temporal stability of the migration rates. According to our predictions, the largest capital cities (Sydney, Melbourne, and Perth) are moving towards more spread-out configurations, while Hobart and Canberra exhibit a more compact structure in equilibrium. The other capitals, Brisbane, Adelaide, and Darwin, are likely to preserve their current configurations in the long run. These results are consistent with those of previous studies predicting the possibility of polycentric transition in Sydney and Melbourne [19,27].
Our predictions are robust with respect to the composition of the migration components, as well as the possible heterogeneity of their relocation patterns, both dynamically and spatially. In particular, we have analytically shown that the long-run equilibrium is independent of the size of each community and their relocation rates. The robustness with respect to the spatial heterogeneity of relocation has been shown numerically through an abstract illustrative example. In this example, the relocation communities have opposite preferences regarding their destinations: Members of the first group prefer central districts, while their counterparts prefer the peripheral ones.
The temporal stability of the migration flows is a crucial element of our long-term analysis. Despite being consistent within the period of observation (2006-2016), they may be affected by multiple factors in future: Human migration is a complicated nonlinear process involving multiple interdependent factors, often leading to various phase transitions and critical phenomena [19,20,22,27,43,56,57]. However, the relocation data may contain some unique features that are not captured in other static human mobility and land-use data (e.g., [13,19,24,27,58]). Thus, we believe that the proposed dynamic framework for intra-urban migration, which enables robust long-term predictions, offers a principled approach to modelling out-of-equilibrium urban development.

Conflicts of Interest:
The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Appendix A. Mathematical Derivations
Appendix A.1. Convergence to the Same Equilibrium in (11) Proposition A1. If, in Equation (11), each component k has a unique equilibrium structure X k,eq , then it is independent of k and has the form X k,eq = α k X eq , where X eq can be found as a solution of the left eigenvector of the matrix H, which corresponds to the unit eigenvalue.
Proof. By the definition of equilibrium X k,eq , we obtain: which is equivalent to where θ = (0, 0, · · · , 0) is a zero row vector. This implies that all X k,eq are eigenvectors of H that correspond to the unit eigenvalue. Exact values of X k,eq can be found from the constraint on total population (each component k has total population α k x). If we choose a vector v satisfying (A3) and whose elements sum to the total population x (∑ N i=1 v i = x), it is easy to see that X k,eq = α k v for each component k, and the total equilibrium structure is given by X eq = ∑ C k=1 X k,eq = v.
Appendix A.2. Derivation of Equation (8) Equation (8) is a consequence of (6), which can be rewritten as where P 1Y and P 5Y are the aggregate one-year and five-year migration rate matrices, respectively. To obtain the first equation in (8), we directly apply decomposition (3) and equate the diagonal elements (without loss of generality, we set h ii = 0). To obtain the second one, we take into account that elements h ij are small and have an order of magnitude of (1/N) 1. This implies that all small powers of H have elements with the same order of magnitude, which can be disregarded. In particular, the elements of H 2 are h ij; 2Y = ∑ N k=1 h ik; 1Y h kj; 1Y ∼ 1/N, the elements of H 3 are h ij; 3Y = ∑ N k=1 h ik; 1Y h kj; 2Y ∼ 1/N, and, hence, the elements of H 5 are h ij; 5Y = ∑ N k=1 h ik; 1Y h kj; 4Y ∼ 1/N (symbol ∼ means "has order of magnitude of"). Hence, the diagonal elements of both matrices P 5 1 and P 5 2 are approximately ≈ (1 − k ) 5 .

Appendix A.3. Migration Model with Memory
We let a i denote the share of people that relocated i years ago (i = 1, 2, . . . τ). We use a 0 to represent the share of people who relocated more than τ years ago or did not relocate at all. Stationary values of a 0 , a 1 , . . . a τ can be found from the following equations: a 1 = a 0 , a 2 = a 1 , . . . a τ = a τ−1 , a 0 + a 1 + · · · + a τ = 1.

(A5)
The solution is: To match the one-year relocation rate, the mobility parameter needs to satisfy the following equation: a 0 (1 − ) + a 1 + · · · + a τ = s 1Y .
The solution can be found numerically for any value of τ, but the case τ = 1 admits an analytical expression for : The model dynamics are now described by an extended transition matrix that describes both migration between suburbs and "mobility states" (from state "migrated i years ago" to state "migrated i + 1 years ago" or state 0, where mobility is unrestricted): MatrixP contains matrices H, (1 − )I and identity matrix I as its N × N blocks. The transition matrix for a five-year period can be calculated using the fifth power ofP. The share of stayers in each suburb k can be calculated as follows: wherep 5 kk;ij is the k-th diagonal element in block i, j of matrixP 5 .  where y actual i is the actual migration outflow from suburb i and y predicted i is its predicted counterpart. The corresponding values are shown in Table A1. In this section, we compare the five-year predictions of the single-component model obtained for different ways of calibration. The share of stayers within a five-year period is calculated based on the fifth power of the one-year migration rate matrix P. We use the following estimates for P: