Philippa Kaur

1.
Chesson, P. Mechanisms of maintenance of species diversity. Annu. Rev. Ecol. Syst. 31, 343–366 (2000).
Article  Google Scholar 
2.
Carson, W. & Schnitzer, S. Tropical Forest Community Ecology (John Wiley & Sons, 2011).

3.
Givnish, T. J. On the causes of gradients in tropical tree diversity. J. Ecol. 87, 193–210 (1999).
Article  Google Scholar 

4.
Condit, R., Engelbrecht, B. M. J., Pino, D., Pérez, R. & Turner, B. L. Species distributions in response to individual soil nutrients and seasonal drought across a community of tropical trees. Proc. Natl Acad. Sci. USA 110, 5064–5068 (2013).
CAS  PubMed  Article  PubMed Central  Google Scholar 

5.
Bagchi, R. et al. Pathogens and insect herbivores drive rainforest plant diversity and composition. Nature 506, 85–88 (2014).
CAS  PubMed  Article  PubMed Central  Google Scholar 

6.
Zalamea, P.-C. et al. Seedling growth responses to phosphorus reflect adult distribution patterns of tropical trees. New Phytol. 212, 400–408 (2016).
CAS  PubMed  Article  PubMed Central  Google Scholar 

7.
Sarmiento, C. et al. Soilborne fungi have host affinity and host-specific effects on seed germination and survival in a lowland tropical forest. Proc. Natl Acad. Sci. USA 114, 11458–11463 (2017).
CAS  PubMed  Article  PubMed Central  Google Scholar 

8.
Ter Steege, H., Pitman, N. & Sabatier, D. A spatial model of tree alpha-diversity and tree density for the Amazon. Biodivers. Conserv. 12, 2255–2277 (2003).
Article  Google Scholar 

9.
Leigh, E. G. Jr. et al. Why do some tropical forests have so many species of trees? Biotropica 36, 447–473 (2004).
Google Scholar 

10.
Rahbek, C. & Graves, G. R. Detection of macro-ecological patterns in South American hummingbirds is affected by spatial scale. Proc. R. Soc. Lond. B. 267, 2259–2265 (2000).
CAS  Article  Google Scholar 

11.
Abrahamczyk, S., Kluge, J., Gareca, Y., Reichle, S. & Kessler, M. The influence of climatic seasonality on the diversity of different tropical pollinator groups. PLoS ONE 6, e27115 (2011).
CAS  PubMed  PubMed Central  Article  Google Scholar 

12.
Tedersoo, L. et al. Global diversity and geography of soil fungi. Science 346, 1256688 (2014).
PubMed  Article  CAS  PubMed Central  Google Scholar 

13.
Tonkin, J. D., Bogan, M. T., Bonada, N., Rios-Touma, B. & Lytle, D. A. Seasonality and predictability shape temporal species diversity. Ecology 98, 1201–1216 (2017).
PubMed  Article  PubMed Central  Google Scholar 

14.
Doležal, J., Lanta, V., Mudrák, O. & Lepš, J. Seasonality promotes grassland diversity: interactions with mowing, fertilization and removal of dominant species. J. Ecol. 107, 203–215 (2019).

15.
Arnold, A. E. et al. Fungal endophytes limit pathogen damage in a tropical tree. Proc. Natl Acad. Sci. USA 100, 15649–15654 (2003).
CAS  PubMed  Article  PubMed Central  Google Scholar 

16.
Arnold, A. E. & Engelbrecht, B. M. J. Fungal endophytes nearly double minimum leaf conductance in seedlings of a neotropical tree species. J. Trop. Ecol. 23, 369–372 (2007).
Article  Google Scholar 

17.
Costa Pinto, L. S., Azevedo, J. L., Pereira, J. O., Carneiro Vieira, M. L. & Labate, C. A. Symptomless infection of banana and maize by endophytic fungi impairs photosynthetic efficiency. New Phytol. 147, 609–615 (2000).
Article  Google Scholar 

18.
U’Ren, J. M. et al. Diversity and evolutionary origins of fungi associated with seeds of a neotropical pioneer tree: a case study for analysing fungal environmental samples. Mycol. Res. 113, 432–449 (2009).
PubMed  Article  CAS  PubMed Central  Google Scholar 

19.
Sanchez-Azofeifa, A., Oki, Y., Wilson Fernandes, G., Ball, R. A. & Gamon, J. Relationships between endophyte diversity and leaf optical properties. Trees 26, 291–299 (2012).
Article  Google Scholar 

20.
Vincent, J. B., Weiblen, G. D. & May, G. Host associations and beta diversity of fungal endophyte communities in New Guinea rainforest trees. Mol. Ecol. 25, 825–841 (2016).
CAS  PubMed  Article  Google Scholar 

21.
Suryanarayanan, T. S., Murali, T. S. & Venkatesan, G. Occurrence and distribution of fungal endophytes in tropical forests across a rainfall gradient. Can. J. Bot. 80, 818–826 (2002).
Article  Google Scholar 

22.
Zimmerman, N. B. & Vitousek, P. M. Fungal endophyte communities reflect environmental structuring across a Hawaiian landscape. Proc. Natl Acad. Sci. USA 109, 13022–13027 (2012).
CAS  PubMed  Article  Google Scholar 

23.
Higgins, K. L., Arnold, A. E., Coley, P. D. & Kursar, T. A. Communities of fungal endophytes in tropical forest grasses: highly diverse host- and habitat generalists characterized by strong spatial structure. Fungal Ecol. 8, 1–11 (2014).
Article  Google Scholar 

24.
Darcy, J. L. et al. Fungal communities living within leaves of native Hawaiian dicots are structured by landscape-scale variables as well as by host plants. Mol. Ecol. 29, 3102–3115 (2020).
Article  Google Scholar 

25.
Arnold, A. E. & Lutzoni, F. Diversity and host range of foliar fungal endophytes: are tropical leaves biodiversity hotspots? Ecology 88, 541–549 (2007).
PubMed  Article  Google Scholar 

26.
Tellez, P. H. Tropical plants and fungal symbionts: Leaf functional traits as drivers of plant-fungal interactions. PhD dissertation (Tulane University, 2019).

27.
U’Ren, J. M. et al. Host availability drives distributions of fungal endophytes in the imperilled boreal realm. Nat. Ecol. Evol. 3, 1430–1437 (2019).
PubMed  Article  Google Scholar 

28.
Arnold, A. E. & Herre, E. A. Canopy cover and leaf age affect colonization by tropical fungal endophytes: ecological pattern and process in Theobroma cacao (Malvaceae). Mycologia 95, 388–398 (2003).
PubMed  Article  PubMed Central  Google Scholar 

29.
Rodríguez-Quiel, E. E., Mendieta-Leiva, G. & Bader, M. Y. Elevational patterns of bryophyte and lichen biomass differ among substrates in the tropical montane forest of Baru Volcano, Panama. J. Bryol. 41, 95–106 (2019).
Article  Google Scholar 

30.
Magill, B., Solomon, J. & Stimmel, H. Tropicos Specimen Data. http://www.tropicos.org (2019).

31.
Arnold, A. E. et al. A phylogenetic estimation of trophic transition networks for ascomycetous fungi: are lichens cradles of symbiotrophic fungal diversification? Syst. Biol. 58, 283–297 (2009).
PubMed  Article  PubMed Central  Google Scholar 

32.
U’Ren, J. M., Lutzoni, F., Miadlikowska, J., Laetsch, A. D. & Arnold, A. E. Host and geographic structure of endophytic and endolichenic fungi at a continental scale. Am. J. Bot. 99, 898–914 (2012).
PubMed  Article  PubMed Central  Google Scholar 

33.
Arnold, A. E., Maynard, Z., Gilbert, G. S., Coley, P. D. & Kursar, T. A. Are tropical fungal endophytes hyperdiverse? Ecol. Lett. 3, 267–274 (2000).
Article  Google Scholar 

34.
Phillips, O. L., Hall, P., Gentry, A. H., Sawyer, S. A. & Vásquez, R. Dynamics and species richness of tropical rain forests. Proc. Natl Acad. Sci. USA 91, 2805–2809 (1994).
CAS  PubMed  Article  PubMed Central  Google Scholar 

35.
Usinowicz, J. et al. Temporal coexistence mechanisms contribute to the latitudinal gradient in forest diversity. Nature 550, 105–108 (2017).
CAS  PubMed  Article  PubMed Central  Google Scholar 

36.
Pianka, E. R. Latitudinal gradients in species diversity: a review of concepts. Am. Nat. 100, 33–46 (1966).
Article  Google Scholar 

37.
Coley, P. D. & Barone, J. A. Herbivory and plant defenses in tropical forests. Annu. Rev. Ecol. Syst. 27, 305–335 (1996).
Article  Google Scholar 

38.
Thrall, P. H., Hochberg, M. E., Burdon, J. J. & Bever, J. D. Coevolution of symbiotic mutualists and parasites in a community context. Trends Ecol. Evol. 22, 120–126 (2007).
PubMed  Article  PubMed Central  Google Scholar 

39.
Poisot, T., Bever, J. D., Nemri, A., Thrall, P. H. & Hochberg, M. E. A conceptual framework for the evolution of ecological specialisation. Ecol. Lett. 14, 841–851 (2011).
PubMed  PubMed Central  Article  Google Scholar 

40.
Van Bael, S., Estrada, C. & Arnold, A. E. Chapter 6: foliar endophyte communities and leaf traits in tropical trees. In The Fungal Community: Its Organization and Role in the Ecosystem. (eds Dighton, J. & White, J. F.) 79–94 (CRC Press, 2017).

41.
Oono, R. et al. Species diversity of fungal endophytes across a stress gradient for plants. New Phytol. 228, 210–225 (2020).

42.
Top, S. M., Preston, C. M., Dukes, J. S. & Tharayil, N. Climate influences the content and chemical composition of foliar tannins in green and senesced tissues of Quercus rubra. Front. Plant Sci. 8, 423 (2017).
PubMed  PubMed Central  Article  Google Scholar 

43.
Higginbotham, S. J. et al. Bioactivity of fungal endophytes as a function of endophyte taxonomy and the taxonomy and distribution of their host plants. PLoS ONE 8, e73192 (2013).
CAS  PubMed  PubMed Central  Article  Google Scholar 

44.
Hijmans, R. J., Cameron, S. E., Parra, J. L., Jones, P. G. & Jarvis, A. Very high resolution interpolated climate surfaces for global land areas. Int. J. Climatol. 25, 1965–1978 (2005).
Article  Google Scholar 

45.
Olson, D. M. et al. Terrestrial ecoregions of the world: a new map of life on earth: a new global map of terrestrial ecoregions provides an innovative tool for conserving biodiversity. Bioscience 51, 933–938 (2001).
Article  Google Scholar 

46.
Prada, C. M. et al. Soils and rainfall drive landscape-scale changes in the diversity and functional composition of tree communities in premontane tropical forest. J. Veg. Sci. 28, 859–870 (2017).
Article  Google Scholar 

47.
Walker, K. Capturing ephemeral forest dynamics with hybrid time-series and composite mapping in the Republic of Panama. Int. J. Appl. Earth Obs. Geoinf. 87, 102029 (2020).
Article  Google Scholar 

48.
Leung, B., Hudgins, E. J., Potapova, A. & Ruiz-Jaen, M. C. A new baseline for countrywide α-diversity and species distributions: illustration using > 6,000 plant species in Panama. Ecol. Appl. 29, e01866 (2019).
PubMed  Article  Google Scholar 

49.
Pyke, C. R., Condit, R., Aguilar, S. & Lao, S. Floristic composition across a climatic gradient in a neotropical lowland forest. J. Veg. Sci. 12, 553–566 (2001).
Article  Google Scholar 

50.
Lieberman, D., Lieberman, M., Peralta, R. & Hartshorn, G. S. Tropical forest structure and composition on a large-scale altitudinal gradient in Costa Rica. J. Ecol. 84, 137–152 (1996).
Article  Google Scholar 

51.
Bowman, E. A. & Arnold, A. E. Distributions of ectomycorrhizal and foliar endophytic fungal communities associated with Pinus ponderosa along a spatially constrained elevation gradient. Am. J. Bot. 105, 687–699 (2018).
PubMed  Article  PubMed Central  Google Scholar 

52.
U’Ren, J. M. et al. Tissue storage and primer selection influence pyrosequencing-based inferences of diversity and community composition of endolichenic and endophytic fungi. Mol. Ecol. Resour. 14, 1032–1048 (2014).
PubMed  PubMed Central  Google Scholar 

53.
U’Ren, J. M. & Arnold, A. E. 96 well DNA extraction protocol for plant and lichen tissue stored in CTAB. protocols.io. https://doi.org/10.17504/protocols.io.fscbnaw (2017).

54.
Daru, B. H., Bowman, E. A., Pfister, D. H. & Arnold, A. E. A novel proof of concept for capturing the diversity of endophytic fungi preserved in herbarium specimens. Philos. Trans. R. Soc. Lond. B Biol. Sci. 374, 20170395 (2018).

55.
Callahan, B. J., McMurdie, P. J. & Holmes, S. P. Exact sequence variants should replace operational taxonomic units in marker-gene data analysis. ISME J. 11, 2639–2643 (2017).
PubMed  PubMed Central  Article  Google Scholar 

56.
Oksanen, J. et al. Vegan: community ecology package, version 2.5-2. https://CRAN.R-project.org/package=vegan (2018).

57.
Schoener, T. W. Food webs from the small to the large: the Robert H. MacArthur award lecture. Ecology 70, 1559–1589 (1989).
Article  Google Scholar 

58.
Apigo, A. & Oono, R. Dimensions of host specificity in foliar fungal endophytes. In Endophytes of Forest Trees: Biology and Applications (eds Pirttilä, A. M. & Frank, A. C.) 15–42 (Springer International Publishing, 2018).

59.
Oita, S. et al. Data from: climate and seasonality drive the richness and composition of tropical fungal endophytes at a landscape scale. figshare https://doi.org/10.6084/m9.figshare.c.5084366.v1 (2020). More

Experiment 1: proof of principle
As a proof of principle, we designed the first experiment to test whether appendages may help to conceal the outline. We created an image of a uniformly light grey coloured circular object with a size of 2950 pixels (px)/250.0 mm circumference and 470 px/39.8 mm radius on a dark grey background using Adobe InDesign CS6 version 8.030. The initial setup started with no appendages added to the outline (Fig. 3a, ‘0’). We then added object-coloured appendages (i.e. lines of 1 Pt/4 px/0.4 mm thickness and 118 px/10.0 mm length) with regular intervals resembling protruding neoptile chick feathers orthogonally to the object outline (‘Basic Scenario’, Fig. 3a). The first image with appendages had 32 appendages added to the outline (Fig. 3a, ‘32′). We then doubled the number of appendages stepwise creating denser spaced appendages to the outline until the extended outline was completely filled (Fig. 3a, ‘full circle’). For the vision of a simulated predator, we used the spatial acuity from humans (Homo sapiens, 72 cycles per degree, cpd)36,45,46 in the basic scenario. The full details for the parameters are provided in Supplementary Table S1 (a–g).
Figure 3

(a) Basic Scenario: Seven stages of the artificial chick setup with varying number of thin, non-transparent appendages having all the same length. (b) Scenario 1: varying appendage thickness applied to the Basic Scenario. (c) Scenario 2: varying appendage transparency applied to the Basic Scenario. (d) Scenario 3: varying appendage length heterogeneity applied to the Basic Scenario. (e) Scenario 4: varying background complexity with chessboard backgrounds. (f) Scenario 5: high, medium and low spatial acuity applied to the Basic Scenario. (a–f) The analysed region of interest (ROI) is highlighted in red for clarification only. The figure was produced in Adobe Photoshop29 and InDesign30.

Full size image

