Выведенное поголовье медоносных пчел обладает устойчивостью к Varroa destructor и связанной с ним вирусной трансмиссией

Аннотация

Эктопаразит Varroa destructor является самой большой угрозой для управляемых колоний медоносных пчел (Apis mellifera) во всем мире. Несмотря на значительные усилия, новые методы борьбы с клещом и переносимыми им патогенами показали ограниченную эффективность, поскольку хозяин остается наивным. Перспективное решение заключается в создании устойчивых к варроа пчел, но недостаток данных по селекции ограничивает их широкое распространение. Здесь мы описываем динамику паразитов и вирусов в устойчивой к варроа популяции медоносных пчел, названной «Pol-line», с помощью крупномасштабного продольного исследования. Результаты демонстрируют заметное снижение уровня варроа в этой популяции, уменьшение титров трех основных вирусов (DWV-A, DWV-B и CBPV) и двукратное увеличение выживаемости. Уровни четвертого вируса, не связанного с варроа — BQCV — не различаются между поголовьями, что свидетельствует о нарушении пути передачи. Кроме того, мы показали, что в отрыве от влияния уровня варроа, титры вирусов не являются сильными независимыми предикторами риска смертности колонии. Эти результаты подчеркивают необходимость переоценки этиологии варроа и предполагают, что производные запасы представляют собой приемлемое решение пандемии варроа.

Введение

Потери колоний медоносных пчел (Apis mellifera) представляют собой серьезную и повсеместную проблему как для отрасли опыления мигрантов, так и для современной сельскохозяйственной безопасности в целом (1,2,3,4). Высокая ежегодная смертность стала повсеместным аспектом коммерческого опыления, несмотря на значительные научные и законодательные усилия по ее снижению (5,6). Только в США, согласно текущим данным, смертность в 2018-2019 годах составляет 37,5% для коммерческих пчеловодческих хозяйств, что значительно выше среднего показателя за 13 лет, составляющего 28,8% (7,8). Причин этого множество, а взаимосвязанный характер посторонних стрессовых факторов, влияющих на здоровье колоний, затрудняет определение основных взаимодействий (9,10). Однако, несмотря на сложность системы, паразитический клещ Varroa destructor неоднократно фигурирует в качестве единственного главного фактора глобальных потерь колоний (1,11,12,13).

Варроа ассоциируется с различными патологиями, включая нарушение развития, иммуносупрессию и поведенческие изменения (14,15,16). Обширный массив исследований показывает, что клещ переносит и усиливает множество вирусных патогенов (17), некоторые из которых сильно связаны с гибелью колоний (12,18,19), и открытие новых типов вирусов продолжается (20). Кроме того, недавние работы показали, что процесс питания варроа сам по себе может нанести значительный ущерб, поскольку эти клещи потребляют жировое тело хозяев медоносной пчелы — орган с многочисленными жизненно важными функциями (12,21). Несмотря на усилия по сдерживанию, варроа продолжает свою 50-летнюю глобальную экспансию в ранее незараженные регионы (1,13,22,23) с соответствующими последствиями для управляемых и диких опылителей (1,24). Одновременно с этим популяции клещей демонстрируют растущую устойчивость к акарицидным обработкам, которые в настоящее время представляют собой последнюю эффективную линию защиты (17,25,26,27,28). Таким образом, угроза, которую представляет собой эта пандемия и ее быстрое распространение в индустриальных сельскохозяйственных системах (29,30), является как непосредственной, так и серьезной.

Многообещающее, устойчивое решение в попытке смягчить последствия варроа заключается в развитии устойчивых к варроа популяций медоносных пчел (17,31,32). Этот подход имеет преимущества перед другими стратегиями управления, поскольку он интегрирован и менее подвержен простым эволюционным адаптациям, которые угрожают химическим методам (17,33,34). В природных популяциях A. mellifera широко известны случаи устойчивости к клещам Варроа (35,36,37). Из-за разрушительного характера инвазии, не поддающейся управлению, естественный отбор колоний, способных сопротивляться или сосуществовать, происходит быстро и интенсивно (38). Механизмы выживания разнообразны, хотя есть некоторые данные, указывающие на конвергентный отбор признаков в изолированных популяциях (17,37). Примечательно, однако, что многие из механизмов, предрасполагающих колонии к устойчивости к варроа в природе, такие как меньшие размеры колоний и повышенная частота роения (38,39,40), противоположны характеристикам, наиболее ценимым в коммерческом пчеловодстве (17). Действительно, усиленный характер динамики варроа-хозяин в управляемых операциях (41,42), в сочетании с эволюционной новизной ассоциации (43), требует механизмов устойчивости, приспособленных к экстремальным условиям. Следовательно, отбор наследственных признаков устойчивости, подходящих для коммерческого пчеловодства, наряду с поддержанием благоприятных характеристик на уровне колонии, является важной задачей для селекционной работы (17,38).

В настоящее время продолжаются крупные программы по созданию устойчивого к варроа товарного поголовья (38) , при этом все большую поддержку оказывает селекция с помощью маркеров (MAS) (44,45) и технологии ‘omics (46). Искусственно отобранные признаки, включая чувствительную к варроа гигиену (VSH), при которой рабочие удаляют зараженный клещами расплод (47,48,49); низкий рост популяции клещей (low MPG), который включает в себя набор поведений, ограничивающих размножение варроа (50,51); и авто-/алло-груминг, процессы удаления и/или повреждения форетических клещей (52,53), показали перспективность в коммерческих сценариях (17,38,54). Однако внедрение полученных линий было ограничено (54), что, предположительно, связано с различной эффективностью в реальных условиях (31,32,55), нехваткой контролируемых крупномасштабных испытаний для подтверждения коммерческой жизнеспособности (38,56,57) и проблемами масштабируемости селекционной работы (58,59). Эти факторы порождают инертность промышленности в отношении эффективной интеграции устойчивых к варроа медоносных пчел, в результате чего глобально управляемая популяция A. mellifera остается адаптивно наивной и, следовательно, остро восприимчивой к паразитизму варроа.

Здесь мы используем контролируемый продольный анализ, чтобы оценить эффективность новой и генетически отличной популяции медоносных пчел, устойчивых к варроа, «Pol-line » (54,60). Это поголовье уникально тем, что оно обладает значительной устойчивостью к клещам, характерной для поведения VSH (61), в сочетании с желательными пчеловодными характеристиками коммерческих итальянских колоний, включая большой размер популяции, значительное производство меда и послушный нрав (54). Мы использовали годичный экспериментальный план, сравнивая Pol-line и коммерческих итальянских медоносных пчел, в котором все колонии отслеживались индивидуально, и в ключевые моменты времени проводились измерения паразитов, патогенов и состояния здоровья. Колонии были интегрированы в миграционную операцию в США, сфокусированную на опылении миндаля и производстве меда, и охватывающую штаты Миссисипи, Южная Дакота и Калифорния. Эта система примечательна тем, что представляет собой один из наиболее интенсивных и стрессовых сценариев в современном пчеловодстве (62,63) и одновременно включает три различных климатических зоны по Кёппен-Гейгеру (64). Таким образом, он обеспечивает строгую проверку функциональности колоний и не ограничивается локальным применением. Дизайн исследования манипулировал применением акарицидов, маршрутом миграции и управлением колониями, чтобы сравнить производительность обоих колоний, одновременно отслеживая уровни варроа и соответствующие титры вирусных патогенов с высокой точностью.

Помимо оценки функциональности медоносных пчел Pol-line, мы использовали нашу экспериментальную систему для изучения относительного влияния варроа и связанных с ним вирусов на смертность в полевых условиях, а также предсказательной полезности того и другого для прогнозирования состояния колоний. В то время как влияние варроа на распространение и усиление вирусов широко изучалось (10,17), конкретные последствия питания варроа как такового, в отличие от влияния передаваемых им вирусов, остаются малоизученными на двухстороннем уровне (12,17,29). Такая скудность информации обусловлена сложностью эффективной изоляции обоих факторов, учитывая их тесно связанную природу, и требуемой стоимостью контролируемых полевых исследований (1,12,65,66,67).

Согласно современным представлениям, РНК-пикорнавирус Deformed wing virus (DWV) — в частности, основные варианты DWV-A (68) и DWV-B (69) — представляет собой главную вирусную угрозу для здоровья колоний (70). Толчком для этого послужила корреляция между DWV и смертностью колоний, а также распространение вируса при переносе варроа (23,71). В присутствии заражения варроа вирус DWV может увеличиться с латентного уровня 6-13% до 100% распространенности в колониях, что сопровождается увеличением титра вируса в миллион раз в течение 3 лет (1,22). В то же время существуют доказательства рекомбинации штаммов DWV, что приводит к появлению более вирулентных гибридных форм, и, что очень важно, отбор на эти формы, когда они передаются в виде вирусной примеси через кормление варроа (68,70,72). Хотя распространенность DWV хорошо установлена, исследования роли вирусного разнообразия и распространения в определении патогенности показали противоречивые результаты (68,73,74,75,76,77,78), и, что удивительно, некоторые исследования, связывающие вирус со смертностью, не учитывают влияние уровня варроа на исход колонии (73,79,80,81). Последнее может быть симптомом концептуализации варроа как вирусного переносчика в первую очередь, и повреждающего агента во вторую; однако это мнение не согласуется ни с биологией клеща (12,17), ни с установленными крупномасштабными данными (67,76,82,83).

Помимо известных вирусов, переносимых варроа, некоторые другие патогены демонстрируют связь с варроа, хотя их точные отношения с клещами неясны (17). Одним из таких случаев является таксономически не идентифицированный РНК-вирус, вирус хронического паралича пчел (CBPV) (84,85). Несмотря на то, что CBPV является одним из первых описанных вирусов медоносной пчелы (86) , он сравнительно мало изучен. Однако известно, что он вызывает смертность взрослых пчел (84), коррелирует с потерей колоний (87,88,89), имеет широкий спектр хозяев (86,90), а его распространенность в мире растет (86,91). Примечательно, что CBPV может существовать на скрытом уровне в здоровых колониях, может увеличиваться в присутствии паразитизма Варроа (92,93,94), и был обнаружен в Варроа (90,95). Однако прямой путь передачи остается неуловимым (84) , а наблюдательные исследования указывают на различную связь (90,96).

Еще одним потенциальным партнером варроа является РНК-пикорнавирус Black queen cell virus (BQCV) (17,97,98). BQCV широко распространен во всем мире (97), вызывает смертельные инфекции у личинок королевы (99), но обычно остается бессимптомным у взрослых рабочих (98,99), демонстрируя ограниченную связь с потерями колоний (100,101). Как и CBPV, BQCV был обнаружен у варроа (95,102), и есть доказательства увеличения передачи вируса другим хозяевам, когда клещи присутствуют (103). Однако, несмотря на то, что некоторые исследования показали, что титры BQCV зависят от зараженности варроа (22,104), другие указывают на то, что вирус не зависит от уровня варроа в колониях (105,106,107), а убедительные доказательства передачи вируса от варроа еще предстоит выяснить (22). Первый пункт очень важен, поскольку текущее понимание предполагает, что хотя передача BQCV, опосредованная варроа, может быть возможной, ее возникновение, по крайней мере, редкость (98). Следовательно, количественная оценка этой несколько загадочной взаимосвязи будет полезна для более широкого изучения динамики варроа-вирусов.

Противоречивые результаты этиологических исследований DWV (66,73,76,81,108), CBPV (90,94,106,109) и BQCV (17,95,105), усугубленные исторически неточной парадигмой исследования питания варроа, которая, по-видимому, недооценивала прямое повреждение жирового тела (12), свидетельствуют о необходимости постоянного изучения смертности колоний, опосредованной вирусом варроа. Пока этот процесс не будет лучше понят, он остается препятствием для эффективного лечения, поскольку нет четкого консенсуса относительно относительной важности ковариативных факторов (1,12,66,110,111).

