#PAGE_PARAMS# #ADS_HEAD_SCRIPTS# #MICRODATA#

Heparan sulfate attachment receptor is a major selection factor for attenuated enterovirus 71 mutants during cell culture adaptation


Authors: Kyousuke Kobayashi aff001;  Katsumi Mizuta aff002;  Satoshi Koike aff001
Authors place of work: Neurovirology Project, Department of Genome Medicine, Tokyo Metropolitan Institute of Medical Science, Tokyo, Japan aff001;  Department of Microbiology, Yamagata Prefectural Institute of Public Health, Yamagata, Japan aff002
Published in the journal: Heparan sulfate attachment receptor is a major selection factor for attenuated enterovirus 71 mutants during cell culture adaptation. PLoS Pathog 16(3): e1008428. doi:10.1371/journal.ppat.1008428
Category: Research Article
doi: https://doi.org/10.1371/journal.ppat.1008428

Summary

Enterovirus 71 (EV71) is a causative agent of hand, foot, and mouth disease (HFMD). However, this infection is sometimes associated with severe neurological complications. Identification of neurovirulence determinants is important to understand the pathogenesis of EV71. One of the problems in evaluating EV71 virulence is that its genome sequence changes rapidly during replication in cultured cells. The factors that induce rapid mutations in the EV71 genome in cultured cells are unclear. Here, we illustrate the population dynamics during adaptation to RD-A cells using EV71 strains isolated from HFMD patients. We identified a reproducible amino acid substitution from glutamic acid (E) to glycine (G) or glutamine (Q) in residue 145 of the VP1 protein (VP1-145) after adaptation to RD-A cells, which was associated with attenuation in human scavenger receptor B2 transgenic (hSCARB2 tg) mice. Because previous reports demonstrated that VP1-145G and Q mutants efficiently infect cultured cells by binding to heparan sulfate (HS), we hypothesized that HS expressed on the cell surface is a major factor for this selection. Supporting this hypothesis, selection of the VP1-145 mutant was prevented by depletion of HS and overexpression of hSCARB2 in RD-A cells. In addition, this mutation promotes the acquisition of secondary amino acid substitutions at various positions of the EV71 capsid to increase its fitness in cultured cells. These results indicate that attachment receptors, especially HS, are important factors for selection of VP1-145 mutants and subsequent capsid mutations. Moreover, we offer an efficient method for isolation and propagation of EV71 virulent strains with minimal selection pressure for attenuation.

Keywords:

Mammalian genomics – RNA viruses – Viral replication – Cell cultures – Viral genomics – Substitution mutation – Microbial mutation – Mutation detection

Introduction

Enterovirus 71 (EV71), which belongs to the genus Enterovirus of the family Picornaviridae, has a single-stranded positive-sense RNA genome encoding a single polyprotein that is flanked by two untranslated regions (UTR), a 5’UTR and a 3’UTR. Polyproteins are comprised of three regions, P1, P2, and P3. The EV71 icosahedral capsid is composed of 60 copies of the capsid proteins (VP1, VP2, VP3, and VP4), which are coded in the P1 region of the viral genome [1]. EV71 is a causative agent of hand, foot, and mouth disease (HFMD), together with other members of the human enterovirus A (EV-A) [2]. HFMD is a mild and self-limiting disease that manifests as vesicular lesions of the hands, feet, and mouth. However, unlike other EV-A serotypes, infection with EV71 is sometimes associated with severe neurological complications. Since the late 1990s, large outbreaks have repeatedly occurred in Asian countries [37]. It is still unclear what viral factors and host factors are critical for EV71 neurovirulence.

It is important to identify neurovirulence determinants in the viral genome to understand EV71 pathogenicity. It has been reported that an amino acid substitution at residue 145 of the VP1 protein (VP1-145) plays an important role in EV71 virulence. Viruses that have glutamic acid (E) at VP1-145 (VP1-145E viruses) are virulent in several animal infection models [815], whereas viruses that have glycine (G) or glutamine (Q) (VP1-145G and VP1-145Q viruses) are avirulent. However, most EV71 strains that are freshly isolated from infected patients have E at this position [16, 17], suggesting that mutation at this site alone is not sufficient to explain the difference in virulence among the circulating EV71 strains, and that additional determinants exist in the EV71 genome. However, other determinants are poorly reported. For example, amino acid mutations at VP1-170 [18], VP1-97 [19], 3D-73, 3D-362 [20], VP1-164, and 2A-68, [21], as well as nucleotide mutations at 5’UTR-158 [21], -272, -488, -700 [22], -485, and 3’UTR-7408 [20], have been reported as virulence determinants of EV71.

One of the main issues in evaluating EV71 virulence determinants is that the genomic sequence of EV71 changes very rapidly during its isolation and propagation in cultured cells, such as RD-A cells and Vero cells [23, 24]. EV71 freshly isolated from human samples grow very poorly in cultured cells. Nucleotide and amino acid substitutions can occur after only a few passages during the isolation of the virus from human specimens or during recovery of the virus from the infectious cDNA clone, which makes it difficult to obtain accurate results even using reverse genetics. The rapid change in the virus sequence suggests that a strong selection pressure is imposed on the viral population in cultured cells. However, it is unclear what factors cause the selection of certain mutations, and detailed information on how the EV71 population changes in culture is lacking. In the case of RNA viruses, including EV71, replication errors in the viral genomic RNA are induced by their low-fidelity RNA-dependent RNA polymerases. This machinery diversifies the RNA virus population, which plays an important role in survival and adaptation to various environments [25]. Selection pressures on the virus population in host animals are thought to be quite different from selection pressures in cultured cells. To replicate in vivo, viruses must overcome various setbacks in a variety of tissues or cells and must escape detection by the host immune system. Through the mutation and selection process, tissue/cell-specific virus populations are generated, as reported for EV71 and poliovirus [19, 26, 27]. Viruses isolated from clinical samples of infected patients should therefore be adapted to various in vivo environments. However, the isolated viruses are usually propagated for characterization in clonal cell cultures, in the absence of an acquired immune system. During this process, cell culture adaptation occurs. It is unclear how this environmental change affects viral populations in vitro. Analyzing the population dynamics of viruses during in vitro adaptation should reveal useful information about the factors that drive adaptive selection of the virus, and how unwanted adaptation in cultured cells could be avoided.

Virus receptors are considered to be one of the selection pressures for virus selection. EV71 infection is initiated by attachment of the virus to the cell surface, followed by its internalization and the release of viral genomic RNA into the cytoplasm of infected cells, a process called uncoating. We previously reported that human scavenger receptor class B, member 2 (hSCARB2) can support these three steps [28]. All EV71 strains can use hSCARB2 as a receptor [29]. hSCARB2 transgenic (tg) mice are susceptible to EV71 infection, and EV71-infected mice show neurological disease [30]. hSCARB2 binds the south rim of the canyon of the EV71 virion [31], and this binding initiates uncoating at a low pH [30]. However, SCARB2 is a lysosomal protein and is not abundantly expressed on the surface of cultured cells. Therefore, this step can be a bottleneck on EV71 replication. Some EV71 strains also use so-called attachment receptors, including P-selectin glycoprotein ligand-1 (PSGL-1) [32], heparan sulfate (HS) [33], annexin II [34], sialic acid [35], nucleolin [36], vimentin [37], and fibronectin [38]. The attachment receptors can bind to the virus at the cell surface and enhance infection, although attachment receptors alone are not sufficient for establishment of infection because they cannot initiate uncoating of the virion. The amino acid residues near the five-fold axis, which includes VP1-145, determine binding specificity to HS and PSGL-1 [13, 24, 33]. The surface of the VP1-145G and VP1-145Q virion around the five-fold axis is rich in positively charged amino acids [39], allowing for electrostatic interaction with HS and highly sulfated PSGL-1. The negative charge of the E residue at VP1-145 neutralizes the positive surface charge, resulting in decreased affinity to HS and PSGL-1 [24, 39]. The binding specificity of EV71 to other attachment receptors has not been elucidated in detail.

We hypothesized that attachment receptors play an important role in the selection of viral populations in vitro. Here, we analyzed the population dynamics of EV71 strains that were initially virulent in vivo during cell culture adaptation. We found that EV71, which acquired a mutation in VP1-145, was effectively selected in cultured cells. This mutation caused attenuation of virulent strains. We hypothesized that HS expressed on the cell surface is a major factor for this selection in RD-A cells. We confirmed this hypothesis using HS-deficient, hSCARB2-overexpressing cells. In addition, this mutation further promotes the acquisition of secondary mutations in the EV71 capsid to increase the fitness of the virus in cultured cells. We propose that attachment receptor usage is a major factor for adaptation of EV71 in vitro.

Results

Viruses carrying VP1-145G and Q mutations selectively replicate during passage of EV71 strains in RD-A cells

