Skip to main content

ORIGINAL RESEARCH article

Front. Physiol., 24 January 2022
Sec. Invertebrate Physiology
This article is part of the Research Topic Advances on the Physiology and Cell Biology of Invertebrate Parasites View all 6 articles

Transcriptome Analysis for Identification of Genes Related to Growth and Development, Digestion and Detoxification, Olfaction in the Litchi Stink Bug Tessaratoma papillosa

\r\nLin Cheng&#x;Lin Cheng1†Shuncai Han&#x;Shuncai Han1†Jingtao JiangJingtao Jiang1Haichao Li*Haichao Li2*Lingfei Peng*\r\nLingfei Peng1*
  • 1Biological Control Research Institute, Fujian Agriculture and Forestry University, China Fruit Fly Research and Control Center of FAO/IAEA, Key Laboratory of Biopesticide and Chemical Biology, Ministry of Education, State Key Laboratory of Ecological Pest Control for Fujian and Taiwan Crops, Fuzhou, China
  • 2Chinese Academy of Sciences Key Laboratory of Insect Developmental and Evolutionary Biology, CAS Center for Excellence in Molecular Plant Sciences, Shanghai Institute of Plant Physiology and Ecology, Chinese Academy of Sciences, Shanghai, China

Tessaratoma papillosa is a major pest of Litchi chinensis and Dimocarpus longan. Adult and nymph secretions are not only harmful to plants but also to humans. At present, there are not a lot of research on T. papillosa, especially omics research. We used high-throughput sequencing technology to sequence the T. papillosa transcriptome and obtained 67,597 unigenes homologous to Halyomorpha halys (88.03%). Subsequently, RNA-SEQ and comparative analyses were performed on the 14 different developmental stages and tissues. A total of 462 unigenes related to growth and development, 1,851 unigenes related to digestion and detoxification, and 70 unigenes related to olfaction were obtained. Moreover, expression analysis showed that the T. papillosa major life activities genes are uniformly expressed across all developmental states. However, the adult midgut gene expression patterns were utterly different from that of the nymphs. Similarly, female fat body genes exhibited distinct expression patterns compared to that of males and nymphs. Thus, different developmental stages and physiological functions affect gene expression patterns. We also found that most of the differential genes were associated with cellular maintenance. This study will help understand the growth and development of litchi stink bugs, their choice of host plants, food digestion and detoxification, and their reproductive behavior. In addition, this result can provide reference information for some target genes in the process of control of T. papillosa.

Introduction

The selection of host plants by herbivorous insects relies on sensory signals, including vision, smell, taste, and touch. These signals enable the insects to distinguish the plants’ physical properties and accurately identify the secondary plant metabolites and subsequent detoxification (Mollo et al., 2017). Furthermore, the insects locate their host plants through chemotaxis. The plants release signal stimuli which are then detected by the insects (Finch and Collier, 2000; Yang et al., 2001). For example, the insect odor binding protein senses and transmits odor substances and regulates oviposition. The larvae of Helicoverpa armigera and Helicoverpa assulta have the habit of feeding on each other. Meanwhile, male moths transmit the OBP10 carrying 1-dodecene to female moths through mating. After mating, the 1-dodecene is then located on the egg surfaces laid on host plants. Consequently, moths detect this compound through the OBP10 in their antennae and avoid the plants containing the eggs; this mechanism also helps prevent cannibalism (Sun et al., 2012) at the same time, to counter the effects of the plant metabolites, insects have evolved various physiological and biochemical adaptation mechanisms. The emergence of a series of detoxification enzymes is one of the defense mechanisms developed by insects against plant metabolites (Oakeshott et al., 2003; Clarkson et al., 2018). A few insects store toxic substances from plants to avoid poisoning and enhance their food utilization efficiency and selection of host plants. Additionally, the intrinsic factors regulating this mechanism mainly depend on the abundance and expression patterns of the insects’ digestive and detoxification enzyme genes (Yang et al., 2001). These genes play a vital role in the growth and development of insects, and they are also important target genes in pest control.

Obtaining the omics data of insects is useful for understanding the growth and development of insects and screening targeted control targets. Rapid development of high-throughput sequencing technologies has greatly improved research on non-model insects, making the acquisition of their genetic, physiological, and ecological information easier. For instance, Wu et al. (2017) used transcriptome sequencing technology to analyze the T. papillosa olfaction gene families. Moreover, targeted acquisition of pest control strategies based on molecular targeting technology offers a viable alternative approach to mitigating pests (Delye et al., 2020; Schmidt et al., 2021).

