Evolutionary history of dog rabies in Brazil

Received 18 August 2010 Accepted 6 October 2010 Center for Information Biology and DNA Data Bank of Japan, National Institute of Genetics, 1111 Yata, Mishima, Shizuoka 411-8540, Japan Nihon University Veterinary Research Center, 1866 Kameino, Fujisawa, Kanagawa 252-8510, Japan Department of Preventive Veterinary Medicine and Animal Health, Faculty of Veterinary Medicine and Zootechny, University of São Paulo, Av. Prof. Dr Orlando Marques de Paiva 87, Cidade Universtiátria, São Paulo 05508-000, Brazil


INTRODUCTION
Rabies virus (RABV), which belongs to genotype 1 of the genus Lyssavirus in the family Rhabdoviridae, is capable of infecting all warm-blooded animals and causes death by lethal encephalitis.Mammals belonging to the order Carnivora, such as dogs, cats, foxes, raccoons and skunks, and the order Chiroptera, bats, are the principal reservoirs of RABV.These two host groups form the basis of seemingly distinct phylogenetic clusters of RABV, that is, the Carnivora-and Chiroptera-related RABV clusters (Badrane & Tordo, 2001).The World Health Organization estimates that approximately 55 000 people die of rabies each year around the world, most of whom become infected with RABV after being bitten by a rabid dog (Knobel et al., 2005).
Rabies is endemic to Brazil, and dogs and vampire bats are the principal RABV transmitters to both humans and livestock (Ito et al., 2001).Although the incidence of dog rabies has decreased in urban areas because of effective vaccination programmes, the frequency of RABV transmission from dogs to other dogs or humans in Brazil is currently the highest in the Americas (Belotto et al., 2005).In addition, RABVs have been isolated from wild canids, particularly foxes in north-eastern Brazil (Bernardi et al., 2005;Carnieli et al., 2008Carnieli et al., , 2009;;Favoretto et al., 2006;Kobayashi et al., 2007).Phylogenetic analyses suggested that inter-species RABV transmission occurred between dogs and foxes in the region, and that the viruses are circulating among the dog and fox populations independently (Bernardi et al., 2005;Carnieli et al., 2008Carnieli et al., , 2009;;Favoretto et al., 2006;Kobayashi et al., 2007).Therefore, dogs are considered to be a critical transmitter for spreading the virus in Brazil, and thus it is important to understand the spatio-temporal characteristics of dog RABV in Brazil in order to develop a strategy for rabies control.
The first report of dog rabies in South America was in Peru in 1803, which was colonized by Spain from 1542 to 1824 (Baer, 2007;Johnson, 1952).Shortly thereafter, in 1806, an outbreak of rabies was observed among hunting dogs that were imported by English officers into Argentina, which was also a colony of Spain from 1516 to 1810 (Baer, 2007;Johnson, 1952).These incidents were followed by other incidents of dog-transmitted rabies in the European colonies of South America.Therefore, it has been suggested that dog rabies was introduced into South America from Europe at the time of European colonization through the The GenBank/EMBL/DDBJ accession numbers for the complete glycoprotein gene sequences of the nine RABV isolates are AB518487-AB518495.importation of dogs infected with RABV (Johnson, 1952;Nadin-Davis & Bingham, 2004).This hypothesis, of rabies having a European origin, is supported by phylogenetic analyses of RABV.Dog RABV isolates collected from Africa and the Americas, both of which were colonized by Europeans, cluster together with isolates from Europe to form what has been referred to as the 'cosmopolitan cluster' (Badrane & Tordo, 2001;Bourhy et al., 2008;Nadin-Davis & Bingham, 2004;Smith et al., 1992).The cosmopolitan cluster also contains carnivore-related RABV isolates from the Middle East, Asia and Russia.It has been estimated that the most recent common ancestor of the cosmopolitan cluster existed 284-504 years ago (Badrane & Tordo, 2001).However, it was unclear when RABV was introduced from Europe to other countries within the cosmopolitan cluster.
From the observation that RABV isolates from different American countries formed distinct clades within the cosmopolitan cluster, dog RABV was suggested to have been introduced into the Americas on multiple occasions (Nadin-Davis & Bingham, 2004).The purpose of the present study was to estimate the time of introduction of dog RABV from Europe into Brazil, as well as that of the spread of dog RABV within Brazil.

