Philippa Kaur

1.
Barnosky, A. D. et al. Has the Earth’s sixth mass extinction already arrived?. Nature 471, 51–57 (2011).
ADS  CAS  PubMed  Article  Google Scholar 
2.
McCauley, D. J. et al. Marine defaunation: Animal loss in the global ocean. Science 347, 1255641–1255641 (2015).
PubMed  Article  CAS  Google Scholar 

3.
Chapin, F. S. III. et al. Consequences of changing biodiversity. Nature 405, 234–242 (2000).
CAS  PubMed  Article  Google Scholar 

4.
Díaz, S., Fargione, J., Chapin, F. S. & Tilman, D. Biodiversity loss threatens human well-being. PLoS Biol. 4, e277 (2006).
PubMed  PubMed Central  Article  CAS  Google Scholar 

5.
Worm, B. et al. Rebuilding global fisheries. Science 325, 578–585 (2009).
ADS  CAS  PubMed  Article  Google Scholar 

6.
Defeo, O. & Castilla, J. C. More than one bag for the world fishery crisis and keys for co-management successes in selected artisanal Latin American shellfisheries. Rev. Fish Biol. Fish. 15, 265–283 (2005).
Article  Google Scholar 

7.
Pauly, D. & Zeller, D. Catch reconstructions reveal that global marine fisheries catches are higher than reported and declining. Nat. Commun. 7, 10244 (2016).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

8.
Defeo, O. et al. Co-management in Latin American small-scale shellfisheries: Assessment from long-term case studies. Fish Fish. 17, 176–192 (2016).
Article  Google Scholar 

9.
Gelcich, S. et al. Fishers’ perceptions on the Chilean coastal TURF system after two decades: Problems, benefits, and emerging needs. Bull. Mar. Sci. 93, 53–67 (2017).
Article  Google Scholar 

10.
Castilla, J. C., Gelcich, S. & Defeo, O. Successes, lessons, and projections from experience in marine benthic invertebrate artisanal fisheries in Chile. In Fisheries Management (eds McClanahan, T. R. & Castilla, J. C.) 23–42 (Blackwell Publishing Ltd, Hoboken, 2007). https://doi.org/10.1002/9780470996072.ch2.
Google Scholar 

11.
Gelcich, S. et al. Navigating transformations in governance of Chilean marine coastal resources. Proc. Natl. Acad. Sci. 107, 16794–16799 (2010).
ADS  CAS  PubMed  Article  Google Scholar 

12.
Kéfi, S. et al. Network structure beyond food webs: Mapping non-trophic and trophic interactions on Chilean rocky shores. Ecology 96, 291–303 (2015).
Article  Google Scholar 

13.
Pérez-Matus, A. et al. Temperate rocky subtidal reef community reveals human impacts across the entire food web. Mar. Ecol. Prog. Ser. 567, 1–16 (2017).
ADS  Article  Google Scholar 

14.
Pérez-Matus, A., Carrasco, S. A., Gelcich, S., Fernandez, M. & Wieters, E. A. Exploring the effects of fishing pressure and upwelling intensity over subtidal kelp forest communities in Central Chile. Ecosphere 8, e01808 (2017).
Article  Google Scholar 

15.
Gelcich, S. et al. Territorial user rights for fisheries as ancillary instruments for marine coastal conservation in Chile: Gelcich et al. Conserv. Biol. 26, 1005–1015 (2012).
PubMed  Article  Google Scholar 

16.
Oyanedel, R., Keim, A., Castilla, J. C. & Gelcich, S. Illegal fishing and territorial user rights in Chile: Illegal fishing. Conserv. Biol. 32, 619–627 (2018).
PubMed  Article  Google Scholar 

17.
Donlan, C. J., Wilcox, C., Luque, G. M. & Gelcich, S. Estimating illegal fishing from enforcement officers. Sci. Rep. 10, 12478 (2020).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

18.
Andreu-Cazenave, M., Subida, M. D. & Fernandez, M. Exploitation rates of two benthic resources across management regimes in central Chile: Evidence of illegal fishing in artisanal fisheries operating in open access areas. PLoS ONE 12, e0180012 (2017).
PubMed  PubMed Central  Article  CAS  Google Scholar 

19.
Castilla, J. C. Coastal marine communities: Trends and perspectives from human-exclusion experiments. Trends Ecol. Evol. 14, 280–283 (1999).
CAS  PubMed  Article  Google Scholar 

20.
Somero, G. N. The physiology of climate change: How potentials for acclimatization and genetic adaptation will determine ‘winners’ and ‘losers’. J. Exp. Biol. 213, 912–920 (2010).
CAS  PubMed  Article  Google Scholar 

21.
Hoegh-Guldberg, O. & Bruno, J. F. The impact of climate change on the world’s marine ecosystems. Science 328, 1523–1528 (2010).
ADS  CAS  PubMed  Article  Google Scholar 

22.
Brose, U. et al. Climate change in size-structured ecosystems. Philos. Trans. R. Soc. B Biol. Sci. 367, 2903–2912 (2012).
Article  Google Scholar 

23.
Doney, S. C. et al. Climate change impacts on marine ecosystems. Annu. Rev. Mar. Sci. 4, 11–37 (2012).
ADS  Article  Google Scholar 

24.
Kwiatkowski, L., Aumont, O. & Bopp, L. Consistent trophic amplification of marine biomass declines under climate change. Glob. Change Biol. 25, 218–229 (2019).
ADS  Article  Google Scholar 

25.
Bakun, A. Coastal ocean upwelling. Science 247, 198–201 (1990).
ADS  CAS  PubMed  Article  Google Scholar 

26.
Bakun, A., Field, D. B., Redondo-Rodriguez, A. & Weeks, S. J. Greenhouse gas, upwelling-favorable winds, and the future of coastal ocean upwelling ecosystems. Glob. Change Biol. 16, 1213–1228 (2010).
ADS  Article  Google Scholar 

27.
Thiel, M. et al. The Humboldt current system of northern and central Chile: Oceanographic processes, ecological interactions and socioeconomic feedback. In Oceanography and Marine Biology Vol. 20074975 (eds Gibson, R. et al.) 195–344 (CRC Press, Boca Raton, 2007).
Google Scholar 

28.
Morales, C., Hormazabal, S., Andrade, I. & Correa-Ramirez, M. Time-space variability of chlorophyll-a and associated physical variables within the region off central-southern Chile. Remote Sens. 5, 5550–5571 (2013).
ADS  Article  Google Scholar 

29.
Aiken, C. M., Navarrete, S. A. & Pelegrí, J. L. Potential changes in larval dispersal and alongshore connectivity on the central Chilean coast due to an altered wind climate. J. Geophys. Res. 116, G04026 (2011).
ADS  Article  Google Scholar 

30.
Blanchard, J. L. et al. Potential consequences of climate change for primary production and fish production in large marine ecosystems. Philos. Trans. R. Soc. B Biol. Sci. 367, 2979–2989 (2012).
Article  Google Scholar 

31.
Testa, G., Masotti, I. & Farías, L. Temporal variability in net primary production in an upwelling area off central Chile (36°S). Front. Mar. Sci. 5, 179 (2018).
Article  Google Scholar 

32.
Batten, S. D. et al. A global plankton diversity monitoring program. Front. Mar. Sci. 6, 321 (2019).
Article  Google Scholar 

33.
Chust, G. et al. Biomass changes and trophic amplification of plankton in a warmer ocean. Glob. Change Biol. 20, 2124–2139 (2014).
ADS  Article  Google Scholar 

34.
Weidberg, N. et al. Spatial shifts in productivity of the coastal ocean over the past two decades induced by migration of the Pacific Anticyclone and Bakun’s effect in the Humboldt Upwelling Ecosystem. Glob. Planet. Change 193, 103259 (2020).
Article  Google Scholar 

35.
Aguirre, C., García-Loyola, S., Testa, G., Silva, D. & Farias, L. Insight into anthropogenic forcing on coastal upwelling off south-central Chile. Elem. Sci. Anth. 6, 59 (2018).
Article  Google Scholar 

36.
Valdovinos, F. S. Mutualistic networks: Moving closer to a predictive theory. Ecol. Lett. 22, 1517–1534 (2019).
PubMed  Article  PubMed Central  Google Scholar 

37.
Pascual, M. & Dunne, J. A. Ecological Networks: Linking Structure to Dynamics in Food Webs (Santa Fe Institute Studies on the Sciences of Complexity) (Oxford University Press, Oxford, 2006).
Google Scholar 

38.
Dunne, J. A., Williams, R. J. & Martinez, N. D. Network structure and biodiversity loss in food webs: Robustness increases with connectance. Ecol. Lett. 5, 558–567 (2002).
Article  Google Scholar 

39.
Curtsdotter, A. et al. Robustness to secondary extinctions: Comparing trait-based sequential deletions in static and dynamic food webs. Basic Appl. Ecol. 12, 571–580 (2011).
Article  Google Scholar 

40.
Ramos-Jiliberto, R., Valdovinos, F. S., Moisset de Espanés, P. & Flores, J. D. Topological plasticity increases robustness of mutualistic networks: Interaction rewiring in mutualistic networks. J. Anim. Ecol. 81, 896–904 (2012).
Article  Google Scholar 

41.
Valdovinos, F. S., Moisset de Espanés, P., Flores, J. D. & Ramos-Jiliberto, R. Adaptive foraging allows the maintenance of biodiversity of pollination networks. Oikos 122, 907–917 (2013).
Article  Google Scholar 

42.
Allesina, S. & Pascual, M. Googling food webs: Can an eigenvector measure species’ importance for coextinctions?. PLoS Comput. Biol. 5, e1000494 (2009).
ADS  MathSciNet  PubMed  PubMed Central  Article  CAS  Google Scholar 

43.
de Santana, C., Rozenfeld, A., Marquet, P. & Duarte, C. Topological properties of polar food webs. Mar. Ecol. Prog. Ser. 474, 15–26 (2013).
ADS  Article  Google Scholar 

44.
Eklöf, A., Tang, S. & Allesina, S. Secondary extinctions in food webs: A Bayesian network approach. Methods Ecol. Evol. 4, 760–770 (2013).
Article  Google Scholar 

45.
Staniczenko, P. P. A., Lewis, O. T., Jones, N. S. & Reed-Tsochas, F. Structural dynamics and robustness of food webs: Structural dynamics and robustness of food webs. Ecol. Lett. 13, 891–899 (2010).
PubMed  Article  Google Scholar 

46.
Albert, R., Jeong, H. & Barabási, A. Error and attack tolerance of complex networks. Nature 406, 378–382 (2000).
ADS  CAS  PubMed  Article  Google Scholar 

47.
Ives, A. R. & Cardinale, B. J. Food–web interactions govern the resistance of communities after non-random extinctions. Nature 429, 174–177 (2004).
ADS  CAS  PubMed  Article  Google Scholar 

48.
Rebolledo, R., Navarrete, S. A., Kéfi, S., Rojas, S. & Marquet, P. A. An open-system approach to complex biological networks. SIAM J. Appl. Math. 79, 619–640 (2019).
MathSciNet  Article  Google Scholar 

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

50.
Glaum, P., Cocco, V. & Valdovinos, F. S. Integrating economic dynamics into ecological networks: The case of fishery sustainability. Sci. Adv. 6, eaaz4891 (2020).
PubMed  PubMed Central  Article  Google Scholar 

51.
Williams, R. J. Network 3D: Visualizing and modelling food webs and other complex networks. Microsoft Res. Camb. UK. http://research.microsoft.com/en-us/um/cambridge/groups/science/tools/network3d/network3d.htm (2010).

52.
Richard, J. W., 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, Berlin, 2006). https://doi.org/10.1007/978-1-4020-5337-5_2.
Google Scholar 

