Transmission dynamics and control of Rickettsia rickettsii in populations of Hydrochoerus hydrochaeris and Amblyomma sculptum

Background Brazilian Spotted Fever (BSF), caused by the bacterium Rickettsia rickettsii, is the tick-borne disease that generates the largest number of human deaths in the world. In Brazil, the current increase of BSF human cases has been associated with the presence and expansion of capybaras Hydrochoerus hydrochaeris, which act as primary hosts for the tick Amblyomma sculptum, vector of the R. rickettsii in this area. Methods We proposed a semi-discrete-time stochastic model to evaluate the role of capybaras in the transmission dynamics of R. rickettsii. Through a sensitivity analysis, we identified the parameters with significant influence on the R. rickettsii establishment. Afterward, we implemented the Gillespie’s algorithm to simulate the impact of potential public health interventions to prevent BSF human cases. Results The introduction of a single infected capybara with at least one infected attached tick is enough to trigger the disease in a non-endemic area. We found that to avoid the formation of new BSF-endemic areas, it is crucial to impede the emigration of capybaras from endemic areas by reducing their birth rate by more than 58%. Model results were corroborated by ex-situ data generated from field studies, and this supports our proposal to prevent BSF human cases by implementing control strategies focused on capybaras. Conclusion The proposed stochastic model illustrates how strategies for the control and prevention of vector-borne infectious diseases can be focused on amplifier hosts management practices. This work provides a basis for future prevention strategies for other neglected vector-borne diseases.


Methods
We proposed a semi-discrete-time stochastic model to evaluate the role of capybaras in the transmission dynamics of R. rickettsii. Through a sensitivity analysis, we identified the parameters with significant influence on the R. rickettsii establishment. Afterward, we implemented the Gillespie's algorithm to simulate the impact of potential public health interventions to prevent BSF human cases.

Results
The introduction of a single infected capybara with at least one infected attached tick is enough to trigger the disease in a non-endemic area. We found that to avoid the formation of new BSF-endemic areas, it is crucial to impede the emigration of capybaras from endemic areas by reducing their birth rate by more than 58%. Model results were corroborated by exsitu data generated from field studies, and this supports our proposal to prevent BSF human cases by implementing control strategies focused on capybaras.

Conclusion
The proposed stochastic model illustrates how strategies for the control and prevention of vector-borne infectious diseases can be focused on amplifier hosts management practices. PLOS

Introduction
Rickettsia rickettsii is the etiological agent of the Brazilian spotted fever (BSF), the deadliest spotted fever in the world. This infection is partially pathogenic to Amblyomma sculptum ticks, main vectors of the R. rickettsii in South America [1,2], generating a drop in the infection rate with each tick generation [3]. In addition, A. sculptum ticks are unable to maintain the R. rickettsii infection in successive generations by transovarial and transstadial transmissions [4]. Hence, the maintenance of R. rickettsii depends on a constant introduction of new susceptible animals (i.e., newborn of vertebrate hosts), which act as amplifier hosts and guarantee the constant creation of new cohorts of infected ticks [4][5][6][7]. This suggests that control strategies focused on vectors [8,9] are not enough for the prevention of this tick-borne infectious disease. In Brazil, the capybara Hydrochoerus hydrochaeris acts as the amplifier host of R. rickettsii infection [10,11]. In southeastern Brazil, both capybaras and BSF occurrences have increased significantly over the last three decades [5,12]. In turn, these occurrences have been spatiotemporally associated with rising production and spatial expansion of sugarcane crops, the main food source of capybaras [13]. In BSF-endemic areas, population densities of capybaras have reached numbers up to 40 times higher than those recorded in natural environments such as the Amazon and Pantanal [14]. However, the effectiveness of control strategies focused on this amplifier host and their impact on the transmission of R. rickettsii are unknown.
The dynamics of complex transmission cycles, such as tick-borne infectious diseases have been broadly analyzed. Hudson et al. published the first deterministic model representing the dynamics of the Louping-ill disease, an acute viral zoonosis which mainly affects sheep [15]. O'Callaghan et al. put forward tick-borne diseases dynamics considering the potential effect of a vaccination program for Ehrlichia ruminantium [16]. Similarly, Rosà et al. developed a model for the tick-borne encephalitis virus including the tick stages and different transmission routes [17,18]. Despite their important role in the understanding of vector-borne diseases dynamics, unfortunately they do not include the seasonal behavior of the vector life cycle, neither the onset nor the extinction of these diseases. As it is well established, stochastic models should be used in phenomena that do not satisfy the law of large numbers such as large communities with minor outbreaks [19]. In fact, the extinction of endemic diseases can only be analyzed with stochastic models, since extinction occurs when the epidemic process deviates from the expected level [19].
In this work, we propose a discrete-state semi-discrete-time stochastic framework to evaluate the role of capybaras in the transmission dynamics of R. rickettsii. We identify the parameters with significant influence on the R. rickettsii establishment and subsequently evaluate the impact of potential public health interventions to prevent BSF in humans.