Tessaratoma papillosa is an important pest of longan and lychee, mainly distributed in Southeast Asia (Thailand and Vietnam), India and southern China (Guangdong, Fujian, Guangxi, and Taiwan) (Li et al., 2014; Tran et al., 2019; Wu et al., 2020). T. papillosa is also an important pest for greening trees. At the same time, when exposed to external stimuli, it releases toxic liquids, irritating human skin and eyes, and causing diseases (Zhang et al., 2009). It is also the vector insect of soot disease and Longan witches broom-associated virus (LWBD) (Chen et al., 2001; Meng et al., 2017). T. papillosa has a relatively long life cycle, and the entire life cycle can be harmful, which brings great difficulties to its prevention and control (Figure 1A). At present, chemical pesticides are not effective in the prevention and control of T. papillosa, and there are fewer natural enemies of T. papillosa, which makes prevention and control difficult. Recent advancements in RNAi technology have, however, made the control of pathogens such as T. papillosa possible (Jain et al., 2021).

FIGURE 1
www.frontiersin.org

Figure 1. Morphology of T. papillosa at different developmental stages and organization Chart. (A) T. papillosa Life cycle. (B) reproductive system of female adult. (C) Digestive tract. (D) Internal organs. (E) reproductive system of male adult.

In this study, the deep transcriptome sequencing technology is used to explore the growth and development-related genes in different developmental stages and tissues of T. papillosa. The findings might be helpful in understanding T. papillosa and other stink bugs survival traits, such as spawning and reproductive cycles, choice of host plants, and food digestion mechanisms.

Materials and Methods

Insects and Samples

Approximately 230 T. papillosa bugs were collected on the wild Lychee trees from four Lychee orchards in Fuzhou, Fujian Province, China. All the bugs were bring back to the lab alive, then the wings, legs, scutellum, and abdominal tergum were cut off with scissors, the abdominal tergum were removed carefully, and the viscera were exposed. The body was fixed with insect pin on a wax plate, and dissected under a stereo-microscope. We collected 14 samples including female and male adults, nymph of different instars, antennae of adults and nymphs, ovaries and testis, fat body (female and male adult, nymphs), midgut (female and male adult, nymphs), and hemolymph. Each samples were mixed with five bugs tissues or organs. After mixing, the samples were immediately frozen in liquid nitrogen for 10 min and then stored at −80°C for RNA extraction. For each samples, experiments were performed in triplicates (Supplementary Table 1).

RNA Extraction, Library Construction, and Sequencing

RNA was extracted from the 14 samples using the TRIzol extraction method and then purified using the RNeasy® MinElute® cleanup kit (Qiagen®). RiboZeroTM Magnetic (Plant leaf) kit (Epicentre®) was used to degrade rRNA from 4 μg of each sample. The removal of rRNA was then verified using a 2100 BioAnalyzer (Agilent Technologies).

After total RNA extraction, mRNA was enriched using oligo dT magnetic beads, and the fragmentation buffer was used to break the mRNA into shorter fragments. Messenger RNA was utilized as a template for the first cDNA strand synthesis using a six-base random primer. Subsequently, dNTP, RNase H, and DNA polymerase I were used for subsequent cDNA strands synthesis. End-repair was then conducted after purification by the QIAquick PCR kit and elution with EB plus foya buffer. Eventually, gel electrophoresis was done for fragment selection, followed by PCR amplification and the constructed library sequenced using BGI-500 (BGI), three replicates per sample.

Transcriptome Assembly, Gene Annotation, and Ontology

We used the short reads assembly software SOAPdenovo (Li et al., 2010) to perform a de novo transcriptome assembly. SOAPdenovo concatenates reads with a certain length of overlap into longer N-free fragments called Contigs. We then compared paired-end reads between the obtained reads and the contigs to determine the variant contigs that originated from the same transcript and the distance between these contigs. Usually, SOAPdenovo links these contigs together, and the intermediate unknown sequences are denoted by N, so that scaffold. The paired-end reads were then used to generate the scaffold with the least N content and shorter extension at both ends. This scaffold was then referred to as unigene. If multiple samples from the same species are sequenced, the unigene can be further sequenced and de-duplicated by the sequence clustering software to obtain the longest and non-redundant consensus unigene. Finally, the unigene sequence was aligned with sequences from the protein databases, including NT, Swiss-Prot, KEGG, and COG databases (e-value < 0.00001), and the best-aligned protein was determined by the Unigene sequence orientation. In case of contradicting comparison results between different libraries, the unigene sequence direction was obtained according to the priority Scores by the NR, Swiss-Prot, KEGG, and COG databases. We used the software ESTScan (Iseli et al., 1999) to predict the coding region and determine the sequence direction for a unigene sequence that was incomparable with the four databases. Moreover, the unigene whose sequence direction was, determined was read in the 5′ to the 3′ direction, whereas the unigene with undetermined sequence direction was read based on the denotations generated by the assembly software.

Sequence Analysis

Differential expression analysis was performed using DESeq R package (1.10.1), which is based on a negative binomial distribution model. DESeq is a statistical program that determines the differential expression in digital gene expression data using a model-based negative binomial distribution (Sun M. et al., 2020), Meanwhile, the P values were adjusted using Benjamini and Hochberg’s approach to account for the false discovery rate (FDR). Genes with a twofold difference (the absolute value of log2 ratio ≥ 1) between two stages and an FDR < 0.01 were categorized as differentially expressed genes (DEGs).

