Skip to main content

Advertisement

We’d like to understand how you use our websites in order to improve them. Register your interest.
Download PDF

Comparative transcriptomic analysis of a wing-dimorphic stonefly reveals candidate wing loss genes

EvoDevo volume 10, Article number: 21 (2019) Cite this article

Abstract

Background

The genetic basis of wing development has been well characterised for model insect species, but remains poorly understood in phylogenetically divergent, non-model taxa. Wing-polymorphic insect species potentially provide ideal systems for unravelling the genetic basis of secondary wing reduction. Stoneflies (Plecoptera) represent an anciently derived insect assemblage for which the genetic basis of wing polymorphism remains unclear. We undertake quantitative RNA-seq of sympatric full-winged versus vestigial-winged nymphs of a widespread wing-dimorphic New Zealand stonefly, Zelandoperla fenestrata, to identify genes potentially involved in wing development and secondary wing loss.

Results

Our analysis reveals substantial differential expression of wing-development genes between full-winged versus vestigial-winged stonefly ecotypes. Specifically, of 23 clusters showing significant similarity to Drosophila wing development-related genes and their pea aphid orthologues, nine were significantly upregulated in full-winged stonefly ecotypes, whereas only one cluster (teashirt) was substantially upregulated in the vestigial-winged ecotype.

Conclusions

These findings suggest remarkable conservation of key wing-development pathways throughout 400 Ma of insect evolution. The finding that two Juvenile Hormone pathway clusters were significantly upregulated in vestigial-winged Zelandoperla supports the hypothesis that Juvenile Hormone may play a key role in modulating insect wing polymorphism, as has previously been suggested for other insect lineages.

Background

The unparalleled diversity of insects is often attributed to the evolution of insect flight, ca. 400 Ma [1,2,3,4,5]. Specifically, wing evolution has carried numerous advantages, including an increased ability for insects to access novel resources and ecosystems, improved predator avoidance, and enhanced mate location [1]. This dispersal capacity has, however, been subsequently lost in numerous insect lineages, across almost all winged orders [6, 7].

The development and maintenance of wings requires significant energy expenditure, and some wing-dimorphic insect species thus exhibit apparent ‘trade-offs’ between flight-ability versus reproductive output [8,9,10,11,12]. Environmental conditions may also contribute to the evolution of wing-reduced lineages [6, 13], with flight loss particularly common in stable ecosystems that lack predators [14, 15]. Wing-reduced lineages are also disproportionally abundant in isolated areas, where there is a high level of mortality in dispersing individuals, and in areas where the energetic costs of flight are high [1, 16]. Likewise, wing loss is particularly common at higher altitudes [17], with recent studies implying that exposure, modulated by the alpine treeline, may be a key driver of insect wing reduction [18, 19].

Wings have been completely lost in many insect species, whereas some species comprise distinct winged and wingless ecotypes, often associated with differing environmental conditions [20,21,22]. As the distinct phenotypes of wing-dimorphic species have very similar genetic backgrounds, such taxa present ideal systems for exploring the molecular basis of wing loss. Previous studies of wing-dimorphic species have suggested that wing polymorphism can be genetically determined [23,24,25,26,27], controlled by environmentally driven gene expression (polyphenism; [11, 28,29,30]), include both genetic and environmental components [31], or be under the control of epigenetic factors [32, 33].

Given the range of environmental, genetic, and epigenetic factors that may contribute to wing polymorphism, elucidating the developmental basis of wing loss can be challenging. While the gene networks underpinning wing development have been well-characterised in model organisms such as Drosophila melanogaster [34, 35] and Tribolium castaneum [36], very little is known about the underlying molecular biology of wing development in non-model insect taxa. Several recent studies, however, have identified the key genes underlying wing development in non-model taxa by comparing gene-expression across winged and wingless morphs of wing polymorphic species (see Additional file 1: Table S1). While some portions of the key insect gene networks appear to be conserved among a variety of holometabolous and hemimetabolous taxa [37,38,39,40], the extent of their phylogenetic conservation across class Insecta remains unclear.

The New Zealand stonefly Zelandoperla fenestrata species group (hereafter simply referred to as Z. fenestrata) contains both full-winged and vestigial-winged ecotypes, with the wing-reduced phenotype particularly common at higher altitudes [20, 41]. The distinct morphotypes were traditionally considered distinct species (full-winged: Zelandoperla fenestrata or Zelandoperla tillyardi, vestigial-winged: Zelandoperla pennulata) [41], but more recent molecular studies demonstrate that vestigial-winged lineages have evolved independently on numerous occasions [20, 21]. In addition, a recent study identified several outlier loci strongly associated with vestigial versus full-winged ecotypes of Z. fenestrata, strongly suggesting a genetic basis for wing-length variation in this species [42].