To further explore the mechanism, we altered appendage characteristics, background and the spatial acuity of the predator. First, we increased appendage thickness to 2 Pt/8 pixels/ 0.7 mm (Scenario 1a) and 3 Pt/12 pixels/1.1 mm (Scenario 1b) resulting in decreased inter-appendage intervals (Fig. 3b and Supplementary Table S1, h–u). Second, we changed appendage transparency to 25% (Scenario 2a) and 50% transparency (Scenario 2b) (Fig. 3c and Supplementary Table S1, v to ai). Third, we varied the appendage length heterogeneity; half of the appendages having 50% of the length (Scenario 3a), and half of the appendages at 25% and one quarter at 50% of the original appendage length (Scenario 3b) (Fig. 3d and Supplementary Table S1, aj to aw). Fourth, we investigated the effect of background complexity on the detectability of the outline. As background, we used a chessboard pattern with large squares (346 pixels/29.3 mm, Scenario 4a) and with small squares (86 pixels/7.3 mm, Scenario 4b) (Fig. 3e and Supplementary Table S1, ax to bk). Fifth, we altered the spatial acuity to test whether or how the visual systems of different predators would affect detectability. We simulated the spatial acuity of a corvid predator (30 cpd, Scenario 5a) and canid predator (10 cpd, Scenario 5b) (Fig. 3f and Supplementary Table S1, bl to by), the two most common predators of ground-nesting plovers16,47,48. This range also covered other potential predators (Supplementary Table S2).
We did not account for differences in colour vision between different predators as the setup mostly consists of greyscale images that predominantly differ in luminance. Note that in many animals, visual acuity is greater for achromatic than chromatic stimuli34,49.
We conducted visual modelling and visual analysis using the Quantitative Colour Pattern Analysis (QCPA) framework27 integrated into the Multispectral Image Analysis and Calibration (MICA) toolbox50 for ImageJ version 1.52a51. We converted the generated images into multispectral images containing the red, green and blue channel in a stack and transformed them further into 32-bits/channel cone-catch images based on the human visual system, which are required by the framework. To create the luminance channel, we averaged the long and medium wave channel, which is thought to be representative of human vision52. We modelled the spatial acuity with Gaussian Acuity Control at a viewing distance of 1300 mm and a minimum resolvable angle (MRA) of 0.01389. To increase biological accuracy, we applied a Receptor Noise Limited (RNL) filter that reduces noise and reconstructs edges in the image. The RNL filter used the Weber fractions “Human 0.05” provided by the framework (longwave 0.05, mediumwave 0.07071, shortwave 0.1657), luminance 0.1, 5 iterations, a radius of 5 pixels and a falloff of 3 pixels as specified in van den Berg et al.27 (Supplementary Fig. S1).
To test for the detectability of the outline, we used LEIA27, which is conceptually similar to the boundary strength analysis34. Boundary strength analysis requires an image with clearly delineated (clustered) colour and luminance pattern elements. However, a large degree of subthreshold details, which may be still perceived by the viewer gets lost in the clustering process. LEIA has the advantage of not requiring such a clustered input and therefore can be directly applied to RNL filtered images. LEIA measures the edge intensity (i.e. the luminance contrast) locally at each position in the image. The output image displays ΔS values in a 32-bit stack of four slices, where each slice shows the values measured in different angles (horizontal, vertical and the two diagonals, for more details, see van den Berg et al.27).
We ran LEIA on the chosen region of interest (ROI) with the same Weber fractions used for the RNL filter. The ROI was a 180 pixel-wide band that included the area of the appendages extended by 30 pixels towards the object inside and towards the outside (Fig. 3a). We log-transformed the ΔS values as recommended for natural scenes53 to make the results comparable to the natural background images used in Experiment 2 (see below). To test whether the size of the ROI affected our results, we ran an additional analysis using a 1500 × 1500 pixel-wide rectangle surrounding the object as the ROI, which included a bigger area of the background and the full object inside (Supplementary Fig. S2).
We extracted the luminance ΔS values from the four slices of the output image stack in ImageJ and stored them in separate matrices for further analysis using R version 3.5.328. ImageJ generally assigned values outside the chosen ROI to zero. Thus, we first discarded all values of zero. We then set all negative values that arose as artefacts in areas without any edges to zero, in order to make them biologically meaningful. We then identified the parallel maximum (R function pmax ()) of the four interrelated direction matrices and transferred this value to a new matrix.
High luminance and colour contrasts imply high conspicuousness34. Consequently, a lower luminance contrast leads to lower conspicuousness and therefore, better camouflage. As the outline is an important cue for predators locating and identifying a prey item7, we assumed that especially low contrasts in the outline of an object improve camouflage. Thus, a reduction of edge intensity in the object outline by the appendages indicates a camouflage improvement. To test whether the object outline became less detectable we compared the edge intensity of the outline pixels in the basic scenario without appendages (Supplementary Table S1, a) with corresponding pixels from other scenarios. The outline pixels were characterised by high edge intensity and constituted a prominent peak. They comprised 1.59% of all pixels in the analysis focused on the contour region (see “Results”, Fig. 1a). For all scenarios, we calculated the mean edge intensity of the high edge intensity pixels (HEI pixels) and identified the changes with parameter variation. Unless otherwise stated we used R28 to produce graphs and panels.
As an alternative mechanism, we tested whether appendages create a transition zone with intermediate luminance around the object (Mean Luminance Comparison (MLC), Supplementary material). We calculated the mean luminance of the object inside up to the border (object region), the area covered by appendages (appendage region) and the background (background region). We predicted that the appendage region would be characterised by intermediate luminance between object and background and therefore provide a luminance transition zone to conceal the object outline.
Experiment 2: chick photographs
Using pictures of young snowy plover chicks hiding when approached by a simulated predator, we tested if protruding neoptile feathers helped to conceal the chicks’ outline and therefore improve their camouflage.
We studied snowy plovers in their natural environment at Bahía de Ceuta, Sinaloa, Mexico. Fieldwork permits were granted by the Secretaría de Medio Ambiente y Recursos Naturales (SEMARNAT). All field activities were performed in accordance with the approved ethical guidelines outlined by SEMARNAT. The breeding site consists of salt flats that are sparsely vegetated and surrounded by mangroves54. The predators of chicks are not well described but likely similar to the egg predator community that includes several mammalian predators such as racoon, opossum, coyote, bob cat, avian predators such as crested caracara Caracara cheriway and reptiles17. General field methodology is provided elsewhere55,56. In 2017, we took photographs of young (one to 3 days old) chicks hiding on the ground, that had already left the nest scrape. To photograph the chicks, two observers approached free-roaming families with two mobile hides within the period one hour after sunrise and one hour before sunset. At a distance of 100–200 m, one observer acted as ‘predator’, left the hide and openly approached the brood while the second observer kept watching the chicks. The chicks responded by crouching to the ground and staying motionless while the parents were alarming. The second observer directed the ‘predator’ to the approximate hiding place. When searching for the chicks, we took great care to reduce the number of steps to avoid modification of the ground through our tracks.
Once the first chick had been found, the second observer joined the ‘predator’ and took the chick photographs. We used a Nikon D7000 camera converted to full spectrum including the UV range (Optic Makario GmbH, Germany) and a Nikkor macro 105 mm lens that allows transmission of light at low wavebands. The equipment was chosen because calibration data were available for this combination50. Each hiding background was photographed with and without the chick using a UV pass filter for the UV spectrum and a UV/IR blocking filter (“IR-Neutralisationsfilter NG”, Optic Makario GmbH, Germany) for the visible spectrum. The camera was set to an aperture of f/8, ISO 400 and the pictures were stored in “RAW” file format. We used exposure bracketing to produce three images to ensure that at least one picture was not over or underexposed. A 25% reflectance standard (Zenith Polymer Diffuse Reflectance Standard provided by SPHEREOPTICS, Germany) placed in the corner of each picture enabled a subsequent standardizing of light conditions.
In total, we took pictures of 32 chicks from 15 families. For 21 chicks we obtained photographs suitable for further analyses with an unobstructed view to the entire chick and only one chick per photograph. Of these, we randomly selected pictures of 15 chicks. Unfortunately, it was not possible to obtain proper alignment of visual and UV pictures in ImageJ as either chick or camera moved slightly in the break between changing filters for the two settings. Therefore, we restricted our analyses to human colour vision and discarded the UV pictures for further analysis.
In each picture, we manually selected the chick outline and the feather-boundary as a basis for the ROIs (Fig. 2a–c). The chick outline included bill, legs, rings and all areas densely covered by feathers without background shining through. We then marked the feather-boundary, i.e., the smoothened line created by the protruding neoptile feather tips. In the next step, we transferred images of chicks with or without protruding feathers, i.e. cropped at feather-boundary or chick outline, respectively, and inserted them into a uniform or the natural background. First, we cropped the chick without protruding feathers and transferred it into a uniform black background. Second, we cropped the chick including all feathers and inserted it into exactly the same hiding spot on the picture of the natural background (Fig. 2b). Third, we cropped the chick excluding the protruding feathers and transferred it into the natural background (Fig. 2c).
We then proceeded with LEIA following the protocol of experiment 1 with the following changes. Again, the selected ROI was the contour region ranging from the chick outline extended by 30 pixels towards the chick inside to the feather-boundary extended by 30 pixels towards the outside. We excluded all areas of the ROI that showed a shadow of the chick as the chicks’ shadow was missing on the empty natural background images to which the cropped chicks were transferred to (Fig. 2a–c). We used the images of the cropped chicks on the black background to determine the threshold of the HEI pixels according to the protocol of experiment 1 for each chick separately. For each cropped chick that was transferred to the picture with the natural background, we compared the mean edge intensity of the HEI pixels provided by LEIA with and without protruding feathers (Fig. 2b,c) using a two-sided paired t-test.
We also calculated mean luminance differences for chick photographs. Details for this MLC are given in the supplementary material. More

1.
Le Quéré, C. et al. Global Carbon Budget 2018. Earth Syst. Sci. Data 10, 2141–2194 (2018).
ADS  Article  Google Scholar 
2.
Cook-Patton, S. C. et al. Mapping carbon accumulation potential from global natural forest regrowth. Nature 585, 545–550 (2020).
ADS  CAS  PubMed  Article  Google Scholar 

3.
Houghton, R. A., Byers, B. & Nassikas, A. A. A role for tropical forests in stabilizing atmospheric CO2. Nat. Clim. Chang. 5, 1022–1023 (2015).
ADS  Article  Google Scholar 

4.
Chazdon, R. L. et al. Carbon sequestration potential of second-growth forest regeneration in the Latin American tropics. Sci. Adv. 2, e1501639 (2016).
ADS  PubMed  PubMed Central  Article  CAS  Google Scholar 

5.
Philipson, C. D. et al. Active restoration accelerates the carbon recovery of human-modified tropical forests. Science 369, 838–841 (2020).
ADS  CAS  PubMed  Article  Google Scholar 

6.
Taubert, F. et al. Global patterns of tropical forest fragmentation. Nature 554, 519–522 (2018).
ADS  CAS  PubMed  Article  Google Scholar 

7.
Tabarelli, M., Lopes, A. V. & Peres, C. A. Edge-effects drive tropical forest fragments towards an early-successional system. Biotropica 40, 657–661 (2008).
Article  Google Scholar 

8.
Arroyo-Rodríguez, V. et al. Multiple successional pathways in human-modified tropical landscapes: new insights from forest succession, forest fragmentation and landscape ecology research. Biol. Rev. Camb. Philos. Soc. 92, 326–340 (2017).
PubMed  Article  Google Scholar 

9.
Collins, C. D. et al. Fragmentation affects plant community composition over time. Ecography 40, 119–130 (2017).
Article  Google Scholar 

10.
Hansen, M. C. et al. The fate of tropical forest fragments. Sci. Adv. 6 eaax8574 (2020).

11.
Qie, L. et al. Long-term carbon sink in Borneo’s forests halted by drought and vulnerable to edge effects. Nat. Commun. 8, 1966 (2017).
ADS  PubMed  PubMed Central  Article  CAS  Google Scholar 

12.
Thirumalai, K., DiNezio, P. N., Okumura, Y. & Deser, C. Extreme temperatures in Southeast Asia caused by El Niño and worsened by global warming. Nat. Commun. 8, 15531 (2017).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

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

14.
Vogel, M. M., Hauser, M. & Seneviratne, S. I. Projected changes in hot, dry and wet extreme events’ clusters in CMIP6 multi-model ensemble. Environ. Res. Lett. 15, 094021 (2020).
ADS  Article  Google Scholar 

15.
McDowell, N. et al. Drivers and mechanisms of tree mortality in moist tropical forests. N. Phytol. https://doi.org/10.1111/[email protected]/(ISSN)1469-8137.DroughtImpactsonTropicalForests (2018).

16.
Walker, A. P. et al. Decadal biomass increment in early secondary succession woody ecosystems is increased by CO2 enrichment. Nat. Commun. 10, 454 (2019).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

17.
Riutta, T. et al. Logging disturbance shifts net primary productivity and its allocation in Bornean tropical forests. Glob. Chang. Biol. 24, 2913–2928 (2018).
ADS  PubMed  Article  Google Scholar 

18.
Both, S. et al. Logging and soil nutrients independently explain plant trait expression in tropical forests. N. Phytol. 221, 1853–1865 (2019).
CAS  Article  Google Scholar 

19.
Swinfield, T. et al. Imaging spectroscopy reveals the effects of topography and logging on the leaf chemistry of tropical forest canopy trees. Glob. Chang. Biol. 26, 989–1002 (2020).
ADS  PubMed  Article  Google Scholar 

20.
Jotan, P., Maycock, C. R., Burslem, D., Berhaman, A. & Both, S. Comparative vessel traits of macaranga gigantea and vatica dulitensis from Malaysian Borneo. J. Trop. Sci. 32, 25–34 (2020).
Google Scholar 

21.
Uriarte, M. et al. Impacts of climate variability on tree demography in second growth tropical forests: the importance of regional context for predicting successional trajectories. Biotropica 48, 780–797 (2016).
Article  Google Scholar 

22.
Park Williams, A. et al. Temperature as a potent driver of regional forest drought stress and tree mortality. Nat. Clim. Chang. 3, 292–297 (2013).
ADS  Article  Google Scholar 

23.
Wolfe, B. T., Sperry, J. S. & Kursar, T. A. Does leaf shedding protect stems from cavitation during seasonal droughts? A test of the hydraulic fuse hypothesis. N. Phytol. 212, 1007–1018 (2016).
Article  Google Scholar 

24.
Slik, J. W. F. El Niño droughts and their effects on tree species composition and diversity in tropical rain forests. Oecologia 141, 114–120 (2004).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

25.
Jucker, T. et al. Topography shapes the structure, composition and function of tropical forest landscapes. Ecol. Lett. 21, 989–1000 (2018).
PubMed  PubMed Central  Article  Google Scholar 

26.
Itoh, A. et al. Importance of topography and soil texture in the spatial distribution of two sympatric dipterocarp trees in a Bornean rainforest. Ecol. Res. 18, 307–320 (2003).
Article  Google Scholar 

27.
Stovall, A. E. L., Shugart, H. H. & Yang, X. Reply to ‘Height-related changes in forest composition explain increasing tree mortality with height during an extreme drought’. Nat. Commun. 11, 3401 (2020).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

28.
Maréchaux, I. et al. Drought tolerance as predicted by leaf water potential at turgor loss point varies strongly across species within an Amazonian forest. Funct. Ecol. 29, 1268–1277 (2015).
Article  Google Scholar 

29.
Cosme, L. H. M., Schietti, J., Costa, F. R. C. & Oliveira, R. S. The importance of hydraulic architecture to the distribution patterns of trees in a central Amazonian forest. N. Phytol. 215, 113–125 (2017).
Article  Google Scholar 

30.
Bittencourt, P. R. L. et al. Amazonia trees have limited capacity to acclimate plant hydraulic properties in response to long-term drought. Glob. Chang. Biol. https://doi.org/10.1111/gcb.15040 (2020).

31.
Luke, S. H. et al. Riparian buffers in tropical agriculture: Scientific support, effectiveness and directions for policy. J. Appl. Ecol. 56, 85–92 (2019).
Article  Google Scholar 

32.
Padfield, R. et al. Co-Producing a Research Agenda for Sustainable Palm Oil. https://doi.org/10.3389/ffgc.2019.00013 (2019).

33.
Zhao, K. et al. Utility of multitemporal lidar for forest and carbon monitoring: Tree growth, biomass dynamics, and carbon flux. Remote Sens. Environ. 204, 883–897 (2018).
ADS  Article  Google Scholar 

