E93 is indispensable for reproduction in ametabolous and hemimetabolous insects Free

The ecdysone-induced protein E93 gene (E93) encodes a transcription factor containing one of two highly conserved helix-turn-helix DNA-binding domains (RHFs) of the Mblk-1 family (Baehrecke and Thummel, 1995; Takayanagi-Kiya et al., 2017). Insect evolution involves transitions from ametaboly to hemimetaboly, and ultimately to holometaboly (Truman and Riddiford, 2019; Belles, 2020); E93 plays a crucial role as the ‘adult specifier’ that triggers metamorphosis in both hemimetabolous and holometabolous insects (Ureña et al., 2014; Belles and Santos, 2014; Chafino et al., 2019). Recently, several studies have shown that E93 is also important for the early embryonic development of ametabolous and hemimetabolous insects (Fernandez-Nicolas et al., 2023), and for female reproduction in hemimetabolous and holometabolous insects. Systemic RNAi knockdown of E93 has been found to decrease the synthesis of vitellogenin (Vg) (vitellogenesis) and the number of eggs laid by the brown planthopper (Mao et al., 2019), red flour beetle (Eid et al., 2020), mosquito (Wang et al., 2021) and migratory locust (Liu et al., 2022). Particularly in adult female mosquitoes, E93 determines reproductive cycle switches by regulating basal metabolism and autophagy in the fat body (FB), which is the major tissue for producing Vg precursors (Wang et al., 2021). However, it remains unknown whether E93 is also essential for reproduction in ametabolous insects and whether E93 directly regulates ovarian development in insects. A comprehensive study of the effect of E93 on ametabolous insect reproduction compared with that of hemimetabolous and holometabolous insects will provide a notable point of interest in evolutionary biology research.

In insects, the maturation of panoistic ovarian follicles, which are composed of oocytes and their surrounding follicular cells, progresses through three major phases: previtellogenesis, vitellogenesis and choriogenesis (Swevers and Iatrou, 2003). The firebrat Thermobia domestica (Zygentoma), the closest extant relative to winged insects (Pterygota), is considered a model organism for ametabolous insects. The gonadotrophic cycles in female T. domestica have been examined previously (Bitsch et al., 1985; Fig. S1A). The German cockroach Blattella germanica (Blattodea) is one of the best-studied hemimetabolous insects with well-characterized gonadotrophic cycles (Martín et al., 1995; Fig. S1B). T. domestica and B. germanica can be used as appropriate model organisms of ametabolous and hemimetabolous insects, respectively, for evolutionary developmental biology (evo-devo) studies (Fernandez-Nicolas et al., 2023), especially those focused on the field of reproductive biology.

Ecdysone receptor (EcR)-transduced 20-hydroxyecdysone (20E) signalling is required for the survival, proliferation and polarity of follicular cells, and thus for ovarian follicle development during previtellogenesis in some hemimetabolous and holometabolous insects (Cruz et al., 2006; Romani et al., 2009; Belles and Piulachs, 2015). 20E signalling is also required for vitellogenesis in some holometabolous insects (Roy et al., 2018), whereas it is indispensable for choriogenesis in most, if not all, insects (Bellés et al., 1993; Papantonis et al., 2015; Zhao et al., 2023). Thus, it would be of great interest to investigate how E93, a 20E primary response gene, regulates insect reproduction from an evolutionary developmental biology perspective.

We have previously reported some specific roles of E93 during metamorphosis in the fruit fly Drosophila melanogaster and the silkworm Bombyx mori (Liu et al., 2014, 2015; Zhang et al., 2023). In an initial investigation, we observed tissue-specific E93 expression in the previtellogenic ovaries (Ov) of T. domestica females, and its expression exhibited a distinct cyclical pattern throughout the gonadotrophic cycle. Subsequently, by integrating mRNA fluorescence in situ hybridization (mRNA-FISH), CRISPR/Cas9, RNAi, RNA sequencing (RNA-seq) and other methods, we demonstrated the indispensable and conserved role of E93 in promoting ovarian follicle development during previtellogenesis in both ametabolous and hemimetabolous insects. In addition, we discovered that E93 acts as a crucial factor in regulating reproductive cycle switches in the adult B. germanica female fat body, which remains conserved in holometabolous insects. By integrating the existing E93 expression and functional data, we present comprehensive insight into the evolution of E93 during insect evolution, which might be helpful for elucidating the mechanisms underlying the origin and evolution of insect metamorphosis.

To understand the roles of E93 in T. domestica, we first examined the spatiotemporal expression pattern of E93 during postembryonic development using quantitative real-time PCR (qPCR) or RNA-seq. To clarify the developmental profile of E93, the abdomen (where most of the reproductive system is located) was sampled, and the firebrats used ranged from young nymphs without sexual traits, including N4D2 (2-day-old in the number 4 instar), N4D4, N6D2 and N6D6 nymphs, to penultimate nymphs (N8D2 and N8D4, when females and males are distinguished by the presence of the ovipositor), through adult stages (N10D2). Surprisingly, the number of E93 transcripts in the abdomen was significantly greater in females than in males from N8 onwards. Moreover, tissue distribution analyses during the first previtellogenesis stage (N9, the last nymphal instar) revealed that the mRNA level of E93 in the Ov was significantly greater than that in other tissues (including the FB) (Fig. 1B). The spatiotemporal expression pattern of E93 in T. domestica suggested that E93 might be involved in the regulation of female reproduction in ametabolous insects.

