Efficient Maximum-Likelihood Inference For The Isolation-With-Initial-Migration Model With Potentially Asymmetric Gene Flow

Abstract

The isolation-with-migration (IM) model is commonly used to make inferences about gene flow during speciation, using polymorphism data. However, Becquet and Przeworski (2009) report that the parameter estimates obtained by fitting the IM model are very sensitive to the model's assumptions -including the assumption of constant gene flow until the present. This paper is concerned with the isolationwith-initial-migration (IIM) model of Wilkinson-Herbots (2012), which drops precisely this assumption. In the IIM model, one ancestral population divides into two descendant subpopulations, between which there is an initial period of gene flow and a subsequent period of isolation. We derive a very fast method of fitting an extended version of the IIM model, which also allows for asymmetric gene flow and unequal population sizes. This is a maximum-likelihood method, applicable to data on the number of segregating sites between pairs of DNA sequences from a large number of independent loci. In addition to obtaining parameter estimates, our method can also be used to distinguish between alternative models representing different evolutionary scenarios, by means of likelihood ratio tests. We illustrate the procedure on pairs of Drosophila sequences from approximately 30,000 loci. The computing time needed to fit the most complex version of the model to this data set is only a couple of minutes. The R code to fit the IIM model can be found in the supplementary files of this paper.

Hybrid asexuality as a primary reproductive barrier: on the interconnection between asexuality and speciation

“…The assembly of mRNA from five C. taenia specimens comprised 20,385 contigs (potential mRNAs) and was used as a reference for mapping of reads from all species (two individuals of each ingroup species and one outgroup). In total, we identified 187,205 SNPs (Appendix S4), which were then analysed by fitting eight coalescent models assuming different scenarios of species divergence and connectivity [36][37][38][39][40]. The models are characterized in Fig 4 and the results of analyses are described in Table 3.…”

“…This model fitted all pairwise species datasets significantly better than both other models of Wilkinson-Herbots [36,37], namely the isolation model I4 assuming no gene flow since the initial split of the species, and the isolation-with-migration model IM4 assuming ongoing hybridization from the initial split until the present (Likelihood Ratio Test; LRT p << 10 -10 for all comparisons involving C. elongatoides, and p < 0.01 for all other species comparisons). To test if the better fit of IIM7 is due to the correct assessment that interspecific gene flow ceased recently, or merely due to allowing an additional size change at time t 1 , we employed five more models that allow an additional population size change and also relax the assumption of equal population sizes during the migration stage (Fig 4; [ 39,40]).…”

“…All models were ranked according to their AIC score and we also calculated the evidence ratio for each model, providing a relative measure of how much less likely a given model is compared to the best-fitting model, given the set of candidate models considered and the data [69]; Table 3. In addition, Likelihood Ratio tests were performed to compare pairs of nested models, where we have assumed that the use of the χ2 distribution with the appropriate number of degrees of freedom is conservative ( [38,39,74]).…”

“…C. strumicae was used as outgroup, whose divergence time from the ingroup was set to 17.7 Mya according to the penalized likelihood time-calibrated tree of [35]. Computation of 95% confidence intervals was based on the profile likelihood as described in [39]. Because of the computational time required, such confidence intervals were provided only for the best-fitting model for each dataset.…”

“…Given that the reference transcriptome was prepared from C. taenia but more distant species were subsequently mapped onto it, our approach introducing 'N' on variant positions minimized the unbalanced mapping of more distant species onto the model transcriptome. Finally, we sorted identified SNPs from each locus into individual-locus-specific matrices.4.2 Detection of interspecific gene flow from SNP data was based on recently introduced coalescent-based maximum-likelihood methods[36][37][38][39][40] which are computationally inexpensive and estimate simultaneously the population sizes and migration rates as well as population splitting times, including some scenarios of time-variable migration rates. The following eight models, which are…”