Multi-scale spatio-temporal analysis of human mobility

The recent availability of digital traces generated by phone calls and online logins has significantly increased the scientific understanding of human mobility. Until now, however, limited data resolution and coverage have hindered a coherent description of human displacements across different spatial and temporal scales. Here, we characterise mobility behaviour across several orders of magnitude by analysing ∼850 individuals’ digital traces sampled every ∼16 seconds for 25 months with ∼10 meters spatial resolution. We show that the distributions of distances and waiting times between consecutive locations are best described by log-normal and gamma distributions, respectively, and that natural time-scales emerge from the regularity of human mobility. We point out that log-normal distributions also characterise the patterns of discovery of new places, implying that they are not a simple consequence of the routine of modern life.


Introduction
Characterising the statistical properties of individual trajectories is necessary to understand the underlying dynamics of human mobility and design reliable predictive models. A trajectory consists of displacements between locations and pauses at locations, where the individual stops and spends time (Fig 1). Thus, the distribution of waiting times (or pause durations), Δt, between movements and the distribution of distances, Δr, travelled between pauses are often used to quantitatively assess the dynamics of human mobility. For example, specific probability distributions of distances and waiting times characterise different types of diffusion processes. Thanks to the recent availability of data used as proxy for human trajectories including mobile phone call records (CDR), location based social networks (LBSN) data, and GPS trajectories of vehicles, the characteristic distributions of distances and waiting times between consecutive locations have been widely investigated. There is no agreement, however, on which distribution best describes these empirical datasets.
Pioneer studies, based on CDR [1,2] and banknote records [3], found that the distribution of displacement Δr is well approximated by a power-law, P(Δr) * Δr −β , (or 'Lévy distribution' [4], as typically 1 < β < 3), and that an exponential cut-off in the distribution may control boundary effects [2]. These findings were confirmed by studies based on GPS trajectories of a1111111111 a1111111111 a1111111111 a1111111111 a1111111111  First, the datasets considered have different spatial resolution and coverage, and few studies have so far considered the whole range of displacements occurring between *10 and 10 7 m (10000 km) (Fig 2). Our analysis suggests that constraining the analysis to a specific distance range may result in different interpretations of the distributions. Another difference concerns the temporal sampling in the datasets analysed so far. Uneven sampling typical of CDR and LBSN data (i) does not allow to distinguish phases of displacement and pause, since individuals could be active also while transiting between locations, and (ii) may fail to capture patterns other than regular ones [31,32], because individuals' voice-call/SMS/data activity may be higher in certain preferred locations. Finally, studies focusing on displacements effectuated using one or several specific transportation modality (private car [24,33], taxi [20], public transportation [34], or walk [7]) capture only a specific aspect of human mobility behaviour.
In this paper, we analyse mobility patterns of *850 individuals involved in the Copenhagen Network Study experiment for over 2 years [35]. Individual trajectories are determined combining GPS and Wi-Fi scans data resulting in a spatial resolution of *10 m, and even sampling every *16 s. Trajectories span more than *10 7 m. Previous studies with comparable spatial coverage (Fig 2) relied on single-transportation modality data [8], unevenly sampled data [16], or small samples (32 individuals in Ref. [5]). To our knowledge, the Copenhagen Network Study data has the best combination of spatio-temporal resolution and sample size among the datasets analysed in the literature to date (see Methods).

Results
We consider an individual to be pausing when he/she spends at least 10 consecutive minutes in the same location, and moving in the complementary case. In the following, we refer to locations as places where individuals pause. The distribution of displacements is robust with respect to variations of the pausing parameter (see S1 File for the results obtained with 15 and 20 minutes pausing).
We start by considering the three distributions most frequently reported in the literature (Table 1), namely • The log-normal distribution of a random variable x, with parameters σ and μ, defined for σ > 0 and x > 0, with probability density function: • The Pareto distribution (i.e. power-law) of a random variable x, with parameter β, defined for x ! 1, and β > 1, with probability density function: • The exponential distribution of a random variable x, with parameter λ, where x ! 0, and λ > 0, with probability density function: In Eq (2) the probability density can be shifted by x 0 and/or scaled by s, as P(x) is identically equivalent to P(y)/s, with y ¼ ðxÀ x 0 Þ s . In Eqs (1) and (3), P(x) is identically equivalent to P(y), with y = (x − x 0 ). In this work, the shift (x 0 ) and scale (s) parameters are considered as additional parameters to take into account the data resolution. With few exceptions, the results presented below hold also imposing no shift, x 0 = 0 (see S1 File). Note also that Pareto distributions with exponential cut-off (or truncated Pareto) are considered below (see also Table 1).

