Lake Sinai virus is a diverse, globally distributed but not emerging multi‐strain honeybee virus

Domesticated honeybees and wild bees are some of the most important beneficial insects for human and environmental health, but infectious diseases pose a serious risk to these pollinators, particularly following the emergence of the ectoparasitic mite Varroa destructor as a viral vector. The acquisition of this novel viral vector from the Asian honeybee Apis ceranae has fundamentally changed viral epidemiology in its new host, the western honeybee A. mellifera. While the recently discovered Lake Sinai Viruses (LSV) have been associated with weak honeybee colonies, they have not been associated with vector‐borne transmission. By combining a large‐scale multi‐year survey of LSV in Chinese A. mellifera and A. cerana honeybee colonies with globally available LSV‐sequence data, we investigate the global epidemiology of this virus. We find that globally distributed LSV is a highly diverse multi‐strain virus, which is predominantly associated with the western honeybee A. mellifera. In contrast to the vector‐borne deformed wing virus, LSV is not an emerging disease. Instead, demographic reconstruction and strong global and local population structure indicates that it is a highly variable multi‐strain virus in a stable association with its main host, the western honeybee. Prevalence patterns in China suggest a potential role for migratory beekeeping in the spread of this pathogen, demonstrating the potential for disease transmission with the man‐made transport of beneficial insects.


| INTRODUC TI ON
Humans rely on livestock and other beneficial animals, but rearing and shipping animals in large numbers not only poses the risk of zoonoses but also of the spread and emergence or re-emergence of animal diseases, as can currently be seen with the spread of avian influenza or African swine fever. This risk is not confined to vertebrates but is also highly relevant to invertebrates, from crustacean and molluscan aquaculture to beneficial insects used for pest control, pollination or increasingly as a human protein source. For example, several bumblebee species, which provide valuable pollination services, have been effectively domesticated over the last decades (Velthuis & van Doorn, 2006). But mass rearing of the North American species Bombus occidentalis was discontinued due to infestations with the previously rare microsporidian Nosema bombi, and this practice was linked to spillover and declines in wild American bumblebee species (Cameron et al., 2016). The sustainable intensification of agriculture will increasingly rely on mass-reared insects, for example, for integrated pest management or for human food and animal feed production. Managed honeybees are amongst the most economically important, common and long-standing domesticated invertebrates, with the earliest records of beekeeping dating back to the Ancient Egyptians. Understanding host-pathogen interactions in managed honeybees can inform us of potential pitfalls in developing mass-rearing operations for other beneficial insects. As is typical for livestock, beekeeping operations are at risk of falling prey to infectious diseases and parasites, with high colony mortality rates and issues such as Colony Collapse Disorder threatening the sustainability of beekeeping (Kulhanek et al., 2017) and poor management resulting in diseased colonies (Jacques et al., 2017). At the same time, emerging honeybee pathogens have been shown to cause spillovers into wild hymenoptera (e.g. Fürst et al., 2014;Loope et al., 2019;Manley et al., 2019Manley et al., , 2020Santamaria et al., 2018).
Because of their economic importance, microbial pathogens and macro-parasites of honeybees have therefore been intensively studied for a long time. With the recent advances in meta-transcriptomics, this field has seen an unprecedented expansion, leading to the discovery of many new pathogens, with over 80 viral pathogens now identified in honeybees (Beaurepaire et al., 2020). Prominent amongst the recently discovered viruses are the Lake Sinai Viruses (LSV), a multi-strain virus (Bigot et al., 2017;Cornman, 2019;Daughenbaugh et al., 2015) originally discovered by a meta-transcriptomic survey in the USA (Runckel et al., 2011). LSV variants have since been detected globally, including in Europe, the Middle East, West and South Africa, South America, Asia, Australia and Pacific Islands using both primer-dependent surveys and meta-transcriptomic approaches (e.g. Amakpe et al., 2016;Brettell et al., 2020;Cavigli et al., 2016;Cepero et al., 2014;D'Alvise et al., 2019;Faurot-Daniels et al., 2020;Gamboa et al., 2015;Malfroy et al., 2016;Ravoet et al., 2013Ravoet et al., , 2015Remnant et al., 2019;Roberts et al., 2017Roberts et al., , 2020Shojaei et al., 2021;Simenc et al., 2020;Thaduri et al., 2018;Tozkar et al., 2015;Traynor et al., 2016;Yang et al., 2020). LSV strains were frequently found to be amongst the most prevalent viruses in these surveys. While the pathology of LSV has not yet been studied experimentally, two studies of Colony Collapse Disorder in the US and Spain found potential links with collapsing honeybee colonies (Cepero et al., 2014;Cornman et al., 2012). However, a systematic study found no association of LSV prevalence with overwinter colony mortality in Belgium (Ravoet et al., 2013). Consistently, studies have found that LSV is associated with weak honeybee colonies: in several surveys from the US, colonies with low population size showed the highest LSV prevalence (Daughenbaugh et al., 2015;Faurot-Daniels et al., 2020;Glenny et al., 2017).
In the western honeybee A. mellifera, host-virus interactions have been dramatically affected by the acquisition of the ectoparasitic mite Varroa destructor as a viral vector. Varroa destructor jumped host from the Asian honeybee A. cerana at the beginning of the last century and has since been spread globally (Oldroyd, 1999).
In particular, this vector acquisition has changed the epidemiology of the important honeybee virus deformed wing virus (DWV) at a global level, increasing virulence and prevalence (reviewed in Martin & Brettell, 2019) with the repeated emergence of epidemic strains (Manley et al., 2019;Paxton et al., 2022;Ryabov et al., 2017;Wilfert et al., 2016) that have spilled-over into other hymenoptera, which themselves are not parasitised by the mite (Fürst et al., 2014;Loope et al., 2019;Manley et al., 2019;Santamaria et al., 2018). Beyond DWV, the prevalence and potentially virulence of viruses such as, for example, slow bee paralysis virus and black queen cell virus also increases in the presence of Varroa (e.g. Carreck et al., 2010;Glenny et al., 2017;Manley et al., 2019;Mondet et al., 2014;Remnant et al., 2019;Santillan-Galicia et al., 2010;Traynor et al., 2016). This pattern however is consistently different for LSV. In a large multiyear survey in the US, Glenny et al. (2017) found no association of LSV and Varroa levels, whereas Traynor et al. (2016) found that LSV levels were highest in colonies with low mite numbers. This pattern has also been verified experimentally: in colonies in which mite levels were artificially lowered or increased, LSV levels were not affected by treatment in contrast to other viruses, notably DWV (Norton et al., 2021). Additionally, while LSV has been detected in Varroa mites (Shojaei et al., 2021), it was found not to replicate in Varroa, in contrast to the detection of active replication in A. mellifera .
While LSV does not appear to be associated with Varroa, it has also been found in several other bee species including in bumblebees in South America (Gamboa et al., 2015) and Europe (Parmentier et al., 2016;Simenc et al., 2020;Toplak et al., 2020) and in solitary bees in Europe , Israel (Daughenbaugh et al., 2021), the US (Dolezal et al., 2016;Iwanowicz et al., 2020) and Australia (Brettell et al., 2020), as well as in ants (Bigot et al., 2017) and in a wasp from China (Yang et al., 2020); additionally it has been detected in butterflies (Pislak Ocepek et al., 2022) and the Chinese mantis Tenodera sinensis (Genbank accession LC646160). In both bumblebees (Parmentier et al., 2016) and solitary bees , replication was verified by strand specific PCR. However, the closest relatives of A. mellifera, the other Asian honeybee species including A. cerana have so far not been systematically examined for LSV prevalence. In Asia, A. cerana is of similar economic importance as A. mellifera and large beekeeping operations for example in China frequently keep both of these cave-dwelling species. While A. cerana is considered to be more robust than A. mellifera, due to stronger hygienic behaviour (e.g. Lin et al., 2016), and generally less affected by infectious diseases (e.g. Ai et al., 2012;Chantawannakul et al., 2016;Chen et al., 2021;Diao et al., 2019), colony mortality is still considerable with an average colony mortality of 12.8% between 2011 and 2014 . The main pathogen of concern in this species is Sacbrood virus, which has led to severe colony losses in outbreaks across Asia, while it is comparatively rare and not of concern in A. mellifera (Chen et al., 2021). This contrasting pattern of prevalence and pathogenicity in A. cerana and A. mellifera highlights the need to study common bee viruses in both domesticated species. While A. cerana does not appear to be the source of the reemerging DWV epidemic (Diao et al., 2019), it is the source of two emerging parasites in A. mellifera: the virus-vectoring V. destructor and the microsporidian Nosema ceranae both jumped to non-native western honeybees in the last century, and have been anthropogenically spread around most of the globe (Klee et al., 2007;Wilfert et al., 2016). In North American A. mellifera populations, N. ceranae infections are linked with LSV (Daughenbaugh et al., 2015;Glenny et al., 2017). This positive association further raises the question of whether LSV may be an important pathogen in A. cerana.
Combining a large-scale multi-year survey of LSV in Chinese A. mellifera and A. cerana colonies with genetic analyses of globally available LSV-sequence data, we investigated firstly whether LSV is specific to A. mellifera based on prevalence and phylogenetic patterns. Secondly, we asked whether LSV shows the signature of a stable association with its host, in contrast to known emerging viruses associated with the acquisition of vector-borne transmission and the ensuing anthropogenic global spread of infected bees. As an emerging disease is typically characterised by strong exponential growth and an absence of geographical population structure, we investigated this question through population genetics and phylogenetic reconstruction. Investigating the results of a honeybee health survey in China, we were also able to investigate potential seasonal and geographical variation in LSV prevalence.  Table S1) as part of routine health screening. No clinical symptoms of virus infection or other pathogens were found in these colonies.

| Field samples
Colonies were kept according to guidelines for beekeeping practice and, in the case of A. mellifera, regularly monitored and treated with acaricides for mites twice a year. There was no effect on colony development reported due to mite infestation during the previous year. Per F I G U R E 1 Collection sites of A. mellifera (a) and A. cerana (b) apiaries in China. At each site, between one and nine apiaries were sampled (indicated by the size of the circle), and within each apiary, three colonies were haphazardly chosen for sampling. LSV prevalence at the location level is indicated by the colour of the circle and was calculated as the number of positive colonies sampled in one location (see Table S2). apiary, in excess of 30-50 adult worker bees were collected from inside each of three haphazardly selected colonies, placed in a small iron wire cage for direct transport or in 50 mL sterile tubes with dry ice for transport to the laboratory, where the samples were stored at −80°C.

| RNA extraction
All individual samples were submerged in liquid nitrogen and ground into a fine powder using a pestle and mortar; this material was used in Trizol-based RNA extractions as described in (Diao et al., 2019).
Total RNA was initially extracted from pooled honey bee samples (30-50 individuals) using TRIzol (Invitrogen) following the manufacturer's instruction in a TissuePrep homogenizer (Gening Scientific).
The obtained RNA was dissolved in 20 μL of sterile water and stored at −80°C prior to analysis. The quantity and purity of the RNA were measured using a Nanodrop spectrophotometer (Thermo Scientific), with samples included for further analysis if OD 260/280 values were in the range of 1.8-2.0.

| RT-PCR amplification
The extracted RNA was used to synthesise first-strand cDNA using an oligo (dT) primer and GoScript Reverse Transcriptase (Promega).
To test for the presence of LSV, we used the non-strain specific primers used as reported in detail in (Bigot et (Bigot et al., 2017) and have been submitted to Genbank under the accession numbers MF667689-MF66702.

| Prevalence measures and statistical analysis
To take account of PCR assay efficiency and sensitivity (set at 95% as indicative values, based on Reiczigel et al., 2010), we used the R library epiR v.0.9-82 (Stevenson et al., 2013) and the function epi.
prev to calculate the true prevalence with 95% confidence based on methods in Blaker (2000). We used Rstudio v. 1.3.1093 with R version 3.5.3 in all analyses. To test whether LSV prevalence at the apiary-level (the number of LSV positive and negative colonies out of three colonies randomly sampled in an apiary at one timepoint, see Table S1) was affected by host species, we conducted a Fisher's exact test, as the sample size for A. cerana did not allow further statistical tests. To test whether geographical location (coded as longitude and latitude) and season affects LSV prevalence at the apiary level in A. mellifera, we used the lme4 package (v1.1-12; Bates et al., 2015) to run generalised linear mixed models (GLMMs) with binomial error distribution and a logit link function. The full model included square-root transformed longitude and latitude as well as season (spring (March-May), summer (June-August) and autumn (September-November) modelled as a continuous trait), with collection year and location as a random factor. The minimum adequate model was identified through removal of non-significant terms and comparison of models using ANOVA.

| Phylogenetic reconstruction and population genetics
For the gene trees and the hypervariable region, we constructed phylogenetic trees in MrBayes 3.2.6 (Ronquist & Huelsenbeck, 2003), using gamma-distributed site-specific general time-reversible models, with parameters estimated from the data during the analysis.
We ran two runs of two chains for 10,000,000 MCMC generations, sampling trees every 10,000 generations. The tree was plotted using the sample size and genetic variation contained in the alignments was not sufficient to include host species as a trait. We report the sampling date of the included sequences, however, spanning a time period of 13 years, the dataset does not contain sufficient temporal information according to analysis in TempEst v 1.5.3. TrNef+I + G was determined as the most suitable substitution model using Jmodeltest (v.2.1) (Darriba et al., 2012;Guindon & Gascuel, 2003).
As the open reading frames of ORF1 and the rdrp-gene overlap in this region, we treated each codon position individually. We used the path sampling maximum likelihood estimator to compare demographic priors (constant population size and exponential growth) and molecular clock rates (strict clock, relaxed clock with lognormal or exponential rate distribution); for all viruses, exponential growth and a relaxed exponential clock were selected. For the exponential population size prior, we used a lognormal distribution with mu = 1 and sigma = 10. All models were checked for convergence in Tracer (v1.6) and run long enough to obtain effective sample sizes >200 for all parameters, with a 10% burn-in, sampling the chain at equal distances to obtain a total of 10,000 trees per analysis. We Tajima's D to test for an excess of rare polymorphisms. To assess the degree of population structure we calculated K ST , a measure of population differentiation based on the proportion of between-population nucleotide differences (Hudson et al., 1992) as well as S NN , which measures how often genetic nearest neighbours come from the same population (Hudson, 2000) between either geographical regions (Asia, Middle East, Europe, North America, South America and Oceania) or host groups. For the latter analysis, we grouped all seven sequences isolated from other Apidae (Andrena vaga (n = 1), Bombus lapidarius (n = 1), Bombus pascuorum (n = 1), Bombus terrestris (n = 2) and Vespa mandarinia (n = 2)) into a group of 'other Apidae' due to the low sample size, with the western honeybee A. mellifera represented by 247 samples and the eastern honeybee A. cerana by 7. Additionally, we calculated the nucleotide diversity π for geographical regions and LSV strains.

| Prevalence of LSV in Chinese apiculture
Sixty-nine Apis cerana colonies from 23 apiaries and 201 A. mellifera colonies from 67 apiaries were sampled across 20 Chinese provinces over a 2.5-year period between March 2015 and September 2017 (see Table S1). Seventy-three of 201 A. mellifera colonies were LSVpositive, whereas this was only the case for a single A. cerana colony from Fujian. Taking into account the uncertainty of PCR-based pathogen detection, assuming a PCR assay sensitivity and specificity of 95%, the estimated true LSV colony-level prevalence, that is, the percentage of LSV-positive colonies, in Chinese A. mellifera populations reaches 34.8% (95% confidence interval (CI) % 27.51-42.51%), whereas the true colony-level prevalence of LSV for A. cerana is estimated only at 0.46% (0%-6.6%), in this case assuming an assay sensitivity and specificity of 99% to account for the low prevalence. At the apiary level, for A. mellifera, 47 of 67 sampled apiaries were positive, while only one of 23 A. cerana apiaries were positive. The difference in prevalence between species is highly significant at both the apiary and colony levels (Fisher's exact test, p < .001 at both levels).
To test whether LSV prevalence (i.e. the number of positive vs. negative colonies in an apiary) in A. mellifera was affected by the geographical location (latitude and longitude; see Figure 1) and the season in which the samples were collected, we analysed prevalence at the apiary level. Location and collection year were accounted for F I G U R E 2 Schematic representation of the LSV-genome, based on the genbank reference genome NC_032433 (black bar). Light grey bars represent the genes/open reading frames and the analysed fragments are shown as dark grey bars, indicating the length and number of samples per alignment. as random factors. Only longitude had a significant effect (χ 2 = 10.45, df = 1, p = .001) on LSV prevalence compared to the null model, with no significant effects of latitude (χ 2 = 2.445, df = 1, p = .118) and no effect of season (χ 2 = 0.468, df = 1, p = .494). To test that the effect of longitude was not caused by the two negative apiaries from Xinjiang in China's northwest (Figures 1 and 3, Figure S1), we tested the effect of longitude when excluding these two data points. The effect of longitude, with prevalence highest in eastern China, persisted when removing these data points (χ 2 = 4.136, df = 1, p = .042).

| Phylogenetic and population genetic analysis
We first reconstructed the LSV phylogeny based on the highly variable 5′-region of the rdrp-gene (263 sequences, 340 bp), excluding any sequences with unknown host origin or from poorly represented geographical region. As Figure 4b shows, a Bayesian phylogeny reconstructs the major LSV groups identified by (Bigot et al., 2017) and and one solitary bee from Europe, two wasps as well as one Varroa and one Chinese mantis from Asia. As shown in Figure 4, these non-A. mellifera isolates are always nested in clades containing A. mellifera-isolates from the same geographical region and study. To allow phylogenetic reconstruction and population genetic analysis of host origin despite the low sample size, we grouped all non-A. mellifera isolates together. We found no evidence for genetic differentiation between LSV sequences isolated from A. mellifera or other bee species. While K ST was non-significant (N other_Apidae = 16, N A.mellifera = 347), the nearest neighbour statistic S NN was significant (p = .013), with a high value of 0.942 consistent with the dominance of A. mellifera as a host species (94% of all isolates). Discrete trait modelling with host species as a discrete trait in BEAST 1.10.4 indicates that A. mellifera is the source population with a root state probability of 91.07% for A. mellifera. All subsequent analyses combine LSV sequences from all host species. In addition to the prevalence analysis within Chinese A. cerana and A. mellifera apiaries, phylogenetic and population genetic analysis thus also point towards LSV being predominantly an A. mellifera virus based on the available data.

LSV was originally isolated by NGS sequencing in the United
States (Runckel et al., 2011) and has since been found in surveys in South America, Europe, Australia and nearby islands (grouped here as Oceania), Asia, South Africa and the Middle East. Based on the sample sizes, we were able to include populations from Asia, Europe, the Middle East, North America, South America, and Oceania for population genetic and phylogenetic analysis. As Figure 4 shows, LSV clades are mostly represented in all six populations, with the major exception of Sister group 2 -only reported from Asia, South America, Oceania as well as South Africa (see Figure 5 and Figure S3 -isolate KY354244). Additionally, strain C (LSV-1) has so far not been found in China and strain B has only been found in Europe, the F I G U R E 3 Lake Sinai viruses (LSV) prevalence per apiary plotted against longitude, showing that LSV was less prevalent in western regions. This relationship remains significant when removing the data from Xinjiang province in the far west (see Figure S1).    Figure S2 in the supplemental material, the accession numbers are indicated alongside the trees.

(a) (b)
Europe to Asia and the Middle East and from the Middle East to South America, as well as from Oceania to Asia, Europe and North America (Bayes Factors >10, posterior probabilities >0.75; see Table S4 for Bayes Factors between all populations), in addition to weaker support for a transition from South America to Oceania (Table S4). Only the latter is not supported by analysing transitions only in sister group 1.
We found no evidence for recent population expansion across LSV. Tajima's D, estimated using a coalescent simulation in DNAsp v.5.10.1, was not significant across the entire phylogeny, nor within clades, with the exception of clade B (Tajima's D = −1.961, p < .05).
While exponential growth was the preferred demographic model in the coalescent sampler BEAST (Suchard et al., 2018), the confidence interval of the exponential growth rate overlapped 0 for the whole phylogeny as well as within the two sister groups.

| DISCUSS ION
Combining a multi-year survey of A. cerana and A. mellifera colonies in China with a global analysis of LSV phylogenetics and population F I G U R E 5 MrBayes tree of partial ORF1, rdrp-and capsid-genes for LSV. The hypervariable overlapping region between ORF1 and the rdrp-gene was removed for these gene trees. The geographical location of samples is colour coded and host origins other than A. mellifera are indicated in the trees with a filled circle. Lake Sinai virus (LSV) strain names used in previous literature are indicated in the tree (LSV-S2: LSV sister group 2; strains A and C have also been referred to as LSV-1 and LSV-2). In Figure S3 in the supplemental material, the accession numbers are indicated alongside the trees. structure, we found that given the currently available data, LSV appears to be predominantly an A. mellifera virus in a stable association with its host. This multi-strain virus shows no signs of recent viral emergence with strong population structure and no evidence of exponential growth, even when analysing patterns within strains.
While LSV was common in A. mellifera with an estimated true prevalence of 34.8%, it was nearly absent in A. cerana (0.46%) in a PCR-based prevalence study of Chinese apiaries, even though it has been found in more distantly related hymenopteran species (e.g. Bigot et al., 2017;Brettell et al., 2020;Daughenbaugh et al., 2021;Dolezal et al., 2016;Gamboa et al., 2015;Iwanowicz et al., 2020;Parmentier et al., 2016;Ravoet et al., 2014;Simenc et al., 2020;Toplak et al., 2020;Yang et al., 2020). This result is consistent with a primerindependent meta-transcriptomic study, where very low read counts for LSV were found in an A. cerana pool from Papua New Guinea compared to read counts two orders of magnitude higher in the A. mellifera pool (~300 vs. ~ 475,000 combined LSV reads; Roberts et al., 2020).
We found that all 16 sequences available from non-A. mellifera hosts cluster with sympatric A. mellifera isolates (see Figure 4) and that there is no effect of population structure by host. Taken together, this suggests that LSV is predominantly an A. mellifera-virus, particularly with regards to the closely related A. cerana, and that there are no host species-specific clades. However, LSV is a multi-strain virus with high genetic diversity and primer-dependent methods are unlikely to uncover the full genetic diversity, and consequently, the true preva- and South America; unusually, even though LSV epidemiology has been well studied in North America (e.g. Daughenbaugh et al., 2015;Faurot-Daniels et al., 2020;Glenny et al., 2017), having been initially discovered in the US (Runckel et al., 2011), North American honeybee populations are also under-represented. Additionally, while we found LSV to be predominantly a virus of the western honeybee A. mellifera, it would be crucial to sample the other honeybee species (e.g. A. andreniformis, florea, dorsata, laboriosa, koschevnikovi and nigrocincta) which are distributed across regions in Asia and, in the case of A. florea, into the Middle East and East Africa.
In stark contrast to the two emerging strains of the Varroavectored DWV (Diao et al., 2019;Manley et al., 2019;Wilfert et al., 2016), LSV shows no evidence of population expansion but strong population structure both globally and locally, as well as within viral strains, which in turn mostly show a wide global distribution.
LSV is thus a highly diverse multi-strain virus, with multiple, genetically diverse strains co-existing globally and locally, again in contrast to the strong genetic bottle necks and loss in genetic diversity experienced by DWV with the spread of vector-borne transmission (Loope et al., 2019;Manley et al., 2019;Martin et al., 2012;Wilfert et al., 2016). LSV shows population genetic patterns indicating that it is a stable and potentially long-established viral infection of A. mellifera, in contrast to DWV but also other viruses such as slow bee paralysis virus that have been affected by the acquisition of vectorborne transmission.
Given the currently available data, no well-supported geographical origin of the entire LSV clade is reconstructed, with similar root probabilities of around 35% for both Oceania and South America. This is driven by sister group 2, which so far has only been found through meta-genomic studies in Oceania, South America and Asia, as well as in a single individual in South Africa, while it has so far not been found in similar meta-transcriptomic studies covering Europe and North America (e.g. Cornman et al., 2012;Ravoet et al., 2014;Remnant et al., 2017;Runckel et al., 2011;Simenc et al., 2020;Thaduri et al., 2018;Tozkar et al., 2015). In contrast,  This raises the potential that, above and beyond a lack of mechanistic association with Varroa such as tissue tropism, seasonal changes in LSV prevalence might be driven by complex patterns of competitive exclusion, for example allowing LSV to attain high prevalences when Varroa-dependent viruses are at their seasonal low. These interactions warrant experimental investigation, particularly given the association of LSV with low colony strength (Daughenbaugh et al., 2015;Faurot-Daniels et al., 2020;Glenny et al., 2017).
Strong patterns of seasonality are common in human and animal diseases, with seasonal flu and respiratory viruses being prominent examples whose transmission depend on climate variation, population density and host behaviour, such as increased contact indoors (e.g. Moriyama et al., 2020). Honeybees are eusocial and highly seasonal organisms, with colonies growing throughout spring and summer with up to around 50,000 individuals per colony at the peak, then reducing reproduction in autumn and overwintering with a rump population of a few thousand bees that are confined to the colony.
Accordingly, viral prevalence in honeybees frequently shows seasonality, including in Chinese populations (Diao et al., 2019), with prevalence lowest in spring in the colonies that survive the winter period and peak prevalences in summer or autumn (e.g. D'Alvise et al., 2019; Diao et al., 2019;Runckel et al., 2011). For LSV in Chinese A. mellifera, we found no effect of season on prevalence in the present study.
Studies in the US (Faurot-Daniels et al., 2020;Traynor et al., 2016) and with other honeybee pathogens, that in turn can have different seasonal patterns, may also contribute to a lack of seasonality.
In contrast to the absence of seasonality, we found a clear effect of longitude on LSV prevalence in Chinese A. mellifera populations, with lower prevalences concentrated towards the central and western regions of China ( Figure 3). This pattern is robust to removing the samples from the far western region of Xinjiang, although it should be noted that sampling density was highest in Eastern regions, particularly around Beijing. This regional pattern in prevalence could potentially be explained by variation in climate or vegetation, but also by apiary density and relative geographical isolation between apiaries, which may be influenced by migratory beekeeping. An estimated 80% of Chinese beekeepers engage in some form of migratory beekeeping (Zheng et al., 2018), with migration from South to North across the year. Migratory beekeeping in China is not well studied, but there is evidence that there are eastern migration routes from Guangdong and Guangxi along the coast or in the east, for example, through Henan and Hubei, separate from more central/western migration routes from Yunnan to Gansu and Xinjiang (in Grillot, 2020).
This raises the potential that migratory routes may drive LSV prevalence as well as other pathogens. As recently reviewed, the impact of migratory beekeeping on disease spread in honeybees and wild bees is an urgent research gap (Martínez-López et al., 2022). Both natural migration and the shipping of diseased animals can lead to the spread of infectious diseases in wildlife, livestock and humans, as illustrated by diseases such as avian influenza or West Nile Virus. With increased domestication and mass-rearing of invertebrates, the importance of man-made transportation for disease spread will likely become more important, particularly for beneficial insects that are used in the field for ecosystem services such as pollination or pest control. For example, hoverflies such as the marmalade hoverfly Episyrphus balteatus are increasingly commercially produced for aphid control, as well as playing an important role in pollination; not only are these insects then transported by humans, but they themselves are migratory between southern and northern Europe (Wotton et al., 2019).
While the precise role of shipping and transhumance or migration on disease spread need to be further studied, the cautionary advice from the role of non-native bumblebees and honeybees in the spread of infectious diseases both within the species itself in commercial or domesticated settings as well as in the spill-over to native wild species is clear: Newly domesticated species that are deployed in the field or cannot be kept in complete confinement, the use should at least be restricted to their native range, taking account of population structure within the insect population but also their pathogens, in order to reduce the risks of population displacement (e.g. Ings et al., 2006) and the spread of disease. Furthermore, the recent discovery of LSV and the continuing addition of genetic variants through metagenomics is a potent reminder that metagenomic, rather than primer-based, surveillance is needed to detect highly variable and evolving pathogens not only in human and livestock populations but also in the increasing array of commercially used mass-reared invertebrates.

ACK N O WLE D G E M ENTS
We would like to express our thanks to beekeepers from all sample sites for assistance in collecting bee samples and for valuable

CO N FLI C T O F I NTE R E S T S TATE M E NT
The authors declare no competing interests.

DATA AVA I L A B I L I T Y S TAT E M E N T
All data on prevalence of LSV are provided as Supporting Information.
All sequences are available on GenBank (accession numbers see Tables S1 and S3).

R E FE R E N C E S
Ai, H. X., Yan, X., & Han, R. C. (2012). Occurrence and prevalence of seven bee viruses in Apis mellifera and Apis cerana apiaries in China.