Ecology

1.Hamza MA, Anderson WK. Soil compaction in cropping systems: a review of the nature, causes and possible solutions. Soil and Tillage Res. 2005;82:121–45.Article 

Google Scholar 
2.Etana A, Larsbo M, Keller T, Arvidsson J, Schjønning P, Forkman J, et al. Persistent subsoil compaction and its effects on preferential flow patterns in a loamy till soil. Geoderma. 2013;192:430–6.Article 

Google Scholar 
3.de Andrade Bonetti J, Anghinoni I, de Moraes MT, Fink JR. Resilience of soils with different texture, mineralogy and organic matter under long-term conservation systems. Soil Tillage Res. 2017;174:104–12.Article 

Google Scholar 
4.FAO. Status of the world’s soil resources – main report. Food and agriculture organization of the United Nations and intergovernmental technical panel on soils. 607 (2015). http://www.fao.org/3/i5199e/I5199E.pdf5.Garrigues E, Corson MS, Angers DA, Van Der Werf HMG, Walter C. Development of a soil compaction indicator in life cycle assessment. Int J Life Cycle Assess. 2013;18:1316–24.Article 

Google Scholar 
6.Schäffer B, Stauber M, Mueller TL, Müller R, Schulin R. Soil and macro-pores under uniaxial compression. I. Mechanical stability of repacked soil and deformation of different types of macro-pores. Geoderma. 2008;146:183–91.Article 

Google Scholar 
7.Pagliai M, Marsili A, Servadio P, Vignozzi N, Pellegrini S. Changes in some physical properties of a clay soil in Central Italy following the passage of rubber tracked and wheeled tractors of medium power. Soil Tillage Res. 2003;73:119–29.Article 

Google Scholar 
8.Gysi M, Ott A, Flühler H. Influence of single passes with high wheel load on a structured, unploughed sandy loam soil. Soil Tillage Res. 1999;52:141–51.Article 

Google Scholar 
9.Czyz EA. Effects of traffic on soil aeration, bulk density and growth of spring barley. Soil Tillage Res. 2004;79:153–66.Article 

Google Scholar 
10.Drewry JJ, Paton RJ, Monaghan RM. Soil compaction and recovery cycle on a Southland dairy farm: Implications for soil monitoring. Aust J Soil Res. 2004;42:851–6.Article 

Google Scholar 
11.Keller T, Colombi T, Ruiz S, Manalili MP, Rek J, Stadelmann V, et al. Long-term soil structure observatory for monitoring post-compaction evolution of soil structure. Vadose Zo. J. 2017;16:1–16.
Google Scholar 
12.Wingate-Hill R, Jakobsen BF. Increased mechanisation and soil damage in forests – a review. New Zeal J For Sci. 1982;12:380–93.
Google Scholar 
13.Fierer N, Schimel JP, Holden PA. Variations in microbial community composition through two soil depth profiles. Soil Biol Biochem. 2003;35:167–76.CAS 
Article 

Google Scholar 
14.Jégou D, Brunotte J, Rogasik H, Capowiez Y, Diestel H, Schrader S, et al. Impact of soil compaction on earthworm burrow systems using X-ray computed tomography: preliminary study. Eur J Soil Biol. 2002;38:329–36.Article 

Google Scholar 
15.Larsen T, Schjønning P, Axelsen J. The impact of soil compaction on euedaphic Collembola. Appl Soil Ecol. 2004;26:273–81.Article 

Google Scholar 
16.Rosolem C, Foloni JS, Tiritan C. Root growth and nutrient accumulation in cover crops as affected by soil compaction. Soil Tillage Res. 2002;65:109–15.Article 

Google Scholar 
17.Arvidsson J, Håkansson I. Response of different crops to soil compaction-Short-term effects in Swedish field experiments. Soil Tillage Res. 2014;138:56–63.Article 

Google Scholar 
18.van der Linden AMA, Jeurisson LJJ, Van Veen JA, Schippers G. Turnover of soil microbial biomass as influence by soil compaction. In: Hansen J., Henriksen K. Editors. Nitrogen in Organic Wastes Applied to Soil, London, UK: Academic Press;1989, pp. 25–36.19.Weisskopf P, Reiser R, Rek J, Oberholzer HR. Effect of different compaction impacts and varying subsequent management practices on soil structure, air regime and microbiological parameters. Soil Tillage Res. 2010;111:65–74.Article 

Google Scholar 
20.Li Q, Lee Allen H, Wollum AG. Microbial biomass and bacterial functional diversity in forest soils: effects of organic matter removal, compaction, and vegetation control. Soil Biol Biochem. 2004;36:571–9.CAS 
Article 

Google Scholar 
21.Tan X, Chang SX. Soil compaction and forest litter amendment affect carbon and net nitrogen mineralization in a boreal forest soil. Soil Tillage Res. 2007;93:77–86.Article 

Google Scholar 
22.Dexter AR. Soil physical quality Part I. Theory, effects of soil texture, density, and organic matter, and effects on root growth. Geoderma. 2004;120:201–14.Article 

Google Scholar 
23.Renault P, Sierra J. Modeling oxygen diffusion in aggregated soils: II. Anaerobiosis in topsoil layers. Soil Sci Soc Am J. 1994;58:1017–23.CAS 
Article 

Google Scholar 
24.Schnurr-Pütz S, Bååth E, Guggenberger G, Drake HL, Küsel K. Compaction of forest soil by logging machinery favours occurrence of prokaryotes. FEMS Microbiol Ecol. 2006;58:503–16.PubMed 
Article 
CAS 

Google Scholar 
25.Hartmann M, Niklaus PA, Zimmermann S, Schmutz S, Kremer J, Abarenkov K, et al. Resistance and resilience of the forest soil microbiome to logging-associated compaction. ISME J. 2014;8:226–44.CAS 
PubMed 
Article 

Google Scholar 
26.Marshall VG. Impacts of forest harvesting on biological processes in northern forest soils. For Ecol Manage. 2000;133:43–60.Article 

Google Scholar 
27.Ponder F, Tadros M. Phospholipid fatty acids in forest soil four years after organic matter removal and Soil compaction. Appl Soil Ecol. 2002;19:173–82.Article 

Google Scholar 
28.Frey B, Kremer J, Rüdt A, Sciacca S, Matthies D, Lüscher P. Compaction of forest soils with heavy logging machinery affects soil bacterial community structure. Eur J Soil Biol. 2009;45:312–20.Article 

Google Scholar 
29.Hartmann M, Howes CG, VanInsberghe D, Yu H, Bachar D, Christen R, et al. Significant and persistent impact of timber harvesting on soil microbial communities in Northern coniferous forests. ISME J. 2012;6:2199–218.CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
30.Calonego JC, Raphael JPA, Rigon JPG, de Oliveira Neto L, Rosolem CA. Soil compaction management and soybean yields with cover crops under no-till and occasional chiseling. Eur J Agron. 2017;85:31–37.Article 

Google Scholar 
31.Flisch R, Sinaj S, Charles R, Richner W. Grundlagen für die Düngung im Acker- und Futterbau (GRUDAF). Agrar Schweiz. 2009;16:6–31.
Google Scholar 
32.Suter D, Rosenberg E, Mosimann E, Frick R. Mélanges standard pour la production fourragère. Recherche agronomique suisse. 2017;8:1–16.33.Bürgmann H, Pesaro M, Widmer F, Zeyer J. A strategy for optimizing quality and quantity of DNA extracted from soil. J Microbiol Methods. 2001;45:7–20.PubMed 
Article 
PubMed Central 

Google Scholar 
34.Frey B, Rime T, Phillips M, Stierli B, Hajdas I, Widmer F, et al. Microbial diversity in European alpine permafrost and active layers. FEMS Microbiol Ecol. 2016;92:1–17.Article 
CAS 

Google Scholar 
35.Tedersoo L, Lindahl B. Fungal identification biases in microbiome projects. Environ Microbiol Rep. 2016;8:774–9.PubMed 
Article 
PubMed Central 

Google Scholar 
36.Caporaso JG, Lauber CL, Walters WA, Berg-Lyons D, Lozupone CA, Turnbaugh PJ, et al. Global patterns of 16S rRNA diversity at a depth of millions of sequences per sample. Proc Natl Acad Sci USA. 2011;108:4516–22.CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
37.Henry S, Baudoin E, López-Gutiérrez JC, Martin-Laurent F, Brauman A, Philippot L. Quantification of denitrifying bacteria in soils by nirK gene targeted real-time PCR. J Microbiol Methods. 2004;59:327–35.CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
38.Kandeler E, Deiglmayr K, Tscherko D, Bru D, Philippot L. Abundance of narG, nirS, nirK, and nosZ genes of denitrifying bacteria during primary successions of a glacier foreland. Appl Environ Microbiol. 2006;72:5957–62.CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
39.Rognes T, Flouri T, Nichols B, Quince C, Mahé F. VSEARCH: a versatile open source tool for metagenomics. PeerJ. 2016;4:e2584.40.Martin M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet J. 2011;17:10–12.Article 

Google Scholar 
41.Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012;9:357–9.CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
42.Edgar RC, Flyvbjerg H. Error filtering, pair assembly and error correction for next-generation sequencing reads. Bioinformatics. 2015;31:3476–82.CAS 
PubMed 
Article 

Google Scholar 
43.Edgar R. UNOISE2: improved error-correction for Illumina 16S and ITS amplicon sequencing. 2016. bioRxiv:081257. Available from: https://doi.org/10.1101/081257.44.Edgar R. UCHIME2: improved chimera prediction for amplicon sequencing. 2016. bioRxiv:074252. Available from: https://doi.org/10.1101/074252.45.Bengtsson-Palme J, Hartmann M, Eriksson KM, Pal C, Thorell K, Larsson DG, et al. metaxa2: improved identification and taxonomic classification of small and large subunit rRNA in metagenomic data. Mol Ecol Resour. 2015;15:1403–14.CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
46.Bengtsson-Palme J, Ryberg M, Hartmann M, Branco S, Wang Z, Godhe A, et al. Improved software detection and extraction of ITS1 and ITS2 from ribosomal ITS sequences of fungi and other eukaryotes for analysis of environmental sequencing data. Methods Ecol. Evol. 2013;4:914–9.
Google Scholar 
47.Edgar R. SINTAX: a simple non-Bayesian taxonomy classifier for 16S and ITS sequences. 2016. bioRxiv:074161. Available from: https://doi.org/10.1101/074161.48.Pruesse E, Quast C, Knittel K, Fuchs BM, Ludwig W, Peplies J, et al. SILVA: a comprehensive online resource for quality checked and aligned ribosomal RNA sequence data compatible with ARB. Nucleic Acids Res. 2007;35:7188–96.CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
49.Abarenkov K, Henrik Nilsson R, Larsson KH, Alexander IJ, Eberhardt U, Erland S, et al. The UNITE database for molecular identification of fungi – recent updates and future perspectives. New Phytol. 2010;186:281–5.PubMed 
Article 
PubMed Central 

Google Scholar 
50.R Core Team. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. 2017. https://www.R-project.org/.51. Dinno A. dunn.test: Dunn’s Test of Multiple Comparisons Using Rank Sums. R  package version 1.3.5. 2017. https://CRAN.R-project.org/package=dunn.test.52.Oksanen J, Blanchet FG, Friendly M, Kindt R, Legendre P, McGlinn D, et al. vegan: Community Ecology Package. R package version 2.5-7. 2020. https://CRAN.R-project.org/package=vegan.53.Martiny JB, Martiny AC, Weihe C, Lu Y, Berlemont R, Brodie EL, et al. Microbial legacies alter decomposition in response to simulated global change. ISME J. 2017;11:490–9.PubMed 
Article 

Google Scholar 
54.Hemkemeyer M, Christensen BT, Tebbe CC, Hartmann M. Taxon-specific fungal preference for distinct soil particle size fractions. Eur J Soil Biol. 2019;94:1–9.Article 

Google Scholar 
55.Anderson MJ. A new method for non-parametric multivariate analysis of variance. Austral Ecol. 2001;26:32–46.
Google Scholar 
56.Maxime Hervé. RVAideMemoire: Testing and Plotting Procedures for Biostatistics. 2020. https://CRAN.R-project.org/package=RVAideMemoire.57.Gower JC. Some distance properties of latent root and vector methods used in multivariate analysis. Biometrika. 1996;53:325–38.Article 

Google Scholar 
58.Anderson MJ, Willis TJ. Canonical analysis of principal coordinates: a useful method of constrained ordination for ecology. Ecology. 2003;84:511–25.Article 

Google Scholar 
59.Kindt R, Coe R. Tree diversity analysis. A manual and software for common statistical methods for ecological and biodiversity studies. World Agroforestry Centre, Nairobi, Kenya (2005) ISBN 92-9059-179-X60.Storey JD, Tibshirani R. Statistical significance for genomewide studies. Proc Natl Acad Sci USA. 2003;100:9440–5.CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
61.Storey JD, Bass A, Dabney A, Robinson, D. qvalue: Q-value estimation for false discovery rate control. R package version 2.22.0; 2020. http://qvalue.princeton.edu.62.Letunic I, Bork P. Interactive Tree of Life (iTOL) v4: recent updates and new developments. Nucleic Acids Res. 2019;47:256–9.Article 
CAS 

Google Scholar 
63.Louca S, Parfrey LW, Doebeli M. Faprotax: decoupling function and taxonomy in the global ocean microbiome. Science. 2016;353:1272–7.CAS 
PubMed 
Article 

Google Scholar 
64.Nguyen NH, Song Z, Bates ST, Branco S, Tedersoo L, Menke J, et al. FUNGuild: an open-annotation database for parsing fungal community datasets by ecological guild. Fungal Ecol. 2016;20:241–8.Article 

Google Scholar 
65.Horn R. Stress-strain effects in structured unsaturated soils on coupled mechanical and hydraulic processes. Geoderma. 2003;116:77–88.Article 

Google Scholar 
66.Reiser R, Stadelmann V, Weisskopf P, Grahm L, Keller T. System for quasi-continuous simultaneous measurement of oxygen diffusion rate and redox potential in soil. J Plant Nutr Soil Sci. 2020;183:316–26.CAS 
Article 

Google Scholar 
67.Keller T, Colombi T, Ruiz S, Schymanski SJ, Weisskopf P, Koestel J, et al. Soil structure recovery following compaction – short-term evolution of soil physical properties in a loamy soil. Soil Sci Soc Am J. 2021;18:1–19.
Google Scholar 
68.Yvan C, Stéphane S, Stéphane C, Pierre B, Guy R, Hubert B. Role of earthworms in regenerating soil structure after compaction in reduced tillage systems. Soil Biol Biochem. 2012;55:93–103.CAS 
Article 

Google Scholar 
69.Jabro JD, Allen BL, Rand T, Dangi SR, Campbell JW. Effect of previous crop roots on soil compaction in 2 yr rotations under a no-tillage system. Land. 2021;202:1–10.
Google Scholar 
70.Batey T. Soil compaction and soil management – a review. Soil Use and Manag. 2009;25:335–45.Article 

Google Scholar 
71.Cambi M, Certini G, Neri F, Marchi E. The impact of heavy traffic on forest soils: a review. For Ecol Manag. 2015;338:124–38.Article 

Google Scholar 
72.Hu W, Tabley F, Beare M, Tregurtha C, Gillespie R, Qiu W, et al. Short-term dynamics of soil physical properties as affected by compaction and tillage in a silt loam soil. Vadose Zo J. 2018;17:1–13.
Google Scholar 
73.Manyiwa T, Dikinya O. Impact of tillage types on compaction and physical properties of soils of Sebele farms in Botswana. Soil Environ. 2014;33:124–32.
Google Scholar 
74.Håkansson I, Reeder RC. Subsoil compaction by vehicles with high axle load-extent, persistence and crop response. Soil Tillage Res. 1994;29:277–304.Article 

