Diapause induction, termination, sequencing, and gene identification
After the completion of induction process, the parasitized host eggs were dissected to verify whether T. dendrolimi entered diapause successfully or not. Diapausing parasitoids should stay at the prepupal stage, noted as Dpre. If diapause induction failed, they would die or continue to develop into pupae or adults. Following the diapause termination process, the parasitized host eggs were transferred to normal development conditions (26°C ± 1°C, 60% ± 5% RH, 16: 8 L: D). If diapause was disrupted, the parasitoids would convert from prepupa to pupa within several days, noted as Dp. Conversely, if diapause was undisrupted, the T. dendrolimi would remain in the prepupa stage. In this study, 99% of the parasitoids entered diapause successfully, and about 95% broke out of the diapause state after the 70 days termination treatment. And the prepupae and pupae of T. dendrolimi developed under normal conditions were also obtained as references noted as NDpre and NDp, respectively.
RNA samples obtained from distinct stages of this species were prepared and sequenced using the Illumina Hiseq2000 sequencing platform. Four cDNA libraries were constructed from the samples described above (Dpre, Dp, NDpre, and NDp). After filtering raw reads (reads containing adaptor, reads containing N larger than 10%, and low-quality reads (Qphred < 20) were removed), clean reads were retained (Table 1). The clean data were assembled by Trinity and Corset with 87,022 unigenes, and an average length of 1,604 nt and an N50 of 3,148. Of the unigenes, 35,231 (40.5%) were longer than 1000 bp (Table 2).
In order to gain information of gene function, transcripts were annotated using BLASTX searches against the non-redundant (NR) sequence database; 39,969 (45.92%) of them displayed homology to known proteins (E < 1e-5; Fig. 1A). Nearly 25,000 annotated unigenes, over 65% of the annotated unigenes, were homologous to T. pretiosum, probably because the genome of T. pretiosum was the only available one in Trichogramma. Fewer unigenes were homologous to Nasonia vitripennis (1,217, 3.1%), Apis florea (26, 0.03%), A. dorsata (12, 0.01%), or A. cerana (15, 0.01%) (Fig. 1B). Less than 40 unigenes were matched to those from Microplitis demolitor. Among all annotated unigenes, 73.0% unigenes had significant homology with an E-value of <10-30 (Fig. 1C), and 52.3% had a similarity higher than 80.0% (Fig. 1D). After filtering and removing redundant sequences, we kept those with significant BLAST hits and constructed a non-redundant dataset containing 7,593 sequences with an average length of 3,351 nt. Based on the annotation, such as gene length, ID, and speculative function, the diapause-related genes and potential genes involved in diapause were sorted out for further analysis.
Identification of DEGs and functional classification
Ten genes were selected for validation with qPCR, and glyceraldehyde phosphate dehydrogenase (GAPDH) was selected as the reference gene after measuring its stable expression level for diapause and non-diapause. The tendencies of these genes’ expression profilings were similar according to RNA-Seq and qPCR (Fig. 2). Among these ten selected genes, all except trehalase (tre) were up-regulated during diapause.
To get more acquainted with diapause-specific transcriptional changes in T. dendrolimi induced by low temperature, we performed pairwise comparisons between different libraries to identify the DEGs. The results showed that there were 5,702 DEGs were identified among four groups. Among these DEGs, there were 3,182 DEGs changed in Dpre compared to NDpre. While, DESeq2 identified 3,251 and 3,442 DEGs exclusively changed in Dp vs NDp and Dpre vs Dp, respectively. In addition, DEGs changed in NDpre vs NDp were 1,511. This group of DEGs may be the genes relative to normal development, namely from prepupa to pupa, not to diapause development. According to the venn diagram, there were 463 genes changed throughout the diapause development process, while in normal development process, the expression of these genes did not change (Fig. 3).
To figure out the potential function of identified DEGs, GO enrichment was performed. From the results, in all combination except for Dpre vs Dp, more genes were up-regulated. But, when compared Dpre to Dp, there was little difference in the number of up- and down-regulated genes. Besides, more DEGs were assigned to the same category among different groups. Regulation of transcription, DNA-templated, oxidation-reduction process, and signal transduction were the top three in all these four groups. Interestingly, the number of DEGs involved in ribosome biogenesis was much higher during diapause development than normal development (Fig. 4). Moreover, the subsequent analyses are based directly upon these results.
Comparative analysis of genes involved in diapause
Diapause is a dynamic process accompanied by a series of physiological transitions. Several studies have focused on the general gene expression pattern of insect diapause without clear elucidation of the diapause mechanism. This is due to the complexity of the diapause process as well as the variations among insect species.
According to the results of GO enrichment, firstly, we focused on the genes enriched in the oxidation-reduction process, regulation of transcription, DNA-templated and signal transduction. In addition, ribosome biogenesis was something that we were looking at.
A total of 342 genes were identified in oxidation-reduction process, 16 of which belong to cytochrome P450s (CYP450s). CYP450s are hemoproteins involved in several physiological processes, such as biosynthesis of hormones and degradation of xenobiotics [36]. There are four clans in P450 supergene family, namely CYP2, CYP3, CYP4, and Mitochondrion (Mito) [37]. In T. dendrolimi transcriptome, 22 CYP450s were identified, and 16 were differentially expressed. These 16 genes belonged to 4 clans. In diapause stages (Dpre), there were ten genes were up-regulated. While, in the pupae after diapause (Dp), five genes highly expressed. And only one gene (CYP9E2) was extremely up-regulated in normal pupae (NDp) (Fig. 5). From these results, we can see that the number of up-regulated genes during diapause in this species was significantly higher than that in other stages. Previous studies showed that CYP450 of Schistosoma mansoni was essential for worm survival and egg development [38]. Besides, CYP4G1 was proved to be related to cuticular hydrocarbon biosynthesis in Drosophila [39]. These genes were up-regulated in diapause individuals, suggesting that when T. dendrolimi entered diapause, the environmental conditions were not suitable for survival. The condition was worse than that under normal condition. In this process, many harmful substances may be produced. Thus, a possible function for these genes is to balance out harmful substances, maintaining cellular homeostasis. Besides, this phenomenon still occurred in pupae after diapause (Dp) compared to that under normal condition (NDp). Given to this result, we could infer that although diapause helped insect survive adverse environments, the process itself also had an effect on insect growth. Diapause development is a tradeoff that insects make in order to survive.
According to the results, there were 36 transcription factors differentially expressed during diapause development. Based on these results, it was speculated that three kinds of transcription factors might be associated with diapause in T. dendrolimi.
Firstly, zinc finger protein. There were three genes encoded zinc finger protein, zinc finger protein 271, zinc finger 184, and zinc finger 544, were identified in the transcriptome. They were all up-regulated that the expressions of these three genes of Dpre were higher than that of Dp. Zinc finger protein gene 271 had an SFP domain. Genes containing this domain were considered to be a putative transcriptional repressor during G2/M (second gap period to mitotic period) transition. One of the most noticeable characteristics of diapause is the blockage of ontogeny, and this blockage always occurs with the cell cycle cessation [40, 41]. In N. vitripennis, the S phase of the cell cycle disappeared in the beginning stage of diapause due to the cells being arrested in the G0/G1 (stop cell division to first gap) and G2 phases[42]. Similarly, in drosophilid fly, Chymomyza costata, the cell cycle of CNS cells was arrested in the G0/G1 (86.6%) and G2 (12.8%) division phases during diapause [43]. The wee1 gene, encoding for a kind of inhibitory kinase, was up-regulated during diapause in this species[43]. We obtained similar results in T. dendrolimi (Fig. 6). Therefore, the wee1 gene might be a molecular marker of diapause. Besides, zinc finger protein 184 contained a GDT1 domain, which is a putative Ca2+/H+ antiporter. Ca2+/H+ antiporter, which maintain the homeostasis, has been studied a lot in plants, but there are few studies in insects.
Second, homeobox domain proteins. In T. dendrolimi transcriptome, 11 homeobox-containing genes were differentially expressed during diapause stage except for pit1, which was significantly up-regulated in the individuals which terminated diapase. Among these genes, homeobox protein homothorax (hth) showed the greatest change in expression, followed by homeotic protein distal-less (dll) and homeobox protein six1 (six1). In Tribolium castaneum, during embryogenesis, hth plays a broad role in the segmentation process and is required for specification of body wall identities in the thorax [44]. In other species, like D. melanogaster, hth has different functions in different tissues. Hth located in head leads to opposite effects on eye and antennal development as a negative regulator of eye development, and it acts with extradenticle (exd) to delimit the eye field and prevent inappropriate eye development [45]. Transcriptome factor dll plays a role in larval and adult appendage development [46]. The gene expression of another homeobox protein six1 was similar to hth, expression of which increased after entering diapause stage. And this gene is postulated to be involved in the regulation of cell proliferation, apoptosis, and embryonic development. In mammals, six1 is essential for early neurogenesis in the development of olfactory epithelium [47]. They were crucial genes regulating the myogenesis and extremely up-regulated during diapause in T. dendrolimi. But the aberrant expression of these genes may cause the cessation of growth. So, we inferred that these genes up-regulated at diapause stage blocked the normal cell cycle in diapause T. dendrolimi.
Third, forkhead box proteins. In previous studies, FOXOs have been identified as promising candidates for the molecular control of embryonic diapause in some species, like Culex pipiens [22, 48]. In T. dendrolimi, three forkhead box proteins, foxo, forkead box protein E3 (foxe3), and forkhead box protein D3 (foxd3) were identified. These genes were up-regulated both in the prepupae and pupae individual in diapause conditions. So, we inferred that they may play a role in diapause development.
In signal transduction term, a variety of enzymes were speculated to play an important role in diapause. Protein phosphatase 2A (PP2A) was one that caught our attention. It is a key serine-threonine protein phosphatase, which regulates several cellular processes, including metabolism, transcription, cell cycle, autophagy, and signal transduction [49, 50]. Cell cycle withdrawal, from G1 to S stage, is negative related to the activity of PP2A [51, 52]. In T. dendrolimi transcriptome, the expression of PP2A was down-regulated during diapause stage. This result was consistent with that obtained in the cotton bollworm, Helicoverpa armigera. Low PP2A expression in diapause-destined individuals contributed the accumulation of p-Akt, and p-Akt leads to H. armigera into diapause [24, 53]. And based on the previous studies, PP2A plays a multi-faceted role in the regulation of a couple of pathways, such as mTOR and wnt signaling pathway, which were related to cell cycle [49, 54]. So, we supposed that PP2A may have a certain function in T. dendrolimi diapause.
Except for these three biological processes, ribosome biogenesis also has an vital role in the control of cell growth and division in eukaryotes and is worth paying attention to [55]. In the present study, ribosome biogenesis involves 31 DEGs, of which 29 genes were up-regulated during diapause in the prepupae stage. Only two genes, 40S ribosomal protein S11 (rpS11) and 28S ribosomal protein S5 (rpS5). And all the 60S ribosomal proteins were up-regulated (Fig. 7). It has been reported that the rate of ribosome synthesis during diapause was lower than that of non-diapause eggs in Bombyx mori, so the up-regulation of ribosomal proteins played an important role in blocking diapause [56, 57]. In mosquitos, C. pipiens, the expression of ribosomal protein S3a (rpS3a) was dramatically reduced for a short time in diapause stage. After the injection of rpS3a dsRNA into non-diapaused females, the development of follicle was arrested, similarly to in the diapause state [58]. Conversely, in T. dendrolimi, diapause prepupae had a higher expression of ribosomal proteins compared to non-diapause wasps. It is noticeable that B. mori and C. pipiens enter diapause as adults, whereas T. dendrolimi enter diapause as prepupae. Although diapause decreased metabolism level, they may still need more energy to maintain the diapause condition. It was the first time to report that the ribosomal protein was related to diapause in T. dendrolimi. But the specific function of these genes involving in diapause still need to be verified.
In addition, some genes, even if they were not involved in the biological process with a considerable number of genes enriched, were still important in diapause development of T. dendrolimi. Taking the following genes for examples.p53 induces growth arrest or apoptosis and has a negative regulation on cell division by controlling an array of genes in this process [59]. In addition, p53 and DNA damage-regulated gene 1 (pdrg1) are involved in multiple cellular processes, such as apoptosis, DNA damage repair, cell cycle. In Artemia sinica, PDRG1 takes an essential role in diapause termination and regulation of cell cycle during early embryonic development [60]. By blocking the expression of pdrg1 in human colon cancer cells, cell growth reduced significantly. Transcriptome based analysis indicated that the expression of p53 significantly increased during diapause, and the contrary result was observed for pdrg1 in T. dendrolimi. Moreover, this result also revealed that apoptosis activity was enhanced although the diapause individual remained in a dormancy state. It is speculated that, the wasps had to face the adverse condition, when more harmful substances could be accumulated during diapause. So, it was necessary to boost apoptosis activity for wasp to survive.
UDP-glucuronosyltransferase (UDPGT) is of central importance in the union and next elimination of toxic xenobiotics and endogenous compounds. Glutathione-S-transferase (GST) belongs to a multifunctional protein family mainly located in the cytoplasm. Both are related to cellular detoxification. The transcriptional expression of GST and UDPGT were up-regulated during T. dendrolimi prepupal diapause. This was different from other species, such as solitary bee Tetrapedia diversipes and Tetranychus urticae in which transcripts are down-regulated during diapause [62]. In other species, the downregulation of GST and UDPGT might be correlated with non-feeding conditions, so the amount of exogenous substances decreases accordingly. However, T. dendrolimi is a kind of endoparasitoid wasp that spends its whole life cycle within the host egg, it means that starvation is a rare situation until the prepupae stage. Therefore, we speculated that the reason for the observed up-regulation of these two genes in diapause T. dendrolimi might be due to the increased resistance under unfavorable environmental conditions.
Lipid metabolism is essential for energy homeostasis. It has been shown that some diapausing insects use lipids as predominant energy stores[63, 64]. During diapause, almost all selected lipid metabolism related genes were up-regulated, coinciding with the mobilization of TAG reserves (Fig. 8). These genes have important influence in the formation of triacylglycerol, which is the main caloric reserve during diapause. In addition, several genes’ expressions significantly increased after diapause termination. This might be due to post-diapause development. In our previous study, diapause T. dendrolimi had higher numbers of parasitized hosts than non-diapause T. dendrolimi [27]. The increasing lipid storage could provide more energy for maintaining activities.