Qualitative assessment of wing color pattern variation

Major wing pattern variation is largely explained by simple bi-allelic inheritance at a handful of loci of large effect in general agreement with previous studies[11],[20],[22],[24]. Overall, natural variation of the size and shape of the bands and bars elements in the fore- and hindwing were consistent with the concept that regions of black scales genetically define the position and size of these elements[25]–[27]. This idea is most easily demonstrated by the forewing band patterns of F1 individuals. In all our crosses, the F1 phenotypes are highly predictable and can be created by the superposition of parental black patterns (Figure 1). This interpretation is very different from the original interpretation of Sheppardet al.[11] who viewed color patterns as being expressed on a black background. Red color elements were also easily understood by allelic variation at a single large effect locus. The phenotypic effects of major loci in specific crosses (Figures 1 and2) are summarized below:

H. himera×H. e. etylus: Allellic variation at two codominant loci largely explains phenotypic variation in ourH. himera×H. e. etylus crosses[11],[20],[22],[24]; see alsoFigure 1A andFigure 2). TheD locus controlled all red and orange variation on both forewings and hindwings. Red and orange pattern elements were codominant and both F1 parents had intermediate phenotypes featuringH. erato derived orange patches at the base of the forewing, and orange hindwing rays blended with red hindwing bands derived fromH. himera. In the F2 offspring, the orange basal forewings and hindwings rays were always inherited together in our crosses and qualitative variation in forewing band size and shape was consistent with bi-allellic inheritance at a single locus[22]. Variation in the size and shape of the forewing band (Figure 3) was controlled by theSd locus[22], which operates by adding patches of melanin in the middle part of the forewing. The almost complete absence of any yellow scales in the forewings of the F1 was due to the activity of alternativeSd alleles that complement each other and induces melaninization either proximally or distally (Figure 3)

H. himera×H. e. erato: Variation in these crosses segregated in a way similar to theH. himera×H. e. etylus crosses. Specifically, we attributed most of the phenotypic variation to allelic variation at two loci, one affecting the presence of red and orange (D) in the forewing and hindwing and the other affecting the size and shape of the forewing band. The gene controlling the broken band phenotype ofH. e. erato was previously calledLy by Sheppardet al.[11] (see alsoFigure 1B andFigure 2). Variation of the orangeH. erato basal forewing and hindwing rays and the redH. himera hindwing bar resembled exactly the segregation observed in theH. himera×H. e. etylus cross for both F1 parents and F2 offspring. Variation in the size and shape of the forewing band again suggested the interaction of alternative complementary melanin alleles. The parental F1 had an intermediate phenotype with a reduced yellow pattern. More precisely there were five yellow spots specific to the parental F1s with their position defined by wing veins: 1) the discal spot, 2) a spot between Cu1b and Cu1a wing veins, 3) a spot between Cu1a and M3, 4) a spot between Sc and R1 wing veins, and 5) a spot between R1 and R2 wing veins (Figure 1B and2). Although the yellow spots in the F1 perfectly reflected the heterozygote intermediatehimera/erato phenotype, the forewing band size and shape variation became more complex in the F2 offspring. In particular, the extent that the forewing band extended towards the apex of the wing was highly variable, suggesting the presence of additional modifier loci (RP and WOM, personal observation).

H. himera×H. erato cyrbia: Most of the phenotypic variation in ourH. himera×H. e. cyrbia crosses could be explained by two major color pattern genes,Cr andD,[14],[20],[23] (see alsoFigure 1C Figure and 2). Red was controlled by theD locus and F1 individuals had red on both the fore and hindwing. In addition a second locus,Cr, was responsible for the yellow ventral line and white fringes of theH. e. cyrbia hindwing and interacted withD to determine the presence/absence of the forewing band. Similar toSd above, theCr locus largely acted by controlling melanic patches on both the forewing and the hindwing, withH. e. cyrbia Cr alleles positioning a large melanic patch across the middle of the forewing and portions of the hindwing; whereas theCr allele inH. himera put melanin across the entire hindwing. The epistatic interaction between theCr^hi andD^cyr alleles created some scoring difficulties forCr. When individuals were homozygous for theH. himera D allele, the forewing band was either completely absent (Cr^cyrCr^cyr), nearly completely covered with melanic scales (Cr^cyrCr^hi), or fully yellow (Cr^hiCr^hi). In contrast individuals heterozygous forD had a red forewing band[20]. There was evidence for additional variation in forewing band size and shape and for variation in the intensity of iridescence (RP and WOM, personal observation), but this variation was not easily scored qualitatively.

H. himera×H. e. notabilis: Qualitative variation in individuals fromH. himera×H. e. notabilis mapping families followed a simple segregation pattern consistent with bi-allelic variation at two major effect loci (Figure 1D andFigure 2). These broods also showed the most quantitative phenotypic variation among our collection of crosses (Figure 1D andFigure 3). The F1 parents of these crosses presented an intermediate phenotype with the red hindwing bar ofH. himera and a partial red forewing band ofH. e. notabilis (Figure 1C andFigure 2), presumably representing alternative alleles at the red locus (D). However, the amount of red present in the forewing band varied extensively in the F2, from nearly fully red to red restricted to the margins of the forewing bands. In addition, yellow pigments were always present in F2 individuals that had theH. himera red hindwing bar, whereas, individuals with red on both the forewing and hindwing, or on the forewing only, had white scale cells in the forewing bands. The geometric shape of the forewing band in ourH. himera andH. e. notabilis crosses was complicated as well. Previous crossing experiments implicated three loci,St,Sd andRound (Ro), as responsible for the distinctive split forewing band ofH. e. notabilis[11]. In our cross, the F1 individuals had an intermediate forewing band betweenH. himera andH. e. notabilis, with black melanin scales phenotypically masking yellow or white scales. However, there were typically some remnants of the upper forewing band in F1 individuals, represented by diffuse pigmentation around the edge of the middle forewing band near the intersection of Cu1a – M3 and M3 – M2 wing veins (Figure 1). F2 individuals showed variable expression of the upper forewing patch (Figure 3D).

Mapping acrossH. erato color pattern radiation

We examined 23 AFLP primer combinations (EcoCN/MseCNN) in several of ourH. erato reference crosses. Using this strategy we identified a number of AFLP markers tightly linked (1–3 cM) to the main color pattern genes (Sd^ety,Sd^not,St^not,D, andLy) segregating in our crosses, and used them as landmarks for cross- comparison of linkage betweenH. erato races and for BAC library screening[15],[28]. For theCr color pattern gene we utilized a gene-based marker linked toCr color pattern gene (“GerTra” -Rab geranylgeranyl transferase beta subunit,βggt-II)[14], which was successfully amplified and mapped in theF2-Not9[15].

To generate a high resolution map of wing color pattern loci across theseH. erato races we examined almost 1500 cumulative AFLP polymorphisms from three independentH. erato×H. himera crosses (F2-Not9.F2-Et2 andF2-He11). In theF2-Not9 cross alone, which was used as the backbone reference linkage map (Figure 4A), we scored 604 AFLP markers representing 245 male informative (MI), 213 backcross informative (BI), and 146 female informative (FI) respectively. After filtering the male and backcross informative AFLP markers for missing genotypes, undetectable phase, and segregation distortion (G test,P<0.05) 20 AFLP loci (5 MI and 16 BI) were discarded and did not enter the final dataset to construct theF2-Not9 linkage map (Figure 4A). In addition to these AFLP loci, we also mapped nine gene-based markers (Table S2). Roughly half of the total AFLP markers scored in our crosses wereH. himera alleles. The segregation ofH. himera alleles in the offspring of our crosses increased the chance to identify co-segregating markers shared between independent crosses ofH. erato (Table S3), and thus the possibility to use them as anchor loci to align the independent linkage analysis and test for color pattern gene homology (Figure 4A). Utilizing gene based markers and major color pattern genes as anchor loci we were also able to compare previous linkage maps inH. erato[22] andH. melpomene[21], which allowed cross-study assignment of linkage groups. As a result, we identified four linkage groups LG 1, LG 10, LG 18 and LG 15 from our analysis that were previously named LG 4, LG 3, LG 6 and LG 2 inH. e. etylus[22] (Figure 4A). For these linkage groups we used the chromosome numbers utilized inH. melpomene[21],[29] since it represent the first fully sequenced and reference genome amongHeliconius butterflies. The remaining linkage groups were ordered by decreasing size in centimorgans (cM)