Google Scholar 
75.Grayston SJ, Wang S, Campbell CD, Edwards AC. Selective influence of plant species on microbial diversity in the rhizosphere. Soil Biol Biochem. 1998;30:369–78.CAS 
Article 

Google Scholar 
76.Degrune F, Theodorakopoulos N, Colinet G, Hiel MP, Bodson B, Taminiau B, et al. Temporal dynamics of soil microbial communities below the seedbed under two contrasting tillage regimes. Front Microbiol. 2017;8:1–15.Article 

Google Scholar 
77.Rousk J, Bååth E, Brookes PC, Lauber CL, Lozupone C, Caporaso JG, et al. Soil bacterial and fungal communities across a pH gradient in an arable soil. ISME J. 2010;4:1340–51.PubMed 
PubMed Central 
Article 

Google Scholar 
78.Bach EM, Baer SG, Meyer CK, Six J. Soil texture affects soil microbial and structural recovery during grassland restoration. Soil Biol Biochem. 2010;42:2182–91.CAS 
Article 

Google Scholar 
79.Evans CA, Coombes PJ, Dunstan RH. Wind, rain and bacteria: the effect of weather on the microbial composition of roof-harvested rainwater. Water Res. 2006;40:37–44.CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
80.He M, Zhang J, Shen L, Xu L, Luo W, Li D, et al. High-throughput sequencing analysis of microbial community diversity in response to indica and japonica bar-transgenic rice paddy soils. PLoS ONE. 2019;14:1–26.
Google Scholar 
81.Gschwend F, Aregger K, Gramlich A, Walter T, Widmer F. Periodic waterlogging consistently shapes agricultural soil microbiomes by promoting specific taxa. Appl Soil Ecol. 2020;155:1–9.Article 

Google Scholar 
82.De Neve S, Hofman G. Influence of soil compaction on carbon and nitrogen mineralization of soil organic matter and crop residues. Biol Fertil Soils. 2000;30:544–9.Article 

Google Scholar 
83.Miransari M, Bahrami HA, Rejali F, Malakouti MJ. Effects of soil compaction and arbuscular mycorrhiza on corn (Zea mays L.) nutrient uptake. Soil Tillage Res. 2009;103:282–90.Article 

Google Scholar 
84.Ruser R, Flessa H, Russow R, Schmidt G, Buegger F, Munch JC. Emission of N2O, N2 and CO2 from soil fertilized with nitrate: effect of compaction, soil moisture and rewetting. Soil Biol Biochem. 2006;38:263–74.CAS 
Article 

Google Scholar 
85.Bao Q, Ju X, Qu Z, Christie P, Lu Y. Response of nitrous oxide and corresponding bacteria to managements in an agricultural soil. Soil Biol Biochem. 2012;76:130–41.CAS 

Google Scholar 
86.Ahemad M, Kibret M. Mechanisms and applications of plant growth promoting rhizobacteria: current perspective. J King Saud Univ Sci. 2014;26:1–20.Article 

Google Scholar 
87.Schwarzott D, Walker C, Schüßler A. Glomus, the largest genus of the arbuscular mycorrhizal fungi (Glomales), is nonmonophyletic. Mol Phylogenet Evol. 2001;21:190–7.CAS 
PubMed 
Article 

Google Scholar 
88.Harman GE, Howell CR, Viterbo A, Chet I, Lorito M. Trichoderma species – opportunistic, avirulent plant symbionts. Nat Rev Microbiol. 2004;2:43–56.CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
89.Waqas M, Khan AL, Hamayun M, Shahzad R, Kang SM, Kim JG, et al. Endophytic fungi promote plant growth and mitigate the adverse effects of stem rot: an example of Penicillium citrinum and Aspergillus terreus. J. Plant Interact. 2015;10:280–7.CAS 
Article 

Google Scholar 
90.Könneke M, Bernhard AE, de la Torre JR, Walker CB, Waterbury JB, Stahl DA. Isolation of an autotrophic ammonia-oxidizing marine archaeon. Nature. 2005;437:543–6.PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
91.Kozlowski TT. Soil compaction and growth of woody plants. Scand J For Res. 1999;14:596–619.Article 

Google Scholar 
92.Haggett KD, Gray PP, Dunn NW. Crystalline cellulose degradation by a strain of Cellulomonas and its mutant derivatives. Eur J Appl Microbiol Biotechnol. 1979;8:183–90.CAS 
Article 

Google Scholar 
93.Rivas R, Trujillo ME, Mateos PF, Martínez-Molina E, Velázquez E. Agromyces ulmi sp. nov., xylanolytic bacterium isolated from Ulmus nigra in Spain. Int J Syst Evol Microbiol. 2004;54:1987–90.CAS 
PubMed 
Article 

Google Scholar 
94.Štursová M, Žifčáková L, Leigh MB, Burgess R, Baldrian P. Cellulose utilization in forest litter and soil: identification of bacterial and fungal decomposers. FEMS Microbiol Ecol. 2012;80:735–46.PubMed 
Article 
CAS 

Google Scholar 
95.Brown AM, Howe DK, Wasala SK, Peetz AB, Zasada IA, Denver DR. Comparative genomics of a plant-parasitic nematode endosymbiont suggest a role in nutritional symbiosis. Genome Biol Evol. 2015;7:2727–46.CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
96.Quandt CA, Beaudet D, Corsaro D, Walochnik J, Michel R, Corradi N, et al. The genome of an intranuclear parasite, Paramicrosporidium saccamoebae, reveals alternative adaptations to obligate intracellular parasitism. Elife. 2017;6:1–19.Article 

Google Scholar 
97.Brussaard L, van Faassen, HG. Effects of compaction on soil biota and soil biological processes. In: Soane BD, van Ouwerkerk C. Soil compaction in crop production. Elsevier; 1994. 11, p. 215–35.98.Shestak CJ, Busse MD. Compaction alters physical but not biological indices of soil health. Soil Sci Soc Am J. 2005;69:236.CAS 
Article 

Google Scholar  More

1.Mendes, R., Garbeva, P. & Raaijmakers, J. M. The rhizosphere microbiome: Significance of plant beneficial, plant pathogenic, and human pathogenic microorganisms. FEMS Microbiol. Rev. 37, 634–663. https://doi.org/10.1111/1574-6976.12028 (2013).CAS 
Article 
PubMed 

Google Scholar 
2.Bhattarai, A., Bhattarai, B. & Pandey, S. Variation of soil microbial population in different soil horizons. J. Microbiol. Exp. 2, 00044. https://doi.org/10.15406/jmen.2015.02.00044 (2015).Article 

Google Scholar 
3.Liu, F. et al. Soil indigenous microbiome and plant genotypes cooperatively modify soybean rhizosphere microbiome assembly. BMC Microbiol. 19, 201. https://doi.org/10.1186/s12866-019-1572-x (2019).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
4.Edwards, J. et al. Structure, variation, and assembly of the root-associated microbiomes of rice. Proc. Natl. Acad. Sci. U.S.A. 112, E911-920. https://doi.org/10.1073/pnas.1414592112 (2015).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
5.Vives-Peris, V., de Ollas, C., Gomez-Cadenas, A. & Perez-Clemente, R. M. Root exudates: From plant to rhizosphere and beyond. Plant Cell Rep. 39, 3–17. https://doi.org/10.1007/s00299-019-02447-5 (2020).CAS 
Article 
PubMed 

Google Scholar 
6.Qu, Q. et al. Rhizosphere microbiome assembly and its impact on plant growth. J. Agric. Food Chem. 68, 5024–5038. https://doi.org/10.1021/acs.jafc.0c00073 (2020).CAS 
Article 
PubMed 

Google Scholar 
7.Schmidt, J. E., Kent, A. D., Brisson, V. L. & Gaudin, A. C. M. Agricultural management and plant selection interactively affect rhizosphere microbial community structure and nitrogen cycling. Microbiome 7, 146. https://doi.org/10.1186/s40168-019-0756-9 (2019).Article 
PubMed 
PubMed Central 

Google Scholar 
8.Cordero, J., de Freitas, J. R. & Germida, J. J. Bacterial microbiome associated with the rhizosphere and root interior of crops in Saskatchewan, Canada. Can. J. Microbiol. 66, 71–85. https://doi.org/10.1139/cjm-2019-0330 (2020).CAS 
Article 
PubMed 

Google Scholar 
9.Bulgarelli, D. et al. Structure and function of the bacterial root microbiota in wild and domesticated barley. Cell Host Microbe 17, 392–403. https://doi.org/10.1016/j.chom.2015.01.011 (2015).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
10.Lundberg, D. S. et al. Defining the core Arabidopsis thaliana root microbiome. Nature 488, 86–90. https://doi.org/10.1038/nature11237 (2012).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
11.Leoni, C. et al. Plant Health and Rhizosphere microbiome: Effects of the bionematicide Aphanocladium album in tomato plants infested by Meloidogyne javanica. Microorganisms https://doi.org/10.3390/microorganisms8121922 (2020).Article 
PubMed 
PubMed Central 

Google Scholar 
12.Vitulo, N. et al. Bark and grape microbiome of vitis vinifera: Influence of geographic patterns and agronomic management on bacterial diversity. Front. Microbiol. 9, 3203. https://doi.org/10.3389/fmicb.2018.03203 (2018).Article 
PubMed 

Google Scholar 
13.Hu, J. et al. Rhizosphere microbiome functional diversity and pathogen invasion resistance build up during plant development. Environ. Microbiol. 22, 5005–5018. https://doi.org/10.1111/1462-2920.15097 (2020).Article 
PubMed 

Google Scholar 
14.Qiao, Q. et al. The variation in the rhizosphere microbiome of cotton with soil type, genotype and developmental stage. Sci. Rep. 7, 3940. https://doi.org/10.1038/s41598-017-04213-7 (2017).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
15.Baudoin, E., Benizri, E. & Guckert, A. Impact of growth stage on the bacterial community structure along maize roots, as determined by metabolic and genetic fingerprinting. Appl. Soil. Ecol. 19, 135–145. https://doi.org/10.1016/S0929-1393(01)00185-8 (2002).Article 

Google Scholar 
16.DeAngelis, K. M. et al. Selective progressive response of soil microbial community to wild oat roots. ISME J. 3, 168–178. https://doi.org/10.1038/ismej.2008.103 (2009).CAS 
Article 
PubMed 

Google Scholar 
17.Ding, L. J. et al. Microbiomes inhabiting rice roots and rhizosphere. FEMS Microbiol. Ecol. https://doi.org/10.1093/femsec/fiz040 (2019).Article 
PubMed 

Google Scholar 
18.Fan, K. et al. Rhizosphere-associated bacterial network structure and spatial distribution differ significantly from bulk soil in wheat crop fields. Soil Biol. Biochem. 113, 275–284. https://doi.org/10.1016/j.soilbio.2017.06.020 (2017).CAS 
Article 

Google Scholar 
19.Jaiswal, S. K., Mohammed, M. & Dakora, F. D. Microbial community structure in the rhizosphere of the orphan legume Kersting’s groundnut [Macrotyloma geocarpum (Harms) Marechal & Baudet]. Mol. Biol. Rep. 46, 4471–4481. https://doi.org/10.1007/s11033-019-04902-8 (2019).CAS 
Article 
PubMed 

Google Scholar 
20.Kuramae, E. E. et al. Soil characteristics more strongly influence soil bacterial communities than land-use type. FEMS Microbiol. Ecol. 79, 12–24. https://doi.org/10.1111/j.1574-6941.2011.01192.x (2012).CAS 
Article 
PubMed 

Google Scholar 
21.Lauber, C. L., Hamady, M., Knight, R. & Fierer, N. Pyrosequencing-based assessment of soil pH as a predictor of soil bacterial community structure at the continental scale. Appl. Environ. Microbiol. 75, 5111–5120. https://doi.org/10.1128/AEM.00335-09 (2009).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
22.Mendes, L. W., Kuramae, E. E., Navarrete, A. A., van Veen, J. A. & Tsai, S. M. Taxonomical and functional microbial community selection in soybean rhizosphere. ISME J. 8, 1577–1587. https://doi.org/10.1038/ismej.2014.17 (2014).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
23.Peiffer, J. A. et al. Diversity and heritability of the maize rhizosphere microbiome under field conditions. Proc. Natl. Acad. Sci. U.S.A. 110, 6548–6553. https://doi.org/10.1073/pnas.1302837110 (2013).ADS 
Article 
PubMed 
PubMed Central 

Google Scholar 
24.Perez-Jaramillo, J. E. et al. Deciphering rhizosphere microbiome assembly of wild and modern common bean (Phaseolus vulgaris) in native and agricultural soils from Colombia. Microbiome 7, 114. https://doi.org/10.1186/s40168-019-0727-1 (2019).Article 
PubMed 
PubMed Central 

Google Scholar 
25.Sugiyama, A., Ueda, Y., Zushi, T., Takase, H. & Yazaki, K. Changes in the bacterial community of soybean rhizospheres during growth in the field. PLoS ONE 9, e100709. https://doi.org/10.1371/journal.pone.0100709 (2014).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
26.Xu, J. et al. The structure and function of the global citrus rhizosphere microbiome. Nat. Commun. 9, 4894. https://doi.org/10.1038/s41467-018-07343-2 (2018).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
27.Haldar, S. & Sengupta, S. Impact of plant development on the rhizobacterial population of Arachis hypogaea: A multifactorial analysis. J. Basic Microbiol. 55, 922–928. https://doi.org/10.1002/jobm.201400683 (2015).CAS 
Article 
PubMed 

Google Scholar 
28.Dai, L. et al. Effect of drought stress and developmental stages on microbial community structure and diversity in peanut rhizosphere soil. Int. J. Mol. Sci. https://doi.org/10.3390/ijms20092265 (2019).Article 
PubMed 
PubMed Central 

Google Scholar 
29.Desmae, H. et al. Genetics, genomics and breeding of groundnut (Arachis hypogaea L.). Plant Breed 138, 425–444. https://doi.org/10.1111/pbr.12645 (2019).Article 
PubMed 

Google Scholar 
30.Pandey, M. K. et al. Translational genomics for achieving higher genetic gains in groundnut. Theor. Appl. Genet. 133, 1679–1702. https://doi.org/10.1007/s00122-020-03592-2 (2020).Article 
PubMed 
PubMed Central 

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

Google Scholar 
32.Parks, D. H. et al. A complete domain-to-species taxonomy for bacteria and archaea. Nat. Biotechnol. 38, 1079–1086. https://doi.org/10.1038/s41587-020-0501-8 (2020).CAS 
Article 
PubMed 

Google Scholar 
33.Parks, D. H. et al. A standardized bacterial taxonomy based on genome phylogeny substantially revises the tree of life. Nat. Biotechnol. 36, 996–1004. https://doi.org/10.1038/nbt.4229 (2018).CAS 
Article 
PubMed 

Google Scholar 
34.Lalucat, J., Mulet, M., Gomila, M. & Garcia-Valdes, E. Genomics in bacterial taxonomy: impact on the genus pseudomonas. Genes https://doi.org/10.3390/genes11020139 (2020).Article 
PubMed 
PubMed Central 

Google Scholar 
35.Correa-Galeote, D., Bedmar, E. J., Fernandez-Gonzalez, A. J., Fernandez-Lopez, M. & Arone, G. J. Bacterial communities in the rhizosphere of Amilaceous Maize (Zea mays L.) as assessed by pyrosequencing. Front. Plant Sci. 7, 1016. https://doi.org/10.3389/fpls.2016.01016 (2016).Article 
PubMed 
PubMed Central 

