Family shapes microbiome differences in Oklahoma salamanders

Introduction: Given the role of microbiomes in promoting host health and homeostasis, understanding the factors shaping skin microbial communities in wild vertebrates has become increasingly important in conservation. This goal is even more pressing for amphibians, for which the skin has multiple critical functions, and pathogens currently decimating populations are linked to significant changes in skin microbiomes. However, because microbiomes are also shaped by environmental and ecological influences, as well as by host phylogeny, it is important to quantify these contributions to microbiome structure in the presence of infection.

Methods: To understand the joint influence of these diverse factors shaping microbiomes, we used 16S rRNA sequencing to characterize the skin microbial communities of six salamander species (families Plethodontidae and Salamandridae) found in Oklahoma and contrasted the effects of infection status, phylogeny, host ecology, and host environment (i.e., climate) on skin microbiomes.

Results: Differences at the level of host family were the main factor influencing microbiome diversity; however, we did not detect a substantial phylogenetic signal. Instead, host ecology and environment were more important in driving microbiome differences among species and genera. Salamanders that tested positive for the skin fungal parasite Batrachochytrium dendrobatidis (Bd) also had slightly less diverse microbiomes than Bd-free animals, but no such differences were associated with the systemic pathogen ranavirus (RV).

Discussion: Together, these results indicate a nuanced relationship between the number and type of microbes present on salamander skin and the factors influencing them. By developing a baseline assessment of the microbiome diversity and richness present on the skin of these focal species, this work also provides a foundation for monitoring and evaluating changes in skin microbiomes as populations continue to experience stressors and diseases.

1 Introduction

Microbiomes, the microorganisms or their genetic material that exist in a particular environment or on a host, are increasingly recognized as providing crucial functions relating to host health, aiding in processes such as digestion and protection from disease (Moloney et al., 2014). There is increasing recognition that even small changes in relative abundance within microbial communities can result in the loss of host innate functions (metabolic, physiological, and immune) – a condition referred to as dysbiosis–with profound negative effects on a host’s overall fitness (e.g., Chimenos-Küstner et al., 2017; Jiménez and Sommer, 2017). Dysbiosis can lead to behavioral and developmental changes as well as decreased immune function. For example, even a slight shift in the abundance or community structure of microbiota can be linked to illness across a wide range of taxa (e.g., humans, amphibians, insects; Hooper and Gordon, 2001; Brucker and Bordenstein, 2012; Kueneman et al., 2016; Chimenos-Küstner et al., 2017; Bletz et al., 2018). Thus, developing a baseline understanding of host-associated microbiomes and identifying the role of microbiomes in shaping susceptibility to infection is critical for many species, especially those of conservation concern (Jiménez and Sommer, 2017).

Amphibians are currently the most threatened vertebrates on the planet, with 41% of amphibian species threatened with extinction due to climate change, habitat loss, and infectious disease (IUCN, 2022). At the end of the twentieth century, two infectious diseases began to be linked to widespread declines in amphibian populations: chytridiomycosis (caused by Batrachochytrium dendrobatidis [Bd] or B. salamandrivorans [Bsal]) and ranavirus disease (caused by infection from one of the several viruses in the genus Ranavirus [family Iridoviridae]; Daszak et al., 1999; Lips et al., 2006; Fisher and Garner, 2020; Scheele et al., 2017; Skerratt et al., 2007). Amphibian susceptibility to Bd infection is influenced by their extremely permeable skin, which allows for cutaneous respiration and osmotic regulation (Campbell et al., 2012; Rollins-Smith, 2009; Smith et al., 2018). Vast microbial diversity is present on the skin of amphibian species, and certain microbial taxa can influence host immunity directly by producing antifungal or antimicrobial compounds (Bletz et al., 2017b; Woodhams et al., 2015). This has been observed with the presence of anti-Bd microbes on the skin of amphibians linking to reduced fungal infection, underscoring the important role beneficial bacteria may play in pathogen susceptibility among hosts (Rebollar et al., 2020; Mutnale et al., 2021; Nava-González et al., 2021). Amphibians vary in their susceptibility to pathogens (DiRenzo et al., 2021); this may be the result of both the host’s immune system and the varied microbial diversity found in the skin microbiome. Thus, it is important to characterize and compare microbial diversity across a range of species.

Complicating our understanding of the role amphibian microbiomes play in inhibiting infection is the fact that infection can change the composition of the skin microbiome (Wilber et al., 2020). Amphibians experience a significant reduction in the skin microbiome's diversity after recovery from chytridiomycosis (Jani and Briggs, 2014; Rebollar et al., 2016; Medina et al., 2017; Bates et al., 2019; Muletz-Wolz et al., 2019; Ruthsatz et al., 2020; Jani et al., 2021). Likewise, the skin microbiomes of amphibians that have experienced ranavirus infection exhibit a shift in composition compared to virus-free individuals or populations (Campbell et al., 2019; Harrison et al., 2019). Interestingly, the skin microbiomes of some amphibian taxa do not appear to be strongly impacted by fungal infection, such as that of the red-backed salamander (Plethodon cinereus; Becker and Harris, 2010; Barnes et al., 2021; Bates et al., 2022). Given the variation among these findings, additional studies documenting the diversity of microbes found in skin microbiomes and their relationship to infection are needed.

Fully understanding the relationship between amphibian infectious disease and host microbiomes requires first determining the baseline composition of the microbiome and the factors driving its diversity. Three main potential contributors to skin microbiomes have been identified—environment, ecology, and phylogeny—each with a unique set of predictions. If the physical environment (i.e., conditions related to local climate and/or biotic community) influences which microbes are able to colonize and thrive on a host, then amphibian species living in the same locality or ecoregion should show similar skin microbial composition. Past studies have found that abiotic aspects of the environment such as temperature, precipitation, and/or proxies for environmental similarity such as geographic distance (Ruthsatz et al., 2020; Walker et al., 2020) are significantly correlated with microbiome diversity. However, if host ecology (e.g., life stage or how a host interacts with the microclimates in its environment) is a greater predictor of skin microbial diversity than the animal’s surrounding environment, species’ skin microbiomes would differ based on a species’ ecological niche, developmental stage, and/or microhabitat preferences. Groups of amphibian species found in similar microhabitats (e.g., arboreal or aquatic) would exhibit similar skin microbiomes, regardless of geographical location. This could be due to lateral transfer among hosts or exposure to the same microhabitats or stressors (Bletz et al., 2017a). Thus, differences in specific microhabitat use, even within the same geographic region, could be significant enough to result in distinct skin microbiomes when compared to species in different niches. There is also evidence for host ecology driving skin microbial structure with differences in life history (Bletz et al., 2017a) and habitat type (Bird et al., 2018; Xu et al., 2020; Smith et al., 2023) shown to correlate significantly with skin microbial diversity.

In contrast to findings that implicate environment or host ecology as driving skin microbial diversity, multiple studies have found that phylogeny–i.e., differences related to either the evolutionary history or, as a proxy, taxonomic grouping reflective of that history–was the strongest predictor of diversity in amphibian skin microbiomes (McKenzie et al., 2012; Prado-Irwin et al., 2017). Certain microbes may be conserved as part of the skin microbiome due to evolutionary differences between amphibian species, and this could dictate which microbes become part of the microbiome. Among the few studies that have investigated skin microbial makeup in amphibians, sampling does not often involve a wide range of taxa, and where it does, phylogeny isn’t often explicitly considered (e.g., Ellison et al., 2019), or is considered by contrasting hosts at different taxonomic levels as a proxy for phylogenetic divergence (e.g., Buttimer et al., 2022). In a study on the salamander genera Batrachoseps, Ensatina, and Taricha, host genus was the main predictor of skin microbial richness (Buttimer et al., 2022). In contrast, a similar study focused only on Batrachoseps and Ensatina found that habitat was the most important factor shaping diversity (Bird et al., 2018). Thus, the taxonomic and phylogenetic sampling of a study, as well as the analytical approaches employed, may influence the inferred importance of factors shaping microbiome diversity. Occasions in which phylogenetic relationships appear to drive skin microbial diversity could also be an example of phylosymbiosis, in which the difference in skin microbiome makeup of two species may correspond to the evolutionary history and phylogenetic distance between the hosts (Lim and Bordenstein, 2020). Previous results are suggestive of the possible presence of phylosymbiotic relationships in amphibians (e.g., Kueneman et al., 2014; Ellison et al., 2019; Buttimer et al., 2022).

Ultimately, multiple drivers are likely to impact the diversity of skin microbiomes simultaneously. For example, a study considering multiple interacting factors in four frog and one newt taxa found that species identity correlated most strongly with microbial community composition differences, but habitat location was a secondary result of significant variation within each species (Kueneman et al., 2014). Similarly, a study of 49 frog (genus Plectrohyla) and 23 salamander (genera Bolitoglossa and Pseudoeurycea) species concluded that family, genus, or species identity (i.e., some aspect of the species’ ancestry; “phylogeny”) was the most important indicator of skin microbial diversity between different orders and families; however, between genera and species, microbial diversity was most impacted by host habitat (Ellison et al., 2019). Despite numerous studies, a consensus on the drivers of amphibian skin microbial diversity has not been reached. Thus, examining multiple interacting drivers simultaneously is necessary to understand how microbiomes assemble and how they intersect with amphibian disease.

In this study, skin microbiome data were collected from six focal, co-distributed species of salamanders representing two families and three genera across four level-III ecoregions (i.e., areas of distinct geography with co-distributed biota) in Oklahoma (Figure 1; Supplementary Table 1). As one of only four states in the country with 12 or more level-III ecoregions (Hoagland and Stoodley, 1998), Oklahoma is particularly well suited for comparisons among differing environments. In terms of host ecology, we sampled two direct-developing species (Plethodon albagula and P. angusticlavius; family Plethodontidae), three species with aquatic larvae and terrestrial or aquatic adults (Eurycea lucifuga, E. longicauda, and E. tynerensis; family Plethodontidae), and one species with an aquatic larval phase, terrestrial juvenile phase, and aquatic adult phase (Notophthalmus viridescens; family Salamandridae). Many of these samples were previously examined for pathogens (Bd and ranavirus), or we evaluated pathogen presence in the context of this study. Our goals were to (1) characterize the skin microbiomes of these six salamander species; (2) evaluate the role of the host’s environment (evaluated by both ecoregional differentiation and directly using climatic variables), host ecology, and evolutionary history on skin microbiome composition and similarity; and (3) determine the relationship between pathogen presence and microbial abundance and diversity. The results provide baseline knowledge about skin microbial composition in salamander species across regions of the United States for which little is known. Additionally, the findings add to a growing body of literature focused on the intersection of infectious disease with other factors structuring the skin microbiomes of vertebrates.

