# Rare interspecific hybridisation between sympatric African four-striped mice, and indications of fine-scale intraspecific spatial structure related to social factors

Guila Ganem, Claire M.S. Dufour, Nico L. Avenant, Karine Berthier, Anne

Loiseau, Pierre Caminade, Neville Pillay, Carine Brouat

► To cite this version:

Guila Ganem, Claire M.S. Dufour, Nico L. Avenant, Karine Berthier, Anne Loiseau, et al. Rare interspecific hybridisation between sympatric African four-striped mice, and indications of fine-scale intraspecific spatial structure related to social factors. Journal of Vertebrate Biology, 2026, 75, pp.25073. 10.25225/jvb.25073. hal-05500308

HAL Id: hal-05500308

https://hal.inrae.fr/hal-05500308v1

Submitted on 9 Feb 2026

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
<figure><img src="image_2.png" /></figure>

# Rare interspecific hybridisation between sympatric African four-striped mice, and indications of fine-scale intraspecific spatial structure related to social factors

Authors: Ganem, Guila, Dufour, Claire M.S., Avenant, Nico L., Berthier, Karine, Loiseau, Anne, et al.

Source: Journal of Vertebrate Biology, 75(25073)

Published By: Institute of Vertebrate Biology, Czech Academy of Sciences

URL: https://doi.org/10.25225/jvb.25073

The BioOne Digital Library (https://bioone.org/) provides worldwide distribution for more than 580 journals and eBooks from BioOne's community of over 150 nonprofit societies, research institutions, and university presses in the biological, ecological, and environmental sciences. The BioOne Digital Library encompasses the flagship aggregation BioOne Complete (https://bioone.org/subscribe), the BioOne Complete Archive (https://bioone.org/archive), and the BioOne eBooks program offerings ESA eBook Collection (https://bioone.org/esa-ebooks) and CSIRO Publishing BioSelect Collection (https://bioone.org/csiro-ebooks).

Your use of this PDF, the BioOne Digital Library, and all posted and associated content indicates your acceptance of BioOne's Terms of Use, available at www.bioone.org/terms-of-use.

Usage of BioOne Digital Library content is strictly limited to personal, educational, and non-commercial use. Commercial inquiries or rights and permissions requests should be directed to the individual publisher as copyright holder.

BioOne is an innovative nonprofit that sees sustainable scholarly publishing as an inherently collaborative enterprise connecting authors, nonprofit publishers, academic institutions, research libraries, and research funders in the common goal of maximizing access to critical research.
# Rare interspecific hybridisation between sympatric African four-striped mice, and indications of fine-scale intraspecific spatial structure related to social factors

Guila GANEM<fnref n="1" />,<fnref n="2" />* ID, Claire M.S. DUFOUR<fnref n="3" /> ID, Nico L. AVENANT<fnref n="4" /> ID, Karine BERTHIER<fnref n="5" /> ID,
Anne LOISEAU<fnref n="6" /> ID, Pierre CAMINADE<fnref n="1" /> ID, Neville PILLAY<fnref n="2" /> ID and Carine BROUAT<fnref n="6" />* ID

► Received 1 August 2025; Accepted 11 December 2025; Published online 30 January 2026

**Abstract.** Spatial and social population structures can influence intra- and interspecific genetic processes. We investigated the spatial and social genetic structure of two sister species and evaluated whether there was evidence of hybridisation between them at a microgeographic scale, in syntopy. Our study models were *Rhabdomys bechuanae* (Rbech) and *Rhabdomys dilectus dilectus* (Rdil). We generated population genetic data using microsatellite markers and cross-referenced some of these results with published behavioural data from a prior study, reanalysed here. We compared genetic diversity and gene flow both in sympatry and allopatry, we tested for hybridisation between the sister species in sympatry and assessed the relationship between kinship and fine-scale spatial and social organisation. Genetic diversity was high within all populations of the two species. Population genetic substructure was more pronounced in Rbech than in Rdil, consistent with stronger group cohesion reported for Rbech. STRUCTURE and NEWHYBRIDS analyses suggested that sympatric individuals, 5% in Rbech, and 2% in Rdil, might have hybrid ancestry. Our results indicate social fences may limit gene flow within and between species. Further, limited hybridisation suggests that, despite syntopy, only rare events of hybridisation may occur, possibly at phases of density troughs and/or that F1 hybrids have a relatively low fitness.

**Key words:** population genetics, microsatellites, species interference, social grouping, kinship, *Rhabdomys bechuanae*, *Rhabdomys dilectus dilectus*, Southern Africa

<footnote n="1">¹ University of Montpellier, ISEM, CNRS, IRD, Montpellier, France; e-mail: guila.ganem@umontpellier.fr, pierre.caminade@umontpellier.fr</footnote>

<footnote n="2">² University of the Witwatersrand, School of Animal, Plant and Environmental Science, Wits, Gauteng, South Africa; e-mail: Neville.Pillay@wits.ac.za</footnote>

<footnote n="3">³ Université Bourgogne Europe, CNRS, Biogéosciences, UMR, Dijon, France; e-mail: claire.dufour@u-bourgogne.fr</footnote>

<footnote n="4">⁴ University of the Free State, Department of Mammalogy, National Museum and Centre for Environmental Management, Bloemfontein, Free State, South Africa; e-mail: navenant@nasmus.co.za</footnote>

<footnote n="5">⁵ INRAE, Pathologie Végétale, Montfavet, France; e-mail: karine.berthier@inrae.fr</footnote>

<footnote n="6">⁶ University of Montpellier, CBGP, IRD, CIRAD, INRAE, Institut Agro, Montpellier, France; e-mail: carine.brouat@ird.fr, anne.loiseau@inrae.fr</footnote>
<figure><img src="image_3.png" /></figure>

Introduction

The spatial and social structures of species may be inferred from behaviour and ecological studies, e.g. home range size, home range overlaps and nest site (NS) sharing (Schradin et al. 2012, Califf et al. 2020, Wauters et al. 2021). Spatial structure and social interactions are reciprocally related, and social organisation may also be shaped or reflect inter-specific competition (e.g. Dufour et al. 2015a, 2019). Population genetics can reveal how intra- and interspecific processes are linked to spatial and social structures by enabling estimates of kinship among individuals, as well as gene flow within and between species (Tung et al. 2008, Morgan et al. 2022). In this study, we investigated how social factors and inter-species competition influences the genetic structure of two species, and evaluated the likelihood of hybridisation between them in sympatry at a microgeographic scale.

Our focal species were two African rodents, *Rhabdomys bechuanae* (hereafter Rbech) and *R. dilectus* (Rdil), whose distributions overlap in central South Africa. Previous studies have mapped the boundaries of these ranges and identified zones of sympatry (Ganem et al. 2012, 2020). Few studies have addressed the population genetics of *Rhabdomys*. The first used allozyme markers in 23 populations across the range of what was then considered a monospecific genus, identifying geographical subgroupings (Mahida et al. 1999). Some of these subgroupings were later confirmed, and geographically distinct mitochondrial lineages were described within the genus (Rambau et al. 2003, du Toit et al. 2012). Another study investigated the impact of climatic variation on population genetic structure in seven populations distributed along a transect from Namibia to South Africa, using microsatellites and MHC genes (Froeschke & Sommer 2014). However, the authors treated populations of two species, *R. bechuanae* (in the more arid Namibian sites) and *R. pumilio* (in the less arid South African sites) as a single taxon. The only microsatellite study that avoided mixing species focused on another species, *R. pumilio*, at a fine spatial scale, testing for spatial structure and sex-biased dispersal (Solmsen et al. 2011).

Phylogenies based on *cytb* and cytochrome oxidase I (COI) showed that Rbech and Rdil diverged approximately 4 million years ago (Rambau et al. 2003, du Toit et al. 2012). They likely evolved in allopatry under distinct environmental conditions: Rbech in arid regions and Rdil in mesic ones (du

Toit et al. 2012). Sympatric populations occurred at the limits of the distribution ranges of the two species, patchily distributed among parapatric populations, in the semi-arid central region of South Africa (Ganem et al. 2020). Both species are diurnal and shelter overnight in nocturnal NS. Rbech has larger home ranges and exploits open habitats more exposed to predation than Rdil (Dufour et al. 2015a, 2019). In parapatry, Rdil is mostly found in riverine vegetation, with dense and more or less continuous cover, while Rbech mostly occupies patchy open scrubland type habitats (Dufour et al. 2015a). In sympatry, Rdil is found in more open microhabitats than in parapatry, although with more cover than those occupied by Rbech (Dufour et al. 2015a, 2019). Finally, Rbech is characterised by more cohesive and stable social groups than Rdil (Dufour et al. 2019), differences that are expected to be reflected in their genetic structures.

A comparison between sympatric and parapatric populations of the two species strongly suggested potential interspecific competition and character displacement between Rdil and Rbech (Dufour et al. 2015a, 2019). Furthermore, smaller-than-expected home range overlaps between the two species in sympatry strongly suggest that they avoid each other (Dufour et al. 2015a). Moreover, Rdil NS are spatially clustered in sympatry as compared to parapatry (Dufour et al. 2019), also suggesting that they may be impacted by species interference. Finally, laboratory trials indicated that in the absence of conspecific potential mates, adults of the two species could engage in mating behaviour with partners of the other species (Dufour et al. 2015b), suggesting that reproductive interference (Groening & Hochkirch 2008) and hybridisation could occur in sympatry.

We had three central aims: 1) to compare genetic diversity, relatedness, and gene flow within populations of Rbech and Rdil in both sympatry and allopatry, 2) to test for hybridisation between Rbech and Rdil in sympatry by assessing admixture, and 3) to integrate genetic data with previously published information on spatial behaviour and social structure of adult individuals and assess the relationship between genetic relatedness and fine-scale spatial and social organisation in both species (Dufour et al. 2019). To achieve these objectives, we generated population genetic data using microsatellite markers and cross-referenced a sub-sample of these results with social and spatial data from a prior telemetry study involving the same individuals (Dufour et al. 2015a, 2019), which we
<figure><img src="image_4.png" /></figure>