Phylogenetic analysis of carnivore-related RABV sequences
A phylogenetic tree, constructed using 122 complete glycoprotein (G)-gene sequences of RABV isolates, is shown in Supplementary Fig. S1.The tree shows how these RABV isolates formed clusters according to their geographical origins (Supplementary Fig. S1).RABV isolates from Brazilian carnivore species formed a monophyletic group at node B and clustered with RABV isolates from carnivore species collected in Europe, the Middle East, China, Africa, North America and Mexico at node A, to form the cosmopolitan cluster (Fig. 1).These results suggest that dog RABVs were introduced into Brazil from European countries during the period between nodes A and B (Fig. 1).
The group consisting of RABV isolates from Brazilian carnivore species was divided into subgroups I and II, which formed clusters at nodes C and D, respectively.Subgroup I consisted mainly of dog and cat RABV isolates collected in the states of Sa ˜o Paulo, Goı ´as, Minas Gerais, Maranha ˜o and Mato Grosso (Fig. 1 and Supplementary Fig. S2, available in JGV Online), suggesting that RABV had occasionally been transmitted from dogs to cats in central to northern Brazil after divergence at node C.
Subgroup II was composed of three lineages (IIa, IIb and IIc), with lineage IIa diverging first from lineages IIb and IIc.Lineage IIa consisted of RABV isolates from dogs collected in the state of Maranha ˜o, and lineages IIb and IIc consisted of RABV isolates collected from foxes (Dusicyon sp.) in the state of Paraı ´ba (Supplementary Fig. S2), suggesting that RABV transmission from dogs to foxes occurred in the northern areas of Brazil on the ancestral branch of lineages IIb and IIc, and that RABV was maintained in foxes thereafter (Fig. 1).

Evolutionary rates and divergence times of Brazilian carnivore RABV
In order to investigate the epidemiological history of carnivore RABV in Brazil, the nucleotide substitution rate of the RABV G gene and the divergence times at nodes A-D in Fig. 1 were estimated using 72 RABV isolates, which formed a single cluster at node A.
The rate constancy among lineages (molecular clock) was rejected for the complete G-gene sequences in the tests using Bayes factor (BF) (relaxed log-normal or exponential model versus strict model: BF.20) (see Methods).However, in the analysis of fourfold-degenerate sites, no significant difference was observed among the three models in the Bayesian framework (BF,10).Thus, in the estimation of the evolutionary rates and divergence times, the relaxed log-normal and exponential models for rate heterogeneity among lineages were used when analysing the complete G-gene sequences, while the strict molecularclock model was used when analysing the fourfolddegenerate sites.Both the constant-size model and exponential-growth model were employed for the population demographic history of each dataset.
In the Bayesian Markov-chain Monte Carlo (MCMC) analysis of the complete G gene, the rate of nucleotide substitution was estimated to be 6.26610 24 -1.3610 23 substitutions per site per year, and the divergence times were estimated to be 1878-1920 at node A, 1909-1961 at node B, 1957-1976 at node C and 1944-1984 at node D (Table 1).Using the fourfold-degenerate sites, the rate of nucleotide substitution was estimated to be 1.73610 23 -2.09610 23 substitutions per site per year, and the divergence times were estimated to be 1885-1906 at node A, 1923-1938 at node B, 1960-1966 at node C and 1952-1962 at node D. The rates estimated by using the exponential population-growth model were slightly higher than those estimated by using the constant-size model.In addition, the rates estimated by using the exponential rate heterogeneity model were approximately twice as great as those estimated by using the relaxed log-normal rate heterogeneity model (Table 1).However, nucleotide substitution rates estimated in the present study were roughly in agreement with the estimates of previous studies (5.5610 24 -1.28610 23 substitutions per site per year) (Badrane & Tordo, 2001;Bourhy et al., 2008;Hanada et al., 2004;Hughes et al., 2004;Talbi et al., 2009), and the 95 % highest probability density (HPD) values of the divergence times largely overlapped among the estimates for each node (Table 1).

DISCUSSION
The Portuguese discovered and started to colonize Brazil in 1500.Since then many immigrants have migrated to the country, but Portugal was the main source of immigrants to Brazil during colonization (1500-1822) (Carvalho-Silva et al., 2001).However, after the Portuguese court opened the seaports of Brazil for trading in 1808 the number of immigrants increased markedly, many immigrants from other countries, mostly Europeans, arrived in Brazil; the influx peaking in approximately 1900 (Carvalho-Silva et al., 2001;IBGE, 2000).
In the present study, it was inferred that dog RABV was introduced into Brazil in the late-19th to early-20th century, which was later than the colonization period but corresponded to the time period of increased immigration from Europe.It should be noted that similar estimates for divergence time were obtained by using different datasets and evolutionary models from those in the present study, which suggests that the results are generally robust.Since dog rabies was most prevalent in Europe in the 19th century (Johnson, 1952), dog RABVs may have been introduced into Brazil from Europe by importing rabid dogs along with immigrants.
In the cosmopolitan clade of dog RABV, Brazilian strains appeared to be located outside of non-Brazilian cosmopolitan strains, suggesting that the latter strains also spread in Fig. 1.The cosmopolitan cluster of carnivorerelated RABV isolates.This enlarged view of the cosmopolitan cluster is also shown within Supplementary Fig. S1 (surrounded by a dashed line) along with additional information.Interior nodes, where the divergence times were estimated, are labelled A-D.Each isolate collected in Brazil is described as follows: name/state isolated in/year isolated.State abbreviations: GO, Goı ´as; MA, Maranha ˜o; MG, Minas Gerais; MT, Mato Grosso; PB, Paraı ´ba; RJ, Rio de Janeiro; and SP, Sa ˜o Paulo.Each isolate collected in another country is described as (name/year isolated).Bootstrap probabilities are provided for interior branches when they are .70%.Bar, nucleotide substitutions per site.the late-19th to early-20th century or later (Nadin-Davis & Bingham, 2004).The importation of dogs from Europe to the Americas was known from the time of the second voyage of Columbus ( 1493), and incidents of dog rabies in South America have been observed since the early 19th century (Baer, 2007).Thus, although it is possible that the first importation of dog RABV to South America occurred before the time estimated in the present study, dog RABVs currently circulating in the countries of South America may have originated from the viruses that were introduced from Europe by multiple importations of rabid dogs since the late-19th century.