Таким образом, используя продольно отслеживаемые колонии в нашей экспериментальной установке, мы проверили влияние уровня варроа, а также титров DWV-A, DWV-B, CBPV и BQCV на смертность колонии. Это позволило оценить вирусы, как тесно (DWV-A и DWV-B), так и слабо (CBPV и BQCV) связанные с паразитизмом варроа. Затем мы контролировали уровень варроа, группируя колонии по уровням заражения клещами, чтобы охарактеризовать относительную аддитивную предсказательную силу каждого вируса в определении смертности колоний, используя эпидемиологический подход. Таким образом, мы стремились оценить эффективность нового поголовья медоносных пчел, устойчивого к варроа, и одновременно исследовать удельный вес факторов, влияющих на смертность колоний от варроа.

Материалы и методы

Колонии

Создание колоний происходило до начала исследования, с марта по май 2017 года, в штате Миссисипи, США. Используя установленные методы, колонии без маток, полученные из крупного коммерческого пчеловодческого хозяйства, были уравнены до среднего расчетного размера популяции ~ 7000 рабочих (112) и размещены в 10-рамочных ульях Лангстрота (Таблица S1). После акклиматизации в течение 24-48 часов каждый из них получал ячейку с маточным расплодом, содержащую маток из одной из двух колод, добавленных к одному и тому же базовому уровню рабочих. Использовались итальянские «коммерческие» матки, выращенные из созданных сотрудниками маток-селекционеров и, таким образом, представляющие промышленный стандарт, и устойчивые к варроа матки «Pol-line » (54). Для обеспечения последовательности все матки выращивались в одних и тех же колониях «клеточного строителя», расположенных в Лаборатории разведения, генетики и физиологии медоносных пчел Министерства сельского хозяйства США в Батон-Руж, штат Луизиана, США. Колонии каждого вида содержались на независимых пасеках, расположенных в 80 км друг от друга для поддержания физической изоляции; для контроля генетической верности девственные королевы были открыто спарены с трутнями того же вида посредством насыщения трутней. Через 14 дней после появления маток колонии были осмотрены, а спаренные матки были помечены краской на грудной клетке для облегчения идентификации: белый цвет соответствовал коммерческим, а синий — Pol-line. Колониям дали акклиматизироваться в течение шести недель до начала отбора образцов, и те колонии, которые не смогли достичь успеха в спаривании или имели неприемлемо высокий уровень варроа [≥ 3,0 «клещей на сто пчел» (MPHB)], были удалены, нормализовав среднюю разницу между колониями по варроа до < 0,1 MPHB с момента начала исследования (Таблица S1). Затем каждой колонии был присвоен случайный идентификационный номер, всего 366 колоний; 193 с коммерческими матками и 173 с матками Pol-line. Все колонии получали дополнительный раствор сахарозы и коммерческие пыльцевые добавки на основе сои ad libitum, в соответствии со стандартной отраслевой практикой.

Схема эксперимента

Для того чтобы оценить влияние миграционных процедур опыления, колонии были случайным образом распределены по одному из трех миграционных маршрутов. Миграционные маршруты были следующими: группа «Калифорния», состоящая из 156 колоний, 80 коммерческих и 76 Pol-line, переехала в Южную Дакоту в мае, затем в Калифорнию в октябре для зимовки и опыления миндаля, и обратно в Миссисипи в конце исследования; группа «Миссисипи», состоящая из 156 колоний; 80 коммерческих и 76 Pol-line, перемещенных в Южную Дакоту в мае и обратно в Миссисипи через Калифорнию в октябре для зимовки; и «стационарная» группа, состоящая из 54 колоний; 32 коммерческих и 22 Pol-line, которые оставались в Миссисипи на протяжении всего исследования (рис. 1). 1). В то время как в Южной Дакоте колонии были распределены по четырем пасекам, в Калифорнии и Миссисипи они содержались на одной пасеке. Такая практика отражает стандартный промышленный протокол и обеспечивает более репрезентативную выборку за счет учета местных различий в климате и кормах.

В дополнение к маршруту миграции, колонии были разделены по частоте применения акарицидной обработки, чтобы лучше выяснить влияние различных слоев инвазии. В рамках каждого миграционного маршрута и поголовья половина колоний получила «высокую» акарицидную обработку, а половина — «низкую». Колонии в группе с высокой степенью обработки обрабатывались дважды, в сентябре и декабре; в то время как колонии в группе с низкой степенью обработки обрабатывались только один раз, в декабре. Все обработки проводились с использованием акарицида амитраза.

Отбор проб

В ходе исследования колонии отбирали пять раз: в мае, июне, сентябре, декабре и феврале. В каждой колонии осматривали царицу, а рамки с рабочими особями, нектаром, пыльцой и расплодом оценивали визуально по процентному покрытию (113). При наблюдении, если у королевы не было существующей метки, на ее грудную клетку наносили одну метку красной краской, чтобы обозначить статус вытеснения. Затем, убедившись, что королева изолирована, чтобы предотвратить непреднамеренный отбор проб, две пробы по ~ 300 рабочих были взяты в герметичные пакеты; одна была заморожена сухим льдом и хранилась при — 80 °C для анализа на патогены, а другая была заморожена обычным льдом для анализа на варроа. Для получения репрезентативной возрастной выборки рабочие особи были извлечены из двух расплодных рамок, содержащих запечатанный и незапечатанный расплод. Рабочих равномерно перемешивали с помощью 20-литрового ведра перед сбором с помощью чашечного мерника и последующим отнесением к одному из двух типов образцов. Затем королева была возвращена в улей вместе с рабочими, не попавшими в выборку.

Анализ на варроа

Анализ на варроа проводился для всех колоний во всех пяти временных точках. Для количественной оценки уровня заражения варроа в колониях использовался метод промывки детергентом (114). Образец пчел сначала помещали в чашку с сетчатыми перегородками и твердым дном, чтобы обеспечить извлечение клещей. Затем добавляли воду до полного погружения пчел и перемешивали моющее средство на основе ПАВ (Procter & Gamble), обеспечивая равномерное покрытие. Крышка чашки была закреплена, и она была помещена в возвратно-поступательный шейкер модели E5850 (Eberbach) для перемешивания в течение 60 минут при 120 об/мин. После перемешивания нижнюю фракцию чашки снимали, а содержимое выливали в неглубокий лоток для подсчета количества клещей, провалившихся сквозь сетку. Этот процесс повторялся до тех пор, пока не были зарегистрированы два последовательных нулевых счета, указывающих на то, что все клещи были удалены. Затем пчел отделяли и подсчитывали с помощью счетчика, чтобы получить точный размер выборки. Наконец, общее количество варроа было преобразовано для получения стандартизированного значения «клещи на сто пчел» (MPHB) следующим образом:

Анализ возбудителей

Были выбраны четыре специфические вирусные мишени: DWV-A, DWV-B, CBPV и BQCV. Анализы проводились на рандомизированном подмножестве 92 колоний из калифорнийской и миссисипской миграционных групп, разделенных равномерно по поголовью и выживаемости. Были выбраны четыре основные временные точки: июнь, сентябрь, декабрь и февраль. Все анализы проводились с помощью RT-qPCR в соответствии с ранее установленными методами (94,115).

Выделение РНК

Для каждой колонии и временной точки из основных образцов ~ 300 пчел случайным образом отбирали пулы из 65 пчел, помещали их в 30 мл пробирки с бисером 19-6358Z (Omni) и хранили при температуре — 80 °C для гомогенизации с помощью Bead Ruptor Elite (Omni). Сразу после гомогенизации в пробирки с образцами добавляли 5 мл раствора для гомогенизации (Promega) и 4 мл стерильного буфера 1xPBS, каждый при температуре 4 °C, перед вихревым перемешиванием в течение 20 с с помощью вортекса Vortex-Genie 2 (Scientific Industries). Затем 1,8 мл каждого образца переносили в пробирки объемом 2 мл (Eppendorf) и центрифугировали в течение 60 с при 5000 об/мин при 4 °C с помощью центрифуги 5430-R (Eppendorf). Затем 400 мкл полученного супернатанта переносили в картриджи Maxwell RSC-48 (Promega), подготовленные в соответствии с инструкциями производителя, и проводили экстракцию РНК с использованием тканевого протокола Maxwell RSC simplyRNA (Promega) (Maxwell RSC simplyRNA Tissue Kit, AS1340, TM416, 2019).

Затем чистоту и выход образца рассчитывали в двух экземплярах с помощью микрообъемного UV-Vis спектрофотометра NanoDrop One (Thermo Fisher Scientific) и проводили соответствующие разведения РНК для обеспечения стандартизации образца до 100 нг/мкл. После стандартизации синтез кДНК проводили с использованием наборов для обратной транскрипции QuantiTect (Qiagen), используя термоциклер T100 (Bio-Rad), выполняя пользовательские протоколы синтеза кДНК QuantiTect (Qiagen) (QuantiTect Reverse Transcription Handbook, 2009). Полученные образцы хранили при температуре — 20 °C до начала анализа методом RT-qPCR.

RT-qPCR

RT-qPCR проводили на 1 мкл аликвотах каждого образца, в трех экземплярах, в общем объеме реакции 10 мкл, используя SsoAdvanced Universal SYBR Green Supermix (Bio-Rad), на 96-луночных оптических ПЦР-планшетах Multiplate (Bio-Rad). Праймеры, использованные для количественного определения DWV-A (116), DWV-B (117), CBPV (96) и BQCV (115), были описаны ранее и проверены на специфичность. Полная информация о последовательности приведена в таблице S2. Для каждой мишени были включены отрицательные контроли, состоящие из 1 мкл свободной от нуклеазы H2O (Promega), которые также проводились в трех экземплярах. Кроме того, для определения среднего количества вирусов были получены стандартные кривые для всех мишеней путем десятикратного серийного разведения известных количеств вирусов, охватывающих 8 порядков величины. Линейность (r2) и эффективность реакции (E) поддерживались на уровне ≥ 0,990 и ≥ 92,5%, соответственно, для всех анализов (Таблица S2). Все анализы проводились на детекторных системах CFX Connect Real-Time PCR (Bio-Rad) с использованием ранее оптимизированных термических протоколов (Таблица S3).

Преобразование данных

Образцы анализировали в трех экземплярах для получения среднего значения Ct. Все трипликаты со стандартным отклонением ≥ 1,0 рассматривались, и расходящиеся реплики удалялись, в результате чего было удалено 26 технических реплик. Если это не приводило к снижению стандартного отклонения до < 1,0, сам образец заменялся и исключался из дальнейшего анализа. Средние значения Ct затем количественно сравнивались со стандартными кривыми для каждой мишени, чтобы рассчитать абсолютный вирусный титр, используя следующее уравнение:

Из-за широкого распределения величин, присущего динамике вирусной репликации (118) , эти значения были затем лог-трансформированы, чтобы получить данные, пригодные для дальнейшего анализа. При наличии необнаруженных значений принималось постоянное среднее значение Ct, равное 45, для последующего преобразования данных (105).

Диапазон факторов оценки

Анализы были разделены на измерения выживаемости колоний, постороннего влияния на них уровней варроа (MPHB) и вирусных титров (log₁₀ среднее количество вируса на мкл кДНК), а также того, как поголовье, обработка клещей и маршрут миграции модулируют эти взаимодействия. Важно отметить, что мы изучили как общее влияние уровня Varroa и вирусных титров на выживаемость колоний, так и предсказательную силу этих показателей в сезонных временных точках. Кроме того, мы определили объяснительную силу вирусных титров, когда колонии были подобраны по уровню варроа, таким образом, выясняя относительное аддитивное влияние вирусов, передающихся через варроа.

Выживание колоний и королев во времени, размер популяции и производство меда

Сначала сравнивали выживаемость колоний и маток с течением времени, размер популяции выживших колоний и производство меда между двумя колониями, используя данные временного ряда, собранные в течение всего исследования. Колонии классифицировались как «мертвые», когда количество рабочих рамок падало до < 1,0, и они подтверждались как свободные при любом последовательном отборе проб. Это определение гибели колонии оставалось неизменным во всех дальнейших анализах. Королевы классифицировались как «мертвые», когда существующая королева не могла быть найдена, и присутствовала новая немаркированная королева или клетки королевы, что указывало на смену (119). Для сравнения численности популяции использовали подсчет рабочих рамок во всех колониях, классифицированных как «живые», то есть имеющих ≥ 1,0 рамок в конце исследования, в феврале. Производство меда оценивалось по чистому весу (кг) меда, извлеченного из каждой колонии в сентябре, исключая колонии в группе стационарной миграции, для которых извлечение меда не регистрировалось.

