A metapopulation model for sylvatic T . cruzi transmission with vector migration

A vector-borne disease of major concern in the Americas, transmitted via insect vectors from the subfamily Triatominae, is Chagas’ disease. Chagas’ disease, discovered in 1909, is widespread in Mexico, Central America, and throughout Latin America. An estimated 8 to 11 million people are currently infected, with many unaware of their infection [56]. However, in the United States, fewer than 10 cases of autochthonous transmission have been reported [5]. Although more attention is being given to Chagas’, incidence of the disease remains underreported, and Chagas’ is classified as a neglected parasitic infection in the United States [10]. Although there have been few human cases in the U.S., the disease remains endemic in sylvatic cycles throughout Mexico and the United States. Sylvatic transmission cycles are vector-host cycles that occur naturally in the wild. In the United States, triatomine vectors are found in 26 states [29], involving approximately 11 triatomine vector species (with 8 of the 11 in Texas) and over 100 mammalian species. In the United States, some of the most common sylvatic hosts include opossums (Didelphis virginiana) and raccoons (Procyon lotor) in the southeastern parts of the country and woodrats (Neotoma micropus) in Texas (extending also into northern parts of Mexico). Other species, such as dogs, armadillos, and skunks, and chickens have also been noted as relevant species in sylvatic settings (with canines part of some domestic and peridomestic cycles [9, 58]). Of the 11 vector species and mammalian species listed here, we identify two primary vector species associated with these hosts in the southeastern U.S., Triatoma sanguisuga, found all along the southeastern Atlantic coast from Florida into central Texas, and Triatoma gerstaeckeri, found mostly from central Texas south into states in northern Mexico [25]. In addition to the complex vector-host cycles, there are also multiple strains of T. cruzi circulating in these populations. There are 6 known strain types of T. cruzi, types I-VI, of which types I and IV are circulating in the United States. There are distinct differences between the strains, from host specificity to virulence. T. cruzi I, associated with Chagas’ disease, is the primary strain circulating in Mexico (also found in hosts in the U.S.), while type IV is almost exclusively found in the United States [44]. Recent work by Crawford and Kribs-Zaleta [14], provides an understanding of how local vector dispersal can be described in terms of global effects, so that we may now consider a model that describes sylvatic cycles of T. cruzi over a large geographic area. Because T. cruzi is maintained in sylvatic cycles, we recognize the need to study the spatial spread of the disease, especially in North America, where risk of Chagas has only recently been studied [25, 22]. Here we investigate several models of T. cruzi, incorporating multiple modes of transmission and multiple patches. More specifically, we wish to focus our efforts on the effects of vector migration on sylvatic T. cruzi strain type IV transmission in two different North American host-vector cycles. To date, the majority of mathematical models for Chagas disease have been studied in humans and vectors, rather than the animal hosts. Velasco-Hernández [53] modeled Chagas in humans using a model structure similar to the Ross-Macdonald malaria model [45, 32], but included another infectious compartment for chronically ill humans. Since infection with T. cruzi is maintained in reservoir (sylvatic) hosts and human transmission cycles cannot be sustained without them [20], recently more attention has been given to the spread of the T. cruzi parasite in Triatoma vectors and associated animal hosts [28, 26]. In each model, Kribs-Zaleta uses a deterministic SI model with one host and one vector to study the effects of alternative transmission modes for


Introduction
A vector-borne disease of major concern in the Americas, transmitted via insect vectors from the subfamily Triatominae, is Chagas' disease.Chagas' disease, discovered in 1909, is widespread in Mexico, Central America, and throughout Latin America.An estimated 8 to 11 million people are currently infected, with many unaware of their infection [56].However, in the United States, fewer than 10 cases of autochthonous transmission have been reported [5].Although more attention is being given to Chagas', incidence of the disease remains underreported, and Chagas' is classified as a neglected parasitic infection in the United States [10].
Although there have been few human cases in the U.S., the disease remains endemic in sylvatic cycles throughout Mexico and the United States.Sylvatic transmission cycles are vector-host cycles that occur naturally in the wild.In the United States, triatomine vectors are found in 26 states [29], involving approximately 11 triatomine vector species (with 8 of the 11 in Texas) and over 100 mammalian species.In the United States, some of the most common sylvatic hosts include opossums (Didelphis virginiana) and raccoons (Procyon lotor) in the southeastern parts of the country and woodrats (Neotoma micropus) in Texas (extending also into northern parts of Mexico).Other species, such as dogs, armadillos, and skunks, and chickens have also been noted as relevant species in sylvatic settings (with canines part of some domestic and peridomestic cycles [9,58]).Of the 11 vector species and mammalian species listed here, we identify two primary vector species associated with these hosts in the southeastern U.S., Triatoma sanguisuga, found all along the southeastern Atlantic coast from Florida into central Texas, and Triatoma gerstaeckeri, found mostly from central Texas south into states in northern Mexico [25].In addition to the complex vector-host cycles, there are also multiple strains of T. cruzi circulating in these populations.There are 6 known strain types of T. cruzi, types I-VI, of which types I and IV are circulating in the United States.There are distinct differences between the strains, from host specificity to virulence.T. cruzi I, associated with Chagas' disease, is the primary strain circulating in Mexico (also found in hosts in the U.S.), while type IV is almost exclusively found in the United States [44].
Recent work by Crawford and Kribs-Zaleta [14], provides an understanding of how local vector dispersal can be described in terms of global effects, so that we may now consider a model that describes sylvatic cycles of T. cruzi over a large geographic area.Because T. cruzi is maintained in sylvatic cycles, we recognize the need to study the spatial spread of the disease, especially in North America, where risk of Chagas has only recently been studied [25,22].Here we investigate several models of T. cruzi, incorporating multiple modes of transmission and multiple patches.More specifically, we wish to focus our efforts on the effects of vector migration on sylvatic T. cruzi strain type IV transmission in two different North American host-vector cycles.
To date, the majority of mathematical models for Chagas disease have been studied in humans and vectors, rather than the animal hosts.Velasco-Hernández [53] modeled Chagas in humans using a model structure similar to the Ross-Macdonald malaria model [45,32], but included another infectious compartment for chronically ill humans.Since infection with T. cruzi is maintained in reservoir (sylvatic) hosts and human transmission cycles cannot be sustained without them [20], recently more attention has been given to the spread of the T. cruzi parasite in Triatoma vectors and associated animal hosts [28,26].In each model, Kribs-Zaleta uses a deterministic SI model with one host and one vector to study the effects of alternative transmission modes for T. cruzi, namely vector consumption by animal hosts and vertical transmission in hosts.Results show that vertical transmission is not enough to maintain the infection cycle alone, but vertical transmission along with even an inefficient host-vector transmission cycle can sustain the T. cruzi infection cycle.Due to the nature of transmission of vector-borne disease, in which vectors and hosts (especially in sylvatic settings) may be easily affected by weather (mainly temperature and humidity) and landscape, spatial spread is a key element in studying a vector-borne disease.
Spatial spread of a disease can be modeled using continuous or discrete space.The majority of mathematical models involving the spatial spread of infectious diseases in continuous time and space are modeled using reactiondiffusion systems taking the form of a system of partial differential equations.Some studies incorporating spatial spread are the spread of rabies in the fox population [38], and the vector borne diseases dengue [33] and West Nile virus [31].The results of such systems are generally described using traveling waves which describe the process of the spread of the disease, most often over a homogeneous landscape.In each model, movement of either hosts or vectors (or both) is considered, with the underlying assumption that the movement is random.
Other types of models incorporating spatial spread include multi-patch metapopulation models in which movement occurs between n patches.Several models have been studied, including a multi-species model by Arino et al. [3], in which analytical results are given for several multi-species, multi-patch models.In the study, a formula is derived for the basic reproductive number R 0 for multiple species and multiple patches and global stability for the disease-free equilibrium is established for R 0 < 1. Allen et al. [2] gives a 3patch model of hantavirus spread in reservoir and spillover species in which the outer patches represent the preferred habitat of the reservoir and spillover species and the middle patch represents the boundary region in which the species overlap.We note here that the overlap region was temporally-and spatially-dependent.The movement here is described in terms of number of visits per year to the boundary region and length of time spent there.Reproductive numbers for each patch were calculated, and it was determined that the greater number of interactions among species caused the reproductive number of the overlap patch to exceed the reproductive numbers of the patches representing preferred habitat, thereby causing a greater possibility of disease persistence.
In this study we consider three different geographical areas we refer to as patches.Each patch is described by a distinct transmission cycle between the hosts and vectors in the model.Patch 1 is defined by the T. gerstaeckeri-woodrat infection cycle in northern Mexico and southern Texas.Patch 3 is the south-eastern United States, including parts of Texas, Louisiana, Mississippi, Alabama, Georgia, and portions of the Carolinas and Florida panhandle and is defined by the T. sanguisuga-raccoon infection cycle, while patch 2 includes the south to southwest parts of Texas and a part of Coahuila, Mexico where we consider both species of vector and host to overlap.We assume the communication between patches and between cycles in patch 2 occurs through the movement of Triatoma vectors.T. sanguisuga migrates between patches 2 and 3, and between raccoons and woodrats in patch 2, while T. gerstaeckeri migrates between patches 1 and 2. We wish to investigate how increased vector migration affects prevalence in the overlap patch compared to the single-cycle patches.
In this study we develop and analyze an S-I metapopulation model with the aforementioned hosts and vectors, in which hosts may exhibit vertical transmission and vectors migrate between patches.We carry out standard analysis techniques, such as calculating R 0 for various sub-models of the larger metapopulation model as well as determining existence of endemic equilibria analytically (when tractable) in order to see effects of vector migration on R 0 .Finally, we perform numerical analysis on the full model to determine effects of migration on prevalence of T. cruzi.
2 Problem formulation and model