34.
Simonson, W., Ruiz-Benito, P., Valladares, F. & Coomes, D. Modelling above-ground carbon dynamics using multi-temporal airborne lidar: insights from a Mediterranean woodland. Biogeosciences 13, 961–973 (2016).
ADS  CAS  Article  Google Scholar 

35.
Leitold, V. et al. El Niño drought increased canopy turnover in Amazon forests. N. Phytol. 219, 959–971 (2018).
CAS  Article  Google Scholar 

36.
Moura, Y. Mde et al. Carbon dynamics in a human-modified tropical forest: a case study using multi-temporal LiDAR data. Remote Sens. 12, 430 (2020).
ADS  Article  Google Scholar 

37.
Simonson, W., Allen, H. & Coomes, D. Effect of tree phenology on LiDAR measurement of mediterranean forest structure. Remote Sens. 10, 659 (2018).
ADS  Article  Google Scholar 

38.
Sullivan, M. J. P. et al. Long-term thermal sensitivity of Earth’s tropical forests. Science 368, 869–874 (2020).
ADS  CAS  PubMed  Article  Google Scholar 

39.
Coomes, D. A. et al. Area-based vs tree-centric approaches to mapping forest carbon in Southeast Asian forests from airborne laser scanning data. Remote Sens. Environ. 194, 77–88 (2017).
ADS  Article  Google Scholar 

40.
Asner, G. P. et al. Mapped aboveground carbon stocks to advance forest conservation and recovery in Malaysian Borneo. Biol. Conserv. 217, 289–310 (2018).
Article  Google Scholar 

41.
Burton, C., Rifai, S. & Malhi, Y. Inter-comparison and assessment of gridded climate products over tropical forests during the 2015/2016 El Niño. Philosophical Transactions of the Royal Society B: Biological Sciences. 373, 1760 (2018).
Article  Google Scholar 

42.
Ordway, E. M. & Asner, G. P. Carbon declines along tropical forest edges correspond to heterogeneous effects on canopy structure and function. Proc. Natl Acad. Sci. USA 117, 7863–7870 (2020).
CAS  PubMed  Article  Google Scholar 

43.
Phillips, O. L. et al. Drought sensitivity of the Amazon rainforest. Science 323, 1344–1347 (2009).
ADS  CAS  PubMed  Article  Google Scholar 

44.
Sevanto, S. Phloem transport and drought. J. Exp. Bot. 65, 1751–1759 (2014).
CAS  PubMed  Article  Google Scholar 

45.
Aleixo, I. et al. Amazonian rainforest tree mortality driven by climate and functional traits. Nat. Clim. Chang. 9, 384–388 (2019).
ADS  Article  Google Scholar 

46.
Woodgate, W. et al. Quantifying the impact of woody material on leaf area index estimation from hemispherical photography using 3D canopy simulations. Agric. Meteorol. 226-227, 1–12 (2016).
Article  Google Scholar 

47.
Nunes, M. H. et al. Changes in leaf functional traits of rainforest canopy trees associated with an El Niño event in Borneo. Environ. Res. Lett. 14, 085005 (2019).
ADS  CAS  Article  Google Scholar 

48.
Tang, H. & Dubayah, R. Light-driven growth in Amazon evergreen forests explained by seasonal variations of vertical canopy structure. Proc. Natl Acad. Sci. USA 114, 2640–2644 (2017).
CAS  PubMed  Article  Google Scholar 

49.
Doughty, C. E. et al. Drought impact on forest carbon dynamics and fluxes in Amazonia. Nature 519, 78–82 (2015).
ADS  CAS  PubMed  Article  Google Scholar 

50.
Marvin, D. C. & Asner, G. P. Branchfall dominates annual carbon flux across lowland Amazonian forests. Environ. Res. Lett. 11, 094027 (2016).
ADS  Article  CAS  Google Scholar 

51.
Roussel, J.-R., Caspersen, J., Béland, M., Thomas, S. & Achim, A. Removing bias from LiDAR-based estimates of canopy height: accounting for the effects of pulse density and footprint size. Remote Sens. Environ. 198, 1–16 (2017).
ADS  Article  Google Scholar 

52.
Sist, P. & Nguyen-Thé, N. Logging damage and the subsequent dynamics of a dipterocarp forest in East Kalimantan (1990–1996). Ecol. Manag. 165, 85–103 (2002).
Article  Google Scholar 

53.
Rutishauser, E. et al. Rapid tree carbon stock recovery in managed Amazonian forests. Curr. Biol. 25, R787–R788 (2015).
CAS  PubMed  Article  PubMed Central  Google Scholar 

54.
Giambelluca, T. W., Ziegler, A. D., Nullet, M. A., Truong, D. M. & Tran, L. T. Transpiration in a small tropical forest patch. Agric. Meteorol. 117, 1–22 (2003).
Article  Google Scholar 

55.
Ewers, R. M. & Banks-Leite, C. Fragmentation impairs the microclimate buffering effect of tropical forests. PLoS ONE 8, e58093 (2013).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

56.
Laurance, W. F. & Curran, T. J. Impacts of wind disturbance on fragmented tropical forests: a review and synthesis. Austral. Ecol. 33, 399–408 (2008).
Article  Google Scholar 

57.
Jucker, T. et al. Canopy structure and topography jointly constrain the microclimate of human-modified tropical landscapes. Glob. Chang. Biol. 24, 5243–5258 (2018).
ADS  PubMed  Article  PubMed Central  Google Scholar 

58.
Laurance, W. F. Theory meets reality: how habitat fragmentation research has transcended island biogeographic theory. Biol. Conserv. 141, 1731–1744 (2008).
Article  Google Scholar 

59.
Muscarella, R., Kolyaie, S., Morton, D. C., Zimmerman, J. K. & Uriarte, M. Effects of topography on tropical forest structure depend on climate context. J. Ecol. 108, 145–159 (2020).
Article  Google Scholar 

60.
Werner, F. A. & Homeier, J. Is tropical montane forest heterogeneity promoted by a resource‐driven feedback cycle? Evidence from nutrient relations, herbivory and litter decomposition along a topographical gradient. Funct. Ecol. 29, 430–440 (2015).
Article  Google Scholar 

61.
Gonzalez‐Akre, E. et al. Patterns of tree mortality in a temperate deciduous forest derived from a large forest dynamics plot. Ecosphere 7, G04014 (2016).
Article  Google Scholar 

62.
Hansen, M. C. et al. High-resolution global maps of 21st-century forest cover change. Science 342, 850–853 (2013).
ADS  CAS  Article  Google Scholar 

63.
Williamson, J. et al. Riparian Buffers Act as Microclimatic Refugia in Oil Palm Landscapes. https://doi.org/10.17863/CAM.57796 (2020).

64.
Hurst, M. D., Mudd, S. M., Walcott, R., Attal, M. & Yoo, K. Using hilltop curvature to derive the spatial distribution of erosion rates: hilltop curvature predicts erosion rates. J. Geophys. Res. 117, F2 (2012).
Google Scholar 

65.
Gaveau, D. L. A. et al. Rapid conversions and avoided deforestation: examining four decades of industrial plantation expansion in Borneo. Sci. Rep. 6, 32017 (2016).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

66.
Ewers, R. M. et al. A large-scale forest fragmentation experiment: the Stability of Altered Forest Ecosystems Project. Philos. Trans. R. Soc. Lond. B Biol. Sci. 366, 3292–3302 (2011).
PubMed  PubMed Central  Article  Google Scholar 

67.
Pfeifer, M. et al. Mapping the structure of Borneo’s tropical forests across a degradation gradient. Remote Sens. Environ. 176, 84–97 (2016).
ADS  Article  Google Scholar 

68.
Reynolds, G., Payne, J., Sinun, W., Mosigil, G. & Walsh, R. P. D. Changes in forest land use and management in Sabah, Malaysian Borneo, 1990-2010, with a focus on the Danum Valley region. Philos. Trans. R. Soc. Lond. B Biol. Sci. 366, 3168–3176 (2011).
PubMed  PubMed Central  Article  Google Scholar 

69.
Clarke, A. Principles of thermal ecology: temperature, energy, and life. Oxford Scholarship Online https://doi.org/10.1093/oso/9780199551668.001.0001 (2017).

70.
Bolton, D. The computation of equivalent potential temperature. Mon. Weather Rev. 108, 1046–1053 (1980).
ADS  Article  Google Scholar 

71.
Asner, G. P. et al. Carnegie airborne observatory-2: increasing science data dimensionality via high-fidelity multi-sensor fusion. Remote Sens. Environ. 124, 454–465 (2012).
ADS  Article  Google Scholar 

72.
Jucker, T. et al. Estimating aboveground carbon density and its uncertainty in Borneo’s structurally complex tropical forests using airborne laser scanning. Biogeosciences, 15, 3811–3830 (2018).
ADS  Article  Google Scholar 

73.
Khosravipour, A., Skidmore, A. K., Isenburg, M., Wang, T. & Hussin, Y. A. Generating pit-free canopy height models from airborne lidar. Photogramm. Eng. Remote Sens. 80, 863–872 (2014).
Article  Google Scholar 

74.
Asner, G. P. & Mascaro, J. Mapping tropical forest carbon: calibrating plot estimates to a simple LiDAR metric. Remote Sens. Environ. 140, 614–624 (2014).
ADS  Article  Google Scholar 

75.
Pearson, T. R. H., Brown, S. & Casarim, F. M. Carbon emissions from tropical forest degradation caused by logging. Environ. Res. Lett. 9, 034017 (2014).
ADS  Article  CAS  Google Scholar 

76.
Gobakken, T. G. & Næsset, E. N. Assessing effects of laser point density, ground sampling intensity, and field sample plot size on biophysical stand properties derived from airborne laser scanner data. Can. J. For. Res. https://doi.org/10.1139/X07-219 (2008).

77.
Gray, C. L., Slade, E. M., Mann, D. J. & Lewis, O. T. Do riparian reserves support dung beetle biodiversity and ecosystem services in oil palm-dominated tropical landscapes? Ecol. Evol. 4, 1049–1060 (2014).
PubMed  PubMed Central  Article  Google Scholar 

78.
Metcalfe, P., Beven, K. & Freer, J. Dynamic TOPMODEL: a new implementation in R and its sensitivity to time and space steps. Environ. Model. Softw. 72, 155–172 (2015).
Article  Google Scholar 

79.
Condit, R. et al. Tropical forest dynamics across a rainfall gradient and the impact of an El Niño Dry Season. J. Trop. Ecol. 20, 51–72 (2004).
Article  Google Scholar 

80.
H’edl, R. et al. A new technique for inventory of permanent plots in tropical forests: a case study from lowland dipterocarp forest in Kuala Belalong, Brunei Darussalam. Blumea J. Biodivers. Evolut. Biogeogr. Plants 54, 124–130 (2009).
Article  Google Scholar 

81.
Kent, R., Lindsell, J. A., Laurin, G. V., Valentini, R. & Coomes, D. A. Airborne LiDAR detects selectively logged tropical forest even in an advanced stage of recovery. Remote Sens. 7, 8348–8367 (2015).
ADS  Article  Google Scholar 

82.
Gonsamo, A., Walter, J.-M. N. & Pellikka, P. Sampling gap fraction and size for estimating leaf area and clumping indices from hemispherical photographs. Can. J. Res. 40, 1588–1603 (2010).
Article  Google Scholar 

83.
Kalacska, M. E. R., Sanchez-Azofeifa, G. A., Calvo-Alvarado, J. C., Rivard, B. & Quesada, M. Effects of season and successional stage on leaf area index and spectral vegetation indices in three mesoamerican tropical dry forests1. Biotropica 37, 486–496 (2005).
Article  Google Scholar 

84.
Thimonier, A., Sedivy, I. & Schleppi, P. Estimating leaf area index in different types of mature forest stands in Switzerland: a comparison of methods. Eur. J. Res. 129, 543–562 (2010).
Article  Google Scholar 

85.
Chen, J. M. & Cihlar, J. Quantifying the effect of canopy architecture on optical measurements of leaf area index using two gap size analysis methods. IEEE Trans. Geosci. Remote Sens. 33, 777–787 (1995).
ADS  Article  Google Scholar 

86.
Schleppi, P., Conedera, M., Sedivy, I. & Thimonier, A. Correcting non-linearity and slope effects in the estimation of the leaf area index of forests from hemispherical photographs. Agric. Meteorol. 144, 236–242 (2007).
Article  Google Scholar 

87.
Harmon, M. E., Whigham, D. F., Sexton, J. & Olmsted, I. Decomposition and mass of woody detritus in the dry tropical forests of the Northeastern Yucatan Peninsula, Mexico. Biotropica 27, 305–316 (1995).
Article  Google Scholar 

88.
Team, R. C. R: A Language and Environment for Statistical Computing. (Foundation for Statistical Computing, Vienna, Austria. 2017). Available online: www.r-project.org (accessed 14 Febuary 2019) (2018).

89.
Ploton, P. et al. Spatial validation reveals poor predictive performance of large-scale ecological mapping models. Nat. Commun. 11, 4540 (2020).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

90.
Wedeux, B. et al. Dynamics of a human-modified tropical peat swamp forest revealed by repeat lidar surveys. Glob. Chang. Biol. https://doi.org/10.1111/gcb.15108 (2020). More

1.
Grant, P. R. & Grant, B. R. How and why Species Multiply: The Radiation of Darwin’s Finches. (Princeton University Press, 2008).
2.
Baldwin, B. G. & Sanderson, M. J. Age and rate of diversification of the Hawaiian silversword alliance (Compositae). Proc. Natl Acad. Sci. USA 95, 9402–9406 (1998).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

3.
Losos, J. B. & Ricklefs, R. E. Adaptation and diversification on islands. Nature 457, 830–836 (2009).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

4.
Macarthur, R. H. & Wilson, E. O. The Theory of Island Biogeography. (Princeton University Press, 1967).

5.
Lewontin, R. C. The organism as the subject and object of evolution. Scientia 77, 65 (1983).
Google Scholar 

6.
Blows, M. W. & Hoffmann, A. A. A reassessment of genetic limits to evolutionary change. Ecology 86, 1371–1384 (2005).
Article  Google Scholar 

7.
Hansen, T. F. & Houle, D. Measuring and comparing evolvability and constraint in multivariate characters. J. Evol. Biol. 21, 1201–1219 (2008).
CAS  PubMed  Article  PubMed Central  Google Scholar 

8.
West-Eberhard, M. J. Developmental Plasticity and Evolution. (Oxford University Press, 2003).

9.
Wagner, G. P. & Altenberg, L. Perspective: complex adaptations and the evolution of evolvability. Evolution 50, 967–976 (1996).
PubMed  Article  PubMed Central  Google Scholar 

10.
Hendrikse, J. L., Parsons, T. E. & Hallgrímsson, B. Evolvability as the proper focus of evolutionary developmental biology. Evol. Dev. 9, 393–401 (2007).
PubMed  Article  PubMed Central  Google Scholar 

11.
Klingenberg, C. P. Studying morphological integration and modularity at multiple levels: concepts and analysis. Philos. Trans. R. Soc. B 369, 20130249 (2014).
Article  Google Scholar 

12.
Jablonski, D. Approaches to macroevolution: 1. General concepts and origin of variation. Evol. Biol. 44, 427–450 (2017).
PubMed  PubMed Central  Article  Google Scholar 

13.
Uller, T., Moczek, A. P., Watson, R. A., Brakefield, P. M. & Laland, K. N. Developmental bias and evolution: a regulatory network perspective. Genetics 209, 949–966 (2018).
PubMed  PubMed Central  Article  Google Scholar 

14.
Hansen, T. F. Is modularity necessary for evolvability? Remarks on the relationship between pleiotropy and evolvability. Biosystems 69, 83–94 (2003).
PubMed  Article  PubMed Central  Google Scholar 

15.
Goswami, A., Binder, W. J., Meachen, J. & O’Keefe, F. R. The fossil record of phenotypic and modularity: a deep-time perspective on developmental and evolutionary dynamics. Proc. Natl Acad. Sci. USA 112, 4891–4896 (2015).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

