Biological Characterization and Evolutionary Dynamics of Pigeon Paramyxovirus Type 1 in China

Pigeon paramyxovirus type 1 (PPMV-1) is considered as an antigenic variant of Newcastle disease virus (NDV) which has an obvious host preference for pigeons and has caused significant economic losses to the global poultry industry. The evolutionary dynamics of PPMV-1 in China, however, are poorly understood. In this study, we characterized seven PPMV-1 isolates from diseased pigeons collected in Jiangsu, Anhui, and Henan provinces during 2020. Phylogenetic analysis revealed that seven isolates belonged to sub-genotype VI.2.1.1.2.2. Biological characterization indicated that seven isolates were mesogenic based on the mean death time (69.6–91.2 h) and intracerebral pathogenicity index (1.19–1.40) and had similar growth kinetics in chicken embryos and CEFs. Furthermore, the four representative viruses (AH/01/20/Pi, JS/06/20/Pi, HN/01/20/Pi, and HN/02/20/Pi) could result in marked cytopathic effects (CPE) in CEFs and induced syncytium formation in Vero cells. Our Bayesian phylogenetic analysis showed that PPMV-1 might first emerge in East China in 1974 and East China had the highest genotypic diversity of PPMV-1. Besides, phylogeographic analysis indicated that East China and South China were probably the major epicenters of dissemination of PPMV-1 in China. Selection pressure analysis and amino acid substitutions analysis revealed that the viral replication complex (NP, P, and L proteins) was likely related with the host preference of PPMV-1. Collectively, this study uncovered the epidemiology and evolutionary dynamics of PPMV-1 circulating in China, emphasizing the importance of strengthening the monitoring of PPMV-1 in East China and South China and providing significant clues for further studies on the molecular mechanism underlying host preference of PPMV-1.


INTRODUCTION
Newcastle disease (ND), caused by virulent Newcastle disease virus (NDV), is an acute and highly contagious disease that has caused significant economic losses to the global poultry industry (1,2). NDV is classified as a member of the genus Orthoavulavirus, family Paramyxoviridae (3). The genome of NDV is ∼15.2 kb in length (4) and encodes six structural proteins, including nucleocapsid (NP) protein, phosphoprotein (P), matrix (M) protein, fusion (F) protein, hemagglutinin-neuraminidase (HN) protein, and large (L) protein (5). Additionally, the P gene encodes two additional non-structural proteins (V and W) (6,7). NDV isolates can be grouped into three pathotypes according to their degree of virulence in chickens: velogenic (highly virulent), mesogenic (moderately virulent) and lentogenic (avirulent or low virulent) (8). NDV isolates are phylogenetically categorized into two distinct classes, class I and class II (9). Class I viruses comprise a single genotype and are frequently isolated from wild birds and domestic waterfowls (10)(11)(12). Class II viruses have emerged as a major cause of ND outbreaks worldwide and can be segregated into 20 genotypes (I to XXI, no XV) (13).
Pigeon paramyxovirus type 1 (PPMV-1) is an antigenic variant of NDV which has a high host preference for pigeons (14,15). PPMV-1 was first identified in the Middle East in the late 1970s and subsequently spread rapidly throughout Europe, resulting in the third ND pandemic during the 1980s (16)(17)(18). Currently, PPMV-1 is still prevalent in several countries worldwide (19)(20)(21)(22). Pigeons are mainly split into three types: meat-type pigeons, homing pigeons, and fancy pigeons, with meat-type pigeons being the most common breed in China. With the rapid development of China's social economy and the continuous improvement in people's quality of life in recent years, the pigeon industry has also developed rapidly, but the epidemic situation of diseases caused by PPMV-1 is not optimistic. Since its first isolation in China in the 1980s, PPMV-1 has long been prevalent in numerous provinces and regions of China, seriously impeding the development of the pigeon industry (23)(24)(25)(26)(27). As a result, it is critical to strengthen the epidemiological surveillance of PPMV-1 in China. Bayesian phylodynamic inference is widely regarded as one of the most employed methods for speculating on the origin of epidemics, tracking the geographic spread and transmission dynamics of infectious pathogens (28). Several previous studies have suggested that PPMV-1 of China might have originated in Europe (14,27,29,30). However, the evolutionary and transmission dynamics of PPMV-1 in China remain not well-understood.
In this study, we identified seven PPMV-1 isolates from diseased pigeons collected in Jiangsu, Anhui and Henan Provinces from January 2020 to December 2020 and determined the biological characteristics of these isolates. We further deeply investigated the prevalence, evolution and transmission dynamics of PPMV-1 prevailing in China. Moreover, we also performed the selection pressure analysis of different genes of PPMV-1 in China. Additionally, we identified specific amino acid substitutions in different proteins of PPMV-1 isolates by comparing the complete genome sequences of these viruses with those of chicken-origin NDVs. This study contributes to a more systematic and comprehensive understanding of PPMV-1 prevailing in China.

