Microbiome Interaction Networks and Community Structure From Laboratory-Reared and Field-Collected Aedes aegypti, Aedes albopictus, and Culex quinquefasciatus Mosquito Vectors

1The microbiome of mosquitoes can be highly variable, both within and between species, and is often dominated by relatively few genera. Understanding the factors that influence this variation is important as microbes drastically alter host biology. For mosquitoes, bacteria can affect a diverse number of traits including immunity, reproduction, survival, and vector competence. These phenotypes have ramifications for the vectorial capacity of pathogens, and as such, microbial-based vector control strategies are under investigation to reduce the burden of arthropod-borne diseases. While our understanding of the contributing factors that affect the composition and abundance of the microbiome is expanding, there are still many unanswered questions regarding microbiome assembly and maintenance within mosquito hosts.

2Exposure to environmental microbes is undoubtedly a major influence on the mosquito microbiome. These effects are particularly pronounced at the aquatic stage as larvae and pupae are immersed in water and can acquire bacteria by filter feeding. Indeed, several studies have shown the larval stages possess a similar microbiome as their larval water environment, and exposure to bacteria at these immature stages has implications for adult traits. Furthermore, the larval habitat can influence the composition of the adult microbiome. Bacteria can be transstadially transmitted to the adult, and newly emerged adults are known to imbibe their larval water, which likely seeds the gut with microbiota.

3Host and bacterial genetics also contribute to microbiome composition and microbial abundance. Mosquitoes can maintain microbiome homeostasis by a variety of different mechanisms. Host pathways and processes known to influence microbiota in mosquitoes include immunity, amino acid metabolism, reactive oxygen species, and calcium transport. Additionally, serial passaging of gut symbionts in mosquitoes selected for isolates that persist in the gut for longer periods of time, indicating that bacterial genetics is also important in shaping the microbiome.

4Adult mosquito feeding behavior also has important implications for microbiome community structure. It is likely that bacteria can be acquired from the nectar of plants, and taking a blood meal alters the microbiome considerably. At 24 h post-blood meal, the bacterial load in the gut drastically increases while species diversity decreases. Culture based assays show that bacterial loads revert to pre-blood fed levels 2–3 days after the blood meal, although other studies have seen high bacterial loads persist for some time and species richness not reverting to the original composition seen prior to the blood meal. Most of these studies either used laboratory-reared mosquitoes to examine culturable bacterial load, or relocated field mosquitoes to the laboratory for experimentation, and as such, there are few studies examining the effect of blood feeding on the microbiome community structure in field populations.

5Recently, it has become evident that a further force affects microbiome composition in mosquitoes – interactions between the microbes themselves. These interactions were first highlighted with the discovery that Wolbachia and Asaia are antagonistic to one another, thereby affecting the vertical transmission of Wolbachia in Anopheles mosquitoes. Further comparisons exploiting 16S rRNA amplicon high throughput sequencing have identified interactions between Wolbachia and other microbes in both Drosophila and mosquitoes. In addition to the specific interactions between Wolbachia and other bacterial taxa, pairwise negative and positive microbial interactions within bacteria or fungi, as well as cross-kingdom interactions (bacterial-fungal) have been reported in the La Crosse virus vectors, Aedes triseriatus, and Ae. japonicus. Taken together, these studies suggest that microbial interactions are important in dictating the composition and abundance of host-associated microbiota. Yet it is unclear how ubiquitous and complex these interactions are within mosquitoes.

616S rRNA amplicon sequencing datasets have been analyzed to create microbial co-occurrence networks for several species and environments, but these networks are lacking for mosquitoes and insects in general. These methods use presence/absence metrics, relative abundance, or both, to examine pairwise interactions to develop a network. Usually, interacting pairs of taxa are designated as having co-occurring or co-exclusionary relationships. Each method used for the identification of co-occurrence/co-exclusion networks has strengths and weaknesses in identifying particular patterns. CoNet, uses an ensemble approach that combines results from a collection of algorithms using presence/absence and relative abundance data to identify statistically significant interactions. Interaction networks provide another methodology to examine the community structure of the microbiome of mosquitoes. Comparing microbiome networks generated from mosquitoes exposed to different conditions may provide insights into factors influencing microbiome structure in mosquitoes and identify pairwise interactions not affected by environmental conditions.

7To expand our understanding of the forces that shape the bacterial microbiome of mosquitoes, we examined the microbial composition and community structure from three major mosquito arboviral vectors, Ae. aegypti, Ae. albopictus, and Cx. quinquefasciatus, collected from the field or reared under uniform insectary conditions. For the field collections, we utilized two trapping methods that primarily attract mosquitoes in different physiological states: host- or oviposition-seeking. Our sampling regime allowed us to examine how factors such as host species, environment, and physiological state in the field influenced the composition of the mosquito microbiome. We also compared the microbiome of mosquitoes that were naturally infected (Ae. albopictus and Cx. quinquefasciatus) and uninfected (Ae. aegypti) with Wolbachia. Furthermore, we developed microbial interaction networks to explore the complexity and nature of microbial interactions in mosquitoes. To demonstrate that microbial interactions influence host colonization, we infected Ae. aegypti larvae either possessing or lacking their native microbiota with a range of bacterial symbionts. Our results highlight the complexities of microbial networks in field collected Ae. aegypti mosquitoes and indicate the native microbiome induces colonization resistance to certain gut microbes.

8Female Ae. albopictus and Cx. quinquefasciatus were collected from an abandoned tire repository in south-eastern Harris County, Houston, TX, United States, while female Ae. aegypti were collected from a separate site in Houston (Supplementary Figure S1). Further details describing the tire repository location were previously reported. All mosquitoes were trapped over a 24 h period with either the Biogents Sentinel (BG) or Harris County gravid (G) traps, which selectively collect host-seeking or gravid female mosquitoes, respectively. Mosquito species were identified using morphological characteristics, surface sterilized (5 min in 70% ethanol followed by three washes in 1× PBS each for 5 min), and stored in ethanol at -20°C while awaiting DNA extraction. 5–7 days old adult sugar fed laboratory-colonized mosquitoes (Ae. aegypti; Galveston strain, Ae. albopictus; Galveston strain, and Cx. quinquefasciatus; Houston strain) were reared under conventional conditions and then processed in the same manner as field samples. All laboratory-reared mosquitoes were reared in the insectary at the University of Texas Medical Branch.

9High-throughput sequencing of the bacterial 16S ribosomal RNA gene was performed using gDNA isolated from each sample. Sequencing libraries for each isolate were generated using universal 16S rRNA V3-V4 region primers in accordance with Illumina 16S rRNA metagenomic sequencing library protocols. The samples were barcoded for multiplexing using Nextera XT Index Kit v2. Sequencing was performed on an Illumina MiSeq instrument using a MiSeq Reagent Kit v2 (500-cycles). The NCBI Bioproject accession number for the raw sequencing data reported here is PRJNA422599.

10To identify the presence of known bacteria, sequences were analyzed using the CLC Genomics Workbench 8.0.1 Microbial Genomics Module1. Reads containing nucleotides below the quality threshold of 0.05 (using the modified Richard Mott algorithm) and those with two or more unknown nucleotides or sequencing adapters were trimmed out. All reads were trimmed to 264 bases for subsequent operational taxonomic unit (OTU) classification. Reference based OTU picking was performed using the SILVA SSU v119 97% database. Sequences present in more than one copy but not clustered to the database were placed into de novo OTUs (97% similarity) and aligned against the reference database with 80% similarity threshold to assign the “closest” taxonomical name where possible. Chimeras were removed from the dataset if the absolute crossover cost was three using a k-mer size of six. Alpha diversity was measured using Shannon entropy (OTU level), rarefaction sampling without replacement, and with 100,000 replicates at each point. Beta diversity was calculated using the Bray-Curtis diversity measure (OTU level). PERmutational Multivariate ANalysis Of VAriance (PERMANOVA) analysis was used to measure effect size and significance on beta diversity for grouping variables. The significance is obtained by a permutation test. For each assessment, a permutation of 99,999 was chosen. Differentially abundant bacteria (genus level, >0.1%) were identified using analysis of composition of microbiomes (ANCOM) with a significance level of P < 0.05, while values quantifying fold change were obtained using the log2 fold change formula.

11Network analysis was conducted independently on each group within a species. OTUs with read counts below 0.1% of total number of reads in all samples were excluded from analysis. The remaining OTUs were combined based on lowest common taxonomy assignments down to genus level, and relative abundance tables were generated by normalizing read counts against total number of reads in the original data. The resulting number of unique entries identified in samples was 33. Interactions (such as co-occurrence and co-exclusion) among these were identified using CoNet app in Cytoscape using the following ensemble of methods: Pearson correlation, Spearman correlation, mutual information, Bray-Curtis dissimilarity, and Kullback-Leibler divergence. Statistical significance of each pair was tested using the row-shuffle randomization option and interactions that scored at the top and bottom 5% of 1000 bootstraps were reported. Resulting statistically significant interactions were categorized by the software into three groups: co-presence, co-exclusion, and unknown. Unknown interactions represent statistically significant patterns that cannot be clearly categorized as co-exclusion or co-occurrence. Since we could not ascribe an interaction pattern, the unknown interactions were excluded from the network. Resulting interaction networks were visualized using Cytoscape software.

12Total bacterial load within each mosquito species or group was assessed by qPCR using gDNA as a template. qPCR was conducted using universal bacterial primers (16rRNA-RTF: TCCTACGGGAGGCAGCAGT and 16rRNA-RTR: GGACTACCAGGGTATCTAATCCTGTT) that amplified the bacterial 16S rRNA gene or primers that specifically amplified Wolbachia (W-Spec-16S-F: CATACCTATTCGAAGGGATA and W-Spec-16s-R: AGCTTCGAGTGAAACCAATTC). Primers specific to S7 gene (aeg-S7-F: ACCGCCG TCTACGATGCCA, aeg-S7-R: ATGGTGGTCTGCTGGTTCTT, Cq-S7-F: CTGGAGATGAACTGGACCT, and Cq-S7-R: CTT GTACACCGACGTGAAGG) or actin gene (alb-act-F: CCCACACAGTCCCCATCTAC, alb-act-F: CGAGTAGCCACGTTCAGTCA) were used to quantify host genomic copies. The PCR was carried out in a 10 μl reaction containing 1 μM of each primer, 1× SYBR Green (Applied Biosystems) and 2 μl of genomic DNA template. Cycling conditions involved an initial denaturation at 95°C for 10 min, 40 cycles of 15 s at 95°C, 1 min at 60°C. Fluorescence readings were taken at 60°C after each cycle before deriving a melting curve (60–95°C) to confirm the identity of the PCR product. The PCR was carried out on the ABI StepOnePlus Real-Time PCR System.

13Relative abundance was calculated by comparing the load of all bacteria or Wolbachia to a single copy mosquito gene. PCRs amplifying the 16S rRNA and wsp genes of Wolbachia were used to screen for Wolbachia in field caught Ae. aegypti, and nematode specific primers were used to screen for nematode infections.

14Aedes aegypti gnotobiotic larvae were generated as previously described. All manipulations of mosquitoes were undertaken in a sterile environment. To synchronize hatching, sterile eggs were transferred to a conical flask and placed under a vacuum for 45 min. To verify sterility, larval water was plated on non-selective LB agar plates. L1 larvae grown without bacteria have slow growth rates and do not reach pupation. Forty-five L1 larvae were transferred to a T75 tissue culture flask and inoculated with 1 × 107 CFU/ml of transgenic symbionts possessing the plasmid the pRAM18dRGA-mCherry that was derived from pRAM18dRGA. Bacterial cultures were quantified with a spectrophotometer (DeNovix DS-11, DeNovix) and validated by plating to determine CFU. For conventional rearing of mosquitoes, eggs (non-sterilized) were vacuum hatched and grown under non-aseptic conditions in a T75 tissue culture flask supplemented with transgenic symbionts at the concentration of 1 × 107 CFU/ml. To feed mosquitoes, ground fish food pellets were sterilized by autoclaving, and mixed with sterile water. The equivalent of 6 mg of fish food was fed to both gnotobiotic and conventionally reared mosquitoes every second day. To quantify the bacterial load of symbionts, surface sterilized L4 larvae were homogenized and plated on a selective media (50 μg/ml Kanamycin) on which only transgenic symbionts grew (pRAM18dRGA-mCherry induces resistance to kanamycin). After incubation at 30 or 37°C (depending on symbiont) for 2–3 days, colonies (expressing mCherry fluorescent protein) were counted. All colonies observed on the kanamycin plate were confirmed to have mCherry fluorescence.

15We sequenced amplicons of the V3/V4 region of the 16S rRNA gene from whole individual adult female mosquitoes either collected from the field or reared in the laboratory. In total, we sequenced 130 adult mosquitoes obtaining 10,668,291 reads (sample size per group and species is reported in Supplementary Table S1). After quality filtering, 7,051,256 reads were assigned to OTUs at 97% identity threshold and on average, there were 54,240 reads per mosquito sample. Rarefaction curve analysis indicated that our sequencing depth was sufficient to observe all OTUs in mosquito samples (Supplementary Figure S2). We identified a total of 4,419 bacterial OTUs in the three mosquito species, but only 58 were present at an infection frequency of over 0.1% within the dataset (Supplementary Table S2). When abundant microbes were classified at higher taxonomic levels, our analysis found 22 families, with Enterobacteriaceae being the most common when disregarding Wolbachia. Bacteria found in mosquitoes were classified into five phyla with bacteria in the phylum Proteobacteria most prevalent in the microbiome, which is consistent with previous studies, although other reports indicated Bacteroidetes or Acinetobacter phyla can be a major component of the microbiome.

16When examining the sequencing data at the genus level, the microbiomes of Cx. quinquefasciatus and Ae. albopictus were dominated by the endosymbiont Wolbachia with 87 and 81% of total reads, respectively (Supplementary Figure S3 and Table S2). While Wolbachia accounted for many of the reads, rarefaction analysis indicated our sampling depth was sufficient to identify rare OTUs. Other highly abundant genera in field-collected mosquitoes included Halomonas, Shewanella, and Asaia in Cx. quinquefasciatus, Halomonas, Pseudomonas, and Zymobacter in Ae. albopictus, and Pseudomonas, Zymobacter, Tatumella, and Enterobacter in Ae. aegypti (Supplementary Figure S3 and Table S2). Similar to a recent finding, we found a small number of Wolbachia reads in a few Ae. aegypti individuals collected from the field in G traps. It is thought that Ae. aegypti are naturally uninfected by Wolbachia, although suggested that some populations may be infected. However, samples that contained Wolbachia reads from our sequencing data could not be independently validated by PCR using several Wolbachia genes (wsp and MLST genes) commonly used to screen for the bacterium, nor were they PCR positive for filarial nematode infection (that carry Wolbachia) when amplifying with primers that detect nematode DNA.

17Given the above finding and since microbiome sequencing can be susceptible to contamination, we examined our data for other possible contamination signatures. While we could not find conclusive evidence of contamination, in this process we observed a possible batch effect for a specific bacterium in field-collected Ae. aegypti. These samples were extracted in two batches (see Supplementary Table S1) and the field-collected samples from the latter extraction were found to have higher loads of Burkholderia compared to those extracted in the first batch. However, all laboratory-reared mosquitoes extracted during this second extraction did not contain high levels of Burkholderia, indicating the prevalence of this microbe in the field-collected samples was not likely due to laboratory contamination. Despite this occurrence, the dataset is suitable for constructing networks as this analysis examines pairwise interactions, and as such, other interacting pairs will not be influenced by Burkholderia. Future studies should consider the use of reagent-only controls and spike controls to help determine if cross-contamination of samples occurs in the sequencing process, although it is questionable if these controls would have been of assistance in this case given that laboratory-reared mosquitoes did not have elevated Burkholderia reads.

18The Shannon diversity index was used to estimate the species richness in mosquitoes. There were significant differences in diversity between the microbial communities of Aedes mosquitoes compared to Cx. quinquefasciatus when compared at the OTU level (Figure 1; Kruskal-Wallis P < 0.0001). When comparing groups within Ae. aegypti, field collected mosquitoes from either the BG or G trap had a significantly greater Shannon diversity index compared to the laboratory mosquitoes. No differences were seen between groups in Ae. albopictus, while Cx. quinquefasciatus caught in the G trap had a significantly lower Shannon diversity index compared to laboratory-reared mosquitoes. Across all species, we found no significant differences in species richness between mosquitoes caught in the BG and G traps.

19Since both Ae. albopictus and Cx. quinquefasciatus were heavily infected with Wolbachia, we examined alpha diversity (OTU level) in these mosquitoes when this endosymbiont was computationally excluded from the dataset (Supplementary Figure S4). In all cases, we observed an increase in Shannon diversity when Wolbachia was excluded. This was significant when analyzed by species (Ae. albopictus P < 0.0001; Cx. quinquefasciatus P < 0.0001), and for all groups with the exception of Ae. albopictus caught in the G traps (Ae. albopictus: BG P < 0.005, Lab P < 0.02; Cx. quinquefasciatus: BG P < 0.04, G P < 0.001, Laboratory P < 0.01). No significant differences were seen when comparing the BG to G groups after removal of Wolbachia (Ae. albopictus: P = 0.23; Cx. quinquefasciatus: P = 0.57).

20We examined how the environmental, physiological state in the field, and host species affected the bacterial community structure using beta diversity analysis by comparing the microbiomes of the three mosquito species or groups. Within each group, distinct microbiome clustering patterns were observed between the three mosquito species and all pairwise comparisons were significantly different (Figure 2; PERMANOVA with Bray-Curtis distance comparison P < 0.05). In general, the microbiome of Ae. aegypti was more divergent compared to the microbiomes of Cx. quinquefasciatus and Ae. albopictus regardless of origin (Figure 2A). When Wolbachia was computationally excluded from the analysis (Figure 2A), the microbiomes of Cx. quinquefasciatus and Ae. albopictus became more divergent in the laboratory, but were not significantly different when considering mosquitoes caught in the BG trap (P = 0.24). Interestingly, there was tight clustering of samples from Cx. quinquefasciatus and Ae. albopictus caught in the G trap (Figures 2A), although these two groups were significantly different (P = 0.00001). When Wolbachia was removed, this clustering became more divergent, yet still was significantly different (P = 0.0052) (Figure 2A).

21To determine how the environment and physiological state influences microbial composition we compared the microbiome of laboratory-reared and field-caught (BG or G) individuals within each mosquito species. For all species, laboratory-reared mosquitoes had a significantly distinct microbiota compared to their field counterparts (Figure 2B). This was most pronounced in Ae. aegypti but less distinct for Ae. albopictus and Cx. quinquefasciatus. No significant differences were observed in the microbiome community structure of mosquitoes caught in the BG or G traps for any of the three mosquito species. No differences were observed when these microbiomes were analyzed with Wolbachia computationally excluded (Figure 2B).

22We examined our data for bacterial genera that were unique to or shared between mosquito species. The majority of bacteria were common between species with the notable exception that Ae. aegypti caught in the G trap possessed 12 genera not present in the other two mosquito species (Figure 3A and Supplementary Table S3). Similarly, when comparing within a species between groups, most bacteria were common to all groups (Laboratory, BG, G). To detect bacteria that may contribute to the observed differences in microbiome community structure of a particular mosquito species, we completed pairwise comparisons to identify bacteria that were differentially abundant (Figure 3B). These comparisons were completed within each group (Laboratory, BG or G). The largest differences in the microbiome were seen when comparing Ae. aegypti to the other two species, which agreed with our beta diversity analysis findings. Several bacterial taxa were differentially abundant between species regardless of group (Laboratory, BG and G traps), suggesting environmental factors are not greatly influential on these specific host-microbe associations. For example, Aeromonas, Serratia, Shewanella, and Wolbachia were less abundant in Ae. aegpyti compared to both Ae. albopictus and Cx. quinquefasciatus regardless of environmental conditions. When comparing differential abundance within a group, we were able to find infection gradients of specific microbes across species. The two examples of these infection clines were for Serratia and Aeromonas in laboratory-reared mosquitoes. These bacteria heavily infected Cx. quinquefasciatus, had moderate infection densities in Ae. albopictus and poorly infected or were absent from Ae. aegypti, despite the fact these three mosquito species were reared under common environmental laboratory conditions (Figure 3B and Supplementary Figure S5).

23Several studies have reported that blood feeding alters the species richness in mosquito guts. While at the community level, the microbiomes of G and BG were not significantly different (Figure 2), we did find specific bacteria that were differentially abundant between these groups (Figure 3C). These changes were mainly seen in Cx. quinquefasciatus. Of the known bacteria, the greatest changes were seen in Acinetobacter, Tatumella, and Pantoea, with the former two being more abundant in the mosquitoes caught in BG traps, while the latter was more abundant in mosquitoes caught in the G trap.

24While high throughput sequencing allows characterization of the composition of the microbiota, it only provides a relative measure of bacterial density. Therefore, to obtain an estimate of the total bacterial load in each vector species, we completed qPCR on mosquitoes with universal eubacterial primers that broadly amplify bacterial species. Culex mosquitoes were seen to have a higher total bacterial load when compared to either of the two Aedes species (Figure 4; Kruskal Wallis P < 0.0001). When comparing within a species between groups, we found that laboratory-reared Ae. aegypti had significantly greater load than those caught in the field (Kruskal Wallis P < 0.0001). This was also the case for Cx. quinquefasciatus (Kruskal Wallis P < 0.0001) while no significant differences were seen between Ae. albopictus groups. As both Ae. albopictus and Cx. quinquefasciatus are infected by Wolbachia, we also quantified Wolbachia by qPCR to determine its relative density in proportion to the total bacterial load of mosquitoes. While it is inappropriate to statistically compare amplicons that have different amplification efficiencies, it is evident that Wolbachia comprises a major component of the microbiome in Ae. albopictus and Cx. quinquefasciatus, which corroborates the high throughput sequencing data.

25The 16S rRNA sequencing data can be analyzed to create microbial interaction networks, providing information on potential interaction patterns of microbes such as co-occurrence and co-exclusion. We created network maps of bacterial interactions using a variety of models that use presence/absence and relative abundance data to identify pairwise relationships. In general, we saw that interaction networks from Ae. aegypti were complex, in that they had more nodes and connections compared to networks from Ae. albopictus and Cx. quinquefasciatus (Figure 5 and Supplementary Table S4). For all mosquito species, both co-occurrence and co-exclusion interaction patterns were observed in all networks (examples of these patterns are shown in Supplementary Figure S6). We were able to identify taxa, or groups of bacteria, within these interactions that appear to be important to the overall structure of the network. For example, in the Ae. aegypti networks from field collected mosquitoes, Enterobacter and Pseudomonas are highly interconnected species having between 6 and 15 interactions in these networks. In the laboratory-reared Ae. aegypti, an Enterobacteriaceae had several interactions with other bacteria. Three-way interactions were also seen in many of the networks, and often these involved Pseudomonas and Enterobacter. Examples of tripartite co-occurrences interactions in Ae. aegypti networks include Pseudomonas, Asaia, and a Clostridium isolate in laboratory-reared mosquitoes, Pseudomonas, Serratia, and Enterobacter and Pseudomonas, Acidovorax, and Enterobacter in the G and BG groups, respectively. Common co-exclusionary interactions were found in the networks generated from field-collected mosquitoes such as Pseudomonas co-excluding Pantoea and Tatumella. Burkholderia co-excluded Enterobacter, Acidovorax, and Escherichia-Shigella, however, these Burkholderia interactions could possibly be an artifact due to extraction batch variation.

26Wolbachia was a highly interconnected taxon in the interaction networks generated from Ae. albopictus and Cx. quinquefasciatus mosquitoes, and often had co-exclusionary relationship with other bacteria. In BG-collected Ae. albopictus, Wolbachia co-excluded six other bacteria including Asaia and Pseudomonas. In other groups, Wolbachia was seen to repeatedly exclude Aeromonas, Serratia, and Shewanella. The three-way interaction of Wolbachia co-excluding the co-occurring Aeromonas and Serratia was observed in both Ae. albopictus and Cx. quinquefasciatus. The Aeromonas and Serratia co-occurrence pattern appears highly robust and independent of environmental factors as this interaction was observed in five of the six Ae. albopictus and Cx. quinquefasciatus groups.

27The microbial interaction networks highlight the multifarious interactions that occur within mosquito systems. We undertook preliminary validation experiments using Ae. aegypti larvae to further demonstrate that microbial interactions can influence microbiome composition and abundance. To this end, we exploited the recently developed gnotobiotic rearing approach where mosquito larvae can be infected with a single bacterial taxon. We compared the density and prevalence of artificially infected symbionts in gnotobiotic lines (infected with the single symbiont) compared to conventionally reared septic mosquitoes (that possessed their resident microbiota) (Figure 6). We completed this with six bacteria isolated from Aedes mosquitoes. Three of the bacteria were isolated from Ae. aegypti (Pantoea, Cedecea [Cedecea-aeg], Asaia) while another three were isolated from Ae. albopictus (Serratia, Cedecea [Cedecea-alb], and Enterobacter). Three bacteria (Serratia; Mann Whitney P < 0.0001, Cedecea-alb; Mann Whitney P < 0.0001, and Cedecea-aeg; Mann Whitney P < 0.0001) were observed to significantly infect Ae. aegypti mosquitoes at higher densities when inoculated into axenic rather than non-axenic mosquitoes. The prevalence (number of individuals infected) of Serratia (Fisher’s exact test, P < 0.0001) and Cedecea-alb (Fisher’s exact test, P < 0.0002) was also significantly higher in gnotobiotic compared to conventionally reared mosquitoes. No change in either the prevalence or density of infection was seen for Enterobacter, Pantoea, or Asaia. Taken together, these data indicate that microbial interactions influence colonization and infection dynamics of specific bacterial species within mosquitoes.

28It is clear that complex factors combine to shape the microbiome of an organism. To further increase our understanding of factors that affect the microbiome of mosquitoes, we sequenced the microbiome of laboratory-reared and field-caught adult mosquitoes exploiting traps that attract host- or oviposition-seeking individuals. Anautogenous mosquitoes require a blood meal to provide the necessary nutrition for egg development and this dramatic influx of blood has been shown to substantially alter the microbiome of laboratory-reared mosquitoes or those caught in the field and blood-fed in the laboratory. Intriguingly, we saw few differences in the microbiome between host- and oviposition-seeking mosquitoes. Across all species, there were no significant differences in alpha or beta diversity for mosquitoes caught in either trap. Similarly, the total bacterial load was similar between mosquitoes caught in the BG and G traps, and this was consistent across all three species. However, when comparing beta diversity between traps, we did see less variation in the groupings of individuals from both Ae. albopictus and Cx. quinquefasciatus in the G trap compared to the BG trap. Overall, our results are in contrast to studies using lab-reared mosquitoes or mosquitoes caught in the field and blood-fed in the laboratory that indicate increased bacterial load, but decreased diversity over a 24–48 h window following a blood meal.

29Several factors could explain the differences between our results examining field-collected samples and those of previous studies. First, mosquitoes often take multiple blood meals, particularly Ae. aegypti. In our collections, the post-blood feeding history of mosquitoes is unknown, and it is possible that mosquitoes caught in either trap may have had a blood meal (or several) prior to being caught. The body of work examining the influence of blood feeding on the microbiome has only examined the effect of a single blood meal, not several, and while it appears the microbiome reverts to a pre-blood fed microbiome several days post-blood meal, it is unknown how quickly this occurs in the field. In field collected Anopheles mosquitoes reared in the laboratory, reversion to a Bacteroidetes-dominated microbiome, which was the dominant phylum seen in sugar fed mosquitoes, was seen 4 days post-blood meal. In our samples collected from G traps, females may have gone several days without a blood meal before finding a suitable oviposition site (i.e., the trap), possibly even relying on nectar-based food sources for sustenance. It is also possible that mosquitoes in search of an oviposition site may have never taken a blood meal. While most mosquitoes usually require a blood meal to develop eggs, autogeny has been reported in these species. Autogeny rates, which are variable and depend on temperature and nutrition, have been reported to range from 3 to 34% for Ae. aegypti, around 5% for Ae. albopictus, and up to 87% in Culex mosquitoes. Little is known regarding the influence of microbiota and autogeny, although in the autogenous mosquito Ae. atropalpus, specific gut taxa have been shown to influence egg production. Finally, while reports indicate that BG and G traps preferential catch mosquitoes in different physiological states, we did not explicitly examine if females were gravid or not. While we contemplated dissecting mosquitoes to examine their parity, we chose not to, as it would increase the potential for contamination of the samples. Although in our work we did not see dramatic differences in the mosquitoes caught in BG or G traps, the earlier points highlight the challenges in undertaking studies on field-derived samples and could explain the disparity between our results and those from studies undertaken in laboratory settings.

30Since both Ae. albopictus and Cx. quinquefasciatus are heavily infected with Wolbachia, we analyzed our data with and without this endosymbiont to garner a better understanding of the other bacterial microbiota in these mosquitoes. Our high throughput sequencing and qPCR results demonstrated that Wolbachia was the most abundant bacterium in Ae. albopictus and Cx. quinquefasciatus. In spite of this, we obtained sufficient sequencing depth to identify other bacterial taxa – a common challenge when characterizing Wolbachia-infected species using amplicon sequencing. Removing Wolbachia from our analysis increased the Shannon diversity index, indicating the remaining microbiota within these mosquitoes is relatively even. For beta diversity, we found that our results were mixed and dependent on group. When comparing within a group, removing Wolbachia made the microbiomes of laboratory-reared Ae. albopictus and Cx. quinquefasciatus more divergent while the microbiomes of field mosquitoes tended to be less different. These findings are consistent with a study by that found the removal of Wolbachia from their analysis led to less distinct differences for field-collected mosquitoes. When comparing groups within a species, exclusion of Wolbachia made the comparison of field and laboratory mosquitoes more distinct, likely due to differences in the gut-associated microbiota. Similar to a recently published study, we found a small number of Wolbachia reads in Ae. aegypti mosquitoes collected in G traps. However, we could not confirm the presence of the bacteria with conventional PCR-based approaches, suggesting these results could be a sequencing or laboratory artifact. A recent PCR-based screening of Ae. aegypti mosquitoes found no evidence for Wolbachia infection, although cryptic low titer infections may be beyond the detection threshold of conventional PCR. Should populations of Ae. aegypti be naturally infected with Wolbachia, this could have important ramifications for biological control strategies being implemented into the field, and as such, further research in this area is warranted.

31Similar to findings in other mosquito species, many bacteria were shared between different mosquito species, and laboratory-reared mosquitoes were seen to have a divergent microbiome from their field counterparts. These common taxa, particularly those that infect at high abundance, could be candidate bacteria for consideration in novel pan-mosquito microbial control strategies as they would likely be compatible for all three vector species. When focusing specifically on bacterial titers, few bacteria were seen to be significantly different between mosquitoes caught in either traps. In Culex mosquitoes, Pantoea was more abundant in individuals caught in the G trap. Similar to this finding, Pantoea has been found to increase in abundance after a blood meal in Anopheles mosquitoes. When considering pairwise comparisons between species, there were more differentially abundant genera when comparing Ae. aegypti to the other two species, which explains the greater divergence of the Ae. aegypti microbiome to the microbiome of the other two species. The Ae. aegypti samples were collected from a different location than the other two species, thus, it is possible that environmental factors could explain these differences. However, we also saw common changes that were consistent regardless of groups (BG, G and Laboratory), such as Aeromonas, Clostridium, Serratia, Shewanella, and Wolbachia, indicating these bacteria were not greatly influenced by the environment and that other factors affected their presence in the particular mosquito species. Notably, we found that in laboratory-reared mosquitoes, Culex harbored significantly higher titers of Serratia and Aeromonas, compared to Ae. albopictus, while these bacteria were at low abundance or absent from Ae. aegypti. These results suggest there are host and/or bacterial related factors that make this particular strain of Ae. aegypti (Galveston) inhospitable for Serratia and Aeromonas as all three mosquito species were subjected to similar uniform environmental conditions when reared in the insectary.

32Current evidence of microbial interactions within mosquitoes is mainly limited to interactions between Wolbachia and other microbiota, or between vertically transmitted symbionts in other arthropod systems. As such, our understanding is generally restricted to inherited symbionts and we have a poor understanding of the scale of interactions between microbes that are associated with insect guts. To address this, we created microbial interaction networks to identify pairwise co-occurrence and co-exclusion patterns. To avoid spurious interactions, which could be due to the presence or absence of a microbe in one environmental condition but not another, we limited our network analysis to within a group for each particular species. Our analysis identified 116 co-occurrence or co-exclusion interactions, substantially increasing the number of bacterial interactions observed in mosquitoes. Bacterial interaction networks generated from Ae. aegypti mosquitoes were more complex than Ae. albopictus or Cx. quinquefasciatus in that they had more nodes and connections. Species richness may explain the differences observed in network structure, as in general, the more complex networks had a greater number of OTUs. Other factors that may have influenced the identification of interacting bacteria are the presence of Wolbachia in Ae. albopictus and Cx. quinquefasciatus, as well as the differences in sample size between mosquito species. Further studies are warranted to determine why some networks are highly interconnected while others are not. Approaches that artificially manipulate the microbiome composition of a mosquito may be particularly useful in this regard.

33Interestingly, in the more complex networks, we saw evidence of hub microbial taxa that were highly interconnected. Pseudomonas and bacteria within the Enterobacteriacaea appear to be important hub taxa. Some of the interactions observed here have been previously reported in Ae. triseriatus and Ae. japonicus, including negative interactions between bacteria within the Burkholderiaceae and Pseudomonas and Acinetobacter as well as Asaia-Enterobacter and Asaia-Ralstonia interactions. The majority of interactions reported in Ae. triseriatus and Ae. japonicus were negative, whereas here we see a mix of both co-occurrence and co-exclusion patterns. Specific hub microbes that are strongly interconnected have been found in plant microbiomes and these taxa have a profound effect on overall microbiome structure. In our data, we also saw three-way interactions. Often these relationships were co-occurring, or were formed by two co-occurring bacteria both co-exclude another bacterium, suggesting the existence of multi-taxa interactions. Further work is required to determine the functionality of these multi-interacting partners and if these interactions represent keystone guilds. The identification of common interaction pairings across several groups indicates that these interactions, and the methods we employed to identify them, are robust and presumably not influenced by physiological state or other environmental conditions.

34It is important to highlight that these network maps represent patterns, and not direct interactions. Many of the observed interactions may be due to microbes sharing a similar ecological niche, and a substantial challenge, particularly for highly interconnected taxa, will be to investigate these interactions further. To undertake initial validation steps and to demonstrate that microbial interactions are an important factor influencing the colonization of gut-associated microbiota, we infected six culturable bacterial taxa into Ae. aegypti larvae that either possessed or lacked their resident microbiota. We used larval infections, rather than adults, as we could exploit the use of the gnotobiotic infection procedure in this life stage. Serratia and Cedecea, which were isolated from Ae. albopictus, poorly infected Ae. aegypti when it possessed its native microbiome. However, when mosquitoes lacked their native microbiota, these bacteria infected at a higher titer. Even more striking was the effect on prevalence (number of individuals infected), which changed from 38 and 50% for Serratia and Cedecea-alb, respectively, in conventionally reared mosquitoes that possessed their resident microbiota, to 100% infection when infected into axenic larvae. This indicated that microbial incompatibility, rather than host or bacterial genetics, was a major factor causing the poor infection rates of mosquitoes. The ability of Serratia and Cedecea-alb to infect some conventionally reared mosquitoes that possessed their native microbiota could be due to variation in the microbiome between individuals, which is likely caused by stochastic processes. In our networks, Serratia and Cedecea have several co-exclusionary relationships with dominant bacterial taxa such as Asaia, Pseudomonas, and Enterobacter, which may explain these results, however specific examination of these interactions in adult mosquitoes is required. Importantly, not all of the bacterial taxa artificially infected into larvae increased in prevalence and density suggesting these enhanced colonization effects are not simply due to mosquitoes lacking their microbiota, but rather are specific for each taxa, likely due to specific microbial interactions. Similar to our findings, it has been reported that antibiotic treatment prior to bacterial supplementation in a sugar meal can increase the prevalence of infection of gut microbes in female Anopheles and Aedes mosquitoes, indicating that resident gut bacteria that are susceptible to antibiotics are antagonistic to the supplemented bacterium. It is important to note that our reinfection study, which exploited the gnotobiotic rearing system, examined interactions in larvae, not adults, and that differences in the gut morphology and function between these two life stages may alter microbial interactions. However, our findings combined with the work of, suggest that microbial interactions between gut-associated bacteria occur within mosquitoes and influence symbiont colonization in aquatic and adult life stages, which likely affects microbiome species richness and evenness. These colonization traits and co-exclusionary associations could offer a possible explanation for the variability seen in the mosquito microbiome between individuals, as bacteria that initially infect the gut may impede colonization by other microbes.

35In this study, we compared microbial interaction networks from field and laboratory mosquitoes to examine the influence of the general environment but it is possible other factors may influence network structure. In particular, it would be interesting to determine if microbial networks differ between tissues within mosquitoes, such as the salivary glands, germline, and the gut. For example, microbial network analysis from the human microbiome project found strong niche specialization in their networks, whereby different body sites had contrasting microbial networks. In Anopheles mosquitoes, salivary glands have a more diverse microbiome compared to the gut, and it is conceivable that elevated species richness would allow for greater network interactions. Furthermore, the germlines of male and female Anopheles mosquitoes share some common taxa but there are also quantitative differences. These differing microbial niches could be exploited to determine the influence of species evenness on microbial networks. Performing network analysis on gut samples of Ae. albopictus and Cx. quinquefasciatus may also overcome any issue with Wolbachia sequestering the majority of the reads, as Wolbachia primarily resides within the germline in these mosquito species. Alternatively, approaches that reduce amplification of Wolbachia could be implemented. Here, we assessed whole mosquitoes to give an initial overall picture of microbial interactions, but analysis of distinct tissues may identify interactions of bacteria that are proximal to one another, and these interactions are more likely to reflect true microbe-microbe interactions, rather than patterns associated with environmental exposure. Within a species we collected mosquitoes from a single site, but future studies examining interaction networks should incorporate diverse sites. Common pairwise interactions identified across sites would indicate robust relationships not influenced by environmental factors.

36In summary, we examined the microbiome of three important mosquito vectors, Ae. aegypti, Ae. albopictus, and Cx. quinquefasciatus. While the overall microbiome structure between host-seeking or ovipositing females was similar, we identified specific bacteria that changed in abundance between mosquito species. Our analysis identified a suite of pairwise interactions used for generating microbial interaction networks, and together with re-infection studies we have demonstrated that microbial interactions affect microbiome composition and abundance of specific bacterial taxa. These findings add to our understanding of microbiome community structure of mosquitoes and factors that influence microbiome acquisition and maintenance in these important disease vectors.

37The NCBI accession number for the raw sequencing data reported here is PRJNA422599.

38SH, KK, LA, JD, EH, YF, and GH designed the experiments. GM, CF, JD, and MD oversaw the fieldwork. SH and MP completed the experiments. KK, MR, LA, SH, EH, and GH undertook analysis. SH, KK, JD, MD, and GH wrote and edited the manuscript. All authors agreed to the final version. YF and GH acquired funding and supervised the work.

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