Skip to main content

Molecular and culture-based assessment of the microbiome in a zebrafish (Danio rerio) housing system during set-up and equilibration



Zebrafish used in research settings are often housed in recirculating aquaculture systems (RAS) which rely on the system microbiome, typically enriched in a biofiltration substrate, to remove the harmful ammonia generated by fish via oxidation. Commercial RAS must be allowed to equilibrate following installation, before fish can be introduced. There is little information available regarding the bacterial community structure in commercial zebrafish housing systems, or the time-point at which the system or biofilter reaches a microbiological equilibrium in RAS in general.


A zebrafish housing system was monitored at multiple different system sites including tank water in six different tanks, pre- and post-particulate filter water, the fluidized bed biofilter substrate, post-carbon filter water, and water leaving the ultra-violet (UV) disinfection unit and entering the tanks. All of these samples were collected in quadruplicate, from prior to population of the system with zebrafish through 18 weeks post-population, and analyzed using both 16S rRNA amplicon sequencing and culture using multiple agars and annotation of isolates via matrix-assisted laser desorption/ionization-time-of-flight (MALDI-TOF) mass spectrometry. Sequencing data were analyzed using traditional methods, network analyses of longitudinal data, and integration of culture and sequence data.


The water microbiome, dominated by Cutibacterium and Staphylococcus spp., reached a relatively stable richness and composition by approximately three to four weeks post-population, but continued to evolve in composition throughout the study duration. The microbiomes of the fluidized bed biofilter and water leaving the UV disinfection unit were distinct from water at all other sites. Core taxa detected using molecular methods comprised 36 amplicon sequence variants, 15 of which represented Proteobacteria including multiple members of the families Burkholderiaceae and Sphingomonadaceae. Culture-based screening yielded 36 distinct isolates, and showed moderate agreement with sequencing data.


The microbiome of commercial RAS used for research zebrafish reaches a relatively stable state by four weeks post-population and would be expected to be suitable for experimental use following that time-point.


The laboratory zebrafish (Danio rerio) are an attractive model species in biomedical and developmental research, owing to a unique combination of qualities that are important for animal model selection including biological characteristics that can be exploited, the availability of advanced imaging and molecular techniques, and financial feasibility [1]. From a biological perspective, zebrafish share many of the same organ systems as higher vertebrates and, in some situations, may represent the optimal model species for translational research [2, 3]. Moreover, the genetic tractability of zebrafish has allowed the targeted testing of gene function as is performed in knock-out and transgenic rodents [4]. From a purely logistical perspective, zebrafish provide several conveniences such as a reduced regulatory burden compared to research rodents, reduced costs associated with housing and husbandry, and the ability to increase throughput in large-scale surveys of compounds.

Zebrafish are also increasingly used in behavioral and neuropsychological research, both in academia and the industrial sector during the pursuit and development of investigational new drugs (INDs) [5, 6]. With comprehensive and well-characterized behavioral ethograms [7, 8], the effects of various stimuli can be assessed [9,10,11,12]. Recent reviews highlight the diversity of behavioral traits and social deficits that have been investigated using both adult and larval zebrafish models [13,14,15].

During the last decade, there has been a growing appreciation of the influence of the gut microbiota (GM) on human and animal health, and in animal models. Specifically, the GM has been identified as a potential source of poor experimental reproducibility in animal models [16], while also providing discoveries regarding the mechanisms by which the GM influences host development and health, both physical and mental [17, 18]. The fact that zebrafish are often housed in aquaria filled with a recirculated filtered water supply brings into question the nature of bacterial populations present at different levels of commercial zebrafish housing systems. Moreover, the establishment of environmental bacterial populations in a recirculating aquaculture system (RAS) is essential to accommodate the removal of nitrogenous waste (i.e., ammonia), the accumulation of which is highly toxic to fish [19, 20]. While it is recommended to allow new zebrafish housing systems to equilibrate for one to two weeks prior to population, and to populate any aquaculture system gradually [21], there are minimal empirical data documenting the time-course at which the microbiological communities within commercial RAS reach a steady state. As the environmental microbiota may influence zebrafish physiology or host-associated microbiota [22], this information is important to ensure robust, reproducible data from zebrafish experiments.

In recent years, there has been an increasing awareness of a reproducibility crisis affecting in vivo experiments. As in other model organisms, lack of reproducibility in zebrafish experiments can be caused by differences in intrinsic factors, such as background genetics [23], or extrinsic factors, which include diet, housing systems and various other aspects of the environment [24]. Environmental microbiota impact human microbiomes and human health [25]. Thus, it is likely that environmental microbiota impact zebrafish microbiomes and zebrafish health, and by extension, experimental reproducibility. To date, the focus on the environment from a reproducibility perspective has been limited to the impact of abiotic factors such as water quality, environmental enrichment, structural complexity, and social factors such as stocking density [26].

The objectives of the current study were therefore to characterize the composition and core taxa present throughout a commercial zebrafish housing system during set-up and population, to characterize the microbial interaction networks during this initial period within the individual system sites and as a whole, and to identify the culturable portion of the system during that time. We therefore sampled six different sites of a new zebrafish housing system at twelve weekly or biweekly time-points, beginning immediately prior to population of several tanks with adult zebrafish (one week post-installation and introduction of circulating water) and ending at 18 weeks post-population. Samples were then subjected to 16S rRNA amplicon sequencing and a culture-based survey of bacterial communities in each system site throughout the study period.


Experimental design and sample collection

The housing system used in the current study was a ZS560 system from Aquaneering, Incorporated (San Diego, CA), equipped with a 6 GPM self-contained filtration system including a reusable solid particle filter, a fluidized bed biological filter, dual carbon filters, and a standard 40-W ultraviolet disinfection unit. Following commercial installation and set-up of the system, deionized (DI) water was introduced and the system was allowed to circulate unpopulated for one week. Initial baseline samples were then collected at time-point 1 (TP1), immediately prior to population of six tanks (labeled A through F) with 6 to 11 adult zebrafish per tank. Samples (1 mL of water or substrate) were collected from six different sites on the system including 1) post-UV disinfection prior to entering the housing tanks, 2) housing tanks, 3) pre-particulate filtration, 4) post-particulate filtration, 5) the fluidized bed biofilter (FBB), and 6) post-Carbon filtration (Fig. 1). Additionally, tank water samples were collected from three previously established tanks from a housing system constructed according to the design of Kim et al. [27], that had housed zebrafish used to populate the new system (three tanks, four replicates/tank). Subsequent samples were then collected from the six aforementioned sites weekly for one month (TP1 through TP4), and then every two weeks for the following three months (TP5 through TP12). Four technical replicates were collected at each site and time-point; additionally, six housing tanks were sampled at each time-point (with four technical replicates per tank).

Fig. 1
figure 1

Schematic diagram of the recirculating system, and the sampled sites including post-UV disinfection water (1), tank water from six different tanks in which fish were introduced following collection of baseline samples, denoted by X (2), pre-particulate filter water (3), post-particulate filter water (4), the fluidized bed biofilter (FBB) (5), and carbon filter water (6) (A), and timeline showing the time-points (TP) at which samples were collected from each site, in quadruplicate (B)

Animals and husbandry

