Experimental bacterial adaptation to the zebrafish gut reveals a primary role for immigration

All animals live in intimate association with microorganisms that profoundly influence their health and development, yet the traits that allow microorganisms to establish and maintain host associations are not well understood. To date, most investigations aimed at identifying traits required for host association have focused on intrahost niches. Consequently, little is known about the relative contribution of extrahost factors such as environmental growth and survival and immigration into hosts from the external environment, as promoters of host association. To address this, we developed a tractable experimental evolution system that investigates both intra- and extrahost factors contributing to bacterial adaptation to the vertebrate gut. We passaged replicate lines of a zebrafish bacterial isolate, Aeromonas veronii, through populations of germ-free larval zebrafish (Danio rerio), each time using gut-associated Aeromonas populations to inoculate the aquatic environment of the next zebrafish population. We observed rapid increased adaptation to the host in all replicate lines. The initial adaptations present in early-evolved isolates did not increase intrahost fitness but rather enhanced both immigration from the environment and interhost transmission. Only in later-evolved isolates did we find evidence for intrahost-specific adaptations, as demonstrated by comparing their competitive fitness in the host genotype to which they evolved to that in a different genotype. Our results show how selection for bacterial transmission between hosts and their environment can shape bacterial-host association. This work illuminates the nature of selective forces present in host–microbe systems and reveals specific mechanisms of increased host association. Furthermore, our findings demonstrate that the entire host–microbe–environment system must be considered when identifying microbial traits that contribute to host adaptation.


Introduction
Animals are intimately associated with highly complex and dynamic communities of microbes that inhabit virtually every surface of their bodies. An explosion of research in recent decades has revealed that both the microbiota and the hosts they colonize have profound influences on the physiology and evolution of one another [1,2]. In vertebrates, it is especially evident that host-microbiota interactions within the gastrointestinal tract are fundamental to host health and development [3]. Wide-ranging surveys of the earth's microbial communities have shown that host-associated microbial communities are distinct from non-host-associated communities, suggesting that there are unique properties of host-microbe systems that drive microbiome assembly [4]. Elucidating the microbial traits that allow host colonization and support host-microbe association could shed light on the ecology of host-associated microbes across the entire symbiotic spectrum, from mutualism (e.g., beneficial microbes) to parasitism (e.g., pathogens). Knowledge of these traits may also provide insights into the development of microbiota-driven diseases and ultimately into innovative ways to treat or prevent them.
Several different approaches have been used in past studies to identify traits that mediate microbial association with vertebrate hosts, including forward genetic screens (e.g., transposon mutagenesis) [5][6][7][8], functional genomics (e.g., transcriptomics) [9,10], and approaches from comparative genomics [11][12][13]. However, these approaches all have substantial limitations; for example, comparative genomics approaches have limited ability to distinguish between differences arising from adaptation to the host and those accumulated over evolutionary time in the non-host environment. Here, we use experimental evolution, which allows for the observation of evolution in real time under controlled laboratory conditions, to identify traits that mediate host association. This approach has the advantage that changes in traits can be directly attributed to the selection imposed by the experimental conditions, and there is no inherent constraint on the spectrum of mutations or traits that can be selected. Furthermore, this approach focuses on a phenotypic-level characterization of the system; adaptive traits may be acquired through many different genetic mechanisms, but their shared phenotypes reveal the selective pressures that drive their evolution.
Historically, experimental evolution has been a powerful strategy to study how microorganisms adapt to in vitro environments [14,15]. In recent years, this approach has been expanded populations undergo dramatic population collapse and regrowth events, on average about once per day [32]. Finally, prior studies showing that this strain is readily outcompeted by other zebrafish gut isolates and that its average gut colonization level is below that of many other strains suggested that it had potential for further host adaptation [32,38].
The amount of genetic variation within a system is an additional parameter that strongly affects the rate of adaptation. Mutator strains carry mutations in the DNA replication or repair machinery, effectively increasing spontaneous mutation rates and, therefore, adaptive potential [39]. To increase the rate of adaptation in our system, we generated a mutator (ΔmutS) strain of Aer01. MutS is a component of the methyl-directed DNA mismatch repair system, and mutS mutants commonly arise in experimental evolution studies [39][40][41] and naturally evolving systems [42,43]. Allelic exchange was used to make an in-frame deletion in the mutS gene of a fluorescently tagged Aer01 strain. A fluctuation assay confirmed that the mutS mutant had an approximately 1,000-fold-increased mutation rate compared to the wild-type (WT) strain (S1 Table).
The trade-off for increasing adaptive potential by using a mutator strain is that background accumulation of mutations confounds identification of adaptive alleles. Since the primary goal of this study was to identify traits that confer host association, rather than the underlying genes, we chose to utilize the mutator phenotype as a way of accelerating evolution and facilitating adaptation in our model system.
Characterization of the colonization cycle. We expanded on previous work investigating Aer01 colonization dynamics by further characterizing key features of the Aer01 colonization cycle in our model system. Larval zebrafish hosts are inoculated by adding Aer01 to the water column, which consists of sterile embryo medium (EM), and thus, colonization begins with Aer01 immigrating into the host. We first determined the robustness of Aer01 colonization as a function of inoculum concentration, which could fluctuate during the evolution experiment. Gut colonization levels were not affected by inoculum size (over >4 log range), consistently reaching about 10 4 colony-forming units per gut (CFU/gut; S1A Fig). However, we found that independent of starting concentration, Aer01 is detected in the EM at about 10 5 CFU/ml by 20 hours post inoculation (S1B Fig). These data show that there is a defined carrying capacity both for populations in the host gut and also in the extrahost environment.
Additional colonization parameters important for adaptation in this system are (1) the size of the initial colonizing population and (2) the frequency of mixing between the extra-and intrahost environments over time. These are important because they influence the pool of genetic variants that will be sampled within this model system, which impacts the rate of adaptation. To estimate the size of the initial (first 24 hours of colonization) bottleneck, we investigated the colonization success of a tagged strain of Aer01 in the gut when introduced at various low frequencies in the inoculum. We found that when a 1:100 mixture of dTomatotagged Aer01 and WT Aer01 was added to the external medium, at 24 hours post inoculation the tagged strain was detected in 80% of fish (n = 54), whereas at a ratio of 1:300 it was only detected in 38% (n = 53) (S2 Fig). We used a binomial sampling model (see Materials and methods) to estimate a 24-hour bottleneck size of 194 ± 140 Aer01 cells, demonstrating that the founder population in a fish gut is composed of much more than a few cells but is still multiple orders of magnitude below the average gut colonization level.
Once the gut population reaches carrying capacity, there may be additional influx as bacteria are expelled from the intestinal tract and new cells migrate in from the environment. The ability of Aer01 to invade already-established populations has previously been reported [32]. We verified this finding and confirmed that members of a differentially tagged Aer01 population, when inoculated into a flask of precolonized Aer01 fish, are present in the gut after 24 hours, at an average ratio of about 2:1 (founder:invader) (S3 Fig). These data show that the larval host is continually sampling from environmental populations over time.
Larval zebrafish undergo rapid development over the 3-day time frame (4-7 days post fertilization [dpf]) during which we conducted our colonization experiments. We wanted to determine, at a broad level, if Aer01 colonization is affected by the stage of host development. To test this, we colonized fish at 4, 5, or 6 dpf and enumerated bacterial load after 24 hours by gut dissection and plating. The average CFU/gut was not significantly different across the three time points, demonstrating that Aer01 colonization level is consistent across these host developmental windows (S4 Fig). These results provide a detailed picture of the Aer01 colonization cycle in our model system ( Fig 1A). This cycle provides multiple opportunities for Aer01 to increase host association (i.e., the chance of being found in the gut). Gut-specific adaptations include increased growth rate, increased abundance (potentially via increased carrying capacity), and increased persistence via tolerance to stresses within the host environment, or resistance to expulsion by peristaltic activity of the gut. Adaptations related to environmental growth, abundance, or persistence could also increase host association because they would increase the probability of being sampled by the host throughout colonization. Finally, changes in the capacity for migration, either increased immigration or decreased emigration, could also improve host association.
Evolution passaging scheme. To identify traits that promote host association in our system, we serially passaged three replicate clonal mutator Aer01 "ancestral" populations through groups of 10-15 WT larval zebrafish (Fig 1B). For each passage, the Aer01 populations were inoculated into flasks of GF 4 dpf larval fish containing EM (the external medium in our system). After 3 days (i.e., 7 dpf), the fish were euthanized and their guts harvested by dissection to specifically select for gut-associated Aer01 populations. The guts were combined and homogenized to release the bacteria; then, an aliquot (approximately half of the total sample) was used as inoculum for the next group of 4 dpf GF larval fish, and the remaining gut sample was cryopreserved. This process was repeated for 22 passages.

