Characterization of Gut Microbiota of Honey Bees in Korea

Shuvo, Md Sarower Hossen; Kim, Sukyung; Jo, Sujin; Rahim, Md Abdur; Barman, Indrajeet; Hossain, Mohammed Solayman; Jeong, Yoonkyoung; Jeong, Hwasik; Kim, Sangrim; Seo, Hoonhee; Song, Ho-Yeon

doi:10.33073/pjm-2025-025

Figures & Tables

Baseline characteristics of the honey bee population studied.
A)Queen bees, worker bees, and drone bees were collected from the Chungnam region of Korea; B) live bees were collected in boxes; C) samples were stored at –80°C; D) body length and weight were manually measured prior to gut extraction. Distinct morphological differences were evident among the castes: drone bees had the most considerable body length and weight, worker bees were the smallest, and queen bees were intermediate in size. Significance among groups was analyzed using one-way ANOVA in GraphPad Prism v9.1.1 (GraphPad Software, USA, www.graphpad.com) (# – p < 0.0001). — Baseline characteristics of the honey bee population studied. A)Queen bees, worker bees, and drone bees were collected from the Chungnam region of Korea; B) live bees were collected in boxes; C) samples were stored at –80°C; D) body length and weight were manually measured prior to gut extraction. Distinct morphological differences were evident among the castes: drone bees had the most considerable body length and weight, worker bees were the smallest, and queen bees were intermediate in size. Significance among groups was analyzed using one-way ANOVA in GraphPad Prism v9.1.1 (GraphPad Software, USA, www.graphpad.com) (# – p < 0.0001).

Mean taxonomic composition of the honey bee gut microbiome.
The taxonomic relative abundances across bee castes were classified at multiple levels: A) phylum, B) class, C) order, D) family, E) genus. Statistical significance among groups (*p < 0.05; **p < 0.01; ***p < 0.001) was analyzed using the Wilcoxon rank-sum test, and relative abundances below 1% are represented as “ETC”. Overall, drone and queen bees exhibited similar levels of taxonomic diversity, with no significant taxa identified between them. In contrast, worker bees showed greater taxonomic diversity compared to the other two castes, with several significant taxa identified. Significant differences compared to worker bees are indicated as “a” for queen bees and “b” for drone bees. — Mean taxonomic composition of the honey bee gut microbiome. The taxonomic relative abundances across bee castes were classified at multiple levels: A) phylum, B) class, C) order, D) family, E) genus. Statistical significance among groups (*p < 0.05; **p < 0.01; ***p < 0.001) was analyzed using the Wilcoxon rank-sum test, and relative abundances below 1% are represented as “ETC”. Overall, drone and queen bees exhibited similar levels of taxonomic diversity, with no significant taxa identified between them. In contrast, worker bees showed greater taxonomic diversity compared to the other two castes, with several significant taxa identified. Significant differences compared to worker bees are indicated as “a” for queen bees and “b” for drone bees.

Box plots of species richness indices among honey bee groups.
Species richness was analyzed using A) Ace, B) Chao1, C) Jackknife indices, D) OTUs (Operational Taxonomic Units) counts. The bars represent the median values, while the hinges indicate the lower and upper quartiles. No statistically significant differences were observed among the three bee castes. — Box plots of species richness indices among honey bee groups. Species richness was analyzed using A) Ace, B) Chao1, C) Jackknife indices, D) OTUs (Operational Taxonomic Units) counts. The bars represent the median values, while the hinges indicate the lower and upper quartiles. No statistically significant differences were observed among the three bee castes.

Box plots of species diversity indices among honey bee groups.
Species diversity was analyzed using A) NPShannon, B) Shannon, C) Simpson, and D) Faith’s phylogenetic diversity indices. The results indicate that worker bees exhibited significantly higher species diversity compared to queen and drone bees in the NPShannon, Shannon, and Simpson analyses. No differences were observed between queen and drone bees. The bars represent the median values, while the hinges indicate the lower and upper quartiles. — Box plots of species diversity indices among honey bee groups. Species diversity was analyzed using A) NPShannon, B) Shannon, C) Simpson, and D) Faith’s phylogenetic diversity indices. The results indicate that worker bees exhibited significantly higher species diversity compared to queen and drone bees in the NPShannon, Shannon, and Simpson analyses. No differences were observed between queen and drone bees. The bars represent the median values, while the hinges indicate the lower and upper quartiles.

Principal coordinate analysis (PCoA) of the gut microbiome across bee castes.
PCoA was conducted using four distance metrics to investigate the microbial community structure among bee castes: A) Bray–Curtis, B) Jensen–Shannon divergence, C) Generalized UniFrac, D) UniFrac. The results are visualized in 3D plots, with blue, green, and red circles representing the worker, queen, and drone bee groups, respectively. The scatter of points illustrates the dissimilarities in microbial community composition among the castes. The analysis reveals that the microbial communities of worker bees are distinct from those of queen and drone bees. — Principal coordinate analysis (PCoA) of the gut microbiome across bee castes. PCoA was conducted using four distance metrics to investigate the microbial community structure among bee castes: A) Bray–Curtis, B) Jensen–Shannon divergence, C) Generalized UniFrac, D) UniFrac. The results are visualized in 3D plots, with blue, green, and red circles representing the worker, queen, and drone bee groups, respectively. The scatter of points illustrates the dissimilarities in microbial community composition among the castes. The analysis reveals that the microbial communities of worker bees are distinct from those of queen and drone bees.

UPGMA clustering analysis of gut microbiome profiles among honey bee castes.
The UPGMA (Unweighted Pair Group Method with Arithmetic Mean) clustering method was used to analyze the gut microbiome profiles of the three honey bee castes: worker, queen, and drone. The dendrogram illustrates the hierarchical grouping of these castes based on the similarity of their gut microbiome compositions. The results highlight the distinctiveness of workers and the relative similarity and closer grouping of queens and drones. — UPGMA clustering analysis of gut microbiome profiles among honey bee castes. The UPGMA (Unweighted Pair Group Method with Arithmetic Mean) clustering method was used to analyze the gut microbiome profiles of the three honey bee castes: worker, queen, and drone. The dendrogram illustrates the hierarchical grouping of these castes based on the similarity of their gut microbiome compositions. The results highlight the distinctiveness of workers and the relative similarity and closer grouping of queens and drones.

Taxonomic biomarker analysis of gut microbiome differences among bee castes using LEfSe (LDA Effect Size).
Based on LEfSe analysis, taxonomic cladograms highlight significant biomarkers among bee castes with a Kruskal–Wallis p-value < 0.05 and an LDA score > 4. The analysis identified distinct microbial features for each caste: Lactobacillus as a biomarker for drones, Gilliamella and Frischella as biomarkers for workers, and Bombella as a biomarker for queens. The taxa presented represent key microbial markers that differentiate the bee groups, providing insights into caste-specific microbial community compositions. — Taxonomic biomarker analysis of gut microbiome differences among bee castes using LEfSe (LDA Effect Size). Based on LEfSe analysis, taxonomic cladograms highlight significant biomarkers among bee castes with a Kruskal–Wallis p-value < 0.05 and an LDA score > 4. The analysis identified distinct microbial features for each caste: Lactobacillus as a biomarker for drones, Gilliamella and Frischella as biomarkers for workers, and Bombella as a biomarker for queens. The taxa presented represent key microbial markers that differentiate the bee groups, providing insights into caste-specific microbial community compositions.

PCR validation of biomarker strains across the three honey bee castes.
The bacterial load of A) Lactobacillus, B) Gilliamella, C) Bombella, D) Frischella in the samples was quantified using qRT-PCR. The results confirmed the biomarker status of each strain: Lactobacillus had the lowest Ct values in drones, Gilliamella and Frischella in workers, and Bombella in queens. Statistical significance among groups was determined using one-way ANOVA (**p < 0.01; ***p < 0.001). — PCR validation of biomarker strains across the three honey bee castes. The bacterial load of A) Lactobacillus, B) Gilliamella, C) Bombella, D) Frischella in the samples was quantified using qRT-PCR. The results confirmed the biomarker status of each strain: Lactobacillus had the lowest Ct values in drones, Gilliamella and Frischella in workers, and Bombella in queens. Statistical significance among groups was determined using one-way ANOVA (**p < 0.01; ***p < 0.001).

Summary of the characterization of caste-specific gut microbiota in honey bees from Korea.
This figure summarizes the study investigating the gut microbiota of honey bees collected in Korea, focusing on worker, queen, and drone castes. The top section outlines a standardized protocol for honey bee microbiome analysis, which includes bee collection, gut dissection, DNA extraction, metagenomic library preparation, and next-generation sequencing (NGS). The bottom section highlights the microbiota composition of the three castes: worker bees exhibit a diverse composition with Lactobacillus (68.3%), Gilliamella (23.0%), and Frischella (4.7%); queen bees are dominated by Lactobacillus (95.4%) with Bombella (4.0%); and drone bees are almost entirely composed of Lactobacillus (98.6%). These results highlight the distinct microbial communities associated with each caste and provide insights into their physiological roles and colony health. — Summary of the characterization of caste-specific gut microbiota in honey bees from Korea. This figure summarizes the study investigating the gut microbiota of honey bees collected in Korea, focusing on worker, queen, and drone castes. The top section outlines a standardized protocol for honey bee microbiome analysis, which includes bee collection, gut dissection, DNA extraction, metagenomic library preparation, and next-generation sequencing (NGS). The bottom section highlights the microbiota composition of the three castes: worker bees exhibit a diverse composition with Lactobacillus (68.3%), Gilliamella (23.0%), and Frischella (4.7%); queen bees are dominated by Lactobacillus (95.4%) with Bombella (4.0%); and drone bees are almost entirely composed of Lactobacillus (98.6%). These results highlight the distinct microbial communities associated with each caste and provide insights into their physiological roles and colony health.

Average taxonomic composition of significantly different taxa among worker, queen, and drone bees_

Phylum

W (%)

Q (%)

D (%)

Class

W (%)

Q (%)

D (%)

Order

W (%)

Q (%)

D (%)

Family

W (%)

Q (%)

D (%)

Genus

W (%)

Q (%)

D (%)

Firmicutes

69.0

95.4 a*

98.6

Bacilli

69.0

95.4 a*

98.6

Lactobacillales

69.0

95.4 a*

98.6

Lactobacillaceae

68.3

95.4 a*

98.6

Lactobacillus

68.3

95.4 a*

98.6

Proteobacteria

30.4

4.3 a***

1.0 b**

Gammaproteobacteria

30.3

0.3 a***

0.3 b**

Orbales

27.7

0.3 a***

0.2 b**

Orbaceae

27.7

0.3 a***

0.2 b**

Gilliamella

23.0

0.2 a***

0.2 b**

Frischella

4.7

0.1 a*

0.0 b*

Alphaproteobacteria

0.1

4.0 a**

0.8

Rhodospirillales

0.1

4.0 a**

0.8

Acetobacteraceae

0.1

4.0 a**

0.8

Bombella

0.1

4.0 a**

0.5

Beta diversity index among the honey bee castes_

Index	Worker-Queen	Worker-Drone	Drone-Queen
Bray–Curtis	NS (p = 0.078)	NS (p = 0.301)	NS (p = 0.313)
Jensen–Shannon	NS (p = 0.438)	NS (p = 0.581)	NS (p = 0.644)
Generalized UniFrac	* (p = 0.031)	* (p = 0.050)	NS (p = 0.495)
UniFrac	** (p = 0.003)	** (p = 0.007)	NS (p = 0.835)

Prediction of functional biomarkers for the three honey bee castes_

Ortholog		Definition	LDA effect size	p-value	Worker (%)	Queen (%)	Drone (%)
Orthology	K11904	Type VI secretion system secreted protein VgrG	3.1	0.000	0.2	0.0	0.0
	K15580	Oligopeptide transport system substrate-binding protein	2.9	0.007	0.1	0.2	0.2
	K07052	Uncharacterized protein	2.9	0.000	0.2	0.3	0.3
	K01439	Succinyl-diaminopimelate desuccinylase	2.9	0.042	0.1	0.2	0.2
	K10947	PadR family transcriptional regulator, regulatory protein PadR	2.9	0.002	0.1	0.3	0.2
Module (PICRUSt)	K02945	Ribosome, bacteria	3.7	0.042	5.7	6.8	6.6
	K11904	Type VI secretion system	3.6	0.000	0.8	0.0	0.0
	K02035, K02031, K02032	Peptides/nickel transport system	3.4	0.000	0.9	1.3	1.3
	K00798, K02226	Cobalamin biosynthesis, cobinamide => cobalamin	3.3	0.060	0.4	0.8	0.7
	K15580	Oligopeptide transport system	3.3	0.002	0.3	0.7	0.6
Module (MinPath)	K17204	Erythritol transport system	3.2	0.000	0.3	0.5	0.6
	K02052	Putative spermidine/putrescine transport system	3.2	0.000	0.3	0.6	0.6
	K20491, K20492	Lantibiotic transport system	3.1	0.000	0.2	0.5	0.5
	K11195	PTS system, fructose-specific II component	3.1	0.000	0.3	0.6	0.6
	K12536	Hemophore/metalloprotease transport system	3.1	0.000	0.2	0.4	0.4
Pathway (PICRUSt)	K00010, K00011, K00015, K00016	Metabolic pathways	3.7	0.001	15.9	16.7	16.7
	K01995, K02031, K02032, K02035	Quorum sensing	3.3	0.000	1.6	2.1	2.1
		Flagellar assembly	3.3	0.000	0.4	0.1	0.0
	K02945	Ribosome	3.2	0.006	2.0	2.4	2.3
	K00010, K00015, K00016, K00033	Microbial metabolism in diverse environments	3.2	0.000	4.8	5.1	5.1
Pathway (MinPath)	K00466	Tryptophan metabolism	3.6	0.000	0.8	1.5	1.6
	K00965, K01198, K02790, K02791	Amino sugar and nucleotide sugar metabolism	3.5	0.000	1.0	1.5	1.6
	K01259	Arginine and proline metabolism	3.5	0.000	1.1	1.8	1.8
	K00011, K00965, K01785, K02744	Galactose metabolism	3.5	0.000	0.7	0.1	0.1
	K13993	Protein processing in endoplasmic reticulum	3.4	0.000	1.2	1.7	1.7

Characterization of Gut Microbiota of Honey Bees in Korea

Figures & Tables

Fig. 1.

Fig. 2.

Fig. 3.

Fig. 4.

Fig. 5.

Fig. 6.

Fig. 7.

Fig. 8.

Fig. 9.

Average taxonomic composition of significantly different taxa among worker, queen, and drone bees_

Beta diversity index among the honey bee castes_

Prediction of functional biomarkers for the three honey bee castes_

Paradigm

My account