Skip to main content
  • Research Article
  • Open access
  • Published:

Monoassociation with bacterial isolates reveals the role of colonization, community complexity and abundance on locomotor behavior in larval zebrafish



Across taxa, animals with depleted intestinal microbiomes show disrupted behavioral phenotypes. Axenic (i.e., microbe-free) mice, zebrafish, and fruit flies exhibit increased locomotor behavior, or hyperactivity. The mechanism through which bacteria interact with host cells to trigger normal neurobehavioral development in larval zebrafish is not well understood. Here, we monoassociated zebrafish with either one of six different zebrafish-associated bacteria, mixtures of these host-associates, or with an environmental bacterial isolate.


As predicted, the axenic cohort was hyperactive. Monoassociation with three different host-associated bacterial species, as well as with the mixtures, resulted in control-like locomotor behavior. Monoassociation with one host-associate and the environmental isolate resulted in the hyperactive phenotype characteristic of axenic larvae, while monoassociation with two other host-associated bacteria partially blocked this phenotype. Furthermore, we found an inverse relationship between the total concentration of bacteria per larvae and locomotor behavior. Lastly, in the axenic and associated cohorts, but not in the larvae with complex communities, we detected unexpected bacteria, some of which may be present as facultative predators.


These data support a growing body of evidence that individual species of bacteria can have different effects on host behavior, potentially related to their success at intestinal colonization. Specific to the zebrafish model, our results suggest that differences in the composition of microbes in fish facilities could affect the results of behavioral assays within pharmacological and toxicological studies.


Microbes have co-evolved with their hosts and are required for host health and development, including neurodevelopment and behavior [1]. Across taxa, animals raised under microbe-free conditions or depleted of their intestinal microbes consistently show alterations in behavior and physiology [1, 2]. For example, both mice and zebrafish lacking microbes show increased locomotor activity, with a critical window in development during which exposure to microbes can result in control-like behavior [3,4,5]. There are several potential mechanisms through which microbes may influence behavioral development. The microbiota-gut-brain axis allows bidirectional communication via activation of the vagus nerve. Another route of influence is through the direct or indirect production of metabolites that can interact with the immune system or penetrate the blood-brain-barrier and exert direct effects on the developing nervous system [2].

While much progress has been made, the specific mechanisms through which microbes influence brain development and behavior remain unknown. In particular, the role of microbial taxonomic and functional diversity is unclear. For example, monoassociation of axenic mice with Bifidobacterium infantis was shown to reduce restraint stress and restore normal behavior [6]. In contrast, only a complex microbiota can restore microglia insufficiency in mice [7], suggesting that effects on microglia are not regulated by bacterial load, but may require either diversity or a suite of particular microbial genes [2]. In some cases, colonization with single bacterial strains, even when closely related, can differentially affect host development [8, 9]. For example, in the fruit fly Drosophila melanogaster, monoassociation with Lactobacillus brevis, but not L. plantarum, restored normal locomotor behavior [9].

In a larval zebrafish model, a powerful system for assessing the function of host-microbe interactions [10, 11], heat-killed bacteria and selected microbe-associated molecular patterns were both insufficient to modify axenic-related hyperactivity [4]. Furthermore, rather than restore normal behavior, monoassociation with either Aeromonas veronii or Vibrio cholerae resulted in a hypoactive phenotype [4]. In a different study, monoassociation of larval zebrafish with L. plantarum resulted in control-like locomotion, but reduced cortisol levels following stress, similar to the microbe-free phenotype [5]. Finally, monoassociation with either Escherichia coli or Bacillus subtilis at different doses only partially reduced axenic-like hyperactivity [12].

Importantly, differences in behavioral effects arising from monoassociation experiments could be mediated, in part, by bacterial abundance, and yet bacterial load is rarely quantified. To our knowledge, confirmation of bacterial identity by sequencing following monoassociation experiments has yet to be applied in zebrafish. A lack of confirmatory sequencing introduces a degree of experimental uncertainty in behavioral manifestations triggered by waterborne inoculation with specific strains of bacteria.

Here, to gain further insight into the mechanisms through which microbes mediate neurobehavioral development and the role of both microbial diversity and absolute abundance, we conventionalized axenic zebrafish with single or mixed host-associated bacterial isolates, as well as with an environmental bacterial isolate. We conducted whole genome sequencing on host-associated isolates and then used 16S rRNA gene sequencing to confirm genus-level identities at the end of the study. Finally, we used droplet digital PCR to quantify the total concentration of bacteria within each host-associated group. Our data suggest the importance of both bacterial identity and intestinal colonization on neurobehavioral development, and further, highlight the importance of understanding host-microbial interactions in animal models used for genetic and toxicity screens.


Individual host-associated isolates differ in the ability to block axenic-related hyperactivity

Zebrafish host-associated isolates were obtained from conventionally colonized zebrafish and used for subsequent association experiments to examine the effects of single and mixed isolates on locomotor behavior (Fig. 1). For host-associated single and mixed isolates, there were significant effects from time (p < 0.0001), colonization status (p < 0.0001), and the interaction between time and status (p < 0.0001) on distance moved (Fig. 2a). Further, there was a significant effect of colonization status on mean distance moved during both the light (p = 0.00048) and dark (p = 5.75 × 10− 13) periods. Because light and dark period behavior were correlated (p = 3.8 × 10− 14, R2 = 0.19) and given our previous observations that the most robust behavioral changes are observed after the light to dark transition [4], we report further behavioral results only for the dark period. As reported previously [4], axenic zebrafish were hyperactive relative to conventionally colonized (CC) or axenic colonized on day 1 (AC1) control groups (Fig. 2b). For host-associated isolates, we found that phenotypically control-like locomotor behavior occurred (i.e., there was no significant difference compared to either CC or AC1 in distance moved) following monoassociation with Aeromonas veronii, Vibrio cholerae, or V. metoecus, or with either defined mixture (Fig. 2b). In contrast, hyperactivity was only partially blocked by monoassociation with Comamonas testosteroni or Delftia tsuruhatensis and was not blocked by monoassociation with Acinetobacter venetianus (Fig. 2b).

Fig. 1
figure 1

A zebrafish experimental system with different colonization statuses. To generate host-associated bacterial isolates, visually pure isolates were grown from homogenized zebrafish larvae and then underwent whole genome sequencing. To examine the effects of colonization status on behavior, larvae were generated to be conventionally colonized (CC), axenic (AX), or axenic colonized on day 1 with either aquaculture facility microbes (AC1) or isolated single or mixed bacterial species (AC1:X). At 10 days post-fertilization, larvae were assayed for behavior, bacterial identity (i.e.16S rRNA gene sequencing), and bacterial concentration (i.e. 16S rDNA droplet digital PCR). Purple hexagons indicate heterogeneous aquaculture facility microbes. Turquoise ovals indicate strains of microbes isolated from zebrafish used for association experiments

Fig. 2
figure 2