Evolutionary adaptation of Aer01 to the larval zebrafish gut occurred quickly and reproducibly across replicate lines
Cryo-archived evolved population samples from passages 4, 8, 13/14, 18, and 22 were streaked for isolation on rich media, and freezer stocks of colony-purified isolates were made. These isolates were assayed for adaptation by competing them against a non-mutator, differentially tagged (green fluorescent protein [GFP]) Aer01 reference strain. For the competitions, ancestral or evolved isolates were mixed with the reference strain and inoculated into the EM of individual flasks of GF fish. Competitions were conducted under the same conditions as the evolution passages, for which strains were competed for 3 days (from 4 to 7 dpf), and all aspects of the colonization cycle are included. At 7 dpf, the fish were dissected and the guts plated to enumerate in vivo abundances of each strain via fluorescence microscopy. Competitive indices (CIs) were calculated by dividing the strain ratio (competitor:reference) at the end of the competition by the starting ratio (Fig 2). Adaptation occurred very quickly (by passage 4) in all three replicate lines, with greater than 10-fold average CIs for the passage 4 isolate from all three lines. This increase in competitive fitness is maintained or increased (approximately 100-to 1,000-fold averages) across isolates from later passages. The increased competitive fitness of the evolved isolates is not due to a general growth advantage, as there are no competitive advantages in vitro, in rich media (S5 Fig).

The competitive advantage of early-evolved isolates is dependent on the mode of inoculation
To determine how the initial adaptations observed in the passage 4 isolates improved fitness, we next investigated the stage of the colonization cycle impacted. To do this, we used an inoculation approach that would isolate where competition takes place, thereby distinguishing fitness specific to the within-host environment from that of extrahost factors. Microgavaging is a technique whereby a blunt needle is used to deliver substances-in this case, bacterial cultures -directly into the gut via the oral route [44]. Competition mixtures of ancestral or passage 4 evolved isolates and reference strain (as above) were inoculated directly into the guts of GF larval fish. Concurrently, the same inocula were introduced into the flasks of GF fish as in the previous competition experiments. For all three lines, the competitive advantage of flask-inoculated evolved isolates was significantly diminished when the strains were introduced directly into the gut (Fig 3A-3C). This loss of competitive advantage was apparent at both 0.5 hours and 5 hours post gavage (Fig 3A-3C). The matched CFU/gut data for the same fish show that Aer01 was gavaged at sufficiently low densities (approximately 5 × 10 2 CFU/gut) to allow for >10-fold increase in abundance during the 5-hour colonization (Fig 3D). If the CIs had increased between 0.5 and 5 hours post gavage, this would have indicated either a higher growth rate of the evolved strains or a better capacity to persist or survive within the gut environment. Since there were no differences in CIs between these time points, neither of these traits can be ascribed to the evolved lines. Instead, the competitive advantage of the evolved isolates is apparent when they are inoculated into the flasks (Fig 3A-3C), with median CIs of 11, 4, and 98 for evolved isolates from lines 1, 2, and 3, respectively. This suggests that the competitive advantage of these isolates is specific to initial colonization events, prior to reaching the gut environment. We wondered whether the competitive advantage was due to the ability to transmit from the mouth into the gut. To investigate this, we gavaged competition mixtures into the mouths of the fish rather than into the gut and then again measured the gut CIs after 5 hours. As with the gut gavages, the evolved strains when introduced into the mouth showed little to no significant advantage over the ancestral strain (Fig 3A-3C). In addition, the variation in these data is high, implying there is a variable-and at times tight-bottleneck imposed between the mouth and the gut such that in some fish one tag or the other dominates the Aeromonas gut population. Combined, these data suggest that the competitive advantage of the early-evolved isolates is driven by initial colonization of the host from the environment.