Distribution of displacements
We start our analysis by investigating the distribution of displacements between consecutive stop-locations P(Δr). First, we consider the overall distribution of the displacements Δr using all available data (851 individuals over 25 months). We find that P(Δr) is best described by a log-normal distribution (Eq (1)) with parameters μ = 6.78 ± 0.07 and σ = 2.45 ± 0.04, which maximises Akaike Information Criterion (see Methods)-among the three models considered -with Akaike weight *1 (Fig 3, see also S1 File).
Second, we investigate if this results holds also for sub-samples of the entire dataset. We bootstrap data 1000 times for samples of 200 and 100 individuals, and we verify that the best distribution is log-normal for all samples, and the average parameters inferred through the bootstrapping procedure are consistent with the parameters found for the entire dataset (see S1 File). In fact, the errors on the value of the parameters reported above are computed by bootstrapping data for samples of 100 randomly selected individuals. This analysis ensures homogeneity within the population considered, and takes into account also that often smaller sample sizes were analysed in previous literature. Third, we zoom in to the individual level. We find that the individual distribution of displacements is best described by a log-normal function for 96.2% of individuals. The best distribution is the Pareto distribution for 1.4%, and exponential for the remaining 2.4%. However, the number of data points per individual tend to be significantly lower in group of individuals exhibiting Pareto or exponential distributions, so that one should be cautious in interpreting the observed deviations from a log-normal distribution. Fig 4 reports the histogram of the individual μ parameters for the 96.2% of the population that is best described by a log-normal distribution, along with three examples of individual distributions.
Finally, we look at large Δr in order to compare our results with precedent studies relying on data with larger spatial resolution. We find that limiting the analysis to large values of Δr results in the selection of a Pareto distribution (Eq (2)). We identify the threshold ΔrÃ = 7420 m as the minimal resolution for which the best fit in ΔrÃ < Δr < 10 7 m is Pareto with coefficient β = 1.81 ± 0.03 and not log-normal. By bootstrapping 1000 times over samples of 100 individuals we find thatDrÃ ¼ 7488:3 AE 328:2 m. Thus, power-law distributions describe mobility behaviour only for large enough distances, while mobility patterns including distances smaller than 7420 m are better described by log-normal distributions.

Distribution of waiting times
We now analyse the distribution of waiting times between displacements. The best model describing the distribution of waiting times over all individuals is the log-normal distribution (Eq (1), Fig 5, see also S1 File), with parameters μ = −0.42 ± 0.04, σ = 2.14 ± 0.02. As above, errors are found by bootstrapping over samples of 100 individuals. Also, by bootstrapping we find that the log-normal distribution is the best descriptor for samples of 200 and 100 randomly selected individuals (see S1 File). As in the case of displacements, we find that restricting the analysis to large values of our observable Δt, and specifically considering only Δt > ΔtÃ = 13 h, results in the selection of the Pareto distribution (Eq (2), see Fig 5), with coefficient β = 1.44 ± 0.01. We find by averaging over 100 samples of 200 individuals that DtÃ ¼ 13:01 AE 0:12. Note that the log-normal distribution is selected as the best model also when the analysis is restricted to Δt < ΔtÃ. The distribution of waiting times shows also the existence of "natural time-scales" of human mobility. We detect local maxima of the distribution at 14.0, 39.3, 64.8, and 89.9 hours. Hence, 14 hours is the typical amount of time that students in the experiment spent home every day, in agreement with previous analyses on human mobility [23,25,26]. Other peaks appear for intervals Δt % 14 + n Á 24, with n = {2, 3. . .}, suggesting individuals spend several days at home. Notice also that the distribution we consider is limited to Δt < 5 days, an interval much shorter than the observation time-window (about 2 years), a fact that guarantees the absence of possible spurious effects [29]. This limit is imposed to control the cases in which students leave their phones home. The upper bound is arbitrarily set to 5 days; however, we have verified that results are consistent with respect to variations of this choice.