Association with select host-associated bacterial species variably blocked the hyperactive axenic behavioral phenotype (n = 4–8 larvae per flask from 3 to 7 flasks per cohort). Distance moved every 2 min in the light and dark periods (a; error bars = standard error of the mean) and mean distance moved every 10 min in the dark period (b). Different letters indicate significant differences (p < 0.05). CC = conventionally colonized; AC1 = axenic colonized on day 1; AX = axenic; AC1:DM1 = A. venetianus + V. metoecus + C. testosteroni + D. tsuruhatensis + A. veronii; AC1:DM2 = A. venetianus + C. testosteroni + D. tsuruhatensis + A. veronii

16S rRNA gene sequencing confirmed bacterial isolates and revealed potential bacteria in axenic cohort

Using 16S rRNA gene sequencing, we assessed the genus level identity of the bacterial isolates after the 10-day study period. For larvae that were monoassociated, we found that the zebrafish-associated genus matched the inoculant species identity with an average of 85% relative read abundance (Fig. 3a; 95% AC1:A. venetianus; 93% AC1:A. veronii 70% AC1:C. testosteroni; 79% AC1:D. tsuruhatensis; 82% AC1:V. cholerae; 89% AC1:V. metoecus). We further examined within-group differences by calculating the distance to centroid on Bray-Curtis dissimilarity scores (Fig. 3b; F = 10.99, p = 1.0 × 10− 4). There was significantly higher heterogeneity within CC, AC1, and AX compared to all association treatment groups (Fig. 3b). To visualize community differences in diversity between CC, AC1, and AX larvae, we created heat trees. These depict the relative abundance of different OTUs using the color and size of the nodes in a taxonomic tree [13]. A visual comparison of the heat trees shows the complex community composition in CC and AC1, while AX larvae are predominantly colonized by two related genera (Fig. 3c).

Fig. 3
figure 3

Results of 16S rRNA gene sequencing. Analysis revealed complex bacterial communities in CC and AC1 larvae and verified genus-level identities of conventionalized groups at the end of each experiment (a). Distance to centroid on a Bray-Curtis dissimilarity matrix (b) revealed high within-group heterogeneity for CC, AC1, and AX, while heterogeneity was low for larvae colonized with isolates. Heat trees (c) show differences in microbial diversity between CC, AC1, and AX larvae. CC = conventionally colonized; AC1 = axenic colonized on day 1; AX = axenic; AC1:DM1 = A. venetianus + V. metoecus + C. testosteroni + D. tsuruhatensis + A. veronii; AC1:DM2 = A. venetianus + C. testosteroni + D. tsuruhatensis + A. veronii

Notably, we detected relatively high read counts for two bacterial operational taxonomic units (OTUs; OTU6 and OTU8) in the AX larvae, identified at the order level as Cytophagales (Fig. 3a, c). These OTUs, as well as an additional Cytophagales (OTU7), were also present in axenic larvae colonized with single isolates or defined mixtures but were not detected in CC or AC1 larvae, nor were they present in the negative controls or fish facility water samples (Supplemental data). To examine whether the presence of these three OTUs may have affected zebrafish locomotor behavior in axenic or monoassociated larvae, we used Bayesian ANCOVA with ‘colonization status’ as the fixed effect and included the three OTUs as separate covariates. The best model was ‘colonization status’ alone, with a Bayes factor of 5.5 × 105, which was 1.29 times better than the second-best model ‘colonization status + OTU8’ (Bayes factor of 4.2 × 105). In an additional analysis, including the three OTUs in the null model resulted in a Bayes factor of 27 for ‘colonization status,’ suggesting strong evidence for the alternative hypothesis that the presence of these OTUs did not affect locomotor behavior.

Bacterial loads differ between groups

Using 16S rRNA gene ddPCR measurements, we normalized relative abundance to estimate the cell-equivalent total bacterial load per larvae, calibrated with the number of 16S rRNA genes per genome. We found that bacterial load differed significantly between groups (Fig. 4a; p = 2.29 × 10− 4) and that this difference related to locomotor behavior. The cell-equivalent concentration of bacteria per larvae was weakly but significantly correlated with mean distance moved across all groups (p = 0.017, R2 = 0.11; Fig. 4b). Furthermore, accurately measuring bacterial concentrations revealed that the two OTUs (OTU6 and OTU8) detected by 16S rRNA gene sequencing were indeed present at a high concentration in AX larvae, at levels that were not significantly different from most treatment groups (Fig. 4a). Lastly, we examined how well the relative abundance data (i.e., 16S rRNA gene sequencing) reflected the absolute (i.e., ddPCR) bacterial load. Across groups, the log sum of the read counts was positively correlated with the log of the 16S rRNA gene concentration data (p = 2.8 × 10− 10, R2 = 0.62) (Fig. S1).

Fig. 4
figure 4

Results of 16S rRNA gene droplet digital PCR. Analysis revealed group differences in bacterial load per larvae (a). These differences correlate with behavior (b; p = 0.017, R2 = 0.11); CC = conventionally colonized; AC1 = axenic colonized on day 1; AX = axenic; AC1:DM1 = A. venetianus + V. metoecus + C. testosteroni + D. tsuruhatensis + A. veronii; AC1:DM2 = A. venetianus + C. testosteroni + D. tsuruhatensis + A. veronii

Monoassociation with environmentally-associated isolate does not block axenic-related hyperactivity

Because we observed variation in the ability of different host-associated bacteria to block hyperactivity, we hypothesized that intestinal colonization may be necessary for neurobehavioral development. To test this, locomotor activity in zebrafish monoassociated with the environmentally-associated isolate Rheinheimera makueensis was compared to CC, AC1, and AX larvae. Overall, significant effects of time (p < 0.0001), colonization status (p < 0.0001), and the interaction between time and status (p < 0.0001) on distance moved were observed (Fig. 5a). Further, for mean distance moved, there was a significant effect of colonization status (p = 8.9 × 10− 7), whereby AX larvae and those exposed to R. makueensis were hyperactive compared to both CC and AC1 groups (Fig. 5b).

Fig. 5
figure 5

Monoassociation with an environmentally-associated bacteria failed to block axenic-like hyperactivity. Distance moved every 2 min in the light and dark periods (a; error bars = standard error of the mean) and mean distance moved every 10 min in the dark period (b). CC = conventionally colonized; AC1 = axenic colonized on day 1; AX = axenic

Isolates that block axenic-related hyperactivity were more abundant than other isolates in defined mixtures

To examine the effect of defined mixtures on bacterial growth patterns, we normalized relative read abundance (i.e., 16S rRNA gene sequencing) by absolute abundance (i.e., 16S rRNA gene ddPCR) for individual species, calibrated with 16S rRNA copies per cell [14]. This resulted in estimates of normalized bacterial loads for each species within both the monoassociation and defined mixture experiments (Table 1). In the monoassociated zebrafish larvae, A. veronii had the highest load followed by C. testosteroni > D. tsuruhatensis > V. metoecus > A. venetianus. For larvae associated with Defined Mixture 1, the estimated bacterial loads for each species were lower than the maximum values of the monoassociated zebrafish, particularly for D. tsuruhatensis (5.66 vs 4.55) and C. testosteroni (5.90 vs 4.15). Ranked by load within Defined Mixture 1, A. veronii was the most abundant, followed by V. metoecus > A. venetianus > D. tsuruhatensis > C. testosteroni. Notably, the two highest-ranked taxa in this mixture, A. veronii and V. metoecus, were the species that fully blocked axenic-related hyperactivity (in addition to V. cholerae) in monoassociation experiments (Fig. 2b). Compared to Defined Mixture 2, the presence of V. metoecus in Defined Mixture 1 was associated with a 0.31 ± 0.05 decrease in the other species (species cell-equivalent log10 mean in Defined Mixture 2 – Defined Mixture 1). For Defined Mixture 2, where V. metoecus was not added, bacterial loads of A. veronii and A. venetianus were estimated at comparable loads to those of the monoassociated zebrafish, while D. tsuruhatensis and C. testosteroni showed significant log reductions, with the ranks as follows: A. veronii > A. venetianus > D. tsuruhatensis > C. testosteroni (Table 1). Generation time and growth rate values for monoassociation and mixture experiments are provided in Supplemental Materials.

