Skin bacterial microbiome of a generalist Puerto Rican frog varies along elevation and land use gradients

Introduction

All animals host a diverse community of microbes. These microbes, along with their genetic contributions, are termed the microbiome, and they provide many important processes for their hosts (McFall-Ngai et al., 2013). For example, the vertebrate gut microbiome helps regulate development of the immune system, as well as playing critical roles in digestion and preventing pathogen colonization (Shreiner, Kao & Young, 2015). Transmission of these microbial symbionts can occur in several ways. In some systems, highly specialized symbionts are vertically transmitted from mother to offspring (Funkhouser & Bordenstein, 2013), while in other systems, many symbionts appear to be obtained and maintained from environmental sources (Bright & Bulgheresi, 2010). For species that obtain symbionts from environmental source pools, variation in environmental conditions that impact those source pools could alter the composition, and ultimately the function, of these microbial communities.

Large-scale differences in environmental conditions contribute to variation in local microbial communities (Fierer & Jackson, 2006), including those that might colonize hosts. For instance, conversion from natural to agricultural land use alters a variety of soil properties, such as temperature and moisture, which in turn impacts the diversity and composition of microorganisms present in soil (Matson et al., 1997). In addition, elevational gradients in soil pH appear to play role in determining the composition of soil bacterial assemblages on mountains (Shen et al., 2013). Despite the known impact of abiotic factors on microbial source pools, relatively little work has focused on understanding variation in the microbiome of wildlife across environmental gradients. Studies of non-human primates and frogs have found changes in the gut microbiome associated with different land use types (e.g., natural versus human-modified habitats), seemingly due to differences in the diet of animals in those habitats (Amato et al., 2013; Chang et al., 2016). In both systems, the gut microbiomes of hosts inhabiting the anthropogenically-modified habitats had characteristics that could negatively impact host health (Amato et al., 2013; Barelli et al., 2015; Gomez et al., 2015; Chang et al., 2016).

The skin microbiome of amphibians is receiving substantial research attention, in part because of its potential role as a barrier against pathogen infection (Woodhams et al., 2016). Variation in amphibian skin microbiome diversity among individuals, populations, and species is well-documented (McKenzie et al., 2012; Kueneman et al., 2014; Walke et al., 2014; Belden et al., 2015; Rebollar et al., 2016). However, the underlying causes and potential consequences of this variation are only just beginning to be understood. Importantly, large-scale factors, such as land use (Krynak, Burke & Benard, 2016) and elevation (Bresciano et al., 2015; Muletz-Wolz et al., 2017), can play a role in determining amphibian skin microbiome diversity. This is not altogether surprising, given that at least a portion of the amphibian skin microbiome is likely obtained from their surroundings and that environmental source pools of microbes appear to be critical for the maintenance of these skin communities (Fitzpatrick & Allison, 2014; Loudon et al., 2014; Michaels, Antwis & Preziosi, 2014; Walke et al., 2014; Rebollar et al., 2016). In addition, skin microbial diversity along these environmental gradients may be influenced by interactions with the chytrid fungus, Batrachochytrium dendrobatidis (Bd; Longcore, Pessier & Nichols, 1999), a widespread amphibian skin pathogen. Bd infection prevalence and intensity vary with elevation (Burrowes, Longo & Rodríguez, 2008; Burrowes et al., 2008; Sapsford, Alford & Schwarzkopf, 2013) and land use (Becker & Zamudio, 2011), and several recent studies have demonstrated a relationship between the presence of Bd and skin bacterial diversity (Jani & Briggs, 2014; Bresciano et al., 2015; Longo et al., 2015; Rebollar et al., 2016; Longo & Zamudio, 2017).

Our aim in this study was to assess variation in the skin microbiome of Coqui frogs, Eleutherodactylus coqui, along land use and elevation gradients. Specifically, we investigated elevational changes in the bacterial communities of the skin along two land use gradients (intact forest transect vs. disturbed forest transect) on the northern slopes of the Sierra de Luquillo in Puerto Rico. Eleutherodactylus coqui is a well-studied, terrestrial, direct-developing frog native to Puerto Rico (Joglar, 1998; Longo & Burrowes, 2010). It is an excellent model species for this study because it is a habitat generalist that can be found in abundance across a range of habitat types and elevations. In addition, E. coqui is vulnerable to Bd, but has persisted in spite of marked declines in the 1990s at high elevation sites in Puerto Rico (Burrowes, Joglar & Green, 2004).

Methods and Materials

Results

In total, we identified 1,087 OTUs from the skin of E. coqui (N = 72 individuals), with a range of 143–818 of those present on individual frogs (mean ± SD = 470 ± 153 OTUs/frog). Most OTUs were from the phyla Proteobacteria (55%), Actinobacteria (17%), Bacteroidetes (13%), and Firmicutes (7%).