53.
Boit, A., Martinez, N. D., Williams, R. J. & Gaedke, U. Mechanistic theory and modelling of complex food-web dynamics in Lake Constance: Mechanistic modelling of complex food web dynamics. Ecol. Lett. 15, 594–602 (2012).
PubMed  Article  Google Scholar 

54.
Kuparinen, A., Boit, A., Valdovinos, F. S., Lassaux, H. & Martinez, N. D. Fishing-induced life-history changes degrade and destabilize harvested ecosystems. Sci. Rep. 6, 22245 (2016).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

55.
Jackson, J. B. C. Historical overfishing and the recent collapse of coastal ecosystems. Science 293, 629–637 (2001).
CAS  PubMed  Article  Google Scholar 

56.
Pauly, D. Fishing down marine food webs. Science 279, 860–863 (1998).
ADS  CAS  PubMed  Article  Google Scholar 

57.
Jordán, F., Okey, T. A., Bauer, B. & Libralato, S. Identifying important species: Linking structure and function in ecological networks. Ecol. Model. 216, 75–80 (2008).
Article  Google Scholar 

58.
Castilla, J. C. & Fernandez, M. Small-scale benthic fisheries in Chile: On co-management and sustainable use of benthic invertebrates. Ecol. Appl. 8, S124–S132 (1998).
Article  Google Scholar 

59.
Allesina, S., Bodini, A. & Pascual, M. Functional links and robustness in food webs. Philos. Trans. R. Soc. B Biol. Sci. 364, 1701–1709 (2009).
Article  Google Scholar 

60.
de Visser, S. N., Freymann, B. P. & Olff, H. The Serengeti food web: Empirical quantification and analysis of topological changes under increasing human impact: Topological changes under human impact. J. Anim. Ecol. 80, 484–494 (2011).
PubMed  Article  Google Scholar 

61.
Srinivasan, U. T., Dunne, J. A., Harte, J. & Martinez, N. D. Response of complex food webs to realistic extinction sequences. Ecology 88, 671–682 (2007).
PubMed  Article  Google Scholar 

62.
Camus, P. A., Arancibia, P. A. & Ávila-Thieme, M. I. A trophic characterization of intertidal consumers on Chilean rocky shores. Rev. Biol. Mar. Oceanogr. 48, 431–450 (2013).
Article  Google Scholar 

63.
Lopez, D. N., Camus, P. A., Valdivia, N. & Estay, S. A. High temporal variability in the occurrence of consumer–resource interactions in ecological networks. Oikos 126, 1699–1707 (2017).
Article  Google Scholar 

64.
Arim, M. & Marquet, P. A. Intraguild predation: A widespread interaction related to species biology: Intraguild predation. Ecol. Lett. 7, 557–564 (2004).
Article  Google Scholar 

65.
Teagle, H., Hawkins, S. J., Moore, P. J. & Smale, D. A. The role of kelp species as biogenic habitat formers in coastal marine ecosystems. J. Exp. Mar. Biol. Ecol. 492, 81–98 (2017).
Article  Google Scholar 

66.
Vásquez, J. A. The brown seaweeds fishery in Chile. In Fisheries and Aquaculture in the Modern World (ed. Mikkola, H.) (InTech, London, 2016). https://doi.org/10.5772/62876.
Google Scholar 

67.
Belmadani, A., Echevin, V., Codron, F., Takahashi, K. & Junquas, C. What dynamics drive future wind scenarios for coastal upwelling off Peru and Chile?. Clim. Dyn. 43, 1893–1914 (2014).
Article  Google Scholar 

68.
Wang, Y., Luo, Y., Lu, J. & Liu, F. Changes in ENSO amplitude under climate warming and cooling. Clim. Dyn. 52, 1871–1882 (2019).
Article  Google Scholar 

69.
Cai, W. et al. Increased variability of eastern Pacific El Niño under greenhouse warming. Nature 564, 201–206 (2018).
ADS  CAS  PubMed  Article  Google Scholar 

70.
Cai, W. et al. Increased frequency of extreme La Niña events under greenhouse warming. Nat. Clim. Change 5, 132–137 (2015).
ADS  Article  Google Scholar 

71.
Fussmann, K. E., Schwarzmüller, F., Brose, U., Jousset, A. & Rall, B. C. Ecological stability in response to warming. Nat. Clim. Change 4, 206–210 (2014).
ADS  Article  Google Scholar 

72.
Hays, G., Richardson, A. & Robinson, C. Climate change and marine plankton. Trends Ecol. Evol. 20, 337–344 (2005).
PubMed  Article  Google Scholar 

73.
Jochum, M., Schneider, F. D., Crowe, T. P., Brose, U. & O’Gorman, E. J. Climate-induced changes in bottom-up and top-down processes independently alter a marine ecosystem. Philos. Trans. R. Soc. B Biol. Sci. 367, 2962–2970 (2012).
Article  Google Scholar 

74.
Hallegraeff, G. M. A review of harmful algal blooms and their apparent global increase. Phycologia 32, 79–99 (1993).
Article  Google Scholar 

75.
He, Q. & Silliman, B. R. Climate change, human impacts, and coastal ecosystems in the anthropocene. Curr. Biol. 29, R1021–R1035 (2019).
CAS  PubMed  Article  Google Scholar 

76.
Brown, C. J., Saunders, M. I., Possingham, H. P. & Richardson, A. J. Interactions between global and local stressors of ecosystems determine management effectiveness in cumulative impact mapping. Divers. Distrib. 20, 538–546 (2014).
Article  Google Scholar 

77.
Crain, C. M., Kroeker, K. & Halpern, B. S. Interactive and cumulative effects of multiple human stressors in marine systems. Ecol. Lett. 11, 1304–1315 (2008).
PubMed  Article  Google Scholar 

78.
Dunne, J. A. et al. The roles and impacts of human hunter-gatherers in North Pacific marine food webs. Sci. Rep. 6, 21179 (2016).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

79.
Hale, K. R. S., Valdovinos, F. S. & Martinez, N. D. Mutualism increases diversity, stability, and function of multiplex networks that integrate pollinators into food webs. Nat. Commun. 11, 2182 (2020).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

80.
Kéfi, S., Miele, V., Wieters, E. A., Navarrete, S. A. & Berlow, E. L. How structured is the entangled bank? The surprisingly simple organization of multiplex ecological networks leads to increased persistence and resilience. PLoS Biol. 14, e1002527 (2016).
PubMed  PubMed Central  Article  CAS  Google Scholar 

81.
Miele, V., Guill, C., Ramos-Jiliberto, R. & Kéfi, S. Non-trophic interactions strengthen the diversity—Functioning relationship in an ecological bioenergetic network model. PLoS Comput. Biol. 15, e1007269 (2019).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

82.
Morgan, S. G., Fisher, J. L., Miller, S. H., McAfee, S. T. & Largier, J. L. Nearshore larval retention in a region of strong upwelling and recruitment limitation. Ecology 90, 3489–3502 (2009).
PubMed  Article  Google Scholar 

83.
Ospina-Alvarez, A., Weidberg, N., Aiken, C. M. & Navarrete, S. A. Larval transport in the upwelling ecosystem of central Chile: The effects of vertical migration, developmental time and coastal topography on recruitment. Prog. Oceanogr. 168, 82–99 (2018).
ADS  Article  Google Scholar 

84.
Sakai, A. K. et al. The population biology of invasive species. Annu. Rev. Ecol. Syst. 32, 305–332 (2001).
Article  Google Scholar 

85.
Thierry, A. et al. Adaptive foraging and the rewiring of size-structured food webs following extinctions. Basic Appl. Ecol. 12, 562–570 (2011).
Article  Google Scholar 

86.
Valdovinos, F. S., Ramos-Jiliberto, R., Garay-Narváez, L., Urbani, P. & Dunne, J. A. Consequences of adaptive behaviour for the structure and dynamics of food webs: Adaptive behaviour in food webs. Ecol. Lett. 13, 1546–1559 (2010).
Article  Google Scholar 

87.
Williams, R. J. Effects of network and dynamical model structure on species persistence in large model food webs. Theor. Ecol. 1, 141–151 (2008).
Article  Google Scholar 

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

89.
Menge, B. A. & Menge, D. N. L. Dynamics of coastal meta-ecosystems: The intermittent upwelling hypothesis and a test in rocky intertidal regions. Ecol. Monogr. 83, 283–310 (2013).
Article  Google Scholar 

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

91.
Berlow, E. L. et al. Simple prediction of interaction strengths in complex food webs. Proc. Natl. Acad. Sci. 106, 187–191 (2009).
ADS  CAS  PubMed  Article  Google Scholar 

92.
Jonsson, T., Kaartinen, R., Jonsson, M. & Bommarco, R. Predictive power of food web models based on body size decreases with trophic complexity. Ecol. Lett. 21, 702–712 (2018).
PubMed  Article  Google Scholar 

93.
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 

94.
Ávila-Thieme, M. I., Corcoran, D., Valdovinos, F. S., Navarrete, S. A. & Marquet, P. A. NetworkExtinction: Extinction Simulation in Food Webs. (R package version 0.1.3., 2018).

95.
Schneider, F. D., Brose, U., Rall, B. C. & Guill, C. Animal diversity and ecosystem functioning in dynamic food webs. Nat. Commun. 7, 12718 (2016).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar  More

