Aim: Understanding patterns and processes of geographic genetic variation within and among closely related species is the essence of phylogeography. Japanese macaques, also called snow monkeys, have been extensively studied, particularly in the fields of sociobiology, ecology, and experimental biology; however, our knowledge of their evolutionary history is relatively limited. In this study we aimed to elucidate the geographic patterns of genetic variation in Japanese macaques and the processes that underlie them.

Location: Japan

Taxon: Japanese macaque, Macaca fuscata; rhesus macaque, M. mulatta; Taiwanese macaque, M. cyclopis

Methods: Double-digest restriction-site associated DNA (RAD) sequencing was used to identify genome-wide single nucleotide variants. We used fineRADstructure, ADMIXTURE, and principal component analyses to estimate the genetic population structure. Phylogenetic relationships were then inferred based on neighbour-net, neighbour-joining, maximum likelihood, and SVDquartets algorithms. We assessed gene flow using demographic inference and ABBA-BABA tests, and estimated past distributions during the Last Glacial Maximum (LGM) using ecological niche modelling.

Results: Japanese macaques show a sister group relationship with a clade comprising Chinese rhesus, Indian rhesus, and Taiwanese macaques. Japanese macaques comprise major northeastern and southwestern clades, with a boundary located near central Japan, and gene flow between the northeastern and southwestern lineages was detected. Refugia during the LGM were estimated to be distributed in limited areas along the south coasts of the Japanese archipelago.

Main Conclusions: Phylogeographic variation of Japanese macaques is likely due mainly to northeast–southwest divergence, which resulted from withdrawal into refugia during the glacial period, and subsequent gene flow.

We used 113 blood- or tissue-derived DNA samples from 10 populations of Macaca fuscata (N = 102), the Chinese and Indian populations of M. mulatta (N = 6), M. cyclopis (N = 2), and M. fascicularis (N = 3). The ddRAD library preparation followed the quaddRAD procedure of Franchini et al. (2017) with minor modifications.

fasta_chr.zip: consensus sequences of each chromosome in fasta format for the calculation of nucleotide diversity (the outputs from Stacks populations function with --fasta-samples option)

fastsimcoal2.zip: data and model settings for fastsimcoal2 demographic analyses
    data: five replicates of multidimensional folded SFS
    models: template (.tpl) and estimation (.est) files of eight demographic models

paup.zip: nexus files for phylogenetic inference using PAUP
    nj_a.nexus: autosomal SNVs for NJ phylogenetic inference
    nj_x.nexus: X-chromosome SNVs for NJ phylogenetic inference
    nj_y.nexus : Y-chromosome SNVs for NJ phylogenetic inference
    svdq_ind.nexus: autosomal SNVs for individual-level SVDquartets phylogenetic inference
    svdq_pop.nexus: autosomal SNVs for population-level SVDquartets phylogenetic inference

populations.haps.radpainter: haplotype data for fineRADstructure analysis

raxml.phy: autosomal SNVs in phylip format for ML phylogenetic inference using RAxML-NG

sex_pop_cluster_list.txt: tab-delimited text file of id, sex (female or male), population, and cluster

splitstree4.dist: uncorrected p-distance data for phylogenetic network analysis using SplitsTree4 (the outputs from PAUP)

treemix.zip: allele frequency data for TreeMix analysis
    treemix_full.frq.gz: full dataset
    treemix_sub.frq.gz: subset with no missing data

vcf.zip: vcf files for population structure analyses (used after only Japanese macaque samples were extracted)
    a.vcf: autosomal SNVs for ADMIXTURE and PCA
    x.vcf: X-chromosome SNVs for ADMIXTURE