Taking advantage of the absence of crossing over in female Lepidoptera[30],[31], we identified ‘chromosome prints’ for the 20 autosomal chromosomes and sex chromosomes using female informative (FI) markers grouped at LOD 8.0[22]. We also used unexpected discrepancies among FI AFLP markers[21],[22], given the lack of recombination in females, to infer the experiment-wide error rate, which was only ≈1%[15]. Amongst the female-informative markers, there were three markers that grouped at >LOD8.0 with sex that were positive only in females, thus identifying the W chromosome (WF-CC-CAC-182, WF-CT-CAG-140, WF-CT-CCT-153). We identified the sex (Z) chromosome from 8 FI markers present in males only. This and the autosomal chromosome prints obtained with FI markers were then used to build the final genetic map by assigning the linkage group for BI markers[21],[22] and using BI and MI markers grouped at LOD 4.0 or higher. Interestingly, we identified 12 BI markers grouping with the Z chromosome print. These 12 markers formed two non-mergeable linkage groups at LOD 4.0 with MI markers. These markers were genotyped in many of the same individuals, eliminating the usual cause for non-merged maps,[32] highlighting the possibility of divergence betweenH. e. notabalis andH. himera Z chromosomes in this study. Finally, our final mapped dataset consisted of 446 loci, 197 censored BI loci[21],[22], 240 MI loci and 9 gene based markers. Out of the total 437 AFLP loci, 171 markers generated 47 haplotypes and therefore reduced the unique breaks in the final genetic map to 313 (47 haplotypes and 266 single markers). Given this level of coverage we identified the maximum likelihood interval containing theD andSd loci to be around 3 cM and 5 cM respectively. If we assume a one-to-one correspondence between physical and recombination size[22],[23] the above interval can be translated to 0.8 Mb and ∼1.4 Mb for theD andSd color pattern genes respectively.

D andWh map to the same region of the genome

We characterized five AFLP loci, all in theF2-Not9 cross, tightly linked toD (Table S4). TheD-linked marker CAC491 also segregated in theF2-Et2 cross and was important for positional cloning theD red locus[15],[28]. Using a suite of nine co-dominant markers tightly linked to CAC491[28], we showedD is responsible for all red pattern variation across our collection ofH. erato crosses (Figure 5A).

Previous researchers[11] speculated that a distinct locus,White (Wh), unlinked to theD locus, controlled the presence of white scales in the forewing ofH. e. notabilis. However, in ourH. himera byH. e. notabilis crosses white and yellow color variation perfectly co-segregated in the two homozygote types: white pigment was never present in individuals that were homozygous for theH. himera D allele, whereas, all the individuals homozygous for theH. e. notabilis D allele showed white scales on their forewings (Figure 3). Heterozygous individuals had mainly white scales in the forewing bands, with the exception of 12 individuals that appeared to have a mix of white and yellow scales. QTL analysis confirmed this pattern and a QTL centered at theD locus largely explained the switch between white and yellow scale cell type in the forewing band (Figure 5C andTable 1). The LOD score around theD locus was so strong compared to any other QTLs that it resulted in the entire linkage group being significant (P<0.01;Figure 5C andFigure 4A). We also observed a second QTL for the last marker at end of LG18, which we believe most likely to be erroneous since it is known that chromosome ends might introduce false positives[33],[34].

Evidence of a major forewing shutter gene:St,Sd, andLy map to the same genomic region acrossH. erato races

A number of distinct loci have been described to explain some aspect of forewing band shape variation inH. erato[11]. To determine if there was evidence for multiple loci segregating in our crosses betweenH. himera and different races ofH. erato, we identified AFLP fragments tightly linked to the qualitative variation in forewing band shape across our collection ofH. himera×H. erato crosses. Among these,AAA145,ACA167 andAAA182 were tightly linked to variation in forewing band shape in theF2-Not9,F2-He11, andF2-Et2 mapping families. All three AFLP loci co-segregating in the three F2 crosses were excised and sequenced. Nucleotide sequences were nearly identical across the different color pattern races suggesting homology (Table S3) and these loci were used as anchors to align the three independent linkage analyses. Remarkably, variation in forewing band shape among the differentH. erato races mapped to the same relative position of LG 10 (seeFigure 6), suggesting that the same genomic interval was responsible for the majority of the forewing band variation observed acrossH. erato races. Moreover, the relative distance between the three AFLP fragments was almost identical in the different mapping families: a) distance betweenLAAA145 andACCA167 was 3 cM inF2-Et2, 2.2 cM inF2-He11, and 2.9 cM inF2-Not9; while betweenACCA167 andLAAA182 the relative distance was 8.8 cM inF2-Et2, 7.5 cM inF2-He11, and 7.1 cM inF2-Not9. The most likely placement of the three color pattern genes in each case was betweenACCA167 andLAAA182 – an interval of ∼2 Mb assuming 276 kb/cM[22]. Given this observation, we have collapsed the followingSd,St, andLy, which were previously thought to map to independent chromosomes or to very distant portion of the same linkage group, into a single relatively small region of the same chromosome. Although this region could easily contain hundreds of genes[15], in the light of recent studies[1],[3],[14], we speculate that our result most likely suggests a single forewing gene, which we callSd for continuity with previous studies.Sd controls all wing color pattern variation in the middle of the forewing.

Quantitative analysis of the forewing band variation in our reference crossF2-Not9 (H. himera×H. e. notabilis) was largely explained by theSd locus. QTL analysis ofH. notabilis band size variation where measured in the large central area of the forewing band characterized proximally by a Big-Spot (BS) and distally by a Not-Spot (NS) (seeFigure S1). These measures captures the complex shape changes of one versus two forewing patches and confirmed the major contribution of theSd locus on LG 10 for both spots of the forewing (Figure 6 andFigure 4). However while the analysis on BS showed a single QTL region centered at theSd locus (Figure 6: QTL region 1), the results for NS displayed two major QTL regions (Figure 6: QTL region 1 and QTL region 2), one of which matchesSd. The second QTL, more towards the bottom end of the LG 10 suggests the presence of an additional locus to delimit the size and shape of the upper forewing band. We hypothesized that this locus could beRo, which has been described as the color pattern gene affecting the distal portion of the upper forewing band (Not-Spot = NS)[11],[35] (see alsoFigure 7).

Finally, theCr locus, which segregates in ourH. himera×H. e. cyrbia crosses, was unlinked toSd and falls on LG 15 on our color pattern linkage map (Figure 4A). Mechanistically, this locus acts similar toSd by altering the distribution of melanic scale cells. InH. erato, it seems to be most important on the hindwing and was responsible for the yellow dorsal and ventral hindwing bar and white fringes. Nonetheless, theH. e. cyrbia Cr allele acted across the middle of the forewing band similar to the action ofSd.

QTL analysis identifies previously uncharacterized modifier loci

Despite the existence of major effect genes, a significant part of the variation observed in our broods had more complicated genetic control, which suggests the contribution of additional genes that may interact with the major effect loci. In order to understand the overall number and effect loci contributing toHeliconius erato wing color pattern variation, we performed a quantitative analysis on one of ourH. himera×H. e. notabilis mapping families. Among our collection of crosses we chose theF2-Not9 cross, because it presented some important advantages compared to the others: a) the most quantitative wing pattern variation, b) the segregation of variation previously ascribed to at least five color pattern genes (D, Sd, St, Ro, Wh), and c) the availability of a very high resolution genetic map. Two aspects to the wing color pattern variation in theF2-Not9 cross emerged from our analysis. First, there was ample variation in the amount of red in the forewing band (Figure 3A andFigure 3B). Our quantitative analysis identified four additional QTLs on different linkage groups (one on LG 4, two on LG 5, and one on LG 12) that helped explain this variation (Figure 4A,Table 1 andTable S5). When all QTLs were considered simultaneously, an epistatic model was supported over a purely additive model (p = 0.002) and explained 88.3% of the variation (Table 1). AfterD, the next two strongest QTLs, one on LG 5 and another on LG 4, respectively explained 15% and 10% of the variation in the amount of red under an epistatic model (Table 1). Second, we similarly found evidence that other regions of the genome interact withD in order to switch between yellow and white scale cell types (Figure 5C,Table 1, andTable S6). A QTL centered inD controls most of the white – yellow variation (Figure 4 andTable 1); however, there was a QTL of moderate effect size on the Z chromosome, which explains over 10% of the variation. In this case, the effect was purely additive, together explaining about 75% of the variation in white – yellow scale types (Figure 1B andTable 1).