There was no indication that any of the individuals we sampled were infected with Bd. Additionally, there were no significant differences between males and females, or between gravid and non-gravid females, in alpha or beta diversity metrics (Tables 2 and 3). Therefore, we included all individuals in subsequent analyses focused on site, elevation, and land use effects.

Alpha diversity metrics varied among sites (Simpson index, Shannon index, phylogenetic diversity: all P < 0.02), but not based on elevation or land use (Table 2, Figs. 1A–1C). Although alpha diversity was variable at all sites, one site had particularly low alpha diversity relative to all other sites: Boca del Yunque (low elevation, intact forest; Fig. 1A).

Individuals displayed community structures ranging from relatively similar to extremely different both within and across sites. Indeed, the range for within site variation was comparable to that for between site variation (range, from the same site, Jaccard: 0.32–0.86, Bray–Curtis: 0.13–0.95, Unifrac: 0.10–0.97; from different sites, Jaccard: 0.33–0.90, Bray–Curtis: 0.18–0.99, Unifrac: 0.13–0.97). Still, the skin bacterial communities of individuals from the same site tended to be more similar in comparison to those from different ponds, at least in terms of OTU presence/absence and relative abundance (Jaccard and Bray–Curtis: P < 0.0001 and P = 0.03, respectively; Unifrac: P = 0.23).

Beta diversity varied across sites, elevations, and land use (PERMANOVA all P < 0.03; Table 3, Figs. 1D–1F). Elevation explained 7–12% of the variation in bacterial community composition; land use explained less than 5% (all indices; Table 3). This may be in part due to relatively subtle differences in OTU relative abundance between land use types, and stronger differences across elevations. For example, 20 OTUs were identified as differing between land use types based on K–S statistics, but K–S measures associated with these OTUs were all relatively low (range 0.25–0.29; Table 4), indicating that the differences in their empirical distribution functions across land use types were quite small. Most of these K–S OTUs that differentiated the intact and disturbed sites had a higher relative abundance in intact, forested habitats (Fig. 2A). By contrast, for elevational differences, K–S measures for the 18 OTUs identified ranged from 0.47 to 0.58 (Table 5). Most of these OTUs displayed the highest relative abundances at mid-elevations (Fig. 2B).

The core microbiome was comprised of 34 OTUs (Table 6). The vast majority (82%) of core OTUs were Proteobacteria; however, at lower taxonomic levels, a variety of bacterial orders (N = 11) and families (N = 14) were represented. Relative abundances of core OTUs varied considerably across individuals, but were generally low (mean relative abundance of 25/34 core OTUs was ≤0.01%; Table 6, Fig. 3). There were only a few OTUs that displayed higher relative abundances on most individuals, most notably X394796 Pseudomonas viridiflava and X4378239 Sanguibacter sp. (relative abundance, mean ± sd, X394796: 10.6 ± 11.2% and X4378239: 7.2 ± 8.2%; Fig. 3). No core OTUs were identified by K–S measure as varying across elevation or land use.

Discussion

Land use change is widely recognized as a major factor impacting biodiversity around the globe. There is increasing evidence linking land use practices to changes in microbiome diversity, especially gut microbiome diversity (Amato et al., 2013; Barelli et al., 2015; Amato et al., 2016; Chang et al., 2016), whereas patterns of microbiome diversity and function along elevation gradients have received less attention (but see Bresciano et al., 2015; Medina et al., 2017; Muletz-Wolz et al., 2017). Our findings reveal subtle, but potentially important, variation in amphibian skin bacterial community structure with site, elevation, and land use.

There are several mechanisms by which these factors could impact skin microbiome diversity. Changes in E. coqui bacterial community structure may reflect direct responses of the skin communities to external environmental conditions experienced by the host, or variation in the environmental microbes available to colonize the skin. Studies comparing the diversity of the amphibian skin bacterial communities to that of free-living bacteria in the surrounding environment suggest that many bacteria that reside on amphibian skin can be also be found in the environment, although relative abundances often differ in environmental samples and amphibian skin samples (Kueneman et al., 2014; Walke et al., 2014; Rebollar et al., 2016; Prado-Irwin et al., 2017). This, along with experimental manipulations demonstrating the importance of environmental source pools to the maintenance of skin bacterial diversity (Loudon et al., 2014; Michaels, Antwis & Preziosi, 2014), provide growing evidence for dispersal of bacteria between amphibians and their habitat. Although assessing the environmental microbiota was beyond the scope of this study, it seems reasonable that variation in the diversity of the environmental microbiota contributes to the variation in E. coqui skin bacterial communities observed in this study.

