Detection and Phylogenetic Analysis of Infectious Bursal Disease Virus Based on both Genome Segments A and B Isolated from Backyard Poultry Punjab, Pakistan

In recent years, the re-emergence of virulent strains of infectious bursal disease virus (IBDV) has resulted in substantial economic losses in Pakistan despite mass and intense vaccination regimens. Infectious bursal disease (IBD) is one of the most important immunosuppressive diseases of the poultry and has been a constraint on sustainable food security around the globe, including Pakistan. This disease damages the bursa of Fabricius (BF) , which causes immunosuppression in birds. A total of 50 tissue samples from backyard chicken flocks presenting suspected symptoms were collected during February 2020 and March 2021 for genetic characterization, followed by phylogentic analysis. A total of 4 isolates were sequenced based on partial VP1 and hyper-variable region of the VP2 genes simultaneously. According to phylogenetic analysis, the study isolates genotype A3B3 were identified as predominant strains in country backyard poultry. Concerning the identity matrix analysis of VP1 representative part, study isolates shared (87-88% nt; 95-96% aa) identity with the vvIBDV and (97-98% nt; 98-99% aa) with non-vvIBDV. While VP2 revealed (99-100% nt; 99-100% aa) identity with previously reported Pakistan vvIBDV strains and (91-95% nt; 91-95% aa) with non vvIBDV. Amino acid alignment analysis of VP1 revealed that current IBDVs have three characteristic aa residues of vvIBDV (287-A, 508-K, and 511-S) and four characteristic aa residues of non-vvIBDV (146-E, 147-G, 242-D and 390-L). While VP2 gene sequence alignment revealed eight characteristic aa residues of vvIBDV (222-A, 242-I, 253-Q, 256-I, 279-D, 284-A, 299-S, and 330-S and a distinct aa 384-I). Based on phylogeny, this is the first identification of IBDV segment reassortants having segment A of A3 (very virulent) and segment B of B3 (early Australian-like) genogroups reported in Pakistan backyard poultry. Further study is required to determine the pathogenicity of the IBDV reassortant and development of new policies for IBDV intervention in the country.


INTRODUCTION
B efore 1960s backyard poultry production was the sole poultry meat source in the country, with flock size of about 20-50 birds of different breeds and ages. It was a low input business without veterinary care facilities (Ahmed et al., 2021). Infectious bursal disease (IBD) is a highly contagious and immunosuppressive disease of young chickens that causes considerable economic losses to the poultry industry globally (Qin and Zheng, 2017;Hussain et al., 2020;Maqbool et al., 2020). The disease is caused by the infectious bursal disease virus (IBDV), a non-enveloped dsRNA virus, which belongs to the genus Avibirnavirus of family Birnaviridae with two distinct serotypes (Müller et al., 2003). Serotype I is pathogenic and induces pathological lesions in chickens, whereas serotype II is non-pathogenic, primarily isolated from turkeys. Pathogenic IBDV is classified into four subtypes, i.e., attenuated, classical virulent, antigenic variants, and very virulent IBDV (Hussain et al., 2020). IBDV targets immature B lymphocytes (Zahid et al., 2018) in the lymphoid organ (bursa of Fabricius), leading to age-

