Understanding how species responded to climate change in the past can help predict the long-term implications of contemporary climate change. The Upper Guinean rainforests of West Africa are a global biodiversity hotspot, and it is well documented that Pleistocene changes in climatic conditions drove the expansion and contraction of rainforest cover in this region. West African slippery frogs, genus Conraua, are rainforest specialists and present an excellent opportunity to study the effects of climate-driven landscape changes on population dynamics and contemporary biogeographic patterns. They exclusively inhabit rainforest streams and show little morphological or ecological disparity. We generated a population-level dataset of genome-wide RAD loci for four species spread across Upper Guinea. The observed phylogeographic structure is consistent with previously postulated macro- and microrefugia. Building on these results, we used demographic modelling to trace demographic trends over time and infer population connectivity patterns. Overall, populations in topographically complex regions, like the Fouta Djallon, showed signatures of long-term local persistence and milder changes in population size. In contrast, we found more dynamic histories of contraction and expansion in the main Upper Guinean rainforest block. Our findings provide insights into regional biodiversity patterns and show large variation in populations' responses to climatic fluctuations. This suggests that local environmental factors have played a key role in shaping population dynamics. Such insights are particularly relevant in relatively understudied biodiversity hotspots, such as the Upper Guinean rainforests of West Africa. Our results have implications for conservation management and prioritization at the species, site, and ecosystem level.

We used tissue samples stored in the collection of the MfN, Germany. We extracted DNA from 221 samples collected across the entire ranges of West African Conraua species (C. kamancamarai: 44 samples; C. alleni: 124; C. sagyimase: 12; C. derooi: 41). DNA was extracted from liver, muscle, or toe clips using the QIAGEN DNeasy Blood and Tissue Kits. Samples were collected from 53 sites, and the average number of samples per site was eight.

To generate genomic data, we applied the ddRADseq protocol described by Peterson et al. (2012). Our library preparation procedure closely resembles those implemented in previous African frog phylogeography studies (Charles et al., 2018; Portik et al., 2017). Briefly, the process involved: (i) digesting 500 ng DNA with restriction enzymes SbfI and MspI, (ii) purifying the digested fragments with magnetic beads and ligating barcoded Illumina adapters, (iii) selecting fragments between 415 and 515 bp after accounting for adapter length, and (iv) quantifying the prepared libraries. The sequencing was conducted on a single Illumina HiSeq4000 lane (50 bp single-end read sequencing with a final read length of 39 bp after trimming 11 bp restriction enzyme overhang and barcode).

After sequencing, we processed the raw Illumina reads using STACKS (V2.55; Catchen et al., 2011). First, we used USTACKS to generate stacks of loci for each individual based on the trimmed reads (39 bp). We built stacks if at least five reads (5X coverage) had a maximum of 2 bp mismatches (5% divergence). Next, we used CSTACKS to build a catalogue of loci, and then we employed SSTACKS to assign loci to each sample. Subsequently, we used POPULATIONS to call SNPs that are present in at least 70% of the individuals using the -r flag. We removed any individuals with missing data for more than 50% of loci. One challenge with RADseq data is allele dropout, which can result from mutations in restriction sites or degraded DNA (Graham et al., 2015; Gautier et al., 2013). To evaluate and mitigate these effects, we ran the STACKS pipeline once for the entire dataset, once for each species, and once for the following two sister species comparisons: (i) C. kamancamarai and C. alleni, and (ii) C. derooi and C. sagyimase. We performed population genetic analyses on the species-specific datasets to decrease the levels of missing data. However, in order to infer the species tree, to investigate interspecific gene flow, and to estimate pairwise F_st values, we employed the combined or full datasets.