<figure><img src="image_5.png" /><figcaption>**Fig. 1.** Geographical locations of the four *Rhabdomys* populations in South Africa. SAN – Sandveld Nature Reserve, TDR – Tussen Die Riviere Game Reserve, BLO – Bloemhof Nature Reserve, LAJ – Lajuma Research Station. Symbols correspond to *R. dilectus dilectus* (yellow triangles), *R. bechuanae* (blue triangle), both species in sympatry (green circle).</figcaption></figure>

reanalysed for the present study. Given the species-
specific differences in social behaviour (Dufour et
al. 2019), Rbech is expected to exhibit more genetic
discontinuities than Rdil. Furthermore, given their
occurrence in syntopy, we expected that, in the
absence of postzygotic isolation, hybridisation may
occur in sympatry.

**Material and Methods**

The project received ethical clearance (French: C34-
264, South African: AESC 2012/13/2A) and local
trapping permits (Free State Province: permits no.
01/15700 and 01/26960, North West Province: permits
no. 01/11262 and HQ10/14-022NW). For Lajuma
Limpopo, samples were obtained from colleagues
with permission from the Lajuma Research Centre.

Study populations

Samples were collected during seven trapping
sessions in South Africa. The two species occurred in
sympatry in one locality: Sandveld Nature Reserve,
Free-State Province (SAN; S 27°43', E 25°45'; sampled

in November 2011, April 2012, November 2012, April
2013). Allopatric localities were: for Rbech, Tussen
Die Riviere Game Reserve, Free-State Province (TDR;
S 30°28', E 26°09'), sampled in November 2012, April
2013, and November 2013; and for Rdil, Bloemhof
Nature Reserve, Northwest Province (BLO; S 27°38',
E 25°40'), sampled in June 2009, November 2012,
and April 2013. In addition, a more distant allopatric
site was included for Rdil, Lajuma Research Centre,
Limpopo Province (LAJ; S 23°02', E 29°26'), sampled
in April 2011 (Fig. 1).

Mice were captured using Sherman and PVC type
live-traps, baited with a mixture of peanut butter,
oats, coarse salt and sunflower seeds. A piece of
cotton wool was placed in each trap to provide
comfort and insulation. Traps were set every 10 m
along transects traversing favourable habitats for the
two species (details in Dufour et al. 2015a). Trapping
geolocations were recorded using a handheld GPS
(Dakota 10, Garmin International, Kansas, USA,
or similar). A small piece of tail (approximately
0.5 cm) was clipped and preserved in 90% ethanol
<figure><img src="image_6.png" /></figure>

for molecular analysis. The mice were genotyped using *COI* and microsatellite markers (see below). Some of the genotyped adult mice were also used in a radiotracking study during which they were equipped with VHF collars (methods and results published in Dufour et al. 2015a, 2019). In the present study, we extracted a subsample of these published data on home range overlap and NS sharing that involved mice that were also genotyped here (SAN: 11 Rdil, 39 Rbech; BLO: 33 Rdil; TDR: 21 Rbech), and for which we identified 61 NS (see below). These published spatial data were reanalysed to assess the relationship between genetic relatedness and spatial/social grouping in both species. Details of the characteristics of all individuals studied are provided in Table S1.

## Mitochondrial identification and microsatellite genotyping

Species identification was based on *COI* sequencing (as described in Ganem et al. 2012). Total genomic DNA was isolated using a standard salt extraction (Bruford et al. 1992). A portion of *COI* was amplified using published primers (NewF and New Rg) in 20 ml reaction volumes, containing 2 µl 10X Taq

Buffer with MgCl₂ and (NH₄)₂SO₄ (Euromedex), 0.6 µl 10 mM primer, 2 µl 10 µM dNTPs, 0.16 µl GoTaq polymerase (Promega, Madison, WI, USA) and 1 µl of DNA extract (ca. 30 ng/µl). PCR annealing temperature was optimised at 61 °C (30 cycles). PCR products were visualized on 1.5% agarose gels containing 2 µl/100 ml ethidium bromide and cycle sequenced using BigDye terminator chemistry (Applied Biosystems, Foster City, CA, USA), purified with Agencourt AMPure XP purification system (Beckman Coulter Genomics, Grenoble, France), and loaded on an ABI Prism 3130xl 16-capillary genetic analyser (Applied Biosystems, Foster City, CA, USA). All sequencing was performed through the technical facilities of the genotyping and sequencing platform GENSEQ. Sequences were aligned by eye in BioEdit (Hall 1999).

Microsatellite genotyping was conducted for ten *loci* (Rh02, Rh09, Rh21, Rh23, Rh24, Rh28, Rh29, Rh31, Rh36, and Rh46; see Table 1), developed specifically for *Rhabdomys* using a high-throughput microsatellite isolation method through 454 GS-FLX Titanium pyrosequencing of enriched DNA libraries (Malausa et al. 2011) by Genoscreen company (Lille,

<table><caption>**Table 1.** Microsatellite *loci* used in genotyping *Rhabdomys* spp. Repeat motifs, primer sequences, cloned allele size and dye labels.</caption><thead><tr><th>Locus</th><th>Repeat motif</th><th>Primer sequence 5'-3'</th><th>Size of cloned allele (bp)</th><th>dye</th></tr></thead><tbody><tr><td>Rh02</td><td>(GATA)12</td><td>F: CATGCAACATATGGTCTAGCAAA<br>R: TGCTTTTTGGATCCTTGACCT</td><td>106</td><td>6FAM</td></tr><tr><td>Rh09</td><td>(AC)21</td><td>F: TGCATCAAAATTCCAACATGA<br>R: TGTTTGTTGGTGGTAGCGAA</td><td>126</td><td>NED</td></tr><tr><td>Rh21</td><td>(TG)12</td><td>F: CGTGTCCACTGAACTGTGCT<br>R: GCTGGAAGATTCCAAAGCAG</td><td>194</td><td>6FAM</td></tr><tr><td>Rh23</td><td>(GT)18</td><td>F: GTTTGGCTTGTGGTTTGGAT<br>R: GAGAACATTCACTGCCACGA</td><td>204</td><td>VIC</td></tr><tr><td>Rh24</td><td>(CA)21</td><td>F: AGTGTCACATGGGGCTAAGG<br>R: AGGCATCACCTTCTGTGAGC</td><td>221</td><td>NED</td></tr><tr><td>Rh28</td><td>(AC)16</td><td>F: GGGCTTCCTCATTGCTATTTC<br>R: TCCTAGGGAGAGCTTGTGGA</td><td>255</td><td>PET</td></tr><tr><td>Rh29</td><td>(CA)19</td><td>F: TTCCTTATTGACCCCTTCCC<br>R: CCCATTTACCATGACCCAAC</td><td>259</td><td>6FAM</td></tr><tr><td>Rh31</td><td>(AC)18</td><td>F: TGGTAGGGAATGTGCCTTCT<br>R: AGACATGTTCTGATACCACAGAAA</td><td>269</td><td>VIC</td></tr><tr><td>Rh36</td><td>(TGTA)11</td><td>F: CCACAGGCATGTGTATCCTCT<br>R: TGCAGGATGTGCTTGCTTAC</td><td>160</td><td>PET</td></tr><tr><td>Rh46</td><td>(CA)12</td><td>F: TGGCTGCAAATATCACCAAA<br>R: CCCTCCAACAGCCTTGTAAA</td><td>261</td><td>NED</td></tr></tbody></table>
<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

France). Genomic DNA was extracted from ethanol-
stored tissues with the DNeasy Blood & Tissue Kit
(Qiagen). Microsatellite loci were amplified in one
multiplex PCR. PCR amplifications were performed
in a final volume of 10 µL containing 5 µL of the
Qiagen multiplex PCR Master Mix (1X), 2 µM of each
primer, and 1 µL of genomic DNA. The PCR protocol
consisted of the following steps: 1) initial denaturation
at 95 °C for 15 min, 2) 30 denaturation cycles at 94 °C
for 30 s, annealing at 57 °C for 90 s, elongation at 72 °C
for 60 s, and 3) final elongation at 60 °C for 30 min.
PCR products were separated and detected with an
ABI 3130 automated sequencer (Applied Biosystems,
Foster City, CA, USA) and analysed with GeneMapper
v.3.7. Every individual successfully genotyped at some
loci but not at some others underwent re-amplification
once by simple PCR (to avoid primer competition) for
each failed locus.

