Ecology

1.Lanza, B., Pastorelli, C., Laghi, P. & Cimmaruta, R. A review of systematics, taxonomy, genetics, biogeography and natural history of the genus Speleomantes Dubois, 1984 (Amphibia Caudata Plethodontidae). Atti Mus civ stor nat Trieste 52, 5–135 (2006).
Google Scholar 
2.Ficetola, G. F. et al. Differences between microhabitat and broad-scale patterns of niche evolution in terrestrial salamanders. Sci Rep 8, 10575, https://doi.org/10.1038/s41598-018-28796-x (2018).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
3.Lunghi, E., Manenti, R. & Ficetola, G. F. Seasonal variation in microhabitat of salamanders: environmental variation or shift of habitat selection? PeerJ 3, e1122, https://doi.org/10.7717/peerj.1122 (2015).Article 
PubMed 
PubMed Central 

Google Scholar 
4.Ficetola, G. F., Lunghi, E. & Manenti, R. Microhabitat analyses support relationships between niche breadth and range size when spatial autocorrelation is strong. Ecography 43, 1–11, https://doi.org/10.1111/ecog.04798 (2020).Article 

Google Scholar 
5.Culver, D. C. & Pipan, T. The biology of caves and other subterranean habitats 2nd edn (Oxford University Press, 2019).6.Bradley, J. G. & Eason, P. K. Predation risk and microhabitat selection by cave salamanders, Eurycea lucifuga (Rafinesque, 1822). Behaviour 155, 841–859, https://doi.org/10.1163/1568539X-00003505 (2019).Article 

Google Scholar 
7.Salvidio, S., Palumbi, G., Romano, A. & Costa, A. Safe caves and dangerous forests? Predation risk may contribute to salamander colonization of subterranean habitats. The Science of Nature 104, 20, https://doi.org/10.1007/s00114-017-1443-y (2017).ADS 
CAS 
Article 
PubMed 

Google Scholar 
8.Manenti, R., Melotto, A., Guillaume, O., Ficetola, G. F. & Lunghi, E. Switching from mesopredator to apex predator: how do responses vary in amphibians adapted to cave living? Behavioral Ecology and Sociobiology 74, 126, https://doi.org/10.1007/s00265-020-02909-x (2020).Article 

Google Scholar 
9.Lunghi, E. et al. Field-recorded data on the diet of six species of European Hydromantes cave salamanders. Sci Data 5, 180083, https://doi.org/10.1038/sdata.2018.83 (2018).Article 
PubMed 
PubMed Central 

Google Scholar 
10.Lunghi, E. & Bruni, G. Long-term reliability of Visual Implant Elastomers in the Italian cave salamander (Hydromantes italicus). Salamandra 54, 283–286 (2018).
Google Scholar 
11.Mace, G. M. & Lande, R. Assessing extinction threats: towards a reevaluation of IUCN threatened species categories. Conservation Biology 5, 148–157 (1991).Article 

Google Scholar 
12.Huey, R. B. et al. Predicting organismal vulnerability to climate warming: roles of behaviour, physiology and adaptation. Philosophical Transaction of the Royal Society B 367, 1665–1679, https://doi.org/10.1098/rstb.2012.0005 (2012).Article 

Google Scholar 
13.Rondinini, C., Battistoni, A., Peronace, V. & Teofili, C. Lista Rossa IUCN dei Vertebrati Italiani. (Comitato Italiano IUCN e Ministero dell’Ambiente e della Tutela del Territorio e del Mare, 2013).14.European Community. Council Directive 92/43/EEC of 21 May 1992 on the conservation of natural habitats and of wild fauna and flora. Official Journal of the European Union L 206/7, 1–44 (1992).
Google Scholar 
15.Régnier, C. et al. Mass extinction in poorly known taxa. Proc Natl Acad Sci USA 112, 7761–7766, https://doi.org/10.1073/pnas.1502350112 (2015).ADS 
CAS 
Article 
PubMed 

Google Scholar 
16.Stuart, S. N. et al. Status and trends of amphibian declines and extinctions worldwide. Science 306, 1783–1786, https://doi.org/10.1126/science.1103538 (2004).ADS 
CAS 
Article 
PubMed 

Google Scholar 
17.Bellard, C., Bertelsmeier, C., Leadley, P., Thuiller, W. & Courchamp, F. Impacts of climate change on the future of biodiversity. Ecol Lett 15, 365–377, https://doi.org/10.1111/j.1461-0248.2011.01736.x (2012).Article 
PubMed 
PubMed Central 

Google Scholar 
18.Connette, G. M., Crawford, J. A. & Peterman, A. E. Climate change and shrinking salamanders: alternative mechanisms for changes in plethodontid salamander body size. Global Change Biology 21, 2834–2843, https://doi.org/10.1111/gcb.12883 (2015).ADS 
Article 
PubMed 

Google Scholar 
19.Heinrichs, J. A., Bender, D. J. & Schumaker, N. H. Habitat degradation and loss as key drivers of regional population extinction. Ecological Modelling 335, 64–73, https://doi.org/10.1016/j.ecolmodel.2016.05.009 (2016).Article 

Google Scholar 
20.Walters, R. J., Blanckenhorn, W. U. & Berger, D. Forecasting extinction risk of ectotherms under climate warming: an evolutionary perspective. Functional Ecology 26, 1324–1338, https://doi.org/10.1111/j.1365-2435.2012.02045.x (2012).Article 

Google Scholar 
21.Zhang, Z. et al. Future climate change will severely reduce habitat suitability of the Critically Endangered Chinese giant salamander. Freshwater Biology 65, 971–980, https://doi.org/10.1111/fwb.13483 (2020).Article 

Google Scholar 
22.Bland, L. M. Global correlates of extinction risk in freshwater crayfish. Animal Conservation 20, 532–542, https://doi.org/10.1111/acv.12350 (2017).Article 

Google Scholar 
23.Lunghi, E. et al. Photographic database of the European cave salamanders, genus Hydromantes. Sci Data 7, 171, https://doi.org/10.1038/s41597-020-0513-8 (2020).Article 
PubMed 
PubMed Central 

Google Scholar 
24.Mammola, S. et al. Continental data on cave-dwelling spider communities across Europe (Arachnida: Araneae). Biodivers Data J 7, e38492, https://doi.org/10.3897/BDJ.7.e38492 (2019).Article 
PubMed 
PubMed Central 

Google Scholar 
25.MacNeil, R. R. & Brcic, J. Coping with the subterranean environment: a thematic content analysis of the narratives of cave explorers. J Hum Perform Environ 13, Article 6, https://doi.org/10.7771/2327-2937.1089 (2017).Article 

Google Scholar 
26.Zagmajster, M., Culver, D. C., Christman, M. C. & Sket, B. Evaluating the sampling bias in pattern of subterranean species richness: combining approaches. Biodivers Conserv 19, 3035–3048, https://doi.org/10.1007/s10531-010-9873-2 (2010).Article 

Google Scholar 
27.Mammola, S. et al. Collecting eco-evolutionary data in the dark: Impediments to subterranean research and how to overcome them. Ecology and Evolution, https://doi.org/10.1002/ece3.7556 (2021).28.Brown, A. W., Kaiser, K. A. & Allison, D. B. Issues with data and analyses: errors, underlying themes, and potential solutions. Proc Natl Acad Sci USA 115, 2563–2570, https://doi.org/10.1073/pnas.1708279115 (2018).CAS 
Article 
PubMed 

Google Scholar 
29.Crovetto, F., Romano, A. & Salvidio, S. Comparison of two non-lethal methods for dietary studies in terrestrial salamanders. Wildlife Research 39, 266–270, https://doi.org/10.1071/WR11103 (2012).Article 

Google Scholar 
30.Lunghi, E. & Veith, M. Are Visual Implant Alpha tags adequate for individually marking European cave salamanders (genus Hydromantes)? Salamandra 53, 541–544 (2017).
Google Scholar 
31.Swanson, J. E., Bailey, L. L., Muths, E. & Funk, W. C. Factors influencing survival and mark retention in postmetamorphic Boreal chorus frogs. Copeia 2013, 670–675, https://doi.org/10.1643/CH-12-129 (2013).Article 

Google Scholar 
32.Sacchi, R. et al. Photographic identification in reptiles: a matter of scales. Amphibia-Reptilia 31, 489–502 (2010).Article 

Google Scholar 
33.Lunghi, E. et al. On the stability of the dorsal pattern of European cave salamanders (genus Hydromantes). Herpetozoa 32, 249–253, https://doi.org/10.3897/herpetozoa.32.e39030 (2019).Article 

Google Scholar 
34.Lunghi, E. et al. What shapes the trophic niche of European plethodontid salamanders? PLoS ONE 13, e0205672, https://doi.org/10.1371/journal.pone.0205672 (2018).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
35.Lunghi, E. et al. The post hoc measurement as a safe and reliable method to age and size plethodontid salamanders. Ecology and Evolution 10, 11111–11116, https://doi.org/10.1002/ece3.6748 (2020).Article 
PubMed 
PubMed Central 

Google Scholar 
36.Hedrick, B. P. et al. Digitization and the future of natural history collections. BioScience 70, 243–251, https://doi.org/10.1093/biosci/biz163 (2020).Article 

Google Scholar 
37.Nelson, G. & Ellis, S. The history and impact of digitization and digital data mobilization on biodiversity research. Philosophical Transactions of the Royal Society B 374, 20170391, https://doi.org/10.1098/rstb.2017.0391 (2019).Article 

Google Scholar 
38.Lunghi, E. et al. Interspecific and inter-population variation in individual diet specialization: do environmental factors have a role? Ecology 101, e03088, https://doi.org/10.1002/ecy.3088 (2020).Article 
PubMed 

Google Scholar 
39.Salvidio, S., Romano, A., Oneto, F., Ottonello, D. & Michelon, R. Different season, different strategies: feeding ecology of two syntopic forest-dwelling salamanders. Acta Oecol 43, 42–50 (2012).ADS 
Article 

Google Scholar 
40.Rosenblatt, A. E. et al. Factors affecting individual foraging specialization and temporal diet stability across the range of a large “generalist” apex predator. Oecologia 178, 5–16, https://doi.org/10.1007/s00442-014-3201-6 (2015).ADS 
Article 
PubMed 

Google Scholar 
41.Lunghi, E. et al. Same diet, different strategies: variability of individual feeding habits across three populations of Ambrosi’s cave salamander (Hydromantes ambrosii). Diversity 12, 180, https://doi.org/10.3390/d12050180 (2020).Article 

Google Scholar 
42.Lunghi, E., Manenti, R. & Ficetola, G. F. Do cave features affect underground habitat exploitation by non-troglobite species? Acta Oecol 55, 29–35, https://doi.org/10.1016/j.actao.2013.11.003 (2014).ADS 
Article 

Google Scholar 
43.Lunghi, E. et al. Cave morphology, microclimate and abundance of five cave predators from the Monte Albo (Sardinia, Italy). Biodivers Data J 8, e48623, https://doi.org/10.3897/BDJ.8.e48623 (2020).Article 
PubMed 
PubMed Central 

Google Scholar 
44.Carvalho-Rocha, V., Cortês, L. B. & Neckel-Oliveira, S. Interindividual patterns of resource use in three subtropical Atlantic Forest frogs. Austral Ecology 43, 150–158, https://doi.org/10.1111/aec.12552 (2018).Article 

Google Scholar 
45.Lunghi, E. et al. Photos and stomach contents of two mainland Italian Speleomantes salamanders: data from summer 2020. figshare https://doi.org/10.6084/m9.figshare.c.5398368 (2021).46.Martel, A. et al. Batrachochytrium salamandrivorans sp. nov. causes lethal chytridiomycosis in amphibians. Proc Natl Acad Sci USA 110, 15325–15329, https://doi.org/10.1073/pnas.1307356110 (2012).ADS 
Article 

Google Scholar 
47.Treilibs, C. E., Pavey, C. R., Hutchinson, M. N. & Bull, C. M. Photographic identification of individuals of a free-ranging, small terrestrial vertebrate. Ecology and Evolution 6, 800–809, https://doi.org/10.1002/ece3.1883 (2016).Article 
PubMed 
PubMed Central 

Google Scholar 
48.Town, C., Marshall, A. & Sethasathien, N. Manta Matcher: automated photographic identification of manta rays using keypoint features. Ecology and Evolution 3, 1902–1914, https://doi.org/10.1002/ece3.587 (2013).Article 
PubMed 
PubMed Central 

Google Scholar 
49.MacCoun, R. & Perlmutter, S. Hide results to seek the truth. Nature 526, 187–189, https://doi.org/10.1038/526187a (2015).ADS 
CAS 
Article 
PubMed 

Google Scholar 
50.Lunghi, E. et al. Thermal equilibrium and temperature differences among body regions in European plethodontid salamanders. J Therm Biol 60, 79–85, https://doi.org/10.1016/j.jtherbio.2016.06.010 (2016).Article 
PubMed 

Google Scholar 
51.Weller, H. I. & Westneat, M. W. Quantitative color profiling of digital images with earth mover’s distance using the R package colordistance. PeerJ 7, e6398, https://doi.org/10.7717/peerj.6398 (2019).Article 
PubMed 
PubMed Central 

Google Scholar 
52.Adams, D., Collyer, M. & Kaliontzopoulou, A. geomorph. Geometric Morphometric Analyses of 2D/3D Landmark Data. R package version 3.2.1, https://github.com/geomorphR/geomorph (2020).53.Bendik, N. F., Morrison, T. A., Gluesenkamp, A. G., Sanders, M. S. & O’Donnell, L. J. Computer-assisted photo identification outperforms visible implant elastomers in an endangered salamander, Eurycea tonkawae. PLoS ONE 8, e59424, https://doi.org/10.1371/journal.pone.0059424 (2013).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
54.Renet, J., Leprêtre, L., Champagnon, J. & Lambret, P. Monitoring amphibian species with complex chromatophore patterns: a non-invasive approach with an evaluation of software effectiveness and reliability. Herpetological Journal 29, 13–22, https://doi.org/10.33256/hj29.1.1322 (2019).Article 

Google Scholar 
55.Allen-Blevins, C. R., You, X., Hinde, K. & Sela, D. A. Handling stress may confound murine gut microbiota studies. PeerJ 5, e2876, https://doi.org/10.7717/peerj.2876 (2017).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
56.Samimi, A. S., Tajik, J., Jarakani, S. & Shojaeepour, S. Evaluation of a five-minute resting period following handling stress on electrocardiogram parameters and cardiac rhythm in sheep. Veterinary Science Development 6, 6481, https://doi.org/10.4081/vsd.2016.6481 (2016).Article 

Google Scholar 
57.Martel, A. et al. Recent introduction of a chytrid fungus endangers Western Palearctic salamanders. Science 346, 630, https://doi.org/10.1126/science.1258268 (2014).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
58.Lunghi, E., Corti, C., Manenti, R. & Ficetola, G. F. Consider species specialism when publishing datasets. Nat Ecol Evol 3, 319, https://doi.org/10.1038/s41559-019-0803-8 (2019).Article 
PubMed 

Google Scholar  More

1.Agarwal, V. et al. Enzymatic halogenation and dehalogenation reactions: pervasive and mechanistically diverse. Chem. Rev. 117, 5619–5674 (2017).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
2.Atashgahi, S., Haggblom, M. M. & Smidt, H. Organohalide respiration in pristine environments: implications for the natural halogen cycle. Environ. Microbiol. 20, 934–948 (2018).CAS 
PubMed 
Article 

