Global Gene Expression in Cotton Fed Upon by Aphis gossypii and Acyrthosiphon gossypii (Hemiptera: Aphididae)

Cotton, Gossypium hirsutum L., is an economically important crop in global agricultural and textile industries. Aphids (Hemiptera: Aphididae) have become important pests of cotton production worldwide and are now considered the dominant pest species of cotton-growing areas in China (Lu and Liang 2016, Lu et al. 2020).

Although insects with piercing-sucking mouthparts cause less mechanical damage to host plants than insects with chewing mouthparts, the damage cycle is longer, resulting in plants exhibiting mild but persistent defense responses to the attack (Moran and Thompson 2001, Thompson and Goggin 2006). Plants possess specialized structures and substances, for example, wax, hairs, spines, glands, and different degrees of ossification or silicification of some tissues on the plant surface, as defense mechanisms (Zhang et al. 2013c). While these specialized tissues and structures help to resist pests, they cannot completely defend the host plant from attack by phytophagous insects.

Insects with piercing-sucking mouthparts can penetrate the plant epidermis to imbibe plant fluids from the phloem and xylem. This may lead to depolarization of the lipid membrane or a disturbance in the transmembrane ion flow in plant cells, causing a change in transmembrane potential across the cell membrane and a change in signal transduction: for example, the calcium ion (Ca²⁺) influx (Bricchi et al. 2012, Luo et al. 2017, Vincent et al. 2017, Yan et al. 2018). Plants also regulate the activity of key proteases to degrade and eliminate reactive oxygen chemical species, including phenolic and quinone compounds (Chen et al. 2011, Liu and Lan 2009, Luo et al. 2008, Tjallingii 2006, Wang et al. 2011, Wu et al. 2015). This defense mechanism likely evolved while resisting the chemicals injected during feeding of insects with piercing-sucking mouthparts (Boyko et al. 2006, Voelckel et al. 2004, Zhou et al. 2009).

Aphis gossypii Glover and Acyrthosiphon gossypii Mordvilko (Hemiptera: Aphididae) are common pests of cotton in China and worldwide. Their feeding may cause leaf curl, stunting of plant growth, and slowing of plant development (Jacobson and Croft 1998). Honeydew produced while feeding may serve as a nutritional substrate for molds that can interfere with light absorption and photosynthetic activity (Hullé et al. 2020). Cotton plants also exhibit a series of physiological and metabolic reactions with stress associated with A. gossypii and Ac. gossypii feeding. These include increased activity of the antioxidative enzymes catalase (CAT), peroxidase (POD), polyphenol oxidase (PPO), lipoxygenase (LOX), and other defense enzymes (Chen et al. 2015; Li et al. 1998a, b).

The alteration in the level of soluble sugars, free proline, and other nutrients due to aphid attack initiates immune defense mechanisms against aphid feeding (Li et al. 2008, Patima et al. 2018). We have previously shown that A. gossypii and Ac. gossypii feeding can cause various defense responses in cotton. For example, A. gossypii was found to cause changes in chlorophyll, soluble protein, proline, malondialdehyde content, and antioxidant enzyme activity in cotton at both the boll and bud stages and, with the extension of stress time, cotton defense ability was enhanced (Deng et al. 2013, Yan et al. 2013). Feeding by Ac. gossypii altered the level of soluble sugar, soluble protein, chlorophyll, carotenoids, malonaldehyde, and the activity of POD in cotton; nutrient metabolism and cell permeability also were altered. At the same time, the activity of related defense enzymes was induced (Zhang et al. 2020). However, the specific gene expression changes that mediate cotton defense responses to A. gossypii and Ac. gossypii attacks remain poorly understood. The aim of this study was, therefore, to use transcriptome sequencing to investigate the differential expression of genes related to biological processes, cell components, and molecular functions in cotton following feeding by A. gossypii and Ac. gossypii.