In this study, we use quantitative RNA-seq to compare gene-expression patterns between the sympatric nymphs of full-winged versus vestigial-winged Z. fenestrata, to identify genes involved in wing development, wing growth and wing differentiation in this dimorphic species. Additionally, we assessed expression differences in candidate genes involved in the juvenile hormone (JH) pathway, as these genes have previously been demonstrated to play a role in insect wing polymorphisms [11, 28]. We identified genes that were significantly differentially expressed between vestigial-winged and full-winged ecotypes, with a particular focus on genes known to be involved in the wing-development networks.

Results

Transcriptome sequencing, assembly, and annotation

Sequencing yielded 204,011,700 125-bp paired reads. After initial adapter trimming and quality filter, 147,396,767 reads remained (Additional file 1: Table S2). These raw RNA-seq reads have been deposited in the NCBI repository under BioProject PRJNA525904. De novo assembly using Trinity generated 552,851 transcripts hierarchically clustered into 442,924 Trinity ‘genes’, with a mean length of 587.6 bp and an N50 of 887 bp (Additional file 1: Table S2). 98.6% of BUSCO proteins were identified as full-length proteins, with only 1.0% fragmented and 0.4% missing (Additional file 1: Table S2). Over 84% of input reads aligned as proper pairs one or more times, and the overall alignment was over 97% (Additional file 1: Table S2). TransDecoder predicted open reading frames in 108,209 transcripts and BLASTp annotated over 70% of these ORFs (UniProt/Swiss-Prot database) (Additional file 1: Table S3).

Corset clustering reduced the number of transcripts to 227,791, which were assigned to 140,592 clusters with a mean transcript length of 1023.5 bp and an N50 of 1668 bp (Additional file 1: Table S2). 98.5% of BUSCOs were found to be complete with 0.7% missing and 0.8% fragmented (Additional file 1: Table S2). 84% of input reads aligned concordantly one or more times, and the overall alignment was over 96% (Additional file 1: Table S2).

Genes differentially expressed across the full-winged and vestigial-winged ecotypes

DE analysis identified 1042 clusters differentially expressed at FDR ≤ 0.01 and LogFC ≥ 1.5 or ≤ − 1.5. Of the 1042 DE clusters, 511 contained TransDecoder-predicted longest ORFs which were annotated by emapper. Of these, 332 clusters were more abundant in full-winged individuals, whilst 179 clusters were more highly expressed in vestigial-winged notum (Fig. 2). Forty-three clusters were annotated with GO:0035220, biological process ‘wing disc development’; 17 of which were more highly expressed in full-winged notum, whilst 26 were more abundant in vestigial-winged notum. All clusters with FDR ≤ 0.01 and LogFC ≥ 1.5 or ≤ − 1.5, including annotations, are provided in Additional file 2: Table S4.

Identification and expression of wing development genes

Twenty-three clusters with significant similarity to Drosophila wing development-related genes and their pea aphid orthologues were identified, including genes in each of the four wing development pathways, as well as homeobox genes (Table 1). As a result of low numbers of reads aligning to the transcripts of one of these clusters, statistics and expression data were available for only 22 clusters. Just over half of these clusters were not significantly differentially expressed across ecotypes (Fig. 3a). However, nine clusters were significantly more expressed in the notum of full-winged individuals (at FDR ≤ 0.01, LogFC ≥ 1.5), including wingless, nubbin, distalless, hedgehog, ventral veins lacking, and engrailed (Fig. 3a). By contrast only a single cluster (teashirt) was upregulated in the vestigial-winged ecotype (Fig. 3a). DE analysis also revealed two clusters with known roles in wing development; optomotor-blind protein and a frizzled receptor family member. Both were more highly expressed in full-winged individuals (Fig. 3b).

Table 1 Wing development genes predicted in the notum of vestigial- and full-winged Zelandoperla fenestrata
Full size table

In addition to wing development genes (above), a juvenile hormone acid methyltransferase cluster was identified as differentially expressed between the two ecotypes, being significantly upregulated in the vestigial-winged samples (Fig. 3b). A juvenile hormone esterase cluster was also significantly differentially expressed between vestigial-winged and full-winged ecotypes, being more abundant in vestigial-winged ecotypes (Fig. 3b).

Discussion

In this study, we compared the expression of genes in the notum of full-winged versus vestigial-winged ecotypes of the wing-dimorphic stonefly Zelandoperla fenestrata. The former specimens had prominent wing-bud development, whereas the latter did not (Fig. 1). We identified over 20 clusters with significant similarity to genes previously implicated in wing development in Drosophila melanogaster. Over 1000 clusters were significantly differentially expressed among the two ecotypes, including numerous clusters annotated as genes involved in key wing development pathways (Figs. 2, 3). We discuss the implications of these findings, particularly in reference to identifying the key gene/s potentially involved in wing loss, below.

Fig. 1
figure1

Dorsal view of Z. fenestrata nymphs sequenced in this study illustrating the extent of wing-bud development in a the vestigial-winged ecotype, and b the full-winged ecotype. P pronotum, M1 mesonotum, M2 metanotum, W1/W2 wing-buds

Full size image
Fig. 2
figure2

