Ecology

1.Millennium Ecosystem Assessment. Ecosystems and Human Well-Being Vol. 5, 563 (Island Press, 2005).
Google Scholar 
2.Parker, I. M. et al. Impact: Toward a framework for understanding the ecological effects of invaders. Biol. Invas. 1(1), 3–19 (1999).Article 

Google Scholar 
3.Crowl, T. A., Crist, T. O., Parmenter, R. R., Belovsky, G. & Lugo, A. E. The spread of invasive species and infectious disease as drivers of ecosystem change. Front. Ecol. Environ. 6(5), 238–246 (2008).Article 

Google Scholar 
4.Mazza, G., Tricarico, E., Genovesi, P. & Gherardi, F. Biological invaders are threats to human health: An overview. Ethol. Ecol. Evol. 26, 112–129 (2014).Article 

Google Scholar 
5.Pimentel, D., Zuniga, R. & Morrison, D. Update on the environmental and economic costs associated with alien invasive species in the United States. Ecol. Econ 52, 273–288 (2005).Article 

Google Scholar 
6.Vilà, M. et al. How well do we understand the impacts of alien species on ecosystem services? A pan-European, cross-taxa assessment. Front. Ecol. Environ. 8(3), 135–144 (2010).Article 

Google Scholar 
7.Lindenmayer, D. B. & Likens, G. E. Adaptive monitoring: A new paradigm for long-term research and monitoring. Trends Ecol. Evol. 24(9), 482–486 (2009).PubMed 
Article 
PubMed Central 

Google Scholar 
8.Lawson Handley, L.-J. et al. Ecological genetics of invasive alien species. Biocontrol 56, 409–428 (2011).Article 

Google Scholar 
9.Fischer, M. L. et al. Multiple founder effects are followed by range expansion and admixture during the invasion process of the raccoon (Procyon lotor) in Europe. Divers. Distrib. 23(4), 409–420 (2017).Article 

Google Scholar 
10.Salgado, I. Is the raccoon (Procyon lotor) out of control in Europe?. Biodivers. Conserv. 27, 2243–2256 (2018).Article 

Google Scholar 
11.Asada, M. “Lag-phase management” as a population management method in low density areas in sika deer (Cervus nippon) and racoon (Procyon lotor). Honyurui Kagaku (Mamm. Sci.) 53(2), 243–255 (2013).
Google Scholar 
12.Gehrt, S. D. & Fritzell, E. K. Duration of familial bonds and dispersal patterns for raccoons in south Texas. J. Mammal. 79(3), 859–872. https://doi.org/10.2307/1383094 (1998).Article 

Google Scholar 
13.Gascoigne, J., Berec, L., Gregory, S. & Courchamp, F. Dangerously few liaisons: A review of mate-finding Allee effects. Popul. Ecol. 51(3), 355–372 (2009).Article 

Google Scholar 
14.Porter, W. F., Mathews, N. E., Underwood, H. B., Sage, R. W. & Behrend, D. F. Social organization in deer: Implications for localized management. Environ. Manage. 15(6), 809–814 (1991).Article 
ADS 

Google Scholar 
15.Long, J. Introduced Mammals of the World: Their History Distribution and Influence (CSIRO Publishing, 2003).Book 

Google Scholar 
16.Gehrt, S. D. Wild mammals of North America. In Raccoons and Allies (eds Feldhamer, G. A. et al.) 611–634 (CABI, 2003).
Google Scholar 
17.Ikeda, T., Asano, M., Matoba, Y. & Abe, G. Present status of invasive alien raccoon and its impact in Japan. Glob. Environ. Res. 8, 125–131 (2004).
Google Scholar 
18.Ministry of the Environment Government of Japan. The Birds and Beasts (Bears) Need Care Habitation Distribution Survey in 2017. Survey Report, Raccoon, Palm civet, Nutria (2018).19.Okuyama, M. W. et al. Genetic population structure of invasive raccoons (Procyon lotor) in Hokkaido, Japan: Unique phenomenon caused by pet escape or abandonment. Sci. Rep. 10(1), 1–10 (2020).Article 
CAS 

Google Scholar 
20.Ochiai, K., Ishii, M. & Furukawa, T. Invasion and distribution of the raccoon, Procyon lotor, in Chiba Prefecture, Central Japan. J. Nat. Hist. Museum Inst. 7, 21–27 (2002).
Google Scholar 
21.Asada, M. Bayesian estimation of population size in raccoon (Procyon lotor) using state-space model based on removal sampling. Honyurui Kagaku (Mamm. Sci.) 54, 207–218 (2014).
Google Scholar 
22.Sugai, T., Matsushima, H. & Mizuno, K. Last 400 ka landform evolution of the Kanto Plain: Under the influence of concurrent glacio-eustatic sea level changes and tectonic activity. J. Geogr. Chigaku Zasshi 122(6), 921–948 (2013).CAS 
Article 

Google Scholar 
23.Bagan, H. & Yamagata, Y. Landsat analysis of urban growth: How Tokyo became the world’s largest megacity during the last 40 years. Remote Sens. Environ. 127, 210–222 (2012).Article 
ADS 

Google Scholar 
24.Japan Meteorology Agency. Search Past Weather Data (2021). http://www.data.jma.go.jp/obd/stats/etrn/view/monthly_h1.php?prec_no=45&block_no=00&year=2020&month=&day=&view=p1. Accessed 28 Feb 2021.25.Geospatial Information Authority of Japan. https://maps.gsi.go.jp/. Accessed 4 Apr 2021.26.Frantz, A. C. et al. Limited mitochondrial DNA diversity is indicative of a small number of founders of the German raccoon (Procyon lotor) population. Eur. J. Wildl. Res. 59, 665–674 (2013).Article 

Google Scholar 
27.Cullingham, C. I., Kyle, C. J., Pond, B. A. & White, B. N. Genetic structure of raccoons in eastern North America based on mtDNA: Implications for subspecies designation and rabies disease dynamics. Can. J. Zool. 86, 947–958 (2008).CAS 
Article 

Google Scholar 
28.Kumar, S., Stecher, G., Li, M., Knyaz, C. & Tamura, K. MEGA X: Molecular evolutionary genetics analysis across computing platforms. Mol. Biol. Evol. 35(6), 1547–1549 (2018).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
29.Cullingham, C. I., Kyle, C. J. & White, B. N. Isolation, characterization and multiplex genotyping of raccoon tetranucleotide microsatellite loci. Mol. Ecol. Notes 6(4), 1030–1032 (2006).CAS 
Article 

Google Scholar 
30.Fike, J. A., Drauch, A. M., Beasley, J. C., Dharmarajan, G. & Rhodes, O. E. Development of 14 multiplexed microsatellite loci for raccoons Procyon lotor: Primer note. Mol. Ecol. Notes 7(3), 525–527 (2007).CAS 
Article 

Google Scholar 
31.Siripunkaw, C. et al. Isolation and characterization of polymorphic microsatellite loci in the raccoon (Procyon lotor). Mol. Ecol. Resour. 8(1), 199–201 (2008).CAS 
PubMed 
Article 

Google Scholar 
32.Koressaar, T. & Remm, M. Enhancements and modifications of primer design program Primer3. Bioinformatics 23(10), 1289–1291 (2007).CAS 
PubMed 
Article 

Google Scholar 
33.Untergasser, A. et al. Primer3-new capabilities and interfaces. Nucleic Acids Res. 40(15), 1–12 (2012).Article 
CAS 

Google Scholar 
34.Brownstein, M. J., Carpten, J. D. & Smith, J. R. Modulation of non-templated nucleotide addition by Taq DNA polymerase: Primer modifications that facilitate genotyping. Biotechniques 20(6), 1004–1010 (1996).CAS 
PubMed 
Article 

Google Scholar 
35.Pompanon, F., Bonin, A., Bellemain, E. & Taberlet, P. Genotyping errors: Causes, consequences and solutions. Nat. Rev. Genet. 6(11), 847–859 (2005).CAS 
PubMed 
Article 

Google Scholar 
36.Holleley, C. E. & Geerts, P. G. Multiplex Manager 1.0: A cross-platform computer program that plans and optimizes multiplex PCR. Biotechniques 46(7), 511–517 (2009).CAS 
PubMed 
Article 

Google Scholar 
37.Yoshida, K., Hirose, M., Hasegawa, M. & Inoue, E. Mitochondrial DNA analyses of invasive raccoons (Procyon lotor) in the Boso Peninsula, Japan. Mamm. Study 45(1), 1–6 (2020).Article 

Google Scholar 
38.Bandelt, H. J., Forster, P. & Rohl, A. Median-joining networks for inferring intraspecific phylogenies. Mol. Biol. Evol. 16, 37–48 (1999).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
39.Cullingham, A., Zalewski, A., Bartoszewicz, M., Okarma, H. & Jędrzejewska, E. The genetic structure of raccoon introduced in Central Europe reflects multiple invasion pathways. Biol. Invas. 16, 1611–1625 (2014).Article 

Google Scholar 
40.Fischer, M. L. et al. Historical invasion records can be misleading: genetic evidence for multiple introductions of invasive raccoons (Procyon lotor) in Germany. PLoS ONE 10(5), e0125441 (2015).PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
41.Alda, F. et al. Genetic evidence for multiple introduction events of raccoons (Procyon lotor) in Spain. Biol. Invas. 15, 687–698 (2013).Article 

Google Scholar 
42.Biedrzycka, A., Zalewski, A., Bartoszewicz, M., Okarma, H. & Jędrzejewska, E. The genetic structure of raccoon introduced in Central Europe reflects multiple invasion pathways. Biol. Invas. 16(8), 1611–1625 (2014).Article 

Google Scholar 
43.Santonastaso, T. T., Dubach, J., Hauver, S. A., Graser, W. H. & Gehrt, S. D. Microsatellite analysis of raccoon (Procyon lotor) population structure across an extensive metropolitan landscape. J. Mammal. 93(2), 447–455 (2012).Article 

Google Scholar 
44.Peakall, R. O. D. & Smouse, P. E. GENALEX 6: Genetic analysis in Excel. Population genetic software for teaching and research. Mol. Ecol. Notes 6(1), 288–295 (2006).Article 