Problem statement and assumptions
To establish stable large-scale demographics for hosts and vectors, we assume that the growth for each species will be logistic, and neither the hosts nor vectors identified in the model exhibit disease-induced mortality.Studies have shown that the T. cruzi infecting raccoons and opossums in areas in the south-southeastern U.S. are not pathogenic, and do not appear to cause any symptoms of Chagas to the usual host [39,57].
Although triatomine bugs feed on many hosts in the wild, we consider woodrats and raccoons based on data that correlates their geographic location very closely with the vectors [40,9,41,24].In literature reviewed here, we have found the only host associated with T. gerstaeckeri is the southern plains woodrat [25,17].Therefore, we consider that the southern plains woodrat is the preferred host for the vector T. gerstaeckeri.
In this model there are several infection rates to be considered.We first distinguish between types of transmission routes.In this study, we consider both horizontal and vertical transmission routes, and we mention here that we group both stercorarian transmission and oral (host predation on infected vectors) transmission into a single horizontal transmission parameter.The infection rates may differ from vector to host and host to vector as well as by patch (geographical region).Biologically, the rate of infection between the hosts and vectors in each patch should be different.However, as seen in section 4 when calculating numerical estimates for the infection rate parameters (using a procedure called back-calculation), it is mathematically necessary to keep some of the rates the same.Thus, we assume that the rate of infection from raccoons to T. sanguisuga is the same in patches 2 and 3 and the rate of infection from T. gerstaeckeri to woodrats is the same in patches 1 and 2.
Literature suggests that T. sanguisuga will feed on other hosts besides raccoons, namely woodrats [17].Thus, in patch 2, T. sanguisuga feeds on raccoons and woodrats.Some proportion of vector-woodrat contacts are made with T. sanguisuga, thus we define q W as the proportion of vector-woodrat contacts made with T. sanguisuga, while 1 − q W is the proportion made with T. gerstaeckeri.Furthermore, since T. sanguisuga feeds on both hosts, we define q S as the proportion of T. sanguisuga contacts made with raccoons, while 1 − q S is the proportion made with woodrats.We denote the per vector infection rate from T. sanguisuga to woodrats as β SW and the per vector infection rate from woodrats to T. sanguisuga to be β W S .
T. cruzi has been confirmed to be transmitted vertically among mice in laboratory conditions [37].However, there is limited data on vertical transmission of T. cruzi in sylvatic hosts.We will assume that T. cruzi can be transmitted vertically in raccoons and woodrats.
Infection contact rates could be limited by the host or vector population.For our model, we assume that the limiting factor for infection from host to vector will be the vector population.The hosts are plentiful enough for vectors to feed as much as desired.Therefore, the vector population density will be the driving force in determining the infection rate, and the contact process saturates more quickly in the vector population than in the hosts.Thus, we consider the infection term from host to vector in the model to be based on several factors.Using similar derivation as in [27], we define the per-vector biting rate as z (in units of contacts per vector per time), and thus the total vector-feeding contact rate as z • N v (with units of bites/time).Thus, to calculate the rate of new vector infections, we multiply the total vector-feeding contact rate by the proportion of contacts that involve uninfected vectors and infected hosts, multiplied by the proportion of contacts that result in an infection (π v ) (units of infected vectors/bite) obtaining To simplify, we write β h = π v z (in units of 1/time).We also need to describe the rate at which vectors infect hosts.Again, assuming that the vector to host infection will be limited by the vector population (vectors feeding as frequently as desired), we multiply the total vector-feeding contact rate by the proportion of contacts involving uninfected hosts and infected vectors, multiplied by the proportion of contacts resulting in an infection (with units infected hosts/bite), In a similar manner, we define β v = π h z.But, we note here that β v is not in units of 1/time, but rather infected hosts per vector per time.We will apply this assumption to the infection terms in the model using the appropriate vector and host subscripts.

The Model
The model presented here is an S-I model incorporating migration and vertical transmission.1).The migration rates represent outgoing rates which must be adjusted by the patch area ratios for incoming rates (see system (1)) For each species we will use the general logistic birth rate b(N ) = rN 1 − N K , where r represents the intrinsic population growth rate, N represents the total population density, and K is the carrying capacity of the population density.We will use this term for each vector and host population using the appropriate subscripts.The logistic birth rate gives a reasonable description of the drop-off in births as population density increases.With respect to the vector population, as intraspecific competition for hosts intensifies, females lay fewer eggs, and the eggs' hatch rate also decreases.We further mention that the hosts in the model exhibit vertical transmission.Therefore, we incorporate this into the model by defining the following functions: f (N, I) and g(N, I) represent the birth rates for the hosts exhibiting vertical transmission.If only a proportion p (0 < p < 1) of infected hosts transmit vertically, then g(N, I) = pI • r 1 − N K , where p is the proportion of the offspring of infected hosts born infected with T. cruzi.Thus f (N, I) = (S + (1 − p)I) • r 1 − N K .We will apply this assumption to both hosts in the model, applying the appropriate subscripts.The natural per host and per vector mortality rates are denoted by µ R , µ W , µ S , and µ G .We distinguish clearly between the logistic birth and linear mortality rate because mortality is assumed to be spread evenly among infected and uninfected vectors and hosts, while births are assumed to contribute only to the susceptible class, except for vertical transmission described previously.
In this model we will assume linear migration based on the idea of local dispersion.We note that there are several different migration rates considered.We assume that hosts and vectors move at different rates, and those rates differ by species, by infection status, and by direction of migration.The migration parameters in the model are denoted by b i or bi for vectors and a i or āi for hosts, designating difference in direction of migration.Each subscript, i, is used to designate the migrating species, R, S, G, or, W. For example, bG represents the migration rate of T. gerstaeckeri from patch 1 to 2, while b G is the rate of T. gerstaeckeri from patch 2 to 1.
Furthermore, we note here that each patch has a different area, and the migration rates derived in [14] are affected by the size of the patch from which migration originates.In the metapopulation model described here, each differential equation represents the change in population density over time, thus the size of patch must also affect the population density of the vectors and hosts in each patch.Since we are considering migration of vectors between patches, the differences in population densities must be accounted for.For example, in order to account for the differing patch sizes, we consider the equation for the absolute number of infected T. sanguisuga vectors in patch 3, where A 3 denotes the area (in m 2 ) of patch 3 and the state variables have units of density.Then, the absolute number of vectors in patch 3 is given by It follows that In the model, we write bS A2 A3 = bS and b S A3 A2 = bS .This notation will be used more generally for northward and southward migration rates, respectively.
Finally, we make the following assumptions about the demographic and migration rates.First, for a given species (whether host or vector), in order for the population not to go extinct in a given patch, we must have r > µ + b (where b is the total per capita migration rate out of the patch)-that is, combined mortality and emigration cannot exceed the maximum growth rate.This seems biologically reasonable since no local vector or host extinction has been reported.Second, we further assume that µ > b; this assumption can be justified by rewriting it as 1/b > 1/µ-that is, the average time before migrating exceeds the average lifetime, which means most individuals will not leave their patch of origin during their lifetimes.Since the patches in this study are large and individual dispersal occurs on a small scale, this is reasonable.Third and last, we make the assumption that the migration rates are low enough not to cause equilibrium population densities in any patch to exceed the region's carrying capacity.If this assumption is violated, and the carrying capacities are exceeded, the model equations would need to be adjusted to distribute the resulting negative logistic term (as a result of additional density-dependent deaths) proportionally between the susceptibles and infectives, since as discussed earlier the equations are written taking logistic terms to measure births only (generally supposed into the susceptible class).All of these assumptions can be viewed as upper bounds on the migration rates, but in practice these bounds do not interfere with the explorations at the heart of this investigation, since (as will be seen in Section 4) estimated migration rates are orders of magnitude smaller than the demographic rates.
We therefore derive the model given below and seen in Figure 1.

