INTRODUCTION

Top Abstract Introduction Materials and Methods Results Discussion References

Analysis of polioviruses isolated either during a given outbreak of poliomyelitis (30, 35, 51) or from sequential fecal samples from infected individuals (29, 31, 34) has revealed rapid changes in the nucleotide sequence of the viral 7.5-kb RNA. Oligonucleotide fingerprinting (51) and genome sequencing (30, 34) suggested that, during the epidemic, the nucleotide substitutions ranged from 1 to 2% per year. The molecular basis of such genetic instability, common to all RNA viruses, resides in a high error rate of the viral RNA-dependent RNA polymerases, estimated to be, for poliovirus, on the order of 10⁴ to 10⁵ substitutions per base per replication (12, 70, 71), and the absence of proofreading mechanisms for the correction of the errors made. The combination of these properties results in a high heterogeneity (the "quasispecies" nature) of all populations of the virus (8, 10, 11, 63).

Less understood, however, are the rules governing the fixation of specific mutations upon passages of a viral population. Generally speaking, a mutation may decrease or increase the level of virus fitness for a particular ecological niche or leave it unchanged. Nucleotide substitutions associated with changes in the "sense" of codons (nonsynonymous mutations) are less likely to be neutral than substitutions resulting in no change of the codon meaning (synonymous mutations). Fixation of mutations conferring a selective advantage is readily understood in the framework of positive Darwinian selection. Adverse mutations may be eliminated by negative selection.

The likelihood of fixation of a mutation depends not only on the associated changes in virus fitness but, to a large extent, on the sizes of propagated populations (48, 52-54). Consecutive passages involving small viral populations, the situation very typical of the natural poliovirus infection, may well result in the accumulation of neutral or even debilitating mutations (5, 7, 13, 18), the phenomenon known as the Muller's ratchet (5, 6, 8). Thus, the consecutive changes in viral lineages are not necessarily adaptive but may be due to random sampling. There are reasons to believe that the majority of generated mutations may in fact be neutral or nearly so (33). Since the same appears to be true of the fixed mutations (60), the rate of nucleotide changes should be approximately constant in time (33). This concept of a "molecular clock" was suggested to be valid for the evolution of RNA viruses (17, 24, 69) and is not infrequently used for estimations of the viral "age" (31, 65, 73).

Here, we undertook an analysis of genomic changes in two evolving poliovirus populations. One was represented by a lineage of wild-type (wt) viruses of T geotype (40) in circulation in the territory of the former Soviet Union (FSU) during 1991 to 1995. The other consisted of the strains successively isolated from an immunodeficient patient (IDP) with paralytic poliomyelitis (31). The obtained sequences revealed basic similarities in the character of accumulation of mutations in the two lineages, even though some quantitative differences were also evident. Different regions of the viral genome have evolved at nonidentical rates, and these rates might have changed during the history of a population. Nevertheless, the major contribution to the evolving populations came from the linear accumulation of predominantly synonymous substitutions. From the perspective of current theories of viral evolution (8, 24), the results obtained suggest a unifying model for the evolution of circulating poliovirus on the one hand and for virus growing in the gut of an individual on the other. This model has important implications for the epidemiology of poliomyelitis.

	MATERIALS AND METHODS

Top Abstract Introduction Materials and Methods Results Discussion References

Virus isolation and typing. Two sets of polioviruses of serotype 1 were analyzed (Table 1; Fig. 1). One contained 24 wt strains isolated in the territory of the FSU, China, and Pakistan from 1991 to 1995. Of these, 19 strains belonged to the T geotype (40). The other set included five vaccine-related strains (GenBank accession numbers for VP1 sequences, AF083933 to AF083937) obtained from an IDP with paralytic poliomyelitis from Missouri in 1981 (31).

View larger version (28K): [in this window] [in a new window]   FIG. 1.   The VP1/2ABC' maximum-likelihood bootstrap consensus tree (PUZZLE, version 4.0) rooted on PV1/169AZB59 showing the relationships of poliovirus type 1 wild and vaccine strains. The scale corresponds to 10% nucleotide divergence. The isolates belonging to the T geotype and IDP lineages are boxed.