qRT-PCR Assays

Total RNA was isolated from the collected samples using the TRIzol reagent (Invitrogen). The RNA served as the template for synthesizing cDNA using the First Strand cDNA Synthesis Kit (Toyobo). The qRT-PCR analysis was performed using the SYBR Green Real-time PCR Master Mix Kit (Toyobo). Details regarding all primers are listed in Supplementary Table 7. Melting-curve analyses were performed for all of the primers. To normalize Ct values obtained for each gene, actin-2 expression levels were used. RT-qPCR was carried out using Mastercycler® ep realplex (Eppendorf). All qRT-PCR assays were repeated three times. qRT-PCR reactions and data were analyzed according to the methods of Livak and Schmittgen (2001) and Bustin et al. (2009). The data were analyzed using a one-way analysis of variance (ANOVA) to look for treatment effects compared with the untreated control.

Results

T. papillosa Midgut Morphology

Mitgut of T. papillosa contains four parts (four stomachs), the first part was the largest one, brown in color and pouched in shape; the second part was short and tubular; the third part was spherical; the fourth part was the longest one, which is up to two times of the body length, and can be separated into two parts, one pouch and followed by a long thin tube (Figure 1B). The fat body was milky white, resembled grapes in shape (Figure 1C). Female adults had a pair of ovaries, each ovary composed of seven ovarioles (Figure 1D) and male adults possessed a pair of red kidney-shaped testicles with deferens ducts (Figure 1E).

T. papillosa Transcriptome Analysis

A total of 7.29 Gb data was obtained. Additionally, the normalized data were regarded as clean if the reads had Q20 of 96.8%, Q30 of 88.11%, and the ratio of 89.19%. Clean reads were assembled, and 67,597 unigenes were obtained after de-redundancy. The total length, average length, N50, and GC content were 74, 010, 888 bp, 1,094 bp, 2,016 bp, and 35.13%, respectively. Moreover, unigene length ranged between 200 and 16,160 bp, with an average of 200–300 bp (16, 754). There were 5,224 unigenes with 3,000 or more base pairs (Supplementary Figure 1A).

Eventually, 32,673 unigenes accounting for 48.33%, 16,901 unigenes accounting for 25.00%, 24,288 unigenes accounting for 35.93%, 23,401 unigenes accounting for 34.62%, 25,225 unigenes accounting for 37.32%, 4,403 unigenes accounting for 6.51%, and 24,097 unigenes accounting for 35.65% were deposited to the NR, NT, SwissProt, KOG, KEGG, GO, and Pfam databases, respectively (Figure 2A). Transdecoder software was then used to detect 27, 031 CDS (Supplementary Figure 1B). Annotated gene species similarity results showed that T. papillosa and Halyomorpha halys had more similar genes (88.05%), followed by Cimex lectularius (2.33%), Riptortus pedestris (1.00%), Nilaparvata lugens (0.59%), and then Zootermopsis nevadensis (0.32%) (Figure 2B).

FIGURE 2
www.frontiersin.org

Figure 2. Overall analysis of transcriptome data of T. papillosa. (A) Five gene function annotation database annotation results Venn diagram. (B) NR annotated species distribution. (C) The abscissa indicates the number of corresponding transcripts, and the ordinate indicates the functional classification of KEGG. (D) The abscissa indicates the number of corresponding transcripts, and the ordinate indicates the GO function classification.

KEGG enrichment analysis of all unigenes and obtained 45 pathways, with the ubiquitin-mediated proteolysis being abundantly represented (Figure 2C). The five pathways with the highest number of unigene enrichment hits were (1) Global and overview maps (3,509), (2) Signal transduction (3,404), (3) Cancers: Overview (2,399), (4) Transport and catabolism (2,055), and (5) Endocrine system (1,915). Additionally, 10,442 unigenes were assigned to the three GO categories: biological processes, cellular components, and molecular functions (Figure 2D). In the biological processes category, level 2 terms cellular processes (1,216) and biological regulation (548) were the most frequently represented, whereas, in the cellular component category, level 2 terms cellular processes (1,145) and membrane (1,100) were the predominant categories. Meanwhile, binding (2,129) and catalytic activities (1,890) were the most abundant processes in the molecular function category. During the comparison and analysis of all unigenes in the Animal TFDB 2.0 database, 4,101 transcription factors were annotated. Six transcription factors families [zf-C2H2(862), RHD(275), TF_bZIP(252), Homeobox(222), HGM(209), and MYB(186)] with the most significant number of unigenes accounted for 48.91% of all transcription factors (Supplementary Figure 2).

Analysis of Differentially Expressed Genes in Different Developmental Stages and Tissues of T. papillosa

