Published online 13 May 2019 Nucleic Acids Research, 2019, Vol. 47, Web Server issue W283–W288 doi: 10.1093/nar/gkz372 iMKT: the integrative McDonald and Kreitman test Jesús Murga-Moreno† , Marta Coronado-Zamora† , Sergi Hervas, Sònia Casillas and Antonio Barbadilla * Institut de Biotecnologia i de Biomedicina and Departament de Genètica i de Microbiologia, Universitat Autònoma de Barcelona, 08193 Bellaterra, Barcelona, Spain Received March 06, 2019; Revised April 18, 2019; Editorial Decision April 29, 2019; Accepted May 03, 2019 ABSTRACT compute the quotient of the number of non-synonymous (dN ) to synonymous (dS ) substitutions, the MKT uses both The McDonald and Kreitman test (MKT) is one of polymorphic and divergence data. Polymorphic data al- the most powerful and widely used methods to de- lows taking into account purifying selection on divergent tect and quantify recurrent natural selection using non-synonymous sites, significantly increasing the detection DNA sequence data. Here we present iMKT (acronym power of recurrent positive selection. The MKT covers the for integrative McDonald and Kreitman test), a novel evolutionary period spanning from the divergence of the web-based service performing four distinct MKT outgroup species to the present. The null model of MKT types. It allows the detection and estimation of is the neutral hypothesis (4,5). Because infrequent adaptive four different selection regimes −adaptive, neutral, mutations fix fast relatively to common neutral mutations, strongly deleterious and weakly deleterious− acting they contribute almost exclusively to divergence and not to on any genomic sequence. iMKT can analyze both polymorphism; therefore, an excess of the divergence ratio relative to the polymorphism ratio is the signal of positive user’s own population genomic data and pre-loaded selection. The fraction of adaptive nonsynonymous substi- Drosophila melanogaster and human sequences of tutions (α) can be estimated from the MKT data (6,7). protein-coding genes obtained from the largest pop- The main drawback of MKT is that it assumes strict neu- ulation genomic datasets to date. Advanced options trality of segregating sites. Because weak negative selection in the website allow testing complex hypotheses abounds in the genomes (1), α estimates are biased down- such as the application example showed here: do ward. Several MKT methodological extensions try to cor- genes located in high recombination regions un- rect the bias. In Appendix 1, four MKT approaches are dergo higher rates of adaptation? We aim that iMKT listed: (i) the standard (original) MKT (2); (ii) the Fay, will become a reference site tool for the study of Wyckoff and Wu correction (FWW MKT) (8); (iii) the ex- evolutionary adaptation in massive population ge- tended MKT (eMKT) (9) and (iv) the asymptotic MKT nomics datasets, especially in Drosophila and hu- (aMKT) (10). Each method has pros and cons as discussed in Appendix 1, and for the comparison of their different mans. iMKT is a free resource online at https://imkt. outputs, it would be very convenient to have a web service uab.cat. to perform at once the four MKTs. Existing web servers compute either the standard MKT (11,12) or more recently the aMKT (13). None of them contains pre-loaded popula- INTRODUCTION tion genomics data of representative species as Drosophila One of the most striking evidence of the power of natural melanogaster or humans. selection is the characteristic footprints that it leaves on the Here we present iMKT (acronym for integrative McDon- patterns of genetic variation. A growing number of statis- ald and Kreitman test), a web-based service performing tical methods to analyze genomic data allows us to detect the four MKT types described in Appendix 1 and Figure and quantify adaptation and other selection regimes in the 1. It allows the detection and estimation of four selection genome at different temporal scales (reviewed in 1). regimes (adaptive, neutral, strongly deleterious and weakly The McDonald and Kreitman test (MKT, 2) is one of deleterious) acting on protein-coding DNA sequences. The the most powerful and robust methods we have to detect benefit of this tool is fourfold. the action of natural selection at the DNA level. MKT tests for the presence of recurrent positive (adaptive) selection on a gene or genome region. Unlike the ω = dN /dS ra- 1. Four MKTs, two of which were not previously available tio (3), which uses only divergence data among species to as open software packages, can be performed at once to * To whom correspondence should be addressed. Tel: +34 935868941; Fax: +34 935812011; Email: antonio.barbadilla@uab.es † The authors wish it to be known that, in their opinion, the first two authors should be regarded as joint First Authors.  C The Author(s) 2019. Published by Oxford University Press on behalf of Nucleic Acids Research. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted reuse, distribution, and reproduction in any medium, provided the original work is properly cited. W284 Nucleic Acids Research, 2019, Vol. 47, Web Server issue analyze user’s own population genomic data in a simple interface offered by a web-based service. 2. It allows the simultaneous comparisons of the results of the different MKTs, which behave differently according to different properties of the data. 3. Taking advantage of the copious information gathered in previous population genome browsers, PopFly (14) and PopHuman (15), it offers a fast tool to estimate the differ- ent selective regimes on thousands of D. melanogaster and human protein-coding genes on several worldwide popu- lations. 4. It allows comparing the selective regimes of a set of cod- ing genes (selected according to the user’s criterion, such as recombination rate bins or chromosome localization) with those of the genome-wide distribution in both hu- mans and D. melanogaster. The incessant accumulation of massive genome data makes this website a timely resource to describe and quan- tify natural selection for any biological species at the genome level. MATERIALS AND METHODS Input data The iMKT server can analyze both user’s own population data and pre-loaded data of D. melanogaster or human protein-coding genes. In the first case, the user can upload as input either poly- morphism and divergence data or aligned multi-FASTA files. For polymorphism and divergence data, the user must upload two files: (i) a tab-delimited file containing the dis- tribution of Derived Allele Frequencies (DAF) (16) of all segregating (polymorphic) variants for two types of sites (putatively under selection and putatively neutral) and (ii) a file containing the counts of divergent positions for the two site types. For aligned multi-FASTA files, the user needs to enter one or more files containing aligned protein-coding sequences for at least two sequences of the same species to estimate polymorphism counts, and one orthologous se- quence from an outgroup species to estimate divergence and ←−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−−− mous and non-synonymous classes for a gene exhibiting an excess of both slightly deleterious and fixed non-synonymous differences with n = 10 sampled chromosomes. (B) The standard MKT for this gene (P-value = 0.09, 2 × 2 Fisher’s exact test). (C) The 2 × 2 table by Fay, Wyckoff and Wu’s correction (24) taking into account only polymorphism found on more than one chromosome (P-value = 0.045, 2 × 2 Fisher’s exact test). (D) Ex- tended MKT (9). The count of segregating sites in non-synonymous sites is partitioned into the number of neutral variants and the number of weakly deleterious variants. PN is substituted with the number of nonsynonymous polymorphisms that is neutral (P-value = 0.042, 2 × 2 Fisher’s exact test). (E) Asymptotic MKT. Example of the result of asymptotic MKT using D. melanogaster 2R chromosome and D. simulans as outgroup. The two vertical lines show the limits of the x cutoff interval used (in the exam- ple [0,0.9]). Black dots indicate the binned values for each DAF category. The solid red curve shows the fitted fit(x). The dashed red line is the final asymptote. The dark gray band indicates the 95% CI around the estima- tion. The blue dashed line shows the estimated using the standard MKT for comparison. For MKT methods definitions, see Appendix 1. Adapted and expanded from 29. Figure 1. Comparison of the four MKT methods implemented in iMKT. (A) The hypothetical derived allele frequency (DAF) spectrum of synony- Nucleic Acids Research, 2019, Vol. 47, Web Server issue W285 infer ancestral alleles. Examples of such files are provided at inbred individuals (22) and non-accessible nucleotides (19) the website. were discarded following the PopHuman methodology (15). For analyzing D. melanogaster or human protein-coding Genome annotations were retrieved from GENCODE (re- genes, the user can use the population genomic data avail- lease 27). Recombination rate values associated with each able in the web server. In this case, the user can either sub- protein-coding gene were obtained from Bhérer (23) and mit a list of protein-coding genes or select them from the correspond to the sex-average estimates. list provided, and select the population(s) and preferred method(s) to analyze the selective regimes on a group of Estimation of the number of synonymous and nonsynonymous protein-coding genes. changes Population genetics pipeline for D. melanogaster and human Inferring the action of natural selection on coding se- data quences relies on the computation of polymorphism and di- vergence data on two distinct types of sites in the genome: We have designed and implemented a custom pipeline for one putatively selected (usually non-synonymous coding analyzing the Drosophila Genome Nexus (17,18) and Hu- sites), and one putatively neutral (usually synonymous cod- man 1000GP Phase III (19) data, which could potentially be ing sites) (2). This implies assigning a selective class for each escalated to any available genomic data source. The pipeline nucleotide site in the genome. This task is not trivial when pre-calculates the DAF and number of divergent synony- different transcripts overlap a genomic region. For example, mous and nonsynonymous sites, which are needed to fur- one nucleotide site can be a non-synonymous site for one ther perform on-the-fly MKTs. A total of 13 753 protein- transcript but a synonymous site for another nested gene coding genes for 16 D. melanogaster populations (17,18) transcript. In these cases, we assign the most selective con- and 20 643 protein-coding genes for 26 human popula- strained class to the nucleotide site. In the example, the site tions of distinct geographical origin (19) were analyzed. is considered non-synonymous. Pre-calculated DAF and divergence values are stored in the server. The complete pipeline is available as a Jupyter Note- Exclusion of low-frequency variants book at https://github.com/BGD-UAB/iMKTData to al- low its reproducibility. Slightly deleterious variants are mainly segregating at low frequency (8,24,25). These rare polymorphisms can be ex- Data retrieval cluded from the analyses by specifying one or several threshold frequency values depending on the FWW MKT, the D. melanogaster population genomic data. Variation data eMKT or the aMKT method. In addition, the aMKT al- generated by the Drosophila Genome Nexus, together with lows removing high-frequency variants that might be due divergence data between D. melanogaster and D. simulans, to polarization errors (10,13). was retrieved from PopFly (14) in FASTA format. Only populations with at least four genome sequences and less Statistical analysis than 20% of missing or ambiguous nucleotides each (after filtering by identity by descent, admixture and heterozygos- For analyses including several protein-coding genes, users ity) were included. DAF spectrum by functional classes was are recommended to select the option Concatenate genes. In estimated by resampling a number of lines with nucleotide this case, iMKT analyzes the selective regimes for the whole information (excluding undetermined sites, N bases) at each gene set instead of for each gene separately and applies a position without replacement. This procedure maximizes statistical test of heterogeneity of the selection acting among the number of informative sites to analyze. The number of the analyzed genes (Cochran-Mantel-Haenszel statistic). In lines resampled for each population was chosen depending addition, the iMKT web server allows performing statistical on the number of lines sequenced and the quality of those enrichment analyses to assess whether a group of genes is ei- sequences (Supplementary Table S1). Positions and genes ther enriched or depleted of positively selected genes when without valid information for at least this defined number compared to the complete genome distribution or to a sec- of lines were discarded for the analysis. The ancestral state ond group of genes submitted by the user. In this case, the of each polymorphic site was inferred from the comparison user should choose also the option Compare against whole- with the outgroup species D. simulans. The genome refer- genome distribution or Compare against a second dataset. ence sequence and annotations correspond to the 5.57 Fly- A resampling 95% confidence interval (CI) is generated by Base release (20). Gene-associated recombination rate for estimating α with the chosen MKT test for 100 bootstrap 100 kb non-overlapping windows were retrieved from Com- replicates by sampling genes with replacement within each eron et al. (21). group. In the asymptotic MKT, 95% CI intervals around the α estimation are already provided in the output. Human population genomic data. Genome variation data and ancestral state of variants generated by the 1000GP Output Phase III (19), together with divergence estimates between humans and chimpanzees, were retrieved from PopHuman The output of iMKT is an extensive report displayed as an (15) in Variant Call Format (VCF). The dataset included HTML page. It contains several sections, starting with a 84.4 million variants detected across 2504 individuals from summary table with the input parameters, a table with de- 26 different populations, which were mapped to the hu- scriptive statistics, and the standard MKT table. Finally, the man reference genome version GRCh37/hg19. Reportedly tests selected by the user are displayed below. W286 Nucleic Acids Research, 2019, Vol. 47, Web Server issue PRACTICAL GUIDE TO THE iMKT WEBSITE tutorial from the main page). The website also contains sam- ple files for each available type of analysis and links to re- The iMKT site allows performing four MK-derived tests as lated resources such as PopFly, PopHuman, and the iMKT a web-based service. The website is divided into different R package. sections, each of which allows performing different types of analyses. APPLICATION EXAMPLE OF THE iMKT WEBSITE MKT analysis The iMKT website is designed to help testing evolutionary This page allows performing diverse MK-derived tests and hypotheses from a population genetics perspective. The on- estimating different selective regimes in your own data. The line tutorial, apart from guiding you in the usage of this input can either be polymorphism and divergence data in resource, contains some worked-out cases that can be ad- two separate files, as described in the Methods section, dressed using iMKT. In the application example developed or protein-coding sequences as aligned multi-FASTA files. here, we want to assess whether recombination rate limits When a multi-FASTA file is uploaded, the server outputs the adaptive potential of protein-coding genes. The specific the DAF spectrum and the divergence calculations, which hypothesis is that genes located in high recombination re- can be downloaded by the user and used in subsequent anal- gions undergo higher rates of adaptation. To test the hy- yses. Note, however, that the former input type gives more pothesis, we start by entering the PopFly data analysis page flexibility to analyze any sort of functional site. As an exam- of iMKT. Next, we use the filtering options below the ta- ple, you might want to test for selection at nonsynonymous ble to select 475 genes having a recombination rate higher coding sites (N) compared to synonymous coding sites (S) than 7 cM/Mb (Min recombination rate: 7). After selecting as the classical MKT was formulated, or to test for selection the genes, we select one or more populations (United States at Conserved Noncoding Sequences (CNS, N) compared to (RAL)) and an MKT test (eMKT). Finally, we choose the non-CNS (S) (26), etc. The choices are unlimited according option Compare against whole-genome distribution, which to the user’s needs. compares the distribution of α for the selected 475 genes located in regions of high recombination against the cor- responding distribution for all D. melanogaster genes. As PopFly/PopHuman data analysis part of an extensive output report, an illustrative box plot If you want to analyze D. melanogaster or human protein- shows a pronounced difference in the level of adaptation coding data, iMKT contains readily available variation data (α) between genes located in regions of high recombina- obtained from the largest genome variation datasets in tion (blue; α mean = 0.602; ±SD = 0.032) and all 13,753 each species (see Methods). The first step is to select the D. melanogaster genes (orange; m α mean = 0.44; ±SD = genes to be analyzed in the table displaying all the avail- 0.055) (Figure 2A). able genes. Genes are identified by either the Gene sym- We can repeat the same procedure for the D. melanogaster bol or the FlyBase/Ensembl ID. Genes in the table can be ancestral population from Zambia (Zambia (ZI)). As pre- sorted/filtered by chromosome and recombination rate, in viously, the output report uncovers a much higher level of addition to the Gene symbol and Flybase/Ensembl ID. In adaptation (α) in genes located in regions of high recombi- case the user needs to analyze a specific list of genes that nation (blue; α mean = 0.633; ±SD = 0.028) compared to cannot be easily filtered from the provided table (e.g. genes the total 13 753 D. melanogaster genes (orange; α mean = related to a specific pathway, as obtained from a search in 0.457; ±SD = 0.053) (Figure 2A). KEGG (27), a list with those genes, were genes are identified Finally, the same analysis in humans for a colonizing by symbol or FlyBase/Ensembl ID, can be uploaded. Sec- population (Utah Residents (CEU)) and an ancestral pop- ond, one or more populations on which to perform the anal- ulation (Yoruba (YRI)) reveals negative α adaptation val- ysis need to be specified. Third, one or more MK-derived ues in most cases and differences between the two groups methods can be chosen. Finally, advanced options are avail- of genes compared (Figure 2B). The results of this straight- able to analyze all the genes as a group instead of ana- forward analysis show that: (i) D. melanogaster undergoes lyzing them separately (option Concatenate genes), and to higher rates of adaptation than humans and (ii) genes lo- compare the results of this gene set against all the genes cated in regions of high recombination undergo higher rates of the genome (option Compare against whole-genome dis- of adaptation in both D. melanogaster and humans. tribution) or against a second group of genes provided by The example application developed here illustrates the the user (option Compare against a second dataset). Poten- power of iMKT to reveal new knowledge about evolu- tial applications include analyzing a single protein-coding tionary processes in Drosophila and humans without the gene or exploring different selective regimes in genes that need for labor-intensive data retrieval and/or processing are expressed tissues, anatomic structures, or developmen- by the user. The wide range of potential queries that tal stages (28)). can be performed using the searching capabilities of the iMKT website remarkably facilitates comprehensive analy- ses of evolutionary adaptation and constraint, even for non- Other sections of the website bioinformaticians. As such, iMKT is a comprehensive ref- The iMKT website includes extensive methodological and erence site for the study of protein adaptation in massive technical documentation (see the section Documentation in population genomics datasets, especially in Drosophila and the website), as well as a complete tutorial on the usage of humans. Finally, we want to emphasize that the flexibility of iMKT, with step-by-step examples (see the section Help and iMKT to input custom data allows analyzing diversity data Nucleic Acids Research, 2019, Vol. 47, Web Server issue W287 ACKNOWLEDGEMENTS We thank Esteve Sanz for help with the informatics infras- tructure in which iMKT is implemented. We also thank two anonymous referees for very helpful comments on the iMKT implementation and manuscript. FUNDING Ministerio de Economı́a y Competitividad (Spain) [CGL2017-89160P to M.S., A.B.]; AGAUR (Generali- tat de Catalunya) [2017SGR-1379 to A.R.]; Secretaria d’Universitats i Recerca del Departament d’Empresa i Coneixement (Generalitat de Catalunya) [FI-DGR2015 to M.C.-Z.]; Departament de Genètica i de Microbiologia (UAB) [PIF to J.M.-M.]; Servei de Genòmica i Bioin- formàtica de la UAB. Funding for open access charge: Ministerio de Economı́a y Competitividad (Spain). Conflict of interest statement. None declared. REFERENCES 1. Casillas,S. and Barbadilla,A. (2017) Molecular population genetics. Genetics, 205, 1003–1035. 2. McDonald,J.H. and Kreitman,M. (1991) Adaptive protein evolution at the Adh locus in drosophila. Nature, 351, 652–654. 3. Kimura,M. (1977) Preponderance of synonymous changes as evidence for the neutral theory of molecular evolution. Nature, 267, 275–276. 4. Kimura,M. (1968) Evolutionary rate at the molecular level. Nature, 217, 624–626. 5. Kimura,M. (1983) The Neutral Theory of Molecular Evolution. Cambridge University Press, Cambridge. 6. Charlesworth,B. (1994) The effect of background selection against deleterious mutations on weakly selected, linked variants. Genet. Res., 63, 213–227. 7. Smith,N.G.C. and Eyre-Walker,A. (2002) Adaptive protein evolution Figure 2. iMKT graphical output of an application example. Sampling dis- in drosophila. Nature, 415, 1022–1024. tribution of α values for protein-coding genes located in regions of high 8. Fay,J.C., Wyckoff,G.J. and Wu,C.-I. (2001) Positive and negative recombination (recombination rate >7 cM/Mb) compared to all protein- selection on the human genome. Genetics, 158, 1227–1234. coding genes in the genome for (A) the D. melanogaster Raleigh (RAL) 9. Mackay,T.F.C., Richards,S., Stone,E.A., Barbadilla,A., Ayroles,J.F., population (blue) and the D. melanogaster Zambia (ZI) population (yel- Zhu,D., Casillas,S., Han,Y., Magwire,M.M., Cridland,J.M. et al. low) and (B) the human Utah Residents (CEU) population (blue) and the (2012) The drosophila melanogaster genetic reference panel. Nature, human Yoruba (YRI) population (yellow). The distribution was calculated 482, 173–178. by randomly sampling 400 genes 100 times from the two lists of genes with 10. Messer,P.W. and Petrov,D.A. (2013) Frequent adaptation and the replacement and estimating α in each bin. Polymorphisms with a frequency McDonald-Kreitman test. Proc. Natl. Acad. Sci. U.S.A., 110, below 0.05 in the analyzed population were discarded (see main text). 8615–8620. 11. Egea,R., Casillas,S. and Barbadilla,A. (2008) Standard and generalized McDonald-Kreitman test: a website to detect selection by comparing different classes of DNA sites. Nucleic Acids Res., 36, outside protein-coding regions. This expands, even more, W157–W162. the hypotheses that can be tested and makes iMKT a key 12. Vos,M., te Beek,T.A.H., van Driel,M.A., Huynen,M.A., tool to test for recurrent adaptation in the genome of any Eyre-Walker,A. and van Passel,M.W.J. (2013) ODoSE: a webserver species. for genome-wide calculation of adaptive divergence in prokaryotes. PLoS One, 8, e62447. 13. Haller,B.C. and Messer,P.W. (2017) asymptoticMK: a web-based tool for the asymptotic McDonald–Kreitman test. G3: Genes Genomes Genetics, 7, 1569–1575. DATA AVAILABILITY 14. Hervas,S., Sanz,E., Casillas,S., Pool,J.E. and Barbadilla,A. (2017) iMKT is a free resource online, open to all users without lo- PopFly: the Drosophila population genomics browser. Bioinformatics, 33, 2779–2780. gin requirement at https://imkt.uab.cat. The corresponding 15. Casillas,S., Mulet,R., Villegas-Mirón,P., Hervas,S., Sanz,E., R package is available for download at https://github.com/ Velasco,D., Bertranpetit,J., Laayouni,H. and Barbadilla,A. (2018) BGD-UAB/iMKT. PopHuman: the human population genomics browser. Nucleic Acids Res., 46, D1003–D1010. 16. Ronen,R., Udpa,N., Halperin,E. and Bafna,V. (2013) Learning natural selection from the site frequency spectrum. Genetics, 195, SUPPLEMENTARY DATA 181–193. 17. Lack,J.B., Cardeno,C.M., Crepeau,M.W., Taylor,W., Supplementary Data are available at NAR Online. Corbett-Detig,R.B., Stevens,K.A., Langley,C.H. and Pool,J.E. (2015) W288 Nucleic Acids Research, 2019, Vol. 47, Web Server issue The Drosophila genome Nexus: A population genomic resource of Fay, Wyckoff and Wu correction (FWW MKT) 623 Drosophila melanogaster genomes, including 197 from a single ancestral range population. Genetics, 199, 1229–1241. Because slightly deleterious variants tend to segregate at 18. Lack,J.B., Lange,J.D., Tang,A.D., Corbett-Detig,R.B. and Pool,J.E. lower frequencies than do neutral mutations, Fay, Wyck- (2016) A thousand fly genomes: An expanded Drosophila genome off and Wu or FWW correction (Figure 1C) propose to re- nexus. Mol. Biol. Evol., 33, 3308–3313. 19. Auton,A., Abecasis,G.R., Altshuler,D.M., Durbin,R.M., move low-frequency polymorphisms from the analysis (8). Abecasis,G.R., Bentley,D.R., Chakravarti,A., Clark,A.G., α is estimated using the standard MKT equation but con- Donnelly,P., Eichler,E.E. et al. (2015) A global reference for human sidering only those polymorphic sites (for both neutral and genetic variation. Nature, 526, 68–74. selected classes) with a frequency above an established cut- 20. Thurmond,J., Goodman,J.L., Strelets,V.B., Attrill,H., Gramates,L.S., off. Charlesworth and Eyre-Walker (25) showed that even Marygold,S.J., Matthews,B.B., Millburn,G., Antonazzo,G., Trovisco,V. et al. (2019) FlyBase 2.0: the next generation. Nucleic removing low-frequency variants, the estimate of α is still Acids Res., 47, D759–D765. downwardly biased. Only these estimates are reasonably ac- 21. Comeron,J.M., Ratnappan,R. and Bailin,S. (2012) The many curate when the rate of adaptive evolution is high and the landscapes of recombination in Drosophila melanogaster. PLoS distribution of fitness effects of slightly deleterious muta- Genet., 8, e1002905. 22. Gazal,S., Sahbatou,M., Babron,M.-C., Génin,E. and tions is leptokurtic (because leptokurtic distributions have a Leutenegger,A.-L. (2015) High level of inbreeding in final phase of smaller proportion of polymorphisms that are slightly dele- 1000 Genomes Project. Sci. Rep., 5, 17453. terious). 23. Bhérer,C., Campbell,C.L. and Auton,A. (2017) Refined genetic maps reveal sexual dimorphism in human meiotic recombination at multiple scales. Nat. Commun., 8, 14994. 24. Templeton,A.R., Hartl,D.L. and Wu,C.-I. (1996) Contingency tests of neutrality using intra/interspecific gene trees: the rejection of Extended MKT (eMKT) neutrality for the evolution of the mitochondrial cytochrome oxidase Mackay et al. (9) proposed the extended MKT (Figure 1D). II gene in the hominoid primates. Genetics, 144, 1263–1270. 25. Charlesworth,J. and Eyre-Walker,A. (2008) The Instead of simply removing low-frequency polymorphism McDonald-Kreitman test and slightly deleterious mutations. Mol. below a given threshold, the count of segregating sites in Biol. Evol., 25, 1007–1015. non-synonymous sites is partitioned in the number of neu- 26. Casillas,S., Barbadilla,A. and Bergman,C.M. (2007) Purifying tral variants (using neutral sites as a proxy) and the num- selection maintains highly conserved noncoding sequences in drosophila. Mol. Biol. Evol., 24, 2222–2234. ber of weakly deleterious variants. This increases the power 27. Kanehisa,M., Sato,Y., Furumichi,M., Morishima,K. and Tanabe,M. of detecting adaptive selection (as it does not remove as (2019) New approach for understanding genome variations in much data as the FWW MKT) and allows the independent KEGG. Nucleic Acids Res., 47, D590–D595. estimation of both adaptive and weakly deleterious sub- 28. Salvador-Martı́nez,I., Coronado-Zamora,M., Castellano,D., stitutions. PN , the count of segregating sites in the non- Barbadilla,A. and Salazar-Ciudad,I. (2018) Mapping selection within drosophila melanogaster Embryo’s anatomy. Mol. Biol. Evol., 35, synonymous class, is discomposed into the number of neu- 66–79. tral variants and the number of weakly deleterious variants, 29. Hahn,M.W. (2018) Molecular Population Genetics. Sinauer PN = PN neutral + PN weakly del. (9). The estimation of both Associates, Sunderland. numbers allows estimating positive (adaptive) and negative selection independently. α is estimated from the standard MKT table discounting weakly deleterious variants: PN is substituted by the expected number of neutral segregating APPENDIX 1 sites, PN neutral . The correct estimate of α is then α = 1 – (PN neutral /PS × DN /DS ). McDonald and Kreitman test (MKT) The standard McDonald and Kreitman test (MKT) (2) was developed to be applied to protein-coding sequences, com- Asymptotic MKT (aMKT) bining both divergence (D) and polymorphism (P) sites, and Messer and Petrov (10) proposed an asymptotic extension categorizing mutations as synonymous (PS , DS ) and non- of MKT that takes slightly deleterious mutations into ac- synonymous (PN , DN ). If all mutations are either strongly count and yields accurate estimates of α (Figure 1E). This deleterious or neutral, then DN /DS is expected to roughly method, named asymptotic MKT, is robust to the presence equal PN /PS . In contrast, if positive selection is operating of selective sweeps (hitchhiking) and to the segregation of in the region, adaptive mutations rapidly reach fixation and slightly deleterious substitutions (BGS). In this method, α thus contribute relatively more to divergence than to poly- is estimated in different frequency intervals (x) and these morphism when compared to neutral mutations, and then values are then adjusted to an exponential function, of the DN /DS > PN /PS (Figure 1A). Assuming that adaptive mu- form: α fit(x) = a + b−cx . The asymptotic α estimate is ob- tations contribute little to polymorphism but substantially tained by extrapolating the value of this function to x = 1: to divergence, the proportion of non-synonymous substitu- α asymptotic = α fit (x = 1) . tions that have been fixed by positive selection can be in- The asymptotic MKT has been extended to estimate ferred as α = 1 – (PN /PS × DN /DS ) (7) (Figure 1B). The both positive (adaptive) and negative selection (Coronado- main limitation of the test is the presence in the popula- Zamora, submitted). aMKT requires a high volume of poly- tion of non-synonymous slightly deleterious variants, bias- morphic data to fit the asymptotic function, being a suitable ing downward the estimates of adaptive evolution (α). Be- method in the case of concatenating numerous variants of low are three proposed methods to correct the bias. multiple genes.