Volcano plot displaying differentially expressed genes in the notum of full-winged and vestigial-winged Zelandoperla fenestrata. The x-axis displays the log fold-expression change, with negative values indicating increased expression in the vestigial-winged ecotype, and positive values indicating increased expression in the full-winged ecotype. The y-axis displays the −log10 adjusted p value. Genes that were significantly more expressed (FDR ≤ 0.01) are coloured (full-winged ecotype = cyan, vestigial-winged ecotype = magenta)

Full size image
Fig. 3
figure3

Comparative expression of a wing development genes and b additional genes of interest in the notum of full-winged versus vestigial-winged ecotypes. Genes in bold are significantly differentially expressed across ecotypes. Positive LogFC values indicate the gene is more expressed in full-winged ecotypes, with the intensity of the cyan shading indicating the extent of increased expression. Negative LogFC values indicate the gene is more expressed in vestigial-winged ecotypes, with the intensity of magenta shading indicating the extent of increased expression. Shading in columns to the right indicate the expression of each gene in individual samples (see key)

Full size image

Some of the clusters with differential expression between full-winged and vestigial-winged forms of Z. fenestrata stoneflies are also critical components of wing patterning in D. melanogaster. While it is unclear how stonefly wing development differs from Drosophila, the key genes appear conserved, and the common ancestry of these insect lineages implies that some aspects of the patterning system are likely conserved, especially those involved in initial patterning. This finding suggests a remarkable conservation of these key developmental pathways, as the ancestors of Plecoptera and Diptera diverged ca. 400 Ma [4].

About half of these candidate wing-development genes were not significantly differentially expressed across the two ecotypes (Fig. 3a). This result is unsurprising, as many of these genes are highly pleiotropic, having additional important developmental roles [40]. Several key wing-development genes were, however, significantly differentially expressed across vestigial-winged and full-winged ecotypes (Fig. 3a). The prominent wing-buds in the full-winged specimens indicate that wing development is already well underway; it, therefore, seems likely that the lower expression level of these wing-development genes in the vestigial-winged ecotype may represent the down-stream consequences of wing reduction, rather than representing the underlying causes of wing reduction.

In Drosophila, three sets of genes act to define the anterior–posterior, the dorsoventral and the proximo-distal axes (reviewed by [43, 44]). The anterior–posterior axis is defined by the expression of engrailed [45]. Engrailed causes the expression of hedgehog [46], which diffuses across the anterior–posterior boundary of the wing, to activate decapentaplegic [47, 48]. Expression of engrailed and hedgehog is reduced in vestigial-winged compared to full-winged stoneflies, as is the expression of optomotor-blind (a gene acting downstream of decapentaplegic), confirming these genes play crucial roles in stonefly wing development.

The dorsoventral axis of the Drosophila wing is patterned by wnt signalling molecules, especially wingless (wg1) [49,50,51,52]. Wingless and three other wg molecules (wg 6, 7, and 16) have reduced expression in vestigial-winged stoneflies, as do their receptors, frizzled. A downstream patterning factor responding to wingless, ventral-veinless (also known as drifter) [53] is also significantly lower in expression in vestigial-winged stoneflies.

Proximo-distal growth of the Drosophila wing is regulated by two key genes, teashirt and nubbin [54]. Nubbin is down-regulated in vestigial-winged stoneflies, while teashirt is up-regulated. Teashirt is known to be highly expressed in the notum and wing-hinge of Drosophila [55], so increased expression was not unexpected in the vestigial-winged nymphs given the slightly higher proportion of notum tissue in the vestigial-winged extractions (see Fig. 1). However, this slightly higher proportion of notum tissue is unlikely to explain the fourfold change in expression we observed (Fig. 3a). While it is not clear what the consequences for wing outgrowth of these changes in gene expression are, they imply mechanisms of proximo-distal axis growth of the wing differ between full-winged and vestigial-winged forms. In Drosophila, repression of teashirt is required for the initiation of wing development [55]. That this gene is upregulated in vestigial-winged stoneflies implies that down-regulation may not occur, potentially providing a proximate mechanism (rather than the ultimate cause) of wing reduction.

Differences in gene expression between full-winged versus vestigial-winged forms of stonefly encompass mechanisms that define and pattern all three wing axes, implying that all aspects of wing growth and patterning may be reduced in vestigial-winged stoneflies. Given the relatively rapid evolution of this trait (see [19]), it seems very unlikely that simultaneous changes across all of these networks underlie the causative changes in the evolution of vestigial-winged forms. We, therefore, suggest that the causative changes may lie upstream of wing patterning and growth, perhaps in juvenile hormone (JH) signalling. JHs regulate several key processes in insect physiology [56], and have previously been associated with wing polymorphism in both Orthoptera [28, 57] and Hemiptera [58, 59].