Факторы, влияющие на выживаемость колоний

Затем мы определили влияние уровня варроа и титров вирусов на выживаемость колоний, используя данные всех пяти и четырех временных точек, соответственно. Влияние каждого фактора, а также поголовья, обработки клещей, маршрута миграции и их согласованных взаимодействий, было оценено количественно. Отслеживая отдельные колонии во времени, в сочетании с их уровнем варроа и титром вируса, мы смогли определить размер влияния факторов и их взаимодействие, приводящее к гибели или выживанию в феврале. Февраль был выбран в качестве интересующей нас временной точки, так как он выпал непосредственно перед опылением миндаля и, таким образом, давал репрезентативный и коммерчески значимый показатель силы колоний, входящих в эту систему. Поскольку измерения вирусного титра проводились для подмножества колоний, для них и для показателей варроа проводился отдельный анализ, чтобы сохранить максимальный размер выборки в каждом случае. Эти данные затем использовались для последующего прогностического анализа на основе уровней варроа и титров вирусов.

Прогностические различия в уровнях варроа и титрах вирусов

Чтобы оценить первоначальные прогностические различия в уровнях варроа и вирусных титрах в зависимости от результатов выживания и запасов, мы связали статус выживания в феврале с двумя прогностическими временными точками: июнь и сентябрь. Примечательно, что первый из них охватывает основной период роста популяции колоний, а второй — период сбора меда, предшествующий зимовке колоний. Эти моменты представляют собой ключевые точки, на которых должны приниматься управленческие решения, и, таким образом, являются практически значимыми в качестве прогностических ориентиров. В анализе колонии классифицировались по поголовью, то есть по их статусу — живые или мертвые в феврале. Затем уровни варроа и титры вирусов сравнивались между группами, что представляет собой базовую проверку полезности данных, полученных в один момент времени, для определения результатов.

Прогностическая сила уровней варроа и титров вирусов

Для количественной оценки эпидемиологической значимости уровня варроа и вирусных титров в определении выживаемости колонии и, следовательно, их прогностической силы, мы провели анализ относительного риска (ОР) для этих факторов, используя гибель колонии к февралю как результат выбора. RR — это общепринятая эпидемиологическая мера, используемая для определения вероятности данного исхода для здоровья при воздействии определенного фактора риска и, таким образом, для сравнения величины риска (120). RR-анализ проводился как для уровней варроа, так и для титров вирусов, опять же для прогностических временных точек июня и сентября. Значимые взаимодействия были дополнительно подтверждены с помощью анализа относительного риска (AR), чтобы обеспечить меру размера эффекта; и дополнительного Байесовского анализа RR, чтобы визуализировать распределение вероятности риска. Анализы были объединены по запасам, чтобы лучше охватить весь спектр страт заражения и инфекции. Кроме того, они были стратифицированы по уровню варроа, чтобы определить, как изменения в тяжести заражения влияют на риск гибели колонии.

Затем мы выделили аддитивный эффект вирусных титров как прогностический фактор, не зависящий от Варроа. Мы посчитали это важным тестом на самостоятельное влияние патогенов, связанных с варроа, на потери колоний, вызванные варроа, поскольку обычно эти факторы трудно разделить. Для этого подгруппа из 60 колоний была сопоставлена по сентябрьскому уровню варроа (средний MPHB, < ± 0,1 между результатами выживания) и оставлена нестратифицированной по титрам вирусов. Сентябрь был выбран в качестве ключевой точки сопоставления для данного анализа, так как на этом этапе колонии показали самые высокие уровни варроа и вирусные титры, и, таким образом, эти значения считались наиболее значимыми для результатов выживания в зимний период. Впоследствии колонии были разделены на живые и мертвые когорты на основании их состояния в феврале, чтобы получить группы, отличающиеся по результатам выживания и титрам вирусов, но с тесно связанными распределениями уровней варроа. Затем мы провели RR-анализ для этой подгруппы, изучая титры вирусов в июне и сентябре как прогностические факторы для выживания колоний, таким образом, представляя собой меру вирусного воздействия, отделенную от варроа.

Модель Варроа

Чтобы продемонстрировать связь между уровнем варроа и смертностью колоний, а также создать модель прогноза лечения для обоих поголовьев, мы объединили все данные по колониям в страты заражения варроа, исходя из уровня варроа в сентябре, и определили процент смертности на основе показателей выживаемости в каждой страте. Затем мы подогнали модели кривых к данным по двум запасам, что позволило предсказать результаты смертности для уровня заражения варроа и запасов, основанные на обширном эмпирическом наборе данных.

Статистический анализ

Для всех данных парных сравнений, в которых анализ предполагал бы нормальное распределение, мы провели тесты Шапиро-Уилка, чтобы проверить нормальность и, следовательно, обосновать применение соответствующих статистических тестов. Во всех случаях данные не были нормально распределены, поэтому были использованы независимые выборочные U-тесты Манна-Уитни. Для этих анализов мы предпочли представить средние ранги, а не медианы, так как на основании описательной статистики мы не предполагали идентичного распределения частот между группами выборки. В связи с большим объемом использованной выборки, мы также посчитали важным включить показатель размера эффекта, что было сделано с помощью пост-специальных тестов Эта-квадрат (η²) (121). Кроме того, для анализов, которые были близки к порогу значимости (0,05-0,06), мы провели перекрестную проверку результатов, применив параметрические методы к преобразованным данным. Это произошло в одном случае, и после преобразования путем отражения квадратного корня значимость теста не изменилась.

Кумулятивные различия в выживаемости колоний и маток между запасами оценивали с помощью моделей пропорциональных рисков Кокса со смешанными эффектами, используя пакеты R ‘coxme’ (122) и ‘multcomp’ (123). Предположение о пропорциональной опасности было проверено в каждом случае путем подтверждения независимости между остатками Шенфельда и временем, как для фиксированных факторов, так и для модели, с помощью пакета R ‘survival’ (124).

Для оценки различий между запасами по размеру популяции колоний и производству меда мы использовали обобщенные линейные модели (ОЛМ). В модели численности популяции в качестве временной точки отклика использовался февраль, а в модели производства меда — сентябрь. Выбор модели был основан на AIC и подтвержден с помощью всесторонних тестов, а также оценкой соотношения отклонений и степеней свободы.

Чтобы определить, как уровни варроа и титры вирусов влияют на выживаемость колоний, и как поголовье, обработка от клещей и маршрут миграции модулируют эти взаимодействия, мы использовали обобщенные линейные смешанные модели (GLMM) с повторными мерами. Для анализа уровня варроа и вирусного титра были созданы отдельные модели, чтобы максимизировать размер выборки и минимизировать сложность модели. Измерения уровня варроа повторяли в пяти временных точках, а измерения вирусного титра — в четырех. Выбор модели был основан на AIC, начиная с полной модели и взаимодействий. Во всех случаях соответствие модели проверялось путем оценки стандартизированных остатков.

При проверке предсказательной силы уровней Varroa и вирусных титров в отдельных временных точках мы использовали анализ относительного риска (RR). Все RR-анализы использовали «аутгруппу» для рассматриваемого эпидемиологического фактора, а затем сравнивали уровень смертности в этой группе с уровнем смертности в группе «воздействия», чтобы определить относительное изменение риска. Группа воздействия во всех анализах вирусных титров определялась как колонии с титрами ≥ 1 × 10⁷, что является предельным значением, указывающим на эпидемиологически тяжелую инфекцию (1,96,98), в то время как «аутгруппа» состояла из колоний с титрами < 1 × 10⁷. Группы воздействия в анализе уровня Varroa были стратифицированы, чтобы охватить диапазон уровней MPHB, присутствующих в колониях: ≥ 1-2,5, > 2,5-5, > 5-7,5, и > 7,5, а аутгруппа поддерживалась на уровне < 1, таким образом, обеспечивая градуированную меру изменения риска при все более тяжелом заражении. Значения RR, превышающие единицу, указывали на соответствующее увеличение риска смертности, в то время как значения, меньшие единицы, указывали на соответствующее снижение. В соответствии со стандартной интерпретацией, значения RR считались незначительными, если полученные 95% доверительные интервалы перекрывали значение единицы (120). Для значимых значений RR мы также рассчитали AR, чтобы продемонстрировать величину абсолютного эффекта, как это рекомендуется для бинарных эпидемиологических исходов (125). Распределения вероятностей значимого риска затем визуализировали с помощью параллельного Байесовского анализа RR, используя пакет R ‘brr’ (126). По всем тестам мы подтвердили необходимый объем выборки для обеспечения минимальной мощности (1 — β) 0,80 при альфа (α) 0,05, используя значения стандартного отклонения (σ) и средней разницы (δ) из данных. Все статистические анализы проводились в PASS (релиз v. 21.0.2), SPSS (релиз v. 25.0.0.0) и R (релиз v. 4.0.2) (127).

Выживание колоний и маток с течением времени

Мы использовали модели пропорциональных рисков Кокса со смешанными эффектами для оценки различий в выживаемости колоний и маток между двумя стадами в течение исследования. В качестве фиксированных факторов-предикторов использовались поголовье, обработка от клещей и маршрут миграции, а в качестве случайного фактора — идентификатор колонии. Поскольку отбор образцов проводился каждые два месяца, кумулятивная выживаемость обновлялась раз в два месяца в течение 10 месяцев.

Размер популяции и производство меда

В ГЛМ, оценивающей влияние поголовья на размер популяции колоний, в качестве гамма-переменной отклика с логарифмической связью использовалось количество выживших колоний в феврале, а в качестве предикторов фиксированных факторов — поголовье, обработка от клещей и маршрут миграции. GLM, оценивающая влияние поголовья на производство меда, использовала чистый вес меда в сентябре как гамма-переменную с логарифмической связью, поголовье и обработку от клещей как предикторы с фиксированным фактором, и опустила маршрут миграции, так как мед был извлечен до транспортировки.

Факторы, влияющие на выживаемость колоний

GLMM, оценивающий влияние уровня варроа на выживаемость колоний, использовал статус выживаемости колоний в феврале как биномиальную переменную отклика с пробит-связью, Varroa MPHB, поголовье, обработку от клещей и маршрут миграции как предикторы с фиксированным фактором, их двустороннее взаимодействие и ID колонии как случайный фактор. В GLMM, оценивающем влияние вирусных титров на выживаемость колоний, статус выживаемости колоний в феврале использовался как биномиальная переменная отклика с логитной связью, log₁₀ DWV-A титр, log10 DWV-B титр, log10 CBPV титр, log10 BQCV титр, поголовье, обработка клещей и маршрут миграции как фиксированные предикторы, их двусторонние взаимодействия и ID колонии как случайный фактор.

Предсказательные различия в уровнях варроа и титрах вирусов

Независимые выборочные U-тесты Манна-Уитни использовались для сравнения различий в средних уровнях варроа и вирусных титрах, основанных на исходе выживания колонии и поголовье, для временных точек июня и сентября.

Прогностическая сила уровней варроа и вирусных титров

Для количественной оценки относительной прогностической силы уровней варроа и вирусных титров в июне и сентябре при определении смертности колоний, а также вирусных титров при сопоставлении с уровнем варроа, мы использовали анализ RR и AR. Для визуализации распределений вероятности риска мы использовали байесовский анализ RR, используя реакцию Пуассона и приоритет Джеффриса в качестве объективного опорного приоритета.

Модель варроа

Для подгонки уравнений кривых к данным модели Варроа использовался анализ оценки кривых, а выбор кривой оценивался по полученным результатам модели, основанным на F-значениях.

Результаты

Выживаемость колоний и маток с течением времени