Experimental treatments. Cotton (New Upland Early Maturity 44 variety) seeds were soaked for 1 h at 55°C, allowed to germinate at room temperature for 24 h, and then planted in vermiculite in plastic basins (12-cm height, 10-cm diameter), which were maintained in environmentally controlled incubators on a 16:8-h photoregime at 24°C from midnight to 8:00 a.m. and 27°C from 8:00 a.m. to 11:59 p.m. Tests with aphids were initiated when cotton seedlings had grown to two true leaves.

Aphis gossypii and Ac. gossypii aphids used in the study were from colonies that had been subcultured on cotton seedlings for more than 30 generations. Aphids from these colonies were transferred individually to cotton seedlings using a fine brush. Nine seedlings were infested with 16 Ac. gossypii per plant; nine were infested with 16 A. gossypii per plant; nine were infested with eight Ac. gossypii and eight A. gossypii per plant; and nine were not infested and served as a check. Once aphids were transferred to the plants, whole plants were covered and placed in an incubator with controlled lighting.

After 3 d, the aphids on the cotton seedlings were removed, and the new plant growth at the top of each plant was excised with scissors and placed individually in 1.5-ml sterile centrifuge tubes. These were immediately placed in liquid nitrogen and transferred to –80°C. Each of the four treatments was replicated three times.

RNA extraction, sequencing, and data analysis. Sample RNA extraction, quality detection, transcriptome sequencing, and statistical analyses were commissioned and performed by Beijing Nuohe Zhiyuan Technology Co., Ltd. (Beijing, China). Briefly, total RNA was extracted from the cotton samples using TRIzol (Tiangen Biotech [Beijing] Co., Ltd., Beijing, China) according to manufacturer's instructions. RNA purity was checked using the NanoPhotometer^® spectrophotometer (Implen, Corston, United Kingdom). RNA concentration was measured using the Qubit^® RNA Assay Kit in Qubit^® 2.0 Flurometer (Life Technologies, Carlsbad, CA), and RNA integrity was assessed with the RNA Nano 6000 Assay Kit of the Bioanalyzer 2100 system (Agilent Technologies, Santa Clara, CA). All samples had an RNA integrity number (RIN) above 6.7. RNA sequencing libraries were prepared with NEBNext^® UltraTM RNA Library Prep Kit for Illumina^® (Illumina, Inc., San Diego, CA) and sequenced on an Illumina Hiseq^™ 2500 platform at an average depth of ∼66 million reads per sample. Raw sequencing reads were quality assessed with FastQ. To pass the initial quality control check, the average Phred score of each base position across all reads had to be at least 30. Reads were further processed by cutting individual low-quality bases and removing adapter and other Illumina-specific sequences with ng-qc using default parameters. HISAT2-2.0.4 was then used to map the trimmed reads to the cotton AD1_ZJU_v2.1 reference genome (Kim et al. 2015, Mortazavi et al. 2008). To quantify gene expression levels, mapped reads were summarized at the gene level using HTSeq version 0.6.0 (Anders 2010). Differential expression analyses were performed with DESeq2 R package (version 1.10.1), and gene ontology (GO) enrichment analyses were conducted using the clusterProfiler R package (Anders and Huber 2012, Wang et al. 2010). The significance threshold used was the adjusted P value of 0.05 and absolute fold change of 2 for the differential expression analysis and adjusted P value less than 0.05 for GO enrichment analysis (Robinson et al. 2010, Young et al. 2010). We used clusterProfiler R package to test the statistical enrichment of differential expression genes in the Kyoto Encyclopedia of Genes and Genomes (KEGG) pathways (Kanehisa et al. 2008).