RD-A cells and Vero cells are commonly used for the isolation of EV71 from clinical samples of HFMD patients. During the isolation process, EV71 replicates very poorly. A cytopathic effect (CPE) does not always appear immediately, but sometimes appears after several blind passages. In addition, virus stocks prepared from a virus with a low passage history tend to have a low virus titer. However, after several passages, the CPE appears much earlier and the viral titer reaches high levels. These results suggested that the cell culture conditions required for viral proliferation are quite different from those in vivo and that virus fitness under cell culture conditions is very low, indicating that adaptation and selection of the adapted virus must occur to overcome this low fitness during this process.

To identify the mutations selected in cultured cells, we analyzed single nucleotide variations (SNVs) occurring in the EV71 genome after passage in RD-A cells. The 2716-Yamagata-03 (2716-Ymg-03) strain, which is classified into subgenogroup B5, was isolated from an HFMD patient using GMK cells, passaged two generations [16], and passaged one generation in RD-A-overexpressing hSCARB2 (RD+hSCARB2) cells. This stock was used as the starting material (passage-0; p-0) for this experiment. The SI/Isehara/Japan/99 (Isehara) strain, which is classified into subgenogroup C2, was isolated from an HFMD patient and passaged several times before we received it. Although the passage history of this virus is not clear, we constructed an infectious cDNA clone for this strain [13], prepared the virus from the cloned cDNA after transfection of the in vitro-transcribed RNA into RD+hSCARB2 cells, and then propagated it by two passages in RD+hSCARB2 cells. This was used as p-0 virus. These viruses were passaged in RD-A cells three times at a low multiplicity of infection (MOI < 0.01). Subsequent generations were designated as p-1, p-2, and p-3. To evaluate mutations arising in each virus population, we performed next-generation sequencing (NGS) on viral RNA from the p-0, p-1, and p-3 virus samples, mapped the reads to the consensus sequences of the p-0 virus, and estimated the abundance of each SNV (Fig 1A and 1B).

SNVs and virulence of EV71 populations passaged in RD-A cells.
Fig. 1. SNVs and virulence of EV71 populations passaged in RD-A cells.
2716-Ymg-03 (A, C–E, G–I) and Isehara (B, F, J) strains were passaged in RD-A cells (A–B, E–F, I–J), Vero cells (C, G), and HeLa+hSCARB2 (D, H) for three generations. (A–D) To detect SNVs in the passaged virus populations, NGS reads were mapped to the consensus sequence of the p-0 virus and analyzed using LoFreq and SnpEff software. The abundance of each annotated SNV was plotted on the corresponding position in the virus genome. Vertical dotted lines indicate borders separating genomic regions (5’UTR, P1, P2, P3, and 3’UTR). Horizontal dotted lines indicate the 1.5 interquartile ranges (IQR) of each SNV mapped to the graph. Arrowheads indicate the position corresponding to VP1-145. (E–H) The abundance of VP1-145 mutants in each virus population before and after passage, estimated from NGS data using ShoRAH software. (I–J) The passaged viruses were used to inoculate hSCARB2 tg mice [106 TCID50/mouse for the 2716-Ymg-03 strain, 105 TCID50/mouse for the Isehara strain; intraperitoneal (ip) injection]. Paralysis and death of the infected mice were observed over 14 days.

In p-1 and p-3 2716-Ymg-03 populations, we identified four mutations with more than 5% abundance at nucleotide positions 186, 2875, 2876, and 6517 of the virus genome (Fig 1A). The nonsynonymous mutations at 2875 and 2876 are located in the same codon, corresponding to VP1-145, where the former is an E-to-Q amino acid substitution and the latter is an E-to-G substitution. In the p-0 population, the VP1-145G mutant was detected at an abundance of approximately 5% and the VP1-145E at 95%, but VP1-145Q was not detected (Fig 1C). After one passage, the VP1-145G virus displaced the VP1-145E population (64.5%), and VP1-145Q appeared (15%). The virus population was completely composed of VP1-145Q virus (100%) after two more passages in RD-A cells (p-3). To exclude the possibility of a double mutation in the same codon (2875 and 2876), which would encode arginine (R), we analyzed the region from 2800 to 2910 by local haplotype analysis (Fig 1C) and detected no E-to-R substitution at VP1-145. The other nonsynonymous mutation at 6517 substituted histidine to tyrosine at amino acid residue 193 of 3D polymerase. This mutation also increased during passage (0%, 8.3%, and 20.9% in p-0, p-1, and p-3, respectively). The deletion of C at position 186 in the 5’UTR was detected in p-0 at an abundance of approximately 12% and was increased in p-1 (59.3%), but it was not detected in p-3.

In the passage experiment with the Isehara strain, two nonsynonymous mutations were detected with more than 5% abundance in the p-1 population, at nucleotide positions 2872 and 4051 (Fig 1B). These mutations lead to amino acid substitutions at VP1-145 and 2B-91, respectively. In the p-3 population, we found four nonsynonymous mutations, at nucleotide positions 1373, 2793, 2872, and 3184, all of which are located in the P1 region and lead to amino acid substitutions at VP2-141, VP1-119, VP1-145, and VP1-249, respectively. The E-to-Q amino acid substitution at VP1-145 was the only mutation detected in both the 2716-Ymg-03 and Isehara populations at a high abundance. In the p-0 population of the Isehara strain, no VP1-145 mutations were detected (Fig 1D). However, after just one passage in RD-A cells, VP1-145Q virus became dominant (92.8% and 99.6% in the p-1 and p-3 populations, respectively). Other mutations also increased during passage but did not reach 100% abundance in the population. Together, VP1-145G and Q mutations were the main mutations in both the 2716-Ymg-03 and Isehara strains, which were reproducibly enriched during passage in RD-A cells. The results suggest that strong selection pressure acts on this site.

To confirm whether an adaptive mutation at VP1-145 and was selected for in other cells, we passaged the 2716-Ymg-03 strain in Vero and HeLa-S3 cells. Because of the low susceptibility of HeLa-S3 cells to EV71, we used HeLa-S3 overexpressing flag-tagged hSCARB2 (HeLa+hSCARB2) for EV71 passages. In the 2716-Ymg-03 strain passaged in Vero cells three times, we detected nearly complete displacement by the VP1-145G mutant (Fig 1C–1G). In the case of HeLa+hSCARB2, the VP1-145Q mutant predominated in the virus population after three passages (Fig 1D–1H). The results suggested that the VP1-145 mutation occurs commonly in other cultured cells as well.

Cell-adapted viruses are attenuated

To evaluate the virulence of EV71 after cell culture adaptation, the parental viruses (p-0) and the passaged viruses (p-1 and p-3) of the 2716-Ymg-03 and Isehara strains in RD-A cells were tested for virulence using hSCARB2 tg mice (Fig 1I and 1J). The parental viruses of the both strains were lethal to a subset of mice under these inoculation conditions (50% and 70%, respectively). However, after just one passage in RD-A cells, neither of the strains were lethal in mice. These data suggest that avirulent mutants were rapidly selected in the cultured cells.

EV71 replicates in HS-deficient and hSCARB2-overexpressing cells without mutation at VP1-145 or attenuation

Surprisingly, the viral population was almost completely replaced by VP1-145G or Q mutants after only a few passages (Fig 1). These results led us to hypothesize that the infection efficiency of VP1-145E in cultured cells is restricted because of the low surface expression of SCARB2. Once adaptive mutants emerge through genome replication by low fidelity RNA-dependent RNA polymerases, they immediately become dominant in the virus population. In addition to this, we hypothesized that HS expressed on the surface of cultured cells plays an important role in the selection of these mutants. It is already known that VP1-145 mutation affects the binding of EV71 virions to HS and PSGL-1 but does not significantly affect the binding affinity to hSCARB2 [13, 24, 33, 39]. Many types of cells, including RD-A, Vero, and HeLa-S3 cells, express HS on the cell surface (Fig 2A), but PSGL-1 expression is not observed in this cell line [32]. Our hypothesis is that HS-binding mutants expressing G or Q at VP1-145 more efficiently infect cultured cells by binding to HS on the cell surface, resulting in the selection of the VP1-145G and VP1-145Q mutants.

HS and SCARB2 expression in genetically modified RD-A cells.
Fig. 2. HS and SCARB2 expression in genetically modified RD-A cells.
(A–B) FACS analysis of HS expression at the cell surface. (C) Western blotting analysis using anti-SCARB2 (top), anti-flag (middle), and anti-ß-actin (bottom) antibodies. The arrowheads indicate hSCARB2.

