Update! Chris Jiggins and Nicola Nadeau have been kind enought to reply with a written response.

Nadeau, et al. address genomic divergence in hybrid zones between subspecies of two mimetic butterfly species in Peru and Ecuador. H. melpomene and H. erato diverged 15-20 million years ago, are presumably reproductively isolated, and today live in sympatry across South America, sharing a host plant. Sympatric subspecies pairs share wing coloration patterns, and the members of each pair were nearly indistinguishable to our group of non-lepidopterists. Previous results suggested that H. erato diversified first, with each subspecies presenting aposematic signals to local predators (birds), and that H. melpomene entered these niches later, facilitating the success of its congeneric as Müllerian mimics, reinforcing the warning coloration pattern provided to their shared predators.

To address four aims - 1. Can you use natural hybrid zones to map traits? 2. How much of the genome is differentiated between subspecies? 3. How much of that differentiation is due to coloration? 4. Are the same regions different between co-mimetic species? - the authors use paired end RAD-seq, map to the H. melpomene reference, and de novo assemble using Stacks (although throwing out the second read and all but 100 bp of the first read!). We felt that the comparison between the de novo and mapping approaches was a little unfair, since a majority of the data was not included in the de novo assembly, making it extremely difficult to establish orthology. The authors discuss that others have first assembled paired-end data into contigs, and mapped to those to call SNPs, and we didn’t see the justification for not using these methods. H. erato shows 15% divergence from H. melpomene, so we feel like a reference-free approach that uses all the data would help, especially in divergent regions of the genome.

Quite a bit is known about the genetic basis of coloration, and most of these major effect loci were recovered using association mapping in most of the hybrid zones. Not surprisingly due to the different methodologies, but still concerning, no SNPs associated with a phenotype overlap between the mapped and the de novo datasets, although similar regions were found. They also used Fst outliers, using only ‘pure’ subspecies as defined by their coloration morphology. They find outlier SNPs that are not associated with coloration loci, and suggest these may be experiencing divergence to the different ecologies found in high and low altitude. This seemed like an untestable claim, since the ‘pure’ color morphs were collected at different altitudes. We were surprised that there is very little differentiation in ‘pure’ subspecies at the edges of the hybrid zone, and it appears that gene flow persists through a region between both morphologically divergent subspecies. We wondered as to the extent of gene flow between these subspecies at the edges of the range, such as where these subspecies co-occur in Brazil.

The language describing the phylogenetic analysis makes it seem that the two species are not included in the same phylogeny, so even barring issues of missing data in H. erato, it seems peculiar to say that the two species have diverged at similar times because they have similar branch lengths. Substitution rates could be very different during the 15-20 million years along lineages leading to these radiations. The supplement is not yet available, so perhaps these questions are addressed there. As written in the main text, claims of the order of divergence of the two species seem poorly supported. We were unconvinced that the order of origin of Müllerian mimicry is conclusively known.

We wondered whether small sample sizes made it difficult to recover truly parallel loci in matched subspecies pairs. The Ro locus, which controls the shape of the distal edge of the forewing band, is found through association mapping and is an Fst outlier in Ecuadorian H. erato. Although the locus is not a significant GWAS result for the trait in Ecuadorian H. melpomone, it looked to us to be a Fst outlier. Lack of significance in Ecuadorian H. melpomone might be due to the small sample size necessitated by the removal of seven individuals belonging to a newly defined population of H. timareta that is co-mimetic with the studied high altitude Ecuadorian Heliconius. QTL results from a previous study in H. melpomone suggested a broad region, now known to be on chromosome 13 that encompasses the Ro locus, further emphasizing the potentially parallel nature of wing coloration evolution in this subspecies pair. With small sample sizes and low statistical power, it seems objectionable to exclude the possibility that a given locus is important in both pairs.

Overall, we felt the study was primarily designed to identify large-effect loci associated with divergence in coloration, and as such fell a bit short on some of the larger claims about addressing all divergence in these populations. This, to us, seemed due to small sample sizes (from neat things like finding a previously undetected co-mimetic species!) and common issues with genomics in a non-model system, such as choosing appropriate thresholds and methods for filtering missing data. The study presents potential evidence for the evolution of aposematic signals alongside ecological specialization and local adaptation, but we feel that future sampling to extricate the signal of wing coloration from altitude and ecology would strengthen these claims.