Additionally, the differences we observed may be a consequence of altered host physiology or immune function under different environmental conditions. For example, in Neotropical amphibians, elevation (Burrowes, Longo & Rodríguez, 2008; Burrowes et al., 2008; Sapsford, Alford & Schwarzkopf, 2013) and land use (Becker & Zamudio, 2011) are important predictors of the pathogenicity of the chytrid fungus, Batrachochytrium dendrobatidis. Bresciano et al. (2015) observed both increased Bd infection prevalence and increased presence of Bd-inhibiting skin bacteria at high elevations relative to low elevations, linking elevational patterns of pathogen presence in the host with changes in skin bacterial diversity and disease-resistance function. This suggests that there may be stronger selection for more Bd-inhibiting bacteria where Bd infection prevalence is greater. On the other hand, physiological differences associated with harboring even mild Bd infections, as occurs in high (>650 m) but not low elevation populations of E. coqui (Burrowes, Joglar & Green, 2004; Burrowes, Longo & Rodríguez, 2008), could manifest as changes to the microbiome. Experimental tests of changes in microbiome diversity point to a role for both direct and host-mediated responses across a wide variety of organisms (Lokmer & Wegner, 2015; Kohl & Yahn, 2016; Zaneveld et al., 2016). Thus, variation in the amphibian skin microbiome likely represents a complex interplay of the effects of both abiotic and biotic conditions on the host and the symbiont community.

One previous study of the microbiome of E. coqui (Longo et al., 2015) provides some additional insight into important drivers of bacterial community structure for this species. That study assessed variation in bacterial community diversity across juveniles and adults infected and not infected with Bd, and determined that developmental stage was a much more important determinant of alpha and beta diversity than Bd infection status (Longo et al., 2015). Based on their findings, Longo et al. (2015) suggested that the environment (what they term “microbial reservoirs”) likely plays a stronger role in shaping the skin microbiota than disease-related factors, such as host immunity. In the present study, no individuals were infected with Bd, so we were unable to assess the effects of immediate disease state on microbiome diversity. The absence of Bd-infected individuals is most likely due to the timing of the study, which was conducted in the warm-wet season. Previous work suggests that both Bd prevalence and infection intensities are lower during this time of the year relative to the cool-dry season (Longo & Burrowes, 2010; Longo et al., 2013; Longo & Zamudio, 2017). However, we expect that Bd could be related to microbiome diversity over longer time scales. Historically, E. coqui experienced die-offs from Bd at higher elevations (Burrowes, Joglar & Green, 2004; Burrowes, Longo & Rodríguez, 2008; Burrowes et al., 2008; Longo et al., 2013), and persisting populations are more vulnerable to the fungus in montane habitats during the dry season (Burrowes et al., 2008; Longo et al., 2013). If these die-off events imposed selection on microbiome diversity, then legacy effects may explain some of the elevational variation in microbiome diversity we observed. Similarly, Bd may contribute to among-site variation in the microbiome of the robber frog Craugastor fitzingeri from Bd endemic versus Bd naïve sites in Panama (Rebollar et al., 2016).

Despite the influence of these varied factors on microbial community structure, E. coqui appear to maintain a core set of bacterial taxa on their skin. Our characterization of E. coqui core taxa is broadly similar to the findings of Longo et al. (2015), who identified five core phylotypes: Actinomycetales, Stenotrophomonas sp., Comamonadaceae, Staphylococcus sp., and Pseudomonadaceae. All of these taxa were represented in our assessment of the core microbiome as well, in addition to many other OTUs. Our understanding of the functional role of the core microbiota of amphibians is limited, but the presence of a core microbiota—which exists in many amphibian species (e.g., Becker et al., 2014; Loudon et al., 2014; Walke et al., 2014)—may imply functional importance. Studies examining the functional role of amphibian skin bacteria have primarily focused on the ability to protect the host against pathogens (e.g., Bletz et al., 2013; Woodhams et al., 2016; Longo & Zamudio, 2017). In this study, many of the core bacteria of E. coqui were in genera containing members with known anti-Bd properties, including Acinetobacter, Janthinobacterium, Pseudomonas, and Stenotrophomonas (Woodhams et al., 2015). However, the skin bacteria may serve other functions. For example, in fish, some skin bacteria may help in the production of the mucosal layer of the skin (Sar & Rosenberg, 1989); this could also be an important function for amphibian skin, especially for terrestrial amphibians in the relatively low humidity conditions typical of disturbed habitat. Additional studies are needed to better understand what role the core and variable components of the microbiota of E. coqui may play in allowing them to persist in a variety of habitats and in presence of Bd. Moreover, it would be interesting to know if the core microbiota are obtained via different mechanisms than the more variable components. For amphibians like E. coqui that provide parental care (Townsend, Stewart & Pough, 1984), there is growing evidence that members of the skin microbiota are shared between parents and eggs (Banning et al., 2008; Walke et al., 2011; Hughey, Delia & Belden, 2017), which could indicate vertical transfer of microbial symbionts.