All zebrafish were maintained in an AAALAC International-accredited facility at the University of Missouri (Columbia, MO) in accordance with the guidelines presented in the Guide for the Care and Use of Laboratory Animals [28]. All husbandry procedures were approved by the University of Missouri Animal Care and Use Committee. The zebrafish were 61 adult, mixed-sex, wild-type fish originally obtained from Aquatica BioTech (Sun City Center, FL). Prior to stocking on the ZS560 system (Aquaneering, Incorporated), zebrafish were housed on an established RAS assembled in-housed based on a previously published design [27]. Zebrafish were maintained on a 14:10-h light:dark cycle at 27 °C. Water quality was tested weekly using Lifegard Test Strips Aquatics 6 way All Purpose Test Kit per the manufacturer instructions (Aquaneering Inc., catalog number STK6), with the following limits of detection: nitrite, 0 to 10 ppm (mg/L); nitrate, 0 to 200 ppm; total hardness, 0 to 300 ppm; total alkalinity, 0 to 300 ppm; and pH, 6.2 to 8.4. Water quality parameters were maintained as follows: pH, 8.08 (at start-up)–7.27; total ammonia nitrogen, 0 ppm; nitrite, 0 ppm; nitrate < 20 ppm; alkalinity 40–80 ppm; and hardness, 80 ppm. System pH was regulated by the addition of salt and sodium bicarbonate to water via an automatic dosing system. Zebrafish were fed once daily with a commercially available feed (TetraMin Plus tropical flake food).

DNA extraction

DNA was extracted from 750 µL of each sample (with no prior filtration or concentration) using QIAamp® DNA PowerFecal® kits (Qiagen®), according to the manufacturer’s instructions with the exception that, rather than performing the initial homogenization of samples using the vortex adapter described in the protocol, samples were homogenized in the provided bead tubes using a TissueLyser II (Qiagen®) for three minutes at 30 Hz/sec, before proceeding according to the protocol and eluting in 100 µL of elution buffer (Qiagen®). DNA yields were quantified via fluorometry (Qubit® 2.0, Invitrogen™) using quant-iT™ BR dsDNA reagent kits (Thermo Fisher Scientific). Routine negative controls consisting of unused reagents reproducibly yield between 0 and 100 sequences. Positive controls for DNA extraction and sequencing consisted of mock community standards containing 10 different microbial taxa (ZymoBIOMICS™, D6300), all of which were detected in the resulting sequencing data with zero contaminating sequences.

16S rRNA library preparation and sequencing

Extracted DNA was processed at the University of Missouri DNA Core Facility (Columbia, MO). Bacterial 16S rRNA amplicon libraries were generated via amplification of the V4 region of the 16S rRNA gene with universal primers (U515F/806R), flanked by Illumina® standard adapter sequences [29, 30]. Primers used for amplification used the TruSeq DUI adapter design, and dual-indexed forward and reverse primers were used in all reactions. PCR was performed in 50 µL reactions containing all available metagenomic DNA concentrated to a uniform volume, primers (0.2 µM each, IDT), dNTPs (200 µM each, NEB), and Phusion™ high-fidelity DNA polymerase (1U, Thermo Fisher Scientific). Amplification parameters were 98 °C(3 min) + [98 °C(15 s) + 50 °C(30 s) + 72 °C(30 s)] × 25 cycles + 72 °C(7 min). Amplicon pools (5 µL/reaction) were combined, thoroughly mixed, and then purified by addition of Axygen® Axyprep Mag™ PCR clean-up beads (Thermo Fisher Scientific) to an equal volume of 50 µL of amplicons and incubated for 15 min at room temperature (RT). Products were washed multiple times with 80% ethanol and the dried pellet resuspended in 32.5 µL EB buffer (Qiagen®), incubated for two minutes at RT, and then placed on the magnetic stand for five minutes. The final amplicon pool was evaluated using the Advanced Analytical Fragment Analyzer automated electrophoresis system, quantified using quant-iT™ HS dsDNA reagent kits (Thermo Fisher Scientific), and diluted according to Illumina’s standard protocol for sequencing on the MiSeq™ instrument using a 2 × 250 bp paired-end read flow cell.


All bioinformatics were performed at the MU Informatics Research Core Facility (Columbia, MO). Cutadapt (version 2.6; was used to remove the primer from the 5' end of the forward read. If found, the reverse complement of the primer to the reverse read was then removed from the forward read as were all bases downstream. Thus, a forward read could be trimmed at both ends if the insert was shorter than the amplicon length. The same approach was used on the reverse read, but with the primers in the opposite roles. Read pairs were rejected if one read or the other did not match a 5' primer, and an error-rate of 0.1 was allowed. Two passes were made over each read to ensure removal of the second primer. A minimal overlap of 3 bp with the 3' end of the primer sequence was required for removal.

The QIIME2 [31] DADA2 [32] plugin (version 1.10.0) was used to denoise, de-replicate, and count amplicon sequence variants (ASVs), incorporating the following parameters: 1) forward and reverse reads were truncated to 150 bases, 2) forward and reverse reads with number of expected errors higher than 2.0 were discarded, and 3) Chimeras were detected using the "consensus" method and removed. R version 3.5.1 [33] and Biom version 2.1.7 were used in QIIME2™. Taxonomies were assigned to final sequences using the Silva.v132 database, using the classify-sklearn procedure.

Total ASV count data obtained via QIIME2 were used to determine detected richness and alpha-diversity, using Past3 software [34]. Rarefaction to a uniform sequence count was not performed due to the low coverage of many samples, lack of a clear threshold in coverage, the fact that low coverage was expected a priori due to low biomass, and multiple reports that rarefaction is neither necessary or advisable in most cases [35, 36]. ASV count data were also used to identify the core taxa at different time-points, i.e., early, mid and late, based on the definitions proposed by Risely [37] and the phyloseq R package [38], designed to be a conservative measure preventing noise and spurious ASVs from being identified as being part of the core. Taxa with a minimum relative abundance of 0.1% that were prevalent in at least 50% of all the samples within each grouping at the different time-points were identified to form the core using the microbiome R package [39]. Additionally, the ggAlluvial package [40] was used to generate plots of taxa found to be overlapping between time-points and within each of the sites. The ASV count information was also used to generate co-occurrence networks within each site, across all time-points, i.e., TP1 through TP12. The log-transformed ASV abundance tables, obtained after processing the 16S rRNA sequences using the DADA2 pipeline implemented into QIIME2, were filtered to remove ASVs with relative abundance greater than or equal to 0.01% of the total number of reads. This level of filtering as based on reports indicating that rare taxa appear in only 1–5% of samples in a metabarcoding dataset [36] and the few, if any, rare taxa are over-represented in microbial networks [41]. Importantly, rare ASVs were removed from network analyses due to their potential to be sequencing artefacts [42] and to reduce false positive results [43]. However, considering a recent report suggesting that low abundance taxa may still influence compositional changes in gut microbiota [44], ASVs below a medium relative abundance of 0.01% were removed to account for any possible role of such taxa in community assembly and composition. Subsequently a weighted conditional-dependence network connecting the ASVs was built using the graphical lasso regression method employed within the Sparse InversE Covariance estimation for Ecological Association and Statistical Inference (Spiec-Easi) package [45]. The graphical lasso method was chosen over a range of independence screening with a nlambda of 30. Manual curation including the removal of self-loops, singleton nodes, and negative edges was performed, followed by generation of network plots using the igraph package [46]. Using this package, the nodes were colored according to their respective taxonomical classifications, while the size of the nodes was adjusted to indicate those with highest degree, centrality, and betweenness. The edges in the network graphs were weighted based on the relative abundance and indicate the correlation between two nodes, i.e., taxa.

Bacterial culture and identification

