Genetic Characterization of Highly Pathogenic Avian Influenza A(H5N8) Virus in Pakistani Live Bird Markets Reveals Rapid Diversification of Clade 2.3.4.4b Viruses

The highly pathogenic (HPAI) avian influenza A(H5N1) viruses have undergone reassortment with multiple non-N1-subtype neuraminidase genes since 2008, leading to the emergence of H5Nx viruses. H5Nx viruses established themselves quickly in birds and disseminated from China to Africa, the Middle East, Europe and North America. Multiple genetic clades have successively evolved through frequent mutations and reassortment, posing a continuous threat to domestic poultry and causing substantial economic losses. Live bird markets are recognized as major sources of avian-to-human infection and for the emergence of zoonotic influenza. In Pakistan, the A(H5N1) virus was first reported in domestic birds in 2007; however, avian influenza surveillance is limited and there is a lack of knowledge on the evolution and transmission of the A(H5) virus in the country. We collected oropharyngeal swabs from domestic poultry and environmental samples from six different live bird markets during 2018–2019. We detected and sequenced HPAI A(H5N8) viruses from two chickens, one quail and one environmental sample in two markets. Temporal phylogenetics indicated that all novel HPAI A(H5N8) viruses belonged to clade 2.3.4.4b, with all eight genes of Pakistan A(H5N8) viruses most closely related to 2017 Saudi Arabia A(H5N8) viruses, which were likely introduced via cross-border transmission from neighboring regions approximately three months prior to virus detection into domestic poultry. Our data further revealed that clade 2.3.4.4b viruses underwent rapid lineage expansion in 2017 and acquired significant amino acid mutations, including mutations associated with increased haemagglutinin affinity to human α-2,6 receptors, prior to the first human A(H5N8) infection in Russian poultry workers in 2020. These results highlight the need for systematic avian influenza surveillance in live bird markets in Pakistan to monitor for potential A(H5Nx) variants that may arise from poultry populations.


Surveillance and Sample Collection
During September 2018-March 2019, avian influenza virus (AIV) surveillance was conducted in six LBMs in three districts (Gujranwala, Lahore and Sheikhupura) of Punjab Province, Pakistan ( Figure 1). These LBMs sell fresh domesticated poultry meat and eggs for consumption. A wide range of bird species are housed, such as backyard and commercially farmed chickens, quails, ducks and pigeons. For each LBM, sample collections were conducted once a month. We collected 449 oropharyngeal swabs from domestic poultry (including chicken and quail) and 696 environmental samples (including chopping board surface, cages, cage drinking water, weighing scale, sewage water and water used for meat processing). All procedures were approved by the Ethical Review Committee at the University of Veterinary and Animal Sciences Lahore in Pakistan (DR/381). Pooled samples (up to n = 5) were cultured in 9-day-old embryonated chicken eggs for 48 h and amnioticallantoic fluids were tested for influenza A virus by hemagglutination assay (HA) using 1% chicken red blood cells [35]. Individual samples from positive pools were subsequently cultured and tested by hemagglutination inhibition (HI) assay for the presence of H9 virus [28]. Hemagglutination-assay-positive but HI-assay-negative samples were then directly submitted for full genome sequencing to identify non-H9 viruses.
Viruses 2021, 13, x FOR PEER REVIEW 3 of 17 surface, cages, cage drinking water, weighing scale, sewage water and water used for meat processing). All procedures were approved by the Ethical Review Committee at the University of Veterinary and Animal Sciences Lahore in Pakistan (DR/381). Pooled samples (up to n = 5) were cultured in 9-day-old embryonated chicken eggs for 48 h and amniotic-allantoic fluids were tested for influenza A virus by hemagglutination assay (HA) using 1% chicken red blood cells [35]. Individual samples from positive pools were subsequently cultured and tested by hemagglutination inhibition (HI) assay for the presence of H9 virus [28]. Hemagglutination-assay-positive but HI-assay-negative samples were then directly submitted for full genome sequencing to identify non-H9 viruses.

