Phylogeography and Genetic Variation of Triatoma dimidiata, the Main Chagas Disease Vector in Central America, and Its Position within the Genus Triatoma

Background Among Chagas disease triatomine vectors, the largest genus, Triatoma, includes species of high public health interest. Triatoma dimidiata, the main vector throughout Central America and up to Ecuador, presents extensive phenotypic, genotypic, and behavioral diversity in sylvatic, peridomestic and domestic habitats, and non-domiciliated populations acting as reinfestation sources. DNA sequence analyses, phylogenetic reconstruction methods, and genetic variation approaches are combined to investigate the haplotype profiling, genetic polymorphism, phylogeography, and evolutionary trends of T. dimidiata and its closest relatives within Triatoma. This is the largest interpopulational analysis performed on a triatomine species so far. Methodology and Findings Triatomines from Mexico, Guatemala, Honduras, Nicaragua, Panama, Cuba, Colombia, Ecuador, and Brazil were used. Triatoma dimidiata populations follow different evolutionary divergences in which geographical isolation appears to have had an important influence. A southern Mexican–northern Guatemalan ancestral form gave rise to two main clades. One clade remained confined to the Yucatan peninsula and northern parts of Chiapas State, Guatemala, and Honduras, with extant descendants deserving specific status. Within the second clade, extant subspecies diversity was shaped by adaptive radiation derived from Guatemalan ancestral populations. Central American populations correspond to subspecies T. d. dimidiata. A southern spread into Panama and Colombia gave the T. d. capitata forms, and a northwestern spread rising from Guatemala into Mexico gave the T. d. maculipennis forms. Triatoma hegneri appears as a subspecific insular form. Conclusions The comparison with very numerous Triatoma species allows us to reach highly supported conclusions not only about T. dimidiata, but also on different, important Triatoma species groupings and their evolution. The very large intraspecific genetic variability found in T. dimidiata sensu lato has never been detected in a triatomine species before. The distinction between the five different taxa furnishes a new frame for future analyses of the different vector transmission capacities and epidemiological characteristics of Chagas disease. Results indicate that T. dimidiata will offer problems for control, although dwelling insecticide spraying might be successful against introduced populations in Ecuador.


