Skip to main content

Stage- and sex-specific transcriptome analyses reveal distinctive sensory gene expression patterns in a butterfly

Abstract

Background

Animal behavior is largely driven by the information that animals are able to extract and process from their environment. However, the function and organization of sensory systems often change throughout ontogeny, particularly in animals that undergo indirect development. As an initial step toward investigating these ontogenetic changes at the molecular level, we characterized the sensory gene repertoire and examined the expression profiles of genes linked to vision and chemosensation in two life stages of an insect that goes through metamorphosis, the butterfly Bicyclus anynana.

Results

Using RNA-seq, we compared gene expression in the heads of late fifth instar larvae and newly eclosed adults that were reared under identical conditions. Over 50 % of all expressed genes were differentially expressed between the two developmental stages, with 4,036 genes upregulated in larval heads and 4,348 genes upregulated in adult heads. In larvae, upregulated vision-related genes were biased toward those involved with eye development, while phototransduction genes dominated the vision genes that were upregulated in adults. Moreover, the majority of the chemosensory genes we identified in the B. anynana genome were differentially expressed between larvae and adults, several of which share homology with genes linked to pheromone detection, host plant recognition, and foraging in other species of Lepidoptera.

Conclusions

These results revealed promising candidates for furthering our understanding of sensory processing and behavior in the disparate developmental stages of butterflies and other animals that undergo metamorphosis.

Peer Review reports

Background

The environment is teeming with information, and the ability to perceive and process this information is critical in shaping the behavior of all animals. Of the various sensory modalities, vision and chemoreception play integral roles in survival and reproduction, including the detection of food sources [1], predator avoidance [2], and locating potential mates [3]. Moreover, both senses are known to drive assortative mating and speciation processes. Visual cues, such as ornaments [4], coloration [5], and courtship displays [6], influence mate choice behaviors and sexual selection in a diverse range of species. Similarly, chemical signals, such as pheromones and cuticular hydrocarbons, have been found to be involved with prezygotic isolation in a wide variety of taxa, ranging from insects [7,8,9] and annelids [10] to mammals [11,12,13].

Despite the significant roles that visual and chemical cues play in animal behavior and sexual selection, considerable morphological differences often exist for the sensory structures that detect these cues throughout ontogeny. This is particularly apparent in animals that undergo metamorphosis from larva to adult life stages, such as holometabolous insects [14], crustaceans [15], and many fishes [16]. For instance, in butterflies, the visual organs of the larval stage typically consist of up to six stemmata per eye, each with a lens and seven photoreceptors that form a tiered rhabdom [17], compared to the much more complex adult compound eyes, which consist of hundreds of tightly packed ommatidia, each containing a facet lens and rhabdom composed of nine photoreceptors [18]. These differences are likely in part due to the different ecological niches that each stage fills; larvae typically reside and forage on host plants until pupation, while adults adopt an aerial lifestyle and are mainly focused on finding a mate and reproducing [14].

Although differences in sensory organ morphology and the behavior of animals that undergo metamorphosis are often readily apparent throughout ontogeny, we still have much to learn about the functional and organizational differences in the sensory systems of pre- and post-metamorphosis life stages, especially at the molecular level. Perhaps one of the most promising taxa in which to dissect these differences is the exceptionally diverse Insecta, which is estimated to consist of 5.5 million species [19]. Indeed, much of what we know about the molecular mechanisms underlying vision and chemosensation has been derived from work on insects, including the common fruit fly, Drosophila melanogaster (see [20] and [21] for review).

Phototransduction in insects is accomplished in the eye through absorption of light by a visual pigment (rhodopsin), which triggers an enzymatic cascade that ultimately leads to depolarization of photoreceptor cells [20]. The perception of different wavelengths of light is dependent upon opsin structure, with peak sensitivities spanning the visible light spectrum and beyond [22]. By contrast, chemosensation in insects occurs at the olfactory sensilla (the sensory structures involved with smell) typically found on head structures, such as the maxillary palps and antennae, and the gustatory sensilla (the sensory structures involved with taste), which are found throughout the insect body, including on the mouthparts, wings, and legs [23, 24]. Odorants are bound by odorant binding proteins (OBPs) or chemosensory proteins (CSPs) and transported through the sensillar lymph to membrane-bound receptors located on the dendrites of olfactory sensory neurons (OSNs) or gustatory sensory neurons (GSNs) [25, 26].

There are three main types of chemoreceptors on chemosensory neurons that are involved with the detection of chemical stimuli from the external environment in insects: odorant receptors (ORs), gustatory receptors (GRs), and ionotropic receptors (IRs). ORs are the foundation of olfaction and are known to selectively detect a diversity of volatile compounds [27,28,29]. In combination with a co-receptor (Orco) and sensory neuron membrane proteins (SNMPs), some ORs have also been found to be involved with the detection of sex pheromones [30, 31]. Insect GRs belong to the same superfamily as ORs [32] but are primarily involved with tasting bitter compounds [33,34,35], sugars [36,37,38], and CO2 [39, 40]. Finally, IRs, which are primitive chemoreceptor proteins that evolved from ionotropic glutamate receptors (iGluRs), are known to be involved with both olfaction and gustation, primarily sensing amines, acids, salt, and pheromones [24, 41,42,43].

Recent work in adult Lepidoptera has focused on elucidating the underpinnings of phototransduction [44] and chemoreception [24, 31, 45, 46], providing a foundation for investigating how sensory systems vary throughout development in an insect order known for its metamorphosis. The squinting bush brown butterfly, Bicyclus anynana, is an ideal model to address this issue, as it has rapidly become a fruitful model system for studying development, evolution, and phenotypic plasticity [47,48,49,50]. Of particular interest, these butterflies rely heavily on visual and chemical cues for mate choice; mates are selected based on the quality of ultraviolet-reflective eyespot pupils and male-specific pheromones [4, 51, 52]. In addition, previous work has identified differences in the visual systems of male and female adults of different seasonal phenotypes, including differences in eye size, facet lens area, facet number per eye, and opsin and eye development gene expression [53, 54]. Importantly, numerous molecular resources are available for this species, including expressed sequence tags [55], transcriptomes [54], and a reference genome assembly [56], making B. anynana amenable to genetic and genomic studies.

Here, we characterized the sensory gene repertoire in B. anynana, comprising genes known to be linked to vision and chemosensation. Specifically, we first identified vision genes in the B. anynana genome involved with phototransduction, eye pigmentation, and eye development, as well as six distinct families of chemosensory genes, consisting of OBPs, CSPs, ORs, IRs, GRs, and SNMPs. We also identified developmental genes that have possibly been co-opted to function as sensory or neural processing genes, such as those known to be involved with wing patterning in B. anynana and other butterflies, which are hypothesized to potentially drive speciation and assortative mating by linking wing pattern traits to preference for those traits [57,58,59]. We then compared the expression patterns of these sensory and developmental genes in the heads of two life stages: late fifth instar larvae and newly eclosed adult butterflies (Fig. 1).

Fig. 1
figure 1

Developmental stages of B. anynana. (A) Late fifth instar larva. (B) Newly eclosed adult. Scale bars are approximate

We predicted that genes directly involved with visual processes (e.g., phototransduction) would be upregulated in the adult phenotype, given the much greater level of complexity of adult compound eyes compared to the relative simplicity of larval stemmata. In addition, we hypothesized that larvae might not express the ultraviolet-sensitive opsin that is critical to eyespot quality evaluation during mate choice in adults. For chemosensory genes, we predicted that genes involved with pheromone and fruit detection would be upregulated in adults, as this stage participates in numerous reproductive behaviors and must locate a food source (i.e., ripe or rotting fruit) separate from the host plant. By contrast, we hypothesized that chemosensory genes linked to host plant recognition and foraging behavior would be upregulated in larvae, given that feeding is the dominant behavior during this stage of development. Furthermore, if genes associated with wing patterning in butterflies are also important for behavioral aspects of assortative mating, we predicted that they would be upregulated in the brains of adults. Finally, we aimed to identify candidate visual and chemosensory genes in the adult and larval phenotypes for future investigation into the sensory ecology of these disparate life stages.

Results