Google Scholar 
36.Xu, Y. et al. Influence of salt stress on the rhizosphere soil bacterial community structure and growth performance of groundnut (Arachis hypogaea L.). Int. Microbiol. 23, 453–465. https://doi.org/10.1007/s10123-020-00118-0 (2020).CAS 
Article 
PubMed 

Google Scholar 
37.Hamonts, K. et al. Field study reveals core plant microbiota and relative importance of their drivers. Environ. Microbiol. 20, 124–140. https://doi.org/10.1111/1462-2920.14031 (2018).CAS 
Article 
PubMed 

Google Scholar 
38.Ansari, F. A. & Ahmad, I. Isolation, functional characterization and efficacy of biofilm-forming rhizobacteria under abiotic stress conditions. Antonie Van Leeuwenhoek 112, 1827–1839. https://doi.org/10.1007/s10482-019-01306-3 (2019).CAS 
Article 
PubMed 

Google Scholar 
39.Singh, T. B. et al. Identification, characterization and evaluation of multifaceted traits of plant growth promoting rhizobacteria from soil for sustainable approach to agriculture. Curr. Microbiol. 77, 3633–3642. https://doi.org/10.1007/s00284-020-02165-2 (2020).CAS 
Article 
PubMed 

Google Scholar 
40.Govindasamy, V. et al. Multi-trait PGP rhizobacterial endophytes alleviate drought stress in a senescent genotype of sorghum [Sorghum bicolor (L.) Moench]. 3 Biotech 10, 13. https://doi.org/10.1007/s13205-019-2001-4 (2020).Article 
PubMed 

Google Scholar 
41.Abedinzadeh, M., Etesami, H. & Alikhani, H. A. Characterization of rhizosphere and endophytic bacteria from roots of maize (Zea mays L.) plant irrigated with wastewater with biotechnological potential in agriculture. Biotechnol. Rep. 21, e00305. https://doi.org/10.1016/j.btre.2019.e00305 (2019).Article 

Google Scholar 
42.Hashem, A., Tabassum, B. & Allah, F. A. Bacillus subtilis: A plant-growth promoting rhizobacterium that also impacts biotic stress. Saudi J. Biol. Sci. 26, 1291–1297. https://doi.org/10.1016/j.sjbs.2019.05.004 (2019).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
43.Parks, D. H. et al. Recovery of nearly 8,000 metagenome-assembled genomes substantially expands the tree of life. Nat. Microbiol. 2, 1533–1542. https://doi.org/10.1038/s41564-017-0012-7 (2017).CAS 
Article 
PubMed 

Google Scholar 
44.Gomez-Lama Cabanas, C. et al. Indigenous Pseudomonas spp. strains from the olive (Olea europaea L.) rhizosphere as effective biocontrol agents against Verticillium dahliae: From the host roots to the bacterial genomes. Front. Microbiol. 9, 277. https://doi.org/10.3389/fmicb.2018.00277 (2018).Article 
PubMed 
PubMed Central 

Google Scholar 
45.Ansari, F. A. & Ahmad, I. Fluorescent pseudomonas-FAP2 and Bacillus licheniformis interact positively in biofilm mode enhancing plant growth and photosynthetic attributes. Sci. Rep. 9, 4547. https://doi.org/10.1038/s41598-019-40864-4 (2019).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
46.Pandey, K. K., Mayilraj, S. & Chakrabarti, T. Pseudomonas indica sp. nov., a novel butane-utilizing species. Int. J. Syst. Evol. Microbiol. 52, 1559–1567. https://doi.org/10.1099/00207713-52-5-1559 (2002).CAS 
Article 
PubMed 

Google Scholar 
47.Shade, A., Jacques, M. A. & Barret, M. Ecological patterns of seed microbiome diversity, transmission, and assembly. Curr. Opin. Microbiol. 37, 15–22. https://doi.org/10.1016/j.mib.2017.03.010 (2017).Article 
PubMed 

Google Scholar 
48.Adam, E., Bernhart, M., Müller, H., Winkler, J. & Berg, G. The Cucurbita pepo seed microbiome: Genotype-specific composition and implications for breeding. Plant Soil 422, 35–49. https://doi.org/10.1007/s11104-016-3113-9 (2018).CAS 
Article 

Google Scholar 
49.Truyens, S., Weyens, N., Cuypers, A. & Vangronsveld, J. Bacterial seed endophytes: Genera, vertical transmission and interaction with plants. Environ. Microbiol. Rep. 7, 40–50. https://doi.org/10.1111/1758-2229.12181 (2015).Article 

Google Scholar 
50.Kong, H. G., Song, G. C. & Ryu, C.-M. Inheritance of seed and rhizosphere microbial communities through plant–soil feedback and soil memory. Environ. Microbiol. Rep. 11, 479–486. https://doi.org/10.1111/1758-2229.12760 (2019).Article 
PubMed 

Google Scholar 
51.Frindte, K., Pape, R., Werner, K., Loffler, J. & Knief, C. Temperature and soil moisture control microbial community composition in an arctic-alpine ecosystem along elevational and micro-topographic gradients. ISME J. 13, 2031–2043. https://doi.org/10.1038/s41396-019-0409-9 (2019).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
52.Cook, R. J. et al. Molecular mechanisms of defense by rhizobacteria against root disease. Proc. Natl. Acad. Sci. U.S.A. 92, 4197–4201. https://doi.org/10.1073/pnas.92.10.4197 (1995).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
53.Chaparro, J. M., Badri, D. V. & Vivanco, J. M. Rhizosphere microbiome assemblage is affected by plant development. ISME J. 8, 790–803. https://doi.org/10.1038/ismej.2013.196 (2014).CAS 
Article 
PubMed 

Google Scholar 
54.Sharma, S. B., Sayyed, R. Z., Trivedi, M. H. & Gobi, T. A. Phosphate solubilizing microbes: Sustainable approach for managing phosphorus deficiency in agricultural soils. Springerplus 2, 587. https://doi.org/10.1186/2193-1801-2-587 (2013).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
55.Kumar, A., Prakash, A. & Johri, B. N. In Bacteria in Agrobiology: Crop Ecosystems (ed. Maheshwari, D. K.) 37–59 (Springer, 2011).Chapter 

Google Scholar 
56.Sachdev, D., Nema, P., Dhakephalkar, P., Zinjarde, S. & Chopade, B. Assessment of 16S rRNA gene-based phylogenetic diversity and promising plant growth-promoting traits of Acinetobacter community from the rhizosphere of wheat. Microbiol. Res. 165, 627–638. https://doi.org/10.1016/j.micres.2009.12.002 (2010).CAS 
Article 
PubMed 

Google Scholar 
57.Lareen, A., Burton, F. & Schafer, P. Plant root-microbe communication in shaping root microbiomes. Plant Mol. Biol. 90, 575–587. https://doi.org/10.1007/s11103-015-0417-8 (2016).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
58.Klindworth, A. et al. Evaluation of general 16S ribosomal RNA gene PCR primers for classical and next-generation sequencing-based diversity studies. Nucleic Acids Res. 41, e1. https://doi.org/10.1093/nar/gks808 (2013).CAS 
Article 
PubMed 

Google Scholar 
59.Callahan, B. J., Sankaran, K., Fukuyama, J. A., McMurdie, P. J. & Holmes, S. P. Bioconductor workflow for microbiome data analysis: From raw reads to community analyses. F1000Res 5, 1492. https://doi.org/10.12688/f1000research.8986.2 (2016).Article 
PubMed 
PubMed Central 

Google Scholar 
60.R Core Team. R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria, https://www.R-project.org/ (2019).61.Alishum, A. DADA2 formatted 16S rRNA gene sequences for both bacteria & archaea. doi: 10.5281/zenodo.2541239 (2019).62.Callahan, B. Silva taxonomic training data formatted for DADA2 (Silva version 132). doi: 10.5281/zenodo.1172783 (2018).63.Callahan, B. RDP taxonomic training data formatted for DADA2 (RDP trainset 16/release 11.5). doi: 10.5281/zenodo.801828 (2017).64.McMurdie, P. J. & Holmes, S. phyloseq: an R package for reproducible interactive analysis and graphics of microbiome census data. PLoS ONE 8, e61217. https://doi.org/10.1371/journal.pone.0061217 (2013).ADS 
CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
65.Kassambara, A. ggpubr: ‘ggplot2’ Based Publication Ready Plots. R package version 0.2.5. https://CRAN.R-project.org/package=ggpubr (2020).66.Lahti, L. & Shetty, S. Tools for microbiome analysis in R Version 2.1.26. http://microbiome.github.com/microbiome (2017).67.Oksanen, J. et al. vegan: Community Ecology Package. R package version 2.5–6. https://CRAN.R-project.org/package=vegan (2019).68.Martinez Arbizu, P. pairwiseAdonis: Pairwise Multilevel Comparison using Adonis. R package version 0.0.1. (2017).69.Love, M. I., Huber, W. & Anders, S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550. https://doi.org/10.1186/s13059-014-0550-8 (2014).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
70.Wickham, H. ggplot2: Elegant Graphics for Data Analysis (Springer, 2016).Book 

Google Scholar 
71.Martin, C. ggConvexHull: Add a convex hull geom to ggplot2. R package version 0.1.0. http://github.com/cmartin/ggConvexHull (2017).72.Campitelli, E. ggnewscale: Multiple Fill and Colour Scales in ‘ggplot2’. R package version 0.4.1. https://CRAN.R-project.org/package=ggnewscale (2020).73.Slowikowski, K. ggrepel: Automatically Position Non-Overlapping Text Labels with ‘ggplot2’. R package version 0.8.2. https://CRAN.R-project.org/package=ggrepel (2020).74.Dowle, M. & Srinivasan, A. data.table: Extension of `data.frame`. R package version 1.12.8. https://CRAN.R-project.org/package=data.table (2019).75.Ammar, R. randomcoloR: Generate Attractive Random Colors. R package version 1.1.0.1. https://CRAN.R-project.org/package=randomcoloR (2019).76.Wickham, H. & Henry, L. tidyr: Tidy Messy Data. R package version 1.0.0. https://CRAN.R-project.org/package=tidyr (2019).77.Wichmann, H. & Seidel, D. scales: Scale Functions for Visualization. R package version 1.1.0. https://CRAN.R-project.org/package=scales (2019).78.Neuwirth, E. RColorBrewer: ColorBrewer Palettes. R package version 1.1–2. https://CRAN.R-project.org/package=RColorBrewer (2014). More

We retrieved 684 Lokiarchaeota and 31 Heimdallarchaeota near-full-length 16S rRNA sequences from sequence libraries generated from sediment sampled at 27 m water depth in 5 cm intervals between 0 and 40 cm below seafloor (cm.b.s.f) in Aarhus Bay (Supplementary Information). The maximum relative read abundance of Lokiarchaeota was 1.6% at 15–20 cm.b.s.f. and 0.1% for Heimdallarchaeota at 10–15 cm.b.s.f. (Fig. 1). The sequences were grouped into 58 Loki- and 3 Heimdallarchaeota operational taxonomic units (OTUs) using a 98% sequence identity threshold and formed three distinct Lokiarchaeota clades and one monophyletic Heimdallarchaeota cluster (Fig. 1). The primer-free sequencing of RNA extracts enabled us to broadly sample the Asgard archaeal diversity in Aarhus Bay sediments and provided a solid database to design oligonucleotide probes for their visualization.Fig. 1: Phylogenetic analysis and depth distribution of Loki- and Heimdallarchaeota 16S rRNA sequences from Aarhus Bay sediments.A Maximum likelihood phylogeny of Loki- and Heimdallarchaeota operational taxonomic units (OTUs) and related sequences selected from the SILVA database (v. 132). Specificities of FISH probes and the number of sequences constituting each OTU are also depicted. TACK archaea were selected as outgroup. Bar: 0.1 substitutions per nucleotide position. B Heatmap and relative abundances of Loki- and Heimdallarchaeota sequences at different sediment depths.Full size imageBased on the newly retrieved full-length sequences, we designed four novel oligonucleotide probes specifically targeting Loki- and Heimdallarchaeota 16S rRNA with high coverage (Fig. 1, Supplementary Table 1). Probe LOK1183 targets almost all sequences in Lokiarchaeota Clade A, which contains 92% of the retrieved Lokiarchaeota sequences from Aarhus Bay sediments, while probe LOK1378 targets 85% of the sequences in all three Lokiarchaeota clades. Probe HEIM329 and HEIM529 each target >97% of the retrieved Heimdallarchaeota sequences. All designed probes cover >89% sequences in their target groups in the SILVA database (v. 132). The two Lokiarchaeota probes match 5 and 10 different non-target sequences in the SILVA database (v. 132), respectively, while the Heimdallarchaeota probes have no match outside their target group. The broad coverage and high specificity suggest that our probes can also be used to detect Loki- and Heimdallarchaeota in other habitats. Furthermore, designing two probes for each phylum enabled us to identify Lokiarchaeota clade A and Heimdallarchaeota cells via double hybridizations with two distinct dyes and thus confidently distinguish true- and false-positive signals (Supplementary Fig. 1). The general archaeal probe ARC915 also targets Lokiarchaeota and thereby provided yet another control for specific hybridization of the two Lokiarchaeota-specific probes, while the non-sense probe NON338 served as the negative control. We also designed competitor probes to minimize the theoretical false-positive hybridizations with the most frequent one and two mismatches [11] in the SILVA database (v. 132) and helper probes to facilitate probe binding [12]. This comprehensive experimental design with appropriate controls enabled reliable detection of low-abundant Loki- and Heimdallarchaeota cells in Aarhus Bay sediments.We used both confocal laser scanning microscopy (CLSM) and three-dimensional super-resolution structured illumination (SR-SIM) microscopy for detailed imaging of dual-labeled Loki- and Heimdallarchaeota signals. Loki- and Heimdallarchaeota cells featured coccoid shapes and often formed clusters (Fig. 2) (Supplementary Fig. 2). Based on SR-SIM imaging, Lokiarchaeota cells (n = 18) were 1.27 ± 0.24 µm in diameter and 1.43 ± 0.25 µm in length, while the width and the length of Heimdallarchaeota cells (n = 11) were 1.30 ± 0.20 µm and 1.37 ± 0.21 µm, respectively (Supplementary Table 2). In addition, we observed a few large ( >3 µm) ovoid and filamentous cells, resembling some of the Lokiarchaeota morphotypes reported from lake sediment [9]; however, we never detected these cell types in double hybridizations with two probes (Supplementary Fig. 1P–R), and therefore consider them false-positives.Fig. 2: Visualization of Loki- and Heimdallarchaeota cells in Aarhus Bay sediments by CARD-FISH.Probe names and the dyes are indicated for each panel. Representative cell morphotypes were imaged in a super-resolution structured illumination microscope (SR-SIM; panels (A), (B), (D), (E)) and confocal laser scanning microscope (CLSM; panels (C) and (F)). For SR-SIM images, single slices from the center of the focal plane are shown. For CLSM images, three-dimensional (3D) surface reconstructions are depicted. All z-stack images taken in CLSM are included in Supplementary Fig. 2. 360° rotation of 3D reconstructed images are also provided in Supplementary Video. Negative and positive controls are shown in Supplementary Fig. 1 together with large ovoid and filamentous false-positive signals. Images are representative of dual labeled Lokiarchaeota (n = 72) and Heimdallarchaeota (n = 70) cells in five individual experiments using two different sediment cores taken from the same sampling site. The scale bar is 1 µm.Full size imageThe DNA stain (4′,6-diamidino-2-phenylindole; DAPI) in the FISH-identified Loki- and Heimdallarchaeota cells was consistently confined to a single spherical central or lateral position (Fig. 2), corroborating the signal pattern suggested for some of the Asgard archaeal cells in lake sediments [9]. Using SR-SIM, we could image a clear gap, which separated the DNA from the ribosome-originated FISH signals with an average width of 0.18 ± 0.07 µm in Heimdallarchaeota and 0.16 ± 0.13 µm in Lokiarchaeota cells (Supplementary Table 2). The spatial separation of DNA and ribosomes in Loki- and Heimdallarchaeota cells represents an unusual observation since DAPI and FISH signals generally overlap partially or completely in prokaryotic cells [13]. Accordingly, SR-SIM imaging of benthic bacteria in Aarhus Bay sediments demonstrated the prevalence of this overlapping signal pattern (Supplementary Fig. 3). Also, the separated DNA signal observed in Loki- and Heimdallarchaeota cells appeared different from the condensed DNA formation previously described, for example, in Escherichia coli cells [14] and the Thaumarcheota Cenarcheum symbiosum [15] and Nitrosopumilus maritimus [16]. To corroborate this, we performed SR-SIM imaging of CARD-FISH-labeled E.coli and N. maritimus cells. Although their DNA was condensed in particular cellular locations, their FISH and DAPI signals always overlapped, indicating that their DNA and ribosomes are partially co-localized and not fully separated (Supplementary Fig. 4).To assess whether the gap between DAPI and FISH signals was indicative of an internal membrane, we tried various dyes to stain membranes of the CARD-FISH-labeled Asgard archaeal cells (Supplementary Information). However, none of these stainings was successful, not even for the outer cell membrane, most likely because cell membranes were disintegrated during the CARD-FISH protocol. We then used wheat germ agglutinin (WGA), a lectin primarily binding to N-acetyl-D-glucosamine but also other glycoconjugates and oligosaccharides [17] to at least be able to visualize the surfaces of Loki- and Heimdallarchaeota cells. WGA consistently decorated a cell surface that enclosed the proximal FISH and DAPI signals, suggesting that both signals originated from the same single cell (Supplementary Fig. 5). The WGA staining also demonstrated extracellular structures connected to some Heimdallarchaeota cells (Supplementary Fig. 5). These structures appear different than the membrane protrusions in the first cultured Lokiarchaeon “Ca. P. syntrophicum”, which has a considerably smaller cell size (550 nm in diameter) and does not possess the separated DNA and ribosome signals [5]. Our observations therefore indicate diverse cellular organizations and morphotypes within Asgard archaea superphylum.Our combined results suggest that genomic material is condensed and spatially distinct from the riboplasm within the detected Loki- and Heimdallarchaeota cells. Considering the anticipated role of Asgard archaea in eukaryogenesis, in particular the presence of ESPs potentially involved in dynamic cytoskeleton formation [18] and membrane remodeling [4, 19], the separation of DNA- and ribosome-derived signals might be indicative of cellular compartmentalization. Alternatively, the observed pattern could be the result of a membrane invagination to form a nucleoid region, similar to membrane organizations for example in Planctomycetes cells [20] or Atribacter laminatus [21].Our study demonstrates the first visualization of diverse Loki- and Heimdallarchaeota cells in the marine environment and provides a protocol for reliable in situ imaging of rare microorganisms in environmental samples. Future research should address the ultrastructure of Asgard archaeal cells using electron microscopy. This would help to elucidate the cell biology of Asgard archaea and provide insights into the emergence of subcellular complexity of the eukaryotic cell. More