Процент выживаемости со временем был значительно выше в колониях Pol-line по сравнению с коммерческими колониями (модель пропорциональных рисков Кокса со смешанными эффектами, влияние поголовья: коэффициент опасности = 0,408, χ2 = 122,800, d.f. = 1, N = 366, P < 0,001, процент выживаемости _Pol-line = 59,884%, процент выживаемости_Commercial = 26. 042%; рис. 2a), и в колониях, получивших высокую обработку от клещей, по сравнению с низкой обработкой от клещей (модель пропорциональных рисков Кокса со смешанными эффектами, эффект обработки от клещей: коэффициент опасности = 0,516, χ2 = 25,192, d.f. = 1, N = 366, P < 0,001, процентная выживаемость_high = 57,051%, процентная выживаемостьal_low = 30,769%). Кроме того, кумулятивная выживаемость значительно различалась между способами миграции (модель пропорциональных рисков Кокса со смешанными эффектами, влияние маршрута миграции: χ2 = 10,854, d.f. = 2, N = 366, P = 0,004, процент выживания_{Стационарная} = 26,415%, процент выживания_{Калифорния} = 39,355%, процент выживания_{Миссисипи} = 50,000%). Эта разница объяснялась более высокой выживаемостью в Калифорнии (модель пропорциональных рисков Кокса со смешанными эффектами, тест Бенджамина-Хохберга: коэффициент опасности = 0,607, z = — 2,854, d.f. = 2, P = 0,004, процентная выживаемость_{Калифорния} = 39,355%, процентная выживаемость_{Стационарная} = 26. 415%), и Миссисипи (модель пропорциональных рисков Кокса со смешанными эффектами, тест Бенджамина-Хохберга: коэффициент опасности = 0,570, z = — 3,206, d.f. = 2, P = 0,002, процентная выживаемость_{Миссисипи} = 50,000%, процентная выживаемость_{Стационарная} = 26,415%), по сравнению с группой стационарной миграции. ID колонии имел значительный случайный эффект (случайный эффект колонии: дисперсия = 0,393, P < 0,001).

Процент выживания маток со временем также был значительно выше в колониях Pol-line по сравнению с коммерческими колониями (модель пропорциональных рисков Кокса со смешанными эффектами, влияние поголовья: коэффициент опасности = 0,355, χ2 = 234,388, d.f. = 1, N = 366, P < 0,001, процент сохранения маток_Pol-line = 88,350%, процент сохранения маток_Commercial = 66. 000%; рис. 2a), в незначительной степени независимые от обработки клещами (модель пропорциональных рисков Кокса со смешанными эффектами, влияние обработки клещами: χ2 = 3,603, d.f. = 1, N = 366, P = 0,058) и маршрута миграции (модель пропорциональных рисков Кокса со смешанными эффектами, влияние маршрута миграции: χ2 = 6,857, d.f. = 2, N = 366, P = 0,052). Идентификатор колонии снова имел значительный случайный эффект (случайный эффект колонии: дисперсия = 0,852, P < 0,001).

Размер популяции и производство меда

Размер популяции выживших колоний в феврале существенно не отличался между двумя стадами (GLM, эффект стада: χ21 = 2,272, N = 153, P = 0,131; рис. 2b) или обработками от клещей (GLM, эффект обработки от клещей: χ²1 = 0,735, N = 153, P = 0. 391), при этом маршрут миграции оказался единственным значимым фактором (GLM, эффект маршрута миграции: χ²2 = 8,935, N = 153, P = 0,011, среднее количество кадров_{Стационар} = 3,240, среднее количество кадров_{Миссисипи} = 4,597, среднее количество кадров_{Калифорния} = 5,614). В частности, эта разница объяснялась более высоким количеством рамок в выживших колониях калифорнийской миграционной группы по сравнению с колониями стационарной миграционной группы (GLM, пост-выборочный тест Холма-Бонферрони: t = 3,154, d.f. = 1, P = 0,005, среднее количество рамок_{Калифорния} = 5,614, среднее количество рамок_{Стационарная} = 3,240). Аналогичным образом, производство меда существенно не различалось между колониями (GLM, эффект колонии: χ²1 = 0,476, N = 278, P = 0,490; рис. 2c) или обработками от клещей (GLM, эффект обработки от клещей: χ²1 = 0,371, N = 278, P = 0,542).

Факторы, влияющие на выживаемость колоний

При изучении влияния уровня варроа на выживаемость в феврале значительно повлияли поголовье колоний, обработка от клещей, маршрут миграции и уровень варроа (рис. 3a,b) (табл. 1). Кроме того, значительные взаимодействия имели место между поголовьем и обработкой от клещей, а также поголовьем и уровнем варроа, что указывает на важность поголовья как модулирующего фактора (рис. 3b; табл. 1). Взаимодействие между обработкой от клещей и уровнем варроа не было значимым (табл. 1), что говорит о том, что хотя уровень варроа был снижен обработкой от клещей, это не изменило его основную связь с выживаемостью колоний.

Тип	Фактор	Статистика	d.f	P
Фиксированный эффект	Запасы	F = 49.233	1, 371	< 0.001
Фиксированный эффект	Обработка от клещей	F = 24.385	1, 384	< 0.001
Фиксированный эффект	Миграционный маршрут	F = 3.337	2, 373	0.037
Фиксированный эффект	Уровень варроа	F = 53.124	1, 373	< 0.001
Фиксированный эффект	поголовье × обработка от клещей	F = 17.502	1, 375	< 0.001
Фиксированный эффект	Запас × уровень варроа	F = 18.912	1, 360	< 0.001
Фиксированный эффект	Обработка от клещей × уровень варроа	F = 1.042	1, 430	0.308
Фиксированный коэффициент	Уровень варроа	Коэффициент = 0,123	373	< 0.001
Фиксированный коэффициент	Пол-линия × уровень варроа	Коэффициент = -0,095	360	< 0.001
Пост-хок тест Холма-Бонферрони	Обработка от клещей: низкий-высокий	t = 5.940	378	< 0.001
Пост-хок тест Холма-Бонферрони	Сток: коммерческий-Pol-line	t = 8.106	381	< 0.001
Пост-хок тест Холма-Бонферрони	Маршрут миграции: Стационарный — Миссисипи	t = 3.563	382	0.001
Пост-хок тест Холма-Бонферрони	Маршрут миграции: Калифорния-Миссисипи	t = 2.603	375	0.019
Пост-хок тест Холма-Бонферрони	Низкая обработка от клещей: коммерческий пол-лайн	t = 7.477	387	< 0.001
Пост-хок тест Холма-Бонферрони	Высокая обработка от клещей: коммерческий Pol-line	t = 1.401	386	0.162
Случайный эффект	Идентификатор колонии	Z = 9.608	—	< 0.001

В частности, при изучении влияния варроа, колонии значительно чаще погибали в течение исследования, если имели более высокий уровень варроа (средний показатель dead = 2,698 MPHB, средний показатель alive = 0,914 MPHB; рис. 3a,b), находились в группе с низким уровнем клещей (средний показатель survival_low = 0,212, средний показатель survival_high = 0. 645), были коммерческими (средняя выживаемость_Commercial = 0,138, средняя выживаемость_Pol-line = 0,666) и относились к группе стационарной или калифорнийской, а не миссисипской миграции (средняя выживаемость_Stationary = 0,184, средняя выживаемость_California = 0,315, средняя выживаемость_{Missississippi} = 0,541) (Таблица 1). Запас повлиял на влияние обработки от клещей на выживаемость, так как коммерческие колонии выиграли от дополнительной обработки от клещей в большей степени, чем колонии Pol-line, что подтверждается значительным снижением выживаемости при низком (средняя выживаемость_Commercial = 0,028, средняя выживаемость_Pol-line = 0,625), но не при высоком уровне обработки от клещей (средняя выживаемость_Commercial = 0,559, средняя выживаемость_Pol-line = 0,724) (Таблица 1). Поголовье также влияло на влияние уровня варроа на выживаемость, так как поголовье Pol-line значительно улучшало прогноз для любого данного уровня варроа (рис. 3b; табл. 1), что приводило к лучшим результатам выживаемости во всех группах заражения. ID колонии имел значительный случайный эффект (Таблица 1), что указывает на необъяснимую межколониальную изменчивость результатов выживания, как и следовало ожидать.

В отличие от уровня варроа, титры DWV-A, DWV-B, CBPV или BQCV не оказали существенного влияния на статус выживания колоний в феврале (рис. 3c,d; табл. 2). Кроме того, при анализе связи между титрами вирусов и выживаемостью не было выявлено влияния поголовья, обработки клещей или маршрута миграции, поскольку колонии отбирались на основе равномерного распределения поголовья и маршрутов миграции по результатам выживаемости, чтобы минимизировать любые сбивающие эффекты на динамику вирусов (Таблица 2). Как и в случае с данными по варроа, идентификатор колонии имел значительный случайный эффект (Таблица 2).

Тип	Фактор	Статистика	d.f	P
Фиксированный эффект	Запас	F = 1.352	1, 64	0.249
Фиксированный эффект	Обработка от клещей	F = 0.216	1, 67	0.644
Фиксированный эффект	Миграционный маршрут	F = 0.135	1, 56	0.715
Фиксированный эффект	титр DWV-A	F = 1.026	1, 66	0.315
Фиксированный эффект	титр DWV-B	F = 0.206	1, 85	0.651
Фиксированный эффект	Титр CBPV	F = 0.516	1, 80	0.475
Фиксированный эффект	Титры BQCV	F = 3.142	1, 69	0.081
Случайный эффект	Идентификатор колонии	Z = 4.185	—	< 0.001

Прогнозируемые различия в уровнях варроа и вирусных титрах

Уровень варроа был значительно выше в колониях, погибших в ходе исследования, при сравнении как в июне (средний ранг погибших = 204,520, средний ранг _живых = 144,590), так и в сентябре (средний ранг _{погибших} = 201,890, средний ранг _живых = 119,550) (рис. 4a; табл. 3). Кроме того, в колониях Pol-line уровень варроа был значительно ниже, чем в коммерческих колониях в июне (средний ранг Pol-line = 140,850, средний ранг коммерческих = 212,530) и сентябре (средний ранг _Pol-line = 106,440, средний ранг _{коммерческих} = 215,200) (рис. 4b; табл. 3), что указывает на прогностическую связь между уровнем варроа, поголовьем и выживаемостью колоний.

Реакция	Фактор	Статистика	N	η2	P
Выживание	Уровень варроа в июне	U = 20,811.000	152, 205	0.082	< 0.001
Выживание	Уровень варроа в сентябре	U = 19,755.000	155, 169	0.193	< 0.001
Запас	Уровень варроа в июне	U = 9,493.500	190, 167	0.120	< 0.001
Запас	Уровень варроа в сентябре	U = 4,308.000	167, 157	0.337	< 0.001

Титры вирусов в июне существенно не отличались между колониями, которые погибли или выжили, для DWV-A (средний ранг _{погибших} = 40,680, средний ранг _{выживших} = 46,790), DWV-B (средний ранг _{погибших} = 49,000, средний ранг _{выживших} = 41,800) или CBPV (средний ранг _{погибших} = 45,440, средний ранг _{выживших} = 43,940) (рис. 4c; табл. 4). Однако титры BQCV в июне были выше в колониях, погибших к февралю (средний ранг = 53,200, средний ранг _{выживших} = 39,280; рис. 4c; табл. 4). В сентябре не было значительных различий между колониями, которые погибли или выжили, ни по одному из измеренных вирусных титров (DWV-A: средний ранг _{погибших} = 43,000, средний ранг _{выживших} = 35,880; DWV-B: средний ранг _{погибших} = 44,840, средний ранг _{выживших} = 34,800; CBPV: средний ранг _{погибших} = 36,960, средний ранг _{выживших} = 39,400; BQCV: средний ранг _{погибших} = 41,140, средний _{ранг выживших} = 36,960; рис. 4d; табл. 4).

Реакция	Фактор	Статистика	N	η2	P
Выживание	Титр DWV-A в июне	U = 781.500	55, 33	0.013	0.277
	Титр DWV-B в июне	U = 1,056.000	55, 33	0.019	0.200
	Титр CBPV в июне	U = 938.500	55, 33	0.001	0.789
	Июньский титр BQCV	U = 1,194.500	55, 33	0.070	0.013
	Титр DWV-A в сентябре	U = 798.000	48, 28	0.021	0.175
	Сентябрь DWV-B титр	U = 849.500	48, 28	0.048	0.056
	√-x (сентябрьский титр DWV-B)	t = 1.951	48, 28	0.051	0.055
	Титры CBPV в сентябре	U = 629.000	48, 28	0.003	0.643
	Сентябрьские титры BQCV	U = 746.000	48, 28	0.008	0.426
