Inferring country-specific import risk of diseases from the world air transportation network

Disease propagation between countries strongly depends on their effective distance, a measure derived from the world air transportation network (WAN). It reduces the complex spreading patterns of a pandemic to a wave-like propagation from the outbreak country, establishing a linear relationship to the arrival time of the unmitigated spread of a disease. However, in the early stages of an outbreak, what concerns decision-makers in countries is understanding the relative risk of active cases arriving in their country—essentially, the likelihood that an active case boarding an airplane at the outbreak location will reach them. While there are data-fitted models available to estimate these risks, accurate mechanistic, parameter-free models are still lacking. Therefore, we introduce the ‘import risk’ model in this study, which defines import probabilities using the effective-distance framework. The model assumes that airline passengers are distributed along the shortest path tree that starts at the outbreak’s origin. In combination with a random walk, we account for all possible paths, thus inferring predominant connecting flights. Our model outperforms other mobility models, such as the radiation and gravity model with varying distance types, and it improves further if additional geographic information is included. The import risk model’s precision increases for countries with stronger connections within the WAN, and it reveals a geographic distance dependence that implies a pull- rather than a push-dynamic in the distribution process.


Introduction
The recent decades have seen a considerable increase in mobility: the world-wide number of passenger cars in use increased between 2006 and 2015 on average by about 4% each year reaching approximately 1 billion in 2015 (1), comparable with the yearly increase in the number of containers shipped over sea (2) and the number of globally scheduled air passenger even increased yearly between 2004 and 2019 by about 6% (3).In other words, the world gets more connected passengerwise on the small (cars) and largest (air traffic) scale, and in terms of the import and export of goods.This increased connectivity eases the distribution of anything related to goods and people, as shown for the distribution of over 400 invasive species by agricultural imports that is best predicted by the global trade network (4).A prime example of unwanted side effects of well-connected regions is the potential of a pandemic, accompanied by death, economic damage and the potential stigmatization of survivors, migrants and minorities (5)(6)(7).Already the first plague pandemic that started 541 in the Nile Delta of Egypt spread in 8 years across the territories (Mediterranean, Northern Europe and Near East) of 2 affected empires because of the intense commerce in the Roman Empire (6).Nowadays, the intensified exchange reduces the time until a pandemic reaches all parts of the world to months as for the 2009 H1N1 virus that spread from Mexico in 5 months to all continents (8,9) or the recent COVID-19 pandemic whose variants spread within a few months across the globe (10)(11)(12)(13).
The connection strength between world regions is only partly explained by their geographic proximity.Instead, due to historic geopolitical relations (14,15) pandemics spread rather along an effective distance that is derived from the world air transportation network (WAN) (16)(17)(18)(19), or, if applied on a smaller scale, also from other means of transportation (16,20).According to the effective distance, region B is closer to region A if the passenger flow from A to B is larger than to other destinations.An interesting extension is the multipath effective distance, that improves the arrival time prediction of a spreading disease by including all paths of a random walker on the WAN (17).The effective distance is regularly used to analyze the impact of mobility on the spread of diseases, as for example for MERS (21), Ebola (22), Zika (23) and most recently COVID-19 (20,(24)(25)(26).While it allows to qualitatively estimate the arrival time of a disease, its use is very limited for the description of import events of infected passengers from a specific source to target.However, these import events are highly relevant for political decisionmakers and to enable modelling predictions.
In this work, we describe these import events via the "import probability" (| ), which is equivalent to the origindestination (OD) matrix whose element    represents the number of trips from A to B, with the difference that the probability is normalized by all trips starting in A, i.e. (| ) =    /  .There exist mobility models that fit the OD matrix, i.e. a reference OD matrix is needed as for the gravity D R A F T model (27)(28)(29)(30)(31)).Yet, it can be extremely difficult to obtain the OD matrix and most often it is estimated by small surveys (32) or alongside a census (33), and even for the air transportation network with a booking system the OD is only an approximation since passengers increasingly book directly at the airlines (in 2015 30% of all Lufthansa flights were booked directly which increased to 52% in 2018 (34)) and not via the big GDS (global distribution systems) from which most ODestimates are derived (35,36).That means to exactly compute the air transportation OD matrix, bookings of all GDSs and about 900 airlines must be purchased/estimated and combined.Thus, models that do not rely on an existing reference OD matrix are important and those either assume an underlying decision process without integrating traffic information as the radiation model (37,38) or they apply a maximum entropy approach to distribute the unknown OD trips along possible routes of a known traffic network (30,39,40).However, none of the above approaches use the effective distance with its validated link to disease propagation and none is based on a mechanistic distribution process on a traffic network.The second point is crucial, since such a process provides us a mechanistic understanding of the observed pattern, and it enables us to study how its modification affects the passenger distribution, e.g.how a containment intervention along the distribution paths reduces the import probability of infected passenger.In this work, we introduce the import risk model that is based on a distribution process following the effective distance shortest path tree of the WAN combined with a random walker that explores all possible paths of the WAN from 2014.As a ground truth, we use the Global Transnational Mobility Dataset from 2014 (35) and investigate the discrepancy to the import risk and alternative mobility models as the gravity (27,31) and radiation model (38) through multiple comparison measures.We find that the import risk model outperforms the alternatives and only marginal improves if geodesic distance information is included.Finally, we analyze the quality of import probability estimation for certain countries and assess if and how the geodesic distance is encoded in the import risk estimate.

