Ecological patterns in the skin microbiota of frogs from tropical Australia

Abstract The microbiota of frog skin can play an important role in protecting against diseases and parasites. The frog skin microbial community represents a complex mix of microbes that are promoted by the chemical environment of the frog skin and influenced by the animal's immediate past environment. The microbial communities of six species of frogs sampled from the campus of Charles Darwin University (CDU) were more similar within species than between species. The microbiota of the introduced cane toad (Rhinella marina) was most dissimilar among the species. Pairwise comparisons showed that the microbial communities of each species were different, except for the terrestrial Litoria nasuta and the arboreal L. rothii. The microbial communities of the six species were not related to ecological habit (arboreal or terrestrial), and neither was the alpha diversity of the microbes. The core microbes (defined as being on ≥90% of individuals of a species or group) were significantly different among all species, although 89 microbial operational taxonomic units (OTUs) were core microbes for all six species at CDU. Two species, Rhinella marina and Litoria rothii, were sampled at additional sites approximately 10 and 30 km from CDU. The microbial communities and the core OTU composition were different among the sites, but there were nevertheless 194 (R. marina) and 181 (L. rothii) core OTUs present at all three sites. Thus, the core microbiota varied with respect to geographic range and sample size.


Biologic factors include transitions in life stages
Despite the perturbations resulting from these environmental and biologic factors, species-specific communities persist and are consistent in different locations (McKenzie et al., 2012;Kueneman et al. 2013;Walke et al., 2014;Belden et al., 2015) and prolonged captivity . This apparent homeostasis is consistent with the notion that several characteristics of the amphibian skin (epidermal structures, skin peptides, and other mucosal components) select and enhance specific bacteria over others (Rollins-Smith & Woodhams, 2012). This dynamic process between a disruptive environment and the homeostatic properties of the skin can explain the seemingly contradictory reports of both variability and consistency of the microbiota of amphibian skins (Kueneman et al., 2013;Adair & Douglas, 2017).
Detailed exploration of host-microbiota associations has resulted in the concept of a "core" community, based on near ubiquity among individuals of a host species. The core microbiota has been variously defined, but is typically defined as being present on >90% of hosts (Loudon et al., 2014;Apprill et al., 2014).
The aim of this study was to enhance our understanding of the ecological context of the skin microbiota by sampling from a previously well-studied group of frogs Young, Christian, Donnellan, Tracy, & Parry, 2005;Tracy, Christian, Betts, & Tracy, 2008;Tracy, Tixier, Le Nöene, & Christian, 2014) from the wet-dry tropics of Australia. We studied six frog species from one site, three of which are arboreal species and three of which are terrestrial. Previous studies of this community of frogs have shown that ecological habit (arboreal versus terrestrial) is associated with some physiological traits in this seasonal tropical environment (see references cited above). To explore the effect of host species, ecological habit, and geographic location on skin microbial community patterns, we further sampled two of the six species from two additional sites at distances of approximately 10 km and 30 km from the main study site.

| Species, study sites, and sampling scheme
Approval to sample frogs was granted by the Charles Darwin University Animal Ethics Committee (project A14012). Three terrestrial frogs (Rhinella marina, Litoria nasuta, and Limnodynastes convexiusculus) and three arboreal frogs (Litoria caerulea (Figure 1), Litoria rubella, and Litoria rothii) species were sampled from on or near the campus of Charles Darwin University (CDU), in Darwin, Northern Territory, Australia. Although some frogs were collected near buildings on campus, others were sampled from a more natural area consisting of native vegetation. Two species, the cane toad, Rhinella marina, and Roth's tree frog, Litoria rothii, were also sampled from Mickett Creek and Howard River, approximately 10 and 30 km, respectively, from CDU. These rural areas primarily consist of native vegetation, with widely dispersed buildings. The cane toad was introduced to Australia approximately 80 years ago (Easteal, 1981), but the other species are native to the area. Twenty individual frogs were sampled from each species and location (total of 200 samples).
New gloves were used for each sample.
After capture by hand, each frog was rinsed twice with 100 mL 0.45 μm filtered high-purity water (Culp, Falkinham, & Belden, 2007;Lauer et al., 2007) before being swabbed with a sterile synthetic swab (MicroRheologics FLOQSwab). Each frog was stroked 30 times to produce a sample, using 10 strokes around body (avoiding cloaca, 4 dorsal, 1 each side, 4 ventral) and 5 strokes on each limb (front F I G U R E 1 Two green tree frogs, Litoria caerulea, on the campus of Charles Darwin University, Darwin, Northern Territory, Australia and back of the foot, front and back of the leg, axial region). Samples were kept on ice while in the field and then frozen at −20°C until the DNA was extracted.

| DNA extraction
DNA was extracted from the samples using the Qiagen DNeasy Blood and Tissue Kit (Qiagen, Valencia, CA), following the manufacturer's protocol for tissue extraction.

| 16S rDNA target sequencing
Two hundred nanograms of DNA was sent to the sequenc-

