Study area

Tissue samples used in this study were distributed from French Guiana and the Brazilian states of Amapá and Pará, situated in the eastern portion of the Guiana region. This region is well-preserved, with large areas of forest protected by a network of conservation units on both sides of Brazil and French Guiana [42, 43, 44]. The main forest areas are found in the western portion of the study area, located on substrates of Precambrian origin [45]. Isolated savannahs are found throughout the study area [46] located between the forests of the interior portion and the coastal region. The coastal zone of southern Amapá is influenced primarily by the adjacent Amazon River, which forms extensive areas of seasonally flooded habitats, with many lakes and islands [47]. The northern coast of the state, north of the Araguari River to French Guiana, encompasses tracts of mangroves interspersed with flooded fields influenced by the Atlantic Ocean (Fig 1). The Brazilian portion of the study area is the most impacted, with extensive deforestation and human colonization along the Atlantic coast and the northern margin of the Amazon River in Pará State [48]. A more recent impact is the ongoing expansion of soybean plantations into the savannah habitats of Amapá State [46, 49]. In French Guiana, the principal threats are mining [50, 51], and deforestation, principally on the coast [43], although forest loss remains below the rates observed in other Amazonian countries [52].

Samples

All tissues samples used in this study were obtained from voucher specimens of the Collection of Amapá Fauna of the Institute for Scientific and Technological Research of the state of Amapá (IEPA) in Macapá; National Institute of Amazon Research (INPA) in Manaus (both in Brazil), and Collection JAGUARS, Association Kwata and Institut Pasteur de la Guyane, Cayenne. All analysed samples were obtained specimens deposited in the aforementioned scientific collections, thus no collecting permits were required. Catalog numbers are given in S1 Appendix. INPA´s Ethics Committee on Animal Use (CEUA) of INPA granted formal waiver of ethics approval due to that tissue of the sequences we generated were obtained from scientific collections.

We sequenced 1100 base pairs (bp) of the CYTB gene in 37 samples of Proechimys cuvieri from 22 localities and 54 samples of P. guyannensis from 20 localities within the study area, between northern French Guiana and Brazil, north of the Amazon River (Amapá and Pará). The two species were sympatric at 11 sites (Fig 1). We also sequenced 800 bp of CYTB in 8 samples collected in Brazil (Amazonas and Roraima states), Venezuela, and Peru (Fig 1) for phylogenetic analyses. Only 25 additional sequences were obtained from GenBank. We analyzed a total of 124 sequences of P. cuvieri and P. guyannensis in 50 localities; details on specimens, localities and GenBank accession numbers, are given in S1 Appendix.

Sequencing

The DNA was extracted from tissue samples, obtain of the muscle and skin byopsy, preserved in 95% ethanol. The DNA was extracted using the NucliSENS EasyMag robot (Biomérieux, Craponne, France) and Wizard Genomic DNA Purification (Promega) kits, according to the manufacturers’ protocols. Three overlapping fragments of the cytb were amplified using internal and external primers. We amplified the first fragment, of 400–500 bps, using the primers H6 (5'TCTTCCATGAGGACAAATATC3') and L15 (5’TCTCCATTTCTGGTTTACAAGAC3'), and the final fragment (also 400–500 bps) using the primers H8 [28] and L2 (5' TACCATGAGGACAAATATC’). The intermediate fragment, with approximately 400 bps, was amplified using the primers PRO 197F (5´TTACACAYATTTGYCGAGAYG3´) and PRO 665R (5´GGGTGRAATGGRATTTTGTCTGA3´), which were designed specifically for this study. The samples were amplified in 40 μL reaction volumes using Platinum Taq (Life Technologies Corp.) containing 4.0 μL of each primer, 5.0 μL of reaction buffer, 3.0 μL of MgCl² 50 mM, 5.0 μL of dNTP mix, 0.5 μL of polymerase, 2.0 μL of the DNA template, and 16,5 μL of H₂O. The polymerase chain reaction (PCR) cycle consisted of initial denaturation at 94°C for 3 min, followed by 30 amplification cycles of denaturation at 94°C for 30 sec, annealing at 45°C for 1 min and extension at 72°C for 2 min, with a final extension at 72°C for 10 min.

The PCR products were sent to Beckman Coulter Genomics (Takeley, UK) for purification and sequencing. The sequences were aligned and checked manually using MEGA 6.0 [53]. We did not observe stop codons, and transition and transversion rates were normal.

Analyses