Table 1 Normalized cell-equivalent log10 bacterial load values for monoassociation and mixture experiments


The interface between host cells and the microbiome is not well understood. While microbe-free animals consistently show altered behavioral phenotypes, including hyperactivity, it is unclear how bacterial colonization interacts with host cells to result in control-like neurobehavioral development. Here, we found variation in the ability of zebrafish-associated microbes to block the hyperactivity phenotype, as measured by distance moved following a light-to-dark transition. Monoassociation with three bacterial isolates resulted in control-like locomotor behavior, two others resulted in an intermediate level of hyperactivity and one failed to block hyperactivity entirely. Association with mixtures of four or five bacterial species resulted in larvae with behavioral phenotypes reflective of those that were conventionalized with a heterogeneous microbiota obtained from the aquaculture facility. In general, these results support those of previous studies in mice and flies in which colonization during development with single strains of bacteria can have different effects on behavior [9, 15]. Monoassociation with an environmentally-associated bacteria did not block the hyperactivity phenotype, but rather these larvae behaved similarly to those in the axenic cohort. This suggests that colonization with microbes that evolved in the host intestinal tract may be required to allow for control-like neurobehavioral development.

We detected two OTUs (OTU6 and OTU8) present in the axenic cohort, as well as a third (OTU7) also present in the monoassociated larvae. These OTUs are in the phylum Bacteriodetes, class Cytophagia, order Cytophagales. Species in this order can be opportunistic and are found mainly in habitats rich in organic material [16]. Cytophagia has previously been detected in the zebrafish intestine [17, 18] and further, was shown to increase in relative abundance in the gut in response to a high fat diet, correlated with markers of inflammation [19]. Because the co-occurrence patterns between OTU6 and OTU8 differed from those of OTU7, we speculate they have different roles. Given that OTU7 co-occurs with the known intestinal bacteria and given its absence in the axenic cohort, we suggest this OTU7 may be a predator in the intestine, as some species within Cytophagales have been identified as facultative predators [20]. They have a gliding motility and biosynthesize metabolites that allow them to extract nutrients by preying on other intestinal bacteria. In contrast, in the axenic larvae, OTU6 and OTU8 were generally the only OTUs detected. These bacteria may be metabolizing organic material produced in the zebrafish gut. Importantly, we did not detect Cytophagales in the larvae with complex microbial communities (CC and AC1). Given that the total read counts were higher in the monoassociated larvae, the lack of detection is unlikely to reflect sequencing limitations. Rather, the complexity and stability of the microbiota within CC and AC1 may serve a protective role and prevent the growth of Cytophagales species. An alternative but not mutually exclusive possibility is that particular species within the communities have adhesion factors that provide protective functions [21].

Given the unexpected presence of the Cytophagales bacteria, we were interested in whether they may have affected larval zebrafish behavior. Using a Bayesian modeling approach, we found a lack of support for the hypothesis that these OTUs affected locomotor behavior in our experiments. Furthermore, the hyperactive behavioral phenotype in the axenic cohort was nearly identical to that observed in our previous studies where the amount of bacteria measured in axenic larvae was significantly lower [4]. The other host-associated strain of bacteria that did not appear to block axenic-related hyperactivity was A. venetianus. The role of this bacteria in the zebrafish is unknown. It may be primarily skin-associated, as A. venetianus has been identified as a skin-associated microbe in juvenile snook, Centropomus undecimalis [22], and as an antifungal skin-associate in the frog, Espadarana prosoblepon [23], though other Acinetobacter have been detected in the zebrafish gut [24, 25]. In our experiments, larvae monoassociated with A. venetianus had the lowest bacterial abundance compared to all other cohorts, and further, A. venetianus was the least impacted by competition in the defined mixtures. In addition, there were very few reads for OTU7 (order Cytophagales) in the A. venetianus-associated cohort. Assuming this OTU is present as a predator in the intestine, as outlined above, these lines of evidence all support the proposition that A. venetianus is primarily inhabiting the skin of larval zebrafish. Lastly, the other bacteria that did not prevent axenic-related hyperactivity was the environmentally-associated R. makueensis, originally isolated from ocean water. However, we note that environmental isolates are not necessarily unable to colonize the zebrafish gut [26] and that additional experiments will be needed to demonstrate to what extent colonization may occur. It may also be worthwhile to examine locomotor behavior in multiple cohorts of axenic larvae serially conventionalized with R. makueensis, as it has been shown that microbes not previously known to colonize the zebrafish intestine can evolve this ability [26].

In contrast to the species that did not appear to affect behavior, evidence suggests that those that did mediate neurobehavioral development colonize host intestinal tracts. The bacteria species that fully blocked axenic dark-phase hyperactivity included V. metoecus, V. cholerae, and A. veronii, all of which are common inhabitants of the zebrafish gut. The colonization dynamics of both Vibrio and A. veronii have been well-studied in the zebrafish intestine [27]. Monocolonization with D. tsuruhatensis and C. testosteroni resulted in an intermediate behavioral phenotype. D. tsuruhatensis has been found in the gut of the grass carp, Ctenopharyngodon idellus [28], and at the genus level in the zebrafish gut [25]; C. testosteroni has previously been identified in the zebrafish intestine [29]. Notably, while colonization of the host intestinal tract may be necessary for neurobehavioral development, critical follow up experiments will be required to determine specific microbially-associated molecules and host responses required to establish the colonized motility phenotype.

In addition to intestinal colonization, there may be a relationship between absolute bacterial abundance and locomotor behavior. We found a weak but significant negative correlation between locomotor behavior and the cell-equivalent bacteria per larvae, whereby larvae with more bacteria show a pattern of reduced hyperactivity. In previous work, axenic zebrafish exposed at 1 day post-fertilization to 100 or 500 cells/mL of A. veronii and 5000 or 25,000 cells/mL of V. cholerae induced hypoactivity compared to colonized controls [4]. In the current study, exposure to 100 cells/mL of these two bacteria resulted in behavior that did not differ significantly from either colonized control group (CC or AC1), though notably exposure to these two strains did result in the lowest activity of all groups. In our previous study, it is possible that relatively high doses of bacteria added to the aquarium water resulted in a rate of bacterial immigration into the host that exceeded emigration, potentially causing sickness behavior [30].

