Skip to main content

Thank you for visiting nature.com. You are using a browser version with limited support for CSS. To obtain the best experience, we recommend you use a more up to date browser (or turn off compatibility mode in Internet Explorer). In the meantime, to ensure continued support, we are displaying the site without styles and JavaScript.

  • Article
  • Published:

North African fox genomes show signatures of repeated introgression and adaptation to life in deserts

Nature Ecology & Evolution volume 7, pages 1267–1286 (2023)Cite this article

Subjects

Abstract

Elucidating the evolutionary process of animal adaptation to deserts is key to understanding adaptive responses to climate change. Here we generated 82 individual whole genomes of four fox species (genus Vulpes) inhabiting the Sahara Desert at different evolutionary times. We show that adaptation of new colonizing species to a hot arid environment has probably been facilitated by introgression and trans-species polymorphisms shared with older desert resident species, including a putatively adaptive 25 Mb genomic region. Scans for signatures of selection implicated genes affecting temperature perception, non-renal water loss and heat production in the recent adaptation of North African red foxes (Vulpes vulpes), after divergence from Eurasian populations approximately 78 thousand years ago. In the extreme desert specialists, Rueppell’s fox (V. rueppellii) and fennec (V. zerda), we identified repeated signatures of selection in genes affecting renal water homeostasis supported by gene expression and physiological differences. Our study provides insights into the mechanisms and genetic underpinnings of a natural experiment of repeated adaptation to extreme conditions.

This is a preview of subscription content, access via your institution

Access options

Buy this article

  • Purchase on Springer Link
  • Instant access to full article PDF
Buy now

Prices may be subject to local taxes which are calculated during checkout

Fig. 1: Species relationships and geographic distribution of North African foxes.
Fig. 2: Genetic structure, diversity and joint evolutionary history of North African red fox and Rueppell’s fox.
Fig. 3: Introgression in North African foxes.
Fig. 4: Demographic and selection history of North African red fox.
Fig. 5: Genome-wide selection scans for Rueppell’s fox and fennec using Ohana based on the structure and tree results from K = 4 highlight unique adaptations supported by expression and physiological differences.

Similar content being viewed by others

Data availability