We implemented Bayesian Inference (BI) and Maximum Likelihood (ML) phylogenetic analysis to reconstruct the relationships within the studied populations of P. guyannensis and P. cuvieri from northeastern Amazonia using the complete dataset. Sequences of the rodent species Trichomys apereoides (AY083332, AY083341) and Hoplomys (HM544128, NC020657) were obtained from GenBank and included in the analyses as outgroups.

The program jModelTest version 3.7 [54] was used to select the HKY+G like the best evolutionary model for the datasets based on Bayesian Information Criterion (BIC). The BI phylogenetic tree was reconstructed in MrBayes 3.2.1 [55], using 10⁷ Markov Chain Monte Carlo (MCMC) generations, with the first 25% being discarded as burn-in. The ML analysis was run in RAxML in the CIPRES Science Gateway [56]. Nodal support was estimated with 100 bootstrap replicates. We estimated the divergence time between the lineages using the relaxed lognormal uncorrelated molecular clock and the Yule Process pure-birth speciation model, run in BEAST 1.6.2 [57]. The substitution rate used in the analysis was the 3.6% divergence per million year, derived from Fabre et al. [58]. We used the following times inferred from Fabre et al. [59]: divergence between [Trichomys] and [Hoplomys + Proechimys spp.] set at 9.4 Ma (95% CI: 8.3–10.9) and divergence between Hoplomys and Proechimys set at 5.5 Ma (95% CI: 4.6–6.7), with these dates being used as calibration for the most recent common ancestor (TMRCA), considering a normal distribution. Since no more recent calibration point was available, we confirmed the absence of substitution saturation on the [Trichomys + Hoplomys + Proechimys] dataset with the Xia test [60] implemented with DAMBE [61]. We generate 10⁷ Markov Chain (MCMC) generations, with the first 25% being discarded as burn-in. We used TRACER 1.6 to evaluate the convergence of the analysis and the parameters estimated, considering an Effective Sample Size (ESS) equal to or greater than 200 [62].

We generated haplotype networks using Network 4.6 and the median-joining (MJ) method [63]. BAPS5 (Bayesian Analysis of Population Structure, [64]) was employed to verify population structure. We ran the mixture model setting the maximum possibility of clusters [64]. The genetic distance between clusters identified in BAPS5 and between clusters composed for sequences with distribution between the rivers (Oiapoque–Maroni; Araguari–Oiapoque, Jari–Araguari, and Jari–Trombetas) (Fig 1). Genetic distance was calculated using MEGA 6.0 [53].

We considered the isolation effect of different river drainages to evaluate gene flow patterns between geographic areas. We used the program Migrate [65] under a Bayesian framework to estimate effective population sizes (Θ) and migration rates M between different populations from north to south in studied area. We ran four chains with thermodynamic integration and a static heating swap scheme (temperatures: 1.0, 1.5, 3.0, 10⁶), sampling at every 100^th step for a total of 20⁶ steps and a burn-in of 50³ steps.

We calculated Fu’s Fs [66], Tajima’s D [67] and R2 [68] to assess potential deviations from neutrality. We used the program DNASP version 5.0 [69] and significance was determined from 10,000 coalescent simulations. For each Proechimys species, we analyzed the complete dataset, the clusters from BAPS analyses, and samples separated by main interfluves (Maroni–Oiapoque; Oiapoque–Araguari; Araguari–Jari; Jari–Trombetas). We used the Mismatch Distribution Analysis [70] with the sum of square deviations (SSD) between the observed and the expected mismatch and the ragggedness index as statistic tests and generated a Bayesian Skyline Plot (BSP) in BEAST 1.6.2 [71] to assess the demographic history of the clusters for both species over time. We only inferred BSP for clusters with more than 10 sequences. We used HKY+G evolutionary model and Yule prior. We ran 10⁷ MCMC generations, excluding the first 25% samples as burn-in. Effective sample sizes, parameter estimates, and convergence of the Markov chains were determined using TRACER 1.6 [62].

Phylogenetic relationships in P. guyannensis and P. cuvieri

