Ecology

1.
Garibaldi, L. The FAO global capture production database: A six-decade effort to catch the trend. Mar. Pol. 36, 760–768 (2012).
Article  Google Scholar 
2.
FAO (Food and Agriculture Organization). The State of World Fisheries and Aquaculture 2018. in Meeting The Sustainable Development Goals. (FAO, Rome, 2018).

3.
Grant, W. S., Jasper, J., Bekkevold, D. & Adkison, M. Responsible genetic approach to stock restoration, sea ranching and stock enhancement of marine fishes and invertebrates. Rev. Fish Biol. Fish. 27, 615–649 (2017).
Article  Google Scholar 

4.
Christie, M. R., Marine, M. L., French, R. A., Waples, R. S. & Blouin, M. S. Effective size of a wild salmonid population is greatly reduced by hatchery supplementation. Heredity 109, 254–260 (2012).
CAS  PubMed  PubMed Central  Article  Google Scholar 

5.
Ryman, N. & Laikre, L. Effects of supportive breeding on the genetically effective population size. Conserv. Biol. 5, 325–329 (1991).
Article  Google Scholar 

6.
Waples, R. S., Hindar, K., Karlsson, S. & Hard, J. J. Evaluating the Ryman–Laikre effect for marine stock enhancement and aquaculture. Curr. Zool. 62, 617–627 (2016).
PubMed  PubMed Central  Article  Google Scholar 

7.
Sun, X. & Hedgecock, D. Temporal genetic change in North American Pacific oyster populations suggests caution in seascape genetics analyses of high gene-flow species. Mar. Ecol. Prog. Ser. 565, 79–93 (2017).
ADS  Article  Google Scholar 

8.
Bacheler, N. M., Wong, R. A. & Buckel, J. A. Movements and mortality rates of striped mullet in North Carolina. N. Am. J. Fish. Manage. 25, 361–373 (2005).
Article  Google Scholar 

9.
Whitfield, A. K., Panfili, J. & Durand, J.-D. A global review of the cosmopolitan flathead mullet Mugil cephalus Linnaeus 1758 (Teleostei: Mugilidae), with emphasis on the biology, genetics, ecology and fisheries aspects of this apparent species complex. Rev. Fish Biol. Fish. 22, 641–681 (2012).
Article  Google Scholar 

10.
Hsu, C.-C., Chang, C.-W., Iizuka, Y. & Tzeng, W.-N. A growth check deposited at estuarine arrival in otoliths of juvenile flathead mullet (Mugil cephalus L.). Zool. Stud. 48(3), 315–323 (2009).
Google Scholar 

11.
Antuofermo, E. et al. First evidence of intersex condition in extensively reared mullets from Sardinian lagoons (central-western Mediterranean, Italy). Ital. J. Anim. Sci. 16, 283–291 (2017).
Article  Google Scholar 

12.
Heras, S., Roldán, M. I. & Castro, M. G. Molecular phylogeny of Mugilidae fishes revised. Rev. Fish Biol. Fish. 19, 217–231 (2009).
Article  Google Scholar 

13.
Heras, S., Maltagliati, F., Fernández, M. V. & Roldán, M. I. Shaken not stirred: A molecular contribution to the systematics of genus Mugil (Teleostei, Mugilidae). Integr. Zool. 11, 263–281 (2016).
PubMed  Article  Google Scholar 

14.
Shen, K.-N., Jamandre, B. W., Hsu, C.-C., Tzeng, W.-N. & Durand, J.-D. Plio-Pleistocene sea level and temperature fluctuations in the northwestern Pacific promoted speciation in the globally-distributed flathead mullet Mugil cephalus. BMC Evol. Biol. 11, 83. https://doi.org/10.1186/1471-2148-11-83 (2011).
CAS  Article  PubMed  PubMed Central  Google Scholar 

15.
Durand, J.-D. et al. Systematics of the grey mullets (Teleostei: Mugiliformes: Mugilidae): Molecular phylogenetic evidence challenges two centuries of morphology-based taxonomy. Mol. Phylogenet. Evol. 64, 73–92 (2012).
PubMed  Article  Google Scholar 

16.
Rossi, A. R., Capula, M., Crosetti, D., Campton, D. E. & Sola, L. Genetic divergence and phylogenetic inferences in five species of Mugilidae (Pisces: Perciformes). Mar. Biol. 131, 213–218 (1998).
CAS  Article  Google Scholar 

17.
Blel, H. et al. Selection footprint at the first intron of the Prl gene in natural populations of the flathead mullet (Mugil cephalus, L. 1758). J. Exp. Mar. Biol. Ecol. 387, 60–67 (2010).
CAS  Article  Google Scholar 

18.
Durand, J., Blel, H., Shen, K., Koutrakis, E. & Guinand, B. Population genetic structure of Mugil cephalus in the Mediterranean and Black Seas: A single mitochondrial clade and many nuclear barriers. Mar. Ecol.-Prog. Ser. 474, 243–261 (2013).
ADS  Article  Google Scholar 

19.
Šegvić-Bubić, T. et al. Range expansion of the non-native oyster Crassostrea gigas in the Adriatic Sea. Acta Adriat. 57(2), 321–330 (2016).
Google Scholar 

20.
Piras, P. et al. A case study on the labeling of bottarga produced in Sardinia from ovaries of grey mullets (Mugil cephalus and Mugil capurrii) caught in Eastern Central Atlantic coasts. Ital. J. Food. Saf. 7(1), 6893. https://doi.org/10.4081/ijfs.2018.6893 (2018).
Article  PubMed  PubMed Central  Google Scholar 

21.
Miggiano, E. et al. Isolation and characterization of microsatellite loci in the striped mullet, Mugil cephalus. Mol. Ecol. Notes 5, 323–326 (2005).
CAS  Article  Google Scholar 

22.
Jinliang, W. A parsimony estimator of the number of populations from a STRUCTURE‐like analysis. Mol. Ecol. Resour. 19, 970–981 (2019).
Article  CAS  Google Scholar 

23.
FAO (Food and Agriculture Organization). Code of Conduct for Responsible Fisheries (FAO, Rome, 1995).
Google Scholar 

24.
Mai, A. C. G. et al. Microsatellite variation and genetic structuring in Mugil liza (Teleostei: Mugilidae) populations from Argentina and Brazil. Estuar. Coast. Shelf Sci. 149, 80–86 (2014).
ADS  Article  Google Scholar 

25.
Pacheco-Almanzar, E., Simons, J., Espinosa-Perez, H., Chiappa-Carrara, X. & Ibanez, A. L. Can the name Mugil cephalus (Pisces: Mugilidae) be used for the species occurring in the north western Atlantic?. Zootaxa 4109, 381–390 (2016).
PubMed  Article  Google Scholar 

26.
Waples, R. S. & Gaggiotti, O. What is a population? An empirical evaluation of some genetic methods for identifying the number of gene pools and their degree of connectivity. Mol. Ecol. 15, 1419–1439 (2006).
CAS  PubMed  Article  Google Scholar 

27.
Hauser, L. & Carvalho, G. R. Paradigm shifts in marine fisheries genetics: Ugly hypotheses slain by beautiful facts. Fish Fish. 9, 333–362 (2008).
Article  Google Scholar 

28.
Waples, R. S. & England, P. R. Estimating contemporary effective population size on the basis of linkage disequilibrium in the face of migration. Genetics 189, 633–644 (2011).
PubMed  PubMed Central  Article  Google Scholar 

29.
Murenu, M., Olita, A., Sabatini, A., Follesa, M. C. & Cau, A. Dystrophy effects on the Liza ramada (Risso, 1826) (Pisces, Mugilidae) population in the Cabras lagoon (Central-Western Sardinia). Chem. Ecol. 20, 425–433 (2004).
Article  Google Scholar 

30.
Ryman, N., Laikre, L. & Hössjer, O. Do estimates of contemporary effective population size tell us what we want to know?. Mol. Ecol. 28, 1904–1918 (2019).
PubMed  PubMed Central  Article  Google Scholar 

31.
Guinand, B. et al. Candidate gene variation in gilthead sea bream reveals complex spatiotemporal selection patterns between marine and lagoon habitats. Mar. Ecol. Prog. Ser. 558, 115–127 (2016).
ADS  CAS  Article  Google Scholar 

32.
Chaoui, L. et al. Microsatellite length variation in candidate genes correlates with habitat in the gilthead sea bream Sparus aurata. Mol. Ecol. 21, 5497–5511 (2012).
CAS  PubMed  Article  Google Scholar 

33.
González-Wangüemert, M. & Pérez-Ruzafa, Á. In two waters: contemporary evolution of lagoonal and marine white seabream (Diplodus sargus) populations. Mar. Ecol. 33, 337–349 (2012).
ADS  Article  Google Scholar 

34.
Cardona, L. Effects of salinity on the habitat selection and growth performance of Mediterranean flathead grey mullet Mugil cephalus (Osteichthyes, Mugilidae). Estuar. Coast. Shelf Sci. 50, 727–737 (2000).
ADS  Article  Google Scholar 

35.
Fortunato, R. C., Galán, A. R., Alonso, I. G., Volpedo, A. & Durà, V. B. Environmental migratory patterns and stock identification of Mugil cephalus in the Spanish Mediterranean Sea, by means of otolith microchemistry. Estuar. Coast. Shelf. Sci. 188, 174–180 (2017).
ADS  Article  CAS  Google Scholar 

36.
Jones, A. G., Small, C. M., Paczolt, K. A. & Ratterman, N. L. A practical guide to methods of parentage analysis. Mol. Ecol. Resour. 10, 6–30 (2010).
PubMed  Article  Google Scholar 

37.
Taylor, H. R. The use and abuse of genetic marker-based estimates of relatedness and inbreeding. Ecol. Evol. 5, 3140–3150 (2015).
PubMed  PubMed Central  Article  Google Scholar 

38.
Coppinger, C. R. et al. Assessing the genetic diversity of catface grouper Epinephelus andersoni in the subtropical Western Indian Ocean. Fish. Res. 218, 186–197 (2019).
Article  Google Scholar 

39.
Cushman, E. L. et al. Development of a standardized molecular tool and estimation of genetic measures for responsible aquaculture-based fisheries enhancement of American Shad in North and South Carolina. Trans. Am. Fish. Soc. 148, 148–162 (2019).
Article  Google Scholar 

40.
Waples, R. S., Punt, A. E. & Cope, J. M. Integrating genetic data into management of marine resources: How can we do it better?. Fish. Fish. 9, 423–449 (2008).
Article  Google Scholar 

41.
Iacchei, M. et al. Combined analyses of kinship and FST suggest potential drivers of chaotic genetic patchiness in high gene-flow populations. Mol. Ecol. 22, 3476–3494 (2013).
PubMed  PubMed Central  Article  Google Scholar 

42.
Bernardi, G., Beldade, R., Holbrook, S. J. & Schmitt, R. J. Full-sibs in cohorts of newly settled coral reef fishes. PLoS ONE 7(e44953), 2012. https://doi.org/10.1371/journal.pone.0044953 (2012).
CAS  Article  Google Scholar 

43.
Como, S., van der Velde, G. & Magni, P. Temporal variation in the trophic levels of secondary consumers in a Mediterranean coastal lagoon (Cabras lagoon, Italy). Estuaries Coasts 41, 218–232 (2018).
CAS  Article  Google Scholar 

44.
Floris, R., Manca, S. & Fois, N. Microbial ecology of intestinal tract of gilthead sea bream (Sparus aurata Linnaeus, 1758) from two coastal lagoons of Sardinia (Italy). Transit. Waters Bullet. 7(2), 4–12. https://doi.org/10.1285/i1825229Xv7n2p4 (2013).
Article  Google Scholar 

45.
Merella, P. & Garippa, G. Metazoan parasites of grey mullets (Teleostea: Mugilidae) from the Mistras lagoon (Sardinia-Western Mediterranean). Sci. Mar. 65, 201–206 (2001).
Article  Google Scholar 

46.
Pitacco, V. et al. Spatial patterns of macrobenthic alpha and beta diversity at different scales in Italian transitional waters (Central Mediterranean). Estuar. Coast. Shelf Sci. 222, 126–138 (2019).
ADS  Article  Google Scholar 

47.
Cioffi, F. & Gallerano, F. From rooted to floating vegetal species in lagoons as a consequence of the increases of external nutrient load: An analysis by model of the species selection mechanism. Appl. Math. Model. 30, 10–37 (2006).
MATH  Article  Google Scholar 

48.
Wasko, A. P., Martins, C., Oliveira, C. & Foresti, F. Non-destructive genetic sampling in fish. An improved method for DNA extraction from fish fins and scales. Hereditas 138, 161–165 (2003).
PubMed  Article  Google Scholar 

49.
Waples, R. S. Testing for Hardy–Weinberg proportions: Have we lost the plot?. J. Hered. 106, 1–19 (2015).
PubMed  Article  Google Scholar 

50.
Rousset, F. genepop’007: A complete re-implementation of the genepop software for Windows and Linux. Mol. Ecol. Resour. 8, 103–106 (2008).
PubMed  Article  Google Scholar 

51.
Kinnison, M. T., Bentzen, P., Unwin, M. J. & Quinn, T. P. Reconstructing recent divergence: Evaluating nonequilibrium population structure in New Zealand chinook salmon. Mol. Ecol. 11, 739–754 (2002).
CAS  PubMed  Article  Google Scholar 