16.
Armbruster, W. S., Pelabon, C., Bolstad, G. H. & Hansen, T. F. Integrated phenotypes: understanding trait covariation in plants and animals. Philos. Trans. R. Soc. B 369, 20130245 (2014).
Article  Google Scholar 

17.
Felice, R. N., Randau, M. & Goswami, A. A fly in a tube: macroevolutionary expectations for integrated phenotypes. Evolution 72, 2580–2594 (2018).
PubMed  PubMed Central  Article  Google Scholar 

18.
Goswami, A., Smaers, J. B., Soligo, C. & Polly, P. D. The macroevolutionary consequences of phenotypic integration: from development to deep time. Philos. Trans. R. Soc. B 369, 20130254 (2014).
CAS  Article  Google Scholar 

19.
Cheverud, J. M. Phenotypic, genetic, and environmental morphological integration in the cranium. Evolution 36, 499–516 (1982).
PubMed  Article  PubMed Central  Google Scholar 

20.
Wagner, G. P., Pavlicev, M. & Cheverud, J. M. The road to modularity. Nat. Rev. Genet. 8, 921–931 (2007).
CAS  PubMed  Article  PubMed Central  Google Scholar 

21.
Melo, D., Porto, A., Cheverud, J. M. & Marroig, G. Modularity: genes, development and evolution. Annu. Rev. Ecol. Evol. Syst. 47, 463–486 (2016).
PubMed  PubMed Central  Article  Google Scholar 

22.
Gerhart, J. & Kirschner, M. The theory of facilitated variation. Proc. Natl Acad. Sci. USA 104, 8582–8589 (2007).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

23.
Villmoare, B., Fish, J. & Jungers, W. Selection, morphological integration, and strepsirrhine locomotor adaptations. Evol. Biol. 38, 88–99 (2011).
Article  Google Scholar 

24.
Navalon, G., Marugan-Lobon, J., Bright, J. A., Cooney, C. R. & Rayfield, E. J. The consequences of craniofacial integration for the adaptive radiations of Darwin’s finches and Hawaiian honeycreepers. Nat. Ecol. Evol. 4, 270–278 (2020).
PubMed  Article  PubMed Central  Google Scholar 

25.
Nicholson, K. E. et al. Mainland colonization by island lizards. J. Biogeogr. 32, 929–938 (2005).
Article  Google Scholar 

26.
Poe, S. et al. A phylogenetic, biogeographic, and taxonomic study of all extant species of Anolis (Squamata; Iguanidae). Syst. Biol. 66, 663–697 (2017).
PubMed  Article  PubMed Central  Google Scholar 

27.
Jackman, T., Losos, J. B., Larson, A. & de Queiroz, K. in Molecular Evolution and Adaptive Radiation (eds Givnish, T. & Systma, K.) 535–557 (Cambridge University Press, 1997).

28.
Underwood, G. The anoles of the Eastern Caribbean (Sauria, Iguanidae). Revisionary notes. Bull. Mus. Comp. Zool., Part III 121, 191–226 (1959).
Google Scholar 

29.
Losos, J. B. Lizards in an Evolutionary Tree: Ecology and Adaptive Radiation of Anoles. Vol. 10 (University of California Press, 2009).

30.
Pinto, G., Mahler, D. L., Harmon, L. J. & Losos, J. B. Testing the island effect in adaptive radiation: rates and patterns of morphological diversification in Caribbean and mainland Anolis lizards. Proc. R. Soc. B 275, 2749–2757 (2008).
PubMed  Article  PubMed Central  Google Scholar 

31.
Poe, S. & Anderson, C. G. The existence and evolution of morphotypes in Anolis lizards: coexistence patterns, not adaptive radiations, distinguish mainland and island faunas. PeerJ 6, e6040 (2019).
PubMed  PubMed Central  Article  Google Scholar 

32.
Irschick, D. J., Vitt, L. J., Zani, P. A. & Losos, J. B. A comparison of evolutionary radiations in mainland and Caribbean Anolis lizards. Ecology 78, 2191–2203 (1997).
Article  Google Scholar 

33.
Macrini, T. E., Irschick, D. J. & Losos, J. B. Ecomorphological differences in toepad characteristics between mainland and island anoles. J. Herpetol. 37, 52–58 (2003).
Article  Google Scholar 

34.
Velasco, J. A. & Herrel, A. Ecomorphology of Anolis lizards of the Choco’ region in Colombia and comparisons with Greater Antillean ecomorphs. Biol. Biol. J. Linn. Soc. 92, 403–403 (2007).
Article  Google Scholar 

35.
Williams, E. E. in Evol. Biol. Vol. 6 (eds Theodosius Dobzhansky, MaxK Hecht, & WilliamC Steere) Ch. 3, 47–89 (Springer US, 1972).

36.
Williams, E. E. in Lizard ecology: studies of a model organism (eds Pianka, E. R., Huey, R. B. & Schoener, T. W.) 326–370 (Harvard University Press, 1983).

37.
Losos, J. B., Jackman, T. R., Larson, A., Queiroz, K. & Rodriguez-Schettino, L. Contingency and determinism in replicated adaptive radiations of island lizards. Science 279, 2115–2118 (1998).
ADS  CAS  PubMed  Article  Google Scholar 

38.
Tinius, A. & Russell, A. P. Geometric morphometric analysis of the breast-shoulder apparatus of lizards: a test case using Jamaican anoles (Squamata: Dactyloidae). Anat. Rec. 297, 410–432 (2014).
Article  Google Scholar 

39.
Tinius, A., Russell, A. P., Jamniczky, H. A. & Anderson, J. S. What is bred in the bone: ecomorphological associations of pelvic girdle form in greater Antillean Anolis lizards. J. Morphol. 279, 1016–1030 (2018).
PubMed  Article  Google Scholar 

40.
Adams, D. C. & Collyer, M. L. Phylogenetic comparative methods and the evolution of multivariate phenotypes. Annu. Rev. Ecol. Evol. Syst. 50, 405–425 (2019).
Article  Google Scholar 

41.
Legendre, P. & Legendre, L. Numerical Ecology. (Elsevier, 2012).

42.
Collyer, M. L., Davis, M. A. & Adams, D. C. Making heads or tails of combined landmark configurations in geometric morphometric data. Evol. Biol. 47, 193–205 (2020).
Article  Google Scholar 

43.
Kumar, S., Stecher, G., Suleski, M. & Hedges, S. B. TimeTree: a resource for timelines, timetrees, and divergence times. Mol. Biol. Evol. 34, 1812–1819 (2017).
CAS  PubMed  Article  PubMed Central  Google Scholar 

44.
Nishimoto, S. & Logan, M. P. O. Subdivision of the lateral plate mesoderm and specification of the forelimb and hindlimb forming domains. Semin. Cell Dev. Biol. 49, 102–108 (2016).
PubMed  Article  PubMed Central  Google Scholar 

45.
Shou, S., Scott, V., Reed, C., Hitzemann, R. & Stadler, H. S. Transcriptome analysis of the murine forelimb and hindlimb autopod. Dev. Dyn. 234, 74–89 (2005).
CAS  PubMed  Article  PubMed Central  Google Scholar 

46.
Margulies, E. H., Kardia, S. L. R. & Innis, J. W. A comparative molecular analysis of developing mouse forelimbs and hindlimbs using Serial Analysis of Gene Expression (SAGE). Genome Res. 11, 1686–1698 (2001).
CAS  PubMed  PubMed Central  Article  Google Scholar 

47.
Adams, D. C. & Collyer, M. L. On the comparison of the strength of morphological integration across morphometric datasets. Evolution 70, 2623–2631 (2016).
PubMed  Article  PubMed Central  Google Scholar 

48.
Adams, D. C. Evaluating modularity in morphometric data: challenges with the RV coefficient and a new test measure. Methods Ecol. Evol. 7, 565–572 (2016).
Article  Google Scholar 

49.
Adams, D. C. & Collyer, M. L. Comparing the strength of modular signal, and evaluating alternative modular hypotheses, using covariance ratio effect sizes with morphometric data. Evolution 73, 2352–2367 (2019).
PubMed  Article  PubMed Central  Google Scholar 

50.
Dellinger, A. S. et al. Modularity increases rate of floral evolution and adaptive success for functionally specialized pollination systems. Commun. Biol. 2, 453 (2019).
PubMed  PubMed Central  Article  Google Scholar 

51.
Venditti, C., Meade, A. & Pagel, M. Multiple routes to mammalian diversity. Nature 479, 393–396 (2011).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

52.
Cooney, C. R. et al. Mega-evolutionary dynamics of the adaptive radiation of birds. Nature 542, 344 (2017).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

53.
Marki, P. Z., Kennedy, J. D., Cooney, C. R., Rahbek, C. & Fjeldsa, J. Adaptive radiation and the evolution of nectarivory in a large songbird clade. Evolution 73, 1226–1240 (2019).
PubMed  Article  PubMed Central  Google Scholar 

54.
Brown, R. L. What evolvability really is. Br. J. Philos. Sci. 65, 549–572 (2013).
MathSciNet  Article  Google Scholar 

55.
Watson, R. A. & Szathmary, E. How can evolution learn? Trends Ecol. Evol. 31, 147–157 (2016).
PubMed  Article  PubMed Central  Google Scholar 

56.
Young, N. M. & Hallgrimsson, B. Serial homology and the evolution of mammalian limb covariation structure. Evolution 59, 2691–2704 (2005).
PubMed  Article  PubMed Central  Google Scholar 

57.
Young, N. M., Wagner, G. P. & Hallgrimsson, B. Development and the evolvability of human limbs. Proc. Natl Acad. Sci. USA 107, 3400–3405 (2010).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

58.
Kelly, E. M. & Sears, K. E. Reduced phenotypic covariation in marsupial limbs and the implications for mammalian evolution. Biol. J. Linn. Soc. 102, 22–36 (2011).
Article  Google Scholar 

59.
Bennett, C. V. & Goswami, A. Does developmental strategy drive limb integration in marsupials and monotremes? Mamm. Biol. 76, 79–83 (2011).
Article  Google Scholar 

60.
Martin-Serra, A. & Benson, R. B. J. Developmental constraints do not influence long-term phenotypic evolution of marsupial forelimbs as revealed by interspecific disparity and integration patterns. Am. Nat. 195, 547–560 (2020).
PubMed  Article  PubMed Central  Google Scholar 

61.
Parter, M., Kashtan, N. & Alon, U. Facilitated variation: how evolution learns from past environments to generalize to new environments. PLoS Comp. Biol. 4, e1000206 (2008).
ADS  Article  CAS  Google Scholar 

62.
Kouvaris, K., Clune, J., Kounios, L., Brede, M. & Watson, R. A. How evolution learns to generalise: Using the principles of learning theory to understand the evolution of developmental organisation. PLoS Comp. Biol. 13, e1005358 (2017).
ADS  Article  CAS  Google Scholar 

63.
Brun-Usan, M., Rago, A., Thies, C., Uller, T. & Watson, R. A. Developmental models reveal the role of phenotypic plasticity in explaining genetic evolvability. bioRxiv https://doi.org/10.1101/2020.06.29.179226 (2020).

64.
Shanahan, T. Phylogenetic inertia and Darwin’s higher law. Stud. Hist. Philos. Sci. Part C 42, 60–68 (2011).
Article  Google Scholar 

65.
Houle, D., Bolstad, G. H., van der Linde, K. & Hansen, T. F. Mutation predicts 40 million years of fly wing evolution. Nature 548, 447–450 (2017).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

66.
Braendle, C., Baer, C. F. & Felix, M. A. Bias and evolution of the mutationally accessible phenotypic space in a developmental system. PLoS Genet. 6, e1000877 (2010).
PubMed  PubMed Central  Article  CAS  Google Scholar 

67.
Haber, A. Phenotypic covariation and morphological diversification in the ruminant skull. Am. Nat. 187, 576–591 (2016).
PubMed  Article  PubMed Central  Google Scholar 

68.
Schluter, D. Adaptive radiation along genetic lines of least resistance. Evolution 50, 1766–1774 (1996).
PubMed  Article  PubMed Central  Google Scholar 

69.
Hanot, P., Herrel, A., Guintard, C. & Cornette, R. The impact of artificial selection on morphological integration in the appendicular skeleton of domestic horses. J. Anat. 232, 657–673 (2018).
PubMed  PubMed Central  Article  Google Scholar 

70.
Penna, A., Melo, D., Bernardi, S., Oyarzabal, M. I. & Marroig, G. The evolution of phenotypic integration: How directional selection reshapes covariation in mice. Evolution 71, 2370–2380 (2017).
PubMed  PubMed Central  Article  Google Scholar 

71.
Watson, R. A., Wagner, G. P., Pavlicev, M., Weinreich, D. M. & Mills, R. The evolution of phenotypic correlations and “developmental memory”. Evolution 68, 1124–1138 (2014).
PubMed  PubMed Central  Article  Google Scholar 

72.
Donihue, C. M. et al. Hurricane effects on Neotropical lizards span geographic and phylogenetic scales. Proc. Natl Acad. Sci. USA 117, 10429–10434 (2020).
CAS  PubMed  Article  PubMed Central  Google Scholar 

73.
Feiner, N., Jackson, I. S. C., Munch, K. L., Radersma, R. & Uller, T. Plasticity and evolutionary convergence in the locomotor skeleton of Greater Antillean Anolis lizards. eLife 9, e57468 (2020).
CAS  PubMed  PubMed Central  Article  Google Scholar 

74.
Vanhooydonck, B. & Irschick, D. in Topics in functional and ecological vertebrate morphology (eds Aerts, P., D’Août, K., Herrel, A. & Van Damme, R.) (Shaker Publishing, 2002).

75.
Schluter, D. The Ecology of Adaptive Radiation. (Oxford: Oxford University Press, 2000).

76.
Roughgarden, J. Anolis Lizards of the Caribbean: Ecology, Evolution, and Plate Tectonics. (Oxford University Press, 1995).

77.
Stacklies, W., Redestig, H., Scholz, M., Walther, D. & Selbig, J. pcaMethods-a bioconductor package providing PCA methods for incomplete data. Bioinformatics 23, 1164–1167 (2007).
CAS  PubMed  Article  PubMed Central  Google Scholar 

78.
Losos, J. B. et al. Evolutionary implications of phenotypic plasticity in the hindlimb of the lizard Anolis sagrei. Evolution 54, 301–305 (2000).
CAS  PubMed  PubMed Central  Google Scholar 

79.
Tinius, A. Geometric morphometric analysis of the breast-shoulder apparatus of Greater Antillean anole ecomorphs PhD thesis, (University of Calgary, 2016).

80.
Cignoni, P. et al. in Eurographics Italian Chapter Conference (eds Scarano, V., De Chiara, R. & Erra, U.) (The Eurographics Association, 2008).

81.
Geomorph: Software for geometric morphometric analyses. R package version 3.1.0. (2019).

82.
Olsen, A. M. & Westneat, M. W. StereoMorph: an R package for the collection of 3D landmarks and curves using a stereo camera set-up. Methods Ecol. Evol. 6, 351–356 (2015).
Article  Google Scholar 

83.
Mahler, D. L., Ingram, T., Revell, L. J. & Losos, J. B. Exceptional convergence on the macroevolutionary landscape in island lizard radiations. Science 341, 292–295 (2013).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

84.
Rohlf, F. J. Shape statistics: procrustes superimpositions and tangent spaces. J. Classif. 16, 197–223 (1999).
MATH  Article  Google Scholar 

85.
Uetz, P., Freed, P. & Hosek, J. The Reptile Database http://www.reptile-database.org (2019).

86.
Pyron, R. A., Burbrink, F. T. & Wiens, J. J. A phylogeny and revised classification of Squamata, including 4161 species of lizards and snakes. BMC Evol. Biol. 13, 93 (2013).
PubMed  PubMed Central  Article  CAS  Google Scholar 

87.
Köhler, G. & Hedges, S. B. A revision of the green anoles of Hispaniola with description of eight new species (Reptilia, Squamata, Dactyloidae). Nov. Carib. 9, 1–135 (2016).
Google Scholar 