Запас	Титр DWV-A в июне	U = 680.500	46, 42	0.065	0.017
	Титр DWV-B в июне	U = 753.000	46, 42	0.036	0.075
	Титр CBPV в июне	U = 954.500	46, 42	< 0.001	0.923
	Июньский титр BQCV	U = 780.500	46, 42	0.027	0.121
	Сентябрь DWV-A титр	U = 124.500	37, 39	0.506	< 0.001
	Титр DWV-B в сентябре	U = 131.500	37, 39	0.495	< 0.001
	Сентябрьские титры CBPV	U = 300.000	37, 39	0.252	< 0.001
	Сентябрьские титры BQCV	U = 597.000	37, 39	0.022	0.196

В июне колонии Pol-line показали значительно более низкие титры DWV-A, чем коммерческие колонии (средний ранг _Pol-line = 37.700, средний ранг_{коммерческий} = 50.710); однако, не было значительных различий между колониями в титрах DWV-B (средний ранг _Pol-line = 39. 430, средний ранг_{коммерческий} = 49,130), CBPV (средний ранг _Pol-line = 44,230, средний ранг_{коммерческий} = 44,750) или титрах BQCV (средний ранг _Pol-line = 40,080, средний ранг_{коммерческий} = 48,530) (рис. 4e; табл. 4). В сентябре титры вирусов были значительно ниже в колониях Pol-line для DWV-A (средний ранг_Pol-line = 23,190, средний ранг_{коммерческий} = 54,640), DWV-B (средний ранг_Pol-line = 23,370, средний ранг_{коммерческий}= 54,450) и CBPV (средний ранг_Pol-line = 27,690, средний ранг_{коммерческий} = 49,890) (рис. 4f; табл. 4), как и следовало ожидать, исходя из их корреляции с уровнем варроа (рис. S1a-c). Однако, напротив, сентябрьские титры BQCV не имели значительных различий между стадами (средний ранг _Pol-line = 35,310, средний ранг_{коммерческий} = 41,860; рис. 4f; табл. 4), и не показали корреляции с уровнем варроа (рис. S1d).

Прогностическая сила уровней варроа и вирусных титров

В июне уровень варроа значительно повышал относительный риск смертности к февралю, если колонии попадали в страты ≥ 1-2,5 (RR, точечная оценка = 1,528, CI_95% 1,256-1,859, AR = 0,238, P < 0,001; рис. 5) или > 2,5-5 (RR, точечная оценка = 1,669, CI_95% 1,330-2,059, AR = 0,303, P < 0,001; рис. 5) зараженных MPHB. Однако малое количество колоний с уровнем Varroa > 5 MPHB не позволило сделать значительный прогноз риска смертности выше этого уровня (рис. 5). В сентябре начальный порог риска был выше, так как уровень варроа начал значительно увеличивать относительный риск смертности только при > 2,5-5 (RR, точечная оценка = 1,594, CI_95% 1,066-2,384, AR = 0,202, P = 0,024; рис. 5). 5), причем риск последовательно возрастал при > 5-7,5 (RR, точечная оценка = 2,326, CI_95% 1,624-3,333, AR = 0,431, P < 0,001; рис. 5) и > 7,5 MPHB (RR, точечная оценка = 2,647, CI_95% 1,916-3,658, AR = 0,536, P < 0,001; рис. 5). Примечательно, что ни один другой фактор существенно не влиял на относительный риск смертности. Список всех факторов и распределение вероятностей риска см. на (рис. 5).

Модель варроа

Связь между уровнем варроа в сентябре и смертностью объяснялась логарифмической функцией в коммерческих колониях (кривая оценки, F1,24 = 12,851, R² = 0,349, P < 0,001; рис. 6) и линейной функцией в колониях Pol-line (кривая оценки, F1,10 = 18,311, R² = 0,647, P = 0,002; рис. 6). При пороговом значении 3 MPHB высокая смертность была предсказана как для коммерческих (кривая оценки, процент смертностиCommercial = 72,618%; рис. 6), так и для колоний Pol-line (кривая оценки, процент смертностиPol-line = 52,210%; рис. 6), хотя в последнем случае она была значительно ниже. Интересно, что две функции пересекаются при 10 MPHB, поскольку, несмотря на значительную смертность, некоторые коммерческие колонии смогли выжить при очень высоких уровнях заражения, что отражается в снижении объяснительной способности коммерческой модели (рис. 6).

Обсуждение

Наши результаты показывают, что в условиях миграционного пчеловодства у пчел Pol-line значительно повышается выживаемость колоний и маток по сравнению со стандартным коммерческим итальянским поголовьем (рис. 2a). Кроме того, продуктивность колоний Pol-line эквивалентна продуктивности коммерческих пчел, о чем свидетельствует сопоставимое производство меда и размер популяции среди выживших колоний (рис. 2b,c). Это улучшение производительности связано с более низким уровнем варроа (рис. 3b, 4b) и одновременным снижением титров DWV-A, DWV-B и CBPV в сентябре (рис. 3d, 4f). Примечательно, что мультипликативный эффект уровня варроа на выживаемость колоний значительно слабее в колониях Pol-line, что приводит к лучшим прогнозам выживаемости во всех, кроме самого высокого уровня инвазии (рис. 3b, 6; табл. 1). Предположительным объяснением этого является постоянное удаление зараженного расплода рабочими Pol-line, что характерно для поведения VSH (48,54), и, таким образом, последующее снижение динамики роста варроа и передачи вируса. Следовательно, колонии Pol-line, получившие только одну обработку от клещей, продемонстрировали уровень выживаемости, эмпирически эквивалентный коммерческим колониям, получившим две обработки. Действительно, польза от дополнительных обработок от клещей была в целом выше для коммерческих колоний, так как в большинстве колоний Pol-line естественным образом поддерживался уровень варроа ниже рекомендуемого порога обработки (114). Это имеет значительные последствия как для повышения выживаемости колоний, так и для снижения растущей потребности в использовании акарицидов (27).

Насколько нам известно, эти данные представляют собой наиболее полную характеристику улучшенной функциональности устойчивых к варроа медоносных пчел на сегодняшний день. Кроме того, поскольку матки, использованные в данном исследовании, были спарены через насыщение трутнями, а полученные колонии сохранили генетику устойчивости к клещу, наши результаты непосредственно применимы для промышленного использования (128). Следует отметить, что настоящее исследование включает данные только за один год, поэтому потребуется дальнейшая работа по проверке производных Pol-line перед коммерческим внедрением. Однако обнадеживает то, что коммерческие полевые испытания необработанных пчел «Hilo «(60) — поголовья, полученного от Pol-line в 2015 году — показали превосходную выживаемость в зимний период, хотя и несколько более низкое производство меда, по сравнению с обработанными акарицидами коммерческими итальянскими и итальянскими × карниоланскими колониями (RGD, личное сообщение) (129). Примечательно, что анализ частоты аллелей помещает пчел Pol-line и Hilo в отдельный генетический кластер, значительно отличающийся от других коммерческих и исследовательских колоний, используемых в настоящее время в США (60). Таким образом, наши результаты служат многообещающей эмпирической характеристикой потенциала выведенных запасов, устойчивых к варроа, для фундаментального улучшения современного пчеловодства.

Сравнительно высокий уровень смертности, наблюдаемый в данном исследовании (7,8) , является следствием двух особенностей экспериментальной схемы. Для того чтобы оценить истинную устойчивость колоний к Варроа, мы обязательно включали минимальную группу обработки от клещей для каждого поголовья, что привело к значительному снижению общего использования акарицидов по сравнению со стандартными коммерческими операциями. Таким образом, уровень варроа был высоким, а его присутствие повсеместным, что имеет соответствующие последствия для распространения клещей и смертности до второй акарицидной обработки — особенно в коммерческих колониях. Кроме того, промышленные опылительные работы создают различные нагрузки на колонии из-за интенсивного транзита, экстремальных климатических условий, высокой плотности колоний и колебаний местных условий (57,63,130). Повышенный уровень смертности, наблюдаемый в группе стационарной миграции, иллюстрирует эту стохастичность, поскольку зимой колонии испытывали голод, что привело к серьезным потерям (Таблица 1). Такой результат подчеркивает важность режима управления и местных стрессовых факторов для влияния на исход колоний (130) , и является актуальным при рассмотрении распределенных систем опыления. Примечательно, что наши результаты не противоречат другим коммерческим сценариям, в которых давление варроа является высоким, что очень часто бывает (82,131). Таким образом, настоящее исследование представляет собой контролируемое, но строгое испытание в условиях крупномасштабной коммерческой деятельности при значительном бремени паразитов. В случае валидации поголовья такая установка предпочтительна, поскольку мы хотели определить постоянную эффективность признаков VSH, полученных в результате ауткросса Pol-line (54) , а для этого необходимо значительное давление варроа.

Результаты нашей оценки варроа, DWV-A, DWV-B, CBPV и BQCV и их влияния на смертность колоний были весьма примечательными. В соответствии с предыдущими исследованиями, уровень варроа был самым сильным предиктором смертности колоний для обеих популяций, расходясь уже в июне и достигая пика в сентябре (рис. 3, 4a,b; табл. 1, 3). Интересно, что влияние как раннего сезона, так и предзимнего уровня варроа на последующую смертность говорит о том, что текущие пороги обработки могут быть недостаточно консервативными, чтобы сдержать значительные потери (рис. 3a, 6). Действительно, дополнительные пробит-модели подтверждают это утверждение для настоящих данных (рис. S2). В отличие от этого, титры DWV-A, DWV-B, CBPV и BQCV не были значимыми независимыми детерминантами выживания колоний (рис. 3c,d, 4c; табл. 2). Единственное значимое вирусное различие наблюдалось в случае июньских титров BQCV, которые были выше в колониях, которые в дальнейшем погибли, однако размер эффекта был промежуточным (Таблица 4) и не был перенесен в многомерный анализ, что указывает на отсутствие причинно-следственной связи.

Хотя сентябрьские титры DWV-A, DWV-B и BQCV были в среднем численно выше в колониях, которые погибли, эта разница не была значительной при контроле ковариаций или при использовании в качестве прогностической меры, а уровни CBPV были фактически ниже в колониях, которые позже погибли (рис. 3c, 4d; табл. 2). DWV-A и DWV-B продемонстрировали в целом схожие модели в зависимости от результатов выживания и запасов, за исключением декабря, когда DWV-B снизился в коммерческих колониях в большей степени, чем DWV-A (рис. 3d). Титры DWV-A, DWV-B и CBPV были значительно снижены в колониях Pol-line в сентябре, а в случае DWV-A — также в июне; однако BQCV не показал различий между запасами (рис. 4e,f). Это может объясняться более низким уровнем варроа в колониях Pol-line, причем наибольшие различия проявляются в сентябре (рис. 3b, 4b; табл. 3), и, следовательно, одновременным снижением уровней прямо и косвенно связанных с варроа патогенов (рис. S1a-c). Примечательно, что такая тенденция указывает на отсутствие значительной ассоциации BQCV с варроа и подтверждается отсутствием прямой корреляции с сентябрьскими уровнями варроа (рис. S1d). Этот вывод важен для этиологического понимания BQCV, поскольку он предполагает, что, несмотря на потенциальную репликацию в клещах (102) и временную корреляцию с инвазией варроа (22) , эмпирические доказательства эпидемиологически значимого пути передачи малочисленны.

В отличие от DWV-A и DWV-B, ассоциация CBPV с варроа была сравнительно слабой (рис. S1c). Таким образом, возможно, что механизмы вирусной резистентности, специфичные для конкретного поголовья (132,133) , или дестабилизация иммунной функции в результате нагрузки варроа (134) , сыграли свою роль в различиях, наблюдаемых между поголовьями (рис. 4f). Генотипическое взаимодействие хозяина и патогена представляет интерес для селекции медоносных пчел (38) и может заслуживать дальнейшего изучения как у пчел Pol-line, так и у других сортов, устойчивых к клещам (66,94). Однако будет важно отделить такие эффекты от оппортунизма патогенов, особенно в тех случаях, когда запасы обеспечивают широкое улучшение здоровья колоний.

