Next Article in Journal
Evaluation of Natural and Factitious Food Sources for Pronematus ubiquitus on Tomato Plants
Next Article in Special Issue
Supercritical Carbon Dioxide Impregnation of Gold Nanoparticles Demonstrates a New Route for the Fabrication of Hybrid Silk Materials
Previous Article in Journal
Effects of Menadione on Survival, Feeding, and Tunneling Activity of the Formosan Subterranean Termite
Previous Article in Special Issue
Products of Sericulture and Their Hypoglycemic Action Evaluated by Using the Silkworm, Bombyx mori (Lepidoptera: Bombycidae), as a Model
 
 
Font Type:
Arial Georgia Verdana
Font Size:
Aa Aa Aa
Line Spacing:
Column Width:
Background:
Article

Transcriptome Analysis Reveals the Gene Expression Changes in the Silkworm (Bombyx mori) in Response to Hydrogen Sulfide Exposure

1
School of Food and Biological Engineering, Hefei University of Technology, Hefei 230009, China
2
School of Biological Science and Engineering, North Minzu University, Yinchuan 750021, China
3
Jiangsu Key Laboratory of Sericultural Biology and Biotechnology, School of Biotechnology, Jiangsu University of Science and Technology, Zhenjiang 212003, China
4
Key Laboratory of Silkworm and Mulberry Genetic Improvement, Ministry of Agriculture, Sericultural Research Institute, Chinese Academy of Agricultural Sciences, Zhenjiang 212018, China
*
Authors to whom correspondence should be addressed.
Co-first author.
Submission received: 3 November 2021 / Revised: 27 November 2021 / Accepted: 9 December 2021 / Published: 13 December 2021
(This article belongs to the Special Issue Silkworm and Silk: Traditional and Innovative Applications)

Abstract

:

Simple Summary

The fat body is one of the most important tissues in the body of insects due to its number of functions. Nowadays the new physiological function of H2S has gained attention as a novel signaling molecule. H2S performs crucial regulatory functions involving growth, the cardiovascular system, oxidative stress, and inflammation in many organisms. In this study, RNA-seq technology was used to investigate the fat body of the silkworm at the transcriptional level after H2S exposure during the 5th larvae stage. A total of 1200 (DEGs) was identified after 7.5 µM H2S treatment, of which 977 DEGs were up-regulated and 223 DEGs were down-regulated. DEGs were mainly involved in the transport pathway, cellular community, carbohydrate metabolism, and immune-associated signal transduction. Present research provides new insights on the gene expression changes in the fat body of silkworms after H2S exposure.

Abstract

Hydrogen sulfide (H2S) has been recognized for its beneficial influence on physiological alterations. The development (body weight) and economic characteristics (cocoon weight, cocoon shell ratio, and cocoon shell weight) of silkworms were increased after continuous 7.5 µM H2S treatment. In the present study, gene expression changes in the fat body of silkworms at the 5th instar larvae in response to the H2S were investigated through comparative transcriptome analysis. Moreover, the expression pattern of significant differentially expressed genes (DEGs) at the 5th instar larvae was confirmed by quantitative real-time PCR (qRT-PCR) after H2S exposure. A total of 1200 (DEGs) was identified, of which 977 DEGs were up-regulated and 223 DEGs were down-regulated. Most of the DEGs were involved in the transport pathway, cellular community, carbohydrate metabolism, and immune-associated signal transduction. The up regulated genes under H2S exposure were involved in endocytosis, glycolysis/gluconeogenesis, the citrate cycle (TCA cycle), and the synthesis of fibroin, while genes related to inflammation were down-regulated, indicating that H2S could promote energy metabolism, the transport pathway, silk synthesis, and inhibit inflammation in the silkworm. In addition, the expression levels of these genes were increased or decreased in a time-dependent manner during the 5th instar larvae. These results provided insight into the effects of H2S on silkworms at the transcriptional level and a substantial foundation for understanding H2S function.

1. Introduction

The silkworm, Bombyx mori, is a pivotal model of Lepidoptera, a holometabolous insect, with a considerable economic importance in the world. The silkworm has been considered a model animal because of its appropriate life cycle, cultivation, and fast reproduction [1,2]. The fat body for insects acts as the dynamic tissue and the metabolic organ, mainly involved in synthesis, nutrient storage, and energy metabolism [3], and it has been widely studied in various studies of silkworm to reflect various biological, physiological, and biochemical processes. The over expression of BmFoxO in the fat body could suppress protein translation in transgenic silkworms [4].
The toxic effects of hydrogen sulfide (H2S) as an environmental pollutant have been studied in the past years. However, nowadays, the new physiological function of H2S has gained attention as a novel signaling molecule. H2S performs crucial regulatory functions involving growth, the cardiovascular system, oxidative stress, and inflammation in many organisms [5,6]. H2S has been accepted as the third transmission gas followed by nitric oxide (NO) and carbon monoxide (CO) in past years with applications to various insects. Miller and Roth (2007) found that treatment with high concentrations of H2S are toxic to Caenorhabditis elegans; However, exposure to low concentrations of H2S increased the thermotolerance and lifespan via increasing the activity of SIR-2.1 [7]. Budde and Roth (2010) found that hif-1 is required when C. elegans were exposed to H2S; in turn, H2S could increase both HIF-1 protein concentration and nuclear localization [8]. H2S is an endogenous regulator of oxidative damage, metabolism, and aging in C. elegans [9]. Wei and Kenyon (2016) found that removing germ cells in C. elegans triggers a level of hydrogen sulfide in nonreproductive somatic tissues. Hydrogen sulfide displayed the protective responses that slow aging [10]. H2S treatment significantly increased the survival and lifespan of Drosophila melanogaster under arid and food-free conditions [11,12].
In our previous research, it was found that the development (body weight) and economic characteristics (cocoon weight, cocoon shell ratio, and cocoon shell weight) of silkworm (Bombyx mori) were increased after continuous treatment with 7.5 μM of H2S. After exposure to H2S, the level of metabolites related to inflammation ((6Z,9Z,12Z)-octadecatrienoic acid, hexadecanoic acid, choline phosphate, and malic acid, etc.) in the hemolymph of silkworms were significantly increased compared to control group [13]. However, when the concentrations of H2S were higher than 7.5 μM, H2S displayed toxic effects to the silkworm, including a decrease in body weight, cocoon weight, cocoon shell weight, etc. RNA-seq technology is widely used in the field of silkworms to investigate physiological and biochemical changes at the transcriptional level under diverse conditions [14]. However, there is no effective report for the effects of H2S on silkworms at the transcriptional level.
The present study aimed to investigate the expression level of genes in the fat body of silkworms after exposure to H2S. The results based on the transcriptomic and bioinformatic analysis showed significantly different responses to H2S and provided useful insights into the role of H2S on silkworms.

