1.1. Biome characterization
The study was performed in two different biomes in Paraíba State, Northeastern Brazil: Atlantic Forest and Caatinga. These biomes are remarkedly different in terms of soil composition, climate characteristics, and vegetation. The Atlantic Forest, a type of altitude tropical forest, extends from Rio Grande do Norte to Rio Grande do Sul States, i.e., along the East coast of th ecountry. It is econpasses a diversity of forest ecosystems with markedly differentiated floristic structures and compositions shaped by the specific climatic characteristics of each region [21]. On the other hand, Caatinga, located in the Northeast of Brazil, is characterized by a semi-arid climate, vegetation adapted to aridity, and a remarkable diversity of species [22]. About 7.31% of Paraiba territory is occupied by the Atlantic Forest, while 92.69% by the Caatinga [23] (IBGE and IBGE, 2019).
1.2. Study design and samplings
Five samples, with an average of 35 Africanized Apis mellifera bee nurses, were collected from five different hives in apiaries located in two contrasting biomes in Brazil. Five samples belonged to the Beekeeping Sector of the Federal University of Paraíba) (Latitude: 6°58'19.77"S, Longitude: 35°43'16.88"W, Altitude 522 m) in the region of Atlantic Forest and the other five samples originated from nested hives in the apiary of São João do Cariri, belonging to the Experimental Station of UFPB) (7° 12' 43'' S, 36° 37' 46'' W; altitude 510 m above sea level) in the Caatinga region (Fig. 1).
The samplings were carried out in the dry period of the year, by the end of November, which is characterized by increased temperatures and absence of rainfall in both regions (Atlantic Forest and Caatinga). PPE equipment necessary for the management of Apis mellifera species bees were used. The bees were placed into sterile tubes containing 70% ethanol and stored at -20ºC until dissection.
1.3. Sample preparation and DNA extraction
The bee intestines were dissected in a sterile environment. Before extraction, the bee samples were left on a sterile filter paper for 10 minutes for thawing and ethanol evaporation. In the dissection process, a cross-section was made in the last segment of the bee abdomen, and then the contents were transferred to sterile microtubes. Each sample represented pooled gut specimens from 35 individual bees. DNA was extracted with a commercial kit (PowerSoil, Qiagen, Hilden, Germany), according to the manufacturer's instructions. DNA integrity was evaluated in agarose gel and DNA concentration was measured in a spectrophotometer (Colibri LB 915, Titertek-Berthold, Germany).
1.4. Library preparation and sequencing
The libraries for 16S rRNA metabarconding were prepared as previously described using a reference protocol (Illumina, 2013). Shortly, the V3-V4 hypervarible region of the 16S rRNA microbial gene was amplified by PCR. The 25 µL master mix consisted of 2.5 µL of DNA template (5 ng/µL), 5 µL of each primer (341F: 5′-TCG TCG GCA GCG TCA GAT GTG TAT AAG AGA CAG CCT ACG GGN GGC WGC A – 3′ and 805R: 5′-GTC TCG TGG GCT CGG AGA TGT GTA TAA GAG ACA GGA CTA CHV GGG TAT CTA ATC C – 3′), and 12 µL 2X KAPA HiFi HotStart ReadyMix (KAPA Biosystems, Wilmington, MA, USA). The thermal cycling condition consisted of an initial denaturation at 95°C for 3 min, followed by 25 cycles at 95°C for 30 s, 55°C for 30 s, and 72°C for 30 s, and a final extension at 72°C for 5 min.
Amplification products were visualized in 1.5% agarose gel before purification using magnetic beads (AMPure XP, Beckman Coulter, United States) to remove excess primer. The dual indices and Illumina sequencing adapters were attached using a Nextera XT Index Kit (Illumina). A second clean up step was then performed using magnetic beads. The purified PCR products were quantified by fluorometry (Qubit 2.0, Life Invitrogen). For sequencing, pooled libraries were denatured with NaOH, diluted with hybridization buffer, then heat denatured. Paired-end sequencing was performed on an Illumina MiSeq with a V2 kit (2 × 250 cycles). At least 5% PhiX DNA was added for sequencing control purposes (Kit PhiX, Illumina).
Paire-end sequencing was performed on a MiSeq platform (Illumina) using a 2x250 cycles V2 kit (Illumina). At least 15% PhiX DNA has been added as a sequencing control (Kit PhiX, Illumina).
1.5. Bioinformatics and statistical analysis
The raw demultiplexed paired-end sequences were processed using QIIME 2-2020.2 [24]. Reads were filtered, denoised, and truncated to a length of 248 base pairs, and then parsed for non-chimeric sequences using DADA2, producing Amplicon Sequence Variants (ASV) [25]. Sequences were aligned using “qiime fragment-insertion sepp” for phylogenetic analysis [26]. Taxonomic composition of the samples were determined rypla pretrained naive Bayes classifier with a 99% sequence similarity threshold for V3-V4 reference sequences (SILVA-132-99-nb-classifier.qza) and the “qiime feature-classifier classify-sklearn.” Negative control samples were examined for potential contaminant taxa. No taxa overlapped between negative control and true samples.
Microbial diversity was quantified using Pielou’s (evenness) and Shannon (richness and abundance) diversity indices. ANOVAs were used to compare diversity between groups in [27] 4.1.0 [28] after testing for normality using a Shapiro-Wilk test.For alpha diversity assessment, we used the phylogenetic diversity indices (Shannon (abundance and richness), Faith (phylogenetic diversity), Pielou (Uniformity of abundance distribution) and Number of observed ASVs in QIIME 2-2020.2 [24]. Differences in alpha diversity between samples from the two biomes were assessed in R version 4.1.0 [27] after normality using Shapiro-Wilk tests; then, as appropriate, used ANOVAs followed by Tukey or Kruskal-Wallis Tests followed by Wilcoxon rank-sum test at 5% probability to compare microbial diversity between honey bee biome using package “laercio” Version 1.0–0.
For qualitative beta diversity analysis, we used Jaccard coefficient to evaluate the dissimilarity across samples regarding the gut microbial structure based on presence and absence of taxa, and Unweighted Unifrac as a phylogenetic-based metric for dissimilarity assessment. For quantitative beta diversity, we used Bray Curtis coefficient and Weighted Unifrac metrics. Bray Curtis was chosen because it combines both abundance and presence/absence data for the compositional dissimilarity assessment, while Weighted Unifrac is a phylogenetic-oriented metric for taxa. Beta-diversity analyses were performed in QIIME 2-2020.2 platform [24]. PERMANOVA was used for each metric as recommended to test putative differences in taxa composition between biomes [29].
The core taxa, i.e., taxa present in all samples, was determined using the "qiime feature-table core-features" command. Analysis of Composition of Microbiomes (ANCOM) was also applied to assess differences in rate abundances between biomes. The relative abundances of ASVs observed in core and ANCOM were assessed. Data were tested for normality using the Shapiro-Wilk test, and an analysis of variance (ANOVA) was applied to observe significant differences under a P value < 0.05. The plugin Emperor 2020.2.0 in QIIME 2-2020.2 [24] was used to visualize beta-diversity outputs. We used R version 4.1.0 [27] to generate alpha-diversity figures.