IP
The phylogenetic analysis in the present study showed that dog RABVs introduced to Brazil were separated into two subgroups (subgroups I and II), and then spread to susceptible animals.Our estimates suggest that RABVs of subgroup I spread rapidly to domestic animals across the country within the past 50 years, while RABVs of subgroup II have been circulating among wild dogs and foxes independently in northern Brazil since the virus was transmitted from dogs to indigenous foxes in the mid-20th century.These estimates correspond to the time when the incidence of dog rabies in Brazil was particularly high (Nadin-Davis & Bingham, 2004).It is well known that the Brazilian economy developed significantly since 1930, with rapid industrial growth and agricultural expansion occurring in this period (Fishlow, 1980).The pioneer settlers, who were immigrants, moved to the virgin territories of Brazil in order to expand their settled land, and as a result of urbanization in Brazil this process accelerated during the 20th century (Fishlow, 1980).Thus, the development and urbanization of Brazil in that time period may have promoted infection by, and spread of, dog RABV in Brazil.

METHODS
Sequences.Complete G-gene sequences (1575 nt) of nine RABV isolates, which were obtained from seven dogs (Canis lupus familiaris) and two cats (Felis silvestris catus) in Brazil and stored at 280 uC (Ito et al., 2001;Kobayashi et al., 2007), were determined using methods described previously (Sato et al., 2006) (GenBank accession nos AB518487-AB518495).In addition, 110 complete RABV-G-gene sequences were retrieved from GenBank: 39 from Brazil and 71 from 14 other countries (Supplementary Tables S1 and S2, available in JGV Online).Since all RABV isolates belonged to the carnivore-related RABV phylogenetic cluster (Supplementary Fig. S1), we considered that they were all derived from carnivores.The complete G-gene sequences of BR-MM1, BR-EF2 and USA7-BT, which belong to the Chiroptera-related RABV cluster (Supplementary Fig. S1), were also obtained from GenBank (AB449216, AB383169 and AF298141).
Estimation of nucleotide substitution rates and divergence times.In order to estimate the rates of nucleotide substitution and RABV divergence time, we analysed 72 complete G-gene sequences isolated at dates ranging from 1981 to 2005 (Supplementary Tables S1  and S2), which formed a single cluster with high bootstrap values and included Brazilian isolates (Fig. 1).Since the molecular clock usually holds for synonymous substitution (Nei & Kumar, 2000), fourfolddegenerate sites of the G gene (194 nt) were also analysed separately using the same set of 72 sequences.
The best-fit models of nucleotide substitution were identified as the transitional model with rate heterogeneity among sites (TIM+! 4 ) (hLRT) and GTR+I+! 4 (AIC) for the complete G-gene sequences, and the three-parameter model with unequal base frequencies (K81uf) (hLRT) and the transversional model (AIC) for the fourfold-degenerate sites.In the present study, GTR+I+! 4 and GTR were used in the analyses of the complete G-gene sequences and the fourfold-degenerate sites, respectively, because these models were more general than the models selected in MODELTEST and were available in the computer program used for estimating the evolutionary rates and divergence times, as indicated below.
The rates of nucleotide substitution (per site per year) and divergence time were estimated simultaneously using the MCMC method available in the BEAST package (Drummond & Rambaut, 2007).The rate heterogeneity among lineages as well as the population demographic history is modelled in this method.The rate of nucleotide substitution was assumed to follow the strict model, which assumes the molecular clock, or the relaxed log-normal model and exponential models, which do not assume a constant rate across lineages.Population size was assumed to be constant (constant model) or to follow an exponential-growth model for each dataset.
The logistic-growth coalescent model, which was also available as a demographic model in BEAST, was not applicable to the data analysed in the present study because the parameter estimates did not converge to reliable values during the computation (data not shown).
For each model, the MCMC analysis was run for 50 000 000 steps.After discarding the initial 5 000 000 steps as burn-in, the parameter estimates were sampled every 1000 steps.In all cases, the convergence of parameter estimates was assessed using TRACER (Suchard et al., 2001), with statistical uncertainly reflected in the interval of the 95 % HPD.The rate-constancy (molecular clock) and rate-heterogeneity models were then compared for each dataset in a Bayesian framework using the BF in TRACER.In general, a BF.20 is considered to be strong support for the favoured model (Drummond & Rambaut, 2007).

Table 1 .
Estimates of nucleotide substitution rate and divergence time for rabies viruses