Study site
The present study was conducted in evergreen broad-leaved subtropical forests in the northern part of Okinawa Island, south-western Japan. Samples were collected in a secondary forest within Yona Experimental Forest of University of the Ryukyus (26°9′ N, 128°5′ E, ca 250–330 m a.s.l.). The mean annual temperature was 20.7 °C and the annual precipitation was 2487 mm. The topography is hilly and dissected. The bedrock is composed of sandstone and slate, and yellow soil has developed. The forest stand was dominated by Castanopsis sieboldii and Schima wallichii ssp. liuliuensis with a maximum height of 20 m28.
Litterbag experiment
Decomposition of leaf litter of six tree species (C. sieboldii, S. wallichii, Daphniphyllum teijsmannii, Persea thunbergii, Distylium racemosum, and Camellia japonica) was studied using a litterbag method, according to the procedure detailed previously13. These six tree species are dominant components of the forest canopy in the study site28. In short, a study plot of 50 m × 10 m (500 m2) was laid out in Yona Experimental Forest and was divided into 125 grids of 2 × 2 m. Freshly fallen leaves of six tree species were collected from the soil surface in March 2008. The leaves were dried in an oven at 40 °C for 1 week to a constant mass. Leaf litter (4.00 g) was placed in litterbags (24 × 18 cm) made of nylon with a mesh size of approximately 2 mm and incubated within the 500 m2 study plot for 18 months from April 2008 to October 2009. Nine litterbags per tree species were retrieved at 3, 6, 9, 12, 15, and 18 months after initiation of the experiment and used for measurement of the remaining mass of whole leaf litter13. In the present study, the bleached leaf area and leaf mass per area (LMA) and chemical compositions of bleached and nonbleached portions were then measured as described below. The LMA indicates the remaining mass of leaf tissues and represents the extent of decomposition in the bleached and nonbleached portions. Leaf litters of S. wallichii and D. teijsmannii collected at 12, 15 and 18 months of decomposition were too fragmented to measure bleached leaf area.
Measurement and chemical analyses
Leaves were pressed between layers of plywood and paper and oven-dried at 40 °C for 1 week. The leaves were photocopied, scanned, and measured for the total leaf area and the proportion of bleached area according to the method described previously29. A 6-mm-diameter cork borer was then used to excise leaf disks, avoiding the primary vein, from the bleached area and surrounding nonbleached area of the same leaves collected for the first 9 months of decomposition. The disks were oven-dried again at 40 °C for 1 week and weighed to calculate LMA. The disks were combined to make 1 sample each of bleached and nonbleached leaf area for each tree species collected at each sampling occasion and used for chemical analyses as described below. Leaf disks could not be excised from leaves collected at 12, 15, and 18 months of decomposition because of fragmentation.
Litter materials were ground in a laboratory mill (0.5-mm screen). The amount of acid unhydrolyzable residue (AUR) and total carbohydrates (TCH) was estimated by means of gravimetry as acid-insoluble residue, using hot sulfuric acid digestion30 and by a phenol–sulfuric acid method31. Total N content was measured by automatic gas chromatography (NC analyzer SUMIGRAPH NC-900, Sumitomo Chemical Co., Osaka, Japan). Details of the methods followed Osono13. The contents of AUR and TCH were expressed in g/g dry litter, and that of total N was in mg/g dry litter. The mass of these components per leaf area was calculated by multiplying the contents and LMA. The AUR fraction contains a mixture of organic compounds in various proportions, including condensed tannins, phenolic and carboxylic compounds, alkyl compounds such as cutins, and true lignin16. No data were available for total N content of the nonbleached portions of D. teijsmannii at 9 months of decomposition because of the small amount of sample.
To analyze the chemical composition of bleached leaf tissues more in detail and to compare it with that of nonbleached portions for multiple tree species, samples of bleached leaf litter were collected during fieldworks in March 2007 and in April 2011. These bleached leaf litter were separated into bleached and nonbleached litter samples to be used for measurement and chemical analyses (Table S1). Bleached leaf litter of 20 tree species was used for measurement of LMA, and the samples of 13 of the 20 tree species were further analyzed for the contents of AUR and TCH, as described above (Table S2). Samples were extracted with alcohol-benzene at room temperature (15–20 °C) to remove extractives (EXT; soluble polyphenols, hydrocarbons, and pigments) and to calculate the content of this fraction.
Solid-state Cross polarization (CP) magic angle spinning (MAS) 13C NMR spectra of bleached and nonbleached litter samples for 12 tree species were obtained with an Alpha 300 FT NMR system (JEOL, Tokyo) operating at 75.45 MHz under the following conditions32: pulse repetition time of 3.1 s, CP contact time of 1 ms, sweep width of 35 kHz, acquisition time of 0.117 s, and MAS of 6 kHz. The finely powdered sample was tightly packed into a high-speed spinning NMR tube (rotor: zirconia, cap: KEL-F, 6-mm i.d., JEOL). Chemical shifts are quoted with respect to tetramethylsilane but were determined by referring to an external sample of adamantane (29.50 ppm). The 13C NMR spectra (Fig. S2) were divided into four chemical shift ranges, as follows: 0 to 45 ppm for alkyl-C (including major C of cutins and suberins), 45 to 110 ppm for O-alkyl-C (oxygen-substituted C in alcohols and ethers, including cellulose, hemicellulose, and other polysaccharides), 110 to 160 ppm for aromatic C (including mainly condensed tannins, hydrolyzable tannins, and lignin), and 160 to 190 ppm for carbonyl C (including carboxylic-C and carbonyl-C)33. The relative area of these chemical shift regions was calculated for each spectrum as the percentage of total area by using computer software ALICE 2 for Windows (JEOL) (Table S3).
Nitrogen attached to leaf litter was determined by extraction and colorimetric analyses of the extractants for 13 tree species. Approximately 100 mg of bleached or nonbleached leaf litter was shaken with 10 ml of 2 M KCl in a 15-ml centrifuge tube on a shaker for 1 h. The suspension was centrifuged at 3000 rpm for 10 min and filtered with glass fiber filters (GF/F, Whatman). The total extractable nitrogen (TEN) in the extractants was measured by the alkali persulfate digestion method34. Ammonium-nitrogen (NH4+-N), nitrate-nitrogen (NO3–-N), and nitrite-nitrogen (NO2–-N) were determined colorimetrically35 for the pre-digested samples. Extractable organic nitrogen (EON) was calculated subtracting these three forms of inorganic nitrogen from TEN (Table S4).
Statistical analyses
Linear relationships between LMA and decomposition time and between contents of AUR, TCH, and total N and accumulated mass loss of leaf tissue were examined separately for bleached and nonbleached portions according to the following equations:

$${text{Ln}},left[ {{text{LMA}},left( {% ,{text{original}},{text{value}}} right)} right] , = a + b times , left( {{text{time}},{text{in}},{text{months}}} right)$$
(1)

$${text{AUR}},{text{TCH}},{text{and}},{text{N}},{text{content }} = a + b times left( {{text{accumulated}},{text{mass}},{text{loss}},{text{of}},{text{leaf}},{text{tissue}}} right)$$
(2)

Accumulated mass loss of leaf tissue of bleached and nonbleached portions after a given period was calculated as the loss of LMA relative to the initial LMA values, expressed as a percentage. Intercepts (a) and slopes (b) of regression equations were calculated for the linear relationships using least-square regression15. The slope of the regression Eq. (1) represented the decomposition constant36. The slopes of the regression Eq. (2) describing AUR and N dynamics represented the N concentration increase rate and the lignin concentration increase rate, respectively15. A paired t-test was used to evaluate the difference between bleached and nonbleached portions in the slopes of regression equations for LMA, AUR, TCH, and N in decomposing leaf litter of 6 tree species and in LMA, contents of proximate organic chemical components, relative area of 13C NMR spectra, and contents of dissolved N in leaf litter of multiple tree species. More

Field observations
Since 1999, the FIA program has operated an extensive, nationally consistent forest inventory designed to monitor changes in forests across all lands in the US61. We used FIA data from 10 states in the continental western US (Washington, Oregon, California, Idaho, Montana, Utah, Nevada, Colorado, Arizona, and New Mexico) to quantify shifts in relative live tree density, excluding Wyoming due to a lack of repeated censuses (Fig. 1). This region spans a wide variety of climatic regimes and forest types, ranging from temperate rain forests of the coastal Pacific Northwest to pinyon-juniper woodlands of the interior southwest62. Although the spatial extent of the FIA plot network represents a large portion of the current range of all species examined in this study (Table 1), substantial portions of some species ranges (e.g., Douglas-fir) extend beyond the study region into Canada and/or Mexico and therefore were not fully addressed here.
The FIA program measures forest attributes on a network of permanent ground plots that are systematically distributed at a rate of ~1 plot per 2428 hectares across the US61. For trees, 12.7 cm DBH and larger, attributes (e.g., species, DBH, live/dead) are measured on a cluster of four 168 m2 subplots61. Trees 2.54–12.7 cm DBH are measured on a microplot (13.5 m2) contained within each subplot, and rare events such as very large trees are measured on an optional macroplot (1012 m2) surrounding each subplot61. In the event a major disturbance (i.e., >1 acre in size, resulting in mortality or damage to >25% of trees) has occurred between measurements on a plot, FIA field crews record the primary disturbance agent (e.g., fire) and estimated year of the event. In the western US, one-tenth of ground plots are measured each year, with remeasurements first occurring in 2011. Please see Data Availability for more information on forest inventory data accessibility.
Forest stability index
Allometric relationships between size and density of live trees make it difficult to interpret many indices of forest change19. Live tree density is expected to decline as trees grow in size, owing to increased individual demand for resources and growing space (i.e., competition)16,23. The expected magnitude of change in tree density, given some change in average tree size, varies considerably across forest communities63, site conditions64, and stand age classes23. Thus, we posit it is useful to contextualize observed changes in live tree density relative to those expected given shifts in average tree size within a stand. To this end, we developed the FSI, a measure of change in relative live tree density that can be applied in stands of any forest community and/or structural type.
To compute the FSI, we first develop a model of maximum size-density relationships for tree populations in our study system (Fig. 6). This model describes the theoretical maximum live tree density (({N}_{max }); in terms of tree number per unit area) attainable in stands as a function of their average tree size ((overline{S})) and will be used as a reference curve to determine the proportionate live tree density of observed stands (i.e., observed density with respect to theoretical maximum density). We use average tree basal area as an index of tree size (one, however, could also use biomass, volume, or other indices of tree size). For stand-type i, the general form of the maximum tree size-density relationship is given by

$${N}_{max }({bar{S}}_{i})={a}_{i}cdot {bar{S}}_{i}^{ {r}_{i}},$$
(1)

where a is a scaling factor that describes the maximum tree density at (bar{S}=1) and r is a negative exponent controlling the decay in maximum tree density with increasing average tree size. Such allometric size-density relationships (i.e., power functions) are widely accepted as quantitative law describing the behavior of even-aged plant populations under self-thinning conditions16,17, and have been used extensively to describe relative stand density in forests23,24. As detailed below, we allow both a and r to vary with stand-type i, as maximum size-density relationships have been shown to vary across forest communities and ecological settings63,65. Allowing a and r to vary by forest community type, for example, allows us to acknowledge that the maximum tree density attainable in a lodgepole pine stand is likely to differ from that of a pinyon-juniper stand with the same average tree size.
Fig. 6: Maximum size-density relationship for an example stand-type.

Individual points represent observed stand-level indices of tree density (N) and average tree size ((overline{S})). Maximum tree density (({N}_{max })) is modeled as a power function of average tree size within a stand. Here, we use quantile regression to estimate ({N}_{max }) as the 99th percentile of N conditional on (overline{S}). The resulting maximum size-density curve can then be used to compute the relative density of observed stands (RD), where relative density is defined as ratio of observed tree density (N) to maximum theoretical density (({N}_{max })), given (overline{S}). Source data are provided as a Source Data file.

Full size image

We next define an index of the relative density of a population of trees j (e.g., species, Pinus edulis) within a stand of type i (e.g., forest community type, pinyon/juniper woodland)

$${{rm{RD}}}_{ij}=sum _{h=1}frac{{N}_{hij}}{{N}_{max }({S}_{hi})},$$
(2)

where N is the density represented by tree h (in terms of tree number per unit area), and S is an index of individual-tree size (e.g., basal area, as used here). The denominator of Equation (2) represents the maximum tree density attainable in a stand of type i with average tree size equal to the size of tree h. We therefore express the relative density of a population j within stand-type i as a sum of the relative densities represented by individual trees within the stand. RD can be interpreted as the proportionate density, or stocking, of a population of trees within stand, where values range from 0 (population j is not present within a stand) to 1 (population j constitutes 100% of a stand and the stand is at maximum theoretical density given its size distribution). As we do in this study, one may apply any range of estimators to summarize the expected relative density of a population of trees j across a range of different stand-types (e.g., estimate the mean and variance of RDj across a broad region containing many different stand-types).
It is important to note that Equation (2) is approximately equal to a simpler method using aggregate indices (i.e., (frac{{sum }_{h = 1}{N}_{hij}}{{N}_{max }(overline{{S}_{i}})})) when tree size-distributions are normally distributed (even age-structures). However, the use of aggregate indices introduces class aggregation bias that results in overestimation of relative density in stands with non-normal size distributions (i.e., uneven age-structures), consistent with other indices of relative tree density66. In contrast, summing tree-level relative densities eliminates such bias and allows RD to accurately compare density conditions across stands in very different structural settings (e.g., even-aged plantation vs. irregularly structured old forest). Furthermore, the partitioning of relative density into tree-level densities allows RD to be accurately summarized within tree size-classes66. That is, it is possible to explicitly estimate the contribution of tree size-classes to overall stand density using RD.
For a given population j within stand-type i, we define the FSI as the average annual change in relative tree density observed between successive measurements of a stand

$${rm{FSI}}=frac{{{Delta }}{rm{RD}}}{{{Delta }}t},$$
(3)

where Δt is the number of years between successive measurement times t1 and t2 and ΔRD is the change in RD over Δt (i.e., ({{rm{RD}}}_{{t}_{2}}-{{rm{RD}}}_{{t}_{1}})). The FSI may also be expressed in units of percent change (%FSI), where average annual change in relative tree density is standardized by previous relative density

$$% {rm{FSI}}=frac{100cdot {rm{FSI}}}{{{rm{RD}}}_{{t}_{1}}}.$$
(4)