The sequencing data (genomic and transcriptomic raw reads) generated in this study are available through NCBI repositories linked to BioProject accession number PRJNA951250. Physiological and morphological data per individual are available from GitHub (https://doi.org/10.5281/zenodo.7826260)132.

Code availability

All code written for this project is available on GitHub (https://doi.org/10.5281/zenodo.7826260)132 archived from repository https://github.com/joanocha/HOTFOXES.

References

  1. Williams, M. Climate Change in Deserts: Past, Present and Future (Cambridge Univ. Press, 2014).

  2. Riddell, E. A. et al. Exposure to climate change drives stability or collapse of desert mammal and bird communities. Science 371, 633–638 (2021).

    Article  CAS  PubMed  Google Scholar 

  3. Liu, Y. & Xue, Y. Expansion of the Sahara desert and shrinking of frozen land of the Arctic. Sci. Rep. 10, 4109 (2020).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  4. Sillero-Zubiri, C. et al. Canids: Foxes, Wolves, Jackals and Dogs (IUCN, 2004).

  5. Brito, J. C. et al. Biogeography and conservation of taxa from remote regions: an application of ecological-niche based models and GIS to North-African canids. Biol. Conserv. 142, 3020–3029 (2009).

    Article  Google Scholar 

  6. Williams, J. B. et al. Energy expenditure and water flux of Rüppell’s foxes in Saudi Arabia. Physiol. Biochem. Zool. 75, 479–488 (2002).

    Article  PubMed  Google Scholar 

  7. Williams, J. B. et al. A phylogenetic analysis of basal metabolism, total evaporative water loss, and life-history among foxes from desert and mesic regions. J. Comp. Physiol. B 174, 29–39 (2004).

    Article  CAS  PubMed  Google Scholar 

  8. Careau, V. et al. Basal metabolic rate of Canidae from hot deserts to cold arctic climates. J. Mammal. 88, 394–400 (2007).

    Article  Google Scholar 

  9. Noll-Banholzer, U. Water balance and kidney structure in the fennec. Comp. Biochem Physiol. A 62, 593–597 (1979).

    Article  Google Scholar 

  10. Maloiy, G. M. O. et al. Thermoregulation and metabolism in a small desert carnivore: the fennec fox (Fennecus zerda) (Mammalia). J. Zool. 198, 279–291 (1982).

    Article  Google Scholar 

  11. Rocha, J. L. et al. Convergent evolution of increased urine-concentrating ability in desert mammals. Mamm. Rev. 51, 482–491 (2021).

    Article  Google Scholar 

  12. Leite, J. V. et al. Differentiation of North African foxes and population genetic dynamics in the desert—insights into the evolutionary history of two sister taxa, Vulpes rueppellii and Vulpes vulpes. Org. Divers Evol. 15, 731–745 (2015).

    Article  Google Scholar 

  13. Sinding, M. H. S. et al. Arctic-adapted dogs emerged at the Pleistocene–Holocene transition. Science 368, 1495–1499 (2020).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  14. Wang, M. S. et al. Ancient hybridization with an unknown population facilitated high-altitude adaptation of canids. Mol. Biol. Evol. 37, 2616–2629 (2020).

    Article  CAS  PubMed  Google Scholar 

  15. von Holdt, B. et al. EPAS1 variants in high altitude Tibetan wolves were selectively introgressed into highland dogs. PeerJ 2017, e3522 (2017).

    Google Scholar 

  16. Li, Y. et al. Population variation revealed high-altitude adaptation of Tibetan mastiffs. Mol. Biol. Evol. 31, 1200–1205 (2014).

    Article  CAS  PubMed  Google Scholar 

  17. Wang, G.-D. et al. Genetic convergence in the adaptation of dogs and humans to the high-altitude environment of the Tibetan plateau. Genome Biol. Evol. 6, 2122–2128 (2014).

    Article  PubMed  PubMed Central  Google Scholar 

  18. Miao, B. et al. Genomic analysis reveals hypoxia adaptation in the Tibetan mastiff by introgression of the gray wolf from the Tibetan plateau. Mol. Biol. Evol. 34, 734–743 (2017).

    CAS  PubMed  Google Scholar 

  19. Signore, A. V. et al. Adaptive changes in hemoglobin function in high-altitude Tibetan canids were derived via gene conversion and introgression. Mol. Biol. Evol. 36, 2227–2237 (2019).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  20. Korneliussen, T. S. et al. ANGSD: analysis of next generation sequencing data. BMC Bioinformatics 15, 356 (2014).

    Article  PubMed  PubMed Central  Google Scholar 

  21. Hasselgren, M. et al. Genomic and fitness consequences of inbreeding in an endangered carnivore. Mol. Ecol. 30, 2790–2799 (2021).

    Article  PubMed  Google Scholar 

  22. Lindblad-Toh, K. et al. Genome sequence, comparative analysis and haplotype structure of the domestic dog. Nature 438, 803–819 (2005).

    Article  CAS  PubMed  Google Scholar 

  23. Statham, M. J. et al. Range-wide multilocus phylogeography of the red fox reveals ancient continental divergence, minimal genomic exchange and distinct demographic histories. Mol. Ecol. 23, 4813–4830 (2014).

    Article  PubMed  Google Scholar 

  24. McDevitt, A.D. et al. Next-generation phylogeography resolves post-glacial colonization patterns in a widespread carnivore, the red fox (Vulpes vulpes), in Europe. Mol. Ecol. 31, 993–1006 (2022).

    Article  PubMed  Google Scholar 

  25. Noskova, E. et al. GADMA: genetic algorithm for inferring demographic history of multiple populations from allele frequency spectrum data. Gigascience 9, giaa005 (2020).

    Article  PubMed  PubMed Central  Google Scholar 

  26. Gutenkunst, R. N. et al. Inferring the joint demographic history of multiple populations from multidimensional SNP frequency data. PLoS Genet. 5, 1000695 (2009).

    Article  Google Scholar 

  27. Bartolini Lucenti, S. & Madurell-Malapeira, J. Unraveling the fossil record of foxes: an updated review on the Plio-Pleistocene Vulpes spp. from Europe. Quat. Sci. Rev. 236, 106296 (2020).

    Article  Google Scholar 

  28. Hoffmann, D. L. et al. Timing and causes of North African wet phases during the last glacial period and implications for modern human migration. Sci. Rep. 6, 36367 (2016).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Drake, N. & Breeze, P. in Africa from MIS 6-2: Population Dynamics and Paleoenvironments (eds Jones, S. C. & Stewart, B. A.) 103–122 (Springer, 2016).

  30. Pausata, F. S. R. et al. The greening of the Sahara: past changes and future implications. One Earth 2, 235–250 (2020).

    Article  Google Scholar 

  31. Osborne, A. H. et al. A humid corridor across the Sahara for the migration of early modern humans out of Africa 120,000 years ago. Proc. Natl Acad. Sci. USA 105, 16444–16447 (2008).

    Article  PubMed  PubMed Central  Google Scholar 

  32. Castañeda, I. S. et al. Wet phases in the Sahara/Sahel region and human migration patterns in North Africa. Proc. Natl Acad. Sci. USA 106, 20159–20163 (2009).

    Article  PubMed  PubMed Central  Google Scholar 

  33. Green, R. E. et al. A draft sequence of the neandertal genome. Science 328, 710–722 (2010).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Pickrell, J. K. & Pritchard, J. K. Inference of population splits and mixtures from genome-wide allele frequency data. PLoS Genet. 8, 1002967 (2012).

    Article  Google Scholar 

  35. Martin, S. H. et al. Evaluating the use of ABBA–BABA statistics to locate introgressed loci. Mol. Biol. Evol. 32, 244–257 (2015).

    Article  CAS  PubMed  Google Scholar 

  36. Ainsworth, H. C. et al. Genetic epidemiology in kidney disease. Nephrol. Dial. Transplant. 32, ii159–ii169 (2017).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  37. Pramod, A. B. et al. SLC6 transporters: structure, function, regulation, disease association and therapeutics. Mol. Asp. Med. 34, 197–219 (2013).

    Article  CAS  Google Scholar 

  38. Hedrick, P. W. Adaptive introgression in animals: examples and comparison to new mutation and standing variation as sources of adaptive variation. Mol. Ecol. 22, 4606–4618 (2013).

    Article  PubMed  Google Scholar 

  39. Charlesworth, D. Balancing selection and its effects on sequences in nearby genome regions. PLoS Genet. 2, e64 (2006).

    Article  PubMed  PubMed Central  Google Scholar 

  40. Makunin, A. I. et al. Genes on B chromosomes of vertebrates. Mol. Cytogenet. 7, 99 (2014).

    Article  PubMed  PubMed Central  Google Scholar 

  41. Graphodatsky, A. S. et al. Phylogenomics of the dog and fox family (Canidae, Carnivora) revealed by chromosome painting. Chromosome Res. 16, 129–143 (2008).

    Article  CAS  PubMed  Google Scholar 

  42. Ilardo, M. A. et al. Physiological and genetic adaptations to diving in sea nomads. Cell 173, 569–580 (2018).

    Article  CAS  PubMed  Google Scholar 

  43. Campbell-Staton, S. C. et al. Ivory poaching and the rapid evolution of tusklessness in African elephants. Science 374, 483–487 (2021).

    Article  CAS  PubMed  Google Scholar 

  44. Stödberg, T. et al. SLC12A2 mutations cause NKCC1 deficiency with encephalopathy and impaired secretory epithelia. Neurol. Genet 6, e478 (2020).

    Article  PubMed  PubMed Central  Google Scholar 

  45. Evans, R. L. et al. Severe impairment of salivation in Na+/K+/2Cl− cotransporter (NKCC1)-deficient mice. J. Biol. Chem. 275, 26720–26726 (2000).

    Article  CAS  PubMed  Google Scholar 

  46. Sinica, V. et al. Human and mouse TRPA1 are heat and cold sensors differentially tuned by voltage. Cells 9, 57 (2019).

    Article  PubMed  PubMed Central  Google Scholar 

  47. Xu, M. et al. TRPV1 and TRPA1 in lung inflammation and airway hyperresponsiveness induced by fine particulate matter (PM2.5). Oxid. Med. Cell Longev. 2019, 7450151 (2019).

    Article  PubMed  PubMed Central  Google Scholar 

  48. Sinica, V. & Vlachová, V. Transient receptor potential ankyrin 1 channel: an evolutionarily tuned thermosensor. Physiol. Res 70, 363 (2021).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Lynch, V. J. et al. Elephantid genomes reveal the molecular bases of woolly mammoth adaptations to the Arctic. Cell Rep. 12, 217–228 (2015).

    Article  CAS  PubMed  Google Scholar 

  50. Key, F. M. et al. Human local adaptation of the TRPM8 cold receptor along a latitudinal cline. PLoS Genet. 14, e1007298 (2018).

    Article  PubMed  PubMed Central  Google Scholar 

  51. Ward, D. The Biology of Deserts 2nd edn (Oxford Univ. Press, 2016).

  52. Moraru, A. et al. THADA regulates the organismal balance between energy storage and heat production. Dev. Cell 41, 72–81 (2017).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Dempersmier, J. et al. Cold-inducible Zfp516 activates UCP1 transcription to promote browning of white fat and development of brown fat. Mol. Cell 57, 235–246 (2015).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Cheng, J. Y. et al. Detecting selection in multiple populations by modeling ancestral admixture components. Mol. Biol. Evol. 39, msab294 (2022).

    Article  CAS  PubMed  Google Scholar 

  55. Kurimoto, J. et al. Deficiency of WFS1 leads to the impairment of AVP secretion under dehydration in male mice. Pituitary 1, 3 (2021).

    Google Scholar 

  56. Bridges, T. E. & James, N. V. The hypothalamo-neurohypophysial system of native Australian desert rodents. The vasopressin and oxytocin contents of Notomys alexis and Pseudomys australis compared with those of the laboratory rat and mouse in different states of water balance. Aust. J. Exp. Biol. Med. Sci. 60, 265–283 (1982).

    Article  CAS  PubMed  Google Scholar 

  57. Donald, J. & Pannabecker, T. L. in Sodium and Water Homeostasis (eds Hyndman, K. & Pannabecker, T.) 191–211 (Springer, 2015).

  58. Schwimmer, H. & Haim, A. Physiological adaptations of small mammals to desert ecosystems. Integr. Zool. 4, 357–366 (2009).

    Article  PubMed  Google Scholar 

  59. Beuchat, C. A. Body size, medullary thickness, and urine concentrating ability in mammals. Am. J. Physiol. Regul. Integr. Comp. Physiol. 258, R298–R308 (1990).

    Article  CAS  Google Scholar 

  60. Solanki, A. K. et al. Mutations in KIRREL1, a slit diaphragm component, cause steroid-resistant nephrotic syndrome. Kidney Int. 96, 883–889 (2019).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  61. Noll-Banholzer, U. Body temperature, oxygen consumption, evaporative water loss and heart rate in the fennec. Comp. Biochem Physiol. A 62, 585–592 (1979).

    Article  Google Scholar 

  62. Ishikawa, A. & Kitano, J. in Thyroid Hormone (ed. Agrawal, N. K.) Ch. 2 (IntechOpen, 2012).

  63. Geffen, E. & Girard, I. in The Swift Fox: Ecology and Conservation of Swift Foxes in a Changing World (eds Sovada, M. A. & Carbyn, L.) 223–229 (Univ. Regina Press, 2003).

  64. Ferreira, M. S. et al. The legacy of recurrent introgression during the radiation of hares. Syst. Biol. 70, 593–607 (2021).

    Article  PubMed  Google Scholar 

  65. Jones, M. R. et al. Adaptive introgression underlies polymorphic seasonal camouflage in snowshoe hares. Science 360, 1355–1358 (2018).

    Article  CAS  PubMed  Google Scholar 

  66. Edelman, N. B. et al. Genomic architecture and introgression shape a butterfly radiation. Science 366, 594–599 (2019).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  67. Oziolor, E. M. et al. Adaptive introgression enables evolutionary rescue from extreme environmental pollution. Science 364, 455–457 (2019).

    Article  CAS  PubMed  Google Scholar 

  68. Witt, K. E. & Huerta-Sánchez, E. Convergent evolution in human and domesticate adaptation to high-altitude environments. Phil. Trans. R. Soc. B 374, 20180235 (2019).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  69. Cao, Y. H. et al. Historical introgression from wild relatives enhanced climatic adaptation and resistance to pneumonia in sheep. Mol. Biol. Evol. 38, 838–855 (2021).

    Article  CAS  PubMed  Google Scholar 

  70. Marques, D. A. et al. A combinatorial view on speciation and adaptive radiation. Trends Ecol. Evol. 34, 531–544 (2019).

    Article  PubMed  Google Scholar 

  71. vonHoldt, B. M. et al. Admixture mapping identifies introgressed genomic regions in North American canids. Mol. Ecol. 25, 2443–2453 (2016).

    Article  PubMed  Google Scholar 

  72. Koepfli, K. P. et al. Genome-wide evidence reveals that African and Eurasian golden jackals are distinct species. Curr. Biol. 25, 2158–2165 (2015).

    Article  CAS  PubMed  Google Scholar 

  73. Gopalakrishnan, S. et al. Interspecific gene flow shaped the evolution of the genus Canis. Curr. Biol. 28, 3441–3449 (2018).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Racimo, F. et al. Archaic adaptive introgression in TBX15/WARS2. Mol. Biol. Evol. 34, msw283 (2016).

    Article  Google Scholar 

  75. Huerta-Sánchez, E. et al. Altitude adaptation in Tibetans caused by introgression of Denisovan-like DNA. Nature 512, 194–197 (2014).

    Article  PubMed  PubMed Central  Google Scholar 

  76. Geffen, E. in The Biology and Conservation of Wild Canids (eds Macdonald, D. W. & Sillero-Zubiri, C.) 199–206 (Oxford Univ. Press, 2004).

  77. Kimura, R. et al. Gene flow and natural selection in oceanic human populations inferred from genome-wide SNP typing. Mol. Biol. Evol. 25, 1750–1761 (2008).

    Article  CAS  PubMed  Google Scholar 

  78. Jamie, G. A. & Meier, J. I. The persistence of polymorphisms across species radiations. Trends Ecol. Evol. 35, 795–808 (2020).

    Article  PubMed  Google Scholar 

  79. Thomas, L. et al. Spatially varying selection between habitats drives physiological shifts and local adaptation in a broadcast spawning coral on a remote atoll in Western Australia. Sci. Adv. 8, 9185 (2022).

    Article  Google Scholar 

  80. Cardona, A. et al. Genome-wide analysis of cold adaptation in indigenous Siberian populations. PLoS ONE 9, e98076 (2014).

    Article  PubMed  PubMed Central  Google Scholar 

  81. Vale, C. G. & Brito, J. C. Desert-adapted species are vulnerable to climate change: insights from the warmest region on Earth. Glob. Ecol. Conserv. 4, 369–379 (2015).

    Article  Google Scholar 

  82. Leclère, D. et al. Bending the curve of terrestrial biodiversity needs an integrated strategy. Nature 585, 551–556 (2020).

    Article  PubMed  Google Scholar 

  83. Rocha, J. L. et al. Life in deserts: the genetic basis of mammalian desert adaptation. Trends Ecol. Evol. 36, 637–650 (2021).

    Article  PubMed  Google Scholar 

  84. Global Aridity and PET Database (CGIAR-CSI, 2009); https://cgiarcsi.community/data/global-aridity-and-pet-database/

  85. Knepper, M. A. et al. Molecular physiology of water balance. N. Engl. J. Med. 372, 1349–1358 (2015).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  86. Sands, J. M. & Layton, H. E. Advances in understanding the urine-concentrating mechanism. Annu. Rev. Physiol. 76, 387–409 (2014).

    Article  CAS  PubMed  Google Scholar 

  87. Xu, M. M. & Wang, D. H. Water deprivation up-regulates urine osmolality and renal aquaporin 2 in Mongolian gerbils (Meriones unguiculatus). Comp. Biochem Physiol. A 194, 37–44 (2016).

    Article  CAS  Google Scholar 

  88. Morgenthaler, N. G. Copeptin: a biomarker of cardiovascular and renal function. Congest. Heart Fail. 16, S37–S44 (2010).

    Article  CAS  PubMed  Google Scholar 

  89. Morgenthaler, N. G. et al. Copeptin: clinical use of a new biomarker. Trends Endocrinol. Metab. 19, 43–49 (2008).

    Article  CAS  PubMed  Google Scholar 

  90. Hou, Z. et al. A cost-effective RNA sequencing protocol for large-scale gene expression studies. Sci. Rep. 5, 9570 (2015).

    Article  PubMed  PubMed Central  Google Scholar 

  91. Martin, M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet J. 17, 10 (2011).

    Article  Google Scholar 

  92. Bolger, A. M. et al. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics 30, 2114–2120 (2014).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  93. Li, H. Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM. Preprint at https://arxiv.org/abs/1303.3997 (2013).

  94. Li, H. et al. The Sequence Alignment/Map format and SAMtools. Bioinformatics 25, 2078–2079 (2009).

    Article  PubMed  PubMed Central  Google Scholar 

  95. Jónsson, H. et al. mapDamage2.0: fast approximate Bayesian estimates of ancient DNA damage parameters. Bioinformatics 29, 1682–1684 (2013).

    Article  PubMed  PubMed Central  Google Scholar 

  96. Kurtz, S. et al. Versatile and open software for comparing large genomes. Genome Biol. 5, R12 (2004).

    Article  PubMed  PubMed Central  Google Scholar 

  97. Katoh, K. & Standley, D. M. MAFFT Multiple Sequence Alignment Software Version 7: improvements in performance and usability. Mol. Biol. Evol. 30, 772–780 (2013).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  98. Gouy, M. et al. SeaView Version 4: a multiplatform graphical user interface for sequence alignment and phylogenetic tree building. Mol. Biol. Evol. 27, 221–224 (2010).

    Article  CAS  PubMed  Google Scholar 

  99. Korneliussen, T. S. & Moltke, I. NgsRelate: a software tool for estimating pairwise relatedness from next-generation sequencing data. Bioinformatics 31, 4009–4011 (2015).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  100. Dobin, A. et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics 29, 15–21 (2013).

    Article  CAS  PubMed  Google Scholar 

  101. Shumate, A. & Salzberg, S. L. Liftoff: accurate mapping of gene annotations. Bioinformatics 37, 1639–1643 (2021).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  102. Okonechnikov, K. et al. Qualimap 2: advanced multi-sample quality control for high-throughput sequencing data. Bioinformatics 32, 292–294 (2016).

    Article  CAS  PubMed  Google Scholar 

  103. Ewels, P. et al. MultiQC: summarize analysis results for multiple tools and samples in a single report. Bioinformatics 32, 3047–3048 (2016).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  104. Li, H. A statistical framework for SNP calling, mutation discovery, association mapping and population genetical parameter estimation from sequencing data. Bioinformatics 27, 2987–2993 (2011).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  105. Danecek, P. et al. The variant call format and VCFtools. Bioinformatics 27, 2156–2158 (2011).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  106. Chang, C. C. et al. Second-generation PLINK: rising to the challenge of larger and richer datasets. Gigascience 4, 7 (2015).

    Article  PubMed  PubMed Central  Google Scholar 

  107. Love, M. I. et al. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 15, 550 (2014).

    Article  PubMed  PubMed Central  Google Scholar 

  108. Kumar, S. et al. MEGA X: molecular evolutionary genetics analysis across computing platforms. Mol. Biol. Evol. 35, 1547–1549 (2018).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  109. Matzke, N. J. Model selection in historical biogeography reveals that founder-event speciation is a crucial process in island clades. Syst. Biol. 63, 951–970 (2014).

    Article  PubMed  Google Scholar 

  110. Matzke, N. J. Probabilistic Historical Biogeography: New Models for Founder-Event Speciation, Imperfect Detection, and Fossils Allow Improved Accuracy and Model-Testing. PhD thesis, Univ. California (2013).

  111. Sacks, B. N. et al. Nuclear genetic analysis of the red fox across its Trans-Pacific range. J. Hered. 109, 573–584 (2018).

    Article  PubMed  Google Scholar 

  112. Drummond, A. J. & Rambaut, A. BEAST: Bayesian evolutionary analysis by sampling trees. BMC Evol. Biol. 7, 214 (2007).

    Article  PubMed  PubMed Central  Google Scholar 

  113. Posada, D. jModelTest: phylogenetic model averaging. Mol. Biol. Evol. 25, 1253–1256 (2008).

    Article  CAS  PubMed  Google Scholar 

  114. Perini, F. A. et al. The evolution of South American endemic canids: a history of rapid diversification and morphological parallelism. J. Evol. Biol. 23, 311–322 (2010).

    Article  CAS  PubMed  Google Scholar 

  115. Kumar, V. et al. Genetic signatures of adaptation revealed from transcriptome sequencing of Arctic and red foxes. BMC Genom. 16, 585 (2015).

    Article  Google Scholar 

  116. Zhao, C. et al. The complete mitochondrial genome of the Tibetan fox (Vulpes ferrilata) and implications for the phylogeny of Canidae. C. R. Biol. 339, 68–77 (2016).

    Article  PubMed  Google Scholar 

  117. Leigh, J. W. & Bryant, D. popart: full-feature software for haplotype network construction. Methods Ecol. Evol. 6, 1110–1116 (2015).

    Article  Google Scholar 

  118. Vieira, F. G. et al. Improving the estimation of genetic distances from next-generation sequencing data. Biol. J. Linn. Soc. 117, 139–149 (2016).

    Article  Google Scholar 

  119. Lefort, V. et al. FastME 2.0: a comprehensive, accurate, and fast distance-based phylogeny inference program. Mol. Biol. Evol. 32, 2798–2800 (2015).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  120. Stamatakis, A. RAxML version 8: a tool for phylogenetic analysis and post-analysis of large phylogenies. Bioinformatics 30, 1312–1313 (2014).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  121. Soraggi, S. et al. Powerful inference with the D-statistic on low-coverage whole-genome data. G3 Genes Genomes Genet. 8, 551–566 (2018).

    Article  Google Scholar 

  122. Koch, E. M. et al. De novo mutation rate estimation in wolves of known pedigree. Mol. Biol. Evol. 36, 2536–2547 (2019).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  123. Kutschera, V. E. et al. A range-wide synthesis and timeline for phylogeographic events in the red fox (Vulpes vulpes). BMC Evol. Biol. 13, 114 (2013).

    Article  PubMed  PubMed Central  Google Scholar 

  124. Lonsinger, R. C. et al. Evaluating effective population size and genetic diversity of a declining kit fox population using contemporary and historical specimens. Ecol. Evol. 8, 12011–12021 (2018).

    Article  PubMed  PubMed Central  Google Scholar 

  125. Yi, X. et al. Sequencing of fifty human exomes reveals adaptations to high altitude. Science 329, 75–78 (2010).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  126. Kent, W. J. et al. The Human Genome Browser at UCSC. Genome Res. 12, 996–1006 (2002).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  127. Malaspinas, A. S. et al. A genomic history of Aboriginal Australia. Nature 538, 207–214 (2016).

    Article  CAS  PubMed  Google Scholar 

  128. Crawford, J. E. & Nielsen, R. Detecting adaptive trait loci in nonmodel systems: divergence or admixture mapping? Mol. Ecol. 22, 6131–6148 (2013).

    Article  PubMed  Google Scholar 

  129. Hudson, R. R. et al. A test of neutral molecular evolution based on nucleotide data. Genetics 116, 153–159 (1987).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  130. Zhai, W. et al. An investigation of the statistical power of neutrality tests based on comparative and population genetic data. Mol. Biol. Evol. 26, 273–283 (2009).

    Article  CAS  PubMed  Google Scholar 

  131. Liu, S. et al. Population genomics reveal recent speciation and rapid evolutionary adaptation in polar bears. Cell 157, 785–794 (2014).

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  132. Rocha, J.L. et al. Data and code from: North-African fox genomes show signatures of repeated introgression and adaptation to life in deserts. GitHub https://doi.org/10.5281/zenodo.7826260 (2023).

Download references

Acknowledgements

We thank M. Samlali, B. Habib and M. Milanese for their assistance in sampling collection in North Africa. We thank J. V. López Bao, H. R. Maior and F. Alvares for support and assistance in sampling collection in Iberia and J. F. Layna (Spain), V. Soeiro (Parque Biológico de Gaia, Portugal), Temara Zoo and Dream Village (Morocco) for access to foxes in captivity. We thank L. Smith (UC Berkeley), P. Ribeiro, S. Mourão and S. Lopes (CIBIO, UP) for support during laboratory work and generation of data. We thank A. Múrias, C. Pacheco, D. Lobo, M. S. Ferreira, D. Brandt, D. Aguillar-Gomez, L. Pipes, H. Wang and T. Linderoth for helpful discussions during data analysis and interpretation of results. We thank L. Dalén and A. Gaur for kindly sharing the arctic fox reference genome and Bengal fox mitogenome, respectively. We acknowledge M. Laranjeira Rocha for the illustrations featured in this work. This project was supported by the Portuguese Foundation for Science and Technology, FCT (PTDC/BIA-EVF/31902/2017 granted to R.G.) and the National Institutes of Health (R01GM138634 granted to R.N.). J.L.R., N.S., M.N., J.C.B. and R.G. worked under scholarships or research contracts from FCT (SFRH/BD/116397/2016, SFRH/BPD/116596/2016, SFRH/BD/144087/2019, CEECINST/00014/2018/CP1512/CT0001 and 2021/00647/CEECIND, respectively).

Author information

Author notes

  1. These authors jointly supervised this work: Rasmus Nielsen, Raquel Godinho.

  2. These authors contributed equally: Rasmus Nielsen, Raquel Godinho.

Authors and Affiliations

  1. CIBIO, Centro de Investigação em Biodiversidade e Recursos Genéticos, InBIO Laboratório Associado, Campus de Vairão, Universidade do Porto, Vairão, Portugal

    Joana L. Rocha, Pedro Silva, Nuno Santos, Mónia Nakamura, Sandra Afonso, Zbyszek Boratynski, José C. Brito & Raquel Godinho

  2. Departamento de Biologia, Faculdade de Ciências, Universidade do Porto, Porto, Portugal

    Joana L. Rocha, Pedro Silva, Mónia Nakamura, José C. Brito & Raquel Godinho

  3. BIOPOLIS Program in Genomics, Biodiversity and Land Planning, CIBIO, Campus de Vairão, Vairão, Portugal

    Joana L. Rocha, Pedro Silva, Nuno Santos, Mónia Nakamura, Sandra Afonso, Zbyszek Boratynski, José C. Brito & Raquel Godinho

  4. Department of Integrative Biology and Department of Statistics, University of California Berkeley, Berkeley, CA, USA

    Joana L. Rocha, Peter H. Sudmant & Rasmus Nielsen

  5. Laboratory of Geophysics and Natural Hazards, Geophysics, Natural Patrimony and Green Chemistry Research Center (GEOPAC), Institut Scientifique, Mohammed V University of Rabat, Rabat, Morocco

    Abdeljebbar Qninba

  6. Center for Computational Biology, University of California, Berkeley, CA, USA

    Peter H. Sudmant & Rasmus Nielsen

  7. Globe Institute, University of Copenhagen, Copenhagen, Denmark

    Rasmus Nielsen

  8. Department of Zoology, University of Johannesburg, Auckland Park, South Africa

    Raquel Godinho

Authors
  1. Joana L. Rocha
    View author publications

    You can also search for this author in PubMed Google Scholar

  2. Pedro Silva
    View author publications

    You can also search for this author in PubMed Google Scholar

  3. Nuno Santos
    View author publications

    You can also search for this author in PubMed Google Scholar

  4. Mónia Nakamura
    View author publications

    You can also search for this author in PubMed Google Scholar

  5. Sandra Afonso
    View author publications

    You can also search for this author in PubMed Google Scholar

  6. Abdeljebbar Qninba
    View author publications

    You can also search for this author in PubMed Google Scholar

  7. Zbyszek Boratynski
    View author publications

    You can also search for this author in PubMed Google Scholar

  8. Peter H. Sudmant
    View author publications

    You can also search for this author in PubMed Google Scholar

  9. José C. Brito
    View author publications

    You can also search for this author in PubMed Google Scholar

  10. Rasmus Nielsen
    View author publications

    You can also search for this author in PubMed Google Scholar

  11. Raquel Godinho
    View author publications

    You can also search for this author in PubMed Google Scholar

Contributions

R.G. and J.L.R. conceived and designed the study with R.N. R.G. and J.L.R. initiated the project. J.L.R., N.S., M.N., R.G., Z.B., J.C.B. and A.Q. collected the samples. J.C.B. and A.Q. provided logistical support in North Africa. J.L.R. performed DNA laboratory work and processed the sequencing data. S.A. performed RNA laboratory work. J.L.R. and N.S. processed samples for physiological analysis. R.N. oversaw all statistical and computational analyses, which were performed by J.L.R. (phenotypic and whole-genome data analyses) and P.S. (transcriptomics). P.H.S., R.G. and R.N. helped with editing and reviewing the paper for submission. J.L.R., R.G. and R.N. wrote the paper with input from all authors.

Corresponding authors

Correspondence to Joana L. Rocha, Rasmus Nielsen or Raquel Godinho.

Ethics declarations

Competing interests

The authors declare no competing interests.

Peer review

Peer review information

Nature Ecology & Evolution thanks the anonymous reviewers for their contribution to the peer review of this work. Peer reviewer reports are available.

Additional information

Publisher’s note Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Extended data

Extended Data Fig. 1 Mitogenomics of Vulpes using the raccoon dog as an outgroup (related to Fig. 1).

(a) Maximum-likelihood tree followed by species distribution range downloaded from IUCN (2021). Fox illustrations by Margarida Laranjeira Rocha. (b) Bayesian phylogenetic tree and time to the most recent common ancestor estimates; timescale in Ma represented bellow (see Supplementary Table 5 for the median time estimates and 95% HPD for each node). Accession numbers of published sequences are available in Supplementary Table 2. North African red foxes are coloured in orange, Eurasian red foxes green and Rueppell’s fox in red. Rueppell’s foxes are outgrouped by a non-coloured clade of 6 red foxes (n = 1 North Africa, n = 6 Middle East), revealing a paraphyletic pattern in red foxes.

Extended Data Fig. 2 ABBA-BABA on windows using the arctic fox (A-D) and the domestic dog (E-H) as reference genomes (related to Fig. 3).

D statistics between Rueppell’s fox (Vulpes rueppellii) and North African red fox (V. vulpes) (a–b, e–f) and between Rueppell’s fox and fennec fox (V. zerda) (c–d, g–h). (a, c, e, g) The pattern of ABBA/BABA in A and E is as follows: red fox Eurasia (H1) and North African red foxes (H2) have a recent common ancestor that diverged from the ancestral of the Rueppell’s fox (H3), using the arctic fox as outgroup (O) in A and the dog as outgroup (O) in E. The pattern of ABBA/BABA in C and G: all red foxes (H1) Rueppell’s foxes (H2) have a recent common ancestor that diverged from the ancestral of the fennec fox (H3), using the pale fox as outgroup (O) in C and the dog as outgroup (O) in E. Illustrations by Margarida Laranjeira Rocha. (b, d, f, h) Manhattan plot of the Z-transformed D statistic for non-overlapping windows of 200 kb across all chromosomes, according to pattern described in A, C, E, G respectively. Each dot represents a window and windows from the same chromosome have the same colour. Line threshold represents regions with Z score values above the 99th percentiles of the empirical distribution.

Extended Data Fig. 3 CRADD gene – top outlier above the 99th percentile of the empirical distribution for ƒd and Z(D) in ABBA-BABA on windows (related to Fig. 3).

a) Maximum-likelihood tree estimated with IQtree for a sequence alignment of pseudohaploids comprising the CRADD gene ( ~ 180 kb) based on sampling of a random base (samples from our dataset are highlighted in the legend of the figure, see Supplementary Table 2 for additional genomes with legend in the tree). (b) PCA with 48 red foxes and 24 Rueppell’s using 23,016 SNPs called from a 1.3 Mbp region containing CRADD; (c) Read depth-based estimate for CRADD per individual (depth for CRADD divided by the genome-wide coverage per individual). There is no evidence for CRADD copy number variation.

Extended Data Fig. 4 Principal component analysis for SLC6A16 (ranked 2nd) and two other genes (MGAM, HPS5) above the 99th percentile of ƒd empirical distribution (related to Fig. 3).

PCA based on sampling of a single read at each site. North African red foxes individuals (n = 30) unambiguously cluster closer to Rueppell’s (n = 24) than other red foxes (n = 18). Arrow highlights hidden North African red fox individuals clustering within Rueppell’s fox (MGAM). See also Figure S9 for gene tree topologies.

Extended Data Fig. 5 Tree topologies and Principal Component Analysis on the ~ 25 Mb outlier region above the 99th percentile of the empirical distribution for ƒd in ABBA-BABA on windows (related to Fig. 3).

Maximum-likelihood tree estimated with IQtree for the affected chromosomal region called with ANGSD (-doHaplo 1) using both (a) the dog genome (NC_006595.3, 38236300–63241923), and (b) the arctic fox as references (scaffolds 47, 126, 155, 200); (c) PCA of the 25 Mb region in 48 red foxes, 24 Rueppell’s foxes and 5 fennecs, based on sampling of a single read at each site. PC1 first separates Rueppell’s and fennec and PC2 separates Rueppell’s and red foxes.

Extended Data Fig. 6 Genetic Divergence for the 25 Mb region above the 99th percentile of the empirical distribution for ƒd on windows (related to Fig. 3).

Estimated genetic divergences (dXY) between different pairs of Vulpes species on sliding windows of 50 kb overlapped every 10 kb for the affected chromosomal region and for all the autosomes (excluding the CRADD gene). Boxes indicate upper and lower quartiles; centre line represents median; whiskers extend to minimum and maximum values within 1.5x interquartile range; points show outliers beyond whiskers. Sample sizes for pairs of species are as follows: n = 24 Rueppell’s fox, n = 48 Red fox.

Extended Data Fig. 7 Evidence for potential introgression of the 25 Mb region supported by KIT (related to Fig. 3).

(a) Tree topology for one of the four scaffolds comprising the 25 Mb region above the 99th percentile of the empirical distribution for ƒd in ABBA-BABA on windows; fennec and Rueppell’s are sister species. (b) Tree topology for KIT; Rueppell’s (n = 24) and red foxes (n = 48) become sister species, as in the average genome-wide tree. (c) Nucleotide diversity for the 25 Mb outlier region in 4 Vulpes species; KIT stands as being highly polymorphic in Rueppell’s fox, suggesting the simultaneous mapping of divergent alleles likely introduced via introgression and alleles located in the B chromosomes which are closer to red foxes. Dotted lines show the 99.95th percentile of the empirical distribution of nucleotide diversity for the region. The fact that KIT is not highly polymorphic in red foxes, which also have autosomal and B alleles mapped to the same region further confirms our hypothesis of shared structural variation between Rueppell’s fox and the fennec.

Extended Data Fig. 8 Gene expression differences between North African and Eurasian red foxes and Rueppell’s foxes (related to Fig. 4).

The plots represent top candidate genes putatively under selection in North African red foxes (CD163, SLC12A2) for which differential expression was significant between North African and Eurasian foxes (p-adjusted CD163 = 9.3e-08, p-adjusted SLC12A2 = 0.001), and/or North African red fox and Eurasian red foxes (p-adjusted CD163 = 2.5e-06; p-adjusted SLC12A2 = 0.079). The arctic fox transcriptome was used as reference.

Extended Data Fig. 9 Selection specific to Eurasian red foxes (related to Fig. 4).

Genome-wide selection scans on windows of 50 kb with a 10 kb slide across autosomal chromosomes, displayed as in the dog genome (x axis): population branch statistics (PBS) and alpha (negative logarithm of alpha) were estimated using Eurasian red foxes as focal group, North African red foxes as closest ingroup and Rueppell’s fox as outgroup; genes in grey are genes that stand out above the 99th percentile of the empirical distribution of f statistics in ABBA-BABA on windows for introgression between Rueppell’s and North African red foxes. The CRADD region and the 25 Mb region identified in Fig. 3 in main text were removed from the analysis. Lines show the 99.95th percentile of the empirical distribution. Names of genes within the highest peaks are shown (non-label peaks represent no overlap with genes). The arctic fox (V. lagopus) was used as reference genome. We note that other genes not displayed in the figure can overlap the outlier regions; a full list can be found in Supplementary Tables S16–S17.

Extended Data Fig. 10 Physiological differences between three species of Vulpes (related to Supplementary Tables 22–23 and Fig. 5).

Boxplots of 12 physiological parameters measured from the serum/plasma of red fox (V. vulpes; n = 10), desert-dwelling Rueppell’s fox (V. rueppellii; n = 13), and fennec fox (V. zerda; n = 10): (a) Copeptin (used as biomarker for vasopressin; pmol/L), (b) Aldosterone (pg/mL), (c) Total T4 (ug/dL), (d) Urea (mg/dL), (e) Uric acid (mg/dL), (f) Creatinine (mg/dL); (g) Albumin (g/dL), (h) Cholesterol (mg/dL), (i) Osmolality (mOSm/Kg), (j) Sodium (mEq/L), (k) Potassium (mEq/L), (l) Chloride (mEq/L); Juveniles were excluded. Boxes indicate upper and lower quartiles; centre line represents median; whiskers extend to minimum and maximum values within 1.5× interquartile range; points show outliers beyond whiskers. Additional pairwise t-tests comparing the three different fox species were performed with Bonferroni correction for multiple comparisons on each parameter, as well as Welch’s two-sample t-test comparing males versus females, North African versus European red foxes, and captive versus free-ranging individuals. All tests are two-sided. No significant differences between species were found except for A-C (pCopeptin = 0.00037; pAldosterone = 0.0034; pT4 = 8.8e-05) and between red fox and fennec (pCopeptin=0.02273; pAldosterone = 0.0027; pT4 = 1.0e-06). No significant differences were found between Rueppell’s and fennec for copeptin (AVP) and Aldosterone (PCopeptin= 0.55044; PAldosterone = 1.00; PT4 = 0.16), nor between captive and wild-caught individuals (PCopeptin=0.9811; PAldosterone = 0.3522; PT4 = 0.3867). After correcting for the potentially confounding factors caused by sex, weight, body mass index and whether the individuals were caught free-ranging or captive as covariates, we still found the same pattern of significance (Supplementary Table 23). This seems to suggest that, at least for these three hormones (A-C), there is no evidence of significant intraspecific differences in red foxes (n = 4 North Africa and n = 6 Iberia), and that the differences observed are not due to individuals being free-ranging or held in captivity, though this should be taken with caution given the small sample size (Supplementary Table 23).

Supplementary information

Supplementary Information

Supplementary Text, Figs. 1–18 and Tables 1–23.

Reporting Summary

Peer Review File

Supplementary Data

Metadata spreadsheet-based table with three tabs, each containing additional information to Supplementary Tables 1 and 3, including geolocation and morphological and physiological data.

Rights and permissions

Springer Nature or its licensor (e.g. a society or other partner) holds exclusive rights to this article under a publishing agreement with the author(s) or other rightsholder(s); author self-archiving of the accepted manuscript version of this article is solely governed by the terms of such publishing agreement and applicable law.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

L. Rocha, J., Silva, P., Santos, N. et al. North African fox genomes show signatures of repeated introgression and adaptation to life in deserts. Nat Ecol Evol 7, 1267–1286 (2023). https://doi.org/10.1038/s41559-023-02094-w

Download citation

  • Received:

  • Accepted:

  • Published:

  • Issue Date:

  • DOI: https://doi.org/10.1038/s41559-023-02094-w

This article is cited by

Search

Advanced search

Quick links

Nature Briefing

Sign up for the Nature Briefing newsletter — what matters in science, free to your inbox daily.

Get the most important science stories of the day, free in your inbox. Sign up for Nature Briefing