Replicate water samples from each site were submitted to the microbiology laboratory at IDEXX BioAnalytics (Columbia, MO) and cultured individually for bacterial growth, except for tank water samples, which were collected from three of the six tanks on the Aquaneering rack sampled for 16S rRNA sequencing and pooled as a single sample for each timepoint. Pooling was accomplished by combining 200 µL from the vortexed water samples from each site into a 1.5 mL sterile microcentrifuge tube. For each water sample, sterile PBS was used to prepare three serial dilutions (1:10, 1:100, and 1:1000). A 100 µL inoculum of each undiluted water sample and 100 µL of each dilution were plated separately onto the following bacterial culture media: BBL™ Trypticase™ Soy Agar with 5% sheep blood (TSA II™; Becton Dickinson), BBL™ CDC 5% Sheep Blood Agar with Phenylethyl Alcohol (PEA; Becton Dickinson), Difco™ Xylose Lysine Deoxycholate Agar (XLD; Becton Dickinson), and Tryptone Yeast Extract Salts (TYES) Agar, which was prepared in-house according to a published formulation. [47] Culture plates were incubated aerobically for 5 days at 22 °C. Morphologically unique colony types were selected and harvested from each plate for proteomic analysis using a direct transfer method as previously described. [48] Transferred bacteria were overlaid with 1 µL of a saturated matrix solution of α-cyano-4-hydroxycinnamic acid in 50% acetonitrile and 2.5% trifluoroacetic acid (HCCA, Bruker Daltronics, Billerica, MA) and analyzed by matrix-assisted laser desorption/ionization time-of-flight mass spectrometry (MALDI-TOF MS) using a mass spectrometer (Microflex™, Bruker Daltronics) and flexControl™ software (Bruker Daltronics). Bacterial identification was achieved using automated analysis by MALDI BioTyper® software (Bruker Daltronics) by comparing the collected spectra with integrated reference spectral databases.


Univariate outcome measures related to sample coverage, richness, and α-diversity were compared using a two-way analysis of variance (ANOVA) in a general linear model, with site and time-point as variables, and implemented in SigmaPlot® 14.0. Differences in β-diversity were tested using permutational multivariate ANOVA (PERMANOVA), with site and time-point as variables, and implemented in Past 3.26 [34]. Differences in β-diversity and dispersion were confirmed using the adonis function from the Vegan package in R, with site and time-point as sequential variables. Comparisons were made using both weighted (Bray–Curtis) and unweighted (Jaccard) similarities.


Sequencing coverage, richness, and α-diversity plateau at four weeks post-population

Samples returned anywhere from 2 to almost 2 million high-quality sequences. Only three samples were removed from further analysis due to low coverage, leaving a total of 537 samples in the study. Several samples from the initial time-point yielded low read counts, however coverage gradually increased at subsequent time-points, most notably in the Fluidized Bed Biofilter (FBB) and Post-UV disinfection samples, until reaching a peak or plateau at approximately TP5 (four weeks post-population) (Additional file 1). Two-way ANOVA revealed significant main effects of time (p < 0.001, F = 26.4) and sample site (p < 0.001, F = 288.0) with significant interactions between factors (p < 0.001, F = 23.2) reflecting that fact that, within sample site, significant time-dependent differences were detected only in FBB samples.

Similarly, there were significant effects of time (p < 0.001, F = 40.2) and sample site (p < 0.001, F = 49.7) on detected richness, which increased gradually in all sample sites between the initial sample and TP4 before reaching a peak or plateau, with the FBB harboring a significantly richer microbiota than other sites at TP4 and all subsequent time-points (p < 0.001, two-way ANOVA) (Fig. 2A). Pairwise comparisons of richness between time-points found significant differences between each of the first four time-points and later time-points, but no significant differences in richness were detected between time-points beyond T4. Likewise, comparison of Shannon diversity detected similar time-dependent (p < 0.001; F = 21.5) and site-dependent (p < 0.001, F = 27.5) differences, and demonstrated a similar pattern among most sites, with the exception of the FBB which experienced a brief decline at TP4, followed by a slow and gradual increase over the remainder of the study duration (Fig. 2B). Pairwise comparisons of Shannon diversity found significant differences between the first three time-points and later time-points, but no pairwise differences in Shannon diversity beyond T3.

Fig. 2
figure 2

Richness and diversity of the system microbiome plateaus by approximately four weeks post-population. Dot plots showing the detected richness (A) and Shannon diversity (B) at each sample site (legend at right), and at each time-point (TP) from TP1 through TP12. p and F values associated with main effects of time and sample site based on two-way analysis of variance (ANOVA)

Dominant taxa established in tanks at four to six weeks post-population

The composition of the water microbiome in the tanks showed substantial inter-tank variability during the early time-points with transient proliferation of specific taxa unique to each tank (e.g., unresolved members of the order Chitinophagales (Ch) in Tank C, Perlucidibaca (Per) in Tank D or Aeromonas (Aer) in Tank E), before ending in a uniform dominance across all tanks of two ASVs representing Staphylococcus sp. and Cutibacterium sp. by TP4 to TP6 (Fig. 3A). This was in contrast to the tank water microbiota found in the pre-existing tanks housing the fish used to populate the new system, dominated by Cetobacterium (Cet), Novosphingobium (Nov), Vibrio (Vib), and Runella (Run) spp. (Fig. 3B). While many of those taxa were detected in the new system (primarily within the FBB), they did not achieve the levels seen in the pre-existing tanks. PCoA plots revealed a gradual time-dependent change in tank water microbiome structure (p < 0.0001, F = 8.6, Fig. 3C) and composition (p < 0.0001, F = 4.6, Fig. 3D), which slowed by TP7 (eight weeks post-population) but continued to shift throughout the entire study duration.

Fig. 3
figure 3

Dominant taxa in tank water are established by two to three weeks post-population. Stacked bar charts showing the progression in microbiome composition in six tanks (A through F), sampled at 12 time-points (4 replicates/time-point) (A) and the pre-existing tanks from which fish originated; ASV = amplicon sequence variant, Aer = Aeromonas, Cet = Cetobacterium, Ch = Chitinophagales, Cut = Cutibacterium, Nov = Novosphigobium, Per = Perlucidibaca, Pl = Plesiomonas, Ru = Runella, Sap = Saprospiraceae, Sta = Staphylococcus, Vib = Vibrio (B). Principal coordinate analysis plots based on Bray–Curtis (C) and Jaccard (D) similarities. TP = time-point, color-bar legend at bottom; taxonomic abbreviations listed at beginning of manuscript. p and F values associated with main effect of time based on one-way permutation multivariate analysis of variance (PERMANOVA)

FBB and post-UV disinfection water harbor distinct microbiomes