Population genetic analyses

Throughout this paper, the term population refers
to all individuals of a given species sampled at a
specific locality. Deviations from Hardy-Weinberg
(HW) equilibrium within *loci* and populations and
genotypic linkage disequilibrium between pairs
of *loci* were assessed using GENEPOP v.4 (Rousset
2008). We corrected for multiple testing by using
the false discovery rate (FDR) approach (Benjamini
& Hochberg 1995) implemented in the QVALUE
package of R software (Storey et al. 2021).

For each population, genetic diversity was estimated
using the allelic richness (a_r), and the observed (H_o)
and expected heterozygosities (H_e) calculated with
FSTAT v. 2.9.3 (Goudet 2001). We also analysed the
subpopulation structure by calculating F IS using
FSTAT, and the kinship coefficient (q) of Loiselle
(Loiselle et al. 1995) between all pairs of individuals in a
population, with SPAGEDI v.1.4 (Hardy & Vekemans
2002); the genotype data for each population
were used as the reference for allelic frequencies.
Genetic differentiation between populations was
summarised by calculating pairwise F ST estimates
(Weir & Cockerham 1984) with FSTAT. The same
software was used to calculate five different indices
for sex-biased dispersal: F IS/ F ST' relatedness, and
mean and variance of assignment indices (I_a), which
quantified how likely an individual's multilocus
genotype is within its sampled population relative
to the average, with positive values indicating likely
residents and negative values indicating potential
immigrants (Goudet et al. 2002). These analyses were
done separately for each species and were limited to
SAN, BLO and TDR, due to small sample size at LAJ.

The dispersal of offspring over limited distances
from their parents results in an increase in genetic
differentiation with geographic distance, through
a process known as isolation by distance (IBD;
Rousset 1997, 2000). We characterised the dispersal
patterns of *Rhabdomys* at the scale of each locality by
conducting IBD analyses for both species separately,
in SAN and TDR for Rbech and in SAN and BLO for
Rdil. IBD was estimated for both sexes combined and
for each sex separately, using the regression method
based on the expected linear relationship between
genetic and geographic distances (Rousset 2000). The
analyses were run in GENEPOP, using the pairwise
genetic differentiation estimator ë calculated between
individuals (Watts et al. 2007) and the logarithms of
Euclidean geographic distances between individuals
for SAN and TDR populations, or the non-
transformed Euclidean geographic distance for BLO,
given the linear distribution of habitats available for
*Rhabdomys* individuals in this locality. Mantel tests
with 10,000 permutations were performed to assess
the correlation between matrices of genetic and
geographic distances.

Factorial correspondence analysis (FCA, Fig. S2) of genotypic data was performed with Genetix v4.01 (Belkhir et al. 1996-2004). The clustering method of Pritchard et al. (2000), implemented in STRUCTURE, was performed to estimate the most likely number of genetic clusters (K) in the data set, and to assign individuals to the clusters. The model probabilistically assigns individuals to source clusters (or jointly to two or more clusters in cases of admixture) based on their genotypes, while simultaneously estimating allelic frequencies. Ten independent runs were carried out for each value of K (1 to 5). For each run, 1,000,000 iterations were carried out after a burn-in period of 600,000 iterations. Analyses were performed without prior assumptions of the species or the locality. Population genetic structure was also investigated through Discriminant Principal Component Analysis (DAPC), which can accommodate the absence of HW equilibrium (Jombart et al. 2010). Analyses were conducted in R software using the *adegenet* and *devtools* packages. The most likely number of genetic clusters was determined using the Bayesian Information Criterion (BIC; Lebarbier & Mary-Huard 2006) using the empirical criterion of a delta-BIC less than six (Kass & Raftery 1995). We performed five independent runs for each K value (1 to 5).

When substructure was detected in a population, it
was further investigated using the spatially explicit
<figure><img src="image_9.png" /></figure>