Cells, Eggs, and Animals
Primary chicken embryo fibroblasts (CEFs) were prepared and cultured as described previously (31). All embryonated SPF

Virus Isolation and Biological Characterization
Seven PPMV-1 strains were isolated from pigeons in Jiangsu, Anhui, and Henan Provinces in China during January-December 2020. Table 1 presents the detailed information about the virus names, hosts, isolation dates, locations, the size of the flock, the clinical outcome for the flock and the age of the pigeons. Hemagglutination (HA), hemagglutination inhibition (HI), and reverse transcription-polymerase chain reaction (RT-PCR) assays were used for detection and identification of PPMV-1 isolates. All PPMV-1 isolates were then purified three times by plaque assay on CEF cells (32) and subsequently propagated in 10-dayold SPF chicken embryos. The virulence of all PPMV-1 isolates was evaluated by mean death time (MDT) in 10-day-old SPF embryonated chicken eggs and intracerebral pathogenicity index (ICPI) in 1-day-old SPF chicks (33): the MDT test in 10-dayold SPF chicken embryos (hours) (>90, lentogenic strains; 60-90, mesogenic strains; <60, velogenic strains), and the ICPI test in 1day-old SPF chicks (<0.7, lentogenic strains; 0.7-1.5, mesogenic strains; >1.5, velogenic strains) (34).
To access the growth kinetics of the PPMV-1 isolates in chicken embryos, the allantoic cavities of 10-day-old SPF embryonated chicken eggs were inoculated with 100 TCID 50 of each virus per embryo. Three chicken embryos were chilled and allantoic fluid containing virus was harvested at 24 and 48 h post-infection (hpi). Additionally, to determine the growth kinetics of the PPMV-1 isolates in chicken embryo fibroblasts (CEFs), primary CEFs grown in six-well plates were infected with the PPMV-1 isolates at a multiplicity of infection (MOI) of 0.01. At 24 and 48 hpi, supernatant samples were taken and replaced with an equal amount of fresh medium. Virus titers were evaluated by performing 50% tissue culture infective dose (TCID 50 ) assays in CEF cells. The TCID 50 was calculated by the Reed-Muench method.
For the cytopathic effect assay, primary CEFs were grown in six-well plates and infected with the four PPMV-1 isolates

Whole-Genome Sequencing
Viral RNA was extracted directly from the collected allantoic fluid with EasyPure R Viral DNA/RNA Kit (TransGen Biotech, Beijing, China). PCR amplification was conducted using 11 pairs of universal primers specific to the complete genomes of PPMV-1. The universal primers used in the present study are shown in Table 2. To obtain the 3 ′ and 5 ′ ends of the viral genomes, the PCR-generated fragments containing the leader and trailer sequences of all 7 PPMV-1 isolates were cloned into the pEASY R -Blunt Zero vector (TransGen, China) for sequencing. The sequencing was carried out by Tsingke Biological Technology (Nanjing, China). All assembled genomes were submitted to GenBank (see Table 1 for accession numbers).

Collection of Sequences
A total of 427 complete F gene sequences of NDV (PPMV-1 and chicken-origin NDV belonging to genotype XX) were collected from GenBank (as of February 2021) and were aligned by using Multiple Alignment with Fast Fourier Transformation (MAFFT v.7.4.50). Additionally, PPMV-1 isolates from China were selected and a dataset containing 186 complete F gene sequences was constructed. Furthermore, 109 complete genomes of NDV (seven strains belonged to genotype XX and the other 102 strains belonged to genotype VI) were compared using AliView software, and the datasets containing the coding  region sequences (CDS) of six genes of NDVs were constructed separately for selection pressure analysis.