The phylogenetic reconstruction presented similar results in BI and ML analyses (Fig 2, see too S2 Appendix and S3 Appendix). Our data for the eastern portion of the Guianas region also formed monophyletic clades for both P. guyannensis (clade B) and P. cuvieri (clade K). The phylogenetic reconstruction identified P. guyannensis clades distributed between the north and southeast of the study area (clade C, Fig 2). The samples from Amazonas, Roraima, and Venezuela formed clade J. In clade C, the samples were organized in four subclades: the haplotypes from French Guiana and northern Amapá (clade D), central Amapá (clade I), south-east Amapá (clade H) and two localities in Pará (clade G). These clades were well-supported statistically either with bootstraps (ML) or posterior probabilities (Bayesian analysis), with the exception of clade D, which was supported only weakly by the ML analysis (Fig 2). Five clades were evidenced in P. cuvieri (clade N, Fig 2), two in northern French Guiana (clades P and R), one that groups haplotypes between central northern French Guiana to the Araguari River in Amapá (clade S), south-east Amapá and Pará, including the sample from Baramita, Guyana (clade V) and haplotypes from the Jari River, associated with those from Santana Island, in the Amazon Estuary (clade U). These clades were supported statistically by the PP or bootstrap analysis (Fig 2).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 2. Molecular dating chronogram using cytochrome b with a HKY + I + G model (ML and BI topology presents similar topology).

Statistical support higher than 0.90 posterior probability (PP) and 70 bootstraps is identified next to the nodes divided in two sides, the left side represents PP values and the right side the bootstrap values; black represent 1 PP and 100 bootstrap, very dark gray represent 0,99 PP and 91–99 bootstrap, dark gray represent values between 0.96–0.98 PP and 80–90 bootstrap, and light gray represents values between 0.95–0.96 PP and 70–80 bootstrap. The letters B, C, D, E, F, G, H, J and I identify clades formed in P. guyannensis. The letters K, L, M, N, O, P, Q, R, S, T, U, V and X identify clades formed in P. cuvieri. Each terminal is identified for sample followed by country of origin. The numbers in the parentheses identified the locality in the Fig 1. We also identified the state of origin in the case of Brazil. Clusters recovered by BAPS are identified after the terminals for P. guyannensis (G1, G2, G3, G4) and for P. cuvieri (C1, C2, C3, C4). Inferred spliting dates are shown on the nodes.

https://doi.org/10.1371/journal.pone.0206660.g002

We infer that P. guyannensis and P. cuvieri diverged in the Pliocene, around 5.03 million years ago (Mya), with 95% confidence of Highest Posterior Density Interval (HPDI) ranging between 3.8–6.2 Mya. The diversification of the clades found in the eastern portion of the Guiana region occurred in late Pliocene, at around 2.85 Mya (95% HPDI: 1.9–3.8 Mya) in P. guyannensis and more recently in early Pleistocene in P. cuvieri, at 2.19 Mya (95% HPDI: 1.4–3.1 Mya).

Genetic structure of the populations

Both Proechimys species had high haplotype diversity and low nucleotide diversity (Table 1). In P. guyannensis, the 67 sequences analyzed provided 43 haplotypes, 155 polymorphic sites, and 101 informative sites for parsimony analysis. In P. cuvieri, the 49 sequences analyzed returned 45 haplotypes, 140 polymorphic sites, and 94 informative sites. Once analyzed with BAPS, the haplotypes of P. guyannensis were recovered in four clusters (G1, G2, G3 and G4). These clusters correspond to the clades recovered in the phylogenetic tree (Fig 2). The haplotypes in the northern part of the study area, including samples from French Guiana and one location in the state of Amapá on the Anotaie River, a tributary of the Oiapoque, formed cluster G1. The haplotypes from the south-east of the study area, from the localities on the Jari and upper Vila Nova rivers formed group G3. Group G4 brings together the samples from two localities between the Jari and Trombetas rivers (Fig 2). The haplotype network recovered the clusters identified by the BAPS, resulting in four clusters separated from each other by more than 20 mutational steps. The geographical layout encompasses clusters distributed from north to south. Significant haplotype sharing was observed in the Cayenne locality (point 12 in Fig 1) in the cluster G1, fifteen samples shared haplotype 1 (H1, Fig 3). Although this pattern was not upheld in other clusters which were characterized by limited haplotype sharing.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 3.

Map (top) showing the locations of the P. guyannensis, with circle colour show the regions within the study area: North, from the Oiapoque river basin to the Maroni River (light blue); Centre, Araguari River basin (red); Southeast, included Jari and Vila Nova river basins (yellow); and Mamiá and Jatuarana villages, state of Pará (blue). Based on cyt-b the haplotype network (bottom) has branch lengths corresponding to the nucleotide substitutions. The number of mutational steps of the largest branches is shown under the branches. Circle sizes correspond to the number of individuals sharing the same haplotype. Circle colors correspond to the colors of the regions identified on the map. Circles in black correspond to the mean vectors. The clusters recovered by BAPS are outlined and identified as G1, G2, G3 and G4 on the map and haplotype network.