Sequencing generated over 387 million high quality single-end (SE) reads (Additional file 1: Table S1). Adapter trimming removed 87,301 reads (0.02 % of the raw sequenced reads) prior to downstream analysis. Approximately 340 million (88 %) of the remaining trimmed reads mapped to the B. anynana genome (v1.2; [56]; http://ensembl.lepbase.org/index.html). Of the 22,642 predicted protein-coding genes in the reference genome, 15,735 (70 %) were overlapped by at least 10 reads across all libraries and used as the expression set for downstream analyses. This gene set corresponded to approximately 12.2 ± 2.3 SD million reads per sample that were used for differential expression analysis (Table S1). Principal components analysis revealed that developmental stage explained 75 % of the variation observed (Additional file 1: Fig. S1). Blast2GO analysis resulted in the functional annotation of 13,498 (60 %) genes in the B. anynana genome, with a total of 40,857 gene ontology (GO) terms assigned to genes in the assembly.

Adult vs. larval heads

A total of 8,384 (53 % of genes in the expression set) genes were differentially expressed between larva and adult stage heads, with 4,348 upregulated in adult heads and 4,036 upregulated in larva heads (FDR < 0.05; Fig. 2; Additional file 2: Table S2). GO enrichment analyses of these upregulated genes found that 255 GO terms were enriched in the heads of adults (Additional file 2: Table S3). When reduced to the most specific terms (i.e., parent functions with a significant child GO term were removed to reduce redundancy), 63 enriched GO terms remained, with the top three terms being oxidation-reduction process (FDR = 1.53 × 10− 18), proton transmembrane transporter activity (FDR = 4.12 × 10− 10), and potassium ion transmembrane transport (FDR = 1.26 × 10− 9) (Table 1; see Additional file 2: Table S4 and Additional file 1: Figs. S2-S4 for full results). By contrast, 212 GO terms were enriched in the heads of larvae (Additional file 2: Table S5). A total of 49 GO terms remained after reduction, with the top three terms being nucleolus (FDR = 4.15 × 10− 11), mRNA splicing, via spliceosome (FDR = 4.26 × 10− 11), and protein folding (FDR = 1.17 × 10− 10) (Table 2; see Additional file 2: Table S6 and Additional file 1: Figs. S5-S7 for full results).

Fig. 2
figure 2

Volcano plot of the false discovery rate (-log10FDR) and expression ratio (log2FC) for each gene in B. anynana adult heads relative to larval heads. Differentially expressed vision, chemosensory, and wing patterning gene homologs are highlighted in blue, red, and orange, respectively. Positive log2FC values indicate higher expression in adults, while negative log2FC values indicate higher expression in larvae. DEGs = Differentially Expressed Genes. This plot was created with ggplot2 v3.3.2 [60] in R v3.6.2 [61]

Table 1 Top 10 most specific GO terms enriched in adult heads. GO = Gene Ontology, BP = Biological Process, MF = Molecular Function, FDR = False Discovery Rate
Table 2 Top 10 most specific GO terms enriched in larval heads. GO = Gene Ontology, BP = Biological Process, MF = Molecular Function, CC = Cellular Component, FDR = False Discovery Rate

Sex-specific differences within each stage

A comparison between male and female adult heads revealed 27 differentially expressed genes, 10 of which were upregulated in male adults, while 17 were upregulated in female adults (FDR < 0.05; Fig. 3A; Additional file 2: Table S7). By contrast, 37 genes were differentially expressed between male and female larval heads, with nine upregulated in male larvae and 28 upregulated in female larvae (FDR < 0.05; Fig. 3B; Additional file 2: Table S8). GO enrichment analyses found no significantly enriched GO terms for the adult or larva differentially expressed gene sets. Furthermore, none of the vision, chemosensory, or wing patterning/development genes identified in the reference genome (see below) were differentially expressed between the sexes of larvae or adults. There were three genes in common between the two stage-specific differentially expressed gene sets, including Putative 115 kDa protein in type-1 retrotransposable element R1DM-like protein (BANY.1.2.g05985) and two copies of neuralized-like protein 4 (BANY.1.2.g11289 and BANY.1.2.g11290), all of which were upregulated in female larvae and adults.

Fig. 3
figure 3

Expression heatmaps of differentially expressed genes for sex-specific comparisons within each developmental stage. (A) Male adults vs. female adults. (B) Male larvae vs. female larvae. Counts were normalized by variance stabilizing transformation, with warmer colors indicating higher expression. Rows denote individual genes, and columns denote samples, both of which are clustered by gene expression. Family indicates the family from which the sample was derived, and Sex indicates the sex of the sample. Heatmaps were created with pheatmap v1.0.12 [62] in R v3.6.2 [61]

Vision genes

Blast hits for 252 of the 274 queried vision genes (92 %) resulted in the identification of 1,555 putative homologs in the B. anynana genome (Additional file 2: Table S9). Of these homologs, 429 were associated with phototransduction, 76 with eye pigment, and 1,050 with eye development. To identify the top homolog candidates for each of the queried vision genes, we collected the best blast hit, resulting in a set of 252 B. anynana vision genes, 250 of which were within the head expression set (Additional file 2: Table S10). Of these top homologs, 165 (65 %) were differentially expressed between larval and adult heads (FDR < 0.05), with 88 (57 phototransduction genes, 10 eye pigment genes, and 21 eye development genes) upregulated in adults and 77 (13 phototransduction genes, 5 eye pigment genes, and 59 eye development genes) upregulated in larvae (Figs. 2 and 4). Only genes associated with phototransduction were significantly enriched in the full differentially expressed gene set (phototransduction, FDR = 1.12 × 10− 5; eye pigment, FDR = 0.13; eye development, FDR = 0.26). Out of the 252 identified vision gene homologs, only five phototransduction genes showed evidence of sex-/stage-specific expression, with the expression of nan_trpv being absent in male larvae and Cib2, pteropsin, santa_maria, and unclassified not showing expression in female larvae (Fig. 5; Additional file 2: Table S11).

Fig. 4
figure 4

Expression heatmaps for differentially expressed top vision homologs in B. anynana. (A) Phototransduction genes. (B) Eye pigment genes. (C) Eye development genes. Counts were normalized by variance stabilizing transformation, with warmer colors indicating higher expression. Rows denote individual genes, and columns denote samples, both of which are clustered by gene expression. Family indicates the family from which the sample was derived, Sex indicates the sex of the sample, and Stage indicates the developmental stage of the sample. Heatmaps were created with pheatmap v1.0.12 [62] in R v3.6.2 [61]

Fig. 5
figure 5

Presence/absence of expression maps for the top vision homologs in B. anynana. Rows denote individual genes, and columns denote sample group. Genes that were expressed are indicated in black, while red indicates genes that were not expressed. Genes colored in grey were identified in the B. anynana genome but were not present in the head expression set (i.e., not expressed in any group). Category indicates the gene family to which each gene belongs. Expression maps were created with pheatmap v1.0.12 [62] in R v3.6.2 [61]

A total of seven opsins (including three visual opsins: UVRh, BRh, and LWRh) were identified in the expression set. While all seven opsins were expressed in both developmental stages, each was significantly upregulated in adults relative to larvae (UVRh, log2FC = 8.72, FDR = 1.99 × 10− 191; BRh, log2FC = 7.84, FDR = 1.58 × 10− 68; LWRh, log2FC = 9.37, FDR < 2.22 × 10− 308; Rh7, log2FC = 1.15, FDR = 5.84 × 10− 4; pteropsin, log2FC = 3.65, FDR = 2.84 × 10− 4; unclassified, log2FC = 14.71, FDR = 1.19 × 10− 52; RGR-like, log2FC = 4.96, FDR = 6.07 × 10− 10; Fig. 6; Additional file 2: Table S10).

Fig. 6
figure 6

Log transformed normalized counts of visual opsin genes in larvae and adults. All opsins were expressed in both sexes of both stages. Horizontal lines within the boxes denote the median. The upper and lower bounds of the boxes indicate the 25th and 75th percentiles, and whiskers extend to the largest count value ≤ 1.5 × the interquartile range. Y-axes are best fit for each gene. Boxplots were created using ggplot2 v3.3.2 [60] in R v3.6.2 [61]

In addition to the Heliconius melpomene and D. melanogaster vision homologs we identified, manual searches of the Blast2GO functional annotation identified an additional 23 vision-related genes, including numerous genes putatively associated with phototransduction and eye development (Additional file 2: Table S12). Of these genes, 14 were differentially expressed between larval and adult heads, with seven upregulated in adults and seven upregulated in larvae (Additional file 1: Fig. S8).

Chemosensory genes

Odorant binding proteins

Blast hits for 28 of the 273 queried OBP genes resulted in the identification of 48 putative homologs in the B. anynana genome (Additional file 2: Table S13). We retained only those containing pfam01395 or smart00708 domains, which resulted in a set of 19 B. anynana OBP genes, 17 being in the head expression set (Additional file 2: Table S14). Of these homologs, 13 (76 %) were differentially expressed (FDR < 0.05), 12 of which were upregulated in adults and one of which was upregulated in larvae (Fig. 7A). OBPs were not significantly enriched in the full differentially expressed gene set (FDR = 0.10). Five OBPs showed stage-/sex-specific expression, including a homolog of Hmel-OBP12 (BANY.1.2.g14367) that was only expressed in the adult stage, and Dple-OBP19 (BANY.1.2.g20356), which was expressed in female but not male adults (Fig. 8; Additional file 2: Table S11). Three other OBPs exhibited sex-specific expression in larvae, with HmelPBP_C (BANY.1.2.g06880) and Dple-PBP-D (BANY.1.2.g06881) showing male-specific expression and Hmel-OBP22 (BANY.1.2.g19953) showing female-specific expression.

Fig. 7
figure 7

Expression heatmaps for differentially expressed top chemosensory homologs in B. anynana. (A) Odorant binding proteins. (B) Chemosensory proteins. (C) Odorant receptors. (D) Ionotropic receptors. (E) Gustatory receptors. (F) Sensory neuron membrane proteins. Counts were normalized by variance stabilizing transformation, with warmer colors indicating higher expression. Rows denote individual genes, and columns denote samples, both of which are clustered by gene expression. Family indicates the family from which the sample was derived, Sex indicates the sex of the sample, and Stage indicates the developmental stage of the sample. Heatmaps were created with pheatmap v1.0.12 [62] in R v3.6.2 [61]

Fig. 8
figure 8

Presence/absence of expression maps for the top chemosensory homologs in B. anynana. Rows denote individual genes, and columns denote sample group. Genes that were expressed are indicated in black, while red indicates genes that were not expressed. Genes colored in grey were identified in the B. anynana genome but were not present in the head expression set (i.e., not expressed in any group). Category indicates the gene family to which each gene belongs. OBP = odorant binding protein, CSP = chemosensory protein, OR = odorant receptor, IR = ionotropic receptor, GR = gustatory receptor, SNMP = sensory neuron membrane protein. Expression maps were created with pheatmap v1.0.12 [62] in R v3.6.2 [61]

Chemosensory proteins

Blast hits for 18 of the 34 queried CSP genes and one gene identified in our manual search resulted in the identification of 27 putative homologs in the B. anynana genome (Additional file 2: Table S15). Twenty-four of these candidate CSP genes contained the pfam03392 domain, 22 of which were in the head expression set (Additional file 2: Table S16). Of these homologs, 19 (86 %) were differentially expressed (FDR < 0.05), with 11 upregulated in adults and eight upregulated in larvae (Fig. 7B). CSPs were significantly enriched in the full differentially expressed gene set (FDR = 3.56 × 10− 3). While there was evidence that all CSPs in the expression set were expressed in both life stages, HmCSP11 (BANY.1.2.g12993) was not expressed in female larvae (Fig. 8; Additional file 2: Table S11).

Odorant receptors

Blast hits for 38 of the 70 OR genes and four genes found via a manual search resulted in the identification of 50 putative homologs in the B. anynana genome (Additional file 2: Table S17). In total, 43 of these were retained as B. anynana OR genes by confirmation of the presence of either the pfam02949 or pfam08395 protein domain (Additional file 2: Table S18). Of these homologs, 38 were in the head expression set, 31 (82 %) of which were differentially expressed between larvae and adults (FDR < 0.05). These differentially expressed genes consisted of 30 that were upregulated in adults and one that was upregulated in larvae (Fig. 7C). Additionally, ORs were significantly enriched in the full differentially expressed gene set (FDR = 1.19 × 10− 3). Thirteen ORs were found to exhibit stage-specific expression, all of which were expressed only in adults (Fig. 8; Additional file 2: Table S11). Moreover, 10 ORs showed sex-specific expression in larvae, with six only expressed in male heads (HmOr13, HmOr19, HmOr36, HmOr42, HmOr54, and HmOr61) and 4 only expressed in female heads (HmOr27, HmOr51, HmOr57, and HmOr70).

Ionotropic receptors

We mapped 24 of the 31 B. anynana IR sequences from [63] to genes in the reference genome, 19 of which were within the head expression set (Additional file 2: Table S19). Of these homologs, 12 (63 %) were differentially expressed (FDR < 0.05), all of which were upregulated in adults (Fig. 7D). However, IRs were not significantly enriched in the full differentially expressed gene set (FDR = 0.34). Three IRs were found to only be expressed in adult heads: BanyIR31a, BanyIR1.2, and BanyIR40a (Fig. 8; Additional file 2: Table S11). In addition, four IRs were expressed in a sex-specific fashion in larvae, with three only expressed in male larvae (BanyIR21a, BanyIR76b, and BanyIR7d.2) and one only expressed in female larvae (BanyIR87a).

Gustatory receptors

Blast hits for 24 of the 73 GR genes resulted in the identification of 39 putative homologs in the B. anynana genome (Additional file 2: Table S20). We retained only those that contained the pfam08395 domain, resulting in a set of 27 B. anynana GR genes, 16 of which were in the head expression set (Additional file 2: Table S21). Of these homologs, six (38 %) were differentially expressed (FDR < 0.05), all of which were upregulated in adults (Fig. 7E). GRs were not significantly enriched in the full differentially expressed gene set (FDR = 0.94). Four GRs, consisting of homologs of HmGr2, HmGr5, and HmGr19, showed adult-specific expression (Fig. 8; Additional file 2: Table S11). Furthermore, one HmGr12 homolog (BANY.1.2.g06465) showed male-specific expression in adults, and six GRs showed sex-specific expression in larvae, with four exhibiting male-specific expression (HmGr27, HmGr44, HmGr52, and HmGr68) and two exhibiting female-specific expression (HmGr57 and HmGr73).

Sensory neuron membrane proteins

We identified 16 putative SNMP homologs, consisting of blast hits for nine of the 33 SNMP query genes (Additional file 2: Table S22). Filtering these putative homologs for genes that contained the pfam01130 domain and were annotated as SNMP genes in the functional annotation resulted in a set of six B. anynana SNMP genes (Additional file 2: Table S23). All of these were within the head expression set, with three (50 %) being differentially expressed between larvae and adults, each of which were upregulated in adults (Fig. 7F). SNMPs, however, were not significantly enriched in the full differentially expressed gene set (FDR = 0.81). Two SNMP homologs showed sex-specific expression in larvae, with SexiSNMP3 (BANY.1.2.g07849) exhibiting male-specific expression and HarmSNMP3 (BANY.1.2.g07846) exhibiting female-specific expression (Fig. 8; Additional file 2: Table S11).

Developmental and wing patterning genes

A total of 52 genes associated with wing patterning in butterflies were found in the head expression set of larval and adult B. anynana (Additional file 2: Table S24; see [64] and [65] for butterfly wing patterning genes). These genes include homologs for al, antp, ap, BarH-1, CD63, Ci, Dll, dpp, dsx, EcR, en, Hh, inv, N, optix, ptc, sal, wg, and several Wnt genes. Of these, 30 (58 %) were differentially expressed between larvae and adults (FDR < 0.05), including homologs for al, BarH-1, CD63, Ci, Dll, dpp, en, Hh, inv, N, ptc, sal, wg, Wnt, and WntA (Figs. 2 and 9). Copies of two genes known to be involved with eyespot development in B. anynana, CD63 (BANY.1.2.g25497) and Ci (BANY.1.2.g11922), were found to exhibit expression specific to adult and larval heads, respectively (Additional file 2: Table S11). Additionally, a copy of BarH-1 (BANY.1.2.g16347) exhibited larva-specific expression, while a Wnt-5 homolog (BANY.1.2.g04762) was expressed in all groups except female larvae.

Fig. 9
figure 9

Expression heatmap of differentially expressed genes linked to wing patterning. Counts were normalized by variance stabilizing transformation, with warmer colors indicating higher expression. Rows denote individual genes, and columns denote samples, both of which are clustered by gene expression. Family indicates the family from which the sample was derived, Sex indicates the sex of the sample, and Stage indicates the developmental stage of the sample. This heatmap was created with pheatmap v1.0.12 [62] in R v3.6.2 [61]

Discussion

Our analysis of the gene expression profiles of larval and adult B. anynana heads revealed considerable differences between the two developmental stages, with > 50 % of all expressed genes showing differential expression. Furthermore, we identified numerous genes involved with vision and chemosensation and elucidated how the expression of these genes, as well as the expression of known wing patterning genes, differs in the head throughout ontogeny and between the sexes. More than 250 B. anynana genes putatively linked to vision-related processes were discovered to be expressed in the head, including genes associated with phototransduction as well as eye pigmentation and development. A total of 143 homologs associated with chemosensation were identified, comprising odorant binding proteins, chemosensory proteins, odorant receptors, ionotropic receptors, gustatory receptors, and sensory neuron membrane proteins. In addition, we found 52 genes previously described as butterfly wing patterning genes that were expressed in larval and/or adult heads, including a WntA homolog, a gene known to play a significant role in wing pattern development across a diversity of nymphalid species [66, 67]. To our knowledge, this study is the first attempt to characterize the sensory gene repertoire of Bicyclus butterflies and provides a promising resource for investigating differences in the sensory biology of larvae and adult butterflies.

Overall expression differences between larval and adult heads

In larval heads, upregulated genes were linked to developmental processes, including multicellular organism development and the Wnt signaling pathway (Table 2 and Additional file 2: Table S6). Wnt signaling is known to be involved with cell differentiation and proliferation in animals [68, 69]. In butterflies, Wnt genes have been found to be involved with wing patterning [67, 70] and are expressed in various B. anynana tissues during embryogenesis, with Wnt7 and Wnt11 both expressed in head tissues [71]. In the current study, we found that numerous Wnt genes are expressed in the heads of fifth instar larvae and adults, including Wnt-1, Wnt-5, Wnt-6, Wnt-10, Wnt-11, and WntA homologs, most of which were upregulated in larval heads.

In addition to developmental processes, genes upregulated in larval heads were also enriched for processes linked to gene expression and protein metabolism. Upregulation of genes involved with these functional categories is possibly in part due to the physiological changes taking place in larval head tissues in preparation for pupation and metamorphosis. As with other holometabolous insects, the larval tissues and organs of butterflies undergo degeneration via autophagy and subsequent remodeling during metamorphosis (see [72] for review). Moreover, in Manduca sexta, metamorphic cell death is associated with a marked drop in protein synthesis [73]. Therefore, the enrichment of GO terms involved with gene expression and protein metabolism might be indicative of a similar decrease in protein synthesis in B. anynana.

In adult heads, chemoreception processes were enriched, with associated GO terms including detection of chemical stimulus involved in sensory perception of smell, olfactory receptor activity, and odorant binding. These results are consistent with elevated chemosensory gene expression in the adult stage of other Lepidoptera [74] and suggest a greater investment in the chemosensory system of adult B. anynana. Indeed, the adult stage partakes in numerous behaviors that larvae do not, including courtship/copulation, oviposition, and foraging for fruit. Furthermore, B. anynana adult females are known to cue in on pheromones produced by males for mate choice [51, 52, 75], and chemical cues appear to play just as important of a role as visual cues in mate choice for this species. Finally, the ability of flight permits these butterflies to perform these behaviors in a greatly expanded three-dimensional space compared to the larval stage, possibly requiring adults to maintain a more sensitive and sophisticated chemosensory system.

Insights into larval and adult phototransduction

A large number of vision-related genes were expressed in B. anynana heads, most of which were differentially expressed between adults and larvae. These results have implications for furthering our understanding of the differences in the visual capabilities and phototransduction signaling cascade for different life stages of lepidopterans and other holometabolous insects.

The primary visual organs of larval and adult butterflies have disparate morphologies, with larvae possessing two simple eyes consisting of up to six stemmata and adults having two compound eyes consisting of hundreds of ommatidia [76]. Therefore, it is likely that at least some of the observed patterns in vision-related gene expression in larval and adult heads are due to substantial differences in structure and cell composition. Moreover, during metamorphosis in holometabolous insects, the larval stemmata migrate to the adult optic lobe and continue to function as extraretinal photoreceptors [17, 76, 77]. Consequently, the presence of both adult and larval visual structures in B. anynana adults might account for a portion of the upregulation observed in vision genes.

Interestingly, the differentially expressed vision genes upregulated in adult heads were dominated by phototransduction genes (65 %), while differentially expressed genes upregulated in larval heads were largely associated with eye development (77 %). A greater emphasis on phototransduction in adults is perhaps not surprising, as a significant proportion of the adult head consists of eye tissue, and optic lobes have been found to comprise nearly 75 % of the butterfly brain [78,79,80]. In comparison, the stemmata of larvae occupy a considerably smaller proportion of the larval head. The upregulation of genes involved with eye development in late fifth instar larvae suggests that compound eye developmental processes have initiated just prior to pupation.

Several phototransduction genes were upregulated in larvae relative to adults, including a copy of norpA, wunen2, and the innexins ogre, inx2, and inx3. norpA encodes the protein phospholipase C (PLC), which is involved with diacylglycerol (DAG) production in the Drosophila phototransduction cascade [81, 82]. Moreover, Macias-Muñoz et al. (2019) hypothesized that wunen plays a similar role in Lepidoptera phototransduction as lazaro in Drosophila, which is involved in DAG level regulation [83]. Lepidoptera have three copies of wunen, and while wunen2 was upregulated in larval B. anynana, the other two (wunen and wunen3) were upregulated in adults. The observed stage-biased expression of specific norpA and wunen copies suggests potential differences in DAG regulation throughout development.

Finally, ogre, inx2, and inx3 form gap junction channels, all of which are critical to visual transmission. Specifically, ogre and inx3 are necessary for visual synaptic transmission in retinal pigment cells in the compound eyes of Drosophila, while inx2 plays an essential role in laminar glial cells [84]. The fact that these genes are upregulated in larvae suggests the possibility that gap junctions might be either more integral to larval phototransduction or present in greater density in larval eyes. Future functional work should explore this possibility.

Intriguingly, none of the identified vision-related genes, including the visual opsins, were differentially expressed between male and female larvae or male and female adults, a result consistent with previous studies with B. anynana [53, 54]. It is also interesting to note that all visual opsins were expressed in larval and adult heads, suggesting that both developmental stages might be capable of perceiving similar wavelengths of light. Future electrophysiological and behavioral studies should explore the spectral sensitivity and behavioral responses of larval and adult B. anynana to different light wavelengths.

From the current study, it is not possible to determine which of these genes are expressed in the eyes, brain, and/or other head structures. Numerous studies have localized phototransduction genes in nonvisual tissues, such as the central nervous system [85, 86]. Therefore, future work should investigate tissue-specific expression and determine the expression patterns of these genes in the eyes. While significant efforts have been made to explore butterfly vision in a number of species, usually focusing on opsins, the results of the current study provide a new set of candidate vision genes for B. anynana and will help to expand our understanding of lepidopteran adult, as well as larval, vision.

Insights into larval and adult chemosensation

In the current study, we identified a total of 143 chemosensory genes, most of which showed differential expression between larval and adult B. anynana heads. Notably, we discovered numerous OBPs, CSPs, ORs, IRs, GRs, and SNMPs with stage-/sex-specific expression or that displayed differential expression between the developmental stages. Many of these genes share homology with chemosensory genes associated with pheromone detection, host plant recognition, and foraging in other species of Lepidoptera. Because the functions and specificity of chemosensory genes in B. anynana are largely unknown, these genes serve as promising targets for further investigation to expand our understanding of chemically mediated behaviors in this species.

A number of the chemosensory genes identified are candidates for pheromone binding and reception. The gene BANY.1.2.g22938 was upregulated in adult heads and is homologous to Hmel-OBP13 (alternative name = HmOBP20; [87]), which was recently found to possibly be involved with species-specific recognition of pheromones in Heliconius butterflies [88]. In addition to the two genes (BANY.1.2.g06880 and BANY.1.2.g06881) sharing sequence homology with H. melpomene and Danaus plexippus pheromone binding proteins (HmelPBP_C and Dple-PBP-D, respectively), BANY.1.2.g22938 might also be involved with pheromone detection in B. anynana. Moreover, two putative pheromone receptor homologs (HmOr3 and HmOr30; [87, 89, 90]) were upregulated in adult heads, suggesting that they might also be involved with chemical-mediated mate choice behaviors. Finally, a homolog of Dple-OBP19 (BANY.1.2.g20356) was expressed in female but not male adult heads, indicating a putative role in male sex pheromone (MSP) detection [75].

In addition to the OBP and OR repertoire, we identified two homologs of SNMP1 (both MsexSNMP1 homologs) that were upregulated in adult heads. SNMP1 a protein that forms a complex with pheromone-detecting ORs and an odorant receptor co-receptor (Orco; HmOR2; identified in the current study as BANY.1.2.g12855 in B. anynana) in insects [30, 31]. SNMP1 is involved with pheromone detection in both D. melanogaster [30, 91] and numerous lepidopteran species [31, 92,93,94,95], and it may also play a role in pheromone detection in B. anynana.

Several genes were identified as putatively involved with host plant recognition and/or foraging behavior. A homolog of Eobl-GOBP2 (BANY.1.2.g06879), a gene involved with the detection of plant volatiles in the moth Ectropis obliqua [96], was upregulated in adult heads. In addition, a homolog of HmOR49 (BANY.1.2.g06204), a putative citral receptor [87, 97], was upregulated in adult heads. Citral is a plant volatile that is present in plant species such as lemongrass and orange [98] and serves as a food attractant for Bombyx mori larvae [99] and an oviposition deterrent in the light brown apple moth, Epiphyas postvittana [100]. Moreover, electroantennography (EAG) recordings indicated that citral evokes a response in B. anynana antennae (Murphy, Joshi, and Westerman, unpublished data). Homologs for two GRs (HmGr9 and HmGr57) that are characterized as putatively being involved with host plant identification via recognition of the plant alkaloid synephrine in H. melpomene [24, 88] were also expressed in B. anynana heads. These genes might also be involved with host plant recognition in B. anynana. Alternatively, they might be involved with detection of B. anynana’s adult food source, ripe/rotting fruit, as synephrine, like citral, is present in citrus fruits [101].

Two IRs upregulated in adult B. anynana, BanyIR1.2 and BanyIR75d, are putatively involved with host plant searching behavior, as IR1.2 and IR75d are upregulated in antennae of mated females of the moth Helicoverpa armigera [63]. In B. anynana, the expression of BanyIR1.2 was specific to adults, consistent with a possible function in oviposition-related behaviors in this species.

Our results also illuminate genes that play potentially important roles in larval chemoreception. One OBP (a homolog for Dple-OBP2), eight CSPs (homologs for HmCSP3, HmCSP7, HmCSP13, HmCSP14, HmCSP16, and HmCSP17), and one OR (a HmOr18 homolog) were upregulated in larval heads relative to adult heads. Previous studies of OBP and CSP gene expression in larvae and adult H. armigera discovered six OBP and four CSP genes that are exclusively expressed in larvae antennae and mouthparts, suggesting that OBP and CSP genes may play a role in larval foraging [102]. Moreover, numerous ORs in the moth Spodoptera littoralis were found to be tuned to plant volatiles and are involved with larval foraging behavior [103, 104]. Thus, it is possible that the chemosensory genes upregulated in B. anynana larvae might play important roles in mediating larval foraging behavior.

Expression of wing patterning genes

Wing patterning genes have been hypothesized to underlie assortative mating behaviors and ultimately speciation in Lepidoptera through associations with preference for the traits they influence [58, 59]. This might occur in two ways: (1) both the trait and preference are controlled by the same gene; or (2) the genes controlling the trait and preference for that trait are separate but maintained in high linkage disequilibrium (i.e., inherited together) [57, 105, 106]. Empirical evidence for either of these hypotheses, or for the genetic basis of assortative mate preference more broadly, is relatively slim. However, numerous wing patterning genes are known to influence sensory organ and neural processes in other insect species, providing a promising set of candidate genes for exploration. For instance, optix, an eye development gene in Drosophila [107], has been co-opted to control red pigmentation in the wings of Heliconius butterflies [108]. Furthermore, engrailed (en), a gene involved with neurogenesis [109], axonal targeting [110], and neuronal cell fate determination [111], is also linked to butterfly eyespot development in Bicyclus [112]. Here, we found that numerous genes known to be involved with wing patterning in butterflies were expressed in B. anynana heads, possibly in the brain, eyes, or both tissues. If these genes, particularly those involved with eyespot development in B. anynana, are linked to preferences for eyespot traits, they might play a role in the great amount of diversity we see in this taxon (80 + species, with many species living in sympatry) [113, 114]. We propose that these wing patterning genes should be investigated as potential drivers of assortative mate preference and speciation in Bicyclus butterflies.

Lack of differential expression of sensory and wing patterning genes between the sexes

It is interesting to note that none of the vision, chemosensory, or wing patterning genes we identified here were differentially expressed between the sexes of either stage. This suggests that developmental stage is likely a larger factor than sex for the expression of these genes in the head. Alternatively, it is also possible that any effects of sex on vision, chemosensory, or wing patterning genes were either too small or too tissue specific to be detected in our data set. For instance, the number of replicates for the sex comparisons was relatively small (Adults: n = 3 males, n = 3 females; Larvae: n = 4 males, n = 2 females); increasing the number of replicates would better capture the biological variability attributed to sex and result in greater power to detect differentially expressed genes. Furthermore, it is possible that increasing tissue specificity (e.g., sequencing the brain neuropils and various sensory tissues separately) might reveal sex-specific gene expression patterns that are obscured when sequencing the whole head.

Conclusions

In this study, we identified the sensory gene repertoire of the butterfly B. anynana and characterized the expression of these genes in larval and adult heads. While visual and chemosensory genes have been explored in many adult Lepidoptera, few studies have investigated the expression of such genes in their larval stages. Our results provide an initial step in elucidating the differences in sensory processing throughout development in butterflies. Moreover, we identified numerous candidate genes for host plant recognition, foraging, and mate choice, including both chemosensory and wing patterning genes expressed in B. anynana heads. Future studies should explore the functions of these candidate genes and determine their tissue specificity.

Methods

Animals

Bicyclus anynana, a nymphalid butterfly native to subtropical Africa, has been maintained in laboratory colonies since 1988. All animals used in this study are descendants of an original population established in Leiden, Netherlands from 80 gravid females that were collected in Malawi [47]. The population at the University of Arkansas was established via the transfer of ~ 1,000 eggs from a population in Singapore to Fayetteville, AR, USA in spring, 2017. Animals were reared in a climate-controlled, USDA-APHIS approved (Permit # P526P-17-00343) greenhouse facility, which was maintained at approximately 27°C, 70 % relative humidity to induce the wet season phenotype in this species [47]. All experiments were conducted between January and March 2019 (sunrise range: 7:06–7:19 am, sunset range: 5:40–7:36 pm) under a 13:11 h light:dark photoperiod. In addition to natural light, the greenhouse was illuminated with full spectrum (including ultraviolet wavelengths) fluorescent lights (lights on: 7:00 am, off: 8:00 pm).

Experimental design and tissue collection

Four unique families were created by pairing one three-day-old virgin male and one three-day-old virgin female together in a small mesh cage (31.8 cm × 31.8 cm × 31.8 cm) at 8:00 am for at least three hours to ensure that copulation occurred. After visual confirmation that the pair had copulated, the female was removed from the mating cage and isolated in a new large mesh cage (39.9 cm × 39.9 cm × 59.9 cm) containing a corn plant (Zea mays) on which to lay eggs and a slice of moistened banana for food. Each female was then given seven days to lay fertilized eggs on the provided corn plant, after which the egg-laden corn plant was transferred to a new small mesh cage (31.8 cm × 31.8 cm × 31.8 cm).

Upon hatching, larvae were reared in their family-specific cages under identical conditions and were fed corn plants ad libitum. To ensure that all four families experienced the same environmental conditions within the greenhouse and to control for any potential unforeseen confounding variables associated with cage location, the physical position of each cage was alternated daily. Upon the morning of reaching the late fifth instar stage, which was determined by the stark change in color from tan/brown to green (Fig. 1A), a subset of the larvae from each family was sacrificed by decapitation with RNase-free scissors (n = 6 larvae total). This stage was chosen to ensure that all larvae were as close as possible in development and because it is the final developmental stage prior to pupation. A second subset from each family was allowed to pupate, and newly eclosed adults (Fig. 1B) were sacrificed by decapitation on the morning of emergence (n = 6 adults total). All decapitations were conducted between 9:30 am–12:00 pm, and heads were immediately transferred into RNase-free, low binding 1.5 ml microcentrifuge tubes (Biotix, San Diego, CA, USA), flash-frozen in liquid nitrogen, and transported to the lab for storage at -80°C until they were processed (Additional file 1: Table S25).

RNA extraction, library preparation, and sequencing

Each frozen head was immersed in pre-chilled RNAlater-ICE (Ambion; Austin, TX, USA) and incubated at -20°C for approximately 16 h prior to tissue processing. After this incubation period, heads were transferred to a dissecting dish filled with RNAlater-ICE, and all residual thoracic tissue was carefully removed with forceps under a dissecting microscope (Zeiss Stemi 508; Jena, Germany), leaving only head tissue. Individual isolated heads (which included antennae and mouthparts) were then disrupted in lysis buffer with an RNase-free, disposable pestle, and small (< 200 nucleotides) and large RNA (> 200 nucleotides) were extracted in separate fractions using the NucleoSpin® miRNA Kit (Macherey-Nagel; Düren, Germany) following the manufacturer’s recommended protocols. RNA purity, concentration, and integrity for each sample were subsequently determined using a NanoDrop 2000 (Thermo Fisher Scientific; Waltham, MA, USA) and TapeStation 2200 (Agilent; Santa Clara, CA, USA).

After confirmation of RNA quality and quantity, a cDNA library for each head was prepared using 500 ng of large RNA as input for the KAPA mRNA HyperPrep Kit (KAPA Biosystems; Wilmington, MA, USA) combined with the KAPA Unique Dual-Indexed Adapter Kit (KAPA Biosystems; Wilmington, MA, USA). The quality of each cDNA library was subsequently verified using a TapeStation 2200 (Agilent; Santa Clara, CA, USA). All libraries (n = 12) were then shipped on dry ice to the University of Chicago Genomics Facility for secondary quality assessment on a 5300 Fragment Analyzer (Agilent; Santa Clara, CA, USA), and 50 base pair (bp) single-end (SE) sequencing was performed on a single lane of a HiSeq 4000 (Illumina; San Diego, CA, USA).

Animal sexing

All adults were sexed based on apparent sexually dimorphic features, specifically the presence/absence of androconia (the male-specific pheromone organ). Because larvae do not display any obvious sexual dimorphism, we extracted DNA from individuals using the KAPA Express Extract Kit (KAPA Biosystems; Wilmington, MA, USA) and amplified a female-specific W-chromosome microsatellite ([115]; Genbank accession no.: AY785080) using the primers from [116]. PCR products were then visualized by gel electrophoresis, and females were identified by the presence of a band at ~ 185 bp, while males lacked this band (Additional file 1: Fig. S9).

Functional annotation

Blast2GO v5.2.5 [117] was used to conduct a de novo functional annotation of all genes in the most current B. anynana reference genome (v1.2; [56]; http://ensembl.lepbase.org/index.html). First, we used BLASTX v2.6.0+ [118] to search the NCBI ‘nr’ protein database (www.ncbi.nlm.nih.gov) and collected the top 10 hits with an E-value < 10− 3. These results were then uploaded into Blast2GO, and further functional classification was performed using the InterProScan [119] function within Blast2GO. Finally, the “Mapping” and “Annotation” steps in Blast2GO were performed using the default parameters, and the resulting functional annotation table was exported.

Differential gene expression analysis

Prior to expression quantification, the quality of the raw reads was assessed using FastQC v0.11.8 (https://www.bioinformatics.babraham.ac.uk/projects/fastqc/), and Illumina adapter sequences were trimmed using Trimmomatic v0.38 [120]. Trimmed reads were then aligned to the B. anynana reference genome (v1.2) using STAR v2.7.1a [121] with the default parameters and the “--twopassMode Basic” option. Reads were quantified using the ‘htseq-count’ script from the HTSeq v0.11.2 Python package [122]. Differential gene expression analysis was conducted using the DESeq2 v1.24.0 package [123] in R v3.6.2 [61].

We conducted two separate differential expression analyses. First, the generalized linear model:

$$y\:\sim\:family\:+\:sex\:+\:stage$$

was fit to each gene using a negative binomial distribution, where y denotes the response variable (expression), family denotes the family to which each individual belongs (family 1–4), sex denotes the sex of each animal (male or female), and stage denotes the life stage of each individual (larva or adult). Using this design enabled us to contrast the effect of stage while controlling for differences in expression associated with lineage and sex. Second, the generalized linear model:

$$y\:\sim\:group$$

was used, where group denotes a grouping variable that combines sex and stage (i.e., male larva, female larva, male adult, and female adult). This design permitted us to contrast the effect of sex for each stage (i.e., male larva vs. female larva and male adult vs. female adult). We note the suboptimal number of replicates (n = 2) for female larvae in the male vs. female larvae analysis; however, we opted to perform this analysis in an attempt to identify any sensory genes that show substantial differential expression between the sexes of this stage. Additionally, because a transcriptome-wide comparison between the sexes of caterpillars has never been conducted (to our knowledge), this analysis will also provide preliminary insights into sex-biased gene expression in lepidopteran larvae.

For both analyses, genes with a total read alignment count < 10 were filtered and not included in the differential expression analysis. Gene expression was calculated as the binary log of the expression fold change (log2FC), and the apeglm method was used for log2FC shrinkage to obtain the most accurate estimates of effect size [124]. Finally, genes with a false discovery rate (FDR; [125]) < 0.05 were retained for downstream analysis.

Gene ontology enrichment analyses

For further characterization, the Fisher’s Exact Test function in Blast2GO was used to test for GO term enrichment. The set of differentially expressed genes identified for the first analysis (y ~ family + sex + stage model) was split into genes that showed increased expression in adults (log2FC > 0) and those that showed increased expression in larvae (log2FC < 0), and each subset was tested separately. For the second analysis (y ~ group model), the differentially expressed gene sets were not split due to the small number of genes in each set. The reference set used for all GO enrichment analyses consisted of all genes in the expression set, and only GO terms with an FDR < 0.05 were considered significantly enriched. The list of enriched GO terms for each analysis was then reduced to the most specific terms for visualization. Additionally, the reduced lists of enriched GO terms were processed using REVIGO (http://revigo.irb.hr/; [126]), which further eliminated redundancy and organized GO terms into treemaps consisting of related superclusters.

Identification of visual genes

To identify genes involved with vision (i.e., phototransduction, eye pigment, and eye development) in B. anynana, we first collected the coding sequences (CDS) of 74 putative H. melpomene phototransduction genes from [44] and the protein sequences of 200 D. melanogaster phototransduction, eye pigment, and eye development genes compiled by [54]. We then used BLASTX and BLASTP (BLAST v2.2.30+; [118]) to query these sequences against the B. anynana reference genome proteins in Lepbase (http://blast.lepbase.org/) and identify homologs. Homologs were determined based on hits with an E-value < 1 × 10− 10, and the top candidate for each query gene was identified as the hit with the lowest E-value. In cases where numerous hits had identical E-values, ties were broken by selecting the hit with the highest bit score. Finally, to identify additional putative vision genes, we manually searched the Blast2GO annotation descriptions, best blast hits, and GO annotations for terms linked to vision, including: “eye,” “ommatidia,” “ommatidium,” “opsin,” “photoreceptor,” “phototransduction,” “retina,” and “visual”. Fisher’s exact tests were conducted to test if the identified vision-related genes were enriched in the differentially expressed gene set, with all expressed genes as the reference set.

To explore stage-/sex-specific expression, we investigated the normalized counts for individuals in each group (i.e., male larva, female larva, male adult, and female adult). Non-zero counts for any individual(s) within a group for a given gene were considered evidence for expression, while zero counts for all individuals in a group were considered evidence that a given gene was not expressed (note: it is possible that these genes are actually lowly expressed but at levels below the detection threshold for the sequencing depth of the current study).

Identification of chemosensory genes

To identify genes involved with chemosensation in B. anynana, we collected 273 lepidopteran and D. melanogaster OBP protein sequences from [45], 34 H. melpomene CSP protein sequences from [87], 70 H. melpomene OR protein sequences from [87], 31 B. anynana IR sequences from [63], 73 H. melpomene GR protein sequences from [24], and 33 lepidopteran SNMP sequences from [31]. These sequences were then queried against the B. anynana reference genome with BLASTX or BLASTP in Lepbase (http://blast.lepbase.org/) to identify putative homologs. With the exception of the previously identified IR sequences in B. anynana, all hits with an E-value < 1 × 10− 10 were further screened for conserved protein domains specific to each gene family using CD-Search [127]. Specifically, sequences with hits for the following domains were retained: OBPs, either pfam01395 (PBP/GOBP family) or smart00708 (Insect pheromone/odorant binding protein domains); CSPs, pfam03392 (Insect pheromone-binding family, A10/OS-D); ORs, either pfam02949 (7tm Odorant receptor) or pfam08395 (7tm Chemosensory receptor); GRs, pfam08395 (7tm Chemosensory receptor); and SNMPs, pfam01130 (CD36 family). Because the CD36 superfamily common to SNMPs consists of three different protein families, only one of which includes SNMPs [128], we filtered the final putative SNMP sequences by only retaining those that were also annotated as SNMPs in our functional annotation.

In addition, we performed a manual search of the Blast2GO functional annotation to identify any additional putative OBP, CSP, OR, IR, GR, and SNMP genes. Specifically, we searched the Blast2GO descriptions and best blast hits for key terms, including: “odorant binding protein,” “pheromone binding protein,” “chemosensory protein,” “ejaculatory bulb-specific protein 3,” “odorant receptor,” “olfactory receptor,” “gustatory receptor,” and “sensory neuron membrane protein” and subjected any putative chemosensory genes to the conserved protein domain filtration described above. Fisher’s exact tests were used to test if the identified chemosensory gene categories were enriched in the differentially expressed gene set, with all expressed genes composing the reference set. P-values resulting from these six tests, along with the p-values resulting from the vision gene enrichment tests, were corrected for multiple comparisons (FDR) [125]. Finally, sex-/stage-specific gene expression was determined as previously described for the vision genes.

Identification of wing patterning genes

To further explore the expression of genes that might be relevant to sensory processing and signaling during mate choice, we manually searched our de novo functional annotation for genes known to be involved with wing patterning, many of which have been hypothesized to also be involved in mate preference as well as other behaviors. These genes included several known B. anynana wing patterning genes, such as those coding the proteins Antennapedia (antp) [129, 130], apterous (ap) [131], CD63 antigen (CD63) [132], Cubitus interruptus (Ci) [133], decapentaplegic (dpp) [134], Distal-less (Dll) [135], doublesex (dsx) [116], Ecdysone Receptor (EcR) [136], Engrailed (en) [137], hedgehog (Hh) [129], Invected (inv) [137], Notch (N) [138], patched (ptc) [138], Spalt (sal) [137], Ultrabithorax (Ubx) [130, 133], and wingless (wg) [132], as well as genes known to be critical for wing patterning in Heliconius and other butterflies, including aristaless [139], BarH-1 [140], cortex [141], optix [108], and Wnt [67, 70].

Availability of data and materials

Raw sequence data associated with this study are accessible through the NCBI Sequence Read Archive (SRA) database under BioProject ID PRJNA730880. All other data presented in this study are available within this manuscript and its additional files.

Abbreviations

bp:

base pair

cDNA:

Complementary DNA

CDS:

Coding sequence

CSP:

Chemosensory protein

DAG:

Diacylglycerol

DEGs:

Differentially expressed genes

EAG:

Electroantennography

FDR:

False discovery rate

GO:

Gene ontology

GR:

Gustatory receptor

GSN:

Gustatory sensory neuron

iGluRs:

Ionotropic glutamate receptors

IR:

Ionotropic receptor

log2FC:

Log2(fold change)

MSP:

Male sex pheromone

OBP:

Odorant binding protein

OR:

Odorant receptor

OSN:

Olfactory sensory neuron

PCR:

Polymerase chain reaction

PLC:

Phospholipase C

RNA-seq:

RNA-sequencing

SD:

Standard deviation

SE:

Single-end

SNMP:

Sensory neuron membrane protein

References

  1. Nevitt GA. Sensory ecology on the high seas: The odor world of the procellariiform seabirds. J Exp Biol. 2008;211:1706–13. doi:https://doi.org/10.1242/jeb.015412.

    Article  PubMed  Google Scholar 

  2. Fischer S, Oberhummer E, Cunha-Saraiva F, Gerber N, Taborsky B. Smell or vision? The use of different sensory modalities in predator discrimination. Behav Ecol Sociobiol. 2017;71:143. doi:https://doi.org/10.1007/s00265-017-2371-8.

    Article  PubMed  PubMed Central  Google Scholar 

  3. Shine R, Webb JK, Lane A, Mason RT. Mate location tactics in garter snakes: effects of rival males, interrupted trails and non-pheromonal cues. Funct Ecol. 2005;19:1017–24. doi:https://doi.org/10.1111/J.1365-2435.2005.01063.X.

    Article  Google Scholar 

  4. Robertson KA, Monteiro A. Female Bicyclus anynana butterflies choose males on the basis of their dorsal UV-reflective eyespot pupils. Proc R Soc B Biol Sci. 2005;272:1541–6. doi:https://doi.org/10.1098/rspb.2005.3142.

    Article  Google Scholar 

  5. Bakker TCM, Mundwiler B. Female mate choice and male red coloration in a natural three-spined stickleback (Gasterosteus aculeatus) population. Behav Ecol. 1994;5:74–80. doi:https://doi.org/10.1093/beheco/5.1.74.

    Article  Google Scholar 

  6. Chapman FM. The courtship of Gould’s manakin (Manacus vitellinus vitellinus) on Barro Colorado Island, Canal Zone. Bull Am Museum Nat Hist. 1935;68:471–525.

    Google Scholar 

  7. Coyne JA, Crittenden AP, Mah K. Genetics of a pheromonal difference contributing to reproductive isolation in Drosophila. Science. 1994;265:1461–4. doi:https://doi.org/10.1126/science.8073292.

    Article  CAS  PubMed  Google Scholar 

  8. Raymond B, Searle JB, Douglas AE. On the processes shaping reproductive isolation in aphids of the Aphis fabae (Scop.) complex (Aphididae: Homoptera). Biol J Linn Soc. 2001;74:205–15. doi:https://doi.org/10.1111/j.1095-8312.2001.tb01387.x.

    Article  Google Scholar 

  9. Pélozuelo L, Meusnier S, Audiot P, Bourguet D, Ponsard S. Assortative mating between European corn borer pheromone races: beyond assortative meeting. PLoS One. 2007;2:e555. doi:https://doi.org/10.1371/journal.pone.0000555.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  10. Sutton R, Bolton E, Bartels-Hardege HD, Eswards M, Reish DJ, Hardege JD. Chemical signal mediated premating reproductive isolation in a marine polychaete, Neanthes acuminata (arenaceodentata). J Chem Ecol. 2005;31:1865–76. doi:https://doi.org/10.1007/s10886-005-5931-8.

    Article  CAS  PubMed  Google Scholar 

  11. Moore RE. Olfactory discrimination as an isolating mechanism between Peromyscus maniculatus and Peromyscus polionotus. Am Midl Nat. 1965;73:85–100. doi: https://doi.org/10.2307/2423324.

    Article  Google Scholar 

  12. Nevo E, Bodmer M, Heth G. Olfactory discrimination as an isolating mechanism in speciating mole rats. Experientia. 1976;32:1511–2. doi:https://doi.org/10.1007/BF01924423.

    Article  CAS  PubMed  Google Scholar 

  13. Laukaitis CM, Critser ES, Karn RC. Salivary androgen-binding protein (abp) mediates sexual isolation in Mus musculus. Evolution. 1997;51:2000–5. doi:https://doi.org/10.1111/j.1558-5646.1997.tb05121.x.

    Article  CAS  PubMed  Google Scholar 

  14. Truman JW, Riddiford LM. The origins of insect metamorphosis. Nature. 1999;401:447–52. doi:https://doi.org/10.1038/46737.

    Article  CAS  PubMed  Google Scholar 

  15. Passano LM. The regulation of crustacean metamorphosis. Am. Zoologist. 1961;1:89–95.

    Article  Google Scholar 

  16. McMenamin SK, Parichy DM. Metamorphosis in teleosts. Curr Top Dev Biol. 2013;103:127–65. doi:https://doi.org/10.1016/B978-0-12-385979-2.00005-8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  17. Gilbert C. Form and function of stemmata in larvae of holometabolous insects. Annu Rev Entomol. 1994;39:323–49. doi:https://doi.org/10.1146/annurev.en.39.010194.001543.

    Article  Google Scholar 

  18. Arikawa K. The eyes and vision of butterflies. J Physiol. 2017;595:5457–64. doi:https://doi.org/10.1113/JP273917.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  19. Stork NE. How many species of insects and other terrestrial arthropods are there on Earth? Annu Rev Entomol. 2018;63:31–45. doi:https://doi.org/10.1146/annurev-ento-020117-043348.

    Article  CAS  PubMed  Google Scholar 

  20. Hardie RC. Phototransduction in Drosophila melanogaster. J Exp Biol. 2001;204:3403–9. doi: https://doi.org/10.1242/jeb.204.20.3403.

    Article  CAS  PubMed  Google Scholar 

  21. Benton R. Chemical sensing in Drosophila. Curr Opin in Neurobiol. 2008;18:357–63. doi:https://doi.org/10.1016/j.conb.2008.08.012.

    Article  CAS  Google Scholar 

  22. Terakita A. The opsins. Genome Biol. 2005;6:213. doi:https://doi.org/10.1186/gb-2005-6-3-213.

    Article  PubMed  PubMed Central  Google Scholar 

  23. Dahanukar A, Hallem EA, Carlson JR. Insect chemoreception. Curr Opin Neurobiol. 2005;15:423–30. doi:https://doi.org/10.1016/j.conb.2005.06.001.

    Article  CAS  PubMed  Google Scholar 

  24. Briscoe AD, Macias-Muñoz A, Kozak KM, Walters JR, Yuan F, Jamie GA, et al. Female behaviour drives expression and evolution of gustatory receptors in butterflies. PLoS Genet. 2013;9:e1003620. doi:https://doi.org/10.1371/journal.pgen.1003620.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Vogt RG, Riddiford LM. Pheromone binding and inactivation by moth antennae. Nature. 1981;293:161–3. doi:https://doi.org/10.1038/293161a0.

    Article  CAS  PubMed  Google Scholar 

  26. Pelosi P, Calvello M, Ban L. Diversity of odorant-binding proteins and chemosensory proteins in insects. Chem Senses; 2005; 30 Suppl 1:i291-2. doi:https://doi.org/10.1093/chemse/bjh229.

    Article  CAS  PubMed  Google Scholar 

  27. Sato K, Pellegrino M, Nakagawa T, Nakagawa T, Vosshall LB, Touhara K. Insect olfactory receptors are heteromeric ligand-gated ion channels. Nature. 2008;452:1002–6. doi:https://doi.org/10.1038/nature06850.

    Article  CAS  PubMed  Google Scholar 

  28. Ha TS, Smith DP. Odorant and pheromone receptors in insects. Front Cell Neurosci. 2009;3:10. doi:https://doi.org/10.3389/neuro.03.010.2009.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  29. Hansson BS, Stensmyr MC. Evolution of insect olfaction. Neuron. 2011;72:698–711. doi:https://doi.org/10.1016/j.neuron.2011.11.003.

    Article  CAS  PubMed  Google Scholar 

  30. Benton R, Vannice KS, Vosshall LB. An essential role for a CD36-related receptor in pheromone detection in Drosophila. Nature. 2007;450:289–93. doi:https://doi.org/10.1038/nature06328.

    Article  CAS  PubMed  Google Scholar 

  31. Zhang HJ, Xu W, Chen QM, Sun LN, Anderson A, Xia QY, et al. A phylogenomics approach to characterizing sensory neuron membrane proteins (SNMPs) in Lepidoptera. Insect Biochem Mol Biol. 2020;118:103313. doi:https://doi.org/10.1016/j.ibmb.2020.103313.

    Article  CAS  PubMed  Google Scholar 

  32. Clyne PJ, Warr CG, Carlson JR. Candidate taste receptors in Drosophila. Science. 2000;287:1830–4. doi:https://doi.org/10.1126/science.287.5459.1830.

    Article  CAS  PubMed  Google Scholar 

  33. Moon SJ, Köttgen M, Jiao Y, Xu H, Montell C. A taste receptor required for the caffeine response in vivo. Curr Biol. 2006;16:1812–7. doi:https://doi.org/10.1016/j.cub.2006.07.024.

    Article  CAS  PubMed  Google Scholar 

  34. Lee Y, Moon SJ, Montell C. Multiple gustatory receptors required for the caffeine response in Drosophila. Proc Natl Acad Sci U S A. 2009;106:4495–500. doi:https://doi.org/10.1073/pnas.0811744106.

    Article  PubMed  PubMed Central  Google Scholar 

  35. Weiss LA, Dahanukar A, Kwon JY, Banerjee D, Carlson JR. The molecular and cellular basis of bitter taste in Drosophila. Neuron. 2011;69:258–72. doi:https://doi.org/10.1016/j.neuron.2011.01.001.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  36. Dahanukar A, Foster K, Van der Goes van Naters WM, Carlson JR. A Gr receptor is required for response to the sugar trehalose in taste neurons of Drosophila. Nat Neurosci. 2001;4:1182–6. doi:https://doi.org/10.1038/nn765.

  37. Chyb S, Dahanukar A, Wickens A, Carlson JR. Drosophila Gr5a encodes a taste receptor tuned to trehalose. Proc Natl Acad Sci U S A. 2003;100 Suppl 2:14526–30. doi:https://doi.org/10.1073/pnas.2135339100.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  38. Slone J, Daniels J, Amrein H. Sugar receptors in Drosophila. Curr Biol. 2007;17:1809–16. doi:https://doi.org/10.1016/j.cub.2007.09.027.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  39. Jones WD, Cayirlioglu P, Kadow IG, Vosshall LB. Two chemosensory receptors together mediate carbon dioxide detection in Drosophila. Nature. 2007;445:86–90. doi:https://doi.org/10.1038/nature05466.

    Article  CAS  PubMed  Google Scholar 

  40. Kwon JY, Dahanukar A, Weiss LA, Carlson JR. The molecular basis of CO2 reception in Drosophila. Proc Natl Acad Sci U S A. 2007;104:3574–8. doi:https://doi.org/10.1073/pnas.0700079104.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  41. Rytz R, Croset V, Benton R. Ionotropic receptors (IRs): chemosensory ionotropic glutamate receptors in Drosophila and beyond. Insect Biochem Mol Biol. 2013;43:888–97. doi:https://doi.org/10.1016/j.ibmb.2013.02.007.

    Article  CAS  PubMed  Google Scholar 

  42. Zhang YV, Ni J, Montell C. The molecular basis for attractive salt-taste coding in Drosophila. Science. 2013;340:1334–8. doi:https://doi.org/10.1126/science.1234133.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  43. Koh TW, He Z, Gorur-Shandilya S, Menuz K, Larter NK, Stewart S, et al. The Drosophila IR20a clade of ionotropic receptors are candidate taste and pheromone receptors. Neuron. 2014;83:850–65. doi:https://doi.org/10.1016/j.neuron.2014.07.012.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  44. Macias-Muñoz A, Rangel Olguin AG, Briscoe AD. Evolution of phototransduction genes in Lepidoptera. Genome Biol Evol. 2019;11:2107–24. doi:https://doi.org/10.1093/gbe/evz150.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  45. Vogt RG, Große-Wilde E, Zhou JJ. The Lepidoptera odorant binding protein gene family: gene gain and loss within the GOBP/PBP complex of moths and butterflies. Insect Biochem Mol Biol. 2015;62:142–53. doi:https://doi.org/10.1016/j.ibmb.2015.03.003.

    Article  CAS  PubMed  Google Scholar 

  46. Liu W, Jiang XC, Cao S, Yang B, Wang GR. Functional studies of sex pheromone receptors in Asian corn borer Ostrinia furnacalis. Front Physiol. 2018;9:591. doi:https://doi.org/10.3389/fphys.2018.00591.

    Article  PubMed  PubMed Central  Google Scholar 

  47. Brakefield PM, Reitsma N. Phenotypic plasticity, seasonal climate and the population biology of Bicyclus butterflies (Satyridae) in Malawi. Ecol Entomol. 1991;16:291–303. doi:https://doi.org/10.1111/j.1365-2311.1991.tb00220.x.

    Article  Google Scholar 

  48. Koch PB, Brakefield PM, Kesbeke F. Ecdysteroids control eyespot size and wing color pattern in the polyphenic butterfly Bicyclus anynana (Lepidoptera: Satyridae). J Insect Physiol. 1996;42:223–30. doi: https://doi.org/10.1016/0022-1910(95)00103-4.

    Article  CAS  Google Scholar 

  49. Kooi RE, Brakefield PM. The critical period for wing pattern induction in the polyphenic tropical butterfly Bicyclus anynana (Satyrinae). J Insect Physiol. 1999;45:201–12. doi:https://doi.org/10.1016/s0022-1910(98)00093-6.

    Article  CAS  PubMed  Google Scholar 

  50. Prudic KL, Jeon C, Cao H, Monteiro A. Developmental plasticity in sexual roles of butterfly species drives mutual sexual ornamentation. Science. 2011;331:73–5. doi:https://doi.org/10.1126/science.1197114.

    Article  CAS  PubMed  Google Scholar 

  51. Costanzo K, Monteiro A. The use of chemical and visual cues in female choice in the butterfly Bicyclus anynana. Proc R Soc B Biol Sci. 2007;274:845–51. doi:https://doi.org/10.1098/rspb.2006.3729.

    Article  Google Scholar 

  52. Westerman EL, Monteiro A. Odour influences whether females learn to prefer or to avoid wing patterns of male butterflies. Anim Behav. 2013;86:1139–45. doi: https://doi.org/10.1016/j.anbehav.2013.09.002.

    Article  Google Scholar 

  53. Everett A, Tong X, Briscoe AD, Monteiro A. Phenotypic plasticity in opsin expression in a butterfly compound eye complements sex role reversal. BMC Evol Biol. 2012;12:232. doi:https://doi.org/10.1186/1471-2148-12-232.

    Article  PubMed  PubMed Central  Google Scholar 

  54. Macias-Munoz A, Smith G, Monteiro A, Briscoe AD. Transcriptome-wide differential gene expression in Bicyclus anynana butterflies: female vision-related genes are more plastic. Mol Biol Evol. 2016;33:79–92. doi:https://doi.org/10.1093/molbev/msv197.

    Article  CAS  PubMed  Google Scholar 

  55. Beldade P, Rudd S, Gruber JD, Long AD. A wing expressed sequence tag resource for Bicyclus anynana butterflies, an evo-devo model. BMC Genomics. 2006;7:130. doi:https://doi.org/10.1186/1471-2164-7-130.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  56. Nowell RW, Elsworth B, Oostra V, Zwaan BJ, Wheat CW, Saastamoinen M, et al. A high-coverage draft genome of the mycalesine butterfly Bicyclus anynana. GigaScience. 2017;6:1–7. doi:https://doi.org/10.1093/gigascience/gix035.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  57. Westerman EL. Searching for the genes driving assortative mating. PLoS Biol. 2019;17: e3000108. doi:https://doi.org/10.1371/journal.pbio.3000108.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Merrill RM, Rastas P, Martin SH, Melo MC, Barker S, Davey J, et al. Genetic dissection of assortative mating behavior. PLoS Biol. 2019;17:e2005902. doi:https://doi.org/10.1371/journal.pbio.2005902.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  59. Kronforst MR, Young LG, Kapan DD, McNeely C, O’Neill RJ, Gilbert LE. Linkage of butterfly mate preference and wing color preference cue at the genomic location of wingless. Proc Natl Acad Sci U S A. 2006;103:6575–80. doi:https://doi.org/10.1073/pnas.0509685103.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  60. Wickham, H. ggplot2: elegant graphics for data analysis. Springer-Verlag New York, 2016.

    Book  Google Scholar 

  61. R Core Team. 2019. R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. URL https://www.R-project.org/.

  62. Kolde, R. 2019. pheatmap: pretty heatmaps. R package version 1.0.12. https://CRAN.R-project.org/package=pheatmap.

  63. Liu NY, Xu W, Dong SL, Zhu JY, Xu YX, Anderson A. Genome-wide analysis of ionotropic receptor gene repertoire in Lepidoptera with an emphasis on its functions of Helicoverpa armigera. Insect Biochem Mol Biol. 2018;99:37–53. doi:https://doi.org/10.1016/j.ibmb.2018.05.005.

    Article  CAS  PubMed  Google Scholar 

  64. Hanly JJ, Wallbank RWR, McMillan WO, Jiggins CD. Conservation and flexibility in the gene regulatory landscape of heliconiine butterfly wings. Evodevo. 2019;10:15. doi:https://doi.org/10.1186/s13227-019-0127-4.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  65. Banerjee TD, Monteiro A. Dissection of larval and pupal wings of Bicyclus anynana butterflies. Methods Protoc. 2020;3:5. doi:https://doi.org/10.3390/mps3010005.

    Article  CAS  PubMed Central  Google Scholar 

  66. Mazo-Vargas A, Concha C, Livraghi L, Massardo D, Wallbank RWR, Zhang L, et al. Macroevolutionary shifts of WntA function potentiate butterfly wing-pattern diversity. Proc Natl Acad Sci U S A. 2017;114:10701–6. doi:https://doi.org/10.1073/pnas.1708149114.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  67. Martin A, Reed RD. Wnt signaling underlies evolution and development of the butterfly wing pattern symmetry systems. Dev Biol. 2014;395:367–78. doi:https://doi.org/10.1016/j.ydbio.2014.08.031.

    Article  CAS  PubMed  Google Scholar 

  68. Logan CY, Nusse R. The Wnt signaling pathway in development and disease. Annu Rev Cell Dev Biol. 2004;20:781–810. doi:https://doi.org/10.1146/annurev.cellbio.20.010403.113126.

    Article  CAS  PubMed  Google Scholar 

  69. Wiese KE, Nusse R, van Amerongen R. Wnt signalling: conquering complexity. Development. 2018;145:dev165902. doi:https://doi.org/10.1242/dev.165902.

    Article  CAS  PubMed  Google Scholar 

  70. Martin A, Reed RD. Wingless and aristaless2 define a developmental ground plan for moth and butterfly wing pattern evolution. Mol Biol Evol. 2010;27:2864–78. doi:https://doi.org/10.1093/molbev/msq173.

    Article  CAS  PubMed  Google Scholar 

  71. Holzem M, Braak N, Brattström O, McGregor AP, Breuker CJ. Wnt gene expression during early embryogenesis in the nymphalid butterfly Bicyclus anynana. Front Ecol Evol. 2019;7:468. doi: https://doi.org/10.3389/fevo.2019.00468.

    Article  Google Scholar 

  72. Tettamanti G, Casartelli M. Cell death during complete metamorphosis. Philos TransR Soc Lond B Biol Sci. 2019;374:20190065. doi:https://doi.org/10.1098/rstb.2019.0065.

    Article  CAS  Google Scholar 

  73. Zakeri Z, Quaglino D, Latham T, Woo K, Lockshin RA. Programmed cell death in the tobacco hornworm, Manduca sexta: alteration in protein synthesis. Microsc Res Tech. 1996;34:192–201. doi:https://doi.org/10.1002/(SICI)1097-0029(19960615)34:3<192::AID-JEMT2>3.0.CO;2-S.

    Article  CAS  PubMed  Google Scholar 

  74. Yang CH, Yang PC, Li J, Yang F, Zhang AB. Transcriptome characterization of Dendrolimus punctatus and expression profiles at different developmental stages. PLoS One. 2016;11:e0161667. doi:https://doi.org/10.1371/journal.pone.0161667.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  75. Nieberding CM, de Vos H, Schneider M V., Lassance JM, Estramil N, Andersson J, et al. The male sex pheromone of the butterfly Bicyclus anynana: towards an evolutionary analysis. PLoS One. 2008;3:e2751. doi:https://doi.org/10.1371/journal.pone.0002751.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  76. Ichikawa T. Integration of colour signals in the medulla of the swallowtail butterfly larva. J Exp Biol. 1991;155:127–145. doi:https://doi.org/10.1242/jeb.155.1.127.

    Article  Google Scholar 

  77. Briscoe AD, White RH. Adult stemmata of the butterfly Vanessa cardui express UV and green opsin mRNAs. Cell Tissue Res. 2005;319:175–9. doi:https://doi.org/10.1007/s00441-004-0994-3.

    Article  CAS  PubMed  Google Scholar 

  78. Ali FA. Structure and metamorphosis of the brain and suboesophageal ganglion of Pieris brassicae (L.) (Lepidoptera: Pieridae). Trans R Entomol Soc London. 1974;125:363–412. doi:https://doi.org/10.1111/j.1365-2311.1974.tb02306.x.

    Article  Google Scholar 

  79. Sivinski J. Mushroom body development in nymphalid butterflies: a correlate of learning? J Insect Behav. 1989;2:277–83. doi:https://doi.org/10.1007/BF01053299.

    Article  Google Scholar 

  80. Heinze S, Reppert SM. Anatomical basis of sun compass navigation I: the general layout of the monarch butterfly brain. J Comp Neurol. 2012;520:1599–628. doi:https://doi.org/10.1002/cne.23054.

    Article  PubMed  Google Scholar 

  81. Bloomquist BT, Shortridge RD, Schneuwly S, Perdew M, Montell C, Steller H, et al. Isolation of a putative phospholipase C gene of Drosophila, norpA, and its role in phototransduction. Cell. 1988;54:723–33. doi:https://doi.org/10.1016/S0092-8674(88)80017-5.

    Article  CAS  PubMed  Google Scholar 

  82. Lee YJ, Shah S, Suzuki E, Zars T, O’Day PM, Hyde DR. The Drosophila dgq gene encodes a G alpha protein that mediates phototransduction. Neuron. 1994;13:1143–57. doi:https://doi.org/10.1016/0896-6273(94)90052-3.

    Article  CAS  PubMed  Google Scholar 

  83. Garcia-Murillas I, Pettitt T, Macdonald E, Okkenhaug H, Georgiev P, Trivedi D, et al. lazaro encodes a lipid phosphate phosphohydrolase that regulates phosphatidylinositol turnover during Drosophila phototransduction. Neuron. 2006;49:533–46. doi:https://doi.org/10.1016/j.neuron.2006.02.001.

    Article  CAS  PubMed  Google Scholar 

  84. Han Y, Xiong L, Xu Y, Tian T, Wang T. The β-alanine transporter BalaT is required for visual neurotransmission in Drosophila. Elife. 2017;6:e29146. doi:https://doi.org/10.7554/eLife.29146.

    Article  PubMed  PubMed Central  Google Scholar 

  85. Kingston ACN, Cronin TW. Short- and long-wavelength-sensitive opsins are involved in photoreception both in the retina and throughout the central nervous system of crayfish. J Comp Physiol A Neuroethol Sens Neural Behav Physiol. 2015;201:1137–45. doi:https://doi.org/10.1007/s00359-015-1043-2.

    Article  CAS  PubMed  Google Scholar 

  86. Donohue MW, Carleton KL, Cronin TW. Opsin expression in the central nervous system of the mantis shrimp Neogonodactylus oerstedii. Biol Bull. 2017;233:58–69. doi:https://doi.org/10.1086/694421.

    Article  CAS  PubMed  Google Scholar 

  87. Heliconius Genome Consortium. Butterfly genome reveals promiscuous exchange of mimicry adaptations among species. Nature. 2012;487:94–8. doi:https://doi.org/10.1038/nature11041.

    Article  CAS  Google Scholar 

  88. van Schooten B, Meléndez-Rosa J, van Belleghem SM, Jiggins CD, Tan JD, McMillan WO, et al. Divergence of chemosensing during the early stages of speciation. Proc Natl Acad Sci U S A. 2020;117:16438–47. doi:https://doi.org/10.1073/pnas.1921318117.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  89. Nakagawa T, Sakurai T, Nishioka T, Touhara K. Insect sex-pheromone signals mediated by specific combinations of olfactory receptors. Science. 2005;307:1638–42. doi:https://doi.org/10.1126/science.1106267.

    Article  CAS  PubMed  Google Scholar 

  90. Wanner KW, Anderson AR, Trowell SC, Theilmann DA, Robertson HM, Newcomb RD. Female-biased expression of odourant receptor genes in the adult antennae of the silkworm, Bombyx mori. Insect Mol Biol. 2007;16:107–19. doi:https://doi.org/10.1111/j.1365-2583.2007.00708.x.

    Article  CAS  PubMed  Google Scholar 

  91. Jin X, Ha TS, Smith DP. SNMP is a signaling component required for pheromone sensitivity in Drosophila. Proc Natl Acad Sci U S A. 2008;105:10996–1001. doi:https://doi.org/10.1073/pnas.0803309105.

    Article  PubMed  PubMed Central  Google Scholar 

  92. Rogers ME, Sun M, Lerner MR, Vogt RG. Snmp-1, a novel membrane protein of olfactory neurons of the silk moth Antheraea polyphemus with homology to the CD36 family of membrane proteins. J Biol Chem. 1997;272:14792–9. doi:https://doi.org/10.1074/jbc.272.23.14792.

    Article  CAS  PubMed  Google Scholar 

  93. Rogers ME, Krieger J, Vogt RG. Antennal SNMPs (sensory neuron membrane proteins) of Lepidoptera define a unique family of invertebrate CD36-like proteins. J Neurobiol. 2001;49:47–61. doi:https://doi.org/10.1002/neu.1065.

    Article  CAS  PubMed  Google Scholar 

  94. Rogers ME, Steinbrecht RA, Vogt RG. Expression of SNMP-1 in olfactory neurons and sensilla of male and female antennae of the silkmoth Antheraea polyphemus. Cell Tissue Res. 2001;303:433–46. doi:https://doi.org/10.1007/s004410000305.

    Article  CAS  PubMed  Google Scholar 

  95. Krieger J, Raming K, Dewer YME, Bette S, Conzelmann S, Breer H. A divergent gene family encoding candidate olfactory receptors of the moth Heliothis virescens. Eur J Neurosci. 2002;16:619–28. doi:https://doi.org/10.1046/j.1460-9568.2002.02109.x.

    Article  PubMed  Google Scholar 

  96. Zhang YL, Fu XB, Cui HC, Zhao L, Yu JZ, Li HL. Functional characteristics, electrophysiological and antennal immunolocalization of general odorant-binding protein 2 in tea geometrid, Ectropis oblique. Int J Mol Sci. 2018;19: 875. doi:https://doi.org/10.3390/ijms19030875.

    Article  CAS  PubMed Central  Google Scholar 

  97. Jordan MD, Anderson A, Begum D, Carraher C, Authier A, Marshall SDG, et al. Odorant receptors from the light brown apple moth (Epiphyas postvittana) recognize important volatile compounds produced by plants. Chem Senses. 2009;34:383–94. doi:https://doi.org/10.1093/chemse/bjp010.

    Article  CAS  PubMed  Google Scholar 

  98. Martins P, Sbaite P, Benites C, Wolf Maciel M. Thermal characterization of orange, lemongrass, and basil essential oils. Chem Eng Trans. 2011;24:463–8. doi:https://doi.org/10.3303/CET1124078.

    Article  Google Scholar 

  99. Hamamura Y, Naito KI. Food selection by silkworm larvæ, Bombyx mori: Citral, linalyl acetate, linalol, and terpinyl acetate as attractants of larvæ. Nature. 1961;190:879–80. doi:https://doi.org/10.1038/190879a0.

    Article  CAS  Google Scholar 

  100. Suckling DM, Karg G, Gibb AR, Bradley SJ. Electroantennogram and oviposition responses of Epiphyas postvittana (Lepidoptera: Tortricidae) to plant volatiles. New Zeal J Crop Hortic Sci. 1996;24:323–33. doi:https://doi.org/10.1080/01140671.1996.9513969.

    Article  Google Scholar 

  101. Stewart I, Newhall WF, Edwards GJ. The isolation and identification of l-synephrine in the leaves and fruit of citrus. J Biol Chem. 1964. 239;930–2. doi:https://doi.org/10.1016/S0021-9258(18)51679-3.

    Article  CAS  Google Scholar 

  102. Chang H, Ai D, Zhang J, Dong S, Liu Y, Wang G. Candidate odorant binding proteins and chemosensory proteins in the larval chemosensory tissues of two closely related noctuidae moths, Helicoverpa armigera and H. assulta. PLoS One. 2017;12: e0179243. doi:https://doi.org/10.1371/journal.pone.0179243.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  103. de Fouchier A, Walker WB 3rd, Montagné N, Steiner C, Binyameen M, Schlyter F, et al. Functional evolution of Lepidoptera olfactory receptors revealed by deorphanization of a moth repertoire. Nat Commun. 2017;8: 15709. doi:https://doi.org/10.1038/ncomms15709.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  104. de Fouchier A, Sun X, Caballero-Vidal G, Travaillard S, Jacquin-Joly E, Montagné N. Behavioral effect of plant volatiles binding to Spodoptera littoralis larval odorant receptors. Front Behav Neurosci. 2018;12:264. doi:https://doi.org/10.3389/fnbeh.2018.00264.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  105. Servedio MR. The role of linkage disequilibrium in the evolution of premating isolation. Heredity (Edinb). 2009;102:51–6. doi:https://doi.org/10.1038/hdy.2008.98.

    Article  CAS  Google Scholar 

  106. Smadja C, Butlin RK. On the scent of speciation: the chemosensory system and its role in premating isolation. Heredity (Edinb). 2009;102:77–97. doi:https://doi.org/10.1038/hdy.2008.55.

    Article  CAS  Google Scholar 

  107. Seimiya M, Gehring WJ. The Drosophila homeobox gene optix is capable of inducing ectopic eyes by an eyeless-independent mechanism. Development. 2000;127:1879–86. doi:https://doi.org/10.1242/dev.127.9.1879.

    Article  CAS  PubMed  Google Scholar 

  108. Reed RD, Papa R, Martin A, Hines HM, Counterman BA, Pardo-Diaz C, et al. optix drives the repeated convergent evolution of butterfly wing pattern mimicry. Science. 2011;333:1137–41. doi:https://doi.org/10.1126/science.1208227.

    Article  CAS  PubMed  Google Scholar 

  109. Patel NH, Martin-Blanco E, Coleman KG, Poole SJ, Ellis MC, Kornberg TB, et al. Expression of engrailed proteins in arthropods, annelids, and chordates. Cell. 1989;58:955–68. doi:https://doi.org/10.1016/0092-8674(89)90947-1.

    Article  CAS  PubMed  Google Scholar 

  110. Whitington PM, Meier T, King P. Segmentation, neurogenesis and formation of early axonal pathways in the centipede, Ethmostigmus rubripes (Brandt). Roux’s Arch Dev Biol. 1991;199:349–63. doi:https://doi.org/10.1007/BF01705928.

    Article  Google Scholar 

  111. Condron BG, Patel NH, Zinn K. Engrailed controls glial/neuronal cell fate decisions at the midline of the central nervous system. Neuron. 1994;13:541–54. doi:https://doi.org/10.1016/0896-6273(94)90024-8.

    Article  CAS  PubMed  Google Scholar 

  112. Keys DN, Lewis DL, Selegue JE, Pearson BJ, Goodrich L V., Johnson RL, et al. Recruitment of a hedgehog regulatory circuit in butterfly eyespot evolution. Science. 1999;283:532–4. doi:https://doi.org/10.1126/science.283.5401.532.

    Article  CAS  PubMed  Google Scholar 

  113. Condamin M. Monographie de genre Bicyclus (Lepidoptera Satyridae). Dakar: IFAN; 1973.

    Google Scholar 

  114. Kodandaramaiah U, Lees DC, Müller CJ, Torres E, Karanth KP, Wahlberg N. Phylogenetics and biogeography of a spectacular Old World radiation of butterflies: the subtribe Mycalesina (Lepidoptera: Nymphalidae: Satyrini). BMC Evol Biol. 2010;10:172. doi:https://doi.org/10.1186/1471-2148-10-172.

    Article  PubMed  PubMed Central  Google Scholar 

  115. Van’t Hof AE, Zwaan BJ, Saccheri IJ, Daly D, Bot ANM, Brakefield PM. Characterization of 28 microsatellite loci for the butterfly Bicyclus anynana. Mol Ecol Notes. 2005;5:169–72. doi:https://doi.org/10.1111/j.1471-8286.2004.00870.x.

    Article  Google Scholar 

  116. Prakash A, Monteiro A. Doublesex mediates the development of sex-specific pheromone organs in Bicyclus butterflies via multiple mechanisms. Mol Biol Evol. 2020;37:1694–707. doi:https://doi.org/10.1093/molbev/msaa039.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  117. Conesa A, Götz S, García-Gómez JM, Terol J, Talón M, Robles M. Blast2GO: A universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics. 2005;21:3674–6. doi:https://doi.org/10.1093/bioinformatics/bti610.

    Article  CAS  PubMed  Google Scholar 

  118. Altschul SF, Gish W, Miller W, Myers EW, Lipman DJ. Basic local alignment search tool. J Mol Biol. 1990;215:403–10. doi:https://doi.org/10.1016/S0022-2836(05)80360-2.

    Article  CAS  PubMed  Google Scholar 

  119. Jones P, Binns D, Chang HY, Fraser M, Li W, McAnulla C, et al. InterProScan 5: genome-scale protein function classification. Bioinformatics. 2014;30:1236–40. doi:https://doi.org/10.1093/bioinformatics/btu031.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  120. Bolger AM, Lohse M, Usadel B. Trimmomatic: a flexible trimmer for Illumina sequence data. Bioinformatics. 2014;30:2114–20. doi:https://doi.org/10.1093/bioinformatics/btu170.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  121. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29:15–21. doi:https://doi.org/10.1093/bioinformatics/bts635.

    Article  CAS  PubMed  Google Scholar 

  122. Anders S, Pyl PT, Huber W. HTSeq–a Python framework to work with high-throughput sequencing data. Bioinformatics. 2015;31:166–9. doi:https://doi.org/10.1093/bioinformatics/btu638.

    Article  CAS  PubMed  Google Scholar 

  123. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2. Genome Biol. 2014;15:550. https://doi.org/10.1186/s13059-014-0550-8.124.

    Article  PubMed  PubMed Central  Google Scholar 

  124. Zhu A, Ibrahim JG, Love MI. Heavy-tailed prior distributions for sequence count data: removing the noise and preserving large differences. Bioinformatics. 2019;35:2084–92. https://doi.org/10.1093/bioinformatics/bty895.

    Article  CAS  PubMed  Google Scholar 

  125. Benjamini Y, Hochberg Y. Controlling the false discovery rate: a practical and powerful approach to multiple testing. J R Stat Soc Ser B. 1995;57:289–300. doi: https://doi.org/10.1111/j.2517-6161.1995.tb02031.x.

    Article  Google Scholar 

  126. Supek F, Bošnjak M, Škunca N, Šmuc T. REVIGO summarizes and visualizes long lists of gene ontology terms. PLoS One. 2011;6:e21800. doi:https://doi.org/10.1371/journal.pone.0021800.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  127. Marchler-Bauer A, Bryant SH. CD-Search: protein domain annotations on the fly. Nucleic Acids Res. 2004;32:W327–31. doi:https://doi.org/10.1093/nar/gkh454.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  128. Vogt RG, Miller NE, Litvack R, Fandino RA, Sparks J, Staples J, et al. The insect SNMP gene family. Insect Biochem Mol Biol. 2009;39:448–56. doi:https://doi.org/10.1016/j.ibmb.2009.03.007.

    Article  CAS  PubMed  Google Scholar 

  129. Saenko SV, Marialva MS, Beldade P. Involvement of the conserved Hox gene Antennapedia in the development and evolution of a novel trait. Evodevo. 2011;2:9. doi:https://doi.org/10.1186/2041-9139-2-9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  130. Matsuoka Y, Monteiro A. Hox genes are essential for the development of eyespots in Bicyclus anynana butterflies. Genetics. 2021;217:iyaa005. doi:https://doi.org/10.1093/genetics/iyaa005

    Article  Google Scholar 

  131. Prakash A, Monteiro A. apterous A specifies dorsal wing patterns and sexual traits in butterflies. Proc R Soc B Biol Sci. 2018;285:20172685. doi:https://doi.org/10.1098/rspb.2017.2685.

    Article  CAS  Google Scholar 

  132. Özsu N, Chan QY, Chen B, Gupta MD, Monteiro A. Wingless is a positive regulator of eyespot color patterns in Bicyclus anynana butterflies. Dev Biol. 2017;429:177–85. doi:https://doi.org/10.1016/j.ydbio.2017.06.030.

    Article  CAS  PubMed  Google Scholar 

  133. Monteiro A, Prudic KL. Multiple approaches to study color pattern evolution in butterflies. Trends Evol Biol. 2010;2:e2. doi: https://doi.org/10.4081/eb.2010.e2.

    Article  Google Scholar 

  134. Connahs H, Tlili S, van Creij J, Loo TYJ, Banerjee TD, Saunders TE, et al. Activation of butterfly eyespots by Distal-less is consistent with a reaction-diffusion process. Development. 2019;146:dev169367. doi:https://doi.org/10.1242/dev.169367.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  135. Monteiro A, Chen B, Ramos DM, Oliver JC, Tong X, Guo M, et al. Distal-less regulates eyespot patterns and melanization in Bicyclus butterflies. J Exp Zool Part B Mol Dev Evol. 2013;320:321–31. doi:https://doi.org/10.1002/jez.b.22503.

    Article  CAS  Google Scholar 

  136. Bhardwaj S, Prudic KL, Bear A, Dasgupta M, Wasik BR, Tong X, et al. Sex differences in 20-hydroxyecdysone hormone levels control sexual dimorphism in Bicyclus anynana wing patterns. Mol Biol Evol. 2018;35:465–72. doi:https://doi.org/10.1093/molbev/msx301.

    Article  CAS  PubMed  Google Scholar 

  137. Monteiro A, Glaser G, Stockslager S, Glansdorp N, Ramos D. Comparative insights into questions of lepidopteran wing pattern homology. BMC Dev Biol. 2006;6:52. doi:https://doi.org/10.1186/1471-213X-6-52.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  138. Beldade P, Peralta CM. Developmental and evolutionary mechanisms shaping butterfly eyespots. Curr Opin Insect Sci. 2017;19:22–9. doi:https://doi.org/10.1016/j.cois.2016.10.006.

    Article  PubMed  Google Scholar 

  139. Westerman EL, VanKuren NW, Massardo D, Tenger-Trolander A, Zhang W, Hill RI, et al. Aristaless controls butterfly wing color variation used in mimicry and mate choice. Curr Biol. 2018;28:3469–-3474.e4. doi:https://doi.org/10.1016/j.cub.2018.08.051.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  140. Woronik A, Tunström K, Perry MW, Neethiraj R, Stefanescu C, Celorio-Mancera MP, et al. A transposable element insertion is associated with an alternative life history strategy. Nat Commun. 2019;10:5757. doi:https://doi.org/10.1038/s41467-019-13596-2.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  141. Nadeau NJ, Pardo-Diaz C, Whibley A, Supple MA, Saenko SV, Wallbank RWR, et al. The gene cortex controls mimicry and crypsis in butterflies and moths. Nature. 2016;534:106–10. doi:https://doi.org/10.1038/nature17961.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We would like to thank Jacob Siebenmorgen for help with starting B. anynana families for this experiment, Elizabeth Ruck, Amanda Scholes, and Tara Stuecker for their advice on molecular techniques, and Grace Hirzel, Matt Murphy, Sushant Potdar, and Nikki Robertson for assistance with butterfly husbandry.

Funding

This research was funded by the University of Arkansas and NSF IOS grant #1937201 to ELW.

Author information

Authors and Affiliations

Authors

Contributions

DAE and ELW designed this experiment and wrote the manuscript. DAE performed the research and analyzed the data. Both authors read and approved the final manuscript.

Corresponding author

Correspondence to David A. Ernst.

Ethics declarations

Ethics approval and consent to participate

No ethics approval was required for this study.

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

Additional file 1:

Supplemental Tables S1 and S25; captions for Supplemental Tables S2-S24; Supplemental Figures S1-S9.

Additional file 2: Supplemental Tables S2-S24

.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons.org/licenses/by/4.0/. 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 in a credit line to the data.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Ernst, D.A., Westerman, E.L. Stage- and sex-specific transcriptome analyses reveal distinctive sensory gene expression patterns in a butterfly. BMC Genomics 22, 584 (2021). https://doi.org/10.1186/s12864-021-07819-4

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12864-021-07819-4

Keywords