Introduction
American trypanosomiasis or Chagas disease is widespread in Latin America from Mexico to Chile and southern Argentina. Although present estimates of 10 to 12 million people infected with the haemoflagellate protozoan species Trypanosoma cruzi represent 6-8 million fewer cases than those reported in the 1980s [1], it remains one of the most serious parasitic diseases of the Americas for its social and economic impact [2]. Although it can also be transmitted by blood transfusion or across the placenta from infected mothers, most human contamination is attributed to insect vectors in poor rural or periurban areas of Central and South America [1].
Chagas disease vectors are haematophagous reduviid (Hemiptera: Heteroptera) insects belonging to the subfamily Triatominae. Species of Triatominae are usually grouped into 17 genera forming five tribes, although other arrangements have been proposed. Of these, Alberproseniini, Bolboderini, Cavernicolini and Rhodniini are considered monophyletic, whereas Triatomini is considered polyphyletic [3]. Among the latter, most of the species (over 70) are included in the genus Triatoma, among which two main clades appear in ribosomal DNA (rDNA) sequence phylogenies, corresponding to species of North and Central America and species of South America separated prior to the closing of the isthmus of Panama about 3 million years ago [4][5][6]. Moreover, Triatoma species are distributed in three main groupings: the Rubrofasciata group (mainly North American and Old World species), the Phyllosoma group (mainly Mesoamerican and Caribbean), and the Infestans group (mainly South American), each including different complexes and subcomplexes in a classification which is progressively updated according to new genetic and morphometric data [7].
A priori, all of the over 130 species currently recognized within Triatominae seem capable of transmitting T. cruzi. Among the species of greatest epidemiological significance as domestic vectors, three belong to the genus Triatoma: T. infestans and T. brasiliensis from South America, and T. dimidiata, distributed in Meso-and Central America from Mexico down to Colombia, Venezuela, Ecuador and northern Peru [3].
Triatoma dimidiata can be found in sylvatic, peridomestic and domestic habitats. Non-domiciliated populations may act as reinfestation sources and become involved in the transmission of the parasite to humans [8,9]. This species includes morphologically variable populations [10,11]. A molecular comparison of Triatominae, including many Central American species of the Phyllosoma complex by means of rDNA second internal transcribed spacer (ITS-2) sequences demonstrated an unusual intraspecific sequence variability in a few T. dimidiata populations studied. This study even revealed differences consistent with a specific status for populations from the Yucatan peninsula, Mexico [4][5][6], thus opening a debate. A large number of recent, multidisciplinary studies using RAPD-PCR, genital structures, morphometrics of head characters, and antennal phenotypes have shown that variation within this species seems much greater than previously considered [8,[12][13][14][15][16]. Morphometric and cuticular hydrocarbon analyses suggest that a sylvatic population from Lanquin, Guatemala, is undergoing a speciation process [13,17]. Chromosomal variation and genome size suggest that T. dimidiata may represent a complex of cryptic species (i.e. morphologically indistinguishable, yet reproductively isolated taxa) [18].
The aim of the present work is to analyze the intraspecific variability, haplotype profiling, phylogeography and genetic polymorphism of populations of the species T. dimidiata, to get a new framework able to facilitate the future understanding of the diferring peculiarities of this crucial vector species throughout its broad geographical distribution. This may also help in understanding the related differences in characteristics of Chagas disease transmission and epidemiology, as well as in responses to control initiatives in the countries concerned. After a deep analysis, it was considered that the most convenient approach would be obtained by using an appropriate marker able to furnish significant information about evolutionary trends of variation on which to construct the new baseline. This new baseline should be, whenever possible, of sufficient weight as to allow its conclusions to be reflected at systematic-taxonomic level.
For this purpose, the rDNA was preferred over mitochondrial DNA (mtDNA) because of its mendelian inheritance, evolutionary rates and overall recognized usefulness in systematics in all metazoan organism groups because of including sequences which allow to distinguish between species and between subspecies units. The better fitting of rDNA for molecular systematics has already been emphasized in large reviews on rDNA/mtDNA marker comparisons in insects [19]. Ribosomal DNA includes excellent genetic markers, because (i) the rDNA operon is tandemly repeated and present in sufficiently high quantities among the genome of an individual thus facilitating sequencing procedures; (ii) the different genes and spacers of the rDNA follow a concerted evolution which, with sufficient time, effectively homologizes the many copies of nuclear rDNA within a genome [20]; this gives rise to a uniformity of their sequences within all individuals of a population and becomes extremely useful from an applied point of view, because it is sufficient to obtain the sequence of only one individual to characterize the local population it belongs to, that is, all other individuals of that population will present the same sequence; (iii) the usefulness of rDNA genes and spacers as genetic markers at different evolutionary levels have already been verified on a large number of very different eukaryotic organism groups including insects, and consequently extensive knowledge on the different rDNA fragments is available [21]. rDNA sequence comparisons offer valuable information about the evolutionary events in triatomine lineages and, by deducing the routes of

Author Summary
Chagas disease is a serious parasitic disease of Latin America. Human contamination in poor rural or periurban areas is mainly attributed to haematophagous triatomine insects. Triatoma includes important vector species, as T. dimidiata in Central and Meso-America. DNA sequences, phylogenetic methods and genetic variation analyses are combined in a large interpopulational approach to investigate T. dimidiata and its closest relatives within Triatoma. The phylogeography of Triatoma indicates two colonization lineages northward and southward of the Panama isthmus during ancient periods, with T. dimidiata presenting a large genetic variability related to evolutionary divergences from a Mexican-Guatemalan origin. spreading of triatomine populations, they may also shed light on the ability of different species to colonize new areas [5].
Within rDNA, ITS-2 was selected as marker because of its wellknown usefulness at species and subspecies levels, including the differentiation of taxa within problematic groups, as is the case of those comprising cryptic or sibling species of other insect groups [22][23][24]. Moreover, the sequences of the ITS-2 have already proved to be a useful tool in the analysis of species, subspecies, hybrids and populations, and for inferring phylogenetic relationships in Triatominae in general [4,5,6,25,26].
In order to be able to assess the ITS-2 evolutionary processes followed by T. dimidiata populations, the ITS-2 sequences of many members of the Phyllosoma, Rubrofasciata and Infestans groups were obtained and analyzed. For this purpose, a large number of rDNA ITS-2 sequences of Triatoma species from numerous geographic origins in Mexico, Guatemala, Honduras, Nicaragua, Panama, Cuba, Colombia, Ecuador, and Brazil was studied. Thus, the nucleotide divergence limits between taxa within the lineage of the genus Triatoma could be established. The present study on T. dimidiata is the largest interpopulational analysis performed on a triatomine species so far.

Triatomine materials
A total of 165 triatomine specimens representing 13 Triatoma species of the Phyllosoma, Rubrofasciata and Infestans groups, among which 137 specimens representing T. dimidiata from 64 different geographic origins, were used for sequencing, genetic variation and phylogenetic analyses (Table 1; Figure 1). The systematic classification recently proposed for the genus Triatoma [7] is used here throughout.

Sequencing of rDNA ITS-2
For DNA extraction, one or two legs fixed in ethanol 70% from each specimen were used and processed individually, as previously described [5,27]. Total DNA was isolated by standard techniques [28] and stored at 220uC until use. The complete ITS-2 fragment was PCR amplified using 4-6 ml of genomic DNA for each 50 ml reaction. Amplifications were generated in a Peltier thermal cycler (MJ Research, Watertown, MA, USA), by 30 cycles of 30 sec at 94uC, 30 sec at 50uC and 1 min at 72uC, preceded by 30 sec at 94uC and followed by 7 min at 72uC. PCR products were purified with Ultra Clean TM PCR Clean-up DNA Purification System (MoBio, Solana Beach, CA, USA) according to the manufacturer's protocol and resuspended in 50 ml of 10 mM TE buffer (10 mM Tris-HCl, 1 mM EDTA, pH 7.6). Sequencing was performed on both strands by the dideoxy chain-termination method, and with the Taq dye-terminator chemistry kit for ABI 3730 and ABI 3700 capillary system (Perkin Elmer, Foster City, CA, USA), using the same amplification PCR primers [6].

Triatomine haplotype code nomenclature
The haplotype (H) terminology used in the present paper follows the nomenclature for composite haplotyping (CH) recently proposed [25]. Accordingly, ITS-2 haplotypes (H) are noted by numbers (Table 1).

Phylogenetic inference
Phylogenies were inferred by maximum-likelihood (ML) using PAUP*4.0b10 [33] and PHYMLv2.4.4 [34]. Maximum-likelihood parameters and the evolutionary model were determined using the hierarchical Likelihood Ratio Test (hLRTs) and the Akaike Information Criterion (AIC) [35,36] implemented in Modeltest 3.7 [37] in conjunction with PAUP*4b10. To assess the reliability of the nodes in the ML tree, a bootstrap analysis using 1000 pseudo-replicates was made with PHYML. Since haplotype sequences for T. dimidiata individuals (populations) are quite similar and potentially subject to homoplasy and recombination, alternative procedures to phylogenetic tree reconstruction revealing their relationships were tested. Therefore, a median-joining network analysis [38] was performed using Network version 4.1.1.2 (available from Fluxus Technology Ltd., http://www. fluxus-engineering.com) with the variable positions in the multiple alignment of the different ITS-2 haplotypes from T. dimidiata populations.
Alternative methods of phylogenetic reconstruction allowing an evaluation of the support for each node were also applied. A distance-based phylogeny using the neighbor-joining (NJ) algorithm [39] with the ML pairwise distances was obtained. Statistical support for the nodes was evaluated with 1000 bootstrap replicates, with and without removal of gapped positions. Finally, a Bayesian phylogeny reconstruction procedure was applied to obtain posterior probabilities (BPP) for the nodes in the ML tree. We used the same evolutionary model as above implemented in MrBayes 3.1 [40] with four chains during 1,000,000 generations and trees were sampled every 100 generations. The last 9,000 trees were used to obtain the consensus tree and posterior probabilities.

Genetic variation studies
Genetic variation within and among populations of T. dimidiata was evaluated using DnaSP version 4 [41] and Arlequin 2000 [42]. Summary parameters include those based on the frequency of variants (haplotype number and diversity) as well as some taking genetic differences among variants into account (gene diversity, polymorphic sites). A hierarchical analysis of molecular variance (AMOVA) was performed using Arlequin. This analysis provides estimates of variance components and F-statistics [43] analogs reflecting the correlation of haplotype diversity at different levels of hierarchical subdivision. Unlike other approaches for partitioning genetic variation based on the analysis of variance of gene frequencies, AMOVA takes into account the genetic relatedness Table 1. Triatoma species and samples studied, including ITS-2 sequence length and AT composition (in percentage). between molecular haplotypes. The hierarchical subdivision was made at three levels. At the top level, different groups were defined on the basis of the phylogenetic relationships for the different T. dimidiata haplotypes obtained. The second level corresponded to countries of sampling within each of these groups, and the third level corresponded to the different haplotypes found in each country within group. AMOVA reports components of variance at the three levels under consideration (among groups, among countries within groups, and within countries within groups) as well as F-statistics analogs. Under the present scheme, F ST is viewed as the correlation of random haplotypes within countries within groups, relative to that of random pairs of haplotypes drawn from the whole species, F CT as the correlation of random haplotypes within groups, relative to that of random pairs of haplotypes drawn from the whole species, and F SC as the correlation of the molecular diversity of random haplotypes within countries within groups, relative to that of random pairs of haplotypes drawn from the corresponding group [44]. Although in the program used (only currently available for molecular variance analysis) the choice for establishing an intermediate level is fully arbitrary and has no influence on the final result of the comparison between units at the higher level, these same analyses were repeated by considering each haplotype, which may encompass several individuals, as a separate group for this intermediate level, because it could be argued that geopolitical country borders was not an appropriate choice despite its interest from the point of view of the control of Chagas disease. The statistical significance of fixation indices was tested using a non-parametric permutation approach [44]. Genetic differentiation between pairs of populations was evaluated by means of F-statistics [43]. Exact tests of population differentiation were performed [45]. Slatkin's linearized F ST 's [46,47] procedure was also followed to obtain estimates of pairwise equilibrium migration rates, both among groups, among countries within groups, and within countries for those cases in which haplotypes from more than one group were present.

Sequence Analyses of Triatoma dimidiata Populations
The 137 ITS-2 sequences revealed the existence of 31 different haplotypes in the T. dimidiata studied (T.dim-H1 to T.dim-H31) (see Tables 1 and 2 for localities and countries). Their length was 489-497 base pairs (bp) (mean, 495.10) with a relative AT-biased nucleotide composition of 75.25-76.85% (75.72%). Sequence similarity analysis of these 31 haplotypes revealed four distinct groupings: grouping 1 (T.dim-H1 to T.dim-H10); grouping 2 (T.dim-H11 to T.dim-H17); grouping 3 (T.dim-H18 to T. dim-H24); and grouping 4 (T. dim-H25 to T. dim-H31) ( Figure 2). These four groupings appear linked to concrete wide geographical areas including neighboring countries and regions. The only exception is Providencia Island, which, although part of Colombia, is located 720 km off the northern coast of Colombia but only 240 km off the western coast of Nicaragua. No haplotype presents a very broad geographical distribution.
The alignment of the 31 T. dimidiata haplotype sequences was 501 bp-long, of which 450 characters were constant and 24 were parsimony-informative. The interrupted microsatellite (AT) [4][5] TTT (AT) [5][6][7] was detected between positions 47 and 73 in all specimens studied. Variability in this microsatellite region and their respective sequence positions are noted in Figure 2.

Sequence Analyses in the Phyllosoma and Rubrofasciata Groups
ITS-2 sequences of T. bassolsae, T. bolivari, T. hegneri, T. mexicana, T. pallidipennis, T. ryckmani, T. flavida, T. nitida, T. gerstaeckeri, and T. rubida, including haplotype length and AT content are listed in Table 1. The comparison analyses which include these ITS-2 sequences and those of the Phyllosoma and Rubrofasciata groups (available in GenBank) provided 48 different haplotypes. Their alignment resulted in a total of 551 characters including gaps, of which 365 sites were constant and 99 parsimony-informative.

Sequence Analyses in the Infestans Group
ITS-2 sequences of T. maculata and T arthurneivai, including haplotype length and AT content are listed in Table 1.
The ITS-2 of T. maculata fits very well within sequences of the Infestans complex species studied in the present work, a total of 6-19 (13.7) mutations, namely 6-11 (7.25) ti and 0-10 (6.5) tv, appearing when comparing the five Infestans complex species in question. The material of Triatoma arthurneivai here analyzed is very

Phylogenetic Analyses
Two different phylogenetic approaches were performed with the 31 T. dimidiata haplotypes, both yielding coincident results. A maximum likelihood tree was reconstructed using the best model of evolution as determined by the lowest AIC, which was GTR+I (2Ln = 887.089), being the proportion of invariable sites (I) of 0.166. Three groups appeared with high support values indicating that their differentiation was not due to random sampling of a low variable sequence (tree not shown). The large group 1 encompassed haplotypes from all the countries, whereas groups 2 (Mexico and Guatemala) and 3 (Mexico, Guatemala and Honduras) were more geographically restricted.
Alternatively, a median-joining network was reconstructed with the 31 different T. dimidiata sequences using the variable sites in the multiple alignment (Figure 3). This network showed the same three groups found in the ML tree. Group 1 occupies a central position in the network and is the most widespread and variable group, so that it most likely corresponds to the ancestral or source set. This is further reinforced by the direct relationship between this group and the two others, more geographically restricted and encompassing fewer variants, group 2 including samples from Mexico and Guatemala, and group 3 including samples from these two countries and Honduras. The group 1 source set would in turn be derived from group 3, which might be interpreted as a geographically restricted relict according to the phylogeographic results. Moreover, sequence variants in group 1 are clustered in two different subgroups, with genetic and geographical borders: subgroup 1A includes sequences from Colombian Providencia island, Ecuador, Guatemala, Honduras, Mexico (only South of Chiapas) and Nicaragua; subgroup 1B encompasses sequences from continental Colombia and Panama. The two closest sequences of each subgroup differ in two sites, which might correspond to haplotypes not found in this sampling.
The relevance of the ITS-2 differences among these T. dimidiata groups and subgroups was assessed by comparison with other Triatoma species. Therefore, a multiple, 562-nucleotide-long alignment was obtained by incorporating 22 additional ITS-2 sequences. This set includes 53 ITS-2 sequences of Triatoma species and, using R. prolixus as outgroup, a ML tree was obtained (2Ln = 2648.5129) using the HKY+G model, according to the AIC results with a gamma distribution shape parameter = 0.58. This tree (Figure 4) shows that: N the 31 T. dimidiata haplotypes appear within a highly supported clade (95/97/100 in ML/NJ/BPP), distributed as follows: a first large subclade, also very well supported (99/97/100), comprising subgroup 1A, subgroup 1B, group 2, and group 3 of the network analysis; subgroup 1A (sequence grouping 1 = T.dim-H1 to T.dim-H10) includes populations from Central America (Honduras, Nicaragua, Guatemala and scattered haplotypes from Mexico, Ecuador and Providence Island); interestingly, the haplotype T.dim-H10 corresponding to phenetically peculiar specimens found in cave-dwellings of Lanquin, Guatemala, appears independent although related to the rest with very high supports; subgroup 1B (sequence grouping 2 = T.dim-H11 to T.dim-H17) comprises populations from continental Colombia and Panama and appears as a monophyletic haplotype cluster; group 2 (sequence grouping 3 = T.dim-H18 to T.dim-H24) shows a well supported branch (91/92/100) and comprises populations from Mexico (Gulf coast, high plains, and Cozumel island) and Guatemala, including the two T. hegneri haplotypes; the second large clade is also highly supported (97/96/100), corresponding to group 3 (sequence grouping 4 = T.dim-H25 to T.dim-H31) and includes populations from the Yucatan peninsula, Holbox and Cozumel islands and northern Chiapas (Mexico), northern Honduras and northern Guatemala;   Triatoma dimidiata groupings appeared well supported, with very high bootstrap proportions (BP.90%) using ML and neighborjoining reconstruction and the highest Bayesian posterior probabilities (BPP = 100%). Similar levels were found for other well established Triatoma species, many of which showed substantially lower support values in the three statistical measurements employed. However, other species presented no ITS-2 nucleotide differences (T. picturata and T. longipennis; T. mazzotti and T. phyllosoma).

Genetic Variation Analyses
The phylogenetic analyses showed that samples from the same country may belong to different clusters. This result, on its own, is not enough to demonstrate the biological distinctiveness of the corresponding populations. Sampled individuals may represent a minor fraction of the total genetic variability in a highly heterogeneous population and the sampling procedure might have resulted, by pure chance, in the observed clustering of some variants. Given that each of these clusters holds some genetic variability of its own, the first task was to evaluate whether the observed groupings were significantly different from each other, in terms of genetic variation, by partitioning the observed genetic variability at three different levels: among groups, among populations (countries) within groups, and within populations. A hierarchical analysis of molecular variance was used to test the null hypothesis of no genetic differentiation among groups considering variation at lower levels. This procedure was first applied to T. dimidiata sequences using three levels as defined above (Table 3a). Most of the genetic variation found was allocated to the among groups level (80.24% of the total variation), with much lower portions of variation assigned to differences among populations within groups level (11.71%) and within populations level (8.05%), although both were still statistically significant after 1000 pseudorandom samples generated for testing. This indicates that, despite genetic variation within and among populations at these three levels, there is a substantial amount of genetic differentiation among them that justifies their consideration as separate groupings for further analysis. The same results were obtained, notwithstanding small numerical differences due to the different numbers of groups, when haplotypes instead of countries were considered at the intermediate level (Table S1). The geographical fitting represents in fact no surprise at all, taking into account that the distribution of T. dimidiata covers different countries which are more or less aligned following a north-south axis because of the relatively slenderness of the Central American bridge. Hence, as any of the two versions of the analyses conveys the same information and leads to the same conclusions, and which one should be reported is simply a matter of opinion, the first considering countries becomes practically more useful because Chagas disease control measures are organized at national level.
The median-joining network reconstructed with the 31 different T. dimidiata ITS-2 sequences revealed the existence of three distinct groups (groups 1, 2 and 3), the first of which further subdivided into two subgroups 1A and 1B. The same AMOVA procedure was applied to ascertain whether these two subgroups could be considered as distinct populations or not. The results (Table 3b) indicate that a significant fraction (60.15%) of the total genetic variation corresponds to differences between these two subgroups which, correspondingly, could be considered as separate populations for the ensuing analyses.
Based on the four groups/subgroups previously described in the median-joining network, a summary of relevant population genetic parameters for T. dimidiata is presented in Table 4. Genetic variation in T. dimidiata populations was quite evenly distributed, with similar levels of nucleotide and haplotype diversities in the four groups/subgroups considered. Nevertheless, for all the parameters studied, subgroup 1A presented higher values than the rest, although significance of the differences was only obtained for haplotype diversity. A similar summary is shown for each country sample within groups in Table S2.
Different estimates of h were obtained based on the expected heterozygosity, the expected number of alleles, the number of polymorphic sites and the nucleotide diversity. The four estimates were quite consistent for the four groups/subgroups and they agreed in assigning a larger value to subgroup 1A. Differences in the genetic composition of the four groups/ subgroups 1A, 1B, 2 and 3 have previously been shown to be statistically significant according to analyses of molecular variance. A further evaluation of this distinctiveness was made (Table 3c), in which the four groups/subgroups were considered for the AMOVA, in correspondence with the previous results. In this case, the amount of among-group variation rose to 86.84% of the total variation, whereas among population within groups and within population levels they were substantially lower, 3.21% and 9.95% respectively.
Genetic differences within and among the ITS-2 locus for T. dimidiata samples were further explored through pairwise comparisons, and estimates of average pairwise differences within and among the four groups/subgroups considered were obtained (Table 5). Subgroup 1A presented the largest value for withingroup pairwise differences. The within-population values were much lower than among-populations comparisons. Among the latter, the smallest number of differences was found between subgroup 1A and 1B, in correspondence with their close phylogenetic relationship. Subgroup 1B was the one with the lowest overall number of pairwise differences, slightly below 1A. On the contrary, the highest value of pairwise differentiation corresponds to group 3, with almost 20 differences (corrected estimate) when compared with any other group.
Within groups genetic differentiation was evaluated by computation of pairwise F ST values for populations defined by country of origin (Table S3). Since all groups/subgroups, with the only exception of subgroup 1A, are characterized by one large (n.10) and several small (n,10) populations, significance values for test of genetic differentiation have to be interpreted cautiously. Hence, there is no apparent differentiation between two populations in subgroup 1B (Colombia2, n = 30, and Panama, n = 4) and similarly in group 2 (Mexico2, n = 23, and Guatemala2, n = 4). The only significant value found in group 3 corresponds to Honduras3 (n = 2) and Guatemala3 (n = 7), for which F ST = 0.529, P,0.05. None of these two populations presented significant differentiation with respect to the largest population in this group, Mexico3 (n = 15). Subgroup 1A includes two large populations, Honduras1 (n = 18) and Guatemala1 (n = 26), which presented a highly significant F ST = 0.193, P,0.001. Although this value, under the assumption of migration-drift equilibrium, corresponds to an estimate of 2.1 migrants per generation between both populations, which would be enough to prevent their complete differentiation, such estimations shall be verified by using larger samples and markers better suited for population genetics analyses. Comparisons between each of these two populations and the smaller ones in subgroup 1A revealed that Honduras1 differed from Mexico1, Guatemala1 was different from Ecuador and Nicaragua, and none of them differed from the only two individuals from Providencia island. Similar comparisons for all pairs of populations assigned to different groups/subgroups resulted in highly significant F ST values (Table S4).

Discussion
Triatoma dimidiata, T. sp. aff. dimidiata and T. hegneri The highest intraspecific ITS-2 variability (absolute nucleotide differences including indels) known in Triatomini members is 2.70% (13/482) in T. infestans specimens collected throughout the very wide geographical distribution of this species [25]. Hence, the result of 10.18% ( = 51/501) detected in T. dimidiata (Figure 2) appears to be pronouncedly outside the limits of the intraspecific variability range known for Triatoma species. Group 3 is the main responsible for such differences (Table 5) and shows a high 2.42% divergence within itself, suggesting an old origin in the light of the relatively reduced geographical area of distribution of these haplotypes in Mexico (Yucatan, Chiapas, Cozumel Island and Holbox Island), Guatemala (Peten) and Honduras (Yoro) only. The time of divergence between group 3 and other T. dimidiata populations was estimated to be of 5.9-10.5 million years ago (Mya) according to a molecular clock analysis based on rDNA evolutionary rates [4].
The divergence of 5.62% shown by the other 24 ITS-2 haplotypes (Figure 2) also appears to be too large, in spite of the wide geographical area they occupy from Mexico down to Ecuador, suggesting a speciation process. However, population average pairwise differences between subgroup 1A, subgroup 1B and group 2 are markedly lower than between these three and group 3 (Table 5), and intragroup differences do fall within the above-mentioned Triatomini range: 2.61% within subgroup 1A, 2.41% within subgroup 1B, and 2.01% within group 2.
Results indicate that several T. dimidiata populations are following different evolutionary divergences in which geographical isolation appears to have had an important influence ( Figure 5). A phenotypic consequence of that process had been observed by other specialists before, who wrote about an assemblage of morphologically variable populations [10]. More recently, significant head shape differences between populations showed a separation between northern, intermediate and southern collections of T. dimidiata and also support an evolutionary divergence of populations within this species [13].
Three subspecies were distinguished on the basis of morphological differences [48,49]: (i) T. d. dimidiata concerns the first description of the species in Peru (no type specimen available; no type locality assigned, but undoubtedly from northern Peru, probably around the locality of Tumbes, near Ecuador) and corresponds to most of the Central American forms; (ii) T. d. maculipennis was proposed for specimens from Mexico (type specimen in Zoologisches Museum Berlin) and corresponds to forms with relatively short heads and large eyes; and (iii) T. d. capitata was proposed for large size specimens typified by longer heads and smaller eyes originally found in Colombia (type specimen in the Academy of Sciences of California). However, these subspecies became later synonymized after results of a morphological re-examination which were interpreted as evidence of a clinal variation along a north-south axis [50,51].
Present ITS-2 sequences and corresponding phylogenetic and genetic variation analyses support the appropriateness to (i) differentiate group 3 as a species of its own (here simply designed as T. sp. aff. dimidiata to avoid further systematic confusion with T. dimidiata, according to taxonomic rules), and (ii) re-assign subspecific status for subgroup 1A, subgroup 1B and group 2. Results of the present study do not support the rise of the abovementioned subspecific taxa to species level for the time being, although it is evident that in the three cases relatively long divergence processes have taken place. Similar genetic studies with other molecular markers may contribute to a more complete assessment of these evolutionary isolation and speciation processes.
The taxon T. sp. aff. dimidiata concerns group 3. This species seems to represent a relatively relict species with a distribution restricted to the Mexican flat areas of the Yucatan peninsula and the northern part of Chiapas state, the northern lowland of Guatemala (and probably also Belize), and only one altitudeadapted haplotype (T.dim-H29) in its most extreme border populations in northern Honduras. The most widely spread haplotype T.dim-H28 is also present in the small island of Holbox and the large island of Cozumel, both near the Yucatan coast, suggesting that this haplotype should be considered the oldest of this species. This species is also of public health importance because of its capacity to transmit Chagas disease [52,53] and the control problems it poses [54,55].
The taxon T. d. dimidiata corresponds to subgroup 1A and populations mainly from Guatemala and Honduras and secondarily Mexico, Nicaragua and Ecuador. The population of the Colombian island of Providence undoubtedly derives from the most widely dispersed haplotype T.dim-H1 on the nearest Caribbean coastal area of Central America and not from continental Colombia. The present populations in Ecuador may derive from introduced specimens originally from the Guatemala-Honduras-Nicaragua region, relatively recently introduced by humans [4], very probably in the period of the early colonialization of northwestern South America by the Spanish 'conquistadores' in which exchange activities between Central American settlements and the Peruvian Tumbes area took place [56]. The type specimens of the original description of the species in northern Peru might also belong to populations derived from such man-made introductions from Central America. The haplotype T.dim-H10 of Lanquin, Alta Verapaz, Guatemala appears in the network analysis as directly derived from an ancestor which gave rise to the subspecies T. d. dimidiata. An isolation phenomenon in caves may explain the albinic characteristics of the specimens presenting this haplotype. These cavernicole specimens from Alta Verapaz have already shown their peculiarity in morphometric and cuticular hydrocarbon studies [13,17].
The taxon T. d. capitata corresponds to subgroup 1B and populations from Colombia and Panama. The isthmus of Panama and the separation/joining process of South and North America towards the end of the Pliocene (3-5 Mya) [57], in a period in which several more or less closely separated islands appeared and evolved up to their fusion into the isthmus, should have played a major role in the isolation and subsequent divergence of these southernmost T. dimidiata populations. The lack of relationship between the haplotypes of Ecuador and those of Colombia is worth mentioning, as the geographical closeness of these two countries could have given rise to the erroneous hypothesis of Colombian forms having derived from Ecuadorian populations. In a recent study of three populations of sylvatic, peridomestic and domestic T. dimidiata from Colombia, the estimated low genetic distances based on RAPD analyses did not discriminate the populations studied, indicating that they maintain the genetic identity of a single recent common ancestor [9].
The taxon T. d. maculipennis corresponds to group 2 and populations mainly from Mexico, but rarely found in Guatemala. According to the network analysis, this subspecies seems to have derived from group 1 probably by isolation in the Mexican part northward from the isthmus of Tehuantepec. Similarly as for other organisms including insects [58], the mountainous Sierra Madre chain throughout southern Mexico and Guatemala areas near the Pacific coast probably played also a role in that isolation process through an area where T. sp. aff. dimidiata did not represent a competition barrier, as T. sp. aff. dimidiata appears to be preferentially a low altitude species in these two countries.
Southern Mexico (including the Yucatan peninsula and Chiapas state) and almost the whole country of Guatemala (at least ten departments) constitute a crucial evolutionary area, where a high number of taxa, including T. d. dimidiata, T. d. maculipennis, and T. sp. aff. dimidiata, overlap. In a morphometric analysis, populations from San Luis Potosi and Veracruz in Mexico were indistinguishable while clearly different from populations from Yucatan in Mexico and Peten in Guatemala [14]. The former correspond to T. d. maculipennis and the latter to T. sp. aff. dimidiata. In Guatemala, a high degree of genetic variation in T. dimidiata sensu lato was shown by RAPD-PCR [12], demonstrating a limited gene flow between different provinces, although barriers between the Atlantic and Pacific drainage slopes did not appear to be significant limiters of a gene flow, according to a hierarchical analysis.
The two haplotypes of T. hegneri differ by only 3 mutations from haplotypes of T. d. maculipennis. This reduced number of nucleotide differences and the location of T. hegneri haplotypes within the clade of T. dimidiata, basal to haplotypes of group 2 (Figure 4), does not support its status as an independent species. The results obtained suggest that it is an insular form of T. d. maculipennis. Originally described from the island of Cozumel [3], a subspecific status T. d. hegneri could be maintained only if morphological characteristics allow a clear differentiation of the insular form, as the phylogenetic analysis somehow separates it in a very close but particular evolutionary line. Triatoma hegneri, although chromatically distinguishable from most forms of T. dimidiata [50], is known to produce fertile hybrids when experimentally crossed with T. dimidiata (R.E. Ryckman, unpublished). Interestingly, the most dispersed haplotypes of both T. d. maculipennis (T.dim-H18) and T. sp. aff. dimidiata (T.dim-H28) are also present on the same island, probably introduced through the intense human transport between the mainland and the island.
The distinction between T. d. dimidiata (subgroup 1A), T. d. capitata (subgroup 1B), T. d. maculipennis (group 2), T. sp. aff. dimidiata (group 3), and T. d. hegneri contributes giving systematic/ taxonomic coherency to present knowledge about morphological and genetic concepts in these taxa. From an ancestral form close to T. sp. aff. dimidiata, it can be postulated that an original diversification focus of T. dimidiata forms took place most probably in Guatemala, with a southern spread into Panama and Colombia to give the capitata forms and a northwestern spread into Mexico to give the maculipennis forms ( Figure 5). Thus, the results of the present paper, obtained from a large amount of samples of T. dimidiata from many different countries covering its whole latitude range, gives rise to a new frame that is different from the previous hypothesis about a clinal variation along a north-south axis, which was formerly suggested to explain both morphological data [50] and preliminary ITS-2 data from a reduced number of samples [6].
Moreover, the distinction between these five entities may facilitate the understanding of different vector transmission capacities and epidemiological characteristics of Chagas disease throughout the very large area where T. dimidiata sensu lato is distributed, from the Mexican northern latitude limit up to the Peruvian southern latitude limit [11]. Recent results obtained by means of a population dynamics model indicate that T. dimidiata in Yucatan, Mexico, is not able to sustain domestic populations, that up to 90% of the individuals found in houses are immigrants, and that consequently Chagas disease control strategies must be adapted to a transmission by non-domiciliated vectors [59]. This might be considered surprising because it does not fit the domiciliation capacity of T. dimidiata in other places, but it appears to be congruent if it is taken into account that in fact the Yucatan vector in question is not T. dimidiata but a different species T. sp. aff. dimidiata.
The results here obtained also suggest that T. d. dimidiata in Ecuador is a good candidate for the design of appropriate vector control intervention, similarly to domestic T. infestans populations in countries such as Uruguay, Chile and Brazil within the successful Southern Cone Initiative [60]. The control and even eradication of T. d. dimidiata in Ecuador by means of insecticidespraying of its domestic habitats might be successful, if it is considered that it is merely an introduced vector species in that area, and a priori it would have difficulties in escaping from the insecticide activity because of its non-adaptativeness to the sylvatic environment in these two countries [61]. Unfortunately, such a control initiative will not be so easy to carry out in Colombia, as results prove that Colombian forms are authochthonous T. d. capitata and not T. d. dimidiata derived from the Ecuadorian introduced form. This fits with the existence of sylvatic populations in Colombia and with the high genetic similarity of sylvatic, peridomestic and domestic populations detected in that country [9]. Similarly to in Colombia, results indicate that T. dimidiata will offer, because of being authochthonous forms, more problems for insecticide-spraying control in Central American countries than introduced T. infestans in Southern Cone countries.

The other Meso-and Central American Triatoma Species
Triatoma bassolsae differs by only one deletion from T. pallidipennis and appears in the branch of the 5 species traditionally included in the Phyllosoma complex: T. longipennis, T. mazzotti, T. picturata, T. pallidipennis and T. phyllosoma. The genetic differences between these taxa are so reduced (sometimes even none at all), that there is no support to maintain them as separated species. Such a low number of nucleotide differences in the ITS is considered as pertaining to organisms able to hybridize [62]. This fully fits the capacity of these taxa to crossbreed and give fertile hybrids [63,64] and agrees with the entomologist conclusion of applying only subspecies level to them [49]. The divergence of members of the phyllosoma complex is estimated at only 0.74-2.28 Mya by the rDNA molecular clock [4], which also seems consistent with a subspecific rank. All further ITS-2 studies have always reached the same conclusion [5,6,65]. By analyzing many interfertility experiments [64], it can be concluded that, in triatomines, morphological differentiation appears to be faster than the installation of reproductive or genetic barriers [66,67]. Rapid morphological changes, associated with ecological adaptation, helps to explain discordance between phenetic and genetic differentiation. Triatomine species with consistent morphological differences would arise through divergent ecological adaptation, a vision which fits with ''evolutionary units'' implying a different evolutionary direction taken by some populations [67]. Until future reproductive isolation thanks to ecological isolation is reached by these morphologically different entities of the Phyllosoma complex, the subspecies concept accurately fits for all these ''evolutionary units'' of the Phyllosoma complex. ITS-2 results indicate that Triatoma bassolsae is one additional taxon to be included in this situation, as has already been suggested [65]. The comparison of the small genetic divergences between these taxa, their distributions exclusively restricted to regions of Mexico, and their different geographical distribution areas slightly overlapping in their bordering zones [3] suggest that genetic exchange might be impeding or delaying definitive divergence processes to reach species level.
Genetic distances between the taxa of the Phyllosoma complex found when analyzing different mtDNA genes proved to be similar to those detected in ITS-2 at the 16S [68], but higher in CytB [65,69], and COI [69]. This agrees with the evolutionary rates of the protein-coding mtDNA genes which are pronouncedly faster than the one of ITS-2. Moreover, aminoacid sequences of the CytB and COI genes show no one difference between the Phyllosoma complex members studied (all are silent mutations or synonymous substitutions) except one aminoacid difference between two populations of the same species T. pallidipennis and one in T. picturata versus the rest [69], which also fit with an intraspecific variability. Additionally, it shall be taken into account that (i) mtDNA becomes monophyletic more rapidly than does a single nuclear gene and far more rapidly than a sample of several nuclear genes, so that mtDNA may make inferences of specieslevel monophyly erroneous [70], and (ii) the known great potential of mtDNA to become monophyletic by selective sweeps can decrease the time to monophyly of a clade and not be reflective of the genealogical processes in the nuclear genome, advantageous mutations occurring on mtDNA causing the entire mitochondrial genome to become monophyletic because of the little or no recombination they have [71]. The crossbreeding capacity and hybrid viability among the Phyllosoma complex taxa in question is well known and, taking into account that their geographical distributions overlap in their border areas and there are no sufficient ecological differences indicating a local spatial separation, it becomes very difficult to support them as separate species from the evolutionary, biogeographical and ecological points of view because there is apparently no barrier for a reproductive isolation. Thus, the results of both ITS-2 and mtDNA genes fit with such an evolutionary, subspecific divergence, when taking into account the peculiarities of both nuclear and mitochondrial markers.
Triatoma mexicana appears to be a good species and its location in the phylogenetic tree fully supports its ascription to the Phyllosoma complex, similarly as suggested by a phylogentic analysis by means of a mtDNA CO1 fragment [69]. Surprisingly, T. gerstaeckeri (Rubrofasciata group) clusters with T. mexicana, suggesting that it should be included in the Phyllosoma complex. All these species, i.e. T. phyllosoma (including its subspecies phyllosoma, longipennis, mazzotti, picturata, pallidipennis and bassolsae), T. dimidiata (with its three subspecies dimidiata, capitata and maculipennis, to which hegneri shall be added), T. sp. aff. dimidiata, T. mexicana and T. gerstaeckeri constitute a well defined clade for which the generic taxon Meccus, proposed long ago [72], afterwards synonymized [50] and recently tentatively revalidated [73], seem to appropriately fit. Previous molecular studies, first with complete ITS-2 sequences [74] and second with partial mtDNA 16S gene sequences [68], also indicate that Meccus might be a valid taxon.
The revalidation of Meccus, as well as that of Nesotriatoma for species of the Flavida complex, has not been accepted because of the close relationship between T. flavida and the Phyllosoma complex [7]. The results of the present study do, however, pose a serious question concerning the inclusion of species as T. bolivari and T. ryckmani in the Phyllosoma complex, as they appear to cluster with T. rubida of the Rubrofasciata group with relatively high support (83 and 96 in ML and BPP, respectively). A T. rubida -T. nitida clade previously detected with weak support under certain conditions in mitochondrial DNA marker analyses [69] does not appear to be supported in the ITS-2 phylogeny.
Although not fully resolved in the tree obtained, the location of the Cuban T. flavida as a species basal to all other North-Central American Triatoma species may be interpreted as a consequence of being a relict insular species close to the ancient first North-Central American Triatoma colonizers. Further studies with other genetic markers are needed to establish the position of T. flavida more adequately.

The South American Triatoma Species
The very scarce ITS-2 sequence differences between T. arthurneivai and T. rubrovaria, a species known in southern Brazil, Uruguay and northern Argentina [75], pose doubts on whether to keep the validity of T. arthurneivai as independent species. Recent genetic and morphometric studies have already raised several questions about T. arthurneivai, indicating that topotypes from Minas Geraes may represent a species different from populations of São Paulo State formerly also referred to T. arthurneivai and suggesting that these São Paulo populations might probably belong to T. wygodzinskyi [76]. This may explain the ITS-2 results, as the two specimens analyzed in the present paper come in fact from Espirito Santo do Pinhal, São Paulo State. Consequently, material of typical T. wygodzinskyi should be sequenced and compared to both true T. arthurneivai from Minas Geraes and T. rubrovaria to ascertain the status of these three taxa.
The South American Triatoma species cluster together with maximum support (100/100/100) and well separated from that of the North and Central American species of the same genus, thus supporting results of previous analyses which indicate an early divergence of about 23-38 Mya between species of the northern (Phyllosoma complex) and southern (T. infestans) continent [4,6].

Supporting Information
Alternative Language Abstract S1 Translation of the abstract into Spanish by S. Mas-Coma.