Two JH pathway clusters were significantly upregulated in the vestigial-winged ecotype, suggesting JH may play a role in modulating wing polymorphism in Z. fenestrata. JH levels in insects can vary both within and across developmental instars, with JH levels typically lower in later instars, but with levels also dropping significantly directly before moulting [56]. Populations of Z. fenestrata may comprise individuals at a variety of developmental stages, as the species has an extended emergence ‘window’ (see [60]). Although the wild-caught nymphs used in this study were all late instars, we were unable to tightly control for their developmental stage, and as a result, some variation in the expression of key JH-related genes might be expected among individuals (e.g. within ecotypes). For instance, the finding that JH acid methyltransferase expression varied substantially among vestigial-winged specimens (Fig. 3b) might potentially reflect minor developmental differences among these specimens. However, there is no evidence for systematic developmental asynchrony among these ecotypes where they co-occur, and indeed there are several localities at which we have simultaneously collected both full-winged and vestigial-winged adults [19, 42]. Systematic differences among ecotypes in the expression of JH-related genes are thus unlikely to result solely from (minor) differences in developmental stage. For instance, JH acid methyltransferase expression is consistently low across all full-winged specimens (despite likely developmental variation within this ecotype, see Fig. 3b), but variable among vestigial-winged individuals. Similarly, JH esterase expression is consistently high across all vestigial-winged specimens, but is highly variable in full-winged individuals. In both cases, the extensive differentiation among ecotypes is unlikely to reflect developmental timing, suggesting that the repression of these genes may play an important role in Z. fenestrata wing development.

In addition to the extensive evidence for effects of JHs on insect wing development, these genes have previously been associated with delayed metamorphosis [61], and variation in several other insect life-history traits (e.g. fecundity; [62]). It is thus conceivable that JH could simultaneously modulate both wing phenotype variation and additional developmental differentiation between high-altitude and low-altitude stonefly lineages (such as differences in body size and emergence timing, see [60, 63]). Future molecular and ecological studies will aim to unravel the potential role of JH in explaining ecotypic variation within and among wing-polymorphic stonefly taxa.

Conclusion

Our study reveals differential expression of wing-development genes between closely related full-winged and vestigial-winged stonefly ecotypes, consistent with previous data from Drosophila, suggesting remarkable conservation of key wing development pathways throughout 400 Ma of insect evolution. The significant gene expression differentiation observed among stonefly ecotypes potentially provides fertile new directions for future evolutionary research on insect wing reduction.

Methods

Sample collection

Late instar Z. fenestrata nymphs were collected by hand from beneath rocks in the stream at an altitude of 670 m in Lug Creek, Rock and Pillar Range, Otago, New Zealand, in July 2017. This site is just above the alpine tree line, and is part of a narrow zone where full-winged and vestigial-winged ecotypes are found in sympatry [19]. Nymphs were stored in containers in water from their natal stream, and kept cool in an ice bath while they were returned to the laboratory. Nymphs were characterised as either full-winged or vestigial-winged based on clear differences in wing-bud development (Fig. 1), as previous work has indicated that these two wing-bud development types are correlated with the distinct adult ecotypes [19, 41, 42].

RNA extraction

Ten nymphs from each of the vestigial-winged and full-winged ecotypes were extracted and sequenced. As sufficient quantities of RNA could not readily be obtained from the wing-buds (alone) of vestigial winged nymphs (as their wing-buds are extremely small, see Fig. 1), we instead extracted RNA from the complete notum (including wing-buds) of both ecotypes. This consistent approach allowed us to directly compare gene expression across full-winged and vestigial-winged ecotypes, although we note the ratio of wing-bud tissue to notum tissue differed across vestigial-winged and full-winged specimens (see Fig. 1). Nymphs were immobilised briefly on dry ice, then pinned through the head and abdomen under a binocular dissection scope before the mesonotum (M1; Fig. 1) and metanotum (M2; Fig. 1) were removed as a single sample. A shallow incision using a razorblade was made across the dorsal surface between the protonotum and M1. A lateral incision was then made on each side separating M1 and M2 from the meso- and metapleuron. A final incision between M2 and the first tergal segment was used to free M1 and M2. Care was taken to remove any fat bodies or viscera that were attached centrally to M1 or M2 to maintain consistency between samples. Once dissected, the notum samples were immediately frozen on dry ice and stored at − 80 °C.

Tissue was homogenised with hard plastic homogenising probes in RTL buffer (supplied as part of the Qiagen RNeasy kit) using an Omni Tissue Homogenizer. RNA isolation was then completed following the protocol supplied with the RNeasy kit (Qiagen), including an on-column DNAse digestion. Once extracted and washed, samples were eluted with 35 μL of RNAse-free water and then stored at − 80 °C.

Library preparation

The 260/230 and 260/280 ratios were examined using a NanoDrop spectrophotometer (Thermo Scientific), and the concentrations were assessed using a Qubit fluorometer (Thermo-Fisher Scientific) with the RNA HS kit. RNA integrity was also assessed on a 1.5% agarose gel after electrophoresis. We used a poly-A capture to select for mRNA, and created a library using a TruSeq Stranded mRNA Library Prep Kit (Illumina), with each individual tagged to separate reads bioinformatically after sequencing. Library quality was assessed using the Agilent 2100 Bioanalyzer and TBS 380 Fluorometer (Turner Biosystems, Sunny-vale, CA, USA), and all libraries were paired-end sequenced (2 × 125 bp) on a single lane of an Illumina HiSeq 2500.