| Processing of sequencing data
Sequence data were processed using a proprietary analysis pipeline (www.mrdnalab.com, Molecular Research DNA, Shallowater, TX).
Sequences were depleted of barcodes, and primers and short sequences <200 bp were removed as well as sequences with ambiguous base calls and homopolymer runs exceeding 6nt. Chimeras were also removed. Sequences were clustered at 3% divergence (97% similarity) to define operational taxonomic units (OTUs), and OTUs with singleton sequences were removed (Dowd et al., 2008;Edgar, 2010;Swanson et al., 2011;Capone, Dowd, Stamatas, & Nikolovski, 2011). OTUs were taxonomically classified using BLASTn against a curated Greengenes database (DeSantis et al., 2006). Furthermore, OTUs were excluded which were not classified as bacteria, occurred in less than 1% of samples (i.e., in less than three samples), or contained fewer sequences than 0.01% of the total sequence abundance. All sequences were subsampled to the lowest common sequence number (2,008 sequences) per sample.
The dataset was subset into two groups: Frogs sampled on or near the CDU campus and the cane toads and Roth's tree frogs from three sites (CDU, Mickett Creek, and Howard River). Bacterial orders that occurred at relative abundances of more than 1% in a frog species or site are shown in taxa plots using phyloseq. A weighted UniFrac distance matrix was created. The distance matrix was visualized using nonmetric multidimensional scaling (nMDS) and a triangle heat map. Alpha diversity and changes in the frog skin microbial communities were analyzed in Primer-7 by permutational MANOVA (PERMANOVA; 9999 permutations) with fixed factors "frog species" and "habitat" for the CDU frog dataset and "frog species" and "site" for the toad and Roth's tree frog dataset. Pairwise tests for frog species, habitat, and site were also conducted if the main test was significant (p < 0.05). Alpha diversity was examined with respect to ecological habit using a linear mixed-effect model in Stata with ecological habit as a fixed effect and frog species as a random effect nested in ecological habit.
We identified the core microbiota for campus frogs, and toad and Roth's tree frogs sampled at the campus, Howard River, and Mickett Creek sites by selecting those OTUs that were present in at least 90% of samples within a group. Core OTUs are shown in Venn diagrams, and PERMANOVA and pairwise analyses were also conducted on these core OTUs to explore significant differences in the core microbiota between frog species and sites.

| Microbial community composition
After data processing, 651 bacterial OTUs were recorded from the six frog species sampled. The composition of frog skin bacteria at the level of order is shown in Figure 2a The alpha diversity, as represented by the Shannon index, is given in Table 1 for both the total complement of OTUs and the core OTUs for the six species on the CDU campus. The diversity of OTUs was significantly higher in Limnodynastes convexiusculus (4.97) than in the other frog species (F 5,119 = 31.3, p < 0.0001), and this was also true for the diversity of core OTUs. The lowest diversity (3.10) was found in Litoria caerulea, and this species also had the lowest diversity of core OTUs. Once the effect of species was taken into account, there was no significant effect of ecological habit on alpha diversity (p = 0.24).