RNA quality and sequencing data results. All indicators of RNA quality (260/280 and 260/230 absorbance ratios and RIN values) were above accepted quality thresholds for all samples, and the extracted RNA could therefore be used for transcriptome sequencing (Table 1). Sequencing data statistics showed that for each sample the sequencing data error rate was ≤0.03%; >97% of sequences had a Phred score of at least Q20, >92% of sequences had a Phred score of at least Q30, and the guanine and cytosine (GC) content was stable around 44%, indicating that the sequencing data had sufficient quality to be used for subsequent analyses (Table 2). Furthermore, the average sequencing depth was >60 million reads per sample, and the filtered sequencing reads had a high alignment with the reference genome, suggesting an appropriate sequencing depth for differential expression analyses (Table 3).

Gene expression analysis. According to the results of gene expression in cotton fed upon by A. gossypii and Ac. gossypii (Fig. 1), the fragments per kilo base of exon per million reads (FPKM) values of gene expression were divided into five levels. Within the range of FPKM values of 0–15, the amount of gene expression in cotton had no significant relationship with the feeding of A. gossypii and Ac. gossypii. In the range of FPKM values of 15–60, A. gossypii feeding had no significant effect on cotton gene expression compared to the control group, while feeding by Ac. gossypii alone and in combination with A. gossypii resulted in significantly reduced cotton gene expression. In the range of FPKM value >60, the amount of gene expression in cotton fed on by Ac. gossypii was significantly lower compared to the other groups (Fig. 1).

Differential gene expression analysis. An overview of the differentially expressed genes in cotton fed upon by A. gossypii and Ac. gossypii is shown in Table 4. A total of 6,565 genes (3,310 genes upregulated and 3,255 genes downregulated) were differentially expressed between cotton fed on by Ac. gossypii and cotton free of aphids (Fig. 2A); 823 genes (470 upregulated and 353 downregulated) were differentially expressed between cotton fed upon by A. gossypii and cotton free of aphids (Fig. 2B); and 2,379 genes (1,003 upregulated and 1,376 downregulated) were differentially expressed between cotton fed upon by A. gossypii and cotton fed upon by Ac. gossypii (Fig. 2C). The upregulated genes in cotton fed upon by Ac. gossypii were mainly concentrated in the photosynthetic metabolism pathway, biosynthesis of secondary metabolites, and the pentose phosphate metabolism pathway. The upregulated genes in cotton fed upon by A. gossypii were mainly concentrated in amino acid metabolism, plant–pathogen interaction, and terpene biosynthesis. When mixed populations of A. gossypii and Ac. gossypii fed on the cotton, the upregulated genes were mainly concentrated in circadian rhythm regulation, photosynthesis, and galactose metabolism pathways (Fig. 3). Among these differentially expressed genes, 280 were common in the first two comparisons, and 96 genes were shared across all three comparisons (Fig. 4).

GO enrichment analysis of differentially expressed genes.Aphis gossypii feeding on cotton resulted in regulation of biological processes related to synthesis and metabolism of tetraterpenoid, carotenoid, and methionine, and molecular functionals related to flavin adenine dinucleotide (FAD) binding and the activity of 5-methyltetrahydropteroyltriglutamate-homocysteine S-methyltransferase and oxidoreductase, among others (Fig. 5A). Acyrthosiphon gossypii feeding on cotton had a modulatory effect on genes involved in photosynthesis, steroid biosynthesis, lipid biosynthesis, and other biological processes, as well as cellular components (e.g., thylakoid, photosystem II oxygen complex, oxidoreductase complex, and cell membrane components in the photosynthetic system). The expression of genes related to molecular functions, such as oxidoreductase activity and steroid dehydrogenase activity, were also affected (Fig. 5B). Feeding by mixed populations of A. gossypii and Ac. gossypii had a significant effect on biological processes related to phosphorus signal transduction in cotton (Fig. 5C).