Transcriptome assembly

The combined vestigial-winged- and full-winged-notum transcriptome was de novo assembled using Trinity v2.8.4 [64, 65]. Reads were quality trimmed using Trimmomatic, run as part of the Trinity package with default settings [66, 67], and normalised in silico prior to assembly. Descriptive statistics were calculated using TrinityStats.pl. Transcriptome completeness was assessed using BUSCO v3.0.2 (an arthropod_odb9 database of 1066 BUSCOs covering 60 species; [68, 69]) and RNA-seq read representation following the trinityrnaseq-wiki (github.com/trinityrnaseq/trinityrnaseq/wiki). Transcript clusters, clustered using Corset (see below), were quality checked in the same way.

The Trinity transcriptome was annotated using Trinotate v3.1.1 [70] and associated packages (TransDecoder v5.5.0, SQLite, blast + v, hmmer v3.2.1, Rnammer v1.2, signalP v4.1, tmhmm v2.0, trinotateR). Annotation of differentially expressed clusters (see below) used emapper-1.0.3 [71] based on eggNOG orthology data [72] with sequence searches performed using HMM and diamond [73, 74]. ‘Contigs of interest’ (i.e. differentially expressed clusters of contigs) were extracted from the Trinity transcriptome fasta file using ‘fetchClusterSeqs.py’ script (github.com/Adamtaranto/Corset-tools) and protein sequences predicted using TransDecoder v 5.5.0. Annotation of the longest ORFs was performed using emapper and two methods: HMMER, using the arthropod protein database artNOG, and diamond, using the eggNOG protein database.

Identification of differentially expressed transcripts

Differential gene expression between vestigial-winged and full-winged individuals was assessed using Corset v1.07 [75] and with read alignment using bowtie2 (v2.3.4.1) with multi-mapping enabled. Statistical analysis used edgeR in R i386 3.5.1 [76]. Analysis followed the Corset wiki (github.com/Oshlack/Corset/wiki).

Genes of interest