2. Materials and Methods

2.1. Silkworm Rearing and H2S Exposure

P50 strains of silkworm were obtained from the Sericulture Research Institute of the Chinese Academy of Agricultural Sciences, Zhenjiang, China. According to the previous study [13], from the 4th instar larvae to the 5th instar larvae, silkworms were divided into two groups (H2S-treated and control) and reared in desiccators using fumigating treatment at 25 ± 1 °C and 12L:12D photoperiod conditions. Silkworms in the H2S-treated group were treated by NaHS (Shandong West Asia Chemical Company, Jinan, China), a donor of H2S. The experimental concentration of H2S was set to 7.5 μM according to the report of Cao et al. (2020) [13]. There were triplicates in each group with 30 silkworms per replicate and the silkworms were fed abundant fresh mulberry leaves twice a day.

2.2. RNA Extraction

The total RNA of the silkworm fat body was extracted on ice from the 1st day (5L1D) to the 5th day (5L5D) at the 5th instar larvae using Trizol reagent (Invitrogen, CA, USA). The concentration and purity of total RNA were analyzed by NanoDrop ND-1000 (NanoDrop, Wilmington, DE, USA) and the ratios of A260/A280 and 260/230 were used to measure the quality of total RNA.

2.3. Sequencing, Data Processing, and Assembly