Google Scholar 
45.Peakall, R. & Smouse, P. E. GenALEx 6.5: Genetic analysis in Excel. Population genetic software for teaching and research-an update. Bioinformatics 28(19), 2537–2539 (2012).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
46.Goudet, J. FSTAT, a Program to Estimate and Test Gene Diversity and Fixation Indices (Version 2.9. 3). http://www2.unil.ch/popgen/softwares/fstat.htm (2001).47.Pritchard, J. K., Stephens, M. & Donnelly, P. Inference of population structure using multilocus genotype data. Genetics 155(2), 945–959 (2000).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
48.Earl, D. A. & von Holdt, B. M. Structure harvester: A website and program for visualizing structure output and implementing the Evanno method. Conserv. Genet. Resour. 4(2), 359–361 (2012).Article 

Google Scholar 
49.Jakobsson, M. & Rosenberg, N. A. CLUMPP: A cluster matching and permutation program for dealing with label switching and multimodality in analysis of population structure. Bioinformatics 23(14), 1801–1806. https://doi.org/10.1093/bioinformatics/btm233 (2007).CAS 
Article 
PubMed 
PubMed Central 

Google Scholar 
50.Rosenberg, N. A. Distruct: A program for the graphical display of population structure. Mol. Ecol. Notes 4(1), 137–138 (2004).Article 

Google Scholar 
51.Ratnayeke, S., Tuskan, G. A. & Pelton, M. R. Genetic relatedness and female spatial organization in a solitary carnivore, the raccoon, Procyon lotor. Mol. Ecol. 11(6), 1115–1124 (2002).CAS 
PubMed 
Article 
PubMed Central 

Google Scholar 
52.Abdelkrim, J., Pascal, M., Calmet, C. & Samadi, S. Importance of assessing population genetic structure before eradication of invasive species: Examples from insular Norway rat populations. Conserv. Biol. 19(5), 1509–1518 (2005).Article 

Google Scholar  More

