Lagrangian particle tracking

Two-dimensional Lagrangian particle tracking was used to make our connectivity calculations15,21,33,37. We used velocity fields from ECCO2 (http://ecco2.org), a high-resolution (1/4) global ocean model that assimilates available satellite and the in situ data26, to advect particles in the surface ocean (Supplementary Fig. 1). ECCO2 is based on a global full-depth ocean and sea-ice configuration of the Massachusetts Institute of Technology general circulation model (MITgcm) and applies an ad-joint approach to generate the physically consistent data assimilations. ECCO2’s resolution is high enough to permit the formation of eddies and other narrow current systems within the ocean.

Particles were advected in the surface ocean using TRACMASS (http://tracmass.org), an off-line particle tracking code that calculates trajectories using Eulerian velocity fields. TRACMASS estimates the trajectory path through each grid cell of every Lagrangian particle, using an analytical solution to a differential equation that depends on the velocities on the grid-box walls. The scheme was originally developed for stationary velocity fields20,38, and thereafter extended for time-dependent fields by solving a linear interpolation of the velocity field both in time and in space over each grid box39. This differs from the Runge-Kutta method, where trajectories are iterated forward in time with short time steps.

Particle seeding and connectivity patches

We seeded six particles in the second depth layer of each ECCO2 grid cell (a total of 4 million particles at each seeding time, or 36 million particles in total over all seeding times). When calculating connectivity, we aggregated the model’s 1/4° × 1/4° grid cells to 11,116 discrete 2° × 2° patches. The size of these connectivity patches was selected as a balance of computational feasibility and biogeographic detail. Each connectivity patch is therefore seeded with 384 particles at each seeding event (9 in total). The second depth layer is between 5 and 20 m depth and was used to avoid potential numerical problems due to how ECCO2 implement a varying sea surface, precipitation, and evaporation. See Supplementary Fig. 2 for the spatial distribution of connectivity patches. Particles were seeded at 9 points in time: 1 January 2001, 1 February 2002, 1 March 2003, 1 April 2004, 1 May 2005, 1 June 2006; 1 July 2007, 1 August 2008; and 1 September 2009 in model years. As a consequence of the multiple seeding times, a total of 3,456 particles were used per patch to estimate connectivity. Particles were then advected using horizontal velocity fields from the second depth layer in ECCO2 so that they were locked in the surface ocean. We looped velocity fields for the years 2000–2010 continuously and advected the particles for 100 years in total. Particle positions were saved every 3 days and used to calculate minimum connection times. No extra diffusivity was added to the movement of the particles. Supplementary Figure 6 shows the relationship between advection time and number of other patches reached. It is clear from this figure that the number of connectivity patches reached saturates after about 12 years. In other words, like the number of particles released, there are diminishing returns to running simulations for longer integration times.

Estimating the timescales of connectivity

The resulting Lagrangian particle trajectories were used to estimate the shortest time taken for water to travel from one patch in the surface ocean to another. This minimum connection time is a variant on the standard measure of ocean distance, which is the expected transit time for water to travel from one patch to another13,14. We use the minimum and not the expected connection time for two reasons. First, the minimum connection time is a more appropriate metric for phytoplankton and bacterial connectivity since asexually reproducing organisms have high reproductive output that attenuates low dispersal probabilities, and only a few individuals are required to ‘connect’ two places, especially in terms of population genetics22. It is therefore possible for such organisms to exploit dispersal routes where the probability to reach a given destination is very low. Second, expected transit times in the global ocean are not properly defined, as water can recirculate for an infinitely long time. There is no limit, therefore, to the distribution of connectivity times over which to calculate expected connection times. Thus, the minimum connection time is a preferable alternative measure of ocean distance for this global application.

Minimum connection times for the global surface ocean are called Min-T, and they are stored in the form of a matrix—the Min-T connectivity matrix, where each i, j element represents the shortest transit time between a given source patch i and destination patch j (Fig. 1a). The raw Min-T matrix, produced from the Lagrangian particle tracking, is highly sparse with most pairs of patches being unconnected.

Network analysis of shortest/quickest paths

Estimation of connection times between all pairs of patches globally using Lagrangian particle simulations alone would require a currently infeasible number of particles23 (see particle density sensitivity test described below). To circumvent this obstacle, a shortest-path algorithm was used to calculate missing values in the raw Min-T connectivity matrix. Here the network is the global ocean, with patches in the ocean as nodes, and minimum connection times as edges connecting the nodes. Applied to this network, shortest path algorithms identify the shortest path between every global ocean patch-pair, accounting for all possible multistep connections (see Supplementary Online Material for details). For example, if there is no direct connection between nodes A and D, then these algorithms identify the multistep connection from A→B→C→D. We use Dijkstra’s algorithm24, which is one of the most commonly used shortest path algorithms, and which fits our specific application. The end result is a modified Min-T connectivity matrix (Fig. 1b), where all possible minimum connection times between patches are calculated.

Each step along these minimum-time routes may be unlikely, and so the conditional multistep probability (of going from A to B to C to D...) can have a very low probability as well. However, we assume that the effect of these low probabilities is attenuated by the large reproductive output of microorganisms drifting with ocean currents. Over the timescales that we are considering, microorganisms moving with water masses can grow by the million25. Hence, there will still be planktonic organisms traveling along the potentially low probability paths identified here. Indeed, if one considers the dispersal of genetic material, then there need only be a small number of individuals traveling along these Min-T routes, to make them evolutionarily relevant17,22,25.

While nodes in a network are usually defined as singular nodes with well-defined distances between them, our ocean patches have relatively large areas and are continuously adjacent to one another. This difference creates a problem when using Dijkstra’s algorithm since a particle seeded next to the boundary of its initial patch can rapidly move to an adjacent patch. However, shortest path algorithms assumes that the travel time across each intermediate node is zero, or at least included in the edge distances. (This phenomenon is also a problem when analysing the speed of tracer transport in General Circulation Models40.) By removing all calculated connectivity times shorter than 1 year before applying the shortest-path algorithm, we limit the effect of not including within-patch crossing times. The 365-day cutoff is based on calculated typical residence times in the patches, which are on the order of weeks. All initial minimum connectivity times are based on travel distances at least an order of magnitude longer than typical patch crossing distances.

The removed connectivity times were added back to the final connectivity matrix, allowing for connection times shorter than 365 days, as shown in Figs 2 and 3. It should be noted that the absolute number of connectivity times shorter than 1 year in Fig. 3 are identical for the Raw and Dijkstra cases. The lack of discontinuities between sub-annual and longer connection times in Fig. 2 (and all other cases we have explored) give us confidence that the resulting connectivity matrix is reasonable and that our approach works.

After applying Dijkstra’s algorithm, we find that the resulting minimum connection time matrices are all connected. However, we do find some areas that are only connected in one direction (that is, there are connection time to, but not from, particular regions). These areas are mainly inland seas—the Baltic and Mediterranean, for example. However, they only account for a small fraction (2%) of the modified Min-T matrix and, consequently, do not impact the general result of the timescales of global surface ocean connectivity.

Particle seeding sensitivity test

Since the number of particles seeded per grid-cell and the seeding times are limited, we have not accounted for all possible Min-T pathways. As a result, our estimates of the timescales of global surface connectivity are conservative, since adding more particles and seeding dates could only lead to shorter Min-T pathways (that is, we look for the shortest connection times over all possibilities including seeding times). Thus, the few seeding dates—although arguably numerically incomplete—strengthen our conclusion that the global surface ocean is well connected over a few decades.

To examine the effect of particle seeding density, we performed a particle sensitivity test. Minimum connection times from a patch in the north pacific to all others were estimatedusing simulations with increasing numbers of seeded particles. Supplementary Figure 3 shows the results of these simulations. It is clear that a larger oceanic extent is reached as the number of particles released increases. However, when we examine only those patches that were reached in all seeding experiments (Supplementary Fig. 4), we can see that increasing the number of particles serves only to decrease the minimum connection times in these patches (Supplementary Fig. 4: with 84 particles some areas are reached after 100 years—the patches in gold, in contrast with 16,660 particles, these same patches are reached after around 20 years—patches now in light red).

Finally, we show the aggregated results of the sensitivity test in Supplementary Fig. 5, where we plot the fraction of patches reached (over the whole ocean) and the median minimum connection time from this study patch. The fraction of patches reached saturates at around 30%, which means that there are diminishing returns (in terms of estimating minimum connection times to new patches) to adding more particles. It also indicates that, to release enough particles to estimate minimum connection times to all patches globally, a currently impossible number of Lagrangian particles would be required. Similarly, the median minimum connection time from this patch saturates at around 8,000 particles released. In our simulations we use 3,456 particles per patch as this achieved a balance of connectivity sampling power and computational efficiency.