Figure 1

Figure 1 Sampling locations and abundance of six salamander host species in this study. Sampling localities (points) occur within the four level III ecoregions (Arkansas Valley, Boston Mountains, Ouachita Mountains, and Ozark Highlands) in eastern Oklahoma, USA. Bar charts refer to the proportion of salamander species found within each ecoregion.

2 Methods

2.1 Sample populations

Here, we used epithelial swab samples and pathogen data, some with associated museum-vouchered specimens, collected by the Herpetology Department of the Sam Noble Oklahoma Museum of Natural History (SNOMNH) from 2015–2021 (Marhanka et al., 2017; Watters et al., 2018; Davis et al., 2019; Watters et al., 2019). We identified 287 samples (Supplementary Data) representing six species of salamanders with sufficient sampling for robust characterization of variation in species-specific skin microbiomes: Eurycea longicauda, E. lucifuga, E. tynerensis, Plethodon albagula and P. angusticlavius, and Notophthalmus viridescens. Samples were selected from four unique Level III ecoregions across eastern Oklahoma: Arkansas Valley, Boston Mountains, Ouachita Mountains, and Ozark Highlands (Figure 1; Hoagland and Stoodley, 1998). Where available, pathogen status (i.e., presence/absence of Bd and ranavirus) was ascertained from prior amphibian surveys containing published pathogen data (Marhanka et al., 2017; Watters et al., 2018; Davis et al., 2019; Watters et al., 2019). Samples that had not been previously assayed for pathogens were analyzed as described below.

2.2 Sample collection

During each collection effort (comprising 24–72 hours each), salamanders were captured in ponds, wetlands, and streams by hand, dipnet, aquatic trap, or seine. Animals were kept individually in sterile plastic bags before being swabbed and either released or euthanized and vouchered in the Herpetology Collection at SNOMNH. Microbiome samples were taken by rubbing a cotton swab tip (Puritan Medical Products, Orlinda, Tennessee, United States) along the ventral and dorsal portions of the trunk, limbs, and toes five times each to collect microbial DNA, and the animal was swabbed in the same way to test for the presence of Bd, using rayon-tipped sterile swabs (Peel Pouch Dryswab Fine Tip [MWE 113], Corsham, Wiltshire, UK), following the methods of Lannoo et al. (2011). Microbial swabs were placed into empty vials and flash-frozen using liquid nitrogen prior to long-term storage at –20°C. Swab samples for Bd were placed in a dry ice cooler prior to long-term storage at –20°C. Animals were sampled for ranavirus by either i) clipping off a portion of the tip of the tail, or ii) euthanization by chlorobutanol or MS-222 solution submersion and then dissection to obtain a tissue sample from the liver, following St-Amour and Lesbarrères (2007) and Gray et al. (2012). Tissue samples were then placed in long-term storage at 20°C. All sampling was conducted in accordance with University of Oklahoma IACUC Protocols R14-026 (2015–2017), R17-031 (2017–2020), and R21-005 (2020–2023).

2.3 Life history stage and host ecology

The life history stage was assessed during sampling or determined from examination of vouchered specimens at the SNOMNH. Individuals were classified as adults if the snout–vent length (SVL) was greater than or equal to 31 mm for E. longicauda, 26 mm for E. tynerensis, 55 mm for P. albagula, and 30 mm for P. angusticlavius (Trauth et al., 2004). Eurycea lucifuga were considered adults if the total length (TL) of the body was 100 mm or greater (Petranka and Scott, 1998). Individuals from either species of Plethodon were classified as juveniles if below the adult threshold SVL values listed above. Members of Eurycea were considered larval if they were below the SVL or TL adult thresholds and had gills; they were classified as juveniles if they had no external gills but were below the length requirement. Adults belonging to E. tynerensis were determined to be paedomorphic if external gills could be seen and they were above the length threshold. All vouchered newts in the SNOMNH Herpetology Collection were adults; they were distinct from larval or eft forms, based upon comparison to published images and descriptions (Trauth et al., 2004).

Samples were labeled as coming from terrestrial or aquatic individuals based on the age and physical characteristics of the animal. If gills were present in the Eurycea, samples were listed as coming from aquatic individuals. Since all newts were determined to be adults, they were listed as aquatic (Trauth et al., 2004). All other animals were categorized as terrestrial.

2.4 Climate and elevation data

Latitude and longitude were recorded for each sample at the time of sample collection using WGS 84. Using these coordinates, elevation was retrieved using an online map tool (The National Map; U. S. Geological Survey, 2017). Abiotic environmental conditions on the local climate—average annual, seasonal, and monthly precipitation and temperature data—associated with sample coordinates and month/year of sampling were assessed for each locality using ClimateNA (Wang et al., 2016). Before use in analyses of alpha and beta diversity, environmental data were centered, with a mean of zero, and scaled, to normalize standard deviations, in R software version 4.2.1 (R Core Team, 2022).

2.5 DNA Extraction and Sequencing

Genomic DNA from epithelial samples was extracted in the Sam Noble Museum Genomics Core Facility using ZymoBIOMICS DNA Miniprep Kits (Zymo Research Products, Irvine, CA, United States). Extraction negatives were run with each extraction, and PCR negative and positive controls were run on each PCR plate. Positive controls consisted of the ZymoBIOMICS Microbial Community Standard and the ZymoBIOMICS Microbial Community DNA Standard (Zymo Research Products, Irvine, CA, United States, RRID : SCR_008968). The V4 subregion of the 16S rRNA gene was amplified using a one-step PCR with the barcoded primers and protocols described in Kozich et al. (2013). We visualized 2 µL of several samples per plate via gel electrophoresis to ensure amplification, after which all samples were cleaned using KAPA Pure Beads (Roche Sequencing Solutions, Pleasanton, CA, United States). DNA concentrations were determined using a Qubit Fluorometer (Thermo Fisher Scientific, Waltham, MA, United States), and samples were normalized to 10 nM of DNA before pooling all samples in a 1.5 mL sterile microcentrifuge tube. The pooled library was then sequenced through a single run on an Illumina MiSeq platform (RRID : SCR_016379) using 2 × 250 bp paired-end sequencing at the University of Oklahoma Consolidated Core Lab.

Pathogen presence data were sourced from previous surveys (Marhanka et al., 2017; Watters et al., 2018; Davis et al., 2019). Recent samples that had not yet been evaluated for pathogens were processed following the methods of Watters et al. (2018). In brief, DNA was extracted using PrepMan Ultra (Bd assays; Life Technologies, Carlsbad, CA, United States; Cheng et al., 2011) and/or a method of high-salt DNA extraction (ranavirus assays; Esselstyn et al., 2008). Pathogens were amplified following previously published qualitative PCR procedures (qPCR; Kerby et al., 2013). Primers targeted the internal transcribed spacer (ITS-1) ribosomal RNA gene (forward primer: ITS1-3 Chytr; reverse primer: 5.8S) for Bd and the major capsid protein (MCP) for ranavirus (Boyle et al., 2004; Forson and Storfer, 2006). All samples were run in triplicate, and any samples with positive values <1.0 gene copies and/or only one well testing positive were run a second time.

2.6 Analysis

Paired-end sequencing reads were trimmed to a minimum length of 30 and aligned using AdapterRemoval v.2 (RRID : SCR_011834, Schubert et al., 2016). Chimeras were removed and non-chimeric sequences were clustered into Amplicon Sequence Variants (ASVs) using VSEARCH v.2.21.1 (Rognes et al., 2016). Amplicon sequence variants were classified against the EzTaxon database (Yoon et al., 2017) and the antifungal isolates database compiled by Woodhams et al. (2015). If an ASV was a 99% or greater match to the EzTaxon database, it was identified to species; 94.5–99% matches were identified to genus and matches that were 88–94.5% identical were identified to family. We included only matches of 99% or above to the antifungal database, following the methods of Barnes et al. (2021). Sequences were aligned with MAFFT (Katoh, 2002). The resulting ASV table, taxa table, metadata, and phylogenetic tree of microbes were imported into R and combined into a phyloseq object using the package phyloseq (RRID : SCR_013080, McMurdie and Holmes, 2013). Samples were rarefied to 200, 500, 1000, and 2000 reads per sample; results listed in the main text refer to rarefaction at 1000 reads (Supplementary Figure 1), a threshold chosen to remove samples with low sequence counts (Weiss et al., 2017).

We calculated three alpha diversity metrics (observed ASVs, Shannon-Wiener Index, and Inverse Simpson’s Index) using the R package phyloseq. These metrics were chosen to provide a range of complimentary insights into the microbiome. For example, observed ASVs provide the number of unique sequence variant clusters in a sample—a direct measure of taxonomic diversity. In contrast, the Shannon-Wiener Index accounts for evenness in addition to the number of distinct taxa, reducing the impact of rare variants. Finally, the Inverse Simpson’s Index accounts for species richness and evenness but weights microbe species by abundance, minimizing the impact of rare taxa. Shannon-Wiener and the Inverse Simpson’s indices are slightly more sensitive to species richness and species evenness, respectively (Johnson and Burnet, 2016). Alpha diversity was analyzed using Kruskal-Wallis tests for explanatory variables with only two qualitative responses and a pairwise Wilcoxon rank-sum test with the Holm p-value correction to compare the remaining qualitative variables.