To confirm this possibility, we generated RD-A cells lacking HS by knocking out the EXT1 or EXT2 genes (RDΔEXT1 or RDΔEXT2), which are involved in elongation of the HS chain. In the knockout cell lines, no expression of HS on the cell surface was observed (Fig 2B). Then, flag-tagged hSCARB2 was introduced into wild-type RD-A (RD+hSCARB2) and HS-deficient RD-A cells (RDΔEXT1+hSCARB2 and RDΔEXT2+hSCARB2). The expression of flag-tagged hSCARB2 was confirmed by western blotting (Fig 2C). These cells are referred to as ΔEXT1+hSCARB2 and ΔEXT2+hSCARB2, respectively, for simplicity. The 2716-Ymg-03 strain was passaged in these cells (Fig 3). Selection of VP1-145G and VP1-145Q mutants in the 2716-Ymg-03 virus population passaged in ΔEXT1 and ΔEXT2 cells was observed (Fig 3A, 3B, 3F and 3G), similar to our observations in RD-A cells (Figs 1A and 3C). The VP1-145G mutant was selected in RD+hSCARB2 cells (Fig 3C and 3H). By contrast, no increase in VP1-145G or VP1-145Q mutants was observed when cells were passaged in the ΔEXT1+hSCARB2 and ΔEXT2+hSCARB2 cells (Fig 3D, 3E, 3I and 3J). This result demonstrates that expressions of HS and hSCARB2 are involved in the rapid selection of VP1-145 mutants in vitro.

The 2716-Ymg-03 strain replicates in HS-deficient and hSCARB2-overexpressing cells without attenuation or mutation in VP1-145.
Fig. 3. The 2716-Ymg-03 strain replicates in HS-deficient and hSCARB2-overexpressing cells without attenuation or mutation in VP1-145.
The 2716-Ymg-03 strain was passaged in genetically modified RD-A cells. (A–E) To detect SNVs in passaged virus populations, NGS reads were mapped to the consensus sequence of the p-0 virus and analyzed using LoFreq and SnpEff software. The abundance of each annotated SNV was plotted on the corresponding position in the virus genome. Vertical dotted lines indicate borders separating genetic regions (5’UTR, P1, P2, P3, and 3’UTR). Horizontal dotted lines indicate the 1.5 IQR of the abundance of each SNV mapped to the graph. Arrowheads indicate the position corresponding to VP1-145. (F–J) The abundance of VP1-145 mutants in the virus population before and after passage, estimated from NGS data using ShoRAH software. (K–O) The passaged 2716-Ymg-03 strains were used to infect hSCARB2 tg mice (106 TCID50, ip). Paralysis and death of the infected mice were observed over 14 days.

The virulence of the passaged 2716-Ymg-03 virus populations was tested using hSCARB2 tg mice (Fig 3K–3O). Although the 2716-Ymg-03 populations passaged in ΔEXT1, ΔEXT2, or RD+hSCARB2 cells became avirulent (Fig 3K–3M), similarly to the viruses passaged in RD-A cells (Fig 1E), the virus populations passaged in ΔEXT1+hSCARB2 or ΔEXT2+hSCARB2 cells remained virulent even after three passages (Fig 3N and 3O).

When the Isehara strain was passaged in these cells, the VP1-145Q mutant was selected in ΔEXT1, ΔEXT2, and RD+hSCARB2 cells (Fig 4A–4C and 4F–4H), as was observed in the RD-A cells (Fig 1B and 1D). However, in ΔEXT1+hSCARB2 and ΔEXT2+hSCARB2 cells, VP1-145E remained dominant even after three passages (Fig 4D, 4E, 4I and 4J). The virus populations passaged in ΔEXT1, ΔEXT2, and RD+hSCARB2 cells became avirulent in mice (Fig 4K–4M), but those passaged in ΔEXT1+hSCARB2 and ΔEXT2+hSCARB2 cells remained virulent (Fig 4N and 4O). These results demonstrate that EV71 replicates in HS-deficient and hSCARB2-overexpressing cells without E-to-G or E-to-Q mutation in VP1-145 or attenuation, but that either HS deletion or hSCARB2 overexpression alone is insufficient to avoid mutation.

The Isehara strain replicates in HS-deficient and hSCARB2-overexpressing cells without attenuation or mutation in VP1-145.
Fig. 4. The Isehara strain replicates in HS-deficient and hSCARB2-overexpressing cells without attenuation or mutation in VP1-145.
The Isehara strain was passaged in genetically modified RD-A cells. (A–E) To detect SNVs in the passaged virus populations, NGS reads were mapped to the consensus sequence of the p-0 virus and analyzed using LoFreq and SnpEff software. The abundance of each annotated SNV was plotted on the corresponding position in the virus genome. Vertical dotted lines indicate borders separating genetic regions (5’UTR, P1, P2, P3, and 3’UTR). Horizontal dotted lines indicate the 1.5 IQR of the abundance of each of the SNVs mapped to the graph. Arrowheads indicate the position corresponding to VP1-145. (F–J) The abundance of VP1-145 mutants in the virus population before and after passage, estimated from NGS data using ShoRAH software. (K–O) Passaged Isehara strains were used to infect hSCARB2 tg mice (105 TCID50, ip). Paralysis and death of the infected mice were observed over 14 days.

VP1-145 mutations occurred and were selected in infectious cDNA-derived virus by passage in RD-A cells

To confirm the hypothesis suggested above, we next investigated how mutations accumulated in the viral genome of a homogenous population. In the experiments described above, we used p-0 viruses that were not completely homogenous. To minimize the effects caused by the heterogeneity of the starting population, experiments were conducted using in vitro-transcribed RNA, which is theoretically quite homogenous. Viral genomic RNA of Isehara strains expressing VP1-145E (Isehara-E), VP1-145G (Isehara-G), or VP1-145Q (Isehara-Q) transcribed from cloned cDNA was transfected into RD-A and ΔEXT1+hSCARB2 cells, and serially passaged three times (Fig 5). For each virus, five independent experiments were performed to assess the reproducibility of selected mutations and the stochastic effect of replication errors in the EV71 genome.

Genome stability and virulence in EV71 populations with or without VP1-145 mutation that were prepared from cDNA clones.
Fig. 5. Genome stability and virulence in EV71 populations with or without VP1-145 mutation that were prepared from cDNA clones.
Viruses recovered from transfection of RD-A or ΔEXT1+hSCARB2 cells with in vitro-transcribed RNA from infectious cDNA of Isehara-E, -G, and -Q were passaged in these cells three times. (A–F) The abundance of VP1-145 mutants in the virus population before and after passage, estimated from NGS data using ShoRAH software. VP1-145E-dominant populations resulting from three passages in ΔEXT1+hSCARB2 cells, VP1-145G-dominant populations passaged in RD-A cells, VP1-145G-dominant populations passaged in ΔEXT1+hSCARB2 cells, VP1-145Q-dominant populations passaged in RD-A cells, and VP1-145Q-dominant populations passaged in ΔEXT1+hSCARB2 (analyzed in Figs 6 and 7 and Table 2) are indicated with *, †, ‡, §, and¶, respectively. (G–L) To detect SNVs in passaged virus populations, NGS reads were mapped to the sequence of the cloned cDNA and analyzed using LoFreq and SnpEff software. The abundance of each annotated SNV was plotted on the corresponding position in the virus genome. Vertical dotted lines indicate borders separating genetic regions (5’UTR, P1, P2, P3, and 3’UTR). Horizontal dotted lines indicate the 1.5 IQR of the abundance of each SNV mapped to the graph. Arrowheads indicate the position corresponding to VP1-145. High frequency (>5%) mutations in the P1 region are boxed with a red dashed line. (M–R) Passaged viruses were used to infect hSCARB2 tg mice (105 TCID50, ip). Paralysis and death of the infected mice were observed over 14 days.

When we prepared Isehara-E p-0 population in RD-A cells, the cells that had probably received viral RNA showed CPE one day after transfection but failed to go completion within seven days in all five independent experiments. As the passage number increased, CPE developed much faster. By contrast, Isehara-G and Isehara-Q showed full CPE two days after transfection to prepare p-0 populations, suggesting that an almost uniform population of VP1-145E virus can hardly infect RD-A cells. When Isehara-E virus was passaged in RD-A cells, VP1-145E virus was diminished by passage in all experiments (Fig 5A). The VP1-145G mutant became dominant in two experiments, while the VP1-145Q mutant became dominant in the other three experiments. In addition, a large number of other nonsynonymous mutations were observed in the P1 region (the area boxed with a red dashed line in Fig 5G). When Isehara-E virus was passaged in ΔEXT1+hSCARB2 cells, CPE appeared within three days after transfection and the following passages, which was much faster than the appearance of CPE in RD-A cells. VP1-145E virus remained dominant in all five experiments (Fig 5B). The abundance of additional nonsynonymous mutations in the P1 region was much lower than in experiments with RD-A cells (the area boxed with a red dashed line in Fig 5H). The result suggested that Isehara-E was unstable in RD-A cells but stable in ΔEXT1+hSCARB2 cells.

When Isehara-G virus was passaged in RD-A cells, the VP1-145G mutation was stable (Fig 5C). In one of five passage experiments with Isehara-G virus using ΔEXT1+hSCARB2 cells, we found an amino acid substitution at VP1-145 from G to E (Fig 5D). In the other four experiments, we found low or undetectable levels of VP1-145E in the virus population, which still predominantly contained the VP1-145G virus. In passage experiments with Isehara-Q virus using RD-A and ΔEXT1+hSCARB2 cells, mutation of VP1-145 was not detected in any experiment (Fig 5E and 5F). In experiments with Isehara-G and Isehara-Q viruses passaged in both cells, nonsynonymous mutations in the P1 region were frequently observed (Fig 5I–5L).