The following experimental workflow was implemented to address the main questions raised in our study (Fig. 1).Fig. 1: The scheme of experimental pipeline used in this study to examine the impact of lytic phage infection on the P. aeruginosa population and the development of phage-resistance.Experiments were conducted as follows: culture preparation (1); biofilm formation (2); phage infection with single or cocktail preparations (3); incubation (4); biofilm and planktonic populations sampling (5); culture plating on TSA agar and isolation of discrete colonies (6); phage typing determination (7); to select isolates with unique patterns (8) for further phenotypic (9) and genome sequencing analyses (10).Full size imageThe P. aeruginosa PAO1 reference strain and four other clinical representatives were infected with distinct lytic phages in a single or different cocktail combination. Randomly picked colonies from the surviving cultures were then tested in terms of susceptibility to inoculated phages as well as to the others from the Pseudomonas phages panel (Table 1). We were interested in exploring the broadest clonal variability developed in phage infected Pseudomonas population. Therefore, the first phase of the study was focused on examining the phenotypic heterogeneity of PAO1 reference mutants (phage typing) within planktonic and biofilm populations. Since the consequences of introducing lytic phages into the bacterial population are difficult to predict, a representative pool of bacterial clones that have survived infection was sampled. A total of 780 P. aeruginosa PAO1 clones were typed with phages (planktonic (320), biofilm populations (400) and 60 control clones). No resistance to phages was observed among the control clones taken from untreated biofilm or plankton. Therefore, three biofilm and three planktonic representatives and the wild-type PAO1 were selected for further genetic and fitness analyses (Table S1). Finally, a pool of 95 isolates has been filtered, representing seventeen different phage susceptibility patterns (Tables S1, 2). This selection was based on the maximum variety of phage-type profiles, without accounting for the origin of the isolate (biofilm/plankton), as the infected planktonic bacteria turned out to be less diverse and all phage types were also present in the biofilm population.Since we did not aim to analyse the differences of planktonic versus sessile cells response to phage infection but rather look for maximum population heterogeneity, we decided to focus the investigation on the biofilm population for the other clinical strains during the second stage of this research. Accordingly, 880 (30 clones from every condition plus 10 control clones for each strain) isolated colonies from A5803, AA43, CHA, and PA biofilm populations were first subjected to phage typing. No phage resistance was observed among clones taken from phage-untreated samples compared to the wild-type strain. Ultimately, 35 phage-treated colonies, three controls, and the wild-type from each strain were selected for further investigation, resulting in a pool of 156 clones in total (39 × 4 strains) representing over twenty different phage susceptibility patterns (Table S1 and S3).Do phages always select for cross-resistance to other phages recognising the same bacterial receptor?The application of monovalent phage against reference PAO1 population generally led to the selection of cross-resistance against phages that recognise the same receptor as the applied one (Table S2). This was observed for 12/17 and 23/24 PAO1 clones isolated after LPS- and T4P-dependent phages treatment, respectively. Similar relation (15/20) was only observed for other clinical cultures infected with phiKZ phage (T4P-dependent) (Table S3). The resistance to both groups of phages was less frequent in monovalent infections (14.5% in PAO1 and 32.5% for other clinical strains) compared to polyvalent infections (61.1%; 33/54) and 51.6% (31/60) for PAO1 and clinical strains, respectively. The use of a cocktail of two phages recognising LPS selected for PAO1 clones resistant only to LPS-dependent phages. In contrast, LPS-dependent phages application was mostly accompanied by the emergence of resistance to phages recognising alternative receptors in clinical strains (28/60 cases).The introduction of a particular phage into the population did not guarantee the isolation of clones resistant to this phage. This event was recorded in the case of single phages, as well as for polyvalent combinations (23 PAO1 mutants). However, the cross-resistance to other phages recognising the same or both receptors did also occur. Interestingly, LUZ7 and KTN6 phages could still infect surviving clinical populations with a frequency of 23/60 and 44/80, respectively. Indicating that the resistance to LPS-dependent phages in clinical strains was more difficult to develop compared to those impaired by giant viruses, with 11/60 and 1/20 still sensitive to phiKZ and PA5oct phages, respectively. Almost all PAO1 (80/95) and clinical (127/140) clones treated with phages developed resistance to phage PA5oct, whereas the resistance to the entire phage panel emerged regardless of the single or cocktails application.To conclude, the selection of cross-resistance to other phages recognising the same bacterial receptor was mostly valid in the PAO1 model, whereas the other clinical strains primarily developed the cross-resistance to T4P-dependent phages.Do phages from different taxonomy groups recognising the same receptor cause the emergence of the same type of resistant mutants? Are the defence response and genome changes correlated with the receptor specificity of infecting phage?To assess the genetic basis of the resistance selected by phages, we performed single nucleotide polymorphisms (SNPs) and mapping analyses of 102 reference PAO1 clones and 156 clones derived from clinical strains (Figs. 2, 3, Table S2–S4). The wild-type P. aeruginosa strains were also re-sequenced with Illumina and PacBio technologies to ascertain their complete genomic background. Missense, nonsense, and frameshift mutation variants were taken into account in the analyses. Mutations that also occurred in control isolates were excluded from further consideration. The remaining mutations were divided into six groups: LPS-related genes, mucoidity-associated genes (EPS production, biofilm formation), T4P-related genes, global regulatory genes, and others (hypothetical or undefined function genes). The comparative analysis showed the presence of point mutations in 64 out of 95 examined PAO1 mutants. The frequency of mutations in PAO1 clones isolated after treatment with single or multiple phages was similar (73% vs 61%, respectively). In most of those isolates, only one gene mutation event was recorded (43%). However, in 23 cases SNPs occurred in two or three genes belonging to different metabolic gene groups. Five PAO1 isolates showed the presence of mutations in two genes from one gene group. The 33 cases of SNPs related to LPS synthesis were found in 29 mutants selected with single LPS-dependent phage preparations or in polyvalent combinations. Among these, the most frequent mutation (21/33 cases) was observed within the wzy gene, encoding the B-band O-antigen polymerase [30]. These frequent mutations in the LPS-biosynthesis cluster confirmed the phage resistance results emerging after LUZ7, KT28, and KTN6 phages propagation. In some cases, the LPS gene modification was accompanied by changes in EPS-related genes, leading to a mucoid phenotype. The T4P-dependent phage treatment also led to the selection of specific mutations in genes responsible for T4P expression, but also alterations in flagella-related genes (flgH, fliN, fliP, flhA). The mutations in global regulatory genes (most frequent yqjG and vfr) and “others” gene groups did not show any correlation to the type of phages used.Fig. 2: Graphical presentation of genetic changes occurring in the population of P. aeruginosa as a result of the infection by selected phages.The colour dots refer to particular gene groups where the point mutations (accumulated results) were recorded within the genomes of examined mutant clones. The lower line contains information on the maximum and minimum size of large deletions (grey bands) and the presence of intact prophages (light blue bands). * means mutation in promoter region of the gene.Full size imageFig. 3: The frequency of genetic changes per clone detected in P. aeruginosa strains.Panel (A) represents the PAO1 clones, and panel (B) represents the clinical strains populations. Populations were selected by specific phages targeting LPS (red dots) or T4P (blue dots) as a single treatment or in cocktails. The colour bars refer to particular gene groups where the point mutations were recorded within the genomes of examined mutant clones. N means the number of analysed clones for each strain.Full size imageApart from point mutations, 23% of phage-resistant PAO1 isolates contained large genomic deletions (23,983 bp–544,729 bp) appeared regardless of the phage-type and cocktail composition used as selective pressure agents. All deletions were located in the same region and despite different starting/ending points, they hold a core element of 19,038 bp. This core element carries the galU gene (responsible for LPS biosynthesis), as well as the hmgA gene, which causes the accumulation of brown pigment in bacterial cells when absent. Besides, the cumulated deletion range contained a total of 706,374 bp, including many key genes involved in the bacterial metabolism.Mutations detected in other clinical phage-resistant clones were classified according to the same criteria as in PAO1 (Figs. 2, 3, Table S3, S4). The genome-driven response to phage infection of A5803 was primarily located in global (71%, cpdA) and other genes (34%, PA2911); of AA43 in other genes (31%, PA2911); of CHA in T4P (34%) and global genes (34%, morA); and of PAK in T4P (25%) and other genes (23%, PA2911). Most of the mutations selected by LPS-dependent phage exposition were found in the global regulatory genes (9–11–25–54%) or “other” genes (17–23–31%), rather than in the LPS biosynthesis locus (0–3–6–17%) depending on the impacted strain (Table S3). That confirmed the phage-typing results where LUZ7 and KTN6 phages remained lytic towards surviving clones. In contrast, the application of phiKZ selected for the cross-resistance to T4P-dependent phages as well as for the genetic modifications in pili-related genes. Mutations in global regulatory and “others” genes show no correlation to the receptor specificity of phages used. Interestingly, a portion of phenotypically phage-resistant clones in each clinical P. aeruginosa population (5-9/35 clones) did not reveal any distinguishable genetic modifications. Consistent with PAO1, large genomic deletions were observed in A5803, AA43, and PAK strains ranging between 92,207 bp and 383,693 bp in size, encompassing the galU region. The MEME analysis of the regions flanking the deletions did not reveal specific motifs that would indicate recombination events. Interestingly, the unique large deletion found in CHA strain (15,126 bp) turned out to be the induced ssDNA filamentous Pf1-like phage.Summarising the analyses one might say that phages from different taxonomy groups recognising the same receptor generally cause the emergence of a similar type of resistance within a particular strain. However, the defence response and genome changes correlated with the receptor specificity of infecting phage differ in a strain-dependent manner.Do different strains of P. aeruginosa react similarly to a specific phage infection?The next step aimed to assess the impact of gaining phage resistance in terms of population growth efficiency as an indicator for bacterial fitness. Three of the examined wild-type strains (PAO1, A5803, and CHA) have a naturally rapid growth rate, while the other two (AA43 and PAK) display moderate growth rates. For this reason, the final results are expressed as the cumulated OD600 (Fig. 4, Table S2, S3). Overall, the majority of PAO1 mutants (61/95; 64%, p  0.001) for the clones resistant to 6–7 phages but only in the PAO1 group. Moreover, only the selection done by phage cocktails gave a statistically significant reduction of bacterial growth (p  > 0.001), while no differences were observed regarding groups treated with single LPS- or T4P- dependent phages. In contrast to the PAO1 reference strain, the statistical analyses conducted in the A5803, AA43, CHA, and PAK strains did not show any differences in terms of phage-typing profile nor phage-type selection pressure versus the population fitness reduction (growth rate).Fig. 4: The population growth efficiency as an indicator for bacterial fitness expressed as the cumulative OD600 values of 18 h culture at 37 °C measured at 20-minute intervals.Dots represent the growth of particular clones: the wild-type and control clones (green dots); mutants selected by LPS-dependent phage (red dots); mutants selected by T4P-dependent phage (blue dots); mutants selected by LPS/T4P-dependent PA5oct phage (orange dots); mutants selected by phage cocktail (black dots). * statistically different cumulative OD value compared to phage-untreated pool (p  More

1.Gardner, M. R. & Ashby, W. R. Connectance of large dynamic (cybernetic) systems: Critical values for stability. Nature 228, 784 (1970).ADS 
CAS 
PubMed 
Article 

Google Scholar 
2.May, R. M. Will a large complex system be stable?. Nature 238, 413–414 (1972).ADS 
CAS 
PubMed 
Article 

Google Scholar 
3.McCann, K. S. The diversity-stability debate. Nature 405, 228–233 (2000).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
4.Dunne, J. A. The network structure of food webs. in Ecological Networks: Linking Structure to Dynamics in Food Webs 27–86 (2006).5.Williams, R. J., Brose, U. & Martinez, N. D. Homage to Yodzis and Innes 1992: Scaling up feeding-based population dynamics to complex ecological networks. in From Energetics to Ecosystems: The Dynamics and Structure of Ecological Systems. 37–51 (Springer, 2007).6.Gross, T., Rudolf, L., Levin, S. A. & Dieckmann, U. Generalized models reveal stabilizing factors in food webs. Science 325, 747–750 (2009).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
7.Fahimipour, A. K., Anderson, K. E. & Williams, R. J. Compensation masks trophic cascades in complex food webs. Theor. Ecol. 10, 245–253 (2017).Article 

Google Scholar 
8.Rooney, N. & McCann, K. S. Integrating food web diversity, structure and stability. Trends Ecol. Evolut. 27, 40–46 (2012).Article 

Google Scholar 
9.Jacquet, C. et al. No complexity-stability relationship in empirical ecosystems. Nat. Commun. 7, 12573 (2016).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
10.Brose, U., Williams, R. J. & Martinez, N. D. Allometric scaling enhances stability in complex food webs. Ecol. Lett. 9, 1228–1236 (2006).PubMed 
Article 
PubMed Central 

Google Scholar 
11.Martinez, N. D. Allometric trophic networks from individuals to socio-ecosystems: Consumer-resource theory of the ecological elephant in the room. Front. Ecol. Evolut. 8, 92 (2020).Article 

Google Scholar 
12.Segel, L. A. & Levin, S. A. Application of nonlinear stability theory to the study of the effects of diffusion on predator-prey interactions. in AIP Conference Proceedings, Vol. 27, 123–152 (American Institute of Physics, 1976).13.Durrett, R. & Levin, S. The importance of being discrete (and spatial). Theor. Popul. Biol. 46, 363–394 (1994).MATH 
Article 

Google Scholar 
14.McCann, K. S., Rasmussen, J. & Umbanhowar, J. The dynamics of spatially coupled food webs. Ecol. Lett. 8, 513–523 (2005).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
15.Fahimipour, A. K. & Hein, A. M. The dynamics of assembling food webs. Ecol. Lett. 17, 606–613 (2014).PubMed 
Article 
PubMed Central 

Google Scholar 
16.Brechtel, A., Gramlich, P., Ritterskamp, D., Drossel, B. & Gross, T. Master stability functions reveal diffusion-driven pattern formation in networks. Phys. Rev. E 97, 032307 (2018).ADS 
CAS 
PubMed 
Article 

Google Scholar 
17.Brechtel, A., Gross, T. & Drossel, B. Far-ranging generalist top predators enhance the stability of meta-foodwebs. Sci. Rep. 9, 1–15 (2019).CAS 
Article 

Google Scholar 
18.Gross, T. & et. al. Modern models of trophic meta-communities. Phil. Trans. R. Soc. B (in press).19.Rooney, N., McCann, K., Gellner, G. & Moore, J. C. Structural asymmetry and the stability of diverse food webs. Nature 442, 265–269 (2006).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
20.Otto, S. B., Rall, B. C. & Brose, U. Allometric degree distributions facilitate food-web stability. Nature 450, 1226–1229 (2007).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
21.Williams, R. J. & Martinez, N. D. Simple rules yield complex food webs. Nature 404, 180–183 (2000).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
22.Cohen, J. E., Pimm, S. L., Yodzis, P. & Saldaña, J. Body sizes of animal predators and animal prey in food webs. J. Anim. Ecol. 67–78 (1993).23.Petchey, O. L., Beckerman, A. P., Riede, J. O. & Warren, P. H. Size, foraging, and food web structure. Proc. Natl. Acad. Sci. 105, 4191–4196 (2008).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
24.Cohen, J. E., Briand, F. & Newman, C. M. Community Food Webs: Data and Theory Vol. 20 (Springer, 2012).MATH 

Google Scholar 
25.Elton, C. S. Animal Ecology (University of Chicago Press, 2001).
Google Scholar 
26.Lindeman, R. L. The trophic-dynamic aspect of ecology. Ecology 23, 399–417 (1942).Article 

Google Scholar 
27.Peters, R. H. & Peters, R. H. The Ecological Implications of Body Size Vol. 2 (Cambridge University Press, 1986).
Google Scholar 
28.Riede, J. O. et al. Stepping in Elton’s footprints: A general scaling model for body masses and trophic levels across ecosystems. Ecol. Lett. 14, 169–178 (2011).29.Kalinkat, G. et al. Body masses, functional responses and predator-prey stability. Ecology letters 16, 1126–1134 (2013).PubMed 
Article 
PubMed Central 

Google Scholar 
30.Costa-Pereira, R., Araújo, M. S., Olivier, R. d. S., Souza, F. L. & Rudolf, V. H. Prey limitation drives variation in allometric scaling of predator-prey interactions. Am. Nat. 192, E139–E149 (2018).31.Guzman, L. M. & Srivastava, D. S. Prey body mass and richness underlie the persistence of a top predator. Proc. R. Soc. B 286, 20190622 (2019).PubMed 
PubMed Central 
Article 

Google Scholar 
32.Brose, U. et al. Consumer-resource body-size relationships in natural food webs. Ecology 87, 2411–2417 (2006).PubMed 
Article 
PubMed Central 

Google Scholar 
33.Barnes, C., Maxwell, D., Reuman, D. C. & Jennings, S. Global patterns in predator-prey size relationships reveal size dependency of trophic transfer efficiency. Ecology 91, 222–232 (2010).PubMed 
Article 
PubMed Central 

Google Scholar 
34.Potapov, A. M., Brose, U., Scheu, S. & Tiunov, A. V. Trophic position of consumers and size structure of food webs across aquatic and terrestrial ecosystems. Am. Nat. 194, 823–839 (2019).PubMed 
Article 
PubMed Central 

Google Scholar 
35.MacArthur, R. H. & Wilson, E. O. The Theory of Island Biogeography Vol. 1 (Princeton University Press, 2001).Book 

Google Scholar 
36.Simberloff, D. S. & Wilson, E. O. Experimental zoogeography of islands: the colonization of empty islands. Ecology 50, 278–296 (1969).Article 

Google Scholar 
37.Brown, J. H. & Kodric-Brown, A. Turnover rates in insular biogeography: effect of immigration on extinction. Ecology 58, 445–449 (1977).Article 

Google Scholar 
38.Levins, R. Some demographic and genetic consequences of environmental heterogeneity for biological control. Am. Entomol. 15, 237–240 (1969).
Google Scholar 
39.Gotelli, N. J. Metapopulation models: The rescue effect, the propagule rain, and the core-satellite hypothesis. Am. Nat. 138, 768–776 (1991).Article 

Google Scholar 
40.Crowley, P. H. Dispersal and the stability of predator-prey interactions. Am. Nat. 118, 673–701 (1981).MathSciNet 
Article 

Google Scholar 
41.Reeve, J. D. Environmental variability, migration, and persistence in host-parasitoid systems. Am. Nat. 132, 810–836 (1988).Article 

Google Scholar 
42.Murdoch, W. W. Population regulation in theory and practice. Ecology 75, 271–287 (1994).Article 

Google Scholar 
43.Briggs, C. J. & Hoopes, M. F. Stabilizing effects in spatial parasitoid-host and predator-prey models: A review. Theor. Popul. Biol. 65, 299–315 (2004).PubMed 
MATH 
Article 
PubMed Central 

Google Scholar 
44.Gravel, D., Massol, F. & Leibold, M. A. Stability and complexity in model meta-ecosystems. Nat. Commun. 7, 1–8 (2016).Article 
CAS 

Google Scholar 
45.Mougi, A. & Kondoh, M. Food-web complexity, meta-community complexity and community stability. Sci. Rep. 6, 24478 (2016).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
46.Domenici, P. The scaling of locomotor performance in predator-prey encounters: from fish to killer whales. Comp. Biochem. Physiol. Part A Mol. Integr. Physiol. 131, 169–182 (2001).CAS 
Article 

Google Scholar 
47.Hirt, M. R., Lauermann, T., Brose, U., Noldus, L. P. & Dell, A. I. The little things that run: a general scaling of invertebrate exploratory speed with body mass. Ecology 98, 2751–2757 (2017).PubMed 
Article 
PubMed Central 

Google Scholar 
48.Hirt, M. R., Jetz, W., Rall, B. C. & Brose, U. A general scaling law reveals why the largest animals are not the fastest. Nat. Ecol. Evolut. 1, 1116–1122 (2017).Article 

Google Scholar 
49.Cloyed, C. S., Grady, J. M., Savage, V. M., Uyeda, J. C. & Dell, A. I. The allometry of locomotion. Ecology e03369 (2021).50.Reiss, M. Scaling of home range size: Body size, metabolic needs and ecology. Trends Ecol. Evolut. 3, 85–86 (1988).CAS 
Article 

Google Scholar 
51.Minns, C. K. Allometry of home range size in lake and river fishes. Can. J. Fish. Aquat. Sci. 52, 1499–1508 (1995).Article 

Google Scholar 
52.Jetz, W., Carbone, C., Fulford, J. & Brown, J. H. The scaling of animal space use. Science 306, 266–268 (2004).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
53.Hendriks, A. J., Willers, B. J., Lenders, H. R. & Leuven, R. S. Towards a coherent allometric framework for individual home ranges, key population patches and geographic ranges. Ecography 32, 929–942 (2009).Article 

Google Scholar 
54.Hein, A. M., Hou, C. & Gillooly, J. F. Energetic and biomechanical constraints on animal migration distance. Ecol. Lett. 15, 104–110 (2012).PubMed 
Article 
PubMed Central 

Google Scholar 
55.Hartfelder, J. et al. The allometry of movement predicts the connectivity of communities. Proc. Natl. Acad. Sci. 117, 22274–22280 (2020).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
56.Vander Zanden, M. J. & Vadeboncoeur, Y. Fishes as integrators of benthic and pelagic food webs in lakes. Ecology 83, 2152–2161 (2002).Article 

Google Scholar 
57.Wolkovich, E. M. et al. Linking the green and brown worlds: The prevalence and effect of multichannel feeding in food webs. Ecology 95, 3376–3386 (2014).Article 

Google Scholar 
58.Lomolino, M. V. Immigrant selection, predation, and the distributions of Microtus pennsylvanicus and Blarina brevicauda on islands. Am. Nat. 123, 468–483 (1984).Article 

Google Scholar 
59.Beisner, B. E., Peres-Neto, P. R., Lindström, E. S., Barnett, A. & Longhi, M. L. The role of environmental and spatial processes in structuring lake communities from bacteria to fish. Ecology 87, 2985–2991 (2006).PubMed 
PubMed Central 
Article 

Google Scholar 
60.De Bie, T. et al. Body size and dispersal mode as key traits determining metacommunity structure of aquatic organisms. Ecol. Lett. 15, 740–747 (2012).PubMed 
Article 

Google Scholar 
61.Kareiva, P. Population dynamics in spatially complex environments: Theory and data. Philos. Trans. R. Soc. Lond. Ser. B: Biol. Sci. 330, 175–190 (1990).ADS 
Article 

Google Scholar 
62.Murray, J. Mathematical Biology II: Spatial Models and Biomedical Applications Vol. 3 (Springer, 2001).
Google Scholar 
63.Rietkerk, M. & Van de Koppel, J. Regular pattern formation in real ecosystems. Trends Ecol. Evolut. 23, 169–175 (2008).Article 

Google Scholar 
64.Pedersen, E. J., Marleau, J. N., Granados, M., Moeller, H. V. & Guichard, F. Nonhierarchical dispersal promotes stability and resilience in a tritrophic metacommunity. Am. Nat. 187, E116–E128 (2016).PubMed 
Article 

Google Scholar 
65.Haegeman, B. & Loreau, M. General relationships between consumer dispersal, resource dispersal and metacommunity diversity. Ecol. Lett. 17, 175–184 (2014).PubMed 
Article 

Google Scholar 
66.Amarasekare, P. Spatial dynamics of foodwebs. Annu. Rev. Ecol. Evol. Syst. 39, 479–500 (2008).Article 

Google Scholar 
67.Fronhofer, E. A., Klecka, J., Melián, C. J. & Altermatt, F. Condition-dependent movement and dispersal in experimental metacommunities. Ecol. Lett. 18, 954–963 (2015).PubMed 
Article 
PubMed Central 

Google Scholar 
68.Toscano, B. J., Gownaris, N. J., Heerhartz, S. M. & Monaco, C. J. Personality, foraging behavior and specialization: integrating behavioral and food web ecology at the individual level. Oecologia 182, 55–69 (2016).ADS 
PubMed 
Article 
PubMed Central 

Google Scholar 
69.Fronhofer, E. A. et al. Bottom-up and top-down control of dispersal across major organismal groups. Nat. Ecol. Evolut. 2, 1859–1863 (2018).Article 

Google Scholar 
70.Gross, T. & Feudel, U. Generalized models as a universal approach to the analysis of nonlinear dynamical systems. Phys. Rev. E 73, 016205 (2006).ADS 
Article 
CAS 

Google Scholar 
71.Yeakel, J. D., Stiefs, D., Novak, M. & Gross, T. Generalized modeling of ecological population dynamics. Theor. Ecol. 4, 179–194 (2011).Article 

Google Scholar 
72.Hirt, M. R. et al. Bridging scales: Allometric random walks link movement and biodiversity research. Trends Ecol. Evolut. 33, 701–712 (2018).Article 

Google Scholar 
73.Othmer, H. G. & Scriven, L. Non-linear aspects of dynamic pattern in cellular networks. J. Theor. Biol. 43, 83–112 (1974).ADS 
CAS 
PubMed 
Article 

Google Scholar 
74.Estes, J. A. et al. Trophic downgrading of planet earth. Science 333, 301–306 (2011).ADS 
CAS 
PubMed 
Article 

Google Scholar 
75.Krause, A. E., Frank, K. A., Mason, D. M., Ulanowicz, R. E. & Taylor, W. W. Compartments revealed in food-web structure. Nature 426, 282–285 (2003).ADS 
CAS 
PubMed 
Article 

Google Scholar 
76.Post, D. M., Conners, M. E. & Goldberg, D. S. Prey preference by a top predator and the stability of linked food chains. Ecology 81, 8–14 (2000).Article 

Google Scholar 
77.Neutel, A.-M., Heesterbeek, J. A. & de Ruiter, P. C. Stability in real food webs: Weak links in long loops. Science 296, 1120–1123 (2002).ADS 
CAS 
PubMed 
Article 

Google Scholar 
78.Leibold, M. A. et al. The metacommunity concept: A framework for multi-scale community ecology. Ecol. Lett. 7, 601–613 (2004).Article 

Google Scholar 
79.Brown, J. H., Gillooly, J. F., Allen, A. P., Savage, V. M. & West, G. B. Toward a metabolic theory of ecology. Ecology 85, 1771–1789 (2004).Article 

Google Scholar 
80.Jenkins, D. G. et al. Does size matter for dispersal distance?. Glob. Ecol. Biogeogr. 16, 415–425 (2007).Article 

Google Scholar 
81.Stevens, V. M. et al. A comparative analysis of dispersal syndromes in terrestrial and semi-terrestrial animals. Ecol. Lett. 17, 1039–1052 (2014).PubMed 
Article 
PubMed Central 

Google Scholar 
82.Guzman, L. M. & Srivastava, D. S. Genomic variation among populations provides insight into the causes of metacommunity survival. Ecology 101, e03182 (2020).PubMed 
Article 
PubMed Central 

Google Scholar 
83.Leitch, K. J., Ponce, F. V., Dickson, W. B., van Breugel, F. & Dickinson, M. H. The long-distance flight behavior of drosophila supports an agent-based model for wind-assisted dispersal in insects. Proc. Natl. Acad. Sci. 118 (2021).84.Bowman, J., Jaeger, J. A. & Fahrig, L. Dispersal distance of mammals is proportional to home range size. Ecology 83, 2049–2055 (2002).Article 

Google Scholar 
85.Shanks, A. L., Grantham, B. A. & Carr, M. H. Propagule dispersal distance and the size and spacing of marine reserves. Ecol. Appl. 13, 159–169 (2003).Article 

Google Scholar 
86.Kartascheff, B., Heckmann, L., Drossel, B. & Guill, C. Why allometric scaling enhances stability in food web models. Theor. Ecol. 3, 195–208 (2010).Article 

Google Scholar 
87.Hudson, L. N. & Reuman, D. C. A cure for the plague of parameters: Constraining models of complex population dynamics with allometries. Proc. R. Soc. B: Biol. Sci. 280, 20131901 (2013).Article 

Google Scholar 
88.Brose, U. et al. Predator traits determine food-web architecture across ecosystems. Nat. Ecol. Evolut. 3, 919–927 (2019).Article 

Google Scholar 
89.Heino, J. et al. Metacommunity organisation, spatial extent and dispersal in aquatic systems: patterns, processes and prospects. Freshw. Biol. 60, 845–869 (2015).Article 

Google Scholar 
90.Siegel, D. et al. The stochastic nature of larval connectivity among nearshore marine populations. Proc. Natl. Acad. Sci. 105, 8974–8979 (2008).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
91.Pillai, P., Loreau, M. & Gonzalez, A. A patch-dynamic framework for food web metacommunities. Theor. Ecol. 3, 223–237 (2010).Article 

Google Scholar 
92.Pillai, P., Gonzalez, A. & Loreau, M. Metacommunity theory explains the emergence of food web complexity. Proc. Natl. Acad. Sci. 108, 19293–19298 (2011).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
93.Plitzko, S. J. & Drossel, B. The effect of dispersal between patches on the stability of large trophic food webs. Theor. Ecol. 8, 233–244 (2015).Article 

Google Scholar 
94.Guichard, F. Recent advances in metacommunities and meta-ecosystem theories. F1000Research 6 (2017).95.Hata, S., Nakao, H. & Mikhailov, A. S. Dispersal-induced destabilization of metapopulations and oscillatory turing patterns in ecological networks. Sci. Rep. 4, 3585 (2014).ADS 
PubMed 
PubMed Central 
Article 

Google Scholar 
96.White, K. & Gilligan, C. Spatial heterogeneity in three species, plant-parasite-hyperparasite, systems. Philos. Trans. R. Soc. Lond. Ser. B: Biol. Sci. 353, 543–557 (1998).Article 

Google Scholar 
97.Gibert, J. P. & Yeakel, J. D. Laplacian matrices and turing bifurcations: Revisiting levin 1974 and the consequences of spatial structure and movement for ecological dynamics. Theor. Ecol. 12, 265–281 (2019).Article 

Google Scholar 
98.Fox, J. W., Vasseur, D., Cotroneo, M., Guan, L. & Simon, F. Population extinctions can increase metapopulation persistence. Nat. Ecol. Evolut. 1, 1271–1278 (2017).Article 

Google Scholar 
99.Hastings, A. Food web theory and stability. Ecology 69, 1665–1668 (1988).Article 

Google Scholar 
100.Anderson, H., Hutson, V. & Law, R. On the conditions for permanence of species in ecological communities. Am. Nat. 139, 663–668 (1992).Article 

Google Scholar 
101.Haydon, D. Pivotal assumptions determining the relationship between stability and complexity: An analytical synthesis of the stability-complexity debate. Am. Nat. 144, 14–29 (1994).Article 

Google Scholar 
102.Chen, X. & Cohen, J. E. Global stability, local stability and permanence in model food webs. J. Theor. Biol. 212, 223–235 (2001).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
103.Bjørnstad, O. N., Ims, R. A. & Lambin, X. Spatial population dynamics: Analyzing patterns and processes of population synchrony. Trends Ecol. Evolut. 14, 427–432 (1999).Article 

Google Scholar 
104.Ims, R. A. & Andreassen, H. P. Spatial synchronization of vole population dynamics by predatory birds. Nature 408, 194–196 (2000).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
105.Sundell, J. et al. Large-scale spatial dynamics of vole populations in Finland revealed by the breeding success of vole-eating avian predators. J. Anim. Ecol. 73, 167–178 (2004).Article 

Google Scholar 
106.Ripple, W. J. et al. Status and ecological effects of the world’s largest carnivores. Science 343 (2014).107.McCauley, D. J. et al. Marine defaunation: Animal loss in the global ocean. Science 347 (2015).108.Parsons, T. The removal of marine predators by fisheries and the impact of trophic structure. Mar. Pollut. Bull. 25, 51–53 (1992).Article 

Google Scholar 
109.Baum, J. K. & Worm, B. Cascading top-down effects of changing oceanic predator abundances. J. Anim. Ecol. 78, 699–714 (2009).PubMed 
Article 
PubMed Central 

Google Scholar 
110.Albert, C. H., Rayfield, B., Dumitru, M. & Gonzalez, A. Applying network theory to prioritize multispecies habitat networks that are robust to climate and land-use change. Conserv. Biol. 31, 1383–1396 (2017).PubMed 
Article 

Google Scholar 
111.Schiesari, L. et al. Towards an applied metaecology. Perspect. Ecol. Conserv. 17, 172–181 (2019).
Google Scholar 
112.Vermaat, J. E., Dunne, J. A. & Gilbert, A. J. Major dimensions in food-web structure properties. Ecology 90, 278–282 (2009).PubMed 
Article 

Google Scholar 
113.White, J. W., Rassweiler, A., Samhouri, J. F., Stier, A. C. & White, C. Ecologists should not use statistical significance tests to interpret simulation model results. Oikos 123, 385–388 (2014).Article 

Google Scholar 
114.R Core Team. R: A Language and Environment for Statistical Computing. https://www.R-project.org/. (R Foundation for Statistical Computing, 2020).115.Aufderheide, H., Rudolf, L., Gross, T. & Lafferty, K. D. How to predict community responses to perturbations in the face of imperfect knowledge and network complexity. Proc. R. Soc. B Biol. Sci. 280, 20132355 (2013).Article 

Google Scholar  More

Evaluating indicatorsAmong the selected parameters the ratio of rural to the urban population directly relates to the per capita renewable water, whereas the population density, internet users, and education index exhibit an inverse relation with the per capita renewable water worldwide. It means the per capita renewable water decreases with decreasing rural to urban population and increasing population density, internet users, and education index. The urban population has increased in developing regions, which feature increasing population density. People’s health is threatened by poor urban sanitary infrastructure leading to disease and social decay. Increasing population density and a reduction in per capita renewable water inflict social harm and disrupt society’s economic growth58. Population density also is positively related to the relative number of elderly and social vulnerability because potential casualties increase with population size40. On the other hand, with the increase of Internet users and education index, the per capita renewable water has increased. As long as the knowledge and awareness of communities improved, the consumption algorithm decreased, leading to a reduction of renewable water per capita. Therefore, the level of literacy and knowledge for a community can be the basis for making the right decisions in agriculture, health, natural resource management, and other activities related to water resources for decision-makers. The latter situation calls for better communication among water users through social media and improved education to learn and develop optimal water management.Evaluating models and developing hydro-social equationsThree soft-computing approaches, namely ANN-LM, ANFIS-SC, and GEP, were applied to develop predictive equations with social indicators worldwide. The ANN-Levenberg–Marquardt (LM) backpropagation algorithm with one hidden layer was applied, and the hidden nodes’ number was determined by trial and error. A hybrid algorithm was combined with the ANFIS-SC models. There is no rule for determining the radii values of the ANFIS-SC models. The final radii values were determined by trial-and-error.The numbers of neurons in the ANN-LM models and the radii values of the ANFIS-SC models are listed in Table 4. The activation functions of the output nodes were linear for all the continents. The activation functions of the hidden nodes of the ANN-LM models for the P1 through P4 indicators were respectively the tangent sigmoid, tangent sigmoid, tangent sigmoid, and logarithm sigmoid for Africa; the activation functions of the proportion of rural to urban population was the tangent sigmoid for all the continents. Table 5 lists the results of the soft computing optimal models’ estimates of the proportion of rural to urban population (PRUP), population density (PD), internet users (IU), and education index (EI), denoted respectively by P1 through P4, during the test period in the world’s continents. Figures 4 and 5 display the characteristics of ANN (the number of neurons and activation functions of hidden and output layers) and ANFIS-SC (radii values) models, respectively. The values of R and RMSE for Africa corresponding to the ANN-LM models were respectively (0.921, 0.981, 0.858, 0.862) and (0.193, 0.058, 0.190, 0.172) associated with the PRUP, PD, IU, and EI parameters, respectively. The values of R and RMSE for Africa corresponding to the ANFIS-SC models equaled respectively (0.933, 0.991, 0.868, 0.891) and (0.130, 0.044, 0.186, 0.156) for the P1 through P4 parameters, respectively. Concerning the GEP models, the root relative squared error (RRSE) was selected as the pressure tree’s fitness function. The values of RMSE for GEP models equaled (0.084, 0.029, 0.178, 0.135), (0.197, 0.056, 0.152, 0.163), (0.151, 0.036, 0.123, 0.210), (0.182, 0.039, 0.148, 0.204) and (0.141, 0.030, 0.226, 0.082) for Africa, America, Asia, Europe, and Oceania, respectively. Table 5 results for the R, RMSE, and MAE values establish the GEP model estimates of PRUP, PD, IU, and EI indicators had the highest R values and the lowest RMSE values. The average R values of the best models (GEP) for all selected social parameters equaled 0.942, 0.909, 0.910, 0.889, and 0.947 for Africa, America, Asia, Europe, and Oceania, respectively. These results indicate the climatic characteristics of the continents influence the performance of the models. The models’ performances for Africa and Oceania associated with the type B dominant Koppen climate classification was the best. The models’ performances for Asia and America that have similar climatic classification were nearly equal. The average model performance for Europe in the type D climate classification was the poorest among the continents.Table 4 The characteristics of ANN (the number of neurons) and ANFIS (radii values) models corresponding to social indicators and continents.Full size tableTable 5 The results of soft computing optimal models corresponding to the testing period in the world’s continents.Full size tableFigure 4The characteristics of optimal ANN models; showing the number of neurons and activation functions of hidden and output layers.Full size imageFigure 5The characteristic of optimal ANFIS-SC model showing the radii values.Full size imageFigures 6, 7, 8, 9 and 10 show the observed and estimated social parameters obtained with the soft-computing models during the test period in Africa, America, Asia, Europe, and Oceania, respectively. Figure 11 compares the R, RMSE, and MAE values from the soft-computing models. The R values for soft-computing models are close to 1, with the quality relations being: RGEP  > RANFIS-SC  > RANN-LM for all social indicators. Figure 11 establishes that the ANFIS-SC model exceeded the ANN-LM models’ performance. Also, the GEP models had better performance than the ANFIS-SC and ANN-LM for estimating the proportion of rural to urban population (PRUP), population density (PD), internet users (IU), and education index (EI) parameters in Africa, America, Asia, Europe, and Oceania.Figure 6Observed and estimated social parameters during the testing period in Africa.Full size imageFigure 7Observed and estimated social parameters during the testing period in America.Full size imageFigure 8Observed and estimated social parameters during the testing period in Asia.Full size imageFigure 9Observed and estimated social parameters during the testing period in Europe.Full size imageFigure 10Observed and estimated social indicators during the testing period in Oceania.Full size imageFigure 11Comparison of R, RMSE and MAE values corresponding to the soft computing methods.Full size imageThe main advantage of the GEP over other soft computing methods (e.g., ANFIS and ANN) is in producing predictive equations. The equations obtained with the optimal models for the social indicators (i.e., the proportion of rural to urban population (PRUP), population density (PD), internet users (IU), and education index (EI) in Africa, America, Asia, Europe, and Oceania) are listed in Table 6. The equations that the GEP model discovers as a structure do not necessarily correspond to reality. The equations listed in Table 6 merely show the optimal equations extracted from the model after the evolution, for all indicators and in all basins (considering renewable water per capita as a decision variable).Table 6 Mathematical equations governing hydro-social indicators.Full size tableThe performance of the GEP models in estimating the social indicators in three ranges of values, namely, 20% of the maximum estimated values (20%max), 60% of median estimated values (60%mid or 20%min to 20%max), and 20% of minimum estimated values (20%min), during the test period for the proportion of rural to urban population (PRUP), population density (PD), internet users (IU) and the education index (EI) parameters of Africa, America, Asia, Europe, and Oceania are listed in Table 7. Table 7’s results indicate there is not a regular rule to determine the best-cited ranges performances. The education index and the population density have the lowest and highest R values among the other parameters in the three different ranges (20%max, 60%mid, and 20%min) in Africa, America, Asia, Europe, and Oceania. Therefore, the results indicate a strong pattern of association between the population density parameter and water resources status in all continents of the world.Table 7 The performance of GEP models with respect to selected ranges.Full size tableFigure 12 depicts the distribution of estimated data values of the social parameters (i = 1, 2, 3, 4) and their comparison through the continents. The box plots are a graphic display integrating multiple numerical relations. One approach to understanding the distribution or dispersion of data is through the box diagram, which is based on the “minimum,” “first quartile-Q1(0.25%)”, “median (0.50%)”, “third quartile-Q3(0.75%)” and “maximum” statistical indicators. Figure 12 shows Oceania and Africa exhibit the smallest and largest values of the rural to urban population, respectively. America has the lowest values of the first to the third quartile. The estimated population density value in Europe has the most values in the third quartile (0.75%). The median values of estimated internet users have the smallest and largest values in Africa and Europe, respectively. America has the lowest values of the first quartile, median, third quartile, and maximum values associated with the estimated education index values among the continents.Figure 12Distribution of estimated data values of social indicators (Pi, i = 1, 2, …, 4).Full size imageThe summary of hydro-social equations performance is listed in Table 8, where it is seen the best models’, performances are such that PD  > PRUP  > EI  > IU, PD  > IU  > EI  > PRUP, PD  > IU  > PRUP  > EI, PD  > PRUP  > IU  > EI and PD  > EI  > IU  > PRUP for Africa, America, Asia, Europe, and Oceania, respectively.Table 8 Summary of hydro-social equations performance.Full size tableThis paper’s results indicate the pattern of association between social parameters and water resources is complex. Renewable water per capita was estimated using social indicators PRUP, PD, IU, and EI based on gene expression programming. The results of GEP to estimate RWPC corresponding to the testing period in the world’s continents as listed in Table 9. The values of RMSE for optimal GEP models equaled 0.089, 0.058, 0.042, 0.049, and 0.036 for Africa, America, Asia, Europe, and Oceania, respectively. Figure 13 displays the observed and estimated RWPC parameter during the test period in the world’s continents. The equations obtained with the optimal models for the renewable water per capita in Africa, America, Asia, Europe, and Oceania are listed in Table 10. The fitted equations can be applied at variable spatial and temporal scales. The derived equations imply that water resources in Africa and Oceania are governed by the PRUP, PD, IU, and EI indicators. Also, the PRUP, PD, and IU indicators in Europe and PD and IU indicators in America and Asia have the most influence on their water resources status. The association between social parameters and water resources in all continents is variable. The linking of these social indicators with the per capita renewable water is a function of the countries’ cultural and economic conditions, thus bearing on the future management and policymaking across continents. This study’s results concerning hydro-social indicators are consistent with the findings by Forouzani et al.2, Carey et al.15, Lima et al.25, Pande et al.7, Diep et al.26, and Diaz et al.22.Table 9 The results of GEP estimating RWPC corresponding to the testing period in the world’s continents.Full size tableFigure 13Observed and estimated RWPC parameters during the test period in the world’s continents.Full size imageTable 10 Mathematical equations governing hydro-social indicators.Full size tableThis paper’s results establish the importance of examining the interactions between climate, the status of water resources, and social indicators. The state and social conditions of a country reflect the status of its water resources. Therefore, this study has shown how significant an impact the management and planning of a country can have on its water resources. Each successful water resources project rests on a successful social setting. More

1.Collins, J. P., & Crump, M. L. Extinction in Our Times: Global Amphibian Decline. (2009).2.Catenazzi, A. State of the world’s amphibians. Annu. Rev. Environ. Resour. 40, 91–119 (2015).Article 

Google Scholar 
3.González-del-Pliego, P. et al. Phylogenetic and trait-based prediction of extinction risk for data deficient amphibians. Curr. Biol. 29, 1557–1563 (2019).PubMed 
Article 
CAS 
PubMed Central 

Google Scholar 
4.Lips, K. R. et al. Emerging infectious disease and the loss of biodiversity in a Neotropical amphibian community. PNAS 102, 3165–3170 (2006).ADS 
Article 
CAS 

Google Scholar 
5.Lips, K. R., Diffendorfer, J., Mendelson, J. R. & Sears, M. W. Riding the wave: Reconciling the roles of disease and climate change in amphibian declines. PLoS Biol. 6, e72 (2008).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
6.Berger, L. et al. Chytridiomycosis causes amphibian mortality associated with population declines in the rain forests of Australia and Central America. PNAS 95, 9031–9036 (1998).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
7.Longcore, J. E., Pessier, A. P. & Nichols, D. K. Batrachochytrium dendrobatidis gen. et sp. nov., a chytrid pathogenic to amphibians. Mycologia 91, 219–227 (1999).Article 

Google Scholar 
8.Berger, L. et al. History and recent progress on chytridiomycosis in amphibians. Fungal Ecol. 19, 89–99 (2016).Article 

Google Scholar 
9.Martel, A. et al. Batrachochytrium salamandrivorans sp. nov. causes lethal chytridiomycosis in amphibians. PNAS 110, 15325–15329 (2013).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
10.Scheele, B. C. et al. Amphibian fungal panzootic causes catastrophic and ongoing loss of biodiversity. Science 363, 1459–1463 (2019).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
11.Lambert, M. R. et al. Comment on “Amphibian fungal panzootic causes catastrophic and ongoing loss of biodiversity”. Science https://doi.org/10.1126/science.aay1838 (2020).Article 
PubMed 
PubMed Central 

Google Scholar 
12.Scheele, B. C. et al. Response to Comment on “Amphibian fungal panzootic causes catastrophic and ongoing loss of biodiversity”. Science https://doi.org/10.1126/science.aay2905 (2020).Article 
PubMed 
PubMed Central 

Google Scholar 
13.Puschendorf, R. et al. Distribution models for the amphibian chytrid Batrachochytrium dendrobatidis in Costa Rica: Proposing climatic refuges as a conservation tool. Divers. Distrib. 15, 401–408 (2009).Article 

Google Scholar 
14.Zumbado-Ulate, H. et al. Endemic infection of Batrachochytrium dendrobatidis in Costa Rica: Implications for amphibian conservation at regional and species level. Diversity 11, 129 (2019).Article 

Google Scholar 
15.Crawford, A. J., Lips, K. R. & Bermingham, E. Epidemic disease decimates amphibian abundance, species diversity, and evolutionary history in the highlands of central Panama. PNAS 107, 13777–13782 (2010).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
16.Woodhams, D. C. et al. Chytridiomycosis and amphibian population declines continue to spread eastward in Panama. EcoHealth 5, 268–274 (2008).PubMed 
Article 
PubMed Central 

Google Scholar 
17.Catenazzi, A., Lehr, E., Rodriguez, L. & Vredenburg, V. Batrachochytrium dendrobatidis and the collapse of anuran species richness and abundance in the upper Manu National Park, southeastern Peru. Conserv. Biol. 25, 382–391 (2011).PubMed 
Article 
PubMed Central 

Google Scholar 
18.Voyles, J. et al. Pathogenesis of chytridiomycosis, a cause of catastrophic amphibian declines. Science 326, 582–585 (2009).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
19.James, T. et al. Disentangling host, pathogen, and environmental determinants of a recently emerged wildlife disease: Lessons from the first 15 years of amphibian chytridiomycosis research. Ecol. Evol. 5, 4079–4097 (2015).PubMed 
PubMed Central 
Article 

Google Scholar 
20.Soto-Azat, C. et al. Xenopus laevis and emerging amphibian pathogens in Chile. EcoHealth 13, 775–783 (2016).PubMed 
Article 
PubMed Central 

Google Scholar 
21.Ron, S. Predicting the distribution of the amphibian pathogen Batrachochytrium dendrobatidis in the new world. Biotropica 37, 209–221 (2005).Article 

Google Scholar 
22.Rödder, D., Kielgast, J. & Lötters, S. Future potential distribution of the emerging amphibian chytrid fungus under anthropogenic climate change. Dis. Aquat. Org. 92, 201–207 (2010).Article 

Google Scholar 
23.Murray, K. A. et al. Assessing spatial patterns of disease risk to biodiversity: Implications for the management of the amphibian pathogen, Batrachochytrium dendrobatidis. J. Appl. Ecol. 48, 163–173 (2011).Article 

Google Scholar 
24.Liu, X., Rohr, J. & Li, Y. Climate, vegetation, introduced hosts and trade shape a global wildlife pandemic. Proc. Biol. Sci. 280, 20122506 (2013).PubMed 
PubMed Central 

Google Scholar 
25.Olson, D. H. et al. Mapping the global emergence of Batrachochytrium dendrobatidis, the amphibian chytrid fungus. PLoS ONE 8, e56802 (2013).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
26.Penner, J. et al. West Africa—A safe haven for frogs? A sub-continental assessment of the chytrid fungus (Batrachochytrium dendrobatidis). PLoS ONE 8, e56236 (2013).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
27.Xie, G. Y., Olson, D. H. & Blaustein, A. R. Projecting the global distribution of the emerging amphibian fungal pathogen, Batrachochytrium dendrobatidis, based on IPCC climate futures. PLoS ONE 11, e0160746 (2016).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
28.Harris, R. N. et al. Skin microbes on frogs prevent morbidity and mortality caused by a lethal skin fungus. ISME J. 3, 818–824 (2009).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
29.Searle, C. L. et al. Differential host susceptibility to Batrachochytrium dendrobatidis, an emerging amphibian pathogen. Conserv. Biol. 25, 965–974 (2011).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
30.Farrer, R. et al. Multiple emergences of genetically diverse amphibian-infecting chytrids include a globalized hypervirulent recombinant lineage. PNAS 108, 18732–18736 (2011).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
31.Lips, K. R. Overview of chytrid emergence and impacts on amphibians. Philos. Trans. R. Soc. B. 371, 20150465 (2016).Article 

Google Scholar 
32.Bolom-Huet, R., Pineda, E., Díaz-Fleischer, F., Muñoz-Alonso, A. L. & Galindo-González, J. Known and estimated distribution in Mexico of Batrachochytrium dendrobatidis, a pathogenic fungus of amphibians. Biotropica 51, 731–746 (2019).Article 

Google Scholar 
33.Zumbado-Ulate, H., García-Rodríguez, A. & Searle, C. L. Species distribution models predict the geographic expansion of an enzootic amphibian pathogen. Biotropica 53, 221–231 (2021).Article 

Google Scholar 
34.Berger, L. et al. Effect of season and temperature on mortality in amphibians due to chytridiomycosis. Aust. Vet. J. 82, 434–439 (2004).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
35.Clare, F. C. et al. Climate forcing of an emerging pathogenic fungus across a montane multi-host community. Philos. Trans. R. Soc. B. 371, 20150454 (2016).Article 
CAS 

Google Scholar 
36.Bacigalupe, L. D., Soto-Azat, C., García-Vera, C., Barría-Oyarzo, I. & Rezende, E. L. Effects of amphibian phylogeny, climate and human impact on the occurrence of the amphibian-killing chytrid fungus. Glob. Change Biol. 23, 3543–3553 (2017).ADS 
Article 

Google Scholar 
37.Raffel, T., Michel, P., Sites, W. & Rohr, J. What drives chytrid infections in newt populations? Associations with substrate, temperature, and shade. EcoHealth 7, 526–536 (2010).PubMed 
Article 
PubMed Central 

Google Scholar 
38.Pounds, J. A. et al. Widespread amphibian extinctions from epidemic disease driven by global warming. Nature 439, 161–167 (2006).ADS 
CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
39.Hudson, M. et al. Reservoir frogs: Seasonality of Batrachochytrium dendrobatidis infection in robber frogs. PeerJ 7, e7021 (2019).PubMed 
PubMed Central 
Article 

Google Scholar 
40.Kriger, K. M. & Hero, J. M. Large-scale seasonal variation in the prevalence and severity of chytridiomycosis. J. Zool. 271, 352–359 (2007).
Google Scholar 
41.Longo, A. V., Burrowes, P. A. & Joglar, R. L. Seasonality of Batrachochytrium dendrobatidis infection in direct-developing frogs suggests a mechanism for persistence. Dis. Aquat. Org. 92, 253–260 (2010).Article 

Google Scholar 
42.Zumbado-Ulate, H., Bolaños, F., Gutiérrez-Espeleta, G. & Puschendorf, R. Extremely low prevalence of Batrachochytrium dendrobatidis in frog populations from Neotropical dry forest of Costa Rica supports the existence of a climatic refuge from disease. EcoHealth 11, 593–602 (2014).PubMed 
Article 
PubMed Central 

Google Scholar 
43.Bacigalupe, L. D. et al. The amphibian-killing fungus in a biodiversity hotspot: Identifying and validating high-risk areas and refugia. Ecosphere. 10, e02724 (2019).Article 

Google Scholar 
44.Flechas, S. V. et al. Current and predicted distribution of the pathogenic fungus Batrachochytrium dendrobatidis in Colombia, a hotspot of amphibian biodiversity. Biotropica 49, 685–694 (2017).Article 

Google Scholar 
45.Lampo, M. et al. Batrachochytrium dendrobatidis in Venezuela. Herpetol. Rev. 39, 449 (2008).
Google Scholar 
46.Valenzuela-Sánchez, A. et al. Genomic epidemiology of the emerging pathogen Batrachochytrium dendrobatidis from native and invasive amphibian species in Chile. Transbound. Emerg. Dis. 65, 309–314 (2018).PubMed 
Article 
PubMed Central 

Google Scholar 
47.O’Hanlon, S. et al. Recent Asian origin of chytrid fungi causing global amphibian declines. Science 360, 621–627 (2018).ADS 
PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
48.Soto-Azat, C. et al. The population decline and extinction of Darwin’s frogs. PLoS ONE 8, e66957 (2013).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
49.Soto-Azat, C. et al. ASG Chile leads update of the extinction risk of Chilean amphibians for the IUCN red list of threatened speciesTM. FrogLog 23, 6–7 (2015).
Google Scholar 
50.Mora, M. et al. High abundance of invasive African clawed frog Xenopus laevis in Chile: Challenges for their control and updated invasive distribution. Manag. Biol. Invasions. 10, 377–388 (2019).Article 

Google Scholar 
51.Solís, R., Penna, M., De la Riva, I., Fisher, M. & Bosch, J. Presence of Batrachochytrium dendrobatdis in anurans from the Andes highlands of northern Chile. Herpetol. J. 24, 55–59 (2015).
Google Scholar 
52.Soto-Azat, C. et al. Is Chytridiomycosis driving Darwin’s frogs to extinction?. PLoS ONE 8, e79862 (2013).ADS 
PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
53.Valenzuela-Sánchez, A. et al. Cryptic disease-induced mortality may cause host extinction in an apparently stable host–parasite system. Proc. Biol. Sci. 284, 20171176 (2017).PubMed 
PubMed Central 

Google Scholar 
54.Lips, K. R., Reeve, J. D. & Witters, L. R. Ecological traits predicting amphibian population declines in Central America. Conserv. Biol. 17, 1078–1088 (2003).Article 

Google Scholar 
55.Hero, J. M., Williams, S. E. & Magnusson, W. E. Ecological traits of declining amphibians in upland areas of eastern Australia. J. Zool. 267(3), 221–232 (2005).Article 

Google Scholar 
56.Kriger, K. M. & Hero, J. M. Altitudinal distribution of chytrid (Batrachochytrium dendrobatidis) infection in subtropical Australian frogs. Austral. Ecol. 33(8), 1022–1032 (2008).Article 

Google Scholar 
57.Skerratt, L. F. et al. Spread of chytridiomycosis has caused the rapid global decline and extinction of frogs. EcoHealth 4, 125 (2007).Article 

Google Scholar 
58.Langwig, K. et al. Context-dependent conservation responses to emerging wildlife diseases. Front. Ecol. Environ. 13, 195–202 (2015).Article 

Google Scholar 
59.Shaw, S. D. et al. The distribution and host range of Batrachochytrium dendrobatidis in New Zealand, 1930–2010. Ecology 94, 2108–2111 (2013).Article 

Google Scholar 
60.Ghirardi, R. et al. Endangered amphibians infected with the chytrid fungus Batrachochytrium dendrobatidis in austral temperate wetlands from Argentina. Herpetol. J. 24, 129–133 (2014).
Google Scholar 
61.Bielby, J., Cooper, N., Cunningham, A., Garner, T. & Purvis, A. Predicting susceptibility to future declines in the world’s frogs. Conserv. Lett. 1, 82–90 (2008).Article 

Google Scholar 
62.Barrionuevo, S. & Mangione, S. Chytridiomycosis in two species of Telmatobius (Anura: Leptodactylidae) from Argentina. Dis. Aquat. Org. 73, 171–174 (2006).Article 

Google Scholar 
63.Seimon, T. A. et al. Upward range extension of Andean anurans and chytridiomycosis to extreme elevations in response to tropical deglaciation. Glob. Change Biol. 13, 288–299 (2007).ADS 
Article 

Google Scholar 
64.Burrowes, P. A. & De la Riva, I. Unraveling the historical prevalence of the invasive chytrid fungus in the Bolivian Andes: Implications in recent amphibian declines. Biol. Invasions. 19, 1781–1794 (2017).Article 

Google Scholar 
65.Vredenburg, V. T., Knapp, R., Tunstall, T. & Briggs, C. Dynamics of an emerging disease drive large-scale amphibian population extinctions. PNAS 107, 9689–9694 (2010).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
66.Azat, C. et al. A flagship for Austral temperate forest conservation: an action plan for Darwin’s frogs bringing together key stakeholders. Oryx 55, 356–363 (2021). Article 

Google Scholar 
67.Pilliod, D. S. et al. Effects of amphibian chytrid fungus on individual survival probability in wild boreal toads. Conserv. Biol. 24, 1259–1267 (2010).PubMed 
Article 
PubMed Central 

Google Scholar 
68.Walker, S. F. et al. Factors driving pathogenicity vs. prevalence of amphibian panzootic chytridiomycosis in Iberia. Ecol. Lett. 13, 372–382 (2010).PubMed 
Article 
PubMed Central 

Google Scholar 
69.Kriger, K. M., Pereoglou, F. & Hero, J. M. Latitudinal variation in the prevalence and intensity of chytrid (Batrachochytrium dendrobatidis) infection in eastern Australia. Conserv. Biol. 21, 1280–1290 (2007).PubMed 
Article 
PubMed Central 

Google Scholar 
70.Petersen, C. E., Lovich, R. E., Phillips, C. A., Dreslik, M. J. & Lannoo, M. J. Prevalence and seasonality of the amphibian chytrid fungus Batrachochytrium dendrobatidis along widely separated longitudes across the United States. EcoHealth 13, 368–382 (2016).PubMed 
Article 
PubMed Central 

Google Scholar 
71.Thorpe, C. J. et al. Climate structuring of Batrachochytrium dendrobatidis infection in the threatened amphibians of the northern Western Ghats, India. R. Soc. Open Sci. 5, 180211 (2018).ADS 
PubMed 
PubMed Central 
Article 

Google Scholar 
72.Sonn, J. M., Utz, R. M. & Richards-Zawacki, C. L. Effects of latitudinal, seasonal, and daily temperature variations on chytrid fungal infections in a North American frog. Ecosphere 10, e02892 (2019).Article 

Google Scholar 
73.Raffel, T. R., Halstead, N. T., McMahon, T. A., Davis, A. K. & Rohr, J. R. Temperature variability and moisture synergistically interact to exacerbate an epizootic disease. Proc. Biol. Sci. 282, 20142039 (2015).PubMed 
PubMed Central 

Google Scholar 
74.Woodhams, D. C. & Alford, R. A. Ecology of chytridiomycosis in rainforest stream frog assemblages of tropical Queensland. Conserv. Biol. 19, 1449–1459 (2005).Article 

Google Scholar 
75.Adams, M. J. et al. Using occupancy models to understand the distribution of an amphibian pathogen, Batrachochytrium dendrobatidis. Ecol. Appl. 20, 289–302 (2010).PubMed 
Article 
PubMed Central 

Google Scholar 
76.Fisher, M., Garner, T. & Walker, S. F. Global emergence of Batrachochytrium dendrobatidis and amphibian chytridiomycosis in space, time, and host. Annu. Rev. Microbiol. 63, 291–310 (2009).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
77.Schloegel, L. M. et al. Novel, panzootic and hybrid genotypes of amphibia chytridiomycosis associated with the bullfrog trade. Mol. Ecol. 21, 5162–5177 (2012).PubMed 
Article 
PubMed Central 

Google Scholar 
78.Wilson, E. A., Briggs, C. J. & Dudley, T. L. Invasive African clawed frogs in California: A reservoir for or predator against the chytrid fungus?. PLoS ONE 13, e0191537 (2018).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
79.Becker, C. G., Longo, A. V., Haddad, C. F. & Zamudio, K. R. Land cover and forest connectivity alter the interactions among host, pathogen and skin microbiome. Proc. Biol. Sci. 284, 20170582 (2017).PubMed 
PubMed Central 

Google Scholar 
80.McCoy, K. A. & Peralta, A. L. Pesticides could alter amphibian skin microbiomes and the effects of Batrachochytrium dendrobatidis. Front. Microbiol. 9, 748 (2018).PubMed 
PubMed Central 
Article 

Google Scholar 
81.Ellis, E. & Ramankutty, N. Putting people in the map: Anthropogenic biomes of the world. Front. Ecol. Environ. 6, 439–447 (2008).Article 

Google Scholar 
82.Rohr, J., Halstead, N. & Raffel, T. Modelling the future distribution of the amphibian chytrid fungus: The influence of climate and human-associated factors. J. Appl. Ecol. 48, 174–176 (2011).Article 

Google Scholar 
83.Myers, N., Mittermeier, R. A., Mittermeier, C. G., Da Fonseca, G. A. & Kent, J. Biodiversity hotspots for conservation priorities. Nature 403, 853–858 (2000).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
84.Echeverria, C., Coomes, D., Hall, M. & Newton, A. Spatially explicit models to analyze forest loss and fragmentation between 1976 and 2020 in southern Chile. Ecol. Model. 212, 439–449 (2008).Article 

Google Scholar 
85.Rodriguez, D., Becker, C., Pupin, C., Haddad, F. & Zamudio, K. Long-term endemism of two highly divergent lineages of the amphibian-killing fungus in the Atlantic Forest of Brazil. Mol. Ecol. 23, 774–787 (2014).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
86.Puschendorf, R., Hodgson, L., Alfors, R. A., Skerrat, L. F. & VanDerWal, J. Underestimated ranges and overlooked refuges from amphibian chytridiomycosis. Divers. Distrib. 19, 1313–1321 (2013).Article 

Google Scholar 
87.Scheele, B. C. et al. After the epidemic: Ongoing declines, stabilizations and recoveries in amphibians afflicted by chytridiomycosis. Biol. Conserv. 206, 37–46 (2017).Article 

Google Scholar 
88.Mendelson, J. R. III., Whitfield, S. M. & Sredl, M. J. A recovery engine strategy for amphibian conservation in the context of disease. Biol. Conserv. 236, 188–191 (2019).Article 

Google Scholar 
89.Van Rooij, P., Martel, A., Haesebrouck, F. & Pasmans, F. Amphibian chytridiomycosis: A review with focus on fungus-host interactions. Vet. Res. 46, 1–22 (2015).Article 
CAS 

Google Scholar 
90.Christie, M. R. & Searle, C. L. Evolutionary rescue in a host–pathogen system results in coexistence not clearance. Evol. Appl. 11, 681–693 (2018).PubMed 
Article 
PubMed Central 

Google Scholar 
91.Bletz, M. C. et al. Mitigating amphibian chytridiomycosis with bioaugmentation: Characteristics of effective probiotics and strategies for their selection and use. Ecol. Lett. 16, 807–820 (2013).PubMed 
Article 
PubMed Central 

Google Scholar 
92.Bosch, J. et al. Successful elimination of a lethal wildlife infectious disease in nature. Biol. Lett. 11, 20150874 (2015).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
93.Olson, D. M. et al. Terrestrial ecoregions of the world: A new map of life on earth. Bioscience 51, 933–938 (2001).Article 

Google Scholar 
94.Pellet, J. & Schmidt, B. R. Monitoring distributions using call surveys: Estimating site occupancy, detection probabilities and inferring absence. Biol. Conserv. 123, 27–35 (2005).Article 

Google Scholar 
95.Drechsler, A. & Bock, D. Ortmann’s funnel trap—A highly efficient tool for monitoring amphibian species. Herpetol. Notes. 3, 13–21 (2010).
Google Scholar 
96.Hudson, M. et al. Dynamics and genetics of a disease-driven species decline to near extinction: Lessons for conservation. Sci. Rep. 6, 1–13 (2016).Article 
CAS 

Google Scholar 
97.R Development Core Team. R: A Language and Environment for Statistical Computing, https://www.R-project.org/ (2019).98.Sanderson, E. W. et al. The human footprint and the last of the wild. Bioscience 52, 891–904 (2002).Article 

Google Scholar 
99.Center for International Earth Science Information Network (CIESIN). Gridded Population of the World, Version 3 (GPWv3). https://doi.org/10.7927/H4639MPP (2005).100.Booth, T. H., Nix, H. A., Busby, J. R. & Hutchinson, M. F. BIOCLIM: The first species distribution modelling package, its early applications and relevance to most current MAXENT studies. Divers. Distrib. 20, 1–9 (2014).Article 

Google Scholar 
101.Center for International Earth Science Information Network (CIESIN). Gridded Species Distribution: Global Amphibian Richness Grids. https://doi.org/10.7927/H4RR1W66 (2015).102.Fick, S. & Hijmans, R. WorldClim 2: New 1-km spatial resolution climate surfaces for global land areas. Int. J. Climatol. 37, 4302–4315 (2017).Article 

Google Scholar 
103.ASTER. ASTER global digital elevation model V003. https://doi.org/10.5067/ASTER/ASTGTM (2018).104.QGIS Development Team. QGIS Geographic Information System. Open Source Geospatial Foundation Project. http://qgis.osgeo.org (2018).105.Zuur, A. F., Ieno, E. N. & Elphick, C. S. A protocol for data exploration to avoid common statistical problems. Methods Ecol. Evol. 1, 3–14 (2010).Article 

Google Scholar 
106.Fox, J. Effect displays in R for generalised linear models. J. Stat. Softw. 8, 1–27 (2003).Article 

Google Scholar 
107.Carpenter, T. E. Methods to investigate spatial and temporal clustering in veterinary epidemiology. Prev. Vet. Med. 48, 303–320 (2001).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
108.Kulldorff, M. A spatial scan statistic. Commun. Stat-Theor. M. 26, 1481–1496 (1997).MathSciNet 
MATH 
Article 

Google Scholar 
109.Kulldorff, M. Information Management Services, Inc. SaTScanTM v.9.4.4: software for the spatial and space-time scan statistics. http://www.satscan.org (2009). More

Photosynthetic characteristicsWater and nitrogen coupling treatment had a significant effect on the photosynthetic characteristics (Fig. 1). Generally, the net photosynthetic rates of the treatments were in the following order: CK, W1N1, W1N3, W3N1, W3N3, W2N1, W1N2, W3N2, W2N3, and W2N2. The treatments with low water and low nitrogen had significantly lower net photosynthetic rates than W2N2. The stomatal conductance and transpiration rate changed in similar patterns. The net photosynthetic rate showed a unimodal trend with the increase of nitrogen application at the same irrigation level. Under the same nitrogen application level, the net photosynthetic rate increased first and then decreased slowly with the increase of irrigation amount, with the highest photosynthetic rates in the order of W2  > W3  > W1. The net photosynthetic rate was the highest, with a mean value was 13.87 μmol m−2 s−1, in treatment W2N2. The results showed that severe water stress and excessive nitrogen were not conducive to the absorption and utilization of water and nutrients by crop roots, which led to the decrease of the photosynthetic rate. The effect of water and nitrogen treatment on the intercellular CO2 concentration was significant (Fig. 1). Under the condition of excessive water or nitrogen, the photosynthesis of Isatis indigotica decreased, and the intercellular CO2 concentration showed a trend opposite to that of the net photosynthetic rate.Compared with N, P, and K deficiency treatments, water–N coupling could increase the Pn of crops, which was the same as that of other fruit trees and vegetables13. Accumulated photoassimilates in the third internode of the upper part of the main stems, as well as in the flag leaf sheath, are mobilized in a higher proportion and can contribute to grain filling in rice plants subjected to water stress in the tillering phase14. The Pn, Gs, and Tr of maize leaves at the seedling stage decreased significantly, while the Ci increased significantly when the nitrogen application rate was low15.The experiments with Isatis indigotica demonstrate that the Pn, Gs, and Tr under the same irrigation level first increased and then decreased with the increase of the nitrogen application rate. The net photosynthetic rate, transpiration rate, and stomatal conductance of Isatis indigotica were improved by rational nitrogen application. Studies have reported similar findings in Isatis indigotica; with the decrease of N level, the net photosynthetic rate, transpiration rate, and stomatal conductance of leaves gradually decreased, while the intercellular CO2 concentration increased16,17. Under reasonable water and nitrogen coordination conditions, the synergistic effect of water and nitrogen increased, which effectively promoted the photosynthesis of Isatis indigotica. Under the condition of too much nitrogen or too little water, the antagonism of water and nitrogen was obvious, and the photosynthesis of Isatis indigotica was inhibited to a certain extent.Yield and water use efficiencyThe Isatis indigotica yield values presented are the average of two consecutive years of water–nitrogen trials (Fig. 2). The I. indigotica yields differed significantly between the water–nitrogen treatments; the W2N2 and W2N3 treatments had the highest yields at 7277.5 and 6820.5 kg hm−2, respectively. The lowest yield of 3264.5 kg hm−2 was recorded in the control treatment. The yields of all treatments were significantly higher than that of the control treatment. The yields of the W2N2 and W2N3 treatments were significantly higher than those of the W1N1 and the W3N1 treatments. With the increase of the nitrogen application rate, the yield first increased and then decreased under the same irrigation conditions.The water use efficiency values of Isatis indigotica presented are the average of 2 consecutive years of water–nitrogen trials (Fig. 2). The water use efficiency of Isatis indigotica differed significantly between the water–nitrogen treatments; the W1N2 and W2N2 treatments had the highest yields at 20.78 and 19.63 kg mm−1 hm−2, respectively. The lowest yield of 13.65 kg mm−1 hm−2 was recorded in the W3N1 treatment. The water use efficiency values of the W1N2 and W2N2 treatments were significantly higher than that of the W3N3 treatment, which was the treatment with excess water and nitrogen fertilizer. The water use efficiency decreased with the increase of irrigation under the same nitrogen application conditions. The water use efficiency first increased and then decreased with the increase in nitrogen application rate under the same irrigation conditions. The W2N2 treatment had the highest yield and water use efficiency. Therefore, the water–nitrogen coupling mode of medium water and medium nitrogen application achieved the highest yield and effectively saved water. This was mainly due to the moderate water and nitrogen to promote the photosynthesis of Isatis indigotica and lead to more dry matter accumulation, so as to increase the yield.Generally, appropriate water deficits can improve crop yield and water use efficiency18,19, and rational fertilization can increase crop yield, such as in fruit trees and vegetables20,21,22. The yield increase in the current experiment was probably related to reasonable water stress and reasonable nitrogen application; the W2N2 treatment had the highest yield and water use efficiency. However, excessive water and nitrogen reduced the yield and water use efficiency of Isatis indigotica. This was consistent with recent research reports23,24. Compared with the local flooding irrigation and excessive nitrogen fertilizer mode, the W2N2 treatment with moderate water and nitrogen application not only obtained a high yield but also significantly improved the water use efficiency. This method could reduce the effect of excessive water and fertilizer application on soil productivity and would be a better water and nitrogen management model for local Isatis indigotica production.QualityThe Isatis indigotica quality values presented are the average of two consecutive years of water–nitrogen trials (Fig. 3). These quality indicators mainly include the following content indicators: indigo, indirubin, (R, S)-goitrin, and polysaccharides. The Isatis indigotica quality indicators differed significantly between the water–nitrogen treatments. The CK treatment had the highest values of all quality indicators. Each quality indicator decreased gradually with the increase of water content under the same nitrogen application conditions. Each quality indicator decreased gradually with the increase of nitrogen application under the same water conditions. The (R, S)-goitrin content of the W2N2 treatment decreased by 6.5% compared with CK and by 3.9% compared with the W1N1 treatment.Water is the medium for improving crop quality. Generally, the crop quality was improved by a suitable water deficit25,26,27 and reasonable fertilization28,29,30. The quality of Isatis indigotica in the current experiment increased gradually with the decrease of water. The water deficit treatment increased the content of effective components and improved the quality of Isatis indigotica. The content of the effective components in all treatments reached the pharmacopoeia standard12. The quality indicator values of each treatment in the current experiment were significantly lower than those of the CK treatment, but there was little difference in the quality indicator values between each treatment. Moreover, the yield of the control treatment was much lower than that of other treatments. Therefore, the effective quality content of the control treatment was lower than other treatments. Excessive water and nitrogen inputs were not conducive to quality, which was not consistent with recent research reports31. Generally, the water-nitrogen coupling type of W2N60 was antagonism basing on the average yield of winter wheat in the 10 years32. Some scholars have studied the irrigation of jujube that WUE and ANUE of jujube cannot reach the maximum at the same time. Different ratio of water and nitrogen will produce coupling and antagonism33. The results showed that total N applications over 200 kg ha−1 did not increase yield or quality. Water deficit treatment could be increased the content of effective components and improve the quality of Isatis indigotica. Due to the high evaporation, moderate water stress and effective use of nitrogen fertilizer, the active components of Isatis indigotica were easier to accumulate in its roots. The synergistic effect of water and nitrogen will lead to the accumulation of active components in Isatis indigotica. More