The present study was carried out in six fjords within New Zealand’s Fiordland system, specifically Breaksea Sound, Chalky Inlet, Doubtful Sound, Dusky Sound, Long Sound, and Wet Jacket Arm, as described in Tobias-Hünefeldt et al.15. Analyses were divided into three categories: (1) a multi-fjord analysis comprising five of the tested fjords (excluding Long Sound), (2) a high-resolution study along Long Sound’s horizontal axis, and (3) a depth profile of Long Sound’s deepest location (at 421 m). These categories were established to identify trends across multiple fjords, and then test the trends using a fjord analysed at a higher resolution. Total community composition (via 16S and 18S rRNA gene sequencing) and metabolic potential did not significantly covary across the five studied fjords (Mantel, r  More

1.Stein, R. & Macdonald, R. W. Organic carbon budget: Arctic Ocean vs. global ocean. In The Organic Carbon Cycle in the Arctic Ocean (eds Stein, R. & Macdonald, R. W.) (Springer, 2004).Chapter 

Google Scholar 
2.Barber, D. G. & Massom, R. A. The role of sea ice in Arctic and Antarctic polynyas. Oceanogr. Ser. 74, 1–54. https://doi.org/10.1016/S0422-9894(06)74001-6 (2007).Article 

Google Scholar 
3.Sheremetevskiy, A. M. Role of meiobenthos of the South Sakhalin shelf, Eastern Kamchatka, and Novosibirsk shallow water area. Issledovaniya Fauny Morei 35, 43 (1987).
Google Scholar 
4.Golikov, A. N. Ecosystems of the New Siberian shoals and fauna of the Laptev Sea and adjacent waters of the Arctic Ocean (in Russian). Explor. Fauna Seas 37, 4 (1990).
Google Scholar 
5.Golikov, A. N. Fauna of the East Siberian Sea. Part III. Explor. Fauna Seas 49, 57 (1994).
Google Scholar 
6.Sirenko, B. I. & Denisenko, S. G. Fauna of the East Siberian Sea, distribution patterns and structure of bottom communities. Explor. Fauna Seas 66, 74 (2010).
Google Scholar 
7.Sirenko, B. I. List of species of free-living invertebrates of Eurasian Arctic seas and adjacent deep waters. Explor. Fauna Seas 51(59), 1–76 (2001).
Google Scholar 
8.Schmidt-Rhaesa, A. Handbook of Zoology: Gastrotricha, Cycloneuralia, Gnathifera Vol. 2, 608 (De Gruyter, 2020).
Google Scholar 
9.Udalov, A. et al. Integrity of benthic assemblages along the arctic estuarine-coastal system. Ecol. Indic. 121, 107115. https://doi.org/10.1016/j.ecolind.2020.107115 (2021).Article 

Google Scholar 
10.Portnova, D., Fedyaeva, M., Udalov, A. & Tchesunov, A. Community structure of nematodes in the Laptev Sea shelf with notes on the lives of ice nematodes. Reg. Stud. Mar. Sci. 31, 100757. https://doi.org/10.1016/j.rsma.2019.100757 (2019).Article 

Google Scholar 
11.Gallucci, F., Moens, T. & Fonseca, G. Small-scale spatial patterns of meiobenthos in the Arctic deep sea. Mar. Biodivers. 39(1), 9–25. https://doi.org/10.1007/s12526-009-0003-x (2009).Article 

Google Scholar 
12.Lei, Y., Stumm, K., Volkenborn, N., Wickham, S. A. & Berninger, U. G. Impact of Arenicola marina (Polychaeta) on the microbial assemblages and meiobenthos in a marine intertidal flat. Mar. Biol. 157(6), 1271–1282. https://doi.org/10.1007/s00227-010-1407-7 (2010).Article 

Google Scholar 
13.Flint, M. V., Poyarkov, S. G. & Rymsky-Korsakov, N. A. Ecosystems of the Siberian Arctic Seas-2017 (Cruise 69 of the R/V Akademik Mstislav Keldysh). Oceanology 58(2), 315–318. https://doi.org/10.1134/S0001437018020042 (2018).ADS 
Article 

Google Scholar 
14.Garlitska, L. A. & Azovsky, A. I. Benthic harpacticoid copepods of the Yenisei Gulf and the adjacent shallow waters of the Kara Sea. J. Nat. Hist. 50, 2941–2959. https://doi.org/10.1080/00222933.2016.1219410 (2016).Article 

Google Scholar 
15.Portnova, D., Garlitska, L., Udalov, A. & Kondar, D. Meiobenthos and nematode community in the Yenisei Bay and adjacent parts of the Kara Sea shelf. Oceanology 57(1), 1–15. https://doi.org/10.1134/S0001437017010155 (2017).Article 

Google Scholar 
16.Carmack, E. et al. Toward quantifying the increasing role of oceanic heat in sea ice loss in the new Arctic. Bull. Am. Meteorol. Soc. 96(12), 2079–2105. https://doi.org/10.1175/BAMS-D-13-00177.1 (2005).ADS 
Article 

Google Scholar 
17.Peterson, B. J. et al. Increasing river discharge to the Arctic Ocean. Science 298(5601), 2171–2173. https://doi.org/10.1126/science.1077445 (2002).ADS 
CAS 
Article 
PubMed 

Google Scholar 
18.Polukhin, A. The role of river runoff in the Kara Sea surface layer acidification and carbonate system changes. ERL 14(10), 105007. https://doi.org/10.1088/1748-9326/ab421e (2019).ADS 
CAS 
Article 

Google Scholar 
19.Lisitzin, A. P. Marginal filter of the oceans. Oceanology 34(5), 735–743 (1994).CAS 

Google Scholar 
20.Moens, T., Braeckman, U., Derycke, S., Fonseca, G., Gallucci, F., Gingold, R., Guilini, Katja, Ingles, J., Leduc, D., Vanaverbeke, J., Van Colen, C., Vanreusel, A, & Vincx, M. Ecology of free-living marine nematodes. In Volume 2 Nematoda, 109–152. De Gruyter (2013)21.Aller, J. Y. & Aller, R. C. General characteristics of benthic faunas on the Amazon inner continental shelf with comparison to the shelf off the Changjiang River, East China Sea. Cont. Shelf Res. 6(1–2), 291–310. https://doi.org/10.1016/0278-4343(86)90065-8 (1986).ADS 
Article 

Google Scholar 
22.Soetaert, K., Vincx, M., Wittoeck, J. & Tulkens, M. Meiobenthic distribution and nematode community structure in five European estuaries. Hydrobiologia 311(1), 185–206. https://doi.org/10.1007/BF00008580 (1995).Article 

Google Scholar 
23.Tank, S. E. et al. The processing and impact of dissolved riverine nitrogen in the Arctic Ocean. Estuaries Coast 35, 401–415. https://doi.org/10.1007/s12237-011-9417-3 (2012).CAS 
Article 

Google Scholar 
24.Galtsova, V. V., Lukina, T. G. & Vladimirov, M. V. Meiobenthos of Chaunskaya Bay, East Siberian Sea. Issledovaniya Fauny Morei 48(56), 67–97 (1994).
Google Scholar 
25.Coull, B. C. Role of meiofauna in estuarine soft‐bottom habitats. Austral Ecol. 24(4), 327–343 (1999).Article 

Google Scholar 
26.Vincx, M., Meire, P., & Heip, C. The distribution of nematodes communities in the Southern Bight of the North Sea. Cah Biol Mar. 31(1), 107–129 (1990).27.Vanaverbeke, J., Gheskiere, T., Steyaert, M., & Vincx, M. Nematode assemblages from subtidal sandbanks in the Southern Bight of the North Sea: effect of small sedimentological differences. J. Sea Res. 48(3), 197–207. https://doi.org/10.1016/S1385-1101(02)00165-X (2002)ADS 
Article 

Google Scholar 
28.Steyaert, M., et al. The importance of fine-scale, vertical profiles in characterising nematode community structure. Estuar Coast Shelf Sci. 58(2), 353–366 (2003).ADS 
Article 

Google Scholar 
29.Alves, A. S., Adão, H., Patrício, J., Neto, J. M., Costa, M. J., & Marques, J. C. Spatial distribution of subtidal meiobenthos along estuarine gradients in two southern European estuaries (Portugal). J. Mar. Biol. Assoc. U. K. 89(8), 1529–1540 (2009).CAS 
Article 

Google Scholar 
30.Garlitska, L. A., Chertoprud, E. S., Portnova, D. A. & Azovsky, A. I. Benthic harpacticoida of the Kara Sea: Species composition and bathymetrically related distribution. Oceanology 59(4), 541–551. https://doi.org/10.1134/S0001437019040064 (2019).ADS 
CAS 
Article 

Google Scholar 
31.Huang, D. et al. Preliminary study on community structures of meiofauna in the middle and eastern Chukchi Sea. Acta Oceanol. Sin. 40(6), 83–91. https://doi.org/10.1007/s13131-021-1777-3 (2021).ADS 
Article 

Google Scholar 
32.Giere, O. Meiobenthology: The Microscopic Motile Fauna in Aquatic Sediments 2nd edn. (Springer, 2009).
Google Scholar 
33.Semiletov, I. et al. The East Siberian Sea as a transition zone between Pacific-derived waters and Arctic shelf waters. Geophys. Res. Lett. https://doi.org/10.1029/2005GL022490 (2005).Article 

Google Scholar 
34.Miroshnikov, A. Y. et al. Ecological state and mineral-geochemical characteristics of the bottom sediments of the East Siberian Sea. Oceanology 60(4), 595–610. https://doi.org/10.31857/S0030157420040152 (2020).Article 

Google Scholar 
35.Frontalini, F. et al. The response of cultured meiofaunal and benthic foraminiferal communities to lead exposure: Results from mesocosm experiments. Environ. Toxicol. Chem. 37(9), 2439–2447. https://doi.org/10.1002/etc.4207 (2018).CAS 
Article 
PubMed 

Google Scholar 
36.Fonseca, G. & Soltwedel, T. Deep-sea meiobenthic communities underneath the marginal ice zone off Eastern Greenland. Polar Biol. 30, 607–618. https://doi.org/10.1007/s00300-006-0220-8 (2007).Article 

Google Scholar 
37.Portnova, D. & Polukhin, A. Meiobenthos of the eastern shelf of the Kara Sea compared with the meiobenthos of other parts of the sea. Reg. Stud. Mar. Sci. 24, 370–378. https://doi.org/10.1016/j.rsma.2018.10.002 (2018).Article 

Google Scholar 
38.Alexeev, D. K., & Galtsova, V. V. Effect of radioactive pollution on the biodiversity of marine benthic ecosystems of the Russian Arctic shelf. Polar Sci. 6(2), 183–195 (2012).ADS 
Article 

Google Scholar 
39.Grzelak, K. & Sørensen, M. V. Diversity and community structure of kinorhynchs around Svalbard: First insights into spatial patterns and environmental drivers. Zool. Anz. 282, 31–43. https://doi.org/10.1016/j.jcz.2019.05.009 (2019).Article 

Google Scholar 
40.Landers, S. C. et al. Kinorhynch communities from Alabama coastal waters. Mar. Biol. Res. 16(6–7), 494–504. https://doi.org/10.1080/17451000.2020.1789660 (2020).Article 

Google Scholar 
41.Holovachov, O. New and known species of the genus Campylaimus Cobb, 1920 (Nematoda: Araeolaimida: Diplopeltidae) from North European marine habitats. Biodivers. Data J. https://doi.org/10.3897/BDJ.7.e46545 (2007).Article 

Google Scholar 
42.Sharma, J. & Bluhm, B. A. Diversity of larger free-living nematodes from macrobenthos ( > 250 μm) in the Arctic deep-sea Canada Basin. Mar. Biodivers. 41(3), 455–465. https://doi.org/10.1007/s12526-010-0060-1 (2010).Article 

Google Scholar 
43.Kotwicki, L., Grzelak, K. & Bełdowski, J. Benthic communities in chemical munitions dumping site areas within the Baltic deeps with special focus on nematodes. Deep Sea Res. II 128, 123–130. https://doi.org/10.1016/j.dsr2.2015.12.012 (2016).CAS 
Article 

Google Scholar 
44.Netto, S. A., Pagliosa, P. R., Colling, A., Fonseca, A. L. & Brauk, K. M. Benthic estuarine assemblages from the Southern Brazilian marine ecoregion. Braz. Estuaries. https://doi.org/10.1007/978-3-319-77779-5_6 (2018).Article 

Google Scholar 
45.Broman, E., et al. Uncovering diversity and metabolic spectrum of animals in dead zone sediments. Commun. Biol. 3(1), 1–12 (2020).46.Zeppilli, D., et al. Characteristics of meiofauna in extreme marine ecosystems: a review. Mar. Biodiver. 48(1), 35–71 (2018).47.Pérez-García, J. A. et al. Nematode diversity of freshwater and anchialine caves of Western Cuba. PBSW 131(1), 144–155. https://doi.org/10.2988/17-00024 (2018).Article 

Google Scholar 
48.Bezzubova, E. M., Seliverstova, A. M., Zamyatin, I. A. & Romanova, N. D. Heterotrophic bacterioplankton of the Laptev and East Siberian Sea shelf affected by freshwater inflow areas. Oceanology 60, 62–73. https://doi.org/10.1134/S0001437020010026 (2020).ADS 
CAS 
Article 

Google Scholar 
49.Vanreusel, A. et al. Meiobenthos of the central Arctic Ocean with special emphasis on the nematode community structure. Deep Sea Res. I 47, 1855–1879. https://doi.org/10.1016/S0967-063728002900007-8 (2000).Article 

Google Scholar 
50.Tahseen, Q. Nematodes in aquatic environments: Adaptations and survival strategies. Biodivers. J. 3(1), 13–40 (2012).
Google Scholar 
51.Williams, W. J. & Carmack, E. C. The ‘interior’ shelves of the Arctic Ocean: Physical oceanographic setting, climatology and effects of sea-ice retreat on cross-shelf exchange. Prog. Ocean 139, 24–41. https://doi.org/10.1016/j.pocean.2015.07.008 (2015).Article 

Google Scholar 
52.Magritsky, D. V. et al. Long-term changes of river water inflow into the seas of the Russian Arctic sector. Polarforschung 87(2), 177–194. https://doi.org/10.2312/polarforschung.87.2.177 (2018).Article 

Google Scholar 
53.Anderson, L. G. et al. East Siberian Sea, an Arctic region of very high biogeochemical activity. Biogeosciences 4, 6. https://doi.org/10.5194/bg-8-1745-2011 (2011).CAS 
Article 

Google Scholar 
54.Dmitrienko, I. A. et al. Impact of the Arctic Ocean Atlantic water layer on Siberian shelf hydrography. J. Geophys. Res. Oceans. https://doi.org/10.1029/2009JC006020 (2010).Article 

Google Scholar 
55.Stein, R. Arctic Ocean Sediments: Processes, PROXIES, and Paleoenvironment (Elsevier, 2008).
Google Scholar 
56.Petrova, V. I., Batova, G. I., Kursheva, A. V. & Litvinenko, I. V. Geochemistry of organic matter of bottom sediments in the rises of the central Arctic Ocean. Russ. Geol. Geophys. 51(1), 88–97. https://doi.org/10.1016/j.rgg.2009.12.008 (2010).ADS 
Article 

Google Scholar 
57.Millero, F. J. Thermodynamics of the carbon dioxide system in oceans. GCA 59(4), 661–677. https://doi.org/10.12691/wjce-3-6-1 (1995).ADS 
CAS 
Article 

Google Scholar 
58.Pavlova, G. Y. et al. Intercalibration of Bruevich’s method to determine the total alkalinity in seawater. Oceanology 48, 438. https://doi.org/10.1134/S0001437008030168 (2008).ADS 
Article 

Google Scholar 
59.Dickson, A. G. & Goyet, C. Handbook of Methods for the Analysis of the Various Parameters of the Carbon Dioxide System in Sea Water. Version 2 (No. ORNL/CDIAC-74) (1994).60.Dickson, A. G., Afghan, J. D. & Anderson, G. C. Reference materials for oceanic CO2 analysis: A method for the certification of total alkalinity. Mar. Chem. 80, 185–197. https://doi.org/10.1016/S0304-4203(02)00133-0 (2003).CAS 
Article 

Google Scholar 
61.Lewis, E. & Wallace, D. W. R. Program Developed for CO2 System Calculations. ORNL/CDIAC-105 (Carbon Dioxide Information Analysis Center, Oak Ridge National Laboratory, U.S. Department of Energy, 1998).Book 

Google Scholar 
62.Shiklomanov, A. I., Holmes, J. W., McClelland, S. E., Tank, R. & Spencer, G.M. Arctic Great Rivers Observatory. Discharge Dataset, Version 20200801 (2020).63.Niemistö, L. A gravity corer for studies of soft sediments. Merentutkimuslait. Julk./Havsforskningsinst. Skr. 238, 33–38 (1974).
Google Scholar 
64.Eleftheriou, A. Methods for the Study of Marine Benthos (Wiley, 2013).Book 

Google Scholar 
65.Wieser, W. Beziehungen zwischen Mundhöhlengestalt, Ernährungsweise und Vorkommen bei freilebenden, marinen Nematoden. Ark. Zool. 2, 439–484 (1953).
Google Scholar 
66.Hammer, Ø., Harper, D. A. T. & Ryan, P. D. PAST paleontological statistics software package for education and data analysis. Palaeontol. Electron. 4, 1–9 (2001).
Google Scholar 
67.Heip, C. & Herman, P. Indices of diversity and evenness. Oceanis 24(4), 61–88 (2001).
Google Scholar  More

1.Grime, J. P. Plant Strategies, Vegetation Processes, and Ecosystem Properties (John Wiley and Sons, 2001).2.Reich, P. B. et al. The evolution of plant functional variation: traits, spectra, and strategies. Int. J. Plant Sci. 164, S143–S164 (2003).Article 

Google Scholar 
3.Díaz, S. et al. The global spectrum of plant form and function. Nature 529, 167–171 (2016).ADS 
PubMed 
Article 
CAS 

Google Scholar 
4.Bergmann, J. et al. The fungal collaboration gradient dominates the root economics space in plants. Sci. Adv. 6, eaba3756 (2020).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
5.Wright, I. J. et al. The worldwide leaf economics spectrum. Nature 428, 821–827 (2004).ADS 
CAS 
PubMed 
Article 

Google Scholar 
6.Kattge, J. et al. TRY plant trait database — enhanced coverage and open access. Glob. Chang. Biol. 26, 119–188 (2020).ADS 
PubMed 
Article 

Google Scholar 
7.Iversen, C. M. et al. A global Fine-Root Ecology Database to address below-ground challenges in plant ecology. New Phytol. 215, 15–26 (2017).PubMed 
Article 

Google Scholar 
8.Guerrero-Ramírez, N. R. et al. Global root traits (GRooT) database. Glob. Ecol. Biogeogr. 30, 25–37 (2021).Article 

Google Scholar 
9.McCormack, M. L. et al. Redefining fine roots improves understanding of below-ground contributions to terrestrial biosphere processes. New Phytol. 207, 505–518 (2015).PubMed 
Article 

Google Scholar 
10.Rasse, D. P., Rumpel, C. & Dignac, M. F. Is soil carbon mostly root carbon? Mechanisms for a specific stabilisation. Plant Soil 269, 341–356 (2005).CAS 
Article 

Google Scholar 
11.Eissenstat, D. M. Costs and benefits of constructing roots of small diameter. J. Plant Nutr. 15, 763–782 (1992).Article 

Google Scholar 
12.Freschet, G. T., Cornelissen, J. H. C., van Logtestijn, R. S. P. & Aerts, R. Evidence of the ‘plant economics spectrum’ in a subarctic flora. J. Ecol. 98, 362–373 (2010).Article 

Google Scholar 
13.Reich, P. B. The world-wide ‘fast–slow’ plant economics spectrum: a traits manifesto. J. Ecol. 102, 275–301 (2014).Article 

Google Scholar 
14.Shen, Y. et al. Linking aboveground traits to root traits and local environment: implications of the plant economics spectrum. Front. Plant Sci. 10, 1412 (2019).Article 

Google Scholar 
15.Kramer-Walter, K. R. et al. Root traits are multidimensional: specific root length is independent from root tissue density and the plant economic spectrum. J. Ecol. 104, 1299–1310 (2016).Article 

Google Scholar 
16.Bergmann, J., Ryo, M., Prati, D., Hempel, S. & Rillig, M. C. Root traits are more than analogues of leaf traits: the case for diaspore mass. New Phytol. 216, 1130–1139 (2017).PubMed 
Article 

Google Scholar 
17.Weemstra, M. et al. Towards a multidimensional root trait framework: a tree root review. New Phytol. 211, 1159–1169 (2016).CAS 
PubMed 
Article 

Google Scholar 
18.Ma, Z. et al. Evolutionary history resolves global organization of root functional traits. Nature 555, 94–97 (2018).ADS 
CAS 
PubMed 
Article 

Google Scholar 
19.de la Riva, E. G. et al. Root traits across environmental gradients in Mediterranean woody communities: are they aligned along the root economics spectrum? Plant Soil 424, 35–48 (2018).Article 
CAS 

Google Scholar 
20.Craine, J. M., Lee, W. G., Bond, W. J., Williams, R. J. & Johnson, L. C. Environmental constraints on a global relationship among leaf and root traits of grasses. Ecology 86, 12–19 (2005).Article 

Google Scholar 
21.Liese, R., Alings, K. & Meier, I. C. Root branching is a leading root trait of the plant economics spectrum in temperate trees. Front. Plant Sci. 8, 315 (2017).PubMed 
PubMed Central 
Article 

Google Scholar 
22.Carmona, C. P. et al. Erosion of global functional diversity across the tree of life. Sci. Adv. 7, eabf2675 (2021).ADS 
PubMed 
PubMed Central 
Article 

Google Scholar 
23.Niklas, K. J. Modelling below- and above-ground biomass for non-woody and woody plants. Ann. Bot. 95, 315–321 (2005).PubMed 
Article 

Google Scholar 
24.Liu, G. et al. Coordinated variation in leaf and root traits across multiple spatial scales in Chinese semi-arid and arid ecosystems. New Phytol. 188, 543–553 (2010).PubMed 
Article 

Google Scholar 
25.Galland, T., Carmona, C. P., Götzenberger, L., Valencia, E. & de Bello, F. Are redundancy indices redundant? An evaluation based on parameterized simulations. Ecol. Indic. 116, 106488 (2020).Article 

Google Scholar 
26.Valverde‐Barrantes, O. J., Maherali, H., Baraloto, C. & Blackwood, C. B. Independent evolutionary changes in fine‐root traits among main clades during the diversification of seed plants. New Phytol. 228, 541–553 (2020).PubMed 
Article 

Google Scholar 
27.Cornwell, W. K. et al. Plant species traits are the predominant control on litter decomposition rates within biomes worldwide. Ecol. Lett. 11, 1065–1071 (2008).PubMed 
Article 

Google Scholar 
28.Freschet, G. T. et al. Climate, soil and plant functional types as drivers of global fine-root trait variation. J. Ecol. 105, 1182–1196 (2017).Article 

Google Scholar 
29.De Deyn, G. B. & Van der Putten, W. H. Linking aboveground and belowground diversity. Trends Ecol. Evol. 20, 625–633 (2005).PubMed 
Article 

Google Scholar 
30.Pausas, J. G. & Bond, W. J. Humboldt and the reinvention of nature. J. Ecol. 107, 1031–1037 (2019).Article 

Google Scholar 
31.Poorter, H. et al. Biomass allocation to leaves, stems and roots: meta-analyses of interspecific variation and environmental control. New Phytol. 193, 30–50 (2012).CAS 
PubMed 
Article 

Google Scholar 
32.Moora, M. Mycorrhizal traits and plant communities: perspectives for integration. J. Veg. Sci. 25, 1126–1132 (2014).Article 

Google Scholar 
33.Freschet, G. T. et al. Root traits as drivers of plant and ecosystem functioning: current understanding, pitfalls and future research needs. New Phytol. https://doi.org/10.1111/nph.17072 (2021).34.McCormack, M. L. & Iversen, C. M. Physical and functional constraints on viable belowground acquisition strategies. Front. Plant Sci. 10, 1215 (2019).PubMed 
PubMed Central 
Article 

Google Scholar 
35.Wells, C. E. & Eissenstat, D. M. Beyond the roots of young seedlings: the influence of age and order on fine root physiology. J. Plant Growth Regul. 21, 324–334 (2002).CAS 
Article 

Google Scholar 
36.Zanne, A. E. et al. Three keys to the radiation of angiosperms into freezing environments. Nature 506, 89–92 (2014).ADS 
CAS 
PubMed 
Article 

Google Scholar 
37.USDA. USDA PLANTS Database (accessed 3rd July 2020); https://plants.sc.egov.usda.gov38.Engemann, K. et al. A plant growth form dataset for the New World. Ecology 97, 3243 (2016).CAS 
PubMed 
Article 

Google Scholar 
39.BGCI. GlobalTreeSearch online database (accessed 3rd July 2020); https://www.bgci.org/globaltree_search.php40.The Plant List. The Plant List (accessed 17th February 2020); http://www.theplantlist.org41.Cayuela, L., Macarro, I., Stein, A. & Oksanen, J. Taxonstand: Taxonomic Standardization of Plant Species Names. R package version 2.2. https://CRAN.R-project.org/package=Taxonstand (2019).42.Stekhoven, D. J. & Buhlmann, P. MissForest—non-parametric missing value imputation for mixed-type data. Bioinformatics 28, 112–118 (2012).CAS 
PubMed 
Article 

Google Scholar 
43.Oliveira, B. F., Sheffers, B. R. & Costa, G. C. Decoupled erosion of amphibians’ phylogenetic and functional diversity due to extinction. Glob. Ecol. Biogeogr. 29, 309–319 (2020).Article 

Google Scholar 
44.Penone, C. et al. Imputation of missing data in life-history trait datasets: which approach performs the best? Methods Ecol. Evol. 5, 961–970 (2014).Article 

Google Scholar 
45.Jin, Y. & Qian, H. V.PhyloMaker: an R package that can generate very large phylogenies for vascular plants. Ecography 42, 1353–1359 (2019).Article 

Google Scholar 
46.Smith, S. A. & Brown, J. W. Constructing a broadly inclusive seed plant phylogeny. Am. J. Bot. 105, 302–314 (2018).PubMed 
Article 

Google Scholar 
47.Whittakker, R. H. Communities and Ecosystems (Macmillan, 1975).48.Stefan, V. & Levin, S. plotbiomes: Plot Whittaker biomes with ggplot2. R package version 0.0.0.9001 https://github.com/valentinitnelav/plotbiomes (2021).49.Ricklefs, R. E. The Economy of Nature (W. H. Freeman and Company, 2008).50.GBIF. GBIF Occurrence Download (accessed 15 December 2019); https://doi.org/10.15468/dl.thlxph51.South, A. rworldmap: a new R package for mapping global data. R J. 3, 35–43 (2011).Article 

Google Scholar 
52.Dinno, A. paran: Horn’s Test of Principal Components/Factors. R package version 1.5.2. https://CRAN.R-project.org/package=paran (2018).53.Dray, S. & Dufour, A.-B. The ade4 package: implementing the duality diagram for ecologists. J. Stat. Softw. https://doi.org/10.18637/jss.v022.i04 (2007).54.Duong, T. ks: kernel density estimation and kernel discriminant analysis for multivariate data in R. J. Stat. Softw. https://doi.org/10.18637/jss.v021.i07 (2015).55.Duong, T. ks: Kernel smoothing. R package version 1.11.5 https://CRAN.R-project.org/package=ks (2019).56.Carmona, C. P., Bello, F., Mason, N. W. H. & Lepš, J. Trait probability density (TPD): measuring functional diversity across scales based on TPD with R. Ecology 100, e02876 (2019).PubMed 
Article 

Google Scholar 
57.Carmona, C. P. TPD: methods for measuring functional diversity based on Trait Probability Density. R package version 1.1.0. https://CRAN.R-project.org/package=TPD (2019).58.Duong, T. & Hazelton, M. L. Plug-in bandwidth matrices for bivariate kernel density estimation. J. Nonparametr. Stat. 15, 17–30 (2003).MathSciNet 
MATH 
Article 

Google Scholar 
59.Carmona, C. P., de Bello, F., Mason, N. W. H. & Lepš, J. Traits without borders: integrating functional diversity across scales. Trends Ecol. Evol. 31, 382–394 (2016).PubMed 
Article 

Google Scholar 
60.Mason, N. W. H., Mouillot, D., Lee, W. G. & Wilson, J. B. Functional richness, functional evenness and functional divergence: the primary components of functional diversity. Oikos 111, 112–118 (2005).Article 

Google Scholar 
61.Villéger, S., Mason, N. W. H. & Mouillot, D. New multidimensional functional diversity indices for a multifaceted framework in functional ecology. Ecology 89, 2290–2301 (2008).PubMed 
Article 

Google Scholar 
62.Oksanen, J. et al. vegan: Community Ecology Package. R package version 2.5-5 https://CRAN.R-project.org/package=vegan (2019).63.Carmona, C. P. et al. Taxonomical and functional diversity turnover in Mediterranean grasslands: interactions between grazing, habitat type and rainfall. J. Appl. Ecol. 49, 1084–1093 (2012).Article 

Google Scholar 
64.Micó, E. et al. Contrasting functional structure of saproxylic beetle assemblages associated to different microhabitats. Sci. Rep. 10, 1520 (2020).ADS 
PubMed 
PubMed Central 
Article 
CAS 

Google Scholar 
65.Blonder, B. et al. New approaches for delineating n-dimensional hypervolumes. Methods Ecol. Evol. 9, 305–319 (2018).Article 

Google Scholar 
66.Carmona, C. P., de Bello, F., Mason, N. W. H. & Lepš, J. The density awakens: a reply to Blonder. Trends Ecol. Evol. 31, 667–669 (2016).PubMed 
Article 

Google Scholar 
67.Mouillot, D. et al. Niche overlap estimates based on quantitative functional traits: a new family of non-parametric indices. Oecologia 145, 345–353 (2005).ADS 
PubMed 
Article 

Google Scholar 
68.de Bello, F., Carmona, C. P., Mason, N. W. H., Sebastià, M.-T. & Lepš, J. Which trait dissimilarity for functional diversity: trait means or trait overlap? J. Veg. Sci. 24, 807–819 (2013).Article 

Google Scholar 
69.Traba, J., Iranzo, E. C., Carmona, C. P. & Malo, J. E. Realised niche changes in a native herbivore assemblage associated with the presence of livestock. Oikos 126, 1400–1409 (2017).Article 

Google Scholar 
70.Cornwell, W. K., Schwilk, D. W. & Ackerly, D. D. A trait-based test for habitat filtering: Convex Hull Volume. Ecology 87, 1465–1471 (2006).PubMed 
Article 

Google Scholar 
71.Blonder, B., Lamanna, C., Violle, C. & Enquist, B. J. The n-dimensional hypervolume. Glob. Ecol. Biogeogr. 23, 595–609 (2014).Article 

Google Scholar 
72.Blonder, B. Hypervolume concepts in niche- and trait-based ecology. Ecography 41, 1441–1455 (2018).Article 

Google Scholar 
73.Ricotta, C. et al. Measuring the functional redundancy of biological communities: a quantitative guide. Methods Ecol. Evol. 7, 1386–1395 (2016).Article 

Google Scholar 
74.Mouillot, D. et al. Functional over-redundancy and high functional vulnerability in global fish faunas on tropical reefs. Proc. Natl. Acad. Sci. USA 111, 13757–13762 (2014).ADS 
CAS 
PubMed 
PubMed Central 
Article 

Google Scholar 
75.Carmona, C. P., de Bello, F., Sasaki, T., Uchida, K. & Pärtel, M. Towards a common toolbox for rarity: a response to Violle et al. Trends Ecol. Evol. 32, 889–891 (2017).PubMed 
Article 

Google Scholar 
76.Violle, C. et al. Functional rarity: the ecology of outliers. Trends Ecol. Evol. 32, 356–367 (2017).PubMed 
PubMed Central 
Article 

Google Scholar 
77.Borcard, D., Legendre, P. & Drapeau, P. Partialling out the spatial component of ecological variation. Ecology 73, 1045–1055 (1992).Article 

Google Scholar 
78.Gower, J. C. General coefficient of similarity and some of its properties. Biometrics 27, 857–871 (1971).Article 

Google Scholar 
79.Carmona, C. P. et al. Agriculture intensification reduces plant taxonomic and functional diversity across European arable systems. Funct. Ecol. 34, 1448–1460 (2020).Article 

Google Scholar 
80.Gherardi, L. A. & Sala, O. E. Global patterns and climatic controls of belowground net carbon fixation. Proc. Natl Acad. Sci. USA 117, 20038–20043 (2020).CAS 
PubMed 
PubMed Central 
Article 

Google Scholar  More

This is an observational study with no intervention on flock and hatchery practices. None of the birds or eggs were exposed to experimental procedures. The study was based mainly on existing data provided by the hatchery company (DanHatch Denmark A/S) from five broiler breeder flocks in Denmark during the period from November 2018 to January 2019 when the breeders were 46 to 62 weeks of age, see details in Table 1. In addition, feed samples from the flock locations and eggs from grocery stores were acquired.Table 1 Flocks and production periods.Full size tableThe average age of breeders was 48–59 weeks (SD from 0.5 to 2.2) ranging from 46–50 weeks to 57–62 weeks (Table 1; Supplementary Fig. S1 online) with observation period ranging from 1.6 to 7.6 weeks in the five flocks. Average laying percent over observation days was 65% (SD = 5.4%) and average hatchability over deliveries was 79% (SD = 5.8%).Feed samplesTwenty-six feed samples were collected for analysis of glyphosate content, 3 to 10 feed samples per flock. The glyphosate concentration related to a given sampling date was assumed representative for the flock from this day and until next sampling. Average duration of the preceding samples were used as duration for the last sampling date within each flock. Glyphosate (N‐(phosphonomethyl) glycine) and the glyphosate degradation product, aminomethylphosphonic acid (AMPA) in feed samples were analysed by the method described by Nørskov et al.4.Production dataData on egg production and hatchability from periods following each feed sampling was obtained from the hatchery company. Daily information was available on laying percent (100% * number of eggs/number of breeders), breeder age (days) and egg weight. For the hatchability, this was calculated as the proportion of eggs placed in incubators from which a viable chicken hatched (but presented as a percentage, i.e. multiplied by 100%). Daily egg weight had been calculated as the average from approx. 30 randomly sampled eggs.Glyphosate concentration of the feed consumed by the breeders during the 10 days prior to laying was the explanatory variable of main interest. The weighted average of glyphosate concentrations across the 10 days of development from follicle to ovulation of egg was used with number of days each glyphosate sample is representative during these 10 days as weights. For hatchability, glyphosate concentrations were aggregated at the level of delivery by weighted averaging using number of hatch eggs as weights.Eggs from grocery storesNo eggs were obtained from the five flocks, however we acquired eight cartons of conventional as well as eight cartons of organic eggs from eight different grocery stores. Three eggs from each carton were selected and egg yolk were analysed for glyphosate by the microLC-MS/MS method as described by Nørskov et al.4 adjusted to the egg yolk matrix.Statistical analysisLaying percent and hatchability were analysed by linear mixed effects models, including a random effect of flock and a first order autoregressive correlation structure to account for the repeated measurements from each flock. Following two covariates were considered for both outcomes: average egg weight (g) and breeder age (decimal weeks). However, since egg weight and breeder age are highly correlated (Pearson’s correlation coefficient ranging from 0.73 to 0.95 in the five flocks; Supplementary Fig. S1 online), only breeder age was included in the models. An important reason for this choice being that average egg weight was missing for 24% and 43% of the days from flock 4 and 5, respectively. In the age range used for this study, laying percent decrease with breeder age (Supplementary Fig. S1 online) as substantiated by a correlation coefficient between − 0.38 and − 0.87. Hatchability also decrease with breeder age (Supplementary Fig. S1 online).In addition, storage time on farm until delivery (1 to 5 days) and storage time at hatchery until incubation starts (1 to 11 days) were included as covariates for hatchability. The incubation start date was determined as date of hatching minus 21 days. For hatchability, covariates obtained from flock production data were aggregated at the level of delivery by weighted averaging; using daily number of eggs as weights for the calculation of average egg weight, number of hatch eggs as weights for average storage time on farm, and current number of breeders as weights for average breeder age. Weighted average storage time on farm until delivery varied from 1.0 to 4.0 and was on average 2.1 days. For storage time at hatchery, deliveries had been split on one to four incubator start dates. Therefore, weighted average of storage days was calculated using number of delivered eggs as weights. Weighted average storage time at hatchery before incubation starts varied from 1.2 to 8.0 days and was on average 4.8 days.Final models were fitted with restricted maximum likelihood estimation using the lme function from the nlme package v. 3.1-152 in R version 4.0.45 and with a significance level of 0.05. Fixed effects were tested by χ2 likelihood ratio tests after maximum likelihood estimation. Model checking was carried out by examination of qq-plots for normality and scatter plots of residuals versus predicted values to look for uncovered trends and variance heterogeneity. More

Based on the extensive camera trap study, the very first information about the occupancy, trapping rate, activity pattern, group size, social structure and vital rates of the critically endangered Western Derby eland in its last refugium, the NKNP in Senegal, is presented here. The first estimation of abundance since 2006 is also provided7.Spatiotemporal behaviour pattern of WDE in the parkThe results of the CT survey in the NKNP highlight the substantially lower occupancy and trapping rate of WDE in comparison to other large ungulates in the park. According to the current results, the WDE occupied less than 5% of the park area during the dry season, being exclusively within the zone of Mont Assirik, and more specifically the Mansa Fara marsh, which can be thus designated as the core area of the WDE distribution. The trapping rate of the roan antelope, which is considered the most abundant antelope species in the NKNP, was 4.04, i.e. more than 11 times higher than that of WDE in the present study. Even the Western hartebeest, which is considered a rare species in the NKNP, had a trapping rate of 0.61, which is ca. twice as high that that of the WDE (see Rabeil et al.8 for further details, and additional ungulate species).In the zone of Mont Assirik, the trapping rate of WDE increased to 2.42, but the trapping rates of other antelope species remained higher still (4.6 for roan antelope and 3.39 for Western hartebeest8). The WDE distribution is therefore strongly localised within an area which seems to also be attractive for the other species, including the incidental records of elephant. The Mont Assirik zone, and more specifically the Mansa Fara marsh area, is therefore the crucial zone within the park for the WDE, and it appears to support a larger number of other antelope species as well. This zone should therefore be considered as a key conservation area, potentially very sensitive to targeted poaching, and thus crucial for efficacy of targeted law enforcement actions.When looking at the diurnal activity pattern, the WDE were active before midnight, approximately 3 h after sunset, in the morning, approximately 2 h after sunrise, and then again in the afternoon, with the peak activity during the hottest part of the day. This activity pattern is different from the typical bimodal activity pattern, which has peaks at dawn and dusk, as reported for most African grazing and browsing herbivores, seen as a behavioural thermoregulation strategy to avoid heat stress41,42,43. Instead, the WDE, being a large body-sized browsing antelope19,44, must stay active throughout the day to seek discretely distributed food, and fulfil foraging requirements by feeding while moving. The WDE appears to be well-adapted to tolerate such high temperatures, similar to kudu45, roan46, and giraffe47. Such behaviour pattern enable the law enforcement patrols, as well to tourists, to detect herds of WDEs and monitor them, thereby increasing their protection against poaching.Individual identification and recapturesThe individual identification of animals was more successful during the daytime, as the light conditions mostly did not allow for the proper visualization of the stripes on the flanks during the night captures (as similarly reported in Jůnek et al.18). When the ID is targeted to be successful during the night (as for the leopards and tigers), the camera traps are often set to the video mode to ensure a higher possibility of identification48. However, the activity of the WDE is not predominantly nocturnal, and the captures were distributed over both the daylight and night hours, and therefore the results are considered representative for the whole period.The AD animals were more likely to be identified in the present study because of their larger body size, resulting in better visibility of their stripes. The higher identification rate of larger individuals also likely contributed to the higher probability of recaptures, which were only recorded for individuals of 2Y and older.Overall, the identification success rate was comparable, maybe even slightly higher, than the previous camera trap study performed on the Eastern Giant eland in Chinko, CAR, specifically in the dataset from the dry season13, which corresponds to the observation period in the present study as well.In the NKNP, recaptures of individuals were recorded, whereas there were none reported in Chinko13. The recapture rate of the WDE in NKNP, with mostly short distances between the capture-recapture sites, even after the long-time gaps between the captures, confirm again that the WDE likely inhabit a relatively limited area of the park.Group size and social structureThe mean group size recorded in the NKNP during the present study was slightly larger than that within Chinko; however, the maximum group size was smaller in NKNP (32 vs. 41 individuals). Mixed herds were the largest in terms of the number of individuals, in both studies. The average group size has been reported as 20–30 individuals49, but Derby elands may form large herds of over 100 individuals in the late dry season14. Similarly, a large herd was reported within NKNP in 2006, having 69 individuals7, and a herd of around 60 WDE was also recently reported by patrols in 2020 (GIE Niokolo, personal communications). It is important to highlight that the results from the present study reflect the number of individuals per event based on visible individuals within the scope of the camera, and that the real group sizes may actually be larger.No adult males were present in the mixed herd in two cases within the present study; however, there were always 2YM and a few unidentified individuals, suggesting that the herd should not be considered as a pure “nursery herd”, as known for sexually dimorphic antelope species50.Calves are born in the NKNP during the period comparable to that of Bandia, Fathala and Chinko, i.e. during the early dry season16. The higher proportion of calves in the dry season corresponds with the nursing period of six months for WDE44. Given a pregnancy length of nine months, the WDE mating season in NKNP peaks in January/February, which also corresponds with the formation of large herds with multiple males, as similarly seen in Chinko and Cameroon13,14.Vital ratesThe sex ratio of the WDE in the NKNP was female-biased. The skewed adult sex ratio reflects the lower survival rate of males in comparison with females, typical for polygynous species51. This result also corresponds with the findings from other Derby eland populations, namely from Chinko, where the bias towards females in the adult sex ratio was even more pronounced (0.67:113). A similar ratio was found in the hunting reserves within Cameroon35, but also in the semi-captive population, without hunting and without predators34. As the ratio in NKNP was less skewed than that within Chinko and Cameroon, a lower or zero selectivity for males by hunters/poachers is expected.The population of WDE in the NKNP showed a lower proportion of adults versus other age categories compared to the demographic structure of the WDE in the semi-captive breeding facilities of the Bandia and Fathala reserve33,52, and to those of the Eastern subspecies of Derby eland in the Central African Republic13 (see Table 3). The data from the present study also showed a surprisingly high breeding rate (likely close to 100%), as well as a high survival rate of yearlings. This combination of demographic characteristics should be highly favourable, and likely to lead to a significant population growth rate; however, this does not seem to be the case of the WDE population in the NKNP (please refer to further discussion about population size).In this context, the population of WDE in the NKNP was explored deeper, to examine possible scenarios of changes within the population structure. The changes in vital rates between two years of monitoring (2017 and 2018) were examined, by taking advantage of the possible recognition of the age category until two years of age, and the knowledge of the life tables of the enclosed, non-predated WDE population in the Bandia reserve34. Life tables were created for each year, and for males (M) and females (F) separately, according to the standard structure2, and based on two scenarios: a) only the observed number of JUV and 1Y (nx), and modelled 2Y (model ‘JUV + 1Y’); b) the observed number of JUV and 2Y (nx) (model ‘JUV + 2Y’). Then, estimations of animals in age categories based on two parameters were calculated: (i) based on the mortality rate (qx) known from the Bandia reserve (Senegal), and (ii) based on the recorded number of animals (NAD), to calculate the estimation of mortality rate (for details, see Additional file 1: Table S2).The resulting values demonstrated that with survival rates comparable to a population without predation and poaching, the number of adults would be twice or three times higher than currently detected in the present study. Yet, considering the recorded number of adult individuals, the annual adult survival rate was considerably low, i.e. 59–69% in males and 67–82% for females. To conclude, the demographic structure of WDE in NKNP showed a high breeding rate, moderate juvenile survival, high survival rate of yearlings, and a low survival rate of adults.Juvenile survival is one of the most fluctuating vital rate parameters, sensitive to population density, stochastic environmental variation, and predation53,54,55. Given the high proportion of juveniles within the population, and the breeding rate higher than that in Cameroon (74%14) and within the captive population (77%34), the juvenile survival rate does not seem to negatively affect the population growth in the NKNP. High breeding rates could be a more robust determinant of population change than AD mortality53, and it is therefore possible that the WDE population size is stable in the NKNP, or even increasing, despite the low adult survival rates. On the other hand, the relatively low numbers of AD individuals in the population indicates low survival rates, which may lead to the decline and final crash of the population54. It is acknowledged that data from two consecutive years was used in the present study, which were not comparable due to different CT settings, and that long-term monitoring, which accounts for variability in vital rates, would be a conservation essential to identify the trend and population change.Based on the present findings of WDE spatiotemporal behaviour and estimates of vital rates, several explanations about multiple processes interacting in the environmental, anthropogenic and conservation context of the park, which inherently affect the small population of WDE, can be inferred. One explanation may suggest that a low proportion of AD WDE and higher JUV survival rates may reflect the influence of growing populations of apex predators in the NKNP, specifically the population of lions56, which may preferentially target the adult individuals57. The age-sex structure also encourages the interpretation that the adult animals are exposed to human-related factors, which prevents them from expanding from the core area of their distribution, exacerbating male-male competition in the limited space34. The poaching activity was also highlighted as an existing threat to WDE populations35. However, law enforcement has been substantially intensified in the core and south-eastern part of the NKNP since 201758, and lion-conservation actions are specifically supported. Thus, the predator populations may have started to grow, which is confirmed by the relative high trapping rate of lions in this core area8. Hence, increased predation may interfere with other environmental factors and consequently affect the WDE population dynamics at the level of AD individuals55,59.A complementary scenario may highlight other factors, specifically, those which maintain the WDE population within a certain spatial extent of the park, i.e. Mont Assirik and Mansa Fara marsh zone. This area can be delimited either ecologically by specific unidentified resources, or by anthropogenic factors, namely a highly frequented trade road crossing the park, wild bushfires, and intensive livestock encroachment in a large band from the borders of the park, inwards (up to 10 km). There is also a vast area in the central part of the park that offers an important space with a supposed carrying capacity for large herbivore populations. This area is, however, outside of the zone of intensified law enforcement, and suffers from inadequate surveillance in the long-term, due to the absence of tracks and therefore being difficult for rangers to access. This area certainly represents an attractive zone for targeted illegal hunting actions. These limiting factors constrain large mammals to concentrate within the zone of Mount Assirik and Mansa Fara marsh, which, in turn, makes animal populations vulnerable to any potential environmental or man-induced incidents, like bush fire.Population sizeThe estimated population size of 195 individuals corresponds with the range of most recent estimates of the WDE population size in the NKNP, i.e. 100–200 (approximately 170) individuals6,7,60. Given the fact that the model contains only the data for AD animals (as no other age category had recapture records), it may be considered that this estimate refers to the number of adult individuals in the population. With regards to Table 3, showing that adults are likely to form 43 to 44% of the whole population, it may be inferred that the actual number of WDE in the NKNP could be higher, even up to 300 individuals, if the data are corrected for the 22% of unidentified individuals. The WDE density estimate of 0.138 individuals/km2 was comparable to densities of Eastern Derby eland in CAR (densities ranging between 0.04 and 0.16 individuals/km2), in Chinko13, and ranging between 0.002 and 0.1 individuals/km2 in the northern CAR61, as well as in Cameroon, with densities ranging between 0.002 and 0.08 individuals/km262. On the other hand, in comparison to other antelope species, the estimated WDE density falls within the range of densities of large herbivores reported from many other sites in African protected areas63, where lower values correspond to the larger areas and are also associated with large browsers, i.e. to the type of diet. Maximum densities of a healthy undisturbed DE population were estimated at about 0.5 individuals/km249, and can reach up to 1.19 individuals/km2 in intensively surveyed hunting zones in Northern CAR61. Thus, the density of WDE in the NKNP could be potentially higher. More