We then conducted a precise time-course transcription analysis of E93 in the Ovs and FBs (as negative controls) of T. domestica females during the last three nymphal instars (N7, N8 and N9) and the first two adult moulting cycles (N10 and N11) using qPCR. Ov and FB samples were collected from N7D4, N8D4, N9D6, N10D1, N10D3, N10D5, N10D7, N10D9, N11D1, N11D4 and N11D8 firebrats. The E93 mRNA abundance in the Ov was much greater than that in the FB, and the ovarian E93 transcript level during previtellogenesis was relatively greater than that during vitellogenesis (N10D1–D5; N11D1-D5) (Fig. 1C). To investigate the responsiveness of E93 to 20E in T. domestica, we performed 20E application experiments at N9D4. The Ovs were isolated for qPCR analysis to detect E93 and E75 (a canonical 20E primary response gene) after 12 h of treatment. 20E significantly induced the expression of these two genes; however, the induction efficiency of E93 was much lower than that of E75 (Fig. 1D).

On this basis, we employed mRNA-FISH to locate E93 mRNAs in ovarioles from the previtellogenic period at N9D3, N9D6, N10D6 and N10D9, and from the vitellogenic period at N10D3. During the first gonadotrophic cycle (N9-N10D5) in T. domestica, the basal 1-3 oocytes mature, and the basal 2-4 oocytes mature during the second gonadotrophic cycle (N10D6-N11D5). We classified and named the ovarian follicles of ovarioles as basal 1st ovarian follicle (BOF1), basal 2nd ovarian follicle (BOF2)…basal 5th ovarian follicle (BOF5), etc., in T. domestica (Fig. 1E). Ovarian follicles that are prepared to mature in a specific gonadotrophic cycle enlarge during previtellogenesis and incorporate Vg during subsequent vitellogenesis. Consistent with the qPCR results, ovarian E93 was cyclically highly expressed during previtellogenesis. However, during vitellogenesis, its transcript was nearly undetectable (Fig. 1E-N). Specifically, E93 mRNAs were detected mainly in the follicular cells surrounding all the oocytes at the early stage of previtellogenesis (N9D3, Fig. 1E′-E⁗; N10D6, Fig. 1K′-K⁗) or in the follicular cells of the oocytes that were not ready to mature at the late stage of previtellogenesis (BOF5; N9D6, Fig. 1G‴-G⁗; N10D9, Fig. 1M‴-M⁗). Notably, there was also an antisense signal in the cytoplasm and nucleus of the oocyte, although to a lesser extent than in the follicular cells. However, at the late stage of previtellogenesis, E93 mRNAs were located mainly in oocytes that were enlarged and nearly able to incorporate Vg (BOF1, BOF2; Fig. 1G′,G″; Fig. 1M′,M″). Overall, when induced by 20E, E93 is predominantly expressed in previtellogenic T. domestica ovarian follicles, suggesting that E93 might play a spatiotemporal-specific role in regulating ovarian development in ametabolous insects.

To evaluate the potential role of E93 in regulating female reproduction in ametabolous insects, we performed E93 loss-of-function experiments in T. domestica. By using CRISPR/Cas9-mediated genome editing, we generated a homozygous E93 mutant strain characterized by a 103 bp deletion within the conserved RHF2 domain (Fig. 2A). In the G0 generation of the E93 mutant, ∼5% of the newly hatched nymphs exhibited abnormal embryonic segmentation (Fig. S2A,B). Because the E93 homozygous T. domestica mutant population yielded almost no offspring, we used the heterozygous population to generate homozygous mutants for all subsequent experiments. In T. domestica, the scales first appear upon ecdysis to the fourth nymphal instar (N4). From a phenotypic perspective, E93 homozygous mutants (E93^−/−) were distinguishable from wild-type (WT) firebrats, and heterozygous mutants (E93^+/−) with smaller body sizes and diminished melanin on scales from N6 onwards [Fig. 2B, also see examples of moulting obstructions (Fig. 2B′)] were confirmed via DNA sequencing. In other words, we were unable to distinguish homozygotes from the offspring of heterozygous populations by phenotype until N6 (Fig. S2C-G′). The hatching rate and survival rate (N1-N5) of the offspring of the heterozygous population were significantly lower than those of the offspring of the WT population (Fig. S2H). Next, the survival rate, body size and instar duration of E93^−/− individuals were compared statistically with those of their E93^+/−/WT counterparts after N6 moulting. As shown in Fig. 2C, the survival rate of the E93^−/− mutants was significantly lower than that of the E93^+/−/WT controls within the ∼60 days after N6 (eight moults experienced by the insects during these 60 days); ∼4% per instar E93^−/− experienced unsuccessful moulting after N6. We continuously monitored the body size of the same batch of individuals from N6 to N11 by taking photographs and calculating their body area on the second day of each instar. Notably, the body size of E93^−/− was significantly reduced compared with that of E93^+/−/WT controls at each instar (Fig. 2D). In addition, nymphal moulting was delayed by a few days in the E93 homozygous mutants compared with the control strains (Fig. 2E). Furthermore, we validated the phenotypes resulting from E93 knockout through RNAi. RNAi experiments showed that sustained dsRNA delivery (injection once per instar, at least twice) yielded a melanin-diminishing phenotype on scales similar to that of the mutants (Fig. S3). The composite data revealed that low-dose E93 consistently plays a specific role in cuticle development and moulting throughout the life cycle in ametabolous T. domestica.

