Synchrony of Bird Migration with Global Dispersal of Avian Influenza Reveals Exposed Bird Orders

Abstract

Highly pathogenic avian influenza virus (HPAIV) A H5, particularly clade 2.3.4.4, has caused worldwide outbreaks in domestic poultry, occasional spillover to humans, and increasing deaths of diverse species of wild birds since 2014. Wild bird migration is currently acknowledged as an important ecological process contributing to the global dispersal of HPAIV H5. However, this mechanism has not been quantified using bird movement data from different species, and the timing and location of exposure of different species is unclear. We sought to explore these questions through phylodynamic analyses based on empirical data of bird movement tracking and virus genome sequences of clade 2.3.4.4 and 2.3.2.1. First, we demonstrate that seasonal bird migration can explain salient features of the global dispersal of clade 2.3.4.4. Second, we detect synchrony between the seasonality of bird annual cycle phases and virus lineage movements. We reveal the differing exposed bird orders at geographical origins and destinations of HPAIV H5 clade 2.3.4.4 lineage movements, including relatively under-discussed orders. Our study provides a phylodynamic framework that links the bird movement ecology and genomic epidemiology of avian influenza; it highlights the importance of integrating bird behavior and life history in avian influenza studies.

Fig. 1. The evolutionary history of HPAIV H5 clade 2.3.4.4 and 2.3.2.1.

A Time-scaled maximum clade credibility (MCC) tree of clade 2.3.4.4 (2010–2023). Sequences sampled from 2010 to 2017 are colored light blue, and those sampled from 2018 to 2023 are darker blue. The phylogeographic reconstructions for the sequences sampled from 2010 to 2017 (B) and for those from 2018 to 2023 (C), shown as time-scaled MCC trees with location annotations to summarize the reconstructions. Abbreviations of locations: NChina (North China), SChina (South China), SEA (South-East Asia), CA (Central Asia), WA (Western Asia), LA (Latin America). Supplementary Dataset 3 lists countries in each aggregated region.

Fig. 2. Contributions of predictors to worldwide diffusion of H5N1 clade 2.3.2.1 and 2.3.4.4.

The virus dispersal is inferred from HA genes by GLM-extended Bayesian phylogeographic inference with heterogeneous evolutionary processes through time. Predictors in the model are shown in Fig. S6. The conditional effect size of each predictor (on a log scale) is presented as mean values with credible intervals. The phylogeographic GLM estimates were obtained from n = 1845, 1844, and 1163 viral sequences sampled over 12, 10, and 11 locations for clade 2.3.4.4 (2018–2023), 2.3.4.4 (2010–2017), and 2.3.2.1, respectively. The light blue, darker blue, and red colors indicate clade 2.3.4.4 (2018–2023), 2.3.4.4 (2010–2017), and 2.3.2.1, respectively.

Fig. 3. The synchrony of virus lineage movements between regions and bird distribution probability at the regions.

A Probability density distribution of the virus lineage migration throughout the year, between locations summarized from the discrete trait phylogeography of HPAIV H5 clade 2.3.4.4 and the Markov jump counts (Section “Discrete trait phylogeography of HPAIVs and counts ofvirus lineage migration”). X axis: Virus lineage migration dates in a year; Y axis: origin region - destination region of the virus lineage migration. The width of the violins represents the virus lineage migration probability density. Non-breeding (blue), migration (yellow) and breeding (red) bird annual cycle phases in general are shown in the south-north migration direction and in the north-south migration direction. Boxes around bird photos show the statistically significant correlation of virus lineage movements and bird order distribution at origin, destination or both regions. We used the following images of bird species from the Macaulay Library at the Cornell Lab of Ornithology: Pandion haliaetus haliaetus (ML523577271) and Cygnus olor (ML72775261). Abbreviations of locations are the same as in Fig. 1. B Schematic diagram of cross-correlation analyses of virus lineage movement between two locations (O: origin, D: destination) and the bird distribution probability at each location. Block bootstrapping was used to calculate confidence intervals and two-tailed p values. See details in Methods “Discrete trait phylogeography of HPAIVs and counts ofvirus lineage migration”.