Here, stability is defined by zero net change in relative tree density over time (i.e., FSI equal to zero), but does not imply zero change in absolute tree density or tree size distributions. For example, a population exhibiting a decrease in absolute tree density (e.g., trees per unit area) may be considered stable if such decline is offset by a compensatory increase in average tree size. However, populations exhibiting expansion (i.e., ({{rm{RD}}}_{{t}_{1}} {,}{{rm{RD}}}_{{t}_{2}})) in relative tree density will be characterized by positive and negative FSI values, respectively.
Statistical analysis
We computed the FSI for all remeasured FIA plots in the western US (N = 24,229). We included plots on both public and private lands and considered all live stems (DBH ≥2.54 cm) in our analysis. As forest management can effect regional shifts in tree density, we excluded plots with evidence of recent (i.e., within 5 years of initial measurement) silvicultural treatment (e.g., harvesting, artificial regeneration, site preparation). All plot measurements occurred from 2001 to 2018, with an average remeasurement interval of 9.78 years (±0.005 years). For brevity, we restricted our analysis to consider the eight most abundant tree species in the western US. We identified the most abundant tree species using the rFIA R package60, defining abundance in terms of estimated total number of trees (DBH ≥ 2.54 cm) in the year 2018. We excluded species that exhibit non-tree growth habits (i.e., shrub, subshrub) across portions of the study region. All statistical analysis was conducted in Program R (4.0.0)67.
We developed a Bayesian quantile regression model to estimate maximum size-density relationships for stand-types observed within our study area. Here, we use TPH as an index of absolute tree density, average tree basal area ((overline{{rm{BA}}}); equivalent to tree basal area per hectare divided by TPH) as an index of average tree size, and forest community type to describe stand-types. We produced stand-level estimates of TPH and (overline{{rm{BA}}}) from the most recent measurements of FIA plots that (1) lack evidence of recent (within remeasurement period or preceding 5 years) disturbance and/or silvicultural treatment and (2) exhibit approximately normal tree diameter distributions (i.e., even-aged). Here we define an approximately normal tree diameter distribution as exhibiting Pearson’s moment coefficient of skewness between −1 and 1.
We transform the nonlinear size-density relationship to a linear function by taking the natural logarithm of TPH and (overline{{rm{BA}}}), and use a linear quantile mixed-effects model to estimate the 99th percentile of TPH conditional on (overline{{rm{BA}}}) (i.e., in log-log space) for all observed forest community types. We allowed both the model intercept and coefficient to vary across observed forest community types (i.e., random slope/intercept model), thereby acknowledging variation in the scaling factor (a) and exponent (r) of the maximum tree size-density relationship across stand-types. We place informative normal priors on the model intercept (μ = 7, σ = 1) and coefficient (μ = 0.8025, σ = 0.1) following the results of decades of previous research in maximum tree size-density relationships16,23,63,65.
The FIA program uses post-stratification to improve precision and reduce non-response bias in estimates of forest variables68, and we used these standard post-stratified estimators to estimate the mean and variance of the FSI for each species across their respective ranges within the study area (see Code Availability for all relevant code). Further, the FIA program uses an annual panel system to estimate current inventories and change, where inventory cycles consist of multiple panels, and individual panels are comprised of mutually exclusive subsets of ground plots measured in the same year within a region. Precision of point and change estimates can often be improved by combining annual panels within an inventory cycle (i.e., by augmenting current data with data collected previously). We used the simple moving average estimator implemented in the rFIA R package60 to compute estimates from a series of eight annual panels (i.e., sets of plots remeasured in the same year) ranging from 2011 to 2018. The simple moving average estimator combines information from annual panels with equal weight (i.e., irrespective of time since remeasurement), thereby allowing us to characterize long-term patterns in relative density shifts. We determine populations to be stable if the 95% confidence intervals for range-averaged FSI included zero. Alternatively, if confidence intervals of range-averaged FSI do not include zero, we determine the population to be expanding when the estimate is positive and declining when the estimate is negative.
To identify changes in species-size distributions, we used the simple moving average estimator to estimate the mean and variance of the FSI by species and size class across the range of each species within our study area. We assign individual trees to size-classes representing 10% quantiles of observed diameter distributions (i.e., diameter at 1.37 m above ground) of each species growing on one of seven site productivity classes (i.e., inherent capacity of a site to grow crops of industrial wood). That is, we allow size class definitions to vary among species and along a gradient of site productivity, thereby acknowledging intra-specific variation in diameter distributions arising from differences in growing conditions. The use of quantiles effectively standardizes absolute size distributions, simplifying both intra-specific and inter-specific comparison of trends in relative density shifts along species-size distributions.
We assessed geographic variation in species relative density shifts at two scales: ecoregion divisions and subsections69. Ecoregion divisions (shown for our study area in Fig. 1) are large geographic units that represent broad-scale patterns in precipitation and temperature across continents. Ecoregion subsections are subclasses of ecoregion divisions, differentiated by variation in climate, vegetation, terrain, and soils at much finer spatial scales than those represented by divisions. We again used the simple moving average estimator to estimate the mean and variance of the FSI by species within each areal unit (i.e., drawing from FIA plots within each areal unit to estimate mean and variance of the FSI). As a direct measure of changes in relative tree density, spatial variation in the FSI is indicative of spatial shifts in species distributions during the remeasurement interval (i.e., range expansion/contraction and/or within-range relative density shifts). That is, the distribution of populations shift toward regions increasing in relative density and away from regions decreasing in relative density during the temporal frame of sampling. We map estimates of the FSI for each areal unit to assess spatial patterns of changes in relative density and identify regions where widespread geographic shifts in species distributions may be underway.
We sought to quantify the average effect of forest disturbance processes on changes in the relative density of top tree species in the western US over the interval 2001–2018. To this end, we developed a hierarchical Bayesian model to determine the average severity and annual probability of disturbances (i.e., wildfire, insect outbreak, and disease) on sites where each species occurs. Average severity was modeled as

$${y}_{jk} sim {rm{normal}}({alpha }_{j}+sum _{l}{beta }_{jl}cdot {x}_{lk}, {varsigma }_{j}^{2}),$$
(5)

where yjk is the FSI of species j on plot k, αj is a species-specific intercept, βjl is a species-specific coefficient corresponding to the binary variable xlk that takes the value of 1 if disturbance l occurred within plot k measurement interval and 0 otherwise. The intercept and regression coefficients each received an uninformative normal prior distribution. The species-specific residual standard deviation ςj received a uninformative uniform prior distribution70.
On average, disturbance will occur at the midpoint of plot remeasurement periods, assuming temporal stationarity in disturbance probability over the study period. As plots in this study are remeasured on 10-year intervals, we assume that tree populations have, on average, 5 years to respond to any disturbance event. Hence, our definition of disturbance severity, βjl’s, cannot be interpreted as the immediate change in relative tree density resulting from disturbance. Rather, disturbance severity (as defined here) includes the immediate effects of disturbance, as well as 5 years of change in relative tree density prior to and following disturbance (where disturbance is assumed to be functionally instantaneous).
Annual probability of disturbance l on plot k was modeled as

$${x}_{lk} sim {rm{binomial}}({{Delta }}{t}_{k},{psi }_{jl}),$$
(6)

where Δtk is the number of years between successive measurements of plot k, viewed here as the number of binomial “trials,” and ψjl is the species-specific probability for disturbance which was assigned a beta(1,1) prior distribution. Hence, annual probability of disturbance is assumed to vary by species j and by disturbance type l.
We estimate the mean effect of forest disturbance processes on changes in species-specific relative tree density by multiplying the posterior distributions of βjl and ψjl. That is, we multiply species-specific disturbance severity by disturbance probability to yield an estimate of the mean change in relative density caused by disturbance over the study period. We then standardize these values across species by dividing by the average relative density of each species at the beginning of the study period. Thus, standardized values can be interpreted as the annual proportionate change in the relative tree density of each species resulting from disturbance over the period 2001–2018.
Reporting summary
Further information on research design is available in the Nature Research Reporting Summary linked to this article. More

Six sponge species, R. odorabile, C. matthewsi, C. foliascens, S. flabelliformis, I. ramosa and C. orientalis (a bioeroding sponge), were selected for metagenomic sequencing (7 ± 0.5 Gbp) as these species represent dominant habitat forming taxa on tropical and temperate Australian reefs and exhibit high intraspecies similarity in their microbiomes. In addition, previously published microbial MAGs from I. ramosa and Aplysina aerophoba were analysed [8, 12], including 62 additional unpublished MAGs from A. aerophoba. The recovered MAGs, averaging 86 ± 12% completeness and 2 ± 2% contamination, made up 72 ± 21% relative abundance of their respective communities (by read mapping) on average and spanned the vast majority of microbial lineages typically seen in marine sponges [45] (Fig. S1 and Table S1), including the bacterial phyla Proteobacteria (331 MAGs), Chloroflexota (242), Actinobacteriota (155), Acidobacteriota (97), Gemmatimonadota (60), Latescibacterota (44; including lineages Anck6, PAUC34 and SAUL), Cyanobacteria (43), Bacteroidota (38), Poribacteria (35), Dadabacteria (22; including SBR1093), Nitrospirota (22), Planctomycetota (15), UBP10 (14), Bdellovibrionota (13), Patescibacteria (9; includes Candidate Phylum Radiation), Spirochaetota (8), Nitrospinota (7), Myxococcota (4), Entotheonella (2) and the archaeal class Nitrososphaeria (21; phylum Crenarchaeota), hereafter referred to by their historical name “Thaumarchaeota” for name recognition. Mapping of the metagenomic reads to the recovered MAGs showed that the communities had high intraspecies similarity across replicates, consistent with previous 16S rRNA gene-based analyses (Fig. S1). In general, taxa present in A. aerophoba, C. foliascens, C. orientalis and S. flabelliformis appeared unique to those sponge species, with only one dominant lineage present in C. orientalis (order Parvibaculales). In contrast, several Actinobacteriota, Acidobacteriota and Cyanobacteria populations were shared across C. matthewsi, R. odorabile and I. ramosa. Further, members of the Thaumarchaeota were detected in all sponge species and were particularly abundant in S. flabelliformis at 12 ± 4% relative abundance (Fig. S1). Addition of these sponge MAGs to genome trees comprising all publicly available sponge symbionts (N = 1188 MAGs) resulted in a phylogenetic gain of 44 and 75% for Bacteria and Archaea, respectively, reflecting substantial novel genomic diversity (Fig. 1).
Comparative genomic analysis of the sponge-derived MAGs provided unique insights into the distribution of metabolic pathways across sponge symbiont taxa. For example, microbial oxidation of ammonia benefits the sponge host by preventing ammonia from accumulating to toxic levels [46], a process thought to be mediated by both symbiotic Bacteria and Archaea (i.e. Thaumarchaeota) [33]. Prior identification of ammonia oxidisers has been based on functional inference from phylogeny (16S rRNA gene amplicon surveys) [47] or homology to specific Pfams (metagenomes) [33]. However, the CuMMO gene family is diverse, encompassing functionally distinct relatives that include amoA, particulate methane monooxygenases and hydrocarbon monooxygenases that cannot be distinguished by homology alone [35]. We used GraftM [32] to recover CuMMO genes from the sponge MAGs and their metagenomic assemblies, as well as previously sequenced metagenomic assemblies from six additional sponge microbiomes where bacterial amoA gene sequences had been identified [33]. Phylogenetic analysis of the recovered CuMMO genes showed that all archaeal homologues came from Thaumarchaeota and fell within the archaeal amoA clade. In contrast, bacterial CuMMO sequences were identified exclusively in MAGs from the phylum UBP10 (formerly unclassified Deltaproteobacteria) and from an unknown taxonomic group in the previous metagenomic assemblies [33]. All recovered bacterial and taxonomically unidentified CuMMO placed within the Deltaproteobacteria/Actinobacteria hmo clade, indicating these genes are specific for hydrocarbons rather than ammonia (Fig. S2). The finding that Thaumarchaeota are the only microbes within any of the surveyed sponge species capable of oxidising ammonia, and their ubiquity across sponges, suggests they are a keystone species for this process.
To further investigate the distribution of functions within the sponge microbiome, a set of highly complete ( >85%) sponge symbiont MAGs were grouped by principal components analysis based on their KEGG and Pfam annotations, as well as orthologous clusters that reflected all gene content. Similar analysis conducted on 37 MAGs from the sponge Aplysina aerophoba suggested the presence of functional guilds, with MAGs from disparate microbial phyla carrying out similar metabolic processes [12] (e.g. carnitine catabolism). Here, we find that MAGs clustered predominately by microbial taxonomy (phylum) rather than function in all three analyses (Fig. S3). While functional guilds could not be identified based on analysis of total genome content, this does not preclude the existence of such guilds based on more specific metabolic pathways.
To identify pathways enriched within the sponge microbiome, sponge-associated MAGs with >85% completeness (N = 798) were compared with a set of coral reef and coastal seawater MAGs (N = 86), 31 derived from published datasets [31] and 55 from this study (Table S1). Seawater MAGs with >85% genome completeness (93 ± 4% completeness and 2 ± 2% contamination; Table S1) spanned the bacterial phyla Proteobacteria (48 MAGs), Bacteroidota (13), Planctomycetota (5), Myxococcota (5), Gemmatimonadota (3), Marinisomatota (3), Actinobacteriota (3), Verrucomicrobiota (2), Cyanobacteriota (2), Bdellovibrionota (1) and the archaeal phylum Nanoarchaeota (1). Comparative analysis revealed that sponge symbionts were enriched in metabolic pathways for carbohydrate metabolism, defence against infection by MGE, amino acid synthesis, eukaryote-like gene repeat proteins (ELRs) and cell–cell attachment (Tables S2–S4).
Genes belonging to GH and carbohydrate esterase (CE) families (Table S2) acting on starch (GH77), arabinose (CAZY families GH127 and GH51), fucose (GH95 and GH29) and xylan polymers (CE7 and CE15), were enriched in sponge-associated lineages, likely reflecting the hosts critical role in catabolising dissolved organic matter (DOM) present in reef seawater (Fig. 2). Microbial GHs from the GH77 family target starch, the main sugar storage compound in marine algae [48], whereas GHs from families 51 and 127 are known to act on plant arabinosaccharides, such as the hydroxyproline-linked arabinosaccharides found in algal extensin glycoproteins [49, 50]. GH127 enzymes are also required for microbial degradation of carrageenan, a complex heteropolysaccharide produced by red algae [51]. Members of the fucosidase GH95 and GH29 enzyme families are known to degrade fucoidan, a complex fucosaccharide prominent in brown algae [50, 52]. Notably, arabino- and fucopolysaccharides also make up a significant proportion of coral mucus, a major component of DOM in coral reefs that sponges have been shown to utilise [53, 54]. Supporting this observation, isotopic investigation of the fate of coral mucus and algal polysaccharides in sponges showed that the microbiome participates in metabolism of these compounds, particularly in sponges with high microbial abundance and diversity [4, 5]. Enzymes from the CE families 15 and 7 have been primarily characterised in terrestrial plants where they act as glucuronyl esterases and acetyl-xylan esterases, degrading lignocellulose and removing acetyl groups from hemicellulose [55] (e.g. xylans). Characterisation of CE15 and CE7 from marine microbes is rare, though activity on xylans, which are a structural component of marine algae, has previously been demonstrated [55,56,57].
Fig. 2: Phylogenetic tree showing the distribution of glycosyl hydrolases and esterases across MAGs with >85% completeness (N = 884).