Model
The R. rickettsii dynamics in populations of H. hydrochaeris and A. sculptum is represented in Fig 1. Our model was adjusted to a semi-discrete time dynamics in order to consider the 1-year life cycle of the tick A. sculptum, which is primarily controlled by the larval behavioral diapause [20]. Thus, larvae exclusively quest and feed from April to July for 110 days, nymphs from July to October for 104 days and adults particularly quest, feed and reproduce from October to March for 151 days, as shown in Fig 1. The tick population is classified according to the life cycle stages as larvae (L), nymphs (N) or adults (A), which could be detached from a capybara or attached to it. When a tick gets infected by an infected capybara, it remains infected until it dies. Thus, each A. sculptum stage is also classified according to whether it is susceptible (S) or infected (I). In this way, the population of each tick stage is represented by three indexes, where the first index denotes the infection status, the second index denotes the tick life cycle stage, and the third index denotes detachment (D) or attachment (A).
Attached ticks at the larval and nymph stages detach at the respective rates of θ L and θ N . Furthermore, detached ticks at larval, nymph and adult stages can attach at rates α L , α N and α A , and can die at rates δ L , δ N and δ A . The production rate ρ of larvae is assumed to be proportional to the total number of adult attached ticks AA.
On the other hand, capybaras are classified as susceptible S C , infected I C and recovered R C , as shown in Fig 1. They reproduce at a constant rate μ C , die at rate δ C and recovered at rate γ. All capybaras have the same susceptibility and there is no increased death rate of infected individuals due to disease. Susceptible capybaras can get infected by an attached larva, nymph or adult tick at a rate of λ L , λ N and λ A , respectively. Once capybaras are infected, they keep the R. rickettsii in the bloodstream for 7 to 10 days [11], during which the infection of new susceptible larvae or nymphs that feed on it can occur at rates β L and β N , respectively. After this period, capybaras become immune to the disease. We do not consider a transmission rate from infected capybaras to susceptible adult ticks β A , because of the time of R. rickettsii infection in ticks is greater than the time of laying. Thus, eggs can be infected only if adult ticks become infected during the larval or nymph stages [4]. Additionally, we consider a vertical transmission and hence infected adult ticks can produce infected detached larvae at rate ρ II . In this way, there are 25 distinct reactions in the stochastic model which are listed in Table 1. The equivalent deterministic equations associated with the stochastic reactions are presented in the S1 File.