In comparison, samples from the effluent water drained from the tanks, analyzed pre- and post-mechanical (particulate) filtration, and post-carbon filtrate returning to the UV disinfection unit revealed apparent similarities to the tank water, particularly at later time-points. While these sites mirrored the composition of the tank water with high relative abundance of Staphylococcus and Cutibacterium spp., the FBB and post-UV disinfection water samples demonstrated certain initial similarities with regard to dominant taxa followed by their own unique site-dependent temporal progressions (Fig. 4A). Post-UV disinfection water shifted dominance between various Alphaproteobacteria (e.g., Sphingopyxis, Novosphingobium, Bradyrhizobium), Betaproteobacteria (e.g., Hydrogenophaga), Gammaproteobacteria (e.g., Pseudomonas), and Actinomycetales (e.g., Cutibacterium, Nocardia), while the FBB started similarly but then showed a steady and gradual increase in evenness among 10 to 12 dominant taxa. PCoA plots reflected a similar progression as seen in tank water, although samples from the post-UV disinfection water and FBB clustered distinctly from the other sites and showed high intra-site similarity at each successive time-point. However, significant time- and site-dependent differences were detected in community structure (time: p < 0.0001, F = 5.2, site: p < 0.0001, F = 6.1, Fig. 4B) and composition (time: p < 0.0001, F = 3.0, site: p < 0.0001, F = 2.2, Fig. 4C). Testing for differences in β-diversity and dispersion according to sample site and time-point using the adonis function in the vegan R package yielded comparable results. In agreement with PERMANOVA, adonis detected significant differences in β-diversity associated with time (p < 0.001, F = 9.6) and site (p < 0.001, F = 8.1). Comparison of multivariate dispersion between groups also detected a significant difference (p < 0.001, F = 12.1). It is worth noting that all sites clustered together at the initial time-point, collected immediately prior to population of the tanks with fish, but began diverging after one week of water circulation through the system, suggesting that the observed progression in community structure over time is largely due to the introduction of zebrafish into the system, as well as the influx of nutrients from feeding them. Line graphs representing sequence numbers in the FBB of dominant taxa, and taxa previously associated with ammonia and nitrite oxidation, suggest weekly log-phase increases beginning almost immediately for Sphingopyxis, other unresolved Sphingomonadaceae, Hydrogenophaga, and Pirellulaceae, followed by similar expansions in Rhodobacteraceae, Pedosphaeraceae, and Blastocatellaceae beginning at later time-points (Additional file 2). Notably, Nitrospira spp. were less abundant than the aforementioned taxa in the FBB by orders of magnitude, and members of the Nitrosomonadaceae family were rare to undetected. Similarly, several unresolved members of the Nitrososphaeraceae (likely ammonia-oxidizing archaea) were detected but at extremely low prevalence and read counts.

Fig. 4
figure 4

Post-UV disinfection water and FBB harbor dynamic communities, distinct from other sites. Stacked bar charts showing the progression in microbiome composition in the post-UV disinfection water entering the tanks (Post-UV disinfection), the pre- and post-particulate filter water (Pre-P and Post-P, respectively), the fluidized bed biofilter (FBB), and post-carbon filter water (Post-C), sampled at 12 time-points (4 replicates/time-point); ASV = amplicon sequence variant, Aer = Aeromonas, Al = Alloiococcus, Aq = Aquabacterium, Bl = Blastocatellaceae, Br = Bradyrhizobium, Bu = Burkholderiaceae, Cet = Cetobacterium, Ch = Chitinophagales, Cut = Cutibacterium, FBB = fluidized bed biofilter, Fer = Ferruginibacter, Hy = Hydrogenophaga, Lim = Limnobacter, Noc = Nocardia, Nov = Novosphingobium, Ped = Pedosphaeraceae, Per = Perlucidibaca, Pl = Plesiomonas, Ps = Pseudomonas, Rh = Rhodobacteraceae, Sap = Saprospiraceae, Spp = Sphingopyxis, Spm = Sphingomonadaceae, Sta = Staphylococcus, Vib = Vibrio (A) and principal coordinate analysis plots of those samples and tank water, based on Bray–Curtis (B) and Jaccard (C) similarities. ASV = amplicon sequence variant, TP = time-point, color-bar and symbol legends at bottom; taxonomic abbreviations listed at beginning of manuscript. p and F values associated with main effects of time and sample site based on two-way permutation multivariate analysis of variance (PERMANOVA)

Core taxa comprise Actinobacteria, Firmicutes, and variable Proteobacteria

To define core taxa during establishment of the system microbiome, amplicon sequence variants (ASVs) were stratified by prevalence and relative abundance. To simplify the analysis and interpretation, time-points were divided categorically into early (TP1 to TP4, at 0-3w post-population), mid (TP5 to TP8 at 4-10w post-population), and late (TP9 to TP12 at 12-18w post-population) time-points (Fig. 5A). Within all three periods, ASV23 and ASV28, representing Cutibacterium and Staphylococcus respectively, were the two dominant taxa system-wide. Other core taxa across all three periods include other members of the Actinobacteria (Lawsonella, Micrococcus, Propionibacterium), Streptococcus, Cetobacterium, and Aeromonas (Additional file 3). While Proteobacteria were commonly identified as core taxa, few were consistently identified as such across all three periods of the study. The distributions of those core taxa within and between each site are shown as alluvial plots (Fig. 5B).

Fig. 5
figure 5

Limited number of taxa comprise core community throughout equilibration. Heatmaps showing prevalence (legend at right) of core ASVs at increasing thresholds of relative abundance at Early (TP1 to TP4), Mid (TP5 to TP8), and Late (TP9 to TP12) time-points (A); alluvial plots showing distribution of core taxa among sample sites within each period of time (legend at right) (B). See Additional file 3 for taxonomic identity of core ASVs

Microbial network analysis indicates central role for Nocardiaceae in FBB

To identify microbial interaction networks associated with the equilibration and stabilization of the tanks, the FBB, and the system as a whole, a graphical lasso regression method was used to infer ecological associations from a sparse matrix such as ASV counts across all time-points. Depending on the site from which the samples were collected, the density of the networks varied. The post-carbon filtration samples had the least density, i.e., connectivity, and were dominated by one or two taxa as indicated by the size of the node. Within each site, smaller clusters were found outside of the largest cluster of taxa indicating the overall niche preferences of the respective taxa and their contributions to the community stability within the system.

The predicted interaction networks in tank water, pre- and post-particulate filter water, and post-carbon filter water over time were relatively sparse in contrast to those detected in the post-UV disinfection water and FBB (Fig. 6A through F). The network in FBB identified Nocardiaceae as a hub taxon between two sub-networks of bacteria, owing to its high level of degree, betweenness, and centrality. A combined network analysis of the entire system again placed Nocardiaceae as pivotal community members alongside Pseudomonadaceae, Nitrosomonadaceae, Rhizobiales incertae sedis, and Rhizobiaceae (Fig. 6G).

Fig. 6
figure 6

Enriched networks present in FBB and post-UV disinfection water dominated by several Proteobacteria. Community network analysis plots based on temporal data for post-UV disinfection unit (A), tanks (B), pre-particulate filter (C), post-particulate filter (D), FBB (E), post-carbon filter (F) samples, and the system as a whole (G). Nodes indicate amplicon sequence variants (ASVs) while the edges in the network graphs are weighted based on the relative abundance and indicate the correlation between two nodes. The nodes are colored based on their respective taxonomical classifications and the size of the nodes indicate the highest degree, betweenness, and centrality within the overall network

Only a fraction of detected ASVs were represented in the culturable fraction

Lastly, to characterize the culturable portion of the system microbiota and aid in identification or resolution of detected taxa, replicates of all samples were serially diluted and plated on four different media selected to grow a broad range of environmental and aquatic bacterial taxa. All cultured isolates were then analyzed via matrix-assisted laser desorption/ionization-time of flight (MALDI-TOF) mass spectrometry and annotated against a protein spectrum database. In total, 269 isolates were recovered, resulting in 33 distinct identifiable isolates (Additional file 4) and three unidentified isolates. The most common isolates were identified as Pseudomonas alcaliphila, Limnobacter thiooxidans, Acidovorax facilis, and Aeromonas veronii.