The total RNA of 5L3D was analyzed in the H2S-treated and control groups by RNA-Sequencing. Total RNA (10 µg) of each sample was subjected to Poly(A) mRNA isolation with poly-T oligo-attached magnetic beads (Invitrogen, CA, USA). The poly(A) RNA was fragmented into small pieces using a Magnesium RNA Fragmentation Module (New England Biolabs (Beijing) LTD., cat. e6150, Beijing, China) under 94 °C for 5–7 min. The cleaved RNA fragments were then reverse-transcribed to create the cDNA by SuperScriptTM II Reverse Transcriptase (Invitrogen, CA, USA), which were next used to synthesize U-labeled second-stranded DNAs. An A-base was then added to the blunt ends of each strand, preparing them for ligation to the indexed adapters. Each adapter contained a T-base overhang for ligating the adapter to the A-tailed fragmented DNA. Single- or dual-index adapters were ligated to the fragments, and size selection was performed with AMPureXP beads. After the heat-labile UDG enzyme (New England Biolabs (Beijing) LTD., cat.m0280, Beijing, China) treatment of the U-labeled second-stranded DNAs, the ligated products were amplified with PCR by the following conditions: initial denaturation at 95 °C for 3 min; 8 cycles of denaturation at 98 °C for 15 s, annealing at 60 °C for 15 s, and extension at 72 °C for 30 s; and final extension at 72 °C for 5 min. The average insert size for the final cDNA library was 300 ± 50 bp. Lastly, we performed the 2 × 150 bp paired-end sequencing (PE150) on an Illumina NovaseqTM 6000 (LC-Bio Technology CO., Ltd., Hangzhou, China) following the vendor’s recommended protocol.
The transcriptome was sequenced on basis of the Illumina paired-end RNA-seq approach. The average insert size for the paired-end libraries was 300 bp (±50 bp). The low-quality reads with sequencing adaptors, sequencing primer, and nucleotide (q quality score (p value corrected by BH algorithm, also called FDR value, q value, p adj value) <20) were removed using Cutadapt software (https://cutadapt.readthedocs.io/en/stable/, accessed on 1 August 2021) [15] before assembly. We used HISAT2 software (https://daehwankimlab.github.io/hisat2/,version:hisat2-2.0.4, accessed on 1 August 2021) [16] to map reads to the genome. Three important parameters, such as Q20 (the percentage of bases with mass value ≥20), Q30 (mass value ≥30), and GC content were verified to evaluate all the reads using the FastQC online tool [17].
The Bombyx mori reference genome [18] was used for the alignment of the sample reads using the HISAT package [19]. StringTie (http://0-ccb-jhu-edu.brum.beds.ac.uk/software/stringtie/, accessed on 1 August 2021) was used for the assembly of the mapped reads of each sample and a comprehensive transcriptome was constructed combining all transcriptomes using Perl scripts [20]. The expression level of mRNAs was used to calculate the fragments per kilobase of exon model per million mapped reads (FPKM) by StringTie (http://www.bioconductor.org/packages/release/bioc/html/ballgown.html, accessed on 1 August 2021) [20]. The FPKM could evaluate the abundance of genes. The DEGs were selected with log2 (fold change) >1 or log2 (fold change) <−1 at p-value < 0.05, and q-value > p-value by R package-edgeR (https://bioconductor.org/packages/release/bioc/html/edgeR.html, accessed on 1 August 2021) [21].

2.4. Functional Annotation and Enrichment Analysis of DEGs

The reference sequences for all the silkworm genes were obtained from the NCBI database and BLASTX was used for annotation. All the sequences were aligned against Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) databases using DIAMOND (0.7.12) [22]. GO function or KEGG pathway significant enrichment analysis first mapped all the DEGs to each GO term in the GO database or each KEGG pathway in the KEGG database, and second calculated the number of genes in each GO term or KEGG pathway and finally applied a hypergeometric test to find a GO term or KEGG pathway that were significantly enriched in DEGs compared to the entire genome background.

2.5. Quantitative Real-Time PCR (qRT-PCR)

The qRT-PCR was used to validate the accuracy of DEGs in L5D3 at the transcriptome level and the expression pattern of DEGs in the 5th instar larvae was investigated. Reverse transcription was performed by HiScript II QRT SuperMix (Vazyme, Nanjing, China) as follows: total RNA (1 μg) (without genome contamination) and 4 × gDNA Wiper were mixed with 5 × HiScript II QRT SuperMix and then incubated at 25 °C for 10 min, 50 °C for 15 min, and 85 °C for 2 min [14]. The cDNA was used for real-time PCR detection using the UltraSYBR Mixture (Cwbio, Beijing, China) and qRT-PCR was executed using the Light Cycler® 480 System (Roche, Basel, Switzerland). The primer sequences for target genes were listed in Table S1 using primer 3.0 software (Premier Biosoft International, Palo Alto, CA, USA). The expression level of action A3, a silkworm housekeeping gene, was regarded as an internal reference in standardization [14]. Before Quantitative real-time PCR, the efficiency and specificity of primers were confirmed. The volume of the reaction system was 20 μL, containing 10 μL UltraSYBR Mixture, 1 μL cDNA, 0.5 μL forward and reverse primers and 8 μL ddH2O. The qRT-PCR program was as follows: initial denaturation at 95 °C for 10 min; followed by 40 times of denaturation at 95 °C for 10 s, annealing at 60 °C for 30 s, and extension at 72 °C for 32 s; a melting curve program was set at 95 °C for 15 s, 60 °C for 1 min, 95 °C for 15 s, and a60 °C for 15 s with a final cooling step at 4 °C for 30 s [23]. The data were analyzed with the LightCycler®96 software (Roche Diagnostics, Indianapolis, IN, USA) using the 2−ΔΔCt method.

2.6. Statistical Analysis

The data were analyzed using one-way ANOVA data analysis. All the data were obtained in triplicate and presented as mean ± standard error (SE) (n = 3). Statistical significance at p < 0.05 was measured using SPSS 20.0 software (IBM, Endicott, NY, USA).

3. Results

3.1. Transcriptome Sequencing and Assembly

The deep sequencing of RNA from the fat body (5L3D) treated with H2S using the illumine sequencing and the data processing results were shown in Table 1. The raw data of H2S-treated and control groups were approximately 7.75 GB and 7.66 GB, respectively. After filtering and removing the low-quality reads, adaptors and poly-N stretches, the valid data and valid reads ratio of each replicate in H2S-treated group were obtained as follows: 7.55 Gb (99.12%), 7.87 Gb (98.91%), and 7.54 Gb (98.39%) and that of the control group were 7.42 Gb (94.00%), 7.37 Gb (98.99%), and 7.50 Gb (98.12%). The quality scores of Q30 exceeded 98% in each sample. Moreover, the mapped ratio in H2S-treated and control groups were 97.35%, 97.43%, and 97.72% and 93.19%, 97.35%, and 97.52%, respectively, which highly matched with the silkworm genome (Table S2). The expression interval of genes was mainly focused on the 0.3~3.57 FPKM interval (FI) and the ratio was about 35% gene of the control group and 39% of H2S-treated group (Table S3). These results revealed a good quality sequencing analysis and the genes were analyzed in the following section.

3.2. DEG Identification and Analysis

DEGs in H2S and control groups were calculated according to the FPKM value of the obtained genes and then a large number of genes were differentially expressed after H2S exposure. The volcano plot of DEGs was shown in Figure 1 and the detailed information of DEGs was listed in Table S4. A total of 1200 DEGs were obtained between the H2S-treated and control groups. The number of up-regulated and down-regulated DEGs was 977 and 223, and the up-regulated genes were higher than down-regulated genes.

3.3. GO Analysis of DEGs

The GO database can provide standard terms for describing the comprehensive properties of genes in organisms as an international classification system. In this study, DEGs were assigned to different GO categories such as biological processes, cellular components, and molecular functions. The major subcategories of DEGs for biological processes were border follicle cell migration, mitotic spindle elongation, myoblast fusion and microtubule-based processes, as well as others. (Figure 2 and Table S5). Through the results of the GO analysis, it was revealed that H2S could influence these multiple biological processes by regulating the expression of genes in the fat body of silkworms at the molecular level.

3.4. KEGG Pathway Analysis

To systematically analyze the biological function and metabolic pathways of DEGs, KEGG pathway annotations were obtained. The 561 genes of 1200 DEGs predicted to encode enzymes were mapped into 225 KEGG pathways. These annotated genes were assigned to six categories including organismal systems (49 pathways), metabolism (74 pathways), human diseases (39 pathways), genetic information processing (19 pathways), environmental information processing (28 pathways), and cellular processes (16 pathways) (Figure 3). In the organismal systems, vascular smooth muscle contraction (7 DEGs) is the most abundant metabolic pathway followed by the insulin signaling pathway (6 DEGs) and toll-like receptor signaling pathway (5 DEGs). The most abundant metabolic category in metabolism is inositol phosphate metabolism (12 DEGs) followed by glycolysis/gluconeogenesis (11 DEGs), glutathione metabolism (10 DEGs), and the citrate cycle (9 DEGs). In human disease, it is mainly related to insulin resistance (17 DEGs). In the category of genetic information processing, it is mainly enriched in ribosomes (33 DEGs) closely related to synthetic genetic material. In environmental information processing, the most abundant pathways are the Hippo signaling pathway (19 DEGs), MAPK signaling pathway (16 DEGs), and phosphatidylinositol signaling system (12 DEGs). In cellular processes, 36 DEGs were mainly enriched in endocytosis followed by focal adhesion (20 DEGs) and regulation of the actin cytoskeleton (15 DEGs).
A significant enrichment analysis was further performed on the metabolic pathways of DEGs, and a total of 26 significantly enriched metabolic pathways were obtained (Table 2). The terms were enriched to multiple molecular pathways such as the transport pathway and cellular community, including endocytosis, focal adhesion, regulation of the actin cytoskeleton, tight junction, and adherens junction; amino acid and carbohydrate metabolisms, such as glycolysis/gluconeogenesis, the citrate cycle (TCA cycle), phenylalanine metabolism, and histidine metabolism; immune-associated signal transduction such VEGF signaling pathway, the TNF signaling pathway, NOD-like receptor signaling, the MAPK signaling pathway, the Hippo signaling pathway-fly, the Toll-like receptor signaling pathway, and the NF-κB signaling pathway (Table S6). The top 10 metabolic pathways affected by H2S were endocytosis, ribosome, focal adhesion, the Hippo signaling pathway, the MAPK signaling pathway, regulation of the actin cytoskeleton, tight junction, inositol phosphate metabolism, glycolysis/gluconeogenesis, and the TCA cycle (Table 2). These annotations indicated that H2S could influence the multiple biological pathways in the fat body of silkworms and provide a new perspective to the effects of H2S on silkworms.

3.5. The Validation of DEGs by qRT-PCR

To validate the accuracy of DEGs obtained from the transcriptome between the H2S-treated and control groups, 28 candidate genes were selected for qRT-PCR analysis. The selected DEGs were mainly annotated to glycolysis/gluconeogenesis, the TCA cycle, focal adhesion, endocytosis, adherens junction, the NF-κB signaling pathway, the MAPK signaling pathway, and the synthesis of fibroin. Good consistency between the results of qRT-PCR and the transcriptome confirmed the accuracy and reliability of the sequencing data and revealed the significant difference of these genes during H2S exposure (Figure 4 and Table S3).

3.6. The Expression Pattern Analysis of DEGs

A total of 11 DEGs was representatively selected to analyze the expression level pattern in the 5th instar larvae using qRT-PCR. Genes (relating to the glycolytic pathway and TCA cycle) such as lactate phosphoglycerate kinase (PGK), phosphoglucomutase (PGM), NADP-dependent isocitrate dehydrogenase (IDH), and pyruvate carboxylase (PC) in the H2S-treated group were up-regulated compared to control group (Figure 4A,B). Meanwhile, the levels of these genes were significantly increased at the 5L3D, 5L4D, and 5L5D comparing with 5L1D in the H2S group (Figure 5A–D). In addition, the expression levels of the members of the Rab family of small GTPases (Rab10) and trehalose transporter 1 (Tret1) in endocytosis were up-regulated after H2S treatment (Figure 4C). The expression levels of Rab10 and Tret1 were significantly increased in a time-dependent manner (Figure 5E,F). On the contrary, transforming growth factor-β activated kinase1 (Tak1) and matrix metalloproteinase-3 (MMP-3) were significantly decreased in the 5th instar larvae after H2S exposure (Figure 5G,H). Moreover, the increased expression of fibroin synthesis genes, such as the fibroin heavy chain (Fib-H), fibroin light chain (Fib-L), and glycoprotein P25 (P25) was observed during the 5th larvae stage after H2S exposure (Figure 4D). At the end of the 5th instar larvae, the expression levels were higher than the early 5th instar larvae (Figure 5I–K).

4. Discussion

H2S has been accepted as the third transmission gas followed by nitric oxide (NO) and carbon monoxide (CO) in recent years with applications to various insects, such as C. elegans [7,8,9,10], D. melanogaster [11,12], and B. mori [13]. RNA-seq technology is widely used in the field of silkworms to investigate the physiological and biochemical changes at the transcriptional level under diverse conditions [4,14,24]. In our previous research, it was found that the development (body weight) and economic characteristics (cocoon weight, cocoon shell ratio, and cocoon shell weight) were increased after continuous treatment with 7.5 μM of H2S. However, there is no effective report for the effects of H2S on the silkworm at the transcriptional level [13]. In the present study, the expression levels of genes in the fat body of silkworms after exposure to H2S was investigated using transcriptomic and bioinformatic analysis. The results showed significantly different responses to H2S and provided useful insights into the role of H2S on silkworms.
Using transcriptome, a total of 1200 DEGs were identified in the fat body of silkworms at the 5th instar larvae. After treatment with H2S, up-regulated DEGs (977) were higher than the down-regulated DEGs (223). These DEGs were mainly involved in glycolysis/gluconeogenesis, endocytosis, and the TCA cycle. The carbohydrate metabolism is one of the important energy metabolism pathways in organisms, incorporating glycolysis/gluconeogenesis, the TCA cycle, oxidative phosphorylation, and the pentose phosphate pathway [25]. Energy metabolism, especially the metabolic control of carbohydrates, is essential for insect growth and development [26]. Glucose is oxidized to pyruvate in glycolysis/gluconeogenesis [25]. Glycolytic activity is regulated by many enzymes in the glycolytic pathway such as PGK and PGM. PGM can convert glucose 1 phosphate into glucose 6 phosphate, providing the raw materials for the glycolytic pathway [27]. PGK is commonly known as a soluble and membrane-bound protein in different organisms [27]. In one of the glycolytic processes, phosphorylation of 1,3-diphosphoglycerate to 3-phosphoglycerate is catalyzed by PGK [28]. Like glycolysis/gluconeogenesis, the TCA cycle also provides the energy and intermediates for diverse biosynthetic pathways [25]. IDH, one of three rate-limiting enzymes in the TCA cycle, catalyzes the TCA cycle reaction wherein isocitrate is oxidized and decarboxylated to produce α-ketoglutarate, NADH, and CO2 [29]. PC catalyzes the carboxylation of pyruvate into oxaloacetate for gluconeogenesis, the urea cycle, lipid synthesis, and other pathways which are the primary anaplerotic pathways for the TCA cycle where the intermediates are replenished by alanine and lactate [30]. In this study, the results showed that the levels of PGM, PGK, IDH, and PC were up-regulated after H2S exposure. According to the expression pattern analysis of these genes, gene levels had gradually increased in a time-dependent manner. Especially for 5L4D and 5L5D, a significant difference was observed compared to 5L1D. These results indicated that the expression levels of the genes involved in carbohydrate metabolism were increased with H2S exposure, suggesting H2S can accelerate the carbohydrate metabolism of silkworms and provide much energy for silkworm metamorphosis.
The uptake and sorting of eukaryotes and the following recycling or degradation of fluid, membranes, membrane proteins, and macromolecules are regulated by endocytosis [31]. Endocytosis often depends on clathrin and occurs at clathrin-coated pits [32]. Utilizing endocytosis with clathrin-dependent or clathrin-independent uptake mechanisms, cargo proteins and adaptor molecules are delivered to early endosomes and sorted [33]. The Rab family of small GTPases can mediate membrane the trafficking processes of intracellular compartmentation, generation, and maintenance [34]. Rab10, which is closely related to clathrin-independent endocytosis is largely expressed in many model animals directing the trafficking of proteins. In Caenorhabditis elegans, Rab10 was required in interneuron postsynaptic glutamate receptor recycling and the transport between early endosomes and recycling endosomes [35]. Rab10 in the H2S-treated group was up-regulated compared to the control group, which indicated that the H2S-treated silkworms were more active in supporting endocytosis for the body. Tret1 participated in the transfer of trehalose, the major hemolymph sugar in most of the insect, which is synthesized in the fat body and released into the hemolymph. The transportation of trehalose with the up-regulated expression level of Tret1 after exposure to TiO2 nanoparticles was accelerated in the silkworm, showing that TiO2 nanoparticles could facilitate carbohydrate and the nutrition metabolism [36]. In this study, the level of Tret1 was up-regulated with the H2S supplement to accelerate the transport rate of trehalose. Moreover, the expression levels of Rab10 and Tret1 were increased in a time-dependent manner during the 5th larvae. These results showed that H2S could accelerate the transmission of mass in the silkworm fat body via up-regulating the expressions of Rab10 and Tret1 to provide material for larva–pupae metamorphosis.
Previous studies reported that H2S exerted anti-inflammatory effects on T lymphocytes. Rats injected with NaHS inhibited inflammatory processes such as leukocyte infiltration and adherence to the vascular endothelium [37]. In this study, DEGs were annotated to various pathways related to inflammation by KEGG analysis, including the TNF signaling pathway, the MAPK signaling pathway, and the NF-κB signaling pathway.
Tak1 is a serine/threonine kinase and a member of the mitogen-activated protein kinase kinase (MAP3K) family which is activated by receptors such as TGF-β, TNF-α, and IL-1. MMP-3 is a matrix-soluble protein containing the binding site of NF-κB and can participate in the degradation of ECM and activate the precursors of MMPs. The content of MMP-3 was increased in the mouse arthritis model and the administration of baicalin can reduce the gene expression level of MMPs and alleviate inflammatory processes [38]. After H2S exposure, levels of Tak1 and MMP-3 were down-regulated compared to the control group. Moreover, the expression levels of Tak1 and MMP-3 at the end of 5th instar larvae were significantly decreased compared to the early 5th instar larvae. These results indicated that H2S might have an anti-inflammatory role in silkworm growth.
In addition, significant changes in the genes related to fibroin synthesis were observed in this study. The domestic silkworm has tremendous economic value and the most important manifestation of its value is silk production [39]. The function of the silk gland is to synthesize silk protein and during the 5th instar larvae it is the biggest organ of B. mori [39]. It was reported that Fib-H, Fib-L, and P25 encoding fibroin were expressed in all the larval instars. However, their highest expression was noticed before the last molting period [40]. The expression levels of Fib-L, P25, and Fib-H were up-regulated in response to H2S compared to the control group. In addition, the expression of Fib-H at the end of the 5th instar larvae (5L4D, 5L5D) was rapidly increased and higher than on the other days. Similarly, the expression levels of P25 and Fib-L were also significantly increased. These results indicated that H2S could promote the fibroin synthesis by up-regulating the expression levels of Fib-L, P25, and Fib-H, which was beneficial to silk production. These results were similar to the previous study which reported the up-regulated expression of Fib-L, P25, and Fib-H after treatment with titanium dioxide nanoparticles (TiO2NPs) [41].

5. Conclusions

To summarize, this is the first-ever report to investigate the fat body of silkworms at the transcriptional level after H2S exposure using the transcriptome analysis during the 5th larvae stage. Transcriptome analysis could detect numerous genes and signaling pathways and most of the DEGs were up-regulated compared to the control group. The expression of key genes related to energy metabolism, transport pathways, and fibroin synthesis in silkworms after H2S exposure was up-regulated; conversely, the genes involved in inflammatory processes were down-regulated, which indicated that H2S may lead to the promotion of biological processes such as energy and nutrition metabolism and inhibition of inflammatory processes. These results provide a novel approach to comprehending the molecular mechanisms underlying the H2S effects in the silkworm fat body.

Supplementary Materials

The following are available online at https://www.mdpi.com/article/10.3390/insects12121110/s1, Tables S1–S6. Table S1. Primer’s sequence of genes used in qRT-PCR analysis. Table S2. Statistical analysis of the RNAs Libraries. Table S3. Distribution of gene expression levels between H2S-treated and the control group and the partial DEGs. Table S4. The list of DEGs comparison between H2S treated group and control group. Table S5. Analysis of the enriched GO terms. Table S6. Analysis of the enriched KEGG pathway.

Author Contributions

Conceptualization, Z.-J.W. and F.H.; methodology, R.Z. and Y.-Y.C.; software, J.D.; validation, R.Z., J.D. and Y.-Y.C.; investigation, R.Z., J.D. and Y.-Y.C.; resources, S.-M.T.; data curation, K.T.; writing—original draft preparation, R.Z., K.T. and Y.-Y.C.; writing—review and editing, Z.-J.W., K.T. and F.H.; supervision, Z.-J.W. and F.H.; project administration, Z.-J.W. and F.H.; funding acquisition, Z.-J.W. All authors have read and agreed to the published version of the manuscript.

Funding

This research was funded by the National Natural Science Foundation of China (31772680, 32072795).

Institutional Review Board Statement

Not applicable.

Informed Consent Statement

Not applicable.

Data Availability Statement

The raw data in present paper was deposited in NCBI with the submission ID SUB10721255 (http://0-www-ncbi-nlm-nih-gov.brum.beds.ac.uk/bioproject/783962, 27 November 2021).

Conflicts of Interest

The authors declare no conflict of interest.

References

  1. Goldsmith, M.R.; Shimada, T.; Abe, H. The genetics and genomics of the silkworm, Bombyx mori. Annu. Rev. Entomol. 2005, 50, 71–100. [Google Scholar] [CrossRef] [PubMed]
  2. Liu, L.; Zhang, P.; Gao, Q.; Feng, X.; Han, L.; Zhang, F.; Bai, Y.; Han, M.; Hu, H.; Dai, F.; et al. Comparative transcriptome analysis reveals bmo-miR-6497-3p regulate circadian clock genes during the embryonic diapause induction process in bivoltine silkworm. Insects 2021, 12, 739. [Google Scholar] [CrossRef] [PubMed]
  3. Zhu, Z.; Tan, Y.; Xiao, S.; Guan, Z.; Zhao, W.; Dai, Z.; Liu, G.; Zhang, Z. Solitary living brings a decreased weight and an increased agility to the domestic silkworm, Bombyx mori. Insects 2021, 12, 809. [Google Scholar] [CrossRef] [PubMed]
  4. Lu, Z.Y.; Meng, Z.; Wen, M.Y.; Kang, X.; Zhang, Y.; Liu, Q.; Zhao, P.; Xia, Q. Overexpression of BmFoxO inhibited larval growth and promoted glucose synthesis and lipolysis in silkworm. Mol. Genet. Genomics. 2019, 294, 1375–1383. [Google Scholar] [CrossRef] [PubMed]
  5. Zheng, S.F.; Jin, X.; Chen, M.H.; Shi, Q.; Zhang, H.; Xu, S. Hydrogen sulfide exposure induces jejunum injury via CYP450s/ROS pathway in broilers. Chemosphere 2019, 214, 25–34. [Google Scholar] [CrossRef]
  6. Magierowski, M.; Magierowska, K.; Hubalewska-Mazgaj, M.; Surmiak, M.; Sliwowski, Z.; Wierdak, M.; Kwiecien, S.; Chmura, A.; Brzozowski, T. Cross-talk between hydrogen sulfide and carbon monoxide in the mechanism of experimental gastric ulcers healing, regulation of gastric blood flow and accompanying inflammation. Biochem. Pharmacol. 2018, 149, 131–142. [Google Scholar] [CrossRef] [PubMed]
  7. Miller, D.L.; Roth, M.B. Hydrogen sulfide increases thermotolerance and lifespan in Caenorhabditis elegans. Proc. Natl. Acad. Sci. USA 2007, 104, 20618–20622. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  8. Budde, M.W.; Roth, M.B. Hydrogen sulfide increases hypoxia-inducible factor-1 activity independently of von Hippel-Lindau tumor suppressor-1 in C. elegans. Mol. Biol. Cell 2010, 21, 212–217. [Google Scholar] [CrossRef] [Green Version]
  9. Qabazard, B.; Li, L.; Gruber, J.; Peh, M.T.; Ng, L.F.; Kumar, S.D.; Rose, P.; Tan, C.H.; Dymock, B.W.; Wei, F.; et al. Hydrogen sulfide is an endogenous regulator of aging in Caenorhabditis elegans. Antioxid. Redox. Signal. 2014, 20, 2621–2630. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  10. Wei, Y.; Kenyon, C. Roles for ROS and hydrogen sulfide in the longevity response to germline loss in Caenorhabditis elegans. Proc. Natl. Acad. Sci. USA 2016, 113, E2832–E2841. [Google Scholar] [CrossRef] [PubMed] [Green Version]
  11. Kabil, H.; Kabil, O.; Banerjee, R.; Harshman, L.G.; Pletcher, S.D. Increased transsulfuration mediates longevity and dietary restriction in Drosophila. Proc. Natl. Acad. Sci. USA 2011, 108, 16831–16836. [Google Scholar] [CrossRef] [Green Version]
  12. Zhong, J.F.; Wang, S.P.; Shi, X.X.; Mu, L.L.; Li, G.Q. Hydrogen sulfide exposure increases desiccation tolerance in Drosophila melanogaster. J. Insect Physiol. 2010, 56, 1777–1782. [Google Scholar] [CrossRef]
  13. Cao, Y.Y.; Peng, L.L.; Jiang, L.; Thakur, K.; Hu, F.; Tang, S.M.; Wei, Z.J. Evaluation of the metabolic effects of hydrogen sulfide on the development of Bombyx mori (Lepidoptera: Bombycidae), using liquid chromatography-mass spectrometry-based metabolomics. J. Insect Sci. 2020, 20, 4. [Google Scholar] [CrossRef]
  14. Jiang, L.; Peng, L.L.; Cao, Y.Y.; Thakur, K.; Tang, S.M.; Hu, F.; Wei, Z.J. Transcriptome analysis reveals gene expression changes of the fat body of silkworm (Bombyx mori L.) in response to selenium treatment. Chemosphere 2020, 245, 125660. [Google Scholar] [CrossRef]
  15. Martin, M. Cutadapt removes adapter sequences from high-throughput sequencing reads. EMBnet. J. 2011, 17, 10–12. [Google Scholar] [CrossRef]
  16. Pertea, M.; Kim, D.; Pertea, G.M.; Leek, J.T.; Salzberg, S.L. Transcript-level expression analysis of RNA-seq experiments with HISAT, StringTie and Ballgown. Nat. Protoc. 2016, 11, 1650–1667. [Google Scholar] [CrossRef]
  17. Beekman, R.; Chapaprieta, V.; Russiñol, N.; Vilarrasa-Blasi, R.; Verdaguer-Dot, N.; Martens, J.H.A.; Duran-Ferrer, M.; Kulis, M.; Serra, F.; Javierre, B.M.; et al. The reference epigenome and regulatory chromatin landscape of chronic lymphocytic leukemia. Nat. Med. 2018, 24, 868–880. [Google Scholar] [CrossRef]
  18. Xia, Q.; Zhou, Z.; Lu, C.; Cheng, D.; Dai, F.; Li, B.; Zhao, P.; Zha, X.; Cheng, T.; Chai, C.; et al. A draft sequence for the genome of the domesticated silkworm (Bombyx mori). Science 2004, 306, 1937–1940. [Google Scholar]
  19. Kim, D.; Langmead, B.; Salzberg, S.L. HISAT: A fast spliced aligner with low memory requirements. Nat. Methods 2015, 12, 357–360. [Google Scholar] [CrossRef] [Green Version]
  20. Pertea, M.; Pertea, G.M.; Antonescu, C.M.; Chang, T.C.; Mendell, J.T.; Salzberg, S.L. StringTie enables improved reconstruction of a transcriptome from RNA-seq reads. Nat. Biotechnol. 2015, 33, 290–295. [Google Scholar] [CrossRef] [Green Version]
  21. Robinson, M.D.; McCarthy, D.J.; Smyth, G.K. edgeR: A Bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics. 2010, 26, 139–140. [Google Scholar] [CrossRef] [Green Version]
  22. Kanehisa, M.; Araki, M.; Goto, S.; Hattori, M.; Hirakawa, M.; Itoh, M.; Katayama, T.; Kawashima, S.; Okuda, S.; Tokimatsu, T.; et al. KEGG for linking genomes to life and the environment. Nucleic Acids Res. 2008, 36, D480–D484. [Google Scholar] [CrossRef]
  23. Zhang, F.; Zhang, Y.Y.; Ma, R.H.; Thakur, K.; Han, J.; Hu, F.; Zhang, J.G.; Wei, Z.J. Multi-omics reveals the anticancer mechanism of asparagus saponin-asparanin A on endometrial cancer Ishikawa cells. Food Funct. 2021, 12, 614–632. [Google Scholar] [CrossRef]
  24. Wang, J.J.; Bai, W.W.; Zhou, W.; Liu, J.; Chen, J.; Liu, X.Y.; Xiang, T.T.; Liu, R.H.; Wang, W.H.; Zhang, B.L.; et al. Transcriptomic analysis of two Beauveria bassiana strains grown on cuticle extracts of the silkworm uncovers their different metabolic response at early infection stage. J. Invertebr. Pathol. 2017, 145, 45–54. [Google Scholar] [CrossRef]
  25. Adam, G.M.; Naama, K. NADH ties One-Carbon metabolism to cellular respiration. Cell Press. 2020, 31, 660–662. [Google Scholar]
  26. Chen, Q.M.; Ma, Z.G.; Wang, X.; Li, Z.; Zhang, Y.; Ma, S.; Zhao, P.; Xia, Q. Comparative proteomic analysis of silkworm fat body after knocking out fibroin heavy chain gene: A novel insight into cross-talk between tissues. Funct. Integr. Genomics. 2015, 15, 611–637. [Google Scholar] [CrossRef]
  27. Long, W.; Wu, J.; Shen, G.; Zhang, H.; Liu, H.; Xu, Y.; Gu, J.; Jia, L.; Lin, Y.; Xia, Q. Estrogen-related receptor participates in regulating glycolysis and influences embryonic development in silkworm Bombyx mori. Inesct Mol. Biol. 2020, 29, 160–169. [Google Scholar] [CrossRef]
  28. Efe, S.; Duvernell, D.D.; Matzkin, M.L.; Duan, Y.H.; Zhu, C.T.; Verrelli, B.C.; Eanes, W.F. Single-Locus latitudinal clines and their relationship to temperate adaptation in metabolic genes and derived alleles in Drosophila melanogaster. Genet. Soc. Am. 2004, 168, 923–931. [Google Scholar]
  29. Long, L.; Wang, C.C.; Zhao, X.Y.; Guan, J.X.; Lei, C.L.; Huang, Q.Y. Isocitrate dehydrogenase-mediated metabolic disorders disrupt active immunization against fungal pathogens in eusocial termites. J. Pest. Sci. 2020, 93, 291–301. [Google Scholar]
  30. Meyer, S.; Gessner, D.K.; Wen, C.P.; Most, E.; Liebisch, G.; Zorn, H.; Ringseis, R.; Eder, K. The Antisteatotic and Hypolipidemic Effect of Insect Meal in Obese Zucker Rats is Accompanied by Profound Changes in Hepatic Phospholipid and 1-Carbon Metabolism. Mol. Nutr. Food Res. 2019, 63, e1801305. [Google Scholar] [CrossRef]
  31. Shi, A.B.; Chen, C.C.; Banerjee, R.; Glodowski, D.; Audhya, A.; Rongo, C.; Grant, B.D. EHBP-1 functions with RAB-10 during endocytic recycling in Caenorhabditis elegans. Mol. Biol. Cell. 2010, 21, 2930–2943. [Google Scholar] [CrossRef] [Green Version]
  32. Le Roy, C.; Wrana, J.L. Clathrin- and non-clathrin-mediated endocytic regulation of cell signaling. Nat. Rev. Mol. Cell Biol. 2005, 6, 112–126. [Google Scholar] [CrossRef]
  33. Nichols, B. Caveosomes and endocytosis of lipid rafts. J. Cell Sci. 2003, 116, 4707–4714. [Google Scholar] [CrossRef] [Green Version]
  34. Stenmark, H. Rab GTPases as coordinators of vesicle traffic. Nat. Rev. Mol. Cell Biol. 2009, 10, 513–525. [Google Scholar] [CrossRef]
  35. Glodowski, D.R.; Chen, C.C.; Schaefer, H.; Grant, B.D.; Rongo, C. RAB-10 regulates glutamate receptor recycling in a cholesterol-dependent endocytosis pathway. Mol. Biol. Cell. 2007, 18, 4387–4396. [Google Scholar] [CrossRef] [Green Version]
  36. Tian, J.H.; Xue, B.; Hu, J.H.; Li, J.X.; Cheng, X.Y.; Hu, J.S.; Li, F.C.; Chen, Y.H.; Li, B. Exogenous substances regulate silkworm fat body protein synthesis through MAPK and PI3K/Akt signaling pathways. Chemosphere 2017, 171, 202–207. [Google Scholar] [CrossRef]
  37. Zanardo, R.C.; Brancaleone, V.; Distrutti, E.; Fiorucci, S.; Cirino, G.; Wallace, J.L. Hydrogen sulfide is an endogenous modulator of leukocyte-mediated inflammation. FASEB J. 2006, 20, 2118–2120. [Google Scholar] [CrossRef]
  38. Jia, Q.; Wang, T.T.; Wang, X.Y.; Xu, H.; Liu, Y.; Wang, Y.; Shi, Q.; Liang, Q. Astragalin suppresses inflammatory responses and bone destruction in mice with collagen-induced arthritis and in human fibroblast-likesynoviocytes. Front. Pharmacol. 2019, 10, 94. [Google Scholar] [CrossRef] [Green Version]
  39. Hu, W.; Wang, X.; Ma, S.; Peng, Z.; Cao, Y.; Xia, Q. CRISPR-Mediated endogenous activation of fibroin heavy chain gene triggers cellular stress responses in Bombyx mori Embryonic Cells. Insects 2021, 12, 552. [Google Scholar] [CrossRef]
  40. Kunz, R.I.; Brancalhao, R.M.; Ribeiro, L.F.; Natali, M.R. Silkworm sericin: Properties and biomedical applications. Biomed. Res. Int. 2016, 2016, 8175701. [Google Scholar] [CrossRef] [Green Version]
  41. Ni, M.; Zhang, H.; Li, F.C.; Wang, B.B.; Xu, K.Z.; Shen, W.D.; Li, B. Nanoparticulate anatase TiO2 (TiO2NPs) upregulates the expression of silkworm (Bombyx mori) neuropeptide receptor and promotes silkworm feeding, growth, and silking. Peptides 2015, 68, 64–71. [Google Scholar] [CrossRef]
Figure 1. The volcano plot of DEGs. The red points represent the up-regulated genes, the blue points represent the down-regulated genes, and the grey points represent the genes without differential expression.
Figure 1. The volcano plot of DEGs. The red points represent the up-regulated genes, the blue points represent the down-regulated genes, and the grey points represent the genes without differential expression.
Insects 12 01110 g001
Figure 2. GO classification of DEGs (p ≤ 0.01).
Figure 2. GO classification of DEGs (p ≤ 0.01).
Insects 12 01110 g002
Figure 3. The KEGG classification of DEGs. The x-axis represents the number of genes annotated into the pathway and the proportion of the number of DEGs annotated to the total number of genes. The y-axis represents the name of the enriched KEGG pathways.
Figure 3. The KEGG classification of DEGs. The x-axis represents the number of genes annotated into the pathway and the proportion of the number of DEGs annotated to the total number of genes. The y-axis represents the name of the enriched KEGG pathways.
Insects 12 01110 g003
Figure 4. The qRT-PCR validation for DEGs expression after exposure to H2S. DGEs related to the glycolysis pathway (A), the TCA cycle (B), endocytosis (C), and inflammation and silk protein synthesis (D) were analyzed. The full names of abbreviations of genes were listed in Table S1.
Figure 4. The qRT-PCR validation for DEGs expression after exposure to H2S. DGEs related to the glycolysis pathway (A), the TCA cycle (B), endocytosis (C), and inflammation and silk protein synthesis (D) were analyzed. The full names of abbreviations of genes were listed in Table S1.
Insects 12 01110 g004
Figure 5. Expression changes of 10 DEGs after H2S exposure at different time points in the 5th instar larvae (5L1D, 5L2D, 5L3D, 5L4D, and 5L5D). (A) PGK; (B) PGM; (C) IDH; (D) PC; (E) Rab-10; (F) Tret1; (G) Tak1; (H) MPP-3; (I) Fib-H; (J) Fib-L; and (K) P25. The lowercase letters (a, b, c, and d) represent the significant difference at p ≤ 0.05 among the expression levels at different times of the H2S treated group; the uppercase letters (A, B, C, and D) represent the significant difference at p ≤ 0.05 among the expression levels at different times of the control group. The full names of abbreviations of genes are listed in Table S1.
Figure 5. Expression changes of 10 DEGs after H2S exposure at different time points in the 5th instar larvae (5L1D, 5L2D, 5L3D, 5L4D, and 5L5D). (A) PGK; (B) PGM; (C) IDH; (D) PC; (E) Rab-10; (F) Tret1; (G) Tak1; (H) MPP-3; (I) Fib-H; (J) Fib-L; and (K) P25. The lowercase letters (a, b, c, and d) represent the significant difference at p ≤ 0.05 among the expression levels at different times of the H2S treated group; the uppercase letters (A, B, C, and D) represent the significant difference at p ≤ 0.05 among the expression levels at different times of the control group. The full names of abbreviations of genes are listed in Table S1.
Insects 12 01110 g005
Table 1. The sequencing data of samples.
Table 1. The sequencing data of samples.
SamplesFB_Control_1FB_Control_2FB_Control_3FB_H2S_1FB_H2S_2FB_H2S_3
Raw reads52,598,69249,647,10450,987,48850,804,42453,068,04251,108,946
(7.89 Gb)(7.45 Gb)(7.65 Gb)(7.62 Gb)(7.96 Gb)(7.67 Gb)
High quality reads49,443,62449,144,22850,029,92650,357,63652,490,30250,286,178
(7.42 Gb)(7.37 Gb)(7.50 Gb)(7.55 Gb)(7.87 Gb)(7.54 Gb)
High quality reads ratio (%)94.0098.9998.1299.1298.9198.39
Q30 (%)98.5399.0899.0798.9398.9899.12
GC content (%)5048494848.5048
Total mapped46,077,04947,843,35448,788,61449,023,73051,140,81749,141,008
Table 2. The main enriched KEGG pathways (p < 0.05) in the H2S-treated and the control groups.
Table 2. The main enriched KEGG pathways (p < 0.05) in the H2S-treated and the control groups.
Pathway IDPathway NameNumber of DEGsp-Value
ko03010Ribosome33<0.01
ko04510Focal adhesion20<0.01
ko05132Salmonella infection7<0.01
ko04370VEGF signaling pathway6<0.01
ko04810Regulation of actin cytoskeleton15<0.01
ko00360Phenylalanine metabolism7<0.01
ko04530Tight junction13<0.01
ko04668TNF signaling pathway40.01
ko00010Glycolysis/Gluconeogenesis110.01
ko04144Endocytosis350.01
ko00562Inositol phosphate metabolism120.01
ko00130Ubiquinone and other terpenoid-quinone biosynthesis60.01
ko04380Osteoclast differentiation30.01
ko04670Leukocyte transendothelial migration30.01
ko00410beta-Alanine metabolism70.01
ko05100Bacterial invasion of epithelial cells40.02
ko04611Platelet activation30.02
ko04520Adherens junction70.02
ko04621NOD-like receptor signaling pathway40.02
ko04391Hippo signaling pathway-fly190.02
ko04010MAPK signaling pathway160.03
ko00340Histidine metabolism40.04
ko05410Hypertrophic cardiomyopathy (HCM)60.04
ko04064NF-kappa B signaling pathway20.04
ko00020Citrate cycle (TCA cycle)90.04
ko04962Vasopressin-regulated water reabsorption40.04
Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Share and Cite

MDPI and ACS Style

Zhang, R.; Cao, Y.-Y.; Du, J.; Thakur, K.; Tang, S.-M.; Hu, F.; Wei, Z.-J. Transcriptome Analysis Reveals the Gene Expression Changes in the Silkworm (Bombyx mori) in Response to Hydrogen Sulfide Exposure. Insects 2021, 12, 1110. https://0-doi-org.brum.beds.ac.uk/10.3390/insects12121110

AMA Style

Zhang R, Cao Y-Y, Du J, Thakur K, Tang S-M, Hu F, Wei Z-J. Transcriptome Analysis Reveals the Gene Expression Changes in the Silkworm (Bombyx mori) in Response to Hydrogen Sulfide Exposure. Insects. 2021; 12(12):1110. https://0-doi-org.brum.beds.ac.uk/10.3390/insects12121110

Chicago/Turabian Style

Zhang, Rui, Yu-Yao Cao, Juan Du, Kiran Thakur, Shun-Ming Tang, Fei Hu, and Zhao-Jun Wei. 2021. "Transcriptome Analysis Reveals the Gene Expression Changes in the Silkworm (Bombyx mori) in Response to Hydrogen Sulfide Exposure" Insects 12, no. 12: 1110. https://0-doi-org.brum.beds.ac.uk/10.3390/insects12121110

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Metrics

Back to TopTop