RNA-SEQ sequencing obtained reads with than 21M (number of bases > 1 Gb) from each set. Both Q20 and Q30 were more than 90% (Supplementary Table 2). The proportion of sequences mapped to the transcriptome data was between 88.92 and 95.57%, and the unique, clean reads were 36.20–56.91% (Supplementary Table 3). Meanwhile, pairwise comparison results showed that the correlation coefficient between the samples was between 0.1907 and 0.9708 (Supplementary Figure 3A). Moreover, PCA analysis showed that the nymph, female and male midguts formed a separate cluster from the female ovary and the fat body cluster. The nymph, male antennae, and testis, however, clustered separately. The results were consistent with the physiological functions of the genes expressed in different body parts, depicting significant differences in the functional structure between different developmental stages tissues (ANOSIM, P < 0.001). In addition, t gene expression had a normal distribution, indicating the occurrence of a normal gene expression distribution state (Supplementary Figure 3C). Considering the dynamic changes associated with gene expression, the co-expression between genes was calculated, and the gene expression regulation model containing essential genes at different developmental stages or under different conditions was searched. Weighted Correlation Network Analysis (WGCNA) results showed that the correlation between genes within the same module was very high but much lower between genes in different modules (Supplementary Figure 3D). Different T. papillosa samples exhibited a direct correlation between the gene expression patterns and the physiological function of tissues during T. papillosa development. According to PCA analysis, the nymph, female and male midguts formed a separate cluster from the female ovaries and fat body. Moreover, the correlation coefficient analysis showed that samples with large correlation variations have different physiological functions. Expression analysis also showed that most genes were expressed across the different developmental states.

Comparison between different tissue parts and developmental stages showed that the commonly expressed genes were more than 50% of all transcripts obtained. There were 40,826 commonly expressed genes among male and female adults, with 4,310 being female specific genes and 15,512 being male genes (Figure 3A). The nymphs had 37,960 widely expressed genes, among which the least specifically expressed genes were only 1,161 (Figure 3D). Moreover, the nymphs′ antennae specifically express genes (8,005) were more than adult genes (3,291) (Figure 3B). Similarly, testes exhibited more expressed genes (55,753) than the ovaries (38,516). The testes and ovaries had the highest number, 32,847 differential genes (up-regulated 18,487, down-regulated 14,360) (Figures 3C, 4), which accounted for 85.28% of ovary unigenes and 58.96% of the number of testes unigenes. This finding indicated that as much as some of the testes and ovaries-related genes maintain their morphological and physiological functions, the expression patterns of other genes are entirely different. Consistently, comparing the midgut and the fat body of L. chinensis at different developmental stages showed that there are fewer specifically expressed genes in females than in nymphs and male adults (midgut 1,867 and Fat body 2,252) (Figures 3E,F) and that the differentially expressed genes in the male and female midguts were the least (4,872). Furthermore, in comparison with nymphs, female adults had 10,493 differentially expressed genes. The differentially expressed genes in adult and larval antennae were 10,953. However, the fat body and nymph comparison indicated that the differentially expressed genes in males were 9,321 and 10,528 in males and females, respectively (Figure 4).

FIGURE 3
www.frontiersin.org

Figure 3. Venn diagram of the number of differential genes between different tissue samples. (A) Comparative analysis of female (TP_FM) and male adult (TP_FM). (B) Comparative analysis of nymph antennae (TP_AA) and adult antennae (TP_FMA). (C) Comparative analysis of ovaries (TP_MO) and testis (TP_FMT). (D) Comparative analysis of nymphs (TP_A), female adult (TP_FM) and male adult (TP_M). (E) Comparative analysis of nymphs midgut (TP_AM), female adult midgut (TP_FMM) and male adult midgut (TP_MM). (F) Comparative analysis of nymphs fat body (TP_AF), female adult fat body (TP_FMF) and male adult fat body (TP_MF).

FIGURE 4
www.frontiersin.org

Figure 4. Analysis of differentially expressed genes (DEGs) at different development stage. The X axis represents the difference comparison scheme for each group, and the Y axis represents the corresponding DEG number. Red represents the number of DEGs that are up-regulated, and blue represents the number of DEGs that are down-regulated.

Differential Analysis of Genes Related to T. papillosa Growth and Development

We compared the expression levels of the obtained genes and found that the midgut genes showed a completely different expression pattern from those in the other tissues. In the adult midgut, there were many highly expressed genes (>50%), which was contrary to those observed for nymphs. However, the Nymphs’ antennae, testis, and ovaries had a relatively large number of highly expressed genes (>10%). In the fat body, the highly expressed genes of the three insect states were about 5% but differed between nymphs and males compared to females (Figure 5).

FIGURE 5
www.frontiersin.org

Figure 5. Cluster analysis of all genes among different tissue samples and different developmental stage.

Upon analysis, 462 growth and development unigenes were obtained, including 58 trypsin, 12 juvenile hormones, 11 ecdysones, 11 chitin, 41 chitinase, 110 cadherin, 7 cuticular proteins, 117 Fatty acdis, and 95 Dynein heavy chains encoding 160 genes (Supplementary Table 4). Additionally, female and male adult worms exhibited noticeable genetic differences between the testes and ovaries. Compared with female adults, nymphs had many silenced genes. Nonetheless, nymph, female and male fat bodies expressed different genes (Figure 6A).