Values represent the copy number of each gene per MAG. Internal branches of the tree are coloured by phylum, while the outer strip is coloured by class. Both are listed clockwise in the order in which they appear. Seawater MAGs are denoted by grey labels with red text.

Full size image

GHs acting on sialic acids (GH33) and glycosaminoglycans (GH88) were also enriched in the sponge-associated MAGs and may act on compounds found within sponge tissue [13] (Fig. 2). In contrast, no genes for the degradation of collagen (collagenases), one of the main structural components of the sponge skeleton were identified. Sialic acid-linked residues are found in the sponge mesohyl [58], and although the impact of cleavage on the host is unknown, analogy can be made to other symbioses. For example, sialidases are common in the commensal bacteria present in the human gut where they are used to cleave and metabolise the sialic acid-containing mucins lining the gut wall [59]. Increased sialidase activity is associated with gut dysbiosis and inflammation [60] and careful control of sialidase-containing commensals is therefore necessary to maintain gut homoeostasis [59]. As glycosaminoglycans are also part of sponge tissue [13, 61], the same may apply to microorganisms encoding GH88 family enzymes. However, these genes are also implicated in the degradation of external sugar compounds, such as ulvans, a major sugar storage compound found in green algae that can make up to 30% of their dry weight [62]. Thus, the ecological role of GH88 family enzymes within the sponge microbiome requires further investigation.
Enrichment of GHs and CEs was largely restricted to the Poribacteria, Latescibacteria (class UBA2968), Spirochaetota, Chloroflexota (classes UBA2235 and Anaerolineae, but not Dehalococcoides) and Acidobacteriota (class Acidobacteriae). These findings corroborate previous targeted genomic characterisations of the Chloroflexota and Poribacteria [13, 14] but show that they are part of a larger set of polysaccharide-degrading lineages. Identification of disparately related microbial taxa across several sponge lineages (Figs. 1 and 2) that encode similar pathways for polysaccharide degradation, and therefore occupy a similar ecological niche, supports the existence of functional guilds within the sponge microbiome when viewed at the level of individual pathways. Given the fundamental role of marine sponges in recycling coral reef DOM, studies targeting these specific guilds are needed to quantify their contribution to reef DOM transformation.
Because sponges filter and retain biomass from an extensive range of reef taxa (eukaryotic algae, bacteria, archaea, etc), they are exposed to a greatly expanded variety of MGEs from these organisms, including viruses, transposable elements and plasmids [33, 63]. For this reason, sponge-associated microorganisms likely require a diverse toolbox of molecular mechanisms for resisting infection. Both RM and CRISPR systems are capable of recognising and cleaving MGEs as part of the bacterial immune repertoire. RM systems are part of the innate immune system of bacteria and archaea and are encoded by a single (Type II) or multiple proteins (Type I, III and IV) that recognise and cleave foreign DNA based on a defined target sequence. In contrast, CRISPR systems are part of the adaptive immune system of some bacteria and archaea and encode a target sequence derived from the genome of a previous infective agent that is used by a CRISPR-associated protein (CAS) to identify and cleave foreign DNA. RM (Fig. S4) and CAS (Fig. S5) genes were both enriched (Table S3) in the sponge-associated MAGs and relatively evenly distributed across taxa, with the exception of the Planctomycetota and Verrucomicrobiota, where they were largely absent. As these MAGS average 93 +/− 5% completeness, this result is not likely due to genome incompleteness. This finding contrasts with comparative investigations of Planctomycetota genomes from other environments [64] and additional research is required to ascertain the mechanisms used by sponge-associated Planctomycetota and Verrucomicrobiota to avoid infection. Although Type III RM genes were enriched in sponge MAGs, they were also present in all seawater MAGs. In contrast, Types I and II RM genes were present almost exclusively in the sponge-associated MAGs. In conjunction with an enrichment in CRISPR systems, this expanded repertoire of defence systems likely reflects the increased burden from MGEs associated with the hosts role in filtering and concentrating diverse sources of reef biomass. Supporting this hypothesis, metagenomic surveys of sponge-associated viruses revealed a more diverse viral population than what could be recovered from the surrounding seawater [63]. Further, we found that genes encoding toxin-antitoxin systems, which are present on MGEs, such as plasmids, were also enriched in sponge-associated MAGs. These observations suggest that RM and CRISPR systems are important features of microbe-sponge symbiosis, allowing the symbionts to colonise and persist within their host by avoiding viral infection or being overtaken by MGEs.
Pathways for the synthesis of amino acids were also enriched in the sponge microbiome. The inability of animals to produce several essential amino acids has been proposed as a primary reason that they harbor microbial symbionts [65,66,67,68] and it has long been thought that sponges acquire at least some of their essential amino acids from their microbiome [69, 70]. Further, gene-centric characterisation of the Xestospongia muta and R. odorabile microbiomes revealed pathways to synthesise and transport essential amino acids [33, 70]. However, these same amino acid pathways are also used catabolically by the microorganisms, and transporters could simply be importing amino acids into the microbial cell. Further, as sponges are almost constantly filter feeding, essential amino acids could be acquired through consumption of microorganisms present in seawater. Comparison of sponge MAGs with those from seawater revealed enrichment of specific pathways for the synthesis of lysine, arginine, histidine, threonine, valine and isoleucine (Table S4). However, visualisation of the distribution of these genes revealed that almost all MAGs in both sponges and seawater produce all amino acids, though specific lineages may use different pathways to achieve this (Fig. S6). The enrichment observed in the sponge MAGs was therefore ascribed to differences in pathway completeness between sponge-associated and seawater microbes, rather than an enhanced ability of sponge symbionts to produce any specific amino acid. In contrast, compounds, such as taurine, carnitine and creatine have also been proposed as important host-derived carbon sources for symbionts [69], but pathways for their catabolism were enriched in seawater rather than sponge-associated MAGs. While these findings do not invalidate the possibility that microbial communities play a role in amino acid provisioning to the host or that they utilise host-derived taurine, carnitine, or creatine, they suggest that these are not key processes mediating microbe-sponge symbiosis.
To form stable symbioses, bacteria must persist within the sponge tissue and avoid phagocytosis by host cells. Microbial proteins containing ELR motifs have been identified in a range of animal and plant-associated microbes and are thought to modulate the host’s intracellular processes to facilitate stable symbiotic associations [71, 72]. For example, ELR-containing proteins from sponge-associated microbes have been shown to confer the ability to evade host phagocytosis when experimentally expressed in E. coli [10, 73]. ELR-containing proteins from the ankyrin (ARP), leucine-rich, tetratricopeptide and HEAT repeat families were enriched in the sponge-associated MAGs. In contrast, WD40 repeats were not found to be enriched but are included here as they have previously been reported as abundant in Poribacteria and symbionts of other marine animals [13, 31]. Most ELRs were present across all taxa but were much more prevalent in specific lineages (Fig. 3). For example, sponge-associated Poribacteria, Latescibacterota and Acidobacteriota encoded a high proportion of all ELR types, while other lineages, such as the Gemmatimonadota (average 0.25% coding genes per sponge-associated MAG versus 0.09% in seawater MAGs), Verrucomicrobiota (2%), Deinococcota (0.85%), Acidobacteriota (0.20%; specifically class Luteitaleia at 0.55%) and Dadabacteria from C. orientalis (0.62%) encoded a comparatively high percentage of ARPs and Nitrospirota encoded a high percentage of HEAT_2 family proteins (0.55% versus 0.05% in seawater MAGs) relative to other taxa. In contrast, ELR abundances were substantially lower, or absent, in the Actinobacteriota, the class Bacteroidia within the phylum Bacteroidota, and the Thaumarchaeota, suggesting these microorganisms utilise alternative mechanisms to maintain their stable associations with the host.
Fig. 3: Phylogenetic tree showing the distribution of eukaryote-like repeat proteins—ankyrin (ARP), leucin-rich (LRR), tetratricopeptide (TPR), HEAT and WD40—across MAGS with >85% completeness (N = 884).

Values represent the percentage of coding genes per MAG devoted to each gene class. Internal branches of the tree are coloured by phylum, while the outer strip is coloured by class, and both are listed clockwise in the order in which they appear. MAGs from seawater are denoted by grey labels with red text.