| Frog species comparisons
An nMDS showed clustering of the microbial composition by frog species (Figure 3a), and PERMANOVA confirmed the clustering of communities according to species (F 5,119 = 36.1, p < 0.001). The microbial communities of the frog species were significantly different from one another in pairwise comparisons (p < 0.001), except not between the terrestrial Litoria nasuta and the arboreal L. rothii.
A triangle heat map of average weighted UniFrac distances between frog species illustrates that the skin microbiota of the introduced cane toads was the most dissimilar when compared to the five native species, and communities were more similar within species than between species (Figure 4). Microbial communities of the three terrestrial species and the three arboreal species were significantly different between the two ecological habits (F 1,119 = 12.5,

| Core OTUs
The number of core OTUs for each frog species is listed in Table 1, and the core OTUs are also expressed as a percentage of the total OTUs for each frog species. Overall, 604 OTUs were core OTUs for at least one frog species. Over 70% of the total OTUs in Limnodynastes convexiusculus were core OTUs, which was a distinctly higher percentage than in the other frogs. Litoria rubella had the lowest percentage as core OTUs (36.6%) among the six species. The composition of the core OTUs was significantly different between frog species at CDU (F 5,1452 = 164.1, with p < 0.001 for all pairwise comparisons), although 89 core OTUs were shared by all six species on the CDU campus. The core OTUs for R. marina and L. rothii were significantly different across three sites (F 2,509 = 9.3, p < 0.0007), but

| D ISCUSS I ON
The six frog species studied had species-specific microbial communities on their skin, despite those at CDU being from the same location. Strikingly, L. caerulea and L. rubella, many of which were particularly those of the genus Pseudomonas, contain several species known to provide protection against pathogenic bacteria and fungi (Chang, Huang, Lin, Huang, & Liao, 2016;Federici et al., 2015). The Enterobacteriales were dominant in this study, and while this order has been reported from the guts of frogs (Chang et al., 2016), we found only two reports from amphibian (toad) skin studies (Hughey et al., 2017;Bataille et al. 2016). Enterobacteria can be found almost everywhere (soil, water, wastewater, animal guts), so frogs could come into contact with them from drains and water bodies (Neave et al., 2014) or via invertebrates they have consumed (Chang et al., 2016).
Taken as a group, the communities of arboreal frogs were significantly different from terrestrial frogs, as are a range of physiological characteristics (Young et al., 2005;Tracy et al., 2014). Nevertheless, the microbial community of the terrestrial frog L. nasuta was not different from that of the arboreal L. rothii, suggesting that the significant difference between the two ecological habits was, in part, due to the general pattern of frog species-specific microbial communities. The similarity between L. nasuta and L. rothii cannot be explained by phylogenetic relationships because the arboreal species L. rothii and L. rubella are closely related, but L. nasuta is more distantly related and groups with other terrestrial species in the genus (Young et al., 2005).
The number of core OTUs for the six frog species on the CDU campus ranged from 256 to 281 (Table 1), and these represented 36.6% to 70.7% of the total microbiota for the frog species, with the core OTUs of Limnodynastes convexiusculus being the highest percentage of the total. For the two species sampled across three sites, the core was variable among the sites, yet there were nevertheless OTUs that were present on all individuals over the 30-km TA B L E 1 (A) The total number of OTUs and core OTUs measured for each of six species sampled from the campus of Charles Darwin University (CDU) and for two species at two additional sites. The %Core is defined as the number of core OTUs expressed as a percentage of the total number of OTUs. The Shannon diversity index is used to represent the total diversity and core diversity for the frogs at CDU. (B) The core is compared using data from three sites for two species. The combined core treats the data from three sites as one group. Thus, combined core OTUs were found in at least 90% of all individuals of a species as opposed to 90% of individuals from a given site F I G U R E 3 Relatedness of the skin microbiota as shown by a nonmetric multidimensional scaling (nMDS) of (a) six frog species sampled from the campus of Charles Darwin University (CDU) and (b) two species (L. rothii and R. marina) sampled from CDU, Mickett Creek (MC), and Howard River (HR). Each nMDS was based on the weighted UniFrac distance matrix of rarefied OTU data and had a stress value <0.14 (a) and <0.10 (b) transect ( Figure 5). If the data from three sites are grouped within a species, the number of OTUs that meet the definition of "core" is substantially smaller than the number of core OTUs from a single site. This raises questions about the definition and determination of the core microbiota. One could argue that the best measurement of the core microbes is represented by the central portion of the Venn diagrams in Figure 5, which is not only based on a larger sample (n = 60), but the larger geographic range of the sample also likely provides a more comprehensive characterization of the microbiota.
However, even if one accepts this argument over a 30-km transect, the definition and determination of the core microbiota become problematic over greater distances because both the environment and the frogs themselves could vary substantially at a larger (i.e., continental) scale. Thus, the optimal microbiota at one site may not be the same at another site with frogs and microbes adapted to each local environment.
The core microbiota is defined on the basis of prevalence (≥90% of individuals), but the unstated assumption is that if it is found in most individuals, then a core microbe is likely to have a role in the microbiota in producing important metabolites that are either important to the frog (i.e., antifungal properties or activation of the immune system ) or in structuring the microbiota. However, near ubiquity in the skin microbiota (i.e., a core microbe) could reflect either functional importance or simply prevalence in the environment. However, the fact that different frog species from the same location have different cores provides indirect support for the notion that the core microbiota is not simply reflecting the ubiquity of microbes in the environment. Analyses that indicate microbial function in the microbiota (Berry & Widder, 2014) are more appropriate for determining the importance of an OTU (i.e., a "keystone") in a microbiota. The relationship between core and keystone OTUs has received little attention. Figure 5 indicates that a core microbe (as determined at a site level) may not be functionally important (given that it is not part of the core at sites 10 or 30 km away), but rather, some environmental microbes may be abundant at some sites but not others. However, is the converse true? That is: Is it possible for a keystone microbe to not be a core microbe?
Understanding the functional importance of microbes within the microbiota is crucial to advancing our understanding of these relationships. Although descriptively satisfying, the concept of a core (as currently defined) may not continue to be a useful construct as our understanding of the structure and function of the frog skin microbiota develops.

ACK N OWLED G M ENTS
Financial support was provided from the Faculty of Engineering, F I G U R E 4 Heat map triangle showing the average weighted UniFrac distances (dissimilarities of the OTU data between species). The circles on the diagonal show the within-species dissimilarities from the 20 individuals of each of the six frog species. Larger circles and numbers represent greater dissimilarities within and between species

CO N FLI C T O F I NTE R E S T
None declared.

AUTH O R CO NTR I B UTI O N S
K.C., K.G., and C.W. designed the project, K.C. and C.W. collected the samples in the field, C.W. extracted the DNA, M.K. and A.R.
analyzed the data and produced the figures, and K.C., with assistance from all other authors, wrote and revised the manuscript.