FIGURE 6
www.frontiersin.org

Figure 6. Cluster analysis of differential genes among different tissue samples and different developmental stage. (A) Genes related to growth and development. (B) Genes of digestion and detoxification enzymes. (C) Genes of olfactory.

We obtained 1,851 unigenes related to digestion and detoxification. These included 108 Acetyltransferase, 11 Alcohol dehydrogenase, 104 Aminopeptidase, 71 Carboxypeptidase, 95 NADH, 285 Cytochrome, 162 Cytochrome P450, 38 glutathione s-transferase, 57 serine protease, 20 carbonic anhydrases (Carbonic anhydrase), 195 ATPase (ATPase), 56 carboxylases (Carboxylase), 8 acyl-CoA dehydrogenase, 74 cathepsins, 560 protein kinase (protein kinase, 3 aldehyde oxidase, and 4 chymotrypsin, encoding 611 genes (Supplementary Table 5). Moreover, we also observed significant genetic differences between male and female adults. Different genes were also reported in the reproductive system, antennae, and fat body (Figure 6B).

Seventy olfactory unigenes, including 28 odorant-binding proteins (OBPs), 6 chemosensory proteins (CSPs), and 4 pheromone-binding proteins (PBPs), 19 odorant receptors (OR), 11 sensory neuron membrane proteins (SNMPs), and 3 aldehyde oxidase (aldehyde oxidase) unigenes were obtained. Each of the 70 unigenes encoded 10 additional unigenes (Supplementary Table 6). Moreover, cluster analysis of these genes revealed that most olfactory genes are expressed in the testis and ovaries, but with different expression patterns indicating that the olfactory gene expression is regulated during the early developmental stages. Similarly, some olfactory-related genes were also found in the insect midgut (Figure 6C).

Based on the differential expression of different tissues in the RNA-SEQ data, we screened 20 genes with significant differences from the digestion and detoxification and olfactory related genes, and used RT-qPCR for gene expression verification (Table 1). The RT-qPCR results were consistent with those obtained through RNA-seq and RNA-SEQ analysis (Figure 7).

TABLE 1
www.frontiersin.org

Table 1. qRT-PCR verification of 20 differentially expressed genes.

FIGURE 7
www.frontiersin.org

Figure 7. Comparison of qRT-PCR expression level and FPKM expression level of 20 genes in T. papillosa. The white column is the qRT-PCR data, and the black column is the RNA-SEQ data.

Discussion

During the feeding process, insects take up some secondary plant metabolites together with the nutrients, and these metabolites eventually affect the normal physiological processes of insects (Pan et al., 2016; Yuan et al., 2020). T. papillosa is a major pest of longan and lychee trees in Southeast Asia. The infestation results in severe damage to the plants involved. The results show that the midgut of T. papillosa has four parts, which is obviously more developed than the midgut system of other insects. At the same time, it also provides the possibility to provide more physical storage. Transcriptome sequencing technology plays a vital role in non-model insect gene discovery and gene function studies. The technology has enabled various species (Bayega et al., 2021; Zhan et al., 2021; Zhao et al., 2021). We performed transcriptome and RNA-SEQ sequencing on the T. papillosa genome at different developmental stages and different tissues. This study show that T. papillosa and H. halys had the highest copies of unigenes. H. halys is a polyphagous pest native to East Asia (Sparks et al., 2020; Sun D. D. et al., 2020). Unlike T. papillosa, H. halys exhibits a wide host range. This may be related to the four stomachs of the T. papillosa.

The insects’ detoxification enzymes mainly include cytochrome P450-dependent monooxygenases (P450s), glutathione-S-transferases (GSTs), and carboxylesterases (COE). There are about 150 members detoxification enzymes in insects, and the P450 family members are twice as many as the members of the other two families (Oakeshott et al., 2003; Yuan et al., 2020). Our results also showed that the number and the expression of Cytochrome P450, carboxypeptidase, protein kinase, and other protein family genes were distinct for different physiological periods, especially in the fat body. These genes were closely related to insects’ tolerance and detoxification (Oakeshott et al., 2003; Clarkson et al., 2018). This phenomenon alludes to that the T. papillosa has high food utilization efficiency, enhanced degradation capacity of exogenous toxic substances, and increased stress tolerance. Moreover, the reduced food spectrum of the T. papillosa also enhances its physiological survival.

Seventy unigenes related to the olfactory-related proteins were obtained, such as OBP and OR, which were lesser than those observed by Wu et al. (2017) in T. papillosa antennae analysis (92) (Wu et al., 2017). In other studies, 110 Apolygus lucorum antennal transcriptome-related genes (An et al., 2016) and 60 H. halys OBPs genes (Paula et al., 2016) were identified. Similarly, the same genes were also highly expressed in different tissues and developmental stages of T. papillosa, thus indicating their significance in T. papillosa developmental process.

Feeding, reproduction, development, and detoxification processes by T. papillosa are physiological responses induced by the complex regulation of multiple genes. Since the insect feeding pattern is regulated by the insect olfactory-related genes, the number of related genes can affect the range of insects’ feeding habits (Li et al., 2013). Concomitantly, stimulation of the environmental factors elicits signals which trigger gene expression, inducing morphological changes of the insects through a series of regulatory processes. The morphology is maintained by a series of genes, which are a result of long-term evolutionary selection. Therefore, most of the genes obtained at different developmental stages through mixed sampling are maintenance-oriented, and only the olfactory, detoxification, and development-related genes of the T. papillosa can be utilized at certain stages for preliminary identification and analysis.

The genes with less expression variation are thought to be housekeeping genes associated with basic life activities. Gene expression patterns in the adult midguts were completely different from that of nymphs. Moreover, female fat body genes also displayed distinct expression patterns from those of males and nymphs. These differences may indicate the various physiological activities occurring at different developmental stages. Therefore, the analysis of the physiological functions of the obtained genes provides insights into the development of alternative pest control systems through the molecular screening basis of target genes.

Conclusion

This study reported the transcriptome data of different developmental stages and different tissues of T. papillosa. Through comparative analysis, some genes related to diet, development and host selection were obtained during the entire development process. Through anatomy, it was found that the midgut of the T. papillosa consists of four parts (four stomachs). At the same time, the digestive and detoxification enzyme genes in the T. papillosa are relatively abundant compared with other types of genes. The findings of these results can provide information for the study of the T. papillosa.

Data Availability Statement

The original contributions presented in the study are publicly available. This data can be found here: PRJNA765492.

Author Contributions

LC: formal analysis, investigation, and data curation. SH: formal analysis and investigation. JJ: investigation. HL: conceptualization, formal analysis, data curation, and writing—review and editing. LP: conceptualization, data curation, and writing—review and editing. All authors read and approved the final manuscript.

Funding

This work was supported by Shanghai Agriculture Applied Technology Development Program (2019-02-08-00-12-F01125), the National Natural Science Foundation of China (31772520), and open fund of Fujian Provincial Key Laboratory of Insect Ecology, Fujian Agriculture and Forestry University (FIE201703).

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We thank MogoEdit (http://www.mogoedit.com/) for editing the English text of a draft of this manuscript.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fphys.2021.774218/full#supplementary-material

Supplementary Figure 1 | Distribution map of unigene transcript length and predicted CDS length. (A) The abscissa represents the transcript length interval, and the ordinate represents the number of corresponding transcripts. (B) The abscissa represents the CDS length interval, and the ordinate represents the number of the corresponding CDS.

Supplementary Figure 2 | The number of genes of each transcription factor families.

Supplementary Figure 3 | (A) Heat map of correlation coefficient between samples. Both X and Y axes represent each sample. The color represents the correlation coefficient (the darker the color, the higher the correlation, and the lighter the color, the lower the correlation). (B) Principal component PCA analysis between samples. (C) Express volume box plot. The X-axis is the sample name, and the Y-axis is log10 (FPKM+1). The box plot of each region corresponds to five statistics (from top to bottom, the upper limit, upper quartile, median, and lower quartile). Number of digits, lower limit, the upper and lower limits do not consider outliers). (D) WGCNA gene module detection. In the figure, a branch of the clustering tree represents a gene, and genes with highly similar expressions are grouped together and belong to the same module (indicated by the same color). The gene correlation between each module is displayed in yellow and red progressive colors. The darker the color (red), the higher the correlation between the two genes.

References

An, X. K., Sun, L., Liu, H. W., Liu, D. F., Ding, Y. X., Li, L. M., et al. (2016). Identification and expression analysis of an olfactory receptor gene family in green plant bug Apolygus lucorum (Meyer-Dür). Sci. Rep. 6:37870. doi: 10.1038/srep37870

PubMed Abstract | CrossRef Full Text | Google Scholar

Bayega, A., Oikonomopoulos, S., Gregoriou, M. E., Tsoumani, K. T., Giakountis, A., Yu, C. W., et al. (2021). Nanopore long-read RNA-seq and absolute quantification delineate transcription dynamics in early embryo development of an insect pest. Sci. Rep. 11:7878. doi: 10.1038/s41598-021-86753-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Bustin, S. A., Benes, V., Garson, J. A., Hellemans, J., Huggett, J., Kubista, M., et al. (2009). The MIQE guidelines: minimum information for publication of quantitative real-time PCR experiments. Clin. Chem. 55, 611–622. doi: 10.1373/clinchem.2008.112797

PubMed Abstract | CrossRef Full Text | Google Scholar

Chen, J. Y., Chen, J. Y., and Xu, X. (2001). Advances in research of longan witches’ broom disease. Acta Hortic. 558, 413–416. doi: 10.17660/ActaHortic.2001.558.66

PubMed Abstract | CrossRef Full Text | Google Scholar

Clarkson, C. S., Temple, H. J., and Miles, A. (2018). The genomics of insecticide resistance: insights from recent studies in African malaria vectors. Curr. Opin. Insect. Sci. 27, 111–115. doi: 10.1016/j.cois.2018.05.017

PubMed Abstract | CrossRef Full Text | Google Scholar

Delye, C., Michel, S., Pernin, F., Gautier, V., Gislard, M., Poncet, C., et al. (2020). Harnessing the power of next-generation sequencing technologies to the purpose of high-throughput pesticide resistance diagnosis. Pest. Manag. Sci. 76, 543–552. doi: 10.1002/ps.5543

PubMed Abstract | CrossRef Full Text | Google Scholar

Finch, S., and Collier, R. H. (2000). Host-plant selection by insects - a theory based on ‘appropriate/inappropriate landings’ by pest insects of cruciferous plants. Entomol. Exp. Appl. 96, 91–102. doi: 10.1046/j.1570-7458.2000.00684.x

CrossRef Full Text | Google Scholar

Iseli, C., Jongeneel, C. V., and Bucher, P. (1999). ESTScan: a program for detecting, evaluating, and reconstructing potential coding regions in EST sequences. Proc. Int. Conf. Intell. Syst. Mol. Biol. 99, 138–148. doi: 10.5555/645634.660818

CrossRef Full Text | Google Scholar

Jain, R. G., Robinson, K. E., Asgari, S., and Mitter, N. (2021). Current scenario of RNAi-based hemipteran control. Pest Manag. Sci. 77, 2188–2196. doi: 10.1002/ps.6153

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, D. S., Liao, C. Y., Zhang, B. X., and Song, Z. W. (2014). Biological control of insect pests in litchi orchards in China. Biol. Control 68, 23–36. doi: 10.1016/j.biocontrol.2013.06.003

CrossRef Full Text | Google Scholar

Li, H. C., Zhang, H., Guan, R. B., and Miao, X. X. (2013). Identification of differential expression genes associated with host selection and adaptation between two sibling insect species by transcriptional profile analysis. BMC Genomics 14:582. doi: 10.1186/1471-2164-14-582

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, R. Q., Zhu, H. M., Ruan, J., Qian, W. B., Fang, X. D., Shi, Z. B., et al. (2010). De novo assembly of human genomes with massively parallel short read sequencing. Genome Res. 20, 265–272. doi: 10.1101/gr.097261.109

PubMed Abstract | CrossRef Full Text | Google Scholar

Livak, K. J., and Schmittgen, T. D. (2001). Analysis of relative gene expression data using real-time quantitative PCR and the 2(T)(-Delta Delta C) method. Methods 25, 402–408. doi: 10.1006/meth.2001.1262

PubMed Abstract | CrossRef Full Text | Google Scholar

Meng, X., Hu, J., and Ouyang, G. (2017). The isolation and identification of pathogenic fungi from Tessaratoma papillosa drury (Hemiptera: Tessaratomidae). PeerJ 5:e3888. doi: 10.7717/peerj.3888

PubMed Abstract | CrossRef Full Text | Google Scholar

Mollo, E., Garson, M. J., Polese, G., Amodeo, P., and Ghiselin, M. T. (2017). Taste and smell in aquatic and terrestrial environments. Nat. Prod. Rep. 34, 496–513. doi: 10.1039/c7np00008a

PubMed Abstract | CrossRef Full Text | Google Scholar

Oakeshott, J. G., Horne, I., Sutherland, T. D., and Russell, R. J. (2003). The genomics of insecticide resistance. Genome Biol. 4:202. doi: 10.1186/gb-2003-4-1-202

PubMed Abstract | CrossRef Full Text | Google Scholar

Pan, L., Ren, L. L., Chen, F., Feng, Y. Q., and Luo, Y. Q. (2016). Antifeedant activity of ginkgo biloba secondary metabolites against hyphantria cunea larvae: mechanisms and applications. PLoS One 11:e0155682. doi: 10.1371/journal.pone.0155682

PubMed Abstract | CrossRef Full Text | Google Scholar

Paula, D. P., Togawa, R. C., Costa, M. M. C., Grynberg, P., Martins, N. F., and Andow, D. A. (2016). Identification and expression profile of odorant-binding proteins in Halyomorpha halys (Hemiptera: Pentatomidae). Insect Mol. Biol. 25, 580–594. doi: 10.1111/imb.12243

PubMed Abstract | CrossRef Full Text | Google Scholar

Schmidt, J. M., Acebes-Doria, A., Blaauw, B., Kheirodin, A., Pandey, S., Lennon, K., et al. (2021). Identifying molecular-based trophic interactions as a resource for advanced integrated pest management. Insects 12:358. doi: 10.3390/insects12040358

PubMed Abstract | CrossRef Full Text | Google Scholar

Sparks, M. E., Bansal, R., Benoit, J. B., Blackburn, M. B., Chao, H., Chen, M. Y., et al. (2020). Brown marmorated stink bug, Halyomorpha halys (Stål), genome: putative underpinnings of polyphagy, insecticide resistance potential and biology of a top worldwide pest. BMC Genomics 21:227. doi: 10.1186/s12864-020-6510-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Sun, D. D., Huang, Y., Qin, Z., Zhan, H., Zhang, J., Liu, Y., et al. (2020). Identification of candidate olfactory genes in the antennal transcriptome of the stink bug Halyomorpha halys. Front. Physiol. 11:876. doi: 10.3389/fphys.2020.00876

PubMed Abstract | CrossRef Full Text | Google Scholar

Sun, M., Huang, D. J., Zhang, A. L., Khan, I., Yan, H. D., Wang, X. S., et al. (2020). Transcriptome analysis of heat stress and drought stress in pearl millet based on Pacbio full-length transcriptome sequencing. BMC Plant Biol. 20:323. doi: 10.1186/s12870-020-02530-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Sun, Y. L., Huang, L. Q., Pelosi, P., and Wang, C. Z. (2012). Expression in antennae and reproductive organs suggests a dual role of an odorant-binding protein in two sibling Helicoverpa species. PLoS One. 7:e30040. doi: 10.1371/journal.pone.0030040

PubMed Abstract | CrossRef Full Text | Google Scholar

Tran, H., Van, H. N., Muniappan, R., Amrine, J., Naidu, R., Gilbertson, R., et al. (2019). Integrated pest management of longan (Sapindales: Sapindaceae) in Vietnam. J. Integr. Pest Manag. 10:18. doi: 10.1093/jipm/pmz016

CrossRef Full Text | Google Scholar

Wu, Y. H., Kamiyama, M. T., Chung, C. C., Tzeng, H. Y., Hsieh, C. H., and Yang, C. S. (2020). Population monitoring, egg parasitoids, and genetic structure of the invasive litchi stink bug, Tessaratoma papillosa in Taiwan. Insects 11:690. doi: 10.3390/insects11100690

PubMed Abstract | CrossRef Full Text | Google Scholar

Wu, Z. Z., Qu, M. Q., Pu, X. H., Cui, Y., Xiao, W. Y., Zhao, H. X., et al. (2017). Transcriptome sequencing of Tessaratoma papillosa antennae to identify and analyze expression patterns of putative olfaction genes. Sci. Rep. 7:3070. doi: 10.1038/s41598-017-03306-7

PubMed Abstract | CrossRef Full Text | Google Scholar

Yang, X. M., Margolies, D. C., Zhu, K. Y., and Buschman, L. L. (2001). Host plant-induced changes in detoxification enzymes and susceptibility to pesticides in the twospotted spider mite (Acari : Tetranychidae). J. Econ. Entomol. 94, 381–387. doi: 10.1603/0022-0493-94.2.381

PubMed Abstract | CrossRef Full Text | Google Scholar

Yuan, Y. F., Li, L. S., Zhao, J. F., and Chen, M. (2020). Effect of tannic acid on nutrition and activities of detoxification enzymes and acetylcholinesterase of the fall webworm (Lepidoptera: Arctiidae). J. Insect. Sci. 20:8. doi: 10.1093/jisesa/ieaa001

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhan, H. X., Dewer, Y., Qu, C., Yang, S. Y., Luo, C., Li, L. J., et al. (2021). Molecular characterization of Donacia provosti (Coleoptera: Chrysomelidae) larval transcriptome by de novo assembly to discover genes associated with underwater environmental adaptations. Insects 12:281. doi: 10.3390/insects12040281

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, Z. M., Wu, W. W., and Li, G. K. (2009). Study of the alarming volatile characteristics of Tessaratoma papillosa using SPME-GC-MS. J. Chromatogr. Sci. 47, 291–296. doi: 10.1093/chromsci/47.4.291

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhao, D. F., Zheng, C. C., Shi, F. M., Xu, Y. B., Zong, S. X., and Tao, J. (2021). Expression analysis of genes related to cold tolerance in Dendroctonus valens. PeerJ 9:e10864. doi: 10.7717/peerj.10864

PubMed Abstract | CrossRef Full Text | Google Scholar

Keywords: Tessaratoma papillosa (Drury), transcriptome, olfactory genes, growth and development, digestion and detoxification

Citation: Cheng L, Han S, Jiang J, Li H and Peng L (2022) Transcriptome Analysis for Identification of Genes Related to Growth and Development, Digestion and Detoxification, Olfaction in the Litchi Stink Bug Tessaratoma papillosa. Front. Physiol. 12:774218. doi: 10.3389/fphys.2021.774218

Received: 11 September 2021; Accepted: 26 October 2021;
Published: 24 January 2022.

Edited by:

Annalisa Grimaldi, University of Insubria, Italy

Reviewed by:

Dandan Wei, Southwest University, China
Xin Yang, Institute of Vegetables and Flowers, Chinese Academy of Agricultural Sciences (CAAS), China

Copyright © 2022 Cheng, Han, Jiang, Li and Peng. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Haichao Li, lihaichao@cemps.ac.cn; Lingfei Peng, lingfeipeng@fafu.edu.cn

These authors have contributed equally to this work and share first authorship

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.