Genetic differentiation in geographically close populations of the water rat Nectomys squamipes ( Rodentia , Sigmodontinae ) from the Brazilian Atlantic Forest

We examined the genetic structure and the effects of a bottleneck in populations of the water rat Nectomys squamipes, a primary host of Schistosoma mansoni. Eight microsatellite loci were studied in 7 populations from the Sumidouro region of the Brazilian state of Rio de Janeiro. Our data, covering a four-year period during which a bottleneck occurred, revealed substantial variation (6-31 alleles per locus) and high levels of both observed (0.718-0.789) and expected (0.748-0.832) heterozygosity. Most populations were in Hardy-Weinberg equilibrium without linkage disequilibrium between loci. Overall average genetic differentiation between populations (estimated with the FST (θ) and RST (ρ) analogues was 0.037 for θ and 0.060 for ρ. There was significant allelic and genotypic differentiation between populations, especially in pairwise comparisons that included the most geographically isolated population. Direct migration estimates showed a low rate of migration, indicating that infected N. squamipes populations had a limited ability to spread S. mansoni. When the preand post-bottleneck populations were compared there was no detectable reduction in heterozygosity or allele number, although a significant excess of heterozygosity was detected in the post-bottleneck population.


Introduction
Water rats (Nectomys squamipes) are semi-aquatic rodents widely distributed along rivers and streams, including those located near the large urban areas of eastern Brazil.These water rats are primary hosts for Schistosoma mansoni (Rey, 1993) and although their fitness is not reduced by infection they are unable to eliminate the parasite (D'Andrea et al., 2000).Dispersion of N. squamipes populations might spread S. mansoni wherever the intermediate host (Biomphalaria sp.) is present, a possibility aggravated by the fact that N. squamipes has a wide geographic distribution which encompasses most areas of Brazil where schistosomiasis is endemic.
The study presented in this paper used eight highly polymorphic microsatellite markers to assess gene flow between seven N. squamipes populations (separated by dis-tances of between 500 m and 6,000 m) inhabiting the Sumidouro region of the Brazilian state of Rio de Janeiro where schistosomiasis is prevalent.Because it has been shown (Gentile et al., 2000) that N. squamipes populations in the study area fluctuate and there are periods when specific populations reach almost local extinction, we investigated the populations both before and after bottlenecking to assess its effects on the genetic diversity of N. squamipes.

Study area
Animals were collected in two adjacent valleys (Pamparrão and Porteira Verde) and in Bela Joana in the Sumidouro region (22°02' S, 40°41' W) of the Brazilian state of Rio de Janeiro (Figure 1).Both Pamparrão and Porteira Verde valleys form part of a common hydrographic system with several small, permanent streams, irrigation channels and flooded areas, the landscape being composed of small rural establishments with plantations, pasture lands and fragments of Atlantic Forest.Distances between capture sites along rivers ranged from 500 to 6,000 m (Table 1).The capture sites in the Pamparrão valley were Sr. Célio (PC), Ponte (PP) and Braço (PB) while the site in the Porteira Verde valley was called Escola (PVE).Two other collecting sites named Banqueta (PVB) and Charco (PVC) were located along an irrigation channel connecting the two valleys.The most isolated site was Bela Joana (BJ), which is located outside these valleys.

Field methods
Capture lines (localities) were established along watercourses and flooded areas with between five and ten wire live traps being set every 13 m.Trapping sessions were conducted for five nights every other month during the mark and recapture phase (September 1995 to November 1996) and the removal phase of the study (April 1997to April 2000), a total of 425 traps being laid per night for each trapping session.Detailed descriptions of this field study have been reported elsewhere (D'Andrea et al., 2000;Gentile et al., 2000).About 1 mL of blood was taken from the hind foot of 23 N. squamipes specimens captured during the mark and recapture phase of the study.Additionally, liver samples were obtained from 139 specimens captured during the removal phase of the study.

DNA isolation and microsatellite analysis
Genomic DNA was extracted from blood samples or liver tissue by the standard proteinase-K/phenolchloroform procedure (Sambrook et al.,1989).Samples were preserved in TE buffer (10 mM Tris, 1 mM EDTA, pH = 7.2) at an approximate concentration of 15 µg/mL.
Fragment analyses were conducted with an automatic DNA sequencer ABI PRISM 377, with standard loading and electrophoresis conditions following a 3 h run.Alleles were sized relative to an internal size standard (ROX GS 500; Applied Biosystems) and analyzed with GENESCAN 2.1 (Applied Biosystems).Size estimates, in base pairs (bp), were rounded to integers differing by 2 bp.

Statistical analyses
Spatial analysis was performed by pooling data from individuals from the same locality.Inter-population differentiation between all population pairs was tested using Fisher's RxC test (Sokal and Rohlf, 1994), a procedure which analyses each locus to determine differences in allele frequencies providing that the populations under investigation are in Hardy-Weinberg (HW) equilibrium and that all alleles are assumed to be independent.Fisher's combined probability test was used as a global test for all loci to determine overall significance.These tests were performed with version 1.3 of the Tools for Population Genetic Analysis (TFPGA) program (Miller, 1997).
Genotypic differentiation between different populations and over all the populations as a whole was tested with the log-likelihood G statistic using genotype permutations rather than allele permutations, a test that does not require   populations to be in HW equilibrium (Goudet et al., 1996).This test was performed with version 3.2a of the GE-NEPOP program (Raymond and Rousset, 1995).Wright's F ST analogue (θ; see Weir and Cockerham, 1984) was used to assess population differentiation, calculations being made using GENEPOP 3.2a.Multilocus values of θ across all populations were calculated with TFPGA 1.3.This estimate relies on the infinite allele model (IAM) which may not be the most appropriate model for analyzing microsatellite data.For this reason, we also calculated population differentiation using the R ST estimator (ρ; see Slatkin, 1995) which is a more accurate estimate when mutations follow the stepwise mutation model (SMM), as may be the case with microsatellites (Valdes et al., 1993).The short time scale of this study may prevent significant departures from the SMM (Goldstein et al., 1995), such deviations normally being caused through factors such as lack of constraints on allele size (Bowcock et al., 1994) and the occurrence of mutations which are independent of allele size (Amos et al., 1996;Levinson and Gut, 1987).Estimates of ρ were calculated with the RstCalc program (Goodman 1997) after standardizing allele length to eliminate bias and with 5,000 bootstrap replicates for calculating 95% confidence intervals, standardization for each locus being performed by dividing allele length by the overall standard deviation of the repeat length (Goodman, 1997).
Isolation by distance was tested with the isolde option of GENEPOP 3.2a, differentiation statistics being linearized based on the stepping stone model (Rousset, 1997) using diff / (1-diff), where diff stands for either θ or ρ values.Geographic distances were not transformed to a logarithmic scale because dispersion patterns of N. squamipes were uni-dimensional (along rivers).The Mantel test (Mantel, 1967) was used to test for a significant relationship between geographic distance and pairwise genetic differentiation.
In spatial analyses, genetic diversity in each population was estimated from the number of alleles per locus (A), heterozygote direct count per locus (H DC ) and mean expected heterozygosity per locus (H E ) under HW equilibrium (Hartl and Clark, 1989).Estimates of expected and observed heterozygosity were calculated with nonparametric statistics using TFPGA 1.3 because most variables did not follow a normal distribution.Significant H E differences between more than two populations were calculated using the Friedman ANOVA program (Zar, 1996).
Linkage disequilibrium between loci was tested with Fisher's exact test using the Markov chain analysis (3,000 de-memorizations; 450,000 iterations) under a null hypothesis postulating lack of genotype association between different loci, the calculations being performed with GENEPOP 3.2a.Each population was tested separately to avoid distortions resulting from the Wahlund effect.
Deviations from HW equilibrium for genotype proportions at each locus were tested by a Markov Chain ap-proximation (3,000 de-memorizations; 450,000 iterations) of unbiased exact P-values (Guo and Thompson, 1992).A global test across loci was performed using Fisher's method, significance levels being calculated per locus, per population for all combined loci and populations (global test) using GENEPOP 3.2a.Wright's F IS analogue (f ; see Weir and Cockerham, 1984) was estimated with TFPGA, 95% confidence intervals being estimated by bootstrapping with 5,000 replications.
Temporal analysis was carried out by pooling animals from different populations with no significant allelic and genotypic differences captured during the same period.For this analysis we used data from a capture period (April to July 1997) before a bottleneck (November 1997 to August 1998, when there were only 2 or 3 captures per collection session) and after this bottleneck (November 1998 to August 1999).
To test for effects on microsatellite diversity (Spencer et al., 2000), we determined A and H DC in pre-and postbottleneck populations and used the Wilcoxon signed-rank test to compare A and H DC values between capture periods; A was compared using the same sample size for both populations (18 animals in each).To detect bottlenecks we used the method described by Cornuet and Luikart (1996) which is based on differences between observed and expected heterozygosity according to allele number in a given population (Nei, 1978), this method relying on the fact that after a bottleneck the number of alleles decreases faster than the expected heterozygosity (Cornuet and Luikart ,1996;Nei et al.,1975).As a negative control, we also tested for heterozygote excess in the pre-bottleneck sample by using a two-phase-model (TPM; 5,000 iterations) with 90% of mutations following the SMM using the BOTTLENECK program (Piry et al., 1999), TPM being an intermediate mutation model between IAM and SMM which is especially recommended for the study of microsatellite loci evolution (Rienzo et al., 1994).A test of significance was performed using the Wilcoxon test option, which provides relatively high power and may be applied with few polymorphic loci and any number of individuals.
To test for multiple sires and as accuracy controls for the microsatellite scores, the microsatellite genotypes of 14 captive-born animals (from 6 captured pregnant females) were compared to the microsatellite genotypes of their mothers.These offspring were excluded from the general population analyses to avoid over-consideration of related alleles.

Differentiation of populations
The number of captured animals in each population is shown in Table 2, which also shows genetic data.We found significant allelic differentiation for all loci and populations (p < 0.001), 15 of 21 population pairs being signifi-cantly different (p < 0.001).When we applied the Bonferroni correction, however, differentiation was not significant for 6 population pairs: PVB/PVE (p = 0.025), PVB/PP (p = 0.027), PVE/PP (p = 0.291), PVE/PC (p = 0.009), PVC/PP (p = 0.089) and PVC/PC (p = 0.013).Analysis of each microsatellite locus showed that Nec15, followed by Nec18 and Nec19, contributed to most of the differentiation between populations while the smallest contribution was due to Nec08 and Nec28.
Significant (p < 0.001) genotypic differentiation was observed for all loci and populations.Each locus showed genotypic differentiation for all populations with the G-likelihood test, with p values ranging from 0.0002 for Nec15 and 0.0000 for all other loci.Genotypic differentia-406 Maroja et al.  (Weir and Cockerham, 1984).Significant single-locus and multi-locus f-values determined using exact P-values calculated by the Markov Chain method with the Bonferroni correction (nominal α = 0.05) are shown in bold type.The significance of the global test across loci was determined using Fisher's method.

Linkage disequilibrium, HW-equilibrium and genetic diversity
All loci were highly polymorphic in all populations, with the mean allele number varying between 11.25 ± 3.88 in the population with the largest sample size (PP; n = 33) to 6.13 ± 1.13 for the population with the second smallest sample size (PVC; n = 8; Table 2).The overall mean allele number was 16.87 ± 7.51, while Nec23 showed the highest allele number (31) followed by Nec15 (24).Mean H DC values varied between 0.789 ± 0.11 for the PC population and 0.718 ± 0.07 for the PP population.Mean H E values varied between 0.832 ± 0.06 for the PP population and 0.748 ± 0.07 for the PVC population.The mean frequency of private alleles over all loci and populations was 0.038, with population BJ having the highest number of private alleles (nine, three of which were in Nec15), and population PVE (which had the lowest sample size) being the only population having no private alleles.Nec15 and Nec23 showed the highest number of private alleles (nine and four, respectively).
When we analyzed the populations separately we found no statistically significant (p > 0.03 with the Bonferroni correction) genotypic linkage disequilibrium but the global test across all populations was significant (p = 0.001) for Nec12 and Nec14, probably because of the pooling of individuals from different populations .When the global test was performed there was a significant (p < 0.001) deviation from HW-equilibrium, although the lower F IS margin was not very high (F IS = 0.1023; 95% CI = 0.0379 -0.1783), probably due to the Wahlund effect.With the Bonferroni correction, HW-equilibrium across loci was significantly (p < 0.001) different from zero for Nec15, Nec19 and Nec23.Nec15 was in disequilibrium in the PVB, PC and BJ populations (p < 0.001), Nec19 was in disequilibrium only in the PP population (p < 0.001) and Nec23 was in HW-equilibrium for each individual population (p > 0.02).When we excluded these three loci from the global test, deviation from HW equilibrium was considerably lower (p < 0.04, F IS = 0.034; 95% CI = 0.067-0.0001).

Effects of the Pamparrão valley population bottleneck
Analysis of the pooled data from the PP, PVB, PVE and PC capture sites did not show any significant difference in the number of alleles per locus or in observed heterozygosity between pre-and post-bottleneck populations.However, heterozygosity excess was significantly higher than zero (p = 0.01) in the post-bottleneck population but was not significant in the pre-bottleneck population (p = 0.84).

Analysis of offspring and mothers
Analysis of 6 litters consisting of 14 captive born offspring showed no inconsistencies in microsatellite allele-scoring.There was no evidence of multiple paternity because we did not find more than 4 alleles per locus among the offspring of each litter.Although the number of families was too small to guarantee that offspring of the same litter were always single-sired, it is likely that copulation plugs observed in recently mated females (L.S. Maroja, field observation) probably hindered insemination by more than one male.

Discussion
The eight microsatellite loci investigated were all highly polymorphic, with a mean observed heterozygosity of 0.747 and mean allele number of 16.9 within populations.These results agree with a study of more geographically separate N. squamipes populations in which the same microsatellite loci showed a mean heterozygosity of 0.83 (Almeida, 2000).This estimate was, however, higher than that reported for small mammals such as the semi-aquatic rodent Arvicola terrestris, described by Steward et al. (1999) as having a mean microsatellite heterozygosity of 0.63, and Mus musculus populations where the mean microsatellite heterozygosity has been calculated to be ≈ 0.61 (Dallas et al., 1995).
Geographically close populations of N. squamipes showed significant differentiation despite the short dis- tances separating them, in agreement with previous studies of small rodents which have shown a significant population structure over small spatial scales (Selander, 1970;Dallas et al., 1995).This significant population differentiation indicates that migration is probably restricted (resulting in limited gene flow), this view being supported by ecological studies showing that N. squamipes dispersion is restricted and that most movement is confined to riverbanks.This finding is relevant in view of the fact that N. squamipes is a reservoir of S. mansoni (Ernest and Mares, 1986) and that its extensive distribution covers most areas where schistosomiasis is endemic (Ribeiro et al., 1998) and where some rodent populations show infection rates higher than 90% (Veiga-Borgeaud et al., 1986, 1987).One of the authors (L.S. Maroja) has noted that direct capture-recapture migration estimates covering the same study area show a low migration rate for N. squamipes, with only 10 (mostly males) out of 162 recaptured animals had moved from their original population (mainly from PVB to PVC) during a six-year period and it is probable that not all of them attained reproductive success.This migration pattern suggests/indicates that N. squamipes may not pose a serious threat of spreading S. mansoni on a large scale.
To assess population differentiation we used both θ as an estimate of F ST under the infinite allele model (IAM) and ρ as an estimate of R ST under the stepwise mutation model (SMM), two estimates being used because both of these mutation models are controversial when applied to microsatellite loci.Estoup et al. (1995) found IAM to be more suitable (although this probably applied to a type of microsatellite repeat motif) whereas other studies have indicated that the SMM is the most adequate model (Valdes et al., 1993), especially when considering large genetic distances (although this was not the case in our study).Although our estimates of ρ were usually higher than θ (mean ρ/θ ratio = 3.40 ± 5.52), both estimators generally agreed and populations with significant genotypic differentiation also showed significant ρ values.Although R ST and F ST values may be expected to vary widely between loci for genealogical reasons (Slatkin and Barton, 1989), a downward bias in mean F ST values (θ) is to be expected in populations that evolved independently by both drift and mutation.This is because F ST does not consider forward and backward mutation in allele size and thus underestimates differentiation between populations (Slatkin, 1995) but also because microsatellites evolve at a much higher mutation rate (Dallas, 1992) which results in a higher number of convergent mutations (homoplasies).Many other studies of natural populations have also reported a downward F ST bias in relation to R ST (Ciofi and Bruford, 1999;Reush et al., 2000;Shaw et al., 1999).
A significant bottleneck was found using the methods of Luikart and Cornuet (1998), despite the fact that allelic number and heterozygosity estimates did not differ be-tween pre-and post-bottleneck populations as was the case for similar studies reported in the literature (Daley, 1992;Ardern et al., 1997).Since the amount of loss of heterozygosity depends on the time a population spends in the bottleneck (Nei et al., 1975), our findings might reflect the fact that the bottleneck encompassed only a few generations and, consequently, did not result in any considerable loss of alleles or of heterozygosity.Three alternative hypotheses might explain these results: firstly, it may have been that alleles were not actually lost and/or any allele loss could not be detected in our small samples; secondly, allele number and heterozygosity could have been re-established by immigration, although this is unlikely to have happened because the populations were analyzed immediately after the bottleneck; thirdly, the bottleneck may not have been as drastic as was indicated by our trap records.The Luikard and Cornuet analysis is preferable to other methods because it does not require genetic data from pre-bottleneck populations (Luikard and Cornuet, 1998), although it is restricted to analyzing relatively recent bottlenecks as was the case for our study.Because of the small number of animals captured we pooled more than one generation, but it appears that this did not affect the results because no excess heterozygosity was detected in the pre-bottleneck population.We may therefore conclude that the genetic diversity of Nectomys squamipes populations remains high despite undergoing regular bottlenecks.Understanding of this process might be valuable for conservation genetics of endangered species that are loosing genetic variability due to a reduction in effective population size.
Effective population size has been found to increase as a result of multiple rather than single paternity (Sugg and Chesser, 1994), multiple paternity having already been reported in several taxa (Avise, 1994) including rodents (Baker et al., 1999).We found no evidence of multiple paternity in the litters studied by us and since only one or two paternal alleles were detected in these litters they were probably sired by a single male.The fact that insemination by other males may have been hindered has been supported by our field observations which have detected copulatory plugs in recently inseminated N. squamipes females.

Table 3 -
Geographically close populations of the water rat Nectomys squamipes 407 Estimates of differentiation between Nectomys squamipes populations.Values above the dashed diagonal are for the F ST estimator (θ) and those below the diagonal the R ST estimator (ρ).