We tested the virulence of these passaged viruses in hSCARB2 tg mice (Fig 5M–5R). Isehara-E viruses passaged in ΔEXT1+hSCARB2 cells, which still had E at VP1-145, showed a virulent phenotype (Fig 5N), whereas those passaged in RD-A cells, which had G or Q at VP1-145, showed an avirulent phenotype (Fig 5M). Almost all the virus populations originating from Isehara-G and Isehara-Q viruses passaged in either cell type induced no symptoms under the same inoculation conditions (Fig 5O–5R). These results support the relationship between the amino acid at VP1-145 and virulence. However, in the fourth passage experiment with Isehara-G virus in ΔEXT1+hSCARB2 cells, although the virus predominantly had E at VP1-145 (93.9%), it was much less virulent (10% paralysis rate and 0% mortality rate) (Fig 5P) than the other VP1-145E virus populations originating from Isehara-E virus (Fig 5N).

Infection efficiency of VP1-145 mutants in RD-A and ΔEXT1+hSCARB2 cells

Selection of certain mutants may result from the emergence of new mutants and competition between the mutants and the parental virus. To understand why VP1-145G and VP1-145Q viruses were rapidly selected in RD-A cells, we estimated the infection efficiency of VP1-145 mutants. In a viral population with a uniform genome sequence, the infection efficiency (E) can be defined using the viral titer (T) and the number of viral particles (P) in the virus sample using the following equation: E = T/P. We determined T and P of the viral populations by microtiter assays and quantitative PCR (qPCR) of the viral RNA genome, and estimated the infection efficiency of the VP1-145E, G, and Q viruses (EE, EG, and EQ, respectively). The virus populations chosen for analysis were the VP1-145E virus passaged in ΔEXT1+hSCARB2 cells (experiment 5), the VP1-145G virus passaged in RD-A cells (experiment 3), and the VP1-145Q virus passaged in RD-A cells (experiment 1), which are shown in Fig 5H–5, 5I–3 and 5K–1, respectively, as representative populations for each virus. EG [0.005 median tissue culture infectious dose (TCID50)/copy] and EQ (0.01 TCID50/copy) were 238 and 489 times higher than EE in RD-A cells (0.000021 TCID50/copy), respectively (Table 1). Thus, once the VP1-145G or Q virus appears in the population during replication, these mutants will infect RD-A cells with more than a 200-fold higher infection efficiency than the parental virus. This result explains, at least in part, the rapid selection of the VP1-145G and Q mutations in RD-A cells. On the other hand, in ΔEXT1+hSCARB2 cells, the infection efficiency of these viruses was notably increased; in particular, EE increased approximately 100-fold. The relative difference between EE and EG or EQ decreased by14–16-fold. This result indicates that selection pressure on the VP1-145G or VP1-145Q mutant is lower in ΔEXT1+hSCARB2 cells than in RD-A cells. Although these estimates might be affected by mutations at other genomic positions that modulate infection efficiency, they are consistent with the population analysis above. According to our previous report [13], the binding efficiency of VP1-145G virus to heparin-agarose is 100–300-fold higher than that of VP1-145E virus, suggesting that the difference in infection efficiency is dependent on the efficiency of viral binding affinity to HS.

Tab. 1. Relative infection efficiency of VP1-145 mutants in RD-A and ΔEXT1+hSCARB2 cells.
Relative infection efficiency of VP1-145 mutants in RD-A and ΔEXT1+hSCARB2 cells.

Secondary mutations were positively selected following E-to-G or E-to-Q mutation at VP1-145

In addition to VP1-145 mutations, we detected nonsynonymous mutations at a high abundance (>5%) that were concentrated in the P1 region (Fig 5G and 5I–5L). These secondary mutations were observed in almost all combinations of viruses and cells except for the Isehara-E virus passaged in ΔEXT1+hSCARB2 cells, where the original amino acid, VP1-145E, was maintained. In contrast to the P1 region, nonsynonymous mutations were not observed in the P2 and P3 regions. These results raised the possibility that once E-to-G or E-to-Q substitutions at VP1-145 occur (spontaneously or artificially by reverse genetics), secondary mutations in the P1 region are introduced and selected. To confirm this possibility, we first classified the p-3 virus populations shown in Fig 5 into three groups according to the VP1-145 amino acid dominant in the p-3 population, and then further classified them according to the cell line used for passage. The populations were classified as VP1-145E-dominant populations passaged in ΔEXT1+hSCARB2 cells, VP1-145G-dominant populations passaged in RD-A cells, VP1-145G-dominant populations passaged in ΔEXT1+hSCARB2 cells, VP1-145Q-dominant populations passaged in RD-A cells, or VP1-145Q-dominant populations passaged in ΔEXT1+hSCARB2 cells, and are indicated by *, †, ‡, §, and¶in Fig 5A–5F.

Nonsynonymous mutations detected in at least two independent experiments shown in Fig 5 at more than 5% abundance are shown in Table 2. Many other minor mutations that appeared only once at a low abundance were not subjected to analysis. We noticed some mutations that appeared preferentially in particular combinations of cells and viruses. In VP1-145G-dominant populations passaged in RD-A cells (indicated by † in Fig 5A and 5C), six mutations, at VP3-144, VP1-164, VP1-169, VP1-224, VP1-241, and VP1-246, fit the above criteria. Of these, three mutations, at VP3-144, VP1-224, and VP1-241, were selected in more than two virus populations of this group. Notably, A-to-V mutation at VP1-224 was detected in five passage experiments at a high abundance (26–87%). In VP1-145Q-dominant populations passaged in RD-A cells (indicated by § in Fig 5A and 5E), two secondary mutations were detected (Table 2). One of them, a P-to-A mutation at VP1-246, was observed in all five independent passage experiments at a high abundance (55–99%). To analyze the spatial distribution of these mutations, we mapped the secondary mutations on the pentameric structure of the EV71 capsid (Fig 6). VP1-224 was located surrounding the pocket factor binding site (Fig 6A). On the other hand, VP1-246 protruded from the virus capsid near the 5-fold axis (Fig 6B). Moreover, in VP1-145G-dominant populations in ΔEXT1+hSCARB2 cells (indicated by ‡ in Fig 5D), L-to-F mutation at VP1-169 was preferentially selected in four independent passage experiments at a high abundance (35–99%), which was different from the results seen in RD-A cells. This residue was not spatially adjacent to VP1-224, where the mutation was selected in RD-A cells (Fig 6C). In contrast to the above passage conditions, eight secondary mutations were observed in VP1-145Q-dominant populations passaged in ΔEXT1+hSCARB2 cells (indicated by ¶ in Fig 5F), but the abundance of these mutations was mostly lower than 50%, suggesting relatively low selection pressure under this passage condition. These results suggested that secondary mutations located at spatially distinct positions on the capsid were selected depending on a combination of the VP1-145 amino acid and the expression of HS and SCARB2 on the host cells, probably due to different selection pressures.

Distribution of secondary mutations induced in VP1-145G and Q mutants.
Fig. 6. Distribution of secondary mutations induced in VP1-145G and Q mutants.
Translucent surface representation of the pentameric structure of the EV71 capsid (Protein Data Bank: 4AED) viewed along the 5-fold axis from the outside. The amino acid residue VP1-145 is indicated in red. The amino acid residues of secondary mutations detected in VP1-145G-dominant populations passaged in RD-A cells (A), VP1-145Q-dominant populations passaged in RD-A cells (B), VP1-145G-dominant populations passaged in ΔEXT1+hSCARB2 cells (C), and VP1-145Q-dominant populations passaged in ΔEXT1+hSCARB2 cells (D), which are marked in Fig 5A–5F with †, ‡, §, and¶, respectively. Residues in which secondary mutations were detected in at least two virus populations in each group are represented by spheres. Pocket factors are indicated in magenta.
Tab. 2. List of secondary mutations detected at greater than 5% abundance in the P1 region of VP1-145G- and VP1-145Q-dominant populations.
List of secondary mutations detected at greater than 5% abundance in the P1 region of VP1-145G- and VP1-145Q-dominant populations.

To confirm whether the secondary mutations were positively selected, we further analyzed the nucleotide diversity (π) of select virus populations using the NGS data and SNPGENIE software [40]. The ratios between the nucleotide diversity at nonsynonymous sites (πN) and synonymous sites (πS) can be used to estimate the balance among neutral selection (πNS = 1), positive selection (πNS > 1), and negative selection (πNS < 1). Positive selection enhances nonsynonymous mutation to select more adaptive variants in an environment. In other words, if the fitness of the virus population is low, variants with nonsynonymous mutations conferring higher fitness will be selected, resulting in an increase in the πNS ratio. In Fig 7, the πNS ratios in the P1 region of VP1-145G- and VP1-145Q-dominant populations passaged in RD-A and ΔEXT1+hSCARB2 cells were higher than those of VP1-145E-dominant populations. This result is evidence of positive selection in the P1 region in association with VP1-145 mutation. On the other hand, we found no apparent correlation between VP1-145 mutation and πNS ratios in the P2 and P3 regions. This evidence supports the induction of secondary nonsynonymous mutations in the P1 region in association with VP1-145G and Q mutations.