Because the factors shaping microbiome diversity likely act simultaneously, we examined the combined effects of the phylogenetic, environmental, and ecological explanatory variables using two complementary approaches. First, we applied a Classification and Regression Tree (CART) framework (R package partykit; Hothorn and Zeileis, 2015), with family, genus, and species included as categorical variables alongside additional predictors. The CART analysis was performed using the conditional inference tree command (ctree) with the number of input variables per tree randomly sampled and no restrictions applied to the number of splits in the tree. Second, we used a phylogenetic generalized linear mixed model (PGLMM; R packages brms, geiger, and PhyloOrchard; O’Meara et al., 2013; Pennell et al., 2014; Bürkner, 2017) to analyze alpha diversity while accounting for phylogenetic non-independence using a phylogeny based upon maximum likelihood inference of 12-genes (12712 bp) for 2871 species from Pyron and Wiens (2011). We used only the species relevant to our study in the model. To prepare data for the PGLMM, continuous explanatory variables were tested for collinearity and removed if they had a variance inflation factor (VIF) greater than four, which represents a conservative metric for minimizing covarying explanatory data (Johnston et al., 2018). Mean annual temperature and mean seasonal precipitation were not included in the model due to VIF > 4. The PGLMM was run with two random effects (phylogeny and species), ten chains, 20,000 iterations (thinned every 10 samples), and a Gaussian response. Phylogenetic signal (h2) was estimated from PGLMMs as the proportional variance explained by phylogeny, using intercept only models (with reduced iterations). The CART and PGLMM analyses were run on each alpha diversity index separately to quantify and contrast the most important drivers of alpha diversity.

To understand differences in abundance and beta diversity microbial community composition was evaluated in the program QIIME2 (RRID : SCR_008249, Bolyen et al., 2019). To identify specific microbes that were more abundant between different families or in the presence/absence of infection, differential abundances were calculated using the packages DEseq2 and ANCOM-BC in R (Love et al., 2014; Lin and Peddada, 2020). To understand differences among microbial communities and the factors driving them, beta diversity metrics were calculated using the UniFrac command to obtain unweighted- and weighted-UniFrac distances in the phyloseq package (Lozupone and Knight, 2005). The unweighted test evaluates the presence and absence of ASVs, whereas the weighted test accounts for the abundance of ASVs. Both tests take phylogenetic distance into account. We generated a phylogeny for our microbial taxa by aligning data with MAFFT (Katoh, 2002) and using the Maximum-Likelihood algorithm implemented in FastTree (v2, Price et al., 2010) with a Jukes-Cantor + CAT model of nucleotide evolution. Beta diversity was analyzed using Permutational Multivariate Analysis of Variance tests on UniFrac results (PERMANOVA; Anderson, 2001) with Holm p-value correction and permutations set to 10,000. To investigate how beta diversity changes by climate while controlling for geographic distance, we conducted partial Mantel tests using the Bray-Curtis matrix of beta diversity values and a geodesic distance matrix and controlled for the geographic distance using the Euclidean distance between sample locations with the R packages vegan and ecodist (RRID : SCR_011950, Goslee and Urban, 2007; Oksanen et al., 2022).

3 Results

The initial dataset consisted of 287 epithelial swab samples (Supplementary Data). After sequencing and quality filtering, there were 280 samples with adequate sequencing depths. In total, 5,256,288 sequences were obtained from these 280 samples, with a total of 13,038 ASVs. The total number of ASVs present in each sample varied widely, from 1 to 429,627, with an average of 11,883 ASVs ± 3,252. This yield is comparable to the number of sequences in other studies on amphibian microbiomes (Muletz-Wolz et al., 2018; Buttimer et al., 2022; 2019). To adjust for differences in ASV sample sizes, samples were rarefied to 1000 reads per sample (Supplementary Figure 1; Supplementary Table 1). After rarefaction, 179 samples remained, with the number of unique ASVs ranging from 22–610 (average 103; Supplementary Data). Of these, 147 animals were tested for Bd presence (59 positive and 88 negative) and 151 were tested for ranavirus (21 positive and 130 negative; 26 were determined from tail clips, 122 from liver samples, one from both liver and tail clip samples, and two from either tail clip or liver samples), resulting in disease status data for both pathogens for 146 animals (see Supplementary Data).

Hereafter, results refer to the final quality filtered dataset. Animals in the family Salamandridae (n = 71) had a skin microbial community composition of 37.3% Comamonadaceae, 10% Enterobacteriaceae, 7% Moraxellaceae, 6% Oxalobacteraceae, 5% Pseudomonadaceae, and 2% Flavobacteriaceae; 34% of the remaining microbes belonged to families that comprised less than 1% of the community (Supplementary Figures 2, 3). The majority of the microbial community on animals in family Plethodontidae (n = 108) was composed of microbial families that made up less than 1% of the community, which combined made up 56% of the microbiome. Plethodontids had a remaining skin microbial composition of 16% Pseudomonadaceae, 9% Enterobacteraceae, 7% Moraxellaceae, 4% Flavobacteraceae, 4% Oxalobacteraceae, and 3% Comamonadaceae (Figure 2; Supplementary Figure 2, 3). There were 79 differentially abundant microbes identified between the Plethodontidae and Salamandridae families by DESeq2 and 79 microbes identified by the ANCOM-BC, with 28 microbial taxa identified by both analyses (Supplementary Tables 3, 4).

Figure 2

Figure 2 Relative compositional differences in skin microbiome communities based on host family. Plethodontids (left) had consistently higher alpha diversity than the Salamandridae (i.e., newts, right) and were distinct from the newts in analyses of beta diversity.

3.1 Alpha diversity

3.1.1 Relationship to phylogenetic divergence among hosts

When contrasting alpha diversity by taxonomic grouping, without considering other explanatory variables, the family Plethodontidae had consistently higher microbial diversity than Salamandridae across all three metrics: observed ASVs, Shannon-Wiener Index, and Inverse Simpson’s Index (significance assessed with Kruskal-Wallis tests; Supplementary Table 4). This same pattern was observed between genera and species, with plethodontid species generally exhibiting greater alpha diversity than newts across all metrics; however, plethodontids were not significantly different from one another at the level of host genus or species (Supplementary Figure 4; Supplementary Tables 4–6).

When host taxa were evaluated in the context of other explanatory variables in our CART analyses, a similar split was evident for both observed ASVs and the Shannon-Wiener Index: family-level divergence was the predominant split in the dataset (Figure 3; Supplementary Figure 5). For the Inverse Simpson’s Index, the primary factor that categorized microbiome diversity was related to species, with E. lucifuga and P. angusticlavius differing from all other taxa, followed by differences between plethodontids and newts (Supplementary Figure 5). However, we found low phylogenetic signal using PGLMMs [ASVs: h² = 0.01 (0.0–0.05 credible interval); Shannon: h² = 0.01 (0.0–0.03 credible interval); and Inverse Simpson: h² = 0.01 (0.0–0.06 credible interval)]. Despite detecting family-level differences in alpha diversity there was no evidence for phylogenetic signal below the level of family (i.e., at the level of genera or species). We found insufficient evidence for a phylogenetic signal for skin microbial diversity, despite detecting family-level differences.

Figure 3

Figure 3 Results from the CART analysis of Observed ASVs. While family is the most important factor influencing the number of unique ASVs, environmental variables are important within family Salamandridae. Similar dominance of differentiation in microbial diversity based upon host phylogeny is observed in the other two metrics (Supplementary Figures 2, 3).

3.1.2 Relationship to environmental differences among hosts

For univariate analyses of the effect of ecoregion on microbiome diversity, we subdivided the data by host family and only contrasted species present in more than one ecoregion. Newts were present in all four ecoregions, but ecoregion was not a significant predictor of skin microbiome diversity for all comparisons in all alpha diversity metrics (Supplementary Table 5). We evaluated the ecoregion in plethodontids by excluding species that were not present in multiple ecoregions, resulting in three species present in two ecoregions (Ozark Highlands and the Boston Mountains; Figure 1; Supplementary Table 1). Across these three species (combined), individuals from the Ozark Highlands had significantly more microbial ASVs than those from the Boston Mountains (Supplementary Table 6).

Because differences among ecoregion may actually be the result of climatic differences among sites, as opposed to the biotic community (both of which are captured by ecoregion), we paired ecoregional comparisons with analyses directly using the underlying climatic data. When contrasting multiple predictor variables using CART analyses, environmental variables pertaining to temperature and precipitation either were not important or were secondary predictors of skin microbial diversity (Figure 3; Supplementary Figures 5, 6). However, when explicitly accounting for phylogeny in the PGLMMs, several abiotic variables had weak but significant positive relationships with alpha diversity: elevation, mean temperature of the month of collection, mean precipitation of season of collection, and mean precipitation of year (Figure 4).

Figure 4

Figure 4 Relative contribution of explanatory variables to alpha diversity based on phylogenetic generalized linear mixed models. The figure contains posterior means and their 95% credible intervals. After accounting for phylogenetic non-independence, we found that ranavirus status has no impact on alpha diversity, but Bd presence has a slight negative impact. Warmer, wetter, and higher elevation sampling sites had slight positive effects on alpha diversity, but terrestrial habitats and some life history stages reduced microbial diversity.

3.1.3 Relationship to ecological differences among hosts

Univariate comparisons of host ecology were limited to contrasts of salamander species that had samples of both aquatic and terrestrial individuals (E. longicauda and E. tynerensis). For these species, habitat was not a significant factor (χ² < 0.76, p > 0.38 for both species across all alpha diversity metrics; Supplementary Table 7). When considered relative to other factors that partition variance in alpha diversity, CART analyses did not identify life stage or habitat as significant drivers of microbial diversity (Figure 3, Supplementary Figures 5, 6). However, when considered jointly with phylogeny and other explanatory factors in the PGLMM, a terrestrial ecology was negatively associated with alpha diversity for the observed ASVs and Inverse Simpson metrics, and skin microbiomes of both paedomorphic adult and larval plethodontids were associated with lower alpha diversity, while those of juvenile plethodontids had slightly increased alpha diversity (Figure 4).

3.2 Beta diversity

3.2.1 Phylogenetic differences among hosts