Reverse transcription, PCR, and sequencing of poliovirus RNAs. RNA was extracted from cell lysates with Ultraspec 3 (Biotexc) and was reverse transcribed using random hexamer primers (Boehringer Ingelheim) with avian myeloblastosis virus reverse transcriptase (Promega). DNA copies of selected genomic portions were amplified by PCR with appropriate primers (the primer sequences are available upon request), and the PCR products were purified with the QIAquick system (Qiagen) and were directly sequenced either automatically (Applied Biosystems) or manually by using the Sequenase kit, version 2.0 (U.S. Biochemicals). The following regions were analyzed in all the 29 strains: the complete VP1 (906 nucleotides [nt]; nt 2500 to 3408), 2A (447 nt; nt 3409 to 3855), and 2B (291 nt; nt 3856 to 4146) genes as well as a portion of the 2C gene (51 nt; nt 4147 to 4197; called 2C'). A portion of the 3D gene (792 nt; nt 6601 to 7392; called 3D') in PV1/169AZB59, -/7TAJ91, -/5794UZB94, -/6484CHE95, -/6486CHE95, -/6490CHE95, -/6405ING95, -/422RUS91, and -/919GEO90 was analyzed. The coordinates are given according to the consensus numbering system of Toyoda et al. (68).

Comparative analysis of nucleotide sequences. Multiple alignment of the determined sequences was performed by the program CLUSTAL W, version 1.74 (67). Phylogenetic analysis was carried out by the PHYLIP, version 3.5c, program package (20), PUZZLE, version 4.0 (64), and METREE (58). Analysis under maximum-likelihood optimality criteria was done with empirical base frequencies and the HKY85 model of substitution (25). The following algorithms have also been used for testing the consistency of branching: maximum parsimony (16), neighbor joining (59), minimum evolution (57), and the unweighted average-linkage method of clustering. Evolutionary distances were estimated using the two-parameter method of Kimura (32) and the Tamura-Nei model (66). The reliability of clustering was tested by 10,000 iterations in the quartet-puzzling method (PUZZLE, version 4.0) and by 1,000 bootstrapping replicates in other methods mentioned above.

Statistical analysis of rare-codon conservation. A codon was defined as rare if its abundance was no more than 30% of that of its most abundant synonym in the poliovirus codon usage table. Five families of codons (encoding Ala, Arg, Pro, Ser, and Thr) contain codons meeting this criterion. The same codons were also scarce in humans. To eliminate the bias due to the close relatedness of some of the strains in our data set, the 6 most-diverged VP1 and 2AB sequences were selected out of 29 sequences available. Upon pairwise comparisons, these sequences exhibited divergences higher than 14 and 18%, respectively. The following strains were used in the analysis of both regions: PV1/169AZB59, -/7TAJ91, -/832GEO85, -/5937RUS94, and -/Day 200; in addition, the strains PV1/6070CHN94 and -/6013TAJ94 were used for the VP1 and 2AB sets, respectively. It was assumed that all of the sites mutated independently and that the mutation rate was uniform across the region in question. The ratios of the respective rare codons to all their synonyms in the whole poliovirus (Sabin 1) RNA were calculated. The Bernoulli formula was used to determine the probability of the occurrence of a certain number of rare codons (1 to 6) in the positions of the six aligned sequences that corresponded to each of the above amino acid families. The positions having significant deviations (P < 0.05) from random codon usage were considered anomalously conserved.

Nucleotide sequence accession numbers. The nucleotide sequence data reported in this paper are available from the GenBank nucleotide sequence database under accession no. AF233098 to AF233222.

	RESULTS

Top Abstract Introduction Materials and Methods Results Discussion References

General characterization of the two poliovirus populations. Two evolving poliovirus populations of serotype 1 were studied (Table 1). One corresponded to 11 wt strains isolated from poliomyelitis patients in the territory of the FSU in 1991 to 1995 and constituted a distinct cluster of the T geotype (Fig. 1). The other set of viruses was represented by five vaccine-related strains successively isolated from an IDP with paralytic poliomyelitis (31). The primary structures of the 150-nt VP1/2A region of the RNA of most of these T geotype isolates and of the VP1 genes of the vaccine-derived viruses have been determined previously (31, 40). Here, the sequenced portions were extended to include the entire 1,694-nt VP1/2ABC' region (nt 2500 to 4197) of all of the viruses and a 792-nt portion of the 3D gene (nt 6601 to 7392) of all IDP isolates and five T geotype strains.

View larger version (16K): [in this window] [in a new window]   FIG. 2.   The crossover region within the 2B gene of the strain Day23M. Dashes, identical bases in all three sequences; framed segment (nt 3973 to 3996), hypothetical recombination region; arrow, beginning of the 2B2 region. Identical nucleotides are boxed.

Time course of genomic evolution of the poliovirus populations. The total number of accumulated mutations in the sequenced portion of the genomes was plotted against time. The time course of mutation accumulation could be well approximated by straight lines (the regression coefficient, r, was 0.95) corresponding to the rates of 1.19 × 10² and 2.62 × 10² substitutions/site/year for the T geotype and IDP lineages, respectively (Fig. 3A and D).

View larger version (20K): [in this window] [in a new window]   FIG. 3.   Mutation fixation rates in the sequenced regions of 11 wt (VP1/2ABC') and 5 IDP (VP1/2ABC'/3D') isolates in substitutions per site per year (r, correlation coefficient) for all (Kt), synonymous (Ks), and nonsynonymous (Ka) sites. y axis, percentage of mutated sites. The position of the first strain in each of the lineages corresponds to x = 0 and y = 0. Dashed lines, 95% confidence intervals.

2

2

2

2

Nonuniform rates of evolution of different parts of the viral genome. The character of mutations accumulated was analyzed for the three genomic regions (VP1, 2AB, and 3D') separately (Fig. 4 and 5). The total and synonymous mutations accumulated linearly with time in all the studied regions of both lineages (with r varying between 0.89 and 0.97). However, in both sets of viruses, different RNA regions varied significantly in their relative mutabilities. The 2AB region evolved 2.5 to 5.5 times faster than the VP1 gene, and again a major contribution to this difference came from synonymous substitutions, with the relevant rate for the 2AB region in the IDP strains being exceptionally high, 18.25 × 10² (Fig. 5E). No significant difference between the rates of evolution of the 2A and 2B genes studied separately was noted (not shown).

View larger version (16K): [in this window] [in a new window]   FIG. 4.   Mutation fixation rates in the T geotype (11 isolates) lineage in all (Kt), synonymous (Ks), and nonsynonymous (Ka) sites. Solid and boldface dashed lines, data for VP1 and 2AB regions, respectively. For other details, see the legend to Fig. 3.

View larger version (32K): [in this window] [in a new window]   FIG. 5.   Mutation fixation rates in the IDP lineage (five isolates) in all (Kt), synonymous (Ks), and nonsynonymous (Ka) sites in the VP1, 2AB, and 3D' regions. For other details, see the legend to Fig. 3.

Character of mutation fixation. In the vaccine-derived viruses, a significant proportion of mutations, once fixed, continued to be present in the genomes of all of the subsequent isolates. This was true of the mutations accumulated both before and after the isolation of the first strain (Fig. 6). Such a character of mutation accumulation was consistent with the assumption that the successive isolates of the IDP lineage were related to each other in the parent-progeny fashion. Such a direct relationship was less evident for the T geotype populations, where a significant proportion of mutations present in a given isolate was absent from some or all subsequent isolates (Fig. 6).

View larger version (42K): [in this window] [in a new window]   FIG. 6.   Characterization of mutation fixation in the wt and IDP lineages. y axis, numbers of mutations that once appeared, became fixed in all (open bars) or at least some (hatched bars) of the subsequent isolates, as well as numbers of solitary mutations (black bars). Insets in the IDP panels, comparison of the first isolate (day 23) to the vaccine sequences. For the wt lineage, the numbers on the x axis correspond to the following isolates: 1, PV1/4TAJ91 (compared to PV1/7TAJ91); 2, PV1/5341UKR92; 3, PV1/5794UZB94; 4, PV1/6011TAJ94; 5, PV1/6064CHE95; 6, PV1/6484CHE95; 7, PV1/6486CHE95; 8, PV1/6490CHE95; and 9, PV1/6405ING95; for the IDP lineage, the numbers correspond to the following isolates: 1, day 48 (compared to day 23); 2, day 126; 3, day 158. The last isolates in the lineages are omitted.

Codon usage and the nonuniform evolution of the viral genome. A possible explanation of different rates of fixation of synonymous mutations along the viral genome may be related to the codon usage bias. Different synonymous codons are known to be used with different frequencies in different organisms, different tissues, and even different genes, contributing thereby to the translational control of gene expression (61, 74). It can be argued that highly preferred codons are favorable for efficient translation and may therefore be less prone to be replaced during evolution. Consequently, the genomic regions with an abundance of such codons might be expected to undergo nucleotide replacements relatively less frequently. No significant heterogeneity in codon usage among the genomic regions studied could be found (not shown). From this analysis, however, a statistically significant conservation of rarely used codons at certain positions of the VP1, but not the 2AB, region was noticed (Fig. 7A). In this analysis, not only the true conservation of a rare codon but also the cases in which a rare codon was replaced by another synonymous rare codon were taken into account. At certain positions, the consecutive conserved rare codons were clustered (Fig. 7A). In the vicinity of 9 out of the 15 conserved rare codons (at a distance of 1 to 3 codons) additional rare codons were also present in at least 75% of the 29 aligned VP1 sequences (not shown).

View larger version (23K): [in this window] [in a new window]   FIG. 7.   Conserved rare codons of the VP1 gene. (A) Sites with anomalous conservation of the rare codons in the six most divergent sequences of the type 1 poliovirus. The codons are numbered from the start of the VP1 gene. The rare codons are underlined. The number of each of the individual codons in a given position of the six sequences is shown in parentheses. (B) VP1 secondary structure wire plot (Protein Data Bank; www.biochem.ucl.ac.uk/bsm/pdbsum/2plv/main.html). The amino acids corresponding to the conserved rare codons in panel A are underlined and in boldface. Straight arrows, beta structural elements; curved arrows, alpha helices. The amino acids are numbered from the start of the VP1 gene.

Does the RNA secondary structure contribute to the nonuniform evolution of the viral genome? It could be hypothesized that local minima in the accumulation of synonymous substitutions were due to constraints imposed by the necessity to maintain specific high-order RNA structures. An analysis of this problem was hampered by the dynamic rather than static nature of the RNA structural organization. RNA regions with relatively stable secondary structure elements are expected to possess a low value of the so-called pairing number (P-num), a quantitative measure of the propensity of a certain base to pair with alternative partners in a collection of suboptimal folds (28, 56). We calculated that the proportions of bases with low P-num values (<100) for the VP1 and 2AB regions of the Mahoney genome were quite similar, being 91 out of 906 nt (10%) and 62 out of 738 nt (8.4%), respectively. This fact failed to support the notion that the two RNA regions differed from one another with regard to possessing strictly conserved secondary structures. No such support was obtained from attempts to derive a consensus secondary structure by comparative sequence analysis of these regions either. However, a higher percentage (17.7% or 140 out of 792 nt) of nucleotides with low P-num values was found in the 3D' region. This observation was indicative of substantial structural constraints, which may be partially responsible for lowering the mutation fixation rate in the sequenced portion of the polymerase gene. The close relatedness of the 13 sequences available did not enable us to search for a consensus secondary structure of the 3D' region.

Deviation from the molecular clock mode of evolution. As shown above, the accumulation of total, synonymous, and nonsynonymous mutations appeared to linearly depend on time in all the RNA regions studied with only one exception (nonlinearity of nonsynonymous substitutions in the T geotype VP1 gene). However, a more serious deviation from linearity was revealed by an analysis of mutation rates throughout the whole history of the Sabin vaccine-derived population. Indeed, the rate of synonymous substitutions in 2AB, during the time interval studied here, was approximately fivefold higher than that in VP1. On the other hand, the proportions of synonymous substitutions (K_s values) in the VP1 and 2AB₁ regions of the earliest IDP isolate, compared to the corresponding values for the Sabin progenitor virus, were similar, 30.4 and 27.9, respectively, suggesting that these regions evolved at comparable rates. This discrepancy became clearly evident when the apparent times of divergence of the relevant IDP lineage genes from the Sabin vaccine were calculated by extrapolations using the estimated evolution rates taken from Fig. 5B and E. For the synonymous nucleotide replacements in the VP1 region (Fig. 8A), this time (9.6 years) corresponded reasonably well to the time (6.9 years) elapsed between the last oral polio vaccine vaccination and the onset of paralysis (31). However, the 2AB₁ region (Fig. 8B) appeared to require only one-fourth of the time (~ 2.5 years) to achieve the observed divergence from the progenitor. Possible reasons for this paradox are discussed below.

View larger version (19K): [in this window] [in a new window]   FIG. 8.   Extrapolation of evolution rates of synonymous sites calculated for the IDP lineage for VP1 (A), 2AB1 (B), and combined VP1/2AB1/3D' (C) regions to zero divergence with Sabin 1. The solid upper portions of the lines correspond to the data taken from Fig. 5B and E and 3E. Open arrowhead, extrapolated age of initial infection as calculated by Kew et al. (31); solid arrowhead, time of the last vaccination with oral polio vaccine; arrow, time of the paralysis onset.

	DISCUSSION

Top Abstract Introduction Materials and Methods Results Discussion References

Basic mechanisms of poliovirus evolution. The rate and character of fixation of mutations generated by an error-prone viral RNA polymerase are controlled by positive selection of more-fit variants, negative selection of variants with decreased fitness, and accumulation of neutral or nearly neutral (33) replacements through random sampling of viral exemplars from the "swarm" of genomic variants (7, 9, 11, 13, 14, 18). Admittedly, the discussion of the nature of evolutionary changes would have been more factual if accompanying alterations in the level of viral fitness were known (8). Unfortunately, no adequate animal model to assess the relative fitness of polioviruses in their natural niche, the human gut, is available. It seems likely, however, that all three of the above-mentioned mechanisms contributed, to a certain extent, to the changes registered in this study. Although the exact contribution of each of them cannot as yet be defined, the preponderance of synonymous replacements and the constant rates of mutation fixation suggest a major role for random sampling as the mechanism of neutral evolution. This appeared to be true of both poliovirus populations studied here, even though they differed from one another in several respects: (i) the starting populations were represented by wt (presumably relatively well adapted) and attenuated (less-fit) viruses; (ii) the wt lineage was interrupted by human-to-human transfers, whereas the other one corresponded to a continuously growing population in a single host; and (iii) accordingly, the wt isolates were hardly related to each other in the direct ancestor-progeny fashion, whereas this was almost certainly true of the IDP isolates.

2

3

Nonuniform rates of evolution of different genomic regions and different genomes. A nonuniform rate of mutation accumulation along the viral RNA is not unique for poliovirus. The bias in nonsynonymous mutation fixation in different genomic regions of different RNA viruses has been observed (24, 27) and can readily be explained by differences in constraints imposed on the structures of different viral proteins. More interesting, in the context of the present study, were unequal rates of accumulation of synonymous substitutions.

View larger version (23K): [in this window] [in a new window]   FIG. 9.   Different topologies of the maximum-likelihood bootstrap consensus trees (PUZZLE, version 4.0) of the poliovirus type 1 wild and vaccine (Sabin 1) strains (rooted on PV1/169AZB59) when based on the VP1 (A) and 2AB (B) sequences. For the sake of brevity, a selected set of strains was used; a comparison of the complete trees also demonstrated the same difference in topology (not shown). Strains of the T geotype lineage are boxed. PV1/422RUS91 belongs to geotype A; PV1/919GEO90 and PV1/717RUS94 correspond to geotype G (40).

Disturbances in the operation of molecular clocks. Although the accumulation of mutations is generally linear, deviation from the molecular clock is possible. Mutations within the 2AB region of the IDP lineage appeared to accumulate nonlinearly with time, the rate of evolution being apparently higher during the last period of the population history.

A hypothetical scenario of viral evolution in the immunodeficient host. The events that occurred in the IDP host between the vaccination and the isolation of the first virus are not documented. Indirect evidence suggests that the lineage studied might have become predominant not long before the first isolate (PV1/2677USA81) was obtained. First, the relevant fecal sample contained two variants, one of which, the minor variant (m), subsequently disappeared (31), suggesting that the major variant (M) exhibited better fitness. This selective advantage of M should have likely precluded the long-term coexistence of the two variants in a single host. Second, the paralytic disorder that developed in the IDP after several years of asymptomatic virus carriage was most likely due to a newly arisen more pathogenic variant.

Implications for epidemiology of poliomyelitis. We propose a unifying model of evolution of both wt poliovirus upon its natural circulation and of vaccine-derived virus replicating in the human gut. According to this model, evolution proceeds essentially through consecutive sampling of small portions of a viral population, ensuring predominant fixation of neutral and possibly adverse mutations. As a result the virus fitness is expected to be maintained at a more or less constant level or to decrease, except for the early stage after the immunization, when loss of some attenuating mutations (47) may result in an enhanced fitness. The possible accumulation of adverse mutations may serve to self-limit viral reproduction in the human gut and during person-to-person circulation and may contribute to the relatively low pathogenicity of poliovirus, which is known to cause overt disease only in a tiny minority of infected nonimmune subjects (44). The fitness decrease may be particularly conspicuous for attenuated viruses, leading to the very low level of their transmissibility. Some of the deleterious mutations may be eliminated by negative selection and by intratypic or intertypic recombination (4, 21, 22, 39, 46). Concurrently, the positive selection of rarely arisen more fit and more neurovirulent variants may facilitate the development of vaccine-associated poliomyelitis in the case of vaccine-derived variants and outbreaks of poliomyelitis in the case of circulating wt viruses.