multi-dice: r package for comparative population genomic inference under hierarchical co-demographic models of independent single-population size changes

Abstract

Population genetic data from multiple taxa can address comparative phylogeographic questions about community-scale response to environmental shifts, and a useful strategy to this end is to employ hierarchical co-demographic models that directly test multi-taxa hypotheses within a single, unified analysis. This approach has been applied to classical phylogeographic data sets such as mitochondrial barcodes as well as reduced-genome polymorphism data sets that can yield 10,000s of SNPs, produced by emergent technologies such as RAD-seq and GBS. A strategy for the latter had been accomplished by adapting the site frequency spectrum to a novel summarization of population genomic data across multiple taxa called the aggregate site frequency spectrum (aSFS), which potentially can be deployed under various inferential frameworks including approximate Bayesian computation, random forest and composite likelihood optimization. Here, we introduce the r package multi-dice, a wrapper program that exploits existing simulation software for flexible execution of hierarchical model-based inference using the aSFS, which is derived from reduced genome data, as well as mitochondrial data. We validate several novel software features such as applying alternative inferential frameworks, enforcing a minimal threshold of time surrounding co-demographic pulses and specifying flexible hyperprior distributions. In sum, multi-dice provides comparative analysis within the familiar R environment while allowing a high degree of user customization, and will thus serve as a tool for comparative phylogeography and population genomics.

Keywords: aggregate site frequency spectrum; approximate Bayesian computation; comparative phylogeography; population genetics software; random forest.

Figure 1

Hierarchical co‐demographic models. (a) Example instantaneous co‐expansion model. (b) Example instantaneous co‐contraction model. Both models are such that eight of the ten taxa are assigned to three synchronous co‐demographic pulses (ѱ = 3; ζ_T = 0.8), with the first pulse containing three taxa (ζ₁ = 0.3), the second pulse containing another two taxa (ζ₂ = 0.2) and the third pulse containing yet another three taxa (ζ₃ = 0.3). Pulse 1 occurs at the most recent time (τ_s,1), pulse 2 occurs at the intermediate time (τ_s,2), and pulse 3 occurs at the most ancient time (τ_s,3). The remaining two taxa are then behaving idiosyncratically in time from all other taxa (τ_i,1 and τ_i,2). Each taxon is allowed nuisance demographic parameter draws independent from each other ({ε₁, …, ε₁₀} and {N ₁, …, N ₁₀})

Figure 2

Flowchart of multi‐dice usage. multi‐dice accomplishes multi‐taxa co‐demographic inference under a hierarchical model through three major steps: model specification, single‐population simulation across multiple taxa and conversion of simulated data to multi‐taxa summary statistics. Hierarchical co‐demographic model specification is conducted across multiple functions in sequence, with preceding functions contained within successive functions. This sequential embedding of functions extends to dice.sims(), allowing the entire model specification process to be performed concurrently with data simulation. Simulated data can then be converted to multi‐taxa summary statistics by either dice.aSFS() or dice.sumstats(), depending on the data type. Additionally, these functions can be applied to empirical data as well. To clarify, only two multi‐dice functions/command lines, dice.sims() and dice.aSFS()/dice.sumstats(), are needed for simplest usage to construct a reference table of multi‐taxa summary statistics under a hierarchical co‐demographic model. This reference table can then be exploited in a downstream software program for hRF or hABC purposes, where appropriate statistical practices should be used to examine robustness and fit. Importantly, exploratory analyses should be performed on the empirical data prior to deploying multi‐dice to better guide its usage, for example, to determine sensible prior distributions and evaluate differences among taxa