В целом, наши данные показывают, что титры DWV-A, DWV-B и, в меньшей степени, CBPV отслеживают уровень варроа, но являются сравнительно более слабыми предикторами выживаемости колоний, чем сам варроа. Анализ относительного и объяснимого риска позволил глубже понять эти результаты. Когда рассматривались только колонии, сгруппированные по сопоставимым уровням варроа, различные титры DWW-A, DWV-B, CBPV и BQCV никогда не являлись значимыми аддитивными предикторами риска смертности (рис. 5). Более того, единственные относительные показатели риска, оказывающие значительное влияние на смертность, относились только к уровню варроа (рис. 5).

Эти результаты имеют значение для понимания DWV-A и DWV-B как возбудителей гибели колоний. Наши данные показывают, что варроа, как фактор, постоянно оказывает наибольшее влияние на смертность, вероятно, поскольку он включает в себя как последствия механического повреждения, через кормление, так и передачу связанных с ним патогенов. Хотя DWV и другие вирусы коррелируют с гибелью колоний и, как было показано, убивают пчел и оказывают на них сублетальное воздействие (65,135), похоже, что когда эти эффекты отделяются от их корреляции с уровнем варроа, их предсказательная сила значительно снижается. Важнейшим методологическим и диагностическим следствием этого является то, что, как минимум, по сравнению с протестированными здесь вирусами, варроа обладает значительной аддитивной силой в прогнозировании исхода колонии. Этот аддитивный эффект может быть объяснен повреждением жирового тела хозяина (12,21) или дополнительной скрытой передачей патогена, не учитываемой мерами, направленными на одну цель (17), и подтверждается предыдущими экспериментальными данными (76,83,101,136). Хотя такая тенденция не отменяет важности вирусов, переносимых варроа, она предполагает, что вред, причиняемый самим кормлением варроа, хронически недооценивается, и подразумевает, что устоявшаяся парадигма потери колоний, опосредованной варроа, может быть неверной. Примечательно, что сравнительно меньший эффект вирусов, наблюдаемый в нашем исследовании, может быть объяснен одновременным использованием контролируемого подхода смешанного моделирования, коммерчески реалистичного размера выборки и измерения только вирусного титра, а не вирусной распространенности. Последний аспект очень важен, поскольку потенциально скрытая этиология DWV-инфекций (137) требует связи между фактическими титрами вирусов и смертностью, особенно если необходимо провести причинные ассоциации. Действительно, вирусы, переносимые варроа, коррелируют с уровнем заражения варроа, но эта связь не является монотонной (1,22), и поэтому будущие исследования, связывающие эти патогены со смертностью, должны включать или контролировать уровень варроа, если они хотят получить значимые выводы.

При рассмотрении потерь колоний, относительный вес DWV или любого из вирусов, переносимых варроа, никогда не был эффективно изолирован от эффекта кормления варроа как такового (12). Следует отметить, что, как показано здесь, другие наборы данных обнаружили различное причинное влияние вирусных титров на здоровье пчел (138) и смертность колоний (81,136). Кроме того, на сегодняшний день ни одно исследование не смогло отделить корреляцию смертности от вируса Варроа от корреляции смертности от Варроа на уровне колоний. Исходя из этого, несколько непонятно, почему переносчики вирусов считаются главными, если не единственными, возбудителями потерь колоний, опосредованных варроа (38,139), особенно учитывая, что способ питания паразита требует повреждения органов (12). Действительно, в некоторых случаях только уровень варроа предсказывает смертность, в то время как титры вирусов этого не делают (83,136). Хорошо известно, что путем удаления клещей варроа или контроля над ними можно косвенно управлять связанными с ними вирусами (111,112), поэтому есть практическая польза в том, чтобы сосредоточиться на варроа, поскольку она является основным эпидемиологическим фактором. Это хорошо иллюстрируется тем, что, хотя поведение VSH может фактически привести к распространению DWV через кукольный каннибализм в лаборатории (140), оно, тем не менее, снижает титры DWV на уровне колоний в поле (рис. 3д, 4ф), поскольку смягчает более мощный путь передачи — кормление варроа. Наряду с недавними достижениями в понимании этиологии варроа (12) , растущий объем исследований также указывает на то, что сами клещи являются не безобидными переносчиками, а инвазивными эпизоотическими паразитами, которые наносят серьезный ущерб своим наивным хозяевам посредством разрушительной липофагии. Наряду с настоящими результатами, это требует серьезной переоценки относительной важности варроа и передаваемых им вирусов в причинении вреда на уровне колонии. Многофакторные подходы к пониманию динамики варроа-вирусов, особенно в отношении весовых коэффициентов факторов, будут важны в будущих эпидемиологических исследованиях, поскольку, возможно, опосредованная варроа гибель колоний все еще требует существенного механистического выяснения.

Наши результаты показывают, что устойчивые к варроа медоносные пчелы имеют большой потенциал для снижения потерь колоний в коммерческом пчеловодстве. Инкорпорирование признака VSH в коммерческое поголовье, по-видимому, обеспечивает как устойчивость к варроа, так и благоприятную продуктивность колоний. Более того, мы обнаружили, что титры DWV-A, DWV-B, CBPV и BQCV хуже предсказывают смертность колоний по сравнению с уровнем варроа. Это означает, что для разработки эффективных методов снижения смертности потребуется дальнейшее изучение относительного ущерба, наносимого питанием варроа как таковым, а также патогенов, передающихся при этом процессе. В целом, наши данные свидетельствуют о том, что стратегия борьбы с варроа, ориентированная на варроа и включающая в себя внедрение устойчивых запасов, в настоящее время представляет собой наиболее эффективное, устойчивое и технологически доступное решение проблемы глобальных потерь медоносных пчел.

Дополнительная информация

Дополнительные рисунки.
Дополнительные таблицы.

Ссылки

1.Martin, S. J. et al. Global honey bee viral landscape altered by a parasitic mite. Science (80–) 336, 1304–1306 (2012).CAS ADS Google Scholar 

2.Neumann, P. & Carreck, N. L. Honey bee colony losses. J. Apic. Res. 49, 1–6 (2010).Google Scholar 

3.Ratnieks, F. L. W. & Carreck, N. L. Clarity on honey bee collapse?. Science (80–) 327, 152–153 (2010).CAS ADS Google Scholar 

4.Potts, S. G. et al. Safeguarding pollinators and their values to human well-being. Nature 540, 220–229 (2016).CAS PubMed ADS Google Scholar 

5.Vanengelsdorp, D. et al. A national survey of managed honey bee 2010–11 winter colony losses in the USA: Results from the bee informed partnership. J. Apic. Res. 51, 115–124 (2012).Google Scholar 

6.Lee, K. V. et al. A national survey of managed honey bee 2013–2014 annual colony losses in the USA. Apidologie 46, 292–305 (2015).
Google Scholar 

7.Kulhanek, K. et al. A national survey of managed honey bee 2015–2016 annual colony losses in the USA. J. Apic. Res. 56, 328–340 (2017).
Google Scholar 

8.Bruckner, S. et al. Honey bee colony losses 2018–2019: Preliminary results. Bee Inf. (2019) https://beeinformed.org/results/honey-bee-colony-losses-2017-2018-preliminary-results/.

9.Nazzi, F. & Pennacchio, F. Disentangling multiple interactions in the hive ecosystem. Trends Parasitol. 30, 556–561 (2014).PubMed Google Scholar 

10.Francis, R. M., Nielsen, S. L. & Kryger, P. Varroa-virus interaction in collapsing honey bee colonies. PLoS One 8, 25 (2013).Google Scholar 

11.Rosenkranz, P., Aumeier, P. & Ziegelmann, B. Biology and control of Varroa destructor. J. Invertebr. Pathol. 103, S96–S119 (2010).
PubMed Google Scholar 

12.Ramsey, S. D. et al. Varroa destructor feeds primarily on honey bee fat body tissue and not hemolymph. Proc. Natl. Acad. Sci. USA 116, 1792–1801 (2019).CAS PubMed PubMed Central Google Scholar 

13.Martin, S. J. The role of Varroa and viral pathogens in the collapse of honeybee colonies: A modelling approach. J. Appl. Ecol. 38, 1082–1093 (2001).Google Scholar 

14.Bowen-Walker, P. L. & Gunn, A. The effect of the ectoparasitic mite, Varroa destructor on adult worker honeybee (Apis mellifera) emergence weights, water, protein, carbohydrate, and lipid levels. Entomol. Exp. Appl. 101, 207–217 (2001).CAS Google Scholar 

15.Yang, X. & Cox-Foster, D. L. Impact of an ectoparasite on the immunity and pathology of an invertebrate: Evidence for host immunosuppression and viral amplification. Proc. Natl. Acad. Sci. USA 102, 7470–7475 (2005).CAS PubMed PubMed Central ADS Google Scholar 

16.Blanken, L. J., van Langevelde, F. & van Dooremalen, C. Interaction between Varroa destructor and imidacloprid reduces flight capacity of honeybees. Proc. R. Soc. B Biol. Sci. 282, 20151738 (2015).Google Scholar 

17.Traynor, K. S. et al. Varroa destructor: A complex parasite, crippling honey bees worldwide. Trends Parasitol. 36, 592–606 (2020).CAS PubMed Google Scholar 

18.Boecking, O. & Genersch, E. Varroosis—the ongoing crisis in bee keeping. J. Verbraucherschutz Leb. 3, 221–228 (2008).Google Scholar 

19.McMenamin, A. J. & Genersch, E. Honey bee colony losses and associated viruses. Curr. Opin. Insect Sci. 8, 121–129 (2015).
PubMed Google Scholar 

20.Levin, S. et al. New viruses from the Ectoparasite Mite Varroa destructor infesting Apis mellifera and Apis cerana. Viruses 11, 94 (2019).CAS PubMed Central Google Scholar 

21.Arrese, E. L. & Soulages, J. L. Insect fat body: Energy, metabolism, and regulation. Annu. Rev. Entomol. 55, 207–225 (2010).
CAS PubMed PubMed Central Google Scholar 

22.Mondet, F., de Miranda, J. R., Kretzschmar, A., Le Conte, Y. & Mercer, A. R. On the front line: Quantitative virus dynamics in honeybee (Apis mellifera L.) colonies along a new expansion front of the parasite Varroa destructor. PLoS Pathog. 10, 25 (2014).Google Scholar 

23.Wilfert, L. et al. Honeybee disease: Deformed wing virus is a recent global epidemic in honeybees driven by Varroa mites. Science (80–) 351, 594–597 (2016).CAS ADS Google Scholar 

24.Santamaria, J. et al. Evidence of Varroa-mediated deformed wing virus spillover in Hawaii. J. Invertebr. Pathol. 151, 126–130 (2017).PubMed Google Scholar 

25.Martin, S. J., Elzen, P. J. & Rubink, W. R. Effect of acaricide resistance on reproductive ability of the honey bee mite Varroa destructor. Exp. Appl. Acarol. 27, 195–207 (2002).CAS PubMed Google Scholar 

26.Evans, J. D. & Cook, S. C. Genetics and physiology of Varroa mites. Curr. Opin. Insect Sci. 26, 130–135 (2018).PubMed Google Scholar 

27.Rinkevich, F. D. Detection of amitraz resistance and reduced treatment efficacy in the Varroa Mite, Varroa destructor, within commercial beekeeping operations. PLoS One 15, 1–12 (2020).Google Scholar 

28.Martin, S. J. Acaricide (pyrethroid) resistance in Varroa destructor. Bee World 85, 67–69 (2004).Google Scholar 

29.Grozinger, C. M. & Flenniken, M. L. Bee viruses: Ecology, pathogenicity, and impacts. Annu. Rev. Entomol. 64, 205–226 (2019).CAS PubMed Google Scholar 