Host-associated microbiota function through the dynamics of community ecology [31]. When microbe-free zebrafish are inoculated with bacterial isolates, initially there is little resource competition, competitive exclusion, or predation in the intestinal tract, thus bacteria can grow rapidly [27, 32]. In general, the larvae associated with isolates and defined mixtures had a higher abundance of bacteria compared to the larvae with complex communities. As suggested by the presence of Cytophagales, this may reflect a state of dysbiosis. In the defined mixture experiments, we found that the two host-associated isolates that fully blocked axenic-related hyperactivity performed better in the mixing experiments (A. veronii and V. metoecus). In other words, while the same number of cells were added to flasks containing axenic zebrafish, these two species yielded the highest normalized abundance and therefore the highest growth rates after 10 days. Competition experiments between A. veronii and V. cholerae have shown intestinal motility can selectively reduce the abundance of A. veronii [33]. However, with mixtures of four or five species, higher-order interactions weaken competition and facilitate diversity [34]. While it may be informative to examine longer-term outcomes in the larvae associated with single isolates or mixtures, these data indicate that species that grow rapidly in isolation outcompeted less robust taxa present in the defined mixtures.


There is now sufficient evidence to suggest that the microbiome modulates the development of brain and behavior across taxa. As such, moving forward it will be critical to consider the microbiome in the design of genetic and toxicity screens for factors that affect zebrafish neurodevelopment. For example, there is a two-way interaction between host microbiota and xenobiotic exposure, whereby chemicals in the environment can modify the composition of the microbiome, and the microbiome (i.e. all of the microorganisms and their encoded genes and associated functions that colonize a host organism), in turn, can bioactivate or inactive xenobiotics [35]. Further, as our data suggest, intramicrobiome interactions are also important to consider. Given that individual bacterial species have different effects on locomotor behavior, these interactions could substantially affect the results of toxicity screens. In particular, there are significant differences between zebrafish laboratory facilities (e.g. strains, water source and treatment methods, and food sources) that harbor the capacity to significantly alter aquaculture microbiota, even between days within the same facility [36]. Thus, consideration should be given not only toward consistently sampling and reporting the bacterial species present in zebrafish experiments using sequencing methods, but also toward accurately measuring their total abundance. While best practices for the determination of zebrafish contamination by culturing flask water have been described [37], our data from the axenic cohort demonstrate that culturing methods and visual inspection of microbial growth can be inadequate. Thus, there is need for standardization in both ontologies and methodological approaches necessary to label an animal “axenic.” In general, it is possible that variation in study results within and across laboratories represents underappreciated manifestations of host-microbial interactions.


Zebrafish husbandry

All experiments involving zebrafish were approved by the Institutional Animal Care and Use Committee at the U.S. EPA National Health and Environmental Effects Research Laboratory and performed in accordance with appropriate guidelines and regulations. A mixed wild type adult zebrafish line (Danio rerio) was used for all studies as previously described [36]. Zebrafish adults were housed on a recirculating system with conditioned reverse osmosis water in 6 L tanks at an approximate density of 8 fish/L, maintained on a 14:10 h light:dark cycle at 28.5 °C (pH = 7.5; conductivity = 800 μS), and fed Gemma Micro (Skretting) once daily and shell free E-Z Egg (Brine Shrimp Direct) twice daily Monday-Friday. On weekends, zebrafish were fed both food sources once daily. For spawning, ~ 80 adults were placed in 10 L angled static breeding tanks overnight. The following morning, adults were transferred to new bottom tanks containing treated reverse osmosis water (60 mg/l sodium bicarbonate and 0.4 g/l Crystal Sea Bioassay Formula Marine Mix) and embryos were collected after 45 min.

Bacterial isolation

