**1. Introduction**

The intestines of animals are occupied by diverse communities of microorganisms that can affect different aspects of host health. The microbiota plays a key role in many

**Citation:** Beribaka, M.; Jeli´c, M.; Tanaskovi´c, M.; Lazi´c, C.; Stamenkovi´c-Radak, M. Life History Traits in Two *Drosophila* Species Differently Affected by Microbiota Diversity under Lead Exposure. *Insects* **2021**, *12*, 1122. https:// doi.org/10.3390/insects12121122

Academic Editor: Nickolas G. Kavallieratos

Received: 8 November 2021 Accepted: 13 December 2021 Published: 15 December 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

aspects of host life, including development, digestion, behavior and the immune system [1–3]. Simple gut microbiota of *Drosophila* species, short generation time and high fecundity are some of the reasons that make them an attractive model for studying the significance of gut microbiome from different aspects. *Drosophila* hosts only a small number of bacterial populations in its gut, but includes species present in the human microbiota as well. *Drosophila* gut microbiota in laboratory is represented by a low-diversity bacterial community [4,5], but it has great implications on its overall health. The gut microbiome of *Drosophila* contributes to a variety of host traits, such as innate immunity [6], lifespan [7–9], nutrition and reproduction [10] and behavior [11,12]. Shifts in microbiota could lead to serious consequences on host physiology, causing even mortality [13]. Thus, it is important to investigate the factors that shape the composition and diversity of microbiota and their possible implications on the host.

One of the greatest problems that animals currently face within the natural environment is pollution. Due to the anthropogenic factors, the presence of pollutants is widespread in the air, soils and water. Heavy metals are pollutants commonly found in nature, with great impact on plants and animals. Recently, it has been reported that heavy metals, such as arsenic, cadmium, lead and mercury, have detrimental effects on the diversity of terrestrial invertebrates at levels below those considered safe for humans [14]. The population-level phenotypic variability of different species depends on the interaction between genetic and environmental variability. The variability in the components of adaptive value (fitness) is an aspect of phenotypic variability. Lead is one of the widespread heavy metals that has previously been reported to have a major negative impact on *Drosophila* fitness [15]. However, it has been shown that its negative impact can be modified depending on the population genetic background [16], genome heterozygosity [17] and genetic variation [18]. Since *Drosophila* species are mainly exposed to lead through food intake, it was suggested that gut microbiota could also have an impact on the species resistance to lead toxicity. Our previous research suggested that bacterial diversity increased in two *Drosophila* species after extended exposure to a lead-saturated substrate. This increase in bacterial diversity underlined certain bacterial genera, such as *Komagataeibacter* and *Acetobacter*, that could be good lead-tolerant members of microbiota [5].

The results obtained in our previous study suggested the difference in shifts of microbiota composition between natural populations of two *Drosophila* species under laboratory conditions on standard and lead-saturated substrate after 13 generations [5]. In the present paper we investigate the microbiota of the same lab-reared *Drosophila* species (*D. melanogaster* and *D. subobscura*), but from each of the two distinct localities on standard and lead-saturated substrate after an additional 22 generations (35 in total). Further, we explore if the changes in life history traits (egg-to-adult viability and developmental time) are associated with microbiota content enabling different responses to lead (Pb) in those experimental groups. The influence of lead exposure was discussed for the population, species, substrate and sex levels for both microbiota and fitness components. The microbiota was analyzed using NGS sequencing of the V3-V4 16S rRNA gene and assessed through diversity indices and taxonomical analysis.

The aim of this study was to determine the impact of prolonged lead exposure in two lab-reared *Drosophila* species, each from two populations on microbiota composition and life history traits, to find potential cause-and-effect relationships between them and to differentiate the response of different origin, species and sex to stress factors.

#### **2. Materials and Methods**

#### *2.1. Sample Collection and Laboratory Maintenance*

In this study, both species *(D. melanogaster* and *D. subobscura*) were sampled from two localities in Serbia: Kalna (43.4217 N, 22.4159 E) and Slankamen (45.1415 N, 20.2586 E). Flies of both species were collected with a sweeping net using fermented apple traps. They were maintained in the laboratory at 19 ± 0.5 ◦C, 12 h:12 h light-dark cycle in 240 mL bottles containing 40 mL of the control substrate (standard-St) or the Pb-acetate-saturated

substrate (labeled as C3). The control substrate consisted of standard molasses corn meal diet (14 g agar, 208 g corn meal, 188 g sugar, 40 g dry active yeast, 5 g Nipagin diluted in 60 mL of 96% ethanol in 2.2 L distilled water) and the Pb-acetate-saturated substrate contained 1000 µg/mL of lead acetate. The standard substrate flies were grown originally for 45 and 30 generations (*D. subobscura* and *D. melanogaster*, respectively). After that, they were maintained for another 35 generations on the standard substrate and the Pb-acetatesaturated substrate. Both species were maintained at 19 ◦C; *D. melanogaster* is successfully reared at 19–25 ◦C, while for *D. subobscura* 19 ◦C is the optimal temperature. Since *D. subobscura* generally does not lay eggs in laboratory conditions without dry yeast powder, the powder was added to the substrate surface in both species in order to maintain equal conditions (detailed description is given in Beribaka et al. [5]). For the same purposes, the Nipagin, which was reported to modify the microbiota associated with the flies to some extent, was added in all substrates in the same amount at the same time [19].

#### *2.2. Experimental Setup*