method MAPI, implemented in the R package 'mapi'
v.1.1.4 (Piry et al. 2016). MAPI is a smoothing-like
procedure applied to pairwise genetic measures
computed between pairs of georeferenced samples.
Pairwise genetic measures are assigned to a network
of ellipses linking the samples. A regular grid of
hexagonal cells is then superimposed on the study
area, and ellipse values are averaged within the
cells they intersect, with a weighting procedure to
limit long-distance effects. The analysis produces a
geographical map representing the spatial variation
in the average level of the genetic similarity or
dissimilarity measures, computed between sample
pairs. The method has low sensitivity to IBD and is
not based on a predefined population genetic model
(Piry et al. 2016). MAPI was run using pairwise
kinship coefficients (Loiselle's Q) computed between
all pairs of individuals within the population to build
the spatial genetic network of ellipses. Following
Piry et al. (2016), ellipse eccentricity was set to
0.975 (i.e. default value), and the resolution of the
grid cell was determined automatically by setting
the beta parameter to 0.25, as recommended for
random sampling. A permutation procedure (10,000
permutations), combined with a false discovery
rate approach (Piry et al. 2016), was applied to
identify major areas of high genetic continuity and
discontinuity.

Identification of hybrids in natural conditions

We conducted analyses to identify potential
interspecific hybrids between Rdil and Rbech
at the scale of the only locality of the dataset
where both species occurred in sympatry (SAN).
First, STRUCTURE was used to infer assignment
probabilities of individuals from SAN to either
Rbech or Rdil. Bayesian assignment probabilities (q_i)
for each species were estimated by constraining the
number of populations to two (K = 2, assuming two
species) and running 500,000 iterations, after a burn-
in of 200,000 iterations, using population information
to assist clustering (USEPOPINFO model) under
both correlated (which assumes that populations
share similar allele frequencies because of their
common origin or gene flow between them) and
independent allele frequency models. The average of
ten independent runs was obtained. High individual
genotype assignment probabilities (q_i > 0.95) indicated
membership to a parental species, and lower scores
(q_i < 0.95) reflected some ancestral admixture.

We also used NEWHYBRIDS v.1.1 (Anderson &
Thompson 2002), which applies Markov chain Monte
Carlo sampling to estimate posterior probabilities

of each individual's membership to BLO and SAN
as a pure parent of Rdil or Rbech, first or second-
generation hybrid (F1 or F2), or an F1 backcrossed
with a pure parent. The analysis was performed
multiple times using different overdispersed starting
values and chain parameters, as recommended
by the developers. Finally, we used Jeffrey's prior
probabilities and default genotype proportions, and
the Markov chain was left to run for 100,000 sweeps
following a burn-in of 100,000 iterations.

Genetic characteristics of behaviourally
identified social groups

These analyses included adult mice that were both
genotyped and radiotracked across four populations
sampled in three localities: SAN (two species), TDR
and BLO (Table S1). An earlier study of SAN and BLO
populations showed that individuals within the same
social group exhibited higher home range overlap
than individuals from different groups (Dufour et al.
2019). We tested whether home range overlap was
related to kinship. Using data from 104 radiotracked
individuals (Dufour et al. 2015a), we assessed the
correlation between home range overlap (measured
as the volume of interaction between two home
ranges: VI) and pairwise kinship (Loiselle's Q) using
Spearman's rank correlation for each population. A
positive correlation would, indirectly, indicate that
related individuals are more likely to form social
groups than unrelated ones.

We then focused on radiotracked mice whose NS
were identified in an earlier study (Dufour et al.
2019). Altogether, we retained 61 NS used by adult
mice for which we also had microsatellite data. Only
one of these was occupied by a sympatric Rdil, so
this category was excluded from further analysis.
We analysed whether NS sharing involved related
individuals in three populations of the two species:
Rbech in Sandveld (18 NS used by genotyped mice)
and TDR (eight NS), and Rdil in Bloemhof (34 NS).
In each case, a minimum of two to five adult mice
shared a NS. We compared the mean kinship among
NS co-occupants and the average kinship across
the population to a reference value of 0 (absence
of relatedness) with a Wilcoxon test. A value
significantly higher than 0 among co-nesters would
indicate that NS sharing occurs preferentially among
related individuals.

Results

Rbech and Rdil genetic structure

Microsatellite genotyping was conducted on 359
<figure><img src="image_10.png" /></figure>

<table><caption>**Table 2.** Estimates of within-population genetic diversity in the two species using ten microsatellite markers. Rbech – *Rhabdomys bechuanae*, Rdil – *R. dilectus dilectus*, SAN – Sandveld Nature Reserve, BLO – Bloemhof Nature Reserve, LAJ – Lajuma Research Station, TDR – Tussen Die Riviere Game Reserve, n – number of genotyped individuals, *a*<sub>r</sub> – allelic richness, estimated for a standard sample size of 17 individuals, *H*<sub>o</sub> – observed heterozygosity, *H*<sub>e</sub> – Nei's unbiased genetic diversity, *F*<sub>is</sub> (10/7) – inbreeding coefficient calculated using ten *loci* or seven *loci* (when excluding those with null alleles).</caption><thead><tr><th>Species</th><th>Locality</th><th>n</th><th>a<sub>r</sub></th><th>H<sub>o</sub></th><th>H<sub>e</sub></th><th>F<sub>is</sub> (10/7)</th></tr></thead><tbody><tr><td>Rdil</td><td>BLO</td><td>47</td><td>9.9 ± 2.5</td><td>0.73 ± 0.21</td><td>0.86 ± 0.05</td><td>0.15/0.02</td></tr><tr><td>Rdil</td><td>LAJ</td><td>19</td><td>10.1 ± 2.8</td><td>0.71 ± 0.21</td><td>0.86 ± 0.07</td><td>0.17/0.06</td></tr><tr><td>Rdil</td><td>SAN</td><td>57</td><td>8.1 ± 2.1</td><td>0.65 ± 0.26</td><td>0.81 ± 0.07</td><td>0.20/0.07</td></tr><tr><td>Rbech</td><td>SAN</td><td>166</td><td>8.4 ± 2.1</td><td>0.55 ± 0.23</td><td>0.80 ± 0.07</td><td>0.30/0.19</td></tr><tr><td>Rbech</td><td>TDR</td><td>70</td><td>8.3 ± 1.8</td><td>0.63 ± 0.24</td><td>0.78 ± 0.09</td><td>0.19/0.05</td></tr></tbody></table>

individuals (see Table S1). Four *loci* (Rh24, Rh29, Rh31, Rh36) were at HW equilibrium. All other *loci* displayed significant heterozygote deficiencies in several populations. Overall, positive *F*<sub>IS</sub> values were obtained for all Rbech and Rdil populations, ranging from 0.15 for Rdil at BLO to 0.3 for Rbech at SAN (Table 2). Linkage disequilibria (LD) between pairs of *loci* were significant for 85 of the 225 tests performed.

HW deviations and significant LD may have resulted from oversampling of closely related individuals. To evaluate this, we calculated the kinship coefficient $\rho$ (Loiselle et al. 1995) between all pairs of individuals within each population (Hardy & Vekemans 2002). The median kinship ranged from -0.018 to -0.010 within populations (Fig. S1). Despite these low median values, very high $\rho$ (> 0.25) were obtained for some pairs of individuals in all populations, indicating the occurrence of close relatives. LD and HW analyses were then re-done on a reduced dataset (167 individuals), excluding 292 individuals involved in pairs with $\rho$ values > 0.25 in each population. All LD tests were non-significant, suggesting that the ten *loci* were genetically independent. Seven *loci* (Rh02, Rh21, Rh24, Rh29, Rh31, Rh36, Rh46) were at HW equilibrium, with the exception of the SAN population of Rbech, which deviated at two *loci* (data not shown). The three other *loci* (Rh09, Rh23, Rh28) exhibited significant heterozygote deficiencies in most populations, associated with some null genotypes (respectively 21, 44, and 11) that suggested the occurrence of null alleles. Estimated null allele frequencies for these three *loci* were similar among populations, with means ranging from 0.21 (Rh09) to 0.26 (Rh28). Given the limited number of markers (ten) considered, we decided to keep these *loci* in the analyses of population genetic structure, while confirming that comparable results were obtained in the analyses performed without these *loci*.

Genetic diversity levels within localities were similar in Rdil and Rbech, using either the full dataset (Table 2) or the reduced dataset (results not shown). Positive *F*<sub>IS</sub> values were obtained for all Rbech and Rdil populations using the full dataset, but only for the Rbech population in SAN when the three *loci* with null alleles were excluded (Table 2). For Rdil, genetic diversity estimates were lower (with or without the three *loci* with null alleles) and were higher in the sympatric locality with Rbech (SAN) compared to the allopatric sites (BLO and LAJ). For Rbech, *F*<sub>IS</sub> was higher in the locality of sympatry with Rdil (SAN) than in allopatry (TDR).

Mean *F*<sub>ST</sub> was 0.11 (95% Confidence Interval = (0.09, 0.13)) over the full dataset. Pairwise *F*<sub>ST</sub> varied from 0.03 (between BLO and LAJ) to 0.07 (between SAN and LAJ) among Rdil populations (mean *F*<sub>ST</sub> = 0.05, 95% CI = (0.03, 0.08)) and equalled 0.11 among Rbech populations (between SAN and TDR).

Different indices for sex-biased dispersal were calculated per species using Fstat. *F*<sub>IS</sub> was significantly higher in males (M) than in females (F) for Rbech (*F*<sub>IS</sub> (M) = 0.33 and *F*<sub>IS</sub> (F) = 0.23, *P* = 0.0002) and Rdil (*F*<sub>IS</sub> (M) = 0.22 and *F*<sub>IS</sub> (F) = 0.14, *P* = 0.0035), which suggested a male-biased dispersal. These patterns are consistent with a Wahlund effect in the dispersing sex due to a hidden substructure resulting from the coexistence of resident and immigrant individuals at each site (Goudet et al. 2002). Mean *I*<sub>a</sub> were significantly negative in males and positive in females for Rbech (*I*<sub>a</sub> (M) = -0.79 and *I*<sub>a</sub> (F) = 0.71, *P* = 0.009) and for Rdil (*I*<sub>a</sub> (M) = -0.91 and *I*<sub>a</sub> (F) = 0.80, *P* = 0.03), which is also expected under male bias in dispersal (Goudet et al. 2002). When the three *loci* with null alleles were excluded, the only significant result was for mean *I*<sub>a</sub> in Rbech, also suggesting male-biased dispersal. IBD analyses were significant (*P* < 0.0001) for both
<figure><img src="image_11.png" /></figure>

<figure><img src="image_12.png" /><figcaption>**Fig. 2.** Best clustering results of *Rhabdomys* individuals according to STRUCTURE and DAPC analyses. Each individual is displayed as a vertical line partitioned into coloured segments, the length of which is proportional to the individual's estimated K cluster membership coefficients. SAN – Sandveld Nature Reserve, TDR – Tussen Die Riviere Game Reserve, BLO – Bloemhof Nature Reserve, LAJ – Lajuma Research Station. A) STRUCTURE analysis, K = 3 and K = 4; B) DAPC analysis, K = 3 and K = 4. The arrow shows the position of individual no. 508 (see text).</figcaption></figure>

Rbech and Rdil within each population for both sexes combined and for each sex separately (SAN and BLO for Rdil, SAN and TDR for Rbech). Slope values, for both species, were similar and significantly positive (Table 3). Similar results were obtained when the three *loci* with null alleles were excluded from the analyses.

FCA, STRUCTURE and DAPC analyses performed on the full dataset provided consistent results. In FCA, the first axis separated the three populations of Rdil from the two populations of Rbech (Fig. S2). Only one

individual (no. 508, then identified as a putative hybrid in dedicated analyses) was misclassified, grouped with Rdil individuals, although it was identified as Rbech using COI. The second axis mainly separated the SAN and TDR populations of Rbech. A similar pattern was obtained on the dataset without *loci* having null alleles. STRUCTURE analysis indicated that the most likely genetic structure for the dataset consisted of three clusters (Fig. 2A). At K = 2, Rdil individuals from all localities grouped together, and Rbech individuals from TDR grouped alternatively with Rdil individuals (in four of the ten runs), or

<table><caption>**Table 3.** Isolation by distance analyses between individuals at a locality scale. Rbech – *R. bechuanae*, Rdil – *R. dilectus dilectus*. Analyses were performed with Genepop, using the pairwise genetic differentiation estimator *ê* calculated between individuals (Watts et al. 2007) and either the logarithms of Euclidean geographic distances between individuals for SAN (Sandveld Nature Reserve) and TDR (Tussen Die Riviere Game Reserve) populations, or the non-transformed Euclidean geographic distance for BLO (Bloemhof Nature Reserve).</caption><thead><tr><th rowspan="2"></th><th rowspan="2"></th><th colspan="2">Sex combined</th><th colspan="2">Males</th><th colspan="2">Females</th></tr><tr><th>Slope</th><th>95% CI</th><th>Slope</th><th>95% CI</th><th>Slope</th><th>95% CI</th></tr></thead><tbody><tr><td rowspan="2">Rbech</td><td>SAN</td><td>0.019</td><td>(0.01, 0.03)</td><td>0.014</td><td>(0.009, 0.027)</td><td>0.024</td><td>(0.015, 0.036)</td></tr><tr><td>TDR</td><td>0.011</td><td>(0.007, 0.015)</td><td>0.013</td><td>(0.007, 0.02)</td><td>0.013</td><td>(0.008, 0.03)</td></tr><tr><td rowspan="2">Rdil</td><td>SAN</td><td>0.019</td><td>(0.013, 0.029)</td><td>0.02</td><td>(0.01, 0.037)</td><td>0.02</td><td>(0.01, 0.03)</td></tr><tr><td>BLO</td><td>2.27e<sup>-6</sup></td><td>(1.96e<sup>-7</sup>, 5.16e<sup>-6</sup>)</td><td>-2.03e<sup>-6</sup></td><td>(-7.90e<sup>-6</sup>, 5.37e<sup>-6</sup>)</td><td>3.24e<sup>-6</sup></td><td>(-4.40e<sup>-7</sup>, 7.13e<sup>-6</sup>)</td></tr></tbody></table>
<figure><img src="image_13.png" /></figure>

<figure><img src="image_14.png" /><figcaption>**Fig. 3.** Spatial genetic structure of *Rhabdomys bechanae* in the locality of SAN, assessed using MAPI. The colour scale represents the spatial variation of the average level of the kinship coefficient of Loiselle, from lowest values in green to highest values in maroon; significant areas of high genetic similarity and discontinuity are delineated by dotted and full black lines, respectively. Sampled individuals are represented by a black dot. Background: Google Earth © 2025 CNES/Airbus.</figcaption></figure>

with Rbech individuals from SAN that formed a distinct group (in six of the ten runs). At K = 3, Rdil individuals from all localities grouped together, Rbech individuals from SAN were largely assigned to a second group (except two individuals: no. 282, which grouped with TDR Rbech individuals, and no. 508, which grouped with Rdil individuals), while Rbech individuals from TDR were assigned to a third group. At K = 4, Rbech individuals from SAN were separated into two different groups. Similar patterns were obtained when the three *loci* having null alleles were excluded from the analyses. DAPC analyses suggested up to four genetic groups, using the criteria of deltaBIC > 6. The pattern inferred was globally similar to that obtained with STRUCTURE (Fig. 2B), including individual no. 508 being classified as Rdil. All these analyses indicated stronger structuration between and within populations of Rbech compared to Rdil. Unfortunately, although individual no. 508 COI type was confirmed, its microsatellite genotype could not be reassessed to confirm that it was indeed a misassigned specimen.

FIS values, DAPC and STRUCTURE analyses suggested substructuring in Rbech at SAN. Rbech sampling was done in three separate areas in SAN, separated by 1 to 2 km, where both genetic clusters identified by STRUCTURE and DAPC were represented. We focused on the central area, where most individuals (133/165) were sampled, excluding individual no. 508, which appeared to be misclassified in the analyses performed on the full dataset. The MAPI analysis run on Loiselle's kinship coefficients clearly identified two areas of high genetic similarity, separated by a significant zone of genetic discontinuity along two road tracks (Fig. 3). This result was consistent with our observation that home range overlap was greater between mice with high genetic similarity than between mice with lower genetic similarity (Table 4).

## Identification of hybrids in SAN

We analysed the assignment probabilities of the 223 individuals, trapped in the sympatric locality (SAN), which were identified using mitochondrial
<figure><img src="image_15.png" /></figure>

<table><caption>**Table 4.** Admixed individuals identified using STRUCTURE and NEWHYBRIDS analyses in the locality of SAN (Sandveld Nature Reserve), where *R. bechuanae* (Rbech) and *R. dilectus dilectus* (Rdil) occurred in sympatry. STRUCTURE analyses were performed using the admixture model for individual ancestry, under both correlated and independent allele frequencies models, and assuming K = 2 (two species).</caption><thead><tr><th rowspan="2">No. Ind.</th><th rowspan="2">Mitochondrial genotype</th><th colspan="3">STRUCTURE analysis</th></tr><tr><th>Correlated allele frequencies</th><th>Independent allele frequencies</th><th>NEWHYBRID analysis – qi values</th></tr></thead><tbody><tr><td>282</td><td>Rbech</td><td>ADM</td><td>ADM</td><td>P<sub>Rbech</sub> = 0.53, F<sub>2</sub> = 0.40</td></tr><tr><td>284</td><td>Rbech</td><td>ADM</td><td>ADM</td><td>Bx<sub>Rbech</sub> = 0.50, F<sub>2</sub> = 0.31</td></tr><tr><td>305</td><td>Rbech</td><td>P<sub>Rbech</sub></td><td>P<sub>Rbech</sub></td><td>P<sub>Rbech</sub> = 0.85, Bx<sub>Rbech</sub> = 0.13</td></tr><tr><td>408</td><td>Rbech</td><td>ADM</td><td>ADM</td><td>F<sub>2</sub> = 0.92</td></tr><tr><td>474</td><td>Rbech</td><td>ADM</td><td>ADM</td><td>F<sub>2</sub> = 0.69, P<sub>Rbech</sub> = 0.29</td></tr><tr><td>486</td><td>Rbech</td><td>ADM</td><td>ADM</td><td>P<sub>Rbech</sub> = 0.81, F<sub>2</sub> = 0.16</td></tr><tr><td>508</td><td>Rbech</td><td>ADM</td><td>ADM</td><td>F<sub>2</sub> = 0.99</td></tr><tr><td>545</td><td>Rbech</td><td>ADM</td><td>ADM</td><td>F<sub>2</sub> = 0.61, Bx<sub>Rbech</sub> = 0.31</td></tr><tr><td>547</td><td>Rbech</td><td>ADM</td><td>ADM</td><td>P<sub>Rbech</sub> = 0.79, F<sub>2</sub> = 0.16</td></tr><tr><td>548</td><td>Rbech</td><td>ADM</td><td>ADM</td><td>P<sub>Rbech</sub> = 0.93</td></tr><tr><td>440</td><td>Rdil</td><td>ADM</td><td>ADM</td><td>F<sub>2</sub> = 0.65, Bx<sub>Rdil</sub> = 0.16</td></tr></tbody></table>

ADM, admixture nature (mean posterior probability of assignment in the given species across ten runs < 0.95); P<sub>Rbech</sub>, pure Rbech, F<sub>2</sub>, second generation hybrid, Bx<sub>Rbech</sub> or Bx<sub>Rdil</sub>, backcross with Rbech or Rdil, respectively.

COI sequences as Rbech or as Rdil, using all ten loci. Using STRUCTURE, we found 157 individuals as 'pure' Rbech, and 56 individuals as 'pure' Rdil under both the correlated allele frequency model and the unlinked allele frequency model (Table 5). Only one individual (no. 440) identified as Rdil with the COI marker was characterised by a mean $q_i < 0.95$, suggesting admixture, although not recent, as it was not assigned to a specific class of recent hybrids. Nine individuals identified as Rbech were characterised by mean $q_i < 0.95$ under both models, but only two of them (no. 284, no. 545) had a high probability of having a Rdil recent ancestor, in their case a grandparent (only no. 545 under the unlinked allele frequency model), and one (no. 508) carried a Rdil genotype ($q_i > 0.8$ under both models).

<table><caption>**Table 5.** Percentage of home range overlaps (measured as the volume of surface interaction, VI) between *R. bechuanae* individuals of the two genetic groups: West (W) and East (E) extracted from the MAPI analysis; see Fig. 3 and main text.</caption><thead><tr><th>Genetic groups of pairs of individuals</th><th>Number of individual pairs</th><th>Mean ± SE overlap Vi (%)</th></tr></thead><tbody><tr><td>W-W (same group)</td><td>36</td><td>23.94 ± 4.05</td></tr><tr><td>W-E</td><td>48</td><td>2.77 ± 0.73</td></tr><tr><td>E-E (same group)</td><td>20</td><td>41.45 ± 6.51</td></tr></tbody></table>

Using NEWHYBRIDS, we identified 156 individuals as 'pure' Rbech, and 56 individuals as 'pure' Rdil ($P > 0.95$). One individual identified with COI as Rdil (no. 440) had a low probability ($P < 0.95$) of being a 'pure' Rdil (Table 3); it was compatible with a F2 crossing type with a probability that was not very strong ($P = 0.65$). Ten individuals identified as Rbech had limited probabilities ($P < 0.95$) to be 'pure' Rbech. Among them, only four were compatible with an F2 crossing type (no. 408, $P = 0.92$; no. 474, $P = 0.67$; no. 508, $P = 0.99$; no. 545, $P = 0.61$) and only one (although with a low probability) with a backcross crossing type (no. 284, $P = 0.50$).

In summary, both methods gave similar results. One individual identified as Rdil (no. 440) and nine individuals identified as Rbech (no. 282, 284, 408, 474, 486, 508, 545, 547 and 548) with the mitochondrial marker, were identified as hybrids based on the nuclear markers by both STRUCTURE and NEWHYBRIDS analyses. One additional hybrid (no. 305) was identified by NEWHYBRIDS analyses only. Only two of these individuals (no. 284 and no. 545) were identified as recent hybrids by both methods, and one (no. 508) as either a recent hybrid (by NEWHYBRIDS) or an Rdil individual with an Rbech mitochondrial haplotype (by STRUCTURE). Spatially, these ten individuals, identified as hybrids using both methods, were relatively dispersed across the sympatric area.
<figure><img src="image_16.png" /></figure>

<figure><img src="image_17.png" /></figure>

## Home range overlap and kinship

A positive correlation between kinship and home range overlap was detected among Rbech adult individuals in SAN (147 pair comparisons: Spearman $r_s = 0.275$, $P = 0.0008$), indicating that related individuals were more likely to share part of their home range surface, even when they did not share an NS/social group, than unrelated individuals. A similar pattern was found for Rbech in TDR, although the trend was not statistically significant (45 pair comparisons: $r_s = 0.282$, $P = 0.06$). For Rdil individuals in BLO for which we had a similar sample size as for Rbech in SAN we did not detect a correlation between overlap and kinship (115 pair comparisons: $r_s = 0.07$, $P = 0.46$), and our sample size was too small to be able to conclude for Rdil when in sympatry in SAN (13 pair comparisons: $r_s = -0.152$, $P = 0.62$).

## Nest site sharing and kinship

Average kinship values, calculated for all pairs in a given locality, varied between $0.001 \pm 0.002$ (Rbech-TDR), $0.002 \pm 0.002$ (Rdil-BLO), and $0.0035 \pm 0.001$ (Rbech-SAN), all means not significantly different from 0 ($P > 0.05$). While within NS, they varied between $0.045 \pm 0.017$ (n = 57, Rdil-BLO), $0.049 \pm 0.02$ (n = 26, Rbech-TDR), and $0.07 \pm 0.02$ (n = 28, Rbech-SAN), all means significantly different from 0 ($P < 0.05$). Adult mice sharing the same NS tended thus to be more related to each other (values > 0), unlike the average calculated for the population.

## Discussion

Investigating the relationships between spatial and social behaviour and population genetics is particularly relevant to the study of the structure and functioning of cryptic species. This study presents the first population genetic analysis of two *Rhabdomys* species and reveals evidence of rare hybridisation occurring between them in sympatry. We also further characterised the social organisation of the two species in light of genetic relatedness, underscoring the role of social fences in limiting gene flow at a microscale. Although this study considered a limited number of populations, it provides new insights into understanding intra- and interspecies interactions.

## Population structure at different geographical scales

The study comprised two populations of Rbech and three populations of Rdil, distributed in four localities, including one site where the two species occurred

in sympatry. Genetic diversity was relatively high within all populations and similar for Rdil and Rbech. Mean genetic diversity estimates ($r, H_o, H_E$) were consistent with those reported in other *Rhabdomys* populations studied at a similar spatial scale (Solmsen et al. 2011, Froeschke & Sommer 2014), as well as in other African rodent species (e.g. Brouat et al. 2007, Gryseels et al. 2016). These relatively high values are expected in rodent populations that are generally assumed to be large. Consistent with findings in *R. pumilio* (Solmsen et al. 2011), we detected evidence of male-biased dispersal in both Rbech and Rdil. This pattern, widespread among rodents (Wolff 2007), likely results from male dispersal and mixing across different 'inbred' demes outside their natal areas.

High $F_{IS}$ values were observed in all populations, but were probably mainly due to null alleles at three *loci*. In the Rbech population of SAN, however, $F_{IS}$ remained high even after excluding these *loci*, suggesting significant substructuring. This likely reflects the presence of different breeding units (or family groups), consistent with the Wahlund effect (Waples 2015). According to this hypothesis, STRUCTURE DAPC analyses revealed the occurrence of two well-differentiated subpopulations of Rbech in SAN. The MAPI analysis suggested genetic discontinuities along two dirt roads that cross the sampled area. However, it is difficult to imagine that these narrow roads act as barriers to Rbech dispersal, as several individuals have been observed crossing them. Genetic discontinuities may instead point to behavioural isolation as a contributing factor. Solmsen et al. (2011) showed that territories in *R. pumilio* acted as social barriers to dispersal. In Rbech, strong social structure and group cohesion, together with nest fidelity, have been reported (Dufour et al. 2019), and may reinforce 'social fences' as suggested by our results. Interestingly, no population sub-structuring was detected in Rdil, and there was no relationship between kinship and home range overlap between Rdil individuals, consistent with earlier findings that Rdil individuals displayed looser group cohesion and NS fidelity compared to Rbech (Dufour et al. 2019), and further indicating that grouping in Rdil probably does not involve related individuals.

The Rdil population at SAN had the lowest genetic diversity observed for the species, suggesting a smaller population size at this site compared to the two other localities. A previous study noted density fluctuations with periods of low density of Rdil compared to Rbech in SAN (Ganem et al. 2020).
<figure><img src="image_18.png" /></figure>

<figure><img src="image_19.png" /></figure>

Because Rdil prefers mesic habitats with dense ground cover (Abu Baker & Brown 2010, Meynard et al. 2012, Dufour et al. 2015a), a smaller population size in the semi-arid SAN and BLO sites is plausible compared to LAJ, which has a subtropical climate and vegetation. Smaller population sizes are also expected at the margin of species distribution, for both species in SAN, where patchy and marginal habitats are expected to exacerbate genetic drift and reduce gene flow (Langin et al. 2017). However, although BLO is located both in a semi-arid area and at the margin of the Rdil distribution, this population did not show lower genetic diversity compared to the sympatric population, as well as the central population of LAJ inhabiting a subtropical mesic habitat. Small population sizes in sympatry have been observed in other rodents and were ascribed to a reduction in habitat quality due to the presence of competitor species (Eccard & Ylönen 2003). An earlier study suggested that species interference resulted in smaller home ranges, restricted individual mobility and NS accessibility as well as forced grouping of individuals in NS of lower quality in Rdil populations in sympatry (SAN) as compared to parapatry (Dufour et al. 2015a, 2019), consistent with the negative impact of Rbech on Rdil population in sympatry. In Rbech, high densities observed in SAN during our main sampling sessions (Dufour et al. 2015a) and similar genetic diversity levels in SAN (where Rbech occurs in sympatry with Rdil) and in TDR (where Rbech does not occur in sympatry with Rdil) suggest the absence of a negative effect of Rdil on Rbech.

Genetic differentiation among populations of Rdil was similar to that found in other *Rhabdomys* populations (Froeschke & Sommer 2014) and other rodent species (e.g. Brouat et al. 2007) at comparable spatial scales. Moreover, pairwise $F_{ST}$ values among the three populations of Rdil were consistent with an IBD effect; the most geographically distant LAJ population was also the most genetically differentiated from the others. Although this pattern should be verified with sampling of more populations, our results suggest that genetic differentiation is primarily driven by limited long-distance dispersal for Rdil (Leblois et al. 2000). In contrast, genetic differentiation between SAN and TDR Rbech populations was surprisingly high and on the same order of magnitude as the genetic differentiation between Rbech and Rdil. Given the similar levels of heterozygosity in Rbech and Rdil populations, and considering that geographical distances were much greater between Rdil populations, we did not expect higher $F_{ST}$ values in Rbech than among populations

of Rdil. The extremely high genetic differentiation observed between TDR and SAN may indicate very limited gene flow, compared to Rdil, potentially reflecting the unique geographical setting of TDR, which is bordered by the Rivers Caledon and Orange and surrounded by unfavourable rocky habitats that could form ecological barriers that may contribute to its isolation from other Rbech populations. However, if this were the case, we would have expected to observe lower genetic diversity in TDR compared to the more connected and less isolated SAN population, unless the TDR population size was particularly large despite its relative isolation. Nevertheless, introgression of Rbech with another unrecognised species of *Rhabdomys* in TDR or the presence in this locality of another cryptic species are hypotheses that warrant further, potentially genomic, investigations.

## Hybridisation

Hybridisation analyses were conducted for populations at the sympatric site (SAN), where trapping and radiotracking indicated that mice of the two species were syntopic, suggesting that interbreeding is possible (Dufour et al. 2015a). In an earlier study, laboratory behavioural trials with individuals from SAN demonstrated that, although potential mates of the two species were capable of mutual discrimination, some adults engaged in inter-specific mating behaviour when conspecific mates were unavailable (Dufour et al. 2015b). This finding suggests that reproductive interference may occur in sympatry (Groening & Hochkirch 2008). However, breeding experiments indicate lower reproductive success of heterospecific as compared to homospecific pairs (Dufour 2014).

STRUCTURE and NEWHYBRIDS analyses suggested that a hybrid ancestry could have occurred in nine out of 166 (5%) Rbech individuals and one out of 57 (2%) Rdil individuals. Most of these putative hybrids had a Rbech mitochondrial genome, suggesting a higher rate of backcrosses with Rbech females, which is consistent with a smaller population size of Rdil compared to Rbech (Dufour et al. 2015a, Ganem et al. 2020). No $F_1$ individuals were detected, and only a few individuals could be confidently assigned as $F_2$ or backcrosses (one identified by both STRUCTURE and NEWHYBRIDS, and two by one of the two methods only). However, hybrid assignments based on microsatellite data may include false positives (Poelstra et al. 2022), and the limited number of *loci* used here may further constrain detection accuracy. Alternatively, the presence of null alleles at three of the ten *loci* may reduce our ability to detect hybrids by
<figure><img src="image_20.png" /></figure>

<figure><img src="image_21.png" /></figure>

misclassifying some of them as homozygotes. A larger
and more reliable genomic dataset would therefore
be required before drawing firm conclusions about
the occurrence and extent of hybridisation between
the two species. For now, our results indicate that
despite syntopy, if hybridisation occurs, it remains
relatively rare – possibly limited to periods when one
of the species experiences a demographic trough.

Finally, the low frequency and likely reduced fitness
of hybrid individuals, combined with evidence
of behavioural discrimination and absence of co-
nesting, suggest the occurrence of both pre- and
postzygotic reproductive barriers. These findings
support the maintenance of reproductive isolation
despite syntopy, suggesting that these sister
species are either in the final stages of speciation or
represent fully differentiated species with infrequent
hybridisation. Our results, therefore, provide
empirical support for the role of, at least, behavioural
barriers in sustaining reproductive boundaries
under conditions of potential gene flow. Similar
patterns have been observed in the vlei rat Otomys
irroratus, where both behavioural and intrinsic (e.g.
karyotypic) barriers contribute to reproductive
isolation between allopatric populations (Pillay et
al. 1995a, b, c), suggesting the occurrence of multiple
isolating mechanisms in the maintenance of species
boundaries among rodents.

Conclusion

Integrating population genetics, spatial behaviour
and social structure provides strong insight into
intra- and interspecific dynamics in natural rodent
populations. Our genetic results are consistent with
behavioural observations of differences between
Rbech and Rdil in social and spatial structuring (i.e.
Rbech being more social than Rdil) and further inform
on the role of kinship in patial structuring in Rbech,
which may also contribute to isolation between the
two species. Overall, microsatellite markers used in

this study revealed little evidence of introgression
between species. However, comprehensive genomic
analyses across a wider range of populations will
be necessary to confirm this pattern. This study
reinforces the importance of integrating behavioural
ecology with population genetics when investigating
the evolutionary dynamics of closely related taxa.

Acknowledgements

We thank the many individuals who contributed to sample
collection, especially J. Watson, J. Du Plessis and J. Perez.
We thank the reserves managers who facilitated our work at
Sandveld Nature Reserve, Bloemhof Nature Reserve, and
Tussen Die Riviere Game Reserve. We also acknowledge the
collaborators who kindly provided samples from Lajuma;
no formal permit documentation was available for these
specimens. We thank the Departments of Environmental
Affairs of the Free State (DESTEA) and the North West
provinces for issuing us the research permits for this
study. Microsatellite genotyping was conducted through
the company Genoscreen, with complementary technical
support from the GenSeq platform at the Institut des
Sciences de l'Evolution de Montpellier (LabEx CeMEB,
ANR 'Investissements d'avenir' program, ANR-10-
LABX-04-01). This work would not have been possible
without specific financial support from ISEM (PEPS
funding) and the CNRS and NRF, which supported the
sampling through PICS (no. 4841, no. 81859). ISEM no.
2025-300.

Author Contributions

G. Ganem: fundraising, field sampling, co-writing,
behavioural data analysis, C.M.S. Dufour: field sampling,
manuscript editing, illustrations, N.L. Avenant: field
sampling, manuscript editing, K. Berthier: genetic
analyses, A. Loiseau: molecular work, P. Caminade: tissue
and DNA samples handling, molecular work, N. Pillay:
co-fundraising, manuscript editing, C. Brouat: genetic
analyses, co-writing.
<figure><img src="image_22.png" /></figure>

<figure><img src="image_23.png" /></figure>

## Literature

Abu Baker M.A. & Brown J.S. 2010: Islands of fear: effects of wooded patches on habitat suitability of the striped mouse in a South African grassland. *Funct. Ecol.* 24: 1313–1322.

Anderson E.C. & Thompson E.A. 2002: A model-based method for identifying species hybrids using multilocus genetic data. *Genetics* 160: 1217–1229.

Belkhir L.C., Raufaste N. & Bonhomme F. 1996–2004: GENETIX 4.05, logiciel sous Windows TM pour la génétique des populations. *Université de Montpellier, Laboratoire Génome, Populations, Interactions, CNRS UMR 5171, France*.

Benjamini Y. & Hochberg Y. 1995: Controlling the false discovery rate: a practical and powerful approach to multiple testing. *J. R. Stat. Soc. B* 57: 289–300.

Brouat C., Loiseau A., Kane M. et al. 2007: Population genetic structure of two ecologically distinct multimammate rats: the commensal *Mastomys natalensis* and the wild *Mastomys erythroleucus* in southeastern Senegal. *Mol. Ecol.* 16: 2985–2997.

Bruford M.W., Hanotte O., Brookfield J.F.Y. et al. 1992: Singlelocus and multilocus DNA fingerprint. In: Ha R. (ed.), Molecular genetic analysis of populations: a practical approach. *IRL Press, Oxford, UK:* 225–270.

Califf K.J., Green D.S., Wagner A.P. et al. 2020: Genetic relatedness and space use in two populations of striped hyenas (*Hyaena hyaena*). *J. Mammal.* 101: 361–372.

du Toit N., van Vuuren B.J., Matthee S. & Matthee C.A. 2012: Biome specificity of distinct genetic lineages within the four-striped mouse *Rhabdomys pumilio* (Rodentia: Muridae) from southern Africa with implications for taxonomy. *Mol. Phylogenet. Evol.* 65: 75–86.

Dufour C.M.S. 2014: Ecologie de la divergence et de la coéxistence: étude empirique chez deux espèces du genre *Rhabdomys*. *PhD thesis, Université Montpellier, France*.

Dufour C.M.S., Meynard C., Watson J. et al. 2015a: Space use variation in co-occurring sister species: response to environmental variation or competition? *PLOS ONE* 10: e0117750.

Dufour C.M.S., Pillay N., Avenant N. et al. 2019: Habitat characteristics and species interference influence space use and nest-site occupancy: implications for social variation in two sister species. *Oikos* 128: 503–516.

Dufour C.M.S., Pillay N. & Ganem G. 2015b: Ventro-ventral copulation in a rodent: a female initiative? *J. Mammal.* 96: 1017–1023.

Eccard J.A. & Ylönen H. 2003: Interspecific competition in small rodents: from populations to individuals. *Evol. Ecol.* 17: 423–440.

Froeschke G. & Sommer S. 2014: Role of selection versus neutral processes determining genetic variation in a small mammal along a climatic gradient in southern Africa. *Evol. Ecol.* 28: 1169–1190.

Ganem G., Dufour C.M.S., Avenant N.L. et al. 2020: An update on the distribution and diversification of *Rhabdomys* sp. (Muridae, Rodentia). *J. Vertebr. Biol.* 69: 20013.

Ganem G., Meynard C.N., Perigault M. et al. 2012: Environmental correlates and co-occurrence of three mitochondrial lineages of striped mice (*Rhabdomys*) in the Free State Province (South Africa). *Acta Oecol.* 42: 30–40.

Goudet J. 2001: FSTAT: a program to estimate and test gene diversities and fixation indices (version 2.9.3.2). https://www2.unil.ch/popgen/software/fstat.htm

Goudet J., Perrin N. & Waser P. 2002: Tests for sex-biased dispersal using bi-parentally inherited genetic markers. *Mol. Ecol.* 11: 1103–1114.

Groening J. & Hochkirch A. 2008: Reproductive interference between animal species. *Q. Rev. Biol.* 83: 257–282.

Gryseels S., Goüy de Bellocq J., Makundi R. et al. 2016: Genetic distinction between contiguous urban and rural multimammate mice in Tanzania despite gene flow. *J. Evol. Biol.* 29: 1952–1967.

Hall T.A. 1999: BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT. *Nucleic Acids Symp. Ser.* 41: 95–98.

Hardy O.J. & Vekemans X. 2002: SPAGEDI: a versatile computer program to analyse spatial genetic structure at the individual or population levels. *Mol. Ecol. Resour.* 2: 618–620.

Jombart T., Devillard S. & Balloux F. 2010: Discriminant analysis of principal components: a new method for the analysis of genetically structured populations. *BMC Genet.* 11: 94.

Kass R.E. & Raftery A.E. 1995: Bayes factors. *J. Am. Statist. Assoc.* 90: 773–795.

Langin K.M., Sillett T.S., Funk W.C. et al. 2017: Partial support for the central-marginal hypothesis within a population: reduced genetic diversity but not increased differentiation at the range edge of an island endemic bird. *Heredity* 119: 8–15.

Lebarbier É. & Mary-Huard T. 2006: Une introduction au critère BIC: fondements théoriques et interpretation. *J. Soc. Fr. Stat.* 147: 39–57.
<figure><img src="image_24.png" /></figure>

<figure><img src="image_25.png" /></figure>

Leblois R., Rousset F., Tikel D. et al. 2000: Absence of evidence for isolation by distance in an expanding cane toad (*Bufo marinus*) population: an individual-based analysis of microsatellite genotypes. *Mol. Ecol.* 9: 1905–1909.

Loiselle B.A., Sork V.L., Nason J. & Graham C. et al. 1995: Spatial genetic structure of a tropical understory shrub, *Psychotria officinalis* (Rubiaceae). *Am. J. Bot.* 82: 1420–1425.

Mahida H., Campbell G.K. & Taylor P.J. 1999: Genetic variation in *Rhabdomys pumilio* (Sparrman 1784) – an allozyme study. *Afr. Zool.* 34: 91–101.

Malausa T., Gilles A., Meglécz E. et al. 2011: High-throughput microsatellite isolation through 454 GS-FLX Titanium pyrosequencing of enriched DNA libraries. *Mol. Ecol. Resour.* 11: 638–644.

Meynard C.N., Pillay N., Perrigault M. et al. 2012: Evidence of environmental niche differentiation in the striped mouse (*Rhabdomys* sp.): inference from its current distribution in southern Africa. *Ecol. Evol.* 2: 1008–1023.

Morgan A.P., Hughes J.J., Didion J.P. et al. 2022: Population structure and inbreeding in wild house mice (*Mus musculus*) at different geographic scales. *Heredity* 129: 183.

Pillay N., Willan K. & Cooke J. 1995c: Evidence of pre-mating reproductive isolation in two populations of the vlei rat *Otomys irroratus*: experiments of intra- and interpopulation male-female encounters. *Z. Säugetierkd.* 60: 352–360.

Pillay N., Willan K. & Meester J. 1995a: Post-zygotic reproductive isolation in two populations of the vlei rat. *Acta Theriol.* 40: 69–76.

Pillay N., Willan K. & Meester J. 1995b: Evidence of pre-mating reproductive isolation in two allopatric populations of the vlei rat (*Otomys irroratus*). *Ethology* 100: 61–71.

Piry S., Chapuis M.P., Gauffre B. et al. 2016: Mapping Averaged Pairwise Information (MAPI): a new exploratory tool to uncover spatial structure. *Methods Ecol. Evol.* 7: 1463–1475.

Poelstra J.W., Montero B.K., Lüdemann J. et al. 2022: RADseq data reveal a lack of admixture in a mouse lemur contact zone contrary to previous microsatellite results. *Proc. Biol. Sci.* 289: 20220596.

Pritchard J.K., Stephens M. & Donnelly P. 2000: Inference of population structure using multilocus genotype data. *Genetics* 155: 945–959.

Rambau R.V., Robinson T.J. & Stanyon R. 2003: Molecular genetics of *Rhabdomys pumilio* subspecies boundaries: mtDNA phylogeography and karyotypic analysis by fluorescence in situ hybridization. *Mol. Phylogenet. Evol.* 28: 564–575.

Rousset F. 1997: Genetic differentiation and estimation of gene flow from F-Statistics under isolation by distance. *Genetics* 145: 1219–1228.

Rousset F. 2000: Genetic differentiation between individuals. *J. Evol. Biol.* 13: 58–62.

Rousset F. 2008: Genepop’007: a complete reimplementation of the Genepop software for Windows and Linux. *Mol. Ecol. Resour.* 8: 103–106.

Schradin C., Lindholm A.K., Johannesen J.E.S. et al. 2012: Social flexibility and social evolution in mammals: a case study of the African striped mouse (*Rhabdomys pumilio*). *Mol. Ecol.* 21: 541–553.

Solmsen N., Johannesen J. & Schradin C. 2011: Highly asymmetric fine-scale genetic structure between sexes of African striped mice and indication for condition dependent alternative male dispersal tactics. *Mol. Ecol.* 20: 1624–1634.

Storey J.D., Bass A.J., Dabney A. et al. 2021: qvalue: Q-value estimation for false discovery rate control. *R package version* 2.24.0. http://github.com/jdstorey/qvalue

Tung J., Charpentier M.J., Garfield D.A. et al. 2008: Genetic evidence reveals temporal change in hybridization patterns in a wild baboon population. *Mol. Ecol.* 17: 1998–2011.

Waples R.S. 2015: Testing for Hardy-Weinberg proportions: have we lost the plot? *J. Hered.* 106: 1–19.

Watts P.C., Rousset F., Saccheri I.J. et al. 2007: Compatible genetic and ecological estimates of dispersal rates in insect (*Coenagrion mercuriale*: Odonata: Zygoptera) populations: analysis of ‘neighbourhood size’ using a more precise estimator. *Mol. Ecol.* 16: 737–751.

Wauters L., Mazzamuto M., Santicchia F. et al. 2021: Personality traits, sex and food abundance shape space use in an arboreal mammal. *Oecologia* 196: 65–76.

Weir B.S. & Cockerham C.C. 1984. Estimating F-statistics for the analysis of population structure. *Evolution* 38: 1358–1370.

Wolff J.O. 2007: Social biology of rodents. *Integr. Zool.* 2: 193–2004.
<figure><img src="image_26.png" /></figure>

<figure><img src="image_27.png" /></figure>

**Supplementary online material**

**Fig. S1.** Distribution and median value of pairwise kinship coefficients of Loiselle et al. (1995) computed between all pairs of individuals from the same population. Rbech – *R. bechuanae*, Rdil – *R. dilectus*. See Fig. 1 for locality codes.

**Fig. S2.** FCA analyses conducted using Genetix v.4.01 on microsatellite data (ten loci) for *Rhabdomys bechuanae* (Rbech) and *Rhabdomys dilectus* (Rdil) individuals from four different localities. Square colours indicated the origin and species identity of individuals: yellow for Rbech and blue for Rdil individuals from SAN (Sandveld Nature Reserve); pink for Rbech individuals from TDR (Tussen Die Riviere Game Reserve); white for Rdil individuals from BLO (Bloemhof Nature Reserve); grey for Rdil individuals from LAJ (Lajuma Research Station). Only the first three axes of the FCA are shown.

(https://www.ivb.cz/wp-content/uploads/JVB-75-2026-GanemG.-et-al.-Fig.-S1-S2.pdf)

**Table S1.** Outline data analysed in this study, and other characteristics of individuals involved (https://www.ivb.cz/wp-content/uploads/JVB-75-2026-GanemG.-et-al.-Table-S1.xlsx).