KEGG enrichment analysis of differentially expressed genes. In the process of cotton differential gene enrichment, porphyrin and chlorophyll metabolism were the most enriched after feeding by Ac. gossypii, with 41 genes in the pathway. The second was monoterpenoid biosynthesis, in which six genes were involved in this pathway, and the least enriched were metabolic pathways (Fig. 6A). Thiamine metabolism was the most enriched in cotton after feeding by A. gossypii, with five genes in the pathway. The second was diterpenoid biosynthesis, which had four genes. Biosynthesis of secondary metabolites was the least enriched (Fig. 6B). Plant circadian rhythm was the most enriched in cotton after feeding by both Ac. gossypii and A. gossypii, with 26 genes in the pathway. The second was thiamine metabolism, which had five genes. The least enrichment was in metabolic pathways (Fig. 6C).

We found that Ac. gossypii feeding on cotton resulted in regulation of 1,221 annotated genes (Table 5), and five KEGG pathways were significantly enriched, including metabolic pathways, porphyrin and chlorophyll metabolism, biosynthesis of secondary metabolites, biosynthesis of carotenoids, and the pentose phosphate pathway (Table 6). The numbers of differentially expressed genes annotated in metabolic pathways and biosynthesis of secondary metabolites were 723 and 403, respectively. Thirty-eight annotated genes were differentially expressed after A. gossypii feeding (Table 5), and four enriched pathways were identified, including thiamine metabolism, glutathione metabolism, plant–pathogen interaction, and sesquiterpene and triterpenoid biosynthesis (Table 6). Among these, the plant–pathogen interaction pathway included the largest number of annotated genes (16). A total of 355 annotated genes were differentially expressed after feeding by mixed populations of A. gossypii and Ac. gossypii (Table 5); six KEGG pathways were significantly enriched, including circadian rhythm regulation, photosynthesis, porphyrin and chlorophyll metabolism, galactose metabolism, and flavonoid biosynthesis (Table 6).

Plants attacked and fed upon by insects with piercing-sucking mouthparts activate related resistance genes (Park et al. 2005), induce an emergency response to the injury at the feeding site, activate the whole-plant defense system to reduce plant damage, and prepare for a rapid defense response to future disturbances (Martinez-Medina et al. 2016). Indeed, host plant defense regulatory mechanisms may differ with insect species, developmental stage, or feeding mechanism (chewing versus piercing-sucking mouthparts), as well as host plant species and characteristics. Yet, all play an important role in host plant defense against pests.

Hettenhausen et al. (2016) demonstrated that feeding by Spodoptera exigua Hübner or Aphis glycines Matsumura increased calcium-dependent protein kinase transcription in soybean, Glycine max (L.) Merrill. Sytykiewicz (2016) reported that feeding by Rhopalosiphum padi (L.) and Sitobion avenae (F.) significantly affected expression of rbohA and rbohD in maize, Zea mays L. In thale cress, Arabidopis thaliana (L.) Heynh., feeding by Bemisia tabaci (Gennadius) biotype B nymphs induced the upstream jasmonic acid (JA) response genes LOX2 and OPR3 and inhibited the downstream JA response gene VSP1, while feeding by adults significantly inhibited the expression of LOX2 and OPR3 (Zhang et al. 2013a, b). In rice (Oryza sativa L.), feeding by either chewing or sucking insects affected the ethylene and JA pathways, and OsHI-LOX was a key gene in JA synthesis (Ma et al. 2020, Zhou et al. 2009). Aphids feeding on tobacco, Nicotiana tabacum L., foliage induced significantly fewer differentially expressed genes compared to feeding by mirids, mealybugs, or lepidopteran larvae (Heidel and Baldwin 2004). Our results showing differential expression of genes and the occurrence of enriched metabolic pathways in cotton after feeding by Ac. gossypii and A. gossypii further support those findings and demonstrate their involvement in the host plant defense response in cotton.