VP1-145 mutation promotes secondary mutations in the P1 region.
Fig. 7. VP1-145 mutation promotes secondary mutations in the P1 region.
The virus populations shown in Fig 5 that predominantly (greater than 50%) contained VP1-145E, VP1-145G, or VP1-145Q mutants were further divided into two groups according to whether they arose from passage in RD-A cells or ΔEXT1+hSCARB2 cells. In each genetic region (P1, P2, and P3), the nucleotide diversity of synonymous and nonsynonymous sites (πS and πN, respectively) was calculated using SNPgenie software. πNS ratios are illustrated by scatter and box plots.

Mutation of VP1-145 exclusively correlated with EV71 virulence

We further tested the possibility that mutations other than VP1-145, including secondary mutations, contribute to EV71 attenuation. In all the codon positions where SNVs were detected (Fig 5G–5L), the Spearman’s correlation coefficient (rho) between the mutation rate and the rate of paralysis or mortality was determined and plotted their significance [-log10 (p-value)] at the corresponding codon position (Fig 8). Since mutations increasing in a virus population transitioning from a virulent to an avirulent phenotype should negatively correlate with virulence (rho < 0), we focused on negatively correlated mutations. We found that mutation at VP1-145 exclusively showed a significant negative correlation with EV71 virulence.

Correlation between SNVs and EV71 virulence.
Fig. 8. Correlation between SNVs and EV71 virulence.
SNVs were detected by mapping of the viral reads acquired in the passage experiments in Fig 5 to the consensus sequence of Isehara-E. Spearman’s correlation coefficients between the SNVs and the paralysis or mortality rate at all codon positions where SNVs were detected were calculated. Negative log-transformed p-values, which were calculated with the null hypothesis that there was no correlation between two factors, were plotted on the corresponding codon positions. Positive and negative correlations in red or blue, respectively. Horizontal dotted lines indicate p = 0.05.

Discussion

EV71 freshly isolated from clinical samples is difficult to grow in cell culture. To overcome this, the virus employs a strategy that allows it to survive under the adverse conditions of cell culture. In this study, we clarified the molecular mechanism involved in virus adaptation to cell culture. HS attachment receptor selected HS-binding viruses very rapidly during cell culture adaptation of EV71. Virus populations in clinical samples directly sequenced without propagation in cell culture are mostly composed of VP1-145E viruses [16, 17]. However, after 1–3 passages in wild-type, HS-deficient, or hSCARB2-overexpressing RD-A cells, VP1-145G or VP1-145Q mutants became dominant in all the populations that we tested. These results suggest that the fitness of the virulent EV71 population is low in the cell culture environment probably because of low surface expression of SCARB2 in cultured cells. Once more adaptive mutants with G or Q at VP1-145 arise through the error-prone RNA replication process, these mutants selectively replicate in cultured cells. Our data indicated that these two mutants were occurred by chance (Fig 5A) and we speculated that a mutant first emerging in the population should be selected as a dominant mutation. VP1-145G or Q mutation eventually becomes fixed in the population, resulting in increased fitness. VP1-145G or VP1-145Q viruses can infect RD-A cells with 200–500-fold higher efficiency than VP1-145E viruses (Table 1). This selection is apparently associated with the balance between the expression of HS and hSCARB2 in cultured cells because this selection was not observed in HS-deficient hSCARB2-overexpressing RD-A cells (Figs 35). We therefore propose that the major selection pressure during propagation of EV71 in cultured cells is the difference in infection efficiency among the VP1-145 mutants due to their different binding affinities to the HS attachment receptor.

During fixation from VP1-145E to VP1-145G or Q, the diversity of the viral population becomes theoretically low due to selective sweeping. However, we found a variety of nonsynonymous mutations in the P1 region that arose after or in parallel with the VP1-145 mutation (Fig 5G). These secondary mutations were also observed in the populations of VP1-145G and Q virus artificially generated by reverse genetics (Fig 5I–5L). Although the actual triggers for secondary mutations are unknown, we propose some possibilities. Some secondary mutations are beneficial for survival in cultured cells by enhancing affinity to receptor(s), altering receptor specificity, enhancing uncoating, stabilizing capsid conformation, and so on. Regardless of the mechanism, these mutations should increase the fitness of viruses expressing VP1-145G or Q in the cultured cells. In summary, the fitness of the virus during cell adaptation progresses through three stages (Fig 9). The VP1-145E virus has the lowest level of fitness (stage 1). After E-to-G or E-to-Q mutation at VP1-145, the mutant viruses have moderate levels of fitness (stage 2). Finally, viruses that acquire secondary mutations will have the highest levels of fitness (stage 3).

A model of the fitness landscape of virulent EV71 populations during cell culture adaptation in RD-A cells.
Fig. 9. A model of the fitness landscape of virulent EV71 populations during cell culture adaptation in RD-A cells.

Several previous studies revealed that VP1-145G and Q mutants replicate in cultured cells more rapidly than the VP1-145E virus by interacting with HS [13, 24]. However, deletion of HS is insufficient to prevent the selection of VP1-145G and Q (Figs 3 and 4), indicating that additional host factor(s) contributing to this selection still remain in HS-deficient cells. If the binding of VP1-145G and Q virus to HS is mediated mainly by electrostatic interactions, other negatively charged molecules expressed on the host cell surface may bind to these mutants and support their infection.

Our study demonstrates that attenuation of EV71 during cell culture adaptation is exclusively associated with the VP1-145 mutation (Fig 8), indicating that VP1-145E is adaptive in an in vivo environment but deleterious in an in vitro environment (stage 1 in Fig 9). Our previous study revealed that VP1-145G is avirulent in hSCARB2 tg mice and cynomolgus monkeys [13, 14]. One possible explanation for these observations is that, because the expression patterns of HS and hSCARB2 are quite different, VP1-145G and Q viruses are adsorbed onto cells with little or no hSCARB2 expression or extracellular matrix via HS, resulting in abortive infection and inefficient dissemination in vivo. We also found that the VP1-145G virus was more easily neutralized with neutralizing antibodies. These two mechanisms likely account for the negative selection of VP1-145G and Q viruses in vivo. Taken together with the findings of this study, these data suggest that HS is a major factor for positive selection of VP1-145G and Q mutants in vitro but for negative selection of these mutants in vivo. These opposing selection processes give rise to different virus populations in vivo and in vitro. This also implies that an unknown attachment receptor(s) specific for VP1-145E virus is expressed in EV71 target cells, such as neuronal cells, but is not expressed in cultured cells, such as RD-A cells.

Similar selections mediated by glycosaminoglycan including HS by cell culture adaptation have been reported in the viruses belonging to the Flaviviridae (Japanese encephalitis virus, Murray Valley encephalitis virus, West Nile virus, and Dengue virus) [4144], Togaviridae (Sindbis virus, Venezuelan equine encephalitis virus, Tick-borne encephalitis virus, and Chikungunya virus) [4548], and Picornaviridae (human Rhinovirus (HRV) C15, HRV89, and foot-and-mouth disease virus) [4951]. These selected viruses were attenuated in vivo. It is probable that the similar mechanisms of selection and attenuation observed in EV71 also occur in these viruses and that this principle is applicable in a wide range of viral infectious diseases.

Attenuation of viruses after tissue culture adaptation is known to contribute greatly to the development of some live attenuated vaccines, including poliovirus, measles virus, and yellow fever virus [52]. This study provides an example of an attenuation mechanism of viruses in tissue culture using EV71. Tsai et al. [53] used this mutation in a live attenuated EV71 vaccine candidate together with a mutation that increased the fidelity of the RNA-dependent RNA polymerase. This neurovirulence determinant may be useful for the development of live attenuated EV71 vaccines if the appearance of the revertant can be completely prevented.

ΔEXT1+hSCARB2 cells, in which HS was depleted and hSCARB2 was overexpressed, failed to support VP1-145G/Q selection. This newly developed cell line allows us to isolate and propagate EV71 strains circulating in human populations with minimal levels of mutation. This cell line will therefore be a useful tool for the identification of neurovirulence determinants of EV71 and for preparation of the challenge virus for efficacy tests of vaccines or antiviral drugs.

Materials and methods

Ethics statements

Experiments using recombinant DNA and pathogens were approved by the Committee for Experiments using Recombinant DNA and Pathogens at the Tokyo Metropolitan Institute of Medical Science. Experiments using mice were approved by the Animal Use and Care Committee and performed in accordance with the Guidelines for the Care and Use of Animals (Tokyo Metropolitan Institute of Medical Science, 2011, approval number: 17036, 18047, and 19063). In the animal infection experiments, mice with severe paralysis in limbs or more than 30% body weight loss were sacrificed by overdose of isoflurane or cervical dislocation.