There was also substantial variation in the size and shape of the forewing bands in offspring of our referenceH. himera×H. e. notabilis F2 brood (Figure 3C and 3B). QTL analysis of the relative size of the large central area of the forewing band (Big-Spot = BS), and the upper band (Not-Spot = NS) indicates the presence of a number of moderate sized (i.e. explaining more than 5% of variation) QTLs across the genome (Figure 4 andTable 1,Table S7 andTable S8). Many of these QTLs were quite large and explained as much variation in the relative size of the forewing band as the majorSd locus. For example, the variation in BS was controlled by four QTLs on four different linkage groups, including one on the chromosome containing theD locus (Figure 4), which acted largely additively to explain 63% variation in the size of this band (Table 1). Similarly, QTLs on three different chromosomes explained the variation in the upperH. e. notabilis forewing band (NS; seeFigure 4A). In this case, QTLs interacted epistatically to explain over 95% of the variation in the size of the HS2 band. Much of this variation was explained by a strong epistatic interaction between the major melanin pattern gene (Sd) on LG 10 and another QTL on the same chromosome. The two QTL peaks were close together (seeFigure 4A) and may reflect the action of a single large effect QTL. These QTLs, together with the interaction between them, explained about 82% of the variation in band size we observed in our crosses. In addition, there was evidence for moderate effect size QTLs on LG 2 and LG 7. Somewhat surprisingly, there was no evidence for any QTL near theCr, which affected the same region of the forewing band in ourH. himera×H.e. cyrbia cross.

Color variation inHeliconius erato

The yellow, orange, and red wing scales that compose the major wing pattern elements inHeliconius are colored by pigments derived from the tryptophan-ommochrome biosynthetic pathway[36]–[38]. Orange and red colors are ommochrome pigments synthesized in the developing scales cells, while the yellow pigment 3-hydroxykynurenine (3-OHK) circulates in the hemolymph and is taken up into scales cells shortly before adult emergence[36]. White scales appear to simply lack any color pigment[26],[27],[35]. White, yellow, orange, and red scales mature in a coordinated manner and share very similar structures. In contrast, melanic black scales develop much later and have distinct cuticle structures compared to non-black scales[39] (Reed pers. obs.). The correlation of developmental timing and cuticle structure between white, 3-OHK, and ommochrome scales implies that they may share a common specification mechanism. This hypothesis is consistent with our finding that inH. erato a major QTL centered on theD locus, is largely responsible for white, yellow, orange and red wing pattern variation (Figure 5C andFigure 4B).

TheD locus has previously been described as a complex of three tightly linked genes, but recent work strongly suggests that the homeobox transcription factoroptix is the sole protein-coding element responsible for the activity ofD[1]. Starting at ∼60 hours after pupationoptix expression precisely prefigures allD-controlled red, orange, and white wing patterns inH. erato, its co-mimicH. melpomene, and otherHeliconius species[1]. Interestingly,optix in situ hybridization matches perfectly the white and red forewing bands inH. melpomene plesseni, which is the co-mimic ofH. e. notabilis, the species used in our QTL study. Although noin situ expression data are currently available forH. e. notabilis, we would expect to observe an identical pattern as inH. m. plesseni given the strong correlation described between the spatial expression ofoptix in all theH. erato andH. melpomene wing tissues examined[1]. Overall, theoptix protein is highly conserved inHeliconius and the switch between red/orange and non-red/orange scales is explained by allelic variation ofoptix cis-regulatory elements[1]. The simplest model for howoptix might regulate color pattern variation imply that differentcis-regulatory elements of this gene read out different aspects of a conserved species-wide pre-pattern during pupal development. Thus, our results together with observations from previous studies suggest thatoptix acts as a selector gene to specify non-black fates of scale precursor cells where it is expressed.

Our study, however, cannot rule out that the large effect of chromosome 18 may fractionate into several QTLs controlling quantitative variation in pattern color composition independently ofoptix. The resolution of a QTL-mapping study depends on the number of offspring and the size of a QTL effects[40]–[42]. Thus, the combination of the small sample size of this study (88 offspring) and the large effect size of theoptix gene restricts our ability to distinguish between a single QTL and a series of tightly linked QTLs[41],[42]. The quantitative variation explained by only theD/optix locus on chromosome 18 under the additive model is significantly large with over 60% for the amount of red and roughly 85% for the proportion of non-red scales that are either white or yellow (Table 1). In spite of these limitations, our QTL study clearly suggest that theD andWh loci, previously speculated to be unlinked inH. erato[11], are actually mapping to a the same region of LG 18. Interestingly, the yellow/white switch inH. cydno is modulated by a different locus, which is linked to thewingless gene[17] and maps our LG 1 (Figure 4). There is no evidence that this linkage group explains any of the variation in white and yellow observed in our crosses, which is nearly completely explained by the large effect QTL in chromosome 18 and a sex-linked QTL (Figure 4B andTable 1). Thus, differentHeliconius species have evolved the ability to modulate the shift between yellow (3-OHK pigment) and white (no pigment) scales cells in at least two distinct ways.

Despite the major control of theD/optix locus on qualitative and quantitative color variation, other regions of the genome play an important role in definingH, erato wing color pattern. Of particular significance is the presence of three QTLs on three different chromosomes. These loci explained nearly 30% of residual forewing variation in the amount of red, which varied from nearly all red to nearly all white/yellow under the epistatic model (Figure 3;Table S5). Two of these QTLs, on LG 5 and LG 4 respectively, accounted for the bulk of this additional variation and each of them interact epistatically withD to reduce the amount of red across the forewing band. As mentioned above, these results are consistent with a recent study demonstrating thatoptix is expressed in the developing wings ofH. m. plesseni (i.e. the comimic ofH. e. notabilis) with a spatial distribution that prefigure the red/white upper and lower forewing bands. Thus, whileoptix acts as a major determinant of non-black wing scale identity, additional modifiers interact with the pathway regulated byoptix. This implies that there are additional loci that repress pigment synthesis in order to temper the final color fate. We speculate these color-repressing modifiers may modulate the expression of transporters or enzymes required for ommochrome pigment synthesis.

Forewing band variation inHeliconius erato

Our crossing experiments included races that represent most of the diversity in forewing band size and shape withinH. erato. Sheppard and colleagues[11] had previously assigned differences in forewing band shape variation inH. erato to a number of distinct loci (Figure 2). However, similar to the red locus, integrative mapping experiments confirm that a single genomic interval, and possibly a single locus centered onSd, is responsible for much of the adaptive variation in forewing band shape, size and position (Figure 6). As discussed for the QTL centered on theD locus, also in this case we cannot exclude the possibility of multiple tightly linked loci aroundSd. In particular, for the relative size of the upper forewing band (NS) there is a peak in LOD scores atSd followed by a curious step rise in values towards the end of the same chromosome (Figure 6). One interpretation of this pattern is that at least two QTLs on LG 10 account for the observed variation: the one centered onSd has an effect on both mid and upper forewing patches whereas the second QTL, at the end of the chromosome effects uniquely the size of the upper forewing band (Figure 6 andFigure 7). As observed for the yellow, red, and white color pattern variation it is likely that the developmental toolkit for altering melanic patterns onHeliconius wings involves the activity of additional loci. At least seven different QTLs across five different linkage groups collectively accounted for differences in the distribution of melanic black pattern elements across the upper and lower forewing.

Recently, high resolution linkage mapping, expression analysis and pharmacological experiments unambiguously suggest theHeliconius homolog ofWntA as theSd locus[19]. The spatial pattern ofWntA expression corresponds to the black pattern[19] and suggests thatWntA/Sd controls the black “shutter” effect produced by pattern expansion from the organizing centers (Fig. 7). Similarly tooptix, variation in thecis regulatory region ofWntA is responsible for the adaptive evolution of black forewing patterns inH. erato[19].WntA is a signaling ligand that creates a morphogen gradient across the developing wing tissue, which possibly induce and maintain the expression of different target genes at distinct concentration thresholds.WntA is transcribed during last-instar imaginal disc development, although pharmacological manipulations suggest thatWntA signaling actually occurs shortly after pupation. Thus,WntA activity precedes color pattern expression ofoptix by a few days. Such observations, leads us to hypothesize thatWntA may be a direct or indirect negative regulator ofoptix, and that the interaction between the two genes is largely responsible for establishing black vs. non-black wing pattern boundaries. Moreover, recent transcriptomic work shows that late-pupal expression of some ommochrome-related genes closely corresponds tooptix expression[43]. Hence, whileWntA possibly controls a variety of molecules during early development to create pattern boundaries that do not directly activate melanin depositionoptix acts directly on the ommochrome pathway switch later in development. We speculate thatWntA andOptix interaction could provide the fundamental mechanisms to control the natural diversity observed inHeliconius wing color patterns.