88.
Hofmann, E. P. & Townsend, J. H. Origins and biogeography of the Anolis crassulus subgroup (Squamata: Dactyloidae) in the highlands of Nuclear Central America. BMC Evol. Biol. 17, 267 (2017).
PubMed  PubMed Central  Article  CAS  Google Scholar 

89.
Mahler, D. L. et al. Discovery of a giant chameleon-like lizard (Anolis) on hispaniola and its significance to understanding replicated adaptive radiations. Am. Nat. 188, 357–364 (2016).
PubMed  Article  PubMed Central  Google Scholar 

90.
Kohler, J., Hahn, M. & Kohler, G. Divergent evolution of hemipenial morphology in two cryptic species of mainland anoles related to Anolis polylepis. Salamandra 48, 1–11 (2012).
Google Scholar 

91.
Kohler, G., Perez, R. G. T., Petersen, C. B. P. & De la Cruz, F. R. M. A revision of the Mexican Anolis (Reptilia, Squamata, Dactyloidae) from the Pacific versant west of the Isthmus de Tehuantepec in the states of Oaxaca, Guerrero, and Puebla, with the description of six new species. Zootaxa 3862, 1 (2014).
PubMed  Article  PubMed Central  Google Scholar 

92.
Revell, L. J. phytools: an R package for phylogenetic comparative biology (and other things). Methods Ecol. Evol. 3, 217–223 (2012).
Article  Google Scholar 

93.
Nicholson, K. E., Crother, B. I., Guyer, C. & Savage, J. M. It is time for a new classification of anoles (Squamata: Dactyloidae). Zootaxa 3477, 1–108 (2012).
Article  Google Scholar 

94.
Goswami, A. & Finarelli, J. A. EMMLi: a maximum likelihood approach to the analysis of modularity. Evolution 70, 1622–1637 (2016).
PubMed  Article  PubMed Central  Google Scholar 

95.
Bookstein, F. L. et al. Cranial integration in Homo: singular warps analysis of the midsagittal plane in ontogeny and evolution. J. Hum. Evol. 44, 167–187 (2003).
PubMed  Article  PubMed Central  Google Scholar 

96.
Adams, D. C. Quantifying and comparing phylogenetic evolutionary rates for shape and other high-dimensional phenotypic data. Syst. Biol. 63, 166–177 (2014).
PubMed  Article  PubMed Central  Google Scholar 

97.
Xie, W. G., Lewis, P. O., Fan, Y., Kuo, L. & Chen, M. H. Improving marginal likelihood estimation for bayesian phylogenetic model selection. Syst. Biol. 60, 150–160 (2011).
PubMed  Article  PubMed Central  Google Scholar 

98.
Brown, M. B. & Forsythe, A. B. Robust tests for the equality of variances. J. Am. Stat. Assoc. 69, 364–367 (1974).
MATH  Article  Google Scholar 

99.
Levene, H. in Contributions to Probability and Statistics (Stanford University Press, 1960).

100.
Kruskal, W. H. & Wallis, W. A. Use of ranks in one-criterion variance analysis. J. Am. Stat. Assoc. 47, 583–621 (1952).
MATH  Article  Google Scholar  More

1.
Edwards, K. J., Becker, K. & Colwell, F. The deep, dark energy biosphere: intraterrestrial life on Earth. Ann. Rev. Earth Planet Sci. 40, 551–568 (2012).
CAS  Article  Google Scholar 
2.
Kallmeyer, J., Pockalny, R., Adhikari, R. R., Smith, D. C. & D’Hondt, S. Global distribution of microbial abundance and biomass in subseafloor sediment. Proc. Nat. Acad. Sci. USA 109, 16213–16216 (2012).
CAS  PubMed  Article  PubMed Central  Google Scholar 

3.
Bar-On, Y. M., Phillips, R. & Milo, R. The biomass distribution on Earth. Proc. Nat. Acad. Sci. USA 115, 6506–6511 (2018).
CAS  PubMed  Article  PubMed Central  Google Scholar 

4.
Magnabosco, C. et al. The biomass and biodiversity of the continental subsurface. Nat. Geosci. 11, 707–717 (2018).
CAS  Article  Google Scholar 

5.
Lau, M. C. Y. et al. An oligotrophic deep-subsurface community dependent on syntrophy is dominated by sulfur-driven autotrophic denitrifiers. Proc. Nat. Acad. Sci. USA 113, E7927–E7936 (2016).
CAS  PubMed  Article  PubMed Central  Google Scholar 

6.
Lopez-Fernandez, M., Broman, E., Simone, D., Bertilsson, S. & Dopson, M. Statistical analysis of community RNA transcripts between organic carbon and ‘geogas’ fed continental deep biosphere groundwaters. mBio 10, e01470–01419 (2019).
CAS  PubMed  PubMed Central  Article  Google Scholar 

7.
Lopez-Fernandez, M. et al. Metatranscriptomes reveal all three domains of life are active, but are dominated by bacteria in the Fennoscandian crystalline granitic continental deep biosphere. mBio 9, e01792–01718 (2018).
PubMed  PubMed Central  Article  Google Scholar 

8.
Borgonie, G. et al. Eukaryotic opportunists dominate the deep-subsurface biosphere in South Africa. Nat. Comm. 6, 8952 (2015).
CAS  Article  Google Scholar 

9.
Wilkins, M. J. et al. Trends and future challenges in sampling the deep terrestrial biosphere. Front. Microbiol. 5, 481 (2014).
PubMed  PubMed Central  Google Scholar 

10.
Guemes, A. G. C. et al. Viruses as winners in the Game of Life. Ann. Rev. Virol. 3, 197–214 (2016).
Article  CAS  Google Scholar 

11.
Dávila-Ramos, S. et al. A review on viral metagenomics in extreme environments. Front. Microbiol. 10, 2403 (2019).
PubMed  PubMed Central  Article  Google Scholar 

12.
Roudnew, B. et al. Bacterial and virus-like particle abundances in purged and unpurged groundwater depth profiles. Ground Water Monit. Remed. 32, 72–77 (2012).
Article  Google Scholar 

13.
Nyyssönen, M. et al. Taxonomically and functionally diverse microbial communities in deep crystalline rocks of the Fennoscandian shield. ISME J. 8, 126–138 (2014).
PubMed  Article  CAS  PubMed Central  Google Scholar 

14.
Daly, R. A. et al. Microbial metabolisms in a 2.5-km-deep ecosystem created by hydraulic fracturing in shales. Nat. Microbiol. 1, 16146 (2016).
CAS  PubMed  Article  PubMed Central  Google Scholar 

15.
Anderson, R. E., Brazelton, W. J. & Baross, J. A. Is the genetic landscape of the deep subsurface biosphere affected by viruses? Front. Microbiol. 2, 219–219 (2011).
PubMed  PubMed Central  Article  Google Scholar 

16.
Anderson, R. E., Brazelton, W. J. & Baross, J. A. The deep viriosphere: assessing the viral impact on microbial community dynamics in the deep subsurface. Rev. Min. Geochem. 75, 649–675 (2013).
CAS  Article  Google Scholar 

17.
Labonté, J. M. et al. Single cell genomics indicates horizontal gene transfer and viral infections in a deep subsurface Firmicutes population. Front. Microbiol. 6, 349–349 (2015).
PubMed  PubMed Central  Google Scholar 

18.
Hallbeck, L. & Pedersen, K. Characterization of microbial processes in deep aquifers of the Fennoscandian Shield. Appl. Geochem. 23, 1796–1819 (2008).
CAS  Article  Google Scholar 

19.
Ström, A., Andersson, J., Skagius, K. & Winberg, A. Site descriptive modelling during characterization for a geological repository for nuclear waste in Sweden. Appl. Geochem. 23, 1747–1760 (2008).
Article  CAS  Google Scholar 

20.
Jägevall, S., Rabe, L. & Pedersen, K. Abundance and diversity of biofilms in natural and artificial aquifers of the Äspö Hard Rock Laboratory, Sweden. Microb. Ecol. 61, 410–422 (2011).
PubMed  Article  Google Scholar 

21.
Pedersen, K. Influence of H2 and O2 on sulphate-reducing activity of a subterranean community and the coupled response in redox potential. FEMS Microbiol. Ecol. 82, 653–665 (2012).
CAS  PubMed  Article  Google Scholar 

22.
Pedersen, K. Metabolic activity of subterranean microbial communities in deep granitic groundwater supplemented with methane and H2. ISME J. 7, 839–849 (2013).
CAS  PubMed  Article  PubMed Central  Google Scholar 

23.
Lopez-Fernandez, M., Broman, E., Wu, X., Bertilsson, S. & Dopson, M. Investigation of viable taxa in the deep terrestrial biosphere suggests high rates of nutrient recycling. FEMS Microbiol. Ecol. 94, fiy121 (2018).

24.
Lopez-Fernandez, M., Åström, M., Bertilsson, S. & Dopson, M. Depth and dissolved organic carbon shape microbial communities in surface influenced but not ancient saline terrestrial aquifers. Front. Microbiol. 9, 2880 (2018).
PubMed  PubMed Central  Article  Google Scholar 

25.
Wu, X. et al. Microbial metagenomes from three aquifers in the Fennoscandian shield terrestrial deep biosphere reveal metabolic partitioning among populations. ISME J. 10, 1192–1203 (2015).
PubMed  PubMed Central  Article  CAS  Google Scholar 

26.
Kyle, J. E., Eydal, H. S., Ferris, F. G. & Pedersen, K. Viruses in granitic groundwater from 69 to 450 m depth of the Äspö hard rock laboratory, Sweden. ISME J. 2, 571–574 (2008).
PubMed  Article  PubMed Central  Google Scholar 

27.
Eydal, H. S., Jagevall, S., Hermansson, M. & Pedersen, K. Bacteriophage lytic to Desulfovibrio aespoeensis isolated from deep groundwater. ISME J. 3, 1139–1147 (2009).
PubMed  Article  PubMed Central  Google Scholar 

28.
Castelle, C. J. et al. Biosynthetic capacity, metabolic variety and unusual biology in the CPR and DPANN radiations. Nat. Rev. Microbiol. 16, 629–645 (2018).
CAS  PubMed  Article  PubMed Central  Google Scholar 

29.
Hurwitz, B. L. & Sullivan, M. B. The Pacific Ocean Virome (POV): A marine viral metagenomic dataset and associated protein clusters for quantitative viral ecology. PLoS ONE 8, e57355 (2013).
CAS  PubMed  PubMed Central  Article  Google Scholar 

30.
Angly, F. E. et al. The marine viromes of four oceanic regions. PLoS Biol. 4, e368 (2006).
PubMed  PubMed Central  Article  CAS  Google Scholar 

31.
Holmfeldt, K. et al. Twelve previously unknown phage genera are ubiquitous in global oceans. Proc. Nat. Acad. Sci. USA 110, 12798–12803 (2013).
CAS  PubMed  Article  PubMed Central  Google Scholar 

32.
Nilsson, E. et al. Genomic and seasonal variations among aquatic phages infecting the Baltic Sea Gammaproteobacterium Rheinheimera sp. strain BAL341. Appl. Environ. Microbiol. 85, e01003-19, https://doi.org/10.1128/aem.01003-19 (2019).
PubMed  PubMed Central  Article  Google Scholar 

33.
Hurwitz, B. L., U’Ren, J. M. & Youens-Clark, K. Computational prospecting the great viral unknown. FEMS Microbiol. Lett. 363, fnw077 (2016).
PubMed  Article  CAS  PubMed Central  Google Scholar 

34.
Bolduc, B. et al. vConTACT: an iVirus tool to classify double-stranded DNA viruses that infect Archaea and Bacteria. PeerJ 5, e3243–e3243 (2017).
PubMed  PubMed Central  Article  Google Scholar 

35.
Lundin, D. & Holmfeldt, K. The deep terrestrial virosphere. Figshare, https://doi.org/10.6084/m6089.figshare.11590494.v11590491 (2020).

36.
Brum, J. R. et al. Patterns and ecological drivers of ocean viral communities. Science 348, 1261498 (2015).
PubMed  Article  CAS  PubMed Central  Google Scholar 

37.
Kadnikov, V. V. et al. Genomes of three bacteriophages from the deep subsurface aquifer. Data Brief. 22, 488–491 (2019).
PubMed  Article  PubMed Central  Google Scholar 

38.
Starnawski, P. et al. Microbial community assembly and evolution in subseafloor sediment. Proc. Nat. Acad. Sci. USA 114, 2940–2945 (2017).
CAS  PubMed  Article  PubMed Central  Google Scholar 

39.
Broman, E., Sjöstedt, J., Pinhassi, J. & Dopson, M. Shifts in coastal sediment oxygenation cause pronounced changes in microbial community composition and associated metabolism. Microbiome 5, 96 (2017).
PubMed  PubMed Central  Article  Google Scholar 

40.
Edwards, R. A., McNair, K., Faust, K., Raes, J. & Dutilh, B. E. Computational approaches to predict bacteriophage-host relationships. FEMS Microbiol. Rev. 40, 258–272 (2016).
CAS  PubMed  Article  PubMed Central  Google Scholar 

41.
Hug, L. A. et al. A new view of the tree of life. Nat. Microbiol. 1, 16048, https://doi.org/10.1038/nmicrobiol.2016.48 (2016).
CAS  PubMed  Article  PubMed Central  Google Scholar 

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

43.
Herrmann, M. et al. Predominance of Cand. Patescibacteria in groundwater is caused by their preferential mobilization from soils and flourishing under oligotrophic conditions. Front. Microbiol. 10, 1407 (2019).
PubMed  PubMed Central  Article  Google Scholar 

44.
Probst, A. J. et al. Differential depth distribution of microbial function and putative symbionts through sediment-hosted aquifers in the deep terrestrial subsurface. Nat. Microbiol. 3, 328–336 (2018).
CAS  PubMed  PubMed Central  Article  Google Scholar 

45.
Anantharaman, K. et al. Thousands of microbial genomes shed light on interconnected biogeochemical processes in an aquifer system. Nat. Comm. 7, https://doi.org/10.1038/ncomms13219 (2016).

46.
Bouvier, T. & del Giorgio, P. A. Key role of selective viral-induced mortality in determining marine bacterial community composition. Environ. Microbiol. 9, 287–297 (2007).
CAS  PubMed  Article  PubMed Central  Google Scholar 

47.
Dinsdale, E. A. et al. Functional metagenomic profiling of nine biomes. Nature 452, 629–632 (2008).
CAS  PubMed  PubMed Central  Article  Google Scholar 

48.
Craig, W. A. & Andes, D. R. in Mandell, Douglas, and Bennett’s Principles and Practice of Infectious Diseases (eds Bennett, J. E., Dolin, R. & Blaser, M. J.) 278–292.e274 (2015).

49.
Hubalek, V. et al. Connectivity to the surface determines diversity patterns in subsurface aquifers of the Fennoscandian shield. ISME J. 10, 2447–2458 (2016).
PubMed  PubMed Central  Article  Google Scholar 

50.
Laaksoharju, M., Gascoyne, M. & Gurban, I. Understanding groundwater chemistry using mixing models. Appl. Geochem. 23, 1921–1940 (2008).
CAS  Article  Google Scholar 

51.
Mathurin, F. A., Astrom, M. E., Laaksoharju, M., Kalinowski, B. E. & Tullborg, E. L. Effect of tunnel excavation on source and mixing of groundwater in a coastal granitoidic fracture network. Environ. Sci. Technol. 46, 12779–12786 (2012).
CAS  PubMed  Article  PubMed Central  Google Scholar 

52.
Smellie, J. A. T., Laaksoharju, M., Wikberg, P. & Äspö, S. E. Sweden—a natural groundwater-flow model derived from hydrogeological observations. J. Hydrol. 172, 147–169 (1995).
CAS  Article  Google Scholar 

53.
John, S. G. et al. A simple and efficient method for concentration of ocean viruses by chemical flocculation. Environ. Microbiol. Rep. 3, 195–202 (2011).
CAS  PubMed  PubMed Central  Article  Google Scholar 

54.
Bolger, A. M., Lohse, M. & Usadel, B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics 30, 2114–2120 (2014).
CAS  PubMed  PubMed Central  Article  Google Scholar 