Evolved isolates have higher rates of immigration into the host gut
To further quantify differences between passage 4 evolved isolates and the ancestral strain, we conducted migration rate experiments. Here, fish were mono-associated (inoculated into the flask) with either the ancestor or a passage 4 evolved isolate. Every 45 minutes, groups of 8-10 fish were dissected and the guts plated to determine (1) the fraction of fish colonized at each time point ( Fig 4A) and (2) the Aeromonas abundance ( Fig 4B). All three evolved isolates were present in a higher proportion of fish at each time point than the ancestor, with close to 100% colonization by the evolved isolates even at the earliest time points, whereas the ancestor was not present in 100% of the fish even at the final sampling (about 5 hours post inoculation). Abundance of the evolved isolates in the guts was also higher than the ancestral strain, with, on average, approximately 10-fold higher CFU/gut at the first time point compared to the ancestor. These data do not reflect a difference in growth or survival in the flask EM, since abundance remained relatively constant for all strains throughout the experiment (S6 Fig).
In addition to assessing bacterial colonization of zebrafish directly from the flask medium, we wanted to test if the increased migration phenotype of the evolved isolates could affect fishto-fish colonization. To do this, we mono-associated fish with either the ancestor or the passage 4 isolate from line 3, washed and added two "donor" fish to a flask of GF "recipient" fish, and dissected and plated the recipient fish 12 hours later. Indeed, we found that the abundance of the evolved isolate was significantly higher than the ancestor in their respective recipient fish (S7 Fig).
We next developed a mathematical stochastic colonization and growth model (see Materials and methods for details) with which we could use the measured distribution of bacterial abundances across fish from the migration rate experiments to estimate the migration and growth rates of each strain. We consider a model in which migration into a fish by an individual bacterium is a stochastic event, with some probability per unit time. The statistics of the colonizers at any time are therefore governed by a Poisson distribution characterized by a mean time to entry, τ. Following entry into the host, the descendants of a colonizer obey logistic growth (initially roughly exponential, constrained by a finite carrying capacity, K), characterized by some growth rate r. We fit the experimental data (i.e., each set of CFU measurements at a particular measurement time t) to the above model of stochastic colonization and growth, determining the best-fit parameters τ and r. The datasets from each measurement time give independent estimates of τ and r; the average and standard error give the final estimates of mean entry time, growth rate, and their uncertainty for each bacterial strain (Table 1). We found that the characteristic time required to enter a larval zebrafish, τ, is approximately 2-fold smaller for the evolved isolates compared to the ancestor. Furthermore, the model output did not show significant differences in in vivo growth rates among these strains.