WntA is exactly the type of molecule predicted by theoretical and developmental models for color pattern formation in butterflies. An important organizing principle of butterfly color patterns is the nymphalid groundplan (NGP), first deduced by Schwanwitsch[44] then integrated and expanded by Nijhout[35]. The NGP proposes that color pattern development is controlled by an array of organizing centers distributed across the wing. Signals from these organizing centers determine the boundaries of the various elements of the color pattern. These organizing centers are distributed in irregular rows that run roughly in an antero-posterior direction across a wing[45]. Among the best-studied of these are the organizing centers for eyespot patterns[46]–[51]. Notch signaling and the expression of distal-less are activated at these centers[46],[52], followed by the expression of spalt and engrailed[48]. Depending on the species, different combinations of these transcriptional regulators determine the black central disk and the outer rings of an eyespot[48],[52],[53].

The bold black regions on the wings ofHeliconius andCharaxes can be shown to come about by the confluence of smaller black pattern elements that expand from organizing centers that are homologous to those that organize eyespot patterns[25],[26],[54]. The locations of these organizing centers, and the manner in which they expand and fuse to form large areas of black pigmentation, are presented inFigure 7. The rows of organizing centers inHeliconius wings are labeled c, d, f and h inFig. 7A, were each letter represents a set of serially homologous centers corresponding to elements of the NGP[25],[35],[55]. Some of the genes that control color patterns inH. erato andH. melpomene have been shown to affect the degree of expansion of these black areas[26]. Gilbert[27] referred to these genes as “shutter genes,” an elegant metaphor that describes how they cause black portions of the pattern to expand and cover up the more brightly pigmented background. This groundplan provides an explanation of the developmental mechanism responsible for the phenotypic effects of theWntA/Sd gene (Fig. 7B). Each of theWntA/Sd alleles appears to affect the degree of distal expansion of black pattern from one (f) or two (d and f) of the serially homologous organizing centers in either the anterior compartment of the wing (anterior to vein M₃[45]), or in both compartments. Thus each allele has a characteristic effect on the same subset of organizing centers. This action is entirely driven by changes cis-regulation ofWntA[19], which effects in which location of the wingWntA expression occurs. Finally, it is noteworthy to mention that there are several additionalWnt family members on the same chromosome (LG 10) nearbyWntA that might explain the two distinctive QTL peaks (Sd/WntA and possibly a color pattern gene namedRo) that we observe in our analysis on such chromosome. These loci become interesting candidates for future developmental analysis.

Comparative mapping and expression experiments suggest that theSd/WntA locus is homologous to theAc locus inH. cydno[17],[19].Ac controls the presence/absence of melanic scales on central portions of the forewing and hindwing inH. cydno[17]. Similar toH. erato, several loci inH. melpomene andH. cydno have been proposed to be responsible for variation in melanin patterns across the forewing and hindwing of different color pattern races of each species. For example, distinct studies list several different loci that have some effect on the shape of the forewing banding patterns inH. melpomene[11],[26],[35] andH. cydno[26],[35]. Similarly, our expectation forH. cydno andH. melpomene is that the number of unique loci involved in black patterning within each adaptive radiation will collapse to a much smaller set of loci, each containing large effect alleles, similar to what we observe inH. erato.

The genetic architecture of adaptation and speciation

By taking advantage of a wide range of crosses generated with distinctH. erato morphs, our study describes the genetic architecture underlying quantitative wing color pattern variation inH. erato at greater detail than has been previously possible. Our work expands upon two previous studies[56],[57] that measured quantitative wing color pattern variation inHeliconius butterflies in order to increase the overall knowledge on the genetic control of this diversity. Although the QTL analysis from inH. melpomene[56] andH. numata[57] used relatively few markers and were largely restricted to chromosome level analysis, their main conclusions are very similar to ours (Table 2). The broadly interpretation of these previous studies and ours is consistent with the presence of a small conserved set of major effect genes that interact with a larger set of “modifiers” in order to modulate wing pattern variation inHeliconius butterflies. For example, Baxter and collaborators[56] identified six chromosomes significantly affecting the size variation of the red forewing band inH. melpomene (Table 2). Among these chromosomes with significant QTLs inH. melpomene, the ones for which we can confidently identify homology withH. erato are LG 18 (D) and LG 10 (Sd) and sex Z chromosome (Table 2). The second quantitative assessment of color pattern inHeliconius wings[57] focused on the variation observed in two broods ofH. numata; a species where a single supergene namedPushmipullyu (P)[58] control almost entirely its natural wing pattern diversity. Unexpectedly, this study was able to demonstrate that other genomic regions beside theP locus, contribute to the overall variation observed in theH. numata wing patterns (Table 2). Similarly to the QTL analysis inH. melpomene[56], Jones and collaborators[57] found that the regions homologous to theD andSd loci ofH. erato had a significant effect onH. numata quantitative wing pattern variation (Table 2). Overall our study adds another level to the analysis, and somewhat paradoxically while simultaneously demonstrating that many previously described gene-loci collapse into alternate alleles of a single locus we have identified new QTL's on other chromosomes that modify the effect of the major wing patterning loci in these species.

In a broader context, the overall genetic architecture ofHeliconius color pattern variation described in this study supports the two-step process of Müllerian mimicry evolution, as proposed by Punnett[59], elaborated by Nicholson[60], and finally presented with a more complete heuristic view by Turner[61]–[63]. This model entails that the first large effect mutation allows the future mimic to cross a fitness valley from a lower adaptive peak to a higher adaptive slope determined by a warningly colored model species. Then, smaller changes will permit climbing, over time, to the highest adaptive peak imposed by the model and generate a perfect mimetic resemblance. Turner adopted this two-step model to explain how pattern evolution inHeliconius butterflies would occur under the selective regime imposed by vertebrate predators. Turner[64],[65] and others[39],[60] recognized that only a finite number of adaptive peaks exist in nature, and that these peaks inHeliconius butterflies are easily visualized as distinctive Müllerian mimicry rings (i.e. groups of distantly related species that share a common warning signal). Each of these rings can contain a large number of distantly related species but only a handful of distinctive rings coexist in any given location[63]. Thus, Turner proposed that mutations of large effect are required to move individuals across phenotypic space between existing mimicry rings. In order to be fixed by selection, such mutations need to be capable of generating a pattern that falls under the protection of an existing mimetic pattern. Then, once a phenotypic shift has occurred, mutations of smaller effect can hone the pattern and perfect mimicry[62],[66]. Several theoretical models have also recently been described in favor of a multi-peaks adaptive landscapes[67]–[70] and field experiment that demonstrated adaptive peaks with valley between them[13]. Moreover, there is a long history field experiments inHeliconius that show strong selection against oddly patterned individuals[13],[71],[72] necessitating initial steps of mimetic resemblance to be ‘large’ and generate a ‘rough & ready’ appearance for any hope of mimetic phenotypes to evolve.

As a final remark, our study provides additional data and insights into the genetic mechanisms driving mimicry evolution and more in general the genetic architecture of adaptive variation. The central debate on mimicry evolution has always focused on the number of steps (few versus many) and their effect size (small versus large) needed in order to resemble the model. Overall, our data are consistent with recent theoretical work modeling the distribution of effect sizes over the course of an adaptive walk towards a local optimum[66],[73]. One of the most robust conclusions from these models is that only a few mutations are expected to account for most phenotypic change. Also, the overall size effect of the mutations that are substituted over entire bouts of adaptation fits an exponential distribution. Such a distribution implies that while the first mutations have large phenotypic effects the subsequent ones decrease exponentially their effects as they accumulate. Under this prediction a population should take few larger steps when far from the optimum, but many small ones when near. A number of empirical studies have been borne out of the above expectations[3],[74]–[77]. This new theoretical understanding coupled with emerging empirical examples underscore a broad conceptual change in our understanding of the genetic basis of adaptation since Fisher's first models over 70 years ago[66],[78]. In particular, one of Fisher's main arguments, is that mutations of small phenotypic effect are more likely to be beneficial than are those of large phenotypic effect. What our data suggests, is that in wing mimicry systems small effect mutations are used to optimize only the characteristics of the phenotype controlled by a major effect gene that needs improvement while preserving, or slightly effecting, the rest. Ultimately, our comprehension of how adaptive novelty arises will become more robust as we integrate the developmental events underpinning the formation of specific morphologies with evolutionary and genetic theory (sensu Stern)[79].