55.
Boisvert, S., Raymond, F., Godzaridis, E., Laviolette, F. & Corbeil, J. Ray Meta: scalable de novo metagenome assembly and profiling. Genome Biol. 13, R122 (2012).
PubMed  PubMed Central  Article  CAS  Google Scholar 

56.
Rodriguez-R, L. M., Gunturu, S., Tiedje, J. M., Cole, J. R. & Konstantinidis, K. T. Nonpareil 3: fast estimation of metagenomic coverage and sequence diversity. mSystems 3, e00039–00018 (2018).
PubMed  PubMed Central  Article  Google Scholar 

57.
Roux, S., Enault, F., Hurwitz, B. L. & Sullivan, M. B. VirSorter: mining viral signal from microbial genomic data. PeerJ 3, e985–e985 (2015).
PubMed  PubMed Central  Article  CAS  Google Scholar 

58.
Brum, J. R. et al. Illuminating structural proteins in viral “dark matter” with metaproteomics. Proc. Nat. Acad. Sci. USA 113, 2436–2441 (2016).
CAS  PubMed  Article  PubMed Central  Google Scholar 

59.
Hugerth, L. W. et al. Metagenome-assembled genomes uncover a global brackish microbiome. Genome Biol. 16, 279 (2015).
PubMed  PubMed Central  Article  CAS  Google Scholar 

60.
Dupont, C. L. et al. Functional tradeoffs underpin salinity-driven divergence in microbial community composition. PloS One 9, e89549 (2014).
PubMed  PubMed Central  Article  CAS  Google Scholar 

61.
Chow, C. E., Winget, D. M., White, R. A., 3rd, Hallam, S. J. & Suttle, C. A. Combining genomic sequencing methods to explore viral diversity and reveal potential virus-host interactions. Front. Microbiol. 6, 265, https://doi.org/10.3389/fmicb.2015.00265 (2015).

62.
Tangherlini, M., Dell’Anno, A., Zeigler Allen, L., Riccioni, G. & Corinaldesi, C. Assessing viral taxonomic composition in benthic marine ecosystems: reliability and efficiency of different bioinformatic tools for viral metagenomic analyses. Sci. Rep. 6, 28428, https://doi.org/10.1038/srep28428 (2016).

63.
Sible, E. et al. Survey of viral populations within Lake Michigan nearshore waters at four Chicago area beaches. Data Brief. 5, 9–12 (2015).
PubMed  PubMed Central  Article  Google Scholar 

64.
Roux, S. et al. Ecogenomics and potential biogeochemical impacts of globally abundant ocean viruses. Nature 537, 689 (2016).
CAS  PubMed  PubMed Central  Article  Google Scholar 

65.
Wu, X. et al. Potential for hydrogen-oxidizing chemolithoautotrophic and diazotrophic populations to initiate biofilm formation in oligotrophic, deep terrestrial subsurface waters. Microbiome 5, 37, https://doi.org/10.1186/s40168-40017-40253-y (2017).

66.
Marcais, G. & Kingsford, C. A fast, lock-free approach for efficient parallel counting of occurrences of k-mers. Bioinformatics 27, 764–770 (2011).
CAS  PubMed  PubMed Central  Article  Google Scholar 

67.
Sullivan, M. J., Petty, N. K. & Beatson, S. A. Easyfig: a genome comparison visualizer. Bioinformatics 27, 1009–1010 (2011).
CAS  PubMed  PubMed Central  Article  Google Scholar 

68.
Simone, D. Domenico-simone/deep-metaviriomes: analysis for paper. Zenodo https://doi.org/10.5281/zenodo.3700451 (2020).
Article  Google Scholar  More

A total number of 3189 lakes have been spatially delimited and analyzed. The delimitation and spatial distribution of the lakes revealed their uneven distribution within the Romanian territory. The largest share of lakes, (41.5%) are distributed within the low plains, located mainly in the South and West of the country. About 36% of the identified lakes are located in the hilly and plateau units, while 11.5% are in the mountain areas and about 11% in the Danube Delta.
The assessment of the state of degradation of the lentic ecosystems by the GIS based extended WRASTIC model revealed that more than half (57%) of the analyzed lakes are classified as semi-degraded. These are mostly distributed in the plains (46%) and in the plateau areas (33%) (Fig. 1).
Figure 1

Distribution of lentic ecosystems in Romania according to the major topography units (this map was created with Arc GIS 10.5 software).

Full size image

The lakes classified as degraded represent 31% of the total, mostly located in the plains (48%) and plateau (36.5%), the least part being located in the Danube Delta (about 0.5%). The lakes in the natural state represent a small share, respectively about 13%, generally located in the Danube Delta (72%) and the mountain areas (21%) (Supplementary Table 1).
Altitudinal stages represent a restrictive factor in the spatial distribution of lakes throughout the country. Most lakes, about 76%, are located at low altitudes, respectively below an altitude of 200 m. The number of lakes identified at higher altitudes, over 800 m, is reduced, representing about 8% of the total number of analyzed lakes (Fig. 2).
Figure 2

Distribution of lentic ecosystems on the Romanian territory according to the altitudinal stages (this map was created with Arc GIS 10.5 software).

Full size image

Regarding the natural lakes, the largest part of them is represented by natural lakes located at altitudes less than 200 m (about 83%). In the same time, about 75% of the lakes evaluated as degraded are located at altitudes less than 200 m.
About 98.5% of the 412 lakes in the natural state are located, partially or totally, in protected areas of national or international interest, such as national and natural parks, Ramsar sites, Biosphere reserves, UNESCO World heritage sites and Natura 2000 sites (Sites of Community Importance and Special Protection Areas) (Fig. 3.).
Figure 3

Distribution of lake categories according to the state of degradation and the relationship with protected areas of national and international interest in Romania (this map was created with Arc GIS 10.5 software).

Full size image

The percentage of semi-degraded and degraded lakes located in protected areas is lower (about 38% of the semi-degraded lakes, respectively 32.5% of the degraded) (Supplementary Table 2). The other lakes are not included in protected areas and do not own any special protection regime.
As a result of applying the methodology based on WRASTIC-HI index, only 9 lakes out of the total of 412 natural lakes are associated with industrial activities and different forms of small-scale exploitation within their hydrographic basins, while the rest of 98% do not present such activities.
The hydrographic basins of the natural lakes do not cover very large areas, 83% of them stretching on area of less than 39 square km, each. In about 97% of cases, agricultural activities in the reception basin represent below 20% of the economic activity, the irrigation percentage being low and the degree of vegetation cover being high.
Considering the degraded lakes, for a large share of them (87%) of them industrial and exploitation activities, such as mines, quarries or dumps are present within the river basin. A high percentage, about 98%, include treatment plants with different types of processing (primary or secondary) within the reception hydrographic basins. For 96% of the degraded lakes, agricultural activities and permanent irrigation activities are covering over 40% of the related basins area.
Spearman correlation showed that there is a good correlation between the degradation state and several components of the WRASTIC-HI index such as Industrial activities, Recreational activities, Wastewater, Ways of transportation, Irrigation and Agricultural activities (Table 1).
Table 1 Correlation between Degradation state and component indices of the WRASTIC-HI index.
Full size table

Between the state of degradation, on the one hand, and the permeability of the soil, the slope and the vegetation cover of the water lily, on the other hand, the analysis showed that there is no correlation. However, there is a weak correlation between the degradation state and Exposition (Table 2).
Table 2 Correlation between degradation state and the HI index variables, which belong to the WRASTIC-HI index.
Full size table

The computational results indicate a negative correlation between the state of degradation, the percentage of the basin included in the protected areas and the percentage of the basin included in the protected areas of the Natura 2000 network (the probability level of 0.0001, being less than 0.05, indicates that there is a statistical significance) (Table 3).
Table 3 The results of the statistical analysis regarding the relationship between the state of degradation and variables within the river basin.
Full size table

The statistical analysis performed shows that there is no correlation between the Degradation State and Altitude, while between the Degradation state and the Relief units there is a weak negative correlation.
The results indicate also a direct correlation between the state of degradation, the number of inhabitants without access to the sewerage and the population density in the lake basin, the correlation coefficient being 0.024, and respectively 0.235.
It is important to study of the state of the lentic ecosystems both regarding the pressures coming from human activity and the ones originating from climatic changes or other natural phenomena. However, the development of an assessment model that can be applied to all aquatic ecosystems, and lentic ecosystem in particular, is a challenge, because available data are not homogenous, each region and lake having its own particularities. Choosing the proper set of indicators that can be useful for the overall characterization of the quality of lentic ecosystems, will lead to a coherent implementation of biodiversity strategies.
To the best of our knowledge, this is the first evaluation of the quality of lentic ecosystems in Romania by multi-criteria analysis. This assessment is part of a national project for implementation of a national policy on biodiversity, as a response to EU requirements.
Over the time, the scientific literature identified various methodologies for the evaluation of the quality of water resources. Such methods include DRASTIC20,21, and methods derived from DRASTIC22,23,24. Other authors used digital surface model (DSM) and a point dataset as the sources of observation and target locations for Geospatial analysis of lake scenery25.
For integrative water quality management, Feng et al.26, proposed a model-based method. Their method integrated three indices derived from three models for assessment of the risk due to nutrient dynamics26.
A multi-attribute value theory to formulate an integrated water quality assessment method was used by Schuwirth27, for aggregation over multiple pollutants and time.
Another category of indexes that are used in evaluating the state of degradation of lacustrine ecosystems aim at analysis of the presence of degradation sources in supplying basin. For example, Mirzaei et al19 calculated WRASTIC index for assessment of pollution risk, respectively degradation sources, from the watershed that feeds a water body. Potential pollutant load index (PPL) was employed by Romanelli et al18 for analyzing the presence and intensity of potential pollution sources from the drainage area of several lakes, with the purpose of establishing degradation classes of the water body. In the same study, the Lake Vulnerability Index highlighted the capacity of the water body to handle the impact generated by degradation sources, taking into account parameters like slope, soil permeability or aspect of slopes18.
According to the characteristics of the study area, modified versions of WRASTIC index are to be found in the literature, and were implemented by using additional criteria or eliminating some parameters28.
Rahimi et al29 used WRASTIC Index for evaluation of wetland water quality. Their results revealed that the activity in adjacent wetland areas exert a large impact on wetland integrity29.
In another study, aquatic ecosystems pollution risk was evaluated by a combined Fuzzy-WRASTIC method. The model was validated by comparison with samples collected from the case study area. The authors concluded that the method has advantages over other methods, as it includes a wide range of drivers and parameters that influence the water quality. The results obtained pointed that areas with high contamination risk are due to the unbalanced arrangement and compact of land uses in the neighborhood of the aquatic ecosystems30. Using analytical survey and experimental studies Mirzaei et al19 investigated the pollution risk for Zayandehrud river, Iran. Agricultural, industrial activities and population centers were the main causes of pollution in the study case area19.
In Romania, the evaluation of aquatic ecosystems was performed by various researchers, either for a specific area or for a hydrographic basin.
For example, Rosca et al31 studied the impact of anthropogenic activity on water quality parameters of glacial lakes from Rodnei mountains. The factors taken in consideration were tourism and livestock. The pollution index was calculated based on three indices, targeted on heavy metal influences, namely, the heavy metal pollution index (indicating the quality of waters related to the heavy metals content), the heavy metal evaluation index (assessment of the quality of water with respect to heavy metals) and the degree of contamination (used to quantify the contamination level with the heavy metal). The physico-chemical parameters pointed a good quality of the study case lakes. The conclusion of the authors was that minor anthropic alteration and a low anthropogenic impact is exerted in these areas. The only anthropic pressure on the aquatic systems in Rodnei Mountains was reported as being exerted by grazing activities31.
Another paper described the assessment of actual water quality and sedimentological conditions of the Corbu lake, Western Black Sea coast. The ecological status of this lake was found to be from good to weak classes for nitrites, ammonium and phosphates, moderate for sulphates and weak for detergents32.
The impact of human interventions and climate changes on the hydro-chemical composition of Techirghiol lake (Romania) was recently investigated by Maftei et al33. The study identified a degradation of this ecosystem between 1970–1998, due to extensive irrigation in the lake region, followed by a major decrease of the lake’s salinity33. Physico-chemical water quality parameters of lake Brăneşti was investigated by Benciu et al34. The water quality parameters for the last 50 years were correlated with the anthropogenic pressure in the region. Analysis of water and soil samples in the vicinity of this lake, revealed that parameters were within legal norms for both water and soil34.
Another study presented by Dumitran et al35 proposed an eutrophication model for describing the ecological behavior of a eutrophic lake. The physical model was mathematically transposed to a set of equations for analysing the selected parameters linked to eutrophication state. The resulted model showed a good correlation with the measured data35.
It is to be noted that most of the available literature is based on the assessment of the water quality, by measurements of physico-chemical parameters, and calculation of pollution indexes. To our knowledge no extensive studies involved the study of the lentic ecosystems with respect to vulnerability and risk of pollution by using multicriterial analysis.
Generally, the precautionary approach is applied by identifying and analyzing the categories of drivers that influence the degradation of lentic ecosystems, especially in the protected areas36,37,38. Three main categories of activities that generate environmental issues have been identified within protected areas included within the Natura 2000 Sites as follows: (a) agricultural activities and forestry practice; (b) sectoral activities (industrial, commercial and tourism sectors); (c) conservation policies (management of protected areas, protection of different species, etc.)39. A cross-sectoral approach is needed in order to resolve medium-term environmental conflicts, thus being be able to extend the assessment towards various categories of protected areas and generating efficient policies for the management of resources40,41.
Identifying and analyzing the categories of conflicts that may be associated with lentic ecosystems provide the possibility of an efficient ecosystem management22,23.
The lakes from this case study comprise both natural lakes (glacial, karst, karst-saline, ponds, lagoons), as well as ponds accumulation lakes, with an important role in ensuring the resources of water for the population and economic activities, as well as the development and maintenance of habitats and species of community interest (birds, amphibians, reptiles, fish, etc.).
Most of the lakes resulting from the analysis as being degraded and semi-degraded are located in the plains, at low altitudes, within areas covered with agricultural lands, industrial facilities and dense transportation routes.
The lentic ecosystems characterized as being in a natural state resulting from the proposed methodology are mainly distributed in the high mountain areas and in the Danube Delta area. Thus, altitude, fragmentation of the relief and accessibility are favorable factors regarding the natural state of water bodies, including lentic ecosystems. The high number of lakes characterized as degraded or semi-degraded compared to that of natural lakes is justified by the existence of a small number of lakes located at high altitudes, over 800 m.
It is to be taken into account that the degradation state classification is directly influenced by the data used in defining WRASTIC indicators, being generally derived data, which may explain the limitation of the method from this point of view.
The lack of correlation or poor correlation resulting from the statistical analysis between the degradation state and the indicators defining the HI index (component part of the WRASTIC-HI index), respectively the permeability, the exposure and the slope, highlight the insignificant role of these parameters in determining the state of lake degradation. However, the processes of erosion and sediment transport on the surface of the basins and their accumulation in lakes can influence the water quality of the lakes, their clogging, their functionality and the services offered, which are important factors in improving the management of the analyzed lakes42.
The status of protected areas offers a high degree of protection by diminishing the anthropic activities and the negative effects on the lakes, being recommended that all economic activities be located outside these protected areas43 the basins being vulnerable to human activities. This aspect is also highlighted by the correlation between the state of degradation obtained and the percentage of the hydrographic basin existing in different categories of protected areas, between which there is a good correlation, as well as by the high number of natural lakes that are included in protected areas (~ 98%), located especially in the Danube Delta Biosphere Reserve.
The local public administrations are directly interested in the management and protection of lentic ecosystems, many lakes being included in different categories of protected areas. Increasing the number of lakes in the natural state involves identifying degraded or semi-degraded lentic ecosystems outside protected areas and carrying out ecological reconstruction activities or diminishing agricultural and industrial activities in their vicinity.
The state of the lakes may also depend on the dynamics of the hydrophilic and hydrophilic vegetation, respectively on the hedge with vegetation cover of the water lily. In our study, the statistical analysis showed that there is no correlation between the state of degradation and the Coverage with vegetation of the water lily. Thus, in this case, the state of degradation of the lentic ecosystems is not influenced by this parameter, although, in some cases, remote sensing analysis revealed the presence of excess algae and aquatic plants in both natural and semi-degraded lakes44.
The lakes located in the low plain areas are also affected by the eutrophication process, amplified by the reduced depth (which ensures the rapid development of algae during the summer), the contribution of nutrients due to the agricultural activities in the vicinity and the development of recreational activities. More