Differences in microbial community composition (Figure 5), as measured by beta diversity, varied significantly in pairwise comparisons between salamander families, genera, and some species in both the weighted- and unweighted-UniFrac analyses (Supplementary Table 7). Notophthalmus viridescens exhibited distinctly different community composition when compared to each of the other salamander species in both analyses (p < 0.01). Similarly, microbial communities on E. tynerensis differed significantly from all other species (unweighted: p ≤ 0.03, r² ≤ 0.05; weighted: p ≤ 0.045, r² ≤ 0.05), except when compared to P. angusticlavius (p > 0.32, r² ≤ 0.03 in both tests). Besides E. tynerensis, none of the other plethodontid microbiomes differed significantly between species (p > 0.1, r² ≤ 0.05 in either metric; Supplementary Table 7).

Figure 5

Figure 5 Principal coordinate analysis of beta diversity values in families Plethodontidae and Salamandridae.

3.2.2 Environmental differences among hosts

To understand the role of the environment in shaping the composition of salamander microbiomes, we compared samples within species present in multiple ecoregions. For N. viridescens, all ecoregion comparisons were significantly different from each other with the exception of Arkansas Valley vs. Boston Mountains (p > 0.05, r² = 0.06) and Boston Mountains vs. Ozark Highlands (p = 0.21, r² = 0.02) in the unweighted-UniFrac comparisons (Supplementary Table 8). In the weighted-UniFrac, microbiomes from all ecoregions were significantly distinct, apart from the Arkansas Valley and Ozark Highlands comparison (p = 0.16, r² = 0.02). When comparing beta diversity across the three plethodontid species present in multiple ecoregions (E. tynerensis, P. albagula, and P. angusticlavius), there was no difference between the Boston Mountains and Ozark Highlands in any comparison (p > 0.1, r² < 0.05; Supplementary Table 9).

We also used Mantel and partial Mantel tests to examine how beta diversity is shaped by climatic differences among sites. As results were similar with and without accounting for geographic distance (Supplementary Tables 10, 11), all following values thus refer to partial Mantel results. Overall, the mean temperature of the month had the strongest effect on beta diversity (all: r = 0.21; newts: r = 0.26; plethodontids: r = 0.17; p < 0.01 in each case); the mean temperature of the season and mean annual temperature had significant but weaker positive relationships with beta diversity (Supplementary Table 10). Elevation had a significant effect on beta diversity across all taxonomic groups but was not significant when considering salamander families separately after correcting for multiple comparisons (all: r = 0.12, p < 0.01; just newts: r = 0.14, p = 0.03; plethodontids: r = 0.11, p = 0.05). In newts, there was a positive relationship between mean precipitation during the month of capture and beta diversity (r = 0.30, p < 0.01), although no other climate variables influenced microbial community structure. In contrast, plethodontids had significant positive relationships with mean annual, seasonal, and monthly precipitation, with the latter being the strongest (annual: r = 0.11, p = 0.01; seasonal: r = 0.17, p < 0.01; monthly: r = 0.26, p < 0.01).

3.2.3 Differences in beta diversity related to host ecology

Microbiomes from hosts found in aquatic and terrestrial habitats differed significantly in both unweighted- and weighted-UniFrac metrics (p < 0.01) when all samples were compared together. However, when testing only species that were sampled in both terrestrial and aquatic habitat types (E. longicauda and E. tynerensis), statistical significance differed based on whether the UniFrac test was weighted (unweighted: r² = 0.02, p = 0.02; weighted: r² = 0.02, p = 0.32; Supplementary Table 9).

To evaluate how life stage impacted beta diversity, we only examined plethodontids, as all newts in this study were adults (Supplementary Table 9). Juvenile plethodontids were significantly different from terrestrial adults, paedomorphic adults, and larvae in the unweighted UniFrac but did not differ from any other life stages in the weighted UniFrac. Adult plethodontids did not differ from paedomorphic adults but did differ from larvae in the unweighted UniFrac but not the weighted UniFrac. Paedomorphic adults did not differ in beta diversity from larvae (Supplementary Table 9).

3.3 Pathogen status and salamander microbiomes

3.3.1 Presence of the fungal pathogen B. dendrobatidis

Animals without Bd present on their skin (n = 88) had skin microbial communities composed of Proteobacteria (59.30%), Firmicutes (12.99%), Bacteroidetes (13.92%), and Actinobacteria (5.40%). While microbial communities of Bd+ animals (n = 59) were composed of the same four groupings, they differed in terms of relative composition: Proteobacteria (81.62%), Firmicutes (6.64%), Bacteroidetes (6.47%), and Actinobacteria (3.06%). Two hundred fifty-nine of the 13,038 ASVs were identified to sequences in the antifungal isolates database. Of these, ten microbe species were differentially abundant between animals with Bd in the DESeq2 analysis and six were differentially abundant between animals without Bd in the ANCOM-BC analysis (Supplementary Tables 12, 13). When considering microbes without antifungal properties, we identified six differentially abundant microbes between animals with and without Bd present on the skin in the DESeq2 analysis and 23 in the ANCOM-BC analysis (Supplementary Tables 14, 15). However, most Bd+ animals were newts (N. viridescens; Supplementary Data). The presence of Bd covaried with family-level classifications; when separated by family, no univariate tests showed differences in alpha diversity by Bd status (Supplementary Tables 5, 6). In contrast, when accounting for phylogenetic dependence and other predictors with a PGLMM, the presence of Bd was associated with a slight decrease in alpha diversity in the observed ASVs, Shannon, and Inverse Simpson metrics (Figures 4, 6).

Figure 6

Figure 6 Alpha diversity differences based on Bd infection status. Alpha diversity of plethodontids and newts (left and right, respectively, within each set of plots) in the presence and absence of Bd. Overall, Bd+ salamanders exhibit a slight but non-significant (Supplementary Tables 4 & 5) reduction in alpha diversity of the skin microbiome.

3.3.2 Presence of ranavirus

Ranavirus was not a significant predictor of alpha diversity, regardless of subdividing data among plethodontids or only newts (p = 1.0 for all univariate tests; Supplementary Tables 5, 6) and was not significant in either combined test of explanatory variables (Figures 3, 4). Similarly, ranavirus status was not significant in either of the beta diversity analyses (r² ≤ 0.01, p ≥ 0.34 for either metric; Supplementary Tables 8, 9).

4 Discussion

We set out to determine the relative importance of potential drivers of skin microbial diversity—phylogenetic distance, host ecology, host environment, and the presence of pathogens—in salamanders. Our results indicate that family-level differences (i.e., differences between Salamandridae and Plethodontidae) are consistently one of the strongest factors shaping the skin microbiomes at this taxonomic and geographic scale. Further analyses indicate these differences do not reflect a broader relationship with phylogenetic distance, but instead suggest that other differences that co-vary at the level of family may be influencing skin microbial diversity. Environmental and ecological factors play a secondary role in shaping both microbial (i.e., alpha) diversity and community differentiation (beta diversity), particularly among species. Although the presence of pathogens was not the most important factor determining microbiome structure or diversity, there was a slight, negative impact of Bd on the alpha diversity of salamander skin microbiomes. Together, these results provide further insight into the importance of multiple factors that simultaneously impact the skin microbiome in the presence of pathogens and underscore the importance of accounting for host family in microbiome studies.

4.1 Family-level differences significantly shape microbiomes

When we examined multiple drivers simultaneously, or accounted for their impacts by subdividing our data, the most important and consistent factor influencing salamander skin microbiomes was the family of the host, with differences in host ecology and climate becoming more important at the interspecific level. However, overall phylogenetic distance did a poor job of explaining these patterns, and most univariate comparisons among species or genera did not indicate increasing microbial dissimilarity with greater host differentiation. These results are generally similar to those of Buttimer et al. (2022), who recovered differences related to both host genera and habitat quality. In our study, family Salamandridae was distinct from Plethodontidae in both alpha and beta diversity metrics–with plethodontids having consistently higher alpha diversity than the newts–as well as overall microbial community membership (as indicated by beta diversity).

Differences in the skin microbiomes that correlate with phylogeny are sometimes attributed to phylosymbiosis, in which the microbes present in or on an organism are a result of long-term natural selection (Brucker and Bordenstein, 2012; Yeoh et al., 2017; Lim and Bordenstein, 2020). Families Plethodontidae and Salamandridae have been distinct for over 150 million years, and the genera Eurycea and Plethodon within Plethodontidae diverged around 75 million years ago (Zhang and Wake, 2009; Martin et al., 2016). This could provide sufficient evolutionary time for phylosymbiosis to develop, potentially driven by fitness consequences for animals that experience changes in the community composition of the skin microbiome (Becker et al., 2015; Bletz et al., 2018). However, while the strongest predictor of microbial abundance and diversity from most analyses in our study was host family (i.e., the greatest level of divergence in our dataset), we did not detect a strong phylogenetic signal at lower taxonomic levels. This suggests thatphylosymbiotic relationships are not the predominant driver of microbiome communities in these salamanders, or we are not able to identify them within our data at this phylogenetic scale. Differences in terms of the magnitude and sign of responses across phylogenetic scales are common (Graham et al., 2018), and are worth further consideration in this system with a broader phylogenetic and geographical dataset.

When comparing microbiomes within plethodontids, we did not find that genus was a significant predictor of microbiome composition, but we found instances where alpha or beta diversity of some host species differed significantly from others However, there was no clear pattern to the variation among hosts. For example, beta diversity of E. tynerensis differed significantly from all other salamander species in the unweighted-UniFrac test (Supplementary Table 7). In terms of alpha diversity, the CART analysis of Inverse Simpson values initially split E. lucifuga and P. angusticlavius from all other species. With no evidence for phylogenetic signals in host identity related to microbiomes, it is possible that unaccounted life history traits, behaviors, or diet are influencing these differences in microbial diversity. Microhabitat preferences also suggest one possibility linking the two species: adult E. lucifuga and P. angusticlavius both dwell in shallow caves and rock ledges (Petranka and Scott, 1998). Identifying and understanding these unaccounted drivers of among-species differences may allow us to better predict how changes in food or habitat availability may drive microbiome shifts in the future.

4.2 Differences in environment and host ecology