As 16S rRNA amplicon sequencing is unable to resolve certain species or genera (due to genetic homogeneity within a given family at the 16S region in question or lack of relevant taxonomies in the database used for annotation), efforts were made to identify culture isolates within the sequencing data via post hoc comparison of poorly resolved ASV sequences against annotations made using the MALDI-TOF mass spectrometry, with the ultimate goal of visualizing the relative abundance of cultured isolates throughout the entire time-course dataset. While not specific to the cultured isolates, candidate ASVs with a 100% nucleotide identity to a culture isolate at the species level were identified for 17 of the 33 culture isolates. Of the remaining 16 isolates, curation of all poorly resolved ASVs in the next higher taxonomic division returned a 99.6% nucleotide identity in four isolates, or was only annotated to the level of genus via both MALDI-TOF and 16S sequencing in 11 isolates. Only one isolate, Tsukamurella sp., could not be matched to a candidate ASV at any taxonomic level. The mean relative abundance in all samples across time of the 32 putative ASVs matching culture isolates is shown in Fig. 7.

Fig. 7
figure 7

Cultivable taxa expand early during equilibration. Heatmap showing cube root-transformed mean relative abundance (legend upper left) at each time point (TP), of amplicon sequence variants (ASVs) matching the taxonomies assigned to culture isolates via MALDI-TOF


Previous studies have suggested that the gut microbiome assemblage of fishes is determined by a combination of extrinsic (e.g., salinity, trophic level) and intrinsic (e.g., host taxonomy) factors [49]. Research focused on the gut microbiome of zebrafish per se suggests that colonization of the zebrafish gut is not a stochastic, or neutral, process, but rather is influenced by active processes including microbe-microbe interactions or host selection [50]. Indeed, Roeselers et al. found evidence of the historical connections between research facilities in the gut microbiome of the fish at those institutions, indicating that zebrafish acquire facility-dependent features within their microbiome over time [51]. While several studies have been performed applying molecular approaches to characterize the microbial community of commercial recirculating aquaculture systems (RAS), we were unable to identify any meaningful surveys of commercial zebrafish housing systems using culture-independent methods. As we were preparing to install and implement a commercial zebrafish housing system for research purposes, we wanted to determine the time point at which the environmental microbiome of a newly installed system had stabilized and would thus be suitable for use in experimental procedures. Based on multiple metrics of community composition, it appears that the tank water and most of the filter water samples have peaked, plateaued, or reached some sort of stable equilibrium by TP4 (3w post-population) to TP6 (6w post-population). This was true of sample coverage across all sites (as a very rough gauge of relative biomass), community richness and α-diversity, and the dominance of the two core taxa, Staphylococcus and Cutibacterium. That is not to say that no change occurs after TP6 as the PCoA plots demonstrated continued drift in β-diversity throughout the entire study duration. Moreover, the FBB and post-UV disinfection samples showed very different patterns throughout the study period, including different community members and kinetics. This is not surprising with the FBB as this is a fluidized fine sand media designed to facilitate bacterial colonization via high surface area, and inoculated with a proprietary mixture of bacteria at the commercial production facility. However, the presence of a distinct microbiome within the effluent from the UV disinfection unit, completely different from the post-carbon filter water entering the UV disinfection unit, was unexpected. For those samples, the small piece of tubing carrying water from the UV disinfection unit to the tanks was temporarily disconnected from the tank and allowed to drain 1 mL of water directly into a sterile container for DNA extraction. The only difference between this and the post-carbon filter water samples is passage through the UV disinfection unit, yet the communities are strikingly different. The ultraviolet exposure achieved in RAS is inadequate to sterilize water, and bacterial taxa vary widely in their susceptibility to ultraviolet irradiation [52, 53]. The disproportionately high coverage of post-UV disinfection water samples at almost every time point (and clear community progression over time) suggests the presence of a stable biofilm upstream of our sampling site. Ultraviolet disinfection is more effective at reducing individual planktonic bacteria suspended in the water column than for bacterial species that aggregate in water [54], or bacteria incorporated into biofilms that are thus poorly penetrated by ultraviolet radiation. Many members of Actinomycetales are hydrophobic, including Nocardia spp. [55], an abundant taxon in the post-UV disinfection water samples, likely reflecting biofilms that became established as the system equilibrated. Alternatively, the bacterial DNA collected from the UV effluent may simply represent a UV-resistant fraction of the community, with reduced fragmentation of the 16S rRNA gene, and subsequent over-representation in the pool of amplicon libraries. These possible explanations are not mutually exclusive.

The high relative abundance of Staphylococcus spp. and Cutibacterium spp. in tank water and other sites is noteworthy due to the role of each genus as dominant members of the human skin microbiome [56]. Moreover, multiple members of each genus are capable of biofilm formation [5761], suggesting they might be particularly well-suited for colonization of RAS following dissemination from the skin of individuals maintaining the system.

That the tank water in the new system never became similar in composition to the tank water from the existing housing system was not entirely unexpected. The older system was not commercially purchased, but rather, was constructed from the necessary material according to a published design [27]. As such, tank size, flow-through rate, and filtration systems all varied between system. Additionally, the existing system had been in use for several years, much longer than the study period for the new system.

Regarding the community within the FBB ostensibly responsible for nitrification of the system via ammonia-oxidizing bacteria (AOB), Nitrosomonas sp. and related taxa were surprisingly rare, being detected in less than 10% of samples and at extremely low relative abundance, and Nitrobacter spp. were not detected at all. In contrast, Nitrospira sp. were present at modest relative abundance, primarily in the FBB. Nitrospira spp. are known to oxidize nitrites in freshwater aquaria [20], and some species are capable of complete nitrification from ammonia to nitrate [62]. Multiple species of Spingomonas and Sphingopyxis also express nitrate reduction machinery [6366], and thus may also play a role in denitrification in the FBB system. In contrast, we note that other dominant taxa identified in the FBB including Hydrogenophaga, Sphingopyxis, Pirellula, and unresolved members of the families Rhodobacteraceae and Sphingomonadaceae are capable of aerobic denitrification. While we were unable to identify any substantive published reports of the microbiome present in RAS used for research purposes, several recent surveys of the water or bioreactors in commercial RAS used for production of fish or shrimp also identified Rhodobacteraceae and Planctomycetes as dominant taxa based on relative abundance of 16S rRNA sequences [6770]. Similarly, multiple studies have implicated Hydrogenophaga as a participant in aerobic denitrification in various closed bioreactor systems [7173]. One limitation of the current study is that nitrogenous compounds were at or below the lower limit of detection, but the lack of any fish mortality during (or after) the study period suggests that nitrification was occurring. The bacteria traditionally associated with nitrification (e.g., Nitrosomonas, Nitrobacter) are extremely slow-growing bacteria with a doubling time of 18–70 hours [7476], with evidence that AOB have a faster generation time than nitrite-oxidizing bacteria (NOB) [75]. Thus, if only these species contributed to denitrification, bacterial reproductive capacity could be a rate-limiting factor in equilibration regardless of how quickly new fish are added. However, in part because the number of prokaryotic taxa that contribute to nitrification is likely larger, it is unclear whether the rate of compositional change would have been different if significantly more fish had been introduced during the initial population of tanks.

Culture-based screening of the system complemented the molecular analysis by demonstrating viability of several core members of the system microbiome, and improving the taxonomic resolution of several of these members. MALDI-TOF mass spectrometry annotations improved the taxonomic resolution relative to putative matches in the 16S rRNA sequencing data in 19 of 33 isolates. The greatest overlap between the core taxa and culturable fraction of the system microbiome was Proteobacteria. Of the 36 taxa defined as core taxa by our criteria (in at least one period of time), only 12 were ostensibly cultured. Of those 12, 10 were Proteobacteria, the only non-Proteobacteria isolates being Micrococcus luteus and Staphylococcus warneri. This is likely, at least partially, due to the aerobic culture approach. However, lack of other requirements in the culture media for growth of certain core taxa is also likely.