We next focused on how the E93 mutation affects female reproduction in T. domestica. We tested four mating combinations (E93^−/− ♀×E93^−/− ♂; E93^−/− ♀×WT ♂; E93^−/− ♂×WT ♀; and WT ♀×WT ♂) to assess the fecundity of E93^−/− across the first three gonadotrophic cycles (N10, N11 and N12). The results showed that the number of eggs was dramatically lower in E93^−/− ♀ firebrats than in their WT ♀ counterparts (Fig. 2F). Notably, the number of eggs laid by WT females crossed with E93^−/− males was also significantly reduced, although the reproductive system of the male was not affected in the E93^−/− mutants (Fig. S4). These statistical findings show that E93 loss of function results in severe fecundity deficiency in female T. domestica, which aligns with E93 spatiotemporal expression.

The analyses of both E93 expression patterns and E93^−/− phenotypes suggest that E93 is indispensable for ovarian development in T. domestica. To elucidate the impact of E93 knockout on ovarian development, we continuously observed morphological changes in Ovs from N9D3 to N12D1, which encompasses the first two gonadotrophic cycles, in both WT and E93^−/−T. domestica individuals. As shown in Figs S1A, S5A and Fig. 3A-E, the first gonadotrophic cycle starts with previtellogenesis (the basal 1-3 follicles in each ovariole enlarge, which prepares for Vg accumulation, characterized by transparent ovarioles) during the last nymphal instar, N9. This process is followed by vitellogenesis (during which the fully enlarged follicles incorporate Vg, sharply swell, and move towards maturation and eventual extrusion, which are marked by milk-white ovarioles), which continues until the egg lays on day 6 of N10. After oviposition, the second gonadotrophic cycle commences, and the basal 2-4 follicles in each of the ovarioles enlarge and become well prepared for Vg accumulation during subsequent vitellogenesis (N10D6-N10D10). After ecdysis, these previtellogenic follicles gradually mature through Vg incorporation until extrusion (N11D1-N11D6), and so forth. However, the ovarian follicles of the E93^−/− did not undergo one successful gonadotrophic cycle during the corresponding time period (Fig. S5B and Fig. 3F-J). The basal follicles were hardly enlarged during previtellogenesis, and only a small proportion of the follicles could partially accumulate yolk during N11 in E93^−/−. As demonstrated by mRNA-FISH analyses, the E93 transcripts were significantly present in the ovarioles of WT individuals but absent in the ovarioles of homozygous mutants (E93^−/−) at both N9D6 (Fig. 3K,L) and N10D8 (Fig. 3N,O), which represent previtellogenesis of the first and second gonadotrophic cycles, respectively. Morphological analysis of Ovs indicated that E93 is essential for previtellogenic follicle development in T. domestica.

We then conducted comparative transcriptome analyses to explore the molecular mechanism underlying the lack of previtellogenic follicle development in E93^−/−. Previtellogenic Ovs were collected from E93^−/− and WT female firebrats at two time points, N9D6 and N10D8, and vitellogenic Ovs were collected from WT female firebrats at N10D3. Through a comparative analysis of N9D6-Ov RNA-seq data from E93^−/− and WT individuals, we identified 284 upregulated genes and 308 downregulated genes (Fig. S6A). KEGG (Kyoto Encyclopedia of Genes and Genomes) enrichment analysis indicated that only one pathway, extracellular matrix (ECM)-receptor interaction, was significantly downregulated in the E93^−/− Ovs compared with the WT Ovs, with an adjusted P-value<0.05 (Fig. 3M). Moreover, the genes upregulated in the E93^−/− Ovs were significantly enriched in pathways responsible for metabolism, with the cytochrome P450 pathway being the most enriched pathway (Fig. S7A). In another comparative RNA-seq analysis of the Ovs at N10D8, 809 genes were significantly upregulated, whereas 989 were downregulated in E93^−/− compared with the WT control (Fig. S6B). The ECM-receptor interaction pathway was also among the pathways enriched in the downregulated genes; in addition, cytoskeleton proteins and cell-adhesion molecules were also enriched (Fig. 3P). Similarly, the cytochrome P450 pathway was the top pathway enriched in the significantly upregulated genes (Fig. S7B). Moreover, the ECM-receptor interaction pathway was again enriched in the downregulated gene set generated from the comparative analysis of the ovarian transcriptome: N10D3 (vitellogenesis) versus N10D8 (previtellogenesis) (Fig. S7D-F). The involved ECM genes, including the synaptic vesicle glycoprotein, laminin, agrin and collagen genes, are labelled on volcano plots and were validated by qPCR (Fig. S6A′,B′). The results of the T. domestica comparative transcriptome analyses suggested that E93 transcriptionally activates ECM-receptor interactions that promote previtellogenic follicle development in ametabolous insects.