Preparation of cell lines

The cells were propagated in Dulbecco’s modified Eagle’s medium (DMEM) containing 5% fetal calf serum (FCS). Human RD-A cells and African green monkey Vero cells were kindly provided by Dr. Hiroyuki Shimizu, National Institute of Infectious Diseases, Japan. In addition to these two cell lines, a suspension culture-adapted HeLa cell line (HeLa-S3) [54] was used for the passage of EV71 strains. ΔEXT1, ΔEXT2, and RD+hSCARB2 cells were described previously [13]. ΔEXT1+hSCARB2 and ΔEXT2+hSCARB2 cells were established using retrovirus vectors as described previously [13].

Viruses

The 2716-Ymg-03 strain was isolated as described previously [16]. Isolated 2716-Ymg-03 was propagated once using RD+hSCARB2 cells. Infectious cDNA from the Isehara-E and Isehara-G strains was described previously [13]. Infectious cDNA from the Isehara-Q strain was established using a previously described method [13]. Briefly, the E-to-Q mutation at amino acid position VP1-145 of the Isehara strain was introduced by amplifying the fragments using a pair of mutated primers (Ise VP1-145Q-2-F: 5’-CTA CCG GGC AAG TTG TCC CGC AAT T-3’, Ise VP1-145Q-2-R: 5’-GGG ACA ACT TGC CCG GTA GGC GTG C-3’) and the pSVA-EV71-Isehara plasmid as a template. The DNA fragment of the original plasmid was replaced with the corresponding mutated fragment. The parental Isehara strain used in Figs 1 and 4 was recovered by transfection of RD+hSCARB2 cells with RNA that had been transcribed in vitro from Isehara-E cDNA using a previously described method [13]. The recovered virus was propagated by two passages in RD+hSCARB2 cells. Parental stocks of the 2716-Ymg-03 and Isehara strains (p-0) were passaged in RD-A cells with or without genetic modifications at a low MOI (<0.01). Infected cells were cultured until approximately more than 80% of infected cells had died. For the passage experiment starting from infectious cDNA, in vitro-transcribed RNAs were transfected into each cell line. Recovered viruses were then passaged in the corresponding cell line. The virus titer was determined by microplate assay using RD-A cells, as previously described [13], and was represented as the TCID50 using the Kaeber method [55].

Animal experiments

Six-to-seven-week-old hSCARB2 tg mice [30] were used for infection experiments. After anesthetization by inhalation of isoflurane, a total of 500 μL of virus solution was inoculated intraperitoneally. A total of ten mice were used per group. The mice were monitored for body weight changes and clinical signs of disease for 2 weeks.

Flowcytometric analysis

Cells were detached by treatment with PBS containing 0.02% ethylenediaminetetraacetic acid (EDTA). Suspended cells were stained with a mouse anti-HS monoclonal IgM antibody (1:500 dilution; F58-10E4; Amsbio) or an isotype control IgM (1:500 dilution; MM-30; BioLegend) for 30 min on ice and then stained with a Cy3-conjugated anti-IgM antibody (1:200 dilution; Jackson ImmunoResearch) for 30 min on ice in the presence of propidium iodide (PI; 1:100 dilution). Stained cells were analyzed on a LSRFortessa X-20 (BD Biosciences).

Western blotting

Cells were lysed in RIPA buffer [25 mM Tris (pH 7.5), 150 mM NaCl, 1% NP-40, 1% sodium deoxycholate, 0.1% sodium dodecyl sulfate (SDS), 1 mM EDTA, and protease inhibitor cocktail (Roche)]. Solubilized proteins were denatured in 2× SDS sample buffer at 90°C for 5 min. Samples were subjected to SDS-polyacrylamide gel electrophoresis (PAGE) and transferred to polyvinylidene fluoride (PVDF) membranes. Proteins were detected using an anti-SCARB2 antibody (AF1966; R&D Systems), an anti-flag antibody (F7425; Sigma-Aldrich), and an anti-β-actin antibody (AC74; Sigma-Aldrich).

Quantitative PCR of the virus genome

Viral RNA was extracted using the QIAamp Viral RNA Mini Kit (Qiagen) followed by reverse transcription using PrimeScript RT Master Mix (Takara). Quantitative PCR was conducted using enterovirus universal primer pairs [56], SYBR Select Master Mix (Applied Biosystems/ThermoFisher Scientific), and a LightCycler 480-II (Roche Life Science).

NGS analysis

RNA was extracted from virus preparations using the QIAamp Viral RNA Extraction Mini Kit (Qiagen). The NEBNext RNA Library Prep Kit for Illumina and the NEBNext Multiplex Oligos for Illumina Dual Index Primer Set 1 (both from New England Biosciences) were used for all library preparation. The average fragment length for the libraries was 300 bp. The quality of the libraries was evaluated on an Agilent 2100 Bioanalyzer (Agilent) using the Agilent DNA 1000 Kit (Agilent) and quantitated using the Kapa Library Quantification Kit (Kapa Biosystems). The libraries were sequenced on a MiSeq instrument (Illumina) using a 150 bp paired-end kit.