Full size image

The mechanisms by which ELRs interact with sponge cells remains largely unknown, although microbes in other host systems are known to deliver ELR-containing effector proteins into host cells via needle-like secretion systems (types III, IV and V) or extracellular contractile injection systems [74, 75], where they interact with the cellular machinery of the host to modify its behaviour. In sponges, it is also possible that ELRs could be secreted into the extracellular space by type I or II secretion systems. Interestingly, although most sponge MAGs encoded eukaryote-like proteins (Fig. 3), few lineages encoded the necessary genes to form secretion systems (Fig. S7). It is therefore unlikely that ELRs are introduced to the sponge host via traditional secretion pathways used in other animal-symbiont systems.
Maintaining stable association with the sponge may also require mechanisms for attachment to the host tissue. For example, cadherin domains are Ca2+-dependent cell–cell adhesion proteins that are abundant in eukaryotes and have been found to serve the same function in bacteria [76]. Similarly, fibronectin III domains mediate cell adhesion in eukaryotes, but also occur in bacteria where they play various roles in carbohydrate binding and biofilm formation [77, 78]. In addition, some bacterial pathogens utilise fibronectin-binding proteins to gain entry into host tissue by binding to host fibronectin [77, 78]. Genes containing cadherin domains were enriched in the sponge-associated MAGs and were identified in most bacterial lineages, but were notably absent in the Cyanobacteriota and Verrucomicrobiota (Fig. 4). Genes containing fibronectin III domains and those for fibronectin-binding proteins were also enriched in sponge-associated MAGs and were distributed across most lineages, though were particularly abundant in the Actinobacteriota and Chloroflexota. However, although fibronectin III-containing genes were taxonomically widespread, those encoding fibronectin-binding proteins were restricted to the phyla Poribacteria, Gemmatimonadota, Latescibacterota, Cyanobacteriota, class Anaerolineae within the Chloroflexota (but not Dehalococcoidia), class Rhodothermia within the Bacteroidota, Spirochaetota, Nitrospirota and the archaeal phylum Thaumarchaeota. Interestingly, the taxonomic distribution of these genes shares significant overlap with lineages encoding the genes for sponge sialic acid and glyosaminoglycans degradation, suggesting that attachment to the host may be necessary for utilisation of these carbohydrates (Fig. 2). However, as the host, bacterial, and archaeal components of the sponge holobiont have fibronectin III domains, symbionts encoding fibronectin-binding proteins may use these to adhere to the host tissue or potentially to form biofilms (bacteria–bacteria attachment). In either case, the enrichment and wide distribution of cadherins, fibronectins and fibronectin-binding proteins in the sponge MAGs suggests that cell–cell adhesion is critical for successful establishment in the sponge niche.
Fig. 4: Phylogenetic tree showing the distribution of cadherins, fibronectins and fibronectin-binding proteins across MAGS with >85% completeness (N = 884).

Values represent the copy number of each gene per MAG. Internal branches of the tree are coloured by phylum while the outer strip is coloured by class. Both are listed clockwise in the order in which they appear. Seawater MAGs are denoted by grey labels with red text.

Full size image

Distribution of genes encoding ELRs, polysaccharide-degrading enzymes (GHs and CEs), cadherins, fibronectins, RMs and CRISPRs across distantly related taxa suggests that they were either acquired from a common ancestor or that they represent more recent LGT events, potentially mediated by MGEs, which are enriched in sponge-associated microbial communities [69]. Here, we identify 4963 LGTs from five sponges for which sufficient sequence data were available ( >100 Mbp total sequence length across all MAGs), as well as 136 LGTs from seawater MAGs, averaging 1.64 and 0.52 LGTs per Mbp sequences, respectively (Fig. 5 and Table S5). Sequence similarity of LGTs from MAGs within a sponge species was higher than between sponge species, indicating relatively recent gene transfers (Fig. S8). A higher frequency (Fig. S9) and lower genetic divergence of LGTs among MAGs derived from the same sponge species likely results from the close physical distance between members of each microbiome, as has been observed in other host-symbiont systems [79]. The identification of lateral transfers between microbes from different sponge species may highlight the horizontal acquisition of these microbes or that a recent ancestor inhabited the same host. Notably, LGTs included a subset of genes that were enriched within the sponge-associated MAGs, such as GH33 (sialidases) and CE7 (acetyl-xylan esterases), attachment proteins (cadherins and fibronectin III), RM and CAS proteins, and members of all ELR families other than WD40 (Figs. 6 and S10). The observation that a significant number of sponge-enriched genes were laterally transferred between disparate microbial lineages suggests that the processes they mediate provide a strong selective advantage within the sponge niche, though further research is required to validate these findings.
Fig. 5: Visualisation of LGTs detected within the MAGs for the five sponges passing the cumulative MAG length criteria ( >100 Mbp).

The inner strip is coloured by phylum while the outer strip is coloured by host sponges. Bands connect donors and recipients, with their colour corresponding to the donors and the width correlating to the number of LGTs.

Full size image

Fig. 6: Visualisation of gene flow among microbial phyla for gene families enriched in sponge-associated MAGs.

The inner ring and band connecting donor and recipient is coloured by protein family of the gene being transferred, with the width of the band correlating to the number of LGTs. Recipient MAGs are shown in grey. The outer ring is coloured by microbial phylum. Representation of RM and CAS gene LGTs can be found in Fig. S10.

Full size image

Sponges are important constituents of coral reef ecosystems because of their critical role in DOM cycling and retention via the sponge-loop. Despite their importance, functional characterisation of sponge symbiont communities has been restricted to just a few lineages of interest, potentially biasing our view of sponge symbiosis. Here we present a comprehensive characterisation of sponge symbiont MAGs spanning the complete range of taxa found in marine sponges (Fig. 7), most of which were previously uncharacterised. We revealed enrichment in glycolytic enzymes (GHs and CEs) reflecting specific functional guilds capable of aiding the sponge in the degradation of reef DOM. Further, we identified several ELRs, CRISPRs and RMs that likely facilitate stable association with the sponge host, showing the specificity of ELR types with individual microbial lineages. We also clarified the role of Thaumarchaeota as a keystone taxon for ammonia oxidation across sponge species and showed that processes previously thought to be important, such as amino acid provisioning and taurine, creatine and carnitine metabolism are unlikely to be central mechanisms mediating sponge-microbe symbiosis. Many of the enriched genes are laterally transferred between microbial lineages, suggesting that LGT plays an important role in conferring a selective advantage to specific sponge-associated microorganisms. Taken together, these data illustrate how evolutionary processes have distributed and partitioned ecological functions across specific sponge symbiont lineages, allowing them to occupy or share specific niches and live symbiotically with their sponge hosts.
Fig. 7: Schematic overview of microbial interactions with the host as inferred from the functional potential encoded by the sponge-associated microbial MAGs.

Fbn fibronectin, cdh cadherins, RM restriction-modification systems, CAS CRISPR-associated proteins, ELP eukaryotic-like repeat proteins, CE7 carbohydrate esterase family 7, GH33 glycosyl hydrolase family 33.

Full size image More

This announcement sent shockwaves around the globe, largely lauding it as an important step in the right direction4. China’s decision is also unprecedented at several levels, which could result in profound and far-reaching impacts for both humans and wildlife.
First, this decision was initiated and adopted by China’s highest legislature — the Standing Committee of the National People’s Congress — with the explicit endorsement of President Xi1,2. In contrast, China’s response to the SARS outbreak in 2003 — a short-lived ban on the trade and consumption of palm civets5 — was initiated by various government agencies at lower legislative levels. In responding to the COVID-19 pandemic, the political will in China has never been stronger or more overt across multiple levels of government6. Within a few months, all 31 provinces in China have published provincial legislation on wildlife farming and consumption. Perhaps more importantly, from May to July, the People’s Congress standing committees at the national and provincial levels conducted nation-wide evaluations on the effectiveness of these policies and their enforcement. The committees concluded that policies were generally well implemented, but there was room for improvement on certain aspects, including finding alternative livelihoods for affected wildlife traders, continuing to revise the protected species lists, and addressing loopholes in wildlife trade monitoring and habitat conservation. In terms of positive outcomes, joint actions and special operations from the government have closed 12,000 wildlife-related businesses, intensified monitoring efforts to include over four million e-commerce platforms, and removed 990,000 online sources of information associated with wildlife trade6.
Second, China’s current decision includes a series of new legislations to build on the achievements of current actions by enhancing the regulation of wildlife farms and markets (Box 1). The revision of the country’s Wildlife Protection Law is expected to bring about long-term and systematic changes to wildlife conservation. Additionally, China is also revising its List of Protected Animals. Species threatened by consumption, such as the pangolin and yellow-breasted bunting, are being promoted to the highest protection level (Class I Wildlife species)7. Furthermore, China’s Ministry of Agriculture and Rural Affairs published an updated Catalogue of Animal Genetic Resource in May 2020. Among wild animals, only species in this catalogue can be farmed or consumed8. There are 64 species of wild animals that are being farmed for consumption, but are not yet included in this catalogue for various reasons (for example, to reduce the risk of sourcing animals from the wild). Nevertheless, the Ministry of Forestry and Grassland has categorized them into two groups: the farming of 45 species (for example, bamboo rat and civet cat) is due to be banned by the end of 2020, and the remaining 19 species (for example, several species of snakes) are allowed to be farmed for non-consumption uses9. Furthermore, the disbursement of government financial compensation to the farmers affected by these new legislations, amounting to over a billion US dollars, is expected to be completed by the end of 2020.
Third, the current ban is likely to galvanize rapid and widespread knock-on actions and impacts. For example, Guangdong has already banned wild vertebrate animals as pets10. The consumption of dogs and cats is banned in the city of Shenzhen11. Pangolin scales are removed as a key ingredient in traditional Chinese medicine, although it is still included as an ingredient in patent medicines in the 2020 Chinese Pharmacopoeia12,13. The pre-COVID wildlife management system in China, especially in its wildlife farming industry, has long been criticized by conservationists to be disordered and outdated. The lack of incentives and capacities had been a major barrier to change for government agencies in the country. The COVID-19 pandemic, at great human and economic costs, has mainstreamed the discourse of wildlife conservation for human well-being, clarified legislations on what species can be farmed, and provided a policy framework for systematic and enforceable wildlife management and conservation. These actions are exactly what scientists have long called for to minimize the risk of zoonotic disease transmission and outbreaks in the future14. More

1.
Ferguson BA, Dreisbach TA, Parks CG, Filip GM, Schmitt CL. Coarse-scale population structure of pathogenic Armillaria species in a mixed-conifer forest in the Blue Mountains of northeast Oregon. Can J Res. 2003;33:612–23.
Article  Google Scholar 
2.
Fricker MD, Heaton LLM, Jones NS, Boddy L. The mycelium as a network. Microbiol Spectr. 2017;5:1–32.
Google Scholar 

3.
Treseder KK, Lennon JT. Fungal traits that drive ecosystem dynamics on land. Microbiol Mol Biol Rev. 2015;79:243–62.
CAS  PubMed  PubMed Central  Article  Google Scholar 

4.
Smith SE, Read DJ. Mycorrhizal symbiosis. 3rd ed. London, UK: Academic Press; 2008.
Google Scholar 

5.
Maron JL, Marler M, Klironomos JN, Cleveland CC. Soil fungal pathogens and the relationship between plant diversity and productivity. Ecol Lett. 2011;14:36–41.
PubMed  Article  Google Scholar 

6.
García-Guzmán G, Heil M. Life histories of hosts and pathogens predict patterns in tropical fungal plant diseases. New Phytol. 2013.

7.
Lin X, Alspaugh JA, Liu H, Harris S. Fungal morphogenesis. Cold Spring Harb Perspect Med. 2014;5:a019679.
PubMed  Article  CAS  PubMed Central  Google Scholar 