Sample collection and crossing strategy

We generated F₂ and backcross mapping families by crossing four different geographic races ofH. erato to the same stock ofH. himera (Figure 1). Races included: 1)H. erato notabilis (north-eastern Ecuador), 2)H. erato etylus (South-Eastern Ecuador), 3)H. erato cyrbia (western Ecuador), and 4)H. erato erato (French Guiana). This design allowed us to scoreH. erato color patterns against the common phenotypic background ofH. himera, which has a fully developed middle forewing yellow band and a red hindwing bar. All crosses were carried out in theHeliconius Insectaries at the University of Puerto Rico from stocks originally collected in the wild. Individuals used to established our stocks were collected from the following locations:H. himera, Vilcabamba, Ecuador (79.13 W, 4.6 S);H. e. notabalis, Puyo, Ecuador (78.0 W, 1.5 S);H. e. cyrbia, Guayquichuma Glen, Ecuador (79.6 W, 3.9 S);H. e. etylus, Zamora River, Ecuador (78.5 W, 3.55 S);H. e. erato, Maripasoula, French Guiana (54.03 W, 3.64 N). Butterflies were cared for as outlined in[80]. Offspring in all mapping families were reared individually and raised until eclosion. Adult butterflies were uniquely identified with a number, their bodies stored at −70°C for genetic analysis, and high-quality standardized pictures of their wings taken with a Nikon Coolpix 995 digital camera and a single Vivitar 2400M external slave flash mounted on a custom arm held at a 30 degree angle with respect to the plane of focus. The camera was fixed in a custom-built copy stand, using constant illumination and camera settings (F10, 1/2000^th shutter speed, 20.2 mm focal length, constant distance from subject and a neutral gray background overlayed with a with 5 mm grid, in sRGB color space) for further morphological analysis.

Scoring wing pattern phenotypes

Genotypes for the color pattern genes segregating in our crosses were scored based on qualitative measures[12],[20]. Although qualitative scoring does not capture all of the segregating variation (Figures 1 andFigure 3), most color pattern difference within individual crosses can be explained by simple allelic variation at a handful of loci. We scored our broods for a total of six previously named color pattern loci (Sd,St,Ly,D,Cr, andWh) that segregated in our collection of crosses (seeFigure 1 andFigure 2). A short description of the action of these genes, based on the work of Sheppard and collaborators[11] and elaborated in several studies[81],[20],[22] is illustrated inFigure 2 and described below.

Action of major color pattern genes

The color pattern variation segregating in our crosses has been previously ascribed to the following hypothesized loci:D,Wh,Sd, andSt, inH. e. notabilis×H. himera;D andSd inH. e. etylus, ×H. himera;D andLy inH. e. erato×H. himera;D andCr inH. e. cyrbia×H. himera (Figure 1 andFigure 2). Brief descriptions of the major phenotypic action of each locus are as follows:

D (modified fromDRY): This locus controls red color patterns inH. erato and was originally proposed by Sheppardet al.[11] to be a “supergene” of three tightly linked lociD(ennis), R(ays), and Y(ellow band).D controls the red/orange at the base of the forewing;R organizes the red/orange rays in the hindwings; andY controls the color of the forewing band, which can be yellow or red. Recently, we presented strong genetic and gene expression evidence that the three loci are cis-regulatory elements of the transcription factoroptix[1]. For simplicity and, in keeping with more recent nomenclature, we refer to theDRY locus asD.

Wh (White): This locus controls the abundance of white scales on the forewing band present in the east Ecuadorian raceH. e. notabilis. Sheppardet al.[11] hypothesized an interaction and/or linkage of theWh locus with another color pattern locus, but they were unable to assign theWh locus to a linkage group.

Sd (Short band) andSt (Split band): Sheppardet al.[11] characterized two tightly linked loci,Sd andSt, in crosses that involved the east Ecuadorian raceH. e. notabilis.Sd is responsible for shortening the forewing band relative to otherH. erato races while theSt gene controls the splitting of the forewing band into two portions (Figure 2). The nameSd was also used to describe the locus controlling forewing band variation inH. e. etylus×H. himera crosses[22] but without evidence of homology withSd as defined in the original cross.

Ly (Broken band): The major phenotypic action of this gene is to alter the forewing band shape, causing it to be fractured relative to the solid forewing band of otherH. erato races. The locus was described from crosses between Suriname hybrids (H. e. amazona×H. e. hydara) and population from Trinidad (H. e. hydara) and was subsequently confirmed by crossing populations ofH. e. amalfreda from Manaus, Brazil andH. e. phyllis from East Brazil[11],[82] (seeFigure 2).

Cr (Cream rectangles): This locus controls a variety of characteristics. It was originally described for its effect of weakening the expression of the forewing yellow line ofH. e. phyllis and the yellow hindwing bar ofH. e. favorinus[11],[81].Cr also controls a series of white marks along the edge of forewings and hindwings (i.e. H. e. cyrbia) and has a pleiotropic effect on the amount of red pigment in the forewing band[20].

Characterization of wing pattern-linked markers

We used Amplified Fragment Length Polymorphisms (AFLPs)[83] to identify a set of markers across theH. erato genome. Genomic DNA was digested withEcoRI andMseI endonucleases and a total of 23 AFLP primer combinations (EcoCN/MseCNN) were analyzed across several crosses using previously described methods[22],[23],[84]. Linkage analyses were performed using previously described methods[21],[22], and AFLP markers of interest were cloned and sequenced[15]. Three crosses -H. himera×H. e. etylus,H. himera×H. e. notabilis, andH. himera×H. e. erato - were screened for co-segregating AFLP loci (i.e. identical fragments segregating in different crosses) on each linkage group containing a color pattern gene. Co-segregating AFLPs linked to pattern genes were excised from gels, cloned size-verified, and sequenced (Table S3). In addition, we also isolated and sequenced an AFLP band (GCAA303) linked toSd^ety identified in a previous mapping work usingH. himera×H. e. etylus broods[22]. Using this strategy we identified a number of AFLP markers tightly linked (1–3 cM) to theSd^ety,Sd^not,St,D, andLy color pattern genes segregating in our crosses and used them as landmarks for linkage comparison betweenH. erato races. A total of 18 AFLP loci within a 3 centimorgan (cM) window from each targeted color pattern gene were identified (Table S4). To supplement the AFLP markers we also designed PCR primers to score several additional anchor loci:D23/24, Has1, Thap, GPCR, optix[1],VanGogh DNAJ, Gn12, Gn47[28],GerTra,EIF3b[14],Wg,RPL3,RPL22,PTC[22], andB303 (SeeTable S2)

Linkage analysis and map construction

Linkage analyses were performed as previously described in Kapanet al. methods[22]. Briefly, we constructed a linkage map using the three-step process first outlined by Jigginset al.[21] and expanded upon by Kapanet al.[22]. Due to a lack of recombination during oogenesis in heterogametic females[30],[31], bands from female-informative (FI) AFLP markers were used to identify linkage groups at LOD 8.0 using JoinMap 3.0[85]. We then used JoinMap 3.0 to compare each “chromosome print” with back-cross informative (BI) AFLPs, as well as codominant anchor loci with segregation ratios of 1∶1∶1∶1 or 1∶2∶1 (see Kapanet al.[22]) using groups that form at LOD 4.0 or higher with each chromosome print. These new groups were then individually examined to verify linkage phase and to identify “forbidden recombinants”[86], which are defined as offspring genotypes that could only appear by crossing-over in females. Presence of forbidden recombinants implies that the grouping is incorrect, unless there are scoring errors. In our analysis we accepted a small amount of scoring error and retained BI AFLP loci with five or fewer forbidden recombinants that have a probability ofp<0.007 or lower of being unlinked under the hypothesis of no scoring error (see Jigginset al.[21]). Finally we extracted the male-informative (MI) component from the BI AFLP scores (i.e. following only AFLP bands inherited from the father) by “censoring” the genotypes from AFLP bands inherited from the mother[21]. The censored BI markers were combined with MI AFLPs and remaining co-dominant anchor markers were re-coded to show only the male-informative allele. For the sex chromosome Z, male offspring are always positive for BI AFLPs and the recombinant analysis is limited to female offspring segregating male-informative markers on Z[21]. For the final map construction linkage groups were assembled at LOD 3.0 utilizing Joinmap 3.0 followed by the use of Mapmaker 3.0 to generate the most likely order and spacing of markers on a chromosome[87].

Quantification of wing color and size band patterns

Heliconius butterflies wing color pattern is not a perfect arrangement of monochromatic areas but it is the result of the ratio between scales of different colors. Such ratio is what determines the final wing color pattern. We used our most phenotypically variable brood as the reference cross (H. himera×H. e. notabilis: F2-Not9), to quantify variation in 1) the amount of red across the central forewing patch, 2) the amount of white versus yellow scales (Figure S2), and the size of the upper and lower forewing bars (Figure S1). The amount of red and white scales was extremely variable depending on the area of the forewing band considered. Wings that appeared purely red had a significant portion of white scales and, vice versa wings that seems all white had many red scales. Thus, to estimate redness of the central forewing patch, we used the program SigmaScan® to calculate the intensity of the blue, green and red emission using the Red Green Blue (RGB) color model across the entire forewing band. Therefore, we measured the amount of red, blue and green pixels and calculated “redness” as red emission divided by blue plus green emission (i.e. redness = red/(blue+green)). This measure provided a reasonably robust quantitative estimate of the amount of red in the whole forewing band. However, the RGB data obtained could not be used to distinguish the relative proportion of white versus yellow scale cells in the forewing band due the limitations of digital photographs to capture differences between yellow and white. To obtain a quantitative measure of the relative proportion of white and yellow scales, we used a spectrophotometer to measure the absorption spectrum (seeFigure S2). The measured region was in the anterior distal region of the middle right forewing patch (seeFigure S2) and was chosen because red scales were mostly absent and thus did not interfere our ability to quantify “whiteness” (or yellowness) of the scale cells in that region. We used the ratio of absorbance at 380 nm and 700 nm as our quantitative measure of “whiteness” (Figure S2). In order to compare the effects of each quantitative trait locus on the different phenotypic measurements (redness, whiteness, yellowness, Big-Spot and Not-Spot) we report the percentage variance explained for each QTL.

Heliconius himera andH. e. notabilis differ notably in the shape of the forewing band (seeFigure 1 andFigure 2). Specifically, the forewing band ofH. e. notabilis is split with a small patch of white/red scale cells displaced distally towards the upper margins of the forewing. In addition, the central forewing band ofH. e. notabilis is smaller than the forewing band ofH. himera. To characterize variation segregating in our focal F2 family, we scored a series of landmarks on standardized wing photos of each offspring (see above) using the computer program DIGIT (courtesy of H. F. Nijhout). DIGIT generates positional (X,Y) coordinates from points on the image. In total, we scored 58 positions across the wing, including 14 “landmark” coordinates chosen to reflect the positions of the major venation patterns and 44 coordinates that reflected the size and position of the proximal and distal forewing bands (seeFigure S1). For our analysis we focused specifically on variation in the relative size of two major forewing patches (Figure S1). The Big-Spot (BS) represented the largest portion of the proximal forewing band and the Not-Spot (NS) represented the bulk of the distal band. Although our method measured the relative size variation at the forewing bands by calculating the area for the two bands in each individual, we argue that this measure represents a close approximation of the forewing band shape as well. Quantitative differences size of the BS and NS were measured as an area of the band standardized by the area of the wing, which could be considered as a biological proxy for shape. The area of each spot was calculated as an area of a polygon, with vertices represented by landmarks (Figure S1).

QTL analysis of color and size data

All QTL analyses were performed in the R statistical package version 2.13.2[88]. We tested the association between the molecular markers, color and band size and shape phenotypes using a linear regression model with results extracted using the ANOVA function. Association was tested by a simple one-way ANOVA with phenotypes treated as dependent variables, and markers as independent variables. The effect of the marker was assessed via the F ratio, which is equivalent to a LOD score in single marker analyses[89]. Statistical significance of the LOD score was assessed via 5,000 matrix permutations[90]. Permutations resulted in the generation of a distribution of LOD scores, and based on this distribution, only those associations with actual LOD scores having less than 0.01 probability of being observed by chance were considered significant. This probability corresponded to LOD score cut-off values of 3.8 or higher depending on all the linkage groups analyzed. Single marker analyses results were confirmed in the program QTL Cartographer[91] (data not reported).

Since the order and position of markers are known in each linkage group, adjacent markers associated with a phenotype were assumed to present the same QTL. After identification of QTLs associated with variation in color and forewing band size, we ran another one-way ANOVA with only one representative marker per QTL. Markers were chosen on the basis of having the highest explanatory power (highest LOD value) within the QTL. For each phenotype, the selected markers were analyzed individually and then together assuming additivity between markers in a linear regression model. The amount of phenotypic variance explained by each marker and each model was estimated as the r-squared value of the linear regression. We then tested if an additive model better explained the phenotypic variation than a model with just one marker (one-way vs. multi-way ANOVA).

The simplest model of genetic architecture assumes that contributions of individual markers to the phenotype are additive (Fisherian-type genetic architecture)[78]. A more complex model allows for epistasis (Wrightian-type genetic architecture)[92], where the relative contribution of each marker to the phenotype could be dependent on the presence of allelic variants at other markers. To test the explanatory power of these two models, we constructed a purely additive model, which included the marker that has the highest explanatory power in that QTL. Sinceoptix is such a strong candidate gene for theD locus[1], we tested whether the phenotypic action of the known major effect locus impacting red/orange coloration interacted with other unknown QTLs. We thus examined whether or not interactions betweenD and other makers were significant, and if the epistatic model better explained the data than a purely additive model (additive vs. interaction ANOVA). Similarly, for variation in forewing band size we compared a purely additive and a more complex epistatic model using one marker per QTL and building our model around the marker with the highest explanatory power.

• 1.Reed RD, Papa R, Martin A, Hines HM, Counterman BA, et al. (2011)optix drives the repeated convergent evolution of butterfly wing pattern mimicry. Science333: 1137–1141. [DOI] [PubMed] [Google Scholar]
• 2.Joron M, Frezal L, Jones RT, Chamberlain NL, Lee SF, et al. (2011) Chromosomal rearrangements maintain a polymorphic supergene controlling butterfly mimicry. Nature477: 203–206. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 3.Papa R, Martin A, Reed RD (2008) Genomic hotspots of adaptation in butterfly wing pattern evolution. Curr Opin Genet Dev18: 559–564. [DOI] [PubMed] [Google Scholar]
• 4.Chamberlain NL, Hill RI, Kapan DD, Gilbert LE, Kronforst MR (2009) Polymorphic butterfly reveals the missing link in ecological speciation. Science326: 847–850. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 5.Quek SP, Counterman BA, Albuquerque de Moura P, Cardoso MZ, Marshall CR, et al. (2010) Dissecting comimetic radiations inHeliconius reveals divergent histories of convergent butterflies. Proc Natl Acad Sci U S A107: 7365–7370. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 6.Hines HM, Counterman BA, Papa R, Albuquerque de Moura P, Cardoso MZ, et al. (2011) Wing patterning gene redefines the mimetic history ofHeliconius butterflies. Proc Natl Acad Sci U S A108: 19666–19671. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 7.Pohl N, Sison-Mangus MP, Yee EN, Liswi SW, Briscoe AD (2009) Impact of duplicate gene copies on phylogenetic analysis and divergence time estimates in butterflies. BMC Evol Biol13: 99. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 8.Cuthill JH, Charleston M (2012) Phylogenetic codivergence supports coevolution of mimeticHeliconius butterflies. PLoS One7: e36464. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 9.Brown KS, Sheppard PM, Turner JRG (1974) Quaternary refugia in tropical America: Evidence from race formationinHeliconius butterflies. Proc R Soc Lond B187: 369–378. [Google Scholar]
• 10.Brown KS (1979) Ecologia Geográfica e Evolução nas Florestas Neotropicais. 2 vols. Universidade Estadual de Campinas, Campinas, Brazil.
• 11.Sheppard PM, Turner JRG, Brown KS, Benson WW, Singer MC (1985) Genetics and the evolution of Müllerian mimicry inHeliconius butterflies. Phil Trans R Soc B308: 433–610. [Google Scholar]
• 12.Mallet J, Barton N (1989) Strong natural selection in a warning-color hybrid zone. Evolution43: 421–431. [DOI] [PubMed] [Google Scholar]
• 13.Kapan DD (2001) Three-butterfly system provides a field test of Müllerian mimicry. Nature409: 338–340. [DOI] [PubMed] [Google Scholar]
• 14.Joron M, Papa R, Beltran M, Chamberlain N, Mavarez J, et al. (2006) A conserved supergene locus controls colour pattern diversity inHeliconius butterflies. PLoS Biol4: e303. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 15.Papa R, Morrison CM, Walters JR, Counterman BA, Chen R, et al. (2008) Highly conserved gene order and numerous novel repetitive elements in genomic regions linked to wing pattern variation inHeliconius butterflies. BMC Genomics9: 345. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 16.Turner JRG, Crane J (1962) The genetics of some polymorphic forms of the butterfliesHeliconius melpomene (Linneaus) andH. erato (Linneaus), I: major genes. Zoologica1962: 141–152. [Google Scholar]
• 17.Kronforst MR, Kapan DD, Gilbert LE (2006) Parallel genetic architecture of parallel adaptive radiations in mimeticHeliconius butterflies. Genetics174: 535–539. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 18.Baxter SW, Papa R, Chamberlain N, Humphray SJ, Joron M, et al. (2008) Convergent evolution in the genetic basis of Müllerian mimicry inHeliconius butterflies. Genetics180: 1567–1577. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 19.Martin A, Papa R, Nadeau NJ, Hill RI, Counterman BA, et al. (2012) Diversification of complex butterfly wing patterns by repeated regulatory evolution of a Wnt ligand. Proc Natl Acad Sci U S A109: 12632–12637. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 20.Jiggins CD, McMillan WO (1997) The genetic basis of an adaptive radiation: warning colour in twoHeliconius species. Proc R Soc Lond B264: 1167–1175. [Google Scholar]
• 21.Jiggins CD, Mavarez J, Beltran M, McMillan WO, Johnston JS, et al. (2005) A genetic linkage map of the mimetic butterflyHeliconius melpomene. Genetics171: 557–570. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 22.Kapan DD, Flanagan NS, Tobler A, Papa R, Reed RD, et al. (2006) Localization of Müllerian mimicry genes on a dense linkage map ofHeliconius erato. Genetics173: 735–757. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 23.Tobler A, Kapan D, Flanagan NS, Gonzalez C, Peterson E, et al. (2005) First-generation linkage map of the warningly colored butterflyHeliconius erato. Heredity94: 408–417. [DOI] [PubMed] [Google Scholar]
• 24.Mallet J, McMillan WO, Jiggins CD (1998) Mimicry and warning color at the boundary between races and species. In Endless forms: species and speciation Howard D.J, Berlocher S.H 1998. Eds New York, NY:Oxford: 390–403.
• 25.Nijhout HF, Wray GA (1988) Homologies in the color patterns of the genusHeliconius (Lepidoptera: Nymphalidae). Biol J Linn Soc33: 345–365. [Google Scholar]
• 26.Nijhout HF, Wray GA, Gilbert LE (1990) An analysis of the phenotypic effects of certain color pattern genes inHeliconius. Biol J Linn Soc40: 357–372. [Google Scholar]
• 27.Gilbert LE (2003) Adaptive novelty through introgression inHeliconius wing patterns: evidence for a shared genetic “tool box” from synthetic hybrid zones and a theory of diversification. In: Boggs CL, Watt WB, Ehrlich PR, editors. Ecology and evolution taking flight: Butterflies as model systems. Chicago: University of Chicago Press. pp. 281–318.
• 28.Counterman BA, Araujo-Perez F, Hines HM, Baxter SW, Morrison CM, et al. (2010) Genomic hotspots for adaptation: the population genetics of Müllerian mimicry inHeliconius erato. PLoS Genet6: e1000796. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 29.Heliconius Genome Consortium (2012) Butterfly genome reveals promiscuous exchange of mimicry adaptations among species. Nature487: 94–98. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 30.Suomalainen E, Cook LM, Turner JRG (1973) Achiasmatic oogenesis in the heliconiine butterflies. Hereditas74: 302–304. [Google Scholar]
• 31.Turner JRG, Smiley J (1975) Absence of crossing-over in female butterflies (Heliconius). Heredity74: 265–269. [DOI] [PubMed] [Google Scholar]
• 32.Van't Hof AE, Marec F, Saccheri IJ, Brakefield PM, Zwaan BJ (2008) Cytogenetic characterization and AFLP-based genetic linkage mapping for the butterfly Bicyclus anynana, covering all 28 karyotyped chromosomes. PLoS One3: e3882. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 33.Crooks L, Nettelblad C, Carlborg Ö (2011) An Improved Method for Estimating Chromosomal Line Origin in QTL Analysis of Crosses Between Outbred Lines. G3: Genes, Genomes, Genetics1: 57–64. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 34.Browman KW, Sen S (2009) A Guide to QTL Mapping with R/qtl. New York: Springer.
• 35.Nijhout HF (1991) The Development and Evolution of Butterfly Wing Patterns. Washington, DC: Smithsonian Institution Press.
• 36.Reed RD, McMillan WO, Nagy LM (2008) Gene expression underlying adaptive variation inHeliconius wing patterns: non-modular regulation of overlapping cinnabar and vermilion prepatterns. Proc Biol Sci B275: 37–45. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 37.Linzen B (1974) The tryptophan–ommochrome pathway in insects. Adv Insect Physiol10: 117–246. [Google Scholar]
• 38.Koch PB (1991) Precursors of pattern specific ommatin in red wing scales of the polyphenic butterfly Araschnia levana L.: haemolymph tryptophan and 3-hydroxykynurenine. Insect Biochem21: 785–794. [Google Scholar]
• 39.Gilbert LE, Forrest HS, Schultz TD, Harvey DJ (1988) Correlations of ultrastructure and pigmentation suggest how genes control development of wing scales of Heliconius butterflies. J Res Lepidoptera26: 141–160. [Google Scholar]
• 40.Falconer DS, Mackay TFC (1996) Introduction to Quantitative Genetics, 4th ed. Harlow, Essex, UK: Longmans Green.
• 41.Beavis WD (1994) The power and deceit of QTL experiments: Lessons from comparative QTL studies. 49th Annual Corn and Sorghum Industry Research Conference ASTA, Washington, DC. pp. 250–266.
• 42.Beavis WD (1998) QTL analyses: power, precision and accuracy. Molecular Analysis of Complex Traits CRC Press, Boca Raton, FL. pp. 145–161.
• 43.Hines HM, Papa R, Ruiz M, Papanicolaou A, Wang C, et al. (2012) Transcriptome analysis reveals novel patterning and pigmentation genes underlyingHeliconius butterfly wing pattern variation. BMC Genomics13: 288. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 44.Schwanwitsch BN (1924) On the ground plan of the wing-pattern in nymphalids and certain other families of rhopalocerous Lepidoptera. Proc Zool Soc Lond B34: 509–528. [Google Scholar]
• 45.Nijhout HF (2001) Elements of butterfly wing patterns. J Exp Zool291: 213–225. [DOI] [PubMed] [Google Scholar]
• 46.Carroll SB, Gates J, Keys DN, Paddock SW, Panganiban GE, et al. (1994) Pattern formation and eyespot determination in butterfly wings. Science265: 109–114. [DOI] [PubMed] [Google Scholar]
• 47.Beldade P, French V, Brakefield PM (2008) Developmental and genetic mechanisms for evolutionary diversification of serial repeats: eyespot size in Bicyclus anynana butterflies. J Exp Zool B Mol Dev Evol310: 191–201. [DOI] [PubMed] [Google Scholar]
• 48.Brunetti CR, Selegue JE, Monteiro A, French V, Brakefield PM, et al. (2001) The generation and diversification of butterfly eyespot color patterns. Curr Biol11: 1578–1585. [DOI] [PubMed] [Google Scholar]
• 49.Nijhout HF (1980) Pattern formation on lepidopteran wings: determination of an eyespot. Dev Biol80: 267–274. [DOI] [PubMed] [Google Scholar]
• 50.Monteiro A, French V, Smit G, Brakefield PM, Metz JA (2001) Butterfly eyespot patterns: evidence for specification by a morphogen diffusion gradient. Acta Biotheor49: 77–88. [DOI] [PubMed] [Google Scholar]
• 51.Monteiro A, Chen B, Scott LC, Vedder L, Prijs HJ, et al. (2007) The combined effect of two mutations that alter serially homologous color pattern elements on the fore and hindwings of a butterfly. BMC Genet8: 22. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 52.Reed RD, Serfas MS (2004) Butterfly wing pattern evolution is associated with changes in a Notch/Distal-less temporal pattern formation process. Curr Biol14: 1159–1166. [DOI] [PubMed] [Google Scholar]
• 53.Nijhout HF (2010) Molecular and physiological basis of colour pattern formation. Adv In Insect Phys38: 219–265. [Google Scholar]
• 54.Nijhout HF, Wray GA (1986) Homologies in the colour patterns of the genus Charaxes (Lepidoptera: Nymphalidae). Biol J Linn Soc28: 387–410. [Google Scholar]
• 55.McMillan WO, Monteiro A, Kapan DD (2002) Development and evolution on the wing. Trends Ecol Evol17: 125–133. [Google Scholar]
• 56.Baxter SW, Johnston SE, Jiggins CD (2009) Butterfly speciation and the distribution of gene effect sizes fixed during adaptation. Heredity102: 57–65. [DOI] [PubMed] [Google Scholar]
• 57.Jones RT, Salazar PA, ffrench-Constant RH, Jiggins CD, Joron M (2012) Evolution of a mimicry supergene from a multilocus architecture. Proc Biol Sci279: 316–325. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 58.Brown KS, Benson WW (1974) Adaptive polymorphism associated with multiple Müllerian mimicry inHeliconius numata. Biotropica6: 205–228. [Google Scholar]
• 59.Punnett RC (1915) Mimicry in butterflies. London: Cambridge University Press.
• 60.Nicholson AJ (1927) A new theory of mimicry in insects. Austral Zool5: 10–104. [Google Scholar]
• 61.Turner JRG (1981) Adaptation and Evolution inHeliconius: a defense of Neo Darwinism. Annual Review of Ecology and Systematics12: 99–121. [Google Scholar]
• 62.Turner JRG (1977) Butterfly mimicry: the genetical evolution of an adaptation. Evol Biol10: 163–206. [Google Scholar]
• 63.Turner JRG (2005) Mimicry. In: Encyclopedia of life sciences. New York: John Wiley and Sons, Ltd. pp. 1–9.
• 64.Turner JRG (1985) Fisher's evolutionary faith and the challenge of mimicry. Oxford Surveys in Evolutionary Biology2: 159–196. [Google Scholar]
• 65.Turner JRG (1988) The evolution of mimicry: a solution to the problem of punctuated equilibrium. Am Nat131: 42–66. [Google Scholar]
• 66.Orr HA (2005) The genetic theory of adaptation: a brief history. Nat Rev Genet6: 119–127. [DOI] [PubMed] [Google Scholar]
• 67.Balogh AC, Leimar O (2005) Müllerian mimicry: an examination of Fisher's theory of gradual evolutionary change. Proc Biol Sci272: 2269–2275. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 68.Balogh AC, Gamberale-Stille G, Tullberg BS, Leimar O (2010) Feature theory and the two-step hypothesis of Müllerian mimicry evolution. Evolution64: 810–822. [DOI] [PubMed] [Google Scholar]
• 69.Franks DW, Sherratt TN (2007) The evolution of multicomponent mimicry. J Theor Biol244: 631–639. [DOI] [PubMed] [Google Scholar]
• 70.Ruxton GD, Franks DW, Balogh AC, Leimar O (2008) Evolutionary implications of the form of predator generalization for aposematic signals and mimicry in prey. Evolution62: 2913–2921. [DOI] [PubMed] [Google Scholar]
• 71.Mallet J, Barton NH (1989) Strong natural selection in a warning color hybrid zone. Evolution43: 421–431. [DOI] [PubMed] [Google Scholar]
• 72.Benson WW (1972) Natural selection for Müllerian mimicryHeliconius erato in Costa Rica. Science176: 936–939. [DOI] [PubMed] [Google Scholar]
• 73.Orr HA, Coyne JA (1992) The genetics of adaptation: a reassessment. Amer Nat140: 725–742. [DOI] [PubMed] [Google Scholar]
• 74.Bradshaw HD, Schemske DW (2003) Allele substitution at a flower colour locus produces a pollinator shift in monkeyflowers. Nature426: 176–178. [DOI] [PubMed] [Google Scholar]
• 75.Colosimo PF, Hosemann KE, Balabhadra S, Villarreal G Jr, Dickson M, et al. (2005) Widespread parallel evolution in sticklebacks by repeated fixation of Ectodysplasin alleles. Science307: 1928–1933. [DOI] [PubMed] [Google Scholar]
• 76.Albertson RC, Streelman JT, Kocher TD, Yelick PC (2005) Integration and evolution of the cichlid mandible: the molecular basis of alternate feeding strategies. Proc Natl Acad Sci U S A102: 16287–16292. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 77.Stinchcombe JR, Weinig C, Heath KD, Brock MT, Schmitt J (2009) Polymorphic genes of major effect: consequences for variation, selection and evolution inArabidopsis thaliana. Genetics182: 911–922. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 78.Fisher RA (1930) The Genetical Theory of Natural Selection. Oxford: Clarendon Press.
• 79.Stern DL, Orgogozo V (2009) Is genetic evolution predictable?Science323: 746–751. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 80.McMillan WO, Jiggins CD, Mallet J (1997) What initiates speciation in passion-vine butterflies?Proc Natl Acad Sci U S A94: 8628–8633. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 81.Mallet J (1989) The genetics of warning colour in Peruvian hybrid zones ofHeliconius erato andH. melpomene. Proc R Soc Lond B236: 163–185. [Google Scholar]
• 82.Sheppard PM (1963) Some genetic studies on Müllerian mimics in butterflies of the genusHeliconius. Zoologica48: 145–154. [Google Scholar]
• 83.Vos P, Hogers R, Bleeker M, Reijans M, van de Lee T, et al. (1995) AFLP: a new technique for DNA fingerprinting. Nucleic Acids Res23: 4407–4414. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 84.Papa R, Troggio M, Ajmone-Marsan P, Nonnis Marzano F (2005) An improved protocol for the production of AFLP markers in complex genomes by means of capillary electrophoresis. J Anim Breed Genet122: 62–68. [DOI] [PubMed] [Google Scholar]
• 85.Van Ooijen JW, Voorrips RE (2001) JoinMap 3.0: software for the calculation of genetic linkage maps. Wageningen, The Netherlands: Plant Research International.
• 86.Shi J, Heckel DG, Goldsmith MR (1995) A genetic linkage map for the domesticated silkworm,Bombyx mori, based on restriction fragment length polymorphisms. Genet Res66: 109–126. [Google Scholar]
• 87.Lincoln SE, Daly MJ, Lander ES (1987) MAPMAKER 3.0. Boston: Whitehead Institute for Biomedical Research/MIT Center for Genome Research.
• 88.Team RDC (2011) R: A language and environment for statistical computing. Vienna, Austria: R Foundation for Statistical Computing.
• 89.Broman KW, Sen S (2009) A Guide to QTL Mapping with R/qtl. Berlin, Germany: Springer.
• 90.Doerge RW, Churchill GA (1996) Permutation tests for multiple loci affecting a quantitative character. Genetics142: 285–294. [DOI] [PMC free article] [PubMed] [Google Scholar]
• 91.Wang S, Batsen CJ, Zheng Z-b (2011) Windows QTL Cartographer 2.5. Raleigh, NC: Department of Statistics, North Carolina State University.
• 92.Wright S (1978) Evolution and the Genetics of Populations. Vol. 4: Variability Within and Among Natural Populations. Chicago: University of Chicago Press.

Supplementary Materials