On 0 days post fertilization (dpf), conventionally colonized zebrafish embryos were loaded into T25 tissue culture flasks (two flasks, n = 15 embryos/flask) containing 10% Hank’s balanced salt solution (HBSS) and incubated at 26 °C. On 1 dpf, an 80% media change was performed, and flasks were housed statically through 6 dpf. From 6 to 9 dpf, media (10% HBSS) was renewed daily and larvae were fed Gemma Micro 75 (Skretting). Ten dpf larvae were anesthetized on ice then homogenized and serial dilutions were prepared (1:10, 1:100, and 1:1000). One hundred uL was cultured on tryptic soy agar (TSA) and Luria-Bertani (LB) plates (n = 4). Plates were wrapped in Parafilm and incubated overnight at 26 °C, then colonies were gently spread to new plates and the process was repeated until visually pure colonies were obtained. Isolates were grown in 5 mL LB broth and incubated overnight at 26 °C with shaking. Optical density (OD) readings were obtained using Nanodrop2000c (OD600 X 10) to identify growth phases and approximate bacterial concentrations (cells/mL) ( Cultures were stored in glycerol stocks at − 80 °C until use.

Whole genome sequencing of host-associated bacterial isolates

We conducted whole genome sequencing on the five isolates of unknown species, which revealed their identity as the following: Acinetobacter venetianus, Vibrio metoecus, Comamonas testosteroni, Delftia tsuruhatensis, and Aeromonas veronii. For detail on genome sequencing, see Keely et al. (in preparation). Briefly, DNA was extracted from bacterial cultures using a ZR-Duet™ DNA/RNA MiniPrep Kit (Zymo Research, #D7003), and samples were stored at − 80 °C. Paired-end indexed libraries were prepared from genomic DNA of each isolate using the Nextera XT Library Prep Kit (Illumina, Inc., San Diego, CA) and Nextera XT Index Kit (Illumina, Inc.) per the manufacturer’s instructions. Libraries were pooled by volume, supplemented with PhiX Control v3 (Illumina, Inc.) to 5% (v/v) and sequenced using the MiSeq System (Illumina, Inc.) with the 500-cycle MiSeq Reagent Kit v2 (Illumina, Inc.) as instructed by the manufacturer. Whole genome sequencing analysis was performed in PATRIC v3.6.2 [38]. Specifically, sequence reads were trimmed and filtered using TrimGalore [39], assembled using SPAdes [40], and classified using MASH [41].

Axenic derivation and association with isolates

The experimental design is summarized in Fig. 1. Conventionally colonized (CC) larvae were generated by collecting embryos in untreated fish facility water at 0 dpf. Axenic (AX) larvae were derived following a series of treatments with antibiotics, poly (vinylpyrrolidone)-iodine (PVP-I) solution, and bleach, as previously described [35, 42]. Briefly, embryos were resuspended in 0.2 μm filter-sterilized 10% Hanks’ Balanced Salt Solution (FS-10% HBSS) containing antibiotics (amphotericin B (0.25 μg/ml), kanamycin (5 μg/ml), ampicillin (100 μg/ml), and gentamycin (0.05 mg/ml) for 4 h at 26 °C. At 5 h post fertilization, embryos were treated with 0.5% PVP-I (CASRN 25655–41-8) for 2 min and 0.05% bleach for 20 min and sorted into sterile T25 tissue culture flasks (25 mL of 10% HBSS).

As a control for the derivation process, a subset of axenic embryos was conventionalized with 10 ml of filtered (5 μm) aquaculture facility water (containing microbiota) at 1 dpf to generate the axenic colonized on day 1 (AC1) cohort. In addition, to increase the similarity between CC and AC1 cohorts, CC flasks also received 10 ml of aquaculture facility water on day 1. This aquaculture facility water was syringe-filtered using a sterile 5 μm filter to remove debris. In addition, a portion of 5 μm filtered fish facility water was syringe-filtered using a sterile 0.2 μm filter to generate microbe-free fish facility water. At 1 dpf, CC and AC1 flasks received an 80% media change consisting of 10 ml of microbe-containing 5 μm filtered fish facility water and 10 ml of sterile 10% Hanks’ Balanced Salt Solution (HBSS). In comparison, AX flasks received 10 ml of sterile 0.2 μm filtered fish facility water and 10 ml of sterile 10% HBSS. All flasks were housed statically through 6 dpf. From 6 to 9 dpf, all flasks underwent a daily 80% media change and received 75 kilogray gamma-irradiated Gemma Micro 75 at a final concentration of 0.04% (vol/vol) to eliminate microbial contributions from the food source.

To test sterility, at 1 and 10 dpf, two tryptic soy agar (TSA) plates (Sigma, #22091) were inoculated with 10 μl of media from each flask and examined for aerobic and anaerobic growth. For the AX-derived cohort, additional sterility testing for aerobic growth was done at dpf 6, 7, 8, and 9. At 10 dpf, the sterility of media from axenic flasks was further tested by inoculating 100 μL of flask media into tubes of Nutrient Broth (Sigma, #70122), Brain Heart Infusion Broth (Sigma, #53286), or Sabouraud Dextrose Broth (Sigma, #S3306). Plates and tubes were incubated at 26 °C under aerobic and anaerobic conditions for at least 7 days. One AX flask was found to be contaminated and was excluded from the study.

To generate larvae that were colonized with single or mixed bacteria isolates, an additional subset of axenic embryos were exposed to 100 cells/ml of A. venetianus, V. metoecus, C. testosteroni, D. tsuruhatensis, A. veronii, V. cholerae (ZWU0020), or one of two defined mixtures (AC1:DM1 = A. venetianus + V. metoecus + C. testosteroni + D. tsuruhatensis + A. veronii; AC1:DM2 = A. venetianus + C. testosteroni + D. tsuruhatensis + A. veronii). Twenty cells/mL per isolate were used for AC1:DM1 and 25 cells/mL per isolate were used for AC1:DM2 such that there were 100 cells/mL total used for all association groups. In a separate experiment following the same derivation methods, axenic larvae were exposed to 100 cells/mL of environmentally-derived Rheinheimera makueensis (AC1:R. makueensis; ATCC, Manassas, VA; sp. Nov.; Strain KH87; batch # 63381256; originally isolated off-shore from Kaneohe Bay, Oahu Hawaii). We selected R. makueensis because Rheinheimera was a genus level taxon present in previous studies [35, 36]. As our strategy to isolate zebrafish-associated microbes did not yield Rheinheimera species, we opted to use a commercially available strain of Rheinheimera to determine whether monoassociation with an environmentally derived strain of Rheinheimera could block axenic hyperactivity.

Behavior testing

At 10 dpf, larvae were assessed by visual inspection. Larvae that appeared morphologically normal were removed from flasks using a sterile transfer pipet and placed into 48-well plates containing 500 μl of FS-10% HBSS per well. Plates were sealed, wrapped in Parafilm, and placed in the dark in a temperature-controlled behavior testing room at 26 °C, for at least 2 h prior to testing. For testing, microtiter plates were placed on a Noldus tracking apparatus and locomotor activity was recorded. The light program consisted of a 20-min acclimation period in the dark (0 lx) followed by a test period comprised of a 10- or 20-min light period (18 lx) and a 10- or 20-min dark period (0 lx) for host- and environmentally-associated isolate trials, respectively [4]. In larval zebrafish, exposure to dark after a light period results in a characteristic increase in locomotor activity (light-seeking behavior) [43, 44]. A luminescence level of 18 lx for periods of 10 or 20 min follows behavioral protocols used previously [4, 42, 45]. Videos were analyzed using Ethovision software v.12 (Noldus Information Technology) as previously described [46].

DNA extraction and 16S rRNA gene sequencing

Microbial DNA was isolated by placing 500 μL of microbial growth collected from LB tubes (i.e. isolates used for monocolonizations) into Red RINO lysis tubes (Next Advance, # REDR5). Drawing from the same flasks from which larvae were taken for behavioral testing, we also examined community identities in all host-associated colonization groups. Whole zebrafish larvae at 10 dpf (3–7 biological replicates, 4–8 larvae per replicate, average of 7.75) were loaded into Red RINO tubes with 500 μL of FS-10% HBSS and anesthetized on ice. There were 3 flasks (i.e., biological replicates) used for CC, AC1, AC1:D. tsuruhatensis, and AC1:V. metoecus; 4 for AC1:A. venetianus, AC1:A. veronii, AC1:C. testosteroni, AC1:V. cholerae, AC1:DM1, and AC1:DM2; and 7 flasks for AX-derived larvae. A higher number of AX flasks was used in order to provide a buffer for flask loss due to contamination [4, 35, 36]. We chose to pool individual larvae to capture the heterogeneity inherent in multiple animals from the same flask (i.e. reared in the same environmental conditions), rather than representing a given flask with only a single larva. Following homogenization of microbial or zebrafish samples using a Bullet Blender, DNA was isolated from each sample using a ZR-Duet™ DNA/RNA MiniPrep Kit (Zymo Research, #D7003), and samples were stored at − 80 °C. Total DNA yield was measured using the Qubit dsDNA High Sensitivity Assay Kit (ThermoFisher, #Q32851) and Qubit 2.0 Fluorometer (ThermoFisher). DNA sequencing of the 16S rRNA gene was done as previously described [4, 36, 42]. Briefly, 250 ng of DNA from each sample was added to PCR reactions along with barcoded primers specific for the V4 region of the 16S rRNA gene [47] and amplified with the Roche FastStart High Fidelity PCR System (Sigma-Aldrich, #4738292001). PCR reactions were run for 2 min at 95 °C, followed by 25 cycles of 30 s at 95 °C, 30 s at 55 °Cn and 1 min at 72 °C, with a 10 min final extension at 72 °C. Triplicate reactions were pooled and products were purified and normalized with the SequalPrep Normalization Plate Kit (ThermoFisher, #A1051001). The DNA library was sequenced on the Illumina MiSeq platform (MiSeq Reagent Kit v2; 500 cycles, Illumina, #MS-102-2003). Positive and negative PCR control reactions were run with every 30 samples and sequenced to assess sequencing error and potential PCR contamination. Positive controls consisted of a mixture of equal concentrations of genomic DNA and negative controls consisted of 10 mM Tris-HCl at pH 8.5, as used in the dilution of DNA extracts.

Analysis of 16S rRNA gene sequences

MiSeq reads were processed using mothur (v1.40) [48], as previously reported [4]. Reads were filtered by Phred quality (Q30 with 50 nucleotide window length) and removed if they failed to form complete contigs. Read-pair contigs were removed if they contained any ambiguous base calls, homopolymers with more than 8 nucleotides, or if they were greater than 275 nucleotides in length. Failure to align to the V4 region of the Silva 16S rRNA gene reference alignment also resulted in removal (v128). The read-pair contigs were denoised using a preclustering algorithm. The UCHIME algorithm of USEARCH software [49] was used to identify and remove chimeras. A Bayesian classifier and the Ribosomal Database Project (RDP, training set version 16) were used to classify the read-pair contigs with a minimum bootstrap of 80% [50]. Read-pair contigs that did not classify at the level of kingdom or that classified as Archaea, Eukaryota, chloroplasts, or mitochondria were removed from further analysis. An operational taxonomic unit (OTU) table was generated with rows and columns representing samples and bacterial taxa counts (binned at 3% dissimilarity) and their taxonomic assignments.

Droplet digital PCR (ddPCR) of 16S rRNA genes

Quantities of 16S rRNA genes in zebrafish samples and potential PCR inhibition were determined by ddPCR using the BioRad QX200 Droplet Digital PCR System as previously described [4]. The exogenous internal amplification control did not show inhibition during ddPCR. For CC, AC1, and AX, the number of bacteria per larvae was determined using the average number of 4.2 copies of 16S rRNA gene per bacterial genome [51]. For monoassociation and mixture experiments, we used the median number of 16S rRNA genes per genome (see below for details).

Statistical analyses

Analyses were done in R version 3.6.1. To analyze behavior data across time points, we used a mixed effects model (lme4 package [52]) with ‘distance moved’ as the dependent variable and ‘colonization status’ and ‘time’ as fixed effects and ‘plate’ as a random effect. ‘Colonization status’ consisted of the 11 different treatment groups. Because flask was not a significant factor in any model, we averaged the distance moved within each light and dark period for individual larvae. ‘Mean distance moved’ and ‘bacteria per larvae’ were non-normally distributed (Shapiro Wilk test), therefore, to examine pairwise differences, the Kruskal-Wallis test and Dunn posthoc analyses with the Benjamini-Hochberg method for p-value adjustment were used to determine significance. To examine within-group similarity, multivariate homogeneity of group dispersions was calculated (i.e., within-group distance-to-centroid) on a Bray-Curtis distance matrix using the function betadisper from the R package ‘vegan’ [53]. To examine significance between groups, Tukey’s ‘Honest Significant Difference’ method was applied [53]. To visually examine differences in community structure from the 16S rRNA gene sequencing data, we used the R library phyloseq [54]. The stacked bar plot was generated using OTUs that contributed at least 1% to any of the sample totals. Heat trees were created using the R package metacoder [13]. To visualize the concentration of ‘bacteria per larvae’, we used a log10 transformation. Using the transformed data, we conducted a linear regression with ‘mean distance moved.’ In addition, we conducted a linear regression of the log transformed ‘bacteria per larvae’ with a sum of the 16S rRNA gene sequencing read counts. To examine the potential role of the three contaminant OTUs on locomotor behavior, we used Bayesian ANCOVA in JASP version 0.12.2 [55]. We set ‘colonization status’ as the fixed effect and included the three OTUs as separate covariates (Cauchy prior scale parameter for fixed effects = 0.5; Cauchy prior scale parameter for covariates = 0.354). CC and AC1 were excluded. In an additional model, we included the three OTUs in the null and examined support for the alternative hypothesis. To estimate the number of cell-equivalent bacterial species for each corresponding OTU, the relative abundance was normalized using the ddPCR measurements and calibrated using the number of 16S rRNA genes per genome, written as:

$$ {B}_{ij}=\left(\frac{otu_{ij}}{\sum {otu}_{ij}}\right)\left(\frac{ddPCR_j}{16{S}_{ij}{b_{ij}}^{-1}}\right) $$

where Bij is the cell-equivalent bacterial quantity estimated by using the relative abundance of the ith otu divided by the total number of otu reads in sample j. The ddPCR quantity of 16S rRNA gene at 10 dpf (ddPCRij) was calibrated by the number of 16S rRNA genes (16Sij) per identified species (bij). Data on the median number of 16S rRNA genes per genome was obtained from Espejo and Plaza [56].

Availability of data and materials

All data generated during this study are included in this published article and its supplementary information files. In addition, the datasets are available by searching for the manuscript title at



Axenic colonized on day 1




Conventionally colonized


Defined mixture 1, containing A. venetianus + V. metoecus + C. testosteroni + D. tsuruhatensis + A. veronii


Defined mixture 2, containing A. venetianus + C. testosteroni + D. tsuruhatensis + A. veronii


Days post-fertilization


  1. Sampson TR, Mazmanian SK. Control of brain development, function, and behavior by the microbiome. Cell Host Microbe. 2015;17:565–76.

    Article  CAS  Google Scholar 

  2. Fung TC, Olson CA, Hsiao EY. Interactions between the microbiota, immune and nervous systems in health and disease. Nat Neurosci. 2017;20:145–55.

    Article  CAS  Google Scholar 

  3. Diaz Heijtz R, Wang S, Anuar F, Qian Y, Bjorkholm B, Samuelsson A, Hibberd ML, Forssberg H, Pettersson S. Normal gut microbiota modulates brain development and behavior. Proc Natl Acad Sci U S A. 2011;108:3047–52.

    Article  Google Scholar 

  4. Phelps D, Brinkman NE, Keely SP, Anneken EM, Catron TR, Betancourt D, Wood CE, Espenschied ST, Rawls JF, Tal T. Microbial colonization is required for normal neurobehavioral development in zebrafish. Sci Rep. 2017;7:11244.

    Article  Google Scholar 

  5. Davis DJ, Bryda EC, Gillespie CH, Ericsson AC. Microbial modulation of behavior and stress responses in zebrafish larvae. Behav Brain Res. 2016;311:219–27.

    Article  Google Scholar 

  6. Sudo N, Chida Y, Aiba Y, Sonoda J, Oyama N, Yu XN, Kubo C, Koga Y. Postnatal microbial colonization programs the hypothalamic-pituitary-adrenal system for stress response in mice. J Physiol-London. 2004;558:263–75.

    Article  CAS  Google Scholar 

  7. Erny D, de Angelis ALH, Jaitin D, Wieghofer P, Staszewski O, David E, Keren-Shaul H, Mahlakoiv T, Jakobshagen K, Buch T, et al. Host microbiota constantly control maturation and function of microglia in the CNS. Nat Neurosci. 2015;18:965.

    Article  CAS  Google Scholar 

  8. Nishino R, Mikami K, Takahashi H, Tomonaga S, Furuse M, Hiramoto T, Aiba Y, Koga Y, Sudo N. Commensal microbiota modulate murine behaviors in a strictly contamination-free environment confirmed by culture-based methods. Neurogastroent Motil. 2013;25.

  9. Schretter CE, Vielmetter J, Bartos I, Marka Z, Marka S, Argade S, Mazmanian SK. A gut microbial factor modulates locomotor behaviour in Drosophila. Nature. 2018;563:402.

    Article  CAS  Google Scholar 

  10. Baumler AJ, Sperandio V. Interactions between the microbiota and pathogenic bacteria in the gut. Nature. 2016;535:85–93.

    Article  CAS  Google Scholar 

  11. Bertotto LB, Catron TR, Tal T. Exploring interactions between xenobiotics, microbiota, and neurotoxicity in zebrafish. Neurotoxicology. 2020;76:235–44.

    Article  CAS  Google Scholar 

  12. Tan F, Limbu SM, Qian Y, Qiao F, Du ZY, Zhang ML. The responses of germ-free Zebrafish (Danio rerio) to varying bacterial concentrations, colonization time points, and exposure duration. Front Microbiol. 2019;10.

  13. Foster ZS, Sharpton TJ, Grunwald NJ. Metacoder: An R package for visualization and manipulation of community taxonomic diversity data. PLoS Comput Biol. 2017;13:e1005404.

    Article  Google Scholar 

  14. Jian C, Luukkonen P, Yki-Jarvinen H, Salonen A, Korpela K. Quantitative PCR provides a simple and accessible method for quantitative microbiota profiling. PLoS One. 2020;15:e0227285.

    Article  CAS  Google Scholar 

  15. Cryan JF, Dinan TG. Mind-altering microorganisms: the impact of the gut microbiota on brain and behaviour. Nat Rev Neurosci. 2012;13:701–12.

    Article  CAS  Google Scholar 

  16. Kampmann K, Ratering S, Kramer I, Schmidt M, Zerr W, Schnell S. Unexpected stability of Bacteroidetes and Firmicutes communities in laboratory biogas reactors fed with different defined substrates. Appl Environ Microb. 2012;78:2106–19.

    Article  CAS  Google Scholar 

  17. Yang HT, Zou SS, Zhai LJ, Wang Y, Zhang FM, An LG, Yang GW. Pathogen invasion changes the intestinal microbiota composition and induces innate immune responses in the zebrafish intestine. Fish Shellfish Immunol. 2017;71:35–42.

    Article  CAS  Google Scholar 

  18. Gaulke C, Martins M, Watral V, Kent M, Sharpton T. Parasitic infection by Pseudocapillaria tomentosa is associated with a longitudinal restructuring of the Zebrafish gut microbiome. BioRxiv. 2016;076596.

  19. Ye L, Mueller O, Bagwell J, Bagnat M, Liddle RA, Rawls JF. High fat diet induces microbiota-dependent silencing of enteroendocrine cells. Elife. 2019;8:e48479.

  20. Linares-Otoya L, Linares-Otoya V, Armas-Mantilla L, Blanco-Olano C, Crusemann M, Ganoza-Yupanqui ML, Campos-Florian J, Konig GM, Schaberle TF. Diversity and antimicrobial potential of predatory Bacteria from the Peruvian coastline. Mar Drugs. 2017;15:308.

  21. Rendueles O, Ferrieres L, Fretaud M, Begaud E, Herbomel P, Levraud JP, Ghigo JM. A new zebrafish model of oro-intestinal pathogen colonization reveals a key role for adhesion in protection by probiotic bacteria. PLoS Pathog. 2012;8:e1002815.

    Article  CAS  Google Scholar 

  22. Tarnecki AM, Brennan NP, Schloesser RW, Rhody NR. Shifts in the skin-associated microbiota of hatchery-reared common Snook Centropomus undecimalis during acclimation to the wild. Microb Ecol. 2019;77:770–81.

    Article  Google Scholar 

  23. Woodhams DC, Alford RA, Antwis RE, Archer H, Becker MH, Belden LK, Bell SC, Bletz M, Daskin JH, Davis LR, et al. Antifungal isolates database of amphibian skin-associated bacteria and function against emerging fungal pathogens. Ecology. 2015;96:595.

    Article  Google Scholar 

  24. Gaulke CA, Martins ML, Watral VG, Humphreys IR, Spagnoli ST, Kent ML, Sharpton TJ. A longitudinal assessment of host-microbe-parasite interactions resolves the zebrafish gut microbiome's link to Pseudocapillaria tomentosa infection and pathology. Microbiome. 2019;7:10.

    Article  Google Scholar 

  25. Roeselers G, Mittge EK, Stephens WZ, Parichy DM, Cavanaugh CM, Guillemin K, Rawls JF. Evidence for a core gut microbiota in the zebrafish. ISME J. 2011;5:1595–608.

    Article  CAS  Google Scholar 

  26. Lebov JF, Schlomann BH, Robinson CD, Bohannan BJM, Comstock LE. Phenotypic parallelism during experimental adaptation of a free-living bacterium to the Zebrafish Gut. mBio. 2020;11:e01519–20.

  27. Stephens WZ, Wiles TJ, Martinez ES, Jemielita M, Burns AR, Parthasarathy R, Bohannan BJM, Guillemin K. Identification of population bottlenecks and colonization factors during assembly of bacterial communities within the Zebrafish Intestine. mBio. 2015;6:e01163–15.

  28. Li H, Zhong QP, Wirth S, Wang WW, Hao YT, Wu SG, Zou H, Li WX, Wang GT. Diversity of autochthonous bacterial communities in the intestinal mucosa of grass carp (Ctenopharyngodon idellus) (Valenciennes) determined by culture-dependent and culture-independent techniques. Aquac Res. 2015;46:2344–59.

    Article  Google Scholar 

  29. Cantas L, Sorby JRT, Alestrom P, Sorum H. Culturable gut microbiota diversity in Zebrafish. Zebrafish. 2012;9:26–37.

    Article  CAS  Google Scholar 

  30. Robinson CD, Klein HS, Murphy KD, Parthasarathy R, Guillemin K, Bohannan BJM. Experimental bacterial adaptation to the zebrafish gut reveals a primary role for immigration. PLoS Biol. 2018;16:e2006893.

  31. Walter J, Maldonado-Gomez MX, Martinez I. To engraft or not to engraft: an ecological framework for gut microbiome modulation with live microbes. Curr Opin Biotechnol. 2018;49:129–39.

    Article  CAS  Google Scholar 

  32. Stephens WZ, Burns AR, Stagaman K, Wong S, Rawls JF, Guillemin K, Bohannan BJ. The composition of the zebrafish intestinal microbial community varies across development. ISME J. 2016;10:644–54.

    Article  Google Scholar 

  33. Wiles TJ, Jemielita M, Baker RP, Schlomann BH, Logan SL, Ganz J, Melancon E, Eisen JS, Guillemin K, Parthasarathy R. Host gut motility promotes competitive exclusion within a model intestinal microbiota. PLoS Biol. 2016;14:e1002517.

  34. Sundarraman D, Hay EA, Martins DM, Shields DS, Pettinari NL, Parthasarathy R. Quantifying multi-species microbial interactions in the larval zebrafish gut. bioRxiv. 2020.

  35. Weitekamp CA, Phelps D, Swank A, McCord J, Sobus J, Catron TR, Keely SP, Brinkman FS, Zurlinden TJ, Wheaton E, et al. Triclosan-selected host associated microbiota perform xenobiotic biotransformations in larval zebrafish. Toxicol Sci. 2019;172:109–22.

    Article  CAS  Google Scholar 

  36. Catron TR, Keely SP, Brinkman NE, Zurlinden TJ, Wood CE, Wright JR, Phelps D, Anneken Wheaton E, Kvasnicka A, Gaballah S, et al. Host developmental toxicity of BPA and BPA alternatives is inversely related to microbiota disruption in zebrafish. Toxicolog Sci. 2018;167:468–83.

  37. Melancon E, Gomez De La Torre Canny S, Sichel S, Kelly M, Wiles TJ, Rawls JF, Eisen JS, Guillemin K: Best practices for germ-free derivation and gnotobiotic zebrafish husbandry. In The Zebrafish - Disease Models and Chemical Screens. 2017: 61–100: Methods in Cell Biology].

  38. Wattam AR, Davis JJ, Assaf R, Boisvert S, Brettin T, Bun C, Conrad N, Dietrich EM, Disz T, Gabbard JL, et al. Improvements to PATRIC, the all-bacterial bioinformatics database and analysis resource center. Nucleic Acids Res. 2017;45:D535–D42.

    Article  CAS  Google Scholar 

  39. Krueger F. Trim Galore: a wrapper tool around Cutadapt and FastQC to consistently apply quality and adapter trimming to FastQ files, with some extra functionality for MspI-digested RRBS-type (Reduced Representation Bisufite-Seq) libraries. 2012.

  40. Bankevich A, Nurk S, Antipov D, Gurevich AA, Dvorkin M, Kulikov AS, Lesin VM, Nikolenko SI, Pham S, Prjibelski AD, et al. SPAdes: a new genome assembly algorithm and its applications to single-cell sequencing. J Comput Biol. 2012;19:455–77.

    Article  CAS  Google Scholar 

  41. Ondov BD, Treangen TJ, Melsted P, Mallonee AB, Bergman NH, Koren S, Phillippy AM. Mash: fast genome and metagenome distance estimation using MinHash. Genome Biol. 2016;17.

  42. Catron TR, Swank A, Wehmas LC, Phelps D, Keely SP, Brinkman NE, McCord J, Singh R, Sobus J, Wood CE, et al. Microbiota alter metabolism and mediate neurodevelopmental toxicity of 17 beta-estradiol. Sci Rep. 2019;9:1–15.

  43. Burgess HA, Granato M. Modulation of locomotor activity in larval zebrafish during light adaptation. J Exp Biol. 2007;210:2526–39.

    Article  Google Scholar 

  44. Fernandes AM, Fero K, Arrenberg AB, Bergeron SA, Driever W, Burgess HA. Deep brain photoreceptors control light-seeking behavior in Zebrafish larvae. Curr Biol. 2012;22:2042–7.

    Article  CAS  Google Scholar 

  45. Gaballah S, Swank A, Sobus JR, Howey XM, Schmid J, Catron T, McCord J, Hines E, Strynar M, Tal T. Evaluation of developmental toxicity, developmental neurotoxicity, and tissue dose in Zebrafish exposed to GenX and other PFAS. Environ Health Perspect. 2020;128:47005.

    Article  CAS  Google Scholar 

  46. Jarema KA, Hunter DL, Shaffer RM, Behl M, Padilla S. Acute and developmental behavioral effects of flame retardants and related chemicals in zebrafish. Neurotoxicol Teratol. 2015;52:194–209.

    Article  CAS  Google Scholar 

  47. Kozich JJ, Westcott SL, Baxter NT, Highlander SK, Schloss PD. Development of a dual-index sequencing strategy and Curation pipeline for analyzing amplicon sequence data on the MiSeq Illumina sequencing platform. Appl Environ Microb. 2013;79:5112–20.

    Article  CAS  Google Scholar 

  48. Schloss PD, Westcott SL, Ryabin T, Hall JR, Hartmann M, Hollister EB, Lesniewski RA, Oakley BB, Parks DH, Robinson CJ, et al. Introducing mothur: open-source, platform-independent, community-supported software for describing and comparing microbial communities. Appl Environ Microb. 2009;75:7537–41.

    Article  CAS  Google Scholar 

  49. Edgar RC. Search and clustering orders of magnitude faster than BLAST. Bioinformatics. 2010;26:2460–1.

    Article  CAS  Google Scholar 

  50. Wang Q, Garrity GM, Tiedje JM, Cole JR. Naïve Bayesian classifier for rapid assignment of rRNA sequences into the new bacterial taxonomy. Appl Environ Microb. 2007;73:5261–7.

    Article  CAS  Google Scholar 

  51. Vetrovsky T, Baldrian P. The variability of the 16S rRNA gene in bacterial genomes and its consequences for bacterial community analyses. PLoS One. 2013;8.

  52. Bates D, Sarkar D, Bates MD, Matrix L. The lme4 package. R Package Version. 2007;2:74.

    Google Scholar 

  53. Oksanen J, Blanchet FG, Kindt R, Legendre P, Minchin PR, O’hara RB, Simpson GL, Solymos P, Stevens MHH, Wagner H, Oksanen MJ: Package ‘vegan’. Community ecology package, version 2.9. 2013;1–295.

  54. Watson M, McMurdie PJ, Holmes S. Phyloseq: an R package for reproducible interactive analysis and graphics of microbiome census data. PloS One. 2013;8:e61217.

  55. Rouder J, Morey R, Speckman P, Province J. Default Bayes factors for ANOVA designs. J Math Psychol. 2012;56:356–74.

    Article  Google Scholar 

  56. Espejo RT, Plaza N. Multiple ribosomal RNA operons in Bacteria; their concerted evolution and potential consequences on the rate of evolution of their 16S rRNA. Front Microbiol. 2018;9:1232.

    Article  Google Scholar 

Download references


We thank Joan Hedge and the US EPA Zebrafish Facility for fish husbandry. We are grateful to David Thomas, Laura Dishaw, and Andrew Hotchkiss for their review of the manuscript. This manuscript has been reviewed by the US EPA and approved for publication. Approval does not signify that the contents reflect the views of the Agency, nor does mention of trade names or commercial products constitute endorsement or recommendation for use.


This work was supported by the US EPA Pathway Innovation Project Award to TT and by the US EPA Office of Research and Development. Open Access funding was enabled and organized by Projekt DEAL.

Author information

Authors and Affiliations



CAW participated in study design, analyzed data, prepared figures, and wrote the paper. AK participated in study design and performed zebrafish experiments. X-MH performed zebrafish experiments, analyzed data, and prepared figures. SK, NB, DP, and TC participated in study design. SK assisted with data analysis and interpretation. NB and EA performed 16S rRNA sequencing. SG, DP, and TC assisted with zebrafish experiments. TT designed the study, assisted with zebrafish experiments, analyzed data, and assisted with paper preparation. All authors reviewed the manuscript. The author(s) read and approved the final manuscript.

Corresponding author

Correspondence to Tamara Tal.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no 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: Fig. S1

. Bacterial abundance (log10) measured by ddPCR and total reads (log10) measured by 16S rRNA gene sequencing show a positive linear relationship (p = 2.8 × 10− 10, R2 = 0.62). Table S1. Generation, generation time, and growth rate for monoassociation and mixture experiments.

Additional file 2:

Data files and R script used to generate figures and conduct statistical analyses.

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

Weitekamp, C.A., Kvasnicka, A., Keely, S.P. et al. Monoassociation with bacterial isolates reveals the role of colonization, community complexity and abundance on locomotor behavior in larval zebrafish. anim microbiome 3, 12 (2021).

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: