The Epidemiological Situation of the Managed Honey Bee (Apis mellifera) Colonies in the Italian Region Emilia-Romagna

Simple Summary The collapse of honey bee colonies is an important phenomenon worldwide. The individual and synergic actions of pathogens are one of the causes of this decline. Monitoring programs are essential to understand and prevent the epidemiological patterns that are involved. The present study aimed to investigate the health status of honey bees in the Emilia–Romagna region (northern Italy) during the year 2021, on workers from 31 apiaries. The prevalence and abundance of DWV, KBV, ABPV, CBPV, Nosema ceranae, and trypanosomatids (Lotmaria passim, Crithidia mellificae, Crithidia bombi) were investigated four times in the year using molecular methods. Trypanosomatids were not found in any of the samples, while DWV, CBPV and N. ceranae were the most prevalent pathogens. Pathogens had different peaks in abundance over the months, showing seasonal trends related to the dynamics of both bee colonies and Varroa destructor infestation. The results of this study suggest that the monitoring program could be useful to understand the dynamics of honey bee pathogens. Abstract The recent decades witnessed the collapse of honey bee colonies at a global level. The major drivers of this collapse include both individual and synergic pathogen actions, threatening the colonies’ survival. The need to define the epidemiological pattern of the pathogens that are involved has led to the establishment of monitoring programs in many countries, Italy included. In this framework, the health status of managed honey bees in the Emilia–Romagna region (northern Italy) was assessed, throughout the year 2021, on workers from 31 apiaries to investigate the presence of major known and emerging honey bee pathogens. The prevalence and abundance of DWV, KBV, ABPV, CBPV, Nosema ceranae, and trypanosomatids (Lotmaria passim, Crithidia mellificae, Crithidia bombi) were assessed by molecular methods. The most prevalent pathogen was DWV, followed by CBPV and N. ceranae. Trypanosomatids were not found in any of the samples. Pathogens had different peaks in abundance over the months, showing seasonal trends that were related to the dynamics of both bee colonies and Varroa destructor infestation. For some of the pathogens, a weak but significant correlation was observed between abundance and geographical longitude. The information obtained in this study increases our understanding of the epidemiological situation of bee colonies in Emilia–Romagna and helps us to implement better disease prevention and improved territorial management of honey bee health.


Introduction
The health status of managed honey bees (Apis mellifera L.) is threatened by several pests and pathogens [1]. The effects of these infections and their interactions with other biotic (e.g., bee genetics, queen behaviour) and abiotic (e.g., climate changes, lack of forage, exposition to agrochemicals, and management practices) factors can lead to colony losses [2][3][4]. Although both the symptoms and effects caused by contact with pesticides have been investigated [3,4], the role of pathogens and their effects in relation to their All the above-reported pathogens are present in Italy; nonetheless, scant information is available for KBV and trypanosomatids. KBV was found in the years 2013 and 2017 in one single apiary of the Latium region [54] and Vespa velutina individuals, respectively [55]. L. passim was reported as present in the Veneto region only [56].
This study was conducted to contribute to our understanding of the epidemiological situation of honey bee diseases in Italy. Samples from 31 apiaries located in the Emilia-Romagna region were collected and analyzed in the year 2021, to quantify both the prevalence and abundance of N. ceranae, DWV, KBV, ABPV, CBPV, L. passim, C. mellificae, and C. bombi. These pathogens were selected for the following reasons: (i) they are not responsible for notifiable diseases in Italy, (ii) molecular methods are available for their assessment, (iii) frozen workers form suitable samples for their quantification, (iv) they are of special practical and scientific interest. The last point was specifically true for L. passim, C. mellificae, and C. bombi, whose distribution in Italy is largely understudied.

Sampling
Thirty-one apiaries distributed in all provinces of the Emilia-Romagna region were investigated ( Figure 1). The apiaries were named after the province code, to which a progressive letter was added (Table 1).   In the study year (2021), the same three colonies were sampled four times for each investigated apiary, namely in April, June, September, and November. Four apiaries (BOC, FCD, RNA, and RNB) missed the first sampling, as the respective owners joined the monitoring plan later.
For each investigated colony, approximately twenty-five worker bees were sampled from the external combs [57], on which older and, therefore, more likely to be infected workers tend to congregate [58,59]. All the samples were stored at −80 • C until analysis.

Extraction of Nucleic Acids
Each sample was processed as a pool. Ten bees were placed in a 2mL microtube with 300 µL of DNA/RNA Shield (Zymo Research, Irvine, CA, USA) and crushed with a TissueLyser II (Qiagen, Hilden, Germany) for 3 min at 30 Hz, as previously reported [60,61]. The obtained suspension was split into two aliquots, from which DNA and RNA were separately extracted. The extraction of the nucleic acids was performed with the Quick DNA Microprep Plus Kit (Zymo Research) and Quick RNA Microprep Plus Kit (Zymo Research). During the process, the modified manufacturer's instructions for solid tissue processing were followed [19,62]. The obtained nucleic acids were eluted in 100 µL of DNAase-RNase-free water and the extracts were stored at −80 • C until the qPCR assays.

Quantitative Real-Time PCR Assays
The extracted DNA and RNA were used in a Real-Time PCR (qPCR) assay to quantify the abundance of each pathogen in the samples. The DNA was used to analyse N. ceranae and trypanosomatids, whereas RNA was used for the viruses. The qPCRs were performed using the primers reported in Table 2. For each target gene, a total reaction volume of 25 µL was prepared as previously described [57,62] using SYBR™ green assays with forward and reverse primers (2 µM), and 3 µL of nucleic acid extract. For the DNA and the RNA, the SYBR PowerUp™ SYBR™ Green Master Mix (ThermoFisher, Waltham, MA, USA) and the Power SYBR™ Green Cells-to-CT™ Kit (ThermoFisher Scientific) were used, respectively. The qPCRs were performed on QuantStudio™ 3 Real-Time PCR System (ThermoFisher Scientific), following the protocols for either gene sequence [63][64][65][66][67]. DNA and RNA, previously extracted from positive honey bee samples, were used as positive controls. Sterile water was used as a negative control in all analytical steps. All the analyses were conducted in duplicate.

Statistical Analysis
The values of pathogens abundance used for statistical analysis referred to the pool of 10 bees for each investigated colony. For each sample, the pathogen abundance was determined by averaging the two technical replicates of each PCR assay. For each apiary, the pathogen abundance was calculated by averaging the data obtained from the three investigated colonies.
A principal component analysis (PCA) was conducted to explore similarities and possible clusters of apiaries according to a regional longitudinal gradient, using the mean number of pathogen copies as a variable. Calculations were carried out over the whole sampling year. After the PCA, the correlation strength between each pathogen's abundance and longitude value was tested with a Pearson correlation test.
Before proceeding with further analysis, the assumption of normal distribution of errors was checked with Shapiro-Wilk test. As the normality test was failed, a nonparametric approach was followed.
To test the effect of the month of sampling on each pathogen's abundance, Kruskal-Wallis tests were used. When the last was found to be significant, Dunn's test was applied as a post-hoc pairwise comparison. The apiary average prevalence of each pathogen was also calculated per month as a percentage ratio between the positive and total samples. A Chi-square test pairwise comparison with Yates' correction was then applied to investigate the pathogen prevalence among the sampling months.
For all statistical tests, a significance threshold at α = 0.05 was assumed. Pathogens whose presence was not detected in any of the samples by real-time PCR were excluded from analyses.
All the statistics were calculated with R studio version 4.1.2. using packages Fac-toMineR and ggplot2.

Results
In total, 360 samples were analyzed. None of the samples were found to be positive for any of the three trypanosomatid species. The highest prevalence, considering the four sampling points, was found for DWV (56.0%), followed by CBPV (48.5%), N. ceranae (42.9%) and ABPV (7.3%). Only three colonies were found to be positive for KBV infection (0.8%) (Table S1).

Geographical Distribution
In the PCA, two components from the pathogen abundance were chosen at each apiary. The first one (PC1) was positively correlated with DWV and N. ceranae. The second one (PC2) was positively correlated with these pathogens, as well as with ABPV and CBPV ( Figure 3 and Table 3).
The principal component analysis showed one cluster along the second principal component, corresponding to a group of easternmost apiaries (FCD, RNA, RNB, FEC), whereas the other apiaries were assembled closer to the chart origin. Three apiaries (MOA, REC and BOA) located in central areas of the region clustered following the second principal There was a positive correlation for N. ceranae (t = 2.192, p-value = 0.0365, cor = 0.377) and longitude value. No significant correlation was found for DWV (t = 1.805, p-value = 0.081, cor = 0.318), ABPV (t = −0.477, p-value = 0.637, cor = −0.088), and CBPV (t = −0.165, p-value = 0.8715, cor = −0.031). No correlation was performed for KBV, because only three apiaries were found to be positive in the total sampling year.
investigate the pathogen prevalence among the sampling months.
For all statistical tests, a significance threshold at α = 0.05 was assumed. Pathogens whose presence was not detected in any of the samples by real-time PCR were excluded from analyses.
All the statistics were calculated with R studio version 4.1.2. using packages FactoMineR and ggplot2.

Results
In total, 360 samples were analyzed. None of the samples were found to be positive for any of the three trypanosomatid species. The highest prevalence, considering the four sampling points, was found for DWV (56.0%), followed by CBPV (48.5%), N. ceranae (42.9%) and ABPV (7.3%). Only three colonies were found to be positive for KBV infection (0.8%) (Table S1).

Geographical Distribution
In the PCA, two components from the pathogen abundance were chosen at each apiary. The first one (PC1) was positively correlated with DWV and N. ceranae. The second one (PC2) was positively correlated with these pathogens, as well as with ABPV and CBPV (Figure 3 and Table 3).

Seasonal Trend
The prevalence of the considered pathogens showed a distinctive seasonal trend (Figure 4). DWV and N. ceranae were most prevalent in April (48.2% and 72.8%, respectively), while their prevalence decreased in June and rose in September and November. On the other hand, ABPV decreased over the months. The prevalence of DWV increased over time, with a peak in November (98.9%). CBPV showed a fluctuating trend, with significant increases in June (65.9%) and November (62.4%). The prevalence of KBV was very low throughout the year and did not significantly fluctuate over time (Figure 4).

Seasonal Trend
The prevalence of the considered pathogens showed a distinctive seasonal trend (Figure 4). DWV and N. ceranae were most prevalent in April (48.2% and 72.8%, respectively), while their prevalence decreased in June and rose in September and November. On the other hand, ABPV decreased over the months. The prevalence of DWV increased over time, with a peak in November (98.9%). CBPV showed a fluctuating trend, with significant increases in June (65.9%) and November (62.4%). The prevalence of KBV was very low throughout the year and did not significantly fluctuate over time ( Figure 4). The results of the Kruskall-Wallis tests conducted on the abundance data of each pathogen are reported in Table 4. Except for KBV, the sampling month resulted in a significant predicting variable. Total pathogen abundance was higher in April and November ( Figure 5). The results of the Kruskall-Wallis tests conducted on the abundance data of each pathogen are reported in Table 4. Except for KBV, the sampling month resulted in a significant predicting variable. Total pathogen abundance was higher in April and November ( Figure 5). The results of the post-hoc tests (Supplementary Table S3) revealed that the majority of significance in the comparisons was found between April and June and between April and November. ABPV and N. ceranae had a higher mean abundance in April, while DWV and CBPV had higher values in November ( Figure 6).   Table S3) revealed that the majority of significance in the comparisons was found between April and June and between April and November. ABPV and N. ceranae had a higher mean abundance in April, while DWV and CBPV had higher values in November ( Figure 6). The results of the post-hoc tests (Supplementary Table S3) revealed that the majority of significance in the comparisons was found between April and June and between April and November. ABPV and N. ceranae had a higher mean abundance in April, while DWV and CBPV had higher values in November ( Figure 6).

Discussion
The frequent losses in honey bee colonies that occurred in recent years in Europe and other continents promoted the establishment of monitoring plans in many countries, to improve our understanding of the underlying factors and adopt the needed countermeasures. This study is part of BeeNet, an Italian multiannual monitoring plan that is conducted nationwide. Despite the large scale of the study, we presented results obtained in apiaries situated in the Emilia Romagna region. We focused on a range of pathogens that were reputed to be implicated in the colony losses [1] detected in the first year of the monitoring action (2021). Five out of eight of the considered pathogens were found to infect the examined apiaries, namely, the viruses DWV, CBPV, ABPV, and KBV, and the microsporidium N. ceranae.
In general, DWV was the most prevalent pathogen, following a previous study carried out in Italy in 2009-2010 where the prevalence value was higher (68%) [2]. A similar prevalence value of DWV-A was recorded in the Veneto region in 2020 [56]. DWV was reported as the most prevalent virus in Europe [68] and the other continents [27]. CBPV was the second most prevalent pathogen, occurring in almost half of the analyzed honey bee samples. The higher prevalence values of this virus (82.2-98.8%) were recently found in Veneto [56]. Since, in 2009-2010, the prevalence of CBPV in Italy was only 8% [2], this study confirms the exponential increase in the incidence of this virus observed over the last decade, as in other countries, including the UK [42], USA [69], and China [70]. ABPV and KBV were the least prevalent viruses in Emilia Romagna, especially KBV, which was only detected in three hives (0.8%). KBV was first identified in Italy in 2013 [54] and has rarely been found since then [2,56], with a predominant distribution in Asia, the USA and Australia [71].
Within the non-viral pathogens, the prevalence of N. ceranae in Emilia Romagna was similar to overall reports for Italy in 2009-2010 [2] but lower than the prevalence recorded in other European countries, including Belgium (93%), Serbia (79-95%) and Spain (66%) [72][73][74]. The three investigated trypanosomatids were not found in any of the colonies. In the last national survey, trypanosomatids were not searched for; therefore, the only data available for comparison are from the neighbouring Veneto region, showing a high prevalence of L. passim (48.8-62.2%) in 2020-2021 [56]. Similar results (62.3%) were obtained in a nine-year survey for trypanosomatids, conducted in Serbia [75]. The only other Italian study investigating trypanosomatids was carried out on honey samples supplied by beekeepers in Northern Italy. A total prevalence of L. passim of 78% was found in honey, with Friuli-Venezia Giulia, Emilia-Romagna and Trentino-Alto Adige scoring the highest values [76]. Data on both the presence and distribution of C. mellificae and C. bombi are even more scattered. C. mellificae was detected in honey bees for the first time in Italy in 2010, in the Latium region [77]; however, after the discovery of L. passim, the sequence was attributed to this new trypanosomatid. In 2020, C. mellificae was found in one single hive in the Veneto region [56]. The presence of trypanosomatids in honey bees is still poorly monitored both in Italy and elsewhere, although they are recognized as among the most prevalent bee parasites, contributing to the increased colony mortality in Europe and the USA [44]. However, our poor understanding of the epidemiologic details of trypanosomatid infections means that comparisons among the results obtained over different years and locations are excessively speculative. This is even more true when the results refer to dissimilar kinds of samples (e.g., worker bees from individual colonies and bulks of honey produced by different hives). Coinfections between L. passim and N. ceranae [75,78] were not confirmed by the present study.
Concerning the abundance of detected pathogens, the highest loads were found for CBPV and DWV (7.60 × 10 10 and 1.36 × 10 10 copies, respectively). A viral load higher than 10 10 pathogen copies is generally related to the presence of disease symptoms for both viruses [29,30,79]. Therefore, the analyzed bees belonging to the monitored colonies did not present any symptoms related to both viral infections. It should be noted that the virulence and subsequent symptomatology of these viruses can vary depending on the means of transmission and the presence of co-infections/infestations [80,81]. ABPV, N. ceranae and KBV were found to have a lower abundance. Given the higher virulence of ABPV and KBV [31,80], apiaries with detected viral loads of around 10 6 should be kept under control.
The geographical distribution of pathogens in the Emilia-Romagna regions highlighted a correlation between the DWV and N. ceranae infection in the East part (Romagna) compared to the west of the regions (Emilia). These differences could be related to the different health management of apiaries by beekeepers, probably due to different traditions or customs. N. ceranae and V. destructor (the main vector of DWV) originate from Asian countries and can replicate in a warm environment [17,82]. The apiary in the east Emilia-Romagna has a closer proximity to the sea, which contributes to a more suitable environment for N. ceranae and DWV infections [83][84][85][86][87]. Additionally, in 2021, a higher average temperature was recorded in Romagna (17 • C) than in Emilia (14 • C), creating a more suitable environment for these organisms (https://simc.arpae.it/dext3r/ (accessed on 13 April 2021)). In the global view of climatic change, this may lead to an increase in cases of exotic pathogens [10,88,89].
In this investigation, the seasonality of pathogen occurrence followed different specific trends for each pathogen. Usually, the prevalence and the abundance of pathogens are strictly related to the bee population dynamics, especially to brood cycles [7,10,14,90].
DWV and ABPV transmission is also linked to V. destructor [26,28]. Their rate of prevalence, abundance and virulence followed the same trends of mite infestation [69,[91][92][93]. The high infection of ABPV and DWW in April may be linked to a constant infection that remains active during the winter. This leads honey bees being infected in the new brood cycle in the spring [94][95][96]. The slight decrease in terms of prevalence and abundance in June could be in line with the increase in bee population, diluting the pathogens among more hosts [97]. The proliferation of both V. destructor and viruses cannot be limited by acaricide treatments. Although mite infestations can successfully be reduced, the virus can persist inside the hosts at low prevalence and loads [10]. In this case, the ABPV cycles decrease over the season, and the summer acaricide treatments seem to be effective at inhibiting its spread within the bee population, as previously reported [94,96]. On the other hand, the summer treatments did not affect the persistence and circulation of DWV, since its infection increased after the summer. Instead, the winter reductions in population and the winter acaricide treatments seem to be useful to limit the damage caused by higher rates of DWV infection.
In recent studies, the CBPV was considered to be strictly linked to spring and related to colony developmentl however, cases of this virus were detected in all beekeeping seasons and during the winter [38,98]. The analyzed CBPV trends confirmed the importance of the re-emerging disease [42]. The virus was present in all seasons in the investigated apiary, as previousy reported [38,41,99]. The higher abundance in April and September reflected the higher prevalence in June and November. These findings could be related to the first and classic aetiology of the disease, where a higher number of cases were found in late spring/early summer [38,39]. The absence of a temporal pattern may be associated with the oral transmission of CBPV without a vector, or due to the genetic diversity of the virus [38,98].
Due to the very limited KBV-positive colony detected in this study, its seasonal dynamics cannot be analyzed. This disease is not present in Italy, except for random findings in some bee colonies in the Latium region and V. velutina sampled in the Liguria region [54,55].
N. ceranae is the most diffused and most frequent microsporidia detected in honey bee colonies worldwide [1,100]. Previous surveys in Europe highlighted the differences and fluctuations in N. ceranae prevalence and abundance depending on the year and the geographical localization of the apiary [101]. In this investigation, its prevalence and abundance were detected in all samples, confirming the lack of seasonality that was previously reported [102,103]. The highest incidence was observed during the colony development in April and in November before the wintering, in line with the life cycles of the microsporidia [104], as previously reported in other countries [100,[105][106][107][108][109].

Conclusions
This investigation focused on the health of managed honey bees, confirming the important role of monitoring surveys on bee diseases. The impact of management practices on bee colonies leads to an intensification of these surveys, not only because the impact of the diseases causes direct damage to beekeeping, but also because the circulation of these pathogens may affect other pollinators [110]. The health of the honey bee is strictly connected to the health of pollinators, and vice versa. A change in perspective is necessary, considering a One-Health Approach to honey bee diseases [111][112][113].
Further studies are needed to better understand the health status of Italian and European honey bee colonies to prevent the diffusion of diseases and to protect the environment.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/vetsci9080437/s1, Table S1. Mean number of copies for investigated apiary in the Emilia-Romagna region during April, June, September, and November 2021. Table S2. The mean annual abundance of DWV, KBV, ABPV, CBPV and N. ceranae in the investigated apiary. Table S3. p values resulting from the post-hoc Dunn's test for pairwise comparisons (ns = non-significant).