52.
Benjamini, Y. & Yekutieli, D. The control of the false discovery rate in multiple testing under dependency. Ann. Stat. 29(4), 1165–1188 (2001).
MathSciNet  MATH  Article  Google Scholar 

53.
R Core Team. R: A Language and Environment for Statistical Computing (R Foundation for Statistical Computing, Vienna, 2018).
Google Scholar 

54.
Cossu, P. et al. Influence of genetic drift on patterns of genetic variation: The footprint of aquaculture practices in Sparus aurata (Teleostei: Sparidae). Mol. Ecol. 28, 3012–3024 (2019).
PubMed  Article  Google Scholar 

55.
Van Oosterhout, C., Hutchinson, W. F., Wills, D. P. M. & Shipley, P. Micro-checker: software for identifying and correcting genotyping errors in microsatellite data. Mol. Ecol. Notes 4, 535–538 (2004).
Article  CAS  Google Scholar 

56.
Chapuis, M.-P. & Estoup, A. Microsatellite null alleles and estimation of population differentiation. Mol. Biol. Evol. 24, 621–631 (2007).
CAS  PubMed  Article  Google Scholar 

57.
Dąbrowski, M. J. et al. Reliability assessment of null allele detection: Inconsistencies between and within different methods. Mol. Ecol. Resour. 14, 361–373 (2014).
PubMed  Article  CAS  Google Scholar 

58.
Foll, M. & Gaggiotti, O. A genome-scan method to identify selected loci appropriate for both dominant and codominant Markers: A Bayesian perspective. Genetics 180, 977–993 (2008).
PubMed  PubMed Central  Article  Google Scholar 

59.
Kauer, M. O., Dieringer, D. & Schlötterer, C. A microsatellite variability screen for positive selection associated with the “Out of Africa” habitat expansion of drosophila melanogaster. Genetics 165, 1137–1148 (2003).
CAS  PubMed  PubMed Central  Google Scholar 

60.
Excoffier, L. & Lischer, H. E. L. Arlequin suite ver 3.5: A new series of programs to perform population genetics analyses under Linux and Windows. Mol. Ecol. Resour. 10, 564–567 (2010).
PubMed  Article  Google Scholar 

61.
Paris, M. et al. Genome scan in the mosquito Aedes rusticus: population structure and detection of positive selection after insecticide treatment. Mol. Ecol. 19, 325–337 (2010).
PubMed  Article  Google Scholar 

62.
Keenan, K., McGinnity, P., Cross, T. F., Crozier, W. W. & Prodöhl, P. A. diveRsity: An R package for the estimation and exploration of population genetics parameters and their associated errors. Methods Ecol. Evol. 4, 782–788 (2013).
Article  Google Scholar 

63.
Piry, S., Luikart, G. & Cornuet, J.-M. Computer note. BOTTLENECK: A computer program for detecting recent reductions in the effective size using allele frequency data. J. Hered. 90, 502–503 (1999).
Article  Google Scholar 

64.
Do, C. et al. NeEstimator v2: Re-implementation of software for the estimation of contemporary effective population size (Ne) from genetic data. Mol. Ecol. Resour. 14, 209–214 (2014).
CAS  PubMed  Article  Google Scholar 

65.
Waples, R. S. & Do, C. Linkage disequilibrium estimates of contemporary Ne using highly variable genetic markers: A largely untapped resource for applied conservation and evolution. Evol. Appl. 3, 244–262 (2010).
PubMed  Article  Google Scholar 

66.
Weir, B. S. & Cockerham, C. C. Estimating F-statistics for the analysis of population structure. Evolution 38, 1358–1370 (1984).
CAS  PubMed  Google Scholar 

67.
Jost, L. GST and its relatives do not measure differentiation. Mol. Ecol. 17, 4015–4026 (2008).
PubMed  Article  Google Scholar 

68.
Ryman, N. & Palm, S. POWSIM: A computer program for assessing statistical power when testing for genetic differentiation. Mol. Ecol. Notes 6, 600–602 (2006).
Article  Google Scholar 

69.
Blouin, M. S., Parsons, M., Lacaille, V. & Lotz, S. Use of microsatellite loci to classify individuals by relatedness. Mol. Ecol. 5, 393–401 (1996).
CAS  PubMed  Article  Google Scholar 

70.
Kraemer, P. & Gerlach, G. Demerelate: Calculating interindividual relatedness for kinship analysis based on codominant diploid genetic markers using R. Mol. Ecol. Resour. 17, 1371–1377 (2017).
CAS  PubMed  Article  Google Scholar 

71.
Kalinowski, S. T., Wagner, A. P. & Taper, M. L. ml-relate: A computer program for maximum likelihood estimation of relatedness and relationship. Mol. Ecol. Notes 6, 576–579 (2006).
CAS  Article  Google Scholar 

72.
Pritchard, J. K., Stephens, M. & Donnelly, P. Inference of population structure using multilocus genotype data. Genetics 155, 945–959 (2000).
CAS  PubMed  PubMed Central  Google Scholar 

73.
Falush, D., Stephens, M. & Pritchard, J. K. Inference of population structure using multilocus genotype data: Linked loci and correlated allele frequencies. Genetics 164, 1567–1587 (2003).
CAS  PubMed  PubMed Central  Google Scholar 

74.
Francis, R. M. pophelper: An R package and web app to analyse and visualize population structure. Mol. Ecol. Resour. 17, 27–32 (2017).
CAS  PubMed  Article  Google Scholar 

75.
Jombart, T., Devillard, S. & Balloux, F. Discriminant analysis of principal components: A new method for the analysis of genetically structured populations. BMC Genet. 11, 94. https://doi.org/10.1186/1471-2156-11-94 (2010).
Article  PubMed  PubMed Central  Google Scholar 

76.
Jombart, T. adegenet: A R package for the multivariate analysis of genetic markers. Bioinformatics 24, 1403–1405 (2008).
CAS  PubMed  Article  Google Scholar  More

All experimental procedures were approved by the University of Toronto animal care committee (Animal Use Protocol number 20012112) and conformed to guidelines prescribed by the Canadian Council on Animal Care.
Animal capture and husbandry
Wild male ruby-throated hummingbirds (Archilochus colubris; (n = 23); mass range during experimental period: 2.59 g to 4.52 g), were caught on the University of Toronto Scarborough Campus (({43.7838}^{circ }hbox {N}), ({79.1875}^{circ }hbox {W})) or the University of Western Ontario campus (({43.0096}^{circ }hbox {N}), ({81.2737}^{circ }hbox {W})) using box traps modified with hook-and-loop fastener tape on a drop door containing hummingbird feeders. Birds were trapped between 06:00 h and 12:00 h during the months of May through September of 2017, 2018, and 2019. Pilot trials were conducted in 03/2018. Birds in the pilot study were on wintering/migratory seasonality with a 12 h daylight schedule. Subsequent trials were conducted between 04/2019 and 01/2020. In 04/2019, birds were under breeding seasonality (14 h daylight) and in 01/2020, birds were under wintering/migratory seasonality (12 h daylight) during experimental trials. The daylight schedule approximated the photoperiod encountered as part of annual migrations to Central America and back. Upon capture, hummingbirds were quickly transported to metal EuroCages (({50.8 times 91.5 times 53.7},hbox {cm}) ((hbox {L}times hbox {W}times hbox {H}))) at the animal care facility where they were housed individually and acclimated to feed from syringe feeders. Birds were provided an 18 % (w/v) Nektar Plus (Guenter Enderle, Tarpon Springs, FL, USA) solution (henceforth referred to as maintenance diet), which was consumed ad libitum, and syringes were replaced daily (range of average consumption of daily maintenance diet during study period was 5.4 mL13.2 mL).
Experimental design
Birds drank solutions of imidacloprid (IMI; Sigma-Aldrich Cat. No. 37894) dissolved in a 20% w/v sucrose solution and were randomly assigned to either control (({0.0},upmu hbox {g g}^{-1}cdot)BW), low (({1.0}upmu hbox {g g}^{-1})), middle (({2.0}upmu hbox {g g}^{-1})), or high dose (({2.5}upmu hbox {g g}^{-1})) groups (n = 7, 4, 8, 4, respectively). Stock solution concentrations were analytically confirmed (low: ({0.32},hbox {gL}^{-1}), middle: ({0.59},hbox {gL}^{-1}), high: ({0.78},hbox {gL}^{-1})) such that a 3 g bird dosed with ({10},upmu hbox {L}) of solution would receive the dosage rate corresponding to either the low, middle or high dose. The volume of imidacloprid stock solution used for dosage was adjusted on a body weight (BW) basis, pipetting from the stock solution into a new nectar syringe and drawing up to a final volume of ({50},upmu hbox {L}) with 20% w/v sucrose solution, ensuring that birds received the same dosage rate throughout the trial. Birds were deprived of their regular nectar solution for 10 min to 15 min consumed the entire small-volume dosing solution within 10 min of being offered the solution. The dose was considered to be delivered when there was no visible solution remaining in the transparent syringe.
Doses were established within a range spanning expected exposure in a bird drinking ({10},hbox {mLd}^{-1}) from contaminated flowers10 up to 10 % of the LD50 in canaries61 (Serinus canaria, LD50: ({25},upmu hbox {g g}^{-1}) to 50 (upmu hbox {g g}^{-1})), similarly small birds with fast metabolic rates to target a sub-lethal concentration expected to produce toxic effects32. When energy demands are high, hummingbirds may consume over three times their body weight in nectar63, therefore ({10},hbox {mLd}^{-1}) is a probable figure for contaminated nectar consumption. Pooled blueberry flower samples collected about 1 year after treatment with imidacloprid contained the neonicotinoids at a concentration of ({5.16},hbox {ng g}^{-1})10. We extrapolated our very low and low dose concentrations based on these data. We stipulate that given the flower sample is a pooled sample, it was collected from flowers long after treatment, and there are different regulations on pesticide use within the ruby-throated hummingbird’s range, these doses were environmentally relevant.
We tested multiple intermediate doses which allowed us to explore dose-response relationships in observed effects64. Pilot experiment data with control, very low (({0.2},upmu hbox {g g}^{-1})), or high dose (({2.5},upmu hbox {g g}^{-1})) (n = 3 per group) are included for cholinesterase activity and toxicokinetic elimination analyses. Other metrics including behaviour and energy expenditure were not included from the pilot study due to differences in data collection protocols and are not strictly comparable. Behavioural data collection, cloacal fluid (CF) collection, and respirometry occurred over 6 days, where pre-dose data were collected for each animal on days 1 through 3, and dosing occurred once per day at 11:00 on days 4 through 6. Body weight measurements were taken daily at 10:00. The body weights of birds on the first day of experimentation ranged from 2.70 g to 4.52 g. For simplicity, 11:00 on days 1-6 is referred to as Dose Time (DT). Terminal sampling and tissue collection occurred 24 h after the third dose was administered. Birds were sacrificed by decapitation following isoflurane overdose, and whole blood, flight muscle, liver, brain, and heart tissues were rapidly excised, flash frozen in liquid nitrogen, and stored at ({-80},^{circ }hbox {C}) until downstream analysis, except in the case of blood which was immediately used for blood smear preparation.
Figure 3

Daily experimental timeline for days 1 through 6 of trials where on days 1 through 3, a control solution (20% w/v sucrose solution) is given in all groups and on days 4 through 6, dosing solutions were administered. Times of data collection are shown relative to Dose Time (DT). Terminal tissue sampling occurred on day 7 at DT, 24 h after the final dose was administered.

Full size image