In addition to identifying differentially expressed genes, expression data for genes with known roles in wing development were sought. Local blast searches were done using Drosophila wing development-related genes and orthologues from the pea aphid, Acyrthosiphon pisum [40], as query terms (taken from the NCBI website: http://www.ncbi.nlm.nih.gov), and using BioEdit software [77]. Searches were done against clusters (Corset clustered transcripts).

Availability of data and materials

The dataset supporting the conclusions of this article is available in the NCBI repository under BioProject PRJNA525904.

References

  1. 1.

    Wagner DL, Liebherr JK. Flightlessness in insects. Trends Ecol Evol. 1992;7:216–20.

  2. 2.

    Kingsolver JG, Koehl M. Selective factors in the evolution of insect wings. Annu Rev Entomol. 1994;39:425–51.

  3. 3.

    Daly HV, Doyen JT, Ehrlich PR. Introduction to insect biology and diversity: McGraw-Hill Book Company.; 1978.

  4. 4.

    Misof B, Liu S, Meusemann K, Peters RS, Donath A, Mayer C, et al. Phylogenomics resolves the timing and pattern of insect evolution. Science. 2014;346:763–7.

  5. 5.

    Medved V, Marden JH, Fescemyer HW, Der JP, Liu J, Mahfooz N, et al. Origin and diversification of wings: insights from a neopteran insect. Proc Natl Acad Sci USA. 2015;112:15946–51.

  6. 6.

    Roff DA. The evolution of flightlessness: is history important? Evol Ecol. 1994;8:639–57.

  7. 7.

    Trautwein MD, Wiegmann BM, Beutel R, Kjer KM, Yeates DK. Advances in insect phylogeny at the dawn of the postgenomic era. Annu Rev Entomol. 2012;57:449–68.

  8. 8.

    Roff DA. The evolution of wing dimorphism in insects. Evolution. 1986;40:1009–20.

  9. 9.

    Roff DA, Bradford MJ. Quantitative genetics of the trade-off between fecundity and wing dimorphism in the cricket Allonemobius socius. Heredity. 1996;76:178–85.

  10. 10.

    Roff DA, Fairbairn DJ. Wing dimorphisms and the evolution of migratory polymorphisms among the Insecta. Am Zool. 1991;31:243–51.

  11. 11.

    Zera AJ, Denno RF. Physiology and ecology of dispersal polymorphism in insects. Annu Rev Entomol. 1997;42:207–30.

  12. 12.

    Vogler AP, Timmermans MJ. Speciation: don’t fly and diversify? Curr Biol. 2012;22:R284–6.

  13. 13.

    Harrison RG. Dispersal polymorphisms in insects. Annu Rev Ecol Syst. 1980;11:95–118.

  14. 14.

    Roff DA. Habitat persistence and the evolution of wing dimorphism in insects. Am Nat. 1994;144:772–98.

  15. 15.

    Denno RF, Roderick GK, Peterson MA, Huberty AF, Dobel HG, Eubanks MD, et al. Habitat persistence underlies intraspecific variation in the dispersal strategies of planthoppers. Ecol Monogr. 1996;66:389–408.

  16. 16.

    Roff DA. The evolution of flightlessness in insects. Ecol Monogr. 1990;60:389–421.

  17. 17.

    Hodkinson ID. Terrestrial insects along elevation gradients: species and community responses to altitude. Biol Rev. 2005;80:489–513.

  18. 18.

    McCulloch GA, Foster BJ, Ingram T, Waters JM. Insect wing loss is tightly linked to the treeline: evidence from a diverse stonefly assemblage. Ecography. 2019;42:811–3.

  19. 19.

    McCulloch GA, Foster BJ, Dutoit L, Ingram T, Hay E, Veale AJ, et al. Ecological gradients drive insect wing loss and speciation: the role of the alpine treeline. Mol Ecol. 2019;28(13):3141–50.

  20. 20.

    McCulloch GA, Wallis GP, Waters JM. Do insects lose flight before they lose their wings? Population genetic structure in subalpine stoneflies. Mol Ecol. 2009;18:4073–87.

  21. 21.

    Dussex N, Chuah A, Waters JM. Genome-wide SNPs reveal fine-scale differentiation among wingless alpine stonefly populations and introgression between winged and wingless forms. Evolution. 2016;70:38–47.

  22. 22.

    Van Belleghem SM, Roelofs D, Hendrickx F. Evolutionary history of a dispersal-associated locus across sympatric and allopatric divergent populations of a wing-polymorphic beetle across Atlantic Europe. Mol Ecol. 2015;24:890–908.

  23. 23.

    Aukema B. Wing-length determination in two wing-dimorphic Calathus species (Coleoptera: Carabidae). Hereditas. 1990;113:189–202.

  24. 24.

    Jackson DJ. XXVII.—The inheritance of long and short wings in the weevil, Sitona hispidula, with a discussion of wing reduction among beetles. Earth Environ Sci. 1928;55:665–735.

  25. 25.

    Lindroth CH. Inheritance of wing dimorphism in Pterostichus anthracinus Ill. Hereditas. 1946;32:37–40.

  26. 26.

    Harrison R. Flight polymorphism in the field cricket Gryllus pennsylvanicus. Oecologia. 1979;40:125–32.

  27. 27.

    Rose D. Dispersal and quality in populations of Cicadulina species (Cicadellidae). J Anim Ecol. 1972;1:589–609.

  28. 28.

    Zera AJ. The endocrine regulation of wing polymorphism in insects: state of the art, recent surprises, and future directions. Integr Comp Biol. 2003;43:607–16.

  29. 29.

    Lin X, Yao Y, Wang B, Emlen DJ, Lavine LC. Ecological trade-offs between migration and reproduction are mediated by the nutrition-sensitive insulin-signaling pathway. Int J Biol Sci. 2016;12:607–16.

  30. 30.

    Xu H-J, Xue J, Lu B, Zhang X-C, Zhuo J-C, He S-F, et al. Two insulin receptors determine alternative wing morphs in planthoppers. Nature. 2015;519:464–7.

  31. 31.

    Honěk A. Factors and consequences of a non-functional alary polymorphism in Pyrrhocoris apterus (Heteroptera: Pyrrhocoridae). Res Popul Ecol. 1995;37:111–8.

  32. 32.

    Zhou X, Chen J, Shike Liang M, Wang F. Differential DNA methylation between two wing phenotypes adults of Sogatella furcifera. Genesis. 2013;51:819–26.

  33. 33.

    Liang S-K, Liang Z-Q, Zhou X-S, Chen J-L, Li G-H, Wang F-H. CpG methylated ribosomal RNA genes in relation to wing polymorphism in the rice pest Sogatella furcifera. J Asia-Pac Entomol. 2015;18:471–5.

  34. 34.

    Williams JA, Carroll SB. The origin, patterning and evolution of insect appendages. BioEssays. 1993;15:567–77.

  35. 35.

    Azpiazu N, Morata G. Function and regulation of homothorax in the wing imaginal disc of Drosophila. Development. 2000;127:2685–93.

  36. 36.

    Consortium TGS. The genome of the model beetle and pest Tribolium castaneum. Nature. 2008;452:949–55.

  37. 37.

    Weatherbee SD, Nijhout HF, Grunert LW, Halder G, Galant R, Selegue J, et al. Ultrabithorax function in butterfly wings and the evolution of insect wing patterns. Curr Biol. 1999;9:109–15.

  38. 38.

    Abouheif E, Wray GA. Evolution of the gene network underlying wing polyphenism in ants. Science. 2002;297:249–52.

  39. 39.

    Van Belleghem SM, Roelofs D, Van Houdt J, Hendrickx F. De novo transcriptome assembly and SNP discovery in the wing polymorphic salt marsh beetle Pogonus chalceus (Coleoptera, Carabidae). PLoS ONE. 2012;7:e42605.

  40. 40.

    Brisson JA, Ishikawa A, Miura T. Wing development genes of the pea aphid and differential gene expression between winged and unwinged morphs. Insect Mol Biol. 2010;19:63–73.

  41. 41.

    McLellan ID. A revision of Zelandoperla Tillyard (Plecoptera: Gripopterygidae: Zelandoperlinae). New Zeal J Zool. 1999;26:199–219.

  42. 42.

    Veale AJ, Foster BJ, Dearden PK, Waters JM. Genotyping-by-sequencing supports a genetic basis for alpine wing-reduction in a New Zealand stonefly. Sci Rep. 2018;8:16275.

  43. 43.

    de Celis JF. Pattern formation in the Drosophila wing: the development of the veins. BioEssays. 2003;25:443–51.

  44. 44.

    Klein T. Wing disc development in the fly: the early stages. Curr Opin Genet Dev. 2001;11:470–5.

  45. 45.

    Bate M, Arias AM. The embryonic origin of imaginal discs in Drosophila. Development. 1991;112:755–61.

  46. 46.

    Guillén I, Mullor JL, Capdevila J, Sánchez-Herrero E, Morata G, Guerrero I. The function of engrailed and the specification of Drosophila wing pattern. Development. 1995;121:3447–56.

  47. 47.

    Tanimoto H, Itoh S, ten Dijke P, Tabata T. Hedgehog creates a gradient of DPP activity in Drosophila wing imaginal discs. Mol Cell. 2000;5:59–71.

  48. 48.

    Zecca M, Basler K, Struhl G. Sequential organizing activities of engrailed, hedgehog and decapentaplegic in the Drosophila wing. Development. 1995;121:2265–78.

  49. 49.

    Cadigan KM, Fish MP, Rulifson EJ, Nusse R. Wingless repression of Drosophila frizzled 2 expression shapes the Wingless morphogen gradient in the wing. Cell. 1998;93:767–77.

  50. 50.

    Couso JP, Bate M, Martinez-Arias A. A wingless-dependent polar coordinate system in Drosophila imaginal discs. Science. 1993;259:484–9.

  51. 51.

    Neumann CJ, Cohen SM. Long-range action of Wingless organizes the dorsal-ventral axis of the Drosophila wing. Development. 1997;124:871–80.

  52. 52.

    Ng M, Diaz-Benjumea FJ, Vincent J-P, Wu J, Cohen SM. Specification of the wing by localized expression of wingless protein. Nature. 1996;381:316–8.

  53. 53.

    Certel K, Hudson A, Carroll SB, Johnson WA. Restricted patterning of vestigial expression in Drosophila wing imaginal discs requires synergistic activation by both Mad and the drifter POU domain transcription factor. Development. 2000;127:3173–83.

  54. 54.

    Zirin JD, Mann RS. Nubbin and Teashirt mark barriers to clonal growth along the proximal–distal axis of the Drosophila wing. Dev Biol. 2007;304:745–58.

  55. 55.

    Wu J, Cohen SM. Repression of Teashirt marks the initiation of wing development. Development. 2002;129:2411–8.

  56. 56.

    Nijhout H. Insect hormones. Princeton: Princeton University Press; 1994. p. 267.

  57. 57.

    Zera AJ. Juvenile Hormone and the endocrine regulation of wing polymorphism in insects: new insights from circadian and functional-genomic studies in Gryllus crickets. Physiol Entomol. 2016;41:313–26.

  58. 58.

    Braendle C, Davis GK, Brisson JA, Stern DL. Wing dimorphism in aphids. Heredity. 2006;97:192–9.

  59. 59.

    Ishikawa A, Gotoh H, Abe T, Miura T. Juvenile hormone titer and wing-morph differentiation in the vetch aphid Megoura crassicauda. J Insect Physiol. 2013;59:444–9.

  60. 60.

    McCulloch GA, Waters JM. Testing for seasonality in alpine streams: how does altitude affect freshwater insect life cycles? Freshw Biol. 2018;63:483–91.

  61. 61.

    Riddiford LM. Juvenile hormone-induced delay of metamorphosis of the viscera of the cecropia silkworm. Biol Bull. 1975;148:429–39.

  62. 62.

    Flatt T, Kawecki TJ. Juvenile hormone as a regulator of the trade-off between reproduction and life span in Drosophila melanogaster. Evolution. 2007;61:1980–91.

  63. 63.

    McCulloch GA, Waters JM. Does wing reduction influence the relationship between altitude and insect body size? A case study using New Zealand’s diverse stonefly fauna. Ecol Evol. 2018;8:953–60.

  64. 64.

    Grabherr MG, Haas BJ, Yassour M, Levin JZ, Thompson DA, Amit I, et al. Trinity: reconstructing a full-length transcriptome without a genome from RNA-Seq data. Nat Biotechnol. 2011;29:644–52.

  65. 65.

    Haas BJ, Papanicolaou A, Yassour M, Grabherr M, Blood PD, Bowden J, et al. De novo transcript sequence reconstruction from RNA-seq using the Trinity platform for reference generation and analysis. Nat Protoc. 2013;8:1494.

  66. 66.

    Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30:2114–20.

  67. 67.

    MacManes MD. On the optimal trimming of high-throughput mRNA sequence data. Front Genet. 2014;5:13.

  68. 68.

    Waterhouse RM, Seppey M, Simão FA, Manni M, Ioannidis P, Klioutchnikov G, et al. BUSCO applications from quality assessments to gene prediction and phylogenomics. Mol Biol Evol. 2017;35:543–8.

  69. 69.

    Simão FA, Waterhouse RM, Ioannidis P, Kriventseva EV, Zdobnov EM. BUSCO: assessing genome assembly and annotation completeness with single-copy orthologs. Bioinformatics. 2015;31:3210–2.

  70. 70.

    Bryant DM, Johnson K, DiTommaso T, Tickle T, Couger MB, Payzin-Dogru D, et al. A tissue-mapped axolotl de novo transcriptome enables identification of limb regeneration factors. Cell Rep. 2017;18:762–76.

  71. 71.

    Huerta-Cepas J, Forslund K, Coelho LP, Szklarczyk D, Jensen LJ, von Mering C, et al. Fast genome-wide functional annotation through orthology assignment by eggNOG-mapper. Mol Biol Evol. 2017;34:2115–22.

  72. 72.

    Huerta-Cepas J, Szklarczyk D, Forslund K, Cook H, Heller D, Walter MC, Rattei T, Mende DR, Sunagawa S, Kuhn M, Jensen LJ. eggNOG 4.5: a hierarchical orthology framework with improved functional annotations for eukaryotic, prokaryotic and viral sequences. Nucleic Acids Res. 2015;20:D286–93.

  73. 73.

    Buchfink B, Xie C, Huson DH. Fast and sensitive protein alignment using DIAMOND. Nat Methods. 2015;12:59–60.

  74. 74.

    Eddy SR. Accelerated profile HMM searches. PLoS Comput Biol. 2011;7:e1002195.

  75. 75.

    Davidson NM, Oshlack A. Corset: enabling differential gene expression analysis for de novo assembled transcriptomes. Genome Biol. 2014;15:410.

  76. 76.

    Robinson MD, McCarthy DJ, Smyth GK. edgeR: a Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010;26:139–40.

  77. 77.

    Hall TA. BioEdit: a user-friendly biological sequence alignment editor and analysis program for Windows 95/98/NT. Nucleic Acids Symp Ser. 1999;41:95–8.

Download references

Acknowledgements

We thank Caroline Beck for assisting in the laboratory, and Brodie Foster for assisting with field collections. In addition, we thank Aaron Jeffs (Otago Genomics & Bioinformatics Facility) for conducting the RNA library preparation.

Funding

This work was supported by Marsden Contract UOO1412 (Royal Society of New Zealand).

Author information

Affiliations

  1. Department of Zoology, University of Otago, P.O. Box 56, Dunedin, 9054, New Zealand
    • Graham A. McCulloch
    • , Andrew J. Veale
    • , Charles W. Ellen
    •  & Jonathan M. Waters
  2. Genomics Aotearoa and Department of Biochemistry, University of Otago, P.O. Box 56, Dunedin, 9054, New Zealand
    • Andrew Oliphant
    •  & Peter K. Dearden
  3. Landcare Research, Private Bag 92170, Auckland Mail Centre, Auckland, 1142, New Zealand
    • Andrew J. Veale
Authors
  1. Graham A. McCulloch
    View author publications
    You can also search for this author in
  2. Andrew Oliphant
    View author publications
    You can also search for this author in
  3. Peter K. Dearden
    View author publications
    You can also search for this author in
  4. Andrew J. Veale
    View author publications
    You can also search for this author in
  5. Charles W. Ellen
    View author publications
    You can also search for this author in
  6. Jonathan M. Waters
    View author publications
    You can also search for this author in

Contributions

JW, PD and AV designed the study. AV and JW conducted the field work. AV and CE extracted and prepared RNA for RNA-seq. AO, PD, and GM analysed the data. GM wrote the initial draft of the manuscript, with significant input from JW, AO, and PD. All authors read and approved the final manuscript.

Corresponding author

Correspondence to Graham A. McCulloch.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Supplementary information

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

McCulloch, G.A., Oliphant, A., Dearden, P.K. et al. Comparative transcriptomic analysis of a wing-dimorphic stonefly reveals candidate wing loss genes. EvoDevo 10, 21 (2019). https://doi.org/10.1186/s13227-019-0135-4

Download citation

Keywords

  • Vestigial winged
  • Wing development
  • Gene expression
  • Zelandoperla fenestrata

EvoDevo

ISSN: 2041-9139