Although we identified environmental and ecological variables as secondary factors in structuring microbial communities in our study, they still played a significant role in explaining microbiome diversity. Although there were no differences by ecoregion, increases in average monthly precipitation, average temperature, and site elevation were associated with increases in alpha diversity. Factors relating to host ecology, such as habitat and life stage, were also significant in the PGLMMs, with juvenile life stage and terrestrial habitat exhibiting lower alpha diversity in comparison to other habitats and life stages. These results align with other studies, where host ecology and environment are important drivers of skin microbial diversity (e.g., Bletz et al., 2017a; Bird et al., 2018; Muletz Wolz et al., 2018; Smith et al., 2023), but these studies differ from one another in terms of the relative importance of predictors.

Although elevation, average precipitation, and average temperature of the month were secondary or tertiary drivers of microbial diversity, all had positive relationships with alpha diversity measures. When considered jointly with phylogeny, local climatic variables predicted alpha diversity; higher temperature or precipitation were associated with greater alpha diversity of microbial communities (Figure 4). None of these explanatory variables related to climate had a strong relationship to community membership alone, but after accounting for geographic distance among sites, partial Mantel tests indicated that microbial communities become more distinct between sites with different temperatures and those with different precipitation in each of the salamander families (Supplementary Table 10). We also compared the effect of ecoregion on the skin microbiome. Our analyses show that while there was no difference in microbial alpha diversity across the ecoregions, each ecoregion significantly differed in community composition from one or more other ecoregions in analyses of beta diversity. Our finding that abiotic factors are more important at the level of family are qualitatively similar to studies that focused on salamander species within the same family or genus (e.g., Kueneman et al., 2014; Walker et al., 2020), further supporting the idea that abiotic environmental variables become more important in structuring the microbiomes of closely related hosts.

When considering explanatory variables relating to host ecology and habitat, we expected to find significant differences in skin microbial diversity between terrestrial and aquatic environments, as well as between different life stages, similar to previous work across habitat types (Xu et al., 2020). The species in our study differ greatly in their life cycles. Adult newts, the only life stage of newts included herein, are aquatic (Petranka and Scott, 1998). In contrast, Plethodon spp. lay eggs on land and spend their entire lives as terrestrial animals; the genus Eurycea lays eggs in water, which hatch into aquatic larvae (Petranka and Scott, 1998). Larval Eurycea either transform into terrestrial juvenile forms which become terrestrial adults or remain in the water to become aquatic adults (Petranka and Scott, 1998). Despite these varied life histories, the difference in alpha diversity between terrestrial and aquatic environments was significant only in the PGLMM, which controls for phylogenetic dependence among hosts rather than directly assessing the effect of taxonomic categories. The difference in alpha diversity between life stages was also only evident in the PGLMM, in which being a juvenile plethodontid had a negative impact on the alpha diversity of the skin microbiome.We did find that juvenile plethodontids have significantly different beta diversity from larval, adult, and paedomorphic adult life stages in the unweighted-UniFrac analyses. Larval plethodontids were also significantly distinct from adults in the unweighted-UniFrac and from paedomorphic adults in the weighted-UniFrac; additionally, the two habitat types did have significantly different beta diversity in the unweighted-UniFrac analysis (Supplementary Table 9). The beta diversity of salamander microbiomes was not significantly different across different habitats when accounting for the phylogenetic distance and abundance among microbes (i.e., using the weighted UniFrac results; Supplementary Table 9). Overall, differences in beta diversity among life stage and habitat are less pronounced when accounting for abundance.

4.3 The presence of infection

We expected to find the presence of Bd on skin to be a substantial driver of the skin microbiome’s community composition, similar to previous work on Anaxyrus boreas toads and the closely related Plethodon cinereus salamanders (Kueneman et al., 2016; Muletz Wolz et al., 2019). After taking host phylogeny into account, we found that animals with Bd present on the skin showed slightly lower alpha diversity (Figures 4, 6). Though only the PGLMM analyses provided statistical support for lower alpha diversity in the presence of disease (e.g., Figure 4), these results fit into our growing understanding of the interaction between chytridiomycosis and amphibian skin microbiomes (e.g., Harris et al., 2009; Bates et al., 2022).

Chytridiomycosis infects the skin of amphibians and disturbs osmoregulation and respiration in its host (Voyles et al., 2007; Rollins-Smith et al., 2011; Campbell et al., 2012; Martel et al., 2014, 2009). Recent studies have documented the significant prevalence of Bd and RV in North America, including among populations of amphibians in Oklahoma (Marhanka et al., 2017; Watters et al., 2018; Davis et al., 2019; Smith et al., 2019, 2019; Watters et al., 2021). Though none of our animals exhibited visible symptoms of disease (e.g., lethargy, increased skin shedding, and skin lesions; Martel et al., 2013, 2014; Voyles et al., 2007, 2009), the presence of Bd on the animal’s skin may indicate an ongoing or recent infection. Varying stages of infection have been shown to have very different impacts on amphibian skin microbiomes; while an overall decline in bacterial abundance has been observed related to Bd infections (Medina et al., 2017), this pattern involves the increase of some microbes and the decrease of others during and after infection (Jani and Briggs, 2014; Jani et al., 2021). Salamanders from our study may be at different stages of infection, which may minimize overall differences among “positive” animals if the microbiome’s alpha and beta diversity change during and again after infection, which could merit further study. Additionally, because we were not able to address changes in the function of the microbes present, there may be large differences in functional microbes but fewer differences in overall alpha or beta diversity (Jani et al., 2021).

Our results also indicate the presence of antifungal microbes in the skin microbiome (Supplementary Table 13), a feature that has been identified in multiple studies (Becker et al., 2015; Woodhams et al., 2015). Laboratory manipulation of antifungal microbes on amphibians’ skin in has been linked to improved chytridiomycosis survival rates (Harris et al., 2009; Bletz et al., 2013; Kueneman et al., 2016; Bletz et al., 2018), which suggests beneficial microbes may naturally increase in abundance during infection in animals that tolerate chytridiomycosis (Campbell et al., 2019). Although there were no significant differences in alpha or beta diversity of antifungal microbes between animals with Bd status in our study, there were multiple differentially abundant antifungal microbe species when contrasting animals with and without Bd (Supplementary Tables 12, 13). However, microbes with antifungal properties may decline in number during infection, suggesting a more nuanced relationship between microbial function and microbiome composition (Medina et al., 2017). More studies about Bd and its effect on the amphibian skin microbiome are needed to aid in answering these questions.

In the case of RV infections, we expected possible virus-driven variation in the composition of the skin microbiome, based on previous work that indicated amphibians with ranavirus infection exhibit changes in microbial community composition (Campbell et al., 2019). However, as infection with RV takes place internally, it is therefore possible that the presence of RV is less likely to influence the skin microbiome than Bd (Gray and Chinchar, 2015). In our study, there was not a significant relationship between ranavirus presence and alpha or beta diversity in any tests. This does not necessarily indicate that there is no relationship between ranavirus infection and changes in the skin microbiome. A recent study has found that the composition of the gut microbiome of Rana sylvatica may influence stress response to ranaviral infection (Hughey et al., 2023). Future studies focused on symptomatic animals or gut microbiomes may provide greater insights into the host-ranavirus-microbiome relationship.

5 Conclusions

Effective conservation planning requires an understanding of not only a target species, but also the impact of its interactions with other species, including microbes (Hird, 2017). This study provides a snapshot of the skin microbiomes of six salamander species and how these microbiomes vary across different habitats, varying environments, and phylogenetic distances. While we found that the most important determinant of skin microbiome composition is the family of the host, the research herein suggests that the skin microbiomes are not more divergent as phylogenetic distance among hosts increases (as would be expected under phylosymbiosis theory). After family-level differences are accounted for, host ecology, host environment, and the presence of disease all influence the diversity and community composition of microbiomes. Furthermore, the presence of differentially abundant microbes linked to antifungal properties on the skin of animals testing positive for Bd suggests a possible interaction occurring between pathogens and skin microbiomes. As infectious diseases continue to decimate amphibian populations, characterizing host microbiomes and examining the function of antifungal microbes is more important than ever.

Data availability statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://www.ncbi.nlm.nih.gov/, PRJNA1002333.

Ethics statement

The animal study was approved by University of Oklahoma Institutional Animal Care and Use Committee. The study was conducted in accordance with the local legislation and institutional requirements.

Author contributions

MK: Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Software, Visualization, Writing – original draft, Writing – review & editing. SS: Investigation, Validation, Writing – review & editing. DB: Methodology, Software, Supervision, Writing – review & editing, Validation. JW: Writing – review & editing, Data curation, Funding acquisition, Investigation, Project administration. KM: Funding acquisition, Supervision, Writing – review & editing, Validation. CS: Funding acquisition, Validation, Writing – review & editing, Conceptualization, Investigation, Methodology, Project administration, Resources, Supervision. HL: Conceptualization, Methodology, Supervision, Validation, Writing – review & editing, Formal analysis, Visualization, Writing – original draft.

Funding

The author(s) declare financial support was received for the research, authorship, and/or publication of this article. This work was financially supported by the Oklahoma Department of Wildlife Conservation OK-WILD Grants #F14F01225 (T-80-1), #F20AF00023 (T-116-R-1), and #F20AF10405 (T-118-R-1); the Sam Noble Oklahoma Museum of Natural History; and the OU Department of Biology.

Acknowledgments

Our sincere thanks to Dr. Krithivasan Sankaranarayanan for help with analyses, and to Katherine Stroh, Gracie Hedgpeth, Alexandria Fulton, Emma Clary, Sam Eliades, Kai Wang, Jordon Henderson, Claudia Goss, and many other Siler Lab members for field sample collection, pathogen DNA extraction, and qPCR screening.

Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

The author(s) declared that they were an editorial board member of Frontiers, at the time of submission. This had no impact on the peer review process and the final decision.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/frmbi.2024.1277645/full#supplementary-material

References

Anderson M. J. (2001). A new method for non-parametric multivariate analysis of variance. Austral. Ecol. 26, 32–46. doi: 10.1111/j.1442-9993.2001.01070.pp.x

CrossRef Full Text | Google Scholar

Barnes E. M., Kutos S., Naghshineh N., Mesko M., You Q., Lewis J. D. (2021). Assembly of the amphibian microbiome is influenced by the effects of land-use change on environmental reservoirs. Environ. Microbiol. 23, 4595–4611. doi: 10.1111/1462-2920.15653