Time series simulations
The stochastic process detailed in Table 1 was simulated in the R programming language employing the Gillespie's algorithm. For each simulation we performed the number of iterations needed to get a stable performance.
We assumed a finite number of individuals distributed over a finite set of discrete states. Given an initial time t 0 and an initial population state X(t 0 ), the Gillespie's algorithm allows us to generate the time-evolution trajectory of the state vector X(t) (X 1 (t), . . ., X N (t)), where X i (t) is the population size of state i at time t and N is the number of states. Changes in the number of individuals in each state occur due to reactions between interacting states. The states interact through M reactions R j , where j = 1, . . ., M denotes the jth reaction. A reaction is defined as any process that instantaneously changes the population size of at least one state [21,22]. The time step to the next reaction was then determined as t ¼ 1 a 0 ðxÞ lnð1=r 1 Þ, where α 0 (x) = Sα j (x) and the index of the next reaction to execute R j is the smallest integer j satisfying j ¼ S j i¼1 a i ðxÞ > r 2 a 0 ðxÞ [21,22]. All parameters were estimated using previously generated data from ex situ field works in southeastern Brazil [4,5,11,14,23]. It is noteworthy that the natural capybara birth rate μ assumes the value of 70% of adults, 64% of females, a litter size mean of 4.2 pups, 1.23 births per female per year and a pregnancy success of 85% [14,23]. If capybaras die at an exponential rate, then δ C is the fraction required to die each day. The birth rate of a susceptible tick was determined assuming a female weight of 500mg [4], CEI (mg egg mass/mg engorged female × 100) of 48.4% [4],18.8 eggs/ 1 mg of eggs [24] and hatching success of 68% [4].
Likewise, the birth rate of an infected tick was determined assuming a transovarial transmission of 42.8% [4], filial infection rate of 50% [4], female weight of 372.20 mg [4], CEI of 39.55% [4], 18.8 eggs/1 mg of eggs, [24] and hatching success of 44.2% [4]. We consider a population of 20 adult female ticks per capybara and the groups of capybaras were restricted to 50 individuals [14,23,25]. A full list of the model's parameters used in the simulations is given in Table 2.

Sensitivity analysis
To quantify the impact of the variation of each parameter on the output of the BSF model, we combined uncertainty through the Latin hypercube sampling (LHS) with the robust Partial rank correlation coefficient (PRCC) method [26,27]. Initially, we obtained a random Table 1. Events and reactions of the tick-capybara-disease stochastic process.