Respirometry
Oxygen consumption and carbon dioxide production rates were measured using open-flow chamber respirometry65. Airflow through three metabolic chambers and one empty reference chamber was maintained at a rate of ({300},hbox {mL min}^{-1}). Excurrent air from the chambers was sub-sampled at ({100},hbox {mL min}^{-1}) sequentially starting with the reference chamber at using a Turbofox-5 (Sable Systems International Las Vegas, NV, USA). Sub-sampled air was passed through a water vapour pressure analyzer, a drying column (Indicating Drierite, W.A. Hammond Drierite, Xenia, Ohio, USA), carbon dioxide meter, and finally an oxygen analyzer (Turbofox-5, Sable Systems International). The oxygen and carbon dioxide analyzers were calibrated according to manufacturer instructions using well-mixed ambient air for the oxygen analyzer, and zero and (0.25 ,%hbox {CO}_{2}) reference gases for the (hbox {CO}_2) analyzer. Respirometry data were recorded at a frequency of 1 Hz using Expedata software (v. 1.84, Sable Systems) for 5 min while sampling from the empty reference chamber, followed by three 7 min recording periods from each of the chambers holding a bird. After this 26 min period, sub-sampling was resumed from the reference chamber for another 5 min followed by another 7 min sub-sampling period from experimental chambers. A final 5 min sampling of the reference chamber concluded the respirometry data collection, and birds were returned to the cloacal fluid collection chambers approximately 60 min after initially being placed in respirometry chambers.
Behavioural data collection and processing
Video recordings of birds were collected for 2 h, starting 4 h after DT (15:00–17:00). At the start of the recording period, birds were returned to their home cages where they could feed ad libitum by hovering and tracking a syringe on a 10 cm arm oscillating through a ({90}^{circ }) range along a lateral arc at a speed of 15 RPM. Video recordings were analyzed for time spent in flight, subdivided into foraging and non-foraging flights. Foraging flights were defined as flights where the bird contacted the hover feeder with their bill. Total consumption of the maintenance diet over this 2 h period was recorded.
Heterophil/lymphocyte ratios
Approximately ({2},upmu hbox {L}) of blood was collected for blood smear preparation immediately following sacrifice. After smearing, slides were left to air dry for a minimum of 3 h before fixing with 100 % methanol and staining with Giemsa–Wright solution (Fisher Scientific Cat. No. 123869). Slides were stained by immersion in eosinophilic dye (5times 1,hbox {s}) followed by (5times {1},hbox {s}) in basophilic dye.
Cholinesterase activity assay
Brain and muscle tissues were homogenized using a sonic dismembrator ((hbox {Fisherbrand}^{mathrm{TM}}) Model 120 Sonic Dismembrator) 1:10 w:v with ice-cold 0.1 M potassium phosphate buffer (pH 7.2). Samples were centrifuged at 10,000 RPM in a Beckman Coulter microfuge 22R centrifuge held at ({4},^{circ }hbox {C}) for 5 min. Total protein concentrations in tissue homogenates were determined by the Bradford assay (Sigma-Aldrich Cat. No. B6916). Cholinesterase activity was measured by the Ellman method adapted for a microplate reader (BioTek Synergy HT)66. Optimal assay conditions were 0.1 M potassium phosphate buffer (pH 7.2), 0.48 mM acetylcholine, 0.64 mM DTNB (Sigma-Aldrich Cat. No. D8130), 1.1 mM sodium bicarbonate. Assays were initiated through the addition of acetylcholine (Sigma-Aldrich, Cat. No. 01480) in a total volume of ({300},upmu hbox {L}). Absorbance was read at 412 nm every 2.5 min for 10 min.
Cloacal fluid
Collection
Cloacal fluid was collected for 1 h at 3 time points each day according to one of two schedules: starting (1) 1 h, 6 h, and 23 h, or (2) 2.5 h, 6.5 h, and 23 h after DT. Cloacal fluid was collected according to schedule (1) in pilot experiments, and (2) in the subsequent trials. A watch glass was placed beneath birds perching in 10 cm W (times) 12 cm H glass cylinder enclosures stopped with 19-gauge galvanized 1 cm hardware mesh openings in order to obtain cloacal fluid. To encourage greater cloacal fluid production, and to simulate the regular feeding behaviour of wild birds, individuals fed ad libitum from a syringe containing a 20% (w/v) sucrose solutions every 5 min to 10 min for the duration of cloacal fluid collection, which took place over 1 h as described under Sect. 4.2. After the collection period, cloacal fluid samples were stored at ({-20},^{circ }hbox {C}) until pooling and refreezing prior to chemical analysis.
Chemical analyses
Cloacal fluid samples and dosing solutions were analyzed for IMI by HPLC-ESI-MS/MS by Laboratory Services, NWRC (National Wildlife Research Centre, Ottawa, ON, Canada). Cloacal fluid samples were pooled by time point across dosing days by individual to reach the necessary minimum volume of ({100},upmu hbox {L}).
Cloacal fluid sample pools from 2018 trials were thawed at room temperature. Each pool was diluted 4 (times) with DI water (({25},upmu hbox {L}) cloacal fluid + ({75},upmu hbox {L}) DI water). The resulting ({100},upmu hbox {L}) diluted samples were then spiked with ({100},upmu hbox {L}) of internal standard (IS) solution. Spiked samples were filtered directly into ({300},upmu hbox {L}) glass inserts using 4 mm PVDF ({0.45},upmu text {m}) Millex filters and ({50},upmu hbox {L}) aliquots were injected. For the 2018 analyses, the minimum detection limit (MDL) and minimum reporting limit (MRL) were ({0.204},hbox {ng}hbox { mL}^{-1}) and ({0.616},hbox {ng}hbox { mL}^{-1}) respectively.
Cloacal fluid sample pools from 2019 trials were thawed at room temperature and ({50},upmu hbox {L}) of IS solution was added to (200,upmu text {L}) of pooled cloacal fluid. In cases where the sample volume was too small, volumes were adjusted: ({100},upmu hbox {L}) cloacal fluid + ({25},upmu hbox {L}) IS or ({80},upmu hbox {L}) cloacal fluid + ({20},upmu hbox {L}) IS as required. In these cases, duplicate injections of ({50},upmu hbox {L}) were not possible. All samples were filtered with 4 mm PVDF ({0.45},upmu text {m}) Millex filters prior to injection. For the 2019 analysis, the MDL and MRL were ({0.051},hbox {ng}hbox { mL}^{-1}) and ({0.154},hbox {ng}hbox { mL}^{-1}) respectively.
Cloacal fluid sample pools and dosing solutions were analyzed according to modifications to the methods of Main et al.9. Briefly, IMI in a cloacal fluid or DI water matrix was quantified by the internal standard method using the API5000 Triple Quadropole Mass Spectrometer (AB Sciex) and the TurboSpray ion source in positive polarity. The calibration curve was constructed from 8 concentrations ranging from ({0.1},hbox {ng}hbox { mL}^{-1}) to ({20},hbox {ng}hbox { mL}^{-1}) yielding an R greater than 0.995 (linear regression, no weighting). Injection cross-contamination was monitored by injecting solvent blanks (water:acetonitrile 80:20) before and after each set of samples. Contamination was also monitored by using a DI water sample blank spiked at ({20},hbox {ng}hbox { mL}^{-1}) IMI. In all cases, no IMI above MDLs was detected. Method precision was evaluated by duplicate injections and/or duplicate dilutions: the RPDs (relative percent differences) were all less than 15 %, demonstrating good method precision. Method accuracy was evaluated by analyzing a ({20},hbox {ng}hbox { mL}^{-1}) QC spike per set: the recoveries ranged between 96 % and 106 %, demonstrating good method accuracy.
Statistical analyses
All statistical analyses were conducted in R version 3.5.267. Birds exhibiting weight loss outside the lower bound of the 95 % confidence interval (CI) of 20 % over the study period ((n=3)) were omitted and data were reanalyzed. Birds exhibiting extreme weight loss across the pre-dose and post-dose conditions were in the control group ((n=2)) and the low dose group ((n=1)) suggesting a adverse response to the experimental period rather than the treatment itself. Data are presented as mean ± standard error. Significance ((p More

1.
Hu, Y. et al. Salinity and nutrient contents of tidal water affects soil respiration and carbon sequestration of high and low tidal flats of Jiuduansha wetlands in different ways. Sci. Total Environ. 565, 637–648 (2016).
ADS  CAS  PubMed  Article  Google Scholar 
2.
Stockmann, U. et al. The knowns, known unknowns and unknowns of sequestration of soil organic carbon. Agr. Ecosyst. Environ. 164, 80–99 (2013).
CAS  Article  Google Scholar 

3.
Yu, J., Wang, Y., Li, Y. & Dong, H. Soil organic carbon storage changes in coastal wetlands of the modern Yellow River Delta from 2000 to 2009. Biogeosciences 9, 2325–2331 (2012).
ADS  CAS  Article  Google Scholar 

4.
Wang, X., Luo, X., Jia, H. & Zheng, H. Dynamic characteristics of soil respiration in Yellow River Delta wetlands, China. Phys. Chem. Earth Parts A/B/C 103, 11–18 (2018).
ADS  Article  Google Scholar 

5.
Pendleton, L. et al. Estimating global “blue carbon” emissions from conversion and degradation of vegetated coastal ecosystems. PLoS ONE 7, 43542–43542 (2012).
ADS  Article  CAS  Google Scholar 

6.
Drake, P. L., Cormick, C. A. M. C. & Smith, M. J. Controls of soil respiration in a salinity-affected ephemeral wetland. Geoderma 221–222, 96–102 (2014).
ADS  Article  Google Scholar 

7.
Martin, D., Lal, T., Sachdev, C. B. & Sharma, J. P. Soil organic carbon storage changes with climate change, landform and land use conditions in Garhwal hills of the Indian Himalayan mountains. Agr. Ecosyst. Environ. 138, 64–73 (2010).
CAS  Article  Google Scholar 

8.
Vågen, T. G., Walsh, M. G. & Shepherd, K. D. Stable isotopes for characterisation of trends in soil carbon following deforestation and land use change in the highlands of Madagascar. Geoderma 135, 133–139 (2006).
ADS  Article  CAS  Google Scholar 

9.
Geissen, V. et al. Effects of land-use change on some properties of tropical soils: an example from Southeast Mexico. Geoderma 151, 87–97 (2009).
ADS  CAS  Article  Google Scholar 

10.
Fernández, S., Santín, C., Marquínez, J. & Álvarez, M. A. Saltmarsh soil evolution after land reclamation in Atlantic estuaries (Bay of Biscay, North coast of Spain). Geomorphology 114, 497–507 (2010).
ADS  Article  Google Scholar 

11.
Zhang, J., Song, C. & Yang, W. Effects of cultivation on soil microbiological properties in a freshwater marsh soil in Northeast China. Soil Tillage Res. 93, 231–235 (2007).
Article  Google Scholar 

12.
Erica, D., Ivana, J. & Luo, Y. On the variability of respiration in terrestrial ecosystems: moving beyond Q10. Glob. Change Biol. 12, 154–164 (2010).
Google Scholar 

13.
Martin, W. et al. Projected loss of soil organic carbon in temperate agricultural soils in the 21stcentury: effects of climate change and carbon input trends. Sci. Rep. 6, 32525 (2016).
ADS  Article  CAS  Google Scholar 

14.
Schulze, E. D. Biological control of the terrestrial carbon sink. Biogeosciences 2, 147–166 (2006).
ADS  Article  Google Scholar 

15.
Bond-Lamberty, B. & Thomson, A. Temperature-associated increases in the global soil respiration record. Nature 464, 579–582 (2010).
ADS  CAS  PubMed  Article  Google Scholar 

16.
Xu, M. & Shang, H. Contribution of soil respiration to the global carbon equation. J. Plant Physiol. 203, 16–28 (2016).
CAS  PubMed  Article  Google Scholar 

17.
Helingerová, M., Frouz, J. & Šantrůčková, H. Microbial activity in reclaimed and unreclaimed post-mining sites near Sokolov (Czech Republic). Ecol. Eng. 36, 768–776 (2010).
Article  Google Scholar 

18.
Kong, D. et al. Evolution of the Yellow River Delta and its relationship with runoff and sediment load from 1983 to 2011. J. Hydrol. 520, 157–167 (2015).
ADS  Article  Google Scholar 

19.
Han, G. et al. Winter soil respiration from different vegetation patches in the Yellow River Delta, China. Environ. Manag. 50, 39–49 (2012).
ADS  Article  Google Scholar 

20.
Scott, D., Baer, S. G. & Blair, J. M. Recovery and relative influence of root, microbial and structural properties of soil on physically sequestered carbon stocks in restored grassland. Soil Sci. Soc. Am. J. 81, 1–11 (2017).
Article  CAS  Google Scholar 

21.
Shi, W.-Y., Yan, M.-J., Zhang, J.-G., Guan, J.-H. & Du, S. Soil CO2 emissions from five different types of land use on the semiarid Loess Plateau of China, with emphasis on the contribution of winter soil respiration. Atmos. Environ. 88, 74–82 (2014).
ADS  CAS  Article  Google Scholar 

22.
Han, G. et al. Agricultural reclamation effects on ecosystem CO2 exchange of a coastal wetland in the Yellow River Delta. Agr. Ecosyst. Environ. 196, 187–198 (2014).
Article  Google Scholar 

23.
Zhao, Q. et al. Effects of water and salinity regulation measures on soil carbon sequestration in coastal wetlands of the Yellow River Delta. Geoderma 319, 219–229 (2018).
ADS  CAS  Article  Google Scholar 

24.
Chen, S. et al. Modeling interannual variability of global soil respiration from climate and soil properties. Agric. For. Meteorol. 150, 590–605 (2010).
ADS  Article  Google Scholar 

25.
Iost, S., Landgraf, D. & Makeschin, F. Chemical soil properties of reclaimed marsh soil from Zhejiang Province P.R. China. Geoderma 142, 245–250 (2007).
ADS  CAS  Article  Google Scholar 

26.
Bu, N.-S. et al. Reclamation of coastal salt marshes promoted carbon loss from previously-sequestered soil carbon pool. Ecol. Eng. 81, 335–339 (2015).
Article  Google Scholar 

27.
Ceschia, E. et al. Management effects on net ecosystem carbon and GHG budgets at European crop sites. Agr. Ecosyst. Environ. 139, 363–383 (2010).
Article  Google Scholar 

28.
Chen, Q., Guo, B., Zhao, C. & Xing, B. Characteristics of CH4 and CO2 emissions and influence of water and salinity in the Yellow River delta wetland China. Environ. Pollut. 239, 289–299 (2018).
CAS  PubMed  Article  Google Scholar 

29.
Melillo, J. M. et al. Long-term pattern and magnitude of soil carbon feedback to the climate system in a warming world. Science 358, 101–105 (2017).
ADS  CAS  PubMed  Article  Google Scholar 

30.
Xie, J., Li, Y., Zhai, C., Li, C. & Lan, Z. CO2 absorption by alkaline soils and its implication to the global carbon cycle. Environ. Geol. 56, 953–961 (2009).
ADS  CAS  Article  Google Scholar 

31.
Wang, C., Yang, J. & Zhang, Q. Soil respiration in six temperate forests in China. Glob. Change Biol. 12, 2103–2114 (2010).
ADS  Article  Google Scholar 

32.
Bahn, M. et al. Soil respiration at mean annual temperature predicts annual total across vegetation types and biomes. Biogeosciences 7, 2147–2157 (2010).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

33.
Liu, X. et al. Diurnal variation in soil respiration under different land uses on Taihang Mountain North China. Atmos. Environ. 125, 283–292 (2016).
ADS  CAS  Article  Google Scholar 

34.
Li, W., Wang, J., Zhang, X., Shi, S. & Cao, W. Effect of degradation and rebuilding of artificial grasslands on soil respiration and carbon and nitrogen pools on an alpine meadow of the Qinghai-Tibetan Plateau. Ecol. Eng. 111, 134–142 (2018).
CAS  Article  Google Scholar 

35.
Kuzyakov, Y. & Gavrichkova, O. Time lag between photosynthesis and carbon dioxide efflux from soil: a review of mechanisms and controls. Glob. Change Biol. 16, 3386–3406 (2010).
ADS  Article  Google Scholar 

36.
Bradford, M. A. et al. Thermal adaptation of soil microbial respiration to elevated temperature. Ecol. Lett. 11, 1316–1327 (2008).
PubMed  Article  Google Scholar 

37.
Jiang, J. et al. Changes in temperature sensitivity of soil respiration in the phases of a three-year crop rotation system. Soil & Tillage Research 150, 139–146 (2015).
Article  Google Scholar 

38.
Baldrian, P. et al. Responses of the extracellular enzyme activities in hardwood forest to soil temperature and seasonality and the potential effects of climate change. Soil Biol. Biochem. 56, 60–68 (2013).
CAS  Article  Google Scholar 

39.
Wang, R. et al. Effects of crop types and nitrogen fertilization on temperature sensitivity of soil respiration in the semi-arid Loess Plateau. Soil & Tillage Research 163, 1–9 (2016).
ADS  Article  Google Scholar 

40.
Janssens, I. A. & Pilegaard, K. Large seasonal changes in Q10 of soil respiration in a beech forest. Glob. Change Biol. 9, 911–918 (2003).
ADS  Article  Google Scholar 

41.
Fang, Q., Wang, G., Liu, T. & Xue, B.-L. Controls of carbon flux in a semi-arid grassland ecosystem experiencing wetland loss: vegetation patterns and environmental variables. Agric. For. Meteorol. 259, 196–210 (2018).
ADS  Article  Google Scholar 

42.
Knowles, J. F., Blanken, P. D. & Williams, M. W. Soil respiration variability across a soil moisture and vegetation community gradient within a snow-scoured alpine meadow. Biogeochemistry 125, 185–202 (2015).
CAS  Article  Google Scholar 

43.
Arredondo, T. et al. Does precipitation affects soil respiration of tropical semiarid grasslands with different plant cover types?. Agr. Ecosyst. Environ. 251, 218–225 (2018).
Article  Google Scholar 

44.
González-Ubierna, S. & Lai, R. Modelling the effects of climate factors on soil respiration across Mediterranean ecosystems. J. Arid Environ. 165, 46–54 (2019).
ADS  Article  Google Scholar 

45.
Curiel Yuste, J., Janssens, I. A., Carrara, A., Meiresonne, L. & Ceulemans, R. Interactive effects of temperature and precipitation on soil respiration in a temperate maritime pine forest. Tree Physiol. 23, 1263–1270 (2003).
CAS  PubMed  Article  Google Scholar 

46.
Palmroth, S. et al. Contrasting responses to drought of forest floor CO2 efflux in a loblolly pine plantation and a nearly oak-hickory forest. Glob. Change Biol. 11, 1–14 (2005).
Article  Google Scholar 

47.
Guo, H. et al. Annual ecosystem respiration of maize was primarily driven by crop growth and soil water conditions. Agr. Ecosyst. Environ. 272, 254–265 (2019).
Article  Google Scholar 

48.
Zhang, L. H., Chen, Y. N., Zhao, R. F. & Li, W. H. Significance of temperature and soil water content on soil respiration in three desert ecosystems in Northwest China. J. Arid Environ. 74, 1200–1211 (2010).
ADS  Article  Google Scholar 

49.
Herbst, M., Tappe, W., Kummer, S. & Vereecken, H. The impact of sieving on heterotrophic respiration response to water content in loamy and sandy topsoils. Geoderma 272, 73–82 (2016).
ADS  Article  Google Scholar 

50.
Tan, Z. H. et al. Soil respiration in an old-growth subtropical forest: Patterns, components, and controls. J. Geophys. Res. Atmos. 118, 2981–2990 (2013).
ADS  Article  Google Scholar 

51.
Lu, D. D. & Shi, D. L. Can rock fragment cover maintain soil and water for saline-sodic soil slopes under coastal reclamation?. CATENA 151, 213–224 (2017).
Article  Google Scholar 

52.
Xie, X. F. et al. Differential effects of various reclamation treatments on soil characteristics: an experimental study of newly reclaimed tidal mudflats on the east China coast. Sci. Total Environ. 768, 144996 (2021).
ADS  CAS  PubMed  Article  Google Scholar 

53.
Almagro, M., López, J., Querejeta, J. I. & Martínez-Mena, M. Temperature dependence of soil CO2 efflux is strongly modulated by seasonal patterns of moisture availability in a Mediterranean ecosystem. Soil Biol. Biochem. 41, 594–605 (2009).
CAS  Article  Google Scholar 

54.
González-Ubierna, S., Teresa de la Cruz, M. & Casermeiro, M. Á. Climate factors mediate soil respiration dynamics in Mediterranean agricultural environments: an empirical approach. Soil Res. 52, 543–553 (2014).
Article  Google Scholar 

55.
Schütt, M., Borken, W., Spott, O., Stange, C. F. & Matzner, E. Temperature sensitivity of C and N mineralization in temperate forest soils at low temperatures. Soil Biol. Biochem. 69, 320–327 (2014).
Article  CAS  Google Scholar 

56.
Zhang, Q., Meng, W. Q. & Li, H. Y. Research of the diurnal soil respiration dynamic in two typical vegetation communities in Tianjin estuarine wetland. Earth Environ. Sci. 69, 1–8 (2016).
Google Scholar 

57.
Eberwein, J. R., Oikawa, P. Y., Allsman, L. A. & Jenerette, G. D. Carbon availability regulates soil respiration response to nitrogen and temperature. Soil Biol. Biochem. 88, 158–164 (2015).
CAS  Article  Google Scholar 

58.
Bujalský, L., Kaneda, S., Dvorščík, P. & Frouz, J. In situ soil respiration at reclaimed and unreclaimed post-mining sites: Responses to temperature and reclamation treatment. Ecol. Eng. 68, 53–59 (2014).
Article  Google Scholar 

59.
Al-Busaidi, K. T. S., Buerkert, A. & Joergensen, R. G. Carbon and nitrogen mineralization at different salinity levels in Omani low organic matter soils. J. Arid Environ. 100–101, 106–110 (2014).
Article  Google Scholar 

60.
Mavi, M. S., Marschner, P., Chittleborough, D. J., Cox, J. W. & Sanderman, J. Salinity and sodicity affect soil respiration and dissolved organic matter dynamics differentially in soils varying in texture. Soil Biol. Biochem. 45, 8–13 (2012).
CAS  Article  Google Scholar 

61.
FAO Soils Portal. Retrieved from http://www.fao.org/soils-portal/soil-management/management-of-some-problem-soils/salt-affected-soils/more-information-on-salt-affected-soils/en/ (2019)

62.
Wu, J., Joergensen, R. G., Pommerening, B., Chaussod, R. & Brookes, P. C. Measurement of soil microbial biomass C by fumigation-extraction: an automated procedure. Soil Biol. Biochem. 22, 1167–1169 (1990).
CAS  Article  Google Scholar 

63.
Bai, J. et al. Spatial distribution characteristics of organic matter and total nitrogen of marsh soils in river marginal wetlands. Geoderma 124, 181–192 (2005).
ADS  CAS  Article  Google Scholar 

64.
Campbell, J. L., Sun, O. J. & Law, B. E. Supply-side controls on soil respiration among Oregon forests. Glob. Change Biol. 10, 1857–1869 (2004).
ADS  Article  Google Scholar 

65.
Zou, J., Tobin, B., Luo, Y. & Osborne, B. Response of soil respiration and its components to experimental warming and water addition in a temperate Sitka spruce forest ecosystem. Agric. For. Meteorol. 260–261, 204–215 (2018).
ADS  Article  Google Scholar 

66.
Li, L. H., Wang, Q. B., Bai, Y. F., Zhou, G. S. & Xing, X. R. Soil respiration of a leymus chinese grassland stand in the Xilin river basin as affected by over-grazing and climate. Acta Phytoecologica Sinica 6, 41–47 (2000).
Google Scholar 

67.
Moinet, G. Y. K. et al. Soil heterotrophic respiration is insensitive to changes in soil water content but related to microbial access to organic matter. Geoderma 274, 68–78 (2016).
ADS  Article  Google Scholar  More

1.
Ineson, P., Leonard, M. A. & Anderson, J. M. Effect of collembolan grazing upon nitrogen and cation leaching from decomposing leaf litter. Soil Biol. Biochem. 14, 601–605. https://doi.org/10.1016/0038-0717(82)90094-3 (1982).
Article  Google Scholar 
2.
Petersen, H. & Luxton, M. A. comparative analysis of soil fauna populations and their role in decomposition processes. Oikos 39, 288–388. https://doi.org/10.1016/j.pedobi.2006.08.006 (1982).
Article  Google Scholar 

3.
Drake, H. L. & Horn, M. A. As the worm turns: The earthworm gut as a transient habitat for soil microbial biomes. Annu. Rev. Microbiol. 61, 169–189. https://doi.org/10.1146/annurev.micro.61.080706.093139 (2007).
CAS  Article  PubMed  Google Scholar 

4.
Liu, Y. et al. Higher soil fauna abundance accelerates litter carbon release across an alpine forest-tundra ecotone. Sci. Rep. 9, 10562. https://doi.org/10.1038/s41598-019-47072-0 (2019).
CAS  Article  Google Scholar 

5.
Hopkin, S. P. Biology of the Springtails (Insecta: Collembola) (Oxford University Press, Oxford, 1997).
Google Scholar 

6.
Maaß, S., Caruso, T. & Rillig, M. C. Functional role of microarthropods in soil aggregation. Pedobiologia 58, 59–63. https://doi.org/10.1016/j.pedobi.2015.03.001 (2015).
Article  Google Scholar 

7.
Bergstrom, D. M., Convey, P. & Huiskes, A. H. L. Trends in Antarctic Terrestrial and Limnetic Ecosystems: Antarctica as a Global Indicator (Springer, Berlin, 2006). .
Google Scholar 

8.
Convey, P. Antarctic terrestrial biodiversity in a changing world. Polar. Biol. 34(11), 1629–1641. https://doi.org/10.1007/s00300-011-1068-0 (2011).
Article  Google Scholar 

9.
Convey, P. et al. The spatial structure of Antarctic biodiversity. Ecol. Monogr. 84(2), 203–244. https://doi.org/10.1890/12-2216.1 (2014).
Article  Google Scholar 

10.
Wauchope, H. S., Shaw, J. D. & Terauds, A. A snapshot of biodiversity protection in Antarctica. Nat. Commun. 10(1), 946. https://doi.org/10.1038/s41467-019-08915-6 (2019).
ADS  CAS  Article  PubMed  PubMed Central  Google Scholar 

11.
Chown, S. L. et al. The changing form of Antarctic biodiversity. Nature 522(7557), 431–438. https://doi.org/10.1038/nature14505 (2015).
ADS  CAS  Article  PubMed  Google Scholar 

12.
Agamennone, V. et al. The microbiome of Folsomia candida: An assessment of bacterial diversity in a Wolbachia-containing animal. FEMS Microbiol. Ecol. 91(11), 1–10. https://doi.org/10.1093/femsec/fiv128 (2015).
CAS  Article  Google Scholar 

13.
Zhu, D. et al. Exposure of soil collembolans to microplastics perturbs their gut microbiota and alters their isotopic composition. Soil Biol. Biochem. 115, 302–310. https://doi.org/10.1016/j.soilbio.2017.10.027 (2018).
CAS  Article  Google Scholar 

14.
Bahrndorff, S. et al. Diversity and metabolic potential of the microbiota associated with a soil arthropod. Sci. Rep. 8(1), 1–8. https://doi.org/10.1038/s41598-018-20967-0 (2018).
CAS  Article  Google Scholar 

15.
Ding, J. et al. Effects of long-term fertilization on the associated microbiota of soil collembolan. Soil Biol. Biochem. 130, 141–149. https://doi.org/10.1016/j.soilbio.2018.12.015 (2019).
CAS  Article  Google Scholar 

16.
Anslan, S., Bahram, M. & Tedersoo, L. Temporal changes in fungal communities associated with guts and appendages of Collembola as based on culturing and high-throughput sequencing. Soil Biol. Biochem. 96, 152–159. https://doi.org/10.1016/j.soilbio.2016.02.006 (2016).
CAS  Article  Google Scholar 

17.
Terauds, A. et al. Conservation biogeography of the Antarctic. Divers Distrib. 18(7), 726–741. https://doi.org/10.1111/j.1472-4642.2012.00925.x (2012).
Article  Google Scholar 

18.
Terauds, A. & Lee, J. R. Antarctic biogeography revisited: Updating the Antarctic Conservation Biogeographic Regions. Divers Distrib. 22(8), 836–840. https://doi.org/10.1111/ddi.12453 (2016).
Article  Google Scholar 

19.
Greenslade, P. An Antarctic biogeographical anomaly resolved: The true identity of a widespread species of Collembola. Polar Biol. 41(5), 969–981. https://doi.org/10.1007/s00300-018-2261-1 (2018).
Article  Google Scholar 

20.
Carapelli, A. et al. Evidence for cryptic diversity in the “pan-Antarctic” springtail Friesea antarctica and the description of two new species. Insects 11, 141. https://doi.org/10.3390/insects11030141 (2020).
Article  PubMed Central  Google Scholar 

21.
Carapelli, A., Convey, P., Frati, F., Spinsanti, G. & Fanciulli, P. P. Population genetics of three sympatric springtail species (Hexapoda: Collembola) from the South Shetland Islands: Evidence for a common biogeographic pattern. Biol. J. Linn. Soc. 120, 788–803. https://doi.org/10.1093/biolinnean/blw004 (2017).
Article  Google Scholar 

22.
Collins, G. E., Hogg, I. D., Convey, P., Barnes, A. D. & McDonald, I. R. Spatial and temporal scales matter when assessing the species and genetic diversity of springtails (Collembola) in Antarctica. Front. Ecol. Evol. 7, 76. https://doi.org/10.3389/fevo.2019.00076 (2019).
Article  Google Scholar 

23.
Collins, G. E. et al. Genetic diversity of soil invertebrates corroborates timing estimates for past collapses of the West Antarctic Ice Sheet. PNAS 117, 22293–22302. https://doi.org/10.1073/pnas.2007925117 (2020).
ADS  CAS  Article  PubMed  Google Scholar 

24.
Holmes, C. J. et al. The Antarctic mite, Alaskozetes antarcticus, shares bacterial microbiome community membership but not abundance between adults and tritonymphs. Polar Biol. 42, 2075–2085. https://doi.org/10.1007/s00300-019-02582-5 (2019).
Article  Google Scholar 

25.
Vecchi, M., Newton, I. L. G., Cesari, M., Rebecchi, L. & Guidetti, R. The microbial community of tardigrades: Environmental influence and species specificity of microbiome structure and composition. Microb. Ecol. 76(2), 467–481. https://doi.org/10.1007/s00248-017-1134-4 (2018).
CAS  Article  PubMed  Google Scholar 

26.
Delgado-Baquerizo, M. et al. Ecological drivers of soil microbial diversity and soil biological networks in the Southern Hemisphere. Ecology 99(3), 583–596. https://doi.org/10.1002/ecy.2137 (2018).
Article  PubMed  Google Scholar 

27.
Chu, H. et al. Soil bacterial diversity in the Arctic is not fundamentally different from that found in other biomes. Environ. Microbiol. 12(11), 2998–3006. https://doi.org/10.1111/j.1462-2920.2010.02277.x (2010).
CAS  Article  PubMed  Google Scholar 

28.
Siciliano, S. D. et al. Soil fertility is associated with fungal and bacterial richness, whereas pH is associated with community composition in polar soil microbial communities. Soil Biol. Biochem. 78, 10–20. https://doi.org/10.1016/j.soilbio.2014.07.005 (2014).
CAS  Article  Google Scholar 

29.
Zouache, K. et al. Composition of bacterial communities associated with natural and laboratory populations of Asobara tabida infected with Wolbachia. Appl. Environ. Microb. 75, 3755–3764. https://doi.org/10.1128/aem.02964-08 (2009).
CAS  Article  Google Scholar 

30.
Potapov, A. A., Semenina, E. E., Korotkevich, A. Y., Kuznetsova, N. A. & Tiunov, A. V. Connecting taxonomy and ecology: Trophic niches of collembolans as related to taxonomic identity and life forms. Soil Biol. Biochem. 101, 20–31. https://doi.org/10.1016/j.soilbio.2016.07.002 (2016).
CAS  Article  Google Scholar 

31.
De Wever, A. et al. Hidden levels of phylodiversity in Antarctic green algae: Further evidence for the existence of glacial refugia. Proc. R. Soc. B 276, 3591–3599. https://doi.org/10.1098/rspb.2009.0994 (2009).
Article  PubMed  Google Scholar 

32.
Vyverman, W. et al. Evidence for widespread endemism among Antarctic micro-organisms. Polar Sci. 4(2), 103–113. https://doi.org/10.1016/j.polar.2010.03.006 (2010).
ADS  Article  Google Scholar 

33.
Finlay, B. J. & Clarke, K. J. Ubiquitous dispersal of microbial species. Nature 400, 828–828. https://doi.org/10.1038/23616 (1999).
ADS  CAS  Article  Google Scholar 

34.
Chown, S. L. & Convey, P. Structure and temporal variability across life’s hierarchies in the terrestrial Antarctic. Philos. Trans. R. Soc. B 362, 2307–23331. https://doi.org/10.1098/rstb.2006.1949 (2007).
Article  Google Scholar 

35.
Convey, P., Biersma, E. M., Casanova-Katny, A. & Maturana, C. S. Refuges of Antarctic diversity. Chapter 10. In Past Antarctica (eds Oliva, M. & Ruiz-Fernández, J.) 181–200 (Academic Press, Burlington, 2020). https://doi.org/10.1016/B978-0-12-817925-3.00010-0.
Google Scholar 

36.
Bolger, A. M., Lohse, M. & Usadel, B. Trimmomatic: A flexible trimmer for Illumina sequence data. Bioinformatics 30(15), 2114–2120. https://doi.org/10.1093/bioinformatics/btu170 (2014).
CAS  Article  PubMed  PubMed Central  Google Scholar 

37.
Bolyen, E. et al. Reproducible, interactive, scalable and extensible microbiome data science using QIIME 2. Nat. Biotechnol. 37, 852–857. https://doi.org/10.1038/s41587-019-0209-9 (2019).
CAS  Article  PubMed  PubMed Central  Google Scholar 

38.
Callahan, B. J. et al. DADA2: High-resolution sample inference from Illumina amplicon data. Nat. Methods 13(7), 581–583. https://doi.org/10.1038/nmeth.3869 (2016).
CAS  Article  PubMed  PubMed Central  Google Scholar 

39.
Rognes, T., Flouri, T., Nichols, B., Quince, C. & Mahé, F. VSEARCH: A versatile open source tool for metagenomics. PeerJ 4, e2584. https://doi.org/10.7717/peerj.2584 (2016).
Article  PubMed  PubMed Central  Google Scholar 

40.
Bokulich, N. A. et al. Optimizing taxonomic classification of marker-gene amplicon sequences with qiime 2’s q2-feature-classifier plugin. Microbiome 6(1), 90. https://doi.org/10.1186/s40168-018-0470-z (2018).
MathSciNet  Article  PubMed  PubMed Central  Google Scholar 

41.
Price, M. N., Dehal, P. S. & Arkin, A. P. Fasttree 2-approximately maximum-likelihood trees for large alignments. PLoS One 5(3), e9490. https://doi.org/10.1371/journal.pone.0009490 (2010).
ADS  CAS  Article  PubMed  PubMed Central  Google Scholar 

42.
Lahti, L. & Shetty, S. Microbiome R package. http://microbiome.github.io (2012–2019).

43.
Ssekagiri, A., Sloan, W. T. & Ijaz, U. Z. microbiomeSeq: An R package for analysis of microbial communities in an environmental context. ISCB Africa ASBCB Conference. http://www.github.com/umerijaz/microbiomeSeq (2017).

44.
McMurdie, P. J. & Holmes, S. phyloseq: An R package for reproducible interactive analysis and graphics of microbiome census data. PLoS One 8(4), e61217. https://doi.org/10.1371/journal.pone.0061217 (2014).
ADS  CAS  Article  Google Scholar 

45.
Oksanen, J., et al. Vegan: Community ecology package. R package version 2.5-6. https://github.com/vegandevs/vegan (2019).

46.
Ritchie, M. E. et al. limma powers differential expression analyses for RNA-sequencing and microarray studies. Nucleic Acids Res. 43(7), e47. https://doi.org/10.1093/nar/gkv007 (2015).
CAS  Article  PubMed  PubMed Central  Google Scholar 

47.
Chen, H. & Boutros, P. C. VennDiagram: A package for the generation of highly-customizable Venn and Euler diagrams in R. BMC Bioinf. 12, 35. https://doi.org/10.1186/1471-2105-12-35 (2011).
Article  Google Scholar 

48.
Wickham, H. ggplot2: Elegant Graphics for Data Analysis. Springer, New York. https://ggplot2.tidyverse.org (2016).

49.
Warnes, G. R., et al. gplots: Various R Programming Tools for Plotting Data. R package version 3.0.1.1. https://CRAN.R-project.org/package=gplots (2019). More

1.
Arrigo, K. R. & Dijken, G. L. Secular trends in Arctic Ocean net primary production. J. Geophys. Res. Oceans. 116, C09011 (2011).
ADS  Google Scholar 
2.
Thomas, D. N. Sea Ice Ch 4 (Wiley Blackwell, Oxford, 2017).
Google Scholar 

3.
Arrigo, K. R. et al. Massive phytoplankton blooms under Arctic sea ice. Science 336, 1408–1408 (2012).
ADS  CAS  Article  Google Scholar 

4.
Assmy, P. et al. Leads in Arctic pack ice enable early phytoplankton blooms below snow-covered sea ice. Sci. Rep. 7, 40850 (2016).
ADS  Article  Google Scholar 

5.
Horvat, C. et al. The frequency and extent of sub-ice phytoplankton bloom in the Arctic Ocean. Sci. Adv. 3, e1601191 (2017).
ADS  Article  Google Scholar 

6.
Ardyna, M. et al. Environmental drivers of under-ice phytoplankton bloom dynamics in the Arctic Ocean. Elem. Sci. Anth. 8, 30 (2020).
Article  Google Scholar 

7.
Ardyna, M. et al. Under-ice phytoplankton blooms: Shedding light on the “invisible” part of Arctic primary production. Front. Mar. Sci. 7, 608032 (2020).
Article  Google Scholar 

8.
Rysgaard, S. & Glud, R. N. Carbon cycling in Arctic marine ecosystems: Case study Young Sound (ed. Rysgaard, S. & Glud, R. N.) 62–94 (Meddelelser om Grønland, Bioscience Vol 58, Copenhagen, Denmark, the Commission for Scientific Research in Greenland, 2007).

9.
Meire, L. et al. Marine-terminating glaciers sustain high productivity in Greenland fjords. Glob. Chang. Biol. 23, 5344–5357 (2017).
ADS  Article  Google Scholar 

10.
Randelhoff, A. et al. Pan-Arctic Ocean primary production constrained by turbulent nitrate fluxes. Front. Mar. Sci. 7, 150 (2020).
Article  Google Scholar 

11.
Holding, J. M. et al. Seasonal and spatial patterns of primary production in a high-latitude fjord affected by Greenland Ice Sheet run-off. Biogeosciences 16, 3777–3792 (2019).
ADS  CAS  Article  Google Scholar 

12.
Juul-Pedersen, T. et al. Seasonal and interannual phytoplankton production in a sub-Arctic tidewater outlet glacier fjord, SW Greenland. Mar. Ecol. Prog. Ser. 524, 27–38 (2015).
ADS  Article  Google Scholar 

13.
Sejr, M. K. et al. Evidence of local and regional freshening of Northeast Greenland coastal waters. Sci. Rep. 7, 13183 (2017).
ADS  Article  Google Scholar 

14.
Boone, W. et al. Circulation and fjord-shelf exchange during the ice-covered period in Young Sound-Tyrolerfjord, Northeast Greenland (74°N). Estuar. Coast. Shelf Sci. 194, 205–216 (2017).
ADS  Article  Google Scholar 

15.
Haine, T. W. N. et al. Arctic freshwater export: Status, mechanisms, and prospects. Glob. Planet Change 125, 13–35 (2015).
ADS  Article  Google Scholar 

16.
Carmack, E. C. et al. Freshwater and its role in the Arctic Marine System: Sources, disposition, storage export, and physical and biogeochemical consequences in the Arctic and global ocean. J. Geophys. Res. Biogeosci. 121, 675–717 (2015).
Article  Google Scholar 

17.
Lund-Hansen, L. C. et al. Will low primary production rates in the Amundsen Basin (Arctic Ocean) remain low in a future ice-free setting, and what governs this production?. J. Mar. Syst. 205, 103287 (2020).
Article  Google Scholar 

18.
Dahl, E., Bagøien, E., Edvardsen, B. & Stenseth, N. C. The dynamics of Chrysochromulina species in the Skagerrak in relation to environmental conditions. J. Sea. Res. 54, 15–24 (2005).
ADS  Article  Google Scholar 

19.
Hansen, P. J., Nielsen, T. G. & Kaas, H. Distribution and growth of protists and mesozooplankton during a bloom of Chrysochromulina spp. (Prymnesiophyceae, Prymnesiales). Phycologia 34, 409–416 (1995).
Article  Google Scholar 

20.
Nielsen, T. G., Kiørboe, T. & Bjørnsen, P. K. Effects of a Chrysochromulina polylepis subsurface bloom on the planktonic community. Mar. Ecol. Prog. Ser. 62, 21–35 (1990).
ADS  Article  Google Scholar 

21.
Hällfors, G. & Niemi, Å. A Chrysochromulina (Haptophyceae) bloom under the ice in the Tvärminne Archipelago, southern coast of Finland. Acta Soc. Fauna Flora Fenn. 50, 89–104 (1974).
Google Scholar 

22.
Manton, I. Chrysochromulina tenuispine sp. nov. from arctic Canada. Br. Phycol. J. 13, 227–234 (1978).
Article  Google Scholar 

23.
Green, J. C. & Leadbeater, B. S. C. The Haptophyte Algae ch. 13 (Systematics Association, London, 1994).
Google Scholar 

24.
Hansen, P. J. & Hjorth, M. Growth and grazing responses of Chrysochromulina ericina (Prymnesiophyceae): The role of irradiance, prey concentration and pH. Mar. Biol. 141, 975–983 (2002).
CAS  Article  Google Scholar 

25.
Anderson, R., Charvet, S. & Hansen, P. J. Mixotrophy in chlorophytes and haptophytes—Effect of irradiance, macronutrient micronutrient and vitamin limitation. Front. Microbiol. 9, 1704 (2018).
Article  Google Scholar 

26.
Anderson, R. & Hansen, P. J. Meteorological conditions induce strong shifts in mixotrophic and heterotrophic flagellate bacterivory over small spatio-temporal scales. Limnol. Oceanogr. 9999, 1–11 (2019).
Google Scholar 

27.
McKie-Krisberg, Z. M., Gast, R. J. & Sanders, R. W. Physiological responses of three species of Antarctic mixotrophic phytoflagellates to changes in light and dissolved nutrients. Microbiol. Ecol 70, 21–29 (2015).
CAS  Article  Google Scholar 

28.
McKie-Krisberg, Z. M., Sanders, R. W. & Gast, R. J. Evaluation of mixotrophy-associated gene expression in two species of polar marine algae. Front. Mar. Sci. 5, 273 (2018).
Article  Google Scholar 

29.
Rysgaard, S., Nielsen, T. G. & Hansen, B. W. Seasonal variation in nutrients, pelagic primary production and grazing in a high-Arctic coastal marine ecosystem, Young Sound, Northeast Greenland. Mar. Ecol. Prog. Ser. 179, 13–25 (1999).
ADS  CAS  Article  Google Scholar 

30.
Bendtsen, J., Mortensen, J. & Rysgaard, S. Seasonal surface layer dynamics and sensitivity to runoff in a high Arctic fjord (Young Sound/Tyrolerfjord, 74°N). J. Geophys. Res. Oceans. 119, 1–18 (2014).
Article  Google Scholar 

31.
Krawczyk, D. W. et al. Spatial and temporal distribution of planktonic protists in the East Greenland fjord and offshore waters. Mar. Ecol. Prog. Ser. 538, 99–116 (2015).
ADS  CAS  Article  Google Scholar 

32.
Søgaard, D. H., Deming, J. W., Meire, L. & Rysgaard, S. Effects of microbial processes and CaCO3 dynamics on inorganic carbon cycling in snow-covered Arctic winter sea ice. Mar. Ecol. Prog. Ser. 611, 31–44 (2019).
ADS  Article  Google Scholar 

33.
Rysgaard, S. et al. Ikaite crystal distribution in winter sea ice and implications for CO2 system dynamics. Cryosphere 7, 707–718 (2013).
ADS  Article  Google Scholar 

34.
Søgaard, D. H. et al. Autotrophic and heterotrophic activity in Arctic first-year sea ice: Seasonal study from Malene Bight, SW Greenland. Mar. Ecol. Prog. Ser. 419, 31–45 (2010).
ADS  Article  Google Scholar 

35.
Grasshoff, K., Kremling, K. & Ehrhardt, M. Methods of Seawater Analysis (WILEY-VCH Verlag GmbH, Weinheim, 1999).
Google Scholar 

36.
Steemann-Nielsen, E. The use of radio-active carbon (C14) for measuring organic production in the sea. ICES J. Mar. Sci. 18, 117–140 (1952).
Article  Google Scholar 

37.
Søgaard, D. H. et al. The relative contributions of biological and abiotic processes to carbon dynamics in subarctic sea ice. Polar Biol. 36, 1761–1777 (2013).
Article  Google Scholar 

38.
Platt, T., Gallegos, C. L. & Harrison, W. G. Photoinhibition of photosynthesis in natural assemblages of marine phytoplankton. J. Mar. Res. 38, 687–701 (1980).
Google Scholar 

39.
Jespersen, A. M. & Christoffersen, K. Measurements of chlorophyll-a from phytoplankton using ethanol as extraction solvent. Arch. Hydrobiol. 109, 445–454 (1987).
CAS  Google Scholar 

40.
Ralph, P. J. & Gademann, R. Rapid light curves: A powerful tool to assess photosynthetic activity. Aquat. Bot. 82, 222–237 (2005).
CAS  Article  Google Scholar 

41.
Jassby, A. D. & Platt, T. Mathematical formulation of the relationship between photosynthesis and light for phytoplankton. Limnol. Oceanogr. 21, 540–547 (1976).
ADS  CAS  Article  Google Scholar  More

Reported SARS-CoV-2 case counts, mortality and testing in SSA as of December 2020
Variables and data sources for reporting data
The numbers of reported cases, deaths and tests for the 48 SSA countries studied (Supplementary Table 1) were sourced from the Africa CDC dashboard on 20 December 2020 (and previously on 23 September and 30 June 2020). The Africa CDC obtains data from the official Africa CDC Regional Collaborating Centre and member state reports. Differences in the timing of reporting by member states results in some variation in the recency of data within the centralized Africa CDC repository, but data should broadly reflect the relative scale of testing and reporting efforts across countries. For Mauritius (https://covid19.mu/) and Rwanda (https://covid19.who.int/region/afro/country/rw), reporting to the Africa CDC was confirmed by comparison to country-specific dashboards.
The countries or member states within SSA in this study follow the United Nations and Africa CDC-listed regions of Southern, Western, Central and Eastern Africa (excluding Sudan). From the Northern Africa region, Mauritania is included in SSA.
For comparison to non-SSA countries, the number of reported cases in other geographical regions were obtained from the Johns Hopkins University Coronavirus Resource Center on 23 September 2020 (https://coronavirus.jhu.edu/map.html).
Case fatality ratios (CFRs) were calculated by dividing the number of reported deaths by the number of reported cases and expressed as a percentage. Positivity was calculated by dividing the number of reported cases by the number of reported tests. Testing and case rates were calculated per 100,000 population using population size estimates for 2020 from the United Nations Population Division (https://population.un.org/wpp/Download/Standard/Population/). Since reported confirmed cases are likely to be an underestimate of the true number of infections, CFRs may be a poor proxy for the IFR, defined as the proportion of infections that result in mortality4.
Variation in testing and mortality rates
Testing rates among SSA countries varied by multiple orders of magnitude as of 30 June and remain highly variable as of 23 September and 20 December 2020. The number of tests completed per 100,000 population ranged from 19.84 in Burundi to 13,508.13 in Mauritius in June 2020; from 65.98 in the Democratic Republic of the Congo to 18,321.83 in Mauritius in September 2020; and from 100.9 in the Democratic Republic of the Congo to 23695.0 in Mauritius in December 2020 (Extended Data Fig. 1a). Tanzania (6.50 tests per 100,000 population) has not reported new tests, cases or deaths to the Africa CDC since April 2020. The number of reported infections (that is, positive tests) was strongly correlated with the number of tests completed in June 2020 (Pearson’s correlation coefficient, r = 0.9667, P 50% coincident with population >50,000) within administrative 2 units61. For some countries, estimates at administrative 2 units were unavailable (Comoros, Cape Verde, Lesotho, Mauritius, Mayotte and Seychelles); estimates at the administrative-1 unit level were used for these cases (these were all island nations, with the exception of Lesotho).
Metapopulation model methods
Once SARS-CoV-2 has been introduced into a country, the degree of spread of the infection within the country is governed by subnational mobility: the pathogen is more likely to be introduced into a location where individuals arrive more frequently than one where incoming travelers are less frequent. Large-scale consistent measures of mobility are rare. However, recently, estimates of accessibility have been produced at a global scale26. Although this is unlikely to perfectly reflect mobility within countries, especially since interventions and travel restrictions are put in place, it provides a starting point for evaluating the role of human mobility in shaping the outbreak pace across SSA. We used the inverse of a measure of the cost of travel between the centroids of administrative level 2 spatial units to describe mobility between locations (estimated by applying the costDistance function in the gdistance package v1.3-6 in R to the friction surfaces supplied in Weiss et al.26). With this, we developed a metapopulation model for each country to develop an overview of the possible range of trajectories of unchecked spread of SARS-CoV-2.
We assumed that the pathogen first arrives in each country in the administrative 2 level unit with the largest population (for example, the largest city) and the population in each administrative 2 level (of size Nj) is entirely susceptible at the time of arrival. We then tracked the spread within and between each of the administrative 2 level units of each country. Within each administrative 2 level unit, dynamics are governed by a discrete time susceptible (S), infected (I) and recovered (R) model with a time step of approximately one week, which is broadly consistent with the serial interval of SARS-CoV-2. Within the spatial unit indexed j, with total size Nj, the number of infected individuals in the next time step is defined by:

$$I_{j,t + 1} = beta I_{j,t}^alpha S_{j,t}/N_j + iota _{j,t}$$

where β captures the magnitude of transmission over the course of one discrete time step; since the discrete time step chosen is set to approximate the serial interval of the virus, this will reflect the R0 of SARS-CoV-2, and is thus set to 2.5; the exponent α = 0.97 is used to capture the effects of discretization62 and Ij,t captures the introduction of new infections into site j at time t. Susceptible and recovered individuals are updated according to:

$$begin{array}{l}S_{j,t + 1} = S_{j,t} + wR_{j,t} – I_{j,t + 1} + b\ R_{j,t + 1} = (1 – w)R_{j,t} + I_{j,t}end{array}$$

where b reflects the introduction of new susceptible individuals resulting from the birth rate, set to reflect the most recent estimates for that country from the World Bank Data (https://data.worldbank.org/indicator/SP.DYN.CBRT.IN), and w reflects the rate of waning of immunity. The population is initiated with Sj,1 = NjRj,1 = 0, and Ij,1 = 0 except for the spatial unit corresponding to the largest population size Nj for each country since this is assumed to be the location of introduction; for this spatial unit, we set Ij,1 = 1.
We made the simplifying assumption that mobility linking locations i and j, denoted as ci,j, scales with the inverse of the cost of travel between sites i and j evaluated according to the friction surface provided in Weiss et al.26. The introduction of an infected individual into location j is then defined by a draw from a Bernouilli distribution following:

$$iota _{j,t} approx {mathrm{Bernouilli}}left( {1 – {mathrm{exp}}left( { – mathop {sum }limits_1^L {c_{i,j}}{I_{i,t}}/{N_i}} right)} right)$$

where L is the total number of administrative 2 units in that country and the rate of introduction is the product of connectivity between the focal location and each other location multiplied by the proportion of population in each other location that is infected.
Some countries show rapid spread between administrative units within the country (for example, a country with parameters that broadly reflect those available for Malawi; Extended Data Fig. 7), while in others (for example, reflecting Madagascar), connectivity may be so low that the outbreak may be over in the administrative unit of the largest size (where it was introduced) before introductions successfully reach other poorly connected administrative units. Where duration of immunity is sufficiently long, the result may be a hump-shaped relationship between the proportion of the population that is infected after five years and the time to the first local extinction of the pathogen (Extended Data Fig. 7, top right). In countries with lower connectivity (for example, resembling Madagascar), local outbreaks can go extinct rapidly before traveling very far; in other countries (for example, resembling Gabon), the pathogen goes extinct rapidly because it travels rapidly and rapidly depletes susceptible individuals everywhere. The U-shaped pattern diminishes as the rate of waning of immunity increases and is replaced by a monotonic negative relationship. With sufficiently rapid waning of immunity, local extinction ceases to occur in the absence of control efforts.
The impact of the pattern of travel between centroids is echoed by the pattern of travel within administrative districts: countries where the pathogen does not reach a large fraction of the administrative 2 units within the country in five years are also those where within-administrative-unit travel is low (Extended Data Fig. 7, right).
These simulations provide a window into qualitative patterns expected for subnational spread of the pandemic virus but there is no clear way of calibrating the absolute rate of travel between regions of relevance for SARS-CoV-2; this is further complicated by the remaining uncertainties around rates of waning of immunity. Thus, the time scales of these simulations should be considered in relative, rather than absolute terms. Variation in lockdown effectiveness, or other changes in mobility for a given country, may also compromise relative comparisons as might large volumes of land border crossings in some settings, which we have not accounted for in this study. Variability in testing and case reporting complicates clarifying this (Extended Data Fig. 7, bottom left and bottom right, respectively) but we have highlighted countries with less connectivity (that is, less synchronous outbreaks expected) relative to the median among SSA countries and with older populations (that is, a greater proportion in higher-risk age groups) (Extended Data Fig. 8).
The University of Oxford’s Blavatnik School of Government generated composite scores of government response, interventions for containment and economic support provided, with each scored from 0 to 100 (Coronavirus Government Response Tracker; https://www.bsg.ox.ac.uk/research/research-projects/coronavirus-government-response-tracker). These data were compared with the day on which ten cases were exceeded in a country according to the Johns Hopkins dashboard data (Johns Hopkins Coronavirus Resource Center; https://coronavirus.jhu.edu/map.html).
While faster waning of immunity will act to increase the rate of spread of the infection, resulting in a higher proportion infected after one year, control efforts will generally act to slow the rate of spread of the infection (Extended Data Fig. 9). Since different countries are likely to have differently effective control efforts (Extended Data Fig. 9), this precludes making country-specific predictions as to the relative impact of control efforts on delay.
Modeling epidemic trajectories in scenarios where transmission rate depends on climate
Climate data sourcing: variation in humidity in SSA
Specific humidity data for selected urban centers comes from the ERA5 using an average climatology (1981–2017)53; we did not consider year-to-year climate variations. Selected cities (n = 56) were chosen to represent the major urban areas in SSA. The largest city in each SSA country was included as well as any additional cities that were among the 25 largest cities or busiest airports in SSA.
Methods for climate-driven modeling of SARS-CoV-2
We used a climate-driven susceptible-infected-recovered-susceptible model to estimate epidemic trajectories (that is, the time of peak incidence) in different cities in 2020, assuming no control measures were in place or a 10 or 20% reduction in R0 beginning 2 weeks after the total reported cases for a country exceeded 10 cases25,63. The model is given by:

$$frac{{mathrm{d}}S}{{mathrm{d}}t} = frac{{N – S – L}}{L} – frac{{beta (t)IS}}{N}$$

$$frac{{mathrm{d}}I}{{mathrm{d}}t} = frac{{beta (t)IS}}{N} – frac{I}{D}$$

where S is the susceptible population, I is the infected population and N is the total population. D is the mean infectious period, set at 5 d following ref. 25.
To investigate the effects on epidemic trajectories of a climate dependency of SARS-CoV-2 on cities with the climate patterns of the selected cities in SSA, we used parameters from the most climate-dependent scenario in ref. 25, based on the endemic betacoronavirus HKU1 in the United States. In this scenario L, the duration of immunity, was 66.25 weeks (that is, >1 year and such that waning immunity did not affect the timing of the epidemic peak). We initially selected a range where R0 declined from R0max = 2.5 to R0min = 1.5 (that is, transmission declined 40% at high humidity) since this exceeds the range observed for influenza and other coronaviruses for which data are available (from the United States). R0max = 2.5 was chosen because 2.5 is often cited as the approximate R0 for SARS-CoV-2. Thus, we initially assumed that the climate dependence of SARS-CoV-2 in SSA would not greatly exceed that of other known coronaviruses from the US context. Then, we explored the effects of different degrees of climate dependency (that is, wider ranges between R0max = 2.5 to R0min = 1.5 and scenarios where R0min approached 1) (Extended Data Fig. 10).
Transmission is governed by β(t), which is related to the basic reproduction number R0 by R0(t) = β(t)D. The basic reproduction number varies based on climate and is related to specific humidity according to the equation:

$$R_0 = {mathrm{exp}}{[a times q(t) + {mathrm{log}}(R_{0{mathrm{max}}} – R_{0{mathrm{min}}})]} + R_{0{mathrm{min}}}$$

where q(t) is specific humidity53 and a is set at −227.5 based on estimated HKU1 parameters25. We assumed the time of introduction for cities to be the date at which the total reported cases for a country exceeded 10 cases.
Sensitivity analysis
Selecting an R0min value of 1, such that epidemic growth stops at high humidities, is likely implausible since simulations indicated no outbreaks would occur in cities such as Antananarivo (countered by the observation that SARS-CoV-2 outbreaks did in fact occur) (Extended Data Fig. 10b; see Supplementary Table 1 for the reported case counts at the country level). Expanding the range between R0min and R0max by increasing R0max resulted in epidemic peaks being reached earlier after outbreak onset but did not increase the difference in timing between cities with different climates (Extended Data Fig. 10c; for example, the difference in timing between peaks in Windhoek and Lomé is similar in 10a and 10c). Finally, we explored scenarios where the R0min was between 1.0 and 1.5. When R0min  > 1.1, epidemic peaks were seen in each SSA city with the difference in timing of the peak growing larger when smaller values of R0min were selected (Extended Data Fig. 10d). However, the difference in timing, even when small values of R0min were selected, was a maximum of 25 weeks and rapidly reduced to only a few weeks when R0min approached 1.5.
Reporting Summary
Further information on research design is available in the Nature Research Reporting Summary linked to this article. More

We selected 10 permanent plots with long-term observations from CERN to include typical forest types of various ages in the East China monsoon forest region, including tropical rainforests, subtropical evergreen coniferous and broad-leaved mixed forests, warm temperate deciduous broad-leaved forests and temperate coniferous and broad-leaved forests, with evident precipitation and temperature gradients from south to north (Fig. 1). The spatial representativeness of the selected 10 sites across the Chinese forest region was evaluated by calculating the Euclidean distance based on various environmental factors. The 10 sites performed well and represented more than 80% of the Chinese forest region (Fig. S1). Of these forests, the Xishuangbanna tropical seasonal rainforest (BNF), Dinghu Mountain subtropical evergreen coniferous and broad-leaved mixed forest (DHF), Ailao Mountain subtropical evergreen broad-leaved forest (ALF), and Changbai Mountain temperate deciduous coniferous and broad-leaved mixed forest (CBF) are mature natural forests; the Shennongjia subtropical evergreen deciduous broad-leaved mixed forest (SNF) and Huitong subtropical evergreen broad-leaved forest (HTF) are natural secondary forests; and the other sites, i.e., the Beijing warm temperate deciduous broad-leaved mixed forest (BJF), Maoxian warm temperate deciduous coniferous mixed forest (MXF), Qianyanzhou subtropical evergreen artificial coniferous mixed forest (QYF), and Heshan subtropical evergreen broad-leaved forest (HSF), are plantations or middle- and young-age forests. All 10 sites are well protected and subject to minimal human activities, thus reflecting the C cycle dynamics under global environmental change, e.g., climate change, increasing CO2 and nitrogen deposition. The detailed characteristics of each plot can be found in their profiles in the CFCCD database.
There are three main steps to create the observation-based basic dataset and assimilated dataset of typical Chinese forests C cycle dynamics:
1.
Observation-based basic data acquisition. An ensemble of daily atmospheric and water data at ten CERN sites were used as forcing datasets for MDF and future scientific analysis; biological and soil data were also collected from CERN and processed by quality control and statistical calculation as benchmark to constrain the model.

2.
Implementation of a multiple data-model fusion framework. The Markov Chain Monte Carlo (MCMC) that integrated the Data Assimilation Linked Ecosystem Carbon (DALEC) model with multiple and dynamic observational data was used to retrieve C-cycle process parameters in a realistic disequilibrium state.

3.
Key process parameters and C function data assimilation. The key parameters of the process-based C cycle model (DALEC) were determined via the model-data fusion method; then the ecosystem C sequestration datasets were simulated by forward running the DALEC model with optimized parameters and then validated based on observational data and other previous studies.

Each step is explained in more detail below.
Observation-based basic data acquisition
Atmospheric and water data
In situ observations of daily air temperature (Ta), photosynthetically active radiation (PAR), relative humidity (RH), precipitation (Precip), and soil moisture (Sw) at the 10 sites from 2005 to 2015 were obtained from the CERN scientific and technological resources service system (http://www.cnern.org.cn/). These atmospheric and water data were mostly observed by an automatic meteorological station at each site. Among them, the PAR was estimated by a LI-COR LI-190SZ Quantum Sensor; Ta and RH were measured by a QMT110 sensor; Sw was estimated by a soil moisture neutron probe or the Time-Domain Reflectometry (TDR) soil moisture probe; the associated saturated soil water capacity (Sc) was measured by the cutting ring method to sample soil in each field campaign and the oven-drying method to measure saturated moisture content after the soil was soaked in water for 48 h; and Precip was artificially observed by CERN staff using an SM1-1 rain gauge. These monitoring data were collected in keeping with CERN’s protocols of observation and quality control29,30.
There were occasional missing data in time-continuous meteorological observations; therefore, the data were processed by standardized gap filling31. Specifically, for Ta, PAR, and RH, which were applied as model driver, we used a linear interpolation method to interpolate continuous missing data with less than three observations; otherwise, we established a regression model using the CERN observations and other observations from adjacent stations of the China Meteorology Administration (756 meteorological stations; http://data.cma.cn/en) to interpolate continuous missing data with more than three observations.
Biological data
Biomass.
At each site, the diameters at breast height (DBHs) and tree heights were measured for each tree in a regular inventory performed at least once every five years. The allometric equations of the DBH and/or tree heights with the biomasses of different plant tissues (i.e., leaves, branches, stems and roots) were developed at each site for various species based on the felled standard trees in the destructive plot. Then, we calculated the biomasses for the ten ecosystems using these allometric equations (FA02 table downloaded from http://www.cnern.org.cn/), which all passed the significance test (0.01 level) and have the R2 most above 0.9 when its estimation compare to observations from standard trees. For some unfelled species under protection, the allometric equations were obtained from Luo et al.32, which were developed based on national inventories and meta analyses from the published literature.
Litterfall.
The aboveground litterfall biomass was measured monthly by ten replicates with 1 m × 1 m baskets during the growing season or once during the nongrowing season. All collected litter was dried at 70 °C for 24 h in the laboratory and then weighed. To avoid the effects of wind on the measurement of litterfall biomass within an individual month, annual litterfall biomass data were finally adopted for each site.
LAI.
The leaf area index (LAI) at each site was measured optically with an LAI-2000 plant canopy analyzer (LI-COR, Lincoln, NE, USA) at least quarterly every year.
Soil data
Soil organic matter (SOM) was measured by the potassium dichromate oxidation titrimetric method. Soil bulk density (SBD) was measured by the cutting ring method in each field campaign at 10 forest sites. Soil particle size (i.e., soil mechanical composition) was measure by the laser particle analyzer. At least three samples were collected from each of the five soil layers (0–10, 10–20, 20–40, 40–60, and 60–100 cm) once every five years.
SOC.
The soil organic C (SOC) content was calculated from SOM, SBD, and volume percentage of gravel with particle size >2 mm at 10 forest sites as follows33:

$$SOC={sum }_{i=1}^{n}0.58times {H}_{i}times {B}_{i}times {O}_{i}times (1-{rm{theta }})times 100$$
(1)

where SOC is the soil organic C density (g C/m2) of all n layers, Hi is the soil thickness (cm), Bi is the soil bulk density (g/cm3), Oi is the SOM content of the ith layer (%), and θ is the volume percentage (%) of gravel with particle size >2 mm. In the absence of soil bulk density or soil organic matter content measurements in some layers, the missing soil measurements corresponding to specific soil depths of theses forest ecosystems were supplemented according to the empirical formulas of the relationships between SOM/soil bulk density and soil depth in different layers, which were developed based on the long-term and across-site CERN soil observations34.
All these raw atmospheric, biological, and soil data mentioned above can be directly download from CERN scientific and technological resources service system (http://www.cnern.org.cn/data/initDRsearch) or obtained after online application via protocol sharing.
Auxiliary flux data
Net ecosystem exchange (NEE).
These data were obtained from ChinaFLUX (http://www.chinaflux.org/), covering CBF, QYF, and BNF. The data were aggregated to the daily time step from half-hourly CO2 flux data measured by the eddy covariance technique and processed with quality control and gap filling procedures35.
Implementation of MDF method
The assimilated data were retrieved from a multiple data-model fusion method (Fig. 2). Specifically, the long-term and dynamic observations of biomass, litterfall, LAI and SOC were used as the model constraint data; Ta, PAR, and RH were used as the meteorological driving data; and the metropolis simulated annealing algorithm, a variation in the MCMC technique36,37, was applied to retrieve the C cycle parameters (e.g., C allocation and C turnover times) against the observations and prior knowledge. Then, we forward-simulated the model to produce the dynamic and time-continuous changes in ecosystem C sequestration function.
Fig. 2

Flowchart of the generation of assimilated datasets in a multiple- and long-term data assimilation framework.

Full size image

Since the dynamic C cycle observations provided an effective solution to constrain the C cycle states without the steady state assumption (SSA), the novelty of our MDF framework involves estimating these C cycle dynamics in better agreement with the actual dynamic disequilibrium state38. Therefore, the uncertainty in allocation and turnover parameters and in C pool states have largely been reduced based on the time-series observations under the non-SSA (NSSA)21,39,40, thereby significantly enhancing the model’s ability to predict the C sequestration function19,41,42.
Carbon cycle process model description
DALEC is a box model of C pools connected via fluxes running at a daily time step and has been applied extensively to the MDF research21,43. Its main structure (i.e., C cycle, C allocation, and turnover process) is generally consistent with state-of-the-art process-based models (Fig. S2; Table S1), with five pools (i.e., foliage (Cf), fine root (Cr), woody (Cw, including branches, stems, and coarse roots), litter (Clit) and SOM (Csom)) for evergreen forests and an additional labile pool (Clab) of stored C that supports leaf flushing for deciduous forests. The C cycle was initiated with the canopy C influx: gross primary productivity (GPP), which was predicted by the Aggregated Canopy Model (ACM)44 (Appendix S1). After GPP is consumed by autotrophic respiration (Ra), the remaining photosynthate (NPP) is allocated to plant tissue pools (Cf, Cr, or Cw). The C exiting from all C reservoirs was based on a first order differential equation with various turnover rates, with temperature and moisture dependency on the turnover from the litter and soil pools. In contrast to the original DALEC model only with temperature scalar fTa, here we added a new function fSw to express soil moisture pressure on litter and soil decomposition processes (Appendix S1). In general, the C pools and fluxes in DALEC were iteratively calculated at a daily time step and determined as a function of the key turnover and allocation parameters. A detailed model description can be found in Williams et al.45 and Fox et al.46.
Multiple data-model fusion at the nonsteady state
In a realistic disequilibrium state, C pools are time-variant, i.e., the C efflux is not equal to the C influx (left(frac{dC}{dt}ne 0right)); thus, the MDF was run via the dynamic and long-term CERN observations to constrain the DALEC model at the non-steady state (Eq. 2). Here, to avoid the uncertainty arising from the spin-up process under SSA, we determined the initial state of the C pools by the initial observations of C stocks or by optimization (i.e., Clab, which cannot be directly observed). Then, the turnover and allocation parameters were retrieved under the disequilibrium state with dynamic environmental forcing. This method avoids the considerable uncertainties when invoking the SSA to estimate the initial state of C pools and the C cycle parameters(e.g., allocation coefficients and turnover rates)39,40,47, which could lead to obvious biases in C sequestration19.

$$left{begin{array}{l}Delta {C}_{i}ne 0\ {C}_{i}left({rm{t}}+1right)={C}_{i}left({rm{t}}right)+{I}_{i}left({rm{t}}right)-{k}_{i}{C}_{i}left({rm{t}}right),{rm{i}}=1,2ldots n\ {C}_{i}left({rm{t}}=0right)={C}_{i}0end{array}right.$$
(2)

where Ci, Ii, and ki represent the size, input and turnover rate of the ith C reservoir, respectively; Ci0 represents the initial state of the ith C reservoir; t represents the specific model-running time step (daily step); and ΔCi represents the ith C pool change between t day and t +1 day when applicable into actual calculation. According to the Bayesian theory, the posterior distributions of the parameters are calculated by maximizing the likelihood function (Eq. 3).

$$L={prod }_{j=1}^{m}{prod }_{i=1}^{{n}_{j}}frac{1}{sqrt{2pi }{sigma }_{j}}{e}^{-{left({x}_{j,i}-{mu }_{j,i}left({boldsymbol{P}}right)right)}^{2}/2{sigma }_{j}^{2}}$$
(3)

where L is the integrated likelihood function; m is the number of multiple data types; n is the number of data points categorized by the jth data type; xj,i is the measured value composed of dynamic C cycle observations; μj,i(P) represents the modeled fluxes and stocks based on parameters under the NSSA (P); and σj is the standard deviation of each data point classified by the jth data type. Moreover, we imposed a sequence of ecological and dynamic constraints on the model parameter inter-relationships and pool dynamics (Appendix S2), which can significantly reduce uncertainty in model parameters and simulations48. The more detailed disequilibrium method can be found in our latest study19.
Key C-cycle process parameters and C sequestration data assimilation
Key process parameter estimation
Here, we mainly focus on how the C input (i.e., the net primary productivity) partitioned into various plant pools (i.e., foliar, wood, and fine roots), i.e., allocation coefficients, which could be directly determined from the optimized parameters (Fig. S3) of the DALEC model after the step 2: MDF method. Another key process parameter, C turnover time, needs further simple statistical calculation based on the model simulations with optimized parameters. Turnover time is commonly estimated by the equation “τ = stock/flux”20,49. Since the C influx is not equal to the C efflux in the realistic dynamic disequilibrium state, the turnover time should be defined as the ratio between the mass of a C pool and its outgoing flux50. Note that with few natural and anthropogenic disturbances in these well-protected CERN sites12,18, the C efflux is approximately equivalent to the Rh from soil and litterfall (mortality) and Ra (growth) from vegetation. Hence, the turnover time for vegetation, soil, and whole ecosystem can be derived as follows:

$${tau }_{veg}=frac{{C}_{live}}{{I}_{live}-Delta {C}_{live}}=frac{{C}_{live}}{litterfall+{R}_{a}}$$
(4)

$${tau }_{soil}=frac{{C}_{dead}}{{I}_{dead}-Delta {C}_{dead}}=frac{{C}_{dead}}{{R}_{h}}$$
(5)

$${tau }_{eco}=frac{{C}_{eco}}{{I}_{eco}-Delta {C}_{eco}}=frac{{C}_{live+}{C}_{dead}}{{R}_{a}+{R}_{h}}$$
(6)

where τveɡ, τsoil, and τeco refer to the biomass, soil and whole-ecosystem turnover times, respectively; Clive, Cdead and Ceco refer to the live biomass C pool size (Cf, Cr, and Cw,), dead organic C pool size (Csoil and Clitter), and the whole-ecosystem C pool size, respectively; Ilive, Idead and Ieco refer to the C input into the live biomass C pool, dead organic C pool, and whole ecosystem C pool, respectively; ΔClive, ΔCdead and ΔCeco refer to the changes in the live biomass C pool, dead organic C pool size, and whole-ecosystem C pool size, respectively; and Ra, Rh and litterfall refer to the autotrophic and heterotrophic respiration, and turnover from all live C pools (i.e., foliage, fine root and woody pools),respectively, as calculated from the DALEC output driven with the estimated parameters during 2005–2015. Since the C reservoirs, fluxes, and turnover times are instantaneous values, we used the averages of the fluxes and reservoirs over multiple years to reflect the average turnover time during a specific period (i.e., 2005–2015).
Time-continuous C sequestration estimation
The optimized parameter values under the NSSA along with the initial observations of the corresponding C pool sizes were used in forward modeling driven by dynamic environmental variables from 2005 to 2015 to obtain the time-continuous soil and vegetation C storage51. The difference between the ecosystem C influx (GPP) and ecosystem respiration (Ra+Rh) is used to examine the ecosystem C sequestration, i.e., net ecosystem productivity (NEP). Similarly, the difference between the ecosystem C influx (GPP) and ecosystem autotrophic respiration (Ra) is used to examine the net primary ecosystem productivity (NPP). More

1.
Alerstam, T., Hedenström, A. & Åkesson, S. Long-distance migration: evolution and determinants. Oikos 103, 247–260 (2003).
Article  Google Scholar 
2.
Gilroy, J. J., Gill, J. A., Butchart, S. H. M., Jones, V. R. & Franco, A. M. A. Migratory diversity predicts population declines in birds. Ecol. Lett. 19, 308–317 (2016).
Article  Google Scholar 

3.
Gill, J. A., Alves, J. A. & Gunnarsson, T. G. Mechanisms driving phenological and range change in migratory species. Philos. Trans. R. Soc. B 374, 20180047 (2019).
Article  Google Scholar 

4.
Méndez, V., Gill, J. A., Alves, J. A., Burton, N. H. K. & Davies, R. G. Consequences of population change for local abundance and site occupancy of wintering waterbirds. Divers. Distrib. 24, 24–35 (2018).
Article  Google Scholar 

5.
Finch, T., Butler, S. J., Franco, A. M. A. & Cresswell, W. Low migratory connectivity is common in long-distance migrant birds. J. Anim. Ecol. 86, 662–673 (2017).
Article  Google Scholar 

6.
Phillips, R. A., Silk, J. R. D., Croxall, J. P., Afanasyev, V. & Bennett, V. J. Summer distribution and migration of nonbreeding albatrosses: Individual consistencies and implications for conservation. Ecology 86, 2386–2396 (2005).
Article  Google Scholar 

7.
Newton, I. The Migration Ecology of Birds (Academic Press, New York, 2008).
Google Scholar 

8.
Grist, H. et al. Site fidelity and individual variation in winter location in partially migratory European shags. PLoS ONE 9, e98562 (2014).
ADS  Article  Google Scholar 

9.
Alves, J. A. et al. Costs, benefits, and fitness consequences of different migratory strategies. Ecology 94, 11–17 (2013).
ADS  Article  Google Scholar 

10.
Evans, D. R. et al. Individual condition, but not fledging phenology, carries over to affect post-fledging survival in a neotropical migratory songbird. Ibis (Lond. 1859) 162, 331–344 (2020).
Article  Google Scholar 

11.
Gill, J. A. et al. Why is timing of bird migration advancing when individuals are not?. Proc. Biol. Sci. 281, 20132161 (2014).
PubMed  PubMed Central  Google Scholar 

12.
Meyburg, B.-U. et al. Orientation of native versus translocated juvenile lesser spotted eagles (Clanga pomarina) on the first autumn migration. J. Exp. Biol. 220, 2765–2776 (2017).
Article  Google Scholar 

13.
Gill, J. A. Does competition really drive population distributions?. Wader Study 126, 166–168 (2019).
Article  Google Scholar 

14.
Berthold, P. Control of Bird Migration (Chapman & Hall, Boca Raton, 1996).
Google Scholar 

15.
Sutherland, W. J. Evidence for flexibility and constraint in migration systems. J. Avian Biol. 29, 441 (1998).
Article  Google Scholar 

16.
Harrison, X. A. et al. Cultural inheritance drives site fidelity and migratory connectivity in a long-distance migrant. Mol. Ecol. 19, 5484–5496 (2010).
Article  Google Scholar 

17.
Piersma, T., Loonstra, A. H. J., Verhoeven, M. A. & Oudman, T. Rethinking classic starling displacement experiments: evidence for innate or for learned migratory directions?. J. Avian Biol. 51, jav.02337 (2020).
Article  Google Scholar 

18.
Þórisson, B. et al. Population size of oystercatchers Haematopus ostralegus wintering in Iceland. Bird Study 65, 274–278 (2018).
Article  Google Scholar 

19.
Méndez, V. et al. Individual variation in migratory behaviour in a sub-Arctic partial migrant shorebird. Behav. Ecol. https://doi.org/10.1093/beheco/araa010 (2020).
Article  Google Scholar 

20.
van de Pol, M. et al. A global assessment of the conservation status of the nominate subspecies of Eurasian oystercatcher Haematopus ostralegus ostralegus. Int. Wader Stud. 20, 47–61 (2014).
Google Scholar 

21.
Méndez, V. et al. Effects of migratory behaviour on breeding phenology and success in a sub-arctic partially migratory shorebird. J. Anim. Ecol. (under review).

22.
Ens, B. J., Safriel, U. N. & Harris, M. P. Divorce in the long-lived and monogamous oystercatcher, Haematopus ostralegus: incompatibility or choosing the better option?. Anim. Behav. 45, 1199–1217 (1993).
Article  Google Scholar 

23.
Ens, B. J., Choudhury, S. & Black, J. M. Mate fidelity and divorce in monogamous birds. In Partnerships in Birds: The Study of Monogamy (ed. Black, J. M.) 344–401 (Oxford University Press, Oxford, 1996).
Google Scholar 

24.
Winger, B. M., Auteri, G. G., Pegan, T. M. & Weeks, B. C. A long winter for the Red Queen: rethinking the evolution of seasonal migration. Biol. Rev. 94, 737–752 (2019).
Article  Google Scholar 

25.
Bulla, M. et al. Unexpected diversity in socially synchronized rhythms of shorebirds. Nature 540, 109–113 (2016).
ADS  CAS  Article  Google Scholar 

26.
Nol, E. Sex roles in the American oystercatcher. Behaviour 95, 232–260 (1985).
Article  Google Scholar 

27.
Reynolds, J. D. & Székely, T. The evolution of parental care in shorebirds: life histories, ecology, and sexual selection. Behav. Ecol. 8, 126–134 (1997).
Article  Google Scholar 

28.
Lazarus, J. The logic of mate desertion. Anim. Behav. 39, 672–684 (1990).
Article  Google Scholar 

29.
Safriel, U. N., Ens, B. J., Kaiser, A. & Goss-Custard, J. D. Rearing to independence. In The Oystercatcher: From Individuals to Populations (ed. Goss-Custard, J. D.) 219–250 (Oxford University Press, Oxford, 1996).
Google Scholar 

30.
Gunnarsson, T. G., Gill, J. A., Sigurbjörnsson, T. & Sutherland, W. J. Pair bonds: arrival synchrony in migratory birds. Nature 431, 646 (2004).
ADS  CAS  Article  Google Scholar 

31.
Gunnarsson, T. G., Gill, J. A., Newton, J., Potts, P. M. & Sutherland, W. J. Seasonal matching of habitat quality and fitness in a migratory bird. Proc. R. Soc. B 272, 2319–2323 (2005).
Article  Google Scholar 

32.
Gunnarsson, T. G. Monitoring wader productivity during autumn passage in Iceland. Wader Study Gr. Bull. 110, 21–29 (2006).
Google Scholar 

33.
Cramp, S. & Simmons, K. E. L. Birds of the Western Palearctic (Oxford University Press, Oxford, 1983).
Google Scholar 

34.
Alves, J. A., Gunnarsson, T. G., Sutherland, W. J., Potts, P. M. & Gill, J. A. Linking warming effects on phenology, demography, and range expansion in a migratory bird population. Ecol. Evol. 9, 2365–2375 (2019).
Article  Google Scholar 

35.
Liebezeit, J. R. et al. Assessing the development of shorebird eggs using the flotation method: species-specific and generalized regression models. Condor 109, 32–47 (2007).
Article  Google Scholar 

36.
Burnham, K. & Anderson, D. Model Selection and Multi-model Inference: A Practical Information-Theoretic Approach (Springer, Berlin, 2002).
Google Scholar 

37.
R Core Team. A Language and Environment for Statistical Computing (R Foundation for Statistical Computing, Vienna, 2020).
Google Scholar  More