Natural intertypic and intratypic recombinants of enterovirus 71 from mainland China during 2009-2018: a complete genome analysis.

Surveillance of recombinant enterovirus 71 (EV71) and subgenotype replacement is vital for preventing and controlling hand, foot, and mouth disease (HFMD) outbreaks. Despite this, data on recombinant variants and phylogeny of circulating EV71 strains in mainland China are limited. In this study, recombinant variants of EV71 were identified in mainland China from 2009 to 2018. Phylogenetic analysis indicated that except for individual strains (CQ2014-86/CQ/CHN/2014 and EV71/Xiamen/2009 (B5)), almost all of the EV71 strains in mainland China belonged to the subgenotype C4a. Analysing complete genome sequences of 196 EV71 isolates, 3 intertypic recombination strains (VR1432, 30-2/2015/BJ, and Guangdong-2009) and 5 intratypic recombination strains (EV71/P1034/2013, VR1432, Henan-ZMD/CHN/2012, Hubei-WH/CHN/2012, and EV71/P868/2013/China) were identified among naturally circulating EV71. The breakpoints of these recombinant strains were located within the P1, P2, and P3 encoding regions. Notably, a double recombinant (VR1432) resulting from recombination between EV71 subgenotype C4a and C4b strain SHZH98 and a CA8 strain Donovan was identified. This study reports these specific intertypic and intratypic recombination events for the first time highlighting the importance of genetic recombination in the emergence of new enterovirus variants.

Introduction

Human enterovirus 71 (EV71) is a small, single-stranded, positive-sense RNA virus with a genome length of approximately 7400 nucleotides. It belongs to the species of Enterovirus A (EV-A), genus Enterovirus, family Picornaviridae [1, 2]. EV71 is the main pathogen responsible for causing hand, foot, and mouth disease (HFMD), an exanthematous illness typically observed in children < 5-years-old. It has been associated with severe clinical symptoms, including encephalitis and poliomyelitis that could lead to death [3]. Based on the VP1 encoding region, EV71 strains have been divided into genogroups A–C [4]. Genogroup A comprises of only one strain, the BrCr prototype, which was initially isolated in California, USA, in 1970. Genogroups B and C are common and consist of sub-genogroups B0–B5 and C1–C5, respectively [5].

Genetic recombination is a common phenomenon essential for the evolution of most RNA viruses [6]. The first natural recombinant RNA virus was reported in the poliovirus in the 1960s. Since then, numerous RNA virus recombinants have been reported, including picornaviruses, coronaviruses, and alphaviruses [7-10]. Previous studies on enterovirus have also suggested extensive recombination between and within enterovirus [11-14]. The recombination events of enterovirus have been observed in many countries, such as Japan, Malaysia, Singapore, Thailand, Netherlands, and China [15-23]. Also, the strains SZ/HK08-5 and SZ/HK08-6 infection of EV71 circulating in China have been speculated to be prevalent due to the emergence of this ‘‘double recombinant’’ strain resulting from the recombination between EV71 genotypes B and C, CA16 strain G-10 and EV71 genotype B [24]. Moreover, it is believed that the emergence of the EV71 triple-recombinant Chinese genotype A strain is due to recombination events among viruses SHZH98, C1, and CVA5 [25].

In mainland China, HFMD was first reported in Shanghai in 1981 [26]. The first severe outbreak of EV71 infection occurred in Linyi, Shandong province in 2007 [15] and then in Fuyang, Anhui province after a year [16]. The HFMD outbreak in Fuyang was caused by a recombinant EV71 strain [15], which indicated that new EV71 variants were potential epidemic strains. Thus, it is possible to generate recombination when at least two viruses infect the same cell under selective pressure from hosts and the environment [5]. The new recombinant variants may have acquired virulence determinants, leading to major fatal outbreaks observed in China over recent years [16, 19, 27]. The naturally occurring recombinant viruses based on homologous recombination are essential for a new epidemic situation. However, in the past decade, the characteristics and the reason for the emergence of new intertypic and intratypic recombinant viruses in the natural population of EV71 in mainland China are yet unknown.

In this study, we systematically investigated the genetic recombination among naturally circulating EV71 isolates in mainland China from 2009–2018, based on the publicly available EV71 genome sequences. The current study provided new clues for the surveillance of various EV71 variants, a prompt understanding of the evolution and epidemic of EV71, and also facilitated the surveillance and control of the virus.