Next, we aimed to determine whether the essential role of E93 observed in the previtellogenic Ovs of firebrats is a shared characteristic among hemimetabolous insects, which evolved from ametabolous ancestors. We used the cockroach B. germanica as a model of hemimetabolous insects. Moreover, the tissue distribution pattern of E93 expression at N6D9 (the mid-late stage of the last nymphal instar, which is the first previtellogenic period, Irles and Piulachs, 2014) was explored by RNA-seq. The results showed that E93 was most highly expressed in the wing, followed by the brain, cercus and Ov (Fig. 4A). Using qPCR, we analysed the developmental profile of E93 expression in female Ovs and FBs (as a control) during N4, N5, N6 and the initial two gonadotrophic cycles. In contrast to firebrats, in which the TdE93 transcript level in FBs was much lower than that in Ovs (Fig. 1C), the relative expression level of BgE93 in FBs was lower but of the same order of magnitude as that in Ovs (Fig. 4B). Importantly, we found that the E93 mRNA abundance in both Ovs and FBs was greater during previtellogenesis but significantly lower during vitellogenesis in B. germanica (Fig. 4B). Similarly, E93 gene expression was significantly upregulated in the Ovs of B. germanica after 20E treatment (Fig. 4C). In B. germanica, only the basal ovarian follicle (BOF) matures in the same gonadotrophic cycle, whereas the sub-basal ovarian follicle (sBOF), which grows only slowly in each gonadotrophic cycle (Alborzi and Piulachs, 2023), can be equivalent to BOF5 during the first two gonadotrophic cycles in T. domestica. mRNA-FISH analysis of the ovarioles from N6D3 (early-), N6D9 (middle-), AD2 (late-previtellogenesis) and AD5 (vitellogenesis) in B. germanica also revealed that the expression level of ovarian E93 was much greater during previtellogenesis than during vitellogenesis (Fig. 4D-K). Additionally, E93 transcripts were present mainly in the follicular cells of whole ovarioles during the entire previtellogenic period in B. germanica (Fig. 4D-D⁗,F-F⁗,H-H⁗), which is generally consistent with the results observed in T. domestica. Nevertheless, there are some differences between the two species: (1) BgE93 mRNAs, which are still located in the follicular cells surrounding the BOF, were not detectable in the oocytes at the late stage of previtellogenesis (Fig. 4H′-H″,Fig. 1G′-G″, Fig. 1M′-M″); and (2) BgE93 mRNAs were unevenly distributed in spaces between the oocyte membrane and follicular epithelia, and the spaces probably resulted from the expansion of the ECM of follicular cells. Our spatiotemporal analyses of E93 expression in the two insect species suggested that the physiological function of E93 might be conserved in the Ov during previtellogenesis but also highlighted a potentially enhancing role for E93 in the FB during the transition from ametaboly to hemimetaboly.

RNA interference and comparative transcriptome experiments were employed to investigate the molecular mechanism by which E93 regulates previtellogenic follicle development in B. germanica. The double-stranded RNA (dsRNA) of E93 was applied at N6D1, and during N6 (previtellogenesis), BOF development was impeded after E93 knockdown, as the average BOF length of RNAi-BgE93 cockroaches was significantly lower than that of RNAi-MusLta controls (Fig. 4L-P). The previtellogenic Ovs were dissected at N6D9 from E93-depleted (RNAi-BgE93) and control (RNAi-MusLta) cockroaches for RNA-seq analysis. Notably, the ECM-receptor interaction and cytochrome P450 pathways were also enriched in the 272 downregulated and 55 upregulated gene sets (Fig. 4Q,R; Fig. S8A), respectively, and the involved ECM genes, such as collagen, were labelled on a volcano plot with validation using qPCR (Fig. S6C,C′). However, many more E93-regulated pathways were found in the Ovs of B. germanica, such as glycan degradation, amino sugar and nucleotide sugar metabolism, and glycosphingolipid biosynthesis (Fig. 4R). In addition, the ECM-receptor interaction pathway was also enriched in the downregulated gene set generated from the comparative transcriptome analysis of AD6 Ovs (vitellogenesis) versus AD1 Ovs (previtellogenesis) (Fig. S8B-D). Based on the above findings, our evolutionary development studies suggest that E93 plays a relatively conserved role in promoting previtellogenic follicle development in both ametabolous and hemimetabolous insects.