Relating the WAN, OD-probability and the effective distance
In this work, we introduce the import risk that estimates the probability of a passenger departing in airport A to end its journey at any other airport world-wide, including airports not directly connected to the origin airport.The estimation is based on the traffic flow of airplanes and the respective maximal passenger capacity between airports, a.k.a. the world air transportation network (WAN), provided by the Official Airline Guide (OAG) (42).This inference-problem is intriguing because it is much easier to monitor the origin and destination of airplanes, than of passengers with possibly multiple connecting flights until their final destination.In our study we use the WAN from 2014 (Fig. 1A) and compare the derived import probabilities to a reference dataset.The reference import probability is based on the Global Transna-tional Mobility Dataset (GTN) from 2014 (35,41), that uses a combination of an origin-final-destination dataset of a major global distribution system (GDS) and a tourism dataset from the World Tourism Organization (Fig. 1B, see Material and Methods for more details on the data).Before introducing the import risk model, we contrast the two datasets, introduce the effective distance (16) and quantify its potential as the base metric for our proposed model.By comparing the world air transportation network (WAN) with the country-specific reference import probability from the GTN (compare Fig. 1A, B), we see that the airports connected via direct links belong to countries that also have a high import probability.Due to physical constraints and logistic optimization, however, not all countries that have nonzero import probability are directly linked to an airport of the source country, but are reached via connecting flights instead.For the import probability, geodesic distance and population of the target country are useful estimates but fail in specific cases, e.g. the import to Italy is about 1.4 times larger than to Germany even if the latter is geographically closer to Canada and has a larger population.The effective distance is an alternative, network-based distance measure that does not rely solely on direct connections and geographic information (16)(17)(18)(19) but instead uses the passenger flow    to  from  and its relation to the outflow   in form of the transition probability    =    /  .Additionally, it increases the distance by a constant  0 for every connecting flight: The effective distance between airports without direct connection is the cumulative distance along the shortest path tree (SPT) derived from  eff , as illustrated for the largest Canadian airport (Toronto Pearson Airport, YYZ) in Fig. 1C.Former studies showed that the countries' arrival time of diseases depends linearly on their effective distance (16)(17)(18)(19).We show that the import probability also correlates with  eff (Fig. 1D), whereby the correlation is higher than for other distance measures (see Fig. S1).In fact, the import probability decays exponentially with effective distance (linear decay on a semi-log scale in Fig. 1D) which can be reproduced in a simplified model for a passenger that travels at a constant effective speed and has constant exit rate.Therefore, the effective distance seems to be a good representation of the underlying distribution process, and is a promising candidate for the base of our proposed import risk model, to directly estimate the import probability.