Analysis
Little information is known about the migration of the vectors.Thus, we will explore several hypotheses.We will first consider that infected vectors are the only species to migrate.In a study done by Añez and East [1] in 1984 on the effect of parasites on the behavior of the vector Rhodnius prolixus, it was shown that that the parasite Trypanosoma rangeli hindered the vector's ability to draw blood, thus causing the bug to bite 25 more times than an uninfected vector.Thus, differential behavior of vectors infected with T. cruzi may affect the transmission of the parasite as well as vector mobility as mentioned in [28].In this context, we will consider the possibility that infected vectors move in only one direction (towards more preferred climates) and the possibility that infected vectors move between patches at different rates for different directions.Furthermore, we consider that uninfected vectors also migrate, but at a rate proportional to that of infected vectors.We will not consider host migration to play a significant role.The hosts in the model are bound by habitat constraints, and thus by definition of the patches, we assume the hosts are not likely to cross patch boundaries.
In order to get a better understanding of the full model, several special cases will be considered.The main identifying characteristics of the model are vertical transmission, migration, and multiple hosts and vectors.

One patch, one host, one vector, no vertical transmission
We begin the analysis of (1) by studying the simple system with one host and one vector.By observing system (1), we see that when the migration terms a i = āi = b j = bj = 0, for i = R, W , j = S, G, the three patches decouple.In this scenario, patch 1 and patch 3 are identical in structure.Thus, analyzing patch 1 and 3 with migration terms set to 0 and p W = p R = 0, we are analyzing the simple one host-one vector system.
Because the vector and host populations always approach an equilibrium, we can consider the limiting system in which N h and N v have reached their positive equilibria, N * h and N * v , where As mentioned previously, for each species, we assume r > µ, which guarantees that all disease-free extinction equilibria are unstable.This assumption will be carried out in this and all of the models hereafter.In system (2), we define βv = . Results by Thieme [48,49] guarantee that the behavior of the full system is asymptotic to the limiting system which is given by This model has been well studied [45,7] and we give results here.The basic reproductive number, calculated using the next-generation matrix [52], is given by which represents the average number of secondary infections caused by an infected individual introduced into a susceptible population.Because of the vector-host dynamics, R 0 represents the geometric mean between the average number of secondary host infections caused by one vector, and the average number of vector infections caused by one host.
In the case when R 0 < 1, the disease will die out and the population will approach the disease free equilibrium.When R 0 > 1, the population will approach a unique endemic state,

Patch 2, no migration, no vertical transmission
We will further analyze the decoupled system (1) by considering the equations representing patch 2 alone with no vertical transmission or migration.Thus, we analyze system (1) with a i = āi = 0, for i = R, W , b j = bj = 0, for j = S, G, and p R = p W = 0.The quantities N R2 , N S2 , N W 2 , and N G2 are asymptotically constant to , and N * G2 , respectively, where Applying the results from Thieme [48,49], we can study the limiting system In analyzing system (4) we observe the disease-free equilibrium is (I * S2 , I * R2 , I * G2 , I * W 2 ) = (0, 0, 0, 0).We find the basic reproductive number R 0 for the system using the next generation matrix method [52].The work can be seen in Appendix A. We determine that where In the terms for R 0 , we see that f 1 represents the T. gerstaeckeri-woodrat transmission cycle, f 2 represents the T. sanguisuga-woodrat cycle, and f 3 represents the T. sanguisuga-raccoon cycle.We further observe that At this point, we wish to observe the importance of the overlap of the transmission cycles between T. sanguisuga and T. gerstaeckeri, and how this overlap affects the ability of the infection to be spread.If there is no overlap, i.e. no T. sanguisuga-woodrat cycle, then f 2 = 0, so then we may define R0 = max{f 1 , f 3 } < R 0 .Since R0 is always less than R 0 , we observe the effect of f 2 is to increase the value of R 0 .It is possible that f 1 < 1 and By investigating the equilibrium conditions, it can be shown that either one or three endemic equilibrium values exist when R 0 > 1.The computations can be seen in Appendix A.

