Multiscale mobility patterns and the restriction of human movement

2.1. The directed graph of baseline UK mobility: quasi-reversibility and commuting travel patterns

By using pre-lockdown mobility data (average of 45 days before 10 March 2020), we construct a strongly connected directed graph G with weighted adjacency matrix A ≠ A^T (figure 1a). The N = 3125 nodes of this graph correspond to geographic tiles (width between 4.8 and 6.1 km, see electronic supplementary material, figure S1a), and the directed edges have weights A_ij corresponding to the average daily number of inter-tile trips from tile i to j (see §4.1 for the notion of ‘trip’ within the Facebook dataset, and some of its caveats). The total average number of daily inter-tile trips is 2 475 527, compared with 10 416 968 intra-tile trips. The matrix A is very sparse, with 99.7% of its entries equal to zero, i.e. there are no direct trips registered between the overwhelming majority of tile pairs. Furthermore, the non-zero edge weights are highly heterogeneous, ranging from 1.4 to 6709 daily trips (average 72.3, coefficient of variation 2.8), underscoring the large variability in trip frequency across the UK.

• Download figure
• Open in new tab
• Download PowerPoint

To assess the directionality of the baseline network, we first compute the pairwise relative asymmetry (PRA) for each pair of tiles ij,

defined for pairs where A_ij + A_ji > 0. The distribution of the PRA_ij (electronic supplementary material, figure S3) shows that 25% of the tile pairs have PRA ≥ 0.23, a substantial asymmetry, including 3226 one-way connections (8.64% of the total) with PRA = 1. It is thus helpful to use network analysis tools that can deal with directed graphs [31].

A natural strategy for the analysis of directed graphs is to employ a diffusive process on the graph to reveal important properties of the network, such as node centrality [32,33] or graph substructures [27], while respecting edge directionality. We consider a discrete-time random walk on graph G defined in terms of the N × N transition probability matrix M:

where $D_{\mathrm{out}}^{+}$ is the pseudo-inverse of D_out = diag(d_out), the diagonal matrix of out-strengths d_out = A 1_N. A key property of the random walk is its stationary distribution π, a 1 × N node vector defined through the equationThe component π_i is a measure of centrality (or importance) of node i; a high value of π_i means that the random walk on G is expected to visit node i often at stationarity [34]. This is equivalent to PageRank without teleportation [35]. As expected, the centralities π_i are highly correlated (R² = 0.97) with another node centrality measure, the out-strengths d_out,i. Interestingly, the centralities π_i are also correlated with the intra-tile mobility (R² = 0.83), a measure that is not used in the computation of π (see electronic supplementary material, figure S2 for details). Therefore, urban areas display high centrality due to the concentration of human mobility in those areas (figure 1a).

A random walk on a directed graph might still not display strong directionality at equilibrium. This is our finding here: the random walk defined by M fulfils approximately the detailed balance condition,

where $\Pi =\text{diag}(\mathit{\pi })$ and ‖ · ‖_F denotes the Frobenius norm (see electronic supplementary material, figure S3 and [27,36] for a more in-depth discussion). The random walk for our mobility graph is therefore close to being time reversible at equilibrium [34], so that the probability of following a particular trajectory from node i to j is almost equal to the probability of going back on the same trajectory from j to i. This property coincides with our intuition that most journeys in the mobility network are linked to commuting travel patterns.

2.2. Unsupervised community detection reveals intrinsic multiscale structure in the baseline mobility data

To extract the inherent scales in the UK mobility data, we apply multiscale community detection to the baseline directed network. We use Markov stability (MS) [25], a methodology that reveals intrinsic, robust graph partitions across all scales through a random walk Laplacian that simulates individual travel on the network based on the observed average daily trips. As random walkers explore the network, they remain contained within small subgraphs (communities) at shorter times and then spill over onto larger communities at longer times. This definition of communities (and partitions) in terms of random walks makes the MS framework generally applicable to a wide range of network topologies, including directed networks, in contrast to standard hierarchical community detection algorithms [37–39]. MS uses an optimization to identify graph communities in which the flow of random walkers is contained over extended periods, uncovering a sequence of robust graph partitions of increasing coarseness (regulated by the Markov scale s). In the context of our mobility network, this set of partitions captures intrinsic scales of human mobility present in the data. See [25–28,40,41] and §4.2 for details of the methodology.