PubMed Abstract | CrossRef Full Text | Google Scholar

Bates K. A., Friesen J., Loyau A., Butler H., Vredenburg V. T., Laufer J., et al. (2022). Environmental and anthropogenic factors shape the skin bacterial communities of a semi-arid amphibian species. Microb. Ecol. 86, 1393–1404. doi: 10.1007/s00248-022-02130-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Bates K. A., Shelton J. M. G., Mercier V. L., Hopkins K. P., Harrison X. A., Petrovan S. O., et al. (2019). Captivity and infection by the fungal pathogen Batrachochytrium salamandrivorans perturb the amphibian skin microbiome. Front. Microbiol. 10. doi: 10.3389/fmicb.2019.01834

CrossRef Full Text | Google Scholar

Becker M. H., Harris R. N. (2010). Cutaneous bacteria of the redback salamander prevent morbidity associated with a lethal disease. PloS One 5, e10957. doi: 10.1371/journal.pone.0010957

PubMed Abstract | CrossRef Full Text | Google Scholar

Becker M. H., Walke J. B., Cikanek S., Savage A. E., Mattheus N., Santiago C. N., et al. (2015). Composition of symbiotic bacteria predicts survival in Panamanian golden frogs infected with a lethal fungus. Proc. R. Soc B. 282, 20142881. doi: 10.1098/rspb.2014.2881

CrossRef Full Text | Google Scholar

Bird A. K., Prado-Irwin S. R., Vredenburg V. T., Zink A. G. (2018). Skin microbiomes of California terrestrial salamanders are influenced by habitat more than host phylogeny. Front. Microbiol. 9. doi: 10.3389/fmicb.2018.00442

CrossRef Full Text | Google Scholar

Bletz M. C., Archer H., Harris R. N., McKenzie V. J., Rabemananjara F. C. E., Rakotoarison A., et al. (2017a). Host ecology rather than host phylogeny drives amphibian skin microbial community structure in the biodiversity hotspot of Madagascar. Front. Microbiol. 8. doi: 10.3389/fmicb.2017.01530

CrossRef Full Text | Google Scholar

Bletz M. C., Kelly M., Sabino-Pinto J., Bales E., Van Praet S., Bert W., et al. (2018). Disruption of skin microbiota contributes to salamander disease. Proc. R. Soc B: Biol. Sci. 285, 20180758. doi: 10.1098/rspb.2018.0758

CrossRef Full Text | Google Scholar

Bletz M. C., Loudon A. H., Becker M. H., Bell S. C., Woodhams D. C., Minbiole K. P. C., et al. (2013). Mitigating amphibian chytridiomycosis with bioaugmentation: characteristics of effective probiotics and strategies for their selection and use. Ecol. Lett. 16, 807–820. doi: 10.1111/ele.12099

PubMed Abstract | CrossRef Full Text | Google Scholar

Bletz M. C., Myers J., Woodhams D. C., Rabemananjara F. C. E., Rakotonirina A., Weldon C., et al. (2017b). Estimating herd immunity to amphibian chytridiomycosis in Madagascar based on the defensive function of amphibian skin bacteria. Front. Microbiol. 8. doi: 10.3389/fmicb.2017.01751

CrossRef Full Text | Google Scholar

Bolyen E., Rideout J. R., Dillon M. R., Bokulich N. A., Abnet C. C., Al-Ghalith G. A., et al. (2019). Reproducible, interactive, scalable and extensible microbiome data science using QIIME 2. Nat. Biotechnol. 37, 852–857. doi: 10.1038/s41587-019-0209-9

PubMed Abstract | CrossRef Full Text | Google Scholar

Boyle D. G., Boyle D. B., Olsen V., Morgan J. A. T., Hyatt A. D. (2004). Rapid quantitative detection of chytridiomycosis (Batrachochytrium dendrobatidis) in amphibian samples using real-time TaqMan PCR assay. Dis. Aquat. Organ. 60, 141–148. doi: 10.3354/dao060141

PubMed Abstract | CrossRef Full Text | Google Scholar

Brucker R. M., Bordenstein S. R. (2012). The roles of host evolutionary relationships (genus: Nasonia) and development in structuring microbial communities: divergence in Nasonia microbial communities. Evolution 66, 349–362. doi: 10.1111/evo.2012.66.issue-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Bürkner P.-C. (2017). brms: an R package for Bayesian multilevel models using Stan. J. Stat. Soft. 80, 1–28. doi: 10.18637/jss.v080.i01

CrossRef Full Text | Google Scholar

Buttimer S., Hernández-Gómez O., Rosenblum E. B. (2022). Skin bacterial metacommunities of San Francisco Bay Area salamanders are structured by host genus and habitat quality. FEMS Microb. Ecol. 97, fiab162. doi: 10.1093/femsec/fiab162

CrossRef Full Text | Google Scholar

Campbell C. R., Voyles J., Cook D. I., Dinudom A. (2012). Frog skin epithelium: electrolyte transport and chytridiomycosis. Int. J. Biochem. Cell Biol. 44, 431-434. doi: 10.1016/j.biocel.2011.12.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Campbell L. J., Garner T. W. J., Hopkins K., Griffiths A. G. F., Harrison X. A. (2019). Outbreaks of an emerging viral disease covary with differences in the composition of the skin microbiome of a wild United Kingdom amphibian. Front. Microbiol. 10. doi: 10.3389/fmicb.2019.01245

CrossRef Full Text | Google Scholar

Cheng T. L., Rovito S. M., Wake D. B., Vredenburg V. T. (2011). Coincident mass extirpation of neotropical amphibians with the emergence of the infectious fungal pathogen Batrachochytrium dendrobatidis. Proc. Natl. Acad. Sci. U.S.A. 108, 9502–9507. doi: 10.1073/pnas.1105538108

PubMed Abstract | CrossRef Full Text | Google Scholar

Chimenos-Küstner E., Giovannoni M. L., Schemel-Suárez M. (2017). Dysbiosis as a determinant factor of systemic and oral pathology: importance of microbiome. Med. Clin. 7, 305. doi: 10.1016/j.medcle.2017.05.023

CrossRef Full Text | Google Scholar

Daszak P., Berger L., Cunningham A. A., Hyatt A. D., Earl Green D., Speare R. (1999). Emerging infectious diseases and amphibian population declines. Emerg. Infect. Dis. 5, 735–748. doi: 10.3201/eid0506.990601

PubMed Abstract | CrossRef Full Text | Google Scholar

Davis D. R., Farkas J. K., Kruisselbrink T. R., Watters J. L., Ellsworth E. D., Kerby J. L., et al. (2019). Prevalence and distribution of ranavirus in amphibians from southeastern Oklahoma, USA. Herpetol. Conserv. Biol. 14, 360–369.

Google Scholar

DiRenzo G. V., Longo A. V., Muletz-Wolz C. R., Pessier A. P., Goodheart J. A., Lips K. R. (2021). Plethodontid salamanders show variable disease dynamics in response to Batrachochytrium salamandrivorans chytridiomycosis. Biol. Invasions 23, 2797–2815. doi: 10.1007/s10530-021-02536-1

CrossRef Full Text | Google Scholar

Ellison S., Rovito S., Parra-Olea G., Vásquez-Almazán C., Flechas S. V., Bi K., et al. (2019). The influence of habitat and phylogeny on the skin microbiome of amphibians in Guatemala and Mexico. Microb. Ecol. 78, 257–267. doi: 10.1007/s00248-018-1288-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Esselstyn J. A., Garcia H. J. D., Saulog M. G., Heaney L. R. (2008). A new species of Desmalopex (Pteropodidae) from the Philippines, with a phylogenetic analysis of the Pteropodini. J. Mammal 89, 815–825. doi: 10.1644/07-MAMM-A-285.1

CrossRef Full Text | Google Scholar

Fisher M., Garner T. (2020). Chytrid fungi and global amphibian declines. Nat. Re. Microbiol. 18, 332–343. doi: 10.1038/s41579-020-0335-x

CrossRef Full Text | Google Scholar

Forson D. D., Storfer A. (2006). Atrazine increases ranavirus susceptibility in the tiger salamander, Ambystoma tigrinum. Ecol. App. 16, 2325–2332. doi: 10.1890/1051-0761(2006)016[2325:AIRSIT]2.0.CO;2

CrossRef Full Text | Google Scholar

Goslee S. C., Urban D. L. (2007). The ecodist package for dissimilarity-based analysis of ecological data. J. Stat. Soft. 22, 1–19. doi: 10.18637/jss.v022.i07

CrossRef Full Text | Google Scholar

Graham C. H., Storch D., Machac A., Isaac N. (2018). Phylogenetic scale in ecology and evolution. Global Ecol. Biogeogr 27, 175–187. doi: 10.1111/geb.12686

CrossRef Full Text | Google Scholar

Gray M. J., Brunner J. L., Earl J. E., Ariel E. (2015). “Design and analysis of ranavirus studies: surveillance and assessing risk,” in Ranaviruses. Eds. Gray M. J., Chinchar V. G. (Springer International Publishing, Cham). doi: 10.1007/978-3-319-13755-1

CrossRef Full Text | Google Scholar

Gray M. J., Miller D. L., Hoverman J. T. (2012). Reliability of nonlethal surveillance methods for detecting ranavirus infection. Dis. Aquat. Org. 99, 1–6. doi: 10.3354/dao02436

CrossRef Full Text | Google Scholar

Harris R. N., Brucker R. M., Walke J. B., Becker M. H., Schwantes C. R., Flaherty D. C., et al. (2009). Skin microbes on frogs prevent morbidity and mortality caused by a lethal skin fungus. ISME J. 3, 818–824. doi: 10.1038/ismej.2009.27

PubMed Abstract | CrossRef Full Text | Google Scholar

Harrison X. A., Price S. J., Hopkins K., Leung W. T. M., Sergeant C., Garner T. W. J. (2019). Diversity-stability dynamics of the amphibian skin microbiome and susceptibility to a lethal viral pathogen. Front. Microbiol. 10. doi: 10.3389/fmicb.2019.02883

CrossRef Full Text | Google Scholar

