Computational protocol: Significant loss of mitochondrial diversity within the last century due to extinction of peripheral populations in eastern gorillas

Similar protocols

Protocol publication

[…] We combined our data with published Grauer’s (N = 7) and mountain (N = 8) gorilla mitochondrial genomes, (Table ). Additionally, we included published mitochondrial genomes from western lowland gorillas (Gorilla gorilla gorilla, N = 43), a subspecies of the western gorillas with a large population size of ~140.000 individuals (Table ). This was done to allow cross-species comparison of genetic diversity. Birth year of individuals from published samples ranged from after 1990 to 2006 (Table ). Geographic coordinates of the sampling locations were obtained from the museum labels (historical samples), field records (fecal samples), or derived from publications (published sequences). All mitochondrial genomes were aligned with Clustal omega 1.2.4. Minimum spanning haplotype networks were constructed in PopART V1.7 (http://popart.otago.ac.nz). We grouped all mitochondrial genomes by sub-species and sample type (historical or modern sample) and calculated haplotype statistics, molecular diversity indices, and neutrality tests in Arlequin3.5. We assessed the effect of unequal sample size by using permutation tests. To this end, we randomly subsampled mitochondrial genomes from each population to match the population with the smallest sample size (modern mountain gorillas, n = 8) and repeatedly (n = 1000) recalculated haplotype and nucleotide diversity.To accurately model eastern gorilla demography, we estimated the gorilla specific mitochondrial mutation rate using BEAST 2.4.6. We used the estimated split date of western and eastern gorillas of ~1 million years ago obtained from nuclear data as calibration point,. To exclude confounding effects of different sample ages, only modern mitochondrial genomes were used. We first estimated genome wide mutation rate using the HKY-site model and enforcing a strict molecular clock. Birth-rate and clock-rate priors were set as gamma distribution with α = 0.001 and β = 1000. The prior for the split time between western and eastern gorillas was set as a log-normal distribution with M = 0.5 and S = 0.1, which corresponds to a divergence time of 0.848–1.18 Mya (95% CI). The Bayesian model was run for an MCMC length of 500 million and we used Tracer 1.6 to confirm run convergence and obtain probability distributions. We then partitioned the aligned mitochondrial genomes into coding genes (splitting triplets into 1st + 2nd and 3rd position), rRNAs, tRNAs, and non-coding regions in Geneious 10.1.2 using the annotated gorilla reference. The same BEAST model was then used to estimate partition-specific mutation rates. The obtained mutation rates were highly similar to those published for the human mitochondrial genome, (Fig. S, Table ).We used the Bayesian skyline model in BEAST 2.4.6 to infer demographic changes within the last centuries in Grauer’s and mountain gorillas. The best partition scheme and the best fitting model for each partition were identified with PartitionFinder 2, only considering models available in the BEAST 2 software package. The skyline analysis was then run for both subspecies separately under a strict molecular clock model using the inferred gorilla mutation rates for the different partitions and MCMC length of 500 million.To account for the possibility of more recent demographic changes that could be bypassed by the skyline model, we also employed an approximate Bayesian computation approach. Since initial exploratory analyses failed to detect the onset of demographic change within the last 2000 years due to the lack of power over this large time period, we narrowed down the time period to the last 110 years. This time frame reflects the temporal distribution of our historical samples and encompasses the documented population size decline within the last four decades in Grauer’s and mountain gorillas,. We divided our dataset into two temporal groups: present (modern samples) and historical (museum specimens).We performed two inferences: (1) a model choice analysis with three competing models in which effective population size was constant, decreased or increased in a single recent event (30–110 years ago); and (2) a parameter estimation of two effective population sizes—present and historical. We ran 30 million simulations (10 million per model) for the model choice analysis and employed 10 acceptance proportions (0.005–0.05%) for assessing consistency of the estimated model likelihoods. For the parameter estimation analysis, we ran 1 billion simulations and applied rejection with a linear regression adjustment over the 0.001% (=10,000) accepted simulations. In a dedicated analysis, we selected the following summary statistics: number of haplotypes, number of private haplotypes and genetic diversity for each of the two temporal groups as well as pairwise differences and FST between the two temporal groups.Finally, to discriminate between negative results that were due to lack of demographic change versus those due to limited information content of the data, we employed pseudo-observed datasets (PODs) to assess the statistical power of our samples and methodology. All analyses were performed separately for Grauer’s and Mountain gorilla mitochondrial genomes in the software BaySICS v1.9 using the estimated gorilla mutation rate (1.28 × 10−8), no gamma parameter, and a Ts/Tv bias of 0.75 and 0.933 in mountain and Grauer’s gorillas, respectively (as identified by the selection model of molecular evolution carried out in MEGA 7). […]

Pipeline specifications