Four culture media, including two permissive (non-selective) media and two selective media were used to isolate a broad array of bacteria. In order to document changing bacterial communities at timepoints occurring before and after the eutrophication associated with the introduction of live zebrafish and zebrafish feed into the system, non-selective media were selected that would facilitate isolation of bacteria that thrive under both oligotrophic and eutrophic nutrient conditions. The BBL™ Trypticase™ Soy Agar with 5% sheep blood (blood agar) was chosen as an enriched differential and non-selective media. Tryptone Yeast Extract Salts (TYES) Agar is a non-selective minimal medium that facilitates the isolation of oligotrophic and fastidious aquatic microorganisms that either will not grow on enriched media or grow very slowly on enriched media, and are thus easily overgrown by rapidly growing bacterial species. Xylose Lysine Deoxycholate (XLD) Agar is a culture medium that selects for a subset of Gram-negative bacteria and is highly differential, aiding in distinguishing between Gram-negative bacteria with similar colony morphologies. Sheep Blood Agar with Phenylethyl Alcohol (PEA) is a selective medium that inhibits Gram-negative bacteria and thus facilitates isolation of Gram-positive bacterial species.

The temporal network data for individual sites revealed that the post-UV disinfection water and FBB substrate showed the highest co-occurrence and dense patterns among all the samples. The other sites, including the tank water and particulate filtration samples demonstrated fewer interactions. It is plausible that these organisms form a multi-cellular biofilm that may be recalcitrant to extreme environments [77]. Interestingly, the overall system network identified Nocardiaceae as a key player in these co-occurring taxa communities in conjunction with Pseudomonadaceae. More work including the use of metagenomics and stable-isotope probing experiments will be needed in the future to delineate the exact mechanisms by which two organisms interact. Additionally, it has been shown that in a nitrifying medium or a nitrate-rich environment, Nitrosomonadaceae and Pseudomonadaceae are thought to be prevalent and mutualistic, especially in biofilms [78, 79]. This may be a plausible reason for the co-occurrence of these taxa, and also for their placement as dominant members of the overall network community. Moreover, Keshvardoust et al. demonstrated that enriching the glucose content of media leads to a shift in dominance from Nitrosomonadaceae to Pseudomonadaceae, indicating a potentially competitive relationship between these taxa within the network. In either case, supplementation, and carefully constructed synthetic community models starting with combinations of specific taxa will be required to validate these findings in the future.

In summary, the data reported here provide a detailed and comprehensive characterization of the prokaryotic communities present at different sites of a research zebrafish RAS during establishment. Collectively, these data suggest that a peak population density occurs at roughly 3 to 4 weeks post-population, although the FBB continued to undergo subtle changes in evenness throughout the 18-week study duration. Moreover, our data strongly suggest the presence of bacterial biofilm communities associated with the UV disinfection unit, representing an unappreciated nidus of bacteria within RAS. Lastly, these data demonstrate the complementary abilities of molecular approaches and traditional culture coupled to MALDI-TOF, to characterize complex microbial communities.

Availability of data and materials

The dataset supporting the conclusions of this article are available in the NCBI Sequence Read Archive (SRA) under BioProject ID: PRJNA674483 (Submission ID: SUB8465482 and SUB8465542).









Amplicon sequence variant














Fluidized bed biofilter




Gut microbiota


