Phylogenetic Analysis and Demographic History
Based on 427 complete F gene sequences of NDV, a maximum likelihood (ML) tree was constructed using the GTR + F + R4 model as suggested by ModelFinder (37) and was implemented in IQ-TREE with 1,000 bootstrap replicates (38). To verify the temporal structure of all PPMV-1 isolates in China, we regressed root-to-tip genetic distances against sampling date using TreeTime v.0.8.1 (39

Discrete Phylogeographic Analysis
To gain insight into the spatial and temporal dynamics of PPMV-1 in China, Bayesian phylogeographic analysis was performed to infer the migration pathways by using BEAST v. SpreaD3 v.0.9.7 was used to visualize the spreading process over time and calculate the Bayes factor (BF) support for each route (45). Significant migration routes were determined based on the established criteria (BF ≥ 3 and the posterior probability >0.5).

Selective Pressure Analysis
To analyze the positive selection sites in different genes of PPMV-1 isolates prevailing in China, the Fixed Effects Likelihood (FEL) method, the Mixed Effects Model of Evolution (MEME) method, the Fast Unconstrained Bayesian AppRoximation (FUBAR) method and Single Likelihood Ancestor Counting (SLAC) method provided on the Datamonkey website (http://www. datamonkey.org/) were used in this study. Significant levels of p < 0.1 in SLAC, p < 0.05 in FEL, p < 0.05 in MEME and posterior probabilities >0.9 in FUBAR were considered as an indication of positive selection.

Amino Acid Substitutions Analysis
Previous studies have shown that PPMV-1 is an antigenic variant of chicken-derived NDV (14,46,47). To clearly elucidate the alteration of specific amino acid residues in all proteins of PPMV-1 (102 strains belonging to genotype VI) compared with chickenorigin NDV (7 strains belonging to genotype XX), this study used AliView to compare the complete genomes of 109 NDV strains.

Statistical Analysis
Statistical analysis was carried out using GraphPad Prism version 8.0 (GraphPad Software Inc., San Diego, CA). Error bars in the data from virus growth kinetics indicate standard errors of the means. P-values were evaluated by one-way ANOVA to determine the significance level of the difference between distinct experimental groups. Statistical significance was taken as follows: * P ≤ 0.05, significant; * * P ≤ 0.01, very significant; * * * P ≤ 0.001, extremely significant.

Isolation and Biological Characterization of PPMV-1 Strains
In the present study, 7 PPMV-1 strains were isolated from diseased pigeons from three different provinces (Anhui, Jiangsu, and Henan) during January 2020-December 2020 ( Figure 1A and Table 1). The clinical signs of the diseased pigeons contained marked depression, green diarrhea, paralysis, neurological symptom, and the gross lesion included severe congestion in the brain (Figure 1A). The HA titers of the PPMV-1 isolates ranged from 5 log2 to 6 log2 ( Table 1). The complete genomes of seven  PPMV-1 strains were determined. The amino acid sequences of the F protein cleavage site ( 112 RRQKRF 117 ) of these isolates were indicative of virulent viruses (Table 1). Nevertheless, the MDT values ( Figure 1B) for all isolates were between 69.6 and 91.2 h and the ICPI values ( Figure 1C) ranged from 1.19 to 1.40, implying that all isolates were moderately virulent (mesogenic). Additionally, the results showed that all isolates exhibited similar growth kinetics in chicken embryos and CEFs at 24 and 48 hpi (Figures 1D-G). HN/01/20/Pi, and HN/02/20/Pi) could lead to syncytium formation at 36 hpi ( Figure 1I). However, these results suggested that there was no significant difference in the fusion activity of the four isolates (AH/01/20/Pi, JS/06/20/Pi, HN/01/20/Pi, and HN/02/20/Pi) ( Figure 1J).

Phylogenetic Analysis
A ML phylogenetic tree was constructed based on the complete F gene sequences. As shown in Figure 2 (Figure 3). The Bayesian analysis also indicated that all PPMV-1 isolates were rooted with the virus strains from East China, and PPMV-   that PPMV-1 has a high host preference for pigeons ( Figure 4A). In terms of the spatial distribution of the 186 PPMV-1 strains (Figure 4B) (Figure 4C).

Demographic History
The historical population dynamics of PPMV-1 circulating in China indicated that the change in the population size might go through six stages (I-VI) (Figure 5). In the first stage (before 1995), the population size of PPMV-1 isolates in China remained relatively constant. Then, the PPMV-1 population

Migration Pattern of PPMV-1 in China
We conducted Bayesian phylogeographic analysis to elucidate the spatial transmission patterns of PPMV-1 in seven geographic regions of China (East China, South China, North China, Northeast China, Southwest China, Northwest China, and Central China). The results showed that the presence of PPMV-1 might limit to East China before 1995 (Figures 6A-C), and then spread from East China to North China, South China, and Southwest China from 2000 to 2005 (Figures 6D,E). PPMV-1 started to spread to other regions from South China by the year 2010 ( Figure 6F). Since then, South China has become another new epicenter besides East China ( Figure 6G). Notably, PPMV-1 has been found in all seven geographic regions of China by 2020 ( Figure 6H).
As exhibited in Figure 7 and Table 3, these outcomes demonstrated that there were 10 migration links in the diffusion of PPMV-1 in China. The highest average transition rate was observed for migration from East China to South China, while the lowest mean transition rate was found for migration from East China to Northwest China (Table 3). Moreover, East and South China might be the two main epicenters responsible for eight migration routes of PPMV-1 (Figure 7). Among the eight migration routes, three transmission routes are decisively

Positive Selection Sites in Different Genes of PPMV-1 in China
To gain a better understanding of the selection pressure acting on different genes of PPMV-1 in China, we conducted a positive selection analysis using four disparate detection methods (FEL, MEME, FUBAR, and SLAC methods). As shown in Table 4, three codons (sites 8, 432 and 486) were found to be under positive selection in the NP gene of PPMV-1 according to one detection method. Furthermore, three positive selection sites (47, 89, and 192) were simultaneously predicted in the P gene by the four methods. One positive selection site (73) was detected in the M gene by the FUBAR method. Within F gene, two positive selection sites (4 and 28) were detected based on all four methods, whereas one positive selection site (13) was predicted using two methods. In addition, four positive selection sites (92, 115, 266, and 407) were identified in the HN gene by using one method.  One positive selection site (1,672) was identified by all four methods, and the other positive selection site (1,365) in the L gene was detected by both MEME and FUBAR methods, Moreover, six positive selection sites (147, 1,759, 1,804, 2,013, 2,029, and 2,064) were found in the L gene of PPMV-1 by one method. Combining the above results, this study demonstrated that the majority of the positive selection sites were located at the NP, P and L genes.

Amino Acid Substitutions Analysis
It has been shown that PPMV-1 might be an antigenic variant of chicken-origin NDV (14,46,47). To investigate the amino acid differences between PPMV-1 and chicken-origin NDV, we compared 109 full-length sequences in this study. The results presented in Table 5 showed the presence of specific amino acid substitutions in all six proteins of PPMV-1 compared to chicken-origin NDV: six different substitutions (G106S, A154V, S405N/K, N463D, T469I, and P477H) in the NP protein; seven disparate substitutions (D10E, S92P, N/I93K/R, E103D/N/G, R163G, A343V/M, and R380K) in the P protein; two distinct substitutions (Q36R and M77T) in the M protein; one different substitution (V50I) in the F protein; six distinct substitutions (V9A, Y112H, G304S, E/K347G, D349E, and A540V/I) in the HN protein; ten disparate substitutions (T129A, T340S, M1079T/I, V1082I, S1099P, N1110A/T, F/S1727P/L, A1757V, K2120R, and D/N2173E/G) in L protein ( Table 5). In summary, PPMV-1 had 23 specific amino acid substitutions in the NP, P, and L proteins and nine specific amino acid substitutions in the M, F and HN proteins compared with chicken-origin NDV.

DISCUSSION
Since the introduction of PPMV-1 from Europe into China, it has long been prevalent among pigeon flocks for a long time, leading to significant economic losses in the pigeon industry nationwide (23,25,26,48). In this study, we identified seven PPMV-1 strains from diseased pigeons in Jiangsu, Anhui, and Henan provinces during 2020. Although the F protein cleavage site sequences of seven PPMV-1 strains were typical of velogenic viruses, these PPMV-1 isolates belonged to mesogenic strains based on MDT and ICPI values which indicated that other structural domains or genes might influence viral virulence in addition to the F protein cleavage site (31,49). Phylogenetic analysis showed that these seven viruses were classified as sub-genotype VI.2.1.1.2.2. Additionally, all PPMV-1 isolates in China were clustered into four sub-genotypes  (26,27). From the host distribution, most of PPMV-1 strains (178/186) were isolated from pigeons, further confirming that PPMV-1 has a significant host preference for pigeons (15). The PPMV-1 strains isolated from East China were located at the root of the MCC tree with the highest posterior probability of 0.69, suggesting that East China was probably the origin of PPMV-1 in China. Moreover, TMRCA of PPMV-1 isolates was dated to the year 1974 (95% HPD: 1956-1986), well before the first genome sequence of the PPMV-1 strain from China was submitted to GenBank. Notably, sample bias, scant sampling and sequencing might potentially influence this conclusion, especially given the scarcity of data in the epidemiological surveillance in the last century (27,50). Over the last 25 years, these viruses transmitted from East China to the rest of the nation and East China and South China were likely the major epicenters responsible for the dissemination of PPMV-1. The suitable climatic conditions, convenient transportation and high population density could account for that East China and South China might play a critical role in the spread of PPMV-1. Hence, strengthening the monitoring of PPMV-1 in East and South China can contribute to understanding the evolutionary dynamics of PPMV-1 in these two regions and predicting the epidemic in other regions in a timely manner, as well as providing a theoretical basis for the development of reasonable PPMV-1 prevention and control strategies. The demographic history indicated that the population size change of PPMV-1 in China might have gone through six different phases. The population size of PPMV-1 in China remained relatively constant during phase I (prior to 1995), which might be attributed to the low demand for the young pigeons and the slow development of the pigeon industry during this period. The population of PPMV-1 increased rapidly at stage II (1995)(1996)(1997)(1998)(1999)(2000)(2001)(2002), probably due to the rapid growth of the pigeon industry in China along with poor feeding management and the lack of a sanitary and epidemic prevention system (50). Then, the PPMV-1 population remained basically constant in the third phase (2003)(2004)(2005)(2006)(2007)(2008), which was likely caused by the fact that the SARS outbreak dealt a huge blow to the pigeon industry in China during 2002-2003, as well as the increasing use of conventional ND vaccines such as La Sota. The rapid increase in PPMV-1 population size at phase IV (2009-2013) could be because of the following two reasons: on the one hand, the pigeon industry entered a new period of rapid expansion with the increasing demand of the young pigeons in China; on the other hand, in the scenario where the circulating PPMV-1 isolates (genotype VI) were sufficiently divergent from the La Sota vaccine strain (genotype II), the vaccine-induced immune response did not confer better protection and prevent virus shedding (51). Interestingly, the population size of PPMV-1 in China decreased slightly at stage V (2013-2014), likely because of the significant impact of H7N9 avian influenza epidemic on the pigeon industry within this period. Subsequently, the PPMV-1 population remained relatively constant at stage VI (after 2014) which might be explained by the fact that although the pigeon breeding industry was able to gradually recover and entered a new stage of development as the demand for suckling pigeons increased, intensive farms with better biosafety developed rapidly with the decrease of individual backyard farms.
Adaptive evolution of genes and genomes determines evolutionary morphology, innovation and species divergence (52,53). Previous studies have shown that PPMV-1 might be a variant of chicken-derived NDV (14,46,47). Besides, our previous study indicated that PPMV-1 changed in virulence and obtained a survival advantage over chicken-origin NDV in pigeons after undergoing adaptive evolution (31). In this study, selection pressure analysis showed that positive selection sites existed in all six genes of PPMV-1, with most of them located at the NP, P, and L genes. In addition, amino acid substitutions analysis indicated that the majority of specific amino acid substitutions (23/32, 71.9%) were situated in the NP, P, and L proteins of PPMV-1 compared to chicken-origin NDV. The aforementioned findings suggested that the viral replication complex (NP, P, and L proteins) might be associated with the host preference of PPMV-1. However, further validation studies utilizing reverse genetics techniques are needed in the future.
All in all, we isolated seven PPMV-1 strains from Jiangsu, Anhui, and Henan provinces in 2020. This study systematically revealed the diffusion pattern of PPMV-1 isolates among different regions of China. Our findings demonstrated that East China and South China might play a critical role in the dissemination of PPMV-1. This work will be of benefit to us for a better understanding of the evolutionary dynamics of PPMV-1 in China and will have major significance for PPMV-1 prevention and control in China, as well as providing important clues for further in-depth studies on the molecular mechanism of host preference of PPMV-1.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

ETHICS STATEMENT
The animal study was reviewed and approved by the Jiangsu Administrative Committee for Laboratory Animals.