We investigated the roles of E93 in the B. germanica FB in view of the obviously cyclical expression pattern of E93 during gonadotrophic cycle switches. We first conducted a time-course transcription analysis of Vg, a marker gene of vitellogenesis, in the FB using qPCR (coupled with E93). FB samples were collected from adult females at various stages during the first gonadotrophic cycle: the previtellogenic phase at AD1, AD2 and AD3; the vitellogenic phase at AD4, AD5, AD6 and AD7; the choriogenesis phase at AD8; and the extrusion phase at AD9 (the beginning of the second gonadotrophic cycle, i.e. previtellogenesis). The Vg mRNA abundance increased slowly from AD1 to AD3 and increased dramatically from AD3 to AD5, reaching a peak at AD5. Subsequently, it decreased by AD8-AD9 to a low level similar to that at AD1. Conversely, the E93 transcript level was relatively higher from AD1 to AD3 and from AD8to AD9, but significantly lower from AD4 to AD7 (Fig. 5A). BgE93 dsRNAs were injected at AD0 (newly emerged adults, i.e. previtellogenesis) or AD6 (vitellogenesis). RNAi-E93, which began at AD0, resulted in a significant reduction in Vg expression in the female adult FB at AD2 and AD6; however, RNAi-E93, which started at AD6, resulted in significant upregulation of Vg expression at AD9 (Fig. 5B). In addition, the cockroaches treated with E93-dsRNA at AD0 exhibited delays in follicle maturation (vitellogenesis) and, eventually, a delay of at least 2 days in oviposition compared with those in control individuals, with no significant E93 knockdown in the Ovs (Fig. 5C-F). The analyses of E93 and Vg gene expression/interaction suggested that E93 is required for the regulation of ovarian maturation during vitellogenesis and plays a key role in determining gonadotrophic cyclicity in the fat body of hemimetabolous insects.

To understand how E93 regulates the gonadotrophic cycle in FBs, we conducted transcriptome analysis. The cockroaches were injected with dsRNAs at AD0 or AD6, and the FBs were dissected from the treated insects at AD6 or AD9 for RNA-seq. We first conducted three groups of differential expression analyses, previtellogenic FB (RNAi-MusLta, AD9) versus vitellogenic FB (RNAi-MusLta, AD6), AD6 FB (RNAi-BgE93 versus RNAi-MusLta) and AD9 FB (RNAi-BgE93 versus RNAi-MusLta), and identified three differentially expressed gene sets (638 downregulated and 560 upregulated genes; 117 downregulated and 331 upregulated genes; 383 downregulated and 264 upregulated genes) (Fig. S9A-C). HeatCluster analysis of the AD6 and AD9 FB samples using the above DEGs revealed that RNAi-E93 was not sufficient to completely switch the gonadotrophic cycle (Fig. 5G). Similarly, KEGG pathway enrichment analysis revealed that most pathways significantly upregulated in the differential analysis of vitellogenic FBs (RNAi-MusLta, AD6) versus previtellogenic FBs (RNAi-MusLta, AD9) were still upregulated in the vitellogenic FBs compared with the RNAi-E93-treated previtellogenic FBs (Fig. S9C-H). However, KEGG analysis also revealed that most of the upregulated pathways in previtellogenic FBs (RNAi-MusLta, AD9) compared with those in vitellogenic FBs (RNAi-MusLta, AD6), which are mainly involved in material metabolism (including drug/xenobiotics metabolism by P450 or other enzymes, ascorbate and aldarate metabolism, pentose and glucoronate interconversions, glycolysis/gluconeogenesis, and the pentose phosphate pathway), were also significantly upregulated in the transcriptomic analysis of RNAi-BgE93 versus RNAi-MusLta during vitellogenesis (AD6). Nevertheless, these pathways were significantly downregulated in the comparative transcriptome analysis of RNAi-BgE93 versus RNAi-MusLta during previtellogenesis (AD9) (Fig. 5H-J). Thus, the RNA-seq results demonstrated that E93 promotes the vitellogenesis-previtellogenesis switch in the FB of B. germanica by regulating genes involved in metabolism.

Despite the relatively low expression level of E93, we conducted an essential RNA-seq analysis on T. domestica FBs. The Ovs of very few E93^−/− individuals could enter the vitellogenic phase. Previtellogenic FBs from E93^−/− and WT female firebrats were collected at N10D6, and vitellogenic FBs were collected from WT female firebrats at N10D3 for RNA-seq. Differential gene expression and KEGG analyses indicated that there was no significant difference in the E93 transcription level between previtellogenic FBs and vitellogenic FBs (Fig. S10A-C). In addition, comparative analysis of the FB transcriptome between E93^−/− and WT (previtellogenesis) strains revealed that the E93 mutation did not affect the expression of Vg; however, loss of function of E93 in the FB resulted in the upregulation of a set of pathways, such as metabolism by P450, retinol metabolism, steroid hormone biosynthesis and bile secretion, which were also significantly upregulated in the vitellogenic FB transcriptome compared with the previtellogenic FB transcriptome (Fig. S10C-F). Notably, there is no commonality in the regulatory mechanism of E93 in FB between T. domestica and B. germanica (Fig. S9B′,D). These results suggest that E93 plays a specific but weak role in regulating gonadotrophic cycles in the FB of ametabolous insects, highlighting the potentially enhancing role of E93 in the FB during the transition from ametaboly to hemimetaboly.