Import risk model
The idea behind the import risk model is a combination of (i) a random walk with exit probability and (ii) a distribution mechanism derived from the  eff SPT (Fig. 2).Applying a random walk is motivated by Iannelli et al. (17) who could improve the arrival-order prediction of  eff by including all possible paths.In the first step, we use the transition network representation of the WAN and let a random walker start at source  0 and after each step it either exits at the current node  with exit probability   or continues to walk.Let us define the walker's probability to continue walking to node  given it was at node  − 1 before and originally started in  0 by with  ,−1 as the transition probability from  − 1 to .Now the probability to walk along a path Γ starting at  0 and exiting at  is the probability to continue walking  ,  along each link (, ) that is part of the path times the exit probability of the final node where we omitted the explicit dependence on the source  0 .Our goal is to describe all possible paths the walker can take from  0 to .We will use the matrix S, whose elements are the probabilities to continue walking  ,  .The element (, ) of the product of the matrix with itself S • S = S 2 sums over all paths of length  = 2 that end in  and start in .We can now specify the probability of a walker to exit at  after taking all paths of length  to Finally, the import risk is the probability to exit at  given all paths of all lengths where we used the convergence of the geometric series with identity matrix I.
In the second step, we approximate the exit probability   ( 0 ) that we used above, but did not specify yet.Thereby, we assume that passengers start in source airport  0 , travel along the SPT and exit at node  with an exit-probability with  () as the population at airport  and Ω(|) as the set of all offspring nodes downstream of  on the SPT centered at source  0 .Thus, the exit probability at  is the ratio of the population at  to all of 's downstream nodes populations, including .We approximate the population at airport  with its outflow on the WAN  () =   and aggregate the import probabilities on country level by summing the targets and applying a weighted average on the source airports with the population as weight.
To clarify how additional information about the geographic distance between nodes influences  ∞ , we explore two variations of the import risk model: In the variation with "geodesic distance weighted" exit probability the populations in Eq. 6 are substituted with N (| 0 ) =  ()/ , 0 , where  , 0 is the geodesic distance between  and  0 .To control for increasing model complexity, we study the "effective distance weighted" exit probability, where N (| 0 ) =  ()/ eff (| 0 ), i.e. no geographic information is used, but the model structure is equivalent.
Alternative models.Many alternative models estimate the OD-matrix from which the import probability can be computed (30,31,37,38,(43)(44)(45)(46).Among those, the gravity (27) and the intervening opportunity (37,38) model are most widely used.A recent version of the latter is the radiation model (38).It is derived from a mechanistic decision process and is in consequence parameter free, and therefore similar to our model.However, it only requires information on the population density and does not integrate flight data.We compare our model to the gravity model with an exponential and power-law distance dependence and the radiation model (see Material and Methods for definitions).These models only use the outflow from the WAN to estimate the node's population and the geographic locations.To incorporate structural information of the WAN, the alternative models are also implemented with the geodesic path distance (the geodesic distance along the SPT) and the effective distance, i.e. there are in total nine alternative models: the radiation model, the gravity model with exponential and with power-law distance decaying function, and each implemented with geodesic, geodesic path and effective distance.
The exponents of the six gravity models are fitted to the reference import probability by assigning the best fitting exponent to each of the six comparison measure (Pearson correlation, root-mean-square error, common part of commuters, Kendalls rank correlation and the correlation and RMSE of the logarithmic measures, all defined in Material and Methods) and taking their mean value (see Figs. S2,S3).As comparison measures, we have chosen three measures that are related to the absolute error and three that are related to the relative error between estimate and reference.
Symmetry by returning visitors.Each of the twelve models provides an estimate for the import probability (| 0 ), from which the OD-matrix T can be computed by a multiplication with the respective source population  ( 0 ).By comparing the symmetry of T with the reference OD-matrix T, we find a much higher and qualitatively different symmetry in the reference data (see Sup.Note 2, Fig. S4).The high symmetry is likely due to visitors (family, business, tourism, etc.) that dominate the international travel.They return to their home-location after a limited period (47) and only the minority of the travelers are migrants, i.e. stay permanently at the destination.Interestingly, the import risk model has the highest symmetry, but is still less symmetric than the reference data by a factor of 4. Thus, before comparing the estimates in detail we correct the import probability estimates by symmetrizing their OD-matrix (by taking the symmetric part and recomputing the import probability, see Material and Methods and Sup.Note 2 for details).This correction can be seen as an alternative version of a doubly constrained model where normally the constraints on in-and out-flow are ensured by an iterative proportionate fitting (31).

Model comparison
In the following, we compare the import probability estimates to the reference data:(i) directly and via their medians to analyze potential systematic errors, (ii) by a collection of six different goodness of fit measures whereby we point out the rank and the relative performance of the single model and (iii) by a classification task that is highly relevant in a pandemic context when the countries at highest import risk needs to be known.
Qualitative Comparison.In Fig. 3 the import probability estimate (| 0 ) of each model is compared to the reference import probability p(| 0 ).The median of the probability estimates is for the gravity models in best agreement with the reference data if the effective distance is used (Fig. 3 first and second column), while the median values of the radiation and import risk model are less affected by the change in distance metric or the weighting by it (first and fourth column).All models overestimate the lowest median import probability (leftmost orange dot in Fig. 3), since a large portion of the observed reference import probability is zero due to a limited number of departing passengers.The overestimation of the median import probability persists up to (| 0 ) ≤ 10 −4 for the gravity and import risk models, but not for the gravity model with exponential distance decaying function and effective distance metric (Fig. 3I) which shows the best agreement of the median with the reference data.The radiation models (first column) systematically overestimates the highest import probabilities ((| 0 ) 10 −1 ) and in consequence underestimates the lower import probabilities.

Goodness of fit by multiple measures.
We compared each model with the reference import probability via the Pearson correlation, the root-mean-square error (RMSE), and the common part of commuters.These measures are more sensitive to strong links, i.e. large import probabilities, which is important when the emphasis is placed on the countries that are most likely to import passengers.However, if the focus is to get a fair comparison including all links, logarithmic versions of the above measures or rank correlations are more appropriate.Thus, we also quantify the agreement by the correlation and the RMSE of the logarithm of the measures and by Kendall's rank correlation.The three import risk model variations outperform the other models in all but one measure, illustrated by their highest ranks, whereby the variation with geodesic distance weighted exit probability performs best (Fig. 4A).The import risk models are followed by the two effective-distance based gravity models, the other models don't show a consistently high rank for all 6 measures, but are rather homogeneously scattered in the lower half.This model categorization also holds for the relative performance of the models (Fig. 4B), with linear scaling of values in between (see Eq. 21).In contrast to the ranks, the median relative performance improves considerably if the effective distance is used in the gravity models, while there is only a marginal difference in between import risk models.The only measure where the import risk models are outperformed by the gravity models with effective distance is the logRMSE (Figs.S5, S6).It is expected from the gravity models' good agreement in median import probability with the reference data over wide ranges and the overestimation of low import probability by the import risk model.This overestimation can be reduced by model-modifications that introduce parameters favoring the exit at nodes with largepopulations (for details, see Sup.Note 3 and Figs.S7,S8).However, we refrain from adding complexity to the model, since its generic nature is its key aspect.
Classification of ten top risk countries.In a pandemic context, it is of specific interest to identify the countries with the highest import probability.We analyzed how well the twelve proxy models can classify, if a country is among the ten countries with highest import probability.Again, the import risk models outperform the other models and the one with geodesic distance-weighted exit probabilities is the top predictor with a sensitivity of 71.1% (Fig. 5D).All effective
e  distance-based models have a high sensitivity ( 65%), including the radiation model with 66.8% that had the lowest relative performance and second-lowest mean rank (Fig. 5I-K).For these high import probabilities, the import risk models now outperform the other models also in terms of RMSE and logRMSE, i.e. the 10 countries at highest risk are not only classified best by the import risk model, but also quantitatively assessed best.

Import risk of countries and regions
Having quantified the performance of the import risk model, we now focus on (i) country specific differences in its prediction quality, (ii) possible limitations due to no concept of administrative units (e.g.countries) whose airports are more interconnected and (iii) how the geodesic distance is encoded in the import risk model, i.e. how a distance dependence emerges from WAN information only.
Country specific performance.In the import risk approach, we assume minimal knowledge of the system, i.e.only the WAN is known.Consequently, we differentiate countries only via their network properties, one of which is the degree of a node, or more precisely the node strength, since the WAN is a weighted network.It is the simplest metric that is also easily adjustable for the country-level perspective.For the country level the node strength directly corre-sponds to the flow out of country This country characteristic represents the potential of a country to dominate the structure of the network, since flows from small-outflow countries are diluted by large-outflow countries.From an ecological point of view, the outflow is strongly correlated with the gross domestic product of a country (Fig. S13).The correlation (logcorr) between the logarithms of the import risk  ∞ and the reference import probability p∞ improves with the outflow of the source country (Fig. 6), as illustrated by Great Britain (GB) as the country with the largest outflow in the WAN and Eritrea (ER) as one of the countries with the lowest outflow.The prediction improvement with the country's outflow suggests that the WAN is dominated by large-outflow countries and therefore predictions worsen for countries with lower WAN outflow.However, the prediction improvement is also present in model alternatives that do not use WAN information at all (e.g.gravity with geodesic distance, Fig. S12).We rule the explanation out that the alternative models show this improvement due to preferential fitting of strong links -and therefore of large-outflow countries -since the models are fitted to the reference data by their import probabilities which ensures equal weighting among countries.It rather suggests that the mobility behavior in low outflow regions is different, also supported by the sudden performance saturation for countries with a WAN outflow of   10 6 (Figs.6,S12).Possibly, their passenger distribution is constrained by additional factors and is limited to the regions in proximity.
There are clear exceptions where the import risk estimation is worse compared to outbreak countries with a similar WAN outflow, as Australia (AU), Israel (IL) and Macao (MO).These countries are connected due to historical relations to specific regions that are either not in their direct neighborhood (European countries for AU and IL) or that are more important than the bare neighborhood would suggest, as Macao that is a special administrative region of China.
For Macao the import risk to China is underestimated which consequently overestimates the import to other countries, and for AU and IL Europe is underestimated which overestimates other regions (Fig. 6).The exceptions AU, IL and MO illustrate that not all information is fully encoded in the WAN and therefore can not be extracted by the import risk model.Another concept that is missing in our methodological approach is the idea of a country or another administrative unit, i.e. we treat airport pairs equally irrespective of their countryaffiliation.Since we know the international flights leaving a specific country from the WAN, we can run a self-consistency analysis, i.e. without the need of reference import probability data.We can estimate the outflow leaving the country  by the import risk model by If we compare it to   the WAN flow out of country  (see Eq. 7), it turns out that the import risk model systematically overestimates the flow out of a country (Fig. S9A).In fact, the relative error increases with the number of airports belonging to the country (Fig. S9B).This overestimation could be due to the lack of a country concept in the import risk model, another explanation is the overestimation of the population in the catchment area of the airport by its outflow (the outflow erroneously counts transit passengers to the population).However, we can easily correct for this overestimation on country-level analysis, by normalizing the airport population such that the WAN country outflow is recovered.

Geodesic distance dependence.
The import risk model estimates import probabilities without explicit geodesicdistance information (excluding the variant with distance weighted exit probability).Since classical models have proven distance to be a good predictor for human mobility, we assume that it is encoded in the WAN structure and by consequence in the import risk estimate (48).For illustrative clarity, we aggregate the import risk on the twenty-two world regions and find that the import risks to a single target decreases with the geodesic distance to the sources in a powerlaw like manner (Fig. 7A,B and Fig. S11).If we switch the viewpoint and analyze the distance-dependence from a single source to all target regions (Fig. 7D,E), the dependence is less in agreement with a power-law fit  ∞ =  •  −   (Fig. 7C).This is surprising, since the import risk is computed via a source-centric view (by computing the exit prob-ability from the shortest path tree originating at each source), which suggest that the distance dependence should be best from one source to its possible targets.A possible explanation is that each target has its own attractiveness independent of the source region, i.e. that the distribution dynamic resembles more a pull-than a push-dynamic.Indeed, we find that the fitted exponent  from the power-law fit decreases the larger the WAN flow out of the target region is, which can serve as a proxy for the attractiveness of a region (Fig. 7F).In other words, the more attractive a region, the larger the import risks from more distant source regions.The fitted exponent  has a high rank correlation with  ( Kendall = 0.89), i.e. also the coefficient is dependent on the attractiveness of the region.

Discussion and Conclusion
Motivated by the import probability's strong dependence on the effective distance, we implemented the import risk model based on the effective distance shortest path tree's exit probability in combination with a random walk on the WAN.Thus, inferring the trip distribution of passengers from the traffic network of their transport vehicle (WAN).By comparing our parameter free model with variations of established mobility models, we find that it outperforms the alternatives in all but a single comparison measure, where the two parameterfitted gravity models with effective distance are performing best.The import risk model is the most accurate in determining countries with the highest import probability, showcasing its importance for epidemic related problem.However, it systematically overestimates low import probabilities and its performance worsens for countries with a passenger outflow below a million per year.Despite the lack of any explicit geodesic distance information, the import risk model recovers a geodesic distance dependence.This is more distinct from all sources to a single target than vice versa, which we connect to a target specific attractiveness estimated by its node strength, i.e. the target's passenger outflow.

D R A F T
The only measure where the gravity models with effective distance outperform the import risk models is the logRMSE.This is likely due to their good agreement over wide ranges of the import probability (Fig. 3I,J).The import risk model performs poorly with respect to the logRMSE due to its systematic overestimation of low import probabilities.Note, that the second parameter free model, the radiation model, systematically underestimates low import probabilities in the same way as import risk does.It is expected since deviation from the assumptions can not be corrected by any parameter adjustment.We identified several ways to reduce the import risk's overestimation of low import probabilities by introducing an additional parameter that scales the population of the respective airport, changes the exit probability along the shortest path tree or only these of the terminal nodes (for details, see Sup.Note 3 and Figs.S7,S8).We conclude that modifications that increases the probability to exit at airports/nodes with large populations reduces the overestimation.However, we leave a possible extension of our model for future studies and highlight that it outperformed the other models in all correlation measures, illustrating its high potential.
We corrected the import probability by the symmetrization of the respective OD-matrices which corresponds to a specific form of a doubly-constrained model.Normally, the constraints only ensure that the out-and inflow of each location corresponds to the observations (31,46,49), in contrast, we assume that both equal each other because of returning visitors and that in consequence the OD-matrix is symmetric.We repeated the model comparison without the correction: it reduced the agreement with the reference data for all but five of the seventy-two model-measures combinations (Fig. S6), which is in agreement with previous studies that report a better performance of doubly constrained models (49).Importantly, the import risk model still outperforms the other models if the import probability estimates are not corrected (compare Figs. 4, S10).We found that without providing any geodesic distance information to the import risk model, a distance dependence is recovered that is stronger for import probabilities to a single target, than from a single source, even if the import probability is computed from a source-centric view.Since the WAN is spatially embedded and has a network dimension of three (48), its connections reflect up to a certain degree the characteristics of the embedding space, explaining that the import risk recovers distance dependence in general.That distance is better predicting in the target-centric view aligns well with a previous study wherein a target-specific humanmobility model can collapse mobility data to multiple targets by assigning each target a specific attractiveness that is proportional to the target's population (45).The import risk model predictions worsen for countries with a small outflow on the WAN, and since the country's WAN outflow is proportional to its gross domestic product, the model performs less good for countries with a lower GDP, i.e. small population and/or low to middle income countries.This is problematic, since our model infers OD information (costly in terms of monitoring) from a low-cost monitoring of traffic flow, and therefore is especially interesting for regions with limited resources.However, we find that the model alternatives (gravity, radiation) also perform poorly for lowoutflow countries and that the passenger distribution of the latter is most likely constrained by the GDP and thus limited to the target-regions in effective proximity.To circumvent this problem, one could aggregate neighboring low-outflow countries until the conglomerate crosses the outflow threshold of   = 10 6 above which we observe a performance saturation (Figs.6,S12).Of course, this compromise comes with a lower spatial resolution and we emphasize the need for future research in this direction.
We have quantified the model performance for the world air transportation network, but it can be applied to other modes D R A F T  Fig. 7. Import risk aggregated on regional level "to target" vs. "from source" and its geodesic distance dependence.The geodesic distance between regions predicts the import risk  ∞ to a single target from all sources (A, B) better than from a single source to all targets (D, E) as can be seen by the p-values (C) of the power law fit  ∞ () =  •  − that is illustrated for each selected examples by a grey line ( A, B, D, E).The fitted exponent  of the import risk to a single target decreases with the respective regional WAN flow out of the target region (F), i.e. the more connected a region, the weaker the import risk decays with distance.The dashed horizontal lines show the average import risk of a single target (A, B) or a single source (D, E).The color of the dots corresponds to the depicted world regions (right).
of transportation, e.g.subway system, cars, buses, trains, and future research will show if there are certain conditions that need to be fulfilled for its application.Also, the crude estimation of the travelling population in the catchment area of an airport by the respective outflow can be improved, since we do not take into consideration the inflated role of hubs and the missing information about transit passengers.The simple framework that only relies on the traffic network is appealing, but in certain scenarios its prediction can be refined by using information about the GDP, Gini-coefficient or population density.We limited our comparison to the parameter-free radiation model and the fitted gravity model, consequently missing their promising variations or other alternatives (30,31,49,50).However, the gravity model is widely applied and has been shown to perform equally well (50) or better than alternatives (49).There are exceptions, e.g. an iterative computation of a gravity-like model outperforms the common gravity model in cases where the complete mobility network is not available (29).Additionally, the radiation model outperforms the gravity model for long-distance connections (49).
Still, the simplicity of the gravity model and its adaptability by parameter adjustment make it a strong counterpart.The model alternatives make use of the WAN-structure information by using the effective distance as done in e.g.Ren et al. (51) where the radiation model with time-distance was better than the travel-distance on the road network to predict the traffic on each link.Analogously, we also found that the effective distance (that relates to the arrival time of diseases) is better than the geodesic path-distance in predicting import probabilities.
The import risk model is fundamentally different from classic approaches that estimate OD trips from traffic data, because the latter find the OD trips that best reproduce the traffic data (28, 30, 39, 40), while our model runs a distribution process on the traffic data network.Thus, our model is mechanistic, while the classic approaches either fit and require the knowledge of the reference trip data (28,30) or are based on the assumption that the trip distribution across the links follows the maximum entropy principle, i.e. the OD trips are considered as most likely that can be realized by the largest number of microstates (39,40).Note that maximum entropy approaches require an estimation of routes and their alternatives between each OD pair, while we allow all routes to be taken by the random walker.Our model is in consequence -to our best knowledge -unique in its mechanistic nature, since it allows studying the modifications of it's underlying distribution process, for example in the form of a containment strategy meant to slow or restrict a pandemic.A straight forward implementation could be the test of a fraction of passengers   ≤ 1 at every transit airport , which corresponds to reducing the probability to continue walking of an infected passenger (Eq.2) to . .] one would allow for a varying testing capacity between the airports.

Material and Methods
Data sources.The WAN provided by OAG (Official Airline Guide) (42) contains the number of flights and the respective maximum seat capacity  ,  between airports  and  aggregated for the year 2014.The reference import probability p ( |) = T / T is based on the "Global Transnational Mobility Dataset" (35,41) that assigns the number of trips in 2014 T from country  to  worldwide by combining the world air transportation origin-final-destination data set from the company SABRE, and cross-boarder visits with an overnight stay from the UNWTO (World Tourism Organization).Thus, p (|) not only represents the mobility via air travel but also via other means (see, road, rail).However, air travel dominates long distance trips which makes it a fair reference set of the air transportation origin-final-destination matrix.For details on how the data sets were combined, see Sup.Note 1.
Alternative models.The gravity model states that the number of trips between regions  and  increase with their population sizes (  and   ) and decrease with distance with  () as a function that grows monotonically with distance , most often chosen as either a power-law  () =   or an exponential  (  ) =   .In the radiation model, the trips from  to  depend on their respective population sizes   ,   (or other measures as job opportunities) and on the number of people   that are in a circle with radius   centered around location  including   and   : The import probability of both models is computed by normalizing the trips with respect to the source-region Trip-symmetrization. We correct the import probability via symmetrizing the OD-matrix by (i) compute the estimated OD-matrix from the import probability estimate, (ii) correct it by computing its symmetric part and (iii) compute the corresponding corrected import probability via By going through these steps, the asymmetry is reduced heavily but still persists.Thus, we repeat steps (i) till (iii) until  (3) ( |), which returns for all models a comparable asymmetry in mean and median to the reference data (see Sup.Note 2 for details).
Comparison measures.We compare the import probability models with the reference data via the Pearson correlation with E [ ] ≡ x as average, the root-mean-square error the common part of commuters ( 50) which is 1 if all links are identical and 0 if none of them agrees.All the above measures are more sensitive to strong links, i.e. large import probabilities.However, if the focus is to get a fair comparison including all links, we are more interested in logarithmic versions of the above measures or rank correlations.Thus, we compare the logarithm of the import probabilities via correlation logcorr( , ) = corr(log( ), log( )) , root-mean-square error and use the Kendall rank correlation coefficient with  and  as the number of concordant and discordant pairs and   and   as ties only in  and , respectively.
To simplify and generalize the comparison we combine the six above defined measures by computing the mean rank of each model, i.e. the best correlating model has the highest (12) and the worst the lowest (0) rank and the mean rank of one model is the average of all six ranks.To quantify the mean difference between the models we define the relative performance of one model  as with  (   ) =  (   , ) as the specific comparison function and best(  (   ), ) and worst(  (   ), ) as the best and worst performing value of all models using this comparison function.Note, that best(. . . ) = max(. . . ) apart for the rmse-measures, where it is min(. . . ) (analog for worst(. . .)).

D R A F T
Grav.

e d
Grav.
e d (geo.path) Grav.The model that agrees best (worst) with the reference import risk according a specific measure has the highest (lowest) rank and a relative performance of one (zero).The relative performance is then a linear interpolation between the best and worst model.The comparison measures are the correlation (corr), correlation between the log-transformed import risks (log_corr), Root-mean-square error (rmse), Root-mean-square error of the log-transformed import risks (log_rmse), Kendall rank correlation (kendalltau) and the common part of commuters (cpc).As exponents of the gravity models the mean optimal parameter is used (horizontal lines in Fig. S3).
e  As exponents of the gravity models the mean optimal parameter is used (horizontal lines in Fig. S3).The colors depict the 4 different models.The solid circles are the models with corrected import probability by symmetrizing their OD-matrix, and the transparent squares are the non-corrected import probabilities of the respective model.  .B: the effective distance offset  0 that penalizes larger hop-distances in the effective distance  eff (| 0 ) =  0 − ln(   ) when creating the shortest path tree (default:  0 = 1).C: the descendant fraction introduced in the shortest path exit probability, where 0.5 is the default value and values larger than 0.5 mean that the exiting at the descendant (or offspring) nodes compared to the current node becomes more likely.D: different weight options introduced for the shortest path tree exit probability.Per default, the node populations are not weighted.The weight is the inverse of either the geodesic or the effective distance.E: manually set shortest path exit probability of leaf nodes (dead-end nodes).Per default, the exit probability is 1.A decrease to 0.9 or 0.8 does not visually change the median.4), however, here the uncorrected model predictions are used, i.e. without symmetrizing the OD-matrix.6 but here instead of the import risk model, the gravity model with power-law distance decaying function using the geodesic  geo (left) or effective  eff (right) distance is applied.Also for these models the logcorr between import probability estimates (| 0 ) and the reference data p(| 0 ) improves for countries with a larger outflow in the WAN.