Physiological metabolic pathways are important regulatory pathways for plants to initiate defense responses. Related metabolic pathways participate in and complement various defense mechanisms such as local defense, systemic defense, and direct defense of plants (Dicke and Poecke 2002, Orians 2005, Vignutelli et al. 1998). For example, Apolygus lucorum Meyer-Dur feeding induced significant changes in flavonoids, phenols, chymotrypsin inhibitors, condensed tannins, and amino acids in grape (Vitis spp.) leaves (Gao et al. 2019). Levels of soluble sugar, soluble protein, and chlorophyll in leaves of Mikania micrantha Kunth increased, while the activity of CAT, superoxide dismutase, and POD decreased after feeding by Pachypeltis sp. (Li et al. 2018). Sitobion avenae feeding on cabbage, Brassica oleracea (L.), and wheat, Triticum aestivum L., could induce increases in PPO, POD, and phenylalanine ammonia lyase activity (Han et al. 2009, Zhang et al. 2005). These studies show that piercing-sucking insect herbivory can cause changes in the reactive oxygen species system, secondary metabolite synthesis, and other physiological metabolic pathways in host plants. Our results in this present study demonstrated that feeding by the aphids Ac. gossypii and A. gossypii also induced changes in multiple physiological metabolic pathways in cotton. These included photosynthetic and secondary metabolic pathways which could improve the cotton plant ability to compensate for damage or loss of photosynthates or other nutrients.

Our findings also showed that feeding by Ac. gossypii and A. gossypii on cotton significantly affected the functional expression of oxidoreductase enzymes in the host plant. The oxidoreductase system, including oxidoreductase lipoxygenase, propylene oxide synthase, propylene oxide cyclase, peroxidase, and polyphenol oxidase, is reported as an important protective enzyme system in defense reactions in cotton (Chung et al. 2013, Si et al. 2020, Ximénez-Embún et al. 2017). We further postulate that cotton initiates oxidoreductase gene expression immediately upon incurring pest damage, thus enhancing the host plant resistance to or tolerance of aphid feeding by increasing the level of protective enzymes (Yan et al. 2013, Zhang et al. 2020).

In addition, we found that when fed upon by either Ac. gossypii or A. gossypii, expression of photosynthesis-related genes was increased in cotton, which supports the findings of Gutsche et al. (2009) that insect feeding can upregulate the expression of photosynthesis-related genes in plants and the conclusion of Kangasjarvi et al. (2012) that photosynthesis is involved in plant defense responses as well as plant physiological functions as a remedy for carbon loss. Furthermore, feeding by combined populations of Ac. gossypii and A. gossypii on cotton significantly affected plant biological processes (e.g., cotton phosphorescence signal transduction). These results provide insight into mechanisms underlying the observed increase in chlorophyll and carotenoid content in cotton leaves when cotton is damaged by Ac. gossypii and A. gossypii (Deng et al. 2013, Zhang et al. 2020).

In our KEGG enrichment analysis, Ac. gossypii increased the expression of biosynthetic pathways of secondary metabolites in cotton, while feeding by A. gossypii increased the expression of sesquiterpenes and triterpenoids. When the two species fed together on the same plant, the expression of flavonoid biosynthesis genes increased. Therefore, under the stress of Ac. gossypii and A. gossypii, cotton initiates defense responses through different pathways involving secondary metabolism. It is known that plant metabolites, including flavonoids, terpenoids, alkaloids, and other secondary metabolites, play an important role in insect feeding induction which, when ingested by the insect, can inhibit digestion, affect feeding, or even kill the insect (Chen et al. 2019, Howe and Jander 2008). Previous studies have shown that Ac. gossypii and A. gossypii feeding increased levels of tannins, flavonoids, total phenols, and other secondary substances in cotton (Liu and Yang 1993, Wu 2020, Zhang 2020) and increased the activity of secondary metabolic enzymes in cotton (Li et al. 1998b, Lu et al. 2017).

