Diversity of avipoxviruses in captive-bred Houbara bustard

Implementation of conservation breeding programs is a key step to ensuring the sustainability of many endangered species. Infectious diseases can be serious threats for the success of such initiatives especially since knowledge on pathogens affecting those species is usually scarce. Houbara bustard species (Chlamydotis undulata and Chlamydotis macqueenii), whose populations have declined over the last decades, have been captive-bred for conservation purposes for more than 15 years. Avipoxviruses are of the highest concern for these species in captivity. Pox lesions were collected from breeding projects in North Africa, the Middle East and Central Asia for 6 years in order to study the diversity of avipoxviruses responsible for clinical infections in Houbara bustard. Molecular and phylogenetic analyses of 113 and 75 DNA sequences for P4b and fpv140 loci respectively, revealed an unexpected wide diversity of viruses affecting Houbara bustard even at a project scale: 17 genotypes equally distributed between fowlpox virus-like and canarypox virus-like have been identified in the present study. This suggests multiple and repeated introductions of virus and questions host specificity and control strategy of avipoxviruses. We also show that the observed high virus burden and co-evolution of diverse avipoxvirus strains at endemic levels may be responsible for the emergence of novel recombinant strains.


Introduction
Avipoxviruses (APV) are large enveloped, double stranded DNA viruses able to naturally infect more than 232 species from 23 orders [1]. They belong to the genus Avipoxvirus, subfamily Chordopoxvirinae and family Poxviridae. To date, 10 APV species have been defined by the International Committee on Taxonomy of Viruses [2]. Phylogenetic analyses of APV have previously been based on a couple of genes. An amplification product of 578 bp from the fpv167 gene, orthologous of vaccinia virus A3L gene encoding a virion core protein 4b, has been widely used since this gene is highly conserved among poxviruses [3]. Other genes have been used in order to validate or improve findings based on fpv167. The fpv140 locus (fpv139, fpv140 and fpv141 genes), has also been used [4] and the phylogenetic analyses based on fpv140, orthologous of vaccinia virus H3L gene encoding a virion envelope protein p35, have provided improved strain discrimination within some subclades [5][6][7]. Phylogenetic analyses have all shown the same distribution of APV into 3 clades: fowlpox virus-like (FWPV-like, clade A), canarypox virus-like (CNPV-like, clade B), and psittacinepox virus-like (clade C).
APV infections can cause significant economic losses in domestic poultry, due to decreased egg production, reduced growth, blindness and increased mortality [8]. In wild birds, the infection can compromise survival and breeding success by decreasing the ability to escape predators [9], to fledge and rear chicks [10], by impairing the pairing success [11] or by increasing mortality [12]. Recently, emergence of APV infections in great tits in Europe and especially in the UK has been described and could represent a threat for the species [13,14]. The impact of the disease can be dramatic for endangered species, especially for island species when APV are introduced in naïve populations as described in Hawaii [15], the Galapagos Islands [16] or Canary islands [17]. The disease is also a concern for the success of conservation programs of endangered species, such as peregrine falcons in Germany [18] or Houbara bustard species (hereafter "Houbara") in the United Arab Emirates (UAE) and Morocco [19,20].
The African Houbara bustard (Chlamydotis undulata) and the Asian Houbara bustard (Chlamydotis macqueenii) inhabit semi-arid areas in North Africa, the Middle East and Central Asia. In the last decades, their populations drastically declined due to over-hunting, habitat degradation and poaching [21] leading the species to be listed as "Vulnerable" by the International Union for Conservation of Nature (IUCN 2012). For more than 15 years, captive-breeding programs of the Houbara have been implemented in North Africa, the Middle East and Central Asia. Their goal is to increase the size of wild populations through release of captive-born individuals, which should substantially decrease extinction risk. APV infections are known to compromise the success of Houbara captive-breeding programs by increasing mortality [19,20]. Persistent presence of infections is observed in captive flocks with morbidity rates of 2-3% in juvenile birds (≤1 year) despite prophylactic programs (G. Le Loc'h, unpublished data) but very little information is available about the APV strains circulating in Houbara.
The present study was aimed at assessing APV diversity within a species as well as virus circulation within captive-breeding programs. Molecular and phylogenetic analyses of APV in Houbara highlight that a large pathogen diversity and burden is possible in an apparent unfavorable environment.

Animal samples
The study was conducted from 2008 through 2013 in 3 Houbara captive-breeding projects: the Emirates Center for Wildlife Propagation in Morocco, the National Avian Research Center and the Sheikh Khalifa Houbara Breeding Center in the UAE and the Emirates Centre for the Conservation of Houbara in Uzbekistan. Birds housed in these projects were either C. macqueenii alone (UAE and Uzbekistan) or C. macqueenii and C. undulata (Morocco). Pox lesions were detected during routine check-up as part of veterinary management of birds. When typical lesions of APV infections, either cutaneous (nodular lesions on nonfeathered areas) or diphtheritic (yellowish lesions on the mucous membranes of the mouth, esophagus, or trachea), were observed [8], they were aseptically removed, transferred into a sterile pot and stored at −80°C until analysis. In total, 169 lesions were sampled (Table 1): 52 in C. undulata and 117 in C. macqueenii, mainly on juvenile birds (≤ 1 year; n = 143, 85%). DNA was extracted from 25 mg of frozen tissue samples with the NucleoSpin Tissue kit (Macherey-Nagel, Düren, Germany) according to the manufacturer's instructions.

Phylogenetic and recombination analyses
After manual editing and excluding primers with BioEdit 7.2.5, P4b and fpv140, DNA sequences were aligned with sequences available on GenBank using ClustalX 2.1. All available sequences on GenBank were initially used, then, one unique sequence of each group was kept for phylogenetic analyses. APV strains collected during the present study were named with the following nomenclature: APV clade (as determined by phylogenetic analysis of P4b: FWPV, CNPV and PSPV)/host species/origin/strain number/collection year. For example, CNPV/Chlamydotis macqueenii/AE/072/2011 is a Canarypox virus-like strain number 072 collected from a Chlamydotis macqueenii in the UAE in 2011. GenBank accession numbers are available in Table 2. Phylogenetic trees were generated using MEGA 5.2.1 by neighbor-joining (NJ) with the Kimura 2parameters model and reliability of trees was tested through 1000 bootstrap replicates. After identification of clades and subclades, within and between subclades mean genetic distances were calculated with the same algorithm.
A recombination analysis was performed with RDP v3.44 in order to detect potential recombination events for CNPV/Chlamydotis undulata/MA/001/2009, which grouped differently when analyzing P4b and fpv140. Concatenated sequences of fpv140 (2750 to 2878 bp) and P4b (487 to 490 bp) loci were treated as linear by RDP, GeneConv, BootScan, MaxChi, Chimaeara, SiScan and 3Seq methods. Settings were adjusted to detect events with a highest acceptable p value equal to 0.0001.
In order to exclude a mix of different sequences in CNPV/Chlamydotis undulata/MA/001/2009, which could explain differences of positioning in phylogenetic trees, primers specific to either B1 (FPV140_03F and FPV1 40_03R) or B2 subclades (FPV140_02F, and FPV140_02R), and targeting the fpv140 locus were used. PCR were performed on CNPV/Chlamydotis undulata/MA/001/ 2009 and one additional specimen of subclades B1 and B2 were used as positive controls.

Phylogeny
Out of the 169 samples, 139 DNA sequences were successfully amplified with at least one set of primers. Discrimination between FWPV-like viruses and CNPV-like viruses based on fpv140 loci fragment size shows a slightly higher proportion of CNPV-like viruses (n = 63, 55%) than FWPVlike viruses (n = 51, 45%).
Thereafter, 113 and 75 DNA sequences were obtained for P4b and fpv140, respectively. P4b sequences were trimmed to a 426 bp length fragment, which was used to carry out phylogenetic analysis. From fpv140 sequences, the 3 open reading frames (ORF) corresponding to the 3 orthologous genes shared by FWPV-like viruses and CNPV-like viruses (fpv139/CNPV184, fpv140/CNPV186 and fpv141/CNPV187) were used to build independent phylogenetic trees.
The NJ tree based on P4b DNA sequences ( Figure 1) provided clear distinction between known APV clades and subclades as previously described [4,23]. The strains identified in the present study clustered in 11 genotypes belonging to clades A and B. Most of the strains of clade A grouped either in subclade A2 (n = 30, genotype represented by FWPV/Chlamydotis macqueenii/AE/066/2012) with 100% identity with sequences obtained from a wide diversity of wild and domestic species, or in subclade A1 (n = 24, genotype represented by FWPV/Chlamydotis macqueenii/AE/070/2011) with 100% identity with sequences from domestic Galliformes (e.g. Gallus gallus). FWPV/ Chlamydotis macqueenii/AE/150/2013 grouped in subclade A3 and was identical to sequences from Columbiformes (e.g. Columba palumbus) and a great bustard (Otis tarda) in Spain. The 8 genotypes of clade B grouped in 3 subclades: the previously described subclades B1 (n = 38) and B2 (n = 19), and a putative new subclade B4 (n = 1).  Figure 1). The closest sequence (90% identity) was isolated from an American robin (Turdus migratorius) in the United States and clustered in subclade B3. The minimal genetic distance between subclade B4 and other B subclades was 10.8% (Table 3). The complete fpv140 gene sequence (933-1020 bp) was used to build an NJ tree ( Figure 2). Its analysis provided the clades/subclades classification as previously described for clade A [4,5]    To exclude a mix of sequences in CNPV/Chlamydotis undulata/MA/001/2009, a PCR was carried out on its DNA showing an amplification product with subclade B2 specific primers but no amplification product with subclade B1 specific primers (Figure 4).