https://doi.org/10.1371/journal.pone.0206660.g003

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 1. Haplotype diversity (HD), nucleotide diversity (π) and the number of individuals (N) in the clusters identified by the BAPS for P. guyannensis (G1, G2, G3, G4) and P. cuvieri (C1, C2, C3, C4), the interfluves of the study region: Maroni-Oiapoque; Oiapoque-Araguari; Araguari-Jari e Jari-Trombetas and the complete dataset for each species.

https://doi.org/10.1371/journal.pone.0206660.t001

The BAPS recovered four P. cuvieri clusters (C1, C2, C3 and C4), matching the clades recovered in the phylogenetic tree (Fig 2). Except the clades P and R (Fig 2) which make up cluster C1 in the BAPS results. The C1 cluster includes the haplotypes from the northern-most area within the study region, between the Maroni and Sinnamary rivers in French Guiana. The haplotypes distributed between central French Guiana and central Amapá, north of the Araguari River, form cluster C2. Haplotypes from the south-east region were recovered in a cluster together with a haplotype from Santana Island, in the Amazon Estuary formed cluster C3. Finally, the C4 cluster linked the samples from the Mamiá, a tributary of the Curuá River, with haplotypes from the Jari River and Baramita, in Guyana (Fig 4). The haplotype network also presented a topology similar to that of the haplotype tree, although the two clades P and R formed a single cluster (C1), in agreement with the BAPS results. The P. cuvieri clusters were separated from each other by fewer mutational steps, between six and 13. There was relatively little sharing of haplotypes, and a large number of unique haplotypes (Fig 4).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 4.

Map (top) showing the locations of P. cuvieri, with circle colour showing the regions within the study area: North, from the Oiapoque River basin to the Maroni River (light blue); Centre, Araguari River basin (red); Southeast, included Jari and Vila Nova River basins (yellow); and Mamiá and Jatuarana villages, state of Pará (blue). Based on cyt-b the haplotype network (bottom) has branch lengths corresponding to the nucleotide substitutions. The number of mutational steps of the largest branches is shown in the figure captions. Circle sizes correspond to the number of individuals sharing the same haplotype. Circle colors correspond to the colors of the regions identified on the map. Circles in black correspond to the mean vectors. The clusters recovered by BAPS are outlined and identified as C1, C2, C3 and C4 on the map and haplotype network.

https://doi.org/10.1371/journal.pone.0206660.g004

We found the largest intraspecific genetic distances for both species when we compared the clusters identified by the BAPS, with p-distance ranging from 3.35% to 4.23% (Table 2) in P. guyannensis. Smaller values (0.84–1.48%) were found when comparing the interfluves (Maroni-Oiapoque, Oiapoque-Araguari, Araguari-Jari, Jari-Trombetas). Similarly, while we found high genetic distance values between interfluvial regions in P. cuvieri (p-distances: 1.2–2.53%), they were still lower than those recorded between clusters (1.84–3.17%; Table 2). Intraspecific genetic distances were smaller in interfluves comparisons than in inter-clusters comparisons. This is well-marked in P. guyannensis (p-values <0.01 in t-tests and Mann-Whitney tests) but less pronounced in P. cuvieri (p-values 0.08 in t-tests; not significant in Mann-Whitney test).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 2. Matrix of the genetic distances (p-distance), based on the sequences of the cytb gene.

We compare values from the clusters defined by the BAPS: P. guyannensis (G1, G2, G3, G4) and P. cuvieri (C1, C2, C3, C4). We compare between the interfluves of the study region: Maroni-Oiapoque; Oiapoque-Araguari; Araguari-Jari e Jari-Trombetas.

https://doi.org/10.1371/journal.pone.0206660.t002

Gene flow between populations

We recorded two-way migrations in all regions analyzed. The largest number of migrants (195.40) in P. guyannensis was recorded from the central region (Araguari basin) to the South region (Jari River, the Amazon and Pará rivers). The second greatest number of migrants (82.60) was recorded between the North region (Oiapoque River and Maroni River), and the South. The lowest values were recorded from the Central region to the North with number of migrants varying between 21.70 and 50.70 (Table 3). In P. cuvieri, the largest numbers of migrants moved from the North to the Central region, and vice versa. The numbers of migrants recorded between the other regions were much smaller, and well balanced around 60.00 (Table 3).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 3. Estimated migration among Proechimys guyannensis and P. cuvieri populations, based on analysis using the Migrate software.

