Transcriptome dynamics during metamorphosis of imaginal discs into wings and thoracic dorsum in Apis mellifera castes

Abstract

Background: Much of the complex anatomy of a holometabolous insect is built from disc-shaped epithelial structures found inside the larva, i.e., the imaginal discs, which undergo a rapid differentiation during metamorphosis. Imaginal discs-derived structures, like wings, are built through the action of genes under precise regulation.

Results: We analyzed 30 honeybee transcriptomes in the search for the gene expression needed for wings and thoracic dorsum construction from the larval wing discs primordia. Analyses were carried out before, during, and after the metamorphic molt and using worker and queen castes. Our RNA-seq libraries revealed 13,202 genes, representing 86.2% of the honeybee annotated genes. Gene Ontology analysis revealed functional terms that were caste-specific or shared by workers and queens. Genes expressed in wing discs and descendant structures showed differential expression profiles dynamics in premetamorphic, metamorphic and postmetamorphic developmental phases, and also between castes. At the metamorphic molt, when ecdysteroids peak, the wing buds of workers showed maximal gene upregulation comparatively to queens, thus underscoring differences in gene expression between castes at the height of the larval-pupal transition. Analysis of small RNA libraries of wing buds allowed us to build miRNA-mRNA interaction networks to predict the regulation of genes expressed during wing discs development.

Conclusion: Together, these data reveal gene expression dynamics leading to wings and thoracic dorsum formation from the wing discs, besides highlighting caste-specific differences during wing discs metamorphosis.

Keywords: Apis mellifera; Caste differential expression; Honeybee; Metamorphosis; RNA-seq; Wing imaginal discs; miRNA.

Fig. 1

Localization of imaginal wing discs in larvae and their descendants in pupae. (1) Imaginal wing discs localization in fifth instar feeding larvae (L5F1). (2) Imaginal wing discs as dissected from an L5F1 larvae. (3) The imaginal wing discs descendants, wings and thoracic dorsum, dissected from a newly-ecdysed pupae (PW phase)

Fig. 2

Multidimensional scaling (MDS) plot showing the relationships (similarity/dissimilarity) between the RNA-seq libraries of wing imaginal discs and their developing descendants, i.e., thoracic dorsum and wings. Worker (W) and queen (Q) samples are indicated for each developmental phase (at the right) from which we extracted the imaginal discs and descendants. L5F1 phase (undifferentiated imaginal discs), L5PP2 phase (metamorphosing imaginal discs), PW phase (pupal wing/thoracic dorsum), PB phase (immature adult wing/thoracic dorsum), PBD phase (mature adult wing/thoracic dorsum). Independent biological samples of a same developmental phase are indicated by _0, _1, and _2

Fig. 3

Dendrogram representing the Euclidean distance analysis grouping the DEGs according the developmental phase (L5F1, L5PP2, PW, PB, or PBD). The letters W and Q ahead of each developmental phase acronym represent workers and queens, respectively. Independent biological samples of a same developmental phase are indicated by _0, _1, and _2. AU cluster supporting values are in red. Branch BP supporting values (10,000 replications) are in green. Edge values are in grey. The small table on the right shows the AU and BP values (corresponding to edges 1, 5, 6, 7, 11, 12, 13, and 21) that could not be indicated in the dendrogram due to insufficient space

Fig. 4

Number of DEGs co-expressed or exclusively expressed in wing discs and their descendants through developmental phases of (A) workers and (B) queens

Fig. 5

Differential gene expression during imaginal disc development to wings/thoracic dorsum in workers and queens. A Volcano graphs represent the upregulated and downregulated DEGs (red dots at the right and left of the volcano graphs, respectively) among the total of genes (blue plus red dots) identified in the RNA-seq analysis. The horizontal line in each volcano graph indicates the adjusted p-value (0.05). B Percentage of upregulated genes in workers or queens. (C) Gene Ontology (GO) functional analysis of DEGs displaying similar expression in workers and queens (see Additional file 5: Supplementary Table 2). (D) GO analysis of DEGs upregulated in the developing wing discs of workers or queens (see Additional file 6: Supplementary Table 3)

Fig. 6

Expression profiles of (A) genes related to ecdysis and (B) genes encoding chitin- and cuticle pigmentation/sclerotization-related enzymes in the wing discs and developing wings/thoracic dorsum of workers and queens in the different developmental phase

Fig. 7

Expression profiles of (A) cuticular protein genes and (B) wing-patterning genes in workers and queens in the different developmental phases

Fig. 8

miRNAs-targets network. (A) Targets of 22 ame-miRNAs whose homologs in other insect species have been related to metamorphosis, larval-to-pupal transition, and ecdysteroid control of metamorphosis. (B) Ten ame-miRNAs among the 22 ones putatively regulating genes with described functional activities in the honeybee. The red circles represent ame-miRNAs (the darker the color, the greater the number of reads). The blue circles represent the downregulated genes in the wing discs of L5PP2 workers