1.
Clavero, M., Brotons, L., Pons, P. & Sol, D. Prominent role of invasive species in avian biodiversity loss. Biol. Conserv. 142, 2043–2049 (2009).
Article  Google Scholar 
2.
Clavero, M. & García-Berthou, E. Homogenization dynamics and introduction routes of invasive freshwater fish in the Iberian Peninsula. Ecol. Appl. 16, 2313–2324 (2006).
PubMed  Article  PubMed Central  Google Scholar 

3.
Tershy, B. R., Shen, K.-W., Newton, K. M., Holmes, N. D. & Croll, D. A. The Importance of islands for the protection of biological and linguistic diversity. Bioscience 65, 592–597 (2015).
Article  Google Scholar 

4.
Jones, H. P. Seabird islands take mere decades to recover following rat eradication. Ecol. Appl. 20, 2075–2080 (2010).
PubMed  Article  PubMed Central  Google Scholar 

5.
Wolf, C. A. et al. Invasive rat eradication strongly impacts plant recruitment on a tropical atoll. PLoS ONE 13, e0200743 (2018).
PubMed  PubMed Central  Article  Google Scholar 

6.
O’Dowd, D. J., Green, P. T. & Lake, P. S. Invasional ‘meltdown’ on an oceanic island. Ecol. Lett. 6, 812–817 (2003).
Article  Google Scholar 

7.
Rogers, H. S. et al. Effects of an invasive predator cascade to plants via mutualism disruption. Nat. Commun. 8, 14557 (2017).
ADS  PubMed  PubMed Central  Article  Google Scholar 

8.
Jones, H. P. et al. Invasive mammal eradication on islands results in substantial conservation gains. Proc. Natl. Acad. Sci. 113, 4033–4038 (2016).
ADS  CAS  PubMed  Article  Google Scholar 

9.
Towns, D. Eradications as reverse invasions: lessons from Pacific rat (Rattus exulans) removals on New Zealand islands. Biol. Invasions 11, 1719–1733 (2008).
Article  Google Scholar 

10.
Donlan, C. J., Croll, D. A. & Tershy, B. R. Islands, exotic herbivores, and invasive plants: their roles in coastal California Restoration. Restor. Ecol. 11, 524–530 (2003).
Article  Google Scholar 

11.
Tabak, M. A., Poncet, S., Passfield, K., Goheen, J. R. & del Rio, C. M. The ghost of invasives past: rat eradication and the community composition and energy flow of island bird communities. Ecosphere 7, e01442 (2016).
Article  Google Scholar 

12.
Kurle, C. M., Croll, D. A. & Tershy, B. R. Introduced rats indirectly change marine rocky intertidal communities from algae- to invertebrate-dominated. Proc. Natl. Acad. Sci. 105, 3800–3804 (2008).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

13.
Thoresen, J. J. et al. Invasive rodents have multiple indirect effects on seabird island invertebrate food web structure. Ecol. Appl. 27, 1190–1198 (2017).
PubMed  Article  PubMed Central  Google Scholar 

14.
Russell, J. Indirect effects of introduced predators on seabird islands. In Seabird Islands: Ecology, Invasion, and Restoration (eds Mulder, C. et al.) (Oxford University Press, Oxford, 2011).
Google Scholar 

15.
Le Corre, M. et al. Seabird recovery and vegetation dynamics after Norway rat eradication at Tromelin Island, western Indian Ocean. Biol. Conserv. 185, 85–94 (2015).
Article  Google Scholar 

16.
Doherty, T. S., Glen, A. S., Nimmo, D. G., Ritchie, E. G. & Dickman, C. R. Invasive predators and global biodiversity loss. Proc. Natl. Acad. Sci. 113, 11261–11265 (2016).
CAS  PubMed  Article  PubMed Central  Google Scholar 

17.
Bellard, C., Genovesi, P. & Jeschke, J. M. Global patterns in threats to vertebrates by biological invasions. Proc. R. Soc. B Biol. Sci. 283, 20152454 (2016).
Article  Google Scholar 

18.
Towns, D. R., Atkinson, I. A. E. & Daugherty, C. H. Have the harmful effects of introduced rats on islands been exaggerated?. Biol. Invasions 8, 863–891 (2006).
Article  Google Scholar 

19.
Jones, H. P. et al. Severity of the effects of invasive rats on seabirds: a global review. Conserv. Biol. 22, 16–26 (2008).
PubMed  Article  PubMed Central  Google Scholar 

20.
Drake, D. R. et al. Direct Impacts of Seabird Predators on Island Biota other than Seabirds. In Seabird Islands: Ecology, Invasion, and Restoration Mulder (eds Anderson, C. P. H. et al.) 91–132 (Oxford University Press, Oxford, 2011). https://doi.org/10.1093/acprof:osobl/9780199735693.003.0004.
Google Scholar 

21.
Towns, D. R. et al. Impacts of Introduced Predators on Seabirds. In Seabird Islands: Ecology, Invasion, and Restoration Mulder (eds Anderson, C. P. H. et al.) 56–90 (Oxford University Press, Oxford, 2011). https://doi.org/10.1093/acprof:osobl/9780199735693.003.0003.
Google Scholar 

22.
Mulder, C. P. H., Anderson, W. B., Towns, D. R. & Bellingham, P. J. Seabird Islands: Ecology, Invasion, and Restoration (Oxford University Press, Oxford, 2011).
Google Scholar 

23.
Croll, D. A., Maron, J. L., Estes, J. A., Danner, E. M. & Byrd, G. V. Introduced predators transform subarctic islands from grassland to Tundra. Science 307, 1959–1961 (2005).
ADS  CAS  PubMed  Article  PubMed Central  Google Scholar 

24.
Aslan, C. E., Zavaleta, E. S., Tershy, B. & Croll, D. Mutualism disruption threatens global plant biodiversity: a systematic review. PLoS ONE 8, e66993 (2013).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

25.
DIISE. The Database of Island Invasive Species Eradications, developed by Island Conservation, Coastal Conservation Action Laboratory UCSC, IUCN SSC Invasive Species Specialist Group, University of Auckland and Landcare Research New Zealand. http://diise.islandconservation.org/ (2018).

26.
Howald, G. et al. Invasive rodent eradication on islands. Conserv. Biol. 21, 1258–1268 (2007).
PubMed  Article  PubMed Central  Google Scholar 

27.
Keitt, B. et al. The Global Islands Invasive Vertebrate Eradication Database: A tool to improve and facilitate restoration of island ecosystems. In Island Invasives: Eradication and Management (eds Veitch, C. et al.) 4 (IUCN, Gland, 2011).
Google Scholar 

28.
Nigro, K. M. et al. Stable isotope analysis as an early monitoring tool for community-scale effects of rat eradication. Restor. Ecol. 25, 1015–1025 (2017).
Article  Google Scholar 

29.
Courchamp, F. et al. Eradication of alien invasive species: surprise effects and conservation successes. In Island Invasives: Eradication and Management (eds Veitch, C. et al.) 285–289 (IUCN, Gland, 2011).
Google Scholar 

30.
Jones, H. & Schmitz, O. Rapid recovery of damaged ecosystems. PLoS ONE 4, e5653 (2009).
ADS  PubMed  PubMed Central  Article  CAS  Google Scholar 

31.
Jones, H. P. et al. Recovery and Restoration on Seabird Islands. In Seabird Islands: Ecology, Invasion, and Restoration (eds Mulder, C. et al.) (Oxford University Press, Oxford, 2011).
Google Scholar 

32.
Buckelew, S., Byrd, V., Howald, G., MacLean, S. & Sheppard, J. Preliminary ecosystem response following invasive Norway rat eradication on Rat Island, Aleutian Islands, Alaska. in Island Invasives: eradicaation and management 5 (IUCN, 2011).

33.
Croll, D. A. et al. Passive recovery of an island bird community after rodent eradication. Biol. Invasions 18, 703–715 (2016).
Article  Google Scholar 

34.
Hanson, K., Goos, M. & Deines, F. G. Introduced arctic fox eradication at Rat Island (Aleutian Islands, Alaska, 1984).
Google Scholar 

35.
ESRI. ESRI ArcMap 10.7.0.10450. (ESRI, 2020).

36.
Lorvelec, O. & Pascal, M. French attempts to eradicate non-indigenous mammals and their consequences for native biota. Biol. Invasions 7, 135–140 (2005).
Article  Google Scholar 

37.
Bellingham, P. J. et al. New Zealand island restoration: seabirds, predators, and the importance of history. N. Z. J. Ecol. 34, 115 (2010).
Google Scholar 

38.
St. Clair, J., Poncet, S., Sheehan, D., Szekely, T. & Hilton, G. Responses of an island endemic invertebrate to rodent invasion and eradication. Anim. Conserv. 14, 66–73 (2011).
Article  Google Scholar 

39.
Monks, J. M., Monks, A. & Towns, D. R. Correlated recovery of five lizard populations following eradication of invasive mammals. Biol. Invasions 16, 167–175 (2014).
Article  Google Scholar 

40.
Whitworth, D. L., Carter, H. R. & Gress, F. Recovery of a threatened seabird after eradication of an introduced predator: Eight years of progress for Scripps’s murrelet at Anacapa Island, California. Biol. Conserv. 162, 52–59 (2013).
Article  Google Scholar 

41.
Brooke, M. L. et al. Seabird population changes following mammal eradications on islands. Anim. Conserv. 21, 3–12 (2017).
Article  Google Scholar 

42.
Bailey, E. P. Introduction of foxes to Alaskan Islands: history, effects on Avifauna, and Eradication. (U.S. Dept. of the Interior, Fish and Wildlife Service ; National Technical Information Service, distributor, 1993).

43.
Byrd G. V., Trapp, J. L., & Zeillemaker, C. F. Removal of Introduced Foxes: A Case Study in Restoration of Native Birds. in vol. 59 317–321 (1994).

44.
Byrd, G. V., Bailey, E. P. & Stahl, W. Restoration of island populations of black oystercatchers and pigeon guillemots by removing introduced foxes. Colon. Waterbirds 20, 253–260 (1997).
Article  Google Scholar 

45.
Ehrenfeld, J. G. Ecosystem consequences of biological invasions. Annu. Rev. Ecol. Evol. Syst. 41, 59–80 (2010).
Article  Google Scholar 

46.
Wootton, J. T. Indirect effects, prey susceptibility, and habitat selection: impacts of birds on limpets and algae. Ecology 73, 981–991 (1992).
Article  Google Scholar 

47.
Ellis, J. C., Chen, W., O’Keefe, B., Shulman, M. J. & Witman, J. D. Predation by gulls on crabs in rocky intertidal and shallow subtidal zones of the Gulf of Maine. J. Exp. Mar. Biol. Ecol. 324, 31–43 (2005).
Article  Google Scholar 

48.
Menge, B. A. Top-down and bottom-up community regulation in marine rocky intertidal habitats. J. Exp. Mar. Biol. Ecol. 250, 257–289 (2000).
CAS  PubMed  Article  PubMed Central  Google Scholar 

49.
Guerry, A. D., Menge, B. A. & Dunmore, R. A. Effects of consumers and enrichment on abundance and diversity of benthic algae in a rocky intertidal community. J. Exp. Mar. Biol. Ecol. 369(2), 155–164 (2009).
Article  Google Scholar 

50.
Wootton, J. T. Effects of birds on sea urchins and algae: a lower-intertidal trophic cascade. Écoscience 2, 321–328 (1995).
Article  Google Scholar 

51.
Ellis, J. C., Shulman, M. J., Wood, M., Witman, J. D. & Lozyniak, S. Regulation of intertidal food webs by avian predators on new england rocky shores. Ecology 88, 853–863 (2007).
PubMed  Article  PubMed Central  Google Scholar 

52.
Freidenburg, T. L., Menge, B. A., Halpin, P. M., Webster, M. & Sutton-Grier, A. Cross-scale variation in top-down and bottom-up control of algal abundance. J. Exp. Mar. Biol. Ecol. 347(1–2), 8–29 (2007).
Article  Google Scholar 

53.
Webster, J. D. Feeding habits of the black oyster-catcher. Condor 43, 175–180 (1941).
Article  Google Scholar 

54.
Trapp, J. L. Variation in summer diet of Glaucous-winged Gulls in the Western Aleutian Islands: an ecological interpretation. Wilson Bull. 91, 412–419 (1979).
Google Scholar 

55.
Irons, D. B., Anthony, R. G. & Estes, J. A. Foraging strategies of Glaucous-winged gulls in a rocky intertidal community. Ecology 67, 1460–1474 (1986).
Article  Google Scholar 

56.
Davis, M. L., Elliott, J. E. & Williams, T. D. Spatial and temporal variation in the dietary ecology of the Glaucous-winged Gull Larus Glaucescens in the Pacific Northwest. Mar. Ornithol. 43, 189–198 (2015).
Google Scholar 

57.
Padilla, D. K. The importance of form: differences in competitive ability, resistance to consumers and environmental stress in an assemblage of coralline algae. J. Exp. Mar. Biol. Ecol. 79, 105–127 (1984).
Article  Google Scholar 

58.
Breitburg, D. Residual effects of grazing – inhibition of competitor recruitment by encrusting coralline algae. Ecology 65, 1136–1143 (1984).
Article  Google Scholar 

59.
Scheibling, R. E. & Hatcher, B. G. Strongylocentrotus droebachiensis. In Developments in Aquaculture and Fisheries Science Vol. 38 (ed. Lawrence, J. M.) 381–412 (Elsevier, Amsterdam, 2013).
Google Scholar 

60.
Estes, J., Tinker, M., Williams, T. & Doak, D. F. Killer whale predation on sea otters linking oceanic and nearshore ecosystems | science. Science 282, 473–476 (1998).
ADS  CAS  PubMed  PubMed Central  Article  Google Scholar 

61.
Estes, J. A., Tinker, M. T. & Bodkin, J. L. Using ecological function to develop recovery criteria for depleted species: sea otters and kelp forests in the aleutian archipelago. Conserv. Biol. 24, 852–860 (2010).
PubMed  Article  PubMed Central  Google Scholar 

62.
Stewart, N. L., Konar, B. & Tinker, M. T. Testing the nutritional-limitation, predator-avoidance, and storm-avoidance hypotheses for restricted sea otter habitat use in the Aleutian Islands, Alaska. Oecologia 177, 645–655 (2015).
ADS  PubMed  Article  PubMed Central  Google Scholar 

63.
Gentemann, C. M., Fewings, M. R. & García-Reyes, M. Satellite sea surface temperatures along the West Coast of the United States during the 2014–2016 northeast Pacific marine heat wave. Geophys Res Lett 44, 312–319 (2017).
ADS  Article  Google Scholar 

64.
Coletti, H. et al. Gulf Watch Alaska: Nearshore Ecosystems in the Gulf of Alaska. Exxon Valdez Oil Spill Restoration Project Annual Report (Restoration Project 18120114-H), Exxon Valdez Oil Spill Trustee Council, Anchorage, Alaska (2019).

65.
Coletti, H. et al. Gulf Watch Alaska: Nearshore Ecosystems in the Gulf of Alaska. Exxon Valdez Oil Spill Restoration Project Annual Report (Restoration Project 18120114-H), Exxon Valdez Oil Spill Trustee Council, Anchorage, Alaska (2020).

66.
Hewson, I. et al. Investigating the Complex Association Between Viral Ecology, Environment, and Northeast Pacific Sea Star Wasting. Front. Mar. Sci. 5, 2018 (2018).
Article  Google Scholar 

67.
Elton, C. S. The Ecology of Invasions by Animals and Plants (Springer, Berlin, 1958).
Google Scholar 

68.
Richardson, D. M. & Pysek, P. Fifty Years of Invasion Ecology: The Legacy of Charles Elton (Blackwell Publishing Ltd., Hoboken, 2008).
Google Scholar 