The role of E93 in reproduction is not well understood, especially in the ovary. Here, we have revealed an important regulatory role for E93 in ovarian follicles during previtellogenesis in T. domestica and B. germanica. The 20E titre in the final nymphal instar or early pupae can trigger previtellogenic development of the Ov in insects (Bitsch et al., 1985; Cruz et al., 2003; Belles and Piulachs, 2015; Ables et al., 2016; Zhang et al., 2017; Rumbo et al., 2023). Follicular cells proliferate intensively during previtellogenesis, reaching the optimal cell number required to provide oocytes with the desired shape and volume before yolk accumulation (Irles and Piulachs, 2014). The proliferation, survival and establishment of polarity in follicular cells require sex-specific 20E signalling (Romani et al., 2009; Morris and Spradling, 2012). In this study, we discovered that the knockout or knockdown of E93, a 20E-response gene, severely inhibited the previtellogenic development of ovarian follicles in T. domestica and B. germanica (Figs 1-4). Moreover, in both T. domestica and B. germanica, elimination or reduction of E93 expression resulted in significant downregulation of genes involved in the ECM-receptor interaction pathway (Figs 3 and 4). The ECM constitutes a crucial component of the oocyte microenvironment, playing a pivotal role in regulating follicular cell development, as it provides physical support to cells and serves as a ligand for cell surface receptors that activate a range of signalling pathways for regulating cell proliferation, differentiation, polarity and morphology (Berkholtz et al., 2006; Hastings et al., 2019; Nagyová et al., 2021). In addition to the ECM-receptor interaction pathway, the MAPK signalling and cytokine/growth factor pathways were enriched in the downregulated gene set when the E93^−/− and WT ovarian transcriptomes were compared; likewise, the PI3K-Akt signalling and cytokine/growth factor pathways were enriched in the downregulated gene set generated from comparison of the ovarian transcriptome at N10D3 (vitellogenesis) and N10D8 (previtellogenesis) (Fig. S7). Upregulated cytokine production is associated with increased expression and clustering of integrin (an ECM receptor) (Lowin and Straub, 2011), which is followed by the activation of signalling pathways such as the MAPK and PI3K-Akt signalling pathways, which promote cell proliferation and survival (Matsubayashi et al., 2004). Given these findings, it is plausible that E93 conservatively promotes follicular cell proliferation via the cytokine-ECM-receptor-MAPK/PI3K-Akt signalling pathway, thus promoting oocyte enlargement during previtellogenesis in both ametabolous and hemimetabolous insects.

E93 is a key factor regulating gonadotrophic cycles in the FB of the adult female mosquito Aedes aegypti (Wang et al., 2021). Here, we have found that E93 plays a conserved role in promoting the vitellogenesis-previtellogenesis switch by regulating metabolism and Vg production in the FBs of hemimetabolous B. germanica. The E93 transcript level was relatively higher after adult emergence (AD1-AD3, i.e. previtellogenesis) but significantly lower during vitellogenesis (AD4-AD7), and the E93 mRNA level was elevated once again at the end of vitellogenesis or at the beginning of the next previtellogenesis (AD8-AD9) (Fig. 5). RNAi-E93 at AD0 resulted in downregulation of Vg at AD2 and AD6. However, RNAi-E93 at AD6 resulted in downregulation of Vg at AD9. These two opposing transcriptional outputs controlled by E93 have been identified in D. melanogaster and A. aegypti (Uyehara et al., 2017; Wang et al., 2021), with E93 regulating its target gene expression by controlling chromatin accessibility. This also suggests that E93 might be involved in the regulation of the gonadotrophic cycle in B. germanica. E93 knockdown at AD0 delayed ovary development without downregulating E93 within the ovaries, indicating that E93 plays a key role in the adult FB (Fig. 5). RNAi-E93 at AD0 did not result in large-scale transcriptional changes at AD2 (105 differentially expressed genes). Our transcriptome-based analysis of the adult female FB at AD6 and AD9 indicated that an E93-dependent decrease/increase in the expression of genes involved in metabolism promotes the end of the reproductive cycle and progression into the next cycle, even if the efficacy is limited (Fig. 5). Nevertheless, the function of E93 in the adult female FB of T. domestica is far weaker than that in the B. germanica FB (Fig. S10).

Taken together (Fig. 6A), the results of this study reveal that E93 conservatively promotes the previtellogenic development of ovarian follicles in both ametabolous and hemimetabolous insects, and that its role in the FB of adult females is gradually enhanced to regulate gonadotrophic cycles during the transition from ametaboly to hemimetaboly.

E93 expression during early embryonic development was found to be high in ametabolous and hemimetabolous insects but low in holometabolous insects, and E93 is crucial for germ band development in B. germanica (Fernandez-Nicolas et al., 2023). We noticed that E93 might regulate embryonic segmentation in T. domestica (Fig. S2) and, as such, is most likely involved in early embryonic development in both ametabolous and hemimetabolous insects; however, it appears that this role of E93 has been weakened or even lost during the transition from hemimetaboly to holometaboly (Fig. 6B).

