RNA-Seq reveals a xenobiotic stress response in the soybean aphid, Aphis glycines, when fed aphid-resistant soybean

Abstract

Background: While much recent research has expanded our understanding of the molecular interactions between aphids and their host plants, it is lacking for the soybean aphid, Aphis glycines. Since its North American invasion, A. glycines has become one of the most damaging insect pests on this important crop. Five soybean genes for host plant resistance to A. glycines have been identified, but populations of A. glycines have already adapted to overcome these resistance genes. Understanding the molecular interactions between resistant soybean and A. glycines can provide clues to its adaptation mechanisms. Here, we used RNA-Sequencing to compare and contrast A. glycines gene expression when fed resistant (Rag1) and susceptible soybean.

Results: Combining results from a previous A. glycines transcriptome, we generated 64,860 high quality transcripts, totaling 41,151,086 bases. Statistical analysis revealed 914 genes with significant differential expression. Most genes with higher expression in A. glycines on resistant plants (N = 352) were related to stress and detoxification such as cytochrome P450s, glutathione-S-transferases, carboxyesterases, and ABC transporters. A total of 562 genes showed lower transcript abundance in A. glycines on resistant plants. From our extensive transcriptome data, we also identified genes encoding for putative salivary effector proteins (N = 73). Among these, 6 effector genes have lower transcript abundance in A. glycines feeding on resistant soybean.

Conclusions: Overall, A. glycines exhibited a pattern typical of xenobiotic challenge, thereby validating antibiosis in Rag1, presumably mediated through toxic secondary metabolites. Additionally, this study identified many A. glycines genes and gene families at the forefront of its molecular interaction with soybean. Further investigation of these genes in other biotypes may reveal adaptation mechanisms to resistant plants.

Figure 1

A. glycines transcriptome annotation and comparative genomics. (A) Length distribution of 64,860 contigs in de novo assembly. Individual contigs are ordered on X-axis based on increasing size. (B) Ortholog hit ratio for transcripts calculated after BLASTx searches to genomes of A. pisum, B. mori, D. melanogaster, N. vitripennis, R. prolixus, and T. casteneum . (C) Venn diagram showing the number of transcript contigs with significant matches (unique and common) to genomes of A. pisum, D. melanogaster, R. prolixus, and T. casteneum. Significant matches (e value <1.0E-3) were calculated after pairwise comparisons (BLASTx) to each individual genome. (D) Comparison of GO term mappings distributions of A. glycines and A. pisum that belong to each of the three top-level GO categories (i.e. biological process, molecular function, and cellular component).

Figure 2

Gene expression changes in A. glycines due to Rag1 -soybean feeding. The expression (log2 fold change) of each gene between insects fed with resistant Rag1-soybean and those fed with susceptible plant is plotted against average expression level of each gene in both treatments. Fold change values for gene expression were considered significant if P values were <0.05. See materials and methods for details.

Figure 3

qRT-PCR validation of RNA-Seq results. Validation of gene expression (14 genes) using Pearson’s correlation (r) between fold changes (log2 scale) observed in qRT-PCR and RNA-Seq results.

Figure 4

Gene expression of effectors in salivary glands of A. glycines . (A) A dissected out salivary gland from an A. glycines adult as viewed through a microscope. The principal salivary gland (PSG), the salivary duct (SD), and the accessory salivary gland (ASG) are indicated. (B) Results of semi quantitative PCR for expression analysis of effector genes in salivary gland and carcass (adult minus salivary gland and developing embryos) are shown. The PCR reactions were run for 35 cycles for all primer pairs except for AyDI2, where it was 40. The effector names and other details are provided in Table 4.

Figure 5

Gene expression comparison among starved and fed A. glycines. Bars represent the relative mRNA levels of different genes in A. glycines using qRT-PCR. The mean (± S.E) expression level is represented for three biological replicates for A. glycines fed with resistant soybean (green bars), susceptible soybean (blue bars), and starved (grey bars). The elongation factor-la (AyEF1α) gene was used as an internal control for cDNA [77]. More details on genes and primer sequences are provided in Additional file 12. (* P < 0.05).