RNA Extraction and Next Generation Sequencing
Total viral RNA was extracted from the cultured isolates using the QIAamp Viral RNA Mini Kit (Qiagen, Hilden, Germany) following the manufacturer's instructions. Fullgenome sequencing of influenza A virus was performed by next-generation sequencing. Briefly, cDNA synthesis and RT-PCR amplification were performed using the SuperScript III Platinum One-step qRT-PCR kit (Thermo Fisher Scientific, Waltham, MA, USA) using the following three primers as previously described [36]: Opti1-F1-5′ GTTACGCGCCAG-CAAAAGCAGG; Opti1-F2-5′ GTTACGCGCCAGCGAAAGCAGG; Opti1-R1-5′ GTTAC-GCGCCAGTAGAAACAAGG (Integrated DNA Technologies Pte. Ltd., Singapore). The cycling conditions were set as follows: 1 cycle of 55

RNA Extraction and Next Generation Sequencing
Total viral RNA was extracted from the cultured isolates using the QIAamp Viral RNA Mini Kit (Qiagen, Hilden, Germany) following the manufacturer's instructions. Full-genome sequencing of influenza A virus was performed by next-generation sequencing. Briefly, cDNA synthesis and RT-PCR amplification were performed using the Su-perScript III Platinum One-step qRT-PCR kit (Thermo Fisher Scientific, Waltham, MA, USA) using the following three primers as previously described [36]: Opti1-F1-5 GT-TACGCGCCAGCAAAAGCAGG; Opti1-F2-5 GTTACGCGCCAGCGAAAGCAGG; Opti1-R1-5 GTTACGCGCCAGTAGAAACAAGG (Integrated DNA Technologies Pte. Ltd., Singapore). The cycling conditions were set as follows: extension of 68 • C for 10 min. The PCR products were measured using the Qubit dsDNA HS Assay Kit (Thermo Fisher Scientific) and diluted to 1 ng/µL of DNA. The libraries were prepared using the Nextera XT DNA Library Preparation Kit (Illumina, Inc. San Diego, CA, USA) and the quality was checked using the Agilent 2100 Bioanalyzer (Agilent Technologies, Palo Alto, CA, USA). The pooled libraries were then run on an Illumina MiSeq (2 × 250 bp, San Diego, CA, USA). The short NGS reads were quality-checked using UGENE v39 [37] and adaptor sequences were trimmed using Trimmomatic v0.39 [38]. De novo assembly of the reads was performed using SPAdes genome assembler v3.13 [39] and the contigs were blasted using BLAST v2.2.30 [40] against a local influenza virus database downloaded from NCBI. The reads were then mapped to the reference influenza genome. We recovered complete genomes of four novel A(H5N8) viruses and the sequences were deposited in the NCBI GenBank database (see Table 1 for accession numbers).

Temporal Phylogenetic Analysis and Molecular Characterization
All available influenza A(H5Nx) sequences were downloaded from NCBI and GI-SAID databases (as of 4 March 2021). Large dataset phylogenies were reconstructed using FastTree in Geneious v7.1.9 (Biomatters Ltd., Auckland, New Zealand) and the datasets were then sub-sampled to include clade 2.3.4.4a-h viruses. For each gene segment, dated phylogenies and changes in relative genetic diversity were estimated using the uncorrelated lognormal relaxed clock with the Gaussian Markov random field (GMRF) Bayesian skyride model and the SRD06 codon position model in BEAST v1.10.4 [41]. At least four independent runs of 100 million generations were used. The convergence of runs was checked using Tracer v.1.7.1 [42] after excluding burn-in values to ensure the effective sampling size values for all parameters were >200. The runs were combined using Log-Combiner and the resulting maximum clade creditability (MCC) trees were summarized using TreeAnnotator. Ancestral amino acid mutations at the nodes were determined for each gene segment phylogeny using the treesub program [43].

Prevalence of AIV in Pakistani Live Bird Markets
During 2018-2019, we detected 77 (17.1%) of 449 samples from chickens and quails as being positive for influenza A virus based on hemagglutination inhibition assay and NGS sequencing (Table 2), including 4 with A(H5N8) virus and 73 with A(H9N2) virus. The positivity rate of influenza A virus in chickens was higher in two LBMs (31.6% in Tollinton and 31.7% in Shahdara) in the Lahore district compared to LBMs in Sheikhupura and Sharaqpur (4.9-13.0%). Quail samples could only be collected from Tollinton market and showed high proportions of influenza A virus (28.6%) ( Table 2). We also tested 696 environmental samples, of which 6 (0.9%) were positive for influenza A virus. Four environmental samples were detected in Lahore LBMs: three samples from chopping boards and one sample from sewage water. Two environmental samples from Sheikhupura LBMs were also positive: one sample from cage drinking water and one sample from a meat weighing scale.