On the first day of the experiment, 30 pairs of flies per bottle were transferred to a fresh substrate for 5 days; they were then removed from the bottles. After eclosion started, virgin males and females of *D. melanogaster* were collected every 12 h and once every 24 h for *D. subobscura* and separated by sex. Bottles with eclosing *D. subobscura* were additionally held in the dark to avoid possible mating. When 150 virgin males and females from each group were collected, 30 females and 30 males per bottle (5 bottles per group) within each group were placed to mate for 3 days. After that, the bottles were covered with Petri dishes that contained substrate and yeast dissolved in distilled water on the substrate surface and turned upside-down to enable flies to lay eggs. The eggs were collected every 8 h from the Petri dishes and transferred into 50 mL vials with 15 mL of the substrate and 1-2 drops of dissolved yeast added to the surface. For each group, 30 vials with 30 eggs were established (additional 5 vials per group were added for larvae collection). The F1 generation eclosion was recorded every 24 h until there were no new individuals in the vial for a period of 72 h; the flies were then stored in EtOH for NGS (separated by sex) at −20 ◦C. The egg-to-adult viability was calculated per vial as the percentage of individuals that emerged from the 30 eggs. The egg-to-adult developmental time was calculated using the formula:

$$DT = \frac{\Sigma n\_d \* d}{\Sigma n\_d} \tag{1}$$

where *n<sup>d</sup>* is the number of flies emerging in *d* days after the eggs were laid.

The third instar larvae were collected from the additional 5 vials per group and after 64 h of incubation they were rinsed with distilled water and stored in EtOH at −20 ◦C for NGS.

The samples were labeled as follows: *D. melanogaster*—Dmel, *D. subobscura*—Dsub, Kalna population—K, Slankamen population—Sl, standard (control) substrate—St, Pbsaturated substrate—C3, males—M, females—F and larvae—L.

#### *2.3. Statistical Analysis of Life History Traits*

Prior to the analyses of life history traits, all data were tested for normality by the Shapiro-Wilk test. Since the data for both egg-to-adult viability and egg-to-adult developmental time were normally distributed, we used the three- and four-factor ANOVA test to ascertain statistical significances. For egg-to-adult developmental time, the analyses were performed for males and females separately. The Bonferroni post-hoc test was used to identify the exact statistically significant differences. All the tests were performed in Statistica, ver. 10.0.228.2 [20].

#### *2.4. Total DNA Extraction and 16S rRNA Sequencing*

Total DNA was extracted from pools of 25 males, 25 females and 40 larvae for each experimental group, with 24 samples in total. DNA isolation was performed according to

the modified protocol by Kapun et al. [21]. The samples were homogenized using handheld motor homogenizer in a 1.5 mL Eppendorf tube in 300 µL of Solution A (1 M Tris HCl, pH 9, 0.5 M EDTA, pH 8 and 1% SDS). Then 4.5 µL of Proteinase K (20 mg/mL) were added, and the samples were vortexed and incubated for 30 min at 56 ◦C, after which they were incubated for another 30 min at 70 ◦C and, finally, for a few minutes at 37 ◦C. Afterwards, 3 µL of RNAse A (10 mg/mL) were added and incubated for 30 min at 37 ◦C. Next, 42 µL 8 M potassium acetate was added and the mixture was kept in the freezer for 30 min. The supernatant was collected after centrifugation at 13,000 rpm for 15 min. One volume of phenol-chloroform-isoamyl alcohol (25:24:1) was added and the mixture was centrifuged for 5 min at 13,000 rpm. After the supernatant was collected, the previous step was repeated with 0.75 volume of pure chloroform. Afterwards, 2.5 volume of 95% ice-cold ethanol was added and centrifuged at 10,000 rpm for 5 min to pellet the DNA. The pellet was washed with 1 mL of 70% ethanol, centrifuged at 13,000 rpm for 5 min and then the supernatant was discarded. The pellet was dried for 30 min and left to resuspend in 50 µL of TE buffer. DNA quality was assessed using the Qubit fluorometer (Thermo Fisher Scientific, Waltham, MA, USA).

The sequencing was performed by Eurofins Genomics Europe Sequencing GmbH (Konstanz, Germany) using the standard procedure InView—Microbiome Profiling 3.0 with MiSeq. The V3-V4 hypervariable region of the 16S rRNA gene was amplified using the forward primer 50 -TACGGGAGGCAGCAG-30 and reverse primer 50 -CCAGGGTATCTAATCC-3 0 [22,23]. The samples passed a quality check and fastq data were delivered for further processing. All sequence data were submitted to the GenBank (SRA) database under accession numbers SRX12591286-SRX12591309 (PRJNA616141 BioProject).

#### *2.5. Amplicon Sequence Variant Inference and Taxonomy Assignment*

