Seasonal gut microbiota dynamics and colony performance under stationary and short-distance migratory beekeeping in Apis mellifera ligustica

Honey bee gut microbiota play important roles in host nutrition, immune regulation, and health, and their composition can vary in response to environmental conditions, diet, seasonality, and beekeeping management. In this two-year field study, we investigated seasonal variation in the gut bacterial communities of Apis mellifera ligustica colonies maintained under stationary or short-distance migratory beekeeping management. Six colonies originating from the same experimental apiary were assigned to stationary or migratory management, with three migratory colonies seasonally transferred to Mount Taburno during summer and three stationary colonies remaining at Castelpoto. Bacterial communities were characterized by high-throughput sequencing of the V1–V3 region of the 16S rRNA gene. The gut microbiota was consistently dominated by characteristic honey bee-associated taxa, particularly Gilliamella and Snodgrassella, although their relative abundances varied substantially across seasons and years. Alpha-diversity analyses based on harmonized genus-level relative abundance profiles revealed significant seasonal variation in both years, whereas no significant overall effect of management system was detected. Similarly, PERMANOVA of Bray–Curtis dissimilarities revealed significant seasonal structuring of bacterial community composition in 2022 and 2023, while stationary versus migratory management did not significantly explain overall community variation. The persistence of dominant gut-associated taxa across sampling periods indicates compositional continuity of major bacterial community members despite temporal variation, but does not by itself demonstrate functional stability or microbial resilience. Colony-performance measurements showed significantly higher annual honey production in migratory than in stationary colonies in 2022, whereas no significant differences in the investigated performance parameters were detected between management systems in 2023. Thus, the productive response was not consistent between years. Although temporal variation in gut bacterial communities occurred alongside differences in colony performance, the present study does not establish a causal or functional relationship between microbiota composition and productivity. Overall, these findings identify seasonality as the predominant factor associated with variation in the gut bacterial communities of the investigated colonies and provide no evidence for a consistent overall effect of short-distance migratory management on bacterial community diversity or structure. Further studies involving larger numbers of independent colonies and apiaries, together with direct functional analyses of the microbiome and characterization of dietary and environmental factors, are required to clarify the consequences of migratory beekeeping for honey bee microbiota and colony performance.

In recent years, honey bee (Apis mellifera) colonies have experienced significant population declines, creating major challenges for beekeepers and threatening agricultural ecosystems that rely on bee pollination (1, 2). Colony mortality and morbidity result from multiple interacting stressors, with Varroa mites as a major driver, while pathogens, pesticide exposure, nutritional deficiencies, and management practices further exacerbate losses (3–5).

The gut microbiota of the honey bee represents a relatively simple but highly specialized community, primarily dominated by five core bacterial lineages: Gilliamella, Snodgrassella, Bombilactobacillus (formerly Lactobacillus Firm-4), Lactobacillus (formerly Firm-5), and Bifidobacterium (6). In addition to these core members, other genera such as Bartonella, Commensalibacter, and Frischella are frequently detected in the gut but are not considered part of the ‘core’ microbiome (7, 8). Both core and non-core bacterial groups play pivotal roles in maintaining bee health by contributing to nutrient metabolism, pathogen resistance, detoxification of dietary components, and regulation of immune activity (9–11). Newly emerged workers acquire their characteristic gut microbiota primarily through social interactions within the hive, including contact with nest-mates and trophallaxis (9, 12–14). Although the honey bee gut microbiota is characterized by a relatively conserved set of core bacterial lineages, differences among individuals may occur in the relative abundance of these taxa and in the presence of transient or opportunistic bacteria (7). The gut microbiota affects both individual and colony-level health, and disturbances leading to dysbiosis can reduce host fitness and resilience (15–18). Given the central role of the gut microbiota in honey bee health, elucidating its composition and the factors driving dysbiosis is essential (19, 20). The composition of the bee gut microbiota is significantly affected by seasonal changes, primarily driven by diet, temperature, pathogenic infection, and ecosystem features (8, 14, 21–31). These changes affect the health of honey bees by strengthening or weakening their immunity, making them more or less susceptible to disease (6, 23). In a recent study, Braglia et al. (32) found that different pollens do not significantly affect the core gut microbiome composition. However, the effect of different diets on the microbiome suggests an indirect impact on the immune system response by modulating and influencing the synthesis of certain immune-related peptides (32).

Migratory beekeeping, practiced worldwide, involves relocating colonies in synchrony with floral cycles to optimize honey yields and support commercial pollination services (33). In the United States, millions of beehives are transported to provide pollination for crops such as almonds, blueberries, and cranberries (34, 35), with similar practices observed in Europe (36–38), and Asia (39–42). While providing essential pollination services, migratory beekeeping can affect biodiversity and honey bee health, potentially increasing or decreasing pollinator diversity depending on agricultural intensity and the health of native populations (43, 44). Short-distance migratory beekeeping is practiced not only in Europe but also in Australia and China at national or regional levels to exploit successive flowering periods, enabling the production of monofloral honeys and supporting colony adaptation to increasingly challenging conditions associated with hotter and drier summers (36). Some studies indicate that migratory beekeeping may negatively affect the health of honey bees (45, 46), while also increasing susceptibility to parasites and infectious diseases (36, 47–52). Cestaro et al. (53) evaluated the presence and prevalence of Paenibacillus larvae, Nosema apis, Nosema ceranae and Varroa destructor in colonies from the east-central region of São Paulo state, Brazil, in two different seasons and under two different management systems (stationary vs. migratory). The results demonstrate that seasonality affects the health of A. mellifera colonies, but this does not depend on the type of management adopted (stationary or migratory) (53). Recent studies have suggested that migratory beekeeping may be associated with changes in honey bee gut microbial diversity and composition (54). Such changes may have biological relevance because the gut microbiota contributes to nutrient metabolism and host defense; however, changes in microbial diversity alone should not be interpreted as evidence of impaired microbiome function or colony health (6). Advances in high-throughput sequencing and bioinformatics have greatly enhanced our understanding of the honeybee microbiome and its role in health (55). However, knowledge remains limited regarding the effects of migratory beekeeping on the gut microbiota of A. mellifera. A recent study by Nishanthini and Kanagarajan (54) on the Indian honeybee (Apis cerana indica) highlighted differences between the microbiomes of bees from different migratory sites. However, A. cerana (Eastern honey bee) and A. mellifera (Western honey bee) exhibit significant phenotypic and productive differences, largely due to their evolutionary adaptation to different environments (56, 57).

A. mellifera ligustica is one of the most widely kept and popular subspecies of A. mellifera in global beekeeping, valued for its docility, high productivity, and adaptability to diverse climates (58–60).

In this study, using Next Generation Sequencing (NGS) technology, we compared seasonal variation in the gut microbiota of honey bee colonies managed under stationary and short-distance migratory beekeeping. To our knowledge, this is the first two-year field study to investigate seasonal gut microbiota dynamics in A. m. ligustica under these two management conditions. We also evaluated colony performance parameters under stationary and migratory management.

Samples collection

The study was conducted over two consecutive years (2022 and 2023) in the Campania region (southern Italy), specifically in the municipality of Castelpoto (41°08′42″ N; 14°42′14″ E; 280 m a.s.l.) and at Piano Melaino on Mount Taburno (41°06′44″ N; 14°35′12″ E; 1350 m a.s.l.), within the Taburno-Camposauro Regional Park. Six A. m. ligustica colonies were initially located in the same experimental apiary in Castelpoto. All colonies were headed by sister queens derived from the same genetic line and were selected to have a comparable initial level of colony development. The use of sister queens allowed partial standardization of the host genetic background across experimental groups, although sister queens are not genetically identical. The original experimental scheme confirms the common genetic lineage and comparable initial colony development.

The six hives were assigned to two management groups. In the stationary beekeeping group, three hives (nos. 3, 4, and 6) remained at the Castelpoto site throughout the experimental period. In the migratory beekeeping group, three hives (nos. 1, 2, and 5) were seasonally transferred from Castelpoto to Piano Melaino on Mount Taburno, where they remained during the summer period and had returned to Castelpoto by the autumn sampling. In 2022, the migratory colonies were transferred to Mount Taburno on 25 June, while hives 3, 4, and 6 remained at Castelpoto. The same experimental management scheme was repeated in 2023. Microbiota data from the three stationary colonies (hives 3, 4, and 6) were previously included in a longitudinal study focused specifically on seasonal gut microbiota dynamics within the single Castelpoto apiary (61). In the present study, these data were incorporated into a broader experimental design together with three short-distance migratory colonies (hives 1, 2, and 5), with the specific aim of comparing seasonal gut microbiota dynamics between stationary and migratory management systems. The present study additionally evaluates colony-performance parameters in relation to the two management conditions.

The two sites were separated by an approximate horizontal geodesic distance of 10.5 km and differed markedly in altitude and surrounding floral resources. During summer, the Castelpoto site was mainly characterized by plants belonging to the genera Helminthotheca, Trifolium, and Beta. In contrast, the higher-altitude Taburno site showed a greater abundance of pollen-producing plants belonging to the genera Echium, Trifolium, Castanea, and Rubus (62, 63).

All honey bee colonies received standardized nutritional supplementation to ensure consistent dietary conditions. During winter, each colony received a single dose of 2.5 kg candied paste composed of sucrose, glucose, and fructose (Fondabee, Beernem, Belgium). At the beginning of spring, colonies were additionally supplemented with 50% (w/v) sucrose syrup, administered weekly for four consecutive weeks.

Microbiological sampling was performed weekly during each season in both 2022 and 2023. At each weekly sampling, fifteen returning adult forager bees were collected at the hive entrance from each colony. Individual worker age within the forager category was not experimentally standardized. After collection, the bees were immediately transported to the reference laboratory under controlled temperature conditions using liquid nitrogen and stored at −80 °C until processing. Bees collected weekly from the same colony during the corresponding season were subsequently combined to generate a single Seasonal Sampling Unit (SSU). For each SSU, all bees were processed and dissected on the same day according to Iorizzo et al. (64). For each bee, the hindgut, comprising the ileum and rectum, was aseptically dissected. The hindguts belonging to bees from the same SSU were pooled in a sterile 1.5-mL centrifuge tube. The resulting SSUs were stored at −80 °C until shipment to the DNA sequencing service provider.

Four SSUs, corresponding to spring, summer, autumn, and winter, were generated for each colony in each year. Thus, 24 SSUs were obtained annually (3 colonies × 2 management systems × 4 seasons), resulting in a total of 48 SSUs over the two-year study period. For each combination of management system, season, and year, the biological replication therefore consisted of three independent colonies (n = 3).

DNA extraction, amplicon library preparation and sequencing

Total genomic DNA was extracted from each Seasonal Sampling Unit (SSU) using the Stool DNA Isolation Kit (Norgen Biotek Corp., Thorold, ON, Canada) according to the manufacturer’s instructions. The bacterial community was characterized by amplification of the V1–V3 hypervariable regions of the 16S rRNA gene using the primer pair 27F (5′-AGAGTTTGATCMTGGCTCAG-3′) and 534R (5′-ATTACCGCGGCTGCTGG-3′) (65, 66). PCR amplicons were purified, quantified, and normalized before library preparation. Equimolar amounts of purified amplicons were pooled to ensure comparable sequencing depth among samples. Sequencing was conducted by Eurofins Genomics Europe Sequencing GmbH (Ebersberg, Germany) on an Illumina MiSeq platform using a paired-end 2 × 300 bp sequencing strategy.

Bioinformatic analysis

Raw paired-end sequences passing the Illumina chastity filter were demultiplexed according to sample-specific indices. Primer sequences were trimmed with no mismatches allowed, and read pairs lacking primer identification were discarded. Forward and reverse reads were merged when sufficient overlap was present; when merging was not possible, only high-quality forward reads were retained. Additional quality filtering removed sequences outside the expected length range and sequences containing ambiguous bases prior to downstream microbiome profiling.

Bioinformatic processing of the quality-filtered sequences was performed using the QIIME pipeline (67), as implemented within the Eurofins analytical workflow. Chimeric sequences were identified and removed using the de novo UCHIME algorithm (68) implemented in VSEARCH (69).

The remaining high-quality sequences were partitioned into operational taxonomic units (OTUs) using Minimum Entropy Decomposition (MED) (70). Taxonomic assignment of OTU representative sequences was performed using DC-MEGABLAST against the NCBI nucleotide (nt) reference database (release 2020-02-03). Reference sequences were considered when they showed at least 70% sequence identity across at least 80% of the OTU representative sequence. The final taxonomic assignment corresponded to the lowest common taxonomic unit shared by the best-matching reference sequences. Bacterial taxonomic abundances were normalized using lineage-specific 16S rRNA gene copy numbers following the approach described by Angly et al. (71). The resulting OTU abundance data were further processed to compute taxonomic relative abundances (%) for each sample. Subsequent data handling and numerical computations were performed using an in-house MATLAB script (version 23.2.0.2485118, R2023b). Graphical and statistical analyses were performed using Python (version 3.11.13) and R (version 4.5.0). To ensure comparability across all sampling seasons and study years, the primary diversity analyses were performed using taxonomic relative abundance profiles harmonized at the genus level (Supplementary Dataset S1). The 0.5% relative-abundance threshold was applied only for graphical presentation of the taxonomic profiles and was not used to filter the abundance tables employed for alpha- or beta-diversity analyses. Diversity analyses were performed using the complete harmonized genus-level relative abundance profiles. This approach was adopted because the 2022 samples were distributed across independently processed OTU datasets, precluding direct combination of the corresponding OTU feature spaces without introducing potential batch-related bias. Alpha diversity was evaluated using the Shannon and Inverse Simpson indices.

Beta diversity was assessed using Bray–Curtis dissimilarities calculated from genus-level relative abundance profiles. Principal Coordinates Analysis (PCoA) was used to visualize differences in bacterial community composition among samples. PERMANOVA was performed separately for each study year to evaluate the effects of season and management system on bacterial community composition. Because the same colonies were sampled repeatedly across seasons, the significance of the seasonal effect was assessed using 9,999 permutations restricted within colony. The effect of management system, which represented a between-colony factor, was evaluated using exact permutations of management assignments at the whole-colony level (20 possible allocations of three colonies to each management group), thereby preserving the repeated-measures structure of the experimental design. PERMDISP was applied to the same Bray–Curtis distance matrices to test for differences in multivariate dispersion among seasons and between management systems.

For 2023, for which all four seasonal sampling periods were represented within a single OTU dataset, alpha- and beta-diversity analyses were additionally performed directly from the OTU abundance table as a confirmatory analysis. OTU-level alpha diversity was evaluated using the Shannon and Inverse Simpson indices and observed OTU richness, whereas beta diversity was assessed using Bray–Curtis dissimilarities followed by PCoA, PERMANOVA, and PERMDISP using the same experimental factors and permutation structure described above. Statistical significance was set at p < 0.05.

Colony performance parameters

To assess colony performance and development, several parameters were evaluated, including colony population size, quantities of hive-stored honey and pollen, and annual honey production (72, 73). Monitoring was conducted between June and October.

Colony population

Colony population was estimated using the Liebefeld method, based on a visual assessment of the number of adult worker bees covering both sides of each comb frame (74–76).

Hive stored honey and pollen

Colonies were housed in Dadant–Blatt hives (47 × 30 cm frames). Each frame was virtually divided into six equal sections per side, resulting in 12 sections per frame. The proportion of each section occupied by sealed (operculated) honey and pollen was visually estimated. Assuming a concentration of 1 g/cm² of stored material, this corresponded to an average of 0.18 kg per one-sixth of a frame, as previously reported (72, 73, 77).

Honey production

Honey productivity was determined by measuring the difference in weight between supers before and after honey extraction, following the standardized COLOSS protocol (74).

Statistical analysis

All statistical analyses of colony-performance parameters were conducted using R (version 4.5.0). The colony was considered the independent experimental unit, with three independent colonies per management system (n = 3). For hive-stored honey, hive-stored pollen, and colony population, the five repeated measurements obtained during each year were first averaged within each colony to obtain a single annual value per colony. Total honey production was determined directly at the colony level as annual production. Stationary and migratory colonies were compared separately for each year using Welch’s two-sample t -test, which does not assume equal variances between groups. Data are reported as mean ± standard deviation (SD) of the three independent colonies. Statistical significance was set at p < 0.05.

Seasonal dynamics of the honey bee gut microbiota under stationary and migratory beekeeping

The present study provides a comprehensive characterization of the seasonal dynamics of the gut microbiota associated with A. m. ligustica under two different beekeeping management systems, namely stationary and short-distance migratory beekeeping. The migratory colonies were transferred from Castelpoto to Mount Taburno during summer and returned to Castelpoto before the autumn sampling period, whereas the stationary colonies remained in Castelpoto throughout the study.

High-throughput sequencing of the V1–V3 region of the bacterial 16S rRNA gene revealed a taxonomically diverse bacterial community comprising 25 phyla, 43 classes, 91 orders, 159 families, and 242 genera. Taxonomic assignments reaching the species level were obtained for 72 bacterial taxa; however, given the intrinsic limitations of partial 16S rRNA gene sequencing for species-level discrimination, these assignments should be interpreted cautiously and are reported primarily as supplementary descriptive information.

Seasonal taxonomic restructuring of the honey bee gut microbiota

At the genus level, the honey bee gut bacterial community was consistently dominated by Gilliamella and Snodgrassella, although their relative abundances varied across seasons, years, and management groups (Figures 1, 2). Gilliamella showed a marked increase during autumn, reaching 70.1% in stationary colonies and 66.3% in migratory colonies in 2023. Conversely, Snodgrassella was generally more abundant during spring and winter, reaching 52.7% in migratory colonies in spring 2022 and 51.3% in stationary colonies in spring 2023. Overall, these two dominant genera displayed contrasting seasonal patterns, with Gilliamella becoming more prominent in autumn and Snodgrassella contributing more strongly to the community during spring and winter.

Among the less abundant genera, Bartonella showed a pronounced seasonal increase, reaching 28.0% in migratory colonies during summer 2022. Fructobacillus also displayed marked seasonal variation, with a pronounced increase in migratory colonies during autumn 2023. Frischella was detected at lower relative abundances and showed variation across seasons and years.

Overall, the genus-level profiles revealed substantial seasonal variation in bacterial community composition in both stationary and migratory colonies, with differences in the magnitude of these changes between management groups and years. Additional taxonomic profiles at the class, order, family, and species levels are provided in the Supplementary Materials (Supplementary Figures S1–S4, respectively), showing the corresponding seasonal patterns under stationary and migratory management in both study years.

Temporal variation in bacterial community composition under stationary and migratory management

Descriptive taxonomic profiles showed temporal variation in bacterial community composition in both stationary and migratory colonies. Differences in the relative abundance of several bacterial genera between the two management groups were observed during individual seasons, particularly during summer and autumn. However, these descriptive differences should not be interpreted as evidence of a consistent overall effect of migratory management on bacterial community structure, as the multivariate analyses did not detect a significant management effect in either study year.

The seasonal relocation of migratory colonies to Mount Taburno exposed them to environmental and floral conditions different from those experienced by stationary colonies at Castelpoto, which may have contributed to some of the taxonomic differences observed during individual sampling periods. Nevertheless, the relative contribution of floral resources, environmental exposure, and other site-associated factors cannot be disentangled from the present experimental design. Therefore, the observed taxonomic variation is more appropriately interpreted within the broader context of temporal and environmental variation rather than as evidence of enhanced microbial turnover or microbiota plasticity associated with migratory management.

Persistence and temporal variation of core gut symbionts

Similar temporal variation in honey bee-associated bacterial communities has been reported in relation to environmental heterogeneity, seasonal changes, and colony translocation (14, 24). In addition, the seasonal microbiota dynamics of the three stationary colonies included in the present study (hives 3, 4, and 6) were previously examined in a longitudinal analysis restricted to the Castelpoto apiary (61). That study documented seasonal variation in the relative abundance of dominant gut-associated taxa, particularly Gilliamella and Snodgrassella. In the present study, data from these stationary colonies were incorporated into a broader comparative design that also included three short-distance migratory colonies, allowing seasonal microbiota patterns to be evaluated under the two management conditions. Despite substantial temporal variation in relative abundance, dominant honey bee-associated bacterial taxa remained detectable across seasons and years. In particular, Gilliamella and Snodgrassella persisted as major members of the gut bacterial community in both stationary and migratory colonies.

Previous studies have associated these core gut symbionts with important functions related to digestion, nutrient metabolism, detoxification of dietary compounds, immune regulation, and protection against pathogens (20, 78, 79). However, microbial functions were not directly measured in the present study. Consequently, the persistence of Gilliamella and Snodgrassella should not be interpreted as evidence that their metabolic or host-associated functions remained unchanged across seasons or management conditions.

Interestingly, Bifidobacterium, commonly recognized as an important member of the honey bee gut microbiota, was detected at very low relative abundance in the present dataset. Previous studies have identified Bifidobacterium as an important symbiotic taxon, particularly associated with the rectum and hindgut compartments of honey bees (7, 79). Its low relative abundance in the present study may reflect several non-mutually exclusive factors, including ecological variation, differences in floral resource availability, gut sampling procedures, and methodological factors associated with 16S rRNA gene amplification and taxonomic profiling. More generally, amplification and taxonomic profiling based on partial 16S rRNA gene sequences may introduce taxon-dependent methodological biases. Therefore, the low relative abundance of Bifidobacterium should be interpreted cautiously and should not be considered evidence of its ecological absence from the investigated colonies. Overall, the persistence of dominant gut-associated taxa alongside substantial temporal changes in their relative abundances indicates that bacterial community variation occurred without the disappearance of its principal members. Nevertheless, taxonomic persistence alone cannot establish microbiome functional stability, resilience, functional redundancy, or the presence or absence of dysbiosis. Direct functional analyses would be required to determine whether the observed compositional changes were accompanied by changes in microbial activity or host-associated functions.

Seasonal patterns in alpha diversity

Alpha diversity was primarily evaluated using the Shannon and Inverse Simpson indices calculated from genus-level relative abundance profiles (Figure 3); individual diversity values and the corresponding statistical outputs are provided in Supplementary Dataset S2. Both indices revealed pronounced temporal variation in bacterial diversity, whereas no significant overall effect of management system was detected. In 2022, season significantly affected both Shannon diversity (p = 0.0001) and Inverse Simpson diversity (p = 0.0002), while differences between stationary and migratory colonies were not significant (p = 0.600 and p = 1.000, respectively). A comparable pattern was observed in 2023, when significant seasonal effects were detected for Shannon (p = 0.0019) and Inverse Simpson (p = 0.0014), whereas management-related differences remained non-significant (p = 0.700 for both indices).

The marked seasonal variation observed in alpha diversity is consistent with previous evidence indicating that honey bee-associated microbial communities can vary in relation to seasonal changes in diet, floral resources, colony physiology, and environmental exposure (80). Nutritional conditions, particularly the diversity and quality of pollen resources, are known to influence honey bee physiology and colony-level responses (81), providing a plausible ecological context for the temporal patterns observed in the present study. However, because diet composition was not directly linked to microbiota variation and microbial functions were not experimentally assessed, the mechanisms underlying these seasonal changes cannot be established from the present data.

Importantly, the absence of a significant overall management effect indicates that short-distance migratory beekeeping was not associated with a consistent increase or decrease in alpha diversity across the two study years. Although individual seasons and colonies exhibited differences in Shannon and Inverse Simpson values, these variations should therefore be interpreted primarily as temporal patterns rather than as evidence of enhanced microbial diversity in migratory colonies. Moreover, alpha diversity alone cannot be used to infer microbiome functional stability, resilience, functional redundancy, or the presence or absence of dysbiosis. Such properties depend on community composition and function and cannot be established solely from taxonomic diversity indices (82, 83).

As a confirmatory analysis, alpha diversity was additionally evaluated directly from the OTU abundance table for the complete 2023 dataset (Supplementary Dataset S3). OTU-level analyses confirmed a significant seasonal effect for both Shannon (p = 0.0113) and Inverse Simpson (p = 0.0415), whereas no significant overall management effect was detected (p = 0.100 and p = 0.200, respectively). Observed OTU richness did not differ significantly according to management system (p = 0.800) and showed a non-significant seasonal pattern (p = 0.0597). Thus, the OTU-level analysis independently supported the principal finding obtained from the harmonized genus-level profiles: temporal variation in bacterial diversity was primarily associated with season rather than with stationary versus migratory management.

Seasonal structuring of bacterial communities revealed by beta diversity

Beta diversity was assessed using Bray–Curtis dissimilarities calculated from genus-level relative abundance profiles. Principal Coordinates Analysis (PCoA) was used to visualize differences in bacterial community composition among samples across seasons (Figure 4). PCoA coordinates and the complete PERMANOVA and PERMDISP outputs are provided in Supplementary Dataset S2. The first two PCoA axes accounted for 74.75% of the represented variation in 2022 and 63.75% in 2023, indicating that a substantial proportion of the variation among samples was represented by the two-dimensional ordinations.

PERMANOVA revealed a significant effect of season on bacterial community composition in both study years. In 2022, season explained 41.0% of the variation in Bray–Curtis dissimilarities (pseudo- F = 4.63, R² = 0.410, p = 0.001), whereas in 2023 it explained 46.5% of the variation (pseudo- F = 5.80, R² = 0.465, p = 0.0001). In contrast, management system did not significantly explain bacterial community variation in either 2022 (pseudo- F = 0.90, R² = 0.039, p = 0.300) or 2023 (pseudo- F = 0.96, R² = 0.042, p = 0.200). These results indicate that seasonality, rather than stationary versus short-distance migratory management, was the predominant factor associated with variation in bacterial community structure.

PERMDISP analyses were performed to determine whether the PERMANOVA results could be influenced by differences in multivariate dispersion. No significant differences in dispersion among seasons were detected in either 2022 (p = 0.767) or 2023 (p = 0.059). Similarly, multivariate dispersion did not differ significantly between stationary and migratory colonies in either 2022 (p = 0.100) or 2023 (p = 0.400). Thus, the significant seasonal effects detected by PERMANOVA were not accompanied by statistically significant differences in multivariate dispersion, although the 2023 result should be interpreted cautiously given its proximity to the conventional significance threshold.

The predominance of seasonal variation is consistent with previous evidence indicating that honey bee gut bacterial communities can vary in relation to temporal changes in diet, floral resource availability, colony physiology, and environmental exposure (80). Seasonal changes in nutritional resources may therefore provide a plausible ecological context for the temporal community patterns observed here. However, the present study did not directly test the contribution of individual dietary or environmental factors to the observed beta-diversity patterns; consequently, these potential mechanisms should be regarded as explanatory hypotheses rather than demonstrated causal relationships. As a confirmatory analysis, beta diversity was additionally evaluated directly from the OTU abundance table for the complete 2023 dataset. OTU-level PERMANOVA confirmed a significant effect of season on bacterial community composition (pseudo- F = 2.02, R² = 0.233, p = 0.0001), whereas the overall effect of management system was not significant (pseudo- F = 0.98, R² = 0.042, p = 0.400). PERMDISP detected no significant differences in multivariate dispersion among seasons (p = 0.1466) or between management systems (p = 0.100). The OTU-level analysis therefore independently supported the main pattern obtained from the harmonized genus-level profiles, showing pronounced seasonal structuring of bacterial community composition without evidence of a consistent overall differentiation between stationary and migratory management systems.

Colony performance under stationary and migratory beekeeping

Colony performance was evaluated over two consecutive years by measuring adult bee population, hive-stored honey, hive-stored pollen, and annual honey production. These parameters provide complementary information on colony development, nutritional reserves, and productive performance. Honey bee colony strength is closely associated with food availability, nutritional status, and environmental conditions, which can influence brood development, colony maintenance, and productivity (84). Nectar collected by foragers is converted into honey, whereas pollen represents an important nutritional resource that is stored within the hive and contributes to colony nutrition during periods of reduced resource availability (85, 86). Individual colony measurements and the corresponding sampling dates, annual colony-level values used for statistical comparisons, and complete statistical outputs are provided in Supplementary Dataset S4, whereas group-level annual values are summarized in Table 1.

Annual honey production was significantly higher in migratory than in stationary colonies in 2022 (20.83 ± 3.04 vs. 12.20 ± 2.14 kg; p = 0.0195). No significant differences between management systems were detected for hive-stored honey, hive-stored pollen, or colony population during the same year (p > 0.05). Thus, although migratory colonies showed numerically higher mean values for hive-stored honey and colony population in 2022, these differences were not statistically supported.

In 2023, none of the investigated colony-performance variables differed significantly between management systems (p > 0.05). Migratory colonies nevertheless showed numerically higher mean values for hive-stored pollen (1.34 ± 0.27 vs. 1.00 ± 0.13 kg), colony population (20,700 ± 328 vs. 19,150 ± 1,489 adult bees), and annual honey production (17.80 ± 1.80 vs. 14.47 ± 5.04 kg), whereas hive-stored honey was essentially comparable between management systems (3.19 ± 0.41 vs. 3.20 ± 0.68 kg). Given the limited number of independent colonies per management system, these numerical differences should be interpreted cautiously and should not be considered evidence of a management effect in the absence of statistical significance.

Overall, the clearest difference in colony performance was observed for annual honey production in 2022, when migratory colonies produced significantly more honey than stationary colonies. However, this difference was not reproduced in 2023, and no consistent significant effect of management was detected for the other colony-performance parameters. These findings therefore indicate that the productive response varied between years and among the parameters investigated. Differences in floral resource availability and other environmental conditions associated with the two management systems may have contributed to the observed patterns; however, the present experimental design does not allow the effects of individual environmental factors to be separated.

The colony-performance results should also be interpreted independently from the observed variation in gut bacterial community composition. The present analyses do not establish a causal relationship between microbiota variation and colony productivity. Similarly, the persistence of dominant core bacterial taxa across sampling periods indicates compositional continuity of major members of the honey bee gut microbiota, but does not by itself demonstrate functional stability, enhanced resilience, or the absence of dysbiosis. Functional properties of the bacterial communities were not directly measured in the present study.

Seasonal variation in bacterial community composition may reflect changes in diet, floral resources, environmental exposure, or other seasonally varying factors. Several dominant honey bee-associated bacterial taxa, including Gilliamella, Snodgrassella, and Bartonella, have previously been associated with important metabolic and host-related functions (20, 78, 87). Nevertheless, their persistence or variation in relative abundance in the present study should not be interpreted as direct evidence of changes in microbial metabolic activity, ecological adaptation, or improved colony performance.

Taken together, the colony-performance data indicate a significantly higher honey production in migratory colonies in 2022, whereas no significant differences in the investigated performance parameters were detected between management systems in 2023. The observed temporal variation in gut bacterial communities occurred alongside these differences in colony performance, but the present study does not establish a functional or causal connection between the two. The higher honey production recorded for migratory colonies in 2022 therefore indicates a productive advantage under the environmental conditions prevailing during that year rather than a consistent effect of migratory management across years. Further studies involving a larger number of independent colonies and direct measurements of microbial functional activity would be required to determine whether, and to what extent, gut microbiota dynamics contribute to colony performance under migratory beekeeping.

The present two-year study shows that the gut bacterial communities of Apis mellifera ligustica undergo marked temporal variation under both stationary and short-distance migratory beekeeping management. Dominant honey bee-associated taxa, particularly Gilliamella and Snodgrassella, persisted across seasons and years, although their relative abundances varied considerably over time and between management conditions. These findings indicate compositional continuity of major gut bacterial members alongside substantial temporal variation in community structure.

Colony-performance responses were not consistent across years. Migratory colonies showed significantly higher annual honey production than stationary colonies in 2022, whereas no significant differences in the investigated performance parameters were detected between management systems in 2023. Thus, the productive advantage observed in 2022 should be interpreted in the context of the environmental conditions prevailing during that year rather than as evidence of a general positive effect of migratory management on colony productivity.

Although temporal variation in gut bacterial communities occurred alongside differences in colony performance, the present study does not establish a causal or functional relationship between these responses. In particular, the persistence of dominant bacterial taxa should not be interpreted as evidence of functional stability, absence of dysbiosis, or enhanced microbiome resilience, as microbial functions were not directly assessed.

Overall, these findings contribute to understanding how honey bee gut bacterial communities vary over time under contrasting beekeeping management conditions and highlight the importance of considering both seasonal and interannual variability when evaluating the effects of short-distance migratory beekeeping. Future studies involving a larger number of independent colonies and apiaries, together with direct functional characterization of the microbiome and detailed assessment of dietary and environmental variables, will be required to clarify the mechanisms linking management practices, gut microbiota dynamics, and colony performance.