Insect metamorphosis is coordinately regulated by 20E and JH, and E93 expression is induced by 20E and inhibited by JH through the MEKRE93 pathway, in which E93 acts as the ‘adult-specifier’ that triggers metamorphosis (Ureña et al., 2014; Belles and Santos, 2014; Chafino et al., 2019). Importantly, E93 transcript levels are strongly upregulated within various metamorphic tissues (such as the wing, cuticle and prothoracic gland) at metamorphic stages during the juvenile period (Fig. 6B; Martín et al., 2021). In contrast, there was no noticeable increase in E93 expression during the transition to adulthood from the final instar nymph stage of T. domestica (Fernandez-Nicolas et al., 2023), corresponding to the fact that firebrats do not undergo typical metamorphosis. However, the phenotypic defects we observed in the E93 mutants or E93 RNAi-knockdown specimens of T. domestica (Fig. 2 and Fig. S3) indicate that E93 might regulate cuticle development and moulting through cumulative effects during post-embryonic development (here, we focus on the function of E93 in regulating reproduction in female firebrats). E93 was previously proposed to have contributed crucially to the origin of metamorphosis (Belles, 2019; Fernandez-Nicolas et al., 2023); therefore, the way in which E93 expression is innovated for regulating metamorphosis is certainly a remarkable challenge for future evolutionary development studies in insects (Fig. 6B).

The role of E93 in reproduction is not well understood, although it would be particularly advantageous to elucidate its potential in ovary development in both ametabolous and hemimetabolous insects. The regulatory mechanism by which E93 is involved in the previtellogenic development of ovarian follicles needs exploration. The role of E93 in the ovary should also be investigated in additional insect species. We speculate that its role is likely widespread. Nevertheless, the role of E93 in the fat body is weak in ametabolous insects and has become crucial during gonadotrophic cycles in the process of insect evolution, as evidenced in both hemimetabolous and holometabolous insects (Wang et al., 2021; Fig. 6B).

From an evolutionary development point of view, we hypothesize that innovations in E93 spatiotemporal expression throughout the life cycle of insects might play a major role in driving the evolution of insect metamorphosis. This hypothesis should apply to other key genes involved in controlling insect metamorphosis, such as Br-C, Kr-h1 and Chinmo.

T. domestica specimens were obtained from a colony reared in the dark at 37±0.5°C and ∼80% relative humidity (RH) in plastic containers and were fed commercial rat food. Under our rearing conditions, firebrat nymphs at nine nymphal stages (N1-2 lasts 3 days; each instar requires 5 days during N3-8, and N9 lasts 7-9 days). Subsequently, the firebrats moult into sexually mature adults (at N10, N11, etc., each requiring approximately 10 days) and pass through 45-60 instars that undergo moulting and reproductive cycles.

B. germanica cockroaches were obtained from a colony fed commercial rat food and sufficient water and reared at 27±1°C and 60-70% RH with a 12 h:12 h (light:dark) photoperiod. The fourth (N4) and fifth (N5) instar nymphs last for 7 and 9 days, respectively, and N6 (the last instar nymph) lasts for 11 days under the above mentioned rearing conditions. Female adults laid eggs enwrapped in ootheca on the ninth day (AD9).

20E (ChemFaces, China) was dissolved in DMSO (100 mg/ml) and diluted 100 times with 0.01 M PBS. 200 nl of the solution was then delivered into the abdomen of the firebrats using a Nanoject II injector (Drummond Scientific), and 1 µl was delivered into German cockroaches using a 10 µl microsyringe (BIOLAND, China). N9D4 firebrats or N6D6 cockroaches were treated with 20E, and an equal volume of 1% DMSO (diluted with 0.01 M PBS) served as a control. Ovs were dissected from the specimens at 12 h after treatment and subjected to qPCR to analyse the expression of E93 and other genes involved in 20E signalling.