Figure 1b summarizes the MS analysis for our network, which was carried out with the PyGenStability Python package [42]. We find nine robust MS partitions H(s_i) at different levels of resolution (s₁, …, s₉) from fine to coarse, which comprise flow communities at different scales of human mobility (see electronic supplementary material, table S1 and figure S4 for further statistics and visualizations). Figure 1c shows that these data-driven flow communities correspond to geographic areas, even though our data only contain relational mobility flows without explicit geographic information. Furthermore, the nine partitions have a strong quasi-hierarchical structure, which is not imposed a priori by our graph partitioning method (see electronic supplementary material, figures S4 and S5). The obtained partitions thus reflect an inherent multiscale structure in the patterns of UK human mobility.

2.3. Comparing the intrinsic mobility scales at baseline with administrative nomenclature of territorial units for statistics regions

Next, we compare the MS partitions with NUTS regions, administrative and geographic regions defined at three hierarchical levels: NUTS1 build upon NUTS2 in turn consisting of NUTS3 regions. In the UK, the 174 NUTS3 regions represent counties and groups of unitary authorities; the 40 NUTS2 regions are groups of counties; and the 12 NUTS1 regions correspond to England regions, plus Scotland, Wales and Northern Ireland as whole nations (see electronic supplementary material, table S2 for further statistics). Our baseline data cover 170 NUTS3 regions, where the missing four are sparsely populated regions in the Scottish Highlands and Islands. The NUTS regions serve as a standard reference point for policy-making, and served to inform regionalized responses to COVID-19 in England (e.g. lockdowns in the North of England were applied to local authorities that form the NUTS2 region of Greater Manchester, Lancashire and West Yorkshire [43]). Comparing the data-driven MS partitions with NUTS regions is thus meant to explore to what degree administrative regions capture the patterns of mobility at the different scales and potential mismatches thereof.

In figure 2, we use the normalized variation of information (NVI; equation (4.3)) to evaluate the similarity of each of the three NUTS levels to the MS partitions at all scales. The best match of each NUTS level (as given by the minimum of NVI) is close to one of the robust partitions: NUTS3 corresponds closely to H(s₁); NUTS2 to H(s₂) and NUTS1 to H(s₄). Hence, these three MS partitions of the mobility network capture the fine, medium and coarse scales in the UK, yet with some significant deviations from the administrative NUTS divisions. For instance, Greater London is separated from the rest of the South East at the level of NUTS1 regions, whereas the whole South East of England forms one flow community in partition H(s₄). Similarly, the south of Wales is connected strongly via flows to the South West of England in partition H(s₄), which is not reflected in the NUTS1 regions. On the other hand, the correspondence between NUTS3 regions and the fine MS partition H(s₁) is strongest (lower value of NVI), with fewer such discrepancies between administrative and flow communities.

• Download figure
• Open in new tab
• Download PowerPoint

To evaluate further how the MS partitions capture the patterns of mobility, we compute two measures: the coverage $\mathcal{C}$ (equation (4.5)) (i.e. the ratio of mobility that remains within communities relative to the total mobility) and the average nodal containment NC (equation (4.7)) (i.e. the ratio of the outflow from each node that remains within its community relative to the total outflow from that node, then averaged over all nodes). High values of these measures (normalized between 0 and 100%) indicate that mobility flows are captured within the boundaries of the communities of the partition.

Table 1 shows that MS partitions are substantially better at reflecting baseline mobility than NUTS divisions since they have higher values for both average $\mathcal{C}$ and NC measures at all scales and especially at the finer scales. We have also evaluated both measures at a local level. Electronic supplementary material, figure S6a shows that the median of the coverage of individual communities ${\mathcal{C}}_k$ (equation (4.4)) is significantly higher for MS partitions (as compared with NUTS) for the fine and medium scales (p < 0.0001, Mann–Whitney). Electronic supplementary material, figure S6b shows that the median of the nodal containment of individual nodes NC_i (equation (4.6)) is also significantly higher for MS partitions relative to NUTS regions at fine, medium and coarse levels (p < 0.001, Mann–Whitney). Indeed, the maps in figure 2b–d show that NC_i(MS) > NC_i(NUTS) in regions where the administrative NUTS boundaries cut through conurbations or closely connected towns or cities. A prominent example is Greater London, where the NUTS2 regions split areas in Central and East London that are tightly linked and thus captured better by the medium MS partition H(s₂), and, similarly, the NUTS1 region of Greater London does not include its wider commuter belt that is naturally captured by the coarse MS partition H(s₄). Similar commuter belt effects are observed, e.g. on the medium level for Birmingham, and on the fine level for Plymouth, which has associated flows across the Cornwall–Devon boundary.

2.4. Comparing the fine mobility scale at baseline with labour-related travel to work areas

We next compare the MS partitions with TTWAs, a different geography that divides the UK into 228 local labour markets computed from 2011 Census data recording place of residency and place of work [29]. Our baseline network has mobility data for 197 of the 228 TTWAs, with missing areas in rural areas of Scotland, Wales and the North of England (see figure 3a and further statistics in electronic supplementary material, table S3).

• Download figure
• Open in new tab
• Download PowerPoint

The TTWAs are intended to reflect local labour markets, and the ensuing commuting between home and place of work, and are thus expected to be linked to small scales. Indeed, as measured by the NVI, the TTWA division is most similar to the NUTS3 level of all NUTS divisions (see electronic supplementary material) and, consistently, most similar to the fine-scale MS partition, H(s₁). Reassuringly, H(s₁) is more similar to TTWA than to NUTS3, since both the TTWA division and our MS partitions are data driven with a basis in mobility patterns. Yet, there are local discrepancies, e.g. H(s₁) combines multiple TTWAs into a single cluster in Cornwall or Northern Ireland, whereas the single TTWA in Greater London corresponds to several smaller communities in H(s₁) (figure 3a).

As mentioned earlier, we evaluate how the TTWA division captures the patterns of mobility. We find that the coverage $\mathcal{C}(\mathrm{TTWA})=95.9\mathrm{\%}$ is higher than for both NUTS3 and H(s₁) (table 1), but the median of the coverage of individual communities ${\mathcal{C}}_k$, a local measure of coverage, is not significantly higher for the TTWAs relatively to H(s₁) (figure 3b). Furthermore, the average nodal containment $\mathrm{NC}(\mathrm{TTWA})=81.3\mathrm{\%}$ is lower than H(s₁) (table 1), and its local version shows that the median of the nodal containment of individual nodes NC_i is significantly lower for TTWA (p < 0.0001, Mann–Whitney, figure 3c). Hence, the fine MS partition H(s₁) captures better the mobility patterns in our baseline data than the TTWA division. This can be explained by potential changes in commuting patterns since the 2011 Census data on which the TTWAs are based, and by the fact that Facebook mobility data also includes trips for leisure, commercial and other activities beyond commuting to work.

2.5. The contraction of the UK mobility under lockdown and its relation to the baseline mobility multiscale network

The first nationwide COVID-19 lockdown in the UK was imposed on 24 March 2020, instructing the British public to stay at home except for limited purposes. Over the following months, restrictions were gradually eased to allow pupils to return to school (1 June 2020 in England but 17 August 2020 in Scotland), businesses to reopen (non-essential shops reopened on 13 June 2020 in England but 13 July in Scotland) and people to travel more freely for leisure purposes (13 May 2020 in England but 8 July 2020 in Scotland) [44]. We have analysed the response to these restrictions using the time-dependent Facebook Movement maps [24] from 10 March to 18 July 2020 (131 days or 18 weeks). We construct mobility networks for each day, G(d), and week, G(w), defined on the same nodes (i.e. tiles) as the baseline network G (see §4.1).

Figure 4a shows the temporal change of the number of trips (intra-tile, inter-tile and total) relative to 10 March 2020, the first day of our study period. It is interesting to note that the decrease in mobility was already taking hold rapidly from 10 March, two weeks before the official enforcement of the lockdown. We find that whilst the total number of trips remained largely unchanged throughout the period, the number of inter-tile trips decreased sharply to approximately 25% of the initial value, followed by a steady increase towards levels of approximately 50% at the end of the study period in mid-July 2020. Conversely, the number of intra-tile trips increased to a maximum of 130% after lockdown before decreasing steadily to approximately 105% by mid-July 2020. Therefore, lockdown induced a redistribution from inter-tile to intra-tile trips as a result of a reduction in commuting and long-distance travel, with mobility reverting to local neighbourhoods.

• Download figure
• Open in new tab
• Download PowerPoint

The observed contraction of human mobility towards local neighbourhoods is consistent with the multiscale structure that was already present in the baseline mobility network pre-lockdown. The coverage $\mathcal{C}$ of all MS partitions increased over lockdown, with larger relative improvement for the finer scales (figure 4b and electronic supplementary material, figure F8). The surge in coverage induced by lockdown, which then decays towards its pre-lockdown value, can be modelled with a simple linear model under an external stimulus αe^−λt. The relative change from the initial value is then given by [45],

from which we estimate the amplitude (α) and characteristic time (1/λ) of the external stimulus, as well as the characteristic recovery time (1/β) of the system towards its pre-stimulus value (see §4 for details). Figure 4b shows the fits of $\mathrm{\Delta }\mathcal{C}(t)$ (dashed lines) with estimated parameters in table 2. The fine MS partition exhibits the largest relative increase $\mathrm{\Delta }\mathcal{C}(t)$ peaking at approximately 6%. The medium and coarse partitions peak at approximately 1 and $0.1$%, respectively. This is also captured by the values of α and indicates that during lockdown people reverted to local mobility neighbourhoods already present in pre-lockdown patterns. The adaptation to the new COVID-19 situation and the pre-announcement of lockdown occurs quickly (over a characteristic time of 1/λ ∼ 2 weeks), signifying that adoption was fast and was already in progress before the official start date of lockdown. Mobility patterns then returned towards pre-lockdown values over longer time scales 1/β between 16.4 weeks (fine scale) and 20.9 weeks (coarse scale) reflecting a slow re-adaptation following the new situation and loosening of restrictions.

To highlight the local differences in the temporal response to the lockdown, figure 5 shows the parameters of the temporal fits for the community coverages $\mathrm{\Delta }{\mathcal{C}}_k$ for all the communities in the fine-scale MS partition. We observe that urban centres like London, Birmingham, Liverpool or Manchester experience the strongest changes in the fine-scale coverage (high values of α) yet with faster recovery times (low values of 1/β). Conversely, rural areas, which were already more constrained to local communities pre-lockdown, exhibit smaller but long-lived effects in the coverage at the local level. Our method also captures divergent trends across the different nations of the UK. For example, Scottish regions show longer time scales of recovery than most regions in England, consistent with the fact that Scotland maintained more stringent lockdown restrictions for a longer time [44], e.g. domestic travel restrictions were eliminated in Scotland only on 8 July 2020 and schools reopened on 17 August 2020, in contrast to 13 May and 1 June 2020 in England, respectively.

• Download figure
• Open in new tab
• Download PowerPoint

4.1. Mobility data and network construction

4.2. Multiscale community detection with Markov stability analysis

Here, we provide a brief outline of the MS framework. For a fuller description, see the electronic supplementary material and in-depth treatments, including extensions to other types of graph processes, in references [25–27,40].

Consider a weighted and directed graph G with adjacency matrix A. Let $L=I-D_{\mathrm{out}}^{+}A$ denote the random walk Laplacian matrix, where I is the identity matrix, and $D_{\mathrm{out}}^{+}$ denotes the pseudo-inverse of the diagonal out-degree matrix. The matrix L defines a continuous-time Markov process on G governed by the diffusive dynamics

where p(r) is a 1 × N node vector of probabilities, and r is the Markov scale. The solution to this equation is given by p(r) = p(0)exp( − Lr), and the matrix exponential defines transition probabilities of the Markov process (see electronic supplementary material). This process converges to a stationary distribution π given by πL = 0.

The goal of MS is to obtain partitions of the graph into c(r) communities such that the probability flow described by (4.1) is optimally contained within the communities as a function of r. MS solves this problem by maximizing the following function:

where $\Pi =\text{diag}(\mathit{\pi })$, and the matrix H(r) is a N × c(r) partition indicator matrix with H(r)_ij = 1 if i is part of community j, and 0 otherwise. We thus obtain a series of optimized partitions over the Markov scales described by the matrices H(r). The scales are more naturally described in log scale, so we redefine the Markov scale as s = log₁₀(r). The optimization (4.2) is carried out using the Louvain algorithm [49] through the implementation in the PyGenStability python package [42].

4.3. Measures of flow containment: coverage and nodal containment

Consider the adjacency matrix A of the mobility graph G and a N × c indicator matrix H for a partition of G into c communities. Let us also define $\tilde{A}$, the adjacency matrix of the graph with self-loops that contains the intra-tile flows on the diagonal. Then $F=H^T\tilde{A}H$ is the c × c lumped adjacency matrix where the element ${(H^T\tilde{A}H)}_{kl}$ corresponds to the mobility flow from community k to community l. The coverage of community k, ${\mathcal{C}}_k(H)$, is defined as follows:

where ${\hat{D}}^{+}$ is the pseudo-inverse of $\hat{D}=\text{diag}(\hat{\mathit{d}}),$ where $\hat{\mathit{d}}=F\hspace{0.17em}{\mathbf{1}}_c$. ${\mathcal{C}}_k(H)$ can be interpreted as the probability of the lumped Markov process to remain in state k; hence, high values of ${\mathcal{C}}_k(H)$ indicate that community k covers well the flows emerging from the community.

The coverage of a partition $\mathcal{C}(H)$ is standard and is defined as the ratio of flows contained within communities by the total amount of flow [55]. It is easy to see that this is given by the weighted average

High values of $\mathcal{C}(H)$ indicate that mobility flows are contained well within the communities of the partition and movement across different communities is limited.

The nodal containment NC_i of node i quantifies the proportion of flow emerging from i that is contained within its community in a partition H,

where C_i is the community of node i and d_i = (A 1_N)_i. Large values of NC_i indicate that the mobility flows emerging from node i are largely contained within its assigned community, indicating a good node assignment. Hence, NC_i measures the containment of flows from a node-centred perspective.

To obtain a partition-level measure, we define the average nodal containment NC(H)

where N is the number of nodes.

4.4. Response to an exponentially decaying shock

The response of a variable x(t) to a shock can be modelled as a linear ordinary differential equation under a stimulus R(t),

where 1/β is the characteristic relaxation time of the system, and we assume here an exponentially decaying external stimulus R(t) : = α e^−λt, with amplitude α ≥ 0 and characteristic decay time 1/λ. The solution of (4.8) is given by [45]We use the Levenberg–Marquardt algorithm [56] implemented in the LMFIT [57] python package to fit the activation response function x(t) to a set of n data points $({\tilde{t}}_i,{\tilde{x}}_i)$ by minimizing the sum of squaresto determine parameter estimates $\hat{\alpha }$, $\hat{\beta }$ and $\hat{\lambda }$. Confidence intervals are obtained from an F-test [58].