Distribution of displacements between discoveries
Log-normal features also characterise patterns of exploration. We consider the temporal sequence of stop-locations that individuals visit for the first time-in our observational window-and characterise the distributions of displacements between these 'discoveries'. We find that the distribution of distances between consecutive discoveries P(Δr) is best described as a log-normal distribution with parameters μ = 6.59 ± 0.02, σ = 1.99 ± 0.01, (Fig 6, see also S1 File). For Δr > 2800 m, the best model fitting the distribution of displacements is the Pareto distribution with coefficient β = 2.07 ± 0.02. This results are verified by bootstrapping (see S1 File).

Correlations between pauses and displacements
We further investigate the properties of individual trajectories by analysing the correlations between the distance Δr and the duration Δt disp characterising a displacement and the time Δt spent at destination. Fig 7A shows a positive correlation between Δr and Δt disp for Δr ≳ 300 m (p < 0.01). As Δr is the distance between the displacement origin and destination, the absence of correlation at short distances could be due to individuals not taking the fastest route. A positive correlation characterises also the distance Δr covered between origin and destination and the waiting time at destination for distances 30 m ≲ Δr ≲ 10 4 m (p < 0.01). Instead, the correlation is negative for distances larger than 5 × 10 4 m (Fig 7B). This could suggest that individuals break long trips with short pauses. We have verified that these results hold also when individuals' most important locations (typically including university and home) are removed from the trajectory, implying that these correlations are not dominated by daily commuting.

Further analysis: Selection of the best model among 68 distributions
In the previous sections we have restricted the analysis of the distributions of displacements and waiting times to the three functional forms that are most frequently found in the literature. We now repeat the selection procedure considering a list of 68 models (see S1 File for the list of distributions) in order to confirm the results described above.
The distributions of displacements and displacements between discoveries are best described by log-normal distributions also when the choice is extended to 68 models, and tails (respectively for Δr > ΔrÃ = 7420 m and Δr > ΔrÃ = 2800 m) are better modelled as generalised Pareto distribution, with form: where ξ is the parameters of the model, such that x ! 0 if ξ ! 0, and 0 x À 1 x if ξ < 0. The best model selected for the whole distribution of waiting time among the 68 models considered is a gamma distribution, defined for x 2 (0, 1), k > 0 and θ > 0 as: Although the gamma distribution is the best model for the distribution of waiting times (see S1 File for the result of the fit), the presence of natural scales could indicate that the whole distribution may be better described as the composition of several models.

Discussion
Using high resolution data we have characterised human mobility patterns across a wide range of scales. We have shown that both the distribution of displacements and waiting times between displacements are best described by a log-normal distribution. We found, however, that powerlaw distributions are selected as the best model when only large spatial or temporal scales are considered, thus explaining (at least partially) the disagreement between previous studies. We also showed that log-normal distributions characterise the distribution of displacements between discoveries, implying that this property is not a simple consequence of the stability of human mobility but a characteristic feature of human behaviour. Finally, we have shown that there exist correlations between displacements' length and the waiting time at destination.
The heavy tailed nature of human mobility has been attributed to various factors, including differences between individual trajectories [36], search optimisation [37][38][39][40], the hierarchical organisation of the streets network [41] and of the transportation system [6,24,42]. On the other hand log-normal distributions can result from multiplicative [43] and additive [44] processes and describe the inter-event time of different human activities such as writing emails, commenting/voting on online content [45] and creating friendship relations on online social networks [46]. Instead, the distribution of inter-event time in mobile-phone call communication activity can be described as the composition of power-laws [47][48][49], a feature attributed to the existence of characteristic scales in communication activity such as the time needed to answer a call, as well as the existence of circadian, weakly and monthly patterns. We also find clear signatures of circadian patterns, which could indicate that the whole distribution may be better described as the composition of several models. However, in our case the best description for times including Δt < Δt Ã is the gamma distribution, which thus is selected both when the whole range of scales is considered and when the analysis is restricted to short times.
Our results come from the analysis of a sample of *850 University students, which of course represent a very specific sample of the whole population. Nevertheless, it is worth noting that many statistical properties of CNS students mobility patterns are consistent with previous results, such as the distribution of the radius of gyration, the Zipf-like behaviour of individual locations frequency-rank plot, and the power-law tail of the distribution of displacements (β = 1.81 ± 0.03 vs. β = 1.75 ± 0.15 of [2]). Details are reported in Supplementary Information of [50].
While identifying the mechanism responsible for the observed mobility patterns is beyond the scope of the present article, we anticipate that a more complete spatio-temporal description of human mobility will help us develop better models of human mobility behaviour [24,51]. Our findings can also help the understanding of phenomena such as the spreading of epidemics at different spatial resolutions, since the nature of heterogeneous waiting times between displacements have a major impact on the spreading of diseases [52].

Data description and pre-processing
The Copenhagen Network Study data collection took place between September 2013 and February 2016 and involved 851 students of Technical University of Denmark (DTU) in Copenhagen. Data collection was approved by the Danish Data Protection Agency. All participants provided informed consent by filling an on-line consent form and all methods were performed in accordance with the relevant guidelines and regulations. Individual trajectories were inferred combining WiFi scans data and GPS scans data recorded on smartphones handed out to all participants. An anthropological field study included in the 2013 deployment of the experiment reported that participants did not alter their habits due to participation in the CNS experiment.
The WiFi scans data provides a time-series of wireless network scans performed by participants' mobile devices. Each record (i, t, SSID, BSSID, RSSI) indicates: • the participant identifier, i • the timestamp in seconds, t • the name of the wireless network scanned, SSID • the unique identifier of the access point (AP) providing access to the wireless network, BSSID • the signal strength in dBm, RSSI.
APs do not have geographical coordinates attached, but their position tend to be fixed. The geographical position of APs is estimated the procedure described in S1 File, which used participants' sequences of GPS scans to obtain APs locations and remove mobile APs. Then, we clustered geo-localised APs to "locations" using a graph-based approach. With our definition, a "location" is a connected component in the graph G d , where a link exists between two APs if their distance is smaller than a threshold d (see [50], SI for more details). Here, we present results obtained for d = 2 m. However, results are robust with respect to the choice of the threshold (see also [50]).
Throughout the experiment, participants' devices scanned for WiFi every Δt seconds. The median time between scans is between Δt M = 16 s and Δt M < 60 s for 90% of the population (see also [50], SI). Data was temporally aggregated in bins of length Δt = 60 s, since we focus here on the pauses between moves. If a participant visits more than one location within a timebin, we assign the location in which they spent the most time to that bin. Given our definition of location and the given time-binning, the median daily time coverage (the fraction of minutes/day that an individual's position is known, where the median is taken across all days) is included between 0.6 and 0.98 for 90% of the population.

Model selection
The best model is selected using Akaike weights [53]. First, we determine the best fit parameters for each of the models via Nelder-Mead numerical Likelihood maximisation [54] (maximisation is considered to fail if convergence with tolerance t = 0.0001 is not reached after 200 Á N iterations, where N is the length of the data). For each model m, we compute the Akaike Information Criterion: where L m is the maximum likelihood for the candidate model m, V m is the number of free parameters in the model, and n is the sample size. The AIC reaches its minimum value AIC min for the model that minimises the expected information loss. Thus, AIC rewards descriptive accuracy via the maximum likelihood and penalises models with large number of parameters. The Akaike w m (AIC) weight of a model m corresponds to its relative likelihood with respect to a set of possible models. Measuring the Akaike weights allows us to compare the descriptive power of several models.
For all distributions considered in this paper, we found one model mÃ such that w m Ã * 1 (which implies all the other models have Akaike weight very close to 0).

Figures
All figures were generated using Matplotlib [55] package (version 1.5.3) for Python.

Related work
We present here more detailed analysis of the literature discussed in the paper.