Recirculating aquaculture system


















  1. Ericsson AC, Crim MJ, Franklin CL. A brief history of animal modeling. Mo Med. 2013;110:201–5.

    PubMed  PubMed Central  Google Scholar 

  2. Vliegenthart AD, Tucker CS, Del Pozo J, et al. Zebrafish as model organisms for studying drug-induced liver injury. Br J Clin Pharmacol. 2014;78:1217–27.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  3. Veldman MB, Lin S. Zebrafish as a developmental model organism for pediatric research. Pediatr Res. 2008;64:470–6.

    Article  PubMed  Google Scholar 

  4. Lin CY, Chiang CY, Tsai HJ. Zebrafish and Medaka: new model organisms for modern biomedical research. J Biomed Sci. 2016;23:19.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Cassar S, Adatto I, Freeman JL, et al. Use of Zebrafish in Drug Discovery Toxicology. Chemical research in toxicology 2020; 33: 95–118. 2019/10/19. DOI:

  6. Wiley DS, Redfield SE, Zon LI. Chemical screening in zebrafish for novel biological and therapeutic discovery. Methods Cell Biol. 2017;138:651–79.

    Article  CAS  PubMed  Google Scholar 

  7. Kalueff AV, Gebhardt M, Stewart AM, et al. Towards a comprehensive catalog of zebrafish behavior 1.0 and beyond. Zebrafish. 2013;10:70–86.

    Article  PubMed  PubMed Central  Google Scholar 

  8. Demin KA, Lakstygal AM, Volgin AD, et al. Cross-species analyses of intra-species behavioral differences in mammals and fish. Neuroscience. 2020;429:33–45.

    Article  CAS  PubMed  Google Scholar 

  9. Demin KA, Lakstygal AM, Chernysh MV, et al. The zebrafish tail immobilization (ZTI) test as a new tool to assess stress-related behavior and a potential screen for drugs affecting despair-like states. J Neurosci Methods. 2020;337:108637.

    Article  CAS  PubMed  Google Scholar 

  10. Golla A, Ostby H, Kermen F. Chronic unpredictable stress induces anxiety-like behaviors in young zebrafish. Scientific reports. 2020;10:10339.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  11. Kermen F, Darnet L, Wiest C, et al. Stimulus-specific behavioral responses of zebrafish to a large range of odors exhibit individual variability. BMC Biol. 2020;18:66.

    Article  PubMed  PubMed Central  Google Scholar 

  12. Spinello C, Yang Y, Macri S, et al. Zebrafish adjust their behavior in response to an interactive robotic predator. Front Robot AI. 2019.

    Article  PubMed  PubMed Central  Google Scholar 

  13. Geng Y, Peterson RT. The zebrafish subcortical social brain as a model for studying social behavior disorders. Disease models & mechanisms. 2019.

    Article  Google Scholar 

  14. Kalueff AV, Stewart AM, Gerlai R. Zebrafish as an emerging model for studying complex brain disorders. Trends Pharmacol Sci. 2014;35:63–75.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Basnet RM, Zizioli D, Taweedet S, et al. Zebrafish Larvae as a Behavioral Model in Neuropharmacology. Biomedicines. 2019.

    Article  PubMed  PubMed Central  Google Scholar 

  16. Franklin CL, Ericsson AC. Microbiota and reproducibility of rodent models. Lab Anim (NY). 2017;46:114–22.

    Article  Google Scholar 

  17. Okazaki F, Zang L, Nakayama H, et al. Microbiome Alteration in Type 2 Diabetes Mellitus Model of Zebrafish. Sci Rep. 2019;9:867.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Catron TR, Swank A, Wehmas LC, et al. Microbiota alter metabolism and mediate neurodevelopmental toxicity of 17beta-estradiol. Sci Rep. 2019;9:7064.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Burrell PC, Phalen CM, Hovanec TA. Identification of bacteria responsible for ammonia oxidation in freshwater aquaria. Appl Environ Microbiol. 2001;67:5791–800.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Hovanec TA, Taylor LT, Blakis A, et al. Nitrospira-like bacteria associated with nitrite oxidation in freshwater aquaria. Appl Environ Microbiol. 1998;64:258–64.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  21. Hammer HS. Water Quality For Zebrafish Culture. In: Cartner SL, Eisen JS, Farmer S, et al. (eds) The Zebrafish in Biomedical Research: Biology, Husbandry, Diseases, and Research Applications. Elsevier, 2020, pp.321–335.

  22. Breen P, Winters AD, Nag D, et al. Internal Versus External Pressures: Effect of Housing Systems on the Zebrafish Microbiome. Zebrafish. 2019;16:388–400.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Crim MJ and Lawrence C. A fish is not a mouse: understanding differences in background genetics is critical for reproducibility. Lab Anim. 2021;50:19–25.

  24. Bauer BA, Besch-Williford C, Livingston RS, et al. Influence of rack design and disease prevalence on detection of rodent pathogens in exhaust debris samples from individually ventilated caging systems. J Am Assoc Lab Anim Sci. 2016;55:782–8.

    PubMed  PubMed Central  Google Scholar 

  25. Trinh P, Zaneveld JR, Safranek S, et al. One Health Relationships Between Human, Animal, and Environmental Microbiomes: A Mini-Review. Front Public Health. 2018;6:235.

    Article  PubMed  PubMed Central  Google Scholar 

  26. Lieggi C, Kalueff AV, Lawrence C, et al. The influence of behavioral, social, and environmental factors on reproducibility and replicability in aquatic animal models. ILAR J. 2020.

    Article  PubMed  PubMed Central  Google Scholar 

  27. Kim S, Carlson R, Zafreen L, et al. Modular, easy-to-assemble, low-cost zebrafish facility. Zebrafish. 2009;6:269–74.

    Article  PubMed  PubMed Central  Google Scholar 

  28. Council NR. Guide for the care and use of laboratory animals. National Academies Press, 2010.

  29. Walters WA, Caporaso JG, Lauber CL, et al. PrimerProspector: de novo design and taxonomic analysis of barcoded polymerase chain reaction primers. Bioinformatics. 2011;27:1159–61.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  30. Caporaso JG, Lauber CL, Walters WA, et al. Global patterns of 16S rRNA diversity at a depth of millions of sequences per sample. Proc Natl Acad Sci USA. 2011;108(Suppl 1):4516–22.

    Article  PubMed  Google Scholar 

  31. Bolyen E, Rideout JR, Dillon MR, et al. Reproducible, interactive, scalable and extensible microbiome data science using QIIME 2. Nat Biotechnol. 2019;37:852–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  32. Callahan BJ, McMurdie PJ, Rosen MJ, et al. DADA2: High-resolution sample inference from Illumina amplicon data. Nat Methods. 2016;13:581–3.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  33. Foster ZS, Sharpton TJ and Grunwald NJ. Metacoder: An R package for visualization and manipulation of community taxonomic diversity data. PLoS computational biology 2017; 13: e1005404. 2017/02/22. DOI:

  34. Hammer O, Harper DAT. PAST: Paleontological statistics software package for education and data analysis. Palaeontol Electron. 2011;4:1–9.

    Google Scholar 

  35. McMurdie PJ, Holmes S. Waste not, want not: why rarefying microbiome data is inadmissible. PLoS Comput Biol. 2014;10:e1003531.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Cao Q, Sun X, Rajesh K, et al. Effects of Rare Microbiome Taxa Filtering on Statistical Analysis. Front Microbiol. 2020;11:607325.

    Article  PubMed  Google Scholar 

  37. Risely A. Applying the core microbiome to understand host-microbe systems. J Anim Ecol. 2020;89:1549–58.

    Article  PubMed  Google Scholar 

  38. McMurdie PJ, Holmes S. phyloseq: an R package for reproducible interactive analysis and graphics of microbiome census data. PLoS ONE. 2013;8: e61217.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  39. Lahti L and Shetty S. Tools for microbiome analysis in R. Version 2.1.26. 2017.

  40. Brunson JC. ggalluvial: Layered grammar for alluvial plots. Journal of Open Source Software 2017; 5. DOI:

  41. Ma B, Wang Y, Ye S, et al. Earth microbial co-occurrence network reveals interconnection pattern across microbiomes. Microbiome 2020; 8: 82. 2020/06/06. DOI:

  42. Martinson VG, Douglas AE, Jaenike J. Community structure of the gut microbiota in sympatric species of wild Drosophila. Ecol Lett. 2017;20:629–39.

    Article  PubMed  Google Scholar 

  43. Weiss S, Van Treuren W, Lozupone C, et al. Correlation detection strategies in microbial data sets vary widely in sensitivity and precision. ISME J. 2016;10:1669–81.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  44. Benjamino J, Lincoln S, Srivastava R, et al. Low-abundant bacteria drive compositional changes in the gut microbiota after dietary alteration. Microbiome. 2018;6:86.

    Article  PubMed  PubMed Central  Google Scholar 

  45. Kurtz ZD, Muller CL, Miraldi ER, et al. Sparse and compositionally robust inference of microbial ecological networks. PLoS Comput Biol. 2015;11:e1004226.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  46. Csardi G and Nepusz T. The igraph software package for complex network research. InterJournal 2006; Complex Systems.

  47. Heil N. National wild fish health survey—laboratory procedures manual. US Fish and Wildlife Service Warm Springs, GA 2009.

  48. Philips BH, Crim MJ, Hankenson FC, et al. Evaluation of presurgical skin preparation agents in African clawed frogs (Xenopus laevis). J Am Assoc Lab Anim Sci. 2015;54:788–98.

    PubMed  PubMed Central  Google Scholar 

  49. Wong S, Rawls JF. Intestinal microbiota composition in fishes is influenced by host ecology and environment. Mol Ecol. 2012;21:3100–2.

    Article  PubMed  PubMed Central  Google Scholar 

  50. Burns AR, Stephens WZ, Stagaman K, et al. Contribution of neutral processes to the assembly of gut microbial communities in the zebrafish over host development. ISME J. 2016;10:655–64.

    Article  CAS  PubMed  Google Scholar 

  51. Roeselers G, Mittge EK, Stephens WZ, et al. Evidence for a core gut microbiota in the zebrafish. ISME J. 2011;5:1595–608.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Lee ES, Yoon TH, Lee MY, et al. Inactivation of environmental mycobacteria by free chlorine and UV. Water Res. 2010;44:1329–34.

    Article  CAS  PubMed  Google Scholar 

  53. Hammer HS. Recirculating aquaculture systems (RAS) for zebrafish culture. In: Cartner SL, Eisen JS, Farmer S, et al. (eds) The Zebrafish in Biomedical Research: Biology, Husbandry, Diseases, and Research Applications. Elsevier, 2020, pp.337–356.

  54. Bohrerova Z, Linden KG. Ultraviolet and chlorine disinfection of mycobacterium in wastewater: effect of aggregation. Water Environ Res. 2006;78:565–71.

    Article  CAS  PubMed  Google Scholar 

  55. Luo Q, Hiessl S, Steinbuchel A. Functional diversity of Nocardia in metabolism. Environ Microbiol. 2014;16:29–48.

    Article  CAS  PubMed  Google Scholar 

  56. Skowron K, Bauza-Kaszewska J, Kraszewska Z, et al. Human skin microbiome: impact of intrinsic and extrinsic factors on skin microbiota. Microorganisms. 2021.

    Article  PubMed  PubMed Central  Google Scholar 

  57. Cucarella C, Solano C, Valle J, et al. Bap, a Staphylococcus aureus surface protein involved in biofilm formation. J Bacteriol. 2001;183:2888–96.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Fredheim EG, Klingenberg C, Rohde H, et al. Biofilm formation by Staphylococcus haemolyticus. J Clin Microbiol. 2009;47:1172–80.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  59. Iwase T, Uehara Y, Shinji H, et al. Staphylococcus epidermidis Esp inhibits Staphylococcus aureus biofilm formation and nasal colonization. Nature. 2010;465:346–9.

    Article  CAS  PubMed  Google Scholar 

  60. Kuehnast T, Cakar F, Weinhaupl T, et al. Comparative analyses of biofilm formation among different Cutibacterium acnes isolates. Int J Med Microbiol. 2018;308:1027–35.

    Article  CAS  PubMed  Google Scholar 

  61. Nakamura K, O’Neill AM, Williams MR, et al. Short chain fatty acids produced by Cutibacterium acnes inhibit biofilm formation by Staphylococcus epidermidis. Sci Rep. 2020;10:21237.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  62. Daims H, Lebedeva EV, Pjevac P, et al. Complete nitrification by Nitrospira bacteria. Nature. 2015;528:504–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  63. Wei S, Wang T, Liu H, et al. Sphingomonas hengshuiensis sp. nov., isolated from lake wetland. Int J Syst Evol Microbiol. 2015;65:4644–9.

    Article  CAS  PubMed  Google Scholar 

  64. Ravintheran SK, Sivaprakasam S, Loke S, et al. Complete genome sequence of Sphingomonas paucimobilis AIMST S2, a xenobiotic-degrading bacterium. Sci Data. 2019;6:280.

    Article  PubMed  PubMed Central  Google Scholar 

  65. Cua LS, Stein LY. Characterization of denitrifying activity by the alphaproteobacterium, Sphingomonas wittichii RW1. Front Microbiol. 2014;5:404.

    Article  PubMed  PubMed Central  Google Scholar 

  66. Garcia-Romero I, Perez-Pulido AJ, Gonzalez-Flores YE, et al. Genomic analysis of the nitrate-respiring Sphingopyxis granuli (formerly Sphingomonas macrogoltabida) strain TFA. BMC Genom. 2016;17:93.

    Article  CAS  Google Scholar 

  67. D'Silva A and Kyndt JA. Metagenomics Analysis - bacterial diversity greatly affects ammonia and overall nitrogen levels in aquabioponics bioflocs systems, based on 16S rRNA gene amplicon metagenomics. Applied Microbiology: Open Access 2020; 6. DOI:

  68. Brailo M, Schreier HJ, McDonald R, et al. Bacterial community analysis of marine recirculating aquaculture system bioreactors for complete nitrogen removal established from a commercial inoculum. Aquaculture. 2019;503:198–206.

    Article  CAS  PubMed  Google Scholar 

  69. Chen Z, Chang Z, Zhang L, et al. Effects of water recirculation rate on the microbial community and water quality in relation to the growth and survival of white shrimp (Litopenaeus vannamei). BMC Microbiol. 2019;19:192.

    Article  PubMed  PubMed Central  Google Scholar 

  70. Huang Z, Wan R, Song X, et al. Metagenomic analysis shows diverse, distinct bacterial communities in biofilters among different marine recirculating aquaculture systems. Aquacult Int. 2016.

    Article  Google Scholar 

  71. Feng C, Huang L, Yu H, et al. Simultaneous phenol removal, nitrification and denitrification using microbial fuel cell technology. Water Res. 2015;76:160–70.

    Article  CAS  PubMed  Google Scholar 

  72. Zhou S-L, Sun Y, Zhang Y-R, et al. Variations in microbial community during nitrogen removal by in situ oxygen-enhanced indigenous nitrogen-removal bacteria. Water Sci Eng. 2018;11:276–87.

    Article  Google Scholar 

  73. Yan Q, Bi Y, Deng Y, et al. Impacts of the three gorges dam on microbial structure and potential function. Sci Rep. 2015;5:8605.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Chen S, Ling J, Blancheton J-P. Nitrification kinetics of biofilm as affected by water quality factors. Aquacult Eng. 2006;34:179–97.

    Article  Google Scholar 

  75. Hagopian DS, Riley JG. A closer look at the bacteriology of nitrification. Aquacult Eng. 1998;18:223–44.

    Article  Google Scholar 

  76. Rurangwa E, Verdegem MC. Microorganisms in recirculating aquaculture systems and their management. Rev Aquac. 2015;7:117–30.

    Article  Google Scholar 

  77. Yin W, Wang Y, Liu L, et al. Biofilms: The Microbial “Protective Clothing” in Extreme Environments. Int J Mol Sci. 2019.

    Article  PubMed  PubMed Central  Google Scholar 

  78. Keshvardoust P, Huron VAA, Clemson M, et al. Biofilm formation inhibition and dispersal of multi-species communities containing ammonia-oxidising bacteria. NPJ Biofilms Microbiomes. 2019;5:22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  79. Petrovich M, Wu CY, Rosenthal A, et al. Nitrosomonas europaea biofilm formation is enhanced by Pseudomonas aeruginosa. FEMS Microbiol Ecol. 2017.

    Article  PubMed  Google Scholar 