O n l i n e F i r s t A r t i c l e
dependent immunosuppression exposing birds to other opportunistic pathogens. Very virulent IBDV (vvIBDV) has been reported in different countries with more than 70% chicken mortality, posing a severe threat to the poultry industry (Dey et al., 2019). IBDV has a bi-segmented genome (segment A and B) of approximately 6.0 kb. Segment A is 3.2 kb long, with two overlapping open reading frames (ORF) encoding four proteins. ORF-1 encodes longer precursor polyprotein, cleaved by autoproteolysis and yields capsid protein VP2, nucleoprotein VP3 and protease VP4. In comparison, ORF-2 encodes only a nonstructural VP5 protein, a significant determinant of viral virulence and host protective antigen (Luque et al., 2009). The VP2 protein is the primary protective antigen containing specific epitopes responsible for inducing neutralizing antibody responses (Reddy et al., 2017). Amino acid substitutions in the VP2 hypervariable region can lead to variations in antibody recognition, antigenicity, virulence, and tissue tropism (Mato et al., 2020). Segment B is 2.8kb long, and it contains only one ORF which encodes RNA-dependent RNA polymerase VP1 protein, involved in virus replication and evolution (Ye et al., 2018). Both VP1 and VP2 proteins have been used to study pathogenicity, virulence, phylogenetic relationship and molecular characterizations .
IBDV may undergo segment reassortment due to the genome's bi-segmented structure. It has been proposed that vvIBDV arose through segment reassortment between segment A of circulating IBDV strains and segment B from an unknown IBDV strain (Hon et al., 2006). Natural reassortant IBDV strains have been reported frequently in recent years from various parts of the world, including reassortment between serotypes 1 and 2 (Soubies et al., 2017;Stoute et al., 2019) and reassortment between different genotypes of serotype 1 IBDVs (Mato et al., 2020).
In 1980 a vvIBDV outbreak was first recorded with high mortality in Punjab Pakistan (Khan et al., 1988) and since then continues to pose a threat by various strains all over the country (Lone et al., 2009;Zahoor et al., 2011). Despite extensive vaccination programs, IBD still prevails at the poultry farms posing significant economic losses in Pakistan, which emphasizes the characterization of field circulating IBDV strains to understand the causes of vaccination failure, hence incorporating local viruses in vaccine development (Shafqat et al., 2017;Sajid et al., 2021). Traditionally, nucleotide sequencing of the segment A coded VP2 hypervariable region is used for genotypic characterization of IBDV strains. Both segments should be characterized to detect genetic variations, influencing strain pathogenicity (Gao et al., 2018). Both genome segments contribute to the virulence of IBDV and segment reassortment plays a role in virus evolution, which is not possible to be evaluated by targeting VP2 gene only. Since VP1 and VP2 proteins are involved in virus evolution, there is still insufficient information to predict the actual dynamics of IBDV in Pakistan. Recently, some authors have characterized the chicken isolated IBDV VP1 and VP2 genes (Hussain et al., , 2020, but the literature on backyard poultry is silent. Therefore, this study was designed for IBDV genotyping and phylogenetic analysis covering both segments from indigenous breed of backyard poultry. The study's primary objective was to detect and genotype IBDV strains circulating in country backyard poultry and add a new piece of knowledge about its characterization to better implement control measures.

Ethical approval
The present work was approved by the ethical review committee of the University of Veterinary and Animal Sciences, Lahore, Pakistan (Permit Number: DR/944).

Study area
The current study was conducted at the Department of Pathology, University of Veterinary and Animal Sciences Lahore, Pakistan.

Clinical samples
A total number of 50 tissue samples (bursa of Fabricius) were collected from backyard chicken flocks for genome extraction from February 2020 to March 2021. The dead bird subjected to virus isolation had typical hemorrhagic lesions on skeletal muscles and mucosal surface of the bursa of Fabricius. The bursal tissues were ground to make homogenate in phosphate buffered saline (pH 7.2) containing penicillin and streptomycin. The homogenate was finally centrifuged at 5 000×g for 5 min at 4°C. Lastly, the upper aquous layer was collected for further detection.

Molecular identification
The virus genome was extracted according to the protocol as described by Qi et al. (2015). Viral RNA was extracted from the bursa homogenate using QIAamp Viral RNA Kit (Qiagen, Germany) as directed by the manufacturer. Viral cDNA was synthesized using RevertAid First Strand cDNA Synthesis Kit (Thermo Fisher Scientific, USA) (Zahid et al., 2016). The selected DNA fragment was amplified through PCR in a 25µL reaction mixture, containing 12.5µL master mix (PrimeSTAR Max DNA polymerase, Catalog Number: R045A), 1µL DNA template, 1µL each forward and reverse primers and 9.5µL O n l i n e

Statistical and phylogenetic analysis
The obtained oligonucleotide sequences were properly trimmed with Chromas Software (version 2.6.6). Multiple sequence alignments were performed on BioEdit Software (version 7.2.5) using the ClustalW (version 1.8) (Thompson et al., 1997). The evolutionary analysis was conducted on MEGA Software (version 6) inferred through the Neighbor-Joining method with the maximum likelihood method. The reliability of the tree was tested by performing 1000 bootstrap replicates (Tamura et al., 2013). The pairwise nucleotide (nt) and amino acid (aa) distance matrix were calculated with UGENE software.

RESULTS AND DISCUSSION
In the present study, the results of RT-PCR revealed that the positive rate of IBDV samples was 32% (16/50). Afterwards, four positive samples were selected randomly to sequence the partial VP1 and VP2 gene fragments and the nucleotide sequences were determined and submitted to NCBI GenBank database ( Table I). The phylogenetic analysis revealed that serotype I was subgrouped into vvIBDV, classical, attenuated, and variant IBDV. All the vvIBDV isolates, including reference strain (NC_004178), were placed in a single clade with a bootstrap value of 100. In contrast, classical, attenuated, and variant IBDV were grouped into another monophyletic clade. In the evolutionary tree based on the 89 nucleotide sequences of VP1 retrieved from GenBank, serotype I strains revealed distinct branches and current study isolates with previously reported Pakistan strains were branched out of the vvIBDV branch and formed a cluster with HLJ0504-like strains, including HLJ0504 and Harbin-1 (Fig. 1). Based on the 72 nucleotide sequences of VP2, the phylogenetic tree was distinctly divided into two major branches, vvIBDV and non-vvIBDv. Present study isolates belonged to vvIBDV branch with already reported isolates (Pakistan reference strains) and formed one subgroup with a bootstrap value of 100 (Fig. 2). Segment A-coded VP2 gene phylogenetic analysis exhibited that our isolate had a close relationship with other Pakistan and Indian isolates. Currently identified isolates acquired both segments (A and B) from different genetic relatives and were further confirmed by sequence alignment and phylogeny. VP2 gene (segment A) of current isolates was clustered into vvIBDV branch while VP1 gene (segment B) branched out of vvIBDV and formed a unique branch with strain like HLJ0504 previously reported in China Qi et al., 2011), Algeria (Nwagbo et al., 2016) and Pakistan . Strains like HLJ0504 and Harbin-1 with unique segment (segment A originated from vvIBDV) and segment B acquired from a unique ancestor. So the current study reported unique segmentreassortant IBDV from backyard poultry. The results agree with (Hussain et al., , 2020, who reported Pakistan originated reassortant strains from commercial poultry. Concerning the identity matrix analysis of VP1 representative part, study isolates shared (87-88% nt; 95-96% aa) homology with the vvIBDV and (97-98% nt; 98-99% aa) with non-vvIBDV. VP1 of isolates have comparatively higher identity of (97-98% nt; 98-99% aa) with previously reported HLJ0504-like strains which are the unique reassortant retrieving segment A from vvIBDV and B from unique ancestor (He et al., 2014;Qi et al., 2011). Concerning the identity matrix analysis of VP2 representative part, present study isolates revealed 99-100% nt and 99-100% aa identity with previously reported Pakistan vvIBDV strains and 91-95% nt and 91-95% aa homology with non vvIBDV.  In VP1 representative part sequence alignment, Pakistan IBDVs has three characteristic aa residues of vvIBDV (287-A, 508-K, and 511-S) and four characteristic aa residues of non-vvIBDV (146-E, 147-G, 242-D and 390-L) (Fig. 3). There are three domains in VP1 of IBDV; N-terminal (1-167aa), the polymerase (168-658aa), and C terminal domain (659-878aa). Pakistan IBDV isolates revealed recombinant characters in the N-terminal domain and the central polymerase domain between very virulent and attenuated strains. Self guanylylation site of RNA dependent RNA polymerase responsible for Protein priming of stains found in the N terminal domain (1-67 aa). Additionally, sequence alignment of VP2 revealed eight characteristic aa residues of vvIBDV (222-A, 242-I, 253-Q, 256-I, 279-D, 284-A, 299-S, and 330-S and a distinct aa 384-I) (Fig. 4) also presented in Pakistan IBDVs previously reported by (Hussain et al., 2020). Likewise the reassortants are also prevailing in the neighboring countries including India and China, it is thought-provoking to discover the relation among them.

O n l i n e F i r s t A r t i c l e
As far as genotype is concerned current isolates belong to A3B3 genotype (Table II) according to the new classification scheme proposed by Wang et al. (2021) and  binomial scheme for genotype characterization suggested by Islam et al. (2021)

CONCLUSION
The naturally occurring reassortant strains of IBDV genotype A3B3, harboring segment A from a very virulent strain and segment B from one unique ancestor, were identified through molecular characterization. To our knowledge, this is the first report on characterization of segment B coding VP1 and segment A-coding VP2 genes in Pakistan backyard poultry. The current work gives insight into a better knowledge of the nature and evolution of IBDV, which will assist in disease prevention and control programs for emerging reassortant strains.

O n l i n e F i r s t A r t i c l e
Phylogenetic Analysis of Infectious Bursal Disease Virus Segment A and B