Hird S. M. (2017). Evolutionary biology needs wild microbiomes. Front. Microbiol. 8. doi: 10.3389/fmicb.2017.00725

CrossRef Full Text | Google Scholar

Hoagland B. W., Stoodley S. (1998). “Ecoregions,” in Riparian Area Management Handbook (Oklahoma Cooperative Extension Service, Oklahoma State University, Stillwater, Oklahoma), 5–18.

Google Scholar

Hooper L. V., Gordon J. I. (2001). Commensal host-bacterial relationships in the gut. Science 5519, 1115–1118. doi: 10.1126/science.1058709

CrossRef Full Text | Google Scholar

Hothorn T., Zeileis A. (2015). partykit: A modular toolkit for recursive partytioning in R. JMLR. 116, 3905–3909.

Google Scholar

Hughey M. C., Warne R., Dulmage A., Reeve R. E., Curtis G. H., Whitfield K., et al. (2023). Diet- and salinity- induced modifications of the gut microbiota are associated with differential physiological responses to ranavirus infection in Rana sylvatica. Phil. Trans. R. Soc B 378, 20220121. doi: 10.1098/rstb.2022.0121

CrossRef Full Text | Google Scholar

IUCN (2022) The IUCN Red List of Threatened Species. Version 2022-2 (Accessed January 17, 2022).

Google Scholar

Jani A. J., Briggs C. J. (2014). The pathogen Batrachochytrium dendrobatidis disturbs the frog skin microbiome during a natural epidemic and experimental infection. Proc. Natl. Acad. Sci. U.S.A. 111, E5049–E5058. doi: 10.1073/pnas.1412752111

PubMed Abstract | CrossRef Full Text | Google Scholar

Jani A. J., Bushell J., Arisdakessian C. G., Belcaid M., Boiano D. M., Brown C., et al. (2021). The amphibian microbiome exhibits poor resilience following pathogen-induced disturbance. ISME J. 15, 1628–1640. doi: 10.1038/s41396-020-00875-w

PubMed Abstract | CrossRef Full Text | Google Scholar

Jiménez R. R., Sommer S. (2017). The amphibian microbiome: natural range of variation, pathogenic dysbiosis, and role in conservation. Biodivers. Conserv. 26, 763–786. doi: 10.1007/s10531-016-1272-x

CrossRef Full Text | Google Scholar

Johnson K. V.-A., Burnet P. W. J. (2016). Microbiome: Should we diversify from diversity? Gut Microbes 7, 455–458. doi: 10.1080/19490976.2016.1241933

PubMed Abstract | CrossRef Full Text | Google Scholar

Johnston R., Jones K., Manley D. (2018). Confounding and collinearity in regression analysis: a cautionary tale and an alternative procedure, illustrated by studies of British voting behaviour. Qual. Quant. 52, 1957–1976. doi: 10.1007/s11135-017-0584-6

PubMed Abstract | CrossRef Full Text | Google Scholar

Katoh K. (2002). MAFFT: a novel method for rapid multiple sequence alignment based on fast Fourier transform. Nucleic Acids Res. 30, 3059–3066. doi: 10.1093/nar/gkf436

PubMed Abstract | CrossRef Full Text | Google Scholar

Kerby J. L., Schieffer A., Brown J. R., Whitfield S. (2013). Utilization of fast qPCR techniques to detect the amphibian chytrid fungus: a cheaper and more efficient alternative method. Methods Ecol. Evol. 4, 162–166. doi: 10.1111/j.2041-210x.2012.00263.x

CrossRef Full Text | Google Scholar

Kozich J. J., Westcott S. L., Baxter N. T., Highlander S. K., Schloss P. D. (2013). Development of a dual-index sequencing strategy and curation pipeline for analyzing amplicon sequence data on the MiSeq Illumina sequencing platform. Appl. Environ. Microbiol. 79, 5112–5120. doi: 10.1128/AEM.01043-13

PubMed Abstract | CrossRef Full Text | Google Scholar

Kueneman J. G., Parfrey L. W., Woodhams D. C., Archer H. M., Knight R., McKenzie V. J. (2014). The amphibian skin-associated microbiome across species, space and life history stages. Mol. Ecol. 23, 1238–1250. doi: 10.1111/mec.12510

PubMed Abstract | CrossRef Full Text | Google Scholar

Kueneman J. G., Woodhams D. C., Van Treuren W., Archer H. M., Knight R., McKenzie V. J. (2016). Inhibitory bacteria reduce fungi on early life stages of endangered Colorado boreal toads (Anaxyrus boreas). ISME J. 10, 934–944. doi: 10.1038/ismej.2015.168

PubMed Abstract | CrossRef Full Text | Google Scholar

Lannoo M. J., Petersen C., Lovich R. E., Nanjappa P., Phillips C., Mitchell J. C., et al. (2011). Do frogs get their kicks on Route 66? Continental U.S. transect reveals spatial and temporal patterns of Batrachochytrium dendrobatidis infection. PloS One 6, e22211. doi: 10.1371/journal.pone.0022211

PubMed Abstract | CrossRef Full Text | Google Scholar

Lim S. J., Bordenstein S. R. (2020). An introduction to phylosymbiosis. Proc. R. Soc B: Biol. Sci. 287, 20192900. doi: 10.1098/rspb.2019.2900

CrossRef Full Text | Google Scholar

Lin H., Peddada S. D. (2020). “Analysis of compositions of microbiomes with bias correction” Nat. Commun. 11, 1–11. doi: 10.1038/s41467-020-17041-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Lips K. R., Brem F., Brenes R., Reeve J. D., Alford R. A., Voyles J., et al. (2006). Emerging infectious disease and the loss of biodiversity in a Neotropical amphibian community. PNAS U.S.A. 103, 3165–3170. doi: 10.1073/pnas.0506889103

CrossRef Full Text | Google Scholar

Love M. I., Huber W., Anders S. (2014). Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550. doi: 10.1186/s13059-014-0550-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Lozupone C., Knight R. (2005). UniFrac: A new phylogenetic method for comparing microbial communities. Appl. Environ. Microbiol. 71, 8228–8235. doi: 10.1128/AEM.71.12.8228-8235.2005

PubMed Abstract | CrossRef Full Text | Google Scholar

Marhanka C. E., Watters L. J., Huron A. N., McMillin L. S., Winfrey C. C., Curtis J. D., et al. (2017). Detection of high prevalence of Batrachochytrium dendrobatidis in amphibians from southern Oklahoma, USA. Herpetol. Rev. 48, 71–74.

Google Scholar

Martel A., Blooi M., Adriaensen C., Van Rooij P., Beukema W., Fisher M. C., et al. (2014). Recent introduction of a chytrid fungus endangers Western Palearctic salamanders. Science 346, 630–631. doi: 10.1126/science.1258268

PubMed Abstract | CrossRef Full Text | Google Scholar

Martel A., Spitzen-van der Sluijs A., Blooi M., Bert W., Ducatelle R., Fisher M. C., et al. (2013). Batrachochytrium salamandrivorans sp. nov. causes lethal chytridiomycosis in amphibians. Proc. Natl. Acad. Sci. U.S.A. 110, 15325–15329. doi: 10.1073/pnas.1307356110

PubMed Abstract | CrossRef Full Text | Google Scholar

Martin S. D., Shepard D. B., Steffen M. A., Phillips J. G., Bonett R. M. (2016). Biogeography and colonization history of plethodontid salamanders from the Interior Highlands of eastern North America. J. Biogeogr. 43, 410–422. doi: 10.1111/jbi.12625

CrossRef Full Text | Google Scholar

McKenzie V. J., Bowers R. M., Fierer N., Knight R., Lauber C. L. (2012). Co-habiting amphibian species harbor unique skin bacterial communities in wild populations. ISME J. 6, 588–596. doi: 10.1038/ismej.2011.129

PubMed Abstract | CrossRef Full Text | Google Scholar

McMurdie P. J., Holmes S. (2013). phyloseq: An R package for reproducible interactive analysis and graphics of microbiome census data. PloS One 8, e61217. doi: 10.1371/journal.pone.0061217

PubMed Abstract | CrossRef Full Text | Google Scholar

Medina D., Hughey M. C., Becker M. H., Walke J. B., Umile T. P., Burzynski E. A., et al. (2017). Variation in metabolite profiles of amphibian skin bacterial communities across elevations in the Neotropics. Microb. Ecol. 74, 227–238. doi: 10.1007/s00248-017-0933-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Moloney R. D., Desbonnet L., Clarke G., Dinan T. G., Cryan J. F. (2014). The microbiome: stress, health, and disease. Mamm. Genome 25, 49–74. doi: 10.1007/s00335-013-9488-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Muletz-Wolz C. R., Fleischer R. C., Lips K. R. (2019). Fungal disease and temperature alter skin microbiome structure in an experimental salamander system. Mol. Ecol. 15122, 2917–2931. doi: 10.1111/mec.15122

CrossRef Full Text | Google Scholar

Muletz Wolz C. R., Yarwood S. A., Campbell Grant E. H., Fleischer R. C., Lips K. R. (2018). Effects of host species and environment on the skin microbiome of plethodontid salamanders. J. Anim. Ecol. 87, 341–353. doi: 10.1111/1365-2656.12726

PubMed Abstract | CrossRef Full Text | Google Scholar

Mutnale M. C., Reddy G. S., Vasudevan K. (2021). Bacterial community in the skin microbiome of frogs in a coldspot of chytridiomycosis infection. Microb. Ecol. 82, 554–558. doi: 10.1007/s00248-020-01669-5

PubMed Abstract | CrossRef Full Text | Google Scholar

Nava-González B., Suazo-Ortuño I., López P. B., Maldonado-López Y., Lopez-Toledo L., Raggi L., et al. (2021). Inhibition of Batrachochytrium dendrobatidis infection by skin bacterial communities in wild amphibian populations. Microb. Ecol. 82, 666–676. doi: 10.1007/s00248-021-01706-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Oksanen J., Simpson G., Blanchet F., Kindt R., Legendre P., Weedon J. (2022). egan: community ecology package. R package version 2.6-4. Available at: https://CRAN.R-project.org/package=vegan.

Google Scholar