Download references


We would like to acknowledge the outstanding support of Nathan Bivens and staff of the MU DNA Core, and Shahida Mohammed Ashraf and Christopher Bottoms at the MU Informatics Research Core.


No funding was received for this study.

Author information

Authors and Affiliations



ACE, DJD, MJC, and ECB conceived and designed the experiments; HN, RAD, GT, and PSO collected, processed, and organized samples, and performed DNA extraction; MC and DE performed bacterial culture and identification via MALDI-TOF mass spectrometry; ACE and SBB analyzed sequencing data; all authors contributed to the writing and review of the manuscript. All authors read and approved the final manuscript.

Corresponding authors

Correspondence to Aaron C. Ericsson or Elizabeth C. Bryda.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

MJC and DCE are employees of IDEXX BioAnalytics, a division of IDEXX Laboratories, Inc., a company that provides veterinary diagnostics, including microbiology services. All other authors declare that they have no other competing interests.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary Information

Additional file 1.

Dot plot showing the number of 16S rRNA amplicon sequences recovered from each sample site (legend at right) at each time-point. Pre-P and Post-P = pre- and post-particulate filter water, FBB = Fluidized bed biofilter substrate, Post-C = post-carbon filter water, TP = time-point. p and F values associated with main effects of time and sample site based on two-way analysis of variance (ANOVA).

Additional file 2.

Line chart showing mean (± SD) sequence number across time in the fluidized bed biofilter (FBB), of dominant taxa and taxa recognized to participate in the oxidation of ammonia and nitrites, or reduction of nitrates and nitrites, on a Log-scale. Sphingomonadaceae NR includes all sequences matched to that family but not resolved to the level of genus, TP = time-point.

Additional file 3.

Core taxa found at a minimum of 0.01% mean relative abundance in 50% of all samples during the Early (TP1 to TP4), Mid (TP5 to TP8), or Late (TP9 to TP12) time-points.

Additional file 4.

Taxonomic identity of 33 different culture isolates listed by time-point at which isolates were cultured (columns TP1-TP3, TP5-TP12), and the sample site (numbers: 1 = post-UV, 2 = tanks, 3 = pre-particulate, 4 = post-particulate, 5 = FBB, 6 = post-carbon). Column labeled ASV_ID and Taxonomic annotation indicate putative matches in the sequencing data, following post hoc BLAST analysis.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Ericsson, A.C., Busi, S.B., Davis, D.J. et al. Molecular and culture-based assessment of the microbiome in a zebrafish (Danio rerio) housing system during set-up and equilibration. anim microbiome 3, 55 (2021).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: