Variable microbiomes between mosquito lines are maintained across different environments

The composition of the microbiome is shaped by both environment and host in most organisms, but in the mosquito Aedes aegypti the role of the host in shaping the microbiome is poorly understood. Previously, we had shown that four lines of Ae. aegypti harbored different microbiomes when reared in the same insectary under identical conditions. To determine whether these lines differed from each other across time and in different environments, we characterized the microbiome of the same four lines of Ae. aegypti reared in the original insectary and at another institution. While it was clear that the environment influenced the microbiomes of these lines, we did still observe distinct differences in the microbiome between lines within each insectary. Clear differences were observed in alpha diversity, beta diversity, and abundance of specific bacterial taxa. To determine if the line specific differences in the microbiome were maintained across environments, pair-wise differential abundances of taxa was compared between insectaries. Lines were most similar to other lines from the same insectary than to the same line reared in a different insectary. Additionally, relatively few differentially abundant taxa identified between pairs of lines were shared across insectaries, indicating that line specific properties of the microbiome are not conserved across environments, or that there were distinct microbiota within each insectary. Overall, these results demonstrate that mosquito lines under the same environmental conditions have different microbiomes across microbially- diverse environments and host by microbe interactions affecting microbiome composition and abundance is dependent on environmentally available bacteria.


Introduction
The microbiome plays a crucial role in the health of various organisms [1,2].Alterations in the relative abundance and overall bacterial community structure can lead to dysbiosis within the organism resulting in disease or morbidity.Factors that shape the acquisition and maintenance of the microbiome vary between organisms and include contributions from both the environment and the host.The role of host genotype in shaping the composition of the microbiome has become apparent in mammalian systems where specific host genomic loci have been associated with specific bacterial taxa [3][4][5].Additionally, the role of host genetic background in invertebrates such as Drosophila and other insects has also been shown to contribute to the composition of the microbiome [6][7][8][9].While the contribution of host genetic background in shaping the microbiome is well studied in mammalian and Drosophila systems, few and contradictory data exist for the role of different mosquito lines, and perhaps different genetic backgrounds, in shaping the microbiome of the mosquito [10][11][12].
Aedes aegypti is the main vector of arthropod-borne viruses (arboviruses) worldwide, such as dengue, Zika, and chikungunya.Examples such as dengue virus (DENV) result in 100-400 million infections annually and remain a major threat to public health [13,14].The microbiome, specifically the presence of distinct isolates, plays a role in the ability of Ae. aegypti to be a successful vector of human pathogens [15,16].The microbiome of Ae. aegypti is largely shaped by the environment [17,18], but the role of the host remains poorly understood.Multiple studies have controlled for environmental variation and characterized the microbiome of diverse lines of Ae. aegypti reared in the same insectary (i.e. the same environmental conditions).These studies have found contradictory results.Two independent studies identified differences in the microbiome from the whole body of Ae. aegypti that were dependent on mosquito line [10,11].A third study did not observe any differences in either the bacterial community structure or the abundance of specific taxa in the midgut from a selection of genetically diverse lines of Ae. aegypti, which represented their worldwide genetic diversity [12].These contradictory results suggest that perhaps the environment is important for detecting line specific differences in the microbiome, or alternatively, there are host and environmental interactions that determine the microbiome composition in Ae. aegypti.
Aedes aegypti occupies a variety of environments worldwide, allowing for the association of numerous bacterial taxa with it.Understanding the relative contributions of the mosquito and the environment in shaping the microbiome composition is important for both teasing apart the role of the microbiome in mosquito vectorial capacity, [15] as well as development of paratransgenic control tools [19][20][21].To understand whether line specific differences in the microbiome is dependent on the environment, we sequenced the 16S ribosomal RNA gene from four lines of Ae. aegypti that previously showed line specific differences in the microbiome [10] and the same four lines after being transferred to a new insectary at another institution.We analyzed the structure of the bacterial community between and within each insectary, as well as the differential abundance of specific taxa.Finally, we identify conserved genera that differ in pairwise comparison between lines at each insectary.Our findings that lines harbor differences in their bacterial composition despite being reared under identical conditions has important implications when comparing phenotypic effects which may be sensitive to the microbiome.

Ae. aegypti Mosquitoes
Colonies of Ae. aegypti used in this study originated from Galveston (USA), Iquitos (Peru), Juchitan (Mexico), and Thailand.At UTMB, the generation of the colonies is not known but the colonies have been in the insectary since 2010 and reared continuously.The colonies were transferred to LSTM in 2018.At each institute the mosquito lines were housed under standard insectary conditions consisting of 28˚C +/-1˚C and 70% +/-10% relative humidity with 12h:12h light dark cycle.Eggs were hatched in deionized water and larvae were fed fine ground tetramin fish food at both insectaries.Adults were held in Bugdorm cages at a density of 300-600 mosquitoes with constant access to 10% sucrose until being harvested.Mosquitoes were not bloodfed for this study.At LSTM colonies are fed fresh human blood from the NHS for egg production.At UTMB colonies are fed fresh defibrinated sheep blood Colorado Serum Company for egg production.

DNA Extractions
Three-five days post emergence, Ae. aegypti were cold anesthetized and females taken for surface sterilization.Individual mosquitoes were surface sterilized in ethanol 70% for 5 minutes followed by 3 washes in sterile PBS.DNA was extracted from 20 female mosquitoes from Galveston, Iquitos, Juchitan and Thailand using the QIAamp DNA Mini Kit (QIAGEN) following the manufacturer protocol with the following modifications: initial volume of 180ul of buffer ATL used for mosquito homogenization and final volume of 100ul nuclease-free water for DNA elution.No RNase A treatment was applied.For the UTMB samples, DNA from mosquitoes from Iquitos and Juchitan was extracted on 25 February 2021 and DNA from mosquitoes from Galveston and Thailand was extracted on 26 February 2021.No-mosquito controls were used for each extraction batch performed at UTMB and sequenced as negative controls.

Library Preparation and sequencing
Sequencing libraries for each isolate were generated using universal 16S rRNA V3-V4 region primers [22] in accordance with Illumina 16S rRNA metagenomic sequencing library protocols.DNA concentrations of each library were determined by Qubit and equal amounts of DNA from each barcoded library were pooled prior to sequencing.The samples were barcoded for multiplexing using Nextera XT Index Kit v2.The pooled libraries were diluted to 4 pM and paired-end 250 bp sequencing was performed on the Illumina Miseq using a MiSeq Reagent Kit v2 (500-cycles).

Data analysis
To identify known bacteria, sequences were analyzed using the CLC Genomics Workbench 21.0.5 Microbial Genomics Module (CLC MGM).Paired end reads were merged as part of the CLC Genomics workflow.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.Reference-based Operational Taxonomic Unit (OTU) picking was performed using the SILVA SSU v132 97% database [23].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 an 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.OTUs with a combined abundance of less than two were removed from the analysis.Low abundance OTUs were removed from the analysis if their combined abundance was below 10 or 0.1% of reads.The number of reads per sample used in the analysis ranged from 2,329-136,685.Only reads that mapped to bacteria were kept.Taxa classified as "Ambiguous Taxa" are reads mapping to bacterial DNA, but that cannot be identified at the taxonomic level.
Sequences obtained from the no-mosquito negative controls were used to filter out contaminants using the Decontam program [24].Read counts of samples and negative controls were analyzed to identify OTUs with a high prevalence in the negative controls, which were marked as contaminants.The default probability threshold of 0.1 was used.A total of six OTUs were removed from the dataset.
Analysis of the OTU table was performed using MicrobiomeAnalyst [25,26].The analysis parameters were set so that OTUs had to have a count of at least 10 in 20% of the samples and above 10% inter-quantile range.Reads were normalized by Total Sum Scaling.Analysis was performed using actual and total sum scale abundances.Alpha diversity was measured using the observed features to identify the community richness using Chao1.Statistical significance was calculated using T-test/ANOVA.Beta diversity was calculated using the Bray-Curtis dissimilarity measure (genus level).Permutational Multivariate Analysis of Variance (PERMA-NOVA) analysis was used to measure effect size and significance on beta diversity for grouping variables [27].Relative abundance analysis was done in MicrobiomeAnalyst at the level of genera.Pairwise differential abundance of specific genera was done in MicrobiomeAnalyst.Statistical significance between groups was determined by T-test/ANOVA and corrected for multiple testing.Presented p-values reflect correction for multiple testing.(S1 Table ).
Rarefaction curves showed that adequate sequencing depth was not achieved in multiple samples (S1 Fig) .Following rarefication of the data, the diversity analysis was re-done, and the beta diversity plot was highly similar to the beta diversity plot generated with non-rarified data (S2 Fig) .Given that the diversity analysis was the same after rarefication, the analysis was completed on the non-rarified data to include rare microbes in an effort to represent the bacterial community as best as possible.

Upset plot
To identify the distribution of bacteria taxa between the different mosquito lines and their rearing insectaries, count data from female mosquitoes were extracted and aggregated to the genus level for each sample.This yielded a maximum of 44 distinct bacterial taxa (including 'uncultured' and 'ambiguous' taxa).This table of count data was then processed using the R package UpsetR [28], stratifying the data based on mosquito strain and the rearing insectary.

Correlation plot
Previous studies have indicated a potential correlation in presence/absence between specific genera of bacteria-namely Cedecea, Enterobacter, Klebsiella, and Serratia.To test for this, counts from the four genera were extracted from the aggregated table described earlier.Using the sum of these counts as the total, percentage relative abundance for each taxa was then counted per individual sample.These results were then passed to the cor() function in R's stats package [29] to obtain a correlation matrix of size four-by-four.This correlation matrix was then visualised using the R package corrplot [30], with the two options method = "color", col = COL2('PiYG').

Between insectary differences
We previously reported that the microbiomes of four Ae.aegypti lines reared in the same insectary were different [10].To confirm whether these differences were conserved across subsequent generations and environments, we characterized the microbiome of these lines in the insectary that they have been continuously reared in (UTMB) and an insectary at a different institution (LSTM) after being transferred and reared for multiple generations.To compare the microbiome diversity and composition between lines we undertook amplicon sequencing on the V3-V4 region of the bacterial 16S rRNA gene from 20 individuals from each line.Out of the 160 individuals sequenced, a total of 63 OTUs were found which represented 44 bacterial genera.
To determine if the microbiome of each Ae.aegypti line differed in diversity of bacterial species present, the within line diversity was determined by calculating the Chao Diversity index (Fig 1A).The richness of the microbiome was greater in mosquitoes reared in the LSTM insectary compared to those reared in the UTMB insectary (ANOVA, p-value < 0.001).To determine if the community structure of any one line was more closely related to the same line reared at a different insectary, or different lines reared within the same insectary, principal component analysis (PCA) was performed on a Bray-Curtis dissimilarity matrix.Overall, individuals from each line shared a more similar bacterial community structure with other lines from the same insectary (Fig 1B ), indicating that the environment and non line-specific factors plays a more influential role in shaping the community structure of the microbiome.

Within insectary differences
To determine whether the lines differed from each other within the same insectary, alpha and beta diversity was compared between the lines from each insectary.Differences in species richness were observed between the lines at both the LSTM insectary ( To measure the extent of overlap of bacteria identified between the two insectaries, the presence of each genera was compared between each line from each insectary (Fig 3).Out of the 44 identified genera, 23 were present in both insectaries.Five genera are specific to the UTMB insectary, and three genera were specific to the LSTM insectary (S2 Table ).We examined the relative abundance of the 20 most abundant genera for each line to determine if the bacterial communities were dominated by the same bacterial genera in each line at both insectaries.Different bacteria genera dominated mosquitoes reared in the LSTM insectary compared to the UTMB insectary (Figs 4 and S4).All four lines reared in the LSTM insectary were dominated by an ambiguous taxa and Perlucidibaca.The four lines reared in the UTMB insectary had a high proportion of Acinetobacter and Asaia.
We previously reported a correlation between Enterobacteriaceae bacteria and Serratia in Ae. aegypti [10].Specifically, if the midgut of Ae. aegypti was colonized with Enterobacteriaceae, Serratia was excluded from colonizing.The robustness of this phenotype after multiple generations and in different environments was measured.We found that this negative correlation between the presence of Enterobacteriaceae and the presence of Serratia was still observable in these lines (S5 Fig).

Environment versus line specific differences
To understand if line specific differences are conserved across environments, pair-wise differences in the abundance of each bacterial genera were calculated between lines within and across insectaries (S1 Table ).Out of the 44 genera identified, a range of 10-37 genera were found in different abundances between lines within each insectary resulting in 25-84% similarity between the lines (Fig 5).To assess whether line specific factors or the environment play a larger role in shaping the microbiome, the percentage of genera differentially abundant was compared for each pairwise comparison within and between insectaries.Lines were more similar to lines from the same insectary than to the matching line in the other insectary (Fig 5), indicating the environment is more important than line specific factors in shaping the microbiome.
To assess if lines differ from each other in the same way across environments, pairwise comparisons of lines were performed within each insectary, and genera that showed statistically differential abundances were identified per comparison.If we identify conservation of differentially abundant genera between the same two lines reared in different environments, this would suggest an interaction between specific bacterial genera and the hosts' background across environments.Overall, a range of 4-14 genera were identified as having conserved differences in their abundance between the same lines in each insectary (Fig 6).The proportion of differentially abundant genera between lines that were shared between institutes is dependent on the specific pairwise comparison (Fig 6 and S3 Table ).The pairwise comparisons including the Thailand line shared fewest conserved differentially abundant genera across insectaries.The other pairwise comparisons that didn't include the Thailand line showed at least half of the differentially abundant bacteria are conserved across insectaries, indicating that the structuring of the microbiome between lines from the same environment is partially conserved across environments.

Discussion
We explored the contribution of mosquito line and environment in shaping the microbiome composition of Ae. aegypti.Having previously observed differences in the microbiome Rows are colored based on mosquito strain (middle-left), and further divided into the rearing Insectary (bottom left).23 taxa are shared across all the different strains between the different insectaries, constituting a potential 'core' set of bacteria.This is then followed by five (Persicitalea, Janthinobacterium, Rahnella, Luteolibacter, Verrucomicrobium) and three taxa (Sphingopyxis, Burkholderia-Caballeronia-Paraburkholderia, Methyloversatilis) that appeared unique to the rearing insectary.All other permutations of intersects contain two or fewer bacterial taxa.https://doi.org/10.1371/journal.pntd.0011306.g003between four lines of Ae. aegypti reared under identical conditions in the same insectary [10], we sought to determine if these lines continue to have distinct microbiomes since our previous work, and if this was consistent after rearing in different insectary environments.In accordance with our previous study [10], we observed differences in alpha and beta diversity between the four lines and interestingly, the microbiome of these lines were still distinct after rearing in separate insectaries.Notably, we found that although the microbiomes of the lines differed from each other within an insectary, any one line was more similar to other lines reared within the same insectary as compared to the same line at the other insectary.Finally, we observed a percentage of differentially abundant genera between lines were shared between insectaries.Together these data demonstrate that while the environment influences the microbiome available to the organisms, different lines of Ae. aegypti can harbor different microbiomes under the same environmental conditions.
A possible source of lines specific differences in the microbiome is the host genetic background.The role of host genotype in determining the microbiome composition has been established in other systems [3][4][5][6][7][8][9].Although the influence of mosquito genotype in shaping the microbiome in the same environment is not well understood, multiple mosquito genes have been identified that influence the composition of the microbiome and gut homeostasis [11,[31][32][33][34]. Genes involved in bloodmeal digestion and immune factors can regulate the abundance of the microbiome as a whole [34], or the abundance of specific taxa [31].Furthermore, the mosquito microbiome can stimulate expression of specific genes to shape the immune status of the mosquito allowing for efficient colonization of specific microbes [33].In addition to immune genes regulating microbiome composition, metabolic signaling and nutrition status of the mosquito can influence the microbiome [11,15].Immunity and nutritional processes are under genetic control and different genotypes of Ae. aegypti could result in differences in levels of micro-and macronutrients, which in turn could influence the ability of specific bacteria to colonize that mosquito line.Additionally, genetic variation in immune genes could influence the success of different bacteria in colonizing the mosquito.
Other possible factors could give rise to different microbiomes between lines reared under identical environmental conditions, such as random events that are selecting for particular microbes over multiple generations.Each mosquito line is bringing along their own microbiome on their egg paper which is subjected to changes during each generation.These random events could enforce drift of the microbial populations between lines or selection of specific microbes in the population with each serial passage of the colony.All of this could be happening independent of the mosquito factors or mosquito genotype.
If mosquito line played a more profound role than the environment in determining which genera colonize the mosquito, we would expect the same line at the two insectaries to share more similarities to each other than the other lines within the same insectary, assuming similar taxa were available to the mosquito across insectaries.Additionally, this would be apparent if pairwise differences in taxa abundance between any two lines within one insectary were conserved in another insectary.We observed the opposite of this where lines within one insectary had more similar abundances of specific genera to each other, compared to the same line reared at a different insectary.We also observed that pairwise differences in the abundance of different genera between lines within one insectary were partially shared in another insectary.The incomplete shared differences in the abundance of specific genera within a line across insectaries suggest that the environment, or the microbes the mosquito has access to, plays a larger role than host specific factors in determining the microbiome.Given that only 23 of the 44 genera identified were shared across insectaries, it is not surprising that line specific abundances of specific taxa are lost across insectaries.This study and previous studies observed differences in the microbiome between the lines when reared under identical environmental conditions in the same insectary [10].This is in contrast to our previous work in which midguts from diverse genotypes of Ae. aegypti harbored the same microbiome when reared under identical conditions [12].The discrepancy between these studies could be a result of what tissues were used for microbiome characterization.Kozlova et al. [10] and this study sequenced the microbiome from the whole mosquito, while Dickson et al. [12] sequenced the midgut microbiome.Additionally, it is possible that the environmental bacteria present in the insectary in Dickson et al. [12] could colonize the midguts of all the lines effectively, while this was not the case in the UTMB and LSTM insectaries.Microbe-microbe interactions are important for determining the microbiome composition [17] and perhaps some of the environmental bacteria present in both the UTMB and LSTM insectaries are incapable of colonizing Ae. aegypti.
Given that the lines of Ae. aegypti were not sequenced in this study we cannot conclude these lines represent different genotypes of Ae. aegypti.The colonies used in the study represent a similar genetic background in relation to global populations of Ae. aegypti [35,36].Without sequencing the different lines to confirm genetic differences, we cannot unequivocally determine if the differences we see in the microbiome composition are correlated with mosquito genotype.Although we expect them to be similar genetically, it is likely that there is genetic divergence between the lines and it is reasonable to conclude that the lines are genetically diverse from each other.
The results from this study provide a unique opportunity to tease apart host control of the microbiome in Ae. aegypti.By exploiting our new cryopreservation method [37] and swapping the microbiome using our recently developed transplantation approach [38] from both lines and insectaries, the role of the host in controlling the microbiome could be further investigated and provide future insight into host genetic mechanisms that underly microbiome composition.Additionally results of this study highlight the need to take into account environmental variation and which mosquito lines are being used when characterizing the microbiome, which could have important implications for developing novel control strategies such as paratransgenesis and phage therapy.
Fig 2A) (ANOVA, pvalue < 0.001) and the UTMB insectary (Fig 2C) (ANOVA, p-value < 0.001).To gauge if the lines differed in the bacterial community structure, PCA was performed on a Bray-Curtis dissimilarity matrix.The PCA plots show that the bacterial communities differ between the lines at both the LSTM insectary (Fig 2B) (PERMANOVA, p-value: 0.001) and the UTMB insectary (Fig 2D) (PERMANOVA, p-value: 0.001).Strikingly, these results indicate that differences in the bacterial community structure of these four lines are maintained across generations and insectaries.To confirm the bacterial community structure differed between lines at each insectary, pairwise beta diversity analysis was performed.The 12 PCA plots from each pairwise comparison show that bacterial community structure differs between each pair of lines at each insectary (S3 Fig)

Fig 1 .
Fig 1. Diversity of the microbiome in individuals reared at the LSTM or UTMB insectaries.Structure of bacterial communities was determined by deep sequencing the V3-V4 region of the bacterial 16S gene in adults from 4 different lines of Ae aegypti (Galveston, Thailand, Iquitos, Juchitan) reared in two different insectaries at UTMB and LSTM.The Bacterial community structure is represented (A) by the species richness index Chao1 and (B) by principal component analysis of Bray-Curtis dissimilarity index.Mosquito lines reared at LSTM are shown in purple, mosquitoes reared to UTMB are shown in turquoise.https://doi.org/10.1371/journal.pntd.0011306.g001

Fig 3 .
Fig 3. Upset plot showing the number of genera that are shared between lines in the LSTM and UTMB insectaries.Each mosquito strain in the different insectary (one per row at the bottom half of the image) is treated as a 'set' with an identified number of bacterial taxa ('Set Size').The various permutations of intersections are denoted by the ball-and-stick diagram at the bottom of the image, and size of these intersections denoted by the bar graph at the top of the image ('Intersection Size').Rows are colored based on mosquito strain (middle-left), and further divided into the rearing Insectary (bottom left).23 taxa are shared across all the different strains between the different insectaries, constituting a potential 'core' set of bacteria.This is then followed by five (Persicitalea, Janthinobacterium, Rahnella, Luteolibacter, Verrucomicrobium) and three taxa (Sphingopyxis, Burkholderia-Caballeronia-Paraburkholderia, Methyloversatilis) that appeared unique to the rearing insectary.All other permutations of intersects contain two or fewer bacterial taxa.

Fig 4 .
Fig 4. Relative abundance of bacteria in each line.The dominant bacterial genera are different between insectaries and between lines.The relative abundance of the 20 most abundant genera in shown for 20 individuals from each line at LSTM and UTMB insectaries.Bacterial genera were assigned to OTUs clustered with a 97% sequence identity cutoff and taxonomically classified with the SILVA database.https://doi.org/10.1371/journal.pntd.0011306.g004

Fig 5 .
Fig 5. Pairwise comparisons of mosquito lines reared in each insectary.The microbiomes of mosquito lines reared in the same insectary are more similar compared to those reared in a different insectary.Results from a pairwise differential abundance analysis are shown for each pair of lines as the percent of genera that are significantly different between the pairs after correcting for multiple comparisons.A light blue color indicates a higher degree of dissimilarity between the lines.https://doi.org/10.1371/journal.pntd.0011306.g005

Fig 6 .
Fig 6.Few genera have conserved differences in abundance between lines at both insectaries.Venn diagrams show the overlap in specific genera that are differentially abundant between each pair of lines (Galveston, Iquitos, Thailand, or Juchitan) at each insectary (LSTM or UTMB).The ID of genera differentially abundant between each pairwise comparison was compared between insectaries.The shared genera represent a genus that is differentially abundant between the lines in both environments.The ID of shared genera can be found in S3 Table. https://doi.org/10.1371/journal.pntd.0011306.g006