From the phytosociological relevés performed in each sampling area it is evident that hives were positioned in grasslands rich in Alpine herbaceous species (Table S1). In fact, among the 169 identified species, 85% were herbaceous species common in meadows (of Arrhenatherion elatioris and Triseto flavescentis-Polygonion bistortae phytosociological alliance) and acidophilus pastures (Siversio-Nardetum). 15% of the species were trees and shrubs (not abundant in the floristic relevés of the apiary areas considered), including some of beekeeping interest such as: Rhododendron ferrugineum, Castanea sativa and Rubus idaeus. From the MDS biplot (Fig. 2) elevation is the main ecological variable that differentiates sampling areas. In particular, the relevés of stations B and F are characterized by a floristic composition which is different from the areas at higher elevation (characterized by a higher presence of microthermal alpine species). This is due to the separation between the sub-montane belt and the high mountain belt vegetation on the 1.300 m a.s.l. line in the study area25.Figure 2MDS of the phytosociological relevés. Capital letters indicate the six sampling areas, the 1.300 m a.s.l. contour line that separates sub-montane belt and high mountain belt vegetation is highlighted in red.Full size imageAlthough the beehives were positioned in mountain grasslands, melissopalynological analysis presented a different picture. The pollen of numerous species detected through the floristic relevés were found in the honey samples via melissopalynological analysis, although the latter did not totally overlap with the floristic characterization of the area, in particular from a “quantitative” point of view. In fact, the floristic relevés showed a relative richness of herbaceous species (Table S1) peculiar of mountain grasslands that would seem promising for the production of wildflower honeys. Conversely, in the melissopalynological analysis the species considered interesting but not predominant in the botanical description were relevant (Fig. 3 and Table S2).Figure 3MDS of the melissopalynological analysis of the six samples (dots) of mountain wildflower honeys produced in the stations considered. The crosses are the pollens found in the honey samples, the most important are indicated.Full size imageThe premises to produce wildflower honey is that the botanical species contributing must be different and sometimes very numerous, without any of them assuming a dominant character. However, this was not fully evident in our research: although it was possible to identify more than seventy species through melissopalynological analysis and even more through the floristic characterization of the areas, most of them were defined as minor or sporadic pollen (Table S2). Even though apiaries were in mountain grasslands, the most relevant role was played by some woody species/shrubs: Rubus (presumably Rubus idaeous L., identified in the floristic relevés) and rhododendron (Rhododendron ferrugineum L.) for the mountain/subalpine belt and Castanea and Ericaceae (heather) in the submountain belt. Following the rules to define ‘‘unifloral honey’’, three of the wildflower honeys could be defined unifloral or bifloral:

Rhododendron unifloral: honey A (Rhododendron 47.18%), and honey C (Rhododendron 62.93%);

Raspberry unifloral: honey B (Rubus 67.12%)

Raspberry and Rhododendron bifloral: honey D (Rhododendron 34.27% and Rubus 34.74%) as well as honey E (Rubus 44.25%, Rhododendron 34.14%).

Honey F, due to the contribution of pollen from Tilia genus (that was detected only in this sample as an important sporadic pollen, 3.5%) Castanea (96.4% in honey F, but it should be noted that chestnut pollen is an overrepresented pollen) and in the second count Ericaceae (32.45%, that was considered a secondary pollen together with Rubus, with a percentage of 38.59% in honey F) differed from the other honeys (Fig. 3).Rubus pollen was anyway present in good amounts in all the samples considered, and was a dominant pollen in honey B, a secondary pollen in honeys C, D, E and F and a minor pollen in honey A. Sorbus and Tilia pollens were detected only in honey F, while no rhododendron was detected in honey F. Honey D was characterized by a percentage higher than the “rare pollen” category of some important alpine essences, such as Liliaceae, Centaurea, Campanulaceae, Anthyllis f., Polygonum bistorta, Lotus alpinus and Potentilla/fragaria (Table S2).Although wildflower honeys are intrinsically characterized by a high variability compared with unifloral honey, this shows the importance of the formal characterization of honey to obtain a product which satisfies consumer expectations, and it was demonstrated that the botanical origin of honey cannot be based on the claims of local beekeepers by considering the predominant flowers surrounding the hive.Although honeybees are considered supergeneralists in their foraging choices, there are certain key species or plant groups that are particularly important in honeybee foraging2, and many were identified in the botanical characterization of the area, including Rubus idaeus L., Calluna vulgaris L., rhododendron and some present in the broad-leaved woods mentioned such as chestnut (Castanea sativa Mill.) or plants of Tilia genus. In the research work by Hawkins et al.2, Rubus fruticosus L. was among the frequently found species and tree pollen belonging to Castanea sativa L. as well as, for example, species of Malus, Salix and Quercus spp, was frequently seen. These kinds of preferences could relate to the ease of availability and abundance of the plant, the quality and abundance of the nectar and pollen and/or specific nutrients or trace elements provided by these species or neurological aspects (as will be discussed further). As referred by beekeepers, over the last decades the production of mountain wildflower honey, that often does not meet the characteristics expected and presents flavours that are reminiscent of other kinds of honey such as rhododendron or linden or chestnut, is becoming more and more critical and this was absolutely confirmed by this study.This could be linked to the fragmentation of an important habitat of the Alps—mountain grasslands (meaning pastures and meadows) for anthropic and climatic reasons8,9. Honeybees from the same colony forage across areas spanning up to several hundred square kilometres, and at linear distances as far as 9 km from the hive41. Onlooker bees are those in charge of finding nectar sources and of giving instructions to the employed bees, the other foraging bees, that communicate the necessity to look for new resources of food to the onlookers through continuous dance communication42. Among the onlookers, there is a difference between the bees that scout for different nectar sources or recruit to well known floral resources43 and there is an optimal ratio of scouts to recruits, for the most effective collective foraging41. However, this balance may change based on the structure of the landscape in which the bees forage for food44,45,46. Theoretical models47,48 and empirical tests49 suggest that when resources are concentrated into a small number of highly rewarding patches, colonies perform best with few scouts and many recruits, while when resource patches are small, evenly distributed, and easy to locate, successful colonies invest more in scouting than in recruitment. This is strictly linked to climate and social changes in the mountains: mountain grasslands are no longer evenly distributed and easily localizable, as they are scattered among expanding areas of shrublands and forests9 and, for the above-mentioned reasons, it is more efficient for the colony to invest in more recruiters than scouters, as recruiters will identify a small number of highly rewarding patches, such as raspberry or rhododendron shrublands or linden and chestnut woods, that are highly rewarding and very different in quality.This overlaps with individual and collective honeybee behaviour driven by proximate physiological mechanisms that involve the tryptophan metabolism via kynurenine pathway that is one of main neuroprotective mechanisms. In this research, many of the differences/similarities among the samples might be attributed to metabolic alterations within this pathway, represented by relative amounts of kynurenic acid. However, different quinoline structures have also been identified (Fig. 4). Neurotransmitters play a central role in several of the biological processes that honeybees require to perform activities such as foraging behaviour50. A considerable amount of literature highlights the involvement of the neuroprotective kynurenine pathway (KP) final product kynurenic acid (KinA) in the regulation of the stress-related hormone dopamine in the honeybee as well as in other animal species51,52. The major known source of dietary KynA are pollen and nectar produced by sweet chestnuts53,54 and it has been verified that this compound is found in high concentrations in chestnut flowers55. This is coherent with the results of this study: chestnut pollen was found in honey F, produced in the lower station where chestnuts also appear in the floristic relevés, and KynA was found to be a dominant compound in honey F. Interestingly, chestnut pollen was found as sporadic pollen in all the other samples, even those produced in the highest apiary stations (Table S2).Figure 4Kinurenic acid and 3-hydroxyquinaldic acid structure and content in the six honey samples, performed in triplicate. The box diagram representing the median with distribution interval between 25 and 75%.Full size imageFurther, KinA may possess positive properties in a number of pathologies of the gastrointestinal tract, especially colitis, colon obstruction or ulceration56,57. It has been proposed that KinA may also possess antioxidative properties56,57,58,59. This was confirmed by this study, since the wildflower honey with a high component of chestnut pollen was the one with the highest antioxidant properties at the FRSA test (66.61 ± 4.77%), even if lower than manuka honey (84.21 ± 1.04%), a dark honey that is a well-known nutraceutical product and has recently attracted attention for its biological properties, especially for its antioxidant and anti-microbial capacities60. Honey A showed the lowest power (22.40 ± 0.28%) while the other honeys ranked around 40% (Fig. 5). Interestingly, metabolomic analysis revealed the presence of 3-hydroxyquinaldic acid (Fig. 4), which is a kynurenic acid isomer and, although its function has not been elucidated in detail, a few literature data indicate its role as a precursor of naturally occurring peptide antibiotics from the quinomycin family61.Figure 5Results of the FRSA test. Capital letters represent the six honey samples considered. Manuka honey was used as a control.Full size imageIn order to evaluate the ability of honey to induce wound closure, a scratch wound assay was performed (Fig. 6)62. Scratch assay creates a gap in confluent keratinocyte monolayer to mimic a wound. It has already been demonstrated that honeys are able to induce wound closure63 to different extents depending on honey origins and properties.Figure 6The scratch wound test in keratinocytes, HaCaT cells, exposed to honeys. (a) The digitalized pictures of scratched cells after 24 h exposure to 0.5% (w/v) of honeys. (b) The closing percentage wound values after 24 h exposure. Statistics on bars indicate differences compared to the control (CTRL) condition determined by a One-Way ANOVA followed by Dunnett’s test (****p  More

Soil and litter samplingMineral soil (0–15 cm) was collected at the Elizabeth Woods site, a 120-year-old deciduous forest in West Virginia, US (39° 32′ 50.6″ N, − 80° 00′ 00.4″ W). Soils were collected from four 20 × 20 m plots dominated by either AM-associated trees (i.e. Liriodendron tulipifera and Acer saccharum), or ECM-associated trees (i.e. Quercus rubra, Quercus velutina and Carya ovata). These sites have been characterized previously as Culleoka-Westmoreland silt loam soils at the AM sites and Dormont and Guernsey silt loams at the ECM sites40. Soils were also characterized by C:N ratios 11.7 and 14.1 for the AM and ECM soils respectively, with a pH of 6.8 for both soils. Soils with the same mycorrhizal status were pooled and homogenized, air-dried at room temperature for ~ 24 h and sieved through 2.0 mm mesh before the initiation of the experiment. Uniformly 13C labeled litter ( > 97 atom % 13C) from Quercus robur (i.e., ECM substrate) and Liriodendron tulipifera (i.e. AM substrate) leaves (Isolife BV, Wageningen, NL) were incubated in soil mesocosms in a factorial design with five replicates for each treatment combination (2 soil types × 2 substrate types), along with five replicate controls (no 13C substrate addition) for each soil type. The 13C enriched substrates were dried and ground to a powder and added in a suspension of 0.5 ml sterile water to 20 g of soil at a concentration of 400 ug 13C g−1 soil. The control soils received 0.5 ml sterile water additions. These incubations were well mixed and kept at 60% water-holding capacity for the 21-day period at room-temperature18. Chemical characteristics of soils and plant substrates are provided in Table S1.DNA processing and qSIPFor quantitative stable isotope probing, DNA was extracted, quantified, ultracentrifuged, fractionated and sequenced as described in18,26. DNA was extracted using a MoBio PowerSoil HTP Kit following the manufacturer’s instructions. For stable isotope probing, 5 ug of DNA was loaded into a 5-ml ultracentrifuge tube with ~ 3.5 ml of a saturated cesium chloride (CsCl) solution and ~ 900 ml gradient buffer (200 mM Tris, 200 mM KCl, 2 mM EDTA). DNA was separated via ultracentrifugation at 127,000g for 72 h using a TLN-100 rotor in an Optima Max bench top ultracentrifuge (Beckman Coulter, Fullerton, CA, USA). Tubes were fractionated into ~ 25 fractions of 150 µl each, and the density of each fraction was measured with a Raichart AR200 digital refractometer. DNA was purified using an isopropanol precipitation method. The 16S rRNA gene was subsequently quantified and sequenced in samples containing DNA, within the density range 1.660–1.735 gml−1 (~ 10 fractions per sample). To quantify the 16S rRNA gene, quantitative PCR was performed in triplicate using a QuantStudio 5 applied biosystems (Thermo Fisher Scientific) and primers 515F (5′-GTGCCAGCMGCCGCGGTAA-3′) and 806R (5′-GGACTACVSGGGTATCTAAT-3′)41. The PCR program used was as follows: 95 °C for 2 min followed by 45 cycles of 95 °C for 30 s, 64.5 °C for 30 s and 72 °C for 1 min. Libraries were sequenced on an Illumina MiSeq instrument (Illumina, Inc., San Diego, CA, USA) using a 300-cycle v2 reagent kit. Fungal 18S rRNA gene copies in each fraction were also quantified using primers 1380F (5′-CCCTGCCHTTTGTACACAC-3′) and 1510R (5′-CCTTCYGCAGGTTCACCTAC-3′). The PCR program used was as follows: 98 °C for 3 min followed by 40 cycles of 98 °C for 45 s, 60 °C for 45 s and 72 °C for 30 s. DNA fractions were amplified for fungal ITS rRNA genes using primers ITS4F (5′-AGCCTCCGCTTATTGATATGCTTAART-3′) and 5.8SF (5′-AACTTTYRRCAAYGGATCWCT-3′)42 and 300-bp paired-end read chemistry on an IlluminaMiSeq (Illumina, Inc., San Diego, CA, USA). The PCR program used was as follows: 95 °C for 6 min followed by 35 cycles of 95 °C for 15 s, 55 °C for 30 s, and 72 °C for 1 min. DNA fractions were then sequenced using a 500 cycle v2 reagent kit.Files came pre-split and joined multiple paired ends that we combined to pick operational taxonomic units (OTU). Open reference OTUs were picked at 97% identity using SILVA 128 release database for Bacteria and RDP database for Fungi. Taxa were analyzed at the ‘OTU’ level from the QIIME L7 table. Calculation of 13C excess atom fraction (EAF) was performed for each taxon as described previously18,19. Briefly, using the CsCl density gradient data, a weighted average density (WAD) was computed for each taxon’s DNA extracted from control soils that did not receive an isotopically enriched substrate. This natural abundance WAD was then compared to the taxon’s WAD following incubation with the 13C enriched material. The change in WAD can be used to quantify the amount of isotope incorporated into the DNA17,18. Preliminary data analysis revealed an effect of ultracentrifuge tube on estimation of phylotype weighted average density, probably a consequence of slight differences in CsCl density gradients between tubes. This technical error was corrected as previously described18,19. In addition to the samples subjected to qSIP analysis we also extracted and analyzed fungal and bacterial OTU’s from control soils where the DNA was extracted prior to incubation.FTICR-MS and lipidomic analysesSoil from substrate-incubated and controls mesocosms were processed and analyzed with Fourier transform ion cyclotron resonance mass spectrometry (FTICR-MS), using a 12 T Bruker SolariX FTICR mass spectrometer at the Environmental Molecular Sciences Laboratory in Richland, WA, as described in Fudyma et al.43. Briefly, 100 mg of dried soil or litter substrate was extracted using an adjusted Folch extraction44. Extraction was performed on each sample by sequentially adding 2 ml MeOH, followed by a 5 s vortex; 4 ml CHCl3, followed by a 5 s vortex; sonication at 25 °C for 1 h (CPX3800 Ultrasonic Bath, Fisherbrand); addition of 1.25 ml of H2O, followed by a slight mix to achieve bi‐layer separation; and incubated at 4 °C overnight. The top, aqueous layer (metabolite—polar) was pipetted off into 1 ml glass vials and stored at − 80 °C until FTICR‐MS. The bottom, chloroform layer was dried down and stored in 50:50 methanol:chloroform until lipidomics analysis.A standard Bruker electrospray ionization (ESI) source was used to generate negatively charged molecular ions in the metabolite fraction. Samples were then introduced directly to the ESI source. The instrument settings were optimized by tuning on a Suwannee River fulvic acid (SRFA) standard, purchased from International Humic Substances Society (IHCC). Blanks (HPLC grade methanol) were analyzed at the beginning and end of the day to monitor potential carry over from one sample to another. The instrument was flushed between samples using a mixture of water and methanol. One hundred and forty‐four individual scans were averaged for each sample and internally calibrated using an organic matter homologous series separated by 14 Da (CH2 groups). The mass measurement accuracy was less than 1 ppm for singly charged ions across a broad m/z range (m/z 300– 800). Data analysis software (Bruker Daltonik version 4.2) was used to convert raw spectra to a list of m/z values, applying the FTMS peak picker module with a signal-to noise ratio (S/N) threshold set to 7 and absolute intensity threshold set to the default value of 100. Chemical formulae were then assigned using in-house software following the compound identification algorithm that was described in Tolić et al.45. Peaks below 200 and above 800 were dropped to select only for calibrated and assigned peaks. Chemical formulae were assigned based on the following criteria: S/N  > 7 and mass measurement error  800 were not detected in our samples. The m/z values represent the molecular mass (in Dalton) of the detected ions since all detected ions were singly charged ions. While our results do not represent a quantitative characterization of OM, the values presented are relative differences and should be representative of the samples. Finally, we would like to acknowledge that we were not able to see any clear evidence of 13C label in our FTICR-MS analysis of the soil samples. The lack of 13C label in our FTICR-MS analysis of the soil samples even though they received labeled substrate could be either due to the fact that most of the labeled substrates produced by microbial activities were of low molecular weight, which cannot be detected by FTICR-MS and/or the leftover labeled substrate was of low abundance compared to the organic compounds previously present in the soil matrix. As such, we used the FTCIR-MS data to identify shifts in the overall composition of the chemical compounds in each soil.Lipids in the chloroform fraction were analyzed by LC‐MS/MS in both positive and negative ESI modes using a linear trap quadropole (LTQ) Orbitrap Velos mass spectrometer (Thermo Fisher Scientific), as described in detail previously46. Lipid species were identified using the LIQUID tool46 followed by manual data inspection. Confidently identified lipid species were quantified using MZmine47 and the peak intensities were normalized by linear regression and central tendency (i.e., identifying a central or typical value for a probability distribution) using InfernoRDN.Statistical analysisAll data analyses were performed using R 3.2.048. To examine the effects of soil type, substrate type and their interaction in the bacterial, fungal and chemical composition of DOM and the lipid pool; Bray–Curtis distance matrices were compared with permutational multivariate analysis of variance (PerMANOVA) and visualized with Principle Coordinate Analysis (PCoA) using vegan package49. PerMANOVA analysis were run on the relative abundance and on the 13C EAF of individual microbial taxa, separately for both bacterial and fungal communities.The analyses for FTICR-MS were performed separately for control and incubated soils using all assigned molecular formulae remaining after quality filtering31. In all cases, we applied a Z-score standardization before calculating Bray–Curtis distance matrices49. We analyzed the results from FTICR-MS as resulting from the decomposition of the added substrates for two reasons. First, this is a fully factorial design where individual soil samples were split to either receive AM poplar or ECM oak litter substrate. Thus, each soil sample starts with the same characteristics and the changes at the end of the incubation period should reflect the processing of litter. Second, we excluded molecular formulae present in the litters and thus, the differences we report in each soil type are derived from this processing (or the lack of it).We calculated aggregated indices that characterize both the composition and the physicochemical properties of the microbial (both bacteria and fungi) and the SOM and lipid pool34,36. For bacterial and fungal communities, we quantified Shannon–Weaver diversity index for each sample H′ = (-{sum }_{i=1}^{S} pi ln(pi)) (where pi is the proportion of species I) using the relative abundance of individual microbial taxa50. To find the percent of substrate assimilation by individual taxa, we calculated the proportion of C assimilated by each group as previously described18,51 as a percent. For SOM and lipid molecular formulae, we separately calculated weighted means of formula-based characteristics (i.e. m/z, Aromaticity Index—AI; H/C, O/C, and Nominal Oxidation State of Carbon-NOSC) as the sum of the product of the single-formula information (i.e. m/zi, AIi, H/Ci and NOSCi) and the relative intensity (Ii) divided by the sum of all intensities (e.g., m/z sample1 = ({sum }_{i=1}^{S})(m/zi ·Ii)/Σ(Ii)). With these metrics we obtained sample-level information related to the molecular size (i.e. m/z), the molecular bioavailability (i.e. higher H/C ratio), the molecular reactiveness (i.e. lower AI) and the energetic rewards from molecular oxidative degradation (i.e. higher NOSC) of the SOM, which allows to infer the potential of decomposition products to form stable SOM12,31,35. Detailed information of the calculated indices can be found in the literature31,35,36.We further tested the effects of soil type, substrate type and their interaction on each index using the “lm” function from the “stats” package. In these analyses, P values were approximated by an F test using Type II ANOVA tests with Kenward-Roger Degrees of Freedom52. When interactions between soil and substrate type were found at P  More