8.
Griffin DH. Fungal physiology. 2nd ed. New York, NY: John Wiley and Sons; 1994.
Google Scholar 

9.
Trewavas A. What is plant behaviour? Plant, Cell Environ. 2009;32:606–16.
Article  Google Scholar 

10.
Karban R. Plant behaviour and communication. Ecol Lett. 2008;11:727–39.
PubMed  Article  Google Scholar 

11.
de Kroon H, Mommer L. Root foraging theory put to the test. Trends Ecol Evol. 2006;21:113–6.
PubMed  Article  Google Scholar 

12.
Novoplansky A. Developmental plasticity in plants: implications of non-cognitive behavior. Evol Ecol. 2002;16:177–88.
Article  Google Scholar 

13.
Lovett Doust L. Population dynamics and local specialization in a clonal perennial (Ranunculus repens): II. The dynamics of leaves, and a reciprocal transplant-replant experiment. J Ecol. 1981;69:757–68.
Article  Google Scholar 

14.
Saiz H, Bittebiere A-K, Benot M-L, Jung V, Mony C. Understanding clonal plant competition for space over time: a fine-scale spatial approach based on experimental communities. J Veg Sci. 2016;27:759–70.
Article  Google Scholar 

15.
Andrews JH. Comparative ecology of microorganisms and macroorganisms. New York, NY: Springer; 1991.
Google Scholar 

16.
Carlile MJ. The success of the hypha and mycelium. In: Gow NAR, Gadd GM, editors. The growing fungus. London: Chapman & Hall; 1995. pp. 3–20.

17.
Boddy L. Saprotrophic cord-forming fungi: meeting the challenge of heterogeneous environments. Mycologia. 1999;91:13.
Article  Google Scholar 

18.
Bielčik M, Aguilar-Trigueros CA, Lakovic M, Jeltsch F, Rillig MC. The role of active movement in fungal ecology and community assembly. Mov Ecol. 2019;7:36.
PubMed  PubMed Central  Article  Google Scholar 

19.
Ritz K, Young IM. Interactions between soil structure and fungi. Mycologist. 2004;18:52–9.
Article  Google Scholar 

20.
Harris K, Young IM, Gilligan CA, Otten W, Ritz K. Effect of bulk density on the spatial organisation of the fungus Rhizoctonia solani in soil. FEMS Microbiol Ecol. 2003;44:45–56.
CAS  PubMed  Article  PubMed Central  Google Scholar 

21.
Otten W, Harris K, Young IM, Ritz K, Gilligan CA. Preferential spread of the pathogenic fungus Rhizoctonia solani through structured soil. Soil Biol Biochem. 2004;36:203–10.
CAS  Article  Google Scholar 

22.
Burges A, Nicholas DP. Use of soil sections in studying the amount of fungal hyphae in soil. Soil Sci. 1961;92:25–9.
Article  Google Scholar 

23.
Dechesne A, Wang G, Gülez G, Or D, Smets BF. Hydration-controlled bacterial motility and dispersal on surfaces. Proc Natl Acad Sci USA. 2010;107:14369–72.
CAS  PubMed  Article  PubMed Central  Google Scholar 

24.
Wolfaardt GM, Hendry MJ, Birkham T, Bressel A, Gardner MN, Sousa AJ, et al. Microbial response to environmental gradients in a ceramic-based diffusion system. Biotechnol Bioeng. 2008;100:141–9.
CAS  PubMed  Article  PubMed Central  Google Scholar 

25.
Otten W, Pajor R, Schmidt S, Baveye PC, Hague R, Falconer RE. Combining X-ray CT and 3D printing technology to produce microcosms with replicable, complex pore geometries. Soil Biol Biochem. 2012;51:53–5.
CAS  Article  Google Scholar 

26.
Deng J, Orner EP, Chau JF, Anderson EM, Kadilak AL, Rubinstein RL, et al. Synergistic effects of soil microstructure and bacterial EPS on drying rate in emulated soil micromodels. Soil Biol Biochem. 2015;83:116–24.
CAS  Article  Google Scholar 

27.
Whitesides GM. The origins and the future of microfluidics. Nature. 2006;442:368–73.
CAS  PubMed  PubMed Central  Article  Google Scholar 

28.
Stanley CE, Grossmann G, Casadevall i Solvas X, DeMello AJ. Soil-on-a-chip: microfluidic platforms for environmental organismal studies. Lab Chip. 2016;16:228–41.
CAS  PubMed  Article  PubMed Central  Google Scholar 

29.
Held M, Kaspar O, Edwards C, Nicolau DV. Intracellular mechanisms of fungal space searching in microenvironments. Proc Natl Acad Sci USA. 2019;116:13543–52.
CAS  PubMed  Article  PubMed Central  Google Scholar 

30.
Soufan R, Delaunay Y, Gonod LV, Shor LM, Garnier P, Otten W, et al. Pore-scale monitoring of the effect of microarchitecture on fungal growth in a two-dimensional soil-like micromodel. Front Environ Sci. 2018;6:68.
Article  Google Scholar 

31.
Schmieder SS, Stanley CE, Rzepiela A, van Swaay D, Sabotič J, Nørrelykke SF, et al. Bidirectional propagation of signals and nutrients in fungal networks via specialized hyphae. Curr Biol. 2019;29:217–28.e4.
CAS  PubMed  Article  Google Scholar 

32.
Aleklett K, Kiers ET, Ohlsson P, Shimizu TS, Caldas VE, Hammer EC. Build your own soil: exploring microfluidics to create microbial habitat structures. ISME J. 2018;12:312–9.
PubMed  Article  Google Scholar 

33.
Veresoglou SD, Wang D, Andrade-Linares DR, Hempel S, Rillig MC. Fungal decision to exploit or explore depends on growth rate. Micro Ecol. 2018;75:289–92.
Article  Google Scholar 

34.
Lange M, Smith, Alexander H. The coprinus ephemerus group. Mycologia. 1953;45:747–80.
Article  Google Scholar 

35.
Hernández-Rodríguez M, Oria-de-Rueda JA, Martín-Pinto P. Post-fire fungal succession in a Mediterranean ecosystem dominated by Cistus ladanifer L. Ecol Manag. 2013;289:48–57.
Article  Google Scholar 

36.
Hughes KW, Petersen RH. Transatlantic disjunction in fleshy fungi III: Gymnopus confluens. MycoKeys. 2015;9:37–63.
Article  Google Scholar 

37.
Garnier-Delcourt M, Reckinger C, Tholl M-T, Turk J. Notes mycologiques luxembourgeoises. IV. Bull Soc Nat Luxemb. 2011;112:39–50.
Google Scholar 

38.
Maynard DS, Bradford MA, Covey KR, Lindner D, Glaeser J, Talbert DA, et al. Consistent trade-offs in fungal trait expression across broad spatial scales. Nat Microbiol. 2019;4:846–53.
CAS  PubMed  Article  Google Scholar 

39.
Held M, Edwards C, Nicolau DV. Probing the growth dynamics of Neurospora crassa with microfluidic structures. Fungal Biol. 2011;115:493–505.
PubMed  Article  Google Scholar 

40.
Dowson CG, Rayner ADM, Boddy L. Spatial dynamics and interactions of the woodland fairy ring fungus, Clitocybe nebularis. N Phytol. 1989;111:699–705.
Article  Google Scholar 

41.
Fukasawa Y, Savoury M, Boddy L. Ecological memory and relocation decisions in fungal mycelial networks: responses to quantity and location of new resources. ISME J. 2020;14:380–8.
PubMed  Article  PubMed Central  Google Scholar 

42.
Held M, Binz M, Edwards C, Nicolau DV. Dynamic behaviour of fungi in microfluidics: a comparative study. In: Proc. SPIE 7182, Imaging, Manipulation, and Analysis of Biomolecules, Cells, and Tissues VII. 2009. pp. 718213. https://doi.org/10.1117/12.822464.

43.
Hanson KL, Nicolau DV, Filipponi L, Wang L, Lee AP, Nicolau DV. Fungi use efficient algorithms for the exploration of microfluidic networks. Small. 2006;2:1212–20.
CAS  PubMed  Article  PubMed Central  Google Scholar 

44.
Falconer RE, Houston AN, Otten W, Baveye PC. Emergent behavior of soil fungal dynamics: influence of soil architecture and water distribution. Soil Sci. 2012;177:111–9.
CAS  Article  Google Scholar 

45.
Schmidt MWI, Torn MS, Abiven S, Dittmar T, Guggenberger G, Janssens IA. et al. Persistence of soil organic matter as an ecosystem property. Nature. 2011;478:49–56.
CAS  PubMed  Article  PubMed Central  Google Scholar 

46.
Nezhad AS, Packirisamy M, Bhat R, Geitmann A. In vitro study of oscillatory growth dynamics of camellia pollen tubes in microfluidic environment. IEEE Trans Biomed Eng. 2013;60:3185–93.
PubMed  Article  PubMed Central  Google Scholar 

47.
Tayagui A, Sun Y, Collings D, Garrill A, Nock V. An elastomeric micropillar platform for the study of protrusive forces in hyphal invasion. Lab Chip. 2017;17:3643–53.
CAS  PubMed  Article  PubMed Central  Google Scholar 

48.
Brand A, Gow NA. Mechanisms of hypha orientation of fungi. Curr Opin Microbiol. 2009;12:350–7.
CAS  PubMed  PubMed Central  Article  Google Scholar 

49.
Pantidou M, Watling R, Gonou Z. Mycelial characters, anamorphs, and teleomorphs in genera and species of various families of Agaricales in culture. Mycotaxon. 1983;17:409–32.
Google Scholar 

50.
Boekhout T, Stalpers J, Verduin SJW, Rademaker J, Noordeloos ME. Experimental taxonomic studies in Psilocybe sect. Psilocybe. Mycol Res. 2002;106:1251–61.
Article  Google Scholar 

51.
Valenzuela E, Garnica S. Pseudohelicomyces, a new anamorph of Psilocybe. Mycol Res. 2000;104:738–41.
Article  Google Scholar 

52.
Gadd GM, Bahri-Esfahani J, Li Q, Rhee YJ, Wei Z, Fomina M, et al. Oxalate production by fungi: Significance in geomycology, biodeterioration and bioremediation. Fungal Biol Rev. 2014;28:36–55.
Article  Google Scholar  More

1.
Gruber, N. et al. The oceanic sink for anthropogenic CO2 from 1994 to 2007. Science 363, 1193–1199 (2019).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 
2.
Friedlingstein, P. et al. Global carbon budget 2019. Earth Syst. Sci. Data 11, 1783–1838 (2019).
ADS  Article  Google Scholar 

3.
Caldeira, K. & Wickett, M. E. Anthropogenic carbon and ocean pH. Nature 425, 365 (2003).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

4.
Feely, R. A. et al. Impact of anthropogenic CO2 on the CaCO3 system in the oceans. Science 305, 362–366 (2004).
ADS  CAS  PubMed  Article  Google Scholar 

5.
Doney, S. C., Fabry, V. J., Feely, R. A. & Kleypas, J. A. Ocean acidification: The other CO2 problem. Ann. Rev. Mar. Sci. 1, 169–192 (2009).
PubMed  Article  PubMed Central  Google Scholar 

6.
Riebesell, U. et al. Reduced calcification of marine plankton in response to increased atmospheric CO2. Nature 407, 364–367 (2000).
ADS  CAS  PubMed  Article  Google Scholar 

7.
Orr, J. C. et al. Anthropogenic ocean acidification over the twenty-first century and its impact on calcifying organisms. Nature 437, 681–686 (2005).
ADS  CAS  PubMed  Article  Google Scholar 

8.
Gazeau, F. et al. Impacts of ocean acidification on marine shelled molluscs. Mar. Biol. 160, 2207–2245 (2013).
CAS  Article  Google Scholar 

