Citation: Garud NR, Messer PW, Buzbas EO, Petrov DA (2015) Recent Selective Sweeps in North American Drosophila melanogaster Show Signatures of Soft Sweeps. PLoS Genet 11(2): e1005004. https://doi.org/10.1371/journal.pgen.1005004
Editor: Gregory P. Copenhaver, The University of North Carolina at Chapel Hill, UNITED STATES
Received: September 15, 2014; Accepted: January 14, 2015; Published: February 23, 2015
Copyright: © 2015 Garud et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited
Data Availability: All relevant data are within the paper and its Supporting Information files.
Funding: This work was supported by the National Institute of Health (www.nih.gov) grants R01 GM100366, R01 GM097415, R01 GM089926 to DAP, and R01 GM081441 to EOB, the National Science Foundation Graduate Research Fellowship (www.nsfgrfp.org) to NRG, and the Human Frontiers Science Program fellowship (www.hfsp.org) to PWM. The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript.
Competing interests: The authors have declared that no competing interests exist.
The ability to identify genomic loci subject to recent positive selection is essential for our efforts to uncover the genetic basis of phenotypic evolution and to understand the overall role of adaptation in molecular evolution. The fruit fly Drosophila melanogaster is one of the classic model organisms for studying the molecular bases and signatures of adaptation. Recent studies have provided evidence for pervasive molecular adaptation in this species, suggesting that approximately 50% of the amino acid changing substitutions, and similarly large proportions of non-coding substitutions, were adaptive [1,2,3,4,5,6,7,8,9]. There is also evidence that at least some of these adaptive events were driven by strong positive selection (~1% or larger), depleting levels of genetic variation on scales of tens of thousands of base pairs in length [10,11].
If adaptation in D. melanogaster is indeed common and often driven by strong selection, it should be possible to detect genomic signatures of recent and strong adaptation [12,13,14]. Three cases of recent and strong adaptation in D. melanogaster are well documented and can inform our intuitions about the expected genomic signatures of such adaptive events. First, resistance to the most commonly used pesticides, carbamates and organophosphates, is known to be largely due to three point mutations at highly conserved sites in the gene Ace, which encodes the neuronal enzyme Acetylcholinesterase [15,16,17]. Second, resistance to DDT evolved via a series of adaptive events that included insertion of an Accord transposon in the 5’ regulatory region of the gene Cyp6g1, duplication of the locus, and additional transposable element insertions into the locus [18,19]. Finally, increased resistance to infection by the sigma virus, as well as resistance to certain organophosphates, has been associated with a transposable element insertion in the protein-coding region of the gene CHKov1 [20,21].
In-depth population genetic studies [17,19,21] of adaptation at these loci revealed that in all three cases adaptation failed to produce classic hard selective sweeps, but instead generated patterns compatible with soft sweeps. In a hard selective sweep, a single adaptive haplotype rises in frequency and removes genetic diversity in the vicinity of the adaptive locus [22,23,24]. In contrast, in a soft sweep multiple adaptive alleles present in the population as standing genetic variation (SGV) or entering as multiple de novo adaptive mutations increase in frequency virtually simultaneously bringing multiple haplotypes to high frequency [25,26,27,28,29]. In the cases of Ace and Cyp6g1, soft sweeps involved multiple de novo mutations [17,19,21] that arose after the introduction of pesticides, whereas in the case of CHKov1, a soft sweep arose in out-of-African populations from standing genetic variation (SGV) [17,19,21] present at low frequencies in the ancestral African population [20,21].
Unfortunately, most scans for selective sweeps in population genomic data have been designed to detect hard selective sweeps (although see ) and focus on such signatures as a dip in neutral diversity around the selected site [22,24,31], an excess of low or high-frequency alleles in the frequency spectrum of polymorphisms surrounding the selected site (i.e. Tajima’s D, Fay and Wu’s H, and Sweepfinder) [32,33,34,35,36], the presence of a single common haplotype , or the observation of a long and unusually frequent haplotype (iHS) [36,38,39,40]. In a soft sweep, however, multiple haplotypes linked to the selected locus can rise to high frequency and levels of diversity and allele frequency spectra should therefore be perturbed to a lesser extent than in a hard sweep. As a result, methods based on the levels and frequency distributions of neutral diversity have low power to detect soft sweeps [13,28,41,42].
Some genomic signatures do have power to detect both hard and soft sweeps. In particular, linkage disequilibrium (LD) measured between pairs of sites or as haplotype homozygosity should be elevated in both hard and soft sweeps. This expectation holds for hard sweeps and for soft sweeps that are not too soft, that is soft sweeps that have such a large number of independent haplotypes bearing adaptive alleles that linkage disequilibrium is no longer elevated beyond neutral expectations [41,43].
Given that none of the described cases of adaptation at Ace, Cyp6g1, and CHKov1 produced hard sweeps, it is possible that additional cases of recent selective sweeps in D. melanogaster remain to be discovered. Here we develop a statistical test based on modified haplotype homozygosity for detecting both hard and soft selective sweeps in population genomic data. We apply this test in a genome-wide scan in a North American population of D. melanogaster using the Drosophila Genetic Reference Panel (DGRP) data set , consisting of 162 fully sequenced isogenic strains from a North Carolina population. Our scan recovers the three known soft sweeps at Ace, Cyp6g1, and CHKov1, and identifies a large number of additional recent and strong selective sweeps. We develop an additional haplotype homozygosity statistic that can distinguish hard from soft sweeps and argue that the haplotype frequency spectra at the top 50 candidate sweeps are best explained by soft selective sweeps.
Slow decay of linkage disequilibrium in the DGRP data
In this paper, we develop a set of new statistics for the detection and characterization of positive selection based on measurements of haplotype homozygosity in a predefined window. Our reasoning in developing these statistics is that haplotype homozygosity, defined as a sum of squares of the frequencies of identical haplotypes in a window, should be a sensitive statistic for the detection of both hard and soft sweeps, as long as the window is large enough that neutral demographic processes are unlikely to elevate haplotype homozygosity by chance [41,43]. At the same time, the window must not be so large that even strong sweeps can no longer generate frequent haplotypes spanning the whole window.
In order to determine an appropriate window length for the measurement of haplotype homozygosity in the DGRP data set, we first assessed the length scale of linkage disequilibrium decay expected in the DGRP data under a range of neutral demographic models for North American D. melanogaster. This length scale should roughly correspond to the window size over which we are unlikely to observe substantial haplotype structure by chance. We considered six demographic models (Fig. 1). The first demographic model is an admixture model of the North American D. melanogaster population proposed by Duchen et al. . In this model, the North American population was co-founded by flies from Africa and Europe 3.05×10–4Ne generations ago (where Ne ≈ 5x106). The second model is a modified admixture model, also proposed by Duchen et al. , in which the founding European population underwent a bottleneck before the admixture event (see S1 Table for complete parameterizations of both admixture models). The third model has a constant effective population size of Ne = 106 , which we considered for its simplicity, computational feasibility and, as we will argue below, its conservativeness for the purposes of detecting selective sweeps using our approach in the DGRP data. The fourth model is a constant Ne = 2.7x106 demographic model fit to Watterson’s θW estimated from short intron autosomal polymorphism data from the DGRP dataset (Methods). Finally, we fit a family of out-of-Africa bottleneck models to short intron regions in the DGRP data set using DaDi  (S2 Table) (Methods). The two bottleneck models we ultimately used are a severe but short bottleneck model (NB = 0.002, TB = 0.0002) and a shallow but long bottleneck model (NB = 0.4, TB = 0.0560), both of which fit the data equally well among a range of other inferred bottleneck models (see S1 Fig. for parameterization). All models except for the constant Ne = 106 model fit the DGRP short intron data in terms of the number of segregating sites (S) and pairwise nucleotide diversity (π) (S3 Table).
Fig 1. Neutral demographic models.
We considered six neutral demographic models for the North American D. melanogaster population: (A) An admixture model as proposed by Duchen et al. . (B) An admixture model with the European population undergoing a bottleneck. This model was also tested by Duchen et al. , but the authors found it to have a poor fit. See S1 Table for parameter estimates and symbol explanations for models A and B. (C) A constant Ne = 106 model. (D) A constant Ne = 2.7x106 model fit to Watterson’s θW measured in short intron autosomal polymorphism data from the DGRP data set. (E) A severe short bottleneck model and (F) a shallow long bottleneck model fit to short intron regions in the DGRP data set using DaDi . See S2 Table for parameter estimates for models E and F. All models except for the constant Ne = 106 model fit the DGRP short intron data in terms of S and π (S3 Table).
We compared the decay in pair-wise LD in the DGRP data at distances from a few base pairs to 10 kb with the expectations under each of the six demographic models using parameters relevant for our subsequent analysis of the DGRP data (Fig. 2). Specifically, we matched the sample depth of the DGRP data set (145 strains after quality control) and assumed a mutation rate (μ) of 10–9 events/bp per generation  and a recombination rate (ρ) of 5×10–7 centimorgans/bp (cM/bp) . In the DGRP data analysis below, we exclude regions with a low recombination rate (ρ < 5x10–7 cM/bp). The use of ρ = 5x10–7 cM/bp should therefore generate higher LD in simulations than in the DGRP data and thus should be conservative for the purposes of defining the expected length scale of LD decay.
Fig 2. Elevated long-range LD in DGRP.
LD in DGRP data is elevated as compared to any neutral demographic model, especially for long distances. Pairwise LD was calculated in DGRP data for regions of the D. melanogaster genome with ρ ≥ 5×10–7 cM/bp. Neutral demographic simulations were generated with ρ = 5×10–7 cM/bp. Pairwise LD was averaged over 3×104 simulations in each neutral demographic scenario.
Fig. 2 shows that LD in the DGRP data is elevated beyond neutral expectations at all length scales (consistent with the observations in ), and dramatically so at the 10 kb length scale. The elevation in LD observed in the data is indicative of either linked positive selection driving haplotypes to high frequency, a lack of fit of current demographic models to the data, or both. Simulations under the most realistic demographic model, admixture , have the fastest decay in LD (S2 Fig.). This is likely because admixture models with two bottlenecks that are fit to diversity statistics generate more haplotypes compared to single bottleneck models, since the same haplotype is unlikely to be sampled independently in both bottlenecked ancestral populations. In contrast, LD under the constant Ne = 106 demographic scenario decays slower than in any other demographic scenario, as expected given that this model has the smallest effective population size.
Fig. 2 suggests that windows of 10 kb are large enough that neutral demography is unlikely to generate high values of LD and elevate haplotype homozygosity by chance, and should thus prevent a high rate of false positives. At the same time, the use of 10 kb windows for the measurement of haplotype homozygosity should still allow us to detect many reasonably strong sweeps, including the known cases of recent adaptation. The footprint of a hard selective sweep extends over approximately s/[log(Nes)ρ] basepairs, where s is the selection strength, Ne the population size, and ρ the recombination rate [22,23,51]. Sweeps with a selection coefficient of s = 0.05% or greater are thus likely to generate sweeps that span 10 kb windows in areas with recombination rate of 5×10–7 cM/bp. As the recombination rate increases, only selective sweeps with s > 0.05% should be observed in the 10 kb windows. Genomic analyses have suggested that adaptation in Drosophila is likely associated with a range of selection strengths, including values of ~1% [7,8,10] or greater as observed at Ace, Cyp6g1, and CHKov1. Our use of 10 kb windows in the rest of the analysis should thus bias the analysis toward detecting the cases of strongest adaptation in Drosophila.
Haplotype spectra expectations under selective sweeps of varying softness
We investigated haplotype spectra in simulations of neutral demography and both hard and soft selective sweeps arising from de novo mutations as well as SGV. For all haplotype spectra and homozygosity analyses in this paper we use windows of 400 SNPs, corresponding roughly to 10 kb in the DGRP data (Fig. 2). Haplotypes within a 400 SNP window are grouped together if they are identical at all SNPs in the window. We fixed the number of SNPs in a window to eliminate variability in the haplotype spectra due to varying numbers of SNPs.
The lower SNP density of the constant Ne = 106 model (S3 Table) effectively increases the size of the analysis window in terms of the number of base pairs when defining the windows in terms of the number of SNPs. Thus, the constant Ne = 106 model should reduce the rate of false positives because the recombination rate under this model is artificially increased. We therefore use the constant Ne = 106 model for the subsequent simulations of neutrality and selective sweeps.
To visualize sample haplotype frequency spectra, we simulated incomplete and complete sweeps with frequencies of the adaptive mutation (PF) at 0.5 or 1 at the time when selection ceased. (Note that below we will investigate a large number of scenarios, focusing on the effects of varying selection strength and the decay of sweep signatures with time). The number of independent haplotypes that rise in frequency simultaneously in soft sweeps—we call this “softness” of a sweep—should increase either (i) when the rate of mutation to de novo adaptive alleles at a locus becomes higher and multiple alleles arise and establish after the onset of selection at a higher rate, or (ii) when adaptation uses SGV with previously neutral or deleterious alleles that are present at higher frequency at the onset of selection [27,29]. More specifically, for sweeps arising from multiple de novo mutations, Pennings and Hermisson  showed that the key population genetic parameter that determines the softness of the sweep is θA = 4NeμA, proportional to the product of Ne, the variance effective population size estimated over the period relevant for adaptation [14,52], and μA, the mutation rate toward adaptive alleles at a locus per individual per generation . The mutation-limited regime with hard sweeps corresponds to θA << 1, whereas θA > 1 specifies the non-mutation-limited regime with primarily soft sweeps. As θA becomes larger, the sweeps become softer as more haplotypes increase in frequency simultaneously . In the case of sweeps arising from SGV, the softness of a sweep is governed by the starting partial frequency of the adaptive allele in the population prior to the onset of selection. For any given rate of recombination, adaptive alleles starting at a higher frequency at the onset of selection should be older and should thus be present on more distinct haplotypes and give rise to softer sweeps .
As can be seen in Fig. 3, most haplotypes in neutral demographic scenarios are unique in our 400 SNP windows, whereas selective sweeps can generate multiple haplotypes at substantial frequencies. Our plot of the haplotype frequency spectra and the expected numbers of adaptive haplotypes show that sweeps arising from de novo mutations become soft with multiple frequent haplotypes in the sample when θA ≥ 1. Sweeps from SGV become soft when the starting partial frequency of the adaptive allele prior to the onset of selection is ≥ 10–4 (100 alleles in the population). In both cases, sweeps become monotonically softer as θA increases or, respectively, the starting partial frequency of the adaptive allele becomes higher. These results conform to the expectations derived in .
Fig 3. Number of adaptive haplotypes in sweeps of varying softness.
The number of origins of adaptive mutations on unique haplotype backgrounds was measured in simulated sweeps of varying softness arising from (A) de novo mutations with θA values ranging from 10–2 to 102 and (D) SGV with starting frequencies ranging from 10–6 to 10–1. Sweeps were simulated under a constant Ne = 106 demographic model with a recombination rate of 5×10–7 cM/bp, selection strength of s = 0.01, partial frequency of the adaptive allele after selection has ceased of PF = 1 and 0.5, and in sample sizes of 145 individuals. 1000 simulations were averaged for each data point. Additionally we show sample haplotype frequency spectra for (B) incomplete and (C) complete sweeps arising from de novo mutations as well as (E) incomplete and (F) complete sweeps arising from SGV. In (G) we show haplotype frequency spectra for a random simulation under the six neutral models considered in this paper. The height of the first bar (light blue) in each frequency spectrum indicates the frequency of the most prevalent haplotype in the sample of 145 individuals, and heights of subsequent colored bars indicate the frequency of the second, third, and so on most frequent haplotypes in a sample. Grey bars indicate singletons. Sweeps generated with a low θA or low starting partial frequency of the adaptive allele prior to the onset of selection have one frequent haplotype in the sample and look hard. In contrast, sweeps look increasingly soft as the θA or starting partial frequency of the adaptive allele prior to the onset of selection increase and there are multiple frequent haplotypes in the sample.
Definitions of haplotype homozygosity statistics H1, H12, and H123
The increase of haplotype population frequencies in both hard and soft sweeps can be captured using haplotype homozygosity [30,39,41]. If pi is the frequency of the ith most common haplotype in a sample, and n is the number of observed haplotypes, then haplotype homozygosity is defined as H1 = Σi = 1, …npi2. We can expect H1 to be particularly high for hard sweeps, with only one adaptive haplotype at high frequency in the sample (Fig. 4A). Thus, H1 is an intuitive candidate for a test of neutrality versus hard sweeps, where the test rejects neutrality for high values of H1. A test based on H1 may also have acceptable power to detect soft sweeps in which only a few haplotypes in the population are present at high frequency. However, as sweeps become softer and the number of sweeping haplotypes increases, the relative contribution of individual haplotypes towards the overall H1 value decreases, and the power of a test based on H1 is expected to decrease.
Fig 4. Haplotype homozygosity statistics.
Depicted are squares of haplotype frequencies for hard (red) and soft (blue) sweeps. Each edge of the square represents haplotype frequencies ranging from 0 to 1. The top row shows incomplete hard sweeps with one prevalent haplotype present in the population at frequency p1, and all other haplotypes present as singletons. The bottom row shows incomplete soft sweeps with one primary haplotype with frequency p1 and a second, less abundant haplotype at frequency p2, with the remaining haplotypes present as singletons. H1 is the sum of the squares of frequencies of each haplotype in a sample and corresponds to the total colored area. Hard sweeps are expected to have a higher H1 value than soft sweeps. In H12, the first and second most abundant haplotype frequencies in a sample are combined into a single combined haplotype frequency and then homozygosity is recalculated using this revised haplotype frequency distribution. By combining the first and second most abundant haplotypes into a single group, H12 should have more similar power to detect hard and soft sweeps than H1. H2 is the haplotype homozygosity calculated after excluding the most abundant haplotype. H2 is expected to be larger for soft sweeps than for hard sweeps. We ultimately use the ratio H2/H1 to differentiate between hard and soft sweeps as we expect this ratio to have even greater discriminatory power than H2 alone.
To have a better ability to detect hard and soft sweeps using homozygosity statistics, we developed a modified homozygosity statistic, H12 = (p1 + p2)2 + Σi>2pi2 = H1 + 2p1p2, in which the frequencies of the first and the second most common haplotype are combined into a single frequency (Fig. 4B). A statistical test based on H12 is expected to be more powerful in detecting soft sweeps than H1 because it combines frequencies of two similarly abundant haplotypes into a single frequency, whereas for hard sweeps the combination of the frequencies of the first and second most abundant haplotypes should not change haplotype homozygosity substantially . We also considered a third test statistic, H123, which combines frequencies of the three most prevalent haplotypes in a sample into a single haplotype and then computes homozygosity. We will primarily employ H12 in subsequent analyses but will consider the effects of using H1 and H123 briefly as well.
Ability of H12 to detect selective sweeps of varying softness
To assess the ability of H12 to detect sweeps of varying softness and to distinguish positive selection from neutrality, we measured H12 in simulated sweeps arising from both de novo mutations and SGV while varying s, PF, and the time since the end of the sweep, TE, measured in units of 4Ne generations in order to model the decay of a sweep through recombination and mutation events over time. We first investigate the behavior of H12 under different selective regimes and then investigate its power in comparison with the popular haplotype statistic iHS.
Fig. 5A shows that for complete and incomplete sweeps with s = 0.01 and TE = 0, H12 monotonically decreases as a function of θA over the interval from 10–2 to 102. When θA ≤ 0.5, many sweeps are hard and H12 values are high. When θA ≈ 1, and practically all sweeps are soft, but not yet extremely soft, H12 retains much of its power. However, for θA > 10, where sweeps are extremely soft, H12 decreases substantially. Similarly, H12 is maximized when the starting frequency of the allele is 10–6 (one copy of the allele in the population generating hard sweeps from SGV) and becomes very small as the frequency of the adaptive allele increases beyond >10-3 (>1000 copies of the allele in the population) (Fig. 5B). Therefore, H12 has reasonable power to detect soft sweeps in samples of hundreds of haplotypes, as long as they are not extremely soft, but remains somewhat biased in favor of detecting hard sweeps.
Fig 5. H12 values in sweeps of varying softness.
H12 values were measured in simulated sweeps arising from (A) de novo mutations with θA values ranging from 10–2 to 102 and (B) SGV with starting frequencies ranging from 10–6 to 10–1. Sweeps were simulated under a constant Ne = 106 demographic model with a recombination rate of 5×10–7 cM/bp, selection strength of s = 0.01, ending partial frequencies of the adaptive allele after selection has ceased, PF = 1 and 0.5, and in samples of 145 individuals. Each data point was averaged over 1000 simulations. H12 values rapidly decline as the softness of a sweep increases and as the ending partial frequency of the adaptive allele decreases. In (C) and (D), s was varied while keeping PF constant at 0.5 for sweeps from de novo mutations and SGV, respectively. H12 values increase as s increases, though for very weak s we observe a ‘hardening’ of sweeps where fewer adaptive alleles reach establishment frequency. In (E) and (F), the time since selection ended (TE) was varied for incomplete (PF = 0.5) and complete (PF = 1) sweeps respectively while keeping s constant at 0.01. As the age of a sweep increases, sweep signatures decay and H12 loses power.
H12 also increases as the ending partial frequency of the adaptive allele after selection ceased (PF) increases from 0.5 to 1 (Fig. 5A and 5B) and as the selection strength increases from 0.001 to 0.1 (Fig. 5C and 5D). We observe that sweeps arising from SGV with low selection coefficients have lower H12 values (Fig. 5D). This is most likely because such weak sweeps are effectively harder: as more of the haplotypes fail to establish, fewer haplotypes end up sweeping in the population leading to higher values of haplotype homozygosity. Fig. 5E and 5F further show that incomplete and complete sweeps decay with time due to recombination and mutation events, resulting in monotonically decreasing values of H12 with time. Overall this analysis demonstrates that H12 has most power to detect recent sweeps driven by strong selection.
We also assessed the ability of H12 to detect selective sweeps as compared to H1 and H123 by calculating the values of H1, H12, and H123 for sweeps generated under the parameters s = 0.01, TE = 0 and PF = 0.5. H12 consistently, albeit modestly, increases the homozygosity for younger soft sweeps as compared to H1 (S3 Fig.). The increase in homozygosity using H123 is marginal relative to homozygosity levels achieved by H12, so we chose not to use this statistic in our study.
Finally, we compared the abilities of H12 and iHS (integrated haplotype score), a haplotype-based statistic designed to detect incomplete hard sweeps [39,40], to detect both hard and soft sweeps. We created receiving operator characteristic (ROC) curves , which plot the true positive rate (TPR) of correctly rejecting neutrality in favor of a sweep (hard or soft) given that a sweep has occurred versus the false positive rate (FPR) of inferring a selective sweep, when in fact a sweep has not occurred.
In our simulations of selective sweeps we used θA = 0.01 as a proxy for scenarios generating almost exclusively hard sweeps, and θA = 10 as a proxy for scenarios generating almost exclusively soft sweeps. We chose θA = 10 for soft sweeps because this is the highest θA value with which H12 can still detect sweeps before substantially losing power given our window size of 400 SNPs and sample size of 145. Note that for soft sweeps with a lower value of θA the power of H12 should be higher. We modeled incomplete sweeps with PF = 0.1, 0.5, and 0.9, with varying times since selection had ceased of TE = 0, 0.001, and 0.01 in units of 4Ne generations. We simulated sweeps under three selection coefficients, s = 0.001, 0.01, and 0.1.
Fig. 6 and S4 Fig. show that the tests based on H12 and iHS have similar power for the detection of hard sweeps, although in the case of old and strong hard sweeps (TE = 0.01, s ≥ 0.01) iHS performs slightly better than H12. On the other hand, H12 substantially outperforms iHS in detecting soft sweeps and has high power when selection is sufficiently strong and the sweeps are sufficiently young. As sweeps become very old, neither statistic can detect them well, as expected.
Fig 6. Power analysis of H12 and iHS under different sweep scenarios.
The plots show ROC curves for H12 and iHS under various sweep scenarios with the specified selection coefficients (s), and the time of the end of selection (TE) in units of 4Ne generations. In all scenarios, the ending partial frequency of the adaptive allele was 0.5. False positive rates (FPR) were calculated by counting the number of neutral simulations that were misclassified as sweeps under a specific cutoff. True positive rates (TPR) were calculated by counting the number of simulations correctly identified as sweeps under the same cutoff. Hard and soft sweeps were simulated from de novo mutations with θA = 0.01 and 10, respectively, under a constant effective population size of Ne = 106, a neutral mutation rate of 10–9 bp/gen, and a recombination rate of 5×10–7 cM/bp. A total of 5000 simulations were conducted for each evolutionary scenario. H12 performs well in identifying recent and strong selective sweeps, and is more powerful than iHS in identifying soft sweeps.
H12 scan of DGRP data
We applied the H12 statistic to DGRP data in sliding windows of 400 SNPs with the centers of each window iterated by 50 SNPs. To classify haplotypes within each analysis window, we assigned the 400 SNP haplotypes into groups according to exact sequence identity. If a haplotype with missing data matched multiple haplotypes at all genotyped sites in the analysis window, then the haplotype was randomly assigned to one of these groups (Methods).
To assess whether the observed H12 values in the DGRP data along the four autosomal arms are unusually high as compared to neutral expectations, we estimated the expected distribution of H12 values under each of the six neutral demographic models. Fig. 7 shows that genome-wide H12 values in DGRP data are substantially elevated as compared to expectations under any of the six neutral demographic models. In addition, there is a long tail of outlier H12 values in the DGRP data suggestive of recent strong selective sweeps.
Fig 7. Elevated H12 values and long-range LD in DGRP data.
(A) Genome-wide H12 values in DGRP data are elevated as compared to expectations under any neutral demographic model tested. Plotted are H12 values for DGRP data reported in analysis windows with ρ ≥ 5×10-7 cM/bp. Red dots overlaid on the distribution of H12 values for DGRP data correspond to the highest H12 values in outlier peaks of the DGRP scan at the 50 top peaks depicted in Fig. 8A. Note that most of the points in the tail of the H12 values calculated in DGRP data are part of the top 50 peaks as well. Neutral demographic simulations were generated with ρ = 5×10–7 cM/bp. Plotted are the result of approximately 1.3x105 simulations under each neutral demographic model, representing ten times the number of analysis windows in DGRP data.
To identify regions of the genome with H12 values significantly higher than expected under neutrality, we calculated critical values (H12o) under each of the six neutral models based on a 1-per-genome false discovery rate (FDR) criterion. Our test rejects neutrality in favor of a selective sweep when H12 > H12o (Methods and S1 Text). The critical H12o values under all neutral demographic models are similar to the median H12 value observed in the DGRP data (Table 1), consistent with the observations of elevated genome-wide haplotype homozygosity and much slower decay in LD at the scale of 10 kb in the DGRP data compared to all neutral expectations (Fig. 2). We focused on the constant Ne = 106 model because it yields a relatively conservative H12o value (Table 1) and preserves the most long-range, pair-wise LD in simulations (Fig. 2).
For our genomic scan we chose to use the 1-per-genome FDR value calculated under the constant Ne = 106 model with a recombination rate of 5×10–7 cM/bp. Note that most H12o values are similar to the genome-wide median H12 value of 0.0155.
In order to call individual sweeps, we first identified all windows with H12 > H12o in the DGRP data set under the constant Ne = 106 model. We then grouped together consecutive windows as belonging to the same ‘peak’ if the H12 values in all of the grouped windows were above H12o for a given model and recombination rate (Methods). We then chose the window with the highest H12 value among all windows in a peak and used this H12 value to represent the entire peak.
We focused on the top 50 peaks with empirically most extreme H12 values, hypothesized to correspond to the strongest and/or most recent selective events (Fig. 8A). The windows with the highest H12 values for each of the top 50 peaks are highlighted in Fig. 8A. The highest H12 values for the top 50 peaks are in the tail of the distribution of H12 values in the DGRP data (Fig. 7) and thus are outliers both compared to the neutral expectations under all six demographic models and the empirical genomic distribution of H12 values. We observed peaks that have H12 values higher than H12o on all chromosomes, but found that there are significantly fewer peaks on 3L (2 peaks) than the approximately 13 out of 50 top peaks expected when assuming a uniform distribution of the top 50 peaks genome-wide (p = 0.00016, two-sided binomial test, Bonferroni corrected).
Fig 8. H12 and iHS scan in DGRP data along the four autosomal arms.
(A) H12 scan. Each data point represents the H12 value calculated over an analysis window of size 400 SNPs centered at the particular genomic position. Grey points indicate regions in the genome with recombination rates lower than 5×10–7 cM/bp we excluded from our analysis. The orange line represents the 1-per-genome FDR line calculated under a neutral demographic model with a constant population size of 106 and a recombination rate of 5×10–7 cM/bp. Red and blue points highlight the top 50 H12 peaks in the DGRP data relative to the 1-per-genome FDR line. Red points indicate the peaks that overlap the top 10% of 100Kb windows with an enrichment of SNPs with |iHS| > 2 in B. We identify three well-characterized cases of selection in D. melanogaster at Ace, CHKov1, and Cyp6g1 as the three highest peaks. (B) iHS scan. Plotted are the number of SNPs in 100Kb windows with |iHS| > 2. Highlighted in red and blue are the top 10%100Kb windows (a total of 95 windows). Red points correspond to those windows that overlap the top 50 peaks in the H12 scan. The positive controls, Ace, CHKov1, and Cyp6g1 are all among the top 10% windows.
The three peaks with the highest observed H12 values correspond to the three known cases of positive selection in D. melanogaster at the genes Ace, Cyp6g1, and CHKov1 [17,19,21], confirming that the H12 scan is capable of identifying previously known cases of adaptation. In S4 Table, we list all genes that overlap with any of the top 50 peaks. Fig. 9A and S5 Fig. show the haplotype frequency spectra observed at the top 50 peaks. In contrast, Fig. 9B shows the frequency spectra observed under the six demographic models with the corresponding critical H12o values.
Fig 9. Haplotype frequency spectra for the top 10 peaks and extreme outliers under neutral demographic scenarios.
(A) Haplotype frequency spectra for the top 10 peaks in the DGRP scan with H12 values ranging from highest to lowest. For each peak, the frequency spectrum corresponding to the analysis window with the highest H12 value is plotted, which should be the “hardest” part of any given peak. At all peaks there are multiple haplotypes present at high frequency, compatible with signatures of soft sweeps shown in Fig. 5. None of the cases have a single haplotype present at high frequency, as would be expected for a hard sweep. (B) In contrast, the haplotype frequency spectra corresponding to the extreme outliers under the six neutral demographic scenarios have critical H120 values that are significantly lower than the H12 values at the top 10 peaks.
We performed several tests to ensure the robustness of the H12 peaks to potential artifacts (S1 Text). We first tested for associations of H12 peaks with inversions in the sample, but did not find any (S1 Text, S5 Table). In addition, we reran the scan in three different data sets of the same population and confirmed that unaccounted population substructure and variability in sequencing quality do not confound our results (S1 Text, S7 Fig.). We also sub-sampled the DGRP data set to 40 strains ten times and plotted the resulting distributions of H12 values. We found that in all subsamples there is an elevation in haplotype homozygosity relative to neutral demographic scenarios, suggesting that the elevation in haplotype homozygosity values is driven by the whole sample and not a particular subset of individuals (S8 Fig.). Finally, to ensure that haplotype homozygosity is not elevated by family structure, we excluded all related individuals and reran the scan, again recovering the majority of our top peaks (S1 Text, S7 Fig.).
We scanned chromosome 3R using H1 and H123 as our test statistics in order to determine the impact of our choice of grouping the two most frequent haplotypes together in our H12 test statistic on the location of the identified peaks (S9 Fig.). We found that the locations of the identified peaks are similar with all three statistics, but that some smaller peaks that cannot be easily identified with H1 are clearly identified with H12 and H123, as expected.
iHS scan of DGRP data
We applied the iHS statistic as described in Voight et al. 2006  to all SNPs in the DGRP data to determine the concordance in the sweep candidates identified by iHS and H12 (Methods). Briefly, we searched for 100 kb windows that have an unusually large number of SNPs with standardized iHS values (|iHS|) > 2. The positive controls Ace, Cyp6g1, and CHKov1 are located within the 95 top 10% iHS 100 kb windows (Fig. 8B), validating this approach.
To determine how often a candidate region identified in the H12 scan is identified in the iHS scan and vice versa, we overlapped the top 50 H12 peaks with the 95 top 10% iHS 100Kb windows. We defined an overlap as the non-empty intersection of the two genomic regions defining the boundaries of a peak in the H12 scan and the non-overlapping 100Kb windows used to calculate enrichment of |iHS| values. We found that 18 H12 peaks overlap 28 |iHS| 100Kb enrichment windows. In contrast, fewer than 5 H12 peaks are expected to overlap approximately 7 iHS 100Kb windows by chance (Methods). The concordance between the two scans confirms that many of the peaks identified in the two scans are likely true selective sweeps and also suggests that the two approaches are not entirely redundant.
Distinguishing hard and soft sweeps based on the statistic H2/H1
Our analysis of H12 haplotype homozygosity and the decay in long range LD in DGRP data suggests that extreme outliers in the H12 DGRP scan are in locations of the genome that may have experienced recent and strong selective sweeps. The visual inspection of the haplotype spectra of the top 10 peaks in Fig. 9A and the remaining 40 peaks in S5 Fig. reveals that they contain many haplotypes at substantial frequency. These spectra do not appear similar to those generated by hard sweeps in Fig. 3 or extreme outliers under neutrality in Fig. 9B, but instead visually resemble incomplete soft sweeps with s = 0.01 and PF = 0.5 either from de novo mutations with θA between 1 and 20 or from SGV starting at partial frequencies of 5x10–5 to 5x10–4 prior to the onset of selection (Fig. 3). The sweeps also appear to become softer as H12 decreases, consistent with our expectation that H12 should lose power for softer sweeps.
In order to gain intuition about whether the haplotype spectra for the top 50 peaks can be more easily generated either by hard or soft sweeps under various evolutionary scenarios, we developed a new haplotype homozygosity statistic, H2/H1, where H2 = Σi>1pi2 = H1—p12 is haplotype homozygosity calculated using all but the most frequent haplotype (Fig. 4C). We expect H2 to be lower for hard sweeps than for soft sweeps because in a hard sweep only one adaptive haplotype is expected to be at very high frequency . The exclusion of the most common haplotype should therefore reduce haplotype homozygosity precipitously. As sweeps get softer, however, multiple haplotypes start appearing at high frequency in the population and the exclusion of the most frequent haplotype should not decrease the haplotype homozygosity to the same extent. Conversely H1, the homozygosity calculated using all haplotypes, is expected to be higher for a hard sweep than for a soft sweep as we described above. The ratio H2/H1 between the two measures should thus increase monotonically as a sweep becomes softer, thereby offering a summary statistic that, in combination with H12, can be used to test whether the observed haplotype patterns are more likely to be generated by hard or soft sweeps. Note that we intend H2/H1 to be measured near the center of the sweep where H12 is the highest. Otherwise, when H2/H1 is estimated further away from the sweep center, mutation and recombination events will decay the haplotype signature and hard and soft sweep signatures can become indistinguishable.
Softness of sweeps at the top 50 H12 peaks
To assess the behavior of H2/H1 as a function of the softness of a sweep, we measured H2/H1 in simulated sweeps of varying softness arising from de novo mutations and SGV with various s, PF, and TE values. Fig. 10 shows that H2/H1 has low values for sweeps with θA ≤ 0.5 or when the starting partial frequency of the adaptive allele prior to the onset of selection is <10–5, i.e., when sweeps are mainly hard. As a sweep becomes softer, H2/H1 values approach one because no single haplotype dominates the haplotype spectrum. In the case of sweeps arising from de novo mutations, H2/H1 values are similar for partial (PF = 0.5) and complete sweeps (PF = 1) and for sweeps of varying strengths (s = 0.001, 0.01, 0.1). However, in the case of sweeps arising from SGV, sweeps with higher selection strengths do have higher H2/H1 values, reflecting the hardening of sweeps for smaller s values as we discussed previously (Fig. 5D). Both sweeps from de novo mutations and SGV have higher H2/H1 values for older sweeps, reflecting the decay of the haplotype frequency spectrum over time.
The beautiful screen on the new Android phones calls out to be touched. But maybe not by your fingers. We have all found, at one time or another, that touching the screen leaves fingerprints and grime on the screen. We have also found that our fingers are fat and clumsy when it comes to painting and writing.
So, we need something with a capacitive end or nib to touch the screen in a far more precise way to really take advantage of some of the newer apps. Enter the SGP Kuel H12 stylus. It's got the look and feel of a pen, minus all the messy ink. For this review, we're going to use apps like Sketchbook Mobile, Note Everything and Handwrite to judge the overall qualities of the stylus.
Read on for our full H12 stylus review!
Criteria for judging a stylus
Think back to the first time you held a fine pen in your hands – not a BIC or a PaperMate pen - but something really fine; think Montblanc, or Waterman. Remember how you picked it up – ever so carefully; how you felt the weight and balance of the instrument. Think about how the ink flowed from the nib onto the paper. A quality pen helps your writing; a quality stylus helps your app experience.
The criteria for judging the stylus will be:
- Appearance and finish
- Nib/tip and sense of flow on the screen
- Handwriting precision
- Drawing/painting capabilities
A fine stylus, like a fine pen, just feels right in the hand. The H12 is a fine stylus. The H12 employs a twist mechanism to open. There are several benefits to this design. First off, you don’t have to remove a cap -- which means you don't have to worry about losing one. Secondly, because you can twist the H12 closed, the nib -- the most fragile part of any stylus -- is protected when it sits in a pocket or purse.
When you hold the H12, it feels nicely balanced. The weight is evenly distributed above and below the middle of the stylus where the twist mechanism lies. The H12 is about the same thickness as a nice ballpoint pen.
The H12 would not be considered “heavy” or is it is light in weight. Rather, it is substantial and solid feeling but not too heavy as to cause fatigue if used for long periods of time. Its length is just a bit shorter than a typical ballpoint pen. The satin finish is easy to grip and there are no sharp edges to worry about with this stylus.
Appearance and finish
The H12 has the look of a premium product, which is no small feat at a very reasonable price of $19.99. My H12 is black with silver accents. While I didn’t notice any flaking or chipping in the finish of the Stylus, the lettering of the Kuel name was starting to wear off after a few weeks of use, but otherwise it maintained a nice satin finish with little signs of wear.
The H12 has a traditional pen clip towards the top, which is both functional and aesthetically pleasing. Sitting in a pocket, the H12 really does look like a nice pen.
Build quality is great for this stylus. The weight really gives it a sold feel and the twist mechanism requires just enough pressure to feel as though it is very well made and designed and will last for many “twistings.”
Nib/Tip and flow on screen
The H12 uses a special silicon tip with a high polymer abrasion resistant coating. The tip is both soft and spongy to the touch. The nib would be considered a moderate size at about 6mm in diameter.
Because the nib of the H12 (and that of most styli) is considerably larger than the tip of a “pen,” there is always the initial feeling that the stylus will somehow not be precise. With the H12, those concerns are allayed with first touch to the screen.
In each of the apps tested with the SGP H12, the feeling was one of precision and confidence when using the stylus for both writing and drawing. Very little pressure was needed to initiate the feeling of ink flowing from the stylus onto the iPad.
The “secret ingredient” of the H12 is that right below the cushioned silicon tip is a harder tip inside. When I pressed just a bit harder on the screen, I could feel this narrower tip make contact with the screen and give me even more precision when writing or drawing thin and straight lines.
There was no lag when writing on the screen with the H12. Perhaps the ultimate compliment to the stylus was when my wife looked at a note I had written in Handwrite and commented: “I can’t read that, it looks exactly like your regular, terrible handwriting!”
I tried a number of tests to examine handwriting precision. From creating numbered lists, to tracing objects, to writing longer sentences – the H12 did not disappoint.
I really had the illusion that I was writing on my screen when using the H12.
If, when drawing a line of text, I paused in the middle of my stroke and then continued, the area in which I paused looked to contain more “ink” – just like writing with a fountain pen. When starting a stroke from top to bottom, as I lifted up the H12 at the bottom of the writing stroke, it looks like more “ink” is at the very bottom – much in the fashion of a real pen. It makes sense, because as you end your stroke, you tend to put a bit more pressure on the pen, which increases the flow of ink.
The overall experience of writing with the H12 fell somewhere between writing with a roller ball and writing with a wider nibbed fountain pen – just in the sweet spot for my tastes.
There are two categories of drawing apps for Android devices; casual drawing apps/games and true artistic expression apps. The H12 excelled at casual drawing apps like Draw Something by OMGPOP.
This is really the perfect stylus for quickly drawing a Pictionary type drawing. I tried each of the “marker” sizes in the game and all the colors and using the H12 was both precise and quick – just what you need in a casual game.
In Sketchbook Mobile, the H12 was certainly acceptable for creating decent works of art. Each of the “tools” – pencil, sketch pen, brush, Outline, Draw and Write tools worked well with the H12. There was some hesitation in using the watercolor-like tool – in that more pressure was needed to really sketch and do things like shading when creating a picture. Most of the tools felt like their real life counterparts in terms of pressure, friction and flow.
Styluses are an emerging category and finding the right criteria on which to judge them can be tricky. Based on the stylus pens I have seen so far, there appear to be three distinct categories; those that excel at writing, those that excel at drawing and those that work well as a multi-use stylus for writing, drawing and basic navigation on the smartphone or tablet.
The Kuel H12 is very good for writing. For apps like Note Everything and Handwrite, this would be a very good stylus to use. Writing is precise and legible and the flow on the screen is quite good.
The H12 is also fine for drawing and painting apps like Sketchbook Mobile. It is certainly passable for the occasional artist, but true hard-core artists may want to look at more of a specialty stylus for artwork.
Where the H12 really excels is in the area that will matter most to the majority of stylus users – this is an excellent multi-purpose stylus. Using the H12 instead of my finger to navigate your Android phone or tablet, activate apps, type on the keyboard in addition to using it for Handwrite and Sketchbook Mobile made this a very versatile stylus.
- Great multipurpose stylus
- Sturdy construction
- Quality nib with nice pointed tip below
- Great for writing and general navigation
- Struggled a bit in drawing/painting apps
- No replaceable nib
The Kuel H12 is certainly among the very best stylus pens on the market today. For general, multi-purpose use you can’t go wrong with this well crafted, nice feeling stylus pen.
Do you use a stylus on your Android phone or tablet? Have a favorite? Let us know in this forum thread.
Other stylus pens worth looking at