30.Owen, R. Role of human action in the spread of honey bee (Hymenoptera: Apidae) pathogens. J. Econ. Entomol. 110, 797–801 (2017).PubMed Google Scholar 

31.Rinderer, T. E., Harris, J. W., Hunt, G. J. & De Guzman, L. I. Breeding for resistance to Varroa destructor in North America. Apidologie 41, 409–424 (2010).Google Scholar 

32.Büchler, R., Berg, S. & Le Conte, Y. Breeding for resistance to Varroa destructor in Europe. Apidologie 41, 393–408 (2010).Google Scholar 

33.Eliash, N. & Mikheyev, A. Varroa mite evolution: A neglected aspect of worldwide bee collapses?. Curr. Opin. Insect Sci. 39, 21–26 (2020).PubMed Google Scholar 

34.Beaurepaire, A. L., Krieger, K. J. & Moritz, R. F. A. Seasonal cycle of inbreeding and recombination of the parasitic mite Varroa destructor in honeybee colonies and its implications for the selection of acaricide resistance. Infect. Genet. Evol. 50, 49–54 (2017).CAS PubMed Google Scholar 

35.Locke, B. Natural Varroa mite-surviving Apis mellifera honeybee populations. Apidologie 47, 467–482 (2016).Google Scholar 

36.Oddie, M. A. Y., Dahle, B. & Neumann, P. Norwegian honey bees surviving Varroa destructor mite infestations by means of natural selection. PeerJ 5, e3956 (2017).PubMed PubMed Central Google Scholar 

37.Oddie, M. et al. Rapid parallel evolution overcomes global honey bee parasite. Sci. Rep. 8, 1–9 (2018).CAS ADS Google Scholar 

38.Mondet, F. et al. Honey bee survival mechanisms against the parasite Varroa destructor: A systematic review of phenotypic and genomic research efforts. Int. J. Parasitol. 50, 433–447 (2020).PubMed Google Scholar 

39.Seeley, T. D. Life-history traits of wild honey bee colonies living in forests around Ithaca, NY, USA. Apidologie 48, 743–754 (2017).Google Scholar 

40.Loftus, J. C., Smith, M. L. & Seeley, T. D. How honey bee colonies survive in the wild: Testing the importance of small nests and frequent swarming. PLoS One 11, 1–11 (2016).Google Scholar 

41.Seeley, T. D. & Smith, M. L. Crowding honeybee colonies in apiaries can increase their vulnerability to the deadly ectoparasite Varroa destructor. Apidologie 46, 716–727 (2015).Google Scholar 

42.Dynes, T. L., Berry, J. A., Delaplane, K. S., Brosi, B. J. & De Roode, J. C. Reduced density and visually complex apiaries reduce parasite load and promote honey production and overwintering survival in honey bees. PLoS One 14, 1–16 (2019).Google Scholar 

43.Di Pasquale, G. et al. Influence of pollen nutrition on honey bee health: Do pollen quality and diversity matter?. PLoS One 8, 1–13 (2013).
Google Scholar 

44.Guarna, M. M. et al. Peptide biomarkers used for the selective breeding of a complex polygenic trait in honey bees. Sci. Rep. 7, 1–10 (2017).CAS Google Scholar 

45.Grozinger, C. M. & Robinson, G. E. The power and promise of applying genomics to honey bee health. Curr. Opin. Insect Sci. 10, 124–132 (2015).PubMed PubMed Central Google Scholar 

46.Techer, M. A. et al. Divergent evolutionary trajectories following speciation in two ectoparasitic honey bee mites. Commun. Biol. 2, 25 (2019).Google Scholar 

47.Harbo, J. R. & Harris, J. W. Responses to Varroa by honey bees with different levels of Varroa sensitive hygiene. J. Apic. Res. 48, 156–161 (2009).Google Scholar 

48.Tsuruda, J. M., Harris, J. W., Bourgeois, L., Danka, R. G. & Hunt, G. J. High-resolution linkage analyses to identify genes that influence varroa sensitive hygiene behavior in honey bees. PLoS One 7, 5 (2012).Google Scholar 

49.Spivak, M. & Danka, R. G. Perspectives on hygienic behavior in Apis mellifera and other social insects. Apidologie https://doi.org/10.1007/s13592-020-00784-z (2020).Article Google Scholar 

50.Rinderer, T. E. et al. The release of ARS Russian honey bees. Am. Bee J. 140, 305–308 (2000).Google Scholar 

51.Guzman, L. I. D. E., Rinderer, T. E. & Frake, A. M. Growth of Varroa destructor (Acari: Varroidae) Populations in Russian Honey Bee (Hymenoptera: Apidae) Colonies. Ecol. Popul. Biol. 100, 187–195 (2007).Google Scholar 

52.Morfin, N., Given, K., Evans, M., Guzman-Novoa, E. & Hunt, G. J. Grooming behavior and gene expression of the Indiana “mite-biter” honey bee stock. Apidologie 51, 1–9 (2019).Google Scholar 

53.Pritchard, D. J. Grooming by honey bees as a component of varroa resistant behavior. J. Apic. Res. 55, 38–48 (2016).Google Scholar 

54.Danka, R. G., Harris, J. W. & Dodds, G. E. Selection of VSH-derived “Pol-line” honey bees and evaluation of their Varroa-resistance characteristics. Apidologie 47, 483–490 (2016).Google Scholar 

55.Guichard, M. et al. Three decades of selecting honey bees that survive infestations by the parasitic mite Varroa destructor: Outcomes, limitations and strategy. Genet. Sel. Evol. https://doi.org/10.1186/s12711-020-00591-1 (2020).Article PubMed PubMed Central Google Scholar 

56.Ibrahim, A., Reuter, G. S. & Spivak, M. Field trial of honey bee colonies bred for mechanisms of resistance against Varroa destructor. Apidologie 38, 67–76 (2007).Google Scholar 

57.Danka, R. G. et al. Functionality of Varroa-resistant honey bees (Hymenoptera: Apidae) when used in migratory beekeeping for crop pollination. J. Econ. Entomol. 105, 313–321 (2012).PubMed Google Scholar 

58.Nielsen, C. & Lund, M. The concept of business model scalability. Soc. Sci. Res. Netw. 20, 1–20 (2015).Google Scholar 

59.Bixby, M. et al. A bio-economic case study of Canadian Honey bee (Hymenoptera: Apidae) colonies: Marker-assisted selection (MAS) in queen breeding affects beekeeper profits. J. Econ. Entomol. 110, 816–825 (2017).PubMed PubMed Central Google Scholar 

60.Saelao, P. et al. Genome-wide patterns of differentiation within and among US commercial honey bee stocks. BMC Genom. 21, 704 (2020).Google Scholar 

61.Harbo, J. R. & Harris, J. W. Resistance to Varroa destructor (Mesostigmata: Varroidae) when mite-resistant queen honey bees (Hymenoptera: Apidae) were free-mated with unselected drones. J. Econ. Entomol. 94, 1319–1323 (2001).CAS PubMed Google Scholar 

62.Free, J. B. Insect Pollination of Crops (Academic Press, 1993).Google Scholar 

63.Delaplane, K. S., Mayer, D. R. & Mayer, D. F. Crop Pollination by Bees (Cabi, 2000).Google Scholar 

64.Beck, H. E. et al. Present and future köppen-geiger climate classification maps at 1-km resolution. Sci. Data 5, 1–12 (2018).ADS Google Scholar 

65.Benaets, K. et al. Covert deformed wing virus infections have long-term deleterious effects on honeybee foraging and survival. Proc. R. Soc. B Biol. Sci. 284, 25 (2017).Google Scholar 

66.Thaduri, S., Stephan, J. G., de Miranda, J. R. & Locke, B. Disentangling host-parasite-pathogen interactions in a varroa-resistant honeybee population reveals virus tolerance as an independent, naturally adapted survival mechanism. Sci. Rep. 9, 1–10 (2019).CAS Google Scholar 

67.Kielmanowicz, M. G. et al. Prospective large-scale field study generates predictive model identifying major contributors to colony losses. PLoS Pathog. 11, 1–20 (2015).CAS Google Scholar 

68.Ryabov, E. V. et al. A virulent strain of deformed wing virus (DWV) of honeybees (Apis mellifera) prevails after Varroa destructor-mediated, or in vitro, transmission. PLoS Pathog. 10, 1004230 (2014).Google Scholar 

69.Ongus, J. R. et al. Complete sequence of a picorna-like virus of the genus Iflavirus replicating in the mite Varroa destructor. J. Gen. Virol. 85, 3747–3755 (2004).CAS PubMed Google Scholar 

70.Dalmon, A. et al. Evidence for positive selection and recombination hotspots in deformed wing virus (DWV). Sci. Rep. 7, 1–12 (2017).Google Scholar 

71.Dainat, B., Evans, J. D., Chen, Y. P., Gauthier, L. & Neumanna, P. Dead or alive: Deformed wing virus and varroa destructor reduce the life span of winter honeybees. Appl. Environ. Microbiol. 78, 981–987 (2012).CAS PubMed PubMed Central ADS Google Scholar 

72.Moore, J. et al. Recombinants between Deformed wing virus and Varroa destructor virus-1 may prevail in Varroa destructor-infested honeybee colonies. J. Gen. Virol. 92, 156–161 (2011).CAS PubMed Google Scholar 

73.Kevill, J. L. et al. DWV-A lethal to honey bees (Apis mellifera): A colony level survey of DWV variants (A, B, and C) in England, Wales, and 32 States across the US. Viruses 11, 25 (2019).Google Scholar 

74.Tehel, A. et al. The two prevalent genotypes of an emerging infectious disease, deformed wing virus, cause equally low pupal mortality and equally high wing deformities in host honey bees. Viruses 11, 25 (2019).Google Scholar 

75.Posada-Florez, F. et al. Deformed wing virus type A, a major honey bee pathogen, is vectored by the mite Varroa destructor in a non-propagative manner. Sci. Rep. 9, 12445 (2019).PubMed PubMed Central ADS Google Scholar 

76.Natsopoulou, M. E. et al. The virulent, emerging genotype B of Deformed wing virus is closely linked to overwinter honeybee worker loss. Sci. Rep. 7, 5242 (2017).PubMed PubMed Central ADS Google Scholar 

77.Ryabov, E. V. et al. Dynamic evolution in the key honey bee pathogen deformed wing virus: Novel insights into virulence and competition using reverse genetics. PLoS Biol. 17, e3000502 (2019).MathSciNet CAS PubMed PubMed Central Google Scholar 

78.Brettell, L. E. & Martin, S. J. Oldest Varroa tolerant honey bee population provides insight into the origins of the global decline of honey bees. Sci. Rep. 7, 1–7 (2017).Google Scholar 

79.Mordecai, G. J. et al. Superinfection exclusion and the long-term survival of honey bees in Varroa-infested colonies. ISME J. 10, 1182–1191 (2016).CAS PubMed Google Scholar 

80.Berthoud, H., Imdorf, A., Haueter, M., Radloff, S. & Neumann, P. Virus infections and winter losses of honey bee colonies (Apis mellifera). J. Apic. Res. 49, 60–65 (2010).Google Scholar 

81.Cirkovic, D. et al. Honey bee viruses in Serbian colonies of different strength. PeerJ 6, e5887 (2018).PubMed PubMed Central Google Scholar 

82.Rinderer, T. E. et al. Functionality of Varroa-resistant honey bees (Hymenoptera: Apidae) when used for Western US Honey production and almond pollination. J. Econ. Entomol. 107, 523–530 (2014).PubMed Google Scholar 

83.Smart, M., Pettis, J., Rice, N., Browning, Z. & Spivak, M. Linking measures of colony and individual honey bee health to survival among apiaries exposed to varying agricultural land use. PLoS One 11, 0152685 (2016).Google Scholar 

84.Seitz, K. et al. A molecular clone of Chronic Bee Paralysis Virus (CBPV) causes mortality in honey bee pupae (Apis mellifera). Sci. Rep. 9, 16274 (2019).PubMed PubMed Central ADS Google Scholar 