Materials and methods

We collected nucleotide sequences of EV71 with the complete genome from NCBI GenBank database (https://www.ncbi.nlm.nih.gov/), and Virus Pathogen Database and Analysis Resource (ViPR) (https://www.viprbrc.org/brc/home.spg?decorator=vipr). Together, we retrieved 196 EV71 strains isolated from different parts of China from 2009 to 2018. These sequences were screened to exclude patented strains and artificial mutants. Consequently, a multiple sequence alignment was performed based on the complete genome sequences and various prototypes, and the oldest available EV-A strains were used as reference sequences for ClustalW2 (http://www.clustal.org/clustal2/) analysis. After removing gaps, a phylogenetic dendrogram was constructed, applying the neighbor-joining method based on the maximum composition likelihood model using MEGA7.0 [28]. RDP5 was used to identify natural recombinant strains within complete genome sequence alignments. RDP5 consists of seven detection methods: RDP (R), GENCONV (G), BootScan (B), MaxChi (M), Chimaera (C), SiScan (S), and 3Seq (T). If at least 3/7 detection methods showed a significant difference (P-value < 0.05), the sequence is considered as a potential recombination sequence [29]. The analysis was performed using SimPlot software package version 3.5.1 to further identify these recombinant strains. Similarity plots and bootscan analyses were used to verify the potential genetic recombination sites in the viral genome [30].

Results and discussion

Phylogenetic analysis of 196 complete genomes showed that strains were mostly clustered in the C4a subgenogroup (Fig. 1), with co-circulation of genotype A [25, 31]. The subgenotype C4b lineage included the oldest virus SHZH98/CHN/1998 strain found in China. The dominant EV71 strains circulating between 2009 and 2018 belonged to the subgenotype C4a. However, individual strains CQ2014-86/CQ/CHN/2014 and EV71/Xiamen/2009 were identified as subgenotype B5 [15]. Strains EV71-Hubei-09-China and EV71/Wuhan/3018/2010 clustered together with genotype A BrCr. In addition, 30-2/2015/BJ clustered into C2, and VR1432 clustered into subgenotype C5 lineage. Subsequently, a phylogenetic tree based on the entire VP1 coding region of EV71 was constructed (Fig. S1), and it was observed that the DL71 strain clustered into the subgenotype C2 lineage. Although the origin of these sporadic individual strains is unknown, it could be speculated that these were imported or resulted from the recombination of strains rather than evolution from the predominant strains of mainland China.

Fig. 1

The phylogenetic tree was based on complete genome sequences of 196 EV71 strains isolated from different parts of China from 2009 to 2018. Various prototypes and the oldest available EV-A strains were used as reference sequences. The phylogenetic dendrogram was constructed by the neighbor-joining method based on the maximum composition likelihood model using MEGA 7.0. The prototype CA16 G-10 strain was used as an outgroup. Bootstrap analysis was performed using 1000 replicates. Bootstrap values < 70% were not shown. The recombinant strains were represented by black dots

The recombination analysis of complete genome sequence alignments using RDP5 detected multiple potential recombinant strains based on 196 EV71 strains and other coxsackievirus (CA) reference strains. As a result, 3 intertypic recombination strains and 5 intratypic recombination strains were identified. The similarity plots and bootscan analyses were applied to confirm these recombinant events. The breakpoint positions and parent strains of the recombinant strains are listed in Table 1. Firstly, 3 intertypic recombinant strains (strain VR1432, 30-2/2015/BJ, and Guangdong-2009) were detected (Fig. 2). The strain VR1432 (accession number KC954664) belongs to the C5 subgenotype and 30-2/2015/BJ (accession number MG214681) belongs to the C2 subgenotype. Both strains have only one representative recombination breakpoint in the 2A nonstructural protein-encoding region. The strain Guangdong-2009 (accession number JF799986) belongs to C4a subgenotype with only one representative and recombinant breakpoint located in the 2C nonstructural protein region [32].

Table 1

Summary of enterovirus A recombination events detected by RDP5

Recombevent no	Breakpoint pos		Recombinant strain	Minor Parent	Major Parent	Detection methods and p-values
	Begin	End				RDP	GENECONV	Bootscan	Maxchi	Chimaera	SiSscan	3Seq
1	19a	3619	VR1432	Henan2-09-China	CA8-Donovan	6.34E−92	3.29E−86	2.45E−93	2.89E−33	3.70E−38	1.46E−77	1.01E−140
3	444	3575	30-2/2015/BJb	Huainan201401	CA8-Donovan	1.78E−61	9.55E−36	3.52E−55	6.08E−38	5.78E−03	3.43E−57	6.88E−86
5	16a	4124	Guangdong_2009	Anhui1-09-China	CA14-G-14	2.92E−80	2.80E−80	2.09E−80	2.33E−19	2.43E−19	2.07E−94	7.43E−06
9	5412	7094	EV71/P1034/2013	B3-MY821-3	Chongqing2-09–China	1.77E−21	NS	7.35E−18	2.59E−22	4.59E−25	5.95E−17	3.66E−50
17	1256	2439	Henan-ZMD/CHN/2012	Hubei-WH/CHN/2012	R615/YN/CHN/2010	2.15E−20	2.14E−21	NS	2.43E−04	7.75E−07	1.88E−10	3.46E−15
22	1596	2567	VR1432	Anhui1-09-China	C4b-SHZH98	2.38E−16	6.48E−21	2.84E−12	8.58E−12	6.73E−10	1.08E−22	9.65E−06
25	1930	3212a	120/EV71/Wenzhou/CHN/2014	SHAPHC5218/SH/CHN/14	Hubei-WH/CHN/2012	9.75E−17	1.28E−07	NS	5.34E−09	1.09E−09	8.47E−08	1.20E−15
49	3762a	5183	EV71/P868/2013/China	EV71/P40/2013/China	Puyang-451–2012	8.03E−09	2.95E−10	8.45E−10	3.69E−04	7.51E−03	6.61E−07	5.53E−07

A complete detailed overview of the recombination events detected by RDP can be found in Supplementary Table S1

Minor Parent Parent contributing the smaller fraction of sequence, Major Parent parent contributing the larger fraction of sequence, NS No significant P-value was recorded for this recombination event using this method

aThe actual breakpoint position is undetermined (it was most likely overprinted by a subsequent recombination event)

bThe recombinant sequence may have been misidentified (one of the identified parents might be the recombinant)

Fig. 2

Identification of intertypic recombination strains of EV71 by similarity plots and bootscan analyses. Similarity plots were conducted in a sliding window size of 200 nucleotides with 20-bp steps using the Kimura 2-parameter distance method; bootscan analyses were performed using the neighbor-joining tree model and the Kimura 2-parameter distance algorithm in a sliding window size of 200 nucleotides with 20-bp steps

For the recombinant strain VR1432 (accession number KC954664), the putative parent strain, Donovan/USA/1949 (accession number AY421766), was first isolated in the USA in 1949 (the strain belongs to the CA8 serotype). Its minor parent was Henan2-09-China (accession number GQ994992) isolated from Henan, China, in 2009 (the strain belongs to the C4a subgenotype). However, the initial breakpoint position of VR1432 has not yet been determined, while the terminal breakpoint position was located at 3619, as detected by RDP5. Simplot verified the breakpoint of the recombinant strain VR1432 located in the 2A gene region (3329–3379 nt), and according to the similarity plot and bootscan analysis, we speculated that the major parent is Henan2-09 -China and the minor parent is CA8.

The recombinant strain 30-2/2015/BJ (accession number MG214681) was isolated in 2015 in Beijing, China. Its minor parent strain was Donovan/USA/1949, derived from the classic strain reference sequence CA8 serotype, while the major parent strain, Huainan201401 (accession number MG773122), was isolated from Huainan, China, in 2014. The recombination breakpoint of strain 30-2/2015/BJ was located in the 2A gene region (3329-3379nt) based on Simplot analysis. The recombination region was located in the P2 and P3 regions. These two intertypic recombinant strains have similar breakpoints in the recombination region. The genetic recombination was observed between CVA8 serotype and C2, C5 subgenotype of EV71 in the nonstructural protein region 2A in this analysis. The breakpoint at 2A region is a recombination hot spot in enterovirus [5]. Based on previous reports [12, 33, 34], we speculated that the 2A region for the CA8 serotype might be prone to recombination events with respect to gene fragment exchange.

The third intertypic recombination, strain Guangdong-2009 (accession number JF799986) belonged to the C4a subgenotype and was isolated in Guangdong, China, in 2009. Its major parent strain, G-14/SOA/1950 (accession number AY421769), belongs to genotype CA14, isolated from South Africa in 1950, while its predicted minor parent, the Anhui1-09-China strain (accession number GQ994988), was isolated from Anhui, China, in 2009, and belongs to the C4a subgenotype. The similarity plot and bootscan analysis revealed the recombination breakpoint at 4200–4300 in the 2C region, which was consistent with the data analyzed by RDP5. Based on the similarity plot and bootscan analysis, we speculated that the major parent is Anhui1-09-China and is the minor parent is CA14 (Fig. 2).

In this study, we identified 5 intratypic recombinant strains (EV71/P1034/2013, Henan-ZMD/CHN/2012, VR1432, 120/EV71/Wenzhou/CHN/2014, and EV71/P868/2013/China) (Fig. 3). For recombinant strain EV71/P1034/2013 (accession number KP289419), the putative major parent was Chongqing2-09-China (accession number GQ994990) isolated from Chongqing, China (the strains belong to a C4a subgenotype), and the minor parent strain was MY821-3/1997 (accession number DQ341367), which belongs to the B3 subgenotype that originated in Singapore in 1997. The initial actual breakpoint position of EV71/P1034/2013 was located at 5412, and the terminal breakpoint position was located at 7094. The breakpoint regions were located at 3C of the P3 nonstructural protein region, in which minor parent MY821-3/1997 shared greater similarity with recombinant strain EV71/P1034/2013. The identical evidence was provided by the bootscan result (Fig. 3).

Fig. 3

Identification of intratypic recombination strains of EV71 by similarity plots and bootscan analyses. The similarity plot was conducted in a sliding window size of 200 nucleotides with 20-bp steps using the Kimura 2-parameter distance method; bootscan analysis was performed using the neighbor-joining tree model and the Kimura 2-parameter distance algorithm in a sliding window size of 200 nucleotides with 20-bp steps, applying the first recombinant strain EV71/P1034/2013. Similarity plots were conducted in a sliding window size of 1000 nucleotides with 20-bp steps using the Kimura 2-parameter distance method; bootscan analyses were performed using the neighbor-joining tree model and the Kimura 2-parameter distance algorithm in a sliding window size of 1000 nucleotides with 20-bp steps, applying the other recombinant strains

For recombinant strain, Henan-ZMD/CHN/2012 (accession KP198624), the major parent R615/YN/CHN/2010(MF662684) was isolated from Yunnan, China, and the strain belongs to the C4a subgenotype, while the minor parent, Hubei-WH/CHN/2012 (accession KP198623), was isolated from Wuhan, China, and belongs to the C4a subgenotype. The recombinant region of the strain was detected in the VP3 region. The initial actual breakpoint position of Henan-ZMD/CHN/2012 was located at 1256, and the terminal breakpoint position was located at 2439 (Fig. 3).

For recombinant strain VR1432, the putative major and minor parent strains were C4b-SHZH98 (accession number AF302996) and Anhui1-09-China (accession number GQ994988), respectively. The recombinant regions were detected in the VP3 and VP1 regions. The actual initial breakpoint position of VR1432 was located at 1596, and the terminal breakpoint position was located at 2567, and the recombinant region was located in VP1–VP3 of the P1 region (Fig. 3). Based on the results of recombination analysis, it was determined that the strain VR1432 is “double recombinant” and involves both intratypic and intertypic recombination events. Interestingly, in the upstream region of the genome, the recombinant strain VR1432 was similar to the subgenotype C4b, while downstream of the genome, a high degree of nucleotide identity with serotype CA8 was established.

For recombinant strain 120/EV71/Wenzhou/CHN/2014 (accession KT345959), the major parent Hubei-WH/CHN/2012 (accession KP198623) strain was isolated from Wuhan, China in 2012 (belongs to C4a subgenotype), and the minor parent was SHAPHC5218/SH/CHN/14 (accession KU936120), isolated from Shanghai, China, in 2014 (belongs to the C4a subgenotype). The recombinant region was determined in the VP3 and VP1 regions. The actual initial breakpoint position of 120/EV71/Wenzhou/CHN/2014 was located in 1930, and the terminal breakpoint position was located at 3212 (Fig. 3).

For recombinant strain EV71/P868/2013/China (accession KP289430), the major parent, Puyang-251-2012 (accession KM211579), was isolated from Puyang, China in 2012, while the minor parent, EV71/P40/2013/China (accession KP289426), was isolated from China in 2013. Both parent strains belong to the C4a subgenotype, and the recombinant region was detected at the P2 region. The initial actual breakpoint position of EV71/P868/2013/China was located at 3762, and the terminal breakpoint position was located at 5183 (Fig. 3).

Recombinant data indicated that using complete genome analysis or sequencing of multiple regions of EV71 for circulating strains is an effective method for viral genetics research. Some publications reported that the high nucleotide sequence similarity in nonstructural regions (P2 and P3) might favor intertypic recombination across homologous regions in poliovirus [35, 36]. Intertypic recombination of CV-A16 and EV-A71 genotype C4 has been reported in Central China [16, 24]. Intertypic recombination enables EV71 to obtain the gene fragment from the other types of viruses; consequently, a highly pathogenic EV71 strain with increased virulence is generated through natural recombination [5, 16, 19, 27, 37–39], which facilitates the recombinant virus to adapt to the new epidemic environment. Notably, the nonstructural protein regions had a high probability of recombination as described previously [16, 22], and that the P1 region encoding structural protein is a vital region for intratypic recombination. Intratypic recombination of the C4a subgenotype was observed in the VP1 region [40]. The VP1 capsid region has been reported to possess several antigenic determinants and induce neutralizing antibodies. The VP1 encodes the capsid protein with the maximal variation, which allows the virus to escape neutralizing antibodies. Thus, the antigenicity was changed to adapt to the subsequent prevalence of HFMD.

The current study indicates that the genome of EV71 undergoes recombination, which might be a key factor in determining the high virulence of some strains. However, the exact mechanism of the recombination events observed in this study is not yet clarified. To understand the recombination process, it will be necessary to develop the experimental model of in vitro selection of recombinant viruses. In addition to recombined strains introduced to China, C4a circulating strains may also be exported to other countries [23, 41]. Local strains and C4a subtypes from other countries might simultaneously go through intra- and inter-genotype recombination events, which produce new pandemic strains. Therefore, the surveillance of the emergence and epidemic of the EV71 recombinant virus might prevent the high epidemic of HFMD.

In conclusion, in this study, we identified 3 intertypic recombinants and 5 intratypic recombinants, including one “double-recombinant” strain based on the analysis of 196 complete genomes of the EV71 in mainland China during 2009–2018. The results indicated that homologous recombination contributes to the emergence of various pathogenic variants. Recombination is not only the evolutionary mechanism of EV71 but also a common pattern of genetic plasticity that leads to viral diversity [42]. Recombination between circulating strains of EV71 and coxsackievirus A of different genotypes may play a role in the evolution of the strains by improving the overall fitness of the virus, such as increased virulence and/or transmissibility [16, 19], which requires further clinical monitoring.

Taken together, this study improves the understanding of the evolution and epidemic of EV71 and facilitates the surveillance and control of the virus. Although there is an EV71 vaccine with a C4 genotype in China, understanding the virus may facilitate the development of a broad-spectrum vaccine that can be used against different genotypes. Therefore, continuous surveillance of EV71 is needed to elucidate its epidemiology and virus evolution.

Below is the link to the electronic supplementary material.

Electronic supplementary material 1 (VCF 0 kb)

Electronic supplementary material 2 (PDF 173 kb) Fig. S1 The phylogenetic tree was based on VP1 sequences of 196 EV71 strains isolated from different parts of China from 2009 to 2018. Various prototypes and the oldest available EV-A strains were used as reference sequences. The phylogenetic dendrogram was constructed by the neighbor-joining method based on the maximum composition likelihood model using MEGA 7.0. The prototype CA16 G-10 strain was used as an outgroup. Bootstrap analysis was performed using 1000 replicates. Bootstrap values <70% were not shown.

Electronic supplementary material 3 (XLSX 26 kb) Table S1 196 EV71 strains were isolated from different parts of China from 2009 to 2018 (Appendix A), and reference strains were used for complete genome analysis (Appendix B).

Electronic supplementary material 4 (XLSX 88 kb) Table S2 The characteristics of detectable enterovirus A recombination events.