The delivered sequences already had the adapters and linkers removed. The primers were removed using cutadapt 3.4 (https://cutadapt.readthedocs.io/en/stable/guide.html, accessed on 2 December 2021) in paired-end mode and reads without detected primers were discarded. Afterwards, dada2 [24] was used to filter and trim the reads. Apart from the default filter and trim options, the forward reads were trimmed to 270 bases while the reverse reads were trimmed to 200 bases and all the reads shorter than 150 bases were removed, as well as all those which had more than two estimated errors for the forward or reverse reads. After error estimation and dereplication, denoising was performed in selfConsist mode (the algorithm alternated between sample inference and error rate estimation until convergence). Sequence pair merger was performed with a minimal overlap of 20 bases without mismatches. Chimeric sequences were removed using default options in remove-BimeraDenovo (https://rdrr.io/bioc/dada2/man/removeBimeraDenovo.html, accessed on 2 December 2021). After the inspection of amplicon sequence variants (ASV) length distribution, all ASV with a length between 400 and 428 were kept. Taxonomy assignment up to the genus level was performed with IDTAXA algorithm [25] using default parameters and the SILVA v138 database (https://www.arb-silva.de/documentation/release-138/, accessed on 2 December 2021). Species level classification was performed with a SILVA species assignment train set (https://zenodo.org/record/3731176, accessed on 2 December 2021) using exact sequence matching without mismatches; multiple matches were allowed (multiple species output). One ASV which had low read counts was not classified at the Kingdom level and was removed from further processing.

#### *2.6. Microbiome and Statistical Analysis*

Sequence diversity within samples (alpha diversity) was estimated using the phyloseq 1.36.0 R package [26] at the ASV level after rarefaction to even depth (sample with the lowest number of reads) and shown through estimators Shannon, Gini-Simpson and invSimpson indices. Observed and estimated richness was determined according to the number of observations (Observed) and the Chao1 index. Comparison of alpha diversity among *D. melanogaster* and *D. subobscura* samples, as well as the impact of the growth

substrate on Chao1 and Shannon indices within each species, was determined using the Wilcoxon Rank-Sum test. For further analysis, prevalence filtering was performed.

For beta diversity and differential abundance estimation between samples, genus level aggregation was used (reads from ASV classified as the same genus were aggregated). Since many low occurrence genera were present in only one or a few samples, prevalence filtering was performed where all genera present in less than four samples were removed. Beta diversity was estimated using Double Principal Coordinate Analysis (DPCoA) [27] of the prevalence filtered data after rarefaction to even depth. Since DPCoA relies also on phylogenetic distances of sequences apart from abundances, a maxim likelihood (GTR+G(4)+I) phylogenetic model was estimated using a multiple sequence alignment of the microbial 16S sequences constructed to take into account RNA secondary structures. To assess if the microbial composition differs among species and among fixed effects (sex and substrate) within species, permutational multivariate analysis of variance (PERMANOVA) using distance matrices via R package vegan 2.5-7 [28] was performed on the DPCoA distance matrices.

To perform differential abundance estimation, a Wilcoxon Rank-Sum test was used to test for differences between the median relative abundances of each genus in each *Drosophila* species. *p*-value adjustment for multiple comparisons was performed according to the Benjamini-Hochberg method [29]. In addition to this, similar tests were performed to check if sex, substrate or population affect abundance of taxa within each species; however, in all cases no taxa were found significant with the nonparametric test. The R package metacoder 0.3.5 [30] was used to create a heat tree which visualizes taxonomic categories significantly differing among groups. The analyses were performed using the R (version 4.1) [31].

#### **3. Results**

#### *3.1. Life History Traits Analysis*

We investigated the influence of the composition of microbiota, different feeding substrates and population origin on two *Drosophila* life history traits, egg-to-adult viability and developmental time. Mean values revealed that egg-to-adult viability was the highest in Dmel\_K\_St (0.88) and the lowest was in Dmel\_Sl\_C3 (0.51), while developmental time was the highest in Dsub\_Sl\_C3\_M (23.06) and the lowest was in Dmel\_K\_St\_F (19.24) (Table 1).

The analysis of different effects showed that population origin, substrate, species and sex were significant for differences in both life history traits, but when it comes to the interactions between these factors, not all interactions were significant (Table 2).

The post-hoc analysis revealed that egg-to-adult viability did not differ significantly between populations in *D. subobscura* on both the standard and the lead-saturated substrate. There were also no significant differences between the substrates within either population for *D. subobscura*. *D. melanogaster* showed no significant differences between the populations on the lead-saturated substrate (Figure 1A). All other interactions showed statistically significant differences (*p* < 0.05) in comparisons.

The results obtained for developmental time included sex as another potential factor that could be differentially affected regarding the population, species and substrate. The post-hoc analysis revealed that, on the standard substrate, there were no significant differences between males and females for both species and both populations. *D. subobscura* showed this pattern on the lead-saturated substrate for both populations as well, while in *D. melanogaster* sex differences were significant on the lead-saturated substrate for both populations. In *D. melanogaster*, both sexes from the same substrate but different populations did not differ significantly (Figure 1B). In *D. subobscura*, this pattern was observed only for lead-saturated substrate (Figure 1C). All other interactions showed statistically significant differences (*p* < 0.05) in comparisons.


**Table 1.** Mean values of egg-to-adult viability and developmental time in all samples.

Dmel—*Drosophila melanogaster*; Dsub—*Drosophila subobscura*; K—Kalna population; Sl—Slankamen population; St—standard (control) substrate; C3—lead-saturated substrate; M—male; F—female.

**Table 2.** Results of the factorial ANOVA test on (a) egg-to-adult viability and (b) developmental time for different effects.


The modified table was presented as a poster presentation at IECE2021 conference.

#### *3.2. Taxonomic Abundance of Microbial Communities* The results obtained for developmental time included sex as another potential factor

statistically significant differences (*p* < 0.05) in comparisons.

*Insects* **2021**, *12*, x FOR PEER REVIEW 7 of 19

A total of 2,706,363 reads in 24 samples were obtained after quality filtering and these were used for ASV inference (Supplementary Table S1). All but one ASV were classified up to the order level (100 ASVs), 81 ASVs were classified up to the genus level and 49 to the species level (Supplementary Tables S2–S4). One ASV could not be assigned to Bacteria, and since it was present in traces only in one sample it was removed from further analysis. that could be differentially affected regarding the population, species and substrate. The post-hoc analysis revealed that, on the standard substrate, there were no significant differences between males and females for both species and both populations. *D. subobscura* showed this pattern on the lead-saturated substrate for both populations as well, while in *D. melanogaster* sex differences were significant on the lead-saturated substrate for both

The post-hoc analysis revealed that egg-to-adult viability did not differ significantly between populations in *D. subobscura* on both the standard and the lead-saturated substrate. There were also no significant differences between the substrates within either population for *D. subobscura*. *D. melanogaster* showed no significant differences between the populations on the lead-saturated substrate (Figure 1A). All other interactions showed

The microbiota of two *Drosophila* species, *D. melanogaster* and *D. subobscura,* was classified into six phyla. At the phylum level, the most abundant phyla in total were Proteobacteria and Firmicutes with 86.6% and 13.3%, respectively; other four phyla accounted for less than 1% in total (Figure 2A). populations. In *D. melanogaster*, both sexes from the same substrate but different populations did not differ significantly (Figure 1B). In *D. subobscura*, this pattern was observed only for lead-saturated substrate (Figure 1C). All other interactions showed statistically significant differences (*p* < 0.05) in comparisons.

**Figure 1.** *Cont*.

**Figure 1.** The LS mean plots of the analyzed traits depict changes that have occurred in different groups and subgroups; (**A**) egg-to-adult viability; (**B**) developmental time in *D. melanogaster;* (**C**) developmental time in *D. subobscura*. The numbers above the graphs indicate the *p*-values of the comparisons. The modified figure (**A**) was presented as a poster presentation at IECE2021 conference. **Figure 1.** The LS mean plots of the analyzed traits depict changes that have occurred in different groups and subgroups; (**A**) egg-to-adult viability; (**B**) developmental time in *D. melanogaster;* (**C**) developmental time in *D. subobscura*. The numbers above the graphs indicate the *p*-values of the comparisons. The modified figure (**A**) was presented as a poster presentation at IECE2021 conference.

*3.2. Taxonomic Abundance of Microbial Communities*  A total of 2,706,363 reads in 24 samples were obtained after quality filtering and these were used for ASV inference (Supplementary Table S1). All but one ASV were classified up to the order level (100 ASVs), 81 ASVs were classified up to the genus level and 49 to the species level (Supplementary Tables S2–S4). One ASV could not be assigned to Bacteria, and since it was present in traces only in one sample it was removed from further analysis. The microbiota of two *Drosophila* species, *D. melanogaster* and *D. subobscura,* was classified into six phyla. At the phylum level, the most abundant phyla in total were Proteobacteria and Firmicutes with 86.6% and 13.3%, respectively; other four phyla accounted for less than 1% in total (Figure 2A). Phylum Proteobacteria was dominant in almost all samples (23 of 24). On the other hand, Firmicutes was highly represented only in *D. melanogaster*, while its prevalence in *D. subobscura* was less than 2% out of total Firmicutes. Phylum Firmicutes was 99.9% represented by the genus *Lactobacillus*. Other most prevalent genera in all the samples were *Wolbachia, Komagataeibacter* and *Acetobacter* (Figure 2B). *Wolbachia* (36.9%) was present in all the samples of *D*. *melanogaster* flies, ranging from 2% to 21.6%, while it was present only in traces in two *D*. *subobscura* samples (<0.01%). The second most represented genus was *Komagataeibacter* (25%), with 95.9% found in *D. subobscura* samples. The increase in abundance of *Komagataeibacter* on the lead-saturated substrate was observed in five of six comparisons. Another highly represented genus from the Acetobacteraceae family was *Acetobacter*, with 65.5% in *D. subobscura* and 33.5% in *D. melanogaster* samples. Within each group, *Acetobacter* was the most prevalent in larvae samples. Another genus represented by >10% in total was *Lactobacillus* (13.3%), with 98.2% found in *D. melanogaster* samples, with higher prevalence in adults. All other genera (such as *Vibrionimonas, Staphylococcus, Sphingomonas* and *Acinetobacter*) were abundant with less than 1% in total. In *D. subobscura* samples, rare genera such as *Staphylococcus, Sphingomonas, Vibrionimonas, Acinetobacter* and a genus from family Xanthobacteraceae were present only in larvae samples. Interestingly, in *D. melanogaster* samples, *Staphylococcus, Sphingomonas* and *Vibrionimonas* were completely absent from larvae samples, but also from all lead-saturated substrate samples.

#### *3.3. Alpha Diversity Analysis*

Alpha diversity was measured using several metrics: Observed, Chao1, Shannon, Gini-Simpson and invSimpson indices (Table 3).

In terms of alpha diversity, the Shannon index did not show statistically significant differences among groups of samples, but the Chao1 estimator was significantly lower in *D. subobscura* samples (Figure 3). When control versus lead substrate samples were compared within each of the species, no significant differences in Chao1 were observed in any of the samples included. However, if larvae samples were excluded, the Wilcoxon Rank-Sum test indicated a significantly lower Chao1 index in lead substrate samples for both species (Figure 3).

**Figure 2.** Relative abundance of the most prevalent (**A**) phyla and (**B**) genera in *Drosophila* species (*D. melanogaster* and *D. subobscura*) from two populations (Kalna and Slankamen) on the control substrate and the lead-saturated substrate in larvae and adult males and females. The modified figure (**B**) was presented as a poster presentation at IECE2021 conference. **Figure 2.** Relative abundance of the most prevalent (**A**) phyla and (**B**) genera in *Drosophila* species (*D. melanogaster* and *D. subobscura*) from two populations (Kalna and Slankamen) on the control substrate and the lead-saturated substrate in larvae and adult males and females. The modified figure (**B**) was presented as a poster presentation at IECE2021 conference.

Phylum Proteobacteria was dominant in almost all samples (23 of 24). On the other hand, Firmicutes was highly represented only in *D. melanogaster*, while its prevalence in *D. subobscura* was less than 2% out of total Firmicutes. Phylum Firmicutes was 99.9% represented by the genus *Lactobacillus*. Other most prevalent genera in all the samples were *Wolbachia, Komagataeibacter* and *Acetobacter* (Figure 2B). *Wolbachia* (36.9%) was present in all the samples of *D*. *melanogaster* flies, ranging from 2% to 21.6%, while it was present only in traces in two *D*. *subobscura* samples (<0.01%). The second most represented genus was *Komagataeibacter* (25%), with 95.9% found in *D. subobscura* samples. The increase in According to the Shannon and Simpson diversity indices, the highest diversity in adults was observed within the *D. melanogaster* female samples from the Kalna (standard) and from Slankamen (standard and lead-saturated substrate). Regarding the larvae samples, the highest diversity was observed in *D. melanogaster*, the Kalna population on lead and both *D. subobscura* samples on lead (Kalna and Slankamen). The lowest diversity in adults was found in both sexes and from both substrates in *D. subobscura* from Kalna population. Larvae samples with the lowest diversity were *D. melanogaster* from Slankamen on lead and *D. subobscura* from the Kalna population on the standard.

abundance of *Komagataeibacter* on the lead-saturated substrate was observed in five of six comparisons. Another highly represented genus from the Acetobacteraceae family was *Acetobacter*, with 65.5% in *D. subobscura* and 33.5% in *D. melanogaster* samples. Within each group, *Acetobacter* was the most prevalent in larvae samples. Another genus represented


**Table 3.** Alpha diversity of ASVs represented in the microbial community of *Drosophila*.

Dmel—*Drosophila melanogaster*; Dsub—*Drosophila subobscura*; K—Kalna population; Sl—Slankamen population; St—standard (control) substrate; C3—lead-saturated substrate; M—male; F—female; L—larvae.

**Figure 3.** Differences in the Chao1 index between the species and the treatments estimated with the Wilcoxon Rank-Sum **Figure 3.** Differences in the Chao1 index between the species and the treatments estimated with the Wilcoxon Rank-Sum test. The number above the data indicates the *p*-value of the comparison.

test. The number above the data indicates the *p*-value of the comparison.

larvae samples.

(Figure 4).

*3.4. Beta Diversity Analysis* 

According to the Shannon and Simpson diversity indices, the highest diversity in adults was observed within the *D. melanogaster* female samples from the Kalna (standard) and from Slankamen (standard and lead-saturated substrate). Regarding the larvae samples, the highest diversity was observed in *D. melanogaster*, the Kalna population on lead and both *D. subobscura* samples on lead (Kalna and Slankamen). The lowest diversity in Similar to the Shannon and Simpson diversity indices, observed richness and Chao1 estimators were highest in adult *D. melanogaster* samples on the standard in the Kalna population, both sexes and the Slankamen population on the standard in males. The lowest richness in adults was observed in seven out of eight *D. subobscura* samples. Larvae richness was higher in *D. subobscura* on lead (both population) than in all *D. melanogaster* larvae samples.

adults was found in both sexes and from both substrates in *D. subobscura* from Kalna population. Larvae samples with the lowest diversity were *D. melanogaster* from Slankamen

population, both sexes and the Slankamen population on the standard in males. The lowest richness in adults was observed in seven out of eight *D. subobscura* samples. Larvae richness was higher in *D. subobscura* on lead (both population) than in all *D. melanogaster* 

Similar to the Shannon and Simpson diversity indices, observed richness and Chao1

Beta diversity was estimated using Double Principal Coordinate Analysis (DPCoA) on prevalence filtered taxa at the genus level aggregation after rarefication to even depth

on lead and *D. subobscura* from the Kalna population on the standard.

#### *3.4. Beta Diversity Analysis*

Beta diversity was estimated using Double Principal Coordinate Analysis (DPCoA) on prevalence filtered taxa at the genus level aggregation after rarefication to even depth (Figure 4). *Insects* **2021**, *12*, x FOR PEER REVIEW 12 of 19

**Figure 4.** Biplot of Double Principal Coordinate Analysis (DPCoA) of bacterial composition in different *Drosophila* samples (*D. melanogaster* and *D. subobscura*; control and lead-saturated substrate indicated in the legend) using prevalence filtered taxa at the genus level aggregation. f—female; m—male; l—larvae; S—Slankamen population; K—Kalna population. **Figure 4.** Biplot of Double Principal Coordinate Analysis (DPCoA) of bacterial composition in different *Drosophila* samples (*D. melanogaster* and *D. subobscura*; control and lead-saturated substrate indicated in the legend) using prevalence filtered taxa at the genus level aggregation. f—female; m—male; l—larvae; S—Slankamen population; K—Kalna population.

Samples belonging to the two *Drosophila* species were separated by the first DPCoA axis. All *D. subobscura* samples clustered together, while *D. melanogaster* samples were dispersed on the biplot. Partial separation between the control substrate and the lead-saturated substrate of *D. subobscura* samples was observed by the second DPCoA axis. In addition, the larvae from both control *D. melanogaster* populations were positioned close to *D. subobscura* samples. Samples belonging to the two *Drosophila* species were separated by the first DPCoA axis. All *D. subobscura* samples clustered together, while *D. melanogaster* samples were dispersed on the biplot. Partial separation between the control substrate and the leadsaturated substrate of *D. subobscura* samples was observed by the second DPCoA axis. In addition, the larvae from both control *D. melanogaster* populations were positioned close to *D. subobscura* samples.

To investigate the source of variation in microbial composition, or more precisely to partition it among fixed effects such as species, sex, substrate and their interaction, PER-MANOVA using DPCoA distances among the samples was performed on all samples and separately for *D. melanogaster* and *D. subobscura* sample subsets. When DPCoA distances among all samples are taken into account, the greatest source of variation was due to species and all interactions were significant (Table 4). When *D. melanogaster* samples were analyzed separately, substrate, sex and sex × substrate interaction were statistically significant (*p* < 0.05), indicating that the substrate had a different effect depending on the sex. For the *D. subobscura* samples DPCoA distance matrix, the substrate, sex and sex × substrate interaction did not have a significant effect. To investigate the source of variation in microbial composition, or more precisely to partition it among fixed effects such as species, sex, substrate and their interaction, PERMANOVA using DPCoA distances among the samples was performed on all samples and separately for *D. melanogaster* and *D. subobscura* sample subsets. When DPCoA distances among all samples are taken into account, the greatest source of variation was due to species and all interactions were significant (Table 4). When *D. melanogaster* samples were analyzed separately, substrate, sex and sex × substrate interaction were statistically significant (*p* < 0.05), indicating that the substrate had a different effect depending on the sex. For the *D. subobscura* samples DPCoA distance matrix, the substrate, sex and sex × substrate interaction did not have a significant effect.

**Table 4.** Results of PERMANOVA using DPCoA distances among all *Drosophila* samples. **Effect df SS MS F R2** *p* Species 1 0.893082 0.893082 132.2773 0.651498 0.001 Differential abundances of the most prevalent genera indicate the dominance of the Acetobacteraceae family in *D. subobscura* samples, whereas the Lactobacillaceae family was predominant in *D. melanogaster* (Figure 5). Substrate and sex had no effect if all samples were taken into account, nor when *D. subobscura* and *D. melanogaster* were analyzed separately.

> Substrate 1 0.046581 0.046581 6.899272 0.033981 0.02 Sex 2 0.100009 0.050004 7.406304 0.072956 0.008

Species × Substrate 1 0.034041 0.034041 5.041911 0.024833 0.029 Species × Sex 2 0.08629 0.043145 6.390341 0.062948 0.013 Substrate × Sex 2 0.065261 0.032631 4.833024 0.047608 0.014 Species × Substrate × Sex 2 0.06453 0.032265 4.778869 0.047074 0.02 Residuals 12 0.081019 0.006752 0.059103


**Table 4.** Results of PERMANOVA using DPCoA distances among all *Drosophila* samples.

**Figure 5.** Differential heat tree of microbial communities' abundances in *Drosophila* species with adjusted *p*-values < 0.05. Orange-colored taxa are more abundant in *D. melanogaster* while blue-colored taxa are more abundant in *D. subobscura*. The size of the circle indicates the number of samples with taxonomic category. **Figure 5.** Differential heat tree of microbial communities' abundances in *Drosophila* species with adjusted *p*-values < 0.05. Orange-colored taxa are more abundant in *D. melanogaster* while bluecolored taxa are more abundant in *D. subobscura*. The size of the circle indicates the number of samples with taxonomic category.

#### **4. Discussion**

**4. Discussion**  In this study we observed that the overall bacterial diversity and richness were higher in adult *D. melanogaster* compared to *D. subobscura* samples. Life history analysis showed significant differences in population origin, substrate, species and sex effects on egg-toadult viability and developmental time. The most prevalent genera in all the samples were *Wolbachia, Komagataeibacter*, *Acetobacter* and *Lactobacillus*. The genus *Lactobacillus* was dominantly abundant in *D. melanogaster* species on the standard substrate, while *Komagataeibacter* genus was dominant in *D. subobscura* on lead-saturated substrate. *Komagataeibacter* genus proved to be a species-specific member of *D. subobscura* microbiota that could In this study we observed that the overall bacterial diversity and richness were higher in adult *D. melanogaster* compared to *D. subobscura* samples. Life history analysis showed significant differences in population origin, substrate, species and sex effects on egg-toadult viability and developmental time. The most prevalent genera in all the samples were *Wolbachia, Komagataeibacter*, *Acetobacter* and *Lactobacillus*. The genus *Lactobacillus* was dominantly abundant in *D. melanogaster* species on the standard substrate, while *Komagataeibacter* genus was dominant in *D. subobscura* on lead-saturated substrate. *Komagataeibacter* genus proved to be a species-specific member of *D. subobscura* microbiota that could be beneficial in overcoming environmental stress.

be beneficial in overcoming environmental stress. Comparing the results obtained from the composition of microbiota and life history traits of two *Drosophila* species reared on different substrates, several potential cause-and-Comparing the results obtained from the composition of microbiota and life history traits of two *Drosophila* species reared on different substrates, several potential cause-and-

effect relationships were discovered. The overall microbial diversity and evenness in *D.* 

study which was done only on the Kalna population [5]. With the addition of the Slankamen population in this research it can be seen that the Kalna population in *D. melanogaster*

effect relationships were discovered. The overall microbial diversity and evenness in *D. subobscura* reared in the laboratory was lower than in *D. melanogaster*, as in our previous study which was done only on the Kalna population [5]. With the addition of the Slankamen population in this research it can be seen that the Kalna population in *D. melanogaster* was mainly advantageous regarding the microbial diversity and richness compared to the Slankamen population. Contrary to that, *D. subobscura* originating from the Slankamen population showed higher diversity than the Kalna population. A similar pattern was observed for both adults and larvae, which indicates that the change of microbiota diversity due to the lead exposure could be population-specific. Microbial richness estimated by the Chao1 estimator was lower in *D. subobscura* in overall analysis, but also in both species on the substrate saturated with lead, when only males and females were included. Beta diversity showed that the differences in bacterial diversity were most expressed on the species level, but also revealed that microbial composition of *D. melanogaster* larvae was more similar to *D. subobscura* adult samples than to *D. melanogaster*. *D. subobscura* samples showed strong clustering and indicated a minor impact of lead exposure to variation in microbial composition. This could be due to the increase in the *Komagataeibacter* genus within the lead-saturated samples, which has been reported as a good probiotic candidate due to its high level of glucose conversion rate and survival rate in the presence of acidic pH and bile salt [32].

Life history traits also varied significantly between the populations and kept a similar pattern as microbiota in *D. melanogaster* on the standard substrate, whereby egg-to-adult viability was higher in the Kalna population compared to Slankamen. *D. subobscura* did not show significant differences in egg-to-adult viability regarding the origin, nor regarding the substrate composition, as was shown by Tanaskovi´c et al. [17]. On the other hand, developmental time results revealed that *D. subobscura* underwent more changes in this trait regarding the microbial composition, population origin and substrate than *D. melanogaster*. Namely, developmental time varied significantly between *D. subobscura* males and females on the standard substrate in both populations, but also within each sex on the standard substrate with different population origin. Males took longer to develop on the standard substrate compared to females, but both males and females originating from the Slankamen population took longer to develop than the Kalna population. This could be due to preadaptation to the polluted environment, as Kalna is probably more polluted area compared to Slankamen. *D. melanogaster* did not show significant differences in developmental time regarding the origin, nor regarding the substrate composition. This indicates that in some species, specific traits could be more susceptible to lead toxicity and changes in microbiota than others.

Microbial composition analysis indicated the dominance of the Acetobacteraceae family in *D. subobscura* samples and the Lactobacillaceae family in *D. melanogaster* species. Thus, *D. melanogaster*, the sample with the highest diversity and richness in microbiota species and the highest representation of genus *Lactobacillus*, showed the highest eggto-adult viability and the shortest developmental time (Dmel\_K\_St, females), while the sample which showed microbial diversity albeit poor richness and domination by the genus *Acetobacter* exhibited lower egg-to-adult viability (Dmel\_Sl\_C3, females). Additionally, it was observed that in larvae of this species, the genus *Acetobacter* was dominant on the standard substrate, while its presence was significantly lower on the substrate with lead. *Lactobacillus* and *Acetobacter* are commonly found in lab-reared *D. melanogaster* [4,33,34]. Both genera can promote growth via different pathways, and in certain conditions, the presence of *Lactobacillus plantarum* helps larval growth and reduces their developmental time, so this could be a reason for short developmental time in Dmel\_K\_St [35–37]. The presence of 98.2% of total *Lactobacillus* ASVs in *D. melanogaster* suggested that *Lactobacillus* could be a species-specific member of *D. melanogaster* gut microbiota. The *Acetobacter* genus was present in larvae samples of both species, with the highest abundance on the standard substrate, but there was a decrease in egg-to-adult viability when *Acetobacter* was accumulated by adults (Dmel\_Sl\_C3 and Dsub\_K\_St). Previous studies have shown

that lead toxicity can drive oxidative stress in many organisms [38–42]. Oxidative stress in *D. melanogaster* can cause various health defects, including reduced lifespan, retarded development, decreased pupation, emergence and survival rates, impaired mobility and reduced egg production [38,43–45]. In addition to the fact that lead has been proven to affect various life history traits, highly reduced viability in *D. melanogaster* on the lead-saturated substrate compared to *D. subobscura* could be a cost of faster development [46,47] and presence of endosymbiotic bacteria (*Wolbachia*), which also have been confirmed to affect the gut microbiota [48,49]. Short developmental time tends to have various fitness costs besides the reduced egg-to-adult viability, such as lower pathogen resistance [50], borderline larval storage of metabolites and reduced adult size [47]. Prolonged developmental time could be a potential mechanism of resistance to heavy metal exposure, providing a higher egg-toadult viability.

Another possible factor that is greatly involved in shaping the microbiota is temperature [51]. *D. melanogaster* is successfully reared at 19–25 ◦C, unlike *D. subobscura,* which has a much tighter temperature range in the wild; 19 ◦C is the optimal rearing temperature in the lab. Heat-stressed *D. subobscura* flies showed changes in bacterial diversity and structure compared to non-stressed flies, and this response demonstrates that the gut microbiota contributes to heat tolerance, which could have important consequences on host fitness [52]. The sub-optimal rearing temperature for *D. melanogaster* could affect the metabolic strategy during the development, but also the growth of the species-specific microbiota. Additional experimental temperature manipulation would probably give a more complete answer in that sense.

Although the *D. subobscura* showed lower richness in microbiota species and lower diversity (Kalna population, both sexes and both substrates), the samples that were predominantly represented by the genus *Komagataeibacter* mainly maintained similar levels of egg-to-adult viability. *D. subobscura* samples from standard substrate, where *Komagataeibacter* genus was highly abundant (Slankamen population), also showed an increase in egg-toadult viability compared to the population with low prevalence of *Komagataeibacter*. This indicates that the high prevalence of the genus *Komagataeibacter* was beneficial for flies' viability in lab-rearing conditions, but also that it could be the key to the higher tolerance to lead exposure in *D. subobscura*. We previously reported the increase of *Komagataeibacter* in lab-reared flies after 13 generations, where its abundance drastically increased on the lead-saturated substrate, pointing to its higher tolerance to this heavy metal if compared to the other members of the microbial community [5]. After 35 generations, its prevalence has been maintained on the lead-saturated substrate, indicating a good heavy metal adaptation of *D. subobscura* species [15]. Measuring of the concentration of lead by inductively coupled plasma optical emission spectrometry (ICP-OES) in *D. melanogaster* and *D. subobscura* flies maintained for more than 30 generations in the control and lead-saturated substrate conditions showed that *D. subobscura* flies on lead-saturated substrate accumulated more lead than *D. melanogaster* (unpublished data). Moreover, the amount of lead accumulated was higher in males than in females in *D. subobscura*, whereas in *D. melanogaster* it was the opposite. The resistance of *D. subobscura* to increased accumulation of lead could be due to the prevalence of the *Komagataeibacter* genus, which could be an example of stable gut-colonizing bacteria in *D. subobscura*, since it has been proven to have a strong anti-oxidant ability in vitro [53] and it is considered to be a good probiotic candidate [32].

Taxonomic analysis revealed that the rare genera (<1%) in *D. subobscura* were present only in larvae, whereby *Staphylococcus* and *Acinetobacter* were present on the lead-saturated substrate, and *Sphingomonas, Vibrionimonas* and a genus from the family Xanthobacteraceae were present in larvae from both substrates. Interestingly, in *D. melanogaster* most of them were almost completely absent in larvae from both substrates, but also in the majority of adults from the lead-saturated samples. These findings suggest a different dynamic of developmental stages, as well as variability in substrate utilization and degradation by larvae in two species. These implicate the modulation in adaptive strategies under different environmental conditions in the two species.

#### **5. Conclusions**

In this study, we observed different patterns of life history traits in accordance with population origin and sex, but also the dominance of different gut microbiota members. The population origin showed a significant influence on life history traits, though each of the traits in the two species was affected differentially. Sex differences were also expressed, but only in *D. subobscura* on the standard substrate, indicating that influence of population origin and sex on life history traits could be species-specific. The presence of the heavy metal caused shifts in developmental time in *D. subobscura*, but maintained the egg-to-adult viability at a similar level. This could be explained by the domination of the *Komagataeibacter* in *D. subobscura* gut microbiota, usually a rare member of the microbiota community. The egg-to-adult viability increased in *D. subobscura* on standard substrate when *Komagataeibacter* was highly abundant, indicating that it could be a valuable member of *D. subobscura* microbiota in overcoming environmental stress. Research of the impact of microbiota on the adaptive response to heavy metals and the potential implications on host fitness is of great importance. Further research could reveal the extent to which species, sex, origin, lead exposure and specific members of microbiota, individually or through interactions, affect the life history traits. It could also help to identify the exact members of the gut microbiota that enable the best possible response to a particular environmental change.

**Supplementary Materials:** The following are available online at http://www.mdpi.com/article/ 10.3390/insects12121122/s1, Table S1: Number of sequences before and after filtering, Table S2: Bacterial phylum taxonomy contingency table of microbiota in all samples, Table S3: Bacterial genus taxonomy contingency table of microbiota in all samples, Table S4: Bacterial species taxonomy contingency table of microbiota in all samples.

**Author Contributions:** Conceptualization, M.S.-R.; methodology M.B. and M.S.-R.; software, M.B.; validation, M.B., M.J., M.T. and M.S.-R.; formal analysis, M.B., M.J. and M.T.; investigation, M.B., M.J. and C.L.; resources, M.B., M.T. and M.S.-R.; data curation, M.B., M.J. and C.L.; writing—original draft preparation, M.B.; writing—review and editing, M.J., M.T. and M.S.-R.; visualization, M.B.; supervision, M.S.-R.; project administration, M.S.-R.; funding acquisition, M.B., M.J., M.T. and M.S.-R. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by Ministry for Scientific and Technological Development, Higher Education and Information Society of the Republic of Srpska, B&H grant number 19/6-020/961- 104/18 for M.B. and C.L.; Ministry of Education, Science and Technological Development of the Republic of Serbia grant number 451-03-9/2021-14/200178 for M.S.-R. and M.J.; Ministry of Education, Science and Technological Development of the Republic of Serbia grant number 451-03-9/2021-14/ 200007 for M.T.

**Institutional Review Board Statement:** Not applicable.

**Data Availability Statement:** The data presented in this study are available in Supplementary Materials.

**Acknowledgments:** We are grateful to Milan Dragiˇcevi´c for his help in bioinformatic analysis of data and Vesna Cvijetinovi´c for proofreading the English language. The research is also supported through European Drosophila Population Genomic Consortium (DrosEU), which is funded by a Special Topic 587 Networks (STN) grant from the European Society for Evolutionary Biology. We are grateful to the reviewers for comments which improved the manuscript. Part of the data from this paper was presented as a poster presentation at 1st International Electronic Conference on Entomology (IECE2021): Beribaka, M.; Jeli´c, M.; Lazi´c, C.; Stamenkovi´c-Radak, M. Microbiota Composition Affects Life History Traits in Drosophila Species, in Proceedings of the 1st International Electronic Conference on Entomology, 1–15 July 2021, MDPI: Basel, Switzerland, doi:10.3390/IECE-10511.

**Conflicts of Interest:** The authors declare no conflict of interest.