When plants are fed upon by insects, they not only synthesize secondary metabolites that are toxic and deterrent, but they also produce changes in primary metabolites such as proteins and soluble sugars (Sulpice and McKeown 2015, Sun et al. 2013). Cotton plants fed upon by Ac. gossypii and A. gossypii respond by increasing soluble protein and sugar content as a defensive mechanism (Deng et al. 2013, Patima et al. 2018, Yan et al. 2013, Zhang et al. 2020), which corresponds to an acceleration of biosynthesis and biological metabolism. Our GO and KEGG enrichment analyses indicated that steroid biosynthesis, lipid biosynthesis, and the pentose phosphate pathway were enhanced in cotton after feeding by Ac. gossypii. Those analyses also showed that cotton on fed by A. gossypii exhibited enhanced methionine metabolism, and that feeding by mixed populations of the aphids enhanced the differential expression of galactose metabolism.

These cotton plant reactions to aphid attack are defense mechanisms. Sterols involved in steroid biosynthesis, lipid biosynthesis, and methionine metabolism play an important role in cell wall formation, cell elongation, and development (Carland et al. 2002, Catterou et al. 2001, Clouse and Sasse 1998, He et al. 2003), while methionine is directly involved in protein biosynthesis (Giovanelli et al. 1985). Metabolites from galactose metabolism can promote the cell wall formation (Atmodjo et al. 2013) and increase the content of soluble sugars in plants (Thoden and Holden 2005). These plant defense mechanisms are energy-consuming processes (Coley et al. 1985, Mooney and Gulmon 1982, Rhoades 1979), and ATP is continuously provided for these processes through the pentose phosphate pathway. Collectively, these biosynthetic processes strongly influence plant morphology, protein and carbohydrate synthesis, and continuous plant defense functions (Limdsey et al. 2003, Schaller 2003).

Molecular studies of plant-pest interactions can reveal crop insect resistance mechanisms. Antibiotic-related substances, such as disease-related proteins, are rapidly produced when plants are fed upon by sucking insects (Park et al. 2005). Furthermore, oxygen-burst reactions occur at injured sites of plants, resulting in accumulation of protin I proteins and injury responses (Kaloshian 2004), activation of mitogen-activated protein kinases, synthesis and interaction of phytohormones, and a series of stress responses in plants (Erb et al. 2012, Zebelo and Maffei 2015). Aphids may also transmit viral plant diseases while feeding (Fereres and Moreno 2009). Cotton will immediately initiate immune factors to resist viral infection (Kørner et al. 2013, Mandadi and Scholthof 2013). In support, our KEGG enrichment analysis showed that the largest number of differentially expressed genes were annotated to the plant–pathogen interaction pathway in cotton fed upon by the aphids.

It should be noted that when A. gossypii and Ac. gossypii feed on cotton in mixed populations, the genes regulating circadian rhythm are differentially expressed, which may be related to the regulation of nutrient homeostasis (Haydon et al. 2015), hormone synthesis and signal transduction (Atamian and Harmer 2016), redox reaction (Zhou et al. 2015), and the changes in levels of some major osmotic regulators (Greenham and McClung 2015). These responses indicate that the biological clock of cotton has a complex regulation when stressed by aphid feeding.

In conclusion, although the gene expression and metabolic pathways of cotton defense responses induced by A. gossypii and Ac. gossypii differ, they all enhance the defense response of cotton through regulating pathways related to photosynthetic substances, oxidoreductase activity, secondary metabolism, and other metabolic activities. This is similar to the defense response pathways induced by most insects with piercing-sucking mouthparts and feeding habits. When the two aphid species damaged the plant simultaneously, the genes regulating cotton photosynthetic phosphorus signal transduction, circadian rhythm regulation, porphyrin and chlorophyll metabolism, photosynthesis, galactose metabolism, flavonoid biosynthesis, and other activities were significantly expressed. Our study thus provides new insights into the complex mechanisms underlying cotton defense responses against aphid attacks. However, in this study, only single omics analysis was used to analyze the mechanisms of cotton defense against aphids. In the future, multigroup analysis should be used to conduct more in-depth analysis at the molecular, metabolic, and protein levels, so as to provide a more comprehensive elucidation of the mechanism of cotton defense against aphids.

This work was supported by the National Natural Science Foundation of China (grant no. 31660519).