Total RNA was prepared from the target tissues with AG RNAex Pro reagent (Accurate Biology, China). The spatial expression profiles of E93 in T. domestica and B. germanica were determined using tissues from the antennas, legs, cercus, paranotum, abdominal tergite, brain, Ovs and FB of firebrats at N9D6, and from cockroaches at N6D9. Previtellogenic Ovs at N9D6 and N10D8 from E93 mutant and WT control T. domestica were collected; correspondingly, previtellogenic Ovs at N6D9 from E93 RNAi and control B. germanica were collected to study the transcriptional regulatory mechanism. In addition, the vitellogenic Ovs of WT T. domestica at N10D3 (compared with N10D8-WT) and the vitellogenic Ovs of WT B. germanica at AD6 (compared with AD1-WT) were sampled for RNA-seq. Furthermore, previtellogenic FBs dissected from E93^−/− and WT female firebrats at N10D6, vitellogenic FBs from WT female firebrats at N10D3, and FBs dissected from adult cockroaches at AD6 (vitellogenesis) and AD9 (previtellogenesis) (comprising RNAi-E93 and RNAi-Lta controls, respectively) were subjected to RNA-seq analysis. Each group contained three independent biological replicates. Library construction, sequencing with an Illumina NovaSeq 6000 and mapping to the genomes with the genomic short-read nucleotide alignment program (GSNAP) were conducted by BMKGENE (Beijing, China). The B. germanica genome (Harrison et al., 2018) sequence and annotation files were downloaded from the NCBI database (https://www.ncbi.nlm.nih.gov/datasets/taxonomy/6973/). Genome sequencing and annotation of T. domestica were conducted in our laboratory (Y.-X.L., Q.-Q. Jia, Y.B., Y.-F.W., Y. Xiao, M.Z., Y.-N.L., D.-W.Y., S.-L. Liu, X. Zhou, S.L., unpublished). Differential gene expression analysis was performed as described previously (Bai et al., 2022), with an adjustment: the absolute value of log2(fold change)>1.

cDNA was synthesized with M-MLV (H-) reverse transcriptase (Vazyme, China) from an aliquot of 1000 ng of total mRNA from each sample. qPCR was performed using 2×Hieff Unicon Universal TaqMan Multiplex qPCR Master Mix (Yeasen, China) and the Applied Biosystems QuantStudio 6 Flex Real-Time PCR System (Thermo Fisher Scientific). The ribosomal protein S26 (RPS26) or Actin5C was chosen as a reference gene for analysing relative expression levels in T. domestica (Bai et al., 2020) or B. germanica, respectively. The sequences of all primers used for qPCR in this study are listed in Table S1.

For in situ hybridization studies, E93 riboprobe templates were prepared from the respective cDNA fragments cloned and inserted into the plasmid. The sequences of the primers used for amplification of the templates are shown in Table S2. The Ovs from N9D3, N9D6, N10D3, N10D6 and N10D9 in T. domestica, and from N6D3, N6D9, AD2 and AD5 in B. germanica were dissected and then fixed in 8% formaldehyde/PBT with 1% EGTA for 1 h at room temperature, washed with 100% methanol and stored at −20°C overnight before staining. Before rehydration, the Ovs were treated with 1:1 (v/v) xylene:ethanol. After rehydration, the Ovs were permeabilized with 80% acetone at −20°C for 10 min and then post-fixed in 8% formaldehyde/PBS for 20 min at room temperature. These steps are detailed in a previously published protocol (Shippy et al., 2009; Tomoyasu et al., 2009), except that a ratio of 1:50 riboprobe:hybridization buffer was used for hybridization and longer washes immediately after hybridization (two 3 h washes). Finally, the Ovs were stained using Fast Red (Sigma, F4648) for fluorescence imaging.

The 20-bp gRNA targets immediately upstream of the protospacer adjacent motif (PAM) were designed with the software sgRNAcas9 (Xie et al., 2014). The mutation strategy that we adopted involved the deletion of a large nucleic acid fragment within the conserved RHF2 of the TdE93 gene with a combination of two gRNAs (marked as gRNA1 and gRNA2). The E93-specific gRNAs were transcribed in vitro from synthesized DNA templates and column-purified using a GeneArt Precision gRNA Synthesis Kit (Invitrogen) according to the manufacturer's instructions, and the primers used for gRNA biosynthesis are listed in Table S3. The synthesized gRNA was stored at −80°C until use. Microinjections were performed using a FemtoJet 4r micromanipulator and FemtoJet 4i microinjector (Eppendorf) under an SZX16 microscope (Olympus). The fertilized eggs were collected within 4 h of oviposition and aligned on double-sided tape on glass slides. The final concentration of each gRNA was 400 ng/µl, and that of Cas9 was 300 ng/µl. The solution was injected into the eggs with glass capillary tubes (Drummond Scientific Company, USA). Injected embryos were incubated at 37°C and 70% RH for 11-12 days until hatching.

After the G0 nymphs were reared until the adult stage, genomic PCR followed by sequencing was carried out to identify the mutant alleles generated by the CRISPR/Cas9 system. Genomic DNA was separately extracted from the leg tips of adults using a TIANcombi DNA Lyse&Det PCR Kit (TIANGEN, China). The G0 individuals were genotyped by PCR product sequencing, and primers designed to detect mutagenesis around the gRNA target sites of E93 are shown in Table S4. The mosaic G0 male was then subjected to mating with WT virgin females to produce F1 offspring. Multiple mating groups were established, and the ratio of males to females was 1:3. Subsequently, we performed cloning sequencing of F1 heterozygous adults to screen the desired gene knockout line, the genotype of which contained a 103 bp deletion between two gRNA sites. Finally, F1 adults with the above genotype were selected for pair mating to generate F2 offspring, one quarter of which were homozygous mutants. Moreover, we identified homozygous mutants through sequencing and phenotypic observation. Owing to the extremely low fecundity of the homozygous mutants, the F2 heterozygous mutants were used to generate the next group of progeny.

Gene-specific primers (Table S5) were designed to amplify a 318 bp fragment for BgE93 and a 329 bp fragment for TdE93. Moreover, a 193 bp fragment of Mus musculus lymphotoxin A (Lta) was also cloned to serve as an unrelated control. The validated plasmid and primers containing the T7 RNA polymerase promoter were subjected to PCR amplification to produce DNA templates that were subsequently used for dsRNA synthesis. dsRNA synthesis and purification were performed using a T7 RiboMAX Express RNAi System (Promega) according to the manufacturer's instructions. In the study of T. domestica E93 function, E93 dsRNA was continuously injected at N7D2, N8D2 and N9D2 (1.5 μg/insect/once). In the study of B. germanica E93 function, E93 dsRNA (20 μg/insect) was applied at N6D1 (for ovary research), AD0 or AD6 (for FB research), and then, qPCR, RNA-seq or phenotypic examination was conducted at the designated time points.

The peer review history is available online at https://journals.biologists.com/dev/lookup/doi/10.1242/dev.202518.reviewer-comments.pdf