Adaptor and low quality sequences were trimmed from raw sequence reads using PRINSEQ software version 0.20.4 (http://prinseq.sourceforge.net/index.html) with the following settings: derep, 14; min_len, 70; min_qual_mean, 25; ns_max_n, 1; lc_threshold, 60; trim_tail_left, 5; trim_tail_right, 5; trim_qual_left, 20; trim_qual_right, 20; trim_qual_window, 5. To remove reads derived from host cells, the remaining reads were mapped to the hg38 human reference genome using Bowtie 2 software version 2.3.4.1 with the default settings [57]. Consensus viral genomic sequences were assembled from the hg38 unmapped reads using IVA software version 1.0.8 with the default settings [58]. Contigs generated by the assembly contained human sequences, indicating the presence of human-derived reads in the viral reads. To diminish the effect of human-derived reads on SNV analysis, non-viral contigs were combined with the reference sequence of the 2716-Ymg-03 strain (GenBank: LC375766.1) or the Isehara strain (LC375764.1). The remaining reads were mapped to the combined multi-FASTA file using Smalt software version 0.7.6 (https://www.sanger.ac.uk/science/tools/smalt-0) and converted to BAM files using Samtools software version 1.9-4-gaelf9d8. Generated BAM files were used for variant calling using Lofreq software version 2.1.2 [59]. SVNs were annotated using SnpEff software version 4.3t [60]. To ensure the accuracy of SNV data, SNVs that met the following criteria were discarded: the abundance ratio was less than 10−3 and the coverage at the position was less than log(1 –p)/log(1 –f)– 1, where the mutation frequency f is detected with a probability p or better [61]. In this study, p was set to 0.95.

Ratios of mutated codons at VP1-145 were calculated based on the abundance of haplotypes estimated using the shorah.py program of ShoRAH software version 1.1.2 [62] with the following settings: analyzed region, 2800–2920; window size, 120; shift, 3. The estimated abundance of the reconstructed haplotypes and their sequences described in the resulting .popl file were used for calculation of the VP1-145 codon abundance.

The synonymous and nonsynonymous nucleotide diversity of each gene was calculated using SNPGENIE software [40].

Statistical analysis

Spearman’s correlation tests (Fig 8) were performed using R version 3.5.1. In multiple comparisons tests, p-values were adjusted by false discovery rate correction (Fig 8).


Zdroje

1. Racaniello V. Picornaviridae: The Viruses and Their Replication. In: Knipe D, Howley PM, editors. Fields Virology. 6th ed ed. Philadelphia: Wolters Kluwer Health/Lippincott Williams & Wilkins; 2013. p. 453–89.

2. Pallansch M, Oberste MS, Whitton JL. Enteroviruses: Polioviruses, Coxackieviruses, Echoviruses, and Newer Enteroviruses. In: Knipe D, Howley PM, editors. Fields Virology. 6th ed ed. Philadelphia: Wolters Kluwer Health/Lippincott Williams & Wilkins; 2013. p. 490–530.

3. Chan LG, Parashar UD, Lye MS, Ong FG, Zaki SR, Alexander JP, et al. Deaths of children during an outbreak of hand, foot, and mouth disease in sarawak, malaysia: clinical and pathological characteristics of the disease. For the Outbreak Study Group. Clin Infect Dis. 2000;31(3):678–83. doi: 10.1086/314032 11017815

4. Khanh TH, Sabanathan S, Thanh TT, Thoa le PK, Thuong TC, Hang V, et al. Enterovirus 71-associated hand, foot, and mouth disease, Southern Vietnam, 2011. Emerging infectious diseases. 2012;18(12):2002–5. doi: 10.3201/eid1812.120929 23194699

5. Ho M, Chen ER, Hsu KH, Twu SJ, Chen KT, Tsai SF, et al. An epidemic of enterovirus 71 infection in Taiwan. Taiwan Enterovirus Epidemic Working Group. The New England journal of medicine. 1999;341(13):929–35. doi: 10.1056/NEJM199909233411301 10498487

6. Zhang Y, Zhu Z, Yang W, Ren J, Tan X, Wang Y, et al. An emerging recombinant human enterovirus 71 responsible for the 2008 outbreak of hand foot and mouth disease in Fuyang city of China. Virol J. 2010;7:94. doi: 10.1186/1743-422X-7-94 20459851

7. Chen KT, Chang HL, Wang ST, Cheng YT, Yang JY. Epidemiologic features of hand-foot-mouth disease and herpangina caused by enterovirus 71 in Taiwan, 1998–2005. Pediatrics. 2007;120(2):e244–52. doi: 10.1542/peds.2006-3331 17671037

8. Chua BH, Phuektes P, Sanders SA, Nicholls PK, McMinn PC. The molecular basis of mouse adaptation by human enterovirus 71. The Journal of general virology. 2008;89(Pt 7):1622–32. doi: 10.1099/vir.0.83676–0

9. Wang YF, Chou CT, Lei HY, Liu CC, Wang SM, Yan JJ, et al. A mouse-adapted enterovirus 71 strain causes neurological disease in mice after oral infection. Journal of virology. 2004;78(15):7916–24. doi: 10.1128/JVI.78.15.7916-7924.2004 15254164

10. Zaini Z, McMinn P. A single mutation in capsid protein VP1 (Q145E) of a genogroup C4 strain of human enterovirus 71 generates a mouse-virulent phenotype. The Journal of general virology. 2012;93(Pt 9):1935–40. doi: 10.1099/vir.0.043893–0

11. Arita M, Ami Y, Wakita T, Shimizu H. Cooperative effect of the attenuation determinants derived from poliovirus sabin 1 strain is essential for attenuation of enterovirus 71 in the NOD/SCID mouse infection model. Journal of virology. 2008;82(4):1787–97. doi: 10.1128/JVI.01798-07 18057246

12. Huang SW, Wang YF, Yu CK, Su IJ, Wang JR. Mutations in VP2 and VP1 capsid proteins increase infectivity and mouse lethality of enterovirus 71 by virus binding and RNA accumulation enhancement. Virology. 2012;422(1):132–43. doi: 10.1016/j.virol.2011.10.015 22078110

13. Kobayashi K, Sudaka Y, Takashino A, Imura A, Fujii K, Koike S. Amino Acid Variation at VP1-145 of Enterovirus 71 Determines Attachment Receptor Usage and Neurovirulence in Human Scavenger Receptor B2 Transgenic Mice. Journal of virology. 2018;92.

14. Fujii K, Sudaka Y, Takashino A, Kobayashi K, Kataoka C, Suzuki T, et al. VP1 Amino Acid Residue 145 of Enterovirus 71 Is a Key Residue for Its Receptor Attachment and Resistance to Neutralizing Antibody during Cynomolgus Monkey Infection. Journal of virology. 2018;92.

15. Kataoka C, Suzuki T, Kotani O, Iwata-Yoshikawa N, Nagata N, Ami Y, et al. The Role of VP1 Amino Acid Residue 145 of Enterovirus 71 in Viral Fitness and Pathogenesis in a Cynomolgus Monkey Model. PLoS Pathog. 2015;11(7):e1005033. doi: 10.1371/journal.ppat.1005033 26181772

16. Mizuta K, Abiko C, Murata T, Matsuzaki Y, Itagaki T, Sanjoh K, et al. Frequent importation of enterovirus 71 from surrounding countries into the local community of Yamagata, Japan, between 1998 and 2003. J Clin Microbiol. 2005;43(12):6171–5. doi: 10.1128/JCM.43.12.6171-6175.2005 16333123

17. Mizuta K, Aoki Y, Suto A, Ootani K, Katsushima N, Itagaki T, et al. Cross-antigenicity among EV71 strains from different genogroups isolated in Yamagata, Japan, between 1990 and 2007. Vaccine. 2009;27(24):3153–8. doi: 10.1016/j.vaccine.2009.03.060 19446185

18. McMinn P, Lindsay K, Perera D, Chan HM, Chan KP, Cardosa MJ. Phylogenetic analysis of enterovirus 71 strains isolated during linked epidemics in Malaysia, Singapore, and Western Australia. Journal of virology. 2001;75(16):7732–8. doi: 10.1128/JVI.75.16.7732-7738.2001 11462047

19. Cordey S, Petty TJ, Schibler M, Martinez Y, Gerlach D, van Belle S, et al. Identification of site-specific adaptations conferring increased neural cell tropism during human enterovirus 71 infection. PLoS Pathog. 2012;8(7):e1002826. doi: 10.1371/journal.ppat.1002826 22910880

20. Arita M, Shimizu H, Nagata N, Ami Y, Suzaki Y, Sata T, et al. Temperature-sensitive mutants of enterovirus 71 show attenuation in cynomolgus monkeys. The Journal of general virology. 2005;86(Pt 5):1391–401. 86/5/1391 [pii] doi: 10.1099/vir.0.80784–0

21. Yeh MT, Wang SW, Yu CK, Lin KH, Lei HY, Su IJ, et al. A Single Nucleotide in Stem Loop II of 5'-Untranslated Region Contributes to Virulence of Enterovirus 71 in Mice. PLoS One. 2011;6(11):e27082. doi: 10.1371/journal.pone.0027082 22069490

22. Li R, Zou Q, Chen L, Zhang H, Wang Y. Molecular analysis of virulent determinants of enterovirus 71. PLoS One. 2011;6(10):e26237. doi: 10.1371/journal.pone.0026237 22039449

23. Chang CK, Wu SR, Chen YC, Lee KJ, Chung NH, Lu YJ, et al. Mutations in VP1 and 5'-UTR affect enterovirus 71 virulence. Scientific reports. 2018;8(1):6688. doi: 10.1038/s41598-018-25091-7 29703921

24. Tan CW, Sam IC, Lee VS, Wong HV, Chan YF. VP1 residues around the five-fold axis of enterovirus A71 mediate heparan sulfate interaction. Virology. 2017;501:79–87. doi: 10.1016/j.virol.2016.11.009 27875780

25. Vignuzzi M, Stone JK, Arnold JJ, Cameron CE, Andino R. Quasispecies diversity determines pathogenesis through cooperative interactions in a viral population. Nature. 2006;439(7074):344–8. doi: 10.1038/nature04388 16327776

26. Huang SW, Huang YH, Tsai HP, Kuo PH, Wang SM, Liu CC, et al. A Selective Bottleneck Shapes the Evolutionary Mutant Spectra of Enterovirus A71 during Viral Dissemination in Humans. Journal of virology. 2017;91(23). doi: 10.1128/jvi.01062-17 28931688

27. Xiao Y, Dolan PT, Goldstein EF, Li M, Farkov M, Brodsky L, et al. Poliovirus intrahost evolution is required to overcome tissue-specific innate immune responses. Nature communications. 2017;8(1):375. doi: 10.1038/s41467-017-00354-5 28851882

28. Yamayoshi S, Yamashita Y, Li J, Hanagata N, Minowa T, Takemura T, et al. Scavenger receptor B2 is a cellular receptor for enterovirus 71. Nature medicine. 2009;15(7):798–801. doi: 10.1038/nm.1992 19543282

29. Yamayoshi S, Iizuka S, Yamashita T, Minagawa H, Mizuta K, Okamoto M, et al. Human SCARB2-dependent infection by coxsackievirus A7, A14, and A16 and enterovirus 71. Journal of virology. 2012;86(10):5686–96. doi: 10.1128/JVI.00020-12 22438546

30. Yamayoshi S, Ohka S, Fujii K, Koike S. Functional comparison of SCARB2 and PSGL1 as receptors for enterovirus 71. Journal of virology. 2013;87(6):3335–47. doi: 10.1128/JVI.02070-12 23302872

31. Zhou D, Zhao Y, Kotecha A, Fry EE, Kelly JT, Wang X, et al. Unexpected mode of engagement between enterovirus 71 and its receptor SCARB2. Nat Microbiol. 2019;4(3):414–9. doi: 10.1038/s41564-018-0319-z 30531980

32. Nishimura Y, Shimojima M, Tano Y, Miyamura T, Wakita T, Shimizu H. Human P-selectin glycoprotein ligand-1 is a functional receptor for enterovirus 71. Nature medicine. 2009;15(7):794–7. nm.1961 [pii] doi: 10.1038/nm.1961 19543284

33. Tan CW, Poh CL, Sam IC, Chan YF. Enterovirus 71 uses cell surface heparan sulfate glycosaminoglycan as an attachment receptor. Journal of virology. 2013;87(1):611–20. doi: 10.1128/JVI.02226-12 23097443

34. Yang SL, Chou YT, Wu CN, Ho MS. Annexin II Binds to Capsid Protein VP1 of Enterovirus 71 and Enhances Viral Infectivity. Journal of virology. 2011;85(22):11809–20. JVI.00297-11 [pii] doi: 10.1128/JVI.00297-11 21900167

35. Su PY, Liu YT, Chang HY, Huang SW, Wang YF, Yu CK, et al. Cell surface sialylation affects binding of enterovirus 71 to rhabdomyosarcoma and neuroblastoma cells. BMC microbiology. 2012;12:162. doi: 10.1186/1471-2180-12-162 22853823

36. Su PY, Wang YF, Huang SW, Lo YC, Wang YH, Wu SR, et al. Cell surface nucleolin facilitates enterovirus 71 binding and infection. Journal of virology. 2015;89(8):4527–38. doi: 10.1128/JVI.03498-14 25673703

37. Du N, Cong H, Tian H, Zhang H, Zhang W, Song L, et al. Cell surface vimentin is an attachment receptor for enterovirus 71. Journal of virology. 2014;88(10):5816–33. doi: 10.1128/JVI.03826-13 24623428

38. He QQ, Ren S, Xia ZC, Cheng ZK, Peng NF, Zhu Y. Fibronectin Facilitates Enterovirus 71 Infection by Mediating Viral Entry. J Virol. 2018;92(9). doi: 10.1128/JVI.02251-17 29467312

39. Nishimura Y, Lee H, Hafenstein S, Kataoka C, Wakita T, Bergelson JM, et al. Enterovirus 71 binding to PSGL-1 on leukocytes: VP1-145 acts as a molecular switch to control receptor interaction. PLoS Pathog. 2013;9(7):e1003511. doi: 10.1371/journal.ppat.1003511 23935488

40. Nelson CW, Moncla LH, Hughes AL. SNPGenie: estimating evolutionary parameters to detect natural selection using pooled next-generation sequencing data. Bioinformatics (Oxford, England). 2015;31(22):3709–11. doi: 10.1093/bioinformatics/btv449 26227143

41. Lee E, Hall RA, Lobigs M. Common E protein determinants for attenuation of glycosaminoglycan-binding variants of Japanese encephalitis and West Nile viruses. Journal of virology. 2004;78(15):8271–80. doi: 10.1128/JVI.78.15.8271-8280.2004 15254199

42. Lee E, Lobigs M. Mechanism of virulence attenuation of glycosaminoglycan-binding variants of Japanese encephalitis virus and Murray Valley encephalitis virus. Journal of virology. 2002;76(10):4901–11. doi: 10.1128/JVI.76.10.4901-4911.2002 11967307

43. Lee E, Wright PJ, Davidson A, Lobigs M. Virulence attenuation of Dengue virus due to augmented glycosaminoglycan-binding affinity and restriction in extraneural dissemination. The Journal of general virology. 2006;87(Pt 10):2791–801. doi: 10.1099/vir.0.82164–0

44. Anez G, Men R, Eckels KH, Lai CJ. Passage of dengue virus type 4 vaccine candidates in fetal rhesus lung cells selects heparin-sensitive variants that result in loss of infectivity and immunogenicity in rhesus macaques. Journal of virology. 2009;83(20):10384–94. doi: 10.1128/JVI.01083-09 19656873

45. Klimstra WB, Ryman KD, Johnston RE. Adaptation of Sindbis virus to BHK cells selects for use of heparan sulfate as an attachment receptor. Journal of virology. 1998;72(9):7357–66. 9696832

46. Bernard KA, Klimstra WB, Johnston RE. Mutations in the E2 glycoprotein of Venezuelan equine encephalitis virus confer heparan sulfate interaction, low morbidity, and rapid clearance from blood of mice. Virology. 2000;276(1):93–103. doi: 10.1006/viro.2000.0546 11021998

47. Mandl CW, Kroschewski H, Allison SL, Kofler R, Holzmann H, Meixner T, et al. Adaptation of tick-borne encephalitis virus to BHK-21 cells results in the formation of multiple heparan sulfate binding sites in the envelope protein and attenuation in vivo. Journal of virology. 2001;75(12):5627–37. doi: 10.1128/JVI.75.12.5627-5637.2001 11356970

48. Gardner CL, Hritz J, Sun C, Vanlandingham DL, Song TY, Ghedin E, et al. Deliberate attenuation of chikungunya virus by adaptation to heparan sulfate-dependent infectivity: a model for rational arboviral vaccine design. PLoS neglected tropical diseases. 2014;8(2):e2719. doi: 10.1371/journal.pntd.0002719 24587470

49. Bochkov YA, Watters K, Basnet S, Sijapati S, Hill M, Palmenberg AC, et al. Mutations in VP1 and 3A proteins improve binding and replication of rhinovirus C15 in HeLa-E8 cells. Virology. 2016;499:350–60. doi: 10.1016/j.virol.2016.09.025 27743961

50. Vlasak M, Goesler I, Blaas D. Human rhinovirus type 89 variants use heparan sulfate proteoglycan for cell attachment. Journal of virology. 2005;79(10):5963–70. doi: 10.1128/JVI.79.10.5963-5970.2005 15857982

51. Sa-Carvalho D, Rieder E, Baxt B, Rodarte R, Tanuri A, Mason PW. Tissue culture adaptation of foot-and-mouth disease virus selects viruses that bind to heparin and are attenuated in cattle. Journal of virology. 1997;71(7):5115–23. 9188578

52. Minor PD. Live attenuated vaccines: Historical successes and current challenges. Virology. 2015;479–480:379–92. doi: 10.1016/j.virol.2015.03.032 25864107

53. Tsai YH, Huang SW, Hsieh WS, Cheng CK, Chang CF, Wang YF, et al. Enterovirus A71 Containing Codon-Deoptimized VP1 and High-Fidelity Polymerase as Next-Generation Vaccine Candidate. J Virol. 2019;93(13). doi: 10.1128/JVI.02308-18 30996087

54. Kohara M, Omata T, Kameda A, Semler BL, Itoh H, Wimmer E, et al. In vitro phenotypic markers of a poliovirus recombinant constructed from infectious cDNA clones of the neurovirulent Mahoney strain and the attenuated Sabin 1 strain. Journal of virology. 1985;53(3):786–92. 2983090

55. Kärber G. Beitrag zur kollektiven Behandlung pharmakologischer Reihenversuche. Naunyn-Schmiedebergs Archiv für Experimentelle Pathologie und Pharmakologie. 1931;162(4):480–3. doi: 10.1007/bf01863914

56. Jonsson N, Gullberg M, Lindberg AM. Real-time polymerase chain reaction as a rapid and efficient alternative to estimation of picornavirus titers by tissue culture infectious dose 50% or plaque forming units. Microbiol Immunol. 2009;53(3):149–54. doi: 10.1111/j.1348-0421.2009.00107.x 19302525

57. Langmead B, Salzberg SL. Fast gapped-read alignment with Bowtie 2. Nat Methods. 2012;9(4):357–9. doi: 10.1038/nmeth.1923 22388286

58. Hunt M, Gall A, Ong SH, Brener J, Ferns B, Goulder P, et al. IVA: accurate de novo assembly of RNA virus genomes. Bioinformatics. 2015;31(14):2374–6. doi: 10.1093/bioinformatics/btv120 25725497

59. Wilm A, Aw PP, Bertrand D, Yeo GH, Ong SH, Wong CH, et al. LoFreq: a sequence-quality aware, ultra-sensitive variant caller for uncovering cell-population heterogeneity from high-throughput sequencing datasets. Nucleic acids research. 2012;40(22):11189–201. doi: 10.1093/nar/gks918 23066108

60. Cingolani P, Platts A, Wang le L, Coon M, Nguyen T, Wang L, et al. A program for annotating and predicting the effects of single nucleotide polymorphisms, SnpEff: SNPs in the genome of Drosophila melanogaster strain w1118; iso-2; iso-3. Fly. 2012;6(2):80–92. doi: 10.4161/fly.19695 22728672

61. Lorenzo-Redondo R, Fryer HR, Bedford T, Kim EY, Archer J, Pond SLK, et al. Persistent HIV-1 replication maintains the tissue reservoir during therapy. Nature. 2016;530(7588):51–6. doi: 10.1038/nature16933 26814962

62. Zagordi O, Bhattacharya A, Eriksson N, Beerenwinkel N. ShoRAH: estimating the genetic diversity of a mixed sample from next-generation sequencing data. BMC Bioinformatics. 2011;12(1):119. doi: 10.1186/1471-2105-12-119 21521499


Článek vyšel v časopise

PLOS Pathogens


2020 Číslo 3
Nejčtenější tento týden
Nejčtenější v tomto čísle
Kurzy Podcasty Doporučená témata Časopisy
Přihlášení
Zapomenuté heslo

Zadejte e-mailovou adresu, se kterou jste vytvářel(a) účet, budou Vám na ni zaslány informace k nastavení nového hesla.

Přihlášení

Nemáte účet?  Registrujte se

#ADS_BOTTOM_SCRIPTS#