Geographic distribution
Out of the 113 P4b sequences obtained in the present study, 57 were collected in Morocco, 46 in the UAE and 10 in Uzbekistan ( Figure 5). Sequences of several subclades were obtained on each site with differences among projects. CNPV-like viruses were predominant in Morocco with 74% of the strains grouping in clade B (47% and 25% in subclades B1 and B2, respectively). FWPV-like viruses were mainly represented by strains of subclade A2 (25%). In the UAE, 87% of the strains were FWPV-like viruses, grouping in subclades A1 and A2 for 50% and 35% of the viruses, respectively. Only CNPV-like viruses were identified in Uzbekistan with 70% of the strains grouping in subclade B1. Strains of subclades A3 and B4 were identified in only one site, the UAE and Morocco, respectively.

Discussion
APV infections in Houbara have been poorly studied and only 2 DNA sequences from these species are published in GenBank. This study contributes to a better comprehension of the disease in Houbara by providing 113 sequences for P4b and 75 sequences for fpv140, collected over 6 years in 3 different countries, which cover the natural distribution of Houbara. Virus isolates have been obtained for a few representative strains (FWPV and CNVP-like viruses: at least 1 isolate per subclade, data not shown).
The NJ phylogenetic tree built from P4b DNA sequences was very similar to those already published [4,5,23] for clade A and distinguished 7 subclades, while CNVP-like sequences from Houbara provided more details in clade B. The use of a second locus gave a better resolution within and between some subclades. The fpv140 locus encompasses an entire ORF (fpv140 gene) ranging from 933 to 1020 bp and despite reports of difficulties to amplify the gene [24], its amplification was successful for 84% of Houbara strains for which a P4b gene sequence had been obtained. The number of APV sequences available has more than doubled thanks to the 75 sequences generated here. A putative new subclade (B4) was identified from an Asian Houbara in Morocco in 2011 (CNPV/Chlamydotis macqueenii/MA/083/2011). Its minimal genetic distance with other B subclades (10.8%) was higher than some between A subclade genetic distances (3.2% to 16.3%, Table 3), which justified considering it as a new subclade. CNPV/Chlamydotis macqueenii/MA/083/2011 seemed to share the same ancestor as subclade B3 isolates, all originating from the USA.
Subclade B1, whose diversity is the highest as compared with other subclades, is usually subdivided into 2 clusters [25]. More recently, a third cluster was described for APV isolates from Magellanic penguins [26]. With Houbara sequences, especially CNPV/Chlamydotis undulata/MA/ 001/2009 and CNPV/Chlamydotis macqueenii/MA/097/ 2012, diversity within subclade B1 increased and distinction of clusters becomes questionable. Also, it questions the definition of this group of virus as a unique subclade since the intra-B1 genetic distance (4.8%) is higher than the genetic distance between subclades A2 and A3 (3.2%). This is the first time that so many sequences are isolated from a single species (n = 84 for C. macqueenii; n = 29 for C. undulata), thus providing a better understanding of the epidemiology of APV for each species. The diversity of APV isolated in Houbara is surprising. Published sequences from Houbara suggested that they are usually infected by CNPV-like viruses [27]. In this study, we show that Houbara viruses grouped in 11 different genotypes and 6 different subclades when looking at P4b (17 genotypes, 7 subclades with fpv140) and almost half of them grouped with FWPV-like viruses. Classically, APV are considered to be host species or order specific and taxonomy was based on this concept. This has been recently questioned by many authors [4][5][6]23] since some taxa like Columbidae and Accipitridae can be infected by a wide diversity of strains, however with different sensitivity depending on the virus involved. This apparent diversity could be explained by some infections occurring as accidental events, especially in zoological collections where many species are housed closely, and it is suggested that such an infection could not lead to sustainable epornitics [23]. However, in the present study, all genotypes identified were responsible for clinical infections in Houbara and lead to nodular lesions on legs, beaks or eyelids without any evidence of phenotype-genotype relationship. CNPV-like viruses as well as FWPV-like viruses were also identified from diphtheritic lesions and could be responsible for outbreaks in Houbara captivebreeding projects as observed in Morocco in 2009 or in the UAE in 2013 with outbreaks of CNPV-like viruses and FWPV-like viruses respectively (G. Le Loc'h, unpublished data). This suggests that the taxonomy of APV may be reviewed and questions host specificity of APV.
The wide diversity of strains affecting Houbara has also been observed within a single breeding project and for a single year since strains of 5 clades were collected in Morocco in 2011 and strains of 4 clades were collected in the UAE in 2013. Previous studies have shown little diversity of strains in geographic areas of the same scale, within one species [26] as well as among different species [28]. This has been explained either by the quick spread of a competent virus in immunologically naïve birds [28] or by a limited number of introductions of strains in insular populations [25,29,30]. The high diversity observed in this study within a single breeding station suggests that multiple and repeated introductions of virus occurred. The origin of introductions can be hypothesized through observation of phylogenetic trees. The gene fpv140 provides a better resolution than the P4b locus for FWP-like viruses since each subclade used to be linked to a host species or order. Moreover the advantage of using fpv140 for virus typing could be attributable to the surface nature of the encoded H3L orthologue ORF. Subclade A1 viruses were exclusively collected from domestic chickens and Houbara strains clustering in this subclade all originated from the UAE where poultry farms  The multiple origin of viruses affecting Houbara questions the role of the species in the APV epidemiological cycle. It could be assumed that this species is especially sensitive to APV viruses regardless of their origin. Management of captive-breeding projects can also provide some clues to understanding our observation. Most of the rearing process takes place in outdoor aviaries [31], where the birds are housed in small groups, allowing potential contact with wild birds and insect vectors. Moreover, although projects are located in quite hostile and isolated landscapes (semi-desert or desert), they are usually built with all accommodations needed to ensure a certain autonomy (alfalfa fields to feed the Houbara, gardens for staff, etc.). Those artificial oases are especially attractive for wild fauna and can concentrate relatively large populations of wild birds and insects. In addition, some projects such as the ones in Morocco are located along migratory pathways. In contrast, due to a high conservation value of each Houbara and despite flock size of thousands of birds in each project, every bird is managed individually and carefully controlled with a high biosecurity level. While this management prevents quick and massive spread of pathogens, it may allow a higher diversity of virus strains entering and co-evolving at an endemic level in the breeding stations by breaking the epidemiological cycle of APV.
This observation questions the prophylactic strategy to control APV infection. Sensitivity of Houbara, as other birds' species, to different APV strains appears higher than expected. Frequent viral challenge even in an apparent unfavorable environment, as shown by the multiple introductions of APV in Houbara breeding stations, is also surprising. This makes the association of different prophylactic measures such as classical biosecurity, prevention of wild bird concentration, control of insect vector populations and vaccination, necessary. Vaccination alone is questionable regarding the wide diversity of strains affecting Houbara and the poorly known cross-protection among APV. Previous studies give contradictory conclusions about APV cross-immunity. While most of them show no cross-protection between clades [32,33] and even inside the same clade [34], a more recent study provides evidence of cross-protection between FWPV-like viruses and CNPV-like viruses in experimentally infected zebra finches [35]. In FWPV-like viruses, the absence of crossprotection has been shown to be linked to the presence of the integrated reticuloendotheliosis virus (REV) sequence in the FWPV genome [34]. In Houbara strains, REV was only detected in the FWPV subclade A1 viruses (data not shown). The latter studies suggest that cross-protection could be related to APV strains but also to host species. Further investigations are warranted to improve our understanding of the determinism of host specificity and pathogenicity among different APV strains.
The presence of several different strains within the same population at the same moment is favorable to the emergence of recombinant strains as already shown in orthopoxviruses [36] and is suspected for APV [23]. In Morocco, strains belonging to subclades B1 and B2 as well as a potential recombinant strain were isolated within the same breeding station in November 2009 and again in July 2011 (G. Le Loc'h, unpublished data), highlighting the high viral load in the field. The analyses confirmed the hypothesis of recombination observed on phylogenetic trees by showing intra and inter-loci recombination events. Interestingly, the parent strains from both subclades B1 and B2 (group of isolates CNPV/Chlamydotis macqueenii/MA/012/2010 and CNPV/Chlamydotis undulata/MA/025/2011, respectively) are the only ones isolated in the 3 Houbara breeding projects and their P4b sequences are identical to those of strains isolated from house sparrows in Morocco [37]. Moreover, the P4b sequence of the subclade B1 parent is identical to a strain distributed worldwide and that has been more especially isolated from great tits in Europe and associated with the emergence of APV infections in this species [13,14]. It can be hypothesized that this is a unique B1 strain that may have the ability to effectively and sustainably infect a wide diversity of avian hosts. However, whole genome sequences would be necessary to conclude on the unicity of the strain. Understanding APV evolution requires better knowledge of viral genomes. Up to now, the genomes of only 4 APV species have been fully sequenced: those of Fowlpox virus [38], Pigeonpox virus, Penguinpox virus [39], and Canarypox virus [40], 288 kbp and 365 kbp long, respectively. Their analyses have shown large genomic rearrangements and suggest significant genomic diversity among APV. Effort should be made to obtain complete genome sequences of several APV in order to better understand their host adaptability and their pathogenicity.