A. cerana from Changbai Mountain in northeast China is adapted to cold temperatures and mountain conditions. Colonies can naturally overwinter in extreme cold conditions at −40°C, and the overwintering period lasts up to six months in Changbai Mountain. It is a cold-resistant honeybee species with economic value (Liu et al., 2022a; Liu et al., 2022b). The honeybee gut microbiota serves as a health indicator in the beekeeping industry (Wu et al., 2020). Thus, studying changes in the structure and composition of the gut microbiota during winter allows understanding of the cold environmental adaptation mechanisms in A. cerana.
As a native species in China, A. cerana along with its gut microbial diversity and function has adapted to the environment through evolution (Zhang et al., 2019). This study analyzed the gut microbial structure and diversity of A. cerana from Changbai Mountain during winter by 16S rDNA high-throughput sequencing.
MATERIALS AND METHODS
Sampling
We selected three healthy overwintering colonies with a similar population in the Jilin Province (128°12′48″E, 41°26′34″N), China. All colonies were checked before the experiments to ensure that no overwintering feed was stored in the colonies. The feed station was filled every three days until the colonies had enough food reserves for overwintering. The insects were fed a mixture of honey and sugar. Five adult worker bees were randomly selected from each hive. The whole intestine of each bee was extracted and merged together. Intestinal samples were immersed in liquid nitrogen and stored at −80°C. They were collected from October 2022 to March 2023 on the same day and time each month.
DNA extraction and 16S rDNA sequencing
Total genomic DNA was extracted using the FastDNA® SPIN Kit for Soil (MP Biomedicals, Santa Ana, CA) according to the manufacturer's instructions. The integrity of genomic DNA was assessed by agarose gel electrophoresis, and the concentration and purity of genomic DNA were determined using a Nanodrop 2000 spectrophotometer and a Qubit 3.0 fluorometer. Amplicon sequencing was performed by Genesky Biotechnologies Inc. (Shanghai, China). The V3–V4 hypervariable regions of the 16S rDNA gene were amplified using the primers 341F (5′-CCTACGGGNGGCWGCAG-3′) and 805R (5′-GACTACHVGGGTATCTAATCC-3′) and were sequenced on an Illumina NovaSeq 6000 platform.
Illumina data processing and analysis
Raw reads were processed with the use of the QIIME 2 software package and associated plugins (Bolyen et al., 2019). Quality control and the identification of amplicon sequence variants (ASVs) were performed with the DADA2 plugin (Callahan et al., 2016). Taxonomic assignments of representative ASV sequences were performed by a pretrained naive Bayes classifier of the Ribosomal Database Project version 11.5 at a confidence threshold of 0.8. Differences in relative abundance at each classification level (p<0.05) were determined by one-way analysis of variance (ANOVA). Differentially abundant operational taxonomic units between groups were identified using the linear discriminant analysis (LDA) effect size (LEfSe) method at a threshold of p<0.05 and log10 LDA scores >2. The community richness indices - observed species, Chao1 and abundance-based coverage estimator [ACE] and the community diversity indices - Shannon, Simpson and coverage were calculated using the R package vegan version 2.5.6. Principal component analysis (PCA) and principal coordinate analysis (PCoA) were performed with the use of R packages ade4 version 1.7.13 and scatter plot3d version 0.3.41, respectively. Gene function prediction was performed with the PICRUSt2 version 2.3.0 in three databases: MetaCyc, Kyoto Encyclopedia of Genes and Genomes (KEGG) and Clusters of Orthologous Groups (COG) (Kanehisa et al., 2012; Langille et al., 2013; Caspi et al., 2016).
RESULTS
Sequencing data and taxonomic classification
Eighteen A. cerana intestinal samples were sequenced, and 1,567,459 high-quality sequences were obtained. The percentage of input passed filter, input merged, and input non-chimeric was 86.7%, 85.3%, and 76.9%, respectively (Tab. S1). The length of most sequences varied from 423 to 428 bp. A total of 313 ASVs were obtained. The number of exclusive ASVs from October 2022 to March 2023 was 35, 47, 36, 40, 36, and 45, respectively, and the number of shared ASVs was 28 (Fig. 1). The most abundant genus shared among samples was Lactobacillus.
Fig. 1.
Distribution of amplicon sequence variants (ASVs) in samples collected from October 2022 (YD_10) to March 2023 (YD_03). Each color represents a group, the number of exclusive ASVs in each month is shown in ellipses, and the number of shared ASVs is shown in circles.
Analysis of microbial community composition
The most abundant groups at different taxonomic levels were Proteobacteria, Firmicutes, Actinobacteria, and Bacteroidetes (phylum); Gammaproteobacteria, Bacilli, Actinobacteria, Alphaproteobacteria, Betaproteobacteria, and Flavobacterium (class); Orbales, Lactobacillales, Bifidobacteriales, Neisseriales, Flavobacteriales, and Rhizobiales (order); Orbaceae, Lactobacillaceae, Bifidobacteriaceae, Neisseriaceae, and Flavobacteriaceae (family); Gilliamella, Lactobacillus, Bifidobacterium, and Snodgrassella (genus); Bifidobacterium asteroides (species) (online Fig. S1).
The ANOVA results showed that several groups were differentially abundant at different taxonomic levels: Firmicutes (phylum) (p<0.05); Bacilli and Betaproteobacteria (class) (p<0.0001, p<0.05); Enterobacteriales, Neisseriales, Burkholderiales, Rhodospirillales, and Lactobacillales (order) (p<0.0001, p<0.0001, p<0.0001, p<0.01, p<0.05); Enterobacteriaceae, Neisseriaceae, Burkholderiaceae, Prevotellaceae, and Acetobacteraceae (family) (all at p<0.0001); Enterobacter, Snodgrassella, Morganella, and Ralstonia (genus) (all at p<0.0001); and Lactobacillus kunkeei (species) (p<0.05). The relative abundance of Lactobacillales fluctuated during the study period and was lowest in January. The relative abundance of Betaproteobacteria, Neisseriales, Neisseriaceae, Snodgrassella, Rhodospirillales and Acetobacteraceae peaked in December. The relative abundance of other groups was low (less than 0.1%) (online Fig. S2). The LEfSe analysis results showed the following:
the greater abundance of Enterobacteriaceae, Enteroacteriales, Burkholderiales and Morganella in November;
the predominance of Neisseriaceae, Beta-proteobacteria, Snodgrassella, Acetobacteraceae, Neisseriales, and Rhodospirillales in December;
Ruminoccaceae more common in January;
the prevelance of Rhizobales, Alpha-proteobacteria and Prevotelaceae in February;
the greater abundance of Ralstonia in March (Fig. 2).
Fig. 2.
Linear discriminant analysis (LDA) effect size analysis of differentially abundant (log10 LDA scores >2) bacterial taxa.
Analysis of community diversity
Alpha diversity: There were no significant differences in the α diversity indexes of observed species, Chao1, ACE and coverage during the overwintering period (p>0.05). In turn, community diversity peaked in December, represented by the highest Shannon index and lowest Simpson index during this month (p<0.05) (Fig. 3).
Fig. 3.
Alpha diversity indexes in samples collected from October 2022 (YD_10) to March 2023 (YD_03).
Beta diversity: PCA and PCoA analyses showed that samples collected in December, January and February clustered together, samples collected in November and March clustered together, while samples collected in October clustered separately (Fig. 4). These results indicate that community composition did not vary significantly from December 2022 to February 2023.
Fig. 4.
Principal component analysis (PCA) and principal coordinate analysis (PCoA). Each point is a sample, and each color represents the samples collected from October 2022 (YD_10) to March 2023 (YD_03). The lines indicate relative distances.
Function prediction
The COG database analysis revealed that microbial genes were enriched in translation, ribosomal structure and biogenesis, amino acid transport and metabolism, carbohydrate transport and metabolism, cell wall/membrane/envelope biogenesis, and transcription (Fig. 5). KEGG analysis indicated that ko00473 (D-alanine metabolism) and ko00471 (D-glutamine and D-glutamate metabolism) were highly enriched in our samples (Fig. 6A). KEGG Orthology (KO) results showed that K01992 (ABC-2 type transport system permease protein [ABC-2.P]), K02004 (putative ABC transport system permease protein [ABC.CD.P]) and K01990 (ABC-2 type transport system ATP-binding protein [ABC-2.A]) were significantly enriched (Fig. 6B). MetaCyc results indicated that PWY-7208 (super pathway of pyrimidine nucleobases salvage) and PWY–7221 (guanosine ribonucleotides de novo biosynthesis) were highly enriched in the samples (Fig. 6C).
Fig. 5.
Functional classification of proteins in a sample collected in December 2022 based on the Clusters of Orthologous Groups database.
Fig. 6.
Relative abundance of Kyoto Encyclopedia of Genes and genomes (KEGG), KEGG Orthology, and MetaCyc pathways. The bars indicate relative abundances.
A. cerana from Changbai Mountain is a cold-tolerant honeybee species. The coldest months in this region are December and January (Liu et al., 2022b), and therefore we focused on the microbial pathways enriched within this period. The ANOVA analysis of KEGG pathways showed that ko00630 (glyoxylate and dicarboxylate metabolism), ko00190 (oxidative phosphorylation), ko00020 (TCA cycle), ko00906 (carotenoid biosynthesis), and ko00860 (porphyrin and chlorophyll metabolism) were significantly enriched in these months. The analysis of KO pathways revealed that K07657 (phoB; two-component system, OmpR family, phosphate regulon response regulator PhoB) and K07132 (MuB; ATP-dependent target DNA activator [EC:3.6.1.3]) were enriched. The analysis of MetaCyc pathways revealed that TCA (TCA cycle I), SULFATE-CYS-PWY (super pathway of sulfate assimilation and cysteine biosynthesis), DTDPRHAMSYN-PWY (dTDP-L-rhamnose biosynthesis I) and REDCITCYC (TCA cycle VIII) were highly enriched in this period (online Fig. S3).
DISCUSSION
Gut microbial communities play an essential role in honeybee physiology, nutrition and health (Hamdi et al., 2011). A. cerana intestinal samples collected from October 2022 to March 2023 were analyzed through 16S rDNA sequencing. The results showed that the most abundant taxa were Proteobacteria, Firmicutes, Actinobacteria and Bacteroidetes. The most predominant genera were Gilliamella, Lactobacillus, Bifidobacterium and Snodgrassella, consistent with the results of previous studies (Guo et al., 2015). The most abundant genus across samples was Lactobacillus, which indicates that it is ubiquitous in the gut of honeybees during winter. This genus protects honeybees from pathogens and improves nectar and pollen utilization (Sabaté et al., 2009; Emery et al., 2017).
The relative abundance of Snodgrassella and Acetobacteraceae peaked in December. Snodgrassella colonizes the hindgut and affects microbial metabolism by reducing pH and redox potential (Martinson et al., 2012; Zheng et al., 2017). Glycolysis is the main metabolic activity in the honeybee gut and its final products of glycolysis vary depending on the type and strain but typically include Lactobacillus and Acetobacter (Kwong et al., 2016; Zheng et al., 2017). Gut bacteria have overlapping metabolic capabilities, and several gut bacteria, Acetobacter and Lactobacillus, metabolize the same food substrates (Bonilla-Rosso et al., 2018). Other non-core low-abundance microbial taxa play an indispensable role by interacting with other gut microorganisms or through their direct or potential impact on the host. For instance, Lactobacillus kunkeei is adapted to acidic and sugar-rich environments, including the insect crop and food stores, and helps preserve bee bread and increase colonization resistance by decreasing pH (Cotter & Hill, 2003; Anderson et al., 2013).
The overwintering period on Changbai Mountain lasts approximately six months from late October to early March, a critical period for the survival of A. cerana colonies. The α-diversity of the gut microbiota of overwintering honeybees increased in early winter from October to December, decreased to the lowest level in January, and increased in late winter from February to March. During overwintering, honeybees consume bee bread because they are unable to obtain fresh pollen from the environment (Liu et al., 2022b). Moreover, honeybees are highly vulnerable to pathogens and multiple stresses. The higher gut community diversity in the early stages of overwintering improves honeybee immunity and adaptation to dietary changes and cold stress. Dietary and metabolic factors affect the gut microbiota of honeybees, and sustained low temperatures and food consumption during mid stages of overwintering decrease gut microbial community diversity. Some gut bacteria, Bifidobacterium and Lactobacillus Firm-5, contain genes that encode enzymes involved in the biosynthesis and utilization of trehalose, a disaccharide used by insects to store energy in the form of nectar (Bottacini et al., 2012; Kwong et al., 2014; Ellegaard et al., 2015). In addition, Gilliamella, Lactobacillus and Bifidobacterium synthesize pectin-degrading enzymes, glycoside hydrolase, polysaccharide lyases and other enzymes that degrade carbohydrates (Engel et al., 2012). In the late stage of overwintering, when the environmental temperature gradually recover, A. cerana starts to collect food, and hypopharyngeal glands' ability to secrete royal jelly is restored. The midgut wall of A. cerana is significantly thicker than that of western honeybees in winter, which indicates that A. cerana maintains midgut development during cold seasons (Liu et al., 2017). These physiological activities can not be separated from the participation of gut microbiota, which may lead to an increase in community diversity.
The analysis of β-diversity showed that gut community composition was similar from December 2022 to February 2023. According to LEfSe analysis, Snodgrassella, Acetobacteraceae and Rhizobiales were significantly enriched in the coldest months from December 2022 to February 2023, which suggests that these taxa are implicated in cold stress tolerance. The greatest abundance of Acetobacteraceae, Rhizobiales and Lactobacillus in the honeybee gut during winter from November to February (Ludvigsen et al., 2015) is in line with our results. Sucrose-rich diets increase the relative abundance of Acetobacteraceae (Taylor et al., 2019). The bee food used in our experiments is a mixture of high-quality granulated sugar and honey, whose main component is sucrose. The increased relative abundance of Rhizobiales may help prevent feces decay in the rectum during winter, when honeybees cannot leave the nest to defecate (D'Alvise et al., 2018). In addition, Neisseriaceae, represented by Snodgrassella, protects the honeybee gut from opportunistic fecal pathogens during winter until they can be excreted in the spring (Bleau et al., 2020).
The functional prediction results suggest that the transport and metabolism of amino acids, carbohydrates and other energy sources and the TCA cycle improve microbial response to cold stress and dietary changes during winter. Liu et al. (2022) studied the full-length transcriptome of A. cerana from Changbai Mountain in winter using Oxford nanopore sequencing. In the functional classification of the COG database, differentially expressed genes and transcripts were significantly enriched in carbohydrate transport and metabolism (Liu et al., 2022b). Li et al. (2022) used genomics, metagenomics and metabolomics to assess changes in the diversity and function of gut bacteria across seasons and indicated that gut bacteria might benefit the host by providing essential amino acids during the period of pollen shortage. In addition, amino acid biosynthesis, the TCA cycle and dicarboxylic acid metabolism were enriched in late winter (Li et al., 2022), which is consistent with our results.
Amplicon data of A. cerana from Changbai Mountain during overwintering were obtained through 16S rDNA sequencing. Our results indicate how the gut microbial community of A. cerana from Changbai Mountain responds to cold stress and provide a basis for further exploration of the diversity of the gut microbiota of A. cerana. In addition, honeybees are a good model for studying the relationship between gut microbiota and longevity. The genotypes of overwintering honeybee, worker honeybee and the queen are identical, but the changes of gut microbiota are closely related to different periods, and the underlying mechanisms need to be further studied. These basic theoretical issues are the directions that we need to consider in future research on the gut microbiota of A. cerana.
, , , , oraz | 22 gru 2023