Event Reaction
Birth of capybara Transmission from an infected capybara to a susceptible larvae Transmission from an infected larvae to a susceptible capybara Stage change from susceptible larvae to detached nymph Stage change from infected larvae to detached nymph Recovery rate of capybara Death of a susceptible capybara Emigration of a susceptible capybara S C ! C S C S C À 1 Death of an infected capybara Emigration of an infected capybara Death of a recovered capybara Emigration of a recovered capybara Engorgement rate of an infected nymph I ND ! a N I ND ðI ND À 1Þ þ ðI NA þ 1Þ Transmission from an infected nymph to a susceptible capybara Transmission from an infected capybara to a susceptible nymph S NA ! parameter distribution divided into one hundred equal probability intervals, which were then sampled. Thus, a LHS matrix was generated with one hundred rows for the number of simulations (sample size) and six columns corresponding to the number of varied parameters (α, μ C , λ, γ C , δ C , C ). The parameters α and λ were varied according to the life cycle season. BSF model solutions were then simulated using each combination of parameter values through ten year simulation. One thousand model outputs were obtained and the parameter and output values were transformed into their ranks. Subsequently, we computed the PRCCs between each parameter and the average infected population size.
Sensitivity analysis of scenario A show that the capybaras birth rate variation had the greatest impact on the average infected population size. In the nymph and adult ticks seasons, when the infection of capybaras is greater, the correlation value between the birth rate of capybaras and the average infected population is PRCC>0.6, being significant in both seasons (Fig 3). Since in BSF-endemic areas of southeastern Brazil the population of capybaras is growing [5,[12][13][14], theoretically, the epidemic will survive forever [31]. Considering μ = δ C + , where μ C represents birth, δ C death and emigration rate of capybaras, the population growth generates a high rate of emigration and consequently the spread of the disease. In this way, to better understand the effect of changes in capybara population, we investigate the host-tick-infection dynamics for different values of μ C , under three additional scenarios: 1) an endemic area with a stable state of capybaras (births equal deaths, no emigration and no importation of the disease from outside); 2) an endemic area with a decrease of 80% and 3) 90% in the birth rate of capybaras.
Contrary to the expected behavior [32], when the population is in a stable state, the disease does not disappear. Indeed, the proportion of infected individuals remains constant over time from year 2 onward as shown in Fig 2B. In this scenario, the average fraction of infected capybaras is 4.6% (95% CI = 0.01%-19.3%), 6.3% (95% CI = 0.01%-23.1%) and 7.4% (95% CI = 0.01%-25.4%) in the larvae, nymphs and adults seasons, respectively. To ensure births equals deaths (μ C = δ C ) and thereby guarantee no emigration ( = 0), it is necessary to control the capybara's birth rate in 58%. Thereby, though the disease does not disappear, this birth rate guarantees that the R. rickettsii will not spread from BSF-endemic areas.
Thus, for the elimination of the R. rickettsii from endemic areas, a decrease in the capybara's birth rate to values lower than 0.0021 is necessary. However, this can lead to a decline in capybaras population over time. When a decrease of 80% (μ C = 0.0011) in the birth rate of capybaras is considered, infected individuals tend to disappear in the fourth year along with a decrease in the total population size, as shown in Fig 2C. Otherwise, when a decrease of of 90% (μ C = 0.0005) in the birth rate of capybaras is considered, infected ticks and capybaras tend to disappear from the second year as shown in Fig 2D. Strategies to reduce the birth rate of capybaras include the reduction of the carrying capacity, their removal, either by euthanasia or regulated hunting, and their reproductive control. As capybaras natality depends primarily on the availability of food sources, as is typically the case for rodents [33], the reduction of the carrying capacity in BSF-endemic areas is a plausible strategy to reduce their birth rate [34]. Polo et al. found a spatiotemporal relationship between the occurrence of BSF and the increment and expansion of sugarcane crops, the main food source of capybaras in southeastern Brazil and the most important agricultural product in the region [13]. Furthermore, in this area, there is a constant availability of water sources, which allow the establishment of capybaras, as this is a semiaquatic vertebrate that depends on water sources for thermic regulation, reproduction, and predator protection [33]. Certainly, controlling these aspects is not feasible.
Additionally, because of the constant increment and abundance of vital resources offered by the environment in BSF-endemic areas of southeastern Brazil, it is important to consider that in response to the removal or elimination of recovering capybaras from theses areas, a reintroduction of susceptible animals can occur [13]. This, along with the long survival of unfed A. sculptum in pastures [35], and the fact that just one infected capybara with a single infected tick attached is sufficient to establish an endemic area, can cause a rapid spread of the disease and consequently an increased risk of transmission to humans.
Reproductive control of capybaras through deferentectomy and ligation of fallopian tubes was previously tested in southeastern Brazil [36]. It was observed that the reproductive management did not negatively influence individual or collective behavioral aspects, with the animals defending their territory and not migrating [36]. The sterilization of capybaras has already been authorized as a way to prevent BSF human cases in a small endemic area of southeastern Brazil [37].
Future studies encompassing field data should be performed to evaluate the role of alternative reservoirs in the dynamics of R. rickettsii. This is a limitation of our work since we considered BSF-endemic areas of southeastern Brazil, where capybaras are the major, though not the exclusive, hosts for either larval, nymphal or adult stages of A. sculptum. In most BSF-endemic areas, the only medium-to large-sized animal species is the capybara; therefore it is the only host for the A. sculptum adult stage, and consequently, the only suitable host species to sustain an A. sculptum population in the area [7,10,38]. However, there have been some previous reports of A. sculptum immature stages (larvae and nymphs) on small animals (wild mice, marsupials, birds) which usually share the same habitat with capybaras [39,40]. Nevertheless, comparing to capybaras, the amount of larvae or nymphs that feed on these animals is minimal; i.e., while hundreds to thousands larvae and nymphs are commonly found feeding on a single capybara [41], we usually find less than 10, or exceptionally a few dozen A. sculptum ticks on individual small animals [39][40][41]. Indeed, this condition is also favored by the fact that all active stages of A. sculptum tend to host-quest on vegetation above 15 cm from the soil, precluding their direct contact with small mammals or ground-feeding passerine birds [20,42].
This work offers an alternative for the planning of prevention strategies for tick-borne neglected infectious diseases. The planning of these type of interventions is usually performed by public entities through heuristic techniques, such as trial and error methods, without obtaining optimal results. The results of our work accurately match with data previously generated from field studies on the R. rickettsii dynamics in southeastern Brazil and will potentially allow the formulation of public policy to prevent BSF human cases based on the control of capybara population growth. Furthermore, this work will provide a basis for the planning of prevention programs for other neglected vector-borne infectious diseases.