O’Meara B. C., Harmon L., Eastman J. (2013). PhyloOrchard: Important and/or useful phylogenetic datasets.

Google Scholar

Pennell M. W., Eastman J. M., Slater G. J., Brown J. W., Uyeda J. C., FitzJohn R. G., et al. (2014). geiger v2.0: An expanded suite of methods for fitting macroevolutionary models to phylogenetic trees. Bioinformatics 30, 2216–2218. doi: 10.1093/bioinformatics/btu181

PubMed Abstract | CrossRef Full Text | Google Scholar

Petranka J. W., Scott D. E. (1998). Salamanders of the United States and Canada. 3rd ed. (Smithsonian Institution) doi: 10.2307/1447629

CrossRef Full Text | Google Scholar

Prado-Irwin S. R., Bird A. K., Zink A. G., Vredenburg V. T. (2017). Intraspecific variation in the skin-associated microbiome of a terrestrial salamander. Microb. Ecol. 74, 745–756. doi: 10.1007/s00248-017-0986-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Price M. N., Dehal P. S., Arkin A. P. (2010). FastTree 2 – Approximately maximum-likelihood trees for large alignments. PloS One 5, e9490. doi: 10.1371/journal.pone.0009490

PubMed Abstract | CrossRef Full Text | Google Scholar

Pyron R. A., Wiens J. J. (2011). A large-scale phylogeny of Amphibia including over 2800 species, and a revised classification of extant frogs, salamanders, and caecilians. Mol. Phylogenet. Evol. 61, 543–583. doi: 10.1016/j.ympev.2011.06.012

PubMed Abstract | CrossRef Full Text | Google Scholar

R Core Team (2022). R: A language and environment for statistical computing (Vienna: R Foundation for Statistical Computing).

Google Scholar

Rebollar E. A., Hughey M. C., Medina D., Harris R. N., Ibáñez R., Belden L. K. (2016). Skin bacterial diversity of Panamanian frogs is associated with host susceptibility and presence of Batrachochytrium dendrobatidis. ISME J. 10, 1682–1695. doi: 10.1038/ismej.2015.234

PubMed Abstract | CrossRef Full Text | Google Scholar

Rebollar E. A., Martínez-Ugalde E., Orta A. H. (2020). The amphibian skin microbiome and its protective role against chytridiomycosis. Herpetologica 76, 167–177. doi: 10.1655/0018-0831-76.2.167

CrossRef Full Text | Google Scholar

Rognes T., Flouri T., Nichols B., Quince C., Mahé F. (2016). VSEARCH: a versatile open source tool for metagenomics. PeerJ 18. doi: 10.7717/peerj.2584

CrossRef Full Text | Google Scholar

Rollins-Smith L. A. (2009). The role of amphibian antimicrobial peptides in protection of amphibians from pathogens linked to global amphibian declines. Biochim. Biophys. Acta 1788, 1593–1599. doi: 10.1016/j.bbamem.2009.03.008

PubMed Abstract | CrossRef Full Text | Google Scholar

Rollins-Smith L. A., Ramsey J. P., Pask J. D., Reinert L. K., Woodhams D. C. (2011). Amphibian immune defenses against chytridiomycosis: Impacts of changing environments. ICB 51, 552–562. doi: 10.1093/icb/icr095

CrossRef Full Text | Google Scholar

Ruthsatz K., Lyra M. L., Lambertini C., Belasen A. M., Jenkinson T. S., da Silva Leite D., et al. (2020). Skin microbiome correlates with bioclimate and Batrachochytrium dendrobatidis infection intensity in Brazil’s Atlantic Forest treefrogs. Sci. Rep. 10, 22311. doi: 10.1038/s41598-020-79130-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Scheele B. C., Skerratt L. F., Grogan L. F., Hunter D. A., Clemann N., McFadden M., et al. (2017). After the epidemic: Ongoing declines, stabilizations and recoveries in amphibians afflicted by chytridiomycosis. Biol. Conserv. 206, 37–46. doi: 10.1016/j.biocon.2016.12.010

CrossRef Full Text | Google Scholar

Schubert M., Lindgreen S., Orlando L. (2016). AdapterRemoval v2: Rapid adapter trimming, identification, and read merging. BMC Res. Notes 9, 88. doi: 10.1186/s13104-016-1900-2

PubMed Abstract | CrossRef Full Text | Google Scholar

Skerratt L. F., Berger L., Speare R., Cashins S., McDonald K. R., Phillott A. D., et al. (2007). Spread of chytridiomycosis has caused the rapid global decline and extinction of frogs. EcoHealth 4, 125. doi: 10.1007/s10393-007-0093-5

CrossRef Full Text | Google Scholar

Smith H. K., Pasmans F., Dhaenens M., Deforce D., Bonte D., Verheyen K., et al. (2018). Skin mucosome activity as an indicator of Batrachochytrium salamandrivorans susceptibility in salamandersSkin mucosome activity as an indicator of Batrachochytrium salamandrivorans susceptibility in salamanders. PLOS One 13, e0199295.

PubMed Abstract | Google Scholar

Smith S. N., Watters J. L., Marhanka E. C., McMillin S. L., Davis D. R., Farkas J. K., et al. (2019). Investigating ranavirus prevalence in Central Oklahoma, USA, Amphibians. Herpetol. Rev. 50, 508–512.

Google Scholar

Smith S. N., Watters J. L., Siler C. D. (2023). Host ecology drives frog skin microbiome diversity across ecotone in South-Central North America. Front. Microbiomes 2. doi: 10.3389/frmbi.2023.1286985

CrossRef Full Text | Google Scholar

St-Amour V., Lesbarrères D. (2007). Genetic evidence of Ranavirus in toe clips: an alternative to lethal sampling methods. Conserv. Genet. 8, 1247–1250. doi: 10.1007/s10592-006-9242-6

CrossRef Full Text | Google Scholar

Trauth S. E., Robison H. W., Plummer M. V. (2004). The Amphibians and Reptiles of Arkansas (Fayetteville, AR, USA: The University of Arkansas Press).

Google Scholar

U. S. Geological Survey (2017) The National Map seamless digital elevation model specifications (Accessed October 3, 2022).

Google Scholar

Voyles J., Berger L., Young S., Speare R., Webb R., Warner J., et al. (2007). Electrolyte depletion and osmotic imbalance in amphibians with chytridiomycosis. Dis. Aquat. Org. 77, 113–118. doi: 10.3354/dao01838

CrossRef Full Text | Google Scholar

Voyles J., Young S., Berger L., Campbell C., Voyles W. F., Dinudom A., et al. (2009). Pathogenesis of chytridiomycosis, a cause of catastrophic amphibian declines. Science 326, 582–585. doi: 10.1126/science.1176765

PubMed Abstract | CrossRef Full Text | Google Scholar

Walker D. M., Hill A. J., Albecker M. A., McCoy M. W., Grisnik M., Romer A., et al. (2020). Variation in the slimy salamander (Plethodon spp.) skin and gut-microbial assemblages is explained by geographic distance and host affinity. Microb. Ecol. 79, 985–997. doi: 10.1007/s00248-019-01456-x

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang T., Hamann A., Spittlehouse D., Carroll C. (2016). Locally downscaled and spatially customizable climate data for historical and future periods for North America. PloS One 11. doi: 10.1371/journal.pone.0156720

CrossRef Full Text | Google Scholar

Watters J. L., Davis D. R., Yuri T., Siler C. D. (2018). Concurrent infection of Batrachochytrium dendrobatidis and ranavirus among native amphibians from northeastern Oklahoma, USA. J. Aquat. Anim. Health 30, 291–301. doi: 10.1002/aah.10041

PubMed Abstract | CrossRef Full Text | Google Scholar

Watters J. L., McMillin S. L., Marhanka E. C., Davis D. R. (2019). Seasonality in Batrachochytrium dendrobatidis detection in amphibians in central Oklahoma, USA. JZWM 50, 492. doi: 10.1638/2018-0103

CrossRef Full Text | Google Scholar

Weiss S., Xu Z. Z., Peddada S., Amir A., Bittinger K., Gonzalez A., et al. (2017). Normalization and microbial differential abundance strategies depend upon data characteristics. Microbiome 5, 27. doi: 10.1186/s40168-017-0237-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Wilber M. Q., Jani A. J., Mihaljevic J. R., Briggs C. J. (2020). Fungal infection alters the selection, dispersal and drift processes structuring the amphibian skin microbiome. Ecol. Lett. 23, 88–98. doi: 10.1111/ele.13414

PubMed Abstract | CrossRef Full Text | Google Scholar

Woodhams D. C., Alford R. A., Antwis R. E., Archer H., Becker M. H., Belden L. K., et al. (2015). Antifungal isolates database of amphibian skin-associated bacteria and function against emerging fungal pathogens. Ecology 96, 595–595. doi: 10.1890/14-1837.1

CrossRef Full Text | Google Scholar

Xu L., Xiang M., Zhu W., Zhang M., Chen H., Huang J., et al. (2020). The behavior of amphibians shapes their symbiotic microbiomes. mSystems 5, e00626-20. doi: 10.1128/mSystems.00626-20

PubMed Abstract | CrossRef Full Text | Google Scholar

Yeoh Y. K., Dennis P. G., Paungfoo-Lonhienne C., Weber L., Brackin R., Ragan M. A., et al. (2017). Evolutionary conservation of a core root microbiome across plant phyla along a tropical soil chronosequence. Nat. Commun. 8, 215. doi: 10.1038/s41467-017-00262-8

PubMed Abstract | CrossRef Full Text | Google Scholar

Yoon S.-H., Ha S.-M., Kwon S., Lim J., Kim Y., Seo H., et al. (2017). Introducing EzBioCloud: A taxonomically united database of 16S rRNA gene sequences and whole-genome assemblies. IJSEM 67, 1613–1617. doi: 10.1099/ijsem.0.001755

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang P., Wake D. B. (2009). Higher-level salamander relationships and divergence dates inferred from complete mitochondrial genomes. Mol. Phylogenet. Evol. 53, 492–508. doi: 10.1016/j.ympev.2009.07.010

PubMed Abstract | CrossRef Full Text | Google Scholar