85.Ribière, M., Olivier, V. & Blanchard, P. Chronic bee paralysis: A disease and a virus like no other?. J. Invertebr. Pathol. 103, S120–S131 (2010).PubMed Google Scholar 

86.Gisder, S. & Genersch, E. Viruses of commercialized insect pollinators. J. Invertebr. Pathol. 147, 51–59 (2017).PubMed Google Scholar 

87.Faucon, J.-P. et al. Honey bee winter mortality in France in 1999 and 2000. Bee World 83, 14–23 (2002).Google Scholar 

88.Antúnez, K., Alessandro, B. D., Corbella, E. & Zunino, P. Detection of Chronic bee paralysis virus and Acute bee paralysis virus in Uruguayan honeybees. J. Invertebr. Pathol. 90, 69–72 (2005).PubMed Google Scholar 

89.Amiri, E., Meixner, M., Büchler, R. & Kryger, P. Chronic bee paralysis virus in honeybee queens: Evaluating susceptibility and infection routes. Viruses 6, 25 (2014).Google Scholar 

90.Celle, O. et al. Detection of Chronic bee paralysis virus (CBPV) genome and its replicative RNA form in various hosts and possible ways of spread. Virus Res. 133, 280–284 (2008).CAS PubMed Google Scholar 

91.Budge, G. E. et al. Chronic bee paralysis as a serious emerging threat to honey bees. Nat. Commun. 11, 2164 (2020).
CAS PubMed PubMed Central ADS Google Scholar 

92.Ball, B. V. & Allen, M. F. The prevalence of pathogens in honey bee (Apis mellifera) colonies infested with the parasitic mite Varroa jacobsoni. Ann. Appl. Biol. 113, 237–244 (1988).Google Scholar 

93.Le Conte, Y., Ellis, M. & Ritter, W. Varroa mites and honey bee health: Can Varroa explain part of the colony losses?. Apidologie 41, 353–363 (2010).Google Scholar 

94.De Guzman, L. I., Simone-Finstrom, M., Frake, A. M. & Tokarz, P. Comb irradiation has limited, interactive effects on colony performance or pathogens in bees, Varroa destructor and wax based on two honey bee stocks. Insects 10, 15 (2019).PubMed Central Google Scholar 

95.Wang, S. et al. Occurrence of multiple honeybee viruses in the ectoparasitic mites Varroa spp. in Apis cerana colonies. J. Invertebr. Pathol. 166, 107225 (2019).PubMed Google Scholar 

96.Blanchard, P. et al. Evaluation of a real-time two-step RT-PCR assay for quantitation of Chronic bee paralysis virus (CBPV) genome in experimentally-infected bee tissues and in life stages of a symptomatic colony. J. Virol. Methods 141, 7–13 (2007).CAS PubMed Google Scholar 

97.Spurny, R. et al. Virion structure of black queen cell virus, a common honeybee pathogen. J. Virol. 91, e02100-e2116 (2017).CAS PubMed PubMed Central Google Scholar 

98.Al Naggar, Y. & Paxton, R. J. Mode of transmission determines the virulence of black queen cell virus in adult honey bees, posing a future threat to bees and apiculture. Viruses 12, 535 (2020).CAS PubMed Central Google Scholar 

99.Murray, E. A. et al. Viral transmission in honey bees and native bees, supported by a global black queen cell virus phylogeny. Environ. Microbiol. 21, 972–983 (2019).PubMed Google Scholar 

100.Desai, S. D. & Currie, R. W. Effects of wintering environment and parasite-pathogen interactions on honey bee colony loss in north temperate regions. PLoS One 11, e0159615 (2016).PubMed PubMed Central Google Scholar 

101.Dainat, B., Evans, J. D., Chen, Y. P., Gauthier, L. & Neumann, P. Predictive markers of honey bee colony collapse. PLoS One 7, 25 (2012).
Google Scholar 

102.Sabahi, Q., Morfin, N., Nehzati-Paghaleh, G. & Guzman-Novoa, E. Detection and replication of deformed wing virus and black queen cell virus in parasitic mites, Varroa destructor, from Iranian honey bee (Apis mellifera) colonies. J. Apic. Res. 59, 211–217 (2020).Google Scholar 

103.Loope, K. J., Baty, J. W., Lester, P. J. & Wilson Rankin, E. E. Pathogen shifts in a honeybee predator following the arrival of the Varroa mite. Proc. R. Soc. B Biol. Sci. 286, 20182499 (2019).CAS Google Scholar 

104.Hamiduzzaman, M. M. et al. Differential responses of Africanized and European honey bees (Apis mellifera) to viral replication following mechanical transmission or Varroa destructor parasitism. J. Invertebr. Pathol. 126, 12–20 (2015).PubMed Google Scholar 

105.Locke, B., Forsgren, E., Fries, I. & de Miranda, J. R. Acaricide treatment affects viral dynamics in Varroa destructor-infested honey bee colonies via both host physiology and mite control. Appl. Environ. Microbiol. 78, 227–235 (2012).CAS PubMed PubMed Central ADS Google Scholar 

106.Tentcheva, D. et al. Prevalence and seasonal variations of six bee viruses in Apis mellifera L. and Varroa destructor mite populations in France. Appl. Environ. Microbiol. 70, 7185–7191 (2004).CAS PubMed PubMed Central ADS Google Scholar 

107.Barroso-Arévalo, S. et al. High load of deformed wing virus and Varroa destructor infestation are related to weakness of honey bee colonies in Southern Spain. Front. Microbiol. 10, 1331 (2019).PubMed PubMed Central Google Scholar 

108.Wu, Y., Liu, Q., Weiss, B., Kaltenpoth, M. & Kadowaki, T. Honey bee suppresses the parasitic mite Vitellogenin by antimicrobial peptide. Front. Microbiol. 11, 25 (2020).Google Scholar 

109.Nguyen, B. K. et al. Effects of honey bee virus prevalence, Varroa destructor load and queen condition on honey bee colony survival over the winter in Belgium. J. Apic. Res. 50, 195–202 (2011).Google Scholar 

110.Kevan, P. G., Hannan, M., Ostiguy, N. & Guzman-novoa, E. A summary of the Varroa-virus disease complex in honey bees. Am. Bee J. 146, 694–697 (2006).Google Scholar 

111.Zhao, Y. et al. The dynamics of deformed wing virus concentration and host defensive gene expression after Varroa mite parasitism in honey bees, Apis mellifera. Insects 10, 16 (2019).PubMed Central Google Scholar 

112.Dietemann, V., Neumann, P. & Ellis, J. D. The COLOSS BEEBOOK, Volume I: Standard methods for Apis mellifera research. J. Apic. Res. 52, 25 (2013).Google Scholar 

113.Rangel, J. & Seeley, T. D. Colony fissioning in honey bees: Size and significance of the swarm fraction. Insectes Soc. 59, 453–462 (2012).
Google Scholar 

114.Dietemann, V. et al. Standard methods for varroa research. J. Apic. Res. 52, 1–54 (2013).Google Scholar 

115.Simone-Finstrom, M., Aronstein, K., Goblirsch, M., Rinkevich, F. & de Guzman, L. Gamma irradiation inactivates honey bee fungal, microsporidian, and viral pathogens and parasites. J. Invertebr. Pathol. 153, 57–64 (2018).PubMed Google Scholar 

116.Boncristiani, H. et al. Direct effect of acaricides on pathogen loads and gene expression levels in honey bees Apis mellifera. J. Insect Physiol. 58, 613–620 (2012).CAS PubMed Google Scholar 

117.Ryabov, E. V. et al. Recent spread of Varroa destructor virus-1, a honey bee pathogen, in the United States. Sci. Rep. 7, 1–10 (2017).
CAS Google Scholar 

118.Brunetto, M. R., Colombatto, P. & Bonino, F. Bio-mathematical models of viral dynamics to tailor antiviral therapy in chronic viral hepatitis. World J. Gastroenterol. 15, 531–537 (2009).CAS PubMed PubMed Central Google Scholar 

119.Tarpy, D. R., Vanengelsdorp, D. & Pettis, J. S. Genetic diversity affects colony survivorship in commercial honey bee colonies. Naturwissenschaften 100, 723–728 (2013).CAS PubMed ADS Google Scholar 

120.Porta, M. A Dictionary of Epidemiology (Oxford University Press, 2014).Google Scholar 

121.Lenhard, W. & Lenhard, A. Calculation of effect sizes. Psychometrica. http://www.psychometrica.de/effect_size.html (2016). https://doi.org/10.13140/RG.2.2.17823.92329.

122.Therneau, T. M. coxme: Mixed effects cox models. R package version 2.2-16. (2020).

123.Hothorn, T., Bretz, F. & Westfall, P. Simultaneous inference in general parametric models. Biometr. J. 50, 346–363 (2008).MathSciNet MATH Google Scholar 

124.Therneau, T. M. A package for survival analysis in R. R package version 3.2-7. (2020).

125.Moher, D. et al. CONSORT 2010 explanation and elaboration: Updated guidelines for reporting parallel group randomised trials. Int. J. Surg. 10, 28–55 (2012).PubMed Google Scholar 

126.Laurent, S. brr: Bayesian inference on the ratio of two poisson rates. R package version 1.0.0. (2015).

127.R Core Team. R: A Language and Environment for Statistical Computing (R Foundation for Statistical Computing, 2020).Google Scholar 

128.Hellmich, R. L. & Waller, G. D. Preparing for Africanized honey bees: Evaluating control in mating apiaries. Am. Bee J. 130, 537–542 (1990).Google Scholar 

129.Hilo Bees. Hilo bee field trials. https://www.hilobees.com/hilo-bee-field-trials.html (2021).

130.Simone-Finstrom, M. et al. Migratory management and environmental conditions affect lifespan and oxidative stress in honey bees. Sci. Rep. 6, 32023 (2016).CAS PubMed PubMed Central ADS Google Scholar 

131.Kefuss, J., Vanpoucke, J., Bolt, M. & Kefuss, C. Selection for resistance to Varroa destructor under commercial beekeeping conditions. J. Apic. Res. 54, 563–576 (2015).Google Scholar 

132.Khongphinitbunjong, K. et al. Responses of Varroa-resistant honey bees (Apis mellifera L.) to Deformed wing virus. J. Asia. Pac. Entomol. 19, 921–927 (2016).Google Scholar 

133.Khongphinitbunjong, K. et al. Differential viral levels and immune gene expression in three stocks of Apis mellifera induced by different numbers of Varroa destructor. J. Insect Physiol. 72, 28–34 (2015).CAS PubMed Google Scholar 

134.Annoscia, D. et al. Haemolymph removal by Varroa mite destabilizes the dynamical interaction between immune effectors and virus in bees, as predicted by Volterra’s model. Proc. R. Soc. B Biol. Sci. 286, 20190331 (2019).CAS Google Scholar 

135.Wells, T. et al. Flight performance of actively foraging honey bees is reduced by a common pathogen. Environ. Microbiol. Rep. 8, 728–737 (2016).PubMed PubMed Central Google Scholar 

136.Clermont, A. et al. Virus status, Varroa levels, and survival of 20 managed honey bee colonies monitored in luxembourg between the summer of 2011 and the spring of 2013. J. Apic. Sci. 59, 59–73 (2015).Google Scholar 

137.McMenamin, A. J. & Flenniken, M. L. Recently identified bee viruses and their impact on bee pollinators. Curr. Opin. Insect Sci. 26, 120–129 (2018).PubMed Google Scholar 

138.Monchanin, C. et al. Hazard of a neonicotinoid insecticide on the homing flight of the honeybee depends on climatic conditions and Varroa infestation. Chemosphere 224, 360–368 (2019).CAS PubMed ADS Google Scholar 

139.Locke, B., Forsgren, E. & De Miranda, J. R. Increased tolerance and resistance to virus infections: A possible factor in the survival of Varroa destructor-resistant honey bees (Apis mellifera). PLoS One 9, e99998 (2014).PubMed PubMed Central ADS Google Scholar 

140.Posada-Florez, F. et al. Pupal cannibalism by worker honey bees contributes to the spread of deformed wing virus. Sci. Rep. 11, 8989 (2021).CAS PubMed PubMed Central ADS Google Scholar