Dengue Virus Type 3 Adaptive Changes during Epidemics in São Jose de Rio Preto, Brazil, 2006–2007

Global dengue virus spread in tropical and sub-tropical regions has become a major international public health concern. It is evident that DENV genetic diversity plays a significant role in the immunopathology of the disease and that the identification of polymorphisms associated with adaptive responses is important for vaccine development. The investigation of naturally occurring genomic variants may play an important role in the comprehension of different adaptive strategies used by these mutants to evade the human immune system. In order to elucidate this role we sequenced the complete polyprotein-coding region of thirty-three DENV-3 isolates to characterize variants circulating under high endemicity in the city of São José de Rio Preto, Brazil, during the onset of the 2006-07 epidemic. By inferring the evolutionary history on a local-scale and estimating rates of synonymous (dS) and nonsynonimous (dN) substitutions, we have documented at least two different introductions of DENV-3 into the city and detected 10 polymorphic codon sites under significant positive selection (dN/dS > 1) and 8 under significant purifying selection (dN/dS < 1). We found several polymorphic amino acid coding sites in the envelope (15), NS1 (17), NS2A (11), and NS5 (24) genes, which suggests that these genes may be experiencing relatively recent adaptive changes. Furthermore, some polymorphisms correlated with changes in the immunogenicity of several epitopes. Our study highlights the existence of significant and informative DENV variability at the spatio-temporal scale of an urban outbreak.


Introduction
Dengue viruses (genus flavivirus, family Flaviviridae) are responsible for dengue fever and severe dengue (previously known as dengue hemorrhagic fever, DHF and dengue shock syndrome, DSS), they are some of the most important arthropod-borne viral diseases worldwide [1][2][3][4]. Increasing incidences of transmission in urban and semi-urban areas in the tropical and sub-tropical countries have promoted dengue to a major public health issue with a hefty economic burden [5][6][7][8][9][10].
Dengue epidemics are largely attributed to three factors: (i) Increased urbanization, (ii) global dissemination of the major mosquito vectors, Aedes aegypti and Aedes albopictus, and (iii) the interaction and evolution of the four DENV serotypes [11]. Inadequate sustained vector control may partially explain unsuccessful disease prevention in developing countries however variables related to the viruses, the human host and the environment should also be considered [12]. Vaccine development is currently underway as a strategy for controlling this worldwide health threat [13][14][15]. Due to the complexity of DENV dynamics, prevailing serotypes and differences in distribution of viruses within different regions of the globe, serotype and lineage replacement events may play an important role in the evolution of DENV viruses and, consequently must be addressed during vaccine development [16].
Gene genealogy-based DENV evolutionary studies helped characterize the genetic diversity and distribution of different serotypes/genotypes in space and time [17][18][19][20][21][22]. Also, they are important for exploring the role of selection on particular dengue proteins [23][24][25][26][27][28][29]. Given that serotype-specific neutralizing antibodies confer limited, if any, protection from subsequent infection by the other three serotypes, the appropriate choice of nucleotide residues in the variable genomic regions is critical for the design of effective tetravalent vaccines [30,31].
Brazil is the South American country with the highest economic impact of dengue and also accounts for the majority of reported cases in the continent as evidenced by a 27 year-long study [8,32].
The highest dengue incidence was among young adults from 2000 to 2007. However, in 2006 the highest incidence rate of severe dengue increased dramatically among children and remained high during severe epidemics in 2007 and 2008 in the state of Rio de Janeiro [32][33][34]. These circumstances highlight the urgent need for continued studies on dengue molecular epidemiology in order to help determine appropriate policies and effective public health campaigns.
The municipality of São José de Rio Preto (SJRP), one of the 645 municipalities within the State of São Paulo according to the Brazilian Institute of Geography and Statistics (Instituto Brasileiro de Geografia e Estatística, IBGE), has a tropical savanna climate (Köppen climate classification), with an annual rainfall of 1410 mm and mean annual temperature of 23.6uC, making it a favorable breeding ground for Aedes mosquitoe vectors. Our previous study on dengue in SJRP [35] using a 399 bp-long portion of the NS5 gene, suggested that DENV-3 was introduced into the city in 2005. Since SJRP is one of the municipalities in São Paulo with the highest number of dengue cases reported yearly according to the Epidemiological Surveillance Center of the state (Centro de Vigilância epidemiológica Alexandre Vranjac, CVE). We chose to further characterize the isolated DENV-3 strains from SJRP from 2006 to 2007 to gain insights into the viral epidemiology that may help control disease.

Area of study
The city of SJRP is located in the northwestern region of the state of São Paulo, Brazil (20u499110 S and 49u229460 W), with a total area of 431 km 2 and an estimated population of 415,769 inhabitants in 2012 (Data obtained from the IBGE). The infestation by Aedes mosquitoes was reported for the first time in 1985 and the first autochthonous cases were reported in 1990, with the introduction of DENV-1 in 1996 (Chiaravalloti-Neto,

Sequencing genomes
A total of 33 plasma samples were collected from dengue patients during the 2006 and 2007 outbreaks as described elsewhere [35]. Viral RNA was extracted directly from plasma using the QIAamp Viral RNA mini kit (Qiagen). cDNA was synthetized in a 20 ml reverse transcription reaction containing: 1 ml of Superscript III Reverse Transcriptase (Invitrogen), 1 ml of random hexamers (50 ng/ml stock), 1 mL of specific reverse primer 59AGAACCTGTTGATTCAACAGCAC39 (10 mM stock) and 5 ml of template RNA. Viral cDNA were diluted to 800 ml in DEPC-treated water and 20 ml were used as template for a 96 specific PCR reactions. Each 10 ml PCR reaction contained: 3 ml of template, 0.03 ml of pfuUltra II polymerase 1 (5 U/ml) (Stratagene), 100 mM dNTPs (Applied Biosystems) and 4 ml of a mixture of forward and reverse primers (0.5 mM stock) (Table S1). Primers were synthesized with M13 sequence tags (forward primers with 59GTAAAACGACGGCCAGT39 and reverse primers with 59CAGGAAACAGCTATGACC39) so that PCR amplicons could be sequenced using universal M13 forward and reverse primers. Each PCR reaction produced 96 overlapping amplicons, each 500-900 nucleotides in length, which were subsequently sequenced bidirectionally using the Big Dye chemistry on ABI3730xl DNA sequencers (Applied Biosystems). The 33 obtained genomic sequences were submitted to the Broad Institute genome project and deposited in GenBank (Table S2 shows all GenBank accession numbers and name of all samples). (http:// www.broadinstitute.org/annotation/viral/Dengue/).

Phylogenetic Inference
The data set includes the obtained 33 genomic sequences from SJRP that were aligned together with 33 DENV-3 worldwide genomic sequences deposited in GenBank (Table 1) using the Muscle 3.8.31 program [37,38]. Visual inspection of the alignment and manual editing were done with the Se-Al v2.0 program (http://tree.bio.ed.ac.uk/software/seal/). A Bayesian maximum clade credibility tree was inferred from a set of plausible trees sampled at the stationary phase of four independent Markov Chain Monte Carlo (MCMC) runs with 40 million generations each using a general time reversible (GTR) model of nucleotide substitution [39], a gamma distributed rate variation and a proportion of invariable sites (GTR + C + I). A relaxed (uncorrelated lognormal) molecular clock [40] and a Bayesian Skyline coalescent tree were used as priors [41]. Analyses were done with Beast v1.6.1 [42]. Convergence of parameters was assessed using Tracer v1.5 program (http://tree.bio.ed.ac.uk/ software/tracer/) until all parameters estimates showed Effective Sample Size (ESS) values over 100.

Selection detection analysis
The SJRP dataset was used in a selection analysis to determine whether amino acid positions were subject to negative or positive selection pressures. Comparisons were made with the dengue reference strain Philippines H87/57, following Barrero and Mistchenko [43] general approach. Both genealogy-based, codon-site models Single Likelihood Ancestor Counting (SLAC) and the Random-Effects (REL) available in the HyPhy package (www. datamonkey.org) were used to estimate the nonsynonymous (dN) and synonymous (dS) rates of substitution [44][45][46]. Polymorphisms were analyzed with the program DNASP v5 [47] and subjected to the Tajima D statistic test, as proposed by Tajima [48] to evaluate deviations from the neutral expectation of molecular evolution [49].

Epitope and immunogenicity prediction
The amino acid sequences coding for the capsid, envelope, NS1, NS2a and NS5 proteins were aligned (http://multalin. toulouse.inra.fr/multalin/multalin.html) and a consensus sequence of each protein was submitted to B (http://www.cbs.dtu. dk/services/BepiPred/) and T (http://www.cbs.dtu.dk/services/ NetCTL/) epitope prediction algorithms. The Allele Frequency Net Database (http://www.allelefrequencies.net/) was used to select the HLA classes most predominant in the Southeast of Brazil and used to set up the NetCTL server. The mean value of the epitope propensity scores for each sequence was then classified and plotted according to its potential immunogenicity.

Results and Discussion
Complete genome analyzes define distinct DENV3 lineages circulating in SJRP The availability of high throughput, parallel DNA sequencing allows for comprehensive information on population genetics at the complete viral genome level. In this study, we determined the sequence of 3390 codon-sites of the complete polyprotein from 33 DENV-3 strains isolated during the 2006 and 2007 outbreaks.
The inferred phylogenetic relationship among DENV-3 strains was summarized in a Bayesian maximum credibility tree (MCC) shown in Figure 1. The majority of Brazilian samples grouped with genotype III, as did most taxa from Latin America isolated since 1994 [20,[50][51][52]. Likewise, the Brazilian strain EF629370 grouped within genotype I, as previously described by Schreiber et al. [53]. DENV-3 isolates, other than genotype III, have been reported in Colombia [54] and in the Brazilian states of Minas Gerais and Rondonia, likely imported from Asia [55,56]. The successive circulation of different genotypes is of public health concern given that distinct combinations could have different epidemic potentials and risk burdens [57].
The samples from SJRP all fell within genotype III. The MCC tree indicates that SJRP sequences group with Brazilian isolates from different regions in highly supported clades (Figure 1). Clade A included most of SJRP isolates and is more closely related to an isolate from Acre in 2004 (EF629367_BR_2004). Clade B includes four SJRP isolates grouping with (i) a sample from Central West Brazil (FJ913015_BR_2001) and (ii) the Brazilian sample Table 2. Non-neutral codons according to the assayed methods (SLAC and REL) at the specified significance levels (p = 0.05 and Bayes factor = 50, respectively).

Succession of DENV-3 genotypes in Brazil
Our data also provides information on the movement of DENV-3 (genotype III) in South America, because sample FJ850079_BR_2002 obtained in the northern part of Brazil grouped with Colombian and Venezuelan sequences suggesting a plausible transmission route.
DENV-3 (genotype III) was isolated for the first time in Brazil from an autochthonous case of DF in December 2000 from Nova Iguaçu, a city located in Rio de Janeiro (RJ) State, which is the second largest metropolitan area in the country, in 2001-02, a large outbreak of DENV-3 occurred in the neighboring city of Rio de Janeiro [60,61]. Following its introduction into RJ, the virus spread to the neighboring State of São Paulo (SP), and eventually to the municipality of SJRP. A previous introduction of DENV-3 took place in Brazil in Rondonia State (Amazon region), but sample EF629370_BR_2002 position in the MCC tree indicates it belongs to genotype V with its sister taxa coming from Asia. Nevertheless, most genotype III isolates from the north of Brazil grouped with the Caribbean strains (samples from Puerto Rico and Trinidad and Tobago), which suggests the existence of preferential routes of spread into Brazil, possibly facilitated by anthropogenic factors [51,62].

The dynamics of viral change
Sequence comparisons using the reference strain Philippines H87/57 (accession no. M93130) as an outgroup, revealed a total of 893 nucleotide changes, of which 132 corresponded to nonsynonymous substitutions. The nucleotide diversity p (that is the average number of nucleotide differences per site between two sequences) [63] within SJRP strains was 0.0052 6 0.001. Nevertheless, this p value was one order of magnitude lower than that for the complete data set used for phylogenetic inference (n = 0.038 6 2610 25 ), which suggests a genetic bottleneck following the introduction into a new human population, with the virus experiencing rapid in situ evolution afterwards. This is a likely scenario given that we sampled from the first documented DENV-3 outbreak in a city where the local human population is expected to be largely susceptible to this serotype. Likewise, the average value of Tajimas' D for all genes was 21.823 (p,0.05) The Tajima's D value near 22 also implies the outcome of negative selection, most likely due to functional constraints imposed on the viral proteins. In agreement with this hypothesis, the overall dN/dS value for the entire polyprotein coding region was 0.068 (estimated 95% CI from 0.057 to 0.080) that implied purifying (stabilizing) selection. Zanotto et al. [64] have argued that at a coarse-grain, the most common pressure acting on DENV is purifying selection with little evidence of recent adaptive evolution when comparing distantly related sequences in space and time, which would be expected due to synonymous change saturation at the serotype level. Moreover, given that arboviruses have adapted to diverse cell types it is generally considered that most mutations produced are likely to be deleterious [24]. Accordingly, the per-site SLAC dN/dS analyses revealed that only eight codons were nonneutral and negatively-selected (p,0.05) ( Table 2). Nevertheless, the presence of low frequency polymorphisms could account for adaptive changes, such as the 10 sites under significant positive selection (Bayes factor .50) that we found with the per-site REL dN/dS analyses.
Selected sites were located in both structural (S) and nonstructural (NS) proteins (Figure 2). Among the four non-neutral sites within the envelope protein, codon 355 coding for A Proline located in the central and dimerization domains, was under negative selection suggesting its structural importance. Likewise, positively-selected codons 404 and 581 in the envelope are near sites associated with disulfide bonds, suggesting adaptive changes near structurally important residues. Interestingly, Barrero and Mistchenko [43] also identified codons 581 and 893 (within the nonstructural NS1 gene) as positively-selected. King et al. [65] detected a statistically significant positively-selected site in the NS1 gene of DENV-3 genotype II. They also detected selective pressure in the C and NS4B genes of genotype I and in the E and NS3 genes of genotype II. We also found additional positivelyselected codon sites; 106 (Endoplasmatic reticulum anchor for the protein Capsid), 2191 (Lumenal domain of NS4A) and 3129 (catalytic domain of the RNA-directed RNA polymerase) that fell within domains that are important for membrane interactions/ rearrangements and virus replication.
Polymorphism analysis of individual genomic regions indicated that the NS5 and the envelope proteins varied the most with 24 and 22 replacement changes respectively, followed by NS1 (17), NS2A (11), NS4B (8) and Capsid (8) ( Table 3). Changes found in the glycoprotein M, the AdoMet-MTase region and the NS3 serine protease were mostly conservative, while replacements in the capsid were less conservative. Twiddy et al. [66] found codon sites 169 and 380 to be under positive selection in DENV-3, with the latter being located within the central and dimerization domains of the glycoprotein, which is involved in cell tropism. Likewise, Rodpothong and Auewarakul [25] reported one codon site in the prM gene and twelve residues in the E gene to be under positive selection. In summary, our findings along with others, suggest that genes other than the envelope may experience adaptive changes following introduction into a susceptible population, which implies a possible role for either host-immunity interaction or vector adaptations that require further study.

Short-term evolution of DENV-3 genotype III
Studies documenting selectively driven evolution in DENV are limited and mostly available for serotype 2 [25,[27][28][29]43,66], further studies to better characterize possible specific adaptive strategies of each serotype are needed. Therefore we looked at allele replacement in DENV 3 (genotype III).
The allele replacement process, from a population genetics standpoint, has an explicit genealogic representation, since frequency increase of an allele will be manifested by appearing  deeper as a synapomorphy in a viral genealogy [67]. Therefore, a synapomorphy indicates a fixed allele in the clade defined by it. To clarify this, we coded each of the ten positively-selected residues as a set of terminal unordered character states, represented as a single capital letter. The most parsimonious reconstruction sets of changes at each internal state in the phylogeny were calculated with MacClade v4.07 [68] using the MCC tree rooted with the reference strain Philippines H87/57. The reconstructed changes are illustrated in Figure 3. They indicate that positively-selected codons (404, 449, 581, 2191, 2879 and 3129) mapped to internal nodes as synapomorphies defining the genotype III. For example, codon 106 Met to Val, (within Protein C, the most nonconservative region in our study) characterized the SJRP clade A and codon 404 Pro to Leu characterized the SJRP clade B. Codon 893 appeared as an autapomorphy for the reference strain H87/57 (M93130). The synapomorphic changes differentiating clades A and B that we observed could reflect some aspects of the intra-genotypic interactions exposed to the human population [22].
Mondini et al. [35] described the pattern of spread of two DENV lineages circulating in SJRP, Brazil, during 2006. Based on 82 NS5 sequences, the phylogenetic analysis indicated that all samples were of DENV-3 and related to strains circulating in Martinique during 2000-2001. We have considerably extended that analysis by further showing that different lineages, closely relate to different Brazilian isolates, were co-circulating while experiencing distinct adaptive changes as indicated by their positively-selected synapomorphies. This finding is of relevance, since in several other viral systems, positively-selected changes are adaptive in nature and indicative of adaptive changes imposed by vectors, animal reservoirs or human hosts (i.e., immune escape and herd immunity) [69][70][71][72][73].
When we observed that DENV-3 sequences, previously clustered in clades A, B and C (Figure 3), we conducted in silico analysis to verify the influence of amino acids substitutions in the immunogenicity of capsid, envelope, NS1, NS2a and NS5 proteins. We found that amino acid substitutions diminished the B/T epitope propensity scores mainly in the DENV-3 sequences clustered in clades B and C (Figure 4). These analyses suggest that amino acid substitutions in sequences grouped in clade A increases its immunogenicity. As a consequence, DENV-3 strains grouped in clade A should expose more antigenic epitopes and more potentially recognized/accessed by the human immune system. We expected that the most immunogenic strains would be negatively selected in a given population. In this particular study we found exactly the opposite. Thus, one can hypothesize that the ADE phenomenon, which would be responsible for an increased pathogenicity and replication in host, can also affect the viral dynamics in a population level. We wonder whether this increased immunogenicity of clade A can induce an antibody dependent enhancement (ADE)-like response in the individual leading to high titers of the virus, which can facilitate the viral transmission to the mosquitoes and, consequently, increase this strain basic reproductive rate. Our findings are of relevance given that ADE may be manifested by challenges imposed by increasing viral immunogenic diversity.
While most studies document the broad geographic expansion of dengue, urban outbreaks dynamics deserves investigation because they contribute to a better understanding of the genetic diversity of strains during local transmission; knowledge that could be exploited for antiviral and vaccine development. In our analyses of 33 complete polyprotein-coding regions we revealed distinct lineages and detected polymorphisms sites that correlated with changes in the immunogenicity of several epitopes. The study transformation(s) in the same character occurs along that branch; X same transformation occurs in another branch (homoplasy) and transformation(s) in the same character occur along that branch. doi:10.1371/journal.pone.0063496.g003 Figure 4. Rank of immunogenicity of DENV-3 sequences. Thirty-three DENV-3 sequences were submitted to B and T epitope prediction using bioinformatics servers and further classified according to their immunogenic potential. Sequences previously clustered in B and C clades have epitopes, in envelope, NS1, NS2a and NS5, putatively less immunogenic when compared to the same epitopes of sequences clustered in clade A, except for sequence ACY70777.1. Substitutions in the capsid were irrelevant. doi:10.1371/journal.pone.0063496.g004 provides fine-grain information about the molecular epidemiology of dengue infection.