The samples were divided into three regions, the North (from the Oiapoque River to Guyana), Central (Araguari basin in Amapá) and South (Amazon Estuary, Jari and Pará rivers). The plus sign (+) after the name of the region indicates a population receiving immigrants. The highest migration rates are marked in bold type.

https://doi.org/10.1371/journal.pone.0206660.t003

Demographic history

We obtained significant results for the R2 (-2.3781), D (-2.3781), and Fs (-4.572) parameters for the P. guyannensis cluster G1. When we analyzed the whole dataset, however, no significant values were obtained, although a significant R2 value was obtained for cluster G4 (Table 4). We obtained similar results when we analyzed the data organized by interfluve, with significant values being recorded for all the parameters for the Oiapoque-Maroni interfluve (Fs = -4.168; D = -2.3200; R2 = -0.0617), which has haplotypes similar to the G1 cluster.

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Table 4. The values of R2 [68], Fs [66], and D [67] for the clusters recuperated by the BAPS (P. guyannensis: G1, G2, G3, and G4; and P. cuvieri: C1, C2, C3, and C4), and for the principal interfluves in the study area (samples between Oiapoque and Maroni rivers, between Oiapoque and Araguari rivers; Araguari and Jari rivers and Jari and Trombetas rivers), and the complete dataset.

https://doi.org/10.1371/journal.pone.0206660.t004

In the complete P. cuvieri dataset, the only deviation from neutrality was found in Fu’s Fs (Table 4). In the analysis of the clusters, we recorded significant Fs values for clusters C2 and C4, and significant R2 values for clusters C1 and C3. In the analysis of the interfluves, we recorded significantly R2 and Fs values for the Araguari-Jari interfluve (Table 4). The Mismatch distribution suggest a significant expansion for G1 only (both SSD and Raggedness tests p<0.05), and to a less extent for G2 (Raggedness p<0.05, SSD p value = 0.08) and G4 (SSD p< 0.01). Bayesian reconstruction of the demographic history revealed an expansion event for both species. In P. guyannensis, the most recent expansion was recorded for cluster G1 and important expansion (revealed by all statistical tests), occurring over the past 30,000 years, whereas in cluster G2, expansion was either quite older and/or of a lower extent (Fig 5). Signals of expansion population were also found in two clusters analyzed of the P. cuvieri, although they benefited from a lower statistical support (Fig 5).

Download:

• PPT
  PowerPoint slide
• PNG
  larger image
• TIFF
  original image

Fig 5.

Bayesian Skyline Plots showing the demographic evolution for the P. guyannensis, G1 (A) and G2 (B) and P. cuvieri C2 (C) and C3 (D).

https://doi.org/10.1371/journal.pone.0206660.g005

Diversification of P. guyannensis and P. cuvieri in northeastern Amazonia

Patton et al. [10] assessed the genetic divergence in P. cuvieri samples from Brazil (Amazonas and Pará), together with localities in Peru, Venezuela and French Guiana, using sequences of the cytb gene and confirmed the monophyly for eastern Amazonia, including Guianas region clade. Our phylogenetics analyses recovered a clade of P. cuvieri in northeastern of the Amazon region which appeared distinct from samples of the other Amazon regions. For P. guyannensis the haplotypes of northeastern region were recovered in a clade separated from the other samples of the distribution area of this species (Amazonas, Roraima and Venezuela). Based on the analysis of cytb and of karyotypes for P. guyannensis samples, Bonvicino et al. [28] recovered three clades in the Guiana region including one north-eastern clade (French Guiana) and two western clades (Roraima and Amazonas/Venezuela).

We inferred a late Pliocene origin for the diversification of the clade from the eastern Guianas region for P. guyannensis at approximately 2.8 Mya and early Pleistocene for P. cuvieri (2.19 Mya). In both species, however, the divergence in four (P. guyannensis) and five (P. cuvieri) lineages occurred during the early Pleistocene (Fig 2). The early Pleistocene is consistent with the intense geomorphological changes occurring in the estuary of the Amazon River, which formed its current course, and the formation of parts of the coast of Amapá [5, 72, 73]. Temperatures also decreased during the Pleistocene [74], leading to the expansion of drier vegetation and the fragmentation of forest blocks [8, 75]. Despite the relatively limited size of our study area probably these phenomena had distinct impacts on its northern and southern portions contributing to the diversification of the populations of both P. cuvieri and P. guyannensis.