Fig. 1 .
Fig. 1.The relation between WAN, OD-probability, SPT and effective distance.A: The world air transportation network (WAN) represents the direct flight connections and maximal seat capacities between airports in 2014, here shown for flights starting from five selected countries.It is based on flight-schedule-data.The lines are bundled and do not represent the specific flight route, but illustrate the links to airports abroad.B: The reference import probability from Canada to all countries, based on the OD matrix (Origin-Destination) of the Global Transnational Mobility Data set (35, 41) in 2014.It combines origin and final-destination trips between countries from the SABRE and the World Tourism Organization (UNWTO).The lines illustrate the connection to the common source country.C: Based on the effective distance  eff =  0 − ln( ) a shortest path tree (SPT) is constructed with the largest Canadian airport as source (YYZ: Toronto Pearson International Airport).The link color and thickness shows the hop distance, i.e. number of connecting flights.D: exponential decay of the reference import probability (as in B but for all countries as source) with the effective distance  eff (derived from the SPT (C) of the WAN (A)).Each dot represents a country-country link, the lines are medians including either all source countries or only from a specific continent.

Fig. 3 .
Fig. 3. Estimates of import probability by the gravity model with exponentially (1st column) and power law (2nd column) decaying distance function, the radiation model (3rd.column) and by the import risk model (4th column).The first three models (1st-3rd column) use as distance the geodesic (1st row), geodesic path (2nd row) and the effective (3rd row) distance.The import risk model is computed from the WAN with the geodesic distance (D) or the effective distance (L) as a weight for the exit probabilities or without weighting (H), i.e. in the last two cases (H, D) only WAN information is used.The orange line depicts the median and the gray line is  =  and illustrates perfect mapping.