69.
Courchamp, F. et al. Invasion biology: specific problems and possible solutions. Trends Ecol. Evol. 32, 13–22 (2017).
PubMed  Article  PubMed Central  Google Scholar 

70.
Cassini, M. H. A review of the critics of invasion biology. Biol. Rev. (2020).

71.
Kurle, C. M. Description of the rocky intertidal communities and Norway rat behavior on Rat Island, Alaska in 2003. 21 (2005).

72.
ESRI. ArcGIS 10.7. (ESRI, 2020).

73.
Simberloff D. Reconstructing the ambiguous: can island ecosystems be restored? in Conservation Sciences Publication (New Zealand). no. 2. (1990). More

Comparison of sampling methods for microbes on the surface of tobacco leaves
According to previous research, two sampling methods for microbes on the surface of tobaccos were selected to separately perform extraction and amplification of genome DNAs after sampling the microbes. As shown in Table 1, as for the first sampling method, the DNAs extracted from two tobacco leaf samples (fresh tobacco leaves and tobacco leaves in the later yellowing period) were both unqualified after subjected to amplification. Therefore, the first method was not suitable for extracting the microbes on the surface of tobacco leaves. The two tobacco samples extracted by using the second method both allowed favorably amplification and their amplification results were both proper. Thus, the second method was applied to sample the microbes on the surface of tobacco leaves subsequently.
Table 1 Comparison of effects of the two methods for sampling microbes on the surface of tobaccos.
Full size table

OTU clustering analysis
To explore the species compositions of various samples, OTUs clustering was carried out on effective Tags of all samples based on 97% of identity; afterwards, species annotation was performed on the OTUs sequences. According to OTUs results obtained through clustering and research requirements, the common and specific OTUs among different samples (groups) were analyzed.
OTU clustering analysis of bacteria in tobacco leaves
The result is shown in Fig. 3. Each petal in the petal diagram represents a group (sample) and different colors mean diverse samples (groups); the number at the core stands for the total number of OTUs in all samples; the number in each petal denotes the number of OTUs specific in the sample (group). It can be seen from the figure that the numbers of the core microbial communities subjected to conventional flue-curing procedure and dry-ball temperature set and wet-ball temperature degradation flue-curing procedure were basically consistent, showing no great change.
Figure 3

Petal diagrams of OTUs in samples flue-cured through conventional procedure and temperature- and humidity-controlled procedure under different sampling stages (SB: the surface bacteria of fresh tobacco leaves; EB: endophytic bacteria of fresh tobacco leaves; CSB: the surface bacteria of tobacco leaves flue-cured using conventional procedure; CEB: endophytic bacteria of tobacco leaves flue-cured using conventional procedure; SSB: the surface bacteria of tobacco leaves flue-cured using temperature- and humidity-controlled procedure; SEB: endophytic bacteria of tobacco leaves flue-cured using temperature- and humidity-controlled procedure; 2–6 represent different sampling stages).

Full size image

OTU clustering analysis of fungi in tobacco leaves
The result is displayed in Fig. 4. As shown in the figure, the core microbial communities flue-cured by conventional procedure and temperature- and humidity-controlled procedure showed basically coincident numbers. The latter was only 5–10 core microbial communities more than the former. Similar to the core bacterial communities, the number of core fungal communities presented no great difference in the flue-curing process.
Figure 4

Petal diagrams of OTUs of fungi in samples flue-cured by using conventional procedure and temperature- and humidity-controlled procedure in different sampling stages (SB: the surface fungi of fresh tobacco leaves; EB: endophytic fungi of fresh tobacco leaves; CSB: the surface fungi of tobacco leaves flue-cured using conventional procedure; CEB: endophytic fungi of tobacco leaves flue-cured using conventional procedure; SSB: the surface fungi of tobacco leaves flue-cured using temperature- and humidity-controlled procedure; SEB: endophytic fungi of tobacco leaves flue-cured using temperature and humidity-controlled procedure; 2–6 denote different sampling stages.

Full size image

Analysis of relative abundances of species
According to the results of species annotations, the species of each sample or group with the maximum abundance ranking the top 10–30 at various classification levels were selected to generate the cumulative histogram of relative abundances of species. The species with a high relative abundance in various samples at different classification levels and their proportions can be intuitively found.
Analysis of relative abundances of bacteria in tobaccos
Based on the results of species annotations, the species of each sample or group with the maximum abundance ranking the top 10–30 at various classification levels were selected to generate the cumulative histogram of relative abundances of species. The species with a high relative abundance in various samples at different classification levels and their proportions can be visualized.
As shown in Fig. 5, at the level of phylum, the bacteria in tobaccos mainly contained Proteobacteria, Actinobacteria, Bacteroidetes, Firmicutes, Planctomycetes, Acidobacteria, Chloroflexi, unidentified bacteria, Thaumarchaeota and Gemmatimonadetes. It can be seen that the surface and endophytic bacterial communities of fresh tobacco leaves slightly differed at the level of phylum. Proteobacteria showed the largest content, followed by Actinobacteria and Bacteroidetes, and the contents of the other bacterial phyla were relatively low. By using conventional flue-curing procedure, the bacterial diversity on the surface of tobacco leaves progressively declined as the flue-curing continued, and the relative content of Proteobacteria rose at first and then reduced; the reduction amplitude of Actinobacteria was relatively stable in the flue-curing process while that of Bacteroidetes was relatively large. For the dry-ball temperature set and wet-ball temperature degradation flue-curing procedure, as the flue-curing proceeded, the relative content of Proteobacteria gradually increased and it did not greatly reduce until reaching the last flue-curing stage. However, its relative content was not significantly different from that in fresh tobacco leaves; similar to Proteobacteria, the relative contents of both Actinobacteria and Bacteroidetes also grew at first and then decreased; in terms of endophytic bacteria of tobaccos, as the flue-curing process continued, the relative contents of Proteobacteria and the other main bacterial communities rapidly dropped while those of the other communities sharply increased.
Figure 5

Histogram of relative abundances of species at the level of phylum.

Full size image

As shown in Fig. 6, at the level of genus, the main dominant bacterial communities in endophytic bacteria of fresh tobacco leaves included Pseudomonas, Sphingomonas, Ralstonia, Methylobacterium, Massilia, Sphingobacterium, Rhizobium, Halomonas, Serratia and Rickettsia.
Figure 6

Column chart of species relative abundance at genus level.

Full size image

When employing conventional flue-curing procedure, the bacterial communities on the surface of tobaccos were relatively marginally changed at the level of genus while the endophytic bacteria varied remarkably. As the flue-curing process proceeded, the relative content of 30 main endophytic bacterial communities found before the flue-curing had dropped to 2% even in the early flue-curing stage (35 ℃). In comparison, the relative content of the bacterial communities in the first two flue-curing stages under dry-ball temperature set and wet-ball temperature degradation flue-curing procedure was higher.
Under conventional procedure, the relative abundance of Pseudomonas on the surface of tobacco leaves increased at first and then decreased, so did that of Sphingomonas. Although no signs of Ralstonia solanacearum were visualized on the surface of the sampled tobacco leaf samples, Ralstonia was found in the analysis of bacterial communities. With the ongoing flue-curing process, the relative content of Ralstonia rapidly reduced; the relative content of Methylobacterium on the surface of fresh tobacco leaves declined to some extent in the flue-curing process, and accounted for a large proportion in bacterial communities on the surface of flue-cured tobacco leaves. The relative contents of the other main bacterial communities were all progressively lowered basically.
When implementing dry-ball temperature set and wet-ball temperature degradation flue-curing procedure, the relative contents of the main bacterial communities in the early flue-curing stage were higher than those in fresh tobacco leaves. The relative contents of them marginally differed from those in fresh tobacco leaves even though flue-curing process was ended; the relative content of Pseudomonas gradually increased in the flue-curing process. By contrast, the relative contents of Sphingomonas and Methylobacterium both grew at first and then declined. The relative contents of the other main bacterial communities relatively slowly varied in the flue-curing process and they did not greatly decrease until the flue-curing process was ended. Moreover, the relative contents of some bacterial genera, including Sphingobacterium and Rickettsia, had remarkably dropped in the early flue-curing stage.
Analysis of relative abundances of fungi in tobaccos
As shown in Fig. 7, at the level of phylum, fungi in tobaccos mainly covered Ascomycota, Basidiomycota, Mortierellomycota, Rozellomycota, Glomeromycota, Chytridiomycota, Kickxellomycota, Mucoromycota and Olpidiomycota. It can be seen from the figure that fungi in tobaccos mainly included Ascomycota and Basidiomycota; the other fungi (phylum) took up a relatively low proportion.
Figure 7

Histogram of relative abundances of species at the level of phylum.

Full size image

When being flue-cured by using conventional procedure, the relative abundance of Ascomycota on the surface of tobacco leaves gradually increased while those of Basidiomycotaand the other fungal phyla gradually decreased with the ongoing flue-curing process; under temperature- and humidity-controlled flue-curing, the evolution law of fungal communities was similar to that using conventional procedure at the level of phylum. To be specific, a trend was shown that the relative abundance of Ascomycota gradually rose while those of Basidiomycota and the other fungal phyla were lowered, which was basically similar to that under conventional flue-curing procedure.
Their proportions were higher than those on the surface of tobacco leaves.The change amplitude of the endophytic fungi of tobacco leaves was less significant than that of fungi on the surface of tobacco leaves. Either under conventional flue-curing or temperature- and humidity-controlled flue-curing, the relative abundance of Ascomycota basically increased at first and then declined while that of Basidiomycota reduced at first, then grew and finally dropped.
As shown in Fig. 8, the change trends of community compositions under conventional flue-curing and temperature- and humidity-controlled flue-curing at the level of genus were similar to those at the level of phylum. There was a great difference only in the relative contents of fungal communities. In terms of fungi on the surface of tobacco leaves, the relative content of Alternaria under conventional flue-curing greatly increased at 38.5 ℃ and 54 ℃, with increases at the same time points under temperature- and humidity-controlled flue-curing; however, the growth amplitude of the relative content was less significant than that under conventional flue-curing. Cladosporium was another main fungal community and its relative content slowly decreased in the later stage of temperature- and humidity-controlled flue-curing. The relative content of Symmetrospora progressively decreased when using conventional flue-curing procedure while its reduction rate slowed down under temperature- and humidity-controlled flue-curing. The relative content of Ophiocordyceps on the surface of tobacco leaves was relatively low and it both gradually reduced when using the two flue-curing technologies. Moreover, the relative contents of the other fungal genera also progressively declined as the flue-curing proceeded.
Figure 8

Histogram of relative abundances of species at the level of genus.

Full size image

As for changes of endophytic fungi of tobaccos when using the two flue-curing technologies, the relative content of Alternaria under conventional flue-curing was higher than that under temperature- and humidity-controlled flue-curing. The result can be found even though the flue-curing process was ended. The relative content of Cladosporium marginally varied under conventional flue-curing and greatly increased at 35 ℃. After completing the flue-curing process, the relative content did not significantly differ from the value in fresh tobacco leaves. However, for dry-ball temperature set and wet-ball temperature degradation flue-curing procedure, the relative content of Cladosporium progressively reduced on the whole and the value after ending the flue-curing process was only about half of that in fresh tobacco leaves. The relative content of Symmetrospora showed a same change trend with Cladosporium under conventional flue-curing and temperature- and humidity-controlled flue-curing. Additionally, the reduction rate of the relative content of Symmetrospora was higher than that of Cladosporium under temperature- and humidity-controlled flue-curing. Relative to fungal communities on the surface of tobaccos, although the relative contents of the other endophytic fungal genera gradually dropped with the flue-curing.
Clustered heat maps of species abundances
According to species annotations and abundances of all samples at the level of genus, genera whose abundances ranked the top 35 were selected. Subsequently, based on the abundances of these genera in each sample, a heat map is drawn by conducting clustering from the two aspects: i.e. species and samples. By doing so, it is convenient to ascertain a species with a high abundance or low content and the sample from which it is found.
Clustered heat map of species abundances of bacteria in tobaccos
The result is displayed in Fig. 9. It can be seen from the clustered heat map that bacterial communities mainly reside on the surface and in the interior of the fresh tobacco leaves at first. As the flue-curing proceeded, the contents of the main bacterial communities were changed to some extent: the main bacterial communities reduced in the content and the tobacco leaves became light in color.
Figure 9

Clustered heat map of species abundances at the level of genus.

Full size image

Clustered heat map of species abundance of fungi in tobaccos
According to species annotations and abundances of all samples at the level of genus, the genera whose abundances ranked the top 35 were selected. Based on the abundances of these genera in each sample, a heat map is drawn by conducting clustering from the two aspects, i.e. species and samples. It can be found from Fig. 10 that the distribution of fungi on the surface of tobaccos greatly differed from that of endophytic fungi. Moreover, the distribution of materials also presented a great difference when using the two flue-curing technologies.
Figure 10

Clustered heat map of species abundances at the level of genus.

Full size image

Correlation with environmental factors
Correlation of bacteria with environmental factors
Temperature and humidity were mainly controlled in the flue-curing process of tobacco leaves; and the main environmental factors involved dry- and wet-bulb temperatures. It can be seen from Fig. 11 that Pantoea, Nesterenkonia, Staphylococcus, Variovorax, Chryseomonas, Rhodococcus, Paracoccus, Massilia, Serratia, Ralstonia and Pseudomonas were more likely to be affected by temperature and humidity. Pantoea and Variovorax exhibited a positive correlation with temperature and humidity; Nesterenkonia, Staphylococcus, Chryseomonas, Rhodococcus, Paracoccus, Serratia and Ralstonia presented a negative correlation with temperature and humidity. Actinomycetospora were negatively correlated with the dry-bulb temperature while positively correlated with the wet-bulb temperature; Stenotrophomonas, Cutibacterium and Sediminibacterium were all negatively correlated with both dry- and wet-bulb temperatures while they presented a higher negative correlation with the dry-bulb temperature.
Figure 11

Heat map of correlation with environmental factors at the level of genus.

Full size image

Correlation of fungi with environmental factors
Similar to the analysis method of environmental factors of bacteria, the correlation between the fungi in tobaccos and environmental factors is displayed in Fig. 12. Temperature and humidity were mainly controlled in the flue-curing process of tobacco leaves; the main environmental factors were dry- and wet-bulb temperatures. As shown in the figure, different from the correlation of bacterial communities in tobaccos with environmental factors, the majority of fungal genera in tobaccos presented a negative correlation with temperature and humidity, for example, Rachicladosporium, Vishniacozyma, Symmetrospora, Sarocladium, Ascochyta, Wallemia, Colletotrichum, Fusarium, Claviceps, Cladosporium, etc. A small number of fungal genera (such as Ustilaginoidea, Septoriella and Alternaria) were positively correlated with temperature and humidity.
Figure 12

Heat map of correlation with environmental factors at the level of genus.

Full size image

Function prediction of bacterial communities in tobaccos
According to the annotation result in the database, the functions with the maximum abundance ranking the top 10 in each sample or group at various layers of annotations were selected to generate the cumulative histogram of relative abundances of functions. Thus, it is convenient to check the functions with a high relative abundance in various samples at different layers of annotations and their proportions.
As shown in Fig. 13, bacterial communities with a half of relative abundances participated in the metabolism and a quarter of bacterial communities took part in the genetic information processing; the rest was engaged in the cellular processes and organismal systems and some bacterial communities were implicated to human diseases.
Figure 13

Relative abundances of function annotations of bacterial communities in tobaccos.

Full size image

Function prediction of fungal communities in tobaccos
Based on amplificon analysis of 16S or ITS, the species classification and abundances of fungi present in the environment can be attained. In many cases, people will also concern what role these species found in the environment play in the ecological environment. By applying FunGuild tool, it is feasible to attain the ecological functions of corresponding fungi based on the species classification of fungi. It can be seen from the Fig. 14 that the fungal communities in tobaccos delivered relatively abundant functions, in which the three fungal communities with the highest relative abundances separately showed the following functions: plant saprophytes, plant pathogens and undefined functional fungal communities; the rest of fungal communities participated in plant parasitism,soil-borne plant pathogen, lichenization, dung saprotroph, etc.
Figure 14

Relative abundances of function annotations of fungal communities in tobaccos.

Full size image More