In spite of the limitation of mitochondrial data, it has already been shown to be robust [76] and has long been considered as a standard for phylogeography and for assessing genetic variation within and between populations [77]. Further, high ration of copies compared to nuclear genome allows for isolation from long dead tissues and collection-preserved samples, and consequently favor work on large dataset, not requiring new field captures for fresh and/or high quality samples. Here we significantly increased the knowledge about the phylogeography of these two species Amazonian Spiny rats [29, 34] confirming the uniqueness of the western portion of the Guyana region. Our conclusions and interpretations should be confirmed by other molecular characters, complex processes of fragmentation, expansion and admixture—as can be inferred from our study—should be confirmed through the use of several and complementary types of molecular markers.

Genetic variability and demography of P. guyannensis and P. cuvieri

Different hydrological factors, including size and course stability [10], may influence the effectiveness of a river to act as a barrier to gene flow. Previous studies suggested ther is no evidence for a major role of rivers in the genetic differentiation of most Amazonian rodents [4, 10]. For example, the spiny rats Proechimys steerei and P. simonsi share haplotypes on both margins of the Rio Juruá which a highly meandering river with intense seasonal flooding in western Amazonia [10, 16]. Similarly our results suggest that the Oiapoque, Araguari and Jari rivers have not acted as barriers to the dispersal of P. guyannensis and P. cuvieri, since haplotypes were shared between both margins of these rivers and included in the same genetic cluster (Figs 3 and 4). In addition, genetic distances for P. guyannensis were significantly smaller in between opposite margins than between clusters (Table 2). Moreover the permeability of riverine barriers is highly dependent on species-specific traits [78]. For example, the Oiapoque river seems to be a more efficient barrier for amphibian species that are leaf-litter dwellers or lack free-living tadpoles [79].

However our results indicate that the northeastern part of the Guianas region seems to have undergone distinct evolutionary events. We found clear genetic structuring in our sample of P. guyannensis, with an association between the haplotype network and the geographic distribution of the four clusters identified as G1, G2, G3 and G4. The diversity of P. guyannensis has been established through karyological [28, 30, 33], molecular [28, 34] and morphological studies [36, 80]. French Guiana is located entirely within the area of the Guianan Shield, this region is an ancient region and previous studies regarded it like a forest refugium of northeastern Amazonian [81]. Van Vuuren et al. [34] suggested a recent expansion based in the high number of shared haplotypes in P. guyannensis from the regions of Cayenne and Petit Saut (see Fig 1), in French Guiana. With the increase of the sampled area in northeastern Amazonia, we identified a cluster (G1) that groups samples from northern of French Guiana to Oiapoque River, with a large number of shared haplotypes. This cluster includes a haplotype (H1) that is shared with the majority of the samples from Cayenne, with one principal haplotype (H35) that forms the center of a star-like network, which may reflect recent colonization and/or are expanding population [34, 82].

In P. cuvieri we recorded also four clades, of which cluster C1 occurs in the north of our study area, in the region between the Maroni and Mana rivers in northern French Guiana. We estimated the diversification at approximately 1.76 Mya for cluster C1 in P. cuvieri and 1.68 Mya for G1 in P. guyannensis (Fig 2). These two datings, almost similar, might be the result of a process of genetic isolation in forest fragments produced by the intense changes in the configuration of the region’s vegetation that occurred during the early Pleistocene and late Pliocene [8, 83]. Cold and warm periods of the Pleistocene, isolating populations and promoting backward and forward migrations, may also have led to competitive interactions between the two species, and facilitated structuration, in a dynamic process comparable to the disturbance-vicariance hypothesis [84].

We also detected a signal of recent expansion indicated by the high degree of haplotype sharing (Fig 3) and demographic changes (Fig 5 and Table 3). The most intense demographic expansion in both P. cuvieri and P. guyannensis occurred over the past 20,000 years (Fig 5), when the climate became more favorable for mesic environments after the Last Glacial Maximum (LGM) (20,000 yr BP) [85, 86]. While there is no clear consensus on the extent to which the region’s forests were reduced or fragmented during the LGM [86], it does seem likely that the forests shrank as temperatures decreased [74, 87]. In north-eastern Amazonia, open habitats, such as dry forests and savannahs, probably expanded during the LGM [88, 89]. It seems likely that the spiny rats, P. cuvieri and P. guyannensis, benefitted from the subsequent increase in temperatures and the expansion of humid forests, as shown with their demographic expansion (Fig 5).