Evolved isolates demonstrate increased motility in swimming assays
In an effort to identify the physiological mechanism underlying increased migration into the host, we investigated the motility behavior of the ancestral and evolved isolates. The role of motility and chemotaxis in host associations is well recognized [5,[45][46][47][48], especially for pathogenic bacteria [49,50]. Furthermore, the importance of motility and chemotaxis in the zebrafish-microbe symbiosis has been previously reported [5,46]. How these traits increase host  colonization has not been characterized extensively, but they could play roles in increasing the rate of microbe-host contact. We conducted low-agar swim plate assays using rich media (tryptic soy broth [TSB]) or media made from fish-conditioned EM (FC-EM), which more closely replicates the conditions experienced by Aer01 when colonizing fish. Motility in the rich media plates was determined via standard protocol, by measuring the diameter of the motility zone. In rich media, population progression away from the site of inoculation is facilitated by growth, motility, and chemotaxis. We observed no significant differences in diameters of the motility zones between the ancestor and evolved isolates, indicating no motility advantage in rich media ( Fig 5A). FC-EM was obtained by filter-sterilizing EM from flasks of hatched 5-6 dpf GF zebrafish larvae. It is essentially a nutrient-poor, dilute salt solution containing only limited fish-derived nutrients (mucus and other secreted or sloughed fish compounds), and therefore, the contribution of growth to progression through this media should be minimal. Cell density in FC-EM swim plates is too low to visually measure a motility zone. Hence, assessment of Aer01's ability to swim in the FC-EM plates was accomplished by sampling an agar plug at a standardized distance away from the site of inoculation ( Fig 5B). The plates were inoculated from rich media liquid cultures; then, a plug was sampled 8 hours and 24 hours post inoculation, homogenized, and plated to enumerate Aer01 cells. In FC-EM swim plates, a significantly higher number of cells were present in the agar plug at the 8 hour postinoculation sampling for each of the evolved isolate populations, compared to the ancestor ( Fig 5C). Replicate FC-EM swim plates sampled at 24 hours show that the ancestor did eventually migrate to the standardized distance given sufficient time ( Fig 5D). Differences in progression of populations through the media in swim plate assays can be attributed to either chemotaxis (i.e., the ability to sense and respond to gradients of chemical attractants and/or repellants) or capacity for motility (e.g., presence of functional flagellar machinery, differences in cellular swim speeds, etc.). To distinguish between these and more precisely characterize the advantage of the evolved isolates, we used microscopy to directly track and measure cellular swimming characteristics of Aer01 ancestral and evolved populations. Ancestral and evolved isolates were incubated in FC-EM for several hours to acclimate and then were imaged. In these experiments, we were specifically testing for differences in general motility, as there are no chemical gradients and therefore no influence of chemotactic behavior. For each strain, 8 movies were recorded (from two independent replicate cultures), which were analyzed with custom tracking software to measure the velocity of thousands of cells each. The swim speed data are plotted in Fig 5E. Each dot represents the mean swim speed across the entire population of cells tracked in an individual movie. The ancestor had a median swim speed of 11.5 μm/second, whereas the evolved isolates from lines 1, 2, and 3 were about 2-fold faster, on average, with median swim speeds of 26.5, 27, and 24.7 μm/second, respectively ( Fig 5E). When we determined the motile fraction of the population from each movie (i.e. the proportion of the cells tracked with >5 μm/second swim speeds), two of the three evolved isolates also had significantly higher motile fractions compared to the ancestor, with the ancestor populations displaying a median motile fraction of 0.13, whereas the medians for evolved isolates from lines 1 and 2 were 0.44 and 0.33, respectively (Fig 5F). Taken together, these data demonstrate that the evolved isolates have a hypermotile phenotype in this nutrient-poor medium, the same condition that Aer01 experiences in a fish flask.

Further-evolved isolates have host genotype-dependent competitive fitness outcomes
The experiments described above were focused specifically on investigating the adaptation mechanisms of the earlier-evolved (passage 4) isolates. Based on the progressive increase in CIs seen in further-evolved isolates (Fig 2), we hypothesized that those strains carried adaptations in addition to the increased immigration phenotype seen in passage 4 isolates. We first verified that further-evolved isolates (passage 18) have the increased-immigration phenotype (S8 Fig; S2 Table). The two passage 18 isolates tested from two of the independent lines both had increased migration rates, of a similar magnitude as the passage 4 isolates, compared to Experimental bacterial adaptation to an animal host the ancestor. This is consistent with the model that they are of the same lineage as the passage 4 isolates, carrying additional adaptive mutations, rather than being independent lineages with more beneficial gain-of-fitness mutations.
We then hypothesized that the further-evolved isolates have adaptations specific to the intrahost environment. Identification of a host genotype in which evolved isolates compete differently compared to the WT host would support our hypothesis that the evolved isolates are adapted to the specific conditions within the WT host.
To test this, we compared CIs when competed in WT hosts (the same host genotype used for the evolution experiment) to a different host genotype, an myd88mutant. MyD88 is a protein adaptor for innate immune receptors, including the Toll-like receptors; hence, mutants are immunodeficient [23,51]. Use of this alternate host genotype is based on the rationale that host immune functions in the vertebrate gut potentially alter, either directly or indirectly, the gut environment [52,53]. Indeed, we have shown previously that myd88 mutant larval zebrafish intestines differ from WT intestines in their paucity of neutrophils, mucus-secreting goblet cells, and proliferating epithelial cells [36,51,54]. We first verified that there were no overt differences in Aeromonas colonization of the myd88fish by showing that ancestral Aer01 competed equally in both WT and myd88hosts (S9 Fig). We then measured the CIs of early-and further-evolved isolates in both host genotypes. As in WT hosts, early-evolved isolates had increased CIs in myd88hosts. The increased migration rate of early-evolved isolates was recapitulated in the myd88host, confirming that the mechanism of this competitive advantage is the same as in WT fish (S10 Fig). For all three lines, the CIs of the passage 18 isolates in WT hosts are significantly higher than the passage 4 isolates in WT hosts (Fig 6). If this advantage was not specific to the host in vivo environment, it would be expected that the CIs of the passage 18 isolates in myd88hosts would also be higher compared to the level observed in the passage 4 isolates in myd88hosts. However, there were no significant differences in CIs in myd88hosts between the passage 4 and passage 18 isolates. Combined, these data support that further-evolved isolates have additional adaptations specific to the intrahost environment, because they confer increased fitness in the WT host but do not provide fitness advantages in a different host genotype.

Discussion
We have developed a new model system in which the roles of both intra-and extrahost environments can be studied simultaneously to investigate the microbial traits that contribute to host Further-evolved isolates have host genotype-dependent competitive fitness outcomes. CIs of evolved isolates and the ancestor, when competed against a differentially tagged ancestral strain in either WT (circles) or myd88 -(triangles) larval zebrafish. Each data point represents a CI from an individual fish. Open symbols represent data points at which the reference strain was undetected on plates; therefore, CIs were calculated by setting the abundance of the reference strain to the limit of detection (5 CFU/gut), thus providing a lower bound for the CI for these points. ���� p < .0001, �� p < .01, � p < .05, ANOVA, log-transformed data, Tukey post hoc analysis. Data combined from two independent experiments for all assays (n � 16/condition). Underlying data are provided in S1 Data. CFU/gut, colony-forming unit per gut; CI, competitive index; ns, not significant; WT, wild type. association. Using this system, we observed that Aer01 quickly and reproducibly increased in zebrafish host association and that the adaptive strategy of the early-evolved isolates was an increased ability to migrate into the host. This was first supported by data showing that competitive fitness is dependent on the mode of inoculation, such that the evolved isolates only outcompeted the ancestral strain when inoculated through the flask water, as opposed to inoculating fish directly by gavage (Fig 3). By investigating the initial colonization dynamics of these isolates at a finer resolution, through sampling mono-associated fish at short time increments, we were then able to observe, quantify, and model this increased immigration (Fig 4; Table 1). This trait contributes to the within-host fitness of the evolved strains by enabling them to constitute a larger portion of the founding population within a host and, consequently, a larger portion of the mature population. Although the predominant impact of early adaptations is on immigration, we cannot rule out that there could be pleiotropic effects of these adaptions on other aspects of host colonization. For example, they could affect gut transit time or expulsion rates. Further investigation exploring the secondary impacts of early adaptions on colonization dynamics could provide additional biological insight into the mechanisms of host colonization.

Immigration impacts multiple facets of colonization dynamics
In addition to playing a role in the initial establishment of the gut community, enhanced immigration could amplify the competitive fitness of the evolved isolates in intrahost populations over time. It has been shown previously that Aer01 goes through major stochastic population collapse events approximately every 24 hours, with decreases in abundance of more than 90% in as little as 30 minutes [32]. These collapse events result in Aer01 expulsion from the fish vent via gut motility and peristaltic activity. Recovery of the gut population from these events can occur from the residual population within the gut, but external Aer01 can also contribute to repopulation [32]. Indeed, we previously established that there is mixing of extrahost with intrahost populations over time ( S3 Fig). Therefore, it is likely that our evolved, hypermigratory Aer01 strains have the potential to out-compete ancestral Aer01 for available open niche space during the repopulation phase that follows each collapse event. This amplification could explain how increased migration evolved so quickly within this system, because it would effectively ratchet up the strength of selection for the immigration phenotype. Moreover, these neutral processes (stochastic collapse events) combined with the deterministic processes (repopulation/regrowth dynamics) may account for the high variability we observe within our competition data. This combination of collapse and regrowth events is not unique to Aer01; these types of dynamics have been reported for other host-associated microbiome constituents -for example, from the human skin [55] and intestinal tract [56]. Whether these cyclic growth dynamics impart a selective pressure for transmission traits and contribute to the evolution of these communities is an open question.
Another way in which the immigration phenotype could impact colonization dynamics within this system is through host-to-host transmission. Our data showing that an evolved isolate had higher abundance than the ancestor in GF recipient fish when introduced from monoassociated donor fish support this hypothesis (S7 Fig). Recent research has demonstrated a role of transmission and dispersal in shaping the gut communities of many animals [23,26,[57][58][59][60]. For example, Burns and colleagues showed that interhost dispersal is a key driver of gut community assembly in zebrafish hosts and can even override the selective forces imposed by the host immune system [23]. This body of literature suggests that selection of traits that improve transmission or dispersal among hosts could be a dominant force in the evolution of host-associated microbes. The work presented here provides support for this hypothesis and has important implications for how we think about the ways in which bacterial symbionts initiate and hone relationships with their hosts. In the future, it will be important to add back complexity to this distilled model system-for example, by increasing microbiota diversity-and determine if the relative importance of immigration is evident. Such studies could reveal trade-offs between colonization of the host and microbe-microbe competitive fitness.
Enhanced transmission is predicted to alter a number of important aspects of host-microbe ecology. For example, it has been proposed that trade-offs exist between transmission and virulence [61] and that mode of transmission influences virulence evolution [62]. However, we do not have any evidence for Aer01 virulence in zebrafish. Previous work has shown that Aer01, in mono-association, does not elicit a host immune response above what is measured in conventional fish and can even attenuate the host response to hyper-immunostimulatory species [38]. Of note, we did not observe any indications of pathogenicity in evolved, hypertransmissible Aer01.
In order to more thoroughly characterize the phenotype of our evolved isolates, we developed a mathematical model to estimate and compare the migration and in vivo growth rates of our ancestral and evolved isolates. Our mathematical model links observed population statistics to essential processes of bacterial colonization, treated stochastically, followed by logistic growth, treated deterministically. We note, however, that additional processes are present and could be elaborated in future models if warranted by relevant data. Bacterial growth, death, and expulsion from the gut [32] are inherently stochastic, and future stochastic models could quantitatively link these processes to measurable population statistics, either through brute force computation or potentially by making use of recently developed analytic tools [63]. Additional processes neglected in this work include the transition time between lag phase and exponential growth upon entering the gut environment (itself stochastic) and interindividual variation in growth rates. Building on the present model to explore how these factors can alter observable features of bacterial abundances will likely lead to further insights.

Bacterial motility is an adaptable trait for host colonization
To further investigate a cellular mechanism to explain the immigration phenotype, we tested the capacity for motility across strains both in classical low-agar swim plate assays and by directly measuring swim velocities of individual cells. Both assays show that the evolved isolates have hypermotile phenotypes compared to the ancestor. It is unclear how this would translate into increased host immigration, but it may act by increasing the likelihood that the evolved isolate will encounter a host. Indeed, once in the host mouth, this phenotype does not promote competitive fitness, since we demonstrated that starting competitions from the host mouth (via mouth gavage) did not confer a competitive advantage (Fig 3). The importance of motility in zebrafish-microbe symbiosis has been previously reported [5,23,46]. For example, Burns and colleagues used ancestral trait reconstruction to infer that genetic pathways associated with motility and chemotaxis were enriched in microbial communities present in the guts of cohoused zebrafish (relative to those in solitary-housed zebrafish), suggesting that motility and chemotaxis may be important to dispersal-adapted gut communities [23]. Stephens and colleagues conducted a transposon-mediated mutational screen to identify genetic determinants for host colonization using two zebrafish gut isolates from different genera, Vibrio and Aeromonas (although a different species than Aer01) [5]. For both species, mutations in chemotaxis and motility genes decreased host-colonization propensity. It was not determined, however, how the impaired motility affected colonization dynamics. Rawls and colleagues previously showed that flagellar motility is important for host colonization of the zebrafish by a human opportunistic pathogen, Pseudomonas aeruginosa, and is required for interaction with the host immune system [46]. The role of motility in other aquatic host-microbe systems has been described as well, such as the Vibrio-squid symbiosis [47]. Additionally, chemotactic motility was shown to be important for maximum virulence of a pathogenic Vibrio anguillarum strain in rainbow trout, and this was shown to be dependent on a natural route of transmission, as motility mutants directly injected into fish did not have decreased virulence [64]. In each of these examples, the stage of colonization impacted by loss of motility functions was not fully resolved. To our knowledge, our study is the first to clearly demonstrate an association between bacterial motility and immigration into the zebrafish host.

Bacterial adaptation to host-microbe metacommunities
This new model system was designed with the specific goal of studying bacterial host association in the context of an entire ecological framework, incorporating multiple hosts and extraand intrahost microbial populations. Although this aquatic model system may not be an intuitive surrogate for understanding terrestrial host-microbe interactions, it is relevant to such systems in many ways. For example, like zebrafish, all animals begin life with very few to no resident microbes. Therefore, there is an initial establishment phase during which immigration of microbes into the naïve host occurs. Likewise, there is growing awareness that individual host microbiomes can be considered as part of a larger metacommunity, connected to the microbiomes of different hosts via transmission and dispersal [22,65]. Compelling evidence for this is the growing body of research showing that cohabiting animals such as mice [66,67] and humans (e.g., [58,68]) have more similar gut microbiomes than those not cohoused. Furthermore, although the mechanism of immigration may be different in different systems (e.g., swimming through the EM in our system versus transmitting between human hosts by attachment to a fomite), the trait (increased immigration) is the same. These points illustrate that the intra-and extrahost factors for host association that we consider in our model are likely broadly applicable to other systems.
We achieved the foremost goal of this study, which was to identify bacterial traits that confer host association via phenotypic characterization of the evolved isolates. Determination of the genetic basis of these adaptations would provide an even deeper understanding of the underlying biological mechanisms at play. A concession of using a mutator strain is that the increased number of mutations could make identification of adaptive ones difficult. In our evolved genomes, initial sequencing and preliminary analysis of passage 4 evolved isolates resulted in an average of 40 mutations per genome, with no obvious adaptive mutation candidates. An important next step is to repeat this study with non-mutator Aer01, in which there will be much fewer mutations, and therefore the adaptive mutations would presumably be easier to identify. This would allow more in-depth investigation of the genetic and physiological mechanisms of host association and adaptation.
Using a highly tractable model system, we were able to demonstrate that an enhanced extrahost process, immigration from the external environment, can increase host association by a bacterium. This finding challenges conventional assumptions about the primary importance of the intrahost environment in shaping host-microbe interactions. Furthermore, it supports a newly emerging paradigm that recognizes that ecological processes such as transmission may play substantive roles in microbiome assembly and dynamics.

Ethics statement
All experiments with zebrafish were done in accordance with protocols approved by the University of Oregon Institutional Animal Care and Use Committee (IACUC); the animal protocol number for this work is 15-98 [69].

Gnotobiotic zebrafish husbandry
All experiments involving zebrafish were conducted following standard protocols and procedures approved by the University of Oregon Institutional Care and Use Committee. For the evolution passaging and bacterial competitions, WT (AB × Tu strain) fish were used. Competition experiments were also conducted using myd88 mutant zebrafish, which were previously generated via CRISPR-Cas9 system and verified to have the expected phenotype of an myd88 KO mutant [23]. Fish were maintained as previously described [69]. Fish were not fed in any of the experiments described here. GF derivation of fish embryos followed protocols previously described [31]. Generally, fish were inoculated with bacterial cultures at 4 dpf. At 7 dpf, fish were euthanized with tricaine (Western Chemical), following approved procedures, and mounted in sterile 4% methylcellulose solution (Fisher), and the intestines were removed by dissection (described in [70]) and used for enumeration of colonizing bacteria or as inoculum for GF fish.

Bacterial strains
The bacterial strain used in this study is the zebrafish isolate Aer01 (A. veronii, PRJNA205571) previously described [33]. Fluorescently tagged (dTomato or superfolder GFP) variants of this strain were generated using a Tn7 transposon-based system as previously described [32,34]. This method results in the integration of a cassette containing the dTomato/GFP gene and a gentamycin resistance gene in the chromosome at a specific target location. Subsequently, mutS deletion mutants (unmarked, clean deletions) were made in the fluorescently tagged genetic background strains following an allelic exchange system [34]. Strains were always grown at 30˚C, with shaking for liquid cultures.

Experiment to estimate bottleneck size of initial colonization
Overnight rich media cultures of untagged and dTomato-tagged Aer01(non-mutator) were pelleted, washed with sterile EM, and mixed (either 1:100 or 1:300, tagged:untagged). These mixed cultures were used as inoculum for flasks of (10-15) 4 dpf larval fish to a starting density of approximately 10 6 CFU/ml. After 24-hour colonization, at 5 dpf, the fish were dissected and the guts transferred into 1.6 ml tubes containing 500 μl sterile EM and approximately 100 μl of 0.5 mm zirconium oxide beads (Next Advance, Averill Park, NY, US) and were then homogenized using a bullet blender tissue homogenizer (Next Advance, Averill Park, NY, US) for 30 seconds at power 4. One hundred microliters of the homogenate was spread-plated on TSA plates, incubated overnight, and screened by fluorescence microscopy for the presence of dTomato-tagged Aer01 colonies. The experiment was repeated twice, on independent days, with two replicate flasks for each ratio on each day. We used a binomial sampling model to estimate a bottleneck size for the first 24 hours of colonization, using the cumulative distribution function (CDF) to determine the number of cells sampled that would result in the probability of colonization observed. The CDF gives the probability of not finding a "success"-in this case, a tagged cell-given a number of trials (the number of cells sampled in 24 hours), with a given success rate (1:100 or 1:300). One minus this probability is our observed frequency of colonized fish. We can then solve for the number of trials that yield the observed colonization rate. This analysis yielded an average bottleneck estimate of 193.7 cells sampled per day, with a standard error of 140.0.

Succession assay to determine intra-and extrahost population mixing
At 4 dpf, GF larval zebrafish were mono-associated with untagged Aer01 (non-mutator) by inoculating to approximately 10 6 CFU/ml. At 5 dpf, dTomato-tagged Aer01 was added to the flask to the same approximately 10 6 CFU/ml density. At 6 dpf (24 hours) and 7 dpf (48 hours), fish were dissected and the guts homogenized and plated as described above. The ratios of untagged (founder) to tagged (invader) were determined by fluorescence microscopy of the colonies.

Fish-to-fish transmission assay
At 4 dpf, flasks of GF larval zebrafish were mono-associated with either the ancestral strain or the passage 4 (line 3) evolved isolate. At 6 dpf, two mono-associated fish ("donors") were washed with sterile EM and then transferred into flasks containing 13 GF larval zebrafish that were 6 dpf ("recipients"). Twelve hours later, all fish in the flasks were dissected and the guts homogenized and plated to determine the colonization level of each strain in the recipient fish.

Evolution experiment
Evolution passaging was initiated using an equal mixture of ΔmutS, attTn7::dTomato and ΔmutS, attTn7::sfGFP strains as a way of tracking gain-of-fitness lineages throughout passaging by monitoring changes in ratios of tagged populations. In all three replicate lines, the ΔmutS, attTn7::sfGFP populations were undetectable by passage 3, suggesting either a fitness defect in this genome or that the emergence of gain-of-fitness mutants arose in the ΔmutS, attTn7::dTomato lineage in all three lines by chance. The first passage was inoculated by pelleting 1 ml TSB overnight cultures of ΔmutS, attTn7::dTomato and ΔmutS, attTn7::sfGFP strains, resuspending them in 1 ml sterile EM, mixing them 1:1, and then adding them to replicate flasks of 4 dpf GF WT fish (10-15 larval fish, 15 ml EM) to a final volume of 10 6 CFU/ml. Inoculated fish were incubated according to IACUC protocol. At 7 dpf, fish were euthanized with tricane, and the intestines were removed by dissection (described in [70]). All guts from a flask of fish were combined into a single 1.6 ml tube containing 500 μl sterile EM and approximately 100 μl of 0.5 mm zirconium oxide beads (Next Advance, Averill Park, NY, US) and were then homogenized using a bullet blender tissue homogenizer (Next Advance, Averill Park, NY, US) for 30 seconds at power 4. Aeromonas populations were monitored by dilution-plating a small aliquot (20 μl) of the combined gut sample and an aliquot of the flask EM on TSA plates. Half (approximately 250 μl) of the homogenate was mixed with 250 μl of sterile 50% glycerol and then stored at −80˚C. The remaining homogenate was stored at 4˚C for 4 days and then used as inocula for the subsequent flasks of 4 dpf GF fish. For subsequent inoculations, all of the 4 C sample (about 200 μl) was added to the next flask of fish (resulting in approximately 10 4 CFU/ml at the beginning of the passage). This process was continued for 22 passages total, for all three lines.

Colony purification of evolved isolates
From freezer stocks of whole populations from selected evolution passages (namely, 4, 8, 14 [13, line 1], 18, and 22), cells were streaked out onto TSA plates for isolation and then incubated at 30˚C for 1 day. Isolated colonies were then picked from the plates into 5 ml TSB cultures and allowed to grow shaking at 30˚C for about 6 hours; then, 25% glycerol freezer stocks were made and stored at −80˚C.

Bacterial competitions in vivo
For in vivo bacterial competitions, to assess fitness, selected strains were grown overnight in TSB from freezer stocks. One milliliter of the overnight cultures was pelleted (8,700 rcf, 2 minutes) and then resuspended in 1 ml sterile EM. Competing strains were mixed at about 1:10 (competitor:reference) and then added to flasks of 4 dpf GF fish (either WT or myd88 mutant zebrafish) at about 10 6 CFU/ml. For all competitions, ancestor (ΔmutS, attTn7::dTomato) or evolved isolates were competed against the WT ancestral strain, Aer01 attTn7::sfGFP. At 7 dpf, fish guts were dissected as described above, and each gut was transferred into a 1.6 ml tube containing 500 μl sterile EM and about 100 μl of bullet beads and then bullet blended as described above. Blended samples were then diluted appropriately in sterile EM, spread-plated on TSA plates, and incubated at 30˚C for 1-2 days, and the colonies were counted. Strains were differentiated by fluorescence microscopy. The limit of detection is 5 CFU/gut.

Bacterial competitions in vitro
For in vitro competitions, selected strains were grown overnight in TSB from freezer stocks. Competing strains were mixed at about 1:4 (competitor:reference); ancestor (ΔmutS, attTn7:: dTomato) or evolved isolates were competed against the WT ancestral strain, Aer01 attTn7:: sfGFP. Mixed cultures were then used to inoculate 5 ml TSB cultures (1:100 back dilution). Cultures were incubated at 30˚C for 24 hours and then diluted and spread-plated on TSA for enumeration via fluorescence microscopy. Data are combined from three independent replicate experiments.

Gavage experiments
Generally, the gavage protocol was followed as previously described [44], with the following modifications. Briefly, gavage needles were made by pulling 3.5-inch (Drummond #3-000-203 GIX) capillaries, microforging (DMF1000, World Precision Instruments) them to an internal diameter of approximately 30 μm, and polishing the ends. Inoculum for gavaging was prepared by pelleting 1 ml of TSB overnight cultures (8,700 rcf, 2 minutes), resuspending them in 1 ml sterile EM, and mixing competing strains at about 1:2 (competitor:reference); ancestor (ΔmutS, attTn7::dTomato) or evolved isolates were competed against the differentially tagged non-mutator ancestral strain, Aer01 attTn7::sfGFP. Culture mixes were then diluted 1:30 in sterile EM for gut gavage (no dilution for mouth gavage). Prepared inocula were incubated at room temperature until gavaging and flask inoculation (approximately 30-60 minutes), allowing time for acclimation to the EM. GF 5 or 6 dpf WT fish were euthanized (120 μg/ml tricaine) and transferred to 3% methylcellulose-coated gavage mold (4% agar). Culture mixes were loaded into gavage needles, and 4.6 nl of culture mix was gavaged directly into the lumen of the gut or the mouth of individual fish using a Nanoject II (Drummond Scientific Company) on the "slow" setting. After gavage, fish were rinsed in sterile EM and then transferred into sterile EM. Immediately after gavaging, flasks of GF fish were inoculated at 10 6 CFU/ml with the same inocula used for gavaging. Approximately 5 hours post gavage, fish were euthanized and dissected, and the guts were plated as described above.

Migration rate experiments
Selected isolates were grown overnight in TSB, at 30˚C, with shaking. Overnight cultures were diluted to 10 −4 in EM, subcultured 1:100 in FC-EM, and grown for an additional 24 hours at 30˚C, with shaking. FC-EM was obtained by filter-sterilizing (0.2 μm) flask EM from flasks of hatched 5-6 dpf GF zebrafish larvae. Flasks of 5 or 6 dpf GF zebrafish were combined and then inoculated with FC-EM cultures to yield approximately 10 4 CFU/ml. The fish were then split into replicate flasks containing 10 fish and 10 ml inoculated flask EM. An EM sample was taken and plated immediately to quantify CFU/ml of inoculating strain. Subsequently, about every 45 minutes (for up to about 350 minutes), all of the fish in a replicate flask were dissected and the guts transferred into 1.6 ml tubes containing 200 μl sterile EM and approximately 100 μl bullet beads. The guts were homogenized as described above and all 200 μl spreadplated on TSA plates. Flask EM samples were also plated to enumerate bacterial CFU/ml. After 24-48 hours incubation, colonies on plates were counted.

Theoretical estimation of growth rate and colonization time
The distribution of bacterial abundances across fish allowed us to estimate the migration and growth rates of each strain. We considered a model in which migration into a fish by an individual bacterium is a stochastic event, with some probability per unit time. The statistics of the colonizers at any time are therefore governed by a Poisson distribution characterized by a mean time to entry, τ. Following entry into the host, the descendants of a colonizer obey logistic growth (initially roughly exponential, constrained by a finite carrying capacity, K), characterized by some growth rate r. For the migration rate experiments shown here, the population size is orders of magnitude less than the carrying capacity, leaving r as the only relevant growth parameter. For each fish gut examined by dissection and plating, we measured the bacterial abundance, N (CFUs). For N from any individual fish, one cannot separately determine τ and r; the same final population can be reached, for example, by rapid colonization followed by slow growth (low τ, high r) or by slow colonization followed by fast growth (high τ, low r). The distribution of N across multiple fish, however, constrains both τ and r. Roughly, high τ and high r gives distributions with large variances and many uncolonized fish (N = 0), being dominated by the stochasticity of rare migration events, whereas low τ and low r gives distributions with smaller variance (S11 Fig).
We fit the experimental data (i.e., each set of CFU measurements at a particular measurement time t) to the above model of stochastic colonization and growth, determining the bestfit parameters τ and r. This is done by varying τ and r across a range of possible values, simulating for each τ and r the model outcome for many (10,000) replicates, and calculating the likelihood that the simulated data match the set of measurements. This gives both the maximum-likelihood parameter values and a measure of the likelihood distribution across the full parameter range. For most datasets (41 of 50), the likelihood p(τ, r | CFU) shows a sharp peak (S12 Fig). For the remainder, it does not, indicating either limitations of this model (e.g., the neglect of exit from the gut) or experimental error or contamination; these datasets were manually discarded. The datasets from each measurement time give independent estimates of τ and r; the average and standard error give the final estimates of mean entry time, growth rate, and their uncertainty for each bacterial strain (Table 1).

Swim plate assays to assess motility
Rich media swim plates were made using TSB containing 0.2% agar (VWR Life Science AMRESCO Agarose). FC-EM swim plates were made by heat-dissolving 100 mg agar in 20 ml sterile EM and then cooling slightly and adding 30 ml sterile FC-EM (0.2% agar). Media (20 ml) were added to 100 mm petri dishes, cooled for 2 hours, and then inoculated. For inoculation, TSB overnight cultures of selected strains were pelleted (8,700 rcf, 2 minutes) and resuspended in one-tenth of the culture volume, and 3 μl was used to inoculate plates. A template for inoculation was made to facilitate insulation and sampling (S6 Fig). Plates were incubated at 30˚C and then sampled and measured at 8 hours (TSB and FC-EM), 24 hours (FC-EM), and 48 hours (FC-EM). To assess motility in TSB swim plates, the diameter of the swim zone was measured in centimeters. To determine motility in FC-EM swim plates, in which cell density was too low to visualize, we sampled agar plugs from the plates at a defined distance from the point of inoculation (Fig 5) using x-tracta Gel Extractor tool (Promega, Madison, WI, US). The agar plug was transferred into a 1.6 ml tube containing 500 ml sterile EM and approximately 100 μl 0.5 mm zirconium oxide beads (Next Advance, Averill Park, NY, US), which were then homogenized using a bullet blender tissue homogenizer (Next Advance, Averill Park, NY, US) for 30 seconds at power 4. Homogenized samples were diluted and spreadplated onto TSA plates and incubated at 30˚C overnight, and colonies were counted to enumerate the number of bacterial cells that migrated to the site of sampling.

Live imaging and measurement of cellular swim speeds
Overnight cultures (TSB) of the ancestral and evolved isolates were washed with sterile EM; then, 50 μl was added to 1 ml sterile FC-EM, and the cultures were incubated at 30˚C for approximately 4 hours to recover and acclimate. The culture was placed in a glass cuvette and imaged in a custom-built light sheet fluorescence microscope, as previously described [73]. The light sheet optically sections bulk samples; thus, the motility of the imaged bacteria is unconstrained by surfaces. Movies in a single optical plane were captured at 30 frames per second for a duration of 20 seconds, with excitation light provided by a 561 nm solid-state laser (Coherent Sapphire 20 mW; all strains express dTomato fluorescent protein). Images were analyzed using standard particle-tracking techniques, with fine localization of bacterial positions determined by a radial symmetry-based algorithm [74], the source code for which is available here: http://pages.uoregon.edu/raghu/particle_tracking.html. For each strain, four movies were recorded from each of two biologically independent replicate cultures (eight total; seven for the ancestor). An average of 6,000 (minimum of 2,300) bacteria were tracked in each movie. In Fig 5, each data point represents average bacterial swimming velocity (μm/sec) across motile cells (velocity > 5 μm/sec) in an individual movie (Fig 5E) or the fraction of the total cells that were considered motile (velocities > 5 μm/sec) in an individual movie (Fig 5F). Mono-associated "donor" fish were incubated with 13 GF "recipient" fish for 12 hours, the guts were dissected, and the abundance of Aer01 in the fish was determined. p-Value = 0.003, determined using two-tailed Student's t test, excluding the donor fish points. Underlying data are provided in S1 Data. GF, germ-free. (TIF)

S8 Fig. Further-evolved (passage 18) isolates demonstrate increased host immigration. (A)
Groups of mono-associated fish were dissected and plated about every 45 minutes, and the fraction of colonized hosts was determined. A higher proportion of fish are colonized at earlier time points for the evolved isolates compared to the ancestor. Mean (± SEM) is plotted for ancestor; ancestor data are the same as those plotted in Fig 4A. (B) CFU/gut (mean ± SEM) from the same samples presented in A show higher gut abundance for the evolved isolates. Ancestor data are the same as those plotted in Fig 4B. (C) EM abundance (CFU/ml) shows strains do not have differences in survival in flask EM. Mean (± SEM) is plotted for ancestor. Data combined from three (ancestor) or one (evolved) independent experiment. Underlying data for A-C are provided in S1 Data. CFU/gut, colony-forming unit per gut; EM, embryo medium.