Fig. 4 .
Fig. 4. Rank and relative performance of import risk estimation models.The different import probability models are compared via their rank (A) and relative performance (B), with the highest values representing the best approach.The rank and relative performance are shown for each (black dots) of the six comparison measures (corr, logcorr, RMSE, logRMSE, cpc,  Kendall ) the box illustrates the interquartile range, the horizontal line the median and the red triangle the mean.The colors of the boxes illustrate the different distance measures in use.The outlier measure of the import risk models (I.R.) is the logRMSE, where the gravity models with effective distance are performing best.See Material and Methods for definitions of comparison measures and Figs.S5, S6 for absolute and detailed relative performance.

Fig. 5 .
Fig. 5. Classification of the 10 countries with the highest import probability by the gravity model with exponentially (1st column) and power law decaying (2nd column) distance function, the radiation model (3rd.column) and by the import risk model (4th column).A true or false positive (T.Pos. or F. Pos.) means that the country is or is not among the 10 countries with the highest reference import probability p.A false negative (F.Neg.) means that it belongs to the reference set but was not detected by the respective model.The pie chart illustrates the sensitivity of the models.

Fig. 6 .
Fig. 6.Source countries' prediction quality and WAN outflow.The correlation between the logarithm of the import risk and the reference import probability logcorr = corr(log(  ∞ ), log( p)) improves with the outflow of the respective source country (top).Examples of source countries with particularly low (ER, Eritrea) and high (GB, Great Britain) outflow and log _  are shown with their import risk and reference import risk to target countries (middle row).Countries with exceptionally low log _  measures compared to source countries with a comparable outflow are either historically linked to specific regions as Australia (AU) and Israel (IL) to European countries (lower right panel) or politically as Macao (MO) as a special administrative region of China.

5 CFig. S1 .Fig. S2 .
Fig. S1.Import probability dependence on the geographic distance (A), the effective distance (B) and the geographic path distance (C).The orange line represents the median and  (, ) is the correlation between the two measures either log-transformed or not.The geographic distance between countries is averaged over all airport pairs.The geographic path distance is the geographic distance along the shortest path derived from the WAN using  eff , i.e. it is a combination of geographic and network information.The axis scale corresponds to the one with the highest correlation, i.e. log-log for distance and path distance (A, C) and y-log for the effective distance (B).

Fig. S3 .
Fig. S3.Mean optimal parameters for gravity models.For each gravity model with exponentially and power law decaying distance function and with one of the three different distance measures (geodesic distance, geodesic path distance and effective distance), the exponent  or  that results in the best fit to the reference import risk is shown.The comparison is quantified via the correlation (corr), correlation between the log-transformed import risks (log_corr), root mean square error (rmse), root mean square error of the log-transformed import risks (log_rmse), Kendall rank correlation (kendalltau) and the common part of commuters (cpc).The mean optimal parameter for each model is marked by a horizontal line and their values are  = [5.68,6.11] * 10 −4 for geographic and geo.path distance and  = 0.863 for the effective distance, and  = [1.84, 1, 87, 5.41] for geo., geo.path, and effective distance, respectively.

Fig. S4 .
Fig. S4.Symmetry check for OD-matrix.Each dot represents the number of passengers that travel between 2 countries and back.The OD-matrix is computed by the radiation model (1st.column), gravity model with exponentially (2nd column) and power law decaying (3rd column) distance function and by the import risk model (4th column).The OD-matrix of the models is computed by multiplying the import probability with the source-outflow.The reference trips and return trips have the highest symmetry (5th column, M).The orange line depicts the median and the gray line is  =  and illustrates perfect symmetry.The mean (AVG(  )) and median (MED(  )) asymmetry of the flows, computed according to Eq. S3, are shown in each panel.The reference trips (M) show the lowest asymmetry, especially for large passenger flows.

Fig. S6 .
Fig. S6.Absolute comparison measures for the import probability estimates.The comparison measures are the correlation (corr), correlation between the log-transformed import risks (log_corr), Root-mean-square error (rmse), Root-mean-square error of the log-transformed import risks (log_rmse), Kendall rank correlation (kendalltau) and the common part of commuters (cpc).As exponents of the gravity models the mean optimal parameter is used (horizontal lines in Fig.S3).The colors depict the 4 different models.The solid circles are the models with corrected import probability by symmetrizing their OD-matrix, and the transparent squares are the non-corrected import probabilities of the respective model.