Google Scholar 
3.Gribble G. W. Naturally Occurring Organohalogen Compounds—A Comprehensive Update: (Wien/New York: Springer, 2010).4.Stringer, R. & Johnston, P. Chlorine and the environment: an overview of the chlorine industry. Environ. Sci. Pollut. Res. 8, 146–159 (2001).Article 

Google Scholar 
5.Falandysz, J., Rose, M. & Fernandes, A. R. Mixed poly-brominated/chlorinated biphenyls (pxbs): widespread food and environmental contaminants. Environ. Int. 44, 118–127 (2012).CAS 
PubMed 
Article 

Google Scholar 
6.Zhang, Z. et al. Halogenated organic pollutants in sediments and organisms from mangrove wetlands of the jiulong river estuary, south china. Environ. Res. 171, 145–152 (2019).CAS 
PubMed 
Article 

Google Scholar 
7.Kunze, C. et al. Cobamide-mediated enzymatic reductive dehalogenation via long-range electron transfer. Nat. Commun. 8, 15858 (2017).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
8.Wang, S. et al. Electron transport chains in organohalide-respiring bacteria and bioremediation implications. Biotechnol. Adv. 36, 1194–1206 (2018).PubMed 
Article 
CAS 

Google Scholar 
9.Atashgahi S., Lu Y., Smidt H. Overview of known organohalide-respiring bacteria—phylogenetic diversity and environmental distribution. In Adrian L., Löffler F. E., editors. Organohalide-Respiring Bacteria. (Springer, Berlin, 2016).10.Fincker, M. & Spormann, A. M. Biochemistry of catabolic reductive dehalogenation. Annu. Rev. Biochem. 86, 357–386 (2017).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
11.Rouzeau-Szynalski, K., Maillard, J. & Holliger, C. Frequent concomitant presence of desulfitobacterium spp. and “dehalococcoides” spp. in chloroethene-dechlorinating microbial communities. Appl. Microbiol. Biotechnol. 90, 361–368 (2011).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
12.Wang, S. & He, J. Dechlorination of commercial pcbs and other multiple halogenated compounds by a sediment-free culture containing dehalococcoides and dehalobacter. Environ. Sci. Technol. 47, 10526–10534 (2013).CAS 
PubMed 
PubMed Central 

Google Scholar 
13.Wang, S. et al. Genomic characterization of three unique dehalococcoides that respire on persistent polychlorinated biphenyls. Proc. Natl. Acad. Sci. U.S.A. 111, 12103–12108 (2014).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
14.Bachmann, H. et al. Availability of public goods shapes the evolution of competing metabolic strategies. Proc. Natl. Acad. Sci. U.S.A. 110, 14302–14307 (2013).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
15.Yan, J., Ritalahti, K. M., Wagner, D. D. & Loffler, F. E. Unexpected specificity of interspecies cobamide transfer from geobacter spp. To organohalide-respiring dehalococcoides mccartyi strains. Appl. Environ. Microbiol. 78, 6630–6636 (2012).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
16.Lawson, C. E. et al. Common principles and best practices for engineering microbiomes. Nat. Rev. Microbiol. 17, 725–741 (2019).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
17.Wang, S., Chen, C., Zhao, S. & He, J. Microbial synergistic interactions for reductive dechlorination of polychlorinated biphenyls. Sci. Total. Environ. 666, 368–376 (2019).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
18.Men, Y. et al. Sustainable syntrophic growth of dehalococcoides ethenogenes strain 195 with desulfovibrio vulgaris hildenborough and methanobacterium congolense: Global transcriptomic and proteomic analyses. ISME J. 6, 410–421 (2012).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
19.Lu, Q. et al. Dehalococcoides as a potential biomarker evidence for uncharacterized organohalides in environmental samples. Front. Microbiol. 8, 1677 (2017).PubMed 
PubMed Central 
Article 

Google Scholar 
20.Xu, G., Lu, Q., Yu, L. & Wang, S. Tetrachloroethene primes reductive dechlorination of polychlorinated biphenyls in a river sediment microcosm. Water Res. 152, 87–95 (2019).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
21.Royer, D. L., Osborne, C. P. & Beerling, D. J. Carbon loss by deciduous trees in a co2-rich ancient polar environment. Nature. 424, 60–62 (2003).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
22.Wang, S. & He, J. Phylogenetically distinct bacteria involve extensive dechlorination of aroclor 1260 in sediment-free cultures. PLoS One 8, e59178 (2013).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
23.Wang, S. et al. Development of an alkaline/acid pre-treatment and anaerobic digestion (apad) process for methane generation from waste activated sludge. Sci. Total Environ. 708, 134564 (2020).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
24.Lu, Q. et al. Inhibitory effects of metal ions on reductive dechlorination of polychlorinated biphenyls and perchloroethene in distinct organohalide-respiring bacteria. Environ. Int. 135, 105373 (2020).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
25.Kumar, S., Stecher, G. & Tamura, K. MEGA7: Molecular evolutionary genetics analysis version 7.0 for bigger datasets. Mol. Biol. Evolut. 33, 1870–1874 (2016).CAS 
Article 

Google Scholar 
26.Cummings, D. E. et al. Diversity of geobacteraceae species inhabiting metal-polluted freshwater lake sediments ascertained by 16S rDNA analyses. Microb. Ecol. 46, 257–269 (2003).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
27.Holmes, V. F., He, J., Lee, P. K. & Alvarez-Cohen, L. Discrimination of multiple dehalococcoides strains in a trichloroethene enrichment by quantification of their reductive dehalogenase genes. Appl. Environ. Microbiol. 72, 5877–5883 (2006).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
28.Joshi N. A., Fass J. N. Sickle: A sliding-window, adaptive, quality-based trimming tool for FastQ files (Version 1.33) [Software]. https://github.com/najoshi/sickle. (2011).29.Bankevich, A. et al. Spades: a new genome assembly algorithm and its applications to single-cell sequencing. J. Comput. Biol. 19, 455–477 (2012).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
30.Langmead, B. & Salzberg, S. L. Fast gapped-read alignment with bowtie 2. Nat. Methods. 9, 357–359 (2012).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
31.Li, H. et al. The sequence alignment/map format and samtools. Bioinformatics. 25, 2078–2079 (2009).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
32.Chen, L. X., Anantharaman, K., Shaiber, A., Eren, A. M. & Banfield, J. F. Accurate and complete genomes from metagenomes. Genome Res. 30, 315–333 (2020).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
33.Parks, D. H., Imelfort, M., Skennerton, C. T., Hugenholtz, P. & Tyson, G. W. Checkm: assessing the quality of microbial genomes recovered from isolates, single cells, and metagenomes. Genome Res. 25, 1043–1055 (2015).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
34.Seemann, T. Prokka: rapid prokaryotic genome annotation. Bioinformatics. 30, 2068–2069 (2014).CAS 
PubMed 
Article 

Google Scholar 
35.Moriya, Y., Itoh, M., Okuda, S., Yoshizawa, A. C. & Kanehisa, M. Kaas: an automatic genome annotation and pathway reconstruction server. Nucleic Acids Res. 35, W182–W185 (2007).PubMed 
PubMed Central 
Article 

Google Scholar 
36.Loffler, F. E. et al. Dehalococcoides mccartyi gen. Nov., sp. Nov., obligately organohalide-respiring anaerobic bacteria relevant to halogen cycling and bioremediation, belong to a novel bacterial class, dehalococcoidia classis nov., order dehalococcoidales ord. Nov. And family dehalococcoidaceae fam. Nov., within the phylum chloroflexi. Int. J. Syst. Evol. Microbiol. 63, 625–635 (2013).CAS 
PubMed 
Article 

Google Scholar 
37.Garcia C.A. Subsurface occurrence and potential source areas of chlorinated Ethenes Identified using concentrations and concentration ratios, Air Force Plant 4 and Naval Air Station-Joint Reserve Base Carswell Field, Fort Worth, Texas. U.S. Geological Survey, (2006).38.Nichols H. Use of electrical resistive heating for the remediation of CVOC and petroleum impacts in soil and groundwater, NEWMOA Conference. New York City (2012).39.Maymó-Gatell, X., Chien, Y., Gossett, J. M. & Zinder, S. H. Isolation of a bacterium that reductively dechlorinates tetrachloroethene to ethene. Science. 276, 1568–1571 (1997).PubMed 
Article 
PubMed Central 

Google Scholar 
40.Sung, Y. et al. Geobacter lovleyi sp. Nov. Strain sz, a novel metal-reducing and tetrachloroethene-dechlorinating bacterium. Appl. Environ. Microbiol. 72, 2775–2782 (2006).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
41.Amos, B. K., Suchomel, E. J., Pennell, K. D. & Löffler, F. E. Spatial and temporal distributions of Geobacter lovleyi and Dehalococcoides spp. during bioenhanced PCE-NAPL dissolution. Environ. Sci. Technol. 43, 1977–1985 (2009).CAS 
PubMed 
Article 

Google Scholar 
42.Lai, Y. & Becker, J. G. Compounded effects of chlorinated ethene inhibition on ecological interactions and population abundance in a Dehalococcoides–Dehalobacter coculture. Environ. Sci. Technol. 47, 1518–1525 (2013).CAS 
PubMed 

Google Scholar 
43.Lovley, D. R. et al. Geobacter: the microbe electric’s physiology, ecology, and practical applications. Adv. Microb. Physiol. 59, 1–100 (2011).CAS 
PubMed 
Article 

Google Scholar 
44.Reguera, G. & Kashefi, K. The electrifying physiology of geobacter bacteria, 30 years on. Adv. Microb. Physiol. 74, 1–96 (2019).PubMed 
Article 

Google Scholar 
45.Kruse, S., Goris, T., Westermann, M., Adrian, L. & Diekert, G. Hydrogen production by sulfurospirillum species enables syntrophic interactions of epsilonproteobacteria. Nat. Commun. 9, 4872 (2018).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
46.Sanford, R. A., Cole, J. R. & Tiedje, J. M. Characterization and description of anaeromyxobacter dehalogenans gen. Nov., sp. Nov., an aryl-halorespiring facultative anaerobic myxobacterium. Appl. Environ. Microbiol. 68, 893–900 (2002).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
47.Sung, Y. et al. Characterization of two tetrachloroethene-reducing, acetate-oxidizing anaerobic bacteria and their description as desulfuromonas michiganensis sp. Nov. Appl. Environ. Microbiol. 69, 2964–2974 (2003).CAS 
PubMed 
Article 

Google Scholar 
48.Peng, P. et al. Organohalide-respiring desulfoluna species isolated from marine environments. ISME J. 14, 815–827 (2020).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
49.Goris, T. et al. Insights into organohalide respiration and the versatile catabolism of sulfurospirillum multivorans gained from comparative genomics and physiological studies. Environ. Microbiol. 16, 3562–3580 (2014).CAS 
PubMed 
Article 

Google Scholar 
50.Goris, T. et al. Proteomics of the organohalide-respiring epsilonproteobacterium sulfurospirillum multivorans adapted to tetrachloroethene and other energy substrates. Sci. Rep. 5, 13794 (2015).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
51.Cross, K. L. et al. Targeted isolation and cultivation of uncultivated bacteria by reverse genomics. Nat. Biotechnol. 37, 1314–1321 (2019).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
52.Lewis, W. H. & Ettema, T. J. G. Culturing the uncultured. Nat. Biotechnol. 37, 1278–1279 (2019).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
53.Koch, H., van Kessel, M. & Lucker, S. Complete nitrification: insights into the ecophysiology of comammox nitrospira. Appl. Microbiol. Biotechnol. 103, 177–189 (2019).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
54.Marcy, Y. et al. Dissecting biological “dark matter” with single-cell genetic analysis of rare and uncultivated TM7 microbes from the human mouth. Proc. Natl. Acad. Sci. U.S.A. 104, 11889–11894 (2007).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
55.Pfeiffer, T., Schuster, S. & Bonhoeffer, S. Cooperation and competition in the evolution of ATP-producing pathways. Science. 292, 504–507 (2001).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
56.Roller, B. R. & Schmidt, T. M. The physiology and ecological implications of efficient growth. ISME J. 9, 1481–1487 (2015).PubMed 
PubMed Central 
Article 

Google Scholar 
57.Lawson, C. E. & Lucker, S. Complete ammonia oxidation: an important control on nitrification in engineered ecosystems? Curr. Opin. Biotechnol. 50, 158–165 (2018).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
58.Kits, K. D. et al. Kinetic analysis of a complete nitrifier reveals an oligotrophic lifestyle. Nature. 549, 269–272 (2017).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar  More