9.
Kroeker, K. J. et al. Impacts of ocean acidification on marine organisms: Quantifying sensitivities and interaction with warming. Glob. Change Biol. 19, 1884–1896 (2013).
ADS  Article  Google Scholar 

10.
Waldbusser, G. G. et al. Saturation-state sensitivity of marine bivalve larvae to ocean acidification. Nat. Clim. Change 5, 273–280 (2015).
ADS  CAS  Article  Google Scholar 

11.
Hoegh-Guldberg, O. et al. Coral reefs under rapid climate change and ocean acidification. Science 318, 1737–1742 (2007).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

12.
Moy, A. D., Howard, W. R., Bray, S. G. & Trull, T. W. Reduced calcification in modern Southern Ocean planktonic foraminifera. Nat. Geosci. 2, 276–280 (2009).
ADS  CAS  Article  Google Scholar 

13.
Bednaršek, N. et al. Extensive dissolution of live pteropods in the Southern Ocean. Nat. Geosci. 5, 881–885 (2012).
ADS  Article  CAS  Google Scholar 

14.
Bednaršek, N. et al. Limacina helicina shell dissolution as an indicator of declining habitat suitability owing to ocean acidification in the California Current Ecosystem. Proc. R. Soc. B Biol. Sci. 281, 20140123 (2014).
Article  CAS  Google Scholar 

15.
Manno, C. et al. Shelled pteropods in peril: Assessing vulnerability in a high CO2 ocean. Earth-Sci. Rev. 169, 132–145 (2017).
ADS  CAS  Article  Google Scholar 

16.
Lischka, S., Büdenbender, J., Boxhammer, T. & Riebesell, U. Impact of ocean acidification and elevated temperatures on early juveniles of the polar shelled pteropod Limacina helicina: Mortality, shell degradation, and shell growth. Biogeosciences 8, 919–932 (2011).
ADS  CAS  Article  Google Scholar 

17.
Bednaršek, N. et al. Exposure history determines pteropod vulnerability to ocean acidification along the US West Coast. Sci. Rep. 7, 1–12 (2017).
Article  CAS  Google Scholar 

18.
Comeau, S. et al. Impact of aragonite saturation state changes on migratory pteropods. Proc. R. Soc. B Biol. Sci. 279, 732–738 (2011).
Article  Google Scholar 

19.
Moya, A. et al. Near-future pH conditions severely impact calcification, metabolism and the nervous system in the pteropod Heliconoides inflatus. Glob. Change Biol. 22, 3888–3900 (2016).
ADS  Article  Google Scholar 

20.
Maas, A., Lawson, G. L., Bergan, A. J. & Tarrant, A. M. Exposure to CO2 influences metabolism, calcification and gene expression of the thecosome pteropod Limacina retroversa. J. Exp. Biol. 221, 164400 (2018).
Article  Google Scholar 

21.
Johnson, K. M. & Hofman, G. E. A transcriptome resource for the Antarctic pteropod Limacina helicina antarctica. Mar. Genom. 28, 25–28 (2016).
Article  Google Scholar 

22.
Feely, R. A. et al. Chemical and biological impacts of ocean acidification along the west coast of North America. Estuar. Coast. Shelf Sci. 183, 260–270 (2016).
ADS  CAS  Article  Google Scholar 

23.
Bednaršek, N. et al. El Niño-related thermal stress coupled with ocean acidification negatively impacts cellular to population-level responses in pteropods along the California Current System with implications for increased bioenergetic costs. Front. Mar. Sci. 5, 486 (2018).
Article  Google Scholar 

24.
Peck, V. L., Tarling, G. A., Manno, C., Harper, E. M. & Tynan, E. Outer organic layer and internal repair mechanism protects pteropod Limacina helicina from ocean acidification. Deep-Sea Res. II 127, 53–56 (2016).
Article  Google Scholar 

25.
Peck, V. L., Oakes, R. L., Harper, E. M., Manno, C. & Tarling, G. A. Pteropods counter mechanical damage and dissolution through extensive shell repair. Nat. Commun. 9, 264 (2018).
ADS  PubMed  PubMed Central  Article  CAS  Google Scholar 

26.
Howes, E. L., Eagle, R. A., Gattuso, J.-P. & Bijma, J. Comparison of Mediterranean pteropod shell biometrics and ultrastructure from historical (1910 and 1921) and present day (2012) samples provides baseline for monitoring effects of global change. PLoS ONE 1, 1–23 (2017).
Google Scholar 

27.
Oakes, R. L. & Sessa, J. A. Determining how biotic and abiotic variables affect the shell condition and parameters of Heliconoides inflatus pteropods from a sediment trap in the Cariaco Basin. Biogeosciences 7, 1975–1990 (2020).
ADS  Article  Google Scholar 

28.
Feely, R. A., Sabine, C. L., Hernandez-Ayon, J. M., Ianson, D. & Hales, B. Evidence for upwelling of corrosive “acidified” water onto the continental shelf. Science 320, 1490–1492 (2008).
ADS  CAS  PubMed  Article  Google Scholar 

29.
Alin, S. R., et al. Dissolved inorganic carbon, total alkalinity, pH on total scale, and other variables collected from profile and discrete sample observations using CTD, Niskin bottle, and other instruments from NOAA Ship Ronald H. Brown in the U.S. West Coast California Current System from 2016-05-08 to 2016-06-06 (NCEI Accession 0169412). Version 1.1. NOAA National Centers for Environmental Information dataset (2017). https://doi.org/10.7289/V5V40SHG.

30.
Northcott, D. et al. Impacts of urban carbon dioxide emissions on sea-air flux and ocean acidification in nearshore waters. PLoS ONE 14, e0214403 (2019).
CAS  PubMed  PubMed Central  Article  Google Scholar 

31.
Wang, K., Hunt, B. P. V., Liang, C., Pauly, D. & Pakhomov, E. A. Reassessment of the life cycle of the pteropod Limacina helicina from a high resolution interannual time series in the temperate North Pacific. ICES J. Mar. Sci. 74, 1906–1920 (2017).
Article  Google Scholar 

32.
Shimizu, K. et al. Phylogeography of the pelagic snail Limacina helicina (Gastropoda: Thecosomata) in the subarctic western North Pacific. J. Mollus. Stud. 84, 30–37 (2017).
Article  Google Scholar 

33.
Sromek, L., Lasota, R. & Wolowicz, M. Impact of glaciations on genetic diversity of pelagic mollusks: Antarctic Limacina Antarctica and Arctic Limacina helicina. Mar. Ecol. Prog. Ser. 525, 143–152 (2015).
ADS  Article  Google Scholar 

34.
Hunt, B. et al. Poles apart: the ‘bipolar’ pteropod species Limacina helicina is genetically distinct between the Arctic and Antarctic Oceans. PLoS ONE 5, e9835 (2010).
ADS  PubMed  PubMed Central  Article  CAS  Google Scholar 

35.
Bednaršek, N. et al. Systematic review and meta-analysis towards synthesis of thresholds of ocean acidification impacts on calcifying pteropods and interactions with warming. Front. Mar. Sci. 6, 227 (2019).
Article  Google Scholar 

36.
Vaquer-Sunyer, R. & Duarte, C. M. Thresholds of hypoxia for marine biodiversity. Proc. Natl. Acad. Sci. 105, 15452–15457 (2008).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

37.
Legaard, K. R. & Thomas, A. C. Spatial patterns in seasonal and interannual variability of chlorophyll and sea surface temperature in the California Current. J. Geophys. Res. 111, C06032 (2006).
ADS  Article  Google Scholar 

38.
Thomsen, J., Casties, I., Pansch, C., Körtzinger, A. & Melzner, F. Food availability outweighs ocean acidification effects in juvenile Mytilus edulis: Laboratory and field experiments. Glob. Change Biol. 19, 1017–1027 (2013).
ADS  Article  Google Scholar 

39.
Maas, A. E., Elder, L. E., Dierssen, H. M. & Seibel, B. A. Metabolic response of Antarctic pteropods (Mollusca: Gastropoda) to food deprivation and regional productivity. Mar. Ecol. Prog. Ser. 441, 129–139 (2011).
ADS  CAS  Article  Google Scholar 

40.
Ramajo, L. et al. Food supply confers calcifiers resistance to ocean acidification. Sci. Rep. 6, 1–6 (2016).
Article  CAS  Google Scholar 

41.
Thomas, A. C. & Strub, P. T. Interannual variability in phytoplankton pigment distribution during the spring transition along the west-coast of North America. J. Geophys. Res. 94, 18095–18117 (1989).
ADS  CAS  Article  Google Scholar 

42.
Bednaršek, N. & Ohman, M. D. Changes in pteropod vertical distribution, abundance and species richness in the California Current System due to ocean acidification. Mar. Ecol. Prog. Ser. 523, 93–103 (2015).
ADS  Article  CAS  Google Scholar 

43.
Lalli, C. M. & Gilmer, R. W. Pelagic Snails: The Biology of Holoplanktonic Gastropod Mollusks (Stanford University Press, Stanford, 1989).
Google Scholar 

44.
Seibel, B. A., Dymowska, A. & Rosenthal, J. Metabolic temperature compensation and coevolution of locomotory performance in pteropod molluscs. Integr. Comp. Biol. 47, 880–891 (2007).
PubMed  Article  PubMed Central  Google Scholar 

45.
Checa, A. G. Physical and biological determinants of the fabrication of molluscan shell microstructures. Front. Mar. Sci. 5, 535 (2018).
Article  Google Scholar 

46.
Marin, F., Le Roy, N. & Marie, B. The formation and mineralization of mollusk shell. Front. Biosci. 4, 1099–1125 (2012).
Article  Google Scholar 

47.
Kroeker, K. J., Kordas, R. L. & Harley, C. D. G. Embracing interactions in ocean acidification research: Confronting multiple stressor scenarios and context dependence. Biol. Lett. 13, 20160802 (2017).
PubMed  PubMed Central  Article  CAS  Google Scholar 

48.
Gruber, N. et al. Rapid progression of ocean acidification in the California Current System. Science 337, 220–223 (2012).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

49.
Buitenhuis, E. T., Le Quéré, C., Bednaršek, N. & Schiebel, R. Large contribution of pteropods to shallow CaCO3 export. Glob. Biogeochem. Cycles 33, 458–468 (2019).
ADS  CAS  Article  Google Scholar 

50.
Mackas, D. L. & Galbraith, M. D. Pteropod time-series from the NE Pacific. ICES J. Mar. Sci. 69, 448–459 (2012).
Article  Google Scholar 

51.
Lueker, T. J., Dickson, A. G. & Keeling, C. D. Ocean pCO2 calculated from dissolved inorganic carbon, alkalinity, and equations for K1 and K2: Validation based on laboratory measurements of CO2 in gas and seawater at equilibrium. Mar. Chem. 70, 105–119 (2000).
CAS  Article  Google Scholar 

52.
Kerney, M. P. & Cameron, R. A. D. A Field Guide to the Land Snails of Britain and North-West Europe (Collins, London, 1979).
Google Scholar 

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

54.
Oksanen, J., et al. Vegan: Community Ecology Package. R package version 2.5-2 (2018).

55.
Wall-Palmer, D. et al. Biogeography and genetic diversity of the atlantid heteropods. Progr. Oceanogr. 160, 1–25 (2018).
ADS  Article  Google Scholar 

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

57.
Barrett, J. C., Fry, B., Maller, J. & Daly, M. J. Haploview: Analysis and visualization of LD and haplotype maps. Bioinformatics 21, 263–265 (2005).
CAS  PubMed  Article  Google Scholar 

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

59.
Darriba, D., Taboada, G. L., Doallo, R. & Posada, D. jModelTest 2: More models, new heuristics and parallel computing. Nat. Methods 9, 772 (2012).
CAS  PubMed  PubMed Central  Article  Google Scholar  More