Fig. S7 .
Fig. S7.Import risk comparison and its deviation from a linear relation.Scatter plot (left) and only median and IQR with an exponential fit (right).

Fig. S8 .
Fig. S8.Variations of the import risk model to investigate how additional parameters influence the relation between the import riskand the reference import risk.A: the flow scaling exponent  that estimates the travelling population  () of the airport  depending on its WAN outflow   via  () =    (default:  = 1).B: the effective distance offset  0 that penalizes larger hop-distances in the effective distance  eff (| 0 ) =  0 − ln(   ) when creating the shortest path tree (default:  0 = 1).C: the descendant fraction introduced in the shortest path exit probability, where 0.5 is the default value and values larger than 0.5 mean that the exiting at the descendant (or

10
Fig. S9.Country outflow reconstruction by import risk.The flow in the WAN leaving a country   is estimated by the import risk model by   = ∈ ∉  ∞ (|)  .Both measures are directly compared (A) and the relative error is computed depending on the number of airports in the respective country     (B).The import risk model does not include the concept of a country which partly explains the overestimation for larger airports.Another explanation is the overestimation of the respective airport population   =   by the WAN outflow for the import risk model (the true population is smaller because of the transit passengers that need to be excluded).Note that the WAN is used here, i.e. we check for self-consistency of the model and no reference data is included.
Fig.S10.Uncorrected models: rank and relative performance.Same analysis as in the main text in Fig.4), however, here the uncorrected model predictions are used, i.e. without symmetrizing the OD-matrix.

Fig. S11 .
Fig. S11.Import risk between world regions to a specific target region.In contrast to its derivation the import risk is displayed in a target-centric view, i.e. each panel displays the import probability to a single target region from all source regions.The distance between world regions is the mean distance between their airport locations.The grey line represents a power-law fit  ∞ =  •  − .The mean import risk is marked for each world region by a horizontal dashed line.The 22 target-world-regions are sorted according to their mean import risk.

Fig. S12 .
Fig. S12.Source countries prediction quality and WAN outflow for two gravity models.Same model-result representation as inFig.6buthere instead of the import risk model, the gravity model with power-law distance decaying function using the geodesic  geo (left) or effective  eff (right) distance is applied.Also for these models the logcorr between import probability estimates (| 0 ) and the reference data p(| 0 ) improves for countries with a larger outflow in the WAN.
Fig. S13.WAN flow out of countries vs. population and GDPThe WAN flow out of a country is best mapped by its gross domestic product (GDP, C) compared to its population (A) or per capita GPD (B).The linear double-logarithmic regression results are shown in the lower part of each panel (r-and p-value).The size of each country corresponds to its population (A) and the color codes its continent.GDP is taken from the World Bank Dataset for the year 2014(54).
• =  (• |★), or continues walking with prob.(1 −  • ).The import risk  ∞ (• |★) (right) is the probability of a walker to exit at node • given it started at node ★ under consideration of all possible paths.