Discussion
Live bird markets (LBMs) are associated with avian influenza transmission and outbreaks, providing environmental sources for virus persistence [57] and amplification, thereby increasing the risk of zoonotic infection. The multispecies composition of LBMs also poses risks for virus mutation and adaptation to different hosts. Through active influenza surveillance of LBMs in Pakistan, we observed varying levels of AIV positivity in domestic poultry, with LBMs in Lahore district having greater AIV prevalence compared to the LBMs in Sheikhupura and Gujranwala districts. Previous surveillance studies have reported the prevalence of AIV in Lahore, where low pathogenic A(H9N2) viruses have frequently been identified in domestic and commercial poultry [27,28,30,31,58]. The circulation of A(H9N2) has led to the generation of novel genotypes, increasing the risks of AIV infection in poultry and related occupational workers such as butchers and vaccinators [29,59,60].
Between 2018 and 2019, we detected and identified four avian HPAI A(H5N8) viruses in chickens and quails from LBMs in Lahore. We also detected environmental contamination with A(H5N8) from the same markets. Our phylogenetic analyses of Pakistan A(H5N8) viruses showed close genetic and phylogenetic similarities, all belonging to clade 2.3.4.4b viruses, showing that LBMs play a significant role in transmission of the virus among domestic birds. We showed that all eight gene segments of Pakistan HPAI A(H5N8) viruses were most closely related to viruses from Saudi Arabia, suggesting they were most likely introduced to Pakistan through the movement of live poultry from neighboring regions, while the similar TMRCA dates across all eight genes indicate that the viruses are likely derived from a single avian source that was introduced approximately three months before virus detection. The presence of viruses on chopping boards suggests that environmental contamination may occur through slaughtering or processing of infected live poultry. A previous study in China indicated that the highest H7 positivity rates came from chopping boards [21]. The prevalence rates of contamination from Pakistani LBMs in this study are also comparable to those of Bangladeshi LBMs [25]. Whilst the A(H9N2) virus has been frequently documented in Pakistan, no H5 or H7 viruses have been reported in poultry, although low levels of H5 antibodies (6.9%) have been reported in backyard poultry collected in 2009 [58].
The emergence of HPAI A(H5N8) clade 2.3.4.4b viruses in Pakistan in chickens is likely associated with poultry production. In Pakistan, poultry production has been growing significantly since the 1960s, with an increase of 126% in total meat production and an increase of 71% in total egg production during 2000-2010, with the main chicken exports being made to neighboring regions, including Afghanistan, Iran and Turkey [59]. A study examining the phylogeographic dynamics of A(H9N2) in Asia highlighted the migration rates from Pakistan to a number of countries, including India, Iran, Israel, Saudi Arabia and the United Arab Emirates [15]. Their data also indicated that live poultry trade and production are important drivers causing the spatial spread of A(H9N2) virus. In addition to domestic poultry, the HAPI A(H5N8) virus has been reported to the World Organization for Animal Health (OIE) by the Pakistan National Laboratory due to H5-positive wild birds (including mallard and swan) and macaws being found in Lahore Zoo, although no genetic data are available from those viruses, meaning cannot rule out introduction via wild birds. More surveillance data is needed to understand the transmission of A(H5N8) among domestic and wild birds in Pakistan and neighboring regions.
In this study, we also showed that clade 2.3.4.4b viruses exhibited markedly rapid diversification that peaked in 2017, with the acquisition of non-synonymous amino acid mutations on most viral genes. Some of these mutations are associated with mammalian adaptation, suggesting domestic poultry and wild birds are playing a significant role in virus evolution and potential cross-species transmission. Environmental contamination in LBMs also provides a persistent source for poultry infection, and effective control strategies to mitigate environmental contamination with AIV in LBMs are needed. The recent infections of poultry workers from Russia and of gray seals [60] indicate a broader host range for A(H5N8) viruses, highlighting their pandemic potential. At present, A(H5N8) clade 2.3.4.4 viruses have evolved into at least 8 distinct lineages, suggesting A(H5N8) variants are emerging at an unprecedented rate and pose a public health concern [61,62]. No sustained human-to-human transmission of A(H5) has been observed; however, the continual evolution and dissemination of A(H5N8) virus across multiple countries highlights the need for active influenza surveillance in domestic poultry and wild birds to detect emerging genetic variants and to assess their pandemic potential. Data Availability Statement: All the avian influenza sequences generated as a part of this study are publicly available in GenBank, individual accession numbers are provided in Table 1.