Database constructionCollection of metagenomic dataBetween March 2012 and May 2016, seawater samples were collected from five different locations around Sendai Bay, Ofunato Bay, and A-line. Samples (N = 142) were collected from the surface and the subsurface chlorophyll-a maximum (SCM) layers (Fig. 1). DNA was extracted from microorganisms trapped in filters (pore sizes, 0.2, 0.8, 5, 20, and 100 μm). Shotgun metagenomic sequence data (N = 454) were acquired (Supplementary Information 1), comprising 3.57 × 109 reads and 3.56 × 1011 bases in total (Supplementary Information 2). 16S rRNA gene sequences in the 0.2-μm fraction were subjected to PCR using universal primers to obtain 111 amplicon metagenomic sequences (6.92 × 106 reads, 1.69 × 109 bases) (Supplementary Information 2).Figure 1Location, changes in water temperature, and changes in chlorophyll-a (Chl-a) concentrations at the sampling points. (a) Sampling points along the Pacific coast of northeastern Japan (Sendai Bay, Ofunato Bay, and A-line). (b) Sampling points C5 and C12 in Sendai Bay. The map was generated using Ocean Data View (https://odv.awi.de) with data imported from the NOAA server (accessed on 22 February 2021). (c) Changes in water temperature and chlorophyll concentrations at Sendai Bay and A-line sampling points. Red circles indicate the depth of the sampled water. X-axis: dates from 2012 to 2014. Y-axis: water depth from the surface.Full size imageTime-series analysis of microbial species compositionWe performed microbial taxonomic assignments through analyses of amplicon and shotgun metagenomic sequence data using the SILVA and NCBI NT databases. For C5 and C12 samples from Sendai Bay and A4 and A21 samples from A-line, sufficient time-series data were available to plot changes in microbiota over time (Fig. 2). By clicking the displayed taxon at the website of Ocean Monitoring Database, the microbiota composition of lower taxa is revealed. For example, Fig. 2 shows that the cyanobacteria community increased in abundance during the summer.Figure 2Time-series analysis of microbial communities along the Pacific coast of northeastern Japan. Each sampling point shows the number of ribosomal sequences normalized to 1000 (excluding no hits). Clicking on the graph at the website of Ocean Monitoring Database exhibits the next taxonomic levels. This figure shows an example of the change in Cyanobacteria communities over time from April 2012 to May 2014. SUF: surface layer (1 m), SCM: the subsurface chlorophyll-a maximum layer.Full size imageWe acquired substantial 16S rRNA gene amplicon sequencing data at the C5 fixed point at Sendai Bay between 2012 and 2014. We generated a 3D graph to simultaneously display the date, water depth, and species composition at this site (Fig. 3). By clicking on the taxon shown in the graph at the website of Ocean Monitoring Database, the composition of the microbiota within a lower taxon is displayed. The 3D display is suitable for presenting a bird’s-eye view of the metagenomic data, which is extremely useful for visualizing and understanding the relationships among microbial communities among sampling points. This innovative function was incorporated into the metagenomic database. Figure 3 shows a contour map of chlorophyll concentrations on the x-axis and the proportion of microbial communities on the z-axis. The proportion of flavobacteria may increase following an increase in chlorophyll concentration during a spring bloom. For example, Buchan et al. reported that the proportion of flavobacteria increase late in a spring bloom16; our results show similar patterns (Fig. 3).Figure 3Three-dimensional (3D) display of microbial communities. 3D display of bacterial communities identified using 16S rRNA gene amplicon analysis of the Sendai Bay C5 samples from 2012 to 2014. The x-axis indicates the date, the y-axis indicates the water depth, and the z-axis indicates the percentage abundance of bacterial genera. The contour plot on the xy plane indicates the chlorophyll concentration. The composition of Flavobacteriaceae is shown as an example.Full size imageDigital DNA chip (DDC) databaseA DDC analysis (DDCA) system is useful for visualizing the characteristics of shotgun metagenomic data as a microarray17 of, for example, filter size, water sampling point, water sampling time, temperature, salinity, and nitrate and phosphate concentrations (Fig. 4a). By mapping sequence data against the probe sets described above, which are associated with environmental factors, we predicted that sequence data would be more enriched and inclusive of environmental information. Figure 4b displays the DDCA shotgun metagenomic data of the 0.2–0.8-μm fraction of the C5 sample collected from Sendai Bay on December 1, 2013. The sample contains a bacterial-fraction DNA marker with filter sizes of 0.2–0.8 µm and a specific DNA marker for December in Sendai Bay. Even if there is only NGS data and no environmental information, just by looking at the digital DNA chip, we can assume this sample is extracted from 0.2 to 0.8 µm fraction and is from Sendai Bay (Fig. 4b).Figure 4Visualization of metagenomics data using digital DNA chips. (a) Overview of in silico probes associated with the environmental factors on a digital DNA chip (See Supplementary Information 8 for details). (b) Digital DNA chip of shotgun metagenomics data of a 0.2–0.8-μm fraction of December 1, 2013, Sendai Bay C5. There are 748 probes, and spots that are positive for digital hybridization are shown in red. Negative spots are black. The hybridization positive probes are an indicator of environmental information of the sequence data.Full size imageDevelopment of a shotgun metagenomic databaseWe assembled the shotgun metagenomic sequence data using Megahit version 1.0.218. There were 57.95 M contigs, with an N50 of 995 bp, a maximum length of 307,212 bp, and a total of 12.39 Gbp (Supplementary Information 3). We calculated the abundance pattern of each contig. Those contigs whose appearance pattern matched with a Pearson correlation coefficient of ≥ 0.95 were clustered into a MAG. We next added the annotation of assembled contigs to the results of the BLAST search of the NCBI NT database and using classification by clustering with Pfam (CCP). This novel annotation method is described below. We developed the database showing the abundance pattern of homologous contigs against a queried sequence by BLAST for each sampling point and filter size (Fig. 5). For example, we found novel PolD families using this database.Figure 5Search for homologous contigs to a query sequence and display of temporal variation patterns. Using nucleotide and amino acid sequences as queries, contigs homologous to the query sequence are identified using BLAST, and the temporal variation patterns and taxonomy information of the hit contigs are displayed.Full size imageDevelopment of a new annotation method for metagenome contigsWe annotated the assembled contigs using BLAST to analyze the NCBI NT database. However, we were unable to annotate  > 50% of the contigs (Fig. 6). Therefore, we developed a novel method, i.e., CCP, to annotate contigs according to their species names. Analysis using a single contig generally does not provide sufficient information for assigning an annotation. However, CCP assigns the appropriate annotation to the sequence because it aggregates the Pfam information of all contigs in a MAG. CCP annotates a MAG by comparing the similarity to the reference genome. Comparing the nucleotide sequences of a MAG directly using blastn to analyze the NCBI NT database shows relatively low homology to de novo virus sequences.Figure 6Comparison between BLAST and CCP annotation results at the super-kingdom level. Comparison of classification results using BLAST to annotate contigs and classification by clustering with Pfam (CCP); the percentage of unknowns was 57% for BLAST and 8% for CCP.Full size imageHowever, a Pfam domain search using HMMER, which employs a different principle19 than BLAST, often detects more informative sequences, even those of viruses. For example, phylogenetic trees constructed according to the type and number of Pfam domains of individual bacterial genomes and those of higher eukaryotes such as humans closely approximate those generated using existing phylogenetic trees20. Thus, genomes with a similar Pfam domain may represent phylogenetically closely related species. We, therefore, searched the Pfam domains for reference genomes of viruses, bacteria, archaea, and eukaryotes included in RefSeq (as of August 31, 2015). We next calculated the number of domains for each species and constructed a CCP database. The types and numbers of Pfam domains contained in the contig obtained from the metagenome were summarized in MAG units. We compared the results using the CCP database and annotated the known genomes with the closest correlation coefficient (Fig. 7).Figure 7Overview of CCP. Flowchart of the search of the Pfam domain against known genomes of viruses, bacteria, archaea, and eukaryotes included in RefSeq to create a Pfam hit database. The Pfam domains were searched in metagenome-assembled genome (MAG) units and the known genomes whose type and number of Pfam domains are closest to the MAG.Full size imageBy annotating the top 10,000 contigs with the highest abundance in our database using CCP,  > 90% of the contigs were explained (Fig. 6). In contrast, the BLAST species search (the existing method) returned  15%, indicating that CCP is a robust method, particularly when applied to the identification of virus annotation. We next compared the agreement between CCP and BLAST annotations using contigs annotated using both CCP and BLAST (Supplementary Information 4). The virus-level agreement between CCP and BLAST was 89.6%, and the kingdom-level agreement was 76.9%. It was difficult to determine whether the contig represented a virus using BLAST; however, CCP showed higher accuracy.Shotgun metagenomic analysisPeriodicity of metagenomic dataFor the bacterial fraction (0.2–0.8 µm) of the shotgun metagenomic data from Sendai Bay, we generated a multidimensional scaling (MDS) plot according to the pattern of the abundance of the assembled contigs (Fig. 8). The MDS plot shows similarities among the samples collected during the same month during different years. Furthermore, the plot reveals that the shotgun metagenomic data exhibit an annual seasonal cycle like the 18S rRNA amplicon data10. However, we did not observe the same annual cycle among all contigs. We, therefore, extracted contigs included in the top 20 highly abundant MAGs from the bacterial fractions of Sendai Bay samples collected from March 2012 to April 2014 and plotted the fluctuation patterns. Only one such contig showed a complete 2-year cycle (Fig. 9a), and four contigs showed an incomplete 2-year cycle with peaks in March 2012 and 2013 but not in 2014 (Fig. 9b). Furthermore, 13 MAGs showed a transient pattern (Fig. 9c), and two MAGs showed peaks with irregular patterns (Fig. 9d). These results suggest that marine microbial communities generally undergo an annual cycle.Figure 8Multidimensional scaling (MDS) plot as a function of the abundance of contigs. MDS plots of bacterial fractions (0.2–0.8 µm) of shotgun metagenomic data from 2012 to 2015 acquired from Sendai Bay according to the pattern of abundance of assembled contigs.Full size imageFigure 9Variation patterns of contigs in the top 20 most abundant metagenome-assembled genomes (MAGs).The top 20 MAGs in the bacterial fractions of Sendai Bay C5 and C12 from March 13, 2012, to April 2, 2014, were classified as follows: (a) complete 1-year cycle for 2.5 years, (b) Incomplete 1-year cycle for 2.5 years, (c) transient peaks, and (d) irregular peaks. A peak within 1 month of ≥ 25% relative to the previous year’s peak was considered cyclical.Full size imageIdentification of repeat sequences in the metagenomesDuring the collection and analysis of the DNA sequencing data, we identified a number of repeat sequences in the metagenomes as follows: (TAG)n, (TGA)n, (GAA)n, and (ACA)n microsatellites. We then determined the frequencies and highest numbers of (TAG)n repeats as a function of filter size. We found that the (TAG)n repeats included up to 7.5% of the 5–20-μm fraction (Supplementary Information 5a,b). To investigate whether this was a characteristic feature of the northeastern coastal region of Japan, we analyzed the shotgun metagenomic sequence data of Tara Oceans21,22. As shown in Supplementary Information 5c,d, Tara Oceans data contained up to 1.9% of TAG repeats.To determine whether these (TAG)n repeats represented artifacts of the NGS method, we performed Southern blot and dot-blot hybridization analyses of the DNA samples extracted from seawater (Fig. 10). The dot-blot hybridization experiment analyzed 13 different samples with various content rates (Fig. 10a). We detected signals from the eight samples containing the (TAG)n that were repeated in  > 0.9% of the labeled d54-mer with the (TAG)18 repeat. In contrast, six samples with a low content ( > 0.2%) were negative (Fig. 10b). To determine whether these repeated sequences originated from a single locus or multiple loci, we performed Southern blot analysis (Fig. 10c, Supplementary Information 6) using two samples with high contents of (TAG)n repeats. A (TAG)n representing a single locus is detectable as a discrete band versus the diffuse bands exhibited by two samples with a high content of (TAG)n repeats. The data (Fig. 10c) suggest that the (TAG)n repeats were derived from multiple loci of distinct genomes. Samples with low numbers of (TAG)n repeats were negative.Figure 10Detection of TAG repeats using Southern blot and dot-blot hybridization analyses. (a) Contents of the TAG repeats of the samples according to next-generation sequencing analysis. (b) Dot-blot analyses. The sample numbers and their amounts, (right side) correspond to the signals of each dot in the left panels. The intact pTV119N plasmid without an insert indicates pTV(0). The calculated contents of TAG repeats (%) are indicated in parentheses. (c) Southern blot analysis of EcoRI-digested samples subjected to 0.8% agarose gel electrophoresis. The plasmid pTV (TAG) (0.35 ng and 1 ng) served as a positive control. E. coli genomic DNA served as a negative control. The calculated contents of TAG repeats (%) are indicated on the bottom of each graph. The length (nt) of the TAG-repeated fragment excised from pTV (TAG) is shown on the right.Full size imageThese results reveal for the first time that such repeat sequences are abundant in the genomes of marine microorganisms. However, their species of origin and functional roles were not identified here. The repeat sequences found in Escherichia coli23, subsequently called CRISPR, led to fundamental discoveries that are essential in the field of genetic engineering24. Thus, understanding the biological significance of trinucleotide repeats in marine microorganisms is of particular importance and may reveal a new research frontier. More

Calculation of the indicatorsFigure 1 shows the distribution of the urban greening benefit indicators computed at European scale under the current scenario, while Fig. 2 shows the cumulative distribution of impervious urban areas by increasing value of each indicator, under the current and future scenarios. It should be stressed that, while the indicators of Eqs. (1–4) are computed for every grid cell, the curves of Fig. 2 reflect also the spatial distribution of impervious urban surfaces, and hence they give more prominence to the values of the indicators in the most densely urbanized areas of the continent.Figure 1Maps of benefits per m2 across Europe for ΔTs (a), ΔT (b), RR/P (c) and CB (d), in the present scenario.Full size imageFigure 2Cumulative curves of urban surfaces versus the indicator ΔTs (a), ΔT (b), RR/P (c) and CB (d). The black line represents present conditions, while lines in color stand each for one climatic scenario. The y-axis is the cumulative surface area of the present European urban areas.Full size imageThe reduction of surface temperature ΔTs (Fig. 1a) is highest in the warmer and not excessively dry climates of Central and Southern Europe, reflecting the patterns of actual evapotranspiration. Most European urban areas would achieve temperature reductions of about 3–3.5 °C (Fig. 2a), slightly increasing with the severity of climate heating under the various scenarios, causing a reduction of sensible heat to the atmosphere, a driver of urban heat island effects, between 20 and 40% (see Appendix 1, Supplementary Material for further details). The highest temperature reduction at the roof surface, ΔT, is mostly perceived in the South of Europe (Fig. 1b), consistent with the pattern of potential evapotranspiration, similarly to the production of dry biomass CB (Fig. 1d). The reduction of temperature at the roof is predicted between 15 and 17 °C for most of Europe under the current scenario, and may increase of about 2 °C under the most severe climate scenario (Fig. 2b). Runoff reduction is significantly higher in areas with moderate precipitation, particularly in the plains, compared to rainier areas such as the Atlantic edge of the continent and high mountain ranges (Fig. 1c).The maximum storage volume, Vmax , calculated by Eq. 6, would allow to reuse 92% of the annual runoff, while Vmin and Vavg would allow to store 77% and 86% of the runoff, respectively, as resulting from a daily balance of the storage volume calculated over the 14 year time series. As the storage volume normalized to the annual runoff Rc is 0.24, 0.36 and 0.51 for Vmin, Vavg and Vmax, respectively (Figure 4b), choosing a storage volume equal to Vmin appears to be the most cost-effective solution. Vmin is mapped as shown in Fig. 3a for the case of constant demand, under the current scenario, while in Fig. 3b the volumes are plotted versus the cumulated areas.Figure 3Storage volume Vmin required to store the runoff in the case of constant demand.Full size imageFigure 4Runoff that could be harvested, and normalized storage volume Vmin versus the annual average runoff (Rc) for the case of constant withdrawal, calculated throughout the 14 year time series.Full size imagePhysical and environmental implicationsThese potential effects of green surfaces at European scale correspond to potential benefits. The total benefits extrapolated for the EU are summarized in Table 1. Results are referred to the impervious surfaces corresponding to building roofs, that are assumed to amount to a total of 26,450 km2 as per Bódis et al.40 . This represents 35% of the European impervious surface. Although it is highly unlikely that the majority of the roofs may support a uniform soil cover of 30 cm, they could still bear patches of that thickness over a part of their surface. Moreover, additional surfaces such as sealed ground could be greened. Overall, having in mind these considerations, we pragmatically regard this 35% of impervious urban areas as a maximum extent that could be greened in Europe. All benefits calculated below would obviously scale proportionally for any reduction of the percentage of area subjected to greening. The quantification of Table 1 is explained below.Table 1 Climatic descriptors and quantification of annual benefits at the European scale in the present and future climatic scenarios, assuming to green all roof surfaces, or 35% of the European impervious surfaces.Full size tableThe reduction of land surface temperatures, ΔTs, reduces the thermal irradiation and convective heat flux from urban surfaces (see Appendix 1 of Supplementary Material), which are the drivers of the heat island effect44. As a first order approximation, the reduction of air temperature at 2 m from the surfaces can be expected to be about a half of ΔTs45 as an average value in summer. The reduction of air temperature would generate economic benefits, like the life cycle extension of electronic material and cars, benefits in the health and transport sectors, reduction of social stress and morbidity, and reduction of damages to trees and animals46,47,48.The reduction of the surface temperature ΔT potentially reduces the cooling demand in summer (Eq. 5) by 92 TWh year−1. This energy saving corresponds to 29.9 Mtons of CO2 for the present scenario, considering emissions of 0.325 kg CO2 equivalent kWh−1 for European electricity39. Our estimate is arguably an upper limit of cooling energy savings. In many cases, underroof spaces of buildings are not cooled and effectively work already as an insulation, hence the reduction in the heat transferred from the roofs to underlying inhabited spaces may be lower than we estimate.The yearly produced biomass CB is a benefit in itself whenever the biomass may be used (e.g. crops from urban agriculture). However, more importantly, it may be appraised in terms of carbon and carbon dioxide sequestration. The carbon dioxide sequestered from the atmosphere through biomass growth is 25.9 Mtons year−1 in the present scenario. This must be summed to the reduction of carbon emission following the expected decrease in cooling energy use for a total of 55.8 Mtons, or about 1.2% of the 4500 Mtons CO2 produced in the EU every year37.It should be stressed that carbon dioxide sequestration by the biomass in green roofs is effective only if residues are not significantly degraded. This may be achieved by removing the biomass periodically before it undergoes respiration and mineralization. One could alternatively employ woody plants with a higher carbon accumulation capacity instead of herbaceous vegetation. Although our calculations are referred to a herbaceous annual crop, the results in terms of dry biomass would not be radically different had we considered a tree or shrub crop, as the dry matter potentially produced per unit surface is relatively independent of the plant49. On the other hand, trees and shrubs may be expected to have higher evapotranspiration, thus enhancing the benefits quantified here for a herbaceous crop.If greening is implemented on about 35% of the impervious urban areas, we expect a reduction of runoff in the order of 17.5% compared to the total. Considering that pollutant loads associated to runoff are estimated in the order of about 30 million population equivalents (PE) in terms of biochemical oxygen demand (BOD), about 18 million PE in terms of total nitrogen and about 6 million PE in terms of total phosphorus 6,35, this can be a sizable contribution to the treatment of pollution from European urban areas. Besides the reduction of runoff volume, greened surfaces may also help reduce the frequency of combined sewer overflows because they buffer runoff and release it more slowly than impervious surfaces. This effect is arguably more important for smaller storm events, and tends to disappear as events cause the saturation of green roof storage.It should be stressed that the above analysis considers a soil thickness of 30 cm on greened surfaces. Using the meta-models proposed in30 for the thickness of 10 cm we obtain a ratio between the indicators for thickness of 10 and 30 cm ranging between 80–97% for the reduction of surface temperatures, 55–57% for roof temperatures, 47–57% for biomass, and 84–86% for runoff. Soil thickness affects in particular the roof temperature, due to the associated thermal insulation effect, and the biomass, because a thicker soil can store a larger amount of water and allows a higher evapotranspiration for vegetation growth, while not impeding root growth. A comparison of different climate scenarios sheds light on the sensitivity of our results to the input climatic predictors (P and ET0). From Table 1, it can be calculated that the range (difference between the maximum and minimum value) of precipitation and potential evapotranspiration, as a percentage of the average value, is 20.3%, and 21.4% respectively. The corresponding ranges are 7% of the average for the cooling reduction, 3.7% for the reduced carbon dioxide emission, and 34% for the runoff reduction. The curves in Fig. 2 visualize the relatively small sensitivity of results to the climatic scenario.Economic implicationsMost of the benefits of green roofs are collective. Only a few (e.g. energy saving in summer, and gardening) have an apparent private nature. The costs of greening roofs, on the contrary, are primarily borne by the private owners50. It has been observed that, in the absence of specific incentives, green roof implementation can be economically convenient only for specific commercial and multifamily buildings25. Therefore, private investments should be encouraged through appropriate fiscal and funding policies if the objective is to facilitate a mainstream uptake of this solution. In this section, an indicative cost-benefit analysis is carried out in order to shed light on the possible financing needs at stake, and considering to green the impervious surfaces covered by roofs.The two main benefits that can be easily monetized are the avoided cost of cooling in summer (based on energy prices) and the reduction of carbon dioxide emissions (based on greenhouse gas emissions market prices). By summing the results of Eq. (5) for all gridcells in Europe where the greened surface is assumed to be 35% of the impervious urban area in the gridcell, cooling savings can reach 18.4 billion Є each year for the current scenario. For comparison, the current expenditure for residential cooling in summer can be assumed to be 78 billion Є year−1, based on an electricity use of 391 TWh51. Therefore, the cooling energy saving is 23.5% (18.4 billion Є/78 billion Є), in agreement with the results of Manso et al.15 for the value of 15% estimated for the hot-summer Mediterranean climate.At the present carbon market price of 22.5 Є tons−1 (Ruf and Mazzoni43), the annual benefit related to the estimated reduction of greenhouse gas emissions corresponds to about 1.26 billion Є. It should be stressed how this is apparently an upper limit of this benefit, because not all greened surfaces may correspond to roofs of cooled building volumes, and because the biomass is likely to undergo at least a partial mineralization if not timely removed from the green surfaces. The benefit associated to the reduction of the heat island effect can also be quantified to some extent on the basis of existing literature studies, although their estimation is very complex and would require ad hoc studies. For example, for the city of Phoenix, this benefit was quantified in 80 € for 1 °C decrease per working resident, considering costs of electronic devices, maintenance of cars and performance of cooling47. In another analysis for the Melbourne area, the annual cost was quantified in 18 € per inhabitant, including health, transport, social distress, electric grid faults and damages to animal and trees48. In Malaysia, the annual cost of hazes, related to the urban heat island, was quantified in 12 € per habitant in 1997, including cost of illness, productivity loss, flight cancellation, tourism reduction, decline in fish landings, fire-fighting, cloud seeding and masks46. Therefore, costs can vary significantly among different contexts. Assuming conservatively a yearly benefit of 20 € for each of the ca. 559.5 million European urban inhabitants living in urban areas (75% of the total52), the Net Present Value (NPV) of this benefit over 40 years would be 221 billion € using a discount rate of 4%.The cost of greening the roofs or other impervious surfaces is more difficult to quantify as it depends on several design details and site-specific conditions. For example, in Finland the cost ranges between 70 and 80 Є m−2, in Germany between 13 and 41 Є m−2, while in Switzerland around 20 Є m−253. Assuming an average unit cost of 50 Є m−2, the costs to turn 26,450 km2 of impervious urban areas in Europe into green surfaces amounts to 1323 billion Є. This corresponds to an annual cost (discount rate 4%, 40 years life) of 63 billion euro. This means a cost of 6.3 € m−3 of annual runoff saved (assuming an average annual runoff saving of 10 km3), which is reasonably in line with an estimate of 9.2 € m−3 for the U.S. context, where the annual runoff volume reduction was 12%54 compared to our estimate of 17.5%.Assuming a lifespan of 40 years55 and a discount rate of 4%50, the NPV of the cost saving of summer cooling over 40 years (18.4 billion Є year−1 in Table 1), that is the main private benefit of a green roof installed in a private building, is 364 billion Є (using a discount rate of 4%). The benefits of CO2 reduction, monetized in an emission trading system, would lead to a NPV of 24.85 ≈ 25 billion Є over 40 years (55.8 Mtons year−1). The NPV of the heat island benefit over 40 years would be 221 billion €. Deducting the sum of these benefits (totalling 610 billion €) from the estimated investment of 1323 billion €, yields a net gap of 713 billion Є, corresponding to an annual cost of about 60 € for each of the 559.5 million European citizens living in urban areas. This estimated annual cost is apparently affected by the uncertainty on green roof costs: it could reduce to 4 Є/year per urban citizen if the cost of the green roof is 25 Є m−2, and 129 Є/year per urban citizen if the cost is 80 Є m−2. An annual cost of 60 Є/year per urban citizen may be in many cases compensated by the additional benefits not quantified here. For example, the average increase of property value (rental prices) was estimated to be 8%15. Other benefits can be associated e.g. to leisure and recreation, socialization, amenity of the urban environment, and the creation of habitat or ecological connections in urban areas, besides the abovementioned positive effects in terms of water pollution and floods. Table 2 summarizes the economic results. Table 2 Summary of benefits and costs of urban greening considered in this study for the European context.Full size tableThe harvesting of runoff is a potential additional benefit, but it also entails costs. These can be quantified as a first approximation considering a cost of the storage volume Cs = 50 € m−3, a lifetime of the storage of 100 years, a discount rate of 4% and annual operation and maintenance costs of 3% of the investment. For a unit greened surface, the runoff potentially harvested equals P-RR and can be computed from Eq. 3, while the required storage volume to harvest it is given by Eq. 6. The cost of harvesting one m3 of runoff (marginal harvesting costs) follows from the abovementioned costing parameters. Figure 5 depicts the cumulate value of runoff as a function of the marginal harvesting cost. It can be seen that about 75% of the runoff can be harvested with marginal costs below 0.7 € m−3, a value compatible with urban water prices usually applied in Europe. Cs may be lower than 50 € m−3 , but often it may also be higher. Hence our calculation can be only regarded as a first indication and is accurate not more than within one order of magnitude. The quality of water from green surface runoff harvesting is arguably adequate for non-potable domestic use, but depends on the type of green roof and vegetation13. Figure 5Cumulate runoff versus the cost of storage per unit of runoff, for a storage cost of 50 € m−3. Different climatic scenarios are shown.Full size image More

1.Carpenter, S. R. et al. Nonpoint pollution of surface waters with phosphorus and nitrogen. Ecol. Appl. 8, 559–568 (1998).Article 

Google Scholar 
2.Galloway, J. N. et al. Transformation of the nitrogen cycle: Recent trends, questions, and potential solutions. Science 320, 889–892 (2008).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
3.LeBauer, DavidS. & Treseder, K. K. Nitrogen limitation of net primary productivity in terrestrial ecosystems is globally distributed. Ecology 89, 371–379 (2008).PubMed 
PubMed Central 
Article 

Google Scholar 
4.Fay, P. A. et al. Grassland productivity limited by multiple nutrients. Nat. Plants 1, 1–5 (2015).Article 
CAS 

Google Scholar 
5.Yue, K. et al. Stimulation of terrestrial ecosystem carbon storage by nitrogen addition: a meta-analysis. Sci. Rep. 6, 1–10 (2016).Article 
CAS 

Google Scholar 
6.Avolio, M. L. et al. Changes in plant community composition, not diversity, during a decade of nitrogen and phosphorus additions drive above-ground productivity in a tallgrass prairie. J. Ecol. 102, 1649–1660 (2014).CAS 
Article 

Google Scholar 
7.Isbell, F. et al. Nutrient enrichment, biodiversity loss, and consequent declines in ecosystem productivity. Proc. Natl Acad. Sci. USA 110, 11911–11916 (2013).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
8.Van der Putten, W. H., Bradford, M. A., Pernilla Brinkman, E., van de Voorde, T. F. J. & Veen, G. F. Where, when and how plant–soil feedback matters in a changing world. Funct. Ecol. 30, 1109–1121 (2016).Article 

Google Scholar 
9.Revillini, D., Gehring, C. A. & Johnson, N. C. The role of locally adapted mycorrhizas and rhizobacteria in plant–soil feedback systems. Funct. Ecol. 30, 1086–1098 (2016).Article 

Google Scholar 
10.Peay, K. G., Kennedy, P. G. & Talbot, J. M. Dimensions of biodiversity in the Earth mycobiome. Nat. Rev. Microbiol. 14, 434–447 (2016).CAS 
PubMed 
Article 

Google Scholar 
11.Semchenko, M. et al. Fungal diversity regulates plant-soil feedbacks in temperate grassland. Sci. Adv. 4, eaau4578 (2018).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
12.Větrovsky, T. et al. A meta-analysis of global fungal distribution reveals climate-driven patterns. Nat. Commun. 10, 1–9 (2019).Article 
CAS 

Google Scholar 
13.Smith, S. E. & Read, D. J. Mycorrhizal Symbiosis. (Academic Press, 2008).14.Johnson, N. C. Resource stoichiometry elucidates the structure and function of arbuscular mycorrhizas across scales. New Phytol. 185, 631–647 (2010).CAS 
PubMed 
Article 

Google Scholar 
15.Velásquez, A. C., Castroverde, C. D. M. & He, S. Y. Plant–pathogen warfare under changing climate conditions. Curr. Biol. 28, R619–R634 (2018).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
16.Mangan, S. A. et al. Negative plant–soil feedback predicts tree-species relative abundance in a tropical forest. Nature 466, 752–755 (2010).ADS 
CAS 
PubMed 
Article 

Google Scholar 
17.Reynolds, H. L., Packer, A., Bever, J. D. & Clay, K. Grassroots ecology: plant-microbe-soil interactions as drivers of plant community structure and dynamics. Ecology 84, 2281–2291 (2003).Article 

Google Scholar 
18.Veresoglou, S. D., Barto, E. K., Menexes, G. & Rillig, M. C. Fertilization affects severity of disease caused by fungal plant pathogens. Plant Pathol. 62, 961–969 (2013).Article 

Google Scholar 
19.Walters, D. R. & Bingham, I. J. Influence of nutrition on disease development caused by fungal pathogens: implications for plant disease control. Ann. Appl. Biol. 151, 307–324 (2007).CAS 
Article 

Google Scholar 
20.Knorr, M., Frey, S. D. & Curtis, P. S. Nitrogen additions and litter decomposition: a meta-analysis. Ecology 86, 3252–3257 (2005).Article 

Google Scholar 
21.Chai, Y. et al. Patterns of taxonomic, phylogenetic diversity during a long-term succession of forest on the Loess Plateau, China: insights into assembly process. Sci. Rep. 6, 27087 (2016).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
22.Crowther, T. W. et al. Sensitivity of global soil carbon stocks to combined nutrient enrichment. Ecol. Lett. 22, 936–945 (2019).CAS 
PubMed 
Article 

Google Scholar 
23.Fogg, K. The effect of added nitrogen on the rate of decomposition of organic matter. Biol. Rev. 63, 433–462 (1988).Article 

Google Scholar 
24.Bonner, M. T. et al. Why does nitrogen addition to forest soils inhibit decomposition? Soil Biol. Biochem. 137, 107570 (2019).CAS 
Article 

Google Scholar 
25.Zak, D. R. et al. Anthropogenic N deposition, fungal gene expression, and an increasing soil carbon sink in the Northern Hemisphere. Ecology 100, 1–8 (2019).Article 

Google Scholar 
26.Hobbie, S. E. et al. Response of decomposing litter and its microbial community to multiple forms of nitrogen enrichment. Ecol. Monogr. 82, 389–405 (2012).Article 

Google Scholar 
27.Leff, J. W. et al. Consistent responses of soil microbial communities to elevated nutrient inputs in grasslands across the globe. Proc. Natl Acad. Sci. USA 112, 10967–10972 (2015).ADS 
CAS 
PubMed 
Article 

Google Scholar 
28.Nguyen, N. H. et al. FUNGuild: an open annotation tool for parsing fungal community datasets by ecological guild. Fungal Ecol. 20, 241–248 (2016).Article 

Google Scholar 
29.Borer, E. T. et al. Finding generality in ecology: a model for globally distributed experiments. Methods Ecol. Evol. 5, 65–73 (2014).Article 

Google Scholar 
30.Barberán, A., Bates, S. T., Casamayor, E. O. & Fierer, N. Using network analysis to explore co-occurrence patterns in soil microbial communities. ISME J. 6, 343–351 (2012).PubMed 
Article 
CAS 
PubMed Central 

Google Scholar 
31.MacKinnon, D. P., Krull, J. L. & Lockwood, C. M. Equivalence of the mediation, confounding and suppression effect. Prev. Sci. 1, 173–181 (2000).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
32.Kiers, E. T. et al. Reciprocal rewards stabilize cooperation in the mycorrhizal symbiosis. Science 333, 880–882 (2011).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
33.Lekberg, Y. et al. Relative importance of competition and plant–soil feedback, their synergy, context dependency and implications for coexistence. Ecol. Lett. 21, 1268–1281 (2018).PubMed 
Article 
PubMed Central 

Google Scholar 
34.Hart, M. M. & Reader, R. J. Taxonomic basis for variation in the colonization strategy of arbuscular mycorrhizal fungi. New Phytol. 153, 335–344 (2002).Article 

Google Scholar 
35.Johnson, N. C. Can fertilization of soil select less mutualistic mycorrhizae? Ecol. Appl. 3, 749–757 (1993).PubMed 
Article 
PubMed Central 

Google Scholar 
36.Weber, S. E. et al. Responses of arbuscular mycorrhizal fungi to multiple coinciding global change drivers. Fungal Ecol. 40, 62–71 (2019).Article 

Google Scholar 
37.Han, Y., Feng, J., Han, M. & Zhu, B. Responses of arbuscular mycorrhizal fungi to nitrogen addition: a meta-analysis. Glob. Change Biol. 26, 7229–7241 (2020).ADS 
Article 

Google Scholar 
38.Treseder, K. K. et al. Arbuscular mycorrhizal fungi as mediators of ecosystem responses to nitrogen deposition: A trait- ­based predictive framework. J. Ecol. 106, 480–489 (2018).39.Sikes, B. A., Cottenie, K. & Klironomos, J. N. Plant and fungal identity determines pathogen protection of plant roots by arbuscular mycorrhizas. J. Ecol. 97, 1274–1280 (2009).Article 

Google Scholar 
40.Maherali, H. & Klironomos, J. N. Influence of phylogeny on fungal community assembly and ecosystem functioning. Science 316, 1746–1748 (2007).ADS 
CAS 
PubMed 
Article 

Google Scholar 
41.Johnson, N. C., Graham, J. H. & Smith, F. A. Functioning of mycorrhizal associations along the mutualism – parasitism continuum. New Phytol. 135, 575–585 (1997).42.Balser, T. C., Treseder, K. K. & Ekenler, M. Using lipid analysis and hyphal length to quantify AM and saprotrophic fungal abundance along a soil chronosequence. Soil Biol. Biochem. 37, 601–604 (2005).CAS 
Article 

Google Scholar 
43.Cappelli, S. L., Pichon, N. A., Kempel, A. & Allan, E. Sick plants in grassland communities: a growth‐defense trade‐off is the main driver of fungal pathogen abundance. Ecol. Lett. 23, 1349–1359 (2020).PubMed 
Article 
PubMed Central 

Google Scholar 
44.Grman, E. Plant species differ in their ability to reduce allocation to non-beneficial arbuscular mycorrhizal fungi. Ecology 93, 711–718 (2012).PubMed 
Article 
PubMed Central 

Google Scholar 
45.Kardol, P., Martijn, Bezemer, T. & van der Putten, W. H. Temporal variation in plant-soil feedback controls succession. Ecol. Lett. 9, 1080–1088 (2006).PubMed 
Article 
PubMed Central 

Google Scholar 
46.Firn, J. et al. Leaf nutrients, not specific leaf area, are consistent indicators of elevated nutrient inputs. Nat. Ecol. Evol. 3, 400–406 (2019).PubMed 
Article 
PubMed Central 

Google Scholar 
47.Cleland, E. E. et al. Belowground biomass response to nutrient enrichment depends on light limitation across globally distributed grasslands. Ecosystems 22, 1466–1477 (2019).CAS 
Article 

Google Scholar 
48.Adler, P. B. et al. Functional traits explain variation in plant life history strategies. Proc. Natl Acad. Sci. USA 111, 10019–10019 (2014).Article 
CAS 

Google Scholar 
49.Kulmatiski, A., Beard, K. H., Stevens, J. R. & Cobbold, S. M. Plant-soil feedbacks: a meta-analytical review. Ecol. Lett. 11, 980–992 (2008).PubMed 
Article 
PubMed Central 

Google Scholar 
50.Tedersoo, L. et al. Global diversity and geography of soil fungi. Science 346, 1052–1053 (2014).Article 
CAS 

Google Scholar 
51.Davison, J. et al. Global assessment of arbuscular mycorrhizal fungus diversity reveals very low endemism. Science 349, 970–973 (2015).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
52.Morriën, E. et al. Soil networks become more connected and take up more carbon as nature restoration progresses. Nat. Commun. 8, 14349 (2017).53.Karimi, B. et al. Microbial diversity and ecological networks as indicators of environmental quality. Environ. Chem. Lett. 15, 265–281 (2017).CAS 
Article 

Google Scholar 
54.Carr, A., Diener, C., Baliga, N. S. & Gibbons, S. M. Use and abuse of correlation analyses in microbial ecology. ISME J. https://doi.org/10.1038/s41396-019-0459-z (2019).55.Prober, S. M. et al. Plant diversity predicts beta but not alpha diversity of soil microbes across grasslands worldwide. Ecol. Lett. 18, 85–95 (2015).PubMed 
Article 

Google Scholar 
56.Steidinger, B. S. et al. Climatic controls of decomposition drive the global biogeography of forest-tree symbioses. Nature 569, 404–408 (2019).ADS 
CAS 
PubMed 
Article 

Google Scholar 
57.Corradi, N. et al. Gene copy number polymorphisms in an arbuscular mycorrhizal fungal population. Appl. Environ. Microbiol. 73, 366–369 (2007).CAS 
PubMed 
Article 

Google Scholar 
58.Tedersoo, L. et al. Response to Comment on “Global diversity and geography of soil fungi. Science 349, 936 (2015).ADS 
CAS 
PubMed 
Article 

Google Scholar 
59.Malcolm, G. M., Kuldau, G. A., Gugino, B. K. & Jiménez-Gasco, M. D. M. Hidden host plant associations of soilborne fungal pathogens: an ecological perspective. Phytopathology 103, 538–544 (2013).CAS 
PubMed 
Article 

Google Scholar 
60.Taylor, D. L. et al. Accurate estimation of fungal diversity and abundance through improved lineage-specific primers optimized for illumina amplicon sequencing. Appl. Environ. Microbiol. 82, 7217–7226 (2016).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
61.Zou, K., Thebault, E., Lacroix, G. & Barot, S. Interactions between the green and brown food web determine ecosystem functioning. Funct. Ecol. 30, 1454–1465 (2016).Article 

Google Scholar 
62.Chen, W. et al. Fertility‐related interplay between fungal guilds underlies plant richness–productivity relationships in natural grasslands. New Phytol. 226, 1129–1143 (2020).PubMed 
Article 

Google Scholar 
63.Busby, P. E., Peay, K. G. & Newcombe, G. Common foliar fungi of Populus trichocarpa modify Melampsora rust disease severity. New Phytol. 209, 1681–1692 (2016).CAS 
PubMed 
Article 

Google Scholar 
64.Li, X., Ding, C., Zhang, T. & Wang, X. Fungal pathogen accumulation at the expense of plant-beneficial fungi as a consequence of consecutive peanut monoculturing. Soil Biol. Biochem. 72, 11–18 (2014).Article 
CAS 

Google Scholar 
65.Lefcheck, J. S. piecewiseSEM: Piecewise structural equation modelling in r for ecology, evolution, and systematics. Methods Ecol. Evol. 7, 573–579 (2016).Article 

Google Scholar 
66.Friedman, J. & Alm, E. J. Inferring correlation networks from genomic survey data. PLoS Comput. Biol. 8, 1–11 (2012).Article 
CAS 

Google Scholar 
67.Kurtz, Z., Mueller, C., Miraldi, E. & Bonneau, R. SpiecEasi: Sparse Inverse Covariance For Ecological Statistical Inference. R package version 1.0.6 (2019).68.Oksanen, J. et al. vegan: Community Ecology Package. R package (2019).69.Wickham, H. ggplot2: Elegant Graphics for Data Analysis. (Springer-Verlag New York, 2016). More

1.Case, T. J. & Gilpin, M. E. Interference competition and niche theory. Proc. Natl. Acad. Sci. 71, 3073–3077 (1974).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
2.Linnell, J. D. & Strand, O. Interference interactions, co-existence and conservation of mammalian carnivores. Divers. Distrib. 6, 169–176 (2000).Article 

Google Scholar 
3.Prugh, L. R. & Sivy, K. J. Enemies with benefits: Integrating positive and negative interactions among terrestrial carnivores. Ecol. Lett. 23, 902–918 (2020).PubMed 
Article 
PubMed Central 

Google Scholar 
4.Polis, G. A., Myers, C. A. & Holt, R. D. The ecology and evolution of intraguild predation: Potential competitors that eat each other. Annu. Rev. Ecol. Syst. 20, 297–330 (1989).Article 

Google Scholar 
5.Belant, J. L., Griffith, B., Zhang, Y., Follmann, E. H. & Adams, L. G. Population-level resource selection by sympatric brown and American black bears in Alaska. Polar Biol. 33, 31–40 (2010).Article 

Google Scholar 
6.Lima, S. L. & Dill, L. M. Behavioral decisions made under the risk of predation: A review and prospectus. Can. J. Zool. 68, 619–640 (1990).Article 

Google Scholar 
7.Laundré, J. W., Hernández, L. & Altendorf, K. B. Wolves, elk, and bison: Reestablishing the “landscape of fear” in Yellowstone National Park, USA. Can. J. Zool. 79, 1401–1409 (2001).Article 

Google Scholar 
8.Moll, R. J. et al. The many faces of fear: A synthesis of the methodological variation in characterizing predation risk. J. Anim. Ecol. 86, 749–765 (2017).PubMed 
Article 
PubMed Central 

Google Scholar 
9.Kohl, M. T. et al. Diel predator activity drives a dynamic landscape of fear. Ecol. Monogr. 88, 638–652 (2018).Article 

Google Scholar 
10.Kuijper, D. P. J. et al. Landscape of fear in Europe: Wolves affect spatial patterns of ungulate browsing in Białowieża Primeval Forest, Poland. Ecography 36, 1263–1275 (2013).Article 

Google Scholar 
11.Smith, J. A., Donadio, E., Pauli, J. N., Sheriff, M. J. & Middleton, A. D. Integrating temporal refugia into landscapes of fear: Prey exploit predator downtimes to forage in risky places. Oecologia 189, 883–890 (2019).ADS 
PubMed 
Article 
PubMed Central 

Google Scholar 
12.Flagel, D. G., Belovsky, G. E. & Beyer, D. E. Natural and experimental tests of trophic cascades: Gray wolves and white-tailed deer in a Great Lakes forest. Oecologia 180, 1183–1194 (2016).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
13.Gaynor, K. M., Brown, J. S., Middleton, A. D., Power, M. E. & Brashares, J. S. Landscapes of fear: Spatial patterns of risk perception and response. Trends Ecol. Evol. 34, 355–368 (2019).PubMed 
Article 
PubMed Central 

Google Scholar 
14.Prugh, L. R. et al. Designing studies of predation risk for improved inference in carnivore-ungulate systems. Biol. Conserv. 232, 194–207 (2019).Article 

Google Scholar 
15.Fisher, J. T., Anholt, B., Bradbury, S., Wheatley, M. & Volpe, J. P. Spatial segregation of sympatric marten and fishers: The influence of landscapes and species-scapes. Ecography 36, 240–248 (2013).Article 

Google Scholar 
16.Manlick, P. J., Woodford, J. E., Zuckerberg, B. & Pauli, J. N. Niche compression intensifies competition between reintroduced American martens (Martes americana) and fishers (Pekania pennanti). J. Mammal. 98, 690–702 (2017).Article 

Google Scholar 
17.Powell, R. A., Buskirk, S. W., & Zielinski, W. J. Fisher and marten. In Wild Mammals of North America: Biology, Management, and Conservation (eds. Feldhamer, G. A et al.), 635–649 (JHU Press, 2003).18.Krohn, W. B., Elowe, K. D. & Boone, R. B. Relations among fishers, snow, and martens: Development and evaluation of two hypotheses. For. Chron. 71, 97–105 (1995).Article 

Google Scholar 
19.Williams, B. W., Gilbert, J. H., & Zollner, P. A. Historical Perspective on the Reintroduction of the Fisher and American Marten in Wisconsin and Michigan, vol. 5. (US Department of Agriculture, Forest Service, Northern Research Station, 2007).20.McCann, N. P., Zollner, P. A. & Gilbert, J. H. Survival of adult martens in northern Wisconsin. J. Wildl. Manag. 74, 1502–1507 (2010).Article 

Google Scholar 
21.Kupferman, C. A. An Expanding Meso-Carnivore: Fisher (Pekania pennanti) Occupancy and Coexistence with Native Mustelids in Southeast Alaska (University of Idaho, 2019).
Google Scholar 
22.Hall, L. K. et al. Vigilance of kit foxes at water sources: A test of competing hypotheses for a solitary carnivore subject to predation. Behav. Proc. 94, 76–82 (2013).Article 

Google Scholar 
23.Chitwood, M. C., Lashley, M. A., Higdon, S. D., DePerno, C. S. & Moorman, C. E. Raccoon vigilance and activity patterns when sympatric with coyotes. Diversity 12, 341 (2020).Article 

Google Scholar 
24.Vanak, A. T., Thaker, M. & Gompper, M. E. Experimental examination of behavioural interactions between free-ranging wild and domestic canids. Behav. Ecol. Sociobiol. 64, 279–287 (2009).Article 

Google Scholar 
25.Croose, E., Bled, F., Fowler, N. L., Beyer, D. E. Jr. & Belant, J. L. American marten and fisher do not segregate in space and time during winter in a mixed-forest system. Ecol. Evol. 9, 4906–4916 (2019).PubMed 
PubMed Central 
Article 

Google Scholar 
26.Gilbert, J. H., Zollner, P. A., Green, A. K., Wright, J. L. & Karasov, W. H. Seasonal field metabolic rates of American martens in Wisconsin. Am. Midl. Nat. 162, 327–334 (2009).Article 

Google Scholar 
27.Hughes, N. K., Price, C. J. & Banks, P. B. Predators are attracted to the olfactory signals of prey. PLoS ONE 5, e13114 (2010).ADS 
PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
28.Bytheway, J. P., Carthey, A. J. & Banks, P. B. Risk vs. reward: How predators and prey respond to aging olfactory cues. Behav. Ecol. Sociobiol. 67, 715–725 (2013).Article 

Google Scholar 
29.Haynes, G. Utilization and skeletal disturbances of North American prey carcasses. Arctic 35, 266–281 (1982).Article 

Google Scholar 
30.Kaufmann, J. H. On the definitions and functions of dominance and territoriality. Biol. Rev. 58, 1–20 (1983).Article 

Google Scholar 
31.Zielinski, W. J., Tucker, J. M. & Rennie, K. M. Niche overlap of competing carnivores across climatic gradients and the conservation implications of climate change at geographic range margins. Biol. Conserv. 209, 533–545 (2017).Article 

Google Scholar 
32.Jensen, P. G. & Humphries, M. M. Abiotic conditions mediate intraguild interactions between mammalian carnivores. J. Anim. Ecol. 88, 1305–1318 (2019).PubMed 
Article 
PubMed Central 

Google Scholar 
33.Manlick, P. J., Windels, S. K., Woodford, J. E. & Pauli, J. N. Can landscape heterogeneity promote carnivore coexistence in human-dominated landscapes?. Landsc. Ecol. 35, 2013–2027 (2020).Article 

Google Scholar 
34.Krohn, W., Hoving, C., Harrison, D., Phillips, D., & Frost, H. Martes foot-loading and snowfall patterns in eastern North America. In Martens and Fishers (Martes) in Human-Altered Environments (eds. Harrison, D. J. et al.) 115–131 (Springer, 2005).35.Hiller, T. L., Etter, D. R., Belant, J. L. & Tyre, A. J. Factors affecting harvests of fishers and American martens in northern Michigan. J. Wildl. Manag. 75, 1399–1405 (2011).Article 

Google Scholar 
36.Childress, M. J. & Lung, M. A. Predation risk, gender and the group size effect: Does elk vigilance depend upon the behaviour of conspecifics?. Anim. Behav. 66, 38–398 (2003).Article 

Google Scholar 
37.Gehr, B. et al. Stay home, stay safe—Site familiarity reduces predation risk in a large herbivore in two contrasting study sites. J. Anim. Ecol. 89, 1329–1339 (2020).PubMed 
Article 
PubMed Central 

Google Scholar 
38.Bull, E. L. & Heater, T. W. Survival, causes of mortality, and reproduction in the American marten in northeastern Oregon. Northwest. Nat. 82, 1–6 (2001).Article 

Google Scholar 
39.White, K. S., Golden, H. N., Hundertmark, K. J. & Lee, G. R. Predation by wolves, Canis lupus, on wolverines, Gulo gulo, and an American marten, Martes americana, Alaska. Can. Field Nat. 116, 132–134 (2002).
Google Scholar 
40.Erb, J., Sampson, B., & Coy, P. Survival and causes of mortality for fisher and marten in Minnesota. Minnesota Department of Natural Resources Summary of Wildlife Research Findings, 2009, 24–31 (2009).41.Wengert, G. M., Gabriel, M. W., Foley, J. E., Kun, T. & Sacks, B. N. Molecular techniques for identifying intraguild predators of fishers and other North American small carnivores. Wildl. Soc. Bull. 37, 659–663 (2013).
Google Scholar 
42.Stricker, H. K. et al. Use of modified snares to estimate bobcat abundance. Wildl. Soc. Bull. 36, 257–263 (2012).Article 

Google Scholar 
43.Kautz, T. M. et al. Predator densities and white-tailed deer fawn survival. J. Wildl. Manag. 83, 1261–1270 (2019).Article 

Google Scholar 
44.Caravaggi, A. et al. A review of camera trapping for conservation behaviour research. Remote Sens. Ecol. Conserv. 3, 109–122 (2017).Article 

Google Scholar 
45.Berger, K. M. & Gese, E. M. Does interference competition with wolves limit the distribution and abundance of coyotes?. J. Anim. Ecol. 76, 1075–1085 (2007).PubMed 
Article 
PubMed Central 

Google Scholar 
46.Merkle, J. A., Stahler, D. R. & Smith, D. W. Interference competition between gray wolves and coyotes in Yellowstone National Park. Can. J. Zool. 87, 56–63 (2009).Article 

Google Scholar 
47.Crimmins, S. M. & Van Deelen, T. R. Limited evidence for mesocarnivore release following wolf recovery in Wisconsin, USA. Wildl. Biol. 2019, 1–7 (2019).Article 

Google Scholar 
48.Petroelje, T. R., Belant, J. L., Beyer, D. E., & Kautz, T. M. Interference competition between wolves and coyotes during variable prey abundance. Ecol. Evol 11, 1413–1431 (2021).
49.Switalski, T. A. Coyote foraging ecology and vigilance in response to gray wolf reintroduction in Yellowstone National Park. Can. J. Zool. 81, 985–993 (2003).Article 

Google Scholar 
50.Hilborn, A. et al. Cheetahs modify their prey handling behavior depending on risks from top predators. Behav. Ecol. Sociobiol. 72, article 74 (2018).Article 

Google Scholar 
51.Elgar, M. A. Predator vigilance and group size in mammals and birds: A critical review of the empirical evidence. Biol. Rev. 64, 13–33 (1989).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
52.Bøving, P. S. & Post, E. Vigilance and foraging behaviour of female caribou in relation to predation risk. Rangifer 17, 55–63 (1997).Article 

Google Scholar 
53.Hunter, L. T. B. & Skinner, J. D. Vigilance behaviour in African ungulates: The role of predation pressure. Behaviour 135, 195–211 (1998).Article 

Google Scholar 
54.Liley, S. & Creel, S. What best explains vigilance in elk: Characteristics of prey, predators, or the environment?. Behav. Ecol. 19, 245–254 (2008).Article 

Google Scholar 
55.Makin, D. F., Chamaillé-Jammes, S. & Shrader, A. M. Herbivores employ a suite of antipredator behaviours to minimize risk from ambush and cursorial predators. Anim. Behav. 127, 225–231 (2017).Article 

Google Scholar 
56.Wikenros, C., Ståhlberg, S. & Sand, H. Feeding under high risk of intraguild predation: Vigilance patterns of two medium-sized generalist predators. J. Mammal. 95, 862–870 (2014).Article 

Google Scholar 
57.Welch, R. J., le Roux, A., Petelle, M. B. & Périquet, S. The influence of environmental and social factors on high-and low-cost vigilance in bat-eared foxes. Behav. Ecol. Sociobiol. 72, article 29 (2018).Article 

Google Scholar 
58.Yang, L. et al. A new generation of the United States National Land Cover Database: Requirements, research priorities, design, and implementation strategies. ISPRS J. Photogramm. Remote. Sens. 146, 108–123 (2018).ADS 
Article 

Google Scholar 
59.Lovallo, M. J. & Anderson, E. M. Bobcat (Lynx rufus) home range size and habitat use in northwest Wisconsin. Am. Midl. Nat. 135, 241–252 (1996).Article 

Google Scholar 
60.Burton, A. C. et al. Wildlife camera trapping: A review and recommendations for linking surveys to ecological processes. J. Appl. Ecol. 52, 675–685 (2015).Article 

Google Scholar 
61.Bates, D., Mächler, M., Bolker, B. & Walker, S. Fitting linear mixed-effects models using lme4. J. Stat. Softw. 67, 1–48 (2015).Article 

Google Scholar 
62.R Core Team. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. http://www.R-project.org/. Accessed Aug 2020.63.National Operational Hydrologic Remote Sensing Center. Snow Data Assimilation System (SNODAS) Data Products at NSIDC, Version 1. Boulder, Colorado USA. https://doi.org/10.7265/N5TB14TC (NSIDC: National Snow and Ice Data Center, 2004).64.Hutchings, M. R. & White, P. C. Mustelid scent-marking in managed ecosystems: Implications for population management. Mammal Rev. 30, 157–169 (2000).Article 

Google Scholar 
65.Mumm, C. A., & Knörnschild, M. Mustelid Communication. In Encyclopedia of Animal Cognition and Behavior (ed. Choe, J.), 1–11 (Springer International, 2018).66.Sullivan, T. P., Nordstrom, L. O. & Sullivan, D. S. Use of predator odors as repellents to reduce feeding damage by herbivores. J. Chem. Ecol. 11, 903–919 (1985).CAS 
PubMed 
Article 

Google Scholar 
67.Rowcliffe, J. M., Kays, R., Kranstauber, B., Carbone, C. & Jansen, P. A. Quantifying levels of animal activity using camera trap data. Methods Ecol. Evol. 5, 1170–1179 (2014).Article 

Google Scholar  More

Tractable systemsTractable systems are models typically of known and reduced complexity that can be operationalized and reproduced over relatively short periods of time. To formalize our study using tractable systems, we consider a (regional) pool ({mathcal{R}} =left{1,2,cdots ,Sright}) of (S) resident species and one non-resident species denoted by “(I)” (Fig. 1A). We denote the resident community by ({mathcal{M}} subseteq {mathcal{R}}), which is the species collection that coexists obtained by assembling all resident species simultaneously. Additionally, let ({ {mathcal{M}} }_{p}) denote the perturbed community formed by the species collection that coexists when assembling all residents and the non-resident species simultaneously assuming the possibility of multiple introductions. Note that this mechanism corresponds to a top-down assembly process.33 Then, the non-resident species is classified as a game-changing species if it changes the number of resident species that coexist: (|{ {mathcal{M}} }_{p}{I}|,ne, | {mathcal{M}} |). Figure 1A illustrates the concept of a game-changing species in a hypothetical microbial community of (S=2) resident species. Here, one possible context for the resident species is that one species excludes the other, e.g., ({mathcal{M}} ={1}) (Fig. 1B). In this case, the non-resident species is a game-changing species if it promotes the establishment of the other resident (Fig. 1C). Note that we do not consider a change when eliminating the current resident species. The other context for the resident species is that both coexist ({mathcal{M}} ={1,2}) (Fig. 1D). In this case, a non-resident species is a game-changing species if it suppresses the establishment of at least one of the residents, e.g., ({ {mathcal{M}} }_{p}={1,I}) (Fig. 1E). Note that a game-changing species can be either colonizer or transient depending on the dynamics.Fig. 1: Game-changing species for a resident microbial community.Illustration of different contexts leading to game-changing and non-game-changing species. Panel A shows a hypothetical microbial community with a pool ({mathcal{R}} ={1,2}) of two resident species (pink and yellow) and one non-resident species “(I)” (green). Panel B shows the context when one species excludes the other, the resident community contains a single resident (({mathcal{M}} ={1})). To change the resident community, the non-resident species needs to promote the establishment of the other species (in this case yellow, but the example can also be done for the pink species). Panel C provides examples of game-changing and non-game-changing non-resident species for the example presented in Panel B (as the outcomes of the perturbed communities ({ {mathcal{M}} }_{p})). Panel D shows another context when the two resident species coexist (({mathcal{M}} ={1,2})). To change the resident community, the non-resident species needs to suppress the establishment of any of the species (green or yellow). Panel E provides potential outcomes of perturbed communities of game-changing and non-game-changing non-resident species for the example presented in Panel D. In this context, the change happens by suppressing the yellow species (the same can be said for the pink species). Note that in all contexts, the non-resident species can be either a colonizer (can become established in the perturbed community) or transient (cannot colonize).Full size imageWe study the generalization of game-changing species under controlled conditions using in vitro experimental soil communities and under changing conditions using in vivo gut microbial communities (see SI for details about these experimental systems). Note that in vitro experiments usually create ad hoc conditions for species by putting them outside of their natural changing habitats, assuring that species survive in monocultures, and by forming interspecific interactions that may not occur otherwise. Instead, in vivo experiments are performed within living systems, resembling much closer the natural habitat of species and their interspecific interactions. Both types of systems are of reduced complexity, allowing the monitoring and reproducibility of experiments. The studied in vitro soil communities are formed by experimental trials of eight interacting heterotrophic soil-dwelling microbes:30 Enterobacter aerogenes, Pseudomonas aurantiaca, Pseudomonas chlororaphis, Pseudomonas citronellolis, Pseudomonas fluorescens, Pseudomonas putida, Pseudomonas veronii, and Serratia marcescens. These experiments were performed by co-inoculating species at different growth–dilution cycles into fresh media. Each species was cultured in isolation. All experiments were carried out in duplicate. The studied in vivo gut communities are formed by experimental trails of five interacting microbes commonly found in the fruit fly Drosophila melanogoster gut microbiota:34 Lactobacillus plantarum, Lactobacillus brevis, Acetobacter pasteurianus, Acetobacter tropicalis, and Acetobacter orientalis. These experiments were performed by co-inoculating species through frequent ingestion in different flies. All experiments were replicated at least 45 times. These two data sets are, to our knowledge, the closest and best described systems of two- and three-species communities currently available describing species coexistence (not just presence/absence records) under two contrasting environmental conditions.Focusing on in vitro communities, we studied all 28 pairs and 56 trios formed by the eight soil species.30 This provided 168 cases, where it is possible to investigate the expected result (({ {mathcal{M}} }_{p})) of assembling a non-resident species together with a resident community (21 cases for each of the 8 studied microbes). The overall competition time was chosen such that species extinctions would have sufficient time to occur, while new mutants would typically not have time to arise and spread. Similarly for the in vivo communities, we studied all 10 pairs and 10 trios formed by the 5 gut species, which provided 30 cases equivalent to the soil experiments.34 Because species extinctions in in vivo communities are harder to establish, we classified as an expected extinction to any species whose relative abundance was less than 10% in at least 71% of all (47–49) replicates, which corresponds to less than 1% of cases under a binomial distribution with (p=0.5) (slightly different thresholds produce qualitatively similar results). Each of these 168 and 30 cases for soil and gut communities, respectively, represents a given resident community (({mathcal{M}})) formed by a pool of two resident species (({mathcal{R}} ={1,2}), ({mathcal{M}} subseteq {mathcal{R}})) where the target for a non-resident species ((I)) can be either to promote or to suppress the establishment of resident species (({mathcal{R}} {cup }left{Iright}), ({mathcal{M}}_{p} setminus {I},ne, {mathcal{M}})). Non-resident species that are expected to survive in ({ {mathcal{M}} }_{p}) are classified as colonizers; otherwise they are classified as transients.Empirical contextsTo investigate the role of context dependency in the game-changing capacity of microbial species, we study the extent to which a given species can be classified as a game-changer regardless of the resident community it interacts with or if it is the resident community that provides the opportunity for a non-resident species to be a game-changer. Specifically, for each species, we calculate the fraction of times such a species changes the resident community conditioned on the type (whether it is a colonizer or a transient) and target (whether promoting or suppressing). Then, we calculate the probability ((p) value) of observing a fraction greater than or equal to the observed fraction under the given type/target (using a one-sided binomial test with mean value given by the empirical frequency within each type/target). High (p) values (e.g., ( > 0.05)) would be indicative of the importance of context-dependency and not of the intrinsic capacity of species.Next, we quantify the average effect of empirical contexts shaping the game-changing capacity of non-resident species. Specifically, we measure the type’s average effect on changing the community using ({E}_{Y}=P(C=1|Y=1)-P(C=1|Y=0)). Here, (C=1) if the non-resident species was a game-changer ((C=0) if it was not), and where (Y=1) if the species was a colonizer ((Y=0) if the species was transient). The non-parametric quantity (P(C|Y)) corresponds to the frequency of observing (C) given (Y). Thus, ({E}_{Y} , > , 0) (resp. ({E}_{Y} , , 0) (resp. (,{E}_{T} , ,0)). Panel H shows that the non-resident species will be able to change the community based on structuralist theory.Full size imageTractable theoretical systems: structuralist theoryWhile mutual invasibility theory has provided key insights regarding population dynamics,29,30,35 it has been shown that it cannot be directly generalized to multispecies communities.37,38,39 Hence, as an alternative potential generalization, we introduce a second heuristic rule based on structuralist theory,32,40,41 Across many areas of biology, the structuralist view has provided a systematic and probabilistic platform for understanding the diversity that we observe in nature,31,42,43 In ecology, structuralist theory assumes that the probability of observing a community is based on the match between the internal constraints established by species interactions (treated as physico-chemical rules of design) within a community and the changing external conditions (treated as unknown conditions).32,44,45 This other premise has also been shown to be as successful as mutual invasibility in predicting the outcome of surviving species in the studied in vitro soil communities,32 but it has not been tested for its generality.Formally, the structuralist framework assumes that the per-capita growth rate of an (i)th species can be approximated by a general phenomenological function ({f}_{i}({N}_{1},cdots ,{N}_{S},{N}_{I};{boldsymbol{theta}})), i.e.,$$frac{d{N}_{i}}{dt}={N}_{i} {f}_{i}({N}_{1},cdots ,{N}_{S},{N}_{I};{mathbf{theta}}),qquad iin {mathcal{R}}{cup}{I}$$
(1)
Above, ({N}_{i}) represents the abundance (or biomass) of species (i). The functions ({f}_{i}) encode the internal constraints of the community dynamics.46 The vector parameter ({mathbf{theta}}) encodes the external (unknown) conditions acting on the community, which can change according to some probability distribution (p({mathbf{theta}})). For a particular value ({mathbf{theta}}={{mathbf{theta}}}^{ast }), a species collection ({mathscr{Z}}subseteq {mathcal{R}} {cup }{I}) is said feasible (potentially observable) for Eq. (1) if there exists equilibrium abundances ({N}_{i}^{ast } , > , 0) for all species (iin {mathscr{Z}}) and ({N}_{i}^{ast }=0) for (i,notin, {mathscr{Z}}) (i.e., ({f}_{i}({N}_{1}^{ast },cdots ,{N}_{S}^{ast },{N}_{I}^{ast };{{mathbf{theta}}}^{ast })=0) for all (i)).40 Then, we can use Eq. (1) to push-forward (p({mathbf{theta}})) and estimate the probability that a randomly chosen species (i) is feasible with ((iin {mathcal{R}} {cup }{I})) and without ((iin {mathcal{R}})) the non-resident species under isotropic changing conditions, respectively. In this form, the effect of a non-resident species (I) on a resident community can be characterized by the expected maximum impact on its feasibility, i.e., ({Delta }_{{{F}}}=p(i| {mathcal{R}} {cup }{I})-p(i| {mathcal{R}} )) (Fig. 2E–H).To make this framework tractable, we leverage on the mathematical properties of the linear Lotka–Volterra (LV) system47 with the per-capita growth rate ({f}_{i}({N}_{1},cdots ,{N}_{S},{N}_{I};{mathbf{theta}})=mathop{sum}limits_{jin {mathcal{R}} {cup }{I}}{a}_{ij}{N}_{j}+{theta}_{i}) for (iin {mathcal{R}} {cup }{I}). While the linear LV system can be interpreted under many different assumptions,47 we follow its most general interpretation as a first-order approximation to Eq. (1).46 In this system, the time-invariant community structure consists of the intraspecific and interspecific species interactions ({bf{A}}=({a}_{ij})in {{mathbb{R}}}^{(S+1)times (S+1)}), and the external factors ({mathbf{theta}}=({theta }_{1},cdots ,{theta }_{S},{theta }_{I})in {{mathbb{R}}}^{S+1}) consist of density-independent intrinsic per-capita growth rates of all species. We assume that (p({boldsymbol{theta }})) is uniform over the positive parameter space (conforming with ergodicity in dynamical systems32) and find analytically the external conditions compatible with the feasibility of a randomly chosen species within a given community ({bf{A}}), i.e., (p(i|{bf{A}})) (see SI). This framework is robust to changes in the system dynamics since (p(i|{bf{A}})) is identical for all systems that are topologically equivalent to the linear LV system32,48 and a lower bound for systems with higher-order terms.41 Note also that while higher-order interactions may impact the dynamics of microbial communities49,50, their incorporation into ecological models as higher-order polynomials rend intractable and super-sensitive systems (no closed-form solutions can be found in terms of radicals)41,51,52,53.To quantify the contribution to feasibility (({Delta }_{{rm{F}}})) of a non-resident species under the structuralist framework defined above, we infer both the resident interaction matrix ({bf{A}}) and the perturbed interaction matrix ({{bf{A}}}_{p}) using only information from experimental monocultures and pairwise cocultures. Interaction matrices were inferred by fitting the linear LV system using the different repetitions of the observed survival data (see SI and Fig. S1). To make the structuralist framework comparable with the mutual invasibility framework, we introduce a heuristic rule based on structuralist theory formalized in a binary variable (({{F}})) that when anticipating the promotion (resp. suppression) of resident species becomes ({{F}}=1) if ({Delta }_{{{F}}} , > , 0); otherwise ({{F}}=0) if ({Delta }_{{{F}}} , More

We successfully sequenced 16S rRNA gene amplicons from 454 samples with universal primers and 283 samples with archaea-specific primers, respectively (see Supplementary Table 2), and recovered 260,266 ASVs in the universal 16S rRNA gene dataset and 28,123 ASVs in the archaea-specific dataset (Supplementary Fig. 2). As samples from the Atacama Trench region included sectioning at higher depth resolution (HR sectioning scheme), we will present and discuss these first.Variability of microbial community composition along the Atacama Trench axisThe sediments of the Atacama Trench showed the shallowest oxygen penetrations found in any hadal trench to date, ranging from 4.1 cm at the southernmost site A6 to 3.1 cm at northernmost A10, and reflected a high organic carbon flux from the Humboldt upwelling system [22, 26]. Correspondingly, nitrate penetration depths ranged from 8 to only 6 cm with dissolved, ferrous iron accumulating below, while hydrogen sulfide was not detected (see Supplementary Fig. 1 and Supplementary Table 1).Comparison of the community structure obtained with universal primers within the Atacama Trench indicated very similar trends with sediment depth for all sites, with a gradual change from the sediment surface to deeper sediment horizons, and with only marginal overlap between individual redox zones (Fig. 2A). Similar patterns were observed using different dissimilarity metrics (Bray Curtis, weighted/unweighted UniFrac), ordination techniques (NMDS/t-SNE), and sectioning schemes (HR/CR; data not shown). The downcore gradient of microbial communities was also evident within individual redox zones. Thus, samples from the same sediment horizon (e.g., 1–2 cm) but from different sites, with geographic distances of up to 430 km, were more similar to each other than to their respective adjacent horizons, above (0–1 cm) or below (2–3 cm) (Supplementary Fig. 3A, B). The horizontal similarity was particularly pronounced in the upper part of the oxic zone, while it decreased toward the bottom of the nitrogenous zone and increased again in the ferruginous zone (Supplementary Fig. 3A, B). The CR sample subset contained triplicates from separate sediment cores originating from two multicorer deployments, and thus included two cores sampled within a distance of 0.1–1 m and one sampled at an estimated distance of 10–100 m from the other two. Triplicate samples from the same sediment horizon were more similar to each other (1 – Bray Curtis dissimilarity ~0.7–0.9) than to samples from the same horizons at other sites (~0.4–0.8; Supplementary Fig. 4A, C). This implies some increase of variability with geographic distance along the trench axis.Fig. 2: Microbial community composition in the Atacama Trench.Principal coordinate analysis (PCoA) of Bray Curtis dissimilarities between hadal samples from the Atacama Trench in the universal 16S rRNA gene (A) and archaea-specific 16S rRNA gene (B) data. The color gradient represents sediment depth and ovals mark the 95% confidence intervals of multivariate normal distributions of the oxic, nitrogenous, and ferruginous zones, respectively. Different symbols correspond to different sites. C Relative read abundance (%) on phylum/class level of the ten most abundant taxonomic groups (universal 16S rRNA gene data) grouped by redox zone and by depth within each zone in hadal samples from the Atacama Trench. Both the color gradient and the number within the squares indicate of the average relative read abundances within the respective sample group.Full size imageThe steep change in community composition with increasing sediment depth yet relatively high similarity of communities from the same sediment depth at different sites was unexpected. Indeed, we expected that the irregular depositional regime of hadal trenches [18, 21], which was indicated in all hadal sediment cores through color layering and by site-specific fluctuations in depth distributions of porosity, TOC content, and cell numbers [22], would leave an imprint in the community. A closer examination of depth trends at individual sites based on Bray Curtis dissimilarity occasionally revealed compositional fluctuations with depth that may be related to depositional events (Supplementary Fig. 5). For instance, at sites A3 and A4 the microbial communities in the oxic zone were more similar to those from the ferruginous zone at around 15–25 and 9–15 cm, respectively, than to the samples of the nitrogenous zone located in between them (Supplementary Fig. 5). This hinted that local depositional events might have entombed parts of the microbial community. We suggest that such events contributed to the enhanced site–site variability in the nitrogenous and ferruginous zones (Supplementary Fig. 4A, C).Compositional changes over sediment depth and redox zonation along the Atacama TrenchThe observed trends in beta diversity were reflected by distinct phylum-level changes that followed redox zonation and sediment depth (Fig. 2C and Supplementary Fig. 6). Conversely, subsurface peaks of microbial abundance [22] were not reflected in the relative abundance patterns of different phyla (Proteobacteria always split to class level yet referred to as phyla for simplicity), which changed steadily with increasing sediment depth. For instance, some of the dominant groups, Gammaproteobacteria (20.9%; mean read abundance per sediment depth), Bacteroidetes (17.1%), and Thaumarchaeota (12.4%), peaked in relative abundance in the oxic zone and then decreased. Alphaproteobacteria (13.7%) on average had the highest relative abundances in the nitrogenous zone and became rather rare in the ferruginous zone. Planctomycetes rose in relative abundance below the oxic zone from around 9 to 15% and remained at this level below. Atribacteria showed the largest change in relative abundance. After being close to detection limit with an average relative abundance of 0.004% in the oxic and nitrogenous zones, their relative abundance increased approximately ten-fold for every centimeter from the transition to the ferruginous zone, until they became the dominant phylum in deeper sediment sections. Aside from Atribacteria, other lineages such as “Candidatus (Ca.) Marinimicrobia” (0.1–4.6%), “Ca. Woesearchaeota” (1.3–11%) and “Ca. Patescibacteria” (0.7–4.7%) increased steadily with increasing sediment depth, while Acidobacteria (~3.4%) and Deltaproteobacteria (~6.1%) showed almost no change. As microbial abundance in the hadal samples of this study fluctuated within less than one order of magnitude with depth [22], these relative abundance patterns resembled absolute abundances estimated by normalizing the data to cell numbers (data not shown).The directional changes in the microbial community composition with sediment depth and associated redox zonation indicated an active community turnover. Similar succession patterns of Gammaproteobacteria, Thaumarchaeota, Planctomycetes and other major microbial taxa have been found in sediments across the entire oceanic depth range from less than a 100 m water depth to the bottom of the Challenger Deep in the Mariana Trench [9, 24, 25]. Our data indicate that some of these directional changes are associated with redox stratification and thus are an inherent characteristic of cohesive marine sediments.Assembly of subsurface phyla in the ferruginous zone in the Atacama TrenchThe high spatial resolution of sampling across redox zones at multiple sites provides new insight into the assembly of deeper microbial communities. For example, combining the relative read abundance of Atribacteria (Fig. 2C) with total cell counts [22] (ignoring potential PCR bias and assuming the same 16S rRNA gene copy numbers in this group as in the community on average), we estimate an absolute increase of Atribacteria from a mean of 5.4 × 103 cells cm−3 at the upper boundary of the ferruginous zone at 6 cm depth to 1.2 × 107 cell cm−3 at 30 cm sediment depth ( >2000-fold increase). Excluding mortality, this can be accomplished in 11 generations, and given an estimated sedimentation rate during periods with no mass depositions of approximately 0.05 cm year−1 (unpublished data) would have occurred over approximately 500 years. Although the number of generations is a minimum estimate, this timeframe appears to leave relatively little opportunity for diversification, as previously concluded for deeper subsurface sediments [11].We further note that bioturbation in hadal sediments is mostly limited to meiofaunal infauna and epibenthic amphipods; hence, sediment mixing is unlikely to affect the depth distribution of microbes below the topmost centimeters [38]. As discussed in previous studies, vertical dispersal of microbes by means of active motility is unlikely to play a role in community assembly in cohesive sediments due to energetic constraints and short-distance chemical gradients [9, 39]. As the sediment in both the Kermadec and in the Atacama trench is cohesive [26], and in accordance with previous studies in deeper redox zones [9, 10], we therefore conclude that selection is likely the dominant force controlling community composition and, e.g., giving Atribacteria their dominant role in the ferruginous zone. They appear to grow from a small seed stock that arrives at the sediment surface and survive burial in an inactive state, until oxygen and nitrate are depleted. Other obligate anaerobes in marine sediments may be subject to similar constraints (see also [10]). This implies that there is little diversification potential for obligate anaerobes in hadal trench sediments. Other anaerobic niches in the hadal zone that might have more diversification potential include hydrothermally active sites and the guts of fauna [40, 41]. However, the conditions in guts and hydrothermally active sediments differ from those of cold deep-sea sediments, and these environments therefore harbor very different microbial communities [42, 43]. Hence, the majority of obligate anaerobes must have originated from the overlying water column and come with the necessary adaptations for the increased hydrostatic pressure and other conditions in hadal sediments. We therefore hypothesize that most obligate anaerobes in hadal sediments tolerate but do not prefer hadal pressures.Frequency distribution of ASVs and taxonomic affiliation of cosmopolitans along the Atacama TrenchDespite the large number of ASVs present in our dataset, only few were found in all samples from a given redox zone, yet these cosmopolitans tended to account for a large fraction of sequencing reads (Supplementary Fig. 7). This was especially pronounced in the rarefied data from the oxic zone where 365 out of 24,844 ASVs occurred across all oxic samples and comprised over 40% of all reads obtained from this zone, thereby contributing substantially to the similarity between cores and sites within the Atacama Trench (Fig. 2). The majority of these ubiquitous reads originated from ASVs belonging Gammaproteobacteria, Thaumarchaeota, Alphaproteobacteria, and Bacteroidetes (Supplementary Fig. 8). The nitrogenous and ferruginous communities were generally more variable, but ubiquitous ASVs accounted for 15% and 10% of all reads, respectively. In both zones most of the reads originated from ASVs classified as Alphaproteobacteria and Bacteroidetes, with ubiquitous ASVs belonging to Ca. Phycisphaerae becoming more abundant in the ferruginous zone.OTUs with high abundances were previously found to be cosmopolitan in deep-sea sediments [6]. Here, we show that this observation does not change when using ASVs and thus a much finer phylogenetic resolution for the formation of ecological units. The decrease of cosmopolitan ASVs in the nitrogenous and ferruginous zones relative to the oxic zone might be due to dispersal barriers in combination with the small seed-stocks of anaerobes in the upper parts of the sediment, which may lead to a higher level of stochasticity in community assembly in deeper sections. In the oxic zone, physical disturbances lead to resuspension of sediment particles and microbes into the water column [44], where they can be transported along trench axes by bottom currents known to ventilate trenches [45]. Therefore, we suggest that dispersal resulted in greater relative read abundances of ubiquitous ASVs in the oxic zone than in the nitrogenous and ferruginous zones.Core microbiomes of each redox zone along the Atacama TrenchTo further analyze the overlaps between abundant community members across redox zones, we performed a core community analysis using the toolset of the ampvis2 R package with adjusted cutoff parameters [32] (see Supplementary Material and Methods). This analysis defines ASVs as part of the core microbiome, when they are above 0.05% relative abundance and within the top 50% of all reads. This classified more than 99% of all ASVs as rare biosphere, while the remaining 441 core ASVs accounted for almost half of all reads (Supplementary Fig. 9A). Each redox zone had a distinct core microbiome, with 196, 66, and 91 core ASVs in the oxic, nitrogenous, and ferruginous zones, respectively, comprising 10.6%, 4.6%, and 8.5% of all reads. The three zones had 17 core ASVs in common that comprised 8.2% of all reads. Aside from these common ASVs, the overlaps between the core microbiomes of the oxic and nitrogenous redox zones were greater than those with the ferruginous zone, with the oxic and nitrogenous sharing an additional 60 ASVs (10.3% of all reads). By contrast, the oxic and nitrogenous zones only shared additional 3 (0.4% of all reads) and 8 (1.5% of all reads) ASVs with the ferruginous zone, respectively. Members of the core microbiome are usually abundant species that are present not merely due to immigration or advection but also through growth, and that are of biogeochemical importance [46, 47]. As the phylum-level composition of the core microbiome in each redox zone was mostly congruent with the overall relative abundance of phyla in each zone (Supplementary Fig. 9B), the distinct shifts in the core microbiome compositions between the zones hint at the potential niche spectra of individual phyla associated with each redox zone (see Supplementary Fig. 9). While many of the core ASVs seemed to thrive in both the oxic and in the nitrogenous zone, the conditions of microbial life seemed to change relatively abruptly when entering the ferruginous zone, resulting in the recruitment of deep-biosphere taxa.Community composition of Archaea along the Atacama TrenchAround 20% of all ASVs in the universal 16S rRNA gene dataset were classified as Archaea. However, due to known mismatches of universal 16S rRNA gene primer sets with archaeal lineages, in particular the phylum Thaumarchaeota, we thus also sequenced archaea-specific 16S rRNA gene amplicons with the same read depth. This primer set recovered approximately three times more thaumarcheotal ASVs than the universal set (4410 vs 1477) and also showed a better coverage over Euryarchaeota (3728 vs 2032), Crenarchaeota (1884 vs 707, including “Ca. Bathyarchaeia”), “Ca. Hydrothermarchaeota” (313 vs 101), and Hadesarchaea (71 vs 26). By contrast, the universal 16S rRNA gene dataset contained 17 times more “Ca. Woesearchaeota” ASVs (42,582 vs 2512), with this phylum even dominating ASV richness over Thaumarchaeota in the universal dataset. Sequencing the HR horizons with this primer set was only successful for only a small number of the samples. As CR horizons were more successful, we focus on this dataset.Thaumarchaeota was the overwhelmingly dominant phylum in the archaeal dataset and drove most of the dissimilarity between individual redox zones (Fig. 2B). In the oxic and nitrogenous zones, they contributed up to 99.3% relative abundance, and other lineages, particularly Crenarchaeota, “Ca. Hydrothermarchaeota,” Euryarchaeota, and “Ca. Asgardaeota,” only increased in relative abundance in the ferruginous zone (Supplementary Fig. 10). Consequently, ordination plots of this dataset only showed a depth gradient to the bottom of the nitrogenous zone, while samples from the ferruginous zone deviated strongly from this gradient (Fig. 2B). Estimates of absolute abundances of the individual archaeal phyla from the universal 16S rRNA gene dataset (Supplementary Fig. 11) showed that the relative depth-wise increase of “Ca. Asgardaeota” (0–9.8%) and Crenarchaeota (0–32.9%) in the archaeal dataset reflects their increase in absolute abundance from 1.4 × 103 to 3.4 × 105 “Ca. Asgardaeota” per ml sediment and 3.7 × 102 to 4.4 × 105 Crenarchaeota per ml sediment. This suggested possible growth of these lineages in hadal sediments.Globally, bacterial lineages dominate over archaeal lineages in marine water columns and surface sediments [48, 49]. The only archaeal phylum in these habitats of comparable abundance is Thaumarchaeota. However, in coastal subsurface sediments, archaeal lineages belonging to “Ca. Lokiarchaeota” and the Miscellaneous Crenarchaeota Group (here referred to as “Ca. Bathyarchaeia”) were found to comprise the majority of intact microbial cells, and Archaea in general contributed significantly to the carbon turnover in these systems [49,50,51,52]. These studies showed that recruitment of archaeal strains occurs in the first few centimeters of these sediments but did not provide a more specific location or connection to biogeochemistry. Our data indicated that the enrichment of Crenarchaeota (including “Ca. Bathyarchaeia”) and “Ca. Asgardaeota” started similarly to that of Atribacteria at the interface of the nitrogenous and ferruginous zones, and was accompanied by an increase in archaeal abundance relative to bacteria. Consequently, both the universal and the archaeal 16S rRNA gene data suggested that the interface between the nitrogenous and ferruginous zones marks the beginning of assembly of subsurface-like microbial communities. This interface also marks the transition from nitrate reduction to iron and/or sulfate reduction as the dominant terminal electron accepting processes. According to existing models, this transition is further associated with a switch in how organic matter is mineralized, with aerobes and denitrifiers being capable of degrading and oxidizing complex organic substrates individually, while a functional division between fermentation and respiration among two sets of organisms is necessary during dissimilatory iron and sulfate reduction [53, 54]. Therefore, we suggest that the distinct differences in microbial communities across the nitrogenous-ferruginous interface are due to the utilization of different electron acceptors and the associated division of labor that causes a rise of fermenters.Community composition across hadal, abyssal, and bathyal sedimentsTo get further insights to factors influencing microbial community composition in hadal sediments, we compared the Atacama Trench to the Kermadec Trench in the less productive western South Pacific, as well as to abyssal and bathyal sites adjacent to these trenches. Kermadec Trench sediments were characterized by deeper oxygen and nitrate penetration depths than in the Atacama Trench (8.5 to >18 cm and 15 to >30 cm, respectively), pushing the ferruginous zone below the sampled sediment horizons at one of the sites (K4 [26], Supplementary Table 1). Consequently, data on the ferruginous zone of the Kermadec Trench were scarce (Supplementary Table 1). In addition, the entire oxic zone was only covered with confidence at site K6, due to potential loss of surface layers at the other stations. Similar to the Atacama Trench, a non-steady state depositional regime in the Kermadec Trench was indicated by fluctuating microbial abundances and visible layering of the sediments.Sediments from the abyssal plains adjacent to both trenches (A7 and K7) showed even deeper oxygen penetration beyond the measured range of the oxygen profiling lander ( >20 cm) and projections indicated that these sediments were oxic across the entire interval analyzed here [26]. In contrast to the trench sites, microbial abundance decayed exponentially with sediment depth and was associated with parallel decreases in TOC content [22]. Conversely, sediment cores from the bathyal (A1) and abyssal (A9) continental slope sites next to the Atacama Trench reached into ferruginous and nitrogenous horizons, respectively, with oxygen penetrating to 1.9 and 6.7 cm, respectively, and nitrate reaching ~6.5 cm at A9. The TOC and microbial abundances showed no clear downcore pattern at these sites [22].At the phylum level, the hadal communities revealed by universal primers were similar in the two trenches, though the drop in Thaumarchaeota abundance was more sharply located at the oxic-nitrogenous interface in the Kermadec Trench than in the Atacama Trench (Fig. 3A and Supplementary Fig. 12A). The Kermadec Trench also exhibited higher relative abundances of “Ca. Woesearchaeota” (18.2% vs 9.2%) in deeper sediment horizons. The archaeal datasets differed more clearly between the two trenches (Supplementary Fig. 12B). While Crenarchaeota reached almost 20% abundance in the Atacama Trench and were detected in the oxic zone, they were essentially absent in the Kermadec Trench. In contrast, “Ca. Diapherotrites” (DPANN) contributed up to 28.6% of relative abundance in the deeper sections of the Kermadec sediments but were almost absent in the Atacama Trench. Similar small-scale differences in the relative abundances of microbial phyla were previously observed between the Japan, Izu-Ogasawara and Mariana trenches [24], Mariana and Mussau trenches [55], as well as between the Mariana and Kermadec trenches [25]. Peoples et al. [25] showed that the Kermadec Trench was enriched in Bacteroidetes, “Ca. Hydrogenedentes” and Planctomycetes in comparison to the Mariana Trench, while the latter had higher relative abundances of “Ca. Marinimicrobia,” Thaumarchaeota, “Ca. Woesearchaeota,” and Chloroflexi. Along with this high similarity between trenches on the phylum level, 58% of all OTUs with ≥97% sequence similarity were shared between the Mariana and Kermadec trenches and these shared OTUs comprised over 95% of all 16S rRNA gene amplicon reads [25]. Thus, they concluded that endemism did not cause the community dissimilarity between the trenches and did not occur on the OTU level. Our ASV-based analysis provided a finer phylogenetic resolution [30], had approximately ten-fold higher sequencing depth, and spanned over three redox zones (up to 40 cm sediment depth) instead of the first 10 cm as in the previous study. A core-to-core comparison from A6 and K6 (the only site with an undisturbed sediment surface in the Kermadec Trench), similar to that of Peoples et al. [25], revealed that the sites shared around 8% of all ASVs, yet these shared ASVs accounted for 62% of all obtained reads from the respective samples. A core microbiome analysis for each redox zone (Fig. 3B) further revealed large overlaps in abundant ASVs between the two trenches. Thus, even with our expansion of the analysis we reach a similar conclusion as Peoples and coworkers that endemism must be relatively rare in hadal sediments. However, the overlap between the trenches decreased from the oxic to the nitrogenous and ferruginous zones. While individual redox zones in these two geographically isolated hadal trenches provide similar ecological niches and are to a large extent inhabited by the same abundant ASVs, this decrease of core microbiome overlaps may be driven by enhanced dispersal barriers, as discussed in the previous paragraph.Fig. 3: Microbial community composition across trenches.A Relative read abundances (%) on phylum/class level of the ten most abundant taxonomic groups (universal 16S rRNA gene data) grouped by individual redox zones in hadal samples of the Kermadec and Atacama trenches. Both the color gradient and the number within the squares indicate the average relative read abundances within the respective sample group. B Number of ASVs and relative fractions of total reads constituted by the core microbiomes of the oxic, nitrogenous, and ferruginous zones, respectively, as unique to the Atacama Trench (blue), unique to the Kermadec Trench (yellow), and shared between trenches (purple).Full size imageUnexpectedly, the phylum-level composition and depth distribution at the continental slope (sites A1 and A9) resembled the results from the hadal zone more closely than those from the abyssal plain, with Thaumarchaeota, Alphaproteobacteria, Gammaproteobacteria, and Bacteroidetes decreasing in relative abundance from the oxic to deeper sediment sections, and Chloroflexi increasing with sediment depth. In contrast, at the abyssal sites, the relative abundance of Thaumarchaeota was almost twice as high (A7: 23% K7: 24%) as in the oxic zone of hadal sediments and did not change significantly over sediment depth (Supplementary Fig. 13A). This was particularly pronounced in the archaea-specific dataset, where Thaumarchaeota comprised more than 99% mean read abundance (Supplementary Fig. 13B). However, in contrast to the hadal and continental slope sediments, where increased relative abundance also indicated growth, estimated absolute abundances at the abyssal plain sites generally decreased with sediment depth (Supplementary Fig. 13C). Hence, we propose that the downcore changes in the abyssal plain sediments may reflect differential persistence and survival capabilities rather than growth, similar to conclusions from subsurface sediments [56]. Thus, the conditions for microbial life differ quite fundamentally between abyssal plains and hadal trenches. The similar directional phylum-level changes in the continental slope and hadal sites further supported biogeochemical forcing as a main driver of microbial community composition on the phylum level and suggested that potential effects associated with oceanic depth or hydrostatic pressure are secondary or mainly apply to finer taxonomic levels.When all sites were compared through PCoA of Bray Curtis dissimilarities, both bathyal and abyssal communities differed from hadal samples (Supplementary Fig. 14A), and ANOSIM confirmed this distinction (p  > 0.001, R = 0.595). The archaea-specific dataset showed similar patterns as the universal dataset except with clearer separation between the two trenches (Supplementary Fig. 14B). This analysis indicated a gradient of microbial community composition with increasing oceanic depth, despite the high similarities of phyla compositions between sediments with similar redox stratifications.Focusing on oxic horizons across sites, and thus removing the strong effect of the redox gradient, the abyssal samples again showed high similarities to each other and to the bathyal site, while the communities of the hadal zone of each trench clustered separately (Supplementary Fig. 14C, D). Still, ANOSIM indicated a clear separation of hadal from shallower samples (p  > 0.001, R = 0.677) with less variation within the individual groups. These patterns were also reflected in the overlaps of core microbiomes of all oxic samples, in which adjacent realms shared more ASVs (Hadal–Abyssal: 37 ASVs; 7% of all reads; Abyssal–Bathyal: 26 ASVs; 3.1%) than the hadal sites and bathyal site (6 ASVs; 0.5%), in addition to the 30 core ASVs (9.2%) found in all realms (Fig. 4B and Supplementary Fig. 15). The relatively small overlap of core ASVs between the hadal and bathyal sites indicates a gradient in community composition of the oxic zone across depth realms with respect to the more abundant members. Factors that could impose such a barrier on core microbiome constituents and benthic microbial communities in general are discussed in the next section.Fig. 4: Microbial community composition across benthic realms.A Principal coordinate analysis (PCoA) of Bray Curtis dissimilarity across samples from the oxic zone (CR sectioning) from the hadal (yellow), abyssal (turquoise), and bathyal (purple) realms in the Kermadec Trench (triangles) and Atacama Trench (circles) based upon the universal 16S rRNA gene dataset. B Number of ASVs and relative fractions of total reads constituted by the core microbiomes of the oxic zones of the hadal (yellow), abyssal (purple), and bathyal (turquoise) realms.Full size imageFactors controlling community composition in hadal vs bathyal and abyssal sedimentsPrevious studies on hadal trench sediments suggested that geochemical factors had a stronger impact on community composition than, for instance, hydrostatic pressure [24]. In this section we aim to test previous hypotheses by determining how well TOC concentration and redox zonation explain variation in the microbial communities. We start by excluding potential confounding effects of oceanic depth and geographic isolation by focusing on the Atacama Trench.Dissimilarity within a trenchIn the Atacama Trench, the ordinations of Bray Curtis dissimilarity already hinted that factors associated with sediment depth and redox zonation were the driving forces of microbial community composition (Fig. 2A, B). Previously, it was shown that variation in TOC concentration (ranging from 0.3 to 1.5%) was one of the best predictors of community variation in abyssal and bathyal sediments [6, 57]. In the Atacama Trench, TOC concentration decreased from the northernmost site A10 (1.44 ± 0.35%; downcore average ± SD) to the southernmost A6 (0.44 ± 0.09%; A6) and this decrease coincided with a decrease in metabolic activity along the trench axis [22, 26]. TOC concentrations also fluctuated with increasing sediment depth at each site, most pronouncedly at sites A2 and A10, where values peaked at around 9 cm depth and varied downcore between 0.3–0.9 and 0.7–2.0%, respectively [22].We delineated the effects of redox zonation from site–site variation and TOC concentration on microbial community composition using variation partitioning on Hellinger-transformed ASV counts (Fig. 5A and Supplementary Fig. 16A). The unique fraction of variation statistically explained by TOC was very low (1%, p = 0.014) yet had a large overlap of 4% with the site-to-site variability. This overlap disappeared completely when A10 was excluded, hinting that much of this trend was driven by high TOC concentrations at this single site. Thus, TOC was a poor predictor of microbial community composition in the Atacama Trench, despite the high downcore fluctuations and the broader concentration range along the trench axis than across all locations of the global dataset on abyssal and bathyal surface sediment of Bienhold et al. [6]. Instead, redox zonation explained the largest unique fraction of variation (24%, p  More