Patches 1 and 2, 1 host 1 vector with vertical transmission and unidirectional migration of infected vectors
In dealing with patches 1 and 2, there are several cases to be considered.We will first consider the scenario with one host and one vector, with vertical transmission, and unidirectional migration of infected vectors.In this case, N R2 = N S2 = 0, q W = q S = 0, and we can apply Thieme's results [48,49] and rewrite I ′ W 1 and I ′ W 2 , passing N W 1 and N W 2 to their limiting values, N * W 1 and N * W 2 .The system therefore simplifies to Cherif et al. [11] studied a similar model for T. cruzi vector transmission dynamics involving two strains (one being more virulent).In their model, a proportion of vectors infected with the more virulent strain migrate to a region in which the less virulent strain is native.However, their model did not include vertical transmission as system (6) does.
We determine the disease-free equilibrium for this system (6) to be (N * G2 , 0, 0, N * G1 , 0, 0), where R 0 can be found via the next generation matrix.For system (6), then the second term of R 0 is larger of the two because bG > 0. In this case max(p W , βGβW To study possible endemic equilibria of system (6), we determine that After substituting these values into the equilibrium conditions for system (6) (seen in Appendix B), we determine existence of an endemic equilibrium when R 0 > 1.We further determine that precisely one endemic equilibria exist in patch 2 alone if and only if R 1 < 1 < R 2 , and in both patches if R 1 > 1 (the patch 2 only endemic equilibrium is unstable in this case).

Patches 1 and 2, 1 host and 1 vector with vertical transmission and unidirectional migration of uninfected and infected vectors
We may also consider the case in which uninfected vectors move at a reduced rate proportional to that of infected vectors.As mentioned previously, if infected vectors exhibit differential behavior causing them to migrate more than uninfected vectors, we would consider the effects of having uninfected vectors migrate as well, but at a reduced rate.We note that this scenario is as far as vertical transmission can be treated analytically regarding R 0 .
The vertical transmission terms in the model do not affect the complexity in computing endemic equilibria since µ is simply replaced with (1−p)µ in the infected host equations.Thus, in computing R 0 , the remaining cases will be done without vertical transmission.We let N R2 = N S2 = 0, q W = q S = 0, a W = āW = 0, and b G = γ bG , where 0 < γ < 1.Since N W 1 and N W 2 are asymptotically constant to we apply the results of Thieme to this system, passing N W 1 and N W 2 to their limiting values, N * W 1 and N * W 2 , and obtaining the following system, Here, we can determine the disease free equilibrium to be (N * G2 , 0, 0, N * G1 , 0, 0), where Observing the terms of N * G2 , we can see that the first term, K G2 1 − µG rG , essentially represents the natural demographic renewal for the population of vectors in patch 2, while the second term represents the population being brought from vectors in patch 1.
R 0 for the system is In general the form of R 0 for system (7) is the same as that of system ( 6), with a different disease free equilibrium for the vector population due to the unidirectional migration of all vectors.We see that the first term of R 0 represents the patch 1 dynamics, while the second term represents patch 2. Similar to the system in section 3.3, if then the second term of R 0 is larger.Determining endemic equilibria for this system is intractable analytically.After a numerical investigation using the parameters estimated in section 4, we determine precisely one unique endemic equilibrium exists when R 0 > 1.

Patch 1 and 2, 1 host 1 vector, no vertical transmission, bidirectional migration of infected vectors
Another scenario we treat in patches 1 and 2 is one host, one vector and bidirectional migration of infected vectors.
We have previously assumed that vectors may have a preferred direction so that our migration is unidirectional.However, we know that vectors will move in every direction (although one direction may be preferred over another); thus we consider bidirectional migration.In this system, p R = p W = 0, q W = q S = 0, a i = āi = 0 for i = R, W .In this system, N W 1 and N W 2 are asymptotically constant to Applying the results by Thieme, we obtain the following system in which N W 1 and N W 2 have reached their limiting values.
The reduced system becomes We determine the disease-free equilibrium for this system (8) to be of similar form to that of (6).R 0 for the system is given as follows: We observe that R 0 for system (8) is of similar form as the R 0 for system (4), with the exception of the migration terms.It is observed that max( We interpret g 1 as the basic reproductive number for patch 1 scaled by the proportion of infected vectors that stay in patch 1, and g 2 is the basic reproductive number for patch 2 scaled by the proportion of infected vectors staying in patch 2.
Determining existence of endemic equilibria is intractable analytically, but after a numerical investigation, we verify the hypothesis that in the case of bidirectional migration, only one endemic equilibria is possible if R 0 > 1.

Patches 1 and 2, 2 hosts 2 vectors, no vertical transmission, unidirectional migration of infected vectors
We end our exploration of two patches by considering the case with 2 patches, 2 hosts, 2 vectors, and unidirectional migration of infected vectors.Then p R = p W = 0, a i = āi = 0 for i = R, W and b S = bS = b G = 0.In this scenario, the quantities N S2 , N R2 , N W 2 , and N W 1 are asymptotically constant to N * S2 , N * R2 , N * W 2 , and N * W 1 , respectively, where The model is For this scenario, we determine the disease-free equilibrium to be (0, 0, N * G2 , 0, 0, N * G1 , 0, 0), where After calculating R 0 , the structure seen is similar to that of R 0 for system (4), and is given by where P and Q are the same expressions as those in (5).We determine existence of endemic equilibria in patch 2 alone if and only if R 1 < 1 and R 2 > 1, and in both patches if and only if R 1 > 1.The computations can be seen in Appendix C.

Patch 1, 2, and 3, 2 hosts 2 vectors, no vertical transmission, unidirectional migration of infected vectors
We finally extend our discussion to all three patches.We now analyze the system represented by patch 1, 2, and 3, with 2 hosts, 2 vectors and unidirectional migration of infected vectors.Thus, we consider p i = 0, In this scenario, the host populations, N R3 , N R2 , N W 2 , and N W 1 are asymptotically constant to the values The reduced system becomes In analysis of system (10), we determine the disease free equilibrium to be (N * S3 , 0, 0, N * S2 , 0, 0, N * G2 , 0, 0, N * G1 , 0, 0), where We determine R 0 for this system to be R 0 = max{R 1 , R 2 , R 3 }, where

R3
and Based on the form of R 0 , we would expect three different scenarios for existence of endemic equilibria.We expect existence of endemic equilibria in all three patches if and only if R 1 > 1, in patch 2 and 3 only if and only if R 2 > 1 and R 1 < 1, and in patch 3 only if and only if R 3 > 1, R 1 < 1, and R 2 < 1. Investigation of these scenarios can be seen in Appendix D, in which we are able to show existence of at least one endemic equilibrium for each of the scenarios mentioned above.

Synthesis
By analyzing many smaller, sub-models of the original system (1), we may make some generalizations regarding the behavior of the full model.We expect that the full system will exhibit classical threshold behavior regarding R 0 , in which we expect a unique endemic equilibrium for R 0 > 1.As stated previously, the form of R 0 in a vector-borne disease is a geometric mean between infections caused by hosts and infections caused by vectors.If the system considers multiple hosts and vectors, the form of R 0 will include separate terms for each transmission cycle considered in the model.As described mathematically in section 3.3 and discussed in [28], we see vertical transmission has an "almost additive" effect on the basic reproductive number.Vertical transmission affects equilibria by effectively rescaling host mortality by a factor of (1 − p), where p is the proportion of births to infected mothers in which vertical transmission occurs.In the case of unidirectional migration of infected vectors, we observe that R 0 consists of as many components as there are patches, and each component for R 0 contains parameters for only one patch.Also, as seen in section 3.7, multiple endemic equilibria are possible depending on the values of the patchwise reproductive numbers, R 1 , R 2 , R 3 .We further note that by examination of ( 7), uninfected vectors migrating in one direction does not complicate the form of R 0 .With bidirectional migration of infected vectors, we determine that the expression for R 0 involves contributions from all patches, rather than having a maximum of several components, which is to be expected since infection is moving in between patches.Based on this determination, we expect the full model to have one component for R 0 , due to bidirectional migration; thus, it will not take on the form of max{R 1 , R 2 , R 3 }.In this case only one endemic equilibrium is possible; either there is no infection in any patch, or infection persists in all patches because all patches are connected by migration of infected vectors.

General demographic parameters
We wish to investigate numerically the results of section 3, as well as investigate the behavior of the full model given by system (1).In order to do this, we will use biological information to estimate the parameters given in our model.Kribs-Zaleta [27] completed a thorough literature study to estimate demographic and T. cruzi infection related parameters regarding hosts and vectors in the United States.We will use the demographic quantities calculated in [27], given in Table 1.We mention here that in our model, each host species has a preferred habitat.That is, the preferred habitat for the raccoons is patch 3; thus we would expect a higher population density of raccoons in patch 3 compared to patch 2. Similarly, the woodrat preferred habitat is prickly pear cactus which predominates in patch 1, with a lower density in patch 2 due to the varying landscapes.Here, we will treat the parameters that differ for each patch.

Patch 1
Kribs-Zaleta [27] obtains woodrat densities for Texas based on several sources [8,43] which estimate woodrat population densities in counties in south and west Texas.These regions, especially the counties in west Texas, are similar to south Texas, dominated by shrub desert, including cactus and honey mesquite.He estimates the woodrat (equilibrium) population density in patch 1 to be 2300 woodrats/km 2 .Because each population is governed by logistic growth with linear per-capita mortality, the populations approach an equilibrium population density, N * , so that the carrying capacity can be back-calculated using N * = K 1 − µ r .Kribs-Zaleta estimates the total vector population density in patch 1 to be 31600 vectors/km 2 .Since we are assuming that the only vector in our model in patch 1 is T. gerstaeckeri, we use this as the density estimate.Using these figures, we obtain patch 1 population density and carrying capacity estimates found in Table 2.

Patch 2
Estimates for southern plains woodrat density in patch 2 vary by geographical location and study.The woodrat density is affected by landscape, climate, and available materials for den construction.Raun [42] determined a positive correlation for population density of woodrats and density of cactus, although he concluded that cactus is not absolutely necessary to support woodrats.Cactus is the preferred material for den construction and food, but woodrats will use other materials to construct dens if cactus is unavailable [50].In general, density of woodrat dens is closely associated with overhead cover.
In an 18 month study in Jim Wells Co., part of the Southern Texas Plains, Merkelz and Kerr [35] record a maximum density of 1.5 wr/ha (during spring 1998) using a 10 ha subplot of a 220 ha study site.Density was calculated by using the number of woodrats captured in the area during a single trapping season.Since they do not give any other density calculations, we use this data as part of our data collection.We note that they did not limit their density calculation to only areas with cactus growth so as to include open areas as part of the normal daily range of woodrats.
Conditt and Ribble [13] estimate a range of 1.6-5.8wr/ha (average 3.7 wr/ha) in Bexar Co located in South Texas.The study was done on a 10 ha area of land with 4 ha dominated by honey mesquite-brush and prickly pear cactus, while the remaining 6 ha dominated by riparian lowland forests.The density was calculated on the 4 ha site due to essentially no woodrats being found on the riparian forested area (cactus-free) of the study site.They mention that the low density may be due to lack of cactus in the region of study and lack of appropriate shelter sites.However, Raymond et al. [43] in 2003 calculate a much higher maximum density of 19.4 wr/ha (with an average of 15.1 wr/ha) in the same county (but a different study site).This study site had limited clumps of prickly pear, but was covered with thick brush and downed trees which served as nest sites.Thus, we see that cactus is not absolutely necessary to maintain a high population density, but rather an abundance of shelter sites.
Although we do not include Oklahoma in patch 2, due to the northern range of T. gerstaeckeri, we refer to a study in Harmon County, Oklahoma in the Mesquite Grass Plains region, in which the estimated woodrat  [51].This region is native to the southern plains woodrat, dominated by Mesquite and prickly-pear cactus.The population density estimate was determined by trapping at 104 of the known 1,129 woodrat houses in the 226 ha study plot, and the density of cactus was not taken into account for the woodrat population density estimate.
Raun [42] estimates a range of 14.8-31.4woodrats/ha (average 23 wr/ha) in San Patricio Co, part of a transitional region between the Southern Texas Plains and Western Gulf Coastal Plain.Thus, the vegetation in this region is diverse, with the major plant communities being Mesquite-Mixedgrass, Chaparral-Mixedgrass, Live Oak-Chaparral, and Prickly Pear-Short grass, with riparian forests along the rivers [19,6].The study site, 9 acres, was reduced to 7.3 acres to estimate the density to eliminate areas that did not support cactus growth.A summary of these results can be found in Table 3.
Although there is a broad range of density estimates, we recall that patch 2 is a region with a diverse landscape, including mesquite, cactus, and savanna regions with areas of tree and prairie grassland.Since some estimates were computed in cactus-free regions, it is important to include each estimate in our computation of the average woodrat density.Thus, we include all 5 estimates in computing the average woodrat density in patch 2, obtaining an average population density of 11.3 wr/ha (1130 wr/km 2 ) in patch 2.
There are relatively few papers regarding raccoon distribution in patch 2. In a 3-year study by Gehrt and Fritzell [19], they estimate an average density of 7.3 raccoons/km 2 in San Patricio Co in southeast Texas.Since this region is not dominated by forest (as mentioned above), we would expect a lower raccoon density than that estimated by Kribs-Zaleta [27] for raccoon density in southeast USA (including patch 3).Using the equilibrium population densities, N * = K 1 − µ r , we calculate the carrying capacities for each species in patch 2, found in Table 4.
In patch 2, the T. cruzi transmission cycles overlap by the association of T. sanguisuga in association with both woodrats and raccoons, while T. gerstaeckeri feeds only on woodrats.The T. sanguisuga move between the woodrat and raccoon populations regularly enough that we will consider the T. sanguisuga a single population.We would like to estimate q W , the proportion of vector-woodrat contacts in patch 2 that are with T. sanguisuga, while 1 − q W is the proportion of vector-woodrat contacts that are T. gerstaeckeri.Also, since T. sanguisuga is associated with both hosts, we must estimate q S , the proportion of T. sanguisuga-host contacts that are raccoons, with 1 − q S the proportion of T. sanguisuga-host contacts made with woodrats.
Eads et al. [17] found 390 vectors from a total of 58 woodrat dens.Of the 390 vectors, 226 were T. sanguisuga, 133 were T. gerstaeckeri, and 31 were T. neotomae.We note here that the proportion of T. neotomae is negligible; thus we will normalize so that the proportions of T. sanguisuga and T. gerstaeckeri sum to 1. Thus, 63% of the vectors found in association with woodrats were T. sanguisuga and 37% were T. gerstaeckeri.Pippin [40] determined that for 85 woodrat dens, of 229 nymph and adult vectors, 58% were T. sanguisuga, while 42% were T. gerstaeckeri.If we pool the data, we determine 61% of the vector-woodrat contacts are with T. sanguisuga, and if we use the weighted average by number of dens excavated, the percentage is 60%.Thus, we estimate q W to be 0.605.
To estimate q S , we will define is the patch 2 woodrat density.V R is the number of vectors per raccoon and V W is the number of vectors per woodrat (scaled by the proportion that are T. sanguisuga).V R can be estimated from the estimates given in [27].We mention here that there are other hosts in patch 3, so not all T. sanguisuga can be found with raccoons.Ideally, we would calculate V R directly to avoid biasing the estimates.Although raccoons are the preferred host of T. sanguisuga, the vector will feed on other hosts including opossums in patch 3. We will estimate the raccoon density equivalent of the opossum density.Based on literature reviews, we determine the population density of opossums to be 10.1/km 2 [27].Thus, we determine the density of opossums is 0.505 times the raccoon density.
We will divide the patch 3 T. sanguisuga population density, N * S3 , estimated in [27] by the raccoon equivalent total host density.Using the estimates for N * S3 and N * R3 from [27], we determine sanguisuga/raccoon.Since both vector species are associated with woodrats, . Kribs-Zaleta [27] estimates that the total Triatoma vector population in patch 2 is 31600 vectors/km 2 .estimated previously, N * W 2 = 1130 woodrats/km 2 .These averages result in an estimate of V W as 16.9 T. sanguisuga/woodrat.Based on this calculation, we arrive at an estimate of q S = 0.286.

Patch 3
We use estimates from [27] for the raccoon and T. sanguisuga density estimates in patch 3, given in Table 5.We note here that the T. sanguisuga population density is based on a single study done by Burkholder et al. [9] regarding population density of Triatoma vectors, which estimates Triatoma density in relation to woodrat nests.Although T. sanguisuga in patch 3 are found with raccoons, we use the same estimate as the total vector population density in patches 1 and 2, due to lack of relevant information on vector population density in patch 3.

Estimation of infection rate parameters
The model here includes 2 modes of host infection: vertical transmission and direct transmission due to biting and to oral transmission via vector consumption.The vertical transmission parameters can be estimated directly via literature.Kribs-Zaleta [27] estimates the vertical transmission proportion to be 0.01.Then p W = 0.01.We estimate p R = 0.1 due to the adaptation of T. cruzi type IV to vertical transmission in raccoons.To estimate the direct infection rate parameters, we utilize the technique outlined in [27] to back-calculate the infection rate parameters, by solving for β R , β S , β SW , β W S , β G , and β W using the equilibrium conditions for model (1) under the assumption that observed prevalence indicates endemic equilibrium.
Since migration is small compared to the demographic processes, we will estimate the infection rate parameters by patch using model (1), assuming no migration.In order to estimate the infection rate parameters, we combine the observed prevalence levels and known demographic parameters in the equilibrium conditions to back-calculate the infection rate parameters, β R , β S , β SW , β W S , β S2 , β W 2 , β G , and β W .
After a thorough literature search, Kribs-Zaleta calculates prevalence levels for each species in each patch.Prevalence levels for T. gerstaeckeri and the woodrat are given for Texas and levels for T. sanguisuga and raccoons are given for Texas and the southeast U.S. To translate these values to this model, we assume that the patch 1 and 2 prevalence levels for T. gerstaeckeri and southern plains woodrat are equivalent to the Texas estimates found in [28].The patch 2 and 3 prevalence levels for T. sanguisuga and the raccoon are the same as the Texas and southeast estimates, respectively.A summary of these values is given in Table 6.We note here that x * G1 =

Species
is the prevalence value for T. gerstaeckeri in patch 1.We utilize a similar notation for the other 3 species in the model.
Using the equilibrium conditions and substituting the population density estimates and prevalence values for patch 1, we obtain estimates for β W and β G , given in Table 7. Substituting the patch 3 population density estimates and prevalence values into the patch 3 equilibrium conditions, we calculate β R and β S , given in Table 7.
Due to the crossover of infection cycles in patch 2, back-calculating the infection rate parameters is complex.After substituting the patch 2 equilibrium population densities and prevalence levels,there are 4 equilibrium conditions remaining We note that there are six β values and four equations; thus the system is underdetermined.Consequently, we will assume that β R and β G have the same values in patch 2 as calculated in patches 3 and 1, respectively.After solving system (11), we determine values for β SW , β W S , β S2 , and β W 2 , given in Table 7.
After an investigation of the effects of the proportions q S and q W on the β W S and β SW with all of the other parameters held fixed as determined above, we determine that any value of q W greater than 0.3 will produce a value of β SW in between 0.11/yr and 0.15/yr, which mitigates any inaccuracy in our estimate of q W = 0.605.We also note that, mathematically, q S needs to be smaller than 0.411 for β W S to be positive, as seen in Figure 2.
Of the vector to host transmission terms, β S and β S2 (vector to raccoon infection rate) are close in value, while β G and β SW (vector to woodrat infection rate) are close in value.But, we note here that there is a factor Table 7: Stercorarian infection rate parameters of 1000 by which woodrats are being infected more rapidly than raccoons are infected by T. sanguisuga.We note this may be attributed partly to the fact that the population densities for woodrats are approximately 100 times as great as raccoon population densities.All of the host to vector transmission parameters are close in value, with the exception of β W S (woodrat infecting T. sanguisuga) which is an order of magnitude less than the raccoon to T. sanguisuga infection rates, which is to say that woodrats are infecting T. sanguisuga at a lower rate than raccoons infecting T. sanguisuga.

Trends in migration rate effects
To investigate the effects of migration numerically, we examine prevalence of T. cruzi as a function of vector migration.In order to get a clear picture, we first look at unidirectional migration of infected vectors.We then investigate unidirectional migration of all vectors as this is the upper bound of the possible unidirectional migration scenarios for vector migration.We then consider bidirectional migration of infected and all vectors.
In each scenario, we consider the effects on T. cruzi prevalence as migration increases.To see the effects of the increase, we consider the effects of factor, k, multiplied by each migration rate.We calculate the prevalence for each patch as k increases from 0 to 10, where k = 0 represents the scenario with no migration and k increasing greater than 1 represents the migration rate increasing past the calculated rate from [14].For each graph given, the host infected prevalence is represented by the lighter colored graph and the vector infected prevalence is represented by the darker colored graph.
In the case of northward migration of infected vectors (Figure 3), we observe a decrease of prevalence in patch 1 since this patch is losing infected vectors and not gaining any, and a slight increase of prevalence in patch 3 as the migration rate increases.The prevalence of T. sanguisuga in patch 2 is much lower than that in patch 3, so although patch 3 is gaining infected vectors, the increase in prevalence at equilibrium is minimal.In this scenario, patch 2 is gaining infected vectors from patch 1(T.gerstaeckeri) and losing infected vectors to patch 3(T.sanguisuga).As seen in Figure 3, the patch 2 prevalence decreases for northward migration.To get a better understanding of why the prevalence decreases in patch 2, we look at the graphs of the prevalence for each vector species in patch 2, seen in Figure 4.It can be observed that the T. sanguisuga infected prevalence decreases by approximately 35% for high migration rates, while the T. gerstaeckeri prevalence increases by approximately 25%.Furthermore, we see a constant decrease in the woodrat prevalence (although T. gerstackeri prevalence increases).We note that in patch 2 N S2 is approximately 50% greater than N G2 .Thus, the export of infected T. sanguisuga patch 2 dominates the import of infected T. gerstaeckeri from patch 1, causing an overall prevalence decrease in patch 2 in the case of northward migration of infected vectors.When we consider northward migration of all vectors, the change in prevalence is undetectable by viewing the graph.Patch 1 prevalence is still reduced since it is losing both infected and uninfected vectors, as seen in Figure 5.One might expect that the decrease should be less when all vectors are migrating, but after a numerical investigation of prevalence, we see that the decrease is slightly greater when all vectors are migrating.Also, we would expect the patch 3 prevalence to increase since this patch is gaining both infected and uninfected vectors from patch 2 with a lower prevalence; however the prevalence decreases slightly.If only infected vectors migrate northward, then R 1 decreases as bG increases, causing prevalence in patch 1 to decrease.When all vectors are migrating northward, R 1 actually decreases more since the vector-host ratio is also decreasing, as the migration rate increases.Thus, the decrease in R 1 is amplified when all vectors are migrating northward.
As we observe southward migration of infected vectors, we again view the prevalence levels for each patch, in Figure 6.Patch 1 is gaining infected vectors from patch 2; thus we see an increase in prevalence.In patch 3, we see a slight decrease (< 1%) in prevalence due to this patch losing infected vectors.In contrast to northward migration, we see a rise in prevalence for patch 2 in the case of southward migration.To understand why prevalence increases even though this patch is losing infected vectors to patch 1, we graph prevalence levels for each species in patch 2 as a function of migration, seen in Figure 7.We observe that prevalence in T. sanguisuga increases, as expected since this vector population is gaining infected vectors from patch 3.An interesting observation is that the prevalence in T. gerstaeckeri rises for small migration rates (0 < k < 1).Then, for higher migration rates, the prevalence for T. gerstaeckeri in patch 2 decreases.As observed in Figure 7, the woodrat T. cruzi prevalence increases in patch 2 for southward migration.Thus, for 0 < k < 1, the sharp rise in prevalence in T. sanguisuga (11% increase) in patch 2 for small migration rates may be enough to increase the prevalence in T. gerstaeckeri through the woodrats.
In the case of southward migration of all vectors, the patch 3 prevalence decreases more than when only infected vectors are moving southward.This can be attributed to the decrease in vector-host ratio as the migration rate increases, thereby causing R 3 to decrease more than if only infected vectors are migrating.Patch 1 prevalence increases more than if only infected vectors are migrating southward, due to the slight increase in vector-host ratio.The patch 2 prevalence increases more than if only infected vectors are migrating southward, which is primarily due to the increase in vector-host ratio causing R 2 to increase more than if only infected vectors are migrating, as seen in Figure 8.
We also wish to investigate bidirectional migration of vectors.We first consider bidirectional migration of infected vectors only.As seen in Figure 9, prevalence decreases in patches 1 and 3, but increases in patch 2. To better understand why prevalence increases in patch 2, we observe the prevalence graphs for patch 2 only.We observe that both T. sanguisuga and T. gerstaeckeri prevalence increases in patch 2 (Figure 10).The T. sanguisuga population in patch 2 is initially at a lower prevalence than patch 3, and since the T. sanguisuga population is gaining vectors from a patch with a higher prevalence, the prevalence increases for this species in patch 2 as seen in Figure 10.The T. gerstaeckeri prevalence also increases, although this population is gaining and losing vectors at the same prevalence.We can most likely attribute the increase in prevalence for T. gerstaeckeri in patch 2 due to the increase in woodrat prevalence (seen in Figure 10) through the T. sanguisuga prevalence increase.
In the case of bidirectional migration of all vectors, the behavior is similar.The patch 1 and 3 prevalence decreases slightly more than in the case of only infected vectors migration, due to the decrease of the vector-host

Calculation of the migration rate
The framework established in [14] describes vector dispersal in terms of three properties: dispersal distance, preferred direction of dispersal, and degree of preference for a particular direction.Since we do not have clear data on a preference for a direction for T. sanguisuga and T. gerstaeckeri, we will consider the simplest case in which the vectors have no preference for direction of dispersal.Here we give the vector migration rates calculated in [14] assuming no preference for direction.We adjust these rates to take into account the area ratios and note that the rates are equal for all vectors (infected and susceptible) in each species.We may also consider that vectors have a preference for direction of migration.Although we do not have clear evidence for the vectors in the model migrating with any clear trend in direction, we consider the possibility of vector migration in a particular direction.In a study on climate change related to Chagas disease distribution, Curto de Casas concludes that higher temperatures may extend the geographical range of the sylvatic vectors of T. cruzi [15].Thus, we may consider vector migration with a northward preference for direction.We note here that when northward is described in section 4.3.1, it is referring to one-directional migration from patches 1 to 2 and 2 to 3, not the actual geographical direction north.In this section we are referring to the actual direction north and migration is bidirectional between patches.
The framework for deriving migration rates with a preference for direction is modeled using a sequence of nested ellipses, in which each ring represents a certain range of dispersal distances.We also assume that vectors have a degree of preference for a direction, which represents the eccentricity of each ellipse.Based on results from [14], we give results for T. cruzi prevalence assuming a northward preference of direction (with moderate degree of preference, e = 0.5) in Table 9.The decrease in prevalence in patch 1 and increase in prevalence for T. gerstaeckeri in patch 2 are consistent with what we expect assuming a northward preference for direction.We note that there is approximately a 1% decrease in prevalence for T. sanguisuga and raccoons in patch 2 if preferred direction is northward, when compared to migration with no preference for direction.There is a more than 3% increase in T. gerstaeckeri prevalence in patch 2 when compared to no migration.Although T. gerstaeckeri feeds only on woodrats, the increase in prevalence for T. gerstaeckeri is not enough to cause a higher increase in woodrat prevalence; thus we only observe only a slight (<1%) increase in prevalence for woodrats in patch 2 when the preferred direction of migration is northward.

Conclusions
Due to the complexity of system (1), several sub-models were analyzed in order to make conclusions regarding the behavior of the full model.The entire system appears to exhibit classical threshold behavior regarding R 0 , and existence of a unique endemic equilibrium when R 0 > 1.In the case of one-directional migration of vectors, R 0 consists of as many components as there are patches with R 0 being the largest value of the patch-specific R 0 values, and multiple endemic equilibria are possible depending on the values of the patch-specific reproductive numbers values.If migration is bidirectional, R 0 will consist of only one component for containing parameters from all patches, and only one endemic equilibrium is possible; thus either the whole system reaches a disease-free state or infection persists in all patches.These results are similar to the multi-patch model results in [3], where it is shown that if patch a is at endemic equilibrium, then the disease is at endemic equilibrium in each patch accessible to patch a. Similar results are observed in the two-patch model in [2], where it is shown that if R 0 < 1, the disease does not persist in either population (patch), and if R 0 > 1, the disease persists in both populations.
In section 4, we obtained several different sets of results with respect to the effects of vector migration on the prevalence of T. cruzi in the geographical region from which we built our model.We examined the effects of one-directional migration, bidirectional migration (at symmetric rates), and bidirectional migration using rates derived from the framework in [14].
In the case of one-directional migration, the prevalence in the outer patches (patches 1 and 3) varies as expected (e.g.patch 1 prevalence decreases for northward migration only).However, the patch 2 dynamics are different for northward vs. southward migration.In the case of northward migration only, the patch 2 prevalence decreases while for southward migration, the prevalence increases.After closer observation of prevalence for each species in patch 2, it is determined that the T. sanguisuga migration has a greater effect on the patch 2 dynamics than the migration of T. gerstaeckeri primarily due to the fact that the difference in population density and prevalence in patches 2 and 3 among T. sanguisuga is greater than the difference in population density and prevalence in T. gerstaeckeri in patches 1 and 2. If vectors are migrating northward only, the export of T. sanguisuga dominates the import of T. gerstaeckeri thereby causing an overall decrease in patch 2 prevalence; we note that the T. sanguisuga population density in patch 2 is approximately 50% greater than the T. gerstaeckeri population density.For southward migration rates, the import of T. sanguisuga dominates the export of T. gerstaeckeri into patch 1.In fact, for lower migration rates, the prevalence for T. gerstaeckeri increases slightly due to the sharp increase in T. sanguisuga through the connection with the woodrat host.The connection between the vectors is the infection cycle with the woodrat host.Thus, the increase in prevalence for T. gerstaeckeri for lower southward migration rates is due to the increase in the prevalence for woodrats through the infection cycle with T. sanguisuga.
If migration is bidirectional, the patch 2 prevalence increases.After analysis of one-directional migration, the reason for the patch 2 prevalence increase is more apparent.We note here that the prevalence for T. sanguisuga is initially at a lower prevalence in patch 2 than in patch 3, so the T. sanguisuga in patch 2 is gaining vectors from a patch with a higher prevalence, thus the increase in T. sanguisuga prevalence.The prevalence for T. gerstaeckeri increases in patch 2 in the case of bidirectional migration due to the increase in woodrat prevalence (again through the infection cycle with T. sanguisuga).We note here that these trends are observed even if migration rates are considered symmetric (independent of patch size).
Because of the differences in patch sizes (especially the large size of patch 3 compared to patches 1 and 2), the vector density in each patch is affected differently by migration.For example, the change in patch 2 vector density will be greater than the change in patch 3 vector density for bidirectional migration.The patch 3 vector density is minimally affected by migration due to the large patch size.
Based upon these results we can conclude that infection dynamics in patch 2 are sensitive to migration, but primarily driven by the T. sanguisuga population.Since the same trends in prevalence change are observed for migration independent of patch size, we should investigate the distinctive transmission characteristics between host and vector.Thus, we consider differences in the T. sanguisuga and T. gerstaeckeri populations.One major difference in the vector populations is the high prevalence for T. sanguisuga in patch 3 compared to patch 2 and the prevalence for T. gerstaeckeri.If we assume no difference in prevalence levels for T. sanguisuga from patch 3 to 2, we can note that for southward migration, the prevalence in patch 2 still increases for small migration rates, but at a much slower rate, while for northward migration, the prevalence decreases but at a much slower rate than if the prevalence in patch 3 is at its current estimated level.These results give a different picture for the patch 2 dynamics, which allow us to see a dampening effect of the T. sanguisuga migration if the prevalence levels for the T. sanguisuga populations in patches 2 and 3 are the same.Thus, we can attribute the majority of the patch 2 dynamics when migration is considered to the higher prevalence in T. sanguisuga in patch 3.
This aim of this study is to determine the effect of vector migration on T. cruzi transmission in the prominent sylvatic cycles ranging from northern Mexico to the southeastern United States.The primary effect of migration is to increase prevalence in the overlap patch where the prevalence is initially at a lower level than the outer (singlecycle) patches.The dominant force is the connection to the large raccoon-T.sanguisuga in the southeastern United States, which is affected little on such a large scale, but which can affect dynamics strongly in the overlap region with woodrats and T. gerstaeckeri.
As in every study, certain limitations must be noted.Due to the very recent awareness of the need to study T. cruzi in the United States, there have been very few studies on the demography of vectors native to the U.S., especially regarding vector population density.To date, there have been virtually no studies on the U.S. Triatoma vectors' dispersal or migration capabilities.Thus, with more studies, we may be able to more accurately describe the dispersal capabilities (especially regarding distance and frequency of vector flights).We acknowledge the need for more experimental or field studies with heavy consideration on the Triatoma vectors native to the United States, especially T. gerstaeckeri and T. sanguisuga which we consider to be the primary vectors in the southeast.
Because of the differing patch sizes, a natural question arises to consider effects of migration for higher spatial resolution where all patches are of uniform size.If the geographical region is broken into smaller, samesize patches, we wish to examine to what extent these results would change.Future work already in progress uses cellular automata to see how migration among smaller patches affects the spread of T. cruzi across a larger geographic region as well as determine a measure for speed of invasion.Although no data exist on large-scale vector migration patterns in this region, this work describes the geographical spread of T. cruzi under different migration scenarios.

A Patch 2 with p = 0 and no migration
We determine the reproductive number, R 0 , for section 3.2, patch 2 with no vertical transmission or migration, using the next generation matrix method [52].
Based on the next generation matrix method, we rewrite system (4), dX dt = f (X) in terms of two vectors: dX dt = F 0 − V 0 .F 0 represents the terms generating new infections, while V 0 consists of the remaining terms.
After computing the derivatives of F 0 and V 0 , we obtain The dominant eigenvalue of where It remains to be shown that P 2 − 4Q is real.Thus, we must show that P 2 − 4Q > 0, as follows: To determine the existence of endemic equilibria, we determine the equilibrium conditions to be where . To simplify, we divide the 1st and 3rd equations by the nonzero values N * S2 and N * G2 respectively and make the substitutions, βG = We solve the 2nd and 3rd equations for x * R2 and x * G2 : x We then solve the 1st equation in order to isolate x * W 2 .We obtain Substituting the result for x * R2 from above, we obtain Multiplying through on both sides by the denominator and combining like terms results in We then divide everything by µ R and µ S to obtain, To simplify, we make the substitutions A = qS βR µS , B = βS2 µR , and C = (1−qS )βW S µS and solve for x * W 2 to obtain We finally solve the 4th equation; substituting x * G2 , we obtain βG Multiplying through by the denominator on both sides, we obtain and dividing all terms by µ W and µ G to get We make the substitutions D = βSW µW , E = βW µG , and F = βG µW and rewrite the previous equation as We can now determine a polynomial in x * S2 , say φ(x * S2 ).The resulting polynomial is of 4th degree, once the disease free equilibrium has been divided out, with constant term φ(0) = C(P − Q − 1), with P, Q as defined in (5), which is positive for R 0 > 1 by the following result R 2 0 > 1 1 2 P + P 2 − 4Q > 1 Therefore, since P 2 − 4Q is real, we have R 0 > 1 ⇔ P − Q > 1.Furthermore, it is observed that φ(1) < 0 which implies that φ(x * S2 ) crosses the x-axis 1 or 3 times between 0 and 1.By inspection of the form of x * R2 and x * G2 , we see that x * R2 is in (0,1) if x * S2 is in (0,1) and x * G2 is in (0,1) if x * W 2 is.Alternatively, we can use the equations ( 14) and ( 16) to obtain a polynomial in x * W 2 and apply the same technique.Thus, there are either 1 or 3 endemic equilibria for this system when R 0 > 1.Thus, if (x * S2 , x * W 2 ) ∈ [0, 1] × [0, 1], then the remaining values, x * R2 and x * G2 are also in [0, 1].
B Patches 1 and 2, 1 host 1 vector, with VT, unidirectional migration of infected vectors To determine the basic reproductive number for section 3.3, we apply the next generation matrix method [52] to (6).We determine the relevant matrices to be After computing F 1 V −1 1 , we obtain the dominant eigenvalue as To study the existence of endemic equilibria, we first determine The remaining equilibrium conditions are given by It can be seen that if 0 ≤ x * W 1 < 1, then 0 ≤ x * G1 < 1.To show 0 < x * W 1 < 1, we substitute this expression for x * G1 into the 4th equation in the equilibrium conditions and expand to obtain a cubic polynomial in x * W 1 .We can factor out the disease free equilibrium (in patch 1), in which x * W 1 = 0 leads us to the simple one host one vector patch 2 only endemic equilibrium.Thus, if R 1 < 1 < R 2 , we have existence of precisely one endemic equilibrium.Otherwise, for x * W 1 > 0, our resulting quadratic polynomial is ζ(x * W 1 ) = a 2 x * W 1 2 + a 1 x * W 1 + a 0 = 0 with ] We first observe that the constant term, ζ(0) = a 0 is positive if and only if R 1 > 1 by the result below: Furthermore ζ(1) < 0. Thus, we can conclude that ζ(x * W 1 ) crosses the x-axis precisely once between 0 and 1.In order to show that 0 ≤ x * G2 < 1 and 0 ≤ x * W 2 < 1, we solve the first equilibrium condition for x * G2 in terms of x * W 2 and x * G1 , We note that if x * W 2 is in (0, 1), < 1 and µ G > bG , then 0 ≤ x * G2 < 1.To show 0 ≤ x * W 2 < 1, we substitute the expression for x * G2 into the 2nd equilibrium condition from (17): which are identical to the equilibrium conditions (13) for patch 2 alone, in which we have determined that if R 2 > 1, either 1 or 3 endemic equilibrium exists, in addition to the disease free equilibrium (which always exists).If instead, R 1 > 1, we assume that x * G1 is the positive equilibrium given in Appendix B for system (6) with p W = 0. Thus, in determining possible endemic equilibrium values for this system, we may solve the 3rd equation of (19) to obtain We substitute x * G2 into the 4th equation obtaining and multiplying through by the denominator we obtain and dividing every term by µ W , µ G , the expression becomes Applying the same substitutions in Appendix B, where the first two equations of ( 19) are identical to the first two equations of (13) in Appendix A and letting X = Finally, we use (23) and (15) from Appendix A to obtain a polynomial in x * S2 , We note that is positive if and only if x * S2 > AB − 1 AB + A , in which case AB > 1 for the expression to be positive (biologically significant).

Figure 2 :
Figure 2: Correspondence of q S and q W with β W S , with units 1/yr, and β SW , with units infected hosts/vector/yr Parameter Value Units β S 0.00025 hosts vector•yr

Figure 3 :
Figure 3: Equilibrium T. cruzi prevalence vs. vector migration rate multiplier for infected vectors migrating north only; dark curve represents vectors, light curve represents hosts

Figure 4 :
Figure 4: Equilibrium T. cruzi prevalence vs. vector migration rate multiplier for infected vectors moving north only (patch 2)

Figure 5 :Figure 6 :
Figure 5: Equilibrium T. cruzi prevalence vs. vector migration rate multiplier for all vectors moving north only; dark curve represents vectors, light curve represents hosts

Figure 7 :
Figure 7: Equilibrium T. cruzi prevalence vs. vector migration rate multiplier for infected vectors moving south only (patch 2)

Figure 8 :Figure 9 :
Figure 8: Equilibrium T. cruzi prevalence vs. vector migration rate multiplier for all vectors moving south only; dark curve represents vectors, light curve represents hosts

Table 4 :
Patch 2 density estimates density was 13 wr/ha where N * R2 is the patch 2 raccoon density and N *

Table 8 :
Migration rates for no preferred direction (units in 1/year)

Table 9 :
Equilibrium prevalence levels for species based on migration (northward preference for direction)